Now the program uses a 2-d torus blocked decomposition.
[charm.git] / examples / charm++ / jacobi2d-iter / jacobi2d.C
1 #include "liveViz.h"
2 #include "jacobi2d.decl.h"
3
4 // See README for documentation
5
6 /*readonly*/ CProxy_Main mainProxy;
7
8 // specify the number of worker chares in each dimension
9 #define num_chare_rows 5
10 #define num_chare_cols 3
11
12 // Each worker chare will process a 4x4 block of elements
13 #define block_width 131
14 #define block_height 61
15
16 // We want to wrap entries around, and because mod operator % sometimes misbehaves on negative values, 
17 // I just wrote these simple wrappers that will make the mod work as expected. -1 maps to the highest value.
18 #define wrap_x(a)  (((a)+num_chare_cols)%num_chare_cols)
19 #define wrap_y(a)  (((a)+num_chare_rows)%num_chare_rows)
20
21 CkArrayID a;
22
23 #define total_iterations 200
24
25 class Main : public CBase_Main
26 {
27 public:
28     int recieve_count;
29     CProxy_Jacobi array;
30     int num_chares;
31     int iterations;
32
33     Main(CkArgMsg* m) {
34         // set iteration counter to zero
35         iterations=0;
36
37         // store the main proxy
38         mainProxy = thisProxy;
39
40         // print info
41         CkPrintf("Running Jacobi on %d processors with (%d,%d) elements\n", CkNumPes(), num_chare_rows, num_chare_cols);
42
43         // Create new array of worker chares
44         array = CProxy_Jacobi::ckNew(num_chare_cols, num_chare_rows);
45
46         // save the total number of worker chares we have in this simulation
47         num_chares = num_chare_rows*num_chare_cols;
48
49         // setup liveviz
50         CkCallback c(CkIndex_Jacobi::requestNextFrame(0),array);
51         liveVizConfig cfg(liveVizConfig::pix_color,true);
52         liveVizInit(cfg,a,c);
53
54         //Start the computation
55         recieve_count = 0;
56         array.begin_iteration();
57     }
58
59     // Each worker reports back to here when it completes an iteration
60     void report(int row, int col) {
61         recieve_count++;
62         if (num_chares == recieve_count) {
63             if (iterations == total_iterations) {
64                 CkPrintf("Completed %d iterations\n", iterations);
65                 CkExit();
66             } else {
67                 CkPrintf("starting new iteration.\n");
68                 recieve_count=0;
69                 iterations++;
70                 // Call begin_iteration on all worker chares in array
71                 array.begin_iteration();
72             }
73         }
74     }
75 };
76
77 class Jacobi: public CBase_Jacobi {
78 public:
79     int messages_due;
80
81     double temperature[block_height][block_width];
82     double left_ghosts[block_height];
83     double right_ghosts[block_height];
84     double top_ghosts[block_width];
85     double bottom_ghosts[block_width];
86
87     // Constructor, initialize values
88     Jacobi() {
89         for(int i=0;i<block_height;++i){
90             for(int j=0;j<block_width;++j){
91                 temperature[i][j] = 0.0;
92             }
93         }
94         BC();
95     }
96
97     // Enforce some boundary conditions
98     void BC(){
99         // Heat left and top edges of each chare's block
100                 for(int i=0;i<block_height;++i)
101             temperature[i][0] = 255.0;
102         for(int j=0;j<block_width;++j)
103             temperature[0][j] = 255.0;
104     }
105
106     // a necessary function which we ignore now
107     // if we were to use load balancing and migration
108     // this function might become useful
109     Jacobi(CkMigrateMessage* m) {}
110
111
112     // Perform one iteration of work
113     // The first step is to send the local state to the neighbors
114     void begin_iteration(void) {
115         messages_due = 4;
116
117         // Copy left column and right column into temporary arrays
118         double left_edge[block_height];
119         double right_edge[block_height];
120                 double top_edge[block_width];
121                 double bottom_edge[block_height];
122
123         for(int i=0;i<block_height;++i){
124             left_edge[i] = temperature[i][0];
125             right_edge[i] = temperature[i][block_width-1];
126         }
127
128         for(int j=0;j<block_width;++j){
129             top_edge[j] = temperature[0][j];
130             bottom_edge[j] = temperature[block_height-1][j];
131         }
132
133         // Send my left edge
134         thisProxy(wrap_x(thisIndex.x-1), thisIndex.y).recieve_neighbor(thisIndex.x, thisIndex.y, block_height, left_edge);
135                 // Send my right edge
136         thisProxy(wrap_x(thisIndex.x+1), thisIndex.y).recieve_neighbor(thisIndex.x, thisIndex.y, block_height, right_edge);
137                 // Send my top edge
138         thisProxy(thisIndex.x, wrap_y(thisIndex.y-1)).recieve_neighbor(thisIndex.x, thisIndex.y, block_width, top_edge);
139                 // Send my bottom edge
140         thisProxy(thisIndex.x, wrap_y(thisIndex.y+1)).recieve_neighbor(thisIndex.x, thisIndex.y, block_width, bottom_edge);
141
142         check_done_iteration();
143     }
144
145     void recieve_neighbor(int neighbor_x, int neighbor_y, int width, double ghost_values[]) {
146         // we have just received a message from worker neighbor_x,neighbor_y with its adjacent
147         // row or column of data values. This worker's index is thisIndex.x, thisIndex.y
148         // We store these in temporary arrays, until all data arrives, then we perform computation
149         // This could be optimized by performing the available computation as soon as the 
150         // required data arrives, but  this example is intentionally simple
151         if(neighbor_x == wrap_x(thisIndex.x-1) && neighbor_y == thisIndex.y){
152             // the ghost data from my LEFT neighbor
153             CkAssert(width == block_height);
154             for(int i=0;i<width;++i){
155                 left_ghosts[i] = ghost_values[i];
156             }
157         } else if(neighbor_x == wrap_x(thisIndex.x+1) && neighbor_y == thisIndex.y){
158             // the ghost data from my RIGHT neighbor
159             CkAssert(width == block_height);
160             for(int i=0;i<width;++i){
161                 right_ghosts[i] = ghost_values[i];
162             }
163         } else  if(neighbor_x == thisIndex.x && neighbor_y == wrap_y(thisIndex.y-1)){
164             // the ghost data from my TOP neighbor
165             CkAssert(width == block_width);
166             for(int i=0;i<width;++i){
167                 top_ghosts[i] = ghost_values[i];
168             }
169         } else if(neighbor_x == thisIndex.x && neighbor_y == wrap_y(thisIndex.y+1)){
170             // the ghost data from my BOTTOM neighbor
171             CkAssert(width == block_width);
172             for(int i=0;i<width;++i){
173                 bottom_ghosts[i] = ghost_values[i];
174             }
175         } else {
176             CkPrintf("Message from non-neighbor chare. I am %d,%d. Message was from %d,%d\n",thisIndex.x,thisIndex.y,neighbor_x,neighbor_y);
177             CkExit();
178         }
179
180         messages_due--;
181         check_done_iteration();
182     }
183
184
185     // Check to see if we have received all neighbor values yet
186     // If all neighbor values have been received, we update our values and proceed
187     void check_done_iteration() {
188         if (messages_due == 0) {
189             // We must create a new array for these values because we don't want to update any of the
190             // the values in temperature[][] array until using them first. Other schemes could be used
191             // to accomplish this same problem. We just put the new values in a temporary array
192             // and write them to temperature[][] after all of the new values are computed.
193             double new_temperature[block_height][block_width];
194     
195             for(int i=0;i<block_height;++i){
196                 for(int j=0;j<block_width;++j){
197                     // first we find the values around the i,j entry in this chare's local block                
198                     double up, down, left, right;
199
200                     if(i==0)
201                         up = top_ghosts[j];
202                     else
203                         up = temperature[i-1][j];
204
205                     if(i==block_height-1)
206                         down = bottom_ghosts[j];
207                     else
208                         down = temperature[i+1][j];
209
210                     if(j==0)
211                         left = left_ghosts[i];
212                     else
213                         left = temperature[i][j-1];
214
215                     if(j==block_width-1)
216                         right = right_ghosts[i];
217                     else
218                         right = temperature[i][j+1];
219
220
221                     // update my value based on the surrounding values
222                                         new_temperature[i][j] = (up+down+left+right+temperature[i][j]) / 5.0;
223
224                 }
225             }
226
227             for(int i=0;i<block_height;++i)
228                 for(int j=0;j<block_width;++j)
229                     temperature[i][j] = new_temperature[i][j];
230
231             // Enforce the boundary conditions again
232             BC();
233
234             mainProxy.report(thisIndex.x, thisIndex.y);
235         }
236     }
237
238
239     // provide my portion of the image to the graphical liveViz client
240     // Currently we just provide some pretty color depending upon the thread id
241     // In a real program we would provide a colored rectangle or pixel that
242     // depends upon the local thread data.
243     void requestNextFrame(liveVizRequestMsg *m){
244                 // These specify the desired total image size requested by the client viewer
245         int wdes = m->req.wid;
246         int hdes = m->req.ht;
247
248         // Deposit a rectangular region to liveViz
249
250         // where to deposit
251         int sx=thisIndex.x*block_width;
252         int sy=thisIndex.y*block_height;
253         int w=block_width,h=block_height; // Size of my rectangular part of the image
254
255         // set the output pixel values for my rectangle
256         // Each component is a char which can have 256 possible values.
257         unsigned char *intensity= new unsigned char[3*w*h];
258         for(int i=0;i<h;++i){
259             for(int j=0;j<w;++j){
260                         intensity[3*(i*w+j)+0] = 255; // RED component
261                         intensity[3*(i*w+j)+1] = 255-temperature[i][j]; // BLUE component
262                         intensity[3*(i*w+j)+2] = 255-temperature[i][j]; // GREEN component
263             }
264         }
265
266         liveVizDeposit(m, sx,sy, w,h, intensity, this);
267         delete[] intensity;
268
269     }
270
271
272
273 };
274
275 #include "jacobi2d.def.h"