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