Skip to content

Fix integer overflows and parallelize the output writers for large meshes - #2924

Open
rois1995 wants to merge 18 commits into
developfrom
fixCGNSOutput
Open

rois1995 wants to merge 18 commits into
developfrom
fixCGNSOutput

Conversation

@rois1995

@rois1995 rois1995 commented Sep 23, 2026 •

Copy link
Copy Markdown
Contributor

Proposed Changes

This PR fixes the output of large meshes and makes the output writing faster. The problems appeared on a mesh with 273M points and 496M cells (mostly prisms) run on 576 ranks: some output files were broken, the CGNS file could not be opened in ParaView, and some formats took very long to write. Each commit is self-contained, and its message gives the details.

1. Integer overflows

Several writers stored sizes, offsets or MPI counts in a 32-bit int. Above the limits below, the run crashed or, worse, wrote a broken file without any error. In all cases, the files of meshes below the limits are unchanged.

  • MPI-IO collective writes (binary restart, PARAVIEW, PARAVIEW_LEGACY): the size in bytes written by one rank was cast to int. This broke when one rank wrote more than 2 GB in one call, e.g. a binary restart with 20 fields and more than ~13M points per rank (common in hybrid MPI+OpenMP runs). The local data is now described with an MPI datatype of 1 MB blocks, so no count can overflow.
  • Sorting of the output data: buffer sizes, offsets and MPI counts were computed as (number of fields) x (number of points) in int, which broke above ~21M points per rank with 100 output fields. The data of each point is now sent as one element of an MPI datatype, and the sizes are size_t.
  • PARAVIEW (.vtu): connectivity and offsets were written as Int32, so meshes with more than 2^31-1 connectivity entries (~358M prisms) gave unreadable files. They are now written as Int64 when needed.
  • PARAVIEW_LEGACY (.vtk): the cell array of the 3.0 format is Int32 only. When it does not fit, the VTK 5.1 layout with Int64 offsets and connectivity is used (readable by ParaView >= 5.8).
  • CGNS gather on the master rank: the MPI counts in bytes were cast to int, which broke above ~33.5M hexahedra per rank. This path is later replaced by the parallel writer (section 4).

2. CGNS: sections above 2^31-1 connectivity entries

The CGNS file of the large mesh is valid, but ParaView crashed when opening it. The VTK CGNS reader stores the connectivity size of each section in an int (reported in https://gitlab.kitware.com/vtk/vtk/-/work_items/20170). Until the reader is fixed, an element type above this limit is written as several sections with consecutive element ranges (<name>_1, <name>_2, ...). Smaller meshes still get one section per element type.

3. CGNS: named boundaries and double precision coordinates

  • Boundaries: the volume CGNS file had only the interior cells, and the surface file merged all markers in one unnamed zone. The volume file now also has one section per marker, named as the marker, with a ZoneBC and a Family that carries the BC type. ParaView shows each marker as a named patch, and SU2 can read the file back as a mesh with the same markers. SURFACE_CGNS writes one zone per marker.
  • Coordinates: they were written in single precision, which moves the points by up to ~1e-7 of the domain size, which can be larger than the smallest cells. When the output was used again as a mesh (a periodic case), the solution changed by up to 5e-5. Coordinates are now always written in double precision; with this, the same test reproduces the solution to ~1e-14.
  • New option WRT_OUTPUT_DOUBLE_PRECISION (default NO): writes the fields of CGNS, PARAVIEW and PARAVIEW_MULTIBLOCK in double precision. The default keeps single precision, so file sizes do not change.

4. Parallel writing

  • CGNS: the master rank gathered all data and wrote the whole file alone. This was slow, needed a lot of memory on the master node, and a large file could take longer than the job time left. All ranks now write their own part with the parallel CGNS API (cgp_*, parallel HDF5, already built by SU2 with MPI). Without MPI, the serial API is used as before. The write bandwidth of the CGNS files is now reported in the output summary (it always showed 0 MB/s).
  • PARAVIEW_ASCII, RESTART_ASCII, TECPLOT_ASCII and the SU2 mesh writers (ASCII and binary): the ranks wrote one after the other, with a barrier after each rank and each block of data, so the time grew with the number of ranks. Each rank now formats its part in memory, and all ranks write at the same time with one collective MPI-IO call. The files are byte-identical to the previous ones, also for multizone meshes.

Testing

The large mesh cannot be shared, and the overflow limits cannot be reached in CI (a single section above 2^31-1 entries needs more than 20 GB of memory). The tests therefore check that the parallel writers give the same files as before:

  • New regression tests (TestCases/output_writers/ascii_output.cfg, parallel regression on 2 ranks): the solution of vandv/rans/flatplate is written without iterating by SU2_CFD (RESTART_ASCII), SU2_SOL (PARAVIEW_ASCII, TECPLOT_ASCII and their surface versions) and SU2_DEF (SU2 mesh file). Each file is compared with a reference file written by the previous writers. The new files are byte-identical to the old ones, on 2 and 4 ranks. The reference files are in Reference files for the parallel ASCII output writer tests TestCases#205.
  • The parallel CGNS writer is run by the existing oneram6 parallel regression test (CGNS and SURFACE_CGNS output).
  • On the large mesh above, the parallel CGNS file is written and read correctly by post-processing scripts.

Temporary: .github/workflows/regression.yml uses the fixCGNSOutput branch of TestCases until su2code/TestCases#205 is merged. It will be set back to develop before merging.

Question for the reviewers

With MPI, the new writers always write in parallel; the old serial path (gather on the master rank, or ranks writing one after the other) is removed. Parallel HDF5 can be slow on some file systems (e.g. NFS). Should we add a config option to write serially (for example only for CGNS), or is the parallel writing enough?

Related Work

PR Checklist

  • I am submitting my contribution to the develop branch.
  • My contribution generates no new compiler warnings (try with --warnlevel=3 when using meson).
  • My contribution is commented and consistent with SU2 style (https://su2code.github.io/docs_v7/Style-Guide/).
  • I used the pre-commit hook to prevent dirty commits and used pre-commit run --all to format old commits.
  • I have added a test case that demonstrates my contribution, if necessary.
  • I have updated appropriate documentation (Tutorials, Docs Page, config_template.cpp), if necessary.

🤖 Generated with Claude Code

rois1995 and others added 14 commits September 18, 2026 11:47
The CGNS writer gathered field and connectivity data on the master rank
with a single MPI_CHAR message per rank, whose byte count was cast to int.
With 64-bit cgsize_t, a rank holding more than ~33.5M hexahedra overflowed
the count, which crashed the run or silently wrote corrupted connectivity.

Send/receive the buffers in chunks of at most 1 GiB, so every MPI count
fits in an int. MPI_CHAR is kept since typed transfers are not portable
across AD builds and platforms where long is 32-bit.

Also accumulate the element counts in CParallelDataSorter with an
unsigned long initial value, to avoid int overflow of the sum.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
WriteMPIBinaryDataAll cast the local size in bytes to int for
MPI_Type_contiguous and MPI_File_write_all. When one rank wrote more than
INT_MAX bytes in one call (e.g. a binary restart with 20 fields and more
than ~13M points per rank, typical of hybrid MPI+OpenMP runs) the write
failed or silently wrote a truncated file.

Describe the local data as 1 MiB blocks plus a remainder in a single
datatype and write it with count 1, so no int count can overflow and the
call stays collective. Output files are unchanged.

Affects the SU2 binary restart and the PARAVIEW / PARAVIEW_LEGACY writers.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The PARAVIEW (.vtu) writer stored connectivity and offsets as Int32, and
the cumulative connectivity counter of CParallelDataSorter was an int.
Meshes with more than 2^31-1 connectivity entries (e.g. more than ~268M
hexahedra or ~358M prisms) got overflowed offsets and unreadable files.

- Write connectivity and offsets as Int64 when the global connectivity
  size does not fit in Int32; smaller meshes keep Int32 (files unchanged).
- Store nElemConn_Send/nElemConn_Cum as unsigned long.
- Write NumberOfPoints/NumberOfCells without casting to int.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The VTK/ParaView CGNS reader stores the connectivity size of each section
in a 32-bit int, so it crashes on sections with more than 2^31-1 entries
(e.g. more than ~358M prisms), although the file is valid.

Write an element type that exceeds the limit as several sections with
consecutive element ranges (<name>_1, <name>_2, ...). Blocks received from
each rank are split at section boundaries. Smaller meshes still get one
section per element type, so their files are unchanged.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The PARAVIEW_LEGACY writer stored all cells in one Int32 array (the number
of nodes followed by the node ids of each cell) and printed the header
counts through int. Meshes with more than 2^31-1 entries in that array
(e.g. more than ~238M hexahedra or ~307M prisms) got broken files.

When the cell array does not fit in Int32, write the legacy 5.1 layout
(readable by VTK >= 9.0 / ParaView >= 5.8) with separate Int64 OFFSETS
and CONNECTIVITY arrays. Smaller meshes keep the 3.0 layout, so their
files are unchanged. Header counts are no longer cast to int.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
SortOutputData and PrepareSendBuffers computed buffer sizes, offsets and
MPI counts as (number of output fields) x (number of points) in int.
When one rank sent or received more than INT_MAX values (e.g. 100 output
fields and more than ~21M points per rank, as in hybrid MPI+OpenMP runs)
the counts overflowed and the output was broken.

Send the data of each point as one element of a contiguous MPI datatype,
so the MPI counts are numbers of points, and use size_t for the buffer
sizes and offsets. Output files are unchanged.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The CGNS output had no boundary information: the volume file only held
the interior cells, and the surface file merged all plotted markers into
one unnamed zone.

Volume file (CGNS):
- After the interior sections, write one boundary section per marker
  (except send-receive), named as the marker. A marker with a single
  element type gets a TRI_3/QUAD_4/BAR_2 section, otherwise MIXED.
  Faces use the volume point numbering, so they share the coordinates
  and fields of the zone.
- Each boundary element is kept on exactly one rank (no halo node and at
  least one owned node) and gathered on the master node.
- Add a ZoneBC with one BC_t per marker (FamilySpecified, PointRange of
  the section, FaceCenter/EdgeCenter) and a Family_t per marker that
  carries the physical BC type (wall, farfield, inflow, ...). The
  ParaView/VTK CGNS reader shows each marker as a named patch, and SU2
  can read the file back as a mesh with the same markers.

Surface file (SURFACE_CGNS):
- Write one zone per plotted marker, named as the marker; the surface
  data is sorted again for each marker.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
…double

The CGNS and Paraview XML writers stored coordinates and fields in single
precision. Coordinates in single precision move the points by up to ~1e-7
of the size of the domain, which can exceed the size of the smallest
cells, so a CGNS output file could not be reused as a mesh without losing
accuracy (in a round trip of a periodic case the solution changed by up
to 5e-5).

- Write the CGNS coordinates always in double precision.
- Add WRT_OUTPUT_DOUBLE_PRECISION (default NO) to write the fields of the
  volume and surface files in double precision, for CGNS, PARAVIEW and
  PARAVIEW_MULTIBLOCK. Single precision remains the default, so files and
  their sizes are unchanged for visualization use.

With double coordinates the same round trip reproduces the coordinates
exactly and the solution to ~1e-14.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The CGNS writer gathered all data on the master node, which wrote the
whole file by itself: the master was an I/O and memory bottleneck, and a
large file could take long enough for the job to end before the file was
closed.

Open the file with cgp_open on all ranks and write each rank's own range
of points, elements and fields with cgp_coord_write_data,
cgp_elements_write_data and cgp_field_write_data. The nodes of the file
(base, zone, sections, solution, fields, boundary conditions, families)
are metadata and are created by all ranks with the same arguments. Ranks
without data take part in the collective calls and write nothing.

Boundary sections are written the same way, so the gather of the boundary
faces on the master node is also gone; markers with more than one element
type use the parallel MIXED API (cgp_poly_*). Without MPI the serial CGNS
API is used, as before. The chunked MPI transfers of the old gather path
are no longer needed and were removed.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
PARAVIEW_ASCII, RESTART_ASCII and TECPLOT_ASCII wrote their files by
letting the ranks take turns: each rank opened the file, appended its
part and closed it, with a barrier after every rank and for every block
of data (points, connectivity, cell types, each field). Only one rank
wrote at any moment, and the number of barriers grows with the number of
ranks, which makes these formats very slow at scale.

Each rank now formats its own data into a string, and all ranks write
their strings at the same time with the collective MPI-IO call of
CFileWriter, at the offset that follows the text of the ranks before it
(new helper CFileWriter::WriteMPIStringAll). The headers are still
written by the master node only.

The files are byte-identical to those of the previous implementation.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
CSU2MeshFileWriter and CSU2MeshBinaryFileWriter also let the ranks take
turns: each rank opened the file, appended its elements or points and
closed it, with a collective call after every rank.

Each rank now builds its own part of the file in memory and all ranks
write it at the same time with the collective MPI-IO call of CFileWriter,
as the other output files do. The global element and point indices come
from the number of elements and points of the ranks before, gathered
once instead of once per rank.

Multizone meshes are written into a single file, so CFileWriter::
OpenMPIFile takes an append flag that continues at the end of an existing
file (and creates it when it does not exist, e.g. for the first zone of
a file of its own).

The files are byte-identical to those of the previous implementation,
also for multizone cases.

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
The CGNS writer never set the bandwidth, so the file writing summary
always printed 0 MB/s for CGNS and SURFACE_CGNS. Time the writing and
read the file size after closing it, as the Tecplot binary writer does.

Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6
Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6
Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6
@rois1995
rois1995 marked this pull request as ready for review September 23, 2026 10:27
The ASCII writers (RESTART_ASCII, PARAVIEW_ASCII, TECPLOT_ASCII and the
SU2 mesh file) now write in parallel. Write the solution of the
vandv/rans/flatplate case on 2 ranks without iterating, with SU2_CFD,
SU2_SOL and SU2_DEF, and compare the files with reference files written
by the previous writers. The reference files are in su2code/TestCases,
folder output_writers.

Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6
Temporary, for the reference files of su2code/TestCases#205. To be
reverted to develop before merging.

Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6
Comment thread SU2_CFD/src/output/filewriter/CTecplotFileWriter.cpp Fixed
Comment thread SU2_CFD/src/output/filewriter/CTecplotFileWriter.cpp Fixed
Use unsigned long for the field counter of the Tecplot ASCII writer, which
is compared with the size of the field names, and rename the offset of the
boundary elements of the SU2 binary mesh writer, which hid the one of the
volume elements.

Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6
Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6
@bigfooted

Copy link
Copy Markdown
Contributor

wasn't this an issue that was raised at a conference a looong time ago? :-)
@rois1995 can we discuss on slack the status and order of importance/merging of all these huge PRs?

This branch has not been deployed

No deployments
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants