rdb_ocean_mle Module

Fox-Kemper, Ferrari & Hallberg (2008) submesoscale mixed-layer-eddy (MLE) restratification. Submesoscale eddies slump lateral buoyancy fronts in the surface mixed layer via an overturning streamfunction Psi(z) = Ce * (H_ml^2 / |f|) * (grad b_bar x z_hat) * mu(z) injected as extra ML-confined per-layer mass transports uhml/vhml added to ms%mass_flux_{x,y}_layer BEFORE the continuity divergence + PPM tracer advection — never touching the velocity fields. The per-layer weights a(k) satisfy sum_k a(k) = mu(0) - mu(-1) = 0 exactly (a closed overturning cell ⇒ mass/tracer conservative).

Vertical convention: k=1 bed, k=nz surface. The ML band walks from the surface down; the FK sigma coordinate is 0 at the surface, -1 at the ML base. H_ml is the EPBL epbl%mld. Timescale forms (use_mom_mixrate): bare Ce/max(|f|,f_floor) (default), or FK11 momentum-mixrate (production-recommended; suppresses restratification under vigorous mixing). Default off ⇒ bit-identical.

Partial-step z-level faces (&vcoord_nml zfixed_closed_faces)

Under z_fixed a face column is not the whole water column: a layer that is an inert FILLER on either side (inside the bed or the ice draft) is a WALL for that layer at that face (metrics%open_u/open_v == 0). Continuity applies that mask to the resolved flux BEFORE mle_fold_x/y adds uhml/vhml, so the overturning must be built on the OPEN part of each face column here. With the knob on:

  • the ML walk (b_bar, htot) at cell centres SKIPS every non-live layer (rdb_vl_is_live), so no ML property is read from a filler and the walk starts at the first live layer from the top (the free surface, or the ice base under a draft);
  • each face marks its open layers ok(k) = open(k) .and. live(h_W(k)) .and. live(h_E(k)) (the set GM’s open column and Redi’s open window use) and walks the FK sigma coordinate over the OPEN face thickness only (hf(k) = 0 off the open set ⇒ sigma is carried unchanged across it ⇒ a(k) = mu(s) - mu(s) = 0): sigma = 0 at the TOP of the open column, so a filler above it (ice draft) carries nothing;
  • H_vel is clamped to the open-column thickness Sum_k hf(k). The sigma walk must reach -1 inside the open column for Sum_k a(k) = mu(0) - mu(-1) = 0; with the ML deeper than a face’s open column (a deep column abutting a shallow step) the unclamped walk stops above -1 and the face transport no longer closes (prototype python_prototypes/mle_zfixed: Sum a = -0.93). The clamped H_vel also sets the H_vel^2 amplitude and the timescale’s H_vel: the cell a face can carry is at most as deep as its open column;
  • uhml/vhml are written EXACTLY zero off the open set.

Hence the FK transport is zero on every closed face-layer and every filler, its face streamfunction vanishes at the top of the open column and below the ML (at the latest at the open column’s bottom), and Sum_k uhml = 0 at every face. The open face thickness is staged in uhml/vhml themselves by a host-gated pre-pass, so the transport kernels never name the (1,1,1) mask placeholders with the knob off. Knob OFF ⇒ byte-identical.

References: Fox-Kemper, Ferrari & Hallberg (2008); Fox-Kemper et al. (2011).


Uses

  • module~~rdb_ocean_mle~~UsesGraph module~rdb_ocean_mle rdb_ocean_mle iso_fortran_env iso_fortran_env module~rdb_ocean_mle->iso_fortran_env module~rdb_constants rdb_constants module~rdb_ocean_mle->module~rdb_constants module~rdb_grid rdb_grid module~rdb_ocean_mle->module~rdb_grid module~rdb_mem_report rdb_mem_report module~rdb_ocean_mle->module~rdb_mem_report module~rdb_multilayer_state rdb_multilayer_state module~rdb_ocean_mle->module~rdb_multilayer_state module~rdb_ocean_boundary_types rdb_ocean_boundary_types module~rdb_ocean_mle->module~rdb_ocean_boundary_types module~rdb_ocean_epbl rdb_ocean_epbl module~rdb_ocean_mle->module~rdb_ocean_epbl module~rdb_ocean_metrics rdb_ocean_metrics module~rdb_ocean_mle->module~rdb_ocean_metrics module~rdb_ocean_surface_stress rdb_ocean_surface_stress module~rdb_ocean_mle->module~rdb_ocean_surface_stress pic_types pic_types module~rdb_constants->pic_types module~rdb_grid->module~rdb_constants module~rdb_mem_report->iso_fortran_env module~rdb_mem_report->module~rdb_constants pic_logger pic_logger module~rdb_mem_report->pic_logger pic_strings pic_strings module~rdb_mem_report->pic_strings module~rdb_multilayer_state->iso_fortran_env module~rdb_multilayer_state->module~rdb_constants module~rdb_multilayer_state->module~rdb_grid module~rdb_multilayer_state->module~rdb_mem_report module~rdb_efp rdb_efp module~rdb_multilayer_state->module~rdb_efp module~rdb_error_ring rdb_error_ring module~rdb_multilayer_state->module~rdb_error_ring module~rdb_tracer rdb_tracer module~rdb_multilayer_state->module~rdb_tracer module~rdb_multilayer_state->pic_logger module~rdb_ocean_boundary_types->iso_fortran_env module~rdb_ocean_boundary_types->module~rdb_constants module~rdb_ocean_boundary_types->module~rdb_grid module~rdb_ocean_boundary_types->module~rdb_mem_report module~rdb_ocean_status rdb_ocean_status module~rdb_ocean_boundary_types->module~rdb_ocean_status module~rdb_ocean_tide_astro rdb_ocean_tide_astro module~rdb_ocean_boundary_types->module~rdb_ocean_tide_astro pic_ascii pic_ascii module~rdb_ocean_boundary_types->pic_ascii module~rdb_ocean_boundary_types->pic_logger module~rdb_ocean_epbl->iso_fortran_env module~rdb_ocean_epbl->module~rdb_constants module~rdb_ocean_epbl->module~rdb_grid module~rdb_ocean_epbl->module~rdb_mem_report module~rdb_ocean_epbl->module~rdb_multilayer_state module~rdb_ocean_epbl->module~rdb_ocean_surface_stress module~rdb_eos rdb_eos module~rdb_ocean_epbl->module~rdb_eos module~rdb_ocean_surface_flux rdb_ocean_surface_flux module~rdb_ocean_epbl->module~rdb_ocean_surface_flux module~rdb_scratch_3d rdb_scratch_3d module~rdb_ocean_epbl->module~rdb_scratch_3d module~rdb_ocean_metrics->iso_fortran_env module~rdb_ocean_metrics->module~rdb_constants module~rdb_ocean_metrics->module~rdb_grid module~rdb_ocean_metrics->module~rdb_mem_report module~rdb_ocean_metrics->module~rdb_error_ring module~rdb_io_netcdf rdb_io_netcdf module~rdb_ocean_metrics->module~rdb_io_netcdf module~rdb_ocean_bipolar rdb_ocean_bipolar module~rdb_ocean_metrics->module~rdb_ocean_bipolar module~rdb_ocean_fold rdb_ocean_fold module~rdb_ocean_metrics->module~rdb_ocean_fold module~rdb_ocean_metrics->module~rdb_ocean_status netcdf netcdf module~rdb_ocean_metrics->netcdf module~rdb_ocean_metrics->pic_logger module~rdb_ocean_metrics->pic_strings module~rdb_ocean_surface_stress->iso_fortran_env module~rdb_ocean_surface_stress->module~rdb_constants module~rdb_ocean_surface_stress->module~rdb_grid module~rdb_ocean_surface_stress->module~rdb_mem_report module~rdb_ocean_surface_stress->module~rdb_multilayer_state module~rdb_ocean_surface_stress->module~rdb_scratch_3d module~rdb_efp->iso_fortran_env ieee_arithmetic ieee_arithmetic module~rdb_efp->ieee_arithmetic module~rdb_eos->module~rdb_constants module~rdb_eos->module~rdb_grid module~rdb_error_ring->pic_logger module~rdb_io_netcdf->iso_fortran_env module~rdb_io_netcdf->module~rdb_constants module~rdb_io_netcdf->module~rdb_error_ring module~rdb_io_netcdf->netcdf module~rdb_io_netcdf->pic_logger module~rdb_io_netcdf->pic_strings iso_c_binding iso_c_binding module~rdb_io_netcdf->iso_c_binding module~rdb_ocean_bipolar->module~rdb_constants module~rdb_ocean_fold->module~rdb_constants module~rdb_ocean_surface_flux->iso_fortran_env module~rdb_ocean_surface_flux->module~rdb_constants module~rdb_ocean_surface_flux->module~rdb_grid module~rdb_ocean_surface_flux->module~rdb_mem_report module~rdb_ocean_surface_flux->module~rdb_multilayer_state module~rdb_ocean_tide_astro->iso_fortran_env module~rdb_ocean_tide_astro->module~rdb_constants module~rdb_scratch_3d->iso_fortran_env module~rdb_scratch_3d->module~rdb_constants module~rdb_scratch_3d->module~rdb_mem_report module~rdb_tracer->iso_fortran_env module~rdb_tracer->module~rdb_constants module~rdb_tracer->module~rdb_grid module~rdb_tracer->module~rdb_mem_report

Used by

  • module~~rdb_ocean_mle~~UsedByGraph module~rdb_ocean_mle rdb_ocean_mle module~rdb_continuity rdb_continuity module~rdb_continuity->module~rdb_ocean_mle module~rdb_ocean_dyn rdb_ocean_dyn module~rdb_ocean_dyn->module~rdb_ocean_mle module~rdb_ocean_dyn->module~rdb_continuity module~rdb_ocean_state rdb_ocean_state module~rdb_ocean_state->module~rdb_ocean_mle module~rdb_ocean_state->module~rdb_continuity module~rdb_ocean_state->module~rdb_ocean_dyn module~rdb_driver rdb_driver module~rdb_driver->module~rdb_ocean_dyn module~rdb_driver->module~rdb_ocean_state module~rdb_ocean_engine rdb_ocean_engine module~rdb_driver->module~rdb_ocean_engine module~rdb_handle rdb_handle module~rdb_handle->module~rdb_ocean_state module~rdb_handle->module~rdb_ocean_engine module~rdb_ice_transport rdb_ice_transport module~rdb_ice_transport->module~rdb_continuity module~rdb_ocean_api rdb_ocean_api module~rdb_ocean_api->module~rdb_ocean_dyn module~rdb_ocean_api->module~rdb_handle module~rdb_ocean_diag_derived rdb_ocean_diag_derived module~rdb_ocean_api->module~rdb_ocean_diag_derived module~rdb_ocean_diag_fills rdb_ocean_diag_fills module~rdb_ocean_api->module~rdb_ocean_diag_fills module~rdb_ocean_api->module~rdb_ocean_engine module~rdb_ocean_diag_derived->module~rdb_ocean_state module~rdb_ocean_diag_derived->module~rdb_ocean_diag_fills module~rdb_ocean_diag_fills->module~rdb_ocean_state module~rdb_ocean_engine->module~rdb_ocean_dyn module~rdb_ocean_engine->module~rdb_ocean_state module~rdb_ocean_engine->module~rdb_ice_transport module~rdb_ocean_engine->module~rdb_ocean_diag_derived module~rdb_ocean_engine->module~rdb_ocean_diag_fills module~rdb_ocean_setup rdb_ocean_setup module~rdb_ocean_engine->module~rdb_ocean_setup module~rdb_ocean_setup->module~rdb_ocean_dyn module~rdb_ocean_setup->module~rdb_ocean_state

Variables

Type Visibility Attributes Name Initial
real(kind=wp), private, parameter :: H_NEGLECT = 1.0e-30_wp

Empty-column guard (spec 3.1 / 3.4).
real(kind=wp), private, parameter :: MLE_H_AVAIL_MIN = 1.0e-6_wp

Floor for the per-layer transport availability cap (m): a donor layer may never be drained below this by the FK overturning. Matched to the continuity / windowed-drain h_min (1e-6 m) so the cap and the drain limiter protect the same minimum thickness.
real(kind=wp), private, parameter :: PI = 3.14159265358979323846_wp
real(kind=wp), private, parameter :: VONKAR = 0.41_wp

von Karman constant for the FK11 momentum-mixrate form.

Derived Types

type, public ::  ocean_mle_t

Fox-Kemper mixed-layer-eddy restratification state. All fields default to the inert (enable=.false.) configuration so an ocean run that never sets &ocean_foxkemper_nml is bit-identical.

Components

Type Visibility Attributes Name Initial
real(kind=wp), public, allocatable :: b_ml(:,:)

ML-averaged buoyancy b_bar (m/s^2), (nx,ny).
real(kind=wp), public :: bodner_mstar = 0.5_wp

Mechanical (u*) weight in w'u' = (mstar·u*³+nstar·w*³)^{2/3}.
real(kind=wp), public :: bodner_nstar = 0.066_wp

Convective (w*) weight; w*³ = max(0,-b0)·mld from epbl%b0.
real(kind=wp), public :: ce = 0.0625_wp

FK08 coefficient Ce (typical 0.06-0.08).
real(kind=wp), public :: cr = 0.0_wp

Bodner efficiency coefficient Cr (nondim).
logical, public :: enable = .false.

Master switch. Off => mle_compute_transports is a no-op and the folds add nothing => bit-identity preserved.
real(kind=wp), public :: f_floor = 1.0e-5_wp

|f| regularisation floor (1/s); ~|f| at 4 degN.
real(kind=wp), public, allocatable :: htot_ml(:,:)

Accumulated ML thickness on the grid (<= mld) (m), (nx,ny).
logical, public :: is_init = .false.

True between init and destroy; tracks GPU attachment.
real(kind=wp), public :: min_wstar2 = 1.0e-24_wp

Floor on w'u' (m²/s²), pure 1/0 armour.
real(kind=wp), public :: mld_decay_time = 0.0_wp

Running-mean MLD filter time-scale (s). 0 ⇒ filter off (instantaneous EPBL MLD, bit-identical). Positive ⇒ the MLD resets instantly to a deeper value but decays over this scale when the diagnosed MLD retreats, bounding Psi ~ MLD² growth under a thinning/oscillating boundary layer.
real(kind=wp), public, allocatable :: mld_filtered(:,:)

Running-mean filtered MLD (m), (nx,ny); persistent across steps. Only used when mld_decay_time > 0. Lazily seeded from the first instantaneous MLD (the < 0 sentinel marks the unseeded state so the running mean starts from the true value rather than a spurious zero — MOM6 seeds from a restart field).
integer, public :: nx_total = 0
integer, public :: ny_total = 0
integer, public :: nz_ml = 0
logical, public :: resolution_taper = .false.

Resolution-function taper hook. Hard config error if on without the resolution function (deferred).
real(kind=wp), public :: tail_dh = 0.0_wp

mu cubic-tail extension below the ML base. Default 0 ⇒ exact FK08 mu (the tail itself is deferred).
real(kind=wp), public, allocatable :: uhml(:,:,:)

FK x-face transport (m^3/s), (nx+1,ny,nz_ml). Filled once per outer step at thermo cadence, then folded into continuity every dynamics call via mle_fold_x. At dt_therm_ratio > 1 the stale values are re-folded on non-thermo steps (over-applies FK; still conservative since sum_k a(k) = 0).
logical, public :: use_bodner = .false.

Bodner et al. (2023) frontogenesis-arrest variant. When true the Ce/|f| timescale is replaced by ts_bod = Cr·Δs·|f|·h/w'u' (a timescale [s], so it drops into the same uDml formula); the frontal-arrest length enters inline as |f|·h/w'u'. Default off.
logical, public :: use_mom_mixrate = .false.

Use the FK11 momentum-mixrate timescale instead of the bare Ce/|f| floor form. Production-recommended; default off.
real(kind=wp), public, allocatable :: vhml(:,:,:)

FK y-face transport (m^3/s), (nx,ny+1,nz_ml); same limitation.

Type-Bound Procedures

procedure, public, non_overridable :: bytes => ocean_mle_bytes
procedure, public :: destroy => ocean_mle_destroy
procedure, public :: enter_data => ocean_mle_enter_data
procedure, public :: exit_data => ocean_mle_exit_data
procedure, public :: init => ocean_mle_init

Functions

public pure function mle_mu_shape(sigma) result(mu)

FK08 second-order vertical structure function mu(sigma). sigma in [-1,0]: 0 = surface interface, -1 = ML base. mu = 0 at sigma=0 (surface) and sigma <= -1 (ML base); peaks near sigma=-0.5.

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: sigma

Return Value real(kind=wp)

private pure function mle_bodner_timescale(f_abs, ustar, h_vel, b0_face, ds, cr, mstar, nstar, min_wstar2) result(ts)

Bodner et al. (2023) frontogenesis-arrest MLE timescale [s]: ts = Cr · ds · |f| · h / w’u’ where the frontal-arrest length enters inline as |f|·h/w'u' and w’u’ = max( (mstar·u³ + nstar·w³)^(2/3), min_wstar2 ), w³ = max(0, -b0)·h (destabilizing buoyancy flux only) ds = sqrt(0.5(dx²+dy²)) is the grid-scale front width. Because the product is dimensionally a TIME, this drops straight into the FK transport form uDml = ts·dyCu·idxCu·db·H² (the swap that turns classic Fox-Kemper into Bodner). w'u' is floored so a quiescent, unforced column (u→0, b0→0) never divides by zero.

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: f_abs
real(kind=wp), intent(in) :: ustar
real(kind=wp), intent(in) :: h_vel
real(kind=wp), intent(in) :: b0_face
real(kind=wp), intent(in) :: ds
real(kind=wp), intent(in) :: cr
real(kind=wp), intent(in) :: mstar
real(kind=wp), intent(in) :: nstar
real(kind=wp), intent(in) :: min_wstar2

Return Value real(kind=wp)

private pure function mle_face_ustar_x(ss, rho0, i, j) result(ustar)

Friction velocity u* = sqrt(|tau|/rho0) at the u-face from the surface wind stress, averaged onto the face. Only used by the FK11 mixrate form. Returns 0 if stress fields are absent.

Arguments

Type IntentOptional Attributes Name
type(ocean_surface_stress_t), intent(in) :: ss
real(kind=wp), intent(in) :: rho0
integer, intent(in) :: i
integer, intent(in) :: j

Return Value real(kind=wp)

private pure function mle_face_ustar_y(ss, rho0, i, j) result(ustar)

u* at the v-face; mirror of mle_face_ustar_x.

Arguments

Type IntentOptional Attributes Name
type(ocean_surface_stress_t), intent(in) :: ss
real(kind=wp), intent(in) :: rho0
integer, intent(in) :: i
integer, intent(in) :: j

Return Value real(kind=wp)

private pure function mle_timescale(f_abs, ustar, h_vel, ce, f_floor, use_mom_mixrate) result(ts)

FK restratification timescale [s].

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: f_abs
real(kind=wp), intent(in) :: ustar
real(kind=wp), intent(in) :: h_vel
real(kind=wp), intent(in) :: ce
real(kind=wp), intent(in) :: f_floor
logical, intent(in) :: use_mom_mixrate

Return Value real(kind=wp)

private pure function ocean_mle_bytes(this) result(nbytes)

Counted allocatable footprint of the MLE / Fox-Kemper slot (0 when unallocated). One arr_bytes term per array — add a term here when a new allocatable joins the type.

Arguments

Type IntentOptional Attributes Name
class(ocean_mle_t), intent(in) :: this

Return Value integer(kind=int64)


Subroutines

public subroutine mle_compute_transports(grid, metrics, this, ms, epbl, ss, dt_limit, bc)

Fill uhml/vhml (m^3/s) with the FK MLE overturning transport. Run once per outer step at thermo cadence, before the continuity divergence. Steps: (1) b_ml + htot_ml at cell centres (surface→bed band to mld, partial-weight the straddling layer); (2) uDml/vDml at faces from grad b_bar, timescale, H_vel²; (2b) optional per-layer availability cap (a scalar shrink keeping sum_k a(k)=0); (3) fold the mu profile a(k) → uhml/vhml. No-op when enable=.false. or the slot / state arrays are absent.

Arguments

Type IntentOptional Attributes Name
type(hgrid_t), intent(in) :: grid
type(ocean_metrics_t), intent(in) :: metrics
type(ocean_mle_t), intent(inout) :: this
type(multilayer_state_t), intent(in) :: ms
type(ocean_epbl_t), intent(in) :: epbl
type(ocean_surface_stress_t), intent(in), optional :: ss
real(kind=wp), intent(in), optional :: dt_limit

Window (s) the FK transport is integrated over (= dt_therm). When present and > 0, the per-layer availability cap bounds the transport so the windowed drain’s hprev reconstruction stays non-negative in thin surface layers. Absent / ≤ 0 ⇒ no cap. PR-8: this cap is UNCONDITIONAL in production — both production call sites (rdb_ocean_dyn.F90) always pass dt_limit=dyn%therm_dt(dt) > 0, so do_limit is always .true. on the live path. The former &ocean_foxkemper_nml apply_cfl_limit knob was deleted as dead/incoherent — it could only ever be set to a thing the code already always does; do not re-add a knob that toggles this cap.
type(ocean_bc_state_t), intent(in), optional :: bc

Per-edge OBC tags. Masks the FK transport on closed (non-periodic) physical wall faces so no MLE overturning crosses a land boundary. Absent ⇒ array-edge zeroing only.

public pure subroutine mle_fold_x(this, mass_flux_x_layer, nx1, ny, nz)

Add the FK x-face transport into the per-layer zonal mass flux, AFTER continuity_zonal_flux fills it and BEFORE the zonal tracer advect / divergence — so the augmented flux transports both h and tracers (conservative; no velocity touched). No-op when disabled.

Arguments

Type IntentOptional Attributes Name
type(ocean_mle_t), intent(in) :: this
real(kind=wp), intent(inout) :: mass_flux_x_layer(nx1,ny,nz)
integer, intent(in) :: nx1
integer, intent(in) :: ny
integer, intent(in) :: nz

public pure subroutine mle_fold_y(this, mass_flux_y_layer, nx, ny1, nz)

Add the FK y-face transport into the per-layer meridional mass flux. Mirror of mle_fold_x. No-op when disabled.

Arguments

Type IntentOptional Attributes Name
type(ocean_mle_t), intent(in) :: this
real(kind=wp), intent(inout) :: mass_flux_y_layer(nx,ny1,nz)
integer, intent(in) :: nx
integer, intent(in) :: ny1
integer, intent(in) :: nz

public pure subroutine mle_layer_weights(h_face, nz, h_vel, a)

Per-layer transport weights a(k) = mu(sigma_top) - mu(sigma_bot), walking surface (k=nz) -> bed (k=1). sum_k a(k) = mu(0)-mu(-1) = 0 (closed cell => conservation). Layers below the ML base get sigma <= -1 => mu=0 on both interfaces => a(k)=0 (ML-confinement).

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: h_face(nz)
integer, intent(in) :: nz
real(kind=wp), intent(in) :: h_vel
real(kind=wp), intent(out) :: a(nz)

private subroutine ocean_mle_destroy(this)

Arguments

Type IntentOptional Attributes Name
class(ocean_mle_t), intent(inout) :: this

private subroutine ocean_mle_enter_data(this)

Arguments

Type IntentOptional Attributes Name
class(ocean_mle_t), intent(inout) :: this

private subroutine ocean_mle_enter_data_impl(this)

Arguments

Type IntentOptional Attributes Name
type(ocean_mle_t), intent(inout) :: this

private subroutine ocean_mle_exit_data(this)

Arguments

Type IntentOptional Attributes Name
class(ocean_mle_t), intent(inout) :: this

private subroutine ocean_mle_exit_data_impl(this)

Arguments

Type IntentOptional Attributes Name
type(ocean_mle_t), intent(inout) :: this

private subroutine ocean_mle_init(this, grid, nz_ml)

Allocate the 2D ML diagnostics + per-layer face transports. Always allocates (configure runs after init); the off-state footprint is two 2D + two face-shaped 3D arrays.

Arguments

Type IntentOptional Attributes Name
class(ocean_mle_t), intent(inout) :: this
type(hgrid_t), intent(in) :: grid
integer, intent(in), optional :: nz_ml