3 #ifndef DUNE_PK3DLOCALBASIS_HH 4 #define DUNE_PK3DLOCALBASIS_HH 8 #include <dune/common/fmatrix.hh> 26 template<
class D,
class R,
unsigned int k>
30 enum {
N = (k+1)*(k+2)*(k+3)/6};
47 std::vector<typename Traits::RangeType>& out)
const 55 for (i[2] = 0; i[2] <= k; ++i[2])
58 for (
unsigned int j = 0; j < i[2]; ++j)
59 factor[2] *= (kx[2]-j) / (i[2]-j);
60 for (i[1] = 0; i[1] <= k - i[2]; ++i[1])
63 for (
unsigned int j = 0; j < i[1]; ++j)
64 factor[1] *= (kx[1]-j) / (i[1]-j);
65 for (i[0] = 0; i[0] <= k - i[1] - i[2]; ++i[0])
68 for (
unsigned int j = 0; j < i[0]; ++j)
69 factor[0] *= (kx[0]-j) / (i[0]-j);
70 i[3] = k - i[0] - i[1] - i[2];
71 D kx3 = k - kx[0] - kx[1] - kx[2];
73 for (
unsigned int j = 0; j < i[3]; ++j)
74 factor[3] *= (kx3-j) / (i[3]-j);
75 out[n++] = factor[0] * factor[1] * factor[2] * factor[3];
84 std::vector<typename Traits::JacobianType>& out)
const 92 for (i[2] = 0; i[2] <= k; ++i[2])
95 for (
unsigned int j = 0; j < i[2]; ++j)
96 factor[2] *= (kx[2]-j) / (i[2]-j);
97 for (i[1] = 0; i[1] <= k - i[2]; ++i[1])
100 for (
unsigned int j = 0; j < i[1]; ++j)
101 factor[1] *= (kx[1]-j) / (i[1]-j);
102 for (i[0] = 0; i[0] <= k - i[1] - i[2]; ++i[0])
105 for (
unsigned int j = 0; j < i[0]; ++j)
106 factor[0] *= (kx[0]-j) / (i[0]-j);
107 i[3] = k - i[0] - i[1] - i[2];
108 D kx3 = k - kx[0] - kx[1] - kx[2];
111 for (
unsigned int j = 0; j < i[3]; ++j)
112 factor[3] /= i[3] - j;
113 R prod_all = factor[0] * factor[1] * factor[2] * factor[3];
114 for (
unsigned int j = 0; j < i[3]; ++j)
117 for (
unsigned int l = 0; l < i[3]; ++l)
124 for (
unsigned int j = 0; j < i[3]; ++j)
125 factor[3] *= kx3 - j;
126 for (
unsigned int m = 0; m < 3; ++m)
129 for (
unsigned int j = 0; j < i[m]; ++j)
132 for (
unsigned int p = 0; p < 3; ++p)
135 for (
unsigned int l = 0; l < i[p]; ++l)
137 prod *= R(k) / (i[p]-l);
139 prod *= (kx[p]-l) / (i[p]-l);
143 out[n][0][m] += prod;
159 std::vector<typename Traits::RangeType>& out)
const 161 auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
162 if (totalOrder == 0) {
165 DUNE_THROW(NotImplemented,
"Desired derivative order is not implemented");
178 template<
class D,
class R>
197 std::vector<typename Traits::RangeType>& out)
const 206 std::vector<typename Traits::JacobianType>& out)
const 221 std::vector<typename Traits::RangeType>& out)
const 223 auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
224 if (totalOrder == 0) {
233 template<
typename E,
typename F,
typename C>
234 void interpolate (
const E& e,
const F& f, std::vector<C>& out)
const R RangeType
range type
Definition: localbasis.hh:55
Lagrange shape functions of arbitrary order on the reference tetrahedron.
Definition: pk3dlocalbasis.hh:27
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 any order of all shape functions.
Definition: pk3dlocalbasis.hh:157
Pk3DLocalBasis()
Standard constructor.
Definition: pk3dlocalbasis.hh:37
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Definition: pk3dlocalbasis.hh:196
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:31
D DomainType
domain type
Definition: localbasis.hh:43
unsigned int order() const
Polynomial order of the shape functions.
Definition: pk3dlocalbasis.hh:170
LocalBasisTraits< D, 3, Dune::FieldVector< D, 3 >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, 3 > > Traits
Definition: pk3dlocalbasis.hh:183
void evaluateJacobian(const typename Traits::DomainType &x, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: pk3dlocalbasis.hh:83
Definition: pk3dlocalbasis.hh:30
void evaluateFunction(const typename Traits::DomainType &x, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: pk3dlocalbasis.hh:46
void interpolate(const E &e, const F &f, std::vector< C > &out) const
Definition: pk3dlocalbasis.hh:234
unsigned int order() const
Definition: pk3dlocalbasis.hh:245
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:15
unsigned int size() const
Definition: pk3dlocalbasis.hh:191
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 any order of all shape functions.
Definition: pk3dlocalbasis.hh:219
unsigned int size() const
number of shape functions
Definition: pk3dlocalbasis.hh:40
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Definition: pk3dlocalbasis.hh:205
LocalBasisTraits< D, 3, Dune::FieldVector< D, 3 >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, 3 > > Traits
Definition: pk3dlocalbasis.hh:34
Definition: pk3dlocalbasis.hh:31