3 #ifndef DUNE_Pk1DLOCALBASIS_HH
4 #define DUNE_Pk1DLOCALBASIS_HH
6 #include <dune/common/fmatrix.hh>
24 template<
class D,
class R,
unsigned int k>
37 Dune::FieldVector<D,1>,
40 Dune::FieldVector<R,1>,
41 Dune::FieldMatrix<R,1,1>
47 for (
unsigned int i=0; i<=k; i++)
48 pos[i] = (1.0*i)/std::max(k,(
unsigned int)1);
59 std::vector<typename Traits::RangeType>& out)
const
63 for (
unsigned int i=0; i<
N; i++)
66 for (
unsigned int alpha=0; alpha<i; alpha++)
67 out[i] *= (x[0]-pos[alpha])/(pos[i]-pos[alpha]);
68 for (
unsigned int gamma=i+1; gamma<=k; gamma++)
69 out[i] *= (x[0]-pos[gamma])/(pos[i]-pos[gamma]);
76 std::vector<typename Traits::JacobianType>& out)
const
80 for (
unsigned int i=0; i<=k; i++) {
85 for (
unsigned int a=0; a<i; a++)
88 for (
unsigned int alpha=0; alpha<i; alpha++)
89 product *= (alpha==a) ? 1.0/(pos[i]-pos[alpha])
90 : (x[0]-pos[alpha])/(pos[i]-pos[alpha]);
91 for (
unsigned int gamma=i+1; gamma<=k; gamma++)
92 product *= (pos[gamma]-x[0])/(pos[gamma]-pos[i]);
93 out[i][0][0] += product;
95 for (
unsigned int c=i+1; c<=k; c++)
98 for (
unsigned int alpha=0; alpha<i; alpha++)
99 product *= (x[0]-pos[alpha])/(pos[i]-pos[alpha]);
100 for (
unsigned int gamma=i+1; gamma<=k; gamma++)
101 product *= (gamma==c) ? -1.0/(pos[gamma]-pos[i])
102 : (pos[gamma]-x[0])/(pos[gamma]-pos[i]);
103 out[i][0][0] += product;
unsigned int order() const
Polynomial order of the shape functions.
Definition: pk1dlocalbasis.hh:110
void evaluateFunction(const typename Traits::DomainType &x, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: pk1dlocalbasis.hh:58
Definition: pk1dlocalbasis.hh:30
Lagrange shape functions of arbitrary order on the 1D reference triangle.
Definition: pk1dlocalbasis.hh:25
Definition: pk1dlocalbasis.hh:33
LocalBasisTraits< D, 1, Dune::FieldVector< D, 1 >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, 1 > > Traits
Definition: pk1dlocalbasis.hh:42
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
unsigned int size() const
number of shape functions
Definition: pk1dlocalbasis.hh:52
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:14
void evaluateJacobian(const typename Traits::DomainType &x, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: pk1dlocalbasis.hh:75
D DomainType
domain type
Definition: localbasis.hh:49
Pk1DLocalBasis()
Standard constructor.
Definition: pk1dlocalbasis.hh:45