Actual source code: plexvtu.c
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) || defined(PETSC_USE_REAL___FP16)
11: /* output in float if single or half precision in memory */
12: static const char precision[] = "Float32";
13: typedef float PetscVTUReal;
14: #define MPIU_VTUREAL MPI_FLOAT
15: #elif defined(PETSC_USE_REAL_DOUBLE) || defined(PETSC_USE_REAL___FLOAT128)
16: /* output in double if double or quad precision in memory */
17: static const char precision[] = "Float64";
18: typedef double PetscVTUReal;
19: #define MPIU_VTUREAL MPI_DOUBLE
20: #else
21: static const char precision[] = "UnknownPrecision";
22: typedef PetscReal PetscVTUReal;
23: #define MPIU_VTUREAL MPIU_REAL
24: #endif
26: static PetscErrorCode TransferWrite(MPI_Comm comm, PetscViewer viewer, FILE *fp, PetscMPIInt srank, PetscMPIInt root, const void *send, void *recv, PetscCount count, MPI_Datatype mpidatatype, PetscMPIInt tag)
27: {
28: PetscMPIInt rank;
30: PetscFunctionBegin;
31: PetscCallMPI(MPI_Comm_rank(comm, &rank));
32: if (rank == srank && rank != root) {
33: PetscCallMPI(MPIU_Send((void *)send, count, mpidatatype, root, tag, comm));
34: } else if (rank == root) {
35: const void *buffer;
36: if (root == srank) { /* self */
37: buffer = send;
38: } else {
39: MPI_Status status;
40: PetscMPIInt nrecv;
41: PetscCallMPI(MPIU_Recv(recv, count, mpidatatype, srank, tag, comm, &status));
42: PetscCallMPI(MPI_Get_count(&status, mpidatatype, &nrecv));
43: PetscCheck(count == nrecv, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Array size mismatch");
44: buffer = recv;
45: }
46: PetscCall(PetscViewerVTKFWrite(viewer, fp, buffer, count, mpidatatype));
47: }
48: PetscFunctionReturn(PETSC_SUCCESS);
49: }
51: static PetscErrorCode DMPlexGetVTKConnectivity(DM dm, PetscBool localized, PieceInfo *piece, PetscVTKInt **oconn, PetscVTKInt **ooffsets, PetscVTKType **otypes)
52: {
53: PetscSection coordSection, cellCoordSection;
54: PetscVTKInt *conn, *offsets;
55: PetscVTKType *types;
56: PetscInt dim, vStart, vEnd, cStart, cEnd, pStart, pEnd, cellHeight, numLabelCells, hasLabel, c, v, countcell, countconn;
58: PetscFunctionBegin;
59: PetscCall(PetscMalloc3(piece->nconn, &conn, piece->ncells, &offsets, piece->ncells, &types));
60: PetscCall(DMGetCoordinateSection(dm, &coordSection));
61: PetscCall(DMGetCellCoordinateSection(dm, &cellCoordSection));
62: PetscCall(DMGetDimension(dm, &dim));
63: PetscCall(DMPlexGetChart(dm, &pStart, &pEnd));
64: PetscCall(DMPlexGetVTKCellHeight(dm, &cellHeight));
65: PetscCall(DMPlexGetHeightStratum(dm, cellHeight, &cStart, &cEnd));
66: PetscCall(DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd));
67: PetscCall(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: PetscCall(DMGetLabelValue(dm, "vtk", c, &value));
79: if (value != 1) continue;
80: }
81: startoffset = countconn;
82: if (localized) PetscCall(PetscSectionGetDof(cellCoordSection, c, &dof));
83: if (!dof) {
84: PetscInt *closure = NULL;
85: PetscInt closureSize;
87: PetscCall(DMPlexGetTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure));
88: for (v = 0; v < closureSize * 2; v += 2) {
89: if ((closure[v] >= vStart) && (closure[v] < vEnd)) {
90: if (!localized) PetscCall(PetscVTKIntCast(closure[v] - vStart, &conn[countconn++]));
91: else PetscCall(PetscVTKIntCast(startoffset + nC, &conn[countconn++]));
92: ++nC;
93: }
94: }
95: PetscCall(DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure));
96: } else {
97: for (nC = 0; nC < dof / dim; nC++) PetscCall(PetscVTKIntCast(startoffset + nC, &conn[countconn++]));
98: }
100: {
101: PetscInt n = PetscMin(nC, 8), s = countconn - nC, i, cone[8];
102: for (i = 0; i < n; ++i) cone[i] = conn[s + i];
103: PetscCall(DMPlexReorderCell(dm, c, cone));
104: for (i = 0; i < n; ++i) PetscCall(PetscVTKIntCast(cone[i], &conn[s + i]));
105: }
106: PetscCall(PetscVTKIntCast(countconn, &offsets[countcell]));
108: nverts = countconn - startoffset;
109: PetscCall(DMPlexVTKGetCellType_Internal(dm, dim, nverts, &celltype));
111: types[countcell] = (PetscVTKType)celltype;
112: countcell++;
113: }
114: PetscCheck(countcell == piece->ncells, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Inconsistent cell count");
115: PetscCheck(countconn == piece->nconn, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Inconsistent connectivity count");
116: *oconn = conn;
117: *ooffsets = offsets;
118: *otypes = types;
119: PetscFunctionReturn(PETSC_SUCCESS);
120: }
122: PETSC_INTERN PetscErrorCode DMPlexGetNonEmptyComm_Private(DM dm, MPI_Comm *comm)
123: {
124: DM_Plex *mesh = (DM_Plex *)dm->data;
126: PetscFunctionBegin;
127: if (mesh->nonempty_comm == MPI_COMM_SELF) { /* Not yet setup */
128: PetscInt cStart, cEnd, cellHeight;
129: MPI_Comm dmcomm = PetscObjectComm((PetscObject)dm);
130: PetscMPIInt color, rank;
132: PetscCall(DMPlexGetVTKCellHeight(dm, &cellHeight));
133: PetscCall(DMPlexGetHeightStratum(dm, cellHeight, &cStart, &cEnd));
134: color = (cStart < cEnd) ? 0 : 1;
135: PetscCallMPI(MPI_Comm_rank(dmcomm, &rank));
136: PetscCallMPI(MPI_Comm_split(dmcomm, color, rank, &mesh->nonempty_comm));
137: if (color == 1) {
138: PetscCallMPI(MPI_Comm_free(&mesh->nonempty_comm));
139: mesh->nonempty_comm = MPI_COMM_NULL;
140: }
141: }
142: *comm = mesh->nonempty_comm;
143: PetscFunctionReturn(PETSC_SUCCESS);
144: }
146: static PetscErrorCode DMGetFieldIfFV_Private(DM dm, PetscInt field, PetscFV *fv)
147: {
148: PetscObject f = NULL;
149: PetscClassId fClass = PETSC_SMALLEST_CLASSID;
150: PetscInt nf;
152: PetscFunctionBegin;
153: *fv = NULL;
154: PetscCall(DMGetNumFields(dm, &nf));
155: if (nf > 0) {
156: PetscCall(DMGetField(dm, field, NULL, &f));
157: PetscCall(PetscObjectGetClassId(f, &fClass));
158: if (fClass == PETSCFV_CLASSID) *fv = (PetscFV)f;
159: }
160: PetscFunctionReturn(PETSC_SUCCESS);
161: }
163: /*
164: Write all fields that have been provided to the viewer
165: Multi-block XML format with binary appended data.
166: */
167: PetscErrorCode DMPlexVTKWriteAll_VTU(DM dm, PetscViewer viewer)
168: {
169: MPI_Comm comm;
170: PetscSection coordSection, cellCoordSection;
171: PetscViewer_VTK *vtk = (PetscViewer_VTK *)viewer->data;
172: PetscViewerVTKObjectLink link;
173: FILE *fp;
174: PetscMPIInt rank, size, tag;
175: PetscInt dimEmbed, cellHeight, cStart, cEnd, vStart, vEnd, numLabelCells, hasLabel, c, v, i;
176: PetscBool localized;
177: PieceInfo piece, *gpiece = NULL;
178: void *buffer = NULL;
179: const char *byte_order = PetscBinaryBigEndian() ? "BigEndian" : "LittleEndian";
180: PetscInt loops_per_scalar;
182: PetscFunctionBegin;
183: PetscCall(DMGetCoordinateDim(dm, &dimEmbed));
184: PetscCall(DMPlexGetVTKCellHeight(dm, &cellHeight));
185: PetscCall(DMPlexGetHeightStratum(dm, cellHeight, &cStart, &cEnd));
186: PetscCall(DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd));
187: PetscCall(DMGetStratumSize(dm, "vtk", 1, &numLabelCells));
188: PetscCall(DMGetCoordinatesLocalized(dm, &localized));
189: PetscCall(DMGetCoordinateSection(dm, &coordSection));
190: PetscCall(DMGetCellCoordinateSection(dm, &cellCoordSection));
191: PetscCall(PetscCommGetNewTag(PetscObjectComm((PetscObject)dm), &tag));
193: PetscCall(DMPlexGetNonEmptyComm_Private(dm, &comm));
194: #if defined(PETSC_USE_COMPLEX)
195: loops_per_scalar = 2;
196: #else
197: loops_per_scalar = 1;
198: #endif
199: if (comm == MPI_COMM_NULL) goto finalize;
200: PetscCallMPI(MPI_Comm_size(comm, &size));
201: PetscCallMPI(MPI_Comm_rank(comm, &rank));
203: PetscCall(PetscFOpen(comm, vtk->filename, "wb", &fp));
204: PetscCall(PetscFPrintf(comm, fp, "<?xml version=\"1.0\"?>\n"));
205: PetscCall(PetscFPrintf(comm, fp, "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"%s\" header_type=\"UInt64\">\n", byte_order));
206: PetscCall(PetscFPrintf(comm, fp, " <UnstructuredGrid>\n"));
208: hasLabel = numLabelCells > 0 ? PETSC_TRUE : PETSC_FALSE;
209: piece.nvertices = 0;
210: piece.ncells = 0;
211: piece.nconn = 0;
212: if (!localized) piece.nvertices = vEnd - vStart;
213: for (c = cStart; c < cEnd; ++c) {
214: PetscInt dof = 0;
215: if (hasLabel) {
216: PetscInt value;
218: PetscCall(DMGetLabelValue(dm, "vtk", c, &value));
219: if (value != 1) continue;
220: }
221: if (localized) PetscCall(PetscSectionGetDof(cellCoordSection, c, &dof));
222: if (!dof) {
223: PetscInt *closure = NULL;
224: PetscInt closureSize;
226: PetscCall(DMPlexGetTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure));
227: for (v = 0; v < closureSize * 2; v += 2) {
228: if ((closure[v] >= vStart) && (closure[v] < vEnd)) {
229: piece.nconn++;
230: if (localized) piece.nvertices++;
231: }
232: }
233: PetscCall(DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure));
234: } else {
235: piece.nvertices += dof / dimEmbed;
236: piece.nconn += dof / dimEmbed;
237: }
238: piece.ncells++;
239: }
240: if (rank == 0) PetscCall(PetscMalloc1(size, &gpiece));
241: PetscCallMPI(MPI_Gather((PetscInt *)&piece, sizeof(piece) / sizeof(PetscInt), MPIU_INT, (PetscInt *)gpiece, sizeof(piece) / sizeof(PetscInt), MPIU_INT, 0, comm));
243: /*
244: * Write file header
245: */
246: if (rank == 0) {
247: PetscInt64 boffset = 0;
249: for (PetscMPIInt r = 0; r < size; r++) {
250: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " <Piece NumberOfPoints=\"%" PetscInt_FMT "\" NumberOfCells=\"%" PetscInt_FMT "\">\n", gpiece[r].nvertices, gpiece[r].ncells));
251: /* Coordinate positions */
252: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " <Points>\n"));
253: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " <DataArray type=\"%s\" Name=\"Position\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, boffset));
254: boffset += gpiece[r].nvertices * 3 * sizeof(PetscVTUReal) + (gpiece[r].nvertices ? sizeof(PetscInt64) : 0);
255: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " </Points>\n"));
256: /* Cell connectivity */
257: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " <Cells>\n"));
258: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " <DataArray type=\"Int32\" Name=\"connectivity\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", boffset));
259: boffset += gpiece[r].nconn * sizeof(PetscVTKInt) + (gpiece[r].nconn ? sizeof(PetscInt64) : 0);
260: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " <DataArray type=\"Int32\" Name=\"offsets\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", boffset));
261: boffset += gpiece[r].ncells * sizeof(PetscVTKInt) + (gpiece[r].ncells ? sizeof(PetscInt64) : 0);
262: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " <DataArray type=\"UInt8\" Name=\"types\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", boffset));
263: boffset += gpiece[r].ncells * sizeof(unsigned char) + (gpiece[r].ncells ? sizeof(PetscInt64) : 0);
264: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " </Cells>\n"));
266: /*
267: * Cell Data headers
268: */
269: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " <CellData>\n"));
270: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " <DataArray type=\"Int32\" Name=\"Rank\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", boffset));
271: boffset += gpiece[r].ncells * sizeof(PetscVTKInt) + (gpiece[r].ncells ? sizeof(PetscInt64) : 0);
272: /* all the vectors */
273: for (link = vtk->link; link; link = link->next) {
274: Vec X = (Vec)link->vec;
275: DM dmX = NULL;
276: PetscInt bs = 1, nfields, field;
277: const char *vecname = "";
278: PetscSection section;
279: if ((link->ft != PETSC_VTK_CELL_FIELD) && (link->ft != PETSC_VTK_CELL_VECTOR_FIELD)) continue;
280: 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. */
281: PetscCall(PetscObjectGetName((PetscObject)X, &vecname));
282: }
283: PetscCall(VecGetDM(X, &dmX));
284: if (!dmX) dmX = dm;
285: PetscCall(PetscObjectQuery(link->vec, "section", (PetscObject *)§ion));
286: if (!section) PetscCall(DMGetLocalSection(dmX, §ion));
287: if (cEnd > cStart) PetscCall(PetscSectionGetDof(section, cStart, &bs));
288: PetscCall(PetscSectionGetNumFields(section, &nfields));
289: field = 0;
290: if (link->field >= 0) {
291: field = link->field;
292: nfields = field + 1;
293: }
294: for (i = 0; field < (nfields ? nfields : 1); field++) {
295: PetscInt fbs, j;
296: PetscFV fv = NULL;
297: const char *fieldname = NULL;
298: char buf[256];
299: PetscBool vector;
301: if (nfields) { /* We have user-defined fields/components */
302: PetscCall(PetscSectionGetFieldDof(section, cStart, field, &fbs));
303: PetscCall(PetscSectionGetFieldName(section, field, &fieldname));
304: } else fbs = bs; /* Say we have one field with 'bs' components */
305: PetscCall(DMGetFieldIfFV_Private(dmX, field, &fv));
306: if (nfields && !fieldname) {
307: PetscCall(PetscSNPrintf(buf, sizeof(buf), "CellField%" PetscInt_FMT, field));
308: fieldname = buf;
309: }
310: vector = PETSC_FALSE;
311: if (link->ft == PETSC_VTK_CELL_VECTOR_FIELD) {
312: vector = PETSC_TRUE;
313: PetscCheck(fbs <= 3, PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_SIZ, "Cell vector fields can have at most 3 components, %" PetscInt_FMT " given", fbs);
314: for (j = 0; j < fbs; j++) {
315: const char *compName = NULL;
316: if (fv) {
317: PetscCall(PetscFVGetComponentName(fv, j, &compName));
318: if (compName) break;
319: }
320: }
321: if (j < fbs) vector = PETSC_FALSE;
322: }
323: if (vector) {
324: #if defined(PETSC_USE_COMPLEX)
325: PetscCall(PetscFPrintf(comm, fp, " <DataArray type=\"%s\" Name=\"%s%s.Re\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, vecname, fieldname, boffset));
326: boffset += gpiece[r].ncells * 3 * sizeof(PetscVTUReal) + (gpiece[r].ncells ? sizeof(PetscInt64) : 0);
327: PetscCall(PetscFPrintf(comm, fp, " <DataArray type=\"%s\" Name=\"%s%s.Im\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, vecname, fieldname, boffset));
328: boffset += gpiece[r].ncells * 3 * sizeof(PetscVTUReal) + (gpiece[r].ncells ? sizeof(PetscInt64) : 0);
329: #else
330: PetscCall(PetscFPrintf(comm, fp, " <DataArray type=\"%s\" Name=\"%s%s\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, vecname, fieldname, boffset));
331: boffset += gpiece[r].ncells * 3 * sizeof(PetscVTUReal) + (gpiece[r].ncells ? sizeof(PetscInt64) : 0);
332: #endif
333: i += fbs;
334: } else {
335: for (j = 0; j < fbs; j++) {
336: const char *compName = NULL;
337: char finalname[256];
338: if (fv) PetscCall(PetscFVGetComponentName(fv, j, &compName));
339: if (compName) {
340: PetscCall(PetscSNPrintf(finalname, 255, "%s%s.%s", vecname, fieldname, compName));
341: } else if (fbs > 1) {
342: PetscCall(PetscSNPrintf(finalname, 255, "%s%s.%" PetscInt_FMT, vecname, fieldname, j));
343: } else {
344: PetscCall(PetscSNPrintf(finalname, 255, "%s%s", vecname, fieldname));
345: }
346: #if defined(PETSC_USE_COMPLEX)
347: PetscCall(PetscFPrintf(comm, fp, " <DataArray type=\"%s\" Name=\"%s.Re\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, finalname, boffset));
348: boffset += gpiece[r].ncells * sizeof(PetscVTUReal) + (gpiece[r].ncells ? sizeof(PetscInt64) : 0);
349: PetscCall(PetscFPrintf(comm, fp, " <DataArray type=\"%s\" Name=\"%s.Im\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, finalname, boffset));
350: boffset += gpiece[r].ncells * sizeof(PetscVTUReal) + (gpiece[r].ncells ? sizeof(PetscInt64) : 0);
351: #else
352: PetscCall(PetscFPrintf(comm, fp, " <DataArray type=\"%s\" Name=\"%s\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, finalname, boffset));
353: boffset += gpiece[r].ncells * sizeof(PetscVTUReal) + (gpiece[r].ncells ? sizeof(PetscInt64) : 0);
354: #endif
355: i++;
356: }
357: }
358: }
359: //PetscCheck(i == bs,comm,PETSC_ERR_PLIB,"Total number of field components %" PetscInt_FMT " != block size %" PetscInt_FMT,i,bs);
360: }
361: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " </CellData>\n"));
363: /*
364: * Point Data headers
365: */
366: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " <PointData>\n"));
367: for (link = vtk->link; link; link = link->next) {
368: Vec X = (Vec)link->vec;
369: DM dmX;
370: PetscInt bs = 1, nfields, field;
371: const char *vecname = "";
372: PetscSection section;
373: if ((link->ft != PETSC_VTK_POINT_FIELD) && (link->ft != PETSC_VTK_POINT_VECTOR_FIELD)) continue;
374: 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. */
375: PetscCall(PetscObjectGetName((PetscObject)X, &vecname));
376: }
377: PetscCall(VecGetDM(X, &dmX));
378: if (!dmX) dmX = dm;
379: PetscCall(PetscObjectQuery(link->vec, "section", (PetscObject *)§ion));
380: if (!section) PetscCall(DMGetLocalSection(dmX, §ion));
381: if (vEnd > vStart) PetscCall(PetscSectionGetDof(section, vStart, &bs));
382: PetscCall(PetscSectionGetNumFields(section, &nfields));
383: field = 0;
384: if (link->field >= 0) {
385: field = link->field;
386: nfields = field + 1;
387: }
388: for (; field < (nfields ? nfields : 1); field++) {
389: PetscInt fbs, j;
390: const char *fieldname = NULL;
391: char buf[256];
392: if (nfields) { /* We have user-defined fields/components */
393: PetscCall(PetscSectionGetFieldDof(section, vStart, field, &fbs));
394: PetscCall(PetscSectionGetFieldName(section, field, &fieldname));
395: } else fbs = bs; /* Say we have one field with 'bs' components */
396: if (nfields && !fieldname) {
397: PetscCall(PetscSNPrintf(buf, sizeof(buf), "PointField%" PetscInt_FMT, field));
398: fieldname = buf;
399: }
400: if (link->ft == PETSC_VTK_POINT_VECTOR_FIELD) {
401: PetscCheck(fbs <= 3, PetscObjectComm((PetscObject)viewer), PETSC_ERR_ARG_SIZ, "Point vector fields can have at most 3 components, %" PetscInt_FMT " given", fbs);
402: #if defined(PETSC_USE_COMPLEX)
403: PetscCall(PetscFPrintf(comm, fp, " <DataArray type=\"%s\" Name=\"%s%s.Re\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, vecname, fieldname, boffset));
404: boffset += gpiece[r].nvertices * 3 * sizeof(PetscVTUReal) + (gpiece[r].nvertices ? sizeof(PetscInt64) : 0);
405: PetscCall(PetscFPrintf(comm, fp, " <DataArray type=\"%s\" Name=\"%s%s.Im\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, vecname, fieldname, boffset));
406: boffset += gpiece[r].nvertices * 3 * sizeof(PetscVTUReal) + (gpiece[r].nvertices ? sizeof(PetscInt64) : 0);
407: #else
408: PetscCall(PetscFPrintf(comm, fp, " <DataArray type=\"%s\" Name=\"%s%s\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, vecname, fieldname, boffset));
409: boffset += gpiece[r].nvertices * 3 * sizeof(PetscVTUReal) + (gpiece[r].nvertices ? sizeof(PetscInt64) : 0);
410: #endif
411: } else {
412: for (j = 0; j < fbs; j++) {
413: const char *compName = NULL;
414: char finalname[256];
415: PetscCall(PetscSectionGetComponentName(section, field, j, &compName));
416: PetscCall(PetscSNPrintf(finalname, 255, "%s%s.%s", vecname, fieldname, compName));
417: #if defined(PETSC_USE_COMPLEX)
418: PetscCall(PetscFPrintf(comm, fp, " <DataArray type=\"%s\" Name=\"%s.Re\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, finalname, boffset));
419: boffset += gpiece[r].nvertices * sizeof(PetscVTUReal) + (gpiece[r].nvertices ? sizeof(PetscInt64) : 0);
420: PetscCall(PetscFPrintf(comm, fp, " <DataArray type=\"%s\" Name=\"%s.Im\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, finalname, boffset));
421: boffset += gpiece[r].nvertices * sizeof(PetscVTUReal) + (gpiece[r].nvertices ? sizeof(PetscInt64) : 0);
422: #else
423: PetscCall(PetscFPrintf(comm, fp, " <DataArray type=\"%s\" Name=\"%s\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%" PetscInt64_FMT "\" />\n", precision, finalname, boffset));
424: boffset += gpiece[r].nvertices * sizeof(PetscVTUReal) + (gpiece[r].nvertices ? sizeof(PetscInt64) : 0);
425: #endif
426: }
427: }
428: }
429: }
430: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " </PointData>\n"));
431: PetscCall(PetscFPrintf(PETSC_COMM_SELF, fp, " </Piece>\n"));
432: }
433: }
435: PetscCall(PetscFPrintf(comm, fp, " </UnstructuredGrid>\n"));
436: PetscCall(PetscFPrintf(comm, fp, " <AppendedData encoding=\"raw\">\n"));
437: PetscCall(PetscFPrintf(comm, fp, "_"));
439: if (rank == 0) {
440: PetscInt maxsize = 0;
441: for (PetscMPIInt r = 0; r < size; r++) {
442: maxsize = PetscMax(maxsize, (PetscInt)(gpiece[r].nvertices * 3 * sizeof(PetscVTUReal)));
443: maxsize = PetscMax(maxsize, (PetscInt)(gpiece[r].ncells * 3 * sizeof(PetscVTUReal)));
444: maxsize = PetscMax(maxsize, (PetscInt)(gpiece[r].nconn * sizeof(PetscVTKInt)));
445: }
446: PetscCall(PetscMalloc(maxsize, &buffer));
447: }
448: for (PetscMPIInt r = 0; r < size; r++) {
449: if (r == rank) {
450: PetscInt nsend;
451: { /* Position */
452: const PetscScalar *x, *cx = NULL;
453: PetscVTUReal *y = NULL;
454: Vec coords, cellCoords;
455: PetscBool copy;
457: PetscCall(DMGetCoordinatesLocal(dm, &coords));
458: PetscCall(VecGetArrayRead(coords, &x));
459: PetscCall(DMGetCellCoordinatesLocal(dm, &cellCoords));
460: if (cellCoords) PetscCall(VecGetArrayRead(cellCoords, &cx));
461: #if defined(PETSC_USE_COMPLEX)
462: copy = PETSC_TRUE;
463: #else
464: copy = (PetscBool)(dimEmbed != 3 || localized || (sizeof(PetscReal) != sizeof(PetscVTUReal)));
465: #endif
466: if (copy) {
467: PetscCall(PetscMalloc1(piece.nvertices * 3, &y));
468: if (localized) {
469: PetscInt cnt;
470: for (c = cStart, cnt = 0; c < cEnd; c++) {
471: PetscInt off, dof;
473: PetscCall(PetscSectionGetDof(cellCoordSection, c, &dof));
474: if (!dof) {
475: PetscInt *closure = NULL;
476: PetscInt closureSize;
478: PetscCall(DMPlexGetTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure));
479: for (v = 0; v < closureSize * 2; v += 2) {
480: if ((closure[v] >= vStart) && (closure[v] < vEnd)) {
481: PetscCall(PetscSectionGetOffset(coordSection, closure[v], &off));
482: if (dimEmbed != 3) {
483: y[cnt * 3 + 0] = (PetscVTUReal)PetscRealPart(x[off + 0]);
484: y[cnt * 3 + 1] = (PetscVTUReal)((dimEmbed > 1) ? PetscRealPart(x[off + 1]) : 0.0);
485: y[cnt * 3 + 2] = (PetscVTUReal)0.0;
486: } else {
487: y[cnt * 3 + 0] = (PetscVTUReal)PetscRealPart(x[off + 0]);
488: y[cnt * 3 + 1] = (PetscVTUReal)PetscRealPart(x[off + 1]);
489: y[cnt * 3 + 2] = (PetscVTUReal)PetscRealPart(x[off + 2]);
490: }
491: cnt++;
492: }
493: }
494: PetscCall(DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure));
495: } else {
496: PetscCall(PetscSectionGetOffset(cellCoordSection, c, &off));
497: if (dimEmbed != 3) {
498: for (i = 0; i < dof / dimEmbed; i++) {
499: y[cnt * 3 + 0] = (PetscVTUReal)PetscRealPart(cx[off + i * dimEmbed + 0]);
500: y[cnt * 3 + 1] = (PetscVTUReal)((dimEmbed > 1) ? PetscRealPart(cx[off + i * dimEmbed + 1]) : 0.0);
501: y[cnt * 3 + 2] = (PetscVTUReal)0.0;
502: cnt++;
503: }
504: } else {
505: for (i = 0; i < dof; i++) y[cnt * 3 + i] = (PetscVTUReal)PetscRealPart(cx[off + i]);
506: cnt += dof / dimEmbed;
507: }
508: }
509: }
510: PetscCheck(cnt == piece.nvertices, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Count does not match");
511: } else {
512: for (i = 0; i < piece.nvertices; i++) {
513: y[i * 3 + 0] = (PetscVTUReal)PetscRealPart(x[i * dimEmbed + 0]);
514: y[i * 3 + 1] = (PetscVTUReal)((dimEmbed > 1) ? PetscRealPart(x[i * dimEmbed + 1]) : 0.);
515: y[i * 3 + 2] = (PetscVTUReal)((dimEmbed > 2) ? PetscRealPart(x[i * dimEmbed + 2]) : 0.);
516: }
517: }
518: }
519: nsend = piece.nvertices * 3;
520: PetscCall(TransferWrite(comm, viewer, fp, r, 0, copy ? (const void *)y : (const void *)x, buffer, nsend, MPIU_VTUREAL, tag));
521: PetscCall(PetscFree(y));
522: PetscCall(VecRestoreArrayRead(coords, &x));
523: if (cellCoords) PetscCall(VecRestoreArrayRead(cellCoords, &cx));
524: }
525: { /* Connectivity, offsets, types */
526: PetscVTKInt *connectivity = NULL, *offsets = NULL;
527: PetscVTKType *types = NULL;
528: PetscCall(DMPlexGetVTKConnectivity(dm, localized, &piece, &connectivity, &offsets, &types));
529: PetscCall(TransferWrite(comm, viewer, fp, r, 0, connectivity, buffer, piece.nconn, MPI_INT, tag));
530: PetscCall(TransferWrite(comm, viewer, fp, r, 0, offsets, buffer, piece.ncells, MPI_INT, tag));
531: PetscCall(TransferWrite(comm, viewer, fp, r, 0, types, buffer, piece.ncells, MPI_CHAR, tag));
532: PetscCall(PetscFree3(connectivity, offsets, types));
533: }
534: { /* Owners (cell data) */
535: PetscVTKInt *owners;
536: PetscMPIInt orank;
538: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &orank));
539: PetscCall(PetscMalloc1(piece.ncells, &owners));
540: for (i = 0; i < piece.ncells; i++) owners[i] = orank;
541: PetscCall(TransferWrite(comm, viewer, fp, r, 0, owners, buffer, piece.ncells, MPI_INT, tag));
542: PetscCall(PetscFree(owners));
543: }
544: /* Cell data */
545: for (link = vtk->link; link; link = link->next) {
546: Vec X = (Vec)link->vec;
547: DM dmX;
548: const PetscScalar *x;
549: PetscVTUReal *y;
550: PetscInt bs = 1, nfields, field;
551: PetscSection section = NULL;
553: if ((link->ft != PETSC_VTK_CELL_FIELD) && (link->ft != PETSC_VTK_CELL_VECTOR_FIELD)) continue;
554: PetscCall(VecGetDM(X, &dmX));
555: if (!dmX) dmX = dm;
556: PetscCall(PetscObjectQuery(link->vec, "section", (PetscObject *)§ion));
557: if (!section) PetscCall(DMGetLocalSection(dmX, §ion));
558: if (cEnd > cStart) PetscCall(PetscSectionGetDof(section, cStart, &bs));
559: PetscCall(PetscSectionGetNumFields(section, &nfields));
560: field = 0;
561: if (link->field >= 0) {
562: field = link->field;
563: nfields = field + 1;
564: }
565: PetscCall(VecGetArrayRead(X, &x));
566: PetscCall(PetscMalloc1(piece.ncells * 3, &y));
567: for (; field < (nfields ? nfields : 1); field++) {
568: PetscInt fbs, j;
569: PetscFV fv = NULL;
570: PetscBool vector;
572: if (nfields && cEnd > cStart) { /* We have user-defined fields/components */
573: PetscCall(PetscSectionGetFieldDof(section, cStart, field, &fbs));
574: } else fbs = bs; /* Say we have one field with 'bs' components */
575: PetscCall(DMGetFieldIfFV_Private(dmX, field, &fv));
576: vector = PETSC_FALSE;
577: if (link->ft == PETSC_VTK_CELL_VECTOR_FIELD) {
578: vector = PETSC_TRUE;
579: for (j = 0; j < fbs; j++) {
580: const char *compName = NULL;
581: if (fv) {
582: PetscCall(PetscFVGetComponentName(fv, j, &compName));
583: if (compName) break;
584: }
585: }
586: if (j < fbs) vector = PETSC_FALSE;
587: }
588: if (vector) {
589: PetscInt cnt, l;
590: for (l = 0; l < loops_per_scalar; l++) {
591: for (c = cStart, cnt = 0; c < cEnd; c++) {
592: const PetscScalar *xpoint;
593: PetscInt off, j;
595: if (hasLabel) { /* Ignore some cells */
596: PetscInt value;
597: PetscCall(DMGetLabelValue(dmX, "vtk", c, &value));
598: if (value != 1) continue;
599: }
600: if (nfields) {
601: PetscCall(PetscSectionGetFieldOffset(section, c, field, &off));
602: } else {
603: PetscCall(PetscSectionGetOffset(section, c, &off));
604: }
605: xpoint = &x[off];
606: for (j = 0; j < fbs; j++) y[cnt++] = (PetscVTUReal)(l ? PetscImaginaryPart(xpoint[j]) : PetscRealPart(xpoint[j]));
607: for (; j < 3; j++) y[cnt++] = 0.;
608: }
609: PetscCheck(cnt == piece.ncells * 3, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Count does not match");
610: PetscCall(TransferWrite(comm, viewer, fp, r, 0, y, buffer, piece.ncells * 3, MPIU_VTUREAL, tag));
611: }
612: } else {
613: for (i = 0; i < fbs; i++) {
614: PetscInt cnt, l;
615: for (l = 0; l < loops_per_scalar; l++) {
616: for (c = cStart, cnt = 0; c < cEnd; c++) {
617: const PetscScalar *xpoint;
618: PetscInt off;
620: if (hasLabel) { /* Ignore some cells */
621: PetscInt value;
622: PetscCall(DMGetLabelValue(dmX, "vtk", c, &value));
623: if (value != 1) continue;
624: }
625: if (nfields) {
626: PetscCall(PetscSectionGetFieldOffset(section, c, field, &off));
627: } else {
628: PetscCall(PetscSectionGetOffset(section, c, &off));
629: }
630: xpoint = &x[off];
631: y[cnt++] = (PetscVTUReal)(l ? PetscImaginaryPart(xpoint[i]) : PetscRealPart(xpoint[i]));
632: }
633: PetscCheck(cnt == piece.ncells, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Count does not match");
634: PetscCall(TransferWrite(comm, viewer, fp, r, 0, y, buffer, piece.ncells, MPIU_VTUREAL, tag));
635: }
636: }
637: }
638: }
639: PetscCall(PetscFree(y));
640: PetscCall(VecRestoreArrayRead(X, &x));
641: }
642: /* point data */
643: for (link = vtk->link; link; link = link->next) {
644: Vec X = (Vec)link->vec;
645: DM dmX;
646: const PetscScalar *x;
647: PetscVTUReal *y;
648: PetscInt bs = 1, nfields, field;
649: PetscSection section = NULL;
651: if ((link->ft != PETSC_VTK_POINT_FIELD) && (link->ft != PETSC_VTK_POINT_VECTOR_FIELD)) continue;
652: PetscCall(VecGetDM(X, &dmX));
653: if (!dmX) dmX = dm;
654: PetscCall(PetscObjectQuery(link->vec, "section", (PetscObject *)§ion));
655: if (!section) PetscCall(DMGetLocalSection(dmX, §ion));
656: if (vEnd > vStart) PetscCall(PetscSectionGetDof(section, vStart, &bs));
657: PetscCall(PetscSectionGetNumFields(section, &nfields));
658: field = 0;
659: if (link->field >= 0) {
660: field = link->field;
661: nfields = field + 1;
662: }
663: PetscCall(VecGetArrayRead(X, &x));
664: PetscCall(PetscMalloc1(piece.nvertices * 3, &y));
665: for (; field < (nfields ? nfields : 1); field++) {
666: PetscInt fbs, j;
667: if (nfields && vEnd > vStart) { /* We have user-defined fields/components */
668: PetscCall(PetscSectionGetFieldDof(section, vStart, field, &fbs));
669: } else fbs = bs; /* Say we have one field with 'bs' components */
670: if (link->ft == PETSC_VTK_POINT_VECTOR_FIELD) {
671: PetscInt cnt, l;
672: for (l = 0; l < loops_per_scalar; l++) {
673: if (!localized) {
674: for (v = vStart, cnt = 0; v < vEnd; v++) {
675: PetscInt off;
676: const PetscScalar *xpoint;
678: if (nfields) {
679: PetscCall(PetscSectionGetFieldOffset(section, v, field, &off));
680: } else {
681: PetscCall(PetscSectionGetOffset(section, v, &off));
682: }
683: xpoint = &x[off];
684: for (j = 0; j < fbs; j++) y[cnt++] = (PetscVTUReal)(l ? PetscImaginaryPart(xpoint[j]) : PetscRealPart(xpoint[j]));
685: for (; j < 3; j++) y[cnt++] = 0.;
686: }
687: } else {
688: for (c = cStart, cnt = 0; c < cEnd; c++) {
689: PetscInt *closure = NULL;
690: PetscInt closureSize, off;
692: PetscCall(DMPlexGetTransitiveClosure(dmX, c, PETSC_TRUE, &closureSize, &closure));
693: for (v = 0, off = 0; v < closureSize * 2; v += 2) {
694: if ((closure[v] >= vStart) && (closure[v] < vEnd)) {
695: PetscInt voff;
696: const PetscScalar *xpoint;
698: if (nfields) {
699: PetscCall(PetscSectionGetFieldOffset(section, closure[v], field, &voff));
700: } else {
701: PetscCall(PetscSectionGetOffset(section, closure[v], &voff));
702: }
703: xpoint = &x[voff];
704: for (j = 0; j < fbs; j++) y[cnt + off++] = (PetscVTUReal)(l ? PetscImaginaryPart(xpoint[j]) : PetscRealPart(xpoint[j]));
705: for (; j < 3; j++) y[cnt + off++] = 0.;
706: }
707: }
708: cnt += off;
709: PetscCall(DMPlexRestoreTransitiveClosure(dmX, c, PETSC_TRUE, &closureSize, &closure));
710: }
711: }
712: PetscCheck(cnt == piece.nvertices * 3, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Count does not match");
713: PetscCall(TransferWrite(comm, viewer, fp, r, 0, y, buffer, piece.nvertices * 3, MPIU_VTUREAL, tag));
714: }
715: } else {
716: for (i = 0; i < fbs; i++) {
717: PetscInt cnt, l;
718: for (l = 0; l < loops_per_scalar; l++) {
719: if (!localized) {
720: for (v = vStart, cnt = 0; v < vEnd; v++) {
721: PetscInt off;
722: const PetscScalar *xpoint;
724: if (nfields) {
725: PetscCall(PetscSectionGetFieldOffset(section, v, field, &off));
726: } else {
727: PetscCall(PetscSectionGetOffset(section, v, &off));
728: }
729: xpoint = &x[off];
730: y[cnt++] = (PetscVTUReal)(l ? PetscImaginaryPart(xpoint[i]) : PetscRealPart(xpoint[i]));
731: }
732: } else {
733: for (c = cStart, cnt = 0; c < cEnd; c++) {
734: PetscInt *closure = NULL;
735: PetscInt closureSize, off;
737: PetscCall(DMPlexGetTransitiveClosure(dmX, c, PETSC_TRUE, &closureSize, &closure));
738: for (v = 0, off = 0; v < closureSize * 2; v += 2) {
739: if ((closure[v] >= vStart) && (closure[v] < vEnd)) {
740: PetscInt voff;
741: const PetscScalar *xpoint;
743: if (nfields) {
744: PetscCall(PetscSectionGetFieldOffset(section, closure[v], field, &voff));
745: } else {
746: PetscCall(PetscSectionGetOffset(section, closure[v], &voff));
747: }
748: xpoint = &x[voff];
749: y[cnt + off++] = (l ? PetscImaginaryPart(xpoint[i]) : PetscRealPart(xpoint[i]));
750: }
751: }
752: cnt += off;
753: PetscCall(DMPlexRestoreTransitiveClosure(dmX, c, PETSC_TRUE, &closureSize, &closure));
754: }
755: }
756: PetscCheck(cnt == piece.nvertices, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Count does not match");
757: PetscCall(TransferWrite(comm, viewer, fp, r, 0, y, buffer, piece.nvertices, MPIU_VTUREAL, tag));
758: }
759: }
760: }
761: }
762: PetscCall(PetscFree(y));
763: PetscCall(VecRestoreArrayRead(X, &x));
764: }
765: } else if (rank == 0) {
766: PetscInt l;
768: PetscCall(TransferWrite(comm, viewer, fp, r, 0, NULL, buffer, gpiece[r].nvertices * 3, MPIU_VTUREAL, tag)); /* positions */
769: PetscCall(TransferWrite(comm, viewer, fp, r, 0, NULL, buffer, gpiece[r].nconn, MPI_INT, tag)); /* connectivity */
770: PetscCall(TransferWrite(comm, viewer, fp, r, 0, NULL, buffer, gpiece[r].ncells, MPI_INT, tag)); /* offsets */
771: PetscCall(TransferWrite(comm, viewer, fp, r, 0, NULL, buffer, gpiece[r].ncells, MPI_CHAR, tag)); /* types */
772: PetscCall(TransferWrite(comm, viewer, fp, r, 0, NULL, buffer, gpiece[r].ncells, MPI_INT, tag)); /* owner rank (cells) */
773: /* all cell data */
774: for (link = vtk->link; link; link = link->next) {
775: Vec X = (Vec)link->vec;
776: PetscInt bs = 1, nfields, field;
777: DM dmX;
778: PetscSection section = NULL;
780: if ((link->ft != PETSC_VTK_CELL_FIELD) && (link->ft != PETSC_VTK_CELL_VECTOR_FIELD)) continue;
781: PetscCall(VecGetDM(X, &dmX));
782: if (!dmX) dmX = dm;
783: PetscCall(PetscObjectQuery(link->vec, "section", (PetscObject *)§ion));
784: if (!section) PetscCall(DMGetLocalSection(dmX, §ion));
785: if (cEnd > cStart) PetscCall(PetscSectionGetDof(section, cStart, &bs));
786: PetscCall(PetscSectionGetNumFields(section, &nfields));
787: field = 0;
788: if (link->field >= 0) {
789: field = link->field;
790: nfields = field + 1;
791: }
792: for (; field < (nfields ? nfields : 1); field++) {
793: PetscInt fbs, j;
794: PetscFV fv = NULL;
795: PetscBool vector;
797: if (nfields && cEnd > cStart) { /* We have user-defined fields/components */
798: PetscCall(PetscSectionGetFieldDof(section, cStart, field, &fbs));
799: } else fbs = bs; /* Say we have one field with 'bs' components */
800: PetscCall(DMGetFieldIfFV_Private(dmX, field, &fv));
801: vector = PETSC_FALSE;
802: if (link->ft == PETSC_VTK_CELL_VECTOR_FIELD) {
803: vector = PETSC_TRUE;
804: for (j = 0; j < fbs; j++) {
805: const char *compName = NULL;
806: if (fv) {
807: PetscCall(PetscFVGetComponentName(fv, j, &compName));
808: if (compName) break;
809: }
810: }
811: if (j < fbs) vector = PETSC_FALSE;
812: }
813: if (vector) {
814: for (l = 0; l < loops_per_scalar; l++) PetscCall(TransferWrite(comm, viewer, fp, r, 0, NULL, buffer, gpiece[r].ncells * 3, MPIU_VTUREAL, tag));
815: } else {
816: for (i = 0; i < fbs; i++) {
817: for (l = 0; l < loops_per_scalar; l++) PetscCall(TransferWrite(comm, viewer, fp, r, 0, NULL, buffer, gpiece[r].ncells, MPIU_VTUREAL, tag));
818: }
819: }
820: }
821: }
822: /* all point data */
823: for (link = vtk->link; link; link = link->next) {
824: Vec X = (Vec)link->vec;
825: DM dmX;
826: PetscInt bs = 1, nfields, field;
827: PetscSection section = NULL;
829: if ((link->ft != PETSC_VTK_POINT_FIELD) && (link->ft != PETSC_VTK_POINT_VECTOR_FIELD)) continue;
830: PetscCall(VecGetDM(X, &dmX));
831: if (!dmX) dmX = dm;
832: PetscCall(PetscObjectQuery(link->vec, "section", (PetscObject *)§ion));
833: if (!section) PetscCall(DMGetLocalSection(dmX, §ion));
834: if (vEnd > vStart) PetscCall(PetscSectionGetDof(section, vStart, &bs));
835: PetscCall(PetscSectionGetNumFields(section, &nfields));
836: field = 0;
837: if (link->field >= 0) {
838: field = link->field;
839: nfields = field + 1;
840: }
841: for (; field < (nfields ? nfields : 1); field++) {
842: PetscInt fbs;
843: if (nfields && vEnd > vStart) { /* We have user-defined fields/components */
844: PetscCall(PetscSectionGetFieldDof(section, vStart, field, &fbs));
845: } else fbs = bs; /* Say we have one field with 'bs' components */
846: if (link->ft == PETSC_VTK_POINT_VECTOR_FIELD) {
847: for (l = 0; l < loops_per_scalar; l++) PetscCall(TransferWrite(comm, viewer, fp, r, 0, NULL, buffer, gpiece[r].nvertices * 3, MPIU_VTUREAL, tag));
848: } else {
849: for (i = 0; i < fbs; i++) {
850: for (l = 0; l < loops_per_scalar; l++) PetscCall(TransferWrite(comm, viewer, fp, r, 0, NULL, buffer, gpiece[r].nvertices, MPIU_VTUREAL, tag));
851: }
852: }
853: }
854: }
855: }
856: }
857: PetscCall(PetscFree(gpiece));
858: PetscCall(PetscFree(buffer));
859: PetscCall(PetscFPrintf(comm, fp, "\n </AppendedData>\n"));
860: PetscCall(PetscFPrintf(comm, fp, "</VTKFile>\n"));
861: PetscCall(PetscFClose(comm, fp));
862: finalize:
863: /* this code sends to rank 0 that writes.
864: It may lead to very unbalanced log_view timings
865: of the next PETSc function logged.
866: Since this call is not performance critical, we
867: issue a barrier here to synchronize the processes */
868: PetscCall(PetscBarrier((PetscObject)viewer));
869: PetscFunctionReturn(PETSC_SUCCESS);
870: }