20 #ifndef LIBMESH_FE_ABSTRACT_H 21 #define LIBMESH_FE_ABSTRACT_H 24 #include "libmesh/reference_counted_object.h" 25 #include "libmesh/point.h" 26 #include "libmesh/vector_value.h" 27 #include "libmesh/fe_type.h" 28 #include "libmesh/fe_map.h" 29 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 30 #include "libmesh/tensor_value.h" 38 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 39 #define virtual_for_inffe virtual 41 #define virtual_for_inffe 50 template <
typename T>
class DenseMatrix;
51 template <
typename T>
class DenseVector;
57 template <
typename T>
class NumericVector;
61 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS 62 class NodeConstraints;
65 #ifdef LIBMESH_ENABLE_PERIODIC 66 class PeriodicBoundaries;
67 class PointLocatorBase;
70 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 71 template <
unsigned int Dim, FEFamily T_radial, InfMapType T_map>
124 static std::unique_ptr<FEAbstract>
build (
const unsigned int dim,
142 const std::vector<Point> *
const pts =
nullptr,
143 const std::vector<Real> *
const weights =
nullptr) = 0;
151 const std::vector<Point> & ,
152 const std::vector<Real> & )
154 libmesh_error_msg(
"Customized dual shape coefficient calculation has not been implemented for this FE type.");
163 libmesh_error_msg(
"Default dual shape coefficient calculation has not been implemented for this FE type.");
175 const unsigned int side,
177 const std::vector<Point> *
const pts =
nullptr,
178 const std::vector<Real> *
const weights =
nullptr) = 0;
189 const unsigned int edge,
191 const std::vector<Point> * pts =
nullptr,
192 const std::vector<Real> * weights =
nullptr) = 0;
200 const unsigned int s,
201 const std::vector<Point> & reference_side_points,
202 std::vector<Point> & reference_points) = 0;
219 std::vector<Point> & nodes);
222 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS 229 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS 231 #ifdef LIBMESH_ENABLE_PERIODIC 233 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS 243 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS 245 #endif // LIBMESH_ENABLE_PERIODIC 276 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 320 {
return _fe_map->get_dxyzdzeta(); }
322 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 454 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 546 virtual void print_phi(std::ostream & os)
const =0;
554 virtual void print_dphi(std::ostream & os)
const =0;
557 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 564 virtual void print_d2phi(std::ostream & os)
const =0;
677 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 761 #endif // LIBMESH_FE_ABSTRACT_H virtual void print_dual_phi(std::ostream &os) const =0
bool calculate_d2phi
Should we calculate shape function hessians?
bool calculate_curl_phi
Should we calculate shape function curls?
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
FEFamily family
The type of finite element.
virtual bool is_hierarchic() const =0
ElemType
Defines an enum for geometric element types.
virtual void print_phi(std::ostream &os) const =0
Prints the value of each shape function at each quadrature point.
void set_calculate_dual(const bool val)
set calculate_dual as needed
unsigned int get_p_level() const
bool calculations_started
Have calculations with this object already been started? Then all get_* functions should already have...
Order
defines an enum for polynomial orders.
bool calculate_phi
Should we calculate shape functions?
virtual void reinit_default_dual_shape_coeffs(const Elem *)
This re-computes the dual shape function coefficients using DEFAULT qrule.
unsigned int _p_level
The p refinement level the current data structures are set up for.
virtual_for_inffe const std::vector< RealGradient > & get_d2xyzdetadzeta() const
bool shapes_on_quadrature
A flag indicating if current data structures correspond to quadrature rule points.
static void compute_node_constraints(NodeConstraints &constraints, const Elem *elem)
Computes the nodal constraint contributions (for non-conforming adapted meshes), using Lagrange geome...
virtual_for_inffe const std::vector< Real > & get_curvatures() const
virtual_for_inffe const std::vector< RealGradient > & get_dxyzdeta() const
unsigned int get_dim() const
ElemType get_type() const
virtual_for_inffe const std::vector< RealGradient > & get_d2xyzdxi2() const
static constexpr Real TOLERANCE
virtual_for_inffe const std::vector< Real > & get_dxidz() const
We're using a class instead of a typedef to allow forward declarations and future flexibility...
unsigned int _elem_p_level
The element p-refinement level the current data structures are set up for.
virtual_for_inffe const std::vector< RealGradient > & get_dxyzdzeta() const
virtual void request_dual_dphi() const =0
This is the base class from which all geometric element types are derived.
virtual unsigned int n_quadrature_points() const =0
FEAbstract(const unsigned int dim, const FEType &fet)
Constructor.
virtual_for_inffe const std::vector< RealGradient > & get_d2xyzdxideta() const
The Node constraint storage format.
virtual void request_dual_phi() const =0
OrderWrapper order
The approximation order of the element.
bool calculate_default_dual_coeff
Are we calculating the coefficient for the dual basis using the default qrule?
The libMesh namespace provides an interface to certain functionality in the library.
virtual_for_inffe const std::vector< std::vector< Point > > & get_tangents() const
bool calculate_div_phi
Should we calculate shape function divergences?
virtual_for_inffe const std::vector< RealGradient > & get_d2xyzdeta2() const
virtual void request_phi() const =0
request phi calculations
static bool on_reference_element(const Point &p, const ElemType t, const Real eps=TOLERANCE)
This is the MeshBase class.
virtual_for_inffe const std::vector< Real > & get_dxidy() const
virtual FEContinuity get_continuity() const =0
virtual void print_dual_d2phi(std::ostream &os) const =0
FEType get_fe_type() const
virtual_for_inffe const std::vector< Real > & get_detadz() const
virtual_for_inffe const std::vector< RealGradient > & get_d2xyzdzeta2() const
virtual unsigned int n_shape_functions() const =0
virtual void print_dphi(std::ostream &os) const =0
Prints the value of each shape function's derivative at each quadrature point.
const unsigned int dim
The dimensionality of the object.
void print_xyz(std::ostream &os) const
Prints the spatial location of each quadrature point (on the physical element).
virtual_for_inffe const std::vector< RealGradient > & get_d2xyzdxidzeta() const
virtual const std::vector< Real > & get_JxWxdecay_sq() const
This function is the variant of get_JxW() for InfFE.
static void get_refspace_nodes(const ElemType t, std::vector< Point > &nodes)
void add_p_level_in_reinit(bool value)
Indicate whether to add p-refinement levels in init/reinit methods.
virtual_for_inffe const std::vector< Real > & get_JxW() const
virtual void reinit(const Elem *elem, const std::vector< Point > *const pts=nullptr, const std::vector< Real > *const weights=nullptr)=0
This is at the core of this class.
virtual bool shapes_need_reinit() const =0
QBase * qrule
A pointer to the quadrature rule employed.
virtual_for_inffe const std::vector< Point > & get_xyz() const
This is the base class for point locators.
virtual_for_inffe const std::vector< Real > & get_detadx() const
This class implements reference counting.
bool calculate_dphiref
Should we calculate reference shape function gradients?
virtual_for_inffe const std::vector< Real > & get_detady() const
void print_JxW(std::ostream &os) const
Prints the Jacobian times the weight for each quadrature point.
void print_info(std::ostream &os) const
Prints all the relevant information about the current element.
virtual void request_dphi() const =0
request dphi calculations
bool _add_p_level_in_reinit
Whether to add p-refinement levels in init/reinit methods.
virtual void edge_reinit(const Elem *elem, const unsigned int edge, const Real tolerance=TOLERANCE, const std::vector< Point > *pts=nullptr, const std::vector< Real > *weights=nullptr)=0
Reinitializes all the physical element-dependent data based on the edge of the element elem...
void set_calculate_default_dual_coeff(const bool val)
set calculate_default_dual_coeff as needed
virtual void compute_shape_functions(const Elem *, const std::vector< Point > &)=0
After having updated the jacobian and the transformation from local to global coordinates in FEMap::c...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual_for_inffe const std::vector< Point > & get_normals() const
bool calculate_dphi
Should we calculate shape function gradients?
bool calculate_map
Are we calculating mapping functions?
FEContinuity
defines an enum for finite element types to libmesh_assert a certain level (or type? Hcurl?) of continuity.
virtual_for_inffe const std::vector< Real > & get_dxidx() const
virtual_for_inffe const std::vector< Real > & get_dzetady() const
virtual_for_inffe const std::vector< Real > & get_dzetadz() const
static std::unique_ptr< FEAbstract > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
FEFamily get_family() const
virtual void print_dual_dphi(std::ostream &os) const =0
friend std::ostream & operator<<(std::ostream &os, const FEAbstract &fe)
Same as above, but allows you to print to a stream.
virtual ~FEAbstract()
Destructor.
virtual void reinit_dual_shape_coeffs(const Elem *, const std::vector< Point > &, const std::vector< Real > &)
This re-computes the dual shape function coefficients using CUSTOMIZED qrule.
void set_fe_order(int new_order)
Sets the base FE order of the finite element.
virtual_for_inffe const std::vector< Real > & get_dzetadx() const
virtual void print_d2phi(std::ostream &os) const =0
Prints the value of each shape function's second derivatives at each quadrature point.
virtual_for_inffe const std::vector< RealGradient > & get_dxyzdxi() const
const FEMap & get_fe_map() const
FEFamily
defines an enum for finite element families.
This class forms the foundation from which generic finite elements may be derived.
FEType fe_type
The finite element type for this object.
std::unique_ptr< FEMap > _fe_map
bool calculate_dual
Are we calculating dual basis?
Class contained in FE that encapsulates mapping (i.e.
ElemType elem_type
The element type the current data structures are set up for.
A Point defines a location in LIBMESH_DIM dimensional Real space.
virtual void attach_quadrature_rule(QBase *q)=0
Provides the class with the quadrature rule.
The QBase class provides the basic functionality from which various quadrature rules can be derived...
void ErrorVector unsigned int
bool calculate_nothing
Are we potentially deliberately calculating nothing?
static void compute_periodic_node_constraints(NodeConstraints &constraints, const PeriodicBoundaries &boundaries, const MeshBase &mesh, const PointLocatorBase *point_locator, const Elem *elem)
Computes the node position constraint equation contributions (for meshes with periodic boundary condi...
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
Computes the reference space quadrature points on the side of an element based on the side quadrature...
bool add_p_level_in_reinit() const
Whether to add p-refinement levels in init/reinit methods.