Adjoints and checkpointing API

The adjoint surface covers cubed-sphere objectives, observations, footprints, preconditioning, and 4D-Var solvers. Tape owns checkpoint schedules and device, pinned-host, and memory-mapped tape storage.

Adjoints

AtmosTransport.Adjoints Module

julia

Adjoints

Prototype adjoint and footprint utilities.

The shipped path is a kernelized cubed-sphere reverse pass for tracer mass with fixed meteorology/air-mass evolution. It accumulates surface-emission footprints (and the initial-condition sensitivity) from an adjoint seed or observation forcing without perturbing each surface cell.

The forward model is the production step by construction: the adjoint replays the production StepPlan (CSStepSchedule: physics cadence, met windows, air-mass resets) with the production operators on a single-tracer state, records its tape through the recorder hooks of the production step (Operators.AbstractStepRecorder) and walks the tape backwards with each operator's adjoint. require_adjoint_support rejects a plan containing an operator without an adjoint before anything is recorded.

AtmosTransport.Adjoints.CELL_AVERAGE_PLACEMENT_CODE Constant

Stencil placement code of a cell average of column records (CSCellAverageObjective, CellAverage): every level of the cell, the inverse-variance-weighted mean of its members' operators.

AtmosTransport.Adjoints.COLUMN_MEAN_PLACEMENT_CODE Constant

Stencil placement code of an observation whose objective is a column mean (CSColumnMeanObjective): every level of the column, no vertical placement. The point objectives use the Output.VERTICAL_PLACEMENT_CODES.

AtmosTransport.Adjoints.FIXED_LEVEL_PLACEMENT_CODE Constant

Stencil placement code of an observation whose objective has a fixed stencil instead of being placed: a CSLayerMeanObjective naming its level, or another fixed operator.

AtmosTransport.Adjoints.RETRIEVAL_COLUMN_PLACEMENT_CODE Constant

Stencil placement code of an observation whose objective is an averaging-kernel column (CSRetrievalColumnObjective): every level of the column, weighted by the retrieval's pressure weights and averaging kernel on the model column stretched onto the retrieval column.

AtmosTransport.Adjoints.AbstractCSOptimType Type

julia

AbstractCSOptimType

Supertype for control-variable parameterizations used by the preconditioner. Concrete subtypes:

LinearOptimType — x = x_b + B^(1/2) χ. Unbounded.
LogNormalOptimType — x = x_b ⊙ exp(B^(1/2) χ). Strictly positive; requires x_b > 0.

AtmosTransport.Adjoints.AbstractCSOptimizer Type

julia

AbstractCSOptimizer

Supertype for CS 4D-Var optimization backends. Concrete subtypes implement cs_surface_flux_4dvar_solve(opt, cost_fn, controls) where cost_fn(controls) -> CS4DVarResult evaluates one cost-and-gradient pass through the 4D-Var driver.

Shipped backends:

CSGradientDescent — the dependency-free backtracking-line-search descent loop.

Additional backends:

CSLBFGS — Optim.jl L-BFGS wrapper.

AtmosTransport.Adjoints.AbstractCSSurfaceFluxCovariance Type

julia

AbstractCSSurfaceFluxCovariance{FT, A}

Supertype for CS surface-flux background-error covariance operators. Concrete subtypes implement apply_B_half! and apply_B_half_adjoint! on a single CSSurfaceFluxControl-shaped tuple NTuple{6, A} of (Nc, Nc) matrices.

AtmosTransport.Adjoints.AbstractCellAveraging Type

julia

AbstractCellAveraging

Whether bind_to_mesh averages column records that fall in one model cell: NoCellAveraging (the default) or CellAverage. The alternative that keeps every record and inflates its model error instead is inflate_cell_model_error.

AtmosTransport.Adjoints.AbstractFlaggedRecordPolicy Type

julia

AbstractFlaggedRecordPolicy

What bind_to_mesh does with records whose flag is not 0 (not to be assimilated): DropFlaggedRecords ("drop", the default) or CosampleFlaggedRecords ("cosample").

AtmosTransport.Adjoints.AbstractObservationKind Type

julia

AbstractObservationKind

What an observation record measures, the kind column of the observation file: PointObservationKind ("point") or ColumnObservationKind ("column"). bind_to_mesh turns a point record into a CSPlacedPointObjective and a column record into a CSColumnMeanObjective.

AtmosTransport.Adjoints.AbstractOutOfRangePolicy Type

julia

AbstractOutOfRangePolicy

Policy applied when an observation maps outside the valid model-step range. Use RejectOutOfRange, SkipOutOfRange, or ClampOutOfRange.

AtmosTransport.Adjoints.CS4DVarResult Type

julia

CS4DVarResult

One 4D-Var cost and gradient evaluation (cs_surface_flux_4dvar): the total, observation and background costs, the simulated values and residuals H(x) − y per observation, the gradient per control (and by control name), the controls and observations evaluated, and stencils, one CSObservationStencil per observation (the steps, weights and levels it was sampled at). An observation whose stencil is not assimilated (a co-sample, or an averaging-kernel column screened out by its surface-pressure mismatch) keeps its simulated value and residual but adds nothing to the cost and the gradient.

AtmosTransport.Adjoints.CSAdjointWorkspace Type

julia

CSAdjointWorkspace(mesh, prototype; column_scratch = false)

Scratch storage for CS adjoint sweeps. lambda_A is one halo-padded panel used as the per-panel transpose output. seam_seeds caches both neighbors' output-seed difference before any panel transpose overwrites those values. column_scratch is the per-column working storage of the FV3 vertical sweep's reverse pass (Nc × Nc × (Nz+1) × 5); it is empty unless built with column_scratch = true, which a tape with FV3-profile vertical sweeps needs (needs_column_scratch(scheme)).

AtmosTransport.Adjoints.CSCellAverageObjective Type

julia

CSCellAverageObjective(panel, i, j, members, weights, ids;
                       measurement_sigma = NaN, model_sigma = NaN, correlation = NaN)

The weighted mean $\sum_{n} w_{n} h_{n}$ of column objectives members of one physical CS interior cell (panel, i, j) (CSRetrievalColumnObjectives or CSColumnMeanObjectives; weights positive, summing to 1; ids the record ids of the members): a cell average of soundings, built by bind_to_mesh with cell_averaging = CellAverage(...), the weights proportional to $1 / σ_{n}^{2}$ . The operators are linear, so its functional is the same weighted sum of the members' functionals: the gradient of the averaged observation is the weighted combination of the members' gradients. Sampled like its members (cs_surface_flux_4dvar(...; column_met)). The surface-pressure screen applies to records, not averages, so an averaging-kernel member must have no surface-pressure limit (max_surface_pressure_mismatch_pa = Inf). measurement_sigma, model_sigma and correlation record how the average's error was built ( $σ_{meas}$ , $σ_{\mod}$ , $ρ$ ; departure files report them).

AtmosTransport.Adjoints.CSColumnMeanObjective Type

julia

CSColumnMeanObjective(panel, i, j)

Scalar objective equal to the final air-mass-weighted column mean mixing ratio at one physical CS interior cell.

AtmosTransport.Adjoints.CSDepartureRecord Type

julia

CSDepartureRecord(; id, tracer, instrument_type, date_components, lat, lon,
                    kind, time_window_s, clamped, panel, i, j,
                    vertical_placement, placement_height_agl_m, placement_pressure_pa,
                    steps, weights, levels,
                    observed_value, simulated_value, departure,
                    value_sigma, model_sigma, normalized_departure,
                    surface_pressure_mismatch_pa = NaN, retrieval_profile = Float64[],
                    profile_model_sigma = 0, assimilated = true, flag = 0, cell_average_id = 0,
                    cell_average_correlation = NaN)

One audit row of the CS departure file (schema v2): the observation record's provenance (id, tracer, instrument_type, date_components, Float64 lat/lon, kind, time_window_s; clamped when its time lies outside the run and ClampOutOfRange kept it), where it was sampled (the cell panel, i, j; the vertical_placement code — Output.VERTICAL_PLACEMENT_CODES for point records, COLUMN_MEAN_PLACEMENT_CODE for column records — with the height above the model surface or the pressure that placed it, NaN when unused; the sampled steps, their time-sampling weights and the levels read at each, k = 1 at the top and 0 for a column mean), and the forward-pass values. An averaging-kernel column row (RETRIEVAL_COLUMN_PLACEMENT_CODE) also carries the model minus retrieval surface pressure surface_pressure_mismatch_pa (Pa), the model profile on the retrieval layers or levels retrieval_profile (surface first, mol/mol) and its profile model error profile_model_sigma ( $\sqrt{\sum_{i} (h_{i} a_{i} σ_{i})^{2}}$ , 0 without a profile uncertainty); other rows have NaN, an empty profile and 0.

assimilated says whether the row itself is a term of the cost. It is false for a co-sampled record (its flag, the record's, is not 0: flagged_records = "cosample"), for a record averaged into a cell average (cell_average_id, the number of the average, 1, 2, … in observation order, is not 0; the average is in the cost, the record is reported on its own), and for an averaging-kernel column the surface-pressure screen left out; a row with a nonzero flag or a cell average is never assimilated.

A cell average itself (CELL_AVERAGE_PLACEMENT_CODE) is a row too, the term of the cost: its cell_average_id (shared with its member rows), an id above every record id of the observation set (the largest record id, at least 0, plus cell_average_id), the observed and simulated averages, value_sigma the measurement error of the average $σ_{meas}$ (with the correlation $ρ$ between its records, cell_average_correlation; NaN on other rows), model_sigma its model error, assimilated, flag 0; tracer, instrument, date, coordinates and window are those of its earliest record. Over the assimilated rows, $\sum {normalized_departure}^{2} = 2 J_{obs}$ . The constructor rejects non-finite values and checks the sign convention departure = simulated_value − observed_value and normalized_departure = departure / σ, $σ = \sqrt{{value_sigma}^{2} + {model_sigma}^{2} + {profile_model_sigma}^{2}}$ .

AtmosTransport.Adjoints.CSDepartureSet Type

julia

CSDepartureSet(records; mesh_Nc, mesh_panel_convention, mesh_cs_definition_tag,
               t_start, dt_seconds, nsteps, run_id = nothing, iteration = nothing)

Collection of CSDepartureRecord rows together with the run-level metadata that lets a downstream consumer reconstruct the forward run that produced the file. mesh_panel_convention and mesh_cs_definition_tag pin the cubed-sphere geometry; t_start, dt_seconds, and nsteps pin the time grid. Optional run_id and iteration are caller-supplied attribution tags.

AtmosTransport.Adjoints.CSFootprintResult Type

julia

CSFootprintResult

Reverse-mode result for one scalar objective. footprints[t] is an NTuple{6} of (Nc, Nc) arrays containing dJ / dE, where E is the per-cell model-storage emission rate the surface-flux operator applies at the palindrome center of model step t. Storage is dry mixing ratio × carrier-air mass; these rates are not physical kg-species per second and are not divided by cell area. The derivative includes the timestep factor. lag_steps[t] == nsteps - t. A step plan without a surface-flux operator (surface_flux = false) has no emission entry point and an empty footprints.

initial_adjoint is dJ/drm at the initial state: an NTuple{6} of halo-padded panels like the initial tracer mass (the halo cells are zero: halos are a function of the interiors).

AtmosTransport.Adjoints.CSGradientDescent Type

julia

CSGradientDescent(; iterations, initial_step, min_step,
                    step_shrink, gradient_tolerance, line_search,
                    log)

Plain backtracking-line-search gradient descent. iterations is the maximum count of accepted descent steps; the loop terminates early when the gradient L2 norm falls below gradient_tolerance or the line search shrinks step below min_step. With line_search = false every candidate step is accepted (constant step size unless the loop exits on tolerance).

log = true captures per-iteration diagnostics into CS4DVarSolveResult.log (see CSIterationLog): cost decomposition (observation vs background), gradient L2 norm, accepted step size, and wall-clock elapsed seconds since the solve started. Default false matches the pre-C3 behavior.

AtmosTransport.Adjoints.CSIterationLog Type

julia

CSIterationLog(entries)

Per-iteration diagnostic log for a CS 4D-Var solve. Populated by optimizers when their log::Bool field is set to true; otherwise the CS4DVarSolveResult.log field is nothing.

The log carries enough information for after-the-fact convergence diagnostics: a monotone cost trajectory check, gradient-norm decay curve, observation-vs-background cost decomposition, and a wall-clock time series.

AtmosTransport.Adjoints.CSIterationLogEntry Type

julia

CSIterationLogEntry(iteration, cost, observation_cost, background_cost,
                    gradient_norm, step_size, elapsed_seconds)

One row of per-iteration solver diagnostics. Captures the cost-decomposition (observation vs background terms), the optimizer's documented gradient norm (L2 for CSGradientDescent, L-infinity for CSLBFGS), the line-search step accepted on this iteration (0 for optimizers that don't expose a step — e.g. L-BFGS), and the wall-clock time since the solve started.

The convention is that an iteration is logged AFTER the descent direction has been accepted; the iteration = 0 row records the initial probe before any descent step.

AtmosTransport.Adjoints.CSLBFGS Type

julia

CSLBFGS(; iterations, gradient_tolerance, m, show_trace, log)

Limited-memory BFGS (L-BFGS) optimizer backed by Optim.jl. m is the history length; the line search is Optim's LBFGS default (Hager–Zhang). iterations is the max outer iteration count; gradient_tolerance is Optim's L-infinity stopping threshold on the χ-space / x-space gradient. The same L-infinity norm is recorded in gradient_norm_history and CSIterationLogEntry.gradient_norm for this backend. show_trace = true forwards verbose progress to Optim's logger. log = true records one CSIterationLogEntry for the initial state (iteration 0) and one per accepted iterate.

Used by the polymorphic dispatch surface cs_surface_flux_4dvar_solve. The cost closure cost_fn(controls) -> CS4DVarResult is handed to Optim as one objective with value, gradient and combined value-and-gradient methods, flattening each control's NTuple{6, Matrix} to a Vector{FT} and back. Each cost_fn call is a forward run plus the adjoint, so all three methods share a cache of the last evaluation: cost_fn runs once per point Optim visits, and the initial point and the returned result reuse those evaluations. (The result costs one more evaluation only when Optim returns an earlier iterate, after a failed line search or a cost increase.)

AtmosTransport.Adjoints.CSLayerMeanObjective Type

julia

CSLayerMeanObjective(panel, i, j, level)

Scalar objective equal to the final-layer mixing ratio rm[panel][i, j, level] / m[panel][i, j, level] on the physical CS interior indices. level = Nz is the surface layer.

AtmosTransport.Adjoints.CSObservation Type

julia

CSObservation(step, objective, value, sigma; id = 0, assimilate = true)
CSObservation(sampling::CSTimeSampling, objective, value, sigma; id = 0, assimilate = true)

Scalar observation of the CS 4D-Var surface-flux path: objective sampled as sampling says (one step, or a weighted window of steps; see CSTimeSampling), with observed value and error standard deviation sigma. step is the last sampled step, after which the observation is complete (the only one for an instantaneous observation). id is the identifier of the observation record it was bound from (bind_to_mesh), 0 for an observation built by hand.

assimilate = false makes it a co-sample: cs_surface_flux_4dvar simulates it on the forward replay and reports its simulated value, residual and stencil (with assimilated = false), but it adds nothing to the cost or the gradient, and does not extend the taped run. bind_to_mesh co-samples records with a nonzero flag under flagged_records = "cosample", and the members of cell averages (CellAverage).

AtmosTransport.Adjoints.CSObservationColumnMet Type

julia

CSObservationColumnMet(; p_top, gravity = STANDARD_GRAVITY,
                         temperature = 280.0, humidity = nothing, mass_basis = DryBasis())

The column meteorology with which CSPlacedPointObjective observations choose their model layer at each sampled state (cs_surface_flux_4dvar(...; column_met)). With the air mass $m_{k}$ (kg per cell, k = 1 at the top) of the sampled state and the cell area $A$ :

interface pressures $p_{1 / 2} = p_{top}$ , $p_{k + 1 / 2} = p_{k - 1 / 2} + g m_{k} / A$ , from p_top (Pa, the model-top interface, the pure-pressure A_ifc[1] of the met model's hybrid coordinate) and gravity (m s⁻²);
hypsometric interface heights above the model surface, $z_{k - 1 / 2} = z_{k + 1 / 2} + (R_{d} T_{v, k} / g) \ln (p_{k + 1 / 2} / p_{k - 1 / 2})$ , with the layer temperature $T_{k}$ : a number (K, every layer; the forward sampler's layer_height_temperature_kelvin, 280 K by default), an Output.AbstractLayerTemperature, per-layer panels (an NTuple{6} of arrays laid out like the air mass, with or without the halo, K, for every step), or a vector of such panels with one entry per model step (the state after that step). With the pressure basis TotalFromHumidity (the default) on a dry basis the pressures are the total ones below and $T_{v, k} = T_{k} (1 + 0.61 q_{k})$ the virtual temperature; with DryInterfaces the dry interfaces above and $T_{v, k} = T_{k}$ ;
for TotalFromHumidity on a dry basis (the heights of every placed point observation, pressure placement, averaging-kernel columns), the total interface pressures $p_{k + 1 / 2} = p_{k - 1 / 2} + g m_{k} / (A (1 - q_{k}))$ with the specific humidity $q_{k}$ (kg kg⁻¹ of moist air; nothing, static panels, or one set of panels per step, as temperature).

mass_basis is the basis of the air mass (DryBasis(), the default of transport binaries, or MoistBasis()).

AtmosTransport.Adjoints.CSObservationRecord Type

julia

CSObservationRecord(id, date_components, lat, lon, kind,
                    elevation_m, intake_height_m, altitude_m, pressure_pa,
                    time_window_s, flag, value, value_sigma, model_sigma,
                    instrument_type, tracer[, retrieval])
CSObservationRecord(; id, date_components, lat, lon, value, value_sigma,
                      instrument_type, tracer, kind = PointObservationKind(),
                      elevation_m = NaN, intake_height_m = NaN, altitude_m = NaN,
                      pressure_pa = NaN, time_window_s = 0, flag = 0, model_sigma = 0,
                      retrieval = nothing)

One observation row of the schema-v2 CSObservationSet.

date_components: (year, month, day, hour, minute, second) UTC, Int16, the centre of the sampling window.
lat, lon: degrees, Float64 (a point within Float32 rounding of a cell edge must bind to the cell that contains it).
kind: PointObservationKind or ColumnObservationKind (the keyword form also takes "point" / "column").
Vertical fields, NaN when the record does not give them: elevation_m (ground, m above sea level), intake_height_m (m above the ground), altitude_m (m above sea level), pressure_pa (total pressure, Pa, positive). Column records give none of intake_height_m, altitude_m, pressure_pa.
time_window_s: length (s, ≥ 0) of the sampling window centred on the record time; 0 samples the model state of one step.
flag: quality flag, 0 = assimilate; bind_to_mesh drops the others.
value, value_sigma: observation and its 1-σ error, in the tracer's native unit (dry mole fraction for trace gases); model_sigma (≥ 0) is a model/representation error added in quadrature (σ² = value_sigma² + model_sigma²).
instrument_type (e.g. "ICOS", "TCCON"), tracer (e.g. "CO2").
retrieval: the retrieval block of an averaging-kernel column record (an Output.RetrievalColumn: pressures from the surface up, pressure weights, column averaging kernel, prior profile in mol/mol, optional model-profile uncertainty, and its mapping: layer means between interfaces, Output.LayerRebin, or values on levels, Output.HatLevelProjection), or nothing. Only column records carry one; such a record is an absolute column average, value in (0, 1] mol/mol with value_sigma and model_sigma not above it. Its model-error term $\sum_{i} (h_{i} a_{i} σ_{i})^{2}$ (Output.retrieval_model_error_variance, zero without a model-profile uncertainty) is added to σ².

AtmosTransport.Adjoints.CSObservationSet Type

julia

CSObservationSet(records, time_origin)

Collection of CSObservationRecord rows, with unique ids, that share an absolute time_origin string (an ISO-8601 timestamp, e.g. "1900-01-01 00:00:00"). The on-disk format always carries the origin even though the per-record dates are stored as absolute components — the origin tags the convention so a downstream consumer that wants to convert to seconds-since-origin can do so without parsing the per- record dates.

AtmosTransport.Adjoints.CSObservationStencil Type

julia

CSObservationStencil(panel, i, j, steps, weights, levels, placement_method,
                     height_agl_m, pressure_pa,
                     surface_pressure_mismatch_pa = NaN, retrieval_profile = Float64[],
                     assimilated = true)

Where an observation was sampled in a 4D-Var evaluation, returned per observation in CS4DVarResult.stencils: the cell (panel, i, j), the sampled steps and their weights (its CSTimeSampling), the model level sampled at each step (levels, k = 1 at the top, resolved from the replayed column at the end of that step; 0 for an operator that reads more than one level, as a column mean), the placement_method (Output.VERTICAL_PLACEMENT_CODES for a CSPlacedPointObjective, COLUMN_MEAN_PLACEMENT_CODE, RETRIEVAL_COLUMN_PLACEMENT_CODE, CELL_AVERAGE_PLACEMENT_CODE or FIXED_LEVEL_PLACEMENT_CODE otherwise), and the height above the model surface (m) or the pressure (Pa) that placed it (NaN when not used). An operator without a single cell has the cell (0, 0, 0).

For an averaging-kernel column (CSRetrievalColumnObjective): surface_pressure_mismatch_pa, the model minus the retrieval surface pressure (Pa; time-weighted over the sampled steps; NaN for other operators); retrieval_profile, the model's dry mixing ratios on the retrieval layers or levels (surface first, the projected profile $P x$ , time-weighted; empty for other operators).

assimilated is false for a co-sample (an observation built with assimilate = false: a flagged record under flagged_records = "cosample", or a member of a cell average) and for an averaging-kernel column whose surface-pressure mismatch exceeded max_surface_pressure_mismatch_pa under SkipUnplaced(): the observation's simulated value is reported, but it is not in the cost.

AtmosTransport.Adjoints.CSPlacedPointObjective Type

julia

CSPlacedPointObjective(panel, i, j, area, placement, pressure_basis)

Mixing ratio rm / m of the model layer that contains an in-situ sample, in the physical CS interior cell (panel, i, j) of area area (m², the cell that contains the sample). placement (an Output.VerticalPlacement) says how the sample is placed in the column — its height above the model surface or its pressure — and is resolved once, by bind_to_mesh with Output.resolve_placement, the rule of the forward observation sampler. The layer itself is chosen every time the objective is sampled, from the replayed air mass of that state and the CSObservationColumnMet (hypsometric layer heights and interface pressures, both on pressure_basis: total pressures and the virtual temperature from the humidity for TotalFromHumidity, the dry interfaces and the layer temperature for DryInterfaces), with Output.sample_column_layer!. Sampling it therefore needs the column meteorology: cs_surface_flux_4dvar(...; column_met); the final-time footprint and Jacobian entry points do not take it.

AtmosTransport.Adjoints.CSRetrievalColumnObjective Type

julia

CSRetrievalColumnObjective(panel, i, j, area, retrieval, pressure_basis;
                           max_surface_pressure_mismatch_pa = Inf,
                           mismatch_policy = ErrorOnUnplaced())

Averaging-kernel column average of the physical CS interior cell (panel, i, j) (area area, m²): a retrieved column (satellite or TCCON sounding) compared with the model through the retrieval's own layers, pressure weights, averaging kernel and prior (retrieval, an Output.RetrievalColumn). At every sampled state the model's interface pressures of pressure_basis (total pressures from the dry air mass and the specific humidity, TotalFromHumidity, or the dry interfaces, DryInterfaces) are stretched onto the retrieval column and the dry mixing ratios $x_{k} = r_{k} / max (m_{k}, ε)$ are mapped onto its layers or levels by the projection $P$ of the retrieval's mapping (the rebin $R$ of TM5-4DVar for a layer column, the hat projection for a level column; Output.retrieval_column_weights!):

X = \sum_{k} H_{k} x_{k} + \sum_{i} h_{i} (1 - a_{i}) x_{a, i}, H_{k} = \sum_{i} h_{i} a_{i} P_{i k} .

max_surface_pressure_mismatch_pa bounds the time-weighted difference between the model and the retrieval surface pressure; a record above it is an error under mismatch_policy = ErrorOnUnplaced() and is left out of the cost (reported, not assimilated) under SkipUnplaced(). Sampling it needs the column meteorology (cs_surface_flux_4dvar(...; column_met)); built by bind_to_mesh from a column record with a retrieval block.

AtmosTransport.Adjoints.CSSeedObjective Type

julia

CSSeedObjective()

Marker objective used when the caller supplies an explicit final adjoint seed (dJ/drm_final) instead of one of the built-in scalar objectives.

AtmosTransport.Adjoints.CSStepSchedule Type

julia

CSStepSchedule(; physics_cadence = nothing, binary_window_contract = false,
               air_mass_reset = :none, steps_per_window = 1, window_dt,
               start_time = 0, window_start_air_mass = nothing,
               window_end_air_mass = nothing)
CSStepSchedule(plan::StepPlan; steps_per_window = 1, window_dt, start_time = 0,
               binary_substeps = false, chemistry = NoChemistry(),
               window_start_air_mass = nothing,
               window_end_air_mass = nothing)

The production step a cubed-sphere adjoint differentiates: the step plan and the met windows it runs on. The keyword form builds the plan as DrivenSimulation does (StepPlan(; physics_cadence, window_contract = binary_window_contract, air_mass_reset)); DrivenSimulation's own schedule is CSStepSchedule(sim). That constructor validates the simulation's operators and retains its chemistry operator; replay never replaces it with a no-op. Unsupported chemistry is rejected by the adjoint setup guard.

steps_per_window — advection steps per met window: an integer (every window) or one entry per window (the run may not go beyond them). Step s of window w advances window_dt / steps_per_window[w], computed as DrivenSimulation computes Δt.
start_time — the run clock [s] at the first step, as DrivenSimulation's start_time; the replay's current_time follows production's clock.
binary_window_contract / binary_substeps — the transport binary's per-window contract: one palindrome subcycle per step (the binary's schedule), and with physics_cadence = :window (the default) convection and chemistry once per met window over window_dt.
air_mass_reset — :none, :preserve_tracer_mass or :preserve_vmr: the air mass is reset to window_start_air_mass[w] on entering window w > 1 and, with the window contract, to window_end_air_mass[w] before the window-end physics (vectors of NTuple{6} halo-padded panels, one per window; a reset refreshes the halos of its target in place, as production does).

The default (CSStepSchedule(window_dt = dt)) is one step per window and the full operator suite at every step, TransportModel.step!'s composition.

AtmosTransport.Adjoints.CSSurfaceFluxControl Type

julia

CSSurfaceFluxControl(window, value; background=nothing, sigma=nothing)

Surface-flux control block for the prototype 4D-Var path. value is an NTuple{6} of (Nc, Nc) rate arrays tied to window. Optional background and sigma add the standard diagonal background term 0.5 * ((value - background) / sigma)^2; sigma can be a scalar or matching panel arrays.

AtmosTransport.Adjoints.CSSurfaceFluxPreconditioner Type

julia

CSSurfaceFluxPreconditioner(covariance, background, optim_type)

Bundle of (covariance, background, optim_type) that performs the change of variables between preconditioned χ-space and physical x-space for one CS surface-flux control. background is the prior x_b shaped as NTuple{6, Matrix{FT}} (matching the covariance's sigma shape and a single CSSurfaceFluxControl.value).

Constructor validates background shapes against the covariance and — for LogNormalOptimType — requires every background entry to be strictly positive (else log(x ./ x_b) is undefined and the log-normal transform is meaningless). It also pre-allocates a panel-tuple scratch buffer reused by apply_preconditioner_inverse! and the LogNormal variant of apply_preconditioner_adjoint!. The preconditioner is therefore not thread-safe — build one instance per worker thread if you intend to call these methods concurrently.

AtmosTransport.Adjoints.CSSurfaceFluxWindow Type

julia

CSSurfaceFluxWindow(name, steps; weights=nothing, normalize=false)

Named surface-flux control window for Jacobian aggregation. steps can be a single step or any iterable of step indices. By default, a window sums per-step surface-emission footprints, which corresponds to perturbing the same rate over the whole window. Pass normalize=true for an average-rate control, or explicit weights for custom temporal basis functions.

AtmosTransport.Adjoints.CSTapeByteEstimate Type

julia

CSTapeByteEstimate

Counts and byte estimate for the CS adjoint tape. Each *_records field is an op count — the number of records of that type the forward pass will push onto the tape — except state_records which is a payload-staging count (full panel tuples written, doubled for nonlinear schemes that stage both panels_m and panels_rm per sweep).

state_bytes is the raw panel-data cost of the tape. It equals state_records * bytes_per_state for split-sweep tapes; Lin-Rood additionally retains face-shaped tuples whose exact byte cost is included in state_bytes. Halo, midpoint and end-of-step records and the diffusion-palindrome's two op records contribute scalar metadata only and are not counted in state_bytes. step_records counts the end-of-step markers (one per step, after convection), where the reverse pass injects observation forcing.

total_records is the op count: sweep + halo + midpoint + diffusion + convection + step. It is not in general equal to state_records + halo + midpoint, because (a) nonlinear schemes have state_records = 2 * sweep_records and (b) the diffusion palindrome contributes two op records per step but stages only one panel tuple.

AtmosTransport.Adjoints.CSTimeSampling Type

julia

CSTimeSampling(steps, weights)
CSTimeSampling(step)

The model states an observation samples: the states after steps (strictly increasing, after convection: the state runtime output sees), averaged with weights (positive, summing to 1). The simulated value is $H = \sum_{s} w_{s} h (r_{s})$ , with $h$ the objective evaluated on the state after step $s$ . CSTimeSampling(step) samples one state (weight 1).

bind_to_mesh builds an instantaneous sampling for a record without a time window (the step whose interval contains the record time) and, for a record with time_window_s > 0, the steps whose intervals overlap the window with weights proportional to the overlap (TM5-4DVAR's Δt/ΣΔt window mean).

AtmosTransport.Adjoints.CellAverage Type

julia

CellAverage(; time_tolerance_s = 60.0, correlation = 0.0)

Average the assimilated column records (flag 0) of one tracer that fall in the same model cell, sample the same model steps with the same weights, and lie within time_tolerance_s of the earliest record of their group (records of one overpass: an OCO-2 orbit crosses a C90 cell in about 15 s, and the next orbit passes 99 minutes later). A group of $n \geq 2$ records $y_{i}$ with measurement errors $σ_{i}$ (value_sigma) becomes one observation with weights $w_{i} = 1 / σ_{i}^{2}$ :

\bar{y} = \frac{\sum_{i} w_{i} y_{i}}{\sum_{i} w_{i}}, \bar{H} = \frac{\sum_{i} w_{i} H_{i}}{\sum_{i} w_{i}}, σ_{meas}^{2} = \frac{\sum_{i} \sum_{j} w_{i} w_{j} ρ_{i j} σ_{i} σ_{j}}{(\sum_{i} w_{i})^{2}}, {\bar{σ}}^{2} = σ_{meas}^{2} + \frac{\sum_{i} w_{i} σ_{\mod, i}^{2}}{\sum_{i} w_{i}},

with $ρ_{i i} = 1$ and $ρ_{i j} =$ correlation between distinct records (0 by default: separate 10-s spans; their within-span correlation is already in $σ_{i}$ ), and $σ_{\mod, i}^{2}$ the record's model error (model_sigma² plus the profile term of an averaging-kernel column), added once, not reduced by the averaging. $\bar{H}$ is a CSCellAverageObjective; the operators are linear, so its gradient is the $w$ -weighted combination of the members'. The members stay bound as co-samples (assimilate = false), so every record keeps its own simulated value. Needs max_surface_pressure_mismatch_pa = Inf (the screen applies to records, not averages).

AtmosTransport.Adjoints.ClampOutOfRange Type

Keep observations outside the valid model-step range: a record whose time window reaches into the run samples the part inside it (truncated and renormalised), any other record the nearest step (first or last).

AtmosTransport.Adjoints.ColumnObservationKind Type

kind = "column": a column-averaged dry mole fraction (TCCON, a satellite sounding). A column record with a retrieval block (its layers, pressure weights, averaging kernel and prior, an Output.RetrievalColumn) is bound to an averaging-kernel column (CSRetrievalColumnObjective); one without it to the air-mass-weighted column mean of its cell (CSColumnMeanObjective). A column record states no intake height, altitude or pressure.

AtmosTransport.Adjoints.CosampleFlaggedRecords Type

flagged_records = "cosample": records with a nonzero flag are bound as co-samples (CSObservation(...; assimilate = false)): simulated on the forward replay and reported with their flag, never in the cost or the gradient. The OCO-2 MIP needs co-samples of every record, the withheld and evaluation-only ones included.

AtmosTransport.Adjoints.DiagonalCSCovariance Type

julia

DiagonalCSCovariance(sigma::NTuple{6, AbstractMatrix})

Pure-diagonal CS surface-flux covariance B = diag(σ)². sigma[p][i, j] is the background standard deviation at panel p, interior cell (i, j). All six panels must share the same (Nc, Nc) shape and every entry must be finite and strictly positive.

This is the cleanest baseline; B^(1/2) is self-adjoint, so apply_B_half_adjoint! reduces to apply_B_half!.

AtmosTransport.Adjoints.DropFlaggedRecords Type

flagged_records = "drop" (default): records with a nonzero flag are not bound.

AtmosTransport.Adjoints.IsotropicGaussianCSCovariance Type

julia

IsotropicGaussianCSCovariance(sigma::NTuple{6, AbstractMatrix},
                              correlation_length_cells::Real)

Panel-local CS surface-flux covariance B = D · C · D with isotropic separable Gaussian correlation of length correlation_length_cells (in interior-cell units).

Each panel's (Nc, Nc) grid is treated as periodic and the correlation matrix C is the circulant whose first row is the wrapped Gaussian. L = C^(1/2) (the symmetric square root) is applied via FFT: in spectral space L is multiplication by sqrt(C̃) where C̃ is the real, non-negative spectral eigenvalues of C. The factorization is B^(1/2) = D · L, so (B^(1/2))^T = L^T · D = L · D (since D and L are both symmetric).

The wrapped-Gaussian normalization pins C[i, i] = 1, hence B[i, i] = σ_i² — the diagonal of B recovers the per-cell variance, independent of correlation length.

To keep B^(-1/2) finite at operational Float32 resolutions, each one-dimensional factor of the separable spectrum is floored at its maximum times sqrt(eps(FT)), then renormalized to preserve the unit diagonal. The corresponding minimum two-dimensional covariance eigenvalue is therefore approximately its maximum times eps(FT). This is a precision-aware regularization of unresolved high-frequency modes.

Limitations:

Cross-panel correlation is dropped. Each panel is smoothed in isolation; the implicit wrap-around at panel boundaries leaves edge artefacts.
The FFT path is CPU-only.
The struct carries mutable fft_buf / fft_scratch scratch buffers reused by every apply_B_half! / _adjoint! / _inverse! call. As a consequence the covariance is not thread-safe — concurrent calls on the same instance race on the scratch. Build one instance per worker thread if needed.

AtmosTransport.Adjoints.LinearOptimType Type

julia

LinearOptimType

Linear additive parameterization: x = x_b + B^(1/2) χ. The physical control x is unbounded; the only nonlinearity is in the cost via the observation operator H(x).

AtmosTransport.Adjoints.LogNormalOptimType Type

julia

LogNormalOptimType

Log-normal multiplicative parameterization: x = x_b ⊙ exp(B^(1/2) χ). The physical control is guaranteed strictly positive for any finite χ, at the cost of the background x_b having to be strictly positive everywhere. Suitable for emissions, mole fractions, and other non-negative atmospheric quantities.

AtmosTransport.Adjoints.NoCellAveraging Type

cell_averaging = NoCellAveraging() (default): every record is its own observation.

AtmosTransport.Adjoints.PointObservationKind Type

kind = "point": an in-situ sample at one place in the column (surface station, tower, flask, aircraft, AirCore). Placed in the model column by its height above the ground or its pressure (intake_height_m, altitude_m with elevation_m or a DEM, pressure_pa); a point record that gives none of the vertical coordinates in use cannot be placed (unplaced_records).

AtmosTransport.Adjoints.RejectOutOfRange Type

Reject the first observation outside the valid model-step range.

AtmosTransport.Adjoints.SkipOutOfRange Type

Drop observations outside the valid model-step range.

AtmosTransport.Adjoints.adjoint_support Method

julia

adjoint_support(x) -> Union{Nothing, String}

nothing if the cubed-sphere adjoint records and transposes x (a step-plan stage, an advection scheme or a diffusion, surface-flux, convection or chemistry operator), otherwise the reason it does not. See require_adjoint_support.

AtmosTransport.Adjoints.apply_B_half! Function

julia

apply_B_half!(y::NTuple{6, AbstractMatrix},
              cov::AbstractCSSurfaceFluxCovariance,
              chi::NTuple{6, AbstractMatrix}) -> y

Apply the (non-symmetric) square root B^(1/2) to chi and write the result into y. Concretely, y = D · L · chi, where D is the diagonal scaling and L is the symmetric correlation square root. y and chi may alias; the implementation writes through a scratch buffer when needed.

AtmosTransport.Adjoints.apply_B_half_adjoint! Function

julia

apply_B_half_adjoint!(g_chi::NTuple{6, AbstractMatrix},
                      cov::AbstractCSSurfaceFluxCovariance,
                      g_phys::NTuple{6, AbstractMatrix}) -> g_chi

Apply the adjoint (B^(1/2))^T to g_phys. With B^(1/2) = D · L and D, L both symmetric, this is g_chi = L · D · g_phys.

AtmosTransport.Adjoints.apply_B_half_inverse! Function

julia

apply_B_half_inverse!(y::NTuple{6, AbstractMatrix},
                      cov::AbstractCSSurfaceFluxCovariance,
                      x::NTuple{6, AbstractMatrix}) -> y

Apply B^(-1/2) to x and write the result into y. With the factorization B^(1/2) = D · L (D diagonal, L symmetric), B^(-1/2) = L^(-1) · D^(-1). Used by the preconditioner to invert the change of variables x = x_b + B^(1/2) χ (Linear) or x = x_b ⊙ exp(B^(1/2) χ) (LogNormal).

The Gaussian inverse is regularized in the high-frequency tail by the precision-scaled spectral floor documented on IsotropicGaussianCSCovariance. It is finite for arbitrary finite inputs, but poorly resolved high-frequency modes remain sensitive to roundoff.

AtmosTransport.Adjoints.apply_preconditioner! Method

julia

apply_preconditioner!(x, prec, chi) -> x

Forward change of variables. With Linear optim_type, x = x_b + B^(1/2) χ. With LogNormal, x = x_b ⊙ exp(B^(1/2) χ). x, chi are NTuple{6, AbstractMatrix}.

AtmosTransport.Adjoints.apply_preconditioner_adjoint! Method

julia

apply_preconditioner_adjoint!(g_chi, prec, x, g_phys) -> g_chi

Adjoint of the tangent at base x. Linear: g_χ = (B^(1/2))^T g_phys. LogNormal: g_χ = (B^(1/2))^T (x ⊙ g_phys).

This is the observation part of the χ-space gradient under the chain rule of T. The caller adds χ from the 0.5 ‖χ‖² term of the preconditioned cost separately.

AtmosTransport.Adjoints.apply_preconditioner_inverse! Method

julia

apply_preconditioner_inverse!(chi, prec, x) -> chi

Inverse change of variables. Linear: χ = B^(-1/2) (x - x_b). LogNormal: χ = B^(-1/2) log(x ⊙ inv(x_b)). LogNormal requires x > 0 everywhere (else log is undefined); the function throws ArgumentError on a non-positive entry.

For a non-trivial Gaussian covariance, B^(-1/2) is numerically unstable for inputs with significant energy in the high-frequency tail of the spectrum where the correlation transfer is tiny. Round-trip T^(-1)(T(χ)) is exact (the FFT factors cancel bit-exactly up to roundoff), but applying the inverse to arbitrary user data may amplify noise. Use the bijection on clean round-trip inputs.

AtmosTransport.Adjoints.apply_preconditioner_tangent! Method

julia

apply_preconditioner_tangent!(delta_x, prec, x, delta_chi) -> delta_x

Tangent-linear of the forward map at the base point x = T(χ). Linear: δx = B^(1/2) δχ. LogNormal: δx = x ⊙ (B^(1/2) δχ). x is ignored for Linear but is part of the uniform API.

AtmosTransport.Adjoints.bind_to_mesh Method

julia

bind_to_mesh(set::CSObservationSet, mesh::CubedSphereMesh, t_start::Dates.DateTime, dt::Real;
             nsteps = nothing, tracer_filter = nothing, out_of_range_policy = RejectOutOfRange(),
             placement = VerticalPlacementSpec(), terrain = ObservationTerrain(),
             max_surface_pressure_mismatch_pa = Inf, flagged_records = DropFlaggedRecords(),
             cell_averaging = NoCellAveraging())
    -> Vector{CSObservation}

Bind each record of set to a CSObservation of the CS 4D-Var path: its time sampling, the cell that contains it, and its objective.

Time. t_start is the model start (UTC) and dt the step length (s): step k covers [t_start + (k−1)dt, t_start + k dt) and is represented by the state at its end, after convection. A record with time_window_s = 0 samples the step whose interval contains its time. A record with time_window_s = W > 0 samples the window [t − W/2, t + W/2] around its time t: every step whose interval overlaps the window, with weight $Δ t_{k} / \sum Δ t$ where $Δ t_{k}$ is the overlap of step k's interval with the window (a step only partly inside counts its overlap; steps fully inside count dt, TM5-4DVAR's Δt/ΣΔt mean). The window is truncated at the run's ends (0 and nsteps·dt) and the weights are renormalised over the part inside the run. See CSTimeSampling.
Out of range. A record whose time falls outside 1:nsteps (only step >= 1 is enforced when nsteps === nothing) is handled by out_of_range_policy: RejectOutOfRange (default, "reject") throws, SkipOutOfRange ("skip") drops it, ClampOutOfRange ("clamp") samples the part of its window inside the run (truncated and renormalised) when the window reaches into the run, else step 1 or nsteps. Skipped and clamped records are counted in an @info. Records with a time window need nsteps.
Ids. Record ids must be unique (departures are joined to records by id).
Cell. (lat, lon) (Float64) → the cell containing the point (Output.locate, the forward observation sampler's locator).
Point records (PointObservationKind) become CSPlacedPointObjectives, placed by the forward sampler's rule (Output.resolve_placement) with placement (an Output.VerticalPlacementSpec: the ordered vertical_coordinates, height above ground and/or pressure; the pressure basis; the altitude option; the unplaced policy) and the static ground fields of terrain (an Output.ObservationTerrain: model surface height, point DEM). A point record that gives none of the coordinates in use is unplaced: an ArgumentError with unplaced_records = "error" (default), dropped and counted with SkipUnplaced().
Column records (ColumnObservationKind) with a retrieval block become CSRetrievalColumnObjectives: the averaging-kernel column on the interface pressures of placement.pressure_basis (total pressures from the humidity by default), with the surface-pressure limit max_surface_pressure_mismatch_pa (Pa, positive; Inf, the default, lets the stretch absorb any mismatch, as TM5) and the unplaced policy of placement as its mismatch policy: a record whose model surface pressure differs from the retrieval's by more than the limit is an error under unplaced_records = "error", and is left out of the cost (and counted) under "skip". The model surface pressure is known only at the sampled states, so the check runs in the 4D-Var's forward replay; it is the surface pressure of the pressure basis, so with DryInterfaces the mismatch also contains the weight of the water vapour (about 0.5 % of the surface pressure, more in the tropics). Column records without a block become CSColumnMeanObjectives.
Records with a nonzero flag are not assimilated: dropped with flagged_records = DropFlaggedRecords() ("drop", the default), bound as co-samples (CSObservation(...; assimilate = false): simulated and reported, never in the cost or the gradient) with CosampleFlaggedRecords ("cosample"). Records of another tracer are dropped when tracer_filter names one.
Cell averaging. With cell_averaging = CellAverage(...), assimilated column records of one cell, step sampling and overpass become one inverse-variance-weighted observation (CellAverage, CSCellAverageObjective), appended after the records, whose own observations become co-samples.
The observation sigma is $\sqrt{σ_{value}^{2} + σ_{model}^{2} + σ_{prof}^{2}}$ (exactly value_sigma without a model error), where $σ_{prof}^{2} = \sum_{i} (h_{i} a_{i} σ_{i})^{2}$ is the model-error term of an averaging-kernel column with a profile uncertainty $σ_{i}$ (Output.retrieval_model_error_variance; fixed by the record, so the cost stays quadratic); the observation carries the record id.

set.time_origin is documentation only.

AtmosTransport.Adjoints.build_departure_set Method

julia

build_departure_set(set::CSObservationSet, result::CS4DVarResult, mesh, t_start, dt, nsteps;
                    run_id = nothing, iteration = nothing) -> CSDepartureSet
build_departure_set(set, observations, simulated, stencils, mesh, t_start, dt, nsteps; ...)

Build a CSDepartureSet from a 4D-Var evaluation: its observations (bound from set by bind_to_mesh), simulated values and stencils (CSObservationStencil), joined to the records of set by record id, one row per observation in observation order. Co-sampled records (flagged_records = "cosample") are rows with their flag and assimilated = false. A cell average (CSCellAverageObjective) is a row (placement code CELL_AVERAGE_PLACEMENT_CODE, its measurement and model sigma and correlation, cell_average_id its number 1, 2, … and a row id above every record id of set), and so are its records, as co-samples with the average's cell_average_id; a hand-built average is checked against its records as bind_to_mesh builds it.

Validation (ArgumentError on any violation):

equal lengths; unique record ids in set; every observation's id names a record of set (hand-built observations have id 0 and no record);
every simulated value is finite;
each observation matches its record as bind_to_mesh binds it: the kind (point → CSPlacedPointObjective, column with a retrieval block → CSRetrievalColumnObjective with that block, column without one → CSColumnMeanObjective), the containing cell of the record's Float64 (lat, lon), the time sampling re-derived from the record date and window on the run grid t_start, dt, nsteps (for a record outside the run, as ClampOutOfRange binds it; such rows are marked clamped), the value, and $σ = \sqrt{{value_sigma}^{2} + {model_sigma}^{2} + {profile_model_sigma}^{2}}$ ;
each stencil matches its observation's cell and time sampling.

Averaging-kernel column rows take the surface-pressure mismatch, the retrieval-level model profile and assimilated from their stencil.

t_start is an ISO-8601 string (Dates.format(t, "yyyy-mm-dd HH:MM:SS")).

AtmosTransport.Adjoints.convert_oco2_mip_observations Method

julia

convert_oco2_mip_observations(path, output; partition = :day, time_range = nothing,
                              data_types = ("LN", "LG", "OG"), assimilate_flags = (1,),
                              model_sigma = :file, mapping = Output.HatLevelProjection(),
                              instrument = nothing, tracer = "CO2") -> Vector{String}

Convert an OCO-2 or OCO-3 v11 MIP 10-second file to observation files of schema v2 (write_observations), one per day (partition = :day, output holds {YYYYMMDD}) or per month (:month, {YYYYMM}), for the records of time_range (all records without one); days or months without records get no file. The records are those of read_oco2_mip_observations with the same options; the files record the source file, the options and the mapping as root attributes. Returns the paths written, in time order.

AtmosTransport.Adjoints.cs_surface_emission_footprint Method

julia

cs_surface_emission_footprint(..., objective; kwargs...) -> CSFootprintResult

Generate reverse-mode footprints for a scalar final-time objective with respect to surface-emission rates at each prior model step.

The tape is recorded from the production step plan of schedule (a CSStepSchedule; default: one step of dt per met window, the full operator suite every step), with a surface-flux operator in the plan (zero rates unless base_emission_rates is given), and walked backwards with each operator's adjoint. require_adjoint_support rejects a plan with an operator the tape cannot transpose before anything is recorded.

Supported CS schemes are the split-sweep UpwindScheme(), SlopesScheme(NoLimiter()), PPMScheme(NoLimiter()) and monotone PPMScheme(), and LinRoodPPMScheme (ORD 5 or 7) with upwind or FV3 kord = 8 vertical transport. The limited PPM paths, including the FV3 vertical sweeps, store tracer branch states from the base trajectory; pass base_emission_rates when differentiating around nonzero surface emissions. Optional ImplicitVerticalDiffusion support transposes the Backward-Euler column solve in kernels on CPU/GPU and uses the runtime's surface-flux coupling at the palindrome center (V(dt/2) -> emissions -> V(dt/2), or emissions -> V(dt) with DiffusiveSurfaceFluxBoundary). Optional CMFMCConvection support transposes the well-mixed sub-cloud, updraft, and tendency passes; optional TM5Convection and CMFMCMatrixConvection support replays the same column matrix and applies the transposed LU solve, at the plan's physics cadence (every step, or once per met window over the window).

The replay uses the caller's convection workspace with the persistent state of a fresh run (its ledger carry and CFL cache are reset when a replay starts) and the diffusion workspaces' factors as scratch; do not pass a live production model's workspaces. cs_adjoint_inputs(sim) captures every input of a DrivenSimulation configuration, with fresh workspaces.

tape_storage selects the per-tape-slot storage policy (:device, :pinned_host, or :mmap). When tape_storage = :mmap, the optional tape_path kwarg points the on-disk tape at a user-owned directory instead of a temporary one — required for tapes that must persist past the call (manual inspection, partial-run debug, multi-session campaigns). tape_path is created on demand and preserved past finalize_tape!; it is rejected for non-:mmap storage. Under StrideCheckpoint each window writes a window_NNNNN/ subdirectory; under RevolveCheckpoint each base-case step writes a step_NNNNN/ subdirectory. The LinRood scheme is pinned to :device and rejects tape_path unconditionally.

AtmosTransport.Adjoints.cs_surface_emission_footprint_from_seed Method

julia

cs_surface_emission_footprint_from_seed(final_adjoint_rm, panels_m0,
                                        panels_am_steps, panels_bm_steps,
                                        panels_cm_steps, mesh; kwargs...)

General reverse-pass entry point. final_adjoint_rm is an NTuple{6} of halo-padded adjoint tracer-mass arrays containing dJ/drm_final for any scalar objective or observation operator. The reverse pass and surface-gradient accumulation use the same CPU/GPU kernels as cs_surface_emission_footprint.

With surface_flux = true (the default) the step plan carries a surface-flux operator and the result holds the emission footprints. surface_flux = false replays a configuration without surface sources (palindrome center V(dt)): the result's footprints are empty and initial_adjoint is the initial-condition sensitivity dJ/drm_0, which every result also carries.

AtmosTransport.Adjoints.cs_surface_flux_4dvar Method

julia

cs_surface_flux_4dvar(..., observations, controls; scheme, kwargs...) -> CS4DVarResult

Evaluate the prototype CS surface-flux 4D-Var cost and gradient. Controls are named CSSurfaceFluxControls over CSSurfaceFluxWindows. Observations are scalar CSObservations sampled after the model steps of their CSTimeSampling (one step, or a weighted time window; after convection, the state runtime output sees). CSPlacedPointObjective observations choose their layer, and CSRetrievalColumnObjective observations (averaging-kernel columns) build their weights, at each sampled state from column_met, a CSObservationColumnMet (required when there are any). An averaging-kernel column whose surface-pressure mismatch exceeds its limit is an error, or, under SkipUnplaced, is reported with its simulated value and residual but left out of the cost and the gradient (its stencil is not assimilated). Co-samples (CSObservation(...; assimilate = false): flagged records bound with flagged_records = "cosample", the members of cell averages) are simulated and reported the same way, never in the cost or the gradient, and do not extend the taped run: adding them leaves the cost and the reverse pass's forcing unchanged, bit for bit.

One cost-and-gradient evaluation is two forward runs and one reverse pass (with FullCheckpoint; the other schedules add their recomputation), whatever the number of observations, as in TM5-4DVAR:

a forward run gives the simulated values Hᵢ(x) = Σₛ ωᵢₛ hᵢ(rₛ), each observation's stencil (CSObservationStencil: steps, weights, levels) and its linear functional (weights from the replayed air mass at each sampled step, times ωᵢₛ);
each observation contributes (Hᵢ(x) − yᵢ)/σᵢ² times its functional's weights to the adjoint forcing of each sampled step;
one taped forward run to the last sampled step and one reverse pass, which adds each step's forcing at that step's end, give ∂J_obs/∂E for every step; each control's window collects its steps.

Both runs replay the production step plan of schedule (a CSStepSchedule; default: one step of dt per met window, the full operator suite every step), with the controls' rates in its surface-flux operator. Steps after the last sampled step do not affect the cost and get a zero gradient. checkpoint (FullCheckpoint(), StrideCheckpoint(k), RevolveCheckpoint()) and tape_storage / tape_path select how the one taped run is stored, as in cs_surface_emission_footprint.

scheme is required: inversion callers must name the advection operator whose forward and adjoint they intend to use; there is no hidden numerical default.

Two modes:

Unconditioned (default, preconditioner = nothing). Each control's .value is the physical-space control directly. Optional diagonal background terms come from .background + .sigma and are added per control via _background_cost_and_gradient!.
Preconditioned (preconditioner non-nothing). Each control's .value is the preconditioned-space variable χ. The function: (1) computes x_k = T(χ_k) per control via the preconditioner; (2) runs the forward simulation with x; (3) reverses through the existing observation-gradient path to get ∂J_obs/∂x_k per control; (4) applies T'(χ_k)^T to get ∂J_obs/∂χ_k and adds χ_k (the derivative of the 0.5 ‖χ‖² regularization). The reported cost is 0.5 ‖χ‖² + observation_cost, the reported gradient is ∂J/∂χ, and the reported controls are the original χ-space inputs. Per-control .background and .sigma are ignored in this mode — the background term comes from 0.5 ‖χ‖², and the background x_b is carried by the preconditioner itself.

preconditioner accepts either a single CSSurfaceFluxPreconditioner (broadcast to every control — the same preconditioner is used for each of them) or a Vector{<:CSSurfaceFluxPreconditioner} aligned 1-to-1 with the control vector. Any other length mismatch throws ArgumentError.

AtmosTransport.Adjoints.cs_surface_flux_4dvar_optimize Method

julia

cs_surface_flux_4dvar_optimize(..., observations, controls;
                                optimizer = nothing, kwargs...)
    -> CS4DVarSolveResult

Run an optimization loop around cs_surface_flux_4dvar.

optimizer::AbstractCSOptimizer (kwarg) selects the backend. When omitted, a CSGradientDescent is constructed from the legacy descent-policy keyword arguments (iterations, initial_step, min_step, step_shrink, gradient_tolerance, line_search, log) so existing call sites keep working unchanged.

log::Bool = false is consumed by the entrypoint, not forwarded to cs_surface_flux_4dvar. It only takes effect on the default CSGradientDescent path; an explicit optimizer = CSLBFGS(..., log = true) (or optimizer = CSGradientDescent(..., log = true)) carries its own setting and overrides the entrypoint kwarg.

Remaining keyword arguments are forwarded to cs_surface_flux_4dvar on every cost evaluation — including preconditioner = ... for the preconditioned-cost path.

AtmosTransport.Adjoints.cs_surface_flux_4dvar_solve Function

julia

cs_surface_flux_4dvar_solve(optimizer::AbstractCSOptimizer,
                             cost_fn,
                             controls) -> CS4DVarSolveResult

Run optimizer against the cost closure cost_fn(controls) -> CS4DVarResult starting from controls. Dispatches on the optimizer's concrete type — additional backends are added by defining a new method here, not by branching inside an existing one.

AtmosTransport.Adjoints.cs_surface_flux_jacobian Method

julia

cs_surface_flux_jacobian(..., objectives, windows; kwargs...)

Compute surface-flux Jacobian maps for several layer/column objectives and named time windows. Each returned footprints[obj, window] entry is an NTuple{6} of (Nc, Nc) arrays. Window aggregation is a weighted sum of per-step emission-rate footprints; use CSSurfaceFluxWindow(...; normalize=true) for average-rate controls or explicit weights for a custom temporal basis. kwargs are those of cs_surface_emission_footprint. Every checkpoint schedule records (and recomputes) the forward run once and walks it once per objective.

AtmosTransport.Adjoints.cs_tape_byte_estimate Method

julia

cs_tape_byte_estimate(panels_m0, panels_am_steps, panels_bm_steps,
                      panels_cm_steps, mesh, scheme;
                      cfl_limit = DEFAULT_CFL_LIMIT,
                      flux_scale = one(eltype(panels_m0[1])),
                      diffusion_op = NoDiffusion(),
                      convection_op = NoConvection(),
                      schedule = nothing, surface_flux = true, dt = nothing)
    -> CSTapeByteEstimate

Size the tape that cs_surface_emission_footprint records for the given problem: per-record-class counts and the total state_bytes, the raw panel-data cost of the tape, independent of storage policy. The transport records (sweeps, halos, Lin-Rood horizontal records, palindrome centers and air-mass resets) are counted by replaying the transport with a counting recorder; convection and diffusion do not change the air mass, so their records are counted from the step plan.

The same state_bytes figure applies to all three storage policies because none of them compress the payload:

tape_storage = :device — state_bytes is the device-resident RAM cost; the reverse loop holds the full tape on the backend.
tape_storage = :pinned_host — state_bytes is the pinned host RAM cost; at most two staged panel sets (one per read role: air mass and tracer mass) are mirrored on the device at a time.
tape_storage = :mmap — state_bytes is the on-disk footprint of records.bin. On a GPU the device read caches hold one NTuple{6, T} per read role and distinct shape signature on top of that.

bytes_per_state is the size of one full-panel state tuple. Lin-Rood's state_bytes also includes the face-panel tuples of its horizontal records, whose dimensions differ from a full state. Halo, midpoint, reset and end-of-step records carry only scalar metadata and are counted in total_records but NOT in state_bytes.

AtmosTransport.Adjoints.inflate_cell_model_error Method

julia

inflate_cell_model_error(set::CSObservationSet, mesh, t_start, dt; nsteps = nothing,
                         time_tolerance_s = 60.0, tracer_filter = nothing,
                         out_of_range_policy = RejectOutOfRange()) -> CSObservationSet

The alternative to CellAverage that keeps every record: the assimilated column records grouped as CellAverage groups them (same tracer, model cell and step sampling as bind_to_mesh would bind them, within time_tolerance_s of the group's earliest record) get their model-error variance multiplied by the number $k$ of records in their group, model_sigma → √k · model_sigma, i.e. model_sigma² times $k$ ; the profile-error term $\sum_{i} (h_{i} a_{i} σ_{i})^{2}$ of an averaging-kernel column (from its retrieval block) is not inflated. The model error of records that sample the same model column is correlated, so $k$ records should not reduce it $k$ -fold; inflating each by $k$ gives the group the weight of one record's model error. Bind the returned set with the same run grid; departures built from it report the inflated model_sigma. Every other record is returned unchanged.

AtmosTransport.Adjoints.read_departures Method

julia

read_departures(path::AbstractString) -> CSDepartureSet

Parse a schema-v2 departures NetCDF written by write_departures and return an in-memory CSDepartureSet.

Validation (fails fast with ArgumentError):

Root attribute cs_departures_schema must equal "v2" (schema-v1 files, with Float32 coordinates and a single step, are rejected).
Every required root attribute must be present: mesh_Nc, mesh_panel_convention, mesh_cs_definition_tag, t_start, dt_seconds, nsteps, departure_sign_convention (equal to "simulated_minus_observed"), and n_clamped (the number of rows marked clamped).
Every required variable must be present, date_components dimensioned (date_component = 6, obs), every obs variable of the obs length, and the sample_counts must add up to the sample length.
The averaging-kernel column block (surface_pressure_mismatch_pa, retrieval_model_profile (retrieval_level, obs), profile_model_sigma, with assimilated) is all or nothing; without it every row has the defaults (NaN, no profile, 0). assimilated, flag and cell_average_id are optional otherwise (defaults: assimilated, flag 0, no cell average). The profile must state ordering = "surface_to_top" and a mole-fraction unit, the mismatch a pressure unit (converted to mol/mol and Pa).
Every row must pass the CSDepartureRecord checks.

AtmosTransport.Adjoints.read_observations Method

julia

read_observations(path::AbstractString) -> CSObservationSet

Parse a schema-v2 observation NetCDF written by write_observations and return an in-memory CSObservationSet.

Validation (fails fast with ArgumentError):

Root attribute cs_observations_schema must equal "v2". Schema-v1 files are rejected: their Float32 coordinates can bind a point near a cell edge to the neighbouring cell, and their single alt cannot place a record by its height above the ground. Rewrite them as v2 with write_observations.
Root attribute time_origin must be present.
Every v2 variable must be present, lat and lon stored as Float64, and date_components dimensioned (date_component = 6, obs).
elevation_m, intake_height_m, altitude_m and pressure_pa read NaN or a declared fill value (_FillValue, e.g. xarray's NaN) as absent; any other variable with a fill value is an error.
The optional retrieval block (see write_observations) is all or nothing: retrieval_pressure_pa, retrieval_pressure_weight, retrieval_averaging_kernel, retrieval_prior_profile and retrieval_surface_pressure_pa (retrieval_profile_sigma, retrieval_mapping and retrieval_prior_column optional), with retrieval_interface = retrieval_level + 1. retrieval_mapping names each block's mapping ("layer_rebin": the pressures are the N + 1 interfaces of N layers; "hat_levels": the N levels themselves); a file without it holds layer records only. Every retrieval-level variable states ordering = "surface_to_top" (index 1 is the retrieval surface); pressures state a pressure unit (converted to Pa), the prior and its uncertainty a mole-fraction unit ("mol mol-1", "1", "ppm", "ppb", …; converted to mol/mol), weights and kernel are dimensionless. A record's values come first and NaN pads the rest (the number of layers or levels is that of its pressure weights); a record whose block entries are all NaN has none. Each block must pass the Output.RetrievalColumn checks (pressures strictly decreasing from the surface, weights ≥ 0 summing to 1, finite values, prior in [0, 1] mol/mol, surface pressure equal to the first pressure), and only column records may carry one.
value, value_sigma and model_sigma without a units attribute are in the tracer's native unit, read as they are. A dry mole-fraction unit ("mol mol-1", "1", "ppm", "ppb", …) is converted to mol/mol; it must then be stated on all three, and a file with averaging-kernel columns must state one on all three.
Every row must pass the CSObservationRecord checks, and record ids must be unique.

AtmosTransport.Adjoints.read_oco2_mip_observations Method

julia

read_oco2_mip_observations(path; time_range = nothing, data_types = ("LN", "LG", "OG"),
                           assimilate_flags = (1,), model_sigma = :file,
                           mapping = Output.HatLevelProjection(), instrument = nothing,
                           tracer = "CO2") -> CSObservationSet

The records of an OCO-2 or OCO-3 v11 MIP 10-second file (read-only) as an observation set of schema v2 (CSObservationRecord), one per record whose time (date, truncated to the second) lies in time_range = (start, stop) (DateTimes, half-open; nothing = all):

record field	from the MIP file
`id`	`sounding_id`
`date_components`	`date[1:6]` (UTC)
`lat`, `lon`	`latitude`, `longitude`
`kind`	`"column"`
`elevation_m`	`altitude` (m above sea level)
`time_window_s`	0 (the model step containing the time)
`value`, `value_sigma`	`xco2_2019_scale`, `xco2_uncertainty` (ppm → mol/mol)
`model_sigma`	`model_error` (`model_sigma = :file`), or a constant in ppm, or 0
`flag`	0 when `assimilate_flag ∈ assimilate_flags` and the data type is in `data_types`; else `1 + assimilate_flag`
`instrument_type`	`"<instrument>_<data type>"`, e.g. `"OCO2_OG"`
`tracer`	`tracer`
`retrieval`	the 20 levels (`sigma_levels × psurf`, hPa → Pa), `pressure_weight`, `xco2_averaging_kernel`, `co2_profile_apriori` (ppm → mol/mol), surface first, and the prior column `xco2_apriori` (ppm → mol/mol; `Output.read_oco2_retrieval_columns`)

data_types name the experiment's data types (LN land nadir, LG land glint, LT land target, LTR land transition, ON ocean nadir, OG ocean glint, OT, OTR, MIX mixed land/water, SAM OCO-3 snapshot area maps; the MIP assimilates LN, LG and OG). assimilate_flags are the MIP assimilate_flag values to assimilate (1; TM5-4DVar's MIP setup also takes the withheld 2). Every record of the range is kept: flag 1 marks assimilate_flag = 0 (sparse, target, transition, mixed scenes; evaluation only), 2 an assimilable record whose data type is outside data_types, 3 a withheld record (assimilate_flag = 2), unless assimilate_flags selects them. Bind with flagged_records = "cosample" to co-sample them.

mapping is Output.HatLevelProjection() (the levels as given, the default) or Output.LayerRebin() (TM5-4DVar's 19 layers, Output.tm5_layer_column; the CLI's --mapping tm5). instrument defaults to the first token of the file's ShortName ("OCO2", "OCO3"). Sounding ids are unique within one file; OCO-2 and OCO-3 ids can coincide, so give each instrument its own set.

AtmosTransport.Adjoints.require_adjoint_support Method

julia

require_adjoint_support(plan::StepPlan, ops::StepOperators...) -> plan

Check that the cubed-sphere adjoint records and transposes every stage of plan and every operator those stages apply (one ops per distinct operator set, e.g. per step when the diffusion operator changes). Throws an ArgumentError naming every unsupported piece. Every adjoint entry point calls it before recording; DrivenSimulation and TransportModel setups can call it directly (require_adjoint_support(sim)).

AtmosTransport.Adjoints.require_adjoint_support Method

julia

require_adjoint_support(scheme) -> scheme

Return the advection scheme if the cubed-sphere reverse pass has an adjoint for it, otherwise throw an ArgumentError that names the missing piece and lists the supported choices. Dispatches on the union of supported schemes that the footprint and 4D-Var entry points also accept, so a configuration is rejected with an explanation before any tape is recorded.

Supported: UpwindScheme(), SlopesScheme(NoLimiter()), PPMScheme with NoLimiter, MonotoneLimiter or CW84Limiter and the default vertical reconstruction, and LinRoodPPMScheme (ORD 5 or 7) with any of its vertical sweeps: upwind or the FV3 kord = 8 profile (vertical = "fv3_kord8" or "fv3_kord8_signed"). Unsupported: limited SlopesScheme, the split-sweep PPMScheme with the FV3 vertical profile, NoAdvection, and every other scheme or limiter.

AtmosTransport.Adjoints.run_cs_footprint_forward Method

julia

run_cs_footprint_forward(..., objective; kwargs...) -> scalar

Run the cubed-sphere production step plan forward and return objective at the final time. Optional emission_rates[t][panel][i, j] entries are surface emission rates in model-storage units per second (dry mixing ratio × carrier air mass per second, integrated over the cell) that the surface-flux operator applies at the palindrome center of step t. Convert physical species-mass inventories and per-area rates before calling this low-level API.

Each step is the production step of schedule (a CSStepSchedule; default: one step of dt per met window and the full operator suite every step, as TransportModel.step!): the transport block — step-start halo refresh, palindrome and its center — and, at the plan's cadence, the convection column operator (CMFMCConvection(), TM5Convection() or CMFMCMatrixConvection) and the air-mass resets. With surface_flux = true (the default) the plan carries a surface-flux operator (zero rates when emission_rates = nothing) and the palindrome center follows the diffusion operator's surface_flux_coupling: V(dt/2) -> emissions -> V(dt/2) (SplitSurfaceFluxCoupling) or emissions -> V(dt) (DiffusiveSurfaceFluxBoundary); surface_flux = false replays a configuration without surface sources, whose center is one solve V(dt). Diffusion requires a panel-native DiffusionWorkspace with filled layer_thickness. flux_scale multiplies all transport fluxes and is shared with the tape/reverse APIs so forward finite differences exercise the same model trajectory.

AtmosTransport.Adjoints.write_departures Method

julia

write_departures(path::AbstractString, set::CSDepartureSet) -> String

Write set to path as a schema-v2 NetCDF departures file (see schemas/cs_departures_v2.toml) and return path; an existing file is overwritten. The sampled steps form a contiguous ragged array on the sample dimension: observation n owns sample_count[n] consecutive entries of sample_step, sample_weight and sample_level. When any row is an averaging-kernel column, the file adds surface_pressure_mismatch_pa, retrieval_model_profile (on a retrieval_level dimension, surface first, NaN-padded) and profile_model_sigma for every row. assimilated is written when any row is an averaging-kernel column or not assimilated, flag when any row has a nonzero flag (co-sampled records), and cell_average_id when any row belongs to a cell average.

AtmosTransport.Adjoints.write_observations Method

julia

write_observations(path::AbstractString, set::CSObservationSet; attributes = Dict()) -> String

Write set to path as a schema-v2 NetCDF observations file (schemas/cs_observations_v2.toml) and return path; an existing file is overwritten. attributes are extra root attributes (provenance; not read back), which may not replace cs_observations_schema or time_origin.

Dimensions: obs (unlimited), date_component (= 6).
Variables: id (i64), date_components (i16), lat, lon (f64), kind (string "point" / "column"), elevation_m, intake_height_m, altitude_m, pressure_pa (f64, NaN = absent), time_window_s (f64), flag (i8), value, value_sigma, model_sigma (f64), instrument_type, tracer (string).
Root attributes: cs_observations_schema = "v2", time_origin = set.time_origin.
When any record has a retrieval block (an averaging-kernel column): the dimensions retrieval_level (the most layers or levels of any record) and retrieval_interface (one more), retrieval_pressure_pa (Pa: the interfaces of a layer record, the levels of a level record), retrieval_pressure_weight, retrieval_averaging_kernel, retrieval_prior_profile (mol mol⁻¹) and, when any record gives one, retrieval_profile_sigma (mol mol⁻¹), each (retrieval_level or retrieval_interface, obs), surface first (ordering = "surface_to_top") and NaN-padded; retrieval_surface_pressure_pa (obs); and retrieval_mapping (obs), the record's mapping ("layer_rebin" or "hat_levels", empty without a block) and, when any record gives one, retrieval_prior_column (obs) (mol mol⁻¹, NaN = absent). Records without a block are NaN throughout. Such a file states units = "mol mol-1" on value, value_sigma and model_sigma: its values are dry mole fractions in mol/mol, the unit of the retrieval priors (every record of the file).

Tape and checkpointing

AtmosTransport.Tape Module

julia

Tape

Tape storage policies (device, pinned-host, mmap on disk), tape record types, and checkpoint schedules for the AtmosTransport.jl cubed-sphere adjoint pipeline.

Tape imports nothing from other AtmosTransport modules. src/AtmosTransport.jl loads it before Adjoints, which imports these names with using ..Tape: .... The forward recorder, reverse loop and checkpoint drivers (src/Footprint/) and the 4D-Var code (src/Inversion/) are files included into the Adjoints module, not separate modules.

Module dependency order: Tape → Adjoints (includes Footprint/ and Inversion/ files)

AtmosTransport.Tape.AbstractHostCSTapeStorage Type

julia

AbstractHostCSTapeStorage <: AbstractCSTapeStorage

Tape storage whose staged slots live off the source device: host memory for PinnedHostCSTapeStorage, a file for MmapCSTapeStorage. Reading a slot copies it into a device-side read cache and returns that cache. There is one cache per (read role, panel shapes) pair, so the two operands that one reverse kernel reads together (the air mass and tracer mass of a limited-PPM sweep record) are always distinct buffers. Subtypes carry the fields read_template (a zero-size array on the source backend, or nothing) and read_caches.

AtmosTransport.Tape.AbstractTapeReadRole Type

julia

AbstractTapeReadRole

The reverse-kernel operand a tape read feeds, passed as _tape_panels(slot, role). Host-side storage keeps one read cache per role, and a read stays valid until the next read with the same role. A reverse kernel that takes two staged states therefore reads them with two different roles. Reading both through one shared cache gave the limited-PPM reverse sweep the same buffer for the air mass and the tracer mass (issue #57).

AtmosTransport.Tape.AirMassRead Type

julia

AirMassRead()

Read role of staged air-mass states (sweep, diffusion and convection records).

AtmosTransport.Tape.DeviceCSTapeStorage Type

julia

DeviceCSTapeStorage()

Tape storage policy that keeps staged adjoint mass states on the same backend as the source panels. This preserves the original in-memory/device-resident tape behavior while making the storage policy explicit.

AtmosTransport.Tape.FullCheckpoint Type

julia

FullCheckpoint()

Default tape schedule. Every forward step contributes records to a single tape; the reverse pass walks the tape once with no recomputation. Peak storage is proportional to nsteps.

AtmosTransport.Tape.MmapCSTapeStorage Type

julia

MmapCSTapeStorage(; dir = mktempdir(prefix="atmostransport-cstape-"),
                   cleanup_on_finalize = (dir == ""))

Tape storage policy that streams staged adjoint mass states onto a single appended raw binary file (records.bin) inside dir and emits a TOML manifest (manifest.toml) at finalisation. dir is created if it does not exist; if omitted, a unique temporary directory is created and deleted when the storage object is garbage-collected.

The policy is backend-agnostic: CPU source panels yield in-place Mmap.mmap views at read time; GPU source panels go through device-side read caches keyed by the read role (AirMassRead, TracerMassRead) and the panel shape signature (NTuple{6, NTuple{3, Int}}). The role keeps the air mass and the tracer mass that one reverse kernel reads together in distinct buffers; the shape lets a reverse loop that alternates between slot shapes reuse each cache instead of reallocating six device buffers on every read.

The on-disk format is the manifest.toml metadata plus appended raw panel payloads in records.bin.

AtmosTransport.Tape.PinnedHostCSTapeStorage Type

julia

PinnedHostCSTapeStorage()

Tape storage policy that stages tape states in host memory and reads them back through device-side read caches during the reverse pass, one per read role (AirMassRead, TracerMassRead) and panel shape, so at most two staged states are mirrored on the device at a time. The generic methods allocate plain Arrays, so the policy also works with CPU panel states; the CUDA extension specialises CuArray panel states to use pinned (CUDA.pin) host memory and a CUDA.synchronize hook.

AtmosTransport.Tape.RevolveCheckpoint Type

julia

RevolveCheckpoint()

Logarithmic-memory checkpoint schedule. The driver bisects the step range recursively: at each level it propagates forward from a saved state to the midpoint, reverses the upper half, restores the saved state, then reverses the lower half. Peak snapshot memory is the recursion depth ~ ceil(log2(nsteps)) state-pair copies; total recomputation work is ~O(nsteps × log nsteps) because each step is re-propagated once per recursion level on the path from root to leaf.

This is a tractable bisection variant of Griewank-Walther Revolve — it does NOT compute optimal binomial splits. For nsteps = 2880 (C180 / 14-day campaign) memory caps at ~12 simultaneous snapshots, which is the production target. A future commit can promote this to the optimal Revolve schedule (RevolveCheckpoint(snapshots) with a user-supplied snapshot cap) without changing the public API for single-call inversions.

Trade-off vs StrideCheckpoint(K):

Scheme	Memory	Forward work
`FullCheckpoint`	O(nsteps)	O(nsteps)
`StrideCheckpoint(K)`	O(nsteps/K)	~2 × nsteps
`RevolveCheckpoint()` (bisect)	O(log nsteps)	O(nsteps × log nsteps)

Parity. Every recorded step starts with production's step-start halo refresh (fill_cs_step_start_halos!), so the halos a step reads are a function of its interior state and a window that restarts at a bisection point records the same states as FullCheckpoint. Footprints agree with FullCheckpoint to rounding (tested at atol = 1e-12) for every supported scheme, including monotone PPM with ImplicitVerticalDiffusion and CMFMCConvection or TM5Convection. An earlier O(1e-7) Revolve drift in that configuration came from the recorders omitting that refresh (the first X sweep after convection read stale tracer halos in the full tape but fresh ones in a restarted window; issue #60), not from limiter flips.

AtmosTransport.Tape.StrideCheckpoint Type

julia

StrideCheckpoint(K)

Save a full panels_m checkpoint every K forward steps and recompute the in-window tape lazily during the reverse pass. K must be a positive integer; K = 1 is degenerate (every step is a checkpoint and the recompute factor matches FullCheckpoint); K >= nsteps collapses to a single window (also degenerate but cheap on tape).

Peak in-memory ops count and per-window tape disk are both ~K times the per-step cost; peak checkpoint memory is cld(nsteps, K) full panels_m copies.

AtmosTransport.Tape.TracerMassRead Type

julia

TracerMassRead()

Read role of the staged tracer-mass state of a limited-PPM or FV3-profile sweep record, which the reverse sweep reads together with the record's air mass.

AtmosTransport.Tape._CSAirMassResetRecord Type

julia

_CSAirMassResetRecord(mode, panels_old_m, panels_new_m)

Forward-recorded air-mass reset to a met target. :preserve_tracer_mass leaves the tracer mass unchanged (no staged air masses); :preserve_vmr scales it by new_m / old_m and stages both air masses.

AtmosTransport.Tape._CSConvectionRecord Type

julia

_CSConvectionRecord{FT, T, C, F}

Forward-recorded convection application. op is the convection operator (CMFMCConvection, TM5Convection or CMFMCMatrixConvection), forcing is the ConvectionForcing it read, panels_m is the air-mass tape slot, dt the interval it integrated over (the advection step or the met window) and substeps the CFL sub-step count the forward call used (CMFMCConvection; 1 for the others).

AtmosTransport.Tape._CSDiffusionRecord Type

julia

_CSDiffusionRecord{FT, T, D, W}

Forward-recorded implicit vertical-diffusion application. op is the diffusion operator (e.g., ImplicitVerticalDiffusion), workspace is the operator's workspace at this step, panels_m is the air-mass tape slot, dt is the (half-)timestep size and time the clock (current_time) at which the forward evaluated the operator's Kz field.

AtmosTransport.Tape._CSHaloRecord Type

julia

_CSHaloRecord(dir)

Forward-recorded cross-panel halo fill (dir ∈ (0, 1, 2): full / X-only / Y-only). The reverse pass applies the adjoint of the halo fill at the same dir.

AtmosTransport.Tape._CSMidpointRecord Type

julia

_CSMidpointRecord(step, dt)

Marker record where the surface-flux operator of model step step adds dt worth of emissions (the palindrome center). The reverse pass accumulates the step's surface footprint dt · λ at this op.

AtmosTransport.Tape._CSStepEndRecord Type

julia

_CSStepEndRecord(step)

Marker record at the end of model step step, after convection: the state that runtime output and the forward observation sampler see. Every recorder pushes one per step. The reverse pass injects the adjoint forcing of the observations sampled after step step here, so the forcing enters before the transpose of that step's convection (TM5-4DVAR adds its forcing just before the adjoint of the operator after which the forward run sampled). Without a forcing the record is a no-op.

AtmosTransport.Tape._CSSweepRecord Type

julia

_CSSweepRecord{FT, T, R, F3, S}

Forward-recorded per-direction advection sweep. direction ∈ (:x, :y, :z), scheme is the advection scheme used, panels_m / panels_rm are the staged air-mass / tracer-mass tape slots (the latter is nothing for linear schemes that don't need a tracer tape), panels_flux is the per-panel face-flux tuple, flux_scale is the scaling applied during the subcycle.

AtmosTransport.Tape._allocate_tape_slot Method

julia

_allocate_tape_slot(storage::MmapCSTapeStorage, panels)

Reserve a span of records.bin for one tape record (six panel arrays) and return a slot that can be staged into and later read from. The disk bookkeeping is backend-agnostic; backend-specific setup (read-cache allocation, synchronisation hooks) happens in _mmap_prepare_for_panels!, which the CUDA extension specialises for NTuple{6, <:CuArray}.

AtmosTransport.Tape._build_window_storage Function

julia

_build_window_storage(tape_storage, tape_path, subdir) -> AbstractCSTapeStorage

Construct a fresh tape storage instance for one window (Stride) or one base-case step (Revolve). When tape_path === nothing this falls through to _tape_storage(tape_storage) and inherits its temp-dir / cleanup behaviour. When tape_path !== nothing, tape_storage must be :mmap; the function builds an MmapCSTapeStorage rooted at joinpath(tape_path, subdir) with cleanup_on_finalize = false so the caller-owned directory is preserved past finalize_tape!.

An empty subdir (the default) is used by the single-tape FullCheckpoint path: the storage is rooted directly at tape_path with no extra nesting. The checkpoint drivers pass per-window / per-step subdirectory names (e.g. "window_00007", "step_00042") so multiple windows can coexist under the same user-supplied tree.

AtmosTransport.Tape._mmap_prepare_for_panels! Method

julia

_mmap_prepare_for_panels!(storage, panels)

Hook called once per slot allocation, before the disk reservation, that lets the storage install per-policy state (the backend of the device-side read caches, a backend synchronize callback).

The default method assumes the source panels are non-Array device-resident buffers (CuArray, MtlArray, …) and records their backend via _register_tape_read_backend!, so the read caches are allocated there on first read. The specialised method for NTuple{6, <:Array} is a no-op so the CPU read path can return mmap views directly without an intermediate copy.

The CUDA extension layers a synchronize = CUDA.synchronize hook on top of the default to amortise asynchronous device transfers; correctness of mmap → device traffic does not depend on that extension because copyto!(::CuArray, ::Array) is synchronous on CUDA.jl.

AtmosTransport.Tape._resolve_tape_path Method

julia

_resolve_tape_path(tape_storage, tape_path) -> Union{Nothing, String}

Validate that tape_path is compatible with tape_storage. Returns nothing when tape_path === nothing (use a temp dir / no disk backing). Returns the path string (after mkpath) when tape_storage is :mmap and a path is supplied. Throws an ArgumentError if a path is supplied with non-:mmap storage (since :device and :pinned_host keep tape state in memory, a path has no meaning), or if tape_storage is a pre-constructed AbstractCSTapeStorage (the storage's own configuration already determines where the tape lives).

AtmosTransport.Tape.checkpoint_window_count Method

julia

checkpoint_window_count(schedule, nsteps) -> Int

Number of forward windows for schedule covering nsteps integration steps. FullCheckpoint is 1 (the existing single tape); StrideCheckpoint(K) is cld(nsteps, K).

AtmosTransport.Tape.checkpoint_window_range Method

julia

checkpoint_window_range(schedule, window_index, nsteps) -> UnitRange{Int}

Step range covered by window_index (1-indexed) for schedule and a run of nsteps integration steps.

AtmosTransport.Tape.finalize_tape! Method

julia

finalize_tape!(storage::MmapCSTapeStorage; quiet = false, strict = false)

Sync any pending mmap state, write the TOML manifest, and close the underlying records.bin handle. Safe to call more than once; subsequent calls are no-ops. If the storage owns its temp directory and cleanup_on_finalize is true, the directory is also removed after manifest emission.

By default (strict = false), failures during manifest emission, IO close, or temp-dir cleanup are caught and surfaced via @warn (or suppressed if quiet = true). This matches the temp-dir use case where leaving a partially-finalised storage behind is benign — the directory gets cleaned up by GC or by the caller's mktempdir do block anyway.

With strict = true, manifest- and close-failures are still caught into local error slots (so the IO handle and device cache are released before throwing), but the first such error is rethrown at the end of the function. This is the right setting whenever the storage is rooted at a caller-supplied tape_path and a missing/corrupt manifest.toml would defeat the purpose of persisting the tape.

AtmosTransport.Tape.get_record Method

julia

get_record(storage::MmapCSTapeStorage, record_id::Integer) -> MmapCSTapeSlot

Return a slot descriptor for an existing record in storage. The slot exposes the same offsets/shapes/eltype fields as one produced by _allocate_tape_slot, so _tape_panels(slot, role) returns mmap views over the recorded panels.

AtmosTransport.Tape.load_mmap_tape Method

julia

load_mmap_tape(dir; readonly = true) -> MmapCSTapeStorage

Reopen a finalised MmapCSTapeStorage directory written by an earlier session. The loader parses manifest.toml, validates the format version and machine endianness, and rebuilds the in-memory record table so callers can fetch slots via get_record and read them through _tape_panels.

The returned storage has finalised = true and cleanup_on_finalize = false; it never touches the on-disk files except through mmap reads on the existing slots. With readonly = true (the default) _bump_cursor!, _allocate_tape_slot, and stage_panels! throw on mutation; readonly = false opens the binary in r+ mode but is intended for future Phase B work (slot reuse), not for forward appends — appends are blocked by the finalised flag.

Validation:

manifest.toml and records.bin must both exist under dir.
Format version must equal v1.
Endianness must match the loading machine.
meta.finalised must be true; a torn tape from an interrupted run is rejected here (Phase B is responsible for repairing it).
records.bin must be at least meta.total_bytes long; the cursor is restored from meta.total_bytes.

AtmosTransport.Tape.stage_panels! Method

julia

stage_panels!(slot::MmapCSTapeSlot, src)

Write the six source panels into the slot's reserved region of records.bin via Mmap.mmap views. Works uniformly for CPU Array and GPU CuArray sources because copyto!(::Array, ::CuArray) performs a device→host transfer.