Take slip normals per boundary, not from the deprecated Gamma_P1 (#538)

src/underworld3/meshing/smoothing/graph.py

@@ -527,33 +527,76 @@ def _min_incident_edge_nd(cells, coords):
 # ``boundary_slip`` orchestration, and ``BoundingSurface`` facet restore.
 # ---------------------------------------------------------------------------
 
-def _slip_normals(mesh, boundary_coords):
-    """Unit outward normals at ``boundary_coords`` from the projected
-    boundary-normal field.
-
-    Re-projects ``mesh._projected_normals`` (``mesh.Gamma_P1``) first so the
-    normals reflect the mesh's *current* coordinates — the projected field is
-    stale after any deform. Returns ``(normals, valid)`` where ``normals`` is
-    ``(k, cdim)`` and ``valid`` is a boolean mask; ``valid`` is ``False`` for
-    nodes with a degenerate (zero / non-finite) normal (e.g. box corners
-    where opposing face normals cancel, or an occasional unlocatable vertex).
-    Such nodes should be pinned, not slipped.
+#: Composite labels that cover every boundary at once. Summing a per-boundary
+#: normal over one of these would give a contribution at every boundary node and
+#: make each of them look like a corner, so they are never iterated.
+_COMPOSITE_BOUNDARY_LABELS = frozenset({"All_Boundaries"})
+
+
+def _slip_normals(mesh, boundary_coords, boundaries=None):
+    """Unit outward normals at ``boundary_coords``, taken PER BOUNDARY.
+
+    Returns ``(normals, valid)``. ``valid`` is ``False`` where no boundary claims
+    the node, where the normal is degenerate, and -- deliberately -- where MORE
+    THAN ONE boundary claims it. Those are corners and 3-D edges: they have no
+    single normal, so they are pinned rather than slipped along a fabricated one.
+
+    This used to evaluate ``mesh.Gamma_P1`` once. That field is deprecated and
+    documented as falling back to a coordinate direction off-kernel; on a unit
+    box it returns the edge TANGENT on Left and Bottom -- ``(0, 1)`` where the
+    outward normal is ``(-1, 0)``. A node whose "normal" lies along its own edge
+    has its tangential motion projected out, which is exactly why
+    ``node_redistribution(..., slip_surfaces=True)`` left those two edges frozen
+    at *exactly* zero displacement while Right and Top moved (#538). Per
+    boundary, ``_assemble_boundary_normal`` is exact on all four edges to 1.6e-11.
+
+    NO NEW MeshVariable IS CREATED. The per-boundary normal is assembled into
+    ``_n_proj`` -- the field ``Gamma_P1`` already owns -- one boundary at a time.
+    That is not tidiness: creating even ONE extra variable on this path
+    restructures the DM, and the parent's static hierarchy stops matching what it
+    recorded (``test_0764::test_redistribute_then_adapt_composition``). Verified
+    by keeping the creation and reverting to the old normals, which still failed;
+    the creation is the trigger, not the normals. ``_n_proj`` is rebuilt by
+    ``_update_projected_normals`` on every call, so overwriting it here is safe.
+
+    ``boundaries`` names the labels to consider; ``None`` uses every
+    non-composite boundary. Prefer passing the slip set: each label costs one
+    evaluation.
     """
     cdim = mesh.cdim
-    n = np.zeros((boundary_coords.shape[0], cdim))
+    count = boundary_coords.shape[0]
+    total = np.zeros((count, cdim))
+    claims = np.zeros(count, dtype=int)
+
+    if boundaries is None:
+        boundaries = [b.name for b in mesh.boundaries
+                      if b.name not in _COMPOSITE_BOUNDARY_LABELS]
+
     try:
-        mesh._update_projected_normals()
-        n = np.asarray(
-            uw.function.evaluate(mesh.Gamma_P1, boundary_coords)
-        ).reshape(-1, cdim)
+        mesh._update_projected_normals()          # creates/refreshes _n_proj
+        scratch = mesh._projected_normals
     except Exception:
-        # Projection unavailable / degenerate on this mesh — fall back to
-        # all-pinned boundaries (valid stays all-False below).
-        n = np.zeros((boundary_coords.shape[0], cdim))
-    mag = np.linalg.norm(n, axis=1)
-    valid = np.isfinite(mag) & (mag > 0.5)
-    out = np.zeros_like(n)
-    out[valid] = n[valid] / mag[valid, None]
+        return np.zeros((count, cdim)), np.zeros(count, dtype=bool)
+
+    for name in boundaries:
+        if name in _COMPOSITE_BOUNDARY_LABELS:
+            continue
+        try:
+            mesh._assemble_boundary_normal(scratch, name)
+            here = np.asarray(
+                uw.function.evaluate(scratch.sym, boundary_coords)
+            ).reshape(-1, cdim)
+        except Exception:
+            continue        # label absent on this rank / this mesh: claims nothing
+        magnitude = np.linalg.norm(here, axis=1)
+        on_it = np.isfinite(magnitude) & (magnitude > 0.5)
+        total[on_it] += here[on_it]
+        claims[on_it] += 1
+
+    magnitude = np.linalg.norm(total, axis=1)
+    valid = (claims == 1) & np.isfinite(magnitude) & (magnitude > 0.5)
+    out = np.zeros_like(total)
+    out[valid] = total[valid] / magnitude[valid, None]
     return out, valid
 
 

tests/test_0764_node_redistribution.py

@@ -222,3 +222,42 @@ def n_cells(m):
         assert n_cells(c_default) == n_cells(c_nvb)
         assert np.array_equal(np.asarray(c_default.X.coords),
                               np.asarray(c_nvb.X.coords))
+
+
+def test_slip_normals_are_outward_on_every_wall():
+    """Every wall gets its OWN outward normal, and corners are pinned (#538).
+
+    `_slip_normals` read the deprecated `mesh.Gamma_P1`, which off-kernel falls
+    back to a coordinate direction rather than a normal: on a unit box it
+    returned the edge TANGENT on Left and Bottom -- (0, 1) where the outward
+    normal is (-1, 0). A node whose "normal" lies along its own edge has its
+    tangential motion projected out, so `node_redistribution(slip_surfaces=True)`
+    froze those two walls at EXACTLY zero displacement while Right and Top moved
+    (left 0.00000 / right 0.18241, bottom 0.00000 / top 0.18427 on a metric band
+    placed symmetrically about the domain centre).
+
+    All four walls are asserted, not one: the defect was confined to two of them,
+    and a test on Right or Top alone passes throughout.
+    """
+    from underworld3.meshing.smoothing.graph import _slip_normals
+
+    mesh = uw.meshing.UnstructuredSimplexBox(
+        minCoords=(0.0, 0.0), maxCoords=(1.0, 1.0), cellSize=1.0 / 8, qdegree=3)
+    along = np.linspace(0.2, 0.8, 4)
+    zeros, ones = np.zeros_like(along), np.ones_like(along)
+
+    for wall, points, expected in (
+            ("Left", np.column_stack([zeros, along]), (-1.0, 0.0)),
+            ("Bottom", np.column_stack([along, zeros]), (0.0, -1.0)),
+            ("Right", np.column_stack([ones, along]), (1.0, 0.0)),
+            ("Top", np.column_stack([along, ones]), (0.0, 1.0))):
+        normals, valid = _slip_normals(mesh, points)
+        assert valid.all(), f"{wall}: {(~valid).sum()} of {len(points)} nodes unusable"
+        assert np.allclose(normals, np.array(expected), atol=1.0e-9), (
+            f"{wall}: got {normals[0]}, expected {expected} — a tangent here "
+            "means the deprecated Gamma_P1 path is back")
+
+    # A corner belongs to two walls and has no single normal: pin, do not slip.
+    corners = np.array([[0.0, 0.0], [1.0, 0.0], [0.0, 1.0], [1.0, 1.0]])
+    _normals, valid = _slip_normals(mesh, corners)
+    assert not valid.any(), "corners must be pinned, not slipped along one wall"