Sunday, March 15, 2009

Adventures in Computing Pi

My Pi Day cake illustrates a method that you could use to compute π. If the circular layer has radius R, then its area is πR2. And since the circle has radius R, then the side of the square is of length 2R, for an area of (2R)2 = 4 R2. Taking the ratio of these areas, we obtain
πR2/(4 R2) = πR2/(4 R2)= π/4.

So if we could figure out a way to measure the areas of the circle and square accurately, we could approximate π.

One way to find the area of a shape is known as Monte Carlo integration. The idea behind Monte Carlo methods is that if you threw darts blindly, the percent of time that you hit inside the shape you're trying to find the area of would be proportional to the area of the shape. In other words, if I had a "blob" shape inside a square, and I threw a bunch of darts "randomly"* at the square, then the relative area of the blob to the square would be the ratio of the number of darts that landed inside the blob to the total number of darts thrown. The more darts we throw, the better accuracy we get.

If you look at the cake, it is decorated with sprinkles.

Ideally, we could count the number of sprinkles on the circular layer and compare that to the total number of sprinkles, and we should get π/4. In reality, due to the imperfections of the cake (e.g., the rounded edges of the square layer, the imperfect shape of the circular layer, the challenge of distributing the sprinkles evenly yet pseudorandomly) the ratio is probably off.

But even if it were a perfect cake with perfectly distributed sprinkles, we would not have a particularly accurate value for π. This is because at most, there are a couple hundred sprinkles on the cake. The error in our calculation of π is proportional to one over the square root of N, where N is the number of sprinkles. So, if we have 100 sprinkles, this means our calculation is right to within ±1/10 of π, or ~0.3. For each additional digit of accuracy, we'll need to increase the number of sprinkles by a factor of 100. So, if we wanted π to three decimal digits of accuracy, we'd need to have on the order of a million sprinkles.

As you can tell, this is not the most efficient way of computing π. In fact, if you needed an accurate value of π to many decimal places, this would definitely not be the way to go. Instead, allow me to suggest the following methods:
  1. Look it up on the internet. Even bigger nerds than me have computed and posted π to many decimal places.
  2. Most computer languages have a constant PI that is accurate to at least 15 decimal places. Generally that's all you need.
  3. You can use a numerical quadrature method to find the integral of a function on an interval which is known to equal π (e.g., the integral of sqrt(1-x2) from 0 to 1 = π/4). You could use a quadrature rule such as Simpson's rule, a higher-order Runge-Kutta method, or a Gaussian quadrature rule, which, in one dimension, have much faster convergence than Monte Carlo integration. (Monte Carlo does not become competitive until we're doing three or higher dimensional integration.)
  4. Finally, there are many formulae that can be used, e.g., the BBP Formula.

* While there's actually no way to throw darts randomly, you can distribute them pseudorandomly. But this is a topic for another post.

No comments: