Skip to content

Add vertical advection#328

Open
brian-oneill wants to merge 43 commits intoE3SM-Project:developfrom
brian-oneill:omega/add-vertical-advection
Open

Add vertical advection#328
brian-oneill wants to merge 43 commits intoE3SM-Project:developfrom
brian-oneill:omega/add-vertical-advection

Conversation

@brian-oneill
Copy link

@brian-oneill brian-oneill commented Jan 11, 2026

Add the VertAdv module to Omega, providing data structures and methods for vertical advection. This module includes:

  • Arrays for storing vertical velocity and tracer fluxes at layer interfaces.
  • A method to compute the vertical velocity from the divergence of the horizontal velocity.
  • Methods for computing the tendencies of thickness, horizontal velocity, and tracers due to vertical advection. The tracer tendency methods include multiple options for order of accuracy of the flux calculation, as well as algorithmic options for standard advection or flux-corrected transport (FCT).

Additional changes in this PR:

  • Extend parallelForOuter with an optional argument to enable the use of scratch arrays in hierarchical parallelism.
  • Add the auxiliary variable ProvThickness, which represents thickness after horizontal thickness flux. Further work required to fully implement provisional thickness.
  • Remove BottomDepth from HorzMesh, this field was previously added to VertCoord.

This PR builds and passes the unit tests successfully on Chrysalis, pm-cpu, pm-gpu, and Frontier (CPU & GPU).

Checklist

  • Documentation:
  • Linting
  • Building
    • CMake build does not produce any new warnings from changes in this PR
  • Testing
    • Add a comment to the PR titled Testing with the following:
      • Which machines CTest unit tests
        have been run on and indicate that are all passing.
      • The Polaris omega_pr test suite
        has passed, using the Polaris e3sm_submodules/Omega baseline
      • Document machine(s), compiler(s), and the build path(s) used for -p for both the baseline (Polaris e3sm_submodules/Omega) and the PR build
      • Indicate "All tests passed" or document failing tests
      • Document testing used to verify the changes including any tests that are added/modified/impacted.
    • New tests:
      • CTest unit tests for new features have been added per the approved design.
      • Polaris tests for new features have been added per the approved design (and included in a test suite)

Fixes #338

@mwarusz
Copy link
Member

mwarusz commented Jan 15, 2026

When running VERTADV_TEST with bounds checking enabled I get out-of-bounds errors. The reason is VertCoord::init(false) gets NVertLayers from OmegaMesh.nc and this value is used to create the member arrays of VertAdv in VertAdv::init(). The value of NVertLayers in VertAdv and MaxLayerCell get overwritten later in the convergence test, but the VertAdv member arrays (e.g. VertFlux) don't get resized.


for (int KVec = 0; KVec < KLen; ++KVec) {
const I4 K = KStart + KVec;
const Real WAvg = 0.5_Real * (LocTotVertVelocity(Cell1, K) +
Copy link

@cbegeman cbegeman Jan 15, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do we want to use

InterpCellToEdge(const HorzMesh *Mesh);

so that we can easily switch from anisotropic to isotropic interpolation for all variables (here and elsewhere)?

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We can’t use InterpCellToEdge here as is, since it takes an Array1DReal as input, while TotalVerticalVelocity is an Array2DReal. Supporting this would require either an overload or a templated version of the interpolator. It may make more sense to handle this outside this PR and identify any other places to implement this refactor.

@brian-oneill brian-oneill marked this pull request as ready for review February 14, 2026 01:01
@brian-oneill
Copy link
Author

Polaris omega_pr suite

  • Baseline workdir: /lcrc/group/e3sm/ac.boneill/polaris/intelTest/baseline_omega_pr/
  • Baseline build: /lcrc/group/e3sm/ac.boneill/dev/intelCPU
  • PR build: /lcrc/group/e3sm/ac.boneill/myFork/intelBuild
  • PR workdir: /lcrc/group/e3sm/ac.boneill/polaris/intelTest/new_omega_pr
  • Machine: chrysalis
  • Partition: compute
  • Compiler: intel
  • Build type: Release
  • Log: /lcrc/group/e3sm/ac.boneill/polaris/intelTest/new_omega_pr/polaris_omega_pr.o1152821
  • Result: All tests passed

Copy link
Collaborator

@mark-petersen mark-petersen left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for your documentation throughout. I fixed a few spelling errors.

@mark-petersen
Copy link
Collaborator

See compiled documentation here:

@hyungyukang
Copy link

Using my time-stepping working branch (https://github.com/hyungyukang/Omega/tree/hkang/omega/merge-vert-vel
), I tested this PR with the merry-go-round test in Polaris. Since the test does not currently support Omega directly, I made several adjustments to run it manually. I will document the procedure for running the test manually in a separate comment later.

As expected, I obtained a convergence rate of 1.971 with VerticalTracerFluxOrder: 2, which matches the convergence rate obtained in MPAS-Ocean:

image

} else if (VAdv->VertAdvChoice == VertAdvOption::FCT) {
ThicknessForVAdv = AuxState->LayerThicknessAux.ProvThickness;
}
VAdv->computeTracerVAdvTend(TracerTend, TracerArray, ThicknessForVAdv,

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
VAdv->computeTracerVAdvTend(TracerTend, TracerArray, ThicknessForVAdv,
VAdv->computeTracerVAdvTend(LocTracerTend, TracerArray, ThicknessForVAdv,

Copy link

@hyungyukang hyungyukang left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I’m approving this PR based on my testing and visual inspection. All the changes I requested have been addressed, and the code works as expected.

@brian-oneill, thank you very much for your hard work on this!

Copy link
Collaborator

@mark-petersen mark-petersen left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Approving based on Hyun's convergence test above, and visual inspection. Successfully passed all cTests on perlmutter CPU and GPU with gnu compilers.

Copy link

@cbegeman cbegeman left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I ran the 5m merry-go-round case with Omega (above) and MPAS-Ocean (below) and the error patterns look very similar. The horizontal velocity for Omega is not yet available, as it relies on reconstruction to the x-direction.

image image

@cbegeman
Copy link

Great work, @brian-oneill! I just wanted to check; has anyone yet verified the 3rd and 4th order convergence options? If not, I'll do this with the polaris merry-go-round test.

@hyungyukang
Copy link

hyungyukang commented Feb 23, 2026

Great work, @brian-oneill! I just wanted to check; has anyone yet verified the 3rd and 4th order convergence options? If not, I'll do this with the polaris merry-go-round test.

@cbegeman , I ran it last time with the 3rd convergence option. However, I was not able to obtain the 3rd order convergence rate for both Omega (1.981) and MPAS-Ocean (1.981). Could you please run it as well just in case I missed something?

@cbegeman
Copy link

Using one processor, I was able to reproduce @hyungyukang's results:
image

Testing whether this is fixed by #316

@cbegeman
Copy link

4th order advection with limiter (MPAS top, Omega bottom):
image
image
image
image

@cbegeman
Copy link

4th order advection without limiter (MPAS top, Omega bottom):
image
image

image image

brian-oneill and others added 26 commits March 6, 2026 08:47
This PR adds ProvThickness to LayerThicknessAux, which is
needed by the tracer equations. PR E3SM-Project#335 removes SSH from
LayerThicknessAux and removes computeVarsOnCells from
computeMomAux.
@sbrus89 sbrus89 force-pushed the omega/add-vertical-advection branch from 10eec48 to fbf582f Compare March 6, 2026 20:01
@sbrus89 sbrus89 self-assigned this Mar 6, 2026
@sbrus89 sbrus89 added the Omega label Mar 6, 2026
@sbrus89
Copy link
Collaborator

sbrus89 commented Mar 9, 2026

Merged locally, passes CTests on pm-gpu and pm-cpu for debug builds. I'm seeing a fail for the barotropic gyre case:

Polaris omega_pr suite

  • Baseline workdir: /pscratch/sd/s/sbrus/polaris_vert_advection_pr_suite_baseline/
  • Baseline build: /global/homes/s/sbrus/community/omega_worktrees/develop/components/omega
  • PR build: /global/homes/s/sbrus/scratch/polaris_vert_advection_pr_suite/build
  • PR workdir: /global/homes/s/sbrus/scratch/polaris_vert_advection_pr_suite
  • Machine: pm-cpu
  • Compiler: gnu
  • Build type: Release
  • Log: not found
  • Result:
    • Diffs (1 of 8):
      • ocean/planar/barotropic_gyre/munk/free-slip

Note, this case is not passing the property check (Nans are introduced somehow) so it is not just a baseline diff.

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

Labels

Projects

None yet

Development

Successfully merging this pull request may close these issues.

BottomDepth should not be read from the horizontal mesh

6 participants