Monte Carlo Method in Python

montecarlo, python

Solution

If you're doing any kind of heavy duty numerical calculation, considering learning `numpy`. Your problem is essentially a one-linear with a numpy setup:

import numpy as np

N   = 10000
pts = np.random.random((N,2))

# Select the points according to your condition
idx = (pts**2).sum(axis=1)  < 1.0
print pts[idx], idx.sum()

Giving:

[[ 0.61255615  0.44319463]
 [ 0.48214768  0.69960483]
 [ 0.04735956  0.18509277]
 ..., 
 [ 0.37543094  0.2858077 ]
 [ 0.43304577  0.45903071]
 [ 0.30838206  0.45977162]], 7854

The last number is count of the number of events that counted, i.e. the count of the points whose radius is less than one.

Problem

I've been attempting to use Python to create a script that lets me generate large numbers of points for use in the Monte Carlo method to calculate an estimate to Pi. The script I have so far is this: ``` import math import random random.seed() n = 10000 for i in range(n): x = random.random() y = random.random() z = (x,y) if x**2+y**2 <= 1: print z else: del z ``` So far, I am able to generate all of the points I need, but what I would like to get is the number of points that are produced when running the script for use in a later calculation. I'm not looking for incredibly precise results, just a good enough estimate. Any suggestions would be greatly appreciated.

Original source