libMesh
fe_lagrange_shape_2D.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 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 // C++ includes
20 
21 // Local includes
22 #include "libmesh/fe.h"
23 #include "libmesh/elem.h"
24 
25 namespace libMesh
26 {
27 
28 
29 
30 
31 template <>
33  const Order order,
34  const unsigned int i,
35  const Point & p)
36 {
37 #if LIBMESH_DIM > 1
38 
39  switch (order)
40  {
41  // linear Lagrange shape functions
42  case FIRST:
43  {
44  switch (type)
45  {
46  case QUAD4:
47  case QUADSHELL4:
48  case QUAD8:
49  case QUADSHELL8:
50  case QUAD9:
51  {
52  // Compute quad shape functions as a tensor-product
53  const Real xi = p(0);
54  const Real eta = p(1);
55 
56  libmesh_assert_less (i, 4);
57 
58  // 0 1 2 3
59  static const unsigned int i0[] = {0, 1, 1, 0};
60  static const unsigned int i1[] = {0, 0, 1, 1};
61 
62  return (FE<1,LAGRANGE>::shape(EDGE2, FIRST, i0[i], xi)*
63  FE<1,LAGRANGE>::shape(EDGE2, FIRST, i1[i], eta));
64  }
65 
66  case TRI3:
67  case TRISHELL3:
68  case TRI6:
69  {
70  const Real zeta1 = p(0);
71  const Real zeta2 = p(1);
72  const Real zeta0 = 1. - zeta1 - zeta2;
73 
74  libmesh_assert_less (i, 3);
75 
76  switch(i)
77  {
78  case 0:
79  return zeta0;
80 
81  case 1:
82  return zeta1;
83 
84  case 2:
85  return zeta2;
86 
87  default:
88  libmesh_error_msg("Invalid shape function index i = " << i);
89  }
90  }
91 
92  default:
93  libmesh_error_msg("ERROR: Unsupported 2D element type: " << type);
94  }
95  }
96 
97 
98  // quadratic Lagrange shape functions
99  case SECOND:
100  {
101  switch (type)
102  {
103  case QUAD8:
104  case QUADSHELL8:
105  {
106  const Real xi = p(0);
107  const Real eta = p(1);
108 
109  libmesh_assert_less (i, 8);
110 
111  switch (i)
112  {
113  case 0:
114  return .25*(1. - xi)*(1. - eta)*(-1. - xi - eta);
115 
116  case 1:
117  return .25*(1. + xi)*(1. - eta)*(-1. + xi - eta);
118 
119  case 2:
120  return .25*(1. + xi)*(1. + eta)*(-1. + xi + eta);
121 
122  case 3:
123  return .25*(1. - xi)*(1. + eta)*(-1. - xi + eta);
124 
125  case 4:
126  return .5*(1. - xi*xi)*(1. - eta);
127 
128  case 5:
129  return .5*(1. + xi)*(1. - eta*eta);
130 
131  case 6:
132  return .5*(1. - xi*xi)*(1. + eta);
133 
134  case 7:
135  return .5*(1. - xi)*(1. - eta*eta);
136 
137  default:
138  libmesh_error_msg("Invalid shape function index i = " << i);
139  }
140  }
141 
142  case QUAD9:
143  {
144  // Compute quad shape functions as a tensor-product
145  const Real xi = p(0);
146  const Real eta = p(1);
147 
148  libmesh_assert_less (i, 9);
149 
150  // 0 1 2 3 4 5 6 7 8
151  static const unsigned int i0[] = {0, 1, 1, 0, 2, 1, 2, 0, 2};
152  static const unsigned int i1[] = {0, 0, 1, 1, 0, 2, 1, 2, 2};
153 
154  return (FE<1,LAGRANGE>::shape(EDGE3, SECOND, i0[i], xi)*
155  FE<1,LAGRANGE>::shape(EDGE3, SECOND, i1[i], eta));
156  }
157 
158  case TRI6:
159  {
160  const Real zeta1 = p(0);
161  const Real zeta2 = p(1);
162  const Real zeta0 = 1. - zeta1 - zeta2;
163 
164  libmesh_assert_less (i, 6);
165 
166  switch(i)
167  {
168  case 0:
169  return 2.*zeta0*(zeta0-0.5);
170 
171  case 1:
172  return 2.*zeta1*(zeta1-0.5);
173 
174  case 2:
175  return 2.*zeta2*(zeta2-0.5);
176 
177  case 3:
178  return 4.*zeta0*zeta1;
179 
180  case 4:
181  return 4.*zeta1*zeta2;
182 
183  case 5:
184  return 4.*zeta2*zeta0;
185 
186  default:
187  libmesh_error_msg("Invalid shape function index i = " << i);
188  }
189  }
190 
191  default:
192  libmesh_error_msg("ERROR: Unsupported 2D element type: " << type);
193  }
194  }
195 
196 
197 
198  // unsupported order
199  default:
200  libmesh_error_msg("ERROR: Unsupported 2D FE order: " << order);
201  }
202 
203  libmesh_error_msg("We'll never get here!");
204 #endif // LIBMESH_DIM > 1
205  return 0.;
206 }
207 
208 
209 
210 template <>
212  const Order order,
213  const unsigned int i,
214  const Point & p)
215 {
216  libmesh_assert(elem);
217 
218  // call the orientation-independent shape functions
219  return FE<2,LAGRANGE>::shape(elem->type(), static_cast<Order>(order + elem->p_level()), i, p);
220 }
221 
222 
223 
224 template <>
226  const Order order,
227  const unsigned int i,
228  const unsigned int j,
229  const Point & p)
230 {
231 #if LIBMESH_DIM > 1
232 
233 
234  libmesh_assert_less (j, 2);
235 
236  switch (order)
237  {
238  // linear Lagrange shape functions
239  case FIRST:
240  {
241  switch (type)
242  {
243  case QUAD4:
244  case QUADSHELL4:
245  case QUAD8:
246  case QUADSHELL8:
247  case QUAD9:
248  {
249  // Compute quad shape functions as a tensor-product
250  const Real xi = p(0);
251  const Real eta = p(1);
252 
253  libmesh_assert_less (i, 4);
254 
255  // 0 1 2 3
256  static const unsigned int i0[] = {0, 1, 1, 0};
257  static const unsigned int i1[] = {0, 0, 1, 1};
258 
259  switch (j)
260  {
261  // d()/dxi
262  case 0:
263  return (FE<1,LAGRANGE>::shape_deriv(EDGE2, FIRST, i0[i], 0, xi)*
264  FE<1,LAGRANGE>::shape (EDGE2, FIRST, i1[i], eta));
265 
266  // d()/deta
267  case 1:
268  return (FE<1,LAGRANGE>::shape (EDGE2, FIRST, i0[i], xi)*
269  FE<1,LAGRANGE>::shape_deriv(EDGE2, FIRST, i1[i], 0, eta));
270 
271  default:
272  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
273  }
274  }
275 
276  case TRI3:
277  case TRISHELL3:
278  case TRI6:
279  {
280  libmesh_assert_less (i, 3);
281 
282  const Real dzeta0dxi = -1.;
283  const Real dzeta1dxi = 1.;
284  const Real dzeta2dxi = 0.;
285 
286  const Real dzeta0deta = -1.;
287  const Real dzeta1deta = 0.;
288  const Real dzeta2deta = 1.;
289 
290  switch (j)
291  {
292  // d()/dxi
293  case 0:
294  {
295  switch(i)
296  {
297  case 0:
298  return dzeta0dxi;
299 
300  case 1:
301  return dzeta1dxi;
302 
303  case 2:
304  return dzeta2dxi;
305 
306  default:
307  libmesh_error_msg("Invalid shape function index i = " << i);
308  }
309  }
310  // d()/deta
311  case 1:
312  {
313  switch(i)
314  {
315  case 0:
316  return dzeta0deta;
317 
318  case 1:
319  return dzeta1deta;
320 
321  case 2:
322  return dzeta2deta;
323 
324  default:
325  libmesh_error_msg("Invalid shape function index i = " << i);
326  }
327  }
328  default:
329  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
330  }
331  }
332 
333  default:
334  libmesh_error_msg("ERROR: Unsupported 2D element type: " << type);
335  }
336  }
337 
338 
339  // quadratic Lagrange shape functions
340  case SECOND:
341  {
342  switch (type)
343  {
344  case QUAD8:
345  case QUADSHELL8:
346  {
347  const Real xi = p(0);
348  const Real eta = p(1);
349 
350  libmesh_assert_less (i, 8);
351 
352  switch (j)
353  {
354  // d/dxi
355  case 0:
356  switch (i)
357  {
358  case 0:
359  return .25*(1. - eta)*((1. - xi)*(-1.) +
360  (-1.)*(-1. - xi - eta));
361 
362  case 1:
363  return .25*(1. - eta)*((1. + xi)*(1.) +
364  (1.)*(-1. + xi - eta));
365 
366  case 2:
367  return .25*(1. + eta)*((1. + xi)*(1.) +
368  (1.)*(-1. + xi + eta));
369 
370  case 3:
371  return .25*(1. + eta)*((1. - xi)*(-1.) +
372  (-1.)*(-1. - xi + eta));
373 
374  case 4:
375  return .5*(-2.*xi)*(1. - eta);
376 
377  case 5:
378  return .5*(1.)*(1. - eta*eta);
379 
380  case 6:
381  return .5*(-2.*xi)*(1. + eta);
382 
383  case 7:
384  return .5*(-1.)*(1. - eta*eta);
385 
386  default:
387  libmesh_error_msg("Invalid shape function index i = " << i);
388  }
389 
390  // d/deta
391  case 1:
392  switch (i)
393  {
394  case 0:
395  return .25*(1. - xi)*((1. - eta)*(-1.) +
396  (-1.)*(-1. - xi - eta));
397 
398  case 1:
399  return .25*(1. + xi)*((1. - eta)*(-1.) +
400  (-1.)*(-1. + xi - eta));
401 
402  case 2:
403  return .25*(1. + xi)*((1. + eta)*(1.) +
404  (1.)*(-1. + xi + eta));
405 
406  case 3:
407  return .25*(1. - xi)*((1. + eta)*(1.) +
408  (1.)*(-1. - xi + eta));
409 
410  case 4:
411  return .5*(1. - xi*xi)*(-1.);
412 
413  case 5:
414  return .5*(1. + xi)*(-2.*eta);
415 
416  case 6:
417  return .5*(1. - xi*xi)*(1.);
418 
419  case 7:
420  return .5*(1. - xi)*(-2.*eta);
421 
422  default:
423  libmesh_error_msg("Invalid shape function index i = " << i);
424  }
425 
426  default:
427  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
428  }
429  }
430 
431  case QUAD9:
432  {
433  // Compute quad shape functions as a tensor-product
434  const Real xi = p(0);
435  const Real eta = p(1);
436 
437  libmesh_assert_less (i, 9);
438 
439  // 0 1 2 3 4 5 6 7 8
440  static const unsigned int i0[] = {0, 1, 1, 0, 2, 1, 2, 0, 2};
441  static const unsigned int i1[] = {0, 0, 1, 1, 0, 2, 1, 2, 2};
442 
443  switch (j)
444  {
445  // d()/dxi
446  case 0:
447  return (FE<1,LAGRANGE>::shape_deriv(EDGE3, SECOND, i0[i], 0, xi)*
448  FE<1,LAGRANGE>::shape (EDGE3, SECOND, i1[i], eta));
449 
450  // d()/deta
451  case 1:
452  return (FE<1,LAGRANGE>::shape (EDGE3, SECOND, i0[i], xi)*
453  FE<1,LAGRANGE>::shape_deriv(EDGE3, SECOND, i1[i], 0, eta));
454 
455  default:
456  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
457  }
458  }
459 
460  case TRI6:
461  {
462  libmesh_assert_less (i, 6);
463 
464  const Real zeta1 = p(0);
465  const Real zeta2 = p(1);
466  const Real zeta0 = 1. - zeta1 - zeta2;
467 
468  const Real dzeta0dxi = -1.;
469  const Real dzeta1dxi = 1.;
470  const Real dzeta2dxi = 0.;
471 
472  const Real dzeta0deta = -1.;
473  const Real dzeta1deta = 0.;
474  const Real dzeta2deta = 1.;
475 
476  switch(j)
477  {
478  case 0:
479  {
480  switch(i)
481  {
482  case 0:
483  return (4.*zeta0-1.)*dzeta0dxi;
484 
485  case 1:
486  return (4.*zeta1-1.)*dzeta1dxi;
487 
488  case 2:
489  return (4.*zeta2-1.)*dzeta2dxi;
490 
491  case 3:
492  return 4.*zeta1*dzeta0dxi + 4.*zeta0*dzeta1dxi;
493 
494  case 4:
495  return 4.*zeta2*dzeta1dxi + 4.*zeta1*dzeta2dxi;
496 
497  case 5:
498  return 4.*zeta2*dzeta0dxi + 4*zeta0*dzeta2dxi;
499 
500  default:
501  libmesh_error_msg("Invalid shape function index i = " << i);
502  }
503  }
504 
505  case 1:
506  {
507  switch(i)
508  {
509  case 0:
510  return (4.*zeta0-1.)*dzeta0deta;
511 
512  case 1:
513  return (4.*zeta1-1.)*dzeta1deta;
514 
515  case 2:
516  return (4.*zeta2-1.)*dzeta2deta;
517 
518  case 3:
519  return 4.*zeta1*dzeta0deta + 4.*zeta0*dzeta1deta;
520 
521  case 4:
522  return 4.*zeta2*dzeta1deta + 4.*zeta1*dzeta2deta;
523 
524  case 5:
525  return 4.*zeta2*dzeta0deta + 4*zeta0*dzeta2deta;
526 
527  default:
528  libmesh_error_msg("Invalid shape function index i = " << i);
529  }
530  }
531  default:
532  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
533  }
534  }
535 
536  default:
537  libmesh_error_msg("ERROR: Unsupported 2D element type: " << type);
538  }
539  }
540 
541  // unsupported order
542  default:
543  libmesh_error_msg("ERROR: Unsupported 2D FE order: " << order);
544  }
545 
546 
547  libmesh_error_msg("We'll never get here!");
548 #endif // LIBMESH_DIM > 1
549  return 0.;
550 }
551 
552 
553 
554 template <>
556  const Order order,
557  const unsigned int i,
558  const unsigned int j,
559  const Point & p)
560 {
561  libmesh_assert(elem);
562 
563 
564  // call the orientation-independent shape functions
565  return FE<2,LAGRANGE>::shape_deriv(elem->type(), static_cast<Order>(order + elem->p_level()), i, j, p);
566 }
567 
568 
569 
570 
571 template <>
573  const Order order,
574  const unsigned int i,
575  const unsigned int j,
576  const Point & p)
577 {
578 #if LIBMESH_DIM > 1
579 
580  // j = 0 ==> d^2 phi / dxi^2
581  // j = 1 ==> d^2 phi / dxi deta
582  // j = 2 ==> d^2 phi / deta^2
583  libmesh_assert_less (j, 3);
584 
585  switch (order)
586  {
587  // linear Lagrange shape functions
588  case FIRST:
589  {
590  switch (type)
591  {
592  case QUAD4:
593  case QUADSHELL4:
594  case QUAD8:
595  case QUADSHELL8:
596  case QUAD9:
597  {
598  // Compute quad shape functions as a tensor-product
599  const Real xi = p(0);
600  const Real eta = p(1);
601 
602  libmesh_assert_less (i, 4);
603 
604  // 0 1 2 3
605  static const unsigned int i0[] = {0, 1, 1, 0};
606  static const unsigned int i1[] = {0, 0, 1, 1};
607 
608  switch (j)
609  {
610  // d^2() / dxi^2
611  case 0:
612  return 0.;
613 
614  // d^2() / dxi deta
615  case 1:
616  return (FE<1,LAGRANGE>::shape_deriv(EDGE2, FIRST, i0[i], 0, xi)*
617  FE<1,LAGRANGE>::shape_deriv(EDGE2, FIRST, i1[i], 0, eta));
618 
619  // d^2() / deta^2
620  case 2:
621  return 0.;
622 
623  default:
624  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
625  }
626  }
627 
628  case TRI3:
629  case TRISHELL3:
630  case TRI6:
631  {
632  // All second derivatives for linear triangles are zero.
633  return 0.;
634  }
635 
636  default:
637  libmesh_error_msg("ERROR: Unsupported 2D element type: " << type);
638 
639  } // end switch (type)
640  } // end case FIRST
641 
642 
643  // quadratic Lagrange shape functions
644  case SECOND:
645  {
646  switch (type)
647  {
648  case QUAD8:
649  case QUADSHELL8:
650  {
651  const Real xi = p(0);
652  const Real eta = p(1);
653 
654  libmesh_assert_less (j, 3);
655 
656  switch (j)
657  {
658  // d^2() / dxi^2
659  case 0:
660  {
661  switch (i)
662  {
663  case 0:
664  case 1:
665  return 0.5*(1.-eta);
666 
667  case 2:
668  case 3:
669  return 0.5*(1.+eta);
670 
671  case 4:
672  return eta - 1.;
673 
674  case 5:
675  case 7:
676  return 0.0;
677 
678  case 6:
679  return -1. - eta;
680 
681  default:
682  libmesh_error_msg("Invalid shape function index i = " << i);
683  }
684  }
685 
686  // d^2() / dxi deta
687  case 1:
688  {
689  switch (i)
690  {
691  case 0:
692  return 0.25*( 1. - 2.*xi - 2.*eta);
693 
694  case 1:
695  return 0.25*(-1. - 2.*xi + 2.*eta);
696 
697  case 2:
698  return 0.25*( 1. + 2.*xi + 2.*eta);
699 
700  case 3:
701  return 0.25*(-1. + 2.*xi - 2.*eta);
702 
703  case 4:
704  return xi;
705 
706  case 5:
707  return -eta;
708 
709  case 6:
710  return -xi;
711 
712  case 7:
713  return eta;
714 
715  default:
716  libmesh_error_msg("Invalid shape function index i = " << i);
717  }
718  }
719 
720  // d^2() / deta^2
721  case 2:
722  {
723  switch (i)
724  {
725  case 0:
726  case 3:
727  return 0.5*(1.-xi);
728 
729  case 1:
730  case 2:
731  return 0.5*(1.+xi);
732 
733  case 4:
734  case 6:
735  return 0.0;
736 
737  case 5:
738  return -1.0 - xi;
739 
740  case 7:
741  return xi - 1.0;
742 
743  default:
744  libmesh_error_msg("Invalid shape function index i = " << i);
745  }
746  }
747 
748  default:
749  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
750  } // end switch (j)
751  } // end case QUAD8
752 
753  case QUAD9:
754  {
755  // Compute QUAD9 second derivatives as tensor product
756  const Real xi = p(0);
757  const Real eta = p(1);
758 
759  libmesh_assert_less (i, 9);
760 
761  // 0 1 2 3 4 5 6 7 8
762  static const unsigned int i0[] = {0, 1, 1, 0, 2, 1, 2, 0, 2};
763  static const unsigned int i1[] = {0, 0, 1, 1, 0, 2, 1, 2, 2};
764 
765  switch (j)
766  {
767  // d^2() / dxi^2
768  case 0:
769  return (FE<1,LAGRANGE>::shape_second_deriv(EDGE3, SECOND, i0[i], 0, xi)*
770  FE<1,LAGRANGE>::shape (EDGE3, SECOND, i1[i], eta));
771 
772  // d^2() / dxi deta
773  case 1:
774  return (FE<1,LAGRANGE>::shape_deriv(EDGE3, SECOND, i0[i], 0, xi)*
775  FE<1,LAGRANGE>::shape_deriv(EDGE3, SECOND, i1[i], 0, eta));
776 
777  // d^2() / deta^2
778  case 2:
779  return (FE<1,LAGRANGE>::shape (EDGE3, SECOND, i0[i], xi)*
781 
782  default:
783  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
784  } // end switch (j)
785  } // end case QUAD9
786 
787  case TRI6:
788  {
789  const Real dzeta0dxi = -1.;
790  const Real dzeta1dxi = 1.;
791  const Real dzeta2dxi = 0.;
792 
793  const Real dzeta0deta = -1.;
794  const Real dzeta1deta = 0.;
795  const Real dzeta2deta = 1.;
796 
797  libmesh_assert_less (j, 3);
798 
799  switch (j)
800  {
801  // d^2() / dxi^2
802  case 0:
803  {
804  switch (i)
805  {
806  case 0:
807  return 4.*dzeta0dxi*dzeta0dxi;
808 
809  case 1:
810  return 4.*dzeta1dxi*dzeta1dxi;
811 
812  case 2:
813  return 4.*dzeta2dxi*dzeta2dxi;
814 
815  case 3:
816  return 8.*dzeta0dxi*dzeta1dxi;
817 
818  case 4:
819  return 8.*dzeta1dxi*dzeta2dxi;
820 
821  case 5:
822  return 8.*dzeta0dxi*dzeta2dxi;
823 
824  default:
825  libmesh_error_msg("Invalid shape function index i = " << i);
826  }
827  }
828 
829  // d^2() / dxi deta
830  case 1:
831  {
832  switch (i)
833  {
834  case 0:
835  return 4.*dzeta0dxi*dzeta0deta;
836 
837  case 1:
838  return 4.*dzeta1dxi*dzeta1deta;
839 
840  case 2:
841  return 4.*dzeta2dxi*dzeta2deta;
842 
843  case 3:
844  return 4.*dzeta1deta*dzeta0dxi + 4.*dzeta0deta*dzeta1dxi;
845 
846  case 4:
847  return 4.*dzeta2deta*dzeta1dxi + 4.*dzeta1deta*dzeta2dxi;
848 
849  case 5:
850  return 4.*dzeta2deta*dzeta0dxi + 4.*dzeta0deta*dzeta2dxi;
851 
852  default:
853  libmesh_error_msg("Invalid shape function index i = " << i);
854  }
855  }
856 
857  // d^2() / deta^2
858  case 2:
859  {
860  switch (i)
861  {
862  case 0:
863  return 4.*dzeta0deta*dzeta0deta;
864 
865  case 1:
866  return 4.*dzeta1deta*dzeta1deta;
867 
868  case 2:
869  return 4.*dzeta2deta*dzeta2deta;
870 
871  case 3:
872  return 8.*dzeta0deta*dzeta1deta;
873 
874  case 4:
875  return 8.*dzeta1deta*dzeta2deta;
876 
877  case 5:
878  return 8.*dzeta0deta*dzeta2deta;
879 
880  default:
881  libmesh_error_msg("Invalid shape function index i = " << i);
882  }
883  }
884 
885  default:
886  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
887  } // end switch (j)
888  } // end case TRI6
889 
890  default:
891  libmesh_error_msg("ERROR: Unsupported 2D element type: " << type);
892  }
893  } // end case SECOND
894 
895 
896 
897  // unsupported order
898  default:
899  libmesh_error_msg("ERROR: Unsupported 2D FE order: " << order);
900 
901  } // end switch (order)
902 
903  libmesh_error_msg("We'll never get here!");
904 #endif // LIBMESH_DIM > 1
905  return 0.;
906 }
907 
908 
909 
910 template <>
912  const Order order,
913  const unsigned int i,
914  const unsigned int j,
915  const Point & p)
916 {
917  libmesh_assert(elem);
918 
919  // call the orientation-independent shape functions
920  return FE<2,LAGRANGE>::shape_second_deriv(elem->type(), static_cast<Order>(order + elem->p_level()), i, j, p);
921 }
922 
923 } // namespace libMesh
static OutputShape shape(const ElemType t, const Order o, const unsigned int i, const Point &p)
virtual ElemType type() const =0
unsigned int p_level() const
Definition: elem.h:2422
ElemType
Defines an enum for geometric element types.
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
static OutputShape shape_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)
The libMesh namespace provides an interface to certain functionality in the library.
libmesh_assert(j)
A specific instantiation of the FEBase class.
Definition: fe.h:89
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Order
defines an enum for polynomial orders.
Definition: enum_order.h:32
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
static OutputShape shape_second_deriv(const ElemType t, const Order o, const unsigned int i, const unsigned int j, const Point &p)