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 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 if (container) { 101 PetscCall(PetscContainerGetPointer(container, (void **)&bytes)); 102 } else SETERRQ(PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Valid to find attached data XDMFViewerContext"); 103 104 PetscCall(PetscObjectTypeCompare((PetscObject)dm, DMSWARM, &isswarm)); 105 PetscCheck(isswarm, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Only valid for DMSwarm"); 106 107 PetscCall(PetscObjectCompose((PetscObject)viewer, "DMSwarm", (PetscObject)dm)); 108 109 PetscCall(PetscViewerASCIIPushTab(viewer)); 110 PetscCall(PetscObjectGetName((PetscObject)dm, &dmname)); 111 if (!dmname) PetscCall(DMGetOptionsPrefix(dm, &dmname)); 112 if (!dmname) { 113 PetscCall(PetscViewerASCIIPrintf(viewer, "<Grid Name=\"DMSwarm\" GridType=\"Uniform\">\n")); 114 } else { 115 PetscCall(PetscViewerASCIIPrintf(viewer, "<Grid Name=\"DMSwarm[%s]\" GridType=\"Uniform\">\n", dmname)); 116 } 117 118 /* create a sub-viewer for topology, geometry and all data fields */ 119 /* name is viewer.name + "_swarm_fields.pbin" */ 120 PetscCall(PetscViewerCreate(PetscObjectComm((PetscObject)viewer), &fviewer)); 121 PetscCall(PetscViewerSetType(fviewer, PETSCVIEWERBINARY)); 122 PetscCall(PetscViewerBinarySetSkipHeader(fviewer, PETSC_TRUE)); 123 PetscCall(PetscViewerBinarySetSkipInfo(fviewer, PETSC_TRUE)); 124 PetscCall(PetscViewerFileSetMode(fviewer, FILE_MODE_WRITE)); 125 126 PetscCall(PetscViewerFileGetName(viewer, &viewername)); 127 PetscCall(private_CreateDataFileNameXDMF(viewername, datafile)); 128 PetscCall(PetscViewerFileSetName(fviewer, datafile)); 129 PetscCall(PetscStrrchr(datafile, '/', &datafilename)); 130 131 PetscCall(DMSwarmGetSize(dm, &ng)); 132 133 /* write topology header */ 134 PetscCall(PetscViewerASCIIPushTab(viewer)); 135 PetscCall(PetscViewerASCIIPrintf(viewer, "<Topology Dimensions=\"%" PetscInt_FMT "\" TopologyType=\"Mixed\">\n", ng)); 136 PetscCall(PetscViewerASCIIPushTab(viewer)); 137 PetscCall(PetscViewerASCIIPrintf(viewer, "<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Int\" Dimensions=\"%" PetscInt_FMT "\" Seek=\"%ld\">\n", ng * 3, bytes[0])); 138 PetscCall(PetscViewerASCIIPushTab(viewer)); 139 PetscCall(PetscViewerASCIIPrintf(viewer, "%s\n", datafilename)); 140 PetscCall(PetscViewerASCIIPopTab(viewer)); 141 PetscCall(PetscViewerASCIIPrintf(viewer, "</DataItem>\n")); 142 PetscCall(PetscViewerASCIIPopTab(viewer)); 143 PetscCall(PetscViewerASCIIPrintf(viewer, "</Topology>\n")); 144 PetscCall(PetscViewerASCIIPopTab(viewer)); 145 146 /* write topology data */ 147 for (k = 0; k < ng; k++) { 148 PetscInt pvertex[3]; 149 150 pvertex[0] = 1; 151 pvertex[1] = 1; 152 pvertex[2] = k; 153 PetscCall(PetscViewerBinaryWrite(fviewer, pvertex, 3, PETSC_INT)); 154 } 155 bytes[0] += sizeof(PetscInt) * ng * 3; 156 157 /* write geometry header */ 158 PetscCall(PetscViewerASCIIPushTab(viewer)); 159 PetscCall(DMGetDimension(dm, &dim)); 160 switch (dim) { 161 case 1: 162 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No support for 1D"); 163 case 2: 164 PetscCall(PetscViewerASCIIPrintf(viewer, "<Geometry Type=\"XY\">\n")); 165 break; 166 case 3: 167 PetscCall(PetscViewerASCIIPrintf(viewer, "<Geometry Type=\"XYZ\">\n")); 168 break; 169 } 170 PetscCall(PetscViewerASCIIPushTab(viewer)); 171 PetscCall(PetscViewerASCIIPrintf(viewer, "<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Float\" Precision=\"8\" Dimensions=\"%" PetscInt_FMT " %" PetscInt_FMT "\" Seek=\"%ld\">\n", ng, dim, bytes[0])); 172 PetscCall(PetscViewerASCIIPushTab(viewer)); 173 PetscCall(PetscViewerASCIIPrintf(viewer, "%s\n", datafilename)); 174 PetscCall(PetscViewerASCIIPopTab(viewer)); 175 PetscCall(PetscViewerASCIIPrintf(viewer, "</DataItem>\n")); 176 PetscCall(PetscViewerASCIIPopTab(viewer)); 177 PetscCall(PetscViewerASCIIPrintf(viewer, "</Geometry>\n")); 178 PetscCall(PetscViewerASCIIPopTab(viewer)); 179 180 PetscCall(DMSwarmGetCellDMActive(dm, &celldm)); 181 PetscCall(DMSwarmCellDMGetCoordinateFields(celldm, &Nfc, &coordFields)); 182 PetscCheck(Nfc == 1, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "We only support a single coordinate field right now, not %" PetscInt_FMT, Nfc); 183 184 /* write geometry data */ 185 PetscCall(DMSwarmCreateGlobalVectorFromField(dm, coordFields[0], &dvec)); 186 PetscCall(VecView(dvec, fviewer)); 187 PetscCall(DMSwarmDestroyGlobalVectorFromField(dm, coordFields[0], &dvec)); 188 bytes[0] += sizeof(PetscReal) * ng * dim; 189 190 PetscCall(PetscViewerDestroy(&fviewer)); 191 PetscFunctionReturn(PETSC_SUCCESS); 192 } 193 194 static PetscErrorCode private_VecView_Swarm_XDMF(Vec x, PetscViewer viewer) 195 { 196 long int *bytes = NULL; 197 PetscContainer container = NULL; 198 const char *viewername; 199 char datafile[PETSC_MAX_PATH_LEN]; 200 char *datafilename; 201 PetscViewer fviewer; 202 PetscInt N, bs; 203 const char *vecname; 204 char fieldname[PETSC_MAX_PATH_LEN]; 205 206 PetscFunctionBegin; 207 PetscCall(PetscObjectQuery((PetscObject)viewer, "XDMFViewerContext", (PetscObject *)&container)); 208 PetscCheck(container, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Unable to find attached data XDMFViewerContext"); 209 PetscCall(PetscContainerGetPointer(container, (void **)&bytes)); 210 PetscCall(PetscViewerFileGetName(viewer, &viewername)); 211 PetscCall(private_CreateDataFileNameXDMF(viewername, datafile)); 212 213 /* re-open a sub-viewer for all data fields */ 214 /* name is viewer.name + "_swarm_fields.pbin" */ 215 PetscCall(PetscViewerCreate(PetscObjectComm((PetscObject)viewer), &fviewer)); 216 PetscCall(PetscViewerSetType(fviewer, PETSCVIEWERBINARY)); 217 PetscCall(PetscViewerBinarySetSkipHeader(fviewer, PETSC_TRUE)); 218 PetscCall(PetscViewerBinarySetSkipInfo(fviewer, PETSC_TRUE)); 219 PetscCall(PetscViewerFileSetMode(fviewer, FILE_MODE_APPEND)); 220 PetscCall(PetscViewerFileSetName(fviewer, datafile)); 221 PetscCall(PetscStrrchr(datafile, '/', &datafilename)); 222 223 PetscCall(VecGetSize(x, &N)); 224 PetscCall(VecGetBlockSize(x, &bs)); 225 N = N / bs; 226 PetscCall(PetscObjectGetName((PetscObject)x, &vecname)); 227 if (!vecname) { 228 PetscCall(PetscSNPrintf(fieldname, PETSC_MAX_PATH_LEN - 1, "swarmfield_%d", ((PetscObject)x)->tag)); 229 } else { 230 PetscCall(PetscSNPrintf(fieldname, PETSC_MAX_PATH_LEN - 1, "%s", vecname)); 231 } 232 233 /* write data header */ 234 PetscCall(PetscViewerASCIIPushTab(viewer)); 235 PetscCall(PetscViewerASCIIPrintf(viewer, "<Attribute Center=\"Node\" Name=\"%s\" Type=\"None\">\n", fieldname)); 236 PetscCall(PetscViewerASCIIPushTab(viewer)); 237 if (bs == 1) { 238 PetscCall(PetscViewerASCIIPrintf(viewer, "<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Float\" Precision=\"8\" Dimensions=\"%" PetscInt_FMT "\" Seek=\"%ld\">\n", N, bytes[0])); 239 } else { 240 PetscCall(PetscViewerASCIIPrintf(viewer, "<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Float\" Precision=\"8\" Dimensions=\"%" PetscInt_FMT " %" PetscInt_FMT "\" Seek=\"%ld\">\n", N, bs, bytes[0])); 241 } 242 PetscCall(PetscViewerASCIIPushTab(viewer)); 243 PetscCall(PetscViewerASCIIPrintf(viewer, "%s\n", datafilename)); 244 PetscCall(PetscViewerASCIIPopTab(viewer)); 245 PetscCall(PetscViewerASCIIPrintf(viewer, "</DataItem>\n")); 246 PetscCall(PetscViewerASCIIPopTab(viewer)); 247 PetscCall(PetscViewerASCIIPrintf(viewer, "</Attribute>\n")); 248 PetscCall(PetscViewerASCIIPopTab(viewer)); 249 250 /* write data */ 251 PetscCall(VecView(x, fviewer)); 252 bytes[0] += sizeof(PetscReal) * N * bs; 253 254 PetscCall(PetscViewerDestroy(&fviewer)); 255 PetscFunctionReturn(PETSC_SUCCESS); 256 } 257 258 static PetscErrorCode private_ISView_Swarm_XDMF(IS is, PetscViewer viewer) 259 { 260 long int *bytes = NULL; 261 PetscContainer container = NULL; 262 const char *viewername; 263 char datafile[PETSC_MAX_PATH_LEN]; 264 char *datafilename; 265 PetscViewer fviewer; 266 PetscInt N, bs; 267 const char *vecname; 268 char fieldname[PETSC_MAX_PATH_LEN]; 269 270 PetscFunctionBegin; 271 PetscCall(PetscObjectQuery((PetscObject)viewer, "XDMFViewerContext", (PetscObject *)&container)); 272 PetscCheck(container, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Unable to find attached data XDMFViewerContext"); 273 PetscCall(PetscContainerGetPointer(container, (void **)&bytes)); 274 PetscCall(PetscViewerFileGetName(viewer, &viewername)); 275 PetscCall(private_CreateDataFileNameXDMF(viewername, datafile)); 276 277 /* re-open a sub-viewer for all data fields */ 278 /* name is viewer.name + "_swarm_fields.pbin" */ 279 PetscCall(PetscViewerCreate(PetscObjectComm((PetscObject)viewer), &fviewer)); 280 PetscCall(PetscViewerSetType(fviewer, PETSCVIEWERBINARY)); 281 PetscCall(PetscViewerBinarySetSkipHeader(fviewer, PETSC_TRUE)); 282 PetscCall(PetscViewerBinarySetSkipInfo(fviewer, PETSC_TRUE)); 283 PetscCall(PetscViewerFileSetMode(fviewer, FILE_MODE_APPEND)); 284 PetscCall(PetscViewerFileSetName(fviewer, datafile)); 285 PetscCall(PetscStrrchr(datafile, '/', &datafilename)); 286 287 PetscCall(ISGetSize(is, &N)); 288 PetscCall(ISGetBlockSize(is, &bs)); 289 N = N / bs; 290 PetscCall(PetscObjectGetName((PetscObject)is, &vecname)); 291 if (!vecname) { 292 PetscCall(PetscSNPrintf(fieldname, PETSC_MAX_PATH_LEN - 1, "swarmfield_%d", ((PetscObject)is)->tag)); 293 } else { 294 PetscCall(PetscSNPrintf(fieldname, PETSC_MAX_PATH_LEN - 1, "%s", vecname)); 295 } 296 297 /* write data header */ 298 PetscCall(PetscViewerASCIIPushTab(viewer)); 299 PetscCall(PetscViewerASCIIPrintf(viewer, "<Attribute Center=\"Node\" Name=\"%s\" Type=\"None\">\n", fieldname)); 300 PetscCall(PetscViewerASCIIPushTab(viewer)); 301 if (bs == 1) { 302 PetscCall(PetscViewerASCIIPrintf(viewer, "<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Int\" Precision=\"4\" Dimensions=\"%" PetscInt_FMT "\" Seek=\"%ld\">\n", N, bytes[0])); 303 } else { 304 PetscCall(PetscViewerASCIIPrintf(viewer, "<DataItem Format=\"Binary\" Endian=\"Big\" DataType=\"Int\" Precision=\"4\" Dimensions=\"%" PetscInt_FMT " %" PetscInt_FMT "\" Seek=\"%ld\">\n", N, bs, bytes[0])); 305 } 306 PetscCall(PetscViewerASCIIPushTab(viewer)); 307 PetscCall(PetscViewerASCIIPrintf(viewer, "%s\n", datafilename)); 308 PetscCall(PetscViewerASCIIPopTab(viewer)); 309 PetscCall(PetscViewerASCIIPrintf(viewer, "</DataItem>\n")); 310 PetscCall(PetscViewerASCIIPopTab(viewer)); 311 PetscCall(PetscViewerASCIIPrintf(viewer, "</Attribute>\n")); 312 PetscCall(PetscViewerASCIIPopTab(viewer)); 313 314 /* write data */ 315 PetscCall(ISView(is, fviewer)); 316 bytes[0] += sizeof(PetscInt) * N * bs; 317 318 PetscCall(PetscViewerDestroy(&fviewer)); 319 PetscFunctionReturn(PETSC_SUCCESS); 320 } 321 322 /*@C 323 DMSwarmViewFieldsXDMF - Write a selection of DMSwarm fields to an XDMF3 file 324 325 Collective 326 327 Input Parameters: 328 + dm - the `DMSWARM` 329 . filename - the file name of the XDMF file (must have the extension .xmf) 330 . nfields - the number of fields to write into the XDMF file 331 - field_name_list - array of length nfields containing the textual name of fields to write 332 333 Level: beginner 334 335 Note: 336 Only fields registered with data type `PETSC_DOUBLE` or `PETSC_INT` can be written into the file 337 338 .seealso: `DM`, `DMSWARM`, `DMSwarmViewXDMF()` 339 @*/ 340 PetscErrorCode DMSwarmViewFieldsXDMF(DM dm, const char filename[], PetscInt nfields, const char *field_name_list[]) 341 { 342 Vec dvec; 343 PetscInt f, N; 344 PetscViewer viewer; 345 346 PetscFunctionBegin; 347 PetscCall(private_PetscViewerCreate_XDMF(PetscObjectComm((PetscObject)dm), filename, &viewer)); 348 PetscCall(private_DMSwarmView_XDMF(dm, viewer)); 349 PetscCall(DMSwarmGetLocalSize(dm, &N)); 350 for (f = 0; f < nfields; f++) { 351 void *data; 352 PetscDataType type; 353 354 PetscCall(DMSwarmGetField(dm, field_name_list[f], NULL, &type, &data)); 355 PetscCall(DMSwarmRestoreField(dm, field_name_list[f], NULL, &type, &data)); 356 if (type == PETSC_DOUBLE) { 357 PetscCall(DMSwarmCreateGlobalVectorFromField(dm, field_name_list[f], &dvec)); 358 PetscCall(PetscObjectSetName((PetscObject)dvec, field_name_list[f])); 359 PetscCall(private_VecView_Swarm_XDMF(dvec, viewer)); 360 PetscCall(DMSwarmDestroyGlobalVectorFromField(dm, field_name_list[f], &dvec)); 361 } else if (type == PETSC_INT) { 362 IS is; 363 const PetscInt *idx; 364 365 PetscCall(DMSwarmGetField(dm, field_name_list[f], NULL, &type, &data)); 366 idx = (const PetscInt *)data; 367 368 PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)dm), N, idx, PETSC_USE_POINTER, &is)); 369 PetscCall(PetscObjectSetName((PetscObject)is, field_name_list[f])); 370 PetscCall(private_ISView_Swarm_XDMF(is, viewer)); 371 PetscCall(ISDestroy(&is)); 372 PetscCall(DMSwarmRestoreField(dm, field_name_list[f], NULL, &type, &data)); 373 } else SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Can only write PETSC_INT and PETSC_DOUBLE"); 374 } 375 PetscCall(private_PetscViewerDestroy_XDMF(&viewer)); 376 PetscFunctionReturn(PETSC_SUCCESS); 377 } 378 379 /*@ 380 DMSwarmViewXDMF - Write `DMSWARM` fields to an XDMF3 file 381 382 Collective 383 384 Input Parameters: 385 + dm - the `DMSWARM` 386 - filename - the file name of the XDMF file (must have the extension .xmf) 387 388 Level: beginner 389 390 Note: 391 Only fields user registered with data type `PETSC_DOUBLE` or `PETSC_INT` will be written into the file 392 393 Developer Note: 394 This should be removed and replaced with the standard use of `PetscViewer` 395 396 .seealso: `DM`, `DMSWARM`, `DMSwarmViewFieldsXDMF()` 397 @*/ 398 PetscErrorCode DMSwarmViewXDMF(DM dm, const char filename[]) 399 { 400 DM_Swarm *swarm = (DM_Swarm *)dm->data; 401 Vec dvec; 402 PetscInt f; 403 PetscViewer viewer; 404 405 PetscFunctionBegin; 406 PetscCall(private_PetscViewerCreate_XDMF(PetscObjectComm((PetscObject)dm), filename, &viewer)); 407 PetscCall(private_DMSwarmView_XDMF(dm, viewer)); 408 for (f = 4; f < swarm->db->nfields; f++) { /* only examine user defined fields - the first 4 are internally created by DMSwarmPIC */ 409 DMSwarmDataField field; 410 411 /* query field type - accept all those of type PETSC_DOUBLE */ 412 field = swarm->db->field[f]; 413 if (field->petsc_type == PETSC_DOUBLE) { 414 PetscCall(DMSwarmCreateGlobalVectorFromField(dm, field->name, &dvec)); 415 PetscCall(PetscObjectSetName((PetscObject)dvec, field->name)); 416 PetscCall(private_VecView_Swarm_XDMF(dvec, viewer)); 417 PetscCall(DMSwarmDestroyGlobalVectorFromField(dm, field->name, &dvec)); 418 } else if (field->petsc_type == PETSC_INT) { 419 IS is; 420 PetscInt N; 421 const PetscInt *idx; 422 void *data; 423 424 PetscCall(DMSwarmGetLocalSize(dm, &N)); 425 PetscCall(DMSwarmGetField(dm, field->name, NULL, NULL, &data)); 426 idx = (const PetscInt *)data; 427 428 PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)dm), N, idx, PETSC_USE_POINTER, &is)); 429 PetscCall(PetscObjectSetName((PetscObject)is, field->name)); 430 PetscCall(private_ISView_Swarm_XDMF(is, viewer)); 431 PetscCall(ISDestroy(&is)); 432 PetscCall(DMSwarmRestoreField(dm, field->name, NULL, NULL, &data)); 433 } 434 } 435 PetscCall(private_PetscViewerDestroy_XDMF(&viewer)); 436 PetscFunctionReturn(PETSC_SUCCESS); 437 } 438