diff --git a/src/underworld3/discretisation/__init__.py b/src/underworld3/discretisation/__init__.py index 30ce1d1a..687b3451 100644 --- a/src/underworld3/discretisation/__init__.py +++ b/src/underworld3/discretisation/__init__.py @@ -23,3 +23,4 @@ from .discretisation_mesh import meshVariable_lookup_by_symbol from .discretisation_mesh import petsc_dm_find_labeled_points_local from .discretisation_mesh import _from_gmsh +from .discretisation_mesh import _gmsh_to_h5 diff --git a/src/underworld3/discretisation/discretisation_mesh.py b/src/underworld3/discretisation/discretisation_mesh.py index 42afa698..5b2afe40 100644 --- a/src/underworld3/discretisation/discretisation_mesh.py +++ b/src/underworld3/discretisation/discretisation_mesh.py @@ -96,18 +96,15 @@ def wrapper(final): @timing.routine_timer_decorator -def _from_gmsh(filename, comm=None, markVertices=False, useRegions=True, useMultipleTags=True): - """Read a Gmsh .msh file from `filename`. - - :kwarg comm: Optional communicator to build the mesh on (defaults to - COMM_WORLD). - """ - - ## NOTE: - this should be smart enough to serialise the msh conversion - ## and then read back in parallel via h5. This is currently done - ## by every gmesh mesh +def _gmsh_to_h5( + filename, + markVertices=False, + useRegions=True, + useMultipleTags=True, +): + """Convert a Gmsh file to PETSc HDF5 without loading the resulting DM.""" - comm = comm or PETSc.COMM_WORLD + h5_filename = filename + ".h5" options = PETSc.Options() options["dm_plex_hash_location"] = None @@ -152,11 +149,11 @@ def _from_gmsh(filename, comm=None, markVertices=False, useRegions=True, useMult # this module, so a top-level import would close a cycle. from underworld3.meshing._mesh_files import _scratch_name - scratch = _scratch_name(filename + ".h5") + scratch = _scratch_name(h5_filename) viewer = PETSc.ViewerHDF5().create(str(scratch), "w", comm=PETSc.COMM_SELF) viewer(plex_0) viewer.destroy() - os.replace(scratch, filename + ".h5") + os.replace(scratch, h5_filename) finally: # The gmsh import options are import-time scratch — meaningful only for # the createFromFile above. Clear the whole namespace so a value set by @@ -166,12 +163,39 @@ def _from_gmsh(filename, comm=None, markVertices=False, useRegions=True, useMult # read as 2-D). Runs on success or failure. _clear_gmsh_import_options() - # Now we have an h5 file and we can hand this to _from_plexh5. The barrier - # is what the atomic write above makes necessary AND sufficient: the other - # ranks must not look for the file before rank 0 has renamed it into place. + # The barrier ensures every rank sees the complete atomically-renamed file. uw.mpi.barrier() - return _from_plexh5(filename + ".h5", comm, return_sf=True) + return filename + ".h5" + + +@timing.routine_timer_decorator +def _from_gmsh( + filename, + comm=None, + markVertices=False, + useRegions=True, + useMultipleTags=True, +): + """Read a Gmsh .msh file from `filename`. + + :kwarg comm: Optional communicator to build the mesh on (defaults to + COMM_WORLD). + """ + + ## NOTE: - this should be smart enough to serialise the msh conversion + ## and then read back in parallel via h5. This is currently done + ## by every gmesh mesh + + comm = comm or PETSc.COMM_WORLD + h5_filename = _gmsh_to_h5( + filename, + markVertices=markVertices, + useRegions=useRegions, + useMultipleTags=useMultipleTags, + ) + + return _from_plexh5(h5_filename, comm, return_sf=True) @timing.routine_timer_decorator diff --git a/src/underworld3/meshing/spherical.py b/src/underworld3/meshing/spherical.py index 0f79a4d4..4f8c97a2 100644 --- a/src/underworld3/meshing/spherical.py +++ b/src/underworld3/meshing/spherical.py @@ -21,6 +21,7 @@ from underworld3 import VarType from underworld3.coordinates import CoordinateSystemType from underworld3.discretisation import _from_gmsh as gmsh2dmplex +from underworld3.discretisation import _gmsh_to_h5 as gmsh2h5 import underworld3.timing as timing import underworld3.cython.petsc_discretisation @@ -509,6 +510,7 @@ def SphericalShellInternalBoundary( refinement=None, gmsh_verbosity=0, verbose=False, + write_mesh_files_only=False, ): """ Generates a spherical shell with an internal boundary using Gmsh. The function creates a 3D mesh of a spherical shell @@ -536,11 +538,19 @@ def SphericalShellInternalBoundary( Gmsh output verbosity (0=quiet). Default is 0. verbose : bool, optional If True, print additional information. Default is False. + write_mesh_files_only : bool, optional + If True, write the Gmsh ``.msh`` and PETSc ``.msh.h5`` files, return + the HDF5 path as a string, and stop without constructing an in-memory + Underworld mesh. This skips HDF5 reloading, coordinate-system setup, + and cell-region classification. If False, return a fully initialized + :class:`underworld3.discretisation.Mesh`, including the ``Inner`` and + ``Outer`` cell-region labels. Default is False. Returns ------- - Mesh - The generated spherical shell mesh with internal boundary. + Mesh or str + The PETSc HDF5 path when ``write_mesh_files_only=True``; otherwise, a + fully initialized spherical-shell mesh. Examples -------- @@ -696,6 +706,14 @@ def bbox_radius(dimtag): write_gmsh(uw_filename) gmsh.finalize() + if write_mesh_files_only: + return gmsh2h5( + uw_filename, + markVertices=True, + useRegions=True, + useMultipleTags=True, + ) + # Ensure boundaries conform (if refined) # This is equivalent to a partial function because it already # knows the configuration of THIS spherical mesh and diff --git a/tests/test_0502_boundary_integrals.py b/tests/test_0502_boundary_integrals.py index fe376474..2dbe5663 100644 --- a/tests/test_0502_boundary_integrals.py +++ b/tests/test_0502_boundary_integrals.py @@ -414,6 +414,61 @@ def test_bd_integral_spherical_internal_boundary_areas(): ) +@pytest.mark.level_2 +@pytest.mark.tier_b +def test_spherical_internal_boundary_mesh_files_only(tmp_path, monkeypatch): + """File generation can bypass Mesh construction and preserve labels.""" + from enum import Enum + + import underworld3.meshing.spherical as spherical + from underworld3.coordinates import CoordinateSystemType + + mesh_file = str(tmp_path / "spherical_internal_mesh_files_only.msh") + + def fail_mesh_construction(*args, **kwargs): + raise AssertionError( + "write_mesh_files_only=True must not construct an Underworld Mesh" + ) + + with monkeypatch.context() as patch: + patch.setattr(spherical, "Mesh", fail_mesh_construction) + h5_file = spherical.SphericalShellInternalBoundary( + radiusOuter=_R_SHELL_OUTER, + radiusInternal=_R_SHELL_INTERNAL, + radiusInner=_R_SHELL_INNER, + cellSize=0.25, + filename=mesh_file, + write_mesh_files_only=True, + ) + + assert h5_file == f"{mesh_file}.h5" + assert (tmp_path / "spherical_internal_mesh_files_only.msh").is_file() + assert (tmp_path / "spherical_internal_mesh_files_only.msh.h5").is_file() + + class Boundaries(Enum): + Centre = 1 + Lower = 11 + Internal = 12 + Upper = 13 + All_Boundaries = 1001 + + reloaded_mesh = uw.discretisation.Mesh( + h5_file, + degree=1, + qdegree=2, + coordinate_system_type=CoordinateSystemType.SPHERICAL, + useMultipleTags=True, + useRegions=True, + markVertices=True, + boundaries=Boundaries, + ) + + for boundary in ("Lower", "Internal", "Upper"): + label = reloaded_mesh.dm.getLabel(boundary) + assert label is not None + assert label.getNumValues() > 0 + + def _build_spherical_shell_for_integrals(): from underworld3.meshing import SphericalShell