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