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