4#ifndef DUNE_PDELAB_FINITEELEMENT_QKDGLAGRANGE_HH
5#define DUNE_PDELAB_FINITEELEMENT_QKDGLAGRANGE_HH
9#include <dune/localfunctions/common/localbasis.hh>
10#include <dune/localfunctions/common/localkey.hh>
11#include <dune/localfunctions/common/localfiniteelementtraits.hh>
12#include <dune/localfunctions/common/localinterpolation.hh>
13#include <dune/localfunctions/common/localtoglobaladaptors.hh>
14#include <dune/localfunctions/common/localinterpolation.hh>
24 template<
int k,
int n>
28 value=(k+1)*QkSize<k,n-1>::value
57 template<
class D,
class R,
int k>
61 for (
int j=0; j<=k; j++)
62 if (j!=i) result *= (k*x-j)/(i-j);
67 template<
class D,
class R,
int k>
72 for (
int j=0; j<=k; j++)
75 R prod( (k*1.0)/(i-j) );
76 for (
int l=0; l<=k; l++)
77 if (l!=i && l!=j) prod *= (k*x-l)/(i-l);
84 template<
class D,
class R,
int k>
89 R p (
int i, D x)
const
92 for (
int j=0; j<=k; j++)
93 if (j!=i) result *= (k*x-j)/(i-j);
98 R dp (
int i, D x)
const
102 for (
int j=0; j<=k; j++)
105 R prod( (k*1.0)/(i-j) );
106 for (
int l=0; l<=k; l++)
107 if (l!=i && l!=j) prod *= (k*x-l)/(i-l);
116 return i/((1.0)*(k+1));
120 template<
int k,
int d>
124 for (
int j=0; j<d; j++)
126 alpha[j] = i % (k+1);
144 template<
class D,
class R,
int k,
int d>
147 enum{ n = QkSize<k,d>::value };
149 typedef LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,
Dune::FieldVector<R,1>,
Dune::FieldMatrix<R,1,d> >
Traits;
154 return QkSize<k,d>::value;
159 std::vector<typename Traits::RangeType>& out)
const
162 for (
size_t i=0; i<
size(); i++)
171 for (
int j=0; j<d; j++)
172 out[i] *= p<D,R,k>(alpha[j],in[j]);
179 std::vector<typename Traits::JacobianType>& out)
const
184 for (
size_t i=0; i<
size(); i++)
190 for (
int j=0; j<d; j++)
194 out[i][0][j] = dp<D,R,k>(alpha[j],in[j]);
197 for (
int l=0; l<d; l++)
199 out[i][0][j] *= p<D,R,k>(alpha[l],in[l]);
207 std::vector<typename Traits::RangeType>& out)
const
210 if (totalOrder == 0) {
230 template<
int k,
int d>
237 for (std::size_t i=0; i<QkSize<k,d>::value; i++)
244 return QkSize<k,d>::value;
254 std::vector<LocalKey> li;
258 template<
int k,
int d,
class LB>
264 template<
typename F,
typename C>
267 typename LB::Traits::DomainType x;
268 auto&& f = Impl::makeFunctionWithCallOperator<typename LB::Traits::DomainType>(ff);
270 out.resize(QkSize<k,d>::value);
272 for (
int i=0; i<QkSize<k,d>::value; i++)
278 for (
int j=0; j<d; j++)
279 x[j] = (1.0*alpha[j])/k;
287 template<
int d,
class LB>
292 template<
typename F,
typename C>
295 typename LB::Traits::DomainType x(0.5);
296 auto&& f = Impl::makeFunctionWithCallOperator<typename LB::Traits::DomainType>(ff);
307 template<
class D,
class R,
int k,
int d>
316 enum{ n = QkStuff::QkSize<k,d>::value };
325 : gt(GeometryTypes::
cube(d))
346 return interpolation;
370 LocalCoefficients coefficients;
371 LocalInterpolation interpolation;
381 template<
class Geometry,
class RF,
int k>
384 QkDGLagrangeLocalFiniteElement<
385 typename Geometry::ctype, RF, k, Geometry::mydimension
395 static const LFE lfe;
402 template<
class Geometry,
class RF,
int k>
403 const typename QkDGFiniteElementFactory<Geometry, RF, k>::LFE
404 QkDGFiniteElementFactory<Geometry, RF, k>::lfe;
A dense n x m matrix.
Definition: fmatrix.hh:117
vector space out of a tensor product of fields.
Definition: fvector.hh:95
Unique label for each type of entities that can occur in DUNE grids.
Definition: type.hh:114
GridImp::ctype ctype
define type used for coordinates in grid module
Definition: geometry.hh:100
static constexpr int mydimension
geometry dimension
Definition: geometry.hh:94
Describe position of one degree of freedom.
Definition: localkey.hh:24
Default exception for dummy implementations.
Definition: exceptions.hh:263
Factory for global-valued QkDG elements.
Definition: qkdglagrange.hh:389
QkDGFiniteElementFactory()
default constructor
Definition: qkdglagrange.hh:399
Definition: qkdglagrange.hh:309
std::size_t size() const
Definition: qkdglagrange.hh:351
const Traits::LocalInterpolationType & localInterpolation() const
Definition: qkdglagrange.hh:344
LocalFiniteElementTraits< LocalBasis, LocalCoefficients, LocalInterpolation > Traits
Definition: qkdglagrange.hh:320
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: qkdglagrange.hh:337
GeometryType type() const
Definition: qkdglagrange.hh:358
QkDGLagrangeLocalFiniteElement()
Definition: qkdglagrange.hh:324
const Traits::LocalBasisType & localBasis() const
Definition: qkdglagrange.hh:330
Lagrange polynomials of degree k at equidistant points as a class.
Definition: qkdglagrange.hh:86
Layout map for Q1 elements.
Definition: qkdglagrange.hh:232
const LocalKey & localKey(std::size_t i) const
get i'th index
Definition: qkdglagrange.hh:248
QkDGLocalCoefficients()
Standard constructor.
Definition: qkdglagrange.hh:235
std::size_t size() const
number of coefficients
Definition: qkdglagrange.hh:242
Lagrange shape functions of order k on the reference cube.
Definition: qkdglagrange.hh:146
unsigned int size() const
number of shape functions
Definition: qkdglagrange.hh:152
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: qkdglagrange.hh:178
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: qkdglagrange.hh:158
void partial(const std::array< unsigned int, Traits::dimDomain > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate partial derivative of all shape functions.
Definition: qkdglagrange.hh:205
unsigned int order() const
Polynomial order of the shape functions.
Definition: qkdglagrange.hh:218
void interpolate(const F &ff, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdglagrange.hh:293
Definition: qkdglagrange.hh:260
void interpolate(const F &ff, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdglagrange.hh:265
Factory for ScalarLocalToGlobalFiniteElementAdaptor objects.
Definition: localtoglobaladaptors.hh:244
#define DUNE_THROW(E, m)
Definition: exceptions.hh:218
constexpr GeometryType cube(unsigned int dim)
Returns a GeometryType representing a hypercube of dimension dim.
Definition: type.hh:462
constexpr T accumulate(Range &&range, T value, F &&f)
Accumulate values.
Definition: hybridutilities.hh:279
Dune namespace.
Definition: alignedallocator.hh:13
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:35
D DomainType
domain type
Definition: localbasis.hh:43
traits helper struct
Definition: localfiniteelementtraits.hh:13
LB LocalBasisType
Definition: localfiniteelementtraits.hh:16
LC LocalCoefficientsType
Definition: localfiniteelementtraits.hh:20
LI LocalInterpolationType
Definition: localfiniteelementtraits.hh:24