20 #include "libmesh/fe.h" 21 #include "libmesh/elem.h" 22 #include "libmesh/fe_interface.h" 23 #include "libmesh/enum_to_string.h" 35 static const Elem * old_elem_ptr =
nullptr;
40 static Real d1xd1x, d2xd2x;
42 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 44 Real clough_raw_shape_second_deriv(
const unsigned int basis_num,
45 const unsigned int deriv_type,
50 Real clough_raw_shape_deriv(
const unsigned int basis_num,
51 const unsigned int deriv_type,
53 Real clough_raw_shape(
const unsigned int basis_num,
58 void clough_compute_coefs(
const Elem * elem)
67 if (elem->
id() == old_elem_id &&
73 const FEType map_fe_type(elem->default_order(), mapping_family);
77 const int n_mapping_shape_functions =
81 std::vector<Point> dofpt;
82 dofpt.push_back(
Point(0));
83 dofpt.push_back(
Point(1));
86 std::vector<Real> dxdxi(2);
87 std::vector<Real> dxidx(2);
92 for (
int p = 0; p != 2; ++p)
94 for (
int i = 0; i != n_mapping_shape_functions; ++i)
96 const Real ddxi = shape_deriv_ptr
97 (map_fe_type, elem, i, 0, dofpt[p],
false);
98 dxdxi[p] += dofpt[p](0) * ddxi;
106 if (elem->id() == old_elem_id &&
107 elem == old_elem_ptr)
109 libmesh_assert_equal_to(d1xd1x, dxdxi[0]);
110 libmesh_assert_equal_to(d2xd2x, dxdxi[1]);
114 old_elem_id = elem->id();
122 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 125 Real clough_raw_shape_second_deriv(
const unsigned int basis_num,
126 const unsigned int deriv_type,
148 libmesh_error_msg(
"Invalid shape function index i = " <<
153 libmesh_error_msg(
"Invalid shape function derivative j = " <<
161 Real clough_raw_shape_deriv(
const unsigned int basis_num,
162 const unsigned int deriv_type,
173 return -6*xi + 6*xi*xi;
175 return 6*xi - 6*xi*xi;
177 return 1 - 4*xi + 3*xi*xi;
179 return -2*xi + 3*xi*xi;
182 libmesh_error_msg(
"Invalid shape function index i = " <<
187 libmesh_error_msg(
"Invalid shape function derivative j = " <<
192 Real clough_raw_shape(
const unsigned int basis_num,
200 return 1 - 3*xi*xi + 2*xi*xi*xi;
202 return 3*xi*xi - 2*xi*xi*xi;
204 return xi - 2*xi*xi + xi*xi*xi;
206 return -xi*xi + xi*xi*xi;
209 libmesh_error_msg(
"Invalid shape function index i = " <<
228 const unsigned int i,
230 const bool add_p_level)
234 clough_compute_coefs(elem);
238 const Order totalorder =
239 static_cast<Order>(order + add_p_level * elem->
p_level());
252 libmesh_assert_less (i, 4);
257 return clough_raw_shape(0, p);
259 return clough_raw_shape(1, p);
261 return d1xd1x * clough_raw_shape(2, p);
263 return d2xd2x * clough_raw_shape(3, p);
265 libmesh_error_msg(
"Invalid shape function index i = " << i);
274 libmesh_error_msg(
"ERROR: Unsupported polynomial order = " << totalorder);
286 libmesh_error_msg(
"Clough-Tocher elements require the real element \nto construct gradient-based degrees of freedom.");
293 const unsigned int i,
295 const bool add_p_level)
310 libmesh_error_msg(
"Clough-Tocher elements require the real element \nto construct gradient-based degrees of freedom.");
319 const unsigned int i,
320 const unsigned int j,
322 const bool add_p_level)
326 clough_compute_coefs(elem);
330 const Order totalorder =
331 static_cast<Order>(order + add_p_level * elem->
p_level());
347 return clough_raw_shape_deriv(0, j, p);
349 return clough_raw_shape_deriv(1, j, p);
351 return d1xd1x * clough_raw_shape_deriv(2, j, p);
353 return d2xd2x * clough_raw_shape_deriv(3, j, p);
355 libmesh_error_msg(
"Invalid shape function index i = " << i);
364 libmesh_error_msg(
"ERROR: Unsupported polynomial order = " << totalorder);
372 const unsigned int i,
373 const unsigned int j,
375 const bool add_p_level)
382 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 387 const unsigned int i,
388 const unsigned int j,
390 const bool add_p_level)
394 clough_compute_coefs(elem);
398 const Order totalorder =
399 static_cast<Order>(order + add_p_level * elem->
p_level());
415 return clough_raw_shape_second_deriv(0, j, p);
417 return clough_raw_shape_second_deriv(1, j, p);
419 return d1xd1x * clough_raw_shape_second_deriv(2, j, p);
421 return d2xd2x * clough_raw_shape_second_deriv(3, j, p);
423 libmesh_error_msg(
"Invalid shape function index i = " << i);
432 libmesh_error_msg(
"ERROR: Unsupported polynomial order = " << totalorder);
444 libmesh_error_msg(
"Clough-Tocher elements require the real element \nto construct gradient-based degrees of freedom.");
451 const unsigned int i,
452 const unsigned int j,
454 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.
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 const dof_id_type invalid_id
An invalid id to distinguish an uninitialized DofObject.
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
FEFamily
defines an enum for finite element families.
virtual ElemType type() const =0
A Point defines a location in LIBMESH_DIM dimensional Real space.
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)