1. /* sdcrpi2.cpp by K.Tsuru */
  2. // function ID 3512 DRADIX, constant pi since version 2.30
  3. /**************************************************************
  4. pi by Ramanujan's second formula
  5. (-1)^n (4n)!(1123+21460n)
  6. 3528/pi = sum_{n=0}^{\inf} ---------------------------
  7. 882^(2n)(4^n n!)^4
  8. B(0)...B(n-1)
  9. = sum_{n=0}^{\inf} --------------- A(n) = s
  10. C(0)...C(n)
  11. pi = 3528*s
  12. See below for A(n),B(k),C(k)
  13. Adding one term the presision rises by eight digits.
  14. See "sdcchpi.cpp" for detail.
  15. Binary Splitting method Nov.14,2016
  16. 49.1(sec) for 4000029 digits
  17. 91.8(sec) 8000029
  18. ***************************************************************/
  19. #ifndef SN_H
  20. #include "sn.h"
  21. #endif
  22. static const SLong D=1123L;
  23. static const SLong Ec=21460L; // To distinguish from E().
  24. static const SLong F=SLong(882L*882L)*SLong(256L);//882^2*4^4
  25. static const SLong G=4L*882L; // =3528
  26. static double upPerTerm = 5.891;
  27. #define BS_RJ_PI2_SLONG 1 // SLong version is faster
  28. #if BS_RJ_PI2_SLONG
  29. ////// SLong version /////// 0.66(sec) for 100000 digits
  30. class BSRamanujanPi2 : public BinarySplitting<SLong> {
  31. public:
  32. /* Constructer : 'L' is upToTerm, 'f' is precision*/
  33. BSRamanujanPi2(long L, long f) : BinarySplitting<SLong>(L, f){}
  34. void setABC(long k, SLong& a, SLong& b, SLong& c) {
  35. // A(k)=(-1)^k*(D+Ek)
  36. if(k&1) a = -D - Ec * k;
  37. else a = D + Ec * k;
  38. // B(k)= 8(k+1)(4k+3)(2k+1)(4k+1)
  39. if(k+1L < LONG_MAX/8L){
  40. long d =8*(k+1L), p = 4L*k+3L, q = 2L*k+1L, r = p-2L;
  41. b= d; b *= p; b *= q; b *= r;
  42. }else{
  43. SLong d =k+1L, p = 4L*k+3L, q = 2L*k+1L, r = p-2L;
  44. b= 8; b *= d; b *= p; b *= q; b *= r;
  45. }
  46. // C(k)=F*k^4
  47. if(k) c = F*Lpow(k, 4);
  48. else c = 1L;
  49. }
  50. SDouble getValue() {
  51. putTogether();
  52. SDouble s = BinarySplitting<SLong>::getValue(true), pi = G*s;
  53. return pi;
  54. }
  55. };
  56. SDouble RamanujanPi2() {
  57. SDouble C;
  58. long f = long(C.EffFig() + C.Hidden()+ 1u)*DFIGURES;
  59. long L = (long)((double)f/upPerTerm);
  60. BSRamanujanPi2 pi(L, f);
  61. return pi.getValue();
  62. }
  63. #else
  64. ////// SDouble version /////// 1.11(sec) for 100000 digits
  65. class BSRamanujanPi2 : public BinarySplitting<SDouble> {
  66. public:
  67. /* Constructer : 'L' is upToTerm, 'f' is precision*/
  68. BSRamanujanPi2(long L, long f) : BinarySplitting<SDouble>(L, f){}
  69. void setABC(long k, SDouble& a, SDouble& b, SDouble& c) {
  70. // A(k)=(-1)^k*(D+Ek)
  71. if(k&1) a = -D - Ec * k;
  72. else a = D + Ec * k;
  73. // B(k)= 8(k+1)(4k+3)(2k+1)(4k+1)
  74. long d =8L*(k+1L), p = 4L*k+3L, q = 2L*k+1L, r = p-2L;
  75. b = d; b *= p; b *= q; b *= r;
  76. // C(k)=F*k^4
  77. if(k) c = F*Lpow(k, 4);
  78. else c = 1L;
  79. }
  80. SDouble getValue() {
  81. putTogether();
  82. SDouble s = BinarySplitting<SDouble>::getValue(true), pi = G*s;
  83. return pi;
  84. }
  85. };
  86. SDouble RamanujanPi2() {
  87. SDouble C;
  88. long f = long( C.EffFig() + C.Hidden() ), L = long(upRate*double(f));
  89. BSRamanujanPi2 pi(L, f);
  90. return pi.getValue();
  91. }
  92. #endif // BS_RJ_PI_SLONG

sdcRjnpi2.cpp : last modifiled at 2017/08/25 15:35:46(3,051 bytes)
created at 2017/10/07 10:21:15
The creation time of this html file is 2017/10/07 10:30:03 (Sat Oct 07 10:30:03 2017).