S Lazy-H
  • Home
  • About
  • Posts
  • Contact
  • Slide Rules
  • A Biker’s Tale

Fun With PI

mathematics
miscellaneous
Author

Sam Hutchins

Published

March 9, 2022

In another post, I posted some observations about Circle Squaring, where people were attempting to square a circle using only a straightedge and a compass. In this post I am concentrating on an “R” program for generating decimal digits of \(\pi\). First, a disclaimer: I am a casual dabbler in R programming.

One of the formulas I showed was a Bailey-Borwein-Plouffe (BBP) Formula:

\[\pi _{n} = \biggl[ \frac{2(-1)^{n+1}(2 n)!}{2^{2 n}B_{2 n}(1-2^{-n)}(1-3^{-n})(1-5^{-n})(1-7^{-n})} \biggr] ^{1/(2n)}\]

and Ramanujan’s 1914 formula:

\[\frac{1}{\pi} = \frac{2 \sqrt{2}}{99^{2}} \sum_{k=0}^{\infty} \frac{(4 k)!}{k!^{4}} \frac{26390 k + 1103}{396^{4k}}\]

I like to work with “R” programming, as it can do most things any other programming language can do, even it was originally designed for statistics and data manipulation. In fact, all the posts on this site are done using rmarkdown(J. Allaire et al. 2024) or bookdown(Xie 2025), and now Quarto(J. J. Allaire et al. 2022) in one form or fashion.

So, I worked up a little “R” program to calculate the first 15 digits of \(\pi\), just for the fun of it. Understand, however, that “R” has \(\pi\) built-in, as do most programming languages, such as Python, Perl, C, C++ and so forth. The formula I am using is the Ramanujan formula shown above.

Code
fac <- function(x) {
    if (x==0) { 
        return(1)
    }
    else{
        r <- x * factorial(x-1)  
        return(r)
    }
}

The above code snippet is a function to execute the factorials in the formula. The next part of the formula in the ramanujan function is the portion \((\frac{2 \sqrt{2}}{99^{2}})\). It is broken out as it is a static value where i is used as the variable, and could be shown as \((\frac{(\sqrt{8})}{9801})\), but is a bit more obscure. The while loop below continues until the femto digit is reached (1e15).

Code
ramanujan <- function() {
    i <- (2*sqrt(2))/99^2
    n <- 0
    tot <- 0
    tmp <- 0
    while (T) {
        tmp <- i * (fac(4*n)/fac(n)^4) * ((26390*n+1103)/396^(4*n))
        tot <- tot+tmp
        
        # stop at digit 15
        if (abs(tmp) <= 1e-15) {
            break
        }
        n <- n+1
    }
    return(1/tot)
}

Increasing the digits calculated serves no purpose as the calculated digits are incorrect after the 15th digit. For example, the correct digits are 3.14159265358979323846, but the above gives 3.14159265358979311599 if extended. Setting the options(digits=22), and changing the if (abs() line to 1e-22, and changing the below line’s format to “%.22f will show the extra incorrect digits.

Code
sprintf(&quot;Calculated PI: %.15f&quot;, ramanujan())

I think this is caused by bit limitations of the underlying operating system (64-bit). However, 64 bits still give a pretty large unsigned number: 264 = 18,446,744,073,709,551,616. There are several programs in the wild that can calculate \(\pi\) to thousands, millions, and more digits accurately. So why bother with this little “R” program? Well, just for the fun of it, mostly. Again, as mentioned above, “R” has its own built-in \(\pi\) value to show the same number of accurate digits. I suspect Ramanujan’s formula can be adjusted to give more accurate results, but I am not sure of the source of some of the constants used in his formula, so cannot do that adjustment myself. I was also getting errors in the if(abs()) line because the formula was producing NaN after so many iterations. I discovered this by printing the intermediate results. To catch the error, I used the below line.

Code
if ( (abs(tmp) <= 1e-15) || tot =='NaN';) { break }

Perhaps more reserch is needed?

I also tried Chudonovsky’s 1987 formula, with no greater success, as it still gives me incorrect digits after the 15th.

\[\frac{1}{\pi} = 12 \sum_{k=0}^{\infty} \frac{(-1)^n (6n)!}{(3n)! (n!)^3} \frac{13591409 + 545140134n}{(640320^3)^{n+\frac{1}{2}}}\]

For example, directly setting a variable to equal the 20 digits noted above (3.14159265358979323846,), and then printing that variable, you will notice it has reverted to the inaccurate value (3.14159265358979311599). As noted above, I suspect the fault lies in the operating system’s available bits being used (64), so I may see what I can do with the packages Brobdingnag and Rmpfr, and setting an option in the program as so: options(scipen=999). More fun with formulae!

In any case, I enjoyed it, and it gave me some practice in “R” that I needed; and also, more insight in structuring math formulas using MathJax. Have a great day, and remember, God is in charge, and can make even bad things work for His good. So Praise God, stay safe, and God Bless you and yours!

References

Allaire, J. J., Charles Teague, Carlos Scheidegger, Yihui Xie, and Christophe Dervieux. 2022. “Quarto.” https://doi.org/10.5281/zenodo.5960048.
Allaire, JJ, Yihui Xie, Christophe Dervieux, Jonathan McPherson, Javier Luraschi, Kevin Ushey, Aron Atkins, et al. 2024. Rmarkdown: Dynamic Documents for r. https://github.com/rstudio/rmarkdown.
Xie, Yihui. 2025. Bookdown: Authoring Books and Technical Documents with r Markdown. https://github.com/rstudio/bookdown.
© S Lazy-H 2019 -