ntl_rr_lanczos.hpp

Boost provides free peer-reviewed portable C++ source libraries. We emphasize libraries that work

         (12870931245150988800uLL),
         (8037811822645051776uLL),
         (3599979517947607200uLL),
         (1206647803780373360uLL),
         (311333643161390640uLL),
         (63030812099294896uLL),
         (10142299865511450uLL),
         (1307535010540395uLL),
         (135585182899530uLL),
         (11310276995381uLL),
         (756111184500uLL),
         (40171771630uLL),
         (1672280820uLL),
         (53327946uLL),
         (1256850uLL),
         (20615uLL),
         (210uLL),
         (1uLL),
      };
      return boost::math::tools::evaluate_rational(num, denom, z, 22);
   }

   template <class T>
   static T lanczos_sum_expG_scaled(const T& z)
   {
      static const T num[22] = {
         boost::lexical_cast<T>("6939996264376682180.277485395074954356211"),
         boost::lexical_cast<T>("6570067992110214451.87201438870245659384"),
         boost::lexical_cast<T>("2961859037444440551.986724631496417064121"),
         boost::lexical_cast<T>("845657339772791245.3541226499766163431651"),
         boost::lexical_cast<T>("171556737035449095.2475716923888737881837"),
         boost::lexical_cast<T>("26296059072490867.7822441885603400926007"),
         boost::lexical_cast<T>("3161305619652108.433798300149816829198706"),
         boost::lexical_cast<T>("305400596026022.4774396904484542582526472"),
         boost::lexical_cast<T>("24094681058862.55120507202622377623528108"),
         boost::lexical_cast<T>("1569055604375.919477574824168939428328839"),
         boost::lexical_cast<T>("84886558909.02047889339710230696942513159"),
         boost::lexical_cast<T>("3827024985.166751989686050643579753162298"),
         boost::lexical_cast<T>("143782298.9273215199098728674282885500522"),
         boost::lexical_cast<T>("4487794.24541641841336786238909171265944"),
         boost::lexical_cast<T>("115618.2025760830513505888216285273541959"),
         boost::lexical_cast<T>("2432.580773108508276957461757328744780439"),
         boost::lexical_cast<T>("41.12782532742893597168530008461874360191"),
         boost::lexical_cast<T>("0.5453771709477689805460179187388702295792"),
         boost::lexical_cast<T>("0.005464211062612080347167337964166505282809"),
         boost::lexical_cast<T>("0.388992321263586767037090706042788910953e-4"),
         boost::lexical_cast<T>("0.1753839324538447655939518484052327068859e-6"),
         boost::lexical_cast<T>("0.3765495513732730583386223384116545391759e-9"),
      };
      static const boost::uint64_t denom[22] = {
         (0uLL),
         (2432902008176640000uLL),
         (8752948036761600000uLL),
         (13803759753640704000uLL),
         (12870931245150988800uLL),
         (8037811822645051776uLL),
         (3599979517947607200uLL),
         (1206647803780373360uLL),
         (311333643161390640uLL),
         (63030812099294896uLL),
         (10142299865511450uLL),
         (1307535010540395uLL),
         (135585182899530uLL),
         (11310276995381uLL),
         (756111184500uLL),
         (40171771630uLL),
         (1672280820uLL),
         (53327946uLL),
         (1256850uLL),
         (20615uLL),
         (210uLL),
         (1uLL),
      };
      return boost::math::tools::evaluate_rational(num, denom, z, 22);
   }


   template<class T>
   static T lanczos_sum_near_1(const T& dz)
   {
      static const T d[21] = {
         boost::lexical_cast<T>("8.318998691953337183034781139546384476554"),
         boost::lexical_cast<T>("-63.15415991415959158214140353299240638675"),
         boost::lexical_cast<T>("217.3108224383632868591462242669081540163"),
         boost::lexical_cast<T>("-448.5134281386108366899784093610397354889"),
         boost::lexical_cast<T>("619.2903759363285456927248474593012711346"),
         boost::lexical_cast<T>("-604.1630177420625418522025080080444177046"),
         boost::lexical_cast<T>("428.8166750424646119935047118287362193314"),
         boost::lexical_cast<T>("-224.6988753721310913866347429589434550302"),
         boost::lexical_cast<T>("87.32181627555510833499451817622786940961"),
         boost::lexical_cast<T>("-25.07866854821128965662498003029199058098"),
         boost::lexical_cast<T>("5.264398125689025351448861011657789005392"),
         boost::lexical_cast<T>("-0.792518936256495243383586076579921559914"),
         boost::lexical_cast<T>("0.08317448364744713773350272460937904691566"),
         boost::lexical_cast<T>("-0.005845345166274053157781068150827567998882"),
         boost::lexical_cast<T>("0.0002599412126352082483326238522490030412391"),
         boost::lexical_cast<T>("-0.6748102079670763884917431338234783496303e-5"),
         boost::lexical_cast<T>("0.908824383434109002762325095643458603605e-7"),
         boost::lexical_cast<T>("-0.5299325929309389890892469299969669579725e-9"),
         boost::lexical_cast<T>("0.994306085859549890267983602248532869362e-12"),
         boost::lexical_cast<T>("-0.3499893692975262747371544905820891835298e-15"),
         boost::lexical_cast<T>("0.7260746353663365145454867069182884694961e-20"),
      };
      T result = 0;
      for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
      {
         result += (-d[k-1]*dz)/(k*dz + k*k);
      }
      return result;
   }

   template<class T>
   static T lanczos_sum_near_2(const T& dz)
   {
      static const T d[21] = {
         boost::lexical_cast<T>("75.39272007105208086018421070699575462226"),
         boost::lexical_cast<T>("-572.3481967049935412452681346759966390319"),
         boost::lexical_cast<T>("1969.426202741555335078065370698955484358"),
         boost::lexical_cast<T>("-4064.74968778032030891520063865996757519"),
         boost::lexical_cast<T>("5612.452614138013929794736248384309574814"),
         boost::lexical_cast<T>("-5475.357667500026172903620177988213902339"),
         boost::lexical_cast<T>("3886.243614216111328329547926490398103492"),
         boost::lexical_cast<T>("-2036.382026072125407192448069428134470564"),
         boost::lexical_cast<T>("791.3727954936062108045551843636692287652"),
         boost::lexical_cast<T>("-227.2808432388436552794021219198885223122"),
         boost::lexical_cast<T>("47.70974355562144229897637024320739257284"),
         boost::lexical_cast<T>("-7.182373807798293545187073539819697141572"),
         boost::lexical_cast<T>("0.7537866989631514559601547530490976100468"),
         boost::lexical_cast<T>("-0.05297470142240154822658739758236594717787"),
         boost::lexical_cast<T>("0.00235577330936380542539812701472320434133"),
         boost::lexical_cast<T>("-0.6115613067659273118098229498679502138802e-4"),
         boost::lexical_cast<T>("0.8236417010170941915758315020695551724181e-6"),
         boost::lexical_cast<T>("-0.4802628430993048190311242611330072198089e-8"),
         boost::lexical_cast<T>("0.9011113376981524418952720279739624707342e-11"),
         boost::lexical_cast<T>("-0.3171854152689711198382455703658589996796e-14"),
         boost::lexical_cast<T>("0.6580207998808093935798753964580596673177e-19"),
      };
      T result = 0;
      T z = dz + 2;
      for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
      {
         result += (-d[k-1]*dz)/(z + k*z + k*k - 1);
      }
      return result;
   }

   static double g(){ return 22.61890999999999962710717227309942245483; }
};

//
// Lanczos Coefficients for N=31 G=32.08067
// Max experimental error (with arbitary precision arithmetic) 0.162e-52
// Generated with compiler: Microsoft Visual C++ version 8.0 on Win32 at May 9 2006
//
struct lanczos31UDT
{
   template <class T>
   static T lanczos_sum(const T& z)
   {
      static const T num[31] = {
         boost::lexical_cast<T>("0.2579646553333513328235723061836959833277e46"),
         boost::lexical_cast<T>("0.2444796504337453845497419271639377138264e46"),
         boost::lexical_cast<T>("0.1119885499016017172212179730662673475329e46"),
         boost::lexical_cast<T>("0.3301983829072723658949204487793889113715e45"),
         boost::lexical_cast<T>("0.7041171040503851585152895336505379417066e44"),
         boost::lexical_cast<T>("0.1156687509001223855125097826246939403504e44"),
         boost::lexical_cast<T>("1522559363393940883866575697565974893306000"),
         boost::lexical_cast<T>("164914363507650839510801418717701057005700"),
         boost::lexical_cast<T>("14978522943127593263654178827041568394060"),
         boost::lexical_cast<T>("1156707153701375383907746879648168666774"),
         boost::lexical_cast<T>("76739431129980851159755403434593664173.2"),
         boost::lexical_cast<T>("4407916278928188620282281495575981079.306"),
         boost::lexical_cast<T>("220487883931812802092792125175269667.3004"),
         boost::lexical_cast<T>("9644828280794966468052381443992828.433924"),
         boost::lexical_cast<T>("369996467042247229310044531282837.6549068"),
         boost::lexical_cast<T>("12468380890717344610932904378961.13494291"),
         boost::lexical_cast<T>("369289245210898235894444657859.0529720075"),
         boost::lexical_cast<T>("9607992460262594951559461829.34885209022"),
         boost::lexical_cast<T>("219225935074853412540086410.981421315799"),
         boost::lexical_cast<T>("4374309943598658046326340.720767382079549"),
         boost::lexical_cast<T>("76008779092264509404014.10530947173485581"),