MFEM  v4.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, Lawrence Livermore National Security, LLC. Produced at
2 // the Lawrence Livermore National Laboratory. LLNL-CODE-443211. All Rights
3 // reserved. See file COPYRIGHT for details.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability see http://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the GNU Lesser General Public License (as published by the Free
10 // Software Foundation) version 2.1 dated February 1999.
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 
81  /// Create from a nonconforming mesh.
82  ParMesh(const ParNCMesh &pncmesh);
83 
84  // Convert the local 'meshgen' to a global one.
85  void ReduceMeshGen();
86 
87  // Determine sedge_ledge and sface_lface.
88  void FinalizeParTopo();
89 
90  // Mark all tets to ensure consistency across MPI tasks; also mark the
91  // shared and boundary triangle faces using the consistently marked tets.
92  virtual void MarkTetMeshForRefinement(DSTable &v_to_v);
93 
94  /// Return a number(0-1) identifying how the given edge has been split
95  int GetEdgeSplittings(Element *edge, const DSTable &v_to_v, int *middle);
96  /// Append codes identifying how the given face has been split to @a codes
97  void GetFaceSplittings(const int *fv, const HashTable<Hashed2> &v_to_v,
98  Array<unsigned> &codes);
99 
100  bool DecodeFaceSplittings(HashTable<Hashed2> &v_to_v, const int *v,
101  const Array<unsigned> &codes, int &pos);
102 
104  int i, IsoparametricTransformation *ElTr);
105 
107  FaceElementTransformations* FETr, Element::Type face_type,
108  Geometry::Type face_geom);
109 
110  /// Update the groups after triangle refinement
111  void RefineGroups(const DSTable &v_to_v, int *middle);
112 
113  /// Update the groups after tetrahedron refinement
114  void RefineGroups(int old_nv, const HashTable<Hashed2> &v_to_v);
115 
116  void UniformRefineGroups2D(int old_nv);
117 
118  // f2qf can be NULL if all faces are quads or there are no quad faces
119  void UniformRefineGroups3D(int old_nv, int old_nedges,
120  const DSTable &old_v_to_v,
121  const STable3D &old_faces,
122  Array<int> *f2qf);
123 
124  void ExchangeFaceNbrData(Table *gr_sface, int *s2l_face);
125 
126  /// Refine a mixed 2D mesh uniformly.
127  virtual void UniformRefinement2D();
128 
129  /// Refine a mixed 3D mesh uniformly.
130  virtual void UniformRefinement3D();
131 
132  virtual void NURBSUniformRefinement();
133 
134  /// This function is not public anymore. Use GeneralRefinement instead.
135  virtual void LocalRefinement(const Array<int> &marked_el, int type = 3);
136 
137  /// This function is not public anymore. Use GeneralRefinement instead.
138  virtual void NonconformingRefinement(const Array<Refinement> &refinements,
139  int nc_limit = 0);
140 
141  virtual bool NonconformingDerefinement(Array<double> &elem_error,
142  double threshold, int nc_limit = 0,
143  int op = 1);
144  void DeleteFaceNbrData();
145 
146  bool WantSkipSharedMaster(const NCMesh::Master &master) const;
147 
148  /// Fills out partitioned Mesh::vertices
149  int BuildLocalVertices(const Mesh& global_mesh, const int *partitioning,
150  Array<int> &vert_global_local);
151 
152  /// Fills out partitioned Mesh::elements
153  int BuildLocalElements(const Mesh& global_mesh, const int *partitioning,
154  const Array<int> &vert_global_local);
155 
156  /// Fills out partitioned Mesh::boundary
157  int BuildLocalBoundary(const Mesh& global_mesh, const int *partitioning,
158  const Array<int> &vert_global_local,
159  Array<bool>& activeBdrElem,
160  Table* &edge_element);
161 
162  void FindSharedFaces(const Mesh &mesh, const int* partition,
163  Array<int>& face_group,
164  ListOfIntegerSets& groups);
165 
166  int FindSharedEdges(const Mesh &mesh, const int* partition,
167  Table* &edge_element, ListOfIntegerSets& groups);
168 
169  int FindSharedVertices(const int *partition, Table* vertex_element,
170  ListOfIntegerSets& groups);
171 
172  void BuildFaceGroup(int ngroups, const Mesh &mesh,
173  const Array<int>& face_group,
174  int &nstria, int &nsquad);
175 
176  void BuildEdgeGroup(int ngroups, const Table& edge_element);
177 
178  void BuildVertexGroup(int ngroups, const Table& vert_element);
179 
180  void BuildSharedFaceElems(int ntri_faces, int nquad_faces,
181  const Mesh &mesh, int *partitioning,
182  const STable3D *faces_tbl,
183  const Array<int> &face_group,
184  const Array<int> &vert_global_local);
185 
186  void BuildSharedEdgeElems(int nedges, Mesh &mesh,
187  const Array<int> &vert_global_local,
188  const Table *edge_element);
189 
190  void BuildSharedVertMapping(int nvert, const Table* vert_element,
191  const Array<int> &vert_global_local);
192 
193 
194 public:
195  /** Copy constructor. Performs a deep copy of (almost) all data, so that the
196  source mesh can be modified (e.g. deleted, refined) without affecting the
197  new mesh. If 'copy_nodes' is false, use a shallow (pointer) copy for the
198  nodes, if present. */
199  explicit ParMesh(const ParMesh &pmesh, bool copy_nodes = true);
200 
201  ParMesh(MPI_Comm comm, Mesh &mesh, int *partitioning_ = NULL,
202  int part_method = 1);
203 
204  /// Read a parallel mesh, each MPI rank from its own file/stream.
205  /** The @a refine parameter is passed to the method Mesh::Finalize(). */
206  ParMesh(MPI_Comm comm, std::istream &input, bool refine = true);
207 
208  /// Create a uniformly refined (by any factor) version of @a orig_mesh.
209  /** @param[in] orig_mesh The starting coarse mesh.
210  @param[in] ref_factor The refinement factor, an integer > 1.
211  @param[in] ref_type Specify the positions of the new vertices. The
212  options are BasisType::ClosedUniform or
213  BasisType::GaussLobatto.
214 
215  The refinement data which can be accessed with GetRefinementTransforms()
216  is set to reflect the performed refinements.
217 
218  @note The constructed ParMesh is linear, i.e. it does not have nodes. */
219  ParMesh(ParMesh *orig_mesh, int ref_factor, int ref_type);
220 
221  virtual void Finalize(bool refine = false, bool fix_orientation = false);
222 
223  MPI_Comm GetComm() const { return MyComm; }
224  int GetNRanks() const { return NRanks; }
225  int GetMyRank() const { return MyRank; }
226 
228 
229  // Face-neighbor elements and vertices
236  // Local face-neighbor elements and vertices ordered by face-neighbor
239 
241 
242  int GetNGroups() const { return gtopo.NGroups(); }
243 
244  ///@{ @name These methods require group > 0
245  int GroupNVertices(int group) { return group_svert.RowSize(group-1); }
246  int GroupNEdges(int group) { return group_sedge.RowSize(group-1); }
247  int GroupNTriangles(int group) { return group_stria.RowSize(group-1); }
248  int GroupNQuadrilaterals(int group) { return group_squad.RowSize(group-1); }
249 
250  int GroupVertex(int group, int i)
251  { return svert_lvert[group_svert.GetRow(group-1)[i]]; }
252  void GroupEdge(int group, int i, int &edge, int &o);
253  void GroupTriangle(int group, int i, int &face, int &o);
254  void GroupQuadrilateral(int group, int i, int &face, int &o);
255  ///@}
256 
257  void GenerateOffsets(int N, HYPRE_Int loc_sizes[],
258  Array<HYPRE_Int> *offsets[]) const;
259 
260  void ExchangeFaceNbrData();
261  void ExchangeFaceNbrNodes();
262 
263  int GetNFaceNeighbors() const { return face_nbr_group.Size(); }
264  int GetFaceNbrGroup(int fn) const { return face_nbr_group[fn]; }
265  int GetFaceNbrRank(int fn) const;
266 
267  /** Similar to Mesh::GetFaceToElementTable with added face-neighbor elements
268  with indices offset by the local number of elements. */
270 
271  /** Get the FaceElementTransformations for the given shared face (edge 2D).
272  In the returned object, 1 and 2 refer to the local and the neighbor
273  elements, respectively. */
275  GetSharedFaceTransformations(int sf, bool fill2 = true);
276 
277  /// Return the number of shared faces (3D), edges (2D), vertices (1D)
278  int GetNSharedFaces() const;
279 
280  /// Return the local face index for the given shared face.
281  int GetSharedFace(int sface) const;
282 
283  /// See the remarks for the serial version in mesh.hpp
284  virtual void ReorientTetMesh();
285 
286  /// Utility function: sum integers from all processors (Allreduce).
287  virtual long ReduceInt(int value) const;
288 
289  /// Load balance the mesh. NC meshes only.
290  void Rebalance();
291 
292  /** Print the part of the mesh in the calling processor adding the interface
293  as boundary (for visualization purposes) using the mfem v1.0 format. */
294  virtual void Print(std::ostream &out = mfem::out) const;
295 
296  /** Print the part of the mesh in the calling processor adding the interface
297  as boundary (for visualization purposes) using Netgen/Truegrid format .*/
298  virtual void PrintXG(std::ostream &out = mfem::out) const;
299 
300  /** Write the mesh to the stream 'out' on Process 0 in a form suitable for
301  visualization: the mesh is written as a disjoint mesh and the shared
302  boundary is added to the actual boundary; both the element and boundary
303  attributes are set to the processor number. */
304  void PrintAsOne(std::ostream &out = mfem::out);
305 
306  /// Old mesh format (Netgen/Truegrid) version of 'PrintAsOne'
307  void PrintAsOneXG(std::ostream &out = mfem::out);
308 
309  /// Returns the minimum and maximum corners of the mesh bounding box. For
310  /// high-order meshes, the geometry is refined first "ref" times.
311  void GetBoundingBox(Vector &p_min, Vector &p_max, int ref = 2);
312 
313  void GetCharacteristics(double &h_min, double &h_max,
314  double &kappa_min, double &kappa_max);
315 
316  /// Print various parallel mesh stats
317  virtual void PrintInfo(std::ostream &out = mfem::out);
318 
319  /// Save the mesh in a parallel mesh format.
320  void ParPrint(std::ostream &out) const;
321 
322  virtual int FindPoints(DenseMatrix& point_mat, Array<int>& elem_ids,
323  Array<IntegrationPoint>& ips, bool warn = true,
324  InverseElementTransformation *inv_trans = NULL);
325 
326  /// Debugging method
327  void PrintSharedEntities(const char *fname_prefix) const;
328 
329  virtual ~ParMesh();
330 
331  friend class ParNCMesh;
332 #ifdef MFEM_USE_PUMI
333  friend class ParPumiMesh;
334 #endif
335 };
336 
337 }
338 
339 #endif // MFEM_USE_MPI
340 
341 #endif
int GetNFaceNeighbors() const
Definition: pmesh.hpp:263
void PrintAsOneXG(std::ostream &out=mfem::out)
Old mesh format (Netgen/Truegrid) version of &#39;PrintAsOne&#39;.
Definition: pmesh.cpp:4291
int Size() const
Logical size of the array.
Definition: array.hpp:118
virtual ~ParMesh()
Definition: pmesh.cpp:5181
void UniformRefineGroups2D(int old_nv)
Definition: pmesh.cpp:3390
int NRanks
Definition: pmesh.hpp:39
int FindSharedVertices(const int *partition, Table *vertex_element, ListOfIntegerSets &groups)
Definition: pmesh.cpp:564
void Set(const int *w)
Definition: pmesh.hpp:58
int BuildLocalVertices(const Mesh &global_mesh, const int *partitioning, Array< int > &vert_global_local)
Fills out partitioned Mesh::vertices.
Definition: pmesh.cpp:281
int GetNGroups() const
Definition: pmesh.hpp:242
virtual void ReorientTetMesh()
See the remarks for the serial version in mesh.hpp.
Definition: pmesh.cpp:2402
Array< Element * > face_nbr_elements
Definition: pmesh.hpp:234
int GetNSharedFaces() const
Return the number of shared faces (3D), edges (2D), vertices (1D)
Definition: pmesh.cpp:2361
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:699
virtual void LocalRefinement(const Array< int > &marked_el, int type=3)
This function is not public anymore. Use GeneralRefinement instead.
Definition: pmesh.cpp:2471
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:4837
void UniformRefineGroups3D(int old_nv, int old_nedges, const DSTable &old_v_to_v, const STable3D &old_faces, Array< int > *f2qf)
Definition: pmesh.cpp:3441
bool have_face_nbr_data
Definition: pmesh.hpp:230
Array< Vert3 > shared_trias
Definition: pmesh.hpp:66
Array< int > face_nbr_vertices_offset
Definition: pmesh.hpp:233
Array< int > face_nbr_group
Definition: pmesh.hpp:231
The inverse transformation of a given ElementTransformation.
Definition: eltrans.hpp:109
void BuildFaceGroup(int ngroups, const Mesh &mesh, const Array< int > &face_group, int &nstria, int &nsquad)
Definition: pmesh.cpp:601
void Set(int v0, int v1, int v2)
Definition: pmesh.hpp:46
virtual void MarkTetMeshForRefinement(DSTable &v_to_v)
Definition: pmesh.cpp:1325
int GroupNQuadrilaterals(int group)
Definition: pmesh.hpp:248
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:1303
bool WantSkipSharedMaster(const NCMesh::Master &master) const
Definition: pmesh.cpp:3850
A parallel extension of the NCMesh class.
Definition: pncmesh.hpp:64
ParNCMesh * pncmesh
Definition: pmesh.hpp:240
Class for PUMI parallel meshes.
Definition: pumi.hpp:72
int FindSharedEdges(const Mesh &mesh, const int *partition, Table *&edge_element, ListOfIntegerSets &groups)
Definition: pmesh.cpp:510
int GetNRanks() const
Definition: pmesh.hpp:224
int GroupVertex(int group, int i)
Definition: pmesh.hpp:250
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:5045
void ExchangeFaceNbrData()
Definition: pmesh.cpp:1655
void BuildSharedVertMapping(int nvert, const Table *vert_element, const Array< int > &vert_global_local)
Definition: pmesh.cpp:813
Vert4(int v0, int v1, int v2, int v3)
Definition: pmesh.hpp:54
void ReduceMeshGen()
Definition: pmesh.cpp:844
virtual void UniformRefinement2D()
Refine a mixed 2D mesh uniformly.
Definition: pmesh.cpp:3590
void Rebalance()
Load balance the mesh. NC meshes only.
Definition: pmesh.cpp:3085
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:362
Array< int > face_nbr_elements_offset
Definition: pmesh.hpp:232
MPI_Comm MyComm
Definition: pmesh.hpp:38
Symmetric 3D Table.
Definition: stable3d.hpp:29
void GetBoundingBox(Vector &p_min, Vector &p_max, int ref=2)
Definition: pmesh.cpp:4821
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:2975
int GetSharedFace(int sface) const
Return the local face index for the given shared face.
Definition: pmesh.cpp:2380
void BuildSharedEdgeElems(int nedges, Mesh &mesh, const Array< int > &vert_global_local, const Table *edge_element)
Definition: pmesh.cpp:776
void PrintAsOne(std::ostream &out=mfem::out)
Definition: pmesh.cpp:3992
FaceElementTransformations * GetSharedFaceTransformations(int sf, bool fill2=true)
Definition: pmesh.cpp:2276
Table send_face_nbr_vertices
Definition: pmesh.hpp:238
Type
Constants for the classes derived from Element.
Definition: element.hpp:41
virtual void Print(std::ostream &out=mfem::out) const
Definition: pmesh.cpp:3863
void Set(int v0, int v1, int v2, int v3)
Definition: pmesh.hpp:56
int GroupNVertices(int group)
Definition: pmesh.hpp:245
Array< Vert4 > shared_quads
Definition: pmesh.hpp:67
int GetFaceNbrGroup(int fn) const
Definition: pmesh.hpp:264
void GroupQuadrilateral(int group, int i, int &face, int &o)
Definition: pmesh.cpp:1314
int GroupNTriangles(int group)
Definition: pmesh.hpp:247
void FindSharedFaces(const Mesh &mesh, const int *partition, Array< int > &face_group, ListOfIntegerSets &groups)
Definition: pmesh.cpp:483
int GetMyRank() const
Definition: pmesh.hpp:225
void RefineGroups(const DSTable &v_to_v, int *middle)
Update the groups after triangle refinement.
Definition: pmesh.cpp:3133
MPI_Comm GetComm() const
Definition: pmesh.hpp:223
virtual bool NonconformingDerefinement(Array< double > &elem_error, double threshold, int nc_limit=0, int op=1)
NC version of GeneralDerefinement.
Definition: pmesh.cpp:3027
Array< Vertex > face_nbr_vertices
Definition: pmesh.hpp:235
void PrintSharedEntities(const char *fname_prefix) const
Debugging method.
Definition: pmesh.cpp:5094
virtual long ReduceInt(int value) const
Utility function: sum integers from all processors (Allreduce).
Definition: pmesh.cpp:4962
void BuildVertexGroup(int ngroups, const Table &vert_element)
Definition: pmesh.cpp:673
void GetFaceNbrElementTransformation(int i, IsoparametricTransformation *ElTr)
Definition: pmesh.cpp:1582
void BuildEdgeGroup(int ngroups, const Table &edge_element)
Definition: pmesh.cpp:647
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:1468
Table group_sedge
Definition: pmesh.hpp:71
void ExchangeFaceNbrNodes()
Definition: pmesh.cpp:2109
int BuildLocalElements(const Mesh &global_mesh, const int *partitioning, const Array< int > &vert_global_local)
Fills out partitioned Mesh::elements.
Definition: pmesh.cpp:329
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:1503
ElementTransformation * GetGhostFaceTransformation(FaceElementTransformations *FETr, Element::Type face_type, Geometry::Type face_geom)
Definition: pmesh.cpp:2244
virtual void UniformRefinement3D()
Refine a mixed 3D mesh uniformly.
Definition: pmesh.cpp:3610
int NGroups() const
virtual void Finalize(bool refine=false, bool fix_orientation=false)
Finalize the construction of a general Mesh.
Definition: pmesh.cpp:1281
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:1453
virtual void PrintXG(std::ostream &out=mfem::out) const
Definition: pmesh.cpp:3649
void GroupEdge(int group, int i, int &edge, int &o)
Definition: pmesh.cpp:1295
Vector data type.
Definition: vector.hpp:48
void GenerateOffsets(int N, HYPRE_Int loc_sizes[], Array< HYPRE_Int > *offsets[]) const
Definition: pmesh.cpp:1539
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:237
virtual void PrintInfo(std::ostream &out=mfem::out)
Print various parallel mesh stats.
Definition: pmesh.cpp:4850
List of integer sets.
Definition: sets.hpp:54
Table * GetFaceToAllElementTable() const
Definition: pmesh.cpp:2186
void FinalizeParTopo()
Definition: pmesh.cpp:850
void DeleteFaceNbrData()
Definition: pmesh.cpp:1634
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:64
GroupTopology gtopo
Definition: pmesh.hpp:227
void ParPrint(std::ostream &out) const
Save the mesh in a parallel mesh format.
Definition: pmesh.cpp:4969
Class for parallel meshes.
Definition: pmesh.hpp:32
Abstract data type element.
Definition: element.hpp:28
int GetFaceNbrRank(int fn) const
Definition: pmesh.cpp:2169
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:246
virtual void NURBSUniformRefinement()
Refine NURBS mesh.
Definition: pmesh.cpp:3641