OpenMP

In this assignment you will use OpenMP to solve part of a real problem.

First some practical details:
Modern compilers can read OpenMP-directives and produce threaded code. So if you have access to e.g. the Intel compilers for C or Fortran just compile by icc -openmp my_prog.c or ifort -openmp my_prog.f90 (plus optimization options).

The gnu-type compilers may or may not support OpenMP, it depends on the version. Our gcc does not yet support OpenMP, but we have a newer gcc, gcc4, which does. So gcc4 -fopenmp my_progr.c and gfortran -fopenmp my_prog.f90 should work.

Now to the problem:

The problem is small and simplified part of the tumour therapy computation I mentioned in the introductory lecture. You are given a file containing the coordinates for a , fairly large, number of cells. Each line in the file contains three floating point numbers, the x-, y- and z-coordinate for the nucleus of that particular cell. The task is to write a program that computes the pairwise distances between all the cells. So, let cj = (xj, yj, zj) and ck = (xk, yk, zk) be coordinates for two cells, then your program should compute

d(cj, ck) = sqrt( (xj - xk)2 + (yj - yk)2 + (zj - zk)2 ), for all j and k, j not equal to k

Due to symmetries in the coordinates, many of the distances are the same, so one should compute all distinct distances and the number of times each distance occurs.

A very small example. Suppose that the cells are located in the corners of a cube, having side one. Looking at one corner there are three cells at a distance of one, three cells have the distance sqrt(2) and one cell is at the distance sqrt(3). Due to symmetry (there are eight corners in the cube), the distance one, and sqrt(2), occur 24 times each. sqrt(3) occurs 8 times. In this example, the program should produce a file containing: 

    1.00000000000000   24
    1.41421356237310   24
    1.73205080756888   8

    n_uniq     3
    sum_uniq  56

 
n_uniq is the number of unique distance values, and  sum_uniq the total number of distance values.

The input file may look like:

0 0 0
1 0 0
1 1 0
0 1 0
0 0 1
1 0 1
1 1 1
0 1 1

In the assignment you may not use the fact that there are symmetries, so you should assume that the coordinates have random values. You should use the fact that two distance values are at least 0.01 apart. If they are closer than this, it is caused by roundoff and we regard such distances as being equal. So if one distance, in the cube, is 1.414215 and the other is 1.414213, they are regarded as being equal. sqrt(2) and sqrt(3) are more than 0.01 apart, so they are different distances.

NOTE: in the real application there would be some 700 000 cells. This makes it impossible to allocate an array to store all the distances. Such an array would require some 0.5 × 7000002  = 2.45 × 1011 elements (0.5, since d(cj, ck) = d(ck, cj)) . So, it is not possible to first compute all the distances and then go through the array and count the number of copies. Since this course deals with real programming, and not how to write toy programs, your code is not allowed to allocate space for all (or half of all) the distances either. In this assignment there are only some 20 000 cells, but you should think of them as being 700 000. Note that you are allowed to store the coordinates of the cells.

Hint: you will get fairly short arrays if you use the fact that two (different) distances are at least 0.01 apart, together with a hashing technique.

Write a program in C or Fortran 77 (if you want to use Omni; Fortran 90 if you use ifort) and parallelize it using OpenMP. You may assume  that we will never use more than four processors. It is not hard to get a good speedup in this application. Test your program on the input file, ~thomas/HPC/cells . To debug your program I suggest that you make a shorter file. Try not to waste memory!

To get you started with handling the file, I have provided examples.

Debug HINTS.

Send me the program by e-mail, but don't send me the cells-file.

Here is simple.c from the lecture and here are the equivalent Fortran programs, simple.f and simple.f90.


Back