Search code examples
rseriespiapproximation

Series vector for approximating pi


I've been set a question about Madhava's approximation of pi. The first part of it is to create a vector which contains the first 20 terms in the series. I know I could just input the first 20 terms into a vector, however that seems like a really long winded way of doing things. I was wondering if there is an easier way to create the vector?

Currently I have the vector

g = c((-3)^(-0)/(2*0+1), (-3)^(-1)/(2*1+1), (-3)^(-2)/(2*2+1), (-3)^(-3)/(2*3+1), (-3)^(-4)/(2*4+1), (-3)^(-5)/(2*5+1), (-3)^(-6)/(2*6+1), (-3)^(-7)/(2*7+1), (-3)^(-8)/(2*8+1), (-3)^(-9)/(2*9+1), (-3)^(-10)/(2*10+1), (-3)^(-11)/(2*11+1), (-3)^(-12)/(2*12+1), (-3)^(-13)/(2*13+1), (-3)^(-14)/(2*14+1), (-3)^(-15)/(2*15+1), (-3)^(-16)/(2*16+1), (-3)^(-17)/(2*17+1), (-3)^(-18)/(2*18+1), (-3)^(-19)/(2*19+1), (-3)^(-20)/(2*20+1))

And

  h = sqrt(12)

So I have done g*h to get the approximation of pi. Surely there's an easier way of doing this? Apologies if this is relatively basic, I am very new to R and still learning how to properly use stack overflow. Thanks.


Solution

  • One of the best features of R is that it is vectorised. This means that we can do operations element-wise on entire vectors rather than having to type out the operation for each element. For example, if you wanted to find the square of the first five natural numbers (starting at one), we can do this:

    (1:5)^2
    

    which results in the output

    [1]  1  4  9 16 25
    

    instead of having to do this:

    c(1^2, 2^2, 3^2, 4^2, 5^2)
    

    which gives the same output.

    Applying this amazing property of R to your situation, instead of having to manually construct the whole vector, we can just do this:

    series <- sqrt(12) * c(1, -1) / 3^(0:19) / seq(from=1, by=2, length.out=20)
    sum(series)
    

    which gives the following output:

    [1] 3.141593
    

    and we can see more decimal places by doing this:

    sprintf("%0.20f", sum(series))
    [1] "3.14159265357140338182"
    

    To explain a little further what I did in that line of code to generate the series:

    • We want to multiply the entire thing by the square root of 12, hence the sqrt(12), which will be applied to every element of the resulting vector
    • We need the signs of the series to alternate, which is accomplished via * c(1, -1); this is because of recycling, where R recycles elements of vectors when doing vector operations. It will multiply the first element by one, the second element by -1, then recycle and multiply the third element by 1, the fourth by -1, etc.
    • We need to divide each element by 1, 3, 9, etc., which is accomplished by / 3^(0:19) which gives / c(3^0, 3^1, ...)
    • Lastly, we also need to divide by 1, 3, 5, 7, etc. which is accomplished by seq(from=1, by=2, length.out=20) (see help(seq))