From 6291cab62b67c4b05548178428ee5e0410247f8e Mon Sep 17 00:00:00 2001 From: lmoresi Date: Tue, 16 Jun 2026 11:39:58 +1000 Subject: [PATCH] Surface.influence_function: respect finite edges (no bleed past the end) MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Surface.influence_function built its influence on `Abs(self.distance.sym[0])`, the Abs of the *signed* distance MeshVariable. The signed distance changes sign across the surface and its zero-contour follows the surface's infinite line/plane, so it extends PAST a finite edge. Any non-nodal evaluation interpolates the signed field through ~0 along that extension, and Abs() of that leaves a spurious high-influence streak well beyond the end of the surface (e.g. a weak fault zone bleeding past the fault tip). The per-node distances are already edge-clamped (the segment/surface distance falls back to the nearest endpoint), so the magnitude is correct at the nodes — the artifact is purely the interpolation of the sign-flipping field. Fix: store a companion UNSIGNED, edge-clamped distance field (`Surface.abs_distance`, populated as |signed distances| in _compute_distance_field) and build influence_function on it. The unsigned field is >= 0 everywhere and, beyond an edge, is just the radial distance to the endpoint, so its interpolant never crosses zero there. Adds tests/test_0851_surface_influence_edge.py: an embedded finite segment whose influence must decay beyond the tip — along the line extension the old basis spikes to ~1.0, the fix keeps it ~0; also checks abs_distance tracks the geometric edge-clamped distance between nodes. Underworld development team with AI support from Claude Code --- src/underworld3/meshing/surfaces.py | 75 ++++++++++++++-- tests/test_0851_surface_influence_edge.py | 103 ++++++++++++++++++++++ 2 files changed, 173 insertions(+), 5 deletions(-) create mode 100644 tests/test_0851_surface_influence_edge.py diff --git a/src/underworld3/meshing/surfaces.py b/src/underworld3/meshing/surfaces.py index 156dfc08..21191979 100644 --- a/src/underworld3/meshing/surfaces.py +++ b/src/underworld3/meshing/surfaces.py @@ -704,6 +704,10 @@ def __init__( # Level 3: Cached proxy MeshVariable for distance self._distance_var: Optional[uw.discretisation.MeshVariable] = None + # Companion UNSIGNED (edge-clamped) distance field. influence_function() + # interpolates THIS rather than Abs() of the signed field — see + # _compute_distance_field for why the signed field is unsafe to abs. + self._abs_distance_var: Optional[uw.discretisation.MeshVariable] = None # Dimension (2 or 3) - detected from mesh or control points self._dim = None @@ -759,6 +763,7 @@ def symbol(self, value: str) -> None: # If distance var exists, it needs to be recreated with new symbol if self._distance_var is not None: self._distance_var = None + self._abs_distance_var = None self._distance_stale = True def _dimensionalise_coords(self, coords: np.ndarray) -> np.ndarray: @@ -1133,6 +1138,37 @@ def distance(self) -> "uw.discretisation.MeshVariable": return self._distance_var + @property + def abs_distance(self) -> "uw.discretisation.MeshVariable": + """Unsigned (edge-clamped) distance from mesh nodes to the surface. + + Unlike ``Abs(self.distance.sym[0])``, this field is safe to interpolate + near a surface edge. The signed ``distance`` field changes sign across + the surface, and its zero-contour extends along the surface's + infinite-line/plane PAST the finite edge; interpolating it (any + non-nodal evaluation) and taking ``Abs`` produces a spurious near-zero + valley beyond the edge. This unsigned field is ``>= 0`` everywhere and, + beyond the edge, is simply the radial distance to the endpoint, so its + interpolant never crosses zero there. Used by :meth:`influence_function`. + + Returns: + MeshVariable with unsigned distance values at each mesh node. + Access ``.sym[0]`` for use in expressions. + """ + if self.mesh is None: + raise RuntimeError( + f"Surface '{self.name}' requires a mesh to compute distance field. " + "Set mesh in constructor or via surface.mesh = mesh" + ) + + self._ensure_discretized() + + if self._distance_stale: + self._compute_distance_field() + self._distance_stale = False + + return self._abs_distance_var + def _compute_distance_field(self) -> None: """Compute signed distance field from mesh nodes to surface. @@ -1187,8 +1223,31 @@ def _compute_distance_field(self) -> None: # Keep signed distance - helpers use sympy.Abs() when needed distances = dist_result.point_data["implicit_distance"] + # Companion UNSIGNED distance field. The per-node distance is already + # edge-clamped (the segment/surface distance falls back to the nearest + # ENDPOINT beyond the edge), so |distances| is the true, edge-aware + # distance to the FINITE surface at every node. We store it separately + # because the *signed* field is unsafe to interpolate near an edge: + # its zero-contour follows the surface's infinite-line/plane and so + # extends PAST the finite edge. A consumer that interpolates the signed + # field (any non-nodal evaluation) and then takes Abs() gets a spurious + # near-zero |distance| valley along that extension — i.e. the influence + # bleeds out beyond the end of the surface. The unsigned field is >= 0 + # and beyond the edge it is just the radial distance to the endpoint, + # so its interpolant never crosses zero there. influence_function() + # therefore uses this field, not Abs(signed). + if self._abs_distance_var is None: + self._abs_distance_var = uw.discretisation.MeshVariable( + f"surf_{self.name}_absdistance", + self.mesh, + 1, + degree=self.mesh.degree, + varsymbol=f"|d_{{{self._symbol}}}|", + ) + with uw.synchronised_array_update(): self._distance_var.data[:, 0] = distances + self._abs_distance_var.data[:, 0] = np.abs(distances) # --- Influence function --- @@ -1204,9 +1263,11 @@ def influence_function( Creates a sympy expression that varies from value_near (at the surface) to value_far (far from the surface) based on the chosen profile. - Uses the absolute value of the signed distance field, so the influence - is symmetric on both sides of the surface. For asymmetric behavior, - access the signed distance directly via ``surface.distance.sym[0]``. + Uses the unsigned, edge-clamped distance field (:attr:`abs_distance`), + so the influence is symmetric on both sides of the surface AND decays + correctly beyond a finite edge (it does not bleed along the surface's + line/plane past the end). For asymmetric behaviour, access the signed + distance directly via ``surface.distance.sym[0]``. Parameters ---------- @@ -1248,8 +1309,12 @@ def influence_function( # Accept quantities and convert to nondimensional mesh coordinates width = _to_nd_length(width) - # Use absolute distance - influence is symmetric about surface - d = sympy.Abs(self.distance.sym[0]) + # Use the UNSIGNED (edge-clamped) distance FIELD, not Abs() of the signed + # field. The signed field's zero-contour extends past a finite edge, so + # Abs() of its interpolant lights up a spurious weak zone beyond the + # surface end; the unsigned field is edge-aware between nodes too. See + # abs_distance / _compute_distance_field. + d = self.abs_distance.sym[0] if profile == "step": return sympy.Piecewise( diff --git a/tests/test_0851_surface_influence_edge.py b/tests/test_0851_surface_influence_edge.py new file mode 100644 index 00000000..a3cdb6f9 --- /dev/null +++ b/tests/test_0851_surface_influence_edge.py @@ -0,0 +1,103 @@ +#!/usr/bin/env python3 +"""Regression test: Surface.influence_function must respect a FINITE edge. + +A surface (here a short 2D line segment) has ends. The influence/weak zone +should decay beyond the segment tip — it must NOT bleed out along the segment's +infinite-line extension. + +The bug (fixed by using the unsigned, edge-clamped distance field rather than +Abs() of the *signed* field): the signed distance changes sign across the line +and its zero-contour runs along the infinite line past the finite tip, so any +non-nodal evaluation interpolates the signed field through ~0 there and Abs() +leaves a spurious high-influence streak beyond the end of the segment. + +We embed a short vertical segment and check the influence at NON-NODAL points +(which is where the interpolation artifact appears). +""" +import numpy as np +import pytest + +import underworld3 as uw + + +def _has(mod): + try: + __import__(mod) + return True + except ImportError: + return False + + +requires_pyvista = pytest.mark.skipif( + not _has("pyvista"), reason="pyvista required for surface distance" +) + +pytestmark = pytest.mark.level_2 + + +@requires_pyvista +def test_influence_decays_beyond_finite_edge(): + mesh = uw.meshing.UnstructuredSimplexBox( + minCoords=(0.0, 0.0), maxCoords=(1.0, 1.0), cellSize=0.05, qdegree=3 + ) + + # Short vertical segment, interior, ends at y = 0.5 (tip). + y = np.linspace(0.3, 0.5, 11) + pts = np.column_stack([np.full_like(y, 0.5), y, np.zeros_like(y)]) + surf = uw.meshing.Surface("seg", mesh, pts, symbol="S") + surf.discretize() + + width = 0.05 + influence = surf.influence_function( + width=width, value_near=1.0, value_far=0.0, profile="gaussian" + ) + + # On the segment (non-nodal): full influence. + v_on = float(np.asarray(uw.function.evaluate(influence, np.array([[0.5, 0.40]]))).reshape(-1)[0]) + assert v_on > 0.9, f"influence on segment should be ~1, got {v_on}" + + # Probe a BAND along the segment's line EXTENSION, well beyond the tip + # (y >= 0.65 is >= 0.15 past the tip at y=0.5, i.e. >= 3*width, so the true + # influence is < 1e-3 everywhere here). The bug spikes the influence to ~1 + # on the mesh edges that cross the extension; the fix keeps it ~0. Use a + # band in x (not just x=0.5) so we actually hit those crossing edges. + xs = np.linspace(0.45, 0.55, 11) + ys = np.linspace(0.65, 0.92, 28) + band = np.array([[x, y] for y in ys for x in xs]) + v_band = np.asarray(uw.function.evaluate(influence, band)).reshape(-1) + assert v_band.max() < 0.05, ( + f"influence bled past the finite edge: max={v_band.max():.3f} in the " + f"extension band >=0.15 beyond the tip (width={width}); expected ~0" + ) + + # Sanity: a perpendicular-far point also decays. + v_perp = float(np.asarray(uw.function.evaluate(influence, np.array([[0.85, 0.40]]))).reshape(-1)[0]) + assert v_perp < 0.05, f"influence far perpendicular should be ~0, got {v_perp}" + + +@requires_pyvista +def test_abs_distance_matches_geometry_beyond_edge(): + """The unsigned distance field is edge-aware between nodes too.""" + from underworld3.utilities.geometry_tools import ( + signed_distance_pointcloud_polyline_2d, + ) + + mesh = uw.meshing.UnstructuredSimplexBox( + minCoords=(0.0, 0.0), maxCoords=(1.0, 1.0), cellSize=0.05, qdegree=3 + ) + y = np.linspace(0.3, 0.5, 11) + xy = np.column_stack([np.full_like(y, 0.5), y]) + surf = uw.meshing.Surface("seg2", mesh, np.column_stack([xy, np.zeros(len(xy))]), symbol="S") + surf.discretize() + + probe = np.array([[0.5, 0.80], [0.5, 0.40], [0.85, 0.40], [0.65, 0.65]]) + d_field = np.asarray( + uw.function.evaluate(surf.abs_distance.sym[0], probe) + ).reshape(-1) + d_true = np.abs(signed_distance_pointcloud_polyline_2d(probe, xy)) + + # Interpolated unsigned field tracks the true edge-clamped distance + # (linear interp of a >=0 cone, so within a cell-size tolerance). + assert np.allclose(d_field, d_true, atol=0.05), ( + f"abs_distance field {d_field} != geometric {d_true}" + )