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

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:05 UTC

Hosted By:
SourceForge.net Logo