MFEM  v4.2.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
pmesh.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2020, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #ifndef MFEM_PMESH
13 #define MFEM_PMESH
14 
15 #include "../config/config.hpp"
16 
17 #ifdef MFEM_USE_MPI
18 
19 #include "../general/communication.hpp"
20 #include "../general/globals.hpp"
21 #include "mesh.hpp"
22 #include "pncmesh.hpp"
23 #include <iostream>
24 
25 namespace mfem
26 {
27 #ifdef MFEM_USE_PUMI
28 class ParPumiMesh;
29 #endif
30 
31 /// Class for parallel meshes
32 class ParMesh : public Mesh
33 {
34 protected:
35  ParMesh() : MyComm(0), NRanks(0), MyRank(-1),
36  have_face_nbr_data(false), pncmesh(NULL) {}
37 
38  MPI_Comm MyComm;
39  int NRanks, MyRank;
40 
41  struct Vert3
42  {
43  int v[3];
44  Vert3() { }
45  Vert3(int v0, int v1, int v2) { v[0] = v0; v[1] = v1; v[2] = v2; }
46  void Set(int v0, int v1, int v2) { v[0] = v0; v[1] = v1; v[2] = v2; }
47  void Set(const int *w) { v[0] = w[0]; v[1] = w[1]; v[2] = w[2]; }
48  };
49 
50  struct Vert4
51  {
52  int v[4];
53  Vert4() { }
54  Vert4(int v0, int v1, int v2, int v3)
55  { v[0] = v0; v[1] = v1; v[2] = v2; v[3] = v3; }
56  void Set(int v0, int v1, int v2, int v3)
57  { v[0] = v0; v[1] = v1; v[2] = v2; v[3] = v3; }
58  void Set(const int *w)
59  { v[0] = w[0]; v[1] = w[1]; v[2] = w[2]; v[3] = w[3]; }
60  };
61 
63  // shared face id 'i' is:
64  // * triangle id 'i', if i < shared_trias.Size()
65  // * quad id 'i-shared_trias.Size()', otherwise
68 
69  /// Shared objects in each group.
72  Table group_stria; // contains shared triangle indices
73  Table group_squad; // contains shared quadrilateral indices
74 
75  /// Shared to local index mapping.
78  // sface ids: all triangles first, then all quads
80 
82 
83  // glob_elem_offset + local element number defines a global element numbering
85  void ComputeGlobalElementOffset() const;
86 
87  /// Create from a nonconforming mesh.
88  ParMesh(const ParNCMesh &pncmesh);
89 
90  // Convert the local 'meshgen' to a global one.
91  void ReduceMeshGen();
92 
93  // Determine sedge_ledge and sface_lface.
94  void FinalizeParTopo();
95 
96  // Mark all tets to ensure consistency across MPI tasks; also mark the
97  // shared and boundary triangle faces using the consistently marked tets.
98  virtual void MarkTetMeshForRefinement(DSTable &v_to_v);
99 
100  /// Return a number(0-1) identifying how the given edge has been split
101  int GetEdgeSplittings(Element *edge, const DSTable &v_to_v, int *middle);
102  /// Append codes identifying how the given face has been split to @a codes
103  void GetFaceSplittings(const int *fv, const HashTable<Hashed2> &v_to_v,
104  Array<unsigned> &codes);
105 
106  bool DecodeFaceSplittings(HashTable<Hashed2> &v_to_v, const int *v,
107  const Array<unsigned> &codes, int &pos);
108 
110  int i, IsoparametricTransformation *ElTr);
111 
113  FaceElementTransformations* FETr, Element::Type face_type,
114  Geometry::Type face_geom);
115 
116  /// Update the groups after triangle refinement
117  void RefineGroups(const DSTable &v_to_v, int *middle);
118 
119  /// Update the groups after tetrahedron refinement
120  void RefineGroups(int old_nv, const HashTable<Hashed2> &v_to_v);
121 
122  void UniformRefineGroups2D(int old_nv);
123 
124  // f2qf can be NULL if all faces are quads or there are no quad faces
125  void UniformRefineGroups3D(int old_nv, int old_nedges,
126  const DSTable &old_v_to_v,
127  const STable3D &old_faces,
128  Array<int> *f2qf);
129 
130  void ExchangeFaceNbrData(Table *gr_sface, int *s2l_face);
131 
132  /// Refine a mixed 2D mesh uniformly.
133  virtual void UniformRefinement2D();
134 
135  /// Refine a mixed 3D mesh uniformly.
136  virtual void UniformRefinement3D();
137 
138  virtual void NURBSUniformRefinement();
139 
140  /// This function is not public anymore. Use GeneralRefinement instead.
141  virtual void LocalRefinement(const Array<int> &marked_el, int type = 3);
142 
143  /// This function is not public anymore. Use GeneralRefinement instead.
144  virtual void NonconformingRefinement(const Array<Refinement> &refinements,
145  int nc_limit = 0);
146 
147  virtual bool NonconformingDerefinement(Array<double> &elem_error,
148  double threshold, int nc_limit = 0,
149  int op = 1);
150 
151  void RebalanceImpl(const Array<int> *partition);
152 
153  void DeleteFaceNbrData();
154 
155  bool WantSkipSharedMaster(const NCMesh::Master &master) const;
156 
157  /// Fills out partitioned Mesh::vertices
158  int BuildLocalVertices(const Mesh& global_mesh, const int *partitioning,
159  Array<int> &vert_global_local);
160 
161  /// Fills out partitioned Mesh::elements
162  int BuildLocalElements(const Mesh& global_mesh, const int *partitioning,
163  const Array<int> &vert_global_local);
164 
165  /// Fills out partitioned Mesh::boundary
166  int BuildLocalBoundary(const Mesh& global_mesh, const int *partitioning,
167  const Array<int> &vert_global_local,
168  Array<bool>& activeBdrElem,
169  Table* &edge_element);
170 
171  void FindSharedFaces(const Mesh &mesh, const int* partition,
172  Array<int>& face_group,
173  ListOfIntegerSets& groups);
174 
175  int FindSharedEdges(const Mesh &mesh, const int* partition,
176  Table* &edge_element, ListOfIntegerSets& groups);
177 
178  int FindSharedVertices(const int *partition, Table* vertex_element,
179  ListOfIntegerSets& groups);
180 
181  void BuildFaceGroup(int ngroups, const Mesh &mesh,
182  const Array<int>& face_group,
183  int &nstria, int &nsquad);
184 
185  void BuildEdgeGroup(int ngroups, const Table& edge_element);
186 
187  void BuildVertexGroup(int ngroups, const Table& vert_element);
188 
189  void BuildSharedFaceElems(int ntri_faces, int nquad_faces,
190  const Mesh &mesh, int *partitioning,
191  const STable3D *faces_tbl,
192  const Array<int> &face_group,
193  const Array<int> &vert_global_local);
194 
195  void BuildSharedEdgeElems(int nedges, Mesh &mesh,
196  const Array<int> &vert_global_local,
197  const Table *edge_element);
198 
199  void BuildSharedVertMapping(int nvert, const Table* vert_element,
200  const Array<int> &vert_global_local);
201 
202  /// Ensure that bdr_attributes and attributes agree across processors
203  void DistributeAttributes(Array<int> &attr);
204 
205 public:
206  /** Copy constructor. Performs a deep copy of (almost) all data, so that the
207  source mesh can be modified (e.g. deleted, refined) without affecting the
208  new mesh. If 'copy_nodes' is false, use a shallow (pointer) copy for the
209  nodes, if present. */
210  explicit ParMesh(const ParMesh &pmesh, bool copy_nodes = true);
211 
212  ParMesh(MPI_Comm comm, Mesh &mesh, int *partitioning_ = NULL,
213  int part_method = 1);
214 
215  /// Read a parallel mesh, each MPI rank from its own file/stream.
216  /** The @a refine parameter is passed to the method Mesh::Finalize(). */
217  ParMesh(MPI_Comm comm, std::istream &input, bool refine = true);
218 
219  /// Create a uniformly refined (by any factor) version of @a orig_mesh.
220  /** @param[in] orig_mesh The starting coarse mesh.
221  @param[in] ref_factor The refinement factor, an integer > 1.
222  @param[in] ref_type Specify the positions of the new vertices. The
223  options are BasisType::ClosedUniform or
224  BasisType::GaussLobatto.
225 
226  The refinement data which can be accessed with GetRefinementTransforms()
227  is set to reflect the performed refinements.
228 
229  @note The constructed ParMesh is linear, i.e. it does not have nodes. */
230  ParMesh(ParMesh *orig_mesh, int ref_factor, int ref_type);
231 
232  virtual void Finalize(bool refine = false, bool fix_orientation = false);
233 
234  virtual void SetAttributes();
235 
236  MPI_Comm GetComm() const { return MyComm; }
237  int GetNRanks() const { return NRanks; }
238  int GetMyRank() const { return MyRank; }
239 
240  /** Map a global element number to a local element number. If the global
241  element is not on this processor, return -1. */
242  int GetLocalElementNum(long global_element_num) const;
243 
244  /// Map a local element number to a global element number.
245  long GetGlobalElementNum(int local_element_num) const;
246 
248 
249  // Face-neighbor elements and vertices
256  // Local face-neighbor elements and vertices ordered by face-neighbor
259 
261 
262  int GetNGroups() const { return gtopo.NGroups(); }
263 
264  ///@{ @name These methods require group > 0
265  int GroupNVertices(int group) { return group_svert.RowSize(group-1); }
266  int GroupNEdges(int group) { return group_sedge.RowSize(group-1); }
267  int GroupNTriangles(int group) { return group_stria.RowSize(group-1); }
268  int GroupNQuadrilaterals(int group) { return group_squad.RowSize(group-1); }
269 
270  int GroupVertex(int group, int i)
271  { return svert_lvert[group_svert.GetRow(group-1)[i]]; }
272  void GroupEdge(int group, int i, int &edge, int &o);
273  void GroupTriangle(int group, int i, int &face, int &o);
274  void GroupQuadrilateral(int group, int i, int &face, int &o);
275  ///@}
276 
277  void GenerateOffsets(int N, HYPRE_Int loc_sizes[],
278  Array<HYPRE_Int> *offsets[]) const;
279 
280  void ExchangeFaceNbrData();
281  void ExchangeFaceNbrNodes();
282 
283  virtual void SetCurvature(int order, bool discont = false, int space_dim = -1,
284  int ordering = 1);
285 
286  int GetNFaceNeighbors() const { return face_nbr_group.Size(); }
287  int GetFaceNbrGroup(int fn) const { return face_nbr_group[fn]; }
288  int GetFaceNbrRank(int fn) const;
289 
290  /** Similar to Mesh::GetFaceToElementTable with added face-neighbor elements
291  with indices offset by the local number of elements. */
293 
294  /** Get the FaceElementTransformations for the given shared face (edge 2D).
295  In the returned object, 1 and 2 refer to the local and the neighbor
296  elements, respectively. */
298  GetSharedFaceTransformations(int sf, bool fill2 = true);
299 
302  {
304 
305  return &FaceNbrTransformation;
306  }
307 
308  /// Return the number of shared faces (3D), edges (2D), vertices (1D)
309  int GetNSharedFaces() const;
310 
311  /// Return the local face index for the given shared face.
312  int GetSharedFace(int sface) const;
313 
314  /// See the remarks for the serial version in mesh.hpp
315  virtual void ReorientTetMesh();
316 
317  /// Utility function: sum integers from all processors (Allreduce).
318  virtual long ReduceInt(int value) const;
319 
320  /** Load balance the mesh by equipartitioning the global space-filling
321  sequence of elements. Works for nonconforming meshes only. */
322  void Rebalance();
323 
324  /** Load balance a nonconforming mesh using a user-defined partition.
325  Each local element 'i' is migrated to processor rank 'partition[i]',
326  for 0 <= i < GetNE(). */
327  void Rebalance(const Array<int> &partition);
328 
329  /** Print the part of the mesh in the calling processor adding the interface
330  as boundary (for visualization purposes) using the mfem v1.0 format. */
331  virtual void Print(std::ostream &out = mfem::out) const;
332 
333 #ifdef MFEM_USE_ADIOS2
334  /** Print the part of the mesh in the calling processor using adios2 bp
335  format. */
336  virtual void Print(adios2stream &out) const;
337 #endif
338 
339  /** Print the part of the mesh in the calling processor adding the interface
340  as boundary (for visualization purposes) using Netgen/Truegrid format .*/
341  virtual void PrintXG(std::ostream &out = mfem::out) const;
342 
343  /** Write the mesh to the stream 'out' on Process 0 in a form suitable for
344  visualization: the mesh is written as a disjoint mesh and the shared
345  boundary is added to the actual boundary; both the element and boundary
346  attributes are set to the processor number. */
347  void PrintAsOne(std::ostream &out = mfem::out);
348 
349  /// Old mesh format (Netgen/Truegrid) version of 'PrintAsOne'
350  void PrintAsOneXG(std::ostream &out = mfem::out);
351 
352  /// Returns the minimum and maximum corners of the mesh bounding box. For
353  /// high-order meshes, the geometry is refined first "ref" times.
354  void GetBoundingBox(Vector &p_min, Vector &p_max, int ref = 2);
355 
356  void GetCharacteristics(double &h_min, double &h_max,
357  double &kappa_min, double &kappa_max);
358 
359  /// Print various parallel mesh stats
360  virtual void PrintInfo(std::ostream &out = mfem::out);
361 
362  /// Save the mesh in a parallel mesh format.
363  void ParPrint(std::ostream &out) const;
364 
365  /** Print the mesh in parallel PVTU format. The PVTU and VTU files will be
366  stored in the directory specified by @a pathname. If the directory does
367  not exist, it will be created. */
368  virtual void PrintVTU(std::string pathname,
370  bool high_order_output=false,
371  int compression_level=0,
372  bool bdr=false);
373 
374  virtual int FindPoints(DenseMatrix& point_mat, Array<int>& elem_ids,
375  Array<IntegrationPoint>& ips, bool warn = true,
376  InverseElementTransformation *inv_trans = NULL);
377 
378  /// Debugging method
379  void PrintSharedEntities(const char *fname_prefix) const;
380 
381  virtual ~ParMesh();
382 
383  friend class ParNCMesh;
384 #ifdef MFEM_USE_PUMI
385  friend class ParPumiMesh;
386 #endif
387 
388 #ifdef MFEM_USE_ADIOS2
389  friend class adios2stream;
390 #endif
391 };
392 
393 }
394 
395 #endif // MFEM_USE_MPI
396 
397 #endif
int GetNFaceNeighbors() const
Definition: pmesh.hpp:286
void PrintAsOneXG(std::ostream &out=mfem::out)
Old mesh format (Netgen/Truegrid) version of &#39;PrintAsOne&#39;.
Definition: pmesh.cpp:4614
int Size() const
Return the logical size of the array.
Definition: array.hpp:124
void GetGhostFaceTransformation(FaceElementTransformations *FETr, Element::Type face_type, Geometry::Type face_geom)
Definition: pmesh.cpp:2368
virtual ~ParMesh()
Definition: pmesh.cpp:5571
long glob_offset_sequence
Definition: pmesh.hpp:84
void UniformRefineGroups2D(int old_nv)
Definition: pmesh.cpp:3694
int NRanks
Definition: pmesh.hpp:39
int FindSharedVertices(const int *partition, Table *vertex_element, ListOfIntegerSets &groups)
Definition: pmesh.cpp:568
ElementTransformation * GetFaceNbrElementTransformation(int i)
Definition: pmesh.hpp:301
void Set(const int *w)
Definition: pmesh.hpp:58
virtual void SetCurvature(int order, bool discont=false, int space_dim=-1, int ordering=1)
Definition: pmesh.cpp:1759
int BuildLocalVertices(const Mesh &global_mesh, const int *partitioning, Array< int > &vert_global_local)
Fills out partitioned Mesh::vertices.
Definition: pmesh.cpp:285
int GetNGroups() const
Definition: pmesh.hpp:262
virtual void ReorientTetMesh()
See the remarks for the serial version in mesh.hpp.
Definition: pmesh.cpp:2575
Array< Element * > face_nbr_elements
Definition: pmesh.hpp:254
int GetNSharedFaces() const
Return the number of shared faces (3D), edges (2D), vertices (1D)
Definition: pmesh.cpp:2508
Array< int > sface_lface
Definition: pmesh.hpp:79
void BuildSharedFaceElems(int ntri_faces, int nquad_faces, const Mesh &mesh, int *partitioning, const STable3D *faces_tbl, const Array< int > &face_group, const Array< int > &vert_global_local)
Definition: pmesh.cpp:703
A specialized ElementTransformation class representing a face and its two neighboring elements...
Definition: eltrans.hpp:467
virtual void LocalRefinement(const Array< int > &marked_el, int type=3)
This function is not public anymore. Use GeneralRefinement instead.
Definition: pmesh.cpp:2765
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
void GetRow(int i, Array< int > &row) const
Return row i in array row (the Table must be finalized)
Definition: table.cpp:191
void GetCharacteristics(double &h_min, double &h_max, double &kappa_min, double &kappa_max)
Definition: pmesh.cpp:5160
long glob_elem_offset
Definition: pmesh.hpp:84
void UniformRefineGroups3D(int old_nv, int old_nedges, const DSTable &old_v_to_v, const STable3D &old_faces, Array< int > *f2qf)
Definition: pmesh.cpp:3745
bool have_face_nbr_data
Definition: pmesh.hpp:250
Array< Vert3 > shared_trias
Definition: pmesh.hpp:66
Array< int > face_nbr_vertices_offset
Definition: pmesh.hpp:253
Array< int > face_nbr_group
Definition: pmesh.hpp:251
virtual void SetAttributes()
Definition: pmesh.cpp:1380
The inverse transformation of a given ElementTransformation.
Definition: eltrans.hpp:172
void RebalanceImpl(const Array< int > *partition)
Definition: pmesh.cpp:3389
void BuildFaceGroup(int ngroups, const Mesh &mesh, const Array< int > &face_group, int &nstria, int &nsquad)
Definition: pmesh.cpp:605
void Set(int v0, int v1, int v2)
Definition: pmesh.hpp:46
virtual void MarkTetMeshForRefinement(DSTable &v_to_v)
Definition: pmesh.cpp:1428
int GroupNQuadrilaterals(int group)
Definition: pmesh.hpp:268
void Set(const int *w)
Definition: pmesh.hpp:47
Array< int > sedge_ledge
Definition: pmesh.hpp:77
Table group_stria
Definition: pmesh.hpp:72
void GroupTriangle(int group, int i, int &face, int &o)
Definition: pmesh.cpp:1406
bool WantSkipSharedMaster(const NCMesh::Master &master) const
Definition: pmesh.cpp:4156
A parallel extension of the NCMesh class.
Definition: pncmesh.hpp:64
ParNCMesh * pncmesh
Definition: pmesh.hpp:260
Class for PUMI parallel meshes.
Definition: pumi.hpp:70
int FindSharedEdges(const Mesh &mesh, const int *partition, Table *&edge_element, ListOfIntegerSets &groups)
Definition: pmesh.cpp:514
int GetNRanks() const
Definition: pmesh.hpp:237
long GetGlobalElementNum(int local_element_num) const
Map a local element number to a global element number.
Definition: pmesh.cpp:1332
int GroupVertex(int group, int i)
Definition: pmesh.hpp:270
virtual int FindPoints(DenseMatrix &point_mat, Array< int > &elem_ids, Array< IntegrationPoint > &ips, bool warn=true, InverseElementTransformation *inv_trans=NULL)
Find the ids of the elements that contain the given points, and their corresponding reference coordin...
Definition: pmesh.cpp:5435
void ExchangeFaceNbrData()
Definition: pmesh.cpp:1779
void BuildSharedVertMapping(int nvert, const Table *vert_element, const Array< int > &vert_global_local)
Definition: pmesh.cpp:817
Vert4(int v0, int v1, int v2, int v3)
Definition: pmesh.hpp:54
void ReduceMeshGen()
Definition: pmesh.cpp:862
virtual void UniformRefinement2D()
Refine a mixed 2D mesh uniformly.
Definition: pmesh.cpp:3894
int GetLocalElementNum(long global_element_num) const
Definition: pmesh.cpp:1324
void Rebalance()
Definition: pmesh.cpp:3379
int BuildLocalBoundary(const Mesh &global_mesh, const int *partitioning, const Array< int > &vert_global_local, Array< bool > &activeBdrElem, Table *&edge_element)
Fills out partitioned Mesh::boundary.
Definition: pmesh.cpp:366
Array< int > face_nbr_elements_offset
Definition: pmesh.hpp:252
MPI_Comm MyComm
Definition: pmesh.hpp:38
Symmetric 3D Table stored as an array of rows each of which has a stack of column, floor, number nodes. The number of the node is assigned by counting the nodes from zero as they are pushed into the table. Diagonals of any kind are not allowed so the row, column and floor must all be different for each node. Only one node is stored for all 6 symmetric entries that are indexable by unique triplets of row, column, and floor.
Definition: stable3d.hpp:34
void GetBoundingBox(Vector &p_min, Vector &p_max, int ref=2)
Definition: pmesh.cpp:5144
Array< Element * > shared_edges
Definition: pmesh.hpp:62
virtual void NonconformingRefinement(const Array< Refinement > &refinements, int nc_limit=0)
This function is not public anymore. Use GeneralRefinement instead.
Definition: pmesh.cpp:3269
int GetSharedFace(int sface) const
Return the local face index for the given shared face.
Definition: pmesh.cpp:2527
void BuildSharedEdgeElems(int nedges, Mesh &mesh, const Array< int > &vert_global_local, const Table *edge_element)
Definition: pmesh.cpp:780
void PrintAsOne(std::ostream &out=mfem::out)
Definition: pmesh.cpp:4305
FaceElementTransformations * GetSharedFaceTransformations(int sf, bool fill2=true)
Definition: pmesh.cpp:2399
Table send_face_nbr_vertices
Definition: pmesh.hpp:258
Type
Constants for the classes derived from Element.
Definition: element.hpp:41
VTKFormat
Definition: vtk.hpp:22
virtual void Print(std::ostream &out=mfem::out) const
Definition: pmesh.cpp:4169
void Set(int v0, int v1, int v2, int v3)
Definition: pmesh.hpp:56
int GroupNVertices(int group)
Definition: pmesh.hpp:265
Array< Vert4 > shared_quads
Definition: pmesh.hpp:67
int GetFaceNbrGroup(int fn) const
Definition: pmesh.hpp:287
void GroupQuadrilateral(int group, int i, int &face, int &o)
Definition: pmesh.cpp:1417
int GroupNTriangles(int group)
Definition: pmesh.hpp:267
void FindSharedFaces(const Mesh &mesh, const int *partition, Array< int > &face_group, ListOfIntegerSets &groups)
Definition: pmesh.cpp:487
int GetMyRank() const
Definition: pmesh.hpp:238
void RefineGroups(const DSTable &v_to_v, int *middle)
Update the groups after triangle refinement.
Definition: pmesh.cpp:3437
MPI_Comm GetComm() const
Definition: pmesh.hpp:236
virtual bool NonconformingDerefinement(Array< double > &elem_error, double threshold, int nc_limit=0, int op=1)
NC version of GeneralDerefinement.
Definition: pmesh.cpp:3321
Array< Vertex > face_nbr_vertices
Definition: pmesh.hpp:255
void PrintSharedEntities(const char *fname_prefix) const
Debugging method.
Definition: pmesh.cpp:5484
virtual long ReduceInt(int value) const
Utility function: sum integers from all processors (Allreduce).
Definition: pmesh.cpp:5285
void BuildVertexGroup(int ngroups, const Table &vert_element)
Definition: pmesh.cpp:677
void GetFaceNbrElementTransformation(int i, IsoparametricTransformation *ElTr)
Definition: pmesh.cpp:1685
void DistributeAttributes(Array< int > &attr)
Ensure that bdr_attributes and attributes agree across processors.
Definition: pmesh.cpp:1338
void BuildEdgeGroup(int ngroups, const Table &edge_element)
Definition: pmesh.cpp:651
void GetFaceSplittings(const int *fv, const HashTable< Hashed2 > &v_to_v, Array< unsigned > &codes)
Append codes identifying how the given face has been split to codes.
Definition: pmesh.cpp:1571
A standard isoparametric element transformation.
Definition: eltrans.hpp:348
Table group_sedge
Definition: pmesh.hpp:71
void ComputeGlobalElementOffset() const
Definition: pmesh.cpp:850
void ExchangeFaceNbrNodes()
Definition: pmesh.cpp:2233
int BuildLocalElements(const Mesh &global_mesh, const int *partitioning, const Array< int > &vert_global_local)
Fills out partitioned Mesh::elements.
Definition: pmesh.cpp:333
Table group_svert
Shared objects in each group.
Definition: pmesh.hpp:70
bool DecodeFaceSplittings(HashTable< Hashed2 > &v_to_v, const int *v, const Array< unsigned > &codes, int &pos)
Definition: pmesh.cpp:1606
virtual void UniformRefinement3D()
Refine a mixed 3D mesh uniformly.
Definition: pmesh.cpp:3918
int NGroups() const
Return the number of groups.
virtual void Finalize(bool refine=false, bool fix_orientation=false)
Finalize the construction of a general Mesh.
Definition: pmesh.cpp:1310
IsoparametricTransformation FaceNbrTransformation
Definition: pmesh.hpp:81
int GetEdgeSplittings(Element *edge, const DSTable &v_to_v, int *middle)
Return a number(0-1) identifying how the given edge has been split.
Definition: pmesh.cpp:1556
virtual void PrintXG(std::ostream &out=mfem::out) const
Definition: pmesh.cpp:3955
void GroupEdge(int group, int i, int &edge, int &o)
Definition: pmesh.cpp:1398
Vector data type.
Definition: vector.hpp:51
void GenerateOffsets(int N, HYPRE_Int loc_sizes[], Array< HYPRE_Int > *offsets[]) const
Definition: pmesh.cpp:1642
Array< int > svert_lvert
Shared to local index mapping.
Definition: pmesh.hpp:76
int MyRank
Definition: pmesh.hpp:39
int RowSize(int i) const
Definition: table.hpp:108
Table send_face_nbr_elements
Definition: pmesh.hpp:257
virtual void PrintInfo(std::ostream &out=mfem::out)
Print various parallel mesh stats.
Definition: pmesh.cpp:5173
List of integer sets.
Definition: sets.hpp:59
Table * GetFaceToAllElementTable() const
Definition: pmesh.cpp:2310
void FinalizeParTopo()
Definition: pmesh.cpp:868
void DeleteFaceNbrData()
Definition: pmesh.cpp:1738
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition: globals.hpp:66
GroupTopology gtopo
Definition: pmesh.hpp:247
void ParPrint(std::ostream &out) const
Save the mesh in a parallel mesh format.
Definition: pmesh.cpp:5292
Class for parallel meshes.
Definition: pmesh.hpp:32
Abstract data type element.
Definition: element.hpp:28
int GetFaceNbrRank(int fn) const
Definition: pmesh.cpp:2293
virtual void PrintVTU(std::string pathname, VTKFormat format=VTKFormat::ASCII, bool high_order_output=false, int compression_level=0, bool bdr=false)
Definition: pmesh.cpp:5368
Table group_squad
Definition: pmesh.hpp:73
Vert3(int v0, int v1, int v2)
Definition: pmesh.hpp:45
int GroupNEdges(int group)
Definition: pmesh.hpp:266
virtual void NURBSUniformRefinement()
Refine NURBS mesh.
Definition: pmesh.cpp:3947