charm examples: clean up and add features to jacobi3d (from pmodels book)
authorJonathan Lifflander <avunda@166.1.168.192.in-addr.arpa>
Thu, 13 Sep 2012 04:10:52 +0000 (23:10 -0500)
committerJonathan Lifflander <avunda@166.1.168.192.in-addr.arpa>
Thu, 13 Sep 2012 04:10:52 +0000 (23:10 -0500)
examples/charm++/jacobi3d-sdag/jacobi3d.C
examples/charm++/jacobi3d-sdag/jacobi3d.ci

index d9263ac79e3d24117b357d07d843981f6c69c81d..a4c4b63344f6e048a61d615787c8302518007a0f 100644 (file)
@@ -1,28 +1,16 @@
-/** \file jacobi3d.C
- *  Author: Abhinav S Bhatele
- *  Date Created: June 01st, 2009
- *
- *     
- *           *****************
- *        *               *  *
- *   ^ *****************     *
- *   | *               *     *
- *   | *               *     *
- *   | *               *     *
- *   Y *               *     *
- *   | *               *     *
- *   | *               *     *
- *   | *               *  * 
- *   ~ *****************    Z
- *     <------ X ------> 
- *
- *   X: left, right --> wrap_x
- *   Y: top, bottom --> wrap_y
- *   Z: front, back --> wrap_z
- */
+#define MAX_ITER               26
+#define LEFT                   1
+#define RIGHT                  2
+#define TOP                    3
+#define BOTTOM                 4
+#define FRONT                  5
+#define BACK                   6
+#define DIVIDEBY7              0.14285714285714285714
+#define DELTA                  0.01
+#define LBPERIOD                       50
+#define CHECKPOINTPERIOD        500
 
 #include "jacobi3d.decl.h"
 
 #include "jacobi3d.decl.h"
-#include "TopoManager.h"
 
 /*readonly*/ CProxy_Main mainProxy;
 /*readonly*/ int arrayDimX;
 
 /*readonly*/ CProxy_Main mainProxy;
 /*readonly*/ int arrayDimX;
 /*readonly*/ int num_chare_y;
 /*readonly*/ int num_chare_z;
 
 /*readonly*/ int num_chare_y;
 /*readonly*/ int num_chare_z;
 
-static unsigned long next = 1;
-
-int myrand(int numpes) {
-  next = next * 1103515245 + 12345;
-  return((unsigned)(next/65536) % numpes);
-}
-
-// We want to wrap entries around, and because mod operator % 
-// sometimes misbehaves on negative values. -1 maps to the highest value.
-#define wrap_x(a)      (((a)+num_chare_x)%num_chare_x)
-#define wrap_y(a)      (((a)+num_chare_y)%num_chare_y)
-#define wrap_z(a)      (((a)+num_chare_z)%num_chare_z)
+#define wrapX(a)       (((a)+num_chare_x)%num_chare_x)
+#define wrapY(a)       (((a)+num_chare_y)%num_chare_y)
+#define wrapZ(a)       (((a)+num_chare_z)%num_chare_z)
 
 #define index(a,b,c)   ((a)+(b)*(blockDimX+2)+(c)*(blockDimX+2)*(blockDimY+2))
 
 
 #define index(a,b,c)   ((a)+(b)*(blockDimX+2)+(c)*(blockDimX+2)*(blockDimY+2))
 
-#define MAX_ITER               26
-#define WARM_ITER              5
-#define LEFT                   1
-#define RIGHT                  2
-#define TOP                    3
-#define BOTTOM                 4
-#define FRONT                  5
-#define BACK                   6
-#define DIVIDEBY7              0.14285714285714285714
-
 double startTime;
 double endTime;
 
 double startTime;
 double endTime;
 
@@ -69,75 +38,67 @@ double endTime;
  *
  */
 class Main : public CBase_Main {
  *
  */
 class Main : public CBase_Main {
-  public:
-    CProxy_Jacobi array;
-    int iterations;
-
-    Main(CkArgMsg* m) {
-      if ( (m->argc != 3) && (m->argc != 7) ) {
-        CkPrintf("%s [array_size] [block_size]\n", m->argv[0]);
-        CkPrintf("OR %s [array_size_X] [array_size_Y] [array_size_Z] [block_size_X] [block_size_Y] [block_size_Z]\n", m->argv[0]);
-        CkAbort("Abort");
-      }
+public:
+  CProxy_Jacobi array;
+  int iterations;
+
+  Main(CkArgMsg* m) {
+    if ( (m->argc != 3) && (m->argc != 7) ) {
+      CkPrintf("%s [array_size] [block_size]\n", m->argv[0]);
+      CkPrintf("OR %s [array_size_X] [array_size_Y] [array_size_Z] [block_size_X] [block_size_Y] [block_size_Z]\n", m->argv[0]);
+      CkAbort("Abort");
+    }
 
 
-      // set iteration counter to zero
-      iterations = 0;
+    // set iteration counter to zero
+    iterations = 0;
 
 
-      // store the main proxy
-      mainProxy = thisProxy;
+    // store the main proxy
+    mainProxy = thisProxy;
        
        
-      if(m->argc == 3) {
-       arrayDimX = arrayDimY = arrayDimZ = atoi(m->argv[1]);
-        blockDimX = blockDimY = blockDimZ = atoi(m->argv[2]); 
-      }
-      else if (m->argc == 7) {
-        arrayDimX = atoi(m->argv[1]);
-       arrayDimY = atoi(m->argv[2]);
-       arrayDimZ = atoi(m->argv[3]);
-        blockDimX = atoi(m->argv[4]); 
-       blockDimY = atoi(m->argv[5]); 
-       blockDimZ = atoi(m->argv[6]);
-      }
-
-      if (arrayDimX < blockDimX || arrayDimX % blockDimX != 0)
-        CkAbort("array_size_X % block_size_X != 0!");
-      if (arrayDimY < blockDimY || arrayDimY % blockDimY != 0)
-        CkAbort("array_size_Y % block_size_Y != 0!");
-      if (arrayDimZ < blockDimZ || arrayDimZ % blockDimZ != 0)
-        CkAbort("array_size_Z % block_size_Z != 0!");
-
-      num_chare_x = arrayDimX / blockDimX;
-      num_chare_y = arrayDimY / blockDimY;
-      num_chare_z = arrayDimZ / blockDimZ;
-
-      // print info
-      CkPrintf("\nSTENCIL COMPUTATION WITH NO BARRIERS\n");
-      CkPrintf("Running Jacobi on %d processors with (%d, %d, %d) chares\n", CkNumPes(), num_chare_x, num_chare_y, num_chare_z);
-      CkPrintf("Array Dimensions: %d %d %d\n", arrayDimX, arrayDimY, arrayDimZ);
-      CkPrintf("Block Dimensions: %d %d %d\n", blockDimX, blockDimY, blockDimZ);
-
-      // Create new array of worker chares
-      array = CProxy_Jacobi::ckNew(num_chare_x, num_chare_y, num_chare_z);
-
-      //Start the computation
-      array.doStep();
+    if(m->argc == 3) {
+      arrayDimX = arrayDimY = arrayDimZ = atoi(m->argv[1]);
+      blockDimX = blockDimY = blockDimZ = atoi(m->argv[2]); 
     }
     }
-
-    // Each worker reports back to here when it completes an iteration
-    void report() {
-      iterations++;
-      if (iterations <= WARM_ITER) {
-       if (iterations == WARM_ITER)
-         startTime = CkWallTimer();
-       array.doStep();
-      }
-      else {
-       CkPrintf("Completed %d iterations\n", MAX_ITER-1);
-       endTime = CkWallTimer();
-       CkPrintf("Time elapsed per iteration: %f\n", (endTime - startTime)/(MAX_ITER-1-WARM_ITER));
-        CkExit();
-      }
+    else if (m->argc == 7) {
+      arrayDimX = atoi(m->argv[1]);
+      arrayDimY = atoi(m->argv[2]);
+      arrayDimZ = atoi(m->argv[3]);
+      blockDimX = atoi(m->argv[4]); 
+      blockDimY = atoi(m->argv[5]); 
+      blockDimZ = atoi(m->argv[6]);
     }
     }
+
+    if (arrayDimX < blockDimX || arrayDimX % blockDimX != 0)
+      CkAbort("array_size_X % block_size_X != 0!");
+    if (arrayDimY < blockDimY || arrayDimY % blockDimY != 0)
+      CkAbort("array_size_Y % block_size_Y != 0!");
+    if (arrayDimZ < blockDimZ || arrayDimZ % blockDimZ != 0)
+      CkAbort("array_size_Z % block_size_Z != 0!");
+
+    num_chare_x = arrayDimX / blockDimX;
+    num_chare_y = arrayDimY / blockDimY;
+    num_chare_z = arrayDimZ / blockDimZ;
+
+    // print info
+    CkPrintf("\nSTENCIL COMPUTATION WITH NO BARRIERS\n");
+    CkPrintf("Running Jacobi on %d processors with (%d, %d, %d) chares\n", CkNumPes(), num_chare_x, num_chare_y, num_chare_z);
+    CkPrintf("Array Dimensions: %d %d %d\n", arrayDimX, arrayDimY, arrayDimZ);
+    CkPrintf("Block Dimensions: %d %d %d\n", blockDimX, blockDimY, blockDimZ);
+
+    // Create new array of worker chares
+    array = CProxy_Jacobi::ckNew(num_chare_x, num_chare_y, num_chare_z);
+
+    //Start the computation
+    array.run();
+    startTime = CkWallTimer();
+  }
+
+  void done(int iterations) {
+    CkPrintf("Completed %d iterations\n", iterations);
+    endTime = CkWallTimer();
+    CkPrintf("Time elapsed per iteration: %f\n", (endTime - startTime) / iterations);
+    CkExit();
+  }
 };
 
 /** \class Jacobi
 };
 
 /** \class Jacobi
@@ -147,210 +108,195 @@ class Main : public CBase_Main {
 class Jacobi: public CBase_Jacobi {
   Jacobi_SDAG_CODE
 
 class Jacobi: public CBase_Jacobi {
   Jacobi_SDAG_CODE
 
-  public:
-    int iterations;
-    int imsg;
+public:
+  int iterations;
+  int remoteCount;
 
 
-    double *temperature;
-    double *new_temperature;
+  double *temperature;
+  double *new_temperature;
+  bool converged;
 
 
-    // Constructor, initialize values
-    Jacobi() {
-      __sdag_init();
+  // Constructor, initialize values
+  Jacobi() {
+    __sdag_init();
 
 
-      int i, j, k;
-      // allocate a three dimensional array
-      temperature = new double[(blockDimX+2) * (blockDimY+2) * (blockDimZ+2)];
-      new_temperature = new double[(blockDimX+2) * (blockDimY+2) * (blockDimZ+2)];
+    usesAtSync = CmiTrue;
+    converged = false;
 
 
-      for(k=0; k<blockDimZ+2; ++k)
-       for(j=0; j<blockDimY+2; ++j)
-         for(i=0; i<blockDimX+2; ++i)
-           temperature[index(i, j, k)] = 0.0;
+    // allocate a three dimensional array
+    temperature = new double[(blockDimX+2) * (blockDimY+2) * (blockDimZ+2)];
+    new_temperature = new double[(blockDimX+2) * (blockDimY+2) * (blockDimZ+2)];
 
 
-      iterations = 0;
-      imsg = 0;
-      constrainBC();
-    }
+    for(int k=0; k<blockDimZ+2; ++k)
+      for(int j=0; j<blockDimY+2; ++j)
+        for(int i=0; i<blockDimX+2; ++i)
+          temperature[index(i, j, k)] = 0.0;
+
+    iterations = 0;
+    constrainBC();
+  }
 
   void pup(PUP::er &p)
   {
     CBase_Jacobi::pup(p);
     __sdag_pup(p);
     p|iterations;
 
   void pup(PUP::er &p)
   {
     CBase_Jacobi::pup(p);
     __sdag_pup(p);
     p|iterations;
-    p|imsg;
+    p|remoteCount;
 
     size_t size = (blockDimX+2) * (blockDimY+2) * (blockDimZ+2);
     if (p.isUnpacking()) {
 
     size_t size = (blockDimX+2) * (blockDimY+2) * (blockDimZ+2);
     if (p.isUnpacking()) {
-       temperature = new double[size];
-       new_temperature = new double[size];
-      }
+      temperature = new double[size];
+      new_temperature = new double[size];
+    }
     p(temperature, size);
     p(new_temperature, size);
   }
 
   Jacobi(CkMigrateMessage* m) {__sdag_init();}
 
     p(temperature, size);
     p(new_temperature, size);
   }
 
   Jacobi(CkMigrateMessage* m) {__sdag_init();}
 
-    ~Jacobi() { 
-      delete [] temperature; 
-      delete [] new_temperature; 
-    }
+  ~Jacobi() { 
+    delete [] temperature; 
+    delete [] new_temperature; 
+  }
 
 
-    // Send ghost faces to the six neighbors
-    void begin_iteration(void) {
-      if (thisIndex.x == 0 && thisIndex.y == 0 && thisIndex.z == 0) {
-          CkPrintf("Start of iteration %d\n", iterations);
-          //BgPrintf("BgPrint> Start of iteration at %f\n");
+  // Send ghost faces to the six neighbors
+  void begin_iteration(void) {
+    // Copy different faces into messages
+    double *leftGhost =  new double[blockDimY*blockDimZ];
+    double *rightGhost =  new double[blockDimY*blockDimZ];
+    double *topGhost =  new double[blockDimX*blockDimZ];
+    double *bottomGhost =  new double[blockDimX*blockDimZ];
+    double *frontGhost =  new double[blockDimX*blockDimY];
+    double *backGhost =  new double[blockDimX*blockDimY];
+
+    for(int k=0; k<blockDimZ; ++k)
+      for(int j=0; j<blockDimY; ++j) {
+        leftGhost[k*blockDimY+j] = temperature[index(1, j+1, k+1)];
+        rightGhost[k*blockDimY+j] = temperature[index(blockDimX, j+1, k+1)];
       }
       }
-      iterations++;
-
-      // Copy different faces into messages
-      double *leftGhost =  new double[blockDimY*blockDimZ];
-      double *rightGhost =  new double[blockDimY*blockDimZ];
-      double *topGhost =  new double[blockDimX*blockDimZ];
-      double *bottomGhost =  new double[blockDimX*blockDimZ];
-      double *frontGhost =  new double[blockDimX*blockDimY];
-      double *backGhost =  new double[blockDimX*blockDimY];
-
-      for(int k=0; k<blockDimZ; ++k)
-       for(int j=0; j<blockDimY; ++j) {
-         leftGhost[k*blockDimY+j] = temperature[index(1, j+1, k+1)];
-         rightGhost[k*blockDimY+j] = temperature[index(blockDimX, j+1, k+1)];
-       }
-
-      for(int k=0; k<blockDimZ; ++k)
-       for(int i=0; i<blockDimX; ++i) {
-         topGhost[k*blockDimX+i] = temperature[index(i+1, 1, k+1)];
-         bottomGhost[k*blockDimX+i] = temperature[index(i+1, blockDimY, k+1)];
-       }
-
-      for(int j=0; j<blockDimY; ++j)
-       for(int i=0; i<blockDimX; ++i) {
-         frontGhost[j*blockDimX+i] = temperature[index(i+1, j+1, 1)];
-         backGhost[j*blockDimX+i] = temperature[index(i+1, j+1, blockDimZ)];
-       }
-
-      // Send my left face
-      thisProxy(wrap_x(thisIndex.x-1), thisIndex.y, thisIndex.z)
-         .receiveGhosts(iterations, RIGHT, blockDimY, blockDimZ, leftGhost);
-      // Send my right face
-      thisProxy(wrap_x(thisIndex.x+1), thisIndex.y, thisIndex.z)
-         .receiveGhosts(iterations, LEFT, blockDimY, blockDimZ, rightGhost);
-      // Send my bottom face
-      thisProxy(thisIndex.x, wrap_y(thisIndex.y-1), thisIndex.z)
-         .receiveGhosts(iterations, TOP, blockDimX, blockDimZ, bottomGhost);
-      // Send my top face
-      thisProxy(thisIndex.x, wrap_y(thisIndex.y+1), thisIndex.z)
-         .receiveGhosts(iterations, BOTTOM, blockDimX, blockDimZ, topGhost);
-      // Send my front face
-      thisProxy(thisIndex.x, thisIndex.y, wrap_z(thisIndex.z-1))
-         .receiveGhosts(iterations, BACK, blockDimX, blockDimY, frontGhost);
-      // Send my back face
-      thisProxy(thisIndex.x, thisIndex.y, wrap_z(thisIndex.z+1))
-         .receiveGhosts(iterations, FRONT, blockDimX, blockDimY, backGhost);
-
-      delete[] leftGhost;
-      delete[] rightGhost;
-      delete[] bottomGhost;
-      delete[] topGhost;
-      delete[] frontGhost;
-      delete[] backGhost;
-    }
 
 
-    void processGhosts(int dir, int height, int width, double gh[]) {
-      switch(dir) {
-       case LEFT:
-         for(int k=0; k<width; ++k)
-           for(int j=0; j<height; ++j) {
-             temperature[index(0, j+1, k+1)] = gh[k*height+j];
-           }
-         break;
-       case RIGHT:
-         for(int k=0; k<width; ++k)
-           for(int j=0; j<height; ++j) {
-             temperature[index(blockDimX+1, j+1, k+1)] = gh[k*height+j];
-           }
-         break;
-       case BOTTOM:
-         for(int k=0; k<width; ++k)
-           for(int i=0; i<height; ++i) {
-             temperature[index(i+1, 0, k+1)] = gh[k*height+i];
-           }
-         break;
-       case TOP:
-         for(int k=0; k<width; ++k)
-           for(int i=0; i<height; ++i) {
-             temperature[index(i+1, blockDimY+1, k+1)] = gh[k*height+i];
-           }
-         break;
-       case FRONT:
-         for(int j=0; j<width; ++j)
-           for(int i=0; i<height; ++i) {
-             temperature[index(i+1, j+1, 0)] = gh[j*height+i];
-           }
-         break;
-       case BACK:
-         for(int j=0; j<width; ++j)
-           for(int i=0; i<height; ++i) {
-             temperature[index(i+1, j+1, blockDimZ+1)] = gh[j*height+i];
-           }
-         break;
-       default:
-          CkAbort("ERROR\n");
+    for(int k=0; k<blockDimZ; ++k)
+      for(int i=0; i<blockDimX; ++i) {
+        topGhost[k*blockDimX+i] = temperature[index(i+1, 1, k+1)];
+        bottomGhost[k*blockDimX+i] = temperature[index(i+1, blockDimY, k+1)];
       }
       }
-    }
 
 
+    for(int j=0; j<blockDimY; ++j)
+      for(int i=0; i<blockDimX; ++i) {
+        frontGhost[j*blockDimX+i] = temperature[index(i+1, j+1, 1)];
+        backGhost[j*blockDimX+i] = temperature[index(i+1, j+1, blockDimZ)];
+      }
 
 
-    void check_and_compute() {
-      compute_kernel();
+    int x = thisIndex.x, y = thisIndex.y, z = thisIndex.z;
+    thisProxy(wrapX(x-1),y,z).updateGhosts(iterations, RIGHT,  blockDimY, blockDimZ, rightGhost);
+    thisProxy(wrapX(x+1),y,z).updateGhosts(iterations, LEFT,   blockDimY, blockDimZ, leftGhost);
+    thisProxy(x,wrapY(y-1),z).updateGhosts(iterations, TOP,    blockDimX, blockDimZ, topGhost);
+    thisProxy(x,wrapY(y+1),z).updateGhosts(iterations, BOTTOM, blockDimX, blockDimZ, bottomGhost);
+    thisProxy(x,y,wrapZ(z-1)).updateGhosts(iterations, BACK,   blockDimX, blockDimY, backGhost);
+    thisProxy(x,y,wrapZ(z+1)).updateGhosts(iterations, FRONT,  blockDimX, blockDimY, frontGhost);
+
+    delete [] leftGhost;
+    delete [] rightGhost;
+    delete [] bottomGhost;
+    delete [] topGhost;
+    delete [] frontGhost;
+    delete [] backGhost;
+  }
 
 
-      // calculate error
-      // not being done right now since we are doing a fixed no. of iterations
+  void updateBoundary(int dir, int height, int width, double* gh) {
+    switch(dir) {
+    case LEFT:
+      for(int k=0; k<width; ++k)
+        for(int j=0; j<height; ++j) {
+          temperature[index(0, j+1, k+1)] = gh[k*height+j];
+        }
+      break;
+    case RIGHT:
+      for(int k=0; k<width; ++k)
+        for(int j=0; j<height; ++j) {
+          temperature[index(blockDimX+1, j+1, k+1)] = gh[k*height+j];
+        }
+      break;
+    case BOTTOM:
+      for(int k=0; k<width; ++k)
+        for(int i=0; i<height; ++i) {
+          temperature[index(i+1, 0, k+1)] = gh[k*height+i];
+        }
+      break;
+    case TOP:
+      for(int k=0; k<width; ++k)
+        for(int i=0; i<height; ++i) {
+          temperature[index(i+1, blockDimY+1, k+1)] = gh[k*height+i];
+        }
+      break;
+    case FRONT:
+      for(int j=0; j<width; ++j)
+        for(int i=0; i<height; ++i) {
+          temperature[index(i+1, j+1, 0)] = gh[j*height+i];
+        }
+      break;
+    case BACK:
+      for(int j=0; j<width; ++j)
+        for(int i=0; i<height; ++i) {
+          temperature[index(i+1, j+1, blockDimZ+1)] = gh[j*height+i];
+        }
+      break;
+    default:
+      CkAbort("ERROR\n");
+    }
+  }
 
 
-      double *tmp;
-      tmp = temperature;
-      temperature = new_temperature;
-      new_temperature = tmp;
+  // Check to see if we have received all neighbor values yet
+  // If all neighbor values have been received, we update our values and proceed
+  double computeKernel() {
+#pragma unroll    
+    for(int k=1; k<blockDimZ+1; ++k)
+      for(int j=1; j<blockDimY+1; ++j)
+        for(int i=1; i<blockDimX+1; ++i) {
+          // update my value based on the surrounding values
+          new_temperature[index(i, j, k)] = (temperature[index(i-1, j, k)] 
+                                             +  temperature[index(i+1, j, k)]
+                                             +  temperature[index(i, j-1, k)]
+                                             +  temperature[index(i, j+1, k)]
+                                             +  temperature[index(i, j, k-1)]
+                                             +  temperature[index(i, j, k+1)]
+                                             +  temperature[index(i, j, k)] ) * DIVIDEBY7;
+        } // end for
+    
+    double error = 0.0, max_error = 0.0;
+
+    for(int k=1; k<blockDimZ+1; ++k)
+      for(int j=1; j<blockDimY+1; ++j)
+        for(int i=1; i<blockDimX+1; ++i) {
+          error = fabs(new_temperature[index(i,j,k)] - temperature[index(i,j,k)]);
+          if (error > max_error) {
+            max_error = error;
+          }
+        }
 
 
-      constrainBC();
+    double *tmp;
+    tmp = temperature;
+    temperature = new_temperature;
+    new_temperature = tmp;
 
 
-      if (iterations <= WARM_ITER || iterations >= MAX_ITER)
-       contribute(0, 0, CkReduction::concat, CkCallback(CkIndex_Main::report(), mainProxy));
-      else
-       doStep();
-    }
+    constrainBC();
 
 
-    // Check to see if we have received all neighbor values yet
-    // If all neighbor values have been received, we update our values and proceed
-    void compute_kernel() {
-#pragma unroll    
-      for(int k=1; k<blockDimZ+1; ++k)
-       for(int j=1; j<blockDimY+1; ++j)
-         for(int i=1; i<blockDimX+1; ++i) {
-           // update my value based on the surrounding values
-           new_temperature[index(i, j, k)] = (temperature[index(i-1, j, k)] 
-                                           +  temperature[index(i+1, j, k)]
-                                           +  temperature[index(i, j-1, k)]
-                                           +  temperature[index(i, j+1, k)]
-                                           +  temperature[index(i, j, k-1)]
-                                           +  temperature[index(i, j, k+1)]
-                                           +  temperature[index(i, j, k)] ) * DIVIDEBY7;
-         } // end for
-    }
+    return max_error;
+  }
 
 
-    // Enforce some boundary conditions
-    void constrainBC() {
-      // Heat left, top and front faces of each chare's block
-      for(int k=1; k<blockDimZ+1; ++k)
-       for(int i=1; i<blockDimX+1; ++i)
-         temperature[index(i, 1, k)] = 255.0;
-      for(int k=1; k<blockDimZ+1; ++k)
-       for(int j=1; j<blockDimY+1; ++j)
-         temperature[index(1, j, k)] = 255.0;
+  // Enforce some boundary conditions
+  void constrainBC() {
+    // Heat left, top and front faces of each chare's block
+    for(int k=1; k<blockDimZ+1; ++k)
+      for(int i=1; i<blockDimX+1; ++i)
+        temperature[index(i, 1, k)] = 255.0;
+    for(int k=1; k<blockDimZ+1; ++k)
       for(int j=1; j<blockDimY+1; ++j)
       for(int j=1; j<blockDimY+1; ++j)
-       for(int i=1; i<blockDimX+1; ++i)
-         temperature[index(i, j, 1)] = 255.0;
-    }
+        temperature[index(1, j, k)] = 255.0;
+    for(int j=1; j<blockDimY+1; ++j)
+      for(int i=1; i<blockDimX+1; ++i)
+        temperature[index(i, j, 1)] = 255.0;
+  }
 
 };
 
 
 };
 
index e744b1e544fdca0286a3c6d9782f5df704219004..10ee4d3820678dfa8a091934021d1e3f10e3851c 100644 (file)
@@ -1,5 +1,4 @@
 mainmodule jacobi3d {
 mainmodule jacobi3d {
-
   readonly CProxy_Main mainProxy;
   readonly int arrayDimX;
   readonly int arrayDimY;
   readonly CProxy_Main mainProxy;
   readonly int arrayDimX;
   readonly int arrayDimY;
@@ -7,39 +6,47 @@ mainmodule jacobi3d {
   readonly int blockDimX;
   readonly int blockDimY;
   readonly int blockDimZ;
   readonly int blockDimX;
   readonly int blockDimY;
   readonly int blockDimZ;
-
   readonly int num_chare_x;
   readonly int num_chare_y;
   readonly int num_chare_z;
 
   mainchare Main {
     entry Main(CkArgMsg *m);
   readonly int num_chare_x;
   readonly int num_chare_y;
   readonly int num_chare_z;
 
   mainchare Main {
     entry Main(CkArgMsg *m);
-    entry void report();
+    entry void done(int iterations);
   };
 
   array [3D] Jacobi {
     entry Jacobi(void);
   };
 
   array [3D] Jacobi {
     entry Jacobi(void);
-    entry void begin_iteration(void);
-    entry void receiveGhosts(int iter, int dir, int height, int width,
-                             double ghosts[height*width]);
-
-    entry void doStep() {
-      atomic "begin_iteration" {
-       begin_iteration();
-      }
-      for(imsg = 0; imsg < 6; imsg++) {
-       // "iterations" keeps track of messages across steps
-       when
-          receiveGhosts[iterations] (int iter, int dir, int height, 
-                                     int width, double ghosts[height*width])
-         atomic "process ghosts" {
-            processGhosts(dir, height, width, ghosts);
+    entry void updateGhosts(int ref, int dir, int w, int h, double gh[w*h]);
+    entry [reductiontarget] void checkConverged(bool result);
+    entry void checkpointFinished();
+    entry void run() {
+      while (!converged) {
+        atomic { begin_iteration(); }
+        for (remoteCount = 0; remoteCount < 6; remoteCount++) {
+          when updateGhosts[iterations](int ref, int dir, int w, int h, double buf[w*h]) atomic {
+            updateBoundary(dir, w, h, buf);
           }
           }
-      }
-      atomic "doWork" {
-       check_and_compute();
+        }
+        atomic {
+          double error = computeKernel();
+          int conv = error < DELTA;
+          if (iterations % 5 == 1) {
+            contribute(sizeof(int), &conv, CkReduction::logical_and, CkCallback(CkReductionTarget(Jacobi, checkConverged), thisProxy));
+          }
+        }
+        if (++iterations % 5 == 0)
+          when checkConverged(bool result)
+            if (result) atomic { mainProxy.done(iterations); converged = true; }
+        if (iterations % LBPERIOD == 0) atomic { AtSync(); }
+        if (iterations % CHECKPOINTPERIOD == 0) {
+          atomic {
+            CkCallback cb(CkIndex_Jacobi::checkpointFinished(), thisProxy);
+            CkStartMemCheckpoint(cb);
+          }
+          when checkpointFinished() { }
+        }
       }
     };
   };
       }
     };
   };
-
 };
 };