# Calculate Π By Throwing Darts

We can approximate π using nothing more than random numbers and some simple geometry: we randomly throw darts at a square board of side r; within the square we inscribe a quadrant of a circle of radius r with its centre at (0, 0). We count all of the 'throws'; if a dart lands within the quadrant, we also count a 'hit'.

For a large number of throws, we see that:

$\frac{\mathrm{hits}}{throws}=\frac{\mathrm{area-of-quadrant}}{area-of-square}$

Some half-remembered geometry tells us that:

$\frac{\mathrm{area-of-quadrant}}{\mathrm{area-of-square}}=\frac{{\mathrm{\pi r}}^{2}/4}{{r}^{2}}=\frac{\pi }{4}$

Or:

$\pi =4\left(\frac{\mathrm{area-of-quadrant}}{\mathrm{area-of-square}}\right)=4\left(\frac{\mathrm{hits}}{\mathrm{throws}}\right)$

I first solved this problem as an undergraduate sometime in 1994 as part of a Computational Physics module. Using FORTRAN.

#!/usr/bin/env python

import random
import math

# class representing a single throw of a dart
class Throw:
def __init__(self):
# generate two random coordinates and work out how far away we are from
# the origin
self.x = random.random()
self.y = random.random()
self.distance = self.distance()

def distance(self):
# the distance from the origin is the hypotenuse of a right-angled
# triangle with sides of length and x and y. Pythagoras told us that:
#    distance = sqrt((x^2) + (y^2))
# which looks like this in python
return math.sqrt(self.x**2 + self.y**2)

# did we land inside the quadrant?
def is_a_hit(self):
return self.distance <= 1.0

# main class
class MonteCarlo:
def __init__(self):
# initialise everything
self.hits = 0
self.throws = 0
self.pi = 0

# this method is called on every throw
def increment(self, throw):
# we always clock up another throw
self.throws += 1
# and accumulate a hit if we scored
if throw.is_a_hit():
self.hits += 1
# then get a new value of pi
self.calculate_pi()

# explanation can be found here: http://icanhaz.com/montecarlo
def calculate_pi(self):
self.pi = 4 * (float(self.hits) / float(self.throws))

# we use this in determining whether to print our status
def divides_by(self, number):
return float(self.throws) % number == 0

# represent the current state as a string
def __repr__(self):
return "Throws: %10d, Hits: %10d, Pi: %10f, Actual Pi: %10f" % (self.throws, self.hits, self.pi, math.pi)

# if we're called on the command line
if __name__ == '__main__':
import sys
try:
step = int(sys.argv[1])
except IndexError:
step = 1000
# construct a new instance
m = MonteCarlo()
# loop forever
while 1:
# keep throwing darts
m.increment(Throw())
# only print on every nth iteration
if m.divides_by(step): print m


There's now a ruby version, too.

Published under a Creative Commons Attribution-NonCommercial licence