4 #ifndef DUNE_QkDG_LOCALFINITEELEMENT_HH 5 #define DUNE_QkDG_LOCALFINITEELEMENT_HH 7 #include <dune/common/fvector.hh> 8 #include <dune/common/deprecated.hh> 10 #include <dune/geometry/type.hh> 12 #include <dune/localfunctions/common/localbasis.hh> 13 #include <dune/localfunctions/common/localfiniteelementtraits.hh> 14 #include <dune/localfunctions/common/localkey.hh> 15 #include <dune/localfunctions/common/localtoglobaladaptors.hh> 28 template<
int k,
int n>
61 template<
class D,
class R,
int k>
65 for (
int j=0; j<=k; j++)
66 if (j!=i) result *= (k*x-j)/(i-j);
71 template<
class D,
class R,
int k>
76 for (
int j=0; j<=k; j++)
79 R prod( (k*1.0)/(i-j) );
80 for (
int l=0; l<=k; l++)
81 if (l!=i && l!=j) prod *= (k*x-l)/(i-l);
88 template<
class D,
class R,
int k>
93 R
p (
int i, D x)
const 96 for (
int j=0; j<=k; j++)
97 if (j!=i) result *= (k*x-j)/(i-j);
102 R
dp (
int i, D x)
const 106 for (
int j=0; j<=k; j++)
109 R prod( (k*1.0)/(i-j) );
110 for (
int l=0; l<=k; l++)
111 if (l!=i && l!=j) prod *= (k*x-l)/(i-l);
120 return i/((1.0)*(k+1));
124 template<
int k,
int d>
127 Dune::FieldVector<int,d> alpha;
128 for (
int j=0; j<d; j++)
130 alpha[j] = i % (k+1);
148 template<
class D,
class R,
int k,
int d>
153 typedef LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,d> >
Traits;
163 std::vector<typename Traits::RangeType>& out)
const 166 for (
size_t i=0; i<size(); i++)
169 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
175 for (
int j=0; j<d; j++)
176 out[i] *= p<D,R,k>(alpha[j],in[j]);
183 std::vector<typename Traits::JacobianType>& out)
const 188 for (
size_t i=0; i<size(); i++)
191 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
194 for (
int j=0; j<d; j++)
198 out[i][0][j] = dp<D,R,k>(alpha[j],in[j]);
201 for (
int l=0; l<d; l++)
203 out[i][0][j] *= p<D,R,k>(alpha[l],in[l]);
221 template<
int k,
int d>
229 li[i] = LocalKey(0,0,i);
245 std::vector<LocalKey> li;
249 template<
int k,
int d,
class LB>
255 template<
typename F,
typename C>
258 typename LB::Traits::DomainType x;
259 typename LB::Traits::RangeType y;
266 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
269 for (
int j=0; j<d; j++)
270 x[j] = (1.0*alpha[j])/k;
272 f.evaluate(x,y); out[i] = y;
278 template<
int d,
class LB>
283 template<
typename F,
typename C>
286 typename LB::Traits::DomainType x(0.5);
287 typename LB::Traits::RangeType y;
298 template<
class D,
class R,
int k,
int d>
311 typedef LocalFiniteElementTraits<LocalBasis,LocalCoefficients,LocalInterpolation>
Traits;
338 return interpolation;
355 LocalCoefficients coefficients;
356 LocalInterpolation interpolation;
366 template<
class Geometry,
class RF,
int k>
368 public ScalarLocalToGlobalFiniteElementAdaptorFactory<
369 QkDGLocalFiniteElement<
370 typename Geometry::ctype, RF, k, Geometry::mydimension
376 typename Geometry::ctype, RF, k, Geometry::mydimension
378 typedef ScalarLocalToGlobalFiniteElementAdaptorFactory<LFE, Geometry> Base;
380 static const LFE lfe;
387 template<
class Geometry,
class RF,
int k>
395 template<
class D,
class R,
int k,
int d>
402 template<
class D,
class R,
int d>
408 class QkCGLocalCoefficients
416 li[i] = LocalKey(i,d,0);
420 std::size_t size ()
const 426 const LocalKey& localKey (std::size_t i)
const 432 std::vector<LocalKey> li;
443 typedef LocalFiniteElementTraits<LocalBasis,QkCGLocalCoefficients,LocalInterpolation>
Traits;
470 return interpolation;
487 QkCGLocalCoefficients coefficients;
488 LocalInterpolation interpolation;
495 template<
class D,
class R>
502 class QkCGLocalCoefficients
510 li[0] = LocalKey(0,2,0);
511 li[1] = LocalKey(2,1,0);
512 li[2] = LocalKey(1,2,0);
513 li[3] = LocalKey(0,1,0);
514 li[4] = LocalKey(0,0,0);
515 li[5] = LocalKey(1,1,0);
516 li[6] = LocalKey(2,2,0);
517 li[7] = LocalKey(3,1,0);
518 li[8] = LocalKey(3,2,0);
522 std::size_t size ()
const 528 const LocalKey& localKey (std::size_t i)
const 534 std::vector<LocalKey> li;
545 typedef LocalFiniteElementTraits<LocalBasis,QkCGLocalCoefficients,LocalInterpolation>
Traits;
572 return interpolation;
589 QkCGLocalCoefficients coefficients;
590 LocalInterpolation interpolation;
598 template<
class D,
class R>
605 class QkCGLocalCoefficients
613 li[24] = LocalKey(6,3,0); li[25] = LocalKey(11,2,0); li[26] = LocalKey(7,3,0);
614 li[21] = LocalKey(8,2,0); li[22] = LocalKey(5,1,0); li[23] = LocalKey(9,2,0);
615 li[18] = LocalKey(4,3,0); li[19] = LocalKey(10,2,0); li[20] = LocalKey(5,3,0);
617 li[15] = LocalKey(2,2,0); li[16] = LocalKey(3,1,0); li[17] = LocalKey(3,2,0);
618 li[12] = LocalKey(0,1,0); li[13] = LocalKey(0,0,0); li[14] = LocalKey(1,1,0);
619 li[ 9] = LocalKey(0,2,0); li[10] = LocalKey(2,1,0); li[11] = LocalKey(1,2,0);
621 li[ 6] = LocalKey(2,3,0); li[ 7] = LocalKey(7,2,0); li[ 8] = LocalKey(3,3,0);
622 li[ 3] = LocalKey(4,2,0); li[ 4] = LocalKey(4,1,0); li[ 5] = LocalKey(5,2,0);
623 li[ 0] = LocalKey(0,3,0); li[ 1] = LocalKey(6,2,0); li[ 2] = LocalKey(1,3,0);
627 std::size_t size ()
const 633 const LocalKey& localKey (std::size_t i)
const 639 std::vector<LocalKey> li;
650 typedef LocalFiniteElementTraits<LocalBasis,QkCGLocalCoefficients,LocalInterpolation>
Traits;
677 return interpolation;
694 QkCGLocalCoefficients coefficients;
695 LocalInterpolation interpolation;
706 template<
class D,
class R,
int k,
int d>
723 std::size_t
size(GeometryType gt)
const 725 if (gt == GeometryType(GeometryType::cube,d))
R dp(int i, D x)
Definition: qkdg.hh:72
QkCGLocalFiniteElement()
Definition: qkdg.hh:654
Dune::FieldVector< int, d > multiindex(int i)
Definition: qkdg.hh:125
QkCGLocalFiniteElement * clone() const
Definition: qkdg.hh:480
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: qkdg.hh:182
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: qkdg.hh:461
static const unsigned int value
Definition: gridfunctionspace/tags.hh:177
QkDGLocalFiniteElement()
Definition: qkdg.hh:315
std::size_t maxLocalSize() const
Definition: qkdg.hh:731
GeometryType type() const
Definition: qkdg.hh:682
LocalFiniteElementTraits< LocalBasis, QkCGLocalCoefficients, LocalInterpolation > Traits
Definition: qkdg.hh:545
GeometryType type() const
Definition: qkdg.hh:475
Definition: adaptivity.hh:27
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: qkdg.hh:329
const Traits::LocalBasisType & localBasis() const
Definition: qkdg.hh:322
R p(int i, D x) const
Definition: qkdg.hh:93
Lagrange shape functions of order k on the reference cube.
Definition: qkdg.hh:149
std::size_t size() const
number of coefficients
Definition: qkdg.hh:233
const Traits::LocalInterpolationType & localInterpolation() const
Definition: qkdg.hh:675
GeometryType type() const
Definition: qkdg.hh:577
R dp(int i, D x) const
Definition: qkdg.hh:102
QkDGLocalFiniteElement * clone() const
Definition: qkdg.hh:348
R p(int i, D x)
Definition: qkdg.hh:62
simple implementation where all entities have the same finite element
Definition: finiteelementmap.hh:107
const Traits::LocalInterpolationType & localInterpolation() const
Definition: qkdg.hh:570
void interpolate(const F &f, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdg.hh:284
QkCGLocalFiniteElement * clone() const
Definition: qkdg.hh:582
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: qkdg.hh:668
QkDGFiniteElementFactory()
default constructor
Definition: qkdg.hh:384
Layout map for Q1 elements.
Definition: qkdg.hh:222
unsigned int order() const
Polynomial order of the shape functions.
Definition: qkdg.hh:209
const Traits::LocalInterpolationType & localInterpolation() const
Definition: qkdg.hh:468
Lagrange polynomials of degree k at equidistant points as a class.
Definition: qkdg.hh:89
R x(int i)
Definition: qkdg.hh:118
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: qkdg.hh:162
LocalFiniteElementTraits< LocalBasis, QkCGLocalCoefficients, LocalInterpolation > Traits
Definition: qkdg.hh:650
void interpolate(const F &f, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdg.hh:256
std::size_t size(GeometryType gt) const
Definition: qkdg.hh:723
QkCGLocalFiniteElement()
Definition: qkdg.hh:447
const Traits::LocalBasisType & localBasis() const
Definition: qkdg.hh:661
const Traits::LocalBasisType & localBasis() const
Definition: qkdg.hh:454
LocalFiniteElementTraits< LocalBasis, QkCGLocalCoefficients, LocalInterpolation > Traits
Definition: qkdg.hh:443
const Traits::LocalInterpolationType & localInterpolation() const
Definition: qkdg.hh:336
QkCGLocalFiniteElement()
Definition: qkdg.hh:549
Factory for global-valued QkDG elements.
Definition: qkdg.hh:367
GeometryType type() const
Definition: qkdg.hh:343
const LocalKey & localKey(std::size_t i) const
get i'th index
Definition: qkdg.hh:239
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: qkdg.hh:563
QkCGLocalFiniteElement * clone() const
Definition: qkdg.hh:687
QkDGLocalCoefficients()
Standard constructor.
Definition: qkdg.hh:226
const Traits::LocalBasisType & localBasis() const
Definition: qkdg.hh:556
LocalFiniteElementTraits< LocalBasis, LocalCoefficients, LocalInterpolation > Traits
Definition: qkdg.hh:311
bool hasDOFs(int codim) const
Definition: qkdg.hh:718
bool fixedSize() const
Definition: qkdg.hh:713
LocalBasisTraits< D, d, Dune::FieldVector< D, d >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, d > > Traits
Definition: qkdg.hh:153
unsigned int size() const
number of shape functions
Definition: qkdg.hh:156