MFEM  v4.5.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
mesh.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2022, 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_MESH
13 #define MFEM_MESH
14 
15 #include "../config/config.hpp"
16 #include "../general/stable3d.hpp"
17 #include "../general/globals.hpp"
18 #include "triangle.hpp"
19 #include "tetrahedron.hpp"
20 #include "vertex.hpp"
21 #include "vtk.hpp"
22 #include "ncmesh.hpp"
23 #include "../fem/eltrans.hpp"
24 #include "../fem/coefficient.hpp"
25 #include "../general/zstr.hpp"
26 #ifdef MFEM_USE_ADIOS2
27 #include "../general/adios2stream.hpp"
28 #endif
29 #include <iostream>
30 
31 namespace mfem
32 {
33 
34 // Data type mesh
35 
36 class GeometricFactors;
37 class FaceGeometricFactors;
38 class KnotVector;
39 class NURBSExtension;
40 class FiniteElementSpace;
41 class GridFunction;
42 struct Refinement;
43 
44 /** An enum type to specify if interior or boundary faces are desired. */
45 enum class FaceType : bool {Interior, Boundary};
46 
47 #ifdef MFEM_USE_MPI
48 class ParMesh;
49 class ParNCMesh;
50 #endif
51 
52 class Mesh
53 {
54 #ifdef MFEM_USE_MPI
55  friend class ParMesh;
56  friend class ParNCMesh;
57 #endif
58  friend class NCMesh;
59  friend class NURBSExtension;
60 
61 #ifdef MFEM_USE_ADIOS2
62  friend class adios2stream;
63 #endif
64 
65 protected:
66  int Dim;
67  int spaceDim;
68 
71  /** These variables store the number of Interior and Boundary faces. Calling
72  fes->GetMesh()->GetNBE() doesn't return the expected value in 3D because
73  periodic meshes in 3D have some of their faces marked as boundary for
74  visualization purpose in GLVis. */
76 
77  int meshgen; // see MeshGenerator()
78  int mesh_geoms; // sum of (1 << geom) for all geom of all dimensions
79 
80  // Counter for Mesh transformations: refinement, derefinement, rebalancing.
81  // Used for checking during Update operations on objects depending on the
82  // Mesh, such as FiniteElementSpace, GridFunction, etc.
83  long sequence;
84 
86  // Vertices are only at the corners of elements, where you would expect them
87  // in the lowest-order mesh. In some cases, e.g. in a Mesh that defines the
88  // patch topology for a NURBS mesh (see LoadPatchTopo()) the vertices may be
89  // empty while NumOfVertices is positive.
93 
94  /** @brief This structure stores the low level information necessary to
95  interpret the configuration of elements on a specific face. This
96  information can be accessed using methods like GetFaceElements(),
97  GetFaceInfos(), FaceIsInterior(), etc.
98 
99  For accessing higher level deciphered information look at
100  Mesh::FaceInformation, and its accessor Mesh::GetFaceInformation().
101 
102  Each face contains information on the indices, local reference faces,
103  orientations, and potential nonconformity for the two neighboring
104  elements on a face.
105  Each face can either be an interior, boundary, or shared interior face.
106  Each interior face is shared by two elements referred as Elem1 and Elem2.
107  For boundary faces only the information on Elem1 is relevant.
108  Shared interior faces correspond to faces where Elem1 and Elem2 are
109  distributed on different MPI ranks.
110  Regarding conformity, three cases are distinguished, conforming faces,
111  nonconforming slave faces, and nonconforming master faces. Master and
112  slave referring to the coarse and fine elements respectively on a
113  nonconforming face.
114  Nonconforming slave faces always have the slave element as Elem1 and
115  the master element as Elem2. On the other side, nonconforming master
116  faces always have the master element as Elem1, and one of the slave
117  element as Elem2. Except for ghost nonconforming slave faces, where
118  Elem1 is the master side and Elem2 is the slave side.
119 
120  The indices of Elem1 and Elem2 can be indirectly extracted from
121  FaceInfo::Elem1No and FaceInfo::Elem2No, read the note below for special
122  cases on the index of Elem2.
123 
124  The local face identifiers are deciphered from FaceInfo::Elem1Inf and
125  FaceInfo::Elem2Inf through the formula: LocalFaceIndex = ElemInf/64,
126  the semantic of the computed local face identifier can be found in
127  fem/geom.cpp. The local face identifier corresponds to an index
128  in the Constants<Geometry>::Edges arrays for 2D element geometries, and
129  to an index in the Constants<Geometry>::FaceVert arrays for 3D element
130  geometries.
131 
132  The orientation of each element relative to a face is obtained through
133  the formula: Orientation = ElemInf%64, the semantic of the orientation
134  can also be found in fem/geom.cpp. The orientation corresponds to
135  an index in the Constants<Geometry>::Orient arrays, providing the
136  sequence of vertices identifying the orientation of an edge/face. By
137  convention the orientation of Elem1 is always set to 0, serving as the
138  reference orientation. The orientation of Elem2 relatively to Elem1 is
139  therefore determined just by using the orientation of Elem2. An important
140  special case is the one of nonconforming faces, the orientation should
141  be composed with the PointMatrix, which also contains orientation
142  information. A special treatment should be done for 2D, the orientation
143  in the PointMatrix is not included, therefore when applying the
144  PointMatrix transformation, the PointMatrix should be flipped, except for
145  shared nonconforming slave faces where the transformation can be applied
146  as is.
147 
148  Another special case is the case of shared nonconforming faces. Ghost
149  faces use a different design based on so called "ghost" faces.
150  Ghost faces, as their name suggest are very well hidden, and they
151  usually have a separate interface from "standard" faces.
152  */
153  struct FaceInfo
154  {
155  // Inf = 64 * LocalFaceIndex + FaceOrientation
157  int NCFace; /* -1 if this is a regular conforming/boundary face;
158  index into 'nc_faces_info' if >= 0. */
159  };
160  // NOTE: in NC meshes, master faces have Elem2No == -1. Slave faces on the
161  // other hand have Elem2No and Elem2Inf set to the master face's element and
162  // its local face number.
163  //
164  // A local face is one generated from a local element and has index i in
165  // faces_info such that i < GetNumFaces(). Also, Elem1No always refers to the
166  // element (slave or master, in the nonconforming case) that generated the
167  // face.
168  // Classification of a local (non-ghost) face based on its FaceInfo:
169  // - Elem2No >= 0 --> local interior face; can be either:
170  // - NCFace == -1 --> conforming face, or
171  // - NCFace >= 0 --> nonconforming slave face; Elem2No is the index of
172  // the master volume element; Elem2Inf%64 is 0, see the note in
173  // Mesh::GenerateNCFaceInfo().
174  // - Elem2No < 0 --> local "boundary" face; can be one of:
175  // - NCFace == -1 --> conforming face; can be either:
176  // - Elem2Inf < 0 --> true boundary face (no element on side 2)
177  // - Elem2Inf >= 0 --> shared face where element 2 is a face-neighbor
178  // element with index -1-Elem2No. This state is initialized by
179  // ParMesh::ExchangeFaceNbrData().
180  // - NCFace >= 0 --> nonconforming face; can be one of:
181  // - Elem2Inf < 0 --> master nonconforming face, interior or shared;
182  // In this case, Elem2No is -1; see GenerateNCFaceInfo().
183  // - Elem2Inf >= 0 --> shared slave nonconforming face where element 2
184  // is the master face-neighbor element with index -1-Elem2No; see
185  // ParNCMesh::GetFaceNeighbors().
186  //
187  // A ghost face is a nonconforming face that is generated by a non-local,
188  // i.e. ghost, element. A ghost face has index i in faces_info such that
189  // i >= GetNumFaces().
190  // Classification of a ghost (non-local) face based on its FaceInfo:
191  // - Elem1No == -1 --> master ghost face? These ghost faces also have:
192  // Elem2No == -1, Elem1Inf == Elem2Inf == -1, and NCFace == -1.
193  // - Elem1No >= 0 --> slave ghost face; Elem1No is the index of the local
194  // master side element, i.e. side 1 IS NOT the side that generated the
195  // face. Elem2No is < 0 and -1-Elem2No is the index of the ghost
196  // face-neighbor element that generated this slave ghost face. In this
197  // case, Elem2Inf >= 0 and NCFace >= 0.
198  // Relevant methods: GenerateFaces(), GenerateNCFaceInfo(),
199  // ParNCMesh::GetFaceNeighbors(),
200  // ParMesh::ExchangeFaceNbrData()
201 
202  struct NCFaceInfo
203  {
204  bool Slave; // true if this is a slave face, false if master face
205  int MasterFace; // if Slave, this is the index of the master face
206  // If not Slave, 'MasterFace' is the local face index of this master face
207  // as a face in the unique adjacent element.
208  const DenseMatrix* PointMatrix; // if Slave, position within master face
209  // (NOTE: PointMatrix points to a matrix owned by NCMesh.)
210 
211  NCFaceInfo() = default;
212 
213  NCFaceInfo(bool slave, int master, const DenseMatrix* pm)
214  : Slave(slave), MasterFace(master), PointMatrix(pm) {}
215  };
216 
219 
224  Table *bel_to_edge; // for 3D
226  mutable Table *face_edge;
227  mutable Table *edge_vertex;
228 
233 
234  // refinement embeddings for forward compatibility with NCMesh
236 
237  // Nodes are only active for higher order meshes, and share locations with
238  // the vertices, plus all the higher- order control points within the
239  // element and along the edges and on the faces.
242 
243  static const int vtk_quadratic_tet[10];
244  static const int vtk_quadratic_pyramid[13];
245  static const int vtk_quadratic_wedge[18];
246  static const int vtk_quadratic_hex[27];
247 
248 #ifdef MFEM_USE_MEMALLOC
249  friend class Tetrahedron;
251 #endif
252 
253  // used during NC mesh initialization only
255 
256 public:
264 
266 
267  /// A list of all unique element attributes used by the Mesh.
269  /// A list of all unique boundary attributes used by the Mesh.
271 
272  NURBSExtension *NURBSext; ///< Optional NURBS mesh extension.
273  NCMesh *ncmesh; ///< Optional nonconforming mesh extension.
274  Array<GeometricFactors*> geom_factors; ///< Optional geometric factors.
276  face_geom_factors; ///< Optional face geometric factors.
277 
278  // Global parameter that can be used to control the removal of unused
279  // vertices performed when reading a mesh in MFEM format. The default value
280  // (true) is set in mesh_readers.cpp.
282 
283 protected:
285 
286  void Init();
287  void InitTables();
288  void SetEmpty(); // Init all data members with empty values
289  void DestroyTables();
291  void DestroyPointers(); // Delete data specifically allocated by class Mesh.
292  void Destroy(); // Delete all owned data.
293  void ResetLazyData();
294 
295  Element *ReadElementWithoutAttr(std::istream &);
296  static void PrintElementWithoutAttr(const Element *, std::ostream &);
297 
298  Element *ReadElement(std::istream &);
299  static void PrintElement(const Element *, std::ostream &);
300 
301  // Readers for different mesh formats, used in the Load() method.
302  // The implementations of these methods are in mesh_readers.cpp.
303  void ReadMFEMMesh(std::istream &input, int version, int &curved);
304  void ReadLineMesh(std::istream &input);
305  void ReadNetgen2DMesh(std::istream &input, int &curved);
306  void ReadNetgen3DMesh(std::istream &input);
307  void ReadTrueGridMesh(std::istream &input);
308  void CreateVTKMesh(const Vector &points, const Array<int> &cell_data,
309  const Array<int> &cell_offsets,
310  const Array<int> &cell_types,
311  const Array<int> &cell_attributes,
312  int &curved, int &read_gf, bool &finalize_topo);
313  void ReadVTKMesh(std::istream &input, int &curved, int &read_gf,
314  bool &finalize_topo);
315  void ReadXML_VTKMesh(std::istream &input, int &curved, int &read_gf,
316  bool &finalize_topo, const std::string &xml_prefix="");
317  void ReadNURBSMesh(std::istream &input, int &curved, int &read_gf);
318  void ReadInlineMesh(std::istream &input, bool generate_edges = false);
319  void ReadGmshMesh(std::istream &input, int &curved, int &read_gf);
320  /* Note NetCDF (optional library) is used for reading cubit files */
321 #ifdef MFEM_USE_NETCDF
322  void ReadCubit(const char *filename, int &curved, int &read_gf);
323 #endif
324 
325  /// Determine the mesh generator bitmask #meshgen, see MeshGenerator().
326  /** Also, initializes #mesh_geoms. */
327  void SetMeshGen();
328 
329  /// Return the length of the segment from node i to node j.
330  double GetLength(int i, int j) const;
331 
332  /** Compute the Jacobian of the transformation from the perfect
333  reference element at the center of the element. */
334  void GetElementJacobian(int i, DenseMatrix &J);
335 
336  void MarkForRefinement();
338  void GetEdgeOrdering(DSTable &v_to_v, Array<int> &order);
339  virtual void MarkTetMeshForRefinement(DSTable &v_to_v);
340 
341  // Methods used to prepare and apply permutation of the mesh nodes assuming
342  // that the mesh elements may be rotated (e.g. to mark triangle or tet edges
343  // for refinement) between the two calls - PrepareNodeReorder() and
344  // DoNodeReorder(). The latter method assumes that the 'faces' have not been
345  // updated after the element rotations.
346  void PrepareNodeReorder(DSTable **old_v_to_v, Table **old_elem_vert);
347  void DoNodeReorder(DSTable *old_v_to_v, Table *old_elem_vert);
348 
350  STable3D *GetElementToFaceTable(int ret_ftbl = 0);
351 
352  /** Red refinement. Element with index i is refined. The default
353  red refinement for now is Uniform. */
354  void RedRefinement(int i, const DSTable &v_to_v,
355  int *edge1, int *edge2, int *middle)
356  { UniformRefinement(i, v_to_v, edge1, edge2, middle); }
357 
358  /** Green refinement. Element with index i is refined. The default
359  refinement for now is Bisection. */
360  void GreenRefinement(int i, const DSTable &v_to_v,
361  int *edge1, int *edge2, int *middle)
362  { Bisection(i, v_to_v, edge1, edge2, middle); }
363 
364  /// Bisect a triangle: element with index @a i is bisected.
365  void Bisection(int i, const DSTable &, int *, int *, int *);
366 
367  /// Bisect a tetrahedron: element with index @a i is bisected.
368  void Bisection(int i, HashTable<Hashed2> &);
369 
370  /// Bisect a boundary triangle: boundary element with index @a i is bisected.
371  void BdrBisection(int i, const HashTable<Hashed2> &);
372 
373  /** Uniform Refinement. Element with index i is refined uniformly. */
374  void UniformRefinement(int i, const DSTable &, int *, int *, int *);
375 
376  /** @brief Averages the vertices with given @a indexes and saves the result
377  in #vertices[result]. */
378  void AverageVertices(const int *indexes, int n, int result);
379 
381  int FindCoarseElement(int i);
382 
383  /// Update the nodes of a curved mesh after refinement
384  void UpdateNodes();
385 
386  /// Helper to set vertex coordinates given a high-order curvature function.
387  void SetVerticesFromNodes(const GridFunction *nodes);
388 
389  void UniformRefinement2D_base(bool update_nodes = true);
390 
391  /// Refine a mixed 2D mesh uniformly.
393 
394  /* If @a f2qf is not NULL, adds all quadrilateral faces to @a f2qf which
395  represents a "face-to-quad-face" index map. When all faces are quads, the
396  array @a f2qf is kept empty since it is not needed. */
397  void UniformRefinement3D_base(Array<int> *f2qf = NULL,
398  DSTable *v_to_v_p = NULL,
399  bool update_nodes = true);
400 
401  /// Refine a mixed 3D mesh uniformly.
403 
404  /// Refine NURBS mesh.
405  virtual void NURBSUniformRefinement();
406 
407  /// This function is not public anymore. Use GeneralRefinement instead.
408  virtual void LocalRefinement(const Array<int> &marked_el, int type = 3);
409 
410  /// This function is not public anymore. Use GeneralRefinement instead.
411  virtual void NonconformingRefinement(const Array<Refinement> &refinements,
412  int nc_limit = 0);
413 
414  /// NC version of GeneralDerefinement.
415  virtual bool NonconformingDerefinement(Array<double> &elem_error,
416  double threshold, int nc_limit = 0,
417  int op = 1);
418  /// Derefinement helper.
419  double AggregateError(const Array<double> &elem_error,
420  const int *fine, int nfine, int op);
421 
422  /// Read NURBS patch/macro-element mesh
423  void LoadPatchTopo(std::istream &input, Array<int> &edge_to_knot);
424 
425  void UpdateNURBS();
426 
427  void PrintTopo(std::ostream &out, const Array<int> &e_to_k) const;
428 
429  /// Used in GetFaceElementTransformations (...)
432  int i);
434  int i);
435  /// Used in GetFaceElementTransformations (...)
437  int i);
438  /// Used in GetFaceElementTransformations (...)
440  int i);
441  /// Used in GetFaceElementTransformations (...)
443  int i);
444  /// Used in GetFaceElementTransformations (...)
446  int i);
447  /// Used in GetFaceElementTransformations (...)
449  int i);
450  /// Used in GetFaceElementTransformations (...)
452  int i);
453 
454  /** Used in GetFaceElementTransformations to account for the fact that a
455  slave face occupies only a portion of its master face. */
457  const FaceInfo &fi, bool is_ghost);
458 
459  bool IsSlaveFace(const FaceInfo &fi) const;
460 
461  /// Returns the orientation of "test" relative to "base"
462  static int GetTriOrientation (const int * base, const int * test);
463  /// Returns the orientation of "test" relative to "base"
464  static int GetQuadOrientation (const int * base, const int * test);
465  /// Returns the orientation of "test" relative to "base"
466  static int GetTetOrientation (const int * base, const int * test);
467 
468  static void GetElementArrayEdgeTable(const Array<Element*> &elem_array,
469  const DSTable &v_to_v,
470  Table &el_to_edge);
471 
472  /** Return element to edge table and the indices for the boundary edges.
473  The entries in the table are ordered according to the order of the
474  nodes in the elements. For example, if T is the element to edge table
475  T(i, 0) gives the index of edge in element i that connects vertex 0
476  to vertex 1, etc. Returns the number of the edges. */
478 
479  /// Used in GenerateFaces()
480  void AddPointFaceElement(int lf, int gf, int el);
481 
482  void AddSegmentFaceElement (int lf, int gf, int el, int v0, int v1);
483 
484  void AddTriangleFaceElement (int lf, int gf, int el,
485  int v0, int v1, int v2);
486 
487  void AddQuadFaceElement (int lf, int gf, int el,
488  int v0, int v1, int v2, int v3);
489  /** For a serial Mesh, return true if the face is interior. For a parallel
490  ParMesh return true if the face is interior or shared. In parallel, this
491  method only works if the face neighbor data is exchanged. */
492  bool FaceIsTrueInterior(int FaceNo) const
493  {
494  return FaceIsInterior(FaceNo) || (faces_info[FaceNo].Elem2Inf >= 0);
495  }
496 
497  void FreeElement(Element *E);
498 
499  void GenerateFaces();
500  void GenerateNCFaceInfo();
501 
502  /// Begin construction of a mesh
503  void InitMesh(int Dim_, int spaceDim_, int NVert, int NElem, int NBdrElem);
504 
505  // Used in the methods FinalizeXXXMesh() and FinalizeTopology()
506  void FinalizeCheck();
507 
508  void Loader(std::istream &input, int generate_edges = 0,
509  std::string parse_tag = "");
510 
511  // If NURBS mesh, write NURBS format. If NCMesh, write mfem v1.1 format.
512  // If section_delimiter is empty, write mfem v1.0 format. Otherwise, write
513  // mfem v1.2 format with the given section_delimiter at the end.
514  void Printer(std::ostream &out = mfem::out,
515  std::string section_delimiter = "") const;
516 
517  /** Creates mesh for the parallelepiped [0,sx]x[0,sy]x[0,sz], divided into
518  nx*ny*nz hexahedra if type=HEXAHEDRON or into 6*nx*ny*nz tetrahedrons if
519  type=TETRAHEDRON. The parameter @a sfc_ordering controls how the elements
520  (when type=HEXAHEDRON) are ordered: true - use space-filling curve
521  ordering, or false - use lexicographic ordering. */
522  void Make3D(int nx, int ny, int nz, Element::Type type,
523  double sx, double sy, double sz, bool sfc_ordering);
524 
525  /** Creates mesh for the rectangle [0,sx]x[0,sy], divided into nx*ny
526  quadrilaterals if type = QUADRILATERAL or into 2*nx*ny triangles if
527  type = TRIANGLE. If generate_edges = 0 (default) edges are not generated,
528  if 1 edges are generated. The parameter @a sfc_ordering controls how the
529  elements (when type=QUADRILATERAL) are ordered: true - use space-filling
530  curve ordering, or false - use lexicographic ordering. */
531  void Make2D(int nx, int ny, Element::Type type, double sx, double sy,
532  bool generate_edges, bool sfc_ordering);
533 
534  /// Creates a 1D mesh for the interval [0,sx] divided into n equal intervals.
535  void Make1D(int n, double sx = 1.0);
536 
537  /// Internal function used in Mesh::MakeRefined
538  void MakeRefined_(Mesh &orig_mesh, const Array<int> ref_factors,
539  int ref_type);
540 
541  /// Initialize vertices/elements/boundary/tables from a nonconforming mesh.
542  void InitFromNCMesh(const NCMesh &ncmesh);
543 
544  /// Create from a nonconforming mesh.
545  explicit Mesh(const NCMesh &ncmesh);
546 
547  // used in GetElementData() and GetBdrElementData()
548  void GetElementData(const Array<Element*> &elem_array, int geom,
549  Array<int> &elem_vtx, Array<int> &attr) const;
550 
551  double GetElementSize(ElementTransformation *T, int type = 0);
552 
553  // Internal helper used in MakeSimplicial (and ParMesh::MakeSimplicial).
554  void MakeSimplicial_(const Mesh &orig_mesh, int *vglobal);
555 
556 public:
557 
558  Mesh() { SetEmpty(); }
559 
560  /** Copy constructor. Performs a deep copy of (almost) all data, so that the
561  source mesh can be modified (e.g. deleted, refined) without affecting the
562  new mesh. If 'copy_nodes' is false, use a shallow (pointer) copy for the
563  nodes, if present. */
564  explicit Mesh(const Mesh &mesh, bool copy_nodes = true);
565 
566  /// Move constructor, useful for using a Mesh as a function return value.
567  Mesh(Mesh &&mesh);
568 
569  /// Move assignment operstor.
570  Mesh& operator=(Mesh &&mesh);
571 
572  /// Explicitly delete the copy assignment operator.
573  Mesh& operator=(const Mesh &mesh) = delete;
574 
575  /** @name Named mesh constructors.
576 
577  Each of these constructors uses the move constructor, and can be used as
578  the right-hand side of an assignment when creating new meshes. */
579  ///@{
580 
581  /** Creates mesh by reading a file in MFEM, Netgen, or VTK format. If
582  generate_edges = 0 (default) edges are not generated, if 1 edges are
583  generated. */
584  static Mesh LoadFromFile(const char *filename,
585  int generate_edges = 0, int refine = 1,
586  bool fix_orientation = true);
587 
588  /** Creates 1D mesh , divided into n equal intervals. */
589  static Mesh MakeCartesian1D(int n, double sx = 1.0);
590 
591  /** Creates mesh for the rectangle [0,sx]x[0,sy], divided into nx*ny
592  quadrilaterals if type = QUADRILATERAL or into 2*nx*ny triangles if
593  type = TRIANGLE. If generate_edges = 0 (default) edges are not generated,
594  if 1 edges are generated. If scf_ordering = true (default), elements are
595  ordered along a space-filling curve, instead of row by row. */
596  static Mesh MakeCartesian2D(
597  int nx, int ny, Element::Type type, bool generate_edges = false,
598  double sx = 1.0, double sy = 1.0, bool sfc_ordering = true);
599 
600  /** Creates mesh for the parallelepiped [0,sx]x[0,sy]x[0,sz], divided into
601  nx*ny*nz hexahedra if type=HEXAHEDRON or into 6*nx*ny*nz tetrahedrons if
602  type=TETRAHEDRON. If sfc_ordering = true (default), elements are ordered
603  along a space-filling curve, instead of row by row and layer by layer. */
604  static Mesh MakeCartesian3D(
605  int nx, int ny, int nz, Element::Type type,
606  double sx = 1.0, double sy = 1.0, double sz = 1.0,
607  bool sfc_ordering = true);
608 
609  /// Create a refined (by any factor) version of @a orig_mesh.
610  /** @param[in] orig_mesh The starting coarse mesh.
611  @param[in] ref_factor The refinement factor, an integer > 1.
612  @param[in] ref_type Specify the positions of the new vertices. The
613  options are BasisType::ClosedUniform or
614  BasisType::GaussLobatto.
615 
616  The refinement data which can be accessed with GetRefinementTransforms()
617  is set to reflect the performed refinements.
618 
619  @note The constructed Mesh is straight-sided. */
620  static Mesh MakeRefined(Mesh &orig_mesh, int ref_factor, int ref_type);
621 
622  /// Create a refined mesh, where each element of the original mesh may be
623  /// refined by a different factor.
624  /** @param[in] orig_mesh The starting coarse mesh.
625  @param[in] ref_factors An array of integers whose size is the number of
626  elements of @a orig_mesh. The @a ith element of
627  @a orig_mesh is refined by refinement factor
628  @a ref_factors[i].
629  @param[in] ref_type Specify the positions of the new vertices. The
630  options are BasisType::ClosedUniform or
631  BasisType::GaussLobatto.
632 
633  The refinement data which can be accessed with GetRefinementTransforms()
634  is set to reflect the performed refinements.
635 
636  @note The constructed Mesh is straight-sided. */
637  /// refined @a ref_factors[i] times in each dimension.
638  static Mesh MakeRefined(Mesh &orig_mesh, const Array<int> &ref_factors,
639  int ref_type);
640 
641  /** Create a mesh by splitting each element of @a orig_mesh into simplices.
642  Quadrilaterals are split into two triangles, prisms are split into
643  3 tetrahedra, and hexahedra are split into either 5 or 6 tetrahedra
644  depending on the configuration.
645  @warning The curvature of the original mesh is not carried over to the
646  new mesh. Periodic meshes are not supported. */
647  static Mesh MakeSimplicial(const Mesh &orig_mesh);
648 
649  /// Create a periodic mesh by identifying vertices of @a orig_mesh.
650  /** Each vertex @a i will be mapped to vertex @a v2v[i], such that all
651  vertices that are coincident under the periodic mapping get mapped to
652  the same index. The mapping @a v2v can be generated from translation
653  vectors using Mesh::CreatePeriodicVertexMapping.
654  @note MFEM requires that each edge of the resulting mesh be uniquely
655  identifiable by a pair of distinct vertices. As a consequence, periodic
656  boundaries must be connected by at least three edges. */
657  static Mesh MakePeriodic(const Mesh &orig_mesh, const std::vector<int> &v2v);
658 
659  ///@}
660 
661  /// @brief Creates a mapping @a v2v from the vertex indices of the mesh such
662  /// that coincident vertices under the given @a translations are identified.
663  /** Each Vector in @a translations should be of size @a sdim (the spatial
664  dimension of the mesh). Two vertices are considered coincident if the
665  translated coordinates of one vertex are within the given tolerance (@a
666  tol, relative to the mesh diameter) of the coordinates of the other
667  vertex.
668  @warning This algorithm does not scale well with the number of boundary
669  vertices in the mesh, and may run slowly on very large meshes. */
670  std::vector<int> CreatePeriodicVertexMapping(
671  const std::vector<Vector> &translations, double tol = 1e-8) const;
672 
673  /// Construct a Mesh from the given primary data.
674  /** The array @a vertices is used as external data, i.e. the Mesh does not
675  copy the data and will not delete the pointer.
676 
677  The data from the other arrays is copied into the internal Mesh data
678  structures.
679 
680  This method calls the method FinalizeTopology(). The method Finalize()
681  may be called after this constructor and after optionally setting the
682  Mesh nodes. */
683  Mesh(double *vertices, int num_vertices,
684  int *element_indices, Geometry::Type element_type,
685  int *element_attributes, int num_elements,
686  int *boundary_indices, Geometry::Type boundary_type,
687  int *boundary_attributes, int num_boundary_elements,
688  int dimension, int space_dimension = -1);
689 
690  /** @anchor mfem_Mesh_init_ctor
691  @brief _Init_ constructor: begin the construction of a Mesh object. */
692  Mesh(int Dim_, int NVert, int NElem, int NBdrElem = 0, int spaceDim_ = -1)
693  {
694  if (spaceDim_ == -1) { spaceDim_ = Dim_; }
695  InitMesh(Dim_, spaceDim_, NVert, NElem, NBdrElem);
696  }
697 
698  /** @name Methods for Mesh construction.
699 
700  These methods are intended to be used with the @ref mfem_Mesh_init_ctor
701  "init constructor". */
702  ///@{
703 
704  Element *NewElement(int geom);
705 
706  int AddVertex(double x, double y = 0.0, double z = 0.0);
707  int AddVertex(const double *coords);
708  /// Mark vertex @a i as nonconforming, with parent vertices @a p1 and @a p2.
709  void AddVertexParents(int i, int p1, int p2);
710 
711  int AddSegment(int v1, int v2, int attr = 1);
712  int AddSegment(const int *vi, int attr = 1);
713 
714  int AddTriangle(int v1, int v2, int v3, int attr = 1);
715  int AddTriangle(const int *vi, int attr = 1);
716  int AddTri(const int *vi, int attr = 1) { return AddTriangle(vi, attr); }
717 
718  int AddQuad(int v1, int v2, int v3, int v4, int attr = 1);
719  int AddQuad(const int *vi, int attr = 1);
720 
721  int AddTet(int v1, int v2, int v3, int v4, int attr = 1);
722  int AddTet(const int *vi, int attr = 1);
723 
724  int AddWedge(int v1, int v2, int v3, int v4, int v5, int v6, int attr = 1);
725  int AddWedge(const int *vi, int attr = 1);
726 
727  int AddPyramid(int v1, int v2, int v3, int v4, int v5, int attr = 1);
728  int AddPyramid(const int *vi, int attr = 1);
729 
730  int AddHex(int v1, int v2, int v3, int v4, int v5, int v6, int v7, int v8,
731  int attr = 1);
732  int AddHex(const int *vi, int attr = 1);
733  void AddHexAsTets(const int *vi, int attr = 1);
734  void AddHexAsWedges(const int *vi, int attr = 1);
735  void AddHexAsPyramids(const int *vi, int attr = 1);
736 
737  /// The parameter @a elem should be allocated using the NewElement() method
738  int AddElement(Element *elem);
739  int AddBdrElement(Element *elem);
740 
741  int AddBdrSegment(int v1, int v2, int attr = 1);
742  int AddBdrSegment(const int *vi, int attr = 1);
743 
744  int AddBdrTriangle(int v1, int v2, int v3, int attr = 1);
745  int AddBdrTriangle(const int *vi, int attr = 1);
746 
747  int AddBdrQuad(int v1, int v2, int v3, int v4, int attr = 1);
748  int AddBdrQuad(const int *vi, int attr = 1);
749  void AddBdrQuadAsTriangles(const int *vi, int attr = 1);
750 
751  int AddBdrPoint(int v, int attr = 1);
752 
754  /// Finalize the construction of a triangular Mesh.
755  void FinalizeTriMesh(int generate_edges = 0, int refine = 0,
756  bool fix_orientation = true);
757  /// Finalize the construction of a quadrilateral Mesh.
758  void FinalizeQuadMesh(int generate_edges = 0, int refine = 0,
759  bool fix_orientation = true);
760  /// Finalize the construction of a tetrahedral Mesh.
761  void FinalizeTetMesh(int generate_edges = 0, int refine = 0,
762  bool fix_orientation = true);
763  /// Finalize the construction of a wedge Mesh.
764  void FinalizeWedgeMesh(int generate_edges = 0, int refine = 0,
765  bool fix_orientation = true);
766  /// Finalize the construction of a hexahedral Mesh.
767  void FinalizeHexMesh(int generate_edges = 0, int refine = 0,
768  bool fix_orientation = true);
769  /// Finalize the construction of any type of Mesh.
770  /** This method calls FinalizeTopology() and Finalize(). */
771  void FinalizeMesh(int refine = 0, bool fix_orientation = true);
772 
773  ///@}
774 
775  /** @brief Finalize the construction of the secondary topology (connectivity)
776  data of a Mesh. */
777  /** This method does not require any actual coordinate data (either vertex
778  coordinates for linear meshes or node coordinates for meshes with nodes)
779  to be available. However, the data generated by this method is generally
780  required by the FiniteElementSpace class.
781 
782  After calling this method, setting the Mesh vertices or nodes, it may be
783  appropriate to call the method Finalize(). */
784  void FinalizeTopology(bool generate_bdr = true);
785 
786  /// Finalize the construction of a general Mesh.
787  /** This method will:
788  - check and optionally fix the orientation of regular elements
789  - check and fix the orientation of boundary elements
790  - assume that #vertices are defined, if #Nodes == NULL
791  - assume that #Nodes are defined, if #Nodes != NULL.
792  @param[in] refine If true, prepare the Mesh for conforming refinement of
793  triangular or tetrahedral meshes.
794  @param[in] fix_orientation
795  If true, fix the orientation of inverted mesh elements
796  by permuting their vertices.
797 
798  Before calling this method, call FinalizeTopology() and ensure that the
799  Mesh vertices or nodes are set. */
800  virtual void Finalize(bool refine = false, bool fix_orientation = false);
801 
802  virtual void SetAttributes();
803 
804  /** This is our integration with the Gecko library. The method finds an
805  element ordering that will increase memory coherency by putting elements
806  that are in physical proximity closer in memory. It can also be used to
807  obtain a space-filling curve ordering for ParNCMesh partitioning.
808  @param[out] ordering Output element ordering.
809  @param iterations Total number of V cycles. The ordering may improve with
810  more iterations. The best iteration is returned at the end.
811  @param window Initial window size. This determines the number of
812  permutations tested at each multigrid level and strongly influences the
813  quality of the result, but the cost of increasing 'window' is exponential.
814  @param period The window size is incremented every 'period' iterations.
815  @param seed Seed for initial random ordering (0 = skip random reorder).
816  @param verbose Print the progress of the optimization to mfem::out.
817  @param time_limit Optional time limit for the optimization, in seconds.
818  When reached, ordering from the best iteration so far is returned
819  (0 = no limit).
820  @return The final edge product cost of the ordering. The function may be
821  called in an external loop with different seeds, and the best ordering can
822  then be retained. */
823  double GetGeckoElementOrdering(Array<int> &ordering,
824  int iterations = 4, int window = 4,
825  int period = 2, int seed = 0,
826  bool verbose = false, double time_limit = 0);
827 
828  /** Return an ordering of the elements that approximately follows the Hilbert
829  curve. The method performs a spatial (Hilbert) sort on the centers of all
830  elements and returns the resulting sequence, which can then be passed to
831  ReorderElements. This is a cheap alternative to GetGeckoElementOrdering.*/
832  void GetHilbertElementOrdering(Array<int> &ordering);
833 
834  /** Rebuilds the mesh with a different order of elements. For each element i,
835  the array ordering[i] contains its desired new index. Note that the method
836  reorders vertices, edges and faces along with the elements. */
837  void ReorderElements(const Array<int> &ordering, bool reorder_vertices = true);
838 
839  /// Deprecated: see @a MakeCartesian3D.
840  MFEM_DEPRECATED
841  Mesh(int nx, int ny, int nz, Element::Type type, bool generate_edges = false,
842  double sx = 1.0, double sy = 1.0, double sz = 1.0,
843  bool sfc_ordering = true)
844  {
845  Make3D(nx, ny, nz, type, sx, sy, sz, sfc_ordering);
846  Finalize(true); // refine = true
847  }
848 
849  /// Deprecated: see @a MakeCartesian2D.
850  MFEM_DEPRECATED
851  Mesh(int nx, int ny, Element::Type type, bool generate_edges = false,
852  double sx = 1.0, double sy = 1.0, bool sfc_ordering = true)
853  {
854  Make2D(nx, ny, type, sx, sy, generate_edges, sfc_ordering);
855  Finalize(true); // refine = true
856  }
857 
858  /// Deprecated: see @a MakeCartesian1D.
859  MFEM_DEPRECATED
860  explicit Mesh(int n, double sx = 1.0)
861  {
862  Make1D(n, sx);
863  // Finalize(); // reminder: not needed
864  }
865 
866  /** Creates mesh by reading a file in MFEM, Netgen, or VTK format. If
867  generate_edges = 0 (default) edges are not generated, if 1 edges are
868  generated. See also @a Mesh::LoadFromFile. */
869  explicit Mesh(const char *filename, int generate_edges = 0, int refine = 1,
870  bool fix_orientation = true);
871 
872  /** Creates mesh by reading data stream in MFEM, Netgen, or VTK format. If
873  generate_edges = 0 (default) edges are not generated, if 1 edges are
874  generated. */
875  explicit Mesh(std::istream &input, int generate_edges = 0, int refine = 1,
876  bool fix_orientation = true);
877 
878  /// Create a disjoint mesh from the given mesh array
879  Mesh(Mesh *mesh_array[], int num_pieces);
880 
881  /// Deprecated: see @a MakeRefined.
882  MFEM_DEPRECATED
883  Mesh(Mesh *orig_mesh, int ref_factor, int ref_type);
884 
885  /** This is similar to the mesh constructor with the same arguments, but here
886  the current mesh is destroyed and another one created based on the data
887  stream again given in MFEM, Netgen, or VTK format. If generate_edges = 0
888  (default) edges are not generated, if 1 edges are generated. */
889  /// \see mfem::ifgzstream() for on-the-fly decompression of compressed ascii
890  /// inputs.
891  virtual void Load(std::istream &input, int generate_edges = 0,
892  int refine = 1, bool fix_orientation = true)
893  {
894  Loader(input, generate_edges);
895  Finalize(refine, fix_orientation);
896  }
897 
898  /// Clear the contents of the Mesh.
899  void Clear() { Destroy(); SetEmpty(); }
900 
901  /** @brief Get the mesh generator/type.
902 
903  The purpose of this is to be able to quickly tell what type of elements
904  one has in the mesh. Examination of this bitmask along with knowledge
905  of the mesh dimension can be used to identify which element types are
906  present.
907 
908  @return A bitmask:
909  - bit 0 - simplices are present in the mesh (triangles, tets),
910  - bit 1 - tensor product elements are present in the mesh (quads, hexes),
911  - bit 2 - the mesh has wedge elements.
912  - bit 3 - the mesh has pyramid elements.
913 
914  In parallel, the result takes into account elements on all processors.
915  */
916  inline int MeshGenerator() { return meshgen; }
917 
918  /** @brief Returns number of vertices. Vertices are only at the corners of
919  elements, where you would expect them in the lowest-order mesh. */
920  inline int GetNV() const { return NumOfVertices; }
921 
922  /// Returns number of elements.
923  inline int GetNE() const { return NumOfElements; }
924 
925  /// Returns number of boundary elements.
926  inline int GetNBE() const { return NumOfBdrElements; }
927 
928  /// Return the number of edges.
929  inline int GetNEdges() const { return NumOfEdges; }
930 
931  /// Return the number of faces in a 3D mesh.
932  inline int GetNFaces() const { return NumOfFaces; }
933 
934  /// Return the number of faces (3D), edges (2D) or vertices (1D).
935  int GetNumFaces() const;
936 
937  /** @brief Return the number of faces (3D), edges (2D) or vertices (1D)
938  including ghost faces. */
939  int GetNumFacesWithGhost() const;
940 
941  /** @brief Returns the number of faces according to the requested type, does
942  not count master nonconforming faces.
943 
944  If type==Boundary returns only the number of true boundary faces
945  contrary to GetNBE() that returns all "boundary" elements which may
946  include actual interior faces.
947  Similarly, if type==Interior, only the true interior faces are counted
948  excluding all master nonconforming faces. */
949  virtual int GetNFbyType(FaceType type) const;
950 
951  /// Utility function: sum integers from all processors (Allreduce).
952  virtual long long ReduceInt(int value) const { return value; }
953 
954  /// Return the total (global) number of elements.
955  long long GetGlobalNE() const { return ReduceInt(NumOfElements); }
956 
957  /** Return vertex to vertex table. The connections stored in the table
958  are from smaller to bigger vertex index, i.e. if i<j and (i, j) is
959  in the table, then (j, i) is not stored. */
960  void GetVertexToVertexTable(DSTable &) const;
961 
962  /** @brief Return the mesh geometric factors corresponding to the given
963  integration rule.
964 
965  The IntegrationRule used with GetGeometricFactors needs to remain valid
966  until the internally stored GeometricFactors objects are destroyed (by
967  either calling Mesh::DeleteGeometricFactors or the Mesh destructor). If
968  the device MemoryType parameter @a d_mt is specified, then the returned
969  object will use that type unless it was previously allocated with a
970  different type. */
972  const IntegrationRule& ir,
973  const int flags,
975 
976  /** @brief Return the mesh geometric factors for the faces corresponding
977  to the given integration rule.
978 
979  The IntegrationRule used with GetFaceGeometricFactors needs to remain
980  valid until the internally stored FaceGeometricFactors objects are
981  destroyed (by either calling Mesh::DeleteGeometricFactors or the Mesh
982  destructor). */
984  const int flags,
985  FaceType type,
987 
988  /// Destroy all GeometricFactors stored by the Mesh.
989  /** This method can be used to force recomputation of the GeometricFactors,
990  for example, after the mesh nodes are modified externally. */
991  void DeleteGeometricFactors();
992 
993  /// @brief This function should be called after the mesh node coordinates
994  /// have changed, e.g. after the mesh has moved.
995  /** It updates internal quantities derived from the node coordinates, such
996  as the GeometricFactors. */
998 
999  /// Equals 1 + num_holes - num_loops
1000  inline int EulerNumber() const
1002  /// Equals 1 - num_holes
1003  inline int EulerNumber2D() const
1004  { return NumOfVertices - NumOfEdges + NumOfElements; }
1005 
1006  int Dimension() const { return Dim; }
1007  int SpaceDimension() const { return spaceDim; }
1008 
1009  /// @brief Return pointer to vertex i's coordinates.
1010  /// @warning For high-order meshes (when #Nodes != NULL) vertices may not be
1011  /// updated and should not be used!
1012  const double *GetVertex(int i) const { return vertices[i](); }
1013 
1014  /// @brief Return pointer to vertex i's coordinates.
1015  /// @warning For high-order meshes (when Nodes != NULL) vertices may not
1016  /// being updated and should not be used!
1017  double *GetVertex(int i) { return vertices[i](); }
1018 
1019  void GetElementData(int geom, Array<int> &elem_vtx, Array<int> &attr) const
1020  { GetElementData(elements, geom, elem_vtx, attr); }
1021 
1022  /// Checks if the mesh has boundary elements
1023  virtual bool HasBoundaryElements() const { return (NumOfBdrElements > 0); }
1024 
1025  void GetBdrElementData(int geom, Array<int> &bdr_elem_vtx,
1026  Array<int> &bdr_attr) const
1027  { GetElementData(boundary, geom, bdr_elem_vtx, bdr_attr); }
1028 
1029  /** @brief Set the internal Vertex array to point to the given @a vertices
1030  array without assuming ownership of the pointer. */
1031  /** If @a zerocopy is `true`, the vertices must be given as an array of 3
1032  doubles per vertex. If @a zerocopy is `false` then the current Vertex
1033  data is first copied to the @a vertices array. */
1034  void ChangeVertexDataOwnership(double *vertices, int len_vertices,
1035  bool zerocopy = false);
1036 
1037  const Element* const *GetElementsArray() const
1038  { return elements.GetData(); }
1039 
1040  const Element *GetElement(int i) const { return elements[i]; }
1041 
1042  Element *GetElement(int i) { return elements[i]; }
1043 
1044  const Element *GetBdrElement(int i) const { return boundary[i]; }
1045 
1046  Element *GetBdrElement(int i) { return boundary[i]; }
1047 
1048  const Element *GetFace(int i) const { return faces[i]; }
1049 
1051  {
1052  return faces[i]->GetGeometryType();
1053  }
1054 
1056  {
1057  return elements[i]->GetGeometryType();
1058  }
1059 
1061  {
1062  return boundary[i]->GetGeometryType();
1063  }
1064 
1065  // deprecated: "base geometry" no longer means anything
1067  { return GetFaceGeometry(i); }
1068 
1070  { return GetElementGeometry(i); }
1071 
1073  { return GetBdrElementGeometry(i); }
1074 
1075  /** @brief Return true iff the given @a geom is encountered in the mesh.
1076  Geometries of dimensions lower than Dimension() are counted as well. */
1077  bool HasGeometry(Geometry::Type geom) const
1078  { return mesh_geoms & (1 << geom); }
1079 
1080  /** @brief Return the number of geometries of the given dimension present in
1081  the mesh. */
1082  /** For a parallel mesh only the local geometries are counted. */
1083  int GetNumGeometries(int dim) const;
1084 
1085  /// Return all element geometries of the given dimension present in the mesh.
1086  /** For a parallel mesh only the local geometries are returned.
1087 
1088  The returned geometries are sorted. */
1089  void GetGeometries(int dim, Array<Geometry::Type> &el_geoms) const;
1090 
1091  /// List of mesh geometries stored as Array<Geometry::Type>.
1092  class GeometryList : public Array<Geometry::Type>
1093  {
1094  protected:
1096  public:
1097  /// Construct a GeometryList of all element geometries in @a mesh.
1099  : Array<Geometry::Type>(geom_buf, Geometry::NumGeom)
1100  { mesh.GetGeometries(mesh.Dimension(), *this); }
1101  /** @brief Construct a GeometryList of all geometries of dimension @a dim
1102  in @a mesh. */
1103  GeometryList(Mesh &mesh, int dim)
1104  : Array<Geometry::Type>(geom_buf, Geometry::NumGeom)
1105  { mesh.GetGeometries(dim, *this); }
1106  };
1107 
1108  /// Returns the indices of the vertices of element i.
1109  void GetElementVertices(int i, Array<int> &v) const
1110  { elements[i]->GetVertices(v); }
1111 
1112  /// Returns the indices of the vertices of boundary element i.
1113  void GetBdrElementVertices(int i, Array<int> &v) const
1114  { boundary[i]->GetVertices(v); }
1115 
1116  /// Return the indices and the orientations of all edges of element i.
1117  void GetElementEdges(int i, Array<int> &edges, Array<int> &cor) const;
1118 
1119  /// Return the indices and the orientations of all edges of bdr element i.
1120  void GetBdrElementEdges(int i, Array<int> &edges, Array<int> &cor) const;
1121 
1122  /** Return the indices and the orientations of all edges of face i.
1123  Works for both 2D (face=edge) and 3D faces. */
1124  void GetFaceEdges(int i, Array<int> &edges, Array<int> &o) const;
1125 
1126  /// Returns the indices of the vertices of face i.
1127  void GetFaceVertices(int i, Array<int> &vert) const
1128  {
1129  if (Dim == 1)
1130  {
1131  vert.SetSize(1); vert[0] = i;
1132  }
1133  else
1134  {
1135  faces[i]->GetVertices(vert);
1136  }
1137  }
1138 
1139  /// Returns the indices of the vertices of edge i.
1140  void GetEdgeVertices(int i, Array<int> &vert) const;
1141 
1142  /// Returns the face-to-edge Table (3D)
1143  Table *GetFaceEdgeTable() const;
1144 
1145  /// Returns the edge-to-vertex Table (3D)
1146  Table *GetEdgeVertexTable() const;
1147 
1148  /// Return the indices and the orientations of all faces of element i.
1149  void GetElementFaces(int i, Array<int> &faces, Array<int> &ori) const;
1150 
1151  /// Return the index and the orientation of the face of bdr element i. (3D)
1152  void GetBdrElementFace(int i, int *f, int *o) const;
1153 
1154  /** Return the vertex index of boundary element i. (1D)
1155  Return the edge index of boundary element i. (2D)
1156  Return the face index of boundary element i. (3D) */
1157  int GetBdrElementEdgeIndex(int i) const;
1158 
1159  /** @brief For the given boundary element, bdr_el, return its adjacent
1160  element and its info, i.e. 64*local_bdr_index+bdr_orientation.
1161 
1162  The returned bdr_orientation is that of the boundary element relative to
1163  the respective face element.
1164 
1165  @sa GetBdrElementAdjacentElement2() */
1166  void GetBdrElementAdjacentElement(int bdr_el, int &el, int &info) const;
1167 
1168  /** @brief For the given boundary element, bdr_el, return its adjacent
1169  element and its info, i.e. 64*local_bdr_index+inverse_bdr_orientation.
1170 
1171  The returned inverse_bdr_orientation is the inverse of the orientation of
1172  the boundary element relative to the respective face element. In other
1173  words this is the orientation of the face element relative to the
1174  boundary element.
1175 
1176  @sa GetBdrElementAdjacentElement() */
1177  void GetBdrElementAdjacentElement2(int bdr_el, int &el, int &info) const;
1178 
1179  /// Returns the type of element i.
1180  Element::Type GetElementType(int i) const;
1181 
1182  /// Returns the type of boundary element i.
1183  Element::Type GetBdrElementType(int i) const;
1184 
1185  /* Return point matrix of element i of dimension Dim X #v, where for every
1186  vertex we give its coordinates in space of dimension Dim. */
1187  void GetPointMatrix(int i, DenseMatrix &pointmat) const;
1188 
1189  /* Return point matrix of boundary element i of dimension Dim X #v, where for
1190  every vertex we give its coordinates in space of dimension Dim. */
1191  void GetBdrPointMatrix(int i, DenseMatrix &pointmat) const;
1192 
1194 
1195  /** Builds the transformation defining the i-th element in the user-defined
1196  variable. */
1198 
1199  /// Returns the transformation defining the i-th element
1201 
1202  /** Return the transformation defining the i-th element assuming
1203  the position of the vertices/nodes are given by 'nodes'. */
1204  void GetElementTransformation(int i, const Vector &nodes,
1206 
1207  /// Returns the transformation defining the i-th boundary element
1210 
1211  /** @brief Returns the transformation defining the given face element in a
1212  user-defined variable. */
1214 
1215  /** @brief A helper method that constructs a transformation from the
1216  reference space of a face to the reference space of an element. */
1217  /** The local index of the face as a face in the element and its orientation
1218  are given by the input parameter @a info, as @a info = 64*loc_face_idx +
1219  loc_face_orientation. */
1220  void GetLocalFaceTransformation(int face_type, int elem_type,
1222  int info);
1223 
1224  /// Returns the transformation defining the given face element
1226 
1227  /** Returns the transformation defining the given edge element.
1228  The transformation is stored in a user-defined variable. */
1230 
1231  /// Returns the transformation defining the given face element
1233 
1234  /// Returns (a pointer to an object containing) the following data:
1235  ///
1236  /// 1) Elem1No - the index of the first element that contains this face this
1237  /// is the element that has the same outward unit normal vector as the
1238  /// face;
1239  ///
1240  /// 2) Elem2No - the index of the second element that contains this face this
1241  /// element has outward unit normal vector as the face multiplied with -1;
1242  ///
1243  /// 3) Elem1, Elem2 - pointers to the ElementTransformation's of the first
1244  /// and the second element respectively;
1245  ///
1246  /// 4) Face - pointer to the ElementTransformation of the face;
1247  ///
1248  /// 5) Loc1, Loc2 - IntegrationPointTransformation's mapping the face
1249  /// coordinate system to the element coordinate system (both in their
1250  /// reference elements). Used to transform IntegrationPoints from face to
1251  /// element. More formally, let:
1252  /// TL1, TL2 be the transformations represented by Loc1, Loc2,
1253  /// TE1, TE2 - the transformations represented by Elem1, Elem2,
1254  /// TF - the transformation represented by Face, then
1255  /// TF(x) = TE1(TL1(x)) = TE2(TL2(x)) for all x in the reference face.
1256  ///
1257  /// 6) FaceGeom - the base geometry for the face.
1258  ///
1259  /// The mask specifies which fields in the structure to return:
1260  /// mask & 1 - Elem1, mask & 2 - Elem2
1261  /// mask & 4 - Loc1, mask & 8 - Loc2, mask & 16 - Face.
1262  /// These mask values are defined in the ConfigMasks enum type as part of the
1263  /// FaceElementTransformations class in fem/eltrans.hpp.
1265  int FaceNo,
1266  int mask = 31);
1267 
1269  {
1270  if (faces_info[FaceNo].Elem2No < 0) { return NULL; }
1271  return GetFaceElementTransformations (FaceNo);
1272  }
1273 
1275 
1276  /// Return the local face index for the given boundary face.
1277  int GetBdrFace(int BdrElemNo) const;
1278 
1279  /// Return true if the given face is interior. @sa FaceIsTrueInterior().
1280  bool FaceIsInterior(int FaceNo) const
1281  {
1282  return (faces_info[FaceNo].Elem2No >= 0);
1283  }
1284 
1285  /** This enumerated type describes the three main face topologies:
1286  - Boundary, for faces on the boundary of the computational domain,
1287  - Conforming, for conforming faces interior to the computational domain,
1288  - Nonconforming, for nonconforming faces interior to the computational
1289  domain. */
1290  enum class FaceTopology { Boundary,
1291  Conforming,
1292  Nonconforming,
1293  NA
1294  };
1295 
1296  /** This enumerated type describes the location of the two elements sharing a
1297  face, Local meaning that the element is local to the MPI rank, FaceNbr
1298  meaning that the element is distributed on a different MPI rank, this
1299  typically means that methods with FaceNbr should be used to access the
1300  relevant information, e.g., ParFiniteElementSpace::GetFaceNbrElementVDofs.
1301  */
1302  enum class ElementLocation { Local, FaceNbr, NA };
1303 
1304  /** This enumerated type describes the topological relation of an element to
1305  a face:
1306  - Coincident meaning that the element's face is topologically equal to
1307  the mesh face.
1308  - Superset meaning that the element's face is topologically coarser than
1309  the mesh face, i.e., the element's face contains the mesh face.
1310  - Subset meaning that the element's face is topologically finer than the
1311  mesh face, i.e., the element's face is contained in the mesh face.
1312  Superset and Subset are only relevant for nonconforming faces.
1313  Master nonconforming faces have a conforming element on one side, and a
1314  fine element on the other side. Slave nonconforming faces have a
1315  conforming element on one side, and a coarse element on the other side.
1316  */
1318 
1319  /** This enumerated type describes the corresponding FaceInfo internal
1320  representation (encoded cases), c.f. FaceInfo's documentation:
1321  Classification of a local (non-ghost) face based on its FaceInfo:
1322  - Elem2No >= 0 --> local interior face; can be either:
1323  - NCFace == -1 --> LocalConforming,
1324  - NCFace >= 0 --> LocalSlaveNonconforming,
1325  - Elem2No < 0 --> local "boundary" face; can be one of:
1326  - NCFace == -1 --> conforming face; can be either:
1327  - Elem2Inf < 0 --> Boundary,
1328  - Elem2Inf >= 0 --> SharedConforming,
1329  - NCFace >= 0 --> nonconforming face; can be one of:
1330  - Elem2Inf < 0 --> MasterNonconforming (shared or not shared),
1331  - Elem2Inf >= 0 --> SharedSlaveNonconforming.
1332  Classification of a ghost (non-local) face based on its FaceInfo:
1333  - Elem1No == -1 --> GhostMaster (includes other unused ghost faces),
1334  - Elem1No >= 0 --> GhostSlave.
1335  */
1336  enum class FaceInfoTag { Boundary,
1342  GhostSlave,
1343  GhostMaster
1344  };
1345 
1346  /** @brief This structure is used as a human readable output format that
1347  decipheres the information contained in Mesh::FaceInfo when using the
1348  Mesh::GetFaceInformation() method.
1349 
1350  The element indices in this structure don't need further processing,
1351  contrary to the ones obtained through Mesh::GetFacesElements and can
1352  directly be used, e.g., Elem1 and Elem2 indices.
1353  Likewise the orientations for Elem1 and Elem2 already take into account
1354  special cases and can be used as is.
1355  */
1357  {
1359 
1360  struct
1361  {
1364  int index;
1367  } element[2];
1368 
1370  int ncface;
1372 
1373  /** @brief Return true if the face is a local interior face which is NOT
1374  a master nonconforming face. */
1375  bool IsLocal() const
1376  {
1377  return element[1].location == Mesh::ElementLocation::Local;
1378  }
1379 
1380  /** @brief Return true if the face is a shared interior face which is NOT
1381  a master nonconforming face. */
1382  bool IsShared() const
1383  {
1384  return element[1].location == Mesh::ElementLocation::FaceNbr;
1385  }
1386 
1387  /** @brief return true if the face is an interior face to the computation
1388  domain, either a local or shared interior face (not a boundary face)
1389  which is NOT a master nonconforming face.
1390  */
1391  bool IsInterior() const
1392  {
1393  return topology == FaceTopology::Conforming ||
1395  }
1396 
1397  /** @brief Return true if the face is a boundary face. */
1398  bool IsBoundary() const
1399  {
1400  return topology == FaceTopology::Boundary;
1401  }
1402 
1403  /// @brief Return true if the face is of the same type as @a type.
1404  bool IsOfFaceType(FaceType type) const
1405  {
1406  switch (type)
1407  {
1408  case FaceType::Interior:
1409  return IsInterior();
1410  case FaceType::Boundary:
1411  return IsBoundary();
1412  default:
1413  return false;
1414  }
1415  }
1416 
1417  /// @brief Return true if the face is a conforming face.
1418  bool IsConforming() const
1419  {
1421  }
1422 
1423  /// @brief Return true if the face is a nonconforming fine face.
1424  bool IsNonconformingFine() const
1425  {
1427  (element[0].conformity == ElementConformity::Superset ||
1428  element[1].conformity == ElementConformity::Superset);
1429  }
1430 
1431  /// @brief Return true if the face is a nonconforming coarse face.
1432  /** Note that ghost nonconforming master faces cannot be clearly
1433  identified as such with the currently available information, so this
1434  method will return false for such faces. */
1436  {
1438  element[1].conformity == ElementConformity::Subset;
1439  }
1440 
1441  /// @brief cast operator from FaceInformation to FaceInfo.
1442  operator Mesh::FaceInfo() const;
1443  };
1444 
1445  /** This method aims to provide face information in a deciphered format, i.e.
1446  Mesh::FaceInformation, compared to the raw encoded information returned
1447  by Mesh::GetFaceElements() and Mesh::GetFaceInfos(). */
1448  FaceInformation GetFaceInformation(int f) const;
1449 
1450  void GetFaceElements (int Face, int *Elem1, int *Elem2) const;
1451  void GetFaceInfos (int Face, int *Inf1, int *Inf2) const;
1452  void GetFaceInfos (int Face, int *Inf1, int *Inf2, int *NCFace) const;
1453 
1454  Geometry::Type GetFaceGeometryType(int Face) const;
1455  Element::Type GetFaceElementType(int Face) const;
1456 
1457  Array<int> GetFaceToBdrElMap() const;
1458 
1459  /// Check (and optionally attempt to fix) the orientation of the elements
1460  /** @param[in] fix_it If `true`, attempt to fix the orientations of some
1461  elements: triangles, quads, and tets.
1462  @return The number of elements with wrong orientation.
1463 
1464  @note For meshes with nodes (e.g. high-order or periodic meshes), fixing
1465  the element orientations may require additional permutation of the nodal
1466  GridFunction of the mesh which is not performed by this method. Instead,
1467  the method Finalize() should be used with the parameter
1468  @a fix_orientation set to `true`.
1469 
1470  @note This method performs a simple check if an element is inverted, e.g.
1471  for most elements types, it checks if the Jacobian of the mapping from
1472  the reference element is non-negative at the center of the element. */
1473  int CheckElementOrientation(bool fix_it = true);
1474 
1475  /// Check the orientation of the boundary elements
1476  /** @return The number of boundary elements with wrong orientation. */
1477  int CheckBdrElementOrientation(bool fix_it = true);
1478 
1479  /// Return the attribute of element i.
1480  int GetAttribute(int i) const { return elements[i]->GetAttribute(); }
1481 
1482  /// Set the attribute of element i.
1483  void SetAttribute(int i, int attr) { elements[i]->SetAttribute(attr); }
1484 
1485  /// Return the attribute of boundary element i.
1486  int GetBdrAttribute(int i) const { return boundary[i]->GetAttribute(); }
1487 
1488  /// Set the attribute of boundary element i.
1489  void SetBdrAttribute(int i, int attr) { boundary[i]->SetAttribute(attr); }
1490 
1491  const Table &ElementToElementTable();
1492 
1493  const Table &ElementToFaceTable() const;
1494 
1495  const Table &ElementToEdgeTable() const;
1496 
1497  /// The returned Table must be destroyed by the caller
1499 
1500  /** Return the "face"-element Table. Here "face" refers to face (3D),
1501  edge (2D), or vertex (1D).
1502  The returned Table must be destroyed by the caller. */
1503  Table *GetFaceToElementTable() const;
1504 
1505  /** This method modifies a tetrahedral mesh so that Nedelec spaces of order
1506  greater than 1 can be defined on the mesh. Specifically, we
1507  1) rotate all tets in the mesh so that the vertices {v0, v1, v2, v3}
1508  satisfy: v0 < v1 < min(v2, v3).
1509  2) rotate all boundary triangles so that the vertices {v0, v1, v2}
1510  satisfy: v0 < min(v1, v2).
1511 
1512  @note Refinement does not work after a call to this method! */
1513  MFEM_DEPRECATED virtual void ReorientTetMesh();
1514 
1515  int *CartesianPartitioning(int nxyz[]);
1516  int *GeneratePartitioning(int nparts, int part_method = 1);
1517  void CheckPartitioning(int *partitioning_);
1518 
1519  void CheckDisplacements(const Vector &displacements, double &tmax);
1520 
1521  // Vertices are only at the corners of elements, where you would expect them
1522  // in the lowest-order mesh.
1523  void MoveVertices(const Vector &displacements);
1524  void GetVertices(Vector &vert_coord) const;
1525  void SetVertices(const Vector &vert_coord);
1526 
1527  // Nodes are only active for higher order meshes, and share locations with
1528  // the vertices, plus all the higher- order control points within the element
1529  // and along the edges and on the faces.
1530  void GetNode(int i, double *coord) const;
1531  void SetNode(int i, const double *coord);
1532 
1533  // Node operations for curved mesh.
1534  // They call the corresponding '...Vertices' method if the
1535  // mesh is not curved (i.e. Nodes == NULL).
1536  void MoveNodes(const Vector &displacements);
1537  void GetNodes(Vector &node_coord) const;
1538  void SetNodes(const Vector &node_coord);
1539 
1540  /// Return a pointer to the internal node GridFunction (may be NULL).
1541  GridFunction *GetNodes() { return Nodes; }
1542  const GridFunction *GetNodes() const { return Nodes; }
1543  /// Return the mesh nodes ownership flag.
1544  bool OwnsNodes() const { return own_nodes; }
1545  /// Set the mesh nodes ownership flag.
1546  void SetNodesOwner(bool nodes_owner) { own_nodes = nodes_owner; }
1547  /// Replace the internal node GridFunction with the given GridFunction.
1548  void NewNodes(GridFunction &nodes, bool make_owner = false);
1549  /** Swap the internal node GridFunction pointer and ownership flag members
1550  with the given ones. */
1551  void SwapNodes(GridFunction *&nodes, int &own_nodes_);
1552 
1553  /// Return the mesh nodes/vertices projected on the given GridFunction.
1554  void GetNodes(GridFunction &nodes) const;
1555  /** Replace the internal node GridFunction with a new GridFunction defined
1556  on the given FiniteElementSpace. The new node coordinates are projected
1557  (derived) from the current nodes/vertices. */
1558  virtual void SetNodalFESpace(FiniteElementSpace *nfes);
1559  /** Replace the internal node GridFunction with the given GridFunction. The
1560  given GridFunction is updated with node coordinates projected (derived)
1561  from the current nodes/vertices. */
1562  void SetNodalGridFunction(GridFunction *nodes, bool make_owner = false);
1563  /** Return the FiniteElementSpace on which the current mesh nodes are
1564  defined or NULL if the mesh does not have nodes. */
1565  const FiniteElementSpace *GetNodalFESpace() const;
1566  /** Make sure that the mesh has valid nodes, i.e. its geometry is described
1567  by a vector finite element grid function (even if it is a low-order mesh
1568  with straight edges). */
1569  void EnsureNodes();
1570 
1571  /** Set the curvature of the mesh nodes using the given polynomial degree,
1572  'order', and optionally: discontinuous or continuous FE space, 'discont',
1573  new space dimension, 'space_dim' (if != -1), and 'ordering' (byVDim by
1574  default). */
1575  virtual void SetCurvature(int order, bool discont = false, int space_dim = -1,
1576  int ordering = 1);
1577 
1578  /// Refine all mesh elements.
1579  /** @param[in] ref_algo %Refinement algorithm. Currently used only for pure
1580  tetrahedral meshes. If set to zero (default), a tet mesh will be refined
1581  using algorithm A, that produces elements with better quality compared to
1582  algorithm B used when the parameter is non-zero.
1583 
1584  For tetrahedral meshes, after using algorithm A, the mesh cannot be
1585  refined locally using methods like GeneralRefinement() unless it is
1586  re-finalized using Finalize() with the parameter @a refine set to true.
1587  Note that calling Finalize() in this way will generally invalidate any
1588  FiniteElementSpace%s and GridFunction%s defined on the mesh. */
1589  void UniformRefinement(int ref_algo = 0);
1590 
1591  /** Refine selected mesh elements. Refinement type can be specified for each
1592  element. The function can do conforming refinement of triangles and
1593  tetrahedra and nonconforming refinement (i.e., with hanging-nodes) of
1594  triangles, quadrilaterals and hexahedra. If 'nonconforming' = -1,
1595  suitable refinement method is selected automatically (namely, conforming
1596  refinement for triangles). Use nonconforming = 0/1 to force the method.
1597  For nonconforming refinements, nc_limit optionally specifies the maximum
1598  level of hanging nodes (unlimited by default). */
1599  void GeneralRefinement(const Array<Refinement> &refinements,
1600  int nonconforming = -1, int nc_limit = 0);
1601 
1602  /** Simplified version of GeneralRefinement taking a simple list of elements
1603  to refine, without refinement types. */
1604  void GeneralRefinement(const Array<int> &el_to_refine,
1605  int nonconforming = -1, int nc_limit = 0);
1606 
1607  /// Refine each element with given probability. Uses GeneralRefinement.
1608  void RandomRefinement(double prob, bool aniso = false,
1609  int nonconforming = -1, int nc_limit = 0);
1610 
1611  /// Refine elements sharing the specified vertex. Uses GeneralRefinement.
1612  void RefineAtVertex(const Vertex& vert,
1613  double eps = 0.0, int nonconforming = -1);
1614 
1615  /** Refine element i if elem_error[i] > threshold, for all i.
1616  Returns true if at least one element was refined, false otherwise. */
1617  bool RefineByError(const Array<double> &elem_error, double threshold,
1618  int nonconforming = -1, int nc_limit = 0);
1619 
1620  /** Refine element i if elem_error(i) > threshold, for all i.
1621  Returns true if at least one element was refined, false otherwise. */
1622  bool RefineByError(const Vector &elem_error, double threshold,
1623  int nonconforming = -1, int nc_limit = 0);
1624 
1625  /** Derefine the mesh based on an error measure associated with each
1626  element. A derefinement is performed if the sum of errors of its fine
1627  elements is smaller than 'threshold'. If 'nc_limit' > 0, derefinements
1628  that would increase the maximum level of hanging nodes of the mesh are
1629  skipped. Returns true if the mesh changed, false otherwise. */
1630  bool DerefineByError(Array<double> &elem_error, double threshold,
1631  int nc_limit = 0, int op = 1);
1632 
1633  /// Same as DerefineByError for an error vector.
1634  bool DerefineByError(const Vector &elem_error, double threshold,
1635  int nc_limit = 0, int op = 1);
1636 
1637  ///@{ @name NURBS mesh refinement methods
1638  void KnotInsert(Array<KnotVector *> &kv);
1639  void KnotInsert(Array<Vector *> &kv);
1640  /* For each knot vector:
1641  new_degree = max(old_degree, min(old_degree + rel_degree, degree)). */
1642  void DegreeElevate(int rel_degree, int degree = 16);
1643  ///@}
1644 
1645  /** Make sure that a quad/hex mesh is considered to be nonconforming (i.e.,
1646  has an associated NCMesh object). Simplex meshes can be both conforming
1647  (default) or nonconforming. */
1648  void EnsureNCMesh(bool simplices_nonconforming = false);
1649 
1650  bool Conforming() const { return ncmesh == NULL; }
1651  bool Nonconforming() const { return ncmesh != NULL; }
1652 
1653  /** Return fine element transformations following a mesh refinement.
1654  Space uses this to construct a global interpolation matrix. */
1656 
1657  /// Return type of last modification of the mesh.
1659 
1660  /** Return update counter. The counter starts at zero and is incremented
1661  each time refinement, derefinement, or rebalancing method is called.
1662  It is used for checking proper sequence of Space:: and GridFunction::
1663  Update() calls. */
1664  long GetSequence() const { return sequence; }
1665 
1666  /// Print the mesh to the given stream using Netgen/Truegrid format.
1667  virtual void PrintXG(std::ostream &os = mfem::out) const;
1668 
1669  /// Print the mesh to the given stream using the default MFEM mesh format.
1670  /// \see mfem::ofgzstream() for on-the-fly compression of ascii outputs
1671  virtual void Print(std::ostream &os = mfem::out) const { Printer(os); }
1672 
1673  /// Save the mesh to a file using Mesh::Print. The given @a precision will be
1674  /// used for ASCII output.
1675  virtual void Save(const char *fname, int precision=16) const;
1676 
1677  /// Print the mesh to the given stream using the adios2 bp format
1678 #ifdef MFEM_USE_ADIOS2
1679  virtual void Print(adios2stream &os) const;
1680 #endif
1681  /// Print the mesh in VTK format (linear and quadratic meshes only).
1682  /// \see mfem::ofgzstream() for on-the-fly compression of ascii outputs
1683  void PrintVTK(std::ostream &os);
1684  /** Print the mesh in VTK format. The parameter ref > 0 specifies an element
1685  subdivision number (useful for high order fields and curved meshes).
1686  If the optional field_data is set, we also add a FIELD section in the
1687  beginning of the file with additional dataset information. */
1688  /// \see mfem::ofgzstream() for on-the-fly compression of ascii outputs
1689  void PrintVTK(std::ostream &os, int ref, int field_data=0);
1690  /** Print the mesh in VTU format. The parameter ref > 0 specifies an element
1691  subdivision number (useful for high order fields and curved meshes).
1692  If @a bdr_elements is true, then output (only) the boundary elements,
1693  otherwise output only the non-boundary elements. */
1694  void PrintVTU(std::ostream &os,
1695  int ref=1,
1696  VTKFormat format=VTKFormat::ASCII,
1697  bool high_order_output=false,
1698  int compression_level=0,
1699  bool bdr_elements=false);
1700  /** Print the mesh in VTU format with file name fname. */
1701  virtual void PrintVTU(std::string fname,
1702  VTKFormat format=VTKFormat::ASCII,
1703  bool high_order_output=false,
1704  int compression_level=0,
1705  bool bdr=false);
1706  /** Print the boundary elements of the mesh in VTU format, and output the
1707  boundary attributes as a data array (useful for boundary conditions). */
1708  void PrintBdrVTU(std::string fname,
1709  VTKFormat format=VTKFormat::ASCII,
1710  bool high_order_output=false,
1711  int compression_level=0);
1712 
1713  void GetElementColoring(Array<int> &colors, int el0 = 0);
1714 
1715  /** @brief Prints the mesh with boundary elements given by the boundary of
1716  the subdomains, so that the boundary of subdomain i has boundary
1717  attribute i+1. */
1718  /// \see mfem::ofgzstream() for on-the-fly compression of ascii outputs
1719  void PrintWithPartitioning (int *partitioning,
1720  std::ostream &os, int elem_attr = 0) const;
1721 
1722  void PrintElementsWithPartitioning (int *partitioning,
1723  std::ostream &out,
1724  int interior_faces = 0);
1725 
1726  /// Print set of disjoint surfaces:
1727  /*!
1728  * If Aface_face(i,j) != 0, print face j as a boundary
1729  * element with attribute i+1.
1730  */
1731  void PrintSurfaces(const Table &Aface_face, std::ostream &out) const;
1732 
1733  void ScaleSubdomains (double sf);
1734  void ScaleElements (double sf);
1735 
1736  void Transform(void (*f)(const Vector&, Vector&));
1737  void Transform(VectorCoefficient &deformation);
1738 
1739  /// Remove unused vertices and rebuild mesh connectivity.
1740  void RemoveUnusedVertices();
1741 
1742  /** Remove boundary elements that lie in the interior of the mesh, i.e. that
1743  have two adjacent faces in 3D, or edges in 2D. */
1744  void RemoveInternalBoundaries();
1745 
1746  /** @brief Get the size of the i-th element relative to the perfect
1747  reference element. */
1748  double GetElementSize(int i, int type = 0);
1749 
1750  double GetElementSize(int i, const Vector &dir);
1751 
1752  double GetElementVolume(int i);
1753 
1754  void GetElementCenter(int i, Vector &center);
1755 
1756  /// Returns the minimum and maximum corners of the mesh bounding box.
1757  /** For high-order meshes, the geometry is first refined @a ref times. */
1758  void GetBoundingBox(Vector &min, Vector &max, int ref = 2);
1759 
1760  void GetCharacteristics(double &h_min, double &h_max,
1761  double &kappa_min, double &kappa_max,
1762  Vector *Vh = NULL, Vector *Vk = NULL);
1763 
1764  /// Auxiliary method used by PrintCharacteristics().
1765  /** It is also used in the `mesh-explorer` miniapp. */
1766  static void PrintElementsByGeometry(int dim,
1767  const Array<int> &num_elems_by_geom,
1768  std::ostream &out);
1769 
1770  /** @brief Compute and print mesh characteristics such as number of vertices,
1771  number of elements, number of boundary elements, minimal and maximal
1772  element sizes, minimal and maximal element aspect ratios, etc. */
1773  /** If @a Vh or @a Vk are not NULL, return the element sizes and aspect
1774  ratios for all elements in the given Vector%s. */
1775  void PrintCharacteristics(Vector *Vh = NULL, Vector *Vk = NULL,
1776  std::ostream &os = mfem::out);
1777 
1778  /** @brief In serial, this method calls PrintCharacteristics(). In parallel,
1779  additional information about the parallel decomposition is also printed.
1780  */
1781  virtual void PrintInfo(std::ostream &os = mfem::out)
1782  {
1783  PrintCharacteristics(NULL, NULL, os);
1784  }
1785 
1786  /** @brief Find the ids of the elements that contain the given points, and
1787  their corresponding reference coordinates.
1788 
1789  The DenseMatrix @a point_mat describes the given points - one point for
1790  each column; it should have SpaceDimension() rows.
1791 
1792  The InverseElementTransformation object, @a inv_trans, is used to attempt
1793  the element transformation inversion. If NULL pointer is given, the
1794  method will use a default constructed InverseElementTransformation. Note
1795  that the algorithms in the base class InverseElementTransformation can be
1796  completely overwritten by deriving custom classes that override the
1797  Transform() method.
1798 
1799  If no element is found for the i-th point, elem_ids[i] is set to -1.
1800 
1801  In the ParMesh implementation, the @a point_mat is expected to be the
1802  same on all ranks. If the i-th point is found by multiple ranks, only one
1803  of them will mark that point as found, i.e. set its elem_ids[i] to a
1804  non-negative number; the other ranks will set their elem_ids[i] to -2 to
1805  indicate that the point was found but assigned to another rank.
1806 
1807  @returns The total number of points that were found.
1808 
1809  @note This method is not 100 percent reliable, i.e. it is not guaranteed
1810  to find a point, even if it lies inside a mesh element. */
1811  virtual int FindPoints(DenseMatrix& point_mat, Array<int>& elem_ids,
1812  Array<IntegrationPoint>& ips, bool warn = true,
1813  InverseElementTransformation *inv_trans = NULL);
1814 
1815  /// Swaps internal data with another mesh. By default, non-geometry members
1816  /// like 'ncmesh' and 'NURBSExt' are only swapped when 'non_geometry' is set.
1817  void Swap(Mesh& other, bool non_geometry);
1818 
1819  /// Destroys Mesh.
1820  virtual ~Mesh() { DestroyPointers(); }
1821 
1822 #ifdef MFEM_DEBUG
1823  /// Output an NCMesh-compatible debug dump.
1824  void DebugDump(std::ostream &out) const;
1825 #endif
1826 };
1827 
1828 /** Overload operator<< for std::ostream and Mesh; valid also for the derived
1829  class ParMesh */
1830 std::ostream &operator<<(std::ostream &out, const Mesh &mesh);
1831 
1832 
1833 /** @brief Structure for storing mesh geometric factors: coordinates, Jacobians,
1834  and determinants of the Jacobians. */
1835 /** Typically objects of this type are constructed and owned by objects of class
1836  Mesh. See Mesh::GetGeometricFactors(). */
1838 {
1839 
1840 private:
1841  void Compute(const GridFunction &nodes,
1843 
1844 public:
1845  const Mesh *mesh;
1848 
1850  {
1851  COORDINATES = 1 << 0,
1852  JACOBIANS = 1 << 1,
1853  DETERMINANTS = 1 << 2,
1854  };
1855 
1856  GeometricFactors(const Mesh *mesh, const IntegrationRule &ir, int flags,
1858 
1859  GeometricFactors(const GridFunction &nodes, const IntegrationRule &ir,
1860  int flags,
1862 
1863  /// Mapped (physical) coordinates of all quadrature points.
1864  /** This array uses a column-major layout with dimensions (NQ x SDIM x NE)
1865  where
1866  - NQ = number of quadrature points per element,
1867  - SDIM = space dimension of the mesh = mesh.SpaceDimension(), and
1868  - NE = number of elements in the mesh. */
1870 
1871  /// Jacobians of the element transformations at all quadrature points.
1872  /** This array uses a column-major layout with dimensions (NQ x SDIM x DIM x
1873  NE) where
1874  - NQ = number of quadrature points per element,
1875  - SDIM = space dimension of the mesh = mesh.SpaceDimension(),
1876  - DIM = dimension of the mesh = mesh.Dimension(), and
1877  - NE = number of elements in the mesh. */
1879 
1880  /// Determinants of the Jacobians at all quadrature points.
1881  /** This array uses a column-major layout with dimensions (NQ x NE) where
1882  - NQ = number of quadrature points per element, and
1883  - NE = number of elements in the mesh. */
1885 };
1886 
1887 /** @brief Structure for storing face geometric factors: coordinates, Jacobians,
1888  determinants of the Jacobians, and normal vectors. */
1889 /** Typically objects of this type are constructed and owned by objects of class
1890  Mesh. See Mesh::GetFaceGeometricFactors(). */
1892 {
1893 public:
1894  const Mesh *mesh;
1898 
1900  {
1901  COORDINATES = 1 << 0,
1902  JACOBIANS = 1 << 1,
1903  DETERMINANTS = 1 << 2,
1904  NORMALS = 1 << 3,
1905  };
1906 
1907  FaceGeometricFactors(const Mesh *mesh, const IntegrationRule &ir, int flags,
1909 
1910  /// Mapped (physical) coordinates of all quadrature points.
1911  /** This array uses a column-major layout with dimensions (NQ x SDIM x NF)
1912  where
1913  - NQ = number of quadrature points per face,
1914  - SDIM = space dimension of the mesh = mesh.SpaceDimension(), and
1915  - NF = number of faces in the mesh. */
1917 
1918  /// Jacobians of the element transformations at all quadrature points.
1919  /** This array uses a column-major layout with dimensions (NQ x SDIM x DIM x
1920  NF) where
1921  - NQ = number of quadrature points per face,
1922  - SDIM = space dimension of the mesh = mesh.SpaceDimension(),
1923  - DIM = dimension of the mesh = mesh.Dimension(), and
1924  - NF = number of faces in the mesh. */
1926 
1927  /// Determinants of the Jacobians at all quadrature points.
1928  /** This array uses a column-major layout with dimensions (NQ x NF) where
1929  - NQ = number of quadrature points per face, and
1930  - NF = number of faces in the mesh. */
1932 
1933  /// Normals at all quadrature points.
1934  /** This array uses a column-major layout with dimensions (NQ x DIM x NF) where
1935  - NQ = number of quadrature points per face,
1936  - SDIM = space dimension of the mesh = mesh.SpaceDimension(), and
1937  - NF = number of faces in the mesh. */
1939 };
1940 
1941 /// Class used to extrude the nodes of a mesh
1943 {
1944 private:
1945  int n, layer;
1946  double p[2], s;
1947  Vector tip;
1948 public:
1949  NodeExtrudeCoefficient(const int dim, const int n_, const double s_);
1950  void SetLayer(const int l) { layer = l; }
1952  virtual void Eval(Vector &V, ElementTransformation &T,
1953  const IntegrationPoint &ip);
1955 };
1956 
1957 
1958 /// Extrude a 1D mesh
1959 Mesh *Extrude1D(Mesh *mesh, const int ny, const double sy,
1960  const bool closed = false);
1961 
1962 /// Extrude a 2D mesh
1963 Mesh *Extrude2D(Mesh *mesh, const int nz, const double sz);
1964 
1965 // shift cyclically 3 integers left-to-right
1966 inline void ShiftRight(int &a, int &b, int &c)
1967 {
1968  int t = a;
1969  a = c; c = b; b = t;
1970 }
1971 
1972 /// @brief Print function for Mesh::FaceInformation.
1973 std::ostream& operator<<(std::ostream& os, const Mesh::FaceInformation& info);
1974 
1975 }
1976 
1977 #endif
static Mesh MakePeriodic(const Mesh &orig_mesh, const std::vector< int > &v2v)
Create a periodic mesh by identifying vertices of orig_mesh.
Definition: mesh.cpp:4891
Abstract class for all finite elements.
Definition: fe_base.hpp:235
const IntegrationRule * IntRule
Definition: mesh.hpp:1895
MFEM_DEPRECATED Mesh(int n, double sx=1.0)
Deprecated: see MakeCartesian1D.
Definition: mesh.hpp:860
void Loader(std::istream &input, int generate_edges=0, std::string parse_tag="")
Definition: mesh.cpp:3966
double GetGeckoElementOrdering(Array< int > &ordering, int iterations=4, int window=4, int period=2, int seed=0, bool verbose=false, double time_limit=0)
Definition: mesh.cpp:2030
Geometry::Constants< Geometry::PYRAMID > pyr_t
Definition: mesh.hpp:263
void PrintSurfaces(const Table &Aface_face, std::ostream &out) const
Print set of disjoint surfaces:
Definition: mesh.cpp:11422
virtual void Save(const char *fname, int precision=16) const
Definition: mesh.cpp:10289
void GetPointMatrix(int i, DenseMatrix &pointmat) const
Definition: mesh.cpp:6260
bool IsNonconformingFine() const
Return true if the face is a nonconforming fine face.
Definition: mesh.hpp:1424
static const int vtk_quadratic_hex[27]
Definition: mesh.hpp:246
int * CartesianPartitioning(int nxyz[])
Definition: mesh.cpp:6989
FaceInformation GetFaceInformation(int f) const
Definition: mesh.cpp:1131
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:90
const DenseMatrix * PointMatrix
Definition: mesh.hpp:208
int GetBdrAttribute(int i) const
Return the attribute of boundary element i.
Definition: mesh.hpp:1486
const double * GetVertex(int i) const
Return pointer to vertex i&#39;s coordinates.
Definition: mesh.hpp:1012
void ScaleElements(double sf)
Definition: mesh.cpp:11559
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:30
Table * GetEdgeVertexTable() const
Returns the edge-to-vertex Table (3D)
Definition: mesh.cpp:6053
void GetBdrElementEdges(int i, Array< int > &edges, Array< int > &cor) const
Return the indices and the orientations of all edges of bdr element i.
Definition: mesh.cpp:5954
void ReadVTKMesh(std::istream &input, int &curved, int &read_gf, bool &finalize_topo)
void FreeElement(Element *E)
Definition: mesh.cpp:11873
int AddQuad(int v1, int v2, int v3, int v4, int attr=1)
Definition: mesh.cpp:1674
int CheckElementOrientation(bool fix_it=true)
Check (and optionally attempt to fix) the orientation of the elements.
Definition: mesh.cpp:5420
static bool remove_unused_vertices
Definition: mesh.hpp:281
void SetVertices(const Vector &vert_coord)
Definition: mesh.cpp:7848
static const int vtk_quadratic_tet[10]
Definition: mesh.hpp:243
void GetEdgeVertices(int i, Array< int > &vert) const
Returns the indices of the vertices of edge i.
Definition: mesh.cpp:6016
void Make2D(int nx, int ny, Element::Type type, double sx, double sy, bool generate_edges, bool sfc_ordering)
Definition: mesh.cpp:3339
bool FaceIsTrueInterior(int FaceNo) const
Definition: mesh.hpp:492
Vector J
Jacobians of the element transformations at all quadrature points.
Definition: mesh.hpp:1925
Vector X
Mapped (physical) coordinates of all quadrature points.
Definition: mesh.hpp:1869
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)
Evaluate the vector coefficient in the element described by T at the point ip, storing the result in ...
Definition: mesh.cpp:12160
static Mesh MakeSimplicial(const Mesh &orig_mesh)
Definition: mesh.cpp:4562
Base class for vector Coefficients that optionally depend on time and space.
void UniformRefinement3D_base(Array< int > *f2qf=NULL, DSTable *v_to_v_p=NULL, bool update_nodes=true)
Definition: mesh.cpp:8165
const GeometricFactors * GetGeometricFactors(const IntegrationRule &ir, const int flags, MemoryType d_mt=MemoryType::DEFAULT)
Return the mesh geometric factors corresponding to the given integration rule.
Definition: mesh.cpp:840
bool IsNonconformingCoarse() const
Return true if the face is a nonconforming coarse face.
Definition: mesh.hpp:1435
static void PrintElementsByGeometry(int dim, const Array< int > &num_elems_by_geom, std::ostream &out)
Auxiliary method used by PrintCharacteristics().
Definition: mesh.cpp:228
void AddHexAsTets(const int *vi, int attr=1)
Definition: mesh.cpp:1753
static const int NumGeom
Definition: geom.hpp:42
void AddHexAsWedges(const int *vi, int attr=1)
Definition: mesh.cpp:1772
virtual void PrintInfo(std::ostream &os=mfem::out)
In serial, this method calls PrintCharacteristics(). In parallel, additional information about the pa...
Definition: mesh.hpp:1781
virtual void Eval(Vector &V, ElementTransformation &T, const IntegrationPoint &ip)=0
Evaluate the vector coefficient in the element described by T at the point ip, storing the result in ...
Array< Element * > boundary
Definition: mesh.hpp:91
Geometry::Constants< Geometry::SQUARE > quad_t
Definition: mesh.hpp:259
int * GeneratePartitioning(int nparts, int part_method=1)
Definition: mesh.cpp:7033
void GetBdrElementFace(int i, int *f, int *o) const
Return the index and the orientation of the face of bdr element i. (3D)
Definition: mesh.cpp:6171
CoarseFineTransformations CoarseFineTr
Definition: mesh.hpp:235
int own_nodes
Definition: mesh.hpp:241
int GetNumGeometries(int dim) const
Return the number of geometries of the given dimension present in the mesh.
Definition: mesh.cpp:5908
virtual bool HasBoundaryElements() const
Checks if the mesh has boundary elements.
Definition: mesh.hpp:1023
bool Conforming() const
Definition: mesh.hpp:1650
void MoveVertices(const Vector &displacements)
Definition: mesh.cpp:7828
int GetNBE() const
Returns number of boundary elements.
Definition: mesh.hpp:926
IsoparametricTransformation Transformation
Definition: mesh.hpp:229
void GetLocalPtToSegTransformation(IsoparametricTransformation &, int)
Used in GetFaceElementTransformations (...)
Definition: mesh.cpp:641
void GetCharacteristics(double &h_min, double &h_max, double &kappa_min, double &kappa_max, Vector *Vh=NULL, Vector *Vk=NULL)
Definition: mesh.cpp:193
static FiniteElement * GetTransformationFEforElementType(Element::Type)
Definition: mesh.cpp:327
FaceGeometricFactors(const Mesh *mesh, const IntegrationRule &ir, int flags, FaceType type, MemoryType d_mt=MemoryType::DEFAULT)
Definition: mesh.cpp:12093
static Mesh MakeCartesian3D(int nx, int ny, int nz, Element::Type type, double sx=1.0, double sy=1.0, double sz=1.0, bool sfc_ordering=true)
Definition: mesh.cpp:3723
void NewNodes(GridFunction &nodes, bool make_owner=false)
Replace the internal node GridFunction with the given GridFunction.
Definition: mesh.cpp:7935
int NumOfEdges
Definition: mesh.hpp:70
virtual void UniformRefinement2D()
Refine a mixed 2D mesh uniformly.
Definition: mesh.hpp:392
void SwapNodes(GridFunction *&nodes, int &own_nodes_)
Definition: mesh.cpp:7957
Element::Type GetElementType(int i) const
Returns the type of element i.
Definition: mesh.cpp:6250
void GetBoundingBox(Vector &min, Vector &max, int ref=2)
Returns the minimum and maximum corners of the mesh bounding box.
Definition: mesh.cpp:129
void ReadNetgen2DMesh(std::istream &input, int &curved)
void ShiftRight(int &a, int &b, int &c)
Definition: mesh.hpp:1966
Element::Type GetBdrElementType(int i) const
Returns the type of boundary element i.
Definition: mesh.cpp:6255
static Mesh MakeRefined(Mesh &orig_mesh, int ref_factor, int ref_type)
Create a refined (by any factor) version of orig_mesh.
Definition: mesh.cpp:3733
int EulerNumber2D() const
Equals 1 - num_holes.
Definition: mesh.hpp:1003
A specialized ElementTransformation class representing a face and its two neighboring elements...
Definition: eltrans.hpp:480
void CreateVTKMesh(const Vector &points, const Array< int > &cell_data, const Array< int > &cell_offsets, const Array< int > &cell_types, const Array< int > &cell_attributes, int &curved, int &read_gf, bool &finalize_topo)
int AddTriangle(int v1, int v2, int v3, int attr=1)
Definition: mesh.cpp:1660
Vector detJ
Determinants of the Jacobians at all quadrature points.
Definition: mesh.hpp:1931
int AddBdrPoint(int v, int attr=1)
Definition: mesh.cpp:1880
void GetGeometries(int dim, Array< Geometry::Type > &el_geoms) const
Return all element geometries of the given dimension present in the mesh.
Definition: mesh.cpp:5919
void DeleteTables()
Definition: mesh.hpp:290
void GetFaceVertices(int i, Array< int > &vert) const
Returns the indices of the vertices of face i.
Definition: mesh.hpp:1127
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
GridFunction * Nodes
Definition: mesh.hpp:240
int NumOfElements
Definition: mesh.hpp:69
Mesh * Extrude1D(Mesh *mesh, const int ny, const double sy, const bool closed)
Extrude a 1D mesh.
Definition: mesh.cpp:12178
void Transform(void(*f)(const Vector &, Vector &))
Definition: mesh.cpp:11629
void FinalizeCheck()
Definition: mesh.cpp:1925
void ApplyLocalSlaveTransformation(FaceElementTransformations &FT, const FaceInfo &fi, bool is_ghost)
Definition: mesh.cpp:1056
int GetBdrElementEdgeIndex(int i) const
Definition: mesh.cpp:6194
void AverageVertices(const int *indexes, int n, int result)
Averages the vertices with given indexes and saves the result in vertices[result].
Definition: mesh.cpp:7970
int GetNE() const
Returns number of elements.
Definition: mesh.hpp:923
const Element * GetFace(int i) const
Definition: mesh.hpp:1048
void FinalizeWedgeMesh(int generate_edges=0, int refine=0, bool fix_orientation=true)
Finalize the construction of a wedge Mesh.
Definition: mesh.cpp:2818
Element * GetElement(int i)
Definition: mesh.hpp:1042
void GenerateBoundaryElements()
Definition: mesh.cpp:1887
Geometry::Type GetFaceBaseGeometry(int i) const
Definition: mesh.hpp:1066
void PrintCharacteristics(Vector *Vh=NULL, Vector *Vk=NULL, std::ostream &os=mfem::out)
Compute and print mesh characteristics such as number of vertices, number of elements, number of boundary elements, minimal and maximal element sizes, minimal and maximal element aspect ratios, etc.
Definition: mesh.cpp:242
Structure for storing mesh geometric factors: coordinates, Jacobians, and determinants of the Jacobia...
Definition: mesh.hpp:1837
NodeExtrudeCoefficient(const int dim, const int n_, const double s_)
Definition: mesh.cpp:12154
Data type for vertex.
Definition: vertex.hpp:22
void SetMeshGen()
Determine the mesh generator bitmask meshgen, see MeshGenerator().
Definition: mesh.cpp:3913
static int GetQuadOrientation(const int *base, const int *test)
Returns the orientation of &quot;test&quot; relative to &quot;base&quot;.
Definition: mesh.cpp:5616
void GetVertices(Vector &vert_coord) const
Definition: mesh.cpp:7837
void ReadNetgen3DMesh(std::istream &input)
Data arrays will be written in ASCII format.
The inverse transformation of a given ElementTransformation.
Definition: eltrans.hpp:185
void ReadInlineMesh(std::istream &input, bool generate_edges=false)
void GetBdrElementAdjacentElement2(int bdr_el, int &el, int &info) const
For the given boundary element, bdr_el, return its adjacent element and its info, i...
Definition: mesh.cpp:6228
FaceElementTransformations * GetInteriorFaceTransformations(int FaceNo)
Definition: mesh.hpp:1268
ElementLocation
Definition: mesh.hpp:1302
void RemoveInternalBoundaries()
Definition: mesh.cpp:11785
ElementConformity
Definition: mesh.hpp:1317
void DebugDump(std::ostream &out) const
Output an NCMesh-compatible debug dump.
Definition: mesh.cpp:12565
Mesh * Extrude2D(Mesh *mesh, const int nz, const double sz)
Extrude a 2D mesh.
Definition: mesh.cpp:12338
virtual FaceElementTransformations * GetFaceElementTransformations(int FaceNo, int mask=31)
Definition: mesh.cpp:961
Geometry::Type GetElementBaseGeometry(int i) const
Definition: mesh.hpp:1069
void SetEmpty()
Definition: mesh.cpp:1472
Array< Element * > faces
Definition: mesh.hpp:92
void FinalizeHexMesh(int generate_edges=0, int refine=0, bool fix_orientation=true)
Finalize the construction of a hexahedral Mesh.
Definition: mesh.cpp:2853
Geometry::Type GetBdrElementGeometry(int i) const
Definition: mesh.hpp:1060
void GetVertexToVertexTable(DSTable &) const
Definition: mesh.cpp:6331
void RedRefinement(int i, const DSTable &v_to_v, int *edge1, int *edge2, int *middle)
Definition: mesh.hpp:354
Operation last_operation
Definition: mesh.hpp:284
void MarkTriMeshForRefinement()
Definition: mesh.cpp:2418
void ReadCubit(const char *filename, int &curved, int &read_gf)
bool HasGeometry(Geometry::Type geom) const
Return true iff the given geom is encountered in the mesh. Geometries of dimensions lower than Dimens...
Definition: mesh.hpp:1077
void InitRefinementTransforms()
Definition: mesh.cpp:9906
IsoparametricTransformation FaceTransformation
Definition: mesh.hpp:231
void UniformRefinement2D_base(bool update_nodes=true)
Definition: mesh.cpp:8006
Array< FaceGeometricFactors * > face_geom_factors
Optional face geometric factors.
Definition: mesh.hpp:276
Array< NCFaceInfo > nc_faces_info
Definition: mesh.hpp:218
Element * ReadElement(std::istream &)
Definition: mesh.cpp:3895
void KnotInsert(Array< KnotVector * > &kv)
Definition: mesh.cpp:5050
void GetLocalTriToPyrTransformation(IsoparametricTransformation &loc, int i)
Used in GetFaceElementTransformations (...)
Definition: mesh.cpp:746
A parallel extension of the NCMesh class.
Definition: pncmesh.hpp:64
std::vector< int > CreatePeriodicVertexMapping(const std::vector< Vector > &translations, double tol=1e-8) const
Creates a mapping v2v from the vertex indices of the mesh such that coincident vertices under the giv...
Definition: mesh.cpp:4925
Element * NewElement(int geom)
Definition: mesh.cpp:3841
FaceElementTransformations FaceElemTr
Definition: mesh.hpp:232
Table * el_to_face
Definition: mesh.hpp:221
void GetEdgeTransformation(int i, IsoparametricTransformation *EdTr)
Definition: mesh.cpp:574
const GridFunction * GetNodes() const
Definition: mesh.hpp:1542
void MarkForRefinement()
Definition: mesh.cpp:2401
Element::Type GetFaceElementType(int Face) const
Definition: mesh.cpp:1432
void SetVerticesFromNodes(const GridFunction *nodes)
Helper to set vertex coordinates given a high-order curvature function.
Definition: mesh.cpp:5362
Structure for storing face geometric factors: coordinates, Jacobians, determinants of the Jacobians...
Definition: mesh.hpp:1891
ElementTransformation * GetBdrElementTransformation(int i)
Returns the transformation defining the i-th boundary element.
Definition: mesh.cpp:431
int AddBdrTriangle(int v1, int v2, int v3, int attr=1)
Definition: mesh.cpp:1837
void PrintTopo(std::ostream &out, const Array< int > &e_to_k) const
Definition: mesh.cpp:10247
int GetNumFaces() const
Return the number of faces (3D), edges (2D) or vertices (1D).
Definition: mesh.cpp:5376
void MakeRefined_(Mesh &orig_mesh, const Array< int > ref_factors, int ref_type)
Internal function used in Mesh::MakeRefined.
Definition: mesh.cpp:4333
Geometry::Constants< Geometry::SEGMENT > seg_t
Definition: mesh.hpp:257
bool IsShared() const
Return true if the face is a shared interior face which is NOT a master nonconforming face...
Definition: mesh.hpp:1382
void GetElementJacobian(int i, DenseMatrix &J)
Definition: mesh.cpp:60
static int GetTetOrientation(const int *base, const int *test)
Returns the orientation of &quot;test&quot; relative to &quot;base&quot;.
Definition: mesh.cpp:5664
long GetSequence() const
Definition: mesh.hpp:1664
MFEM_DEPRECATED Mesh(int nx, int ny, int nz, Element::Type type, bool generate_edges=false, double sx=1.0, double sy=1.0, double sz=1.0, bool sfc_ordering=true)
Deprecated: see MakeCartesian3D.
Definition: mesh.hpp:841
int AddTri(const int *vi, int attr=1)
Definition: mesh.hpp:716
int AddVertex(double x, double y=0.0, double z=0.0)
Definition: mesh.cpp:1613
bool Nonconforming() const
Definition: mesh.hpp:1651
Vector J
Jacobians of the element transformations at all quadrature points.
Definition: mesh.hpp:1878
FaceType
Definition: mesh.hpp:45
void MoveNodes(const Vector &displacements)
Definition: mesh.cpp:7896
void GetEdgeOrdering(DSTable &v_to_v, Array< int > &order)
Definition: mesh.cpp:2433
void GetElementData(int geom, Array< int > &elem_vtx, Array< int > &attr) const
Definition: mesh.hpp:1019
void PrintWithPartitioning(int *partitioning, std::ostream &os, int elem_attr=0) const
Prints the mesh with boundary elements given by the boundary of the subdomains, so that the boundary ...
Definition: mesh.cpp:10930
bool IsInterior() const
return true if the face is an interior face to the computation domain, either a local or shared inter...
Definition: mesh.hpp:1391
NCFaceInfo(bool slave, int master, const DenseMatrix *pm)
Definition: mesh.hpp:213
void FinalizeTriMesh(int generate_edges=0, int refine=0, bool fix_orientation=true)
Finalize the construction of a triangular Mesh.
Definition: mesh.cpp:1939
double f(const Vector &xvec)
Definition: lor_mms.hpp:32
double GetElementSize(ElementTransformation *T, int type=0)
Definition: mesh.cpp:76
void UpdateNURBS()
Definition: mesh.cpp:5124
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 PrintVTK(std::ostream &os)
Definition: mesh.cpp:10303
void FinalizeTetMesh(int generate_edges=0, int refine=0, bool fix_orientation=true)
Finalize the construction of a tetrahedral Mesh.
Definition: mesh.cpp:2777
void AddBdrQuadAsTriangles(const int *vi, int attr=1)
Definition: mesh.cpp:1865
void Make3D(int nx, int ny, int nz, Element::Type type, double sx, double sy, double sz, bool sfc_ordering)
Definition: mesh.cpp:3058
void AddSegmentFaceElement(int lf, int gf, int el, int v0, int v1)
Definition: mesh.cpp:6481
virtual void SetAttributes()
Definition: mesh.cpp:1559
void EnsureNCMesh(bool simplices_nonconforming=false)
Definition: mesh.cpp:9473
GeometryList(Mesh &mesh, int dim)
Construct a GeometryList of all geometries of dimension dim in mesh.
Definition: mesh.hpp:1103
Vector detJ
Determinants of the Jacobians at all quadrature points.
Definition: mesh.hpp:1884
void AddHexAsPyramids(const int *vi, int attr=1)
Definition: mesh.cpp:1790
int FindCoarseElement(int i)
Definition: mesh.cpp:9918
GeometryList(Mesh &mesh)
Construct a GeometryList of all element geometries in mesh.
Definition: mesh.hpp:1098
double b
Definition: lissajous.cpp:42
int AddTet(int v1, int v2, int v3, int v4, int attr=1)
Definition: mesh.cpp:1688
void GetFaceTransformation(int i, IsoparametricTransformation *FTr)
Returns the transformation defining the given face element in a user-defined variable.
Definition: mesh.cpp:497
void GetLocalSegToTriTransformation(IsoparametricTransformation &loc, int i)
Definition: mesh.cpp:656
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
Definition: mesh.cpp:9816
void AddQuadFaceElement(int lf, int gf, int el, int v0, int v1, int v2, int v3)
Definition: mesh.cpp:6546
void CheckPartitioning(int *partitioning_)
Definition: mesh.cpp:7433
void LoadPatchTopo(std::istream &input, Array< int > &edge_to_knot)
Read NURBS patch/macro-element mesh.
Definition: mesh.cpp:5189
void ReadGmshMesh(std::istream &input, int &curved, int &read_gf)
STable3D * GetElementToFaceTable(int ret_ftbl=0)
Definition: mesh.cpp:6798
virtual void SetCurvature(int order, bool discont=false, int space_dim=-1, int ordering=1)
Definition: mesh.cpp:5343
int MeshGenerator()
Get the mesh generator/type.
Definition: mesh.hpp:916
Type
Constants for the classes derived from Element.
Definition: element.hpp:41
void GetBdrPointMatrix(int i, DenseMatrix &pointmat) const
Definition: mesh.cpp:6278
void InitTables()
Definition: mesh.cpp:1466
Geometry::Type GetFaceGeometry(int i) const
Definition: mesh.hpp:1050
void CheckDisplacements(const Vector &displacements, double &tmax)
Definition: mesh.cpp:7751
This structure stores the low level information necessary to interpret the configuration of elements ...
Definition: mesh.hpp:153
VTKFormat
Data array format for VTK and VTU files.
Definition: vtk.hpp:96
int mesh_geoms
Definition: mesh.hpp:78
int AddBdrSegment(int v1, int v2, int attr=1)
Definition: mesh.cpp:1823
Geometry::Constants< Geometry::TETRAHEDRON > tet_t
Definition: mesh.hpp:260
int AddElement(Element *elem)
The parameter elem should be allocated using the NewElement() method.
Definition: mesh.cpp:1809
void ReadLineMesh(std::istream &input)
virtual void Load(std::istream &input, int generate_edges=0, int refine=1, bool fix_orientation=true)
Definition: mesh.hpp:891
void SetNodesOwner(bool nodes_owner)
Set the mesh nodes ownership flag.
Definition: mesh.hpp:1546
void SetLayer(const int l)
Definition: mesh.hpp:1950
int GetElementToEdgeTable(Table &, Array< int > &)
Definition: mesh.cpp:6356
void GetBdrElementAdjacentElement(int bdr_el, int &el, int &info) const
For the given boundary element, bdr_el, return its adjacent element and its info, i...
Definition: mesh.cpp:6206
bool IsBoundary() const
Return true if the face is a boundary face.
Definition: mesh.hpp:1398
void GetHilbertElementOrdering(Array< int > &ordering)
Definition: mesh.cpp:2197
int nbBoundaryFaces
Definition: mesh.hpp:75
void RandomRefinement(double prob, bool aniso=false, int nonconforming=-1, int nc_limit=0)
Refine each element with given probability. Uses GeneralRefinement.
Definition: mesh.cpp:9497
const Element * GetElement(int i) const
Definition: mesh.hpp:1040
void GetLocalQuadToHexTransformation(IsoparametricTransformation &loc, int i)
Used in GetFaceElementTransformations (...)
Definition: mesh.cpp:771
void ResetLazyData()
Definition: mesh.cpp:1549
static const int vtk_quadratic_pyramid[13]
Definition: mesh.hpp:244
void GetElementCenter(int i, Vector &center)
Definition: mesh.cpp:68
IsoparametricTransformation BdrTransformation
Definition: mesh.hpp:230
void GetLocalFaceTransformation(int face_type, int elem_type, IsoparametricTransformation &Transf, int info)
A helper method that constructs a transformation from the reference space of a face to the reference ...
Definition: mesh.cpp:896
IsoparametricTransformation Transformation2
Definition: mesh.hpp:229
void Init()
Definition: mesh.cpp:1448
virtual void SetNodalFESpace(FiniteElementSpace *nfes)
Definition: mesh.cpp:5285
Element * GetBdrElement(int i)
Definition: mesh.hpp:1046
Vector normal
Normals at all quadrature points.
Definition: mesh.hpp:1938
A class for non-conforming AMR. The class is not used directly by the user, rather it is an extension...
Definition: ncmesh.hpp:121
FaceElementTransformations * GetBdrFaceTransformations(int BdrElemNo)
Definition: mesh.cpp:1095
Geometry::Constants< Geometry::PRISM > pri_t
Definition: mesh.hpp:262
int Dimension() const
Definition: mesh.hpp:1006
struct mfem::Mesh::FaceInformation::@13 element[2]
virtual MFEM_DEPRECATED void ReorientTetMesh()
Definition: mesh.cpp:6927
prob_type prob
Definition: ex25.cpp:153
int NumOfBdrElements
Definition: mesh.hpp:69
Table * el_to_edge
Definition: mesh.hpp:220
void ReadMFEMMesh(std::istream &input, int version, int &curved)
virtual void LocalRefinement(const Array< int > &marked_el, int type=3)
This function is not public anymore. Use GeneralRefinement instead.
Definition: mesh.cpp:8853
Data type tetrahedron element.
Definition: tetrahedron.hpp:22
virtual void MarkTetMeshForRefinement(DSTable &v_to_v)
Definition: mesh.cpp:2458
Array< Triple< int, int, int > > tmp_vertex_parents
Definition: mesh.hpp:254
bool IsOfFaceType(FaceType type) const
Return true if the face is of the same type as type.
Definition: mesh.hpp:1404
MFEM_DEPRECATED Mesh(int nx, int ny, Element::Type type, bool generate_edges=false, double sx=1.0, double sy=1.0, bool sfc_ordering=true)
Deprecated: see MakeCartesian2D.
Definition: mesh.hpp:851
List of mesh geometries stored as Array&lt;Geometry::Type&gt;.
Definition: mesh.hpp:1092
void Destroy()
Definition: mesh.cpp:1520
void PrintVTU(std::ostream &os, int ref=1, VTKFormat format=VTKFormat::ASCII, bool high_order_output=false, int compression_level=0, bool bdr_elements=false)
Definition: mesh.cpp:10503
This structure is used as a human readable output format that decipheres the information contained in...
Definition: mesh.hpp:1356
int CheckBdrElementOrientation(bool fix_it=true)
Check the orientation of the boundary elements.
Definition: mesh.cpp:5797
int SpaceDimension() const
Definition: mesh.hpp:1007
GeometricFactors(const Mesh *mesh, const IntegrationRule &ir, int flags, MemoryType d_mt=MemoryType::DEFAULT)
Definition: mesh.cpp:12013
void PrintElementsWithPartitioning(int *partitioning, std::ostream &out, int interior_faces=0)
Definition: mesh.cpp:11049
bool FaceIsInterior(int FaceNo) const
Return true if the given face is interior.
Definition: mesh.hpp:1280
void GetLocalQuadToWdgTransformation(IsoparametricTransformation &loc, int i)
Used in GetFaceElementTransformations (...)
Definition: mesh.cpp:793
void AddPointFaceElement(int lf, int gf, int el)
Used in GenerateFaces()
Definition: mesh.cpp:6449
void GetFaceElements(int Face, int *Elem1, int *Elem2) const
Definition: mesh.cpp:1394
void Make1D(int n, double sx=1.0)
Creates a 1D mesh for the interval [0,sx] divided into n equal intervals.
Definition: mesh.cpp:3517
const Element *const * GetElementsArray() const
Definition: mesh.hpp:1037
Array< int > bdr_attributes
A list of all unique boundary attributes used by the Mesh.
Definition: mesh.hpp:270
virtual void NonconformingRefinement(const Array< Refinement > &refinements, int nc_limit=0)
This function is not public anymore. Use GeneralRefinement instead.
Definition: mesh.cpp:9095
virtual ~Mesh()
Destroys Mesh.
Definition: mesh.hpp:1820
void FinalizeMesh(int refine=0, bool fix_orientation=true)
Finalize the construction of any type of Mesh.
Definition: mesh.cpp:2883
Table * el_to_el
Definition: mesh.hpp:222
void GetFaceInfos(int Face, int *Inf1, int *Inf2) const
Definition: mesh.cpp:1400
int AddSegment(int v1, int v2, int attr=1)
Definition: mesh.cpp:1646
Geometry::Constants< Geometry::CUBE > hex_t
Definition: mesh.hpp:261
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:96
Geometry::Type geom_buf[Geometry::NumGeom]
Definition: mesh.hpp:1095
static Mesh LoadFromFile(const char *filename, int generate_edges=0, int refine=1, bool fix_orientation=true)
Definition: mesh.cpp:3695
void RemoveUnusedVertices()
Remove unused vertices and rebuild mesh connectivity.
Definition: mesh.cpp:11678
static void PrintElement(const Element *, std::ostream &)
Definition: mesh.cpp:3907
virtual void NURBSUniformRefinement()
Refine NURBS mesh.
Definition: mesh.cpp:5094
int AddHex(int v1, int v2, int v3, int v4, int v5, int v6, int v7, int v8, int attr=1)
Definition: mesh.cpp:1737
const IntegrationRule * IntRule
Definition: mesh.hpp:1846
MemoryType
Memory types supported by MFEM.
Definition: mem_manager.hpp:31
Vector X
Mapped (physical) coordinates of all quadrature points.
Definition: mesh.hpp:1916
bool IsLocal() const
Return true if the face is a local interior face which is NOT a master nonconforming face...
Definition: mesh.hpp:1375
void FinalizeTopology(bool generate_bdr=true)
Finalize the construction of the secondary topology (connectivity) data of a Mesh.
Definition: mesh.cpp:2890
static void GetElementArrayEdgeTable(const Array< Element * > &elem_array, const DSTable &v_to_v, Table &el_to_edge)
Definition: mesh.cpp:6309
GridFunction * GetNodes()
Return a pointer to the internal node GridFunction (may be NULL).
Definition: mesh.hpp:1541
void SetSize(int nsize)
Change the logical size of the array, keep existing entries.
Definition: array.hpp:679
void PrepareNodeReorder(DSTable **old_v_to_v, Table **old_elem_vert)
Definition: mesh.cpp:2481
int AddBdrQuad(int v1, int v2, int v3, int v4, int attr=1)
Definition: mesh.cpp:1851
bool IsSlaveFace(const FaceInfo &fi) const
Definition: mesh.cpp:1051
Array< Vertex > vertices
Definition: mesh.hpp:90
static Mesh MakeCartesian2D(int nx, int ny, Element::Type type, bool generate_edges=false, double sx=1.0, double sy=1.0, bool sfc_ordering=true)
Definition: mesh.cpp:3713
void InitMesh(int Dim_, int spaceDim_, int NVert, int NElem, int NBdrElem)
Begin construction of a mesh.
Definition: mesh.cpp:1590
void SetBdrAttribute(int i, int attr)
Set the attribute of boundary element i.
Definition: mesh.hpp:1489
void ReadXML_VTKMesh(std::istream &input, int &curved, int &read_gf, bool &finalize_topo, const std::string &xml_prefix="")
Geometry::Type GetElementGeometry(int i) const
Definition: mesh.hpp:1055
void RefineAtVertex(const Vertex &vert, double eps=0.0, int nonconforming=-1)
Refine elements sharing the specified vertex. Uses GeneralRefinement.
Definition: mesh.cpp:9516
static Mesh MakeCartesian1D(int n, double sx=1.0)
Definition: mesh.cpp:3705
Array< Element * > elements
Definition: mesh.hpp:85
const Table & ElementToElementTable()
Definition: mesh.cpp:6395
double GetLength(int i, int j) const
Return the length of the segment from node i to node j.
Definition: mesh.cpp:6294
int meshgen
Definition: mesh.hpp:77
const CoarseFineTransformations & GetRefinementTransforms()
Definition: mesh.cpp:9928
bool RefineByError(const Array< double > &elem_error, double threshold, int nonconforming=-1, int nc_limit=0)
Definition: mesh.cpp:9542
Table * bel_to_edge
Definition: mesh.hpp:224
Operation GetLastOperation() const
Return type of last modification of the mesh.
Definition: mesh.hpp:1658
Table * GetFaceToElementTable() const
Definition: mesh.cpp:6114
void GetLocalTriToTetTransformation(IsoparametricTransformation &loc, int i)
Used in GetFaceElementTransformations (...)
Definition: mesh.cpp:696
double a
Definition: lissajous.cpp:41
int EulerNumber() const
Equals 1 + num_holes - num_loops.
Definition: mesh.hpp:1000
NURBSExtension * NURBSext
Optional NURBS mesh extension.
Definition: mesh.hpp:272
int GetNV() const
Returns number of vertices. Vertices are only at the corners of elements, where you would expect them...
Definition: mesh.hpp:920
int nbInteriorFaces
Definition: mesh.hpp:75
void PrintBdrVTU(std::string fname, VTKFormat format=VTKFormat::ASCII, bool high_order_output=false, int compression_level=0)
Definition: mesh.cpp:10495
void GetNode(int i, double *coord) const
Definition: mesh.cpp:7857
Array< int > be_to_edge
Definition: mesh.hpp:223
const Table & ElementToFaceTable() const
Definition: mesh.cpp:6431
Class for integration point with weight.
Definition: intrules.hpp:25
void Swap(Mesh &other, bool non_geometry)
Definition: mesh.cpp:9285
Element * ReadElementWithoutAttr(std::istream &)
Definition: mesh.cpp:3865
int GetNumFacesWithGhost() const
Return the number of faces (3D), edges (2D) or vertices (1D) including ghost faces.
Definition: mesh.cpp:5387
A standard isoparametric element transformation.
Definition: eltrans.hpp:361
void GetElementEdges(int i, Array< int > &edges, Array< int > &cor) const
Return the indices and the orientations of all edges of element i.
Definition: mesh.cpp:5932
static const int vtk_quadratic_wedge[18]
Definition: mesh.hpp:245
bool OwnsNodes() const
Return the mesh nodes ownership flag.
Definition: mesh.hpp:1544
void GetElementTransformation(int i, IsoparametricTransformation *ElTr)
Definition: mesh.cpp:348
void GetLocalSegToQuadTransformation(IsoparametricTransformation &loc, int i)
Definition: mesh.cpp:676
const FiniteElementSpace * GetNodalFESpace() const
Definition: mesh.cpp:5338
virtual void Print(std::ostream &os=mfem::out) const
Definition: mesh.hpp:1671
Table * face_edge
Definition: mesh.hpp:226
void FinalizeQuadMesh(int generate_edges=0, int refine=0, bool fix_orientation=true)
Finalize the construction of a quadrilateral Mesh.
Definition: mesh.cpp:1968
void GetBdrElementData(int geom, Array< int > &bdr_elem_vtx, Array< int > &bdr_attr) const
Definition: mesh.hpp:1025
virtual void Finalize(bool refine=false, bool fix_orientation=false)
Finalize the construction of a general Mesh.
Definition: mesh.cpp:2983
bool DerefineByError(Array< double > &elem_error, double threshold, int nc_limit=0, int op=1)
Definition: mesh.cpp:9212
void ReorderElements(const Array< int > &ordering, bool reorder_vertices=true)
Definition: mesh.cpp:2249
void GetFaceEdges(int i, Array< int > &edges, Array< int > &o) const
Definition: mesh.cpp:5986
void GetLocalQuadToPyrTransformation(IsoparametricTransformation &loc, int i)
Used in GetFaceElementTransformations (...)
Definition: mesh.cpp:817
long long GetGlobalNE() const
Return the total (global) number of elements.
Definition: mesh.hpp:955
void DegreeElevate(int rel_degree, int degree=16)
Definition: mesh.cpp:5107
int dim
Definition: ex24.cpp:53
Table * edge_vertex
Definition: mesh.hpp:227
long sequence
Definition: mesh.hpp:83
virtual ~NodeExtrudeCoefficient()
Definition: mesh.hpp:1954
constexpr int dimension
This example only works in 3D. Kernels for 2D are not implemented.
Definition: hooke.cpp:45
const Mesh * mesh
Definition: mesh.hpp:1845
void MakeSimplicial_(const Mesh &orig_mesh, int *vglobal)
Definition: mesh.cpp:4569
void DestroyPointers()
Definition: mesh.cpp:1494
Array< FaceInfo > faces_info
Definition: mesh.hpp:217
STable3D * GetFacesTable()
Definition: mesh.cpp:6735
std::ostream & operator<<(std::ostream &os, SparseMatrix const &mat)
Definition: sparsemat.hpp:710
NCMesh * ncmesh
Optional nonconforming mesh extension.
Definition: mesh.hpp:273
Geometry::Type GetBdrElementBaseGeometry(int i) const
Definition: mesh.hpp:1072
virtual void UniformRefinement3D()
Refine a mixed 3D mesh uniformly.
Definition: mesh.hpp:402
int AddBdrElement(Element *elem)
Definition: mesh.cpp:1816
RefCoord t[3]
static void PrintElementWithoutAttr(const Element *, std::ostream &)
Definition: mesh.cpp:3883
Mesh(int Dim_, int NVert, int NElem, int NBdrElem=0, int spaceDim_=-1)
Init constructor: begin the construction of a Mesh object.
Definition: mesh.hpp:692
int AddWedge(int v1, int v2, int v3, int v4, int v5, int v6, int attr=1)
Definition: mesh.cpp:1709
virtual bool NonconformingDerefinement(Array< double > &elem_error, double threshold, int nc_limit=0, int op=1)
NC version of GeneralDerefinement.
Definition: mesh.cpp:9164
int Dim
Definition: mesh.hpp:66
Geometry::Type GetFaceGeometryType(int Face) const
Definition: mesh.cpp:1413
void Clear()
Clear the contents of the Mesh.
Definition: mesh.hpp:899
int GetNEdges() const
Return the number of edges.
Definition: mesh.hpp:929
ElementConformity conformity
Definition: mesh.hpp:1363
Geometry::Constants< Geometry::TRIANGLE > tri_t
Definition: mesh.hpp:258
double AggregateError(const Array< double > &elem_error, const int *fine, int nfine, int op)
Derefinement helper.
Definition: mesh.cpp:9145
double * GetVertex(int i)
Return pointer to vertex i&#39;s coordinates.
Definition: mesh.hpp:1017
int GetNFaces() const
Return the number of faces in a 3D mesh.
Definition: mesh.hpp:932
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: mesh.cpp:11898
void SetAttribute(int i, int attr)
Set the attribute of element i.
Definition: mesh.hpp:1483
Mesh & operator=(Mesh &&mesh)
Move assignment operstor.
Definition: mesh.cpp:3689
virtual void PrintXG(std::ostream &os=mfem::out) const
Print the mesh to the given stream using Netgen/Truegrid format.
Definition: mesh.cpp:9991
Vector data type.
Definition: vector.hpp:60
void ReadTrueGridMesh(std::istream &input)
const FaceGeometricFactors * GetFaceGeometricFactors(const IntegrationRule &ir, const int flags, FaceType type, MemoryType d_mt=MemoryType::DEFAULT)
Return the mesh geometric factors for the faces corresponding to the given integration rule...
Definition: mesh.cpp:860
Array< int > GetFaceToBdrElMap() const
Definition: mesh.cpp:1437
virtual int GetNFbyType(FaceType type) const
Returns the number of faces according to the requested type, does not count master nonconforming face...
Definition: mesh.cpp:5392
static int GetTriOrientation(const int *base, const int *test)
Returns the orientation of &quot;test&quot; relative to &quot;base&quot;.
Definition: mesh.cpp:5568
void SetNode(int i, const double *coord)
Definition: mesh.cpp:7876
virtual long long ReduceInt(int value) const
Utility function: sum integers from all processors (Allreduce).
Definition: mesh.hpp:952
int GetBdrFace(int BdrElemNo) const
Return the local face index for the given boundary face.
Definition: mesh.cpp:1113
void DestroyTables()
Definition: mesh.cpp:1478
void GetLocalTriToWdgTransformation(IsoparametricTransformation &loc, int i)
Used in GetFaceElementTransformations (...)
Definition: mesh.cpp:720
void GenerateFaces()
Definition: mesh.cpp:6573
void GetElementColoring(Array< int > &colors, int el0=0)
Definition: mesh.cpp:10855
void ChangeVertexDataOwnership(double *vertices, int len_vertices, bool zerocopy=false)
Set the internal Vertex array to point to the given vertices array without assuming ownership of the ...
Definition: mesh.cpp:3775
void DeleteGeometricFactors()
Destroy all GeometricFactors stored by the Mesh.
Definition: mesh.cpp:882
int spaceDim
Definition: mesh.hpp:67
IsoparametricTransformation EdgeTransformation
Definition: mesh.hpp:231
Defines the coarse-fine transformations of all fine elements.
Definition: ncmesh.hpp:70
MemAlloc< Tetrahedron, 1024 > TetMemory
Definition: mesh.hpp:250
Table * GetFaceEdgeTable() const
Returns the face-to-edge Table (3D)
Definition: mesh.cpp:6025
void DoNodeReorder(DSTable *old_v_to_v, Table *old_elem_vert)
Definition: mesh.cpp:2547
void SetNodes(const Vector &node_coord)
Definition: mesh.cpp:7920
void SetNodalGridFunction(GridFunction *nodes, bool make_owner=false)
Definition: mesh.cpp:5332
void ReadNURBSMesh(std::istream &input, int &curved, int &read_gf)
int NumOfVertices
Definition: mesh.hpp:69
void ScaleSubdomains(double sf)
Definition: mesh.cpp:11489
Array< int > be_to_face
Definition: mesh.hpp:225
int AddPyramid(int v1, int v2, int v3, int v4, int v5, int attr=1)
Definition: mesh.cpp:1723
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
void InitFromNCMesh(const NCMesh &ncmesh)
Initialize vertices/elements/boundary/tables from a nonconforming mesh.
Definition: mesh.cpp:9241
void Bisection(int i, const DSTable &, int *, int *, int *)
Bisect a triangle: element with index i is bisected.
Definition: mesh.cpp:9571
void GeneralRefinement(const Array< Refinement > &refinements, int nonconforming=-1, int nc_limit=0)
Definition: mesh.cpp:9405
ElementLocation location
Definition: mesh.hpp:1362
Class for parallel meshes.
Definition: pmesh.hpp:32
Abstract data type element.
Definition: element.hpp:28
int GetAttribute(int i) const
Return the attribute of element i.
Definition: mesh.hpp:1480
void AddTriangleFaceElement(int lf, int gf, int el, int v0, int v1, int v2)
Definition: mesh.cpp:6518
const Table & ElementToEdgeTable() const
Definition: mesh.cpp:6440
Array< GeometricFactors * > geom_factors
Optional geometric factors.
Definition: mesh.hpp:274
void BdrBisection(int i, const HashTable< Hashed2 > &)
Bisect a boundary triangle: boundary element with index i is bisected.
Definition: mesh.cpp:9779
void GetBdrElementVertices(int i, Array< int > &v) const
Returns the indices of the vertices of boundary element i.
Definition: mesh.hpp:1113
void GenerateNCFaceInfo()
Definition: mesh.cpp:6674
void GetElementVertices(int i, Array< int > &v) const
Returns the indices of the vertices of element i.
Definition: mesh.hpp:1109
double GetElementVolume(int i)
Definition: mesh.cpp:112
Array< int > attributes
A list of all unique element attributes used by the Mesh.
Definition: mesh.hpp:268
Table * GetVertexToElementTable()
The returned Table must be destroyed by the caller.
Definition: mesh.cpp:6079
const Element * GetBdrElement(int i) const
Definition: mesh.hpp:1044
void EnsureNodes()
Definition: mesh.cpp:5291
void GetElementData(const Array< Element * > &elem_array, int geom, Array< int > &elem_vtx, Array< int > &attr) const
Definition: mesh.cpp:9341
const DenseMatrix * point_matrix
Definition: mesh.hpp:1371
void UpdateNodes()
Update the nodes of a curved mesh after refinement.
Definition: mesh.cpp:7991
void GetElementFaces(int i, Array< int > &faces, Array< int > &ori) const
Return the indices and the orientations of all faces of element i.
Definition: mesh.cpp:6148
void GreenRefinement(int i, const DSTable &v_to_v, int *edge1, int *edge2, int *middle)
Definition: mesh.hpp:360
bool IsConforming() const
Return true if the face is a conforming face.
Definition: mesh.hpp:1418
void Printer(std::ostream &out=mfem::out, std::string section_delimiter="") const
Definition: mesh.cpp:10155
void NodesUpdated()
This function should be called after the mesh node coordinates have changed, e.g. after the mesh has ...
Definition: mesh.hpp:997
Class used to extrude the nodes of a mesh.
Definition: mesh.hpp:1942
void AddVertexParents(int i, int p1, int p2)
Mark vertex i as nonconforming, with parent vertices p1 and p2.
Definition: mesh.cpp:1630
int NumOfFaces
Definition: mesh.hpp:70