dune-localfunctions  2.3.1-rc1
prismp2localbasis.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_PRISM_P2_LOCALBASIS_HH
4 #define DUNE_PRISM_P2_LOCALBASIS_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 18;
34  }
35 
37  inline void evaluateFunction (const typename Traits::DomainType& in,
38  std::vector<typename Traits::RangeType>& out) const
39  {
40  out.resize(18);
41 
42 
43  int coeff;
44  R a[2], b[2], c[2], a1d, b1d, c1d;
45 
46 
47  // lower triangle:
48  coeff= 2;
49  a[0] = 1;
50  a[1] = 0.5;
51  b[0] = -1;
52  b[1] = -1;
53  c[0] = -1;
54  c[1] = -1;
55  a1d = 1;
56  b1d = -3;
57  c1d = 2;
58  out[0] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
59 
60 
61  coeff= 2;
62  a[0] = 0;
63  a[1] = -0.5;
64  b[0] = 1;
65  b[1] = 0;
66  c[0] = 1;
67  c[1] = 0;
68  a1d = 1;
69  b1d = -3;
70  c1d = 2;
71  out[1] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
72 
73 
74  coeff= 2;
75  a[0] = 0;
76  a[1] = -0.5;
77  b[0] = 0;
78  b[1] = 1;
79  c[0] = 0;
80  c[1] = 1;
81  a1d = 1;
82  b1d = -3;
83  c1d = 2;
84  out[2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
85 
86  //upper triangle
87  coeff= 2;
88  a[0] = 1;
89  a[1] = 0.5;
90  b[0] = -1;
91  b[1] = -1;
92  c[0] = -1;
93  c[1] = -1;
94  a1d = 0;
95  b1d = -1;
96  c1d = 2;
97  out[3] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
98 
99 
100  coeff= 2;
101  a[0] = 0;
102  a[1] = -0.5;
103  b[0] = 1;
104  b[1] = 0;
105  c[0] = 1;
106  c[1] = 0;
107  a1d = 0;
108  b1d = -1;
109  c1d = 2;
110  out[4] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
111 
112 
113  coeff= 2;
114  a[0] = 0;
115  a[1] = -0.5;
116  b[0] = 0;
117  b[1] = 1;
118  c[0] = 0;
119  c[1] = 1;
120  a1d = 0;
121  b1d = -1;
122  c1d = 2;
123  out[5] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
124 
125  // vertical edges
126  coeff= 2;
127  a[0] = 1;
128  a[1] = 0.5;
129  b[0] = -1;
130  b[1] = -1;
131  c[0] = -1;
132  c[1] = -1;
133  a1d = 0;
134  b1d = 4;
135  c1d = -4;
136  out[6] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
137 
138 
139  coeff= 2;
140  a[0] = 0;
141  a[1] = -0.5;
142  b[0] = 1;
143  b[1] = 0;
144  c[0] = 1;
145  c[1] = 0;
146  a1d = 0;
147  b1d = 4;
148  c1d = -4;
149  out[7] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
150 
151 
152  coeff= 2;
153  a[0] = 0;
154  a[1] = -0.5;
155  b[0] = 0;
156  b[1] = 1;
157  c[0] = 0;
158  c[1] = 1;
159  a1d = 0;
160  b1d = 4;
161  c1d = -4;
162  out[8] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
163 
164  // lower triangle edges
165  coeff= 4;
166  a[0] = 0;
167  a[1] = 1;
168  b[0] = 1;
169  b[1] = 0;
170  c[0] = -1;
171  c[1] = -1;
172  a1d = 1;
173  b1d = -3;
174  c1d = 2;
175  out[9] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
176 
177 
178  coeff= 4;
179  a[0] = 0;
180  a[1] = 1;
181  b[0] = 0;
182  b[1] = 1;
183  c[0] = -1;
184  c[1] = -1;
185  a1d = 1;
186  b1d = -3;
187  c1d = 2;
188  out[10] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
189 
190 
191  coeff= 4;
192  a[0] = 0;
193  a[1] = 0;
194  b[0] = 1;
195  b[1] = 0;
196  c[0] = 0;
197  c[1] = 1;
198  a1d = 1;
199  b1d = -3;
200  c1d = 2;
201  out[11] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
202 
203  // upper triangle edges
204  coeff= 4;
205  a[0] = 0;
206  a[1] = 1;
207  b[0] = 1;
208  b[1] = 0;
209  c[0] = -1;
210  c[1] = -1;
211  a1d = 0;
212  b1d = -1;
213  c1d = 2;
214  out[12] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
215 
216 
217  coeff= 4;
218  a[0] = 0;
219  a[1] = 1;
220  b[0] = 0;
221  b[1] = 1;
222  c[0] = -1;
223  c[1] = -1;
224  a1d = 0;
225  b1d = -1;
226  c1d = 2;
227  out[13] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
228 
229 
230  coeff= 4;
231  a[0] = 0;
232  a[1] = 0;
233  b[0] = 1;
234  b[1] = 0;
235  c[0] = 0;
236  c[1] = 1;
237  a1d = 0;
238  b1d = -1;
239  c1d = 2;
240  out[14] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
241 
242  // quadrilateral sides
243  coeff= 4;
244  a[0] = 0;
245  a[1] = 1;
246  b[0] = 1;
247  b[1] = 0;
248  c[0] = -1;
249  c[1] = -1;
250  a1d = 0;
251  b1d = 4;
252  c1d = -4;
253  out[15] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
254 
255 
256  coeff= 4;
257  a[0] = 0;
258  a[1] = 1;
259  b[0] = 0;
260  b[1] = 1;
261  c[0] = -1;
262  c[1] = -1;
263  a1d = 0;
264  b1d = 4;
265  c1d = -4;
266  out[16] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
267 
268 
269  coeff= 4;
270  a[0] = 0;
271  a[1] = 0;
272  b[0] = 1;
273  b[1] = 0;
274  c[0] = 0;
275  c[1] = 1;
276  a1d = 0;
277  b1d = 4;
278  c1d = -4;
279  out[17] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1])*(a1d + in[2]*b1d + in[2]*in[2]*c1d);
280  }
281 
283  inline void
284  evaluateJacobian (const typename Traits::DomainType& in, // position
285  std::vector<typename Traits::JacobianType>& out) const // return value
286  {
287  out.resize(18);
288 
289 
290 
291  int coeff;
292  R a[2], b[2], c[2], aa[2], bb[2][2], a1d, b1d, c1d;
293 
294 
295  // lower triangle:
296  coeff= 2;
297  a[0] = 1;
298  a[1] = 0.5;
299  b[0] = -1;
300  b[1] = -1;
301  c[0] = -1;
302  c[1] = -1;
303  a1d = 1;
304  b1d = -3;
305  c1d = 2;
306  aa[0] = -3;
307  aa[1] = -3;
308  bb[0][0] = 4;
309  bb[0][1] = 4;
310  bb[1][0] = 4;
311  bb[1][1] = 4;
312  out[0][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
313  out[0][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
314  out[0][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
315 
316 
317 
318  coeff= 2;
319  a[0] = 0;
320  a[1] = -0.5;
321  b[0] = 1;
322  b[1] = 0;
323  c[0] = 1;
324  c[1] = 0;
325  a1d = 1;
326  b1d = -3;
327  c1d = 2;
328  aa[0] = -1;
329  aa[1] = 0;
330  bb[0][0] = 4;
331  bb[0][1] = 0;
332  bb[1][0] = 0;
333  bb[1][1] = 0;
334  out[1][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
335  out[1][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
336  out[1][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
337 
338 
339  coeff= 2;
340  a[0] = 0;
341  a[1] = -0.5;
342  b[0] = 0;
343  b[1] = 1;
344  c[0] = 0;
345  c[1] = 1;
346  a1d = 1;
347  b1d = -3;
348  c1d = 2;
349  aa[0] = 0;
350  aa[1] = -1;
351  bb[0][0] = 0;
352  bb[0][1] = 0;
353  bb[1][0] = 0;
354  bb[1][1] = 4;
355  out[2][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
356  out[2][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
357  out[2][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
358 
359 
360  //upper triangle
361  coeff= 2;
362  a[0] = 1;
363  a[1] = 0.5;
364  b[0] = -1;
365  b[1] = -1;
366  c[0] = -1;
367  c[1] = -1;
368  a1d = 0;
369  b1d = -1;
370  c1d = 2;
371  aa[0] = -3;
372  aa[1] = -3;
373  bb[0][0] = 4;
374  bb[0][1] = 4;
375  bb[1][0] = 4;
376  bb[1][1] = 4;
377  out[3][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
378  out[3][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
379  out[3][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
380 
381 
382 
383  coeff= 2;
384  a[0] = 0;
385  a[1] = -0.5;
386  b[0] = 1;
387  b[1] = 0;
388  c[0] = 1;
389  c[1] = 0;
390  a1d = 0;
391  b1d = -1;
392  c1d = 2;
393  aa[0] = -1;
394  aa[1] = 0;
395  bb[0][0] = 4;
396  bb[0][1] = 0;
397  bb[1][0] = 0;
398  bb[1][1] = 0;
399  out[4][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
400  out[4][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
401  out[4][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
402 
403 
404 
405  coeff= 2;
406  a[0] = 0;
407  a[1] = -0.5;
408  b[0] = 0;
409  b[1] = 1;
410  c[0] = 0;
411  c[1] = 1;
412  a1d = 0;
413  b1d = -1;
414  c1d = 2;
415  aa[0] = 0;
416  aa[1] = -1;
417  bb[0][0] = 0;
418  bb[0][1] = 0;
419  bb[1][0] = 0;
420  bb[1][1] = 4;
421  out[5][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
422  out[5][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
423  out[5][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
424 
425 
426 
427  // vertical edges
428  coeff= 2;
429  a[0] = 1;
430  a[1] = 0.5;
431  b[0] = -1;
432  b[1] = -1;
433  c[0] = -1;
434  c[1] = -1;
435  a1d = 0;
436  b1d = 4;
437  c1d = -4;
438  aa[0] = -3;
439  aa[1] = -3;
440  bb[0][0] = 4;
441  bb[0][1] = 4;
442  bb[1][0] = 4;
443  bb[1][1] = 4;
444  out[6][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
445  out[6][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
446  out[6][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
447 
448 
449 
450  coeff= 2;
451  a[0] = 0;
452  a[1] = -0.5;
453  b[0] = 1;
454  b[1] = 0;
455  c[0] = 1;
456  c[1] = 0;
457  a1d = 0;
458  b1d = 4;
459  c1d = -4;
460  aa[0] = -1;
461  aa[1] = 0;
462  bb[0][0] = 4;
463  bb[0][1] = 0;
464  bb[1][0] = 0;
465  bb[1][1] = 0;
466  out[7][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
467  out[7][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
468  out[7][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
469 
470 
471 
472  coeff= 2;
473  a[0] = 0;
474  a[1] = -0.5;
475  b[0] = 0;
476  b[1] = 1;
477  c[0] = 0;
478  c[1] = 1;
479  a1d = 0;
480  b1d = 4;
481  c1d = -4;
482  aa[0] = 0;
483  aa[1] = -1;
484  bb[0][0] = 0;
485  bb[0][1] = 0;
486  bb[1][0] = 0;
487  bb[1][1] = 4;
488  out[8][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
489  out[8][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
490  out[8][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
491 
492 
493 
494 
495  // lower triangle edges
496  coeff= 4;
497  a[0] = 0;
498  a[1] = 1;
499  b[0] = 1;
500  b[1] = 0;
501  c[0] = -1;
502  c[1] = -1;
503  a1d = 1;
504  b1d = -3;
505  c1d = 2;
506  aa[0] = 4;
507  aa[1] = 0;
508  bb[0][0] = -8;
509  bb[0][1] = -4;
510  bb[1][0] = -4;
511  bb[1][1] = 0;
512  out[9][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
513  out[9][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
514  out[9][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
515 
516 
517 
518 
519  coeff= 4;
520  a[0] = 0;
521  a[1] = 1;
522  b[0] = 0;
523  b[1] = 1;
524  c[0] = -1;
525  c[1] = -1;
526  a1d = 1;
527  b1d = -3;
528  c1d = 2;
529  aa[0] = 0;
530  aa[1] = 4; //changed from zero to 4
531  bb[0][0] = 0;
532  bb[0][1] = -4;
533  bb[1][0] = -4;
534  bb[1][1] = -8;
535  out[10][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
536  out[10][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
537  out[10][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
538 
539 
540 
541  coeff= 4;
542  a[0] = 0;
543  a[1] = 0;
544  b[0] = 1;
545  b[1] = 0;
546  c[0] = 0;
547  c[1] = 1;
548  a1d = 1;
549  b1d = -3;
550  c1d = 2;
551  aa[0] = 0;
552  aa[1] = 0;
553  bb[0][0] = 0;
554  bb[0][1] = 4;
555  bb[1][0] = 4;
556  bb[1][1] = 0;
557  out[11][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
558  out[11][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
559  out[11][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
560 
561 
562 
563  // upper triangle edges
564  coeff= 4;
565  a[0] = 0;
566  a[1] = 1;
567  b[0] = 1;
568  b[1] = 0;
569  c[0] = -1;
570  c[1] = -1;
571  a1d = 0;
572  b1d = -1;
573  c1d = 2;
574  aa[0] = 4;
575  aa[1] = 0;
576  bb[0][0] = -8;
577  bb[0][1] = -4;
578  bb[1][0] = -4;
579  bb[1][1] = 0;
580  out[12][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
581  out[12][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
582  out[12][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
583 
584 
585 
586  coeff= 4;
587  a[0] = 0;
588  a[1] = 1;
589  b[0] = 0;
590  b[1] = 1;
591  c[0] = -1;
592  c[1] = -1;
593  a1d = 0;
594  b1d = -1;
595  c1d = 2;
596  aa[0] = 0;
597  aa[1] = 4;
598  bb[0][0] = 0;
599  bb[0][1] = -4;
600  bb[1][0] = -4;
601  bb[1][1] = -8;
602  out[13][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
603  out[13][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
604  out[13][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
605 
606 
607 
608  coeff= 4;
609  a[0] = 0;
610  a[1] = 0;
611  b[0] = 1;
612  b[1] = 0;
613  c[0] = 0;
614  c[1] = 1;
615  a1d = 0;
616  b1d = -1;
617  c1d = 2;
618  aa[0] = 0;
619  aa[1] = 0;
620  bb[0][0] = 0;
621  bb[0][1] = 4;
622  bb[1][0] = 4;
623  bb[1][1] = 0;
624  out[14][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
625  out[14][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
626  out[14][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
627 
628 
629 
630  // quadrilateral sides
631  coeff= 4;
632  a[0] = 0;
633  a[1] = 1;
634  b[0] = 1;
635  b[1] = 0;
636  c[0] = -1;
637  c[1] = -1;
638  a1d = 0;
639  b1d = 4;
640  c1d = -4;
641  aa[0] = 4;
642  aa[1] = 0;
643  bb[0][0] = -8;
644  bb[0][1] = -4;
645  bb[1][0] = -4;
646  bb[1][1] = 0;
647  out[15][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
648  out[15][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
649  out[15][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
650 
651 
652 
653 
654 
655  coeff= 4;
656  a[0] = 0;
657  a[1] = 1;
658  b[0] = 0;
659  b[1] = 1;
660  c[0] = -1;
661  c[1] = -1;
662  a1d = 0;
663  b1d = 4;
664  c1d = -4;
665  aa[0] = 0;
666  aa[1] = 4;
667  bb[0][0] = 0;
668  bb[0][1] = -4;
669  bb[1][0] = -4;
670  bb[1][1] = -8;
671  out[16][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
672  out[16][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
673  out[16][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
674 
675 
676 
677 
678  coeff= 4;
679  a[0] = 0;
680  a[1] = 0;
681  b[0] = 1;
682  b[1] = 0;
683  c[0] = 0;
684  c[1] = 1;
685  a1d = 0;
686  b1d = 4;
687  c1d = -4;
688  aa[0] = 0;
689  aa[1] = 0;
690  bb[0][0] = 0;
691  bb[0][1] = 4;
692  bb[1][0] = 4;
693  bb[1][1] = 0;
694  out[17][0][0] = (aa[0] + bb[0][0]*in[0] + bb[1][0]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
695  out[17][0][1] = (aa[1] + bb[0][1]*in[0] + bb[1][1]*in[1]) * (a1d + in[2]*b1d + in[2]*in[2]*c1d);
696  out[17][0][2] = coeff * (a[0] + b[0]*in[0] + b[1]*in[1]) * (a[1] + c[0]*in[0] + c[1]*in[1]) * (b1d + 2*c1d*in[2]);
697 
698  }
699 
701  unsigned int order () const
702  {
703  return 2;
704  }
705  };
706 }
707 #endif
Type traits for LocalBasisVirtualInterface.
Definition: localbasis.hh:39
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: prismp2localbasis.hh:28
void evaluateJacobian(const typename Traits::DomainType &in, std::vector< typename Traits::JacobianType > &out) const
Evaluate Jacobian of all shape functions.
Definition: prismp2localbasis.hh:284
unsigned int size() const
number of shape functions
Definition: prismp2localbasis.hh:31
void evaluateFunction(const typename Traits::DomainType &in, std::vector< typename Traits::RangeType > &out) const
Evaluate all shape functions.
Definition: prismp2localbasis.hh:37
D DomainType
domain type
Definition: localbasis.hh:51
unsigned int order() const
Polynomial order of the shape functions.
Definition: prismp2localbasis.hh:701
Quadratic Lagrange shape functions on the prism.
Definition: prismp2localbasis.hh:23