# A Monte Carlo Algorithm for Computing pi;

## To illustrate running a job under Condor, you will use a Monte Carlo method for computing π

This is an example of an embarrassingly parallel problem that can easily be run on the Campus

Condor Pool.

To illustrate running a job under Condor, you will use a Monte Carlo method for computing . This is an example of an embarrassingly parallel problem that can easily be run on the Campus Condor Pool.

##### Monte Carlo Algorithm

Consider the familiar example of finding the area of a circle. Figure~?? shows a circle of radius inscribed within a square. The area of the circle is , and the area of the square is . The ratio of the area of the circle to the area of the square is:

If we could compute the ratio , then we could multiple it by four to obtain the value of . One particularly simple way to do this is to pick regularly spaced lattice points in the square and count how many of them lie inside the circle. Suppose that we pick the points defined by:

Then we find that there are points inside the circle and points outside the circle. The ratio is then and so our estimate of is . This is obviously not a very good approximation.

Rather than picking a regularly spaced lattice, we can use Monte Carlo to compute . Monte Carlo methods can be thought of as statistical simulation methods that utilize a sequences of random numbers to perform the simulation. In a typical Monte Carlo process, we compute the number of points in a set that lie inside a region . The ratio of the number of points that fall inside to the total number of points tried is equal to the ratio of the two areas (or volume in three dimensions). The accuracy of the ratio depends on the number of points used, with more points leading to a more accurate value. Monte Carlo methods like this are a good way of computing the values of complicated higher-dimensional integrals.

A simple Monte Carlo simulation to approximate the value of is to randomly select points in the unit square and determine the ratio , where is number of points that satisfy .

Write a C program that computes using this Monte Carlo method. Your code should take two command line arguments: the first should specify an integer number of points to use in the Monte Carlo, and the second should be a string which specifies the name of a file that your program should write its output to. Your program should print out to *the specified file* (not to the screen) columns containing the number of points inside the circle , the total number of points , the ratio , the value of , and the absolute error in the value of . To compute the error, compare your computed value of to the C standard library value stored in the constant **M_PI**. You should use double precision numbers for these calculations and make sure that the columns in the file contain sufficient precision using the appropriate format statement in your **fprintf**.

To compute the random numbers that you will need in this Monte Carlo, you should use the C standard library functions **rand()** and **srand()**. A digital computer cannot generate truly random numbers unless it is connected to a physically random source, however the C library provides a way of generating pseudo-random numbers that will be sufficient for our purpose.

Random-number generation is a complex topic. The book Numerical Recipes in C provides an excellent discussion of practical random number generation issues in Chapter 7. For the this exercise, you can generate a random number between zero and one with the following code

x = (double) rand() / (double) RAND_MAX;

You will need to **#include <stdlib.h>** to use **rand()**. Type

man 3 rand

for more information on the function.

Pseudo random number generators need a *seed* to start the random number generation. If you do not specify a seed to **rand()** it will use the number 1 and generate the same sequence of “random” numbers each time. Sometimes this is desirable, but not in our Monte Carlo simulation. Before you call **rand()** for the first time, you should call **srand()** to seed the random number generator. **srand()** takes an integer as its only argument. A different integer will generate a different sequence of “random'” numbers.

There are a variety of ways to seed a pseudo random number. The best way is to read an integer from the system file **/dev/urandom** using the usual **fopen()**, **fscanf()**, and **fclose()** calls.

**/dev/urandom** is a special system file that looks at read-world entropy (time between keystrokes, packets on the network, etc.) to generate random data. If you cannot use **/dev/urandom**, then you can use the system call **time(NULL)** to return the time in seconds since 00:00:00 UTC, January 1, 1970. If you use **time()** you will need to **#include <time.h>**. Make sure your program uses one of these methods to seed the random number generator in your Monte Carlo code.

Use the **time** command to time how long your program takes to run for , and points (use the value of **user** time, which gives the amount of time your code spent on the CPU performing calculations). Make a plot that shows the error in your computed value of as a function of the number of points. Your report should state how quickly your code is converging on the value of (i.e. what is the order of convergence as a function of )? A log-log plot will make it easier to see this.

##### A Large-Scale Monte Carlo on Condor

Computing a very accurate value of in the way described above will take some time, as the convergence as a function of is quite slow. When we encounter a problem such as this, high throughput computing is a good way of solving the computational bottleneck.

First recompile your Monte Carlo code, but add the commands **-static** and **-static-libgcc** to your compile line. For example, if your program is called **computepi.c** then you should compile it with

gcc -static -static-libgcc -Wall -o computepi computepi.c -lm

This makes your code more portable for running on the campus Condor pool. Check that this worked by running the command **dd** on your executable. If it prints

not a dynamic executable

then your code compile correctly. If it prints something like

linux-vdso.so.1 => (0x00007fff8edff000)

libc.so.6 => /lib64/libc.so.6 (0x00000038cba00000)

/lib64/ld-linux-x86-64.so.2 (0x00000038cb600000)

then you need to re-compile with the static options. Use **scp** to copy your executable from

**phy607.phy.syr.edu** to **its-condor-submit.syr.edu**

Write a Condor submit file that executes your program on the Campus condor pool for points. To specify an argument to your program, you need to add the line with the \verb|arguments| keyword in your submit file. For example, to run the program **computepi** with 100 points and write to the file **pi_output.txt**, the Condor submit file would begin

universe = vanilla

executable = computepi

arguments = 100 pi_output.txt

Don’t forget you need the rest of the submit file described previously to get your job to run. I suggest copying and modifying **hello.sub** for this purpose. Submit this job and check that it works. If you are having trouble running a Condor job, see me.

Using your timing data from the previous exercise, estimate how many points it would take for your program to execute for approximately two hours. Modify the **arguments** line in your Condor submit file to use this value of . To parallelize the task, use the **queue** line in your submit file to submit 1000 jobs to the pool. We now have a problem: if all of your jobs write to the file **output.txt** they will overwrite each others results.

Fortunately, we can use the Condor (process) variables to ensure that each job writes to a different file. To do this, your Condor submit file should read

arguments = 100 pi_output.(process).out

Once you have done this, submit the jobs to the Condor pool. This may take some time to run, so leave yourself plenty of time to do this.

Write a second program that reads in the values of and from the output files of **computepi** and used *all* of the measures values to compute . How accurate is your value now? Add a point to your graph showing the accuracy as a function of . Is it consistent with your predictions? How long would it take you to compute to an accuracy of , and using this method?

Your report should include all your C code and your Condor submit files.