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/libmesh_config.h"
22 : #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
23 :
24 : #include "libmesh/elem.h"
25 : #include "libmesh/enum_to_string.h"
26 : #include "libmesh/fe.h"
27 : #include "libmesh/fe_interface.h"
28 : #include "libmesh/fe_macro.h"
29 :
30 : namespace libMesh
31 : {
32 :
33 : // ------------------------------------------------------------
34 : // Szabo-Babuska-specific implementations, Steffen Petersen 2004
35 :
36 : // Anonymous namespace for local helper functions
37 : namespace {
38 :
39 576 : void szabab_nodal_soln(const Elem * elem,
40 : const Order order,
41 : const std::vector<Number> & elem_soln,
42 : std::vector<Number> & nodal_soln,
43 : const bool add_p_level)
44 : {
45 576 : const unsigned int n_nodes = elem->n_nodes();
46 :
47 576 : const ElemType elem_type = elem->type();
48 :
49 576 : nodal_soln.resize(n_nodes);
50 :
51 624 : const Order totalorder = order + add_p_level*elem->p_level();
52 :
53 : // FEType object to be passed to various FEInterface functions below.
54 48 : FEType fe_type(order, SZABAB);
55 :
56 576 : switch (totalorder)
57 : {
58 : // Constant shape functions
59 0 : case CONSTANT:
60 : {
61 0 : libmesh_assert_equal_to (elem_soln.size(), 1);
62 :
63 0 : std::fill(nodal_soln.begin(), nodal_soln.end(), elem_soln[0]);
64 :
65 0 : return;
66 : }
67 :
68 :
69 : // For other bases do interpolation at the nodes
70 : // explicitly.
71 576 : case FIRST:
72 : case SECOND:
73 : case THIRD:
74 : case FOURTH:
75 : case FIFTH:
76 : case SIXTH:
77 : case SEVENTH:
78 : {
79 : const unsigned int n_sf =
80 576 : FEInterface::n_shape_functions(fe_type, elem);
81 :
82 96 : std::vector<Point> refspace_nodes;
83 576 : FEBase::get_refspace_nodes(elem_type,refspace_nodes);
84 48 : libmesh_assert_equal_to (refspace_nodes.size(), n_nodes);
85 48 : libmesh_assert_equal_to (elem_soln.size(), n_sf);
86 :
87 : // Zero before summation
88 48 : std::fill(nodal_soln.begin(), nodal_soln.end(), 0);
89 :
90 4032 : for (unsigned int n=0; n<n_nodes; n++)
91 : // u_i = Sum (alpha_i phi_i)
92 24192 : for (unsigned int i=0; i<n_sf; i++)
93 22464 : nodal_soln[n] += elem_soln[i] *
94 22464 : FEInterface::shape(fe_type, elem, i, refspace_nodes[n]);
95 :
96 48 : return;
97 : }
98 :
99 0 : default:
100 0 : libmesh_error_msg("ERROR: Invalid total order " << totalorder);
101 : }
102 : } // szabab_nodal_soln()
103 :
104 :
105 :
106 58966 : unsigned int szabab_n_dofs(const ElemType t, const Order o)
107 : {
108 58966 : switch (o)
109 : {
110 : // Szabo-Babuska 1st-order polynomials.
111 8526 : case FIRST:
112 : {
113 8526 : switch (t)
114 : {
115 0 : case NODEELEM:
116 0 : return 1;
117 :
118 168 : case EDGE2:
119 : case EDGE3:
120 168 : return 2;
121 :
122 930 : case TRI3:
123 : case TRISHELL3:
124 : case TRI6:
125 : case TRI7:
126 930 : return 3;
127 :
128 764 : case QUAD4:
129 : case QUADSHELL4:
130 : case QUAD8:
131 : case QUADSHELL8:
132 : case QUAD9:
133 : case QUADSHELL9:
134 764 : return 4;
135 :
136 0 : case INVALID_ELEM:
137 0 : return 0;
138 :
139 0 : default:
140 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
141 : }
142 : }
143 :
144 :
145 : // Szabo-Babuska 2nd-order polynomials.
146 17484 : case SECOND:
147 : {
148 17484 : switch (t)
149 : {
150 0 : case NODEELEM:
151 0 : return 1;
152 :
153 220 : case EDGE2:
154 : case EDGE3:
155 220 : return 3;
156 :
157 2928 : case TRI6:
158 : case TRI7:
159 2928 : return 6;
160 :
161 0 : case QUAD8:
162 : case QUADSHELL8:
163 0 : return 8;
164 1024 : case QUAD9:
165 : case QUADSHELL9:
166 1024 : return 9;
167 :
168 0 : case INVALID_ELEM:
169 0 : return 0;
170 :
171 0 : default:
172 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
173 : }
174 : }
175 :
176 :
177 : // Szabo-Babuska 3rd-order polynomials.
178 16478 : case THIRD:
179 : {
180 16478 : switch (t)
181 : {
182 0 : case NODEELEM:
183 0 : return 1;
184 :
185 220 : case EDGE2:
186 : case EDGE3:
187 220 : return 4;
188 :
189 2718 : case TRI6:
190 : case TRI7:
191 2718 : return 10;
192 :
193 4844 : case QUAD8:
194 : case QUADSHELL8:
195 : case QUAD9:
196 : case QUADSHELL9:
197 4844 : return 16;
198 :
199 0 : case INVALID_ELEM:
200 0 : return 0;
201 :
202 0 : default:
203 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
204 : }
205 : }
206 :
207 :
208 : // Szabo-Babuska 4th-order polynomials.
209 16478 : case FOURTH:
210 : {
211 16478 : switch (t)
212 : {
213 0 : case NODEELEM:
214 0 : return 1;
215 :
216 220 : case EDGE2:
217 : case EDGE3:
218 220 : return 5;
219 :
220 10834 : case TRI6:
221 : case TRI7:
222 10834 : return 15;
223 :
224 4844 : case QUAD8:
225 : case QUADSHELL8:
226 : case QUAD9:
227 : case QUADSHELL9:
228 4844 : return 25;
229 :
230 0 : case INVALID_ELEM:
231 0 : return 0;
232 :
233 0 : default:
234 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
235 : }
236 : }
237 :
238 :
239 : // Szabo-Babuska 5th-order polynomials.
240 0 : case FIFTH:
241 : {
242 0 : switch (t)
243 : {
244 0 : case NODEELEM:
245 0 : return 1;
246 :
247 0 : case EDGE2:
248 : case EDGE3:
249 0 : return 6;
250 :
251 0 : case TRI6:
252 : case TRI7:
253 0 : return 21;
254 :
255 0 : case QUAD8:
256 : case QUADSHELL8:
257 : case QUAD9:
258 : case QUADSHELL9:
259 0 : return 36;
260 :
261 0 : case INVALID_ELEM:
262 0 : return 0;
263 :
264 0 : default:
265 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
266 : }
267 : }
268 :
269 :
270 : // Szabo-Babuska 6th-order polynomials.
271 0 : case SIXTH:
272 : {
273 0 : switch (t)
274 : {
275 0 : case NODEELEM:
276 0 : return 1;
277 :
278 0 : case EDGE2:
279 : case EDGE3:
280 0 : return 7;
281 :
282 0 : case TRI6:
283 : case TRI7:
284 0 : return 28;
285 :
286 0 : case QUAD8:
287 : case QUADSHELL8:
288 : case QUAD9:
289 : case QUADSHELL9:
290 0 : return 49;
291 :
292 0 : case INVALID_ELEM:
293 0 : return 0;
294 :
295 0 : default:
296 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
297 : }
298 : }
299 :
300 : // Szabo-Babuska 7th-order polynomials.
301 0 : case SEVENTH:
302 : {
303 0 : switch (t)
304 : {
305 0 : case NODEELEM:
306 0 : return 1;
307 :
308 0 : case EDGE2:
309 : case EDGE3:
310 0 : return 8;
311 :
312 0 : case TRI6:
313 : case TRI7:
314 0 : return 36;
315 :
316 0 : case QUAD8:
317 : case QUADSHELL8:
318 : case QUAD9:
319 : case QUADSHELL9:
320 0 : return 64;
321 :
322 0 : case INVALID_ELEM:
323 0 : return 0;
324 :
325 0 : default:
326 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
327 : }
328 : }
329 :
330 :
331 0 : default:
332 0 : libmesh_error_msg("ERROR: Invalid Order " << Utility::enum_to_string(o) << " selected for SZABAB FE family!");
333 : }
334 : } // szabab_n_dofs()
335 :
336 :
337 :
338 54247 : unsigned int szabab_n_dofs(const Elem * e, const Order o)
339 : {
340 9239 : libmesh_assert(e);
341 58966 : return szabab_n_dofs(e->type(), o);
342 : }
343 :
344 :
345 :
346 452466 : unsigned int szabab_n_dofs_at_node(const ElemType t,
347 : const Order o,
348 : const unsigned int n)
349 : {
350 452466 : switch (o)
351 : {
352 : // The first-order Szabo-Babuska shape functions
353 23114 : case FIRST:
354 : {
355 23114 : switch (t)
356 : {
357 0 : case NODEELEM:
358 0 : return 1;
359 :
360 : // The 1D Szabo-Babuska defined on a three-noded edge
361 2756 : case EDGE2:
362 : case EDGE3:
363 : {
364 2756 : switch (n)
365 : {
366 273 : case 0:
367 : case 1:
368 273 : return 1;
369 :
370 0 : case 2:
371 0 : return 0;
372 :
373 0 : default:
374 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for EDGE2/3!");
375 : }
376 : }
377 :
378 :
379 : // The 2D Szabo-Babuska defined on a 6-noded triangle
380 12506 : case TRI3:
381 : case TRISHELL3:
382 : case TRI6:
383 : case TRI7:
384 : {
385 11362 : switch (n)
386 : {
387 1144 : case 0:
388 : case 1:
389 : case 2:
390 1144 : return 1;
391 :
392 0 : case 3:
393 : case 4:
394 : case 5:
395 0 : return 0;
396 :
397 0 : case 6:
398 0 : libmesh_assert(t == TRI7);
399 0 : return 0;
400 :
401 0 : default:
402 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for TRI3/6/7!");
403 : }
404 : }
405 :
406 :
407 : // The 2D tensor-product Szabo-Babuska defined on a
408 : // nine-noded quadrilateral.
409 7852 : case QUAD4:
410 : case QUADSHELL4:
411 : case QUAD8:
412 : case QUADSHELL8:
413 : case QUAD9:
414 : case QUADSHELL9:
415 : {
416 7111 : switch (n)
417 : {
418 741 : case 0:
419 : case 1:
420 : case 2:
421 : case 3:
422 741 : return 1;
423 :
424 0 : case 4:
425 : case 5:
426 : case 6:
427 : case 7:
428 0 : return 0;
429 :
430 0 : case 8:
431 0 : return 0;
432 :
433 0 : default:
434 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD4/8/9!");
435 : }
436 : }
437 :
438 0 : case INVALID_ELEM:
439 0 : return 0;
440 :
441 0 : default:
442 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
443 : }
444 : }
445 :
446 :
447 :
448 : // The second-order Szabo-Babuska shape functions
449 155520 : case SECOND:
450 : {
451 155520 : switch (t)
452 : {
453 0 : case NODEELEM:
454 0 : return 1;
455 :
456 : // The 1D Szabo-Babuska defined on a three-noded edge
457 5109 : case EDGE2:
458 : case EDGE3:
459 : {
460 5109 : switch (n)
461 : {
462 325 : case 0:
463 : case 1:
464 325 : return 1;
465 :
466 156 : case 2:
467 156 : return 1;
468 :
469 0 : default:
470 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for EDGE2/3!");
471 : }
472 : }
473 :
474 :
475 : // The 2D Szabo-Babuska defined on a 6-noded triangle
476 113946 : case TRI6:
477 : case TRI7:
478 : {
479 113946 : switch (n)
480 : {
481 4626 : case 0:
482 : case 1:
483 : case 2:
484 4626 : return 1;
485 :
486 4190 : case 3:
487 : case 4:
488 : case 5:
489 4190 : return 1;
490 :
491 1170 : case 6:
492 585 : libmesh_assert(t == TRI7);
493 1170 : return 0;
494 :
495 0 : default:
496 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for TRI6/7!");
497 : }
498 : }
499 :
500 :
501 : // The 2D tensor-product Szabo-Babuska defined on a
502 : // nine-noded quadrilateral.
503 36465 : case QUAD8:
504 : case QUADSHELL8:
505 : case QUAD9:
506 : case QUADSHELL9:
507 : {
508 36465 : switch (n)
509 : {
510 1469 : case 0:
511 : case 1:
512 : case 2:
513 : case 3:
514 1469 : return 1;
515 :
516 1365 : case 4:
517 : case 5:
518 : case 6:
519 : case 7:
520 1365 : return 1;
521 :
522 351 : case 8:
523 351 : return 1;
524 :
525 0 : default:
526 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD8/9!");
527 : }
528 : }
529 :
530 0 : case INVALID_ELEM:
531 0 : return 0;
532 :
533 0 : default:
534 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
535 : }
536 : }
537 :
538 :
539 :
540 : // The third-order Szabo-Babuska shape functions
541 136916 : case THIRD:
542 : {
543 136916 : switch (t)
544 : {
545 0 : case NODEELEM:
546 0 : return 1;
547 :
548 : // The 1D Szabo-Babuska defined on a three-noded edge
549 5109 : case EDGE2:
550 : case EDGE3:
551 : {
552 5109 : switch (n)
553 : {
554 325 : case 0:
555 : case 1:
556 325 : return 1;
557 :
558 312 : case 2:
559 312 : return 2;
560 :
561 0 : default:
562 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for EDGE2/3!");
563 : }
564 : }
565 :
566 :
567 : // The 2D Szabo-Babuska defined on a 6-noded triangle
568 51272 : case TRI7:
569 51272 : if (n == 6)
570 637 : return 1;
571 : libmesh_fallthrough();
572 : case TRI6:
573 : {
574 88140 : switch (n)
575 : {
576 4082 : case 0:
577 : case 1:
578 : case 2:
579 4082 : return 1;
580 :
581 7488 : case 3:
582 : case 4:
583 : case 5:
584 7488 : return 2;
585 :
586 0 : default:
587 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for TRI6/7!");
588 : }
589 : }
590 :
591 :
592 : // The 2D tensor-product Szabo-Babuska defined on a
593 : // nine-noded quadrilateral.
594 36465 : case QUAD8:
595 : case QUADSHELL8:
596 : case QUAD9:
597 : case QUADSHELL9:
598 : {
599 33280 : switch (n)
600 : {
601 1469 : case 0:
602 : case 1:
603 : case 2:
604 : case 3:
605 1469 : return 1;
606 :
607 1365 : case 4:
608 : case 5:
609 : case 6:
610 : case 7:
611 1365 : return 2;
612 :
613 351 : case 8:
614 351 : return 4;
615 :
616 0 : default:
617 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD8/9!");
618 : }
619 : }
620 :
621 0 : case INVALID_ELEM:
622 0 : return 0;
623 :
624 0 : default:
625 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
626 : }
627 : }
628 :
629 :
630 :
631 : // The fourth-order Szabo-Babuska shape functions
632 136916 : case FOURTH:
633 : {
634 136916 : switch (t)
635 : {
636 0 : case NODEELEM:
637 0 : return 1;
638 :
639 : // The 1D Szabo-Babuska defined on a three-noded edge
640 5109 : case EDGE2:
641 : case EDGE3:
642 : {
643 5109 : switch (n)
644 : {
645 325 : case 0:
646 : case 1:
647 325 : return 1;
648 :
649 312 : case 2:
650 312 : return 3;
651 :
652 0 : default:
653 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for EDGE2/3!");
654 : }
655 : }
656 :
657 :
658 : // The 2D Szabo-Babuska defined on a 6/7-noded triangle
659 51272 : case TRI7:
660 51272 : if (n == 6)
661 637 : return 3;
662 : libmesh_fallthrough();
663 : case TRI6:
664 : {
665 88140 : switch (n)
666 : {
667 4082 : case 0:
668 : case 1:
669 : case 2:
670 4082 : return 1;
671 :
672 7488 : case 3:
673 : case 4:
674 : case 5:
675 7488 : return 3;
676 :
677 0 : default:
678 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for TRI6!");
679 : }
680 : }
681 :
682 :
683 : // The 2D tensor-product Szabo-Babuska defined on a
684 : // nine-noded quadrilateral.
685 36465 : case QUAD8:
686 : case QUADSHELL8:
687 : case QUAD9:
688 : case QUADSHELL9:
689 : {
690 33280 : switch (n)
691 : {
692 1469 : case 0:
693 : case 1:
694 : case 2:
695 : case 3:
696 1469 : return 1;
697 :
698 1365 : case 4:
699 : case 5:
700 : case 6:
701 : case 7:
702 1365 : return 3;
703 :
704 351 : case 8:
705 351 : return 9;
706 :
707 0 : default:
708 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD8/9!");
709 : }
710 : }
711 :
712 0 : case INVALID_ELEM:
713 0 : return 0;
714 :
715 0 : default:
716 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
717 : }
718 : }
719 :
720 :
721 :
722 : // The fifth-order Szabo-Babuska shape functions
723 0 : case FIFTH:
724 : {
725 0 : switch (t)
726 : {
727 0 : case NODEELEM:
728 0 : return 1;
729 :
730 : // The 1D Szabo-Babuska defined on a three-noded edge
731 0 : case EDGE2:
732 : case EDGE3:
733 : {
734 0 : switch (n)
735 : {
736 0 : case 0:
737 : case 1:
738 0 : return 1;
739 :
740 0 : case 2:
741 0 : return 4;
742 :
743 0 : default:
744 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for EDGE2/3!");
745 : }
746 : }
747 :
748 :
749 : // The 2D Szabo-Babuska defined on a 6/7-noded triangle
750 0 : case TRI7:
751 0 : if (n == 6)
752 0 : return 6;
753 : libmesh_fallthrough();
754 : case TRI6:
755 : {
756 0 : switch (n)
757 : {
758 0 : case 0:
759 : case 1:
760 : case 2:
761 0 : return 1;
762 :
763 0 : case 3:
764 : case 4:
765 : case 5:
766 0 : return 4;
767 :
768 0 : default:
769 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for TRI6!");
770 : }
771 : }
772 :
773 :
774 : // The 2D tensor-product Szabo-Babuska defined on a
775 : // nine-noded quadrilateral.
776 0 : case QUAD8:
777 : case QUADSHELL8:
778 : case QUAD9:
779 : case QUADSHELL9:
780 : {
781 0 : switch (n)
782 : {
783 0 : case 0:
784 : case 1:
785 : case 2:
786 : case 3:
787 0 : return 1;
788 :
789 0 : case 4:
790 : case 5:
791 : case 6:
792 : case 7:
793 0 : return 4;
794 :
795 0 : case 8:
796 0 : return 16;
797 :
798 0 : default:
799 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD8/9!");
800 : }
801 : }
802 :
803 0 : case INVALID_ELEM:
804 0 : return 0;
805 :
806 0 : default:
807 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
808 : }
809 : }
810 :
811 :
812 :
813 : // The sixth-order Szabo-Babuska shape functions
814 0 : case SIXTH:
815 : {
816 0 : switch (t)
817 : {
818 0 : case NODEELEM:
819 0 : return 1;
820 :
821 : // The 1D Szabo-Babuska defined on a three-noded edge
822 0 : case EDGE2:
823 : case EDGE3:
824 : {
825 0 : switch (n)
826 : {
827 0 : case 0:
828 : case 1:
829 0 : return 1;
830 :
831 0 : case 2:
832 0 : return 5;
833 :
834 0 : default:
835 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for EDGE2/3!");
836 : }
837 : }
838 :
839 :
840 : // The 2D Szabo-Babuska defined on a 6/7-noded triangle
841 0 : case TRI7:
842 0 : if (n == 6)
843 0 : return 10;
844 : libmesh_fallthrough();
845 : case TRI6:
846 : {
847 0 : switch (n)
848 : {
849 0 : case 0:
850 : case 1:
851 : case 2:
852 0 : return 1;
853 :
854 0 : case 3:
855 : case 4:
856 : case 5:
857 0 : return 5;
858 :
859 0 : default:
860 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for TRI6!");
861 : }
862 : }
863 :
864 :
865 : // The 2D tensor-product Szabo-Babuska defined on a
866 : // nine-noded quadrilateral.
867 0 : case QUAD8:
868 : case QUADSHELL8:
869 : case QUAD9:
870 : case QUADSHELL9:
871 : {
872 0 : switch (n)
873 : {
874 0 : case 0:
875 : case 1:
876 : case 2:
877 : case 3:
878 0 : return 1;
879 :
880 0 : case 4:
881 : case 5:
882 : case 6:
883 : case 7:
884 0 : return 5;
885 :
886 0 : case 8:
887 0 : return 25;
888 :
889 0 : default:
890 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD8/9!");
891 : }
892 : }
893 :
894 0 : case INVALID_ELEM:
895 0 : return 0;
896 :
897 0 : default:
898 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
899 : }
900 : }
901 :
902 :
903 : // The seventh-order Szabo-Babuska shape functions
904 0 : case SEVENTH:
905 : {
906 0 : switch (t)
907 : {
908 0 : case NODEELEM:
909 0 : return 1;
910 :
911 : // The 1D Szabo-Babuska defined on a three-noded edge
912 0 : case EDGE2:
913 : case EDGE3:
914 : {
915 0 : switch (n)
916 : {
917 0 : case 0:
918 : case 1:
919 0 : return 1;
920 :
921 0 : case 2:
922 0 : return 6;
923 :
924 0 : default:
925 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for EDGE2/3!");
926 : }
927 : }
928 :
929 :
930 : // The 2D Szabo-Babuska defined on a 6/7-noded triangle
931 0 : case TRI7:
932 0 : if (n == 6)
933 0 : return 15;
934 : libmesh_fallthrough();
935 : case TRI6:
936 : {
937 0 : switch (n)
938 : {
939 0 : case 0:
940 : case 1:
941 : case 2:
942 0 : return 1;
943 :
944 0 : case 3:
945 : case 4:
946 : case 5:
947 0 : return 6;
948 :
949 0 : default:
950 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for TRI6!");
951 : }
952 : }
953 :
954 :
955 : // The 2D tensor-product Szabo-Babuska defined on a
956 : // nine-noded quadrilateral.
957 0 : case QUAD8:
958 : case QUADSHELL8:
959 : case QUAD9:
960 : case QUADSHELL9:
961 : {
962 0 : switch (n)
963 : {
964 0 : case 0:
965 : case 1:
966 : case 2:
967 : case 3:
968 0 : return 1;
969 :
970 0 : case 4:
971 : case 5:
972 : case 6:
973 : case 7:
974 0 : return 6;
975 :
976 0 : case 8:
977 0 : return 36;
978 :
979 0 : default:
980 0 : libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD8/9!");
981 : }
982 : }
983 :
984 0 : case INVALID_ELEM:
985 0 : return 0;
986 :
987 0 : default:
988 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
989 : }
990 : }
991 :
992 :
993 0 : default:
994 0 : libmesh_error_msg("ERROR: Invalid Order " << Utility::enum_to_string(o) << " selected for SZABAB FE family!");
995 : }
996 : } // szabab_n_dofs_at_node()
997 :
998 :
999 :
1000 111228 : unsigned int szabab_n_dofs_at_node(const Elem & e,
1001 : const Order o,
1002 : const unsigned int n)
1003 : {
1004 122148 : return szabab_n_dofs_at_node(e.type(), o, n);
1005 : }
1006 :
1007 :
1008 :
1009 62098 : unsigned int szabab_n_dofs_per_elem(const ElemType t, const Order o)
1010 : {
1011 62098 : switch (o)
1012 : {
1013 : // The first-order Szabo-Babuska shape functions
1014 6864 : case FIRST:
1015 : {
1016 6240 : switch (t)
1017 : {
1018 0 : case NODEELEM:
1019 0 : return 0;
1020 :
1021 : // The 1D Szabo-Babuska defined on an edge
1022 117 : case EDGE2:
1023 : case EDGE3:
1024 117 : return 0;
1025 :
1026 : // The 2D Szabo-Babuska defined on a triangle
1027 338 : case TRI3:
1028 : case TRISHELL3:
1029 : case TRI6:
1030 : case TRI7:
1031 338 : return 0;
1032 :
1033 : // The 2D tensor-product Szabo-Babuska defined on a
1034 : // quadrilateral.
1035 169 : case QUAD4:
1036 : case QUADSHELL4:
1037 : case QUAD8:
1038 : case QUADSHELL8:
1039 : case QUAD9:
1040 : case QUADSHELL9:
1041 169 : return 0;
1042 :
1043 0 : case INVALID_ELEM:
1044 0 : return 0;
1045 :
1046 0 : default:
1047 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
1048 : }
1049 : }
1050 :
1051 :
1052 :
1053 : // The second-order Szabo-Babuska shape functions
1054 20030 : case SECOND:
1055 : {
1056 18348 : switch (t)
1057 : {
1058 0 : case NODEELEM:
1059 0 : return 0;
1060 :
1061 : // The 1D Szabo-Babuska defined on a two-noded edge
1062 0 : case EDGE2:
1063 0 : return 1;
1064 :
1065 : // The 1D Szabo-Babuska defined on a three-noded edge
1066 143 : case EDGE3:
1067 143 : return 0;
1068 :
1069 : // The 2D Szabo-Babuska defined on a 6-noded triangle
1070 1240 : case TRI6:
1071 : case TRI7:
1072 1240 : return 0;
1073 :
1074 : // The 2D tensor-product Szabo-Babuska defined on a
1075 : // eight-noded quadrilateral.
1076 0 : case QUAD8:
1077 : case QUADSHELL8:
1078 0 : return 0;
1079 :
1080 : // The 2D tensor-product Szabo-Babuska defined on a
1081 : // nine-noded quadrilateral.
1082 299 : case QUAD9:
1083 : case QUADSHELL9:
1084 299 : return 0;
1085 :
1086 0 : case INVALID_ELEM:
1087 0 : return 0;
1088 :
1089 0 : default:
1090 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
1091 : }
1092 : }
1093 :
1094 :
1095 :
1096 : // The third-order Szabo-Babuska shape functions
1097 17602 : case THIRD:
1098 : {
1099 17602 : switch (t)
1100 : {
1101 0 : case NODEELEM:
1102 0 : return 0;
1103 :
1104 : // The 1D Szabo-Babuska defined on a two-noded edge
1105 0 : case EDGE2:
1106 0 : return 2;
1107 :
1108 : // The 1D Szabo-Babuska defined on a three-noded edge
1109 143 : case EDGE3:
1110 143 : return 0;
1111 :
1112 : // The 2D Szabo-Babuska defined on a 6-noded triangle
1113 1118 : case TRI6:
1114 1118 : return 1;
1115 :
1116 559 : case TRI7:
1117 559 : return 0;
1118 :
1119 : // The 2D tensor-product Szabo-Babuska defined on a
1120 : // eight-noded quadrilateral.
1121 0 : case QUAD8:
1122 : case QUADSHELL8:
1123 0 : return 4;
1124 :
1125 : // The 2D tensor-product Szabo-Babuska defined on a
1126 : // nine-noded quadrilateral.
1127 299 : case QUAD9:
1128 : case QUADSHELL9:
1129 299 : return 0;
1130 :
1131 0 : case INVALID_ELEM:
1132 0 : return 0;
1133 :
1134 0 : default:
1135 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
1136 : }
1137 : }
1138 :
1139 :
1140 :
1141 : // The fourth-order Szabo-Babuska shape functions
1142 17602 : case FOURTH:
1143 : {
1144 17602 : switch (t)
1145 : {
1146 0 : case NODEELEM:
1147 0 : return 0;
1148 :
1149 : // The 1D Szabo-Babuska defined on a two-noded edge
1150 0 : case EDGE2:
1151 0 : return 3;
1152 :
1153 : // The 1D Szabo-Babuska defined on a three-noded edge
1154 143 : case EDGE3:
1155 143 : return 0;
1156 :
1157 : // The 2D Szabo-Babuska defined on a 6-noded triangle
1158 1118 : case TRI6:
1159 1118 : return 3;
1160 :
1161 559 : case TRI7:
1162 559 : return 0;
1163 :
1164 : // The 2D tensor-product Szabo-Babuska defined on a
1165 : // eight-noded quadrilateral.
1166 0 : case QUAD8:
1167 : case QUADSHELL8:
1168 0 : return 9;
1169 :
1170 : // The 2D tensor-product Szabo-Babuska defined on a
1171 : // nine-noded quadrilateral.
1172 299 : case QUAD9:
1173 : case QUADSHELL9:
1174 299 : return 0;
1175 :
1176 0 : case INVALID_ELEM:
1177 0 : return 0;
1178 :
1179 0 : default:
1180 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
1181 : }
1182 : }
1183 :
1184 :
1185 :
1186 : // The fifth-order Szabo-Babuska shape functions
1187 0 : case FIFTH:
1188 : {
1189 0 : switch (t)
1190 : {
1191 0 : case NODEELEM:
1192 0 : return 0;
1193 :
1194 : // The 1D Szabo-Babuska defined on a two-noded edge
1195 0 : case EDGE2:
1196 0 : return 4;
1197 :
1198 : // The 1D Szabo-Babuska defined on a three-noded edge
1199 0 : case EDGE3:
1200 0 : return 0;
1201 :
1202 : // The 2D Szabo-Babuska defined on a 6-noded triangle
1203 0 : case TRI6:
1204 0 : return 6;
1205 :
1206 0 : case TRI7:
1207 0 : return 0;
1208 :
1209 : // The 2D tensor-product Szabo-Babuska defined on a
1210 : // eight-noded quadrilateral.
1211 0 : case QUAD8:
1212 : case QUADSHELL8:
1213 0 : return 16;
1214 :
1215 : // The 2D tensor-product Szabo-Babuska defined on a
1216 : // nine-noded quadrilateral.
1217 0 : case QUAD9:
1218 : case QUADSHELL9:
1219 0 : return 0;
1220 :
1221 0 : case INVALID_ELEM:
1222 0 : return 0;
1223 :
1224 0 : default:
1225 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
1226 : }
1227 : }
1228 :
1229 :
1230 : // The sixth-order Szabo-Babuska shape functions
1231 0 : case SIXTH:
1232 : {
1233 0 : switch (t)
1234 : {
1235 0 : case NODEELEM:
1236 0 : return 0;
1237 :
1238 : // The 1D Szabo-Babuska defined on a two-noded edge
1239 0 : case EDGE2:
1240 0 : return 5;
1241 :
1242 : // The 1D Szabo-Babuska defined on a three-noded edge
1243 0 : case EDGE3:
1244 0 : return 0;
1245 :
1246 : // The 2D Szabo-Babuska defined on a 6-noded triangle
1247 0 : case TRI6:
1248 0 : return 10;
1249 :
1250 0 : case TRI7:
1251 0 : return 0;
1252 :
1253 : // The 2D tensor-product Szabo-Babuska defined on a
1254 : // eight-noded quadrilateral.
1255 0 : case QUAD8:
1256 : case QUADSHELL8:
1257 0 : return 25;
1258 :
1259 : // The 2D tensor-product Szabo-Babuska defined on a
1260 : // nine-noded quadrilateral.
1261 0 : case QUAD9:
1262 : case QUADSHELL9:
1263 0 : return 0;
1264 :
1265 0 : case INVALID_ELEM:
1266 0 : return 0;
1267 :
1268 0 : default:
1269 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
1270 : }
1271 : }
1272 :
1273 :
1274 : // The seventh-order Szabo-Babuska shape functions
1275 0 : case SEVENTH:
1276 : {
1277 0 : switch (t)
1278 : {
1279 0 : case NODEELEM:
1280 0 : return 0;
1281 :
1282 : // The 1D Szabo-Babuska defined on a two-noded edge
1283 0 : case EDGE2:
1284 0 : return 6;
1285 :
1286 : // The 1D Szabo-Babuska defined on a three-noded edge
1287 0 : case EDGE3:
1288 0 : return 0;
1289 :
1290 : // The 2D Szabo-Babuska defined on a 6-noded triangle
1291 0 : case TRI6:
1292 0 : return 15;
1293 :
1294 0 : case TRI7:
1295 0 : return 0;
1296 :
1297 : // The 2D tensor-product Szabo-Babuska defined on a
1298 : // eight-noded quadrilateral.
1299 0 : case QUAD8:
1300 : case QUADSHELL8:
1301 0 : return 36;
1302 :
1303 : // The 2D tensor-product Szabo-Babuska defined on a
1304 : // nine-noded quadrilateral.
1305 0 : case QUAD9:
1306 : case QUADSHELL9:
1307 0 : return 0;
1308 :
1309 0 : case INVALID_ELEM:
1310 0 : return 0;
1311 :
1312 0 : default:
1313 0 : libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
1314 : }
1315 : }
1316 :
1317 :
1318 : // Otherwise no DOFS per element
1319 0 : default:
1320 0 : return 0;
1321 : }
1322 : } // szabab_n_dofs_per_elem
1323 :
1324 :
1325 :
1326 56672 : unsigned int szabab_n_dofs_per_elem(const Elem & e, const Order o)
1327 : {
1328 62098 : return szabab_n_dofs_per_elem(e.type(), o);
1329 : }
1330 :
1331 : } // anonymous namespace
1332 :
1333 :
1334 : // Instantiate (side_) nodal_soln() function for every dimension
1335 576 : LIBMESH_FE_NODAL_SOLN(SZABAB, szabab_nodal_soln)
1336 384 : LIBMESH_FE_SIDE_NODAL_SOLN(SZABAB)
1337 :
1338 :
1339 : // Full specialization of n_dofs() function for every dimension
1340 0 : template <> unsigned int FE<0,SZABAB>::n_dofs(const ElemType t, const Order o) { return szabab_n_dofs(t, o); }
1341 0 : template <> unsigned int FE<1,SZABAB>::n_dofs(const ElemType t, const Order o) { return szabab_n_dofs(t, o); }
1342 0 : template <> unsigned int FE<2,SZABAB>::n_dofs(const ElemType t, const Order o) { return szabab_n_dofs(t, o); }
1343 0 : template <> unsigned int FE<3,SZABAB>::n_dofs(const ElemType t, const Order o) { return szabab_n_dofs(t, o); }
1344 :
1345 0 : template <> unsigned int FE<0,SZABAB>::n_dofs(const Elem * e, const Order o) { return szabab_n_dofs(e, o); }
1346 3018 : template <> unsigned int FE<1,SZABAB>::n_dofs(const Elem * e, const Order o) { return szabab_n_dofs(e, o); }
1347 55948 : template <> unsigned int FE<2,SZABAB>::n_dofs(const Elem * e, const Order o) { return szabab_n_dofs(e, o); }
1348 0 : template <> unsigned int FE<3,SZABAB>::n_dofs(const Elem * e, const Order o) { return szabab_n_dofs(e, o); }
1349 :
1350 : // Full specialization of n_dofs_at_node() function for every dimension.
1351 0 : template <> unsigned int FE<0,SZABAB>::n_dofs_at_node(const ElemType t, const Order o, const unsigned int n) { return szabab_n_dofs_at_node(t, o, n); }
1352 12818 : template <> unsigned int FE<1,SZABAB>::n_dofs_at_node(const ElemType t, const Order o, const unsigned int n) { return szabab_n_dofs_at_node(t, o, n); }
1353 317500 : template <> unsigned int FE<2,SZABAB>::n_dofs_at_node(const ElemType t, const Order o, const unsigned int n) { return szabab_n_dofs_at_node(t, o, n); }
1354 0 : template <> unsigned int FE<3,SZABAB>::n_dofs_at_node(const ElemType t, const Order o, const unsigned int n) { return szabab_n_dofs_at_node(t, o, n); }
1355 :
1356 0 : template <> unsigned int FE<0,SZABAB>::n_dofs_at_node(const Elem & e, const Order o, const unsigned int n) { return szabab_n_dofs_at_node(e, o, n); }
1357 5265 : template <> unsigned int FE<1,SZABAB>::n_dofs_at_node(const Elem & e, const Order o, const unsigned int n) { return szabab_n_dofs_at_node(e, o, n); }
1358 116883 : template <> unsigned int FE<2,SZABAB>::n_dofs_at_node(const Elem & e, const Order o, const unsigned int n) { return szabab_n_dofs_at_node(e, o, n); }
1359 0 : template <> unsigned int FE<3,SZABAB>::n_dofs_at_node(const Elem & e, const Order o, const unsigned int n) { return szabab_n_dofs_at_node(e, o, n); }
1360 :
1361 : // Full specialization of n_dofs_per_elem() function for every dimension.
1362 0 : template <> unsigned int FE<0,SZABAB>::n_dofs_per_elem(const ElemType t, const Order o) { return szabab_n_dofs_per_elem(t, o); }
1363 0 : template <> unsigned int FE<1,SZABAB>::n_dofs_per_elem(const ElemType t, const Order o) { return szabab_n_dofs_per_elem(t, o); }
1364 0 : template <> unsigned int FE<2,SZABAB>::n_dofs_per_elem(const ElemType t, const Order o) { return szabab_n_dofs_per_elem(t, o); }
1365 0 : template <> unsigned int FE<3,SZABAB>::n_dofs_per_elem(const ElemType t, const Order o) { return szabab_n_dofs_per_elem(t, o); }
1366 :
1367 0 : template <> unsigned int FE<0,SZABAB>::n_dofs_per_elem(const Elem & e, const Order o) { return szabab_n_dofs_per_elem(e, o); }
1368 5980 : template <> unsigned int FE<1,SZABAB>::n_dofs_per_elem(const Elem & e, const Order o) { return szabab_n_dofs_per_elem(e, o); }
1369 56118 : template <> unsigned int FE<2,SZABAB>::n_dofs_per_elem(const Elem & e, const Order o) { return szabab_n_dofs_per_elem(e, o); }
1370 0 : template <> unsigned int FE<3,SZABAB>::n_dofs_per_elem(const Elem & e, const Order o) { return szabab_n_dofs_per_elem(e, o); }
1371 :
1372 : // Szabab FEMs are C^0 continuous
1373 0 : template <> FEContinuity FE<0,SZABAB>::get_continuity() const { return C_ZERO; }
1374 18976 : template <> FEContinuity FE<1,SZABAB>::get_continuity() const { return C_ZERO; }
1375 64215 : template <> FEContinuity FE<2,SZABAB>::get_continuity() const { return C_ZERO; }
1376 0 : template <> FEContinuity FE<3,SZABAB>::get_continuity() const { return C_ZERO; }
1377 :
1378 : // Szabab FEMs are hierarchic
1379 0 : template <> bool FE<0,SZABAB>::is_hierarchic() const { return true; }
1380 48 : template <> bool FE<1,SZABAB>::is_hierarchic() const { return true; }
1381 202 : template <> bool FE<2,SZABAB>::is_hierarchic() const { return true; }
1382 0 : template <> bool FE<3,SZABAB>::is_hierarchic() const { return true; }
1383 :
1384 : #ifdef LIBMESH_ENABLE_AMR
1385 : // compute_constraints() specializations are only needed for 2 and 3D
1386 : template <>
1387 3000 : void FE<2,SZABAB>::compute_constraints (DofConstraints & constraints,
1388 : DofMap & dof_map,
1389 : const unsigned int variable_number,
1390 : const Elem * elem)
1391 3000 : { compute_proj_constraints(constraints, dof_map, variable_number, elem); }
1392 :
1393 : template <>
1394 0 : void FE<3,SZABAB>::compute_constraints (DofConstraints & constraints,
1395 : DofMap & dof_map,
1396 : const unsigned int variable_number,
1397 : const Elem * elem)
1398 0 : { compute_proj_constraints(constraints, dof_map, variable_number, elem); }
1399 : #endif // #ifdef LIBMESH_ENABLE_AMR
1400 :
1401 : // Szabab shapes need reinit only for approximation orders >= 3,
1402 : // but we might reach that with p refinement
1403 0 : template <> bool FE<0,SZABAB>::shapes_need_reinit() const { return true; }
1404 612 : template <> bool FE<1,SZABAB>::shapes_need_reinit() const { return true; }
1405 8444 : template <> bool FE<2,SZABAB>::shapes_need_reinit() const { return true; }
1406 0 : template <> bool FE<3,SZABAB>::shapes_need_reinit() const { return true; }
1407 :
1408 : } // namespace libMesh
1409 :
1410 : #endif //LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
|
