Introduction to Parallel Computing

For (i), I showed slides numbers 2, 3, 4 as information about myself, slides 13, 14, 17, 20, 21 about a prototypical PDE problem, slides 60-67 to introduce how a parallel computing cluster looks like and what the issues are with distributed-memory clusters, and slides 78, 79, 81, 82, 83 to motivate the following conlusion: The most important purpose of parallel computing is to pool the memory of several compute nodes to allow the solution of larger problems than possible on a serial node. The second most important purpose is to speed up calculations to make calculations feasible within reasonable amoung of time.

For (ii), I used my iPad2 with NoteTakerHD to write and project hand-written lectures, which will be posted as PDF transcripts here, I annotated the slides with red color (slides posted without annotations here), and used iSSH to log in briefly to the cluster hrz-cs400.hrz.uni-kassel.de that we will use in this class.

Problem 1 introduces several things: We logged in to the cluster hrz-cs400.hrz.uni-kassel.de; when doing this and the following, we cross-checked against information on the webpage of the IT Servicezentrum. I demonstrated basic Unix/Linux commands such as pwd, ls -l, mkdir, cd, rm. For Problem 1, I showed how to program in the language C, by writing a "Hello, world!" program from scratch; see attached file hello_serial.c. To do this, we had to learn the basics of the vi editor, namely to enter it, simply say "vi filename" which puts you in command mode of the editor; type "i" to enter input mode; then you can enter text. Use ESC to escape the input mode and return to command mode. In this mode, you can use the keyboard keys "h", "j", "k", "l" (or apparently the arrow keys) to go left, up, down, and right, respectively. Finally, to quit the editor and save ('write to file') the result, type ":wq" and hit return.
We used the compile command for parallel MPI code, namely mpicc to compile hello_serial.c. The default name for the resulting executable is a.out; to give it a more meaningful name, use the "-o" option, so we used
mpicc hello_serial.c -o hello_serial
You can run it by saying "./hello_serial".

Problem 3 asks you to download two posted files. The easiest method at this point turned out to be to click on the links of the posted file to display them in the web browser, highlight all contents, and copy-and-paste it into the file by the desired name.
I then showed what should happen when you run this in parallel by submitting the job to the scheduler user
qsub mpich2_pgi.qsub
You should follow carefully the steps in Problem 3, because this should not work the first time around.

Regarding use of this cluster, we used qsub to submit a job, then qstat to see the status of all currently running or queued jobs (recommended options are -a and -n) and the Ganglia software, which you access by pointing a web browser to
hrz-cs409.hrz.uni-kassel.de/ganglia
(Notice the 409 [not 400] in this URL!).
Then we read the hello_parallel.c code and realized it is nearly the most basic MPI command, only containing MPI_Init and MPI_Finalize at beginning and end, and MPI_Comm_size and MPI_Comm_rank to determine the size of the given communicator and the rank of the process that is executing the code. See Chapter 3 in Pacheco and see the man pages for each of these commands. Notice that we delayed explaining some details of C programming (command-line arguments argc and argv to main program and MPI_Init function, role of ampersand & in &np, pointers as address locations of variables, ambiguity of pointer vs. array, string really being array of characters, and more); please understand that I am applying a tactic that teaches you the necessary C skills by observing supplied code. We will return to more explanations later.

We discussed some ideas for how to improve this code. Two specific issues right at the beginning are that the compilation with mpicc gives a warning (this should be fixed) and the output is poor. I pointed out that the output raises the question what mathematical problem we are solving anyhow; looking through the code for a, b, and f(x) gives the information, and this allows us to obtain the true integral value I=1/3. Using this, I wanted to teach the reasonable habit to output all available information such as, at least, the number of subintervals n, the step size h, true value I, the approximation In (which is the variable total in Pacheco's code), and the error I-In (coded as I-total); in fact, since the latter one should have magnitude h^2, you may also want to print out h*h, as well. I spent some time to demonstrate some pitfalls of C here, such as coding I=1/3 just like that, which gives integer division 1/3=0 and hence float I=0.0. Moreover, this outputs as I=0.000000 for format conversion %f, which is what Pacheco uses in several places. Once we fixed the computation as I as I=1.0/3.0, we still a strange thing like error=I-total=0.000000, which is odd, since there should be some error on the order of h^2=1e-7. The issue is that the output format %f simply does not show enough digits. I suggested highly to always use %24.16e, which is in exponential format, normalized to one non-zero decimal digit before the decimal point, and with 16 decimal digits after the decimal point, then multiplication with powers of 10 as indicated by notation e-12 at the end. Single-precision numbers have about 8 significant decimal digits of accuracy, so the above is too many clearly (%15.8e would suffice). But for very large n values, we will have reason to switch to using double-precision numbers anyways, so I anticipate this with the %24.16e. The 24 indicates the fieldwidth of the output field. We need to allow for 8 additional characters here, since double-precision numbers range from 1e-308 to 1e+308 approximately, hence there can be up to 3 digits necessary in the exponent of 10 as in e-123.

Then we introduced MPI_Gather, MPI_Scatter, and MPI_Allgather from Sec. 5.7, but using strictly vectors x and l_x as examples (not matrices like in Pacheco).

Finally, I lectured a little on the Unix make utility with an example of C code split into two files file1.c and file2.c. The make utility makes this more efficient then by only re-compiling one of these files, if only one of them has changed since the last compilation. These date/time stamp comparisons are the contents of the conditional lines (those with a colon ":") in the Makefile. To show this, I logged in to the cluster and showed the files supplied for HW 4 and 5, which include a Makefile and code split into files main.c and memory.c. I demonstrated "make" as well as "make clean". We will look at this more as part of the lab.

MISTAKE ALERT: The file mpich2_pgi_power.qsub may have a mistake in it (depending on which version I posted): The mpijob line should read
mpijob="./power 8192 1.0e-12 50"
with reference to executable power in the current directory by "./power"; it may say "../../power" in the supplied file.

We covered the IEEE Standard 754 for binary floating-point arithmetic. More specifically, we covered its portion on number representation. It is implied that operations among such numbers also result in accuracies, which cause round-off errors of the same type as that caused by the representation itself.

• Use the webpage www.netlib.org to "Browse" to "blas", then (in the appropriate list of BLAS1, BLAS2, BLAS3) find the desired function. Example: A matrix-vector product is a BLAS2 operation. Looking through the list of BLAS2 functions, we find dgemv, the double-precision version of the matrix-vector product.
• For real-life purposes, where we just want to use the function, we will not usually read the code. Rather, I would like you to know that there are man pages, so try "man functioname" on the hrz-cs400 cluster. Example: From the documentation in "man dgemv", we see that the function actually implements y := alpha*A*x + beta*y with matrix A and y := alpha*A'*x + beta* with the transpose of matrix A. This choice is controlled by the first input argument, which is either 'N' (= 'no-transpose') or 'T' (= 'transpose'). One very important thing to note is that the output y is also an input. That is very typical of BLAS functions, many of which are designed to accumulate results. So, to avoid unintended accumulation, we should set beta=0.0.
• However, we might find that this documents the Fortran interface. It does give the meaning of the arguments, but we also need to know the C interface of the function. To this end, I want you know that the file acml.h has useful comments at the beginning and that you can search it for the function name to see its C interface. To find acml.h, use the command "locate acml.h", then use any editor to open and read the file; if in doubt, look at the 64-bit version of it (the one with linux86-64 in the directory). Example: I said "vi /opt/pgi-7.0-5/linux86-64/7.0-5/include/acml.h" to edit the file acml.h and searched it for "dgemv". This confirms the exact C syntax of the function interface. While many of the arguments were as expected int for the integers and double for alpha and beta, this confirms that all array are defined as double*, even the matrix A; this is the one-dimensional storage of A that we talked about. Also, we can see here that the first argument is a char; notice it is a char itself, not a char*; this is exactly what a single character 'N' or 'T' is in C.

• Focus on integers when debugging; design test output of integers and print them before using them (even comment out their use in MPI commands to avoid hanging code until after testing them).
• As part of the printf in a parallel code, print out the id number as in printf("[%3d] ...\n", id, ...); this is always a good habit, so you know where a message comes from.
• Initialize all defined variables to something, for instance some odd value like -69, so that you can recognize in output if a variable was never set at all.
• Parallel code might not run if np=1. So, put id(np>1) as condition around MPI calls; to make the code correct mathematically, you often need an else case where date is set. For example, inside a parallel dot product you might have computed the local dot product in variable l_dot; inside the if(np>1) you would call an MPI_Allreduce to obtain the global variable dot; the else case for a single-process (serial) run needs to be dot=l_dot then.