nth digit of PI
BlitzMax Forums/BlitzMax Beginners Area/nth digit of PI
| ||
Hi, I'm trying to program the formula to find the Nth digit of PI. Here is one version of the formula, but I'm doing something wrong and can't make it work... Quoted from http://www.math.hmc.edu/funfacts/ffiles/20010.5.shtml Here is a very interesting formula for pi, discovered by David Bailey, Peter Borwein, and Simon Plouffe in 1995: Pi = SUMk=0 to infinity 16-k [ 4/(8k+1) - 2/(8k+4) - 1/(8k+5) - 1/(8k+6) ]. ... Now Wolfram says something different: http://mathworld.wolfram.com/PiDigits.html Xn = frac( 16Xn-1 + 120n^2 - 89n +16 / 512n^4 - 1024n^3 + 712n^2 - 206n - 21) I can't get either of these to work. It's got to be my interpretation of the formulas. My non-working attempt... '4/(8k+1) - 2/(8k+4) - 1/(8k+5) - 1/(8k+6) Local k:Double Local j:Double Local x:Double For k=0 To 4 j=16^-k*( (4/(8*k+1)) - (2/(8*k+4)) - (1/(8*k+5)) - (1/(8*k+6))) Print j Next 'Xn = frac( 16Xn-1 + (120n^2 - 89n +16 / 512n^4 - 1024n^3 + 712n^2 - 206n - 21)) For k=0 To 4 x= 16*(k-1) + (((120*n^2) - (89*n +16)) / ((512*n^4) - (1024*n^3) + (712*n^2) - (206*n) - 21)) Print x Next |
| ||
you forgot to sum the first one:'4/(8k+1) - 2/(8k+4) - 1/(8k+5) - 1/(8k+6) Local k:Double Local j:Double Local x:Double For k=0 To 24 j:+16^-k*( (4/(8*k+1)) - (2/(8*k+4)) - (1/(8*k+5)) - (1/(8*k+6))) Print j Next |
| ||
Thanks skidracer, but I guess I'm still confused. If I want to find the 29th digit of PI, I should be able to make k=29, send it to the formula and have 7 as the answer. Do I need to change to a different print routine to make it work correctly? The whole point of this formula is to make it so you don't have to calculate all the previous digits. I just stuck it in a loop to see what answers I'd get. Maybe I'm just not smart enough to figure this out.... ho hummmm. :( |
| ||
There's a lot of detail missing from that description. And it ultimately gives base 16 digits. Pi in base 10 is 3.14159..., which means 3 + 1/10 + 4/100 + 1/1000 etc. Pi in base 16 is 3.243F6..., which is 3 + 2/16 + 4/256 + 3/4096 + 15/65536 etc. First, here's a demonstration that the infinite sum gives the value of Pi. Then we pick out the first few hexadecimal digits. Local k:Double Local j:Double Local x:Double Local Pie:Double For k=0 To 15 ' not quite infinity, but close enough. Pie :+ 16^-k*( (4/(8*k+1)) - (2/(8*k+4)) - (1/(8*k+5)) - (1/(8*k+6))) Next Print Print " 0 to 15 sum: " + Pie Print " Value of Pi: " + Pi Print ' Now let's see some hexadecimal digits of Pi. Local n:Int, d:Int, t:String For n = 0 To 10 d = Int(Pie) Pie :- d Pie :* 16 t :+ "0123456789ABCDEF"[d..d+1] Next Print "Pi in base 16: " + t[0..1] + "." + t[1..] + "..." The point of all this is that the initial sum to compute Pi contains powers of 1/16. At first it looks like you need to do the entire sum. That's what I do in the sample code above. But with careful estimates you can figure out what a particular digit is with much less computation. As I said, there are a lot of details to fill in. Also, you really can't get very far with double precision floats. They give 53-bit accuracy. Calculating Pi out to the 29th hexadecimal digit needs far more as each digit is four bits. Allowing for some extra precision to ensure rounding to the proper digit you would probably need about 130-bit accuracy. Last edited 2010 |
| ||
This is interesting. It looks like that MathWorld formula has a typo. I implemented it and the resulting digits looked nothing like reality. Notice that the terms in the numerator go + - +. The denominator goes + - + - -. We can reasonably speculate that should be + - + - +. Now the result agrees with the code I posted earlier: 3.243F6A8885 Local x:Double, n:Double, d:Int For n = 1 To 10 x = 16*x + (120*n^2 - 89*n +16) / (512*n^4 - 1024*n^3 + 712*n^2 - 206*n + 21) x :- Floor(x) ' x = Frac(x) d = Int(16*x) ' nth hexadecimal digit after 3. Print d NextNote this can't be continued for long. Every time through the loop we multiply by 16 and discard the integer part. This throws away 4 bits. So with 53-bit arithmetic we can get twelve digits, maybe thirteen. After that the results are certainly junk. |