* Usage: java -Dpj.np=K edu.rit.clu.heat.HotSpotClu imagefile
* H W rl1 cl1 ru1 cu1 temp1 [
* rl2 cl2 ru2 cu2 temp2 . . . ]
*
K = Number of parallel processes
*
imagefile = Output image file name
*
H = Mesh height in pixels (H >= 1)
*
W = Mesh width in pixels (W >= 1)
*
rl1 = First hot spot lower row (1 <= rl1 <=
* W)
*
cl1 = First hot spot lower column (1 <= cl1 <=
* W)
*
ru1 = First hot spot upper row (1 <= ru1 <=
* W)
*
cu1 = First hot spot upper column (1 <= cu1 <=
* W)
*
temp1 = First hot spot temperature (0.0 <= temp1 <=
* 100.0)
*
* The program sets up a mesh of equally-spaced points with H+2 rows and * W+2 columns. The temperature of each boundary point [r,c], * where r = 0, r = H+1, c = 0, or c = * W+1, is fixed at 0.0. The temperatures at certain interior points, * known as "hot spots," are fixed at certain values greater than 0.0. * Specifically, the temperature of each point in the rectangle from * [rl1,cl1] to [ru1,cu1] inclusive is fixed at temp1, the * temperature of each point in the rectangle from [rl2,cl2] to * [ru2,cu2] inclusive is fixed at temp2, and so on. The program * calculates the temperature at every interior point (other than the hot spots) * using successive overrelaxation with Chebyshev acceleration and red-black * updating. The program outputs a PJG color image H+2 pixels high and * W+2 pixels wide. Each pixel's hue depends on the corresponding mesh * point's temperature. A temperature of 0.0 is blue; a temperature of 100.0 is * red; intermediate temperatures are intermediate hues. Each process writes its * own PJG image file with a slice of the image. If imagefile is * specified as "out.pjg", for example, then process 0 writes file * "out_0.pjg", process 1 writes file "out_1.pjg", and so on. *
* The computation is performed in parallel in multiple processors. The program
* measures the computation's running time.
*
* @author Alan Kaminsky
* @version 11-Apr-2008
*/
public class HotSpotClu
{
// Prevent construction.
private HotSpotClu()
{
}
// Hidden constants.
private static final double MIN_TEMP = 0.0;
private static final double MAX_TEMP = 100.0;
private static final double DELTA_TEMP = MAX_TEMP - MIN_TEMP;
private static final double MIN_HUE = 4.0/6.0;
private static final double MAX_HUE = 0.0;
private static final double DELTA_HUE = MAX_HUE - MIN_HUE;
private static final double EPS = 1.0e-3;
private static final int FIRST = 0;
private static final int MIDDLE = 1;
private static final int LAST = 2;
private static final int SINGLE = 3;
// Hidden variables.
// World communicator.
static Comm world;
static int size;
static int rank;
static int position;
static int predRank;
static int succRank;
// Command line arguments.
static File imagefile;
static int H;
static int W;
// Row slice index ranges.
static Range[] slices;
static Range mySlice;
static int myLb;
static int myUb;
static int myLen;
// Temperature grid.
static double[][] h;
// Mesh of hot spot locations.
static boolean[][] hotspot;
// Variables for total absolute residual.
static double EPS_initialTotalAbsXi;
static double totalAbsXi;
// Other variables used in the successive overrrelaxation algorithm.
static int MAXITER;
static double rho_s_sqr;
static double omega_over_4;
static int iterations;
// Communication buffers.
static DoubleItemBuf xibuf;
static DoubleBuf hbuf_pred_red;
static DoubleBuf hbuf_pred_black;
static DoubleBuf hbuf_top_red;
static DoubleBuf hbuf_top_black;
static DoubleBuf hbuf_bottom_red;
static DoubleBuf hbuf_bottom_black;
static DoubleBuf hbuf_succ_red;
static DoubleBuf hbuf_succ_black;
// Main program.
/**
* Main program.
*/
public static void main
(String[] args)
throws Exception
{
// Start timing.
long t1 = System.currentTimeMillis();
// Initialize world communicator.
Comm.init (args);
world = Comm.world();
size = world.size();
rank = world.rank();
if (size == 1) position = SINGLE;
else if (rank == 0) position = FIRST;
else if (rank < size-1) position = MIDDLE;
else position = LAST;
predRank = rank - 1;
succRank = rank + 1;
// Parse command line arguments.
if (args.length < 8 || (args.length % 5) != 3) usage();
imagefile = new File (args[0]);
H = Integer.parseInt (args[1]);
W = Integer.parseInt (args[2]);
if (H < 1) usage();
if (W < 1) usage();
// Determine row slice index ranges.
slices = new Range (1, H) .subranges (size);
mySlice = slices[rank];
myLb = mySlice.lb();
myUb = mySlice.ub();
myLen = mySlice.length();
// Initialize temperature and hot spot meshes.
h = new double [H+2] [];
Arrays.allocate (h, new Range (myLb-1, myUb+1), W+2);
hotspot = new boolean [H+2] [];
Arrays.allocate (hotspot, new Range (myLb-1, myUb+1), W+2);
// Record hot spot coordinates and temperatures.
int n = (args.length - 3) / 5;
for (int i = 0; i < n; ++ i)
{
int rl = Integer.parseInt (args[3+5*i]);
int cl = Integer.parseInt (args[4+5*i]);
int ru = Integer.parseInt (args[5+5*i]);
int cu = Integer.parseInt (args[6+5*i]);
double temp = Double.parseDouble (args[7+5*i]);
if (1 > rl || rl > W) usage();
if (1 > cl || cl > H) usage();
if (1 > ru || ru > W) usage();
if (1 > cu || cu > H) usage();
if (MIN_TEMP > temp || temp > MAX_TEMP) usage();
for (int r = rl; r <= ru; ++ r)
{
double[] h_r = h[r];
boolean[] hotspot_r = hotspot[r];
if (h_r != null)
{
for (int c = cl; c <= cu; ++ c)
{
h_r[c] = temp;
hotspot_r[c] = true;
}
}
}
}
// Initialize communication buffers.
xibuf = DoubleBuf.buffer();
hbuf_pred_red = redBuffer (myLb-1);
hbuf_pred_black = blackBuffer (myLb-1);
hbuf_top_red = redBuffer (myLb);
hbuf_top_black = blackBuffer (myLb);
hbuf_bottom_red = redBuffer (myUb);
hbuf_bottom_black = blackBuffer (myUb);
hbuf_succ_red = redBuffer (myUb+1);
hbuf_succ_black = blackBuffer (myUb+1);
// Compute initial total absolute residual, then multiply by EPS.
totalAbsXi = 0.0;
double xi;
for (int r = myLb; r <= myUb; ++ r)
{
double[] h_rm1 = h[r-1];
double[] h_r = h[r];
double[] h_rp1 = h[r+1];
boolean[] hotspot_r = hotspot[r];
for (int c = 1; c <= W; ++ c)
{
xi =
hotspot_r[c] ?
0.0 :
h_rm1[c]+h_rp1[c]+h_r[c-1]+h_r[c+1]-4.0*h_r[c];
totalAbsXi += Math.abs (xi);
}
}
xibuf.item = totalAbsXi;
world.allReduce (xibuf, DoubleOp.SUM);
totalAbsXi = xibuf.item;
EPS_initialTotalAbsXi = EPS * totalAbsXi;
// Initialize other variables.
MAXITER = 2 * (W + H);
rho_s_sqr = 0.5 * (Math.cos (Math.PI/W) + Math.cos (Math.PI/H));
rho_s_sqr = rho_s_sqr * rho_s_sqr;
omega_over_4 = 0.25;
iterations = 0;
long t2 = System.currentTimeMillis();
// Perform successive overrelaxation.
do
{
totalAbsXi = 0.0;
// Red half-sweep.
for (int r = myLb; r <= myUb; ++ r)
{
double[] h_rm1 = h[r-1];
double[] h_r = h[r];
double[] h_rp1 = h[r+1];
boolean[] hotspot_r = hotspot[r];
for (int c = 1 + (r&1); c <= W; c += 2)
{
xi =
hotspot_r[c] ?
0.0 :
h_rm1[c]+h_rp1[c]+h_r[c-1]+h_r[c+1]-4.0*h_r[c];
totalAbsXi += Math.abs (xi);
h_r[c] += omega_over_4 * xi;
}
}
omega_over_4 = 0.25 /
(1.0 - rho_s_sqr * (iterations == 0 ? 0.5 : omega_over_4));
// Exchange boundary row red cells with neighboring processes.
switch (position)
{
case FIRST:
world.send (succRank, hbuf_bottom_red);
world.receive (succRank, hbuf_succ_red);
break;
case MIDDLE:
world.sendReceive
(succRank, hbuf_bottom_red,
predRank, hbuf_pred_red);
world.sendReceive
(predRank, hbuf_top_red,
succRank, hbuf_succ_red);
break;
case LAST:
world.receive (predRank, hbuf_pred_red);
world.send (predRank, hbuf_top_red);
break;
}
// Black half-sweep.
for (int r = myLb; r <= myUb; ++ r)
{
double[] h_rm1 = h[r-1];
double[] h_r = h[r];
double[] h_rp1 = h[r+1];
boolean[] hotspot_r = hotspot[r];
for (int c = 2 - (r&1); c <= W; c += 2)
{
xi =
hotspot_r[c] ?
0.0 :
h_rm1[c]+h_rp1[c]+h_r[c-1]+h_r[c+1]-4.0*h_r[c];
totalAbsXi += Math.abs (xi);
h_r[c] += omega_over_4 * xi;
}
}
omega_over_4 = 0.25 / (1.0 - rho_s_sqr * omega_over_4);
// Exchange boundary row black cells with neighboring processes.
switch (position)
{
case FIRST:
world.send (succRank, hbuf_bottom_black);
world.receive (succRank, hbuf_succ_black);
break;
case MIDDLE:
world.sendReceive
(succRank, hbuf_bottom_black,
predRank, hbuf_pred_black);
world.sendReceive
(predRank, hbuf_top_black,
succRank, hbuf_succ_black);
break;
case LAST:
world.receive (predRank, hbuf_pred_black);
world.send (predRank, hbuf_top_black);
break;
}
// Determine total absolute residual from all processes.
xibuf.item = totalAbsXi;
world.allReduce (xibuf, DoubleOp.SUM);
totalAbsXi = xibuf.item;
++ iterations;
}
while (iterations < MAXITER && totalAbsXi >= EPS_initialTotalAbsXi);
// Check for convergence.
if (iterations == MAXITER)
{
System.err.println ("HotSpotClu: Did not converge");
System.exit (1);
}
long t3 = System.currentTimeMillis();
// Generate image.
int[][] matrix = new int [H+2] [];
int rlb = rank == 0 ? myLb-1 : myLb;
int rub = rank == size-1 ? myUb+1 : myUb;
Arrays.allocate (matrix, new Range (rlb, rub), W+2);
ColorImageRow matrix_r = new ColorImageRow (matrix[rlb]);
for (int r = rlb; r <= rub; ++ r)
{
double[] h_r = h[r];
matrix_r.setArray (matrix[r]);
for (int c = 0; c <= W+1; ++ c)
{
matrix_r.setPixelHSB
(/*c */ c,
/*hue*/ (float)
((h_r[c]-MIN_TEMP)/DELTA_TEMP*DELTA_HUE+MIN_HUE),
/*sat*/ 1.0f,
/*bri*/ 1.0f);
}
}
PJGHueImage image = new PJGHueImage (H+2, W+2, matrix);
PJGImage.Writer writer =
image.prepareToWrite
(new BufferedOutputStream
(new FileOutputStream
(Files.fileForRank (imagefile, rank))));
writer.writeRowSlice (new Range (rlb, rub));
writer.close();
// Stop timing.
long t4 = System.currentTimeMillis();
System.out.println (iterations + " iterations " + rank);
System.out.println ((t2-t1) + " msec pre " + rank);
System.out.println ((t3-t2) + " msec calc " + rank);
System.out.println ((t4-t3) + " msec post " + rank);
System.out.println ((t4-t1) + " msec total " + rank);
}
// Hidden operations.
/**
* Returns a communication buffer for the red columns of the given row of
* the h matrix.
*
* @param r Row index.
*
* @return Communication buffer.
*/
private static DoubleBuf redBuffer
(int r)
{
return DoubleBuf.sliceBuffer (h[r], new Range (1 + (r&1), W, 2));
}
/**
* Returns a communication buffer for the black columns of the given row of
* the h matrix.
*
* @param r Row index.
*
* @return Communication buffer.
*/
private static DoubleBuf blackBuffer
(int r)
{
return DoubleBuf.sliceBuffer (h[r], new Range (2 - (r&1), W, 2));
}
/**
* Print a usage message and exit.
*/
private static void usage()
{
System.err.println ("Usage: java -Dpj.np=