5 #ifndef DUNE_DUAL_P1_LOCALBASIS_HH
6 #define DUNE_DUAL_P1_LOCALBASIS_HH
10 #include <dune/common/fvector.hh>
11 #include <dune/common/fmatrix.hh>
32 template<
class D,
class R,
int dim,
bool faceDualT=false>
50 std::vector<typename Traits::RangeType>& out)
const
53 std::vector<typename Traits::RangeType> p1Values(
size());
57 for (
int i=0; i<dim; i++) {
59 p1Values[i+1] = in[i];
65 for (
int i=0; i<=dim; i++) {
66 out[i] = (dim+!
faceDual)*p1Values[i];
67 for (
int j=0; j<i; j++)
68 out[i] -= p1Values[j];
70 for (
int j=i+1; j<=dim; j++)
71 out[i] -= p1Values[j];
78 std::vector<typename Traits::JacobianType>& out)
const
81 std::vector<typename Traits::JacobianType> p1Jacs(
size());
83 for (
int i=0; i<dim; i++)
86 for (
int i=0; i<dim; i++)
87 for (
int j=0; j<dim; j++)
88 p1Jacs[i+1][0][j] = (i==j);
93 for (
size_t i=0; i<=dim; i++) {
95 out[i][0].axpy(dim+!
faceDual,p1Jacs[i][0]);
97 for (
size_t j=0; j<i; j++)
98 out[i][0] -= p1Jacs[j][0];
100 for (
int j=i+1; j<=dim; j++)
101 out[i][0] -= p1Jacs[j][0];
108 std::vector<typename Traits::RangeType>& out)
const
110 auto totalOrder = std::accumulate(
order.begin(),
order.end(), 0);
111 if (totalOrder == 0) {
114 DUNE_THROW(NotImplemented,
"Desired derivative order is not implemented");
Definition: bdfmcube.hh:18
Type traits for LocalBasisVirtualInterface.
Definition: common/localbasis.hh:34
D DomainType
domain type
Definition: common/localbasis.hh:42
Dual Lagrange shape functions on the simplex.
Definition: dualp1localbasis.hh:34
unsigned int order() const
Polynomial order of the shape functions.
Definition: dualp1localbasis.hh:119
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: dualp1localbasis.hh:49
static const bool faceDual
Determines if the basis is only biorthogonal on adjacent faces.
Definition: dualp1localbasis.hh:37
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 all shape functions.
Definition: dualp1localbasis.hh:106
unsigned int size() const
number of shape functions
Definition: dualp1localbasis.hh:43
LocalBasisTraits< D, dim, Dune::FieldVector< D, dim >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, dim > > Traits
export type traits for function signature
Definition: dualp1localbasis.hh:40
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: dualp1localbasis.hh:77