4 #ifndef DUNE_LOCALFUNCTIONS_QKLOCALBASIS_HH 5 #define DUNE_LOCALFUNCTIONS_QKLOCALBASIS_HH 9 #include <dune/common/deprecated.hh> 10 #include <dune/common/fvector.hh> 11 #include <dune/common/fmatrix.hh> 12 #include <dune/common/power.hh> 14 #include <dune/geometry/type.hh> 34 template<
class D,
class R,
int k,
int d>
37 enum { n = StaticPower<k+1,d>::power };
40 static R p (
int i, D x)
43 for (
int j=0; j<=k; j++)
44 if (j!=i) result *= (k*x-j)/(i-j);
49 static R dp (
int i, D x)
53 for (
int j=0; j<=k; j++)
56 R prod( (k*1.0)/(i-j) );
57 for (
int l=0; l<=k; l++)
59 prod *= (k*x-l)/(i-l);
66 static Dune::FieldVector<int,d> multiindex (
int i)
68 Dune::FieldVector<int,d> alpha;
69 for (
int j=0; j<d; j++)
78 typedef LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,d>, 1>
Traits;
83 return StaticPower<k+1,d>::power;
88 std::vector<typename Traits::RangeType>& out)
const 91 for (
size_t i=0; i<
size(); i++)
94 Dune::FieldVector<int,d> alpha(multiindex(i));
100 for (
int j=0; j<d; j++)
101 out[i] *= p(alpha[j],in[j]);
111 std::vector<typename Traits::JacobianType>& out)
const 116 for (
size_t i=0; i<
size(); i++)
119 Dune::FieldVector<int,d> alpha(multiindex(i));
122 for (
int j=0; j<d; j++)
126 out[i][0][j] = dp(alpha[j],in[j]);
129 for (
int l=0; l<d; l++)
131 out[i][0][j] *= p(alpha[l],in[l]);
143 std::vector<typename Traits::RangeType>& out)
const 145 auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
157 for (
size_t i=0; i<
size(); i++)
160 Dune::FieldVector<int,d> alpha(multiindex(i));
166 for (std::size_t l=0; l<d; l++)
167 out[i][0] *= (order[l]) ? dp(alpha[l],in[l]) : p(alpha[l],in[l]);
172 DUNE_THROW(NotImplemented,
"Desired derivative order is not implemented");
181 template<
int diffOrder>
182 inline void DUNE_DEPRECATED_MSG(
"Use method 'partial' instead!")
184 const
std::array<
int,1>& direction,
185 const typename Traits::DomainType& in,
186 std::vector<typename Traits::RangeType>& out)
const 188 static_assert(diffOrder == 1,
"We only can compute first derivatives");
192 for (
size_t i=0; i<
size(); i++)
195 Dune::FieldVector<int,d> alpha(multiindex(i));
198 std::size_t j = direction[0];
202 out[i][0] = dp(alpha[j],in[j]);
205 for (std::size_t l=0; l<d; l++)
207 out[i][0] *= p(alpha[l],in[l]);
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: qklocalbasis.hh:87
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: qklocalbasis.hh:110
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:15
void partial(const std::array< unsigned int, d > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate partial derivatives of any order of all shape functions.
Definition: qklocalbasis.hh:141
void evaluate(const std::array< int, 1 > &direction, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate derivative in a given direction.
Definition: qklocalbasis.hh:183
D DomainType
domain type
Definition: localbasis.hh:49
unsigned int size() const
number of shape functions
Definition: qklocalbasis.hh:81
LocalBasisTraits< D, d, Dune::FieldVector< D, d >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, d >, 1 > Traits
Definition: qklocalbasis.hh:78
unsigned int order() const
Polynomial order of the shape functions.
Definition: qklocalbasis.hh:212
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
Lagrange shape functions of order k on the reference cube.
Definition: qklocalbasis.hh:35