Computer Science 400 Parallel Processing Siena College Fall 2008

Your task this week is to write a C or C++ program using OpenMP that solves Laplace's equation on a two-dimensional, uniform, square grid, using Jacobi iteration. Don't worry if none of those terms make any sense - this document tells you what little you need to know about the math and physics.

Laplace's equation is an elliptic partial differential equation that governs physical phenomena such as heat. In two dimensions, it can be written

Phi_xx + Phi_yy = 0.

Given a spatial region and values for points on the boundaries of the region, the goal is to approximate the steady-state solution for points in the interior. We do this by covering the region with an evenly-spaced grid of points. A grid of 8 ×8 would look like this:

* * * * * * * * * *
* . . . . . . . . *
* . . . . . . . . *
* . . . . . . . . *
* . . . . . . . . *
* . . . . . . . . *
* . . . . . . . . *
* . . . . . . . . *
* . . . . . . . . *
* * * * * * * * * *

There are several iterative methods for solving Laplace's equation. Your program is to use Jacobi iteration, which is the simplest and easily parallelizable, though certainly not the most efficient in terms of convergence rate.

In Jacobi iteration, the new value for each grid point in the interior is set to the average of the old values of the four points left, right, above, and below it. This process is repeated until the program terminates. Note that some of the values used for the average will be boundary points.

• To avoid special cases at the boundaries, you should allocate your grid for an n ×n simulation to be (n+2) ×(n+2), so you can use row and column 0 and row and column n+1 to store the boundary values.
• You will need two copies of the grid, one to store the "current" solution values and one to store the "next" solution values. Do not copy the "next" solution to the "current" solution at each step. Instead, always do two iterations. The first uses one grid as "current" and the other as "next," then the second swaps their roles. This is an example of a simple loop unrolling.
• Each grid cell computation looks something like this:

    nextgrid[i][j] = (grid[i-1][j] + grid[i+1][j] +
                      grid[i][j-1] + grid[i][j+1]) * 0.25;
  

  Doing * 0.25 is usually faster than / 4 since multiplication is an easier operation for a computer than division.
• After each pair of iterations, the corresponding values in each cell of the two grids are compared, and the maximum such value is computed. If this value is less than eps, the computation terminates. Otherwise, it continues.
• A maximum number of iterations is also specified, after which the computation stops even if the error tolerance has not been reached.
• Your program should take two command-line arguments: the maximum number of iterations and the error tolerance eps.
• Initialization of the boundary and interior will determine the exact problem being solved. I suggest initializing your interior to all 0's and the boundary to have two sides set to 1, two sides set to 0, as follows:

  1 1 1 1 1 1 1 1 1 1
  1 . . . . . . . . 0
  1 . . . . . . . . 0
  1 . . . . . . . . 0
  1 . . . . . . . . 0
  1 . . . . . . . . 0
  1 . . . . . . . . 0
  1 . . . . . . . . 0
  1 . . . . . . . . 0
  1 0 0 0 0 0 0 0 0 0
  

  This keeps the initialization simple, but still allows for some work during the solution phase.
• At the end of the computation, print out the total number of iterations needed, and the time taken to achieve the solution.
• Get the serial version working first. Choose values for the size of the grid, maximum number of iterations, and eps so that your program will run for a few minutes.
• Parallelize your program. Your parallel program should create and destroy the OpenMP threads only once.
• You may break down the computation any way you'd like. Justify your choice in the README file that you include in your submission.
• You will need a number of OpenMP directives. Explain each directive (what you needed to do and why you chose that particular directive) in detailed comments in your source code.
• Run your program using 1, 2, and 4 threads on the four-processor Solaris nodes of the cluster (wetteland or rivera). Include the timings and your analysis of those timings in the README file that you include with your submission.
• For a small amount of extra credit, you can extend your program to include extra features, such as support for more interesting initial conditions and boundary conditions. Describe these features in your README file.

When you are finished, submit using the turnin utility as lab5. Your submission should include an appropriate Makefile, your C or C++ source code (including the timer code from class, if you choose to use it), a brief README file expaining how to run your program and with your timings and analysis. Please do not include object files or your executable.

determined by correctness, design, documentation, and style, as well as the presentation of your timing results.