6 #ifndef DUNE_PDELAB_FINITEELEMENT_QKDGLEGENDRE_HH 7 #define DUNE_PDELAB_FINITEELEMENT_QKDGLEGENDRE_HH 11 #include <dune/common/fvector.hh> 12 #include <dune/common/deprecated.hh> 14 #include <dune/geometry/type.hh> 15 #include <dune/geometry/quadraturerules.hh> 17 #include <dune/localfunctions/common/localbasis.hh> 18 #include <dune/localfunctions/common/localfiniteelementtraits.hh> 19 #include <dune/localfunctions/common/localkey.hh> 20 #include <dune/localfunctions/common/localtoglobaladaptors.hh> 25 namespace LegendreStuff
30 template<
int k,
int n>
62 template<
int k,
int d>
65 Dune::FieldVector<int,d> alpha;
66 for (
int j=0; j<d; j++)
75 template<
class D,
class R,
int k>
80 void p (D x, std::vector<R>&
value)
const 85 for (
int n=2; n<=k; n++)
86 value[n] = ((2*n-1)*(2*x-1)*value[n-1]-(n-1)*value[n-2])/n;
90 R
p (
int i, D x)
const 92 std::vector<R> v(k+1);
98 void dp (D x, std::vector<R>& derivative)
const 101 derivative.resize(k+1);
106 for (
int n=2; n<=k; n++)
108 value[n] = ((2*n-1)*(2*x-1)*value[n-1]-(n-1)*value[n-2])/n;
109 derivative[n] = (2*x-1)*derivative[n-1] + 2*n*value[n-1];
114 void pdp (D x, std::vector<R>&
value, std::vector<R>& derivative)
const 117 derivative.resize(k+1);
122 for (
int n=2; n<=k; n++)
124 value[n] = ((2*n-1)*(2*x-1)*value[n-1]-(n-1)*value[n-2])/n;
125 derivative[n] = (2*x-1)*derivative[n-1] + 2*n*value[n-1];
130 R
dp (
int i, D x)
const 132 std::vector<R> v(k+1);
138 template<
class D,
class R>
143 void p (D x, std::vector<R>&
value)
const 150 R
p (
int i, D x)
const 156 void dp (D x, std::vector<R>& derivative)
const 158 derivative.resize(1);
163 R
dp (
int i, D x)
const 169 void pdp (D x, std::vector<R>&
value, std::vector<R>& derivative)
const 172 derivative.resize(1);
178 template<
class D,
class R>
183 void p (D x, std::vector<R>&
value)
const 191 R
p (
int i, D x)
const 193 return (1-i) + i*(2*x-1);
197 void dp (D x, std::vector<R>& derivative)
const 199 derivative.resize(2);
205 R
dp (
int i, D x)
const 207 return (1-i)*0 + i*(2);
211 void pdp (D x, std::vector<R>&
value, std::vector<R>& derivative)
const 214 derivative.resize(2);
234 template<
class D,
class R,
int k,
int d>
239 mutable std::vector<std::vector<R> > v;
240 mutable std::vector<std::vector<R> > a;
243 typedef LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,d> >
Traits;
256 std::vector<typename Traits::RangeType>&
value)
const 262 for (
size_t j=0; j<d; j++) poly.
p(x[j],v[j]);
265 for (
size_t i=0; i<n; i++)
268 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
274 for (
int j=0; j<d; j++) value[i] *= v[j][alpha[j]];
281 std::vector<typename Traits::JacobianType>&
value)
const 284 value.resize(size());
287 for (
size_t j=0; j<d; j++) poly.
pdp(x[j],v[j],a[j]);
290 for (
size_t i=0; i<n; i++)
293 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
296 for (
int j=0; j<d; j++)
299 value[i][0][j] = a[j][alpha[j]];
302 for (
int l=0; l<d; l++)
304 value[i][0][j] *= v[l][alpha[l]];
318 template<
class D,
class R,
int k,
int d>
324 Dune::GeometryType gt;
332 template<
typename F,
typename C>
336 typedef typename LB::Traits::RangeType RangeType;
337 const Dune::QuadratureRule<R,d>&
338 rule = Dune::QuadratureRules<R,d>::rule(gt,2*k);
342 std::vector<R> diagonal(n);
343 for (
int i=0; i<n; i++) { out[i] = 0.0; diagonal[i] = 0.0; }
346 for (
typename Dune::QuadratureRule<R,d>::const_iterator
347 it=rule.begin(); it!=rule.end(); ++it)
350 typename LB::Traits::DomainType x;
352 for (
int i=0; i<d; i++) x[i] = it->position()[i];
356 std::vector<RangeType> phi(n);
360 for (
int i=0; i<n; i++) {
361 out[i] += y*phi[i]*it->weight();
362 diagonal[i] += phi[i]*phi[i]*it->weight();
365 for (
int i=0; i<n; i++) out[i] /= diagonal[i];
375 template<
int k,
int d>
384 for (std::size_t i=0; i<n; i++)
385 li[i] = LocalKey(0,0,i);
401 std::vector<LocalKey> li;
408 template<
class D,
class R,
int k,
int d>
421 typedef LocalFiniteElementTraits<LocalBasis,LocalCoefficients,LocalInterpolation>
Traits;
448 return interpolation;
465 LocalCoefficients coefficients;
466 LocalInterpolation interpolation;
476 template<
class Geometry,
class RF,
int k>
478 public ScalarLocalToGlobalFiniteElementAdaptorFactory<
479 QkDGLegendreLocalFiniteElement<
480 typename Geometry::ctype, RF, k, Geometry::mydimension
486 typename Geometry::ctype, RF, k, Geometry::mydimension
488 typedef ScalarLocalToGlobalFiniteElementAdaptorFactory<LFE, Geometry> Base;
490 static const LFE lfe;
497 template<
class Geometry,
class RF,
int k>
503 #endif // DUNE_PDELAB_FINITEELEMENT_QKDGLEGENDRE_HH QkDGLegendreLocalFiniteElement()
Definition: qkdglegendre.hh:425
const LocalKey & localKey(std::size_t i) const
get i'th index
Definition: qkdglegendre.hh:395
Lagrange shape functions of order k on the reference cube.
Definition: qkdglegendre.hh:235
unsigned int size() const
number of shape functions
Definition: qkdglegendre.hh:249
The 1d Legendre Polynomials (k=0,1 are specialized below)
Definition: qkdglegendre.hh:76
Definition: qkdglegendre.hh:409
DGLegendreLocalCoefficients()
Standard constructor.
Definition: qkdglegendre.hh:382
void pdp(D x, std::vector< R > &value, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:169
GeometryType type() const
Definition: qkdglegendre.hh:453
void p(D x, std::vector< R > &value) const
evaluate all polynomials at point x
Definition: qkdglegendre.hh:183
R dp(int i, D x) const
Definition: qkdglegendre.hh:130
R dp(int i, D x) const
Definition: qkdglegendre.hh:163
void dp(D x, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:156
DGLegendreLocalBasis()
Definition: qkdglegendre.hh:245
void dp(D x, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:98
R p(int i, D x) const
Definition: qkdglegendre.hh:90
R dp(int i, D x)
Definition: qkdglagrange.hh:64
void pdp(D x, std::vector< R > &value, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:114
DGLegendreLocalInterpolation()
Definition: qkdglegendre.hh:328
Definition: qkdglegendre.hh:31
const Traits::LocalBasisType & localBasis() const
Definition: qkdglegendre.hh:432
LocalFiniteElementTraits< LocalBasis, LocalCoefficients, LocalInterpolation > Traits
Definition: qkdglegendre.hh:421
const Traits::LocalInterpolationType & localInterpolation() const
Definition: qkdglegendre.hh:446
determine degrees of freedom
Definition: qkdglegendre.hh:319
LocalBasisTraits< D, d, Dune::FieldVector< D, d >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, d > > Traits
Definition: qkdglegendre.hh:243
void interpolate(const F &f, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdglegendre.hh:333
const P & p
Definition: constraints.hh:147
void pdp(D x, std::vector< R > &value, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:211
Dune::FieldVector< int, d > multiindex(int i)
Definition: qkdglegendre.hh:63
void dp(D x, std::vector< R > &derivative) const
Definition: qkdglegendre.hh:197
QkDGLegendreFiniteElementFactory()
default constructor
Definition: qkdglegendre.hh:494
Definition: qkdglegendre.hh:34
R p(int i, D x) const
Definition: qkdglegendre.hh:191
void evaluateFunction(const typename Traits::DomainType &x, std::vector< typename Traits::RangeType > &value) const
Evaluate all shape functions.
Definition: qkdglegendre.hh:255
void evaluateJacobian(const typename Traits::DomainType &x, std::vector< typename Traits::JacobianType > &value) const
Evaluate Jacobian of all shape functions.
Definition: qkdglegendre.hh:280
R dp(int i, D x) const
Definition: qkdglegendre.hh:205
Layout map for Q1 elements.
Definition: qkdglegendre.hh:376
R p(int i, D x) const
Definition: qkdglegendre.hh:150
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
unsigned int order() const
Polynomial order of the shape functions.
Definition: qkdglegendre.hh:310
void p(D x, std::vector< R > &value) const
evaluate all polynomials at point x
Definition: qkdglegendre.hh:143
void p(D x, std::vector< R > &value) const
evaluate all polynomials at point x
Definition: qkdglegendre.hh:80
QkDGLegendreLocalFiniteElement * clone() const
Definition: qkdglegendre.hh:458
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: qkdglegendre.hh:439
std::size_t size() const
number of coefficients
Definition: qkdglegendre.hh:389
Factory for global-valued DGLegendre elements.
Definition: qkdglegendre.hh:477