Calculating pi by Recursion

[danbaron]

Hi Gang!

Maybe everyone else knew it, but I didn't.

EB can do recursion, i.e., a function can call itself.

Here is a contrived little program, that uses a recursive function, to calculate pi.

(I realize that in this case, and in almost all cases, recursion is less efficient than iteration.
But, I guess, sometimes, recursion is much clearer than iteration (but again, not here).)

Dan.

Code Select Expand


'----------------------------------------------------------------------------------------

'FILE = PICALC.EBA
'PI CALCULATION USING RECURSION
'EB PROGRAM
'RUN AS "CONSOLE EXE".

'----------------------------------------------------------------------------------------

CONST C1 =     8.0
CONST C2 =  9801.0
CONST C3 =  1103.0
CONST C4 = 26390.0
CONST C5 =   396.0

'----------------------------------------------------------------------------------------

DOUBLE EPS = 1.0E-17
UINT64 FACTPROD
DOUBLE PISUM = 0
DOUBLE PI
STRING S

'----------------------------------------------------------------------------------------

PI = CALCPI(0)

'VALUE OF PI, CALCULATED USING FUNCTION "CALCPI".
PRINT USING("##.################", PI)

'VALUE OF PI, CALCULATED USING TRIGONOMETRY, TO CHECK "CALCPI".
PRINT USING("##.################", 4 * ATAN(1.0))

INPUT S
END

'----------------------------------------------------------------------------------------

SUB CALCPI(INT I), DOUBLE

'PI FUNCTION

'PI FORMULA FROM, "NUMERICAL RECIPES IN C", 2ND EDITION, P.924.
'FORMULA DISCOVERED BY, S. RAMANUJAN, CIRCA, 1920.
'FORMULA IS IN THE FORM OF AN INFINITE SERIES.
'FOR EXACT VALUE OF PI, I = 0,1,2,3,4, ... , TO INFINITY.
'FOR THIS CALCULATION, QUIT WHEN THE NEXT TERM IS SMALLER THAN THE MACHINE EPSILON,
'              "EPS", 1.0E-17 (<-- VALUE, I CALCULATED FOR DOUBLEs, ON MY MACHINE).
'IN OTHER WORDS, WHEN A TERM IS SMALLER THAN "EPS", ITS ADDITION TO A SUM, DOES NOT CHANGE THE SUM.
'SO, FOR THIS CALCULATION, IT TURNS OUT THAT, I = 0,1,2; ONLY THREE TERMS.

DOUBLE PITERM = SQRT(C1) / C2 * FACT(4 * I) * (C3 + C4 * I) / (FACT(I) * C5 ^ I) ^ 4
IF PITERM < EPS
GOTO CP
ELSE
PISUM += PITERM
CALCPI(I + 1)
ENDIF
LABEL CP
RETURN 1.0 / PISUM
ENDSUB

'----------------------------------------------------------------------------------------

SUB FACT(INT F), UINT64

'FACTORIAL FUNCTION

IF F = 0 THEN RETURN 1
IF F = 1
RETURN F
ELSE
FACTPROD = F * FACT(F - 1)
ENDIF
RETURN FACTPROD
ENDSUB

'----------------------------------------------------------------------------------------