libMesh
Public Member Functions | Public Attributes | List of all members
libMesh::FEType Class Reference

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

#include <fe_type.h>

Public Member Functions

 FEType (const int o=1, const FEFamily f=LAGRANGE)
 Constructor. More...
 
 FEType (const int o=1, const FEFamily f=LAGRANGE, const int ro=THIRD, const FEFamily rf=JACOBI_20_00, const InfMapType im=CARTESIAN)
 Constructor. More...
 
bool operator== (const FEType &f2) const
 Tests equality. More...
 
bool operator!= (const FEType &f2) const
 Tests inequality. More...
 
bool operator< (const FEType &f2) const
 An ordering to make FEType useful as a std::map key. More...
 
Order default_quadrature_order () const
 
std::unique_ptr< QBasedefault_quadrature_rule (const unsigned int dim, const int extraorder=0) const
 
Order unweighted_quadrature_order () const
 
std::unique_ptr< QBaseunweighted_quadrature_rule (const unsigned int dim, const int extraorder=0) const
 

Public Attributes

OrderWrapper order
 The approximation order of the element. More...
 
FEFamily family
 The type of finite element. More...
 
OrderWrapper radial_order
 The approximation order in radial direction of the infinite element. More...
 
FEFamily radial_family
 The type of approximation in radial direction. More...
 
InfMapType inf_map
 The coordinate mapping type of the infinite element. More...
 

Detailed Description

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

Author
Benjamin S. Kirk
Date
2002 Manages the family, order, etc. parameters for a given FE.

Definition at line 196 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 206 of file fe_type.h.

207  :
208  order(o),
209  family(f)
210  {}
FEFamily family
The type of finite element.
Definition: fe_type.h:221
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:215

◆ 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 234 of file fe_type.h.

238  :
239  order(o),
240  radial_order(ro),
241  family(f),
242  radial_family(rf),
243  inf_map(im)
244  {}
FEFamily family
The type of finite element.
Definition: fe_type.h:221
OrderWrapper radial_order
The approximation order in radial direction of the infinite element.
Definition: fe_type.h:254
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:215
InfMapType inf_map
The coordinate mapping type of the infinite element.
Definition: fe_type.h:275
FEFamily radial_family
The type of approximation in radial direction.
Definition: fe_type.h:267

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 on affine elements.

Definition at line 371 of file fe_type.h.

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

Referenced by alternative_fe_assembly(), LinearElasticity::assemble(), AssembleOptimization::assemble_A_and_F(), assemble_divgrad(), assemble_elasticity(), assemble_ellipticdg(), assemble_graddiv(), libMesh::ClawSystem::assemble_jump_coupling_matrix(), libMesh::ClawSystem::assemble_mass_matrix(), assemble_poisson(), assemble_shell(), assemble_stokes(), compute_jacobian(), libMesh::FEGenericBase< FEOutputType< T >::type >::compute_periodic_constraints(), libMesh::FEGenericBase< FEOutputType< T >::type >::compute_proj_constraints(), compute_residual(), compute_stresses(), LinearElasticityWithContact::compute_stresses(), LinearElasticity::compute_stresses(), LargeDeformationElasticity::compute_stresses(), default_quadrature_rule(), fe_assembly(), form_functionA(), form_functionB(), form_matrixA(), LargeDeformationElasticity::jacobian(), OverlappingCouplingFunctor::operator()(), periodic_bc_test_poisson(), libMesh::FE< Dim, LAGRANGE_VEC >::reinit_default_dual_shape_coeffs(), LargeDeformationElasticity::residual(), LinearElasticityWithContact::residual_and_jacobian(), DualShapeTest::setUp(), FETestBase< order, family, elem_type, 1 >::setUp(), InfFERadialTest::testInfQuants(), InfFERadialTest::testSides(), InfFERadialTest::testSingleOrder(), libMesh::Elem::true_centroid(), and libMesh::Elem::volume().

372 {
373  return static_cast<Order>(2*static_cast<unsigned int>(order.get_order()) + 1);
374 }
Order
defines an enum for polynomial orders.
Definition: enum_order.h:40
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:215
int get_order() const
Explicitly request the order as an int.
Definition: fe_type.h:80

◆ 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, with an additional power on the basis order to help account for nonlinearities and/or nonuniform coefficients. Higher or lower degree rules can be chosen by changing the extraorder parameter.

Definition at line 34 of file fe_type.C.

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

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::RBConstruction::add_scaled_matrix_and_vector(), assemble_biharmonic(), assemble_func(), assemble_laplace(), assemble_shell(), libMesh::System::calculate_norm(), libMesh::FEGenericBase< FEOutputType< T >::type >::coarsened_dof_values(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::BoundaryProjectSolution::operator()(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectEdges::operator()(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectSides::operator()(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectInteriors::operator()(), Biharmonic::JR::residual_and_jacobian(), libMesh::HPCoarsenTest::select_refinement(), FETest< order, family, elem_type >::testCustomReinit(), and libMesh::FEMContext::use_default_quadrature_rules().

36 {
37  // Clough elements have at least piecewise cubic functions
38  if (family == CLOUGH)
39  {
40  Order o = static_cast<Order>(std::max(static_cast<unsigned int>(this->default_quadrature_order()),
41  static_cast<unsigned int>(7 + extraorder)));
42  return std::make_unique<QClough>(dim, o);
43  }
44 
45  if (family == SUBDIVISION)
46  return std::make_unique<QGauss>(dim, static_cast<Order>(1 + extraorder));
47 
48  return std::make_unique<QGauss>(dim, static_cast<Order>(this->default_quadrature_order() + extraorder));
49 }
FEFamily family
The type of finite element.
Definition: fe_type.h:221
Order
defines an enum for polynomial orders.
Definition: enum_order.h:40
unsigned int dim
Order default_quadrature_order() const
Definition: fe_type.h:371

◆ operator!=()

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

Tests inequality.

Definition at line 297 of file fe_type.h.

298  {
299  return !(*this == f2);
300  }

◆ operator<()

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

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

Definition at line 305 of file fe_type.h.

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

306  {
307  if (order != f2.order)
308  return (order < f2.order);
309  if (family != f2.family)
310  return (family < f2.family);
311 
312 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
313  if (radial_order != f2.radial_order)
314  return (radial_order < f2.radial_order);
315  if (radial_family != f2.radial_family)
316  return (radial_family < f2.radial_family);
317  if (inf_map != f2.inf_map)
318  return (inf_map < f2.inf_map);
319 #endif // ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
320  return false;
321  }
FEFamily family
The type of finite element.
Definition: fe_type.h:221
OrderWrapper radial_order
The approximation order in radial direction of the infinite element.
Definition: fe_type.h:254
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:215
InfMapType inf_map
The coordinate mapping type of the infinite element.
Definition: fe_type.h:275
FEFamily radial_family
The type of approximation in radial direction.
Definition: fe_type.h:267

◆ operator==()

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

Tests equality.

Definition at line 282 of file fe_type.h.

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

283  {
284  return (order == f2.order
285  && family == f2.family
286 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
287  && radial_order == f2.radial_order
288  && radial_family == f2.radial_family
289  && inf_map == f2.inf_map
290 #endif // ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
291  );
292  }
FEFamily family
The type of finite element.
Definition: fe_type.h:221
OrderWrapper radial_order
The approximation order in radial direction of the infinite element.
Definition: fe_type.h:254
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:215
InfMapType inf_map
The coordinate mapping type of the infinite element.
Definition: fe_type.h:275
FEFamily radial_family
The type of approximation in radial direction.
Definition: fe_type.h:267

◆ unweighted_quadrature_order()

Order libMesh::FEType::unweighted_quadrature_order ( ) const
inline
Returns
The default quadrature order for integrating unweighted basis functions of this FEType. The unweighted quadrature order is calculated assuming a polynomial of degree order and is based on integrating the shape functions for such an element exactly on affine elements.

Definition at line 377 of file fe_type.h.

References order.

Referenced by unweighted_quadrature_rule().

378 {
379  return order;
380 }
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:215

◆ unweighted_quadrature_rule()

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

Definition at line 53 of file fe_type.C.

References libMesh::CLOUGH, dim, family, libMesh::SUBDIVISION, and unweighted_quadrature_order().

Referenced by libMesh::JumpErrorEstimator::estimate_error(), and libMesh::FEMContext::use_unweighted_quadrature_rules().

55 {
56  // Clough elements have at least piecewise cubic functions
57  if (family == CLOUGH)
58  {
59  Order o = static_cast<Order>(std::max(static_cast<unsigned int>(this->unweighted_quadrature_order()),
60  static_cast<unsigned int>(3 + extraorder)));
61  return std::make_unique<QClough>(dim, o);
62  }
63 
64  if (family == SUBDIVISION)
65  return std::make_unique<QGauss>(dim, static_cast<Order>(1 + extraorder));
66 
67  return std::make_unique<QGauss>(dim, static_cast<Order>(this->unweighted_quadrature_order() + extraorder));
68 }
FEFamily family
The type of finite element.
Definition: fe_type.h:221
Order
defines an enum for polynomial orders.
Definition: enum_order.h:40
Order unweighted_quadrature_order() const
Definition: fe_type.h:377
unsigned int dim

Member Data Documentation

◆ family

FEFamily libMesh::FEType::family

The type of finite element.

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 221 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_local_dofs_node_major(), libMesh::DofMap::distribute_local_dofs_var_major(), libMesh::DofMap::distribute_scalar_dofs(), libMesh::DofMap::dof_indices(), libMesh::FEInterface::extra_hanging_dofs(), TIMPI::OpFunction< libMesh::FEType >::fetype_max(), TIMPI::OpFunction< libMesh::FEType >::fetype_min(), libMesh::FEInterface::field_type(), libMesh::FEMContext::find_hardest_fe_type(), libMesh::FEInterface::get_continuity(), libMesh::FEAbstract::get_family(), libMesh::System::get_info(), libMesh::InfFE< Dim, T_radial, T_map >::inf_compute_constraints(), libMesh::FEMContext::init_internal_data(), libMesh::FEInterface::is_hierarchic(), main(), libMesh::FEInterface::max_order(), libMesh::Variable::n_components(), libMesh::FEInterface::n_vec_dim(), libMesh::DifferentiablePhysics::nonlocal_mass_residual(), libMesh::DofMap::old_dof_indices(), libMesh::BoundaryProjectSolution::operator()(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SortAndCopy::operator()(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectVertices::operator()(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectEdges::operator()(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectSides::operator()(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectInteriors::operator()(), std::hash< libMesh::FEType >::operator()(), operator<(), operator==(), libMesh::DofMap::process_mesh_constraint_rows(), libMesh::System::project_vector(), libMesh::System::read_header(), libMesh::System::read_parallel_data(), libMesh::StaticCondensationDofMap::reinit(), libMesh::DofMap::reinit(), TIMPI::Attributes< libMesh::FEType >::set_highest(), TIMPI::Attributes< libMesh::FEType >::set_lowest(), libMesh::PetscDMWrapper::set_point_range_in_section(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SubFunctor::SubFunctor(), libMesh::MeshFunctionSolutionTransfer::transfer(), unweighted_quadrature_rule(), libMesh::System::write_header(), libMesh::Nemesis_IO_Helper::write_nodal_solution(), libMesh::System::write_parallel_data(), and libMesh::System::write_serialized_vector().

◆ 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 275 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 the base of the infinite element.

Definition at line 215 of file fe_type.h.

Referenced by libMesh::PetscDMWrapper::add_dofs_to_section(), libMesh::FEMContext::build_new_fe(), libMesh::FEGenericBase< FEOutputType< T >::type >::coarsened_dof_values(), libMesh::FEAbstract::compute_node_constraints(), libMesh::FEGenericBase< FEOutputType< T >::type >::compute_proj_constraints(), libMesh::GMVIO::copy_nodal_solution(), default_quadrature_order(), libMesh::DofMap::distribute_dofs(), libMesh::DofMap::distribute_scalar_dofs(), libMesh::DofMap::dof_indices(), libMesh::FEInterface::dofs_on_edge(), libMesh::FEInterface::dofs_on_side(), TIMPI::OpFunction< libMesh::FEType >::fetype_max(), TIMPI::OpFunction< libMesh::FEType >::fetype_min(), libMesh::FEMContext::find_hardest_fe_type(), libMesh::System::get_info(), libMesh::FEAbstract::get_order(), main(), 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::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::SortAndCopy::operator()(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectEdges::operator()(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectSides::operator()(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::ProjectInteriors::operator()(), std::hash< libMesh::FEType >::operator()(), operator<(), operator==(), libMesh::System::read_header(), libMesh::System::read_SCALAR_dofs(), libMesh::DofMap::reinit(), libMesh::DofMap::SCALAR_dof_indices(), libMesh::FEAbstract::set_fe_order(), TIMPI::Attributes< libMesh::FEType >::set_highest(), TIMPI::Attributes< libMesh::FEType >::set_lowest(), libMesh::FE< Dim, LAGRANGE_VEC >::shape(), libMesh::InfFE< Dim, T_radial, T_map >::shape(), libMesh::FEInterface::shape(), libMesh::InfFE< Dim, T_radial, T_map >::shape_deriv(), libMesh::FE< Dim, LAGRANGE_VEC >::shape_deriv(), libMesh::FEInterface::shape_deriv(), libMesh::FE< Dim, LAGRANGE_VEC >::shape_second_deriv(), libMesh::FEInterface::shape_second_deriv(), libMesh::FEInterface::shapes(), libMesh::FEInterface::side_nodal_soln(), unweighted_quadrature_order(), and libMesh::System::write_header().

◆ radial_family

FEFamily libMesh::FEType::radial_family

◆ radial_order

OrderWrapper libMesh::FEType::radial_order

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