6 #ifndef DUNE_PDELAB_FINITEELEMENT_PK1D_HH 7 #define DUNE_PDELAB_FINITEELEMENT_PK1D_HH 11 #include <dune/common/fmatrix.hh> 12 #include <dune/geometry/type.hh> 14 #include<dune/localfunctions/common/localbasis.hh> 15 #include<dune/localfunctions/common/localkey.hh> 16 #include<dune/localfunctions/common/localfiniteelementtraits.hh> 25 template<
class D,
class R>
31 Dune::GeometryType gt;
36 typedef Dune::LocalBasisTraits<D,1,Dune::FieldVector<D,1>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,1>, 1>
Traits;
39 Pk1dLocalBasis (std::size_t k_) : gt(Dune::GeometryType::cube,1), k(k_), n(k_+1),
s(n)
41 for (std::size_t i=0; i<=k; i++) s[i] = (1.0*i)/k;
45 std::size_t size ()
const {
return n; }
48 inline void evaluateFunction (
const typename Traits::DomainType& in,
49 std::vector<typename Traits::RangeType>& out)
const {
51 for (std::size_t i=0; i<=k; i++)
54 for (std::size_t j=0; j<=k; j++)
55 if (i!=j) out[i] *= (in[0]-s[j])/(s[i]-s[j]);
61 evaluateJacobian (
const typename Traits::DomainType& in,
62 std::vector<typename Traits::JacobianType>& out)
const {
64 for (std::size_t i=0; i<=k; i++)
69 for (std::size_t j=0; j<=k; j++)
72 denominator *= s[i]-s[j];
74 for (std::size_t l=j+1; l<=k; l++)
79 out[i][0][0] += factor*a;
82 out[i][0][0] /= denominator;
87 unsigned int order ()
const {
92 Dune::GeometryType
type ()
const {
return gt; }
96 class Pk1dLocalCoefficients
99 Pk1dLocalCoefficients (std::size_t k_) : k(k_), n(k_+1), li(k_+1) {
100 li[0] = Dune::LocalKey(0,1,0);
101 for (
int i=1; i<int(k); i++) li[i] = Dune::LocalKey(0,0,i-1);
102 li[k] = Dune::LocalKey(1,1,0);
106 std::size_t size ()
const {
return n; }
109 const Dune::LocalKey& localKey (
int i)
const {
116 std::vector<Dune::LocalKey> li;
120 template<
typename LB>
121 class Pk1dLocalInterpolation
124 Pk1dLocalInterpolation (std::size_t k_) : k(k_), n(k_+1) {}
127 template<
typename F,
typename C>
128 void interpolate (
const F& f, std::vector<C>& out)
const 131 typename LB::Traits::DomainType x;
132 typename LB::Traits::RangeType y;
133 for (
int i=0; i<=int(k); i++)
145 Dune::GeometryType gt;
146 Pk1dLocalBasis basis;
147 Pk1dLocalCoefficients coefficients;
148 Pk1dLocalInterpolation<Pk1dLocalBasis> interpolation;
151 typedef Dune::LocalFiniteElementTraits<Pk1dLocalBasis,
152 Pk1dLocalCoefficients,
153 Pk1dLocalInterpolation<Pk1dLocalBasis> >
Traits;
156 : gt(
Dune::GeometryType::cube,1), basis(k), coefficients(k), interpolation(k)
171 return interpolation;
174 Dune::GeometryType
type ()
const {
return gt; }
182 #endif // DUNE_PDELAB_FINITEELEMENT_PK1D_HH Define the Pk Lagrange basis functions in 1d on the reference interval.
Definition: pk1d.hh:26
const std::string s
Definition: function.hh:1101
void interpolate(const F &f, const GFS &gfs, XG &xg)
interpolation from a given grid function
Definition: interpolate.hh:191
const Traits::LocalCoefficientsType & localCoefficients() const
Definition: pk1d.hh:164
Dune::LocalFiniteElementTraits< Pk1dLocalBasis, Pk1dLocalCoefficients, Pk1dLocalInterpolation< Pk1dLocalBasis > > Traits
Definition: pk1d.hh:153
Pk1dLocalFiniteElement * clone() const
Definition: pk1d.hh:176
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
const Traits::LocalBasisType & localBasis() const
Definition: pk1d.hh:159
const Traits::LocalInterpolationType & localInterpolation() const
Definition: pk1d.hh:169
Pk1dLocalFiniteElement(std::size_t k)
Definition: pk1d.hh:155
Dune::GeometryType type() const
Definition: pk1d.hh:174