$extrastylesheet
00001 // The libMesh Finite Element Library. 00002 // Copyright (C) 2002-2014 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<1,LAGRANGE>::shape(const ElemType, 00033 const Order order, 00034 const unsigned int i, 00035 const Point& p) 00036 { 00037 const Real xi = p(0); 00038 00039 00040 switch (order) 00041 { 00042 // Lagrange linears 00043 case FIRST: 00044 { 00045 libmesh_assert_less (i, 2); 00046 00047 switch (i) 00048 { 00049 case 0: 00050 return .5*(1. - xi); 00051 00052 case 1: 00053 return .5*(1. + xi); 00054 00055 default: 00056 libmesh_error_msg("Invalid shape function index i = " << i); 00057 } 00058 } 00059 00060 00061 00062 // Lagrange quadratics 00063 case SECOND: 00064 { 00065 libmesh_assert_less (i, 3); 00066 00067 switch (i) 00068 { 00069 case 0: 00070 return .5*xi*(xi - 1.); 00071 00072 case 1: 00073 return .5*xi*(xi + 1); 00074 00075 case 2: 00076 return (1. - xi*xi); 00077 00078 default: 00079 libmesh_error_msg("Invalid shape function index i = " << i); 00080 } 00081 } 00082 00083 00084 00085 // Lagrange cubics 00086 case THIRD: 00087 { 00088 libmesh_assert_less (i, 4); 00089 00090 switch (i) 00091 { 00092 case 0: 00093 return 9./16.*(1./9.-xi*xi)*(xi-1.); 00094 00095 case 1: 00096 return -9./16.*(1./9.-xi*xi)*(xi+1.); 00097 00098 case 2: 00099 return 27./16.*(1.-xi*xi)*(1./3.-xi); 00100 00101 case 3: 00102 return 27./16.*(1.-xi*xi)*(1./3.+xi); 00103 00104 default: 00105 libmesh_error_msg("Invalid shape function index i = " << i); 00106 } 00107 } 00108 00109 default: 00110 libmesh_error_msg("ERROR: Unsupported polynomial order = " << order); 00111 } 00112 00113 libmesh_error_msg("We'll never get here!"); 00114 return 0.; 00115 } 00116 00117 00118 00119 template <> 00120 Real FE<1,LAGRANGE>::shape(const Elem* elem, 00121 const Order order, 00122 const unsigned int i, 00123 const Point& p) 00124 { 00125 libmesh_assert(elem); 00126 00127 return FE<1,LAGRANGE>::shape(elem->type(), static_cast<Order>(order + elem->p_level()), i, p); 00128 } 00129 00130 00131 00132 template <> 00133 Real FE<1,LAGRANGE>::shape_deriv(const ElemType, 00134 const Order order, 00135 const unsigned int i, 00136 const unsigned int libmesh_dbg_var(j), 00137 const Point& p) 00138 { 00139 // only d()/dxi in 1D! 00140 00141 libmesh_assert_equal_to (j, 0); 00142 00143 const Real xi = p(0); 00144 00145 00146 switch (order) 00147 { 00148 // Lagrange linear shape function derivatives 00149 case FIRST: 00150 { 00151 libmesh_assert_less (i, 2); 00152 00153 switch (i) 00154 { 00155 case 0: 00156 return -.5; 00157 00158 case 1: 00159 return .5; 00160 00161 default: 00162 libmesh_error_msg("Invalid shape function index i = " << i); 00163 } 00164 } 00165 00166 00167 // Lagrange quadratic shape function derivatives 00168 case SECOND: 00169 { 00170 libmesh_assert_less (i, 3); 00171 00172 switch (i) 00173 { 00174 case 0: 00175 return xi-.5; 00176 00177 case 1: 00178 return xi+.5; 00179 00180 case 2: 00181 return -2.*xi; 00182 00183 default: 00184 libmesh_error_msg("Invalid shape function index i = " << i); 00185 } 00186 } 00187 00188 00189 // Lagrange cubic shape function derivatives 00190 case THIRD: 00191 { 00192 libmesh_assert_less (i, 4); 00193 00194 switch (i) 00195 { 00196 case 0: 00197 return -9./16.*(3.*xi*xi-2.*xi-1./9.); 00198 00199 case 1: 00200 return -9./16.*(-3.*xi*xi-2.*xi+1./9.); 00201 00202 case 2: 00203 return 27./16.*(3.*xi*xi-2./3.*xi-1.); 00204 00205 case 3: 00206 return 27./16.*(-3.*xi*xi-2./3.*xi+1.); 00207 00208 default: 00209 libmesh_error_msg("Invalid shape function index i = " << i); 00210 } 00211 } 00212 00213 00214 default: 00215 libmesh_error_msg("ERROR: Unsupported polynomial order = " << order); 00216 } 00217 00218 libmesh_error_msg("We'll never get here!"); 00219 return 0.; 00220 } 00221 00222 00223 00224 template <> 00225 Real FE<1,LAGRANGE>::shape_deriv(const Elem* elem, 00226 const Order order, 00227 const unsigned int i, 00228 const unsigned int j, 00229 const Point& p) 00230 { 00231 libmesh_assert(elem); 00232 00233 return FE<1,LAGRANGE>::shape_deriv(elem->type(), 00234 static_cast<Order>(order + elem->p_level()), i, j, p); 00235 } 00236 00237 00238 00239 00240 template <> 00241 Real FE<1,LAGRANGE>::shape_second_deriv(const ElemType, 00242 const Order order, 00243 const unsigned int i, 00244 const unsigned int libmesh_dbg_var(j), 00245 const Point& p) 00246 { 00247 // Don't need to switch on j. 1D shape functions 00248 // depend on xi only! 00249 00250 const Real xi = p(0); 00251 libmesh_assert_equal_to (j, 0); 00252 00253 switch (order) 00254 { 00255 // linear Lagrange shape functions 00256 case FIRST: 00257 { 00258 // All second derivatives of linears are zero.... 00259 return 0.; 00260 } 00261 00262 // quadratic Lagrange shape functions 00263 case SECOND: 00264 { 00265 switch (i) 00266 { 00267 case 0: 00268 return 1.; 00269 00270 case 1: 00271 return 1.; 00272 00273 case 2: 00274 return -2.; 00275 00276 default: 00277 libmesh_error_msg("Invalid shape function index i = " << i); 00278 } 00279 } // end case SECOND 00280 00281 case THIRD: 00282 { 00283 switch (i) 00284 { 00285 case 0: 00286 return -9./16.*(6.*xi-2); 00287 00288 case 1: 00289 return -9./16.*(-6*xi-2.); 00290 00291 case 2: 00292 return 27./16.*(6*xi-2./3.); 00293 00294 case 3: 00295 return 27./16.*(-6*xi-2./3.); 00296 00297 default: 00298 libmesh_error_msg("Invalid shape function index i = " << i); 00299 } 00300 } // end case THIRD 00301 00302 00303 default: 00304 libmesh_error_msg("ERROR: Unsupported polynomial order = " << order); 00305 } // end switch (order) 00306 00307 libmesh_error_msg("We'll never get here!"); 00308 return 0.; 00309 } 00310 00311 00312 00313 template <> 00314 Real FE<1,LAGRANGE>::shape_second_deriv(const Elem* elem, 00315 const Order order, 00316 const unsigned int i, 00317 const unsigned int j, 00318 const Point& p) 00319 { 00320 libmesh_assert(elem); 00321 00322 return FE<1,LAGRANGE>::shape_second_deriv(elem->type(), 00323 static_cast<Order>(order + elem->p_level()), i, j, p); 00324 } 00325 00326 } // namespace libMesh