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 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 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, (void **)&bytes)); 58 PetscCall(PetscFree(bytes)); 59 PetscCall(PetscContainerDestroy(&container)); 60 } 61 PetscCall(PetscViewerDestroy(&viewer)); 62 *v = NULL; 63 PetscFunctionReturn(PETSC_SUCCESS); 64 } 65 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 84 static PetscErrorCode private_DMSwarmView_XDMF(DM dm, PetscViewer viewer) 85 { 86 PetscBool isswarm = PETSC_FALSE; 87 const char *viewername; 88 char datafile[PETSC_MAX_PATH_LEN]; 89 char *datafilename; 90 PetscViewer fviewer; 91 PetscInt k, ng, dim; 92 Vec dvec; 93 long int *bytes = NULL; 94 PetscContainer container = NULL; 95 const char *dmname; 96 97 PetscFunctionBegin; 98 PetscCall(PetscObjectQuery((PetscObject)viewer, "XDMFViewerContext", (PetscObject *)&container)); 99 if (container) { 100 PetscCall(PetscContainerGetPointer(container, (void **)&bytes)); 101 } else SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Valid to find attached data XDMFViewerContext"); 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 /* write geometry data */ 180 PetscCall(DMSwarmCreateGlobalVectorFromField(dm, DMSwarmPICField_coor, &dvec)); 181 PetscCall(VecView(dvec, fviewer)); 182 PetscCall(DMSwarmDestroyGlobalVectorFromField(dm, DMSwarmPICField_coor, &dvec)); 183 bytes[0] += sizeof(PetscReal) * ng * dim; 184 185 PetscCall(PetscViewerDestroy(&fviewer)); 186 PetscFunctionReturn(PETSC_SUCCESS); 187 } 188 189 static PetscErrorCode private_VecView_Swarm_XDMF(Vec x, PetscViewer viewer) 190 { 191 long int *bytes = NULL; 192 PetscContainer container = NULL; 193 const char *viewername; 194 char datafile[PETSC_MAX_PATH_LEN]; 195 char *datafilename; 196 PetscViewer fviewer; 197 PetscInt N, bs; 198 const char *vecname; 199 char fieldname[PETSC_MAX_PATH_LEN]; 200 201 PetscFunctionBegin; 202 PetscCall(PetscObjectQuery((PetscObject)viewer, "XDMFViewerContext", (PetscObject *)&container)); 203 PetscCheck(container, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Unable to find attached data XDMFViewerContext"); 204 PetscCall(PetscContainerGetPointer(container, (void **)&bytes)); 205 PetscCall(PetscViewerFileGetName(viewer, &viewername)); 206 PetscCall(private_CreateDataFileNameXDMF(viewername, datafile)); 207 208 /* re-open a sub-viewer for all data fields */ 209 /* name is viewer.name + "_swarm_fields.pbin" */ 210 PetscCall(PetscViewerCreate(PetscObjectComm((PetscObject)viewer), &fviewer)); 211 PetscCall(PetscViewerSetType(fviewer, PETSCVIEWERBINARY)); 212 PetscCall(PetscViewerBinarySetSkipHeader(fviewer, PETSC_TRUE)); 213 PetscCall(PetscViewerBinarySetSkipInfo(fviewer, PETSC_TRUE)); 214 PetscCall(PetscViewerFileSetMode(fviewer, FILE_MODE_APPEND)); 215 PetscCall(PetscViewerFileSetName(fviewer, datafile)); 216 PetscCall(PetscStrrchr(datafile, '/', &datafilename)); 217 218 PetscCall(VecGetSize(x, &N)); 219 PetscCall(VecGetBlockSize(x, &bs)); 220 N = N / bs; 221 PetscCall(PetscObjectGetName((PetscObject)x, &vecname)); 222 if (!vecname) { 223 PetscCall(PetscSNPrintf(fieldname, PETSC_MAX_PATH_LEN - 1, "swarmfield_%d", ((PetscObject)x)->tag)); 224 } else { 225 PetscCall(PetscSNPrintf(fieldname, PETSC_MAX_PATH_LEN - 1, "%s", vecname)); 226 } 227 228 /* write data header */ 229 PetscCall(PetscViewerASCIIPushTab(viewer)); 230 PetscCall(PetscViewerASCIIPrintf(viewer, "<Attribute Center=\"Node\" Name=\"%s\" Type=\"None\">\n", fieldname)); 231 PetscCall(PetscViewerASCIIPushTab(viewer)); 232 if (bs == 1) { 233 PetscCall(PetscViewerASCIIPrintf(viewer, "<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Float\" Precision=\"8\" Dimensions=\"%" PetscInt_FMT "\" Seek=\"%ld\">\n", N, bytes[0])); 234 } else { 235 PetscCall(PetscViewerASCIIPrintf(viewer, "<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Float\" Precision=\"8\" Dimensions=\"%" PetscInt_FMT " %" PetscInt_FMT "\" Seek=\"%ld\">\n", N, bs, bytes[0])); 236 } 237 PetscCall(PetscViewerASCIIPushTab(viewer)); 238 PetscCall(PetscViewerASCIIPrintf(viewer, "%s\n", datafilename)); 239 PetscCall(PetscViewerASCIIPopTab(viewer)); 240 PetscCall(PetscViewerASCIIPrintf(viewer, "</DataItem>\n")); 241 PetscCall(PetscViewerASCIIPopTab(viewer)); 242 PetscCall(PetscViewerASCIIPrintf(viewer, "</Attribute>\n")); 243 PetscCall(PetscViewerASCIIPopTab(viewer)); 244 245 /* write data */ 246 PetscCall(VecView(x, fviewer)); 247 bytes[0] += sizeof(PetscReal) * N * bs; 248 249 PetscCall(PetscViewerDestroy(&fviewer)); 250 PetscFunctionReturn(PETSC_SUCCESS); 251 } 252 253 static PetscErrorCode private_ISView_Swarm_XDMF(IS is, PetscViewer viewer) 254 { 255 long int *bytes = NULL; 256 PetscContainer container = NULL; 257 const char *viewername; 258 char datafile[PETSC_MAX_PATH_LEN]; 259 char *datafilename; 260 PetscViewer fviewer; 261 PetscInt N, bs; 262 const char *vecname; 263 char fieldname[PETSC_MAX_PATH_LEN]; 264 265 PetscFunctionBegin; 266 PetscCall(PetscObjectQuery((PetscObject)viewer, "XDMFViewerContext", (PetscObject *)&container)); 267 PetscCheck(container, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Unable to find attached data XDMFViewerContext"); 268 PetscCall(PetscContainerGetPointer(container, (void **)&bytes)); 269 PetscCall(PetscViewerFileGetName(viewer, &viewername)); 270 PetscCall(private_CreateDataFileNameXDMF(viewername, datafile)); 271 272 /* re-open a sub-viewer for all data fields */ 273 /* name is viewer.name + "_swarm_fields.pbin" */ 274 PetscCall(PetscViewerCreate(PetscObjectComm((PetscObject)viewer), &fviewer)); 275 PetscCall(PetscViewerSetType(fviewer, PETSCVIEWERBINARY)); 276 PetscCall(PetscViewerBinarySetSkipHeader(fviewer, PETSC_TRUE)); 277 PetscCall(PetscViewerBinarySetSkipInfo(fviewer, PETSC_TRUE)); 278 PetscCall(PetscViewerFileSetMode(fviewer, FILE_MODE_APPEND)); 279 PetscCall(PetscViewerFileSetName(fviewer, datafile)); 280 PetscCall(PetscStrrchr(datafile, '/', &datafilename)); 281 282 PetscCall(ISGetSize(is, &N)); 283 PetscCall(ISGetBlockSize(is, &bs)); 284 N = N / bs; 285 PetscCall(PetscObjectGetName((PetscObject)is, &vecname)); 286 if (!vecname) { 287 PetscCall(PetscSNPrintf(fieldname, PETSC_MAX_PATH_LEN - 1, "swarmfield_%d", ((PetscObject)is)->tag)); 288 } else { 289 PetscCall(PetscSNPrintf(fieldname, PETSC_MAX_PATH_LEN - 1, "%s", vecname)); 290 } 291 292 /* write data header */ 293 PetscCall(PetscViewerASCIIPushTab(viewer)); 294 PetscCall(PetscViewerASCIIPrintf(viewer, "<Attribute Center=\"Node\" Name=\"%s\" Type=\"None\">\n", fieldname)); 295 PetscCall(PetscViewerASCIIPushTab(viewer)); 296 if (bs == 1) { 297 PetscCall(PetscViewerASCIIPrintf(viewer, "<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Int\" Precision=\"4\" Dimensions=\"%" PetscInt_FMT "\" Seek=\"%ld\">\n", N, bytes[0])); 298 } else { 299 PetscCall(PetscViewerASCIIPrintf(viewer, "<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Int\" Precision=\"4\" Dimensions=\"%" PetscInt_FMT " %" PetscInt_FMT "\" Seek=\"%ld\">\n", N, bs, bytes[0])); 300 } 301 PetscCall(PetscViewerASCIIPushTab(viewer)); 302 PetscCall(PetscViewerASCIIPrintf(viewer, "%s\n", datafilename)); 303 PetscCall(PetscViewerASCIIPopTab(viewer)); 304 PetscCall(PetscViewerASCIIPrintf(viewer, "</DataItem>\n")); 305 PetscCall(PetscViewerASCIIPopTab(viewer)); 306 PetscCall(PetscViewerASCIIPrintf(viewer, "</Attribute>\n")); 307 PetscCall(PetscViewerASCIIPopTab(viewer)); 308 309 /* write data */ 310 PetscCall(ISView(is, fviewer)); 311 bytes[0] += sizeof(PetscInt) * N * bs; 312 313 PetscCall(PetscViewerDestroy(&fviewer)); 314 PetscFunctionReturn(PETSC_SUCCESS); 315 } 316 317 /*@C 318 DMSwarmViewFieldsXDMF - Write a selection of DMSwarm fields to an XDMF3 file 319 320 Collective 321 322 Input Parameters: 323 + dm - the `DMSWARM` 324 . filename - the file name of the XDMF file (must have the extension .xmf) 325 . nfields - the number of fields to write into the XDMF file 326 - field_name_list - array of length nfields containing the textual name of fields to write 327 328 Level: beginner 329 330 Note: 331 Only fields registered with data type `PETSC_DOUBLE` or `PETSC_INT` can be written into the file 332 333 .seealso: `DM`, `DMSWARM`, `DMSwarmViewXDMF()` 334 @*/ 335 PetscErrorCode DMSwarmViewFieldsXDMF(DM dm, const char filename[], PetscInt nfields, const char *field_name_list[]) 336 { 337 Vec dvec; 338 PetscInt f, N; 339 PetscViewer viewer; 340 341 PetscFunctionBegin; 342 PetscCall(private_PetscViewerCreate_XDMF(PetscObjectComm((PetscObject)dm), filename, &viewer)); 343 PetscCall(private_DMSwarmView_XDMF(dm, viewer)); 344 PetscCall(DMSwarmGetLocalSize(dm, &N)); 345 for (f = 0; f < nfields; f++) { 346 void *data; 347 PetscDataType type; 348 349 PetscCall(DMSwarmGetField(dm, field_name_list[f], NULL, &type, &data)); 350 PetscCall(DMSwarmRestoreField(dm, field_name_list[f], NULL, &type, &data)); 351 if (type == PETSC_DOUBLE) { 352 PetscCall(DMSwarmCreateGlobalVectorFromField(dm, field_name_list[f], &dvec)); 353 PetscCall(PetscObjectSetName((PetscObject)dvec, field_name_list[f])); 354 PetscCall(private_VecView_Swarm_XDMF(dvec, viewer)); 355 PetscCall(DMSwarmDestroyGlobalVectorFromField(dm, field_name_list[f], &dvec)); 356 } else if (type == PETSC_INT) { 357 IS is; 358 const PetscInt *idx; 359 360 PetscCall(DMSwarmGetField(dm, field_name_list[f], NULL, &type, &data)); 361 idx = (const PetscInt *)data; 362 363 PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)dm), N, idx, PETSC_USE_POINTER, &is)); 364 PetscCall(PetscObjectSetName((PetscObject)is, field_name_list[f])); 365 PetscCall(private_ISView_Swarm_XDMF(is, viewer)); 366 PetscCall(ISDestroy(&is)); 367 PetscCall(DMSwarmRestoreField(dm, field_name_list[f], NULL, &type, &data)); 368 } else SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Can only write PETSC_INT and PETSC_DOUBLE"); 369 } 370 PetscCall(private_PetscViewerDestroy_XDMF(&viewer)); 371 PetscFunctionReturn(PETSC_SUCCESS); 372 } 373 374 /*@ 375 DMSwarmViewXDMF - Write `DMSWARM` fields to an XDMF3 file 376 377 Collective 378 379 Input Parameters: 380 + dm - the `DMSWARM` 381 - filename - the file name of the XDMF file (must have the extension .xmf) 382 383 Level: beginner 384 385 Note: 386 Only fields user registered with data type `PETSC_DOUBLE` or `PETSC_INT` will be written into the file 387 388 Developer Note: 389 This should be removed and replaced with the standard use of `PetscViewer` 390 391 .seealso: `DM`, `DMSWARM`, `DMSwarmViewFieldsXDMF()` 392 @*/ 393 PetscErrorCode DMSwarmViewXDMF(DM dm, const char filename[]) 394 { 395 DM_Swarm *swarm = (DM_Swarm *)dm->data; 396 Vec dvec; 397 PetscInt f; 398 PetscViewer viewer; 399 400 PetscFunctionBegin; 401 PetscCall(private_PetscViewerCreate_XDMF(PetscObjectComm((PetscObject)dm), filename, &viewer)); 402 PetscCall(private_DMSwarmView_XDMF(dm, viewer)); 403 for (f = 4; f < swarm->db->nfields; f++) { /* only examine user defined fields - the first 4 are internally created by DMSwarmPIC */ 404 DMSwarmDataField field; 405 406 /* query field type - accept all those of type PETSC_DOUBLE */ 407 field = swarm->db->field[f]; 408 if (field->petsc_type == PETSC_DOUBLE) { 409 PetscCall(DMSwarmCreateGlobalVectorFromField(dm, field->name, &dvec)); 410 PetscCall(PetscObjectSetName((PetscObject)dvec, field->name)); 411 PetscCall(private_VecView_Swarm_XDMF(dvec, viewer)); 412 PetscCall(DMSwarmDestroyGlobalVectorFromField(dm, field->name, &dvec)); 413 } else if (field->petsc_type == PETSC_INT) { 414 IS is; 415 PetscInt N; 416 const PetscInt *idx; 417 void *data; 418 419 PetscCall(DMSwarmGetLocalSize(dm, &N)); 420 PetscCall(DMSwarmGetField(dm, field->name, NULL, NULL, &data)); 421 idx = (const PetscInt *)data; 422 423 PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)dm), N, idx, PETSC_USE_POINTER, &is)); 424 PetscCall(PetscObjectSetName((PetscObject)is, field->name)); 425 PetscCall(private_ISView_Swarm_XDMF(is, viewer)); 426 PetscCall(ISDestroy(&is)); 427 PetscCall(DMSwarmRestoreField(dm, field->name, NULL, NULL, &data)); 428 } 429 } 430 PetscCall(private_PetscViewerDestroy_XDMF(&viewer)); 431 PetscFunctionReturn(PETSC_SUCCESS); 432 } 433