Author Topic: Estimating Pi  (Read 154 times)

FellippeHeitor

• QB64 Developer
• Forum Resident
• Posts: 2209
• LET IT = BE
Estimating Pi
« on: August 01, 2020, 09:11:11 PM »
Of course we already have _PI, this is just a cool way to visualize it.

Code: QB64: [Select]
1. SCREEN _NEWIMAGE(600, 600, 32)
2.
3. DIM inCircle AS _UNSIGNED LONG, total AS _UNSIGNED LONG
4.
5.     px = RND * _WIDTH
6.     py = RND * _HEIGHT
7.
8.     total = total + 1
9.     IF dist(_WIDTH / 2, _HEIGHT / 2, px, py) <= 300 THEN
10.         inCircle = inCircle + 1
11.         PSET (px, py), _RGB32(100)
12.         PSET (px, py), _RGB32(70)
13.     pi = (inCircle * 4) / total
14.
15.     LOCATE 1, 1
16.     PRINT pi
17.
18.     CIRCLE (_WIDTH / 2, _HEIGHT / 2), 299, _RGB32(216, 144, 22)
19.
20. FUNCTION dist! (x1!, y1!, x2!, y2!)
21.     dist! = _HYPOT((x2! - x1!), (y2! - y1!))

As seen on:
« Last Edit: August 01, 2020, 09:12:53 PM by FellippeHeitor »

_vince

• Seasoned Forum Regular
• Posts: 251
Re: Estimating Pi
« Reply #1 on: August 02, 2020, 12:45:23 AM »
It is a great simulation! Another one is called "buffon's needles", I wanted to write a program that demonstrates the experiment graphically and then numerically computes the results to iteratively approximated pi as more needles are dropped, open challenge I guess

ie
(https://en.wikipedia.org/wiki/Buffon%27s_needle_problem#/media/File:Buffon_needle_experiment_compressed.gif)
« Last Edit: August 02, 2020, 12:46:39 AM by _vince »

SierraKen

• Forum Resident
• Posts: 931
Re: Estimating Pi
« Reply #2 on: August 02, 2020, 12:55:12 AM »
Really awesome Felippe! I have no idea how you and others find Pi, was just ages ago since school. But I remembered a year ago (which turned out to be exactly a year ago), I found code to output as many digits as you want. I posted it in another thread if anyone is interested. Just funny that it's been one year, around the Sun, using Pi LOL. Oh, I also updated it tonight to make it show how many digits it has left to calculate before it shows a Notepad of it.
« Last Edit: August 02, 2020, 01:06:16 AM by SierraKen »

FellippeHeitor

• QB64 Developer
• Forum Resident
• Posts: 2209
• LET IT = BE
Re: Estimating Pi
« Reply #3 on: August 02, 2020, 01:04:52 AM »
Well, I guess that only proves that life is a circle, @SierraKen 😂

@_vince looking forward to your project!

SierraKen

• Forum Resident
• Posts: 931
Re: Estimating Pi
« Reply #4 on: August 02, 2020, 01:05:25 AM »
LOL true.

bplus

• Forum Resident
• Posts: 4548
• B+ nots
Re: Estimating Pi
« Reply #5 on: August 02, 2020, 12:49:32 PM »
It is a great simulation! Another one is called "buffon's needles", I wanted to write a program that demonstrates the experiment graphically and then numerically computes the results to iteratively approximated pi as more needles are dropped, open challenge I guess

ie
(https://en.wikipedia.org/wiki/Buffon%27s_needle_problem#/media/File:Buffon_needle_experiment_compressed.gif)

A very slow road to Pi:
Code: QB64: [Select]
1. _TITLE "Buffon Pi Calc " 'b+ 2020-08-2
2. '2020-08-02 try this with Buffon Pi Calc  pi ~ 2 * nTrials / nCross  from
3. ' pi ~ 2*L*N/(T*H) L=length of needle, T=spacing between lines N = Number of trials (or pins) H = number that cross line
4.
5. CONST xmax = 800, ymax = 600
6. SCREEN _NEWIMAGE(xmax, ymax, 32)
7. _DELAY .25
8. L = 100 ' both pin length and line spacing
9. N = 0 '   number of trials
10. C = 0 '   number of pins that cross the line
11. SS = 0 '  pin head and point lay in Same Space
12.
13.     'redraw lines
14.     FOR y = 0 TO _HEIGHT - 1 STEP L
15.         LINE (0, y)-(_WIDTH - 1, y), &HFFFFFFFF
16.
17.     'drop a random pin
18.     px1 = RND * (_WIDTH - 200) + 100
19.     py1 = RND * (_HEIGHT - 200) + 100
20.     ra = RND * _PI(2)
21.     px2 = px1 + L * COS(ra)
22.     py2 = py1 + L * SIN(ra)
23.     LINE (px1, py1)-(px2, py2), _RGB32(50 + RND * 150)
24.     N = N + 1 '                             another trial pin
25.
26.     'did our pin cross a line
27.     FOR y = 0 TO _HEIGHT - 1 STEP L
28.         IF py1 >= y AND py1 < y + L THEN ' one side of pin is between y and y + L
29.             IF py2 >= y AND py2 < y + L THEN ' the pin did NOT cross a line both the head and point lay in same space
30.                 SS = SS + 1
31.     C = N - SS
32.     'Print results of trials
33.     LOCATE 1, 1: PRINT "Trials:"; N
34.     LOCATE 2, 1: PRINT "Crossings:"; C
35.     LOCATE 3, 1: PRINT "Pi est:"; 2 * N / C
36.     IF N = 1 THEN
37.     ELSEIF N < 5 THEN
38.     ELSEIF N < 10 THEN
39.         _DELAY .5
40.     ELSEIF N < 50 THEN
41.         _DELAY .25
42.     ELSEIF N < 1000 THEN
43.         _DELAY .05
44.     ELSEIF N <= 10000 THEN
45.         IF N MOD 1000 = 0 THEN _DELAY .25
46.         IF N MOD 10000 = 0 THEN _DELAY .05
47.
48.
49.

_vince

• Seasoned Forum Regular
• Posts: 251
Re: Estimating Pi
« Reply #6 on: August 02, 2020, 01:10:45 PM »
Yes, that is it. I've uncommented CLS to watch the needles pile up. Nice work, bplus

bplus

• Forum Resident
• Posts: 4548
• B+ nots
Re: Estimating Pi
« Reply #7 on: August 02, 2020, 01:25:36 PM »
Thanks @_vince  this might be slightly more efficient for check needle crossing the line, skipping a FOR loop and maybe an extra IF THEN check:
Code: QB64: [Select]
1. _TITLE "Buffon Pi Calc " 'b+ 2020-08-2
2. '2020-08-02 try this with Buffon Pi Calc  pi ~ 2 * nTrials / nCross  from
3. ' pi ~ 2*L*N/(T*H) L=length of needle, T=spacing between lines N = Number of trials (or pins) H = number that cross line
4.
5. CONST xmax = 800, ymax = 600
6. SCREEN _NEWIMAGE(xmax, ymax, 32)
7. _DELAY .25
8. L = 100 ' both pin length and line spacing
9. N = 0 '   number of trials
10. C = 0 '   number of pins that cross the line
11. 'SS = 0 '  pin head and point lay in Same Space
12.
13.     'redraw lines
14.     FOR y = 0 TO _HEIGHT - 1 STEP L
15.         LINE (0, y)-(_WIDTH - 1, y), &HFFFFFFFF
16.
17.     'drop a random pin
18.     px1 = RND * (_WIDTH - 200) + 100
19.     py1 = RND * (_HEIGHT - 200) + 100
20.     ra = RND * _PI(2)
21.     px2 = px1 + L * COS(ra)
22.     py2 = py1 + L * SIN(ra)
23.     LINE (px1, py1)-(px2, py2), _RGB32(50 + RND * 150)
24.     N = N + 1 '                             another trial pin
25.
26.     'did our pin cross a line
27.     IF py1 < py2 THEN
28.         IF py2 >= 100 * INT(py1 / 100) + L THEN C = C + 1
29.         IF py1 >= 100 * INT(py2 / 100) + L THEN C = C + 1
30.     'FOR y = 0 TO _HEIGHT - 1 STEP L
31.     '    IF py1 >= y AND py1 < y + L THEN ' one side of pin is between y and y + L
32.     '        IF py2 >= y AND py2 < y + L THEN ' the pin did NOT cross a line both the head and point lay in same space
33.     '            SS = SS + 1
34.     '            EXIT FOR ' speed up a little
35.     '        END IF
36.     '    END IF
37.     'NEXT
38.     'C = N - SS
39.     'Print results of trials
40.     LOCATE 1, 1: PRINT "Trials:"; N
41.     LOCATE 2, 1: PRINT "Crossings:"; C
42.     LOCATE 3, 1: PRINT "Pi est:"; 2 * N / C
43.     IF N = 1 THEN
44.     ELSEIF N < 5 THEN
45.     ELSEIF N < 10 THEN
46.     ELSEIF N < 50 THEN
47.         _DELAY .25
48.     ELSEIF N < 1000 THEN
49.         _DELAY .05
50.     ELSEIF N <= 10000 THEN
51.         IF N MOD 1000 = 0 THEN _DELAY .25
52.         IF N MOD 10000 = 0 THEN _DELAY .05
53.
54.
55.

bplus

• Forum Resident
• Posts: 4548
• B+ nots
Re: Estimating Pi
« Reply #8 on: August 02, 2020, 02:21:40 PM »
There might be a problem with the code or the Random Number Generator as I have a definite pattern of "random" needle color "waves" along the bottom and right borders of the screen.

Anyone have an idea of the error, blunder or glitch?

Both versions produce this pattern.
« Last Edit: August 02, 2020, 02:23:22 PM by bplus »

_vince

• Seasoned Forum Regular
• Posts: 251
Re: Estimating Pi
« Reply #9 on: August 02, 2020, 02:33:57 PM »
looks like the bounds of the needle coordinates are too far inside the screen

bplus

• Forum Resident
• Posts: 4548
• B+ nots
Re: Estimating Pi
« Reply #10 on: August 02, 2020, 04:01:50 PM »
looks like the bounds of the needle coordinates are too far inside the screen

No I don't think it was that, here is the fix of the "pattern" and even speedier calc's seems 3.140 is limit of precision.

Code: QB64: [Select]
1. _TITLE "Buffon Pi Calc " 'b+ 2020-08-2
2. '2020-08-02 try this with Buffon Pi Calc  pi ~ 2 * nTrials / nCross  from
3. ' pi ~ 2*L*N/(T*H) L=length of needle, T=spacing between lines N = Number of trials (or pins) H = number that cross line
4.
5. CONST xmax = 800, ymax = 600
6. SCREEN _NEWIMAGE(xmax, ymax, 32)
7. _DELAY .25
8.
9. L = 100 ' both pin length and line spacing
10. N = 0 '   number of trials
11. C = 0 '   number of pins that cross the line
12.     'CLS
13.     'redraw lines
14.     'drop a random pin
15.     px1 = RND * (_WIDTH - 201) + 100
16.     py1 = RND * (_HEIGHT - 201) + 100
17.     ra = RND * _PI(2)
18.     IF RND < .5 THEN dir = -1 ELSE dir = 1 ' << this fixed the pattern developing along right and bottom border of screen
19.     px2 = px1 + dir * L * COS(ra)
20.     py2 = py1 + dir * L * SIN(ra)
21.
22.     N = N + 1 '                             another trial pin
23.
24.     'did our pin cross a line
25.     IF py1 < py2 THEN
26.         IF py2 >= 100 * INT(py1 / 100) + L THEN C = C + 1
27.         IF py1 >= 100 * INT(py2 / 100) + L THEN C = C + 1
28.     IF N < 100000 THEN 'display
29.         FOR y = 0 TO _HEIGHT - 1 STEP L
30.             LINE (0, y)-(_WIDTH - 1, y), &HFFFFFFFF
31.         LINE (px1, py1)-(px2, py2), _RGB32(50 + RND * 150)
32.
33.         'Print results of trials
34.         LOCATE 1, 1: PRINT "Trials:"; N
35.         LOCATE 2, 1: PRINT "Crossings:"; C
36.         LOCATE 3, 1: PRINT "Pi est:"; 2 * N / C
37.         IF N = 1 THEN
38.         ELSEIF N < 5 THEN
39.         ELSEIF N < 10 THEN
40.         ELSEIF N < 50 THEN
41.             _DELAY .25
42.         ELSEIF N < 1000 THEN
43.             _DELAY .05
44.         ELSEIF N <= 10000 THEN
45.             IF N MOD 1000 = 0 THEN _DELAY .25
46.             IF N MOD 10000 = 0 THEN _DELAY .05
47.         IF N MOD 1000000 = 0 THEN
48.             LOCATE 1, 1: PRINT "Trials:"; N
49.             LOCATE 2, 1: PRINT "Crossings:"; C
50.             LOCATE 3, 1: PRINT "Pi est:"; 2 * N / C
51.