26 #include "libmesh/libmesh_common.h" 27 #include "libmesh/fe.h" 28 #include "libmesh/fe_interface.h" 29 #include "libmesh/quadrature.h" 30 #include "libmesh/elem.h" 31 #include "libmesh/libmesh_logging.h" 32 #include "libmesh/tensor_value.h" 40 #define REINIT_ERROR(_dim, _type, _func) \ 42 void FE<_dim,_type>::_func(const Elem *, \ 45 const std::vector<Point> * const, \ 46 const std::vector<Real> * const) 48 #define SIDEMAP_ERROR(_dim, _type, _func) \ 50 void FE<_dim,_type>::_func(const Elem *, \ 53 const std::vector<Point> &, \ 56 #define FACE_EDGE_SHAPE_ERROR(_dim, _func) \ 58 void FEMap::_func<_dim>(const std::vector<Point> &, \ 61 libmesh_error_msg("ERROR: This method makes no sense for low-D elements!"); \ 65 #define LIBMESH_ERRORS_IN_LOW_D(_type) \ 66 REINIT_ERROR(0, _type, reinit) { libmesh_error_msg("ERROR: Cannot reinit 0D " #_type " elements!"); } \ 67 REINIT_ERROR(0, _type, edge_reinit) { libmesh_error_msg("ERROR: Cannot edge_reinit 0D " #_type " elements!"); } \ 68 SIDEMAP_ERROR(0, _type, side_map) { libmesh_error_msg("ERROR: Cannot side_map 0D " #_type " elements!"); } \ 69 SIDEMAP_ERROR(0, _type, edge_map) { libmesh_error_msg("ERROR: Cannot edge_map 0D " #_type " elements!"); } \ 70 REINIT_ERROR(1, _type, edge_reinit) { libmesh_error_msg("ERROR: Cannot edge_reinit 1D " #_type " elements!"); } \ 71 SIDEMAP_ERROR(1, _type, edge_map) { libmesh_error_msg("ERROR: Cannot edge_map 1D " #_type " elements!"); } 92 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES 108 template <
unsigned int Dim, FEFamily T>
110 const unsigned int s,
112 const std::vector<Point> *
const pts,
113 const std::vector<Real> *
const weights)
123 if (this->calculating_nothing())
125 this->calculations_started =
true;
132 this->_fe_map->add_calculations();
133 this->_fe_map->get_JxW();
134 this->_fe_map->get_xyz();
135 this->determine_calculations();
142 unsigned int side_p_level = elem->
p_level();
151 this->shapes_on_quadrature =
false;
154 this->_fe_map->template init_face_shape_functions<Dim>(*pts, side.get());
157 if (weights !=
nullptr)
159 this->_fe_map->compute_face_map (Dim, *weights, side.get());
163 std::vector<Real> dummy_weights (pts->size(), 1.);
164 this->_fe_map->compute_face_map (Dim, dummy_weights, side.get());
172 this->qrule->init(side->type(), side_p_level);
174 if (this->qrule->shapes_need_reinit())
175 this->shapes_on_quadrature =
false;
180 if ((this->get_type() != elem->
type()) ||
181 (side->type() != last_side) ||
182 (this->_elem_p_level != side_p_level) ||
183 this->shapes_need_reinit() ||
184 !this->shapes_on_quadrature)
187 this->elem_type = elem->
type();
188 this->_elem_p_level = side_p_level;
191 last_side = side->type();
194 this->_p_level = this->_add_p_level_in_reinit * side_p_level;
197 this->_fe_map->template init_face_shape_functions<Dim>(this->qrule->get_points(), side.get());
201 this->_fe_map->compute_face_map (Dim, this->qrule->get_weights(), side.get());
204 this->shapes_on_quadrature =
true;
208 const std::vector<Real> JxW_int(this->_fe_map->get_JxW());
211 bool shapes_on_quadrature_side = this->shapes_on_quadrature;
215 const std::vector<Point> * ref_qp;
219 ref_qp = &this->qrule->get_points();
221 std::vector<Point> qp;
222 this->side_map(elem, side.get(), s, *ref_qp, qp);
226 this->reinit (elem, &qp);
228 this->shapes_on_quadrature = shapes_on_quadrature_side;
231 this->_fe_map->get_JxW() = JxW_int;
236 template <
unsigned int Dim, FEFamily T>
238 const unsigned int e,
240 const std::vector<Point> *
const pts,
241 const std::vector<Real> *
const weights)
246 libmesh_assert_not_equal_to (Dim, 1);
251 this->_fe_map->add_calculations();
252 this->_fe_map->get_JxW();
253 this->_fe_map->get_xyz();
254 this->determine_calculations();
264 this->shapes_on_quadrature =
false;
267 this->_fe_map->template init_edge_shape_functions<Dim> (*pts, edge.get());
270 if (weights !=
nullptr)
272 this->_fe_map->compute_edge_map (Dim, *weights, edge.get());
276 std::vector<Real> dummy_weights (pts->size(), 1.);
277 this->_fe_map->compute_edge_map (Dim, dummy_weights, edge.get());
285 this->qrule->init(edge->type(), elem->
p_level());
287 if (this->qrule->shapes_need_reinit())
288 this->shapes_on_quadrature =
false;
291 if ((this->get_type() != elem->
type()) ||
292 (edge->type() != last_edge) ||
293 this->shapes_need_reinit() ||
294 !this->shapes_on_quadrature)
297 this->elem_type = elem->
type();
300 last_edge = edge->type();
303 this->_fe_map->template init_edge_shape_functions<Dim> (this->qrule->get_points(), edge.get());
307 this->_fe_map->compute_edge_map (Dim, this->qrule->get_weights(), edge.get());
310 this->shapes_on_quadrature =
true;
314 const std::vector<Real> JxW_int(this->_fe_map->get_JxW());
318 const std::vector<Point> * ref_qp;
322 ref_qp = & this->qrule->get_points();
324 std::vector<Point> qp;
325 this->edge_map(elem, edge.get(), e, *ref_qp, qp);
329 this->reinit (elem, &qp);
332 this->_fe_map->get_JxW() = JxW_int;
335 template <
unsigned int Dim, FEFamily T>
338 const unsigned int s,
339 const std::vector<Point> & reference_side_points,
340 std::vector<Point> & reference_points)
343 this->calculate_phi =
true;
344 this->determine_calculations();
346 unsigned int side_p_level = elem->
p_level();
350 if (side->
type() != last_side ||
351 side_p_level != this->_elem_p_level ||
352 !this->shapes_on_quadrature)
355 this->elem_type = elem->
type();
356 this->_elem_p_level = side_p_level;
357 this->_p_level = this->_add_p_level_in_reinit * side_p_level;
360 last_side = side->
type();
363 this->_fe_map->template init_face_shape_functions<Dim>(reference_side_points, side);
366 const unsigned int n_points =
367 cast_int<unsigned int>(reference_side_points.size());
368 reference_points.resize(n_points);
369 for (
unsigned int i = 0; i < n_points; i++)
370 reference_points[i].
zero();
372 std::vector<Point> refspace_nodes;
373 this->get_refspace_nodes(elem->
type(), refspace_nodes);
375 const std::vector<std::vector<Real>> & psi_map = this->_fe_map->get_psi();
381 for (
unsigned int p=0; p<n_points; p++)
382 reference_points[p].add_scaled (side_node, psi_map[i][p]);
386 template <
unsigned int Dim, FEFamily T>
389 const unsigned int e,
390 const std::vector<Point> & reference_edge_points,
391 std::vector<Point> & reference_points)
394 this->calculate_phi =
true;
395 this->determine_calculations();
397 unsigned int edge_p_level = elem->
p_level();
399 if (edge->
type() != last_edge ||
400 edge_p_level != this->_elem_p_level ||
401 !this->shapes_on_quadrature)
404 this->elem_type = elem->
type();
405 this->_elem_p_level = edge_p_level;
406 this->_p_level = this->_add_p_level_in_reinit * edge_p_level;
409 last_edge = edge->
type();
412 this->_fe_map->template init_edge_shape_functions<Dim>(reference_edge_points, edge);
415 const unsigned int n_points =
416 cast_int<unsigned int>(reference_edge_points.size());
417 reference_points.resize(n_points);
418 for (
unsigned int i = 0; i < n_points; i++)
419 reference_points[i].
zero();
421 std::vector<Point> refspace_nodes;
422 this->get_refspace_nodes(elem->
type(), refspace_nodes);
424 const std::vector<std::vector<Real>> & psi_map = this->_fe_map->get_psi();
430 for (
unsigned int p=0; p<n_points; p++)
431 reference_points[p].add_scaled (edge_node, psi_map[i][p]);
436 template<
unsigned int Dim>
441 LOG_SCOPE(
"init_face_shape_functions()",
"FEMap");
454 const unsigned int n_qp = cast_int<unsigned int>(qp.size());
457 const unsigned int n_mapping_shape_functions =
463 this->
psi_map.resize (n_mapping_shape_functions);
468 this->
dpsidxi_map.resize (n_mapping_shape_functions);
469 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 479 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 494 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 499 for (
unsigned int i=0; i<n_mapping_shape_functions; i++)
504 this->
psi_map[i].resize (n_qp);
509 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 518 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 530 for (
unsigned int p=0; p<n_qp; p++)
538 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 552 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 564 template<
unsigned int Dim>
569 LOG_SCOPE(
"init_edge_shape_functions()",
"FEMap");
582 const unsigned int n_qp = cast_int<unsigned int>(qp.size());
585 const unsigned int n_mapping_shape_functions =
591 this->
psi_map.resize (n_mapping_shape_functions);
593 this->
dpsidxi_map.resize (n_mapping_shape_functions);
594 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 605 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 608 #endif // LIBMESH_ENABLE_SECOND_DERIVATIVES 610 for (
unsigned int i=0; i<n_mapping_shape_functions; i++)
615 this->
psi_map[i].resize (n_qp);
618 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 625 for (
unsigned int p=0; p<n_qp; p++)
631 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 649 LOG_SCOPE(
"compute_face_map()",
"FEMap");
652 const unsigned int n_qp = cast_int<unsigned int>(qw.size());
659 const unsigned int n_mapping_shape_functions =
673 this->
xyz.resize(n_qp);
678 this->
JxW.resize(n_qp);
695 libmesh_assert_equal_to (side->
node_id(0),
703 libmesh_assert_equal_to (this->
psi_map.size(), 1);
706 for (
unsigned int p=0; p<n_qp; p++)
716 this->
JxW[p] = 1.0*qw[p];
732 this->
xyz.resize(n_qp);
739 this->
JxW.resize(n_qp);
742 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 753 for (
unsigned int p=0; p<n_qp; p++)
759 this->
tangents[p].resize(LIBMESH_DIM-1);
762 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 769 for (
unsigned int i=0; i<n_mapping_shape_functions; i++)
773 for (
unsigned int p=0; p<n_qp; p++)
776 this->
xyz[p].add_scaled (side_point, this->
psi_map[i][p]);
779 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 789 for (
unsigned int p=0; p<n_qp; p++)
799 #elif LIBMESH_DIM == 3 822 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 836 libmesh_assert_not_equal_to (denominator, 0);
843 for (
unsigned int p=0; p<n_qp; p++)
847 libmesh_assert_greater (the_jac, 0.);
849 this->
JxW[p] = the_jac*qw[p];
865 this->
xyz.resize(n_qp);
872 this->
JxW.resize(n_qp);
874 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 886 for (
unsigned int p=0; p<n_qp; p++)
892 this->
tangents[p].resize(LIBMESH_DIM-1);
896 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 907 for (
unsigned int i=0; i<n_mapping_shape_functions; i++)
911 for (
unsigned int p=0; p<n_qp; p++)
914 this->
xyz[p].add_scaled (side_point, this->
psi_map[i][p]);
920 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 934 for (
unsigned int p=0; p<n_qp; p++)
941 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 955 const Real G = this->dxyzdeta_map[p].norm_sq();
957 const Real numerator = E*N -2.*F*M + G*L;
958 const Real denominator = E*G - F*F;
959 libmesh_assert_not_equal_to (denominator, 0.);
967 for (
unsigned int p=0; p<n_qp; p++)
977 const Real g21 = g12;
986 libmesh_assert_greater (the_jac, 0.);
988 this->
JxW[p] = the_jac*qw[p];
998 libmesh_error_msg(
"Invalid dimension dim = " <<
dim);
1006 const std::vector<Real> & qw,
1020 libmesh_assert_equal_to (
dim, 3);
1022 LOG_SCOPE(
"compute_edge_map()",
"FEMap");
1028 const unsigned int n_qp = cast_int<unsigned int>(qw.size());
1032 this->
xyz.resize(n_qp);
1039 this->
JxW.resize(n_qp);
1041 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 1052 for (
unsigned int p=0; p<n_qp; p++)
1055 this->
xyz[p].zero();
1062 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 1073 for (
unsigned int i=0,
1074 psi_map_size=cast_int<unsigned int>(
psi_map.size());
1075 i != psi_map_size; i++)
1079 for (
unsigned int p=0; p<n_qp; p++)
1082 this->
xyz[p].add_scaled (edge_point, this->
psi_map[i][p]);
1085 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 1095 for (
unsigned int p=0; p<n_qp; p++)
1105 libmesh_assert_greater (the_jac, 0.);
1107 this->
JxW[p] = the_jac*qw[p];
1113 FACE_EDGE_SHAPE_ERROR(0,init_face_shape_functions)
1114 template LIBMESH_EXPORT
void FEMap::init_face_shape_functions<1>(
const std::vector<Point> &,
const Elem *);
1115 template LIBMESH_EXPORT
void FEMap::init_face_shape_functions<2>(
const std::vector<Point> &,
const Elem *);
1116 template LIBMESH_EXPORT
void FEMap::init_face_shape_functions<3>(
const std::vector<Point> &,
const Elem *);
1118 FACE_EDGE_SHAPE_ERROR(0,init_edge_shape_functions)
1119 template LIBMESH_EXPORT
void FEMap::init_edge_shape_functions<1>(
const std::vector<Point> &,
const Elem *);
1120 template LIBMESH_EXPORT
void FEMap::init_edge_shape_functions<2>(
const std::vector<Point> &,
const Elem *);
1121 template LIBMESH_EXPORT
void FEMap::init_edge_shape_functions<3>(
const std::vector<Point> &,
const Elem *);
1125 #define REINIT_AND_SIDE_MAPS(_type) \ 1126 template LIBMESH_EXPORT void FE<1,_type>::reinit(Elem const *, unsigned int, Real, const std::vector<Point> * const, const std::vector<Real> * const); \ 1127 template LIBMESH_EXPORT void FE<1,_type>::side_map(Elem const *, Elem const *, const unsigned int, const std::vector<Point> &, std::vector<Point> &); \ 1128 template LIBMESH_EXPORT void FE<2,_type>::reinit(Elem const *, unsigned int, Real, const std::vector<Point> * const, const std::vector<Real> * const); \ 1129 template LIBMESH_EXPORT void FE<2,_type>::side_map(Elem const *, Elem const *, const unsigned int, const std::vector<Point> &, std::vector<Point> &); \ 1130 template LIBMESH_EXPORT void FE<2,_type>::edge_map(Elem const *, Elem const *, const unsigned int, const std::vector<Point> &, std::vector<Point> &); \ 1131 template LIBMESH_EXPORT void FE<2,_type>::edge_reinit(Elem const *, unsigned int, Real, const std::vector<Point> * const, const std::vector<Real> * const); \ 1132 template LIBMESH_EXPORT void FE<3,_type>::reinit(Elem const *, unsigned int, Real, const std::vector<Point> * const, const std::vector<Real> * const); \ 1133 template LIBMESH_EXPORT void FE<3,_type>::side_map(Elem const *, Elem const *, const unsigned int, const std::vector<Point> &, std::vector<Point> &); \ 1134 template LIBMESH_EXPORT void FE<3,_type>::edge_map(Elem const *, Elem const *, const unsigned int, const std::vector<Point> &, std::vector<Point> &); \ 1135 template LIBMESH_EXPORT void FE<3,_type>::edge_reinit(Elem const *, unsigned int, Real, const std::vector<Point> * const, const std::vector<Real> * const) 1146 REINIT_AND_SIDE_MAPS(
CLOUGH);
1147 REINIT_AND_SIDE_MAPS(
HERMITE);
1150 REINIT_AND_SIDE_MAPS(
SCALAR);
1151 REINIT_AND_SIDE_MAPS(
XYZ);
1153 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES 1155 REINIT_AND_SIDE_MAPS(
SZABAB);
Real(* shape_deriv_ptr)(const FEType fet, const Elem *elem, const unsigned int i, const unsigned int j, const Point &p, const bool add_p_level)
Typedef for pointer to a function that returns FE shape function derivative values.
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
void init_edge_shape_functions(const std::vector< Point > &qp, const Elem *edge)
Same as before, but for an edge.
Real dzdeta_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
std::vector< std::vector< Real > > dpsideta_map
Map for the derivative of the side function, d(psi)/d(eta).
std::vector< std::vector< Real > > psi_map
Map for the side shape functions, psi.
Order
defines an enum for polynomial orders.
std::vector< std::vector< Real > > d2psideta2_map
Map for the second derivatives (in eta) of the side shape functions.
virtual void side_map(const Elem *elem, const Elem *side, const unsigned int s, const std::vector< Point > &reference_side_points, std::vector< Point > &reference_points) override
Computes the reference space quadrature points on the side of an element based on the side quadrature...
static shape_ptr shape_function(const unsigned int dim, const FEType &fe_t, const ElemType t)
const Elem * interior_parent() const
bool calculate_dxyz
Should we calculate mapping gradients?
void compute_edge_map(int dim, const std::vector< Real > &qw, const Elem *side)
Same as before, but for an edge.
static Point inverse_map(const unsigned int dim, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true, const bool extra_checks=true)
Real dxdeta_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
Real dzdxi_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
This is the base class from which all geometric element types are derived.
unsigned int p_level() const
ADRealEigenVector< T, D, asd > sqrt(const ADRealEigenVector< T, D, asd > &)
The libMesh namespace provides an interface to certain functionality in the library.
virtual std::unique_ptr< Elem > build_side_ptr(const unsigned int i, bool proxy=false)=0
std::vector< std::vector< Real > > d2psidxideta_map
Map for the second (cross) derivatives in xi, eta of the side shape functions.
std::vector< RealGradient > d2xyzdxideta_map
Vector of mixed second partial derivatives in xi-eta: d^2(x)/d(xi)d(eta) d^2(y)/d(xi)d(eta) d^2(z)/d(...
virtual void reinit(const Elem *elem, const std::vector< Point > *const pts=nullptr, const std::vector< Real > *const weights=nullptr) override
This is at the core of this class.
std::vector< RealGradient > dxyzdxi_map
Vector of partial derivatives: d(x)/d(xi), d(y)/d(xi), d(z)/d(xi)
TypeVector< T > unit() const
std::vector< RealGradient > d2xyzdeta2_map
Vector of second partial derivatives in eta: d^2(x)/d(eta)^2.
REINIT_ERROR(1, RAVIART_THOMAS, reinit)
SIDEMAP_ERROR(1, NEDELEC_ONE, side_map)
static unsigned int n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
std::vector< RealGradient > d2xyzdxi2_map
Vector of second partial derivatives in xi: d^2(x)/d(xi)^2, d^2(y)/d(xi)^2, d^2(z)/d(xi)^2.
virtual unsigned int local_side_node(unsigned int side, unsigned int side_node) const =0
LIBMESH_ERRORS_IN_LOW_D(CLOUGH)
bool calculate_d2xyz
Should we calculate mapping hessians?
virtual unsigned int local_edge_node(unsigned int edge, unsigned int edge_node) const =0
Similar to Elem::local_side_node(), but instead of a side id, takes an edge id and a node id on that ...
bool calculate_xyz
Should we calculate physical point locations?
TypeVector< typename CompareTypes< T, T2 >::supertype > cross(const TypeVector< T2 > &v) const
std::vector< Real > curvatures
The mean curvature (= one half the sum of the principal curvatures) on the boundary at the quadrature...
std::vector< std::vector< Real > > d2psidxi2_map
Map for the second derivatives (in xi) of the side shape functions.
static shape_deriv_ptr shape_deriv_function(const unsigned int dim, const FEType &fe_t, const ElemType t)
std::vector< Real > JxW
Jacobian*Weight values at quadrature points.
std::vector< Point > normals
Normal vectors on boundary at quadrature points.
std::vector< Point > xyz
The spatial locations of the quadrature points.
const Elem * neighbor_ptr(unsigned int i) const
Real(* shape_ptr)(const FEType fe_t, const Elem *elem, const unsigned int i, const Point &p, const bool add_p_level)
Typedef for pointer to a function that returns FE shape function values.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Real dydeta_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
static shape_second_deriv_ptr shape_second_deriv_function(const unsigned int dim, const FEType &fe_t, const ElemType t)
virtual void edge_map(const Elem *elem, const Elem *edge, const unsigned int e, const std::vector< Point > &reference_edge_points, std::vector< Point > &reference_points)
Computes the reference space quadrature points on the side of an element based on the edge quadrature...
std::vector< std::vector< Real > > dpsidxi_map
Map for the derivative of the side functions, d(psi)/d(xi).
std::vector< std::vector< Point > > tangents
Tangent vectors on boundary at quadrature points.
std::vector< RealGradient > dxyzdeta_map
Vector of partial derivatives: d(x)/d(eta), d(y)/d(eta), d(z)/d(eta)
static Point map_deriv(const unsigned int dim, const Elem *elem, const unsigned int j, const Point &reference_point)
Real(* shape_second_deriv_ptr)(const FEType fet, const Elem *elem, const unsigned int i, const unsigned int j, const Point &p, const bool add_p_level)
Typedef for pointer to a function that returns FE shape function second derivative values...
Real dxdxi_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
void determine_calculations()
Determine which values are to be calculated.
Real dydxi_map(const unsigned int p) const
Used in FEMap::compute_map(), which should be be usable in derived classes, and therefore protected...
virtual void compute_face_map(int dim, const std::vector< Real > &qw, const Elem *side)
Same as compute_map, but for a side.
FEFamily
defines an enum for finite element families.
virtual std::unique_ptr< Elem > build_edge_ptr(const unsigned int i)=0
virtual Order default_order() const =0
virtual void edge_reinit(const Elem *elem, const unsigned int edge, const Real tolerance=TOLERANCE, const std::vector< Point > *const pts=nullptr, const std::vector< Real > *const weights=nullptr) override
Reinitializes all the physical element-dependent data based on the edge.
virtual ElemType type() const =0
A Point defines a location in LIBMESH_DIM dimensional Real space.
dof_id_type node_id(const unsigned int i) const
const Point & point(const unsigned int i) const
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
void init_face_shape_functions(const std::vector< Point > &qp, const Elem *side)
Initializes the reference to physical element map for a side.
static FEFamily map_fe_type(const Elem &elem)