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>
28 template<
class D,
class R>
34 Dune::GeometryType gt;
39 typedef Dune::LocalBasisTraits<D,1,Dune::FieldVector<D,1>,R,1,Dune::FieldVector<R,1>,Dune::FieldMatrix<R,1,1>, 1>
Traits;
42 Pk1dLocalBasis (std::size_t k_) : gt(Dune::GeometryType::cube,1), k(k_), n(k_+1),
s(n)
44 for (std::size_t i=0; i<=k; i++)
s[i] = (1.0*i)/k;
48 std::size_t size ()
const {
return n; }
51 inline void evaluateFunction (
const typename Traits::DomainType& in,
52 std::vector<typename Traits::RangeType>& out)
const {
54 for (std::size_t i=0; i<=k; i++)
57 for (std::size_t j=0; j<=k; j++)
58 if (i!=j) out[i] *= (in[0]-
s[j])/(
s[i]-
s[j]);
64 evaluateJacobian (
const typename Traits::DomainType& in,
65 std::vector<typename Traits::JacobianType>& out)
const {
67 for (std::size_t i=0; i<=k; i++)
72 for (std::size_t j=0; j<=k; j++)
75 denominator *=
s[i]-
s[j];
77 for (std::size_t l=j+1; l<=k; l++)
82 out[i][0][0] += factor*a;
85 out[i][0][0] /= denominator;
90 unsigned int order ()
const {
95 Dune::GeometryType
type ()
const {
return gt; }
99 class Pk1dLocalCoefficients
102 Pk1dLocalCoefficients (std::size_t k_) : k(k_), n(k_+1), li(k_+1) {
103 li[0] = Dune::LocalKey(0,1,0);
104 for (
int i=1; i<int(k); i++) li[i] = Dune::LocalKey(0,0,i-1);
105 li[k] = Dune::LocalKey(1,1,0);
109 std::size_t size ()
const {
return n; }
112 const Dune::LocalKey& localKey (
int i)
const {
119 std::vector<Dune::LocalKey> li;
123 template<
typename LB>
124 class Pk1dLocalInterpolation
127 Pk1dLocalInterpolation (std::size_t k_) : k(k_), n(k_+1) {}
130 template<
typename F,
typename C>
131 void interpolate (
const F& f, std::vector<C>& out)
const
134 typename LB::Traits::DomainType x;
135 typename LB::Traits::RangeType y;
136 for (
int i=0; i<=int(k); i++)
148 Dune::GeometryType gt;
149 Pk1dLocalBasis basis;
150 Pk1dLocalCoefficients coefficients;
151 Pk1dLocalInterpolation<Pk1dLocalBasis> interpolation;
154 typedef Dune::LocalFiniteElementTraits<Pk1dLocalBasis,
155 Pk1dLocalCoefficients,
156 Pk1dLocalInterpolation<Pk1dLocalBasis> >
Traits;
159 : gt(
Dune::GeometryType::cube,1), basis(k), coefficients(k), interpolation(k)
174 return interpolation;
177 Dune::GeometryType
type ()
const {
return gt; }
185 #endif // DUNE_PDELAB_FINITEELEMENT_PK1D_HH