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