Line data Source code
1 : // The libMesh Finite Element Library.
2 : // Copyright (C) 2002-2026 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 :
4 : // This library is free software; you can redistribute it and/or
5 : // modify it under the terms of the GNU Lesser General Public
6 : // License as published by the Free Software Foundation; either
7 : // version 2.1 of the License, or (at your option) any later version.
8 :
9 : // This library is distributed in the hope that it will be useful,
10 : // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 : // Lesser General Public License for more details.
13 :
14 : // You should have received a copy of the GNU Lesser General Public
15 : // License along with this library; if not, write to the Free Software
16 : // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 :
18 :
19 :
20 : // Local includes
21 : #include "libmesh/fe_interface.h"
22 : #include "libmesh/fe_lagrange_gll_shape_1D.h"
23 :
24 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
25 : #include "libmesh/fe_interface_macros.h"
26 : #include "libmesh/inf_fe.h"
27 : #endif
28 :
29 : #include "libmesh/elem.h"
30 : #include "libmesh/fe.h"
31 : #include "libmesh/fe_compute_data.h"
32 : #include "libmesh/dof_map.h"
33 : #include "libmesh/enum_fe_family.h"
34 : #include "libmesh/enum_order.h"
35 : #include "libmesh/enum_elem_type.h"
36 : #include "libmesh/enum_to_string.h"
37 :
38 : namespace libMesh
39 : {
40 :
41 : //------------------------------------------------------------
42 : //FEInterface class members
43 0 : FEInterface::FEInterface()
44 : {
45 0 : libmesh_error_msg("ERROR: Do not define an object of this type.");
46 : }
47 :
48 :
49 :
50 : #ifndef LIBMESH_ENABLE_INFINITE_ELEMENTS
51 :
52 : bool
53 0 : FEInterface::is_InfFE_elem(const ElemType)
54 : {
55 0 : return false;
56 : }
57 :
58 : #else
59 :
60 : bool
61 1166008072 : FEInterface::is_InfFE_elem(const ElemType et)
62 : {
63 :
64 1166008072 : switch (et)
65 : {
66 184371 : case INFEDGE2:
67 : case INFQUAD4:
68 : case INFQUAD6:
69 : case INFHEX8:
70 : case INFHEX16:
71 : case INFHEX18:
72 : case INFPRISM6:
73 : case INFPRISM12:
74 184371 : return true;
75 :
76 1165500961 : default:
77 1165500961 : return false;
78 : }
79 : }
80 :
81 : #endif //ifndef LIBMESH_ENABLE_INFINITE_ELEMENTS
82 :
83 : #define fe_family_case_func(family, dim, func_and_args, prefix, suffix) \
84 : case family: \
85 : prefix FE<dim,family>::func_and_args suffix \
86 :
87 : #define fe_family_case(family) \
88 : case family:
89 :
90 : #define fe_family_scalar_case_func(dim, func_and_args, prefix, suffix) \
91 : fe_family_case_func(CLOUGH, dim, func_and_args, prefix, suffix) \
92 : fe_family_case_func(HERMITE, dim, func_and_args, prefix, suffix) \
93 : fe_family_case_func(HIERARCHIC, dim, func_and_args, prefix, suffix) \
94 : fe_family_case_func(L2_HIERARCHIC, dim, func_and_args, prefix, suffix) \
95 : fe_family_case_func(SIDE_HIERARCHIC, dim, func_and_args, prefix, suffix) \
96 : fe_family_case_func(LAGRANGE, dim, func_and_args, prefix, suffix) \
97 : fe_family_case_func(L2_LAGRANGE, dim, func_and_args, prefix, suffix) \
98 : fe_family_case_func(L2_LAGRANGE_GLL, dim, func_and_args, prefix, suffix) \
99 : fe_family_case_func(MONOMIAL, dim, func_and_args, prefix, suffix) \
100 : fe_family_case_func(SCALAR, dim, func_and_args, prefix, suffix) \
101 : fe_family_case_func(XYZ, dim, func_and_args, prefix, suffix) \
102 : fe_family_case_func(SUBDIVISION, 2, func_and_args, \
103 : libmesh_assert_equal_to (dim, 2); prefix, suffix)
104 :
105 : #define fe_family_scalar_case() \
106 : fe_family_case(CLOUGH) \
107 : fe_family_case(HERMITE) \
108 : fe_family_case(HIERARCHIC) \
109 : fe_family_case(L2_HIERARCHIC) \
110 : fe_family_case(SIDE_HIERARCHIC) \
111 : fe_family_case(LAGRANGE) \
112 : fe_family_case(L2_LAGRANGE) \
113 : fe_family_case(L2_LAGRANGE_GLL) \
114 : fe_family_case(MONOMIAL) \
115 : fe_family_case(SCALAR) \
116 : fe_family_case(XYZ) \
117 : fe_family_case(SUBDIVISION) \
118 :
119 : #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
120 :
121 : #define fe_family_horder_case_func(dim, func_and_args, prefix, suffix) \
122 : fe_family_case_func(BERNSTEIN, dim, func_and_args, prefix, suffix) \
123 : fe_family_case_func(SZABAB, dim, func_and_args, prefix, suffix) \
124 : fe_family_case_func(RATIONAL_BERNSTEIN, dim, func_and_args, prefix, suffix) \
125 :
126 : #define fe_family_horder_case() \
127 : fe_family_case(BERNSTEIN) \
128 : fe_family_case(SZABAB) \
129 : fe_family_case(RATIONAL_BERNSTEIN) \
130 :
131 : #else
132 :
133 : #define fe_family_horder_case_func(dim, func_and_args, prefix, suffix)
134 : #define fe_family_horder_case()
135 :
136 : #endif
137 :
138 : #define fe_family_vector_case_func(dim, func_and_args, prefix, suffix) \
139 : fe_family_case_func(HIERARCHIC_VEC, dim, func_and_args, prefix, suffix) \
140 : fe_family_case_func(L2_HIERARCHIC_VEC, dim, func_and_args, prefix, suffix) \
141 : fe_family_case_func(LAGRANGE_VEC, dim, func_and_args, prefix, suffix) \
142 : fe_family_case_func(L2_LAGRANGE_VEC, dim, func_and_args, prefix, suffix) \
143 : fe_family_case_func(MONOMIAL_VEC, dim, func_and_args, prefix, suffix) \
144 : fe_family_case_func(NEDELEC_ONE, dim, func_and_args, prefix, suffix) \
145 : fe_family_case_func(RAVIART_THOMAS, dim, func_and_args, prefix, suffix) \
146 : fe_family_case_func(L2_RAVIART_THOMAS, dim, func_and_args, prefix, suffix) \
147 :
148 : #define fe_family_vector_case() \
149 : fe_family_case(HIERARCHIC_VEC) \
150 : fe_family_case(L2_HIERARCHIC_VEC) \
151 : fe_family_case(LAGRANGE_VEC) \
152 : fe_family_case(L2_LAGRANGE_VEC) \
153 : fe_family_case(MONOMIAL_VEC) \
154 : fe_family_case(NEDELEC_ONE) \
155 : fe_family_case(RAVIART_THOMAS) \
156 : fe_family_case(L2_RAVIART_THOMAS) \
157 :
158 : #define fe_family_switch(dim, func_and_args, prefix, suffix) \
159 : do { \
160 : switch (fe_t.family) \
161 : { \
162 : fe_family_scalar_case_func(dim, func_and_args, prefix, suffix) \
163 : fe_family_horder_case_func(dim, func_and_args, prefix, suffix) \
164 : default: \
165 : libmesh_error_msg("Unsupported family = " << Utility::enum_to_string(fe_t.family)); \
166 : } \
167 : } while (0)
168 :
169 : #define fe_family_with_vec_switch(dim, func_and_args, prefix, suffix) \
170 : do { \
171 : switch (fe_t.family) \
172 : { \
173 : fe_family_scalar_case_func(dim, func_and_args, prefix, suffix) \
174 : fe_family_horder_case_func(dim, func_and_args, prefix, suffix) \
175 : fe_family_vector_case_func(dim, func_and_args, prefix, suffix) \
176 : default: \
177 : libmesh_error_msg("Unsupported family = " << Utility::enum_to_string(fe_t.family)); \
178 : } \
179 : } while (0)
180 :
181 : #define fe_scalar_vec_error_switch(dim, func_and_args, prefix, suffix) \
182 : do { \
183 : switch (fe_t.family) \
184 : { \
185 : fe_family_scalar_case_func(dim, func_and_args, prefix, suffix) \
186 : fe_family_horder_case_func(dim, func_and_args, prefix, suffix) \
187 : fe_family_vector_case() \
188 : libmesh_error_msg("Error: Can only request scalar valued elements for Real FEInterface::func_and_args"); \
189 : default: \
190 : libmesh_error_msg("Unsupported family = " << Utility::enum_to_string(fe_t.family)); \
191 : } \
192 : } while (0)
193 :
194 :
195 : #define fe_vector_scalar_error_switch(dim, func_and_args, prefix, suffix) \
196 : do { \
197 : switch (fe_t.family) \
198 : { \
199 : fe_family_vector_case_func(dim, func_and_args, prefix, suffix) \
200 : fe_family_scalar_case() \
201 : fe_family_horder_case() \
202 : libmesh_error_msg("Error: Can only request vector valued elements for RealGradient FEInterface::func_and_args"); \
203 : default: \
204 : libmesh_error_msg("Unsupported family = " << Utility::enum_to_string(fe_t.family)); \
205 : } \
206 : } while (0)
207 :
208 :
209 : #define fe_switch(func_and_args) \
210 : do { \
211 : switch (dim) \
212 : { \
213 : /* 0D */ \
214 : case 0: \
215 : fe_family_switch (0, func_and_args, return, ;); \
216 : /* 1D */ \
217 : case 1: \
218 : fe_family_switch (1, func_and_args, return, ;); \
219 : /* 2D */ \
220 : case 2: \
221 : fe_family_switch (2, func_and_args, return, ;); \
222 : /* 3D */ \
223 : case 3: \
224 : fe_family_switch (3, func_and_args, return, ;); \
225 : default: \
226 : libmesh_error_msg("Invalid dim = " << dim); \
227 : } \
228 : } while (0)
229 :
230 : #define fe_with_vec_switch(func_and_args) \
231 : do { \
232 : switch (dim) \
233 : { \
234 : /* 0D */ \
235 : case 0: \
236 : fe_family_with_vec_switch (0, func_and_args, return, ;); \
237 : /* 1D */ \
238 : case 1: \
239 : fe_family_with_vec_switch (1, func_and_args, return, ;); \
240 : /* 2D */ \
241 : case 2: \
242 : fe_family_with_vec_switch (2, func_and_args, return, ;); \
243 : /* 3D */ \
244 : case 3: \
245 : fe_family_with_vec_switch (3, func_and_args, return, ;); \
246 : default: \
247 : libmesh_error_msg("Invalid dim = " << dim); \
248 : } \
249 : } while (0)
250 :
251 : #define void_fe_with_vec_switch(func_and_args) \
252 : do { \
253 : switch (dim) \
254 : { \
255 : /* 0D */ \
256 : case 0: \
257 : fe_family_with_vec_switch (0, func_and_args, , ; return;); \
258 : /* 1D */ \
259 : case 1: \
260 : fe_family_with_vec_switch (1, func_and_args, , ; return;); \
261 : /* 2D */ \
262 : case 2: \
263 : fe_family_with_vec_switch (2, func_and_args, , ; return;); \
264 : /* 3D */ \
265 : case 3: \
266 : fe_family_with_vec_switch (3, func_and_args, , ; return;); \
267 : default: \
268 : libmesh_error_msg("Invalid dim = " << dim); \
269 : } \
270 : } while (0)
271 :
272 :
273 :
274 : #ifdef LIBMESH_ENABLE_DEPRECATED
275 0 : unsigned int FEInterface::n_shape_functions(const unsigned int dim,
276 : const FEType & fe_t,
277 : const ElemType t)
278 : {
279 : libmesh_deprecated();
280 :
281 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
282 :
283 : // We no longer support calling this function with an infinite
284 : // ElemType, you must call the FeInterface::n_shape_functions()
285 : // taking an Elem* instead.
286 0 : if (is_InfFE_elem(t))
287 0 : libmesh_not_implemented();
288 :
289 : #endif
290 :
291 0 : const Order o = fe_t.order;
292 :
293 0 : fe_with_vec_switch(n_shape_functions(t, o));
294 : }
295 : #endif // LIBMESH_ENABLE_DEPRECATED
296 :
297 :
298 :
299 : unsigned int
300 8435330 : FEInterface::n_shape_functions(const FEType & fe_t,
301 : const Elem * elem,
302 : const bool add_p_level)
303 : {
304 : // dim is required by the fe_with_vec_switch macro
305 8435330 : auto dim = elem->dim();
306 :
307 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
308 : /*
309 : * Since the FEType, stored in DofMap/(some System child), has to
310 : * be the _same_ for InfFE and FE, we have to catch calls
311 : * to infinite elements through the element type.
312 : */
313 :
314 1861874 : if (is_InfFE_elem(elem->type()))
315 0 : return ifem_n_shape_functions(fe_t, elem);
316 :
317 : #endif
318 :
319 : // Account for Elem::p_level() when computing total_order
320 9198435 : auto total_order = fe_t.order + add_p_level*elem->p_level();
321 :
322 8435330 : fe_with_vec_switch(n_dofs(elem, total_order));
323 : }
324 :
325 :
326 :
327 : unsigned int
328 1663885730 : FEInterface::n_shape_functions(const FEType & fe_t,
329 : const int extra_order,
330 : const Elem * elem)
331 : {
332 : // dim is required by the fe_with_vec_switch macro
333 1663885730 : auto dim = elem->dim();
334 :
335 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
336 : /*
337 : * Since the FEType, stored in DofMap/(some System child), has to
338 : * be the _same_ for InfFE and FE, we have to catch calls
339 : * to infinite elements through the element type.
340 : */
341 :
342 307533488 : if (is_InfFE_elem(elem->type()))
343 0 : return ifem_n_shape_functions(fe_t, elem);
344 :
345 : #endif
346 :
347 : // Ignore Elem::p_level() and instead use extra_order to compute total_order.
348 228182102 : auto total_order = fe_t.order + extra_order;
349 :
350 1663885730 : fe_with_vec_switch(n_dofs(elem, total_order));
351 : }
352 :
353 :
354 :
355 : #ifdef LIBMESH_ENABLE_DEPRECATED
356 0 : unsigned int FEInterface::n_dofs(const unsigned int dim,
357 : const FEType & fe_t,
358 : const ElemType t)
359 : {
360 : libmesh_deprecated();
361 :
362 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
363 :
364 : // We no longer support calling this function with an infinite
365 : // ElemType, you must call the FeInterface::n_dofs() taking an Elem*
366 : // instead.
367 0 : if (is_InfFE_elem(t))
368 0 : libmesh_not_implemented();
369 :
370 : #endif
371 :
372 0 : const Order o = fe_t.order;
373 :
374 0 : fe_with_vec_switch(n_dofs(t, o));
375 : }
376 :
377 :
378 :
379 :
380 : unsigned int
381 0 : FEInterface::n_dofs (const unsigned int dim,
382 : const FEType & fe_t,
383 : const Elem * elem)
384 : {
385 : libmesh_deprecated();
386 :
387 0 : fe_with_vec_switch(n_dofs(elem, fe_t.order + elem->p_level()));
388 : }
389 : #endif // LIBMESH_ENABLE_DEPRECATED
390 :
391 :
392 :
393 : unsigned int
394 9100199 : FEInterface::n_dofs(const FEType & fe_t,
395 : const Elem * elem,
396 : const bool add_p_level)
397 : {
398 : // dim is required by the fe_with_vec_switch macro
399 9100199 : auto dim = elem->dim();
400 :
401 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
402 :
403 : // InfElems currently don't support p_level()
404 1589879 : if (is_InfFE_elem(elem->type()))
405 223092 : return ifem_n_dofs(fe_t, elem);
406 :
407 : #endif
408 :
409 : // Account for Elem::p_level() when computing total_order
410 8877107 : fe_with_vec_switch(n_dofs(elem, fe_t.order + add_p_level*elem->p_level()));
411 : }
412 :
413 :
414 :
415 : unsigned int
416 24751485 : FEInterface::n_dofs(const FEType & fe_t,
417 : int extra_order,
418 : const Elem * elem)
419 : {
420 : // dim is required by the fe_with_vec_switch macro
421 24751485 : auto dim = elem->dim();
422 :
423 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
424 :
425 : // InfElems currently don't support p_level()
426 24750797 : if (is_InfFE_elem(elem->type()))
427 4257 : return ifem_n_dofs(fe_t, elem);
428 :
429 : #endif
430 :
431 : // Elem::p_level() is ignored, extra_order is used instead.
432 24746463 : auto total_order = fe_t.order + extra_order;
433 :
434 24747228 : fe_with_vec_switch(n_dofs(elem, total_order));
435 : }
436 :
437 :
438 :
439 : #ifdef LIBMESH_ENABLE_DEPRECATED
440 0 : unsigned int FEInterface::n_dofs_at_node(const unsigned int dim,
441 : const FEType & fe_t,
442 : const ElemType t,
443 : const unsigned int n)
444 : {
445 : libmesh_deprecated();
446 :
447 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
448 :
449 0 : if (is_InfFE_elem(t))
450 0 : return ifem_n_dofs_at_node(dim, fe_t, t, n);
451 :
452 : #endif
453 :
454 0 : const Order o = fe_t.order;
455 :
456 0 : fe_with_vec_switch(n_dofs_at_node(t, o, n));
457 : }
458 :
459 :
460 :
461 : FEInterface::n_dofs_at_node_ptr
462 0 : FEInterface::n_dofs_at_node_function(const unsigned int dim,
463 : const FEType & fe_t)
464 : {
465 : libmesh_deprecated();
466 :
467 0 : fe_with_vec_switch(n_dofs_at_node);
468 : }
469 : #endif // LIBMESH_ENABLE_DEPRECATED
470 :
471 :
472 :
473 : FEInterface::n_dofs_at_node_ptr
474 309480063 : FEInterface::n_dofs_at_node_function(const FEType & fe_t,
475 : const Elem * elem)
476 : {
477 : // dim is required by the fe_with_vec_switch macro
478 309480063 : auto dim = elem->dim();
479 :
480 309480063 : fe_with_vec_switch(n_dofs_at_node);
481 : }
482 :
483 :
484 :
485 : unsigned int
486 192578215 : FEInterface::n_dofs_at_node(const FEType & fe_t,
487 : const Elem * elem,
488 : const unsigned int n,
489 : const bool add_p_level)
490 : {
491 : // dim is required by the fe_with_vec_switch macro
492 192578215 : auto dim = elem->dim();
493 :
494 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
495 :
496 30778274 : if (is_InfFE_elem(elem->type()))
497 99078 : return ifem_n_dofs_at_node(fe_t, elem, n);
498 :
499 : #endif
500 :
501 : // Account for Elem::p_level() when computing total_order
502 204887426 : auto total_order = fe_t.order + add_p_level*elem->p_level();
503 :
504 192479137 : fe_with_vec_switch(n_dofs_at_node(*elem, total_order, n));
505 : }
506 :
507 :
508 :
509 : unsigned int
510 159660 : FEInterface::n_dofs_at_node(const FEType & fe_t,
511 : const int extra_order,
512 : const Elem * elem,
513 : const unsigned int n)
514 : {
515 : // dim is required by the fe_with_vec_switch macro
516 159660 : auto dim = elem->dim();
517 :
518 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
519 :
520 158364 : if (is_InfFE_elem(elem->type()))
521 157932 : return ifem_n_dofs_at_node(fe_t, elem, n);
522 :
523 : #endif
524 :
525 : // Ignore Elem::p_level() and instead use extra_order to compute total_order.
526 288 : auto total_order = fe_t.order + extra_order;
527 :
528 1728 : fe_with_vec_switch(n_dofs_at_node(*elem, total_order, n));
529 : }
530 :
531 :
532 :
533 : #ifdef LIBMESH_ENABLE_DEPRECATED
534 0 : unsigned int FEInterface::n_dofs_per_elem(const unsigned int dim,
535 : const FEType & fe_t,
536 : const ElemType t)
537 : {
538 : libmesh_deprecated();
539 :
540 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
541 :
542 0 : if (is_InfFE_elem(t))
543 0 : return ifem_n_dofs_per_elem(dim, fe_t, t);
544 :
545 : #endif
546 :
547 0 : const Order o = fe_t.order;
548 :
549 0 : fe_with_vec_switch(n_dofs_per_elem(t, o));
550 : }
551 : #endif // LIBMESH_ENABLE_DEPRECATED
552 :
553 :
554 :
555 : unsigned int
556 24007961 : FEInterface::n_dofs_per_elem(const FEType & fe_t,
557 : const Elem * elem,
558 : const bool add_p_level)
559 : {
560 : // dim is required by the fe_with_vec_switch macro
561 24007961 : auto dim = elem->dim();
562 :
563 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
564 :
565 4626167 : if (is_InfFE_elem(elem->type()))
566 7932 : return ifem_n_dofs_per_elem(fe_t, elem);
567 :
568 : #endif
569 :
570 : // Account for Elem::p_level() when computing total_order
571 25764675 : auto total_order = fe_t.order + add_p_level*elem->p_level();
572 :
573 24000029 : fe_with_vec_switch(n_dofs_per_elem(*elem, total_order));
574 : }
575 :
576 :
577 :
578 : unsigned int
579 309480063 : FEInterface::n_dofs_per_elem(const FEType & fe_t,
580 : const int extra_order,
581 : const Elem * elem)
582 : {
583 : // dim is required by the fe_with_vec_switch macro
584 309480063 : auto dim = elem->dim();
585 :
586 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
587 :
588 58075789 : if (is_InfFE_elem(elem->type()))
589 14820 : return ifem_n_dofs_per_elem(fe_t, elem);
590 :
591 : #endif
592 :
593 : // Ignore Elem::p_level() and instead use extra_order to compute total_order.
594 50125122 : auto total_order = fe_t.order + extra_order;
595 :
596 309465243 : fe_with_vec_switch(n_dofs_per_elem(*elem, total_order));
597 : }
598 :
599 :
600 :
601 930161 : void FEInterface::dofs_on_side(const Elem * const elem,
602 : const unsigned int dim,
603 : const FEType & fe_t,
604 : unsigned int s,
605 : std::vector<unsigned int> & di,
606 : const bool add_p_level)
607 : {
608 157230 : const Order o = fe_t.order;
609 :
610 930161 : void_fe_with_vec_switch(dofs_on_side(elem, o, s, di, add_p_level));
611 : }
612 :
613 :
614 :
615 148298 : void FEInterface::dofs_on_edge(const Elem * const elem,
616 : const unsigned int dim,
617 : const FEType & fe_t,
618 : unsigned int e,
619 : std::vector<unsigned int> & di,
620 : const bool add_p_level)
621 : {
622 21578 : const Order o = fe_t.order;
623 :
624 148298 : void_fe_with_vec_switch(dofs_on_edge(elem, o, e, di, add_p_level));
625 : }
626 :
627 :
628 :
629 15889262 : void FEInterface::nodal_soln(const unsigned int dim,
630 : const FEType & fe_t,
631 : const Elem * elem,
632 : const std::vector<Number> & elem_soln,
633 : std::vector<Number> & nodal_soln,
634 : const bool add_p_level,
635 : const unsigned int vdim)
636 : {
637 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
638 :
639 3469759 : if (is_InfFE_elem(elem->type()))
640 : {
641 0 : ifem_nodal_soln(dim, fe_t, elem, elem_soln, nodal_soln);
642 0 : return;
643 : }
644 :
645 : #endif
646 :
647 2884090 : const Order order = fe_t.order;
648 :
649 15889262 : void_fe_with_vec_switch(nodal_soln(elem, order, elem_soln, nodal_soln, add_p_level, vdim));
650 : }
651 :
652 :
653 :
654 495912 : void FEInterface::side_nodal_soln(const FEType & fe_t,
655 : const Elem * elem,
656 : const unsigned int side,
657 : const std::vector<Number> & elem_soln,
658 : std::vector<Number> & nodal_soln,
659 : const bool add_p_level,
660 : const unsigned int vdim)
661 : {
662 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
663 :
664 123402 : if (is_InfFE_elem(elem->type()))
665 : {
666 0 : libmesh_not_implemented();
667 : return;
668 : }
669 :
670 : #endif
671 :
672 82780 : const Order order = fe_t.order;
673 495912 : const unsigned int dim = elem->dim();
674 :
675 495912 : void_fe_with_vec_switch(side_nodal_soln(elem, order, side, elem_soln, nodal_soln, add_p_level, vdim));
676 : }
677 :
678 :
679 :
680 : #ifdef LIBMESH_ENABLE_DEPRECATED
681 0 : Point FEInterface::map(unsigned int dim,
682 : const FEType & fe_t,
683 : const Elem * elem,
684 : const Point & p)
685 : {
686 : libmesh_deprecated();
687 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
688 0 : if (is_InfFE_elem(elem->type()))
689 0 : return ifem_map(dim, fe_t, elem, p);
690 : #endif
691 0 : fe_with_vec_switch(map(elem, p));
692 : }
693 :
694 :
695 :
696 0 : Point FEInterface::inverse_map (const unsigned int dim,
697 : const FEType & fe_t,
698 : const Elem * elem,
699 : const Point & p,
700 : const Real tolerance,
701 : const bool secure)
702 : {
703 : libmesh_deprecated();
704 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
705 :
706 0 : if (is_InfFE_elem(elem->type()))
707 0 : return ifem_inverse_map(dim, fe_t, elem, p,tolerance, secure);
708 :
709 : #endif
710 :
711 0 : fe_with_vec_switch(inverse_map(elem, p, tolerance, secure));
712 : }
713 :
714 :
715 :
716 0 : void FEInterface::inverse_map (const unsigned int dim,
717 : const FEType & fe_t,
718 : const Elem * elem,
719 : const std::vector<Point> & physical_points,
720 : std::vector<Point> & reference_points,
721 : const Real tolerance,
722 : const bool secure)
723 : {
724 : libmesh_deprecated();
725 :
726 0 : const std::size_t n_pts = physical_points.size();
727 :
728 : // Resize the vector
729 0 : reference_points.resize(n_pts);
730 :
731 0 : if (n_pts == 0)
732 : {
733 0 : libMesh::err << "WARNING: empty vector physical_points!"
734 0 : << std::endl;
735 0 : libmesh_here();
736 0 : return;
737 : }
738 :
739 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
740 :
741 0 : if (is_InfFE_elem(elem->type()))
742 : {
743 0 : ifem_inverse_map(dim, fe_t, elem, physical_points, reference_points, tolerance, secure);
744 0 : return;
745 : // libmesh_not_implemented();
746 : }
747 :
748 : #endif
749 :
750 0 : void_fe_with_vec_switch(inverse_map(elem, physical_points, reference_points, tolerance, secure));
751 : }
752 :
753 :
754 :
755 0 : bool FEInterface::on_reference_element(const Point & p,
756 : const ElemType t,
757 : const Real eps)
758 : {
759 0 : return FEBase::on_reference_element(p,t,eps);
760 : }
761 :
762 :
763 :
764 0 : Real FEInterface::shape(const unsigned int dim,
765 : const FEType & fe_t,
766 : const ElemType t,
767 : const unsigned int i,
768 : const Point & p)
769 : {
770 : libmesh_deprecated();
771 :
772 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
773 :
774 0 : if (is_InfFE_elem(t))
775 0 : return ifem_shape(dim, fe_t, t, i, p);
776 :
777 : #endif
778 :
779 0 : const Order o = fe_t.order;
780 :
781 0 : fe_switch(shape(t,o,i,p));
782 : }
783 :
784 0 : Real FEInterface::shape(const unsigned int dim,
785 : const FEType & fe_t,
786 : const Elem * elem,
787 : const unsigned int i,
788 : const Point & p)
789 : {
790 : libmesh_deprecated();
791 :
792 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
793 :
794 0 : if (elem && is_InfFE_elem(elem->type()))
795 0 : return ifem_shape(fe_t, elem, i, p);
796 :
797 : #endif
798 :
799 0 : const Order o = fe_t.order;
800 :
801 0 : fe_switch(shape(elem,o,i,p));
802 : }
803 : #endif // LIBMESH_ENABLE_DEPRECATED
804 :
805 :
806 :
807 : Real
808 223878891 : FEInterface::shape(const FEType & fe_t,
809 : const Elem * elem,
810 : const unsigned int i,
811 : const Point & p,
812 : const bool add_p_level)
813 : {
814 : // dim is required by the fe_switch macro
815 223878891 : auto dim = elem->dim();
816 :
817 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
818 :
819 51887423 : if (elem && is_InfFE_elem(elem->type()))
820 0 : return ifem_shape(fe_t, elem, i, p);
821 :
822 : #endif
823 :
824 : // We are calling
825 : //
826 : // FE<X,Y>::shape(Elem *, Order, unsigned, Point, true)
827 : //
828 : // See fe.h for more details.
829 223878891 : fe_switch(shape(elem, fe_t.order, i, p, add_p_level));
830 : }
831 :
832 :
833 :
834 : Real
835 1333151491 : FEInterface::shape(const FEType & fe_t,
836 : int extra_order,
837 : const Elem * elem,
838 : const unsigned int i,
839 : const Point & p)
840 : {
841 : // dim is required by the fe_switch macro
842 1333151491 : auto dim = elem->dim();
843 :
844 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
845 :
846 305113357 : if (elem && is_InfFE_elem(elem->type()))
847 0 : return ifem_shape(fe_t, elem, i, p);
848 :
849 : #endif
850 :
851 : // We are calling
852 : //
853 : // FE<X,Y>::shape(Elem *, Order, unsigned, Point, false)
854 : //
855 : // with the last parameter set to "false" so that the
856 : // Elem::p_level() is not used internally and the "total_order" that
857 : // we compute is used instead. See fe.h for more details.
858 217810650 : auto total_order = fe_t.order + extra_order;
859 :
860 1333151491 : fe_switch(shape(elem, total_order, i, p, false));
861 : }
862 :
863 :
864 :
865 : #ifdef LIBMESH_ENABLE_DEPRECATED
866 : template<>
867 0 : void FEInterface::shape<Real>(const unsigned int dim,
868 : const FEType & fe_t,
869 : const ElemType t,
870 : const unsigned int i,
871 : const Point & p,
872 : Real & phi)
873 : {
874 : libmesh_deprecated();
875 :
876 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
877 :
878 0 : if (is_InfFE_elem(t))
879 : {
880 0 : phi = ifem_shape(dim, fe_t, t, i, p);
881 0 : return;
882 : }
883 :
884 : #endif
885 :
886 0 : const Order o = fe_t.order;
887 :
888 0 : switch(dim)
889 : {
890 0 : case 0:
891 0 : fe_scalar_vec_error_switch(0, shape(t,o,i,p), phi = , ; break;);
892 0 : break;
893 0 : case 1:
894 0 : fe_scalar_vec_error_switch(1, shape(t,o,i,p), phi = , ; break;);
895 0 : break;
896 0 : case 2:
897 0 : fe_scalar_vec_error_switch(2, shape(t,o,i,p), phi = , ; break;);
898 0 : break;
899 0 : case 3:
900 0 : fe_scalar_vec_error_switch(3, shape(t,o,i,p), phi = , ; break;);
901 0 : break;
902 0 : default:
903 0 : libmesh_error_msg("Invalid dimension = " << dim);
904 : }
905 :
906 0 : return;
907 : }
908 :
909 :
910 :
911 : template<>
912 0 : void FEInterface::shape<Real>(const unsigned int dim,
913 : const FEType & fe_t,
914 : const Elem * elem,
915 : const unsigned int i,
916 : const Point & p,
917 : Real & phi)
918 : {
919 : libmesh_deprecated();
920 :
921 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
922 :
923 0 : if (elem && is_InfFE_elem(elem->type()))
924 : {
925 0 : phi = ifem_shape(fe_t, elem, i, p);
926 0 : return;
927 : }
928 : #endif
929 :
930 0 : const Order o = fe_t.order;
931 :
932 0 : switch(dim)
933 : {
934 0 : case 0:
935 0 : fe_scalar_vec_error_switch(0, shape(elem,o,i,p), phi = , ; break;);
936 0 : break;
937 0 : case 1:
938 0 : fe_scalar_vec_error_switch(1, shape(elem,o,i,p), phi = , ; break;);
939 0 : break;
940 0 : case 2:
941 0 : fe_scalar_vec_error_switch(2, shape(elem,o,i,p), phi = , ; break;);
942 0 : break;
943 0 : case 3:
944 0 : fe_scalar_vec_error_switch(3, shape(elem,o,i,p), phi = , ; break;);
945 0 : break;
946 0 : default:
947 0 : libmesh_error_msg("Invalid dimension = " << dim);
948 : }
949 :
950 0 : return;
951 : }
952 : #endif // LIBMESH_ENABLE_DEPRECATED
953 :
954 :
955 :
956 : template<>
957 0 : void FEInterface::shape<Real>(const FEType & fe_t,
958 : const Elem * elem,
959 : const unsigned int i,
960 : const Point & p,
961 : Real & phi)
962 : {
963 : // dim is required by the fe_switch macro
964 0 : auto dim = elem->dim();
965 :
966 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
967 :
968 0 : if (is_InfFE_elem(elem->type()))
969 : {
970 0 : phi = ifem_shape(fe_t, elem, i, p);
971 0 : return;
972 : }
973 :
974 : #endif
975 :
976 : // Below we call FE<X,Y>::shape(Elem *, Order, unsigned, Point, true)
977 : // so that the Elem::p_level() is accounted for internally.
978 0 : switch(dim)
979 : {
980 0 : case 0:
981 0 : fe_scalar_vec_error_switch(0, shape(elem, fe_t.order, i, p, true), phi = , ; break;);
982 0 : break;
983 0 : case 1:
984 0 : fe_scalar_vec_error_switch(1, shape(elem, fe_t.order, i, p, true), phi = , ; break;);
985 0 : break;
986 0 : case 2:
987 0 : fe_scalar_vec_error_switch(2, shape(elem, fe_t.order, i, p, true), phi = , ; break;);
988 0 : break;
989 0 : case 3:
990 0 : fe_scalar_vec_error_switch(3, shape(elem, fe_t.order, i, p, true), phi = , ; break;);
991 0 : break;
992 0 : default:
993 0 : libmesh_error_msg("Invalid dimension = " << dim);
994 : }
995 0 : }
996 :
997 :
998 :
999 : template<>
1000 0 : void FEInterface::shape<Real>(const FEType & fe_t,
1001 : int extra_order,
1002 : const Elem * elem,
1003 : const unsigned int i,
1004 : const Point & p,
1005 : Real & phi)
1006 : {
1007 : // dim is required by the fe_switch macro
1008 0 : auto dim = elem->dim();
1009 :
1010 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1011 :
1012 0 : if (is_InfFE_elem(elem->type()))
1013 : {
1014 0 : phi = ifem_shape(fe_t, elem, i, p);
1015 0 : return;
1016 : }
1017 :
1018 : #endif
1019 :
1020 : // Ignore Elem::p_level() and instead use extra_order to compute total_order
1021 0 : auto total_order = fe_t.order + extra_order;
1022 :
1023 : // Below we call
1024 : //
1025 : // FE<X,Y>::shape(Elem *, Order, unsigned, Point, false)
1026 : //
1027 : // so that the Elem::p_level() is ignored and the total_order that
1028 : // we compute is used instead.
1029 0 : switch(dim)
1030 : {
1031 0 : case 0:
1032 0 : fe_scalar_vec_error_switch(0, shape(elem, total_order, i, p, false), phi = , ; break;);
1033 0 : break;
1034 0 : case 1:
1035 0 : fe_scalar_vec_error_switch(1, shape(elem, total_order, i, p, false), phi = , ; break;);
1036 0 : break;
1037 0 : case 2:
1038 0 : fe_scalar_vec_error_switch(2, shape(elem, total_order, i, p, false), phi = , ; break;);
1039 0 : break;
1040 0 : case 3:
1041 0 : fe_scalar_vec_error_switch(3, shape(elem, total_order, i, p, false), phi = , ; break;);
1042 0 : break;
1043 0 : default:
1044 0 : libmesh_error_msg("Invalid dimension = " << dim);
1045 : }
1046 0 : }
1047 :
1048 :
1049 :
1050 : template<>
1051 5485828 : void FEInterface::shapes<Real>(const unsigned int dim,
1052 : const FEType & fe_t,
1053 : const Elem * elem,
1054 : const unsigned int i,
1055 : const std::vector<Point> & p,
1056 : std::vector<Real> & phi,
1057 : const bool add_p_level)
1058 : {
1059 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1060 :
1061 1318792 : if (elem && is_InfFE_elem(elem->type()))
1062 : {
1063 0 : FEType elevated = fe_t;
1064 0 : elevated.order = fe_t.order + add_p_level * elem->p_level();
1065 0 : for (auto qpi : index_range(p))
1066 0 : phi[qpi] = ifem_shape(elevated, elem, i, p[qpi]);
1067 0 : return;
1068 : }
1069 : #endif
1070 :
1071 937786 : const Order o = fe_t.order;
1072 :
1073 5485828 : switch(dim)
1074 : {
1075 0 : case 0:
1076 0 : fe_scalar_vec_error_switch(0, shapes(elem,o,i,p,phi,add_p_level), , ; return;);
1077 : break;
1078 4188 : case 1:
1079 4188 : fe_scalar_vec_error_switch(1, shapes(elem,o,i,p,phi,add_p_level), , ; return;);
1080 : break;
1081 2490672 : case 2:
1082 2490672 : fe_scalar_vec_error_switch(2, shapes(elem,o,i,p,phi,add_p_level), , ; return;);
1083 : break;
1084 2990968 : case 3:
1085 2990968 : fe_scalar_vec_error_switch(3, shapes(elem,o,i,p,phi,add_p_level), , ; return;);
1086 : break;
1087 0 : default:
1088 0 : libmesh_error_msg("Invalid dimension = " << dim);
1089 : }
1090 :
1091 : return;
1092 : }
1093 :
1094 :
1095 : template<>
1096 134715035 : void FEInterface::all_shapes<Real>(const unsigned int dim,
1097 : const FEType & fe_t,
1098 : const Elem * elem,
1099 : const std::vector<Point> & p,
1100 : std::vector<std::vector<Real>> & phi,
1101 : const bool add_p_level)
1102 : {
1103 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1104 :
1105 30684347 : if (elem && is_InfFE_elem(elem->type()))
1106 : {
1107 0 : for (auto i : index_range(phi))
1108 0 : FEInterface::shapes<Real>(dim, fe_t, elem, i, p, phi[i], add_p_level);
1109 0 : return;
1110 : }
1111 : #endif
1112 :
1113 24213036 : const Order o = fe_t.order;
1114 :
1115 134715035 : switch(dim)
1116 : {
1117 198780 : case 0:
1118 198780 : fe_scalar_vec_error_switch(0, all_shapes(elem,o,p,phi,add_p_level), , ; return;);
1119 : break;
1120 2822426 : case 1:
1121 2822426 : fe_scalar_vec_error_switch(1, all_shapes(elem,o,p,phi,add_p_level), , ; return;);
1122 : break;
1123 98406346 : case 2:
1124 98406346 : fe_scalar_vec_error_switch(2, all_shapes(elem,o,p,phi,add_p_level), , ; return;);
1125 : break;
1126 33287483 : case 3:
1127 33287483 : fe_scalar_vec_error_switch(3, all_shapes(elem,o,p,phi,add_p_level), , ; return;);
1128 : break;
1129 0 : default:
1130 0 : libmesh_error_msg("Invalid dimension = " << dim);
1131 : }
1132 :
1133 : return;
1134 : }
1135 :
1136 :
1137 :
1138 : #ifdef LIBMESH_ENABLE_DEPRECATED
1139 : template<>
1140 0 : void FEInterface::shape<RealGradient>(const unsigned int dim,
1141 : const FEType & fe_t,
1142 : const ElemType t,
1143 : const unsigned int i,
1144 : const Point & p,
1145 : RealGradient & phi)
1146 : {
1147 : libmesh_deprecated();
1148 :
1149 : // This API does not currently support infinite elements.
1150 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1151 0 : if (is_InfFE_elem(t))
1152 : {
1153 0 : libmesh_not_implemented();
1154 : }
1155 : #endif
1156 0 : const Order o = fe_t.order;
1157 :
1158 0 : switch(dim)
1159 : {
1160 0 : case 0:
1161 0 : fe_vector_scalar_error_switch(0, shape(t,o,i,p), phi = , ; break;);
1162 0 : break;
1163 0 : case 1:
1164 0 : fe_vector_scalar_error_switch(1, shape(t,o,i,p), phi = , ; break;);
1165 0 : break;
1166 0 : case 2:
1167 0 : fe_vector_scalar_error_switch(2, shape(t,o,i,p), phi = , ; break;);
1168 0 : break;
1169 0 : case 3:
1170 0 : fe_vector_scalar_error_switch(3, shape(t,o,i,p), phi = , ; break;);
1171 0 : break;
1172 0 : default:
1173 0 : libmesh_error_msg("Invalid dimension = " << dim);
1174 : }
1175 0 : }
1176 : #endif // LIBMESH_ENABLE_DEPRECATED
1177 :
1178 :
1179 :
1180 : template<>
1181 0 : void FEInterface::shape<RealGradient>(const FEType & fe_t,
1182 : const Elem * elem,
1183 : const unsigned int i,
1184 : const Point & p,
1185 : RealGradient & phi)
1186 : {
1187 : // This API does not currently support infinite elements.
1188 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1189 0 : if (is_InfFE_elem(elem->type()))
1190 : {
1191 0 : libmesh_not_implemented();
1192 : }
1193 : #endif
1194 :
1195 0 : auto dim = elem->dim();
1196 :
1197 : // We are calling
1198 : //
1199 : // FE<X,Y>::shape(Elem *, Order, unsigned, Point, true)
1200 : //
1201 : // with the last parameter set to "true" so that the Elem::p_level()
1202 : // is accounted for internally. See fe.h for more details.
1203 :
1204 0 : switch(dim)
1205 : {
1206 0 : case 0:
1207 0 : fe_vector_scalar_error_switch(0, shape(elem, fe_t.order, i, p, true), phi = , ; break;);
1208 0 : break;
1209 0 : case 1:
1210 0 : fe_vector_scalar_error_switch(1, shape(elem, fe_t.order, i, p, true), phi = , ; break;);
1211 0 : break;
1212 0 : case 2:
1213 0 : fe_vector_scalar_error_switch(2, shape(elem, fe_t.order, i, p, true), phi = , ; break;);
1214 0 : break;
1215 0 : case 3:
1216 0 : fe_vector_scalar_error_switch(3, shape(elem, fe_t.order, i, p, true), phi = , ; break;);
1217 0 : break;
1218 0 : default:
1219 0 : libmesh_error_msg("Invalid dimension = " << dim);
1220 : }
1221 0 : }
1222 :
1223 :
1224 :
1225 : template<>
1226 166085800 : void FEInterface::shape<RealGradient>(const FEType & fe_t,
1227 : int extra_order,
1228 : const Elem * elem,
1229 : const unsigned int i,
1230 : const Point & p,
1231 : RealGradient & phi)
1232 : {
1233 : // This API does not currently support infinite elements.
1234 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1235 33003240 : if (is_InfFE_elem(elem->type()))
1236 : {
1237 0 : libmesh_not_implemented();
1238 : }
1239 : #endif
1240 :
1241 166085800 : auto dim = elem->dim();
1242 :
1243 : // We are calling
1244 : //
1245 : // FE<X,Y>::shape(Elem *, Order, unsigned, Point, false)
1246 : //
1247 : // with the last parameter set to "false" so that the
1248 : // Elem::p_level() is not used internally and the "total_order" that
1249 : // we compute is used instead. See fe.h for more details.
1250 22221960 : auto total_order = fe_t.order + extra_order;
1251 :
1252 166085800 : switch(dim)
1253 : {
1254 0 : case 0:
1255 0 : fe_vector_scalar_error_switch(0, shape(elem, total_order, i, p, false), phi = , ; break;);
1256 0 : break;
1257 0 : case 1:
1258 0 : fe_vector_scalar_error_switch(1, shape(elem, total_order, i, p, false), phi = , ; break;);
1259 0 : break;
1260 46325814 : case 2:
1261 46325814 : fe_vector_scalar_error_switch(2, shape(elem, total_order, i, p, false), phi = , ; break;);
1262 3641840 : break;
1263 119759986 : case 3:
1264 119759986 : fe_vector_scalar_error_switch(3, shape(elem, total_order, i, p, false), phi = , ; break;);
1265 7469140 : break;
1266 0 : default:
1267 0 : libmesh_error_msg("Invalid dimension = " << dim);
1268 : }
1269 166085800 : }
1270 :
1271 :
1272 :
1273 : template<>
1274 0 : void FEInterface::shapes<RealGradient>(const unsigned int dim,
1275 : const FEType & fe_t,
1276 : const Elem * elem,
1277 : const unsigned int i,
1278 : const std::vector<Point> & p,
1279 : std::vector<RealGradient> & phi,
1280 : const bool add_p_level)
1281 : {
1282 : // This is actually an issue for infinite elements: They require type 'Gradient'!
1283 0 : if (elem->infinite())
1284 0 : libmesh_not_implemented();
1285 :
1286 0 : const Order o = fe_t.order;
1287 :
1288 0 : switch(dim)
1289 : {
1290 0 : case 0:
1291 0 : fe_vector_scalar_error_switch(0, shapes(elem,o,i,p,phi,add_p_level), , ; return;);
1292 : break;
1293 0 : case 1:
1294 0 : fe_vector_scalar_error_switch(1, shapes(elem,o,i,p,phi,add_p_level), , ; return;);
1295 : break;
1296 0 : case 2:
1297 0 : fe_vector_scalar_error_switch(2, shapes(elem,o,i,p,phi,add_p_level), , ; return;);
1298 : break;
1299 0 : case 3:
1300 0 : fe_vector_scalar_error_switch(3, shapes(elem,o,i,p,phi,add_p_level), , ; return;);
1301 : break;
1302 0 : default:
1303 0 : libmesh_error_msg("Invalid dimension = " << dim);
1304 : }
1305 :
1306 : return;
1307 : }
1308 :
1309 :
1310 :
1311 : template<>
1312 6445360 : void FEInterface::all_shapes<RealGradient>(const unsigned int dim,
1313 : const FEType & fe_t,
1314 : const Elem * elem,
1315 : const std::vector<Point> & p,
1316 : std::vector<std::vector<RealGradient>> & phi,
1317 : const bool add_p_level)
1318 : {
1319 : // This is actually an issue for infinite elements: They require type 'Gradient'!
1320 1597436 : if (elem->infinite())
1321 0 : libmesh_not_implemented();
1322 :
1323 1075052 : const Order o = fe_t.order;
1324 :
1325 6445360 : switch(dim)
1326 : {
1327 0 : case 0:
1328 0 : fe_vector_scalar_error_switch(0, all_shapes(elem,o,p,phi,add_p_level), , ; return;);
1329 : break;
1330 0 : case 1:
1331 0 : fe_vector_scalar_error_switch(1, all_shapes(elem,o,p,phi,add_p_level), , ; return;);
1332 : break;
1333 1994299 : case 2:
1334 1994299 : fe_vector_scalar_error_switch(2, all_shapes(elem,o,p,phi,add_p_level), , ; return;);
1335 : break;
1336 4451061 : case 3:
1337 4451061 : fe_vector_scalar_error_switch(3, all_shapes(elem,o,p,phi,add_p_level), , ; return;);
1338 : break;
1339 0 : default:
1340 0 : libmesh_error_msg("Invalid dimension = " << dim);
1341 : }
1342 :
1343 : return;
1344 : }
1345 :
1346 :
1347 : #ifdef LIBMESH_ENABLE_DEPRECATED
1348 : FEInterface::shape_ptr
1349 0 : FEInterface::shape_function(const unsigned int dim,
1350 : const FEType & fe_t,
1351 : const ElemType libmesh_inf_var(t))
1352 : {
1353 : libmesh_deprecated();
1354 :
1355 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1356 0 : if (is_InfFE_elem(t))
1357 : {
1358 0 : inf_fe_switch(shape);
1359 : }
1360 : #endif
1361 0 : fe_switch(shape);
1362 : }
1363 : #endif // LIBMESH_ENABLE_DEPRECATED
1364 :
1365 :
1366 :
1367 : FEInterface::shape_ptr
1368 358084951 : FEInterface::shape_function(const FEType & fe_t,
1369 : const Elem * elem)
1370 : {
1371 : // dim is needed by the fe_switch macros below
1372 358084951 : auto dim = elem->dim();
1373 :
1374 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1375 66893619 : if (is_InfFE_elem(elem->type()))
1376 0 : inf_fe_switch(shape);
1377 : #endif
1378 358084951 : fe_switch(shape);
1379 : }
1380 :
1381 :
1382 :
1383 : #ifdef LIBMESH_ENABLE_DEPRECATED
1384 0 : Real FEInterface::shape_deriv(const unsigned int dim,
1385 : const FEType & fe_t,
1386 : const ElemType t,
1387 : const unsigned int i,
1388 : const unsigned int j,
1389 : const Point & p)
1390 : {
1391 : libmesh_deprecated();
1392 :
1393 0 : libmesh_assert_greater (dim,j);
1394 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1395 :
1396 0 : if (is_InfFE_elem(t)){
1397 0 : return ifem_shape_deriv(dim, fe_t, t, i, j, p);
1398 : }
1399 :
1400 : #endif
1401 :
1402 0 : const Order o = fe_t.order;
1403 0 : fe_switch(shape_deriv(t,o,i,j,p));
1404 : return 0;
1405 : }
1406 :
1407 :
1408 0 : Real FEInterface::shape_deriv(const unsigned int dim,
1409 : const FEType & fe_t,
1410 : const Elem * elem,
1411 : const unsigned int i,
1412 : const unsigned int j,
1413 : const Point & p)
1414 : {
1415 : libmesh_deprecated();
1416 :
1417 0 : libmesh_assert_greater (dim,j);
1418 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1419 :
1420 0 : if (elem->infinite()){
1421 0 : return ifem_shape_deriv(fe_t, elem, i, j, p);
1422 : }
1423 :
1424 : #endif
1425 :
1426 0 : const Order o = fe_t.order;
1427 :
1428 0 : switch(dim)
1429 : {
1430 0 : case 0:
1431 0 : fe_family_switch (0, shape_deriv(elem, o, i, j, p), return , ;);
1432 0 : case 1:
1433 0 : fe_family_switch (1, shape_deriv(elem, o, i, j, p), return , ;);
1434 0 : case 2:
1435 0 : fe_family_switch (2, shape_deriv(elem, o, i, j, p), return , ;);
1436 0 : case 3:
1437 0 : fe_family_switch (3, shape_deriv(elem, o, i, j, p), return , ;);
1438 0 : default:
1439 0 : libmesh_error_msg("Invalid dimension = " << dim);
1440 : }
1441 : return 0;
1442 : }
1443 : #endif // LIBMESH_ENABLE_DEPRECATED
1444 :
1445 :
1446 :
1447 0 : Real FEInterface::shape_deriv(const FEType & fe_t,
1448 : const Elem * elem,
1449 : const unsigned int i,
1450 : const unsigned int j,
1451 : const Point & p)
1452 : {
1453 0 : auto dim = elem->dim();
1454 :
1455 0 : libmesh_assert_greater (dim, j);
1456 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1457 :
1458 0 : if (is_InfFE_elem(elem->type())){
1459 0 : return ifem_shape_deriv(fe_t, elem, i, j, p);
1460 : }
1461 :
1462 : #endif
1463 :
1464 : // Account for Elem::p_level() when computing total order. Note: we are calling
1465 : // FE::shape_deriv() with the final argument == true so that the Elem::p_level()
1466 : // is accounted for automatically internally.
1467 0 : fe_switch(shape_deriv(elem, fe_t.order, i, j, p, true));
1468 :
1469 : // We'll never get here
1470 : return 0;
1471 : }
1472 :
1473 :
1474 :
1475 1181419464 : Real FEInterface::shape_deriv(const FEType & fe_t,
1476 : int extra_order,
1477 : const Elem * elem,
1478 : const unsigned int i,
1479 : const unsigned int j,
1480 : const Point & p)
1481 : {
1482 1181419464 : auto dim = elem->dim();
1483 :
1484 95937904 : libmesh_assert_greater (dim, j);
1485 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1486 :
1487 274572966 : if (elem->infinite()){
1488 0 : return ifem_shape_deriv(fe_t, elem, i, j, p);
1489 : }
1490 :
1491 : #endif
1492 :
1493 : // Ignore Elem::p_level() when computing total order, use
1494 : // extra_order instead.
1495 191875808 : auto total_order = fe_t.order + extra_order;
1496 :
1497 : // We call shape_deriv() with the final argument == false so that
1498 : // the Elem::p_level() is ignored internally.
1499 1181419464 : switch(dim)
1500 : {
1501 0 : case 0:
1502 0 : fe_family_switch (0, shape_deriv(elem, total_order, i, j, p, false), return , ;);
1503 1369218 : case 1:
1504 1369218 : fe_family_switch (1, shape_deriv(elem, total_order, i, j, p, false), return , ;);
1505 294841782 : case 2:
1506 294841782 : fe_family_switch (2, shape_deriv(elem, total_order, i, j, p, false), return , ;);
1507 885208464 : case 3:
1508 885208464 : fe_family_switch (3, shape_deriv(elem, total_order, i, j, p, false), return , ;);
1509 0 : default:
1510 0 : libmesh_error_msg("Invalid dimension = " << dim);
1511 : }
1512 :
1513 : // We'll never get here
1514 : return 0;
1515 : }
1516 :
1517 :
1518 :
1519 : template<>
1520 16140090 : void FEInterface::shape_derivs<Real>(const FEType & fe_t,
1521 : const Elem * elem,
1522 : const unsigned int i,
1523 : const unsigned int j,
1524 : const std::vector<Point> & p,
1525 : std::vector<Real> & dphi,
1526 : const bool add_p_level)
1527 : {
1528 2767348 : const Order o = fe_t.order;
1529 :
1530 16140090 : switch(elem->dim())
1531 : {
1532 0 : case 0:
1533 0 : fe_scalar_vec_error_switch(0, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ; break;);
1534 0 : break;
1535 4800 : case 1:
1536 4800 : fe_scalar_vec_error_switch(1, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ; break;);
1537 401 : break;
1538 7182564 : case 2:
1539 7182564 : fe_scalar_vec_error_switch(2, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ; break;);
1540 676944 : break;
1541 8952726 : case 3:
1542 8952726 : fe_scalar_vec_error_switch(3, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ; break;);
1543 706329 : break;
1544 0 : default:
1545 0 : libmesh_error_msg("Invalid dimension = " << elem->dim());
1546 : }
1547 :
1548 16140090 : return;
1549 : }
1550 :
1551 :
1552 :
1553 : // Only instantiating Real because we don't expect vector-valued
1554 : // mapping functions, and that's what's using this interface.
1555 : template<>
1556 38822226 : void FEInterface::all_shape_derivs<Real>(const unsigned int dim,
1557 : const FEType & fe_t,
1558 : const Elem * elem,
1559 : const std::vector<Point> & p,
1560 : std::vector<std::vector<Real>> * comps[3],
1561 : const bool add_p_level)
1562 : {
1563 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1564 :
1565 9877300 : if (elem && is_InfFE_elem(elem->type()))
1566 : {
1567 0 : for (auto j : make_range(dim))
1568 0 : for (auto i : index_range(*comps[j]))
1569 0 : FEInterface::shape_derivs<Real>(fe_t, elem, i, j, p, (*comps[j])[i], add_p_level);
1570 0 : return;
1571 : }
1572 : #endif
1573 :
1574 6837356 : const Order o = fe_t.order;
1575 :
1576 38822226 : switch(dim)
1577 : {
1578 0 : case 0:
1579 0 : fe_scalar_vec_error_switch(0, all_shape_derivs(elem,o,p,comps,add_p_level), , ; return;);
1580 : break;
1581 145248 : case 1:
1582 145248 : fe_scalar_vec_error_switch(1, all_shape_derivs(elem,o,p,comps,add_p_level), , ; return;);
1583 : break;
1584 14592427 : case 2:
1585 14592427 : fe_scalar_vec_error_switch(2, all_shape_derivs(elem,o,p,comps,add_p_level), , ; return;);
1586 : break;
1587 24084551 : case 3:
1588 24084551 : fe_scalar_vec_error_switch(3, all_shape_derivs(elem,o,p,comps,add_p_level), , ; return;);
1589 : break;
1590 0 : default:
1591 0 : libmesh_error_msg("Invalid dimension = " << dim);
1592 : }
1593 :
1594 : return;
1595 : }
1596 :
1597 :
1598 :
1599 : template<>
1600 0 : void FEInterface::shape_derivs<RealGradient>(const FEType & fe_t,
1601 : const Elem * elem,
1602 : const unsigned int i,
1603 : const unsigned int j,
1604 : const std::vector<Point> & p,
1605 : std::vector<RealGradient> & dphi,
1606 : const bool add_p_level)
1607 : {
1608 0 : const Order o = fe_t.order;
1609 :
1610 0 : switch(elem->dim())
1611 : {
1612 0 : case 0:
1613 0 : fe_vector_scalar_error_switch(0, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ; break;);
1614 0 : break;
1615 0 : case 1:
1616 0 : fe_vector_scalar_error_switch(1, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ; break;);
1617 0 : break;
1618 0 : case 2:
1619 0 : fe_vector_scalar_error_switch(2, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ; break;);
1620 0 : break;
1621 0 : case 3:
1622 0 : fe_vector_scalar_error_switch(3, shape_derivs(elem,o,i,j,p,dphi,add_p_level), , ; break;);
1623 0 : break;
1624 0 : default:
1625 0 : libmesh_error_msg("Invalid dimension = " << elem->dim());
1626 : }
1627 :
1628 0 : return;
1629 : }
1630 :
1631 :
1632 :
1633 : #ifdef LIBMESH_ENABLE_DEPRECATED
1634 : FEInterface::shape_deriv_ptr
1635 0 : FEInterface::shape_deriv_function(const unsigned int dim,
1636 : const FEType & fe_t,
1637 : const ElemType libmesh_inf_var(t))
1638 : {
1639 : libmesh_deprecated();
1640 :
1641 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1642 0 : if (is_InfFE_elem(t))
1643 : {
1644 0 : inf_fe_switch(shape_deriv);
1645 : }
1646 : #endif
1647 0 : fe_switch(shape_deriv);
1648 : }
1649 : #endif // LIBMESH_ENABLE_DEPRECATED
1650 :
1651 :
1652 :
1653 : FEInterface::shape_deriv_ptr
1654 1203091436 : FEInterface::shape_deriv_function(const FEType & fe_t,
1655 : const Elem * elem)
1656 : {
1657 : // dim is needed by the fe_switch macros below
1658 1203091436 : auto dim = elem->dim();
1659 :
1660 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1661 216396481 : if (is_InfFE_elem(elem->type()))
1662 : {
1663 0 : inf_fe_switch(shape_deriv);
1664 : }
1665 : #endif
1666 1203091436 : fe_switch(shape_deriv);
1667 : }
1668 :
1669 :
1670 : #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
1671 : #ifdef LIBMESH_ENABLE_DEPRECATED
1672 0 : Real FEInterface::shape_second_deriv(const unsigned int dim,
1673 : const FEType & fe_t,
1674 : const ElemType t,
1675 : const unsigned int i,
1676 : const unsigned int j,
1677 : const Point & p)
1678 : {
1679 : libmesh_deprecated();
1680 :
1681 0 : libmesh_assert_greater_equal (dim*(dim-1),j);
1682 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1683 0 : if (is_InfFE_elem(t))
1684 0 : libmesh_not_implemented();
1685 : #endif
1686 :
1687 0 : const Order o = fe_t.order;
1688 :
1689 0 : switch(dim)
1690 : {
1691 0 : case 0:
1692 0 : fe_family_switch (0, shape_second_deriv(t, o, i, j, p), return , ;);
1693 0 : case 1:
1694 0 : fe_family_switch (1, shape_second_deriv(t, o, i, j, p), return , ;);
1695 0 : case 2:
1696 0 : fe_family_switch (2, shape_second_deriv(t, o, i, j, p), return , ;);
1697 0 : case 3:
1698 0 : fe_family_switch (3, shape_second_deriv(t, o, i, j, p), return , ;);
1699 0 : default:
1700 0 : libmesh_error_msg("Invalid dimension = " << dim);
1701 : }
1702 : return 0;
1703 : }
1704 :
1705 :
1706 0 : Real FEInterface::shape_second_deriv(const unsigned int dim,
1707 : const FEType & fe_t,
1708 : const Elem * elem,
1709 : const unsigned int i,
1710 : const unsigned int j,
1711 : const Point & p)
1712 : {
1713 : libmesh_deprecated();
1714 :
1715 0 : libmesh_assert_greater_equal (dim*(dim-1),j);
1716 :
1717 0 : if (elem->infinite())
1718 0 : libmesh_not_implemented();
1719 :
1720 0 : const Order o = fe_t.order;
1721 :
1722 0 : switch(dim)
1723 : {
1724 0 : case 0:
1725 0 : fe_family_switch (0, shape_second_deriv(elem, o, i, j, p), return , ;);
1726 0 : case 1:
1727 0 : fe_family_switch (1, shape_second_deriv(elem, o, i, j, p), return , ;);
1728 0 : case 2:
1729 0 : fe_family_switch (2, shape_second_deriv(elem, o, i, j, p), return , ;);
1730 0 : case 3:
1731 0 : fe_family_switch (3, shape_second_deriv(elem, o, i, j, p), return , ;);
1732 0 : default:
1733 0 : libmesh_error_msg("Invalid dimension = " << dim);
1734 : }
1735 : return 0;
1736 : }
1737 : #endif // LIBMESH_ENABLE_DEPRECATED
1738 :
1739 :
1740 :
1741 0 : Real FEInterface::shape_second_deriv(const FEType & fe_t,
1742 : const Elem * elem,
1743 : const unsigned int i,
1744 : const unsigned int j,
1745 : const Point & p)
1746 : {
1747 0 : auto dim = elem->dim();
1748 :
1749 0 : libmesh_assert_greater_equal (dim*(dim-1),j);
1750 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1751 0 : if (is_InfFE_elem(elem->type()))
1752 0 : libmesh_not_implemented();
1753 : #endif
1754 :
1755 : // We are calling FE::shape_second_deriv() with the final argument
1756 : // == true so that the Elem::p_level() is accounted for
1757 : // automatically internally.
1758 0 : switch(dim)
1759 : {
1760 0 : case 0:
1761 0 : fe_family_switch (0, shape_second_deriv(elem, fe_t.order, i, j, p, true), return , ;);
1762 0 : case 1:
1763 0 : fe_family_switch (1, shape_second_deriv(elem, fe_t.order, i, j, p, true), return , ;);
1764 0 : case 2:
1765 0 : fe_family_switch (2, shape_second_deriv(elem, fe_t.order, i, j, p, true), return , ;);
1766 0 : case 3:
1767 0 : fe_family_switch (3, shape_second_deriv(elem, fe_t.order, i, j, p, true), return , ;);
1768 0 : default:
1769 0 : libmesh_error_msg("Invalid dimension = " << dim);
1770 : }
1771 :
1772 : // We'll never get here
1773 : return 0;
1774 : }
1775 :
1776 :
1777 :
1778 317113182 : Real FEInterface::shape_second_deriv(const FEType & fe_t,
1779 : int extra_order,
1780 : const Elem * elem,
1781 : const unsigned int i,
1782 : const unsigned int j,
1783 : const Point & p)
1784 : {
1785 317113182 : auto dim = elem->dim();
1786 :
1787 26273056 : libmesh_assert_greater_equal (dim*(dim-1),j);
1788 :
1789 78761391 : if (elem->infinite())
1790 0 : libmesh_not_implemented();
1791 :
1792 : // Ignore Elem::p_level() when computing total order, use
1793 : // extra_order instead.
1794 52546112 : auto total_order = fe_t.order + extra_order;
1795 :
1796 : // We are calling FE::shape_second_deriv() with the final argument
1797 : // == false so that the Elem::p_level() is ignored and the
1798 : // total_order we compute is used instead.
1799 317113182 : switch(dim)
1800 : {
1801 0 : case 0:
1802 0 : fe_family_switch (0, shape_second_deriv(elem, total_order, i, j, p, false), return , ;);
1803 10308 : case 1:
1804 10308 : fe_family_switch (1, shape_second_deriv(elem, total_order, i, j, p, false), return , ;);
1805 2246994 : case 2:
1806 2246994 : fe_family_switch (2, shape_second_deriv(elem, total_order, i, j, p, false), return , ;);
1807 314855880 : case 3:
1808 314855880 : fe_family_switch (3, shape_second_deriv(elem, total_order, i, j, p, false), return , ;);
1809 0 : default:
1810 0 : libmesh_error_msg("Invalid dimension = " << dim);
1811 : }
1812 :
1813 : // We'll never get here
1814 : return 0;
1815 : }
1816 :
1817 :
1818 :
1819 : #ifdef LIBMESH_ENABLE_DEPRECATED
1820 : FEInterface::shape_second_deriv_ptr
1821 0 : FEInterface::shape_second_deriv_function(const unsigned int dim,
1822 : const FEType & fe_t,
1823 : const ElemType libmesh_inf_var(t))
1824 : {
1825 : libmesh_deprecated();
1826 :
1827 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1828 0 : if (is_InfFE_elem(t))
1829 0 : libmesh_not_implemented();
1830 : #endif
1831 0 : fe_switch(shape_second_deriv);
1832 : }
1833 : #endif // LIBMESH_ENABLE_DEPRECATED
1834 :
1835 :
1836 :
1837 : FEInterface::shape_second_deriv_ptr
1838 70982232 : FEInterface::shape_second_deriv_function(const FEType & fe_t,
1839 : const Elem * elem)
1840 : {
1841 70982232 : auto dim = elem->dim();
1842 :
1843 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1844 17144266 : if (is_InfFE_elem(elem->type()))
1845 0 : libmesh_not_implemented();
1846 : #endif
1847 70982232 : fe_switch(shape_second_deriv);
1848 : }
1849 :
1850 : #endif //LIBMESH_ENABLE_SECOND_DERIVATIVES
1851 :
1852 :
1853 : #ifdef LIBMESH_ENABLE_DEPRECATED
1854 : template<>
1855 0 : void FEInterface::shape<RealGradient>(const unsigned int dim,
1856 : const FEType & fe_t,
1857 : const Elem * elem,
1858 : const unsigned int i,
1859 : const Point & p,
1860 : RealGradient & phi)
1861 : {
1862 : libmesh_deprecated();
1863 :
1864 : // This is actually an issue for infinite elements: They require type 'Gradient'!
1865 0 : if (elem->infinite())
1866 0 : libmesh_not_implemented();
1867 :
1868 0 : const Order o = fe_t.order;
1869 :
1870 0 : switch(dim)
1871 : {
1872 0 : case 0:
1873 0 : fe_vector_scalar_error_switch(0, shape(elem,o,i,p), phi = , ; break;);
1874 0 : break;
1875 0 : case 1:
1876 0 : fe_vector_scalar_error_switch(1, shape(elem,o,i,p), phi = , ; break;);
1877 0 : break;
1878 0 : case 2:
1879 0 : fe_vector_scalar_error_switch(2, shape(elem,o,i,p), phi = , ; break;);
1880 0 : break;
1881 0 : case 3:
1882 0 : fe_vector_scalar_error_switch(3, shape(elem,o,i,p), phi = , ; break;);
1883 0 : break;
1884 0 : default:
1885 0 : libmesh_error_msg("Invalid dimension = " << dim);
1886 : }
1887 :
1888 0 : return;
1889 : }
1890 : #endif // LIBMESH_ENABLE_DEPRECATED
1891 :
1892 :
1893 7793931 : void FEInterface::compute_data(const unsigned int dim,
1894 : const FEType & fe_t,
1895 : const Elem * elem,
1896 : FEComputeData & data)
1897 : {
1898 : #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
1899 :
1900 721454 : if (elem && is_InfFE_elem(elem->type()))
1901 : {
1902 0 : data.init();
1903 0 : ifem_compute_data(dim, fe_t, elem, data);
1904 0 : return;
1905 : }
1906 :
1907 : #endif
1908 :
1909 7793931 : const unsigned int n_dof = n_dofs (fe_t, elem);
1910 7793931 : const Point & p = data.p;
1911 7793931 : data.shape.resize(n_dof);
1912 :
1913 7793931 : if (data.need_derivative())
1914 : {
1915 0 : data.dshape.resize(n_dof);
1916 0 : data.local_transform.resize(dim);
1917 :
1918 0 : for (unsigned int d=0; d<dim; d++)
1919 0 : data.local_transform[d].resize(dim);
1920 :
1921 0 : auto fe = FEBase::build(dim, fe_t);
1922 0 : std::vector<Point> pt = {p};
1923 0 : fe->get_dphideta(); // to compute the map
1924 0 : fe->reinit(elem, &pt);
1925 :
1926 : // compute the reference->physical map.
1927 0 : data.local_transform[0][0] = fe->get_dxidx()[0];
1928 0 : if (dim > 1)
1929 : {
1930 0 : data.local_transform[1][0] = fe->get_detadx()[0];
1931 0 : data.local_transform[1][1] = fe->get_detady()[0];
1932 0 : data.local_transform[0][1] = fe->get_dxidy()[0];
1933 0 : if (dim > 2)
1934 : {
1935 0 : data.local_transform[2][0] = fe->get_dzetadx()[0];
1936 0 : data.local_transform[2][1] = fe->get_dzetady()[0];
1937 0 : data.local_transform[2][2] = fe->get_dzetadz()[0];
1938 0 : data.local_transform[1][2] = fe->get_detadz()[0];
1939 0 : data.local_transform[0][2] = fe->get_dxidz()[0];
1940 : }
1941 : }
1942 0 : }
1943 :
1944 : // set default values for all the output fields
1945 7793931 : data.init();
1946 :
1947 49189544 : for (unsigned int n=0; n<n_dof; n++)
1948 : {
1949 : // Here we pass the original fe_t object. Additional p-levels
1950 : // (if any) are handled internally by the shape() and
1951 : // shape_deriv() functions since they have access to the elem
1952 : // pointer. Note that we are already using the n_dof value
1953 : // appropriate to the elevated p-level.
1954 41395613 : data.shape[n] = shape(fe_t, elem, n, p);
1955 41395613 : if (data.need_derivative())
1956 : {
1957 0 : for (unsigned int j=0; j<dim; j++)
1958 0 : data.dshape[n](j) = shape_deriv(fe_t, elem, n, j, p);
1959 : }
1960 : }
1961 7072477 : }
1962 :
1963 :
1964 :
1965 : #ifdef LIBMESH_ENABLE_AMR
1966 :
1967 7863979 : void FEInterface::compute_constraints (DofConstraints & constraints,
1968 : DofMap & dof_map,
1969 : const unsigned int variable_number,
1970 : const Elem * elem)
1971 : {
1972 698537 : libmesh_assert(elem);
1973 :
1974 698537 : const FEType & fe_t = dof_map.variable_type(variable_number);
1975 :
1976 7863979 : switch (elem->dim())
1977 : {
1978 13356 : case 0:
1979 : case 1:
1980 : {
1981 : // No constraints in 0D/1D.
1982 13356 : return;
1983 : }
1984 :
1985 :
1986 4295201 : case 2:
1987 : {
1988 4295201 : switch (fe_t.family)
1989 : {
1990 38660 : case CLOUGH:
1991 38660 : FE<2,CLOUGH>::compute_constraints (constraints,
1992 : dof_map,
1993 : variable_number,
1994 38660 : elem); return;
1995 :
1996 18240 : case HERMITE:
1997 18240 : FE<2,HERMITE>::compute_constraints (constraints,
1998 : dof_map,
1999 : variable_number,
2000 18240 : elem); return;
2001 :
2002 3209263 : case LAGRANGE:
2003 3209263 : FE<2,LAGRANGE>::compute_constraints (constraints,
2004 : dof_map,
2005 : variable_number,
2006 3209263 : elem); return;
2007 :
2008 292412 : case HIERARCHIC:
2009 292412 : FE<2,HIERARCHIC>::compute_constraints (constraints,
2010 : dof_map,
2011 : variable_number,
2012 292412 : elem); return;
2013 :
2014 7860 : case HIERARCHIC_VEC:
2015 7860 : FE<2,HIERARCHIC_VEC>::compute_constraints (constraints,
2016 : dof_map,
2017 : variable_number,
2018 7860 : elem); return;
2019 :
2020 65168 : case SIDE_HIERARCHIC:
2021 65168 : FE<2,SIDE_HIERARCHIC>::compute_constraints (constraints,
2022 : dof_map,
2023 : variable_number,
2024 65168 : elem); return;
2025 :
2026 7256 : case LAGRANGE_VEC:
2027 7256 : FE<2,LAGRANGE_VEC>::compute_constraints (constraints,
2028 : dof_map,
2029 : variable_number,
2030 7256 : elem); return;
2031 :
2032 58726 : fe_family_horder_case_func(2, compute_constraints (constraints,
2033 : dof_map,
2034 : variable_number,
2035 : elem), , ; return;)
2036 49713 : default:
2037 49713 : return;
2038 : }
2039 : }
2040 :
2041 :
2042 3412974 : case 3:
2043 : {
2044 3412974 : switch (fe_t.family)
2045 : {
2046 156 : case HERMITE:
2047 156 : FE<3,HERMITE>::compute_constraints (constraints,
2048 : dof_map,
2049 : variable_number,
2050 156 : elem); return;
2051 :
2052 1310021 : case LAGRANGE:
2053 1310021 : FE<3,LAGRANGE>::compute_constraints (constraints,
2054 : dof_map,
2055 : variable_number,
2056 1310021 : elem); return;
2057 :
2058 30012 : case HIERARCHIC:
2059 30012 : FE<3,HIERARCHIC>::compute_constraints (constraints,
2060 : dof_map,
2061 : variable_number,
2062 30012 : elem); return;
2063 :
2064 145248 : case SIDE_HIERARCHIC:
2065 145248 : FE<3,SIDE_HIERARCHIC>::compute_constraints (constraints,
2066 : dof_map,
2067 : variable_number,
2068 145248 : elem); return;
2069 :
2070 11400 : case LAGRANGE_VEC:
2071 11400 : FE<3,LAGRANGE_VEC>::compute_constraints (constraints,
2072 : dof_map,
2073 : variable_number,
2074 11400 : elem); return;
2075 :
2076 0 : case HIERARCHIC_VEC:
2077 0 : FE<3,HIERARCHIC_VEC>::compute_constraints (constraints,
2078 : dof_map,
2079 : variable_number,
2080 0 : elem); return;
2081 :
2082 16120 : fe_family_horder_case_func(3, compute_constraints (constraints,
2083 : dof_map,
2084 : variable_number,
2085 : elem), , ; return;)
2086 154914 : default:
2087 154914 : return;
2088 : }
2089 : }
2090 :
2091 :
2092 0 : default:
2093 0 : libmesh_error_msg("Invalid dimension = " << elem->dim());
2094 : }
2095 : }
2096 :
2097 : #endif // #ifdef LIBMESH_ENABLE_AMR
2098 :
2099 :
2100 :
2101 : #ifdef LIBMESH_ENABLE_PERIODIC
2102 :
2103 264164 : void FEInterface::compute_periodic_constraints (DofConstraints & constraints,
2104 : DofMap & dof_map,
2105 : const PeriodicBoundaries & boundaries,
2106 : const MeshBase & mesh,
2107 : const PointLocatorBase * point_locator,
2108 : const unsigned int variable_number,
2109 : const Elem * elem)
2110 : {
2111 67082 : const FEType & fe_t = dof_map.variable_type(variable_number);
2112 : // No element-specific optimizations currently exist, although
2113 : // we do have to select the right compute_periodic_constraints
2114 : // for the OutputType of this FEType.
2115 264164 : switch (field_type(fe_t))
2116 : {
2117 264164 : case TYPE_SCALAR:
2118 264164 : FEBase::compute_periodic_constraints(
2119 : constraints, dof_map, boundaries, mesh, point_locator, variable_number, elem);
2120 264164 : break;
2121 0 : case TYPE_VECTOR:
2122 0 : FEVectorBase::compute_periodic_constraints(
2123 : constraints, dof_map, boundaries, mesh, point_locator, variable_number, elem);
2124 0 : break;
2125 0 : default:
2126 0 : libmesh_error_msg(
2127 : "compute_periodic_constraints only set up for vector or scalar FEFieldTypes");
2128 : }
2129 264164 : }
2130 :
2131 : #endif // #ifdef LIBMESH_ENABLE_PERIODIC
2132 :
2133 :
2134 :
2135 31314218 : unsigned int FEInterface::max_order(const FEType & fe_t,
2136 : const ElemType & el_t)
2137 : {
2138 : // Yeah, I know, infinity is much larger than 11, but our
2139 : // solvers don't seem to like high degree polynomials, and our
2140 : // quadrature rules and number_lookups tables
2141 : // need to go up higher.
2142 1917518 : const unsigned int unlimited = 11;
2143 :
2144 : // If we used 0 as a default, then elements missing from this
2145 : // table (e.g. infinite elements) would be considered broken.
2146 1917518 : const unsigned int unknown = unlimited;
2147 :
2148 31314218 : switch (fe_t.family)
2149 : {
2150 18733672 : case LAGRANGE:
2151 : case LAGRANGE_VEC:
2152 18733672 : switch (el_t)
2153 : {
2154 121838 : case EDGE2:
2155 : case EDGE3:
2156 : case EDGE4:
2157 121838 : return 3;
2158 229940 : case TRI3:
2159 : case TRISHELL3:
2160 : case C0POLYGON:
2161 229940 : return 1;
2162 7896 : case TRI6:
2163 7896 : return 2;
2164 856 : case TRI7:
2165 856 : return 3;
2166 458168 : case QUAD4:
2167 : case QUADSHELL4:
2168 458168 : return 1;
2169 159952 : case QUAD8:
2170 : case QUADSHELL8:
2171 : case QUAD9:
2172 : case QUADSHELL9:
2173 159952 : return 2;
2174 78044 : case TET4:
2175 78044 : return 1;
2176 2976 : case TET10:
2177 2976 : return 2;
2178 3648 : case TET14:
2179 3648 : return 3;
2180 129848 : case HEX8:
2181 129848 : return 1;
2182 8452 : case HEX20:
2183 : case HEX27:
2184 8452 : return 2;
2185 4328 : case PRISM6:
2186 4328 : return 1;
2187 4672 : case PRISM15:
2188 : case PRISM18:
2189 4672 : return 2;
2190 2376 : case PRISM20:
2191 : case PRISM21:
2192 2376 : return 3;
2193 1848 : case PYRAMID5:
2194 1848 : return 1;
2195 1008 : case PYRAMID13:
2196 : case PYRAMID14:
2197 1008 : return 2;
2198 672 : case PYRAMID18:
2199 672 : return 3;
2200 52 : case C0POLYHEDRON:
2201 52 : return 1;
2202 4716 : default:
2203 4716 : return unknown;
2204 : }
2205 : break;
2206 14040 : case L2_LAGRANGE_GLL:
2207 14040 : switch (el_t)
2208 : {
2209 1872 : case EDGE2:
2210 : case EDGE3:
2211 : case EDGE4:
2212 : case QUAD4:
2213 : case QUADSHELL4:
2214 : case QUAD8:
2215 : case QUADSHELL8:
2216 : case QUAD9:
2217 : case QUADSHELL9:
2218 : case HEX8:
2219 : case HEX20:
2220 : case HEX27:
2221 : // A basis of degree p needs p+1 Gauss-Lobatto points, and QGaussLobatto tabulates
2222 : // at most fe_lagrange_gll_max_points_1D
2223 1872 : return fe_lagrange_gll_max_points_1D - 1;
2224 0 : default:
2225 0 : return unknown;
2226 : }
2227 : break;
2228 660808 : case MONOMIAL:
2229 : case L2_LAGRANGE:
2230 : case L2_LAGRANGE_VEC:
2231 : case L2_HIERARCHIC:
2232 : case L2_HIERARCHIC_VEC:
2233 : case MONOMIAL_VEC:
2234 330404 : switch (el_t)
2235 : {
2236 330396 : case EDGE2:
2237 : case EDGE3:
2238 : case EDGE4:
2239 : case C0POLYGON:
2240 : case TRI3:
2241 : case TRISHELL3:
2242 : case TRI6:
2243 : case TRI7:
2244 : case QUAD4:
2245 : case QUADSHELL4:
2246 : case QUAD8:
2247 : case QUADSHELL8:
2248 : case QUAD9:
2249 : case QUADSHELL9:
2250 : case TET4:
2251 : case TET10:
2252 : case TET14:
2253 : case HEX8:
2254 : case HEX20:
2255 : case HEX27:
2256 : case PRISM6:
2257 : case PRISM15:
2258 : case PRISM18:
2259 : case PRISM20:
2260 : case PRISM21:
2261 : case PYRAMID5:
2262 : case PYRAMID13:
2263 : case PYRAMID14:
2264 : case PYRAMID18:
2265 : case C0POLYHEDRON:
2266 330396 : return unlimited;
2267 8 : default:
2268 8 : return unknown;
2269 : }
2270 : break;
2271 : #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
2272 383476 : case BERNSTEIN:
2273 : case RATIONAL_BERNSTEIN:
2274 383476 : switch (el_t)
2275 : {
2276 312 : case EDGE2:
2277 : case EDGE3:
2278 : case EDGE4:
2279 312 : return unlimited;
2280 208 : case TRI3:
2281 : case TRISHELL3:
2282 208 : return 1;
2283 624 : case TRI6:
2284 : case TRI7:
2285 624 : return 6;
2286 104 : case QUAD4:
2287 : case QUADSHELL4:
2288 104 : return 1;
2289 7392 : case QUAD8:
2290 : case QUADSHELL8:
2291 : case QUAD9:
2292 : case QUADSHELL9:
2293 7392 : return unlimited;
2294 2496 : case TET4:
2295 2496 : return 1;
2296 2496 : case TET10:
2297 : case TET14:
2298 2496 : return 2;
2299 104 : case HEX8:
2300 104 : return 1;
2301 0 : case HEX20:
2302 0 : return 2;
2303 2920 : case HEX27:
2304 2920 : return 4;
2305 0 : case C0POLYGON:
2306 : case PRISM6:
2307 : case PRISM15:
2308 : case PRISM18:
2309 : case PRISM20:
2310 : case PRISM21:
2311 : case PYRAMID5:
2312 : case PYRAMID13:
2313 : case PYRAMID14:
2314 : case PYRAMID18:
2315 : case C0POLYHEDRON:
2316 0 : return 0;
2317 9668 : default:
2318 9668 : return unknown;
2319 : }
2320 : break;
2321 12588 : case SZABAB:
2322 12588 : switch (el_t)
2323 : {
2324 52 : case EDGE2:
2325 52 : return 1;
2326 156 : case EDGE3:
2327 : case EDGE4:
2328 156 : return 7;
2329 104 : case TRI3:
2330 : case TRISHELL3:
2331 104 : return 1;
2332 688 : case TRI6:
2333 : case TRI7:
2334 688 : return 7;
2335 52 : case QUAD4:
2336 : case QUADSHELL4:
2337 52 : return 1;
2338 156 : case QUAD8:
2339 : case QUADSHELL8:
2340 : case QUAD9:
2341 : case QUADSHELL9:
2342 156 : return 7;
2343 0 : case C0POLYGON:
2344 : case TET4:
2345 : case TET10:
2346 : case TET14:
2347 : case HEX8:
2348 : case HEX20:
2349 : case HEX27:
2350 : case PRISM6:
2351 : case PRISM15:
2352 : case PRISM18:
2353 : case PRISM20:
2354 : case PRISM21:
2355 : case PYRAMID5:
2356 : case PYRAMID13:
2357 : case PYRAMID14:
2358 : case PYRAMID18:
2359 : case C0POLYHEDRON:
2360 0 : return 0;
2361 0 : default:
2362 0 : return unknown;
2363 : }
2364 : break;
2365 : #endif
2366 50720 : case XYZ:
2367 25360 : switch (el_t)
2368 : {
2369 25360 : case EDGE2:
2370 : case EDGE3:
2371 : case EDGE4:
2372 : case C0POLYGON:
2373 : case TRI3:
2374 : case TRISHELL3:
2375 : case TRI6:
2376 : case TRI7:
2377 : case QUAD4:
2378 : case QUADSHELL4:
2379 : case QUAD8:
2380 : case QUADSHELL8:
2381 : case QUAD9:
2382 : case QUADSHELL9:
2383 : case TET4:
2384 : case TET10:
2385 : case TET14:
2386 : case HEX8:
2387 : case HEX20:
2388 : case HEX27:
2389 : case PRISM6:
2390 : case PRISM15:
2391 : case PRISM18:
2392 : case PRISM20:
2393 : case PRISM21:
2394 : case PYRAMID5:
2395 : case PYRAMID13:
2396 : case PYRAMID14:
2397 : case PYRAMID18:
2398 : case C0POLYHEDRON:
2399 25360 : return unlimited;
2400 0 : default:
2401 0 : return unknown;
2402 : }
2403 : break;
2404 112492 : case CLOUGH:
2405 112492 : switch (el_t)
2406 : {
2407 0 : case EDGE2:
2408 : case EDGE3:
2409 0 : return 3;
2410 0 : case EDGE4:
2411 : case TRI3:
2412 : case TRISHELL3:
2413 0 : return 0;
2414 9900 : case TRI6:
2415 : case TRI7:
2416 9900 : return 3;
2417 0 : case C0POLYGON:
2418 : case QUAD4:
2419 : case QUADSHELL4:
2420 : case QUAD8:
2421 : case QUADSHELL8:
2422 : case QUAD9:
2423 : case QUADSHELL9:
2424 : case TET4:
2425 : case TET10:
2426 : case TET14:
2427 : case HEX8:
2428 : case HEX20:
2429 : case HEX27:
2430 : case PRISM6:
2431 : case PRISM15:
2432 : case PRISM18:
2433 : case PRISM20:
2434 : case PRISM21:
2435 : case PYRAMID5:
2436 : case PYRAMID13:
2437 : case PYRAMID14:
2438 : case PYRAMID18:
2439 : case C0POLYHEDRON:
2440 0 : return 0;
2441 0 : default:
2442 0 : return unknown;
2443 : }
2444 : break;
2445 667232 : case HERMITE:
2446 667232 : switch (el_t)
2447 : {
2448 21152 : case EDGE2:
2449 : case EDGE3:
2450 21152 : return unlimited;
2451 0 : case EDGE4:
2452 : case TRI3:
2453 : case TRISHELL3:
2454 : case TRI6:
2455 : case TRI7:
2456 0 : return 0;
2457 0 : case QUAD4:
2458 : case QUADSHELL4:
2459 0 : return 3;
2460 4596 : case QUAD8:
2461 : case QUADSHELL8:
2462 : case QUAD9:
2463 : case QUADSHELL9:
2464 4596 : return unlimited;
2465 0 : case TET4:
2466 : case TET10:
2467 : case TET14:
2468 0 : return 0;
2469 0 : case HEX8:
2470 0 : return 3;
2471 52 : case HEX20:
2472 : case HEX27:
2473 52 : return unlimited;
2474 0 : case C0POLYGON:
2475 : case PRISM6:
2476 : case PRISM15:
2477 : case PRISM18:
2478 : case PRISM20:
2479 : case PRISM21:
2480 : case PYRAMID5:
2481 : case PYRAMID13:
2482 : case PYRAMID14:
2483 : case PYRAMID18:
2484 : case C0POLYHEDRON:
2485 0 : return 0;
2486 0 : default:
2487 0 : return unknown;
2488 : }
2489 : break;
2490 2843630 : case HIERARCHIC:
2491 : case HIERARCHIC_VEC:
2492 2843630 : switch (el_t)
2493 : {
2494 856 : case EDGE2:
2495 : case EDGE3:
2496 : case EDGE4:
2497 856 : return unlimited;
2498 104 : case TRI3:
2499 : case TRISHELL3:
2500 104 : return 1;
2501 27576 : case TRI6:
2502 : case TRI7:
2503 27576 : return unlimited;
2504 52 : case QUAD4:
2505 : case QUADSHELL4:
2506 52 : return 1;
2507 42068 : case QUAD8:
2508 : case QUADSHELL8:
2509 : case QUAD9:
2510 : case QUADSHELL9:
2511 42068 : return unlimited;
2512 1248 : case TET4:
2513 1248 : return 1;
2514 1248 : case TET10:
2515 1248 : return 2;
2516 4040 : case TET14:
2517 4040 : return unlimited;
2518 52 : case HEX8:
2519 : case HEX20:
2520 52 : return 1;
2521 2308 : case HEX27:
2522 2308 : return unlimited;
2523 1080 : case PRISM6:
2524 : case PRISM15:
2525 : case PRISM18:
2526 : case PRISM20:
2527 : case PRISM21:
2528 1080 : return unlimited;
2529 0 : case C0POLYGON:
2530 : case PYRAMID5:
2531 : case PYRAMID13:
2532 : case PYRAMID14:
2533 : case PYRAMID18:
2534 : case C0POLYHEDRON:
2535 0 : return 0;
2536 0 : default:
2537 0 : return unknown;
2538 : }
2539 : break;
2540 1157596 : case SIDE_HIERARCHIC:
2541 1157596 : switch (el_t)
2542 : {
2543 264 : case EDGE2:
2544 : case EDGE3:
2545 : case EDGE4:
2546 264 : return unlimited; // although it's all the same as 0...
2547 0 : case TRI3:
2548 : case TRISHELL3:
2549 0 : return 0;
2550 18568 : case TRI6:
2551 : case TRI7:
2552 18568 : return unlimited;
2553 0 : case C0POLYGON:
2554 : case QUAD4:
2555 : case QUADSHELL4:
2556 0 : return 0;
2557 2728 : case QUAD8:
2558 : case QUADSHELL8:
2559 : case QUAD9:
2560 : case QUADSHELL9:
2561 2728 : return unlimited;
2562 0 : case TET4:
2563 : case TET10:
2564 0 : return 0;
2565 43776 : case TET14:
2566 43776 : return unlimited;
2567 0 : case HEX8:
2568 : case HEX20:
2569 0 : return 0;
2570 1440 : case HEX27:
2571 1440 : return unlimited;
2572 0 : case PRISM6:
2573 : case PRISM15:
2574 : case PRISM18:
2575 0 : return 0;
2576 3200 : case PRISM20:
2577 : case PRISM21:
2578 3200 : return unlimited;
2579 0 : case PYRAMID5:
2580 : case PYRAMID13:
2581 : case PYRAMID14:
2582 : case PYRAMID18:
2583 : case C0POLYHEDRON:
2584 0 : return 0;
2585 0 : default:
2586 0 : return unknown;
2587 : }
2588 : break;
2589 2592 : case SUBDIVISION:
2590 1296 : switch (el_t)
2591 : {
2592 1296 : case TRI3SUBDIVISION:
2593 1296 : return unlimited;
2594 0 : default:
2595 0 : return unknown;
2596 : }
2597 : break;
2598 386618 : case NEDELEC_ONE:
2599 386618 : switch (el_t)
2600 : {
2601 11632 : case TRI6:
2602 : case TRI7:
2603 : case QUAD8:
2604 : case QUAD9:
2605 11632 : return 5;
2606 6400 : case TET10:
2607 : case TET14:
2608 : case HEX20:
2609 : case HEX27:
2610 6400 : return 1;
2611 0 : default:
2612 0 : return 0;
2613 : }
2614 : break;
2615 1514630 : case RAVIART_THOMAS:
2616 : case L2_RAVIART_THOMAS:
2617 1514630 : switch (el_t)
2618 : {
2619 33200 : case TRI6:
2620 : case TRI7:
2621 : case QUAD8:
2622 : case QUAD9:
2623 33200 : return 5;
2624 74872 : case TET14:
2625 : case HEX27:
2626 74872 : return 1;
2627 0 : default:
2628 0 : return 0;
2629 : }
2630 : break;
2631 0 : default:
2632 0 : return 0;
2633 : break;
2634 : }
2635 : }
2636 :
2637 :
2638 :
2639 320718417 : bool FEInterface::extra_hanging_dofs(const FEType & fe_t)
2640 : {
2641 320718417 : switch (fe_t.family)
2642 : {
2643 19662019 : case LAGRANGE:
2644 : case L2_LAGRANGE:
2645 : case L2_LAGRANGE_GLL:
2646 : case MONOMIAL:
2647 : case MONOMIAL_VEC:
2648 : case L2_HIERARCHIC:
2649 : case SIDE_HIERARCHIC:
2650 : case XYZ:
2651 : case SUBDIVISION:
2652 : case LAGRANGE_VEC:
2653 : case L2_LAGRANGE_VEC:
2654 : case NEDELEC_ONE:
2655 : case RAVIART_THOMAS:
2656 : case L2_RAVIART_THOMAS:
2657 19662019 : return false;
2658 92347059 : case CLOUGH:
2659 : case HERMITE:
2660 : case HIERARCHIC:
2661 : case HIERARCHIC_VEC:
2662 : case L2_HIERARCHIC_VEC:
2663 : fe_family_horder_case()
2664 : default:
2665 92347059 : return true;
2666 : }
2667 : }
2668 :
2669 3333426 : FEFieldType FEInterface::field_type(const FEType & fe_type)
2670 : {
2671 3333426 : return FEInterface::field_type(fe_type.family);
2672 : }
2673 :
2674 36469023 : FEFieldType FEInterface::field_type (const FEFamily & fe_family)
2675 : {
2676 36469023 : switch (fe_family)
2677 : {
2678 8317 : fe_family_vector_case()
2679 8317 : return TYPE_VECTOR;
2680 36182936 : default:
2681 36182936 : return TYPE_SCALAR;
2682 : }
2683 : }
2684 :
2685 123975797 : bool FEInterface::orientation_dependent (const FEFamily & fe_family)
2686 : {
2687 123975797 : switch (fe_family)
2688 : {
2689 604848 : case HIERARCHIC:
2690 : case L2_HIERARCHIC:
2691 : case HIERARCHIC_VEC:
2692 : case L2_HIERARCHIC_VEC:
2693 : case BERNSTEIN:
2694 : case RATIONAL_BERNSTEIN:
2695 : case SZABAB:
2696 : case NEDELEC_ONE:
2697 : case RAVIART_THOMAS:
2698 604848 : return true;
2699 106718257 : default:
2700 106718257 : return false;
2701 : }
2702 : }
2703 :
2704 32962742 : unsigned int FEInterface::n_vec_dim (const MeshBase & mesh,
2705 : const FEType & fe_type)
2706 : {
2707 : // We assume the number of vector components is the mesh spatial dimension.
2708 32962742 : return field_type(fe_type.family) == TYPE_VECTOR ? mesh.spatial_dimension() : 1;
2709 : }
2710 :
2711 :
2712 :
2713 :
2714 : bool
2715 6850 : FEInterface::is_hierarchic (const FEType & fe_type)
2716 : {
2717 6850 : switch (fe_type.family)
2718 : {
2719 284 : case HERMITE:
2720 : case HIERARCHIC:
2721 : case HIERARCHIC_VEC:
2722 : case L2_HIERARCHIC:
2723 : case L2_HIERARCHIC_VEC:
2724 : case MONOMIAL:
2725 : case MONOMIAL_VEC:
2726 : case SIDE_HIERARCHIC:
2727 : case SZABAB:
2728 : case XYZ:
2729 284 : return true;
2730 :
2731 2726 : case BERNSTEIN:
2732 : case CLOUGH: // maybe some day?
2733 : case LAGRANGE:
2734 : case LAGRANGE_VEC:
2735 : case L2_LAGRANGE:
2736 : case L2_LAGRANGE_GLL:
2737 : case L2_LAGRANGE_VEC:
2738 : case NEDELEC_ONE:
2739 : case RATIONAL_BERNSTEIN:
2740 : case RAVIART_THOMAS:
2741 : case L2_RAVIART_THOMAS:
2742 : case SCALAR:
2743 : case SUBDIVISION:
2744 2726 : return false;
2745 :
2746 0 : default:
2747 0 : libmesh_error_msg("Unknown FE Family " << Utility::enum_to_string(fe_type.family));
2748 : }
2749 : }
2750 :
2751 :
2752 :
2753 15112033 : FEContinuity FEInterface::get_continuity(const FEType & fe_type)
2754 : {
2755 15112033 : switch (fe_type.family)
2756 : {
2757 : // Discontinuous elements
2758 12784 : case MONOMIAL:
2759 : case MONOMIAL_VEC:
2760 : case L2_HIERARCHIC:
2761 : case L2_LAGRANGE:
2762 : case L2_LAGRANGE_GLL:
2763 : case XYZ:
2764 : case SCALAR:
2765 : case L2_RAVIART_THOMAS:
2766 : case L2_HIERARCHIC_VEC:
2767 : case L2_LAGRANGE_VEC:
2768 12784 : return DISCONTINUOUS;
2769 :
2770 : // C0 elements
2771 14580256 : case LAGRANGE:
2772 : case HIERARCHIC:
2773 : case BERNSTEIN:
2774 : case SZABAB:
2775 : case RATIONAL_BERNSTEIN:
2776 : case INFINITE_MAP:
2777 : case JACOBI_20_00:
2778 : case JACOBI_30_00:
2779 : case LEGENDRE:
2780 : case LAGRANGE_VEC:
2781 : case HIERARCHIC_VEC:
2782 14580256 : return C_ZERO;
2783 :
2784 : // C1 elements
2785 31988 : case CLOUGH:
2786 : case HERMITE:
2787 : case SUBDIVISION:
2788 31988 : return C_ONE;
2789 :
2790 3536 : case NEDELEC_ONE:
2791 3536 : return H_CURL;
2792 :
2793 6720 : case RAVIART_THOMAS:
2794 6720 : return H_DIV;
2795 :
2796 : // Side elements
2797 47478 : case SIDE_HIERARCHIC:
2798 47478 : return SIDE_DISCONTINUOUS;
2799 :
2800 0 : default:
2801 0 : libmesh_error_msg("Unknown FE Family " << Utility::enum_to_string(fe_type.family));
2802 : }
2803 : }
2804 :
2805 :
2806 : } // namespace libMesh
|
