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

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

Hosted By:
SourceForge.net Logo