Line data Source code
1 : // The libMesh Finite Element Library.
2 : // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 :
4 : // This library is free software; you can redistribute it and/or
5 : // modify it under the terms of the GNU Lesser General Public
6 : // License as published by the Free Software Foundation; either
7 : // version 2.1 of the License, or (at your option) any later version.
8 :
9 : // This library is distributed in the hope that it will be useful,
10 : // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 : // Lesser General Public License for more details.
13 :
14 : // You should have received a copy of the GNU Lesser General Public
15 : // License along with this library; if not, write to the Free Software
16 : // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 :
18 :
19 :
20 : #ifndef LIBMESH_ELEM_H
21 : #define LIBMESH_ELEM_H
22 :
23 : // Local includes
24 : #include "libmesh/libmesh_common.h"
25 : #include "libmesh/bounding_box.h"
26 : #include "libmesh/dof_object.h"
27 : #include "libmesh/id_types.h"
28 : #include "libmesh/reference_counted_object.h"
29 : #include "libmesh/node.h"
30 : #include "libmesh/enum_elem_type.h" // INVALID_ELEM
31 : #include "libmesh/multi_predicates.h"
32 : #include "libmesh/pointer_to_pointer_iter.h"
33 : #include "libmesh/int_range.h"
34 : #include "libmesh/simple_range.h"
35 : #include "libmesh/variant_filter_iterator.h"
36 : #include "libmesh/hashword.h" // Used in compute_key() functions
37 :
38 : // C++ includes
39 : #include <algorithm>
40 : #include <cstddef>
41 : #include <iostream>
42 : #include <limits.h> // CHAR_BIT
43 : #include <set>
44 : #include <vector>
45 : #include <memory>
46 : #include <array>
47 :
48 : namespace libMesh
49 : {
50 :
51 : // Forward declarations
52 : class BoundaryInfo;
53 : class Elem;
54 : class MeshBase;
55 : class MeshRefinement;
56 : #ifdef LIBMESH_ENABLE_PERIODIC
57 : class PeriodicBoundaries;
58 : class PointLocatorBase;
59 : #endif
60 : template <class SideType, class ParentType>
61 : class Side;
62 : enum ElemQuality : int;
63 : enum IOPackage : int;
64 : enum Order : int;
65 :
66 :
67 : /**
68 : * This is the base class from which all geometric element types are
69 : * derived. The \p Elem class provides standard information such as
70 : * the number of nodes, edges, faces, vertices, children, and
71 : * neighbors it has, as well as access to (or the ability to
72 : * construct) these entities.
73 : *
74 : * An \p Elem has pointers to its \p Node objects. Some of these
75 : * nodes live at the vertices of the element, while others may live on
76 : * edges (and faces in 3D), or interior to the element. The number of
77 : * nodes in a given element, \p n_nodes(), is encoded into the name of
78 : * the class. For example, a \p Tri3 has three nodes which correspond
79 : * to the vertices, while a \p Tri6 has six nodes, three of which are
80 : * located at vertices, and three which are located at the midpoint of
81 : * each edge. Nodes on edges, faces, and element interiors are called
82 : * second-order nodes.
83 : *
84 : * A 1D Elem is an \p Edge, a 2D Elem is a \p Face, and a 3D Elem is a
85 : * \p Cell. An \p Elem is composed of a number of sides, which can
86 : * be accessed as dim-1 dimensional \p Elem types. For example, a \p
87 : * Hex8 is a 3D hexahedral element. A \p Hex8 has 6 sides, which are
88 : * \p Faces of type Quad4.
89 : *
90 : * \author Benjamin S. Kirk
91 : * \date 2002-2007
92 : * \brief The base class for all geometric element types.
93 : */
94 : class Elem : public ReferenceCountedObject<Elem>,
95 : public DofObject
96 : {
97 : protected:
98 :
99 : /**
100 : * Constructor. Creates an element with \p n_nodes nodes,
101 : * \p n_sides sides, \p n_children possible children, and
102 : * parent \p p. The constructor allocates the memory necessary
103 : * to support this data.
104 : */
105 : Elem (const unsigned int n_nodes,
106 : const unsigned int n_sides,
107 : Elem * parent,
108 : Elem ** elemlinkdata,
109 : Node ** nodelinkdata);
110 :
111 : public:
112 :
113 : /**
114 : * Elems are responsible for allocating and deleting space for
115 : * storing pointers to their children during refinement, so they
116 : * cannot currently be (default) copy-constructed or copy-
117 : * assigned. We therefore explicitly delete these operations. In
118 : * addition, the DofObject base class currently has private copy
119 : * construction and assignment operators, so that prevents us from
120 : * copying Elems as well.
121 : */
122 : Elem (Elem &&) = delete;
123 : Elem (const Elem &) = delete;
124 : Elem & operator= (const Elem &) = delete;
125 : Elem & operator= (Elem &&) = delete;
126 :
127 : /**
128 : * Destructor.
129 : */
130 4293329305 : virtual ~Elem() = default;
131 :
132 : /**
133 : * \returns The \p Point associated with local \p Node \p i.
134 : */
135 : const Point & point (const unsigned int i) const;
136 :
137 : /**
138 : * \returns The \p Point associated with local \p Node \p i
139 : * as a writable reference.
140 : */
141 : Point & point (const unsigned int i);
142 :
143 : /**
144 : * \returns The \p Point associated with local \p Node \p i,
145 : * in master element rather than physical coordinates.
146 : */
147 : virtual Point master_point (const unsigned int i) const = 0;
148 :
149 : /**
150 : * \returns The global id number of local \p Node \p i.
151 : */
152 : dof_id_type node_id (const unsigned int i) const;
153 :
154 : /**
155 : * \returns The local id number of global \p Node id \p i,
156 : * or \p invalid_uint if Node id \p i is not local.
157 : */
158 : unsigned int local_node (const dof_id_type i) const;
159 :
160 : /**
161 : * \returns The local index for the \p Node pointer \p node_ptr,
162 : * or \p invalid_uint if \p node_ptr is not a local node.
163 : */
164 : unsigned int get_node_index (const Node * node_ptr) const;
165 :
166 : /**
167 : * \returns A pointer to an array of local node pointers.
168 : */
169 : const Node * const * get_nodes () const;
170 :
171 : /**
172 : * \returns A const pointer to local \p Node \p i.
173 : */
174 : const Node * node_ptr (const unsigned int i) const;
175 :
176 : /**
177 : * \returns A non-const pointer to local \p Node \p i.
178 : */
179 : Node * node_ptr (const unsigned int i);
180 :
181 : /**
182 : * \returns A const reference to local \p Node \p i.
183 : */
184 : const Node & node_ref (const unsigned int i) const;
185 :
186 : /**
187 : * \returns A writable reference to local \p Node \p i.
188 : */
189 : Node & node_ref (const unsigned int i);
190 :
191 : #ifdef LIBMESH_ENABLE_DEPRECATED
192 : /**
193 : * \returns The pointer to the \p Node with local number \p i as a
194 : * writable reference.
195 : *
196 : * \deprecated This setter cannot update the multiple node pointers
197 : * used in a general polyhedron; use the \p set_node overload that
198 : * takes an argument.
199 : */
200 : virtual Node * & set_node (const unsigned int i);
201 : #endif // LIBMESH_ENABLE_DEPRECATED
202 :
203 : /**
204 : * Sets local \p Node \p i to refer to \p node.
205 : */
206 : virtual void set_node (const unsigned int i,
207 : Node * node);
208 :
209 : /**
210 : * Nested classes for use iterating over all nodes of an element.
211 : */
212 : class NodeRefIter;
213 : class ConstNodeRefIter;
214 :
215 : /**
216 : * Returns a range with all nodes of an element, usable in
217 : * range-based for loops. The exact type of the return value here
218 : * may be subject to change in future libMesh releases, but the
219 : * iterators will always dereference to produce a reference to a
220 : * Node.
221 : */
222 : SimpleRange<NodeRefIter> node_ref_range();
223 :
224 : SimpleRange<ConstNodeRefIter> node_ref_range() const;
225 :
226 : /**
227 : * \returns The subdomain that this element belongs to.
228 : */
229 : subdomain_id_type subdomain_id () const;
230 :
231 : /**
232 : * \returns The subdomain that this element belongs to as a
233 : * writable reference.
234 : */
235 : subdomain_id_type & subdomain_id ();
236 :
237 : /**
238 : * A static integral constant representing an invalid subdomain id.
239 : * See also DofObject::{invalid_id, invalid_unique_id, invalid_processor_id}.
240 : *
241 : * \note We don't use the static_cast(-1) trick here since
242 : * \p subdomain_id_type is sometimes a *signed* integer for
243 : * compatibility reasons (see libmesh/id_types.h).
244 : *
245 : * \note Normally you can declare static const integral types
246 : * directly in the header file (C++ standard, 9.4.2/4) but
247 : * std::numeric_limits<T>::max() is not considered a "constant
248 : * expression". This one is therefore defined in elem.C.
249 : * http://stackoverflow.com/questions/2738435/using-numeric-limitsmax-in-constant-expressions
250 : */
251 : static const subdomain_id_type invalid_subdomain_id;
252 :
253 : /**
254 : * \returns true iff this element type can vary in topology (e.g.
255 : * have different numbers of sides and/or nodes) at runtime. For
256 : * such general polygons or polyhedra, APIs which assume a fixed
257 : * topology are not safe to use.
258 : */
259 285566633 : virtual bool runtime_topology() const { return false; }
260 :
261 : /**
262 : * \returns A pointer to the "reference element" associated
263 : * with this element. The reference element is the image of this
264 : * element in reference parametric space. Importantly, it is *not*
265 : * an actual element in the mesh, but rather a Singleton-type
266 : * object, so for example all \p Quad4 elements share the same
267 : * \p reference_elem().
268 : *
269 : * If the element is of a type that can admit multiple topologies,
270 : * such as a Polygon subtype, then there is no reference element;
271 : * for such types this method should not be used.
272 : */
273 : const Elem * reference_elem () const;
274 :
275 : /**
276 : * \returns An id associated with the \p s side of this element.
277 : * The id is not necessarily unique, but should be close.
278 : */
279 : virtual dof_id_type key (const unsigned int s) const = 0;
280 :
281 : /**
282 : * \returns An id associated with the \p s side of this element, as
283 : * defined solely by element vertices. The id is not necessarily
284 : * unique, but should be close. This is particularly useful in the
285 : * \p MeshBase::find_neighbors() routine.
286 : */
287 : virtual dof_id_type low_order_key (const unsigned int s) const = 0;
288 :
289 : /**
290 : * \returns An id associated with the global node ids of this
291 : * element. The id is not necessarily unique, but should be
292 : * close. Uses the same hash as the key(s) function, so for example
293 : * if "tri3" is side 0 of "tet4", then tri3->key()==tet4->key(0).
294 : */
295 : virtual dof_id_type key () const;
296 :
297 : /**
298 : * \returns \p true if two elements are equivalent, false otherwise.
299 : * This is true if the elements are connected to identical global
300 : * nodes, regardless of how those nodes might be numbered local
301 : * to the elements.
302 : */
303 : bool operator == (const Elem & rhs) const;
304 :
305 : /**
306 : * \returns \p true if two elements have equal topologies, false
307 : * otherwise.
308 : * This is true if the elements connect to nodes of the same id in
309 : * the same order, and neighbors of the same id on each side, the
310 : * same id on any parent and/or interior_parent link, etc.
311 : */
312 : bool topologically_equal (const Elem & rhs) const;
313 :
314 : /**
315 : * \returns A const pointer to the \f$ i^{th} \f$ neighbor of this
316 : * element, or \p nullptr if \p MeshBase::find_neighbors() has not been
317 : * called.
318 : *
319 : * \note If \p MeshBase::find_neighbors() has been called and this
320 : * function still returns \p nullptr, then the side is on a boundary of
321 : * the domain.
322 : */
323 : const Elem * neighbor_ptr (unsigned int i) const;
324 :
325 : /**
326 : * \returns A non-const pointer to the \f$ i^{th} \f$ neighbor of this element.
327 : */
328 : Elem * neighbor_ptr (unsigned int i);
329 :
330 : /**
331 : * Nested "classes" for use iterating over all neighbors of an element.
332 : */
333 : typedef Elem * const * NeighborPtrIter;
334 : typedef const Elem * const * ConstNeighborPtrIter;
335 :
336 : /**
337 : * Returns a range with all neighbors of an element, usable in
338 : * range-based for loops. The exact type of the return value here
339 : * may be subject to change in future libMesh releases, but the
340 : * iterators will always dereference to produce a pointer to a
341 : * neighbor element (or a null pointer, for sides which have no
342 : * neighbors).
343 : */
344 : SimpleRange<NeighborPtrIter> neighbor_ptr_range();
345 :
346 : SimpleRange<ConstNeighborPtrIter> neighbor_ptr_range() const;
347 :
348 : #ifdef LIBMESH_ENABLE_PERIODIC
349 : /**
350 : * \returns A pointer to the \f$ i^{th} \f$ neighbor of this element
351 : * for interior elements. If an element is on a periodic
352 : * boundary, it will return a corresponding element on the opposite
353 : * side.
354 : */
355 : const Elem * topological_neighbor (const unsigned int i,
356 : const MeshBase & mesh,
357 : const PointLocatorBase & point_locator,
358 : const PeriodicBoundaries * pb) const;
359 :
360 : /**
361 : * \returns A writable pointer to the \f$ i^{th} \f$ neighbor of
362 : * this element for interior elements. If an element is on a
363 : * periodic boundary, it will return a corresponding element on the
364 : * opposite side.
365 : */
366 : Elem * topological_neighbor (const unsigned int i,
367 : MeshBase & mesh,
368 : const PointLocatorBase & point_locator,
369 : const PeriodicBoundaries * pb);
370 :
371 : /**
372 : * \returns \p true if the element \p elem in question is a neighbor or
373 : * topological neighbor of this element, \p false otherwise.
374 : */
375 : bool has_topological_neighbor (const Elem * elem,
376 : const MeshBase & mesh,
377 : const PointLocatorBase & point_locator,
378 : const PeriodicBoundaries * pb) const;
379 : #endif
380 :
381 : /**
382 : * Assigns \p n as the \f$ i^{th} \f$ neighbor.
383 : */
384 : void set_neighbor (const unsigned int i, Elem * n);
385 :
386 : /**
387 : * \returns \p true if the element \p elem in question is a neighbor
388 : * of this element, \p false otherwise.
389 : */
390 : bool has_neighbor (const Elem * elem) const;
391 :
392 : /**
393 : * \returns If \p elem is a neighbor of a child of this element, a
394 : * pointer to that child, otherwise \p nullptr.
395 : */
396 : Elem * child_neighbor (Elem * elem);
397 :
398 : /**
399 : * \returns If \p elem is a neighbor of a child of this element, a
400 : * pointer to that child, otherwise \p nullptr.
401 : */
402 : const Elem * child_neighbor (const Elem * elem) const;
403 :
404 : /**
405 : * \returns \p true if this element has a side coincident
406 : * with a boundary (indicated by a \p nullptr neighbor), \p false
407 : * otherwise.
408 : */
409 : bool on_boundary () const;
410 :
411 : /**
412 : * \returns \p true if this element is "semilocal" to the calling
413 : * processor, which must specify its rank.
414 : *
415 : * This method is discouraged, as it uses the *old* definition of
416 : * semilocal (elements which are not local but which are point
417 : * neighbors of something local) rather than any of the new
418 : * definitions discussed in ghosting_functor.h
419 : */
420 : bool is_semilocal (const processor_id_type my_pid) const;
421 :
422 : /**
423 : * This function tells you which neighbor \p e is.
424 : * I.e. if s = a->which_neighbor_am_i(e); then
425 : * a->neighbor(s) will be an ancestor of e.
426 : */
427 : unsigned int which_neighbor_am_i(const Elem * e) const;
428 :
429 : /**
430 : * This function tells you which side the boundary element \p e is.
431 : * I.e. if e = a->build_side_ptr(s) or e = a->side_ptr(s); then
432 : * a->which_side_am_i(e) will be s.
433 : *
434 : * \note An \e exact floating point comparison of the nodal
435 : * positions of \p e is made with the nodal positions of \p this in
436 : * order to perform this test. The idea is that the test will return
437 : * a valid side id if \p e either directly shares Node pointers with
438 : * \p this, or was created by exactly copying some of the nodes of
439 : * \p this (e.g. through BoundaryMesh::sync()). In these
440 : * circumstances, non-fuzzy floating point equality is expected.
441 : *
442 : * \returns The side of \p this the element which \p e is, otherwise
443 : * \p invalid_uint.
444 : */
445 : unsigned int which_side_am_i(const Elem * e) const;
446 :
447 : /**
448 : * \returns The local node id for node \p side_node on side \p side of
449 : * this Elem. Simply relies on the \p side_nodes_map for each of the
450 : * derived types. For example,
451 : * Tri3::local_side_node(0, 0) -> 0
452 : * Tri3::local_side_node(0, 1) -> 1
453 : * Tri3::local_side_node(1, 0) -> 1
454 : * Tri3::local_side_node(1, 1) -> 2
455 : * etc...
456 : */
457 : virtual unsigned int local_side_node(unsigned int side,
458 : unsigned int side_node) const = 0;
459 :
460 : /**
461 : * Similar to Elem::local_side_node(), but instead of a side id, takes
462 : * an edge id and a node id on that edge and returns a local node number
463 : * for the Elem. The implementation relies on the "edge_nodes_map" tables
464 : * for 3D elements. For 2D elements, calls local_side_node(). Throws an
465 : * error if called on 1D elements.
466 : */
467 : virtual unsigned int local_edge_node(unsigned int edge,
468 : unsigned int edge_node) const = 0;
469 :
470 : #ifdef LIBMESH_ENABLE_DEPRECATED
471 : /**
472 : * This function is deprecated, call local_side_node(side, side_node) instead.
473 : */
474 : unsigned int which_node_am_i(unsigned int side,
475 : unsigned int side_node) const;
476 : #endif // LIBMESH_ENABLE_DEPRECATED
477 :
478 : /**
479 : * \returns \p true if a vertex of \p e is contained
480 : * in this element. If \p mesh_connection is true, looks
481 : * specifically for containment possibilities of an element \p e
482 : * that is connected to \p this via membership in the same manifold
483 : * of the same mesh.
484 : */
485 : bool contains_vertex_of(const Elem * e, bool mesh_connection=false) const;
486 :
487 : /**
488 : * \returns \p true if an edge of \p e is contained in
489 : * this element. (Internally, this is done by checking whether at
490 : * least two vertices of \p e are contained in this element).
491 : */
492 : bool contains_edge_of(const Elem * e) const;
493 :
494 : /**
495 : * This function finds all active elements (including this one)
496 : * which are in the same manifold as this element and which touch
497 : * the current active element at the specified point, which should
498 : * be a point in the current element.
499 : *
500 : * Elements which are not "in the same manifold" (e.g. the
501 : * interior_parent of a boundary element) will not be found with
502 : * this method.
503 : *
504 : * Elements which overlap the specified point but which are only
505 : * connected to the current element via elements which do not
506 : * overlap that point (e.g. in a folded or tangled mesh) are not
507 : * considered to "touch" the current element and will not be found
508 : * with this method.
509 : */
510 : void find_point_neighbors(const Point & p,
511 : std::set<const Elem *> & neighbor_set) const;
512 :
513 : /**
514 : * This function finds all active elements (including this one) in
515 : * the same manifold as this element which touch this active element
516 : * at any point.
517 : */
518 : void find_point_neighbors(std::set<const Elem *> & neighbor_set) const;
519 :
520 : /**
521 : * This function finds all active elements (including this one) in
522 : * the same manifold as start_elem (which must be active and must
523 : * touch this element) which touch this element at any point.
524 : */
525 : void find_point_neighbors(std::set<const Elem *> & neighbor_set,
526 : const Elem * start_elem) const;
527 :
528 : /**
529 : * Non-const version of function above. Fills a set of non-const Elem pointers.
530 : */
531 : void find_point_neighbors(std::set<Elem *> & neighbor_set,
532 : Elem * start_elem);
533 :
534 : /**
535 : * This function finds all active elements in the same manifold as
536 : * this element which touch the current active element along the
537 : * whole edge defined by the two points \p p1 and \p p2.
538 : */
539 : void find_edge_neighbors(const Point & p1,
540 : const Point & p2,
541 : std::set<const Elem *> & neighbor_set) const;
542 :
543 : /**
544 : * This function finds all active elements in the same manifold as
545 : * this element which touch the current active element along any
546 : * edge (more precisely, at at least two points).
547 : *
548 : * In this case, elements are included even if they do not touch a
549 : * *whole* edge of this element.
550 : */
551 : void find_edge_neighbors(std::set<const Elem *> & neighbor_set) const;
552 :
553 : /**
554 : * This function finds all active elements (*not* including this
555 : * one) in the parent manifold of this element whose intersection
556 : * with this element has non-zero measure.
557 : */
558 : void find_interior_neighbors(std::set<const Elem *> & neighbor_set) const;
559 :
560 : /**
561 : * Non-const version of function above that fills up a vector of
562 : * non-const Elem pointers instead.
563 : */
564 : void find_interior_neighbors(std::set<Elem *> & neighbor_set);
565 :
566 : /**
567 : * Resets this element's neighbors' appropriate neighbor pointers
568 : * and its parent's and children's appropriate pointers
569 : * to point to null instead of to this.
570 : *
571 : * To be used before an element is deleted from a mesh.
572 : */
573 : void remove_links_to_me ();
574 :
575 : /**
576 : * Resets this element's neighbors' appropriate neighbor pointers
577 : * and its parent's and children's appropriate pointers
578 : * to point to the global remote_elem instead of this.
579 : * Used by the library before an element becomes remote on the
580 : * local processor.
581 : */
582 : void make_links_to_me_remote ();
583 :
584 : /**
585 : * Resets the \p neighbor_side pointers of our nth neighbor (and
586 : * its descendants, if appropriate) to point to this Elem instead of
587 : * to the global remote_elem. Used by the library when a formerly
588 : * remote element is being added to the local processor.
589 : */
590 : void make_links_to_me_local (unsigned int n, unsigned int neighbor_side);
591 :
592 : /**
593 : * \returns \p true if this element is remote, false otherwise.
594 : *
595 : * A remote element (see \p RemoteElem) is a syntactic convenience --
596 : * it is a placeholder for an element which exists on some other
597 : * processor. Local elements are required to have valid neighbors,
598 : * and these ghost elements may have remote neighbors for data
599 : * structure consistency. The use of remote elements helps ensure
600 : * that any element we may access has a \p nullptr neighbor only if it
601 : * lies on the physical boundary of the domain.
602 : */
603 7280126 : virtual bool is_remote () const
604 7280126 : { return false; }
605 :
606 : /**
607 : * \returns The connectivity for this element in a specific
608 : * format, which is specified by the IOPackage tag.
609 : */
610 : virtual void connectivity(const unsigned int sc,
611 : const IOPackage iop,
612 : std::vector<dof_id_type> & conn) const = 0;
613 :
614 : /**
615 : * Writes the element connectivity for various IO packages
616 : * to the passed ostream "out". Not virtual, since it is
617 : * implemented in the base class.
618 : */
619 : void write_connectivity (std::ostream & out,
620 : const IOPackage iop) const;
621 :
622 : /**
623 : * \returns The type of element that has been derived from this
624 : * base class.
625 : */
626 : virtual ElemType type () const = 0;
627 :
628 : /**
629 : * This array maps the integer representation of the \p ElemType enum
630 : * to the geometric dimension of the element.
631 : *
632 : * This is currently usable even for complicated subclasses with
633 : * runtime-varying topology.
634 : */
635 : static const unsigned int type_to_dim_map[INVALID_ELEM];
636 :
637 : /**
638 : * \returns The dimensionality of the object.
639 : */
640 : virtual unsigned short dim () const = 0;
641 :
642 : /**
643 : * This array maps the integer representation of the \p ElemType enum
644 : * to the number of nodes in the element.
645 : *
646 : * This is only usable for simple types for which the node number
647 : * is fixed; for more general types like Polygon subclasses an actual
648 : * instantiated Elem must be queried.
649 : */
650 : static const unsigned int type_to_n_nodes_map[INVALID_ELEM];
651 :
652 : /**
653 : * \returns The number of nodes this element contains.
654 : */
655 : virtual unsigned int n_nodes () const = 0;
656 :
657 : /**
658 : * The maximum number of nodes *any* element can contain.
659 : * This is useful for replacing heap vectors with stack arrays.
660 : */
661 : static const unsigned int max_n_nodes = 27;
662 :
663 : /**
664 : * \returns An integer range from 0 up to (but not including)
665 : * the number of nodes this element contains.
666 : */
667 : IntRange<unsigned short> node_index_range () const;
668 :
669 : /**
670 : * \returns The number of nodes the given child of this element
671 : * contains. Except in odd cases like pyramid refinement this will
672 : * be the same as the number of nodes in the parent element.
673 : */
674 18635148 : virtual unsigned int n_nodes_in_child (unsigned int /*c*/) const
675 18635148 : { return this->n_nodes(); }
676 :
677 : /**
678 : * This array maps the integer representation of the \p ElemType enum
679 : * to the number of sides on the element.
680 : *
681 : * This is only usable for simple types for which the node number
682 : * is fixed; for more general types like Polygon subclasses an actual
683 : * instantiated Elem must be queried.
684 : */
685 : static const unsigned int type_to_n_sides_map[INVALID_ELEM];
686 :
687 : /**
688 : * \returns The number of sides the element that has been derived
689 : * from this class has. In 2D the number of sides is the number
690 : * of edges, in 3D the number of sides is the number of faces.
691 : */
692 : virtual unsigned int n_sides () const = 0;
693 :
694 : /**
695 : * \returns The type of element for side \p s.
696 : */
697 : virtual ElemType side_type (const unsigned int s) const = 0;
698 :
699 : /**
700 : * \returns An integer range from 0 up to (but not including)
701 : * the number of sides this element has.
702 : */
703 : IntRange<unsigned short> side_index_range () const;
704 :
705 : /**
706 : * \returns The number of neighbors the element that has been derived
707 : * from this class has.
708 : *
709 : * Only face (or edge in 2D) neighbors are stored, so this method
710 : * returns n_sides(). At one point we intended to allow derived
711 : * classes to override this, but too much current libMesh code
712 : * assumes n_neighbors==n_sides.
713 : */
714 85301839 : unsigned int n_neighbors () const
715 152196476 : { return this->n_sides(); }
716 :
717 : /**
718 : * \returns The number of vertices the element that has been derived
719 : * from this class has.
720 : */
721 : virtual unsigned int n_vertices () const = 0;
722 :
723 : /**
724 : * \returns The number of edges the element that has been derived
725 : * from this class has.
726 : */
727 : virtual unsigned int n_edges () const = 0;
728 :
729 : /**
730 : * \returns An integer range from 0 up to (but not including)
731 : * the number of edges this element has.
732 : */
733 : IntRange<unsigned short> edge_index_range () const;
734 :
735 : /**
736 : * This array maps the integer representation of the \p ElemType enum
737 : * to the number of edges on the element.
738 : *
739 : * This is only usable for simple types for which the node number
740 : * is fixed; for more general types like Polygon subclasses an actual
741 : * instantiated Elem must be queried.
742 : */
743 : static const unsigned int type_to_n_edges_map[INVALID_ELEM];
744 :
745 : /**
746 : * \returns The number of faces the element that has been derived
747 : * from this class has.
748 : */
749 : virtual unsigned int n_faces () const = 0;
750 :
751 : /**
752 : * \returns An integer range from 0 up to (but not including)
753 : * the number of faces this element has.
754 : */
755 : IntRange<unsigned short> face_index_range () const;
756 :
757 : /**
758 : * \returns The number of children the element that has been derived
759 : * from this class may have.
760 : */
761 : virtual unsigned int n_children () const = 0;
762 :
763 : /**
764 : * \returns \p true if the specified (local) node number is a vertex node.
765 : */
766 : virtual bool is_vertex(const unsigned int i) const = 0;
767 :
768 : /**
769 : * \returns \p true if the specified child has a vertex at the
770 : * specified (child-local) node number.
771 : * Except in odd cases like pyramid refinement the child will have
772 : * the same local structure as the parent element.
773 : */
774 169263 : virtual bool is_vertex_on_child (unsigned int /*c*/,
775 : unsigned int n) const
776 169263 : { return this->is_vertex(n); }
777 :
778 : /**
779 : * \returns \p true if this element has a vertex at the specified
780 : * (child-local) node number \p n of the specified child \p c.
781 : */
782 : virtual bool is_vertex_on_parent(unsigned int c,
783 : unsigned int n) const;
784 :
785 : /**
786 : * \returns \p true if the specified (local) node number is an edge node.
787 : * For 1D elements, is_edge() is equivalent to is_internal().
788 : */
789 : virtual bool is_edge(const unsigned int i) const = 0;
790 :
791 : /**
792 : * \returns \p true if the specified (local) node number is a face node.
793 : * For 2D elements, is_face() is equivalent to is_internal().
794 : * For 1D elements, is_face() == false.
795 : */
796 : virtual bool is_face(const unsigned int i) const = 0;
797 :
798 : /**
799 : * \returns \p true if the specified (local) node number is an internal node.
800 : */
801 : bool is_internal(const unsigned int i) const;
802 :
803 : /**
804 : * \returns \p true if the specified (local) node number is on the
805 : * specified side.
806 : */
807 : virtual bool is_node_on_side(const unsigned int n,
808 : const unsigned int s) const = 0;
809 :
810 : /**
811 : * \returns the (local) node numbers on the specified side
812 : */
813 : virtual std::vector<unsigned int> nodes_on_side(const unsigned int /*s*/) const = 0;
814 :
815 : /**
816 : * \returns the (local) node numbers on the specified edge
817 : */
818 : virtual std::vector<unsigned int> nodes_on_edge(const unsigned int /*e*/) const = 0;
819 :
820 : /**
821 : * \returns the (local) side numbers that touch the specified edge
822 : */
823 : virtual std::vector<unsigned int> sides_on_edge(const unsigned int /*e*/) const = 0;
824 :
825 : /**
826 : * \returns the (local) edge numbers that touch the specified node
827 : */
828 : virtual std::vector<unsigned int> edges_adjacent_to_node(const unsigned int /*n*/) const = 0;
829 :
830 : /**
831 : * \returns \p true if the specified (local) node number is on the
832 : * specified edge.
833 : */
834 : virtual bool is_node_on_edge(const unsigned int n,
835 : const unsigned int e) const = 0;
836 :
837 : /**
838 : * \returns \p true if the specified edge is on the specified side.
839 : */
840 : virtual bool is_edge_on_side(const unsigned int e,
841 : const unsigned int s) const = 0;
842 :
843 : /**
844 : * \returns The side number opposite to \p s (for a tensor product
845 : * element), or throws an error otherwise.
846 : */
847 : virtual unsigned int opposite_side(const unsigned int s) const;
848 :
849 : /**
850 : * \returns The local node number for the node opposite to node n
851 : * on side \p opposite_side(s) (for a tensor product element), or
852 : * throws an error otherwise.
853 : */
854 : virtual unsigned int opposite_node(const unsigned int n,
855 : const unsigned int s) const;
856 :
857 : /**
858 : * \returns The number of sub-elements this element may be broken
859 : * down into for visualization purposes. For example, 1 for a
860 : * linear triangle, 4 for a quadratic (6-noded) triangle, etc...
861 : */
862 : virtual unsigned int n_sub_elem () const = 0;
863 :
864 : /**
865 : * \returns A temporary element coincident with side \p i.
866 : *
867 : * This method returns the _minimum_ element necessary to uniquely
868 : * identify the side. For example, the side of a hexahedron is
869 : * always returned as a 4-noded quadrilateral, regardless of what
870 : * type of hex you are dealing with. Important data like subdomain
871 : * id, p level, or mapping type may be omitted from the temporary
872 : * element. If you want a first-class full-ordered face (i.e. a
873 : * 9-noded quad face for a 27-noded hexahedron), use the
874 : * build_side_ptr method.
875 : *
876 : * \note The const version of this function is non-virtual; it
877 : * simply calls the virtual non-const version and const_casts the
878 : * return type.
879 : */
880 : virtual std::unique_ptr<Elem> side_ptr (unsigned int i) = 0;
881 : std::unique_ptr<const Elem> side_ptr (unsigned int i) const;
882 :
883 : /**
884 : * Resets the loose element \p side, which may currently point to a
885 : * different side than \p i or even a different element than \p
886 : * this, to point to side \p i on \p this. If \p side is currently
887 : * an element of the wrong type, it will be freed and a new element
888 : * allocated; otherwise no memory allocation will occur.
889 : *
890 : * This will cause \p side to be a minimum-ordered element, even if
891 : * it is handed a higher-ordered element that must be replaced.
892 : *
893 : * The const version of this function is non-virtual; it simply
894 : * calls the virtual non-const version and const_casts the return
895 : * type.
896 : */
897 : virtual void side_ptr (std::unique_ptr<Elem> & side, const unsigned int i) = 0;
898 : void side_ptr (std::unique_ptr<const Elem> & side, const unsigned int i) const;
899 :
900 : /**
901 : * \returns An temporary element coincident with side \p i wrapped
902 : * in a smart pointer.
903 : *
904 : * The element returned is full-ordered and full-featured, in
905 : * contrast to the side method. For example, calling
906 : * build_side_ptr(0) on a 20-noded hex in subdomain 5 will build a
907 : * 8-noded quadrilateral coincident with face 0, assign it subdomain
908 : * id 5, and pass back the pointer.
909 : *
910 : * The side element's id() is undefined; it is a temporary element
911 : * not added to any mesh.
912 : *
913 : * A \p std::unique_ptr<Elem> is returned to prevent a memory leak.
914 : * This way the user need not remember to delete the object.
915 : *
916 : * The const version of this function is non-virtual; it simply
917 : * calls the virtual non-const version and const_casts the return
918 : * type.
919 : */
920 : virtual std::unique_ptr<Elem> build_side_ptr (const unsigned int i) = 0;
921 : std::unique_ptr<const Elem> build_side_ptr (const unsigned int i) const;
922 :
923 : #ifdef LIBMESH_ENABLE_DEPRECATED
924 : /*
925 : * Older versions of libMesh supported a "proxy" option here.
926 : */
927 0 : virtual std::unique_ptr<Elem> build_side_ptr (const unsigned int i, bool proxy)
928 0 : { if (proxy) libmesh_error(); libmesh_deprecated(); return this->build_side_ptr(i); }
929 :
930 : std::unique_ptr<const Elem> build_side_ptr (const unsigned int i, bool proxy) const
931 : { if (proxy) libmesh_error(); libmesh_deprecated(); return this->build_side_ptr(i); }
932 : #endif
933 :
934 : /**
935 : * Resets the loose element \p side, which may currently point to a
936 : * different side than \p i or even a different element than \p
937 : * this, to point to side \p i on \p this. If \p side is currently
938 : * an element of the wrong type, it will be freed and a new element
939 : * allocated; otherwise no memory allocation will occur.
940 : *
941 : * This will cause \p side to be a full-ordered element, even if it
942 : * is handed a lower-ordered element that must be replaced.
943 : *
944 : * The const version of this function is non-virtual; it simply
945 : * calls the virtual non-const version and const_casts the return
946 : * type.
947 : */
948 : virtual void build_side_ptr (std::unique_ptr<Elem> & side, const unsigned int i) = 0;
949 : void build_side_ptr (std::unique_ptr<const Elem> & side, const unsigned int i) const;
950 :
951 : /**
952 : * \returns An element coincident with edge \p i wrapped in a smart pointer.
953 : *
954 : * The element returned is full-ordered. For example, calling
955 : * build_edge_ptr(0) on a 20-noded hex will build a 3-noded edge
956 : * coincident with edge 0 and pass back the pointer. A \p
957 : * std::unique_ptr<Elem> is returned to prevent a memory leak. This way
958 : * the user need not remember to delete the object.
959 : *
960 : * The const version of this function is non-virtual; it simply
961 : * calls the virtual non-const version and const_casts the return
962 : * type.
963 : */
964 : virtual std::unique_ptr<Elem> build_edge_ptr (const unsigned int i) = 0;
965 : std::unique_ptr<const Elem> build_edge_ptr (const unsigned int i) const;
966 :
967 : /**
968 : * Resets the loose element \p edge, which may currently point to a
969 : * different edge than \p i or even a different element than \p
970 : * this, to point to edge \p i on \p this. If \p edge is currently
971 : * an element of the wrong type, it will be freed and a new element
972 : * allocated; otherwise no memory allocation will occur.
973 : *
974 : * This will cause \p edge to be a full-ordered element, even if it
975 : * is handed a lower-ordered element that must be replaced.
976 : *
977 : * The const version of this function is non-virtual; it simply
978 : * calls the virtual non-const version and const_casts the return
979 : * type.
980 : */
981 : virtual void build_edge_ptr (std::unique_ptr<Elem> & edge, const unsigned int i) = 0;
982 : void build_edge_ptr (std::unique_ptr<const Elem> & edge, const unsigned int i) const;
983 :
984 : /**
985 : * This array maps the integer representation of the \p ElemType enum
986 : * to the default approximation order of elements of that type.
987 : *
988 : * This is currently usable even for complicated subclasses with
989 : * runtime-varying topology.
990 : */
991 : static const Order type_to_default_order_map[INVALID_ELEM];
992 :
993 : /**
994 : * \returns The default approximation order for this element. This
995 : * is the order that will be used to compute the map to the
996 : * reference element.
997 : */
998 : virtual Order default_order () const = 0;
999 :
1000 : /**
1001 : * \returns The maximum supported approximation order for nodal
1002 : * (Lagrange or Rational Bezier-Bernstein) variables on this element
1003 : * type. This is usually the same as the default order.
1004 : */
1005 38570746 : virtual Order supported_nodal_order() const { return default_order(); }
1006 :
1007 : /**
1008 : * \returns The default approximation order for side elements of
1009 : * this element type. This may be lower for elements with 'bubble
1010 : * functions' in the Lagrange basis.
1011 : */
1012 481188 : virtual Order default_side_order () const { return default_order(); }
1013 :
1014 : #ifdef LIBMESH_ENABLE_DEPRECATED
1015 : /**
1016 : * Calls Elem::vertex_average() for backwards compatibility.
1017 : *
1018 : * \deprecated This method has now been deprecated, so it will be
1019 : * removed at some point in the future. Calls to Elem::centroid()
1020 : * in user code should be updated to either call
1021 : * Elem::vertex_average() or Elem::true_centroid() on a case by case
1022 : * basis.
1023 : */
1024 : virtual Point centroid () const;
1025 : #endif // LIBMESH_ENABLE_DEPRECATED
1026 :
1027 : /**
1028 : * \returns The "true" geometric centroid of the element, c=(cx, cy,
1029 : * cz), where:
1030 : *
1031 : * [cx] [\int x dV]
1032 : * [cy] := (1/V) * [\int y dV]
1033 : * [cz] [\int z dV]
1034 : *
1035 : * This method is virtual since some derived elements might want to
1036 : * use shortcuts to compute their centroid. For most element types,
1037 : * this method is more expensive than calling vertex_average(), so
1038 : * if you only need a point which is located "somewhere" in the
1039 : * interior of the element, consider calling vertex_average() instead.
1040 : */
1041 : virtual Point true_centroid () const;
1042 :
1043 : /**
1044 : * \returns A Point at the average of the elment's vertices.
1045 : *
1046 : * \note This used to be the base class centroid() implementation, but
1047 : * the centroid is only equal to the vertex average in some special cases.
1048 : * The centroid() implementation now returns the "true" centroid of the
1049 : * element (up to quadrature error).
1050 : */
1051 : Point vertex_average () const;
1052 :
1053 : /**
1054 : * \returns The "circumcenter of mass" (area-weighted average of
1055 : * triangulation circumcenters) of the element.
1056 : *
1057 : * Not implemented for infinite elements, not currently implemented
1058 : * for 3D elements, currently ignores curvature of element edges.
1059 : */
1060 0 : virtual Point quasicircumcenter () const
1061 0 : { libmesh_not_implemented(); }
1062 :
1063 : /**
1064 : * \returns The minimum vertex separation for the element.
1065 : */
1066 : virtual Real hmin () const;
1067 :
1068 : /**
1069 : * \returns The maximum vertex separation for the element.
1070 : */
1071 : virtual Real hmax () const;
1072 :
1073 : /**
1074 : * \returns The (length/area/volume) of the geometric element.
1075 : *
1076 : * If the element is twisted or inverted such that the mapping
1077 : * Jacobian is singular at any point, implementations of this method
1078 : * may return a "net" volume or may simply return NaN.
1079 : */
1080 : virtual Real volume () const;
1081 :
1082 : /**
1083 : * \returns A bounding box (not necessarily the minimal bounding box)
1084 : * containing the geometric element.
1085 : *
1086 : * The base class implementation determines a bounding box for the
1087 : * element *nodes*, which should be sufficient for first order
1088 : * finite elements. Higher order geometric elements will need to
1089 : * override with an implementation which takes curved elements into
1090 : * account.
1091 : */
1092 : virtual BoundingBox loose_bounding_box () const;
1093 :
1094 : /**
1095 : * \returns A quantitative assessment of element quality based on
1096 : * the quality metric \p q specified by the user. Not all ElemQuality
1097 : * metrics are supported for all Elem types; consult the Elem::quality()
1098 : * overrides for specific Elem types to determine which quality metrics
1099 : * are supported. The ElemQuality metrics with generic support for all
1100 : * Elems with dimension > 1 are:
1101 : * .) EDGE_LENGTH_RATIO - ratio of maximum to minimum edge (in 2D,
1102 : * side) length, where the min/max is taken over all Elem edges.
1103 : * .) MIN,MAX_ANGLE - The minimum (respectively maximum) angle
1104 : * between all pairs of adjacent Elem edges, in degrees. Note
1105 : * that, in 3D, these are *not* the dihedral angles (angle
1106 : * between adjacent planar faces of the element), which we plan to
1107 : * add support for in the future. In 2D, we compute the angle
1108 : * between adjacent sides for this metric.
1109 : */
1110 : virtual Real quality (const ElemQuality q) const;
1111 :
1112 : /**
1113 : * \returns The suggested quality bounds for the Elem based on
1114 : * quality measure \p q.
1115 : *
1116 : * These are the values suggested by the CUBIT User's Manual. Since
1117 : * this function can have no possible meaning for an abstract Elem,
1118 : * it is an error in the base class.
1119 : */
1120 0 : virtual std::pair<Real,Real> qual_bounds (const ElemQuality) const
1121 0 : { libmesh_not_implemented(); return std::make_pair(0.,0.); }
1122 :
1123 : /**
1124 : * \returns \p true if the physical point p is contained in this
1125 : * element, false otherwise.
1126 : *
1127 : * For linear elements, performs an initial tight bounding box check
1128 : * (as an optimization step) and (if that passes) then uses the
1129 : * user-defined tolerance "tol" in a call to inverse_map() to actually
1130 : * test if the point is in the element. For quadratic elements, the
1131 : * bounding box optimization is skipped, and only the inverse_map()
1132 : * steps are performed.
1133 : *
1134 : * \note This routine should not be used to determine if a point
1135 : * is merely "nearby" an element to within some tolerance. For that,
1136 : * use Elem::close_to_point() instead.
1137 : */
1138 : virtual bool contains_point (const Point & p, Real tol=TOLERANCE) const;
1139 :
1140 : /**
1141 : * \returns \p true if the master-space point p is contained in the
1142 : * reference element corresponding to this element, false otherwise.
1143 : *
1144 : * Since we are doing floating point comparisons here the parameter
1145 : * \p eps can be specified to indicate a tolerance. For example,
1146 : * \f$ x \le 1 \f$ becomes \f$ x \le 1 + \epsilon \f$.
1147 : */
1148 : virtual bool on_reference_element(const Point & p,
1149 : const Real eps = TOLERANCE) const = 0;
1150 :
1151 : /**
1152 : * \returns \p true if this element is "close" to the point p, where
1153 : * "close" is determined by the tolerance tol.
1154 : */
1155 : virtual bool close_to_point(const Point & p, Real tol) const;
1156 :
1157 : /**
1158 : * \returns \p true if edge \p i is positively oriented. An edge is
1159 : * positively oriented iff its first vertex (i.e. zeroth node) is
1160 : * lexicographically greater than its second vertex (i.e. first node).
1161 : */
1162 : bool positive_edge_orientation(const unsigned int i) const;
1163 :
1164 : /**
1165 : * \returns \p true if face \p i is positively oriented. A face is
1166 : * positively oriented iff the triangle defined by the lexicographically
1167 : * least vertex and its two adjacent vertices on the same face is
1168 : * positively oriented. Said triangle is positively oriented iff its
1169 : * vertices are an odd permutation of their lexicographic ordering.
1170 : */
1171 : bool positive_face_orientation(const unsigned int i) const;
1172 :
1173 :
1174 : /**
1175 : * A helper function for copying generic element data (mapping,
1176 : * subdomain, processor) from an element to a derived (child, side,
1177 : * edge) element. Useful for forwards compatibility when new data
1178 : * is added.
1179 : */
1180 : void inherit_data_from(const Elem & src);
1181 :
1182 :
1183 : private:
1184 : /**
1185 : * Shared private implementation used by the contains_point()
1186 : * and close_to_point() routines. The box_tol tolerance is
1187 : * used in the bounding box optimization, the map_tol tolerance is used
1188 : * in the calls to inverse_map() and on_reference_element().
1189 : */
1190 : bool point_test(const Point & p, Real box_tol, Real map_tol) const;
1191 :
1192 : public:
1193 : /**
1194 : * \returns \p true if the element map is definitely affine (i.e. the same at
1195 : * every quadrature point) within numerical tolerances.
1196 : */
1197 0 : virtual bool has_affine_map () const { return false; }
1198 :
1199 : /**
1200 : * \returns \p true if the element map is invertible everywhere on
1201 : * the element, to within a user-specified tolerance. The tolerance
1202 : * is generally used in comparisons against zero, so it should be an
1203 : * absolute rather than a relative tolerance. Throws a
1204 : * libmesh_not_implemented() error unless specialized by derived
1205 : * classes.
1206 : */
1207 : virtual bool has_invertible_map(Real tol = TOLERANCE*TOLERANCE) const;
1208 :
1209 : /**
1210 : * \returns \p true if the Lagrange shape functions on this element
1211 : * are linear.
1212 : */
1213 0 : virtual bool is_linear () const { return false; }
1214 :
1215 : /**
1216 : * Prints relevant information about the element.
1217 : */
1218 : void print_info (std::ostream & os=libMesh::out) const;
1219 :
1220 : /**
1221 : * Prints relevant information about the element to a string.
1222 : */
1223 : std::string get_info () const;
1224 :
1225 : /**
1226 : * \returns \p true if the element is active (i.e. has no active
1227 : * descendants) or AMR is disabled, \p false otherwise.
1228 : *
1229 : * \note It suffices to check the first child only.
1230 : */
1231 : bool active () const;
1232 :
1233 : /**
1234 : * \returns \p true if the element is an ancestor (i.e. has an
1235 : * active child or ancestor child), \p false otherwise or when AMR
1236 : * is disabled.
1237 : */
1238 : bool ancestor () const;
1239 :
1240 : /**
1241 : * \returns \p true if the element is subactive (i.e. has no active
1242 : * descendants), \p false otherwise or if AMR is disabled.
1243 : */
1244 : bool subactive () const;
1245 :
1246 : /**
1247 : * \returns \p true if the element has any children (active or not),
1248 : * \p false otherwise, or if AMR is disabled.
1249 : */
1250 : bool has_children () const;
1251 :
1252 : /**
1253 : * \returns \p true if the element has any descendants other than
1254 : * its immediate children, \p false otherwise, or if AMR is disabled.
1255 : */
1256 : bool has_ancestor_children () const;
1257 :
1258 : /**
1259 : * \returns \p true if \p descendant is a child of \p this, or a
1260 : * child of a child of \p this, etc., \p false otherwise or if AMR
1261 : * is disabled.
1262 : */
1263 : bool is_ancestor_of(const Elem * descendant) const;
1264 :
1265 : /**
1266 : * \returns A const pointer to the element's parent, or \p nullptr if
1267 : * the element was not created via refinement.
1268 : */
1269 : const Elem * parent () const;
1270 :
1271 : /**
1272 : * \returns A pointer to the element's parent, or \p nullptr if
1273 : * the element was not created via refinement.
1274 : */
1275 : Elem * parent ();
1276 :
1277 : /**
1278 : * Sets the pointer to the element's parent.
1279 : * Dangerous! Only use this if you know what you are doing!
1280 : */
1281 : void set_parent (Elem * p);
1282 :
1283 : /**
1284 : * \returns A pointer to the element's top-most (i.e. level-0) parent.
1285 : *
1286 : * That is, \p this if this is a level-0 element, this element's parent
1287 : * if this is a level-1 element, this element's grandparent if this is
1288 : * a level-2 element, etc...
1289 : */
1290 : const Elem * top_parent () const;
1291 :
1292 : /**
1293 : * \returns The higher-dimensional Elem for which this Elem is a face.
1294 : *
1295 : * In some cases it is desirable to extract the boundary (or a subset thereof)
1296 : * of a D-dimensional mesh as a (D-1)-dimensional manifold. In this case
1297 : * we may want to know the 'parent' element from which the manifold elements
1298 : * were extracted. We can easily do that for the level-0 manifold elements
1299 : * by storing the D-dimensional parent. This method provides access to that
1300 : * element.
1301 : *
1302 : * This method returns nullptr if this->dim() == LIBMESH_DIM; in
1303 : * such cases no data storage for an interior parent pointer has
1304 : * been allocated.
1305 : */
1306 : const Elem * interior_parent () const;
1307 :
1308 : Elem * interior_parent ();
1309 :
1310 : /**
1311 : * Sets the pointer to the element's interior_parent.
1312 : * Dangerous! Only use this if you know what you are doing!
1313 : */
1314 : void set_interior_parent (Elem * p);
1315 :
1316 : /**
1317 : * \returns The distance between nodes n1 and n2.
1318 : *
1319 : * Useful for computing the lengths of the sides of elements.
1320 : */
1321 : Real length (const unsigned int n1,
1322 : const unsigned int n2) const;
1323 :
1324 : /**
1325 : * \returns The number of adjacent vertices that uniquely define the
1326 : * location of the \f$ n^{th} \f$ second-order node, or 0 for linear
1327 : * elements.
1328 : *
1329 : * This method is useful when converting linear elements to quadratic
1330 : * elements.
1331 : *
1332 : * \note \p n has to be greater than or equal to \p this->n_vertices().
1333 : */
1334 : virtual unsigned int n_second_order_adjacent_vertices (const unsigned int n) const;
1335 :
1336 : /**
1337 : * \returns The element-local number of the \f$ v^{th} \f$ vertex
1338 : * that defines the \f$ n^{th} \f$ second-order node, or 0 for
1339 : * linear elements.
1340 : *
1341 : * \note The value is always less than \p this->n_vertices(), while
1342 : * \p n has to be greater than or equal to \p this->n_vertices().
1343 : */
1344 : virtual unsigned short int second_order_adjacent_vertex (const unsigned int n,
1345 : const unsigned int v) const;
1346 :
1347 : /**
1348 : * \returns A pair (c,v), where
1349 : * c == child index, and
1350 : * v == element-local index of the \p \f$ n^{th} \f$
1351 : * second-order node on the parent element.
1352 : * For linear elements, (0,0) is returned.
1353 : *
1354 : * \note The return values are always less than \p this->n_children()
1355 : * and \p this->child_ptr(c)->n_vertices().
1356 : *
1357 : * \note \p n has to be greater than or equal to \p this->n_vertices().
1358 : *
1359 : * \note On refined second-order elements, the return value will
1360 : * satisfy \p this->node_ptr(n) == this->child_ptr(c)->node_ptr(v).
1361 : */
1362 : virtual std::pair<unsigned short int, unsigned short int>
1363 : second_order_child_vertex (const unsigned int n) const;
1364 :
1365 : /**
1366 : * \returns The ElemType of the associated second-order element
1367 : * (which will be the same as the input if the input is already a
1368 : * second-order ElemType) or INVALID_ELEM for elements that cannot be
1369 : * converted into higher order equivalents.
1370 : *
1371 : * For example, when \p this is a \p TET4, then \p TET10 is returned.
1372 : *
1373 : * For some elements, there exist two second-order equivalents, e.g.
1374 : * for \p Quad4 there is \p Quad8 and \p Quad9. When the optional
1375 : * \p full_ordered is \p true, then \p QUAD9 is returned. When
1376 : * \p full_ordered is \p false, then \p QUAD8 is returned.
1377 : */
1378 : static ElemType second_order_equivalent_type (const ElemType et,
1379 : const bool full_ordered=true);
1380 :
1381 : /**
1382 : * \returns The element type of the associated first-order element,
1383 : * or \p INVALID_ELEM for first-order or other elements that cannot be
1384 : * converted into lower order equivalents.
1385 : *
1386 : * For example, when \p this is a \p TET10, then \p TET4 is returned.
1387 : */
1388 : static ElemType first_order_equivalent_type (const ElemType et);
1389 :
1390 : /**
1391 : * \returns The ElemType of the associated "complete" order element
1392 : * (which will be the same as the input if the input is already a
1393 : * complete-order ElemType), or INVALID_ELEM for elements that cannot be
1394 : * converted into complete-order equivalents.
1395 : *
1396 : * The "complete" version of an element is an element which can
1397 : * represent the same geometry but which has nodes available to
1398 : * restore degrees of freedom on any vertex, edge, or face.
1399 : *
1400 : * For example, when \p this is a \p TET4, then \p TET14 is returned.
1401 : */
1402 : static ElemType complete_order_equivalent_type (const ElemType et);
1403 :
1404 : /**
1405 : * \returns The refinement level of the current element.
1406 : *
1407 : * If the element's parent is \p nullptr then by convention it is at
1408 : * level 0, otherwise it is simply at one level greater than its
1409 : * parent.
1410 : */
1411 : unsigned int level () const;
1412 :
1413 : /**
1414 : * \returns The value of the p refinement level of an active
1415 : * element, or the minimum value of the p refinement levels
1416 : * of an ancestor element's descendants.
1417 : */
1418 : unsigned int p_level () const;
1419 :
1420 : /**
1421 : * \returns \p true if the specified child is on the specified side.
1422 : */
1423 : virtual bool is_child_on_side(const unsigned int c,
1424 : const unsigned int s) const = 0;
1425 :
1426 : /**
1427 : * \returns The value of the mapping type for the element.
1428 : */
1429 : ElemMappingType mapping_type () const;
1430 :
1431 : /**
1432 : * Sets the value of the mapping type for the element.
1433 : */
1434 : void set_mapping_type (const ElemMappingType type);
1435 :
1436 : /**
1437 : * \returns The value of the mapping data for the element.
1438 : */
1439 : unsigned char mapping_data () const;
1440 :
1441 : /**
1442 : * Sets the value of the mapping data for the element.
1443 : */
1444 : void set_mapping_data (const unsigned char data);
1445 :
1446 :
1447 : #ifdef LIBMESH_ENABLE_AMR
1448 :
1449 : /**
1450 : * Enumeration of possible element refinement states.
1451 : */
1452 : enum RefinementState { COARSEN = 0,
1453 : DO_NOTHING,
1454 : REFINE,
1455 : JUST_REFINED,
1456 : JUST_COARSENED,
1457 : INACTIVE,
1458 : COARSEN_INACTIVE,
1459 : INVALID_REFINEMENTSTATE };
1460 :
1461 : /**
1462 : * \returns A constant pointer to the \f$ i^{th} \f$ child for this element.
1463 : * For internal use only - skips assertions about null pointers.
1464 : */
1465 : const Elem * raw_child_ptr (unsigned int i) const;
1466 :
1467 : /**
1468 : * \returns A constant pointer to the \f$ i^{th} \f$ child for this element.
1469 : * Do not call if this element has no children, i.e. is active.
1470 : */
1471 : const Elem * child_ptr (unsigned int i) const;
1472 :
1473 : /**
1474 : * \returns A non-constant pointer to the \f$ i^{th} \f$ child for this element.
1475 : * Do not call if this element has no children, i.e. is active.
1476 : */
1477 : Elem * child_ptr (unsigned int i);
1478 :
1479 : /**
1480 : * Nested classes for use iterating over all children of a parent
1481 : * element.
1482 : */
1483 : class ChildRefIter;
1484 : class ConstChildRefIter;
1485 :
1486 : /**
1487 : * Returns a range with all children of a parent element, usable in
1488 : * range-based for loops. The exact type of the return value here
1489 : * may be subject to change in future libMesh releases, but the
1490 : * iterators will always dereference to produce a reference to a
1491 : * child element.
1492 : */
1493 : SimpleRange<ChildRefIter> child_ref_range();
1494 :
1495 : SimpleRange<ConstChildRefIter> child_ref_range() const;
1496 :
1497 : private:
1498 : /**
1499 : * Sets the pointer to the \f$ i^{th} \f$ child for this element.
1500 : * Do not call if this element has no children, i.e. is active.
1501 : */
1502 : void set_child (unsigned int c, Elem * elem);
1503 :
1504 : public:
1505 : /**
1506 : * \returns The child index which \p e corresponds to.
1507 : *
1508 : * I.e. if c = a->which_child_am_i(e); then a->child_ptr(c) will be
1509 : * e.
1510 : */
1511 : unsigned int which_child_am_i(const Elem * e) const;
1512 :
1513 : /**
1514 : * \returns \p true if the specified child is on the specified edge.
1515 : */
1516 : virtual bool is_child_on_edge(const unsigned int c,
1517 : const unsigned int e) const;
1518 :
1519 : /**
1520 : * Adds a child pointer to the array of children of this element.
1521 : * If this is the first child to be added, this method allocates
1522 : * memory in the parent's _children array, otherwise, it just sets
1523 : * the pointer.
1524 : */
1525 : void add_child (Elem * elem);
1526 :
1527 : /**
1528 : * Adds a new child pointer to the specified index in the array of
1529 : * children of this element. If this is the first child to be added,
1530 : * this method allocates memory in the parent's _children array,
1531 : * otherwise, it just sets the pointer.
1532 : */
1533 : void add_child (Elem * elem, unsigned int c);
1534 :
1535 : /**
1536 : * Replaces the child pointer at the specified index in the child array.
1537 : */
1538 : void replace_child (Elem * elem, unsigned int c);
1539 :
1540 : /**
1541 : * Fills the vector \p family with the children of this element,
1542 : * recursively. Calling this method on a twice-refined element
1543 : * will give you the element itself, its direct children, and their
1544 : * children, etc... When the optional parameter \p reset is
1545 : * true, the vector will be cleared before the element and its
1546 : * descendants are added.
1547 : *
1548 : * The family tree only includes ancestor and active elements. To
1549 : * include subactive elements as well, use total_family_tree().
1550 : */
1551 : void family_tree (std::vector<const Elem *> & family,
1552 : bool reset = true) const;
1553 :
1554 : /**
1555 : * Non-const version of function above; fills a vector of non-const pointers.
1556 : */
1557 : void family_tree (std::vector<Elem *> & family,
1558 : bool reset = true);
1559 :
1560 : /**
1561 : * Same as the \p family_tree() member, but also adds any subactive
1562 : * descendants.
1563 : */
1564 : void total_family_tree (std::vector<const Elem *> & family,
1565 : bool reset = true) const;
1566 :
1567 : /**
1568 : * Non-const version of function above; fills a vector of non-const pointers.
1569 : */
1570 : void total_family_tree (std::vector<Elem *> & family,
1571 : bool reset = true);
1572 :
1573 : /**
1574 : * Same as the \p family_tree() member, but only adds the active
1575 : * children. Can be thought of as removing all the inactive
1576 : * elements from the vector created by \p family_tree, but is
1577 : * implemented more efficiently.
1578 : */
1579 : void active_family_tree (std::vector<const Elem *> & active_family,
1580 : bool reset = true) const;
1581 :
1582 : /**
1583 : * Non-const version of function above; fills a vector of non-const pointers.
1584 : */
1585 : void active_family_tree (std::vector<Elem *> & active_family,
1586 : bool reset = true);
1587 :
1588 : /**
1589 : * Same as the \p family_tree() member, but only adds elements
1590 : * which are next to \p side.
1591 : */
1592 : void family_tree_by_side (std::vector<const Elem *> & family,
1593 : unsigned int side,
1594 : bool reset = true) const;
1595 :
1596 : /**
1597 : * Non-const version of function above; fills a vector of non-const pointers.
1598 : */
1599 : void family_tree_by_side (std::vector<Elem *> & family,
1600 : unsigned int side,
1601 : bool reset = true);
1602 :
1603 : /**
1604 : * Same as the \p active_family_tree() member, but only adds elements
1605 : * which are next to \p side.
1606 : */
1607 : void active_family_tree_by_side (std::vector<const Elem *> & family,
1608 : unsigned int side,
1609 : bool reset = true) const;
1610 :
1611 : /**
1612 : * Non-const version of function above; fills a vector of non-const pointers.
1613 : */
1614 : void active_family_tree_by_side (std::vector<Elem *> & family,
1615 : unsigned int side,
1616 : bool reset = true);
1617 :
1618 : /**
1619 : * Same as the \p family_tree() member, but only adds elements
1620 : * which are next to \p neighbor.
1621 : */
1622 : void family_tree_by_neighbor (std::vector<const Elem *> & family,
1623 : const Elem * neighbor,
1624 : bool reset = true) const;
1625 :
1626 : /**
1627 : * Non-const version of function above; fills a vector of non-const pointers.
1628 : */
1629 : void family_tree_by_neighbor (std::vector<Elem *> & family,
1630 : Elem * neighbor,
1631 : bool reset = true);
1632 :
1633 : /**
1634 : * Same as the \p family_tree_by_neighbor() member, but also adds
1635 : * any subactive descendants.
1636 : */
1637 : void total_family_tree_by_neighbor (std::vector<const Elem *> & family,
1638 : const Elem * neighbor,
1639 : bool reset = true) const;
1640 :
1641 : /**
1642 : * Non-const version of function above; fills a vector of non-const pointers.
1643 : */
1644 : void total_family_tree_by_neighbor (std::vector<Elem *> & family,
1645 : Elem * neighbor,
1646 : bool reset = true);
1647 :
1648 : /**
1649 : * Same as the \p family_tree() member, but only adds elements
1650 : * which are next to \p subneighbor. Only applicable when
1651 : * \p this->has_neighbor(neighbor) and
1652 : * \p neighbor->is_ancestor(subneighbor)
1653 : */
1654 : void family_tree_by_subneighbor (std::vector<const Elem *> & family,
1655 : const Elem * neighbor,
1656 : const Elem * subneighbor,
1657 : bool reset = true) const;
1658 :
1659 : /**
1660 : * Non-const version of function above; fills a vector of non-const pointers.
1661 : */
1662 : void family_tree_by_subneighbor (std::vector<Elem *> & family,
1663 : Elem * neighbor,
1664 : Elem * subneighbor,
1665 : bool reset = true);
1666 :
1667 : /**
1668 : * Same as the \p family_tree_by_subneighbor() member, but also adds
1669 : * any subactive descendants.
1670 : */
1671 : void total_family_tree_by_subneighbor (std::vector<const Elem *> & family,
1672 : const Elem * neighbor,
1673 : const Elem * subneighbor,
1674 : bool reset = true) const;
1675 :
1676 : /**
1677 : * Non-const version of function above; fills a vector of non-const pointers.
1678 : */
1679 : void total_family_tree_by_subneighbor (std::vector<Elem *> & family,
1680 : Elem * neighbor,
1681 : Elem * subneighbor,
1682 : bool reset = true);
1683 :
1684 : /**
1685 : * Same as the \p active_family_tree() member, but only adds elements
1686 : * which are next to \p neighbor.
1687 : */
1688 : void active_family_tree_by_neighbor (std::vector<const Elem *> & family,
1689 : const Elem * neighbor,
1690 : bool reset = true) const;
1691 :
1692 : /**
1693 : * Non-const version of function above; fills a vector of non-const pointers.
1694 : */
1695 : void active_family_tree_by_neighbor (std::vector<Elem *> & family,
1696 : Elem * neighbor,
1697 : bool reset = true);
1698 :
1699 : /**
1700 : * Same as the \p active_family_tree_by_neighbor() member, but the
1701 : * \p neighbor here may be a topological (e.g. periodic boundary
1702 : * condition) neighbor, not just a local neighbor.
1703 : */
1704 : void active_family_tree_by_topological_neighbor (std::vector<const Elem *> & family,
1705 : const Elem * neighbor,
1706 : const MeshBase & mesh,
1707 : const PointLocatorBase & point_locator,
1708 : const PeriodicBoundaries * pb,
1709 : bool reset = true) const;
1710 :
1711 : /**
1712 : * Non-const version of function above; fills a vector of non-const pointers.
1713 : */
1714 : void active_family_tree_by_topological_neighbor (std::vector<Elem *> & family,
1715 : Elem * neighbor,
1716 : const MeshBase & mesh,
1717 : const PointLocatorBase & point_locator,
1718 : const PeriodicBoundaries * pb,
1719 : bool reset = true);
1720 :
1721 : /**
1722 : * \returns The value of the refinement flag for the element.
1723 : */
1724 : RefinementState refinement_flag () const;
1725 :
1726 : /**
1727 : * Sets the value of the refinement flag for the element.
1728 : */
1729 : void set_refinement_flag (const RefinementState rflag);
1730 :
1731 : /**
1732 : * \returns The value of the p-refinement flag for the element.
1733 : */
1734 : RefinementState p_refinement_flag () const;
1735 :
1736 : /**
1737 : * Sets the value of the p-refinement flag for the element.
1738 : */
1739 : void set_p_refinement_flag (const RefinementState pflag);
1740 :
1741 : /**
1742 : * \returns The maximum value of the p-refinement levels of
1743 : * an ancestor element's descendants.
1744 : */
1745 : unsigned int max_descendant_p_level () const;
1746 :
1747 : /**
1748 : * \returns The minimum p-refinement level of elements which are
1749 : * descended from this element, and which share a side with the
1750 : * active \p neighbor.
1751 : */
1752 : unsigned int min_p_level_by_neighbor (const Elem * neighbor,
1753 : unsigned int current_min) const;
1754 :
1755 : /**
1756 : * \returns The minimum new p-refinement level (i.e. after refinement
1757 : * and coarsening is done) of elements which are descended from this
1758 : * element and which share a side with the active \p neighbor.
1759 : */
1760 : unsigned int min_new_p_level_by_neighbor (const Elem * neighbor,
1761 : unsigned int current_min) const;
1762 :
1763 : /**
1764 : * Sets the value of the p-refinement level for the element.
1765 : *
1766 : * \note The maximum p-refinement level is currently 255.
1767 : */
1768 : void set_p_level (const unsigned int p);
1769 :
1770 : /**
1771 : * Sets the value of the p-refinement level for the element
1772 : * without altering the p-level of its ancestors
1773 : */
1774 : void hack_p_level (const unsigned int p);
1775 :
1776 : /**
1777 : * Sets the value of the p-refinement level for the element
1778 : * without altering the p-level of its ancestors; also sets the
1779 : * p_refinement_flag, simultaneously so that they can be safely
1780 : * checked for mutual consistency
1781 : */
1782 : void hack_p_level_and_refinement_flag (const unsigned int p,
1783 : RefinementState pflag);
1784 :
1785 : /**
1786 : * Refine the element.
1787 : */
1788 : virtual void refine (MeshRefinement & mesh_refinement);
1789 :
1790 : /**
1791 : * Coarsen the element. This function is non-virtual since it is the same
1792 : * for all element types.
1793 : */
1794 : void coarsen ();
1795 :
1796 : /**
1797 : * Contract an active element, i.e. remove pointers to any
1798 : * subactive children. This should only be called via
1799 : * MeshRefinement::contract, which will also remove subactive
1800 : * children from the mesh.
1801 : */
1802 : void contract ();
1803 :
1804 : #endif
1805 :
1806 : #ifndef NDEBUG
1807 : /**
1808 : * Checks for consistent neighbor links on this element.
1809 : */
1810 : void libmesh_assert_valid_neighbors() const;
1811 :
1812 : /**
1813 : * Checks for a valid id and pointers to nodes with valid ids on
1814 : * this element.
1815 : */
1816 : void libmesh_assert_valid_node_pointers() const;
1817 : #endif // !NDEBUG
1818 :
1819 : /**
1820 : * \returns The local node index of the given point IF said node
1821 : * has a singular Jacobian for this element. If the given point
1822 : * is not a node or is a node and does not have a singular Jacobian,
1823 : * this will return invalid_uint.
1824 : *
1825 : * The intention is for this to be overridden in derived element
1826 : * classes that do have nodes that have singular Jacobians. When
1827 : * mapping failures are caught, we can check this to see if the
1828 : * failed physical point is actually a singular point and
1829 : * return the correct master point.
1830 : */
1831 0 : virtual unsigned int local_singular_node(const Point & /* p */, const Real /* tol */ = TOLERANCE*TOLERANCE) const
1832 0 : { return invalid_uint; }
1833 :
1834 : /**
1835 : * \returns true iff the node at the given index has a singular
1836 : * mapping; i.e. is the degree-4 node on a Pyramid.
1837 : */
1838 0 : virtual bool is_singular_node(unsigned int /* node_i */) const { return false; }
1839 :
1840 : /**
1841 : * \returns The local index of the center node on the side \p side.
1842 : *
1843 : * A center node is a node that is located at the centroid of the given side.
1844 : * If the given side does not have a center node, this will return invalid_uint.
1845 : */
1846 : virtual unsigned int center_node_on_side(const unsigned short side) const;
1847 :
1848 : protected:
1849 :
1850 : /**
1851 : * The protected nested SideIter class is used to iterate over the
1852 : * sides of this Elem. It is a specially-designed class since
1853 : * no sides are actually stored by the element. This iterator-like
1854 : * class has to provide the following three operations
1855 : * 1) operator*
1856 : * 2) operator++
1857 : * 3) operator==
1858 : * The definition can be found at the end of this header file.
1859 : */
1860 : class SideIter;
1861 :
1862 : public:
1863 : /**
1864 : * Useful iterator typedefs
1865 : */
1866 : typedef Predicates::multi_predicate Predicate;
1867 :
1868 : /**
1869 : * Data structure for iterating over sides. Defined at the end of
1870 : * this header file.
1871 : */
1872 : struct side_iterator;
1873 :
1874 : /**
1875 : * Iterator accessor functions
1876 : */
1877 : side_iterator boundary_sides_begin();
1878 : side_iterator boundary_sides_end();
1879 :
1880 : private:
1881 : /**
1882 : * Side iterator helper functions. Used to replace the begin()
1883 : * and end() functions of the STL containers.
1884 : */
1885 : SideIter _first_side();
1886 : SideIter _last_side();
1887 :
1888 : public:
1889 :
1890 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1891 :
1892 : /**
1893 : * \returns \p true if the element is an infinite element,
1894 : * \p false otherwise.
1895 : */
1896 : virtual bool infinite () const = 0;
1897 :
1898 : /**
1899 : * \returns \p true if the specified (local) node number is a
1900 : * "mid-edge" node on an infinite element edge.
1901 : *
1902 : * This is false for all nodes on non-infinite elements, so we won't
1903 : * make it pure virtual, to simplify their code.
1904 : */
1905 0 : virtual bool is_mid_infinite_edge_node(const unsigned int /* n */) const
1906 0 : { libmesh_assert (!this->infinite()); return false; }
1907 :
1908 : /**
1909 : * \returns The origin for an infinite element.
1910 : *
1911 : * Currently, all infinite elements used in a mesh share the same
1912 : * origin. Override this in infinite element classes.
1913 : */
1914 0 : virtual Point origin () const { libmesh_not_implemented(); return Point(); }
1915 :
1916 : #else
1917 :
1918 : static constexpr bool infinite () { return false; }
1919 :
1920 : #endif
1921 :
1922 : /**
1923 : * \returns An Elem of type \p type wrapped in a smart pointer.
1924 : */
1925 : static std::unique_ptr<Elem> build (const ElemType type,
1926 : Elem * p=nullptr);
1927 :
1928 : /**
1929 : * Calls the build() method above with a nullptr parent, and
1930 : * additionally sets the newly-created Elem's id. This can be useful
1931 : * when adding pre-numbered Elems to a Mesh via add_elem() calls.
1932 : */
1933 : static std::unique_ptr<Elem> build_with_id (const ElemType type,
1934 : dof_id_type id);
1935 :
1936 : /**
1937 : * \returns An Elem of the same type as \p this, wrapped in a smart
1938 : * pointer.
1939 : *
1940 : * This is not a complete clone() method (since e.g. it does not set
1941 : * node pointers; the standard use case reassigns node pointers from
1942 : * a different mesh), but it is necessary to use this instead of
1943 : * build() for runtime-polymorphic elements like Polygon subtypes
1944 : * whose "type" depends on more than their type(), and it is useful
1945 : * to use this for elements whose id, unique_id, extra integers,
1946 : * etc. should be preserved in the near-clone.
1947 : */
1948 : virtual std::unique_ptr<Elem> disconnected_clone () const;
1949 :
1950 : /**
1951 : * Returns the number of independent permutations of element nodes -
1952 : * e.g. a cube can be reoriented to put side 0 where side N is (for
1953 : * 0 <= N < 6) and then rotated in one of four ways, giving 24
1954 : * possible permutations.
1955 : *
1956 : * Permutations which change the mapping Jacobian of an element
1957 : * (i.e. flipping the element) are not allowed in this definition.
1958 : */
1959 : virtual unsigned int n_permutations() const = 0;
1960 :
1961 : /**
1962 : * Permutes the element (by swapping node and neighbor pointers)
1963 : * according to the specified index.
1964 : *
1965 : * This is useful for regression testing, by making it easy to make
1966 : * a structured mesh behave more like an arbitrarily unstructured
1967 : * mesh.
1968 : *
1969 : * This is so far *only* used for regression testing, so we do
1970 : * not currently provide a way to permute any boundary side/edge ids
1971 : * along with the element permutation.
1972 : */
1973 : virtual void permute(unsigned int perm_num) = 0;
1974 :
1975 : /**
1976 : * Flips the element (by swapping node and neighbor pointers) to
1977 : * have a mapping Jacobian of opposite sign.
1978 : *
1979 : * This is useful for automatically fixing up elements that have
1980 : * been newly created (e.g. from extrusions) with a negative
1981 : * Jacobian.
1982 : *
1983 : * If \p boundary_info is not null, swap boundary side/edge ids
1984 : * consistently.
1985 : */
1986 : virtual void flip(BoundaryInfo * boundary_info) = 0;
1987 :
1988 : /**
1989 : * \returns Whether the element is flipped compared to standard
1990 : * libMesh (e.g. clockwise for 2D elements) node orientations.
1991 : *
1992 : * Always returns \p false if a 2D element is not in the XY plane or
1993 : * a 1D element is not on the X axis; user code designed to work for
1994 : * embedded manifolds should handle any consistent orientation, and
1995 : * determining whether an orientation is consistent is not a local
1996 : * operation.
1997 : */
1998 : virtual bool is_flipped() const = 0;
1999 :
2000 : /**
2001 : * Flips the element (by swapping node and neighbor pointers) to
2002 : * have a mapping Jacobian of opposite sign, iff we find a negative
2003 : * orientation. This only fixes flipped elements; for tangled
2004 : * elements the only fixes possible are non-local.
2005 : */
2006 : void orient(BoundaryInfo * boundary_info);
2007 :
2008 : #ifdef LIBMESH_ENABLE_AMR
2009 :
2010 : /**
2011 : * \returns The local node id on the parent which corresponds to node
2012 : * \p n of child \p c, or \p invalid_uint if no such parent
2013 : * node exists.
2014 : */
2015 : virtual unsigned int as_parent_node (unsigned int c,
2016 : unsigned int n) const;
2017 :
2018 : /**
2019 : * \returns All the pairs of nodes (indexed by local node id) which
2020 : * should bracket node \p n of child \p c.
2021 : */
2022 : virtual
2023 : const std::vector<std::pair<unsigned char, unsigned char>> &
2024 : parent_bracketing_nodes(unsigned int c,
2025 : unsigned int n) const;
2026 :
2027 : /**
2028 : * \returns All the pairs of nodes (indexed by global node id) which
2029 : * should bracket node \p n of child \p c.
2030 : */
2031 : virtual
2032 : const std::vector<std::pair<dof_id_type, dof_id_type>>
2033 : bracketing_nodes(unsigned int c,
2034 : unsigned int n) const;
2035 :
2036 :
2037 : /**
2038 : * \returns The embedding matrix entry for the requested child.
2039 : */
2040 : virtual Real embedding_matrix (const unsigned int child_num,
2041 : const unsigned int child_node_num,
2042 : const unsigned int parent_node_num) const = 0;
2043 :
2044 : /**
2045 : * \returns A "version number" that identifies which embedding
2046 : * matrix is in use.
2047 : *
2048 : * Some element types may use a different embedding matrix depending
2049 : * on their geometric characteristics.
2050 : */
2051 0 : virtual unsigned int embedding_matrix_version () const { return 0; }
2052 :
2053 : #endif // LIBMESH_ENABLE_AMR
2054 :
2055 :
2056 : protected:
2057 :
2058 : /**
2059 : * Default tolerance to use in has_affine_map().
2060 : */
2061 : static constexpr Real affine_tol = TOLERANCE*TOLERANCE;
2062 :
2063 : /**
2064 : * \returns A hash key computed from a single node id.
2065 : */
2066 : static dof_id_type compute_key (dof_id_type n0);
2067 :
2068 : /**
2069 : * \returns A hash key computed from two node ids.
2070 : */
2071 : static dof_id_type compute_key (dof_id_type n0,
2072 : dof_id_type n1);
2073 :
2074 : /**
2075 : * \returns A hash key computed from three node ids.
2076 : */
2077 : static dof_id_type compute_key (dof_id_type n0,
2078 : dof_id_type n1,
2079 : dof_id_type n2);
2080 :
2081 : /**
2082 : * \returns A hash key computed from four node ids.
2083 : */
2084 : static dof_id_type compute_key (dof_id_type n0,
2085 : dof_id_type n1,
2086 : dof_id_type n2,
2087 : dof_id_type n3);
2088 :
2089 : /**
2090 : * Swaps two node_ptrs
2091 : */
2092 28794629 : void swap2nodes(unsigned int n1, unsigned int n2)
2093 : {
2094 4334448 : Node * temp = this->node_ptr(n1);
2095 30961853 : this->set_node(n1, this->node_ptr(n2));
2096 28794629 : this->set_node(n2, temp);
2097 28794629 : }
2098 :
2099 : /**
2100 : * Swaps two neighbor_ptrs
2101 : */
2102 6786404 : void swap2neighbors(unsigned int n1, unsigned int n2)
2103 : {
2104 793040 : Elem * temp = this->neighbor_ptr(n1);
2105 548184 : this->set_neighbor(n1, this->neighbor_ptr(n2));
2106 548184 : this->set_neighbor(n2, temp);
2107 6861068 : }
2108 :
2109 : /**
2110 : * Swaps two sides in \p boundary_info, if it is non-null.
2111 : */
2112 : void swap2boundarysides(unsigned short s1, unsigned short s2,
2113 : BoundaryInfo * boundary_info) const;
2114 :
2115 : /**
2116 : * Swaps two edges in \p boundary_info, if it is non-null.
2117 : */
2118 : void swap2boundaryedges(unsigned short e1, unsigned short e2,
2119 : BoundaryInfo * boundary_info) const;
2120 :
2121 : /**
2122 : * Swaps three node_ptrs, "rotating" them.
2123 : */
2124 11312324 : void swap3nodes(unsigned int n1, unsigned int n2, unsigned int n3)
2125 : {
2126 12189036 : swap2nodes(n1, n2);
2127 12189036 : swap2nodes(n2, n3);
2128 11312324 : }
2129 :
2130 : /**
2131 : * Swaps three neighbor_ptrs, "rotating" them.
2132 : */
2133 3071806 : void swap3neighbors(unsigned int n1, unsigned int n2,
2134 : unsigned int n3)
2135 : {
2136 3071806 : swap2neighbors(n1, n2);
2137 3071806 : swap2neighbors(n2, n3);
2138 3071806 : }
2139 :
2140 : /**
2141 : * Swaps four node_ptrs, "rotating" them.
2142 : */
2143 3278918 : void swap4nodes(unsigned int n1, unsigned int n2, unsigned int n3,
2144 : unsigned int n4)
2145 : {
2146 2952118 : swap3nodes(n1, n2, n3);
2147 3278918 : swap2nodes(n3, n4);
2148 3278918 : }
2149 :
2150 : /**
2151 : * Swaps four neighbor_ptrs, "rotating" them.
2152 : */
2153 618012 : void swap4neighbors(unsigned int n1, unsigned int n2,
2154 : unsigned int n3, unsigned int n4)
2155 : {
2156 618012 : swap3neighbors(n1, n2, n3);
2157 618012 : swap2neighbors(n3, n4);
2158 618012 : }
2159 :
2160 :
2161 : /**
2162 : * An implementation for simple (all sides equal) elements
2163 : */
2164 : template <typename Sideclass, typename Subclass>
2165 : std::unique_ptr<Elem>
2166 : simple_build_side_ptr(const unsigned int i);
2167 :
2168 : /**
2169 : * An implementation for simple (all sides equal) elements
2170 : */
2171 : template <typename Subclass>
2172 : void simple_build_side_ptr(std::unique_ptr<Elem> & side,
2173 : const unsigned int i,
2174 : ElemType sidetype);
2175 :
2176 : /**
2177 : * An implementation for simple (all sides equal) elements
2178 : */
2179 : template <typename Subclass, typename Mapclass>
2180 : void simple_side_ptr(std::unique_ptr<Elem> & side,
2181 : const unsigned int i,
2182 : ElemType sidetype);
2183 :
2184 : /**
2185 : * An implementation for simple (all edges equal) elements
2186 : */
2187 : template <typename Edgeclass, typename Subclass>
2188 : std::unique_ptr<Elem>
2189 : simple_build_edge_ptr(const unsigned int i);
2190 :
2191 : /**
2192 : * An implementation for simple (all edges equal) elements
2193 : */
2194 : template <typename Subclass>
2195 : void simple_build_edge_ptr(std::unique_ptr<Elem> & edge,
2196 : const unsigned int i,
2197 : ElemType edgetype);
2198 :
2199 :
2200 : #ifdef LIBMESH_ENABLE_AMR
2201 :
2202 : /**
2203 : * Elem subclasses which don't do their own bracketing node
2204 : * calculations will need to supply a static cache, since the
2205 : * default calculation is slow.
2206 : */
2207 : virtual
2208 : std::vector<std::vector<std::vector<std::vector<std::pair<unsigned char, unsigned char>>>>> &
2209 0 : _get_bracketing_node_cache() const
2210 : {
2211 0 : static std::vector<std::vector<std::vector<std::vector<std::pair<unsigned char, unsigned char>>>>> c;
2212 0 : libmesh_error();
2213 : return c;
2214 : }
2215 :
2216 : /**
2217 : * Elem subclasses which don't do their own child-to-parent node
2218 : * calculations will need to supply a static cache, since the
2219 : * default calculation is slow.
2220 : */
2221 : virtual
2222 : std::vector<std::vector<std::vector<signed char>>> &
2223 0 : _get_parent_indices_cache() const
2224 : {
2225 0 : static std::vector<std::vector<std::vector<signed char>>> c;
2226 0 : libmesh_error();
2227 : return c;
2228 : }
2229 :
2230 : #endif // LIBMESH_ENABLE_AMR
2231 :
2232 : public:
2233 :
2234 : /**
2235 : * Replaces this element with \p nullptr for all of its neighbors.
2236 : * This is useful when deleting an element.
2237 : */
2238 : void nullify_neighbors ();
2239 :
2240 : protected:
2241 :
2242 : /**
2243 : * Pointers to the nodes we are connected to.
2244 : */
2245 : Node ** _nodes;
2246 :
2247 : /**
2248 : * Pointers to this element's parent and neighbors, and for
2249 : * lower-dimensional elements' interior_parent.
2250 : */
2251 : Elem ** _elemlinks;
2252 :
2253 : #ifdef LIBMESH_ENABLE_AMR
2254 : /**
2255 : * unique_ptr to array of this element's children.
2256 : *
2257 : * A Mesh ultimately owns the child Elems so we are not responsible
2258 : * for deleting them, but we are responsible for cleaning up the
2259 : * array allocated to hold those Elems, hence the unique_ptr.
2260 : */
2261 : std::unique_ptr<Elem *[]> _children;
2262 : #endif
2263 :
2264 : /**
2265 : * The subdomain to which this element belongs.
2266 : */
2267 : subdomain_id_type _sbd_id;
2268 :
2269 : #ifdef LIBMESH_ENABLE_AMR
2270 : /**
2271 : * h refinement flag. This is stored as an unsigned char
2272 : * to save space.
2273 : */
2274 : unsigned char _rflag;
2275 :
2276 : /**
2277 : * p refinement flag. This is stored as an unsigned char
2278 : * to save space.
2279 : */
2280 : unsigned char _pflag;
2281 :
2282 : /**
2283 : * p refinement level - the difference between the
2284 : * polynomial degree on this element and the minimum
2285 : * polynomial degree on the mesh.
2286 : * This is stored as an unsigned char to save space.
2287 : * In theory, these last four bytes might have
2288 : * been padding anyway.
2289 : */
2290 : unsigned char _p_level;
2291 : #endif
2292 :
2293 : /**
2294 : * Mapping function type; currently either 0 (LAGRANGE) or 1
2295 : * (RATIONAL_BERNSTEIN).
2296 : */
2297 : unsigned char _map_type;
2298 :
2299 : /**
2300 : * Mapping function data; currently used when needed to store the
2301 : * RATIONAL_BERNSTEIN nodal weight data index.
2302 : */
2303 : unsigned char _map_data;
2304 : };
2305 :
2306 :
2307 :
2308 : // ------------------------------------------------------------
2309 : // Elem helper classes
2310 : //
2311 : class
2312 : Elem::NodeRefIter : public PointerToPointerIter<Node>
2313 : {
2314 : public:
2315 22574732 : NodeRefIter (Node * const * nodepp) : PointerToPointerIter<Node>(nodepp) {}
2316 : };
2317 :
2318 :
2319 : class
2320 : Elem::ConstNodeRefIter : public PointerToPointerIter<const Node>
2321 : {
2322 : public:
2323 17714944 : ConstNodeRefIter (const Node * const * nodepp) : PointerToPointerIter<const Node>(nodepp) {}
2324 : };
2325 :
2326 :
2327 : #ifdef LIBMESH_ENABLE_AMR
2328 : class
2329 : Elem::ChildRefIter : public PointerToPointerIter<Elem>
2330 : {
2331 : public:
2332 6438024 : ChildRefIter (Elem * const * childpp) : PointerToPointerIter<Elem>(childpp) {}
2333 : };
2334 :
2335 :
2336 : class
2337 : Elem::ConstChildRefIter : public PointerToPointerIter<const Elem>
2338 : {
2339 : public:
2340 2483426 : ConstChildRefIter (const Elem * const * childpp) : PointerToPointerIter<const Elem>(childpp) {}
2341 : };
2342 :
2343 :
2344 :
2345 : inline
2346 57316340 : SimpleRange<Elem::ChildRefIter> Elem::child_ref_range()
2347 : {
2348 3219012 : libmesh_assert(_children);
2349 59997033 : return {_children.get(), _children.get() + this->n_children()};
2350 : }
2351 :
2352 :
2353 : inline
2354 6674206 : SimpleRange<Elem::ConstChildRefIter> Elem::child_ref_range() const
2355 : {
2356 1241713 : libmesh_assert(_children);
2357 6901529 : return {_children.get(), _children.get() + this->n_children()};
2358 : }
2359 : #endif // LIBMESH_ENABLE_AMR
2360 :
2361 :
2362 :
2363 :
2364 : // ------------------------------------------------------------
2365 : // global Elem functions
2366 :
2367 : inline
2368 0 : std::ostream & operator << (std::ostream & os, const Elem & e)
2369 : {
2370 0 : e.print_info(os);
2371 0 : return os;
2372 : }
2373 :
2374 :
2375 : // ------------------------------------------------------------
2376 : // Elem class member functions
2377 : inline
2378 1131182165 : Elem::Elem(const unsigned int nn,
2379 : const unsigned int ns,
2380 : Elem * p,
2381 : Elem ** elemlinkdata,
2382 1131182165 : Node ** nodelinkdata) :
2383 658400884 : _nodes(nodelinkdata),
2384 658400884 : _elemlinks(elemlinkdata),
2385 658400884 : _sbd_id(0),
2386 : #ifdef LIBMESH_ENABLE_AMR
2387 658400884 : _rflag(Elem::DO_NOTHING),
2388 658400884 : _pflag(Elem::DO_NOTHING),
2389 658400884 : _p_level(0),
2390 : #endif
2391 658982800 : _map_type(p ? p->mapping_type() : 0),
2392 1159984618 : _map_data(p ? p->mapping_data() : 0)
2393 : {
2394 1131182165 : this->processor_id() = DofObject::invalid_processor_id;
2395 :
2396 : // If this ever legitimately fails we need to increase max_n_nodes
2397 57048105 : libmesh_assert_less_equal(nn, max_n_nodes);
2398 :
2399 : // We currently only support refinement of elements into child
2400 : // elements of the same type. We can't test elem->type() here,
2401 : // because that's virtual and we're still in the base class
2402 : // constructor, but we can at least usually verify constency with
2403 : // the arguments we were handed.
2404 : #ifndef NDEBUG
2405 57048105 : if (p && !p->runtime_topology())
2406 : {
2407 287216 : libmesh_assert_equal_to(nn, p->n_nodes());
2408 287216 : libmesh_assert_equal_to(ns, p->n_sides());
2409 : }
2410 : #endif
2411 :
2412 : // Initialize the nodes data structure if we're given a pointer to
2413 : // memory for it.
2414 1131182165 : if (_nodes)
2415 : {
2416 5364286509 : for (unsigned int n=0; n<nn; n++)
2417 4233152891 : _nodes[n] = nullptr;
2418 : }
2419 :
2420 : // Initialize the neighbors/parent data structure
2421 : // _elemlinks = new Elem *[ns+1];
2422 :
2423 : // Initialize the elements data structure if we're given a pointer
2424 : // to memory for it. If we *weren't* given memory for it, e.g.
2425 : // because a subclass like an arbitrary Polygon needs to
2426 : // heap-allocate this memory, then that subclass will have to handle
2427 : // this initialization too.
2428 1131182165 : if (_elemlinks)
2429 : {
2430 1131149997 : _elemlinks[0] = p;
2431 :
2432 5171663553 : for (unsigned int n=1; n<ns+1; n++)
2433 4040513556 : _elemlinks[n] = nullptr;
2434 :
2435 : // Optionally initialize data from the parent
2436 1131149997 : if (this->parent())
2437 : {
2438 28802453 : this->subdomain_id() = this->parent()->subdomain_id();
2439 28802453 : this->processor_id() = this->parent()->processor_id();
2440 28802453 : _map_type = this->parent()->_map_type;
2441 28802453 : _map_data = this->parent()->_map_data;
2442 :
2443 : #ifdef LIBMESH_ENABLE_AMR
2444 29097153 : this->set_p_level(this->parent()->p_level());
2445 : #endif
2446 : }
2447 : }
2448 1131182165 : }
2449 :
2450 :
2451 :
2452 : inline
2453 3417382241 : const Point & Elem::point (const unsigned int i) const
2454 : {
2455 3417382241 : libmesh_assert_less (i, this->n_nodes());
2456 3417382241 : libmesh_assert(_nodes[i]);
2457 3417382241 : libmesh_assert_not_equal_to (_nodes[i]->id(), Node::invalid_id);
2458 :
2459 47420710097 : return *_nodes[i];
2460 : }
2461 :
2462 :
2463 :
2464 : inline
2465 2660622 : Point & Elem::point (const unsigned int i)
2466 : {
2467 2660622 : libmesh_assert_less (i, this->n_nodes());
2468 :
2469 167008525 : return *_nodes[i];
2470 : }
2471 :
2472 :
2473 :
2474 : inline
2475 87767806 : dof_id_type Elem::node_id (const unsigned int i) const
2476 : {
2477 87767806 : libmesh_assert_less (i, this->n_nodes());
2478 87767806 : libmesh_assert(_nodes[i]);
2479 87767806 : libmesh_assert_not_equal_to (_nodes[i]->id(), Node::invalid_id);
2480 :
2481 978150922 : return _nodes[i]->id();
2482 : }
2483 :
2484 :
2485 :
2486 : inline
2487 489210 : unsigned int Elem::local_node (const dof_id_type i) const
2488 : {
2489 4385954 : for (auto n : make_range(this->n_nodes()))
2490 4385954 : if (this->node_id(n) == i)
2491 14180 : return n;
2492 :
2493 0 : return libMesh::invalid_uint;
2494 : }
2495 :
2496 :
2497 :
2498 : inline
2499 23413305 : const Node * const * Elem::get_nodes () const
2500 : {
2501 130815390 : return _nodes;
2502 : }
2503 :
2504 :
2505 :
2506 : inline
2507 113358762 : const Node * Elem::node_ptr (const unsigned int i) const
2508 : {
2509 113358762 : libmesh_assert_less (i, this->n_nodes());
2510 113358762 : libmesh_assert(_nodes[i]);
2511 :
2512 24228745030 : return _nodes[i];
2513 : }
2514 :
2515 :
2516 :
2517 : inline
2518 206788501 : Node * Elem::node_ptr (const unsigned int i)
2519 : {
2520 206788501 : libmesh_assert_less (i, this->n_nodes());
2521 206788501 : libmesh_assert(_nodes[i]);
2522 :
2523 4714648114 : return _nodes[i];
2524 : }
2525 :
2526 :
2527 :
2528 : inline
2529 19108430 : const Node & Elem::node_ref (const unsigned int i) const
2530 : {
2531 175165868 : return *this->node_ptr(i);
2532 : }
2533 :
2534 :
2535 :
2536 : inline
2537 26066614 : Node & Elem::node_ref (const unsigned int i)
2538 : {
2539 54599782 : return *this->node_ptr(i);
2540 : }
2541 :
2542 :
2543 :
2544 : inline
2545 1582548 : unsigned int Elem::get_node_index (const Node * node_ptr) const
2546 : {
2547 11605332 : for (auto n : make_range(this->n_nodes()))
2548 11593420 : if (this->_nodes[n] == node_ptr)
2549 160618 : return n;
2550 :
2551 0 : return libMesh::invalid_uint;
2552 : }
2553 :
2554 :
2555 :
2556 : #ifdef LIBMESH_ENABLE_DEPRECATED
2557 : inline
2558 0 : Node * & Elem::set_node (const unsigned int i)
2559 : {
2560 0 : libmesh_assert_less (i, this->n_nodes());
2561 :
2562 : libmesh_deprecated();
2563 :
2564 0 : return _nodes[i];
2565 : }
2566 : #endif // LIBMESH_ENABLE_DEPRECATED
2567 :
2568 :
2569 :
2570 : inline
2571 6243932322 : void Elem::set_node (const unsigned int i,
2572 : Node * node)
2573 : {
2574 177702 : libmesh_assert_less (i, this->n_nodes());
2575 :
2576 6532121837 : _nodes[i] = node;
2577 5880014622 : }
2578 :
2579 :
2580 :
2581 : inline
2582 57627387 : subdomain_id_type Elem::subdomain_id () const
2583 : {
2584 530259858 : return _sbd_id;
2585 : }
2586 :
2587 :
2588 :
2589 : inline
2590 42932469 : subdomain_id_type & Elem::subdomain_id ()
2591 : {
2592 124484898 : return _sbd_id;
2593 : }
2594 :
2595 :
2596 :
2597 : inline
2598 98444635 : const Elem * Elem::neighbor_ptr (unsigned int i) const
2599 : {
2600 98444635 : libmesh_assert_less (i, this->n_neighbors());
2601 :
2602 467770564 : return _elemlinks[i+1];
2603 : }
2604 :
2605 :
2606 :
2607 : inline
2608 18202813 : Elem * Elem::neighbor_ptr (unsigned int i)
2609 : {
2610 18202813 : libmesh_assert_less (i, this->n_neighbors());
2611 :
2612 1218359443 : return _elemlinks[i+1];
2613 : }
2614 :
2615 :
2616 :
2617 : inline
2618 9747370 : void Elem::set_neighbor (const unsigned int i, Elem * n)
2619 : {
2620 9747370 : libmesh_assert_less (i, this->n_neighbors());
2621 :
2622 405737878 : _elemlinks[i+1] = n;
2623 508185802 : }
2624 :
2625 :
2626 :
2627 : inline
2628 91038344 : bool Elem::has_neighbor (const Elem * elem) const
2629 : {
2630 272326604 : for (auto n : this->neighbor_ptr_range())
2631 261607843 : if (n == elem)
2632 80622214 : return true;
2633 :
2634 10416130 : return false;
2635 : }
2636 :
2637 :
2638 :
2639 : inline
2640 : Elem * Elem::child_neighbor (Elem * elem)
2641 : {
2642 : for (auto n : elem->neighbor_ptr_range())
2643 : if (n && n->parent() == this)
2644 : return n;
2645 :
2646 : return nullptr;
2647 : }
2648 :
2649 :
2650 :
2651 : inline
2652 : const Elem * Elem::child_neighbor (const Elem * elem) const
2653 : {
2654 : for (auto n : elem->neighbor_ptr_range())
2655 : if (n && n->parent() == this)
2656 : return n;
2657 :
2658 : return nullptr;
2659 : }
2660 :
2661 :
2662 :
2663 : inline
2664 : SimpleRange<Elem::NodeRefIter>
2665 192233172 : Elem::node_ref_range()
2666 : {
2667 201726102 : return {_nodes, _nodes+this->n_nodes()};
2668 : }
2669 :
2670 :
2671 :
2672 : inline
2673 : SimpleRange<Elem::ConstNodeRefIter>
2674 261972210 : Elem::node_ref_range() const
2675 : {
2676 279197531 : return {_nodes, _nodes+this->n_nodes()};
2677 : }
2678 :
2679 :
2680 :
2681 : inline
2682 : IntRange<unsigned short>
2683 60954974 : Elem::node_index_range() const
2684 : {
2685 1406919005 : return {0, cast_int<unsigned short>(this->n_nodes())};
2686 : }
2687 :
2688 :
2689 :
2690 : inline
2691 : IntRange<unsigned short>
2692 410691 : Elem::edge_index_range() const
2693 : {
2694 41942485 : return {0, cast_int<unsigned short>(this->n_edges())};
2695 : }
2696 :
2697 :
2698 :
2699 : inline
2700 : IntRange<unsigned short>
2701 308674 : Elem::face_index_range() const
2702 : {
2703 3655740 : return {0, cast_int<unsigned short>(this->n_faces())};
2704 : }
2705 :
2706 :
2707 :
2708 : inline
2709 : IntRange<unsigned short>
2710 8891330 : Elem::side_index_range() const
2711 : {
2712 357402518 : return {0, cast_int<unsigned short>(this->n_sides())};
2713 : }
2714 :
2715 :
2716 :
2717 :
2718 : inline
2719 826716870 : std::unique_ptr<const Elem> Elem::side_ptr (unsigned int i) const
2720 : {
2721 : // Call the non-const version of this function, return the result as
2722 : // a std::unique_ptr<const Elem>.
2723 75970818 : Elem * me = const_cast<Elem *>(this);
2724 902687688 : return me->side_ptr(i);
2725 : }
2726 :
2727 :
2728 :
2729 : inline
2730 : void
2731 15336 : Elem::side_ptr (std::unique_ptr<const Elem> & elem,
2732 : const unsigned int i) const
2733 : {
2734 : // Hand off to the non-const version of this function
2735 432 : Elem * me = const_cast<Elem *>(this);
2736 864 : std::unique_ptr<Elem> e {const_cast<Elem *>(elem.release())};
2737 15336 : me->side_ptr(e, i);
2738 14904 : elem = std::move(e);
2739 15336 : }
2740 :
2741 :
2742 :
2743 : inline
2744 : std::unique_ptr<const Elem>
2745 105832738 : Elem::build_side_ptr (const unsigned int i) const
2746 : {
2747 : // Call the non-const version of this function, return the result as
2748 : // a std::unique_ptr<const Elem>.
2749 2828158 : Elem * me = const_cast<Elem *>(this);
2750 115975310 : return me->build_side_ptr(i);
2751 : }
2752 :
2753 :
2754 :
2755 : inline
2756 : void
2757 38454834 : Elem::build_side_ptr (std::unique_ptr<const Elem> & elem,
2758 : const unsigned int i) const
2759 : {
2760 : // Hand off to the non-const version of this function
2761 827644 : Elem * me = const_cast<Elem *>(this);
2762 1655288 : std::unique_ptr<Elem> e {const_cast<Elem *>(elem.release())};
2763 38454834 : me->build_side_ptr(e, i);
2764 36801655 : elem = std::move(e);
2765 38454834 : }
2766 :
2767 :
2768 :
2769 : template <typename Sideclass, typename Subclass>
2770 : inline
2771 : std::unique_ptr<Elem>
2772 155767396 : Elem::simple_build_side_ptr (const unsigned int i)
2773 : {
2774 33123278 : libmesh_assert_less (i, this->n_sides());
2775 :
2776 155767396 : std::unique_ptr<Elem> face = std::make_unique<Sideclass>();
2777 1200819175 : for (auto n : face->node_index_range())
2778 1045051779 : face->set_node(n, this->node_ptr(Subclass::side_nodes_map[i][n]));
2779 :
2780 155767396 : face->set_interior_parent(this);
2781 144786270 : face->inherit_data_from(*this);
2782 :
2783 155767396 : return face;
2784 0 : }
2785 :
2786 :
2787 :
2788 : template <typename Subclass>
2789 : inline
2790 : void
2791 38136889 : Elem::simple_build_side_ptr (std::unique_ptr<Elem> & side,
2792 : const unsigned int i,
2793 : ElemType sidetype)
2794 : {
2795 1665316 : libmesh_assert_less (i, this->n_sides());
2796 :
2797 38136889 : if (!side.get() || side->type() != sidetype)
2798 : {
2799 147137 : Subclass & real_me = cast_ref<Subclass&>(*this);
2800 1448595 : side = real_me.Subclass::build_side_ptr(i);
2801 : }
2802 : else
2803 : {
2804 37339023 : side->set_interior_parent(this);
2805 35820839 : side->inherit_data_from(*this);
2806 265840699 : for (auto n : side->node_index_range())
2807 228501676 : side->set_node(n, this->node_ptr(Subclass::side_nodes_map[i][n]));
2808 : }
2809 38136889 : }
2810 :
2811 :
2812 :
2813 : template <typename Subclass, typename Mapclass>
2814 : inline
2815 : void
2816 237434443 : Elem::simple_side_ptr (std::unique_ptr<Elem> & side,
2817 : const unsigned int i,
2818 : ElemType sidetype)
2819 : {
2820 4771950 : libmesh_assert_less (i, this->n_sides());
2821 :
2822 237434443 : if (!side.get() || side->type() != sidetype)
2823 : {
2824 16110 : Subclass & real_me = cast_ref<Subclass&>(*this);
2825 1165816 : side = real_me.Subclass::side_ptr(i);
2826 : }
2827 : else
2828 : {
2829 241556912 : side->subdomain_id() = this->subdomain_id();
2830 :
2831 806745786 : for (auto n : side->node_index_range())
2832 569902306 : side->set_node(n, this->node_ptr(Mapclass::side_nodes_map[i][n]));
2833 : }
2834 237434443 : }
2835 :
2836 :
2837 :
2838 : inline
2839 : std::unique_ptr<const Elem>
2840 46699008 : Elem::build_edge_ptr (const unsigned int i) const
2841 : {
2842 : // Call the non-const version of this function, return the result as
2843 : // a std::unique_ptr<const Elem>.
2844 27212782 : Elem * me = const_cast<Elem *>(this);
2845 48808049 : return me->build_edge_ptr(i);
2846 : }
2847 :
2848 :
2849 :
2850 : inline
2851 : void
2852 95756 : Elem::build_edge_ptr (std::unique_ptr<const Elem> & elem,
2853 : const unsigned int i) const
2854 : {
2855 : // Hand off to the non-const version of this function
2856 2924 : Elem * me = const_cast<Elem *>(this);
2857 5848 : std::unique_ptr<Elem> e {const_cast<Elem *>(elem.release())};
2858 95756 : me->build_edge_ptr(e, i);
2859 92832 : elem = std::move(e);
2860 95756 : }
2861 :
2862 :
2863 : template <typename Edgeclass, typename Subclass>
2864 : inline
2865 : std::unique_ptr<Elem>
2866 27878431 : Elem::simple_build_edge_ptr (const unsigned int i)
2867 : {
2868 9812036 : libmesh_assert_less (i, this->n_edges());
2869 :
2870 27878431 : std::unique_ptr<Elem> edge = std::make_unique<Edgeclass>();
2871 :
2872 97685766 : for (auto n : edge->node_index_range())
2873 69807335 : edge->set_node(n, this->node_ptr(Subclass::edge_nodes_map[i][n]));
2874 :
2875 27878431 : edge->set_interior_parent(this);
2876 26092187 : edge->inherit_data_from(*this);
2877 :
2878 27878431 : return edge;
2879 0 : }
2880 :
2881 :
2882 :
2883 :
2884 : template <typename Subclass>
2885 : inline
2886 : void
2887 137744 : Elem::simple_build_edge_ptr (std::unique_ptr<Elem> & edge,
2888 : const unsigned int i,
2889 : ElemType edgetype)
2890 : {
2891 4196 : libmesh_assert_less (i, this->n_edges());
2892 :
2893 137744 : if (!edge.get() || edge->type() != edgetype)
2894 : {
2895 25 : Subclass & real_me = cast_ref<Subclass&>(*this);
2896 655 : edge = real_me.Subclass::build_edge_ptr(i);
2897 : }
2898 : else
2899 : {
2900 133233 : edge->inherit_data_from(*this);
2901 414504 : for (auto n : edge->node_index_range())
2902 277100 : edge->set_node(n, this->node_ptr(Subclass::edge_nodes_map[i][n]));
2903 : }
2904 137744 : }
2905 :
2906 :
2907 :
2908 : inline
2909 39 : bool Elem::on_boundary () const
2910 : {
2911 : // By convention, the element is on the boundary
2912 : // if it has a nullptr neighbor.
2913 351 : return this->has_neighbor(nullptr);
2914 : }
2915 :
2916 :
2917 :
2918 : inline
2919 155938134 : unsigned int Elem::which_neighbor_am_i (const Elem * e) const
2920 : {
2921 8215728 : libmesh_assert(e);
2922 :
2923 8215728 : const Elem * eparent = e;
2924 :
2925 156729608 : while (eparent->level() > this->level())
2926 : {
2927 263020 : eparent = eparent->parent();
2928 260944 : libmesh_assert(eparent);
2929 : }
2930 :
2931 408466446 : for (auto s : make_range(this->n_sides()))
2932 408458999 : if (this->neighbor_ptr(s) == eparent)
2933 8215728 : return s;
2934 :
2935 0 : return libMesh::invalid_uint;
2936 : }
2937 :
2938 :
2939 :
2940 : inline
2941 180562955 : bool Elem::active() const
2942 : {
2943 : #ifdef LIBMESH_ENABLE_AMR
2944 3758463266 : if ((this->refinement_flag() == INACTIVE) ||
2945 126867886 : (this->refinement_flag() == COARSEN_INACTIVE))
2946 53803755 : return false;
2947 : else
2948 126759200 : return true;
2949 : #else
2950 : return true;
2951 : #endif
2952 : }
2953 :
2954 :
2955 :
2956 :
2957 :
2958 : inline
2959 417253691 : bool Elem::subactive() const
2960 : {
2961 : #ifdef LIBMESH_ENABLE_AMR
2962 46093962 : if (this->active())
2963 34007967 : return false;
2964 12085995 : if (!this->has_children())
2965 3662996 : return true;
2966 85393454 : for (const Elem * my_ancestor = this->parent();
2967 297531560 : my_ancestor != nullptr;
2968 28662050 : my_ancestor = my_ancestor->parent())
2969 24850667 : if (my_ancestor->active())
2970 10516 : return true;
2971 : #endif
2972 :
2973 8412483 : return false;
2974 : }
2975 :
2976 :
2977 :
2978 : inline
2979 46255328 : bool Elem::has_children() const
2980 : {
2981 : #ifdef LIBMESH_ENABLE_AMR
2982 358472208 : if (!_children)
2983 9120301 : return false;
2984 : else
2985 37135027 : return true;
2986 : #else
2987 : return false;
2988 : #endif
2989 : }
2990 :
2991 :
2992 : inline
2993 : bool Elem::has_ancestor_children() const
2994 : {
2995 : #ifdef LIBMESH_ENABLE_AMR
2996 : if (!_children)
2997 : return false;
2998 : else
2999 : for (auto & c : child_ref_range())
3000 : if (c.has_children())
3001 : return true;
3002 : #endif
3003 : return false;
3004 : }
3005 :
3006 :
3007 :
3008 : inline
3009 0 : bool Elem::is_ancestor_of(const Elem *
3010 : #ifdef LIBMESH_ENABLE_AMR
3011 : descendant
3012 : #endif
3013 : ) const
3014 : {
3015 : #ifdef LIBMESH_ENABLE_AMR
3016 0 : const Elem * e = descendant;
3017 10626 : while (e)
3018 : {
3019 8638 : if (this == e)
3020 0 : return true;
3021 0 : e = e->parent();
3022 : }
3023 : #endif
3024 0 : return false;
3025 : }
3026 :
3027 :
3028 :
3029 : inline
3030 1176063936 : const Elem * Elem::parent () const
3031 : {
3032 2484017581 : return _elemlinks[0];
3033 : }
3034 :
3035 :
3036 :
3037 : inline
3038 98477048 : Elem * Elem::parent ()
3039 : {
3040 1584386004 : return _elemlinks[0];
3041 : }
3042 :
3043 :
3044 :
3045 : inline
3046 342794 : void Elem::set_parent (Elem * p)
3047 : {
3048 : // We no longer support using parent() as interior_parent()
3049 342794 : libmesh_assert_equal_to(this->dim(), p ? p->dim() : this->dim());
3050 10050050 : _elemlinks[0] = p;
3051 525761 : }
3052 :
3053 :
3054 :
3055 : inline
3056 824530 : const Elem * Elem::top_parent () const
3057 : {
3058 824530 : const Elem * tp = this;
3059 :
3060 : // Keep getting the element's parent
3061 : // until that parent is at level-0
3062 4370355 : while (tp->parent() != nullptr)
3063 2045469 : tp = tp->parent();
3064 :
3065 824530 : libmesh_assert(tp);
3066 824530 : libmesh_assert_equal_to (tp->level(), 0);
3067 :
3068 824530 : return tp;
3069 : }
3070 :
3071 :
3072 :
3073 : inline
3074 7220288398 : unsigned int Elem::level() const
3075 : {
3076 : #ifdef LIBMESH_ENABLE_AMR
3077 :
3078 : // if I don't have a parent I was
3079 : // created directly from file
3080 : // or by the user, so I am a
3081 : // level-0 element
3082 19962356362 : if (this->parent() == nullptr)
3083 112104713 : return 0;
3084 :
3085 : // if the parent and this element are of different
3086 : // dimensionality we are at the same level as
3087 : // the parent (e.g. we are the 2D side of a
3088 : // 3D element)
3089 13041110071 : if (this->dim() != this->parent()->dim())
3090 0 : return this->parent()->level();
3091 :
3092 : // otherwise we are at a level one
3093 : // higher than our parent
3094 13164294483 : return (this->parent()->level() + 1);
3095 :
3096 : #else
3097 :
3098 : // Without AMR all elements are
3099 : // at level 0.
3100 : return 0;
3101 :
3102 : #endif
3103 : }
3104 :
3105 :
3106 :
3107 : inline
3108 2172897227 : unsigned int Elem::p_level() const
3109 : {
3110 : #ifdef LIBMESH_ENABLE_AMR
3111 53656004728 : return _p_level;
3112 : #else
3113 : return 0;
3114 : #endif
3115 : }
3116 :
3117 :
3118 :
3119 : inline
3120 100875910 : ElemMappingType Elem::mapping_type () const
3121 : {
3122 1260769642 : return static_cast<ElemMappingType>(_map_type);
3123 : }
3124 :
3125 :
3126 :
3127 : inline
3128 46752459 : void Elem::set_mapping_type(const ElemMappingType type)
3129 : {
3130 301183312 : _map_type = cast_int<unsigned char>(type);
3131 46752459 : }
3132 :
3133 :
3134 :
3135 : inline
3136 33328490 : unsigned char Elem::mapping_data () const
3137 : {
3138 138115602 : return _map_data;
3139 : }
3140 :
3141 :
3142 :
3143 : inline
3144 46752452 : void Elem::set_mapping_data(const unsigned char data)
3145 : {
3146 300920538 : _map_data = data;
3147 47015559 : }
3148 :
3149 :
3150 :
3151 : #ifdef LIBMESH_ENABLE_AMR
3152 :
3153 : inline
3154 0 : const Elem * Elem::raw_child_ptr (unsigned int i) const
3155 : {
3156 4240 : if (!_children)
3157 0 : return nullptr;
3158 :
3159 2880 : return _children[i];
3160 : }
3161 :
3162 : inline
3163 90033288 : const Elem * Elem::child_ptr (unsigned int i) const
3164 : {
3165 90033288 : libmesh_assert(_children);
3166 90033288 : libmesh_assert(_children[i]);
3167 :
3168 337242775 : return _children[i];
3169 : }
3170 :
3171 : inline
3172 2095948 : Elem * Elem::child_ptr (unsigned int i)
3173 : {
3174 2095948 : libmesh_assert(_children);
3175 2095948 : libmesh_assert(_children[i]);
3176 :
3177 23606432 : return _children[i];
3178 : }
3179 :
3180 :
3181 : inline
3182 133432 : void Elem::set_child (unsigned int c, Elem * elem)
3183 : {
3184 133432 : libmesh_assert (this->has_children());
3185 :
3186 52216507 : _children[c] = elem;
3187 27078529 : }
3188 :
3189 :
3190 :
3191 : inline
3192 109568436 : unsigned int Elem::which_child_am_i (const Elem * e) const
3193 : {
3194 27330114 : libmesh_assert(e);
3195 27330114 : libmesh_assert (this->has_children());
3196 :
3197 109568436 : unsigned int nc = this->n_children();
3198 292061801 : for (unsigned int c=0; c != nc; c++)
3199 292061801 : if (this->child_ptr(c) == e)
3200 27330114 : return c;
3201 :
3202 0 : libmesh_error_msg("ERROR: which_child_am_i() was called with a non-child!");
3203 :
3204 : return libMesh::invalid_uint;
3205 : }
3206 :
3207 :
3208 :
3209 : inline
3210 352189607 : Elem::RefinementState Elem::refinement_flag () const
3211 : {
3212 4600349316 : return static_cast<RefinementState>(_rflag);
3213 : }
3214 :
3215 :
3216 :
3217 : inline
3218 3563222 : void Elem::set_refinement_flag(RefinementState rflag)
3219 : {
3220 60230522 : _rflag = cast_int<unsigned char>(rflag);
3221 11435464 : }
3222 :
3223 :
3224 :
3225 : inline
3226 86954956 : Elem::RefinementState Elem::p_refinement_flag () const
3227 : {
3228 234303843 : return static_cast<RefinementState>(_pflag);
3229 : }
3230 :
3231 :
3232 :
3233 : inline
3234 2520938 : void Elem::set_p_refinement_flag(RefinementState pflag)
3235 : {
3236 2520590 : if (this->p_level() == 0)
3237 2515304 : libmesh_assert_not_equal_to
3238 : (pflag, Elem::JUST_REFINED);
3239 :
3240 40625417 : _pflag = cast_int<unsigned char>(pflag);
3241 25173199 : }
3242 :
3243 :
3244 :
3245 : inline
3246 0 : unsigned int Elem::max_descendant_p_level () const
3247 : {
3248 : // This is undefined for subactive elements,
3249 : // which have no active descendants
3250 0 : libmesh_assert (!this->subactive());
3251 0 : if (this->active())
3252 0 : return this->p_level();
3253 :
3254 0 : unsigned int max_p_level = _p_level;
3255 0 : for (auto & c : child_ref_range())
3256 0 : max_p_level = std::max(max_p_level,
3257 0 : c.max_descendant_p_level());
3258 0 : return max_p_level;
3259 : }
3260 :
3261 :
3262 :
3263 : inline
3264 46719312 : void Elem::hack_p_level(unsigned int p)
3265 : {
3266 46719312 : if (p == 0)
3267 46635093 : libmesh_assert_not_equal_to
3268 : (this->p_refinement_flag(), Elem::JUST_REFINED);
3269 :
3270 356170601 : _p_level = cast_int<unsigned char>(p);
3271 295499477 : }
3272 :
3273 :
3274 : inline
3275 3602 : void Elem::hack_p_level_and_refinement_flag (unsigned int p,
3276 : RefinementState pflag)
3277 : {
3278 34520678 : _pflag = cast_int<unsigned char>(pflag);
3279 3602 : this->hack_p_level(p);
3280 34501851 : }
3281 :
3282 : #endif // ifdef LIBMESH_ENABLE_AMR
3283 :
3284 :
3285 : inline
3286 290889 : void Elem::orient(BoundaryInfo * boundary_info)
3287 : {
3288 304149 : if (this->is_flipped())
3289 138978 : this->flip(boundary_info);
3290 290889 : }
3291 :
3292 :
3293 : inline
3294 43734 : dof_id_type Elem::compute_key (dof_id_type n0)
3295 : {
3296 43734 : return n0;
3297 : }
3298 :
3299 :
3300 :
3301 : inline
3302 6482525 : dof_id_type Elem::compute_key (dof_id_type n0,
3303 : dof_id_type n1)
3304 : {
3305 : // Order the two so that n0 < n1
3306 177644014 : if (n0 > n1) std::swap (n0, n1);
3307 :
3308 177644014 : return Utility::hashword2(n0, n1);
3309 : }
3310 :
3311 :
3312 :
3313 : inline
3314 49007988 : dof_id_type Elem::compute_key (dof_id_type n0,
3315 : dof_id_type n1,
3316 : dof_id_type n2)
3317 : {
3318 49007988 : std::array<dof_id_type, 3> array = {{n0, n1, n2}};
3319 1765184 : std::sort(array.begin(), array.end());
3320 50773172 : return Utility::hashword(array);
3321 : }
3322 :
3323 :
3324 :
3325 : inline
3326 35834259 : dof_id_type Elem::compute_key (dof_id_type n0,
3327 : dof_id_type n1,
3328 : dof_id_type n2,
3329 : dof_id_type n3)
3330 : {
3331 35834259 : std::array<dof_id_type, 4> array = {{n0, n1, n2, n3}};
3332 1075098 : std::sort(array.begin(), array.end());
3333 36909357 : return Utility::hashword(array);
3334 : }
3335 :
3336 :
3337 :
3338 : inline
3339 227319879 : void Elem::inherit_data_from (const Elem & src)
3340 : {
3341 60780466 : this->set_mapping_type(src.mapping_type());
3342 60780466 : this->set_mapping_data(src.mapping_data());
3343 242618350 : this->subdomain_id() = src.subdomain_id();
3344 60780466 : this->processor_id(src.processor_id());
3345 : #ifdef LIBMESH_ENABLE_AMR
3346 242618350 : this->set_p_level(src.p_level());
3347 : #endif
3348 227319879 : }
3349 :
3350 :
3351 :
3352 : /**
3353 : * The definition of the protected nested SideIter class.
3354 : */
3355 0 : class Elem::SideIter
3356 : {
3357 : public:
3358 : // Constructor with arguments.
3359 0 : SideIter(const unsigned int side_number,
3360 : Elem * parent)
3361 0 : : _side(),
3362 0 : _side_ptr(nullptr),
3363 0 : _parent(parent),
3364 0 : _side_number(side_number)
3365 0 : {}
3366 :
3367 :
3368 : // Empty constructor.
3369 : SideIter()
3370 : : _side(),
3371 : _side_ptr(nullptr),
3372 : _parent(nullptr),
3373 : _side_number(libMesh::invalid_uint)
3374 : {}
3375 :
3376 :
3377 : // Copy constructor
3378 0 : SideIter(const SideIter & other)
3379 0 : : _side(),
3380 0 : _side_ptr(nullptr),
3381 0 : _parent(other._parent),
3382 0 : _side_number(other._side_number)
3383 0 : {}
3384 :
3385 :
3386 : // op=
3387 : SideIter & operator=(const SideIter & other)
3388 : {
3389 : this->_parent = other._parent;
3390 : this->_side_number = other._side_number;
3391 : return *this;
3392 : }
3393 :
3394 : // unary op*
3395 0 : Elem *& operator*() const
3396 : {
3397 : // Set the std::unique_ptr
3398 0 : this->_update_side_ptr();
3399 :
3400 : // Return a reference to _side_ptr
3401 0 : return this->_side_ptr;
3402 : }
3403 :
3404 : // op++
3405 0 : SideIter & operator++()
3406 : {
3407 0 : ++_side_number;
3408 0 : return *this;
3409 : }
3410 :
3411 : // op== Two side iterators are equal if they have
3412 : // the same side number and the same parent element.
3413 0 : bool operator == (const SideIter & other) const
3414 : {
3415 0 : return (this->_side_number == other._side_number &&
3416 0 : this->_parent == other._parent);
3417 : }
3418 :
3419 :
3420 : // Consults the parent Elem to determine if the side
3421 : // is a boundary side. Note: currently side N is a
3422 : // boundary side if neighbor N is nullptr. Be careful,
3423 : // this could possibly change in the future?
3424 0 : bool side_on_boundary() const
3425 : {
3426 0 : return this->_parent->neighbor_ptr(_side_number) == nullptr;
3427 : }
3428 :
3429 : private:
3430 : // Update the _side pointer by building the correct side.
3431 : // This has to be called before dereferencing.
3432 0 : void _update_side_ptr() const
3433 : {
3434 : // Construct new side, store in std::unique_ptr
3435 0 : this->_side = this->_parent->build_side_ptr(this->_side_number);
3436 :
3437 : // Also set our internal naked pointer. Memory is still owned
3438 : // by the std::unique_ptr.
3439 0 : this->_side_ptr = _side.get();
3440 0 : }
3441 :
3442 : // std::unique_ptr to the actual side, handles memory management for
3443 : // the sides which are created during the course of iteration.
3444 : mutable std::unique_ptr<Elem> _side;
3445 :
3446 : // Raw pointer needed to facilitate passing back to the user a
3447 : // reference to a non-temporary raw pointer in order to conform to
3448 : // the variant_filter_iterator interface. It points to the same
3449 : // thing the std::unique_ptr "_side" above holds. What happens if the user
3450 : // calls delete on the pointer passed back? Well, this is an issue
3451 : // which is not addressed by the iterators in libMesh. Basically it
3452 : // is a bad idea to ever call delete on an iterator from the library.
3453 : mutable Elem * _side_ptr;
3454 :
3455 : // Pointer to the parent Elem class which generated this iterator
3456 : Elem * _parent;
3457 :
3458 : // A counter variable which keeps track of the side number
3459 : unsigned int _side_number;
3460 : };
3461 :
3462 :
3463 :
3464 :
3465 :
3466 :
3467 : // Private implementation functions in the Elem class for the side iterators.
3468 : // They have to come after the definition of the SideIter class.
3469 : inline
3470 0 : Elem::SideIter Elem::_first_side()
3471 : {
3472 0 : return SideIter(0, this);
3473 : }
3474 :
3475 :
3476 :
3477 : inline
3478 0 : Elem::SideIter Elem::_last_side()
3479 : {
3480 0 : return SideIter(this->n_neighbors(), this);
3481 : }
3482 :
3483 :
3484 :
3485 :
3486 : /**
3487 : * The definition of the struct used for iterating over sides.
3488 : */
3489 : struct
3490 : Elem::side_iterator : variant_filter_iterator<Elem::Predicate, Elem *>
3491 : {
3492 : // Templated forwarding ctor -- forwards to appropriate variant_filter_iterator ctor
3493 : template <typename PredType, typename IterType>
3494 0 : side_iterator (const IterType & d,
3495 : const IterType & e,
3496 : const PredType & p ) :
3497 0 : variant_filter_iterator<Elem::Predicate, Elem *>(d,e,p) {}
3498 : };
3499 :
3500 :
3501 :
3502 : inline
3503 77818706 : SimpleRange<Elem::NeighborPtrIter> Elem::neighbor_ptr_range()
3504 : {
3505 80891005 : return {_elemlinks+1, _elemlinks + 1 + this->n_neighbors()};
3506 : }
3507 :
3508 :
3509 : inline
3510 386731604 : SimpleRange<Elem::ConstNeighborPtrIter> Elem::neighbor_ptr_range() const
3511 : {
3512 388151768 : return {_elemlinks+1, _elemlinks + 1 + this->n_neighbors()};
3513 : }
3514 :
3515 : } // namespace libMesh
3516 :
3517 :
3518 : // Helper function for default caches in Elem subclasses
3519 :
3520 : #define LIBMESH_ENABLE_TOPOLOGY_CACHES \
3521 : virtual \
3522 : std::vector<std::vector<std::vector<std::vector<std::pair<unsigned char, unsigned char>>>>> & \
3523 : _get_bracketing_node_cache() const override \
3524 : { \
3525 : static std::vector<std::vector<std::vector<std::vector<std::pair<unsigned char, unsigned char>>>>> c; \
3526 : return c; \
3527 : } \
3528 : \
3529 : virtual \
3530 : std::vector<std::vector<std::vector<signed char>>> & \
3531 : _get_parent_indices_cache() const override \
3532 : { \
3533 : static std::vector<std::vector<std::vector<signed char>>> c; \
3534 : return c; \
3535 : }
3536 :
3537 :
3538 :
3539 :
3540 :
3541 :
3542 : #endif // LIBMESH_ELEM_H
|
