Actual source code: swarm_ex1.c
1: static char help[] = "Tests DMSwarm\n\n";
3: #include <petscdm.h>
4: #include <petscdmda.h>
5: #include <petscdmswarm.h>
7: PETSC_EXTERN PetscErrorCode DMSwarmCollect_General(DM, PetscErrorCode (*)(DM, void *, PetscInt *, PetscInt **), size_t, void *, PetscInt *);
8: PETSC_EXTERN PetscErrorCode DMSwarmCollect_DMDABoundingBox(DM, PetscInt *);
10: PetscErrorCode ex1_1(void)
11: {
12: DM dms;
13: Vec x;
14: PetscMPIInt rank, size;
15: PetscInt p;
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");
22: PetscCall(DMCreate(PETSC_COMM_WORLD, &dms));
23: PetscCall(DMSetType(dms, DMSWARM));
24: PetscCall(PetscObjectSetName((PetscObject)dms, "Particles"));
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));
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: }
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: }
47: PetscCall(DMSwarmCreateGlobalVectorFromField(dms, "viscosity", &x));
48: PetscCall(DMSwarmDestroyGlobalVectorFromField(dms, "viscosity", &x));
50: PetscCall(DMSwarmVectorDefineField(dms, "strain"));
51: PetscCall(DMCreateGlobalVector(dms, &x));
52: PetscCall(VecDestroy(&x));
54: {
55: PetscInt *rankval;
56: PetscInt npoints[2], npoints_orig[2];
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: }
86: PetscCall(DMDestroy(&dms));
87: PetscFunctionReturn(PETSC_SUCCESS);
88: }
90: PetscErrorCode ex1_2(void)
91: {
92: DM dms;
93: Vec x;
94: PetscMPIInt rank, size;
95: PetscInt p;
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));
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];
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));
133: PetscCall(DMSwarmGetField(dms, "DMSwarm_rank", NULL, NULL, (void **)&rankval));
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));
150: PetscCall(DMSwarmCreateGlobalVectorFromField(dms, "viscosity", &x));
151: PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD));
152: PetscCall(DMSwarmDestroyGlobalVectorFromField(dms, "viscosity", &x));
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));
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: }
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;
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));
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));
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;
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;
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;
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));
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];
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: }
287: typedef struct {
288: PetscReal cx[2];
289: PetscReal radius;
290: } CollectZoneCtx;
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;
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));
306: r2 = zone->radius * zone->radius;
307: p2collect = 0;
308: for (p = 0; p < npoints; p++) {
309: PetscReal sep2;
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: }
316: PetscCall(PetscMalloc1(p2collect + 1, &plist));
317: p2collect = 0;
318: for (p = 0; p < npoints; p++) {
319: PetscReal sep2;
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));
331: *nfound = p2collect;
332: *foundlist = plist;
333: PetscFunctionReturn(PETSC_SUCCESS);
334: }
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;
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"));
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));
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;
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;
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;
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];
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: }
478: int main(int argc, char **argv)
479: {
480: PetscInt test_mode = 4;
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: }
498: /*TEST
500: build:
501: requires: !complex double
503: test:
504: args: -test_mode 1
505: filter: grep -v atomic
506: filter_output: grep -v atomic
508: test:
509: suffix: 2
510: args: -test_mode 2
511: filter: grep -v atomic
512: filter_output: grep -v atomic
514: test:
515: suffix: 3
516: args: -test_mode 3
517: filter: grep -v atomic
518: filter_output: grep -v atomic
519: TODO: broken
521: test:
522: suffix: 4
523: args: -test_mode 4
524: filter: grep -v atomic
525: filter_output: grep -v atomic
527: test:
528: suffix: 5
529: nsize: 4
530: args: -test_mode 1
531: filter: grep -v atomic
532: filter_output: grep -v atomic
534: test:
535: suffix: 6
536: nsize: 4
537: args: -test_mode 2
538: filter: grep -v atomic
539: filter_output: grep -v atomic
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
549: test:
550: suffix: 8
551: nsize: 4
552: args: -test_mode 4
553: filter: grep -v atomic
554: filter_output: grep -v atomic
556: TEST*/