Stokes_Constrained: selfp Schur PC default, viscosity-scaled penalty, drop redundant saddle_preconditioner, and nullspace re-setup fix by lmoresi · Pull Request #229 · underworldcode/underworld3

lmoresi · 2026-06-10T02:28:18Z

Architectural Review Submission

Review Type: Code Implementation
Priority: MEDIUM
Review Period: 2026-06

📄 Review Document

File: docs/developer/design/CONSTRAINED_FREESLIP_MULTIPLIER.md
Tracking Index: N/A

📋 Executive Summary

This PR improves robustness of uw.systems.Stokes_Constrained for strong viscosity contrast by defaulting the constrained Schur preconditioner to PETSc selfp, removes a redundant constrained-solver knob, and corrects viscosity scaling for optional incompressibility penalty terms in base Stokes (penalty != 0).
A follow-up fix also ensures petsc_use_pressure_nullspace toggles force solver re-setup (is_setup = False) and replaces a local-path code comment reference with a stable in-repo design-doc link.

🎯 Scope

Systems Affected:

Core solvers (Stokes, Poisson, Advection-Diffusion)
Variable system (Mesh/Swarm variables, arrays)
Units & non-dimensionalization
Parallel operations (MPI, collective operations)
Testing infrastructure
Documentation & examples
Other: constrained free-slip multiplier/topography workflow

Files Changed:

New files: 0
Modified files: 5
Total LOC: ~moderate, solver+tests+design-doc focused

✅ Testing Status

All tests passing: pixi run underworld-test
New tests added: constrained coverage extensions in test_1061
Test coverage: N/A
Regression tests validated
Performance benchmarks run (if applicable)

Test Results:

# Command used
pixi run -e default pytest tests/test_1061_constrained_freeslip.py tests/test_1062_constrained_solcx.py tests/test_1010_stokesCart.py tests/test_1015_analytic_solcx.py -v

# Results summary
Targeted constrained/base-Stokes paths green (including test_1061 parameterized set = 15 cases)

📊 Key Metrics & Impact

Performance:

Before: constrained Schur with a11 degraded with high viscosity contrast (iteration growth/divergence cases)
After: constrained Schur defaults to selfp, with stable low iteration counts across large contrast in tested annulus setup
Improvement: materially improved robustness and convergence behavior under strong contrast

Quality:

Test coverage: improved for constrained free-slip contrast/default-options path
Code complexity: unchanged to slightly increased (small targeted solver-path updates)
Documentation: design note expanded with Schur/augmentation/topography guidance

🔍 Review Checklist (for Reviewers)

Design & Architecture:

Design rationale is clear and well-justified
Trade-offs are documented with alternatives considered
System architecture is comprehensible
Integration points are identified

Implementation:

Implementation matches documented design
Code quality meets project standards
Breaking changes are identified and justified
Backward compatibility is addressed

Testing & Validation:

Testing strategy is adequate
Test coverage is sufficient
Edge cases are covered
Performance impact is assessed

Documentation:

Known limitations are clearly documented
Benefits are quantified
User-facing changes are documented
Migration guide provided (if needed)

⚠️ Known Limitations

Constrained solver path remains serial-only (consistent with existing implementation).
Base Stokes penalty scaling change only affects non-default penalty != 0 workflows; users relying on prior unscaled behavior should review.

🔗 Dependencies & Related Work

Depends on:

Review #XXX (must be approved first)

Blocks:

Review #YYY (waiting on this review)

Related:

Discussion: N/A
Prior art: docs/developer/design/CONSTRAINED_FREESLIP_MULTIPLIER.md

🖊️ Sign-Off

This PR merge represents formal approval of the architectural review.

Reviewers: Please review the full document at docs/developer/design/CONSTRAINED_FREESLIP_MULTIPLIER.md and approve this PR only when satisfied with the design, implementation, and documentation.

Role	GitHub Handle	Sign-Off Date	Status
Author	@lmoresi	2026-06-10	✅ Submitted
Primary Reviewer	@reviewer1		⏳ Pending
Secondary Reviewer	@reviewer2		⏳ Pending
Project Lead	@lead		⏳ Pending

📝 Additional Context

Follow-up feedback fixes included in latest branch state:

petsc_use_pressure_nullspace setter now marks is_setup = False so post-setup toggles correctly re-register nullspace during _setup_solver.
Local-path comment reference removed and replaced with stable in-repo doc link.
These do not change the core PR goals; they tighten correctness/maintainability around the introduced solver behavior.

Underworld development team with AI support from Claude Code

- Default pc_fieldsplit_schur_precondition to `selfp` for Stokes_Constrained (base Stokes stays `a11`). selfp builds S ~ A11 - A10 diag(A00)^-1 A01 from the actual operator blocks; its lambda-lambda corner is the true constraint Schur -C diag(A)^-1 C^T. It cracks strong boundary-viscosity-contrast walls that the bare `a11` mass cannot, makes the augmentation optional on contrast problems that have a velocity anchor, and keeps the recovered multiplier (dynamic topography) clean at small augmentation. Override via petsc_options if needed. - Add experimental `multiplier_schur_pc` flag (default off): injects a 1/mu mass into the lambda block's preconditioner only (Pmat, via PetscDSSetJacobianPreconditioner, reusing pressure's _pp_G0) while keeping the operator (lambda,lambda) block exact. Superseded by selfp for production but kept as an escape hatch. - Extend docs/developer/design/CONSTRAINED_FREESLIP_MULTIPLIER.md: selfp; the augmentation's roles (Schur conditioning AND rigid-body null-space regularisation of the inner velocity block); the mixed / penalty / augmented topography triad (lambda is to u.n=g as p is to div u=0) and how large augmentation contaminates the recovered topography; and a warning that the inner velocity-block rigid-body singularity must NOT be removed via the coupled petsc_velocity_nullspace_basis on a constrained problem (it corrupts the solution -- rotation carries real circulation, translations hit the Euclidean vs L2 inner-product mismatch); the augmentation is the working remedy. - test_1061: contrast-annulus selfp section + topography recovery. Underworld development team with AI support from Claude Code

The incompressibility penalty term in F1 was `penalty * div_u * I` with `penalty` a bare constant. Under spatially-variable viscosity a constant over-stiffens low-viscosity regions (penalty/mu huge there) and is negligible in high-viscosity regions, locking the velocity on the soft side -- e.g. a constant penalty of 100 gives SolCx (1e6 contrast) velocity error 0.1. Scale the penalty by the local viscosity (constitutive_model.K) so the term is `penalty * mu * div_u * I` and the ratio penalty/mu is uniform; `penalty` is now a dimensionless O(1) base. Verified: penalty=1 (mu-scaled) gives SolCx velocity error 7.9e-6 (vs 0.1 for the old constant 100). This matches how the boundary augmented-Lagrangian is already scaled (augmentation_base * mu). Inherited by Stokes_Constrained. Default penalty=0 is unchanged (term vanishes) -- no regression (test_1061 + test_1010 green). Document the consequent pressure correction: with the penalty on, the recovered `p` is the multiplier, and the mechanical pressure is `p_mech = p - penalty*mu*div_u` (they agree at convergence; differ pointwise by ~a few percent at penalty=O(1)). Use p_mech for pressure-dependent constitutive laws. The penalty is usually unnecessary -- the 1/mu Schur preconditioner already conditions the pressure block. Underworld development team with AI support from Claude Code

…erty The pressure/constraint Schur preconditioner is built automatically: `selfp` forms it from the operator blocks, and the pressure mass it needs is the 1/viscosity (1/constitutive_model.K) term supplied by the existing fallback when saddle_preconditioner is unset. Setting `saddle_preconditioner = 1/viscosity` is therefore exactly the automatic value (verified bit-identical, rel-diff 0.0) -- it adds nothing but a chance to supply an inconsistent viscosity. - Stokes_Constrained: override the property so it is inert (getter returns None) and assigning to it raises with a clear message. The internal `_pp_G0` uses the 1/K fallback as before. - Strip the redundant `saddle_preconditioner = 1/viscosity` lines from the constrained-solver tests (test_1061, test_1062); they stay green on the fallback. - Base SNES_Stokes keeps `saddle_preconditioner` as an advanced override; its docstring now states it is usually unset (the automatic 1/K is the default) and is only for non-standard K / deliberate tuning. Underworld development team with AI support from Claude Code

Copilot

Pull request overview

This PR improves robustness and usability of the constrained free-slip Stokes solver (uw.systems.Stokes_Constrained) under strong boundary viscosity contrast, adjusts augmented-Lagrangian incompressibility penalty scaling to be viscosity-weighted, and removes a redundant constrained-solver knob (saddle_preconditioner). It also updates tests and design documentation accordingly.

Changes:

Default Stokes_Constrained grouped Schur preconditioner to PETSc selfp, and make constrained saddle_preconditioner inert (raises on set).
Viscosity-scale the optional grad-div (augmented Lagrangian) incompressibility penalty term (penalty * mu * div(u) * I) and document the resulting pressure correction.
Extend constrained free-slip tests (including a strong-contrast annulus case) and expand the developer design note; add an opt-in multiplier-block Schur Pmat option (multiplier_schur_pc).

Reviewed changes

Copilot reviewed 5 out of 5 changed files in this pull request and generated 2 comments.

Show a summary per file

File	Description
`tests/test_1062_constrained_solcx.py`	Removes redundant constrained-solver `saddle_preconditioner` assignment in SolCx regression test.
`tests/test_1061_constrained_freeslip.py`	Removes redundant `saddle_preconditioner` assignments; adds strong-contrast annulus coverage and asserts default `selfp`.
`src/underworld3/systems/solvers.py`	Viscosity-scales incompressibility penalty in `F1`; documents penalty/pressure implications; sets `selfp` default for constrained solver; disables constrained `saddle_preconditioner`.
`src/underworld3/cython/petsc_generic_snes_solvers.pyx`	Adds opt-in `multiplier_schur_pc` and wires Pmat-only multiplier preconditioner support.
`docs/developer/design/CONSTRAINED_FREESLIP_MULTIPLIER.md`	Expands design/conditioning notes (selfp Schur, augmentation roles, nullspace cautions) and documents penalty scaling rationale.

💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.

    @petsc_use_pressure_nullspace.setter
    def petsc_use_pressure_nullspace(self, value):
        self._petsc_use_pressure_nullspace = bool(value)
        self._reset_stokes_nullspace()
-        self.is_setup = False
+


+        # makes the augmentation optional. Default OFF (opt-in): a monolithic lu solve
+        # factorizes the Pmat (pc_use_amat is a no-op here), so this Pmat term is NOT
+        # inert for direct solves — flipping the default on regresses the direct-lu
+        # constraint enforcement. See ~/.claude/plans for the load-bearing-Pmat note.


…ce toggle - Replace a `~/.claude/plans` reference in a code comment with the in-repo design note path (CONSTRAINED_FREESLIP_MULTIPLIER.md). - `petsc_use_pressure_nullspace` setter now also sets `is_setup = False` (matching `petsc_use_nullspace`): `_attach_stokes_nullspace` runs inside the is_setup-gated `_setup_solver`, so toggling the flag after setup previously would not re-register the changed nullspace with PETSc. No effect on the usual set-before-solve path. Underworld development team with AI support from Claude Code

lmoresi · 2026-06-10T04:47:32Z

Thanks @copilot — both points addressed in 19fd5b8:

Removed the ~/.claude/plans reference from the code comment; it now points to docs/developer/design/CONSTRAINED_FREESLIP_MULTIPLIER.md.
petsc_use_pressure_nullspace setter now also sets is_setup = False (matching petsc_use_nullspace), so toggling it after setup re-registers the nullspace — _attach_stokes_nullspace runs inside the is_setup-gated _setup_solver. No effect on the usual set-before-solve path; constrained suite stays green (test_1061 15, test_1062).

lmoresi · 2026-06-10T05:44:13Z

This seems good to me because I wrote it. But this is also largely orthogonal to existing solvers.

lmoresi added 3 commits June 10, 2026 07:32

Copilot AI review requested due to automatic review settings June 10, 2026 02:28

Copilot started reviewing on behalf of lmoresi June 10, 2026 02:28 View session

Copilot AI reviewed Jun 10, 2026

View reviewed changes

Copilot started work on behalf of lmoresi June 10, 2026 04:47 View session

Copilot AI changed the title ~~Stokes_Constrained: selfp Schur PC default, viscosity-scaled penalty, drop redundant saddle_preconditioner~~ Stokes_Constrained: selfp Schur PC default, viscosity-scaled penalty, drop redundant saddle_preconditioner, and nullspace re-setup fix Jun 10, 2026

Copilot finished work on behalf of lmoresi June 10, 2026 04:49

lmoresi merged commit 87d778d into development Jun 10, 2026
2 checks passed

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Stokes_Constrained: selfp Schur PC default, viscosity-scaled penalty, drop redundant saddle_preconditioner, and nullspace re-setup fix#229

Stokes_Constrained: selfp Schur PC default, viscosity-scaled penalty, drop redundant saddle_preconditioner, and nullspace re-setup fix#229
lmoresi merged 4 commits into
developmentfrom
feature/block-multiplier-preconditioner

lmoresi commented Jun 10, 2026 •

edited by Copilot AI

Loading

Uh oh!

Copilot AI left a comment

Uh oh!

lmoresi commented Jun 10, 2026

Uh oh!

lmoresi commented Jun 10, 2026

Uh oh!

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

2 participants

Conversation

lmoresi commented Jun 10, 2026 • edited by Copilot AI Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Architectural Review Submission

📄 Review Document

📋 Executive Summary

🎯 Scope

✅ Testing Status

📊 Key Metrics & Impact

🔍 Review Checklist (for Reviewers)

⚠️ Known Limitations

🔗 Dependencies & Related Work

🖊️ Sign-Off

📝 Additional Context

Uh oh!

Copilot AI left a comment

Choose a reason for hiding this comment

Pull request overview

Reviewed changes

Uh oh!

lmoresi commented Jun 10, 2026

Uh oh!

lmoresi commented Jun 10, 2026

Uh oh!

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

2 participants

lmoresi commented Jun 10, 2026 •

edited by Copilot AI

Loading