GEOS
VTKUtilities.hpp
Go to the documentation of this file.
1 /*
2  * ------------------------------------------------------------------------------------------------------------
3  * SPDX-License-Identifier: LGPL-2.1-only
4  *
5  * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC
6  * Copyright (c) 2018-2024 TotalEnergies
7  * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University
8  * Copyright (c) 2023-2024 Chevron
9  * Copyright (c) 2019- GEOS/GEOSX Contributors
10  * All rights reserved
11  *
12  * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details.
13  * ------------------------------------------------------------------------------------------------------------
14  */
15 
20 #ifndef GEOS_MESH_GENERATORS_VTKUTILITIES_HPP
21 #define GEOS_MESH_GENERATORS_VTKUTILITIES_HPP
22 
23 #include "common/DataTypes.hpp"
24 #include "common/MpiWrapper.hpp"
26 
27 #include <vtkDataSet.h>
28 #include <vtkMultiProcessController.h>
29 #include <vtkSmartPointer.h>
30 
31 #include <numeric>
32 #include <unordered_set>
33 
34 namespace geos
35 {
36 namespace vtk
37 {
38 
43 {
44  parmetis,
45  ptscotch,
46 };
47 
50  "parmetis",
51  "ptscotch" );
52 
59 using CellMapType = std::map< ElementType, std::unordered_map< int, stdVector< vtkIdType > > >;
60 
65 vtkSmartPointer< vtkMultiProcessController > getController();
66 
70 class AllMeshes
71 {
72 public:
73  AllMeshes() = default;
74 
80  AllMeshes( vtkSmartPointer< vtkDataSet > const & main,
81  std::map< string, vtkSmartPointer< vtkDataSet > > const & faceBlocks )
82  : m_main( main ),
83  m_faceBlocks( faceBlocks )
84  { }
85 
89  vtkSmartPointer< vtkDataSet > getMainMesh()
90  {
91  return m_main;
92  }
93 
97  std::map< string, vtkSmartPointer< vtkDataSet > > & getFaceBlocks()
98  {
99  return m_faceBlocks;
100  }
101 
106  void setMainMesh( vtkSmartPointer< vtkDataSet > main )
107  {
108  m_main = main;
109  }
110 
115  void setFaceBlocks( std::map< string, vtkSmartPointer< vtkDataSet > > const & faceBlocks )
116  {
117  m_faceBlocks = faceBlocks;
118  }
119 
120 private:
122  vtkSmartPointer< vtkDataSet > m_main;
123 
125  std::map< string, vtkSmartPointer< vtkDataSet > > m_faceBlocks;
126 };
127 
135 AllMeshes loadAllMeshes( Path const & filePath,
136  string const & mainBlockName,
137  string_array const & faceBlockNames );
138 
146 
160 AllMeshes
161 redistributeMeshes( integer const logLevel,
162  vtkSmartPointer< vtkDataSet > loadedMesh,
163  std::map< string, vtkSmartPointer< vtkDataSet > > & namesToFractures,
164  MPI_Comm const comm,
165  PartitionMethod const method,
166  int const partitionRefinement,
167  int const useGlobalIds,
168  string const & structuredIndexAttributeName,
169  int const numPartZ );
170 
179 CellMapType buildCellMap( vtkDataSet & mesh,
180  string const & attributeName );
181 
189 void printMeshStatistics( vtkDataSet & mesh,
190  CellMapType const & cellMap,
191  MPI_Comm const comm );
192 
199 vtkDataArray * findArrayForImport( vtkDataSet & mesh, string const & sourceName );
200 
208 bool hasArray( vtkDataSet & mesh, string const & sourceName );
209 
216 string buildCellBlockName( ElementType const type, int const regionId );
217 
225  vtkDataArray * vtkArray,
226  dataRepository::WrapperBase & wrapper );
227 
235  vtkDataArray * vtkArray,
236  dataRepository::WrapperBase & wrapper );
237 
243 void importRegularField( vtkDataArray * vtkArray,
244  dataRepository::WrapperBase & wrapper );
245 
246 
247 } // namespace vtk
248 
258 void writeNodes( integer const logLevel,
259  vtkDataSet & mesh,
260  string_array & nodesetNames,
261  CellBlockManager & cellBlockManager,
262  const geos::R1Tensor & translate,
263  const geos::R1Tensor & scale );
264 
273 void writeCells( integer const logLevel,
274  vtkDataSet & mesh,
275  vtk::CellMapType const & cellMap,
276  string const & structuredIndexAttributeName,
277  CellBlockManager & cellBlockManager );
278 
288 void writeSurfaces( integer const logLevel,
289  vtkDataSet & mesh,
290  const geos::vtk::CellMapType & cellMap,
291  CellBlockManager & cellBlockManager );
292 
298 std::pair< real64, real64 > getGlobalLengthAndOffset( vtkDataSet & mesh );
299 
300 } // namespace geos
301 
302 #endif /* GEOS_MESH_GENERATORS_VTKUTILITIES_HPP */
AllMeshes loadAllMeshes(Path const &filePath, string const &mainBlockName, string_array const &faceBlockNames)
Load the VTK file into the VTK data structure.
void importMaterialField(stdVector< vtkIdType > const &cellIds, vtkDataArray *vtkArray, dataRepository::WrapperBase &wrapper)
Imports 2d and 3d arrays from vtkArray to wrapper, only for cellIds.
vtkDataArray * findArrayForImport(vtkDataSet &mesh, string const &sourceName)
Collect the data to be imported.
std::map< ElementType, std::unordered_map< int, stdVector< vtkIdType > > > CellMapType
Type of map used to store cell lists.
vtkSmartPointer< vtkMultiProcessController > getController()
Return a VTK controller for multiprocessing.
void printMeshStatistics(vtkDataSet &mesh, CellMapType const &cellMap, MPI_Comm const comm)
Print statistics for a vtk mesh.
void importRegularField(stdVector< vtkIdType > const &cellIds, vtkDataArray *vtkArray, dataRepository::WrapperBase &wrapper)
Imports 1d and 2d arrays from vtkArray to wrapper, only for cellIds.
string buildCellBlockName(ElementType const type, int const regionId)
build cell block name from regionId and cellType
CellMapType buildCellMap(vtkDataSet &mesh, string const &attributeName)
Collect lists of VTK cell indices organized by type and attribute value.
bool hasArray(vtkDataSet &mesh, string const &sourceName)
Check if the vtk mesh as a cell data field of name sourceName.
PartitionMethod
Choice of advanced mesh partitioner.
@ ptscotch
Use PTScotch library.
@ parmetis
Use ParMETIS library.
stdVector< int > findNeighborRanks(stdVector< vtkBoundingBox > boundingBoxes)
Compute the rank neighbor candidate list.
AllMeshes redistributeMeshes(integer const logLevel, vtkSmartPointer< vtkDataSet > loadedMesh, std::map< string, vtkSmartPointer< vtkDataSet > > &namesToFractures, MPI_Comm const comm, PartitionMethod const method, int const partitionRefinement, int const useGlobalIds, string const &structuredIndexAttributeName, int const numPartZ)
Generate global point/cell IDs and redistribute the mesh among MPI ranks.
The CellBlockManager class provides an interface to ObjectManagerBase in order to manage CellBlock da...
Class describing a file Path.
Definition: Path.hpp:35
Base class for all wrappers containing common operations.
Definition: WrapperBase.hpp:56
Gathers all the vtk meshes together.
std::map< string, vtkSmartPointer< vtkDataSet > > & getFaceBlocks()
void setMainMesh(vtkSmartPointer< vtkDataSet > main)
Defines the main 3d mesh for the simulation.
vtkSmartPointer< vtkDataSet > getMainMesh()
void setFaceBlocks(std::map< string, vtkSmartPointer< vtkDataSet > > const &faceBlocks)
Defines the face blocks/fractures.
AllMeshes(vtkSmartPointer< vtkDataSet > const &main, std::map< string, vtkSmartPointer< vtkDataSet > > const &faceBlocks)
Builds the compound from values.
stdVector< string > string_array
A 1-dimensional array of geos::string types.
Definition: DataTypes.hpp:361
std::pair< real64, real64 > getGlobalLengthAndOffset(vtkDataSet &mesh)
Compute the global length of the mesh and its offset.
void writeNodes(integer const logLevel, vtkDataSet &mesh, string_array &nodesetNames, CellBlockManager &cellBlockManager, const geos::R1Tensor &translate, const geos::R1Tensor &scale)
Build all the vertex blocks.
void writeCells(integer const logLevel, vtkDataSet &mesh, vtk::CellMapType const &cellMap, string const &structuredIndexAttributeName, CellBlockManager &cellBlockManager)
Build all the cell blocks.
mapBase< TKEY, TVAL, std::integral_constant< bool, true > > map
Ordered map type.
Definition: DataTypes.hpp:339
void writeSurfaces(integer const logLevel, vtkDataSet &mesh, const geos::vtk::CellMapType &cellMap, CellBlockManager &cellBlockManager)
Build the "surface" node sets from the surface information.
ElementType
Denotes type of cell/element shape.
Definition: ElementType.hpp:32
int integer
Signed integer type.
Definition: DataTypes.hpp:81
ENUM_STRINGS(LinearSolverParameters::SolverType, "direct", "cg", "gmres", "fgmres", "bicgstab", "richardson", "preconditioner")
Declare strings associated with enumeration values.
internal::StdVectorWrapper< T, Allocator, USE_STD_CONTAINER_BOUNDS_CHECKING > stdVector