Programming Lab #3: Parallel Gaussian Elimination

Due: Wednesday, November 20th, 2013 at 6pm


Changelog

Date Description
30-Oct-13 Original posting
9-Nov-13 Added language to permit verify.cpp to be modified
9-Nov-13 Added clarification about fastest and slowest varying dimensions


In this assignment you'll explore multithreaded implementations of Gaussian Elimination using C++11 threads. The procedure for obtaining the starter code is similar to the previous assignments. For instructions on how to obtain starter code, visit the following web page:

http://www-cse.ucsd.edu/classes/fa13/cse160-a/HW/A3/StarterCode.html

Operation of the Starter code

The organization of the code is described in a README file. Be sure to look over this file so you are familiar with the various code modules. To ensure that your program is correctly graded, do not change any of the provided code, except for the following two files where lu.cpp, is the driver program elim.cpp performs Gaussian Elimination with partial pivoting, and verify.cpp verifies the result.

lu.cp      elim.cpp       verify.cpp

The code takes 4 command line options: the input size (-n), the number of threads (-t), and a matrix input code (-h or -l)

./lu  [-n < matrix size >] [-t < # threads >] [-[h|l]]

If you don't specify the input code, you get a matrix which has the digit two (2) in all diagonal elements. This corresponds to the matrix 2*I, where I is the identity matrix. This matrix can help you catch simple errors, since all elements are already zero below the main diagonal. However, it will obviously miss certain errors, so we have provided more interesting test cases. The -l option gives you a more useful test case, though it does not require partial pivoting. Use this matrix to test everything except partial pivoting. For this input, the -n parameter must have an integral square root, i.e. n=1024 has an integer square root (32), but n=2048 does not (45.24, not an integer).

Testing your partial pivoting code

The next case tests partial pivoting. The matrix is called the Hadamard matrix, named after the French Mathematician Hadamard. A Hadamard matrix has entries that are +1 or -1 only, and N must be an integral power of two, i.e. N=512 = 29, 1024=210, and so on. (The code will report an error if you us a non-integral power of two, such as 1536). Here is a Hadamard matrix of order 4 (-h -n 4).

     1     1     1     1
     1     1    -1    -1
     1    -1     1    -1
     1    -1    -1     1

You can use this simple matrix when debugging the partial pivoting code. If you don't turn on partial pivoting, the provided code will report an error: a value called Not a Number (NaN), which gets outputted as the value nan (We will talk about NaN later in the course).

$ ./lu -h
N = 16, nt = 1
Using a 'Hadamard' Matrix
NOT Using partial pivoting
Non-zero residual = -nan
N=16: 0.000003 seconds (0.954204 gflops)

    N, NT, tP, gflops, partial Pivoting? verified?, mtxCode
    @ 16 1 0.000003 0.954204 N N H

Command line arguments are parsed in cmdLine.cpp.

The Assignment

Parallelize the code with C++11 threads. Most your effort will focus on parallelizing the Gaussian elimination routine contained in elim.cpp, but you'll also need to modify code in lu.cpp. No changes need be made to any other files, though in the course of developing your code, you may decide to make temporary modifications. But be sure to restore those files before you do final measurement and testing, using 'svn update.' (You may need to delete the modfied code or change its name before doing an update). When grading your turnin, we will use the clean copies of the Provided files to build your code.

The code stores the matrices in row major order such that entries with the same first subscript lay in the same row, and row members with consecutive second subscripts follow one another contiguously in memory. Thus, we say that the column dimension varies more quickly in memory than the row dimension.

Implement 1-dimensional partitioning and split only the slowest varying dimension in memory. Why do we partition the slowest varying dimension rather than the fastest varying dimension? It is not necessary to use cyclic decomposition. Your code must correctly handle the case when the number of threads doesn't evenly divide the size of the matrix. In particular, the amount of work assigned to the most heavily-loaded and most lightly-loaded threads should differ by not more than one column of the matrix

You may assume that

n >  2 × the number of threads

Try to make the code go as fast as you can. In the spirit of friendly competition, post your performance results to the A3 Moodle forum.

Testing your code

Test your code by comparing against the provided implementation. The provided code includes a routine to verify your results, though you'll need to modify it to run on multiple threads.

The code verifies your results by computing a measure called the residual. This residual will generally be a small positive number ε which we have set to 106 in the provided code. When your code verifies correctly, you'll see the following message

LU decomposition verified OK [eps < 1e-06]

If there is an error, you will get a message like this

Non-zero residual = 7.63808

The residual relates to a solution to linear system of equations. If the values are tiny, but larger than 10-6--for example on the order of 10-5-- you have a subtle error. The above error is a not-so-subtle error.

To help with debugging, the provided code will output the input matrix and the outcome, if sufficiently small. A small matrix is currently set at N<10 (see lu.c) and you can change this if you like. Note that the entries below the main diagonal are the "L" matrix in the LU decomposition, but that the diagonal entries have been overwritten by the "U" matrix. In the actual L matrix, the diagonal entries would all be one (1).

Test against various combinations of the -n and -t command line parameters, along with the various test matrices specified with the -h and -l flags. Verify that your code produces correct results that are independent of the number of threads and that it runs correctly whether or not -nt divides -n.

While the application will not likely attain perfect scaling (100% efficiency) it should deliver a reasonable speedup to at least 4 cores. Since the running time of LU decomposition varies as O(N3), use care when scaling the -n parameter. Do you notice any change in scaling as you increase n?

Performance Measurement

Use the provided timer to take timings. Run the program at least 3 times -- for each input size or any other variable that affects the running time -- until you see a consistent timing emerge, that is, where a majority of the times agree to within 10%. If you do not observe a consensus, try doubling the number of measurements. Exclude outlying timings (but note how large they are).

Out of courtesy to other users, be sure that your jobs-- including those comprising several runs-- do not last for more than a few minutes. You can always repeat your runs again if you aren't getting consistent results. Keeping your runs short will increase the likelihood that your results will be consistent. There is no need for your program to run for more than about 10 seconds on 1 core. This corresponds to an input size in the neighborhood of 2048 (2K).

As in the previous assignments, we have provided a high resolution timers based on gettimeofday(). This timer is found in Timer.cpp. When run on a multithreaded program, the timer enables you to time individual threads. Recall that the running time of the application is determined by the last core to finish. Consequently, you must provide synchronization to ensure that you collect the correct time which will enable you to compute the parallel running time, Tp.

You'll notice that at the end of the run, you'll get line of output that looks like the following. We will use this terse line of output to test your code. If you are doing several runs, you can quickly search for the appearance of the string "@ " to locate the performance data which includes the input size (-n), the number of threads (-t) running time, and the Gigaflops rate. Following that are flags telling you whether or not you used partial pivoting (P/n), the resulted verified to be correct (Y/N) and which input matrix you used (D/L/H, "D" refers to the default matrix)

@ 2048 8 3.788521 1.511572 N Y D


Turnin

The instructions for turning in this assignment are similar to the previous one, except for some assignment specific details. There is also one modification to the submitted report: it must be in pdf format, and have the following name:

A3-report-GROUP-ID-Student AID's.pdf

GROUP-ID is the name of your group. If your group ID was G-123-456, and your AID numbers are 12345678 and 12345679, then file name of your report should be

A3-report-123-456_12345678-12345679.pdf

If you working alone, then a legal file name would be, for example,

A3-report-123_12345678.pdf

For the remaining instructions on what and how to turnin, visit the following web page

http://www-cse.ucsd.edu/classes/fa13/cse160-a/HW/A3/WhatToTurnin.html

Grading

Your grade will based on results obtained from our testing and examination of your code, as well as documentation.

We will test your code on bang. Be sure to do your final testing on that machine, to avoid any surprises at turnin time. Our test script will build your code against fresh copies of the provided source files mentioned above. Then, it will run your code against various inputs to measure your code's performance, and test for correctness. We will also examine your code manually, since some correctness errors cannot be tested automatically.

Documentation includes the report, comments in the code, and svn log entries. The grade for the code will depend on two factors: performance and correctness. The grade for your documentation will be based on clear explanations, and in the case of the report, insight into the results and presentation of experimental data. Remember: the report is brief so you need not go into long explanations. Focus on what's important. When you commit files with svn, be sure to add brief, but descriptive comments, suitable for understanding your development process. If you aren't sure, be sure to ask us.

To get full credit on the programming portion of the assignment your code must

  1. Be correct.
  2. Run not less than 10% slower than the provided reference, for each combination that we test, of the matrix size -n, the number of threads -t, for each matrix input type -h and -l. Some partial credit will be given for lower performance, but not below about 50% of what the reference can achieve.
  3. Your code must correctly handle the case where the number of threads does not divide the matrix size -n, and should fairly balance the workloads. As described previously, "fairly" means that the number of columns assigned to the threads will differ by at most a single column.

Frequently asked questions

We will add questions to this section as the need arises.


Copyright © 2013 Scott B. Baden Fri   [Wed Oct 30 18:39:11 PDT 2013]