Actual source code: swarmpic_view.c
petsc-3.10.5 2019-03-28
1: #include <petscdm.h>
2: #include <petscdmda.h>
3: #include <petscdmswarm.h>
4: #include <petsc/private/dmswarmimpl.h>
5: #include "../src/dm/impls/swarm/data_bucket.h"
7: PetscErrorCode private_PetscViewerCreate_XDMF(MPI_Comm comm,const char filename[],PetscViewer *v)
8: {
9: long int *bytes;
10: PetscContainer container;
11: PetscViewer viewer;
13:
15: PetscViewerCreate(comm,&viewer);
16: PetscViewerSetType(viewer,PETSCVIEWERASCII);
17: PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
18: PetscViewerFileSetName(viewer,filename);
19:
20: PetscMalloc1(1,&bytes);
21: bytes[0] = 0;
22: PetscContainerCreate(comm,&container);
23: PetscContainerSetPointer(container,(void*)bytes);
24: PetscObjectCompose((PetscObject)viewer,"XDMFViewerContext",(PetscObject)container);
25:
26: /* write xdmf header */
27: PetscViewerASCIIPrintf(viewer,"<?xml version=\"1.0\" encoding=\"utf-8\"?>\n");
28: PetscViewerASCIIPrintf(viewer,"<Xdmf xmlns:xi=\"http://www.w3.org/2001/XInclude\" Version=\"2.99\">\n");
29: /* write xdmf domain */
30: PetscViewerASCIIPushTab(viewer);
31: PetscViewerASCIIPrintf(viewer,"<Domain>\n");
32: *v = viewer;
33: return(0);
34: }
36: PetscErrorCode private_PetscViewerDestroy_XDMF(PetscViewer *v)
37: {
38: PetscViewer viewer;
39: DM dm = NULL;
40: long int *bytes;
41: PetscContainer container = NULL;
43:
45: if (!v) return(0);
46: viewer = *v;
47:
48: PetscObjectQuery((PetscObject)viewer,"DMSwarm",(PetscObject*)&dm);
49: if (dm) {
50: PetscViewerASCIIPrintf(viewer,"</Grid>\n");
51: PetscViewerASCIIPopTab(viewer);
52: }
53:
54: /* close xdmf header */
55: PetscViewerASCIIPrintf(viewer,"</Domain>\n");
56: PetscViewerASCIIPopTab(viewer);
57: PetscViewerASCIIPrintf(viewer,"</Xdmf>\n");
58:
59: PetscObjectQuery((PetscObject)viewer,"XDMFViewerContext",(PetscObject*)&container);
60: if (container) {
61: PetscContainerGetPointer(container,(void**)&bytes);
62: PetscFree(bytes);
63: PetscContainerDestroy(&container);
64: }
65:
66: PetscViewerDestroy(&viewer);
67: *v = NULL;
68: return(0);
69: }
71: PetscErrorCode private_CreateDataFileNameXDMF(const char filename[],char dfilename[])
72: {
73: char *ext;
74: PetscBool flg;
76:
78: PetscStrrchr(filename,'.',&ext);
79: PetscStrcmp("xmf",ext,&flg);
80: if (flg) {
81: size_t len;
82: char viewername_minus_ext[PETSC_MAX_PATH_LEN];
83:
84: PetscStrlen(filename,&len);
85: PetscStrncpy(viewername_minus_ext,filename,len-2);
86: PetscSNPrintf(dfilename,PETSC_MAX_PATH_LEN-1,"%s_swarm_fields.pbin",viewername_minus_ext);
87: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"File extension must by .xmf");
89: return(0);
90: }
92: PetscErrorCode private_DMSwarmView_XDMF(DM dm,PetscViewer viewer)
93: {
94: PetscBool isswarm = PETSC_FALSE;
95: const char *viewername;
96: char datafile[PETSC_MAX_PATH_LEN];
97: PetscViewer fviewer;
98: PetscInt k,ng,dim;
99: Vec dvec;
100: long int *bytes = NULL;
101: PetscContainer container = NULL;
102: const char *dmname;
104:
106: PetscObjectQuery((PetscObject)viewer,"XDMFViewerContext",(PetscObject*)&container);
107: if (container) {
108: PetscContainerGetPointer(container,(void**)&bytes);
109: } else SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Valid to find attached data XDMFViewerContext");
110:
111: PetscObjectTypeCompare((PetscObject)dm,DMSWARM,&isswarm);
112: if (!isswarm) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Only valid for DMSwarm");
113:
114: PetscObjectCompose((PetscObject)viewer,"DMSwarm",(PetscObject)dm);
115:
116: PetscViewerASCIIPushTab(viewer);
117: PetscObjectGetName((PetscObject)dm,&dmname);
118: if (!dmname) {
119: DMGetOptionsPrefix(dm,&dmname);
120: }
121: if (!dmname) {
122: PetscViewerASCIIPrintf(viewer,"<Grid Name=\"DMSwarm\" GridType=\"Uniform\">\n");
123: } else {
124: PetscViewerASCIIPrintf(viewer,"<Grid Name=\"DMSwarm[%s]\" GridType=\"Uniform\">\n",dmname);
125: }
126:
127: /* create a sub-viewer for topology, geometry and all data fields */
128: /* name is viewer.name + "_swarm_fields.pbin" */
129: PetscViewerCreate(PetscObjectComm((PetscObject)viewer),&fviewer);
130: PetscViewerSetType(fviewer,PETSCVIEWERBINARY);
131: PetscViewerBinarySetSkipHeader(fviewer,PETSC_TRUE);
132: PetscViewerBinarySetSkipInfo(fviewer,PETSC_TRUE);
133: PetscViewerFileSetMode(fviewer,FILE_MODE_WRITE);
134:
135: PetscViewerFileGetName(viewer,&viewername);
136: private_CreateDataFileNameXDMF(viewername,datafile);
137: PetscViewerFileSetName(fviewer,datafile);
138:
139: DMSwarmGetSize(dm,&ng);
140:
141: /* write topology header */
142: PetscViewerASCIIPushTab(viewer);
143: PetscViewerASCIIPrintf(viewer,"<Topology Dimensions=\"%D\" TopologyType=\"Mixed\">\n",ng);
144: PetscViewerASCIIPushTab(viewer);
145: PetscViewerASCIIPrintf(viewer,"<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Int\" Dimensions=\"%D\" Seek=\"%D\">\n",ng*3,bytes[0]);
146: PetscViewerASCIIPushTab(viewer);
147: PetscViewerASCIIPrintf(viewer,"%s\n",datafile);
148: PetscViewerASCIIPopTab(viewer);
149: PetscViewerASCIIPrintf(viewer,"</DataItem>\n");
150: PetscViewerASCIIPopTab(viewer);
151: PetscViewerASCIIPrintf(viewer,"</Topology>\n");
152: PetscViewerASCIIPopTab(viewer);
153:
154: /* write topology data */
155: for (k=0; k<ng; k++) {
156: PetscInt pvertex[3];
157:
158: pvertex[0] = 1;
159: pvertex[1] = 1;
160: pvertex[2] = k;
161: PetscViewerBinaryWrite(fviewer,pvertex,3,PETSC_INT,PETSC_FALSE);
162: }
163: bytes[0] += sizeof(PetscInt) * ng * 3;
164:
165: /* write geometry header */
166: PetscViewerASCIIPushTab(viewer);
167: DMGetDimension(dm,&dim);
168: switch (dim) {
169: case 1:
170: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No support for 1D");
171: break;
172: case 2:
173: PetscViewerASCIIPrintf(viewer,"<Geometry Type=\"XY\">\n");
174: break;
175: case 3:
176: PetscViewerASCIIPrintf(viewer,"<Geometry Type=\"XYZ\">\n");
177: break;
178: }
179: PetscViewerASCIIPushTab(viewer);
180: PetscViewerASCIIPrintf(viewer,"<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Float\" Precision=\"8\" Dimensions=\"%D %D\" Seek=\"%D\">\n",ng,dim,bytes[0]);
181: PetscViewerASCIIPushTab(viewer);
182: PetscViewerASCIIPrintf(viewer,"%s\n",datafile);
183: PetscViewerASCIIPopTab(viewer);
184: PetscViewerASCIIPrintf(viewer,"</DataItem>\n");
185: PetscViewerASCIIPopTab(viewer);
186: PetscViewerASCIIPrintf(viewer,"</Geometry>\n");
187: PetscViewerASCIIPopTab(viewer);
188:
189: /* write geometry data */
190: DMSwarmCreateGlobalVectorFromField(dm,DMSwarmPICField_coor,&dvec);
191: VecView(dvec,fviewer);
192: DMSwarmDestroyGlobalVectorFromField(dm,DMSwarmPICField_coor,&dvec);
193: bytes[0] += sizeof(PetscReal) * ng * dim;
194:
195: PetscViewerDestroy(&fviewer);
196:
197: return(0);
198: }
200: PetscErrorCode private_VecView_Swarm_XDMF(Vec x,PetscViewer viewer)
201: {
202: long int *bytes = NULL;
203: PetscContainer container = NULL;
204: const char *viewername;
205: char datafile[PETSC_MAX_PATH_LEN];
206: PetscViewer fviewer;
207: PetscInt N,bs;
208: const char *vecname;
209: char fieldname[PETSC_MAX_PATH_LEN];
211:
213: PetscObjectQuery((PetscObject)viewer,"XDMFViewerContext",(PetscObject*)&container);
214: if (container) {
215: PetscContainerGetPointer(container,(void**)&bytes);
216: } else SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Valid to find attached data XDMFViewerContext");
217:
218: PetscViewerFileGetName(viewer,&viewername);
219: private_CreateDataFileNameXDMF(viewername,datafile);
220:
221: /* re-open a sub-viewer for all data fields */
222: /* name is viewer.name + "_swarm_fields.pbin" */
223: PetscViewerCreate(PetscObjectComm((PetscObject)viewer),&fviewer);
224: PetscViewerSetType(fviewer,PETSCVIEWERBINARY);
225: PetscViewerBinarySetSkipHeader(fviewer,PETSC_TRUE);
226: PetscViewerBinarySetSkipInfo(fviewer,PETSC_TRUE);
227: PetscViewerFileSetMode(fviewer,FILE_MODE_APPEND);
228: PetscViewerFileSetName(fviewer,datafile);
229:
230: VecGetSize(x,&N);
231: VecGetBlockSize(x,&bs);
232: N = N/bs;
233: PetscObjectGetName((PetscObject)x,&vecname);
234: if (!vecname) {
235: PetscSNPrintf(fieldname,PETSC_MAX_PATH_LEN-1,"swarmfield_%D",((PetscObject)x)->tag);
236: } else {
237: PetscSNPrintf(fieldname,PETSC_MAX_PATH_LEN-1,"%s",vecname);
238: }
239:
240: /* write data header */
241: PetscViewerASCIIPushTab(viewer);
242: PetscViewerASCIIPrintf(viewer,"<Attribute Center=\"Node\" Name=\"%s\" Type=\"None\">\n",fieldname);
243: PetscViewerASCIIPushTab(viewer);
244: if (bs == 1) {
245: PetscViewerASCIIPrintf(viewer,"<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Float\" Precision=\"8\" Dimensions=\"%D\" Seek=\"%D\">\n",N,bytes[0]);
246: } else {
247: PetscViewerASCIIPrintf(viewer,"<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Float\" Precision=\"8\" Dimensions=\"%D %D\" Seek=\"%D\">\n",N,bs,bytes[0]);
248: }
249: PetscViewerASCIIPushTab(viewer);
250: PetscViewerASCIIPrintf(viewer,"%s\n",datafile);
251: PetscViewerASCIIPopTab(viewer);
252: PetscViewerASCIIPrintf(viewer,"</DataItem>\n");
253: PetscViewerASCIIPopTab(viewer);
254: PetscViewerASCIIPrintf(viewer,"</Attribute>\n");
255: PetscViewerASCIIPopTab(viewer);
256:
257: /* write data */
258: VecView(x,fviewer);
259: bytes[0] += sizeof(PetscReal) * N * bs;
260:
261: PetscViewerDestroy(&fviewer);
262:
263: return(0);
264: }
266: PetscErrorCode private_ISView_Swarm_XDMF(IS is,PetscViewer viewer)
267: {
268: long int *bytes = NULL;
269: PetscContainer container = NULL;
270: const char *viewername;
271: char datafile[PETSC_MAX_PATH_LEN];
272: PetscViewer fviewer;
273: PetscInt N,bs;
274: const char *vecname;
275: char fieldname[PETSC_MAX_PATH_LEN];
277:
279: PetscObjectQuery((PetscObject)viewer,"XDMFViewerContext",(PetscObject*)&container);
280: if (container) {
281: PetscContainerGetPointer(container,(void**)&bytes);
282: } else SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Valid to find attached data XDMFViewerContext");
283:
284: PetscViewerFileGetName(viewer,&viewername);
285: private_CreateDataFileNameXDMF(viewername,datafile);
286:
287: /* re-open a sub-viewer for all data fields */
288: /* name is viewer.name + "_swarm_fields.pbin" */
289: PetscViewerCreate(PetscObjectComm((PetscObject)viewer),&fviewer);
290: PetscViewerSetType(fviewer,PETSCVIEWERBINARY);
291: PetscViewerBinarySetSkipHeader(fviewer,PETSC_TRUE);
292: PetscViewerBinarySetSkipInfo(fviewer,PETSC_TRUE);
293: PetscViewerFileSetMode(fviewer,FILE_MODE_APPEND);
294: PetscViewerFileSetName(fviewer,datafile);
295:
296: ISGetSize(is,&N);
297: ISGetBlockSize(is,&bs);
298: N = N/bs;
299: PetscObjectGetName((PetscObject)is,&vecname);
300: if (!vecname) {
301: PetscSNPrintf(fieldname,PETSC_MAX_PATH_LEN-1,"swarmfield_%D",((PetscObject)is)->tag);
302: } else {
303: PetscSNPrintf(fieldname,PETSC_MAX_PATH_LEN-1,"%s",vecname);
304: }
305:
306: /* write data header */
307: PetscViewerASCIIPushTab(viewer);
308: PetscViewerASCIIPrintf(viewer,"<Attribute Center=\"Node\" Name=\"%s\" Type=\"None\">\n",fieldname);
309: PetscViewerASCIIPushTab(viewer);
310: if (bs == 1) {
311: PetscViewerASCIIPrintf(viewer,"<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Int\" Precision=\"4\" Dimensions=\"%D\" Seek=\"%D\">\n",N,bytes[0]);
312: } else {
313: PetscViewerASCIIPrintf(viewer,"<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Int\" Precision=\"4\" Dimensions=\"%D %D\" Seek=\"%D\">\n",N,bs,bytes[0]);
314: }
315: PetscViewerASCIIPushTab(viewer);
316: PetscViewerASCIIPrintf(viewer,"%s\n",datafile);
317: PetscViewerASCIIPopTab(viewer);
318: PetscViewerASCIIPrintf(viewer,"</DataItem>\n");
319: PetscViewerASCIIPopTab(viewer);
320: PetscViewerASCIIPrintf(viewer,"</Attribute>\n");
321: PetscViewerASCIIPopTab(viewer);
322:
323: /* write data */
324: ISView(is,fviewer);
325: bytes[0] += sizeof(PetscInt) * N * bs;
326:
327: PetscViewerDestroy(&fviewer);
328:
329: return(0);
330: }
332: /*@C
333: DMSwarmViewFieldsXDMF - Write a selection of DMSwarm fields to an XDMF3 file
334:
335: Collective on DM
336:
337: Input parameters:
338: + dm - the DMSwarm
339: . filename - the file name of the XDMF file (must have the extension .xmf)
340: . nfields - the number of fields to write into the XDMF file
341: - field_name_list - array of length nfields containing the textual name of fields to write
342:
343: Level: beginner
345: Notes:
346: Only fields registered with data type PETSC_DOUBLE or PETSC_INT can be written into the file
347:
348: .seealso: DMSwarmViewXDMF()
349: @*/
350: PETSC_EXTERN PetscErrorCode DMSwarmViewFieldsXDMF(DM dm,const char filename[],PetscInt nfields,const char *field_name_list[])
351: {
353: Vec dvec;
354: PetscInt f,N;
355: PetscViewer viewer;
356:
358: private_PetscViewerCreate_XDMF(PetscObjectComm((PetscObject)dm),filename,&viewer);
359: private_DMSwarmView_XDMF(dm,viewer);
360: DMSwarmGetLocalSize(dm,&N);
361: for (f=0; f<nfields; f++) {
362: void *data;
363: PetscDataType type;
365: DMSwarmGetField(dm,field_name_list[f],NULL,&type,&data);
366: DMSwarmRestoreField(dm,field_name_list[f],NULL,&type,&data);
368: if (type == PETSC_DOUBLE) {
369: DMSwarmCreateGlobalVectorFromField(dm,field_name_list[f],&dvec);
370: PetscObjectSetName((PetscObject)dvec,field_name_list[f]);
371: private_VecView_Swarm_XDMF(dvec,viewer);
372: DMSwarmDestroyGlobalVectorFromField(dm,field_name_list[f],&dvec);
373: } else if (type == PETSC_INT) {
374: IS is;
375: const PetscInt *idx;
377: DMSwarmGetField(dm,field_name_list[f],NULL,&type,&data);
378: idx = (const PetscInt*)data;
379:
380: ISCreateGeneral(PetscObjectComm((PetscObject)dm),N,idx,PETSC_USE_POINTER,&is);
381: PetscObjectSetName((PetscObject)is,field_name_list[f]);
382: private_ISView_Swarm_XDMF(is,viewer);
383: ISDestroy(&is);
384: DMSwarmRestoreField(dm,field_name_list[f],NULL,&type,&data);
385: } else SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Can only write PETSC_INT and PETSC_DOUBLE");
387: }
388: private_PetscViewerDestroy_XDMF(&viewer);
389: return(0);
390: }
392: /*@C
393: DMSwarmViewXDMF - Write DMSwarm fields to an XDMF3 file
394:
395: Collective on DM
396:
397: Input parameters:
398: + dm - the DMSwarm
399: - filename - the file name of the XDMF file (must have the extension .xmf)
400:
401: Level: beginner
403: Notes:
404: Only fields user registered with data type PETSC_DOUBLE or PETSC_INT will be written into the file
406: .seealso: DMSwarmViewFieldsXDMF()
407: @*/
408: PETSC_EXTERN PetscErrorCode DMSwarmViewXDMF(DM dm,const char filename[])
409: {
410: DM_Swarm *swarm = (DM_Swarm*)dm->data;
412: Vec dvec;
413: PetscInt f;
414: PetscViewer viewer;
415:
417: private_PetscViewerCreate_XDMF(PetscObjectComm((PetscObject)dm),filename,&viewer);
418: private_DMSwarmView_XDMF(dm,viewer);
419: for (f=4; f<swarm->db->nfields; f++) { /* only examine user defined fields - the first 4 are internally created by DMSwarmPIC */
420: DMSwarmDataField field;
421:
422: /* query field type - accept all those of type PETSC_DOUBLE */
423: field = swarm->db->field[f];
424: if (field->petsc_type == PETSC_DOUBLE) {
425: DMSwarmCreateGlobalVectorFromField(dm,field->name,&dvec);
426: PetscObjectSetName((PetscObject)dvec,field->name);
427: private_VecView_Swarm_XDMF(dvec,viewer);
428: DMSwarmDestroyGlobalVectorFromField(dm,field->name,&dvec);
429: } else if (field->petsc_type == PETSC_INT) {
430: IS is;
431: PetscInt N;
432: const PetscInt *idx;
433: void *data;
434:
435: DMSwarmGetLocalSize(dm,&N);
436: DMSwarmGetField(dm,field->name,NULL,NULL,&data);
437: idx = (const PetscInt*)data;
438:
439: ISCreateGeneral(PetscObjectComm((PetscObject)dm),N,idx,PETSC_USE_POINTER,&is);
440: PetscObjectSetName((PetscObject)is,field->name);
441: private_ISView_Swarm_XDMF(is,viewer);
442: ISDestroy(&is);
443: DMSwarmRestoreField(dm,field->name,NULL,NULL,&data);
444: }
445: }
446: private_PetscViewerDestroy_XDMF(&viewer);
447: return(0);
448: }