Go to the documentation of this file.
4 #ifndef DUNE_PDELAB_FINITEELEMENT_QKDGLAGRANGE_HH
5 #define DUNE_PDELAB_FINITEELEMENT_QKDGLAGRANGE_HH
7 #include <dune/localfunctions/common/localbasis.hh>
8 #include <dune/localfunctions/common/localkey.hh>
9 #include <dune/localfunctions/common/localfiniteelementtraits.hh>
10 #include <dune/localfunctions/common/localtoglobaladaptors.hh>
20 template<
int k,
int n>
53 template<
class D,
class R,
int k>
57 for (
int j=0; j<=k; j++)
58 if (j!=i) result *= (k*x-j)/(i-j);
63 template<
class D,
class R,
int k>
68 for (
int j=0; j<=k; j++)
71 R prod( (k*1.0)/(i-j) );
72 for (
int l=0; l<=k; l++)
73 if (l!=i && l!=j) prod *= (k*x-l)/(i-l);
80 template<
class D,
class R,
int k>
85 R
p (
int i, D
x)
const
88 for (
int j=0; j<=k; j++)
89 if (j!=i) result *= (k*
x-j)/(i-j);
94 R
dp (
int i, D
x)
const
98 for (
int j=0; j<=k; j++)
101 R prod( (k*1.0)/(i-j) );
102 for (
int l=0; l<=k; l++)
103 if (l!=i && l!=j) prod *= (k*
x-l)/(i-l);
112 return i/((1.0)*(k+1));
116 template<
int k,
int d>
119 Dune::FieldVector<int,d> alpha;
120 for (
int j=0; j<d; j++)
122 alpha[j] = i % (k+1);
140 template<
class D,
class R,
int k,
int d>
145 typedef LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,d> >
Traits;
155 std::vector<typename Traits::RangeType>& out)
const
158 for (
size_t i=0; i<
size(); i++)
161 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
167 for (
int j=0; j<d; j++)
168 out[i] *= p<D,R,k>(alpha[j],in[j]);
175 std::vector<typename Traits::JacobianType>& out)
const
180 for (
size_t i=0; i<
size(); i++)
183 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
186 for (
int j=0; j<d; j++)
190 out[i][0][j] = dp<D,R,k>(alpha[j],in[j]);
193 for (
int l=0; l<d; l++)
195 out[i][0][j] *= p<D,R,k>(alpha[l],in[l]);
213 template<
int k,
int d>
221 li[i] = LocalKey(0,0,i);
237 std::vector<LocalKey> li;
241 template<
int k,
int d,
class LB>
247 template<
typename F,
typename C>
250 typename LB::Traits::DomainType x;
251 typename LB::Traits::RangeType y;
258 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
261 for (
int j=0; j<d; j++)
262 x[j] = (1.0*alpha[j])/k;
264 f.evaluate(x,y); out[i] = y;
270 template<
int d,
class LB>
275 template<
typename F,
typename C>
278 typename LB::Traits::DomainType x(0.5);
279 typename LB::Traits::RangeType y;
290 template<
class D,
class R,
int k,
int d>
303 typedef LocalFiniteElementTraits<LocalBasis,LocalCoefficients,LocalInterpolation>
Traits;
308 : gt(GeometryTypes::cube(d))
329 return interpolation;
346 LocalCoefficients coefficients;
347 LocalInterpolation interpolation;
357 template<
class Geometry,
class RF,
int k>
359 public ScalarLocalToGlobalFiniteElementAdaptorFactory<
360 QkDGLagrangeLocalFiniteElement<
361 typename Geometry::ctype, RF, k, Geometry::mydimension
367 typename Geometry::ctype, RF, k, Geometry::mydimension
369 typedef ScalarLocalToGlobalFiniteElementAdaptorFactory<LFE, Geometry> Base;
371 static const LFE lfe;
378 template<
class Geometry,
class RF,
int k>
379 const typename QkDGFiniteElementFactory<Geometry, RF, k>::LFE
380 QkDGFiniteElementFactory<Geometry, RF, k>::lfe;
385 #endif // DUNE_PDELAB_FINITEELEMENT_QKDGLAGRANGE_HH
QkDGLagrangeLocalFiniteElement()
Definition: qkdglagrange.hh:307
Layout map for Q1 elements.
Definition: qkdglagrange.hh:214
R x(int i)
Definition: qkdglagrange.hh:110
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: qkdglagrange.hh:174
Factory for global-valued QkDG elements.
Definition: qkdglagrange.hh:358
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
QkDGLagrangeLocalFiniteElement * clone() const
Definition: qkdglagrange.hh:339
@ n
Definition: qkdglagrange.hh:299
Definition: qkdglagrange.hh:291
R p(int i, D x) const
Definition: qkdglagrange.hh:85
const LocalKey & localKey(std::size_t i) const
get i'th index
Definition: qkdglagrange.hh:231
@ value
Definition: qkdglagrange.hh:24
QkDGLocalCoefficients()
Standard constructor.
Definition: qkdglagrange.hh:218
Definition: qkdglagrange.hh:242
const Traits::LocalInterpolationType & localInterpolation() const
Definition: qkdglagrange.hh:327
Lagrange shape functions of order k on the reference cube.
Definition: qkdglagrange.hh:141
Dune::FieldVector< int, d > multiindex(int i)
Definition: qkdglagrange.hh:117
LocalFiniteElementTraits< LocalBasis, LocalCoefficients, LocalInterpolation > Traits
Definition: qkdglagrange.hh:303
R dp(int i, D x)
Definition: qkdglagrange.hh:64
void interpolate(const F &f, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdglagrange.hh:248
Definition: qkdglagrange.hh:21
Lagrange polynomials of degree k at equidistant points as a class.
Definition: qkdglagrange.hh:81
R dp(int i, D x) const
Definition: qkdglagrange.hh:94
unsigned int size() const
number of shape functions
Definition: qkdglagrange.hh:148
LocalBasisTraits< D, d, Dune::FieldVector< D, d >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, d > > Traits
Definition: qkdglagrange.hh:145
const Traits::LocalBasisType & localBasis() const
Definition: qkdglagrange.hh:313
unsigned int order() const
Polynomial order of the shape functions.
Definition: qkdglagrange.hh:201
GeometryType type() const
Definition: qkdglagrange.hh:334
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
std::size_t size() const
number of coefficients
Definition: qkdglagrange.hh:225
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: qkdglagrange.hh:154
QkDGFiniteElementFactory()
default constructor
Definition: qkdglagrange.hh:375
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: qkdglagrange.hh:320
void interpolate(const F &f, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdglagrange.hh:276
R p(int i, D x)
Definition: qkdglagrange.hh:54