Line data Source code
1 : /*
2 : * Copyright (c) 2010-2015: G-CSC, Goethe University Frankfurt
3 : * Author: Andreas Vogel
4 : *
5 : * This file is part of UG4.
6 : *
7 : * UG4 is free software: you can redistribute it and/or modify it under the
8 : * terms of the GNU Lesser General Public License version 3 (as published by the
9 : * Free Software Foundation) with the following additional attribution
10 : * requirements (according to LGPL/GPL v3 §7):
11 : *
12 : * (1) The following notice must be displayed in the Appropriate Legal Notices
13 : * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
14 : *
15 : * (2) The following notice must be displayed at a prominent place in the
16 : * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
17 : *
18 : * (3) The following bibliography is recommended for citation and must be
19 : * preserved in all covered files:
20 : * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
21 : * parallel geometric multigrid solver on hierarchically distributed grids.
22 : * Computing and visualization in science 16, 4 (2013), 151-164"
23 : * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
24 : * flexible software system for simulating pde based models on high performance
25 : * computers. Computing and visualization in science 16, 4 (2013), 165-179"
26 : *
27 : * This program is distributed in the hope that it will be useful,
28 : * but WITHOUT ANY WARRANTY; without even the implied warranty of
29 : * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 : * GNU Lesser General Public License for more details.
31 : */
32 :
33 :
34 : #include "./lagrangep1.h"
35 :
36 : namespace ug{
37 :
38 : /// \cond HIDDEN_SYMBOLS
39 :
40 : ///////////////////////////////////////
41 : // ReferenceVertex
42 : ///////////////////////////////////////
43 :
44 : template<>
45 : LagrangeP1<ReferenceVertex>::shape_type
46 0 : LagrangeP1<ReferenceVertex>::
47 : shape(size_t i, const MathVector<dim>& x) const
48 : {
49 0 : return 1.0;
50 : };
51 :
52 : template<>
53 : void
54 0 : LagrangeP1<ReferenceVertex>::
55 : grad(grad_type& value, size_t i, const MathVector<dim>& x) const
56 : {
57 0 : }
58 :
59 : template<>
60 0 : bool LagrangeP1<ReferenceVertex>::
61 : position(size_t i, MathVector<dim>& value) const
62 : {
63 0 : return true;
64 : }
65 :
66 : ///////////////////////////////////////
67 : // ReferenceEdge
68 : ///////////////////////////////////////
69 :
70 : template<>
71 : LagrangeP1<ReferenceEdge>::shape_type
72 0 : LagrangeP1<ReferenceEdge>::
73 : shape(size_t i, const MathVector<dim>& x) const
74 : {
75 0 : switch(i)
76 : {
77 0 : case 0: return (1.-x[0]);
78 0 : case 1: return x[0];
79 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
80 : }
81 : };
82 :
83 : template<>
84 : void
85 0 : LagrangeP1<ReferenceEdge>::
86 : grad(grad_type& value, size_t i, const MathVector<dim>& x) const
87 : {
88 0 : switch(i)
89 : {
90 0 : case 0: value[0] = -1.0; break;
91 0 : case 1: value[0] = 1.0; break;
92 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
93 : }
94 0 : }
95 :
96 : template<>
97 0 : bool LagrangeP1<ReferenceEdge>::
98 : position(size_t i, MathVector<dim>& value) const
99 : {
100 0 : switch(i)
101 : {
102 0 : case 0: value[0] = 0.0; return true;
103 0 : case 1: value[0] = 1.0; return true;
104 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
105 : }
106 : return true;
107 : }
108 :
109 : ///////////////////////////////////////
110 : // ReferenceTriangle
111 : ///////////////////////////////////////
112 :
113 : template<>
114 : LagrangeP1<ReferenceTriangle>::shape_type
115 0 : LagrangeP1<ReferenceTriangle>::
116 : shape(size_t i, const MathVector<dim>& x) const
117 : {
118 0 : switch (i)
119 : {
120 0 : case 0: return(1.0-x[0]-x[1]);
121 0 : case 1: return(x[0]);
122 0 : case 2: return(x[1]);
123 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
124 : }
125 : };
126 :
127 : template<>
128 : void
129 0 : LagrangeP1<ReferenceTriangle>::
130 : grad(grad_type& value, size_t i, const MathVector<dim>& x) const
131 : {
132 0 : switch(i)
133 : {
134 0 : case 0: value[0] = -1.0;
135 0 : value[1] = -1.0; break;
136 0 : case 1: value[0] = 1.0;
137 0 : value[1] = 0.0; break;
138 0 : case 2: value[0] = 0.0;
139 0 : value[1] = 1.0; break;
140 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
141 : }
142 0 : }
143 :
144 : template<>
145 0 : bool LagrangeP1<ReferenceTriangle>::
146 : position(size_t i, MathVector<dim>& value) const
147 : {
148 0 : switch(i)
149 : {
150 0 : case 0: value[0] = 0.0;
151 0 : value[1] = 0.0; return true;
152 0 : case 1: value[0] = 1.0;
153 0 : value[1] = 0.0; return true;
154 0 : case 2: value[0] = 0.0;
155 0 : value[1] = 1.0; return true;
156 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
157 : }
158 : return true;
159 : }
160 :
161 : ///////////////////////////////////////
162 : // ReferenceQuadrilateral
163 : ///////////////////////////////////////
164 :
165 : template<>
166 : LagrangeP1<ReferenceQuadrilateral>::shape_type
167 0 : LagrangeP1<ReferenceQuadrilateral>::
168 : shape(size_t i, const MathVector<dim>& x) const
169 : {
170 0 : switch(i)
171 : {
172 0 : case 0: return((1.0-x[0])*(1.0-x[1]));
173 0 : case 1: return(x[0]*(1.0-x[1]));
174 0 : case 2: return(x[0]*x[1]);
175 0 : case 3: return((1.0-x[0])*x[1]);
176 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
177 : }
178 : };
179 :
180 : template<>
181 : void
182 0 : LagrangeP1<ReferenceQuadrilateral>::
183 : grad(grad_type& value, size_t i, const MathVector<dim>& x) const
184 : {
185 0 : switch(i)
186 : {
187 0 : case 0: value[0] = -1.0 + x[1];
188 0 : value[1] = -1.0 + x[0]; break;
189 0 : case 1: value[0] = 1.0 - x[1];
190 0 : value[1] = - x[0]; break;
191 0 : case 2: value[0] = x[1];
192 0 : value[1] = x[0]; break;
193 0 : case 3: value[0] = - x[1];
194 0 : value[1] = 1.0 - x[0]; break;
195 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
196 : }
197 0 : }
198 :
199 : template<>
200 : bool
201 0 : LagrangeP1<ReferenceQuadrilateral>::
202 : position(size_t i, MathVector<dim>& value) const
203 : {
204 0 : switch(i)
205 : {
206 0 : case 0: value[0] = 0.0;
207 0 : value[1] = 0.0; return true;
208 0 : case 1: value[0] = 1.0;
209 0 : value[1] = 0.0; return true;
210 0 : case 2: value[0] = 1.0;
211 0 : value[1] = 1.0; return true;
212 0 : case 3: value[0] = 0.0;
213 0 : value[1] = 1.0; return true;
214 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
215 : }
216 : return true;
217 : }
218 :
219 : ///////////////////////////////////////
220 : // ReferenceTetrahedron
221 : ///////////////////////////////////////
222 :
223 : template<>
224 : LagrangeP1<ReferenceTetrahedron>::shape_type
225 0 : LagrangeP1<ReferenceTetrahedron>::
226 : shape(size_t i, const MathVector<dim>& x) const
227 : {
228 0 : switch(i)
229 : {
230 0 : case 0: return(1.0-x[0]-x[1]-x[2]);
231 0 : case 1: return(x[0]);
232 0 : case 2: return(x[1]);
233 0 : case 3: return(x[2]);
234 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
235 : }
236 : };
237 :
238 : template<>
239 : void
240 0 : LagrangeP1<ReferenceTetrahedron>::
241 : grad(grad_type& value, size_t i, const MathVector<dim>& x) const
242 : {
243 0 : switch(i)
244 : {
245 0 : case 0: value[0] = -1.0;
246 0 : value[1] = -1.0;
247 0 : value[2] = -1.0; break;
248 0 : case 1: value[0] = 1.0;
249 0 : value[1] = 0.0;
250 0 : value[2] = 0.0; break;
251 0 : case 2: value[0] = 0.0;
252 0 : value[1] = 1.0;
253 0 : value[2] = 0.0; break;
254 0 : case 3: value[0] = 0.0;
255 0 : value[1] = 0.0;
256 0 : value[2] = 1.0; break;
257 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
258 : }
259 0 : }
260 :
261 : template<>
262 0 : bool LagrangeP1<ReferenceTetrahedron>::
263 : position(size_t i, MathVector<dim>& value) const
264 : {
265 0 : switch(i)
266 : {
267 0 : case 0: value[0] = 0.0;
268 0 : value[1] = 0.0;
269 0 : value[2] = 0.0; return true;
270 0 : case 1: value[0] = 1.0;
271 0 : value[1] = 0.0;
272 0 : value[2] = 0.0; return true;
273 0 : case 2: value[0] = 0.0;
274 0 : value[1] = 1.0;
275 0 : value[2] = 0.0; return true;
276 0 : case 3: value[0] = 0.0;
277 0 : value[1] = 0.0;
278 0 : value[2] = 1.0; return true;
279 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
280 : }
281 : return true;
282 : }
283 :
284 : ///////////////////////////////////////
285 : // ReferencePyramid
286 : ///////////////////////////////////////
287 :
288 : template<>
289 : LagrangeP1<ReferencePyramid>::shape_type
290 0 : LagrangeP1<ReferencePyramid>::
291 : shape(size_t i, const MathVector<dim>& x) const
292 : {
293 0 : switch(i)
294 : {
295 : case 0 :
296 0 : if (x[0] > x[1])
297 0 : return((1.0-x[0])*(1.0-x[1]) + x[2]*(x[1]-1.0));
298 : else
299 0 : return((1.0-x[0])*(1.0-x[1]) + x[2]*(x[0]-1.0));
300 : case 1 :
301 0 : if (x[0] > x[1])
302 0 : return(x[0]*(1.0-x[1]) - x[2]*x[1]);
303 : else
304 0 : return(x[0]*(1.0-x[1]) - x[2]*x[0]);
305 : case 2 :
306 0 : if (x[0] > x[1])
307 0 : return(x[0]*x[1] + x[2]*x[1]);
308 : else
309 0 : return(x[0]*x[1] + x[2]*x[0]);
310 : case 3 :
311 0 : if (x[0] > x[1])
312 0 : return((1.0-x[0])*x[1] - x[2]*x[1]);
313 : else
314 0 : return((1.0-x[0])*x[1] - x[2]*x[0]);
315 0 : case 4 : return(x[2]);
316 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
317 : }
318 : };
319 :
320 : template<>
321 : void
322 0 : LagrangeP1<ReferencePyramid>::
323 : grad(grad_type& value, size_t i, const MathVector<dim>& x) const
324 : {
325 0 : switch(i)
326 : {
327 : case 0:
328 0 : if (x[0] > x[1])
329 : {
330 0 : value[0] = -(1.0-x[1]);
331 0 : value[1] = -(1.0-x[0]) + x[2];
332 0 : value[2] = -(1.0-x[1]);
333 0 : break;
334 : }
335 : else
336 : {
337 0 : value[0] = -(1.0-x[1]) + x[2];
338 0 : value[1] = -(1.0-x[0]);
339 0 : value[2] = -(1.0-x[0]);
340 0 : break;
341 : }
342 : case 1:
343 0 : if (x[0] > x[1])
344 : {
345 0 : value[0] = (1.0-x[1]);
346 0 : value[1] = -x[0] - x[2];
347 0 : value[2] = -x[1];
348 0 : break;
349 : }
350 : else
351 : {
352 0 : value[0] = (1.0-x[1]) - x[2];
353 0 : value[1] = -x[0];
354 0 : value[2] = -x[0];
355 0 : break;
356 : }
357 : case 2:
358 0 : if (x[0] > x[1])
359 : {
360 0 : value[0] = x[1];
361 0 : value[1] = x[0] + x[2];
362 0 : value[2] = x[1];
363 0 : break;
364 : }
365 : else
366 : {
367 0 : value[0] = x[1] + x[2];
368 0 : value[1] = x[0];
369 0 : value[2] = x[0];
370 0 : break;
371 : }
372 : case 3:
373 0 : if (x[0] > x[1])
374 : {
375 0 : value[0] = -x[1];
376 0 : value[1] = 1.0-x[0] - x[2];
377 0 : value[2] = -x[1];
378 0 : break;
379 : }
380 : else
381 : {
382 0 : value[0] = -x[1] - x[2];
383 0 : value[1] = 1.0-x[0];
384 0 : value[2] = -x[0];
385 0 : break;
386 : }
387 : case 4:
388 0 : value[0] = 0.0;
389 0 : value[1] = 0.0;
390 0 : value[2] = 1.0;
391 0 : break;
392 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
393 : }
394 0 : }
395 :
396 : template<>
397 0 : bool LagrangeP1<ReferencePyramid>::
398 : position(size_t i, MathVector<dim>& value) const
399 : {
400 : static const DimReferenceElement<3>& refElem
401 0 : = ReferenceElementProvider::get<3>(ROID_PYRAMID);
402 :
403 0 : value = refElem.corner(i);
404 0 : return true;
405 : }
406 :
407 : ///////////////////////////////////////
408 : // ReferencePrism
409 : ///////////////////////////////////////
410 :
411 : template<>
412 : LagrangeP1<ReferencePrism>::shape_type
413 0 : LagrangeP1<ReferencePrism>::
414 : shape(size_t i, const MathVector<dim>& x) const
415 : {
416 0 : switch(i)
417 : {
418 0 : case 0 : return((1.0-x[0]-x[1])*(1.0-x[2]));
419 0 : case 1 : return(x[0]*(1.0-x[2]));
420 0 : case 2 : return(x[1]*(1.0-x[2]));
421 0 : case 3 : return((1.0-x[0]-x[1])*x[2]);
422 0 : case 4 : return(x[0]*x[2]);
423 0 : case 5 : return(x[1]*x[2]);
424 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
425 : }
426 : };
427 :
428 : template<>
429 : void
430 0 : LagrangeP1<ReferencePrism>::
431 : grad(grad_type& value, size_t i, const MathVector<dim>& x) const
432 : {
433 0 : switch(i)
434 : {
435 : case 0:
436 0 : value[0] = -(1.0-x[2]);
437 0 : value[1] = -(1.0-x[2]);
438 0 : value[2] = -(1.0-x[0]-x[1]);
439 0 : break;
440 : case 1:
441 0 : value[0] = (1.0-x[2]);
442 0 : value[1] = 0.0;
443 0 : value[2] = -x[0];
444 0 : break;
445 : case 2:
446 0 : value[0] = 0.0;
447 0 : value[1] = (1.0-x[2]);
448 0 : value[2] = -x[1];
449 0 : break;
450 : case 3:
451 0 : value[0] = -x[2];
452 0 : value[1] = -x[2];
453 0 : value[2] = 1.0-x[0]-x[1];
454 0 : break;
455 : case 4:
456 0 : value[0] = x[2];
457 0 : value[1] = 0.0;
458 0 : value[2] = x[0];
459 0 : break;
460 : case 5:
461 0 : value[0] = 0.0;
462 0 : value[1] = x[2];
463 0 : value[2] = x[1];
464 0 : break;
465 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
466 : }
467 0 : }
468 :
469 : template<>
470 0 : bool LagrangeP1<ReferencePrism>::
471 : position(size_t i, MathVector<dim>& value) const
472 : {
473 : static const DimReferenceElement<3>& refElem
474 0 : = ReferenceElementProvider::get<3>(ROID_PRISM);
475 :
476 0 : value = refElem.corner(i);
477 0 : return true;
478 : }
479 :
480 : ///////////////////////////////////////
481 : // ReferenceHexahedron
482 : ///////////////////////////////////////
483 :
484 : template<>
485 : LagrangeP1<ReferenceHexahedron>::shape_type
486 0 : LagrangeP1<ReferenceHexahedron>::
487 : shape(size_t i, const MathVector<dim>& x) const
488 : {
489 0 : switch(i)
490 : {
491 0 : case 0: return((1.0-x[0])*(1.0-x[1])*(1.0-x[2]));
492 0 : case 1: return((x[0])*(1.0-x[1])*(1.0-x[2]));
493 0 : case 2: return((x[0])*(x[1])*(1.0-x[2]));
494 0 : case 3: return((1.0-x[0])*(x[1])*(1.0-x[2]));
495 0 : case 4: return((1.0-x[0])*(1.0-x[1])*(x[2]));
496 0 : case 5: return((x[0])*(1.0-x[1])*(x[2]));
497 0 : case 6: return((x[0])*(x[1])*(x[2]));
498 0 : case 7: return((1.0-x[0])*(x[1])*(x[2]));
499 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
500 : }
501 : };
502 :
503 : template<>
504 : void
505 0 : LagrangeP1<ReferenceHexahedron>::
506 : grad(grad_type& value, size_t i, const MathVector<dim>& x) const
507 : {
508 0 : switch(i)
509 : {
510 : case 0:
511 0 : value[0] = -(1.0-x[1])*(1.0-x[2]);
512 0 : value[1] = -(1.0-x[0])*(1.0-x[2]);
513 0 : value[2] = -(1.0-x[0])*(1.0-x[1]);
514 0 : break;
515 : case 1:
516 0 : value[0] = (1.0-x[1])*(1.0-x[2]);
517 0 : value[1] = -(x[0])*(1.0-x[2]);
518 0 : value[2] = -(x[0])*(1.0-x[1]);
519 0 : break;
520 : case 2:
521 0 : value[0] = (x[1])*(1.0-x[2]);
522 0 : value[1] = (x[0])*(1.0-x[2]);
523 0 : value[2] = -x[0]*x[1];
524 0 : break;
525 : case 3:
526 0 : value[0] = -(x[1])*(1.0-x[2]);
527 0 : value[1] = (1.0-x[0])*(1.0-x[2]);
528 0 : value[2] = -(1.0-x[0])*(x[1]);
529 0 : break;
530 : case 4:
531 0 : value[0] = -(1.0-x[1])*(x[2]);
532 0 : value[1] = -(1.0-x[0])*(x[2]);
533 0 : value[2] = (1.0-x[0])*(1.0-x[1]);
534 0 : break;
535 : case 5:
536 0 : value[0] = (1.0-x[1])*x[2];
537 0 : value[1] = -(x[0])*x[2];
538 0 : value[2] = (x[0])*(1.0-x[1]);
539 0 : break;
540 : case 6:
541 0 : value[0] = (x[1])*x[2];
542 0 : value[1] = (x[0])*x[2];
543 0 : value[2] = x[0]*x[1];
544 0 : break;
545 : case 7:
546 0 : value[0] = -(x[1])*x[2];
547 0 : value[1] = (1.0-x[0])*x[2];
548 0 : value[2] = (1.0-x[0])*x[1];
549 0 : break;
550 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
551 : }
552 0 : }
553 :
554 : template<>
555 0 : bool LagrangeP1<ReferenceHexahedron>::
556 : position(size_t i, MathVector<dim>& value) const
557 : {
558 : static const DimReferenceElement<3>& refElem
559 0 : = ReferenceElementProvider::get<3>(ROID_HEXAHEDRON);
560 :
561 0 : value = refElem.corner(i);
562 0 : return true;
563 : }
564 :
565 : ///////////////////////////////////////
566 : // ReferenceOctahedron
567 : ///////////////////////////////////////
568 :
569 : template<>
570 : LagrangeP1<ReferenceOctahedron>::shape_type
571 0 : LagrangeP1<ReferenceOctahedron>::
572 : shape(size_t i, const MathVector<dim>& x) const
573 : {
574 : // the octahedral shape functions correspond to the tetrahedral ones,
575 : // but locally piecewise linear on a subdivision of the octahedron
576 : // into 4 sub-tetrahedrons.
577 0 : const number z_sgn = (x[2] < 0) ? -1.0 : 1.0;
578 :
579 0 : switch(i)
580 : {
581 : case 0 :
582 0 : if (x[2] < 0)
583 0 : return(-1.0*x[2]);
584 : else
585 : return(0.0);
586 : case 1 :
587 0 : if (x[0] > x[1])
588 0 : return(1.0-x[0]-z_sgn*x[2]);
589 : else
590 0 : return(1.0-x[1]-z_sgn*x[2]);
591 : case 2 :
592 0 : if (x[0] > x[1])
593 0 : return(x[0]-x[1]);
594 : else
595 : return(0.0);
596 : case 3 :
597 0 : if (x[0] > x[1])
598 : return(x[1]);
599 : else
600 0 : return(x[0]);
601 : case 4 :
602 0 : if (x[0] > x[1])
603 : return(0.0);
604 : else
605 0 : return(x[1]-x[0]);
606 : case 5 :
607 0 : if (x[2] < 0)
608 0 : return(0.0);
609 : else
610 : return(x[2]);
611 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
612 : }
613 : };
614 :
615 : template<>
616 : void
617 0 : LagrangeP1<ReferenceOctahedron>::
618 : grad(grad_type& value, size_t i, const MathVector<dim>& x) const
619 : {
620 : // the octahedral shape functions correspond to the tetrahedral ones,
621 : // but locally piecewise linear on a subdivision of the octahedron
622 : // into 4 sub-tetrahedrons.
623 0 : const number z_sgn = (x[2] < 0) ? -1.0 : 1.0;
624 :
625 0 : switch(i)
626 : {
627 : case 0:
628 0 : if (x[2] < 0.0)
629 : {
630 0 : value[0] = 0.0;
631 0 : value[1] = 0.0;
632 0 : value[2] = -1.0;
633 0 : break;
634 : }
635 : else
636 : {
637 0 : value[0] = 0.0;
638 0 : value[1] = 0.0;
639 0 : value[2] = 0.0;
640 0 : break;
641 : }
642 : case 1:
643 0 : if (x[0] > x[1])
644 : {
645 0 : value[0] = -1.0;
646 0 : value[1] = 0.0;
647 0 : value[2] = -z_sgn;
648 0 : break;
649 : }
650 : else
651 : {
652 0 : value[0] = 0.0;
653 0 : value[1] = -1.0;
654 0 : value[2] = -z_sgn;
655 0 : break;
656 : }
657 : case 2:
658 0 : if (x[0] > x[1])
659 : {
660 0 : value[0] = 1.0;
661 0 : value[1] = -1.0;
662 0 : value[2] = 0.0;
663 0 : break;
664 : }
665 : else
666 : {
667 0 : value[0] = 0.0;
668 0 : value[1] = 0.0;
669 0 : value[2] = 0.0;
670 0 : break;
671 : }
672 : case 3:
673 0 : if (x[0] > x[1])
674 : {
675 0 : value[0] = 0.0;
676 0 : value[1] = 1.0;
677 0 : value[2] = 0.0;
678 0 : break;
679 : }
680 : else
681 : {
682 0 : value[0] = 1.0;
683 0 : value[1] = 0.0;
684 0 : value[2] = 0.0;
685 0 : break;
686 : }
687 : case 4:
688 0 : if (x[0] > x[1])
689 : {
690 0 : value[0] = 0.0;
691 0 : value[1] = 0.0;
692 0 : value[2] = 0.0;
693 0 : break;
694 : }
695 : else
696 : {
697 0 : value[0] = -1.0;
698 0 : value[1] = 1.0;
699 0 : value[2] = 0.0;
700 0 : break;
701 : }
702 : case 5:
703 0 : if (x[2] < 0.0)
704 : {
705 0 : value[0] = 0.0;
706 0 : value[1] = 0.0;
707 0 : value[2] = 0.0;
708 0 : break;
709 : }
710 : else
711 : {
712 0 : value[0] = 0.0;
713 0 : value[1] = 0.0;
714 0 : value[2] = 1.0;
715 0 : break;
716 : }
717 0 : default: UG_THROW("LagrangeP1: Invalid shape fct index: "<<i);
718 : }
719 0 : }
720 :
721 : template<>
722 0 : bool LagrangeP1<ReferenceOctahedron>::
723 : position(size_t i, MathVector<dim>& value) const
724 : {
725 : static const DimReferenceElement<3>& refElem
726 0 : = ReferenceElementProvider::get<3>(ROID_OCTAHEDRON);
727 :
728 0 : value = refElem.corner(i);
729 0 : return true;
730 : }
731 :
732 : /// \endcond
733 :
734 : }
735 :
|
