Faster and parallel reading; velocity-aware GMSH preprocessing - #95
Merged
Merged
Conversation
The new CTest suite runs pumgen with 1 to 4 ranks on small gmsh-generated fixtures and compares the written datasets with h5diff against reference PUML files. fixtures/generate.py computes the references from the gmsh model directly, i.e. independently of PUMGen. The element and vertex counts of the layered fixture do not divide evenly among 3 or 4 ranks. CI now runs the tests for the Release and the Debug (AddressSanitizer) build. AI-generated. Model: Claude Opus 5.5
Vertices and elements were distributed in ceil-sized chunks, while SerialMeshFile allocates the chunks given by getChunksize. Whenever the counts do not divide evenly (e.g. 815206 elements on 3 ranks), ranks received more data than their buffers hold; pumgen then aborted or hung. Group numbers and boundary conditions were routed to their elements with the same mismatch, and readBoundaries only read the first nBoundaries/P faces in many cases, so with 2 ranks half of the boundary conditions were silently dropped. All four readers now use getChunksize/getChunksum, a new getChunkOwner maps elements to ranks, and groups and boundaries are read completely and ended by an explicit terminator message. The serial reader no longer reads past a group section that ends exactly at a chunk boundary or past the last section. AI-generated. Model: Claude Opus 5.5
--filter-enable failed for every mesh: the chunk shape of the one-dimensional datasets (group, boundary) was set with rank 2 and a zero second dimension, which HDF5 rejects. The chunk rows were also derived from the local size, although the collective dataset creation needs the same shape on all ranks. The chunk shape now uses the dataset rank and the global size. With more ranks than cells, the empty ranks dereferenced the end iterator of their empty insphere list. They now contribute +inf to the minimum and take part in the collective writes with an empty selection. AI-generated. Model: Claude Opus 5.5
The write chunk size was read from --filter-chunksize, so --chunksize had no effect and e.g. --filter-chunksize 4096 shrank the write buffer to 4 KiB. --order took an optional argument, so "-o 2" left the order empty (read as 0, reported as "Order too high: 0") and consumed the input file name as a positional argument; the order is now a required argument and has to be at least 1. --compactify-datatypes is a flag, but required an argument and thus swallowed the next command line argument. AI-generated. Model: Claude Opus 5.5
Writing with compact types aborted at the geometry: the datatype was closed whenever the option was set, including the immutable predefined float type. Besides, the reduced integer size was derived from the number of cells rather than from the stored values. The bit-packed boundary codes then did not fit and HDF5 clipped them (306 of 1597 cells in the layered test mesh), and vertex indices could exceed the reduced range as well. The size is now the smallest one that holds the global value range of the dataset, including a sign bit for signed types; data with negative values is not reduced. Only copied types are closed. AI-generated. Model: Claude Opus 5.5
Parse errors were only appended to a message and parsing continued with the value 0, e.g. a coordinate "abc" was silently read as 0 and pumgen succeeded; an erroneous node tag became 0 and wrapped around when shifted to a 0-based index. Errors now throw, and parseFile returns false with the first error message. Further messages and checks: - binary MSH files are rejected with a hint instead of "Expected 0" - MSH 4.0 is rejected, since the parser implements the 4.1 layout - unknown element types are reported instead of indexing past the node count table (which aborted with std::out_of_range) - missing $Entities/$Nodes/$Elements sections are reported The node buffer is sized by the largest element (up to 1000 nodes) instead of the number of element types (136), which overflowed the stack in the msh2 parser, and msh2 elements without tags get tag 0 instead of the tag of the previous element. AI-generated. Model: Claude Opus 5.5
Only rank 0 parses the GMSH file, but every rank asked its own (empty) builder whether the mesh has a vertex identification. Rank 0 therefore scattered and wrote the identification alone, and pumgen hung in the collective HDF5 write. The flag is now broadcast. Adds a periodic fixture (a box, periodic in x) whose reference identification is computed from gmsh's periodic node pairs; the tests compare the identify dataset as well. AI-generated. Model: Claude Opus 5.5
The msh4 parser used tag - 1 as the vertex index and only checked it with an assert, so any file whose tags do not form 1..N wrote out of bounds, e.g. dense tags starting at 1001. Element and periodic node references were not checked in either parser, and duplicate tags left vertices undefined. Now the node tags have to be contiguous; in msh4 files they may start at any positive tag. Duplicate tags and references to unknown nodes are reported, and $Nodes has to precede $Elements and $Periodic. Gmsh writes contiguous tags also for models with volumes outside of any physical group, so non-contiguous tags are rejected instead of remapped. Parametric node coordinates (Mesh.SaveParametric = 1) are skipped instead of being read as the next node, and facet node indices use std::size_t like the cell node indices. AI-generated. Model: Claude Opus 5.5
The msh4 parser picked the physical tag of an element from the surface map only for element type 2 (linear triangles) and from the volume map for all other types. The boundary conditions of higher-order meshes were thus looked up with surface entity tags in the volume map (208 of 320 cells wrong in the order-2 test mesh). The map is now chosen by the entity dimension of the element block, and lookups no longer insert entries. Elements of surfaces or volumes without a physical group still get 0, but pumgen now warns about them. Partitioned MSH files, whose physical groups live in $PartitionedEntities, are rejected instead of being read with all groups and boundary conditions 0. The minimum insphere assumed 4 nodes per cell and reported 0 for order 2; it now uses the four vertices of each cell. AI-generated. Model: Claude Opus 5.5
The lexer read the file character by character through std::istream::get and converted numbers with strtol/strtod; on a mesh with 3.6M cells (170 MB MSH 4.1, 196 MB MSH 2.2) parsing took 2.18 s and 2.82 s. The new MshInput reads the file in 16 MiB blocks and parses the numbers in place with std::from_chars (with a strtod fallback for standard libraries without floating-point from_chars), which brings parsing down to 0.41 s and 0.54 s (5.3x); the results are identical. The parsers become part of src/meshreader and read the values they expect instead of generic tokens. Unneeded sections are skipped as a whole. Errors are still reported with line and column; binary MSH 2 files get a message of their own. A new unit test parses the fixtures with a buffer of 256 bytes, so that tokens and section markers straddle the buffer boundaries. AI-generated. Model: Claude Opus 5.5
Binary files are parsed as the ASCII ones, but with the values read as raw 4-byte ints, doubles and data-size integers; node tags, coordinates and element nodes are read in bulk. Other byte orders (detected by the marker in $MeshFormat) and a data size of 4 are supported as well. On the mesh with 3.6M cells, parsing the binary file takes 0.18 s instead of 0.43 s for the ASCII file (and 2.18 s for the ASCII file before the buffered reader). Errors in binary data are reported with their byte offset. The binary fixtures are written by gmsh; fixtures/generate.py derives big-endian and 4-byte-size variants and first checks that its converter reproduces gmsh's file byte for byte. The references of the binary tests are compared with a tolerance, since they are built from the ASCII files (%.16g). AI-generated. Model: Claude Opus 5.5
In one of about 90 runs, a test expecting an error did not see pumgen's error message: mpiexec may drop the output of a process which ends with MPI_Abort. The tests expecting an error on one rank now start pumgen as MPI singleton, which also shortens them from about 2.3 s to 0.35 s each. AI-generated. Model: Claude Opus 5.5
The matching kept a list of adjacent surface facets for every vertex and intersected these lists for each face of each cell. It now marks the vertices of the surface mesh and looks up only the faces of cells with at least three such vertices, in the surface faces sorted by their vertices. On the mesh with 3.6M cells (151016 boundary faces), 161886 faces remain to be looked up; the matching takes 0.105 s instead of 0.21 s, and the peak memory of reading the mesh drops from 230 MB to 215 MB. The results are identical. Adds a test for a surface face given twice, which remains an error. AI-generated. Model: Claude Opus 5.5
Rank 0 kept the complete parsed mesh until pumgen ended, i.e. also while computing the inspheres and writing. Each part (vertices, elements, groups, boundary conditions, identification) is now released right after it has been scattered, and the surface mesh right after the boundary conditions have been matched. On the mesh with 3.6M cells and one rank, the memory after reading drops from 392 MB to 192 MB and the peak of the whole run from 534 MB to 388 MB; the peak is now reached while the elements are scattered. The scatter no longer dereferences the data of empty vectors on the other ranks. AI-generated. Model: Claude Opus 5.5
The settings of velocity-aware meshing (refinement cuboids, frequency, easi file, elements per wavelength) and their XML parsing move from MeshAttributes to src/sizing, as does the mesh size computation from EasiMeshSize, which now evaluates the easi model for many points at once. EasiMeshSize keeps only the SimModSuite part, i.e. finding the region of a point, and passes single points to the shared computation; the settings and the size formula stay the same. easi becomes a dependency of its own (option EASI, implied by SIMMETRIX), so that tools without SimModSuite can use the same settings and sizes. AI-generated. Model: Claude Opus 5.5
CellVertices fetches the coordinates of the vertices of the local cells which are owned by other ranks once and provides the vertices per cell, so that further per-cell computations can share them with the insphere calculation. AI-generated. Model: Claude Opus 5.5
pumgen-sizefield reads the VelocityAwareMeshing element of a mesh attributes file and writes, for each refinement cuboid, a binary grid for gmsh's Structured field with the mesh size of the cuboid's group (given by bypassFindRegionAndUseGroup, which is required here). The easi model is evaluated for all grid nodes at once, distributed over the ranks in x slabs; each rank writes its slab with MPI-IO. Each node takes the minimum of its 3x3x3 neighbourhood, so that gmsh's trilinear interpolation never exceeds the smallest size sampled at the corners of a grid cell; the grid covers the cuboid with a margin of two spacings, so that all grid cells which intersect the cuboid have defined corners. A .geo file restricts each field to the physical volume of its group and uses their minimum as background field; it refers to the grids relative to its own directory. On the layered test model (20 x 20 x 10 km, 2 km sediment), gmsh's HXT meshes with these fields have the same edge length distribution as meshes from the exact per-point size callback (816945 vs. 815206 cells). AI-generated. Model: Claude Opus 5.5
With --velocity-check <mesh attributes>, pumgen compares the edge lengths of the cells in the refinement regions of the VelocityAwareMeshing settings with the mesh size at their barycentre, for the material of the cell's group. It logs the distribution of the mean and of the longest edge length over the mesh size (quantiles from a histogram with 0.23 % wide bins, the maximum, and the number of cells above 1). This works for meshes from any source, so that meshes of different generators can be compared by the same measure. For the layered test model, meshed by gmsh with a size callback, the check reproduces the independently computed quantiles (mean edge: p50 1.37, p90 1.57, p99 1.71, max 2.14). AI-generated. Model: Claude Opus 5.5
For binary MSH 4.1 files, rank 0 now only reads the structure of the file (entities, block headers, periodic nodes) and skips the node and element data. With the resulting offsets, every rank reads a share of the nodes, cells and boundary faces. The nodes are sent to the ranks owning them in the chunk distribution; the cells are read directly by their owners, in file order as before. The boundary faces go to the ranks given by a hash of their vertices, their vertices are marked in a bit set shared by all ranks, and only cell faces whose vertices all lie on the surface mesh are looked up. Periodic nodes are joined into classes on rank 0 and sent to the owners. On the mesh with 3.6M cells, the output is identical to the one of the serial reader for 1, 2, 4 and 8 ranks. The highest peak memory of a rank drops from 255 MB to 161 MB with 4 ranks and from 235 MB to 129 MB with 8 ranks; with one rank, time and memory stay about the same (1.07 s vs. 1.11 s, 404 MB vs. 378 MB). --gmsh-reader=serial keeps reading on rank 0; ASCII files are always read there. AI-generated. Model: Claude Opus 5.5
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.
Also add a lot of tests.