Skip to content

HDF5 mesh parser as an alternative for SeisSol 3D Meshes#97

Draft
piyushkarki wants to merge 67 commits intopiyush/thomas/dockerfrom
piyush/HDF5-mesh-parser-alternative
Draft

HDF5 mesh parser as an alternative for SeisSol 3D Meshes#97
piyushkarki wants to merge 67 commits intopiyush/thomas/dockerfrom
piyush/HDF5-mesh-parser-alternative

Conversation

@piyushkarki
Copy link
Copy Markdown
Collaborator

@piyushkarki piyushkarki commented Mar 14, 2026

Closes Issue: N/A

Summary:

This PR introduces support for parsing HDF5 meshes (specifically targeting SeisSol 3D Meshes) as an alternative to the existing GMSH format. To achieve this cleanly, the mesh parsing logic has been refactored to use a polymorphic base class meshParser, with GMSHParser and H5Parser as specific implementations. HDF5 is introduced as an optional dependency, ensuring the core application can still build without it.

Description: Changes

  • Mesh Parser Abstraction: Introduced a new meshParser base class and refactored the existing GMSH parsing logic into a derived GMSHParser class.

  • HDF5 Implementation: Added the H5Parser implementation to read .h5 files.

  • Dynamic Format Resolution: Added runtime/static checks (meshInGMSHFile(), meshInH5File()) to dynamically instantiate the correct parser based on the mesh file format.

  • Optional HDF5 Dependency (CMake): Based on feedback, HDF5 is not a hard dependency. Added the ENABLE_HDF5 CMake flag.

    • If ENABLE_HDF5=ON, the HDF5 libraries are linked.
    • If ENABLE_HDF5=OFF, the H5Parser class remains a valid compilation unit, but relies on a stub implementation that safely throws a std::runtime_error if a user attempts to parse an .h5 file.
  • Docker/CI Integration: Updated Dockerfiles and GitHub Actions to handle HDF5 dependencies dynamically based on hdf5_version.

Tests

  • CI Matrix Expansion: Updated .github/workflows/build-and-test-tandem.yml to include an enable_hdf5: ["ON", "OFF"] matrix configuration. All 8 configurations (gcc/clang, 2D/3D, ON/OFF) are successfully passing.
  • H5Parser Tests:
    • Check if each face has correct tags based on the 32 bit encoding

    • Deduplicate faces if a face is shared by two elements

    • Check if correctly parsed if only some faces of an element are tagged

    • Correct error handling of non existent file

    • No error before any parse attempt

    • Throw correct error in case no HDF5 support available (ENABLE_HDF5=OFF)

      3D Restriction: Note that HDF5 parsing is currently explicitly restricted to 3D problems; using it in 2D will safely throw an error.

Additional Info

  • Merge Conflicts: There is currently a merge conflict in test/CMakeLists.txt that I will resolve before moving this out of Draft.

Note: The Docker images for dependencies have been updated/pushed to include the hdf5_version tags.

tldr:

Compile with -DENABLE_HDF5=ON and include the mesh file with the correct extension in the toml script mesh_file = "filename.h5" to run simulation with HDF5 parser. Only works with 3D meshes.

@piyushkarki piyushkarki changed the base branch from release-v1.2.0 to piyush/thomas/docker March 14, 2026 16:21
@piyushkarki piyushkarki self-assigned this Mar 17, 2026
@hpc4geo
Copy link
Copy Markdown
Collaborator

hpc4geo commented Mar 17, 2026

@piyushkarki In general, I would like HDF5 to not become a hard dependency. Could you please add support in CMAKE to allow HDF5 to be optional (maybe with a flag such as TANDEM_HAVE_HDF5 for example). My idea would be to have H5Parser.cpp remain as a compilation unit whether TANDEM_HAVE_HDF5 is TRUE or FALSE. Inside H5Parser.cpp you would guard all required H5 calls and headers. Then you could also guard the entire h5 parser class. In the event that TANDEM_HAVE_HDF5 = FALSE, the alternate class could be instantiated but would return a meaningful runtime error if you ever tried to perform a parse of a H5 class. The other more intrusive option would be to throw an exception if the H5 parser class was instantiated if TANDEM_HAVE_HDF5 = FALSE. I say this is intrusive because it would also require you add TANDEM_HAVE_HDF5 guards in static.cpp and tandem.cpp. I think I prefer the design that you can "create" a H5 parser object independent of the value of TANDEM_HAVE_HDF5 , and that it should throw an exception if you call any of its methods if TANDEM_HAVE_HDF5 = FALSE.

What are your thoughts?

@piyushkarki
Copy link
Copy Markdown
Collaborator Author

@hpc4geo I agree with the assessment. I will be pushing the strategy you described to always create an H5 parser and throwing an exception if it is used without the proper CMake variable. In this regard, I have also modified tests for the particular HDF5 parser class to test for the parsing features if the build is built with ENABLE_HDF5 options ON and test for the correct error throw in case the build has it turned off.

@piyushkarki piyushkarki force-pushed the piyush/HDF5-mesh-parser-alternative branch from a47a22a to 1186419 Compare March 24, 2026 17:38
@yohaimagen
Copy link
Copy Markdown
Collaborator

Hey @piyushkarki, one small request regarding the tandem-base image.

Currently the image is only tagged with the PETSc version (e.g. ghcr.io/tear-erc/tandem/tandem-base:3.22.5). Could you also push a :latest tag pointing to the same image whenever the base is built?
That way, the release image workflow in #96 can always reference ghcr.io/tear-erc/tandem/tandem-base:latest without hardcoding a specific PETSc version.

In build-dependencies.yml, it would just be adding latest alongside the versioned tag in the tags field of the Build and push base image step:

tags: |
${{ inputs.base_image }}:${{ inputs.petsc_version }}
${{ inputs.base_image }}:latest


// The boundary condition for this tet is packed as four 8-bit tags,
// one per face: bits [7:0]=face0, [15:8]=face1, [23:16]=face2, [31:24]=face3.
const uint32_t boundaryCondition = boundaryData[i];
Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hi, SeisSol dev here.

One short comment w.r.t. the boundary formats; there are actually three of them (with this 8-bit format being the most common for now). It's effectively denoted in all modern PUML files by an attribute (if not, the 8-bit format can be assumed in almost all cases). The same is also true for periodicity by now.

Cf. https://github.com/SeisSol/SeisSol/blob/bb700decaa218162de1e157ce49a1da00488cc61/src/Geometry/PUMLReader.h#L27-L44
and https://seissol.readthedocs.io/en/latest/puml-mesh-format.html#puml-mesh-format .

It's mostly a corner case (as the larger formats are only used for fault networks with many different tags to my knowledge); but I wanted to mention it just in case—maybe it's even more relevant for tandem than SeisSol.

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I see. i will document this in the comment. Thank you!


// Extract the tag for this face. A tag of 0 means interior (skip).
const uint8_t faceTag = (boundaryCondition >> (8 * face)) & 0xFF;
if (faceTag == 0) {
Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It's better to check for !(faceTag == 3 || faceTag > 64) instead (or rather, only execute the core for these faces); those are guaranteed DR faces.

Depending on which faces you scan here (I presume, all in the mesh?), you'll need to also exclude all boundary conditions. Also face tag 6 is another, rarely-used interior face tag.

Cf. https://github.com/SeisSol/SeisSol/blob/master/src/Geometry/PUMLReader.cpp#L83-L109 for a list of all face types. (the tag > 64 are also DR faces; I'm gonna clarify that in a future PR)

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I see. But this portion is aggregating all tags of interest that go into the DGOperator and not just DR faces (so 1,3 and 5). Also, I think currently we do not want to explicitly support all SeisSol mesh BC tags. We would rather want the user to modify the mesh so as to translate the tags to the ones supported by Tandem. So for us, the tags of interest would rather be 1,3 and 5 for the time being but fault tagging will soon change that.

Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ah I see; just looked through the corresponding builder files. ... (talking about https://github.com/TEAR-ERC/tandem/blob/main/src/form/BC.h )

To ask: what do the Natural and Dirichlet BCs correspond or "do" to in tandem? And out of interest; why can the None boundary condition faces be ignored here?

Maybe, to raise the further question: does it even make sense to keep the same boundary tag arrays for SeisSol and tandem—or could something happen differently numerically/physically?
Because then the tags might be just as well be interpreted as bare numbers and we'll need some extra mapping file (or an extra boundary array) for either SeisSol or tandem anyways.
But if not, then maybe it's nice if a SeisSol mesh file could be used out of the box in tandem (given that all BCs are supported by both codes)—or what would you think about that?

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Natural corresponds to a free-surface boundary condition (zero traction), while Dirichlet is used to prescribe displacement fields. Both enter through the boundary integrals in the weak form. You are correct to question if None can be ignored. It can't! It represents interior boundaries that are not assigned a specific BC but still need to be included in certain weak form terms. This likely comes from an earlier implementation and will be fixed. This did not get caught when I tested with a mesh because the MeshData class when distributing the BC's assigns any face without a BC to be BC::None. So that is why I wasn't noticing this. Thanks!

Regarding tagging: in Tandem we currently use 1 (free surface), 3 (fault), and 5 (Dirichlet), with 7 (traction) and 9 (free slip) planned. Facet tags are also coming. Since these lists are evolving and the physics differ somewhat from SeisSol, maintaining identical tag schemes may not be practical, though there will be significant overlap. Alignment is definitely worth discussing between the teams.

Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks alot; it does look like the mapping already diverged by now I guess. :)

For reference, here's how the current SeisSol face tags should map to tandem (given that I didn't miss any description):

  • 0, 6 -> 0 (none)
  • 1 -> 1 (natural)
  • 3, 65, 66, ..., 2**31-1 -> 3 (DR)
  • 4, 7 -> 5 (Dirichlet—for SeisSol/7 you'd need time dependent support)
  • 5 -> not supported

The new BCs sound interesting; maybe those will need to go to SeisSol as well (though they will have different tags probably).

Two ad hoc ideas:

We could have a "common" tag list or an option to use a common tag list. Or one tag list adopts the other on demand (parameter option).
Though the question might be—will tandem and SeisSol operations need different tags? Or could some boundary condition even split into two different ones when switching codes for the same scenario?

To suggest; one possible intermediate solution might be to introduce an extra tag map file for SeisSol and/or tandem. I.e. something like a type-to-ID-maps as in

- freeSurface: 0, 1
- dynamicRupture: 3, 4, 65
- regular: 100

(implementing that on the SeisSol side is pretty straight-forward)

That way, we could keep the total amount of boundary arrays to 1.

Copy link
Copy Markdown
Collaborator

@Thomas-Ulrich Thomas-Ulrich left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

LGTM but see comments

cxx: clang++-18
dimension: [2, 3]
degree: [3]
enable_hdf5: ["ON", "OFF"]
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I guess enable_hdf5 is not necessary if the image as the hdf5 version in its name?
(suggesting is it always enabled)

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes HDF5 is always available in the image but we would want to test for when someone HDF5 is off too. And this way we get to use single image and still test scenarios flexibly because the compilation does not know whether HDF5 is there or not unless told explicitly by this option.

name: logs-${{ matrix.dimension }}d-p${{ matrix.degree }}-hdf5-${{ matrix.enable_hdf5 }}-${{ matrix.compiler.cc }}
path: |
${{ github.workspace }}/build_${{ matrix.dimension }}d_p${{ matrix.degree }}/Testing/Temporary/*.log
${{ github.workspace }}/build_${{ matrix.dimension }}d_p${{ matrix.degree }}_hdf5_${{ matrix.enable_hdf5 }}/Testing/Temporary/*
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

keeping *.log would make things clearer

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Done

id: check-base
run: |
if docker manifest inspect "${{ inputs.base_image }}:${{ inputs.petsc_version }}" > /dev/null 2>&1; then
# Added the HDF5 version to the tag check
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

remove comment

std::optional<ProbeOutputConfig> domain_probe_output;
std::optional<GfCheckpointConfig> gf_checkpoint_config;
TsCheckpointConfig ts_checkpoint_config;
bool meshInGMSHFile() {
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

meshInGMSHFormat
or isGMSHFile

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Not relevant now. Changed the logic.

@@ -2,6 +2,8 @@
#define GLOBALSIMPLEXMESHBUILDER_20200901_H

#include "io/GMSHParser.h"
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

can be removed?

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes removed.

app/static.cpp Outdated
if (rank == 0) {
GMSHParser parser(&builder);
ok = parser.parseFile(*cfg->mesh_file);
std::unique_ptr<meshParser> parser; // Pointer to the base class
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe all of that could be part of meshParser to avoid code duplication?

--enable-parallel \
--enable-shared \
CC=mpicc CXX=mpicxx && \
make -j$(nproc) && \
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

ninja?

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Will change this at the end after approval to avoid creating a new image but yes.

void GlobalSimplexMeshBuilder<D>::addElement(long type, long tag, long* node,
std::size_t numNodes) {
if (is_gmsh_simplex<D>(type)) {
if (is_gmsh_simplex<D>(type)) { // Higher order simplex parsing (tetrahedra for example)
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Parsing D-dimensional simplices (e.g., Tetrahedra if D=3)
?

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

added detail. Thanks!

H5Dclose(dataset);
}
} else {
std::cout << "Invalid dataset name" << std::endl;
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Not considered volume tagging and the mesh example I got did not have group (or had no significantly important values in them - can't recall). But added this now.

return true;
}

bool H5Parser::parseBoundary(hid_t file) {
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

do you really need this function? It is not doing what its names suggests.

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We do need the function. But I agree about the name. Changed it to readBoundaryData

@piyushkarki piyushkarki force-pushed the piyush/HDF5-mesh-parser-alternative branch 2 times, most recently from 17322c6 to 0608ee5 Compare April 8, 2026 22:49
@piyushkarki piyushkarki force-pushed the piyush/HDF5-mesh-parser-alternative branch from e227a57 to 0b0b3d8 Compare April 11, 2026 02:20
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

5 participants