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