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 { 34 PetscFunctionBegin; 35 if (p_init++) PetscFunctionReturn(PETSC_SUCCESS); 36 37 PetscCallMPI(MPI_Comm_size(MPI_COMM_WORLD, &PCTFS_num_nodes)); 38 PetscCallMPI(MPI_Comm_rank(MPI_COMM_WORLD, &PCTFS_my_id)); 39 40 PetscCheck(PCTFS_num_nodes <= (INT_MAX >> 1), PETSC_COMM_SELF, PETSC_ERR_PLIB, "Can't have more then MAX_INT/2 nodes!!!"); 41 42 PetscCall(PCTFS_ivec_zero((PetscInt *)edge_node, sizeof(PetscInt) * 32)); 43 44 PCTFS_floor_num_nodes = 1; 45 PCTFS_i_log2_num_nodes = modfl_num_nodes = 0; 46 while (PCTFS_floor_num_nodes <= PCTFS_num_nodes) { 47 edge_node[PCTFS_i_log2_num_nodes] = PCTFS_my_id ^ PCTFS_floor_num_nodes; 48 PCTFS_floor_num_nodes <<= 1; 49 PCTFS_i_log2_num_nodes++; 50 } 51 52 PCTFS_i_log2_num_nodes--; 53 PCTFS_floor_num_nodes >>= 1; 54 modfl_num_nodes = (PCTFS_num_nodes - PCTFS_floor_num_nodes); 55 56 if ((PCTFS_my_id > 0) && (PCTFS_my_id <= modfl_num_nodes)) edge_not_pow_2 = ((PCTFS_my_id | PCTFS_floor_num_nodes) - 1); 57 else if (PCTFS_my_id >= PCTFS_floor_num_nodes) edge_not_pow_2 = ((PCTFS_my_id ^ PCTFS_floor_num_nodes) + 1); 58 else edge_not_pow_2 = 0; 59 PetscFunctionReturn(PETSC_SUCCESS); 60 } 61 62 /***********************************comm.c*************************************/ 63 PetscErrorCode PCTFS_giop(PetscInt *vals, PetscInt *work, PetscInt n, PetscInt *oprs) 64 { 65 PetscInt mask, edge; 66 PetscInt type, dest; 67 vfp fp; 68 MPI_Status status; 69 70 PetscFunctionBegin; 71 /* ok ... should have some data, work, and operator(s) */ 72 PetscCheck(vals && work && oprs, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop() :: vals=%p, work=%p, oprs=%p", (void *)vals, (void *)work, (void *)oprs); 73 74 /* non-uniform should have at least two entries */ 75 PetscCheck(!(oprs[0] == NON_UNIFORM) || !(n < 2), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop() :: non_uniform and n=0,1?"); 76 77 /* check to make sure comm package has been initialized */ 78 if (!p_init) PetscCall(PCTFS_comm_init()); 79 80 /* if there's nothing to do return */ 81 if ((PCTFS_num_nodes < 2) || (!n)) PetscFunctionReturn(PETSC_SUCCESS); 82 83 /* a negative number if items to send ==> fatal */ 84 PetscCheck(n >= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop() :: n=%" PetscInt_FMT "<0?", n); 85 86 /* advance to list of n operations for custom */ 87 if ((type = oprs[0]) == NON_UNIFORM) oprs++; 88 89 /* major league hack */ 90 PetscCheck(fp = (vfp)PCTFS_ivec_fct_addr(type), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop() :: Could not retrieve function pointer!"); 91 92 /* all msgs will be of the same length */ 93 /* if not a hypercube must colapse partial dim */ 94 if (edge_not_pow_2) { 95 if (PCTFS_my_id >= PCTFS_floor_num_nodes) { 96 PetscCallMPI(MPI_Send(vals, n, MPIU_INT, edge_not_pow_2, MSGTAG0 + PCTFS_my_id, MPI_COMM_WORLD)); 97 } else { 98 PetscCallMPI(MPI_Recv(work, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG0 + edge_not_pow_2, MPI_COMM_WORLD, &status)); 99 PetscCall((*fp)(vals, work, n, oprs)); 100 } 101 } 102 103 /* implement the mesh fan in/out exchange algorithm */ 104 if (PCTFS_my_id < PCTFS_floor_num_nodes) { 105 for (mask = 1, edge = 0; edge < PCTFS_i_log2_num_nodes; edge++, mask <<= 1) { 106 dest = PCTFS_my_id ^ mask; 107 if (PCTFS_my_id > dest) { 108 PetscCallMPI(MPI_Send(vals, n, MPIU_INT, dest, MSGTAG2 + PCTFS_my_id, MPI_COMM_WORLD)); 109 } else { 110 PetscCallMPI(MPI_Recv(work, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG2 + dest, MPI_COMM_WORLD, &status)); 111 PetscCall((*fp)(vals, work, n, oprs)); 112 } 113 } 114 115 mask = PCTFS_floor_num_nodes >> 1; 116 for (edge = 0; edge < PCTFS_i_log2_num_nodes; edge++, mask >>= 1) { 117 if (PCTFS_my_id % mask) continue; 118 119 dest = PCTFS_my_id ^ mask; 120 if (PCTFS_my_id < dest) { 121 PetscCallMPI(MPI_Send(vals, n, MPIU_INT, dest, MSGTAG4 + PCTFS_my_id, MPI_COMM_WORLD)); 122 } else { 123 PetscCallMPI(MPI_Recv(vals, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG4 + dest, MPI_COMM_WORLD, &status)); 124 } 125 } 126 } 127 128 /* if not a hypercube must expand to partial dim */ 129 if (edge_not_pow_2) { 130 if (PCTFS_my_id >= PCTFS_floor_num_nodes) { 131 PetscCallMPI(MPI_Recv(vals, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG5 + edge_not_pow_2, MPI_COMM_WORLD, &status)); 132 } else { 133 PetscCallMPI(MPI_Send(vals, n, MPIU_INT, edge_not_pow_2, MSGTAG5 + PCTFS_my_id, MPI_COMM_WORLD)); 134 } 135 } 136 PetscFunctionReturn(PETSC_SUCCESS); 137 } 138 139 /***********************************comm.c*************************************/ 140 PetscErrorCode PCTFS_grop(PetscScalar *vals, PetscScalar *work, PetscInt n, PetscInt *oprs) 141 { 142 PetscInt mask, edge; 143 PetscInt type, dest; 144 vfp fp; 145 MPI_Status status; 146 147 PetscFunctionBegin; 148 /* ok ... should have some data, work, and operator(s) */ 149 PetscCheck(vals && work && oprs, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop() :: vals=%p, work=%p, oprs=%p", (void *)vals, (void *)work, (void *)oprs); 150 151 /* non-uniform should have at least two entries */ 152 PetscCheck(!(oprs[0] == NON_UNIFORM) || !(n < 2), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop() :: non_uniform and n=0,1?"); 153 154 /* check to make sure comm package has been initialized */ 155 if (!p_init) PetscCall(PCTFS_comm_init()); 156 157 /* if there's nothing to do return */ 158 if ((PCTFS_num_nodes < 2) || (!n)) PetscFunctionReturn(PETSC_SUCCESS); 159 160 /* a negative number of items to send ==> fatal */ 161 PetscCheck(n >= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "gdop() :: n=%" PetscInt_FMT "<0?", n); 162 163 /* advance to list of n operations for custom */ 164 if ((type = oprs[0]) == NON_UNIFORM) oprs++; 165 166 PetscCheck((fp = (vfp)PCTFS_rvec_fct_addr(type)), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop() :: Could not retrieve function pointer!"); 167 168 /* all msgs will be of the same length */ 169 /* if not a hypercube must colapse partial dim */ 170 if (edge_not_pow_2) { 171 if (PCTFS_my_id >= PCTFS_floor_num_nodes) { 172 PetscCallMPI(MPI_Send(vals, n, MPIU_SCALAR, edge_not_pow_2, MSGTAG0 + PCTFS_my_id, MPI_COMM_WORLD)); 173 } else { 174 PetscCallMPI(MPI_Recv(work, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG0 + edge_not_pow_2, MPI_COMM_WORLD, &status)); 175 PetscCall((*fp)(vals, work, n, oprs)); 176 } 177 } 178 179 /* implement the mesh fan in/out exchange algorithm */ 180 if (PCTFS_my_id < PCTFS_floor_num_nodes) { 181 for (mask = 1, edge = 0; edge < PCTFS_i_log2_num_nodes; edge++, mask <<= 1) { 182 dest = PCTFS_my_id ^ mask; 183 if (PCTFS_my_id > dest) { 184 PetscCallMPI(MPI_Send(vals, n, MPIU_SCALAR, dest, MSGTAG2 + PCTFS_my_id, MPI_COMM_WORLD)); 185 } else { 186 PetscCallMPI(MPI_Recv(work, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG2 + dest, MPI_COMM_WORLD, &status)); 187 PetscCall((*fp)(vals, work, n, oprs)); 188 } 189 } 190 191 mask = PCTFS_floor_num_nodes >> 1; 192 for (edge = 0; edge < PCTFS_i_log2_num_nodes; edge++, mask >>= 1) { 193 if (PCTFS_my_id % mask) continue; 194 195 dest = PCTFS_my_id ^ mask; 196 if (PCTFS_my_id < dest) { 197 PetscCallMPI(MPI_Send(vals, n, MPIU_SCALAR, dest, MSGTAG4 + PCTFS_my_id, MPI_COMM_WORLD)); 198 } else { 199 PetscCallMPI(MPI_Recv(vals, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG4 + dest, MPI_COMM_WORLD, &status)); 200 } 201 } 202 } 203 204 /* if not a hypercube must expand to partial dim */ 205 if (edge_not_pow_2) { 206 if (PCTFS_my_id >= PCTFS_floor_num_nodes) { 207 PetscCallMPI(MPI_Recv(vals, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG5 + edge_not_pow_2, MPI_COMM_WORLD, &status)); 208 } else { 209 PetscCallMPI(MPI_Send(vals, n, MPIU_SCALAR, edge_not_pow_2, MSGTAG5 + PCTFS_my_id, MPI_COMM_WORLD)); 210 } 211 } 212 PetscFunctionReturn(PETSC_SUCCESS); 213 } 214 215 /***********************************comm.c*************************************/ 216 PetscErrorCode PCTFS_grop_hc(PetscScalar *vals, PetscScalar *work, PetscInt n, PetscInt *oprs, PetscInt dim) 217 { 218 PetscInt mask, edge; 219 PetscInt type, dest; 220 vfp fp; 221 MPI_Status status; 222 223 PetscFunctionBegin; 224 /* ok ... should have some data, work, and operator(s) */ 225 PetscCheck(vals && work && oprs, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop_hc() :: vals=%p, work=%p, oprs=%p", (void *)vals, (void *)work, (void *)oprs); 226 227 /* non-uniform should have at least two entries */ 228 PetscCheck(!(oprs[0] == NON_UNIFORM) || !(n < 2), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop_hc() :: non_uniform and n=0,1?"); 229 230 /* check to make sure comm package has been initialized */ 231 if (!p_init) PetscCall(PCTFS_comm_init()); 232 233 /* if there's nothing to do return */ 234 if ((PCTFS_num_nodes < 2) || (!n) || (dim <= 0)) PetscFunctionReturn(PETSC_SUCCESS); 235 236 /* the error msg says it all!!! */ 237 PetscCheck(!modfl_num_nodes, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop_hc() :: PCTFS_num_nodes not a power of 2!?!"); 238 239 /* a negative number of items to send ==> fatal */ 240 PetscCheck(n >= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop_hc() :: n=%" PetscInt_FMT "<0?", n); 241 242 /* can't do more dimensions then exist */ 243 dim = PetscMin(dim, PCTFS_i_log2_num_nodes); 244 245 /* advance to list of n operations for custom */ 246 if ((type = oprs[0]) == NON_UNIFORM) oprs++; 247 248 PetscCheck(fp = (vfp)PCTFS_rvec_fct_addr(type), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_grop_hc() :: Could not retrieve function pointer!"); 249 250 for (mask = 1, edge = 0; edge < dim; edge++, mask <<= 1) { 251 dest = PCTFS_my_id ^ mask; 252 if (PCTFS_my_id > dest) { 253 PetscCallMPI(MPI_Send(vals, n, MPIU_SCALAR, dest, MSGTAG2 + PCTFS_my_id, MPI_COMM_WORLD)); 254 } else { 255 PetscCallMPI(MPI_Recv(work, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG2 + dest, MPI_COMM_WORLD, &status)); 256 PetscCall((*fp)(vals, work, n, oprs)); 257 } 258 } 259 260 if (edge == dim) mask >>= 1; 261 else { 262 while (++edge < dim) mask <<= 1; 263 } 264 265 for (edge = 0; edge < dim; edge++, mask >>= 1) { 266 if (PCTFS_my_id % mask) continue; 267 268 dest = PCTFS_my_id ^ mask; 269 if (PCTFS_my_id < dest) { 270 PetscCallMPI(MPI_Send(vals, n, MPIU_SCALAR, dest, MSGTAG4 + PCTFS_my_id, MPI_COMM_WORLD)); 271 } else { 272 PetscCallMPI(MPI_Recv(vals, n, MPIU_SCALAR, MPI_ANY_SOURCE, MSGTAG4 + dest, MPI_COMM_WORLD, &status)); 273 } 274 } 275 PetscFunctionReturn(PETSC_SUCCESS); 276 } 277 278 /******************************************************************************/ 279 PetscErrorCode PCTFS_ssgl_radd(PetscScalar *vals, PetscScalar *work, PetscInt level, PetscInt *segs) 280 { 281 PetscInt edge, type, dest, mask; 282 PetscInt stage_n; 283 MPI_Status status; 284 PetscMPIInt *maxval, flg; 285 286 PetscFunctionBegin; 287 /* check to make sure comm package has been initialized */ 288 if (!p_init) PetscCall(PCTFS_comm_init()); 289 290 /* all msgs are *NOT* the same length */ 291 /* implement the mesh fan in/out exchange algorithm */ 292 for (mask = 0, edge = 0; edge < level; edge++, mask++) { 293 stage_n = (segs[level] - segs[edge]); 294 if (stage_n && !(PCTFS_my_id & mask)) { 295 dest = edge_node[edge]; 296 type = MSGTAG3 + PCTFS_my_id + (PCTFS_num_nodes * edge); 297 if (PCTFS_my_id > dest) { 298 PetscCallMPI(MPI_Send(vals + segs[edge], stage_n, MPIU_SCALAR, dest, type, MPI_COMM_WORLD)); 299 } else { 300 type = type - PCTFS_my_id + dest; 301 PetscCallMPI(MPI_Recv(work, stage_n, MPIU_SCALAR, MPI_ANY_SOURCE, type, MPI_COMM_WORLD, &status)); 302 PetscCall(PCTFS_rvec_add(vals + segs[edge], work, stage_n)); 303 } 304 } 305 mask <<= 1; 306 } 307 mask >>= 1; 308 for (edge = 0; edge < level; edge++) { 309 stage_n = (segs[level] - segs[level - 1 - edge]); 310 if (stage_n && !(PCTFS_my_id & mask)) { 311 dest = edge_node[level - edge - 1]; 312 type = MSGTAG6 + PCTFS_my_id + (PCTFS_num_nodes * edge); 313 PetscCallMPI(MPI_Comm_get_attr(MPI_COMM_WORLD, MPI_TAG_UB, &maxval, &flg)); 314 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_LIB, "MPI error: MPI_Comm_get_attr() is not returning a MPI_TAG_UB"); 315 PetscCheck(*maxval > type, PETSC_COMM_SELF, PETSC_ERR_PLIB, "MPI_TAG_UB for your current MPI implementation is not large enough to use PCTFS"); 316 if (PCTFS_my_id < dest) { 317 PetscCallMPI(MPI_Send(vals + segs[level - 1 - edge], stage_n, MPIU_SCALAR, dest, type, MPI_COMM_WORLD)); 318 } else { 319 type = type - PCTFS_my_id + dest; 320 PetscCallMPI(MPI_Recv(vals + segs[level - 1 - edge], stage_n, MPIU_SCALAR, MPI_ANY_SOURCE, type, MPI_COMM_WORLD, &status)); 321 } 322 } 323 mask >>= 1; 324 } 325 PetscFunctionReturn(PETSC_SUCCESS); 326 } 327 328 /***********************************comm.c*************************************/ 329 PetscErrorCode PCTFS_giop_hc(PetscInt *vals, PetscInt *work, PetscInt n, PetscInt *oprs, PetscInt dim) 330 { 331 PetscInt mask, edge; 332 PetscInt type, dest; 333 vfp fp; 334 MPI_Status status; 335 336 PetscFunctionBegin; 337 /* ok ... should have some data, work, and operator(s) */ 338 PetscCheck(vals && work && oprs, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop_hc() :: vals=%p, work=%p, oprs=%p", (void *)vals, (void *)work, (void *)oprs); 339 340 /* non-uniform should have at least two entries */ 341 PetscCheck(!(oprs[0] == NON_UNIFORM) || !(n < 2), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop_hc() :: non_uniform and n=0,1?"); 342 343 /* check to make sure comm package has been initialized */ 344 if (!p_init) PetscCall(PCTFS_comm_init()); 345 346 /* if there's nothing to do return */ 347 if ((PCTFS_num_nodes < 2) || (!n) || (dim <= 0)) PetscFunctionReturn(PETSC_SUCCESS); 348 349 /* the error msg says it all!!! */ 350 PetscCheck(!modfl_num_nodes, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop_hc() :: PCTFS_num_nodes not a power of 2!?!"); 351 352 /* a negative number of items to send ==> fatal */ 353 PetscCheck(n >= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop_hc() :: n=%" PetscInt_FMT "<0?", n); 354 355 /* can't do more dimensions then exist */ 356 dim = PetscMin(dim, PCTFS_i_log2_num_nodes); 357 358 /* advance to list of n operations for custom */ 359 if ((type = oprs[0]) == NON_UNIFORM) oprs++; 360 361 PetscCheck(fp = (vfp)PCTFS_ivec_fct_addr(type), PETSC_COMM_SELF, PETSC_ERR_PLIB, "PCTFS_giop_hc() :: Could not retrieve function pointer!"); 362 363 for (mask = 1, edge = 0; edge < dim; edge++, mask <<= 1) { 364 dest = PCTFS_my_id ^ mask; 365 if (PCTFS_my_id > dest) { 366 PetscCallMPI(MPI_Send(vals, n, MPIU_INT, dest, MSGTAG2 + PCTFS_my_id, MPI_COMM_WORLD)); 367 } else { 368 PetscCallMPI(MPI_Recv(work, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG2 + dest, MPI_COMM_WORLD, &status)); 369 PetscCall((*fp)(vals, work, n, oprs)); 370 } 371 } 372 373 if (edge == dim) mask >>= 1; 374 else { 375 while (++edge < dim) mask <<= 1; 376 } 377 378 for (edge = 0; edge < dim; edge++, mask >>= 1) { 379 if (PCTFS_my_id % mask) continue; 380 381 dest = PCTFS_my_id ^ mask; 382 if (PCTFS_my_id < dest) { 383 PetscCallMPI(MPI_Send(vals, n, MPIU_INT, dest, MSGTAG4 + PCTFS_my_id, MPI_COMM_WORLD)); 384 } else { 385 PetscCallMPI(MPI_Recv(vals, n, MPIU_INT, MPI_ANY_SOURCE, MSGTAG4 + dest, MPI_COMM_WORLD, &status)); 386 } 387 } 388 PetscFunctionReturn(PETSC_SUCCESS); 389 } 390