libMesh::FEInterface Class Reference

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

#include <fe_interface.h>

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 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 62 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 }

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 844 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::libmesh_assert(), libMesh::SZABAB, and libMesh::DofMap::variable_type().

Referenced by ~FEInterface().

848 {
849  libmesh_assert(elem);
850 
851  const FEType & fe_t = dof_map.variable_type(variable_number);
852 
853  switch (elem->dim())
854  {
855  case 0:
856  case 1:
857  {
858  // No constraints in 0D/1D.
859  return;
860  }
861 
862 
863  case 2:
864  {
865  switch (fe_t.family)
866  {
867  case CLOUGH:
869  dof_map,
870  variable_number,
871  elem); return;
872 
873  case HERMITE:
875  dof_map,
876  variable_number,
877  elem); return;
878 
879  case LAGRANGE:
881  dof_map,
882  variable_number,
883  elem); return;
884 
885  case HIERARCHIC:
887  dof_map,
888  variable_number,
889  elem); return;
890 
891  case L2_HIERARCHIC:
893  dof_map,
894  variable_number,
895  elem); return;
896 
897  case LAGRANGE_VEC:
899  dof_map,
900  variable_number,
901  elem); return;
902 
903 
904 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
905  case SZABAB:
907  dof_map,
908  variable_number,
909  elem); return;
910 
911  case BERNSTEIN:
913  dof_map,
914  variable_number,
915  elem); return;
916 
917 #endif
918  default:
919  return;
920  }
921  }
922 
923 
924  case 3:
925  {
926  switch (fe_t.family)
927  {
928  case HERMITE:
930  dof_map,
931  variable_number,
932  elem); return;
933 
934  case LAGRANGE:
936  dof_map,
937  variable_number,
938  elem); return;
939 
940  case HIERARCHIC:
942  dof_map,
943  variable_number,
944  elem); return;
945 
946  case L2_HIERARCHIC:
948  dof_map,
949  variable_number,
950  elem); return;
951 
952  case LAGRANGE_VEC:
954  dof_map,
955  variable_number,
956  elem); return;
957 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
958  case SZABAB:
960  dof_map,
961  variable_number,
962  elem); return;
963 
964  case BERNSTEIN:
966  dof_map,
967  variable_number,
968  elem); return;
969 
970 #endif
971  default:
972  return;
973  }
974  }
975 
976 
977  default:
978  libmesh_error_msg("Invalid dimension = " << elem->dim());
979  }
980 }
libmesh_assert(j)
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 returns the values also through data. See this as a generalization of shape(). Currently, with disabled infinite elements, returns a vector of all shape functions of elem evaluated ap p.

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

Definition at line 808 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 ~FEInterface().

812 {
813 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
814 
815  if ( elem && is_InfFE_elem(elem->type()) )
816  {
817  data.init();
818  ifem_compute_data(dim, fe_t, elem, data);
819  return;
820  }
821 
822 #endif
823 
824  FEType p_refined = fe_t;
825  p_refined.order = static_cast<Order>(p_refined.order + elem->p_level());
826 
827  const unsigned int n_dof = n_dofs (dim, p_refined, elem->type());
828  const Point & p = data.p;
829  data.shape.resize(n_dof);
830 
831  // set default values for all the output fields
832  data.init();
833 
834  for (unsigned int n=0; n<n_dof; n++)
835  data.shape[n] = shape(dim, p_refined, elem, n, p);
836 
837  return;
838 }
static unsigned int n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:414
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:628
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:453
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 988 of file fe_interface.C.

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

Referenced by ~FEInterface().

995 {
996  // No element-specific optimizations currently exist
998  dof_map,
999  boundaries,
1000  mesh,
1001  point_locator,
1002  variable_number,
1003  elem);
1004 }
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:1658
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 497 of file fe_interface.C.

References libMesh::FEType::order.

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

502 {
503  const Order o = fe_t.order;
504 
505  void_fe_with_vec_switch(dofs_on_edge(elem, o, e, di));
506 
507  libmesh_error_msg("We'll never get here!");
508 }
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:497
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 482 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()(), libMesh::BoundaryProjectSolution::operator()(), and ~FEInterface().

487 {
488  const Order o = fe_t.order;
489 
490  void_fe_with_vec_switch(dofs_on_side(elem, o, s, di));
491 
492  libmesh_error_msg("We'll never get here!");
493 }
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:482
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 1374 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(), libMesh::DofMap::reinit(), and ~FEInterface().

1375 {
1376  switch (fe_t.family)
1377  {
1378  case LAGRANGE:
1379  case L2_LAGRANGE:
1380  case MONOMIAL:
1381  case L2_HIERARCHIC:
1382  case XYZ:
1383  case SUBDIVISION:
1384  case LAGRANGE_VEC:
1385  case NEDELEC_ONE:
1386  return false;
1387  case CLOUGH:
1388  case HERMITE:
1389  case HIERARCHIC:
1390 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
1391  case BERNSTEIN:
1392  case SZABAB:
1393 #endif
1394  default:
1395  return true;
1396  }
1397 }
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 1404 of file fe_interface.C.

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

1405 {
1406  switch (fe_family)
1407  {
1408  case LAGRANGE_VEC:
1409  case NEDELEC_ONE:
1410  return TYPE_VECTOR;
1411  default:
1412  return TYPE_SCALAR;
1413  }
1414 }
void libMesh::FEInterface::ifem_compute_data ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
FEComputeData data 
)
staticprivate

Definition at line 830 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(), and ~FEInterface().

834 {
835  switch (dim)
836  {
837  // 1D
838  case 1:
839  {
840  switch (fe_t.radial_family)
841  {
842  /*
843  * For no derivatives (and local coordinates, as
844  * given in \p p) the infinite element shapes
845  * are independent of mapping type
846  */
847  case INFINITE_MAP:
849  break;
850 
851  case JACOBI_20_00:
853  break;
854 
855  case JACOBI_30_00:
857  break;
858 
859  case LEGENDRE:
861  break;
862 
863  case LAGRANGE:
865  break;
866 
867  default:
868  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
869  }
870 
871  break;
872  }
873 
874 
875  // 2D
876  case 2:
877  {
878  switch (fe_t.radial_family)
879  {
880  case INFINITE_MAP:
882  break;
883 
884  case JACOBI_20_00:
886  break;
887 
888  case JACOBI_30_00:
890  break;
891 
892  case LEGENDRE:
894  break;
895 
896  case LAGRANGE:
898  break;
899 
900  default:
901  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
902  }
903 
904  break;
905  }
906 
907 
908  // 3D
909  case 3:
910  {
911  switch (fe_t.radial_family)
912  {
913  case INFINITE_MAP:
915  break;
916 
917  case JACOBI_20_00:
919  break;
920 
921  case JACOBI_30_00:
923  break;
924 
925  case LEGENDRE:
927  break;
928 
929  case LAGRANGE:
931  break;
932 
933  default:
934  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
935  }
936 
937  break;
938  }
939 
940 
941  default:
942  libmesh_error_msg("Invalid dim = " << dim);
943  break;
944  }
945 }
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 459 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(), and ~FEInterface().

465 {
466  switch (dim)
467  {
468  // 1D
469  case 1:
470  {
471  switch (fe_t.inf_map)
472  {
473  case CARTESIAN:
474  return InfFE<1,JACOBI_20_00,CARTESIAN>::inverse_map(elem, p, tolerance, secure);
475 
476  case SPHERICAL:
477  case ELLIPSOIDAL:
478  libmesh_not_implemented_msg("ERROR: Spherical and Ellipsoidal IFEMs not (yet) implemented.");
479 
480  /*
481  case SPHERICAL:
482  return InfFE<1,JACOBI_20_00,SPHERICAL>::inverse_map(elem, p, tolerance);
483 
484  case ELLIPSOIDAL:
485  return InfFE<1,JACOBI_20_00,ELLIPSOIDAL>::inverse_map(elem, p, tolerance);
486  */
487 
488  default:
489  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
490  }
491  }
492 
493 
494  // 2D
495  case 2:
496  {
497  switch (fe_t.inf_map)
498  {
499  case CARTESIAN:
500  return InfFE<2,JACOBI_20_00,CARTESIAN>::inverse_map(elem, p, tolerance, secure);
501 
502  case SPHERICAL:
503  case ELLIPSOIDAL:
504  libmesh_not_implemented_msg("ERROR: Spherical and Ellipsoidal IFEMs not (yet) implemented.");
505 
506  /*
507  case SPHERICAL:
508  return InfFE<2,JACOBI_20_00,SPHERICAL>::inverse_map(elem, p, tolerance);
509 
510  case ELLIPSOIDAL:
511  return InfFE<2,JACOBI_20_00,ELLIPSOIDAL>::inverse_map(elem, p, tolerance);
512  */
513 
514  default:
515  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
516  }
517  }
518 
519 
520  // 3D
521  case 3:
522  {
523  switch (fe_t.inf_map)
524  {
525  case CARTESIAN:
526  return InfFE<3,JACOBI_20_00,CARTESIAN>::inverse_map(elem, p, tolerance, secure);
527 
528  case SPHERICAL:
529  case ELLIPSOIDAL:
530  libmesh_not_implemented_msg("ERROR: Spherical and Ellipsoidal IFEMs not (yet) implemented.");
531 
532  /*
533  case SPHERICAL:
534  return InfFE<3,JACOBI_20_00,SPHERICAL>::inverse_map(elem, p, tolerance);
535 
536  case ELLIPSOIDAL:
537  return InfFE<3,JACOBI_20_00,ELLIPSOIDAL>::inverse_map(elem, p, tolerance);
538  */
539 
540  default:
541  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
542  }
543  }
544 
545  default:
546  libmesh_error_msg("Invalid dim = " << dim);
547  }
548 
549  libmesh_error_msg("We'll never get here!");
550  Point pt;
551  return pt;
552 }
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 556 of file fe_interface_inf_fe.C.

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

563 {
564  switch (dim)
565  {
566  // 1D
567  case 1:
568  {
569  switch (fe_t.inf_map)
570  {
571  case CARTESIAN:
572  InfFE<1,JACOBI_20_00,CARTESIAN>::inverse_map(elem, physical_points, reference_points, tolerance, secure);
573  return;
574 
575  default:
576  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
577  }
578  }
579 
580 
581  // 2D
582  case 2:
583  {
584  switch (fe_t.inf_map)
585  {
586  case CARTESIAN:
587  InfFE<2,JACOBI_20_00,CARTESIAN>::inverse_map(elem, physical_points, reference_points, tolerance, secure);
588  return;
589 
590  default:
591  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
592  }
593  }
594 
595 
596  // 3D
597  case 3:
598  {
599  switch (fe_t.inf_map)
600  {
601  case CARTESIAN:
602  InfFE<3,JACOBI_20_00,CARTESIAN>::inverse_map(elem, physical_points, reference_points, tolerance, secure);
603  return;
604 
605  default:
606  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
607  }
608  }
609 
610  default:
611  libmesh_error_msg("Invalid dim = " << dim);
612  }
613 }
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 428 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(), and ~FEInterface().

432 {
433  switch (fe_t.inf_map)
434  {
435  case CARTESIAN:
436  {
437  switch (dim)
438  {
439  case 1:
440  return InfFE<1,JACOBI_20_00,CARTESIAN>::map(elem, p);
441  case 2:
442  return InfFE<2,JACOBI_20_00,CARTESIAN>::map(elem, p);
443  case 3:
444  return InfFE<3,JACOBI_20_00,CARTESIAN>::map(elem, p);
445  default:
446  libmesh_error_msg("Invalid dim = " << dim);
447  }
448  }
449  case SPHERICAL:
450  case ELLIPSOIDAL:
451  libmesh_not_implemented_msg("ERROR: Spherical and Ellipsoidal IFEMs not (yet) implemented.");
452  default:
453  libmesh_error_msg("Invalid map = " << fe_t.inf_map);
454  }
455 }
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 74 of file fe_interface_inf_fe.C.

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

Referenced by n_dofs(), and ~FEInterface().

77 {
78  switch (dim)
79  {
80  // 1D
81  case 1:
82  /*
83  * Since InfFE<Dim,T_radial,T_map>::n_dofs(...)
84  * is actually independent of T_radial and T_map, we can use
85  * just any T_radial and T_map
86  */
88 
89  // 2D
90  case 2:
92 
93  // 3D
94  case 3:
96 
97  default:
98  libmesh_error_msg("Unsupported dim = " << dim);
99  }
100 
101  libmesh_error_msg("We'll never get here!");
102  return 0;
103 }
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 108 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(), and ~FEInterface().

112 {
113  switch (dim)
114  {
115  // 1D
116  case 1:
117  /*
118  * Since InfFE<Dim,T_radial,T_map>::n_dofs_at_node(...)
119  * is actually independent of T_radial and T_map, we can use
120  * just any T_radial and T_map
121  */
123 
124  // 2D
125  case 2:
127 
128  // 3D
129  case 3:
131 
132  default:
133  libmesh_error_msg("Unsupported dim = " << dim);
134  }
135 
136  libmesh_error_msg("We'll never get here!");
137  return 0;
138 }
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 144 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(), and ~FEInterface().

147 {
148  switch (dim)
149  {
150  // 1D
151  case 1:
152  /*
153  * Since InfFE<Dim,T_radial,T_map>::n_dofs(...)
154  * is actually independent of T_radial and T_map, we can use
155  * just any T_radial and T_map
156  */
158 
159  // 2D
160  case 2:
162 
163  // 3D
164  case 3:
166 
167  default:
168  libmesh_error_msg("Unsupported dim = " << dim);
169  }
170 
171  libmesh_error_msg("We'll never get here!");
172  return 0;
173 }
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(), and ~FEInterface().

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 
66  libmesh_error_msg("We'll never get here!");
67  return 0;
68 }
virtual unsigned int n_shape_functions() const libmesh_override
Definition: inf_fe.h:453
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 178 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(), and ~FEInterface().

183 {
184  switch (dim)
185  {
186 
187  // 1D
188  case 1:
189  {
190  switch (fe_t.radial_family)
191  {
192  case INFINITE_MAP:
193  libmesh_error_msg("ERROR: INFINTE_MAP is not a valid shape family for radial approximation.");
194 
195  case JACOBI_20_00:
196  {
197  switch (fe_t.inf_map)
198  {
199  case CARTESIAN:
200  {
202  break;
203  }
204  default:
205  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
206  }
207  break;
208  }
209 
210  case JACOBI_30_00:
211  {
212  switch (fe_t.inf_map)
213  {
214  case CARTESIAN:
215  {
217  break;
218  }
219  default:
220  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
221  }
222  break;
223  }
224 
225  case LEGENDRE:
226  {
227  switch (fe_t.inf_map)
228  {
229  case CARTESIAN:
230  {
231  InfFE<1,LEGENDRE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
232  break;
233  }
234  default:
235  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
236  }
237  break;
238  }
239 
240  case LAGRANGE:
241  {
242  switch (fe_t.inf_map)
243  {
244  case CARTESIAN:
245  {
246  InfFE<1,LAGRANGE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
247  break;
248  }
249  default:
250  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
251  }
252  break;
253  }
254 
255  default:
256  libmesh_error_msg("ERROR: Bad FEType.radial_family= " << fe_t.radial_family);
257  }
258 
259  break;
260  }
261 
262 
263 
264 
265  // 2D
266  case 2:
267  {
268  switch (fe_t.radial_family)
269  {
270  case INFINITE_MAP:
271  libmesh_error_msg("ERROR: INFINTE_MAP is not a valid shape family for radial approximation.");
272 
273  case JACOBI_20_00:
274  {
275  switch (fe_t.inf_map)
276  {
277  case CARTESIAN:
278  {
280  break;
281  }
282  default:
283  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
284  }
285  break;
286  }
287 
288  case JACOBI_30_00:
289  {
290  switch (fe_t.inf_map)
291  {
292  case CARTESIAN:
293  {
295  break;
296  }
297  default:
298  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
299  }
300  break;
301  }
302 
303  case LEGENDRE:
304  {
305  switch (fe_t.inf_map)
306  {
307  case CARTESIAN:
308  {
309  InfFE<2,LEGENDRE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
310  break;
311  }
312  default:
313  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
314  }
315  break;
316  }
317 
318  case LAGRANGE:
319  {
320  switch (fe_t.inf_map)
321  {
322  case CARTESIAN:
323  {
324  InfFE<2,LAGRANGE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
325  break;
326  }
327  default:
328  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
329  }
330  break;
331  }
332 
333  default:
334  libmesh_error_msg("ERROR: Bad FEType.radial_family= " << fe_t.radial_family);
335  }
336 
337  break;
338  }
339 
340 
341 
342 
343  // 3D
344  case 3:
345  {
346  switch (fe_t.radial_family)
347  {
348  case INFINITE_MAP:
349  libmesh_error_msg("ERROR: INFINTE_MAP is not a valid shape family for radial approximation.");
350 
351  case JACOBI_20_00:
352  {
353  switch (fe_t.inf_map)
354  {
355  case CARTESIAN:
356  {
358  break;
359  }
360  default:
361  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
362  }
363  break;
364  }
365 
366  case JACOBI_30_00:
367  {
368  switch (fe_t.inf_map)
369  {
370  case CARTESIAN:
371  {
373  break;
374  }
375  default:
376  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
377  }
378  break;
379  }
380 
381  case LEGENDRE:
382  {
383  switch (fe_t.inf_map)
384  {
385  case CARTESIAN:
386  {
387  InfFE<3,LEGENDRE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
388  break;
389  }
390  default:
391  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
392  }
393  break;
394  }
395 
396  case LAGRANGE:
397  {
398  switch (fe_t.inf_map)
399  {
400  case CARTESIAN:
401  {
402  InfFE<3,LAGRANGE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
403  break;
404  }
405  default:
406  libmesh_error_msg("ERROR: Spherical & Ellipsoidal IFEMs not implemented.");
407  }
408  break;
409  }
410 
411 
412 
413  default:
414  libmesh_error_msg("ERROR: Bad FEType.radial_family= " << fe_t.radial_family);
415  }
416 
417  break;
418  }
419 
420  default:
421  libmesh_error_msg("Invalid dim = " << dim);
422  }
423 }
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:513
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 618 of file fe_interface_inf_fe.C.

References libMesh::FEAbstract::on_reference_element().

Referenced by ~FEInterface().

621 {
622  return FEBase::on_reference_element(p,t,eps);
623 }
static bool on_reference_element(const Point &p, const ElemType t, const Real eps=TOLERANCE)
Definition: fe_abstract.C:556
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 628 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(), and ~FEInterface().

633 {
634  switch (dim)
635  {
636  // 1D
637  case 1:
638  {
639  switch (fe_t.radial_family)
640  {
641  /*
642  * For no derivatives (and local coordinates, as
643  * given in \p p) the infinite element shapes
644  * are independent of mapping type
645  */
646  case INFINITE_MAP:
647  return InfFE<1,INFINITE_MAP,CARTESIAN>::shape(fe_t, t, i, p);
648 
649  case JACOBI_20_00:
650  return InfFE<1,JACOBI_20_00,CARTESIAN>::shape(fe_t, t, i, p);
651 
652  case JACOBI_30_00:
653  return InfFE<1,JACOBI_30_00,CARTESIAN>::shape(fe_t, t, i, p);
654 
655  case LEGENDRE:
656  return InfFE<1,LEGENDRE,CARTESIAN>::shape(fe_t, t, i, p);
657 
658  case LAGRANGE:
659  return InfFE<1,LAGRANGE,CARTESIAN>::shape(fe_t, t, i, p);
660 
661  default:
662  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
663  }
664  }
665 
666 
667  // 2D
668  case 2:
669  {
670  switch (fe_t.radial_family)
671  {
672  case INFINITE_MAP:
673  return InfFE<2,INFINITE_MAP,CARTESIAN>::shape(fe_t, t, i, p);
674 
675  case JACOBI_20_00:
676  return InfFE<2,JACOBI_20_00,CARTESIAN>::shape(fe_t, t, i, p);
677 
678  case JACOBI_30_00:
679  return InfFE<2,JACOBI_30_00,CARTESIAN>::shape(fe_t, t, i, p);
680 
681  case LEGENDRE:
682  return InfFE<2,LEGENDRE,CARTESIAN>::shape(fe_t, t, i, p);
683 
684  case LAGRANGE:
685  return InfFE<2,LAGRANGE,CARTESIAN>::shape(fe_t, t, i, p);
686 
687  default:
688  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
689  }
690  }
691 
692 
693  // 3D
694  case 3:
695  {
696  switch (fe_t.radial_family)
697  {
698  case INFINITE_MAP:
699  return InfFE<3,INFINITE_MAP,CARTESIAN>::shape(fe_t, t, i, p);
700 
701  case JACOBI_20_00:
702  return InfFE<3,JACOBI_20_00,CARTESIAN>::shape(fe_t, t, i, p);
703 
704  case JACOBI_30_00:
705  return InfFE<3,JACOBI_30_00,CARTESIAN>::shape(fe_t, t, i, p);
706 
707  case LEGENDRE:
708  return InfFE<3,LEGENDRE,CARTESIAN>::shape(fe_t, t, i, p);
709 
710  case LAGRANGE:
711  return InfFE<3,LAGRANGE,CARTESIAN>::shape(fe_t, t, i, p);
712 
713  default:
714  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
715  }
716  }
717 
718  default:
719  libmesh_error_msg("Invalid dim = " << dim);
720  }
721 
722  libmesh_error_msg("We'll never get here!");
723  return 0.;
724 }
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 729 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().

734 {
735  switch (dim)
736  {
737  // 1D
738  case 1:
739  {
740  switch (fe_t.radial_family)
741  {
742  /*
743  * For no derivatives (and local coordinates, as
744  * given in \p p) the infinite element shapes
745  * are independent of mapping type
746  */
747  case INFINITE_MAP:
748  return InfFE<1,INFINITE_MAP,CARTESIAN>::shape(fe_t, elem, i, p);
749 
750  case JACOBI_20_00:
751  return InfFE<1,JACOBI_20_00,CARTESIAN>::shape(fe_t, elem, i, p);
752 
753  case JACOBI_30_00:
754  return InfFE<1,JACOBI_30_00,CARTESIAN>::shape(fe_t, elem, i, p);
755 
756  case LEGENDRE:
757  return InfFE<1,LEGENDRE,CARTESIAN>::shape(fe_t, elem, i, p);
758 
759  case LAGRANGE:
760  return InfFE<1,LAGRANGE,CARTESIAN>::shape(fe_t, elem, i, p);
761 
762  default:
763  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
764  }
765  }
766 
767 
768  // 2D
769  case 2:
770  {
771  switch (fe_t.radial_family)
772  {
773  case INFINITE_MAP:
774  return InfFE<2,INFINITE_MAP,CARTESIAN>::shape(fe_t, elem, i, p);
775 
776  case JACOBI_20_00:
777  return InfFE<2,JACOBI_20_00,CARTESIAN>::shape(fe_t, elem, i, p);
778 
779  case JACOBI_30_00:
780  return InfFE<2,JACOBI_30_00,CARTESIAN>::shape(fe_t, elem, i, p);
781 
782  case LEGENDRE:
783  return InfFE<2,LEGENDRE,CARTESIAN>::shape(fe_t, elem, i, p);
784 
785  case LAGRANGE:
786  return InfFE<2,LAGRANGE,CARTESIAN>::shape(fe_t, elem, i, p);
787 
788  default:
789  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
790  }
791  }
792 
793 
794  // 3D
795  case 3:
796  {
797  switch (fe_t.radial_family)
798  {
799  case INFINITE_MAP:
800  return InfFE<3,INFINITE_MAP,CARTESIAN>::shape(fe_t, elem, i, p);
801 
802  case JACOBI_20_00:
803  return InfFE<3,JACOBI_20_00,CARTESIAN>::shape(fe_t, elem, i, p);
804 
805  case JACOBI_30_00:
806  return InfFE<3,JACOBI_30_00,CARTESIAN>::shape(fe_t, elem, i, p);
807 
808  case LEGENDRE:
809  return InfFE<3,LEGENDRE,CARTESIAN>::shape(fe_t, elem, i, p);
810 
811  case LAGRANGE:
812  return InfFE<3,LAGRANGE,CARTESIAN>::shape(fe_t, elem, i, p);
813 
814  default:
815  libmesh_error_msg("Invalid radial family = " << fe_t.radial_family);
816  }
817  }
818 
819  default:
820  libmesh_error_msg("Invalid dim = " << dim);
821  }
822 
823  libmesh_error_msg("We'll never get here!");
824  return 0.;
825 }
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 556 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::InfPrism6::contains_point(), libMesh::InfHex8::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(), libMesh::HPCoarsenTest::select_refinement(), and ~FEInterface().

562 {
563 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
564 
565  if ( is_InfFE_elem(elem->type()) )
566  return ifem_inverse_map(dim, fe_t, elem, p,tolerance, secure);
567 
568 #endif
569 
570  fe_with_vec_switch(inverse_map(elem, p, tolerance, secure));
571 
572  libmesh_error_msg("We'll never get here!");
573  return Point();
574 }
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:556
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:453
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 579 of file fe_interface.C.

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

586 {
587  const std::size_t n_pts = physical_points.size();
588 
589  // Resize the vector
590  reference_points.resize(n_pts);
591 
592  if (n_pts == 0)
593  {
594  libMesh::err << "WARNING: empty vector physical_points!"
595  << std::endl;
596  libmesh_here();
597  return;
598  }
599 
600 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
601 
602  if ( is_InfFE_elem(elem->type()) )
603  {
604  ifem_inverse_map(dim, fe_t, elem, physical_points, reference_points, tolerance, secure);
605  return;
606  // libmesh_not_implemented();
607  }
608 
609 #endif
610 
611  void_fe_with_vec_switch(inverse_map(elem, physical_points, reference_points, tolerance, secure));
612 
613  libmesh_error_msg("We'll never get here!");
614 }
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:556
OStreamProxy err(std::cerr)
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:453
bool libMesh::FEInterface::is_InfFE_elem ( const ElemType  et)
inlinestaticprivate
Returns
true if et is an element to be processed by class InfFE. Otherwise, it returns false. For compatibility with disabled infinite elements it always returns false.

Definition at line 453 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(), shape(), and ~FEInterface().

454 {
455  return false;
456 }
Point libMesh::FEInterface::map ( unsigned int  dim,
const FEType fe_t,
const Elem elem,
const Point p 
)
static

Returns the point in physical space of the reference point refpoint which is passed in.

Definition at line 537 of file fe_interface.C.

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

Referenced by libMesh::Elem::point_test(), and ~FEInterface().

541 {
542 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
543  if (is_InfFE_elem(elem->type()))
544  return ifem_map(dim, fe_t, elem, p);
545 #endif
546  fe_with_vec_switch(map(elem, p));
547 
548  libmesh_error_msg("We'll never get here!");
549  return Point();
550 }
static Point map(unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p)
Definition: fe_interface.C:537
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:453
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 1010 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::SUBDIVISION, libMesh::SZABAB, libMesh::TET10, libMesh::TET4, libMesh::TRI3, libMesh::TRI3SUBDIVISION, libMesh::TRI6, libMesh::TRISHELL3, and libMesh::XYZ.

Referenced by libMesh::DofMap::reinit(), and ~FEInterface().

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

417 {
418 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
419 
420  if ( is_InfFE_elem(t) )
421  return ifem_n_dofs(dim, fe_t, t);
422 
423 #endif
424 
425  const Order o = fe_t.order;
426 
427  fe_with_vec_switch(n_dofs(t, o));
428 
429  libmesh_error_msg("We'll never get here!");
430  return 0;
431 }
static unsigned int n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:414
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:453
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 436 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(), libMesh::DofMap::old_dof_indices(), libMesh::GenericProjector< FFunctor, GFunctor, FValue, ProjectionAction >::operator()(), libMesh::BoundaryProjectSolution::operator()(), libMesh::DofMap::reinit(), and ~FEInterface().

440 {
441 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
442 
443  if ( is_InfFE_elem(t) )
444  return ifem_n_dofs_at_node(dim, fe_t, t, n);
445 
446 #endif
447 
448  const Order o = fe_t.order;
449 
450  fe_with_vec_switch(n_dofs_at_node(t, o, n));
451 
452  libmesh_error_msg("We'll never get here!");
453  return 0;
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:436
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:453
static unsigned int ifem_n_dofs_at_node(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
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 460 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(), libMesh::DofMap::reinit(), and ~FEInterface().

463 {
464 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
465 
466  if ( is_InfFE_elem(t) )
467  return ifem_n_dofs_per_elem(dim, fe_t, t);
468 
469 #endif
470 
471  const Order o = fe_t.order;
472 
473  fe_with_vec_switch(n_dofs_per_elem(t, o));
474 
475  libmesh_error_msg("We'll never get here!");
476  return 0;
477 }
static unsigned int n_dofs_per_elem(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:460
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:453
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 
406  libmesh_error_msg("We'll never get here!");
407  return 0;
408 }
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:453
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 1416 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(), libMesh::EquationSystems::build_variable_names(), and ~FEInterface().

1418 {
1419  switch (fe_type.family)
1420  {
1421  //FIXME: We currently assume that the number of vector components is tied
1422  // to the mesh dimension. This will break for mixed-dimension meshes.
1423  case LAGRANGE_VEC:
1424  case NEDELEC_ONE:
1425  return mesh.mesh_dimension();
1426  default:
1427  return 1;
1428  }
1429 }
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 513 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(), libMesh::EnsightIO::write_vector_ascii(), and ~FEInterface().

518 {
519 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
520 
521  if ( is_InfFE_elem(elem->type()) )
522  {
523  ifem_nodal_soln(dim, fe_t, elem, elem_soln, nodal_soln);
524  return;
525  }
526 
527 #endif
528 
529  const Order order = fe_t.order;
530 
531  void_fe_with_vec_switch(nodal_soln(elem, order, elem_soln, nodal_soln));
532 }
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:453
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:513
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 here the parameter eps can be specified to indicate a tolerance. For example, $ \xi \le 1 $ becomes $ \xi \le 1 + \epsilon $.

Definition at line 618 of file fe_interface.C.

References libMesh::FEAbstract::on_reference_element().

Referenced by libMesh::InfQuad4::contains_point(), libMesh::InfPrism6::contains_point(), libMesh::InfHex8::contains_point(), libMesh::Elem::point_test(), and ~FEInterface().

621 {
622  return FEBase::on_reference_element(p,t,eps);
623 }
static bool on_reference_element(const Point &p, const ElemType t, const Real eps=TOLERANCE)
Definition: fe_abstract.C:556
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.

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

Definition at line 628 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(), libMesh::InfFE< Dim, T_radial, T_map >::shape(), and ~FEInterface().

633 {
634 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
635 
636  if ( is_InfFE_elem(t) )
637  return ifem_shape(dim, fe_t, t, i, p);
638 
639 #endif
640 
641  const Order o = fe_t.order;
642 
643  fe_switch(shape(t,o,i,p));
644 
645  libmesh_error_msg("We'll never get here!");
646  return 0.;
647 }
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:628
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:453
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.

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

Definition at line 649 of file fe_interface.C.

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

654 {
655 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
656 
657  if ( elem && is_InfFE_elem(elem->type()) )
658  return ifem_shape(dim, fe_t, elem, i, p);
659 
660 #endif
661 
662  const Order o = fe_t.order;
663 
664  fe_switch(shape(elem,o,i,p));
665 
666  libmesh_error_msg("We'll never get here!");
667  return 0.;
668 }
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:628
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:453
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.

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.

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 671 of file fe_interface.C.

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

677 {
678 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
679 
680  if ( is_InfFE_elem(t) )
681  phi = ifem_shape(dim, fe_t, t, i, p);
682 
683 #endif
684 
685  const Order o = fe_t.order;
686 
687  switch(dim)
688  {
689  case 0:
690  fe_scalar_vec_error_switch(0, shape(t,o,i,p), phi = , ; break;);
691  break;
692  case 1:
693  fe_scalar_vec_error_switch(1, shape(t,o,i,p), phi = , ; break;);
694  break;
695  case 2:
696  fe_scalar_vec_error_switch(2, shape(t,o,i,p), phi = , ; break;);
697  break;
698  case 3:
699  fe_scalar_vec_error_switch(3, shape(t,o,i,p), phi = , ; break;);
700  break;
701  default:
702  libmesh_error_msg("Invalid dimension = " << dim);
703  }
704 
705  return;
706 }
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:628
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:453
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 709 of file fe_interface.C.

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

715 {
716 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
717 
718  if ( elem && is_InfFE_elem(elem->type()) )
719  phi = ifem_shape(dim, fe_t, elem, i, p);
720 
721 #endif
722 
723  const Order o = fe_t.order;
724 
725  switch(dim)
726  {
727  case 0:
728  fe_scalar_vec_error_switch(0, shape(elem,o,i,p), phi = , ; break;);
729  break;
730  case 1:
731  fe_scalar_vec_error_switch(1, shape(elem,o,i,p), phi = , ; break;);
732  break;
733  case 2:
734  fe_scalar_vec_error_switch(2, shape(elem,o,i,p), phi = , ; break;);
735  break;
736  case 3:
737  fe_scalar_vec_error_switch(3, shape(elem,o,i,p), phi = , ; break;);
738  break;
739  default:
740  libmesh_error_msg("Invalid dimension = " << dim);
741  }
742 
743  return;
744 }
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:628
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:453
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 747 of file fe_interface.C.

References shape().

753 {
754  const Order o = fe_t.order;
755 
756  switch(dim)
757  {
758  case 0:
759  fe_vector_scalar_error_switch(0, shape(t,o,i,p), phi = , ; break;);
760  break;
761  case 1:
762  fe_vector_scalar_error_switch(1, shape(t,o,i,p), phi = , ; break;);
763  break;
764  case 2:
765  fe_vector_scalar_error_switch(2, shape(t,o,i,p), phi = , ; break;);
766  break;
767  case 3:
768  fe_vector_scalar_error_switch(3, shape(t,o,i,p), phi = , ; break;);
769  break;
770  default:
771  libmesh_error_msg("Invalid dimension = " << dim);
772  }
773 
774  return;
775 }
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:628
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 778 of file fe_interface.C.

References shape().

784 {
785  const Order o = fe_t.order;
786 
787  switch(dim)
788  {
789  case 0:
790  fe_vector_scalar_error_switch(0, shape(elem,o,i,p), phi = , ; break;);
791  break;
792  case 1:
793  fe_vector_scalar_error_switch(1, shape(elem,o,i,p), phi = , ; break;);
794  break;
795  case 2:
796  fe_vector_scalar_error_switch(2, shape(elem,o,i,p), phi = , ; break;);
797  break;
798  case 3:
799  fe_vector_scalar_error_switch(3, shape(elem,o,i,p), phi = , ; break;);
800  break;
801  default:
802  libmesh_error_msg("Invalid dimension = " << dim);
803  }
804 
805  return;
806 }
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:628

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