[OT] 355/113 and pi



#2

A friend pointed me at this article about 355/113 and it's remarkable closeness to pi. I thought some of you might find it interesting.

Dave


#3

That was a very interesting read! Instead of memorizing a fraction that is "close enough" to pi, I just went and memorized 3.1415926535897932384626433 and yes that is all from memory... :)


#4

Quote:
That was a very interesting read! Instead of memorizing a fraction that is "close enough" to pi, I just went and memorized 3.1415926535897932384626433 and yes that is all from memory... :)

...83279 (from memory as well - which shows you what happens to bored teenagers in high school English classes. At least this one. Couldn't diagram a sentence to save my life, though.)

But back OT: as you noted below the approximation 355/113 gives more than six digits accuracy and is arguably easier to remember:

- First three odd integers, doubled: 113355

- Split in the middle and divide

One of my math professors had worked for NASA and brought this topic up in a numerical analysis class. He said they were computing trajectories for the Apollo 8 flight using two programs but kept getting different results. The problem turned out to be that "some &#$! engineer" had put 22/7 in for pi in one of them. "At 240,000 miles it makes a difference." I've always remembered that!


Edited: 26 July 2012, 12:20 a.m.

#5

Quote:
That was a very interesting read! Instead of memorizing a fraction that is "close enough" to pi, I just went and memorized 3.1415926535897932384626433 and yes that is all from memory... :)

I teach at a K-12 school. A couple of months ago at the elementary school talent show, a kid got up and recited pi to 50 digits. It was both impressive and entertaining.

-wes

#6

There is an answer to this issue (that that article doesn't have).

The continued fraction expansion for pi is

3; 7, 15, 1, 292...

which works out to

3, 22/7, 333/106, 355/113, 103993/33102,...

There's no mystery where 355/113 comes from; it's one of the convergents of the continued fraction.

http://en.wikipedia.org/wiki/Convergent_(continued_fraction)


#7

Excellent catch. In case you don't wanna bother with wiki, I copied the pertinent info here

Every real number can be expressed as a regular continued fraction in canonical form. Each convergent of that continued fraction is in a sense the best possible rational approximation to that real number, for a given number of digits. Such a convergent is usually about as accurate as a finite decimal expansion having as many digits as the total number of digits in the nth numerator and nth denominator. For example, the third convergent 333/106 for π (Pi) is roughly 3.1415094, which is not quite as accurate as the 6-digit 3.14159; the fourth convergent 355/113 = 3.14159292 is more accurate than the 6-digit decimal.


Edited: 25 July 2012, 6:49 p.m.


#8

I wrote some FORTRAN programs for an IBM 360 mainframe in the 1970s that used 355/113 for pi, since it was accurate to within the machine accuracy. Some of these programs are still in use today on PCs, and 355/113 is a sufficiently accurate representation of pi for the purposes of the computations.


#9

Yep, I never needed more accuracy either actually. I learned to code in FORTRAN on a DEC VAX 11/780 in the early eighties myself.

#10

There was a post here a while ago (thus, it's in the archives somewhere) about pi to hundreds of thousands or millions of decimal places. Within the thread, someone mentioned a website which contains a searchable and browsable representation of pi. The image is a colour-coded representation, with each digit represented by a different coloured pixel. You could also move your mouse around this picture and the hypertext would display the digits of pi where you were browsing.

#11

pi and other rational-number approximations I've found for use in scaled-integer math: http://wilsonminesco.com/16bitMathTables/RationalApprox.html


#12

In 1976 I found this while toying around with my beloved HP25:
Pi~(401699/133)^(1/7).
HP fans do strange things, don't they ;)


#13

That's pretty close to 6317674554 / 2010978267 = 3.14159265551....

And a resounding "NO! We're not strange...the rest of the world is!"

Edited: 26 July 2012, 2:36 p.m.

#14

Quote:
HP fans do strange things, don't they ;)

I would rephrase this as "Pi fans do strange things, don't they?" :-)

The other day I presented the following 12-digit approximation involving ln(2) and e:

Jim Markovitch, whom I also showed it to, not only improved it to 16 decimal places but rewrote it in a very interesting way:

#15

Interesting...

The Bailey-Borwein-Plouffe (BBP) formula is also interesting :

On the 50G

« 0 DUP ROT FOR k
`(1/(16^k))*((((4/((8*k)+1))-(2/((8*k)+4)))-(1/((8*k)+5)))-(1/((8*k)+6)))` +
NEXT
»
'~PI' STO

5 ~PI

'47/15+53/6552+829/5026560+79/15590400+857/4561108992+1901/244729774080'

->NUM 3.14159265323

I computed 2000 digits of PI with this formula on my 50G

The new HP39gII give for PI with the embedded a b/c function :
1146408/364913

Same result with ->Q on the 50G wich convert any real to a close rationnal number


Edited: 26 July 2012, 7:51 a.m.


#16

I've read the BBP formula is used to compute individual hexadecimal digits of pi without calculation all the previous ones. There is a Mathcad implementation of the algorithm here, but I have no idea on how to port it to RPL.

Quote:
I computed 2000 digits of PI with this formula on my 50G

How long did that take? Not too long, I presume.

One year ago I computed pi to one thousand places on the HP 50g in slightly less than three hours, using a modified Archimedes algorithm and a third-party extended precision libray. That's a lot of time, but the original method would have required five times as much!

http://www.hpmuseum.org/cgi-sys/cgiwrap/hpmuseum/archv020.cgi?read=188443

Regards,

Gerson.


#17

Quote:
One year ago I computed pi to one thousand places on the HP 50g in slightly less than three hours, using a modified Archimedes algorithm and a third-party extended precision libray. That's a lot of time, but the original method would have required five times as much!

My 50g pies: http://sense.net/~egan/hpgcc/#Example:%20%20π%20Shootout

Fastest uses FFTs. 8109 digits in 36 seconds. Of course I used C.


#18

Thanks for this !

It's very very interesting. Not only for PI but for yours explanations about "C" and HP50G !

I have the idea than (in theory !) it's possible to calculte more digits in RPL than in C but it will take more more time (probably too much to be really usable)!! My 5000 RPL PI digits is running ... Wait and see ;)

Edited: 26 July 2012, 5:28 p.m.


#19

Quote:
Thanks for this !

It's very very interesting. Not only for PI but for yours explanations about "C" and HP50G !

I have the idea than (in theory !) it's possible to calculte more digits in RPL than in C but it will take more more time (probably too much to be really usable)!! My 5000 RPL PI digits is running ... Wait and see ;)


Not with a single set of batteries. :-) Connect up your USB port.

If you read my page you'll see that I did 200,000 digits in ~17.5 hrs. In C @ 75MHz.


#20

Yes I saw that : really impressive on the calc at 75MHz!

I've got my 5000 digits this night (but with emu48). Less than 8 hours but I don't know exactly (I runned with EVAL instead of TEVAL)

I've an idea for a much smaller (!) and perhaps quicker program (User RPL)

Gilles


Edited: 27 July 2012, 2:37 a.m.


#21

Lighter and faster USR RPL :

Here for 2000 digits :

«
0 @ initial Value of PI
2000 ALOG @ Number of digits to find
1 @ 16^k initial value = 16^0

0 1701 FOR k @ Number of iteration = 2+log(8n)+n*log(16)
k 1. DISP @ Display iterations
OVER [ 120 151 47 ] k PEVAL * @ Nominator * 2000 ALOG
OVER 16 * SWAP ROT @ 16^k and 16^(k+1) on the stack
[ 512 1024 712 194 15 ] k PEVAL * @ Denominator
IQUOT 4. ROLL + UNROT @ Integer quotient and add to ~PI
NEXT

DROP2
»

181 Bytes and 371 sec with EMU48 for 2000 digits on my old laptop

It use Plouffe formula

'(1/(16^k))*((((4/((8*k)+1))-(2/((8*k)+4)))-(1/((8*k)+5)))-(1/((8*k)+6)))'

EVAL EXPAND

->

'(k^2*120+151*k+47)/((512*k^4+1024*k^3+712*k^2+194*k+15)*2^(4*k))'

EDIT 1 : I just realise that 1701 is less than 2+log(8n)+n*log(16). I found this formula on the web after some tests for the value of the iterations. With 1701 the 4 or 5 lasts digits are false.

EDIT 2 : 5000 digit in 7500 sec (~ 2 hours) with EMU48. The last four digits are wrong again- curious

EDIT 3 : 10.000 in progress ;) But 200.000 will takes years or centuries ! :O :D IQUOT use lot of time with so big numbers !
note that 10000 ALOG -> 'integer too large', but 5000 ALOG DUP * seems OK. The program use very few memory (few KBytes)


Edited: 27 July 2012, 7:29 a.m.

#22

Quote:
Fastest uses FFTs. 8109 digits in 36 seconds. Of course I used C.

That's really very impressive! Of course I wasn't trying to break any record with a 2,500 year old algorithm :-)

Also very interesting is your Pi Day Rematch: Apple II vs. HP-41C last year. Were I brave enough, I would make my 8-bit MSX computer join the fight :-)

#23

Hi Gerson,

here is the program for 2000 decimals.

«
2000 ALOG
0 DUP 1701 FOR k
k 1. DISP
k DUPDUP 151 * SWAP SQ DUP 120 * ROT + 47 + 5. PICK *
UNROT
DUP SQ 512 * UNROT DUP2 * 1024 * UNROT 712 * SWAP 194 * 15 + + + + 16 k ^ *
IQUOT +
NEXT
»

It's only 207 Bytes. It requires no extension precision library. It only use the 'infinite' integer possibility of the 50G. In this case, it works with a 2000 digits integer number.

The idea is :

10^2000 * * 10^2000

I evaluate the Plouffe formula with the 50G and transform it in full RPN (but it is not obvious it's faster that the algebraic, I dont try, but there is sometimes surprises with this !)

Must run in exact mode, uncheck LASTARG

With full speed emu48, it takes 18 minutes. Much more on a real 50G (many hours ?) (the emulator is much faster that the real calc with long integer, more than for other calulations...)

The number of iterations (here 1701 for 2000 décimals) is calculated with :
2+log(8n)+n*log(16)

For example, if you want 200 decimals to test, just change :

2000 ALOG -> 200 ALOG
1701 -> 247 (numbers of iterations)

So it possible to modify the program to ask the number or decimal you want... I tried with 5000 but it seems that a 5000 digits integer is too big for the 50G (Just try 5000 ALOG, you will get an error) ! It is very 'economic' in memory because 1 digit is only a nibble ( 1 octet for 2 digits !)

PS : I deleted some library on my EMU48 and try again with 5000 digits.... My computer will work this night ;)

PS : Gerson in your link, there is a PASCAL version. Seems very easy to use it for the 39Gii (copy and paste and few change !)

What does the :nn here

WriteLn(n:4,a:22:17,b:22:17);

?


Edited: 26 July 2012, 5:22 p.m.


#24

Thanks for your explanations!

Quote:
What does the :nn here

WriteLn(n:4,a:22:17,b:22:17);

?


That's just number formatting:

n:4 -> 4 printing positions for integer variable n

a:22:17 -> 22 positions for real variable a, 17 of which reserved for the fractional part.

----------------------------------------

P.S.: Here is a C++ version, for clarity:

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

int main(int argc, char *argv[])
{
char i, m;
int n = 6;
double a, b, c, d, p;
scanf("%d", &m);
system("cls");
b = 2*sqrt(3);
a = 3;
printf("\ n-gon\t a\t b\t\t p\n\n");
for (i = 1; i < m; i++)
{
printf("%12d %.16f %.16f\n",n,a,b);
b *= 2*a/(a+b);
a = sqrt(a*b);
n += n;
}
d = 2*a*b/(a+b);
c = sqrt(a*d);
p = (640*c*(28*c-139*a)/(c-4*a)+4*b*(3737+14256*b/(d-4*b))+1056*d-2568*a-6453*pow(a*a*b,1/3.))/11655;
printf("%12d %.16f %.16f %.16f\n",n,a,b,p);
printf("\n");
system("pause");
return 0;
}

n-gon a b p

6 3.0000000000000000 3.4641016151377544
12 3.1058285412302493 3.2153903091734723
24 3.1326286132812382 3.1596599420975005
48 3.1393502030468672 3.1460862151314348
96 3.1410319508905093 3.1427145996453683 3.1415926535897927

Pressione qualquer tecla para continuar. . .


Edited: 26 July 2012, 9:12 p.m.

#25

Hi

here is a 50G UsrRPL program to calculate the Nth digits of PI in hexadecimal. It is bases on the french explanation of wikipedia

« -> n a « 0 n `n+15` FOR k `16^(n-k)/(8*k+a)` + NEXT » ->NUM » 'B' STO

« -> n a « 0 0 `n-1` FOR k k 1. DISP MEM 2. DISP `(8*k+a)` DUP MODSTO 16 n k - POWMOD SWAP / + NEXT » ->NUM » 'A' STO

«  -> n
'4*(A(n,1)+B(n,1))-2*(A(n,4)+B(n,4))-(A(n,5)+B(n,5))-(A(n,6)+B(n,6))'
DUP SIGN * FP 16 4 ^ * HEX R->B
» 'PIn' STO

Just type

n PIn to get the hexa digits for n-1 to n+3

For example 200 PIn


Edited: 28 July 2012, 6:13 p.m.

#26

How did you get the extra 1900+ digits when the 50G only calculates to less than 16 digits of precision?


#27

50G is only limited by memory for integers.


#28

Still not clear. Are you saying that the 1900+ extra digits were calculated by one of the above algorithms one at a time?


#29

Hi Matt,

see my previous post and new algorithm.

The idea is (for example ) :

1/3 ~= 0.3333333333

but we want to work in integer because the 50G can work with very big numbers in integer (in fact I dont know the limite....)

Suppose we want 5 decimals for 1/3 but we want to work in integer only :

100000*(1/3) = 100000/3 ~= 100000*0.333333 -> 33333(.333)

With the 50G :

(100000 * 1) IQUOT 3 @ IQUOT = Integer quotient on the 50G

-> 33333 (the 5 first decimals )

<< 1000 ALOG 3 IQUOT >> -> 3333...3333 (one hundred'3')

And the HP50 can handle numbers like 5000 ALOG (1 followed by 5000 zero , 10^5000, which is much more thans the number a atoms in the Universe 10^80 (as far we know !)

PS :I cannot test now, but

'\GS(k=0,50,IQUOT(ALOG(64)*PEVAL([ 120 151 47 ],k), PEVAL([ 512 1024 712 194 15 ],k)*16^k)))'

could give then 64 first decimal of PI. Not sure it will work


Edited: 27 July 2012, 5:37 a.m.

#30

Quote:
Still not clear. Are you saying that the 1900+ extra digits were calculated by one of the above algorithms one at a time?

Yes. Arbitrary precision math often uses large integers or better stated, an array of integers.

Possibly Related Threads...
Thread Author Replies Views Last Post
  [OT] Mathematica free for Raspberry PI BruceH 32 2,283 11-23-2013, 05:24 AM
Last Post: Nick_S
  Computing pi with the PC-1300S Kiyoshi Akima 0 226 11-17-2013, 12:24 AM
Last Post: Kiyoshi Akima
  Calculating Pi LHH 9 709 09-27-2013, 10:50 PM
Last Post: Gerson W. Barbosa
  Visualization of pi Bruce Bergman 13 953 08-17-2013, 05:00 PM
Last Post: Howard Owen
  OT: Happy Pi Day! Eddie W. Shore 13 981 03-22-2013, 10:44 AM
Last Post: Les Koller
  Totally OT ... Pi Day for my car Maximilian Hohmann 18 1,206 03-10-2013, 01:15 PM
Last Post: chris smith
  [WP34S] A funny bug in Pi (prod) Eduardo Duenez 3 442 01-28-2013, 03:41 AM
Last Post: Walter B
  28S Pi Functionality Matt Fegenbush 3 401 10-17-2012, 02:15 AM
Last Post: Nick_S
  e^pi - pi + 9/10^4 + 1/(10^4*ln(2) + sqrt(10)/6)^2 Gerson W. Barbosa 47 3,195 08-08-2012, 10:58 PM
Last Post: Les Koller
  Calculating Pi Andrew Davie 6 496 07-06-2012, 03:30 AM
Last Post: Gjermund Skailand

Forum Jump: