3 #ifndef DUNE_DUAL_Q1_LOCALBASIS_HH 4 #define DUNE_DUAL_Q1_LOCALBASIS_HH 9 #include <dune/common/fvector.hh> 10 #include <dune/common/fmatrix.hh> 25 template<
class D,
class R,
int dim>
32 void setCoefficients(
const std::array<Dune::FieldVector<R, (1<<dim)> ,(1<<dim)>& coefficients)
34 coefficients_ = coefficients;
45 std::vector<typename Traits::RangeType>& out)
const 48 std::vector<typename Traits::RangeType> q1Values(
size());
50 for (
size_t i=0; i<
size(); i++) {
54 for (
int j=0; j<dim; j++)
56 q1Values[i] *= (i & (1<<j)) ? in[j] : 1-in[j];
62 for (
size_t i=0; i<
size(); i++)
65 for (
size_t i=0; i<
size(); i++)
66 for (
size_t j=0; j<
size(); j++)
67 out[i] += coefficients_[i][j]*q1Values[j];
75 std::vector<typename Traits::JacobianType>& out)
const 78 std::vector<typename Traits::JacobianType> q1Jacs(
size());
81 for (
size_t i=0; i<
size(); i++) {
84 for (
int j=0; j<dim; j++) {
88 q1Jacs[i][0][j] = (i & (1<<j)) ? 1 : -1;
90 for (
int k=0; k<dim; k++) {
94 q1Jacs[i][0][j] *= (i & (1<<k)) ? in[k] : 1-in[k];
104 for (
size_t i=0; i<
size(); i++)
107 for (
size_t i=0; i<
size(); i++)
108 for (
size_t j=0; j<
size(); j++)
109 out[i].axpy(coefficients_[i][j],q1Jacs[j]);
116 std::vector<typename Traits::RangeType>& out)
const 118 auto totalOrder = std::accumulate(
order.begin(),
order.end(), 0);
119 if (totalOrder == 0) {
122 DUNE_THROW(NotImplemented,
"Desired derivative order is not implemented");
133 std::array<Dune::FieldVector<R, (1<<dim)> ,(1<<dim)> coefficients_;
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: dualq1localbasis.hh:74
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: dualq1localbasis.hh:44
void setCoefficients(const std::array< Dune::FieldVector< R,(1<< dim)>,(1<< dim)> &coefficients)
Definition: dualq1localbasis.hh:32
Dual Lagrange shape functions of order 1 on the reference cube.
Definition: dualq1localbasis.hh:26
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:31
unsigned int order() const
Polynomial order of the shape functions.
Definition: dualq1localbasis.hh:127
D DomainType
domain type
Definition: localbasis.hh:43
unsigned int size() const
number of shape functions
Definition: dualq1localbasis.hh:38
LocalBasisTraits< D, dim, Dune::FieldVector< D, dim >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, dim > > Traits
Definition: dualq1localbasis.hh:30
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:15
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: dualq1localbasis.hh:114