dune-localfunctions  2.3.0
p23dlocalbasis.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_P2_3DLOCALBASIS_HH
4 #define DUNE_P2_3DLOCALBASIS_HH
5 
6 #include <dune/common/fmatrix.hh>
7 
9 
10 namespace Dune
11 {
22  template<class D, class R>
24  {
25  public:
27  typedef LocalBasisTraits<D,3,Dune::FieldVector<D,3>,R,1,Dune::FieldVector<R,1>,
28  Dune::FieldMatrix<R,1,3> > Traits;
29 
31  unsigned int size () const
32  {
33  return 10;
34  }
35 
37  inline void evaluateFunction (const typename Traits::DomainType& in,
38  std::vector<typename Traits::RangeType>& out) const
39  {
40  out.resize(10);
41 
42  int coeff;
43  R a[2], b[3], c[3];
44 
45  // case 0:
46  coeff=2;
47  a[0]=1.0;
48  a[1]=0.5;
49  b[0]=-1.0;
50  b[1]=-1.0;
51  b[2]=-1.0;
52  c[0]=-1.0;
53  c[1]=-1.0;
54  c[2]=-1.0;
55 
56  out[0] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1] + b[2]*in[2]) * (a[1] + c[0]*in[0] + c[1]*in[1] + c[2]*in[2]);
57 
58  // case 1:
59  coeff=2;
60  a[0]=0.0;
61  a[1]=-0.5;
62  b[0]=1.0;
63  b[1]=0.0;
64  b[2]=0.0;
65  c[0]=1.0;
66  c[1]=0.0;
67  c[2]=0.0;
68 
69  out[1] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1] + b[2]*in[2]) * (a[1] + c[0]*in[0] + c[1]*in[1] + c[2]*in[2]);
70 
71  // case 2:
72  coeff=2;
73  a[0]=0.0;
74  a[1]=-0.5;
75  b[0]=0.0;
76  b[1]=1.0;
77  b[2]=0.0;
78  c[0]=0.0;
79  c[1]=1.0;
80  c[2]=0.0;
81 
82  out[2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1] + b[2]*in[2]) * (a[1] + c[0]*in[0] + c[1]*in[1] + c[2]*in[2]);
83 
84  // case 3:
85  coeff=2;
86  a[0]=0.0;
87  a[1]=-0.5;
88  b[0]=0.0;
89  b[1]=0.0;
90  b[2]=1.0;
91  c[0]=0.0;
92  c[1]=0.0;
93  c[2]=1.0;
94 
95  out[3] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1] + b[2]*in[2]) * (a[1] + c[0]*in[0] + c[1]*in[1] + c[2]*in[2]);
96 
97  // case 4:
98  coeff=4;
99  a[0]=0.0;
100  a[1]=1.0;
101  b[0]=1.0;
102  b[1]=0.0;
103  b[2]=0.0;
104  c[0]=-1.0;
105  c[1]=-1.0;
106  c[2]=-1.0;
107 
108  out[4] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1] + b[2]*in[2]) * (a[1] + c[0]*in[0] + c[1]*in[1] + c[2]*in[2]);
109 
110  // case 5:
111  coeff=4;
112  a[0]=0.0;
113  a[1]=0.0;
114  b[0]=1.0;
115  b[1]=0.0;
116  b[2]=0.0;
117  c[0]=0.0;
118  c[1]=1.0;
119  c[2]=0.0;
120 
121  out[5] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1] + b[2]*in[2]) * (a[1] + c[0]*in[0] + c[1]*in[1] + c[2]*in[2]);
122 
123  // case 6:
124  coeff=4;
125  a[0]=0.0;
126  a[1]=1.0;
127  b[0]=0.0;
128  b[1]=1.0;
129  b[2]=0.0;
130  c[0]=-1.0;
131  c[1]=-1.0;
132  c[2]=-1.0;
133 
134  out[6] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1] + b[2]*in[2]) * (a[1] + c[0]*in[0] + c[1]*in[1] + c[2]*in[2]);
135 
136  // case 7:
137  coeff=4;
138  a[0]=0.0;
139  a[1]=1.0;
140  b[0]=0.0;
141  b[1]=0.0;
142  b[2]=1.0;
143  c[0]=-1.0;
144  c[1]=-1.0;
145  c[2]=-1.0;
146 
147  out[7] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1] + b[2]*in[2]) * (a[1] + c[0]*in[0] + c[1]*in[1] + c[2]*in[2]);
148 
149  // case 8:
150  coeff=4;
151  a[0]=0.0;
152  a[1]=0.0;
153  b[0]=1.0;
154  b[1]=0.0;
155  b[2]=0.0;
156  c[0]=0.0;
157  c[1]=0.0;
158  c[2]=1.0;
159 
160  out[8] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1] + b[2]*in[2]) * (a[1] + c[0]*in[0] + c[1]*in[1] + c[2]*in[2]);
161 
162  // case 9:
163  coeff=4;
164  a[0]=0.0;
165  a[1]=0.0;
166  b[0]=0.0;
167  b[1]=1.0;
168  b[2]=0.0;
169  c[0]=0.0;
170  c[1]=0.0;
171  c[2]=1.0;
172 
173  out[9] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1] + b[2]*in[2]) * (a[1] + c[0]*in[0] + c[1]*in[1] + c[2]*in[2]);
174 
175  }
176 
178  inline void
179  evaluateJacobian (const typename Traits::DomainType& in, // position
180  std::vector<typename Traits::JacobianType>& out) const // return value
181  {
182  out.resize(10);
183 
184  R aa[3][3], bb[3][3];
185  // case 0:
186  //x derivative
187  aa[0][0]=-3.0;
188  bb[0][0]=4.0;
189  bb[1][0]=4.0;
190  bb[2][0]=4.0;
191  //y derivative
192  aa[0][1]=-3.0;
193  bb[0][1]=4.0;
194  bb[1][1]=4.0;
195  bb[2][1]=4.0;
196  // z derivative
197  aa[0][2]=-3.0;
198  bb[0][2]=4.0;
199  bb[1][2]=4.0;
200  bb[2][2]=4.0;
201 
202  out[0][0][0] = aa[0][0] + bb[0][0]*in[0] + bb[1][0]*in[1] + bb[2][0]*in[2];
203  out[0][0][1] = aa[0][1] + bb[0][1]*in[0] + bb[1][1]*in[1] + bb[2][1]*in[2];
204  out[0][0][2] = aa[0][2] + bb[0][2]*in[0] + bb[1][2]*in[1] + bb[2][2]*in[2];
205 
206  //case 1:
207  //x derivative
208  aa[0][0]=-1.0;
209  bb[0][0]=4.0;
210  bb[1][0]=0.0;
211  bb[2][0]=0.0;
212  //y derivative
213  aa[0][1]=0.0;
214  bb[0][1]=0.0;
215  bb[1][1]=0.0;
216  bb[2][1]=0.0;
217  // z derivative
218  aa[0][2]=0.0;
219  bb[0][2]=0.0;
220  bb[1][2]=0.0;
221  bb[2][2]=0.0;
222 
223  out[1][0][0] = aa[0][0] + bb[0][0]*in[0] + bb[1][0]*in[1] + bb[2][0]*in[2];
224  out[1][0][1] = aa[0][1] + bb[0][1]*in[0] + bb[1][1]*in[1] + bb[2][1]*in[2];
225  out[1][0][2] = aa[0][2] + bb[0][2]*in[0] + bb[1][2]*in[1] + bb[2][2]*in[2];
226 
227  // case 2:
228  //x derivative
229  aa[0][0]=0.0;
230  bb[0][0]=0.0;
231  bb[1][0]=0.0;
232  bb[2][0]=0.0;
233  //y derivative
234  aa[0][1]=-1.0;
235  bb[0][1]=0.0;
236  bb[1][1]=4.0;
237  bb[2][1]=0.0;
238  // z derivative
239  aa[0][2]=0.0;
240  bb[0][2]=0.0;
241  bb[1][2]=0.0;
242  bb[2][2]=0.0;
243 
244  out[2][0][0] = aa[0][0] + bb[0][0]*in[0] + bb[1][0]*in[1] + bb[2][0]*in[2];
245  out[2][0][1] = aa[0][1] + bb[0][1]*in[0] + bb[1][1]*in[1] + bb[2][1]*in[2];
246  out[2][0][2] = aa[0][2] + bb[0][2]*in[0] + bb[1][2]*in[1] + bb[2][2]*in[2];
247 
248  // case 3:
249  //x derivative
250  aa[0][0]=0.0;
251  bb[0][0]=0.0;
252  bb[1][0]=0.0;
253  bb[2][0]=0.0;
254  //y derivative
255  aa[0][1]=0.0;
256  bb[0][1]=0.0;
257  bb[1][1]=0.0;
258  bb[2][1]=0.0;
259  // z derivative
260  aa[0][2]=-1.0;
261  bb[0][2]=0.0;
262  bb[1][2]=0.0;
263  bb[2][2]=4.0;
264 
265  out[3][0][0] = aa[0][0] + bb[0][0]*in[0] + bb[1][0]*in[1] + bb[2][0]*in[2];
266  out[3][0][1] = aa[0][1] + bb[0][1]*in[0] + bb[1][1]*in[1] + bb[2][1]*in[2];
267  out[3][0][2] = aa[0][2] + bb[0][2]*in[0] + bb[1][2]*in[1] + bb[2][2]*in[2];
268 
269  // case 4:
270  //x derivative
271  aa[0][0]=4.0;
272  bb[0][0]=-8.0;
273  bb[1][0]=-4.0;
274  bb[2][0]=-4.0;
275  //y derivative
276  aa[0][1]=0.0;
277  bb[0][1]=-4.0;
278  bb[1][1]=0.0;
279  bb[2][1]=0.0;
280  // z derivative
281  aa[0][2]=0.0;
282  bb[0][2]=-4.0;
283  bb[1][2]=0.0;
284  bb[2][2]=0.0;
285 
286  out[4][0][0] = aa[0][0] + bb[0][0]*in[0] + bb[1][0]*in[1] + bb[2][0]*in[2];
287  out[4][0][1] = aa[0][1] + bb[0][1]*in[0] + bb[1][1]*in[1] + bb[2][1]*in[2];
288  out[4][0][2] = aa[0][2] + bb[0][2]*in[0] + bb[1][2]*in[1] + bb[2][2]*in[2];
289 
290  // case 5:
291  //x derivative
292  aa[0][0]=0.0;
293  bb[0][0]=0.0;
294  bb[1][0]=4.0;
295  bb[2][0]=0.0;
296  //y derivative
297  aa[0][1]=0.0;
298  bb[0][1]=4.0;
299  bb[1][1]=0.0;
300  bb[2][1]=0.0;
301  // z derivative
302  aa[0][2]=0.0;
303  bb[0][2]=0.0;
304  bb[1][2]=0.0;
305  bb[2][2]=0.0;
306 
307  out[5][0][0] = aa[0][0] + bb[0][0]*in[0] + bb[1][0]*in[1] + bb[2][0]*in[2];
308  out[5][0][1] = aa[0][1] + bb[0][1]*in[0] + bb[1][1]*in[1] + bb[2][1]*in[2];
309  out[5][0][2] = aa[0][2] + bb[0][2]*in[0] + bb[1][2]*in[1] + bb[2][2]*in[2];
310 
311  // case 6:
312  //x derivative
313  aa[0][0]=0.0;
314  bb[0][0]=0.0;
315  bb[1][0]=-4.0;
316  bb[2][0]=0.0;
317  //y derivative
318  aa[0][1]=4.0;
319  bb[0][1]=-4.0;
320  bb[1][1]=-8.0;
321  bb[2][1]=-4.0;
322  // z derivative
323  aa[0][2]=0.0;
324  bb[0][2]=0.0;
325  bb[1][2]=-4.0;
326  bb[2][2]=0.0;
327 
328  out[6][0][0] = aa[0][0] + bb[0][0]*in[0] + bb[1][0]*in[1] + bb[2][0]*in[2];
329  out[6][0][1] = aa[0][1] + bb[0][1]*in[0] + bb[1][1]*in[1] + bb[2][1]*in[2];
330  out[6][0][2] = aa[0][2] + bb[0][2]*in[0] + bb[1][2]*in[1] + bb[2][2]*in[2];
331 
332  // case 7:
333  //x derivative
334  aa[0][0]=0.0;
335  bb[0][0]=0.0;
336  bb[1][0]=0.0;
337  bb[2][0]=-4.0;
338  //y derivative
339  aa[0][1]=0.0;
340  bb[0][1]=0.0;
341  bb[1][1]=0.0;
342  bb[2][1]=-4.0;
343  // z derivative
344  aa[0][2]=4.0;
345  bb[0][2]=-4.0;
346  bb[1][2]=-4.0;
347  bb[2][2]=-8.0;
348 
349  out[7][0][0] = aa[0][0] + bb[0][0]*in[0] + bb[1][0]*in[1] + bb[2][0]*in[2];
350  out[7][0][1] = aa[0][1] + bb[0][1]*in[0] + bb[1][1]*in[1] + bb[2][1]*in[2];
351  out[7][0][2] = aa[0][2] + bb[0][2]*in[0] + bb[1][2]*in[1] + bb[2][2]*in[2];
352 
353  //case 8:
354  //x derivative
355  aa[0][0]=0.0;
356  bb[0][0]=0.0;
357  bb[1][0]=0.0;
358  bb[2][0]=4.0;
359  //y derivative
360  aa[0][1]=0.0;
361  bb[0][1]=0.0;
362  bb[1][1]=0.0;
363  bb[2][1]=0.0;
364  // z derivative
365  aa[0][2]=0.0;
366  bb[0][2]=4.0;
367  bb[1][2]=0.0;
368  bb[2][2]=0.0;
369 
370  out[8][0][0] = aa[0][0] + bb[0][0]*in[0] + bb[1][0]*in[1] + bb[2][0]*in[2];
371  out[8][0][1] = aa[0][1] + bb[0][1]*in[0] + bb[1][1]*in[1] + bb[2][1]*in[2];
372  out[8][0][2] = aa[0][2] + bb[0][2]*in[0] + bb[1][2]*in[1] + bb[2][2]*in[2];
373 
374  // case 9:
375  //x derivative
376  aa[0][0]=0.0;
377  bb[0][0]=0.0;
378  bb[1][0]=0.0;
379  bb[2][0]=0.0;
380  //y derivative
381  aa[0][1]=0.0;
382  bb[0][1]=0.0;
383  bb[1][1]=0.0;
384  bb[2][1]=4.0;
385  // z derivative
386  aa[0][2]=0.0;
387  bb[0][2]=0.0;
388  bb[1][2]=4.0;
389  bb[2][2]=0.0;
390 
391  out[9][0][0] = aa[0][0] + bb[0][0]*in[0] + bb[1][0]*in[1] + bb[2][0]*in[2];
392  out[9][0][1] = aa[0][1] + bb[0][1]*in[0] + bb[1][1]*in[1] + bb[2][1]*in[2];
393  out[9][0][2] = aa[0][2] + bb[0][2]*in[0] + bb[1][2]*in[1] + bb[2][2]*in[2];
394 
395  }
396 
398  unsigned int order () const
399  {
400  return 2;
401  }
402  };
403 }
404 #endif
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:39
unsigned int order() const
Polynomial order of the shape functions.
Definition: p23dlocalbasis.hh:398
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: p23dlocalbasis.hh:37
unsigned int size() const
number of shape functions
Definition: p23dlocalbasis.hh:31
LocalBasisTraits< D, 3, Dune::FieldVector< D, 3 >, R, 1, Dune::FieldVector< R, 1 >, Dune::FieldMatrix< R, 1, 3 > > Traits
export type traits for function signature
Definition: p23dlocalbasis.hh:28
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: p23dlocalbasis.hh:179
Quadratic Lagrange shape functions on the tetrahedron.
Definition: p23dlocalbasis.hh:23
D DomainType
domain type
Definition: localbasis.hh:51