Require the PCD Earth radius for spherical meshes by xylar · Pull Request #547 · E3SM-Project/Omega

xylar · 2026-09-11T21:54:22Z

Require the sphere radius of a spherical input mesh to match the Physical Constants Dictionary Earth radius to a relative 1e-8 instead of 1e-4, now that the icosahedral CTest mesh has been regenerated on that sphere. Fix a spherical test field that was singular at the poles, add a check in HorzMeshTest that the mesh coordinates, the sphere_radius attribute and REarth agree, and document the requirement.

Fixes #346.

Please review specifically:

The layer thickness field scalarB in the spherical tracer diffusion test of TendencyTermsTest (from Add tracer tendency terms and auxiliary variables #140, @brian-oneill) has no limit at the poles, so the exact tendency at a polar cell depended on that cell's arbitrary longitude; the regenerated mesh exposed this (LInf error 0.049 to 0.235, all of it in the two polar pentagons). The new field is smooth on the sphere and the error is 0.016 on either mesh (c25835d):
```
old:  scalarB = 2 + cos(lon) sin(lat)
new:  scalarB = 2 + cos(lon) cos(lat) sin(lat)   (= 2 + x z / R^2)
```
The radius check compares the file value with pcd::mean_radius in double precision before the cast to Real, so a single-precision build applies the same tolerance.

The CTests need cosine_bell_icos480.omega_vars.260911.nc as OmegaSphereMesh.nc; the companion Polaris PR E3SM-Project/polaris#779 points omega_ctest.py at it. The quick start also now links the planar mesh PlanarPeriodic48x48.omega_vars.260825.nc (the one with vector reconstruction weights) that omega_ctest.py has downloaded since E3SM-Project/polaris#721; the Omega docs were not updated at the time. The expected error norms of the spherical tests are updated for the new mesh, which has the same connectivity as the old one but slightly different cell positions.

Checklist

🤖 Generated with Claude Code

xylar · 2026-09-11T21:54:34Z

Testing

CTests on Chrysalis with Intel and OpenMPI, Release, built from 6accad7772 with OmegaSphereMesh.nc linked to cosine_bell_icos480.omega_vars.260911.nc:

/gpfs/fs1/home/ac.xylar/e3sm_work/polaris/update-omega-ctest-sphere-mesh/build_omega/build_chrysalis_intel

All 49 tests pass. With the previous 260807 mesh the five spherical tests abort on the radius check, as intended. The Polaris omega_pr suite has not been run yet.

Posted by Claude Code on @xylar's behalf. The testing, analysis and wording above are AI-authored; please check them accordingly.

xylar · 2026-09-12T07:20:59Z

Testing

I ran the Polaris omega_pr suite on Chrysalis with Intel and OpenMPI from E3SM-Project/polaris#779, building Omega from 6accad7772, and compared against a baseline built from this PR's base d086d4984d with the same Polaris branch:

/lcrc/group/e3sm/ac.xylar/polaris_1.1/chrysalis/test_20260912/omega_pr_baseline_develop_d086d49
/lcrc/group/e3sm/ac.xylar/polaris_1.1/chrysalis/test_20260912/omega_pr_earth_radius_6accad7

23 of 24 tasks pass and are bit for bit with the baseline, including the spherical cosine_bell and QU.240km tasks. The remaining task, barotropic_channel/short, fails its baseline comparison for a Polaris reason unrelated to this PR (E3SM-Project/polaris#780); its other variables are bit for bit.

Posted by Claude Code on @xylar's behalf. The testing, analysis and wording above are AI-authored; please check them accordingly.

xylar · 2026-09-12T07:26:15Z

@brian-oneill, it looks like you are the author of scalarB. Are you okay with redefining it so it's smooth at the pole?

mwarusz

Looks good, but I have a few comments.

mwarusz · 2026-09-14T21:00:48Z

                          cos(LatEdge(IEdge)) * cos(Pi * TSim / Tau);
-             const R8 normalVel = REarth * (u * cos(AngleEdge(IEdge)) +
-                                            v * sin(AngleEdge(IEdge)));
+             const R8 normalVel = SphereRadius * (u * cos(AngleEdge(IEdge)) +


If we require REarth to be very close to HorzMesh.SphereRadius why not do the opposite, i.e. use REarth everywhere in the code. REarth is simpler to get, since it is a global constant.

Agreed. With the check at 1e-8 the mesh value adds nothing, so cc31158 uses REarth in the operators, the prescribed velocities and the test setups. The coupled driver still reads HMesh->SphereRadius; I left it since it is outside this PR.

Posted by Claude Code on @xylar's behalf. The testing, analysis and wording above are AI-authored; please check them accordingly.

The coupled driver still reads HMesh->SphereRadius; I left it since it is outside this PR.

@andrewdnolan, perhaps consider whether this should be changed as well.

mwarusz · 2026-09-14T21:10:25Z

+   // scalarB is differentiated by the tracer diffusion test, so it must be
+   // smooth on the sphere: cos(Lon)*cos(Lat)*sin(Lat) = x*z/R^2 vanishes at
+   // the poles, whereas cos(Lon)*sin(Lat) has no limit there and makes the
+   // exact diffusion at a polar cell depend on the arbitrary longitude of
+   // that cell's center.


LLM comment talking about what the code used to do. Please rephrase or remove it. Every test field defined in this file was always supposed to be smooth.

Reworded in 0bd7df2 to state only the requirement.

Posted by Claude Code on @xylar's behalf. The testing, analysis and wording above are AI-authored; please check them accordingly.

mwarusz · 2026-09-14T21:16:54Z

+      if (std::abs((pcd::mean_radius - TmpRadius) / pcd::mean_radius) >
+          RadiusRelTol)
+         ABORT_ERROR("Input mesh sphere radius {} does not match the Earth "
+                     "radius {} from the Physical Constants Dictionary "
+                     "(relative tolerance {})",
+                     TmpRadius, pcd::mean_radius, RadiusRelTol);


This is a very minor suggestion, but I think we should get in the habit of stating requirements using OMEGA_REQUIRE.

Suggested change

if (std::abs((pcd::mean_radius - TmpRadius) / pcd::mean_radius) >

RadiusRelTol)

ABORT_ERROR("Input mesh sphere radius {} does not match the Earth "

"radius {} from the Physical Constants Dictionary "

"(relative tolerance {})",

TmpRadius, pcd::mean_radius, RadiusRelTol);

OMEGA_REQUIRE(std::abs((pcd::mean_radius - TmpRadius) / pcd::mean_radius) <

RadiusRelTol,

"Input mesh sphere radius {} does not match the Earth "

"radius {} from the Physical Constants Dictionary "

"(relative tolerance {})",

TmpRadius, pcd::mean_radius, RadiusRelTol);

Done in 427f0cf.

Posted by Claude Code on @xylar's behalf. The testing, analysis and wording above are AI-authored; please check them accordingly.

xylar · 2026-09-15T10:04:06Z

Testing

I ran the Polaris omega_pr suite on Chrysalis with gnu and OpenMPI, using an Omega build from 427f0cf0e8 and Polaris main (ea15edbdbd):

/lcrc/group/e3sm/ac.xylar/polaris_1.1/chrysalis/test_20260915/earth_radius_omega_pr

All 23 tasks pass and are bit-for-bit with the main Omega baseline (ee6397a627) at /lcrc/group/e3sm/ac.xylar/polaris_1.1/chrysalis/test_20260911/omega_pr_baseline_main_ea15edb. The property-check failures in column/vmix_unstable and column/ekman are the same in the baseline run.

Posted by Claude Code on @xylar's behalf. The testing, analysis and wording above are AI-authored; please check them accordingly.

xylar · 2026-09-16T16:21:17Z

And @mwarusz, thanks very much for your review!

mwarusz

Thanks for addressing my comments @xylar . Approving based on successful unit testing on frontier with craygnu and craygnu-mphipcc. I have one more comment that may be useful to other people looking at his PR.

mwarusz · 2026-09-16T19:53:56Z

- [Global Mesh](https://web.lcrc.anl.gov/public/e3sm/polaris/ocean/omega_ctest/cosine_bell_icos480.omega_vars.260807.nc)
- [Planar Mesh](https://web.lcrc.anl.gov/public/e3sm/polaris/ocean/omega_ctest/PlanarPeriodic48x48.omega_vars.260720.nc)
+- [Global Mesh](https://web.lcrc.anl.gov/public/e3sm/polaris/ocean/omega_ctest/cosine_bell_icos480.omega_vars.260911.nc)
+- [Planar Mesh](https://web.lcrc.anl.gov/public/e3sm/polaris/ocean/omega_ctest/PlanarPeriodic48x48.omega_vars.260825.nc)


It seems that this PR also updates the planar test mesh with vector reconstruction weights. I'm not opposed to the change, but wanted to highlight it for others, since it isn't mentioned anywhere in the PR description. I was confused why this PR needs to update this mesh.

The update here is merely that the update (already made and documented previously in Polaris) was somehow missing from the Omega. I will follow up with details.

Good catch; the description now says so. The planar mesh itself is not new to this PR: omega_ctest.py has downloaded PlanarPeriodic48x48.omega_vars.260825.nc since E3SM-Project/polaris#721 (2026-08-26), which added the reconstruction weights, but the Omega quick start was not updated then. This PR only brings the link into line with what the utility already uses.

Posted by Claude Code on @xylar's behalf. The testing, analysis and wording above are AI-authored; please check them accordingly.

brian-oneill · 2026-09-23T13:51:13Z

Testing

I ran the ctests successfully (using the new OmegaSphereMesh.nc) on Frontier with craygnu and craygnu-mphipcc, on pm-cpu with intel, and on pm-gpu with gnugpu.

brian-oneill

The function I had initial used in TendencyTermsTest was chosen arbitrarily, so I'm fine with a more robust choice of test function. Everything here looks good to me. Approved.

xylar · 2026-09-23T14:45:41Z

Thanks @brian-oneill! I'll rebase to fix the conflict right now.

SecondDerivativeOnCell and the prescribed nondivergent and divergent velocities used the REarth constant to scale spherical geometry, while the rest of the mesh geometry comes from the input file. Take the radius from HorzMesh::SphereRadius instead, so a single value describes the sphere the mesh lies on; HorzMesh checks it against REarth on read. DetermineSphericalPatchGeometry now receives the radius as an argument, and the unit test that projects a synthetic hexagonal patch onto a sphere passes REarth explicitly, since it builds its own sphere and also runs on the planar mesh. Also correct the comment on the divergent-flow period (12 days, not 14). Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

The spherical test setups evaluated their analytic fields with REarth while the mesh they compared against could lie on a slightly different sphere; the CTest icosahedral mesh was 3.3e-5 larger, which contributed a systematic error of that size to every spherical operator test. Read the radius from the default HorzMesh instead, which resolves the long-standing TODO in three of the setups. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

HorzMeshTest verified that cell centers are equidistant from the origin but never compared that distance with the SphereRadius attribute read from the file, nor with the REarth constant the model assumes. Add both checks, with a tolerance that follows the build precision since the coordinates are stored as Real. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

The relative tolerance on the input mesh's sphere radius was 1e-4 so that the CTest icosahedral mesh, generated with an older Earth radius, would still load. That mesh has been regenerated on the sphere of the Physical Constants Dictionary, so require agreement to 1e-8. Compare the file value with pcd::mean_radius in double precision, before the cast to Real, so that a single-precision build applies the same test. Closes E3SM-Project#346 Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

TestSetupSphere::scalarB was 2 + cos(Lon) sin(Lat), which has no limit at the poles: its value there depends on the direction of approach. The tracer diffusion test differentiates it, so the exact tendency at a polar cell depended on the arbitrary longitude assigned to that cell's center, through a term 2 cos^3(Lon)/R. On the old CTest mesh the polar pentagons happened to sit near Lon = 58 and 257 degrees, where the term is small (LInf 0.049); on the regenerated mesh they sit near 20 and 191 degrees and the LInf error rose to 0.235, all of it in those two cells. Use 2 + cos(Lon) cos(Lat) sin(Lat) = 2 + x z / R^2 instead, which is smooth on the sphere, and update the exact divergence accordingly. The error is then 0.0156 on either mesh, with its maximum at 30 degrees latitude as ordinary truncation error. The other uses of scalarB on the sphere are pointwise at edges, away from the poles, and unaffected. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

The icosahedral CTest mesh, cosine_bell_icos480.omega_vars.260911.nc, was regenerated by Polaris on the sphere with the Physical Constants Dictionary radius. The mesh has the same connectivity as before but JIGSAW places the cells up to 7e-4 (relative) differently, and the analytic test fields are now evaluated with the mesh radius, so the error norms of the spherical operator, auxiliary variable, tendency and forcing tests move by up to a few tenths of a percent (2 percent for those with a radius-scaled analytic solution). Record the new values. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

The MPAS mesh page said there was no requirement on the sphere radius, which has not been true since the model started checking it. State that spherical meshes must use the Physical Constants Dictionary Earth radius, and point to the Polaris converter for older MPAS-Ocean meshes. Mention the check on the horizontal mesh page, and update the quick start links to the regenerated icosahedral mesh and the current planar mesh. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

With the input mesh now required to lie on the sphere of radius REarth to 1e-8, carrying HorzMesh::SphereRadius into SecondDerivativeOnCell, the prescribed velocities and the test setups adds plumbing without adding information. Use the REarth constant directly, as suggested in review. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

The FCT low-order flux test built its reference solution from a single hardcoded edge length per geometry, 0.012028130608117 on the plane and 277371 m on the sphere. The spherical value was the mean dvEdge of an earlier icosahedral CTest mesh, so the test measured how far each edge departs from that mean and its expected errors had to be re-derived whenever the mesh changed. Use Mesh->DvEdge instead. The low-order flux of a uniform tracer is then exactly DvEdge times the normal velocity, and the expected errors are zero on both meshes. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>

The spherical FCT test advects with a constant lon/lat velocity, so the provisional layer thickness error is dt |div u| = 3 dt tan(lat) / R and scales with 1/R. On the regenerated icosahedral CTest mesh, with the Physical Constants Dictionary radius, the L2 errors of the provisional thickness and its inverse rise by 2.6e-5 relative (3.3e-5 from the radius alone) and the LInf errors by 1e-5 to 7e-5. Record the new values. The FCTHNewInv check still yields inf on every mesh; that is separate from the mesh change and left alone here. Co-Authored-By: Claude Opus 5.5 (1M context) <noreply@anthropic.com>

xylar · 2026-09-23T16:18:10Z

Rebased onto develop to pick up FCT (#447). Two new commits: the spherical FCT low-order flux test now takes edge lengths from the mesh instead of a hardcoded value from the old mesh, and the spherical FCT provisional-thickness errors are updated for the regenerated mesh.

Testing

All 50 CTests pass on Chrysalis with Intel and the 260911 mesh. That build needed a local workaround for #572. The FCT test hides some failures that are unrelated to this PR; see #573.

Posted by Claude Code on @xylar's behalf. The testing, analysis and wording above are AI-authored; please check them accordingly.

sbrus89 · 2026-09-24T22:17:59Z

I did a local merge and I'm seeing a bunch of CTest fails with gnugpu on pm-gpu:

6: aborting job:
6: application called MPI_Abort(MPI_COMM_WORLD, -1) - process 6
6: Kokkos::Cuda ERROR: Failed to call Kokkos::Cuda::finalize()

I don't see how these could be related to this PR.

sbrus89 · 2026-09-24T22:21:16Z

CTests on pm-cpu (gnu) pass and the omega_pr suite passes with both gnu and gnugpu on pm-cpu/pm-gpu.

Polaris `omega_pr` suite

Baseline workdir: /global/homes/s/sbrus/scratch/polaris_baseline_omega_pr_omega_gnu_mpich/
Baseline build: /global/homes/s/sbrus/scratch/polaris_baseline_omega_pr_omega_gnu_mpich/ocean_build
PR build: /global/homes/s/sbrus/scratch/polaris_earth_radius_omega_pr_omega_gnu_mpich/ocean_build
PR workdir: /global/homes/s/sbrus/scratch/polaris_earth_radius_omega_pr_omega_gnu_mpich
Machine: pm-cpu
Compiler: gnu
Build type: Release
Log: /global/homes/s/sbrus/scratch/polaris_earth_radius_omega_pr_omega_gnu_mpich/polaris_omega_pr.o58832465
Result: All tests passed

Polaris `omega_pr` suite

Baseline workdir: /global/homes/s/sbrus/scratch/polaris_baseline_omega_pr_omega_gnugpu_mpich/
Baseline build: /global/homes/s/sbrus/scratch/polaris_baseline_omega_pr_omega_gnugpu_mpich/ocean_build
PR build: /global/homes/s/sbrus/scratch/polaris_earth_radius_omega_pr_omega_gnugpu_mpich/ocean_build
PR workdir: /global/homes/s/sbrus/scratch/polaris_earth_radius_omega_pr_omega_gnugpu_mpich
Machine: pm-gpu
Compiler: gnugpu
Build type: Release
Log: /global/homes/s/sbrus/scratch/polaris_earth_radius_omega_pr_omega_gnugpu_mpich/polaris_omega_pr.o58832467
Result: All tests passed

xylar · 2026-09-25T07:13:27Z

@sbrus89 That abort is a mesh file Omega could not read, not the radius check this PR tightens.

MPI_Abort(MPI_COMM_WORLD, -1) comes from SCORPIO's piodie. Omega's Error::abort uses ErrorCode::Critical, which is 3. Further up your log there should be a PIO: FATAL ERROR: Aborting... Opening file (...) failed line naming the file — could you post that, and ls -lL test/Omega*.nc from your build? My guess is stale or dangling mesh symlinks from an earlier build directory.

You do separately need E3SM-Project/polaris#779 for OmegaSphereMesh.nc. Without it the CTests link the old icos480 mesh and exactly five tests fail, but with a self-explanatory message rather than a bare abort, so that is not what you hit:

[critical] [HorzMesh.cpp:195] Input mesh sphere radius 6371220 does not match the Earth radius 6371008.7714 from the Physical Constants Dictionary (relative tolerance 1e-08)

Worth checking both, so that relinking the mesh does not just leave you on the IO failure.

pm-gpu, gnugpu, mpich

With polaris#779 and this branch at 1a564eb:

CTests:   100% tests passed, 0 tests failed out of 50
omega_pr: 25 tasks, PASS: All passed successfully! (0:14:07)

/pscratch/sd/x/xylar/polaris_1.1/pm-gpu/earth_radius_20260924

Posted by Claude Code on @xylar's behalf. The testing, analysis and wording above are AI-authored; please check them accordingly.

sbrus89 · 2026-09-25T14:54:04Z

Yep, I was testing multiple things at a time and forgot to link the meshes for this one. CTests pass on pm-gpu (gnugpu) Thanks @xylar.

sbrus89 · 2026-09-25T17:17:52Z

Passes e3sm_omega_developer on pm-cpu with gnu.

xylar mentioned this pull request Sep 12, 2026

Use the Icos480 CTest mesh regenerated on the PCD sphere E3SM-Project/polaris#779

Merged

1 task done

xylar requested review from brian-oneill and mwarusz September 12, 2026 06:52

mwarusz reviewed Sep 14, 2026

View reviewed changes

mwarusz approved these changes Sep 16, 2026

View reviewed changes

xylar mentioned this pull request Sep 21, 2026

CTests fail for current HEAD on chrysalis #565

Closed

brian-oneill approved these changes Sep 23, 2026

View reviewed changes

xylar assigned sbrus89 Sep 23, 2026

xylar and others added 12 commits September 23, 2026 09:54

State the smoothness requirement on scalarB without history

d8fecb8

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

Express the sphere radius check with OMEGA_REQUIRE

630d17f

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>

xylar force-pushed the omega/earth-radius branch from 427f0cf to 1a564eb Compare September 23, 2026 16:18

sbrus89 merged commit 272e03d into E3SM-Project:develop Sep 25, 2026
1 check passed

andrewdnolan mentioned this pull request Sep 25, 2026

Add MOAB coupler support for Omega #523

Merged

4 of 26 tasks

xylar deleted the omega/earth-radius branch September 26, 2026 07:50

This was referenced Sep 28, 2026

Add MOC stream function analysis capability #481

Merged

Update Omega submodule E3SM-Project/polaris#828

Merged

Conversation

xylar commented Sep 11, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

xylar commented Sep 11, 2026

Testing

Uh oh!

xylar commented Sep 12, 2026

Testing

Uh oh!

xylar commented Sep 12, 2026

Uh oh!

mwarusz left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

xylar commented Sep 15, 2026

Testing

Uh oh!

xylar commented Sep 16, 2026

Uh oh!

mwarusz left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

brian-oneill commented Sep 23, 2026

Testing

Uh oh!

brian-oneill left a comment

Choose a reason for hiding this comment

Uh oh!

xylar commented Sep 23, 2026

Uh oh!

xylar commented Sep 23, 2026

Testing

Uh oh!

sbrus89 commented Sep 24, 2026

Uh oh!

sbrus89 commented Sep 24, 2026

Polaris omega_pr suite

Polaris omega_pr suite

Uh oh!

xylar commented Sep 25, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

pm-gpu, gnugpu, mpich

Uh oh!

sbrus89 commented Sep 25, 2026

Uh oh!

sbrus89 commented Sep 25, 2026

Uh oh!

Uh oh!

Reviewers

Uh oh!

Assignees

Labels

Projects

Milestone

Development

Uh oh!

4 participants

xylar commented Sep 11, 2026 •

edited

Loading

Polaris `omega_pr` suite

Polaris `omega_pr` suite

xylar commented Sep 25, 2026 •

edited

Loading