Actual source code: vtkv.c
petsc-3.13.6 2020-09-29
1: #include <../src/sys/classes/viewer/impls/vtk/vtkvimpl.h>
3: /*MC
4: PetscViewerVTKWriteFunction - functional form used to provide writer to the PetscViewerVTK
6: Synopsis:
7: #include <petscviewer.h>
8: PetscViewerVTKWriteFunction(PetscObject object,PetscViewer viewer)
10: Input Parameters:
11: + object - the PETSc object to be written
12: - viewer - viewer it is to be written to
14: Level: developer
16: .seealso: PetscViewerVTKAddField()
17: M*/
19: /*@C
20: PetscViewerVTKAddField - Add a field to the viewer
22: Collective
24: Input Arguments:
25: + viewer - VTK viewer
26: . dm - DM on which Vec lives
27: . PetscViewerVTKWriteFunction - function to write this Vec
28: . fieldnum - which field of the DM to write (PETSC_DEFAULT if the whle vector should be written)
29: . fieldtype - Either PETSC_VTK_POINT_FIELD or PETSC_VTK_CELL_FIELD
30: . checkdm - whether to check for identical dm arguments as fields are added
31: - vec - Vec from which to write
33: Note:
34: This routine keeps exclusive ownership of the Vec. The caller should not use or destroy the Vec after adding it.
36: Level: developer
38: .seealso: PetscViewerVTKOpen(), DMDAVTKWriteAll(), PetscViewerVTKWriteFunction, PetscViewerVTKGetDM()
39: @*/
40: PetscErrorCode PetscViewerVTKAddField(PetscViewer viewer,PetscObject dm,PetscErrorCode (*PetscViewerVTKWriteFunction)(PetscObject,PetscViewer),PetscInt fieldnum,PetscViewerVTKFieldType fieldtype,PetscBool checkdm,PetscObject vec)
41: {
48: PetscUseMethod(viewer,"PetscViewerVTKAddField_C",(PetscViewer,PetscObject,PetscErrorCode (*)(PetscObject,PetscViewer),PetscInt,PetscViewerVTKFieldType,PetscBool,PetscObject),(viewer,dm,PetscViewerVTKWriteFunction,fieldnum,fieldtype,checkdm,vec));
49: return(0);
50: }
52: /*@C
53: PetscViewerVTKGetDM - get the DM associated with the viewer
55: Collective
57: Input Arguments:
58: + viewer - VTK viewer
59: - dm - DM associated with the viewer (as PetscObject)
61: Level: developer
63: .seealso: PetscViewerVTKOpen(), DMDAVTKWriteAll(), PetscViewerVTKWriteFunction, PetscViewerVTKAddField()
64: @*/
65: PetscErrorCode PetscViewerVTKGetDM(PetscViewer viewer,PetscObject *dm)
66: {
71: PetscUseMethod(viewer,"PetscViewerVTKGetDM_C",(PetscViewer,PetscObject*),(viewer,dm));
72: return(0);
73: }
75: static PetscErrorCode PetscViewerDestroy_VTK(PetscViewer viewer)
76: {
77: PetscViewer_VTK *vtk = (PetscViewer_VTK*)viewer->data;
78: PetscErrorCode ierr;
81: PetscFree(vtk->filename);
82: PetscFree(vtk);
83: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileSetName_C",NULL);
84: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileGetName_C",NULL);
85: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileSetMode_C",NULL);
86: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerFileGetMode_C",NULL);
87: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerVTKAddField_C",NULL);
88: PetscObjectComposeFunction((PetscObject)viewer,"PetscViewerVTKGetDM_C",NULL);
89: return(0);
90: }
92: static PetscErrorCode PetscViewerFlush_VTK(PetscViewer viewer)
93: {
94: PetscViewer_VTK *vtk = (PetscViewer_VTK*)viewer->data;
95: PetscErrorCode ierr;
96: PetscViewerVTKObjectLink link,next;
99: if (vtk->link && (!vtk->dm || !vtk->write)) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_ARG_WRONGSTATE,"No fields or no grid");
100: if (vtk->write) {(*vtk->write)(vtk->dm,viewer);}
101: for (link=vtk->link; link; link=next) {
102: next = link->next;
103: PetscObjectDestroy(&link->vec);
104: PetscFree(link);
105: }
106: PetscObjectDestroy(&vtk->dm);
107: vtk->write = NULL;
108: vtk->link = NULL;
109: return(0);
110: }
112: PetscErrorCode PetscViewerFileSetName_VTK(PetscViewer viewer,const char name[])
113: {
114: PetscViewer_VTK *vtk = (PetscViewer_VTK*)viewer->data;
115: PetscErrorCode ierr;
116: PetscBool isvtk,isvts,isvtu,isvtr;
117: size_t len;
120: PetscViewerFlush(viewer);
121: PetscFree(vtk->filename);
122: PetscStrlen(name,&len);
123: if (!len) {
124: isvtk = PETSC_TRUE;
125: } else {
126: PetscStrcasecmp(name+len-4,".vtk",&isvtk);
127: PetscStrcasecmp(name+len-4,".vts",&isvts);
128: PetscStrcasecmp(name+len-4,".vtu",&isvtu);
129: PetscStrcasecmp(name+len-4,".vtr",&isvtr);
130: }
131: if (isvtk) {
132: if (viewer->format == PETSC_VIEWER_DEFAULT) viewer->format = PETSC_VIEWER_ASCII_VTK;
133: if (viewer->format != PETSC_VIEWER_ASCII_VTK) SETERRQ2(PetscObjectComm((PetscObject)viewer),PETSC_ERR_ARG_INCOMP,"Cannot use file '%s' with format %s, should have '.vtk' extension",name,PetscViewerFormats[viewer->format]);
134: } else if (isvts) {
135: if (viewer->format == PETSC_VIEWER_DEFAULT) viewer->format = PETSC_VIEWER_VTK_VTS;
136: if (viewer->format != PETSC_VIEWER_VTK_VTS) SETERRQ2(PetscObjectComm((PetscObject)viewer),PETSC_ERR_ARG_INCOMP,"Cannot use file '%s' with format %s, should have '.vts' extension",name,PetscViewerFormats[viewer->format]);
137: } else if (isvtu) {
138: if (viewer->format == PETSC_VIEWER_DEFAULT) viewer->format = PETSC_VIEWER_VTK_VTU;
139: if (viewer->format != PETSC_VIEWER_VTK_VTU) SETERRQ2(PetscObjectComm((PetscObject)viewer),PETSC_ERR_ARG_INCOMP,"Cannot use file '%s' with format %s, should have '.vtu' extension",name,PetscViewerFormats[viewer->format]);
140: } else if (isvtr) {
141: if (viewer->format == PETSC_VIEWER_DEFAULT) viewer->format = PETSC_VIEWER_VTK_VTR;
142: if (viewer->format != PETSC_VIEWER_VTK_VTR) SETERRQ2(PetscObjectComm((PetscObject)viewer),PETSC_ERR_ARG_INCOMP,"Cannot use file '%s' with format %s, should have '.vtr' extension",name,PetscViewerFormats[viewer->format]);
143: } else SETERRQ1(PetscObjectComm((PetscObject)viewer),PETSC_ERR_ARG_UNKNOWN_TYPE,"File '%s' has unrecognized extension",name);
144: PetscStrallocpy(len ? name : "stdout",&vtk->filename);
145: return(0);
146: }
148: PetscErrorCode PetscViewerFileGetName_VTK(PetscViewer viewer,const char **name)
149: {
150: PetscViewer_VTK *vtk = (PetscViewer_VTK*)viewer->data;
152: *name = vtk->filename;
153: return(0);
154: }
156: PetscErrorCode PetscViewerFileSetMode_VTK(PetscViewer viewer,PetscFileMode type)
157: {
158: PetscViewer_VTK *vtk = (PetscViewer_VTK*)viewer->data;
161: vtk->btype = type;
162: return(0);
163: }
165: PetscErrorCode PetscViewerFileGetMode_VTK(PetscViewer viewer,PetscFileMode *type)
166: {
167: PetscViewer_VTK *vtk = (PetscViewer_VTK*)viewer->data;
170: *type = vtk->btype;
171: return(0);
172: }
174: PetscErrorCode PetscViewerVTKAddField_VTK(PetscViewer viewer,PetscObject dm,PetscErrorCode (*PetscViewerVTKWriteFunction)(PetscObject,PetscViewer),PetscInt fieldnum,PetscViewerVTKFieldType fieldtype,PetscBool checkdm,PetscObject vec)
175: {
176: PetscViewer_VTK *vtk = (PetscViewer_VTK*)viewer->data;
177: PetscViewerVTKObjectLink link, tail = vtk->link;
178: PetscErrorCode ierr;
181: if (vtk->dm) {
182: if (checkdm && dm != vtk->dm) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_ARG_INCOMP,"Refusing to write a field from more than one grid to the same VTK file. Set checkdm = PETSC_FALSE to skip this check.");
183: } else {
184: PetscObjectReference(dm);
185: vtk->dm = dm;
186: }
187: vtk->write = PetscViewerVTKWriteFunction;
188: PetscNew(&link);
189: link->ft = fieldtype;
190: link->vec = vec;
191: link->field = fieldnum;
192: link->next = NULL;
193: /* Append to list */
194: if (tail) {
195: while (tail->next) tail = tail->next;
196: tail->next = link;
197: } else vtk->link = link;
198: return(0);
199: }
201: PetscErrorCode PetscViewerVTKGetDM_VTK(PetscViewer viewer,PetscObject *dm)
202: {
203: PetscViewer_VTK *vtk = (PetscViewer_VTK*)viewer->data;
206: *dm = vtk->dm;
207: return(0);
208: }
210: /*MC
211: PETSCVIEWERVTK - A viewer that writes to a VTK file
214: .seealso: PetscViewerVTKOpen(), PetscViewerHDF5Open(), PetscViewerStringSPrintf(), PetscViewerSocketOpen(), PetscViewerDrawOpen(), PETSCVIEWERSOCKET,
215: PetscViewerCreate(), PetscViewerASCIIOpen(), PetscViewerBinaryOpen(), PETSCVIEWERBINARY, PETSCVIEWERDRAW, PETSCVIEWERSTRING,
216: PetscViewerMatlabOpen(), VecView(), DMView(), PetscViewerMatlabPutArray(), PETSCVIEWERASCII, PETSCVIEWERMATLAB,
217: PetscViewerFileSetName(), PetscViewerFileSetMode(), PetscViewerFormat, PetscViewerType, PetscViewerSetType()
219: Level: beginner
220: M*/
222: PETSC_EXTERN PetscErrorCode PetscViewerCreate_VTK(PetscViewer v)
223: {
224: PetscViewer_VTK *vtk;
225: PetscErrorCode ierr;
228: PetscNewLog(v,&vtk);
230: v->data = (void*)vtk;
231: v->ops->destroy = PetscViewerDestroy_VTK;
232: v->ops->flush = PetscViewerFlush_VTK;
233: vtk->btype = (PetscFileMode) -1;
234: vtk->filename = NULL;
236: PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileSetName_C",PetscViewerFileSetName_VTK);
237: PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileGetName_C",PetscViewerFileGetName_VTK);
238: PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileSetMode_C",PetscViewerFileSetMode_VTK);
239: PetscObjectComposeFunction((PetscObject)v,"PetscViewerFileGetMode_C",PetscViewerFileGetMode_VTK);
240: PetscObjectComposeFunction((PetscObject)v,"PetscViewerVTKAddField_C",PetscViewerVTKAddField_VTK);
241: PetscObjectComposeFunction((PetscObject)v,"PetscViewerVTKGetDM_C",PetscViewerVTKGetDM_VTK);
242: return(0);
243: }
245: /*@C
246: PetscViewerVTKOpen - Opens a file for VTK output.
248: Collective
250: Input Parameters:
251: + comm - MPI communicator
252: . name - name of file
253: - type - type of file
254: $ FILE_MODE_WRITE - create new file for binary output
255: $ FILE_MODE_READ - open existing file for binary input (not currently supported)
256: $ FILE_MODE_APPEND - open existing file for binary output (not currently supported)
258: Output Parameter:
259: . vtk - PetscViewer for VTK input/output to use with the specified file
261: Level: beginner
263: Note:
264: This PetscViewer should be destroyed with PetscViewerDestroy().
267: .seealso: PetscViewerASCIIOpen(), PetscViewerPushFormat(), PetscViewerDestroy(),
268: VecView(), MatView(), VecLoad(), MatLoad(),
269: PetscFileMode, PetscViewer
270: @*/
271: PetscErrorCode PetscViewerVTKOpen(MPI_Comm comm,const char name[],PetscFileMode type,PetscViewer *vtk)
272: {
276: PetscViewerCreate(comm,vtk);
277: PetscViewerSetType(*vtk,PETSCVIEWERVTK);
278: PetscViewerFileSetMode(*vtk,type);
279: PetscViewerFileSetName(*vtk,name);
280: return(0);
281: }
283: /*@C
284: PetscViewerVTKFWrite - write binary data preceded by 32-bit int length (in bytes), does not do byte swapping.
286: Logically collective on PetscViewer
288: Input Parameters:
289: + viewer - logically collective viewer, data written from rank 0
290: . fp - file pointer valid on rank 0
291: . data - data pointer valid on rank 0
292: . n - number of data items
293: - dtype - data type
295: Level: developer
297: Notes:
298: If PetscScalar is __float128 then the binary files are written in double precision
301: .seealso: DMDAVTKWriteAll(), DMComplexVTKWriteAll(), PetscViewerPushFormat(), PetscViewerVTKOpen(), PetscBinaryWrite()
302: @*/
303: PetscErrorCode PetscViewerVTKFWrite(PetscViewer viewer,FILE *fp,const void *data,PetscInt n,MPI_Datatype dtype)
304: {
306: PetscMPIInt rank;
307: MPI_Datatype vdtype=dtype;
308: #if defined(PETSC_USE_REAL___FLOAT128)
309: double *tmp;
310: PetscInt i;
311: PetscReal *ttmp = (PetscReal*)data;
312: #endif
315: if (n < 0) SETERRQ1(PetscObjectComm((PetscObject)viewer),PETSC_ERR_ARG_OUTOFRANGE,"Trying to write a negative amount of data %D",n);
316: if (!n) return(0);
317: MPI_Comm_rank(PetscObjectComm((PetscObject)viewer),&rank);
318: if (!rank) {
319: size_t count;
320: PetscMPIInt dsize;
321: PetscVTKInt bytes;
323: #if defined(PETSC_USE_REAL___FLOAT128)
324: if (dtype == MPIU___FLOAT128) {
325: PetscMalloc1(n,&tmp);
326: for (i=0; i<n; i++) tmp[i] = ttmp[i];
327: data = (void*) tmp;
328: vdtype = MPI_DOUBLE;
329: }
330: #endif
331: MPI_Type_size(vdtype,&dsize);
332: bytes = PetscVTKIntCast(dsize*n);
334: count = fwrite(&bytes,sizeof(int),1,fp);
335: if (count != 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FILE_WRITE,"Error writing byte count");
336: count = fwrite(data,dsize,(size_t)n,fp);
337: if ((PetscInt)count != n) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_FILE_WRITE,"Wrote %D/%D array members of size %d",(PetscInt)count,n,dsize);
338: #if defined(PETSC_USE_REAL___FLOAT128)
339: if (dtype == MPIU___FLOAT128) {
340: PetscFree(tmp);
341: }
342: #endif
343: }
344: return(0);
345: }