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