MFEM  v4.6.0
Finite element discretization library
ncmesh.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, 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_NCMESH
13 #define MFEM_NCMESH
14 
15 #include "../config/config.hpp"
16 #include "../general/hash.hpp"
17 #include "../general/globals.hpp"
18 #include "../general/sort_pairs.hpp"
19 #include "../linalg/densemat.hpp"
20 #include "element.hpp"
21 #include "vertex.hpp"
22 #include "../fem/geom.hpp"
23 
24 #include <vector>
25 #include <map>
26 #include <iostream>
27 
28 namespace mfem
29 {
30 
31 /** Represents the index of an element to refine, plus a refinement type.
32  The refinement type is needed for anisotropic refinement of quads and hexes.
33  Bits 0,1 and 2 of 'ref_type' specify whether the element should be split
34  in the X, Y and Z directions, respectively (Z is ignored for quads). */
35 struct Refinement
36 {
37  enum : char { X = 1, Y = 2, Z = 4, XY = 3, XZ = 5, YZ = 6, XYZ = 7 };
38  int index; ///< Mesh element number
39  char ref_type; ///< refinement XYZ bit mask (7 = full isotropic)
40 
41  Refinement() = default;
42 
43  Refinement(int index, int type = Refinement::XYZ)
44  : index(index), ref_type(type) {}
45 };
46 
47 
48 /// Defines the position of a fine element within a coarse element.
49 struct Embedding
50 {
51  /// Coarse %Element index in the coarse mesh.
52  int parent;
53 
54  /** The (geom, matrix) pair determines the sub-element transformation for the
55  fine element: CoarseFineTransformations::point_matrices[geom](matrix) is
56  the point matrix of the region within the coarse element reference domain.*/
57  unsigned geom : 4;
58  unsigned matrix : 27;
59 
60  /// For internal use: 0 if regular fine element, 1 if parallel ghost element.
61  unsigned ghost : 1;
62 
63  Embedding() = default;
64  Embedding(int elem, Geometry::Type geom, int matrix = 0, bool ghost = false)
65  : parent(elem), geom(geom), matrix(matrix), ghost(ghost) {}
66 };
67 
68 
69 /// Defines the coarse-fine transformations of all fine elements.
71 {
72  /// Fine element positions in their parents.
74 
75  /** A "dictionary" of matrices for IsoparametricTransformation. Use
76  Embedding::{geom,matrix} to access a fine element point matrix. */
78 
79  /** Invert the 'embeddings' array: create a Table with coarse elements as
80  rows and fine elements as columns. If 'want_ghosts' is false, parallel
81  ghost fine elements are not included in the table. */
82  void MakeCoarseToFineTable(Table &coarse_to_fine,
83  bool want_ghosts = false) const;
84 
85  void Clear();
86  bool IsInitialized() const;
87  long MemoryUsage() const;
88 
89  MFEM_DEPRECATED
90  void GetCoarseToFineMap(const Mesh &fine_mesh, Table &coarse_to_fine) const
91  { MakeCoarseToFineTable(coarse_to_fine, true); (void) fine_mesh; }
92 };
93 
94 void Swap(CoarseFineTransformations &a, CoarseFineTransformations &b);
95 
96 struct MatrixMap; // for internal use
97 
98 
99 /** \brief A class for non-conforming AMR. The class is not used directly
100  * by the user, rather it is an extension of the Mesh class.
101  *
102  * In general, the class is used by MFEM as follows:
103  *
104  * 1. NCMesh is constructed from elements of an existing Mesh. The elements
105  * are copied and become roots of the refinement hierarchy.
106  *
107  * 2. Some elements are refined with the Refine() method. Both isotropic and
108  * anisotropic refinements of quads/hexes are supported.
109  *
110  * 3. A new Mesh is created from NCMesh containing the leaf elements.
111  * This new Mesh may have non-conforming (hanging) edges and faces and
112  * is the one seen by the user.
113  *
114  * 4. FiniteElementSpace asks NCMesh for a list of conforming, master and
115  * slave edges/faces and creates the conforming interpolation matrix P.
116  *
117  * 5. A continuous/conforming solution is obtained by solving P'*A*P x = P'*b.
118  *
119  * 6. Repeat from step 2.
120  */
121 class NCMesh
122 {
123 public:
124  //// Initialize with elements from an existing 'mesh'.
125  explicit NCMesh(const Mesh *mesh);
126 
127  /** Load from a stream. The id header is assumed to have been read already
128  from \param[in] input . \param[in] version is 10 for the v1.0 NC format,
129  or 1 for the legacy v1.1 format. \param[out] curved is set to 1 if the
130  curvature GridFunction follows after mesh data. \param[out] is_nc (again
131  treated as a boolean) is set to 0 if the legacy v1.1 format in fact
132  defines a conforming mesh. See Mesh::Loader for details. */
133  NCMesh(std::istream &input, int version, int &curved, int &is_nc);
134 
135  /// Deep copy of another instance.
136  NCMesh(const NCMesh &other);
137 
138  /// Copy assignment not supported
139  NCMesh& operator=(NCMesh&) = delete;
140 
141  virtual ~NCMesh();
142 
143  /// Return the dimension of the NCMesh.
144  int Dimension() const { return Dim; }
145  /// Return the space dimension of the NCMesh.
146  int SpaceDimension() const { return spaceDim; }
147 
148  /// Return the number of vertices in the NCMesh.
149  int GetNVertices() const { return NVertices; }
150  /// Return the number of edges in the NCMesh.
151  int GetNEdges() const { return NEdges; }
152  /// Return the number of (2D) faces in the NCMesh.
153  int GetNFaces() const { return NFaces; }
154  virtual int GetNGhostElements() const { return 0; }
155 
156  /** Perform the given batch of refinements. Please note that in the presence
157  of anisotropic splits additional refinements may be necessary to keep
158  the mesh consistent. However, the function always performs at least the
159  requested refinements. */
160  virtual void Refine(const Array<Refinement> &refinements);
161 
162  /** Check the mesh and potentially refine some elements so that the maximum
163  difference of refinement levels between adjacent elements is not greater
164  than 'max_nc_level'. */
165  virtual void LimitNCLevel(int max_nc_level);
166 
167  /** Return a list of derefinement opportunities. Each row of the table
168  contains Mesh indices of existing elements that can be derefined to form
169  a single new coarse element. Row numbers are then passed to Derefine.
170  This function works both in serial and parallel. */
171  const Table &GetDerefinementTable();
172 
173  /** Check derefinements returned by GetDerefinementTable and mark those that
174  can be done safely so that the maximum NC level condition is not violated.
175  On return, level_ok.Size() == deref_table.Size() and contains 0/1s. */
176  virtual void CheckDerefinementNCLevel(const Table &deref_table,
177  Array<int> &level_ok, int max_nc_level);
178 
179  /** Perform a subset of the possible derefinements (see GetDerefinementTable).
180  Note that if anisotropic refinements are present in the mesh, some of the
181  derefinements may have to be skipped to preserve mesh consistency. */
182  virtual void Derefine(const Array<int> &derefs);
183 
184  // master/slave lists
185 
186  /// Identifies a vertex/edge/face in both Mesh and NCMesh.
187  struct MeshId
188  {
189  int index; ///< Mesh number
190  int element; ///< NCMesh::Element containing this vertex/edge/face
191  signed char local; ///< local number within 'element'
192  signed char geom; ///< Geometry::Type (faces only) (char to save RAM)
193 
194  Geometry::Type Geom() const { return Geometry::Type(geom); }
195 
196  MeshId() = default;
197  MeshId(int index, int element, int local, int geom = -1)
199  };
200 
201  /** Nonconforming edge/face that has more than one neighbor. The neighbors
202  are stored in NCList::slaves[i], slaves_begin <= i < slaves_end. */
203  struct Master : public MeshId
204  {
205  int slaves_begin, slaves_end; ///< slave faces
206 
207  Master() = default;
208  Master(int index, int element, int local, int geom, int sb, int se)
210  , slaves_begin(sb), slaves_end(se) {}
211  };
212 
213  /// Nonconforming edge/face within a bigger edge/face.
214  struct Slave : public MeshId
215  {
216  int master; ///< master number (in Mesh numbering)
217  unsigned matrix : 24; ///< index into NCList::point_matrices[geom]
218  unsigned edge_flags : 8; ///< orientation flags, see OrientedPointMatrix
219 
220  Slave() = default;
221  Slave(int index, int element, int local, int geom)
223  , master(-1), matrix(0), edge_flags(0) {}
224  };
225 
226  /// Lists all edges/faces in the nonconforming mesh.
227  struct NCList
228  {
232 
233  /// List of unique point matrices for each slave geometry.
235 
236  /// Return the point matrix oriented according to the master and slave edges
237  void OrientedPointMatrix(const Slave &slave,
238  DenseMatrix &oriented_matrix) const;
239 
240  void Clear();
241  bool Empty() const { return !conforming.Size() && !masters.Size(); }
242  long TotalSize() const;
243  long MemoryUsage() const;
244 
245  const MeshId& LookUp(int index, int *type = NULL) const;
246 
247  ~NCList() { Clear(); }
248  private:
249  mutable Array<int> inv_index;
250  };
251 
252  /// Return the current list of conforming and nonconforming faces.
254  {
255  if (face_list.Empty()) { BuildFaceList(); }
256  return face_list;
257  }
258 
259  /// Return the current list of conforming and nonconforming edges.
261  {
262  if (edge_list.Empty()) { BuildEdgeList(); }
263  return edge_list;
264  }
265 
266  /** Return a list of vertices (in 'conforming'); this function is provided
267  for uniformity/completeness. Needed in ParNCMesh/ParFESpace. */
269  {
270  if (vertex_list.Empty()) { BuildVertexList(); }
271  return vertex_list;
272  }
273 
274  /// Return vertex/edge/face list (entity = 0/1/2, respectively).
275  const NCList& GetNCList(int entity)
276  {
277  switch (entity)
278  {
279  case 0: return GetVertexList();
280  case 1: return GetEdgeList();
281  default: return GetFaceList();
282  }
283  }
284 
285 
286  // coarse/fine transforms
287 
288  /** Remember the current layer of leaf elements before the mesh is refined.
289  Needed by GetRefinementTransforms(), must be called before Refine(). */
290  void MarkCoarseLevel();
291 
292  /** After refinement, calculate the relation of each fine element to its
293  parent coarse element. Note that Refine() or LimitNCLevel() can be called
294  multiple times between MarkCoarseLevel() and this function. */
296 
297  /** After derefinement, calculate the relations of previous fine elements
298  (some of which may no longer exist) to the current leaf elements.
299  Unlike for refinement, Derefine() may only be called once before this
300  function so there is no MarkFineLevel(). */
302 
303  /// Free all internal data created by the above three functions.
304  void ClearTransforms();
305 
306 
307  // grid ordering
308 
309  /** Return a space filling curve for a rectangular grid of elements.
310  Implemented is a generalized Hilbert curve for arbitrary grid dimensions.
311  If the width is odd, height should be odd too, otherwise one diagonal
312  (vertex-neighbor) step cannot be avoided in the curve. Even dimensions
313  are recommended. */
314  static void GridSfcOrdering2D(int width, int height,
315  Array<int> &coords);
316 
317  /** Return a space filling curve for a 3D rectangular grid of elements.
318  The Hilbert-curve-like algorithm works well for even dimensions. For odd
319  width/height/depth it tends to produce some diagonal (edge-neighbor)
320  steps. Even dimensions are recommended. */
321  static void GridSfcOrdering3D(int width, int height, int depth,
322  Array<int> &coords);
323 
324 
325  // utility
326 
327  /// Return Mesh vertex indices of an edge identified by 'edge_id'.
328  void GetEdgeVertices(const MeshId &edge_id, int vert_index[2],
329  bool oriented = true) const;
330 
331  /** Return "NC" orientation of an edge. As opposed to standard Mesh edge
332  orientation based on vertex IDs, "NC" edge orientation follows the local
333  edge orientation within the element 'edge_id.element' and is thus
334  processor independent. TODO: this seems only partially true? */
335  int GetEdgeNCOrientation(const MeshId &edge_id) const;
336 
337  /** Return Mesh vertex and edge indices of a face identified by 'face_id'.
338  The return value is the number of face vertices. */
339  int GetFaceVerticesEdges(const MeshId &face_id,
340  int vert_index[4], int edge_index[4],
341  int edge_orientation[4]) const;
342 
343  /** Given an edge (by its vertex indices v1 and v2) return the first
344  (geometric) parent edge that exists in the Mesh or -1 if there is no such
345  parent. */
346  int GetEdgeMaster(int v1, int v2) const;
347 
348  /** Get a list of vertices (2D/3D) and edges (3D) that coincide with boundary
349  elements with the specified attributes (marked in 'bdr_attr_is_ess').
350  In 3D this function also reveals "hidden" boundary edges. In parallel it
351  helps identifying boundary vertices/edges affected by non-local boundary
352  elements. */
353  virtual void GetBoundaryClosure(const Array<int> &bdr_attr_is_ess,
354  Array<int> &bdr_vertices,
355  Array<int> &bdr_edges);
356 
357  /// Return element geometry type. @a index is the Mesh element number.
359  { return elements[leaf_elements[index]].Geom(); }
360 
361  /// Return face geometry type. @a index is the Mesh face number.
363  { return Geometry::Type(face_geom[index]); }
364 
365  /// Return the number of root elements.
366  int GetNumRootElements() { return root_state.Size(); }
367 
368  /// Return the distance of leaf 'i' from the root.
369  int GetElementDepth(int i) const;
370 
371  /** Return the size reduction compared to the root element (ignoring local
372  stretching and curvature). */
373  int GetElementSizeReduction(int i) const;
374 
375  /// Return the faces and face attributes of leaf element 'i'.
377  Array<int> &fattr) const;
378 
379 
380  /// I/O: Print the mesh in "MFEM NC mesh v1.0" format.
381  void Print(std::ostream &out) const;
382 
383  /// I/O: Return true if the mesh was loaded from the legacy v1.1 format.
384  bool IsLegacyLoaded() const { return Legacy; }
385 
386  /// I/O: Return a map from old (v1.1) vertex indices to new vertex indices.
387  void LegacyToNewVertexOrdering(Array<int> &order) const;
388 
389  /// Save memory by releasing all non-essential and cached data.
390  virtual void Trim();
391 
392  /// Return total number of bytes allocated.
393  long MemoryUsage() const;
394 
395  int PrintMemoryDetail() const;
396 
397  typedef std::int64_t RefCoord;
398 
399 
400 protected: // non-public interface for the Mesh class
401 
402  friend class Mesh;
403 
404  /// Fill Mesh::{vertices,elements,boundary} for the current finest level.
405  void GetMeshComponents(Mesh &mesh) const;
406 
407  /** Get edge and face numbering from 'mesh' (i.e., set all Edge::index and
408  Face::index) after a new mesh was created from us. */
409  void OnMeshUpdated(Mesh *mesh);
410 
411  /** Delete top-level vertex coordinates if the Mesh became curved, e.g.,
412  by calling Mesh::SetCurvature or otherwise setting the Nodes. */
414 
415 
416 protected: // implementation
417 
418  int Dim, spaceDim; ///< dimensions of the elements and the vertex coordinates
419  int MyRank; ///< used in parallel, or when loading a parallel file in serial
420  bool Iso; ///< true if the mesh only contains isotropic refinements
421  int Geoms; ///< bit mask of element geometries present, see InitGeomFlags()
422  bool Legacy; ///< true if the mesh was loaded from the legacy v1.1 format
423 
424  static const int MaxElemNodes =
425  8; ///< Number of nodes of an element can have
426  static const int MaxElemEdges =
427  12; ///< Number of edges of an element can have
428  static const int MaxElemFaces =
429  6; ///< Number of faces of an element can have
430  static const int MaxElemChildren =
431  10; ///< Number of children of an element can have
432 
433  /** A Node can hold a vertex, an edge, or both. Elements directly point to
434  their corner nodes, but edge nodes also exist and can be accessed using
435  a hash-table given their two end-point node IDs. All nodes can be
436  accessed in this way, with the exception of top-level vertex nodes.
437  When an element is being refined, the mid-edge nodes are readily
438  available with this mechanism. The new elements "sign in" to the nodes
439  by increasing the reference counts of their vertices and edges. The
440  parent element "signs off" its nodes by decrementing the ref counts. */
441  struct Node : public Hashed2
442  {
445 
446  Node() : vert_refc(0), edge_refc(0), vert_index(-1), edge_index(-1) {}
447  ~Node();
448 
449  bool HasVertex() const { return vert_refc > 0; }
450  bool HasEdge() const { return edge_refc > 0; }
451 
452  // decrease vertex/edge ref count, return false if Node should be deleted
453  bool UnrefVertex() { --vert_refc; return vert_refc || edge_refc; }
454  bool UnrefEdge() { --edge_refc; return vert_refc || edge_refc; }
455  };
456 
457  /** Similarly to nodes, faces can be accessed by hashing their four vertex
458  node IDs. A face knows about the one or two elements that are using it.
459  A face that is not on the boundary and only has one element referencing
460  it is either a master or a slave face. */
461  struct Face : public Hashed4
462  {
463  int attribute; ///< boundary element attribute, -1 if internal face
464  int index; ///< face number in the Mesh
465  int elem[2]; ///< up to 2 elements sharing the face
466 
467  Face() : attribute(-1), index(-1) { elem[0] = elem[1] = -1; }
468 
469  bool Boundary() const { return attribute >= 0; }
470  bool Unused() const { return elem[0] < 0 && elem[1] < 0; }
471 
472  // add or remove an element from the 'elem[2]' array
473  void RegisterElement(int e);
474  void ForgetElement(int e);
475 
476  /// Return one of elem[0] or elem[1] and make sure the other is -1.
477  int GetSingleElement() const;
478  };
479 
480  /** This is an element in the refinement hierarchy. Each element has
481  either been refined and points to its children, or is a leaf and points
482  to its vertex nodes. */
483  struct Element
484  {
485  char geom; ///< Geometry::Type of the element (char for storage only)
486  char ref_type; ///< bit mask of X,Y,Z refinements (bits 0,1,2 respectively)
487  char tet_type; ///< tetrahedron split type, currently always 0
488  char flag; ///< generic flag/marker, can be used by algorithms
489  int index; ///< element number in the Mesh, -1 if refined
490  int rank; ///< processor number (ParNCMesh), -1 if undefined/unknown
492  union
493  {
494  int node[MaxElemNodes]; ///< element corners (if ref_type == 0)
495  int child[MaxElemChildren]; ///< 2-10 children (if ref_type != 0)
496  };
497  int parent; ///< parent element, -1 if this is a root element, -2 if free'd
498 
499  Element(Geometry::Type geom, int attr);
500 
501  Geometry::Type Geom() const { return Geometry::Type(geom); }
502  bool IsLeaf() const { return !ref_type && (parent != -2); }
503  };
504 
505 
506  // primary data
507 
508  HashTable<Node> nodes; // associative container holding all Nodes
509  HashTable<Face> faces; // associative container holding all Faces
510 
511  BlockArray<Element> elements; // storage for all Elements
512  Array<int> free_element_ids; // unused element ids - indices into 'elements'
513 
514  /** Initial traversal state (~ element orientation) for each root element
515  NOTE: M = root_state.Size() is the number of root elements.
516  NOTE: the first M items of 'elements' is the coarse mesh. */
518 
519  /** Coordinates of top-level vertices (organized as triples). If empty,
520  the Mesh is curved (Nodes != NULL) and NCMesh is topology-only. */
522 
523 
524  // secondary data
525 
526  /** Apart from the primary data structure, which is the element/node/face
527  hierarchy, there is secondary data that is derived from the primary
528  data and needs to be updated when the primary data changes. Update()
529  takes care of that and needs to be called after each refinement and
530  derefinement. */
531  virtual void Update();
532 
533  // set by UpdateLeafElements, UpdateVertices and OnMeshUpdated
535 
536  // NOTE: the serial code understands the bare minimum about ghost elements and
537  // other ghost entities in order to be able to load parallel partial meshes
539 
540  Array<int> leaf_elements; ///< finest elements, in Mesh ordering (+ ghosts)
541  Array<int> leaf_sfc_index; ///< natural tree ordering of leaf elements
542  Array<int> vertex_nodeId; ///< vertex-index to node-id map, see UpdateVertices
543 
544  NCList face_list; ///< lazy-initialized list of faces, see GetFaceList
545  NCList edge_list; ///< lazy-initialized list of edges, see GetEdgeList
546  NCList vertex_list; ///< lazy-initialized list of vertices, see GetVertexList
547 
548  Array<int> boundary_faces; ///< subset of all faces, set by BuildFaceList
549  Array<char> face_geom; ///< face geometry by face index, set by OnMeshUpdated
550 
551  Table element_vertex; ///< leaf-element to vertex table, see FindSetNeighbors
552 
553 
554  /// Update the leaf elements indices in leaf_elements
555  void UpdateLeafElements();
556 
557  /** @brief This method assigns indices to vertices (Node::vert_index) that
558  will be seen by the Mesh class and the rest of MFEM.
559 
560  We must be careful to:
561  1. Stay compatible with the conforming code, which expects top-level
562  (original) vertices to be indexed first, otherwise GridFunctions
563  defined on a conforming mesh would no longer be valid when the
564  mesh is converted to an NC mesh.
565 
566  2. Make sure serial NCMesh is compatible with the parallel ParNCMesh,
567  so it is possible to read parallel partial solutions in serial code
568  (e.g., serial GLVis). This means handling ghost elements, if present.
569 
570  3. Assign vertices in a globally consistent order for parallel meshes:
571  if two vertices i,j are shared by two ranks r1,r2, and i<j on r1,
572  then i<j on r2 as well. This is true for top-level vertices but also
573  for the remaining shared vertices thanks to the globally consistent
574  SFC ordering of the leaf elements. This property reduces communication
575  and simplifies ParNCMesh. */
576  void UpdateVertices(); ///< update Vertex::index and vertex_nodeId
577 
578  /** Collect the leaf elements in leaf_elements, and the ghost elements in
579  ghosts. Compute and set the element indices of @a elements. On quad and
580  hex refined elements tries to order leaf elements along a space-filling
581  curve according to the given @a state variable. */
582  void CollectLeafElements(int elem, int state, Array<int> &ghosts,
583  int &counter);
584 
585  /** Try to find a space-filling curve friendly orientation of the root
586  elements: set 'root_state' based on the ordering of coarse elements.
587  Note that the coarse mesh itself must be ordered as an SFC by e.g.
588  Mesh::GetGeckoElementOrdering. */
589  void InitRootState(int root_count);
590 
591  /** Compute the Geometry::Type present in the root elements (coarse elements)
592  and set @a Geoms bitmask accordingly. */
593  void InitGeomFlags();
594 
595  /// Return true if the mesh contains prism elements.
596  bool HavePrisms() const { return Geoms & (1 << Geometry::PRISM); }
597 
598  /// Return true if the mesh contains pyramid elements.
599  bool HavePyramids() const { return Geoms & (1 << Geometry::PYRAMID); }
600 
601  /// Return true if the mesh contains tetrahedral elements.
602  bool HaveTets() const { return Geoms & (1 << Geometry::TETRAHEDRON); }
603 
604  /// Return true if the Element @a el is a ghost element.
605  bool IsGhost(const Element &el) const { return el.rank != MyRank; }
606 
607 
608  // refinement/derefinement
609 
610  Array<Refinement> ref_stack; ///< stack of scheduled refinements (temporary)
611  HashTable<Node> shadow; ///< temporary storage for reparented nodes
612  Array<Triple<int, int, int> > reparents; ///< scheduled node reparents (tmp)
613 
614  Table derefinements; ///< possible derefinements, see GetDerefinementTable
615 
616  /** Refine the element @a elem with the refinement @a ref_type
617  (c.f. Refinement::enum) */
618  void RefineElement(int elem, char ref_type);
619 
620  /// Derefine the element @a elem, does nothing on leaf elements.
621  void DerefineElement(int elem);
622 
623  // Add an Element @a el to the NCMesh, optimized to reuse freed elements.
624  int AddElement(const Element &el)
625  {
626  if (free_element_ids.Size())
627  {
628  int idx = free_element_ids.Last();
630  elements[idx] = el;
631  return idx;
632  }
633  return elements.Append(el);
634  }
635 
636  // Free the element with index @a id.
637  void FreeElement(int id)
638  {
640  elements[id].ref_type = 0;
641  elements[id].parent = -2; // mark the element as free
642  }
643 
644  int NewHexahedron(int n0, int n1, int n2, int n3,
645  int n4, int n5, int n6, int n7, int attr,
646  int fattr0, int fattr1, int fattr2,
647  int fattr3, int fattr4, int fattr5);
648 
649  int NewWedge(int n0, int n1, int n2,
650  int n3, int n4, int n5, int attr,
651  int fattr0, int fattr1,
652  int fattr2, int fattr3, int fattr4);
653 
654  int NewTetrahedron(int n0, int n1, int n2, int n3, int attr,
655  int fattr0, int fattr1, int fattr2, int fattr3);
656 
657  int NewPyramid(int n0, int n1, int n2, int n3, int n4, int attr,
658  int fattr0, int fattr1, int fattr2, int fattr3,
659  int fattr4);
660 
661  int NewQuadrilateral(int n0, int n1, int n2, int n3, int attr,
662  int eattr0, int eattr1, int eattr2, int eattr3);
663 
664  int NewTriangle(int n0, int n1, int n2,
665  int attr, int eattr0, int eattr1, int eattr2);
666 
667  int NewSegment(int n0, int n1, int attr, int vattr1, int vattr2);
668 
669  mfem::Element* NewMeshElement(int geom) const;
670 
671  int QuadFaceSplitType(int v1, int v2, int v3, int v4, int mid[5]
672  = NULL /*optional output of mid-edge nodes*/) const;
673 
674  bool TriFaceSplit(int v1, int v2, int v3, int mid[3] = NULL) const;
675 
676  void ForceRefinement(int vn1, int vn2, int vn3, int vn4);
677 
678  void FindEdgeElements(int vn1, int vn2, int vn3, int vn4,
679  Array<MeshId> &prisms) const;
680 
681  void CheckAnisoPrism(int vn1, int vn2, int vn3, int vn4,
682  const Refinement *refs, int nref);
683 
684  void CheckAnisoFace(int vn1, int vn2, int vn3, int vn4,
685  int mid12, int mid34, int level = 0);
686 
687  void CheckIsoFace(int vn1, int vn2, int vn3, int vn4,
688  int en1, int en2, int en3, int en4, int midf);
689 
690  void ReparentNode(int node, int new_p1, int new_p2);
691 
692  int FindMidEdgeNode(int node1, int node2) const;
693  int GetMidEdgeNode(int node1, int node2);
694 
695  int GetMidFaceNode(int en1, int en2, int en3, int en4);
696 
697  void ReferenceElement(int elem);
698  void UnreferenceElement(int elem, Array<int> &elemFaces);
699 
700  Face* GetFace(Element &elem, int face_no);
701  void RegisterFaces(int elem, int *fattr = NULL);
702  void DeleteUnusedFaces(const Array<int> &elemFaces);
703 
704  void CollectDerefinements(int elem, Array<Connection> &list);
705 
706  /// Return el.node[index] correctly, even if the element is refined.
707  int RetrieveNode(const Element &el, int index);
708 
709  /// Extended version of find_node: works if 'el' is refined.
710  int FindNodeExt(const Element &el, int node, bool abort = true);
711 
712 
713  // face/edge lists
714 
715  static int find_node(const Element &el, int node);
716  static int find_element_edge(const Element &el, int vn0, int vn1,
717  bool abort = true);
718  static int find_local_face(int geom, int a, int b, int c);
719 
720  struct Point;
721  struct PointMatrix;
722 
723  int ReorderFacePointMat(int v0, int v1, int v2, int v3,
724  int elem, const PointMatrix &pm,
725  PointMatrix &reordered) const;
726 
727  void TraverseQuadFace(int vn0, int vn1, int vn2, int vn3,
728  const PointMatrix& pm, int level, Face* eface[4],
729  MatrixMap &matrix_map);
730  bool TraverseTriFace(int vn0, int vn1, int vn2,
731  const PointMatrix& pm, int level,
732  MatrixMap &matrix_map);
733  void TraverseTetEdge(int vn0, int vn1, const Point &p0, const Point &p1,
734  MatrixMap &matrix_map);
735  void TraverseEdge(int vn0, int vn1, double t0, double t1, int flags,
736  int level, MatrixMap &matrix_map);
737 
738  virtual void BuildFaceList();
739  virtual void BuildEdgeList();
740  virtual void BuildVertexList();
741 
742  virtual void ElementSharesFace(int elem, int local, int face) {} // ParNCMesh
743  virtual void ElementSharesEdge(int elem, int local, int enode) {} // ParNCMesh
744  virtual void ElementSharesVertex(int elem, int local, int vnode) {} // ParNCMesh
745 
746 
747  // neighbors / element_vertex table
748 
749  /** Return all vertex-, edge- and face-neighbors of a set of elements.
750  The neighbors are returned as a list (neighbors != NULL), as a set
751  (neighbor_set != NULL), or both. The sizes of the set arrays must match
752  that of leaf_elements. The function is intended to be used for large
753  sets of elements and its complexity is linear in the number of leaf
754  elements in the mesh. */
755  void FindSetNeighbors(const Array<char> &elem_set,
756  Array<int> *neighbors, /* append */
757  Array<char> *neighbor_set = NULL);
758 
759  /** Return all vertex-, edge- and face-neighbors of a single element.
760  You can limit the number of elements being checked using 'search_set'.
761  The complexity of the function is linear in the size of the search set.*/
762  void FindNeighbors(int elem,
763  Array<int> &neighbors, /* append */
764  const Array<int> *search_set = NULL);
765 
766  /** Expand a set of elements by all vertex-, edge- and face-neighbors.
767  The output array 'expanded' will contain all items from 'elems'
768  (provided they are in 'search_set') plus their neighbors. The neighbor
769  search can be limited to the optional search set. The complexity is
770  linear in the sum of the sizes of 'elems' and 'search_set'. */
771  void NeighborExpand(const Array<int> &elems,
772  Array<int> &expanded,
773  const Array<int> *search_set = NULL);
774 
775 
776  void CollectEdgeVertices(int v0, int v1, Array<int> &indices);
777  void CollectTriFaceVertices(int v0, int v1, int v2, Array<int> &indices);
778  void CollectQuadFaceVertices(int v0, int v1, int v2, int v3,
779  Array<int> &indices);
781 
783  {
785  }
786 
787  int GetVertexRootCoord(int elem, RefCoord coord[3]) const;
788  void CollectIncidentElements(int elem, const RefCoord coord[3],
789  Array<int> &list) const;
790 
791  /** Return elements neighboring to a local vertex of element 'elem'. Only
792  elements from within the same refinement tree ('cousins') are returned.
793  Complexity is proportional to the depth of elem's refinement tree. */
794  void FindVertexCousins(int elem, int local, Array<int> &cousins) const;
795 
796 
797  // coarse/fine transformations
798 
799  struct Point
800  {
801  int dim;
802  double coord[3];
803 
804  Point() { dim = 0; }
805 
806  Point(const Point &) = default;
807 
808  Point(double x)
809  { dim = 1; coord[0] = x; }
810 
811  Point(double x, double y)
812  { dim = 2; coord[0] = x; coord[1] = y; }
813 
814  Point(double x, double y, double z)
815  { dim = 3; coord[0] = x; coord[1] = y; coord[2] = z; }
816 
817  Point(const Point& p0, const Point& p1)
818  {
819  dim = p0.dim;
820  for (int i = 0; i < dim; i++)
821  {
822  coord[i] = (p0.coord[i] + p1.coord[i]) * 0.5;
823  }
824  }
825 
826  Point(const Point& p0, const Point& p1, const Point& p2, const Point& p3)
827  {
828  dim = p0.dim;
829  MFEM_ASSERT(p1.dim == dim && p2.dim == dim && p3.dim == dim, "");
830  for (int i = 0; i < dim; i++)
831  {
832  coord[i] = (p0.coord[i] + p1.coord[i] + p2.coord[i] + p3.coord[i])
833  * 0.25;
834  }
835  }
836 
837  Point& operator=(const Point& src)
838  {
839  dim = src.dim;
840  for (int i = 0; i < dim; i++) { coord[i] = src.coord[i]; }
841  return *this;
842  }
843  };
844 
845  /** @brief The PointMatrix stores the coordinates of the slave face using the
846  master face coordinate as reference.
847 
848  In 2D, the point matrix has the orientation of the parent
849  edge, so its columns need to be flipped when applying it, see
850  ApplyLocalSlaveTransformation.
851 
852  In 3D, the orientation part of Elem2Inf is encoded in the point
853  matrix.
854 
855  The following transformation gives the relation between the
856  reference quad face coordinates (xi, eta) in [0,1]^2, and the fine quad
857  face coordinates (x, y):
858  x = a0*(1-xi)*(1-eta) + a1*xi*(1-eta) + a2*xi*eta + a3*(1-xi)*eta
859  y = b0*(1-xi)*(1-eta) + b1*xi*(1-eta) + b2*xi*eta + b3*(1-xi)*eta
860  */
861  struct PointMatrix
862  {
863  int np;
865 
866  PointMatrix() : np(0) {}
867 
868  PointMatrix(const Point& p0, const Point& p1)
869  { np = 2; points[0] = p0; points[1] = p1; }
870 
871  PointMatrix(const Point& p0, const Point& p1, const Point& p2)
872  { np = 3; points[0] = p0; points[1] = p1; points[2] = p2; }
873 
874  PointMatrix(const Point& p0, const Point& p1, const Point& p2, const Point& p3)
875  { np = 4; points[0] = p0; points[1] = p1; points[2] = p2; points[3] = p3; }
876 
877  PointMatrix(const Point& p0, const Point& p1, const Point& p2,
878  const Point& p3, const Point& p4)
879  {
880  np = 5;
881  points[0] = p0; points[1] = p1; points[2] = p2;
882  points[3] = p3; points[4] = p4;
883  }
884  PointMatrix(const Point& p0, const Point& p1, const Point& p2,
885  const Point& p3, const Point& p4, const Point& p5)
886  {
887  np = 6;
888  points[0] = p0; points[1] = p1; points[2] = p2;
889  points[3] = p3; points[4] = p4; points[5] = p5;
890  }
891  PointMatrix(const Point& p0, const Point& p1, const Point& p2,
892  const Point& p3, const Point& p4, const Point& p5,
893  const Point& p6, const Point& p7)
894  {
895  np = 8;
896  points[0] = p0; points[1] = p1; points[2] = p2; points[3] = p3;
897  points[4] = p4; points[5] = p5; points[6] = p6; points[7] = p7;
898  }
899 
900  Point& operator()(int i) { return points[i]; }
901  const Point& operator()(int i) const { return points[i]; }
902 
903  bool operator==(const PointMatrix &pm) const;
904 
905  void GetMatrix(DenseMatrix& point_matrix) const;
906  };
907 
915 
916  static const PointMatrix& GetGeomIdentity(Geometry::Type geom);
917 
918  void GetPointMatrix(Geometry::Type geom, const char* ref_path,
919  DenseMatrix& matrix);
920 
921  typedef std::map<std::string, int> RefPathMap;
922 
923  void TraverseRefinements(int elem, int coarse_index,
924  std::string &ref_path, RefPathMap &map);
925 
926  /// storage for data returned by Get[De]RefinementTransforms()
928 
929  /// state of leaf_elements before Refine(), set by MarkCoarseLevel()
931 
932  void InitDerefTransforms();
933  void SetDerefMatrixCodes(int parent, Array<int> &fine_coarse);
934 
935 
936  // vertex temporary data, used by GetMeshComponents
937 
938  struct TmpVertex
939  {
940  bool valid, visited;
941  double pos[3];
942  TmpVertex() : valid(false), visited(false) {}
943  };
944 
946 
947  const double *CalcVertexPos(int node) const;
948 
949 
950  // utility
951 
952  int GetEdgeMaster(int node) const;
953 
954  void FindFaceNodes(int face, int node[4]);
955 
956  int EdgeSplitLevel(int vn1, int vn2) const;
957  int TriFaceSplitLevel(int vn1, int vn2, int vn3) const;
958  void QuadFaceSplitLevel(int vn1, int vn2, int vn3, int vn4,
959  int& h_level, int& v_level) const;
960 
961  void CountSplits(int elem, int splits[3]) const;
962  void GetLimitRefinements(Array<Refinement> &refinements, int max_level);
963 
964 
965  // I/O
966 
967  /// Print the "vertex_parents" section of the mesh file.
968  int PrintVertexParents(std::ostream *out) const;
969  /// Load the vertex parent hierarchy from a mesh file.
970  void LoadVertexParents(std::istream &input);
971 
972  /** Print the "boundary" section of the mesh file.
973  If out == NULL, only return the number of boundary elements. */
974  int PrintBoundary(std::ostream *out) const;
975  /// Load the "boundary" section of the mesh file.
976  void LoadBoundary(std::istream &input);
977 
978  /// Print the "coordinates" section of the mesh file.
979  void PrintCoordinates(std::ostream &out) const;
980  /// Load the "coordinates" section of the mesh file.
981  void LoadCoordinates(std::istream &input);
982 
983  /// Count root elements and initialize root_state.
984  void InitRootElements();
985  /// Return the index of the last top-level node plus one.
986  int CountTopLevelNodes() const;
987  /// Return true if all root_states are zero.
988  bool ZeroRootStates() const;
989 
990  /// Load the element refinement hierarchy from a legacy mesh file.
991  void LoadCoarseElements(std::istream &input);
992  void CopyElements(int elem, const BlockArray<Element> &tmp_elements);
993  /// Load the deprecated MFEM mesh v1.1 format for backward compatibility.
994  void LoadLegacyFormat(std::istream &input, int &curved, int &is_nc);
995 
996 
997  // geometry
998 
999  /// This holds in one place the constants about the geometries we support
1000  struct GeomInfo
1001  {
1002  int nv, ne, nf; // number of: vertices, edges, faces
1003  int edges[MaxElemEdges][2]; // edge vertices (up to 12 edges)
1004  int faces[MaxElemFaces][4]; // face vertices (up to 6 faces)
1005  int nfv[MaxElemFaces]; // number of face vertices
1006 
1008  GeomInfo() : initialized(false) {}
1009  void InitGeom(Geometry::Type geom);
1010  };
1011 
1013 
1014 #ifdef MFEM_DEBUG
1015 public:
1016  void DebugLeafOrder(std::ostream &out) const;
1017  void DebugDump(std::ostream &out) const;
1018 #endif
1019 
1020  friend class ParNCMesh; // for ParNCMesh::ElementSet
1021  friend struct MatrixMap;
1022  friend struct PointMatrixHash;
1023 };
1024 
1025 }
1026 
1027 #endif
NCList face_list
lazy-initialized list of faces, see GetFaceList
Definition: ncmesh.hpp:544
NCList edge_list
lazy-initialized list of edges, see GetEdgeList
Definition: ncmesh.hpp:545
The PointMatrix stores the coordinates of the slave face using the master face coordinate as referenc...
Definition: ncmesh.hpp:861
void DebugDump(std::ostream &out) const
Definition: ncmesh.cpp:6329
bool TriFaceSplit(int v1, int v2, int v3, int mid[3]=NULL) const
Definition: ncmesh.cpp:2731
int PrintBoundary(std::ostream *out) const
Definition: ncmesh.cpp:5507
void LoadLegacyFormat(std::istream &input, int &curved, int &is_nc)
Load the deprecated MFEM mesh v1.1 format for backward compatibility.
Definition: ncmesh.cpp:6027
static void GridSfcOrdering3D(int width, int height, int depth, Array< int > &coords)
Definition: ncmesh.cpp:4975
int faces[MaxElemFaces][4]
Definition: ncmesh.hpp:1004
virtual int GetNGhostElements() const
Definition: ncmesh.hpp:154
bool IsLegacyLoaded() const
I/O: Return true if the mesh was loaded from the legacy v1.1 format.
Definition: ncmesh.hpp:384
void CheckAnisoPrism(int vn1, int vn2, int vn3, int vn4, const Refinement *refs, int nref)
Definition: ncmesh.cpp:787
void Print(std::ostream &out) const
I/O: Print the mesh in "MFEM NC mesh v1.0" format.
Definition: ncmesh.cpp:5634
CoarseFineTransformations transforms
storage for data returned by Get[De]RefinementTransforms()
Definition: ncmesh.hpp:927
void DebugLeafOrder(std::ostream &out) const
Definition: ncmesh.cpp:6304
static void GridSfcOrdering2D(int width, int height, Array< int > &coords)
Definition: ncmesh.cpp:4960
int NewQuadrilateral(int n0, int n1, int n2, int n3, int attr, int eattr0, int eattr1, int eattr2, int eattr3)
Definition: ncmesh.cpp:590
static PointMatrix pm_seg_identity
Definition: ncmesh.hpp:908
void QuadFaceSplitLevel(int vn1, int vn2, int vn3, int vn4, int &h_level, int &v_level) const
Definition: ncmesh.cpp:5285
int elem[2]
up to 2 elements sharing the face
Definition: ncmesh.hpp:465
char ref_type
refinement XYZ bit mask (7 = full isotropic)
Definition: ncmesh.hpp:39
bool ZeroRootStates() const
Return true if all root_states are zero.
Definition: ncmesh.cpp:5625
TmpVertex * tmp_vertex
Definition: ncmesh.hpp:945
int GetEdgeMaster(int v1, int v2) const
Definition: ncmesh.cpp:5108
void SetDerefMatrixCodes(int parent, Array< int > &fine_coarse)
Definition: ncmesh.cpp:2036
const CoarseFineTransformations & GetDerefinementTransforms()
Definition: ncmesh.cpp:4623
int child[MaxElemChildren]
2-10 children (if ref_type != 0)
Definition: ncmesh.hpp:495
char flag
generic flag/marker, can be used by algorithms
Definition: ncmesh.hpp:488
void FindVertexCousins(int elem, int local, Array< int > &cousins) const
Definition: ncmesh.cpp:3957
void FindEdgeElements(int vn1, int vn2, int vn3, int vn4, Array< MeshId > &prisms) const
Definition: ncmesh.cpp:740
Refinement(int index, int type=Refinement::XYZ)
Definition: ncmesh.hpp:43
void OrientedPointMatrix(const Slave &slave, DenseMatrix &oriented_matrix) const
Return the point matrix oriented according to the master and slave edges.
Definition: ncmesh.cpp:3368
PointMatrix(const Point &p0, const Point &p1, const Point &p2)
Definition: ncmesh.hpp:871
This holds in one place the constants about the geometries we support.
Definition: ncmesh.hpp:1000
int GetVertexRootCoord(int elem, RefCoord coord[3]) const
Definition: ncmesh.cpp:3881
Array< Triple< int, int, int > > reparents
scheduled node reparents (tmp)
Definition: ncmesh.hpp:612
bool HavePyramids() const
Return true if the mesh contains pyramid elements.
Definition: ncmesh.hpp:599
Point(const Point &p0, const Point &p1)
Definition: ncmesh.hpp:817
static const int NumGeom
Definition: geom.hpp:42
Array< Slave > slaves
Definition: ncmesh.hpp:231
void CollectTriFaceVertices(int v0, int v1, int v2, Array< int > &indices)
Definition: ncmesh.cpp:3490
char tet_type
tetrahedron split type, currently always 0
Definition: ncmesh.hpp:487
int NewSegment(int n0, int n1, int attr, int vattr1, int vattr2)
Definition: ncmesh.cpp:642
virtual void LimitNCLevel(int max_nc_level)
Definition: ncmesh.cpp:5437
virtual void Trim()
Save memory by releasing all non-essential and cached data.
Definition: ncmesh.cpp:6210
int index
element number in the Mesh, -1 if refined
Definition: ncmesh.hpp:489
PointMatrix(const Point &p0, const Point &p1, const Point &p2, const Point &p3)
Definition: ncmesh.hpp:874
std::int64_t RefCoord
Definition: ncmesh.hpp:397
void LegacyToNewVertexOrdering(Array< int > &order) const
I/O: Return a map from old (v1.1) vertex indices to new vertex indices.
Definition: ncmesh.cpp:6185
Embedding()=default
int Dimension() const
Return the dimension of the NCMesh.
Definition: ncmesh.hpp:144
void ForceRefinement(int vn1, int vn2, int vn3, int vn4)
Definition: ncmesh.cpp:682
PointMatrix(const Point &p0, const Point &p1, const Point &p2, const Point &p3, const Point &p4, const Point &p5, const Point &p6, const Point &p7)
Definition: ncmesh.hpp:891
Lists all edges/faces in the nonconforming mesh.
Definition: ncmesh.hpp:227
long MemoryUsage() const
Return total number of bytes allocated.
Definition: ncmesh.cpp:6253
bool Legacy
true if the mesh was loaded from the legacy v1.1 format
Definition: ncmesh.hpp:422
unsigned matrix
index into NCList::point_matrices[geom]
Definition: ncmesh.hpp:217
int NVertices
Definition: ncmesh.hpp:534
void ClearTransforms()
Free all internal data created by the above three functions.
Definition: ncmesh.cpp:4703
int ReorderFacePointMat(int v0, int v1, int v2, int v3, int elem, const PointMatrix &pm, PointMatrix &reordered) const
Definition: ncmesh.cpp:2864
Table derefinements
possible derefinements, see GetDerefinementTable
Definition: ncmesh.hpp:614
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
NCList vertex_list
lazy-initialized list of vertices, see GetVertexList
Definition: ncmesh.hpp:546
void InitDerefTransforms()
Definition: ncmesh.cpp:2017
static PointMatrix pm_prism_identity
Definition: ncmesh.hpp:912
char geom
Geometry::Type of the element (char for storage only)
Definition: ncmesh.hpp:485
void InitRootState(int root_count)
Definition: ncmesh.cpp:2331
int NewTetrahedron(int n0, int n1, int n2, int n3, int attr, int fattr0, int fattr1, int fattr2, int fattr3)
Definition: ncmesh.cpp:535
PointMatrix(const Point &p0, const Point &p1, const Point &p2, const Point &p3, const Point &p4)
Definition: ncmesh.hpp:877
void CheckAnisoFace(int vn1, int vn2, int vn3, int vn4, int mid12, int mid34, int level=0)
Definition: ncmesh.cpp:812
static PointMatrix pm_tri_identity
Definition: ncmesh.hpp:909
bool HasEdge() const
Definition: ncmesh.hpp:450
Geometry::Type Geom() const
Definition: ncmesh.hpp:194
void FreeElement(int id)
Definition: ncmesh.hpp:637
bool HasVertex() const
Definition: ncmesh.hpp:449
bool IsGhost(const Element &el) const
Return true if the Element el is a ghost element.
Definition: ncmesh.hpp:605
void CollectDerefinements(int elem, Array< Connection > &list)
Definition: ncmesh.cpp:1904
Point & operator()(int i)
Definition: ncmesh.hpp:900
void GetLimitRefinements(Array< Refinement > &refinements, int max_level)
Definition: ncmesh.cpp:5407
void GetMatrix(DenseMatrix &point_matrix) const
Definition: ncmesh.cpp:3988
void TraverseEdge(int vn0, int vn1, double t0, double t1, int flags, int level, MatrixMap &matrix_map)
Definition: ncmesh.cpp:3202
static const int MaxElemEdges
Number of edges of an element can have.
Definition: ncmesh.hpp:426
double coord[3]
Definition: ncmesh.hpp:802
bool UnrefVertex()
Definition: ncmesh.hpp:453
int GetNEdges() const
Return the number of edges in the NCMesh.
Definition: ncmesh.hpp:151
int spaceDim
dimensions of the elements and the vertex coordinates
Definition: ncmesh.hpp:418
Array< int > vertex_nodeId
vertex-index to node-id map, see UpdateVertices
Definition: ncmesh.hpp:542
int attribute
boundary element attribute, -1 if internal face
Definition: ncmesh.hpp:463
void FindFaceNodes(int face, int node[4])
Definition: ncmesh.cpp:5164
Array< int > root_state
Definition: ncmesh.hpp:517
static const PointMatrix & GetGeomIdentity(Geometry::Type geom)
Definition: ncmesh.cpp:4025
bool HavePrisms() const
Return true if the mesh contains prism elements.
Definition: ncmesh.hpp:596
void DeleteAll()
Delete the whole array.
Definition: array.hpp:854
static int find_element_edge(const Element &el, int vn0, int vn1, bool abort=true)
Definition: ncmesh.cpp:2768
int index
Mesh element number.
Definition: ncmesh.hpp:38
int master
master number (in Mesh numbering)
Definition: ncmesh.hpp:216
const NCList & GetFaceList()
Return the current list of conforming and nonconforming faces.
Definition: ncmesh.hpp:253
void OnMeshUpdated(Mesh *mesh)
Definition: ncmesh.cpp:2540
A parallel extension of the NCMesh class.
Definition: pncmesh.hpp:64
int NewWedge(int n0, int n1, int n2, int n3, int n4, int n5, int attr, int fattr0, int fattr1, int fattr2, int fattr3, int fattr4)
Definition: ncmesh.cpp:503
void CollectEdgeVertices(int v0, int v1, Array< int > &indices)
Definition: ncmesh.cpp:3478
void CollectLeafElements(int elem, int state, Array< int > &ghosts, int &counter)
Definition: ncmesh.cpp:2055
static int find_local_face(int geom, int a, int b, int c)
Definition: ncmesh.cpp:2786
Array< Refinement > ref_stack
stack of scheduled refinements (temporary)
Definition: ncmesh.hpp:610
Point(const Point &p0, const Point &p1, const Point &p2, const Point &p3)
Definition: ncmesh.hpp:826
void GetMeshComponents(Mesh &mesh) const
Fill Mesh::{vertices,elements,boundary} for the current finest level.
Definition: ncmesh.cpp:2435
void ReferenceElement(int elem)
Definition: ncmesh.cpp:322
int GetElementDepth(int i) const
Return the distance of leaf &#39;i&#39; from the root.
Definition: ncmesh.cpp:5117
void GetPointMatrix(Geometry::Type geom, const char *ref_path, DenseMatrix &matrix)
Definition: ncmesh.cpp:4042
void CheckIsoFace(int vn1, int vn2, int vn3, int vn4, int en1, int en2, int en3, int en4, int midf)
Definition: ncmesh.cpp:907
void InitGeomFlags()
Definition: ncmesh.cpp:229
Face * GetFace(Element &elem, int face_no)
Definition: ncmesh.cpp:435
long MemoryUsage() const
Definition: ncmesh.cpp:6225
void TraverseRefinements(int elem, int coarse_index, std::string &ref_path, RefPathMap &map)
Definition: ncmesh.cpp:4538
int GetMidEdgeNode(int node1, int node2)
Definition: ncmesh.cpp:307
Element(Geometry::Type geom, int attr)
Definition: ncmesh.cpp:460
Geometry::Type Geom() const
Definition: ncmesh.hpp:501
int index
face number in the Mesh
Definition: ncmesh.hpp:464
const Table & GetDerefinementTable()
Definition: ncmesh.cpp:1937
Table element_vertex
leaf-element to vertex table, see FindSetNeighbors
Definition: ncmesh.hpp:551
int GetNFaces() const
Return the number of (2D) faces in the NCMesh.
Definition: ncmesh.hpp:153
DenseTensor point_matrices[Geometry::NumGeom]
Definition: ncmesh.hpp:77
unsigned matrix
Definition: ncmesh.hpp:58
int Append(const T &el)
Append element &#39;el&#39; to array, resize if necessary.
Definition: array.hpp:759
const double * CalcVertexPos(int node) const
Definition: ncmesh.cpp:2410
int AddElement(const Element &el)
Definition: ncmesh.hpp:624
int SpaceDimension() const
Return the space dimension of the NCMesh.
Definition: ncmesh.hpp:146
void GetEdgeVertices(const MeshId &edge_id, int vert_index[2], bool oriented=true) const
Return Mesh vertex indices of an edge identified by &#39;edge_id&#39;.
Definition: ncmesh.cpp:5007
virtual void Derefine(const Array< int > &derefs)
Definition: ncmesh.cpp:1981
double b
Definition: lissajous.cpp:42
void FindNeighbors(int elem, Array< int > &neighbors, const Array< int > *search_set=NULL)
Definition: ncmesh.cpp:3739
void InitGeom(Geometry::Type geom)
Definition: ncmesh.cpp:27
friend struct MatrixMap
Definition: ncmesh.hpp:1021
Array< int > coarse_elements
state of leaf_elements before Refine(), set by MarkCoarseLevel()
Definition: ncmesh.hpp:930
void UpdateVertices()
This method assigns indices to vertices (Node::vert_index) that will be seen by the Mesh class and th...
Definition: ncmesh.cpp:2151
int FindMidEdgeNode(int node1, int node2) const
Definition: ncmesh.cpp:291
int NGhostVertices
Definition: ncmesh.hpp:538
Master(int index, int element, int local, int geom, int sb, int se)
Definition: ncmesh.hpp:208
bool Empty() const
Definition: ncmesh.hpp:241
Array< DenseMatrix * > point_matrices[Geometry::NumGeom]
List of unique point matrices for each slave geometry.
Definition: ncmesh.hpp:234
void RegisterFaces(int elem, int *fattr=NULL)
Definition: ncmesh.cpp:396
static PointMatrix pm_pyramid_identity
Definition: ncmesh.hpp:913
Array< Embedding > embeddings
Fine element positions in their parents.
Definition: ncmesh.hpp:73
virtual void BuildEdgeList()
Definition: ncmesh.cpp:3230
int GetElementSizeReduction(int i) const
Definition: ncmesh.cpp:5129
const CoarseFineTransformations & GetRefinementTransforms()
Definition: ncmesh.cpp:4572
void LoadCoarseElements(std::istream &input)
Load the element refinement hierarchy from a legacy mesh file.
Definition: ncmesh.cpp:5956
Identifies a vertex/edge/face in both Mesh and NCMesh.
Definition: ncmesh.hpp:187
int parent
parent element, -1 if this is a root element, -2 if free&#39;d
Definition: ncmesh.hpp:497
void MarkCoarseLevel()
Definition: ncmesh.cpp:4524
signed char local
local number within &#39;element&#39;
Definition: ncmesh.hpp:191
virtual void ElementSharesFace(int elem, int local, int face)
Definition: ncmesh.hpp:742
int GetNumRootElements()
Return the number of root elements.
Definition: ncmesh.hpp:366
bool HaveTets() const
Return true if the mesh contains tetrahedral elements.
Definition: ncmesh.hpp:602
void LoadCoordinates(std::istream &input)
Load the "coordinates" section of the mesh file.
Definition: ncmesh.cpp:5604
virtual void ElementSharesVertex(int elem, int local, int vnode)
Definition: ncmesh.hpp:744
int MyRank
used in parallel, or when loading a parallel file in serial
Definition: ncmesh.hpp:419
A class for non-conforming AMR. The class is not used directly by the user, rather it is an extension...
Definition: ncmesh.hpp:121
int EdgeSplitLevel(int vn1, int vn2) const
Definition: ncmesh.cpp:5261
unsigned geom
Definition: ncmesh.hpp:57
virtual ~NCMesh()
Definition: ncmesh.cpp:250
int CountTopLevelNodes() const
Return the index of the last top-level node plus one.
Definition: ncmesh.cpp:5754
friend struct PointMatrixHash
Definition: ncmesh.hpp:1022
virtual void GetBoundaryClosure(const Array< int > &bdr_attr_is_ess, Array< int > &bdr_vertices, Array< int > &bdr_edges)
Definition: ncmesh.cpp:5188
int element
NCMesh::Element containing this vertex/edge/face.
Definition: ncmesh.hpp:190
static GeomInfo GI[Geometry::NumGeom]
Definition: ncmesh.hpp:1012
int slaves_end
slave faces
Definition: ncmesh.hpp:205
void MakeTopologyOnly()
Definition: ncmesh.hpp:413
void ReparentNode(int node, int new_p1, int new_p2)
Definition: ncmesh.cpp:275
static const int MaxElemNodes
Number of nodes of an element can have.
Definition: ncmesh.hpp:424
int PrintMemoryDetail() const
Definition: ncmesh.cpp:6276
int NewTriangle(int n0, int n1, int n2, int attr, int eattr0, int eattr1, int eattr2)
Definition: ncmesh.cpp:616
void DeleteUnusedFaces(const Array< int > &elemFaces)
Definition: ncmesh.cpp:410
Geometry::Type GetFaceGeometry(int index) const
Return face geometry type. index is the Mesh face number.
Definition: ncmesh.hpp:362
bool Iso
true if the mesh only contains isotropic refinements
Definition: ncmesh.hpp:420
Point(double x, double y, double z)
Definition: ncmesh.hpp:814
bool IsLeaf() const
Definition: ncmesh.hpp:502
void Swap(Array< T > &, Array< T > &)
Definition: array.hpp:638
std::map< std::string, int > RefPathMap
Definition: ncmesh.hpp:921
int QuadFaceSplitType(int v1, int v2, int v3, int v4, int mid[5]=NULL) const
Definition: ncmesh.cpp:2684
int GetFaceVerticesEdges(const MeshId &face_id, int vert_index[4], int edge_index[4], int edge_orientation[4]) const
Definition: ncmesh.cpp:5038
unsigned ghost
For internal use: 0 if regular fine element, 1 if parallel ghost element.
Definition: ncmesh.hpp:61
Geometry::Type GetElementGeometry(int index) const
Return element geometry type. index is the Mesh element number.
Definition: ncmesh.hpp:358
virtual void BuildFaceList()
Definition: ncmesh.cpp:3110
Nonconforming edge/face within a bigger edge/face.
Definition: ncmesh.hpp:214
void DerefineElement(int elem)
Derefine the element elem, does nothing on leaf elements.
Definition: ncmesh.cpp:1710
int TriFaceSplitLevel(int vn1, int vn2, int vn3) const
Definition: ncmesh.cpp:5268
const NCList & GetVertexList()
Definition: ncmesh.hpp:268
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
Point(double x)
Definition: ncmesh.hpp:808
void RegisterElement(int e)
Definition: ncmesh.cpp:421
NCMesh & operator=(NCMesh &)=delete
Copy assignment not supported.
Embedding(int elem, Geometry::Type geom, int matrix=0, bool ghost=false)
Definition: ncmesh.hpp:64
long TotalSize() const
Definition: ncmesh.cpp:3408
Array< int > leaf_sfc_index
natural tree ordering of leaf elements
Definition: ncmesh.hpp:541
const Point & operator()(int i) const
Definition: ncmesh.hpp:901
PointMatrix(const Point &p0, const Point &p1, const Point &p2, const Point &p3, const Point &p4, const Point &p5)
Definition: ncmesh.hpp:884
void DeleteLast()
Delete the last entry of the array.
Definition: array.hpp:196
void UpdateLeafElements()
Update the leaf elements indices in leaf_elements.
Definition: ncmesh.cpp:2123
int PrintVertexParents(std::ostream *out) const
Print the "vertex_parents" section of the mesh file.
Definition: ncmesh.cpp:5455
int RetrieveNode(const Element &el, int index)
Return el.node[index] correctly, even if the element is refined.
Definition: ncmesh.cpp:1668
int Geoms
bit mask of element geometries present, see InitGeomFlags()
Definition: ncmesh.hpp:421
PointMatrix(const Point &p0, const Point &p1)
Definition: ncmesh.hpp:868
int NGhostElements
Definition: ncmesh.hpp:538
void MakeCoarseToFineTable(Table &coarse_to_fine, bool want_ghosts=false) const
Definition: ncmesh.cpp:4681
void BuildElementToVertexTable()
Definition: ncmesh.cpp:3548
HashTable< Node > nodes
Definition: ncmesh.hpp:508
double a
Definition: lissajous.cpp:41
void ForgetElement(int e)
Definition: ncmesh.cpp:428
static PointMatrix pm_tet_identity
Definition: ncmesh.hpp:911
int GetMidFaceNode(int en1, int en2, int en3, int en4)
Definition: ncmesh.cpp:314
void RefineElement(int elem, char ref_type)
Definition: ncmesh.cpp:924
static const int MaxElemChildren
Number of children of an element can have.
Definition: ncmesh.hpp:430
int NewHexahedron(int n0, int n1, int n2, int n3, int n4, int n5, int n6, int n7, int attr, int fattr0, int fattr1, int fattr2, int fattr3, int fattr4, int fattr5)
Definition: ncmesh.cpp:473
void TraverseTetEdge(int vn0, int vn1, const Point &p0, const Point &p1, MatrixMap &matrix_map)
Definition: ncmesh.cpp:3014
Array< int > boundary_faces
subset of all faces, set by BuildFaceList
Definition: ncmesh.hpp:548
int NGhostFaces
Definition: ncmesh.hpp:538
Array< double > coordinates
Definition: ncmesh.hpp:521
const NCList & GetNCList(int entity)
Return vertex/edge/face list (entity = 0/1/2, respectively).
Definition: ncmesh.hpp:275
int Size() const
Returns the number of TYPE I elements.
Definition: table.hpp:92
int edges[MaxElemEdges][2]
Definition: ncmesh.hpp:1003
Array< int > leaf_elements
finest elements, in Mesh ordering (+ ghosts)
Definition: ncmesh.hpp:540
static PointMatrix pm_quad_identity
Definition: ncmesh.hpp:910
Array< Master > masters
Definition: ncmesh.hpp:230
int node[MaxElemNodes]
element corners (if ref_type == 0)
Definition: ncmesh.hpp:494
Array< MeshId > conforming
Definition: ncmesh.hpp:229
int nfv[MaxElemFaces]
Definition: ncmesh.hpp:1005
unsigned edge_flags
orientation flags, see OrientedPointMatrix
Definition: ncmesh.hpp:218
MFEM_DEPRECATED void GetCoarseToFineMap(const Mesh &fine_mesh, Table &coarse_to_fine) const
Definition: ncmesh.hpp:90
signed char geom
Geometry::Type (faces only) (char to save RAM)
Definition: ncmesh.hpp:192
Array< int > free_element_ids
Definition: ncmesh.hpp:512
void TraverseQuadFace(int vn0, int vn1, int vn2, int vn3, const PointMatrix &pm, int level, Face *eface[4], MatrixMap &matrix_map)
Definition: ncmesh.cpp:2895
int index(int i, int j, int nx, int ny)
Definition: life.cpp:235
bool operator==(const PointMatrix &pm) const
Definition: ncmesh.cpp:3974
int NGhostEdges
Definition: ncmesh.hpp:538
void GetElementFacesAttributes(int i, Array< int > &faces, Array< int > &fattr) const
Return the faces and face attributes of leaf element &#39;i&#39;.
Definition: ncmesh.cpp:5143
NCMesh(const Mesh *mesh)
Definition: ncmesh.cpp:105
int parent
Coarse Element index in the coarse mesh.
Definition: ncmesh.hpp:52
Point(double x, double y)
Definition: ncmesh.hpp:811
static PointMatrix pm_hex_identity
Definition: ncmesh.hpp:914
const MeshId & LookUp(int index, int *type=NULL) const
Definition: ncmesh.cpp:3413
virtual void CheckDerefinementNCLevel(const Table &deref_table, Array< int > &level_ok, int max_nc_level)
Definition: ncmesh.cpp:1952
int Size() const
Return the logical size of the array.
Definition: array.hpp:141
T & Last()
Return the last element in the array.
Definition: array.hpp:792
bool Boundary() const
Definition: ncmesh.hpp:469
void NeighborExpand(const Array< int > &elems, Array< int > &expanded, const Array< int > *search_set=NULL)
Definition: ncmesh.cpp:3822
BlockArray< Element > elements
Definition: ncmesh.hpp:511
void FindSetNeighbors(const Array< char > &elem_set, Array< int > *neighbors, Array< char > *neighbor_set=NULL)
Definition: ncmesh.cpp:3628
virtual void Update()
Definition: ncmesh.cpp:238
Slave(int index, int element, int local, int geom)
Definition: ncmesh.hpp:221
bool Unused() const
Definition: ncmesh.hpp:470
mfem::Element * NewMeshElement(int geom) const
Definition: ncmesh.cpp:2395
Data type point element.
Definition: point.hpp:22
void InitRootElements()
Count root elements and initialize root_state.
Definition: ncmesh.cpp:5714
int GetEdgeNCOrientation(const MeshId &edge_id) const
Definition: ncmesh.cpp:5026
HashTable< Node > shadow
temporary storage for reparented nodes
Definition: ncmesh.hpp:611
Defines the coarse-fine transformations of all fine elements.
Definition: ncmesh.hpp:70
virtual void BuildVertexList()
Definition: ncmesh.cpp:3333
const NCList & GetEdgeList()
Return the current list of conforming and nonconforming edges.
Definition: ncmesh.hpp:260
MeshId(int index, int element, int local, int geom=-1)
Definition: ncmesh.hpp:197
void UpdateElementToVertexTable()
Definition: ncmesh.hpp:782
void LoadBoundary(std::istream &input)
Load the "boundary" section of the mesh file.
Definition: ncmesh.cpp:5546
int NewPyramid(int n0, int n1, int n2, int n3, int n4, int attr, int fattr0, int fattr1, int fattr2, int fattr3, int fattr4)
Definition: ncmesh.cpp:560
void PrintCoordinates(std::ostream &out) const
Print the "coordinates" section of the mesh file.
Definition: ncmesh.cpp:5586
void CopyElements(int elem, const BlockArray< Element > &tmp_elements)
Definition: ncmesh.cpp:5938
virtual void ElementSharesEdge(int elem, int local, int enode)
Definition: ncmesh.hpp:743
int FindNodeExt(const Element &el, int node, bool abort=true)
Extended version of find_node: works if &#39;el&#39; is refined.
Definition: ncmesh.cpp:2758
int GetNVertices() const
Return the number of vertices in the NCMesh.
Definition: ncmesh.hpp:149
char ref_type
bit mask of X,Y,Z refinements (bits 0,1,2 respectively)
Definition: ncmesh.hpp:486
static const int MaxElemFaces
Number of faces of an element can have.
Definition: ncmesh.hpp:428
int NElements
Definition: ncmesh.hpp:534
void CollectQuadFaceVertices(int v0, int v1, int v2, int v3, Array< int > &indices)
Definition: ncmesh.cpp:3514
static int find_node(const Element &el, int node)
Definition: ncmesh.cpp:2748
void UnreferenceElement(int elem, Array< int > &elemFaces)
Definition: ncmesh.cpp:353
int index
Mesh number.
Definition: ncmesh.hpp:189
void LoadVertexParents(std::istream &input)
Load the vertex parent hierarchy from a mesh file.
Definition: ncmesh.cpp:5484
Point points[MaxElemNodes]
Definition: ncmesh.hpp:864
Rank 3 tensor (array of matrices)
Definition: densemat.hpp:1096
Abstract data type element.
Definition: element.hpp:28
void CollectIncidentElements(int elem, const RefCoord coord[3], Array< int > &list) const
Definition: ncmesh.cpp:3934
void CountSplits(int elem, int splits[3]) const
Definition: ncmesh.cpp:5312
Point & operator=(const Point &src)
Definition: ncmesh.hpp:837
Refinement()=default
int GetSingleElement() const
Return one of elem[0] or elem[1] and make sure the other is -1.
Definition: ncmesh.cpp:443
int rank
processor number (ParNCMesh), -1 if undefined/unknown
Definition: ncmesh.hpp:490
bool TraverseTriFace(int vn0, int vn1, int vn2, const PointMatrix &pm, int level, MatrixMap &matrix_map)
Definition: ncmesh.cpp:3051
HashTable< Face > faces
Definition: ncmesh.hpp:509
Defines the position of a fine element within a coarse element.
Definition: ncmesh.hpp:49
virtual void Refine(const Array< Refinement > &refinements)
Definition: ncmesh.cpp:1620
Array< char > face_geom
face geometry by face index, set by OnMeshUpdated
Definition: ncmesh.hpp:549