Skip to content

C19 round 3: RDLS on the engine; RANS2P defaults to Schur + BoomerAMG (ASM without hypre) - #212

Draft
cekees wants to merge 7 commits into
devfrom
c19-r3
Draft

C19 round 3: RDLS on the engine; RANS2P defaults to Schur + BoomerAMG (ASM without hypre)#212
cekees wants to merge 7 commits into
devfrom
c19-r3

Conversation

@cekees

@cekees cekees commented Nov 5, 2026 •
edited
Loading

Copy link
Copy Markdown
Owner

C19 round 3 (stream D3). Each commit message has its numbers.

commit what numbers
1994e8f RDLS on the engine's scalar SUPG kernel: the third two-phase slice and the first nonlinear Hamiltonian, checked against proteus's RDLS on damBreak residual 4.1e-14 of its max (4.5e-17 absolute), Jacobian 3.9e-17, next Newton iterate 1.9e-16; the VOF/NCLS fixtures are unchanged
6095102 3D level-set dam break deck. The perf harness now compares solutions across partitions (node coordinates) and can set the velocity split's block size partition noise floor (ASM, he=0.15, 4 ranks): velocity 1.3e-5, pressure 3.6e-7, water volume 2.1e-9
6c0b4f1, be4b929 RANS2P default flip to 'auto' (fieldsplit with hypre, ASM without), then reverted after the 96-rank check: the default fieldsplit has 3.2-3.3e-4 differences at wall nodes (z = 0, y = 0.53-0.57) that recur on a different partition; the fieldsplit stays available as preconditioner = 'schur_boomeramg' 24 ranks 7.5x faster on the solve (55 s vs 416 s), 96 ranks 1.8x (44 s vs 80 s), KSP its flat (2.8); see below
752005b Coverage floor NCLS3P.py 76.2 -> 75.6, the same hunk as #209/#211 main's 76f3e12 (scatterDOFArray de-duplication) took NCLS3P to 75.7%; this PR touches no NCLS3P path
32b650e TwoPhaseFlow's iterative options PETSc never read: VOF/NCLS/MCorr pc_pc_hypre_type, RDLS pc_pc_asm_type and max_it. Only RDLS's ASM type changes behaviour (restrict -> the basic it asked for) damBreakLevelSet 2D, 2 ranks, -options_left: 5 unused options -> none; RDLS KSP its 2.4 -> 2.1, Newton unchanged; solution moves 1.9e-6 (RDLS, NCLS), 1.2e-6 (velocity), water volume 2.6e-14; a rerun is bitwise
bfa1a68 MCorr on the engine's generic element kernel (not the SUPG template), with a bindPoint(q) hook for point fields; checked against proteus's MCorr on damBreak residual 6.9e-15, Jacobian 3.3e-16, q_r 5.4e-15, next Newton iterate 2.5e-14, residual there 3.2e-15 of the first, q_u/q_r/q_n there 1.1e-14; the Richards tests print the same digits as before the hook

The 3D check (damBreakLevelSet3D, 110 k tets, 4 ranks, M4 Max at load 55-100; times are the min of 2 interleaved runs):

RANS2P Newton its/step KSP its/Newton RANS2P/Newton solve velocity vs tight ref
'asm' 3.0 6.0-7.1 3333 ms 60.0 s 1.5-1.7e-5
'schur_boomeramg' 3.0 2.8-3.0 516 ms 14.4 s 9.9e-6
+ velocity bs 3 3.0 3.0-3.1 614 ms 19.5 s 3.9e-6
+ bs 3 + nodal coarsening 3.0 3.0 1.1e-4 (rejected)

The floor is 1.5e-5 (the two tight references differ by that much), so everything except nodal coarsening is within it. The block size buys no iterations, so proteus keeps bs = 1.

ceds-model, 24 and 96 ranks (247 k nodes, 2 steps, quiet machine, baseline §9.2):

ranks RANS2P Newton its/step KSP its/Newton RANS2P/Newton solve
24 'asm' 4.0 13.4 23742 ms 415.8 s
24 'schur_boomeramg' 3.0 2.8 2530 ms 55.2 s (7.5x)
96 'asm' 4.0 10.8 3619 ms 79.9 s
96 'schur_boomeramg' 3.0 2.8 1953 ms 44.4 s (1.8x)

The fieldsplit's KSP count does not grow with ranks. Solution: five of six runs agree to rms 1.4e-6; ASM at 96 ranks is the outlier (rms 1.2e-5 from all of them, uniformly). The fieldsplit at 96 ranks matches them except at 7 floor nodes (z = 0, y = 0.53-0.57, up to 3.3e-4 of max |vel|, against 6e-5 between the other runs); the other 96-rank fieldsplit run has none. The repeat of that run (different partition: 7 of 247 482 nodes with the same owner) has the same worst node, 3.2e-4 at (1.22, 0.547, 0), and the same second node; the difference fields correlate at 0.56. So not partition noise, and the default flip is reverted (be4b929). Whether it is the fieldsplit or ASM's looser convergence at 96 ranks needs a tight ASM reference there; the flip can return after.

The hypre survey: ON_NODE_PERFORMANCE_BASELINE.md §9.3 (local).

  • Have hypre: conda-forge real builds (all four platforms), proteus's feedstock and CI envs, the PETSc fork's configure scripts, the colonized prefixes, the HPCMP and LSU module builds, Spack py-proteus.
  • No hypre (the opt-in 'schur_boomeramg' would be refused there; the default is ASM since the revert): Homebrew's petsc, Spack's petsc without +hypre, and a PyPI petsc built without options.

Tests: test/TwoPhaseFlow/test_rans2p_preconditioner.py (5 passed on the Mac); experimental/adrh vof/ncls/rdls/mcorr tests and richards_layers/richards_nd_test pass with adrh_proto on the Mac (the Kokkos tests need adrh_kokkos, CI only).

Found and not acted on (C0's territory): RDLS.h writes q[('grad(u)',0)] from an uninitialized local (dir).

🤖 Generated with Claude Code

cekees and others added 3 commits November 5, 2026 22:00
@cekees @claude
C19: RDLS on the engine's scalar SUPG kernel, against proteus's RDLS …
1994e8f 
…on damBreak

The third model of the two-phase level-set family, done the way VOF (#182) and
NCLS (#198) were. RDLS (mprans/RDLS.h, ELLIPTIC_REDISTANCING=0, freezeLevelSet=1,
useExact=0, no time integration, shock capturing and subgrid error not lagged: the
damBreakLevelSet configuration) is the first nonlinear Hamiltonian on the template:

  S = H_eps(phi) - (1 - H_eps(phi)), eps = epsFact h_phi,
  H = S |grad u|, dH = S grad u/(|grad u| + 1e-12), r = -S.

- vof.h: three additions, each behind a functor trait or a parameter, so VOF and
  NCLS compile as before (vof_test.py and ncls_test.py give the same numbers):
  a reaction r (hasR) in the weak form and the strong residual; RDLS's tau
  1/(sqrt(dH.G dH) + 1e-8) (tauForm = 1; VOF's form would differ at 1e-8 with
  dm_t = 0); shock capturing at the current u (SUPGParameters::lag_sc = 0) plus a
  background diffusion the functor returns (nuExtra).
- rdls.h: RDLSCoefficients (fields: NCLS's level set, h_phi, h_e), the interface
  pin rdlsPinElement (rows of nodes with |u_last| < eps replaced by
  (u - u_last) wDF h_e; eps is the element's last quadrature point's, as RDLS.h
  reuses the loop variable), and rdlsTraceFace (ebqe_u, ebqe_n; RDLS has no
  boundary terms). The Jacobian is proteus's: exact in the H term, dH and tau
  frozen in the subgrid term, nu not differentiated.
- vof_capture.py model=rdls; the CSR structure comes from the matrix (RDLS's
  kernel only gets element-local offsets). rdls_test.py + rdls_damBreak_he0.05.npz
  (733 triangles, t = 0.062, first Newton iteration of a redistancing solve,
  30 pinned nodes / 165 pinned element rows); ctest proto_rdls_test.

Engine vs proteus (max relative difference):
  residual 4.1e-14 of its max (4.5e-17 absolute, 3e-15 of a node's H and r
  terms, which cancel where |grad u| = 1), Jacobian 3.9e-17, shock-capturing
  viscosity 4.7e-15, CFL 3.2e-16, boundary trace bitwise, normalized boundary
  gradient 2.4e-15, next Newton iterate 1.9e-16, residual at the next iterate
  (pinned rows nonzero there) 6.9e-14 of its max.

Finding, not acted on (C0's territory): RDLS.h writes q[('grad(u)',0)] from an
uninitialized local (`dir`; its assignment is commented out) in
calculateResidual and calculateJacobian.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
@cekees @claude
TwoPhaseFlow: a 3D level-set dam break; perf harness matches solution…
6095102 
…s across partitions

The 3D RANS2P case the round-2 handoff asked for (no TwoPhaseFlow deck ran RANS2P in
3D) and the two harness changes its check needs.

- test/TwoPhaseFlow/damBreakLevelSet3D.py: RANS2P, VOF, NCLS, RDLS and MCorr
  (useDefaultModels(0, 0)) in a 3.22 x 1.0 x 1.8 tank, the water column a block in
  one corner (x < 1.2, y < 0.5, z < 0.6), so the flow is 3D from the first step;
  TetGen mesh from he at load (poly written by rank 0 only).
- model_baseline_mpi.py solution=PREFIX also writes each model's owned node
  coordinates (MODEL_x), so runs on different partitions can be compared node by
  node. Needed because neither ParMETIS nor PT-Scotch gives the same partition twice
  here (4 ranks: different owned counts run to run), and the 'current' partition
  keeps TetGen's node order: 78 k elements per rank with ghosts for 22 k nodes,
  against 33 k with PT-Scotch.
- model_baseline_mpi.py velocity_block_size=1 gives the Schur fieldsplit's velocity
  split block size nSpace, so BoomerAMG sees 3 functions per node (systems AMG, and
  what -..._pc_hypre_boomeramg_nodal_coarsen needs; proteus leaves bs = 1, which
  ksp_view confirms: rows=11430, bs=3 with the option, no bs without it).

Partition noise floor (ASM, he = 0.15, 4 ranks, two PT-Scotch partitions): velocity
1.3e-5, pressure 3.6e-7 relative, water volume 2.1e-9.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
@cekees @claude
TwoPhaseFlow: RANS2P defaults to the Schur fieldsplit + BoomerAMG, AS…
6c0b4f1 
…M without hypre (behaviour change)

SystemPhysics.modelDict['flow'].preconditioner now defaults to 'auto':
'schur_boomeramg' when PETSc has hypre, 'asm' (with a log line) when it does not.
An explicit 'schur_boomeramg' without hypre is still refused at setup; 'asm' is
still available on request. Applies to TwoPhaseFlow's iterative RANS2P solve
(useSuperlu = False, the SimpleNavierStokes3D smoother); the test decks'
useSuperlu = True direct solves are unchanged. cekees's two conditions (11-05):

1. The 3D RANS2P check. damBreakLevelSet3D (6095102), he = 0.08 (110 k tets,
   22 k nodes, 89 k RANS2P dofs), 4 ranks (PT-Scotch), T = 0.03, M4 Max at load
   55-100, so times are min of 2 interleaved runs and indicative:

   | RANS2P preconditioner | Newton its/step | KSP its/Newton | RANS2P/Newton | solve |
   |---|---|---|---|---|
   | 'asm' (GMRES + ASM, PETSc LU) | 3.0 | 6.0-7.1 | 3333 ms | 60.0 s |
   | 'schur_boomeramg' | 3.0 | 2.8-3.0 | 516 ms | 14.4 s |
   | + velocity block size 3 (BoomerAMG systems) | 3.0 | 3.0-3.1 | 614 ms | 19.5 s |
   | + bs 3 + nodal coarsening 1 | 3.0 | 3.0 | | |

   ASM's cost in 3D is the subdomain LU (2.4-3.3 s of PCSetUpOnBlocks per Newton
   iteration). Velocity, max relative difference from a tight reference (ASM,
   Newton 1e-8): 'asm' 1.5-1.7e-5, 'schur_boomeramg' 9.9e-6, + bs 3 3.9e-6, a tight
   'schur_boomeramg' reference 1.5e-5; that is the floor here (every run gets its
   own PT-Scotch partition and RDLS/MCorr keep their loose tolerances). Water volume
   within 4e-8. Nodal coarsening is out: 1.1e-4, 7x the floor. The block size buys
   no iterations, so proteus keeps bs = 1.
2. The hypre survey (ON_NODE_PERFORMANCE_BASELINE.md §9.3). conda-forge's real
   builds on all four platforms, proteus's feedstock and CI envs, the PETSc fork's
   configure scripts, the colonized prefixes, the HPCMP (Carpenter, Narwhal,
   Blueback, Barfoot) and LSU (Mike, QB) module builds and Spack's py-proteus all
   have hypre. Homebrew's petsc, Spack's petsc without +hypre and a PyPI petsc
   built without options do not: they get 'asm'.

2D (round 2, #198, ceds-model): one damBreakLevelSet step 2.78 s at
24 ranks and 1.26 s at 96 with the fieldsplit, 10.4 s and 4.92 s with 'asm'.

test_rans2p_preconditioner.py: the default is the fieldsplit with hypre and ASM
without (petsc_has_hypre patched), 'asm' on request, an explicit
'schur_boomeramg' without hypre and an unknown name are refused; 5 passed (Mac).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
@cekees cekees changed the title C19 round 3: RDLS on the engine; the 3D RANS2P check for the schur_boomeramg default C19 round 3: RDLS on the engine; RANS2P defaults to Schur + BoomerAMG (ASM without hypre) Nov 5, 2026
cekees added a commit that referenced this pull request Nov 5, 2026
@cekees @claude
Coverage floors: NCLS3P.py 76.2 -> 75.6 (round 5 moved 15 covered sta…
7c2a5d5 
…tements to Transport)

The floor was recorded on main ce43b73 (33f30b1): 522 of 685 statements, 76.2%. #199 (C0 round 5,
76f3e12) replaced NCLS3P's own scatterDOFArray (15 statements, all covered) by the one on OneLevelTransport, so main
c09280c has 507 of 670, 75.67%: the same lines untested, a smaller file. Every pull request on c09280c without #204
(whose staticmethod lines add covered statements to the file) is below the floor: this one's coverage job and #212's
both report NCLS3P 75.7%. Line by line, the NCLS3P lines hit in this branch's run and in #204's run (mapped through
#204's edits) are the same set. The new floor is the measured value rounded down to a tenth.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
@cekees cekees mentioned this pull request Nov 5, 2026
Draft
cekees and others added 2 commits November 5, 2026 23:33
@cekees @claude
Coverage floors: NCLS3P.py 76.2 -> 75.6 for main's scatterDOFArray de…
752005b 
…-duplication (76f3e12)

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
@cekees @claude
TwoPhaseFlow: VOF/NCLS/MCorr/RDLS iterative solver options PETSc neve…
32b650e 
…r read

With useSuperlu = False, five options TwoPhaseFlow sets were left unused
(-options_left): VOF, NCLS and MCorr set `pc_pc_hypre_type` (for
`pc_hypre_type`), RDLS set `pc_pc_asm_type` (for `pc_asm_type`) and `max_it`
(for `ksp_max_it`).

Behaviour: hypre's default type is already boomeramg and no solve reaches
2000 (or PETSc's default 10000) iterations, so only RDLS's ASM type changes:
from PETSc's default 'restrict' to the 'basic' the code always asked for.

Measured: test/TwoPhaseFlow/damBreakLevelSet.py, 2 ranks ('current'
partition), final_time 0.05 (5 steps), superlu=0, model_baseline_mpi.py with
petsc='-options_left', M4 Max:

|                         | before | after |
|---|---|---|
| unused options          | 5 (the five above) | none |
| RDLS KSP its/solve      | 2.4 (2-3) | 2.1 (2-2) |
| RDLS Newton its/step    | 3.0 | 3.0 |
| other models' KSP/Newton | unchanged | unchanged |

Solution, max relative difference after vs before (a rerun of the same code
is identical to 0): RDLS 1.9e-6, NCLS 1.9e-6, VOF 1.3e-6, velocity 1.2e-6,
pressure 1.3e-7, MCorr 1.1e-4 of its max (a correction near zero); water
volume 2.6e-14. All within RDLS's own tolerance (nl_atol 0.1 he). No CI test
runs these options: the test decks set useSuperlu = True.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
cekees added a commit that referenced this pull request Nov 6, 2026
@cekees @claude
Coverage floors: NCLS3P.py 76.2 -> 75.6 (round 5 moved 15 covered sta…
33ad04f 
…tements to Transport)

The floor was recorded on main ce43b73 (33f30b1): 522 of 685 statements, 76.2%. #199 (C0 round 5,
76f3e12) replaced NCLS3P's own scatterDOFArray (15 statements, all covered) by the one on OneLevelTransport, so main
c09280c has 507 of 670, 75.67%: the same lines untested, a smaller file. Every pull request on c09280c without #204
(whose staticmethod lines add covered statements to the file) is below the floor: this PR's coverage job (run
42113031204, NCLS3P the only file under its floor) and #212's both report NCLS3P 75.7%. The new floor is the measured
value rounded down to a tenth. Same hunk as 7c2a5d5 on c0-r6, so the train merges both cleanly.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
@cekees @claude
C19: MCorr on the engine's generic element kernel, against proteus's …
bfa1a68 
…MCorr on damBreak

The fourth model of the two-phase level-set family, and the first on the generic
element kernel (engine.h elementKernelSpace, Richards' path) rather than the
VOF/NCLS/RDLS SUPG template. MCorr (mprans/MCorr.h, useExact = 0, useExact_s = 0:
the damBreakLevelSet configuration) is a steady nonlinear reaction-diffusion problem:

  R_i = sum_q H_s [ r w_i + eps_diff grad u . grad w_i ] dV,
  r = porosity (H_epsH(phi + u) - H_vof), dr = porosity delta_epsD(phi + u),
  eps_X = epsFact_X (useMetrics h_phi + (1 - useMetrics) h_e), H_s = H_0(phi_s).

- engine.h: one hook, bindPoint(q) beside bindElement(eN). A functor with point
  data (here NCLS's phi, VOF's H, the porosity, h_phi, h_e, phi_s) is evaluated as
  C.bindElement(eN).bindPoint(q); a functor without it compiles as before.
  richards_layers.py and richards_nd_test.py print identical output with the
  module built before and after; vof/ncls/rdls_test.py pass with their numbers.
- mcorr.h: MCorrCoefficients (masks: no mass, advection or Hamiltonian; eps_diff
  independent of u; exact Jacobian) and mcorrPointOutputs (q_u, q_r without H_s,
  q_n = grad u / sqrt(1e-8 + |grad u|^2), what MCorr hands the other models).
  With useExact_s = 0 every element is active, so MCorr.py's masking of the
  residual and Jacobian rows is the identity; the cut-cell path is not reproduced.
- vof_capture.py model=mcorr (also q_u, q_r, q_n at the next iterate, since MCorr
  starts from u = 0); mcorr_test.py + mcorr_damBreak_he0.05.npz (733 triangles,
  t = 0.063, 125 elements within eps of the interface); ctest proto_mcorr_test.

Engine vs proteus (max relative difference):
  residual 6.9e-15, Jacobian 3.3e-16, q_r 5.4e-15 (q_u = q_n = 0 at u = 0 in both),
  every DOF and element active in both, next Newton iterate 2.5e-14 (MCorr has no
  line search), residual at the next iterate 3.2e-15 of the first residual
  (4.6e-19 absolute; Newton cut it 1e4, to 2.6e-8), q_u/q_r/q_n there 1.1e-14.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
cekees added a commit that referenced this pull request Nov 6, 2026
@cekees @claude
Coverage floors: NCLS3P.py 76.2 -> 75.6 (round 5 moved 15 covered sta…
9d52853 
…tements to Transport)

The floor was recorded on main ce43b73 (33f30b1): 522 of 685 statements, 76.2%. #199 (C0 round 5,
76f3e12) replaced NCLS3P's own scatterDOFArray (15 statements, all covered) by the one on OneLevelTransport, so main
c09280c has 507 of 670, 75.67%: the same lines untested, a smaller file. Every pull request on c09280c without #204
(whose staticmethod lines add covered statements to the file) is below the floor: this one's coverage job and #212's
both report NCLS3P 75.7%. Line by line, the NCLS3P lines hit in this branch's run and in #204's run (mapped through
#204's edits) are the same set. The new floor is the measured value rounded down to a tenth.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
(cherry picked from commit 7c2a5d5)
@cekees @claude
Revert "TwoPhaseFlow: RANS2P defaults to the Schur fieldsplit + Boome…
be4b929 
…rAMG, ASM without hypre"

This reverts 6c0b4f1. The 96-rank check on ceds-model (247 k nodes, 2 steps,
quiet machine) found the default fieldsplit above the solution floor at a few
wall nodes, and the repeat run shows it is not partition noise.

Measured (damBreakLevelSet3D, he 0.035, T 0.015; max velocity difference relative
to max |vel| = 0.643, node by node, against the other fieldsplit run at 96 ranks
with bs 3 + nodal coarsening, which agrees with the 24-rank runs to 5e-5):

|                           | worst node          | nodes > 1e-4 | rms   |
|---|---|---|---|
| fieldsplit, 96 ranks, run 1 | 3.3e-4 at (1.22, 0.547, 0) | 7 | 1.4e-6 |
| fieldsplit, 96 ranks, repeat | 3.2e-4 at (1.22, 0.547, 0) | 5 | 1.2e-6 |
| ASM, 96 ranks            | 6.1e-5 (uniform rms 1.2e-5 from all others) | 0 | |
| fieldsplit, 24 ranks, tight reference | 5.0e-5 | 0 | 4.8e-7 |

The two 96-rank runs had different partitions (7 of 247 482 nodes with the same
owner), yet the worst node and the second (1.236, 0.534, 0) recur at the same size
(3.3/3.2e-4, 1.7e-4) and the difference fields correlate at 0.56. All on the floor
z = 0, y = 0.53-0.57. The Mac's floor was 1.5e-5 (4 ranks, he 0.08), 6e-5 here
between the other runs. The default 'asm' is not changed by this revert.

What the revert gives up: the fieldsplit solves 7.5x faster at 24 ranks (55 s vs
416 s) and 1.8x at 96 (44 s vs 80 s), with KSP its flat in ranks (2.8); it stays
available as preconditioner = 'schur_boomeramg'. The numbers and raw data:
ON_NODE_PERFORMANCE_BASELINE.md §9.2, c19_r3_data/ceds/r3/. Whether the 96-rank
wall-node difference comes from the fieldsplit or from ASM's looser convergence
at 96 ranks needs a tight ASM reference there; the default flip can return after.

test_rans2p_preconditioner.py is back to its previous 4 tests (4 passed, Mac).
The option-typo commit 32b650e is kept.

Co-Authored-By: Claude Sonnet 5.5 <noreply@anthropic.com>
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Reviewers

No reviews

Assignees

No one assigned

Labels

None yet

Projects

None yet

Milestone

No milestone

Development

Successfully merging this pull request may close these issues.

1 participant

@cekees