$extrastylesheet
libMesh::FEType Class Reference

#include <fe_type.h>

List of all members.

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
UniquePtr< 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.

                                      :
    order(o),
    family(f)
  {}
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.

                                          :
    order(o),
    radial_order(ro),
    family(f),
    radial_family(rf),
    inf_map(im)
  {}

Member Function Documentation

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< OutputType >::compute_periodic_constraints(), libMesh::FEGenericBase< OutputType >::compute_proj_constraints(), default_quadrature_rule(), libMesh::JumpErrorEstimator::estimate_error(), and libMesh::Elem::volume().

{
  return static_cast<Order>(2*static_cast<unsigned int>(order) + 1);
}
UniquePtr< 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< OutputType >::coarsened_dof_values(), libMesh::ExactErrorEstimator::estimate_error(), libMesh::FEMContext::FEMContext(), libMesh::WeightedPatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::PatchRecoveryErrorEstimator::EstimateError::operator()(), libMesh::ProjectSolution::operator()(), libMesh::BoundaryProjectSolution::operator()(), and libMesh::HPCoarsenTest::select_refinement().

{

  // Clough elements have at least piecewise cubic functions
  if (family == CLOUGH)
    {
      // this seems ridiculous but for some reason gcc 3.3.5 wants
      // this when using complex numbers (spetersen 04/20/06)
      const unsigned int seven = 7;

      return UniquePtr<QBase>
        (new QClough(dim,
                     static_cast<Order>
                     (std::max(static_cast<unsigned int>
                               (this->default_quadrature_order()), seven + extraorder))));
    }

  if (family == SUBDIVISION)
    return UniquePtr<QBase>
      (new QGauss(dim, static_cast<Order>(1 + extraorder)));

  return UniquePtr<QBase>
    (new QGauss(dim, static_cast<Order>(this->default_quadrature_order()
                                        + extraorder)));
}
bool libMesh::FEType::operator!= ( const FEType f2) const [inline]

Tests inequality

Definition at line 147 of file fe_type.h.

  {
    return !(*this == f2);
  }
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.

  {
    if (order != f2.order)
      return (order < f2.order);
    if (family != f2.family)
      return (family < f2.family);

#ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
    if (radial_order != f2.radial_order)
      return (radial_order < f2.radial_order);
    if (radial_family != f2.radial_family)
      return (radial_family < f2.radial_family);
    if (inf_map != f2.inf_map)
      return (inf_map < f2.inf_map);
#endif // ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
    return false;
  }
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.

  {
    return (order == f2.order
            && family == f2.family
#ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
            && radial_order == f2.radial_order
            && radial_family == f2.radial_family
            && inf_map == f2.inf_map
#endif // ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
            );
  }

Member Data Documentation

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

Referenced by libMesh::DofMap::_dof_indices(), libMesh::DofMap::add_neighbors_to_send_list(), libMesh::FEMSystem::assembly(), libMesh::FEMap::build(), libMesh::FETransformationBase< OutputShape >::build(), libMesh::FEAbstract::build(), libMesh::FEGenericBase< OutputType >::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::FEInterface::extra_hanging_dofs(), libMesh::FEMContext::FEMContext(), libMesh::FEInterface::field_type(), libMesh::FEAbstract::get_family(), libMesh::System::get_info(), libMesh::FEInterface::max_order(), libMesh::Variable::n_components(), libMesh::FEInterface::n_vec_dim(), libMesh::DofMap::old_dof_indices(), libMesh::ProjectSolution::operator()(), libMesh::ProjectFEMSolution::operator()(), libMesh::BoundaryProjectSolution::operator()(), operator<(), operator==(), libMesh::System::project_vector(), libMesh::System::read_header(), libMesh::System::read_parallel_data(), libMesh::System::read_serialized_vectors(), libMesh::DofMap::reinit(), libMesh::System::write_header(), libMesh::System::write_parallel_data(), libMesh::System::write_serialized_vector(), and libMesh::System::write_serialized_vectors().

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< OutputType >::build_InfFE(), libMesh::System::get_info(), libMesh::FEInterface::ifem_inverse_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().


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