lanczos.hpp 76 KB

12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576777879808182838485868788899091929394959697989910010110210310410510610710810911011111211311411511611711811912012112212312412512612712812913013113213313413513613713813914014114214314414514614714814915015115215315415515615715815916016116216316416516616716816917017117217317417517617717817918018118218318418518618718818919019119219319419519619719819920020120220320420520620720820921021121221321421521621721821922022122222322422522622722822923023123223323423523623723823924024124224324424524624724824925025125225325425525625725825926026126226326426526626726826927027127227327427527627727827928028128228328428528628728828929029129229329429529629729829930030130230330430530630730830931031131231331431531631731831932032132232332432532632732832933033133233333433533633733833934034134234334434534634734834935035135235335435535635735835936036136236336436536636736836937037137237337437537637737837938038138238338438538638738838939039139239339439539639739839940040140240340440540640740840941041141241341441541641741841942042142242342442542642742842943043143243343443543643743843944044144244344444544644744844945045145245345445545645745845946046146246346446546646746846947047147247347447547647747847948048148248348448548648748848949049149249349449549649749849950050150250350450550650750850951051151251351451551651751851952052152252352452552652752852953053153253353453553653753853954054154254354454554654754854955055155255355455555655755855956056156256356456556656756856957057157257357457557657757857958058158258358458558658758858959059159259359459559659759859960060160260360460560660760860961061161261361461561661761861962062162262362462562662762862963063163263363463563663763863964064164264364464564664764864965065165265365465565665765865966066166266366466566666766866967067167267367467567667767867968068168268368468568668768868969069169269369469569669769869970070170270370470570670770870971071171271371471571671771871972072172272372472572672772872973073173273373473573673773873974074174274374474574674774874975075175275375475575675775875976076176276376476576676776876977077177277377477577677777877978078178278378478578678778878979079179279379479579679779879980080180280380480580680780880981081181281381481581681781881982082182282382482582682782882983083183283383483583683783883984084184284384484584684784884985085185285385485585685785885986086186286386486586686786886987087187287387487587687787887988088188288388488588688788888989089189289389489589689789889990090190290390490590690790890991091191291391491591691791891992092192292392492592692792892993093193293393493593693793893994094194294394494594694794894995095195295395495595695795895996096196296396496596696796896997097197297397497597697797897998098198298398498598698798898999099199299399499599699799899910001001100210031004100510061007100810091010101110121013101410151016101710181019102010211022102310241025102610271028102910301031103210331034103510361037103810391040104110421043104410451046104710481049105010511052105310541055105610571058105910601061106210631064106510661067106810691070107110721073107410751076107710781079108010811082108310841085108610871088108910901091109210931094109510961097109810991100110111021103110411051106110711081109111011111112111311141115111611171118111911201121112211231124112511261127112811291130113111321133113411351136113711381139114011411142114311441145114611471148114911501151115211531154115511561157115811591160116111621163116411651166116711681169117011711172117311741175117611771178117911801181118211831184118511861187118811891190119111921193119411951196119711981199120012011202120312041205120612071208120912101211121212131214121512161217121812191220122112221223122412251226122712281229123012311232123312341235123612371238123912401241124212431244124512461247124812491250125112521253125412551256125712581259126012611262126312641265126612671268126912701271127212731274127512761277127812791280128112821283128412851286128712881289129012911292129312941295129612971298129913001301130213031304130513061307
  1. // (C) Copyright John Maddock 2006.
  2. // Use, modification and distribution are subject to the
  3. // Boost Software License, Version 1.0. (See accompanying file
  4. // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
  5. #ifndef BOOST_MATH_SPECIAL_FUNCTIONS_LANCZOS
  6. #define BOOST_MATH_SPECIAL_FUNCTIONS_LANCZOS
  7. #ifdef _MSC_VER
  8. #pragma once
  9. #endif
  10. #include <boost/config.hpp>
  11. #include <boost/math/tools/big_constant.hpp>
  12. #include <boost/mpl/if.hpp>
  13. #include <boost/limits.hpp>
  14. #include <boost/cstdint.hpp>
  15. #include <boost/math/tools/rational.hpp>
  16. #include <boost/math/policies/policy.hpp>
  17. #include <boost/mpl/less_equal.hpp>
  18. #include <limits.h>
  19. #if defined(__GNUC__) && defined(BOOST_MATH_USE_FLOAT128)
  20. //
  21. // This is the only way we can avoid
  22. // warning: non-standard suffix on floating constant [-Wpedantic]
  23. // when building with -Wall -pedantic. Neither __extension__
  24. // nor #pragma dianostic ignored work :(
  25. //
  26. #pragma GCC system_header
  27. #endif
  28. namespace boost{ namespace math{ namespace lanczos{
  29. //
  30. // Individual lanczos approximations start here.
  31. //
  32. // Optimal values for G for each N are taken from
  33. // http://web.mala.bc.ca/pughg/phdThesis/phdThesis.pdf,
  34. // as are the theoretical error bounds.
  35. //
  36. // Constants calculated using the method described by Godfrey
  37. // http://my.fit.edu/~gabdo/gamma.txt and elaborated by Toth at
  38. // http://www.rskey.org/gamma.htm using NTL::RR at 1000 bit precision.
  39. //
  40. // Begin with a small helper to force initialization of constants prior
  41. // to main. This makes the constant initialization thread safe, even
  42. // when called with a user-defined number type.
  43. //
  44. template <class Lanczos, class T>
  45. struct lanczos_initializer
  46. {
  47. struct init
  48. {
  49. init()
  50. {
  51. T t(1);
  52. Lanczos::lanczos_sum(t);
  53. Lanczos::lanczos_sum_expG_scaled(t);
  54. Lanczos::lanczos_sum_near_1(t);
  55. Lanczos::lanczos_sum_near_2(t);
  56. Lanczos::g();
  57. }
  58. void force_instantiate()const{}
  59. };
  60. static const init initializer;
  61. static void force_instantiate()
  62. {
  63. initializer.force_instantiate();
  64. }
  65. };
  66. template <class Lanczos, class T>
  67. typename lanczos_initializer<Lanczos, T>::init const lanczos_initializer<Lanczos, T>::initializer;
  68. //
  69. // Lanczos Coefficients for N=6 G=5.581
  70. // Max experimental error (with arbitary precision arithmetic) 9.516e-12
  71. // Generated with compiler: Microsoft Visual C++ version 8.0 on Win32 at Mar 23 2006
  72. //
  73. struct lanczos6 : public mpl::int_<35>
  74. {
  75. //
  76. // Produces slightly better than float precision when evaluated at
  77. // double precision:
  78. //
  79. template <class T>
  80. static T lanczos_sum(const T& z)
  81. {
  82. lanczos_initializer<lanczos6, T>::force_instantiate(); // Ensure our constants get initialized before main()
  83. static const T num[6] = {
  84. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 8706.349592549009182288174442774377925882)),
  85. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 8523.650341121874633477483696775067709735)),
  86. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 3338.029219476423550899999750161289306564)),
  87. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 653.6424994294008795995653541449610986791)),
  88. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 63.99951844938187085666201263218840287667)),
  89. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 2.506628274631006311133031631822390264407))
  90. };
  91. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint16_t) denom[6] = {
  92. static_cast<boost::uint16_t>(0u),
  93. static_cast<boost::uint16_t>(24u),
  94. static_cast<boost::uint16_t>(50u),
  95. static_cast<boost::uint16_t>(35u),
  96. static_cast<boost::uint16_t>(10u),
  97. static_cast<boost::uint16_t>(1u)
  98. };
  99. return boost::math::tools::evaluate_rational(num, denom, z);
  100. }
  101. template <class T>
  102. static T lanczos_sum_expG_scaled(const T& z)
  103. {
  104. lanczos_initializer<lanczos6, T>::force_instantiate(); // Ensure our constants get initialized before main()
  105. static const T num[6] = {
  106. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 32.81244541029783471623665933780748627823)),
  107. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 32.12388941444332003446077108933558534361)),
  108. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 12.58034729455216106950851080138931470954)),
  109. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 2.463444478353241423633780693218408889251)),
  110. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 0.2412010548258800231126240760264822486599)),
  111. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 0.009446967704539249494420221613134244048319))
  112. };
  113. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint16_t) denom[6] = {
  114. static_cast<boost::uint16_t>(0u),
  115. static_cast<boost::uint16_t>(24u),
  116. static_cast<boost::uint16_t>(50u),
  117. static_cast<boost::uint16_t>(35u),
  118. static_cast<boost::uint16_t>(10u),
  119. static_cast<boost::uint16_t>(1u)
  120. };
  121. return boost::math::tools::evaluate_rational(num, denom, z);
  122. }
  123. template<class T>
  124. static T lanczos_sum_near_1(const T& dz)
  125. {
  126. lanczos_initializer<lanczos6, T>::force_instantiate(); // Ensure our constants get initialized before main()
  127. static const T d[5] = {
  128. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 2.044879010930422922760429926121241330235)),
  129. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, -2.751366405578505366591317846728753993668)),
  130. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 1.02282965224225004296750609604264824677)),
  131. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, -0.09786124911582813985028889636665335893627)),
  132. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 0.0009829742267506615183144364420540766510112)),
  133. };
  134. T result = 0;
  135. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  136. {
  137. result += (-d[k-1]*dz)/(k*dz + k*k);
  138. }
  139. return result;
  140. }
  141. template<class T>
  142. static T lanczos_sum_near_2(const T& dz)
  143. {
  144. lanczos_initializer<lanczos6, T>::force_instantiate(); // Ensure our constants get initialized before main()
  145. static const T d[5] = {
  146. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 5.748142489536043490764289256167080091892)),
  147. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, -7.734074268282457156081021756682138251825)),
  148. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 2.875167944990511006997713242805893543947)),
  149. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, -0.2750873773533504542306766137703788781776)),
  150. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 35, 0.002763134585812698552178368447708846850353)),
  151. };
  152. T result = 0;
  153. T z = dz + 2;
  154. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  155. {
  156. result += (-d[k-1]*dz)/(z + k*z + k*k - 1);
  157. }
  158. return result;
  159. }
  160. static double g(){ return 5.581000000000000405009359383257105946541; }
  161. };
  162. //
  163. // Lanczos Coefficients for N=11 G=10.900511
  164. // Max experimental error (with arbitary precision arithmetic) 2.16676e-19
  165. // Generated with compiler: Microsoft Visual C++ version 8.0 on Win32 at Mar 23 2006
  166. //
  167. struct lanczos11 : public mpl::int_<60>
  168. {
  169. //
  170. // Produces slightly better than double precision when evaluated at
  171. // extended-double precision:
  172. //
  173. template <class T>
  174. static T lanczos_sum(const T& z)
  175. {
  176. lanczos_initializer<lanczos11, T>::force_instantiate(); // Ensure our constants get initialized before main()
  177. static const T num[11] = {
  178. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 38474670393.31776828316099004518914832218)),
  179. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 36857665043.51950660081971227404959150474)),
  180. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 15889202453.72942008945006665994637853242)),
  181. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 4059208354.298834770194507810788393801607)),
  182. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 680547661.1834733286087695557084801366446)),
  183. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 78239755.00312005289816041245285376206263)),
  184. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 6246580.776401795264013335510453568106366)),
  185. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 341986.3488721347032223777872763188768288)),
  186. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 12287.19451182455120096222044424100527629)),
  187. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 261.6140441641668190791708576058805625502)),
  188. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 2.506628274631000502415573855452633787834))
  189. };
  190. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint32_t) denom[11] = {
  191. static_cast<boost::uint32_t>(0u),
  192. static_cast<boost::uint32_t>(362880u),
  193. static_cast<boost::uint32_t>(1026576u),
  194. static_cast<boost::uint32_t>(1172700u),
  195. static_cast<boost::uint32_t>(723680u),
  196. static_cast<boost::uint32_t>(269325u),
  197. static_cast<boost::uint32_t>(63273u),
  198. static_cast<boost::uint32_t>(9450u),
  199. static_cast<boost::uint32_t>(870u),
  200. static_cast<boost::uint32_t>(45u),
  201. static_cast<boost::uint32_t>(1u)
  202. };
  203. return boost::math::tools::evaluate_rational(num, denom, z);
  204. }
  205. template <class T>
  206. static T lanczos_sum_expG_scaled(const T& z)
  207. {
  208. lanczos_initializer<lanczos11, T>::force_instantiate(); // Ensure our constants get initialized before main()
  209. static const T num[11] = {
  210. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 709811.662581657956893540610814842699825)),
  211. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 679979.847415722640161734319823103390728)),
  212. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 293136.785721159725251629480984140341656)),
  213. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 74887.5403291467179935942448101441897121)),
  214. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 12555.29058241386295096255111537516768137)),
  215. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 1443.42992444170669746078056942194198252)),
  216. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 115.2419459613734722083208906727972935065)),
  217. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 6.30923920573262762719523981992008976989)),
  218. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 0.2266840463022436475495508977579735223818)),
  219. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 0.004826466289237661857584712046231435101741)),
  220. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 0.4624429436045378766270459638520555557321e-4))
  221. };
  222. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint32_t) denom[11] = {
  223. static_cast<boost::uint32_t>(0u),
  224. static_cast<boost::uint32_t>(362880u),
  225. static_cast<boost::uint32_t>(1026576u),
  226. static_cast<boost::uint32_t>(1172700u),
  227. static_cast<boost::uint32_t>(723680u),
  228. static_cast<boost::uint32_t>(269325u),
  229. static_cast<boost::uint32_t>(63273u),
  230. static_cast<boost::uint32_t>(9450u),
  231. static_cast<boost::uint32_t>(870u),
  232. static_cast<boost::uint32_t>(45u),
  233. static_cast<boost::uint32_t>(1u)
  234. };
  235. return boost::math::tools::evaluate_rational(num, denom, z);
  236. }
  237. template<class T>
  238. static T lanczos_sum_near_1(const T& dz)
  239. {
  240. lanczos_initializer<lanczos11, T>::force_instantiate(); // Ensure our constants get initialized before main()
  241. static const T d[10] = {
  242. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 4.005853070677940377969080796551266387954)),
  243. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, -13.17044315127646469834125159673527183164)),
  244. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 17.19146865350790353683895137079288129318)),
  245. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, -11.36446409067666626185701599196274701126)),
  246. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 4.024801119349323770107694133829772634737)),
  247. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, -0.7445703262078094128346501724255463005006)),
  248. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 0.06513861351917497265045550019547857713172)),
  249. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, -0.00217899958561830354633560009312512312758)),
  250. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 0.17655204574495137651670832229571934738e-4)),
  251. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, -0.1036282091079938047775645941885460820853e-7)),
  252. };
  253. T result = 0;
  254. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  255. {
  256. result += (-d[k-1]*dz)/(k*dz + k*k);
  257. }
  258. return result;
  259. }
  260. template<class T>
  261. static T lanczos_sum_near_2(const T& dz)
  262. {
  263. lanczos_initializer<lanczos11, T>::force_instantiate(); // Ensure our constants get initialized before main()
  264. static const T d[10] = {
  265. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 19.05889633808148715159575716844556056056)),
  266. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, -62.66183664701721716960978577959655644762)),
  267. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 81.7929198065004751699057192860287512027)),
  268. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, -54.06941772964234828416072865069196553015)),
  269. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 19.14904664790693019642068229478769661515)),
  270. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, -3.542488556926667589704590409095331790317)),
  271. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 0.3099140334815639910894627700232804503017)),
  272. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, -0.01036716187296241640634252431913030440825)),
  273. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, 0.8399926504443119927673843789048514017761e-4)),
  274. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 60, -0.493038376656195010308610694048822561263e-7)),
  275. };
  276. T result = 0;
  277. T z = dz + 2;
  278. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  279. {
  280. result += (-d[k-1]*dz)/(z + k*z + k*k - 1);
  281. }
  282. return result;
  283. }
  284. static double g(){ return 10.90051099999999983936049829935654997826; }
  285. };
  286. //
  287. // Lanczos Coefficients for N=13 G=13.144565
  288. // Max experimental error (with arbitary precision arithmetic) 9.2213e-23
  289. // Generated with compiler: Microsoft Visual C++ version 8.0 on Win32 at Mar 23 2006
  290. //
  291. struct lanczos13 : public mpl::int_<72>
  292. {
  293. //
  294. // Produces slightly better than extended-double precision when evaluated at
  295. // higher precision:
  296. //
  297. template <class T>
  298. static T lanczos_sum(const T& z)
  299. {
  300. lanczos_initializer<lanczos13, T>::force_instantiate(); // Ensure our constants get initialized before main()
  301. static const T num[13] = {
  302. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 44012138428004.60895436261759919070125699)),
  303. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 41590453358593.20051581730723108131357995)),
  304. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 18013842787117.99677796276038389462742949)),
  305. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 4728736263475.388896889723995205703970787)),
  306. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 837910083628.4046470415724300225777912264)),
  307. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 105583707273.4299344907359855510105321192)),
  308. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 9701363618.494999493386608345339104922694)),
  309. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 654914397.5482052641016767125048538245644)),
  310. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 32238322.94213356530668889463945849409184)),
  311. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 1128514.219497091438040721811544858643121)),
  312. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 26665.79378459858944762533958798805525125)),
  313. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 381.8801248632926870394389468349331394196)),
  314. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 2.506628274631000502415763426076722427007))
  315. };
  316. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint32_t) denom[13] = {
  317. static_cast<boost::uint32_t>(0u),
  318. static_cast<boost::uint32_t>(39916800u),
  319. static_cast<boost::uint32_t>(120543840u),
  320. static_cast<boost::uint32_t>(150917976u),
  321. static_cast<boost::uint32_t>(105258076u),
  322. static_cast<boost::uint32_t>(45995730u),
  323. static_cast<boost::uint32_t>(13339535u),
  324. static_cast<boost::uint32_t>(2637558u),
  325. static_cast<boost::uint32_t>(357423u),
  326. static_cast<boost::uint32_t>(32670u),
  327. static_cast<boost::uint32_t>(1925u),
  328. static_cast<boost::uint32_t>(66u),
  329. static_cast<boost::uint32_t>(1u)
  330. };
  331. return boost::math::tools::evaluate_rational(num, denom, z);
  332. }
  333. template <class T>
  334. static T lanczos_sum_expG_scaled(const T& z)
  335. {
  336. lanczos_initializer<lanczos13, T>::force_instantiate(); // Ensure our constants get initialized before main()
  337. static const T num[13] = {
  338. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 86091529.53418537217994842267760536134841)),
  339. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 81354505.17858011242874285785316135398567)),
  340. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 35236626.38815461910817650960734605416521)),
  341. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 9249814.988024471294683815872977672237195)),
  342. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 1639024.216687146960253839656643518985826)),
  343. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 206530.8157641225032631778026076868855623)),
  344. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 18976.70193530288915698282139308582105936)),
  345. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 1281.068909912559479885759622791374106059)),
  346. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 63.06093343420234536146194868906771599354)),
  347. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 2.207470909792527638222674678171050209691)),
  348. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 0.05216058694613505427476207805814960742102)),
  349. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 0.0007469903808915448316510079585999893674101)),
  350. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 0.4903180573459871862552197089738373164184e-5))
  351. };
  352. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint32_t) denom[13] = {
  353. static_cast<boost::uint32_t>(0u),
  354. static_cast<boost::uint32_t>(39916800u),
  355. static_cast<boost::uint32_t>(120543840u),
  356. static_cast<boost::uint32_t>(150917976u),
  357. static_cast<boost::uint32_t>(105258076u),
  358. static_cast<boost::uint32_t>(45995730u),
  359. static_cast<boost::uint32_t>(13339535u),
  360. static_cast<boost::uint32_t>(2637558u),
  361. static_cast<boost::uint32_t>(357423u),
  362. static_cast<boost::uint32_t>(32670u),
  363. static_cast<boost::uint32_t>(1925u),
  364. static_cast<boost::uint32_t>(66u),
  365. static_cast<boost::uint32_t>(1u)
  366. };
  367. return boost::math::tools::evaluate_rational(num, denom, z);
  368. }
  369. template<class T>
  370. static T lanczos_sum_near_1(const T& dz)
  371. {
  372. lanczos_initializer<lanczos13, T>::force_instantiate(); // Ensure our constants get initialized before main()
  373. static const T d[12] = {
  374. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 4.832115561461656947793029596285626840312)),
  375. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, -19.86441536140337740383120735104359034688)),
  376. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 33.9927422807443239927197864963170585331)),
  377. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, -31.41520692249765980987427413991250886138)),
  378. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 17.0270866009599345679868972409543597821)),
  379. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, -5.5077216950865501362506920516723682167)),
  380. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 1.037811741948214855286817963800439373362)),
  381. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, -0.106640468537356182313660880481398642811)),
  382. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 0.005276450526660653288757565778182586742831)),
  383. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, -0.0001000935625597121545867453746252064770029)),
  384. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 0.462590910138598083940803704521211569234e-6)),
  385. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, -0.1735307814426389420248044907765671743012e-9)),
  386. };
  387. T result = 0;
  388. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  389. {
  390. result += (-d[k-1]*dz)/(k*dz + k*k);
  391. }
  392. return result;
  393. }
  394. template<class T>
  395. static T lanczos_sum_near_2(const T& dz)
  396. {
  397. lanczos_initializer<lanczos13, T>::force_instantiate(); // Ensure our constants get initialized before main()
  398. static const T d[12] = {
  399. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 26.96979819614830698367887026728396466395)),
  400. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, -110.8705424709385114023884328797900204863)),
  401. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 189.7258846119231466417015694690434770085)),
  402. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, -175.3397202971107486383321670769397356553)),
  403. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 95.03437648691551457087250340903980824948)),
  404. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, -30.7406022781665264273675797983497141978)),
  405. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 5.792405601630517993355102578874590410552)),
  406. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, -0.5951993240669148697377539518639997795831)),
  407. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 0.02944979359164017509944724739946255067671)),
  408. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, -0.0005586586555377030921194246330399163602684)),
  409. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, 0.2581888478270733025288922038673392636029e-5)),
  410. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 72, -0.9685385411006641478305219367315965391289e-9)),
  411. };
  412. T result = 0;
  413. T z = dz + 2;
  414. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  415. {
  416. result += (-d[k-1]*dz)/(z + k*z + k*k - 1);
  417. }
  418. return result;
  419. }
  420. static double g(){ return 13.1445650000000000545696821063756942749; }
  421. };
  422. //
  423. // Lanczos Coefficients for N=22 G=22.61891
  424. // Max experimental error (with arbitary precision arithmetic) 2.9524e-38
  425. // Generated with compiler: Microsoft Visual C++ version 8.0 on Win32 at Mar 23 2006
  426. //
  427. struct lanczos22 : public mpl::int_<120>
  428. {
  429. //
  430. // Produces slightly better than 128-bit long-double precision when
  431. // evaluated at higher precision:
  432. //
  433. template <class T>
  434. static T lanczos_sum(const T& z)
  435. {
  436. lanczos_initializer<lanczos22, T>::force_instantiate(); // Ensure our constants get initialized before main()
  437. static const T num[22] = {
  438. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 46198410803245094237463011094.12173081986)),
  439. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 43735859291852324413622037436.321513777)),
  440. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 19716607234435171720534556386.97481377748)),
  441. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 5629401471315018442177955161.245623932129)),
  442. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 1142024910634417138386281569.245580222392)),
  443. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 175048529315951173131586747.695329230778)),
  444. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 21044290245653709191654675.41581372963167)),
  445. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 2033001410561031998451380.335553678782601)),
  446. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 160394318862140953773928.8736211601848891)),
  447. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 10444944438396359705707.48957290388740896)),
  448. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 565075825801617290121.1466393747967538948)),
  449. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 25475874292116227538.99448534450411942597)),
  450. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 957135055846602154.6720835535232270205725)),
  451. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 29874506304047462.23662392445173880821515)),
  452. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 769651310384737.2749087590725764959689181)),
  453. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 16193289100889.15989633624378404096011797)),
  454. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 273781151680.6807433264462376754578933261)),
  455. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 3630485900.32917021712188739762161583295)),
  456. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 36374352.05577334277856865691538582936484)),
  457. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 258945.7742115532455441786924971194951043)),
  458. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 1167.501919472435718934219997431551246996)),
  459. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 2.50662827463100050241576528481104525333))
  460. };
  461. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint64_t) denom[22] = {
  462. BOOST_MATH_INT_VALUE_SUFFIX(0, uLL),
  463. BOOST_MATH_INT_VALUE_SUFFIX(2432902008176640000, uLL),
  464. BOOST_MATH_INT_VALUE_SUFFIX(8752948036761600000, uLL),
  465. BOOST_MATH_INT_VALUE_SUFFIX(13803759753640704000, uLL),
  466. BOOST_MATH_INT_VALUE_SUFFIX(12870931245150988800, uLL),
  467. BOOST_MATH_INT_VALUE_SUFFIX(8037811822645051776, uLL),
  468. BOOST_MATH_INT_VALUE_SUFFIX(3599979517947607200, uLL),
  469. BOOST_MATH_INT_VALUE_SUFFIX(1206647803780373360, uLL),
  470. BOOST_MATH_INT_VALUE_SUFFIX(311333643161390640, uLL),
  471. BOOST_MATH_INT_VALUE_SUFFIX(63030812099294896, uLL),
  472. BOOST_MATH_INT_VALUE_SUFFIX(10142299865511450, uLL),
  473. BOOST_MATH_INT_VALUE_SUFFIX(1307535010540395, uLL),
  474. BOOST_MATH_INT_VALUE_SUFFIX(135585182899530, uLL),
  475. BOOST_MATH_INT_VALUE_SUFFIX(11310276995381, uLL),
  476. BOOST_MATH_INT_VALUE_SUFFIX(756111184500, uLL),
  477. BOOST_MATH_INT_VALUE_SUFFIX(40171771630, uLL),
  478. BOOST_MATH_INT_VALUE_SUFFIX(1672280820, uLL),
  479. BOOST_MATH_INT_VALUE_SUFFIX(53327946, uLL),
  480. BOOST_MATH_INT_VALUE_SUFFIX(1256850, uLL),
  481. BOOST_MATH_INT_VALUE_SUFFIX(20615, uLL),
  482. BOOST_MATH_INT_VALUE_SUFFIX(210, uLL),
  483. BOOST_MATH_INT_VALUE_SUFFIX(1, uLL)
  484. };
  485. return boost::math::tools::evaluate_rational(num, denom, z);
  486. }
  487. template <class T>
  488. static T lanczos_sum_expG_scaled(const T& z)
  489. {
  490. lanczos_initializer<lanczos22, T>::force_instantiate(); // Ensure our constants get initialized before main()
  491. static const T num[22] = {
  492. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 6939996264376682180.277485395074954356211)),
  493. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 6570067992110214451.87201438870245659384)),
  494. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 2961859037444440551.986724631496417064121)),
  495. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 845657339772791245.3541226499766163431651)),
  496. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 171556737035449095.2475716923888737881837)),
  497. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 26296059072490867.7822441885603400926007)),
  498. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 3161305619652108.433798300149816829198706)),
  499. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 305400596026022.4774396904484542582526472)),
  500. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 24094681058862.55120507202622377623528108)),
  501. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 1569055604375.919477574824168939428328839)),
  502. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 84886558909.02047889339710230696942513159)),
  503. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 3827024985.166751989686050643579753162298)),
  504. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 143782298.9273215199098728674282885500522)),
  505. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 4487794.24541641841336786238909171265944)),
  506. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 115618.2025760830513505888216285273541959)),
  507. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 2432.580773108508276957461757328744780439)),
  508. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 41.12782532742893597168530008461874360191)),
  509. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.5453771709477689805460179187388702295792)),
  510. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.005464211062612080347167337964166505282809)),
  511. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.388992321263586767037090706042788910953e-4)),
  512. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.1753839324538447655939518484052327068859e-6)),
  513. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.3765495513732730583386223384116545391759e-9))
  514. };
  515. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint64_t) denom[22] = {
  516. BOOST_MATH_INT_VALUE_SUFFIX(0, uLL),
  517. BOOST_MATH_INT_VALUE_SUFFIX(2432902008176640000, uLL),
  518. BOOST_MATH_INT_VALUE_SUFFIX(8752948036761600000, uLL),
  519. BOOST_MATH_INT_VALUE_SUFFIX(13803759753640704000, uLL),
  520. BOOST_MATH_INT_VALUE_SUFFIX(12870931245150988800, uLL),
  521. BOOST_MATH_INT_VALUE_SUFFIX(8037811822645051776, uLL),
  522. BOOST_MATH_INT_VALUE_SUFFIX(3599979517947607200, uLL),
  523. BOOST_MATH_INT_VALUE_SUFFIX(1206647803780373360, uLL),
  524. BOOST_MATH_INT_VALUE_SUFFIX(311333643161390640, uLL),
  525. BOOST_MATH_INT_VALUE_SUFFIX(63030812099294896, uLL),
  526. BOOST_MATH_INT_VALUE_SUFFIX(10142299865511450, uLL),
  527. BOOST_MATH_INT_VALUE_SUFFIX(1307535010540395, uLL),
  528. BOOST_MATH_INT_VALUE_SUFFIX(135585182899530, uLL),
  529. BOOST_MATH_INT_VALUE_SUFFIX(11310276995381, uLL),
  530. BOOST_MATH_INT_VALUE_SUFFIX(756111184500, uLL),
  531. BOOST_MATH_INT_VALUE_SUFFIX(40171771630, uLL),
  532. BOOST_MATH_INT_VALUE_SUFFIX(1672280820, uLL),
  533. BOOST_MATH_INT_VALUE_SUFFIX(53327946, uLL),
  534. BOOST_MATH_INT_VALUE_SUFFIX(1256850, uLL),
  535. BOOST_MATH_INT_VALUE_SUFFIX(20615, uLL),
  536. BOOST_MATH_INT_VALUE_SUFFIX(210, uLL),
  537. BOOST_MATH_INT_VALUE_SUFFIX(1, uLL)
  538. };
  539. return boost::math::tools::evaluate_rational(num, denom, z);
  540. }
  541. template<class T>
  542. static T lanczos_sum_near_1(const T& dz)
  543. {
  544. lanczos_initializer<lanczos22, T>::force_instantiate(); // Ensure our constants get initialized before main()
  545. static const T d[21] = {
  546. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 8.318998691953337183034781139546384476554)),
  547. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -63.15415991415959158214140353299240638675)),
  548. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 217.3108224383632868591462242669081540163)),
  549. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -448.5134281386108366899784093610397354889)),
  550. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 619.2903759363285456927248474593012711346)),
  551. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -604.1630177420625418522025080080444177046)),
  552. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 428.8166750424646119935047118287362193314)),
  553. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -224.6988753721310913866347429589434550302)),
  554. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 87.32181627555510833499451817622786940961)),
  555. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -25.07866854821128965662498003029199058098)),
  556. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 5.264398125689025351448861011657789005392)),
  557. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -0.792518936256495243383586076579921559914)),
  558. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.08317448364744713773350272460937904691566)),
  559. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -0.005845345166274053157781068150827567998882)),
  560. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.0002599412126352082483326238522490030412391)),
  561. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -0.6748102079670763884917431338234783496303e-5)),
  562. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.908824383434109002762325095643458603605e-7)),
  563. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -0.5299325929309389890892469299969669579725e-9)),
  564. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.994306085859549890267983602248532869362e-12)),
  565. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -0.3499893692975262747371544905820891835298e-15)),
  566. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.7260746353663365145454867069182884694961e-20)),
  567. };
  568. T result = 0;
  569. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  570. {
  571. result += (-d[k-1]*dz)/(k*dz + k*k);
  572. }
  573. return result;
  574. }
  575. template<class T>
  576. static T lanczos_sum_near_2(const T& dz)
  577. {
  578. lanczos_initializer<lanczos22, T>::force_instantiate(); // Ensure our constants get initialized before main()
  579. static const T d[21] = {
  580. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 75.39272007105208086018421070699575462226)),
  581. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -572.3481967049935412452681346759966390319)),
  582. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 1969.426202741555335078065370698955484358)),
  583. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -4064.74968778032030891520063865996757519)),
  584. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 5612.452614138013929794736248384309574814)),
  585. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -5475.357667500026172903620177988213902339)),
  586. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 3886.243614216111328329547926490398103492)),
  587. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -2036.382026072125407192448069428134470564)),
  588. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 791.3727954936062108045551843636692287652)),
  589. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -227.2808432388436552794021219198885223122)),
  590. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 47.70974355562144229897637024320739257284)),
  591. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -7.182373807798293545187073539819697141572)),
  592. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.7537866989631514559601547530490976100468)),
  593. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -0.05297470142240154822658739758236594717787)),
  594. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.00235577330936380542539812701472320434133)),
  595. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -0.6115613067659273118098229498679502138802e-4)),
  596. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.8236417010170941915758315020695551724181e-6)),
  597. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -0.4802628430993048190311242611330072198089e-8)),
  598. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.9011113376981524418952720279739624707342e-11)),
  599. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, -0.3171854152689711198382455703658589996796e-14)),
  600. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 120, 0.6580207998808093935798753964580596673177e-19)),
  601. };
  602. T result = 0;
  603. T z = dz + 2;
  604. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  605. {
  606. result += (-d[k-1]*dz)/(z + k*z + k*k - 1);
  607. }
  608. return result;
  609. }
  610. static double g(){ return 22.61890999999999962710717227309942245483; }
  611. };
  612. //
  613. // Lanczos Coefficients for N=6 G=1.428456135094165802001953125
  614. // Max experimental error (with arbitary precision arithmetic) 8.111667e-8
  615. // Generated with compiler: Microsoft Visual C++ version 8.0 on Win32 at Mar 23 2006
  616. //
  617. struct lanczos6m24 : public mpl::int_<24>
  618. {
  619. //
  620. // Use for float precision, when evaluated as a float:
  621. //
  622. template <class T>
  623. static T lanczos_sum(const T& z)
  624. {
  625. static const T num[6] = {
  626. static_cast<T>(58.52061591769095910314047740215847630266L),
  627. static_cast<T>(182.5248962595894264831189414768236280862L),
  628. static_cast<T>(211.0971093028510041839168287718170827259L),
  629. static_cast<T>(112.2526547883668146736465390902227161763L),
  630. static_cast<T>(27.5192015197455403062503721613097825345L),
  631. static_cast<T>(2.50662858515256974113978724717473206342L)
  632. };
  633. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint16_t) denom[6] = {
  634. static_cast<boost::uint16_t>(0u),
  635. static_cast<boost::uint16_t>(24u),
  636. static_cast<boost::uint16_t>(50u),
  637. static_cast<boost::uint16_t>(35u),
  638. static_cast<boost::uint16_t>(10u),
  639. static_cast<boost::uint16_t>(1u)
  640. };
  641. return boost::math::tools::evaluate_rational(num, denom, z);
  642. }
  643. template <class T>
  644. static T lanczos_sum_expG_scaled(const T& z)
  645. {
  646. static const T num[6] = {
  647. static_cast<T>(14.0261432874996476619570577285003839357L),
  648. static_cast<T>(43.74732405540314316089531289293124360129L),
  649. static_cast<T>(50.59547402616588964511581430025589038612L),
  650. static_cast<T>(26.90456680562548195593733429204228910299L),
  651. static_cast<T>(6.595765571169314946316366571954421695196L),
  652. static_cast<T>(0.6007854010515290065101128585795542383721L)
  653. };
  654. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint16_t) denom[6] = {
  655. static_cast<boost::uint16_t>(0u),
  656. static_cast<boost::uint16_t>(24u),
  657. static_cast<boost::uint16_t>(50u),
  658. static_cast<boost::uint16_t>(35u),
  659. static_cast<boost::uint16_t>(10u),
  660. static_cast<boost::uint16_t>(1u)
  661. };
  662. return boost::math::tools::evaluate_rational(num, denom, z);
  663. }
  664. template<class T>
  665. static T lanczos_sum_near_1(const T& dz)
  666. {
  667. static const T d[5] = {
  668. static_cast<T>(0.4922488055204602807654354732674868442106L),
  669. static_cast<T>(0.004954497451132152436631238060933905650346L),
  670. static_cast<T>(-0.003374784572167105840686977985330859371848L),
  671. static_cast<T>(0.001924276018962061937026396537786414831385L),
  672. static_cast<T>(-0.00056533046336427583708166383712907694434L),
  673. };
  674. T result = 0;
  675. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  676. {
  677. result += (-d[k-1]*dz)/(k*dz + k*k);
  678. }
  679. return result;
  680. }
  681. template<class T>
  682. static T lanczos_sum_near_2(const T& dz)
  683. {
  684. static const T d[5] = {
  685. static_cast<T>(0.6534966888520080645505805298901130485464L),
  686. static_cast<T>(0.006577461728560758362509168026049182707101L),
  687. static_cast<T>(-0.004480276069269967207178373559014835978161L),
  688. static_cast<T>(0.00255461870648818292376982818026706528842L),
  689. static_cast<T>(-0.000750517993690428370380996157470900204524L),
  690. };
  691. T result = 0;
  692. T z = dz + 2;
  693. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  694. {
  695. result += (-d[k-1]*dz)/(z + k*z + k*k - 1);
  696. }
  697. return result;
  698. }
  699. static double g(){ return 1.428456135094165802001953125; }
  700. };
  701. //
  702. // Lanczos Coefficients for N=13 G=6.024680040776729583740234375
  703. // Max experimental error (with arbitary precision arithmetic) 1.196214e-17
  704. // Generated with compiler: Microsoft Visual C++ version 8.0 on Win32 at Mar 23 2006
  705. //
  706. struct lanczos13m53 : public mpl::int_<53>
  707. {
  708. //
  709. // Use for double precision, when evaluated as a double:
  710. //
  711. template <class T>
  712. static T lanczos_sum(const T& z)
  713. {
  714. static const T num[13] = {
  715. static_cast<T>(23531376880.41075968857200767445163675473L),
  716. static_cast<T>(42919803642.64909876895789904700198885093L),
  717. static_cast<T>(35711959237.35566804944018545154716670596L),
  718. static_cast<T>(17921034426.03720969991975575445893111267L),
  719. static_cast<T>(6039542586.35202800506429164430729792107L),
  720. static_cast<T>(1439720407.311721673663223072794912393972L),
  721. static_cast<T>(248874557.8620541565114603864132294232163L),
  722. static_cast<T>(31426415.58540019438061423162831820536287L),
  723. static_cast<T>(2876370.628935372441225409051620849613599L),
  724. static_cast<T>(186056.2653952234950402949897160456992822L),
  725. static_cast<T>(8071.672002365816210638002902272250613822L),
  726. static_cast<T>(210.8242777515793458725097339207133627117L),
  727. static_cast<T>(2.506628274631000270164908177133837338626L)
  728. };
  729. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint32_t) denom[13] = {
  730. static_cast<boost::uint32_t>(0u),
  731. static_cast<boost::uint32_t>(39916800u),
  732. static_cast<boost::uint32_t>(120543840u),
  733. static_cast<boost::uint32_t>(150917976u),
  734. static_cast<boost::uint32_t>(105258076u),
  735. static_cast<boost::uint32_t>(45995730u),
  736. static_cast<boost::uint32_t>(13339535u),
  737. static_cast<boost::uint32_t>(2637558u),
  738. static_cast<boost::uint32_t>(357423u),
  739. static_cast<boost::uint32_t>(32670u),
  740. static_cast<boost::uint32_t>(1925u),
  741. static_cast<boost::uint32_t>(66u),
  742. static_cast<boost::uint32_t>(1u)
  743. };
  744. return boost::math::tools::evaluate_rational(num, denom, z);
  745. }
  746. template <class T>
  747. static T lanczos_sum_expG_scaled(const T& z)
  748. {
  749. static const T num[13] = {
  750. static_cast<T>(56906521.91347156388090791033559122686859L),
  751. static_cast<T>(103794043.1163445451906271053616070238554L),
  752. static_cast<T>(86363131.28813859145546927288977868422342L),
  753. static_cast<T>(43338889.32467613834773723740590533316085L),
  754. static_cast<T>(14605578.08768506808414169982791359218571L),
  755. static_cast<T>(3481712.15498064590882071018964774556468L),
  756. static_cast<T>(601859.6171681098786670226533699352302507L),
  757. static_cast<T>(75999.29304014542649875303443598909137092L),
  758. static_cast<T>(6955.999602515376140356310115515198987526L),
  759. static_cast<T>(449.9445569063168119446858607650988409623L),
  760. static_cast<T>(19.51992788247617482847860966235652136208L),
  761. static_cast<T>(0.5098416655656676188125178644804694509993L),
  762. static_cast<T>(0.006061842346248906525783753964555936883222L)
  763. };
  764. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint32_t) denom[13] = {
  765. static_cast<boost::uint32_t>(0u),
  766. static_cast<boost::uint32_t>(39916800u),
  767. static_cast<boost::uint32_t>(120543840u),
  768. static_cast<boost::uint32_t>(150917976u),
  769. static_cast<boost::uint32_t>(105258076u),
  770. static_cast<boost::uint32_t>(45995730u),
  771. static_cast<boost::uint32_t>(13339535u),
  772. static_cast<boost::uint32_t>(2637558u),
  773. static_cast<boost::uint32_t>(357423u),
  774. static_cast<boost::uint32_t>(32670u),
  775. static_cast<boost::uint32_t>(1925u),
  776. static_cast<boost::uint32_t>(66u),
  777. static_cast<boost::uint32_t>(1u)
  778. };
  779. return boost::math::tools::evaluate_rational(num, denom, z);
  780. }
  781. template<class T>
  782. static T lanczos_sum_near_1(const T& dz)
  783. {
  784. static const T d[12] = {
  785. static_cast<T>(2.208709979316623790862569924861841433016L),
  786. static_cast<T>(-3.327150580651624233553677113928873034916L),
  787. static_cast<T>(1.483082862367253753040442933770164111678L),
  788. static_cast<T>(-0.1993758927614728757314233026257810172008L),
  789. static_cast<T>(0.004785200610085071473880915854204301886437L),
  790. static_cast<T>(-0.1515973019871092388943437623825208095123e-5L),
  791. static_cast<T>(-0.2752907702903126466004207345038327818713e-7L),
  792. static_cast<T>(0.3075580174791348492737947340039992829546e-7L),
  793. static_cast<T>(-0.1933117898880828348692541394841204288047e-7L),
  794. static_cast<T>(0.8690926181038057039526127422002498960172e-8L),
  795. static_cast<T>(-0.2499505151487868335680273909354071938387e-8L),
  796. static_cast<T>(0.3394643171893132535170101292240837927725e-9L),
  797. };
  798. T result = 0;
  799. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  800. {
  801. result += (-d[k-1]*dz)/(k*dz + k*k);
  802. }
  803. return result;
  804. }
  805. template<class T>
  806. static T lanczos_sum_near_2(const T& dz)
  807. {
  808. static const T d[12] = {
  809. static_cast<T>(6.565936202082889535528455955485877361223L),
  810. static_cast<T>(-9.8907772644920670589288081640128194231L),
  811. static_cast<T>(4.408830289125943377923077727900630927902L),
  812. static_cast<T>(-0.5926941084905061794445733628891024027949L),
  813. static_cast<T>(0.01422519127192419234315002746252160965831L),
  814. static_cast<T>(-0.4506604409707170077136555010018549819192e-5L),
  815. static_cast<T>(-0.8183698410724358930823737982119474130069e-7L),
  816. static_cast<T>(0.9142922068165324132060550591210267992072e-7L),
  817. static_cast<T>(-0.5746670642147041587497159649318454348117e-7L),
  818. static_cast<T>(0.2583592566524439230844378948704262291927e-7L),
  819. static_cast<T>(-0.7430396708998719707642735577238449585822e-8L),
  820. static_cast<T>(0.1009141566987569892221439918230042368112e-8L),
  821. };
  822. T result = 0;
  823. T z = dz + 2;
  824. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  825. {
  826. result += (-d[k-1]*dz)/(z + k*z + k*k - 1);
  827. }
  828. return result;
  829. }
  830. static double g(){ return 6.024680040776729583740234375; }
  831. };
  832. //
  833. // Lanczos Coefficients for N=17 G=12.2252227365970611572265625
  834. // Max experimental error (with arbitary precision arithmetic) 2.7699e-26
  835. // Generated with compiler: Microsoft Visual C++ version 8.0 on Win32 at Mar 23 2006
  836. //
  837. struct lanczos17m64 : public mpl::int_<64>
  838. {
  839. //
  840. // Use for extended-double precision, when evaluated as an extended-double:
  841. //
  842. template <class T>
  843. static T lanczos_sum(const T& z)
  844. {
  845. lanczos_initializer<lanczos17m64, T>::force_instantiate(); // Ensure our constants get initialized before main()
  846. static const T num[17] = {
  847. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 553681095419291969.2230556393350368550504)),
  848. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 731918863887667017.2511276782146694632234)),
  849. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 453393234285807339.4627124634539085143364)),
  850. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 174701893724452790.3546219631779712198035)),
  851. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 46866125995234723.82897281620357050883077)),
  852. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 9281280675933215.169109622777099699054272)),
  853. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 1403600894156674.551057997617468721789536)),
  854. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 165345984157572.7305349809894046783973837)),
  855. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 15333629842677.31531822808737907246817024)),
  856. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 1123152927963.956626161137169462874517318)),
  857. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 64763127437.92329018717775593533620578237)),
  858. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 2908830362.657527782848828237106640944457)),
  859. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 99764700.56999856729959383751710026787811)),
  860. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 2525791.604886139959837791244686290089331)),
  861. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 44516.94034970167828580039370201346554872)),
  862. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 488.0063567520005730476791712814838113252)),
  863. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 2.50662827463100050241576877135758834683))
  864. };
  865. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint64_t) denom[17] = {
  866. BOOST_MATH_INT_VALUE_SUFFIX(0, uLL),
  867. BOOST_MATH_INT_VALUE_SUFFIX(1307674368000, uLL),
  868. BOOST_MATH_INT_VALUE_SUFFIX(4339163001600, uLL),
  869. BOOST_MATH_INT_VALUE_SUFFIX(6165817614720, uLL),
  870. BOOST_MATH_INT_VALUE_SUFFIX(5056995703824, uLL),
  871. BOOST_MATH_INT_VALUE_SUFFIX(2706813345600, uLL),
  872. BOOST_MATH_INT_VALUE_SUFFIX(1009672107080, uLL),
  873. BOOST_MATH_INT_VALUE_SUFFIX(272803210680, uLL),
  874. BOOST_MATH_INT_VALUE_SUFFIX(54631129553, uLL),
  875. BOOST_MATH_INT_VALUE_SUFFIX(8207628000, uLL),
  876. BOOST_MATH_INT_VALUE_SUFFIX(928095740, uLL),
  877. BOOST_MATH_INT_VALUE_SUFFIX(78558480, uLL),
  878. BOOST_MATH_INT_VALUE_SUFFIX(4899622, uLL),
  879. BOOST_MATH_INT_VALUE_SUFFIX(218400, uLL),
  880. BOOST_MATH_INT_VALUE_SUFFIX(6580, uLL),
  881. BOOST_MATH_INT_VALUE_SUFFIX(120, uLL),
  882. BOOST_MATH_INT_VALUE_SUFFIX(1, uLL)
  883. };
  884. return boost::math::tools::evaluate_rational(num, denom, z);
  885. }
  886. template <class T>
  887. static T lanczos_sum_expG_scaled(const T& z)
  888. {
  889. lanczos_initializer<lanczos17m64, T>::force_instantiate(); // Ensure our constants get initialized before main()
  890. static const T num[17] = {
  891. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 2715894658327.717377557655133124376674911)),
  892. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 3590179526097.912105038525528721129550434)),
  893. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 2223966599737.814969312127353235818710172)),
  894. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 856940834518.9562481809925866825485883417)),
  895. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 229885871668.749072933597446453399395469)),
  896. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 45526171687.54610815813502794395753410032)),
  897. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 6884887713.165178784550917647709216424823)),
  898. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 811048596.1407531864760282453852372777439)),
  899. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 75213915.96540822314499613623119501704812)),
  900. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 5509245.417224265151697527957954952830126)),
  901. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 317673.5368435419126714931842182369574221)),
  902. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 14268.27989845035520147014373320337523596)),
  903. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 489.3618720403263670213909083601787814792)),
  904. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 12.38941330038454449295883217865458609584)),
  905. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 0.2183627389504614963941574507281683147897)),
  906. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 0.002393749522058449186690627996063983095463)),
  907. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 0.1229541408909435212800785616808830746135e-4))
  908. };
  909. static const BOOST_MATH_INT_TABLE_TYPE(T, boost::uint64_t) denom[17] = {
  910. BOOST_MATH_INT_VALUE_SUFFIX(0, uLL),
  911. BOOST_MATH_INT_VALUE_SUFFIX(1307674368000, uLL),
  912. BOOST_MATH_INT_VALUE_SUFFIX(4339163001600, uLL),
  913. BOOST_MATH_INT_VALUE_SUFFIX(6165817614720, uLL),
  914. BOOST_MATH_INT_VALUE_SUFFIX(5056995703824, uLL),
  915. BOOST_MATH_INT_VALUE_SUFFIX(2706813345600, uLL),
  916. BOOST_MATH_INT_VALUE_SUFFIX(1009672107080, uLL),
  917. BOOST_MATH_INT_VALUE_SUFFIX(272803210680, uLL),
  918. BOOST_MATH_INT_VALUE_SUFFIX(54631129553, uLL),
  919. BOOST_MATH_INT_VALUE_SUFFIX(8207628000, uLL),
  920. BOOST_MATH_INT_VALUE_SUFFIX(928095740, uLL),
  921. BOOST_MATH_INT_VALUE_SUFFIX(78558480, uLL),
  922. BOOST_MATH_INT_VALUE_SUFFIX(4899622, uLL),
  923. BOOST_MATH_INT_VALUE_SUFFIX(218400, uLL),
  924. BOOST_MATH_INT_VALUE_SUFFIX(6580, uLL),
  925. BOOST_MATH_INT_VALUE_SUFFIX(120, uLL),
  926. BOOST_MATH_INT_VALUE_SUFFIX(1, uLL)
  927. };
  928. return boost::math::tools::evaluate_rational(num, denom, z);
  929. }
  930. template<class T>
  931. static T lanczos_sum_near_1(const T& dz)
  932. {
  933. lanczos_initializer<lanczos17m64, T>::force_instantiate(); // Ensure our constants get initialized before main()
  934. static const T d[16] = {
  935. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 4.493645054286536365763334986866616581265)),
  936. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -16.95716370392468543800733966378143997694)),
  937. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 26.19196892983737527836811770970479846644)),
  938. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -21.3659076437988814488356323758179283908)),
  939. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 9.913992596774556590710751047594507535764)),
  940. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -2.62888300018780199210536267080940382158)),
  941. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 0.3807056693542503606384861890663080735588)),
  942. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -0.02714647489697685807340312061034730486958)),
  943. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 0.0007815484715461206757220527133967191796747)),
  944. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -0.6108630817371501052576880554048972272435e-5)),
  945. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 0.5037380238864836824167713635482801545086e-8)),
  946. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -0.1483232144262638814568926925964858237006e-13)),
  947. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 0.1346609158752142460943888149156716841693e-14)),
  948. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -0.660492688923978805315914918995410340796e-15)),
  949. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 0.1472114697343266749193617793755763792681e-15)),
  950. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -0.1410901942033374651613542904678399264447e-16)),
  951. };
  952. T result = 0;
  953. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  954. {
  955. result += (-d[k-1]*dz)/(k*dz + k*k);
  956. }
  957. return result;
  958. }
  959. template<class T>
  960. static T lanczos_sum_near_2(const T& dz)
  961. {
  962. lanczos_initializer<lanczos17m64, T>::force_instantiate(); // Ensure our constants get initialized before main()
  963. static const T d[16] = {
  964. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 23.56409085052261327114594781581930373708)),
  965. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -88.92116338946308797946237246006238652361)),
  966. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 137.3472822086847596961177383569603988797)),
  967. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -112.0400438263562152489272966461114852861)),
  968. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 51.98768915202973863076166956576777843805)),
  969. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -13.78552090862799358221343319574970124948)),
  970. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 1.996371068830872830250406773917646121742)),
  971. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -0.1423525874909934506274738563671862576161)),
  972. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 0.004098338646046865122459664947239111298524)),
  973. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -0.3203286637326511000882086573060433529094e-4)),
  974. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 0.2641536751640138646146395939004587594407e-7)),
  975. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -0.7777876663062235617693516558976641009819e-13)),
  976. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 0.7061443477097101636871806229515157914789e-14)),
  977. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -0.3463537849537988455590834887691613484813e-14)),
  978. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 0.7719578215795234036320348283011129450595e-15)),
  979. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -0.7398586479708476329563577384044188912075e-16)),
  980. };
  981. T result = 0;
  982. T z = dz + 2;
  983. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  984. {
  985. result += (-d[k-1]*dz)/(z + k*z + k*k - 1);
  986. }
  987. return result;
  988. }
  989. static double g(){ return 12.2252227365970611572265625; }
  990. };
  991. //
  992. // Lanczos Coefficients for N=24 G=20.3209821879863739013671875
  993. // Max experimental error (with arbitary precision arithmetic) 1.0541e-38
  994. // Generated with compiler: Microsoft Visual C++ version 8.0 on Win32 at Mar 23 2006
  995. //
  996. struct lanczos24m113 : public mpl::int_<113>
  997. {
  998. //
  999. // Use for long-double precision, when evaluated as an long-double:
  1000. //
  1001. template <class T>
  1002. static T lanczos_sum(const T& z)
  1003. {
  1004. lanczos_initializer<lanczos24m113, T>::force_instantiate(); // Ensure our constants get initialized before main()
  1005. static const T num[24] = {
  1006. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 2029889364934367661624137213253.22102954656825019111612712252027267955023987678816620961507)),
  1007. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 2338599599286656537526273232565.2727349714338768161421882478417543004440597874814359063158)),
  1008. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1288527989493833400335117708406.3953711906175960449186720680201425446299360322830739180195)),
  1009. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 451779745834728745064649902914.550539158066332484594436145043388809847364393288132164411521)),
  1010. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 113141284461097964029239556815.291212318665536114012605167994061291631013303788706545334708)),
  1011. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 21533689802794625866812941616.7509064680880468667055339259146063256555368135236149614592432)),
  1012. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 3235510315314840089932120340.71494940111731241353655381919722177496659303550321056514776757)),
  1013. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 393537392344185475704891959.081297108513472083749083165179784098220158201055270548272414314)),
  1014. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 39418265082950435024868801.5005452240816902251477336582325944930252142622315101857742955673)),
  1015. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 3290158764187118871697791.05850632319194734270969161036889516414516566453884272345518372696)),
  1016. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 230677110449632078321772.618245845856640677845629174549731890660612368500786684333975350954)),
  1017. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 13652233645509183190158.5916189185218250859402806777406323001463296297553612462737044693697)),
  1018. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 683661466754325350495.216655026531202476397782296585200982429378069417193575896602446904762)),
  1019. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 28967871782219334117.0122379171041074970463982134039409352925258212207710168851968215545064)),
  1020. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1036104088560167006.2022834098572346459442601718514554488352117620272232373622553429728555)),
  1021. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 31128490785613152.8380102669349814751268126141105475287632676569913936040772990253369753962)),
  1022. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 779327504127342.536207878988196814811198475410572992436243686674896894543126229424358472541)),
  1023. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 16067543181294.643350688789124777020407337133926174150582333950666044399234540521336771876)),
  1024. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 268161795520.300916569439413185778557212729611517883948634711190170998896514639936969855484)),
  1025. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 3533216359.10528191668842486732408440112703691790824611391987708562111396961696753452085068)),
  1026. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 35378979.5479656110614685178752543826919239614088343789329169535932709470588426584501652577)),
  1027. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 253034.881362204346444503097491737872930637147096453940375713745904094735506180552724766444)),
  1028. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1151.61895453463992438325318456328526085882924197763140514450975619271382783957699017875304)),
  1029. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 2.50662827463100050241576528481104515966515623051532908941425544355490413900497467936202516))
  1030. };
  1031. static const T denom[24] = {
  1032. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.0)),
  1033. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.112400072777760768e22)),
  1034. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.414847677933545472e22)),
  1035. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 6756146673770930688000.0)),
  1036. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 6548684852703068697600.0)),
  1037. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 4280722865357147142912.0)),
  1038. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 2021687376910682741568.0)),
  1039. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 720308216440924653696.0)),
  1040. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 199321978221066137360.0)),
  1041. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 43714229649594412832.0)),
  1042. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 7707401101297361068.0)),
  1043. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1103230881185949736.0)),
  1044. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 129006659818331295.0)),
  1045. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 12363045847086207.0)),
  1046. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 971250460939913.0)),
  1047. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 62382416421941.0)),
  1048. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 3256091103430.0)),
  1049. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 136717357942.0)),
  1050. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 4546047198.0)),
  1051. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 116896626.0)),
  1052. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 2240315.0)),
  1053. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 30107.0)),
  1054. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 253.0)),
  1055. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1.0))
  1056. };
  1057. return boost::math::tools::evaluate_rational(num, denom, z);
  1058. }
  1059. template <class T>
  1060. static T lanczos_sum_expG_scaled(const T& z)
  1061. {
  1062. lanczos_initializer<lanczos24m113, T>::force_instantiate(); // Ensure our constants get initialized before main()
  1063. static const T num[24] = {
  1064. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 3035162425359883494754.02878223286972654682199012688209026810841953293372712802258398358538)),
  1065. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 3496756894406430103600.16057175075063458536101374170860226963245118484234495645518505519827)),
  1066. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1926652656689320888654.01954015145958293168365236755537645929361841917596501251362171653478)),
  1067. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 675517066488272766316.083023742440619929434602223726894748181327187670231286180156444871912)),
  1068. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 169172853104918752780.086262749564831660238912144573032141700464995906149421555926000038492)),
  1069. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 32197935167225605785.6444116302160245528783954573163541751756353183343357329404208062043808)),
  1070. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 4837849542714083249.37587447454818124327561966323276633775195138872820542242539845253171632)),
  1071. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 588431038090493242.308438203986649553459461798968819276505178004064031201740043314534404158)),
  1072. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 58939585141634058.6206417889192563007809470547755357240808035714047014324843817783741669733)),
  1073. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 4919561837722192.82991866530802080996138070630296720420704876654726991998309206256077395868)),
  1074. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 344916580244240.407442753122831512004021081677987651622305356145640394384006997569631719101)),
  1075. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 20413302960687.8250598845969238472629322716685686993835561234733641729957841485003560103066)),
  1076. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1022234822943.78400752460970689311934727763870970686747383486600540378889311406851534545789)),
  1077. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 43313787191.9821354846952908076307094286897439975815501673706144217246093900159173598852503)),
  1078. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1549219505.59667418528481770869280437577581951167003505825834192510436144666564648361001914)),
  1079. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 46544421.1998761919380541579358096705925369145324466147390364674998568485110045455014967149)),
  1080. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1165278.06807504975090675074910052763026564833951579556132777702952882101173607903881127542)),
  1081. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 24024.759267256769471083727721827405338569868270177779485912486668586611981795179894572115)),
  1082. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 400.965008113421955824358063769761286758463521789765880962939528760888853281920872064838918)),
  1083. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 5.28299015654478269617039029170846385138134929147421558771949982217659507918482272439717603)),
  1084. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.0528999024412510102409256676599360516359062802002483877724963720047531347449011629466149805)),
  1085. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.000378346710654740685454266569593414561162134092347356968516522170279688139165340746957511115)),
  1086. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.172194142179211139195966608011235161516824700287310869949928393345257114743230967204370963e-5)),
  1087. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.374799931707148855771381263542708435935402853962736029347951399323367765509988401336565436e-8))
  1088. };
  1089. static const T denom[24] = {
  1090. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.0)),
  1091. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.112400072777760768e22)),
  1092. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.414847677933545472e22)),
  1093. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 6756146673770930688000.0)),
  1094. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 6548684852703068697600.0)),
  1095. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 4280722865357147142912.0)),
  1096. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 2021687376910682741568.0)),
  1097. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 720308216440924653696.0)),
  1098. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 199321978221066137360.0)),
  1099. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 43714229649594412832.0)),
  1100. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 7707401101297361068.0)),
  1101. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1103230881185949736.0)),
  1102. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 129006659818331295.0)),
  1103. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 12363045847086207.0)),
  1104. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 971250460939913.0)),
  1105. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 62382416421941.0)),
  1106. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 3256091103430.0)),
  1107. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 136717357942.0)),
  1108. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 4546047198.0)),
  1109. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 116896626.0)),
  1110. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 2240315.0)),
  1111. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 30107.0)),
  1112. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 253.0)),
  1113. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1.0))
  1114. };
  1115. return boost::math::tools::evaluate_rational(num, denom, z);
  1116. }
  1117. template<class T>
  1118. static T lanczos_sum_near_1(const T& dz)
  1119. {
  1120. lanczos_initializer<lanczos24m113, T>::force_instantiate(); // Ensure our constants get initialized before main()
  1121. static const T d[23] = {
  1122. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 7.4734083002469026177867421609938203388868806387315406134072298925733950040583068760685908)),
  1123. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -50.4225805042247530267317342133388132970816607563062253708655085754357843064134941138154171)),
  1124. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 152.288200621747008570784082624444625293884063492396162110698238568311211546361189979357019)),
  1125. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -271.894959539150384169327513139846971255640842175739337449692360299099322742181325023644769)),
  1126. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 319.240102980202312307047586791116902719088581839891008532114107693294261542869734803906793)),
  1127. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -259.493144143048088289689500935518073716201741349569864988870534417890269467336454358361499)),
  1128. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 149.747518319689708813209645403067832020714660918583227716408482877303972685262557460145835)),
  1129. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -61.9261301009341333289187201425188698128684426428003249782448828881580630606817104372760037)),
  1130. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 18.3077524177286961563937379403377462608113523887554047531153187277072451294845795496072365)),
  1131. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -3.82011322251948043097070160584761236869363471824695092089556195047949392738162970152230254)),
  1132. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.549382685505691522516705902336780999493262538301283190963770663549981309645795228539620711)),
  1133. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -0.0524814679715180697633723771076668718265358076235229045603747927518423453658004287459638024)),
  1134. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.00315392664003333528534120626687784812050217700942910879712808180705014754163256855643360698)),
  1135. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -0.000110098373127648510519799564665442121339511198561008748083409549601095293123407080388658329)),
  1136. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.19809382866681658224945717689377373458866950897791116315219376038432014207446832310901893e-5)),
  1137. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -0.152278977408600291408265615203504153130482270424202400677280558181047344681214058227949755e-7)),
  1138. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.364344768076106268872239259083188037615571711218395765792787047015406264051536972018235217e-10)),
  1139. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -0.148897510480440424971521542520683536298361220674662555578951242811522959610991621951203526e-13)),
  1140. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.261199241161582662426512749820666625442516059622425213340053324061794752786482115387573582e-18)),
  1141. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -0.780072664167099103420998436901014795601783313858454665485256897090476089641613851903791529e-24)),
  1142. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.303465867587106629530056603454807425512962762653755513440561256044986695349304176849392735e-24)),
  1143. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -0.615420597971283870342083342286977366161772327800327789325710571275345878439656918541092056e-25)),
  1144. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.499641233843540749369110053005439398774706583601830828776209650445427083113181961630763702e-26)),
  1145. };
  1146. T result = 0;
  1147. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  1148. {
  1149. result += (-d[k-1]*dz)/(k*dz + k*k);
  1150. }
  1151. return result;
  1152. }
  1153. template<class T>
  1154. static T lanczos_sum_near_2(const T& dz)
  1155. {
  1156. lanczos_initializer<lanczos24m113, T>::force_instantiate(); // Ensure our constants get initialized before main()
  1157. static const T d[23] = {
  1158. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 61.4165001061101455341808888883960361969557848005400286332291451422461117307237198559485365)),
  1159. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -414.372973678657049667308134761613915623353625332248315105320470271523320700386200587519147)),
  1160. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1251.50505818554680171298972755376376836161706773644771875668053742215217922228357204561873)),
  1161. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -2234.43389421602399514176336175766511311493214354568097811220122848998413358085613880612158)),
  1162. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 2623.51647746991904821899989145639147785427273427135380151752779100215839537090464785708684)),
  1163. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -2132.51572435428751962745870184529534443305617818870214348386131243463614597272260797772423)),
  1164. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 1230.62572059218405766499842067263311220019173335523810725664442147670956427061920234820189)),
  1165. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -508.90919151163744999377586956023909888833335885805154492270846381061182696305011395981929)),
  1166. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 150.453184562246579758706538566480316921938628645961177699894388251635886834047343195475395)),
  1167. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -31.3937061525822497422230490071156186113405446381476081565548185848237169870395131828731397)),
  1168. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 4.51482916590287954234936829724231512565732528859217337795452389161322923867318809206313688)),
  1169. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -0.431292919341108177524462194102701868233551186625103849565527515201492276412231365776131952)),
  1170. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.0259189820815586225636729971503340447445001375909094681698918294680345547092233915092128323)),
  1171. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -0.000904788882557558697594884691337532557729219389814315972435534723829065673966567231504429712)),
  1172. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.162793589759218213439218473348810982422449144393340433592232065020562974405674317564164312e-4)),
  1173. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -0.125142926178202562426432039899709511761368233479483128438847484617555752948755923647214487e-6)),
  1174. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.299418680048132583204152682950097239197934281178261879500770485862852229898797687301941982e-9)),
  1175. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -0.122364035267809278675627784883078206654408225276233049012165202996967011873995261617995421e-12)),
  1176. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.21465364366598631597052073538883430194257709353929022544344097235100199405814005393447785e-17)),
  1177. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -0.641064035802907518396608051803921688237330857546406669209280666066685733941549058513986818e-23)),
  1178. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.249388374622173329690271566855185869111237201309011956145463506483151054813346819490278951e-23)),
  1179. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, -0.505752900177513489906064295001851463338022055787536494321532352380960774349054239257683149e-24)),
  1180. static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 113, 0.410605371184590959139968810080063542546949719163227555918846829816144878123034347778284006e-25)),
  1181. };
  1182. T result = 0;
  1183. T z = dz + 2;
  1184. for(unsigned k = 1; k <= sizeof(d)/sizeof(d[0]); ++k)
  1185. {
  1186. result += (-d[k-1]*dz)/(z + k*z + k*k - 1);
  1187. }
  1188. return result;
  1189. }
  1190. static double g(){ return 20.3209821879863739013671875; }
  1191. };
  1192. //
  1193. // placeholder for no lanczos info available:
  1194. //
  1195. struct undefined_lanczos : public mpl::int_<INT_MAX - 1> { };
  1196. #if 0
  1197. #ifndef BOOST_NO_LIMITS_COMPILE_TIME_CONSTANTS
  1198. #define BOOST_MATH_FLT_DIGITS ::std::numeric_limits<float>::digits
  1199. #define BOOST_MATH_DBL_DIGITS ::std::numeric_limits<double>::digits
  1200. #define BOOST_MATH_LDBL_DIGITS ::std::numeric_limits<long double>::digits
  1201. #else
  1202. #define BOOST_MATH_FLT_DIGITS FLT_MANT_DIG
  1203. #define BOOST_MATH_DBL_DIGITS DBL_MANT_DIG
  1204. #define BOOST_MATH_LDBL_DIGITS LDBL_MANT_DIG
  1205. #endif
  1206. #endif
  1207. typedef mpl::list<
  1208. lanczos6m24,
  1209. /* lanczos6, */
  1210. lanczos13m53,
  1211. /* lanczos13, */
  1212. lanczos17m64,
  1213. lanczos24m113,
  1214. lanczos22,
  1215. undefined_lanczos> lanczos_list;
  1216. template <class Real, class Policy>
  1217. struct lanczos
  1218. {
  1219. typedef typename mpl::if_<
  1220. typename mpl::less_equal<
  1221. typename policies::precision<Real, Policy>::type,
  1222. mpl::int_<0>
  1223. >::type,
  1224. mpl::int_<INT_MAX - 2>,
  1225. typename policies::precision<Real, Policy>::type
  1226. >::type target_precision;
  1227. typedef typename mpl::deref<typename mpl::find_if<
  1228. lanczos_list,
  1229. mpl::less_equal<target_precision, mpl::_1> >::type>::type type;
  1230. };
  1231. } // namespace lanczos
  1232. } // namespace math
  1233. } // namespace boost
  1234. #if !defined(_CRAYC) && !defined(__CUDACC__) && (!defined(__GNUC__) || (__GNUC__ > 3) || ((__GNUC__ == 3) && (__GNUC_MINOR__ > 3)))
  1235. #if ((defined(_M_IX86_FP) && (_M_IX86_FP >= 2)) || defined(__SSE2__) || defined(_M_AMD64) || defined(_M_X64)) && !defined(_MANAGED)
  1236. #include <boost/math/special_functions/detail/lanczos_sse2.hpp>
  1237. #endif
  1238. #endif
  1239. #endif // BOOST_MATH_SPECIAL_FUNCTIONS_LANCZOS