Many computational methods in science rely on the generation of random numbers. We will briefly discuss how to generate `pseudo-random numbers' on a computer, and how to use such a sequence of numbers to perform an integration.
The essential issue in the generation of a sequence of random numbers is to prevent
corrolations among parts of the sequence.
Roughly speaking, most generators use a linear congruent method to generate
a sequence of integers
, where m is a large integer - perhaps the maximum
integer represented on the machine. Dividing by m will produce a
random number between 0 and 1. If - and this is a big if - a and c are
chosen well, this sequence will have no corrolations, up to order m.
But if a and c are badly chosen, there will be corrolations far earlier.
Consider further the use of say a triple of random numbers to represent
a point in 3-space. Pairs of points define a plane, and there are only
distinct such planes - thus a much worse
sampling of space than one might think.
Numerical Recipes gives three routines for generating a random sequence. One way shuffles the numbers generated by a generator resident in most machines (ususally called ran or rand). The other two are full-fledged generators themselves. BE WARY - the codes in the first edition are not correct.
Now consider the difficulties in determining the area of a complicated shaped region in the plane. Traditional integration methods might run into trouble - even though you might be able to tell when a given point is inside or outside the region, writing down formulaes to perform the multiple integrals could be difficult. However, if you could imbed the region in a regularly shaped region Monte Carlo methods would give you the desired area. Choose a sequence of pairs of random numbers, with the pair thought of as a point in the plane, and test to see if the point is inside or outside the complicated region. The proportion of points `inside' is the proportion of the total area of the regular region occupied by the complicated region.
As a simle example, consider the problem of finding the area of
a quarter-circle of unit radius, embedded in the unit square [0,1]X[0,1].
Choose a pair of random numbers, and check to see if the point is
less than one from the origin. If yes, add one to a running tally of
`points inside'. The fraction of points inside is an estimate of
.
Unfortunately, this method has a very low order of accuracy -
the error in the estimate decreases only like
,
where N is the number of trial points.
Run the following code and determine the number of sample points needed to obtain 4 digits of accuracy.