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/localtoglobaladaptors.hh>
22 template<
int k,
int n>
55 template<
class D,
class R,
int k>
59 for (
int j=0; j<=k; j++)
60 if (j!=i) result *= (k*x-j)/(i-j);
65 template<
class D,
class R,
int k>
70 for (
int j=0; j<=k; j++)
73 R prod( (k*1.0)/(i-j) );
74 for (
int l=0; l<=k; l++)
75 if (l!=i && l!=j) prod *= (k*x-l)/(i-l);
82 template<
class D,
class R,
int k>
87 R
p (
int i, D
x)
const
90 for (
int j=0; j<=k; j++)
91 if (j!=i) result *= (k*
x-j)/(i-j);
96 R
dp (
int i, D
x)
const
100 for (
int j=0; j<=k; j++)
103 R prod( (k*1.0)/(i-j) );
104 for (
int l=0; l<=k; l++)
105 if (l!=i && l!=j) prod *= (k*
x-l)/(i-l);
114 return i/((1.0)*(k+1));
118 template<
int k,
int d>
121 Dune::FieldVector<int,d> alpha;
122 for (
int j=0; j<d; j++)
124 alpha[j] = i % (k+1);
142 template<
class D,
class R,
int k,
int d>
147 typedef LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,d> >
Traits;
157 std::vector<typename Traits::RangeType>& out)
const
160 for (
size_t i=0; i<
size(); i++)
163 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
169 for (
int j=0; j<d; j++)
170 out[i] *= p<D,R,k>(alpha[j],in[j]);
177 std::vector<typename Traits::JacobianType>& out)
const
182 for (
size_t i=0; i<
size(); i++)
185 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
188 for (
int j=0; j<d; j++)
192 out[i][0][j] = dp<D,R,k>(alpha[j],in[j]);
195 for (
int l=0; l<d; l++)
197 out[i][0][j] *= p<D,R,k>(alpha[l],in[l]);
203 void partial(
const std::array<unsigned int,Traits::dimDomain>& DUNE_UNUSED(
order),
204 const typename Traits::DomainType& DUNE_UNUSED(in),
205 std::vector<typename Traits::RangeType>& DUNE_UNUSED(out))
const
207 auto totalOrder = std::accumulate(
order.begin(),
order.end(), 0);
208 if (totalOrder == 0) {
211 DUNE_THROW(NotImplemented,
"Desired derivative order is not implemented");
228 template<
int k,
int d>
236 li[i] = LocalKey(0,0,i);
252 std::vector<LocalKey> li;
256 template<
int k,
int d,
class LB>
262 template<
typename F,
typename C>
265 typename LB::Traits::DomainType x;
272 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
275 for (
int j=0; j<d; j++)
276 x[j] = (1.0*alpha[j])/k;
284 template<
int d,
class LB>
289 template<
typename F,
typename C>
292 typename LB::Traits::DomainType x(0.5);
303 template<
class D,
class R,
int k,
int d>
316 typedef LocalFiniteElementTraits<LocalBasis,LocalCoefficients,LocalInterpolation>
Traits;
321 : gt(GeometryTypes::cube(d))
342 return interpolation;
366 LocalCoefficients coefficients;
367 LocalInterpolation interpolation;
377 template<
class Geometry,
class RF,
int k>
379 public ScalarLocalToGlobalFiniteElementAdaptorFactory<
380 QkDGLagrangeLocalFiniteElement<
381 typename Geometry::ctype, RF, k, Geometry::mydimension
387 typename Geometry::ctype, RF, k, Geometry::mydimension
389 typedef ScalarLocalToGlobalFiniteElementAdaptorFactory<LFE, Geometry> Base;
391 static const LFE lfe;
398 template<
class Geometry,
class RF,
int k>
399 const typename QkDGFiniteElementFactory<Geometry, RF, k>::LFE
400 QkDGFiniteElementFactory<Geometry, RF, k>::lfe;
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::FieldVector< int, d > multiindex(int i)
Definition: qkdglagrange.hh:119
R p(int i, D x)
Definition: qkdglagrange.hh:56
R dp(int i, D x)
Definition: qkdglagrange.hh:66
Definition: qkdglagrange.hh:24
@ value
Definition: qkdglagrange.hh:26
Lagrange polynomials of degree k at equidistant points as a class.
Definition: qkdglagrange.hh:84
R x(int i)
Definition: qkdglagrange.hh:112
R dp(int i, D x) const
Definition: qkdglagrange.hh:96
R p(int i, D x) const
Definition: qkdglagrange.hh:87
Lagrange shape functions of order k on the reference cube.
Definition: qkdglagrange.hh:144
unsigned int size() const
number of shape functions
Definition: qkdglagrange.hh:150
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: qkdglagrange.hh:176
void partial(const std::array< unsigned int, Traits::dimDomain > &DUNE_UNUSED(order), const typename Traits::DomainType &DUNE_UNUSED(in), std::vector< typename Traits::RangeType > &DUNE_UNUSED(out)) const
Evaluate partial derivative of all shape functions.
Definition: qkdglagrange.hh:203
LocalBasisTraits< D, d, Dune::FieldVector< D, d >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, d > > Traits
Definition: qkdglagrange.hh:147
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: qkdglagrange.hh:156
unsigned int order() const
Polynomial order of the shape functions.
Definition: qkdglagrange.hh:216
Layout map for Q1 elements.
Definition: qkdglagrange.hh:230
const LocalKey & localKey(std::size_t i) const
get i'th index
Definition: qkdglagrange.hh:246
QkDGLocalCoefficients()
Standard constructor.
Definition: qkdglagrange.hh:233
std::size_t size() const
number of coefficients
Definition: qkdglagrange.hh:240
Definition: qkdglagrange.hh:258
void interpolate(const F &f, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdglagrange.hh:263
void interpolate(const F &f, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdglagrange.hh:290
Definition: qkdglagrange.hh:305
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: qkdglagrange.hh:333
const Traits::LocalBasisType & localBasis() const
Definition: qkdglagrange.hh:326
const Traits::LocalInterpolationType & localInterpolation() const
Definition: qkdglagrange.hh:340
std::size_t size() const
Definition: qkdglagrange.hh:347
QkDGLagrangeLocalFiniteElement * clone() const
Definition: qkdglagrange.hh:359
LocalFiniteElementTraits< LocalBasis, LocalCoefficients, LocalInterpolation > Traits
Definition: qkdglagrange.hh:316
GeometryType type() const
Definition: qkdglagrange.hh:354
QkDGLagrangeLocalFiniteElement()
Definition: qkdglagrange.hh:320
@ n
Definition: qkdglagrange.hh:312
Factory for global-valued QkDG elements.
Definition: qkdglagrange.hh:385
QkDGFiniteElementFactory()
default constructor
Definition: qkdglagrange.hh:395
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139