dune-localfunctions  2.2.1
q2localbasis.hh
Go to the documentation of this file.
1 // -*- tab-width: 2; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=2 sw=2 sts=2:
3 #ifndef DUNE_Q2_LOCALBASIS_HH
4 #define DUNE_Q2_LOCALBASIS_HH
5 
6 #include <dune/common/fmatrix.hh>
7 
9 
10 namespace Dune
11 {
23 template<class D, class R, int dim>
25 {
26 public:
27  typedef LocalBasisTraits<D,dim,Dune::FieldVector<D,dim>,R,1,Dune::FieldVector<R,1>,
28  Dune::FieldMatrix<R,1,dim> > Traits;
29 
31  unsigned int size () const
32  {
33  int size = 1;
34  for (int i=0; i<dim; i++)
35  size *= 3;
36  return size;
37  }
38 
40  inline void evaluateFunction (const typename Traits::DomainType& in,
41  std::vector<typename Traits::RangeType>& out) const
42  {
43  out.resize(size());
44 
45  // Evaluate the Lagrange functions
46  array<array<R,3>, dim> X;
47 
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];
52  }
53 
54 
55  // Compute function values: they are products of 1d Lagrange function values
56  for (size_t i=0; i<out.size(); i++) {
57 
58  out[i] = 1;
59 
60  // Construct the i-th Lagrange point
61  size_t ternary = i;
62  for (int j=0; j<dim; j++) {
63 
64  int digit = ternary%3;
65  ternary /= 3;
66 
67  // Multiply the 1d Lagrange shape functions together
68  out[i] *= X[j][digit];
69 
70  }
71 
72  }
73 
74  }
75 
77  inline void
78  evaluateJacobian (const typename Traits::DomainType& in, // position
79  std::vector<typename Traits::JacobianType>& out) const // return value
80  {
81  out.resize(size());
82 
83  // Evaluate the 1d Lagrange functions and their derivatives
84  array<array<R,3>, dim> X, DX;
85 
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];
90 
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);
94  }
95 
96 
97  // Compute the derivatives by deriving the products of 1d Lagrange functions
98  for (size_t i=0; i<out.size(); i++) {
99 
100  // Computing the j-th partial derivative
101  for (int j=0; j<dim; j++) {
102 
103  out[i][0][j] = 1;
104 
105  // Loop over the 'dim' terms in the product rule
106  size_t ternary = i;
107  for (int k=0; k<dim; k++) {
108 
109  int digit = ternary%3;
110  ternary /= 3;
111 
112  out[i][0][j] *= (k==j) ? DX[k][digit] : X[k][digit];
113 
114  }
115 
116  }
117 
118  }
119 
120 
121 
122  }
123 
125  unsigned int order () const
126  {
127  return 2;
128  }
129 };
130 }
131 #endif