Actual source code: comm.c


  2: /***********************************comm.c*************************************

  4: Author: Henry M. Tufo III

  6: e-mail: hmt@cs.brown.edu

  8: snail-mail:
  9: Division of Applied Mathematics
 10: Brown University
 11: Providence, RI 02912

 13: Last Modification:
 14: 11.21.97
 15: ***********************************comm.c*************************************/
 16: #include <../src/ksp/pc/impls/tfs/tfs.h>

 18: /* global program control variables - explicitly exported */
 19: PetscMPIInt PCTFS_my_id            = 0;
 20: PetscMPIInt PCTFS_num_nodes        = 1;
 21: PetscMPIInt PCTFS_floor_num_nodes  = 0;
 22: PetscMPIInt PCTFS_i_log2_num_nodes = 0;

 24: /* global program control variables */
 25: static PetscInt p_init = 0;
 26: static PetscInt modfl_num_nodes;
 27: static PetscInt edge_not_pow_2;

 29: static PetscInt edge_node[sizeof(PetscInt) * 32];

 31: /***********************************comm.c*************************************/
 32: PetscErrorCode PCTFS_comm_init(void)
 33: {
 34:   if (p_init++) return 0;

 36:   MPI_Comm_size(MPI_COMM_WORLD, &PCTFS_num_nodes);
 37:   MPI_Comm_rank(MPI_COMM_WORLD, &PCTFS_my_id);


 41:   PCTFS_ivec_zero((PetscInt *)edge_node, sizeof(PetscInt) * 32);

 43:   PCTFS_floor_num_nodes  = 1;
 44:   PCTFS_i_log2_num_nodes = modfl_num_nodes = 0;
 45:   while (PCTFS_floor_num_nodes <= PCTFS_num_nodes) {
 46:     edge_node[PCTFS_i_log2_num_nodes] = PCTFS_my_id ^ PCTFS_floor_num_nodes;
 47:     PCTFS_floor_num_nodes <<= 1;
 48:     PCTFS_i_log2_num_nodes++;
 49:   }

 51:   PCTFS_i_log2_num_nodes--;
 52:   PCTFS_floor_num_nodes >>= 1;
 53:   modfl_num_nodes = (PCTFS_num_nodes - PCTFS_floor_num_nodes);

 55:   if ((PCTFS_my_id > 0) && (PCTFS_my_id <= modfl_num_nodes)) edge_not_pow_2 = ((PCTFS_my_id | PCTFS_floor_num_nodes) - 1);
 56:   else if (PCTFS_my_id >= PCTFS_floor_num_nodes) edge_not_pow_2 = ((PCTFS_my_id ^ PCTFS_floor_num_nodes) + 1);
 57:   else edge_not_pow_2 = 0;
 58:   return 0;
 59: }

 61: /***********************************comm.c*************************************/
 62: PetscErrorCode PCTFS_giop(PetscInt *vals, PetscInt *work, PetscInt n, PetscInt *oprs)
 63: {
 64:   PetscInt   mask, edge;
 65:   PetscInt   type, dest;
 66:   vfp        fp;
 67:   MPI_Status status;

 69:   /* ok ... should have some data, work, and operator(s) */

 72:   /* non-uniform should have at least two entries */

 75:   /* check to make sure comm package has been initialized */
 76:   if (!p_init) PCTFS_comm_init();

 78:   /* if there's nothing to do return */
 79:   if ((PCTFS_num_nodes < 2) || (!n)) return 0;

 81:   /* a negative number if items to send ==> fatal */

 84:   /* advance to list of n operations for custom */
 85:   if ((type = oprs[0]) == NON_UNIFORM) oprs++;

 87:   /* major league hack */

 90:   /* all msgs will be of the same length */
 91:   /* if not a hypercube must colapse partial dim */
 92:   if (edge_not_pow_2) {
 93:     if (PCTFS_my_id >= PCTFS_floor_num_nodes) {
 94:       MPI_Send(vals, n, MPIU_INT, edge_not_pow_2, MSGTAG0 + PCTFS_my_id, MPI_COMM_WORLD);
 95:     } else {
 96:       MPI_Recv(work, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG0 + edge_not_pow_2, MPI_COMM_WORLD, &status);
 97:       (*fp)(vals, work, n, oprs);
 98:     }
 99:   }

101:   /* implement the mesh fan in/out exchange algorithm */
102:   if (PCTFS_my_id < PCTFS_floor_num_nodes) {
103:     for (mask = 1, edge = 0; edge < PCTFS_i_log2_num_nodes; edge++, mask <<= 1) {
104:       dest = PCTFS_my_id ^ mask;
105:       if (PCTFS_my_id > dest) {
106:         MPI_Send(vals, n, MPIU_INT, dest, MSGTAG2 + PCTFS_my_id, MPI_COMM_WORLD);
107:       } else {
108:         MPI_Recv(work, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG2 + dest, MPI_COMM_WORLD, &status);
109:         (*fp)(vals, work, n, oprs);
110:       }
111:     }

113:     mask = PCTFS_floor_num_nodes >> 1;
114:     for (edge = 0; edge < PCTFS_i_log2_num_nodes; edge++, mask >>= 1) {
115:       if (PCTFS_my_id % mask) continue;

117:       dest = PCTFS_my_id ^ mask;
118:       if (PCTFS_my_id < dest) {
119:         MPI_Send(vals, n, MPIU_INT, dest, MSGTAG4 + PCTFS_my_id, MPI_COMM_WORLD);
120:       } else {
121:         MPI_Recv(vals, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG4 + dest, MPI_COMM_WORLD, &status);
122:       }
123:     }
124:   }

126:   /* if not a hypercube must expand to partial dim */
127:   if (edge_not_pow_2) {
128:     if (PCTFS_my_id >= PCTFS_floor_num_nodes) {
129:       MPI_Recv(vals, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG5 + edge_not_pow_2, MPI_COMM_WORLD, &status);
130:     } else {
131:       MPI_Send(vals, n, MPIU_INT, edge_not_pow_2, MSGTAG5 + PCTFS_my_id, MPI_COMM_WORLD);
132:     }
133:   }
134:   return 0;
135: }

137: /***********************************comm.c*************************************/
138: PetscErrorCode PCTFS_grop(PetscScalar *vals, PetscScalar *work, PetscInt n, PetscInt *oprs)
139: {
140:   PetscInt   mask, edge;
141:   PetscInt   type, dest;
142:   vfp        fp;
143:   MPI_Status status;

145:   /* ok ... should have some data, work, and operator(s) */

148:   /* non-uniform should have at least two entries */

151:   /* check to make sure comm package has been initialized */
152:   if (!p_init) PCTFS_comm_init();

154:   /* if there's nothing to do return */
155:   if ((PCTFS_num_nodes < 2) || (!n)) return 0;

157:   /* a negative number of items to send ==> fatal */

160:   /* advance to list of n operations for custom */
161:   if ((type = oprs[0]) == NON_UNIFORM) oprs++;


165:   /* all msgs will be of the same length */
166:   /* if not a hypercube must colapse partial dim */
167:   if (edge_not_pow_2) {
168:     if (PCTFS_my_id >= PCTFS_floor_num_nodes) {
169:       MPI_Send(vals, n, MPIU_SCALAR, edge_not_pow_2, MSGTAG0 + PCTFS_my_id, MPI_COMM_WORLD);
170:     } else {
171:       MPI_Recv(work, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG0 + edge_not_pow_2, MPI_COMM_WORLD, &status);
172:       (*fp)(vals, work, n, oprs);
173:     }
174:   }

176:   /* implement the mesh fan in/out exchange algorithm */
177:   if (PCTFS_my_id < PCTFS_floor_num_nodes) {
178:     for (mask = 1, edge = 0; edge < PCTFS_i_log2_num_nodes; edge++, mask <<= 1) {
179:       dest = PCTFS_my_id ^ mask;
180:       if (PCTFS_my_id > dest) {
181:         MPI_Send(vals, n, MPIU_SCALAR, dest, MSGTAG2 + PCTFS_my_id, MPI_COMM_WORLD);
182:       } else {
183:         MPI_Recv(work, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG2 + dest, MPI_COMM_WORLD, &status);
184:         (*fp)(vals, work, n, oprs);
185:       }
186:     }

188:     mask = PCTFS_floor_num_nodes >> 1;
189:     for (edge = 0; edge < PCTFS_i_log2_num_nodes; edge++, mask >>= 1) {
190:       if (PCTFS_my_id % mask) continue;

192:       dest = PCTFS_my_id ^ mask;
193:       if (PCTFS_my_id < dest) {
194:         MPI_Send(vals, n, MPIU_SCALAR, dest, MSGTAG4 + PCTFS_my_id, MPI_COMM_WORLD);
195:       } else {
196:         MPI_Recv(vals, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG4 + dest, MPI_COMM_WORLD, &status);
197:       }
198:     }
199:   }

201:   /* if not a hypercube must expand to partial dim */
202:   if (edge_not_pow_2) {
203:     if (PCTFS_my_id >= PCTFS_floor_num_nodes) {
204:       MPI_Recv(vals, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG5 + edge_not_pow_2, MPI_COMM_WORLD, &status);
205:     } else {
206:       MPI_Send(vals, n, MPIU_SCALAR, edge_not_pow_2, MSGTAG5 + PCTFS_my_id, MPI_COMM_WORLD);
207:     }
208:   }
209:   return 0;
210: }

212: /***********************************comm.c*************************************/
213: PetscErrorCode PCTFS_grop_hc(PetscScalar *vals, PetscScalar *work, PetscInt n, PetscInt *oprs, PetscInt dim)
214: {
215:   PetscInt   mask, edge;
216:   PetscInt   type, dest;
217:   vfp        fp;
218:   MPI_Status status;

220:   /* ok ... should have some data, work, and operator(s) */

223:   /* non-uniform should have at least two entries */

226:   /* check to make sure comm package has been initialized */
227:   if (!p_init) PCTFS_comm_init();

229:   /* if there's nothing to do return */
230:   if ((PCTFS_num_nodes < 2) || (!n) || (dim <= 0)) return 0;

232:   /* the error msg says it all!!! */

235:   /* a negative number of items to send ==> fatal */

238:   /* can't do more dimensions then exist */
239:   dim = PetscMin(dim, PCTFS_i_log2_num_nodes);

241:   /* advance to list of n operations for custom */
242:   if ((type = oprs[0]) == NON_UNIFORM) oprs++;


246:   for (mask = 1, edge = 0; edge < dim; edge++, mask <<= 1) {
247:     dest = PCTFS_my_id ^ mask;
248:     if (PCTFS_my_id > dest) {
249:       MPI_Send(vals, n, MPIU_SCALAR, dest, MSGTAG2 + PCTFS_my_id, MPI_COMM_WORLD);
250:     } else {
251:       MPI_Recv(work, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG2 + dest, MPI_COMM_WORLD, &status);
252:       (*fp)(vals, work, n, oprs);
253:     }
254:   }

256:   if (edge == dim) mask >>= 1;
257:   else {
258:     while (++edge < dim) mask <<= 1;
259:   }

261:   for (edge = 0; edge < dim; edge++, mask >>= 1) {
262:     if (PCTFS_my_id % mask) continue;

264:     dest = PCTFS_my_id ^ mask;
265:     if (PCTFS_my_id < dest) {
266:       MPI_Send(vals, n, MPIU_SCALAR, dest, MSGTAG4 + PCTFS_my_id, MPI_COMM_WORLD);
267:     } else {
268:       MPI_Recv(vals, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG4 + dest, MPI_COMM_WORLD, &status);
269:     }
270:   }
271:   return 0;
272: }

274: /******************************************************************************/
275: PetscErrorCode PCTFS_ssgl_radd(PetscScalar *vals, PetscScalar *work, PetscInt level, PetscInt *segs)
276: {
277:   PetscInt     edge, type, dest, mask;
278:   PetscInt     stage_n;
279:   MPI_Status   status;
280:   PetscMPIInt *maxval, flg;

282:   /* check to make sure comm package has been initialized */
283:   if (!p_init) PCTFS_comm_init();

285:   /* all msgs are *NOT* the same length */
286:   /* implement the mesh fan in/out exchange algorithm */
287:   for (mask = 0, edge = 0; edge < level; edge++, mask++) {
288:     stage_n = (segs[level] - segs[edge]);
289:     if (stage_n && !(PCTFS_my_id & mask)) {
290:       dest = edge_node[edge];
291:       type = MSGTAG3 + PCTFS_my_id + (PCTFS_num_nodes * edge);
292:       if (PCTFS_my_id > dest) {
293:         MPI_Send(vals + segs[edge], stage_n, MPIU_SCALAR, dest, type, MPI_COMM_WORLD);
294:       } else {
295:         type = type - PCTFS_my_id + dest;
296:         MPI_Recv(work, stage_n, MPIU_SCALAR, MPI_ANY_SOURCE, type, MPI_COMM_WORLD, &status);
297:         PCTFS_rvec_add(vals + segs[edge], work, stage_n);
298:       }
299:     }
300:     mask <<= 1;
301:   }
302:   mask >>= 1;
303:   for (edge = 0; edge < level; edge++) {
304:     stage_n = (segs[level] - segs[level - 1 - edge]);
305:     if (stage_n && !(PCTFS_my_id & mask)) {
306:       dest = edge_node[level - edge - 1];
307:       type = MSGTAG6 + PCTFS_my_id + (PCTFS_num_nodes * edge);
308:       MPI_Comm_get_attr(MPI_COMM_WORLD, MPI_TAG_UB, &maxval, &flg);
311:       if (PCTFS_my_id < dest) {
312:         MPI_Send(vals + segs[level - 1 - edge], stage_n, MPIU_SCALAR, dest, type, MPI_COMM_WORLD);
313:       } else {
314:         type = type - PCTFS_my_id + dest;
315:         MPI_Recv(vals + segs[level - 1 - edge], stage_n, MPIU_SCALAR, MPI_ANY_SOURCE, type, MPI_COMM_WORLD, &status);
316:       }
317:     }
318:     mask >>= 1;
319:   }
320:   return 0;
321: }

323: /***********************************comm.c*************************************/
324: PetscErrorCode PCTFS_giop_hc(PetscInt *vals, PetscInt *work, PetscInt n, PetscInt *oprs, PetscInt dim)
325: {
326:   PetscInt   mask, edge;
327:   PetscInt   type, dest;
328:   vfp        fp;
329:   MPI_Status status;

331:   /* ok ... should have some data, work, and operator(s) */

334:   /* non-uniform should have at least two entries */

337:   /* check to make sure comm package has been initialized */
338:   if (!p_init) PCTFS_comm_init();

340:   /* if there's nothing to do return */
341:   if ((PCTFS_num_nodes < 2) || (!n) || (dim <= 0)) return 0;

343:   /* the error msg says it all!!! */

346:   /* a negative number of items to send ==> fatal */

349:   /* can't do more dimensions then exist */
350:   dim = PetscMin(dim, PCTFS_i_log2_num_nodes);

352:   /* advance to list of n operations for custom */
353:   if ((type = oprs[0]) == NON_UNIFORM) oprs++;


357:   for (mask = 1, edge = 0; edge < dim; edge++, mask <<= 1) {
358:     dest = PCTFS_my_id ^ mask;
359:     if (PCTFS_my_id > dest) {
360:       MPI_Send(vals, n, MPIU_INT, dest, MSGTAG2 + PCTFS_my_id, MPI_COMM_WORLD);
361:     } else {
362:       MPI_Recv(work, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG2 + dest, MPI_COMM_WORLD, &status);
363:       (*fp)(vals, work, n, oprs);
364:     }
365:   }

367:   if (edge == dim) mask >>= 1;
368:   else {
369:     while (++edge < dim) mask <<= 1;
370:   }

372:   for (edge = 0; edge < dim; edge++, mask >>= 1) {
373:     if (PCTFS_my_id % mask) continue;

375:     dest = PCTFS_my_id ^ mask;
376:     if (PCTFS_my_id < dest) {
377:       MPI_Send(vals, n, MPIU_INT, dest, MSGTAG4 + PCTFS_my_id, MPI_COMM_WORLD);
378:     } else {
379:       MPI_Recv(vals, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG4 + dest, MPI_COMM_WORLD, &status);
380:     }
381:   }
382:   return 0;
383: }