Distributed Array Type Mapping

  

This contains listings of two files. The first, mandel.c, is an OpenMP program that creates a Mandelbrot set. The second, cyclic_array.tcl, implements a type mapping that allows TotalView to display the distributed array created by mandel.c in one Variable Window.

The mandel.c Program

The following program is a simple Mandelbrot program that uses a distributed array. This program is not a good Mandelbrot program as its purpose is to demonstrate a method for handling a cyclically distributed array and how TotalView can reconstruct the global array from the distributed components. Filling the array with the Mandelbrot set simply makes it easier to see that everything is working.

#include <stdio.h>
#include <stdlib.h>
#include <mpi.h>
    
/*
 * This represents a square two dimensional array, where the 
 * first index is cyclically distributed over the processors.
 */
struct cyclic_array 
{
    int * local_elements;
    int local_count;
    int global_count;
    int numprocs;
    int myproc;
};
    
/*
 * Set up the array.
 */
void setup_array (struct cyclic_array *a, int count, int np, int 
me)
{
    a->numprocs = np;
    a->myproc   = me;
    
    a->local_count = count/np + 1;
    a->global_count= count;

    a->local_elements = malloc (sizeof(int) 
                * count * a->local_count);
}
    
/*
 * Compute the owning processor.
 */
int owner_of (struct cyclic_array * a, int ix)
{
    return ix%(a->numprocs);
}
    
/*
 * Compute the address of the local column given the global
 * column's index.
 */
int * local_column (struct cyclic_array * a, int ix)
{
    int owner = owner_of (a, ix);
    
    if (owner != a->myproc)
        return (int *)0;

    return a->local_elements + (ix/a->numprocs) 
                * a->global_count;
}
    
/*
 * Functions for handling complex numbers 
 */
struct complex 
{
    float re;
    float im;
};
    
float cabs_squared (struct complex * z)
{
    return z->re*z->re + z->im*z->im;
}
    
struct complex 
ctimes (struct complex * z0, struct complex * z1)
{
    struct complex res;

    res.re = z0->re * z1->re - z0->im * z1->im;
    res.im = z0->re * z1->im + z0->im * z1->re;
    
    return res;
}
    
struct complex 
cplus (struct complex * z0, struct complex * z1)
{
    struct complex res;
    
    res.re = z0->re + z1->re;
    res.im = z0->im + z1->im;
    
    return res;
}
    
/*
 * Compute the number of iterations needed to determine that 
 * the value is outside the Mandelbrot set.
 */
int mandel (struct complex z0)
{
    struct complex z = z0;
    int i;
    
    for (i=0; i<31; i++)
    {
        if (cabs_squared (&z) >= 4.0)
            break;

         /* The Mandelbrot step, z = z*z + z0 */
        z = ctimes (&z, &z);
        z = cplus (&z, &z0);
    }
    
    return i;
}
    
/* Convert an index into a scaled position in -2..2 */
float position (struct cyclic_array * a, int idx)
{
    idx = idx-(a->global_count/2);
    
    return 4*((0.5 + (float) idx)/a->global_count);
}
    
/*
 * Collect the distributed array and print it.
 * May be rotated...
 */
void print_array (struct cyclic_array * a)
{
    int ix;
    int * buffer = malloc (sizeof (int) * a->global_count);
    
    for (ix = 0; ix<a->global_count; ix++)
    {
        int * col = local_column (a, ix);
        MPI_Status stat;
    
        if (a->myproc == 0)
        { /* First process does the printing */
            int iy;
    
            if (col == 0)
            {
                MPI_Recv (buffer, a->global_count, MPI_INT,
                        MPI_ANY_SOURCE, ix,
                        MPI_COMM_WORLD, &stat);
                col = buffer;
            }
    
            for (iy=0; iy<a->global_count; iy++)
            {
                printf ("%2d ", col[iy]);
            }
            printf ("\n");
        }
        else
        { /* All others send the data */
            if (col != 0)
            {
                MPI_Send (col, a->global_count, MPI_INT, 0,
                        ix, MPI_COMM_WORLD);
            }
        }
    }
    
    MPI_Barrier (MPI_COMM_WORLD);
    free (buffer);
}
    
/* We could read this from the user if we wanted to */
#define GLOBAL_SIZE 40
    
/*
 * Finally we can do the business !
 */
int main (int argc, char ** argv)
{
    int myrank;
    int nprocs;
    struct cyclic_array a;
    int ix;

    MPI_Init(&argc, &argv);
    MPI_Comm_rank(MPI_COMM_WORLD, &myrank);
    MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
    
    setup_array (&a, GLOBAL_SIZE, nprocs, myrank);
    
    /* For each column we own calculate the results */
    for (ix = 0; ix < GLOBAL_SIZE; ix++)
    {
        int * lp = local_column (&a, ix);
    
        if (lp != 0)
        {
            struct complex c;
            int iy;
    
            c.re = position (&a, ix);
    
            for (iy = 0; iy < GLOBAL_SIZE; iy++)
            {
                c.im = - position (&a, iy);        /* - to flip it vertically */
    
                lp [iy] = mandel (c);
            }
        }
    }
    
    print_array (&a);
}