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
 
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 UniquePtr< 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

UniquePtr< 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 that the amount of virtual members is kept to a minimum, and the sophisticated template scheme of FE is quite likely to offer acceptably fast code.

All calls to static members of the FE classes should be requested through the FEInterface. This interface class offers sort-of runtime polymorphism for the templated finite element classes. Even internal library classes, like DofMap, request the number of dof's through this interface class. Note that this 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 110 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 600 of file fe_abstract.h.

601  :
602  _fe_map( FEMap::build(fet) ),
603  dim(d),
604  calculations_started(false),
605  calculate_phi(false),
606  calculate_dphi(false),
607  calculate_d2phi(false),
608  calculate_curl_phi(false),
609  calculate_div_phi(false),
610  calculate_dphiref(false),
611  fe_type(fet),
613  _p_level(0),
615  shapes_on_quadrature(false)
616 {
617 }
unsigned int _p_level
Definition: fe_abstract.h:571
const class libmesh_nullptr_t libmesh_nullptr
const unsigned int dim
Definition: fe_abstract.h:516
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
static UniquePtr< FEMap > build(FEType fe_type)
Definition: fe_map.C:50
libMesh::FEAbstract::~FEAbstract ( )
inlinevirtual

Destructor.

Definition at line 621 of file fe_abstract.h.

622 {
623 }

Member Function Documentation

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

Builds a specific finite element type. A UniquePtr<FEAbstract> is returned to prevent a memory leak. This way the user need not remember to delete the object.

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::FEMContext(), 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 UniquePtr<FEAbstract>(new FE<0,CLOUGH>(fet));
56 
57  case HERMITE:
58  return UniquePtr<FEAbstract>(new FE<0,HERMITE>(fet));
59 
60  case LAGRANGE:
61  return UniquePtr<FEAbstract>(new FE<0,LAGRANGE>(fet));
62 
63  case LAGRANGE_VEC:
64  return UniquePtr<FEAbstract>(new FE<0,LAGRANGE_VEC>(fet));
65 
66  case L2_LAGRANGE:
67  return UniquePtr<FEAbstract>(new FE<0,L2_LAGRANGE>(fet));
68 
69  case HIERARCHIC:
70  return UniquePtr<FEAbstract>(new FE<0,HIERARCHIC>(fet));
71 
72  case L2_HIERARCHIC:
73  return UniquePtr<FEAbstract>(new FE<0,L2_HIERARCHIC>(fet));
74 
75  case MONOMIAL:
76  return UniquePtr<FEAbstract>(new FE<0,MONOMIAL>(fet));
77 
78 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
79  case SZABAB:
80  return UniquePtr<FEAbstract>(new FE<0,SZABAB>(fet));
81 
82  case BERNSTEIN:
83  return UniquePtr<FEAbstract>(new FE<0,BERNSTEIN>(fet));
84 #endif
85 
86  case XYZ:
87  return UniquePtr<FEAbstract>(new FEXYZ<0>(fet));
88 
89  case SCALAR:
90  return UniquePtr<FEAbstract>(new 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 UniquePtr<FEAbstract>(new FE<1,CLOUGH>(fet));
103 
104  case HERMITE:
105  return UniquePtr<FEAbstract>(new FE<1,HERMITE>(fet));
106 
107  case LAGRANGE:
108  return UniquePtr<FEAbstract>(new FE<1,LAGRANGE>(fet));
109 
110  case LAGRANGE_VEC:
111  return UniquePtr<FEAbstract>(new FE<1,LAGRANGE_VEC>(fet));
112 
113  case L2_LAGRANGE:
114  return UniquePtr<FEAbstract>(new FE<1,L2_LAGRANGE>(fet));
115 
116  case HIERARCHIC:
117  return UniquePtr<FEAbstract>(new FE<1,HIERARCHIC>(fet));
118 
119  case L2_HIERARCHIC:
120  return UniquePtr<FEAbstract>(new FE<1,L2_HIERARCHIC>(fet));
121 
122  case MONOMIAL:
123  return UniquePtr<FEAbstract>(new FE<1,MONOMIAL>(fet));
124 
125 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
126  case SZABAB:
127  return UniquePtr<FEAbstract>(new FE<1,SZABAB>(fet));
128 
129  case BERNSTEIN:
130  return UniquePtr<FEAbstract>(new FE<1,BERNSTEIN>(fet));
131 #endif
132 
133  case XYZ:
134  return UniquePtr<FEAbstract>(new FEXYZ<1>(fet));
135 
136  case SCALAR:
137  return UniquePtr<FEAbstract>(new 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 UniquePtr<FEAbstract>(new FE<2,CLOUGH>(fet));
152 
153  case HERMITE:
154  return UniquePtr<FEAbstract>(new FE<2,HERMITE>(fet));
155 
156  case LAGRANGE:
157  return UniquePtr<FEAbstract>(new FE<2,LAGRANGE>(fet));
158 
159  case LAGRANGE_VEC:
160  return UniquePtr<FEAbstract>(new FE<2,LAGRANGE_VEC>(fet));
161 
162  case L2_LAGRANGE:
163  return UniquePtr<FEAbstract>(new FE<2,L2_LAGRANGE>(fet));
164 
165  case HIERARCHIC:
166  return UniquePtr<FEAbstract>(new FE<2,HIERARCHIC>(fet));
167 
168  case L2_HIERARCHIC:
169  return UniquePtr<FEAbstract>(new FE<2,L2_HIERARCHIC>(fet));
170 
171  case MONOMIAL:
172  return UniquePtr<FEAbstract>(new FE<2,MONOMIAL>(fet));
173 
174 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
175  case SZABAB:
176  return UniquePtr<FEAbstract>(new FE<2,SZABAB>(fet));
177 
178  case BERNSTEIN:
179  return UniquePtr<FEAbstract>(new FE<2,BERNSTEIN>(fet));
180 #endif
181 
182  case XYZ:
183  return UniquePtr<FEAbstract>(new FEXYZ<2>(fet));
184 
185  case SCALAR:
186  return UniquePtr<FEAbstract>(new FEScalar<2>(fet));
187 
188  case NEDELEC_ONE:
189  return UniquePtr<FEAbstract>(new FENedelecOne<2>(fet));
190 
191  case SUBDIVISION:
192  return UniquePtr<FEAbstract>(new 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 UniquePtr<FEAbstract>(new FE<3,HERMITE>(fet));
210 
211  case LAGRANGE:
212  return UniquePtr<FEAbstract>(new FE<3,LAGRANGE>(fet));
213 
214  case LAGRANGE_VEC:
215  return UniquePtr<FEAbstract>(new FE<3,LAGRANGE_VEC>(fet));
216 
217  case L2_LAGRANGE:
218  return UniquePtr<FEAbstract>(new FE<3,L2_LAGRANGE>(fet));
219 
220  case HIERARCHIC:
221  return UniquePtr<FEAbstract>(new FE<3,HIERARCHIC>(fet));
222 
223  case L2_HIERARCHIC:
224  return UniquePtr<FEAbstract>(new FE<3,L2_HIERARCHIC>(fet));
225 
226  case MONOMIAL:
227  return UniquePtr<FEAbstract>(new FE<3,MONOMIAL>(fet));
228 
229 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
230  case SZABAB:
231  return UniquePtr<FEAbstract>(new FE<3,SZABAB>(fet));
232 
233  case BERNSTEIN:
234  return UniquePtr<FEAbstract>(new FE<3,BERNSTEIN>(fet));
235 #endif
236 
237  case XYZ:
238  return UniquePtr<FEAbstract>(new FEXYZ<3>(fet));
239 
240  case SCALAR:
241  return UniquePtr<FEAbstract>(new FEScalar<3>(fet));
242 
243  case NEDELEC_ONE:
244  return UniquePtr<FEAbstract>(new 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 
255  libmesh_error_msg("We'll never get here!");
256  return UniquePtr<FEAbstract>();
257 }
const unsigned int dim
Definition: fe_abstract.h:516
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 791 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::libmesh_assert(), libmesh_nullptr, libMesh::FEInterface::n_dofs(), libMesh::Elem::n_sides(), libMesh::Elem::neighbor_ptr(), libMesh::Elem::parent(), libMesh::Real, libMesh::remote_elem, libMesh::FEInterface::shape(), libMesh::Threads::spin_mtx, and libMesh::Elem::subactive().

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

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

References libMesh::ReferenceCounter::_enable_print_counter.

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

101 {
102  _enable_print_counter = true;
103  return;
104 }
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 376 of file fe_abstract.h.

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

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

Definition at line 263 of file fe_abstract.h.

References _fe_map.

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

Definition at line 293 of file fe_abstract.h.

References _fe_map.

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

Definition at line 257 of file fe_abstract.h.

References _fe_map.

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

Definition at line 279 of file fe_abstract.h.

References _fe_map.

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

Definition at line 287 of file fe_abstract.h.

References _fe_map.

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

Definition at line 271 of file fe_abstract.h.

References _fe_map.

272  { return this->_fe_map->get_d2xyzdzeta2(); }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
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 323 of file fe_abstract.h.

References _fe_map.

324  { return this->_fe_map->get_detadx(); }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
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 330 of file fe_abstract.h.

References _fe_map.

331  { return this->_fe_map->get_detady(); }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
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 337 of file fe_abstract.h.

References _fe_map.

338  { return this->_fe_map->get_detadz(); }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
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 302 of file fe_abstract.h.

References _fe_map.

303  { return this->_fe_map->get_dxidx(); }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
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 309 of file fe_abstract.h.

References _fe_map.

310  { return this->_fe_map->get_dxidy(); }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
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 316 of file fe_abstract.h.

References _fe_map.

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

Definition at line 244 of file fe_abstract.h.

References _fe_map.

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

Definition at line 237 of file fe_abstract.h.

References _fe_map.

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

Definition at line 251 of file fe_abstract.h.

References _fe_map.

252  { return _fe_map->get_dxyzdzeta(); }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
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 344 of file fe_abstract.h.

References _fe_map.

345  { return this->_fe_map->get_dzetadx(); }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
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 351 of file fe_abstract.h.

References _fe_map.

352  { return this->_fe_map->get_dzetady(); }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
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 358 of file fe_abstract.h.

References _fe_map.

359  { return this->_fe_map->get_dzetadz(); }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
FEFamily libMesh::FEAbstract::get_family ( ) const
inline
Returns
the finite element family of this element.

Definition at line 440 of file fe_abstract.h.

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

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

440 { return fe_type.family; }
FEFamily family
Definition: fe_type.h:206
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 419 of file fe_abstract.h.

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

419 { return static_cast<Order>(fe_type.order + _p_level); }
unsigned int _p_level
Definition: fe_abstract.h:571
OrderWrapper order
Definition: fe_type.h:200
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 409 of file fe_abstract.h.

References _p_level.

409 { return _p_level; }
unsigned int _p_level
Definition: fe_abstract.h:571
void libMesh::FEAbstract::get_refspace_nodes ( const ElemType  t,
std::vector< Point > &  nodes 
)
static

returns the reference space nodes coordinates given the element type

Definition at line 259 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::TET10, libMesh::TET4, libMesh::TRI3, libMesh::TRI6, and libMesh::TRISHELL3.

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

Definition at line 364 of file fe_abstract.h.

References _fe_map.

365  { return this->_fe_map->get_tangents(); }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
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 403 of file fe_abstract.h.

References elem_type.

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

403 { 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 160 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().

161 {
162  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
163  std::pair<unsigned int, unsigned int> & p = _counts[name];
164 
165  p.first++;
166 }
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 173 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().

174 {
175  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
176  std::pair<unsigned int, unsigned int> & p = _counts[name];
177 
178  p.second++;
179 }
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 556 of file fe_abstract.C.

References libMesh::EDGE2, libMesh::EDGE3, libMesh::EDGE4, libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::INFHEX8, 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::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().

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

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

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

765 {
766  os << "phi[i][j]: Shape function i at quadrature pt. j" << std::endl;
767  this->print_phi(os);
768 
769  os << "dphi[i][j]: Shape function i's gradient at quadrature pt. j" << std::endl;
770  this->print_dphi(os);
771 
772  os << "XYZ locations of the quadrature pts." << std::endl;
773  this->print_xyz(os);
774 
775  os << "Values of JxW at the quadrature pts." << std::endl;
776  this->print_JxW(os);
777 }
void print_JxW(std::ostream &os) const
Definition: fe_abstract.C:751
void print_xyz(std::ostream &os) const
Definition: fe_abstract.C:758
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 751 of file fe_abstract.C.

References _fe_map.

Referenced by get_fe_map(), and print_info().

752 {
753  this->_fe_map->print_JxW(os);
754 }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
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 758 of file fe_abstract.C.

References _fe_map.

Referenced by get_fe_map(), and print_info().

759 {
760  this->_fe_map->print_xyz(os);
761 }
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:510
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 that 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 paremeter 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 424 of file fe_abstract.h.

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

424 { fe_type.order = new_order; }
OrderWrapper order
Definition: fe_type.h:200
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 780 of file fe_abstract.C.

Referenced by get_fe_map().

781 {
782  fe.print_info(os);
783  return os;
784 }

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 134 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 128 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 123 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 571 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 542 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 547 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 552 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 527 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 522 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 565 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 576 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 582 of file fe_abstract.h.

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


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