dune-localfunctions  2.2.1
brezzidouglasmarini1q2dlocalbasis.hh
Go to the documentation of this file.
1 #ifndef DUNE_LOCALFUNCTIONS_BREZZIDOUGLASMARINI1Q2DLOCALBASIS_HH
2 #define DUNE_LOCALFUNCTIONS_BREZZIDOUGLASMARINI1Q2DLOCALBASIS_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;
28 
31  {
32  sign0 = sign1 = sign2 = sign3 = 1.0;
33  }
34 
40  BDM1Q2DLocalBasis (unsigned int s)
41  {
42  sign0 = sign1 = sign2 = sign3 = 1.0;
43  if (s & 1)
44  {
45  sign0 = -1.0;
46  }
47  if (s & 2)
48  {
49  sign1 = -1.0;
50  }
51  if (s & 4)
52  {
53  sign2 = -1.0;
54  }
55  if (s & 8)
56  {
57  sign3 = -1.0;
58  }
59  }
60 
62  unsigned int size () const
63  {
64  return 8;
65  }
66 
73  inline void evaluateFunction (const typename Traits::DomainType& in,
74  std::vector<typename Traits::RangeType>& out) const
75  {
76  out.resize(8);
77 
78  out[0][0] = sign0*(in[0] - 1.0);
79  out[0][1] = 0.0;
80  out[1][0] = 6.0*in[0]*in[1] - 3.0*in[0]-6*in[1] + 3.0;
81  out[1][1] = -3.0*in[1]*in[1] + 3.0*in[1];
82  out[2][0] = sign1*(in[0]);
83  out[2][1] = 0.0;
84  out[3][0] = -6.0*in[0]*in[1] + 3.0*in[0];
85  out[3][1] = 3.0*in[1]*in[1] - 3.0*in[1];
86  out[4][0] = 0.0;
87  out[4][1] = sign2*(in[1] - 1.0);
88  out[5][0] = 3.0*in[0]*in[0] - 3.0*in[0];
89  out[5][1] = -6.0*in[0]*in[1] + 6.0*in[0] + 3.0*in[1] - 3.0;
90  out[6][0] = 0.0;
91  out[6][1] = sign3*(in[1]);
92  out[7][0] = -3.0*in[0]*in[0] + 3.0*in[0];
93  out[7][1] = 6.0*in[0]*in[1] - 3.0*in[1];
94  }
95 
102  inline void evaluateJacobian (const typename Traits::DomainType& in,
103  std::vector<typename Traits::JacobianType>& out) const
104  {
105  out.resize(8);
106 
107  out[0][0][0] = sign0;
108  out[0][0][1] = 0.0;
109  out[0][1][0] = 0.0;
110  out[0][1][1] = 0.0;
111 
112  out[1][0][0] = 6.0*in[1] - 3.0;
113  out[1][0][1] = 6.0*in[0] - 6.0;
114  out[1][1][0] = 0.0;
115  out[1][1][1] = -6.0*in[1] + 3.0;
116 
117  out[2][0][0] = sign1;
118  out[2][0][1] = 0.0;
119  out[2][1][0] = 0.0;
120  out[2][1][1] = 0.0;
121 
122  out[3][0][0] = -6.0*in[1] + 3.0;
123  out[3][0][1] = -6.0*in[0];
124  out[3][1][0] = 0.0;
125  out[3][1][1] = 6.0*in[1] - 3.0;
126 
127  out[4][0][0] = 0.0;
128  out[4][0][1] = 0.0;
129  out[4][1][0] = 0.0;
130  out[4][1][1] = sign2;
131 
132  out[5][0][0] = 6.0*in[0] - 3.0;
133  out[5][0][1] = 0.0;
134  out[5][1][0] = -6.0*in[1] + 6.0;
135  out[5][1][1] = -6.0*in[0] + 3.0;
136 
137  out[6][0][0] = 0.0;
138  out[6][0][1] = 0.0;
139  out[6][1][0] = 0.0;
140  out[6][1][1] = sign3;
141 
142  out[7][0][0] = -6.0*in[0] + 3.0;
143  out[7][0][1] = 0.0;
144  out[7][1][0] = 6.0*in[1];
145  out[7][1][1] = 6.0*in[0] - 3.0;
146  }
147 
149  unsigned int order () const
150  {
151  return 2;
152  }
153 
154 private:
155  R sign0, sign1, sign2, sign3;
156  };
157 }
158 #endif // DUNE_LOCALFUNCTIONS_BREZZIDOUGLASMARINI1Q2DLOCALBASIS_HH