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/localinterpolation.hh>
21#include <dune/localfunctions/common/localkey.hh>
22#include <dune/localfunctions/common/localtoglobaladaptors.hh>
23#include <dune/localfunctions/common/localinterpolation.hh>
28 namespace LegendreStuff
33 template<
int k,
int n>
65 template<
int k,
int d>
68 Dune::FieldVector<int,d> alpha;
69 for (
int j=0; j<d; j++)
78 template<
class D,
class R,
int k>
83 void p (D x, std::vector<R>&
value)
const
88 for (
int n=2; n<=k; n++)
93 R
p (
int i, D x)
const
95 std::vector<R> v(k+1);
101 void dp (D x, std::vector<R>& derivative)
const
104 derivative.resize(k+1);
109 for (
int n=2; n<=k; n++)
112 derivative[n] = (2*x-1)*derivative[n-1] + 2*n*
value[n-1];
117 void pdp (D x, std::vector<R>&
value, std::vector<R>& derivative)
const
120 derivative.resize(k+1);
125 for (
int n=2; n<=k; n++)
128 derivative[n] = (2*x-1)*derivative[n-1] + 2*n*
value[n-1];
133 R
dp (
int i, D x)
const
135 std::vector<R> v(k+1);
141 template<
class D,
class R>
146 void p (D x, std::vector<R>&
value)
const
153 R
p (
int i, D x)
const
159 void dp (D x, std::vector<R>& derivative)
const
161 derivative.resize(1);
166 R
dp (
int i, D x)
const
172 void pdp (D x, std::vector<R>&
value, std::vector<R>& derivative)
const
175 derivative.resize(1);
181 template<
class D,
class R>
186 void p (D x, std::vector<R>&
value)
const
194 R
p (
int i, D x)
const
196 return (1-i) + i*(2*x-1);
200 void dp (D x, std::vector<R>& derivative)
const
202 derivative.resize(2);
208 R
dp (
int i, D x)
const
210 return (1-i)*0 + i*(2);
214 void pdp (D x, std::vector<R>&
value, std::vector<R>& derivative)
const
217 derivative.resize(2);
237 template<
class D,
class R,
int k,
int d>
242 mutable std::vector<std::vector<R> > v;
243 mutable std::vector<std::vector<R> > a;
246 typedef LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,d> >
Traits;
259 std::vector<typename Traits::RangeType>&
value)
const
265 for (
size_t j=0; j<d; j++) poly.
p(x[j],v[j]);
268 for (
size_t i=0; i<n; i++)
271 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
277 for (
int j=0; j<d; j++)
value[i] *= v[j][alpha[j]];
284 std::vector<typename Traits::JacobianType>&
value)
const
290 for (
size_t j=0; j<d; j++) poly.
pdp(x[j],v[j],a[j]);
293 for (
size_t i=0; i<n; i++)
296 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
299 for (
int j=0; j<d; j++)
302 value[i][0][j] = a[j][alpha[j]];
305 for (
int l=0; l<d; l++)
307 value[i][0][j] *= v[l][alpha[l]];
314 const typename Traits::DomainType& in,
315 std::vector<typename Traits::RangeType>& out)
const {
316 auto totalOrder = std::accumulate(
order.begin(),
order.end(), 0);
317 if (totalOrder == 0) {
320 DUNE_THROW(NotImplemented,
"Desired derivative order is not implemented");
333 template<
class D,
class R,
int k,
int d>
339 Dune::GeometryType gt;
347 template<
typename F,
typename C>
351 typedef typename LB::Traits::RangeType RangeType;
353 auto&& f = Impl::makeFunctionWithCallOperator<typename LB::Traits::DomainType>(ff);
354 const Dune::QuadratureRule<R,d>&
355 rule = Dune::QuadratureRules<R,d>::rule(gt,2*k);
359 std::vector<R> diagonal(n);
360 for (
int i=0; i<n; i++) { out[i] = 0.0; diagonal[i] = 0.0; }
363 for (
typename Dune::QuadratureRule<R,d>::const_iterator
364 it=rule.begin(); it!=rule.end(); ++it)
367 typename LB::Traits::DomainType x;
369 for (
int i=0; i<d; i++) x[i] = it->position()[i];
373 std::vector<RangeType> phi(n);
377 for (
int i=0; i<n; i++) {
378 out[i] += y*phi[i]*it->weight();
379 diagonal[i] += phi[i]*phi[i]*it->weight();
382 for (
int i=0; i<n; i++) out[i] /= diagonal[i];
392 template<
int k,
int d>
401 for (std::size_t i=0; i<n; i++)
402 li[i] = LocalKey(0,0,i);
418 std::vector<LocalKey> li;
425 template<
class D,
class R,
int k,
int d>
438 typedef LocalFiniteElementTraits<LocalBasis,LocalCoefficients,LocalInterpolation>
Traits;
443 : gt(GeometryTypes::cube(d))
464 return interpolation;
488 LocalCoefficients coefficients;
489 LocalInterpolation interpolation;
499 template<
class Geometry,
class RF,
int k>
501 public ScalarLocalToGlobalFiniteElementAdaptorFactory<
502 QkDGLegendreLocalFiniteElement<
503 typename Geometry::ctype, RF, k, Geometry::mydimension
509 typename Geometry::ctype, RF, k, Geometry::mydimension
511 typedef ScalarLocalToGlobalFiniteElementAdaptorFactory<LFE, Geometry> Base;
513 static const LFE lfe;
520 template<
class Geometry,
class RF,
int k>
521 const typename QkDGLegendreFiniteElementFactory<Geometry, RF, k>::LFE
522 QkDGLegendreFiniteElementFactory<Geometry, RF, k>::lfe;
const P & p
Definition: constraints.hh:148
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Dune::FieldVector< int, d > multiindex(int i)
Definition: qkdglegendre.hh:66
Definition: qkdglegendre.hh:35
@ value
Definition: qkdglegendre.hh:37
The 1d Legendre Polynomials (k=0,1 are specialized below)
Definition: qkdglegendre.hh:80
R p(int i, D x) const
Definition: qkdglegendre.hh:93
void pdp(D x, std::vector< R > &value, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:117
void p(D x, std::vector< R > &value) const
evaluate all polynomials at point x
Definition: qkdglegendre.hh:83
void dp(D x, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:101
R dp(int i, D x) const
Definition: qkdglegendre.hh:133
void pdp(D x, std::vector< R > &value, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:172
void p(D x, std::vector< R > &value) const
evaluate all polynomials at point x
Definition: qkdglegendre.hh:146
void dp(D x, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:159
R p(int i, D x) const
Definition: qkdglegendre.hh:153
R dp(int i, D x) const
Definition: qkdglegendre.hh:166
void dp(D x, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:200
R p(int i, D x) const
Definition: qkdglegendre.hh:194
void p(D x, std::vector< R > &value) const
evaluate all polynomials at point x
Definition: qkdglegendre.hh:186
R dp(int i, D x) const
Definition: qkdglegendre.hh:208
void pdp(D x, std::vector< R > &value, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:214
Lagrange shape functions of order k on the reference cube.
Definition: qkdglegendre.hh:239
unsigned int size() const
number of shape functions
Definition: qkdglegendre.hh:252
void evaluateFunction(const typename Traits::DomainType &x, std::vector< typename Traits::RangeType > &value) const
Evaluate all shape functions.
Definition: qkdglegendre.hh:258
void evaluateJacobian(const typename Traits::DomainType &x, std::vector< typename Traits::JacobianType > &value) const
Evaluate Jacobian of all shape functions.
Definition: qkdglegendre.hh:283
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: qkdglegendre.hh:313
DGLegendreLocalBasis()
Definition: qkdglegendre.hh:248
unsigned int order() const
Polynomial order of the shape functions.
Definition: qkdglegendre.hh:325
LocalBasisTraits< D, d, Dune::FieldVector< D, d >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, d > > Traits
Definition: qkdglegendre.hh:246
determine degrees of freedom
Definition: qkdglegendre.hh:335
DGLegendreLocalInterpolation()
Definition: qkdglegendre.hh:343
void interpolate(const F &ff, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdglegendre.hh:348
Layout map for Q1 elements.
Definition: qkdglegendre.hh:394
DGLegendreLocalCoefficients()
Standard constructor.
Definition: qkdglegendre.hh:399
std::size_t size() const
number of coefficients
Definition: qkdglegendre.hh:406
const LocalKey & localKey(std::size_t i) const
get i'th index
Definition: qkdglegendre.hh:412
Definition: qkdglegendre.hh:427
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: qkdglegendre.hh:455
LocalFiniteElementTraits< LocalBasis, LocalCoefficients, LocalInterpolation > Traits
Definition: qkdglegendre.hh:438
const Traits::LocalInterpolationType & localInterpolation() const
Definition: qkdglegendre.hh:462
QkDGLegendreLocalFiniteElement * clone() const
Definition: qkdglegendre.hh:481
GeometryType type() const
Definition: qkdglegendre.hh:476
std::size_t size() const
Definition: qkdglegendre.hh:469
const Traits::LocalBasisType & localBasis() const
Definition: qkdglegendre.hh:448
@ n
Definition: qkdglegendre.hh:434
QkDGLegendreLocalFiniteElement()
Definition: qkdglegendre.hh:442
Factory for global-valued DGLegendre elements.
Definition: qkdglegendre.hh:507
QkDGLegendreFiniteElementFactory()
default constructor
Definition: qkdglegendre.hh:517
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139