Skip to content

Add reaction-based spherical-shell geoid projection - #686

Open
gthyagi wants to merge 1 commit into
underworldcode:developmentfrom
gthyagi:feature/spherical-reaction-projection
Open

Add reaction-based spherical-shell geoid projection#686
gthyagi wants to merge 1 commit into
underworldcode:developmentfrom
gthyagi:feature/spherical-reaction-projection

Conversation

@gthyagi

@gthyagi gthyagi commented Sep 5, 2026

Copy link
Copy Markdown
Contributor

Summary

Supersedes #646 with a single, narrowly scoped commit from current
development (9da04b72). The original PR accumulated unrelated mantle
convection development; this replacement contains only the outstanding
spherical-shell reaction-projection changes.

The generic Stokes.boundary_normal_traction_integral() implementation is
already upstream through #648. This PR reuses that API rather than adding a
second implementation or changing shared solver/Cython code.

No SUPG/SLCN, checkpoint, evaluator, scalar heat-flux, or solver-configuration
changes are included. The diff is five files: one postprocessing module, two
test files, the subsystem documentation, and the changelog.

Why and How

A nodal reaction represents an integrated force, not a pointwise traction.
When only a harmonic coefficient is needed, recovering pointwise P2 values and
then reconstructing a surface triangulation is unnecessary. It also encounters
the curved-boundary P2 issues discussed in #414, #431, and #633/#649.

For the scalar normal-reaction load R and nodal samples of the harmonic
phi = P_l(cos(theta)), the existing upstream integral API evaluates:

I(phi) = phi^T R - (1^T R / A) * integral_boundary(phi)

The adapter divides by BdIntegral(phi**2) and applies its existing boundary
sign and buoyancy scale. Numerator and denominator therefore use the same
faceted finite-element boundary geometry. Reaction ownership and reductions
remain the responsibility of the upstream API.

The reaction path needs no global boundary-sample gather, rank-zero
ConvexHull, or pointwise traction recovery. Consumers needing a nodal
topography field continue to use the existing pointwise APIs.

Fresh Validation of This Branch

Tested commit: fd0a3bcc. A wheel was built from this clean branch and installed
into an isolated target using the existing Pixi dependencies. PYTHONPATH
selected that installation for the tests; its imported UW3 had the upstream
reaction API and no AdvDiffusionSUPG class. The installed geoid.py checksum
matches the branch source.

Environment: macOS/arm64, Python 3.12, PETSc 3.25.0, Open MPI 5.0.10,
MPI-enabled h5py.

Test Result Time
test_1070_postprocessing_geoid.py 13 passed, including explicit reaction projection and invalid-selector coverage Included below
test_1120_SLVectorCartesian.py All 3 cases passed, including the case previously flagged in #646 27.11 s for both serial files
parallel/test_1071_spherical_shell_geoid_parallel.py Passed on all 4 MPI ranks; eight response coefficients checked against Zhong Table 2 Approximately 5.18 s

The MPI comparison allows near-machine-precision reduction differences
(rtol=1e-13, atol=1e-14) rather than requiring bitwise equality.
Serial output contained three non-interactive plotting warnings and no failures.

Separately, the synchronized mantle integration branch passed strict
Level 1 selection: 1676 passed, 36 skipped, 2 expected failures in 1053.89 s.
That broader check is not presented as a full Level 1 run of this standalone
branch. Linux CI remains to be verified.

Build, extraction, and validation report.

Historical Production Evidence

These are the controlled Gadi results already reported in #646, not new
measurements of this extracted commit
. The layered Zhong case used degree 5,
load depth 0.5, cellsize=1/64, a fitted 10^4 viscosity lid, P2P1, rotated
free slip, 192 ranks, and Stokes tolerance 1e-5.

Self-gravity response Centroid error Reaction error
Surface topography -2.62% +0.08%
CMB topography +0.91% <0.01% in magnitude
Surface geoid -4.01% +0.12%
CMB geoid +1.19% <0.01% in magnitude
Stage Centroid Reaction
Stokes solve 784.55 s 772.01 s
Response integrals, topography, and geoid 120.14 s 36.92 s

The response stage was 3.25 times faster. It includes common velocity/divergence
integrals, not geoid evaluation alone. The small Stokes timing change is
run-to-run variation; no Stokes speedup or whole-job memory reduction is claimed.
Velocity is unchanged by projection, including the approximately 8.46% residual
error in this case's small surface velocity.

Production benchmark
and historical results.

Review

@lmoresi: this isolates the spherical adapter from the unrelated changes you
flagged in #646. Please review its use of the already-merged reaction integral
and matching boundary norm. No generic reaction-recovery implementation is
being replaced. The mantle branch remains the development branch and is not
deleted when the superseded PR is closed.

…ntegral

Extract only the outstanding spherical postprocessing changes from PR underworldcode#646 onto upstream development at 9da04b7. Reuse boundary_normal_traction_integral already provided by underworldcode#648; no solver, Cython, SLCN/SUPG, checkpoint, or evaluation changes are included.

Add opt-in projection=reaction to contract the normal reaction with the axisymmetric harmonic and normalize by its matching finite-element boundary inner product. Preserve centroid as the compatibility default.

Carry the projection-validation and Zhong serial/MPI regressions, including near-machine-precision rank comparisons, and update the postprocessing documentation and changelog. Post-extraction build and regression results are recorded separately in the benchmark repository.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant