libMesh::FEInterface Class Reference

#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_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

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, 2002-2007

Definition at line 64 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.

References libMesh::err.

33 {
34  libMesh::err << "ERROR: Do not define an object of this type."
35  << std::endl;
36  libmesh_error();
37 }
virtual libMesh::FEInterface::~FEInterface ( )
inlinevirtual

Destructor.

Definition at line 78 of file fe_interface.h.

78 {}

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

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

830 {
831  libmesh_assert(elem);
832 
833  const FEType& fe_t = dof_map.variable_type(variable_number);
834 
835  switch (elem->dim())
836  {
837  case 0:
838  case 1:
839  {
840  // No constraints in 0D/1D.
841  return;
842  }
843 
844 
845  case 2:
846  {
847  switch (fe_t.family)
848  {
849  case CLOUGH:
851  dof_map,
852  variable_number,
853  elem); return;
854 
855  case HERMITE:
857  dof_map,
858  variable_number,
859  elem); return;
860 
861  case LAGRANGE:
863  dof_map,
864  variable_number,
865  elem); return;
866 
867  case HIERARCHIC:
869  dof_map,
870  variable_number,
871  elem); return;
872 
873  case L2_HIERARCHIC:
875  dof_map,
876  variable_number,
877  elem); return;
878 
879  case LAGRANGE_VEC:
881  dof_map,
882  variable_number,
883  elem); return;
884 
885 
886 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
887  case SZABAB:
889  dof_map,
890  variable_number,
891  elem); return;
892 
893  case BERNSTEIN:
895  dof_map,
896  variable_number,
897  elem); return;
898 
899 #endif
900  default:
901  return;
902  }
903  }
904 
905 
906  case 3:
907  {
908  switch (fe_t.family)
909  {
910  case HERMITE:
912  dof_map,
913  variable_number,
914  elem); return;
915 
916  case LAGRANGE:
918  dof_map,
919  variable_number,
920  elem); return;
921 
922  case HIERARCHIC:
924  dof_map,
925  variable_number,
926  elem); return;
927 
928  case L2_HIERARCHIC:
930  dof_map,
931  variable_number,
932  elem); return;
933 
934  case LAGRANGE_VEC:
936  dof_map,
937  variable_number,
938  elem); return;
939 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
940  case SZABAB:
942  dof_map,
943  variable_number,
944  elem); return;
945 
946  case BERNSTEIN:
948  dof_map,
949  variable_number,
950  elem); return;
951 
952 #endif
953  default:
954  return;
955  }
956  }
957 
958 
959  default:
960  libmesh_error();
961  }
962 }
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 790 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::operator()().

794 {
795 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
796 
797  if ( is_InfFE_elem(elem->type()) )
798  {
799  data.init();
800  ifem_compute_data(dim, fe_t, elem, data);
801  return;
802  }
803 
804 #endif
805 
806  FEType p_refined = fe_t;
807  p_refined.order = static_cast<Order>(p_refined.order + elem->p_level());
808 
809  const unsigned int n_dof = n_dofs (dim, p_refined, elem->type());
810  const Point& p = data.p;
811  data.shape.resize(n_dof);
812 
813  // set default values for all the output fields
814  data.init();
815 
816  for (unsigned int n=0; n<n_dof; n++)
817  data.shape[n] = shape(dim, p_refined, elem, n, p);
818 
819  return;
820 }
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 970 of file fe_interface.C.

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

977 {
978  // No element-specific optimizations currently exist
980  dof_map,
981  boundaries,
982  mesh,
983  point_locator,
984  variable_number,
985  elem);
986 }
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 487 of file fe_interface.C.

References libMesh::FEType::order.

Referenced by libMesh::FEGenericBase< T >::coarsened_dof_values(), libMesh::ProjectSolution::operator()(), libMesh::ProjectFEMSolution::operator()(), and libMesh::BoundaryProjectSolution::operator()().

492 {
493  const Order o = fe_t.order;
494 
495  void_fe_with_vec_switch(dofs_on_edge(elem, o, e, di));
496 
497  libmesh_error();
498 }
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 472 of file fe_interface.C.

References libMesh::FEType::order.

Referenced by libMesh::FEGenericBase< T >::coarsened_dof_values(), libMesh::FEGenericBase< T >::compute_periodic_constraints(), libMesh::FEGenericBase< T >::compute_proj_constraints(), libMesh::ProjectSolution::operator()(), libMesh::ProjectFEMSolution::operator()(), and libMesh::BoundaryProjectSolution::operator()().

477 {
478  const Order o = fe_t.order;
479 
480  void_fe_with_vec_switch(dofs_on_side(elem, o, s, di));
481 
482  libmesh_error();
483 }
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 1319 of file fe_interface.C.

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

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

1320 {
1321  switch (fe_t.family)
1322  {
1323  case LAGRANGE:
1324  case L2_LAGRANGE:
1325  case MONOMIAL:
1326  case L2_HIERARCHIC:
1327  case XYZ:
1328  case LAGRANGE_VEC:
1329  case NEDELEC_ONE:
1330  return false;
1331  case CLOUGH:
1332  case HERMITE:
1333  case HIERARCHIC:
1334 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
1335  case BERNSTEIN:
1336  case SZABAB:
1337 #endif
1338  default:
1339  return true;
1340  }
1341 }
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 1348 of file fe_interface.C.

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

1349 {
1350  switch (fe_family)
1351  {
1352  case LAGRANGE_VEC:
1353  case NEDELEC_ONE:
1354  return TYPE_VECTOR;
1355  default:
1356  return TYPE_SCALAR;
1357  }
1358 }
void libMesh::FEInterface::ifem_compute_data ( const unsigned int  dim,
const FEType fe_t,
const Elem elem,
FEComputeData data 
)
staticprivate

Definition at line 871 of file fe_interface_inf_fe.C.

References libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::compute_data(), libMeshEnums::INFINITE_MAP, libMeshEnums::JACOBI_20_00, libMeshEnums::JACOBI_30_00, libMeshEnums::LAGRANGE, libMeshEnums::LEGENDRE, and libMesh::FEType::radial_family.

Referenced by compute_data().

875 {
876  switch (dim)
877  {
878  // 1D
879  case 1:
880  {
881  switch (fe_t.radial_family)
882  {
883  /*
884  * For no derivatives (and local coordinates, as
885  * given in \p p) the infinite element shapes
886  * are independent of mapping type
887  */
888  case INFINITE_MAP:
890  break;
891 
892  case JACOBI_20_00:
894  break;
895 
896  case JACOBI_30_00:
898  break;
899 
900  case LEGENDRE:
902  break;
903 
904  case LAGRANGE:
906  break;
907 
908  default:
909  libmesh_error();
910  }
911 
912  break;
913  }
914 
915 
916  // 2D
917  case 2:
918  {
919  switch (fe_t.radial_family)
920  {
921  case INFINITE_MAP:
923  break;
924 
925  case JACOBI_20_00:
927  break;
928 
929  case JACOBI_30_00:
931  break;
932 
933  case LEGENDRE:
935  break;
936 
937  case LAGRANGE:
939  break;
940 
941  default:
942  libmesh_error();
943  }
944 
945  break;
946  }
947 
948 
949  // 3D
950  case 3:
951  {
952  switch (fe_t.radial_family)
953  {
954  case INFINITE_MAP:
956  break;
957 
958  case JACOBI_20_00:
960  break;
961 
962  case JACOBI_30_00:
964  break;
965 
966  case LEGENDRE:
968  break;
969 
970  case LAGRANGE:
972  break;
973 
974  default:
975  libmesh_error();
976  }
977 
978  break;
979  }
980 
981 
982  default:
983  libmesh_error();
984  break;
985  }
986 
987  return;
988 }
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 470 of file fe_interface_inf_fe.C.

References libMeshEnums::CARTESIAN, libMeshEnums::ELLIPSOIDAL, libMesh::err, libMesh::FEType::inf_map, libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::inverse_map(), and libMeshEnums::SPHERICAL.

Referenced by inverse_map().

476 {
477  switch (dim)
478  {
479  // 1D
480  case 1:
481  {
482  switch (fe_t.inf_map)
483  {
484  case CARTESIAN:
485  return InfFE<1,JACOBI_20_00,CARTESIAN>::inverse_map(elem, p, tolerance, secure);
486 
487  case SPHERICAL:
488  case ELLIPSOIDAL:
489  {
490  libMesh::err << "ERROR: Spherical and Ellipsoidal IFEMs not (yet) " << std::endl
491  << "implemented." << std::endl;
492  libmesh_error();
493  }
494 
495 /*
496  case SPHERICAL:
497  return InfFE<1,JACOBI_20_00,SPHERICAL>::inverse_map(elem, p, tolerance);
498 
499  case ELLIPSOIDAL:
500  return InfFE<1,JACOBI_20_00,ELLIPSOIDAL>::inverse_map(elem, p, tolerance);
501 */
502 
503  default:
504  libmesh_error();
505  }
506  }
507 
508 
509  // 2D
510  case 2:
511  {
512  switch (fe_t.inf_map)
513  {
514  case CARTESIAN:
515  return InfFE<2,JACOBI_20_00,CARTESIAN>::inverse_map(elem, p, tolerance, secure);
516 
517  case SPHERICAL:
518  case ELLIPSOIDAL:
519  {
520  libMesh::err << "ERROR: Spherical and Ellipsoidal IFEMs not (yet) " << std::endl
521  << "implemented." << std::endl;
522  libmesh_error();
523  }
524 
525 /*
526  case SPHERICAL:
527  return InfFE<2,JACOBI_20_00,SPHERICAL>::inverse_map(elem, p, tolerance);
528 
529  case ELLIPSOIDAL:
530  return InfFE<2,JACOBI_20_00,ELLIPSOIDAL>::inverse_map(elem, p, tolerance);
531 */
532 
533  default:
534  libmesh_error();
535  }
536 
537  }
538 
539 
540  // 3D
541  case 3:
542  {
543  switch (fe_t.inf_map)
544  {
545  case CARTESIAN:
546  return InfFE<3,JACOBI_20_00,CARTESIAN>::inverse_map(elem, p, tolerance, secure);
547 
548  case SPHERICAL:
549  case ELLIPSOIDAL:
550  {
551  libMesh::err << "ERROR: Spherical and Ellipsoidal IFEMs not (yet) " << std::endl
552  << "implemented." << std::endl;
553  libmesh_error();
554  }
555 
556 /*
557  case SPHERICAL:
558  return InfFE<3,JACOBI_20_00,SPHERICAL>::inverse_map(elem, p, tolerance);
559 
560  case ELLIPSOIDAL:
561  return InfFE<3,JACOBI_20_00,ELLIPSOIDAL>::inverse_map(elem, p, tolerance);
562 */
563 
564  default:
565  libmesh_error();
566  }
567 
568  }
569 
570 
571  default:
572  libmesh_error();
573  }
574 
575 
576  libmesh_error();
577  Point pt;
578  return pt;
579 }
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 583 of file fe_interface_inf_fe.C.

References libMeshEnums::CARTESIAN, libMesh::FEType::inf_map, and libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::inverse_map().

590 {
591  switch (dim)
592  {
593  // 1D
594  case 1:
595  {
596  switch (fe_t.inf_map)
597  {
598  case CARTESIAN:
599  InfFE<1,JACOBI_20_00,CARTESIAN>::inverse_map(elem, physical_points, reference_points, tolerance, secure);
600  return;
601 
602  default:
603  libmesh_error();
604  }
605  }
606 
607 
608  // 2D
609  case 2:
610  {
611  switch (fe_t.inf_map)
612  {
613  case CARTESIAN:
614  InfFE<2,JACOBI_20_00,CARTESIAN>::inverse_map(elem, physical_points, reference_points, tolerance, secure);
615  return;
616 
617  default:
618  libmesh_error();
619  }
620 
621  }
622 
623 
624  // 3D
625  case 3:
626  {
627  switch (fe_t.inf_map)
628  {
629  case CARTESIAN:
630  InfFE<3,JACOBI_20_00,CARTESIAN>::inverse_map(elem, physical_points, reference_points, tolerance, secure);
631  return;
632 
633  default:
634  libmesh_error();
635  }
636 
637  }
638 
639 
640  default:
641  libmesh_error();
642  }
643 
644  libmesh_error();
645  return;
646 }
unsigned int libMesh::FEInterface::ifem_n_dofs ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t 
)
staticprivate

Definition at line 75 of file fe_interface_inf_fe.C.

References libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::n_dofs().

Referenced by n_dofs().

78 {
79  switch (dim)
80  {
81  // 1D
82  case 1:
83  /*
84  * Since InfFE<Dim,T_radial,T_map>::n_dofs(...)
85  * is actually independent of T_radial and T_map, we can use
86  * just any T_radial and T_map
87  */
89 
90  // 2D
91  case 2:
93 
94  // 3D
95  case 3:
97 
98  default:
99  libmesh_error();
100  }
101 
102 
103  libmesh_error();
104  return 0;
105 }
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 110 of file fe_interface_inf_fe.C.

References libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::n_dofs_at_node().

Referenced by n_dofs_at_node().

114 {
115  switch (dim)
116  {
117  // 1D
118  case 1:
119  /*
120  * Since InfFE<Dim,T_radial,T_map>::n_dofs_at_node(...)
121  * is actually independent of T_radial and T_map, we can use
122  * just any T_radial and T_map
123  */
125 
126  // 2D
127  case 2:
129 
130  // 3D
131  case 3:
133 
134  default:
135  libmesh_error();
136  }
137 
138 
139  libmesh_error();
140  return 0;
141 }
unsigned int libMesh::FEInterface::ifem_n_dofs_per_elem ( const unsigned int  dim,
const FEType fe_t,
const ElemType  t 
)
staticprivate

Definition at line 147 of file fe_interface_inf_fe.C.

References libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::n_dofs_per_elem().

Referenced by n_dofs_per_elem().

150 {
151  switch (dim)
152  {
153  // 1D
154  case 1:
155  /*
156  * Since InfFE<Dim,T_radial,T_map>::n_dofs(...)
157  * is actually independent of T_radial and T_map, we can use
158  * just any T_radial and T_map
159  */
161 
162  // 2D
163  case 2:
165 
166  // 3D
167  case 3:
169 
170  default:
171  libmesh_error();
172  }
173 
174 
175  libmesh_error();
176  return 0;
177 }
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< friend_Dim, friend_T_radial, friend_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();
64  }
65 
66 
67  libmesh_error();
68  return 0;
69 }
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 182 of file fe_interface_inf_fe.C.

References libMeshEnums::CARTESIAN, libMesh::err, libMesh::FEType::inf_map, libMeshEnums::INFINITE_MAP, libMeshEnums::JACOBI_20_00, libMeshEnums::JACOBI_30_00, libMeshEnums::LAGRANGE, libMeshEnums::LEGENDRE, libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::nodal_soln(), and libMesh::FEType::radial_family.

Referenced by nodal_soln().

187 {
188  switch (dim)
189  {
190 
191  // 1D
192  case 1:
193  {
194  switch (fe_t.radial_family)
195  {
196  case INFINITE_MAP:
197  {
198  libMesh::err << "ERROR: INFINTE_MAP is not a valid shape family for radial approximation." << std::endl;
199  libmesh_error();
200  break;
201  }
202 
203  case JACOBI_20_00:
204  {
205  switch (fe_t.inf_map)
206  {
207  case CARTESIAN:
208  {
210  break;
211  }
212  default:
213  libMesh::err << "ERROR: Spherical & Ellipsoidal IFEMs not implemented." << std::endl;
214  libmesh_error();
215  }
216  break;
217  }
218 
219  case JACOBI_30_00:
220  {
221  switch (fe_t.inf_map)
222  {
223  case CARTESIAN:
224  {
226  break;
227  }
228  default:
229  libMesh::err << "ERROR: Spherical & Ellipsoidal IFEMs not implemented." << std::endl;
230  libmesh_error();
231  }
232  break;
233  }
234 
235  case LEGENDRE:
236  {
237  switch (fe_t.inf_map)
238  {
239  case CARTESIAN:
240  {
241  InfFE<1,LEGENDRE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
242  break;
243  }
244  default:
245  libMesh::err << "ERROR: Spherical & Ellipsoidal IFEMs not implemented." << std::endl;
246  libmesh_error();
247  }
248  break;
249  }
250 
251  case LAGRANGE:
252  {
253  switch (fe_t.inf_map)
254  {
255  case CARTESIAN:
256  {
257  InfFE<1,LAGRANGE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
258  break;
259  }
260  default:
261  libMesh::err << "ERROR: Spherical & Ellipsoidal IFEMs not implemented." << std::endl;
262  libmesh_error();
263  }
264  break;
265  }
266 
267 
268 
269  default:
270  libMesh::err << "ERROR: Bad FEType.radial_family= " << fe_t.radial_family << std::endl;
271  libmesh_error();
272  break;
273  }
274 
275  break;
276  }
277 
278 
279 
280 
281  // 2D
282  case 2:
283  {
284  switch (fe_t.radial_family)
285  {
286  case INFINITE_MAP:
287  {
288  libMesh::err << "ERROR: INFINTE_MAP is not a valid shape family for radial approximation." << std::endl;
289  libmesh_error();
290  break;
291  }
292 
293  case JACOBI_20_00:
294  {
295  switch (fe_t.inf_map)
296  {
297  case CARTESIAN:
298  {
300  break;
301  }
302  default:
303  libMesh::err << "ERROR: Spherical & Ellipsoidal IFEMs not implemented." << std::endl;
304  libmesh_error();
305  }
306  break;
307  }
308 
309  case JACOBI_30_00:
310  {
311  switch (fe_t.inf_map)
312  {
313  case CARTESIAN:
314  {
316  break;
317  }
318  default:
319  libMesh::err << "ERROR: Spherical & Ellipsoidal IFEMs not implemented." << std::endl;
320  libmesh_error();
321  }
322  break;
323  }
324 
325  case LEGENDRE:
326  {
327  switch (fe_t.inf_map)
328  {
329  case CARTESIAN:
330  {
331  InfFE<2,LEGENDRE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
332  break;
333  }
334  default:
335  libMesh::err << "ERROR: Spherical & Ellipsoidal IFEMs not implemented." << std::endl;
336  libmesh_error();
337  }
338  break;
339  }
340 
341  case LAGRANGE:
342  {
343  switch (fe_t.inf_map)
344  {
345  case CARTESIAN:
346  {
347  InfFE<2,LAGRANGE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
348  break;
349  }
350  default:
351  libMesh::err << "ERROR: Spherical & Ellipsoidal IFEMs not implemented." << std::endl;
352  libmesh_error();
353  }
354  break;
355  }
356 
357 
358 
359  default:
360  libMesh::err << "ERROR: Bad FEType.radial_family= " << fe_t.radial_family << std::endl;
361  libmesh_error();
362  break;
363  }
364 
365  break;
366  }
367 
368 
369 
370 
371  // 3D
372  case 3:
373  {
374  switch (fe_t.radial_family)
375  {
376  case INFINITE_MAP:
377  {
378  libMesh::err << "ERROR: INFINTE_MAP is not a valid shape family for radial approximation." << std::endl;
379  libmesh_error();
380  break;
381  }
382 
383  case JACOBI_20_00:
384  {
385  switch (fe_t.inf_map)
386  {
387  case CARTESIAN:
388  {
390  break;
391  }
392  default:
393  libMesh::err << "ERROR: Spherical & Ellipsoidal IFEMs not implemented." << std::endl;
394  libmesh_error();
395  }
396  break;
397  }
398 
399  case JACOBI_30_00:
400  {
401  switch (fe_t.inf_map)
402  {
403  case CARTESIAN:
404  {
406  break;
407  }
408  default:
409  libMesh::err << "ERROR: Spherical & Ellipsoidal IFEMs not implemented." << std::endl;
410  libmesh_error();
411  }
412  break;
413  }
414 
415  case LEGENDRE:
416  {
417  switch (fe_t.inf_map)
418  {
419  case CARTESIAN:
420  {
421  InfFE<3,LEGENDRE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
422  break;
423  }
424  default:
425  libMesh::err << "ERROR: Spherical & Ellipsoidal IFEMs not implemented." << std::endl;
426  libmesh_error();
427  }
428  break;
429  }
430 
431  case LAGRANGE:
432  {
433  switch (fe_t.inf_map)
434  {
435  case CARTESIAN:
436  {
437  InfFE<3,LAGRANGE,CARTESIAN>::nodal_soln(fe_t, elem, elem_soln, nodal_soln);
438  break;
439  }
440  default:
441  libMesh::err << "ERROR: Spherical & Ellipsoidal IFEMs not implemented." << std::endl;
442  libmesh_error();
443  }
444  break;
445  }
446 
447 
448 
449  default:
450  libMesh::err << "ERROR: Bad FEType.radial_family= " << fe_t.radial_family << std::endl;
451  libmesh_error();
452  break;
453  }
454 
455  break;
456  }
457 
458  default:
459  libmesh_error();
460  }
461  return;
462 }
bool libMesh::FEInterface::ifem_on_reference_element ( const Point p,
const ElemType  t,
const Real  eps 
)
staticprivate

Definition at line 651 of file fe_interface_inf_fe.C.

References libMesh::FEAbstract::on_reference_element().

654 {
655  return FEBase::on_reference_element(p,t,eps);
656 }
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 661 of file fe_interface_inf_fe.C.

References libMeshEnums::INFINITE_MAP, libMeshEnums::JACOBI_20_00, libMeshEnums::JACOBI_30_00, libMeshEnums::LAGRANGE, libMeshEnums::LEGENDRE, libMesh::FEType::radial_family, and libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::shape().

Referenced by shape().

666 {
667  switch (dim)
668  {
669  // 1D
670  case 1:
671  {
672  switch (fe_t.radial_family)
673  {
674  /*
675  * For no derivatives (and local coordinates, as
676  * given in \p p) the infinite element shapes
677  * are independent of mapping type
678  */
679  case INFINITE_MAP:
680  return InfFE<1,INFINITE_MAP,CARTESIAN>::shape(fe_t, t, i, p);
681 
682  case JACOBI_20_00:
683  return InfFE<1,JACOBI_20_00,CARTESIAN>::shape(fe_t, t, i, p);
684 
685  case JACOBI_30_00:
686  return InfFE<1,JACOBI_30_00,CARTESIAN>::shape(fe_t, t, i, p);
687 
688  case LEGENDRE:
689  return InfFE<1,LEGENDRE,CARTESIAN>::shape(fe_t, t, i, p);
690 
691  case LAGRANGE:
692  return InfFE<1,LAGRANGE,CARTESIAN>::shape(fe_t, t, i, p);
693 
694  default:
695  libmesh_error();
696  }
697  }
698 
699 
700  // 2D
701  case 2:
702  {
703  switch (fe_t.radial_family)
704  {
705  case INFINITE_MAP:
706  return InfFE<2,INFINITE_MAP,CARTESIAN>::shape(fe_t, t, i, p);
707 
708  case JACOBI_20_00:
709  return InfFE<2,JACOBI_20_00,CARTESIAN>::shape(fe_t, t, i, p);
710 
711  case JACOBI_30_00:
712  return InfFE<2,JACOBI_30_00,CARTESIAN>::shape(fe_t, t, i, p);
713 
714  case LEGENDRE:
715  return InfFE<2,LEGENDRE,CARTESIAN>::shape(fe_t, t, i, p);
716 
717  case LAGRANGE:
718  return InfFE<2,LAGRANGE,CARTESIAN>::shape(fe_t, t, i, p);
719 
720  default:
721  libmesh_error();
722  }
723 
724  }
725 
726 
727  // 3D
728  case 3:
729  {
730  switch (fe_t.radial_family)
731  {
732  case INFINITE_MAP:
733  return InfFE<3,INFINITE_MAP,CARTESIAN>::shape(fe_t, t, i, p);
734 
735  case JACOBI_20_00:
736  return InfFE<3,JACOBI_20_00,CARTESIAN>::shape(fe_t, t, i, p);
737 
738  case JACOBI_30_00:
739  return InfFE<3,JACOBI_30_00,CARTESIAN>::shape(fe_t, t, i, p);
740 
741  case LEGENDRE:
742  return InfFE<3,LEGENDRE,CARTESIAN>::shape(fe_t, t, i, p);
743 
744  case LAGRANGE:
745  return InfFE<3,LAGRANGE,CARTESIAN>::shape(fe_t, t, i, p);
746 
747  default:
748  libmesh_error();
749  }
750 
751  }
752 
753 
754  default:
755  libmesh_error();
756  }
757 
758 
759  libmesh_error();
760  return 0.;
761 }
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 766 of file fe_interface_inf_fe.C.

References libMeshEnums::INFINITE_MAP, libMeshEnums::JACOBI_20_00, libMeshEnums::JACOBI_30_00, libMeshEnums::LAGRANGE, libMeshEnums::LEGENDRE, libMesh::FEType::radial_family, and libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::shape().

771 {
772  switch (dim)
773  {
774  // 1D
775  case 1:
776  {
777  switch (fe_t.radial_family)
778  {
779  /*
780  * For no derivatives (and local coordinates, as
781  * given in \p p) the infinite element shapes
782  * are independent of mapping type
783  */
784  case INFINITE_MAP:
785  return InfFE<1,INFINITE_MAP,CARTESIAN>::shape(fe_t, elem, i, p);
786 
787  case JACOBI_20_00:
788  return InfFE<1,JACOBI_20_00,CARTESIAN>::shape(fe_t, elem, i, p);
789 
790  case JACOBI_30_00:
791  return InfFE<1,JACOBI_30_00,CARTESIAN>::shape(fe_t, elem, i, p);
792 
793  case LEGENDRE:
794  return InfFE<1,LEGENDRE,CARTESIAN>::shape(fe_t, elem, i, p);
795 
796  case LAGRANGE:
797  return InfFE<1,LAGRANGE,CARTESIAN>::shape(fe_t, elem, i, p);
798 
799  default:
800  libmesh_error();
801  }
802  }
803 
804 
805  // 2D
806  case 2:
807  {
808  switch (fe_t.radial_family)
809  {
810  case INFINITE_MAP:
811  return InfFE<2,INFINITE_MAP,CARTESIAN>::shape(fe_t, elem, i, p);
812 
813  case JACOBI_20_00:
814  return InfFE<2,JACOBI_20_00,CARTESIAN>::shape(fe_t, elem, i, p);
815 
816  case JACOBI_30_00:
817  return InfFE<2,JACOBI_30_00,CARTESIAN>::shape(fe_t, elem, i, p);
818 
819  case LEGENDRE:
820  return InfFE<2,LEGENDRE,CARTESIAN>::shape(fe_t, elem, i, p);
821 
822  case LAGRANGE:
823  return InfFE<2,LAGRANGE,CARTESIAN>::shape(fe_t, elem, i, p);
824 
825  default:
826  libmesh_error();
827  }
828 
829  }
830 
831 
832  // 3D
833  case 3:
834  {
835  switch (fe_t.radial_family)
836  {
837  case INFINITE_MAP:
838  return InfFE<3,INFINITE_MAP,CARTESIAN>::shape(fe_t, elem, i, p);
839 
840  case JACOBI_20_00:
841  return InfFE<3,JACOBI_20_00,CARTESIAN>::shape(fe_t, elem, i, p);
842 
843  case JACOBI_30_00:
844  return InfFE<3,JACOBI_30_00,CARTESIAN>::shape(fe_t, elem, i, p);
845 
846  case LEGENDRE:
847  return InfFE<3,LEGENDRE,CARTESIAN>::shape(fe_t, elem, i, p);
848 
849  case LAGRANGE:
850  return InfFE<3,LAGRANGE,CARTESIAN>::shape(fe_t, elem, i, p);
851 
852  default:
853  libmesh_error();
854  }
855 
856  }
857 
858 
859  default:
860  libmesh_error();
861  }
862 
863 
864  libmesh_error();
865  return 0.;
866 }
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 542 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< T >::coarsened_dof_values(), libMesh::FEAbstract::compute_node_constraints(), libMesh::FEGenericBase< T >::compute_periodic_constraints(), libMesh::FEAbstract::compute_periodic_node_constraints(), libMesh::FEGenericBase< T >::compute_proj_constraints(), libMesh::InfQuad4::contains_point(), libMesh::InfPrism6::contains_point(), libMesh::InfHex8::contains_point(), libMesh::MeshFunction::gradient(), libMesh::MeshFunction::hessian(), inverse_map(), libMesh::DGFEMContext::neighbor_side_fe_reinit(), 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().

548 {
549 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
550 
551  if ( is_InfFE_elem(elem->type()) )
552  return ifem_inverse_map(dim, fe_t, elem, p,tolerance, secure);
553 
554 #endif
555 
556  fe_with_vec_switch(inverse_map(elem, p, tolerance, secure));
557 
558  libmesh_error();
559  return Point();
560 }
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 565 of file fe_interface.C.

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

572 {
573  const std::size_t n_pts = physical_points.size();
574 
575  // Resize the vector
576  reference_points.resize(n_pts);
577 
578  if (n_pts == 0)
579  {
580  libMesh::err << "WARNING: empty vector physical_points!"
581  << std::endl;
582  libmesh_here();
583  return;
584  }
585 
586 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
587 
588  if ( is_InfFE_elem(elem->type()) )
589  {
590  ifem_inverse_map(dim, fe_t, elem, physical_points, reference_points, tolerance, secure);
591  return;
592 
593 // libMesh::err << "ERROR: Not implemented!"
594 // << std::endl;
595 // libmesh_error();
596  }
597 
598 #endif
599 
600  void_fe_with_vec_switch(inverse_map(elem, physical_points, reference_points, tolerance, secure));
601 
602  libmesh_error();
603  return;
604 }
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 449 of file fe_interface.h.

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

450 {
451  return false;
452 }
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 527 of file fe_interface.C.

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

531 {
532  fe_with_vec_switch(map(elem, p));
533 
534  libmesh_error();
535  return Point();
536 }
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 992 of file fe_interface.C.

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

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

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

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

Referenced by libMesh::FEGenericBase< T >::coarsened_dof_values(), compute_data(), libMesh::FEAbstract::compute_node_constraints(), libMesh::FEAbstract::compute_periodic_node_constraints(), libMesh::DofMap::dof_indices(), libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::n_dofs(), libMesh::DofMap::old_dof_indices(), and libMesh::HPCoarsenTest::select_refinement().

407 {
408 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
409 
410  if ( is_InfFE_elem(t) )
411  return ifem_n_dofs(dim, fe_t, t);
412 
413 #endif
414 
415  const Order o = fe_t.order;
416 
417  fe_with_vec_switch(n_dofs(t, o));
418 
419  libmesh_error();
420  return 0;
421 }
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 426 of file fe_interface.C.

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

Referenced by libMesh::FEGenericBase< T >::coarsened_dof_values(), libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::compute_shape_indices(), libMesh::DofMap::dof_indices(), libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::n_dofs_at_node(), libMesh::DofMap::old_dof_indices(), libMesh::ProjectSolution::operator()(), libMesh::ProjectFEMSolution::operator()(), libMesh::BoundaryProjectSolution::operator()(), and libMesh::DofMap::reinit().

430 {
431 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
432 
433  if ( is_InfFE_elem(t) )
434  return ifem_n_dofs_at_node(dim, fe_t, t, n);
435 
436 #endif
437 
438  const Order o = fe_t.order;
439 
440  fe_with_vec_switch(n_dofs_at_node(t, o, n));
441 
442  libmesh_error();
443  return 0;
444 }
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 450 of file fe_interface.C.

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

Referenced by libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::compute_shape_indices(), libMesh::DofMap::dof_indices(), libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::n_dofs_per_elem(), libMesh::DofMap::old_dof_indices(), and libMesh::DofMap::reinit().

453 {
454 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
455 
456  if ( is_InfFE_elem(t) )
457  return ifem_n_dofs_per_elem(dim, fe_t, t);
458 
459 #endif
460 
461  const Order o = fe_t.order;
462 
463  fe_with_vec_switch(n_dofs_per_elem(t, o));
464 
465  libmesh_error();
466  return 0;
467 }
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 375 of file fe_interface.C.

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

378 {
379 
380 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
381  /*
382  * Since the FEType, stored in DofMap/(some System child), has to
383  * be the _same_ for InfFE and FE, we have to catch calls
384  * to infinite elements through the element type.
385  */
386 
387  if ( is_InfFE_elem(t) )
388  return ifem_n_shape_functions(dim, fe_t, t);
389 
390 #endif
391 
392  const Order o = fe_t.order;
393 
394  fe_with_vec_switch(n_shape_functions(t, o));
395 
396  libmesh_error();
397  return 0;
398 }
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 1360 of file fe_interface.C.

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

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

1361 {
1362  switch (fe_type.family)
1363  {
1364  //FIXME: We currently assume that the number of vector components is tied
1365  // to the mesh dimension. This will break for mixed-dimension meshes.
1366  case LAGRANGE_VEC:
1367  case NEDELEC_ONE:
1368  return mesh.mesh_dimension();
1369  default:
1370  return 1;
1371  }
1372 }
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 503 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_solution_vector(), libMesh::EnsightIO::write_scalar_ascii(), and libMesh::EnsightIO::write_vector_ascii().

508 {
509 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
510 
511  if ( is_InfFE_elem(elem->type()) )
512  {
513  ifem_nodal_soln(dim, fe_t, elem, elem_soln, nodal_soln);
514  return;
515  }
516 
517 #endif
518 
519  const Order order = fe_t.order;
520 
521  void_fe_with_vec_switch(nodal_soln(elem, order, elem_soln, nodal_soln));
522 }
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 608 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(), and libMesh::Elem::point_test().

611 {
612  return FEBase::on_reference_element(p,t,eps);
613 }
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 618 of file fe_interface.C.

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

Referenced by compute_data(), libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::compute_data(), libMesh::FEAbstract::compute_node_constraints(), libMesh::FEAbstract::compute_periodic_node_constraints(), shape(), and libMesh::InfFE< friend_Dim, friend_T_radial, friend_T_map >::shape().

623 {
624 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
625 
626  if ( is_InfFE_elem(t) )
627  return ifem_shape(dim, fe_t, t, i, p);
628 
629 #endif
630 
631  const Order o = fe_t.order;
632 
633  fe_switch(shape(t,o,i,p));
634 
635  libmesh_error();
636  return 0.;
637 }
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 639 of file fe_interface.C.

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

644 {
645 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
646 
647  if ( is_InfFE_elem(elem->type()) )
648  return ifem_shape(dim, fe_t, elem, i, p);
649 
650 #endif
651 
652  const Order o = fe_t.order;
653 
654  fe_switch(shape(elem,o,i,p));
655 
656  libmesh_error();
657  return 0.;
658 }
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 661 of file fe_interface.C.

References libMesh::dim.

667 {
668 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
669 
670  if ( is_InfFE_elem(t) )
671  phi = ifem_shape(dim, fe_t, t, i, p);
672 
673 #endif
674 
675  const Order o = fe_t.order;
676 
677  switch(dim)
678  {
679  case 0:
680  fe_scalar_vec_error_switch(0, shape(t,o,i,p), phi = , ; break;);
681  break;
682  case 1:
683  fe_scalar_vec_error_switch(1, shape(t,o,i,p), phi = , ; break;);
684  break;
685  case 2:
686  fe_scalar_vec_error_switch(2, shape(t,o,i,p), phi = , ; break;);
687  break;
688  case 3:
689  fe_scalar_vec_error_switch(3, shape(t,o,i,p), phi = , ; break;);
690  break;
691  }
692 
693  return;
694 }
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 697 of file fe_interface.C.

References libMesh::dim.

703 {
704 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
705 
706  if ( is_InfFE_elem(elem->type()) )
707  phi = ifem_shape(dim, fe_t, elem, i, p);
708 
709 #endif
710 
711  const Order o = fe_t.order;
712 
713  switch(dim)
714  {
715  case 0:
716  fe_scalar_vec_error_switch(0, shape(elem,o,i,p), phi = , ; break;);
717  break;
718  case 1:
719  fe_scalar_vec_error_switch(1, shape(elem,o,i,p), phi = , ; break;);
720  break;
721  case 2:
722  fe_scalar_vec_error_switch(2, shape(elem,o,i,p), phi = , ; break;);
723  break;
724  case 3:
725  fe_scalar_vec_error_switch(3, shape(elem,o,i,p), phi = , ; break;);
726  break;
727  }
728 
729  return;
730 }
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 733 of file fe_interface.C.

References libMesh::dim.

739 {
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  }
757 
758  return;
759 }
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 762 of file fe_interface.C.

References libMesh::dim.

768 {
769  const Order o = fe_t.order;
770 
771  switch(dim)
772  {
773  case 0:
774  fe_vector_scalar_error_switch(0, shape(elem,o,i,p), phi = , ; break;);
775  break;
776  case 1:
777  fe_vector_scalar_error_switch(1, shape(elem,o,i,p), phi = , ; break;);
778  break;
779  case 2:
780  fe_vector_scalar_error_switch(2, shape(elem,o,i,p), phi = , ; break;);
781  break;
782  case 3:
783  fe_vector_scalar_error_switch(3, shape(elem,o,i,p), phi = , ; break;);
784  break;
785  }
786 
787  return;
788 }

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

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

Hosted By:
SourceForge.net Logo