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