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