1 2 /***********************************comm.c************************************* 3 4 Author: Henry M. Tufo III 5 6 e-mail: hmt@cs.brown.edu 7 8 snail-mail: 9 Division of Applied Mathematics 10 Brown University 11 Providence, RI 02912 12 13 Last Modification: 14 11.21.97 15 ***********************************comm.c*************************************/ 16 #include <../src/ksp/pc/impls/tfs/tfs.h> 17 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; 23 24 /* global program control variables */ 25 static PetscInt p_init = 0; 26 static PetscInt modfl_num_nodes; 27 static PetscInt edge_not_pow_2; 28 29 static PetscInt edge_node[sizeof(PetscInt) * 32]; 30 31 /***********************************comm.c*************************************/ 32 PetscErrorCode PCTFS_comm_init(void) { 33 PetscFunctionBegin; 34 if (p_init++) PetscFunctionReturn(0); 35 36 MPI_Comm_size(MPI_COMM_WORLD, &PCTFS_num_nodes); 37 MPI_Comm_rank(MPI_COMM_WORLD, &PCTFS_my_id); 38 39 PetscCheck(PCTFS_num_nodes <= (INT_MAX >> 1), PETSC_COMM_SELF, PETSC_ERR_PLIB, "Can't have more then MAX_INT/2 nodes!!!"); 40 41 PCTFS_ivec_zero((PetscInt *)edge_node, sizeof(PetscInt) * 32); 42 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 } 50 51 PCTFS_i_log2_num_nodes--; 52 PCTFS_floor_num_nodes >>= 1; 53 modfl_num_nodes = (PCTFS_num_nodes - PCTFS_floor_num_nodes); 54 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 PetscFunctionReturn(0); 59 } 60 61 /***********************************comm.c*************************************/ 62 PetscErrorCode PCTFS_giop(PetscInt *vals, PetscInt *work, PetscInt n, PetscInt *oprs) { 63 PetscInt mask, edge; 64 PetscInt type, dest; 65 vfp fp; 66 MPI_Status status; 67 68 PetscFunctionBegin; 69 /* ok ... should have some data, work, and operator(s) */ 70 PetscCheck(vals && work && oprs, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop() :: vals=%p, work=%p, oprs=%p", vals, work, oprs); 71 72 /* non-uniform should have at least two entries */ 73 PetscCheck(!(oprs[0] == NON_UNIFORM) || !(n < 2), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop() :: non_uniform and n=0,1?"); 74 75 /* check to make sure comm package has been initialized */ 76 if (!p_init) PCTFS_comm_init(); 77 78 /* if there's nothing to do return */ 79 if ((PCTFS_num_nodes < 2) || (!n)) PetscFunctionReturn(0); 80 81 /* a negative number if items to send ==> fatal */ 82 PetscCheck(n >= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop() :: n=%" PetscInt_FMT "<0?", n); 83 84 /* advance to list of n operations for custom */ 85 if ((type = oprs[0]) == NON_UNIFORM) oprs++; 86 87 /* major league hack */ 88 PetscCheck(fp = (vfp)PCTFS_ivec_fct_addr(type), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop() :: Could not retrieve function pointer!"); 89 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 PetscCallMPI(MPI_Send(vals, n, MPIU_INT, edge_not_pow_2, MSGTAG0 + PCTFS_my_id, MPI_COMM_WORLD)); 95 } else { 96 PetscCallMPI(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 } 100 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 PetscCallMPI(MPI_Send(vals, n, MPIU_INT, dest, MSGTAG2 + PCTFS_my_id, MPI_COMM_WORLD)); 107 } else { 108 PetscCallMPI(MPI_Recv(work, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG2 + dest, MPI_COMM_WORLD, &status)); 109 (*fp)(vals, work, n, oprs); 110 } 111 } 112 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; 116 117 dest = PCTFS_my_id ^ mask; 118 if (PCTFS_my_id < dest) { 119 PetscCallMPI(MPI_Send(vals, n, MPIU_INT, dest, MSGTAG4 + PCTFS_my_id, MPI_COMM_WORLD)); 120 } else { 121 PetscCallMPI(MPI_Recv(vals, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG4 + dest, MPI_COMM_WORLD, &status)); 122 } 123 } 124 } 125 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 PetscCallMPI(MPI_Recv(vals, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG5 + edge_not_pow_2, MPI_COMM_WORLD, &status)); 130 } else { 131 PetscCallMPI(MPI_Send(vals, n, MPIU_INT, edge_not_pow_2, MSGTAG5 + PCTFS_my_id, MPI_COMM_WORLD)); 132 } 133 } 134 PetscFunctionReturn(0); 135 } 136 137 /***********************************comm.c*************************************/ 138 PetscErrorCode PCTFS_grop(PetscScalar *vals, PetscScalar *work, PetscInt n, PetscInt *oprs) { 139 PetscInt mask, edge; 140 PetscInt type, dest; 141 vfp fp; 142 MPI_Status status; 143 144 PetscFunctionBegin; 145 /* ok ... should have some data, work, and operator(s) */ 146 PetscCheck(vals && work && oprs, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop() :: vals=%p, work=%p, oprs=%p", vals, work, oprs); 147 148 /* non-uniform should have at least two entries */ 149 PetscCheck(!(oprs[0] == NON_UNIFORM) || !(n < 2), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop() :: non_uniform and n=0,1?"); 150 151 /* check to make sure comm package has been initialized */ 152 if (!p_init) PCTFS_comm_init(); 153 154 /* if there's nothing to do return */ 155 if ((PCTFS_num_nodes < 2) || (!n)) PetscFunctionReturn(0); 156 157 /* a negative number of items to send ==> fatal */ 158 PetscCheck(n >= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "gdop() :: n=%" PetscInt_FMT "<0?", n); 159 160 /* advance to list of n operations for custom */ 161 if ((type = oprs[0]) == NON_UNIFORM) oprs++; 162 163 PetscCheck(fp = (vfp)PCTFS_rvec_fct_addr(type), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop() :: Could not retrieve function pointer!"); 164 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 PetscCallMPI(MPI_Send(vals, n, MPIU_SCALAR, edge_not_pow_2, MSGTAG0 + PCTFS_my_id, MPI_COMM_WORLD)); 170 } else { 171 PetscCallMPI(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 } 175 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 PetscCallMPI(MPI_Send(vals, n, MPIU_SCALAR, dest, MSGTAG2 + PCTFS_my_id, MPI_COMM_WORLD)); 182 } else { 183 PetscCallMPI(MPI_Recv(work, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG2 + dest, MPI_COMM_WORLD, &status)); 184 (*fp)(vals, work, n, oprs); 185 } 186 } 187 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; 191 192 dest = PCTFS_my_id ^ mask; 193 if (PCTFS_my_id < dest) { 194 PetscCallMPI(MPI_Send(vals, n, MPIU_SCALAR, dest, MSGTAG4 + PCTFS_my_id, MPI_COMM_WORLD)); 195 } else { 196 PetscCallMPI(MPI_Recv(vals, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG4 + dest, MPI_COMM_WORLD, &status)); 197 } 198 } 199 } 200 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 PetscCallMPI(MPI_Recv(vals, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG5 + edge_not_pow_2, MPI_COMM_WORLD, &status)); 205 } else { 206 PetscCallMPI(MPI_Send(vals, n, MPIU_SCALAR, edge_not_pow_2, MSGTAG5 + PCTFS_my_id, MPI_COMM_WORLD)); 207 } 208 } 209 PetscFunctionReturn(0); 210 } 211 212 /***********************************comm.c*************************************/ 213 PetscErrorCode PCTFS_grop_hc(PetscScalar *vals, PetscScalar *work, PetscInt n, PetscInt *oprs, PetscInt dim) { 214 PetscInt mask, edge; 215 PetscInt type, dest; 216 vfp fp; 217 MPI_Status status; 218 219 PetscFunctionBegin; 220 /* ok ... should have some data, work, and operator(s) */ 221 PetscCheck(vals && work && oprs, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop_hc() :: vals=%p, work=%p, oprs=%p", vals, work, oprs); 222 223 /* non-uniform should have at least two entries */ 224 PetscCheck(!(oprs[0] == NON_UNIFORM) || !(n < 2), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop_hc() :: non_uniform and n=0,1?"); 225 226 /* check to make sure comm package has been initialized */ 227 if (!p_init) PCTFS_comm_init(); 228 229 /* if there's nothing to do return */ 230 if ((PCTFS_num_nodes < 2) || (!n) || (dim <= 0)) PetscFunctionReturn(0); 231 232 /* the error msg says it all!!! */ 233 PetscCheck(!modfl_num_nodes, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop_hc() :: PCTFS_num_nodes not a power of 2!?!"); 234 235 /* a negative number of items to send ==> fatal */ 236 PetscCheck(n >= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop_hc() :: n=%" PetscInt_FMT "<0?", n); 237 238 /* can't do more dimensions then exist */ 239 dim = PetscMin(dim, PCTFS_i_log2_num_nodes); 240 241 /* advance to list of n operations for custom */ 242 if ((type = oprs[0]) == NON_UNIFORM) oprs++; 243 244 PetscCheck(fp = (vfp)PCTFS_rvec_fct_addr(type), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop_hc() :: Could not retrieve function pointer!"); 245 246 for (mask = 1, edge = 0; edge < dim; edge++, mask <<= 1) { 247 dest = PCTFS_my_id ^ mask; 248 if (PCTFS_my_id > dest) { 249 PetscCallMPI(MPI_Send(vals, n, MPIU_SCALAR, dest, MSGTAG2 + PCTFS_my_id, MPI_COMM_WORLD)); 250 } else { 251 PetscCallMPI(MPI_Recv(work, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG2 + dest, MPI_COMM_WORLD, &status)); 252 (*fp)(vals, work, n, oprs); 253 } 254 } 255 256 if (edge == dim) mask >>= 1; 257 else { 258 while (++edge < dim) mask <<= 1; 259 } 260 261 for (edge = 0; edge < dim; edge++, mask >>= 1) { 262 if (PCTFS_my_id % mask) continue; 263 264 dest = PCTFS_my_id ^ mask; 265 if (PCTFS_my_id < dest) { 266 PetscCallMPI(MPI_Send(vals, n, MPIU_SCALAR, dest, MSGTAG4 + PCTFS_my_id, MPI_COMM_WORLD)); 267 } else { 268 PetscCallMPI(MPI_Recv(vals, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG4 + dest, MPI_COMM_WORLD, &status)); 269 } 270 } 271 PetscFunctionReturn(0); 272 } 273 274 /******************************************************************************/ 275 PetscErrorCode PCTFS_ssgl_radd(PetscScalar *vals, PetscScalar *work, PetscInt level, PetscInt *segs) { 276 PetscInt edge, type, dest, mask; 277 PetscInt stage_n; 278 MPI_Status status; 279 PetscMPIInt *maxval, flg; 280 281 PetscFunctionBegin; 282 /* check to make sure comm package has been initialized */ 283 if (!p_init) PCTFS_comm_init(); 284 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 PetscCallMPI(MPI_Send(vals + segs[edge], stage_n, MPIU_SCALAR, dest, type, MPI_COMM_WORLD)); 294 } else { 295 type = type - PCTFS_my_id + dest; 296 PetscCallMPI(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 PetscCallMPI(MPI_Comm_get_attr(MPI_COMM_WORLD, MPI_TAG_UB, &maxval, &flg)); 309 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_LIB, "MPI error: MPI_Comm_get_attr() is not returning a MPI_TAG_UB"); 310 PetscCheck(*maxval > type, PETSC_COMM_SELF, PETSC_ERR_PLIB, "MPI_TAG_UB for your current MPI implementation is not large enough to use PCTFS"); 311 if (PCTFS_my_id < dest) { 312 PetscCallMPI(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 PetscCallMPI(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 PetscFunctionReturn(0); 321 } 322 323 /***********************************comm.c*************************************/ 324 PetscErrorCode PCTFS_giop_hc(PetscInt *vals, PetscInt *work, PetscInt n, PetscInt *oprs, PetscInt dim) { 325 PetscInt mask, edge; 326 PetscInt type, dest; 327 vfp fp; 328 MPI_Status status; 329 330 PetscFunctionBegin; 331 /* ok ... should have some data, work, and operator(s) */ 332 PetscCheck(vals && work && oprs, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop_hc() :: vals=%p, work=%p, oprs=%p", vals, work, oprs); 333 334 /* non-uniform should have at least two entries */ 335 PetscCheck(!(oprs[0] == NON_UNIFORM) || !(n < 2), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop_hc() :: non_uniform and n=0,1?"); 336 337 /* check to make sure comm package has been initialized */ 338 if (!p_init) PCTFS_comm_init(); 339 340 /* if there's nothing to do return */ 341 if ((PCTFS_num_nodes < 2) || (!n) || (dim <= 0)) PetscFunctionReturn(0); 342 343 /* the error msg says it all!!! */ 344 PetscCheck(!modfl_num_nodes, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop_hc() :: PCTFS_num_nodes not a power of 2!?!"); 345 346 /* a negative number of items to send ==> fatal */ 347 PetscCheck(n >= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop_hc() :: n=%" PetscInt_FMT "<0?", n); 348 349 /* can't do more dimensions then exist */ 350 dim = PetscMin(dim, PCTFS_i_log2_num_nodes); 351 352 /* advance to list of n operations for custom */ 353 if ((type = oprs[0]) == NON_UNIFORM) oprs++; 354 355 PetscCheck(fp = (vfp)PCTFS_ivec_fct_addr(type), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop_hc() :: Could not retrieve function pointer!"); 356 357 for (mask = 1, edge = 0; edge < dim; edge++, mask <<= 1) { 358 dest = PCTFS_my_id ^ mask; 359 if (PCTFS_my_id > dest) { 360 PetscCallMPI(MPI_Send(vals, n, MPIU_INT, dest, MSGTAG2 + PCTFS_my_id, MPI_COMM_WORLD)); 361 } else { 362 PetscCallMPI(MPI_Recv(work, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG2 + dest, MPI_COMM_WORLD, &status)); 363 (*fp)(vals, work, n, oprs); 364 } 365 } 366 367 if (edge == dim) mask >>= 1; 368 else { 369 while (++edge < dim) mask <<= 1; 370 } 371 372 for (edge = 0; edge < dim; edge++, mask >>= 1) { 373 if (PCTFS_my_id % mask) continue; 374 375 dest = PCTFS_my_id ^ mask; 376 if (PCTFS_my_id < dest) { 377 PetscCallMPI(MPI_Send(vals, n, MPIU_INT, dest, MSGTAG4 + PCTFS_my_id, MPI_COMM_WORLD)); 378 } else { 379 PetscCallMPI(MPI_Recv(vals, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG4 + dest, MPI_COMM_WORLD, &status)); 380 } 381 } 382 PetscFunctionReturn(0); 383 } 384