libMesh::FEInterface Class Reference

Interface class which provides access to FE functions. More...

#include <fe_interface.h>

Public Types

typedef unsigned int(* n_dofs_at_node_ptr) (const ElemType, const Order, const unsigned int)
 

Public Member Functions

virtual ~FEInterface ()
 
template<>
void shape (const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p, Real &phi)
 
template<>
void shape (const unsigned int dim, const FEType &fe_t, const Elem *elem, const unsigned int i, const Point &p, Real &phi)
 
template<>
void shape (const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p, RealGradient &phi)
 
template<>
void shape (const unsigned int dim, const FEType &fe_t, const Elem *elem, const unsigned int i, const Point &p, RealGradient &phi)
 

Static Public Member Functions

static unsigned int n_shape_functions (const unsigned int dim, const FEType &fe_t, const ElemType t)
 
static unsigned int n_dofs (const unsigned int dim, const FEType &fe_t, const ElemType t)
 
static unsigned int n_dofs_at_node (const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
 
static n_dofs_at_node_ptr n_dofs_at_node_function (const unsigned int dim, const FEType &fe_t)
 
static unsigned int n_dofs_per_elem (const unsigned int dim, const FEType &fe_t, const ElemType t)
 
static void dofs_on_side (const Elem *const elem, const unsigned int dim, const FEType &fe_t, unsigned int s, std::vector< unsigned int > &di)
 
static void dofs_on_edge (const Elem *const elem, const unsigned int dim, const FEType &fe_t, unsigned int e, std::vector< unsigned int > &di)
 
static void nodal_soln (const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
 
static Point map (unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p)
 
static Point inverse_map (const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
 
static void inverse_map (const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Point > &physical_points, std::vector< Point > &reference_points, const Real tolerance=TOLERANCE, const bool secure=true)
 
static bool on_reference_element (const Point &p, const ElemType t, const Real eps=TOLERANCE)
 
static Real shape (const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
 
static Real shape (const unsigned int dim, const FEType &fe_t, const Elem *elem, const unsigned int i, const Point &p)
 
template<typename OutputType >
static void shape (const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p, OutputType &phi)
 
template<typename OutputType >
static void shape (const unsigned int dim, const FEType &fe_t, const Elem *elem, const unsigned int i, const Point &p, OutputType &phi)
 
static void compute_data (const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
 
static void compute_constraints (DofConstraints &constraints, DofMap &dof_map, const unsigned int variable_number, const Elem *elem)
 
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)
 
static unsigned int max_order (const FEType &fe_t, const ElemType &el_t)
 
static bool extra_hanging_dofs (const FEType &fe_t)
 
static FEFieldType field_type (const FEType &fe_type)
 
static FEFieldType field_type (const FEFamily &fe_family)
 
static unsigned int n_vec_dim (const MeshBase &mesh, const FEType &fe_type)
 

Private Member Functions

 FEInterface ()
 

Static Private Member Functions

static bool is_InfFE_elem (const ElemType et)
 
static unsigned int ifem_n_shape_functions (const unsigned int dim, const FEType &fe_t, const ElemType t)
 
static unsigned int ifem_n_dofs (const unsigned int dim, const FEType &fe_t, const ElemType t)
 
static unsigned int ifem_n_dofs_at_node (const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
 
static unsigned int ifem_n_dofs_per_elem (const unsigned int dim, const FEType &fe_t, const ElemType t)
 
static void ifem_nodal_soln (const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
 
static Point ifem_map (const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p)
 
static Point ifem_inverse_map (const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
 
static void ifem_inverse_map (const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Point > &physical_points, std::vector< Point > &reference_points, const Real tolerance=TOLERANCE, const bool secure=true)
 
static bool ifem_on_reference_element (const Point &p, const ElemType t, const Real eps)
 
static Real ifem_shape (const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
 
static Real ifem_shape (const unsigned int dim, const FEType &fe_t, const Elem *elem, const unsigned int i, const Point &p)
 
static void ifem_compute_data (const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
 

Detailed Description

Interface class which provides access to FE functions.

This class provides an encapsulated access to all static public member functions of finite element classes. Using this class, one need not worry about the correct finite element class.

Author
Daniel Dreyer
Date
2002-2007

Definition at line 63 of file fe_interface.h.

Member Typedef Documentation

typedef unsigned int(* libMesh::FEInterface::n_dofs_at_node_ptr) (const ElemType, const Order, const unsigned int)

Definition at line 113 of file fe_interface.h.

Constructor & Destructor Documentation

libMesh::FEInterface::FEInterface ( )
private

Empty constructor. Do not create an object of this type.

Definition at line 32 of file fe_interface.C.

33 {
34  libmesh_error_msg("ERROR: Do not define an object of this type.");
35 }
virtual libMesh::FEInterface::~FEInterface ( )
inlinevirtual

Destructor.

Definition at line 77 of file fe_interface.h.

References n_dofs(), n_dofs_at_node(), and n_shape_functions().

77 {}

Member Function Documentation

void libMesh::FEInterface::compute_constraints ( DofConstraints constraints,
DofMap dof_map,
const unsigned int  variable_number,
const Elem elem 
)
static

Computes the constraint matrix contributions (for non-conforming adapted meshes) corresponding to variable number var_number.

Definition at line 830 of file fe_interface.C.

References libMesh::BERNSTEIN, libMesh::CLOUGH, libMesh::FE< Dim, T >::compute_constraints(), libMesh::Elem::dim(), libMesh::FEType::family, libMesh::HERMITE, libMesh::HIERARCHIC, libMesh::L2_HIERARCHIC, libMesh::LAGRANGE, libMesh::LAGRANGE_VEC, libMesh::SZABAB, and libMesh::DofMap::variable_type().

834 {
835  libmesh_assert(elem);
836 
837  const FEType & fe_t = dof_map.variable_type(variable_number);
838 
839  switch (elem->dim())
840  {
841  case 0:
842  case 1:
843  {
844  // No constraints in 0D/1D.
845  return;
846  }
847 
848 
849  case 2:
850  {
851  switch (fe_t.family)
852  {
853  case CLOUGH:
855  dof_map,
856  variable_number,
857  elem); return;
858 
859  case HERMITE:
861  dof_map,
862  variable_number,
863  elem); return;
864 
865  case LAGRANGE:
867  dof_map,
868  variable_number,
869  elem); return;
870 
871  case HIERARCHIC:
873  dof_map,
874  variable_number,
875  elem); return;
876 
877  case L2_HIERARCHIC:
879  dof_map,
880  variable_number,
881  elem); return;
882 
883  case LAGRANGE_VEC:
885  dof_map,
886  variable_number,
887  elem); return;
888 
889 
890 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
891  case SZABAB:
893  dof_map,
894  variable_number,
895  elem); return;
896 
897  case BERNSTEIN:
899  dof_map,
900  variable_number,
901  elem); return;
902 
903 #endif
904  default:
905  return;
906  }
907  }
908 
909 
910  case 3:
911  {
912  switch (fe_t.family)
913  {
914  case HERMITE:
916  dof_map,
917  variable_number,
918  elem); return;
919 
920  case LAGRANGE:
922  dof_map,
923  variable_number,
924  elem); return;
925 
926  case HIERARCHIC:
928  dof_map,
929  variable_number,
930  elem); return;
931 
932  case L2_HIERARCHIC:
934  dof_map,
935  variable_number,
936  elem); return;
937 
938  case LAGRANGE_VEC:
940  dof_map,
941  variable_number,
942  elem); return;
943 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
944  case SZABAB:
946  dof_map,
947  variable_number,
948  elem); return;
949 
950  case BERNSTEIN:
952  dof_map,
953  variable_number,
954  elem); return;
955 
956 #endif
957  default:
958  return;
959  }
960  }
961 
962 
963  default:
964  libmesh_error_msg("Invalid dimension = " << elem->dim());
965  }
966 }
static void compute_constraints(DofConstraints &constraints, DofMap &dof_map, const unsigned int variable_number, const Elem *elem)
void libMesh::FEInterface::compute_data ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
FEComputeData data 
)
static

Lets the appropriate child of FEBase compute the requested data for the input specified in data, and sets the values in data. See this as a generalization of shape(). With infinite elements disabled, computes values for all shape functions of elem evaluated at p.

Note
On a p-refined element, fe_t.order should be the base order of the element.

Definition at line 794 of file fe_interface.C.

References ifem_compute_data(), libMesh::FEComputeData::init(), is_InfFE_elem(), n_dofs(), libMesh::FEType::order, libMesh::FEComputeData::p, libMesh::Elem::p_level(), libMesh::FEComputeData::shape, shape(), and libMesh::Elem::type().

Referenced by libMesh::MeshFunction::discontinuous_value(), libMesh::MeshFunction::operator()(), and libMesh::System::point_value().

798 {
799 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
800 
801  if (elem && is_InfFE_elem(elem->type()))
802  {
803  data.init();
804  ifem_compute_data(dim, fe_t, elem, data);
805  return;
806  }
807 
808 #endif
809 
810  FEType p_refined = fe_t;
811  p_refined.order = static_cast<Order>(p_refined.order + elem->p_level());
812 
813  const unsigned int n_dof = n_dofs (dim, p_refined, elem->type());
814  const Point & p = data.p;
815  data.shape.resize(n_dof);
816 
817  // set default values for all the output fields
818  data.init();
819 
820  for (unsigned int n=0; n<n_dof; n++)
821  data.shape[n] = shape(dim, p_refined, elem, n, p);
822 
823  return;
824 }
static unsigned int n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:411
static void ifem_compute_data(const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition: fe_interface.C:614
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
IterBase * data
void libMesh::FEInterface::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 
)
static

Computes the constraint matrix contributions (for periodic boundary conditions) corresponding to variable number var_number.

Definition at line 974 of file fe_interface.C.

References libMesh::FEGenericBase< OutputType >::compute_periodic_constraints().

981 {
982  // No element-specific optimizations currently exist
984  dof_map,
985  boundaries,
986  mesh,
987  point_locator,
988  variable_number,
989  elem);
990 }
MeshBase & mesh
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)
Definition: fe_base.C:1647
void libMesh::FEInterface::dofs_on_edge ( const Elem *const  elem,
const unsigned int  dim,
const FEType fe_t,
unsigned int  e,
std::vector< unsigned int > &  di 
)
static

Fills the vector di with the local degree of freedom indices associated with edge e of element elem Automatically decides which finite element class to use.

On a p-refined element, fe_t.order should be the base order of the element.

Definition at line 493 of file fe_interface.C.

References libMesh::FEType::order.

Referenced by libMesh::FEGenericBase< OutputType >::coarsened_dof_values(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::operator()(), and libMesh::BoundaryProjectSolution::operator()().

498 {
499  const Order o = fe_t.order;
500 
501  void_fe_with_vec_switch(dofs_on_edge(elem, o, e, di));
502 }
static void dofs_on_edge(const Elem *const elem, const unsigned int dim, const FEType &fe_t, unsigned int e, std::vector< unsigned int > &di)
Definition: fe_interface.C:493
void libMesh::FEInterface::dofs_on_side ( const Elem *const  elem,
const unsigned int  dim,
const FEType fe_t,
unsigned int  s,
std::vector< unsigned int > &  di 
)
static

Fills the vector di with the local degree of freedom indices associated with side s of element elem Automatically decides which finite element class to use.

On a p-refined element, fe_t.order should be the base order of the element.

Definition at line 480 of file fe_interface.C.

References libMesh::FEType::order.

Referenced by libMesh::FEGenericBase< OutputType >::coarsened_dof_values(), libMesh::FEGenericBase< OutputType >::compute_periodic_constraints(), libMesh::FEGenericBase< OutputType >::compute_proj_constraints(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::operator()(), and libMesh::BoundaryProjectSolution::operator()().

485 {
486  const Order o = fe_t.order;
487 
488  void_fe_with_vec_switch(dofs_on_side(elem, o, s, di));
489 }
static void dofs_on_side(const Elem *const elem, const unsigned int dim, const FEType &fe_t, unsigned int s, std::vector< unsigned int > &di)
Definition: fe_interface.C:480
bool libMesh::FEInterface::extra_hanging_dofs ( const FEType fe_t)
static
Returns
true if separate degrees of freedom must be allocated for vertex DoFs and edge/face DoFs at a hanging node.

Definition at line 1369 of file fe_interface.C.

References libMesh::BERNSTEIN, libMesh::CLOUGH, libMesh::FEType::family, libMesh::HERMITE, libMesh::HIERARCHIC, libMesh::L2_HIERARCHIC, libMesh::L2_LAGRANGE, libMesh::LAGRANGE, libMesh::LAGRANGE_VEC, libMesh::MONOMIAL, libMesh::NEDELEC_ONE, libMesh::SUBDIVISION, libMesh::SZABAB, and libMesh::XYZ.

Referenced by libMesh::DofMap::_dof_indices(), libMesh::DofMap::old_dof_indices(), and libMesh::DofMap::reinit().

1370 {
1371  switch (fe_t.family)
1372  {
1373  case LAGRANGE:
1374  case L2_LAGRANGE:
1375  case MONOMIAL:
1376  case L2_HIERARCHIC:
1377  case XYZ:
1378  case SUBDIVISION:
1379  case LAGRANGE_VEC:
1380  case NEDELEC_ONE:
1381  return false;
1382  case CLOUGH:
1383  case HERMITE:
1384  case HIERARCHIC:
1385 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
1386  case BERNSTEIN:
1387  case SZABAB:
1388 #endif
1389  default:
1390  return true;
1391  }
1392 }
FEFieldType libMesh::FEInterface::field_type ( const FEType fe_type)
static
FEFieldType libMesh::FEInterface::field_type ( const FEFamily fe_family)
static
Returns
The number of components of a vector-valued element. Scalar-valued elements return 1.

Definition at line 1399 of file fe_interface.C.

References libMesh::LAGRANGE_VEC, libMesh::NEDELEC_ONE, libMesh::TYPE_SCALAR, and libMesh::TYPE_VECTOR.

1400 {
1401  switch (fe_family)
1402  {
1403  case LAGRANGE_VEC:
1404  case NEDELEC_ONE:
1405  return TYPE_VECTOR;
1406  default:
1407  return TYPE_SCALAR;
1408  }
1409 }
void libMesh::FEInterface::ifem_compute_data ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
FEComputeData data 
)
staticprivate

Definition at line 808 of file fe_interface_inf_fe.C.

References libMesh::InfFE< Dim, T_radial, T_map >::compute_data(), libMesh::INFINITE_MAP, libMesh::JACOBI_20_00, libMesh::JACOBI_30_00, libMesh::LAGRANGE, libMesh::LEGENDRE, and libMesh::FEType::radial_family.

Referenced by compute_data().

812 {
813  switch (dim)
814  {
815  // 1D
816  case 1:
817  {
818  switch (fe_t.radial_family)
819  {
820  /*
821  * For no derivatives (and local coordinates, as
822  * given in \p p) the infinite element shapes
823  * are independent of mapping type
824  */
825  case INFINITE_MAP:
827  break;
828 
829  case JACOBI_20_00:
831  break;
832 
833  case JACOBI_30_00:
835  break;
836 
837  case LEGENDRE:
839  break;
840 
841  case LAGRANGE:
843  break;
844 
845  default:
846  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
847  }
848 
849  break;
850  }
851 
852 
853  // 2D
854  case 2:
855  {
856  switch (fe_t.radial_family)
857  {
858  case INFINITE_MAP:
860  break;
861 
862  case JACOBI_20_00:
864  break;
865 
866  case JACOBI_30_00:
868  break;
869 
870  case LEGENDRE:
872  break;
873 
874  case LAGRANGE:
876  break;
877 
878  default:
879  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
880  }
881 
882  break;
883  }
884 
885 
886  // 3D
887  case 3:
888  {
889  switch (fe_t.radial_family)
890  {
891  case INFINITE_MAP:
893  break;
894 
895  case JACOBI_20_00:
897  break;
898 
899  case JACOBI_30_00:
901  break;
902 
903  case LEGENDRE:
905  break;
906 
907  case LAGRANGE:
909  break;
910 
911  default:
912  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
913  }
914 
915  break;
916  }
917 
918 
919  default:
920  libmesh_error_msg("Invalid dim = " << dim);
921  break;
922  }
923 }
static void compute_data(const FEType &fe_t, const Elem *inf_elem, FEComputeData &data)
IterBase * data
Point libMesh::FEInterface::ifem_inverse_map ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
const Point p,
const Real  tolerance = TOLERANCE,
const bool  secure = true 
)
staticprivate

Definition at line 447 of file fe_interface_inf_fe.C.

References libMesh::CARTESIAN, libMesh::ELLIPSOIDAL, libMesh::FEType::inf_map, libMesh::InfFE< Dim, T_radial, T_map >::inverse_map(), and libMesh::SPHERICAL.

Referenced by inverse_map().

453 {
454  switch (dim)
455  {
456  // 1D
457  case 1:
458  {
459  switch (fe_t.inf_map)
460  {
461  case CARTESIAN:
462  return InfFE<1,JACOBI_20_00,CARTESIAN>::inverse_map(elem, p, tolerance, secure);
463 
464  case SPHERICAL:
465  case ELLIPSOIDAL:
466  libmesh_not_implemented_msg("ERROR: Spherical and Ellipsoidal IFEMs not (yet) implemented.");
467 
468  /*
469  case SPHERICAL:
470  return InfFE<1,JACOBI_20_00,SPHERICAL>::inverse_map(elem, p, tolerance);
471 
472  case ELLIPSOIDAL:
473  return InfFE<1,JACOBI_20_00,ELLIPSOIDAL>::inverse_map(elem, p, tolerance);
474  */
475 
476  default:
477  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
478  }
479  }
480 
481 
482  // 2D
483  case 2:
484  {
485  switch (fe_t.inf_map)
486  {
487  case CARTESIAN:
488  return InfFE<2,JACOBI_20_00,CARTESIAN>::inverse_map(elem, p, tolerance, secure);
489 
490  case SPHERICAL:
491  case ELLIPSOIDAL:
492  libmesh_not_implemented_msg("ERROR: Spherical and Ellipsoidal IFEMs not (yet) implemented.");
493 
494  /*
495  case SPHERICAL:
496  return InfFE<2,JACOBI_20_00,SPHERICAL>::inverse_map(elem, p, tolerance);
497 
498  case ELLIPSOIDAL:
499  return InfFE<2,JACOBI_20_00,ELLIPSOIDAL>::inverse_map(elem, p, tolerance);
500  */
501 
502  default:
503  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
504  }
505  }
506 
507 
508  // 3D
509  case 3:
510  {
511  switch (fe_t.inf_map)
512  {
513  case CARTESIAN:
514  return InfFE<3,JACOBI_20_00,CARTESIAN>::inverse_map(elem, p, tolerance, secure);
515 
516  case SPHERICAL:
517  case ELLIPSOIDAL:
518  libmesh_not_implemented_msg("ERROR: Spherical and Ellipsoidal IFEMs not (yet) implemented.");
519 
520  /*
521  case SPHERICAL:
522  return InfFE<3,JACOBI_20_00,SPHERICAL>::inverse_map(elem, p, tolerance);
523 
524  case ELLIPSOIDAL:
525  return InfFE<3,JACOBI_20_00,ELLIPSOIDAL>::inverse_map(elem, p, tolerance);
526  */
527 
528  default:
529  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
530  }
531  }
532 
533  default:
534  libmesh_error_msg("Invalid dim = " << dim);
535  }
536 }
static Point inverse_map(const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: inf_fe_map.C:89
void libMesh::FEInterface::ifem_inverse_map ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
const std::vector< Point > &  physical_points,
std::vector< Point > &  reference_points,
const Real  tolerance = TOLERANCE,
const bool  secure = true 
)
staticprivate

Definition at line 540 of file fe_interface_inf_fe.C.

References libMesh::CARTESIAN, libMesh::FEType::inf_map, and libMesh::InfFE< Dim, T_radial, T_map >::inverse_map().

547 {
548  switch (dim)
549  {
550  // 1D
551  case 1:
552  {
553  switch (fe_t.inf_map)
554  {
555  case CARTESIAN:
556  InfFE<1,JACOBI_20_00,CARTESIAN>::inverse_map(elem, physical_points, reference_points, tolerance, secure);
557  return;
558 
559  default:
560  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
561  }
562  }
563 
564 
565  // 2D
566  case 2:
567  {
568  switch (fe_t.inf_map)
569  {
570  case CARTESIAN:
571  InfFE<2,JACOBI_20_00,CARTESIAN>::inverse_map(elem, physical_points, reference_points, tolerance, secure);
572  return;
573 
574  default:
575  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
576  }
577  }
578 
579 
580  // 3D
581  case 3:
582  {
583  switch (fe_t.inf_map)
584  {
585  case CARTESIAN:
586  InfFE<3,JACOBI_20_00,CARTESIAN>::inverse_map(elem, physical_points, reference_points, tolerance, secure);
587  return;
588 
589  default:
590  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
591  }
592  }
593 
594  default:
595  libmesh_error_msg("Invalid dim = " << dim);
596  }
597 }
static Point inverse_map(const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: inf_fe_map.C:89
Point libMesh::FEInterface::ifem_map ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
const Point p 
)
staticprivate

Definition at line 416 of file fe_interface_inf_fe.C.

References libMesh::CARTESIAN, libMesh::ELLIPSOIDAL, libMesh::FEType::inf_map, libMesh::InfFE< Dim, T_radial, T_map >::map(), and libMesh::SPHERICAL.

Referenced by map().

420 {
421  switch (fe_t.inf_map)
422  {
423  case CARTESIAN:
424  {
425  switch (dim)
426  {
427  case 1:
428  return InfFE<1,JACOBI_20_00,CARTESIAN>::map(elem, p);
429  case 2:
430  return InfFE<2,JACOBI_20_00,CARTESIAN>::map(elem, p);
431  case 3:
432  return InfFE<3,JACOBI_20_00,CARTESIAN>::map(elem, p);
433  default:
434  libmesh_error_msg("Invalid dim = " << dim);
435  }
436  }
437  case SPHERICAL:
438  case ELLIPSOIDAL:
439  libmesh_not_implemented_msg("ERROR: Spherical and Ellipsoidal IFEMs not (yet) implemented.");
440  default:
441  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
442  }
443 }
static Point map(const Elem *inf_elem, const Point &reference_point)
Definition: inf_fe_map.C:40
unsigned int libMesh::FEInterface::ifem_n_dofs ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t 
)
staticprivate

Definition at line 71 of file fe_interface_inf_fe.C.

References libMesh::InfFE< Dim, T_radial, T_map >::n_dofs().

Referenced by n_dofs().

74 {
75  switch (dim)
76  {
77  // 1D
78  case 1:
79  /*
80  * Since InfFE<Dim,T_radial,T_map>::n_dofs(...)
81  * is actually independent of T_radial and T_map, we can use
82  * just any T_radial and T_map
83  */
85 
86  // 2D
87  case 2:
89 
90  // 3D
91  case 3:
93 
94  default:
95  libmesh_error_msg("Unsupported dim = " << dim);
96  }
97 }
static unsigned int n_dofs(const FEType &fet, const ElemType inf_elem_type)
Definition: inf_fe_static.C:55
unsigned int libMesh::FEInterface::ifem_n_dofs_at_node ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t,
const unsigned int  n 
)
staticprivate

Definition at line 102 of file fe_interface_inf_fe.C.

References libMesh::InfFE< Dim, T_radial, T_map >::n_dofs_at_node().

Referenced by n_dofs_at_node().

106 {
107  switch (dim)
108  {
109  // 1D
110  case 1:
111  /*
112  * Since InfFE<Dim,T_radial,T_map>::n_dofs_at_node(...)
113  * is actually independent of T_radial and T_map, we can use
114  * just any T_radial and T_map
115  */
117 
118  // 2D
119  case 2:
121 
122  // 3D
123  case 3:
125 
126  default:
127  libmesh_error_msg("Unsupported dim = " << dim);
128  }
129 }
static unsigned int n_dofs_at_node(const FEType &fet, const ElemType inf_elem_type, const unsigned int n)
Definition: inf_fe_static.C:72
unsigned int libMesh::FEInterface::ifem_n_dofs_per_elem ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t 
)
staticprivate

Definition at line 135 of file fe_interface_inf_fe.C.

References libMesh::InfFE< Dim, T_radial, T_map >::n_dofs_per_elem().

Referenced by n_dofs_per_elem().

138 {
139  switch (dim)
140  {
141  // 1D
142  case 1:
143  /*
144  * Since InfFE<Dim,T_radial,T_map>::n_dofs(...)
145  * is actually independent of T_radial and T_map, we can use
146  * just any T_radial and T_map
147  */
149 
150  // 2D
151  case 2:
153 
154  // 3D
155  case 3:
157 
158  default:
159  libmesh_error_msg("Unsupported dim = " << dim);
160  }
161 }
static unsigned int n_dofs_per_elem(const FEType &fet, const ElemType inf_elem_type)
unsigned int libMesh::FEInterface::ifem_n_shape_functions ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t 
)
staticprivate

Definition at line 39 of file fe_interface_inf_fe.C.

References libMesh::InfFE< Dim, T_radial, T_map >::n_shape_functions().

Referenced by n_shape_functions().

42 {
43  switch (dim)
44  {
45  // 1D
46  case 1:
47  /*
48  * Since InfFE<Dim,T_radial,T_map>::n_shape_functions(...)
49  * is actually independent of T_radial and T_map, we can use
50  * just any T_radial and T_map
51  */
53 
54  // 2D
55  case 2:
57 
58  // 3D
59  case 3:
61 
62  default:
63  libmesh_error_msg("Unsupported dim = " << dim);
64  }
65 }
virtual unsigned int n_shape_functions() const libmesh_override
Definition: inf_fe.h:459
void libMesh::FEInterface::ifem_nodal_soln ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
const std::vector< Number > &  elem_soln,
std::vector< Number > &  nodal_soln 
)
staticprivate

Definition at line 166 of file fe_interface_inf_fe.C.

References libMesh::CARTESIAN, libMesh::FEType::inf_map, libMesh::INFINITE_MAP, libMesh::JACOBI_20_00, libMesh::JACOBI_30_00, libMesh::LAGRANGE, libMesh::LEGENDRE, libMesh::InfFE< Dim, T_radial, T_map >::nodal_soln(), and libMesh::FEType::radial_family.

Referenced by nodal_soln().

171 {
172  switch (dim)
173  {
174 
175  // 1D
176  case 1:
177  {
178  switch (fe_t.radial_family)
179  {
180  case INFINITE_MAP:
181  libmesh_error_msg("ERROR: INFINITE_MAP is not a valid shape family for radial approximation.");
182 
183  case JACOBI_20_00:
184  {
185  switch (fe_t.inf_map)
186  {
187  case CARTESIAN:
188  {
190  break;
191  }
192  default:
193  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
194  }
195  break;
196  }
197 
198  case JACOBI_30_00:
199  {
200  switch (fe_t.inf_map)
201  {
202  case CARTESIAN:
203  {
205  break;
206  }
207  default:
208  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
209  }
210  break;
211  }
212 
213  case LEGENDRE:
214  {
215  switch (fe_t.inf_map)
216  {
217  case CARTESIAN:
218  {
219  InfFE<1,LEGENDRE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
220  break;
221  }
222  default:
223  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
224  }
225  break;
226  }
227 
228  case LAGRANGE:
229  {
230  switch (fe_t.inf_map)
231  {
232  case CARTESIAN:
233  {
234  InfFE<1,LAGRANGE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
235  break;
236  }
237  default:
238  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
239  }
240  break;
241  }
242 
243  default:
244  libmesh_error_msg("ERROR: Bad FEType.radial_family= " << fe_t.radial_family);
245  }
246 
247  break;
248  }
249 
250 
251 
252 
253  // 2D
254  case 2:
255  {
256  switch (fe_t.radial_family)
257  {
258  case INFINITE_MAP:
259  libmesh_error_msg("ERROR: INFINITE_MAP is not a valid shape family for radial approximation.");
260 
261  case JACOBI_20_00:
262  {
263  switch (fe_t.inf_map)
264  {
265  case CARTESIAN:
266  {
268  break;
269  }
270  default:
271  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
272  }
273  break;
274  }
275 
276  case JACOBI_30_00:
277  {
278  switch (fe_t.inf_map)
279  {
280  case CARTESIAN:
281  {
283  break;
284  }
285  default:
286  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
287  }
288  break;
289  }
290 
291  case LEGENDRE:
292  {
293  switch (fe_t.inf_map)
294  {
295  case CARTESIAN:
296  {
297  InfFE<2,LEGENDRE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
298  break;
299  }
300  default:
301  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
302  }
303  break;
304  }
305 
306  case LAGRANGE:
307  {
308  switch (fe_t.inf_map)
309  {
310  case CARTESIAN:
311  {
312  InfFE<2,LAGRANGE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
313  break;
314  }
315  default:
316  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
317  }
318  break;
319  }
320 
321  default:
322  libmesh_error_msg("ERROR: Bad FEType.radial_family= " << fe_t.radial_family);
323  }
324 
325  break;
326  }
327 
328 
329 
330 
331  // 3D
332  case 3:
333  {
334  switch (fe_t.radial_family)
335  {
336  case INFINITE_MAP:
337  libmesh_error_msg("ERROR: INFINITE_MAP is not a valid shape family for radial approximation.");
338 
339  case JACOBI_20_00:
340  {
341  switch (fe_t.inf_map)
342  {
343  case CARTESIAN:
344  {
346  break;
347  }
348  default:
349  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
350  }
351  break;
352  }
353 
354  case JACOBI_30_00:
355  {
356  switch (fe_t.inf_map)
357  {
358  case CARTESIAN:
359  {
361  break;
362  }
363  default:
364  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
365  }
366  break;
367  }
368 
369  case LEGENDRE:
370  {
371  switch (fe_t.inf_map)
372  {
373  case CARTESIAN:
374  {
375  InfFE<3,LEGENDRE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
376  break;
377  }
378  default:
379  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
380  }
381  break;
382  }
383 
384  case LAGRANGE:
385  {
386  switch (fe_t.inf_map)
387  {
388  case CARTESIAN:
389  {
390  InfFE<3,LAGRANGE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
391  break;
392  }
393  default:
394  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
395  }
396  break;
397  }
398 
399 
400 
401  default:
402  libmesh_error_msg("ERROR: Bad FEType.radial_family= " << fe_t.radial_family);
403  }
404 
405  break;
406  }
407 
408  default:
409  libmesh_error_msg("Invalid dim = " << dim);
410  }
411 }
static void nodal_soln(const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
Definition: fe_interface.C:507
static void nodal_soln(const FEType &fet, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
bool libMesh::FEInterface::ifem_on_reference_element ( const Point p,
const ElemType  t,
const Real  eps 
)
staticprivate

Definition at line 602 of file fe_interface_inf_fe.C.

References libMesh::FEAbstract::on_reference_element().

605 {
606  return FEBase::on_reference_element(p,t,eps);
607 }
static bool on_reference_element(const Point &p, const ElemType t, const Real eps=TOLERANCE)
Definition: fe_abstract.C:554
Real libMesh::FEInterface::ifem_shape ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t,
const unsigned int  i,
const Point p 
)
staticprivate

Definition at line 612 of file fe_interface_inf_fe.C.

References libMesh::INFINITE_MAP, libMesh::JACOBI_20_00, libMesh::JACOBI_30_00, libMesh::LAGRANGE, libMesh::LEGENDRE, libMesh::FEType::radial_family, and libMesh::InfFE< Dim, T_radial, T_map >::shape().

Referenced by shape().

617 {
618  switch (dim)
619  {
620  // 1D
621  case 1:
622  {
623  switch (fe_t.radial_family)
624  {
625  /*
626  * For no derivatives (and local coordinates, as
627  * given in \p p) the infinite element shapes
628  * are independent of mapping type
629  */
630  case INFINITE_MAP:
631  return InfFE<1,INFINITE_MAP,CARTESIAN>::shape(fe_t, t, i, p);
632 
633  case JACOBI_20_00:
634  return InfFE<1,JACOBI_20_00,CARTESIAN>::shape(fe_t, t, i, p);
635 
636  case JACOBI_30_00:
637  return InfFE<1,JACOBI_30_00,CARTESIAN>::shape(fe_t, t, i, p);
638 
639  case LEGENDRE:
640  return InfFE<1,LEGENDRE,CARTESIAN>::shape(fe_t, t, i, p);
641 
642  case LAGRANGE:
643  return InfFE<1,LAGRANGE,CARTESIAN>::shape(fe_t, t, i, p);
644 
645  default:
646  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
647  }
648  }
649 
650 
651  // 2D
652  case 2:
653  {
654  switch (fe_t.radial_family)
655  {
656  case INFINITE_MAP:
657  return InfFE<2,INFINITE_MAP,CARTESIAN>::shape(fe_t, t, i, p);
658 
659  case JACOBI_20_00:
660  return InfFE<2,JACOBI_20_00,CARTESIAN>::shape(fe_t, t, i, p);
661 
662  case JACOBI_30_00:
663  return InfFE<2,JACOBI_30_00,CARTESIAN>::shape(fe_t, t, i, p);
664 
665  case LEGENDRE:
666  return InfFE<2,LEGENDRE,CARTESIAN>::shape(fe_t, t, i, p);
667 
668  case LAGRANGE:
669  return InfFE<2,LAGRANGE,CARTESIAN>::shape(fe_t, t, i, p);
670 
671  default:
672  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
673  }
674  }
675 
676 
677  // 3D
678  case 3:
679  {
680  switch (fe_t.radial_family)
681  {
682  case INFINITE_MAP:
683  return InfFE<3,INFINITE_MAP,CARTESIAN>::shape(fe_t, t, i, p);
684 
685  case JACOBI_20_00:
686  return InfFE<3,JACOBI_20_00,CARTESIAN>::shape(fe_t, t, i, p);
687 
688  case JACOBI_30_00:
689  return InfFE<3,JACOBI_30_00,CARTESIAN>::shape(fe_t, t, i, p);
690 
691  case LEGENDRE:
692  return InfFE<3,LEGENDRE,CARTESIAN>::shape(fe_t, t, i, p);
693 
694  case LAGRANGE:
695  return InfFE<3,LAGRANGE,CARTESIAN>::shape(fe_t, t, i, p);
696 
697  default:
698  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
699  }
700  }
701 
702  default:
703  libmesh_error_msg("Invalid dim = " << dim);
704  }
705 }
static Real shape(const FEType &fet, const ElemType t, const unsigned int i, const Point &p)
Real libMesh::FEInterface::ifem_shape ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
const unsigned int  i,
const Point p 
)
staticprivate

Definition at line 710 of file fe_interface_inf_fe.C.

References libMesh::INFINITE_MAP, libMesh::JACOBI_20_00, libMesh::JACOBI_30_00, libMesh::LAGRANGE, libMesh::LEGENDRE, libMesh::FEType::radial_family, and libMesh::InfFE< Dim, T_radial, T_map >::shape().

715 {
716  switch (dim)
717  {
718  // 1D
719  case 1:
720  {
721  switch (fe_t.radial_family)
722  {
723  /*
724  * For no derivatives (and local coordinates, as
725  * given in \p p) the infinite element shapes
726  * are independent of mapping type
727  */
728  case INFINITE_MAP:
729  return InfFE<1,INFINITE_MAP,CARTESIAN>::shape(fe_t, elem, i, p);
730 
731  case JACOBI_20_00:
732  return InfFE<1,JACOBI_20_00,CARTESIAN>::shape(fe_t, elem, i, p);
733 
734  case JACOBI_30_00:
735  return InfFE<1,JACOBI_30_00,CARTESIAN>::shape(fe_t, elem, i, p);
736 
737  case LEGENDRE:
738  return InfFE<1,LEGENDRE,CARTESIAN>::shape(fe_t, elem, i, p);
739 
740  case LAGRANGE:
741  return InfFE<1,LAGRANGE,CARTESIAN>::shape(fe_t, elem, i, p);
742 
743  default:
744  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
745  }
746  }
747 
748 
749  // 2D
750  case 2:
751  {
752  switch (fe_t.radial_family)
753  {
754  case INFINITE_MAP:
755  return InfFE<2,INFINITE_MAP,CARTESIAN>::shape(fe_t, elem, i, p);
756 
757  case JACOBI_20_00:
758  return InfFE<2,JACOBI_20_00,CARTESIAN>::shape(fe_t, elem, i, p);
759 
760  case JACOBI_30_00:
761  return InfFE<2,JACOBI_30_00,CARTESIAN>::shape(fe_t, elem, i, p);
762 
763  case LEGENDRE:
764  return InfFE<2,LEGENDRE,CARTESIAN>::shape(fe_t, elem, i, p);
765 
766  case LAGRANGE:
767  return InfFE<2,LAGRANGE,CARTESIAN>::shape(fe_t, elem, i, p);
768 
769  default:
770  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
771  }
772  }
773 
774 
775  // 3D
776  case 3:
777  {
778  switch (fe_t.radial_family)
779  {
780  case INFINITE_MAP:
781  return InfFE<3,INFINITE_MAP,CARTESIAN>::shape(fe_t, elem, i, p);
782 
783  case JACOBI_20_00:
784  return InfFE<3,JACOBI_20_00,CARTESIAN>::shape(fe_t, elem, i, p);
785 
786  case JACOBI_30_00:
787  return InfFE<3,JACOBI_30_00,CARTESIAN>::shape(fe_t, elem, i, p);
788 
789  case LEGENDRE:
790  return InfFE<3,LEGENDRE,CARTESIAN>::shape(fe_t, elem, i, p);
791 
792  case LAGRANGE:
793  return InfFE<3,LAGRANGE,CARTESIAN>::shape(fe_t, elem, i, p);
794 
795  default:
796  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
797  }
798  }
799 
800  default:
801  libmesh_error_msg("Invalid dim = " << dim);
802  }
803 }
static Real shape(const FEType &fet, const ElemType t, const unsigned int i, const Point &p)
Point libMesh::FEInterface::inverse_map ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
const Point p,
const Real  tolerance = TOLERANCE,
const bool  secure = true 
)
static
Returns
The location (on the reference element) of the point p located in physical space. This function requires inverting the (probably nonlinear) transformation map, so it is not trivial. The optional parameter tolerance defines how close is "good enough." The map inversion iteration computes the sequence $ \{ p_n \} $, and the iteration is terminated when $ \|p - p_n\| < \mbox{\texttt{tolerance}} $

Definition at line 547 of file fe_interface.C.

References ifem_inverse_map(), is_InfFE_elem(), and libMesh::Elem::type().

Referenced by libMesh::HPCoarsenTest::add_projection(), libMesh::FEMContext::build_new_fe(), libMesh::FEGenericBase< OutputType >::coarsened_dof_values(), libMesh::FEAbstract::compute_node_constraints(), libMesh::FEGenericBase< OutputType >::compute_periodic_constraints(), libMesh::FEAbstract::compute_periodic_node_constraints(), libMesh::FEGenericBase< OutputType >::compute_proj_constraints(), libMesh::InfQuad4::contains_point(), libMesh::InfPrism::contains_point(), libMesh::InfHex::contains_point(), libMesh::MeshFunction::discontinuous_gradient(), libMesh::MeshFunction::discontinuous_value(), libMesh::MeshFunction::gradient(), libMesh::MeshFunction::hessian(), inverse_map(), libMesh::DGFEMContext::neighbor_side_fe_reinit(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::operator()(), libMesh::MeshFunction::operator()(), libMesh::System::point_gradient(), libMesh::System::point_hessian(), libMesh::Elem::point_test(), libMesh::System::point_value(), libMesh::JumpErrorEstimator::reinit_sides(), and libMesh::HPCoarsenTest::select_refinement().

553 {
554 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
555 
556  if (is_InfFE_elem(elem->type()))
557  return ifem_inverse_map(dim, fe_t, elem, p,tolerance, secure);
558 
559 #endif
560 
561  fe_with_vec_switch(inverse_map(elem, p, tolerance, secure));
562 }
static Point ifem_inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_interface.C:547
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
void libMesh::FEInterface::inverse_map ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
const std::vector< Point > &  physical_points,
std::vector< Point > &  reference_points,
const Real  tolerance = TOLERANCE,
const bool  secure = true 
)
static
Returns
The location (on the reference element) of the points physical_points located in physical space. This function requires inverting the (probably nonlinear) transformation map, so it is not trivial. The location of each point on the reference element is returned in the vector reference_points. The optional parameter tolerance defines how close is "good enough." The map inversion iteration computes the sequence $ \{ p_n \} $, and the iteration is terminated when $ \|p - p_n\| < \mbox{\texttt{tolerance}} $

Definition at line 567 of file fe_interface.C.

References libMesh::err, ifem_inverse_map(), inverse_map(), is_InfFE_elem(), libMesh::TypeVector< T >::size(), and libMesh::Elem::type().

574 {
575  const std::size_t n_pts = physical_points.size();
576 
577  // Resize the vector
578  reference_points.resize(n_pts);
579 
580  if (n_pts == 0)
581  {
582  libMesh::err << "WARNING: empty vector physical_points!"
583  << std::endl;
584  libmesh_here();
585  return;
586  }
587 
588 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
589 
590  if (is_InfFE_elem(elem->type()))
591  {
592  ifem_inverse_map(dim, fe_t, elem, physical_points, reference_points, tolerance, secure);
593  return;
594  // libmesh_not_implemented();
595  }
596 
597 #endif
598 
599  void_fe_with_vec_switch(inverse_map(elem, physical_points, reference_points, tolerance, secure));
600 }
static Point ifem_inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_interface.C:547
OStreamProxy err(std::cerr)
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
bool libMesh::FEInterface::is_InfFE_elem ( const ElemType  et)
inlinestaticprivate
Returns
true if et is an element to be processed by class InfFE, false otherwise or when !LIBMESH_ENABLE_INFINITE_ELEMENTS.

Definition at line 468 of file fe_interface.h.

References libMesh::INFEDGE2, libMesh::INFHEX16, libMesh::INFHEX18, libMesh::INFHEX8, libMesh::INFPRISM12, libMesh::INFPRISM6, libMesh::INFQUAD4, and libMesh::INFQUAD6.

Referenced by compute_data(), inverse_map(), map(), n_dofs(), n_dofs_at_node(), n_dofs_per_elem(), n_shape_functions(), nodal_soln(), and shape().

469 {
470  return false;
471 }
Point libMesh::FEInterface::map ( unsigned int  dim,
const FEType fe_t,
const Elem elem,
const Point p 
)
static
Returns
The point in physical space corresponding to the reference point p which is passed in.

Definition at line 531 of file fe_interface.C.

References ifem_map(), is_InfFE_elem(), and libMesh::Elem::type().

Referenced by libMesh::Elem::point_test().

535 {
536 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
537  if (is_InfFE_elem(elem->type()))
538  return ifem_map(dim, fe_t, elem, p);
539 #endif
540  fe_with_vec_switch(map(elem, p));
541 }
static Point map(unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p)
Definition: fe_interface.C:531
static Point ifem_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p)
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
unsigned int libMesh::FEInterface::max_order ( const FEType fe_t,
const ElemType el_t 
)
static
Returns
The maximum polynomial degree that the given finite element family can support on the given geometric element.

Definition at line 996 of file fe_interface.C.

References libMesh::BERNSTEIN, libMesh::CLOUGH, libMesh::EDGE2, libMesh::EDGE3, libMesh::EDGE4, libMesh::FEType::family, libMesh::HERMITE, libMesh::HEX20, libMesh::HEX27, libMesh::HEX8, libMesh::HIERARCHIC, libMesh::L2_HIERARCHIC, libMesh::L2_LAGRANGE, libMesh::LAGRANGE, libMesh::LAGRANGE_VEC, libMesh::MONOMIAL, libMesh::NEDELEC_ONE, libMesh::PRISM15, libMesh::PRISM18, libMesh::PRISM6, libMesh::PYRAMID13, libMesh::PYRAMID14, libMesh::PYRAMID5, libMesh::QUAD4, libMesh::QUAD8, libMesh::QUAD9, libMesh::QUADSHELL4, libMesh::QUADSHELL8, libMesh::SUBDIVISION, libMesh::SZABAB, libMesh::TET10, libMesh::TET4, libMesh::TRI3, libMesh::TRI3SUBDIVISION, libMesh::TRI6, libMesh::TRISHELL3, and libMesh::XYZ.

Referenced by libMesh::DofMap::reinit().

998 {
999  // Yeah, I know, infinity is much larger than 11, but our
1000  // solvers don't seem to like high degree polynomials, and our
1001  // quadrature rules and number_lookups tables
1002  // need to go up higher.
1003  const unsigned int unlimited = 11;
1004 
1005  // If we used 0 as a default, then elements missing from this
1006  // table (e.g. infinite elements) would be considered broken.
1007  const unsigned int unknown = unlimited;
1008 
1009  switch (fe_t.family)
1010  {
1011  case LAGRANGE:
1012  case L2_LAGRANGE: // TODO: L2_LAGRANGE can have higher "max_order" than LAGRANGE
1013  case LAGRANGE_VEC:
1014  switch (el_t)
1015  {
1016  case EDGE2:
1017  case EDGE3:
1018  case EDGE4:
1019  return 3;
1020  case TRI3:
1021  case TRISHELL3:
1022  return 1;
1023  case TRI6:
1024  return 2;
1025  case QUAD4:
1026  case QUADSHELL4:
1027  return 1;
1028  case QUAD8:
1029  case QUADSHELL8:
1030  case QUAD9:
1031  return 2;
1032  case TET4:
1033  return 1;
1034  case TET10:
1035  return 2;
1036  case HEX8:
1037  return 1;
1038  case HEX20:
1039  case HEX27:
1040  return 2;
1041  case PRISM6:
1042  return 1;
1043  case PRISM15:
1044  case PRISM18:
1045  return 2;
1046  case PYRAMID5:
1047  return 1;
1048  case PYRAMID13:
1049  case PYRAMID14:
1050  return 2;
1051  default:
1052  return unknown;
1053  }
1054  break;
1055  case MONOMIAL:
1056  switch (el_t)
1057  {
1058  case EDGE2:
1059  case EDGE3:
1060  case EDGE4:
1061  case TRI3:
1062  case TRISHELL3:
1063  case TRI6:
1064  case QUAD4:
1065  case QUADSHELL4:
1066  case QUAD8:
1067  case QUADSHELL8:
1068  case QUAD9:
1069  case TET4:
1070  case TET10:
1071  case HEX8:
1072  case HEX20:
1073  case HEX27:
1074  case PRISM6:
1075  case PRISM15:
1076  case PRISM18:
1077  case PYRAMID5:
1078  case PYRAMID13:
1079  case PYRAMID14:
1080  return unlimited;
1081  default:
1082  return unknown;
1083  }
1084  break;
1085 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
1086  case BERNSTEIN:
1087  switch (el_t)
1088  {
1089  case EDGE2:
1090  case EDGE3:
1091  case EDGE4:
1092  return unlimited;
1093  case TRI3:
1094  case TRISHELL3:
1095  return 1;
1096  case TRI6:
1097  return 6;
1098  case QUAD4:
1099  case QUADSHELL4:
1100  return 1;
1101  case QUAD8:
1102  case QUADSHELL8:
1103  case QUAD9:
1104  return unlimited;
1105  case TET4:
1106  return 1;
1107  case TET10:
1108  return 2;
1109  case HEX8:
1110  return 1;
1111  case HEX20:
1112  return 2;
1113  case HEX27:
1114  return 4;
1115  case PRISM6:
1116  case PRISM15:
1117  case PRISM18:
1118  case PYRAMID5:
1119  case PYRAMID13:
1120  case PYRAMID14:
1121  return 0;
1122  default:
1123  return unknown;
1124  }
1125  break;
1126  case SZABAB:
1127  switch (el_t)
1128  {
1129  case EDGE2:
1130  return 1;
1131  case EDGE3:
1132  case EDGE4:
1133  return 7;
1134  case TRI3:
1135  case TRISHELL3:
1136  return 1;
1137  case TRI6:
1138  return 7;
1139  case QUAD4:
1140  case QUADSHELL4:
1141  return 1;
1142  case QUAD8:
1143  case QUADSHELL8:
1144  case QUAD9:
1145  return 7;
1146  case TET4:
1147  case TET10:
1148  case HEX8:
1149  case HEX20:
1150  case HEX27:
1151  case PRISM6:
1152  case PRISM15:
1153  case PRISM18:
1154  case PYRAMID5:
1155  case PYRAMID13:
1156  case PYRAMID14:
1157  return 0;
1158  default:
1159  return unknown;
1160  }
1161  break;
1162 #endif
1163  case XYZ:
1164  switch (el_t)
1165  {
1166  case EDGE2:
1167  case EDGE3:
1168  case EDGE4:
1169  case TRI3:
1170  case TRISHELL3:
1171  case TRI6:
1172  case QUAD4:
1173  case QUADSHELL4:
1174  case QUAD8:
1175  case QUADSHELL8:
1176  case QUAD9:
1177  case TET4:
1178  case TET10:
1179  case HEX8:
1180  case HEX20:
1181  case HEX27:
1182  case PRISM6:
1183  case PRISM15:
1184  case PRISM18:
1185  case PYRAMID5:
1186  case PYRAMID13:
1187  case PYRAMID14:
1188  return unlimited;
1189  default:
1190  return unknown;
1191  }
1192  break;
1193  case CLOUGH:
1194  switch (el_t)
1195  {
1196  case EDGE2:
1197  case EDGE3:
1198  return 3;
1199  case EDGE4:
1200  case TRI3:
1201  case TRISHELL3:
1202  return 0;
1203  case TRI6:
1204  return 3;
1205  case QUAD4:
1206  case QUADSHELL4:
1207  case QUAD8:
1208  case QUADSHELL8:
1209  case QUAD9:
1210  case TET4:
1211  case TET10:
1212  case HEX8:
1213  case HEX20:
1214  case HEX27:
1215  case PRISM6:
1216  case PRISM15:
1217  case PRISM18:
1218  case PYRAMID5:
1219  case PYRAMID13:
1220  case PYRAMID14:
1221  return 0;
1222  default:
1223  return unknown;
1224  }
1225  break;
1226  case HERMITE:
1227  switch (el_t)
1228  {
1229  case EDGE2:
1230  case EDGE3:
1231  return unlimited;
1232  case EDGE4:
1233  case TRI3:
1234  case TRISHELL3:
1235  case TRI6:
1236  return 0;
1237  case QUAD4:
1238  case QUADSHELL4:
1239  return 3;
1240  case QUAD8:
1241  case QUADSHELL8:
1242  case QUAD9:
1243  return unlimited;
1244  case TET4:
1245  case TET10:
1246  return 0;
1247  case HEX8:
1248  return 3;
1249  case HEX20:
1250  case HEX27:
1251  return unlimited;
1252  case PRISM6:
1253  case PRISM15:
1254  case PRISM18:
1255  case PYRAMID5:
1256  case PYRAMID13:
1257  case PYRAMID14:
1258  return 0;
1259  default:
1260  return unknown;
1261  }
1262  break;
1263  case HIERARCHIC:
1264  switch (el_t)
1265  {
1266  case EDGE2:
1267  case EDGE3:
1268  case EDGE4:
1269  return unlimited;
1270  case TRI3:
1271  case TRISHELL3:
1272  return 1;
1273  case TRI6:
1274  return unlimited;
1275  case QUAD4:
1276  case QUADSHELL4:
1277  return 1;
1278  case QUAD8:
1279  case QUADSHELL8:
1280  case QUAD9:
1281  return unlimited;
1282  case TET4:
1283  case TET10:
1284  return 0;
1285  case HEX8:
1286  case HEX20:
1287  return 1;
1288  case HEX27:
1289  return unlimited;
1290  case PRISM6:
1291  case PRISM15:
1292  case PRISM18:
1293  case PYRAMID5:
1294  case PYRAMID13:
1295  case PYRAMID14:
1296  return 0;
1297  default:
1298  return unknown;
1299  }
1300  break;
1301  case L2_HIERARCHIC:
1302  switch (el_t)
1303  {
1304  case EDGE2:
1305  case EDGE3:
1306  case EDGE4:
1307  return unlimited;
1308  case TRI3:
1309  case TRISHELL3:
1310  return 1;
1311  case TRI6:
1312  return unlimited;
1313  case QUAD4:
1314  case QUADSHELL4:
1315  return 1;
1316  case QUAD8:
1317  case QUADSHELL8:
1318  case QUAD9:
1319  return unlimited;
1320  case TET4:
1321  case TET10:
1322  return 0;
1323  case HEX8:
1324  case HEX20:
1325  return 1;
1326  case HEX27:
1327  return unlimited;
1328  case PRISM6:
1329  case PRISM15:
1330  case PRISM18:
1331  case PYRAMID5:
1332  case PYRAMID13:
1333  case PYRAMID14:
1334  return 0;
1335  default:
1336  return unknown;
1337  }
1338  break;
1339  case SUBDIVISION:
1340  switch (el_t)
1341  {
1342  case TRI3SUBDIVISION:
1343  return unlimited;
1344  default:
1345  return unknown;
1346  }
1347  break;
1348  case NEDELEC_ONE:
1349  switch (el_t)
1350  {
1351  case TRI6:
1352  case QUAD8:
1353  case QUAD9:
1354  case HEX20:
1355  case HEX27:
1356  return 1;
1357  default:
1358  return 0;
1359  }
1360  break;
1361  default:
1362  return 0;
1363  break;
1364  }
1365 }
unsigned int libMesh::FEInterface::n_dofs ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t 
)
static
Returns
The number of shape functions associated with this finite element. Automatically decides which finite element class to use.

On a p-refined element, fe_t.order should be the total order of the element.

Definition at line 411 of file fe_interface.C.

References ifem_n_dofs(), is_InfFE_elem(), and libMesh::FEType::order.

Referenced by libMesh::DofMap::_dof_indices(), libMesh::FEGenericBase< OutputType >::coarsened_dof_values(), compute_data(), libMesh::FEAbstract::compute_node_constraints(), libMesh::FEAbstract::compute_periodic_node_constraints(), libMesh::InfFE< Dim, T_radial, T_map >::n_dofs(), libMesh::HPCoarsenTest::select_refinement(), and ~FEInterface().

414 {
415 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
416 
417  if (is_InfFE_elem(t))
418  return ifem_n_dofs(dim, fe_t, t);
419 
420 #endif
421 
422  const Order o = fe_t.order;
423 
424  fe_with_vec_switch(n_dofs(t, o));
425 }
static unsigned int n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:411
static unsigned int ifem_n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
unsigned int libMesh::FEInterface::n_dofs_at_node ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t,
const unsigned int  n 
)
static
Returns
The number of dofs at node n for a finite element of type fe_t. Automatically decides which finite element class to use.

On a p-refined element, fe_t.order should be the total order of the element.

Definition at line 430 of file fe_interface.C.

References ifem_n_dofs_at_node(), is_InfFE_elem(), and libMesh::FEType::order.

Referenced by libMesh::DofMap::_dof_indices(), libMesh::FEGenericBase< OutputType >::coarsened_dof_values(), libMesh::InfFE< Dim, T_radial, T_map >::compute_shape_indices(), libMesh::InfFE< Dim, T_radial, T_map >::n_dofs_at_node(), n_dofs_at_node_function(), libMesh::DofMap::old_dof_indices(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::operator()(), libMesh::BoundaryProjectSolution::operator()(), libMesh::DofMap::reinit(), and ~FEInterface().

434 {
435 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
436 
437  if (is_InfFE_elem(t))
438  return ifem_n_dofs_at_node(dim, fe_t, t, n);
439 
440 #endif
441 
442  const Order o = fe_t.order;
443 
444  fe_with_vec_switch(n_dofs_at_node(t, o, n));
445 }
static unsigned int n_dofs_at_node(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
Definition: fe_interface.C:430
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
static unsigned int ifem_n_dofs_at_node(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
FEInterface::n_dofs_at_node_ptr libMesh::FEInterface::n_dofs_at_node_function ( const unsigned int  dim,
const FEType fe_t 
)
static
Returns
A function which evaluates n_dofs_at_node for the requested FE type and dimension.

Definition at line 450 of file fe_interface.C.

References n_dofs_at_node().

Referenced by libMesh::DofMap::_dof_indices(), and libMesh::DofMap::old_dof_indices().

452 {
453  fe_with_vec_switch(n_dofs_at_node);
454 }
static unsigned int n_dofs_at_node(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
Definition: fe_interface.C:430
unsigned int libMesh::FEInterface::n_dofs_per_elem ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t 
)
static
Returns
The number of dofs interior to the element, not associated with any interior nodes. Automatically decides which finite element class to use.

On a p-refined element, fe_t.order should be the total order of the element.

Definition at line 461 of file fe_interface.C.

References ifem_n_dofs_per_elem(), is_InfFE_elem(), and libMesh::FEType::order.

Referenced by libMesh::DofMap::_dof_indices(), libMesh::InfFE< Dim, T_radial, T_map >::compute_shape_indices(), libMesh::InfFE< Dim, T_radial, T_map >::n_dofs_per_elem(), libMesh::DofMap::old_dof_indices(), and libMesh::DofMap::reinit().

464 {
465 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
466 
467  if (is_InfFE_elem(t))
468  return ifem_n_dofs_per_elem(dim, fe_t, t);
469 
470 #endif
471 
472  const Order o = fe_t.order;
473 
474  fe_with_vec_switch(n_dofs_per_elem(t, o));
475 }
static unsigned int n_dofs_per_elem(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:461
static unsigned int ifem_n_dofs_per_elem(const unsigned int dim, const FEType &fe_t, const ElemType t)
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
unsigned int libMesh::FEInterface::n_shape_functions ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t 
)
static
Returns
The number of shape functions associated with this finite element of type fe_t. Automatically decides which finite element class to use.

On a p-refined element, fe_t.order should be the total order of the element.

Definition at line 385 of file fe_interface.C.

References ifem_n_shape_functions(), is_InfFE_elem(), and libMesh::FEType::order.

Referenced by ~FEInterface().

388 {
389 
390 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
391  /*
392  * Since the FEType, stored in DofMap/(some System child), has to
393  * be the _same_ for InfFE and FE, we have to catch calls
394  * to infinite elements through the element type.
395  */
396 
397  if (is_InfFE_elem(t))
398  return ifem_n_shape_functions(dim, fe_t, t);
399 
400 #endif
401 
402  const Order o = fe_t.order;
403 
404  fe_with_vec_switch(n_shape_functions(t, o));
405 }
static unsigned int ifem_n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
static unsigned int n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:385
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
unsigned int libMesh::FEInterface::n_vec_dim ( const MeshBase mesh,
const FEType fe_type 
)
static
Returns
The number of components of a vector-valued element. Scalar-valued elements return 1.

Definition at line 1411 of file fe_interface.C.

References libMesh::FEType::family, libMesh::LAGRANGE_VEC, libMesh::MeshBase::mesh_dimension(), and libMesh::NEDELEC_ONE.

Referenced by libMesh::ExactSolution::_compute_error(), libMesh::EquationSystems::build_parallel_solution_vector(), and libMesh::EquationSystems::build_variable_names().

1413 {
1414  switch (fe_type.family)
1415  {
1416  //FIXME: We currently assume that the number of vector components is tied
1417  // to the mesh dimension. This will break for mixed-dimension meshes.
1418  case LAGRANGE_VEC:
1419  case NEDELEC_ONE:
1420  return mesh.mesh_dimension();
1421  default:
1422  return 1;
1423  }
1424 }
MeshBase & mesh
void libMesh::FEInterface::nodal_soln ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
const std::vector< Number > &  elem_soln,
std::vector< Number > &  nodal_soln 
)
static

Build the nodal soln from the element soln. This is the solution that will be plotted. Automatically passes the request to the appropriate finite element class member. To indicate that results from this specific implementation of nodal_soln should not be used, the vector nodal_soln is returned empty.

On a p-refined element, fe_t.order should be the base order of the element.

Definition at line 507 of file fe_interface.C.

References ifem_nodal_soln(), is_InfFE_elem(), libMesh::FEType::order, and libMesh::Elem::type().

Referenced by libMesh::EquationSystems::build_discontinuous_solution_vector(), libMesh::EquationSystems::build_parallel_solution_vector(), libMesh::EnsightIO::write_scalar_ascii(), and libMesh::EnsightIO::write_vector_ascii().

512 {
513 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
514 
515  if (is_InfFE_elem(elem->type()))
516  {
517  ifem_nodal_soln(dim, fe_t, elem, elem_soln, nodal_soln);
518  return;
519  }
520 
521 #endif
522 
523  const Order order = fe_t.order;
524 
525  void_fe_with_vec_switch(nodal_soln(elem, order, elem_soln, nodal_soln));
526 }
static void ifem_nodal_soln(const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
static void nodal_soln(const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
Definition: fe_interface.C:507
bool libMesh::FEInterface::on_reference_element ( const Point p,
const ElemType  t,
const Real  eps = TOLERANCE 
)
static
Returns
true if the point p is located on the reference element for element type t, false otherwise.

Since we are doing floating point comparisons, the parameter eps can be specified to indicate a tolerance. For example, $ \xi \le 1 $ becomes $ \xi \le 1 + \epsilon $.

Definition at line 604 of file fe_interface.C.

References libMesh::FEAbstract::on_reference_element().

Referenced by libMesh::InfQuad4::contains_point(), libMesh::InfPrism::contains_point(), libMesh::InfHex::contains_point(), and libMesh::Elem::point_test().

607 {
608  return FEBase::on_reference_element(p,t,eps);
609 }
static bool on_reference_element(const Point &p, const ElemType t, const Real eps=TOLERANCE)
Definition: fe_abstract.C:554
Real libMesh::FEInterface::shape ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t,
const unsigned int  i,
const Point p 
)
static
Returns
The value of the $ i^{th} $ shape function at point p. This method allows you to specify the dimension, element type, and order directly. Automatically passes the request to the appropriate finite element class member.
Note
On a p-refined element, fe_t.order should be the total order of the element.

Definition at line 614 of file fe_interface.C.

References ifem_shape(), is_InfFE_elem(), and libMesh::FEType::order.

Referenced by libMesh::InfFE< Dim, T_radial, T_map >::compute_data(), compute_data(), libMesh::FEAbstract::compute_node_constraints(), libMesh::FEAbstract::compute_periodic_node_constraints(), shape(), and libMesh::InfFE< Dim, T_radial, T_map >::shape().

619 {
620 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
621 
622  if (is_InfFE_elem(t))
623  return ifem_shape(dim, fe_t, t, i, p);
624 
625 #endif
626 
627  const Order o = fe_t.order;
628 
629  fe_switch(shape(t,o,i,p));
630 }
static Real ifem_shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition: fe_interface.C:614
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
Real libMesh::FEInterface::shape ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
const unsigned int  i,
const Point p 
)
static
Returns
The value of the $ i^{th} $ shape function at point p. This method allows you to specify the dimension, element type, and order directly. Automatically passes the request to the appropriate finite element class member.
Note
On a p-refined element, fe_t.order should be the base order of the element.

Definition at line 632 of file fe_interface.C.

References ifem_shape(), is_InfFE_elem(), libMesh::FEType::order, shape(), and libMesh::Elem::type().

637 {
638 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
639 
640  if (elem && is_InfFE_elem(elem->type()))
641  return ifem_shape(dim, fe_t, elem, i, p);
642 
643 #endif
644 
645  const Order o = fe_t.order;
646 
647  fe_switch(shape(elem,o,i,p));
648 }
static Real ifem_shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition: fe_interface.C:614
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
template<typename OutputType >
static void libMesh::FEInterface::shape ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t,
const unsigned int  i,
const Point p,
OutputType &  phi 
)
static
Returns
The value of the $ i^{th} $ shape function at point p. This method allows you to specify the dimension, element type, and order directly. Automatically passes the request to the appropriate scalar finite element class member.
Note
On a p-refined element, fe_t.order should be the total order of the element.
template<typename OutputType >
static void libMesh::FEInterface::shape ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
const unsigned int  i,
const Point p,
OutputType &  phi 
)
static
Returns
The value of the $ i^{th} $ shape function at point p. This method allows you to specify the dimension, element type, and order directly. Automatically passes the request to the appropriate scalar finite element class member.
Note
On a p-refined element, fe_t.order should be the total order of the element.
template<>
void libMesh::FEInterface::shape ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t,
const unsigned int  i,
const Point p,
Real phi 
)

Definition at line 651 of file fe_interface.C.

References ifem_shape(), is_InfFE_elem(), and shape().

657 {
658 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
659 
660  if (is_InfFE_elem(t))
661  {
662  phi = ifem_shape(dim, fe_t, t, i, p);
663  return;
664  }
665 
666 #endif
667 
668  const Order o = fe_t.order;
669 
670  switch(dim)
671  {
672  case 0:
673  fe_scalar_vec_error_switch(0, shape(t,o,i,p), phi = , ; break;);
674  break;
675  case 1:
676  fe_scalar_vec_error_switch(1, shape(t,o,i,p), phi = , ; break;);
677  break;
678  case 2:
679  fe_scalar_vec_error_switch(2, shape(t,o,i,p), phi = , ; break;);
680  break;
681  case 3:
682  fe_scalar_vec_error_switch(3, shape(t,o,i,p), phi = , ; break;);
683  break;
684  default:
685  libmesh_error_msg("Invalid dimension = " << dim);
686  }
687 
688  return;
689 }
static Real ifem_shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition: fe_interface.C:614
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
template<>
void libMesh::FEInterface::shape ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
const unsigned int  i,
const Point p,
Real phi 
)

Definition at line 692 of file fe_interface.C.

References ifem_shape(), is_InfFE_elem(), and shape().

698 {
699 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
700 
701  if (elem && is_InfFE_elem(elem->type()))
702  {
703  phi = ifem_shape(dim, fe_t, elem, i, p);
704  return;
705  }
706 #endif
707 
708  const Order o = fe_t.order;
709 
710  switch(dim)
711  {
712  case 0:
713  fe_scalar_vec_error_switch(0, shape(elem,o,i,p), phi = , ; break;);
714  break;
715  case 1:
716  fe_scalar_vec_error_switch(1, shape(elem,o,i,p), phi = , ; break;);
717  break;
718  case 2:
719  fe_scalar_vec_error_switch(2, shape(elem,o,i,p), phi = , ; break;);
720  break;
721  case 3:
722  fe_scalar_vec_error_switch(3, shape(elem,o,i,p), phi = , ; break;);
723  break;
724  default:
725  libmesh_error_msg("Invalid dimension = " << dim);
726  }
727 
728  return;
729 }
static Real ifem_shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition: fe_interface.C:614
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
template<>
void libMesh::FEInterface::shape ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t,
const unsigned int  i,
const Point p,
RealGradient phi 
)

Definition at line 732 of file fe_interface.C.

References shape().

738 {
739  // This even does not handle infinite elements at all!?
740  const Order o = fe_t.order;
741 
742  switch(dim)
743  {
744  case 0:
745  fe_vector_scalar_error_switch(0, shape(t,o,i,p), phi = , ; break;);
746  break;
747  case 1:
748  fe_vector_scalar_error_switch(1, shape(t,o,i,p), phi = , ; break;);
749  break;
750  case 2:
751  fe_vector_scalar_error_switch(2, shape(t,o,i,p), phi = , ; break;);
752  break;
753  case 3:
754  fe_vector_scalar_error_switch(3, shape(t,o,i,p), phi = , ; break;);
755  break;
756  default:
757  libmesh_error_msg("Invalid dimension = " << dim);
758  }
759 
760  return;
761 }
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition: fe_interface.C:614
template<>
void libMesh::FEInterface::shape ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
const unsigned int  i,
const Point p,
RealGradient phi 
)

Definition at line 764 of file fe_interface.C.

References shape().

770 {
771  const Order o = fe_t.order;
772 
773  switch(dim)
774  {
775  case 0:
776  fe_vector_scalar_error_switch(0, shape(elem,o,i,p), phi = , ; break;);
777  break;
778  case 1:
779  fe_vector_scalar_error_switch(1, shape(elem,o,i,p), phi = , ; break;);
780  break;
781  case 2:
782  fe_vector_scalar_error_switch(2, shape(elem,o,i,p), phi = , ; break;);
783  break;
784  case 3:
785  fe_vector_scalar_error_switch(3, shape(elem,o,i,p), phi = , ; break;);
786  break;
787  default:
788  libmesh_error_msg("Invalid dimension = " << dim);
789  }
790 
791  return;
792 }
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition: fe_interface.C:614

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