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