Design docs for Omega analysis by xylar · Pull Request #706 · E3SM-Project/polaris

xylar · 2026-08-11T12:59:13Z

Design docs for ocean analysis in Polaris

Two design documents, no code. They are the design the six implementation branches below them are written against, and they are still expected to change as we discuss them.

ocean_analysis.md is the umbrella document: the long-term direction for analyzing Omega and MPAS-Ocean simulations in Polaris, why this is written from scratch rather than ported from MPAS-Analysis, where the line falls between what Omega computes in situ and what Polaris computes offline, the organizing principles the implementation follows, and a phased roadmap. It also carries the requirements that outlive Phase 1 — presentation and provenance, pruning and archiving, concurrency, and workflow integration — and the list of items Phase 1 knowingly defers.

ocean_analysis_initial.md designs the concrete near-term deliverable: the analysis capabilities the E3SM Ocean Team has committed to providing for Omega's initial coupled runs by September 15, 2026. Climatology maps, ocean heat content as maps and as a time series, global statistics time series, a MOC plot, and the publication layer that turns them into a browsable gallery, all produced by a Polaris suite pointed at a completed simulation.

Where review is most useful

The requirements, and the open questions collected at the end of the second document — those are where the decisions we still need to make as a team are gathered. Two in particular are worth an opinion: whether the gallery's presentation layer is the right shape given that the LCRC portal throttles, and what a product means for archiving as opposed to publication, which the design now raises but deliberately does not settle.

This is one of seven PRs, reviewed in order

The work is a chain in dependency order. GitHub cannot base a cross-fork PR on another fork branch, so every one of these is based on main and its diff contains everything below it. Reviewing them in order is what keeps each diff small — review only the commits this branch adds.

order	PR	branch	what it adds
1	#706	`add-analysis-designs`	the two design documents
2	#735	`add-omega-analysis-elevation`	the vertical reduction kernel
3	#726	`add-omega-analysis-scaffolding`	the suite, tasks and step structure
4	#736	`add-omega-analysis-publisher`	manifests, the staging tree, the gallery
5	#729	`add-omega-analysis-global-stats`	global statistics time series
6	#730	`add-omega-analysis-climatology-maps`	the ncclimo climatology and its maps
7	#737	`add-omega-analysis-heat-content`	ocean heat content, maps and time series

A gallery generated by the whole chain against the QU240 mock-up is published here, which is the quickest way to see what the suite produces: https://web.lcrc.anl.gov/public/e3sm/diagnostic_output/xasaydavis/omega_analysis_qu240_mockup_20260830/

Checklist

Documentation has been built locally and changes look as expected

cbegeman · 2026-08-24T17:08:18Z

+Whether this can be delivered in Omega by September 15 is a question for the
+team rather than something this design can settle.  **If it cannot**, Polaris
+shall fall back to computing mixed-layer depth offline from the monthly-mean
+conservative temperature and absolute salinity, using the same density
+threshold criterion, and shall label the resulting plots to make clear that
+they are computed from monthly means.  The fallback is described in the
+mixed-layer-depth algorithm design, along with what it does and does not
+capture.


@vanroekel Do you have thoughts about this?

Leaving this open for @vanroekel. The design carries the offline fallback with its caveats stated, and lists replacing it with Omega's in-situ diagnostic as the first follow-up if the fallback is what ships.

Sorry for the slow reply on this one. MLD will come in with the submesoscale PR which needs further testing. My guess is that doesn't get there by 9/15 but it should be discussed at the meeting tomorrow (unfortunately I'm out tomorrow). Your fallback plan seems quite reasonable to me.

we could try extract the MLD portion of the submeso PR for merging quicker potentially as well.

cbegeman · 2026-08-24T17:21:44Z

+Reconstruction is a linear operator on the edge field, so reconstructing from a
+climatology of normal velocity gives exactly the climatology of the
+reconstructed velocity.  There is therefore no accuracy argument for doing the
+reconstruction in the model, and doing it offline avoids adding a dependency on
+Omega work.  If Omega later writes reconstructed velocities directly, the step
+uses them in preference.


I agree that we should allow for offline reconstruction and use the online-reconstructed variables when present. I think we should plan to have zonal, meridional velocity be the default omega monthly output for memory reasons.

Agreed and flipped. The omega-monthly-means requirement now asks for reconstructed zonal and meridional velocity at cell centers rather than normalVelocity.

cbegeman · 2026-08-24T17:28:18Z

+so that a user can experiment with the TEOS-10 value without a code change.
+
+In the long run we do want $c_p^0$: because Omega carries conservative
+temperature, $\rho_0 c_p^0 \Theta$ is the definition of heat content rather


Suggested change

temperature, $\rho_0 c_p^0 \Theta$ is the definition of heat content rather

temperature, $\rho c_p^0 \Theta$ is the definition of heat content rather

Taken, though phrased per unit mass rather than per unit volume, since that is the form TEOS-10 actually defines: $c_p^0 \Theta$ is potential enthalpy per unit mass by definition, so the mass-weighted integral above is the heat content rather than an approximation to it. Same point, and it now reads consistently with the change in the previous thread.

xylar · 2026-08-25T11:46:26Z

@cbegeman, I revised the design quite a lot based on your feedback and the thinking that it invoked and also because of task-parallelism issues that arose (see #717). Please give it another look when you have time.

As before, feel free to comment as you have been doing. But also feel free to push changes if that's more efficient. This design can be as collaborative as you like.

cbegeman · 2026-08-25T14:53:30Z

If not already in ocean_analysis.md, can we plan to have the default behavior be a user supplying the full filepath to the omega.yml file for the simulation and polaris reading the format of output file names, replacing cfg options like monthly_mean_template?

Adds docs/design_docs/ocean_analysis.md, the umbrella document for analysis of Omega and MPAS-Ocean simulations in Polaris. It records the long-term direction: why the capability is written from scratch with Polaris and Omega in mind rather than ported from MPAS-Analysis, where the line falls between diagnostics Omega computes in situ and those Polaris computes offline, and a three-phase roadmap. Most requirements are deliberately placeholders to be filled in as each piece is designed. The exceptions are the conventions that the rest of the family depends on: MPAS-Ocean variable and dimension names as the Polaris standard, observations remapped onto the MPAS mesh rather than the model onto a comparison grid, and per-year decomposition of expensive work. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

Adds docs/design_docs/ocean_analysis_initial.md, designing the set of analysis capabilities the E3SM Ocean Team has committed to delivering for Omega's initial coupled runs by September 15, 2026. The deliverable is an omega_analysis suite, pointed at a completed simulation through a user-supplied config file, that produces map-view climatologies at configurable elevations, ocean heat content over elevation ranges as both maps and a global time series, time series from Omega's GlobalStats output, and a latitude-elevation plot of the global MOC, each with the netCDF behind it. The document also covers the Omega-side dependencies it rests on (monthly means, the geometric vertical coordinate, and a mixed-layer depth diagnostic with an offline fallback), and how repeated analysis over different date ranges reuses per-year work through range-keyed and year-keyed shared steps. Out of scope and stated as such: zppy integration, comparison with observations, regional analysis, and a global MOC time series. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

Addresses cbegeman's review comment asking whether analysis should use PseudoThickness directly rather than converting it to a thickness. Omega prognoses pseudo-height, z-tilde = -p / (rho0 g), so PseudoThickness is a normalized pressure increment in meters, and rho0 * h-tilde = rho * h = dp / g is the layer's mass per unit area, exactly, by hydrostatic balance. The geometric thickness is the derived quantity, h = rho0 * alpha * h-tilde, which is why it needs the equation of state and cannot be recovered offline. Record this in the conventions as a deliberate exception to the MPAS-Ocean-names rule: PseudoThickness is not translated to layerThickness, because the two mean the same thing for a mass-weighted integral and different things for a geometric one. Analysis asks for geometry via zMid and zInterface, and for mass via a get_layer_mass helper that reads each model's own mass-like thickness variable. This matches what mpaso_to_omega.yaml already does for the state variables. Also correct the vertical-geometry section, which was missing the factor of rho0 in Omega's accumulation of geometric elevation. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

Addresses cbegeman's suggestion to write the heat content integral with the full density and with pseudo-thickness weights. Heat content becomes a change of variable rather than an approximation: Q = cp0 * int(rho * Theta dz) = rho0 * cp0 * int(Theta dz-tilde) ~ rho0 * cp0 * sum_k Theta_k * w-tilde_k Because rho0 * h-tilde_k is the layer's mass per unit area exactly, a range covering whole layers carries no reference-density error at all. This removes the in-situ-versus-reference density error of the MPAS-Analysis formulation, which weighted by geometric thickness. The only quadrature error left is treating Theta as uniform within a layer. The geometric coordinate now enters only through the partial layers at a range boundary, as w-tilde_k = (w_k / h_k) * h-tilde_k. That fraction is exact with respect to Omega's own discretization, since Omega uses a single specific volume per layer, making z linear in z-tilde across the layer. Ranges stay geometric, and the doc now says why: "0 to 700 m" means 700 geometric meters in MPAS-Analysis and in the observational products this is compared against. Pseudo-depth ranges would be more natural for a mass-conserving model, and the trade-off is stated so the choice can be revisited. rho0 is no longer a config option. Unlike cp0 it is not a modeling choice we are making: it is the constant that defines pseudo-height, so it must be the model's, read from the same PCD constant that pstar.py uses to build the coordinate. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

Follows the algorithm change in the previous commit through to the implementation and testing sections. elevation_range_weights now takes the layer mass and returns mass per unit area within the range, computing the geometric overlap only to form the partial-layer fraction; heat_content drops its density argument. Passing interface differences as the layer mass recovers the purely geometric weights, which is what the elevation-slice utilities want. The one place that knows how each model spells its mass-like thickness stays get_layer_mass. The unit tests are rewritten in pseudo-height. Every synthetic column is built with h-tilde != h, so a formulation that confuses the two coordinates fails rather than coincidentally passing, and the whole-column result is checked to be invariant under perturbing the geometric interfaces while holding pseudo-thickness fixed. The previous tests leaned on a uniform-density column, in which the two coordinates agree and the distinction being tested disappears. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

This has now been worked out twice from scratch -- once for the tranche release, once for the publisher -- and both times the bump landed in the middle of a series and had to be moved. Record the tension rather than a single rule, because both halves matter: add the bump the moment a dependency changes, since a forgotten one silently leaves everyone on a stale environment, and rebase it to the end before merge, since a release marker under later work claims a version for a tree that was never released. Also note that the number is provisional while several branches are in flight, since each will expect the same next one and only the first to merge is right. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

The publication section said the publish step reads every fragment but never said what makes it run after the steps that write them. Nothing in the suite ordering guarantees it, so say that it declares each of those steps as a dependency through Step.add_dependency(), which exists for the case where what a step consumes is not known at setup. Those dependencies cross task boundaries, so record the two costs that come with that -- a pickle per dependency and a symlink per dependency -- and the one question they raise, which is what should happen when the gallery is asked for after only some products have been computed. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

Phase 1 publishes into a staging tree and stops, so getting results onto a machine's web portal is a manual copy and nothing makes them readable once they are there. Record it with the other deferred items, along with what is already available and what is not. mache carries base_path and base_url per machine and Polaris's merged config already exposes them, so reporting a URL is the same substitution MPAS-Analysis makes. What neither of those gives us is world-readable permissions, and publishing products by symlink into a private work tree is not something a web server can follow, so one of the two has to give. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

The previous version justified add_dependency() by saying the publish step cannot know what it consumes at setup. That is not true of the fragments: their filename is a constant and each step's path is fixed when the suite is built, so they can be declared as ordinary input files pointing into the steps that write them, which is what the convergence and baroclinic_channel viz steps already do. Say that, and say why walking the work directory for anything named manifest.json is the wrong shape: a step that finds its own inputs is invisible to the checking Polaris does at setup and to any scheduler that orders steps by what they read. Keep add_dependency() for what declaration cannot reach -- attributes only known after a step runs, and reporting a dependency that did not run -- and note that neither mechanism reorders anything, so publish still goes last. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

The design left open whether a fragment is mandatory, and the implementation chose optional on the grounds that it lets a developer publish after running only some products. It does not: add_dependency() makes each dependency's pickle an input of the publish step, so a step that did not run fails the publish step before a missing fragment is ever reached. The only case optionality covers is a step that ran and made nothing, and an empty product list covers that. So say that every step writes a fragment even when it produced no products, which is what lets the publish step declare each one as an input and have the framework check it. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

The snippet for the bounding-layer search counted valid midpoints at or above z and subtracted one, which is the index only when the column starts at level zero. z_mid is NaN above the valid range as well as below it, so the count is k1 - k_min + 1, and without the k_min term the index lands that many layers too high. Every column outside an ice-shelf cavity has k_min = 0, which is why this survived: it is invisible in the common case. In a cavity it inverts a clamp -- for k_min = 3, k_max = 5 and a z below every valid midpoint, the index clips to k_min and clamps to the top of the column where the rule above it says the answer is the bottom. While here, write down the two implementation details that follow from the rules rather than adding to them: clipping the weight into [0, 1] reproduces both clamping cases without a branch, and a column with one valid layer needs its own, since the clip range is empty there. Reported by the agent implementing the reduction, before writing any code. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

The word was doing two jobs. A manifest Product is a figure with the data behind it, and group, gallery and title are gallery vocabulary; the design also calls climatologies and reduced heat content "intermediate products", which are governed by reuse rather than by publication and never reach the manifest. Nothing said so plainly, so the answer had to be read out of Product.__init__. State it in the fragment section, and raise the question it exposes under pruning and archiving: what counts as a product depends on who is asking. An analysis that outlives its simulation cannot recompute the intermediates that requirement treats as discardable, so the archiving set is not the published set. The likely shape is a second collector reading the same fragments rather than a second list, since the two differ in what they select and in the vocabulary they need, not in how they learn what a step made. None of it is settled and Phase 1 builds none of it. What Phase 1 should avoid is foreclosing it, which is one more reason for the base class rather than the publisher's consumers to write a step's fragment. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

The work items were ordered by what could be on screen soonest, which put the vertical reduction last. That was right for delivery and wrong for review: it had the climatology maps document elevations as asked for and reported as skipped, and heat content introduce the vocabulary for elevation ranges while producing only the whole column, so a later branch had to retract prose and API that earlier branches had written. Mandatory fragments sharpened it, since a partial product set is now a visibly incomplete gallery. Put the vertical reduction second, right after the design. It imports nothing from polaris.Step, is unit tested against synthetic columns, and depends on nothing else in the list, so landing it first costs little and saves every branch above it from shipping a limitation it will have to un-say. Publication moves ahead of the products for the same reason it was already cut below them, now that every step that makes products must leave a fragment. Replace the branching rule with the one this implies. Independent items branching from a common parent was tried: siblings cannot see each other, so a helper on one arm was out of reach on the other; the integrated tree existed only as a throwaway merge branch whose documentation conflicts were re-resolved every time anything moved; and a rebase low in the tree forked it silently more than once. The earlier order was Carolyn's request, and what changed is the reason for it rather than anyone's mind about heat content drift mattering early. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

The list did not ask for a column with a single valid layer, where the clip range is empty, nor for an elevation below the second midpoint of a cavity column, which is the one case that tells apart a bounding-layer search that adds k_min back from one that does not. Both are now in the code and both would have failed against the formula the design carried until f6e0d79. Reported by the agent implementing the reduction, which verified that the targeted case and the linear-in-z test both fail when the term is reverted. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

Two sections disagreed. Reuse comes from seeded accumulators said the remaining months are spread over a process pool inside the step, in the present tense; One step per accumulator in Phase 1 says they are processed one after another and that Phase 1 runs serially, and the umbrella document defers process pools inside steps with the rest of the concurrency work. Keep the design intent -- when this is parallelized it belongs in a pool inside the step, not in more steps -- and say plainly that Phase 1 does not do it, in the same shape the replotting section already uses for the same distinction. Found by the agent implementing the accumulator, which implemented the serial reading rather than picking one and moving on. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

The design said every step leaves a fragment but not who writes it, and the one hint it gave -- "a product is declared once", pointing at add_produced_file() -- names a helper that lives on a different branch. That left the guarantee existing only in code. Record it: AnalysisStep.runtime_setup() writes the empty fragment and describing a product rewrites it, so no step author has an end-of-run call to forget. The step that would forget is the one that made nothing, which is the case the empty fragment exists for. The fragment is a declared output too, so a step that leaves none names itself instead of appearing as a missing input of publish. Also say why makes_products gates it without making the base class gallery-aware: an archiver would want a step to write a fragment without being a gallery step, and keeping the two questions separate is what leaves room for that. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

Seven sections were edited after midnight and still carried 2026/08/29, or in two cases a date from before this week. The session began on the 29th, which is where the wrong date came from, but the header should agree with the commit timestamps rather than with when someone started working. Updated: climatology-maps and its testing section, repeated-analysis, publication, order of work, pruning and archiving, and code organization. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

xylar · 2026-08-30T13:03:34Z

Status update, 2026-08-30

The whole analysis chain now runs end to end and publishes a browsable gallery. Three things changed today that affect every PR in the series.

The branches were relinearized. They were a tree with independent arms; they are now a line in dependency order, and each PR should be reviewed against the one below it:

order	PR	branch
1	#706	`add-analysis-designs`
2	#735	`add-omega-analysis-elevation`
3	#726	`add-omega-analysis-scaffolding`
4	#736	`add-omega-analysis-publisher`
5	#729	`add-omega-analysis-global-stats`
6	#730	`add-omega-analysis-climatology-maps`
7	#737	`add-omega-analysis-heat-content`

GitHub cannot base a cross-fork PR on another fork branch, so all seven are based on main and each diff contains everything below it. Review only the commits the branch itself adds.

Two things moved as part of this. The vertical reduction was item 11 and is now item 2, at the bottom: it is a dependency-light leaf module that everything vertically reduced needs, and landing it first means no branch above it ships a limitation that a later branch retracts — the climatology maps now support every elevation from the start, and heat content produces every elevation range rather than the whole column only. Publication moved ahead of the products for the same reason, since every step that makes products must leave a manifest fragment for the publish step to declare as an input.

Every product is now described in a manifest, so the gallery is populated rather than empty. Products are grouped by kind — climatology_maps and time_series — with galleries inside each: one per field for the maps, and one holding all of the global statistics, since a page per statistic would be a page with one plot on it. These groups are deliberately coarse for Phase 1; finer-grained grouping closer to MPAS-Analysis's belongs with the deferred work on how results are best displayed.

All seven branches were force-pushed after the relinearization, so local copies of any of them are stale.

What it produces

A gallery generated by the whole chain against the QU240 one-year mock-up, which is the quickest way to see what the suite does:

https://web.lcrc.anl.gov/public/e3sm/diagnostic_output/xasaydavis/omega_analysis_qu240_mockup_20260830/

Three steps publish nothing, and all three are waiting on Omega rather than on Polaris: velocityZonal/velocityMeridional, mixedLayerDepth, and the MOC. Each is reported as skipped in the step's log rather than failing the suite, so the products appear with no Polaris change once Omega writes them.

xylar · 2026-08-30T13:18:45Z

Testing

On add-analysis-designs at 05f680b48f, on chrysalis with intel/openmpi:

pytest tests --- 460 passed
pre-commit run --all-files --- clean, no issues found in 401 source files

The counts grow down the chain, since each branch's diff contains the branches below it: 460 at the design branch through 715 at the top.

End to end

The whole chain was run as the omega_analysis suite against a one-year QU240 Omega mock-up, using the derived monthly means rather than the raw history snapshots. All five tasks pass in about 1 minute 48 seconds, and the publish step reports:

published 69 products from 9 manifests, rendering 69 thumbnails
3 steps published nothing:
  analysis/climatology_maps/0001-0001/velocity
  analysis/climatology_maps/0001-0001/mixed_layer_depth
  analysis/moc/0001-0001

The resulting gallery is published here: https://web.lcrc.anl.gov/public/e3sm/diagnostic_output/xasaydavis/omega_analysis_qu240_mockup_20260830/

Two checks on the numbers rather than on the plumbing, both from the mock-up: each numeric elevation agrees with its nearest layer index to the offset between them --- temperature at -100 m is 0.17 K colder in the mean than k9 at -95.0 m, and at -2000 m only 0.02 K colder than k45, shrinking with depth as stratification weakens --- and the three heat content bands sum to the whole column to a relative 2e-14. The twelve monthly heat content values are equal to round-off, which is what an unforced mock-up should give and is a real check that the mass-weighted integral conserves.

Inheritance was exercised separately, since one simulated year cannot: a two-year simulation was faked with symlinks, the 1--2 range inherited year 1 from the completed 1--1 step bit-identically and computed only year 2, and re-running with a different specific heat capacity correctly rejected both the sibling seed and the step's own stale cache and recomputed all 24 months.

xylar · 2026-08-30T17:52:18Z

@cbegeman, do you want to give this another look? It is the first in the chain of PRs to get reviewed and merge. I know it's a ton of work so please let me know what I can to to make it as painless as possible.

@vanroekel, we have an question about MLD above that it would be good to get your take on, basically what you view as the best way forward by the 15th.

cbegeman · 2026-08-31T14:33:29Z

@xylar I gave it a pretty thorough look the first time around but there seems to have been a rebase that makes the "commits since last review" feature report all commits. Is there any change in particular I should take a look at because I'm inclined to approve on the basis of my previous review.

xylar · 2026-08-31T14:41:22Z

@cbegeman, there was a lot of rewriting that happened, in part based on your last round of commits. It was necessary to rebase to pull in a bunch of changes on main.

I'm fine with you approving based on the last round if you are okay with that. We can always revise the design as implementation details get finalized. The main thing would be if there are major choices that you want to revisit, but with 6 PRs waiting downstream, I think maybe that ship has sailed for this initial implementation at least.

xylar · 2026-08-31T14:41:48Z

And, sorry, I meant for this to be out of draft already.

cbegeman

@xylar Thanks for all your work on this. Apologies that I don't have the bandwidth to give this another thorough read. I am comfortable with going back and revising as we change the implementation. My general comment about these docs is that I think claude tends to over-explain which makes it difficult for me to maintain my attention for the whole document. It might be worth considering instructing claude to be concise and to assume prior knowledge of things that are written in Omega documentation, e.g., the bit about pseudo-thickness.

xylar · 2026-08-31T18:10:19Z

@cbegeman, yes sorry about that. Believe it or not, I do try to rein in Claude's tendency to over explain but I should push harder.

cbegeman · 2026-08-31T18:57:14Z

@xylar I believe it. And I should emphasize how much appreciate the awesome progress on this and the follow-on PRs. I'm very impressed

xylar self-assigned this Aug 11, 2026

xylar added documentation Improvements or additions to documentation ocean Related to the ocean component labels Aug 11, 2026

cbegeman reviewed Aug 11, 2026

View reviewed changes