libMesh::FEType Class Reference

Manages the family, order, etc. parameters for a given FE. More...

#include <fe_type.h>

Public Member Functions

 FEType (const int o=1, const FEFamily f=LAGRANGE)
 
 FEType (const int o=1, const FEFamily f=LAGRANGE, const int ro=THIRD, const FEFamily rf=JACOBI_20_00, const InfMapType im=CARTESIAN)
 
bool operator== (const FEType &f2) const
 
bool operator!= (const FEType &f2) const
 
bool operator< (const FEType &f2) const
 
Order default_quadrature_order () const
 
std::unique_ptr< QBasedefault_quadrature_rule (const unsigned int dim, const int extraorder=0) const
 

Public Attributes

OrderWrapper order
 
FEFamily family
 
OrderWrapper radial_order
 
FEFamily radial_family
 
InfMapType inf_map
 

Detailed Description

Manages the family, order, etc. parameters for a given FE.

class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialized finite element families.

Author
Benjamin S. Kirk
Date
2002

Definition at line 179 of file fe_type.h.

Constructor & Destructor Documentation

◆ FEType() [1/2]

libMesh::FEType::FEType ( const int  o = 1,
const FEFamily  f = LAGRANGE 
)
inline

Constructor. Optionally takes the approximation Order and the finite element family FEFamily

Definition at line 189 of file fe_type.h.

190  :
191  order(o),
192  family(f)
193  {}
FEFamily family
Definition: fe_type.h:204
OrderWrapper order
Definition: fe_type.h:198

◆ FEType() [2/2]

libMesh::FEType::FEType ( const int  o = 1,
const FEFamily  f = LAGRANGE,
const int  ro = THIRD,
const FEFamily  rf = JACOBI_20_00,
const InfMapType  im = CARTESIAN 
)
inline

Constructor. Optionally takes the approximation Order and the finite element family FEFamily.

Note
For non-infinite elements, the order and base order are the same, as with the family and base_family. It must be so, otherwise what we switch on would change when infinite elements are not compiled in.

Definition at line 217 of file fe_type.h.

221  :
222  order(o),
223  radial_order(ro),
224  family(f),
225  radial_family(rf),
226  inf_map(im)
227  {}
FEFamily family
Definition: fe_type.h:204
OrderWrapper radial_order
Definition: fe_type.h:237
OrderWrapper order
Definition: fe_type.h:198
InfMapType inf_map
Definition: fe_type.h:258
FEFamily radial_family
Definition: fe_type.h:250

Member Function Documentation

◆ default_quadrature_order()

Order libMesh::FEType::default_quadrature_order ( ) const
inline
Returns
The default quadrature order for this FEType. The default quadrature order is calculated assuming a polynomial of degree order and is based on integrating the mass matrix for such an element exactly.

Definition at line 333 of file fe_type.h.

References libMesh::OrderWrapper::get_order(), and order.

Referenced by libMesh::FEGenericBase< FEOutputType< T >::type >::compute_periodic_constraints(), libMesh::FEGenericBase< FEOutputType< T >::type >::compute_proj_constraints(), default_quadrature_rule(), and libMesh::Elem::volume().

334 {
335  return static_cast<Order>(2*static_cast<unsigned int>(order.get_order()) + 1);
336 }
OrderWrapper order
Definition: fe_type.h:198
int get_order() const
Definition: fe_type.h:78

◆ default_quadrature_rule()

std::unique_ptr< QBase > libMesh::FEType::default_quadrature_rule ( const unsigned int  dim,
const int  extraorder = 0 
) const
Returns
A quadrature rule of appropriate type and order for this FEType. The default quadrature rule is based on integrating the mass matrix for such an element exactly. Higher or lower degree rules can be chosen by changing the extraorder parameter.

Definition at line 31 of file fe_type.C.

References libMesh::CLOUGH, default_quadrature_order(), family, std::max(), and libMesh::SUBDIVISION.

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::System::calculate_norm(), libMesh::FEGenericBase< FEOutputType< T >::type >::coarsened_dof_values(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::FEMContext::init_internal_data(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::operator()(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::BoundaryProjectSolution::operator()(), and libMesh::HPCoarsenTest::select_refinement().

33 {
34  // Clough elements have at least piecewise cubic functions
35  if (family == CLOUGH)
36  {
37  Order o = static_cast<Order>(std::max(static_cast<unsigned int>(this->default_quadrature_order()),
38  static_cast<unsigned int>(7 + extraorder)));
39  return libmesh_make_unique<QClough>(dim, o);
40  }
41 
42  if (family == SUBDIVISION)
43  return libmesh_make_unique<QGauss>(dim, static_cast<Order>(1 + extraorder));
44 
45  return libmesh_make_unique<QGauss>(dim, static_cast<Order>(this->default_quadrature_order() + extraorder));
46 }
FEFamily family
Definition: fe_type.h:204
Order default_quadrature_order() const
Definition: fe_type.h:333
long double max(long double a, double b)

◆ operator!=()

bool libMesh::FEType::operator!= ( const FEType f2) const
inline

Tests inequality

Definition at line 280 of file fe_type.h.

281  {
282  return !(*this == f2);
283  }

◆ operator<()

bool libMesh::FEType::operator< ( const FEType f2) const
inline

An ordering to make FEType useful as a std::map key

Definition at line 288 of file fe_type.h.

References family, inf_map, order, radial_family, and radial_order.

289  {
290  if (order != f2.order)
291  return (order < f2.order);
292  if (family != f2.family)
293  return (family < f2.family);
294 
295 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
296  if (radial_order != f2.radial_order)
297  return (radial_order < f2.radial_order);
298  if (radial_family != f2.radial_family)
299  return (radial_family < f2.radial_family);
300  if (inf_map != f2.inf_map)
301  return (inf_map < f2.inf_map);
302 #endif // ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
303  return false;
304  }
FEFamily family
Definition: fe_type.h:204
OrderWrapper radial_order
Definition: fe_type.h:237
OrderWrapper order
Definition: fe_type.h:198
InfMapType inf_map
Definition: fe_type.h:258
FEFamily radial_family
Definition: fe_type.h:250

◆ operator==()

bool libMesh::FEType::operator== ( const FEType f2) const
inline

Tests equality

Definition at line 265 of file fe_type.h.

References family, inf_map, order, radial_family, and radial_order.

266  {
267  return (order == f2.order
268  && family == f2.family
269 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
270  && radial_order == f2.radial_order
271  && radial_family == f2.radial_family
272  && inf_map == f2.inf_map
273 #endif // ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
274  );
275  }
FEFamily family
Definition: fe_type.h:204
OrderWrapper radial_order
Definition: fe_type.h:237
OrderWrapper order
Definition: fe_type.h:198
InfMapType inf_map
Definition: fe_type.h:258
FEFamily radial_family
Definition: fe_type.h:250

Member Data Documentation

◆ family

FEFamily libMesh::FEType::family

The type of finite element. Valid types are LAGRANGE, HIERARCHIC, etc...

The type of approximation in radial direction. Valid types are JACOBI_20_00, JACOBI_30_00, etc...

Definition at line 204 of file fe_type.h.

Referenced by libMesh::DofMap::_dof_indices(), libMesh::FEMSystem::assembly(), libMesh::FETransformationBase< OutputShape >::build(), libMesh::FEMap::build(), libMesh::FEAbstract::build(), libMesh::FEGenericBase< FEOutputType< T >::type >::build(), libMesh::FEMContext::build_new_fe(), libMesh::FEInterface::compute_constraints(), libMesh::GMVIO::copy_nodal_solution(), default_quadrature_rule(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DofMap::dof_indices(), libMesh::JumpErrorEstimator::estimate_error(), libMesh::FEInterface::extra_hanging_dofs(), libMesh::FEInterface::field_type(), libMesh::FEAbstract::get_family(), libMesh::System::get_info(), libMesh::FEMContext::init_internal_data(), libMesh::FEInterface::max_order(), libMesh::Variable::n_components(), libMesh::FEInterface::n_vec_dim(), libMesh::DofMap::old_dof_indices(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::operator()(), libMesh::BoundaryProjectSolution::operator()(), operator<(), operator==(), libMesh::System::read_header(), libMesh::System::read_parallel_data(), libMesh::System::read_serialized_vectors(), libMesh::DofMap::reinit(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::System::write_header(), libMesh::System::write_parallel_data(), libMesh::System::write_serialized_vector(), and libMesh::System::write_serialized_vectors().

◆ inf_map

InfMapType libMesh::FEType::inf_map

The coordinate mapping type of the infinite element. When the infinite elements are defined over a surface with a separable coordinate system (sphere, spheroid, ellipsoid), the infinite elements may take advantage of this fact.

Definition at line 258 of file fe_type.h.

Referenced by libMesh::FEGenericBase< FEOutputType< T >::type >::build_InfFE(), libMesh::System::get_info(), libMesh::FEInterface::ifem_inverse_map(), libMesh::FEInterface::ifem_map(), libMesh::FEInterface::ifem_nodal_soln(), libMesh::InfFE< Dim, T_radial, T_map >::InfFE(), operator<(), operator==(), libMesh::System::read_header(), and libMesh::System::write_header().

◆ order

OrderWrapper libMesh::FEType::order

The approximation order of the element.

The approximation order in radial direction of the infinite element.

Definition at line 198 of file fe_type.h.

Referenced by libMesh::DofMap::_dof_indices(), libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::FEMContext::build_new_fe(), libMesh::FEGenericBase< FEOutputType< T >::type >::coarsened_dof_values(), libMesh::FEInterface::compute_data(), libMesh::FEGenericBase< FEOutputType< T >::type >::compute_proj_constraints(), libMesh::DofMap::constrain_p_dofs(), libMesh::GMVIO::copy_nodal_solution(), default_quadrature_order(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DofMap::dof_indices(), libMesh::FEInterface::dofs_on_edge(), libMesh::FEInterface::dofs_on_side(), libMesh::System::get_info(), libMesh::FEAbstract::get_order(), libMesh::FEMContext::init_internal_data(), libMesh::FESubdivision::init_shape_functions(), libMesh::Variable::n_components(), libMesh::FEInterface::n_dofs(), libMesh::FEInterface::n_dofs_at_node(), libMesh::FEInterface::n_dofs_per_elem(), libMesh::FEInterface::n_shape_functions(), libMesh::FEInterface::nodal_soln(), libMesh::DofMap::old_dof_indices(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::operator()(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), operator<(), operator==(), libMesh::System::read_header(), libMesh::System::read_SCALAR_dofs(), libMesh::DofMap::reinit(), libMesh::DofMap::SCALAR_dof_indices(), libMesh::HPCoarsenTest::select_refinement(), libMesh::FEAbstract::set_fe_order(), libMesh::FEInterface::shape(), and libMesh::System::write_header().

◆ radial_family

FEFamily libMesh::FEType::radial_family

For InfFE, family contains the radial shape family, while base_family contains the approximation type in circumferential direction. Valid types are LAGRANGE, HIERARCHIC, etc...

Definition at line 250 of file fe_type.h.

Referenced by libMesh::FEGenericBase< FEOutputType< T >::type >::build_InfFE(), libMesh::System::get_info(), libMesh::FEInterface::ifem_compute_data(), libMesh::FEInterface::ifem_nodal_soln(), libMesh::FEInterface::ifem_shape(), libMesh::InfFE< Dim, T_radial, T_map >::InfFE(), operator<(), operator==(), libMesh::System::read_header(), and libMesh::System::write_header().

◆ radial_order


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