$extrastylesheet
fe_lagrange_shape_1D.C
Go to the documentation of this file.
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