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