libMesh
fe_l2_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,L2_LAGRANGE>::shape(EDGE2, FIRST, i0[i], xi)*
63  FE<1,L2_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,L2_LAGRANGE>::shape(EDGE3, SECOND, i0[i], xi)*
155  FE<1,L2_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  return 0.;
205 
206 #endif
207 }
208 
209 
210 
211 template <>
213  const Order order,
214  const unsigned int i,
215  const Point & p)
216 {
217  libmesh_assert(elem);
218 
219  // call the orientation-independent shape functions
220  return FE<2,L2_LAGRANGE>::shape(elem->type(), static_cast<Order>(order + elem->p_level()), i, p);
221 }
222 
223 
224 
225 template <>
227  const Order order,
228  const unsigned int i,
229  const unsigned int j,
230  const Point & p)
231 {
232 #if LIBMESH_DIM > 1
233 
234 
235  libmesh_assert_less (j, 2);
236 
237  switch (order)
238  {
239  // linear Lagrange shape functions
240  case FIRST:
241  {
242  switch (type)
243  {
244  case QUAD4:
245  case QUADSHELL4:
246  case QUAD8:
247  case QUADSHELL8:
248  case QUAD9:
249  {
250  // Compute quad shape functions as a tensor-product
251  const Real xi = p(0);
252  const Real eta = p(1);
253 
254  libmesh_assert_less (i, 4);
255 
256  // 0 1 2 3
257  static const unsigned int i0[] = {0, 1, 1, 0};
258  static const unsigned int i1[] = {0, 0, 1, 1};
259 
260  switch (j)
261  {
262  // d()/dxi
263  case 0:
264  return (FE<1,L2_LAGRANGE>::shape_deriv(EDGE2, FIRST, i0[i], 0, xi)*
265  FE<1,L2_LAGRANGE>::shape (EDGE2, FIRST, i1[i], eta));
266 
267  // d()/deta
268  case 1:
269  return (FE<1,L2_LAGRANGE>::shape (EDGE2, FIRST, i0[i], xi)*
270  FE<1,L2_LAGRANGE>::shape_deriv(EDGE2, FIRST, i1[i], 0, eta));
271 
272  default:
273  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
274  }
275  }
276 
277  case TRI3:
278  case TRISHELL3:
279  case TRI6:
280  {
281  libmesh_assert_less (i, 3);
282 
283  const Real dzeta0dxi = -1.;
284  const Real dzeta1dxi = 1.;
285  const Real dzeta2dxi = 0.;
286 
287  const Real dzeta0deta = -1.;
288  const Real dzeta1deta = 0.;
289  const Real dzeta2deta = 1.;
290 
291  switch (j)
292  {
293  // d()/dxi
294  case 0:
295  {
296  switch(i)
297  {
298  case 0:
299  return dzeta0dxi;
300 
301  case 1:
302  return dzeta1dxi;
303 
304  case 2:
305  return dzeta2dxi;
306 
307  default:
308  libmesh_error_msg("Invalid shape function index i = " << i);
309  }
310  }
311  // d()/deta
312  case 1:
313  {
314  switch(i)
315  {
316  case 0:
317  return dzeta0deta;
318 
319  case 1:
320  return dzeta1deta;
321 
322  case 2:
323  return dzeta2deta;
324 
325  default:
326  libmesh_error_msg("Invalid shape function index i = " << i);
327  }
328  }
329  default:
330  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
331  }
332  }
333 
334  default:
335  libmesh_error_msg("ERROR: Unsupported 2D element type: " << type);
336  }
337  }
338 
339 
340  // quadratic Lagrange shape functions
341  case SECOND:
342  {
343  switch (type)
344  {
345  case QUAD8:
346  case QUADSHELL8:
347  {
348  const Real xi = p(0);
349  const Real eta = p(1);
350 
351  libmesh_assert_less (i, 8);
352 
353  switch (j)
354  {
355  // d/dxi
356  case 0:
357  switch (i)
358  {
359  case 0:
360  return .25*(1. - eta)*((1. - xi)*(-1.) +
361  (-1.)*(-1. - xi - eta));
362 
363  case 1:
364  return .25*(1. - eta)*((1. + xi)*(1.) +
365  (1.)*(-1. + xi - eta));
366 
367  case 2:
368  return .25*(1. + eta)*((1. + xi)*(1.) +
369  (1.)*(-1. + xi + eta));
370 
371  case 3:
372  return .25*(1. + eta)*((1. - xi)*(-1.) +
373  (-1.)*(-1. - xi + eta));
374 
375  case 4:
376  return .5*(-2.*xi)*(1. - eta);
377 
378  case 5:
379  return .5*(1.)*(1. - eta*eta);
380 
381  case 6:
382  return .5*(-2.*xi)*(1. + eta);
383 
384  case 7:
385  return .5*(-1.)*(1. - eta*eta);
386 
387  default:
388  libmesh_error_msg("Invalid shape function index i = " << i);
389  }
390 
391  // d/deta
392  case 1:
393  switch (i)
394  {
395  case 0:
396  return .25*(1. - xi)*((1. - eta)*(-1.) +
397  (-1.)*(-1. - xi - eta));
398 
399  case 1:
400  return .25*(1. + xi)*((1. - eta)*(-1.) +
401  (-1.)*(-1. + xi - eta));
402 
403  case 2:
404  return .25*(1. + xi)*((1. + eta)*(1.) +
405  (1.)*(-1. + xi + eta));
406 
407  case 3:
408  return .25*(1. - xi)*((1. + eta)*(1.) +
409  (1.)*(-1. - xi + eta));
410 
411  case 4:
412  return .5*(1. - xi*xi)*(-1.);
413 
414  case 5:
415  return .5*(1. + xi)*(-2.*eta);
416 
417  case 6:
418  return .5*(1. - xi*xi)*(1.);
419 
420  case 7:
421  return .5*(1. - xi)*(-2.*eta);
422 
423  default:
424  libmesh_error_msg("Invalid shape function index i = " << i);
425  }
426 
427  default:
428  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
429  }
430  }
431 
432  case QUAD9:
433  {
434  // Compute quad shape functions as a tensor-product
435  const Real xi = p(0);
436  const Real eta = p(1);
437 
438  libmesh_assert_less (i, 9);
439 
440  // 0 1 2 3 4 5 6 7 8
441  static const unsigned int i0[] = {0, 1, 1, 0, 2, 1, 2, 0, 2};
442  static const unsigned int i1[] = {0, 0, 1, 1, 0, 2, 1, 2, 2};
443 
444  switch (j)
445  {
446  // d()/dxi
447  case 0:
448  return (FE<1,L2_LAGRANGE>::shape_deriv(EDGE3, SECOND, i0[i], 0, xi)*
449  FE<1,L2_LAGRANGE>::shape (EDGE3, SECOND, i1[i], eta));
450 
451  // d()/deta
452  case 1:
453  return (FE<1,L2_LAGRANGE>::shape (EDGE3, SECOND, i0[i], xi)*
454  FE<1,L2_LAGRANGE>::shape_deriv(EDGE3, SECOND, i1[i], 0, eta));
455 
456  default:
457  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
458  }
459  }
460 
461  case TRI6:
462  {
463  libmesh_assert_less (i, 6);
464 
465  const Real zeta1 = p(0);
466  const Real zeta2 = p(1);
467  const Real zeta0 = 1. - zeta1 - zeta2;
468 
469  const Real dzeta0dxi = -1.;
470  const Real dzeta1dxi = 1.;
471  const Real dzeta2dxi = 0.;
472 
473  const Real dzeta0deta = -1.;
474  const Real dzeta1deta = 0.;
475  const Real dzeta2deta = 1.;
476 
477  switch(j)
478  {
479  case 0:
480  {
481  switch(i)
482  {
483  case 0:
484  return (4.*zeta0-1.)*dzeta0dxi;
485 
486  case 1:
487  return (4.*zeta1-1.)*dzeta1dxi;
488 
489  case 2:
490  return (4.*zeta2-1.)*dzeta2dxi;
491 
492  case 3:
493  return 4.*zeta1*dzeta0dxi + 4.*zeta0*dzeta1dxi;
494 
495  case 4:
496  return 4.*zeta2*dzeta1dxi + 4.*zeta1*dzeta2dxi;
497 
498  case 5:
499  return 4.*zeta2*dzeta0dxi + 4*zeta0*dzeta2dxi;
500 
501  default:
502  libmesh_error_msg("Invalid shape function index i = " << i);
503  }
504  }
505 
506  case 1:
507  {
508  switch(i)
509  {
510  case 0:
511  return (4.*zeta0-1.)*dzeta0deta;
512 
513  case 1:
514  return (4.*zeta1-1.)*dzeta1deta;
515 
516  case 2:
517  return (4.*zeta2-1.)*dzeta2deta;
518 
519  case 3:
520  return 4.*zeta1*dzeta0deta + 4.*zeta0*dzeta1deta;
521 
522  case 4:
523  return 4.*zeta2*dzeta1deta + 4.*zeta1*dzeta2deta;
524 
525  case 5:
526  return 4.*zeta2*dzeta0deta + 4*zeta0*dzeta2deta;
527 
528  default:
529  libmesh_error_msg("Invalid shape function index i = " << i);
530  }
531  }
532  default:
533  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
534  }
535  }
536 
537  default:
538  libmesh_error_msg("ERROR: Unsupported 2D element type: " << type);
539  }
540  }
541 
542 
543 
544  // unsupported order
545  default:
546  libmesh_error_msg("ERROR: Unsupported 2D FE order: " << order);
547  }
548 
549 
550  libmesh_error_msg("We'll never get here!");
551  return 0.;
552 
553 #endif
554 }
555 
556 
557 
558 template <>
560  const Order order,
561  const unsigned int i,
562  const unsigned int j,
563  const Point & p)
564 {
565  libmesh_assert(elem);
566 
567 
568  // call the orientation-independent shape functions
569  return FE<2,L2_LAGRANGE>::shape_deriv(elem->type(), static_cast<Order>(order + elem->p_level()), i, j, p);
570 }
571 
572 
573 
574 
575 template <>
577  const Order order,
578  const unsigned int i,
579  const unsigned int j,
580  const Point & p)
581 {
582 #if LIBMESH_DIM > 1
583 
584  // j = 0 ==> d^2 phi / dxi^2
585  // j = 1 ==> d^2 phi / dxi deta
586  // j = 2 ==> d^2 phi / deta^2
587  libmesh_assert_less (j, 3);
588 
589  switch (order)
590  {
591  // linear Lagrange shape functions
592  case FIRST:
593  {
594  switch (type)
595  {
596  case QUAD4:
597  case QUADSHELL4:
598  case QUAD8:
599  case QUADSHELL8:
600  case QUAD9:
601  {
602  // Compute quad shape functions as a tensor-product
603  const Real xi = p(0);
604  const Real eta = p(1);
605 
606  libmesh_assert_less (i, 4);
607 
608  // 0 1 2 3
609  static const unsigned int i0[] = {0, 1, 1, 0};
610  static const unsigned int i1[] = {0, 0, 1, 1};
611 
612  switch (j)
613  {
614  // d^2() / dxi^2
615  case 0:
616  return 0.;
617 
618  // d^2() / dxi deta
619  case 1:
620  return (FE<1,L2_LAGRANGE>::shape_deriv(EDGE2, FIRST, i0[i], 0, xi)*
621  FE<1,L2_LAGRANGE>::shape_deriv(EDGE2, FIRST, i1[i], 0, eta));
622 
623  // d^2() / deta^2
624  case 2:
625  return 0.;
626 
627  default:
628  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
629  }
630  }
631 
632  case TRI3:
633  case TRISHELL3:
634  case TRI6:
635  {
636  // All second derivatives for linear triangles are zero.
637  return 0.;
638  }
639 
640  default:
641  libmesh_error_msg("ERROR: Unsupported 2D element type: " << type);
642 
643  } // end switch (type)
644  } // end case FIRST
645 
646 
647  // quadratic Lagrange shape functions
648  case SECOND:
649  {
650  switch (type)
651  {
652  case QUAD8:
653  {
654  const Real xi = p(0);
655  const Real eta = p(1);
656 
657  libmesh_assert_less (j, 3);
658 
659  switch (j)
660  {
661  // d^2() / dxi^2
662  case 0:
663  {
664  switch (i)
665  {
666  case 0:
667  case 1:
668  return 0.5*(1.-eta);
669 
670  case 2:
671  case 3:
672  return 0.5*(1.+eta);
673 
674  case 4:
675  return eta - 1.;
676 
677  case 5:
678  case 7:
679  return 0.0;
680 
681  case 6:
682  return -1. - eta;
683 
684  default:
685  libmesh_error_msg("Invalid shape function index i = " << i);
686  }
687  }
688 
689  // d^2() / dxi deta
690  case 1:
691  {
692  switch (i)
693  {
694  case 0:
695  return 0.25*( 1. - 2.*xi - 2.*eta);
696 
697  case 1:
698  return 0.25*(-1. - 2.*xi + 2.*eta);
699 
700  case 2:
701  return 0.25*( 1. + 2.*xi + 2.*eta);
702 
703  case 3:
704  return 0.25*(-1. + 2.*xi - 2.*eta);
705 
706  case 4:
707  return xi;
708 
709  case 5:
710  return -eta;
711 
712  case 6:
713  return -xi;
714 
715  case 7:
716  return eta;
717 
718  default:
719  libmesh_error_msg("Invalid shape function index i = " << i);
720  }
721  }
722 
723  // d^2() / deta^2
724  case 2:
725  {
726  switch (i)
727  {
728  case 0:
729  case 3:
730  return 0.5*(1.-xi);
731 
732  case 1:
733  case 2:
734  return 0.5*(1.+xi);
735 
736  case 4:
737  case 6:
738  return 0.0;
739 
740  case 5:
741  return -1.0 - xi;
742 
743  case 7:
744  return xi - 1.0;
745 
746  default:
747  libmesh_error_msg("Invalid shape function index i = " << i);
748  }
749  }
750 
751 
752  default:
753  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
754  } // end switch (j)
755  } // end case QUAD8
756 
757  case QUAD9:
758  {
759  // Compute QUAD9 second derivatives as tensor product
760  const Real xi = p(0);
761  const Real eta = p(1);
762 
763  libmesh_assert_less (i, 9);
764 
765  // 0 1 2 3 4 5 6 7 8
766  static const unsigned int i0[] = {0, 1, 1, 0, 2, 1, 2, 0, 2};
767  static const unsigned int i1[] = {0, 0, 1, 1, 0, 2, 1, 2, 2};
768 
769  switch (j)
770  {
771  // d^2() / dxi^2
772  case 0:
773  return (FE<1,L2_LAGRANGE>::shape_second_deriv(EDGE3, SECOND, i0[i], 0, xi)*
774  FE<1,L2_LAGRANGE>::shape (EDGE3, SECOND, i1[i], eta));
775 
776  // d^2() / dxi deta
777  case 1:
778  return (FE<1,L2_LAGRANGE>::shape_deriv(EDGE3, SECOND, i0[i], 0, xi)*
779  FE<1,L2_LAGRANGE>::shape_deriv(EDGE3, SECOND, i1[i], 0, eta));
780 
781  // d^2() / deta^2
782  case 2:
783  return (FE<1,L2_LAGRANGE>::shape (EDGE3, SECOND, i0[i], xi)*
785 
786  default:
787  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
788  } // end switch (j)
789  } // end case QUAD9
790 
791  case TRI6:
792  {
793  const Real dzeta0dxi = -1.;
794  const Real dzeta1dxi = 1.;
795  const Real dzeta2dxi = 0.;
796 
797  const Real dzeta0deta = -1.;
798  const Real dzeta1deta = 0.;
799  const Real dzeta2deta = 1.;
800 
801  libmesh_assert_less (j, 3);
802 
803  switch (j)
804  {
805  // d^2() / dxi^2
806  case 0:
807  {
808  switch (i)
809  {
810  case 0:
811  return 4.*dzeta0dxi*dzeta0dxi;
812 
813  case 1:
814  return 4.*dzeta1dxi*dzeta1dxi;
815 
816  case 2:
817  return 4.*dzeta2dxi*dzeta2dxi;
818 
819  case 3:
820  return 8.*dzeta0dxi*dzeta1dxi;
821 
822  case 4:
823  return 8.*dzeta1dxi*dzeta2dxi;
824 
825  case 5:
826  return 8.*dzeta0dxi*dzeta2dxi;
827 
828  default:
829  libmesh_error_msg("Invalid shape function index i = " << i);
830  }
831  }
832 
833  // d^2() / dxi deta
834  case 1:
835  {
836  switch (i)
837  {
838  case 0:
839  return 4.*dzeta0dxi*dzeta0deta;
840 
841  case 1:
842  return 4.*dzeta1dxi*dzeta1deta;
843 
844  case 2:
845  return 4.*dzeta2dxi*dzeta2deta;
846 
847  case 3:
848  return 4.*dzeta1deta*dzeta0dxi + 4.*dzeta0deta*dzeta1dxi;
849 
850  case 4:
851  return 4.*dzeta2deta*dzeta1dxi + 4.*dzeta1deta*dzeta2dxi;
852 
853  case 5:
854  return 4.*dzeta2deta*dzeta0dxi + 4.*dzeta0deta*dzeta2dxi;
855 
856  default:
857  libmesh_error_msg("Invalid shape function index i = " << i);
858  }
859  }
860 
861  // d^2() / deta^2
862  case 2:
863  {
864  switch (i)
865  {
866  case 0:
867  return 4.*dzeta0deta*dzeta0deta;
868 
869  case 1:
870  return 4.*dzeta1deta*dzeta1deta;
871 
872  case 2:
873  return 4.*dzeta2deta*dzeta2deta;
874 
875  case 3:
876  return 8.*dzeta0deta*dzeta1deta;
877 
878  case 4:
879  return 8.*dzeta1deta*dzeta2deta;
880 
881  case 5:
882  return 8.*dzeta0deta*dzeta2deta;
883 
884  default:
885  libmesh_error_msg("Invalid shape function index i = " << i);
886  }
887  }
888 
889  default:
890  libmesh_error_msg("ERROR: Invalid derivative index j = " << j);
891  } // end switch (j)
892  } // end case TRI6
893 
894  default:
895  libmesh_error_msg("ERROR: Unsupported 2D element type: " << type);
896  }
897  } // end case SECOND
898 
899 
900 
901  // unsupported order
902  default:
903  libmesh_error_msg("ERROR: Unsupported 2D FE order: " << order);
904 
905  } // end switch (order)
906 
907 
908  libmesh_error_msg("We'll never get here!");
909  return 0.;
910 #endif // LIBMESH_DIM > 1
911 }
912 
913 
914 
915 template <>
917  const Order order,
918  const unsigned int i,
919  const unsigned int j,
920  const Point & p)
921 {
922  libmesh_assert(elem);
923 
924  // call the orientation-independent shape functions
925  return FE<2,L2_LAGRANGE>::shape_second_deriv(elem->type(), static_cast<Order>(order + elem->p_level()), i, j, p);
926 }
927 
928 } // 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)