Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
75 changes: 70 additions & 5 deletions src/underworld3/meshing/surfaces.py
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down Expand Up @@ -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:
Expand Down Expand Up @@ -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.

Expand Down Expand Up @@ -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 ---

Expand All @@ -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
----------
Expand Down Expand Up @@ -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(
Expand Down
103 changes: 103 additions & 0 deletions tests/test_0851_surface_influence_edge.py
Original file line number Diff line number Diff line change
@@ -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}"
)
Loading