Kodeclik Blog
How to make a Monte Carlo Simulator in Python
Many problems—like estimating chances, measuring areas, or predicting outcomes—are too complex or tedious to solve with exact math alone. That’s where the power of computers and randomness comes in.
What is Monte Carlo Simulation?
Monte Carlo simulation is a mathematical technique that uses random sampling to solve problems that might be deterministic in principle but are difficult to analyze directly.
Monte Carlo simulation is named after the famous casino city. We will show how it can be especially valuable for demonstrating concepts like probability and statistical estimation, as it makes abstract ideas tangible and visual.
Monte Carlo simulation lets us experiment with thousands or millions of random trials, making it possible to estimate solutions and see patterns emerge, even when the underlying system is complicated. Using Python, this process becomes accessible and engaging, allowing anyone to visualize probability and uncertainty in action before learning the theory behind it.
Example 1: Monte Carlo Simulation to find the area of a circle
Here is a simple Monte Carlo simulation in Python to estimate the area of a circle. This approach randomly generates points within a square and checks how many fall inside the circle inscribed within the square. The ratio gives an estimate of the circle’s area.
import numpy as np
def estimate_circle_area(num_samples, radius=1):
x = np.random.uniform(-radius, radius, num_samples)
y = np.random.uniform(-radius, radius, num_samples)
inside_circle = (x**2 + y**2) <= radius**2
estimated_area = (np.sum(inside_circle) / num_samples) * (2 * radius) ** 2
print(f"{num_samples}, {estimated_area:.5f}")
if __name__ == "__main__":
for samples in [10, 100, 1000, 10000, 100000, 1000000, 10000000]:
estimate_circle_area(samples)In the above program, the estimate_circle_area function generates a given number of random points within a square that fully contains the circle. It then counts how many points fall inside the circle by checking if each point's distance from the origin is less than or equal to the radius.
The estimated area is calculated as the proportion of points inside the circle multiplied by the area of the square.
In the main part of the program, this function is called with increasing numbers of random samples (from 10 up to 1,000,000).
The output will be:
100, 3.32000
1000, 3.18000
10000, 3.16400
100000, 3.14608
1000000, 3.13944
10000000, 3.14140Note that the accuracy gradually improves as we use more samples. This demonstrates that as the sample size grows, the estimated area converges toward the actual value.
Example 2: Monte Carlo Simulation to estimate heads
A second fun and simple example of a Monte Carlo simulation is estimating the probability of getting heads when flipping coins. This project uses the same basic idea of random sampling but where each run or simulation denotes the act of flipping a coin.
import numpy as np
def simulate_coin_flips(num_flips):
flips = np.random.choice(['Heads', 'Tails'], size=num_flips)
heads_count = np.sum(flips == 'Heads')
probability_heads = heads_count / num_flips
print(f"Flips: {num_flips}, Estimated Probability of Heads: {probability_heads:.3f}")
for flips in [10, 100, 1000, 10000]:
simulate_coin_flips(flips)This program randomly “flips” a coin many times, counts how often it lands on heads, and calculates the estimated probability. The function simulate_coin_flips takes a number of coin flips as input, then uses numpy’s random.choice() to randomly generate that many outcomes, each either 'Heads' or 'Tails'. It counts how many times 'Heads' appears using numpy, and divides this by the total number of flips to estimate the probability of getting heads. The program then prints both the number of flips and the estimated probability (rounded to three decimals).
In the main loop, this function is called with increasing numbers of flips (10, 100, 1000, 10000) so students can see how the estimated probability gets closer to 0.5—which is the expected probability for a fair coin—as the sample size grows.
The true probability is 0.5, which you will observe as the number of flips increases. For instance, in our case, the output would be something like:
Flips: 10, Estimated Probability of Heads: 0.200
Flips: 100, Estimated Probability of Heads: 0.490
Flips: 1000, Estimated Probability of Heads: 0.476
Flips: 10000, Estimated Probability of Heads: 0.500Example 3: Monte Carlo Simulation to estimate doubles
A third example of a Monte Carlo simulation is estimating the probability of rolling doubles with two dice. To determine the probability, Imagine that you roll two dice many times (or simulate it on a computer), and then we simply count how often both die show the same number, and then divide that count by the total number of rolls to estimate the probability of rolling doubles.
import numpy as np
def estimate_doubles(num_rolls):
die1 = np.random.randint(1, 7, num_rolls)
die2 = np.random.randint(1, 7, num_rolls)
doubles = die1 == die2
probability = np.sum(doubles) / num_rolls
print(f"Rolls: {num_rolls}, Estimated Probability of Doubles: {probability:.4f}")
for rolls in [10, 100, 1000, 10000]:
estimate_doubles(rolls)Here, the function estimate_doubles takes a number of rolls as input and, for each roll, randomly generates two numbers between 1 and 6 (simulating two dice). It counts how many times both dice show the same number (a double) by comparing the two arrays. The estimated probability is then calculated as the proportion of double rolls out of the total rolls, and the result is printed with the sample size and the probability.
In the main loop, the function runs for increasing sample sizes (10, 100, 1000, 10000), letting students see how the estimated probability (which should be close to 1/6 or about 0.1667 for fair dice) gets more accurate as more rolls are simulated.
In our case, the program’s output will be (something like):
Rolls: 10, Estimated Probability of Doubles: 0.1000
Rolls: 100, Estimated Probability of Doubles: 0.1800
Rolls: 1000, Estimated Probability of Doubles: 0.1620
Rolls: 10000, Estimated Probability of Doubles: 0.1631Summary
Monte Carlo simulations use random sampling to estimate outcomes or probabilities for complex processes where analytical solutions may be difficult or impossible. In the examples above, one program estimates the area of a circle by randomly generating points inside a square and checking how many fall within the circle; another estimates the probability of getting heads in repeated coin flips; and a third simulates rolling two dice to find the chance of rolling doubles.
Each simulation demonstrates that with more samples, the estimated probability or measurement becomes more precise, showing that randomness can be quite useful.
Enjoy this blogpost? Want to learn Python with us? Sign up for 1:1 or small group classes.
