libMesh::FEType Class Reference

#include <fe_type.h>

Public Member Functions

 FEType (const Order o=FIRST, const FEFamily f=LAGRANGE)
 
 FEType (const Order o=FIRST, const FEFamily f=LAGRANGE, const Order 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
 
AutoPtr< QBasedefault_quadrature_rule (const unsigned int dim, const int extraorder=0) const
 

Public Attributes

Order order
 
FEFamily family
 
Order radial_order
 
FEFamily radial_family
 
InfMapType inf_map
 

Detailed Description

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

Definition at line 43 of file fe_type.h.

Constructor & Destructor Documentation

libMesh::FEType::FEType ( const Order  o = FIRST,
const FEFamily  f = LAGRANGE 
)
inline

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

Definition at line 53 of file fe_type.h.

54  :
55  order(o),
56  family(f)
57  {}
libMesh::FEType::FEType ( const Order  o = FIRST,
const FEFamily  f = LAGRANGE,
const Order  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 that 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 83 of file fe_type.h.

87  :
88  order(o),
89  radial_order(ro),
90  family(f),
91  radial_family(rf),
92  inf_map(im)
93  {}

Member Function Documentation

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

References order.

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

201 {
202  return static_cast<Order>(2*static_cast<unsigned int>(order) + 1);
203 }
AutoPtr< 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 libMeshEnums::CLOUGH, default_quadrature_order(), family, and std::max().

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::UniformRefinementEstimator::_estimate_error(), libMesh::System::calculate_norm(), libMesh::FEGenericBase< T >::coarsened_dof_values(), libMesh::FEMContext::FEMContext(), libMesh::for(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::ProjectSolution::operator()(), libMesh::BoundaryProjectSolution::operator()(), and libMesh::HPCoarsenTest::select_refinement().

33 {
34 
35  // Clough elements have at least piecewise cubic functions
36  if (family == CLOUGH)
37  {
38  // this seems ridiculous but for some reason gcc 3.3.5 wants
39  // this when using complex numbers (spetersen 04/20/06)
40  const unsigned int seven = 7;
41 
42  return AutoPtr<QBase>
43  (new QClough(dim,
44  static_cast<Order>
45  (std::max(static_cast<unsigned int>
46  (this->default_quadrature_order()), seven + extraorder))));
47  }
48 
49  return AutoPtr<QBase>
50  (new QGauss(dim, static_cast<Order>(this->default_quadrature_order()
51  + extraorder)));
52 }
bool libMesh::FEType::operator!= ( const FEType f2) const
inline

Tests inequality

Definition at line 147 of file fe_type.h.

148  {
149  return !(*this == f2);
150  }
bool libMesh::FEType::operator< ( const FEType f2) const
inline

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

Definition at line 155 of file fe_type.h.

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

156  {
157  if (order != f2.order)
158  return (order < f2.order);
159  if (family != f2.family)
160  return (family < f2.family);
161 
162 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
163  if (radial_order != f2.radial_order)
164  return (radial_order < f2.radial_order);
165  if (radial_family != f2.radial_family)
166  return (radial_family < f2.radial_family);
167  if (inf_map != f2.inf_map)
168  return (inf_map < f2.inf_map);
169 #endif // ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
170  return false;
171  }
bool libMesh::FEType::operator== ( const FEType f2) const
inline

Tests equality

Definition at line 132 of file fe_type.h.

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

133  {
134  return (order == f2.order
135  && family == f2.family
136 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
137  && radial_order == f2.radial_order
138  && radial_family == f2.radial_family
139  && inf_map == f2.inf_map
140 #endif // ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
141  );
142  }

Member Data Documentation

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

Referenced by libMesh::FEGenericBase< T >::build_InfFE(), libMesh::System::get_info(), libMesh::FEInterface::ifem_inverse_map(), libMesh::FEInterface::ifem_nodal_soln(), libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::InfFE(), operator<(), operator==(), libMesh::System::read_header(), and libMesh::System::write_header().

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

Referenced by libMesh::FEGenericBase< T >::build_InfFE(), libMesh::System::get_info(), libMesh::FEInterface::ifem_compute_data(), libMesh::FEInterface::ifem_nodal_soln(), libMesh::FEInterface::ifem_shape(), libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::InfFE(), operator<(), operator==(), libMesh::System::read_header(), and libMesh::System::write_header().


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

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:24 UTC

Hosted By:
SourceForge.net Logo