Nested Davies relaxation via SpecificForcing on a PrescribedAtmosphere parent (fixes boundary cooling)

[ewquon]

What

For a PrescribedAtmosphere (reference-driving) parent, the nested-atmosphere interior Davies relaxation now nudges the intensive fields — specific θ, u, v — via Breeze's SpecificForcing, instead of the density-weighted ρθ/ρu/ρv. ρᵈ (mass) and the moisture density stay density-weighted; the lateral BCs are unchanged.

Closes #457.

Why: nonphysical boundary cooling

The relaxation nudged the density-weighted ρθ toward the parent. Under CompressibleDynamics, ρθ is the compressible pressure (p ∝ (ρθ)^γ), so relaxing it toward the smooth parent bled the near-wall acoustic restoring pressure — mass piled into ρᵈ while ρθ was held down, and since T = θ = ρθ/ρᵈ, the pile-up surfaced as a growing near-wall cold rim (24-h CONUS boundary cells reaching ~170 K).

Supplied under the specific key θ, Breeze wraps it in SpecificForcing, which multiplies by the live child ρᵈ at kernel time — the effective tendency is rate·mask·(ρᵈ_local·θ_parent − ρθ), so θ → θ_parent with no pressure bled and ρᵈ/momentum untouched. u/v get the same treatment for consistency. On the 24-h CONUS run this takes the worst near-wall surface cell from 171 K → 257 K and the ring-averaged near-wall dip from −13 K → −1.4 K; a minimal steady reproducer goes from −18.7 K → +2.9 K.

Generalized for telescoping nests

The intensive-vs-extensive choice is dispatched on the parent kind (davies_forcing_variables):

• ::PrescribedAtmosphere → intensive (ρᵈ, θ, u, v)
• live prognostic parent → density-weighted (ρᵈ, ρθ, ρu, ρv)

The extensive branch is a forward-looking hook for telescoping nests: nested_atmosphere_model still dispatches on ::PrescribedAtmosphere, so every call currently takes the intensive path, but once the entry point accepts a live parent nest the density-weighted relaxation activates automatically — a live parent's ρu/ρθ are a genuine prognostic solution to match at the boundary, not merely a reference.

Implementation

The state exchanger already computed the parent's intensive θ (to form ρθ = ρᵈ·θ) and reads the parent velocities u/v; it now carries all three as extra window-managed members alongside the density-weighted ones. Adds a davies_forcing_variables dispatch unit test.

(The exchanger.prognostic field naming is intentionally left as-is here; renaming it is discussed separately in #459.)

🤖 Generated with Claude Code

[codecov]

Codecov Report

✅ All modified and coverable lines are covered by tests.

📢 Thoughts on this report? Let us know!

[ewquon]

For reference, see the surface temperature deviation (first column)

before:

after:

[ewquon]

Follow-up on the #457 ocean MWE with this branch: this PR closes #457. The density-weighted relaxation was the dominant cause of the reported drift, and the small residual that remains after the fix is now fully attributed — it is expected free-interior physics of the bare-nest configuration, not a nesting defect.

Reproduction and effect of this PR

Same MWE (N = 60, 12 h, cfl = 0.5, single GPU). A density-weighted baseline reproduces #457 (surface Δθᵛ = −3.59 K vs −3.64 K reported); this branch cuts it to −0.62 K (~83 %) and collapses the vertical dipole.

The residual −0.62 K is not an LBC problem

• Temperature, not moisture: Δθ-part −0.61 K, moisture-part −0.01 K.
• ERA5 itself is not steady over the window: box-mean 1000 hPa θᵛ falls −0.24 K (00–12 UTC = local afternoon → night), so ~40 % is the child correctly tracking the parent.
• Not microphysics: microphysics = nothing → −0.627 vs −0.622 (identical drift, map, extremes).
• Not the relaxation: relaxation_rate = nothing (open BCs only) → same interior cooling (−0.84 vs −0.87).
• Not frame-origin: child-minus-parent excess at 30 m is +0.09 in the relaxation frame vs −0.59 uniformly across the free interior, including the domain center that 12 h advection cannot reach — the inverse of a boundary-defect footprint.
• Not numerics: Δt 10 → 5 s gives −0.618; WENO 5 → 9 on momentum and scalars gives −0.637 (±2 %).

What it actually is

The MWE removes every process that maintains ERA5's marine boundary layer (surface drag, surface fluxes, radiation, SGS closure). Two resolved responses follow:

1. Frictionless BL spin-up: released from ERA5's drag-balanced −5.3 m s⁻¹, the BL wind executes an inertial oscillation about the geostrophic wind (box-mean ug = −6.9 m s⁻¹): ū(30 m) reaches −9.7 m s⁻¹ by h12, matching the ≈ −9.4 m s⁻¹ inertial extreme.
2. Mixed-layer restratification: with nothing keeping the ML mixed, the SST-tied horizontal θᵛ contrasts drive thermally-direct overturning (corr(w′, θᵛ′) = +0.35; ⟨w′θᵛ′⟩ peaks +2.2×10⁻³ K m s⁻¹ ≈ 2.7 W m⁻² near 450 m). Its divergence cools the 30 m level by −0.37 K/12 h, and the initially neutral lowest ~150 m becomes stably stratified (bottom cools fastest) — the classic slumping signature.

Budget at 30 m: parent diurnal −0.24 + restratification flux −0.37 ≈ −0.61 vs −0.62 observed. The mid-troposphere warming is the retained (real) trade subsidence no longer balanced by radiative cooling; it saturates on a 4–6 h timescale as the inversion sharpens.

Positive control: re-enabling just bulk_drag (z₀ = 10⁻⁴ m) removes ~25 % of the interior drift — the spin-up-coupled part; the remainder persists because drag restores no stirring or surface heating. Full restoration needs the full BL budget, i.e. the real downscaling configuration.

For reference, a survey of LBC practice (WRF, MPAS-A, ERF, FV3-LAM, ICON-LAM, COSMO, ALADIN) shows the scheme in this PR matches established practice — closest to ICON-LAM (intensive θᵥ/u/v, density relaxed) — while the pre-PR density-weighted-θ-without-comoving-mass configuration has no analog in any surveyed model.

Bare-nest takeaway: with no BL physics, the interior should not hold ERA5's flux-maintained near-surface state; −0.6 K/12 h is its natural inviscid restratification plus ERA5's own diurnal cycle. bulk_drag (+ closure/radiation as available) restores the budget.

🤖 Generated with Claude Code