Actual source code: plexvtu.c
petsc-3.10.5 2019-03-28
1: #include <petsc/private/dmpleximpl.h>
2: #include <../src/sys/classes/viewer/impls/vtk/vtkvimpl.h>
4: typedef struct {
5: PetscInt nvertices;
6: PetscInt ncells;
7: PetscInt nconn; /* number of entries in cell->vertex connectivity array */
8: } PieceInfo;
10: #if defined(PETSC_USE_REAL_SINGLE)
11: static const char precision[] = "Float32";
12: #elif defined(PETSC_USE_REAL_DOUBLE)
13: static const char precision[] = "Float64";
14: #else
15: static const char precision[] = "UnknownPrecision";
16: #endif
18: static PetscErrorCode TransferWrite(PetscViewer viewer,FILE *fp,PetscMPIInt srank,PetscMPIInt root,const void *send,void *recv,PetscMPIInt count,MPI_Datatype mpidatatype,PetscMPIInt tag)
19: {
20: PetscMPIInt rank;
22: MPI_Comm comm;
25: PetscObjectGetComm((PetscObject)viewer,&comm);
26: MPI_Comm_rank(comm,&rank);
28: if (rank == srank && rank != root) {
29: MPI_Send((void*)send,count,mpidatatype,root,tag,comm);
30: } else if (rank == root) {
31: const void *buffer;
32: if (root == srank) { /* self */
33: buffer = send;
34: } else {
35: MPI_Status status;
36: PetscMPIInt nrecv;
37: MPI_Recv(recv,count,mpidatatype,srank,tag,comm,&status);
38: MPI_Get_count(&status,mpidatatype,&nrecv);
39: if (count != nrecv) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Array size mismatch");
40: buffer = recv;
41: }
42: PetscViewerVTKFWrite(viewer,fp,buffer,count,mpidatatype);
43: }
44: return(0);
45: }
47: static PetscErrorCode DMPlexGetVTKConnectivity(DM dm, PetscBool localized, PieceInfo *piece,PetscVTKInt **oconn,PetscVTKInt **ooffsets,PetscVTKType **otypes)
48: {
50: PetscSection coordSection;
51: PetscVTKInt *conn,*offsets;
52: PetscVTKType *types;
53: PetscInt dim,vStart,vEnd,cStart,cEnd,pStart,pEnd,cellHeight,numLabelCells,hasLabel,c,v,countcell,countconn;
56: PetscMalloc3(piece->nconn,&conn,piece->ncells,&offsets,piece->ncells,&types);
57: DMGetCoordinateSection(dm, &coordSection);
58: DMGetDimension(dm,&dim);
59: DMPlexGetChart(dm,&pStart,&pEnd);
60: DMPlexGetVTKCellHeight(dm, &cellHeight);
61: DMPlexGetHeightStratum(dm, cellHeight, &cStart, &cEnd);
62: DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd);
63: DMGetStratumSize(dm, "vtk", 1, &numLabelCells);
64: hasLabel = numLabelCells > 0 ? PETSC_TRUE : PETSC_FALSE;
66: countcell = 0;
67: countconn = 0;
68: for (c = cStart; c < cEnd; ++c) {
69: PetscInt nverts,dof = 0,celltype,startoffset,nC=0;
71: if (hasLabel) {
72: PetscInt value;
74: DMGetLabelValue(dm, "vtk", c, &value);
75: if (value != 1) continue;
76: }
77: startoffset = countconn;
78: if (localized) {
79: PetscSectionGetDof(coordSection, c, &dof);
80: }
81: if (!dof) {
82: PetscInt *closure = NULL;
83: PetscInt closureSize;
85: DMPlexGetTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);
86: for (v = 0; v < closureSize*2; v += 2) {
87: if ((closure[v] >= vStart) && (closure[v] < vEnd)) {
88: if (!localized) conn[countconn++] = closure[v] - vStart;
89: else conn[countconn++] = startoffset + nC;
90: ++nC;
91: }
92: }
93: DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);
94: } else {
95: for (nC = 0; nC < dof/dim; nC++) conn[countconn++] = startoffset + nC;
96: }
97: DMPlexInvertCell(dim, nC, &conn[countconn-nC]);
99: offsets[countcell] = countconn;
101: nverts = countconn - startoffset;
102: DMPlexVTKGetCellType_Internal(dm,dim,nverts,&celltype);
104: types[countcell] = celltype;
105: countcell++;
106: }
107: if (countcell != piece->ncells) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Inconsistent cell count");
108: if (countconn != piece->nconn) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Inconsistent connectivity count");
109: *oconn = conn;
110: *ooffsets = offsets;
111: *otypes = types;
112: return(0);
113: }
115: /*
116: Write all fields that have been provided to the viewer
117: Multi-block XML format with binary appended data.
118: */
119: PetscErrorCode DMPlexVTKWriteAll_VTU(DM dm,PetscViewer viewer)
120: {
121: MPI_Comm comm;
122: PetscSection coordSection;
123: PetscViewer_VTK *vtk = (PetscViewer_VTK*)viewer->data;
124: PetscViewerVTKObjectLink link;
125: FILE *fp;
126: PetscMPIInt rank,size,tag;
127: PetscErrorCode ierr;
128: PetscInt dimEmbed,cellHeight,cStart,cEnd,vStart,vEnd,numLabelCells,hasLabel,c,v,r,i;
129: PetscBool localized;
130: PieceInfo piece,*gpiece = NULL;
131: void *buffer = NULL;
134: PetscObjectGetComm((PetscObject)dm,&comm);
135: #if defined(PETSC_USE_COMPLEX)
136: SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_SUP,"Complex values not supported");
137: #endif
138: MPI_Comm_size(comm,&size);
139: MPI_Comm_rank(comm,&rank);
140: tag = ((PetscObject)viewer)->tag;
142: PetscFOpen(comm,vtk->filename,"wb",&fp);
143: PetscFPrintf(comm,fp,"<?xml version=\"1.0\"?>\n");
144: #if defined(PETSC_WORDS_BIGENDIAN)
145: PetscFPrintf(comm,fp,"<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"BigEndian\">\n");
146: #else
147: PetscFPrintf(comm,fp,"<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">\n");
148: #endif
149: PetscFPrintf(comm,fp," <UnstructuredGrid>\n");
151: DMGetCoordinateDim(dm, &dimEmbed);
152: DMPlexGetVTKCellHeight(dm, &cellHeight);
153: DMPlexGetHeightStratum(dm, cellHeight, &cStart, &cEnd);
154: DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd);
155: DMGetStratumSize(dm, "vtk", 1, &numLabelCells);
156: DMGetCoordinatesLocalized(dm, &localized);
157: DMGetCoordinateSection(dm, &coordSection);
159: hasLabel = numLabelCells > 0 ? PETSC_TRUE : PETSC_FALSE;
160: piece.nvertices = 0;
161: piece.ncells = 0;
162: piece.nconn = 0;
163: if (!localized) piece.nvertices = vEnd - vStart;
164: for (c = cStart; c < cEnd; ++c) {
165: PetscInt dof = 0;
166: if (hasLabel) {
167: PetscInt value;
169: DMGetLabelValue(dm, "vtk", c, &value);
170: if (value != 1) continue;
171: }
172: if (localized) {
173: PetscSectionGetDof(coordSection, c, &dof);
174: }
175: if (!dof) {
176: PetscInt *closure = NULL;
177: PetscInt closureSize;
179: DMPlexGetTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);
180: for (v = 0; v < closureSize*2; v += 2) {
181: if ((closure[v] >= vStart) && (closure[v] < vEnd)) {
182: piece.nconn++;
183: if (localized) piece.nvertices++;
184: }
185: }
186: DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);
187: } else {
188: piece.nvertices += dof/dimEmbed;
189: piece.nconn += dof/dimEmbed;
190: }
191: piece.ncells++;
192: }
193: if (!rank) {PetscMalloc1(size,&gpiece);}
194: MPI_Gather((PetscInt*)&piece,sizeof(piece)/sizeof(PetscInt),MPIU_INT,(PetscInt*)gpiece,sizeof(piece)/sizeof(PetscInt),MPIU_INT,0,comm);
196: /*
197: * Write file header
198: */
199: if (!rank) {
200: PetscInt boffset = 0;
202: for (r=0; r<size; r++) {
203: PetscFPrintf(PETSC_COMM_SELF,fp," <Piece NumberOfPoints=\"%D\" NumberOfCells=\"%D\">\n",gpiece[r].nvertices,gpiece[r].ncells);
204: /* Coordinate positions */
205: PetscFPrintf(PETSC_COMM_SELF,fp," <Points>\n");
206: PetscFPrintf(PETSC_COMM_SELF,fp," <DataArray type=\"%s\" Name=\"Position\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%D\" />\n",precision,boffset);
207: boffset += gpiece[r].nvertices*3*sizeof(PetscScalar) + sizeof(int);
208: PetscFPrintf(PETSC_COMM_SELF,fp," </Points>\n");
209: /* Cell connectivity */
210: PetscFPrintf(PETSC_COMM_SELF,fp," <Cells>\n");
211: PetscFPrintf(PETSC_COMM_SELF,fp," <DataArray type=\"Int32\" Name=\"connectivity\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",boffset);
212: boffset += gpiece[r].nconn*sizeof(PetscInt) + sizeof(int);
213: PetscFPrintf(PETSC_COMM_SELF,fp," <DataArray type=\"Int32\" Name=\"offsets\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",boffset);
214: boffset += gpiece[r].ncells*sizeof(PetscInt) + sizeof(int);
215: PetscFPrintf(PETSC_COMM_SELF,fp," <DataArray type=\"UInt8\" Name=\"types\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",boffset);
216: boffset += gpiece[r].ncells*sizeof(unsigned char) + sizeof(int);
217: PetscFPrintf(PETSC_COMM_SELF,fp," </Cells>\n");
219: /*
220: * Cell Data headers
221: */
222: PetscFPrintf(PETSC_COMM_SELF,fp," <CellData>\n");
223: PetscFPrintf(PETSC_COMM_SELF,fp," <DataArray type=\"Int32\" Name=\"Rank\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",boffset);
224: boffset += gpiece[r].ncells*sizeof(int) + sizeof(int);
225: /* all the vectors */
226: for (link=vtk->link; link; link=link->next) {
227: Vec X = (Vec)link->vec;
228: PetscInt bs,nfields,field;
229: const char *vecname = "";
230: if ((link->ft != PETSC_VTK_CELL_FIELD) && (link->ft != PETSC_VTK_CELL_VECTOR_FIELD)) continue;
231: if (((PetscObject)X)->name || link != vtk->link) { /* If the object is already named, use it. If it is past the first link, name it to disambiguate. */
232: PetscObjectGetName((PetscObject)X,&vecname);
233: }
234: PetscSectionGetDof(dm->defaultSection,cStart,&bs);
235: PetscSectionGetNumFields(dm->defaultSection,&nfields);
236: for (field=0,i=0; field<(nfields?nfields:1); field++) {
237: PetscInt fbs,j;
238: PetscFV fv = NULL;
239: PetscObject f;
240: PetscClassId fClass;
241: const char *fieldname = NULL;
242: char buf[256];
243: if (nfields) { /* We have user-defined fields/components */
244: PetscSectionGetFieldDof(dm->defaultSection,cStart,field,&fbs);
245: PetscSectionGetFieldName(dm->defaultSection,field,&fieldname);
246: } else fbs = bs; /* Say we have one field with 'bs' components */
247: DMGetField(dm,field,&f);
248: PetscObjectGetClassId(f,&fClass);
249: if (fClass == PETSCFV_CLASSID) {
250: fv = (PetscFV) f;
251: }
252: if (!fieldname) {
253: PetscSNPrintf(buf,sizeof(buf),"CellField%D",field);
254: fieldname = buf;
255: }
256: for (j=0; j<fbs; j++) {
257: const char *compName = NULL;
258: if (fv) {
259: PetscFVGetComponentName(fv,j,&compName);
260: }
261: if (compName) {
262: PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s%s.%s\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",precision,vecname,fieldname,compName,boffset);
263: }
264: else {
265: PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s%s.%D\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",precision,vecname,fieldname,j,boffset);
266: }
267: boffset += gpiece[r].ncells*sizeof(PetscScalar) + sizeof(int);
268: i++;
269: }
270: }
271: if (i != bs) SETERRQ2(comm,PETSC_ERR_PLIB,"Total number of field components %D != block size %D",i,bs);
272: }
273: PetscFPrintf(PETSC_COMM_SELF,fp," </CellData>\n");
275: /*
276: * Point Data headers
277: */
278: PetscFPrintf(PETSC_COMM_SELF,fp," <PointData>\n");
279: for (link=vtk->link; link; link=link->next) {
280: Vec X = (Vec)link->vec;
281: PetscInt bs,nfields,field;
282: const char *vecname = "";
283: if ((link->ft != PETSC_VTK_POINT_FIELD) && (link->ft != PETSC_VTK_POINT_VECTOR_FIELD)) continue;
284: if (((PetscObject)X)->name || link != vtk->link) { /* If the object is already named, use it. If it is past the first link, name it to disambiguate. */
285: PetscObjectGetName((PetscObject)X,&vecname);
286: }
287: PetscSectionGetDof(dm->defaultSection,vStart,&bs);
288: PetscSectionGetNumFields(dm->defaultSection,&nfields);
289: for (field=0,i=0; field<(nfields?nfields:1); field++) {
290: PetscInt fbs,j;
291: const char *fieldname = NULL;
292: char buf[256];
293: if (nfields) { /* We have user-defined fields/components */
294: PetscSectionGetFieldDof(dm->defaultSection,vStart,field,&fbs);
295: PetscSectionGetFieldName(dm->defaultSection,field,&fieldname);
296: } else fbs = bs; /* Say we have one field with 'bs' components */
297: if (!fieldname) {
298: PetscSNPrintf(buf,sizeof(buf),"PointField%D",field);
299: fieldname = buf;
300: }
301: for (j=0; j<fbs; j++) {
302: PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s%s.%D\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",precision,vecname,fieldname,j,boffset);
303: boffset += gpiece[r].nvertices*sizeof(PetscScalar) + sizeof(int);
304: }
305: }
306: }
307: PetscFPrintf(PETSC_COMM_SELF,fp," </PointData>\n");
309: PetscFPrintf(PETSC_COMM_SELF,fp," </Piece>\n");
310: }
311: }
313: PetscFPrintf(comm,fp," </UnstructuredGrid>\n");
314: PetscFPrintf(comm,fp," <AppendedData encoding=\"raw\">\n");
315: PetscFPrintf(comm,fp,"_");
317: if (!rank) {
318: PetscInt maxsize = 0;
319: for (r=0; r<size; r++) {
320: maxsize = PetscMax(maxsize, (PetscInt) (gpiece[r].nvertices*3*sizeof(PetscScalar)));
321: maxsize = PetscMax(maxsize, (PetscInt) (gpiece[r].ncells*sizeof(PetscScalar)));
322: maxsize = PetscMax(maxsize, (PetscInt) (gpiece[r].nconn*sizeof(PetscVTKInt)));
323: }
324: PetscMalloc(maxsize,&buffer);
325: }
326: for (r=0; r<size; r++) {
327: if (r == rank) {
328: PetscInt nsend;
329: { /* Position */
330: const PetscScalar *x;
331: PetscScalar *y = NULL;
332: Vec coords;
334: DMGetCoordinatesLocal(dm,&coords);
335: VecGetArrayRead(coords,&x);
336: if (dimEmbed != 3 || localized) {
337: PetscMalloc1(piece.nvertices*3,&y);
338: if (localized) {
339: PetscInt cnt;
340: for (c=cStart,cnt=0; c<cEnd; c++) {
341: PetscInt off, dof;
343: PetscSectionGetDof(coordSection, c, &dof);
344: if (!dof) {
345: PetscInt *closure = NULL;
346: PetscInt closureSize;
348: DMPlexGetTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);
349: for (v = 0; v < closureSize*2; v += 2) {
350: if ((closure[v] >= vStart) && (closure[v] < vEnd)) {
351: PetscSectionGetOffset(coordSection, closure[v], &off);
352: if (dimEmbed != 3) {
353: y[cnt*3+0] = x[off+0];
354: y[cnt*3+1] = (dimEmbed > 1) ? x[off+1] : 0.0;
355: y[cnt*3+2] = 0.0;
356: } else {
357: y[cnt*3+0] = x[off+0];
358: y[cnt*3+1] = x[off+1];
359: y[cnt*3+2] = x[off+2];
360: }
361: cnt++;
362: }
363: }
364: DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);
365: } else {
366: PetscSectionGetOffset(coordSection, c, &off);
367: if (dimEmbed != 3) {
368: for (i=0; i<dof/dimEmbed; i++) {
369: y[cnt*3+0] = x[off + i*dimEmbed + 0];
370: y[cnt*3+1] = (dimEmbed > 1) ? x[off + i*dimEmbed + 1] : 0.0;
371: y[cnt*3+2] = 0.0;
372: cnt++;
373: }
374: } else {
375: for (i=0; i<dof; i ++) {
376: y[cnt*3+i] = x[off + i];
377: }
378: cnt += dof/dimEmbed;
379: }
380: }
381: }
382: if (cnt != piece.nvertices) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Count does not match");
383: } else {
384: for (i=0; i<piece.nvertices; i++) {
385: y[i*3+0] = x[i*dimEmbed+0];
386: y[i*3+1] = (dimEmbed > 1) ? x[i*dimEmbed+1] : 0;
387: y[i*3+2] = 0.0;
388: }
389: }
390: }
391: nsend = piece.nvertices*3;
392: TransferWrite(viewer,fp,r,0,y ? y : x,buffer,nsend,MPIU_SCALAR,tag);
393: PetscFree(y);
394: VecRestoreArrayRead(coords,&x);
395: }
396: { /* Connectivity, offsets, types */
397: PetscVTKInt *connectivity = NULL, *offsets = NULL;
398: PetscVTKType *types = NULL;
399: DMPlexGetVTKConnectivity(dm,localized,&piece,&connectivity,&offsets,&types);
400: TransferWrite(viewer,fp,r,0,connectivity,buffer,piece.nconn,MPI_INT,tag);
401: TransferWrite(viewer,fp,r,0,offsets,buffer,piece.ncells,MPI_INT,tag);
402: TransferWrite(viewer,fp,r,0,types,buffer,piece.ncells,MPI_CHAR,tag);
403: PetscFree3(connectivity,offsets,types);
404: }
405: { /* Owners (cell data) */
406: PetscVTKInt *owners;
407: PetscMalloc1(piece.ncells,&owners);
408: for (i=0; i<piece.ncells; i++) owners[i] = rank;
409: TransferWrite(viewer,fp,r,0,owners,buffer,piece.ncells,MPI_INT,tag);
410: PetscFree(owners);
411: }
412: /* Cell data */
413: for (link=vtk->link; link; link=link->next) {
414: Vec X = (Vec)link->vec;
415: const PetscScalar *x;
416: PetscScalar *y;
417: PetscInt bs;
418: if ((link->ft != PETSC_VTK_CELL_FIELD) && (link->ft != PETSC_VTK_CELL_VECTOR_FIELD)) continue;
419: PetscSectionGetDof(dm->defaultSection,cStart,&bs);
420: VecGetArrayRead(X,&x);
421: PetscMalloc1(piece.ncells,&y);
422: for (i=0; i<bs; i++) {
423: PetscInt cnt;
424: for (c=cStart,cnt=0; c<cEnd; c++) {
425: PetscScalar *xpoint;
426: if (hasLabel) { /* Ignore some cells */
427: PetscInt value;
428: DMGetLabelValue(dm, "vtk", c, &value);
429: if (value != 1) continue;
430: }
431: DMPlexPointLocalRead(dm,c,x,&xpoint);
432: y[cnt++] = xpoint[i];
433: }
434: if (cnt != piece.ncells) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Count does not match");
435: TransferWrite(viewer,fp,r,0,y,buffer,piece.ncells,MPIU_SCALAR,tag);
436: }
437: PetscFree(y);
438: VecRestoreArrayRead(X,&x);
439: }
441: for (link=vtk->link; link; link=link->next) {
442: Vec X = (Vec)link->vec;
443: const PetscScalar *x;
444: PetscScalar *y;
445: PetscInt bs;
446: if ((link->ft != PETSC_VTK_POINT_FIELD) && (link->ft != PETSC_VTK_POINT_VECTOR_FIELD)) continue;
447: PetscSectionGetDof(dm->defaultSection,vStart,&bs);
448: VecGetArrayRead(X,&x);
449: PetscMalloc1(piece.nvertices,&y);
450: for (i=0; i<bs; i++) {
451: PetscInt cnt;
452: if (!localized) {
453: for (v=vStart,cnt=0; v<vEnd; v++) {
454: PetscScalar *xpoint;
455: DMPlexPointLocalRead(dm,v,x,&xpoint);
456: y[cnt++] = xpoint[i];
457: }
458: } else {
459: for (c=cStart,cnt=0; c<cEnd; c++) {
460: PetscInt *closure = NULL;
461: PetscInt closureSize, off;
463: DMPlexGetTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);
464: for (v = 0, off = 0; v < closureSize*2; v += 2) {
465: if ((closure[v] >= vStart) && (closure[v] < vEnd)) {
466: PetscScalar *xpoint;
468: DMPlexPointLocalRead(dm,closure[v],x,&xpoint);
469: y[cnt + off++] = xpoint[i];
470: }
471: }
472: cnt += off;
473: DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);
474: }
475: }
476: if (cnt != piece.nvertices) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Count does not match");
477: TransferWrite(viewer,fp,r,0,y,buffer,piece.nvertices,MPIU_SCALAR,tag);
478: }
479: PetscFree(y);
480: VecRestoreArrayRead(X,&x);
481: }
482: } else if (!rank) {
483: TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].nvertices*3,MPIU_SCALAR,tag); /* positions */
484: TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].nconn,MPI_INT,tag); /* connectivity */
485: TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].ncells,MPI_INT,tag); /* offsets */
486: TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].ncells,MPI_CHAR,tag); /* types */
487: TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].ncells,MPI_INT,tag); /* owner rank (cells) */
488: /* all cell data */
489: for (link=vtk->link; link; link=link->next) {
490: PetscInt bs;
491: if ((link->ft != PETSC_VTK_CELL_FIELD) && (link->ft != PETSC_VTK_CELL_VECTOR_FIELD)) continue;
492: PetscSectionGetDof(dm->defaultSection,cStart,&bs);
493: for (i=0; i<bs; i++) {
494: TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].ncells,MPIU_SCALAR,tag);
495: }
496: }
497: /* all point data */
498: for (link=vtk->link; link; link=link->next) {
499: PetscInt bs;
500: if ((link->ft != PETSC_VTK_POINT_FIELD) && (link->ft != PETSC_VTK_POINT_VECTOR_FIELD)) continue;
501: PetscSectionGetDof(dm->defaultSection,vStart,&bs);
502: for (i=0; i<bs; i++) {
503: TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].nvertices,MPIU_SCALAR,tag);
504: }
505: }
506: }
507: }
508: PetscFree(gpiece);
509: PetscFree(buffer);
510: PetscFPrintf(comm,fp,"\n </AppendedData>\n");
511: PetscFPrintf(comm,fp,"</VTKFile>\n");
512: PetscFClose(comm,fp);
513: return(0);
514: }