Pythonmatics

Solving mathematical problems with Python

Monte Carlo Pie


The Monte Carlo algorithm uses randomness to provide an estimated result for a mathematical problem instead of computing the exact result, which might be complicated. For example, the probability of getting the sum 6 when two dice are thrown is 5/36 (approximately 0.1389). This code estimates this probability by simulating 1,000,000 throws and calculating the event ratio:

from random import *

draws = 1000000
results = [randrange(1,7) + randrange(1,7) for _ in range(draws)]
print(results.count(6)/draws)

Output:

0.139124

In most cases, such as the example above, the Monte Carlo algorithm is used to provide an estimation of the probability of a certain event occurring. However, the more intriguing application of the Monte Carlo randomized algorithm involves problems that are not intuitively related to probability.

One such example is using a Monte Carlo-style code to estimate the value of π. Here is a way to connect the value of π to a probability problem: Imagine a square-shaped target with a circle inscribed within it.

Assuming that the radius of the circle is 1 unit (meaning the edges of the square are 2 units), then the ratio between the circle’s area and the square’s area is (π·12)/(22), or π/4. Statistically, if a shooter randomly fires shots at the target, then the ratio between the number of hits within the circle and the total number of shots should be around π/4. Alternatively, the product of the measured ratio by 4 should be around π.

Of course, there’s no need for a human shooter or to waste ammunition. That’s what computer programming is made for.

It’s Python time

from random import *
from math import *

rounds = 1000000

shots = [(uniform(-1,1),uniform(-1,1))
                        for _ in range(rounds)]  # simulate random coordinates between (0,0) and (1,1)
dists = [sqrt(x*x+y*y) for (x,y) in shots]       # calculate the distance of each shot from the origin
hits = len([i for i in dists if i < 1])          # count the number of hits within the circle

print(hits/rounds * 4)

Output

3.141444

Discover more from Pythonmatics

Subscribe now to keep reading and get access to the full archive.

Continue reading