libMesh::FEAbstract Class Referenceabstract

#include <fe_abstract.h>

Inheritance diagram for libMesh::FEAbstract:

Public Member Functions

virtual ~FEAbstract ()
 
virtual void reinit (const Elem *elem, const std::vector< Point > *const pts=libmesh_nullptr, const std::vector< Real > *const weights=libmesh_nullptr)=0
 
virtual void reinit (const Elem *elem, const unsigned int side, const Real tolerance=TOLERANCE, const std::vector< Point > *const pts=libmesh_nullptr, const std::vector< Real > *const weights=libmesh_nullptr)=0
 
virtual void edge_reinit (const Elem *elem, const unsigned int edge, const Real tolerance=TOLERANCE, const std::vector< Point > *pts=libmesh_nullptr, const std::vector< Real > *weights=libmesh_nullptr)=0
 
virtual void side_map (const Elem *elem, const Elem *side, const unsigned int s, const std::vector< Point > &reference_side_points, std::vector< Point > &reference_points)=0
 
unsigned int get_dim () const
 
const std::vector< Point > & get_xyz () const
 
const std::vector< Real > & get_JxW () const
 
const std::vector< RealGradient > & get_dxyzdxi () const
 
const std::vector< RealGradient > & get_dxyzdeta () const
 
const std::vector< RealGradient > & get_dxyzdzeta () const
 
const std::vector< RealGradient > & get_d2xyzdxi2 () const
 
const std::vector< RealGradient > & get_d2xyzdeta2 () const
 
const std::vector< RealGradient > & get_d2xyzdzeta2 () const
 
const std::vector< RealGradient > & get_d2xyzdxideta () const
 
const std::vector< RealGradient > & get_d2xyzdxidzeta () const
 
const std::vector< RealGradient > & get_d2xyzdetadzeta () const
 
const std::vector< Real > & get_dxidx () const
 
const std::vector< Real > & get_dxidy () const
 
const std::vector< Real > & get_dxidz () const
 
const std::vector< Real > & get_detadx () const
 
const std::vector< Real > & get_detady () const
 
const std::vector< Real > & get_detadz () const
 
const std::vector< Real > & get_dzetadx () const
 
const std::vector< Real > & get_dzetady () const
 
const std::vector< Real > & get_dzetadz () const
 
const std::vector< std::vector< Point > > & get_tangents () const
 
const std::vector< Point > & get_normals () const
 
const std::vector< Real > & get_curvatures () const
 
virtual void attach_quadrature_rule (QBase *q)=0
 
virtual unsigned int n_shape_functions () const =0
 
virtual unsigned int n_quadrature_points () const =0
 
ElemType get_type () const
 
unsigned int get_p_level () const
 
FEType get_fe_type () const
 
Order get_order () const
 
void set_fe_order (int new_order)
 
virtual FEContinuity get_continuity () const =0
 
virtual bool is_hierarchic () const =0
 
FEFamily get_family () const
 
const FEMapget_fe_map () const
 
void print_JxW (std::ostream &os) const
 
virtual void print_phi (std::ostream &os) const =0
 
virtual void print_dphi (std::ostream &os) const =0
 
virtual void print_d2phi (std::ostream &os) const =0
 
void print_xyz (std::ostream &os) const
 
void print_info (std::ostream &os) const
 

Static Public Member Functions

static std::unique_ptr< FEAbstractbuild (const unsigned int dim, const FEType &type)
 
static bool on_reference_element (const Point &p, const ElemType t, const Real eps=TOLERANCE)
 
static void get_refspace_nodes (const ElemType t, std::vector< Point > &nodes)
 
static void compute_node_constraints (NodeConstraints &constraints, const Elem *elem)
 
static void compute_periodic_node_constraints (NodeConstraints &constraints, const PeriodicBoundaries &boundaries, const MeshBase &mesh, const PointLocatorBase *point_locator, const Elem *elem)
 
static std::string get_info ()
 
static void print_info (std::ostream &out=libMesh::out)
 
static unsigned int n_objects ()
 
static void enable_print_counter_info ()
 
static void disable_print_counter_info ()
 

Protected Types

typedef std::map< std::string, std::pair< unsigned int, unsigned int > > Counts
 

Protected Member Functions

 FEAbstract (const unsigned int dim, const FEType &fet)
 
virtual void compute_shape_functions (const Elem *, const std::vector< Point > &)=0
 
virtual bool shapes_need_reinit () const =0
 
void increment_constructor_count (const std::string &name)
 
void increment_destructor_count (const std::string &name)
 

Protected Attributes

std::unique_ptr< FEMap_fe_map
 
const unsigned int dim
 
bool calculations_started
 
bool calculate_phi
 
bool calculate_dphi
 
bool calculate_d2phi
 
bool calculate_curl_phi
 
bool calculate_div_phi
 
bool calculate_dphiref
 
FEType fe_type
 
ElemType elem_type
 
unsigned int _p_level
 
QBaseqrule
 
bool shapes_on_quadrature
 

Static Protected Attributes

static Counts _counts
 
static Threads::atomic< unsigned int > _n_objects
 
static Threads::spin_mutex _mutex
 
static bool _enable_print_counter = true
 

Friends

std::ostream & operator<< (std::ostream &os, const FEAbstract &fe)
 

Detailed Description

This class forms the foundation from which generic finite elements may be derived. In the current implementation, the templated derived class FE offers a wide variety of commonly used finite element concepts. Check there for details. Use the FEAbstract::build() method to create an object of any of the derived classes.

Note
In the present design, the number of virtual members is kept to a minimum for performance reasons, although this is not based on rigorous profiling.

All calls to static members of the FE classes should be requested through the FEInterface. This interface class approximates runtime polymorphism for the templated finite element classes. Even internal library classes, like DofMap, request the number of DOFs through this interface class. This approach also enables the co-existence of various element-based schemes.

Author
Benjamin S. Kirk
Date
2002

Definition at line 93 of file fe_abstract.h.

Member Typedef Documentation

typedef std::map<std::string, std::pair<unsigned int, unsigned int> > libMesh::ReferenceCounter::Counts
protectedinherited

Data structure to log the information. The log is identified by the class name.

Definition at line 119 of file reference_counter.h.

Constructor & Destructor Documentation

libMesh::FEAbstract::FEAbstract ( const unsigned int  dim,
const FEType fet 
)
inlineprotected

Constructor. Optionally initializes required data structures. Protected so that this base class cannot be explicitly instantiated.

Definition at line 608 of file fe_abstract.h.

609  :
610  _fe_map( FEMap::build(fet) ),
611  dim(d),
612  calculations_started(false),
613  calculate_phi(false),
614  calculate_dphi(false),
615  calculate_d2phi(false),
616  calculate_curl_phi(false),
617  calculate_div_phi(false),
618  calculate_dphiref(false),
619  fe_type(fet),
621  _p_level(0),
623  shapes_on_quadrature(false)
624 {
625 }
unsigned int _p_level
Definition: fe_abstract.h:579
static std::unique_ptr< FEMap > build(FEType fe_type)
Definition: fe_map.C:51
const class libmesh_nullptr_t libmesh_nullptr
const unsigned int dim
Definition: fe_abstract.h:523
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
libMesh::FEAbstract::~FEAbstract ( )
inlinevirtual

Destructor.

Definition at line 629 of file fe_abstract.h.

630 {
631 }

Member Function Documentation

std::unique_ptr< FEAbstract > libMesh::FEAbstract::build ( const unsigned int  dim,
const FEType type 
)
static

Builds a specific finite element type.

Returns
A std::unique_ptr<FEAbstract> to the FE object to prevent memory leaks.

Definition at line 44 of file fe_abstract.C.

References libMesh::BERNSTEIN, libMesh::CLOUGH, libMesh::FEType::family, libMesh::HERMITE, libMesh::HIERARCHIC, libMesh::L2_HIERARCHIC, libMesh::L2_LAGRANGE, libMesh::LAGRANGE, libMesh::LAGRANGE_VEC, libMesh::MONOMIAL, libMesh::NEDELEC_ONE, libMesh::SCALAR, libMesh::SUBDIVISION, libMesh::SZABAB, and libMesh::XYZ.

Referenced by libMesh::DGFEMContext::DGFEMContext(), libMesh::FEMContext::init_internal_data(), and libMesh::DofMap::use_coupled_neighbor_dofs().

46 {
47  switch (dim)
48  {
49  // 0D
50  case 0:
51  {
52  switch (fet.family)
53  {
54  case CLOUGH:
55  return libmesh_make_unique<FE<0,CLOUGH>>(fet);
56 
57  case HERMITE:
58  return libmesh_make_unique<FE<0,HERMITE>>(fet);
59 
60  case LAGRANGE:
61  return libmesh_make_unique<FE<0,LAGRANGE>>(fet);
62 
63  case LAGRANGE_VEC:
64  return libmesh_make_unique<FE<0,LAGRANGE_VEC>>(fet);
65 
66  case L2_LAGRANGE:
67  return libmesh_make_unique<FE<0,L2_LAGRANGE>>(fet);
68 
69  case HIERARCHIC:
70  return libmesh_make_unique<FE<0,HIERARCHIC>>(fet);
71 
72  case L2_HIERARCHIC:
73  return libmesh_make_unique<FE<0,L2_HIERARCHIC>>(fet);
74 
75  case MONOMIAL:
76  return libmesh_make_unique<FE<0,MONOMIAL>>(fet);
77 
78 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
79  case SZABAB:
80  return libmesh_make_unique<FE<0,SZABAB>>(fet);
81 
82  case BERNSTEIN:
83  return libmesh_make_unique<FE<0,BERNSTEIN>>(fet);
84 #endif
85 
86  case XYZ:
87  return libmesh_make_unique<FEXYZ<0>>(fet);
88 
89  case SCALAR:
90  return libmesh_make_unique<FEScalar<0>>(fet);
91 
92  default:
93  libmesh_error_msg("ERROR: Bad FEType.family= " << fet.family);
94  }
95  }
96  // 1D
97  case 1:
98  {
99  switch (fet.family)
100  {
101  case CLOUGH:
102  return libmesh_make_unique<FE<1,CLOUGH>>(fet);
103 
104  case HERMITE:
105  return libmesh_make_unique<FE<1,HERMITE>>(fet);
106 
107  case LAGRANGE:
108  return libmesh_make_unique<FE<1,LAGRANGE>>(fet);
109 
110  case LAGRANGE_VEC:
111  return libmesh_make_unique<FE<1,LAGRANGE_VEC>>(fet);
112 
113  case L2_LAGRANGE:
114  return libmesh_make_unique<FE<1,L2_LAGRANGE>>(fet);
115 
116  case HIERARCHIC:
117  return libmesh_make_unique<FE<1,HIERARCHIC>>(fet);
118 
119  case L2_HIERARCHIC:
120  return libmesh_make_unique<FE<1,L2_HIERARCHIC>>(fet);
121 
122  case MONOMIAL:
123  return libmesh_make_unique<FE<1,MONOMIAL>>(fet);
124 
125 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
126  case SZABAB:
127  return libmesh_make_unique<FE<1,SZABAB>>(fet);
128 
129  case BERNSTEIN:
130  return libmesh_make_unique<FE<1,BERNSTEIN>>(fet);
131 #endif
132 
133  case XYZ:
134  return libmesh_make_unique<FEXYZ<1>>(fet);
135 
136  case SCALAR:
137  return libmesh_make_unique<FEScalar<1>>(fet);
138 
139  default:
140  libmesh_error_msg("ERROR: Bad FEType.family= " << fet.family);
141  }
142  }
143 
144 
145  // 2D
146  case 2:
147  {
148  switch (fet.family)
149  {
150  case CLOUGH:
151  return libmesh_make_unique<FE<2,CLOUGH>>(fet);
152 
153  case HERMITE:
154  return libmesh_make_unique<FE<2,HERMITE>>(fet);
155 
156  case LAGRANGE:
157  return libmesh_make_unique<FE<2,LAGRANGE>>(fet);
158 
159  case LAGRANGE_VEC:
160  return libmesh_make_unique<FE<2,LAGRANGE_VEC>>(fet);
161 
162  case L2_LAGRANGE:
163  return libmesh_make_unique<FE<2,L2_LAGRANGE>>(fet);
164 
165  case HIERARCHIC:
166  return libmesh_make_unique<FE<2,HIERARCHIC>>(fet);
167 
168  case L2_HIERARCHIC:
169  return libmesh_make_unique<FE<2,L2_HIERARCHIC>>(fet);
170 
171  case MONOMIAL:
172  return libmesh_make_unique<FE<2,MONOMIAL>>(fet);
173 
174 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
175  case SZABAB:
176  return libmesh_make_unique<FE<2,SZABAB>>(fet);
177 
178  case BERNSTEIN:
179  return libmesh_make_unique<FE<2,BERNSTEIN>>(fet);
180 #endif
181 
182  case XYZ:
183  return libmesh_make_unique<FEXYZ<2>>(fet);
184 
185  case SCALAR:
186  return libmesh_make_unique<FEScalar<2>>(fet);
187 
188  case NEDELEC_ONE:
189  return libmesh_make_unique<FENedelecOne<2>>(fet);
190 
191  case SUBDIVISION:
192  return libmesh_make_unique<FESubdivision>(fet);
193 
194  default:
195  libmesh_error_msg("ERROR: Bad FEType.family= " << fet.family);
196  }
197  }
198 
199 
200  // 3D
201  case 3:
202  {
203  switch (fet.family)
204  {
205  case CLOUGH:
206  libmesh_error_msg("ERROR: Clough-Tocher elements currently only support 1D and 2D");
207 
208  case HERMITE:
209  return libmesh_make_unique<FE<3,HERMITE>>(fet);
210 
211  case LAGRANGE:
212  return libmesh_make_unique<FE<3,LAGRANGE>>(fet);
213 
214  case LAGRANGE_VEC:
215  return libmesh_make_unique<FE<3,LAGRANGE_VEC>>(fet);
216 
217  case L2_LAGRANGE:
218  return libmesh_make_unique<FE<3,L2_LAGRANGE>>(fet);
219 
220  case HIERARCHIC:
221  return libmesh_make_unique<FE<3,HIERARCHIC>>(fet);
222 
223  case L2_HIERARCHIC:
224  return libmesh_make_unique<FE<3,L2_HIERARCHIC>>(fet);
225 
226  case MONOMIAL:
227  return libmesh_make_unique<FE<3,MONOMIAL>>(fet);
228 
229 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
230  case SZABAB:
231  return libmesh_make_unique<FE<3,SZABAB>>(fet);
232 
233  case BERNSTEIN:
234  return libmesh_make_unique<FE<3,BERNSTEIN>>(fet);
235 #endif
236 
237  case XYZ:
238  return libmesh_make_unique<FEXYZ<3>>(fet);
239 
240  case SCALAR:
241  return libmesh_make_unique<FEScalar<3>>(fet);
242 
243  case NEDELEC_ONE:
244  return libmesh_make_unique<FENedelecOne<3>>(fet);
245 
246  default:
247  libmesh_error_msg("ERROR: Bad FEType.family= " << fet.family);
248  }
249  }
250 
251  default:
252  libmesh_error_msg("Invalid dimension dim = " << dim);
253  }
254 }
const unsigned int dim
Definition: fe_abstract.h:523
void libMesh::FEAbstract::compute_node_constraints ( NodeConstraints constraints,
const Elem elem 
)
static

Computes the nodal constraint contributions (for non-conforming adapted meshes), using Lagrange geometry

Definition at line 793 of file fe_abstract.C.

References std::abs(), libMesh::Elem::build_side_ptr(), libMesh::Elem::default_order(), libMesh::Elem::dim(), fe_type, libMesh::FEInterface::inverse_map(), libMesh::LAGRANGE, libMesh::Elem::level(), libmesh_nullptr, libMesh::FEInterface::n_dofs(), libMesh::Elem::neighbor_ptr(), libMesh::Elem::parent(), libMesh::Real, libMesh::remote_elem, libMesh::FEInterface::shape(), libMesh::Elem::side_index_range(), libMesh::Threads::spin_mtx, and libMesh::Elem::subactive().

795 {
796  libmesh_assert(elem);
797 
798  const unsigned int Dim = elem->dim();
799 
800  // Only constrain elements in 2,3D.
801  if (Dim == 1)
802  return;
803 
804  // Only constrain active and ancestor elements
805  if (elem->subactive())
806  return;
807 
808  // We currently always use LAGRANGE mappings for geometry
809  const FEType fe_type(elem->default_order(), LAGRANGE);
810 
811  std::vector<const Node *> my_nodes, parent_nodes;
812 
813  // Look at the element faces. Check to see if we need to
814  // build constraints.
815  for (auto s : elem->side_index_range())
816  if (elem->neighbor_ptr(s) != libmesh_nullptr &&
817  elem->neighbor_ptr(s) != remote_elem)
818  if (elem->neighbor_ptr(s)->level() < elem->level()) // constrain dofs shared between
819  { // this element and ones coarser
820  // than this element.
821  // Get pointers to the elements of interest and its parent.
822  const Elem * parent = elem->parent();
823 
824  // This can't happen... Only level-0 elements have NULL
825  // parents, and no level-0 elements can be at a higher
826  // level than their neighbors!
827  libmesh_assert(parent);
828 
829  const std::unique_ptr<const Elem> my_side (elem->build_side_ptr(s));
830  const std::unique_ptr<const Elem> parent_side (parent->build_side_ptr(s));
831 
832  const unsigned int n_side_nodes = my_side->n_nodes();
833 
834  my_nodes.clear();
835  my_nodes.reserve (n_side_nodes);
836  parent_nodes.clear();
837  parent_nodes.reserve (n_side_nodes);
838 
839  for (unsigned int n=0; n != n_side_nodes; ++n)
840  my_nodes.push_back(my_side->node_ptr(n));
841 
842  for (unsigned int n=0; n != n_side_nodes; ++n)
843  parent_nodes.push_back(parent_side->node_ptr(n));
844 
845  for (unsigned int my_side_n=0;
846  my_side_n < n_side_nodes;
847  my_side_n++)
848  {
849  libmesh_assert_less (my_side_n, FEInterface::n_dofs(Dim-1, fe_type, my_side->type()));
850 
851  const Node * my_node = my_nodes[my_side_n];
852 
853  // The support point of the DOF
854  const Point & support_point = *my_node;
855 
856  // Figure out where my node lies on their reference element.
857  const Point mapped_point = FEInterface::inverse_map(Dim-1, fe_type,
858  parent_side.get(),
859  support_point);
860 
861  // Compute the parent's side shape function values.
862  for (unsigned int their_side_n=0;
863  their_side_n < n_side_nodes;
864  their_side_n++)
865  {
866  libmesh_assert_less (their_side_n, FEInterface::n_dofs(Dim-1, fe_type, parent_side->type()));
867 
868  const Node * their_node = parent_nodes[their_side_n];
869  libmesh_assert(their_node);
870 
871  const Real their_value = FEInterface::shape(Dim-1,
872  fe_type,
873  parent_side->type(),
874  their_side_n,
875  mapped_point);
876 
877  const Real their_mag = std::abs(their_value);
878 #ifdef DEBUG
879  // Protect for the case u_i ~= u_j,
880  // in which case i better equal j.
881  if (their_mag > 0.999)
882  {
883  libmesh_assert_equal_to (my_node, their_node);
884  libmesh_assert_less (std::abs(their_value - 1.), 0.001);
885  }
886  else
887 #endif
888  // To make nodal constraints useful for constructing
889  // sparsity patterns faster, we need to get EVERY
890  // POSSIBLE constraint coupling identified, even if
891  // there is no coupling in the isoparametric
892  // Lagrange case.
893  if (their_mag < 1.e-5)
894  {
895  // since we may be running this method concurrently
896  // on multiple threads we need to acquire a lock
897  // before modifying the shared constraint_row object.
898  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
899 
900  // A reference to the constraint row.
901  NodeConstraintRow & constraint_row = constraints[my_node].first;
902 
903  constraint_row.insert(std::make_pair (their_node,
904  0.));
905  }
906  // To get nodal coordinate constraints right, only
907  // add non-zero and non-identity values for Lagrange
908  // basis functions.
909  else // (1.e-5 <= their_mag <= .999)
910  {
911  // since we may be running this method concurrently
912  // on multiple threads we need to acquire a lock
913  // before modifying the shared constraint_row object.
914  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
915 
916  // A reference to the constraint row.
917  NodeConstraintRow & constraint_row = constraints[my_node].first;
918 
919  constraint_row.insert(std::make_pair (their_node,
920  their_value));
921  }
922  }
923  }
924  }
925 }
double abs(double a)
static unsigned int n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:411
const class libmesh_nullptr_t libmesh_nullptr
spin_mutex spin_mtx
Definition: threads.C:29
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition: fe_interface.C:614
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_interface.C:547
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::map< const Node *, Real, std::less< const Node * >, Threads::scalable_allocator< std::pair< const Node *const, Real > > > NodeConstraintRow
Definition: dof_map.h:137
const RemoteElem * remote_elem
Definition: remote_elem.C:57
void libMesh::FEAbstract::compute_periodic_node_constraints ( NodeConstraints constraints,
const PeriodicBoundaries boundaries,
const MeshBase mesh,
const PointLocatorBase point_locator,
const Elem elem 
)
static

Computes the node position constraint equation contributions (for meshes with periodic boundary conditions)

Definition at line 936 of file fe_abstract.C.

References libMesh::Elem::active(), libMesh::PeriodicBoundaries::boundary(), libMesh::BoundaryInfo::boundary_ids(), libMesh::Elem::build_side_ptr(), libMesh::Elem::default_order(), libMesh::Elem::dim(), fe_type, libMesh::MeshBase::get_boundary_info(), libMesh::PeriodicBoundaryBase::get_corresponding_pos(), libMesh::invalid_uint, libMesh::FEInterface::inverse_map(), libMesh::LAGRANGE, libMesh::Elem::level(), libMesh::FEInterface::n_dofs(), libMesh::PeriodicBoundaries::neighbor(), libMesh::Elem::neighbor_ptr(), libMesh::PeriodicBoundaryBase::pairedboundary, libMesh::Real, libMesh::FEInterface::shape(), libMesh::Elem::side_index_range(), libMesh::BoundaryInfo::side_with_boundary_id(), and libMesh::Threads::spin_mtx.

941 {
942  // Only bother if we truly have periodic boundaries
943  if (boundaries.empty())
944  return;
945 
946  libmesh_assert(elem);
947 
948  // Only constrain active elements with this method
949  if (!elem->active())
950  return;
951 
952  const unsigned int Dim = elem->dim();
953 
954  // We currently always use LAGRANGE mappings for geometry
955  const FEType fe_type(elem->default_order(), LAGRANGE);
956 
957  std::vector<const Node *> my_nodes, neigh_nodes;
958 
959  // Look at the element faces. Check to see if we need to
960  // build constraints.
961  std::vector<boundary_id_type> bc_ids;
962  for (auto s : elem->side_index_range())
963  {
964  if (elem->neighbor_ptr(s))
965  continue;
966 
967  mesh.get_boundary_info().boundary_ids (elem, s, bc_ids);
968  for (std::vector<boundary_id_type>::const_iterator id_it=bc_ids.begin(); id_it!=bc_ids.end(); ++id_it)
969  {
970  const boundary_id_type boundary_id = *id_it;
971  const PeriodicBoundaryBase * periodic = boundaries.boundary(boundary_id);
972  if (periodic)
973  {
974  libmesh_assert(point_locator);
975 
976  // Get pointers to the element's neighbor.
977  const Elem * neigh = boundaries.neighbor(boundary_id, *point_locator, elem, s);
978 
979  // h refinement constraints:
980  // constrain dofs shared between
981  // this element and ones as coarse
982  // as or coarser than this element.
983  if (neigh->level() <= elem->level())
984  {
985  unsigned int s_neigh =
986  mesh.get_boundary_info().side_with_boundary_id(neigh, periodic->pairedboundary);
987  libmesh_assert_not_equal_to (s_neigh, libMesh::invalid_uint);
988 
989 #ifdef LIBMESH_ENABLE_AMR
990  libmesh_assert(neigh->active());
991 #endif // #ifdef LIBMESH_ENABLE_AMR
992 
993  const std::unique_ptr<const Elem> my_side (elem->build_side_ptr(s));
994  const std::unique_ptr<const Elem> neigh_side (neigh->build_side_ptr(s_neigh));
995 
996  const unsigned int n_side_nodes = my_side->n_nodes();
997 
998  my_nodes.clear();
999  my_nodes.reserve (n_side_nodes);
1000  neigh_nodes.clear();
1001  neigh_nodes.reserve (n_side_nodes);
1002 
1003  for (unsigned int n=0; n != n_side_nodes; ++n)
1004  my_nodes.push_back(my_side->node_ptr(n));
1005 
1006  for (unsigned int n=0; n != n_side_nodes; ++n)
1007  neigh_nodes.push_back(neigh_side->node_ptr(n));
1008 
1009  // Make sure we're not adding recursive constraints
1010  // due to the redundancy in the way we add periodic
1011  // boundary constraints, or adding constraints to
1012  // nodes that already have AMR constraints
1013  std::vector<bool> skip_constraint(n_side_nodes, false);
1014 
1015  for (unsigned int my_side_n=0;
1016  my_side_n < n_side_nodes;
1017  my_side_n++)
1018  {
1019  libmesh_assert_less (my_side_n, FEInterface::n_dofs(Dim-1, fe_type, my_side->type()));
1020 
1021  const Node * my_node = my_nodes[my_side_n];
1022 
1023  // Figure out where my node lies on their reference element.
1024  const Point neigh_point = periodic->get_corresponding_pos(*my_node);
1025 
1026  const Point mapped_point = FEInterface::inverse_map(Dim-1, fe_type,
1027  neigh_side.get(),
1028  neigh_point);
1029 
1030  // If we've already got a constraint on this
1031  // node, then the periodic constraint is
1032  // redundant
1033  {
1034  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1035 
1036  if (constraints.count(my_node))
1037  {
1038  skip_constraint[my_side_n] = true;
1039  continue;
1040  }
1041  }
1042 
1043  // Compute the neighbors's side shape function values.
1044  for (unsigned int their_side_n=0;
1045  their_side_n < n_side_nodes;
1046  their_side_n++)
1047  {
1048  libmesh_assert_less (their_side_n, FEInterface::n_dofs(Dim-1, fe_type, neigh_side->type()));
1049 
1050  const Node * their_node = neigh_nodes[their_side_n];
1051 
1052  // If there's a constraint on an opposing node,
1053  // we need to see if it's constrained by
1054  // *our side* making any periodic constraint
1055  // on us recursive
1056  {
1057  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1058 
1059  if (!constraints.count(their_node))
1060  continue;
1061 
1062  const NodeConstraintRow & their_constraint_row =
1063  constraints[their_node].first;
1064 
1065  for (unsigned int orig_side_n=0;
1066  orig_side_n < n_side_nodes;
1067  orig_side_n++)
1068  {
1069  libmesh_assert_less (orig_side_n, FEInterface::n_dofs(Dim-1, fe_type, my_side->type()));
1070 
1071  const Node * orig_node = my_nodes[orig_side_n];
1072 
1073  if (their_constraint_row.count(orig_node))
1074  skip_constraint[orig_side_n] = true;
1075  }
1076  }
1077  }
1078  }
1079  for (unsigned int my_side_n=0;
1080  my_side_n < n_side_nodes;
1081  my_side_n++)
1082  {
1083  libmesh_assert_less (my_side_n, FEInterface::n_dofs(Dim-1, fe_type, my_side->type()));
1084 
1085  if (skip_constraint[my_side_n])
1086  continue;
1087 
1088  const Node * my_node = my_nodes[my_side_n];
1089 
1090  // Figure out where my node lies on their reference element.
1091  const Point neigh_point = periodic->get_corresponding_pos(*my_node);
1092 
1093  // Figure out where my node lies on their reference element.
1094  const Point mapped_point = FEInterface::inverse_map(Dim-1, fe_type,
1095  neigh_side.get(),
1096  neigh_point);
1097 
1098  for (unsigned int their_side_n=0;
1099  their_side_n < n_side_nodes;
1100  their_side_n++)
1101  {
1102  libmesh_assert_less (their_side_n, FEInterface::n_dofs(Dim-1, fe_type, neigh_side->type()));
1103 
1104  const Node * their_node = neigh_nodes[their_side_n];
1105  libmesh_assert(their_node);
1106 
1107  const Real their_value = FEInterface::shape(Dim-1,
1108  fe_type,
1109  neigh_side->type(),
1110  their_side_n,
1111  mapped_point);
1112 
1113  // since we may be running this method concurrently
1114  // on multiple threads we need to acquire a lock
1115  // before modifying the shared constraint_row object.
1116  {
1117  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1118 
1119  NodeConstraintRow & constraint_row =
1120  constraints[my_node].first;
1121 
1122  constraint_row.insert(std::make_pair(their_node,
1123  their_value));
1124  }
1125  }
1126  }
1127  }
1128  }
1129  }
1130  }
1131 }
static unsigned int n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:411
const unsigned int invalid_uint
Definition: libmesh.h:184
MeshBase & mesh
spin_mutex spin_mtx
Definition: threads.C:29
int8_t boundary_id_type
Definition: id_types.h:51
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition: fe_interface.C:614
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_interface.C:547
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::map< const Node *, Real, std::less< const Node * >, Threads::scalable_allocator< std::pair< const Node *const, Real > > > NodeConstraintRow
Definition: dof_map.h:137
virtual void libMesh::FEAbstract::compute_shape_functions ( const Elem ,
const std::vector< Point > &   
)
protectedpure virtual

After having updated the jacobian and the transformation from local to global coordinates in FEMap::compute_map(), the first derivatives of the shape functions are transformed to global coordinates, giving dphi, dphidx, dphidy, and dphidz. This method should rarely be re-defined in derived classes, but still should be usable for children. Therefore, keep it protected. This needs to be implemented in the derived class since this function depends on whether the shape functions are vector-valued or not.

Implemented in libMesh::FEXYZ< Dim >, libMesh::InfFE< Dim, T_radial, T_map >, libMesh::FEGenericBase< OutputType >, and libMesh::FEGenericBase< FEOutputType< T >::type >.

Referenced by get_fe_map().

void libMesh::ReferenceCounter::disable_print_counter_info ( )
staticinherited
virtual void libMesh::FEAbstract::edge_reinit ( const Elem elem,
const unsigned int  edge,
const Real  tolerance = TOLERANCE,
const std::vector< Point > *  pts = libmesh_nullptr,
const std::vector< Real > *  weights = libmesh_nullptr 
)
pure virtual

Reinitializes all the physical element-dependent data based on the edge of the element elem. The tolerance parameter is passed to the involved call to inverse_map(). By default the element data are computed at the quadrature points specified by the quadrature rule qrule, but any set of points on the reference edge element may be specified in the optional argument pts.

Implemented in libMesh::InfFE< Dim, T_radial, T_map >, libMesh::FE< Dim, T >, libMesh::FE< 2, SUBDIVISION >, libMesh::FE< Dim, HIERARCHIC >, libMesh::FE< Dim, SCALAR >, libMesh::FE< Dim, L2_LAGRANGE >, libMesh::FE< Dim, NEDELEC_ONE >, libMesh::FE< Dim, HERMITE >, libMesh::FE< Dim, CLOUGH >, libMesh::FE< Dim, MONOMIAL >, libMesh::FE< Dim, XYZ >, libMesh::FE< Dim, LAGRANGE >, libMesh::FE< Dim, L2_HIERARCHIC >, and libMesh::FE< Dim, LAGRANGE_VEC >.

void libMesh::ReferenceCounter::enable_print_counter_info ( )
staticinherited

Methods to enable/disable the reference counter output from print_info()

Definition at line 101 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter.

Referenced by libMesh::ReferenceCounter::n_objects().

102 {
103  _enable_print_counter = true;
104  return;
105 }
virtual FEContinuity libMesh::FEAbstract::get_continuity ( ) const
pure virtual
Returns
The continuity level of the finite element.

Implemented in libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::InfFE< Dim, T_radial, T_map >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< 2, SUBDIVISION >, libMesh::FE< Dim, HIERARCHIC >, libMesh::FE< Dim, SCALAR >, libMesh::FE< Dim, L2_LAGRANGE >, libMesh::FE< Dim, NEDELEC_ONE >, libMesh::FE< Dim, HERMITE >, libMesh::FE< Dim, CLOUGH >, libMesh::FE< Dim, MONOMIAL >, libMesh::FE< Dim, XYZ >, libMesh::FE< Dim, LAGRANGE >, libMesh::FE< Dim, L2_HIERARCHIC >, libMesh::FE< Dim, LAGRANGE_VEC >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, and libMesh::FE< Dim, T >.

Referenced by libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::operator()(), and set_fe_order().

const std::vector<Real>& libMesh::FEAbstract::get_curvatures ( ) const
inline
Returns
The curvatures for use in face integration.

Definition at line 383 of file fe_abstract.h.

References _fe_map, attach_quadrature_rule(), n_quadrature_points(), and n_shape_functions().

384  { return this->_fe_map->get_curvatures();}
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<RealGradient>& libMesh::FEAbstract::get_d2xyzdeta2 ( ) const
inline
Returns
The second partial derivatives in eta.

Definition at line 270 of file fe_abstract.h.

References _fe_map.

271  { return this->_fe_map->get_d2xyzdeta2(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<RealGradient>& libMesh::FEAbstract::get_d2xyzdetadzeta ( ) const
inline
Returns
The second partial derivatives in eta-zeta.

Definition at line 300 of file fe_abstract.h.

References _fe_map.

301  { return this->_fe_map->get_d2xyzdetadzeta(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<RealGradient>& libMesh::FEAbstract::get_d2xyzdxi2 ( ) const
inline
Returns
The second partial derivatives in xi.

Definition at line 264 of file fe_abstract.h.

References _fe_map.

265  { return this->_fe_map->get_d2xyzdxi2(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<RealGradient>& libMesh::FEAbstract::get_d2xyzdxideta ( ) const
inline
Returns
The second partial derivatives in xi-eta.

Definition at line 286 of file fe_abstract.h.

References _fe_map.

287  { return this->_fe_map->get_d2xyzdxideta(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<RealGradient>& libMesh::FEAbstract::get_d2xyzdxidzeta ( ) const
inline
Returns
The second partial derivatives in xi-zeta.

Definition at line 294 of file fe_abstract.h.

References _fe_map.

295  { return this->_fe_map->get_d2xyzdxidzeta(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<RealGradient>& libMesh::FEAbstract::get_d2xyzdzeta2 ( ) const
inline
Returns
The second partial derivatives in zeta.

Definition at line 278 of file fe_abstract.h.

References _fe_map.

279  { return this->_fe_map->get_d2xyzdzeta2(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<Real>& libMesh::FEAbstract::get_detadx ( ) const
inline
Returns
The deta/dx entry in the transformation matrix from physical to local coordinates.

Definition at line 330 of file fe_abstract.h.

References _fe_map.

331  { return this->_fe_map->get_detadx(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<Real>& libMesh::FEAbstract::get_detady ( ) const
inline
Returns
The deta/dy entry in the transformation matrix from physical to local coordinates.

Definition at line 337 of file fe_abstract.h.

References _fe_map.

338  { return this->_fe_map->get_detady(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<Real>& libMesh::FEAbstract::get_detadz ( ) const
inline
Returns
The deta/dz entry in the transformation matrix from physical to local coordinates.

Definition at line 344 of file fe_abstract.h.

References _fe_map.

345  { return this->_fe_map->get_detadz(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
unsigned int libMesh::FEAbstract::get_dim ( ) const
inline
Returns
the dimension of this FE

Definition at line 223 of file fe_abstract.h.

References dim.

224  { return dim; }
const unsigned int dim
Definition: fe_abstract.h:523
const std::vector<Real>& libMesh::FEAbstract::get_dxidx ( ) const
inline
Returns
The dxi/dx entry in the transformation matrix from physical to local coordinates.

Definition at line 309 of file fe_abstract.h.

References _fe_map.

310  { return this->_fe_map->get_dxidx(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<Real>& libMesh::FEAbstract::get_dxidy ( ) const
inline
Returns
The dxi/dy entry in the transformation matrix from physical to local coordinates.

Definition at line 316 of file fe_abstract.h.

References _fe_map.

317  { return this->_fe_map->get_dxidy(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<Real>& libMesh::FEAbstract::get_dxidz ( ) const
inline
Returns
The dxi/dz entry in the transformation matrix from physical to local coordinates.

Definition at line 323 of file fe_abstract.h.

References _fe_map.

324  { return this->_fe_map->get_dxidz(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<RealGradient>& libMesh::FEAbstract::get_dxyzdeta ( ) const
inline
Returns
The element tangents in eta-direction at the quadrature points.

Definition at line 251 of file fe_abstract.h.

References _fe_map.

252  { return this->_fe_map->get_dxyzdeta(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<RealGradient>& libMesh::FEAbstract::get_dxyzdxi ( ) const
inline
Returns
The element tangents in xi-direction at the quadrature points.

Definition at line 244 of file fe_abstract.h.

References _fe_map.

245  { return this->_fe_map->get_dxyzdxi(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<RealGradient>& libMesh::FEAbstract::get_dxyzdzeta ( ) const
inline
Returns
The element tangents in zeta-direction at the quadrature points.

Definition at line 258 of file fe_abstract.h.

References _fe_map.

259  { return _fe_map->get_dxyzdzeta(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<Real>& libMesh::FEAbstract::get_dzetadx ( ) const
inline
Returns
The dzeta/dx entry in the transformation matrix from physical to local coordinates.

Definition at line 351 of file fe_abstract.h.

References _fe_map.

352  { return this->_fe_map->get_dzetadx(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<Real>& libMesh::FEAbstract::get_dzetady ( ) const
inline
Returns
The dzeta/dy entry in the transformation matrix from physical to local coordinates.

Definition at line 358 of file fe_abstract.h.

References _fe_map.

359  { return this->_fe_map->get_dzetady(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
const std::vector<Real>& libMesh::FEAbstract::get_dzetadz ( ) const
inline
Returns
The dzeta/dz entry in the transformation matrix from physical to local coordinates.

Definition at line 365 of file fe_abstract.h.

References _fe_map.

366  { return this->_fe_map->get_dzetadz(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
FEFamily libMesh::FEAbstract::get_family ( ) const
inline
Returns
The finite element family of this element.

Definition at line 447 of file fe_abstract.h.

References libMesh::FEType::family, and fe_type.

Referenced by libMesh::FE< Dim, T >::FE().

447 { return fe_type.family; }
FEFamily family
Definition: fe_type.h:204
FEType libMesh::FEAbstract::get_fe_type ( ) const
inline
std::string libMesh::ReferenceCounter::get_info ( )
staticinherited

Gets a string containing the reference information.

Definition at line 47 of file reference_counter.C.

References libMesh::ReferenceCounter::_counts, and libMesh::Quality::name().

Referenced by libMesh::ReferenceCounter::print_info().

48 {
49 #if defined(LIBMESH_ENABLE_REFERENCE_COUNTING) && defined(DEBUG)
50 
51  std::ostringstream oss;
52 
53  oss << '\n'
54  << " ---------------------------------------------------------------------------- \n"
55  << "| Reference count information |\n"
56  << " ---------------------------------------------------------------------------- \n";
57 
58  for (Counts::iterator it = _counts.begin();
59  it != _counts.end(); ++it)
60  {
61  const std::string name(it->first);
62  const unsigned int creations = it->second.first;
63  const unsigned int destructions = it->second.second;
64 
65  oss << "| " << name << " reference count information:\n"
66  << "| Creations: " << creations << '\n'
67  << "| Destructions: " << destructions << '\n';
68  }
69 
70  oss << " ---------------------------------------------------------------------------- \n";
71 
72  return oss.str();
73 
74 #else
75 
76  return "";
77 
78 #endif
79 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
const std::vector<Point>& libMesh::FEAbstract::get_normals ( ) const
inline
Order libMesh::FEAbstract::get_order ( ) const
inline
Returns
The approximation order of the finite element.

Definition at line 426 of file fe_abstract.h.

References _p_level, fe_type, and libMesh::FEType::order.

426 { return static_cast<Order>(fe_type.order + _p_level); }
unsigned int _p_level
Definition: fe_abstract.h:579
OrderWrapper order
Definition: fe_type.h:198
unsigned int libMesh::FEAbstract::get_p_level ( ) const
inline
Returns
The p refinement level that the current shape functions have been calculated for.

Definition at line 416 of file fe_abstract.h.

References _p_level.

416 { return _p_level; }
unsigned int _p_level
Definition: fe_abstract.h:579
void libMesh::FEAbstract::get_refspace_nodes ( const ElemType  t,
std::vector< Point > &  nodes 
)
static
Returns
The reference space coordinates of nodes based on the element type.

Definition at line 256 of file fe_abstract.C.

References libMesh::EDGE2, libMesh::EDGE3, libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::PRISM15, libMesh::PRISM18, libMesh::PRISM6, libMesh::PYRAMID13, libMesh::PYRAMID14, libMesh::PYRAMID5, libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::QUADSHELL4, libMesh::QUADSHELL8, libMesh::TET10, libMesh::TET4, libMesh::TRI3, libMesh::TRI6, and libMesh::TRISHELL3.

257 {
258  switch(itemType)
259  {
260  case EDGE2:
261  {
262  nodes.resize(2);
263  nodes[0] = Point (-1.,0.,0.);
264  nodes[1] = Point (1.,0.,0.);
265  return;
266  }
267  case EDGE3:
268  {
269  nodes.resize(3);
270  nodes[0] = Point (-1.,0.,0.);
271  nodes[1] = Point (1.,0.,0.);
272  nodes[2] = Point (0.,0.,0.);
273  return;
274  }
275  case TRI3:
276  case TRISHELL3:
277  {
278  nodes.resize(3);
279  nodes[0] = Point (0.,0.,0.);
280  nodes[1] = Point (1.,0.,0.);
281  nodes[2] = Point (0.,1.,0.);
282  return;
283  }
284  case TRI6:
285  {
286  nodes.resize(6);
287  nodes[0] = Point (0.,0.,0.);
288  nodes[1] = Point (1.,0.,0.);
289  nodes[2] = Point (0.,1.,0.);
290  nodes[3] = Point (.5,0.,0.);
291  nodes[4] = Point (.5,.5,0.);
292  nodes[5] = Point (0.,.5,0.);
293  return;
294  }
295  case QUAD4:
296  case QUADSHELL4:
297  {
298  nodes.resize(4);
299  nodes[0] = Point (-1.,-1.,0.);
300  nodes[1] = Point (1.,-1.,0.);
301  nodes[2] = Point (1.,1.,0.);
302  nodes[3] = Point (-1.,1.,0.);
303  return;
304  }
305  case QUAD8:
306  case QUADSHELL8:
307  {
308  nodes.resize(8);
309  nodes[0] = Point (-1.,-1.,0.);
310  nodes[1] = Point (1.,-1.,0.);
311  nodes[2] = Point (1.,1.,0.);
312  nodes[3] = Point (-1.,1.,0.);
313  nodes[4] = Point (0.,-1.,0.);
314  nodes[5] = Point (1.,0.,0.);
315  nodes[6] = Point (0.,1.,0.);
316  nodes[7] = Point (-1.,0.,0.);
317  return;
318  }
319  case QUAD9:
320  {
321  nodes.resize(9);
322  nodes[0] = Point (-1.,-1.,0.);
323  nodes[1] = Point (1.,-1.,0.);
324  nodes[2] = Point (1.,1.,0.);
325  nodes[3] = Point (-1.,1.,0.);
326  nodes[4] = Point (0.,-1.,0.);
327  nodes[5] = Point (1.,0.,0.);
328  nodes[6] = Point (0.,1.,0.);
329  nodes[7] = Point (-1.,0.,0.);
330  nodes[8] = Point (0.,0.,0.);
331  return;
332  }
333  case TET4:
334  {
335  nodes.resize(4);
336  nodes[0] = Point (0.,0.,0.);
337  nodes[1] = Point (1.,0.,0.);
338  nodes[2] = Point (0.,1.,0.);
339  nodes[3] = Point (0.,0.,1.);
340  return;
341  }
342  case TET10:
343  {
344  nodes.resize(10);
345  nodes[0] = Point (0.,0.,0.);
346  nodes[1] = Point (1.,0.,0.);
347  nodes[2] = Point (0.,1.,0.);
348  nodes[3] = Point (0.,0.,1.);
349  nodes[4] = Point (.5,0.,0.);
350  nodes[5] = Point (.5,.5,0.);
351  nodes[6] = Point (0.,.5,0.);
352  nodes[7] = Point (0.,0.,.5);
353  nodes[8] = Point (.5,0.,.5);
354  nodes[9] = Point (0.,.5,.5);
355  return;
356  }
357  case HEX8:
358  {
359  nodes.resize(8);
360  nodes[0] = Point (-1.,-1.,-1.);
361  nodes[1] = Point (1.,-1.,-1.);
362  nodes[2] = Point (1.,1.,-1.);
363  nodes[3] = Point (-1.,1.,-1.);
364  nodes[4] = Point (-1.,-1.,1.);
365  nodes[5] = Point (1.,-1.,1.);
366  nodes[6] = Point (1.,1.,1.);
367  nodes[7] = Point (-1.,1.,1.);
368  return;
369  }
370  case HEX20:
371  {
372  nodes.resize(20);
373  nodes[0] = Point (-1.,-1.,-1.);
374  nodes[1] = Point (1.,-1.,-1.);
375  nodes[2] = Point (1.,1.,-1.);
376  nodes[3] = Point (-1.,1.,-1.);
377  nodes[4] = Point (-1.,-1.,1.);
378  nodes[5] = Point (1.,-1.,1.);
379  nodes[6] = Point (1.,1.,1.);
380  nodes[7] = Point (-1.,1.,1.);
381  nodes[8] = Point (0.,-1.,-1.);
382  nodes[9] = Point (1.,0.,-1.);
383  nodes[10] = Point (0.,1.,-1.);
384  nodes[11] = Point (-1.,0.,-1.);
385  nodes[12] = Point (-1.,-1.,0.);
386  nodes[13] = Point (1.,-1.,0.);
387  nodes[14] = Point (1.,1.,0.);
388  nodes[15] = Point (-1.,1.,0.);
389  nodes[16] = Point (0.,-1.,1.);
390  nodes[17] = Point (1.,0.,1.);
391  nodes[18] = Point (0.,1.,1.);
392  nodes[19] = Point (-1.,0.,1.);
393  return;
394  }
395  case HEX27:
396  {
397  nodes.resize(27);
398  nodes[0] = Point (-1.,-1.,-1.);
399  nodes[1] = Point (1.,-1.,-1.);
400  nodes[2] = Point (1.,1.,-1.);
401  nodes[3] = Point (-1.,1.,-1.);
402  nodes[4] = Point (-1.,-1.,1.);
403  nodes[5] = Point (1.,-1.,1.);
404  nodes[6] = Point (1.,1.,1.);
405  nodes[7] = Point (-1.,1.,1.);
406  nodes[8] = Point (0.,-1.,-1.);
407  nodes[9] = Point (1.,0.,-1.);
408  nodes[10] = Point (0.,1.,-1.);
409  nodes[11] = Point (-1.,0.,-1.);
410  nodes[12] = Point (-1.,-1.,0.);
411  nodes[13] = Point (1.,-1.,0.);
412  nodes[14] = Point (1.,1.,0.);
413  nodes[15] = Point (-1.,1.,0.);
414  nodes[16] = Point (0.,-1.,1.);
415  nodes[17] = Point (1.,0.,1.);
416  nodes[18] = Point (0.,1.,1.);
417  nodes[19] = Point (-1.,0.,1.);
418  nodes[20] = Point (0.,0.,-1.);
419  nodes[21] = Point (0.,-1.,0.);
420  nodes[22] = Point (1.,0.,0.);
421  nodes[23] = Point (0.,1.,0.);
422  nodes[24] = Point (-1.,0.,0.);
423  nodes[25] = Point (0.,0.,1.);
424  nodes[26] = Point (0.,0.,0.);
425  return;
426  }
427  case PRISM6:
428  {
429  nodes.resize(6);
430  nodes[0] = Point (0.,0.,-1.);
431  nodes[1] = Point (1.,0.,-1.);
432  nodes[2] = Point (0.,1.,-1.);
433  nodes[3] = Point (0.,0.,1.);
434  nodes[4] = Point (1.,0.,1.);
435  nodes[5] = Point (0.,1.,1.);
436  return;
437  }
438  case PRISM15:
439  {
440  nodes.resize(15);
441  nodes[0] = Point (0.,0.,-1.);
442  nodes[1] = Point (1.,0.,-1.);
443  nodes[2] = Point (0.,1.,-1.);
444  nodes[3] = Point (0.,0.,1.);
445  nodes[4] = Point (1.,0.,1.);
446  nodes[5] = Point (0.,1.,1.);
447  nodes[6] = Point (.5,0.,-1.);
448  nodes[7] = Point (.5,.5,-1.);
449  nodes[8] = Point (0.,.5,-1.);
450  nodes[9] = Point (0.,0.,0.);
451  nodes[10] = Point (1.,0.,0.);
452  nodes[11] = Point (0.,1.,0.);
453  nodes[12] = Point (.5,0.,1.);
454  nodes[13] = Point (.5,.5,1.);
455  nodes[14] = Point (0.,.5,1.);
456  return;
457  }
458  case PRISM18:
459  {
460  nodes.resize(18);
461  nodes[0] = Point (0.,0.,-1.);
462  nodes[1] = Point (1.,0.,-1.);
463  nodes[2] = Point (0.,1.,-1.);
464  nodes[3] = Point (0.,0.,1.);
465  nodes[4] = Point (1.,0.,1.);
466  nodes[5] = Point (0.,1.,1.);
467  nodes[6] = Point (.5,0.,-1.);
468  nodes[7] = Point (.5,.5,-1.);
469  nodes[8] = Point (0.,.5,-1.);
470  nodes[9] = Point (0.,0.,0.);
471  nodes[10] = Point (1.,0.,0.);
472  nodes[11] = Point (0.,1.,0.);
473  nodes[12] = Point (.5,0.,1.);
474  nodes[13] = Point (.5,.5,1.);
475  nodes[14] = Point (0.,.5,1.);
476  nodes[15] = Point (.5,0.,0.);
477  nodes[16] = Point (.5,.5,0.);
478  nodes[17] = Point (0.,.5,0.);
479  return;
480  }
481  case PYRAMID5:
482  {
483  nodes.resize(5);
484  nodes[0] = Point (-1.,-1.,0.);
485  nodes[1] = Point (1.,-1.,0.);
486  nodes[2] = Point (1.,1.,0.);
487  nodes[3] = Point (-1.,1.,0.);
488  nodes[4] = Point (0.,0.,1.);
489  return;
490  }
491  case PYRAMID13:
492  {
493  nodes.resize(13);
494 
495  // base corners
496  nodes[0] = Point (-1.,-1.,0.);
497  nodes[1] = Point (1.,-1.,0.);
498  nodes[2] = Point (1.,1.,0.);
499  nodes[3] = Point (-1.,1.,0.);
500 
501  // apex
502  nodes[4] = Point (0.,0.,1.);
503 
504  // base midedge
505  nodes[5] = Point (0.,-1.,0.);
506  nodes[6] = Point (1.,0.,0.);
507  nodes[7] = Point (0.,1.,0.);
508  nodes[8] = Point (-1,0.,0.);
509 
510  // lateral midedge
511  nodes[9] = Point (-.5,-.5,.5);
512  nodes[10] = Point (.5,-.5,.5);
513  nodes[11] = Point (.5,.5,.5);
514  nodes[12] = Point (-.5,.5,.5);
515 
516  return;
517  }
518  case PYRAMID14:
519  {
520  nodes.resize(14);
521 
522  // base corners
523  nodes[0] = Point (-1.,-1.,0.);
524  nodes[1] = Point (1.,-1.,0.);
525  nodes[2] = Point (1.,1.,0.);
526  nodes[3] = Point (-1.,1.,0.);
527 
528  // apex
529  nodes[4] = Point (0.,0.,1.);
530 
531  // base midedge
532  nodes[5] = Point (0.,-1.,0.);
533  nodes[6] = Point (1.,0.,0.);
534  nodes[7] = Point (0.,1.,0.);
535  nodes[8] = Point (-1,0.,0.);
536 
537  // lateral midedge
538  nodes[9] = Point (-.5,-.5,.5);
539  nodes[10] = Point (.5,-.5,.5);
540  nodes[11] = Point (.5,.5,.5);
541  nodes[12] = Point (-.5,.5,.5);
542 
543  // base center
544  nodes[13] = Point (0.,0.,0.);
545 
546  return;
547  }
548 
549  default:
550  libmesh_error_msg("ERROR: Unknown element type " << itemType);
551  }
552 }
const std::vector<std::vector<Point> >& libMesh::FEAbstract::get_tangents ( ) const
inline
Returns
The tangent vectors for face integration.

Definition at line 371 of file fe_abstract.h.

References _fe_map.

372  { return this->_fe_map->get_tangents(); }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
ElemType libMesh::FEAbstract::get_type ( ) const
inline
Returns
The element type that the current shape functions have been calculated for. Useful in determining when shape functions must be recomputed.

Definition at line 410 of file fe_abstract.h.

References elem_type.

Referenced by libMesh::InfFE< Dim, T_radial, T_map >::reinit().

410 { return elem_type; }
void libMesh::ReferenceCounter::increment_constructor_count ( const std::string &  name)
inlineprotectedinherited

Increments the construction counter. Should be called in the constructor of any derived class that will be reference counted.

Definition at line 185 of file reference_counter.h.

References libMesh::ReferenceCounter::_counts, libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCounter::n_objects(), and libMesh::ReferenceCountedObject< RBParametrized >::ReferenceCountedObject().

186 {
187  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
188  std::pair<unsigned int, unsigned int> & p = _counts[name];
189 
190  p.first++;
191 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
spin_mutex spin_mtx
Definition: threads.C:29
void libMesh::ReferenceCounter::increment_destructor_count ( const std::string &  name)
inlineprotectedinherited

Increments the destruction counter. Should be called in the destructor of any derived class that will be reference counted.

Definition at line 198 of file reference_counter.h.

References libMesh::ReferenceCounter::_counts, libMesh::Quality::name(), and libMesh::Threads::spin_mtx.

Referenced by libMesh::ReferenceCounter::n_objects(), and libMesh::ReferenceCountedObject< RBParametrized >::~ReferenceCountedObject().

199 {
200  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
201  std::pair<unsigned int, unsigned int> & p = _counts[name];
202 
203  p.second++;
204 }
std::string name(const ElemQuality q)
Definition: elem_quality.C:39
spin_mutex spin_mtx
Definition: threads.C:29
virtual bool libMesh::FEAbstract::is_hierarchic ( ) const
pure virtual
Returns
true if the finite element's higher order shape functions are hierarchic

Implemented in libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::InfFE< Dim, T_radial, T_map >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< 2, SUBDIVISION >, libMesh::FE< Dim, HIERARCHIC >, libMesh::FE< Dim, SCALAR >, libMesh::FE< Dim, L2_LAGRANGE >, libMesh::FE< Dim, NEDELEC_ONE >, libMesh::FE< Dim, HERMITE >, libMesh::FE< Dim, CLOUGH >, libMesh::FE< Dim, MONOMIAL >, libMesh::FE< Dim, XYZ >, libMesh::FE< Dim, LAGRANGE >, libMesh::FE< Dim, L2_HIERARCHIC >, libMesh::FE< Dim, LAGRANGE_VEC >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, and libMesh::FE< Dim, T >.

Referenced by set_fe_order().

static unsigned int libMesh::ReferenceCounter::n_objects ( )
inlinestaticinherited
virtual unsigned int libMesh::FEAbstract::n_shape_functions ( ) const
pure virtual
bool libMesh::FEAbstract::on_reference_element ( const Point p,
const ElemType  t,
const Real  eps = TOLERANCE 
)
static
Returns
true if the point p is located on the reference element for element type t, false otherwise. Since we are doing floating point comparisons here the parameter eps can be specified to indicate a tolerance. For example, $ x \le 1 $ becomes $ x \le 1 + \epsilon $.

Definition at line 554 of file fe_abstract.C.

References libMesh::EDGE2, libMesh::EDGE3, libMesh::EDGE4, libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::INFHEX16, libMesh::INFHEX18, libMesh::INFHEX8, libMesh::INFPRISM12, libMesh::INFPRISM6, libMesh::NODEELEM, libMesh::PRISM15, libMesh::PRISM18, libMesh::PRISM6, libMesh::PYRAMID13, libMesh::PYRAMID14, libMesh::PYRAMID5, libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::QUADSHELL4, libMesh::QUADSHELL8, libMesh::Real, libMesh::TET10, libMesh::TET4, libMesh::TRI3, libMesh::TRI6, and libMesh::TRISHELL3.

Referenced by libMesh::FEInterface::ifem_on_reference_element(), libMesh::FE< Dim, T >::inverse_map(), and libMesh::FEInterface::on_reference_element().

555 {
556  libmesh_assert_greater_equal (eps, 0.);
557 
558  const Real xi = p(0);
559 #if LIBMESH_DIM > 1
560  const Real eta = p(1);
561 #else
562  const Real eta = 0.;
563 #endif
564 #if LIBMESH_DIM > 2
565  const Real zeta = p(2);
566 #else
567  const Real zeta = 0.;
568 #endif
569 
570  switch (t)
571  {
572  case NODEELEM:
573  {
574  return (!xi && !eta && !zeta);
575  }
576  case EDGE2:
577  case EDGE3:
578  case EDGE4:
579  {
580  // The reference 1D element is [-1,1].
581  if ((xi >= -1.-eps) &&
582  (xi <= 1.+eps))
583  return true;
584 
585  return false;
586  }
587 
588 
589  case TRI3:
590  case TRISHELL3:
591  case TRI6:
592  {
593  // The reference triangle is isosceles
594  // and is bound by xi=0, eta=0, and xi+eta=1.
595  if ((xi >= 0.-eps) &&
596  (eta >= 0.-eps) &&
597  ((xi + eta) <= 1.+eps))
598  return true;
599 
600  return false;
601  }
602 
603 
604  case QUAD4:
605  case QUADSHELL4:
606  case QUAD8:
607  case QUADSHELL8:
608  case QUAD9:
609  {
610  // The reference quadrilateral element is [-1,1]^2.
611  if ((xi >= -1.-eps) &&
612  (xi <= 1.+eps) &&
613  (eta >= -1.-eps) &&
614  (eta <= 1.+eps))
615  return true;
616 
617  return false;
618  }
619 
620 
621  case TET4:
622  case TET10:
623  {
624  // The reference tetrahedral is isosceles
625  // and is bound by xi=0, eta=0, zeta=0,
626  // and xi+eta+zeta=1.
627  if ((xi >= 0.-eps) &&
628  (eta >= 0.-eps) &&
629  (zeta >= 0.-eps) &&
630  ((xi + eta + zeta) <= 1.+eps))
631  return true;
632 
633  return false;
634  }
635 
636 
637  case HEX8:
638  case HEX20:
639  case HEX27:
640  {
641  /*
642  if ((xi >= -1.) &&
643  (xi <= 1.) &&
644  (eta >= -1.) &&
645  (eta <= 1.) &&
646  (zeta >= -1.) &&
647  (zeta <= 1.))
648  return true;
649  */
650 
651  // The reference hexahedral element is [-1,1]^3.
652  if ((xi >= -1.-eps) &&
653  (xi <= 1.+eps) &&
654  (eta >= -1.-eps) &&
655  (eta <= 1.+eps) &&
656  (zeta >= -1.-eps) &&
657  (zeta <= 1.+eps))
658  {
659  // libMesh::out << "Strange Point:\n";
660  // p.print();
661  return true;
662  }
663 
664  return false;
665  }
666 
667  case PRISM6:
668  case PRISM15:
669  case PRISM18:
670  {
671  // Figure this one out...
672  // inside the reference triangle with zeta in [-1,1]
673  if ((xi >= 0.-eps) &&
674  (eta >= 0.-eps) &&
675  (zeta >= -1.-eps) &&
676  (zeta <= 1.+eps) &&
677  ((xi + eta) <= 1.+eps))
678  return true;
679 
680  return false;
681  }
682 
683 
684  case PYRAMID5:
685  case PYRAMID13:
686  case PYRAMID14:
687  {
688  // Check that the point is on the same side of all the faces
689  // by testing whether:
690  //
691  // n_i.(x - x_i) <= 0
692  //
693  // for each i, where:
694  // n_i is the outward normal of face i,
695  // x_i is a point on face i.
696  if ((-eta - 1. + zeta <= 0.+eps) &&
697  ( xi - 1. + zeta <= 0.+eps) &&
698  ( eta - 1. + zeta <= 0.+eps) &&
699  ( -xi - 1. + zeta <= 0.+eps) &&
700  ( zeta >= 0.-eps))
701  return true;
702 
703  return false;
704  }
705 
706 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
707  case INFHEX8:
708  case INFHEX16:
709  case INFHEX18:
710  {
711  // The reference infhex8 is a [-1,1]^3.
712  if ((xi >= -1.-eps) &&
713  (xi <= 1.+eps) &&
714  (eta >= -1.-eps) &&
715  (eta <= 1.+eps) &&
716  (zeta >= -1.-eps) &&
717  (zeta <= 1.+eps))
718  {
719  return true;
720  }
721  return false;
722  }
723 
724  case INFPRISM6:
725  case INFPRISM12:
726  {
727  // inside the reference triangle with zeta in [-1,1]
728  if ((xi >= 0.-eps) &&
729  (eta >= 0.-eps) &&
730  (zeta >= -1.-eps) &&
731  (zeta <= 1.+eps) &&
732  ((xi + eta) <= 1.+eps))
733  {
734  return true;
735  }
736 
737  return false;
738  }
739 #endif
740 
741  default:
742  libmesh_error_msg("ERROR: Unknown element type " << t);
743  }
744 
745  // If we get here then the point is _not_ in the
746  // reference element. Better return false.
747 
748  return false;
749 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void libMesh::FEAbstract::print_d2phi ( std::ostream &  os) const
pure virtual

Prints the value of each shape function's second derivatives at each quadrature point. Implement in derived class since this depends on whether the element is vector-valued or not.

Implemented in libMesh::FEGenericBase< OutputType >, and libMesh::FEGenericBase< FEOutputType< T >::type >.

Referenced by get_fe_map().

virtual void libMesh::FEAbstract::print_dphi ( std::ostream &  os) const
pure virtual

Prints the value of each shape function's derivative at each quadrature point. Implement in derived class since this depends on whether the element is vector-valued or not.

Implemented in libMesh::FEGenericBase< OutputType >, and libMesh::FEGenericBase< FEOutputType< T >::type >.

Referenced by get_fe_map(), and print_info().

void libMesh::ReferenceCounter::print_info ( std::ostream &  out = libMesh::out)
staticinherited

Prints the reference information, by default to libMesh::out.

Definition at line 88 of file reference_counter.C.

References libMesh::ReferenceCounter::_enable_print_counter, and libMesh::ReferenceCounter::get_info().

Referenced by libMesh::LibMeshInit::LibMeshInit().

89 {
91  out_stream << ReferenceCounter::get_info();
92 }
static std::string get_info()
void libMesh::FEAbstract::print_info ( std::ostream &  os) const

Prints all the relevant information about the current element.

Definition at line 766 of file fe_abstract.C.

References print_dphi(), print_JxW(), print_phi(), and print_xyz().

Referenced by get_fe_map(), and libMesh::operator<<().

767 {
768  os << "phi[i][j]: Shape function i at quadrature pt. j" << std::endl;
769  this->print_phi(os);
770 
771  os << "dphi[i][j]: Shape function i's gradient at quadrature pt. j" << std::endl;
772  this->print_dphi(os);
773 
774  os << "XYZ locations of the quadrature pts." << std::endl;
775  this->print_xyz(os);
776 
777  os << "Values of JxW at the quadrature pts." << std::endl;
778  this->print_JxW(os);
779 }
void print_JxW(std::ostream &os) const
Definition: fe_abstract.C:753
void print_xyz(std::ostream &os) const
Definition: fe_abstract.C:760
virtual void print_phi(std::ostream &os) const =0
virtual void print_dphi(std::ostream &os) const =0
void libMesh::FEAbstract::print_JxW ( std::ostream &  os) const

Prints the Jacobian times the weight for each quadrature point.

Definition at line 753 of file fe_abstract.C.

References _fe_map.

Referenced by get_fe_map(), and print_info().

754 {
755  this->_fe_map->print_JxW(os);
756 }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
virtual void libMesh::FEAbstract::print_phi ( std::ostream &  os) const
pure virtual

Prints the value of each shape function at each quadrature point. Implement in derived class since this depends on whether the element is vector-valued or not.

Implemented in libMesh::FEGenericBase< OutputType >, and libMesh::FEGenericBase< FEOutputType< T >::type >.

Referenced by get_fe_map(), and print_info().

void libMesh::FEAbstract::print_xyz ( std::ostream &  os) const

Prints the spatial location of each quadrature point (on the physical element).

Definition at line 760 of file fe_abstract.C.

References _fe_map.

Referenced by get_fe_map(), and print_info().

761 {
762  this->_fe_map->print_xyz(os);
763 }
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:517
virtual void libMesh::FEAbstract::reinit ( const Elem elem,
const std::vector< Point > *const  pts = libmesh_nullptr,
const std::vector< Real > *const  weights = libmesh_nullptr 
)
pure virtual

This is at the core of this class. Use this for each new element in the mesh. Reinitializes the requested physical element-dependent data based on the current element elem. By default the element data are computed at the quadrature points specified by the quadrature rule qrule, but any set of points on the reference element may be specified in the optional argument pts.

Note
The FE classes decide which data to initialize based on which accessor functions such as get_phi() or get_d2phi() have been called, so all such accessors should be called before the first reinit().

Implemented in libMesh::FEXYZ< Dim >, libMesh::FESubdivision, libMesh::InfFE< Dim, T_radial, T_map >, libMesh::FE< Dim, T >, libMesh::FE< 2, SUBDIVISION >, libMesh::FE< Dim, HIERARCHIC >, libMesh::FE< Dim, SCALAR >, libMesh::FE< Dim, L2_LAGRANGE >, libMesh::FE< Dim, NEDELEC_ONE >, libMesh::FE< Dim, HERMITE >, libMesh::FE< Dim, CLOUGH >, libMesh::FE< Dim, MONOMIAL >, libMesh::FE< Dim, XYZ >, libMesh::FE< Dim, LAGRANGE >, libMesh::FE< Dim, L2_HIERARCHIC >, and libMesh::FE< Dim, LAGRANGE_VEC >.

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::FEMContext::build_new_fe(), libMesh::System::calculate_norm(), libMesh::ExactErrorEstimator::find_squared_element_error(), and libMesh::JumpErrorEstimator::reinit_sides().

virtual void libMesh::FEAbstract::reinit ( const Elem elem,
const unsigned int  side,
const Real  tolerance = TOLERANCE,
const std::vector< Point > *const  pts = libmesh_nullptr,
const std::vector< Real > *const  weights = libmesh_nullptr 
)
pure virtual

Reinitializes all the physical element-dependent data based on the side of the element elem. The tolerance parameter is passed to the involved call to inverse_map(). By default the element data are computed at the quadrature points specified by the quadrature rule qrule, but any set of points on the reference side element may be specified in the optional argument pts.

Implemented in libMesh::FEXYZ< Dim >, libMesh::FESubdivision, libMesh::InfFE< Dim, T_radial, T_map >, libMesh::FE< Dim, T >, libMesh::FE< 2, SUBDIVISION >, libMesh::FE< Dim, HIERARCHIC >, libMesh::FE< Dim, SCALAR >, libMesh::FE< Dim, L2_LAGRANGE >, libMesh::FE< Dim, NEDELEC_ONE >, libMesh::FE< Dim, HERMITE >, libMesh::FE< Dim, CLOUGH >, libMesh::FE< Dim, MONOMIAL >, libMesh::FE< Dim, XYZ >, libMesh::FE< Dim, LAGRANGE >, libMesh::FE< Dim, L2_HIERARCHIC >, libMesh::FE< Dim, LAGRANGE_VEC >, libMesh::FEXYZ< Dim >, and libMesh::FEXYZ< Dim >.

void libMesh::FEAbstract::set_fe_order ( int  new_order)
inline

Sets the base FE order of the finite element.

Definition at line 431 of file fe_abstract.h.

References fe_type, get_continuity(), is_hierarchic(), and libMesh::FEType::order.

431 { fe_type.order = new_order; }
OrderWrapper order
Definition: fe_type.h:198
virtual bool libMesh::FEAbstract::shapes_need_reinit ( ) const
protectedpure virtual
Returns
true when the shape functions (for this FEFamily) depend on the particular element, and therefore needs to be re-initialized for each new element. false otherwise.

Implemented in libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::InfFE< Dim, T_radial, T_map >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< 2, SUBDIVISION >, libMesh::FE< Dim, HIERARCHIC >, libMesh::FE< Dim, SCALAR >, libMesh::FE< Dim, L2_LAGRANGE >, libMesh::FE< Dim, NEDELEC_ONE >, libMesh::FE< Dim, HERMITE >, libMesh::FE< Dim, CLOUGH >, libMesh::FE< Dim, MONOMIAL >, libMesh::FE< Dim, XYZ >, libMesh::FE< Dim, LAGRANGE >, libMesh::FE< Dim, L2_HIERARCHIC >, libMesh::FE< Dim, LAGRANGE_VEC >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, libMesh::FE< Dim, T >, and libMesh::FE< Dim, T >.

virtual void libMesh::FEAbstract::side_map ( const Elem elem,
const Elem side,
const unsigned int  s,
const std::vector< Point > &  reference_side_points,
std::vector< Point > &  reference_points 
)
pure virtual

Friends And Related Function Documentation

std::ostream& operator<< ( std::ostream &  os,
const FEAbstract fe 
)
friend

Same as above, but allows you to print to a stream.

Definition at line 782 of file fe_abstract.C.

Referenced by get_fe_map().

783 {
784  fe.print_info(os);
785  return os;
786 }

Member Data Documentation

ReferenceCounter::Counts libMesh::ReferenceCounter::_counts
staticprotectedinherited
bool libMesh::ReferenceCounter::_enable_print_counter = true
staticprotectedinherited

Flag to control whether reference count information is printed when print_info is called.

Definition at line 143 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::disable_print_counter_info(), libMesh::ReferenceCounter::enable_print_counter_info(), and libMesh::ReferenceCounter::print_info().

Threads::spin_mutex libMesh::ReferenceCounter::_mutex
staticprotectedinherited

Mutual exclusion object to enable thread-safe reference counting.

Definition at line 137 of file reference_counter.h.

Threads::atomic< unsigned int > libMesh::ReferenceCounter::_n_objects
staticprotectedinherited

The number of objects. Print the reference count information when the number returns to 0.

Definition at line 132 of file reference_counter.h.

Referenced by libMesh::ReferenceCounter::n_objects(), libMesh::ReferenceCounter::ReferenceCounter(), and libMesh::ReferenceCounter::~ReferenceCounter().

unsigned int libMesh::FEAbstract::_p_level
protected

The p refinement level the current data structures are set up for.

Definition at line 579 of file fe_abstract.h.

Referenced by get_order(), and get_p_level().

bool libMesh::FEAbstract::calculate_curl_phi
mutableprotected

Should we calculate shape function curls?

Definition at line 549 of file fe_abstract.h.

Referenced by libMesh::FEGenericBase< FEOutputType< T >::type >::get_curl_phi().

bool libMesh::FEAbstract::calculate_div_phi
mutableprotected

Should we calculate shape function divergences?

Definition at line 554 of file fe_abstract.h.

Referenced by libMesh::FEGenericBase< FEOutputType< T >::type >::get_div_phi().

bool libMesh::FEAbstract::calculate_dphiref
mutableprotected

Should we calculate reference shape function gradients?

Definition at line 559 of file fe_abstract.h.

Referenced by libMesh::FEGenericBase< FEOutputType< T >::type >::get_curl_phi(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phi(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phideta2(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidetadzeta(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidx2(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidxdy(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidxdz(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidxi2(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidxideta(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidxidzeta(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidy2(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidydz(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidz2(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidzeta2(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_div_phi(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphi(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphideta(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphidx(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphidxi(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphidy(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphidz(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphidzeta(), and libMesh::InfFE< Dim, T_radial, T_map >::reinit().

bool libMesh::FEAbstract::calculate_phi
mutableprotected

Should we calculate shape functions?

Definition at line 534 of file fe_abstract.h.

Referenced by libMesh::FEGenericBase< FEOutputType< T >::type >::get_phi(), and libMesh::InfFE< Dim, T_radial, T_map >::reinit().

bool libMesh::FEAbstract::calculations_started
mutableprotected

Have calculations with this object already been started? Then all get_* functions should already have been called.

Definition at line 529 of file fe_abstract.h.

Referenced by libMesh::FEGenericBase< FEOutputType< T >::type >::get_curl_phi(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phi(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phideta2(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidetadzeta(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidx2(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidxdy(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidxdz(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidxi2(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidxideta(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidxidzeta(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidy2(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidydz(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidz2(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_d2phidzeta2(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_div_phi(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphi(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphideta(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphidx(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphidxi(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphidy(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphidz(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_dphidzeta(), libMesh::FEGenericBase< FEOutputType< T >::type >::get_phi(), and libMesh::FESubdivision::init_shape_functions().

const unsigned int libMesh::FEAbstract::dim
protected
ElemType libMesh::FEAbstract::elem_type
protected

The element type the current data structures are set up for.

Definition at line 573 of file fe_abstract.h.

Referenced by libMesh::FESubdivision::attach_quadrature_rule(), get_type(), and libMesh::InfFE< Dim, T_radial, T_map >::reinit().

QBase* libMesh::FEAbstract::qrule
protected

A pointer to the quadrature rule employed

Definition at line 584 of file fe_abstract.h.

Referenced by libMesh::InfFE< Dim, T_radial, T_map >::attach_quadrature_rule(), and libMesh::FESubdivision::attach_quadrature_rule().

bool libMesh::FEAbstract::shapes_on_quadrature
protected

A flag indicating if current data structures correspond to quadrature rule points

Definition at line 590 of file fe_abstract.h.

Referenced by libMesh::FESubdivision::attach_quadrature_rule().


The documentation for this class was generated from the following files: