From f6b609a4f6441251479e0833855d2bddbc5ba221 Mon Sep 17 00:00:00 2001
From: lmoresi <louis.moresi@anu.edu.au>
Date: Sun, 16 Aug 2026 08:56:27 +1000
Subject: [PATCH 1/4] Analytic solutions declare their own boundary conditions
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

The suite carried two mixins, FreeSlipWalls and FixedWalls, each looping
over every boundary of a solution and applying one condition to all of
them. That encodes "every boundary is a wall of the same kind" — true of
the classical box benchmarks, false of what the suite has to serve next:
a spherical shell or annulus with different conditions on the two radii,
a faulted disc, a channel driven at one end. FreeSlipWalls would have had
to grow a parameter for each of those, and inheritance advertises the
choice as if it were part of what the solution is.

Both mixins are gone. Every solution now writes apply_boundary_conditions
itself and composes three module-level helpers, each taking the
boundaries it applies to rather than reading self.boundaries:

    free_slip(solver, boundaries, normal=None)
    prescribed_velocity(solver, boundaries, velocity)
    prescribed_scalar(solver, boundaries, field)

A solution needing two kinds of condition calls two of them. The scalar
helper also removes the third copy of the same Dirichlet loop, which the
transport and Richards families each carried.

The pressure nullspace is now stated by the solution rather than set as a
side effect of a wall type. It is a property of the domain — enclosed, so
the pressure is fixed only up to a constant — and not of any one
boundary's condition, and two different mixins both setting it hid that.

free_slip takes normal= so a curved-geometry solution can ask for an
analytic normal. CylindricalStokes, the one curved case here, keeps the
geometric default deliberately: it is the direction the straight-facet
boundary integral actually sees, and an analytic X/|X| would trade a
machine-precision pressure gauge for a consistency error that grows with
facet non-uniformity (#561).

This is a pure refactor and it was verified rather than asserted. A
recording mock solver captured the ordered sequence of calls and
attribute sets that apply_boundary_conditions produces for every
registered solution — sympy arguments by srepr, floats by their exact hex
form — on the merge base and on this branch. The two recordings are
byte-identical across all 21 cases (20 registered solutions plus both
CylindricalStokes boundary cases).

The equivalence check did expose one latent defect, left in place and
marked with a TODO(BUG): CylindricalStokes with boundary="zero" does not
remove the pressure nullspace, where every other enclosed case does. An
annulus held at zero velocity on both arcs is enclosed, so its pressure
is determined only up to a constant. Fixing it would change behaviour and
belongs under its own test.

Underworld development team with AI support from Claude Code
---
 src/underworld3/analytic/__init__.py  |  12 ++-
 src/underworld3/analytic/_base.py     | 119 ++++++++++++++++++++------
 src/underworld3/analytic/inclusion.py |  16 +++-
 src/underworld3/analytic/kramer.py    |  28 ++++--
 src/underworld3/analytic/richards.py  |  12 ++-
 src/underworld3/analytic/transport.py |  18 ++--
 src/underworld3/analytic/velic.py     |  86 ++++++++++++++++---
 tests/test_1016_analytic_contract.py  |  83 ++++++++++++++----
 8 files changed, 296 insertions(+), 78 deletions(-)

diff --git a/src/underworld3/analytic/__init__.py b/src/underworld3/analytic/__init__.py
index f7a45858..6d7690c9 100644
--- a/src/underworld3/analytic/__init__.py
+++ b/src/underworld3/analytic/__init__.py
@@ -26,7 +26,12 @@
 underworld3.analytic._base : the contract each solution satisfies.
 """
 
-from ._base import AnalyticSolution, FreeSlipWalls, FixedWalls
+from ._base import (
+    AnalyticSolution,
+    free_slip,
+    prescribed_scalar,
+    prescribed_velocity,
+)
 
 from .inclusion import EllipticalInclusion
 from .kramer import CylindricalStokes
@@ -36,8 +41,9 @@
 
 __all__ = [
     "AnalyticSolution",
-    "FreeSlipWalls",
-    "FixedWalls",
+    "free_slip",
+    "prescribed_scalar",
+    "prescribed_velocity",
     "AdvectedFront",
     "EllipticalInclusion",
     "ErfcDiffusion",
diff --git a/src/underworld3/analytic/_base.py b/src/underworld3/analytic/_base.py
index 8bd5ee77..b029e9eb 100644
--- a/src/underworld3/analytic/_base.py
+++ b/src/underworld3/analytic/_base.py
@@ -522,11 +522,20 @@ def set_richards_field(self, solution, conductivity, capacity):
         self.fn_capacity = sympy.sympify(capacity)
 
     def apply_boundary_conditions(self, solver):
-        """Impose the boundary conditions this solution is posed under."""
+        """Impose the boundary conditions this solution is posed under.
+
+        Every solution writes this itself, because the conditions are part of
+        the problem it answers and are not always the same on every boundary.
+        :func:`free_slip`, :func:`prescribed_velocity` and
+        :func:`prescribed_scalar` are here to be composed inside it — each takes
+        the boundaries it applies to, so a solution with a slipping wall and a
+        driven one calls both.
+        """
 
         raise NotImplementedError(
             f"{type(self).__name__} does not declare its boundary conditions. "
-            f"Mix in FreeSlipWalls or FixedWalls, or override this method."
+            f"Write apply_boundary_conditions, composing free_slip, "
+            f"prescribed_velocity or prescribed_scalar as the problem requires."
         )
 
     def _object_viewer(self):
@@ -547,39 +556,93 @@ def _object_viewer(self):
                 display(Markdown(rf"{label}: $\displaystyle {equation}$"))
 
 
-class FreeSlipWalls:
-    r"""Mixin: the solution is posed with free slip on every wall.
+# The boundary conditions solutions build from. Functions rather than mixins:
+# a mixin can only say "every boundary of this solution gets the same
+# condition", which is true of the classical box benchmarks and false of most
+# of what the suite has to serve next — a spherical shell or annulus with
+# different conditions on the two radii, a faulted disc, a channel driven at one
+# end. Composed inside a solution's own `apply_boundary_conditions`, these say
+# what is imposed WHERE, and a solution that needs two kinds calls two of them.
+#
+# The pressure nullspace is deliberately NOT set here. It is a property of the
+# domain — enclosed, so the pressure is determined only up to a constant — and
+# not of any one wall's condition, so each solution states it for itself:
+#
+#     solver.petsc_use_pressure_nullspace = True
+#
+# Leaving it out on an enclosed domain is the failure this whole suite exists to
+# catch: a direct solve on the singular saddle returns a quiet, wrong answer
+# that only an exact solution exposes.
+
+
+def free_slip(solver, boundaries, normal=None):
+    r"""Impose free slip, :math:`\mathbf u\cdot\hat{\mathbf n}=0`, on *boundaries*.
+
+    Applied as a strong *rotated* constraint (satisfied to machine precision)
+    rather than by zeroing a velocity component. On an axis-aligned box the two
+    agree; on a curved, tilted or adapted boundary only the rotated form is
+    correct, so it is the one that still holds when a solution is used to
+    validate an adapted mesh.
+
+    Parameters
+    ----------
+    solver : Stokes
+        The solver to configure.
+    boundaries : sequence of str
+        The boundary labels to constrain — not necessarily all of them.
+    normal : sympy 1 x dim Matrix, optional
+        Analytic outward normal. The default takes the solver's geometric facet
+        normal, which is consistent with the straight-facet integral the
+        assembler evaluates and is normally the right choice even on a curved
+        wall. Supply one only when the constraint must follow the TRUE surface
+        rather than the mesh; see the "Which normal to use" section of
+        ``docs/developer/subsystems/rotated-freeslip.md``.
+    """
+
+    # Passed only when the solution has actually chosen a normal, so that
+    # otherwise the solver's own geometric default governs.
+    choice = {} if normal is None else {"normal": normal}
+
+    for boundary in boundaries:
+        solver.add_rotated_freeslip_bc(0.0, boundary, **choice)
+
 
-    Free slip is imposed as a strong *rotated* constraint
-    (:math:`\mathbf u\cdot\hat{\mathbf n}=0` to machine precision) rather than by
-    zeroing a velocity component. On an axis-aligned box the two agree; on a
-    curved, tilted or adapted boundary only the rotated form is correct, so it is
-    the one that still holds when a solution is used to validate an adapted mesh.
+def prescribed_velocity(solver, boundaries, velocity):
+    """Impose a velocity on *boundaries*, component by component.
 
-    The domain is enclosed, so the pressure carries a constant nullspace and the
-    solver is told to remove it. Leaving that out is the failure this whole suite
-    exists to catch: a direct solve on the singular saddle returns a quiet, wrong
-    answer that only an exact solution exposes.
+    For a solution driven through its boundaries rather than by a body force —
+    a far-field shear, say — and for any solution whose exact velocity is not
+    tangential to the domain.
+
+    Parameters
+    ----------
+    solver : Stokes
+        The solver to configure.
+    boundaries : sequence of str
+        The boundary labels to constrain.
+    velocity : sympy Matrix, shape (1, dim)
+        Usually the solution's own ``fn_velocity``.
     """
 
-    def apply_boundary_conditions(self, solver):
-        for boundary in self.boundaries:
-            solver.add_rotated_freeslip_bc(0.0, boundary)
+    for boundary in boundaries:
+        solver.add_dirichlet_bc(velocity, boundary)
 
-        solver.petsc_use_pressure_nullspace = True
 
+def prescribed_scalar(solver, boundaries, field):
+    """Impose a scalar field on *boundaries*.
 
-class FixedWalls:
-    """Mixin: velocity is prescribed on every wall, from the exact solution.
+    The scalar counterpart of :func:`prescribed_velocity`, for the transport and
+    Richards families: temperature, pressure head, hydraulic head.
 
-    For solutions driven by their boundaries rather than by a body force — a
-    far-field shear, say — and for manufactured solutions whose exact velocity is
-    not tangential to the domain. The domain is again enclosed, so the pressure
-    nullspace is removed.
+    Parameters
+    ----------
+    solver : Poisson or Darcy
+        The solver to configure.
+    boundaries : sequence of str
+        The boundary labels to constrain.
+    field : sympy expression
+        Usually the solution's own ``fn_solution``.
     """
 
-    def apply_boundary_conditions(self, solver):
-        for boundary in self.boundaries:
-            solver.add_dirichlet_bc(self.fn_velocity, boundary)
-
-        solver.petsc_use_pressure_nullspace = True
+    for boundary in boundaries:
+        solver.add_dirichlet_bc([field], boundary)
diff --git a/src/underworld3/analytic/inclusion.py b/src/underworld3/analytic/inclusion.py
index b21500a1..515de208 100644
--- a/src/underworld3/analytic/inclusion.py
+++ b/src/underworld3/analytic/inclusion.py
@@ -36,7 +36,7 @@
 
 import sympy
 
-from ._base import AnalyticSolution, FixedWalls
+from ._base import AnalyticSolution, prescribed_velocity
 
 
 def _conjugate(expression):
@@ -172,7 +172,7 @@ def _potentials(zeta, viscosity_ratio, aspect_ratio, alpha, pure_shear, simple_s
     }
 
 
-class EllipticalInclusion(FixedWalls, AnalyticSolution):
+class EllipticalInclusion(AnalyticSolution):
     r"""A viscous elliptical inclusion in a matrix under general shear.
 
     Parameters
@@ -376,6 +376,18 @@ def __init__(
         self._potentials = potentials
         self._centre = tuple(float(c) for c in centre)
 
+    def apply_boundary_conditions(self, solver):
+        """The exact velocity on every wall; the enclosed box has a pressure nullspace.
+
+        The flow is driven by the far field, so the domain boundary is where that
+        far field enters: the walls carry the exact solution rather than a slip
+        condition, and cutting the box closer to the inclusion changes only how
+        much of the analytic far field is imposed rather than solved for.
+        """
+
+        prescribed_velocity(solver, self.boundaries, self.fn_velocity)
+        solver.petsc_use_pressure_nullspace = True
+
     @property
     def semi_axes(self):
         """Long and short semi-axes of the inclusion, in physical units."""
diff --git a/src/underworld3/analytic/kramer.py b/src/underworld3/analytic/kramer.py
index 620fbe3c..f386170b 100644
--- a/src/underworld3/analytic/kramer.py
+++ b/src/underworld3/analytic/kramer.py
@@ -31,7 +31,7 @@
 
 import numpy as np
 
-from ._base import AnalyticSolution
+from ._base import AnalyticSolution, free_slip, prescribed_velocity
 
 
 _INSTALL_MESSAGE = (
@@ -254,14 +254,28 @@ def error(self, field, meshvar, norm="l2"):
         return float(np.sqrt(difference / magnitude))
 
     def apply_boundary_conditions(self, solver):
-        """Free-slip or zero-slip on both walls, as the case declares."""
+        """Free-slip or zero-slip on both arcs, as the case declares.
+
+        The two boundaries are curved, and the constraint uses the solver's own
+        geometric normal rather than the analytic radial one. That is deliberate:
+        the geometric normal is the direction the straight-facet boundary
+        integral actually sees, so the constant pressure stays a null vector to
+        machine precision, where ``X/|X|`` leaves a consistency error that grows
+        with facet non-uniformity. Reach for ``normal=`` only when the constraint
+        must follow the true circle rather than the mesh — see the "Which normal
+        to use" section of ``docs/developer/subsystems/rotated-freeslip.md``.
+        """
 
         if self.boundary == "zero":
-            for wall in self.boundaries:
-                solver.add_dirichlet_bc((0.0, 0.0), wall)
+            # TODO(BUG): the zero-slip case leaves the pressure nullspace in
+            # place, where every other enclosed case in the suite removes it. An
+            # annulus with both arcs held at zero velocity is enclosed, so its
+            # pressure is determined only up to a constant and a direct solve on
+            # the singular saddle can return a quiet, wrong answer. Preserved
+            # here because this refactor is behaviour-preserving by contract;
+            # fix it under its own test.
+            prescribed_velocity(solver, self.boundaries, (0.0, 0.0))
             return
 
-        for wall in self.boundaries:
-            solver.add_rotated_freeslip_bc(0.0, wall)
-
+        free_slip(solver, self.boundaries)
         solver.petsc_use_pressure_nullspace = True
diff --git a/src/underworld3/analytic/richards.py b/src/underworld3/analytic/richards.py
index 2556e70f..af2cb9ab 100644
--- a/src/underworld3/analytic/richards.py
+++ b/src/underworld3/analytic/richards.py
@@ -26,7 +26,7 @@
 
 import sympy
 
-from ._base import AnalyticSolution
+from ._base import AnalyticSolution, prescribed_scalar
 
 
 def gardner_steady_saturation(y, psi_0, psi_L, L, alpha):
@@ -91,10 +91,14 @@ class _Gardner(AnalyticSolution):
     )
 
     def apply_boundary_conditions(self, solver):
-        """Prescribe the exact head on every wall."""
+        """Prescribe the exact head on every wall.
 
-        for boundary in self.boundaries:
-            solver.add_dirichlet_bc([self.fn_solution], boundary)
+        Both Gardner solutions are column problems posed with their own head on
+        the boundary. A member needing a flux condition on the base instead
+        writes its own method and composes the helpers per boundary.
+        """
+
+        prescribed_scalar(solver, self.boundaries, self.fn_solution)
 
     def _gardner_material(self, head):
         r"""Conductivity and capacity for a given head expression.
diff --git a/src/underworld3/analytic/transport.py b/src/underworld3/analytic/transport.py
index ce2c0054..15304440 100644
--- a/src/underworld3/analytic/transport.py
+++ b/src/underworld3/analytic/transport.py
@@ -16,25 +16,29 @@
 
 import sympy
 
-from ._base import AnalyticSolution
+from ._base import AnalyticSolution, prescribed_scalar
 
 
 class _Transport(AnalyticSolution):
     """Shared metadata and boundary conditions for the scalar solutions.
 
-    These prescribe the field itself on the boundary, so they cannot reuse the
-    Stokes mixins — `FixedWalls` applies a velocity, which a scalar solution does
-    not have.
+    The whole family is posed with its own exact field on every wall, which is
+    why the condition is stated once here. A member that needs something else on
+    one boundary — a no-flux base, an inflow face — writes its own
+    ``apply_boundary_conditions`` and composes the helpers per boundary.
     """
 
     solves = "transport"
     dim = 2
 
     def apply_boundary_conditions(self, solver):
-        """Prescribe the exact field on every wall."""
+        """Prescribe the exact field on every wall.
 
-        for boundary in self.boundaries:
-            solver.add_dirichlet_bc([self.fn_solution], boundary)
+        There is no pressure nullspace to remove: the field is pinned by its own
+        Dirichlet data everywhere on the boundary.
+        """
+
+        prescribed_scalar(solver, self.boundaries, self.fn_solution)
 
 
 class Poisson1D(_Transport):
diff --git a/src/underworld3/analytic/velic.py b/src/underworld3/analytic/velic.py
index c97205e1..cc0e9ec3 100644
--- a/src/underworld3/analytic/velic.py
+++ b/src/underworld3/analytic/velic.py
@@ -17,7 +17,7 @@
 
 import sympy
 
-from ._base import AnalyticSolution, FixedWalls, FreeSlipWalls
+from ._base import AnalyticSolution, free_slip, prescribed_velocity
 from ._transcribe import (
     CSource,
     evaluate_block,
@@ -113,7 +113,7 @@ def _solcx_kernel(variant="_solCx_A"):
     }
 
 
-class SolCx(FreeSlipWalls, AnalyticSolution):
+class SolCx(AnalyticSolution):
     r"""Isoviscous-column Stokes flow with a viscosity step — the SolCx benchmark.
 
     Viscosity jumps from :math:`\eta_A` to :math:`\eta_B` at :math:`x = x_c` on
@@ -234,6 +234,12 @@ def __init__(self, mesh, eta_A=1.0, eta_B=1.0e6, x_c=0.5, n=1, reference=False):
         if reference:
             self._use_reference_kernel()
 
+    def apply_boundary_conditions(self, solver):
+        """Free slip on all four walls; the enclosed box has a pressure nullspace."""
+
+        free_slip(solver, self.boundaries)
+        solver.petsc_use_pressure_nullspace = True
+
     def velocity_error(self, velocity_var):
         """Global relative L2 velocity error. Equivalent to ``error("velocity", ...)``."""
 
@@ -344,7 +350,7 @@ def block(name, signature):
     }
 
 
-class SolNL(FixedWalls, AnalyticSolution):
+class SolNL(AnalyticSolution):
     r"""Power-law viscous flow — the SolNL nonlinear benchmark.
 
     A manufactured solution for a shear-thinning fluid: the viscosity depends on
@@ -447,6 +453,12 @@ def __init__(self, mesh, eta_0=1.0, n=1, r=1.5, reference=False):
         if reference:
             self._use_reference_kernel()
 
+    def apply_boundary_conditions(self, solver):
+        """The exact velocity on all four walls — SolNL is driven by its boundaries."""
+
+        prescribed_velocity(solver, self.boundaries, self.fn_velocity)
+        solver.petsc_use_pressure_nullspace = True
+
     def _use_reference_kernel(self):
         """Point-evaluation only: opaque to the JIT, so no solver can use it."""
 
@@ -523,7 +535,7 @@ def _solkx_kernel():
     }
 
 
-class SolKx(FreeSlipWalls, AnalyticSolution):
+class SolKx(AnalyticSolution):
     r"""Stokes flow with an exponentially varying viscosity — the SolKx benchmark.
 
     Viscosity :math:`\eta = e^{2Bx}` on the unit box, driven by the density
@@ -624,6 +636,12 @@ def __init__(self, mesh, B=2.302585092994046, n=3, m=2):
             ),
         )
 
+    def apply_boundary_conditions(self, solver):
+        """Free slip on all four walls; the enclosed box has a pressure nullspace."""
+
+        free_slip(solver, self.boundaries)
+        solver.petsc_use_pressure_nullspace = True
+
 
 _Y = sympy.Symbol("y")
 _BETA = sympy.Symbol("Beta")
@@ -677,11 +695,21 @@ def block(method):
     return fields
 
 
-class _SolDB(FixedWalls, AnalyticSolution):
+class _SolDB(AnalyticSolution):
     """Shared assembly for the Dohrmann–Bochev manufactured solutions."""
 
     stress_is_deviatoric = True
 
+    def apply_boundary_conditions(self, solver):
+        """The exact velocity on every wall; the enclosed box has a pressure nullspace.
+
+        These are manufactured solutions whose velocity is not tangential to the
+        box, so free slip would be a different problem.
+        """
+
+        prescribed_velocity(solver, self.boundaries, self.fn_velocity)
+        solver.petsc_use_pressure_nullspace = True
+
     def _assemble(self, mesh, values, names):
         kernel = {
             field: expression.subs(values)
@@ -826,7 +854,7 @@ def _solkz_kernel():
     }
 
 
-class SolKz(FreeSlipWalls, AnalyticSolution):
+class SolKz(AnalyticSolution):
     r"""Stokes flow with a depth-dependent viscosity — the SolKz benchmark.
 
     Viscosity :math:`\eta = e^{2Bz}` on the unit box, free slip everywhere,
@@ -909,6 +937,12 @@ def __init__(self, mesh, B=2.302585092994046, n=3, m=2):
             ),
         )
 
+    def apply_boundary_conditions(self, solver):
+        """Free slip on all four walls; the enclosed box has a pressure nullspace."""
+
+        free_slip(solver, self.boundaries)
+        solver.petsc_use_pressure_nullspace = True
+
 
 _SIGMA = sympy.Symbol("sigma")
 
@@ -950,7 +984,7 @@ def _solab_kernel(name):
     }
 
 
-class _SolAB(FreeSlipWalls, AnalyticSolution):
+class _SolAB(AnalyticSolution):
     """Shared assembly for the two isoviscous Velic solutions."""
 
     dim = 2
@@ -1032,6 +1066,12 @@ def __init__(self, mesh, sigma=1.0, eta=1.0, n=3, m=2):
             ),
         )
 
+    def apply_boundary_conditions(self, solver):
+        """Free slip on all four walls; the enclosed box has a pressure nullspace."""
+
+        free_slip(solver, self.boundaries)
+        solver.petsc_use_pressure_nullspace = True
+
 
 class SolA(_SolAB):
     r"""Isoviscous Stokes flow with a sinusoidal body force — the SolA benchmark.
@@ -1157,7 +1197,7 @@ def block(method):
     }
 
 
-class SolM(FreeSlipWalls, AnalyticSolution):
+class SolM(AnalyticSolution):
     r"""Stokes flow with a laterally oscillating viscosity — the SolM benchmark.
 
     Viscosity :math:`\eta = 1 + \eta_0(1 + \cos(r\pi x))` on the unit box, free
@@ -1250,6 +1290,12 @@ def __init__(self, mesh, eta_0=1.0, n=3, m=2, r=4.0):
             ),
         )
 
+    def apply_boundary_conditions(self, solver):
+        """Free slip on all four walls; the enclosed box has a pressure nullspace."""
+
+        free_slip(solver, self.boundaries)
+        solver.petsc_use_pressure_nullspace = True
+
 
 _XC_C = sympy.Symbol("xc_solc")
 
@@ -1321,7 +1367,7 @@ def _solc_kernel(modes):
     return totals
 
 
-class SolC(FreeSlipWalls, AnalyticSolution):
+class SolC(AnalyticSolution):
     r"""Isoviscous flow driven by a dense column — the SolC benchmark.
 
     Constant viscosity on the unit box, free slip everywhere, driven by a density
@@ -1410,6 +1456,12 @@ def __init__(self, mesh, sigma=1.0, eta=1.0, x_c=0.5, modes=40):
             ),
         )
 
+    def apply_boundary_conditions(self, solver):
+        """Free slip on all four walls; the enclosed box has a pressure nullspace."""
+
+        free_slip(solver, self.boundaries)
+        solver.petsc_use_pressure_nullspace = True
+
 
 _ZC, _DX, _X0 = sympy.symbols("zc dx x0")
 
@@ -1508,7 +1560,7 @@ def across(symbol):
     return totals
 
 
-class SolDA(FreeSlipWalls, AnalyticSolution):
+class SolDA(AnalyticSolution):
     r"""A dense column in a layered fluid — the SolDA benchmark.
 
     A rectangular density anomaly of width *dx* centred at *x_0*, in a fluid
@@ -1620,6 +1672,12 @@ def __init__(
             ),
         )
 
+    def apply_boundary_conditions(self, solver):
+        """Free slip on all four walls; the enclosed box has a pressure nullspace."""
+
+        free_slip(solver, self.boundaries)
+        solver.petsc_use_pressure_nullspace = True
+
 
 _DY = sympy.Symbol("dy")
 
@@ -1695,7 +1753,7 @@ def _solh_kernel(modes):
     return totals
 
 
-class SolH(FreeSlipWalls, AnalyticSolution):
+class SolH(AnalyticSolution):
     r"""A dense block in three dimensions — the SolH benchmark.
 
     Isoviscous flow in the unit cube, free slip everywhere, driven by a
@@ -1783,3 +1841,9 @@ def __init__(self, mesh, sigma=1.0, eta=1.0, dx=0.5, dy=0.5, modes=6):
                 (kernel["stress_xz"], kernel["stress_yz"], kernel["stress_zz"]),
             ),
         )
+
+    def apply_boundary_conditions(self, solver):
+        """Free slip on all six faces; the enclosed cube has a pressure nullspace."""
+
+        free_slip(solver, self.boundaries)
+        solver.petsc_use_pressure_nullspace = True
diff --git a/tests/test_1016_analytic_contract.py b/tests/test_1016_analytic_contract.py
index 6b714c0b..1ca28dbd 100644
--- a/tests/test_1016_analytic_contract.py
+++ b/tests/test_1016_analytic_contract.py
@@ -2,9 +2,9 @@
 
 `uw.analytic` is the namespace for exact solutions used to validate a solve.
 These tests fix the contract itself — that the base class exposes the fields and
-error norms every solution promises, that the boundary-condition mixins configure
-a solver rather than returning something the caller has to apply, and that a
-solution reached through the new namespace is the *same object* as the one
+error norms every solution promises, that the boundary-condition helpers
+configure a solver rather than returning something the caller has to apply, and
+that a solution reached through the new namespace is the *same object* as the one
 reached through the old one.
 
 That last check matters more than it looks: `uw.function.analytic` is a compiled
@@ -31,7 +31,7 @@ def mesh():
     )
 
 
-class Quadratic(uw.analytic.FixedWalls, uw.analytic.AnalyticSolution):
+class Quadratic(uw.analytic.AnalyticSolution):
     r"""A manufactured Stokes solution used to exercise the contract.
 
     Divergence-free by construction: :math:`\mathbf u = (y, -x)` is a rigid
@@ -56,10 +56,23 @@ def __init__(self, mesh):
         self.fn_strainrate = sympy.Matrix([[0, 0], [0, 0]])
         self.fn_stress = sympy.Matrix([[0, 0], [0, 0]])
 
+    def apply_boundary_conditions(self, solver):
+        uw.analytic.prescribed_velocity(solver, self.boundaries, self.fn_velocity)
+        solver.petsc_use_pressure_nullspace = True
+
 
 def test_contract_is_exported():
-    """The contract and the mixins are reachable from the package namespace."""
-    for name in ("AnalyticSolution", "FreeSlipWalls", "FixedWalls"):
+    """The contract and the BC helpers are reachable from the package namespace.
+
+    A solution written outside this package has to be able to reach the helpers
+    it composes, so they are public API rather than a private convenience.
+    """
+    for name in (
+        "AnalyticSolution",
+        "free_slip",
+        "prescribed_scalar",
+        "prescribed_velocity",
+    ):
         assert name in uw.analytic.__all__
         assert hasattr(uw.analytic, name)
 
@@ -266,8 +279,8 @@ class NoBoundaryConditions(uw.analytic.AnalyticSolution):
         sol.apply_boundary_conditions(stokes)
 
 
-def test_fixed_walls_configures_the_solver(mesh):
-    """FixedWalls imposes the exact velocity on every wall and kills the nullspace."""
+def test_prescribed_velocity_configures_the_solver(mesh):
+    """The velocity helper imposes the exact velocity on the boundaries it is given."""
     sol = Quadratic(mesh)
     stokes = uw.systems.Stokes(mesh)
 
@@ -277,8 +290,8 @@ def test_fixed_walls_configures_the_solver(mesh):
     assert {bc.boundary for bc in stokes.essential_bcs} == set(sol.boundaries)
 
 
-def test_free_slip_walls_uses_the_rotated_constraint(mesh):
-    """FreeSlipWalls imposes u.n = 0 by rotation, not by masking a component.
+def test_free_slip_uses_the_rotated_constraint(mesh):
+    """free_slip imposes u.n = 0 by rotation, not by masking a component.
 
     Component masking is only equivalent on an axis-aligned box; the rotated form
     is what still holds when a solution is used to validate a curved or adapted
@@ -286,20 +299,58 @@ def test_free_slip_walls_uses_the_rotated_constraint(mesh):
     is the only signal that distinguishes the two paths.
     """
 
-    class FreeSlip(uw.analytic.FreeSlipWalls, Quadratic):
-        pass
-
-    sol = FreeSlip(mesh)
+    sol = Quadratic(mesh)
     stokes = uw.systems.Stokes(mesh)
 
-    sol.apply_boundary_conditions(stokes)
+    uw.analytic.free_slip(stokes, sol.boundaries)
 
-    assert stokes.petsc_use_pressure_nullspace
     registered = {boundary for boundary, _ in stokes._rotated_freeslip_bcs}
     assert registered == set(sol.boundaries)
     assert stokes.essential_bcs == []
 
 
+def test_boundaries_can_carry_different_conditions(mesh):
+    """The reason the helpers take a boundary list rather than reading self.
+
+    A shell, an annulus or a channel is not "walls of one kind": the condition
+    differs from one boundary to the next. A solution says so by calling more
+    than one helper, each with the boundaries it applies to.
+    """
+
+    sol = Quadratic(mesh)
+    stokes = uw.systems.Stokes(mesh)
+
+    uw.analytic.free_slip(stokes, ["Left", "Right"])
+    uw.analytic.prescribed_velocity(stokes, ["Bottom", "Top"], sol.fn_velocity)
+
+    slipping = {boundary for boundary, _ in stokes._rotated_freeslip_bcs}
+    driven = {bc.boundary for bc in stokes.essential_bcs}
+
+    assert slipping == {"Left", "Right"}
+    assert driven == {"Bottom", "Top"}
+
+
+def test_free_slip_passes_an_analytic_normal_through(mesh):
+    """A curved-boundary solution can choose the normal; the default is geometric.
+
+    Both calls must register, and the analytic one must carry the normal it was
+    given — otherwise a spherical or annulus solution silently gets the facet
+    normal it explicitly declined.
+    """
+
+    stokes = uw.systems.Stokes(mesh)
+    x, y = mesh.X
+    radial = sympy.Matrix([[x, y]]) / sympy.sqrt(x**2 + y**2)
+
+    uw.analytic.free_slip(stokes, ["Left"])
+    uw.analytic.free_slip(stokes, ["Right"], normal=radial)
+
+    registered = dict(stokes._rotated_freeslip_bcs)
+
+    assert registered["Left"] is None
+    assert registered["Right"] == radial
+
+
 def test_boundaries_follow_the_mesh_dimension(mesh):
     sol = Quadratic(mesh)
     assert sol.boundaries == ["Left", "Right", "Bottom", "Top"]

From 4643069e119d2e3cf21a06a35c1d47c230a18bc5 Mon Sep 17 00:00:00 2001
From: lmoresi <louis.moresi@anu.edu.au>
Date: Sun, 16 Aug 2026 08:56:56 +1000
Subject: [PATCH 2/4] Document how an analytic solution states its boundary
 conditions

The subsystem note described the two mixins. Replace that with the
composition the code now uses: a section on the three helpers, what each
takes, why functions rather than mixins, why the pressure nullspace is
the solution's statement and not a wall type's side effect, and which
normal to use on a curved boundary.

Also corrects the one stale cross-reference left in the conformance
suite's comments, which still pointed at "the mixins".

Underworld development team with AI support from Claude Code
---
 .../subsystems/analytic-solutions.md          | 67 +++++++++++++++++--
 tests/test_1024_analytic_conformance.py       |  2 +-
 2 files changed, 63 insertions(+), 6 deletions(-)

diff --git a/docs/developer/subsystems/analytic-solutions.md b/docs/developer/subsystems/analytic-solutions.md
index 2fe4ae6b..9f7f079b 100644
--- a/docs/developer/subsystems/analytic-solutions.md
+++ b/docs/developer/subsystems/analytic-solutions.md
@@ -255,7 +255,7 @@ than by being exempted:
 |---|---|---|
 | `sample_points` | unit box or cube, plus faces and corners | the elliptical inclusion is not box-filling, and the conformal map is singular at its foci — a generic sampler lands on both |
 | `boundaries` | box wall labels | a curved geometry has its own |
-| `apply_boundary_conditions` | free slip or Dirichlet mixin | — |
+| `apply_boundary_conditions` | refuses, with the helpers named | every solution writes it — see below |
 | `stress_is_deviatoric` | `False` | the source publishes $\tau$ |
 
 One consequence worth knowing: `_validation.sample` returns real values, and if
@@ -265,6 +265,62 @@ real-valued without SymPy being able to prove it — but a *genuinely* complex
 result would mean the construction is wrong, and silently taking the real part
 would hide exactly that.
 
+## Boundary conditions: composed functions, not inherited wall types
+
+Every solution writes its own `apply_boundary_conditions`. It composes three
+module-level helpers, each of which takes **the boundaries it applies to**:
+
+```python
+from underworld3.analytic import free_slip, prescribed_velocity, prescribed_scalar
+
+free_slip(solver, boundaries, normal=None)          # strong rotated u.n = 0
+prescribed_velocity(solver, boundaries, velocity)   # Dirichlet velocity
+prescribed_scalar(solver, boundaries, field)        # Dirichlet scalar
+```
+
+A typical Velic solution is then three lines:
+
+```python
+    def apply_boundary_conditions(self, solver):
+        """Free slip on all four walls; the enclosed box has a pressure nullspace."""
+
+        free_slip(solver, self.boundaries)
+        solver.petsc_use_pressure_nullspace = True
+```
+
+**Why not mixins.** This started as `FreeSlipWalls` and `FixedWalls` — mixed in,
+looping over every boundary, applying one condition. That encodes "every boundary
+is a wall of the same kind", which is true of the classical box benchmarks and
+false of most of what the suite has to serve next: a spherical shell or annulus
+with different conditions on the two radii, a faulted disc, a channel driven at
+one end. A mixin cannot express those without growing a parameter for each, and
+inheritance advertises the choice as if it were part of what the solution *is*.
+Functions taking an explicit boundary list say what is imposed *where*, and a
+solution needing two kinds calls two of them:
+
+```python
+    def apply_boundary_conditions(self, solver):
+        free_slip(solver, ["Upper"])
+        prescribed_velocity(solver, ["Lower"], self.fn_velocity)
+```
+
+**The pressure nullspace is stated by the solution, not by a wall type.** It is a
+property of the domain — enclosed, so the pressure is determined only up to a
+constant — and not of any one boundary's condition. Both mixins used to set it,
+which hid that. Leaving it out on an enclosed domain is the failure this whole
+suite exists to catch: a direct solve on the singular saddle returns a quiet,
+wrong answer that only an exact solution exposes.
+
+**Curved boundaries.** `free_slip` takes an optional `normal=`. Leave it out —
+the solver's geometric facet normal is measure-weighted to match the
+straight-facet integral the assembler evaluates, and keeps the constant pressure
+a null vector to machine precision. Pass an analytic normal such as `X/|X|` only
+when the constraint must follow the *true* surface rather than the mesh; it is
+exact for the geometry but keeps a consistency error that grows with facet
+non-uniformity. `CylindricalStokes`, the one curved-geometry solution here, uses
+the default deliberately. See
+[rotated-freeslip.md](rotated-freeslip.md) ("Which normal to use").
+
 ## The stress convention is not uniform across the family
 
 **Check which stress a kernel publishes. Do not read it off the variable name.**
@@ -545,8 +601,9 @@ Starting a comparison from a smooth profile at $t > 0$ rather than from the step
 itself is the point of using these: the step is not representable on the mesh,
 which is what makes `test_1100`'s current comparison fragile.
 
-They cannot reuse the Stokes boundary-condition mixins, which apply a *velocity*.
-`_Transport` prescribes `fn_solution` on every wall instead.
+Their boundary condition prescribes the field itself, not a velocity, so it uses
+`prescribed_scalar`. `_Transport` applies `fn_solution` on every wall, and there
+is no pressure nullspace to remove.
 
 ### Transient solutions are singular at t = 0, and that changes how you use them
 
@@ -802,8 +859,8 @@ a spatial singularity (the inclusion's foci, SolCx's isoviscous limit) needs a
 integer wavenumbers; the Gardner solutions require $\psi < 0$. Ranges become
 constructor validation, which is where they stop being folklore.
 
-1. Subclass `AnalyticSolution` and one of the boundary-condition mixins
-   (`FreeSlipWalls`, `FixedWalls`) — or `_Transport` for a scalar solution.
+1. Subclass `AnalyticSolution` — or `_Transport` for a scalar solution — and
+   write `apply_boundary_conditions`, composing the helpers below.
 2. Build the exact fields on `mesh.X` in `__init__`; set `dim`, `reference`, and
    the `eqn_*` LaTeX strings that document the *problem*. Set the fields through
    `set_fields` (Stokes) or `set_scalar_field` (transport) rather than assigning
diff --git a/tests/test_1024_analytic_conformance.py b/tests/test_1024_analytic_conformance.py
index 209be181..a323ef03 100644
--- a/tests/test_1024_analytic_conformance.py
+++ b/tests/test_1024_analytic_conformance.py
@@ -270,7 +270,7 @@ def test_stress_and_strain_rate_agree(name, built):
             assert np.abs(_validation.sample(sol, difference, points)).max() / scale < 1.0e-8
 
 
-# The boundary-condition mixins are exercised in test_1016_analytic_contract.py.
+# The boundary-condition helpers are exercised in test_1016_analytic_contract.py.
 # Building a Stokes solver per solution here as well was not worth what it cost:
 # it dominated the runtime of this file without checking anything the contract
 # tests do not already cover.

From 51e0b09b309b4c25d11f4c7f0ef37118bbda174c Mon Sep 17 00:00:00 2001
From: lmoresi <louis.moresi@anu.edu.au>
Date: Sun, 16 Aug 2026 09:27:56 +1000
Subject: [PATCH 3/4] Name issue #577 at the CylindricalStokes nullspace TODO

Underworld development team with AI support from Claude Code
---
 src/underworld3/analytic/kramer.py | 14 +++++++-------
 1 file changed, 7 insertions(+), 7 deletions(-)

diff --git a/src/underworld3/analytic/kramer.py b/src/underworld3/analytic/kramer.py
index f386170b..13330513 100644
--- a/src/underworld3/analytic/kramer.py
+++ b/src/underworld3/analytic/kramer.py
@@ -267,13 +267,13 @@ def apply_boundary_conditions(self, solver):
         """
 
         if self.boundary == "zero":
-            # TODO(BUG): the zero-slip case leaves the pressure nullspace in
-            # place, where every other enclosed case in the suite removes it. An
-            # annulus with both arcs held at zero velocity is enclosed, so its
-            # pressure is determined only up to a constant and a direct solve on
-            # the singular saddle can return a quiet, wrong answer. Preserved
-            # here because this refactor is behaviour-preserving by contract;
-            # fix it under its own test.
+            # TODO(BUG): issue #577 — the zero-slip case leaves the pressure
+            # nullspace in place, where every other enclosed case in the suite
+            # removes it. An annulus with both arcs held at zero velocity is
+            # enclosed, so its pressure is determined only up to a constant and
+            # a direct solve on the singular saddle can return a quiet, wrong
+            # answer. Preserved here because this refactor is
+            # behaviour-preserving by contract; fix it under its own test.
             prescribed_velocity(solver, self.boundaries, (0.0, 0.0))
             return
 

From df1b7b27af69dbddc0b23a177c580375a156b85a Mon Sep 17 00:00:00 2001
From: lmoresi <louis.moresi@anu.edu.au>
Date: Sun, 16 Aug 2026 09:32:24 +1000
Subject: [PATCH 4/4] Review fixes: test the free_slip call, not the solver's
 default

The normal test asserted `registered["Left"] is None`, which passes because
add_rotated_freeslip_bc's own default happens to be None. The claim being made
is that free_slip does not pass `normal` when none was chosen, so record the
call instead; a separate test confirms the recorded call is one a real Stokes
solver accepts.

Also widen prescribed_velocity's documented argument type, which excluded the
one in-tree caller that passes a tuple, and say why prescribed_scalar brackets
its argument where prescribed_velocity does not.

Underworld development team with AI support from Claude Code
---
 pixi.lock                            | 39 +++++-----------------------
 src/underworld3/analytic/_base.py    | 12 ++++++---
 tests/test_1016_analytic_contract.py | 38 ++++++++++++++++++++-------
 3 files changed, 44 insertions(+), 45 deletions(-)

diff --git a/pixi.lock b/pixi.lock
index c2601b56..6e023a1a 100644
--- a/pixi.lock
+++ b/pixi.lock
@@ -5,8 +5,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -1332,8 +1330,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -2422,6 +2418,7 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstandard-0.25.0-py312hd41f8a7_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstd-1.5.7-h85ac4a6_6.conda
       - pypi: https://files.pythonhosted.org/packages/63/af/0f5db57b9397a0f3b7fc204cbef143401a7cadaf982330f97f1ce3d39f34/anthropic-0.120.2-py3-none-any.whl
+      - pypi: https://files.pythonhosted.org/packages/db/56/df7622ef4c96b61a5563a6dd999f10c19ef4da113906a1718f4a0b856a2a/assess-1.4-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/a7/5f/ed01f9a3cdffbd5a008556fc7b2a08ddb1cc6ace7effa7340604b1d16699/docstring_parser-0.18.0-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/88/6c/59bd309cab4460c54cf1079f3eb7fe7af6a4c895c5c957a53378693bad2b/jiter-0.16.0-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl
       - pypi: https://files.pythonhosted.org/packages/ef/81/6809b33b9ff88d75164d9740efe362b80201709557995bbf55cf82e597e8/pygmsh-7.1.12-py3-none-any.whl
@@ -2928,8 +2925,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -4018,6 +4013,7 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstandard-0.25.0-py312hd41f8a7_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstd-1.5.7-h85ac4a6_6.conda
       - pypi: https://files.pythonhosted.org/packages/63/af/0f5db57b9397a0f3b7fc204cbef143401a7cadaf982330f97f1ce3d39f34/anthropic-0.120.2-py3-none-any.whl
+      - pypi: https://files.pythonhosted.org/packages/db/56/df7622ef4c96b61a5563a6dd999f10c19ef4da113906a1718f4a0b856a2a/assess-1.4-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/a7/5f/ed01f9a3cdffbd5a008556fc7b2a08ddb1cc6ace7effa7340604b1d16699/docstring_parser-0.18.0-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/88/6c/59bd309cab4460c54cf1079f3eb7fe7af6a4c895c5c957a53378693bad2b/jiter-0.16.0-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl
       - pypi: https://files.pythonhosted.org/packages/ef/81/6809b33b9ff88d75164d9740efe362b80201709557995bbf55cf82e597e8/pygmsh-7.1.12-py3-none-any.whl
@@ -4524,8 +4520,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -5851,8 +5845,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -6942,6 +6934,7 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstandard-0.25.0-py312hd41f8a7_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstd-1.5.7-h85ac4a6_6.conda
       - pypi: https://files.pythonhosted.org/packages/63/af/0f5db57b9397a0f3b7fc204cbef143401a7cadaf982330f97f1ce3d39f34/anthropic-0.120.2-py3-none-any.whl
+      - pypi: https://files.pythonhosted.org/packages/db/56/df7622ef4c96b61a5563a6dd999f10c19ef4da113906a1718f4a0b856a2a/assess-1.4-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/a7/5f/ed01f9a3cdffbd5a008556fc7b2a08ddb1cc6ace7effa7340604b1d16699/docstring_parser-0.18.0-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/88/6c/59bd309cab4460c54cf1079f3eb7fe7af6a4c895c5c957a53378693bad2b/jiter-0.16.0-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl
       - pypi: https://files.pythonhosted.org/packages/ef/81/6809b33b9ff88d75164d9740efe362b80201709557995bbf55cf82e597e8/pygmsh-7.1.12-py3-none-any.whl
@@ -7449,8 +7442,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -8784,8 +8775,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -9881,6 +9870,7 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstandard-0.25.0-py312hd41f8a7_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstd-1.5.7-h85ac4a6_6.conda
       - pypi: https://files.pythonhosted.org/packages/63/af/0f5db57b9397a0f3b7fc204cbef143401a7cadaf982330f97f1ce3d39f34/anthropic-0.120.2-py3-none-any.whl
+      - pypi: https://files.pythonhosted.org/packages/db/56/df7622ef4c96b61a5563a6dd999f10c19ef4da113906a1718f4a0b856a2a/assess-1.4-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/a7/5f/ed01f9a3cdffbd5a008556fc7b2a08ddb1cc6ace7effa7340604b1d16699/docstring_parser-0.18.0-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/88/6c/59bd309cab4460c54cf1079f3eb7fe7af6a4c895c5c957a53378693bad2b/jiter-0.16.0-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl
       - pypi: https://files.pythonhosted.org/packages/ef/81/6809b33b9ff88d75164d9740efe362b80201709557995bbf55cf82e597e8/pygmsh-7.1.12-py3-none-any.whl
@@ -10390,8 +10380,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -11793,8 +11781,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -13176,8 +13162,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -14306,6 +14290,7 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstandard-0.25.0-py312hd41f8a7_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstd-1.5.7-h85ac4a6_6.conda
       - pypi: https://files.pythonhosted.org/packages/63/af/0f5db57b9397a0f3b7fc204cbef143401a7cadaf982330f97f1ce3d39f34/anthropic-0.120.2-py3-none-any.whl
+      - pypi: https://files.pythonhosted.org/packages/db/56/df7622ef4c96b61a5563a6dd999f10c19ef4da113906a1718f4a0b856a2a/assess-1.4-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/a7/5f/ed01f9a3cdffbd5a008556fc7b2a08ddb1cc6ace7effa7340604b1d16699/docstring_parser-0.18.0-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/88/6c/59bd309cab4460c54cf1079f3eb7fe7af6a4c895c5c957a53378693bad2b/jiter-0.16.0-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl
       - pypi: https://files.pythonhosted.org/packages/ef/81/6809b33b9ff88d75164d9740efe362b80201709557995bbf55cf82e597e8/pygmsh-7.1.12-py3-none-any.whl
@@ -14831,8 +14816,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-20_gnu.conda
@@ -15197,8 +15180,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -16579,8 +16560,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -17708,6 +17687,7 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstandard-0.25.0-py312hd41f8a7_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstd-1.5.7-h85ac4a6_6.conda
       - pypi: https://files.pythonhosted.org/packages/63/af/0f5db57b9397a0f3b7fc204cbef143401a7cadaf982330f97f1ce3d39f34/anthropic-0.120.2-py3-none-any.whl
+      - pypi: https://files.pythonhosted.org/packages/db/56/df7622ef4c96b61a5563a6dd999f10c19ef4da113906a1718f4a0b856a2a/assess-1.4-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/a7/5f/ed01f9a3cdffbd5a008556fc7b2a08ddb1cc6ace7effa7340604b1d16699/docstring_parser-0.18.0-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/88/6c/59bd309cab4460c54cf1079f3eb7fe7af6a4c895c5c957a53378693bad2b/jiter-0.16.0-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl
       - pypi: https://files.pythonhosted.org/packages/ef/81/6809b33b9ff88d75164d9740efe362b80201709557995bbf55cf82e597e8/pygmsh-7.1.12-py3-none-any.whl
@@ -18235,8 +18215,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -19624,8 +19602,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
@@ -20759,6 +20735,7 @@ environments:
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstandard-0.25.0-py312hd41f8a7_1.conda
       - conda: https://conda.anaconda.org/conda-forge/linux-aarch64/zstd-1.5.7-h85ac4a6_6.conda
       - pypi: https://files.pythonhosted.org/packages/63/af/0f5db57b9397a0f3b7fc204cbef143401a7cadaf982330f97f1ce3d39f34/anthropic-0.120.2-py3-none-any.whl
+      - pypi: https://files.pythonhosted.org/packages/db/56/df7622ef4c96b61a5563a6dd999f10c19ef4da113906a1718f4a0b856a2a/assess-1.4-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/a7/5f/ed01f9a3cdffbd5a008556fc7b2a08ddb1cc6ace7effa7340604b1d16699/docstring_parser-0.18.0-py3-none-any.whl
       - pypi: https://files.pythonhosted.org/packages/88/6c/59bd309cab4460c54cf1079f3eb7fe7af6a4c895c5c957a53378693bad2b/jiter-0.16.0-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl
       - pypi: https://files.pythonhosted.org/packages/ef/81/6809b33b9ff88d75164d9740efe362b80201709557995bbf55cf82e597e8/pygmsh-7.1.12-py3-none-any.whl
@@ -21287,8 +21264,6 @@ environments:
     - url: https://conda.anaconda.org/conda-forge/
     indexes:
     - https://pypi.org/simple
-    options:
-      pypi-prerelease-mode: if-necessary-or-explicit
     packages:
       linux-64:
       - conda: https://conda.anaconda.org/conda-forge/linux-64/_openmp_mutex-4.5-7_kmp_llvm.conda
diff --git a/src/underworld3/analytic/_base.py b/src/underworld3/analytic/_base.py
index b029e9eb..fbec7ec3 100644
--- a/src/underworld3/analytic/_base.py
+++ b/src/underworld3/analytic/_base.py
@@ -620,8 +620,10 @@ def prescribed_velocity(solver, boundaries, velocity):
         The solver to configure.
     boundaries : sequence of str
         The boundary labels to constrain.
-    velocity : sympy Matrix, shape (1, dim)
-        Usually the solution's own ``fn_velocity``.
+    velocity : sequence of length dim
+        Usually the solution's own ``fn_velocity``, a 1 x dim sympy Matrix, but
+        any dim-long sequence will do — a solution holding a wall at rest passes
+        a tuple of zeros.
     """
 
     for boundary in boundaries:
@@ -641,7 +643,11 @@ def prescribed_scalar(solver, boundaries, field):
     boundaries : sequence of str
         The boundary labels to constrain.
     field : sympy expression
-        Usually the solution's own ``fn_solution``.
+        Usually the solution's own ``fn_solution``. A single expression, not a
+        sequence — the bracket below is this function's job, because
+        ``add_dirichlet_bc`` takes one entry per component and a scalar solution
+        has exactly one. That is the only difference from
+        :func:`prescribed_velocity`, whose argument already is the sequence.
     """
 
     for boundary in boundaries:
diff --git a/tests/test_1016_analytic_contract.py b/tests/test_1016_analytic_contract.py
index 1ca28dbd..7ebdac4d 100644
--- a/tests/test_1016_analytic_contract.py
+++ b/tests/test_1016_analytic_contract.py
@@ -331,24 +331,42 @@ def test_boundaries_can_carry_different_conditions(mesh):
 
 
 def test_free_slip_passes_an_analytic_normal_through(mesh):
-    """A curved-boundary solution can choose the normal; the default is geometric.
+    """A curved-boundary solution can choose the normal; otherwise it says nothing.
 
-    Both calls must register, and the analytic one must carry the normal it was
-    given — otherwise a spherical or annulus solution silently gets the facet
-    normal it explicitly declined.
+    The distinction is between passing `normal=None` and not passing `normal` at
+    all: the first overrides the solver's default with nothing, the second leaves
+    the solver to choose. So this records the *call*, rather than reading the
+    value back off the solver — reading it back would only confirm what the
+    solver's own default happens to be.
     """
 
-    stokes = uw.systems.Stokes(mesh)
     x, y = mesh.X
     radial = sympy.Matrix([[x, y]]) / sympy.sqrt(x**2 + y**2)
 
-    uw.analytic.free_slip(stokes, ["Left"])
-    uw.analytic.free_slip(stokes, ["Right"], normal=radial)
+    class Recorder:
+        def __init__(self):
+            self.calls = []
+
+        def add_rotated_freeslip_bc(self, value, boundary, **kwargs):
+            self.calls.append((value, boundary, kwargs))
+
+    solver = Recorder()
+    uw.analytic.free_slip(solver, ["Left"])
+    uw.analytic.free_slip(solver, ["Right"], normal=radial)
+
+    assert solver.calls == [
+        (0.0, "Left", {}),
+        (0.0, "Right", {"normal": radial}),
+    ]
 
-    registered = dict(stokes._rotated_freeslip_bcs)
 
-    assert registered["Left"] is None
-    assert registered["Right"] == radial
+def test_free_slip_reaches_a_real_solver(mesh):
+    """The recorded call is the one a Stokes solver actually accepts."""
+
+    stokes = uw.systems.Stokes(mesh)
+    uw.analytic.free_slip(stokes, ["Left"])
+
+    assert {boundary for boundary, _ in stokes._rotated_freeslip_bcs} == {"Left"}
 
 
 def test_boundaries_follow_the_mesh_dimension(mesh):