# Random Numbers and Monte Carlo Simulations

Document Sample

```					Random Numbers and Monte Carlo Simulations

Monte Carlo calculations are a class of calculations that barely existed

E.g., a particle passing through a gas, several different outcomes are
possible with known probabilities. This could be formulated as a problem
in mathematical physics.

An alternative easier approach is to simulate the fate of a large number of
particles entering the detector using random numbers to determine what
happens to them and ﬁnally count how many are detected.
Random Numbers and Monte Carlo Simulations

Monte Carlo calculations are a class of calculations that barely existed

E.g., a particle passing through a gas, several different outcomes are
possible with known probabilities. This could be formulated as a problem
in mathematical physics.

An alternative easier approach is to simulate the fate of a large number of
particles entering the detector using random numbers to determine what
happens to them and ﬁnally count how many are detected.

Quantum Mechanics gives an intrinsecally probabilistic description of
physics phenomena. Monte Carlo can give a probabilistic simulations of
the physics!
Monte Carlo: why the name? A bit of history

Monte Carlo methods were known as ”statistical sampling”.

The name ”Monte Carlo” was due to physicists Ulam, Fermi, von
Neumann and Metropolis: the name is a reference to Monaco’s casino
which Ulam’s uncle would borrow money to gamble at !

The most famous early use was by Fermi in 1930, he used a random
method to calculate the properties of the newly-discovered neutron.

Monte Carlo methods were central to the simulations required for the
Manhattan Project, though limited by the computational tools at the time.

It was only after electronic computers were ﬁrst built (from 1945 on) that
Monte Carlo methods began to be studied in depth.

In the 1950s they were efﬁciently used at Los Alamos for early work
relating to the development of the hydrogen bomb.
The problem posed by the Conte de Buffon in 1777 provides a simple
illustration of a simulation.

“Consider parallel lines of separation L to be drawn on the ﬂoor. A
needle of length L is thrown ‘at random’; for what fraction of the
throws does it land crossing a line?”

The position of any thrown needle is then represented as a point in the
(x, φ) plane. If the needle is thrown ‘at random’ any value of (x, φ) is
equally likely, or the ‘density’ of points in the (x, φ) plane is uniform.

You can show the needle crosses a line if x < 1 L sin φ or L − x < 1 L sin φ.
2                    2
So the curves x = 1 L sin φ and L − x = 1 L sin φ divide the plane into areas
2                     2
where the needle does or does not cross a line.

x
L

hit

L/2          miss               miss

hit

0                                     φ
π/2          π
0

N
Fraction crossing a line is ratio of the ‘hit’ area to total = N hit
total
So the curves x = 1 L sin φ and L − x = 1 L sin φ divide the plane into areas
2                     2
where the needle does or does not cross a line.

x
L

hit

L/2          miss                               miss

hit

0                                                     φ
π/2                π
0

N
Fraction crossing a line is ratio of the ‘hit’ area to total = N hit
total
Rπ 1
2   0 2 L sin φ dφ      2
=
Lπ               π
So the curves x = 1 L sin φ and L − x = 1 L sin φ divide the plane into areas
2                     2
where the needle does or does not cross a line.

x
L

hit

L/2          miss                               miss

hit

0                                                     φ
π/2                π
0

N
Fraction crossing a line is ratio of the ‘hit’ area to total = N hit
total
Rπ 1
2   0 2 L sin φ dφ      2
=
Lπ               π

A statistical estimate of π can be obtained!
Simple computer simulation: pairs of numbers x and φ are generated
using a pseudo-random number generator covering the ranges 0 < x < 1
2
and 0 < φ < π.

(See lecture notes for a few examples of ‘random’ number generation
algorithms and exercises with these.)
1
The number of pairs that satisfy x < 2 sin φ is counted: i.e., a simple hit or
miss Monte Carlo. Ten repetitions each involving 104 trials gave the
0.6364   0.6361     0.6271    0.6370    0.6343
following results
0.6385   0.6443     0.6339    0.6352    0.6324

Average is 0.6355 ± 0.0014 compared with theoretical result 0.63662.

High accuracy is impossible. This used 105 random numbers.

It would take a 100 times as many to get another ﬁgure in the average.
√
It can be shown that the error of N trials is 2.37/ N .
In this example the simulation is estimating an area.

Most Monte Carlo calculations can be regarded as evaluating an area,
volume or ‘hyper-volume’ by counting whether random points are inside or
outside.

They are, in some sense, performing an integration, although this may not
be explicit.

Later we turn this argument round to evaluate integrals.

One can also devise adaptive Monte Carlo algorithms to perform
integrations of rapidly varying integrand functions, using the so-called
importance sampling technique (more on this later).
Computer Exercises

Extension of the Buffon problem
Program computes the probability of the needle crossing a line.

Python Code
n=100000
for k in range(1,6):
c=0
for i in range(0,n):
x=uniform(1)
phi=pi*uniform(1)
sx=0.5*sin(phi)
if ((x<sx) or (1-x<sx)):
c+=1

print k, float(c)/n, 2/pi
Simple extension of the Buffon needle problem: needle falling on a ﬂoor
with square tiles.
The Problem of Random Flights.

One of a variety of problems known as ‘drunkard’s walks’ (proposed by
Pearson in 1906).

A drunkard starts from the origin and walks a distance L in a straight line;
he then turns through a random angle and walks a distance L in a second
straight line. He repeats this process n times.

What is the probability that he is a distance between r and r + δr from the
starting point?

Exercise involves a computer simulation of the three-dimensional
generalisation of this problem.
Generation of random directions.
In two dimensions generate a random angle φ between 0 and 2π.

(x, y) → (x + L cos φ, y + L sin φ).

In three dimensions use polar angles θ and φ

θ

φ

(x, y, z) → (x + L sin θ cos φ, y + L sin θ sin φ, z + L cos θ).
Complication: in order to get an equal probability of all directions it is cos θ
not θ that must have a uniform distribution.
Python code to compute distance for one drunkard

def distance(n):
"Computes distance moved after n random steps"
x=y=z=0                            # Start at origin
for i in range(0,n):             # Count steps
c=2.0*uniform(1)-1.0       # c=cos(theta)
s=sqrt(1.0-c*c)             # s=sin(theta)
phi=2.0*pi*uniform(1)      # phi is random azimuth
x+=s*cos(phi)                # calculate new position
y+=s*sin(phi)
z+=c

return sqrt(x*x+y*y+z*z)    # final distance from origin
Displaying the results as a histogram
Four examples, with 15 bins 500 walks, 15 bins 10000 walks 50 bins 10000 walks
50 bins 100000 walks

1400
60                                                  1200

1000
40                                                  800

600
20                                                  400

200

0   2        4   6    8    10   12   14               0   2        4   6    8    10   12   14

4000
400

3000
300

200                                                  2000

100                                                  1000

0       10       20       30    40        50          0       10       20       30    40        50
nbins=20
jmax=nbins-1                       # Size of array of distances
n=3                                  # Number of steps
nwalk=500
deltaR=(1.1*n)/jmax
m=[0 for i in range(0,jmax)]     # Zero array
for k in range(0,nwalk):         # Count walks
R=distance(n)                # Determine distance from start
j=int(floor(R/deltaR+0.5))   # Convert to bin number
if j<jmax: m[j]+=1            # Increment count of distance bin

histogram(m)
Remember that the function P(r) is presumably a smooth curve.

For keen students only:

The NFB and C programs contain a deliberate mistake!
This illustrates the difﬁculties of debugging Monte Carlo programs so
perfectly that it is worth trying the exercise in these languages as well.
It is not possible to make the mistake in Python.

```
DOCUMENT INFO
Shared By:
Categories:
Stats:
 views: 8 posted: 6/1/2010 language: English pages: 18
How are you planning on using Docstoc?