1 #ifndef DUNE_NPDE_PK1D_HH
2 #define DUNE_NPDE_PK1D_HH
6 #include<dune/common/fvector.hh>
7 #include<dune/common/fmatrix.hh>
8 #include<dune/common/exceptions.hh>
10 #include<dune/geometry/type.hh>
11 #include<dune/geometry/referenceelements.hh>
12 #include<dune/geometry/quadraturerules.hh>
14 #include<dune/localfunctions/common/localbasis.hh>
15 #include<dune/localfunctions/common/localkey.hh>
16 #include<dune/localfunctions/common/localfiniteelementtraits.hh>
27 template<
class D,
class R>
33 Dune::GeometryType gt;
38 typedef Dune::LocalBasisTraits<D,1,Dune::FieldVector<D,1>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,1>, 1>
Traits;
41 Pk1dLocalBasis (std::size_t k_) : gt(Dune::GeometryType::cube,1), k(k_), n(k_+1),
s(n)
43 for (std::size_t i=0; i<=k; i++) s[i] = (1.0*i)/k;
47 std::size_t size ()
const {
return n; }
50 inline void evaluateFunction (
const typename Traits::DomainType& in,
51 std::vector<typename Traits::RangeType>& out)
const {
53 for (std::size_t i=0; i<=k; i++)
56 for (std::size_t j=0; j<=k; j++)
57 if (i!=j) out[i] *= (in[0]-s[j])/(s[i]-s[j]);
63 evaluateJacobian (
const typename Traits::DomainType& in,
64 std::vector<typename Traits::JacobianType>& out)
const {
66 for (std::size_t i=0; i<=k; i++)
71 for (std::size_t j=0; j<=k; j++)
74 denominator *= s[i]-s[j];
76 for (std::size_t l=j+1; l<=k; l++)
81 out[i][0][0] += factor*a;
84 out[i][0][0] /= denominator;
89 unsigned int order ()
const {
94 Dune::GeometryType
type ()
const {
return gt; }
98 class Pk1dLocalCoefficients
101 Pk1dLocalCoefficients (std::size_t k_) : k(k_), n(k_+1), li(k_+1) {
102 li[0] = Dune::LocalKey(0,1,0);
103 for (
int i=1; i<int(k); i++) li[i] = Dune::LocalKey(0,0,i-1);
104 li[k] = Dune::LocalKey(1,1,0);
108 std::size_t size ()
const {
return n; }
111 const Dune::LocalKey& localKey (
int i)
const {
118 std::vector<Dune::LocalKey> li;
122 template<
typename LB>
123 class Pk1dLocalInterpolation
126 Pk1dLocalInterpolation (std::size_t k_) : k(k_), n(k_+1) {}
129 template<
typename F,
typename C>
130 void interpolate (
const F& f, std::vector<C>& out)
const
133 typename LB::Traits::DomainType x;
134 typename LB::Traits::RangeType y;
135 for (
int i=0; i<=int(k); i++)
147 Dune::GeometryType gt;
148 Pk1dLocalBasis basis;
149 Pk1dLocalCoefficients coefficients;
150 Pk1dLocalInterpolation<Pk1dLocalBasis> interpolation;
153 typedef Dune::LocalFiniteElementTraits<Pk1dLocalBasis,
154 Pk1dLocalCoefficients,
155 Pk1dLocalInterpolation<Pk1dLocalBasis> >
Traits;
158 : gt(
Dune::GeometryType::cube,1), basis(k), coefficients(k), interpolation(k)
173 return interpolation;
176 Dune::GeometryType
type ()
const {
return gt; }
190 template<
class D,
class R>
205 std::size_t
size(GeometryType gt)
const
208 return _k > 0 ? 1 : 0;
210 return _k > 0 ? _k - 1 : 1;
220 const std::size_t _k;
FiniteElementMap for the Pk basis in 1d.
Definition: pk1dbasis.hh:191
Pk1dLocalFiniteElementMap(std::size_t k)
Definition: pk1dbasis.hh:195
Pk1dLocalFiniteElement * clone() const
Definition: pk1dbasis.hh:178
Dune::LocalFiniteElementTraits< Pk1dLocalBasis, Pk1dLocalCoefficients, Pk1dLocalInterpolation< Pk1dLocalBasis > > Traits
Definition: pk1dbasis.hh:155
std::size_t maxLocalSize() const
Definition: pk1dbasis.hh:214
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: pk1dbasis.hh:166
Dune::GeometryType type() const
Definition: pk1dbasis.hh:176
const Traits::LocalBasisType & localBasis() const
Definition: pk1dbasis.hh:161
std::size_t size(GeometryType gt) const
Definition: pk1dbasis.hh:205
Define the Pk Lagrange basis functions in 1d on the reference interval.
Definition: pk1dbasis.hh:28
bool fixedSize() const
Definition: pk1dbasis.hh:200
Pk1dLocalFiniteElement(std::size_t k)
Definition: pk1dbasis.hh:157
Definition: adaptivity.hh:27
simple implementation where all entities have the same finite element
Definition: finiteelementmap.hh:107
const Traits::LocalInterpolationType & localInterpolation() const
Definition: pk1dbasis.hh:171
const std::string s
Definition: function.hh:1102
void interpolate(const F &f, const GFS &gfs, XG &xg)
interpolation from a given grid function
Definition: interpolate.hh:191