6 #ifndef DUNE_PDELAB_FINITEELEMENT_QKDGLEGENDRE_HH
7 #define DUNE_PDELAB_FINITEELEMENT_QKDGLEGENDRE_HH
12 #include <dune/common/fvector.hh>
13 #include <dune/common/deprecated.hh>
15 #include <dune/geometry/type.hh>
16 #include <dune/geometry/quadraturerules.hh>
18 #include <dune/localfunctions/common/localbasis.hh>
19 #include <dune/localfunctions/common/localfiniteelementtraits.hh>
20 #include <dune/localfunctions/common/localkey.hh>
21 #include <dune/localfunctions/common/localtoglobaladaptors.hh>
26 namespace LegendreStuff
31 template<
int k,
int n>
63 template<
int k,
int d>
66 Dune::FieldVector<int,d> alpha;
67 for (
int j=0; j<d; j++)
76 template<
class D,
class R,
int k>
81 void p (D x, std::vector<R>&
value)
const
86 for (
int n=2; n<=k; n++)
91 R
p (
int i, D x)
const
93 std::vector<R> v(k+1);
99 void dp (D x, std::vector<R>& derivative)
const
102 derivative.resize(k+1);
107 for (
int n=2; n<=k; n++)
110 derivative[n] = (2*x-1)*derivative[n-1] + 2*n*
value[n-1];
115 void pdp (D x, std::vector<R>&
value, std::vector<R>& derivative)
const
118 derivative.resize(k+1);
123 for (
int n=2; n<=k; n++)
126 derivative[n] = (2*x-1)*derivative[n-1] + 2*n*
value[n-1];
131 R
dp (
int i, D x)
const
133 std::vector<R> v(k+1);
139 template<
class D,
class R>
144 void p (D x, std::vector<R>&
value)
const
151 R
p (
int i, D x)
const
157 void dp (D x, std::vector<R>& derivative)
const
159 derivative.resize(1);
164 R
dp (
int i, D x)
const
170 void pdp (D x, std::vector<R>&
value, std::vector<R>& derivative)
const
173 derivative.resize(1);
179 template<
class D,
class R>
184 void p (D x, std::vector<R>&
value)
const
192 R
p (
int i, D x)
const
194 return (1-i) + i*(2*x-1);
198 void dp (D x, std::vector<R>& derivative)
const
200 derivative.resize(2);
206 R
dp (
int i, D x)
const
208 return (1-i)*0 + i*(2);
212 void pdp (D x, std::vector<R>&
value, std::vector<R>& derivative)
const
215 derivative.resize(2);
235 template<
class D,
class R,
int k,
int d>
240 mutable std::vector<std::vector<R> > v;
241 mutable std::vector<std::vector<R> > a;
244 typedef LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,d> >
Traits;
257 std::vector<typename Traits::RangeType>&
value)
const
263 for (
size_t j=0; j<d; j++) poly.
p(x[j],v[j]);
266 for (
size_t i=0; i<n; i++)
269 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
275 for (
int j=0; j<d; j++)
value[i] *= v[j][alpha[j]];
282 std::vector<typename Traits::JacobianType>&
value)
const
288 for (
size_t j=0; j<d; j++) poly.
pdp(x[j],v[j],a[j]);
291 for (
size_t i=0; i<n; i++)
294 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
297 for (
int j=0; j<d; j++)
300 value[i][0][j] = a[j][alpha[j]];
303 for (
int l=0; l<d; l++)
305 value[i][0][j] *= v[l][alpha[l]];
311 void partial(
const std::array<unsigned int, Traits::dimDomain>& DUNE_UNUSED(
order),
312 const typename Traits::DomainType& DUNE_UNUSED(in),
313 std::vector<typename Traits::RangeType>& DUNE_UNUSED(out))
const {
314 auto totalOrder = std::accumulate(
order.begin(),
order.end(), 0);
315 if (totalOrder == 0) {
318 DUNE_THROW(NotImplemented,
"Desired derivative order is not implemented");
331 template<
class D,
class R,
int k,
int d>
337 Dune::GeometryType gt;
345 template<
typename F,
typename C>
349 typedef typename LB::Traits::RangeType RangeType;
350 const Dune::QuadratureRule<R,d>&
351 rule = Dune::QuadratureRules<R,d>::rule(gt,2*k);
355 std::vector<R> diagonal(n);
356 for (
int i=0; i<n; i++) { out[i] = 0.0; diagonal[i] = 0.0; }
359 for (
typename Dune::QuadratureRule<R,d>::const_iterator
360 it=rule.begin(); it!=rule.end(); ++it)
363 typename LB::Traits::DomainType x;
365 for (
int i=0; i<d; i++) x[i] = it->position()[i];
369 std::vector<RangeType> phi(n);
373 for (
int i=0; i<n; i++) {
374 out[i] += y*phi[i]*it->weight();
375 diagonal[i] += phi[i]*phi[i]*it->weight();
378 for (
int i=0; i<n; i++) out[i] /= diagonal[i];
388 template<
int k,
int d>
397 for (std::size_t i=0; i<n; i++)
398 li[i] = LocalKey(0,0,i);
414 std::vector<LocalKey> li;
421 template<
class D,
class R,
int k,
int d>
434 typedef LocalFiniteElementTraits<LocalBasis,LocalCoefficients,LocalInterpolation>
Traits;
439 : gt(GeometryTypes::cube(d))
460 return interpolation;
484 LocalCoefficients coefficients;
485 LocalInterpolation interpolation;
495 template<
class Geometry,
class RF,
int k>
497 public ScalarLocalToGlobalFiniteElementAdaptorFactory<
498 QkDGLegendreLocalFiniteElement<
499 typename Geometry::ctype, RF, k, Geometry::mydimension
505 typename Geometry::ctype, RF, k, Geometry::mydimension
507 typedef ScalarLocalToGlobalFiniteElementAdaptorFactory<LFE, Geometry> Base;
509 static const LFE lfe;
516 template<
class Geometry,
class RF,
int k>
517 const typename QkDGLegendreFiniteElementFactory<Geometry, RF, k>::LFE
518 QkDGLegendreFiniteElementFactory<Geometry, RF, k>::lfe;
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::FieldVector< int, d > multiindex(int i)
Definition: qkdglegendre.hh:64
Definition: qkdglegendre.hh:33
@ value
Definition: qkdglegendre.hh:35
The 1d Legendre Polynomials (k=0,1 are specialized below)
Definition: qkdglegendre.hh:78
R p(int i, D x) const
Definition: qkdglegendre.hh:91
void pdp(D x, std::vector< R > &value, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:115
void p(D x, std::vector< R > &value) const
evaluate all polynomials at point x
Definition: qkdglegendre.hh:81
void dp(D x, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:99
R dp(int i, D x) const
Definition: qkdglegendre.hh:131
void pdp(D x, std::vector< R > &value, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:170
void p(D x, std::vector< R > &value) const
evaluate all polynomials at point x
Definition: qkdglegendre.hh:144
void dp(D x, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:157
R p(int i, D x) const
Definition: qkdglegendre.hh:151
R dp(int i, D x) const
Definition: qkdglegendre.hh:164
void dp(D x, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:198
R p(int i, D x) const
Definition: qkdglegendre.hh:192
void p(D x, std::vector< R > &value) const
evaluate all polynomials at point x
Definition: qkdglegendre.hh:184
R dp(int i, D x) const
Definition: qkdglegendre.hh:206
void pdp(D x, std::vector< R > &value, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:212
Lagrange shape functions of order k on the reference cube.
Definition: qkdglegendre.hh:237
unsigned int size() const
number of shape functions
Definition: qkdglegendre.hh:250
void evaluateFunction(const typename Traits::DomainType &x, std::vector< typename Traits::RangeType > &value) const
Evaluate all shape functions.
Definition: qkdglegendre.hh:256
void evaluateJacobian(const typename Traits::DomainType &x, std::vector< typename Traits::JacobianType > &value) const
Evaluate Jacobian of all shape functions.
Definition: qkdglegendre.hh:281
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: qkdglegendre.hh:311
DGLegendreLocalBasis()
Definition: qkdglegendre.hh:246
unsigned int order() const
Polynomial order of the shape functions.
Definition: qkdglegendre.hh:323
LocalBasisTraits< D, d, Dune::FieldVector< D, d >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, d > > Traits
Definition: qkdglegendre.hh:244
determine degrees of freedom
Definition: qkdglegendre.hh:333
void interpolate(const F &f, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdglegendre.hh:346
DGLegendreLocalInterpolation()
Definition: qkdglegendre.hh:341
Layout map for Q1 elements.
Definition: qkdglegendre.hh:390
const LocalKey & localKey(std::size_t i) const
get i'th index
Definition: qkdglegendre.hh:408
DGLegendreLocalCoefficients()
Standard constructor.
Definition: qkdglegendre.hh:395
std::size_t size() const
number of coefficients
Definition: qkdglegendre.hh:402
Definition: qkdglegendre.hh:423
const Traits::LocalBasisType & localBasis() const
Definition: qkdglegendre.hh:444
LocalFiniteElementTraits< LocalBasis, LocalCoefficients, LocalInterpolation > Traits
Definition: qkdglegendre.hh:434
QkDGLegendreLocalFiniteElement * clone() const
Definition: qkdglegendre.hh:477
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: qkdglegendre.hh:451
GeometryType type() const
Definition: qkdglegendre.hh:472
std::size_t size() const
Definition: qkdglegendre.hh:465
const Traits::LocalInterpolationType & localInterpolation() const
Definition: qkdglegendre.hh:458
@ n
Definition: qkdglegendre.hh:430
QkDGLegendreLocalFiniteElement()
Definition: qkdglegendre.hh:438
Factory for global-valued DGLegendre elements.
Definition: qkdglegendre.hh:503
QkDGLegendreFiniteElementFactory()
default constructor
Definition: qkdglegendre.hh:513
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139