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 19 /* global program control variables - explicitly exported */ 20 PetscMPIInt my_id = 0; 21 PetscMPIInt num_nodes = 1; 22 PetscMPIInt floor_num_nodes = 0; 23 PetscMPIInt i_log2_num_nodes = 0; 24 25 /* global program control variables */ 26 static PetscInt p_init = 0; 27 static PetscInt modfl_num_nodes; 28 static PetscInt edge_not_pow_2; 29 30 static PetscInt edge_node[sizeof(PetscInt)*32]; 31 32 /***********************************comm.c*************************************/ 33 PetscErrorCode comm_init (void) 34 { 35 36 if (p_init++) PetscFunctionReturn(0); 37 38 MPI_Comm_size(MPI_COMM_WORLD,&num_nodes); 39 MPI_Comm_rank(MPI_COMM_WORLD,&my_id); 40 41 if (num_nodes> (INT_MAX >> 1)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Can't have more then MAX_INT/2 nodes!!!"); 42 43 ivec_zero((PetscInt*)edge_node,sizeof(PetscInt)*32); 44 45 floor_num_nodes = 1; 46 i_log2_num_nodes = modfl_num_nodes = 0; 47 while (floor_num_nodes <= num_nodes) 48 { 49 edge_node[i_log2_num_nodes] = my_id ^ floor_num_nodes; 50 floor_num_nodes <<= 1; 51 i_log2_num_nodes++; 52 } 53 54 i_log2_num_nodes--; 55 floor_num_nodes >>= 1; 56 modfl_num_nodes = (num_nodes - floor_num_nodes); 57 58 if ((my_id > 0) && (my_id <= modfl_num_nodes)) 59 {edge_not_pow_2=((my_id|floor_num_nodes)-1);} 60 else if (my_id >= floor_num_nodes) 61 {edge_not_pow_2=((my_id^floor_num_nodes)+1); 62 } 63 else 64 {edge_not_pow_2 = 0;} 65 PetscFunctionReturn(0); 66 } 67 68 /***********************************comm.c*************************************/ 69 PetscErrorCode giop(PetscInt *vals, PetscInt *work, PetscInt n, PetscInt *oprs) 70 { 71 PetscInt mask, edge; 72 PetscInt type, dest; 73 vfp fp; 74 MPI_Status status; 75 PetscInt ierr; 76 77 PetscFunctionBegin; 78 /* ok ... should have some data, work, and operator(s) */ 79 if (!vals||!work||!oprs) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_PLIB,"giop() :: vals=%D, work=%D, oprs=%D",vals,work,oprs); 80 81 /* non-uniform should have at least two entries */ 82 if ((oprs[0] == NON_UNIFORM)&&(n<2)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"giop() :: non_uniform and n=0,1?"); 83 84 /* check to make sure comm package has been initialized */ 85 if (!p_init) 86 {comm_init();} 87 88 /* if there's nothing to do return */ 89 if ((num_nodes<2)||(!n)) 90 { 91 PetscFunctionReturn(0); 92 } 93 94 95 /* a negative number if items to send ==> fatal */ 96 if (n<0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_PLIB,"giop() :: n=%D<0?",n); 97 98 /* advance to list of n operations for custom */ 99 if ((type=oprs[0])==NON_UNIFORM) 100 {oprs++;} 101 102 /* major league hack */ 103 if (!(fp = (vfp) ivec_fct_addr(type))) { 104 ierr = PetscInfo(0,"giop() :: hope you passed in a rbfp!\n");CHKERRQ(ierr); 105 fp = (vfp) oprs; 106 } 107 108 /* all msgs will be of the same length */ 109 /* if not a hypercube must colapse partial dim */ 110 if (edge_not_pow_2) 111 { 112 if (my_id >= floor_num_nodes) 113 {ierr = MPI_Send(vals,n,MPIU_INT,edge_not_pow_2,MSGTAG0+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 114 else 115 { 116 ierr = MPI_Recv(work,n,MPIU_INT,MPI_ANY_SOURCE,MSGTAG0+edge_not_pow_2, MPI_COMM_WORLD,&status);CHKERRQ(ierr); 117 (*fp)(vals,work,n,oprs); 118 } 119 } 120 121 /* implement the mesh fan in/out exchange algorithm */ 122 if (my_id<floor_num_nodes) 123 { 124 for (mask=1,edge=0; edge<i_log2_num_nodes; edge++,mask<<=1) 125 { 126 dest = my_id^mask; 127 if (my_id > dest) 128 {ierr = MPI_Send(vals,n,MPIU_INT,dest,MSGTAG2+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 129 else 130 { 131 ierr = MPI_Recv(work,n,MPIU_INT,MPI_ANY_SOURCE,MSGTAG2+dest,MPI_COMM_WORLD, &status);CHKERRQ(ierr); 132 (*fp)(vals, work, n, oprs); 133 } 134 } 135 136 mask=floor_num_nodes>>1; 137 for (edge=0; edge<i_log2_num_nodes; edge++,mask>>=1) 138 { 139 if (my_id%mask) 140 {continue;} 141 142 dest = my_id^mask; 143 if (my_id < dest) 144 {ierr = MPI_Send(vals,n,MPIU_INT,dest,MSGTAG4+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 145 else 146 { 147 ierr = MPI_Recv(vals,n,MPIU_INT,MPI_ANY_SOURCE,MSGTAG4+dest,MPI_COMM_WORLD, &status);CHKERRQ(ierr); 148 } 149 } 150 } 151 152 /* if not a hypercube must expand to partial dim */ 153 if (edge_not_pow_2) 154 { 155 if (my_id >= floor_num_nodes) 156 { 157 ierr = MPI_Recv(vals,n,MPIU_INT,MPI_ANY_SOURCE,MSGTAG5+edge_not_pow_2,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 158 } 159 else 160 {ierr = MPI_Send(vals,n,MPIU_INT,edge_not_pow_2,MSGTAG5+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 161 } 162 PetscFunctionReturn(0); 163 } 164 165 /***********************************comm.c*************************************/ 166 PetscErrorCode grop(PetscScalar *vals, PetscScalar *work, PetscInt n, PetscInt *oprs) 167 { 168 PetscInt mask, edge; 169 PetscInt type, dest; 170 vfp fp; 171 MPI_Status status; 172 PetscErrorCode ierr; 173 174 PetscFunctionBegin; 175 /* ok ... should have some data, work, and operator(s) */ 176 if (!vals||!work||!oprs) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_PLIB,"grop() :: vals=%D, work=%D, oprs=%D",vals,work,oprs); 177 178 /* non-uniform should have at least two entries */ 179 if ((oprs[0] == NON_UNIFORM)&&(n<2)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"grop() :: non_uniform and n=0,1?"); 180 181 /* check to make sure comm package has been initialized */ 182 if (!p_init) 183 {comm_init();} 184 185 /* if there's nothing to do return */ 186 if ((num_nodes<2)||(!n)) 187 { PetscFunctionReturn(0);} 188 189 /* a negative number of items to send ==> fatal */ 190 if (n<0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_PLIB,"gdop() :: n=%D<0?",n); 191 192 /* advance to list of n operations for custom */ 193 if ((type=oprs[0])==NON_UNIFORM) 194 {oprs++;} 195 196 if (!(fp = (vfp) rvec_fct_addr(type))) { 197 ierr = PetscInfo(0,"grop() :: hope you passed in a rbfp!\n");CHKERRQ(ierr); 198 fp = (vfp) oprs; 199 } 200 201 /* all msgs will be of the same length */ 202 /* if not a hypercube must colapse partial dim */ 203 if (edge_not_pow_2) 204 { 205 if (my_id >= floor_num_nodes) 206 {ierr = MPI_Send(vals,n,MPIU_SCALAR,edge_not_pow_2,MSGTAG0+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 207 else 208 { 209 ierr = MPI_Recv(work,n,MPIU_SCALAR,MPI_ANY_SOURCE,MSGTAG0+edge_not_pow_2,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 210 (*fp)(vals,work,n,oprs); 211 } 212 } 213 214 /* implement the mesh fan in/out exchange algorithm */ 215 if (my_id<floor_num_nodes) 216 { 217 for (mask=1,edge=0; edge<i_log2_num_nodes; edge++,mask<<=1) 218 { 219 dest = my_id^mask; 220 if (my_id > dest) 221 {ierr = MPI_Send(vals,n,MPIU_SCALAR,dest,MSGTAG2+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 222 else 223 { 224 ierr = MPI_Recv(work,n,MPIU_SCALAR,MPI_ANY_SOURCE,MSGTAG2+dest,MPI_COMM_WORLD, &status);CHKERRQ(ierr); 225 (*fp)(vals, work, n, oprs); 226 } 227 } 228 229 mask=floor_num_nodes>>1; 230 for (edge=0; edge<i_log2_num_nodes; edge++,mask>>=1) 231 { 232 if (my_id%mask) 233 {continue;} 234 235 dest = my_id^mask; 236 if (my_id < dest) 237 {ierr = MPI_Send(vals,n,MPIU_SCALAR,dest,MSGTAG4+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 238 else 239 { 240 ierr = MPI_Recv(vals,n,MPIU_SCALAR,MPI_ANY_SOURCE,MSGTAG4+dest,MPI_COMM_WORLD, &status);CHKERRQ(ierr); 241 } 242 } 243 } 244 245 /* if not a hypercube must expand to partial dim */ 246 if (edge_not_pow_2) 247 { 248 if (my_id >= floor_num_nodes) 249 { 250 ierr = MPI_Recv(vals,n,MPIU_SCALAR,MPI_ANY_SOURCE,MSGTAG5+edge_not_pow_2, MPI_COMM_WORLD,&status);CHKERRQ(ierr); 251 } 252 else 253 {ierr = MPI_Send(vals,n,MPIU_SCALAR,edge_not_pow_2,MSGTAG5+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 254 } 255 PetscFunctionReturn(0); 256 } 257 258 /***********************************comm.c*************************************/ 259 PetscErrorCode grop_hc(PetscScalar *vals, PetscScalar *work, PetscInt n, PetscInt *oprs, PetscInt dim) 260 { 261 PetscInt mask, edge; 262 PetscInt type, dest; 263 vfp fp; 264 MPI_Status status; 265 PetscErrorCode ierr; 266 267 PetscFunctionBegin; 268 /* ok ... should have some data, work, and operator(s) */ 269 if (!vals||!work||!oprs) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_PLIB,"grop_hc() :: vals=%D, work=%D, oprs=%D",vals,work,oprs); 270 271 /* non-uniform should have at least two entries */ 272 if ((oprs[0] == NON_UNIFORM)&&(n<2)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"grop_hc() :: non_uniform and n=0,1?"); 273 274 /* check to make sure comm package has been initialized */ 275 if (!p_init) 276 {comm_init();} 277 278 /* if there's nothing to do return */ 279 if ((num_nodes<2)||(!n)||(dim<=0)) 280 {PetscFunctionReturn(0);} 281 282 /* the error msg says it all!!! */ 283 if (modfl_num_nodes) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"grop_hc() :: num_nodes not a power of 2!?!"); 284 285 /* a negative number of items to send ==> fatal */ 286 if (n<0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_PLIB,"grop_hc() :: n=%D<0?",n); 287 288 /* can't do more dimensions then exist */ 289 dim = PetscMin(dim,i_log2_num_nodes); 290 291 /* advance to list of n operations for custom */ 292 if ((type=oprs[0])==NON_UNIFORM) 293 {oprs++;} 294 295 if (!(fp = (vfp) rvec_fct_addr(type))) { 296 ierr = PetscInfo(0,"grop_hc() :: hope you passed in a rbfp!\n");CHKERRQ(ierr); 297 fp = (vfp) oprs; 298 } 299 300 for (mask=1,edge=0; edge<dim; edge++,mask<<=1) 301 { 302 dest = my_id^mask; 303 if (my_id > dest) 304 {ierr = MPI_Send(vals,n,MPIU_SCALAR,dest,MSGTAG2+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 305 else 306 { 307 ierr = MPI_Recv(work,n,MPIU_SCALAR,MPI_ANY_SOURCE,MSGTAG2+dest,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 308 (*fp)(vals, work, n, oprs); 309 } 310 } 311 312 if (edge==dim) 313 {mask>>=1;} 314 else 315 {while (++edge<dim) {mask<<=1;}} 316 317 for (edge=0; edge<dim; edge++,mask>>=1) 318 { 319 if (my_id%mask) 320 {continue;} 321 322 dest = my_id^mask; 323 if (my_id < dest) 324 {ierr = MPI_Send(vals,n,MPIU_SCALAR,dest,MSGTAG4+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 325 else 326 { 327 ierr = MPI_Recv(vals,n,MPIU_SCALAR,MPI_ANY_SOURCE,MSGTAG4+dest,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 328 } 329 } 330 PetscFunctionReturn(0); 331 } 332 333 /******************************************************************************/ 334 PetscErrorCode ssgl_radd( PetscScalar *vals, PetscScalar *work, PetscInt level, PetscInt *segs) 335 { 336 PetscInt edge, type, dest, mask; 337 PetscInt stage_n; 338 MPI_Status status; 339 PetscErrorCode ierr; 340 341 PetscFunctionBegin; 342 /* check to make sure comm package has been initialized */ 343 if (!p_init) 344 {comm_init();} 345 346 347 /* all msgs are *NOT* the same length */ 348 /* implement the mesh fan in/out exchange algorithm */ 349 for (mask=0, edge=0; edge<level; edge++, mask++) 350 { 351 stage_n = (segs[level] - segs[edge]); 352 if (stage_n && !(my_id & mask)) 353 { 354 dest = edge_node[edge]; 355 type = MSGTAG3 + my_id + (num_nodes*edge); 356 if (my_id>dest) 357 {ierr = MPI_Send(vals+segs[edge],stage_n,MPIU_SCALAR,dest,type, MPI_COMM_WORLD);CHKERRQ(ierr);} 358 else 359 { 360 type = type - my_id + dest; 361 ierr = MPI_Recv(work,stage_n,MPIU_SCALAR,MPI_ANY_SOURCE,type,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 362 rvec_add(vals+segs[edge], work, stage_n); 363 } 364 } 365 mask <<= 1; 366 } 367 mask>>=1; 368 for (edge=0; edge<level; edge++) 369 { 370 stage_n = (segs[level] - segs[level-1-edge]); 371 if (stage_n && !(my_id & mask)) 372 { 373 dest = edge_node[level-edge-1]; 374 type = MSGTAG6 + my_id + (num_nodes*edge); 375 if (my_id<dest) 376 {ierr = MPI_Send(vals+segs[level-1-edge],stage_n,MPIU_SCALAR,dest,type,MPI_COMM_WORLD);CHKERRQ(ierr);} 377 else 378 { 379 type = type - my_id + dest; 380 ierr = MPI_Recv(vals+segs[level-1-edge],stage_n,MPIU_SCALAR, MPI_ANY_SOURCE,type,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 381 } 382 } 383 mask >>= 1; 384 } 385 PetscFunctionReturn(0); 386 } 387 388 /******************************************************************************/ 389 PetscErrorCode new_ssgl_radd( PetscScalar *vals, PetscScalar *work, PetscInt level, PetscInt *segs) 390 { 391 PetscInt edge, type, dest, mask; 392 PetscInt stage_n; 393 MPI_Status status; 394 PetscErrorCode ierr; 395 396 PetscFunctionBegin; 397 /* check to make sure comm package has been initialized */ 398 if (!p_init) 399 {comm_init();} 400 401 /* all msgs are *NOT* the same length */ 402 /* implement the mesh fan in/out exchange algorithm */ 403 for (mask=0, edge=0; edge<level; edge++, mask++) 404 { 405 stage_n = (segs[level] - segs[edge]); 406 if (stage_n && !(my_id & mask)) 407 { 408 dest = edge_node[edge]; 409 type = MSGTAG3 + my_id + (num_nodes*edge); 410 if (my_id>dest) 411 {ierr = MPI_Send(vals+segs[edge],stage_n,MPIU_SCALAR,dest,type, MPI_COMM_WORLD);CHKERRQ(ierr);} 412 else 413 { 414 type = type - my_id + dest; 415 ierr = MPI_Recv(work,stage_n,MPIU_SCALAR,MPI_ANY_SOURCE,type, MPI_COMM_WORLD,&status);CHKERRQ(ierr); 416 rvec_add(vals+segs[edge], work, stage_n); 417 } 418 } 419 mask <<= 1; 420 } 421 mask>>=1; 422 for (edge=0; edge<level; edge++) 423 { 424 stage_n = (segs[level] - segs[level-1-edge]); 425 if (stage_n && !(my_id & mask)) 426 { 427 dest = edge_node[level-edge-1]; 428 type = MSGTAG6 + my_id + (num_nodes*edge); 429 if (my_id<dest) 430 {ierr = MPI_Send(vals+segs[level-1-edge],stage_n,MPIU_SCALAR,dest,type,MPI_COMM_WORLD);CHKERRQ(ierr);} 431 else 432 { 433 type = type - my_id + dest; 434 ierr = MPI_Recv(vals+segs[level-1-edge],stage_n,MPIU_SCALAR, MPI_ANY_SOURCE,type,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 435 } 436 } 437 mask >>= 1; 438 } 439 PetscFunctionReturn(0); 440 } 441 442 /***********************************comm.c*************************************/ 443 PetscErrorCode giop_hc(PetscInt *vals, PetscInt *work, PetscInt n, PetscInt *oprs, PetscInt dim) 444 { 445 PetscInt mask, edge; 446 PetscInt type, dest; 447 vfp fp; 448 MPI_Status status; 449 PetscErrorCode ierr; 450 451 PetscFunctionBegin; 452 /* ok ... should have some data, work, and operator(s) */ 453 if (!vals||!work||!oprs) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_PLIB,"giop_hc() :: vals=%D, work=%D, oprs=%D",vals,work,oprs); 454 455 /* non-uniform should have at least two entries */ 456 if ((oprs[0] == NON_UNIFORM)&&(n<2)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"giop_hc() :: non_uniform and n=0,1?"); 457 458 /* check to make sure comm package has been initialized */ 459 if (!p_init) 460 {comm_init();} 461 462 /* if there's nothing to do return */ 463 if ((num_nodes<2)||(!n)||(dim<=0)) 464 { PetscFunctionReturn(0);} 465 466 /* the error msg says it all!!! */ 467 if (modfl_num_nodes) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"giop_hc() :: num_nodes not a power of 2!?!"); 468 469 /* a negative number of items to send ==> fatal */ 470 if (n<0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_PLIB,"giop_hc() :: n=%D<0?",n); 471 472 /* can't do more dimensions then exist */ 473 dim = PetscMin(dim,i_log2_num_nodes); 474 475 /* advance to list of n operations for custom */ 476 if ((type=oprs[0])==NON_UNIFORM) 477 {oprs++;} 478 479 if (!(fp = (vfp) ivec_fct_addr(type))){ 480 ierr = PetscInfo(0,"giop_hc() :: hope you passed in a rbfp!\n");CHKERRQ(ierr); 481 fp = (vfp) oprs; 482 } 483 484 for (mask=1,edge=0; edge<dim; edge++,mask<<=1) 485 { 486 dest = my_id^mask; 487 if (my_id > dest) 488 {ierr = MPI_Send(vals,n,MPIU_INT,dest,MSGTAG2+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 489 else 490 { 491 ierr = MPI_Recv(work,n,MPIU_INT,MPI_ANY_SOURCE,MSGTAG2+dest,MPI_COMM_WORLD, &status);CHKERRQ(ierr); 492 (*fp)(vals, work, n, oprs); 493 } 494 } 495 496 if (edge==dim) 497 {mask>>=1;} 498 else 499 {while (++edge<dim) {mask<<=1;}} 500 501 for (edge=0; edge<dim; edge++,mask>>=1) 502 { 503 if (my_id%mask) 504 {continue;} 505 506 dest = my_id^mask; 507 if (my_id < dest) 508 {ierr = MPI_Send(vals,n,MPIU_INT,dest,MSGTAG4+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 509 else 510 { 511 ierr = MPI_Recv(vals,n,MPIU_INT,MPI_ANY_SOURCE,MSGTAG4+dest,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 512 } 513 } 514 PetscFunctionReturn(0); 515 } 516