Skip to content

Nitsche penalty scale regressed with mesh.cell_size(): the Zhong spherical-shell response went from 0.2% to 2.4-5.7% off #734

Description

@lmoresi

add_nitsche_bc scales its penalty as gamma*mu/h and takes h from
mesh.cell_size() by default (local_h=True). #692 and #687 changed what
cell_size() returns. The Zhong et al. (2008) spherical-shell response
benchmark in tests/test_1064_constrained_spherical_shell_response.py has
degraded by an order of magnitude as a result, and nobody saw it because that
test matched no glob in scripts/test.sh and had never run in CI (#721).

Measured

Same case, same mesh, same code — only the Nitsche h differs:

Nitsche h surface error vs Zhong CMB error
local_h=True — current default, mesh.cell_size() 2.36% 0.94%
local_h=False — legacy global get_min_radius() 0.20% 0.22%

(SphericalShell r 0.55→1.0, internal boundary 0.775, cellSize 1/8, degree-2
harmonic load, NITSCHE_GAMMA as the test sets it. Both converged,
reason=5.)

The legacy global h reproduces the published response essentially exactly. The
current default does not.

Why the scale moved

On this mesh the two measures differ by a factor of 1.745:

mesh.cell_size()   mean 0.09282   min 0.06625   max 0.11970
get_min_radius()        0.05319

so the penalty gamma*mu/h is 0.573x — 43% weaker — than it was. A weaker
Nitsche penalty enforces v.n = 0 less tightly, which is exactly the quantity
this benchmark measures.

Note the docstring on add_nitsche_bc:

Set local_h=False to restore the legacy single global-minimum scalar
(mesh.get_min_radius()); on a uniform mesh the two coincide.

They do not coincide here, on a near-uniform spherical shell — 1.745x comparing
mean to global, and still 1.25x comparing minimum to minimum. The two are
measuring different things. This is the "three live definitions of cell radius"
problem (#694) landing somewhere that has a published answer to check against.

Timeline

date event
2026-06-24 test calibrated against Zhong (#265); local_h=True becomes the Nitsche default the same day (#275)
2026-09-07 #692 changes what cell_size() returns
2026-09-08 #687 touches it again
2026-09-12 un-deferring test_106* runs the test for the first time; it fails

local_h is not the culprit — it was already the default when the reference was
calibrated. What moved is the value cell_size() returns, and the test was dark
throughout.

It is now platform-decided

CI measured 5.67% and failed the 5% tolerance; the same case on a developer
Mac gives 2.36% and passes. gmsh triangulates differently across platforms
(this repo already records it — the serial_reference skips note cells=434 vs
430). At 0.2% that difference was irrelevant; at 2-6% the platform decides the
outcome.

What this is not

It is not a tolerance to relax. add_nitsche_bc is used across the free-slip
and fault work, so the question is which penalty scale is correct in general:

  1. the old scale was right, and the new cell_size() makes Nitsche too weak
    everywhere; or
  2. the new scale is right, and the Zhong tolerance was calibrated against an
    over-stiff penalty — in which case the 0.2% agreement it used to show was
    partly luck.

The first is the more likely reading given a 0.2% match to a published result,
but it needs the resolution ladder to settle rather than one data point.

Suggested first steps

  • Run the response at cellSize 1/8, 1/12, 1/16 under both h choices: if the
    local-h error converges away it is a scale issue, if it plateaus it is a
    consistency issue.
  • Decide what h the Nitsche facet term should actually use — the adjacent
    cell's insphere radius is the usual choice for an interior-penalty method, and
    it is not obvious that either current measure is that.
  • Correct or remove the "on a uniform mesh the two coincide" claim in the
    add_nitsche_bc docstring; it is false as measured.

Found during the 2026-09 audit, by turning on tests that had never run.
Related: #721 (the dark tests), #694 (cell radius definitions).

Underworld development team with AI support from Claude Code

Activity

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

      Milestone

      No milestone

      Relationships

      None yet

      Development

      No branches or pull requests

      Issue actions