21 #include "libmesh/libmesh_config.h" 22 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES 27 #include "libmesh/libmesh_common.h" 28 #include "libmesh/fe.h" 29 #include "libmesh/elem.h" 30 #include "libmesh/utility.h" 60 libmesh_error_msg(
"Invalid shape function index i = " << i);
68 return (1./4.)*pow<2>(1.-xi);
70 return (1./4.)*pow<2>(1.+xi);
72 return (1./2.)*(1.-xi)*(1.+xi);
74 libmesh_error_msg(
"Invalid shape function index i = " << i);
82 return (1./8.)*pow<3>(1.-xi);
84 return (1./8.)*pow<3>(1.+xi);
86 return (3./8.)*(1.+xi)*pow<2>(1.-xi);
88 return (3./8.)*pow<2>(1.+xi)*(1.-xi);
90 libmesh_error_msg(
"Invalid shape function index i = " << i);
98 return (1./16.)*pow<4>(1.-xi);
100 return (1./16.)*pow<4>(1.+xi);
102 return (1./ 4.)*(1.+xi)*pow<3>(1.-xi);
104 return (3./ 8.)*pow<2>(1.+xi)*pow<2>(1.-xi);
106 return (1./ 4.)*pow<3>(1.+xi)*(1.-xi);
108 libmesh_error_msg(
"Invalid shape function index i = " << i);
117 return (1./32.)*pow<5>(1.-xi);
119 return (1./32.)*pow<5>(1.+xi);
121 return (5./32.)*(1.+xi)*pow<4>(1.-xi);
123 return (5./16.)*pow<2>(1.+xi)*pow<3>(1.-xi);
125 return (5./16.)*pow<3>(1.+xi)*pow<2>(1.-xi);
127 return (5./32.)*pow<4>(1.+xi)*(1.-xi);
129 libmesh_error_msg(
"Invalid shape function index i = " << i);
138 return ( 1./64.)*pow<6>(1.-xi);
140 return ( 1./64.)*pow<6>(1.+xi);
142 return ( 3./32.)*(1.+xi)*pow<5>(1.-xi);
144 return (15./64.)*pow<2>(1.+xi)*pow<4>(1.-xi);
146 return ( 5./16.)*pow<3>(1.+xi)*pow<3>(1.-xi);
148 return (15./64.)*pow<4>(1.+xi)*pow<2>(1.-xi);
150 return ( 3./32.)*pow<5>(1.+xi)*(1.-xi);
152 libmesh_error_msg(
"Invalid shape function index i = " << i);
161 const int p_order =
static_cast<int>(order);
162 const int m = p_order-i+1;
165 Real binomial_p_i = 1;
172 static_cast<unsigned long>(n)));
177 return binomial_p_i *
std::pow((1-xi)/2, p_order);
179 return binomial_p_i *
std::pow((1+xi)/2, p_order);
182 return binomial_p_i *
std::pow((1+xi)/2,n)
195 const unsigned int i,
197 const bool add_p_level)
203 static_cast<Order>(order + add_p_level*elem->
p_level()), i, p);
210 const unsigned int i,
212 const bool add_p_level)
224 const unsigned int i,
225 const unsigned int libmesh_dbg_var(j),
230 libmesh_assert_equal_to (j, 0);
232 const Real xi = p(0);
247 libmesh_error_msg(
"Invalid shape function index i = " << i);
261 libmesh_error_msg(
"Invalid shape function index i = " << i);
269 return -0.375*pow<2>(1.-xi);
271 return 0.375*pow<2>(1.+xi);
273 return -0.375 -.75*xi +1.125*pow<2>(xi);
275 return 0.375 -.75*xi -1.125*pow<2>(xi);
277 libmesh_error_msg(
"Invalid shape function index i = " << i);
285 return -0.25*pow<3>(1.-xi);
287 return 0.25*pow<3>(1.+xi);
289 return -0.5 +1.5*pow<2>(xi)-pow<3>(xi);
291 return 1.5*(pow<3>(xi)-xi);
293 return 0.5 -1.5*pow<2>(xi)-pow<3>(xi);
295 libmesh_error_msg(
"Invalid shape function index i = " << i);
303 return -(5./32.)*pow<4>(xi-1.);
305 return (5./32.)*pow<4>(xi+1.);
307 return (5./32.)*pow<4>(1.-xi) -(5./8.)*(1.+xi)*pow<3>(1.-xi);
309 return (5./ 8.)*(1.+xi)*pow<3>(1.-xi) -(15./16.)*pow<2>(1.+xi)*pow<2>(1.-xi);
311 return -(5./ 8.)*pow<3>(1.+xi)*(1.-xi) +(15./16.)*pow<2>(1.+xi)*pow<2>(1.-xi);
313 return (5./ 8.)*pow<3>(1.+xi)*(1.-xi) -(5./32.)*pow<4>(1.+xi);
315 libmesh_error_msg(
"Invalid shape function index i = " << i);
323 return -( 3./32.)*pow<5>(1.-xi);
325 return ( 3./32.)*pow<5>(1.+xi);
327 return ( 3./32.)*pow<5>(1.-xi)-(15./32.)*(1.+xi)*pow<4>(1.-xi);
329 return (15./32.)*(1.+xi)*pow<4>(1.-xi)-(15./16.)*pow<2>(1.+xi)*pow<3>(1.-xi);
331 return -(15./ 8.)*xi +(15./4.)*pow<3>(xi)-(15./8.)*pow<5>(xi);
333 return -(15./32.)*(1.-xi)*pow<4>(1.+xi)+(15./16.)*pow<2>(1.-xi)*pow<3>(1.+xi);
335 return (15./32.)*pow<4>(1.+xi)*(1.-xi)-(3./32.)*pow<5>(1.+xi);
337 libmesh_error_msg(
"Invalid shape function index i = " << i);
346 const int p_order =
static_cast<int>(order);
347 const int m = p_order-(i-1);
350 Real binomial_p_i = 1;
357 static_cast<unsigned long>(n)));
362 return binomial_p_i * (-1./2.) * p_order *
std::pow((1-xi)/2, p_order-1);
364 return binomial_p_i * ( 1./2.) * p_order *
std::pow((1+xi)/2, p_order-1);
368 return binomial_p_i * (1./2. * n *
std::pow((1+xi)/2,n-1) *
std::pow((1-xi)/2,m)
382 const unsigned int i,
383 const unsigned int j,
385 const bool add_p_level)
391 static_cast<Order>(order + add_p_level*elem->
p_level()), i, j, p);
397 const unsigned int i,
398 const unsigned int j,
400 const bool add_p_level)
409 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES 414 const unsigned int i,
415 const unsigned int libmesh_dbg_var(j),
420 libmesh_assert_equal_to (j, 0);
422 const Real xi = p(0);
436 libmesh_error_msg(
"Invalid shape function index i = " << i);
449 libmesh_error_msg(
"Invalid shape function index i = " << i);
461 return -.75 + 2.25*xi;
463 return -.75 - 2.25*xi;
465 libmesh_error_msg(
"Invalid shape function index i = " << i);
473 return 0.75*pow<2>(1.-xi);
475 return 0.75*pow<2>(1.+xi);
477 return 3*(xi - pow<2>(xi));
479 return 1.5*(3*pow<2>(xi)-1);
481 return -3*xi-3*pow<2>(xi);
483 libmesh_error_msg(
"Invalid shape function index i = " << i);
491 return -(5./8.)*pow<3>(xi-1.);
493 return (5./8.)*pow<3>(xi+1.);
495 return -(5./4.)*pow<3>(1.-xi) + (15./8.)*(1.+xi)*pow<2>(1.-xi);
497 return -(15./ 4.)*(1.+xi)*pow<2>(1.-xi) + (5./ 8.)*pow<3>(1.-xi)
498 + (15./8.)*pow<2>(1.+xi)*(1.-xi);
500 return (5./ 8.)*pow<3>(1.+xi) - (15./ 4.)*pow<2>(1.+xi)*(1.-xi)
501 +(15./8.)*(1.+xi)*pow<2>(1.-xi);
503 return -(5./ 8.)*pow<3>(1.+xi) + (15./ 8.)*pow<2>(1.+xi)*(1.-xi)
504 -(5./8.)*pow<3>(1.+xi);
506 libmesh_error_msg(
"Invalid shape function index i = " << i);
514 return ( 15./32.)*pow<4>(1.-xi);
516 return ( 15./32.)*pow<4>(1.+xi);
518 return -( 15./8.)*pow<4>(1.-xi) +
519 ( 15./8.)*(1.+xi)*pow<3>(1.-xi);
521 return -(15./4.)*(1.+xi)*pow<3>(1.-xi)
522 + (15./32.)*pow<4>(1.-xi)
523 + (45./16.)*pow<2>(1.+xi)*pow<2>(1.-xi);
525 return -(15./ 8.) +(45./4.)*pow<2>(xi) - (75./8.)*pow<4>(xi);
527 return -(15./4.)*(1.-xi)*pow<3>(1.+xi)
528 + (15./32.)*pow<4>(1.+xi)
529 + (45./16.)*pow<2>(1.-xi)*pow<2>(1.+xi);
531 return -(15./16.)*pow<4>(1.+xi)
532 + (15./8.)*pow<3>(1.+xi)*(1.-xi);
534 libmesh_error_msg(
"Invalid shape function index i = " << i);
543 const int p_order =
static_cast<int>(order);
544 const int m = p_order-(i-1);
547 Real binomial_p_i = 1;
554 static_cast<unsigned long>(n)));
559 return binomial_p_i * (1./4.) * p_order * (p_order-1) *
std::pow((1-xi)/2, p_order-2);
561 return binomial_p_i * (1./4.) * p_order * (p_order-1) *
std::pow((1+xi)/2, p_order-2);
569 binomial_p_i * (-1./4. * m *
std::pow((1-xi)/2,m-1));
572 binomial_p_i * (-1./4. * n * m *
std::pow((1+xi)/2,n-1) *
std::pow((1-xi)/2,m-1) +
576 val += binomial_p_i * (-1./4. * n *
std::pow((1+xi)/2,n-1));
579 binomial_p_i * (1./4. * m * (m-1) *
std::pow((1+xi)/2,n) *
std::pow((1-xi)/2,m-2)
595 const unsigned int i,
596 const unsigned int j,
598 const bool add_p_level)
604 static_cast<Order>(order + add_p_level*elem->
p_level()), i, j, p);
610 const unsigned int i,
611 const unsigned int j,
613 const bool add_p_level)
626 #endif //LIBMESH_ENABLE_HIGHER_ORDER_SHAPES 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)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
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)