diff --git a/src/underworld3/discretisation/discretisation_mesh.py b/src/underworld3/discretisation/discretisation_mesh.py index b2527343..bc83ad18 100644 --- a/src/underworld3/discretisation/discretisation_mesh.py +++ b/src/underworld3/discretisation/discretisation_mesh.py @@ -3705,6 +3705,7 @@ def get_min_radius_old(self) -> float: return self._min_radius + @uw.collective_operation def get_min_radius(self) -> float: """ This method returns the global minimum distance from any cell centroid to a face. @@ -3721,6 +3722,7 @@ def get_min_radius(self) -> float: return all_min_radii.min() + @uw.collective_operation def get_max_radius(self) -> float: """ This method returns the global maximum distance from any cell centroid to a face. @@ -3735,6 +3737,35 @@ def get_max_radius(self) -> float: return all_max_radii.max() + @uw.collective_operation + def get_mean_radius(self) -> float: + """ + Global mean of the characteristic cell length scale + (``volume^(1/dim)``, i.e. the equivalent radius derived from each + cell's volume — the same quantity averaged by ``get_min_radius`` + and ``get_max_radius`` to obtain global min/max). Parallel-safe + via MPI allreduce of the local sum and count. + + Together with :meth:`get_min_radius` / :meth:`get_max_radius` + this is the canonical "mesh length" API. Use this anywhere you + need a representative h0 (smoothing-length defaults, diffusion- + stability heuristics, problem-scale normalisation) rather than + reaching for the rank-local ``self._radii`` array, which gives + different answers on different MPI ranks and leaks downstream + (e.g. into JIT C source via per-rank pointwise-function inputs). + """ + + import numpy as np + from mpi4py import MPI + + radii = np.asarray(self._radii) + local_sum = float(radii.sum()) + local_n = int(radii.size) + if uw.mpi.size > 1: + local_sum = uw.mpi.comm.allreduce(local_sum, op=MPI.SUM) + local_n = uw.mpi.comm.allreduce(local_n, op=MPI.SUM) + return local_sum / max(local_n, 1) + # This should be deprecated in favour of using integrals def stats(self, uw_function, uw_meshVariable, basis=None): """ diff --git a/tests/parallel/ptest_0008_mesh_radii_accessors.py b/tests/parallel/ptest_0008_mesh_radii_accessors.py new file mode 100644 index 00000000..d56b40db --- /dev/null +++ b/tests/parallel/ptest_0008_mesh_radii_accessors.py @@ -0,0 +1,46 @@ +"""MPI smoke test for the parallel-safe mesh radius accessors. + +Run via mpi_runner.sh (mpirun -np N python ptest_0008_*.py). +Asserts: + - min/mean/max return positive, ordered values on every rank + - every rank sees IDENTICAL min, mean, max (parallel-consistent) + - the methods are tagged as @uw.collective_operation +""" +import underworld3 as uw + +mesh = uw.meshing.Annulus( + radiusOuter=1.0, radiusInner=0.5, cellSize=1.0 / 16, qdegree=3) + +rmin = mesh.get_min_radius() +rmean = mesh.get_mean_radius() +rmax = mesh.get_max_radius() + +assert rmin > 0.0, f"rank {uw.mpi.rank}: min must be positive, got {rmin}" +assert rmin <= rmean <= rmax, ( + f"rank {uw.mpi.rank}: expected min <= mean <= max, got " + f"{rmin} <= {rmean} <= {rmax}") + +# All ranks must agree to within float-round-off (the values come from +# global allreduce / gather, so any disagreement would indicate the +# accessor leaked rank-local state). +all_min = uw.mpi.comm.allgather(rmin) +all_mean = uw.mpi.comm.allgather(rmean) +all_max = uw.mpi.comm.allgather(rmax) + +tol = 1.0e-12 +assert max(all_min) - min(all_min) < tol, ( + f"min_radius differs across ranks: {all_min}") +assert max(all_mean) - min(all_mean) < tol, ( + f"mean_radius differs across ranks: {all_mean}") +assert max(all_max) - min(all_max) < tol, ( + f"max_radius differs across ranks: {all_max}") + +# Decorator must be present so any selective_ranks() misuse fails fast +for name in ("get_min_radius", "get_mean_radius", "get_max_radius"): + fn = getattr(mesh, name) + assert getattr(fn, "_is_collective", False), ( + f"rank {uw.mpi.rank}: {name} must be @uw.collective_operation") + +if uw.mpi.rank == 0: + print(f"OK: min={rmin:.6f} mean={rmean:.6f} max={rmax:.6f} " + f"agree across {uw.mpi.size} ranks", flush=True) diff --git a/tests/test_0008_mesh_radii_accessors.py b/tests/test_0008_mesh_radii_accessors.py new file mode 100644 index 00000000..7c79378c --- /dev/null +++ b/tests/test_0008_mesh_radii_accessors.py @@ -0,0 +1,50 @@ +"""Smoke tests for the parallel-safe mesh radius accessors. + +Verifies that :meth:`Mesh.get_min_radius`, :meth:`Mesh.get_mean_radius`, +and :meth:`Mesh.get_max_radius` return sensible values in serial. +A matching MPI smoke test lives in tests/parallel/. +""" +import pytest + +pytestmark = pytest.mark.level_1 + + +def _check_radii_ordering(mesh): + rmin = mesh.get_min_radius() + rmean = mesh.get_mean_radius() + rmax = mesh.get_max_radius() + assert rmin > 0.0, f"min radius must be positive, got {rmin}" + assert rmin <= rmean <= rmax, ( + f"expected min <= mean <= max, got " + f"min={rmin}, mean={rmean}, max={rmax}") + + +def test_mesh_radii_accessors_annulus(): + import underworld3 as uw + mesh = uw.meshing.Annulus( + radiusOuter=1.0, radiusInner=0.5, + cellSize=1.0 / 16, qdegree=3) + _check_radii_ordering(mesh) + + +def test_mesh_radii_accessors_box(): + import underworld3 as uw + mesh = uw.meshing.UnstructuredSimplexBox( + minCoords=(0.0, 0.0), maxCoords=(1.0, 1.0), + cellSize=1.0 / 16) + _check_radii_ordering(mesh) + + +def test_mesh_radii_accessors_are_collective(): + """The @uw.collective_operation decorator must be present so a + selective_ranks() misuse fails fast rather than deadlocking.""" + import underworld3 as uw + mesh = uw.meshing.Annulus( + radiusOuter=1.0, radiusInner=0.5, + cellSize=1.0 / 16, qdegree=3) + for name in ("get_min_radius", "get_mean_radius", + "get_max_radius"): + fn = getattr(mesh, name) + assert getattr(fn, "_is_collective", False), ( + f"{name} must be decorated with @uw.collective_operation " + f"so misuse inside selective_ranks() is caught early.")