Sunday, January 20, 2008

Pi, Probability, and Python

π and probability theory intersect more often than you would expect. The mix gets more interesting when you combine that with a quick look at math history, and toss in a short program as well.

Buffon's Problem
George Louis Leclerc, Comte de Buffon (1707-1788) was a general scientist, mathematician, and gambler, as well as a man of means with time to spare on interesting math problems. One which he put forth and subsequently solved was this: Given a needle of length L, to be tossed upon a horizontal planed with parallel straight lines separated uniformly by a distance d, what is the probability that the needle will intersect with one of the lines?

Assume that the position and orientation of the needle are two variables that are random and independent. Let the distance from the needle's center to the nearest line be x, with the orientation of the needle being the angle φ created by the needle intersecting the line in question. Consider x and φ per the needle's intersection with a single line, since all intersections with any other line will be the same.

Now from the figure above, it can be seen that the needle will intersect a line if and only if

x < ½ L sinφ

Therefore, Buffon's problem can be resolved by finding the following probability:
P(x < ½ L sinφ)

In order to find this probability, we can visually represent this by graphing with the rectangular coordinates x and φ, where
range: 0 < x < d/2
domain: 0° < φ < 180°

Given these intervals for x and φ, any point that lies within the parallelogram described by points A and B correspond to one and only one position and orientation of the needle thrown. Furthermore, since all such combinations are equally probable, the area of AB represents the sum of all probable outcomes for the thrown needle's position and orientation, with the area under the curve x < ½ L sinφ representing the sum of all possibilities where the needle intersects a line.

Hence, the probability that the thrown needle intersects a line is the ratio of the the sum of possible outcomes where the needle intersects a line, to the sum total of all possible outcomes. This can be expressed by:


Laplace and the Monte Carlo Method
Pierre Simon Laplace (1749-1827) saw Buffon's solution from a different perspective, which gave rise to a new method in computing. From that last expression above, Laplace generalized it to:

The end result was a new method whereby π could now be calculated by determining the probability P. The calculation is simplified if one sets d equal to L. But the key point is that Laplace discovered a method whereby a numerical value is found by realizing a random event many times and observing the outcomes. This method is known as the Monte Carlo method.

So where does the Python come in? I wrote a method in Python for computing π in Monte Carlo fashion:

from math import *
from random import *

def execute():
count = 0.0
step = 1000.0
for i in range(1, 10001):
for j in range(0, step):
x = 0.5 * random()
phi = randrange(0, 180, 1)
if x < 0.5 * sin(radians(phi)):
count += 1
P = count/(i * step)
print "%d: P = %0.25f, pi = %0.25f" % (i * step, P, 2.0/P)
print "all pau!"


This runs a Monte Carlo simulation of 10 million throws of the needle, computing the probability of intersection every 1000 steps of the way. One execution yielded the following output:

...
9995000: P = 0.6366198099049524827819369, pi = 3.1415924683503022585284725
9996000: P = 0.6366196478591437113436768, pi = 3.1415932680144256217147358
9997000: P = 0.6366192857857356779405222, pi = 3.1415950547767912404140134
9998000: P = 0.6366201240248049453285262, pi = 3.1415909182319108339243030
9999000: P = 0.6366206620662065995830403, pi = 3.1415882631091953669510985
10000000: P = 0.6366190000000000459579041, pi = 3.1415964650756573739442956
all pau!


Hm... only correct to 5 places after 10 million trials? I think I may have to see about calculating the probably of computing π to n places in t trials. While it may not be all that efficient for calculating π, the Monte Carlo method is nevertheless a powerful approach for resolving problems with incredibly large sets of computer-simulated trials.

With inspiration from
A History Of Pi, P. Beckmann
History of Mathematics, D.E. Smith

No comments: