6 #ifndef DUNE_QkDGGL_LOCALFINITEELEMENT_HH 7 #define DUNE_QkDGGL_LOCALFINITEELEMENT_HH 9 #include <dune/common/fvector.hh> 10 #include <dune/common/deprecated.hh> 12 #include <dune/geometry/type.hh> 13 #include <dune/geometry/quadraturerules.hh> 15 #include <dune/localfunctions/common/localbasis.hh> 16 #include <dune/localfunctions/common/localfiniteelementtraits.hh> 17 #include <dune/localfunctions/common/localkey.hh> 18 #include <dune/localfunctions/common/localtoglobaladaptors.hh> 29 template<
class D,
class R,
int k>
40 for (
int order=1; order<=40; order++)
42 const Dune::QuadratureRule<D,1>& rule = Dune::QuadratureRules<D,1>::rule(Dune::GeometryType::cube,order,Dune::QuadratureType::GaussLobatto);
45 matched_order = order;
46 matched_size = rule.size();
51 if (matched_order<0) DUNE_THROW(Dune::Exception,
"could not find Gauss Lobatto rule of appropriate size");
52 const Dune::QuadratureRule<D,1>& rule = Dune::QuadratureRules<D,1>::rule(Dune::GeometryType::cube,matched_order,Dune::QuadratureType::GaussLobatto);
54 for (
typename Dune::QuadratureRule<D,1>::const_iterator it=rule.begin(); it!=rule.end(); ++it)
57 size_t member=count%2;
59 if (member==1) newj=group;
else newj=k-group;
60 xi_gl[newj] = it->position()[0];
61 w_gl[newj] = it->weight();
64 for (
size_t j=0; j<matched_size/2; j++)
68 xi_gl[j] = xi_gl[k-j];
82 R
p (
int i, D
x)
const 85 for (
int j=0; j<=k; j++)
86 if (j!=i) result *= (x-xi_gl[j])/(xi_gl[i]-xi_gl[j]);
91 R
dp (
int i, D
x)
const 95 for (
int j=0; j<=k; j++)
98 R prod( 1.0/(xi_gl[i]-xi_gl[j]) );
99 for (
int l=0; l<=k; l++)
100 if (l!=i && l!=j) prod *= (x-xi_gl[l])/(xi_gl[i]-xi_gl[l]);
131 template<
class D,
class R,
int k,
int d>
138 typedef LocalBasisTraits<D,d,Dune::FieldVector<D,d>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,d> >
Traits;
148 std::vector<typename Traits::RangeType>& out)
const 151 for (
size_t i=0; i<size(); i++)
154 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
160 for (
int j=0; j<d; j++)
161 out[i] *= poly.
p(alpha[j],in[j]);
168 std::vector<typename Traits::JacobianType>& out)
const 173 for (
size_t i=0; i<size(); i++)
176 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
179 for (
int j=0; j<d; j++)
183 out[i][0][j] = poly.
dp(alpha[j],in[j]);
186 for (
int l=0; l<d; l++)
188 out[i][0][j] *= poly.
p(alpha[l],in[l]);
201 template<
int k,
int d,
class LB>
209 template<
typename F,
typename C>
212 typename LB::Traits::DomainType
x;
213 typename LB::Traits::RangeType y;
220 Dune::FieldVector<int,d> alpha(multiindex<k,d>(i));
223 for (
int j=0; j<d; j++)
224 x[j] = poly.
x(alpha[j]);
226 f.evaluate(x,y); out[i] = y;
232 template<
int d,
class LB>
237 template<
typename F,
typename C>
240 typename LB::Traits::DomainType
x(0.5);
241 typename LB::Traits::RangeType y;
252 template<
class D,
class R,
int k,
int d>
265 typedef LocalFiniteElementTraits<LocalBasis,LocalCoefficients,LocalInterpolation>
Traits;
292 return interpolation;
309 LocalCoefficients coefficients;
310 LocalInterpolation interpolation;
320 template<
class Geometry,
class RF,
int k>
322 public ScalarLocalToGlobalFiniteElementAdaptorFactory<
323 QkDGGLLocalFiniteElement<
324 typename Geometry::ctype, RF, k, Geometry::mydimension
330 typename Geometry::ctype, RF, k, Geometry::mydimension
332 typedef ScalarLocalToGlobalFiniteElementAdaptorFactory<LFE, Geometry> Base;
334 static const LFE lfe;
341 template<
class Geometry,
class RF,
int k>
352 template<
class D,
class R,
int k,
int d>
369 std::size_t
size(GeometryType gt)
const 371 if (gt == GeometryType(GeometryType::cube,d))
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: qkdggl.hh:147
Definition: qkdggl.hh:253
LocalFiniteElementTraits< LocalBasis, LocalCoefficients, LocalInterpolation > Traits
Definition: qkdggl.hh:265
Lagrange shape functions of order k on the reference cube.
Definition: qkdggl.hh:132
static const unsigned int value
Definition: gridfunctionspace/tags.hh:177
unsigned int size() const
number of shape functions
Definition: qkdggl.hh:141
R p(int i, D x) const
Definition: qkdggl.hh:82
GaussLobattoLagrangePolynomials()
Definition: qkdggl.hh:36
Definition: adaptivity.hh:27
const Traits::LocalInterpolationType & localInterpolation() const
Definition: qkdggl.hh:290
simple implementation where all entities have the same finite element
Definition: finiteelementmap.hh:107
Factory for global-valued QkDG elements.
Definition: qkdggl.hh:321
Layout map for Q1 elements.
Definition: qkdg.hh:222
std::size_t size(GeometryType gt) const
Definition: qkdggl.hh:369
LocalBasisTraits< D, d, Dune::FieldVector< D, d >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, d > > Traits
Definition: qkdggl.hh:138
Lagrange polynomials at Gauss-Lobatto points.
Definition: qkdggl.hh:30
QkDGGLLocalFiniteElement * clone() const
Definition: qkdggl.hh:302
const Traits::LocalBasisType & localBasis() const
Definition: qkdggl.hh:276
void interpolate(const F &f, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdggl.hh:210
R dp(int i, D x) const
Definition: qkdggl.hh:91
bool fixedSize() const
Definition: qkdggl.hh:359
QkDGGLFiniteElementFactory()
default constructor
Definition: qkdggl.hh:338
Definition: qkdggl.hh:353
GeometryType type() const
Definition: qkdggl.hh:297
bool hasDOFs(int codim) const
Definition: qkdggl.hh:364
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: qkdggl.hh:167
R w(int i) const
Definition: qkdggl.hh:113
QkDGGLLocalFiniteElement()
Definition: qkdggl.hh:269
R x(int i) const
Definition: qkdggl.hh:107
void interpolate(const F &f, std::vector< C > &out) const
Local interpolation of a function.
Definition: qkdggl.hh:238
Definition: qkdggl.hh:202
std::size_t maxLocalSize() const
Definition: qkdggl.hh:377
unsigned int order() const
Polynomial order of the shape functions.
Definition: qkdggl.hh:194
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: qkdggl.hh:283