dune-localfunctions  2.2.1
raviartthomas12dlocalbasis.hh
Go to the documentation of this file.
1 #ifndef DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS12DLOCALBASIS_HH
2 #define DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS12DLOCALBASIS_HH
3 
4 #include <vector>
5 
6 #include <dune/common/fmatrix.hh>
7 
8 #include "../../common/localbasis.hh"
9 
10 namespace Dune
11 {
21  template<class D, class R>
23  {
24 
25 public:
26  typedef LocalBasisTraits<D,2,Dune::FieldVector<D,2>,R,2,Dune::FieldVector<R,2>,
27  Dune::FieldMatrix<R,2,2> > Traits;
30  {
31  sign0 = sign1 = sign2 = 1.0;
32  }
33 
39  RT12DLocalBasis (unsigned int s)
40  {
41  sign0 = sign1 = sign2 = 1.0;
42  if (s & 1)
43  {
44  sign0 = -1.0;
45  }
46  if (s & 2)
47  {
48  sign1 = -1.0;
49  }
50  if (s & 4)
51  {
52  sign2 = -1.0;
53  }
54  }
55 
57  unsigned int size () const
58  {
59  return 8;
60  }
61 
68  inline void evaluateFunction (const typename Traits::DomainType& in,
69  std::vector<typename Traits::RangeType>& out) const
70  {
71  out.resize(8);
72  out[0][0] = sign0*(in[0] - 4.0*in[0]*in[1]);
73  out[0][1] = sign0*(-1.0 + 5.0*in[1] - 4.0*in[1]*in[1]);
74  out[1][0] = sign1*(-1.0 + 5.0*in[0] - 4.0*in[0]*in[0]);
75  out[1][1] = sign1*(in[1] - 4.0*in[0]*in[1]);
76  out[2][0] = sign2*(-3.0*in[0] + 4.0*in[0]*in[0] + 4.0*in[1]*in[0]);
77  out[2][1] = sign2*(-3.0*in[1] + 4.0*in[0]*in[1] + 4.0*in[1]*in[1]);
78  out[3][0] = -5.0*in[0] + 8.0*in[0]*in[0] + 4.0*in[1]*in[0];
79  out[3][1] = 3.0 - 6.0*in[0] - 7.0*in[1] + 8.0*in[0]*in[1] + 4.0*in[1]*in[1];
80  out[4][0] = -3.0 + 7.0*in[0] + 6.0*in[1] - 4.0*in[0]*in[0] - 8.0*in[1]*in[0];
81  out[4][1] = 5.0*in[1] - 4.0*in[0]*in[1] - 8.0*in[1]*in[1];
82  out[5][0] = in[0] - 4.0*in[0]*in[0] + 4.0*in[1]*in[0];
83  out[5][1] = -1.0*in[1] - 4.0*in[0]*in[1] + 4.0*in[1]*in[1];
84  out[6][0] = 16.0*in[0] - 16.0*in[0]*in[0] - 8.0*in[1]*in[0];
85  out[6][1] = 8.0*in[1] - 16.0*in[0]*in[1] - 8.0*in[1]*in[1];
86  out[7][0] = 8.0*in[0] - 8.0*in[0]*in[0] - 16.0*in[1]*in[0];
87  out[7][1] = 16.0*in[1] - 8.0*in[0]*in[1] - 16.0*in[1]*in[1];
88  }
89 
96  inline void evaluateJacobian (const typename Traits::DomainType& in,
97  std::vector<typename Traits::JacobianType>& out) const
98  {
99  out.resize(8);
100 
101  out[0][0][0] = sign0*(1.0 - 4.0*in[1]);
102  out[0][0][1] = sign0*(-4.0*in[0]);
103  out[0][1][0] = 0.0;
104  out[0][1][1] = sign0*(5.0 - 8.0*in[1]);
105 
106  out[1][0][0] = sign1*(5.0 - 8.0*in[0]);
107  out[1][0][1] = 0.0;
108  out[1][1][0] = sign1*(-4.0*in[1]);
109  out[1][1][1] = sign1*(1.0 - 4.0*in[0]);
110 
111  out[2][0][0] = sign2*(-3.0 + 8.0*in[0] + 4.0*in[1]);
112  out[2][0][1] = sign2*(4.0*in[0]);
113  out[2][1][0] = sign2*(4.0*in[1]);
114  out[2][1][1] = sign2*(-3.0 + 4.0*in[0] + 8.0*in[1]);
115 
116  out[3][0][0] = -5.0 + 16.0*in[0] + 4.0*in[1];
117  out[3][0][1] = 4.0*in[0];
118  out[3][1][0] = -6.0 + 8.0*in[1];
119  out[3][1][1] = -7.0 + 8.0*in[0] + 8.0*in[1];
120 
121  out[4][0][0] = 7.0 - 8.0*in[0] - 8.0*in[1];
122  out[4][0][1] = 6.0 - 8.0*in[0];
123  out[4][1][0] = -4.0*in[1];
124  out[4][1][1] = 5.0 - 4.0*in[0] - 16.0*in[1];
125 
126  out[5][0][0] = 1.0 - 8.0*in[0] + 4*in[1];
127  out[5][0][1] = 4.0*in[0];
128  out[5][1][0] = -4.0*in[1];
129  out[5][1][1] = -1.0 - 4.0*in[0] + 8.0*in[1];
130 
131  out[6][0][0] = 16.0 - 32.0*in[0] - 8.0*in[1];
132  out[6][0][1] = -8.0*in[0];
133  out[6][1][0] = -16.0*in[1];
134  out[6][1][1] = 8.0 - 16.0*in[0] - 16.0*in[1];
135 
136  out[7][0][0] = 8.0 - 16.0*in[0] - 16.0*in[1];
137  out[7][0][1] = -16.0*in[0];
138  out[7][1][0] = -8.0*in[1];
139  out[7][1][1] = 16.0 - 8.0*in[0] - 32.0*in[1];
140  }
141 
143  unsigned int order () const
144  {
145  return 2;
146  }
147 
148 private:
149  R sign0, sign1, sign2;
150  };
151 }
152 #endif // DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS12DLOCALBASIS_HH