3 #ifndef DUNE_Q1_LOCALBASIS_HH 4 #define DUNE_Q1_LOCALBASIS_HH 6 #include <dune/common/fmatrix.hh> 23 template<
class D,
class R,
int dim>
38 std::vector<typename Traits::RangeType>& out)
const 42 for (
size_t i=0; i<
size(); i++) {
46 for (
int j=0; j<dim; j++)
48 out[i] *= (i & (1<<j)) ? in[j] : 1-in[j];
57 std::vector<typename Traits::JacobianType>& out)
const 62 for (
size_t i=0; i<
size(); i++) {
65 for (
int j=0; j<dim; j++) {
69 out[i][0][j] = (i & (1<<j)) ? 1 : -1;
71 for (
int k=0; k<dim; k++) {
75 out[i][0][j] *= (i & (1<<k)) ? in[k] : 1-in[k];
unsigned int order() const
Polynomial order of the shape functions.
Definition: q1localbasis.hh:86
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:37
unsigned int size() const
number of shape functions
Definition: q1localbasis.hh:31
LocalBasisTraits< D, dim, Dune::FieldVector< D, dim >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, dim > > Traits
Definition: q1localbasis.hh:28
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:14
Lagrange shape functions of order 1 on the reference cube.
Definition: q1localbasis.hh:24
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: q1localbasis.hh:56
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: q1localbasis.hh:37