1 2 static char help[] = "Tests DMSwarm\n\n"; 3 4 #include <petscdm.h> 5 #include <petscdmda.h> 6 #include <petscdmswarm.h> 7 8 PETSC_EXTERN PetscErrorCode DMSwarmCollect_General(DM,PetscErrorCode (*)(DM,void*,PetscInt*,PetscInt**),size_t,void*,PetscInt*); 9 PETSC_EXTERN PetscErrorCode DMSwarmCollect_DMDABoundingBox(DM,PetscInt*); 10 11 PetscErrorCode ex1_1(void) 12 { 13 DM dms; 14 Vec x; 15 PetscMPIInt rank,size; 16 PetscInt p; 17 18 PetscFunctionBegin; 19 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 20 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 21 PetscCheck(!(size > 1) || !(size != 4),PETSC_COMM_WORLD,PETSC_ERR_SUP,"Must be run wuth 4 MPI ranks"); 22 23 PetscCall(DMCreate(PETSC_COMM_WORLD,&dms)); 24 PetscCall(DMSetType(dms,DMSWARM)); 25 PetscCall(PetscObjectSetName((PetscObject) dms, "Particles")); 26 27 PetscCall(DMSwarmInitializeFieldRegister(dms)); 28 PetscCall(DMSwarmRegisterPetscDatatypeField(dms,"viscosity",1,PETSC_REAL)); 29 PetscCall(DMSwarmRegisterPetscDatatypeField(dms,"strain",1,PETSC_REAL)); 30 PetscCall(DMSwarmFinalizeFieldRegister(dms)); 31 PetscCall(DMSwarmSetLocalSizes(dms,5+rank,4)); 32 PetscCall(DMView(dms,PETSC_VIEWER_STDOUT_WORLD)); 33 34 { 35 PetscReal *array; 36 PetscCall(DMSwarmGetField(dms,"viscosity",NULL,NULL,(void**)&array)); 37 for (p=0; p<5+rank; p++) { 38 array[p] = 11.1 + p*0.1 + rank*100.0; 39 } 40 PetscCall(DMSwarmRestoreField(dms,"viscosity",NULL,NULL,(void**)&array)); 41 } 42 43 { 44 PetscReal *array; 45 PetscCall(DMSwarmGetField(dms,"strain",NULL,NULL,(void**)&array)); 46 for (p=0; p<5+rank; p++) { 47 array[p] = 2.0e-2 + p*0.001 + rank*1.0; 48 } 49 PetscCall(DMSwarmRestoreField(dms,"strain",NULL,NULL,(void**)&array)); 50 } 51 52 PetscCall(DMSwarmCreateGlobalVectorFromField(dms,"viscosity",&x)); 53 PetscCall(DMSwarmDestroyGlobalVectorFromField(dms,"viscosity",&x)); 54 55 PetscCall(DMSwarmVectorDefineField(dms,"strain")); 56 PetscCall(DMCreateGlobalVector(dms,&x)); 57 PetscCall(VecDestroy(&x)); 58 59 { 60 PetscInt *rankval; 61 PetscInt npoints[2],npoints_orig[2]; 62 63 PetscCall(DMSwarmGetLocalSize(dms,&npoints_orig[0])); 64 PetscCall(DMSwarmGetSize(dms,&npoints_orig[1])); 65 PetscCall(DMSwarmGetField(dms,"DMSwarm_rank",NULL,NULL,(void**)&rankval)); 66 if ((rank == 0) && (size > 1)) { 67 rankval[0] = 1; 68 rankval[3] = 1; 69 } 70 if (rank == 3) { 71 rankval[2] = 1; 72 } 73 PetscCall(DMSwarmRestoreField(dms,"DMSwarm_rank",NULL,NULL,(void**)&rankval)); 74 PetscCall(DMSwarmMigrate(dms,PETSC_TRUE)); 75 PetscCall(DMSwarmGetLocalSize(dms,&npoints[0])); 76 PetscCall(DMSwarmGetSize(dms,&npoints[1])); 77 PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 78 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"rank[%d] before(%" PetscInt_FMT ",%" PetscInt_FMT ") after(%" PetscInt_FMT ",%" PetscInt_FMT ")\n",rank,npoints_orig[0],npoints_orig[1],npoints[0],npoints[1])); 79 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 80 PetscCall(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 81 } 82 { 83 PetscCall(DMSwarmCreateGlobalVectorFromField(dms,"viscosity",&x)); 84 PetscCall(VecView(x,PETSC_VIEWER_STDOUT_WORLD)); 85 PetscCall(DMSwarmDestroyGlobalVectorFromField(dms,"viscosity",&x)); 86 } 87 { 88 PetscCall(DMSwarmCreateGlobalVectorFromField(dms,"strain",&x)); 89 PetscCall(VecView(x,PETSC_VIEWER_STDOUT_WORLD)); 90 PetscCall(DMSwarmDestroyGlobalVectorFromField(dms,"strain",&x)); 91 } 92 93 PetscCall(DMDestroy(&dms)); 94 PetscFunctionReturn(0); 95 } 96 97 PetscErrorCode ex1_2(void) 98 { 99 DM dms; 100 Vec x; 101 PetscMPIInt rank,size; 102 PetscInt p; 103 104 PetscFunctionBegin; 105 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 106 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 107 PetscCall(DMCreate(PETSC_COMM_WORLD,&dms)); 108 PetscCall(DMSetType(dms,DMSWARM)); 109 PetscCall(PetscObjectSetName((PetscObject) dms, "Particles")); 110 PetscCall(DMSwarmInitializeFieldRegister(dms)); 111 112 PetscCall(DMSwarmRegisterPetscDatatypeField(dms,"viscosity",1,PETSC_REAL)); 113 PetscCall(DMSwarmRegisterPetscDatatypeField(dms,"strain",1,PETSC_REAL)); 114 PetscCall(DMSwarmFinalizeFieldRegister(dms)); 115 PetscCall(DMSwarmSetLocalSizes(dms,5+rank,4)); 116 PetscCall(DMView(dms,PETSC_VIEWER_STDOUT_WORLD)); 117 { 118 PetscReal *array; 119 PetscCall(DMSwarmGetField(dms,"viscosity",NULL,NULL,(void**)&array)); 120 for (p=0; p<5+rank; p++) { 121 array[p] = 11.1 + p*0.1 + rank*100.0; 122 } 123 PetscCall(DMSwarmRestoreField(dms,"viscosity",NULL,NULL,(void**)&array)); 124 } 125 { 126 PetscReal *array; 127 PetscCall(DMSwarmGetField(dms,"strain",NULL,NULL,(void**)&array)); 128 for (p=0; p<5+rank; p++) { 129 array[p] = 2.0e-2 + p*0.001 + rank*1.0; 130 } 131 PetscCall(DMSwarmRestoreField(dms,"strain",NULL,NULL,(void**)&array)); 132 } 133 { 134 PetscInt *rankval; 135 PetscInt npoints[2],npoints_orig[2]; 136 137 PetscCall(DMSwarmGetLocalSize(dms,&npoints_orig[0])); 138 PetscCall(DMSwarmGetSize(dms,&npoints_orig[1])); 139 PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 140 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"rank[%d] before(%" PetscInt_FMT ",%" PetscInt_FMT ")\n",rank,npoints_orig[0],npoints_orig[1])); 141 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 142 PetscCall(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 143 144 PetscCall(DMSwarmGetField(dms,"DMSwarm_rank",NULL,NULL,(void**)&rankval)); 145 146 if (rank == 1) { 147 rankval[0] = -1; 148 } 149 if (rank == 2) { 150 rankval[1] = -1; 151 } 152 if (rank == 3) { 153 rankval[3] = -1; 154 rankval[4] = -1; 155 } 156 PetscCall(DMSwarmRestoreField(dms,"DMSwarm_rank",NULL,NULL,(void**)&rankval)); 157 PetscCall(DMSwarmCollectViewCreate(dms)); 158 PetscCall(DMSwarmGetLocalSize(dms,&npoints[0])); 159 PetscCall(DMSwarmGetSize(dms,&npoints[1])); 160 PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 161 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"rank[%d] after(%" PetscInt_FMT ",%" PetscInt_FMT ")\n",rank,npoints[0],npoints[1])); 162 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 163 PetscCall(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 164 165 PetscCall(DMSwarmCreateGlobalVectorFromField(dms,"viscosity",&x)); 166 PetscCall(VecView(x,PETSC_VIEWER_STDOUT_WORLD)); 167 PetscCall(DMSwarmDestroyGlobalVectorFromField(dms,"viscosity",&x)); 168 169 PetscCall(DMSwarmCollectViewDestroy(dms)); 170 PetscCall(DMSwarmGetLocalSize(dms,&npoints[0])); 171 PetscCall(DMSwarmGetSize(dms,&npoints[1])); 172 PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 173 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"rank[%d] after_v(%" PetscInt_FMT ",%" PetscInt_FMT ")\n",rank,npoints[0],npoints[1])); 174 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 175 PetscCall(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 176 177 PetscCall(DMSwarmCreateGlobalVectorFromField(dms,"viscosity",&x)); 178 PetscCall(VecView(x,PETSC_VIEWER_STDOUT_WORLD)); 179 PetscCall(DMSwarmDestroyGlobalVectorFromField(dms,"viscosity",&x)); 180 } 181 PetscCall(DMDestroy(&dms)); 182 PetscFunctionReturn(0); 183 } 184 185 /* 186 splot "c-rank0.gp","c-rank1.gp","c-rank2.gp","c-rank3.gp" 187 */ 188 PetscErrorCode ex1_3(void) 189 { 190 DM dms; 191 PetscMPIInt rank,size; 192 PetscInt is,js,ni,nj,overlap; 193 DM dmcell; 194 195 PetscFunctionBegin; 196 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 197 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 198 overlap = 2; 199 PetscCall(DMDACreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,13,13,PETSC_DECIDE,PETSC_DECIDE,1,overlap,NULL,NULL,&dmcell)); 200 PetscCall(DMSetFromOptions(dmcell)); 201 PetscCall(DMSetUp(dmcell)); 202 PetscCall(DMDASetUniformCoordinates(dmcell,-1.0,1.0,-1.0,1.0,-1.0,1.0)); 203 PetscCall(DMDAGetCorners(dmcell,&is,&js,NULL,&ni,&nj,NULL)); 204 PetscCall(DMCreate(PETSC_COMM_WORLD,&dms)); 205 PetscCall(DMSetType(dms,DMSWARM)); 206 PetscCall(PetscObjectSetName((PetscObject) dms, "Particles")); 207 PetscCall(DMSwarmSetCellDM(dms,dmcell)); 208 209 /* load in data types */ 210 PetscCall(DMSwarmInitializeFieldRegister(dms)); 211 PetscCall(DMSwarmRegisterPetscDatatypeField(dms,"viscosity",1,PETSC_REAL)); 212 PetscCall(DMSwarmRegisterPetscDatatypeField(dms,"coorx",1,PETSC_REAL)); 213 PetscCall(DMSwarmRegisterPetscDatatypeField(dms,"coory",1,PETSC_REAL)); 214 PetscCall(DMSwarmFinalizeFieldRegister(dms)); 215 PetscCall(DMSwarmSetLocalSizes(dms,ni*nj*4,4)); 216 PetscCall(DMView(dms,PETSC_VIEWER_STDOUT_WORLD)); 217 218 /* set values within the swarm */ 219 { 220 PetscReal *array_x,*array_y; 221 PetscInt npoints,i,j,cnt; 222 DMDACoor2d **LA_coor; 223 Vec coor; 224 DM dmcellcdm; 225 226 PetscCall(DMGetCoordinateDM(dmcell,&dmcellcdm)); 227 PetscCall(DMGetCoordinates(dmcell,&coor)); 228 PetscCall(DMDAVecGetArray(dmcellcdm,coor,&LA_coor)); 229 PetscCall(DMSwarmGetLocalSize(dms,&npoints)); 230 PetscCall(DMSwarmGetField(dms,"coorx",NULL,NULL,(void**)&array_x)); 231 PetscCall(DMSwarmGetField(dms,"coory",NULL,NULL,(void**)&array_y)); 232 cnt = 0; 233 for (j=js; j<js+nj; j++) { 234 for (i=is; i<is+ni; i++) { 235 PetscReal xp,yp; 236 xp = PetscRealPart(LA_coor[j][i].x); 237 yp = PetscRealPart(LA_coor[j][i].y); 238 array_x[4*cnt+0] = xp - 0.05; if (array_x[4*cnt+0] < -1.0) { array_x[4*cnt+0] = -1.0+1.0e-12; } 239 array_x[4*cnt+1] = xp + 0.05; if (array_x[4*cnt+1] > 1.0) { array_x[4*cnt+1] = 1.0-1.0e-12; } 240 array_x[4*cnt+2] = xp - 0.05; if (array_x[4*cnt+2] < -1.0) { array_x[4*cnt+2] = -1.0+1.0e-12; } 241 array_x[4*cnt+3] = xp + 0.05; if (array_x[4*cnt+3] > 1.0) { array_x[4*cnt+3] = 1.0-1.0e-12; } 242 243 array_y[4*cnt+0] = yp - 0.05; if (array_y[4*cnt+0] < -1.0) { array_y[4*cnt+0] = -1.0+1.0e-12; } 244 array_y[4*cnt+1] = yp - 0.05; if (array_y[4*cnt+1] < -1.0) { array_y[4*cnt+1] = -1.0+1.0e-12; } 245 array_y[4*cnt+2] = yp + 0.05; if (array_y[4*cnt+2] > 1.0) { array_y[4*cnt+2] = 1.0-1.0e-12; } 246 array_y[4*cnt+3] = yp + 0.05; if (array_y[4*cnt+3] > 1.0) { array_y[4*cnt+3] = 1.0-1.0e-12; } 247 cnt++; 248 } 249 } 250 PetscCall(DMSwarmRestoreField(dms,"coory",NULL,NULL,(void**)&array_y)); 251 PetscCall(DMSwarmRestoreField(dms,"coorx",NULL,NULL,(void**)&array_x)); 252 PetscCall(DMDAVecRestoreArray(dmcellcdm,coor,&LA_coor)); 253 } 254 { 255 PetscInt npoints[2],npoints_orig[2],ng; 256 257 PetscCall(DMSwarmGetLocalSize(dms,&npoints_orig[0])); 258 PetscCall(DMSwarmGetSize(dms,&npoints_orig[1])); 259 PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 260 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"rank[%d] before(%" PetscInt_FMT ",%" PetscInt_FMT ")\n",rank,npoints_orig[0],npoints_orig[1])); 261 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 262 PetscCall(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 263 PetscCall(DMSwarmCollect_DMDABoundingBox(dms,&ng)); 264 265 PetscCall(DMSwarmGetLocalSize(dms,&npoints[0])); 266 PetscCall(DMSwarmGetSize(dms,&npoints[1])); 267 PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 268 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"rank[%d] after(%" PetscInt_FMT ",%" PetscInt_FMT ")\n",rank,npoints[0],npoints[1])); 269 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 270 PetscCall(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 271 } 272 { 273 PetscReal *array_x,*array_y; 274 PetscInt npoints,p; 275 FILE *fp = NULL; 276 char name[PETSC_MAX_PATH_LEN]; 277 278 PetscCall(PetscSNPrintf(name,PETSC_MAX_PATH_LEN-1,"c-rank%d.gp",rank)); 279 fp = fopen(name,"w"); 280 PetscCheck(fp,PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN,"Cannot open file %s",name); 281 PetscCall(DMSwarmGetLocalSize(dms,&npoints)); 282 PetscCall(DMSwarmGetField(dms,"coorx",NULL,NULL,(void**)&array_x)); 283 PetscCall(DMSwarmGetField(dms,"coory",NULL,NULL,(void**)&array_y)); 284 for (p=0; p<npoints; p++) { 285 fprintf(fp,"%+1.4e %+1.4e %1.4e\n",array_x[p],array_y[p],(double)rank); 286 } 287 PetscCall(DMSwarmRestoreField(dms,"coory",NULL,NULL,(void**)&array_y)); 288 PetscCall(DMSwarmRestoreField(dms,"coorx",NULL,NULL,(void**)&array_x)); 289 fclose(fp); 290 } 291 PetscCall(DMDestroy(&dmcell)); 292 PetscCall(DMDestroy(&dms)); 293 PetscFunctionReturn(0); 294 } 295 296 typedef struct { 297 PetscReal cx[2]; 298 PetscReal radius; 299 } CollectZoneCtx; 300 301 PetscErrorCode collect_zone(DM dm,void *ctx,PetscInt *nfound,PetscInt **foundlist) 302 { 303 CollectZoneCtx *zone = (CollectZoneCtx*)ctx; 304 PetscInt p,npoints; 305 PetscReal *array_x,*array_y,r2; 306 PetscInt p2collect,*plist; 307 PetscMPIInt rank; 308 309 PetscFunctionBegin; 310 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 311 PetscCall(DMSwarmGetLocalSize(dm,&npoints)); 312 PetscCall(DMSwarmGetField(dm,"coorx",NULL,NULL,(void**)&array_x)); 313 PetscCall(DMSwarmGetField(dm,"coory",NULL,NULL,(void**)&array_y)); 314 315 r2 = zone->radius * zone->radius; 316 p2collect = 0; 317 for (p=0; p<npoints; p++) { 318 PetscReal sep2; 319 320 sep2 = (array_x[p] - zone->cx[0])*(array_x[p] - zone->cx[0]); 321 sep2 += (array_y[p] - zone->cx[1])*(array_y[p] - zone->cx[1]); 322 if (sep2 < r2) { 323 p2collect++; 324 } 325 } 326 327 PetscCall(PetscMalloc1(p2collect+1,&plist)); 328 p2collect = 0; 329 for (p=0; p<npoints; p++) { 330 PetscReal sep2; 331 332 sep2 = (array_x[p] - zone->cx[0])*(array_x[p] - zone->cx[0]); 333 sep2 += (array_y[p] - zone->cx[1])*(array_y[p] - zone->cx[1]); 334 if (sep2 < r2) { 335 plist[p2collect] = p; 336 p2collect++; 337 } 338 } 339 PetscCall(DMSwarmRestoreField(dm,"coory",NULL,NULL,(void**)&array_y)); 340 PetscCall(DMSwarmRestoreField(dm,"coorx",NULL,NULL,(void**)&array_x)); 341 342 *nfound = p2collect; 343 *foundlist = plist; 344 PetscFunctionReturn(0); 345 } 346 347 PetscErrorCode ex1_4(void) 348 { 349 DM dms; 350 PetscMPIInt rank,size; 351 PetscInt is,js,ni,nj,overlap,nn; 352 DM dmcell; 353 CollectZoneCtx *zone; 354 PetscReal dx; 355 356 PetscFunctionBegin; 357 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 358 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 359 nn = 101; 360 dx = 2.0/ (PetscReal)(nn-1); 361 overlap = 0; 362 PetscCall(DMDACreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,nn,nn,PETSC_DECIDE,PETSC_DECIDE,1,overlap,NULL,NULL,&dmcell)); 363 PetscCall(DMSetFromOptions(dmcell)); 364 PetscCall(DMSetUp(dmcell)); 365 PetscCall(DMDASetUniformCoordinates(dmcell,-1.0,1.0,-1.0,1.0,-1.0,1.0)); 366 PetscCall(DMDAGetCorners(dmcell,&is,&js,NULL,&ni,&nj,NULL)); 367 PetscCall(DMCreate(PETSC_COMM_WORLD,&dms)); 368 PetscCall(DMSetType(dms,DMSWARM)); 369 PetscCall(PetscObjectSetName((PetscObject) dms, "Particles")); 370 371 /* load in data types */ 372 PetscCall(DMSwarmInitializeFieldRegister(dms)); 373 PetscCall(DMSwarmRegisterPetscDatatypeField(dms,"viscosity",1,PETSC_REAL)); 374 PetscCall(DMSwarmRegisterPetscDatatypeField(dms,"coorx",1,PETSC_REAL)); 375 PetscCall(DMSwarmRegisterPetscDatatypeField(dms,"coory",1,PETSC_REAL)); 376 PetscCall(DMSwarmFinalizeFieldRegister(dms)); 377 PetscCall(DMSwarmSetLocalSizes(dms,ni*nj*4,4)); 378 PetscCall(DMView(dms,PETSC_VIEWER_STDOUT_WORLD)); 379 380 /* set values within the swarm */ 381 { 382 PetscReal *array_x,*array_y; 383 PetscInt npoints,i,j,cnt; 384 DMDACoor2d **LA_coor; 385 Vec coor; 386 DM dmcellcdm; 387 388 PetscCall(DMGetCoordinateDM(dmcell,&dmcellcdm)); 389 PetscCall(DMGetCoordinates(dmcell,&coor)); 390 PetscCall(DMDAVecGetArray(dmcellcdm,coor,&LA_coor)); 391 PetscCall(DMSwarmGetLocalSize(dms,&npoints)); 392 PetscCall(DMSwarmGetField(dms,"coorx",NULL,NULL,(void**)&array_x)); 393 PetscCall(DMSwarmGetField(dms,"coory",NULL,NULL,(void**)&array_y)); 394 cnt = 0; 395 for (j=js; j<js+nj; j++) { 396 for (i=is; i<is+ni; i++) { 397 PetscReal xp,yp; 398 399 xp = PetscRealPart(LA_coor[j][i].x); 400 yp = PetscRealPart(LA_coor[j][i].y); 401 array_x[4*cnt+0] = xp - dx*0.1; /*if (array_x[4*cnt+0] < -1.0) { array_x[4*cnt+0] = -1.0+1.0e-12; }*/ 402 array_x[4*cnt+1] = xp + dx*0.1; /*if (array_x[4*cnt+1] > 1.0) { array_x[4*cnt+1] = 1.0-1.0e-12; }*/ 403 array_x[4*cnt+2] = xp - dx*0.1; /*if (array_x[4*cnt+2] < -1.0) { array_x[4*cnt+2] = -1.0+1.0e-12; }*/ 404 array_x[4*cnt+3] = xp + dx*0.1; /*if (array_x[4*cnt+3] > 1.0) { array_x[4*cnt+3] = 1.0-1.0e-12; }*/ 405 array_y[4*cnt+0] = yp - dx*0.1; /*if (array_y[4*cnt+0] < -1.0) { array_y[4*cnt+0] = -1.0+1.0e-12; }*/ 406 array_y[4*cnt+1] = yp - dx*0.1; /*if (array_y[4*cnt+1] < -1.0) { array_y[4*cnt+1] = -1.0+1.0e-12; }*/ 407 array_y[4*cnt+2] = yp + dx*0.1; /*if (array_y[4*cnt+2] > 1.0) { array_y[4*cnt+2] = 1.0-1.0e-12; }*/ 408 array_y[4*cnt+3] = yp + dx*0.1; /*if (array_y[4*cnt+3] > 1.0) { array_y[4*cnt+3] = 1.0-1.0e-12; }*/ 409 cnt++; 410 } 411 } 412 PetscCall(DMSwarmRestoreField(dms,"coory",NULL,NULL,(void**)&array_y)); 413 PetscCall(DMSwarmRestoreField(dms,"coorx",NULL,NULL,(void**)&array_x)); 414 PetscCall(DMDAVecRestoreArray(dmcellcdm,coor,&LA_coor)); 415 } 416 PetscCall(PetscMalloc1(1,&zone)); 417 if (size == 4) { 418 if (rank == 0) { 419 zone->cx[0] = 0.5; 420 zone->cx[1] = 0.5; 421 zone->radius = 0.3; 422 } 423 if (rank == 1) { 424 zone->cx[0] = -0.5; 425 zone->cx[1] = 0.5; 426 zone->radius = 0.25; 427 } 428 if (rank == 2) { 429 zone->cx[0] = 0.5; 430 zone->cx[1] = -0.5; 431 zone->radius = 0.2; 432 } 433 if (rank == 3) { 434 zone->cx[0] = -0.5; 435 zone->cx[1] = -0.5; 436 zone->radius = 0.1; 437 } 438 } else { 439 if (rank == 0) { 440 zone->cx[0] = 0.5; 441 zone->cx[1] = 0.5; 442 zone->radius = 0.8; 443 } else { 444 zone->cx[0] = 10.0; 445 zone->cx[1] = 10.0; 446 zone->radius = 0.0; 447 } 448 } 449 { 450 PetscInt npoints[2],npoints_orig[2],ng; 451 452 PetscCall(DMSwarmGetLocalSize(dms,&npoints_orig[0])); 453 PetscCall(DMSwarmGetSize(dms,&npoints_orig[1])); 454 PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 455 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"rank[%d] before(%" PetscInt_FMT ",%" PetscInt_FMT ")\n",rank,npoints_orig[0],npoints_orig[1])); 456 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 457 PetscCall(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 458 PetscCall(DMSwarmCollect_General(dms,collect_zone,sizeof(CollectZoneCtx),zone,&ng)); 459 PetscCall(DMSwarmGetLocalSize(dms,&npoints[0])); 460 PetscCall(DMSwarmGetSize(dms,&npoints[1])); 461 PetscCall(PetscViewerASCIIPushSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 462 PetscCall(PetscViewerASCIISynchronizedPrintf(PETSC_VIEWER_STDOUT_WORLD,"rank[%d] after(%" PetscInt_FMT ",%" PetscInt_FMT ")\n",rank,npoints[0],npoints[1])); 463 PetscCall(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 464 PetscCall(PetscViewerASCIIPopSynchronized(PETSC_VIEWER_STDOUT_WORLD)); 465 } 466 { 467 PetscReal *array_x,*array_y; 468 PetscInt npoints,p; 469 FILE *fp = NULL; 470 char name[PETSC_MAX_PATH_LEN]; 471 472 PetscCall(PetscSNPrintf(name,PETSC_MAX_PATH_LEN-1,"c-rank%d.gp",rank)); 473 fp = fopen(name,"w"); 474 PetscCheck(fp,PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN,"Cannot open file %s",name); 475 PetscCall(DMSwarmGetLocalSize(dms,&npoints)); 476 PetscCall(DMSwarmGetField(dms,"coorx",NULL,NULL,(void**)&array_x)); 477 PetscCall(DMSwarmGetField(dms,"coory",NULL,NULL,(void**)&array_y)); 478 for (p=0; p<npoints; p++) { 479 fprintf(fp,"%+1.4e %+1.4e %1.4e\n",array_x[p],array_y[p],(double)rank); 480 } 481 PetscCall(DMSwarmRestoreField(dms,"coory",NULL,NULL,(void**)&array_y)); 482 PetscCall(DMSwarmRestoreField(dms,"coorx",NULL,NULL,(void**)&array_x)); 483 fclose(fp); 484 } 485 PetscCall(DMDestroy(&dmcell)); 486 PetscCall(DMDestroy(&dms)); 487 PetscCall(PetscFree(zone)); 488 PetscFunctionReturn(0); 489 } 490 491 int main(int argc,char **argv) 492 { 493 PetscInt test_mode = 4; 494 495 PetscCall(PetscInitialize(&argc,&argv,(char*)0,help)); 496 PetscCall(PetscOptionsGetInt(NULL,NULL,"-test_mode",&test_mode,NULL)); 497 if (test_mode == 1) { 498 PetscCall(ex1_1()); 499 } else if (test_mode == 2) { 500 PetscCall(ex1_2()); 501 } else if (test_mode == 3) { 502 PetscCall(ex1_3()); 503 } else if (test_mode == 4) { 504 PetscCall(ex1_4()); 505 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_USER,"Unknown test_mode value, should be 1,2,3,4"); 506 PetscCall(PetscFinalize()); 507 return 0; 508 } 509 510 /*TEST 511 512 build: 513 requires: !complex double 514 515 test: 516 args: -test_mode 1 517 filter: grep -v atomic 518 filter_output: grep -v atomic 519 520 test: 521 suffix: 2 522 args: -test_mode 2 523 filter: grep -v atomic 524 filter_output: grep -v atomic 525 526 test: 527 suffix: 3 528 args: -test_mode 3 529 filter: grep -v atomic 530 filter_output: grep -v atomic 531 TODO: broken 532 533 test: 534 suffix: 4 535 args: -test_mode 4 536 filter: grep -v atomic 537 filter_output: grep -v atomic 538 539 test: 540 suffix: 5 541 nsize: 4 542 args: -test_mode 1 543 filter: grep -v atomic 544 filter_output: grep -v atomic 545 546 test: 547 suffix: 6 548 nsize: 4 549 args: -test_mode 2 550 filter: grep -v atomic 551 filter_output: grep -v atomic 552 553 test: 554 suffix: 7 555 nsize: 4 556 args: -test_mode 3 557 filter: grep -v atomic 558 filter_output: grep -v atomic 559 TODO: broken 560 561 test: 562 suffix: 8 563 nsize: 4 564 args: -test_mode 4 565 filter: grep -v atomic 566 filter_output: grep -v atomic 567 568 TEST*/ 569