Fix XDMF output to explicitly write cell-to-vertex topology and ensure strict ParaView compatibility#205
Open
tiannh7 wants to merge 5 commits into
Open
Conversation
Contributor
There was a problem hiding this comment.
Pull request overview
Note
Copilot was unable to run its full agentic suite in this review.
Adds ParaView-compatible “/viz” HDF5 geometry/topology outputs for meshes written via write_timestep, and validates that the generated XDMF references these groups.
Changes:
- Write
/viz/geometry/verticesand/viz/topology/cellsinto the mesh HDF5 whencreate_xdmfis enabled. - Add MPI gather-based assembly of global vertex lists and dense connectivity for visualization.
- Add a regression test ensuring
/vizdatasets exist and that the XDMF points to them.
Reviewed changes
Copilot reviewed 2 out of 2 changed files in this pull request and generated 9 comments.
| File | Description |
|---|---|
| tests/test_0005_xdmf_compat.py | Adds a regression test validating /viz datasets and XDMF references. |
| src/underworld3/discretisation/discretisation_mesh.py | Implements /viz mesh group writer, hooks it into write_timestep, and adds extra XDMF warnings/validation. |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
| f"to DMPlex vertex count ({n_local_vertices}) for {mesh_h5_path}." | ||
| ) | ||
|
|
||
| local_vertices = numpy.column_stack((vertex_gids, coords_local)) |
Comment on lines
+4120
to
+4131
| pStart, _ = dm.getDepthStratum(0) | ||
| cStart, cEnd = dm.getHeightStratum(0) | ||
| vertex_numbering = dm.getVertexNumbering().getIndices() | ||
| vertex_gids = _petsc_numbering_to_global_ids(vertex_numbering) | ||
| cell_num_points = mesh.element.entities[mesh.dim] | ||
|
|
||
| cell_points_list = [] | ||
| for cell_id in range(cStart, cEnd): | ||
| points = numpy.asarray( | ||
| dm.getTransitiveClosure(cell_id)[0][-cell_num_points:], | ||
| dtype=numpy.int64, | ||
| ) |
Comment on lines
+4169
to
+4176
| gathered_vertices = uw.mpi.comm.gather(local_vertices, root=0) | ||
| gathered_cells = uw.mpi.comm.gather(local_cells, root=0) | ||
| uw.mpi.barrier() | ||
|
|
||
| if uw.mpi.rank == 0: | ||
| import h5py | ||
|
|
||
| vertex_blocks = [block for block in gathered_vertices if block is not None and block.size > 0] |
| cell_blocks = [block for block in gathered_cells if block is not None and block.size > 0] | ||
|
|
||
| if vertex_blocks: | ||
| all_vertices = numpy.vstack(vertex_blocks) |
Comment on lines
+4189
to
+4206
| vertex_lookup = {} | ||
| for row in all_vertices: | ||
| gid = int(row[0]) | ||
| if gid not in vertex_lookup: | ||
| vertex_lookup[gid] = numpy.asarray(row[1:], dtype=numpy.float64) | ||
|
|
||
| ordered_gids = numpy.array(sorted(vertex_lookup.keys()), dtype=numpy.int64) | ||
| if ordered_gids.size: | ||
| ordered_vertices = numpy.vstack([vertex_lookup[int(gid)] for gid in ordered_gids]) | ||
| else: | ||
| ordered_vertices = numpy.empty((0, mesh.dim), dtype=numpy.float64) | ||
|
|
||
| dense_lookup = {int(gid): idx for idx, gid in enumerate(ordered_gids.tolist())} | ||
| if all_cells.size: | ||
| dense_cells = numpy.asarray( | ||
| [[dense_lookup[int(gid)] for gid in row] for row in all_cells], | ||
| dtype=numpy.int64, | ||
| ) |
Comment on lines
+4189
to
+4206
| vertex_lookup = {} | ||
| for row in all_vertices: | ||
| gid = int(row[0]) | ||
| if gid not in vertex_lookup: | ||
| vertex_lookup[gid] = numpy.asarray(row[1:], dtype=numpy.float64) | ||
|
|
||
| ordered_gids = numpy.array(sorted(vertex_lookup.keys()), dtype=numpy.int64) | ||
| if ordered_gids.size: | ||
| ordered_vertices = numpy.vstack([vertex_lookup[int(gid)] for gid in ordered_gids]) | ||
| else: | ||
| ordered_vertices = numpy.empty((0, mesh.dim), dtype=numpy.float64) | ||
|
|
||
| dense_lookup = {int(gid): idx for idx, gid in enumerate(ordered_gids.tolist())} | ||
| if all_cells.size: | ||
| dense_cells = numpy.asarray( | ||
| [[dense_lookup[int(gid)] for gid in row] for row in all_cells], | ||
| dtype=numpy.int64, | ||
| ) |
Comment on lines
+4326
to
+4334
| cells_data = cells[...] | ||
| c_min, c_max = cells_data.min(), cells_data.max() | ||
| if c_min < 0 or c_max >= numVertices: | ||
| warnings.warn( | ||
| f"XDMF connectivity is invalid! cells max {c_max} >= " | ||
| f"numVertices {numVertices} or min {c_min} < 0. ParaView will likely crash. " | ||
| f"Ensure cell-to-vertex connectivity is written." | ||
| ) | ||
|
|
Comment on lines
+4321
to
+4324
| warnings.warn( | ||
| "Using raw '/topology/cells' for XDMF. This may not be Paraview-compatible. " | ||
| "Expected '/viz/topology/cells'." | ||
| ) |
Comment on lines
+4137
to
+4142
| if mesh.dim == 3: | ||
| if dm.isSimplex(): | ||
| reorder = [0, 2, 1, 3] | ||
| else: | ||
| reorder = [0, 3, 2, 1, 4, 5, 6, 7] | ||
| cell_points_list = [pts[reorder] for pts in cell_points_list] |
Member
|
@lmoresi I'm also experiencing ParaView crashes when viewing mesh xdmf output, which this PR fixes. I've tested and confirmed working. commit 22dd9ab on feature/annulus-spherical-benchmark-debug takes a different approach (adds |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
This PR resolves a critical issue where generating XDMF files could silently output PETSc DMPlex's internal point numbering instead of zero-indexed vertex connectivity. This malformed connectivity causes abrupt crashes in ParaView 5.13+ and renders corrupted meshes in older versions.
Key changes:
/vizgeneration: Ported the_write_mesh_viz_groupsgeneration logic intowrite_timestep. It now correctly gathers and constructs true/viz/geometry/verticesand/viz/topology/cellsacross MPI ranks. This restores strict XDMF schema compliance and is 100% backward-compatible with all older versions of ParaView.checkpoint_xdmfto explicitly warn users if the connectivity goes out of bounds (cells.max() >= numVertices), ifnumCornersmismatches the topology type, or if Field attribute lengths (Node/Cellcenter) are malformed./vizgroups are consistently generated and XDMF reference paths correctly point to them without triggering validation warnings.While this PR hardens our existing XDMF/HDF5 pipeline and ensures stability, it might be worth exploring alternative or complementary formats (such as native
.vtu/.vtk) down the road. Integrating tools likepyvistaormeshiocould potentially offer a more integrated, single-file visualization workflow. This is just a thought for the future.