19 #include "libmesh/fe.h" 20 #include "libmesh/elem.h" 21 #include "libmesh/fe_interface.h" 22 #include "libmesh/utility.h" 23 #include "libmesh/enum_to_string.h" 30 void hermite_compute_coefs(
const Elem * elem,
Real & d1xd1x,
Real & d2xd2x)
37 const int n_mapping_shape_functions =
44 std::vector<Real> dxdxi(2);
45 std::vector<Real> dxidx(2);
50 for (
int p = 0; p != 2; ++p)
53 for (
int i = 0; i != n_mapping_shape_functions; ++i)
55 const Real ddxi = shape_deriv_ptr
56 (map_fe_type, elem, i, 0, dofpt[p],
false);
57 dxdxi[p] += elem->
point(i)(0) * ddxi;
77 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 91 return 0.5 * (-1. + 3.*xi);
93 return 0.5 * (1. + 3.*xi);
95 return (8.*xi*xi + 4.*(xi*xi-1.))/24.;
97 return (8.*xi*xi*xi + 12.*xi*(xi*xi-1.))/120.;
102 Real denominator = 720., xipower = 1.;
103 for (
unsigned n=6; n != i; ++n)
106 denominator *= (n+1);
108 return (8.*pow<4>(xi)*xipower +
109 (8.*(i-4)+4.)*xi*xi*xipower*(xi*xi-1.) +
110 (i-4)*(i-5)*xipower*(xi*xi-1.)*(xi*xi-1.))/denominator;
123 return 0.75 * (-1. + xi*xi);
125 return 0.75 * (1. - xi*xi);
127 return 0.25 * (-1. - 2.*xi + 3.*xi*xi);
129 return 0.25 * (-1. + 2.*xi + 3.*xi*xi);
131 return 4.*xi * (xi*xi-1.)/24.;
133 return (4*xi*xi*(xi*xi-1.) + (xi*xi-1.)*(xi*xi-1.))/120.;
137 Real denominator = 720., xipower = 1.;
138 for (
unsigned n=6; n != i; ++n)
141 denominator *= (n+1);
143 return (4*xi*xi*xi*xipower*(xi*xi-1.) +
144 (i-4)*xi*xipower*(xi*xi-1.)*(xi*xi-1.))/denominator;
154 return 0.25 * (2. - 3.*xi + xi*xi*xi);
156 return 0.25 * (2. + 3.*xi - xi*xi*xi);
158 return 0.25 * (1. - xi - xi*xi + xi*xi*xi);
160 return 0.25 * (-1. - xi + xi*xi + xi*xi*xi);
163 return (xi*xi-1.) * (xi*xi-1.)/24.;
165 return xi * (xi*xi-1.) * (xi*xi-1.)/120.;
169 Real denominator = 720., xipower = 1.;
170 for (
unsigned n=6; n != i; ++n)
173 denominator *= (n+1);
175 return (xi*xi*xipower*(xi*xi-1.)*(xi*xi-1.))/denominator;
183 const Order libmesh_dbg_var(order),
184 const unsigned int i,
186 const bool libmesh_dbg_var(add_p_level))
195 hermite_compute_coefs(elem, d1xd1x, d2xd2x);
200 const unsigned int totalorder =
201 order + add_p_level * elem->
p_level();
210 libmesh_assert_less (i, totalorder+1);
239 libmesh_error_msg(
"Hermite elements require the real element \nto construct gradient-based degrees of freedom.");
247 const unsigned int i,
249 const bool add_p_level)
257 const Order libmesh_dbg_var(order),
258 const unsigned int i,
261 const bool libmesh_dbg_var(add_p_level))
270 hermite_compute_coefs(elem, d1xd1x, d2xd2x);
275 const unsigned int totalorder =
276 order + add_p_level * elem->
p_level();
285 libmesh_assert_less (i, totalorder+1);
314 libmesh_error_msg(
"Hermite elements require the real element \nto construct gradient-based degrees of freedom.");
321 const unsigned int i,
322 const unsigned int j,
324 const bool add_p_level)
329 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 334 const Order libmesh_dbg_var(order),
335 const unsigned int i,
338 const bool libmesh_dbg_var(add_p_level))
347 hermite_compute_coefs(elem, d1xd1x, d2xd2x);
352 const unsigned int totalorder =
353 order + add_p_level * elem->
p_level();
362 libmesh_assert_less (i, totalorder+1);
391 libmesh_error_msg(
"Hermite elements require the real element \nto construct gradient-based degrees of freedom.");
398 const unsigned int i,
399 const unsigned int j,
401 const bool add_p_level)
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...
ElemType
Defines an enum for geometric element types.
Order
defines an enum for polynomial orders.
static OutputShape shape(const ElemType t, const Order o, const unsigned int i, const Point &p)
This is the base class from which all geometric element types are derived.
static OutputShape shape_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
unsigned int p_level() const
OrderWrapper order
The approximation order of the element.
The libMesh namespace provides an interface to certain functionality in the library.
static Real hermite_raw_shape_second_deriv(const unsigned int basis_num, const Real xi)
1D hermite functions on unit interval
static Real hermite_raw_shape(const unsigned int basis_num, const Real xi)
LIBMESH_DEFAULT_VECTORIZED_FE(template<>Real FE< 0, BERNSTEIN)
static unsigned int n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
static shape_deriv_ptr shape_deriv_function(const unsigned int dim, const FEType &fe_t, const ElemType t)
std::string enum_to_string(const T e)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static Real hermite_raw_shape_deriv(const unsigned int basis_num, const Real xi)
FEFamily
defines an enum for finite element families.
virtual Order default_order() const =0
virtual ElemType type() const =0
A Point defines a location in LIBMESH_DIM dimensional Real space.
const Point & point(const unsigned int i) const
static OutputShape shape_second_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
static FEFamily map_fe_type(const Elem &elem)