Just finished calculating pi to 30 trillion places.

[Kernelpanic]

Pete, Kama's formula with "n = 1" should lead to this; see screenshot.

If I take it over directly, it also works, but not in the program.

It doesn't work because I don't understand the formula properly, so I don't implement it correctly. That's the problem.
Where is the error in the program?

[Jack]

Kernelpanic
Ramanujan’s Formula for Pi would go something like this, consider it pseudo-code

Code: (Select All)

pie=0
digits=16 'number digits
m=digits/8 'number of terms

for n=0 to m
    tmp=factorial(n)
    tmp=factorial(4*n)/(tmp*tmp*tmp*tmp) 'tmp^4
    pie=pie+tmp*(26390*n+1103)/(396^(4*n))
next
pie=pie*sqr(8)/9801
pie=1/pie
print pie

function factorial(n)
    f=1
    if n=0 then
        factorial=f
        exit function
    end if
    for i=1 to n
        f=f*i
    next
    factorial=f
end function

[Pete]

Kernelpanic
Ramanujan’s Formula for Pi would go something like this, consider it pseudo-code

Code: (Select All)

pie=0
digits=16 'number digits
m=digits/8 'number of terms

for n=0 to m
    tmp=factorial(n)
    tmp=factorial(4*n)/(tmp*tmp*tmp*tmp) 'tmp^4
    pie=pie+tmp*(26390*n+1103)/(396^(4*n))
next
pie=pie*sqr(8)/9801
pie=1/pie
print pie

function factorial(n)
    f=1
    if n=0 then
        factorial=f
        exit function
    end if
    for i=1 to n
        f=f*i
    next
    factorial=f
end function

Nice work at applying the formula, Jack. It works up to 8 iterations before the binary computer math breaks it. Makes me with I had powers and square roots coded for my routine.

I'll hang on to this for some time in the future.

Pete

[Jack]

the simplest and fast way to calculate pi is to use the Gauss–Legendre algorithm https://en.wikipedia.org/wiki/Gauss%E2%8..._algorithm
the Pi - Chudnovsky with binary splitting is much faster here's a version in FreeBasic https://www.freebasic.net/forum/viewtopi...85#p178885

all calculations need to be performed to full precision

Code: (Select All)

   digits=100
   pow2=1
    c0=0: ak = c0: bk = c0: ab = c0: asq = c0
    c1=1: a = c1: ck = c1
    c2=2: b = c2
    c05=.5: sum = c05

    b=sqr(b)
    b=c1/b
    for k=0 to fix(log(digits)*1.44269504088896)
        ak=a+b
        ak=c05*ak
        ab=a*b
        bk=sqr(ab)
        asq=ak*ak
        ck=asq-ab
        pow2=pow2*2
        tmp=ck*pow2
        sum=sum-tmp
        a = ak: b = bk
    next
    tmp=asq/sum
    pie=c2*tmp
    print pie

[JRace]

@KernelPanic: the first constant:

0.0002885855652225477 is being stored into pi as:
0.0002885855792555958




Edited to add: also, (Sqr(8) / 9801) = 0.0002885855792555958 according to QB64.
WolframAlpha says it's: 0.00028858556522254775

I haven't programmed square roots for string math yet, but with a precision calculator, sqrt(8) = 2.8284271247461900976

so with string math: 2.8284271247461900976 / 9801=0.00028858556522254775

0.00028858556522254775 (String math)

0.0002885855792555958  (QB64 reported by JRice)
0.00028858556522254775 (WolframAlpha)

So I would concur WolframAlpha is correct.

Pete

Yeah, looks like a precision issue, but the location of the discrepancy strikes me as odd.

[Pete]

the simplest and fast way to calculate pi is to use the Gauss–Legendre algorithm https://en.wikipedia.org/wiki/Gauss%E2%8..._algorithm
the Pi - Chudnovsky with binary splitting is much faster here's a version in FreeBasic https://www.freebasic.net/forum/viewtopi...85#p178885

all calculations need to be performed to full precision

Code: (Select All)

   digits=100
   pow2=1
    c0=0: ak = c0: bk = c0: ab = c0: asq = c0
    c1=1: a = c1: ck = c1
    c2=2: b = c2
    c05=.5: sum = c05

    b=sqr(b)
    b=c1/b
    for k=0 to fix(log(digits)*1.44269504088896)
        ak=a+b
        ak=c05*ak
        ab=a*b
        bk=sqr(ab)
        asq=ak*ak
        ck=asq-ab
        pow2=pow2*2
        tmp=ck*pow2
        sum=sum-tmp
        a = ak: b = bk
    next
    tmp=asq/sum
    pie=c2*tmp
    print pie

@Jack

I can't say I'm a fan of this one. Digits variable at 100 runs 6 iterations. I tried other iterations with less accurate or completely skewed results. So the...

Code: (Select All)

fix(log(digits)*1.44269504088896)

...is lost on me. at digits = 100 the output was: 3.141697. Is there any way to get that more accurate and produce more digits?

Pete

[Jack]

Pete, that should work but you need to perform all calculations to digits accuracy
there's a quartic version of the Gauss–Legendre algorithm for Pi, that is the accuracy quadruples on each iteration, but I find that's a bit slower https://en.wikipedia.org/wiki/Approximat...algorithms
the quartic version would be something like the following

Code: (Select All)

    y=2 : y=sqr(y)
    a=6 : a=a-4*y : y=y-1
    p=2
    for k=0 to int(log(NUMBER_OF_DIGITS)*0.72134752044448) '1/ln(4)
        y=y*y : y=y*y 'y^4
        y=1-y : y=sqr(y) : y=sqr(y)
        y=(1-y)/(1+y)
        p=p*4
        yk=y+1
        yk=yk*yk : yk=yk*yk
        a=a*yk-p*y*(1+y+y*y)
    next
    pie = 1/a
    print pie

[mdijkens]

A really nice and quick one is this one:

Code: (Select All)

$Console:Only
dig& = Val(Command$)
If dig& <= 0 Then Input "how many digits? ", dig&
Print calcPI$(dig&)

Function calcPI$ (dig&)
  arrSize& = (dig& + 1) Or 1: If arrSize& > 200 Then ten% = 10 Else ten% = 100: arrSize& = (arrSize& + 1) / 2
  Dim POWER(arrSize&), TERM(arrSize&), RESULT(arrSize&)
 
  For pass% = 1 To 2
    For L& = 0 To arrSize&
      POWER(L&) = 0: TERM(L&) = 0: If pass% = 1 Then RESULT(L&) = 0
    Next L&
    POWER(0) = 16 / pass% ^ 2
    mode% = 0: If pass% = 1 Then divisor% = 5 Else divisor% = 239
    GoSub divide
    xponent% = 1: sign% = 3 - 2 * pass%
    Do
      For L& = 0 To arrSize&: TERM(L&) = POWER(L&): Next L& '
      mode% = 1: divisor% = xponent%: GoSub divide
      mode% = (sign% < 0) + Abs(sign% > 0)
      carry% = 0
      For l1& = arrSize& To 0 Step -1
        sum% = RESULT(l1&) + TERM(l1&) * mode% + carry% * mode%: carry% = 0
        If (mode% = 1) * (sum% < ten%) + (mode% = -1) * (sum% >= 0) Then
          RESULT(l1&) = sum%
        Else
          RESULT(l1&) = sum% + mode% * -ten%: carry% = 1
        End If
        'Locate 20, (l1& + 1) * 2: Print _Trim$(Str$(RESULT(l1&)));: _Delay .001
      Next l1&
      xponent% = xponent% + 2: sign% = -sign%
      mode% = 0
      If pass% = 1 Then divisor% = 25 Else divisor% = 239
      GoSub divide
      If pass% = 2 Then GoSub divide
    Loop While zero% <> 0
  Next pass%
  result$ = _Trim$(Str$(RESULT(0))) + "."
  For L& = 1 To arrSize&
    result$ = result$ + _Trim$(Str$(RESULT(L&)))
  Next L&
  calcPI$ = result$: Exit Function

  divide:
  digit% = 0: zero% = 0
  For l1& = 0 To arrSize&
    digit% = digit% + TERM(l1&) * mode% + POWER(l1&) - POWER(l1&) * mode%
    quotient% = Int(digit% / divisor%)
    residue% = digit% Mod divisor%
    zero% = zero% Or (quotient% + residue%)
    If mode% Then TERM(l1&) = quotient% Else POWER(l1&) = quotient%
    digit% = residue% * ten%
  Next l1&
  mode% = 0
  Return
End Function

[Jack]

Pete, you could adapt this code https://qb64forum.alephc.xyz/index.php?P...#msg133377

[Kernelpanic]

Kernelpanic

Ramanujan’s Formula for Pi would go something like this, consider it pseudo-code

Thanks! This works, but only once. "n" and "m" apparently have no meaning.
Actually, 8 new correct digits should appear with each run.