20 #ifndef LIBMESH_FE_BASE_H 21 #define LIBMESH_FE_BASE_H 24 #include "libmesh/libmesh_common.h" 25 #include "libmesh/compare_types.h" 26 #include "libmesh/fe_abstract.h" 27 #include "libmesh/fe_transformation_base.h" 28 #include "libmesh/point.h" 29 #include "libmesh/reference_counted_object.h" 30 #include "libmesh/tensor_tools.h" 31 #include "libmesh/type_n_tensor.h" 32 #include "libmesh/vector_value.h" 33 #include "libmesh/dense_matrix.h" 45 template <
typename T>
class DenseMatrix;
46 template <
typename T>
class DenseVector;
52 template <
typename T>
class NumericVector;
57 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS 61 #ifdef LIBMESH_ENABLE_PERIODIC 66 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 67 template <
unsigned int Dim, FEFamily T_radial, InfMapType T_map>
84 template <
typename OutputType>
112 static std::unique_ptr<FEGenericBase>
build (
const unsigned int dim,
130 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 140 static std::unique_ptr<FEGenericBase>
build_InfFE (
const unsigned int dim,
145 #ifdef LIBMESH_ENABLE_AMR 155 const unsigned int variable_number,
167 const Elem * coarse_elem,
169 const unsigned int var,
170 const bool use_old_dof_indices =
false);
180 const Elem * coarse_elem,
182 const bool use_old_dof_indices =
false);
184 #endif // #ifdef LIBMESH_ENABLE_AMR 186 #ifdef LIBMESH_ENABLE_PERIODIC 198 const unsigned int variable_number,
201 #endif // LIBMESH_ENABLE_PERIODIC 207 const std::vector<std::vector<OutputShape>> &
get_phi()
const 230 const std::vector<std::vector<OutputGradient>> &
get_dphi()
const 261 const std::vector<std::vector<OutputDivergence>> &
get_div_phi()
const 269 const std::vector<std::vector<OutputShape>> &
get_dphidx()
const 277 const std::vector<std::vector<OutputShape>> &
get_dphidy()
const 285 const std::vector<std::vector<OutputShape>> &
get_dphidz()
const 313 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 319 const std::vector<std::vector<OutputTensor>> &
get_d2phi()
const 423 #endif //LIBMESH_ENABLE_SECOND_DERIVATIVES 425 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 519 virtual void print_phi(std::ostream & os)
const override;
526 virtual void print_dphi(std::ostream & os)
const override;
529 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 535 virtual void print_d2phi(std::ostream & os)
const override;
545 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 572 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 609 std::unique_ptr<FETransformationBase<OutputType>>
_fe_trans;
614 std::vector<std::vector<OutputShape>>
phi;
620 std::vector<std::vector<OutputGradient>>
dphi;
636 std::vector<std::vector<OutputDivergence>>
div_phi;
641 std::vector<std::vector<OutputShape>>
dphidxi;
656 std::vector<std::vector<OutputShape>>
dphidx;
661 std::vector<std::vector<OutputShape>>
dphidy;
666 std::vector<std::vector<OutputShape>>
dphidz;
669 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 674 std::vector<std::vector<OutputTensor>>
d2phi;
740 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 773 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS 780 template <
unsigned int friend_Dim, FEFamily friend_T_radial, InfMapType friend_T_map>
791 template <
typename OutputType>
795 "Computation of dual shape functions for vector finite element " 796 "families is not currently implemented");
799 template <
typename OutputType>
803 "Computation of dual shape functions for vector finite element " 804 "families is not currently implemented");
825 template <
typename OutputType>
843 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
859 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
869 template <
typename OutputType>
877 #endif // LIBMESH_FE_BASE_H FEGenericBase(const unsigned int dim, const FEType &fet)
Constructor.
bool calculate_d2phi
Should we calculate shape function hessians?
bool calculate_curl_phi
Should we calculate shape function curls?
std::vector< std::vector< OutputTensor > > d2phi
Shape function second derivative values.
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
std::vector< std::vector< OutputShape > > dphidxi
Shape function derivatives in the xi direction.
std::vector< std::vector< OutputShape > > d2phidxdz
Shape function second derivatives in the x-z direction.
const std::vector< std::vector< OutputShape > > & get_d2phidxdy() const
virtual_for_inffe const std::vector< std::vector< OutputDivergence > > & get_div_phi() const
std::vector< std::vector< OutputShape > > dphidzeta
Shape function derivatives in the zeta direction.
bool calculations_started
Have calculations with this object already been started? Then all get_* functions should already have...
bool calculating_nothing() const
std::vector< std::vector< OutputShape > > d2phidydz
Shape function second derivatives in the y-z direction.
bool calculate_phi
Should we calculate shape functions?
A specific instantiation of the FEBase class.
const std::vector< std::vector< OutputShape > > & get_dual_phi() const
std::unique_ptr< FETransformationBase< OutputType > > _fe_trans
Object that handles computing shape function values, gradients, etc in the physical domain...
TensorTools::DecrementRank< OutputNumber >::type OutputNumberDivergence
virtual const std::vector< Real > & get_Sobolev_weightxR_sq() const
virtual void print_dphi(std::ostream &os) const override
Prints the value of each shape function's derivative at each quadrature point.
void compute_dual_shape_functions()
Compute dual_phi, dual_dphi, dual_d2phi It is only valid for this to be called after reinit has occur...
virtual void request_dphi() const override
request dphi calculations
const std::vector< std::vector< OutputShape > > & get_d2phideta2() const
virtual void init_base_shape_functions(const std::vector< Point > &qp, const Elem *e)=0
Initialize the data fields for the base of an an infinite element.
We're using a class instead of a typedef to allow forward declarations and future flexibility...
const std::vector< std::vector< OutputShape > > & get_d2phidxideta() const
std::vector< std::vector< OutputShape > > d2phidxideta
Shape function second derivatives in the xi-eta direction.
std::vector< std::vector< OutputTensor > > dual_d2phi
const std::vector< std::vector< OutputShape > > & get_d2phidydz() const
This is the base class from which all geometric element types are derived.
TensorTools::IncrementRank< OutputNumberGradient >::type OutputNumberTensor
virtual const std::vector< Real > & get_Sobolev_weight() const
virtual void request_phi() const override
request phi calculations
std::vector< Real > weight
Used for certain infinite element families: the additional radial weight in local coordinates...
std::vector< std::vector< OutputGradient > > dual_dphi
TensorTools::IncrementRank< OutputNumber >::type OutputNumberGradient
The Node constraint storage format.
This class defines a vector in LIBMESH_DIM dimensional Real or Complex space.
TensorTools::IncrementRank< OutputShape >::type OutputGradient
std::vector< std::vector< OutputShape > > d2phidx2
Shape function second derivatives in the x direction.
std::vector< std::vector< OutputShape > > curl_phi
Shape function curl values.
The libMesh namespace provides an interface to certain functionality in the library.
const std::vector< std::vector< OutputShape > > & get_dphideta() const
DenseMatrix< Real > dual_coeff
Coefficient matrix for the dual basis.
std::vector< std::vector< OutputShape > > dphidy
Shape function derivatives in the y direction.
FEGenericBase< RealGradient > FEVectorBase
TensorTools::DecrementRank< OutputShape >::type OutputDivergence
std::vector< std::vector< OutputShape > > d2phidy2
Shape function second derivatives in the y direction.
virtual void request_dual_phi() const override
bool calculate_div_phi
Should we calculate shape function divergences?
This is the MeshBase class.
const std::vector< std::vector< OutputShape > > & get_dphidzeta() const
const std::vector< OutputGradient > & get_dphase() const
const DenseMatrix< Real > & get_dual_coeff() const
std::vector< std::vector< OutputShape > > d2phidetadzeta
Shape function second derivatives in the eta-zeta direction.
const std::vector< std::vector< OutputShape > > & get_dphidx() const
static void coarsened_dof_values(const NumericVector< Number > &global_vector, const DofMap &dof_map, const Elem *coarse_elem, DenseVector< Number > &coarse_dofs, const unsigned int var, const bool use_old_dof_indices=false)
Creates a local projection on coarse_elem, based on the DoF values in global_vector for it's children...
std::vector< std::vector< OutputShape > > d2phidxidzeta
Shape function second derivatives in the xi-zeta direction.
This class handles the numbering of degrees of freedom on a mesh.
std::vector< std::vector< OutputShape > > d2phidxdy
Shape function second derivatives in the x-y direction.
std::vector< std::vector< OutputShape > > dphidx
Shape function derivatives in the x direction.
const std::vector< std::vector< OutputGradient > > & get_dphi() const
std::vector< std::vector< OutputShape > > phi
Shape function values.
void compute_dual_shape_coeffs(const std::vector< Real > &JxW, const std::vector< std::vector< OutputShape >> &phi)
Compute the dual basis coefficients dual_coeff we rely on the JxW (or weights) and the phi values...
std::vector< std::vector< OutputShape > > d2phideta2
Shape function second derivatives in the eta direction.
const unsigned int dim
The dimensionality of the object.
std::vector< std::vector< OutputShape > > dual_phi
const std::vector< std::vector< OutputGradient > > & get_dual_dphi() const
TensorTools::IncrementRank< OutputGradient >::type OutputTensor
virtual const std::vector< std::vector< OutputShape > > & get_phi_over_decayxR() const
const std::vector< std::vector< OutputShape > > & get_d2phidx2() const
virtual void print_phi(std::ostream &os) const override
Prints the value of each shape function at each quadrature point.
virtual const std::vector< std::vector< OutputGradient > > & get_dphi_over_decay() const
std::vector< OutputGradient > dphase
Used for certain infinite element families: the first derivatives of the phase term in global coordin...
static std::unique_ptr< FEGenericBase > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
virtual void print_dual_dphi(std::ostream &os) const override
virtual ~FEGenericBase()
Destructor.
const std::vector< std::vector< OutputShape > > & get_dphidy() const
static void compute_proj_constraints(DofConstraints &constraints, DofMap &dof_map, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for non-conforming adapted meshes) corresponding to var...
const std::vector< std::vector< OutputShape > > & get_d2phidy2() const
This is the base class for point locators.
const std::vector< std::vector< OutputTensor > > & get_d2phi() const
OutputType OutputShape
Convenient typedefs for gradients of output, hessians of output, and potentially-complex-valued versi...
TensorTools::MakeNumber< OutputShape >::type OutputNumber
std::vector< std::vector< OutputDivergence > > div_phi
Shape function divergence values.
const std::vector< std::vector< OutputShape > > & get_d2phidxi2() const
FEGenericBase< Real > FEBase
virtual_for_inffe void determine_calculations()
Determine which values are to be calculated, for both the FE itself and for the FEMap.
const std::vector< std::vector< OutputShape > > & get_d2phidxdz() const
bool calculate_dphiref
Should we calculate reference shape function gradients?
std::vector< std::vector< OutputGradient > > dphi
Shape function derivative values.
virtual void request_dual_dphi() const override
static void compute_periodic_constraints(DofConstraints &constraints, DofMap &dof_map, const PeriodicBoundaries &boundaries, const MeshBase &mesh, const PointLocatorBase *point_locator, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for meshes with periodic boundary conditions) correspon...
const std::vector< std::vector< OutputShape > > & get_d2phidxidzeta() const
const std::vector< std::vector< OutputShape > > & get_d2phidetadzeta() const
bool calculate_dphi
Should we calculate shape function gradients?
bool calculate_map
Are we calculating mapping functions?
virtual const std::vector< RealGradient > & get_Sobolev_dweight() const
virtual void compute_shape_functions(const Elem *elem, const std::vector< Point > &qp) override
After having updated the jacobian and the transformation from local to global coordinates in FEMap::c...
std::vector< std::vector< OutputShape > > d2phidz2
Shape function second derivatives in the z direction.
virtual void print_dual_phi(std::ostream &os) const override
const std::vector< std::vector< OutputShape > > & get_d2phidzeta2() const
static std::unique_ptr< FEGenericBase > build_InfFE(const unsigned int dim, const FEType &type)
Builds a specific infinite element type.
virtual const std::vector< RealGradient > & get_Sobolev_dweightxR_sq() const
std::vector< std::vector< OutputShape > > d2phidxi2
Shape function second derivatives in the xi direction.
const std::vector< std::vector< OutputShape > > & get_dphidxi() const
This class forms the foundation from which generic finite elements may be derived.
virtual void print_d2phi(std::ostream &os) const override
Prints the value of each shape function's second derivatives at each quadrature point.
bool calculate_dual
Are we calculating dual basis?
virtual void print_dual_d2phi(std::ostream &os) const override
The constraint matrix storage format.
bool calculate_nothing
Are we potentially deliberately calculating nothing?
This class forms the foundation from which generic finite elements may be derived.
std::vector< std::vector< OutputShape > > d2phidzeta2
Shape function second derivatives in the zeta direction.
virtual const std::vector< std::vector< OutputGradient > > & get_dphi_over_decayxR() const
std::vector< std::vector< OutputShape > > dphidz
Shape function derivatives in the z direction.
const std::vector< std::vector< OutputShape > > & get_d2phidz2() const
virtual_for_inffe const std::vector< std::vector< OutputShape > > & get_curl_phi() const
const std::vector< std::vector< OutputShape > > & get_phi() const
std::vector< std::vector< OutputShape > > dphideta
Shape function derivatives in the eta direction.
std::vector< RealGradient > dweight
Used for certain infinite element families: the global derivative of the additional radial weight ...
const std::vector< std::vector< OutputTensor > > & get_dual_d2phi() const
const std::vector< std::vector< OutputShape > > & get_dphidz() const