Concurrent Monte Carlo Estimate of Pi in Python

Last Updated on September 12, 2022

The ProcessPoolExecutor class in Python can be used to estimate Pi by performing multiple Monte Carlo simulations at the same time.

This can dramatically speed-up your program compared to using a single CPU core to run simulations.

In this tutorial, you will discover how to estimate pi using a concurrent implementation of the Monte Carlo method.

After completing this tutorial, you will know:

How to estimate the value of Pi using the Monte Carlo method and how to implement it in Python.
How to use the ProcessPoolExecutor to manage a pool of worker processes.
How to use the ProcessPoolExecutor to improve the estimate of Pi using the Monte Carlo method.

Let’s dive in.

Table of Contents

Monte Carlo Estimate Pi Sequentially (slowly)

Pi is a mathematical constant and has an approximate value of 3.14159.

It is defined by the ratio of a circle’s circumference (around the outside) to its diameter (across the middle), and it is used all throughout geometry.

The precision or number of decimal points for Pi does not end, making it a fascinating number.

There are many ways to estimate the value of Pi, but the Monte Carlo method is a particularly interesting method that we can use as the basis for exploring concurrent programming in Python.

Monte Carlo Estimate for Pi

Python has a value of Pi that we can use directly as the Pi constant in the math module.

math — Mathematical functions

Nevertheless, we can calculate a value for Pi.

The Monte Carlo method is a statistical technique for estimating quantities using simulations.

There is a standard procedure for estimating Pi using the Monte Carlo method, as follows:

Draw a unit square, e.g. width and height are 1.0.
Draw a quarter circle, a quadrant, in the square from the top left to the bottom right points.
Uniformly scatter points in the square.
Count the number of points in the quadrant.
Estimate Pi as the ratio of points in the quadrant over the total points multiplied by four.

Or put simply:

Pi = (points in quadrant / total points) * 4

Monte Carlo Estimate for Pi in Python

We can implement the Monte Carlo estimate of Pi in Python.

First, we need to generate uniformly random points in the unit square.

This can be achieved by using the random.random() function.

random — Generate pseudo-random numbers

We can call the random() function for both the x and y coordinates of each point. Each point can be defined as a tuple of floating point values and we can generate the desired number of tuples in a list comprehension.

For example:

...

# generate uniformly random points in [0,1)

points = [(random(),random()) for _ in range(n)]

We can then check the number of points that are within the quadrant of the circle.

This can be achieved by checking if the distance between the point and the middle of the circle is less than the radius of the circle.

The center of the circle is at the coordinate (0,0), therefore the distance can be calculated as the sum of the square of each coordinate, for example:

distance = x^2 + y^2

If this distance is less than 1, we know the point is within the circle.

For example:

# count all points in the circle

count = 0

for x,y in points:

# euclidean distance from origin (0,0)

distance = (x**2 + y**2)

# check if within radius

if distance <= 1:

count += 1

This is a little bloated, we can simplify the code into a list comprehension that creates a list of 1 values for all points in the circle that we can then sum, for example:

...

# count the number of points in the circle

result = sum([1 for x,y in points if (x**2 + y**2) < 1])

We can tie this together and define a function task() that generates the desired number of points n and returns the total number of points that fall within the quadrant.

# monte carlo simulation of points in a quadrant within unit square

def task(n):

# generate random points

points = [(random(),random()) for _ in range(n)]

# count the number of points in the circle

return sum([1 for x,y in points if (x**2 + y**2) < 1])

Now it is just a matter of calling this function and estimating the value of pi.

First, we can define the number of points to generate, such as 100 million, then count the number of points in the quadrant.

...

# define number of points per process

n = 1000000 * 100

# count total points in the quadrant

total_in_quadrant = task(n)

We can then estimate the value of Pi using this count and report the value.

...

# estimate pi

pi_estimate = total_in_quadrant / n * 4

print(f'pi (estimate):\t{pi_estimate}')

It might also be interesting to report the actual value of Pi in the math module as a point of comparison.

1 2	... print(f'pi (actual):\t{pi}')

That’s it.

Tying this together, the complete example of estimating Pi using the Monte Carlo method is listed below.

# SuperFastPython.com

# estimate pi using the monte carlo method

from random import random

from math import pi

# monte carlo simulation of points in a quadrant within unit square

def task(n):

# generate random points

points = [(random(),random()) for _ in range(n)]

# count the number of points in the circle

return sum([1 for x,y in points if (x**2 + y**2) < 1])

# entry point

if __name__ == '__main__':

# define number of points per process

n = 1000000 * 100

# count total points in the quadrant

total_in_quadrant = task(n)

# estimate pi

pi_estimate = total_in_quadrant / n * 4

print(f'pi (estimate):\t{pi_estimate}')

print(f'pi (actual):\t{pi}')

Running the example generates the points and counts the number of points that land within the quadrant then uses this to estimate a value for pi.

We can see that the estimate is a rough approximation of the actual value down to 3 decimal places (3.141).

It is also slow.

On my system this example took about 36.1 seconds. That is a lot of work for a very rough approximation of pi.

How long did it take on your system?
Let me know in the comments below.

1 2	pi (estimate): 3.14138952 pi (actual): 3.141592653589793

Next, let’s look at how we can adapt our program for estimating Pi to make use of a process pool to improve the estimate.

Run loops using all CPUs, download your FREE book to learn how.

Monte Carlo Estimate Pi Concurrently

We can update our program for estimating Pi to use multiple CPU cores with very little change.

The task() function can be called for each CPU and we can sum the counts and divide the result by the sum of all points generated.

First, let’s define the number of process workers we want in the pool.

This should match the number of physical CPUs in your system. We could use the number of logical CPUs (after hyperthreading) but this would likely lead to worse performance (test and see).

I have four physical CPU cores, so we will use this many workers (change for your system).

...

# define the number of workers

workers = 4

We can then create the ProcessPoolExecutor with the desired number of workers using the context manager.

...

# create the process pool

with ProcessPoolExecutor(workers) as executor:

# ...

Next we can submit the tasks.

We can use the map() function to call the task() function with the desired number of points (100 million) for each physical CPU in the system (e.g. 4).

This can be achieved by defining a list with 4 elements, where each element is 100 million defined by n.

For example:

...

# define number of points per process

n = 1000000 * 100

# input to the map

inputs = [n, n, n, n]

Or more simply:

...

# define number of points per process

n = 1000000 * 100

# input to the map

inputs = [n] * workers

The call to map() will return an iterable over the number of points in the quadrant.

...

# submit the tasks

counts = executor.map(task, [n]*workers)

We just want the sum of these values, therefore we can sum the return value from map() directly:

...

# submit tasks

total_in_quadrant = sum(executor.map(task, [n]*workers))

Then, when we estimate the value of Pi, we divide the sum of the points in quadrants by the total number of workers multiplied by the number of points generated per worker.

...

# estimate pi

pi_estimate = total_in_quadrant / (workers * n) * 4

And that’s it.

Tying this together, the complete example is listed below.

Our expectation is that that program takes about the same amount of time to complete, but will generate 4 times the number of points.

# SuperFastPython.com

# estimate pi using the monte carlo method concurrently

from random import random

from math import pi

from concurrent.futures import ProcessPoolExecutor

# monte carlo simulation of points in a quadrant within unit square

def task(n):

# generate random points

points = [(random(),random()) for _ in range(n)]

# count the number of points in the circle

return sum([1 for x,y in points if (x**2 + y**2) < 1])

# entry point

if __name__ == '__main__':

# define the number of workers

workers = 4

# define number of points per process

n = 1000000 * 100

# create the process pool

with ProcessPoolExecutor(workers) as executor:

# submit tasks

total_in_quadrant = sum(executor.map(task, [n]*workers))

# estimate pi

pi_estimate = total_in_quadrant / (workers * n) * 4

print(f’pi (estimate):\t{pi_estimate}')

print(f'pi (actual):\t{pi}')

Running the program estimates the value of pi using the Monte Carlo method as before.

In this case, we generate 4 times the number of points to use in the estimate, or 400 million points.

The program runs for approximately the same amount of time on my system, e.g. 39.7 seconds compared to 36.1 seconds with the serial case, or about 3.5 seconds longer.

We can see that this improves the estimate of pi slightly in this case, adding one more decimal point of precision (3.1415).

Note, your specific estimate of Pi will differ given the use of randomness in the program.

1 2	pi (estimate): 3.14157306 pi (actual): 3.141592653589793

Download Now: Free ProcessPoolExecutor PDF Cheat Sheet

Extensions

This section lists ideas for extending the tutorial.

Increase Workers to Match Logical CPUs. Update the example to use the number of logical CPUs in your system and compare the execution time to using the number of physical CPUs.
Increase Number of Points Per Worker. Increase the number of points generated per worker to five or ten times the number and compare the precision in the estimate of Pi.
Increase Number of Tasks. Update the example to submit smaller blocks of numbers to each task and perhaps use submit() and as_completed() to update the estimate of pi as each task is completed.

Share your extensions in the comments below, it would be great to see what you come up with.

Free Python ProcessPoolExecutor Course

Download your FREE ProcessPoolExecutor PDF cheat sheet and get BONUS access to my free 7-day crash course on the ProcessPoolExecutor API.

Discover how to use the ProcessPoolExecutor class including how to configure the number of workers and how to execute tasks asynchronously.

Learn more