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