rdb_ocean_vdiff Module

Constant momentum vertical viscosity (m^2/s). Zero is a no-op.

Constant tracer vertical diffusivity (m^2/s). Zero is a no-op.

Sub-diagonal for the tracer solve. Shape (nx, ny, nz).

East-face (u) tridiagonal scratch.

North-face (v) tridiagonal scratch.

Main diagonal.

East-face (u) tridiagonal scratch.

North-face (v) tridiagonal scratch.

MOM6 DRAG_BG_VEL (&ocean_bdrag_nml bg_vel).

MOM6 CDRAG. Quadratic: &ocean_bdrag_nml cd. Linear: the equivalent r·hbbl/bg_vel, so that CDRAG·DRAG_BG_VEL = r·hbbl — the stress of the HBBL-distributed linear drag, r·hbbl·u.

Cell salinity concentration, as bbl_conc_t.

Cell temperature CONCENTRATION (nx, ny, nz) for the BBL stratification measure, read through the I1′ column rule (rdb_vl_column_conc: a filler reads its donor). Workspace of vdiff_set_viscous_bbl; (1,1,1) unless bbl_per_face.

Drag law of the BBL (BBL_FORM_*).

MOM6 bottomdraglaw parity for the interface coupling (&ocean_vdiff_nml bbl_glue, PGF_BUG.md §9). MOM6 holds a layered basin at rest NOT by computing a clean PGF — its PFu over grounded/sliver layers is identical to ours — but by absorbing it viscously every step (du_dt_visc = −PFu to machine precision, measured): find_coupling_coef raises the interface coupling to a_cpl = kv_bbl/h_shear within the bottom boundary layer (botfn weight, z_i accumulated from HARMONIC thicknesses so grounded stacks sit at z≈0), and the BED coupling is the piston kv_bbl/(h₁/2) — unbounded as the bottom layer thins (measured a_cpl 3e7–1.3e8 m/s on the rest-state reproducer). This knob ports both halves into the momentum tridiagonal: * interface kv → kv + (kv_bbl − KV)·botfn (KV = kv_bbl_bg), with h_shear capped toward bbl_thick under the same botfn weight; * bed row → dt·kv_bbl/(hf₁·(min(hvel₁/2, bbl_thick))) — THE bed sink whenever the glue is on: it replaces the dt·λ_bot Rayleigh fold, and the driver skips the explicit drag apply (the explicit du_drag still feeds the barotropic F_slow, the split implicit_drag takes). kv_bbl and bbl_thick are PER FACE: vdiff_set_viscous_bbl (MOM6 set_viscous_BBL, quadratic or linear drag law, KW99 rotation/stratification-limited thickness) fills them once per outer step when bbl_per_face (configure: a bottom drag is configured). A hand-built slot without bbl_per_face gets the historical constants bbl_piston·hbbl_visc / hbbl_visc. Requires hvel_mom6 (supplies the height-above-bed bookkeeping), enforced at configure. With no bottom drag configured the configure latch turns it off (nothing to glue).

MOM6 HBBL (&ocean_bdrag_nml hbbl): the thickness over which the near-bottom velocity is averaged for the drag law.

The MOM6 bottom boundary layer (set_viscous_BBL) is computed per face once per outer step by vdiff_set_viscous_bbl into kv_bbl_u/v and bbl_thick_u/v, which the glue then reads. Set at configure when bbl_glue is on AND a bottom drag is configured with &ocean_bdrag_nml hbbl > 0 (MOM6 requires HBBL). .false. ⇒ the glue uses the historical constants bbl_piston·hbbl_visc / hbbl_visc, refilled each call (the hand-built test path).

BBL drag piston velocity u* (m/s) for bbl_glue — MOM6 CDRAG·DRAG_BG_VEL (0.003·0.1 with reference defaults). kv_bbl = bbl_piston·hbbl_visc.

Boussinesq reference density of the BBL stratification measure.

MOM6 RiNo_mix (= kappa-shear on): the BBL is capped at 0.5·HBBL.

MOM6 BBL_THICK_MIN (&ocean_bdrag_nml bbl_thick_min; MOM6 default 0).

BBL thickness (m) at u-faces, (nx+1, ny).

BBL thickness (m) at v-faces, (nx, ny+1).

Super-diagonal (overwritten by c_prime).

East-face (u) tridiagonal scratch.

North-face (v) tridiagonal scratch.

Bottom-boundary-layer scale for the botfn blend (MOM6 HBBL, 10.0 m in the double_gyre reference) when the BBL glue is off. Under the per-face glue the blend is normalised by each face’s bbl_thick instead (MOM6 I_Hbbl = 1/bbl_thick), and this is only the HBBL fallback for a drag configured with &ocean_bdrag_nml hbbl = 0 (MOM6 has ONE HBBL).

Which MOM6 face-thickness branch hvel_mom6 builds. .false. (default) is MOM6’s DEFAULT, HARMONIC_VISC = False (vertvisc_coef): the face thickness is the ARITHMETIC mean, blended toward the harmonic mean near the bed when the flow runs from the thin side to the thick one, and the height above the bed is max(zh, z_clear) — z_clear the height of the higher of the two cells’ interfaces above the SHALLOWER of the two beds. So every face layer below the shallower bed of a step sits at z_i ≈ 0, inside the bottom boundary layer, and the BBL glue couples it. .true. is MOM6 HARMONIC_VISC = True: harmonic face thickness, blended toward arithmetic when the flow runs thick -> thin, height above bed from harmonic thicknesses alone (the historical roundabout hvel_mom6, gated by hvel_upwind).

MOM6 HARMONIC_VISC = True parity for the MOMENTUM face thickness (hvel). Roundabout’s h_u is the ARITHMETIC mean unconditionally while MOM6’s is the HARMONIC mean blended back toward arithmetic near the bed when the flow runs thick->thin – and MOM6’s h_shear is then the ARITHMETIC mean of those hvels. The two means are effectively SWAPPED relative to MOM6. For a grounded sliver harm(1e-10, 20) ~ 2e-10 against arith ~ 10: eleven orders, and that gap IS the ~1e19 suppression that renders the spurious sliver PGF harmless in MOM6.

Near-bed upwind (arithmetic-donor) blend inside the hvel_mom6 face-thickness build. .true. = MOM6 parity (the historical hvel_mom6 behaviour, bit-identical). .false. = pure harmonic hvel: the blend’s u-sign test flip-flops on roundoff velocities at (near-)rest and collapses the BBL glue at flipped faces (PGF_BUG.md §9.8) — the rest-state envelope runs with it off.

Fold the quadratic / linear bottom drag into the backward-Euler vertical-friction tridiagonal as a stress bottom-BC diagonal coupling (&ocean_vdiff_nml implicit_drag). Roundabout bed is k = 1: the drag rate λ_bot adds to the b_diag(:, :, 1) diagonal only (a positive add ⇒ unconditionally stable on thin bottom layers where the explicit u·(1 − dt·c_d|U|/h) reverses sign once dt·c_d|U|/h > 2). λ_bot is the SAME Rayleigh rate (c_d·|U_bbl|/h_1 quadratic, r linear, |U| frozen at uⁿ) the bottom-drag slot already forms; the row is already normalized by h_1 so the diagonal increment is dt·λ_bot directly. Mutually exclusive with &ocean_bdrag_nml implicit (configure fails loud); the explicit drag apply is gated off when on. Default .false. ⇒ bit-identical.

Fold the surface wind stress into the backward-Euler vertical- friction tridiagonal as a Neumann top-BC right-hand-side term (&ocean_vdiff_nml implicit_stress). Roundabout is bottom-up, so the SURFACE is k = nz: the kinematic stress τ/ρ₀ adds to the rhs(:, :, nz) row only (the diagonal is unchanged), filtered through the implicit operator together with the interior shear. When .true. the explicit surface-stress pre-solve apply in the dyn run-stage is suppressed (the driver gates it) so the forcing is not double-counted. Default .false. ⇒ matrix + RHS built exactly as before ⇒ bit-identical to the explicit path.

Fold the ICE-SHELF TOP drag into the k = nz DIAGONAL (&ocean_vdiff_nml implicit_top_drag) instead of the explicit pre-solve add, and MASK the wind-stress RHS off on the faces the ice covers. The mirror of implicit_drag at the other end of the column: a drag is a diagonal term and the wind is an RHS term, so the two share the surface row without either being approximated, and the cover mask is what makes the sharing physical (there is no atmosphere under a shelf). λ_top is the SAME Rayleigh rate (C_d·|U|/h_nz quadratic, r linear, |U| frozen at uⁿ) the top-drag slot forms, and it arrives already cover- and wet-masked. Requires &ocean_tdrag_nml enable; mutually exclusive with &ocean_tdrag_nml implicit and with htbl > 0; all fail loud at configure. The explicit top-drag apply is gated off when on. Default .false. ⇒ bit-identical.

True between init and destroy.

MOM6 KV: the background viscosity the BBL viscosity REPLACES within the boundary layer, Kv_tot = Kv_tot + (kv_bbl − KV)·botfn (find_coupling_coef) — the shear / boundary-layer contributions are kept. Set from &ocean_vmix_nml kv_bg.

BBL viscosity kv_bbl (m²/s) at u-faces, (nx+1, ny).

BBL viscosity at v-faces, (nx, ny+1).

Diffusivity workspace used when no 3D source is provided.

Right-hand side / solution.

East-face (u) tridiagonal scratch.

North-face (v) tridiagonal scratch.

MOM6 HARMONIC_VISC analogue. When .true. the implicit vdiff solver uses the harmonic mean of adjacent layer thicknesses in the face-thickness denominator instead of the arithmetic mean. Better-conditioned at thin / vanishing layers — harm = 2·h_a·h_b / (h_a + h_b) stays small when one neighbour is small, whereas arith = 0.5·(h_a + h_b) is dominated by the thicker neighbour and produces stiff coefficients that the Thomas solve can’t handle gracefully. Default .false. preserves the arithmetic-mean behaviour bit-identically.

&vcoord_nml zfixed_closed_faces — z-level partial steps. Latched by configure_ocean_closed_faces, forwarded as a plain scalar to diffuse_velocity_columns_impl, where it switches the face column to min(h_L, h_R) and CUTS the tridiagonal coupling across any interface touching a closed (inert-filler-on-one-side) layer. See that routine’s zlevel_faces docstring for the full argument. Default .false. ⇒ bit-identical.

Boussinesq reference density for the implicit surface-stress fold, dt*(tau/rho0)/h_nz. READ ONLY when that fold is active (do_stress below); every other path never reaches it. There is deliberately NO defaulted value: the one rho0 of record is &ocean_ic_nml rho_0 -> eos%rho0, which the production call site hands over as surface_stress%rho0 (configure_ocean_reference_density). Omitting it WITH the fold active is a wiring bug and fails loud — it used to fall back to a silent literal 1035, which is how a run configured at a different rho_0 could fold its wind stress in on the wrong one.

Explicit-shape for the reason spelled out on the argument list above.

.true. = build the matrices and (re)fill visc_rem_u/v WITHOUT solving for / modifying the velocities — the pre-substep visc_rem refresh (PGF_BUG.md §9; MOM6 computes vertvisc_coef before btstep, so its BT weights never lag). Requires both visc_rem_u/v present. Default .false..

Corner-staggered interface viscosity source, (nx+1, ny+1, nz+1) — see the corner add-on note above. Typically the vertex kappa-shear kd_corner carrier (device-resident). Explicit-shape for the reason spelled out on the argument list above.

Scale applied to kv_corner_source (Kv = Pr·Kd; the vertex kappa-shear prandtl_turb). Default 1.

Temperature diffusivity, (nx, ny, nz+1), interface-located.

Salinity + passive-tracer diffusivity, same shape. MUST be accompanied by kt_source (case 2 fails loud).

Cell-centred extents (ghosts included) of the workspace.

Cell-centred extents (ghosts included) of the workspace.

Cell-centred extents (ghosts included) of the workspace.

BBL_FORM_LINEAR or BBL_FORM_QUADRATIC.

Quadratic drag coefficient (dimensionless).

Linear drag rate (1/s) of the HBBL-distributed linear drag.

MOM6 HBBL (m).

MOM6 DRAG_BG_VEL (m/s).

MOM6 BBL_THICK_MIN (m).

MOM6 RiNo_mix (kappa-shear on).

Boussinesq reference density (kg/m^3).

MOM6 KV, the background viscosity (m^2/s).

The run’s equation of state, for ∂ρ/∂T and ∂ρ/∂S.

Coriolis parameter at cell SW corners (1/s), cor%f_corner.

Cell T/S concentrations; read only when use_eos (then (nx, ny, nz); a (1,1,1) placeholder otherwise).

Cell-centred wet/dry land mask (1 wet, 0 land). The wind-stress fold multiplies tau_face by the face mask min of the two bounding cells — matching the explicit surface-stress kernel — so no stress is injected at a no-normal-flow land face. All-wet (mask ≡ 1) ⇒ bit-identical. Only read when do_stress.

Fold the surface wind stress into the k = nz RHS row.

Fold the bottom drag into the k = k_bot_face diagonal.

Boussinesq reference density for the τ/ρ₀ stress conversion.

Wind stress (N/m²) on this face. Present iff do_stress.

Bottom-drag Rayleigh rate λ (1/s) on this face. Present iff do_drag.

Fold the ice-shelf top drag into the k = nz diagonal AND mask the wind-stress RHS by (1 − cover_face).

Top-drag Rayleigh rate λ (1/s) on this face. Present iff do_top.

Face ice-cover mask (0 open, 1 under ice), the OR of the two abutting cells (rdb_ocean_top_drag). Present iff do_top.

.false. = remnant-only mode: build the matrix and compute visc_rem_out WITHOUT touching u_face (the pre-substep visc_rem refresh, PGF_BUG.md §9 — MOM6 computes vertvisc_coef before btstep every stage, so its BT weights never lag; the stage-end producer alone leaves visc_rem ≡ 1 for the whole first stage, and the Δu corrector then deposits the spurious column-mean into wet layers). .true. = the normal solve.

Fill visc_rem_out with the viscous remnant γ_k — the sensitivity of the post-friction layer velocity to a uniform barotropic acceleration, γ_k ≡ (1/Δt)·∂u_k^{n+1}/∂Ā. γ solves the SAME tridiagonal system as the momentum solve (same a_diag/b_diag, and the already-factorized c_diag left behind by the momentum forward sweep) with RHS ≡ 1 — Roundabout’s rows are pre-normalized by h_k, so the remnant RHS is the unit vector, NOT h_u(k) (MOM6’s un-normalized convention). Reusing the surviving factorization means γ cannot drift from the operator actually solved for momentum.

Output γ_k, bottom-up (k = 1 bed → γ smallest, k = nz surface → γ → 1). Clamped to min(γ, 1.0) at production (cheap FP insurance on a quantity the maximum principle already bounds to (0, 1]). Present iff do_remnant.

MOM6 HARMONIC_VISC parity for h_u + h_shear (see the slot-type docstring). .false. => the historical arithmetic-h_u / face_thick-dz pair, bit-identical.

Bottom-layer scale for the botfn blend (MOM6 HBBL).

MOM6 bottomdraglaw coupling parity (see the slot-type docstring). Configure guarantees hvel_mom6 when this is on; the bed piston then applies whether or not do_drag. .false. ⇒ bit-identical.

.false. = skip the near-bed upwind blend (pure harmonic hvel). See the slot-type docstring. hvel_harmonic branch only.

MOM6 HARMONIC_VISC: .false. = MOM6’s default z_clear face-thickness branch, .true. = the harmonic branch (slot-type docstring). Only read when hvel_mom6.

BBL viscosity kv_bbl (m²/s) of this face (MOM6 visc%Kv_bbl_u). Only read when bbl_glue.

BBL thickness (m) of this face (MOM6 visc%bbl_thick_u); also the height-above-bed normaliser 1/I_Hbbl under the glue, as in MOM6’s bottomdraglaw branch. Only read when bbl_glue.

MOM6 KV: the background viscosity the BBL viscosity replaces within botfn reach of the bed, Kv_tot + (kv_bbl − KV)·botfn.

Add the corner-staggered viscosity to every face interface. Gated INSIDE the single DC (no split loop — NVHPC penalty); kv_corner is guaranteed present when do_corner.

Scale on kv_corner (Kv = Pr·Kd). Only read when do_corner.

Corner-staggered interface viscosity (SW-corner convention: corner (i,j) is the SW corner of cell (i,j)). A face reads the 2-point average of its two END corners: u-face (i,j) → corners (i,j)/(i,j+1); v-face (i,j) → corners (i,j)/(i+1,j) — the direct corner→face route, never via a tracer point. Present iff do_corner.

&vcoord_nml zfixed_closed_faces — z-level partial steps.

ms%k_top_u / k_top_v — the index of the first layer that is LIVE on BOTH sides of this face, nz wherever nothing vanishes against the top. The surface row of the column solve, the ice-shelf top-drag diagonal fold and the wind-stress RHS fold all sit on THIS row instead of on nz: under a quasi-geopotential coordinate beneath a shelf the rows above it are inert fillers, decoupled to the identity by the zlevel_faces gates, so a dt*lambda_top added at nz would damp a row that is not coupled to the ocean at all — and the drag rate lambda_top is itself captured at k_top_face by rdb_ocean_top_drag, so the two MUST agree. With k_top_face ≡ nz every expression below is the one it replaced, bit for bit.

ms%k_bot_u / k_bot_v — the first layer LIVE on BOTH sides of this face counting UP from the bed (max of its two columns), 1 wherever nothing vanishes against the bed. The bed row of the column solve (no-flux-below BC, the implicit bottom-drag diagonal fold and the bbl_glue piston) sits on THIS row instead of on k = 1, and the rows below it — inert z_fixed bed fillers — are the identity. The drag rate lambda_bot is itself captured at k_bot_face by ocean_bottom_drag_compute_tendencies, so the two MUST agree. The height-above-bed stack zint (BBL glue, upwind blend) also starts accumulating here. With k_bot_face ≡ 1 every expression below is the one it replaced, bit for bit.