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