3 #ifndef DUNE_P1_LOCALBASIS_HH 4 #define DUNE_P1_LOCALBASIS_HH 9 #include <dune/common/deprecated.hh> 10 #include <dune/common/fmatrix.hh> 27 template<
class D,
class R,
int dim>
33 Dune::FieldMatrix<R,1,dim>, 2>
Traits;
43 std::vector<typename Traits::RangeType>& out)
const 47 for (
size_t i=0; i<dim; i++) {
56 std::vector<typename Traits::JacobianType>& out)
const 60 for (
int i=0; i<dim; i++)
63 for (
int i=0; i<dim; i++)
64 for (
int j=0; j<dim; j++)
65 out[i+1][0][j] = (i==j);
76 std::vector<typename Traits::RangeType>& out)
const 78 auto totalOrder = std::accumulate(order.begin(), order.end(), 0);
82 else if (totalOrder==1)
84 auto direction = std::find(order.begin(), order.end(), 1);
88 for (
int i=0; i<dim; i++)
89 out[i+1] = (i==(direction-order.begin()));
95 for (
int i=0; i<dim+1; i++)
101 template<
unsigned int k>
102 inline void DUNE_DEPRECATED_MSG(
"Use method 'partial' instead!")
104 const typename Traits::DomainType& in,
105 std::vector<typename Traits::RangeType>& out)
const 114 for (
int i=0; i<dim; i++)
115 out[i+1] = (i==directions[0]);
121 for (
int i=0; i<dim+1; i++)
unsigned int order() const
Polynomial order of the shape functions.
Definition: p1localbasis.hh:127
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:15
Linear Lagrange shape functions on the simplex.
Definition: p1localbasis.hh:28
unsigned int size() const
number of shape functions
Definition: p1localbasis.hh:36
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: p1localbasis.hh:42
D DomainType
domain type
Definition: localbasis.hh:49
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: p1localbasis.hh:55
LocalBasisTraits< D, dim, Dune::FieldVector< D, dim >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, dim >, 2 > Traits
export type traits for function signature
Definition: p1localbasis.hh:33
void partial(const std::array< unsigned int, dim > &order, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate partial derivatives of any order of all shape functions.
Definition: p1localbasis.hh:74
void evaluate(const typename std::array< int, k > &directions, const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: p1localbasis.hh:103
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37