3 #ifndef DUNE_PRISM_P2_LOCALBASIS_HH 4 #define DUNE_PRISM_P2_LOCALBASIS_HH 8 #include <dune/common/fmatrix.hh> 24 template<
class D,
class R>
40 std::vector<typename Traits::RangeType>& out)
const 46 R a[2], b[2], c[2], a1d, b1d, c1d;
60 out[0] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
73 out[1] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
86 out[2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
99 out[3] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
112 out[4] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
125 out[5] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
138 out[6] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
151 out[7] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
164 out[8] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
177 out[9] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
190 out[10] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
203 out[11] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
216 out[12] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
229 out[13] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
242 out[14] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
255 out[15] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
268 out[16] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
281 out[17] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
287 std::vector<typename Traits::JacobianType>& out)
const 294 R a[2], b[2], c[2], aa[2], bb[2][2], a1d, b1d, c1d;
314 out[0][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
315 out[0][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
316 out[0][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
336 out[1][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
337 out[1][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
338 out[1][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
357 out[2][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
358 out[2][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
359 out[2][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
379 out[3][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
380 out[3][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
381 out[3][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
401 out[4][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
402 out[4][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
403 out[4][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
423 out[5][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
424 out[5][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
425 out[5][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
446 out[6][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
447 out[6][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
448 out[6][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
468 out[7][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
469 out[7][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
470 out[7][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
490 out[8][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
491 out[8][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
492 out[8][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
514 out[9][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
515 out[9][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
516 out[9][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
537 out[10][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
538 out[10][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
539 out[10][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
559 out[11][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
560 out[11][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
561 out[11][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
582 out[12][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
583 out[12][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
584 out[12][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
604 out[13][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
605 out[13][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
606 out[13][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
626 out[14][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
627 out[14][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
628 out[14][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
649 out[15][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
650 out[15][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
651 out[15][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
673 out[16][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
674 out[16][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
675 out[16][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
696 out[17][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
697 out[17][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
698 out[17][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
705 std::vector<typename Traits::RangeType>& out)
const 707 auto totalOrder = std::accumulate(
order.begin(),
order.end(), 0);
708 if (totalOrder == 0) {
710 }
else if (totalOrder == 1) {
713 auto const direction = std::distance(
order.begin(), std::find(
order.begin(),
order.end(), 1));
716 out[0] = (-3 + 4*(in[0] + in[1])) * (1 - 3*in[2] + 2*in[2]*in[2]);
717 out[1] = (-1 + 4*in[0]) * (1 - 3*in[2] + 2*in[2]*in[2]);
719 out[3] = (-3 + 4*(in[0] + in[1])) * (-in[2] + 2*in[2]*in[2]);
720 out[4] = (-1 + 4*in[0]) * (-in[2] + 2*in[2]*in[2]);
722 out[6] = (-3 + 4*(in[0] + in[1])) * 4*(in[2] - in[2]*in[2]);
723 out[7] = (-1 + 4*in[0]) * 4*(in[2] - in[2]*in[2]);
725 out[9] = (4 - 8*in[0] - 4*in[1]) * (1 - 3*in[2] + 2*in[2]*in[2]);
726 out[10] = -4*in[1] * (1 - 3*in[2] + 2*in[2]*in[2]);
727 out[11] = 4*in[1] * (1 - 3*in[2] + 2*in[2]*in[2]);
728 out[12] = (4 - 8*in[0] - 4*in[1]) * (-in[2] + 2*in[2]*in[2]);
729 out[13] = -4*in[1] * (-in[2] + 2*in[2]*in[2]);
730 out[14] = 4*in[1] * (-in[2] + 2*in[2]*in[2]);
731 out[15] = (4 - 8*in[0] - 4*in[1]) * (4*in[2] - 4*in[2]*in[2]);
732 out[16] = -4*in[1] * (4*in[2] - 4*in[2]*in[2]);
733 out[17] = 4*in[1] * (4*in[2] - 4*in[2]*in[2]);
737 out[0] = (-3 + 4*(in[0] + in[1])) * (1 - 3*in[2] + 2*in[2]*in[2]);
739 out[2] = (-1 + 4*in[1]) * (1 - 3*in[2] + 2*in[2]*in[2]);
740 out[3] = (-3 + 4*(in[0] + in[1])) * (-in[2] + 2*in[2]*in[2]);
742 out[5] = (-1 + 4*in[1]) * (-in[2] + 2*in[2]*in[2]);
743 out[6] = (-3 + 4*(in[0] + in[1])) * 4*(in[2] - in[2]*in[2]);
745 out[8] = (-1 + 4*in[1]) * 4*(in[2] - in[2]*in[2]);
746 out[9] = -4*in[0] * (1 - 3*in[2] + 2*in[2]*in[2]);
747 out[10] = (4 - 4*in[0] - 8*in[1]) * (1 - 3*in[2] + 2*in[2]*in[2]);
748 out[11] = 4*in[0] * (1 - 3*in[2] + 2*in[2]*in[2]);
749 out[12] = -4*in[0] * (-in[2] + 2*in[2]*in[2]);
750 out[13] = (4 - 4*in[0] - 8*in[1]) * (-in[2] + 2*in[2]*in[2]);
751 out[14] = 4*in[0] * (-in[2] + 2*in[2]*in[2]);
752 out[15] = -4*in[0] * 4*(in[2] - in[2]*in[2]);
753 out[16] = (4 - 4*in[0] - 8*in[1]) * (4*in[2] - 4*in[2]*in[2]);
754 out[17] = 4*in[0] * (4*in[2] - 4*in[2]*in[2]);
758 out[0] = 2 * (1 - in[0] - in[1]) * (0.5 - in[0] - in[1]) * (-3 + 4*in[2]);
759 out[1] = 2 * in[0] * (-0.5 + in[0]) * (-3 + 4*in[2]);
760 out[2] = 2 * in[1] * (-0.5 + in[1]) * (-3 + 4*in[2]);
761 out[3] = 2 * (1 - in[0] - in[1]) * (0.5 - in[0] - in[1]) * (-1 + 4*in[2]);
762 out[4] = 2 * in[0] * (-0.5 + in[0]) * (-1 + 4*in[2]);
763 out[5] = 2 * in[1] * (-0.5 + in[1]) * (-1 + 4*in[2]);
764 out[6] = 2 * (1 - in[0] - in[1]) * (0.5 - in[0] - in[1]) * (4 - 8*in[2]);
765 out[7] = 2 * in[0] * (-0.5 + in[0]) * (4 - 8*in[2]);
766 out[8] = 2*in[1] * (-0.5 + in[1]) * (4 - 8*in[2]);
767 out[9] = 4*in[0] * (1 - in[0] - in[1]) * (-3 + 4*in[2]);
768 out[10] = 4*in[1] * (1 - in[0] - in[1]) * (-3 + 4*in[2]);
769 out[11] = 4*in[0]*in[1] * (-3 + 4*in[2]);
770 out[12] = 4 * in[0] * (1 - in[0] - in[1]) * (-1 + 4*in[2]);
771 out[13] = 4 * in[1] * (1 - in[0] - in[1]) * (-1 + 4*in[2]);
772 out[14] = 4*in[0]*in[1] * (-1 + 4*in[2]);
773 out[15] = 4 * in[0] * (1 - in[0] - in[1]) * (4 - 8*in[2]);
774 out[16] = 4 * in[1] * (1 - in[0] - in[1]) * (4 - 8*in[2]);
775 out[17] = 4*in[0]*in[1] * (4 - 8*in[2]);
779 DUNE_THROW(RangeError,
"Component out of range.");
782 DUNE_THROW(NotImplemented,
"Desired derivative order is not implemented");
LocalBasisTraits< D, 3, Dune::FieldVector< D, 3 >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, 3 > > Traits
export type traits for function signature
Definition: prismp2localbasis.hh:30
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:31
unsigned int size() const
number of shape functions
Definition: prismp2localbasis.hh:33
void partial(const std::array< unsigned int, 3 > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate partial derivatives of all shape functions.
Definition: prismp2localbasis.hh:703
D DomainType
domain type
Definition: localbasis.hh:43
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: prismp2localbasis.hh:286
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:15
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: prismp2localbasis.hh:39
Quadratic Lagrange shape functions on the prism.
Definition: prismp2localbasis.hh:25
unsigned int order() const
Polynomial order of the shape functions.
Definition: prismp2localbasis.hh:787