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