The principle of the code is as follows: we consider the unit square with diagonally opposite vertices at (0,0) and (1,1). We scatter points in that square at random, and count the number that fall into the positive quadrant of the unit circle. The ratio of this number to the total number of points approximates Pi/4.
main()
{
int i, j;
float pi, piEst, realI, realJ, xSquared, ySquared;
double x,y;
FILE *outdata;
pi = 4.0 * atan(1.0);
srand48(1.0);
if ((outdata=fopen("monty3","w")) == NULL)
{
printf("failed to open the write file\n");
exit(1);
}
j = 0;
for (i=1; i<= 1000000; i++)
{
x = drand48();
y = drand48();
xSquared = x * x;
ySquared = y * y;
if ((xSquared + ySquared) < 1.0)j++;
realI = i;
realJ = j;
piEst = 4.0 * realJ/realI;
if (i%1000==0)fprintf(outdata,"%d %f %f \n",i,pi,piEst);
}
}
You can copy this code to your own machine; to run it, you may need to replace the calls to srand48 and drand48 with a call to your local random-number generator.
To apply this method to an economics problem, you would replace x and y with inflation rate, selling price, or other variables. The current version of the code will give you a random distribution of values between 0 and 1, which will probably not match the possible values of the variables. If the actual range for selling price of a product is $4.00 to $6.00, the expression
price = 4.00 + 2 * drand48();
will give you a uniform distribution of prices over this range. But what if you know that the distribution is not uniform, as it probably isn't? This can be handled by approximating the probability function with a histogram, and mapping the interval (0,1) to the cumulative value of the histogram. Thus, for example, if there is a 10% probability that the price will be in the range ($4.00, 4.50), you can use the expression
if (0 < drand48() < 0.1) price = $4.25
and so on for the rest of the range.