fe_l2_lagrange_shape_2D.C
Go to the documentation of this file.00001 // The libMesh Finite Element Library. 00002 // Copyright (C) 2002-2012 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner 00003 00004 // This library is free software; you can redistribute it and/or 00005 // modify it under the terms of the GNU Lesser General Public 00006 // License as published by the Free Software Foundation; either 00007 // version 2.1 of the License, or (at your option) any later version. 00008 00009 // This library is distributed in the hope that it will be useful, 00010 // but WITHOUT ANY WARRANTY; without even the implied warranty of 00011 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00012 // Lesser General Public License for more details. 00013 00014 // You should have received a copy of the GNU Lesser General Public 00015 // License along with this library; if not, write to the Free Software 00016 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA 00017 00018 00019 // C++ inlcludes 00020 00021 // Local includes 00022 #include "libmesh/fe.h" 00023 #include "libmesh/elem.h" 00024 00025 namespace libMesh 00026 { 00027 00028 00029 00030 00031 template <> 00032 Real FE<2,L2_LAGRANGE>::shape(const ElemType type, 00033 const Order order, 00034 const unsigned int i, 00035 const Point& p) 00036 { 00037 #if LIBMESH_DIM > 1 00038 00039 switch (order) 00040 { 00041 // linear Lagrange shape functions 00042 case FIRST: 00043 { 00044 switch (type) 00045 { 00046 case QUAD4: 00047 case QUAD8: 00048 case QUAD9: 00049 { 00050 // Compute quad shape functions as a tensor-product 00051 const Real xi = p(0); 00052 const Real eta = p(1); 00053 00054 libmesh_assert_less (i, 4); 00055 00056 // 0 1 2 3 00057 static const unsigned int i0[] = {0, 1, 1, 0}; 00058 static const unsigned int i1[] = {0, 0, 1, 1}; 00059 00060 return (FE<1,L2_LAGRANGE>::shape(EDGE2, FIRST, i0[i], xi)* 00061 FE<1,L2_LAGRANGE>::shape(EDGE2, FIRST, i1[i], eta)); 00062 } 00063 00064 case TRI3: 00065 case TRI6: 00066 { 00067 const Real zeta1 = p(0); 00068 const Real zeta2 = p(1); 00069 const Real zeta0 = 1. - zeta1 - zeta2; 00070 00071 libmesh_assert_less (i, 3); 00072 00073 switch(i) 00074 { 00075 case 0: 00076 return zeta0; 00077 00078 case 1: 00079 return zeta1; 00080 00081 case 2: 00082 return zeta2; 00083 00084 default: 00085 libmesh_error(); 00086 00087 } 00088 } 00089 00090 default: 00091 { 00092 libMesh::err << "ERROR: Unsupported 2D element type!: " << type 00093 << std::endl; 00094 libmesh_error(); 00095 } 00096 } 00097 } 00098 00099 00100 // quadratic Lagrange shape functions 00101 case SECOND: 00102 { 00103 switch (type) 00104 { 00105 case QUAD8: 00106 { 00107 const Real xi = p(0); 00108 const Real eta = p(1); 00109 00110 libmesh_assert_less (i, 8); 00111 00112 switch (i) 00113 { 00114 case 0: 00115 return .25*(1. - xi)*(1. - eta)*(-1. - xi - eta); 00116 00117 case 1: 00118 return .25*(1. + xi)*(1. - eta)*(-1. + xi - eta); 00119 00120 case 2: 00121 return .25*(1. + xi)*(1. + eta)*(-1. + xi + eta); 00122 00123 case 3: 00124 return .25*(1. - xi)*(1. + eta)*(-1. - xi + eta); 00125 00126 case 4: 00127 return .5*(1. - xi*xi)*(1. - eta); 00128 00129 case 5: 00130 return .5*(1. + xi)*(1. - eta*eta); 00131 00132 case 6: 00133 return .5*(1. - xi*xi)*(1. + eta); 00134 00135 case 7: 00136 return .5*(1. - xi)*(1. - eta*eta); 00137 00138 default: 00139 libmesh_error(); 00140 } 00141 } 00142 00143 case QUAD9: 00144 { 00145 // Compute quad shape functions as a tensor-product 00146 const Real xi = p(0); 00147 const Real eta = p(1); 00148 00149 libmesh_assert_less (i, 9); 00150 00151 // 0 1 2 3 4 5 6 7 8 00152 static const unsigned int i0[] = {0, 1, 1, 0, 2, 1, 2, 0, 2}; 00153 static const unsigned int i1[] = {0, 0, 1, 1, 0, 2, 1, 2, 2}; 00154 00155 return (FE<1,L2_LAGRANGE>::shape(EDGE3, SECOND, i0[i], xi)* 00156 FE<1,L2_LAGRANGE>::shape(EDGE3, SECOND, i1[i], eta)); 00157 } 00158 00159 case TRI6: 00160 { 00161 const Real zeta1 = p(0); 00162 const Real zeta2 = p(1); 00163 const Real zeta0 = 1. - zeta1 - zeta2; 00164 00165 libmesh_assert_less (i, 6); 00166 00167 switch(i) 00168 { 00169 case 0: 00170 return 2.*zeta0*(zeta0-0.5); 00171 00172 case 1: 00173 return 2.*zeta1*(zeta1-0.5); 00174 00175 case 2: 00176 return 2.*zeta2*(zeta2-0.5); 00177 00178 case 3: 00179 return 4.*zeta0*zeta1; 00180 00181 case 4: 00182 return 4.*zeta1*zeta2; 00183 00184 case 5: 00185 return 4.*zeta2*zeta0; 00186 00187 default: 00188 libmesh_error(); 00189 00190 } 00191 } 00192 00193 default: 00194 { 00195 libMesh::err << "ERROR: Unsupported 2D element type!: " << type 00196 << std::endl; 00197 libmesh_error(); 00198 } 00199 } 00200 } 00201 00202 00203 00204 // unsupported order 00205 default: 00206 { 00207 libMesh::err << "ERROR: Unsupported 2D FE order!: " << order 00208 << std::endl; 00209 libmesh_error(); 00210 } 00211 } 00212 00213 libmesh_error(); 00214 return 0.; 00215 00216 #endif 00217 } 00218 00219 00220 00221 template <> 00222 Real FE<2,L2_LAGRANGE>::shape(const Elem* elem, 00223 const Order order, 00224 const unsigned int i, 00225 const Point& p) 00226 { 00227 libmesh_assert(elem); 00228 00229 // call the orientation-independent shape functions 00230 return FE<2,L2_LAGRANGE>::shape(elem->type(), static_cast<Order>(order + elem->p_level()), i, p); 00231 } 00232 00233 00234 00235 template <> 00236 Real FE<2,L2_LAGRANGE>::shape_deriv(const ElemType type, 00237 const Order order, 00238 const unsigned int i, 00239 const unsigned int j, 00240 const Point& p) 00241 { 00242 #if LIBMESH_DIM > 1 00243 00244 00245 libmesh_assert_less (j, 2); 00246 00247 switch (order) 00248 { 00249 // linear Lagrange shape functions 00250 case FIRST: 00251 { 00252 switch (type) 00253 { 00254 case QUAD4: 00255 case QUAD8: 00256 case QUAD9: 00257 { 00258 // Compute quad shape functions as a tensor-product 00259 const Real xi = p(0); 00260 const Real eta = p(1); 00261 00262 libmesh_assert_less (i, 4); 00263 00264 // 0 1 2 3 00265 static const unsigned int i0[] = {0, 1, 1, 0}; 00266 static const unsigned int i1[] = {0, 0, 1, 1}; 00267 00268 switch (j) 00269 { 00270 // d()/dxi 00271 case 0: 00272 return (FE<1,L2_LAGRANGE>::shape_deriv(EDGE2, FIRST, i0[i], 0, xi)* 00273 FE<1,L2_LAGRANGE>::shape (EDGE2, FIRST, i1[i], eta)); 00274 00275 // d()/deta 00276 case 1: 00277 return (FE<1,L2_LAGRANGE>::shape (EDGE2, FIRST, i0[i], xi)* 00278 FE<1,L2_LAGRANGE>::shape_deriv(EDGE2, FIRST, i1[i], 0, eta)); 00279 00280 default: 00281 libmesh_error(); 00282 } 00283 } 00284 00285 case TRI3: 00286 case TRI6: 00287 { 00288 libmesh_assert_less (i, 3); 00289 00290 const Real dzeta0dxi = -1.; 00291 const Real dzeta1dxi = 1.; 00292 const Real dzeta2dxi = 0.; 00293 00294 const Real dzeta0deta = -1.; 00295 const Real dzeta1deta = 0.; 00296 const Real dzeta2deta = 1.; 00297 00298 switch (j) 00299 { 00300 // d()/dxi 00301 case 0: 00302 { 00303 switch(i) 00304 { 00305 case 0: 00306 return dzeta0dxi; 00307 00308 case 1: 00309 return dzeta1dxi; 00310 00311 case 2: 00312 return dzeta2dxi; 00313 00314 default: 00315 libmesh_error(); 00316 } 00317 } 00318 // d()/deta 00319 case 1: 00320 { 00321 switch(i) 00322 { 00323 case 0: 00324 return dzeta0deta; 00325 00326 case 1: 00327 return dzeta1deta; 00328 00329 case 2: 00330 return dzeta2deta; 00331 00332 default: 00333 libmesh_error(); 00334 00335 } 00336 } 00337 default: 00338 libmesh_error(); 00339 } 00340 } 00341 00342 default: 00343 { 00344 libMesh::err << "ERROR: Unsupported 2D element type!: " << type 00345 << std::endl; 00346 libmesh_error(); 00347 } 00348 } 00349 } 00350 00351 00352 // quadratic Lagrange shape functions 00353 case SECOND: 00354 { 00355 switch (type) 00356 { 00357 case QUAD8: 00358 { 00359 const Real xi = p(0); 00360 const Real eta = p(1); 00361 00362 libmesh_assert_less (i, 8); 00363 00364 switch (j) 00365 { 00366 // d/dxi 00367 case 0: 00368 switch (i) 00369 { 00370 case 0: 00371 return .25*(1. - eta)*((1. - xi)*(-1.) + 00372 (-1.)*(-1. - xi - eta)); 00373 00374 case 1: 00375 return .25*(1. - eta)*((1. + xi)*(1.) + 00376 (1.)*(-1. + xi - eta)); 00377 00378 case 2: 00379 return .25*(1. + eta)*((1. + xi)*(1.) + 00380 (1.)*(-1. + xi + eta)); 00381 00382 case 3: 00383 return .25*(1. + eta)*((1. - xi)*(-1.) + 00384 (-1.)*(-1. - xi + eta)); 00385 00386 case 4: 00387 return .5*(-2.*xi)*(1. - eta); 00388 00389 case 5: 00390 return .5*(1.)*(1. - eta*eta); 00391 00392 case 6: 00393 return .5*(-2.*xi)*(1. + eta); 00394 00395 case 7: 00396 return .5*(-1.)*(1. - eta*eta); 00397 00398 default: 00399 libmesh_error(); 00400 } 00401 00402 // d/deta 00403 case 1: 00404 switch (i) 00405 { 00406 case 0: 00407 return .25*(1. - xi)*((1. - eta)*(-1.) + 00408 (-1.)*(-1. - xi - eta)); 00409 00410 case 1: 00411 return .25*(1. + xi)*((1. - eta)*(-1.) + 00412 (-1.)*(-1. + xi - eta)); 00413 00414 case 2: 00415 return .25*(1. + xi)*((1. + eta)*(1.) + 00416 (1.)*(-1. + xi + eta)); 00417 00418 case 3: 00419 return .25*(1. - xi)*((1. + eta)*(1.) + 00420 (1.)*(-1. - xi + eta)); 00421 00422 case 4: 00423 return .5*(1. - xi*xi)*(-1.); 00424 00425 case 5: 00426 return .5*(1. + xi)*(-2.*eta); 00427 00428 case 6: 00429 return .5*(1. - xi*xi)*(1.); 00430 00431 case 7: 00432 return .5*(1. - xi)*(-2.*eta); 00433 00434 default: 00435 libmesh_error(); 00436 } 00437 00438 default: 00439 libmesh_error(); 00440 } 00441 } 00442 00443 case QUAD9: 00444 { 00445 // Compute quad shape functions as a tensor-product 00446 const Real xi = p(0); 00447 const Real eta = p(1); 00448 00449 libmesh_assert_less (i, 9); 00450 00451 // 0 1 2 3 4 5 6 7 8 00452 static const unsigned int i0[] = {0, 1, 1, 0, 2, 1, 2, 0, 2}; 00453 static const unsigned int i1[] = {0, 0, 1, 1, 0, 2, 1, 2, 2}; 00454 00455 switch (j) 00456 { 00457 // d()/dxi 00458 case 0: 00459 return (FE<1,L2_LAGRANGE>::shape_deriv(EDGE3, SECOND, i0[i], 0, xi)* 00460 FE<1,L2_LAGRANGE>::shape (EDGE3, SECOND, i1[i], eta)); 00461 00462 // d()/deta 00463 case 1: 00464 return (FE<1,L2_LAGRANGE>::shape (EDGE3, SECOND, i0[i], xi)* 00465 FE<1,L2_LAGRANGE>::shape_deriv(EDGE3, SECOND, i1[i], 0, eta)); 00466 00467 default: 00468 libmesh_error(); 00469 } 00470 } 00471 00472 case TRI6: 00473 { 00474 libmesh_assert_less (i, 6); 00475 00476 const Real zeta1 = p(0); 00477 const Real zeta2 = p(1); 00478 const Real zeta0 = 1. - zeta1 - zeta2; 00479 00480 const Real dzeta0dxi = -1.; 00481 const Real dzeta1dxi = 1.; 00482 const Real dzeta2dxi = 0.; 00483 00484 const Real dzeta0deta = -1.; 00485 const Real dzeta1deta = 0.; 00486 const Real dzeta2deta = 1.; 00487 00488 switch(j) 00489 { 00490 case 0: 00491 { 00492 switch(i) 00493 { 00494 case 0: 00495 return (4.*zeta0-1.)*dzeta0dxi; 00496 00497 case 1: 00498 return (4.*zeta1-1.)*dzeta1dxi; 00499 00500 case 2: 00501 return (4.*zeta2-1.)*dzeta2dxi; 00502 00503 case 3: 00504 return 4.*zeta1*dzeta0dxi + 4.*zeta0*dzeta1dxi; 00505 00506 case 4: 00507 return 4.*zeta2*dzeta1dxi + 4.*zeta1*dzeta2dxi; 00508 00509 case 5: 00510 return 4.*zeta2*dzeta0dxi + 4*zeta0*dzeta2dxi; 00511 00512 default: 00513 libmesh_error(); 00514 } 00515 } 00516 00517 case 1: 00518 { 00519 switch(i) 00520 { 00521 case 0: 00522 return (4.*zeta0-1.)*dzeta0deta; 00523 00524 case 1: 00525 return (4.*zeta1-1.)*dzeta1deta; 00526 00527 case 2: 00528 return (4.*zeta2-1.)*dzeta2deta; 00529 00530 case 3: 00531 return 4.*zeta1*dzeta0deta + 4.*zeta0*dzeta1deta; 00532 00533 case 4: 00534 return 4.*zeta2*dzeta1deta + 4.*zeta1*dzeta2deta; 00535 00536 case 5: 00537 return 4.*zeta2*dzeta0deta + 4*zeta0*dzeta2deta; 00538 00539 default: 00540 libmesh_error(); 00541 } 00542 } 00543 default: 00544 libmesh_error(); 00545 } 00546 } 00547 00548 default: 00549 { 00550 libMesh::err << "ERROR: Unsupported 2D element type!: " << type 00551 << std::endl; 00552 libmesh_error(); 00553 } 00554 } 00555 } 00556 00557 00558 00559 // unsupported order 00560 default: 00561 { 00562 libMesh::err << "ERROR: Unsupported 2D FE order!: " << order 00563 << std::endl; 00564 libmesh_error(); 00565 } 00566 } 00567 00568 00569 libmesh_error(); 00570 return 0.; 00571 00572 #endif 00573 } 00574 00575 00576 00577 template <> 00578 Real FE<2,L2_LAGRANGE>::shape_deriv(const Elem* elem, 00579 const Order order, 00580 const unsigned int i, 00581 const unsigned int j, 00582 const Point& p) 00583 { 00584 libmesh_assert(elem); 00585 00586 00587 // call the orientation-independent shape functions 00588 return FE<2,L2_LAGRANGE>::shape_deriv(elem->type(), static_cast<Order>(order + elem->p_level()), i, j, p); 00589 } 00590 00591 00592 00593 00594 template <> 00595 Real FE<2,L2_LAGRANGE>::shape_second_deriv(const ElemType type, 00596 const Order order, 00597 const unsigned int i, 00598 const unsigned int j, 00599 const Point& p) 00600 { 00601 #if LIBMESH_DIM > 1 00602 00603 // j = 0 ==> d^2 phi / dxi^2 00604 // j = 1 ==> d^2 phi / dxi deta 00605 // j = 2 ==> d^2 phi / deta^2 00606 libmesh_assert_less (j, 3); 00607 00608 switch (order) 00609 { 00610 // linear Lagrange shape functions 00611 case FIRST: 00612 { 00613 switch (type) 00614 { 00615 case QUAD4: 00616 case QUAD8: 00617 case QUAD9: 00618 { 00619 // Compute quad shape functions as a tensor-product 00620 const Real xi = p(0); 00621 const Real eta = p(1); 00622 00623 libmesh_assert_less (i, 4); 00624 00625 // 0 1 2 3 00626 static const unsigned int i0[] = {0, 1, 1, 0}; 00627 static const unsigned int i1[] = {0, 0, 1, 1}; 00628 00629 switch (j) 00630 { 00631 // d^2() / dxi^2 00632 case 0: 00633 return 0.; 00634 00635 // d^2() / dxi deta 00636 case 1: 00637 return (FE<1,L2_LAGRANGE>::shape_deriv(EDGE2, FIRST, i0[i], 0, xi)* 00638 FE<1,L2_LAGRANGE>::shape_deriv(EDGE2, FIRST, i1[1], 0, eta)); 00639 00640 // d^2() / deta^2 00641 case 2: 00642 return 0.; 00643 00644 default: 00645 { 00646 libMesh::err << "ERROR: Invalid derivative requested! " 00647 << std::endl; 00648 libmesh_error(); 00649 } 00650 } 00651 } 00652 00653 case TRI3: 00654 case TRI6: 00655 { 00656 // All second derivatives for linear triangles are zero. 00657 return 0.; 00658 } 00659 00660 default: 00661 { 00662 libMesh::err << "ERROR: Unsupported 2D element type!: " << type 00663 << std::endl; 00664 libmesh_error(); 00665 } 00666 00667 } // end switch (type) 00668 } // end case FIRST 00669 00670 00671 // quadratic Lagrange shape functions 00672 case SECOND: 00673 { 00674 switch (type) 00675 { 00676 case QUAD8: 00677 { 00678 const Real xi = p(0); 00679 const Real eta = p(1); 00680 00681 libmesh_assert_less (j, 3); 00682 00683 switch (j) 00684 { 00685 // d^2() / dxi^2 00686 case 0: 00687 { 00688 switch (i) 00689 { 00690 case 0: 00691 case 1: 00692 return 0.5*(1.-eta); 00693 00694 case 2: 00695 case 3: 00696 return 0.5*(1.+eta); 00697 00698 case 4: 00699 return eta - 1.; 00700 00701 case 5: 00702 case 7: 00703 return 0.0; 00704 00705 case 6: 00706 return -1. - eta; 00707 00708 default: 00709 { 00710 libMesh::err << "Invalid shape function index requested!" 00711 << std::endl; 00712 libmesh_error(); 00713 } 00714 } 00715 } 00716 00717 // d^2() / dxi deta 00718 case 1: 00719 { 00720 switch (i) 00721 { 00722 case 0: 00723 return 0.25*( 1. - 2.*xi - 2.*eta); 00724 00725 case 1: 00726 return 0.25*(-1. - 2.*xi + 2.*eta); 00727 00728 case 2: 00729 return 0.25*( 1. + 2.*xi + 2.*eta); 00730 00731 case 3: 00732 return 0.25*(-1. + 2.*xi - 2.*eta); 00733 00734 case 4: 00735 return xi; 00736 00737 case 5: 00738 return -eta; 00739 00740 case 6: 00741 return -xi; 00742 00743 case 7: 00744 return eta; 00745 00746 default: 00747 { 00748 libMesh::err << "Invalid shape function index requested!" 00749 << std::endl; 00750 libmesh_error(); 00751 } 00752 } 00753 } 00754 00755 // d^2() / deta^2 00756 case 2: 00757 { 00758 switch (i) 00759 { 00760 case 0: 00761 case 3: 00762 return 0.5*(1.-xi); 00763 00764 case 1: 00765 case 2: 00766 return 0.5*(1.+xi); 00767 00768 case 4: 00769 case 6: 00770 return 0.0; 00771 00772 case 5: 00773 return -1.0 - xi; 00774 00775 case 7: 00776 return xi - 1.0; 00777 00778 default: 00779 { 00780 libMesh::err << "Invalid shape function index requested!" 00781 << std::endl; 00782 libmesh_error(); 00783 } 00784 } 00785 } 00786 00787 00788 default: 00789 { 00790 libMesh::err << "ERROR: Invalid derivative requested! " 00791 << std::endl; 00792 libmesh_error(); 00793 } 00794 } // end switch (j) 00795 } // end case QUAD8 00796 00797 case QUAD9: 00798 { 00799 // Compute QUAD9 second derivatives as tensor product 00800 const Real xi = p(0); 00801 const Real eta = p(1); 00802 00803 libmesh_assert_less (i, 9); 00804 00805 // 0 1 2 3 4 5 6 7 8 00806 static const unsigned int i0[] = {0, 1, 1, 0, 2, 1, 2, 0, 2}; 00807 static const unsigned int i1[] = {0, 0, 1, 1, 0, 2, 1, 2, 2}; 00808 00809 switch (j) 00810 { 00811 // d^2() / dxi^2 00812 case 0: 00813 return (FE<1,L2_LAGRANGE>::shape_second_deriv(EDGE3, SECOND, i0[i], 0, xi)* 00814 FE<1,L2_LAGRANGE>::shape (EDGE3, SECOND, i1[i], eta)); 00815 00816 // d^2() / dxi deta 00817 case 1: 00818 return (FE<1,L2_LAGRANGE>::shape_deriv(EDGE3, SECOND, i0[i], 0, xi)* 00819 FE<1,L2_LAGRANGE>::shape_deriv(EDGE3, SECOND, i1[i], 0, eta)); 00820 00821 // d^2() / deta^2 00822 case 2: 00823 return (FE<1,L2_LAGRANGE>::shape (EDGE3, SECOND, i0[i], xi)* 00824 FE<1,L2_LAGRANGE>::shape_second_deriv(EDGE3, SECOND, i1[i], 0, eta)); 00825 00826 default: 00827 { 00828 libMesh::err << "ERROR: Invalid derivative requested! " 00829 << std::endl; 00830 libmesh_error(); 00831 } 00832 } // end switch (j) 00833 } // end case QUAD9 00834 00835 case TRI6: 00836 { 00837 const Real dzeta0dxi = -1.; 00838 const Real dzeta1dxi = 1.; 00839 const Real dzeta2dxi = 0.; 00840 00841 const Real dzeta0deta = -1.; 00842 const Real dzeta1deta = 0.; 00843 const Real dzeta2deta = 1.; 00844 00845 libmesh_assert_less (j, 3); 00846 00847 switch (j) 00848 { 00849 // d^2() / dxi^2 00850 case 0: 00851 { 00852 switch (i) 00853 { 00854 case 0: 00855 return 4.*dzeta0dxi*dzeta0dxi; 00856 00857 case 1: 00858 return 4.*dzeta1dxi*dzeta1dxi; 00859 00860 case 2: 00861 return 4.*dzeta2dxi*dzeta2dxi; 00862 00863 case 3: 00864 return 8.*dzeta0dxi*dzeta1dxi; 00865 00866 case 4: 00867 return 8.*dzeta1dxi*dzeta2dxi; 00868 00869 case 5: 00870 return 8.*dzeta0dxi*dzeta2dxi; 00871 00872 default: 00873 { 00874 libMesh::err << "Invalid shape function index requested!" 00875 << std::endl; 00876 libmesh_error(); 00877 } 00878 } 00879 } 00880 00881 // d^2() / dxi deta 00882 case 1: 00883 { 00884 switch (i) 00885 { 00886 case 0: 00887 return 4.*dzeta0dxi*dzeta0deta; 00888 00889 case 1: 00890 return 4.*dzeta1dxi*dzeta1deta; 00891 00892 case 2: 00893 return 4.*dzeta2dxi*dzeta2deta; 00894 00895 case 3: 00896 return 4.*dzeta1deta*dzeta0dxi + 4.*dzeta0deta*dzeta1dxi; 00897 00898 case 4: 00899 return 4.*dzeta2deta*dzeta1dxi + 4.*dzeta1deta*dzeta2dxi; 00900 00901 case 5: 00902 return 4.*dzeta2deta*dzeta0dxi + 4.*dzeta0deta*dzeta2dxi; 00903 00904 default: 00905 { 00906 libMesh::err << "Invalid shape function index requested!" 00907 << std::endl; 00908 libmesh_error(); 00909 } 00910 } 00911 } 00912 00913 // d^2() / deta^2 00914 case 2: 00915 { 00916 switch (i) 00917 { 00918 case 0: 00919 return 4.*dzeta0deta*dzeta0deta; 00920 00921 case 1: 00922 return 4.*dzeta1deta*dzeta1deta; 00923 00924 case 2: 00925 return 4.*dzeta2deta*dzeta2deta; 00926 00927 case 3: 00928 return 8.*dzeta0deta*dzeta1deta; 00929 00930 case 4: 00931 return 8.*dzeta1deta*dzeta2deta; 00932 00933 case 5: 00934 return 8.*dzeta0deta*dzeta2deta; 00935 00936 default: 00937 { 00938 libMesh::err << "Invalid shape function index requested!" 00939 << std::endl; 00940 libmesh_error(); 00941 } 00942 } 00943 } 00944 00945 default: 00946 { 00947 libMesh::err << "ERROR: Invalid derivative requested! " 00948 << std::endl; 00949 libmesh_error(); 00950 } 00951 } // end switch (j) 00952 } // end case TRI6 00953 00954 default: 00955 { 00956 libMesh::err << "ERROR: Unsupported 2D element type!: " << type 00957 << std::endl; 00958 libmesh_error(); 00959 } 00960 } 00961 } // end case SECOND 00962 00963 00964 00965 // unsupported order 00966 default: 00967 { 00968 libMesh::err << "ERROR: Unsupported 2D FE order!: " << order 00969 << std::endl; 00970 libmesh_error(); 00971 } 00972 00973 } // end switch (order) 00974 00975 00976 libmesh_error(); 00977 return 0.; 00978 #endif // LIBMESH_DIM > 1 00979 } 00980 00981 00982 00983 template <> 00984 Real FE<2,L2_LAGRANGE>::shape_second_deriv(const Elem* elem, 00985 const Order order, 00986 const unsigned int i, 00987 const unsigned int j, 00988 const Point& p) 00989 { 00990 libmesh_assert(elem); 00991 00992 // call the orientation-independent shape functions 00993 return FE<2,L2_LAGRANGE>::shape_second_deriv(elem->type(), static_cast<Order>(order + elem->p_level()), i, j, p); 00994 } 00995 00996 } // namespace libMesh
Site Created By: libMesh Developers
Last modified: February 05 2013 19:54:46 UTC
Hosted By: