Calculating pi with Monte Carlo integration

[mjordan2nd]

Hello!

I am taking a computational physics class this semester and just got out of a test. One of the questions provided us with 10 random numbers and asked us to approximate pi. The method I came up with was realizing that

[tex] \int_0^1 \sqrt{1-x^2} dx = \frac{\pi}{4},[/tex]

and then performing this integral numerically. To perform the integral numerically, I plugged in the random values in place of x and evaluated the integrand. I then summed up the integrand values for each of the 10 random numbers and divided this value by 10 to give me pi/4. This method game me a value of pi of around 3.49. I figured this was good enough since I only had 10 values. However, I decided to come back and actually code this into the computer. Using java, I wrote the following code:

import java.io.*;
import java.util.*;
import java.math.*;

public class Test
{
	public static void main(String [] args)
	{
		double sum = 0;
		for(int i = 0; i< 100000; i++)
		{
			sum+=Math.sqrt(1-(Math.random() * Math.random()));
		}
		sum = (sum / 100000) * 4;
		System.out.println(sum);
	}
}

However, even with a 100000 random numbers I'm still only getting values of around 3.415. There seems to be some kind of systematic error in my logic, however I can't pinpoint what it is. Any help clearing up this matter would be appreciated.

Thanks.

 

[mjordan2nd]

Never mind, I found my error. The math.random() * math.random() should read Math.pow(Math.random(),2). This fixed it, and gave me a value of 3.14158

 

[Mark44]

First, I don't see any point in generating two random numbers per iteration. Just generate one and square it.

Second, I generally try to avoid using pow() when all I'm doing is squaring or cubing a number. After all, x*x produces the same result as pow(x, 2) without the overhead of a function call.

Third, it's not clear to me what you are doing, but it doesn't seem like a Monte Carlo simulation. In a Monte Carlo simulation of the value of ##\pi## you would be simulating throwing darts at a quarter circle of radius 1 that sits inside a unit square. In the simulation you generate x and y values at random (in the range 0 to 1). If the point lies inside the circle, count it. If it lies outside the circle, don't count it. With a large number of darts, the ratio of points inside the circle to the total number of darts thrown should be close to ##\pi/4##.

Whatever it is that your code is doing, here's how I would modify your code, with an eye to speeding it up.

for(int i = 0; i< 100000; i++)
{
   x = Math.random();
   sum += Math.sqrt(1.0 - x*x); 
}

 

[D H]

Whatever it is that your code is doing ...

What his code is doing is using random sampling to estimate the expected value of the random variable Y, where [itex]Y = \sqrt{1-X^2}[/itex] and X is U(0,1). Since the expected value is an integral (in this case, [itex]E(Y) = \int_0^1 \sqrt{1-x^2}\,dx[/itex]), I would say this qualifies as a monte carlo integration technique -- and it is more efficient than the rejection technique you described.

 

[Mark44]

What his code is doing is using random sampling to estimate the expected value of the random variable Y, where [itex]Y = \sqrt{1-X^2}[/itex] and X is U(0,1). Since the expected value is an integral (in this case, [itex]E(Y) = \int_0^1 \sqrt{1-x^2}\,dx[/itex]), I would say this qualifies as a monte carlo integration technique -- and it is more efficient than the rejection technique you described.

OK, good to know. It bothered me that the result was not too far off from ##\pi##.