OGS
VtuInterface-impl.h
Go to the documentation of this file.
1 
19 #include <vtkNew.h>
20 #include <vtkSmartPointer.h>
21 #include <vtkUnstructuredGrid.h>
22 
23 #include "BaseLib/Logging.h"
25 #include "VtuInterface.h"
26 
27 #ifdef USE_PETSC
28 #include <mpi.h>
29 #endif
30 
31 class vtkXMLPUnstructuredGridWriter;
32 
33 namespace MeshLib
34 {
35 namespace IO
36 {
37 template <typename UnstructuredGridWriter>
38 bool VtuInterface::writeVTU(std::string const& file_name,
39  [[maybe_unused]] const int num_partitions,
40  [[maybe_unused]] const int rank)
41 {
42  if (!_mesh)
43  {
44  ERR("VtuInterface::write(): No mesh specified.");
45  return false;
46  }
47 
48  vtkNew<MeshLib::VtkMappedMeshSource> vtkSource;
49  vtkSource->SetMesh(_mesh);
50 
51  vtkSmartPointer<UnstructuredGridWriter> vtuWriter =
52  vtkSmartPointer<UnstructuredGridWriter>::New();
53 
54  vtkSource->Update();
55  vtuWriter->SetInputData(vtkSource->GetOutput());
56 
57  if (_use_compressor)
58  {
59  vtuWriter->SetCompressorTypeToZLib();
60  }
61  else
62  {
63  vtuWriter->SetCompressorTypeToNone();
64  }
65 
66  vtuWriter->SetDataMode(_data_mode);
67  if (_data_mode == vtkXMLWriter::Appended)
68  {
69  vtuWriter->SetEncodeAppendedData(1);
70  }
71  if (_data_mode == vtkXMLWriter::Ascii)
72  {
73  vtkSource->Update();
74  vtkSmartPointer<vtkUnstructuredGrid> tempGrid =
75  vtkSmartPointer<vtkUnstructuredGrid>::New();
76  tempGrid->DeepCopy(vtkSource->GetOutput());
77  vtuWriter->SetInputDataObject(tempGrid);
78  }
79 
80  vtuWriter->SetFileName(file_name.c_str());
81 
82  if constexpr (std::is_same_v<UnstructuredGridWriter,
83  vtkXMLPUnstructuredGridWriter>)
84  {
85  vtuWriter->SetGhostLevel(1);
86  vtuWriter->SetNumberOfPieces(num_partitions);
87  vtuWriter->SetStartPiece(rank);
88  vtuWriter->SetEndPiece(rank);
89  }
90 
91 #ifdef VTK_USE_64BIT_IDS
92  vtuWriter->SetHeaderTypeToUInt64();
93  // set SetIdTypeToInt64() as well?
94 #endif
95 
96  return (vtuWriter->Write() > 0);
97 }
98 
99 } // end namespace IO
100 } // end namespace MeshLib
void ERR(fmt::format_string< Args... > fmt, Args &&... args)
Definition: Logging.h:44
VtkMappedMeshSource is a source class to transform OGS meshes into complete vtkUnstructuredGrids....
Implementation of the VtuInterface class.
bool writeVTU(std::string const &file_name, const int num_partitions=1, const int rank=1)
const MeshLib::Mesh * _mesh
Definition: VtuInterface.h:68