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