Parallel Craps in Python

Code along video

Learning objectives

Summing the elements of a list in parallel was a nice first example, but it never beat the serial approach. In this code along, we return to the game of craps from Chapter 2 and estimate its house edge in parallel, splitting the trials among several worker processes that each report back through a queue. We will then time the parallel and serial approaches against each other and see that, once we run enough trials, parallel programming finally pays off.

Code along summary

Setup

Create a new folder parallel_craps inside your python/src directory. Our project will contain three files: serial.py, which holds the serial craps simulator that we already wrote; parallel.py, where we will write the parallel version; and main.py, which runs and times the two approaches.

Add the following code to serial.py. These functions are the ones that we wrote in Chapter 2 when we simulated craps with Monte Carlo simulation; the only change is that we have renamed the highest level function compute_craps_house_edge_serial() so that we can tell it apart from its parallel cousin.

"""
serial.py

Serial craps simulation: roll dice, play one game, compute the
house edge over many trials.
"""

import random


def roll_die() -> int:
    """Simulate rolling a single six-sided die."""
    return random.randint(1, 6)


def sum_dice(num_dice: int) -> int:
    """Simulate rolling num_dice six-sided dice and summing the results."""
    total = 0
    for _ in range(num_dice):
        total += roll_die()
    return total


def play_craps_once() -> bool:
    """
    Simulate a single game of craps.
    Return True if player wins, False if loses.
    """
    first_roll = sum_dice(2)
    if first_roll == 7 or first_roll == 11:
        return True
    elif first_roll == 2 or first_roll == 3 or first_roll == 12:
        return False
    else:
        while True:
            next_roll = sum_dice(2)
            if next_roll == first_roll:
                return True
            elif next_roll == 7:
                return False


def compute_craps_house_edge_serial(num_trials: int) -> float:
    """
    Simulate num_trials games of craps in serial.
    """
    win_total = 0
    for _ in range(num_trials):
        if play_craps_once():
            win_total += 1
        else:
            win_total -= 1
    return win_total / num_trials

Next, add the following starter code to parallel.py. It imports multiprocessing and borrows play_craps_once() from serial.py; we will see why shortly. It also contains the signatures of the two functions that we will write.

"""
parallel.py

Parallel craps simulation using multiprocessing.Process + Queue.

Divide num_trials evenly across num_procs worker processes, run
play_craps_once in each, and combine the per-worker totals to get the
house edge.
"""

import multiprocessing

from serial import play_craps_once


def compute_craps_house_edge_parallel(num_trials: int, num_procs: int) -> float:
    """
    Simulate num_trials trials of craps in parallel, divided over
    num_procs worker processes, and return the house edge of the game
    as the amount won/lost (+/-) divided by num_trials.
    """
    pass


def total_win_one_proc(num_trials_one_proc: int,
                       result_queue: multiprocessing.Queue) -> None:
    """
    One worker bee playing craps num_trials_one_proc times.
    It places its amount won/lost (+/-) into the input queue.
    """
    pass

Finally, add the following code to main.py. This is our old friend the timing code from the previous code along: it times the parallel approach, then times the serial approach, and prints the speedup, which is the ratio of the serial run time to the parallel run time. If the parallel approach is three times faster, then the serial approach takes three times longer, and the speedup is 3.

import time
import multiprocessing

from serial import compute_craps_house_edge_serial
from parallel import compute_craps_house_edge_parallel


def main() -> None:
    print("Parallel programming and craps.")

    num_procs = multiprocessing.cpu_count()
    num_trials = 1000

    # time the parallel approach
    start_parallel = time.time()
    compute_craps_house_edge_parallel(num_trials, num_procs)
    end_parallel = time.time()
    total_parallel = end_parallel - start_parallel
    print(f"Time taken for parallel craps house edge computation of {num_trials} trials using {num_procs} processes is: {total_parallel}")

    # time the serial approach
    start_serial = time.time()
    compute_craps_house_edge_serial(num_trials)
    end_serial = time.time()
    total_serial = end_serial - start_serial
    print(f"Time taken for serial craps house edge computation of {num_trials} trials is: {total_serial}")

    print(f"Speedup provided by parallel approach: {total_serial / total_parallel}")


if __name__ == "__main__":
    main()

A quick refresher on our craps simulator

Our craps simulator was our first introduction to modular, top-down design, and that design is about to pay off. The highest level function, compute_craps_house_edge_serial(), plays num_trials games, assuming that we wager a single unit of currency each time. If we win a game, our win total increases by one; if we lose, it decreases by one. At the end, we divide the win total by the number of trials, so that a negative result means that we lost money on average, and a positive result means that we won.

All of the hard work is passed to play_craps_once(), which returns True if we win and False if we lose. This is where the rules of craps live: we roll two dice, and a 7 or 11 on the first roll wins, whereas a 2, 3, or 12 loses. Otherwise, we keep rolling until we either match our first roll (a win) or roll a 7, which is now a losing number. Below play_craps_once(), sum_dice() simulates rolling any number of dice by calling roll_die() repeatedly, and roll_die() returns a random integer between 1 and 6.

Here is the key point: to parallelize this simulation, we do not need to know anything about the rules of craps. We only need to split the work of playing many games among several workers.

A five-step recipe for parallel programs

When we summed a list in parallel, we followed an order of operations that will apply in many parallel programs. Let’s write it down as comments inside compute_craps_house_edge_parallel(), replacing pass.

# 1. create a queue to store our results
# 2. initialize our worker processes
# 3. start the processes, which will put their return values into the queue
# 4. range over all processes and "join" them (wait until all are finished before continuing)
# 5. grab all results from the queue and combine

Why a queue? When many pieces of our code are running at the same time, it becomes hard to say what it means for a function to “return” a value. A queue is accessible to every worker process, so it gives us a shared place to collect results: when a worker finishes, it simply drops whatever it would have returned into the queue.

Before writing any parallel code, let’s fill in what the function does at a high level. Picture yourself sitting in the middle of the casino while a crowd of players each plays some number of games and then comes back to tell you how much they won or lost. You keep a running win_total over all of the workers, and at the end, you divide it by the total number of trials, just as in the serial case.

# win_total stores total amount won/lost over all workers
win_total = 0

return win_total / num_trials

Step 1 is to create a queue to store our results, along with an empty list that will hold our worker processes. Add these lines after the declaration of win_total.

# create our queue and an empty collection of processes
result_queue = multiprocessing.Queue()
processes = []

How many trials should each worker get?

Before we initialize the workers, we need to decide how much work each one receives.

With 100,000 trials and 100 workers, each worker should play 1,000 games. In general, each worker gets num_trials divided by num_procs trials. Because this needs to be an integer, we use integer division (//).

# how many trials should each worker be given?
trials_one_proc = num_trials // num_procs

# be careful: num_trials = 1003, num_procs = 10
# trials_one_proc = 1003 // 10 = 100
# final worker should get 3 extra trials

As with summing a list, we need to be careful. If num_trials is 1,003 and num_procs is 10, then trials_one_proc is 1003 // 10, which is 100. If we give each of 10 workers 100 trials, we only play 1,000 games and come up three short. So the first num_procs - 1 workers will each get trials_one_proc trials, and the final worker will pick up the leftovers.

Starting the first workers

Steps 2 and 3 are to initialize and start the first num_procs - 1 workers. We don’t care about the loop index, so we use an underscore. Each worker is a multiprocessing.Process whose target is the function it runs and whose args is a tuple of that function’s arguments. The function that each worker runs is total_win_one_proc(), which takes the number of trials that the worker plays and the queue into which it places its result.

# start first num_procs - 1 workers
for _ in range(num_procs - 1):
    p = multiprocessing.Process(target=total_win_one_proc, args=(trials_one_proc, result_queue))
    processes.append(p)
    p.start()

We append each process to processes because we will need to join them all later, and p.start() begins the current worker’s job. Note that the for loop does not wait for a worker to finish before creating the next one. Each worker starts its work, and the loop carries on. In fact, if we give the workers enough to do, our code could reach the return statement while the workers are still playing craps, which is why we will need a separate step to wait for them.

Next, we start the final worker. How many trials does it get?

The extra trials are the remainder when we divide num_trials by num_procs, which we compute with the modulus operator, %. In our example, 1003 % 10 is 3. The final worker gets a little extra work, but not much, since the remainder is always smaller than num_procs. Once we know how many trials the final worker receives, we create its process in the same way as before.

# start final worker's job too
# how many extra trials does the final worker get?
extra_trials = num_trials % num_procs
trials_final_proc = trials_one_proc + extra_trials
p = multiprocessing.Process(target=total_win_one_proc, args=(trials_final_proc, result_queue))
processes.append(p)
p.start()

Waiting on the workers and combining their results

Step 4 is to wait for all of the workers to finish by ranging over the processes and calling join() on each one.

# we need to wait on all the workers to finish
for p in processes:
    p.join()
    # two ways of thinking of this:
    # 1. every time through the loop, the current process will not allow the code to continue executing unless this process has finished
    # 2. these two lines of code wait for every worker to finish

There are two ways to think about these two lines. First, each time through the loop, p.join() will not let our code continue executing until process p has finished. Second, the only way to make it all the way through the loop is if every worker is done. Either way, once we pass this loop, every worker has placed its result into the queue.

Step 5 is to collect the results. Each value in the queue is the total amount won or lost by one worker. We have no control over the order in which the workers report back, but we don’t need it, since all we do is add the values together. There is one value per worker, so we call result_queue.get() num_procs times; each call grabs whichever value is at the front of the line.

# collect results from queue
# each worker is "reporting back" with amount won/lost
for _ in range(num_procs):
    win_total += result_queue.get()

Putting everything together, our function is complete.

def compute_craps_house_edge_parallel(num_trials: int, num_procs: int) -> float:
    """
    Simulate num_trials trials of craps in parallel, divided over
    num_procs worker processes, and return the house edge of the game
    as the amount won/lost (+/-) divided by num_trials.
    """
    # 1. create a queue to store our results
    # 2. initialize our worker processes
    # 3. start the processes, which will put their return values into the queue
    # 4. range over all processes and "join" them (wait until all are finished before continuing)
    # 5. grab all results from the queue and combine

    # win_total stores total amount won/lost over all workers
    win_total = 0

    # create our queue and an empty collection of processes
    result_queue = multiprocessing.Queue()
    processes = []

    # how many trials should each worker be given?
    trials_one_proc = num_trials // num_procs

    # be careful: num_trials = 1003, num_procs = 10
    # trials_one_proc = 1003 // 10 = 100
    # final worker should get 3 extra trials

    # start first num_procs - 1 workers
    for _ in range(num_procs - 1):
        p = multiprocessing.Process(target=total_win_one_proc, args=(trials_one_proc, result_queue))
        processes.append(p)
        p.start()

    # start final worker's job too
    # how many extra trials does the final worker get?
    extra_trials = num_trials % num_procs
    trials_final_proc = trials_one_proc + extra_trials
    p = multiprocessing.Process(target=total_win_one_proc, args=(trials_final_proc, result_queue))
    processes.append(p)
    p.start()

    # we need to wait on all the workers to finish
    for p in processes:
        p.join()
        # two ways of thinking of this:
        # 1. every time through the loop, the current process will not allow the code to continue executing unless this process has finished
        # 2. these two lines of code wait for every worker to finish

    # collect results from queue
    # each worker is "reporting back" with amount won/lost
    for _ in range(num_procs):
        win_total += result_queue.get()

    return win_total / num_trials

The worker bee

We have done the hard work, but total_win_one_proc() is still empty. What does one worker bee do? It plays craps num_trials_one_proc times and then places the result into the queue that it was given. This function is easier than you might think.

def total_win_one_proc(num_trials_one_proc: int,
                       result_queue: multiprocessing.Queue) -> None:
    """
    One worker bee playing craps num_trials_one_proc times.
    It places its amount won/lost (+/-) into the input queue.
    """
    win_total = 0

    # range over the number of trials this worker gets, and simulate the game each time, adding in the result to win_total
    for _ in range(num_trials_one_proc):
        if play_craps_once():
            # win!
            win_total += 1
        else:
            # lost :(
            win_total -= 1

    result_queue.put(win_total)

There is nothing parallel about this function, and it looks almost exactly like compute_craps_house_edge_serial(). It differs in only two ways. First, it does not divide by the number of trials; the function that combines the results takes care of that. Second, instead of returning its win total, it calls result_queue.put() to place the win total into the queue.

And that’s it! We don’t need to write anything else, because play_craps_once() is imported from serial.py, and its subroutines sum_dice() and roll_die() were written long ago. Once we have a modular serial implementation, a parallel implementation is quick to write and plugs right into the lower level subroutines that we already have.

Checking our house edge

Before timing anything, let’s make sure that our parallel function produces a sensible answer. Add the following two lines to main() after the declarations of num_procs and num_trials. While we check the edge, we don’t want the timing code to run, so temporarily switch it off. Phillip does this by placing a triple quote (""") on the line before the comment # time the parallel approach and another after the final print(), which turns the timing code into a string that Python ignores; commenting out each line with # works just as well.

edge = compute_craps_house_edge_parallel(num_trials, num_procs)
print("House edge:", edge)

Navigate into python/src/parallel_craps and run python3 main.py (or python main.py on Windows). The edge is the average amount won or lost per game, so a negative value means that the house is winning. With 1,000 trials, the estimate bounces around from run to run, and it may even come out positive. Your numbers will differ.

Parallel programming and craps.
House edge: -0.06

As we saw when we first introduced Monte Carlo simulation, if we raise the number of trials to 100,000 or a million, the law of large numbers takes over, and the estimates tighten around the true house edge of about 1.4 cents lost per dollar wagered. Once you are convinced that the function works, switch the timing code back on and switch the two lines that you just added off (Phillip moves the triple quotes so that they surround these two lines instead).

Timing the speedup

Now for the part that we really care about. We start with 1,000 trials. Phillip’s laptop has 10 processors, so multiprocessing.cpu_count() returns 10, and each of the 10 workers plays 100 games.

When we run our code, the result is not encouraging. On our computer, which has 8 cores, we obtained the following (your times will differ).

Parallel programming and craps.
Time taken for parallel craps house edge computation of 1000 trials using 8 processes is: 1.083463191986084
Time taken for serial craps house edge computation of 1000 trials is: 0.003910064697265625
Speedup provided by parallel approach: 0.0036088578976994413

In the video, the parallel approach takes 0.18 seconds while the serial approach takes 0.0013 seconds, for a speedup of just 0.007. A speedup below 1 means that the parallel approach is slower, and a speedup close to 0 means that it is much slower. The problem is that creating and starting processes carries overhead, and 100 games of craps per worker is not enough work to justify that overhead.

So let’s give the workers more to do. Change num_trials to 10,000, so that each of Phillip’s workers plays 1,000 games, and run the code again. The speedup jumps from 0.007 to 0.09. Still not worth it, but a big jump. At 100,000 trials, the speedup climbs to 0.72, and we have almost caught up. At a million trials, the parallel approach finally wins: it finishes in 0.47 seconds while the serial approach takes 1.34, a speedup of 2.87. Here is what we saw on our 8-core machine at a million trials.

Parallel programming and craps.
Time taken for parallel craps house edge computation of 1000000 trials using 8 processes is: 3.188666343688965
Time taken for serial craps house edge computation of 1000000 trials is: 10.620057821273804
Speedup provided by parallel approach: 3.3305641533467782

Finally, at 5 million trials, Phillip’s parallel simulation finishes in 1.4 seconds while the serial one takes 6.7 seconds, for a speedup of 4.76. The figure below plots the speedups from the video.

The best speedup that we could ever hope for on a 10-core machine is 10. But if this Monte Carlo simulation keeps improving as the input grows, we could move it to a machine with 1,000 cores, and a speedup of 500 would make writing our code in parallel more than worth it. This is where we begin to distinguish parallel programming on a single computer from distributed programming across a machine with many, many cores, where the real benefits of parallelism arise.

The lesson is a simple but deep one: parallelism is not worth it for small inputs, but for large inputs, and with many cores, it is very much worth it. Since our datasets keep growing and our computers keep gaining cores, parallel programming matters more today than ever.

Looking ahead

Our parallel craps simulator gave us our first real speedup, and it did so by reusing the serial subroutines that we already had. We now turn toward a more scientific simulation. In the next code along, we will simulate a single particle diffusing through space as a random walk, which is the first step toward building a simulator of particles that diffuse and react with each other.