inf_fe_legendre_eval.c

一个用来实现偏微分方程中网格的计算库

// $Id: inf_fe_legendre_eval.C 2789 2008-04-13 02:24:40Z roystgnr $

// The libMesh Finite Element Library.
// Copyright (C) 2002-2007  Benjamin S. Kirk, John W. Peterson

// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.

// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// Lesser General Public License for more details.

// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA


// Local Includes
#include "libmesh_config.h"
#ifdef ENABLE_INFINITE_ELEMENTS
#include "inf_fe.h"
#include "inf_fe_macro.h"



template <unsigned int Dim, FEFamily T_radial, InfMapType T_map>
Real InfFE<Dim,T_radial,T_map>::eval(const Real v,
                                     const Order, /* not used */
                                     const unsigned int i)
{
  libmesh_assert (-1.-1.e-5 <= v && v < 1.);

  switch (i)
    {
    case 0:
      return 1.;

    case 1:
      return v+1.;

    case 2:
      return 1.5*v*v-1.5;

    case 3:
      return 1.+(-1.5+2.5*v*v)*v;

    case 4:
      return -.625+(-3.75+4.375*v*v)*v*v;

    case 5:
      return 1.+(1.875+(-8.75+7.875*v*v)*v*v)*v;

    case 6:
      return -1.3125+(6.5625+(-19.6875+14.4375*v*v)*v*v)*v*v;

    case 7:
      return 1.+(-2.1875+(19.6875+(-43.3125+26.8125*v*v)*v*v)*v*v)*v;

    case 8:
      return -.7265625+(-9.84375+(54.140625+(-93.84375+50.2734375*v*v)*v*v)*v*v)*v*v;

    case 9:
      return 1.+(2.4609375+(-36.09375+(140.765625+(-201.09375+94.9609375*v*v)*v*v)*v*v)*v*v)*v;

    case 10:
      return -1.24609375+(13.53515625+(-117.3046875+(351.9140625+(-427.32421875+180.42578125*v*v)*v*v)*v*v)*v*v)*v*v;

    case 11:
      return 1.+(-2.70703125+(58.65234375+(-351.9140625+(854.6484375+(-902.12890625+344.44921875*v*v)*v*v)*v*v)*v*v)*v*v)*v;

    case 12:
      return -.7744140625+(-17.595703125+(219.9462890625+(-997.08984375+(2029.7900390625+(-1894.470703125+660.1943359375*v*v)*v*v)*v*v)*v*v)*v*v)*v*v;

    case 13:
      return 1.+(2.9326171875+(-87.978515625+(747.8173828125+(-2706.38671875+(4736.1767578125+(-3961.166015625+1269.6044921875*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v;

    case 14:
      return -1.20947265625+(21.99462890625+(-373.90869140625+(2368.08837890625+(-7104.26513671875+(10893.20654296875+(-8252.42919921875+2448.52294921875*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v;

    case 15:
      return 1.+(-3.14208984375+(124.63623046875+(-1420.85302734375+(7104.26513671875+(-18155.34423828125+(24757.28759765625+(-17139.66064453125+4733.81103515625*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v;

    case 16:
      return -.803619384765625+(-26.707763671875+(592.0220947265625+(-4972.985595703125+(20424.76226806641+(-45388.36059570313+(55703.89709472656+(-35503.58276367188+9171.758880615234*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v;

    case 17:
      return 1.+(3.338470458984375+(-169.149169921875+(2486.492797851563+(-16339.80981445313+(56735.45074462891+(-111407.7941894531+(124262.5396728516+(-73374.07104492188+17804.00253295898*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v;

    case 18:
      return -1.185470581054688+(31.71546936035156+(-888.0331420898438+(9531.555725097656+(-51061.90567016602+(153185.717010498+(-269235.5026245117+(275152.766418457+(-151334.0215301514+34618.89381408691*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v;


    default:
      std::cerr << "bad index i = " << i << std::endl;
      libmesh_error();

    }

  // we never end up here.
  return 0.;
}



template <unsigned int Dim, FEFamily T_radial, InfMapType T_map>
Real InfFE<Dim,T_radial,T_map>::eval_deriv(const Real v,
                                           const Order, /* not used */
                                           const unsigned int i)
{
  libmesh_assert (-1.-1.e-5 <= v && v < 1.);

  switch (i)
    {
    case 0:
      return 0.;

    case 1:
      return 1.;

    case 2:
      return 3.*v;

    case 3:
      return 7.5*v*v-1.5;

    case 4:
      return (-7.5+17.5*v*v)*v;

    case 5:
      return 1.875+(-26.25+39.375*v*v)*v*v;

    case 6:
      return (13.125+(-78.75+86.625*v*v)*v*v)*v;

    case 7:
      return -2.1875+(59.0625+(-216.5625+187.6875*v*v)*v*v)*v*v;

    case 8:
      return (-19.6875+(216.5625+(-563.0625+402.1875*v*v)*v*v)*v*v)*v;

    case 9:
      return 2.4609375+(-108.28125+(703.828125+(-1407.65625+854.6484375*v*v)*v*v)*v*v)*v*v;

    case 10:
      return (27.0703125+(-469.21875+(2111.484375+(-3418.59375+1804.2578125*v*v)*v*v)*v*v)*v*v)*v;

    case 11:
      return -2.70703125+(175.95703125+(-1759.5703125+(5982.5390625+(-8119.16015625+3788.94140625*v*v)*v*v)*v*v)*v*v)*v*v;

    case 12:
      return (-35.19140625+(879.78515625+(-5982.5390625+(16238.3203125+(-18944.70703125+7922.33203125*v*v)*v*v)*v*v)*v*v)*v*v)*v;

    case 13:
      return 2.9326171875+(-263.935546875+(3739.0869140625+(-18944.70703125+(42625.5908203125+(-43572.826171875+16504.8583984375*v*v)*v*v)*v*v)*v*v)*v*v)*v*v;

    case 14:
      return (43.9892578125+(-1495.634765625+(14208.5302734375+(-56834.12109375+(108932.0654296875+(-99029.150390625+34279.3212890625*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v;

    case 15:
      return -3.14208984375+(373.90869140625+(-7104.26513671875+(49729.85595703125+(-163398.0981445313+(272330.1635742188+(-222815.5883789063+71007.16552734375*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v;

    case 16:
      return (-53.41552734375+(2368.08837890625+(-29837.91357421875+(163398.0981445313+(-453883.6059570313+(668446.7651367188+(-497050.1586914063+146748.1420898438*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v;

    case 17:
      return 3.338470458984375+(-507.447509765625+(12432.46398925781+(-114378.6687011719+(510619.0567016602+(-1225485.736083984+(1615413.01574707+(-1100611.065673828+302668.0430603027*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v;

    case 18:
      return (63.43093872070313+(-3552.132568359375+(57189.33435058594+(-408495.2453613281+(1531857.17010498+(-3230826.031494141+(3852138.729858398+(-2421344.344482422+623140.0886535645*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v*v)*v;


    default:
      std::cerr << "bad index i = " << i << std::endl;
      libmesh_error();

    }

  // we never end up here.
  return 0.;
}



//-----------------------------------------------------------------
// Explicit instantiations
// template class InfFE<1,LEGENDRE,CARTESIAN>;
// template class InfFE<2,LEGENDRE,CARTESIAN>;
// template class InfFE<3,LEGENDRE,CARTESIAN>;

// template class InfFE<1,LEGENDRE,SPHERICAL>;
// template class InfFE<2,LEGENDRE,SPHERICAL>;
// template class InfFE<3,LEGENDRE,SPHERICAL>;

// template class InfFE<1,LEGENDRE,ELLIPSOIDAL>;
// template class InfFE<2,LEGENDRE,ELLIPSOIDAL>;
// template class InfFE<3,LEGENDRE,ELLIPSOIDAL>;

template Real InfFE<1,LEGENDRE,CARTESIAN>::eval(const Real,const Order,const unsigned int);
template Real InfFE<1,LEGENDRE,CARTESIAN>::eval_deriv(const Real,const Order,const unsigned int);
template Real InfFE<2,LEGENDRE,CARTESIAN>::eval(const Real,const Order,const unsigned int);
template Real InfFE<2,LEGENDRE,CARTESIAN>::eval_deriv(const Real,const Order,const unsigned int);
template Real InfFE<3,LEGENDRE,CARTESIAN>::eval(const Real,const Order,const unsigned int);
template Real InfFE<3,LEGENDRE,CARTESIAN>::eval_deriv(const Real,const Order,const unsigned int);


#endif //ifdef ENABLE_INFINITE_ELEMENTS