libMesh
Loading...
Searching...
No Matches
fe_szabab.C
Go to the documentation of this file.
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
30namespace libMesh
31{
32
33// ------------------------------------------------------------
34// Szabo-Babuska-specific implementations, Steffen Petersen 2004
35
36// Anonymous namespace for local helper functions
37namespace {
38
39void 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 const unsigned int n_nodes = elem->n_nodes();
46
47 const ElemType elem_type = elem->type();
48
49 nodal_soln.resize(n_nodes);
50
51 const Order totalorder = order + add_p_level*elem->p_level();
52
53 // FEType object to be passed to various FEInterface functions below.
54 FEType fe_type(order, SZABAB);
55
56 switch (totalorder)
57 {
58 // Constant shape functions
59 case CONSTANT:
60 {
61 libmesh_assert_equal_to (elem_soln.size(), 1);
62
63 std::fill(nodal_soln.begin(), nodal_soln.end(), elem_soln[0]);
64
65 return;
66 }
67
68
69 // For other bases do interpolation at the nodes
70 // explicitly.
71 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 FEInterface::n_shape_functions(fe_type, elem);
81
82 std::vector<Point> refspace_nodes;
83 FEBase::get_refspace_nodes(elem_type,refspace_nodes);
84 libmesh_assert_equal_to (refspace_nodes.size(), n_nodes);
85 libmesh_assert_equal_to (elem_soln.size(), n_sf);
86
87 // Zero before summation
88 std::fill(nodal_soln.begin(), nodal_soln.end(), 0);
89
90 for (unsigned int n=0; n<n_nodes; n++)
91 // u_i = Sum (alpha_i phi_i)
92 for (unsigned int i=0; i<n_sf; i++)
93 nodal_soln[n] += elem_soln[i] *
94 FEInterface::shape(fe_type, elem, i, refspace_nodes[n]);
95
96 return;
97 }
98
99 default:
100 libmesh_error_msg("ERROR: Invalid total order " << totalorder);
101 }
102} // szabab_nodal_soln()
103
104
105
106unsigned int szabab_n_dofs(const ElemType t, const Order o)
107{
108 switch (o)
109 {
110 // Szabo-Babuska 1st-order polynomials.
111 case FIRST:
112 {
113 switch (t)
114 {
115 case NODEELEM:
116 return 1;
117
118 case EDGE2:
119 case EDGE3:
120 return 2;
121
122 case TRI3:
123 case TRISHELL3:
124 case TRI6:
125 case TRI7:
126 return 3;
127
128 case QUAD4:
129 case QUADSHELL4:
130 case QUAD8:
131 case QUADSHELL8:
132 case QUAD9:
133 case QUADSHELL9:
134 return 4;
135
136 case INVALID_ELEM:
137 return 0;
138
139 default:
140 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 case SECOND:
147 {
148 switch (t)
149 {
150 case NODEELEM:
151 return 1;
152
153 case EDGE2:
154 case EDGE3:
155 return 3;
156
157 case TRI6:
158 case TRI7:
159 return 6;
160
161 case QUAD8:
162 case QUADSHELL8:
163 return 8;
164 case QUAD9:
165 case QUADSHELL9:
166 return 9;
167
168 case INVALID_ELEM:
169 return 0;
170
171 default:
172 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 case THIRD:
179 {
180 switch (t)
181 {
182 case NODEELEM:
183 return 1;
184
185 case EDGE2:
186 case EDGE3:
187 return 4;
188
189 case TRI6:
190 case TRI7:
191 return 10;
192
193 case QUAD8:
194 case QUADSHELL8:
195 case QUAD9:
196 case QUADSHELL9:
197 return 16;
198
199 case INVALID_ELEM:
200 return 0;
201
202 default:
203 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 case FOURTH:
210 {
211 switch (t)
212 {
213 case NODEELEM:
214 return 1;
215
216 case EDGE2:
217 case EDGE3:
218 return 5;
219
220 case TRI6:
221 case TRI7:
222 return 15;
223
224 case QUAD8:
225 case QUADSHELL8:
226 case QUAD9:
227 case QUADSHELL9:
228 return 25;
229
230 case INVALID_ELEM:
231 return 0;
232
233 default:
234 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 case FIFTH:
241 {
242 switch (t)
243 {
244 case NODEELEM:
245 return 1;
246
247 case EDGE2:
248 case EDGE3:
249 return 6;
250
251 case TRI6:
252 case TRI7:
253 return 21;
254
255 case QUAD8:
256 case QUADSHELL8:
257 case QUAD9:
258 case QUADSHELL9:
259 return 36;
260
261 case INVALID_ELEM:
262 return 0;
263
264 default:
265 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 case SIXTH:
272 {
273 switch (t)
274 {
275 case NODEELEM:
276 return 1;
277
278 case EDGE2:
279 case EDGE3:
280 return 7;
281
282 case TRI6:
283 case TRI7:
284 return 28;
285
286 case QUAD8:
287 case QUADSHELL8:
288 case QUAD9:
289 case QUADSHELL9:
290 return 49;
291
292 case INVALID_ELEM:
293 return 0;
294
295 default:
296 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 case SEVENTH:
302 {
303 switch (t)
304 {
305 case NODEELEM:
306 return 1;
307
308 case EDGE2:
309 case EDGE3:
310 return 8;
311
312 case TRI6:
313 case TRI7:
314 return 36;
315
316 case QUAD8:
317 case QUADSHELL8:
318 case QUAD9:
319 case QUADSHELL9:
320 return 64;
321
322 case INVALID_ELEM:
323 return 0;
324
325 default:
326 libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
327 }
328 }
329
330
331 default:
332 libmesh_error_msg("ERROR: Invalid Order " << Utility::enum_to_string(o) << " selected for SZABAB FE family!");
333 }
334} // szabab_n_dofs()
335
336
337
338unsigned int szabab_n_dofs(const Elem * e, const Order o)
339{
341 return szabab_n_dofs(e->type(), o);
342}
343
344
345
346unsigned int szabab_n_dofs_at_node(const ElemType t,
347 const Order o,
348 const unsigned int n)
349{
350 switch (o)
351 {
352 // The first-order Szabo-Babuska shape functions
353 case FIRST:
354 {
355 switch (t)
356 {
357 case NODEELEM:
358 return 1;
359
360 // The 1D Szabo-Babuska defined on a three-noded edge
361 case EDGE2:
362 case EDGE3:
363 {
364 switch (n)
365 {
366 case 0:
367 case 1:
368 return 1;
369
370 case 2:
371 return 0;
372
373 default:
374 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 case TRI3:
381 case TRISHELL3:
382 case TRI6:
383 case TRI7:
384 {
385 switch (n)
386 {
387 case 0:
388 case 1:
389 case 2:
390 return 1;
391
392 case 3:
393 case 4:
394 case 5:
395 return 0;
396
397 case 6:
398 libmesh_assert(t == TRI7);
399 return 0;
400
401 default:
402 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 case QUAD4:
410 case QUADSHELL4:
411 case QUAD8:
412 case QUADSHELL8:
413 case QUAD9:
414 case QUADSHELL9:
415 {
416 switch (n)
417 {
418 case 0:
419 case 1:
420 case 2:
421 case 3:
422 return 1;
423
424 case 4:
425 case 5:
426 case 6:
427 case 7:
428 return 0;
429
430 case 8:
431 return 0;
432
433 default:
434 libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD4/8/9!");
435 }
436 }
437
438 case INVALID_ELEM:
439 return 0;
440
441 default:
442 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 case SECOND:
450 {
451 switch (t)
452 {
453 case NODEELEM:
454 return 1;
455
456 // The 1D Szabo-Babuska defined on a three-noded edge
457 case EDGE2:
458 case EDGE3:
459 {
460 switch (n)
461 {
462 case 0:
463 case 1:
464 return 1;
465
466 case 2:
467 return 1;
468
469 default:
470 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 case TRI6:
477 case TRI7:
478 {
479 switch (n)
480 {
481 case 0:
482 case 1:
483 case 2:
484 return 1;
485
486 case 3:
487 case 4:
488 case 5:
489 return 1;
490
491 case 6:
492 libmesh_assert(t == TRI7);
493 return 0;
494
495 default:
496 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 case QUAD8:
504 case QUADSHELL8:
505 case QUAD9:
506 case QUADSHELL9:
507 {
508 switch (n)
509 {
510 case 0:
511 case 1:
512 case 2:
513 case 3:
514 return 1;
515
516 case 4:
517 case 5:
518 case 6:
519 case 7:
520 return 1;
521
522 case 8:
523 return 1;
524
525 default:
526 libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD8/9!");
527 }
528 }
529
530 case INVALID_ELEM:
531 return 0;
532
533 default:
534 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 case THIRD:
542 {
543 switch (t)
544 {
545 case NODEELEM:
546 return 1;
547
548 // The 1D Szabo-Babuska defined on a three-noded edge
549 case EDGE2:
550 case EDGE3:
551 {
552 switch (n)
553 {
554 case 0:
555 case 1:
556 return 1;
557
558 case 2:
559 return 2;
560
561 default:
562 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 case TRI7:
569 if (n == 6)
570 return 1;
571 libmesh_fallthrough();
572 case TRI6:
573 {
574 switch (n)
575 {
576 case 0:
577 case 1:
578 case 2:
579 return 1;
580
581 case 3:
582 case 4:
583 case 5:
584 return 2;
585
586 default:
587 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 case QUAD8:
595 case QUADSHELL8:
596 case QUAD9:
597 case QUADSHELL9:
598 {
599 switch (n)
600 {
601 case 0:
602 case 1:
603 case 2:
604 case 3:
605 return 1;
606
607 case 4:
608 case 5:
609 case 6:
610 case 7:
611 return 2;
612
613 case 8:
614 return 4;
615
616 default:
617 libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD8/9!");
618 }
619 }
620
621 case INVALID_ELEM:
622 return 0;
623
624 default:
625 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 case FOURTH:
633 {
634 switch (t)
635 {
636 case NODEELEM:
637 return 1;
638
639 // The 1D Szabo-Babuska defined on a three-noded edge
640 case EDGE2:
641 case EDGE3:
642 {
643 switch (n)
644 {
645 case 0:
646 case 1:
647 return 1;
648
649 case 2:
650 return 3;
651
652 default:
653 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 case TRI7:
660 if (n == 6)
661 return 3;
662 libmesh_fallthrough();
663 case TRI6:
664 {
665 switch (n)
666 {
667 case 0:
668 case 1:
669 case 2:
670 return 1;
671
672 case 3:
673 case 4:
674 case 5:
675 return 3;
676
677 default:
678 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 case QUAD8:
686 case QUADSHELL8:
687 case QUAD9:
688 case QUADSHELL9:
689 {
690 switch (n)
691 {
692 case 0:
693 case 1:
694 case 2:
695 case 3:
696 return 1;
697
698 case 4:
699 case 5:
700 case 6:
701 case 7:
702 return 3;
703
704 case 8:
705 return 9;
706
707 default:
708 libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD8/9!");
709 }
710 }
711
712 case INVALID_ELEM:
713 return 0;
714
715 default:
716 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 case FIFTH:
724 {
725 switch (t)
726 {
727 case NODEELEM:
728 return 1;
729
730 // The 1D Szabo-Babuska defined on a three-noded edge
731 case EDGE2:
732 case EDGE3:
733 {
734 switch (n)
735 {
736 case 0:
737 case 1:
738 return 1;
739
740 case 2:
741 return 4;
742
743 default:
744 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 case TRI7:
751 if (n == 6)
752 return 6;
753 libmesh_fallthrough();
754 case TRI6:
755 {
756 switch (n)
757 {
758 case 0:
759 case 1:
760 case 2:
761 return 1;
762
763 case 3:
764 case 4:
765 case 5:
766 return 4;
767
768 default:
769 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 case QUAD8:
777 case QUADSHELL8:
778 case QUAD9:
779 case QUADSHELL9:
780 {
781 switch (n)
782 {
783 case 0:
784 case 1:
785 case 2:
786 case 3:
787 return 1;
788
789 case 4:
790 case 5:
791 case 6:
792 case 7:
793 return 4;
794
795 case 8:
796 return 16;
797
798 default:
799 libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD8/9!");
800 }
801 }
802
803 case INVALID_ELEM:
804 return 0;
805
806 default:
807 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 case SIXTH:
815 {
816 switch (t)
817 {
818 case NODEELEM:
819 return 1;
820
821 // The 1D Szabo-Babuska defined on a three-noded edge
822 case EDGE2:
823 case EDGE3:
824 {
825 switch (n)
826 {
827 case 0:
828 case 1:
829 return 1;
830
831 case 2:
832 return 5;
833
834 default:
835 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 case TRI7:
842 if (n == 6)
843 return 10;
844 libmesh_fallthrough();
845 case TRI6:
846 {
847 switch (n)
848 {
849 case 0:
850 case 1:
851 case 2:
852 return 1;
853
854 case 3:
855 case 4:
856 case 5:
857 return 5;
858
859 default:
860 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 case QUAD8:
868 case QUADSHELL8:
869 case QUAD9:
870 case QUADSHELL9:
871 {
872 switch (n)
873 {
874 case 0:
875 case 1:
876 case 2:
877 case 3:
878 return 1;
879
880 case 4:
881 case 5:
882 case 6:
883 case 7:
884 return 5;
885
886 case 8:
887 return 25;
888
889 default:
890 libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD8/9!");
891 }
892 }
893
894 case INVALID_ELEM:
895 return 0;
896
897 default:
898 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 case SEVENTH:
905 {
906 switch (t)
907 {
908 case NODEELEM:
909 return 1;
910
911 // The 1D Szabo-Babuska defined on a three-noded edge
912 case EDGE2:
913 case EDGE3:
914 {
915 switch (n)
916 {
917 case 0:
918 case 1:
919 return 1;
920
921 case 2:
922 return 6;
923
924 default:
925 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 case TRI7:
932 if (n == 6)
933 return 15;
934 libmesh_fallthrough();
935 case TRI6:
936 {
937 switch (n)
938 {
939 case 0:
940 case 1:
941 case 2:
942 return 1;
943
944 case 3:
945 case 4:
946 case 5:
947 return 6;
948
949 default:
950 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 case QUAD8:
958 case QUADSHELL8:
959 case QUAD9:
960 case QUADSHELL9:
961 {
962 switch (n)
963 {
964 case 0:
965 case 1:
966 case 2:
967 case 3:
968 return 1;
969
970 case 4:
971 case 5:
972 case 6:
973 case 7:
974 return 6;
975
976 case 8:
977 return 36;
978
979 default:
980 libmesh_error_msg("ERROR: Invalid node ID " << n << " selected for QUAD8/9!");
981 }
982 }
983
984 case INVALID_ELEM:
985 return 0;
986
987 default:
988 libmesh_error_msg("ERROR: Invalid ElemType " << Utility::enum_to_string(t) << " selected for SZABAB FE family!");
989 }
990 }
991
992
993 default:
994 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
1000unsigned int szabab_n_dofs_at_node(const Elem & e,
1001 const Order o,
1002 const unsigned int n)
1003{
1004 return szabab_n_dofs_at_node(e.type(), o, n);
1005}
1006
1007
1008
1009unsigned int szabab_n_dofs_per_elem(const ElemType t, const Order o)
1010{
1011 switch (o)
1012 {
1013 // The first-order Szabo-Babuska shape functions
1014 case FIRST:
1015 {
1016 switch (t)
1017 {
1018 case NODEELEM:
1019 return 0;
1020
1021 // The 1D Szabo-Babuska defined on an edge
1022 case EDGE2:
1023 case EDGE3:
1024 return 0;
1025
1026 // The 2D Szabo-Babuska defined on a triangle
1027 case TRI3:
1028 case TRISHELL3:
1029 case TRI6:
1030 case TRI7:
1031 return 0;
1032
1033 // The 2D tensor-product Szabo-Babuska defined on a
1034 // quadrilateral.
1035 case QUAD4:
1036 case QUADSHELL4:
1037 case QUAD8:
1038 case QUADSHELL8:
1039 case QUAD9:
1040 case QUADSHELL9:
1041 return 0;
1042
1043 case INVALID_ELEM:
1044 return 0;
1045
1046 default:
1047 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 case SECOND:
1055 {
1056 switch (t)
1057 {
1058 case NODEELEM:
1059 return 0;
1060
1061 // The 1D Szabo-Babuska defined on a two-noded edge
1062 case EDGE2:
1063 return 1;
1064
1065 // The 1D Szabo-Babuska defined on a three-noded edge
1066 case EDGE3:
1067 return 0;
1068
1069 // The 2D Szabo-Babuska defined on a 6-noded triangle
1070 case TRI6:
1071 case TRI7:
1072 return 0;
1073
1074 // The 2D tensor-product Szabo-Babuska defined on a
1075 // eight-noded quadrilateral.
1076 case QUAD8:
1077 case QUADSHELL8:
1078 return 0;
1079
1080 // The 2D tensor-product Szabo-Babuska defined on a
1081 // nine-noded quadrilateral.
1082 case QUAD9:
1083 case QUADSHELL9:
1084 return 0;
1085
1086 case INVALID_ELEM:
1087 return 0;
1088
1089 default:
1090 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 case THIRD:
1098 {
1099 switch (t)
1100 {
1101 case NODEELEM:
1102 return 0;
1103
1104 // The 1D Szabo-Babuska defined on a two-noded edge
1105 case EDGE2:
1106 return 2;
1107
1108 // The 1D Szabo-Babuska defined on a three-noded edge
1109 case EDGE3:
1110 return 0;
1111
1112 // The 2D Szabo-Babuska defined on a 6-noded triangle
1113 case TRI6:
1114 return 1;
1115
1116 case TRI7:
1117 return 0;
1118
1119 // The 2D tensor-product Szabo-Babuska defined on a
1120 // eight-noded quadrilateral.
1121 case QUAD8:
1122 case QUADSHELL8:
1123 return 4;
1124
1125 // The 2D tensor-product Szabo-Babuska defined on a
1126 // nine-noded quadrilateral.
1127 case QUAD9:
1128 case QUADSHELL9:
1129 return 0;
1130
1131 case INVALID_ELEM:
1132 return 0;
1133
1134 default:
1135 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 case FOURTH:
1143 {
1144 switch (t)
1145 {
1146 case NODEELEM:
1147 return 0;
1148
1149 // The 1D Szabo-Babuska defined on a two-noded edge
1150 case EDGE2:
1151 return 3;
1152
1153 // The 1D Szabo-Babuska defined on a three-noded edge
1154 case EDGE3:
1155 return 0;
1156
1157 // The 2D Szabo-Babuska defined on a 6-noded triangle
1158 case TRI6:
1159 return 3;
1160
1161 case TRI7:
1162 return 0;
1163
1164 // The 2D tensor-product Szabo-Babuska defined on a
1165 // eight-noded quadrilateral.
1166 case QUAD8:
1167 case QUADSHELL8:
1168 return 9;
1169
1170 // The 2D tensor-product Szabo-Babuska defined on a
1171 // nine-noded quadrilateral.
1172 case QUAD9:
1173 case QUADSHELL9:
1174 return 0;
1175
1176 case INVALID_ELEM:
1177 return 0;
1178
1179 default:
1180 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 case FIFTH:
1188 {
1189 switch (t)
1190 {
1191 case NODEELEM:
1192 return 0;
1193
1194 // The 1D Szabo-Babuska defined on a two-noded edge
1195 case EDGE2:
1196 return 4;
1197
1198 // The 1D Szabo-Babuska defined on a three-noded edge
1199 case EDGE3:
1200 return 0;
1201
1202 // The 2D Szabo-Babuska defined on a 6-noded triangle
1203 case TRI6:
1204 return 6;
1205
1206 case TRI7:
1207 return 0;
1208
1209 // The 2D tensor-product Szabo-Babuska defined on a
1210 // eight-noded quadrilateral.
1211 case QUAD8:
1212 case QUADSHELL8:
1213 return 16;
1214
1215 // The 2D tensor-product Szabo-Babuska defined on a
1216 // nine-noded quadrilateral.
1217 case QUAD9:
1218 case QUADSHELL9:
1219 return 0;
1220
1221 case INVALID_ELEM:
1222 return 0;
1223
1224 default:
1225 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 case SIXTH:
1232 {
1233 switch (t)
1234 {
1235 case NODEELEM:
1236 return 0;
1237
1238 // The 1D Szabo-Babuska defined on a two-noded edge
1239 case EDGE2:
1240 return 5;
1241
1242 // The 1D Szabo-Babuska defined on a three-noded edge
1243 case EDGE3:
1244 return 0;
1245
1246 // The 2D Szabo-Babuska defined on a 6-noded triangle
1247 case TRI6:
1248 return 10;
1249
1250 case TRI7:
1251 return 0;
1252
1253 // The 2D tensor-product Szabo-Babuska defined on a
1254 // eight-noded quadrilateral.
1255 case QUAD8:
1256 case QUADSHELL8:
1257 return 25;
1258
1259 // The 2D tensor-product Szabo-Babuska defined on a
1260 // nine-noded quadrilateral.
1261 case QUAD9:
1262 case QUADSHELL9:
1263 return 0;
1264
1265 case INVALID_ELEM:
1266 return 0;
1267
1268 default:
1269 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 case SEVENTH:
1276 {
1277 switch (t)
1278 {
1279 case NODEELEM:
1280 return 0;
1281
1282 // The 1D Szabo-Babuska defined on a two-noded edge
1283 case EDGE2:
1284 return 6;
1285
1286 // The 1D Szabo-Babuska defined on a three-noded edge
1287 case EDGE3:
1288 return 0;
1289
1290 // The 2D Szabo-Babuska defined on a 6-noded triangle
1291 case TRI6:
1292 return 15;
1293
1294 case TRI7:
1295 return 0;
1296
1297 // The 2D tensor-product Szabo-Babuska defined on a
1298 // eight-noded quadrilateral.
1299 case QUAD8:
1300 case QUADSHELL8:
1301 return 36;
1302
1303 // The 2D tensor-product Szabo-Babuska defined on a
1304 // nine-noded quadrilateral.
1305 case QUAD9:
1306 case QUADSHELL9:
1307 return 0;
1308
1309 case INVALID_ELEM:
1310 return 0;
1311
1312 default:
1313 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 default:
1320 return 0;
1321 }
1322} // szabab_n_dofs_per_elem
1323
1324
1325
1326unsigned int szabab_n_dofs_per_elem(const Elem & e, const Order o)
1327{
1328 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
1337
1338
1339// Full specialization of n_dofs() function for every dimension
1340template <> unsigned int FE<0,SZABAB>::n_dofs(const ElemType t, const Order o) { return szabab_n_dofs(t, o); }
1341template <> unsigned int FE<1,SZABAB>::n_dofs(const ElemType t, const Order o) { return szabab_n_dofs(t, o); }
1342template <> unsigned int FE<2,SZABAB>::n_dofs(const ElemType t, const Order o) { return szabab_n_dofs(t, o); }
1343template <> unsigned int FE<3,SZABAB>::n_dofs(const ElemType t, const Order o) { return szabab_n_dofs(t, o); }
1344
1345template <> unsigned int FE<0,SZABAB>::n_dofs(const Elem * e, const Order o) { return szabab_n_dofs(e, o); }
1346template <> unsigned int FE<1,SZABAB>::n_dofs(const Elem * e, const Order o) { return szabab_n_dofs(e, o); }
1347template <> unsigned int FE<2,SZABAB>::n_dofs(const Elem * e, const Order o) { return szabab_n_dofs(e, o); }
1348template <> 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.
1351template <> 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); }
1352template <> 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); }
1353template <> 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); }
1354template <> 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
1356template <> 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); }
1357template <> 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); }
1358template <> 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); }
1359template <> 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.
1362template <> unsigned int FE<0,SZABAB>::n_dofs_per_elem(const ElemType t, const Order o) { return szabab_n_dofs_per_elem(t, o); }
1363template <> unsigned int FE<1,SZABAB>::n_dofs_per_elem(const ElemType t, const Order o) { return szabab_n_dofs_per_elem(t, o); }
1364template <> unsigned int FE<2,SZABAB>::n_dofs_per_elem(const ElemType t, const Order o) { return szabab_n_dofs_per_elem(t, o); }
1365template <> unsigned int FE<3,SZABAB>::n_dofs_per_elem(const ElemType t, const Order o) { return szabab_n_dofs_per_elem(t, o); }
1366
1367template <> unsigned int FE<0,SZABAB>::n_dofs_per_elem(const Elem & e, const Order o) { return szabab_n_dofs_per_elem(e, o); }
1368template <> unsigned int FE<1,SZABAB>::n_dofs_per_elem(const Elem & e, const Order o) { return szabab_n_dofs_per_elem(e, o); }
1369template <> unsigned int FE<2,SZABAB>::n_dofs_per_elem(const Elem & e, const Order o) { return szabab_n_dofs_per_elem(e, o); }
1370template <> 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
1373template <> FEContinuity FE<0,SZABAB>::get_continuity() const { return C_ZERO; }
1374template <> FEContinuity FE<1,SZABAB>::get_continuity() const { return C_ZERO; }
1375template <> FEContinuity FE<2,SZABAB>::get_continuity() const { return C_ZERO; }
1376template <> FEContinuity FE<3,SZABAB>::get_continuity() const { return C_ZERO; }
1377
1378// Szabab FEMs are hierarchic
1379template <> bool FE<0,SZABAB>::is_hierarchic() const { return true; }
1380template <> bool FE<1,SZABAB>::is_hierarchic() const { return true; }
1381template <> bool FE<2,SZABAB>::is_hierarchic() const { return true; }
1382template <> 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
1386template <>
1388 DofMap & dof_map,
1389 const unsigned int variable_number,
1390 const Elem * elem)
1391{ compute_proj_constraints(constraints, dof_map, variable_number, elem); }
1392
1393template <>
1395 DofMap & dof_map,
1396 const unsigned int variable_number,
1397 const Elem * elem)
1398{ 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
1403template <> bool FE<0,SZABAB>::shapes_need_reinit() const { return true; }
1404template <> bool FE<1,SZABAB>::shapes_need_reinit() const { return true; }
1405template <> bool FE<2,SZABAB>::shapes_need_reinit() const { return true; }
1406template <> bool FE<3,SZABAB>::shapes_need_reinit() const { return true; }
1407
1408} // namespace libMesh
1409
1410#endif //LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
The constraint matrix storage format.
Definition dof_map.h:112
This class handles the numbering of degrees of freedom on a mesh.
Definition dof_map.h:181
This is the base class from which all geometric element types are derived.
Definition elem.h:96
static void get_refspace_nodes(const ElemType t, std::vector< Point > &nodes)
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
static unsigned int n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
virtual bool is_hierarchic() const override
static unsigned int n_dofs_per_elem(const ElemType t, const Order o)
virtual FEContinuity get_continuity() const override
static void compute_constraints(DofConstraints &constraints, DofMap &dof_map, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for non-conforming adapted meshes) corresponding to var...
static unsigned int n_dofs(const ElemType t, const Order o)
static unsigned int n_dofs_at_node(const ElemType t, const Order o, const unsigned int n)
virtual bool shapes_need_reinit() const override
std::string enum_to_string(const T e)
The libMesh namespace provides an interface to certain functionality in the library.
LIBMESH_FE_NODAL_SOLN(BERNSTEIN, bernstein_nodal_soln)
ElemType
Defines an enum for geometric element types.
libmesh_assert(ctx)
LIBMESH_FE_SIDE_NODAL_SOLN(HIERARCHIC_VEC)
const dof_id_type n_nodes
Definition tecplot_io.C:67