The adapt level loop stops on a collective fact, not a rank-local one

src/underworld3/discretisation/discretisation_mesh.py

@@ -8226,6 +8226,34 @@ def eval_metric(centroids):
             def eval_metric(centroids):
                 return numpy.asarray(_sampler(metric_sym, centroids)).reshape(-1)
 
+        def marking_metric(centroids):
+            """``eval_metric``, clipped, with "nobody has cells" decided together.
+
+            Returns ``None`` when NO rank owns cells, which the marking loops
+            read as "mark nothing". The point of routing the emptiness question
+            through a reduction is that the loops below then have a collective
+            fact to stop on, instead of each rank deciding for itself whether
+            to leave — which is the defect this replaced (#512).
+
+            Note on what this does NOT fix. #512 describes the rank-local
+            ``if cur_h.size:`` guards as skipping a collective, on the grounds
+            that ``eval_metric`` is ``global_evaluate`` for a field or
+            expression metric. Measured at np=2, that is not so: with one rank
+            asleep for 10 s, the other's ``global_evaluate`` returned in 0.06 s,
+            both for in-domain points and for points that strand outside the
+            mesh entirely. It does not wait for its peers on these shapes, so a
+            cell-less rank skipping it does not hang. The guards are routed
+            through here for uniformity with the stop below, not because they
+            were hanging.
+            """
+
+            # mpi4py allreduce defaults to MPI.SUM, as at the collective stop
+            # further down.
+            if uw.mpi.comm.allreduce(int(centroids.shape[0])) == 0:
+                return None
+
+            return numpy.clip(eval_metric(centroids), 1e-30, None)
+
         def cell_geometry(dm):
             """Per-cell centroid and size for every cell of a simplicial DM.
 
@@ -8458,8 +8486,8 @@ def _relax_generation(engine_obj, carry, rcarry):
             n_gen = dim * max_levels
             for level in range(n_gen):
                 centroids, cur_h, cs = cell_geometry(current_dm)
-                if cur_h.size:
-                    M = numpy.clip(eval_metric(centroids), 1e-30, None)
+                M = marking_metric(centroids)
+                if M is not None and cur_h.size:
                     h_target = 1.0 / numpy.sqrt(M)
                     sel = numpy.where(cur_h > h_target)[0]
                     if node_budget is not None and sel.size > node_budget:
@@ -8543,8 +8571,8 @@ def _relax_generation(engine_obj, carry, rcarry):
             current_dm = base_finest
             for level in range(n_pass):
                 centroids, _proxy_h, cs = cell_geometry(current_dm)
-                if centroids.shape[0]:
-                    M = numpy.clip(eval_metric(centroids), 1e-30, None)
+                M = marking_metric(centroids)
+                if M is not None and centroids.shape[0]:
                     h_target = 1.0 / numpy.sqrt(M)
                     diameter = edge_split.cell_diameters(current_dm)
                     sel = numpy.where(diameter > h_target)[0]
@@ -8639,6 +8667,13 @@ def _relax_generation(engine_obj, carry, rcarry):
                                   regions=rcarry)
             n_gen = dim * max_levels
             for level in range(n_gen):
+                # Rank-local marking and a rank-local break, deliberately: this
+                # is the pure-Python cell-list engine, reached only when the
+                # native uwnvb transform is absent, and that combination raises
+                # NotImplementedError at np>1 further up. It is serial by
+                # construction, so the collective discipline the other engines
+                # need here (#512) has nothing to protect. Port it before
+                # porting this loop.
                 centroids, cur_h, cids = nvb.centroids_h()
                 M = numpy.clip(eval_metric(centroids), 1e-30, None)
                 h_target = 1.0 / numpy.sqrt(M)
@@ -8696,17 +8731,30 @@ def _relax_generation(engine_obj, carry, rcarry):
             current_dm = base_finest             # base finest, current geometry
             for level in range(max_levels):
                 centroids, cur_h, cs = cell_geometry(current_dm)
-                if cur_h.size == 0:
-                    break
 
                 # Metric M = 1/h_target² at the cell centroids (parent field).
                 # A callable is evaluated directly on the centroids; a field/expr
                 # goes through global_evaluate (parallel) or evaluate (serial).
-                M = numpy.clip(eval_metric(centroids), 1e-30, None)
-                h_target = 1.0 / numpy.sqrt(M)
+                M = marking_metric(centroids)
+                if M is None:
+                    break                     # no rank has cells: leaving together
 
-                refine = numpy.where(cur_h > h_target)[0]
-                if refine.size == 0:
+                if cur_h.size:
+                    h_target = 1.0 / numpy.sqrt(M)
+                    refine = numpy.where(cur_h > h_target)[0]
+                else:
+                    refine = numpy.empty(0, dtype=int)
+
+                # Collective stop. `refine.size == 0` is a rank-local fact: this
+                # rank's cells are all fine enough, which says nothing about its
+                # peers'. Breaking on it alone takes the rank out of the loop for
+                # good while the others go round again into the DM refinement,
+                # which IS collective — so the job hangs, and not latently.
+                # Measured at np=2 with a callable metric demanding h=0.01 inside
+                # r<0.15 of one corner and nothing elsewhere: the rank holding no
+                # corner cells left at the first level and `mpirun` reached its
+                # timeout with neither rank returning. See test_0873.
+                if uw.mpi.comm.allreduce(int(refine.size)) == 0:
                     if verbose:
                         uw.pprint(0, f"[adapt] level {level}: no cells need refinement")
                     break

src/underworld3/utilities/reconnect.py

@@ -833,20 +833,37 @@ def _rebuild_point_sf(new, dm, point_map, nroots):
     # petsc4py will not narrow an index array for us — the graph must arrive as
     # PETSc's own integer type or `setGraph` raises an unsafe-cast TypeError.
     new_sf = PETSc.SF().create(comm=dm.comm)
-    if ilocal is None or not len(ilocal):
+
+    # A rank that shares nothing has no leaves to renumber, but it still has to
+    # reach the verdict below with its peers — which is why the emptiness test
+    # is read here and acted on after, rather than returning straight out of it.
+    shares_nothing = ilocal is None or not len(ilocal)
+    if shares_nothing:
+        leaves = local = remote_index = None
+        deleted_here = False
+    else:
+        leaves = np.asarray(ilocal, dtype=np.int64)
+        local = point_map[leaves - pStart]
+        remote_index = leaf_new[leaves - pStart]
+        deleted_here = bool((local < 0).any() or (remote_index < 0).any())
+
+    # Assert-class: fires only on an internal contract violation. It is still
+    # raised collectively — one rank raising while its peers carry on into the
+    # next collective is a hang on top of the error, and the error is the thing
+    # worth seeing (#512).
+    offenders = [rank for rank, bad
+                 in enumerate(dm.comm.tompi4py().allgather(deleted_here)) if bad]
+    if offenders:
+        raise RuntimeError(
+            f"reconnect: a shared point was deleted on rank(s) {offenders}. "
+            "The removal pass must freeze the seam; see remove_vertices.")
+
+    if shares_nothing:
         new_sf.setGraph(nroots, np.zeros(0, dtype=PETSc.IntType),
                         np.zeros(0, dtype=PETSc.IntType))
         _install_point_sf(new, new_sf)
         return
 
-    leaves = np.asarray(ilocal, dtype=np.int64)
-    local = point_map[leaves - pStart]
-    remote_index = leaf_new[leaves - pStart]
-    if (local < 0).any() or (remote_index < 0).any():
-        raise RuntimeError(
-            "reconnect: a shared point was deleted. The removal pass must "
-            "freeze the seam; see remove_vertices.")
-
     remote = np.empty((len(leaves), 2), dtype=PETSc.IntType)
     remote[:, 0] = np.asarray(iremote).reshape(-1, 2)[:, 0]
     remote[:, 1] = remote_index

tests/parallel/ptest_0845_relax_pinned_band_parallel.py

@@ -61,7 +61,21 @@ def _pinned_indices(mesh, name):
     return np.asarray(iset.getIndices(), dtype=np.int64) - vS
 
 
-def test_pinned_set_is_partition_independent():
+def test_every_rank_pins_every_band_vertex_it_holds():
+    """A vertex in the band is pinned on EVERY rank that holds it.
+
+    The previous form of this test allgathered the pinned coordinates, took the
+    union, and compared ``len(union)`` with ``allreduce(len(union), MAX)``. The
+    union is built from the same gathered list on every rank, so those two
+    numbers are the same by construction and the assertion could not fail
+    whatever the pinning did (#512).
+
+    What can fail is this: a rank-local pin misses a vertex its peer pinned, so
+    the owner is free to move a vertex the seam expects to stay. Every rank
+    therefore has to pin every band vertex present in its own coordinate array,
+    not merely some of them.
+    """
+
     mesh, surface = _fixture()
     name = mesh.label_interface_band(surface, offset=0.0, halo=1)
     X = _coords(mesh)
@@ -70,13 +84,21 @@ def test_pinned_set_is_partition_independent():
     # Compare the pinned COORDINATES, not counts: a shared vertex is held by
     # every rank on the seam, so a count double-counts it and would mask exactly
     # the defect this test exists to catch.
-    local = {(round(float(x), 12), round(float(y), 12)) for x, y in X[idx]}
-    gathered = uw.mpi.comm.allgather(local)
-    union = set().union(*gathered)
+    def key(x, y):
+        return (round(float(x), 12), round(float(y), 12))
+
+    local = {key(x, y) for x, y in X[idx]}
+    union = set().union(*uw.mpi.comm.allgather(local))
+
+    assert union, "nothing pinned; the fixture is not exercising the band"
+
+    held_here = {key(x, y) for x, y in X}
+    missing = (union & held_here) - local
 
-    total = uw.mpi.comm.allreduce(len(union), op=MPI.MAX)
-    assert len(union) == total
-    assert total > 0, "nothing pinned; the fixture is not exercising the band"
+    assert not missing, (
+        f"rank {uw.mpi.rank} holds {len(missing)} band vertex/vertices that a "
+        f"peer pinned and it did not: {sorted(missing)[:5]} — a rank-local pin "
+        "lets the owner move a vertex the seam expects to stay")
 
 
 def test_pinned_vertices_including_shared_ones_do_not_move():

tests/parallel/test_0873_adapt_collective_stop_mpi.py

@@ -0,0 +1,104 @@
+"""The adapt level loop must stop on a collective fact, not a rank-local one (#512).
+
+Each level of the marking loop ends with a DM refinement, which is collective.
+Deciding to leave that loop is therefore a decision every rank has to take
+together. The sbr path used to break on ``refine.size == 0`` — "this rank's own
+cells are all fine enough" — which says nothing about anyone else's. The rank
+that finished first left for good and its peers went round again into the
+refinement and waited.
+
+It is not latent. Measured at np=2 on `development` before the fix: a callable
+metric demanding h = 0.01 within r < 0.15 of one corner and nothing elsewhere,
+on a 132-cell-per-rank base, hung with neither rank returning; `mpirun` reached
+its timeout.
+
+Run under MPI::
+
+    mpirun -np 2 python -m pytest --with-mpi \
+        tests/parallel/test_0873_adapt_collective_stop_mpi.py
+
+A callable metric is used deliberately, so that nothing here depends on how
+``global_evaluate`` behaves — it is evaluated rank-locally. #512 attributes the
+hazard to ``eval_metric`` being collective; measured at np=2 that is not so,
+with one rank asleep for 10 s the other's ``global_evaluate`` returned in
+0.06 s, for in-domain and for stranded points alike. The collective in this
+loop is the refinement, not the metric.
+"""
+
+import numpy as np
+import pytest
+
+import underworld3 as uw
+
+pytestmark = [
+    pytest.mark.mpi(min_size=2),
+    pytest.mark.timeout(300),
+    pytest.mark.level_2,
+    pytest.mark.tier_a,
+]
+
+ENGINES = ["sbr", "edge_split", "nvb"]
+
+
+def _corner_metric(points):
+    """Fine in one corner, already satisfied everywhere else.
+
+    This is what splits the ranks: whichever rank holds no corner cells has
+    nothing to refine at the first level, while its peer has work for several.
+    """
+
+    p = np.asarray(points)
+    r = np.sqrt((p[:, 0] - 0.02) ** 2 + (p[:, 1] - 0.02) ** 2)
+
+    return 1.0 / np.where(r < 0.15, 0.01, 1.0) ** 2
+
+
+def _base_mesh():
+    return uw.meshing.UnstructuredSimplexBox(
+        minCoords=(0.0, 0.0),
+        maxCoords=(1.0, 1.0),
+        cellSize=0.2,
+        regular=False,
+        qdegree=2,
+        refinement=1,
+    )
+
+
+def test_premise_the_metric_splits_the_ranks():
+    """Test premise: without a rank that runs out of work, nothing here is tested."""
+
+    mesh = _base_mesh()
+    cs, ce = mesh.dm.getHeightStratum(0)
+
+    centroid = np.asarray(
+        [mesh.dm.computeCellGeometryFVM(c)[1] for c in range(cs, ce)]
+    ).reshape(ce - cs, -1)[:, : mesh.dim]
+
+    wants_refining = int((_corner_metric(centroid) > 1.0).sum())
+    per_rank = uw.mpi.comm.allgather(wants_refining)
+
+    assert sum(per_rank) > 0, f"no rank has anything to refine: {per_rank}"
+    assert min(per_rank) == 0, (
+        f"every rank has work at np={uw.mpi.size} ({per_rank}) — this test "
+        "cannot detect the rank-local stop it exists for"
+    )
+
+
+@pytest.mark.parametrize("engine", ENGINES)
+def test_adapt_returns_when_one_rank_runs_out_of_work(engine):
+    """The whole assertion is that this returns. The failure mode is a hang."""
+
+    mesh = _base_mesh()
+
+    child = mesh.adapt(_corner_metric, max_levels=4, engine=engine)
+
+    c0, c1 = child.dm.getHeightStratum(0)
+    total = uw.mpi.comm.allreduce(c1 - c0)
+    base = uw.mpi.comm.allreduce(
+        mesh.dm.getHeightStratum(0)[1] - mesh.dm.getHeightStratum(0)[0]
+    )
+
+    assert total > base, (
+        f"{engine}: the corner metric demands refinement, but the child has "
+        f"{total} cells against the base's {base}"
+    )