Wave Equation Implementation (Version 1)
Wave Equation Implementation (Version 1)
Take a look at the code below and note the MPI functions.
Click on the name of each MPI routine to read a
short description of that routine's purpose and syntax.
/*-------------------------------------------------------------------
* This program implements the concurrent wave equation described
* in Chapter 5 of Fox et al., 1988, Solving Problems on Concurrent
* Processors, vol 1.
*
* A vibrating string is decomposed into points. Each process is
* responsible for updating the amplitude of a number of points over
* time.
*
* At each iteration, each process exchanges boundary points with
* nearest neighbors. This version uses low level sends and receives
* to exchange boundary points.
*
* Revised: 11/27/97 by Saleh Elmohamed
*----------------------------------------------------------------- */
#include
#include
#include
#include
#include
#define MASTER 0
#define MAXPOINTS 1000 /* recommended max number of points */
#define MAXSTEPS 10000 /* recommended max number of steps */
#define PI 3.14159265
int E_RtoL = 10;
int E_LtoR = 20;
int E_OUT1 = 30;
int E_OUT2 = 40;
void init_master(void);
void init_workers(void);
void init_line(void);
void update (int left, int right);
void output_master(void);
void output_workers(void);
MPI_Status status;
MPI_Request request;
int procid, /* processor ID */
nproc, /* number of processes */
nsteps, /* number of time steps */
tpoints, /* total points along string */
npoints, /* number of points handled by this processor */
first, /* index of 1st point handled by this processor */
rcode; /* generic return code */
double values[MAXPOINTS+2], /* values at time t */
oldval[MAXPOINTS+2], /* values at time (t-dt) */
newval[MAXPOINTS+2]; /* values at time (t+dt) */
double begin_time, stop_time;
/* time_t t1, t2; */
/*------------------------------------------------------
* Master obtains input values from user
*------------------------------------------------------*/
void init_master(void)
{
int buffer[2];
/*==============================================*/
/* Set number of points, number of iterations */
/* These numbers can be set by the user. It is */
/* recommended that tpoints not to exceed 1000 */
/* and nsteps not to exceed 10000. */
/*==============================================*/
tpoints = 100;
nsteps = 1000;
#if 0
while ((tpoints < nproc) || (tpoints > MAXPOINTS))
{
printf("Enter number of points along vibrating string\n");
scanf("%s", tchar);
tpoints = atoi(tchar);
if ((tpoints < nproc) || (tpoints > MAXPOINTS))
printf("enter value between %d and %d\n", nproc, MAXPOINTS);
}
while ((nsteps < 1) || (nsteps > MAXSTEPS))
{
printf("Enter number of time steps\n");
scanf("%s", tchar);
nsteps = atoi(tchar);
if ((nsteps < 1) || (nsteps > MAXSTEPS))
printf("enter value between 1 and %d\n", MAXSTEPS);
}
#endif
printf("Processor %d: points = %d, steps = %d\n",
procid, tpoints, nsteps);
/* then broadcast time advance parameter, total points, time steps */
buffer[0] = tpoints;
buffer[1] = nsteps;
MPI_Bcast(&buffer, 2, MPI_INT, MASTER, MPI_COMM_WORLD);
}
/*-------------------------------------------------------
* Workers receive input values from master
*------------------------------------------------------*/
void init_workers(void)
{
int buffer[2];
/* then receive total points, time steps */
MPI_Bcast(&buffer, 2, MPI_INT, MASTER, MPI_COMM_WORLD);
tpoints = buffer[0];
nsteps = buffer[1];
}
/*---------------------------------------------------------
* All processes initialize points on line
*------------------------------------------------------- */
void init_line(void)
{
int nmin, nleft, npts, i, j, k;
double x, fac;
/* calculate initial values based on sine curve */
nmin = tpoints/nproc;
nleft = tpoints%nproc;
fac = 2.0 * PI;
for (i = 0, k = 0; i < nproc; i++)
{
npts = (i < nleft) ? nmin + 1 : nmin;
if (procid == i)
{
first = k + 1;
npoints = npts;
printf ("Processor %d: first = %d, npoints = %d\n",
procid, first, npts);
for (j = 1; j <= npts; j++, k++)
{
x = (double)k/(double)(tpoints - 1);
values[j] = sin (fac * x);
}
}
else
k += npts;
}
for (i = 1; i <= npoints; i++)
oldval[i] = values[i];
}
/*------------------------------------------------------------
* Calculate new values using wave equation
*-----------------------------------------------------------*/
void do_math(int i)
{
const double dtime = 0.3;
const double c = 1.0;
const double dx = 1.0;
double tau, sqtau;
tau = (c * dtime / dx);
sqtau = tau * tau;
newval[i] = (2.0 * values[i]) - oldval[i]
+ (sqtau * (values[i-1] - (2.0 * values[i]) + values[i+1]));
}
/*--------------------------------------------------------------------
* All processes update their points a specified number of times
*-------------------------------------------------------------------*/
void update(int left, int right)
{
int i, j, id_rtol, id_ltor;
/* update values for each point along string */
for (i = 1; i <= nsteps; i++)
{
/* exchange data with "left-hand" neighbor */
if (first != 1)
{
MPI_Isend(&values[1], 1,
MPI_DOUBLE, left, E_RtoL, MPI_COMM_WORLD, &request);
MPI_Wait(&request, &status);
MPI_Recv(&values[0], 1,
MPI_DOUBLE, left, E_LtoR, MPI_COMM_WORLD, &status);
}
/* exchange data with "right-hand" neighbor */
if (first + npoints -1 != tpoints)
{
MPI_Isend(&values[npoints], 1,
MPI_DOUBLE, right, E_LtoR, MPI_COMM_WORLD, &request);
MPI_Wait(&request, &status);
MPI_Recv(&values[npoints+1],
1, MPI_DOUBLE, right, E_RtoL, MPI_COMM_WORLD, &status);
}
/* update points along line */
for (j = 1; j <= npoints; j++)
{
/* global endpoints */
if ((first + j - 1 == 1) || (first + j - 1 == tpoints))
newval[j] = 0.0;
else
do_math(j);
}
for (j = 1; j <= npoints; j++)
{
oldval[j] = values[j];
values[j] = newval[j];
}
}
}
/*-----------------------------------------------------------
* Master receives results from workers and prints
*---------------------------------------------------------- */
void output_master(void)
{
int i, start, npts, tpts, buffer[2];
double results[MAXPOINTS];
/* store worker's results in results array */
for (i = 1; i < nproc; i++)
{
/* receive number of points and first point */
MPI_Recv(&buffer, 2, MPI_INT,
MPI_ANY_SOURCE, E_OUT1, MPI_COMM_WORLD, &status);
start = buffer[0];
npts = buffer[1];
/* receive results */
MPI_Recv(&results[start-1],
npts, MPI_DOUBLE, MPI_ANY_SOURCE, E_OUT2, MPI_COMM_WORLD, &status);
}
/* store master's results in results array */
for (i = first; i < first + npoints; i++)
results[i-1] = values[i];
tpts = (tpoints < 10) ? tpoints: 10;
printf("The first %d points (just for validation):\n", tpts);
for (i = 0; i < tpts; i++)
printf("%4.2f ", results[i]);
printf("\n");
}
/*-----------------------------------------------------------
* Workers send the updated values to the master
*----------------------------------------------------------*/
void output_workers(void)
{
int buffer[2];
/* send first point and number of points handled to master */
buffer[0] = first;
buffer[1] = npoints;
MPI_Isend(&buffer, 2, MPI_INT, 0,
E_OUT1, MPI_COMM_WORLD, &request);
MPI_Wait(&request, &status);
/* send results to master */
MPI_Isend(&values[1], npoints,
MPI_DOUBLE, 0, E_OUT2, MPI_COMM_WORLD, &request);
MPI_Wait(&request, &status);
}
/* -----------------------------------------------
* Main program
* --------------------------------------------- */
main(int argc, char **argv)
{
int left, right;
/* learn number of tasks and task ID */
rcode = MPI_Init(&argc, &argv);
begin_time = MPI_Wtime();
MPI_Comm_rank(MPI_COMM_WORLD,
&procid);
MPI_Comm_size(MPI_COMM_WORLD,
&nproc);
if (rcode >= 0)
printf ("Processor %d: Wave Program running ... \n",
procid);
/* determine left and right neighbors */
if (procid == nproc-1)
right = 0;
else
right = procid + 1;
if (procid == 0)
left = nproc - 1;
else
left = procid - 1;
/* get program parameters and initialize wave values */
if (procid == MASTER)
init_master();
else
init_workers();
init_line();
/* t1 = time(NULL); */
/* update values along the line for nstep time steps */
update(left, right);
/* collect results and print */
if (procid == MASTER)
output_master();
else
output_workers();
stop_time = (MPI_Wtime() - begin_time);
printf("... it took %6.3f seconds\n", stop_time);
/* t2 = time(NULL); */
/*
printf("The whole thing took about %4.2f seconds\n",
difftime(t2, t1));
*/
MPI_Finalize();
}
A run on 4 processors will result in the following output:
Processor 1: Wave Program running ...
Processor 0: Wave Program running ...
Processor 3: Wave Program running ...
Processor 2: Wave Program running ...
Processor 0: points = 100, steps = 1000
Processor 0: first = 1, npoints = 25
Processor 1: first = 26, npoints = 25
Processor 2: first = 51, npoints = 25
Processor 3: first = 76, npoints = 25
... it took 2.035 seconds
... it took 2.035 seconds
... it took 2.035 seconds
The first 10 points (just for validation):
0.00 0.06 0.12 0.19 0.25 0.31 0.36 0.42 0.48 0.53
... it took 2.039 seconds