# Pi Day - Calculate Pi on the TI

## Recommended Posts

You guys are so much better at math than the person who wrote this scene (doctored with music and footage from a movie)

• 1
• 3

##### Share on other sites
2 hours ago, Vorticon said:

I let the MLC version run for a while, and it seems that the Pi value is getting stuck around the 1.427 value +/- 0.005 or so after over 119 million iterations. I wonder if this is a limitation of the RND algorithm.

The 1 million iterations trial with the XB version is still running, and I'll be curious what kind of result I will get (should be complete within the next 48hrs).

my RND runs from 0 to 4095 but the algorithm considers a radius of 4096. So, more points than expected are included into the 1/4 of circle.

So PI is overestimated.

A change is needed!

Guillaume.

##### Share on other sites

So 1 million iterations in plain XB yielded 3.1434, clearly a very very slow conversion. The MLC version gives about the same answer at 1 million iterations.

OK so now I'm going to do a trial run for 1 billion under Classic99 to see if I can gain another decimal. This will take a while  (just kidding. This will take about 3 years to complete! I need a supercomputer cluster. On the other hand, with overdrive on in Classic 99, I might be able to bring it down to about 3.5 months. That's feasible).

##### Share on other sites
3 hours ago, Vorticon said:

I let the MLC version run for a while, and it seems that the Pi value is getting stuck around the 1.427 value +/- 0.005 or so after over 119 million iterations. I wonder if this is a limitation of the RND algorithm.

The 1 million iterations trial with the XB version is still running, and I'll be curious what kind of result I will get (should be complete within the next 48hrs).

Thanks to your remark, I made two changes in the MLC program:

1) the correct divisor to tranform 4095^2 into 32768 is... 4095^2/32768. This one is sent via CALL LINK to the INITSQ routine

2) the internal routine in MLC that turns a float into an integer does only work from -32768 to +32767. But the last value was +32768 leading to a bad conversion. So the last value of the square table is written directly and not computed.

It looks like the PI value is closer to the one expected now. But, I didn't let the program run up to 119 millions.... Just a few millions..

Here is the new ZIP file with the fixed program.

Guillaume.

##### Share on other sites
5 hours ago, Vorticon said:

I let the MLC version run for a while, and it seems that the Pi value is getting stuck around the 1.427 value +/- 0.005 or so after over 119 million iterations. I wonder if this is a limitation of the RND algorithm.

a program that tests the distribution of one million random numbers:

into a table from test(0) to test(4095) it counts how many times one number is returned

If this was a perfect random generator, we should have:

• table(i)=1.000.000/4096 = 244
• mean = 2047.5
• standard deviation = 1182.4

Here is what I get after millions and millions tries (if we look at the last run):

• all table(i) are included into [178 ; 296]
• all values are reached
• mean is 2048
• standard deviation is 1183

It appears that this generator is not so bad.

Guillaume.

##### Share on other sites
Posted (edited)

Hi again,

I thought that the PI calculation requires two consecutive random numbers to represent X and Y.

What about the differences between two consecutive terms?

**** UPDATE ****

I forgot in a first idea that the differences are not equivalent in probability.

To get 4095 : only one solution 4095-0.

But to get zero : 0-0, or 1-1, or 2-2 etc...

So I change the rest of the post

***************

Picking numbers from 0 to 4095 leads to:

• differences are from -4095 to + 4095
• the perfect mean is zero
• the standard deviation 1672
• the perfect mean of absolute differences should be 1365

Looking at the result, the generator is correct even when thinking about the independency between two consecutive terms.

Guillaume.

Edited by moulinaie

##### Share on other sites

Much better ! But I still don't get the 5 of 3.1415...

##### Share on other sites
3 hours ago, moulinaie said:

Hi again,

I thought that the PI calculation requires two consecutive random numbers to represent X and Y.

What about the differences between two consecutive terms?

**** UPDATE ****

I forgot in a first idea that the differences are not equivalent in probability.

To get 4095 : only one solution 4095-0.

But to get zero : 0-0, or 1-1, or 2-2 etc...

So I change the rest of the post

***************

Picking numbers from 0 to 4095 leads to:

• differences are from -4095 to + 4095
• the perfect mean is zero
• the standard deviation 1672
• the perfect mean of absolute differences should be 1365

Looking at the result, the generator is correct even when thinking about the independency between two consecutive terms.

Guillaume.

Impressive results!

I went ahead and did a quick and dirty evaluation of the XB RND function generating 3000 numbers from 1-10, the maximum array size I could fit in XB. This is the ideal result:

and this is the actual result:

Really not bad at all for only 3000 data points.

##### Share on other sites
1 hour ago, moulinaie said:

Much better ! But I still don't get the 5 of 3.1415...

Excellent! According to Dewdney's article, it took about 4 billion iterations to get the first 4 decimals correct! Keep it running

##### Share on other sites
23 minutes ago, Vorticon said:

Impressive results!

I went ahead and did a quick and dirty evaluation of the XB RND function generating 3000 numbers from 1-10, the maximum array size I could fit in XB. This is the ideal result:

and this is the actual result:

Really not bad at all for only 3000 data points.

Sure ! It's very close to perfection !

##### Share on other sites

Here's a listing of Nilakantha's series, which converges pretty fast:

```100 PI=3
110 N=2
120 PI=PI+4/(N*(N+1)*(N+2))
130 N=N+2
140 PRINT PI
150 PI=PI-4/(N*(N+1)*(N+2))
160 N=N+2
170 PRINT PI
180 GOTO 120```

##### Share on other sites

Now run it in Forth, it should be pi times faster..

##### Share on other sites
6 hours ago, Vorticon said:

Impressive results!

I went ahead and did a quick and dirty evaluation of the XB RND function generating 3000 numbers from 1-10, the maximum array size I could fit in XB. This is the ideal result:

and this is the actual result:

Really not bad at all for only 3000 data points.

Could you do this for RXB or TI BASIC? It'd be interesting to see how they compare

##### Share on other sites

Hi everybody,

I modified the program : now you can stop and save the data on disk to run it again later. On and on... up to 4 billions?

While the computer is calculating, if you press the SPACE BAR (until a new line is displayed):

And when you run it later, you can reload from file:

That's it !

Guillaume.

##### Share on other sites
7 hours ago, senior_falcon said:

Could you do this for RXB or TI BASIC? It'd be interesting to see how they compare

Here you go, with the caveat that the array size is limited to 1500 in TI Basic:

```10 CALL CLEAR
20 RANDOMIZE
30 OPTION BASE 1
40 DIM N(1500)
50 PRINT "PHASE 1:"
60 FOR I=1 TO 1500
70 N(I)=(RND*10)+1
80 T=T+N(I)
90 NEXT I
100 PRINT "PHASE 2:"
110 MEAN=T/1500
120 FOR I=1 TO 1500
130 T1=T1+(N(I)-MEAN)^2
140 PRINT I
150 NEXT I
160 SD=SQR(T1/1500)
170 PRINT "MEAN:";MEAN:"SD:";SD```

Using the same array size in XB gives:

XB is still better, but not my much, although I do expect the gap to widen with a larger number of iterations.

##### Share on other sites
18 minutes ago, Vorticon said:

Here you go, with the caveat that the array size is limited to 1500 in TI Basic:

Hi,

If you change the organization of your DATA, there is no, limitation.

You can use a table, let's say form T(0) to T(9) that count how many times a 0, or a 1 , or ... a 9 has appeared with RND.

This way, you just have to reserve memory for ten numbers.

Then, there is a simplified formula for standard deviation :

SD = SQR ( [sum of squares] / N - [square of mean] )

```100 CALL CLEAR
110 DIM T(9)
120 RANDOMIZE
125 INPUT "How many ":N
130 FOR I=1 TO N
140 X=INT(RND*10)
150 T(X)=T(X)+1
160 NEXT I
170 S=0
180 S2=0
190 FOR I=0 TO 9
200 S=S+T(I)*I
210 S2=S2+T(I)*I*I
220 NEXT I
230 S=S/N
240 S2=S2/N
250 SD=SQR(S2-S*S)
260 PRINT "Mean=";S
270 PRINT "S.Dev=";SD```

I prefer 0 to 9 than 1 to 10, so I skip the "+1" at each random!

The expected values are:

MEAN = 4.5

SD = 2.87

Works both in TI Basic and XB.

Guillaume.

##### Share on other sites
33 minutes ago, moulinaie said:

Hi,

If you change the organization of your DATA, there is no, limitation.

You can use a table, let's say form T(0) to T(9) that count how many times a 0, or a 1 , or ... a 9 has appeared with RND.

This way, you just have to reserve memory for ten numbers.

Then, there is a simplified formula for standard deviation :

SD = SQR ( [sum of squares] / N - [square of mean] )

```100 CALL CLEAR
110 DIM T(9)
120 RANDOMIZE
125 INPUT "How many ":N
130 FOR I=1 TO N
140 X=INT(RND*10)
150 T(X)=T(X)+1
160 NEXT I
170 S=0
180 S2=0
190 FOR I=0 TO 9
200 S=S+T(I)*I
210 S2=S2+T(I)*I*I
220 NEXT I
230 S=S/N
240 S2=S2/N
250 SD=SQR(S2-S*S)
260 PRINT "Mean=";S
270 PRINT "S.Dev=";SD```

I prefer 0 to 9 than 1 to 10, so I skip the "+1" at each random!

The expected values are:

MEAN = 4.5

SD = 2.87

Works both in TI Basic and XB.

Guillaume.

Yes of course! Did not think of it that way. Well done.

##### Share on other sites
2 hours ago, Vorticon said:

Here you go, with the caveat that the array size is limited to 1500 in TI Basic:

```10 CALL CLEAR
20 RANDOMIZE
30 OPTION BASE 1
40 DIM N(1500)
50 PRINT "PHASE 1:"
60 FOR I=1 TO 1500
70 N(I)=(RND*10)+1
80 T=T+N(I)
90 NEXT I
100 PRINT "PHASE 2:"
110 MEAN=T/1500
120 FOR I=1 TO 1500
130 T1=T1+(N(I)-MEAN)^2
140 PRINT I
150 NEXT I
160 SD=SQR(T1/1500)
170 PRINT "MEAN:";MEAN:"SD:";SD```

Using the same array size in XB gives:

XB is still better, but not my much, although I do expect the gap to widen with a larger number of iterations.

RXB got:

MEAN: 5.923564439

SD:     2.929720887

That seems odd as second run says:

MEAN: 5.965401107

SD:     2.921225948

##### Share on other sites
38 minutes ago, RXB said:

RXB got:

MEAN: 5.923564439

SD:     2.929720887

That seems odd as second run says:

MEAN: 5.965401107

SD:     2.921225948

Try using Guillaume's version with a larger n, perhaps 5000, and see what you get.

##### Share on other sites

Here is the RXB versus XB RND function comparison using Guillaume's program for an N of 5000. Again expected values are mean= 4.5 and SD=2.87

RXB:

XB:

Not a huge difference for all practical purposes unless a very large set of random numbers are needed such as in the PI approximation program in which case the XB RND will be superior.

##### Share on other sites
23 minutes ago, Vorticon said:

Here is the RXB versus XB RND function comparison using Guillaume's program for an N of 5000. Again expected values are mean= 4.5 and SD=2.87

RXB:

XB:

Not a huge difference for all practical purposes unless a very large set of random numbers are needed such as in the PI approximation program in which case the XB RND will be superior.

Yea well I swapped out the XB RND for TI Basic RND for the speed difference not accuracy. (Bonus was a lot of more GRAM space.)

##### Share on other sites

Some 4.2 billion iterations and we have the first 3 decimals of PI=3.141 using Guillaume's program! It took about 3 days to get this far under Classic99 in Overdrive mode and TI XB.

I'm letting it run for a while longer and see if we can score a 4th decimal

##### Share on other sites

I have been working on a Forth version, which I could post now, but I want to be able to discuss it, and that will not happen until I get home on Thursday—too much going on. A day or so after that, I should probably be able to get back to it. We shall see.

Just a teaser: I first used floating point math in the loop, but that approached XB’s pace. Then I worked out a way with 32-bit integers, which is very much faster. I was as disappointed as you, Walid, to discover how long it was likely to take for a reasonable value.

TTFN...

..lee

• 3
• 1

##### Share on other sites

Well it seems that the random shot method has its limits. After 14.3 billion iterations, PI has converged to 3.141076. I guess we're not going to get that 4th decimal after all. Interesting...

##### Share on other sites

The standard calculation of pi produces a number that extends to infinity but we don't really need that.  Back in 1706 a good approximation calculation was produced by John Machin and we can program that to produce quite a few accurate decimal places fairly quickly in Extended Basic.

I wrote the program below back in March 1991.  It uses the Extended Basic "logical operator" AND. You can see what this does with FOR T=1 TO 20 :: PRINT T;T AND 2 :: NEXT T

However- calculating PI:
1 ! from W F DOSSETT,
AUSTIN, TEXAS,
25 feb 88
2 ! REC NEWSLETTER March 1988 Vol 3 #2
3 ! for ti99/4a s shaw 3/91
4 ! from John Machin, 1706
7 ! n=7 for 13 digit accuracy of ti99/4a
8 ! n>7 if language goes beyond 13 digits!
9 ! eg n=35 for 53 digits
10 !
100 N=7 :: K=-1
110 K=K+1
120 I=-((2 AND 2*K)-1)
130 U=U+4*I/((2*K+1)*5^(2*K+1))
140 V=V+I/((2*K+1)*239^(2*K+1))
150 IF K<=N THEN 110
160 W=U-V :: PI2=4*W
170 ! DONE
180 PRINT PI2,PI
181 ! calculated number then onboard number
190 PRINT (PI2-PI)*1000,! COMPARE with on board variable
191 ! difference times 1000 makes it easier to see
200 PRINT (PI2-(4*ATN(1)))*1000 ! compare with TI Basic's calculation
210 PRINT " ";(PI2-3.14)*100;:: DISPLAY AT(24,2)SIZE(4):"3.14" ! displays ALL the digits actually in the variable!!!=13 incl dec point.
220 PRINT " 3.14159265358979 -etc" ! more accurate actual value for comparison!
230 END
see https://en.wikipedia.org/wiki/John_Machin

"The benefit of the new formula, a variation on the Gregory/Leibniz series (Pi/4 = arctan 1), was that it had a significantly increased rate of convergence, which made it a much more practical method of calculation."

s

Edited by blackbox
references

## Join the conversation

You can post now and register later. If you have an account, sign in now to post with your account.

×   Pasted as rich text.   Paste as plain text instead

Only 75 emoji are allowed.