3 #ifndef DUNE_Q2_LOCALBASIS_HH
4 #define DUNE_Q2_LOCALBASIS_HH
6 #include <dune/common/fmatrix.hh>
23 template<
class D,
class R,
int dim>
34 for (
int i=0; i<dim; i++)
41 std::vector<typename Traits::RangeType>& out)
const
46 array<array<R,3>, dim> X;
48 for (
size_t i=0; i<dim; i++) {
49 X[i][0] = R(2)*in[i]*in[i] - R(3)*in[i]+R(1);
50 X[i][1] = -R(4)*in[i]*in[i] + R(4)*in[i];
51 X[i][2] = R(2)*in[i]*in[i] - in[i];
56 for (
size_t i=0; i<out.size(); i++) {
62 for (
int j=0; j<dim; j++) {
64 int digit = ternary%3;
68 out[i] *= X[j][digit];
79 std::vector<typename Traits::JacobianType>& out)
const
84 array<array<R,3>, dim> X, DX;
86 for (
size_t i=0; i<dim; i++) {
87 X[i][0] = R(2)*in[i]*in[i] - R(3)*in[i]+R(1);
88 X[i][1] = -R(4)*in[i]*in[i] + R(4)*in[i];
89 X[i][2] = R(2)*in[i]*in[i] - in[i];
91 DX[i][0] = R(4)*in[i] - R(3);
92 DX[i][1] = -R(8)*in[i] + R(4);
93 DX[i][2] = R(4)*in[i] - R(1);
98 for (
size_t i=0; i<out.size(); i++) {
101 for (
int j=0; j<dim; j++) {
107 for (
int k=0; k<dim; k++) {
109 int digit = ternary%3;
112 out[i][0][j] *= (k==j) ? DX[k][digit] : X[k][digit];
unsigned int size() const
number of shape functions
Definition: q2localbasis.hh:31
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:39
LocalBasisTraits< D, dim, Dune::FieldVector< D, dim >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, dim > > Traits
Definition: q2localbasis.hh:28
unsigned int order() const
Polynomial order of the shape functions.
Definition: q2localbasis.hh:125
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: q2localbasis.hh:78
D DomainType
domain type
Definition: localbasis.hh:51
Lagrange shape functions of order 2 on the reference cube.
Definition: q2localbasis.hh:24
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: q2localbasis.hh:40