1 #define PETSCKSP_DLL 2 3 /***********************************comm.c************************************* 4 5 Author: Henry M. Tufo III 6 7 e-mail: hmt@cs.brown.edu 8 9 snail-mail: 10 Division of Applied Mathematics 11 Brown University 12 Providence, RI 02912 13 14 Last Modification: 15 11.21.97 16 ***********************************comm.c*************************************/ 17 #include "src/ksp/pc/impls/tfs/tfs.h" 18 19 20 /* global program control variables - explicitly exported */ 21 PetscMPIInt my_id = 0; 22 PetscMPIInt num_nodes = 1; 23 PetscInt floor_num_nodes = 0; 24 PetscInt i_log2_num_nodes = 0; 25 26 /* global program control variables */ 27 static PetscInt p_init = 0; 28 static PetscInt modfl_num_nodes; 29 static PetscInt edge_not_pow_2; 30 31 static PetscInt edge_node[sizeof(PetscInt)*32]; 32 33 /***********************************comm.c*************************************/ 34 PetscErrorCode comm_init (void) 35 { 36 37 if (p_init++) PetscFunctionReturn(0); 38 39 MPI_Comm_size(MPI_COMM_WORLD,&num_nodes); 40 MPI_Comm_rank(MPI_COMM_WORLD,&my_id); 41 42 if (num_nodes> (INT_MAX >> 1)) 43 {error_msg_fatal("Can't have more then MAX_INT/2 nodes!!!");} 44 45 ivec_zero((PetscInt*)edge_node,sizeof(PetscInt)*32); 46 47 floor_num_nodes = 1; 48 i_log2_num_nodes = modfl_num_nodes = 0; 49 while (floor_num_nodes <= num_nodes) 50 { 51 edge_node[i_log2_num_nodes] = my_id ^ floor_num_nodes; 52 floor_num_nodes <<= 1; 53 i_log2_num_nodes++; 54 } 55 56 i_log2_num_nodes--; 57 floor_num_nodes >>= 1; 58 modfl_num_nodes = (num_nodes - floor_num_nodes); 59 60 if ((my_id > 0) && (my_id <= modfl_num_nodes)) 61 {edge_not_pow_2=((my_id|floor_num_nodes)-1);} 62 else if (my_id >= floor_num_nodes) 63 {edge_not_pow_2=((my_id^floor_num_nodes)+1); 64 } 65 else 66 {edge_not_pow_2 = 0;} 67 PetscFunctionReturn(0); 68 } 69 70 /***********************************comm.c*************************************/ 71 PetscErrorCode giop(PetscInt *vals, PetscInt *work, PetscInt n, PetscInt *oprs) 72 { 73 PetscInt mask, edge; 74 PetscInt type, dest; 75 vfp fp; 76 MPI_Status status; 77 PetscInt ierr; 78 79 PetscFunctionBegin; 80 /* ok ... should have some data, work, and operator(s) */ 81 if (!vals||!work||!oprs) 82 {error_msg_fatal("giop() :: vals=%D, work=%D, oprs=%D",vals,work,oprs);} 83 84 /* non-uniform should have at least two entries */ 85 if ((oprs[0] == NON_UNIFORM)&&(n<2)) 86 {error_msg_fatal("giop() :: non_uniform and n=0,1?");} 87 88 /* check to make sure comm package has been initialized */ 89 if (!p_init) 90 {comm_init();} 91 92 /* if there's nothing to do return */ 93 if ((num_nodes<2)||(!n)) 94 { 95 PetscFunctionReturn(0); 96 } 97 98 /* a negative number if items to send ==> fatal */ 99 if (n<0) 100 {error_msg_fatal("giop() :: n=%D<0?",n);} 101 102 /* advance to list of n operations for custom */ 103 if ((type=oprs[0])==NON_UNIFORM) 104 {oprs++;} 105 106 /* major league hack */ 107 if (!(fp = (vfp) ivec_fct_addr(type))) { 108 error_msg_warning("giop() :: hope you passed in a rbfp!\n"); 109 fp = (vfp) oprs; 110 } 111 112 /* all msgs will be of the same length */ 113 /* if not a hypercube must colapse partial dim */ 114 if (edge_not_pow_2) 115 { 116 if (my_id >= floor_num_nodes) 117 {ierr = MPI_Send(vals,n,MPIU_INT,edge_not_pow_2,MSGTAG0+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 118 else 119 { 120 ierr = MPI_Recv(work,n,MPIU_INT,MPI_ANY_SOURCE,MSGTAG0+edge_not_pow_2, MPI_COMM_WORLD,&status);CHKERRQ(ierr); 121 (*fp)(vals,work,n,oprs); 122 } 123 } 124 125 /* implement the mesh fan in/out exchange algorithm */ 126 if (my_id<floor_num_nodes) 127 { 128 for (mask=1,edge=0; edge<i_log2_num_nodes; edge++,mask<<=1) 129 { 130 dest = my_id^mask; 131 if (my_id > dest) 132 {ierr = MPI_Send(vals,n,MPIU_INT,dest,MSGTAG2+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 133 else 134 { 135 ierr = MPI_Recv(work,n,MPIU_INT,MPI_ANY_SOURCE,MSGTAG2+dest,MPI_COMM_WORLD, &status);CHKERRQ(ierr); 136 (*fp)(vals, work, n, oprs); 137 } 138 } 139 140 mask=floor_num_nodes>>1; 141 for (edge=0; edge<i_log2_num_nodes; edge++,mask>>=1) 142 { 143 if (my_id%mask) 144 {continue;} 145 146 dest = my_id^mask; 147 if (my_id < dest) 148 {ierr = MPI_Send(vals,n,MPIU_INT,dest,MSGTAG4+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 149 else 150 { 151 ierr = MPI_Recv(vals,n,MPIU_INT,MPI_ANY_SOURCE,MSGTAG4+dest,MPI_COMM_WORLD, &status);CHKERRQ(ierr); 152 } 153 } 154 } 155 156 /* if not a hypercube must expand to partial dim */ 157 if (edge_not_pow_2) 158 { 159 if (my_id >= floor_num_nodes) 160 { 161 ierr = MPI_Recv(vals,n,MPIU_INT,MPI_ANY_SOURCE,MSGTAG5+edge_not_pow_2,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 162 } 163 else 164 {ierr = MPI_Send(vals,n,MPIU_INT,edge_not_pow_2,MSGTAG5+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 165 } 166 PetscFunctionReturn(0); 167 } 168 169 /***********************************comm.c*************************************/ 170 PetscErrorCode grop(PetscScalar *vals, PetscScalar *work, PetscInt n, int *oprs) 171 { 172 PetscInt mask, edge; 173 PetscInt type, dest; 174 vfp fp; 175 MPI_Status status; 176 PetscErrorCode ierr; 177 178 PetscFunctionBegin; 179 /* ok ... should have some data, work, and operator(s) */ 180 if (!vals||!work||!oprs) 181 {error_msg_fatal("grop() :: vals=%D, work=%D, oprs=%D",vals,work,oprs);} 182 183 /* non-uniform should have at least two entries */ 184 if ((oprs[0] == NON_UNIFORM)&&(n<2)) 185 {error_msg_fatal("grop() :: non_uniform and n=0,1?");} 186 187 /* check to make sure comm package has been initialized */ 188 if (!p_init) 189 {comm_init();} 190 191 /* if there's nothing to do return */ 192 if ((num_nodes<2)||(!n)) 193 { PetscFunctionReturn(0);} 194 195 /* a negative number of items to send ==> fatal */ 196 if (n<0) 197 {error_msg_fatal("gdop() :: n=%D<0?",n);} 198 199 /* advance to list of n operations for custom */ 200 if ((type=oprs[0])==NON_UNIFORM) 201 {oprs++;} 202 203 if (!(fp = (vfp) rvec_fct_addr(type))) { 204 error_msg_warning("grop() :: hope you passed in a rbfp!\n"); 205 fp = (vfp) oprs; 206 } 207 208 /* all msgs will be of the same length */ 209 /* if not a hypercube must colapse partial dim */ 210 if (edge_not_pow_2) 211 { 212 if (my_id >= floor_num_nodes) 213 {ierr = MPI_Send(vals,n,MPIU_SCALAR,edge_not_pow_2,MSGTAG0+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 214 else 215 { 216 ierr = MPI_Recv(work,n,MPIU_SCALAR,MPI_ANY_SOURCE,MSGTAG0+edge_not_pow_2,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 217 (*fp)(vals,work,n,oprs); 218 } 219 } 220 221 /* implement the mesh fan in/out exchange algorithm */ 222 if (my_id<floor_num_nodes) 223 { 224 for (mask=1,edge=0; edge<i_log2_num_nodes; edge++,mask<<=1) 225 { 226 dest = my_id^mask; 227 if (my_id > dest) 228 {ierr = MPI_Send(vals,n,MPIU_SCALAR,dest,MSGTAG2+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 229 else 230 { 231 ierr = MPI_Recv(work,n,MPIU_SCALAR,MPI_ANY_SOURCE,MSGTAG2+dest,MPI_COMM_WORLD, &status);CHKERRQ(ierr); 232 (*fp)(vals, work, n, oprs); 233 } 234 } 235 236 mask=floor_num_nodes>>1; 237 for (edge=0; edge<i_log2_num_nodes; edge++,mask>>=1) 238 { 239 if (my_id%mask) 240 {continue;} 241 242 dest = my_id^mask; 243 if (my_id < dest) 244 {ierr = MPI_Send(vals,n,MPIU_SCALAR,dest,MSGTAG4+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 245 else 246 { 247 ierr = MPI_Recv(vals,n,MPIU_SCALAR,MPI_ANY_SOURCE,MSGTAG4+dest,MPI_COMM_WORLD, &status);CHKERRQ(ierr); 248 } 249 } 250 } 251 252 /* if not a hypercube must expand to partial dim */ 253 if (edge_not_pow_2) 254 { 255 if (my_id >= floor_num_nodes) 256 { 257 ierr = MPI_Recv(vals,n,MPIU_SCALAR,MPI_ANY_SOURCE,MSGTAG5+edge_not_pow_2, MPI_COMM_WORLD,&status);CHKERRQ(ierr); 258 } 259 else 260 {ierr = MPI_Send(vals,n,MPIU_SCALAR,edge_not_pow_2,MSGTAG5+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 261 } 262 PetscFunctionReturn(0); 263 } 264 265 /***********************************comm.c*************************************/ 266 PetscErrorCode grop_hc(PetscScalar *vals, PetscScalar *work, PetscInt n, PetscInt *oprs, PetscInt dim) 267 { 268 PetscInt mask, edge; 269 PetscInt type, dest; 270 vfp fp; 271 MPI_Status status; 272 PetscErrorCode ierr; 273 274 PetscFunctionBegin; 275 /* ok ... should have some data, work, and operator(s) */ 276 if (!vals||!work||!oprs) 277 {error_msg_fatal("grop_hc() :: vals=%D, work=%D, oprs=%D",vals,work,oprs);} 278 279 /* non-uniform should have at least two entries */ 280 if ((oprs[0] == NON_UNIFORM)&&(n<2)) 281 {error_msg_fatal("grop_hc() :: non_uniform and n=0,1?");} 282 283 /* check to make sure comm package has been initialized */ 284 if (!p_init) 285 {comm_init();} 286 287 /* if there's nothing to do return */ 288 if ((num_nodes<2)||(!n)||(dim<=0)) 289 {PetscFunctionReturn(0);} 290 291 /* the error msg says it all!!! */ 292 if (modfl_num_nodes) 293 {error_msg_fatal("grop_hc() :: num_nodes not a power of 2!?!");} 294 295 /* a negative number of items to send ==> fatal */ 296 if (n<0) 297 {error_msg_fatal("grop_hc() :: n=%D<0?",n);} 298 299 /* can't do more dimensions then exist */ 300 dim = PetscMin(dim,i_log2_num_nodes); 301 302 /* advance to list of n operations for custom */ 303 if ((type=oprs[0])==NON_UNIFORM) 304 {oprs++;} 305 306 if (!(fp = (vfp) rvec_fct_addr(type))) { 307 error_msg_warning("grop_hc() :: hope you passed in a rbfp!\n"); 308 fp = (vfp) oprs; 309 } 310 311 for (mask=1,edge=0; edge<dim; edge++,mask<<=1) 312 { 313 dest = my_id^mask; 314 if (my_id > dest) 315 {ierr = MPI_Send(vals,n,MPIU_SCALAR,dest,MSGTAG2+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 316 else 317 { 318 ierr = MPI_Recv(work,n,MPIU_SCALAR,MPI_ANY_SOURCE,MSGTAG2+dest,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 319 (*fp)(vals, work, n, oprs); 320 } 321 } 322 323 if (edge==dim) 324 {mask>>=1;} 325 else 326 {while (++edge<dim) {mask<<=1;}} 327 328 for (edge=0; edge<dim; edge++,mask>>=1) 329 { 330 if (my_id%mask) 331 {continue;} 332 333 dest = my_id^mask; 334 if (my_id < dest) 335 {ierr = MPI_Send(vals,n,MPIU_SCALAR,dest,MSGTAG4+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 336 else 337 { 338 ierr = MPI_Recv(vals,n,MPIU_SCALAR,MPI_ANY_SOURCE,MSGTAG4+dest,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 339 } 340 } 341 PetscFunctionReturn(0); 342 } 343 344 /***********************************comm.c*************************************/ 345 PetscErrorCode gfop(void *vals, void *work, PetscInt n, vbfp fp, MPI_Datatype dt, PetscInt comm_type) 346 { 347 PetscInt mask, edge; 348 PetscInt dest; 349 MPI_Status status; 350 MPI_Op op; 351 PetscErrorCode ierr; 352 353 PetscFunctionBegin; 354 /* check to make sure comm package has been initialized */ 355 if (!p_init) 356 {comm_init();} 357 358 /* ok ... should have some data, work, and operator(s) */ 359 if (!vals||!work||!fp) 360 {error_msg_fatal("gop() :: v=%D, w=%D, f=%D",vals,work,fp);} 361 362 /* if there's nothing to do return */ 363 if ((num_nodes<2)||(!n)) 364 {PetscFunctionReturn(0);} 365 366 /* a negative number of items to send ==> fatal */ 367 if (n<0) 368 {error_msg_fatal("gop() :: n=%D<0?",n);} 369 370 if (comm_type==MPI) 371 { 372 ierr = MPI_Op_create(fp,TRUE,&op);CHKERRQ(ierr); 373 ierr = MPI_Allreduce (vals, work, n, dt, op, MPI_COMM_WORLD);CHKERRQ(ierr); 374 ierr = MPI_Op_free(&op);CHKERRQ(ierr); 375 } 376 377 /* if not a hypercube must colapse partial dim */ 378 if (edge_not_pow_2) 379 { 380 if (my_id >= floor_num_nodes) 381 {ierr = MPI_Send(vals,n,dt,edge_not_pow_2,MSGTAG0+my_id, MPI_COMM_WORLD);CHKERRQ(ierr);} 382 else 383 { 384 ierr = MPI_Recv(work,n,dt,MPI_ANY_SOURCE,MSGTAG0+edge_not_pow_2,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 385 (*fp)(vals,work,&n,&dt); 386 } 387 } 388 389 /* implement the mesh fan in/out exchange algorithm */ 390 if (my_id<floor_num_nodes) 391 { 392 for (mask=1,edge=0; edge<i_log2_num_nodes; edge++,mask<<=1) 393 { 394 dest = my_id^mask; 395 if (my_id > dest) 396 {ierr = MPI_Send(vals,n,dt,dest,MSGTAG2+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 397 else 398 { 399 ierr = MPI_Recv(work,n,dt,MPI_ANY_SOURCE,MSGTAG2+dest,MPI_COMM_WORLD, &status);CHKERRQ(ierr); 400 (*fp)(vals, work, &n, &dt); 401 } 402 } 403 404 mask=floor_num_nodes>>1; 405 for (edge=0; edge<i_log2_num_nodes; edge++,mask>>=1) 406 { 407 if (my_id%mask) 408 {continue;} 409 410 dest = my_id^mask; 411 if (my_id < dest) 412 {ierr = MPI_Send(vals,n,dt,dest,MSGTAG4+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 413 else 414 { 415 ierr = MPI_Recv(vals,n,dt,MPI_ANY_SOURCE,MSGTAG4+dest, MPI_COMM_WORLD, &status);CHKERRQ(ierr); 416 } 417 } 418 } 419 /* if not a hypercube must expand to partial dim */ 420 if (edge_not_pow_2) 421 { 422 if (my_id >= floor_num_nodes) 423 { 424 ierr = MPI_Recv(vals,n,dt,MPI_ANY_SOURCE,MSGTAG5+edge_not_pow_2, MPI_COMM_WORLD,&status);CHKERRQ(ierr); 425 } 426 else 427 {ierr = MPI_Send(vals,n,dt,edge_not_pow_2,MSGTAG5+my_id, MPI_COMM_WORLD);CHKERRQ(ierr);} 428 } 429 PetscFunctionReturn(0); 430 } 431 432 /******************************************************************************/ 433 PetscErrorCode ssgl_radd( PetscScalar *vals, PetscScalar *work, PetscInt level, PetscInt *segs) 434 { 435 PetscInt edge, type, dest, mask; 436 PetscInt stage_n; 437 MPI_Status status; 438 PetscErrorCode ierr; 439 440 PetscFunctionBegin; 441 /* check to make sure comm package has been initialized */ 442 if (!p_init) 443 {comm_init();} 444 445 446 /* all msgs are *NOT* the same length */ 447 /* implement the mesh fan in/out exchange algorithm */ 448 for (mask=0, edge=0; edge<level; edge++, mask++) 449 { 450 stage_n = (segs[level] - segs[edge]); 451 if (stage_n && !(my_id & mask)) 452 { 453 dest = edge_node[edge]; 454 type = MSGTAG3 + my_id + (num_nodes*edge); 455 if (my_id>dest) 456 {ierr = MPI_Send(vals+segs[edge],stage_n,MPIU_SCALAR,dest,type, MPI_COMM_WORLD);CHKERRQ(ierr);} 457 else 458 { 459 type = type - my_id + dest; 460 ierr = MPI_Recv(work,stage_n,MPIU_SCALAR,MPI_ANY_SOURCE,type,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 461 rvec_add(vals+segs[edge], work, stage_n); 462 } 463 } 464 mask <<= 1; 465 } 466 mask>>=1; 467 for (edge=0; edge<level; edge++) 468 { 469 stage_n = (segs[level] - segs[level-1-edge]); 470 if (stage_n && !(my_id & mask)) 471 { 472 dest = edge_node[level-edge-1]; 473 type = MSGTAG6 + my_id + (num_nodes*edge); 474 if (my_id<dest) 475 {ierr = MPI_Send(vals+segs[level-1-edge],stage_n,MPIU_SCALAR,dest,type,MPI_COMM_WORLD);CHKERRQ(ierr);} 476 else 477 { 478 type = type - my_id + dest; 479 ierr = MPI_Recv(vals+segs[level-1-edge],stage_n,MPIU_SCALAR, MPI_ANY_SOURCE,type,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 480 } 481 } 482 mask >>= 1; 483 } 484 PetscFunctionReturn(0); 485 } 486 487 /******************************************************************************/ 488 PetscErrorCode new_ssgl_radd( PetscScalar *vals, PetscScalar *work, PetscInt level, PetscInt *segs) 489 { 490 PetscInt edge, type, dest, mask; 491 PetscInt stage_n; 492 MPI_Status status; 493 PetscErrorCode ierr; 494 495 PetscFunctionBegin; 496 /* check to make sure comm package has been initialized */ 497 if (!p_init) 498 {comm_init();} 499 500 /* all msgs are *NOT* the same length */ 501 /* implement the mesh fan in/out exchange algorithm */ 502 for (mask=0, edge=0; edge<level; edge++, mask++) 503 { 504 stage_n = (segs[level] - segs[edge]); 505 if (stage_n && !(my_id & mask)) 506 { 507 dest = edge_node[edge]; 508 type = MSGTAG3 + my_id + (num_nodes*edge); 509 if (my_id>dest) 510 {ierr = MPI_Send(vals+segs[edge],stage_n,MPIU_SCALAR,dest,type, MPI_COMM_WORLD);CHKERRQ(ierr);} 511 else 512 { 513 type = type - my_id + dest; 514 ierr = MPI_Recv(work,stage_n,MPIU_SCALAR,MPI_ANY_SOURCE,type, MPI_COMM_WORLD,&status);CHKERRQ(ierr); 515 rvec_add(vals+segs[edge], work, stage_n); 516 } 517 } 518 mask <<= 1; 519 } 520 mask>>=1; 521 for (edge=0; edge<level; edge++) 522 { 523 stage_n = (segs[level] - segs[level-1-edge]); 524 if (stage_n && !(my_id & mask)) 525 { 526 dest = edge_node[level-edge-1]; 527 type = MSGTAG6 + my_id + (num_nodes*edge); 528 if (my_id<dest) 529 {ierr = MPI_Send(vals+segs[level-1-edge],stage_n,MPIU_SCALAR,dest,type,MPI_COMM_WORLD);CHKERRQ(ierr);} 530 else 531 { 532 type = type - my_id + dest; 533 ierr = MPI_Recv(vals+segs[level-1-edge],stage_n,MPIU_SCALAR, MPI_ANY_SOURCE,type,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 534 } 535 } 536 mask >>= 1; 537 } 538 PetscFunctionReturn(0); 539 } 540 541 /***********************************comm.c*************************************/ 542 PetscErrorCode giop_hc(PetscInt *vals, PetscInt *work, PetscInt n, PetscInt *oprs, PetscInt dim) 543 { 544 PetscInt mask, edge; 545 PetscInt type, dest; 546 vfp fp; 547 MPI_Status status; 548 PetscErrorCode ierr; 549 550 PetscFunctionBegin; 551 /* ok ... should have some data, work, and operator(s) */ 552 if (!vals||!work||!oprs) 553 {error_msg_fatal("giop_hc() :: vals=%D, work=%D, oprs=%D",vals,work,oprs);} 554 555 /* non-uniform should have at least two entries */ 556 if ((oprs[0] == NON_UNIFORM)&&(n<2)) 557 {error_msg_fatal("giop_hc() :: non_uniform and n=0,1?");} 558 559 /* check to make sure comm package has been initialized */ 560 if (!p_init) 561 {comm_init();} 562 563 /* if there's nothing to do return */ 564 if ((num_nodes<2)||(!n)||(dim<=0)) 565 { PetscFunctionReturn(0);} 566 567 /* the error msg says it all!!! */ 568 if (modfl_num_nodes) 569 {error_msg_fatal("giop_hc() :: num_nodes not a power of 2!?!");} 570 571 /* a negative number of items to send ==> fatal */ 572 if (n<0) 573 {error_msg_fatal("giop_hc() :: n=%D<0?",n);} 574 575 /* can't do more dimensions then exist */ 576 dim = PetscMin(dim,i_log2_num_nodes); 577 578 /* advance to list of n operations for custom */ 579 if ((type=oprs[0])==NON_UNIFORM) 580 {oprs++;} 581 582 if (!(fp = (vfp) ivec_fct_addr(type))){ 583 error_msg_warning("giop_hc() :: hope you passed in a rbfp!\n"); 584 fp = (vfp) oprs; 585 } 586 587 for (mask=1,edge=0; edge<dim; edge++,mask<<=1) 588 { 589 dest = my_id^mask; 590 if (my_id > dest) 591 {ierr = MPI_Send(vals,n,MPIU_INT,dest,MSGTAG2+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 592 else 593 { 594 ierr = MPI_Recv(work,n,MPIU_INT,MPI_ANY_SOURCE,MSGTAG2+dest,MPI_COMM_WORLD, &status);CHKERRQ(ierr); 595 (*fp)(vals, work, n, oprs); 596 } 597 } 598 599 if (edge==dim) 600 {mask>>=1;} 601 else 602 {while (++edge<dim) {mask<<=1;}} 603 604 for (edge=0; edge<dim; edge++,mask>>=1) 605 { 606 if (my_id%mask) 607 {continue;} 608 609 dest = my_id^mask; 610 if (my_id < dest) 611 {ierr = MPI_Send(vals,n,MPIU_INT,dest,MSGTAG4+my_id,MPI_COMM_WORLD);CHKERRQ(ierr);} 612 else 613 { 614 ierr = MPI_Recv(vals,n,MPIU_INT,MPI_ANY_SOURCE,MSGTAG4+dest,MPI_COMM_WORLD,&status);CHKERRQ(ierr); 615 } 616 } 617 PetscFunctionReturn(0); 618 } 619