css_wp
css_wp

Reputation: 203

Plotting Pi using Monte Carlo Method

I can evaluate the value of pi using different data points by Python. But for each repeat I want to plot the scatter plot like this: enter image description here

My python code for finding pi using monte carlo method is :

from random import *
from math import sqrt
inside=0
n=10**6
for i in range(0,n):
    x=random()
    y=random()
    if sqrt(x*x+y*y)<=1:
        inside+=1
pi=4*inside/n
print (pi)

Upvotes: 2

Views: 13538

Answers (4)

Carlo Carandang
Carlo Carandang

Reputation: 225

Here is a variation on hiro protagonist's code, using random.uniform() to allow for random numbers between -1.0 and 1.0, allowing all the points to be plotted, and not just 1/4 of it (not the most elegant code, but it is spelled-out to learn the basics of the Monte Carlo Simulation):

import matplotlib.pyplot as plt
import random

inside = 0
n = 10**3

x_inside = []
y_inside = []
x_outside = []
y_outside = []

for _ in range(n):
    x = random.uniform(-1.0,1.0)
    y = random.uniform(-1.0,1.0)
    if x**2+y**2 <= 1:
        inside += 1
        x_inside.append(x)
        y_inside.append(y)
    else:
        x_outside.append(x)
        y_outside.append(y)

To estimate pi, the points in the circle correspond to the area of the circle enclosing it (pi*radius^2) and the total points correspond to the area of the square enclosing it (2*radius)^2. So this translates into:

(points in the circle)/(total points) = (pi*radius^2)/(2*radius)^2

Solving for pi, the equation becomes:

pi=4*(points in the circle)/(total points)

pi = 4*inside/n
print(pi)

Plot the points inside and outside the circle:

fig, ax = plt.subplots()
ax.set_aspect('equal')
ax.scatter(x_inside, y_inside, color='g', marker='s')
ax.scatter(x_outside, y_outside, color='r', marker='s')
fig.show()

Plot of points inside and outside the circle

Upvotes: 0

hiro protagonist
hiro protagonist

Reputation: 46839

building from your code, this may get you started:

import matplotlib.pyplot as plt

from random import random

inside = 0
n = 10**3

x_inside = []
y_inside = []
x_outside = []
y_outside = []

for _ in range(n):
    x = random()
    y = random()
    if x**2+y**2 <= 1:
        inside += 1
        x_inside.append(x)
        y_inside.append(y)
    else:
        x_outside.append(x)
        y_outside.append(y)

pi = 4*inside/n
print(pi)

fig, ax = plt.subplots()
ax.set_aspect('equal')
ax.scatter(x_inside, y_inside, color='g', marker='s')
ax.scatter(x_outside, y_outside, color='r', marker='s')
fig.show()

although i prefer this answer that uses numpy from the start...

Upvotes: 0

DYZ
DYZ

Reputation: 57033

To further elaborate Robbie's code:

import numpy as np
import matplotlib.pyplot as plt

n = 1000
xy = np.random.uniform(-1, 1, 2 * n).reshape((2, n))
in_marker = xy[0]**2 + xy[1]**2 <= 1
pi = np.sum(in_marker) / n * 4
in_xy = xy[:, in_marker]
out_xy = xy[:, ~in_marker]


fig, ax = plt.subplots(1)
ax.scatter(*in_xy,c='b')
ax.scatter(*out_xy,c='r')
ax.set_aspect('equal')
fig.show()

Upvotes: 3

Robbie
Robbie

Reputation: 4872

If you get errors about the backend use this:

import matplotlib as mp
mp.use('Tkagg')

Which will set the backend to TkAgg, which uses the Tkinter user interface toolkit.

import numpy as np
import matplotlib.pyplot as plt

n=1e3
x = 1-2*np.random.random(int(n))
y = 1-2.*np.random.random(int(n))
insideX,  insideY  = x[(x*x+y*y)<=1],y[(x*x+y*y)<=1]
outsideX, outsideY = x[(x*x+y*y)>1],y[(x*x+y*y)>1]

fig, ax = plt.subplots(1)
ax.scatter(insideX, insideY, c='b', alpha=0.8, edgecolor=None)
ax.scatter(outsideX, outsideY, c='r', alpha=0.8, edgecolor=None)
ax.set_aspect('equal')
fig.show()

enter image description here

Upvotes: 5

Related Questions