ocean_kappa_shear_t Derived Type

type, public :: ocean_kappa_shear_t


Inherits

type~~ocean_kappa_shear_t~~InheritsGraph type~ocean_kappa_shear_t ocean_kappa_shear_t type~eos_t eos_t type~ocean_kappa_shear_t->type~eos_t eos

Inherited by

type~~ocean_kappa_shear_t~~InheritedByGraph type~ocean_kappa_shear_t ocean_kappa_shear_t type~ocean_state_t ocean_state_t type~ocean_state_t->type~ocean_kappa_shear_t kshear type~ocean_engine_t ocean_engine_t type~ocean_engine_t->type~ocean_state_t state type~ocean_handle_t ocean_handle_t type~ocean_handle_t->type~ocean_state_t state type~ocean_handle_t->type~ocean_engine_t engine

Components

Type Visibility Attributes Name Initial
logical, public :: at_vertex = .false.

Solve the JHL08 columns at C-grid CORNERS (vorticity points) instead of tracer points, then average the corner Kd back to tracer points (MOM6 VERTEX_SHEAR; the OM5-class production setting). The corner column sees the native face velocities without the u_h/v_h centre average, so resolved shear is not damped before the solve. Default off — bit-identical. Kd: corner solve -> Pass-C average -> kd_int at tracer points. Kv: routed corner->face (MOM6 Kv_shear_Bu consumed in vertvisc) — kd_corner feeds vdiff_apply_momentum’s kv_corner_source seam scaled by prandtl_turb, and the cell-centred kv merge is suppressed (no corner->centre->face smoothing, no double-count). Corner TKE is not carried (tke_int is zeroed in vertex mode).
real(kind=wp), public :: c_n = 0.24_wp

TKE decay vs N (TKE_N_DECAY_CONST).
real(kind=wp), public :: c_s = 0.14_wp

TKE decay vs shear (TKE_SHEAR_DECAY_CONST).
logical, public :: enable = .false.

Master switch. Default off — existing namelists and tests stay bit-identical. Requires vmix%use_closure + thermodynamics (validated at configure).
type(eos_t), public :: eos
real(kind=wp), public, allocatable :: f_centre(:,:)

|f| at cell centres (1/s); filled by set_f_centre.
real(kind=wp), public, allocatable :: f_corner(:,:)

SIGNED Coriolis f at C-grid corners (1/s), (nx+1,ny+1); corner (i,j) is the SW corner of cell (i,j). The kernel squares it (MOM6 vertex form takes f^2 straight at the corner, no 4-point average). Filled by set_f_corner (beta-plane) or fill_coriolis_corner at configure.
real(kind=wp), public :: fri_curvature = -0.97_wp

Ri-function curvature (FRI_CURVATURE).
logical, public :: is_init = .false.

True between init and destroy.
real(kind=wp), public :: kappa_0 = 1.0e-7_wp

Background diffusivity (m^2/s); also the pre-step kappa.
real(kind=wp), public :: kappa_seed = 1.0_wp

Iteration seed diffusivity (m^2/s).
real(kind=wp), public :: kappa_trunc = 1.0e-9_wp

Diffusivity below this -> 0 (m^2/s).
real(kind=wp), public, allocatable :: kd_corner(:,:,:)

Corner diffusivity at interfaces (m^2/s), (nx+1,ny+1,nz+1), global bottom-up. The Pass-B -> Pass-C carrier (MOM6 kappa_vertex): the corner solve writes it, the scatter kernel averages it to tracer points, and the momentum vdiff reads it as the corner Kv source (kv_corner_source, scaled by prandtl_turb — the corner->face viscosity seam, MOM6 Kv_shear_Bu). A REQUIRED snapshot — fusing solve+scatter would be a read-neighbour/write-own do concurrent race. Ring corners (ic=1, ic=nx+1, jc=1, jc=ny+1) are never solved and stay exactly 0.
real(kind=wp), public, allocatable :: kd_int(:,:,:)

Kappa-shear diffusivity at interfaces (m^2/s), (nx,ny,nz+1), global bottom-up: zero at bed (K=1) and surface (K=nz+1).
real(kind=wp), public :: lambda = 0.82_wp

Buoyancy length-scale coefficient (KAPPA_BUOY_SCALE_COEF).
real(kind=wp), public :: lz_rescale = 1.0_wp

Boundary-distance length-scale rescale (LZ_RESCALE).
logical, public :: massless_merge = .false.

is merged onto its massive sub-grid (rdb_massless), solved on nzc <= nz layers, and the kappa/TKE interface fields are interpolated back — replacing the blunt max(h, H_VANISHED) gather floor. Default off (bit-identical). A per-column any(h < H_VANISHED) precheck makes healthy columns bypass the merge machinery entirely, so knob-on stays bit-identical on healthy envelopes (invariant I1).
integer, public :: max_inner_it = 50

Inner Picard iteration cap (MAX_RINO_IT).
integer, public :: max_substep_it = 13

Outer adaptive substep cap (MAX_KAPPA_SHEAR_IT).
real(kind=wp), public :: prandtl_turb = 1.0_wp

Kv = prandtl_turb * Kd into the momentum solve.
real(kind=wp), public :: rho0 = 1035.0_wp

Boussinesq reference density (kg/m^3).
real(kind=wp), public :: ri_crit = 0.25_wp

Critical Richardson number (MOM6 RINO_CRIT).
real(kind=wp), public :: shearmix_rate = 0.089_wp

Source-rate coefficient (SHEARMIX_RATE).
real(kind=wp), public :: src_max_chg = 10.0_wp

Adaptive-dt source-change tolerance band.
real(kind=wp), public :: tke_bg = 0.0_wp

Background TKE (m^2/s^2); Q is a denominator, floored.
real(kind=wp), public, allocatable :: tke_int(:,:,:)

Time-mean TKE at interfaces (m^2/s^2), same shape/convention — diagnostic (currently filled to 0; reserved for the TKE budget diag). Vertex mode zeroes it (corner TKE not carried).
real(kind=wp), public :: tol_err = 0.1_wp

Picard convergence tolerance (KAPPA_SHEAR_TOL_ERR).
real(kind=wp), public :: vel_underflow = 0.0_wp

Velocity snap-to-zero magnitude (m/s) in the projection.
real(kind=wp), public :: vertex_geomean_kdmin = 0.0_wp

Floor (m^2/s) applied to each corner Kd BEFORE the geometric mean (MOM6 VERTEX_SHEAR_GEOMETRIC_MEAN_KDMIN; inert unless vertex_geometric_mean). With 0 the geometric mean hard- zeros Kd at every shear-zone edge; OM5 configs use 1e-9.
logical, public :: vertex_geometric_mean = .false.

Corner->centre averaging: geometric mean of the 4 corner Kd (MOM6 VERTEX_SHEAR_GEOMETRIC_MEAN) instead of the plain arithmetic mean. A geometric mean is 0 if ANY corner is 0 — pair with vertex_geomean_kdmin (see below).

Type-Bound Procedures

procedure, public, non_overridable :: bytes => ocean_kappa_shear_bytes

  • private pure function ocean_kappa_shear_bytes(this) result(nbytes)

    Counted allocatable footprint of the kappa-shear 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_kappa_shear_t), intent(in) :: this

    Return Value integer(kind=int64)

procedure, public, non_overridable :: destroy => ocean_kappa_shear_destroy

procedure, public, non_overridable :: enter_data => ocean_kappa_shear_enter_data

procedure, public, non_overridable :: exit_data => ocean_kappa_shear_exit_data

procedure, public, non_overridable :: init => ocean_kappa_shear_init

  • private subroutine ocean_kappa_shear_init(this, grid, nz_ml)

    Allocate the persistent fields. Always allocates (configure runs after init, so enable is not known yet); the off-cost is the f_centre + two interface fields.

    Arguments

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

procedure, public, non_overridable :: init_vertex => ocean_kappa_shear_init_vertex

  • private subroutine ocean_kappa_shear_init_vertex(this, grid, nz_ml)

    Allocate the vertex-mode corner fields and set at_vertex. Called at CONFIGURE time (after init, before enter_data) — deliberately NOT from init, so the (nx+1,ny+1,nz+1) corner carrier is only ever allocated when the vertex form is actually selected (~327 MB at 1000x800x50).

    Arguments

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

procedure, public, non_overridable :: set_f_centre => ocean_kappa_shear_set_f_centre

  • private subroutine ocean_kappa_shear_set_f_centre(this, grid, f_0, beta, y_ref)

    Fill f_centre with the beta-plane Coriolis magnitude at cell centres: |f_0 + beta*(y - y_ref)|. Mirrors EPBL’s set_f_centre. Call after init, before enter_data.

    Arguments

    Type IntentOptional Attributes Name
    class(ocean_kappa_shear_t), intent(inout) :: this
    type(hgrid_t), intent(in) :: grid
    real(kind=wp), intent(in) :: f_0
    real(kind=wp), intent(in) :: beta
    real(kind=wp), intent(in) :: y_ref

procedure, public, non_overridable :: set_f_corner => ocean_kappa_shear_set_f_corner

  • private subroutine ocean_kappa_shear_set_f_corner(this, grid, f_0, beta, y_ref)

    Fill f_corner with the SIGNED beta-plane Coriolis at C-grid corners: f_0 + beta(y - y_ref), corner row j at y = (j-1-nghost)dy (half a cell below centre row j — corner (i,j) is the SW corner of cell (i,j)). Bit-identical to metrics_fill_coriolis’s beta-plane corner fill. Call after init_vertex, before enter_data.

    Arguments

    Type IntentOptional Attributes Name
    class(ocean_kappa_shear_t), intent(inout) :: this
    type(hgrid_t), intent(in) :: grid
    real(kind=wp), intent(in) :: f_0
    real(kind=wp), intent(in) :: beta
    real(kind=wp), intent(in) :: y_ref

Source Code

   type :: ocean_kappa_shear_t
      logical :: is_init = .false.
         !! True between `init` and `destroy`.

      ! ---- Scheme selection + master switch ----
      logical :: enable = .false.
         !! Master switch.  Default off — existing namelists and tests
         !! stay bit-identical.  Requires `vmix%use_closure` +
         !! thermodynamics (validated at configure).

      ! ---- JHL08 knobs (defaults = paper / OM4 production) ----
      real(wp) :: ri_crit = 0.25_wp
         !! Critical Richardson number (MOM6 RINO_CRIT).
      real(wp) :: shearmix_rate = 0.089_wp
         !! Source-rate coefficient (SHEARMIX_RATE).
      real(wp) :: fri_curvature = -0.97_wp
         !! Ri-function curvature (FRI_CURVATURE).
      real(wp) :: c_n = 0.24_wp
         !! TKE decay vs N (TKE_N_DECAY_CONST).
      real(wp) :: c_s = 0.14_wp
         !! TKE decay vs shear (TKE_SHEAR_DECAY_CONST).
      real(wp) :: lambda = 0.82_wp
         !! Buoyancy length-scale coefficient (KAPPA_BUOY_SCALE_COEF).
      real(wp) :: lz_rescale = 1.0_wp
         !! Boundary-distance length-scale rescale (LZ_RESCALE).
      real(wp) :: kappa_0 = 1.0e-7_wp
         !! Background diffusivity (m^2/s); also the pre-step kappa.
      real(wp) :: kappa_seed = 1.0_wp
         !! Iteration seed diffusivity (m^2/s).
      real(wp) :: kappa_trunc = 1.0e-9_wp
         !! Diffusivity below this -> 0 (m^2/s).
      real(wp) :: tke_bg = 0.0_wp
         !! Background TKE (m^2/s^2); Q is a denominator, floored.
      real(wp) :: tol_err = 0.1_wp
         !! Picard convergence tolerance (KAPPA_SHEAR_TOL_ERR).
      integer :: max_inner_it = 50
         !! Inner Picard iteration cap (MAX_RINO_IT).
      integer :: max_substep_it = 13
         !! Outer adaptive substep cap (MAX_KAPPA_SHEAR_IT).
      real(wp) :: src_max_chg = 10.0_wp
         !! Adaptive-dt source-change tolerance band.
      real(wp) :: prandtl_turb = 1.0_wp
         !! Kv = prandtl_turb * Kd into the momentum solve.
      real(wp) :: vel_underflow = 0.0_wp
         !! Velocity snap-to-zero magnitude (m/s) in the projection.
      logical :: massless_merge = .false.
         !! D4: when on, a column carrying vanished (< H_VANISHED) layers
         !! is merged onto its massive sub-grid (rdb_massless), solved on
         !! `nzc <= nz` layers, and the kappa/TKE interface fields are
         !! interpolated back — replacing the blunt `max(h, H_VANISHED)`
         !! gather floor.  Default off (bit-identical).  A per-column
         !! `any(h < H_VANISHED)` precheck makes healthy columns bypass
         !! the merge machinery entirely, so knob-on stays bit-identical
         !! on healthy envelopes (invariant I1).

      ! ---- Vertex (corner) form — MOM6 VERTEX_SHEAR ----
      logical :: at_vertex = .false.
         !! Solve the JHL08 columns at C-grid CORNERS (vorticity points)
         !! instead of tracer points, then average the corner Kd back to
         !! tracer points (MOM6 `VERTEX_SHEAR`; the OM5-class production
         !! setting).  The corner column sees the native face velocities
         !! without the u_h/v_h centre average, so resolved shear is not
         !! damped before the solve.  Default off — bit-identical.
         !! Kd: corner solve -> Pass-C average -> `kd_int` at tracer
         !! points.  Kv: routed corner->face (MOM6 `Kv_shear_Bu`
         !! consumed in vertvisc) — `kd_corner` feeds
         !! `vdiff_apply_momentum`'s `kv_corner_source` seam scaled by
         !! `prandtl_turb`, and the cell-centred kv merge is suppressed
         !! (no corner->centre->face smoothing, no double-count).
         !! Corner TKE is not carried (`tke_int` is zeroed in vertex
         !! mode).
      logical :: vertex_geometric_mean = .false.
         !! Corner->centre averaging: geometric mean of the 4 corner Kd
         !! (MOM6 `VERTEX_SHEAR_GEOMETRIC_MEAN`) instead of the plain
         !! arithmetic mean.  A geometric mean is 0 if ANY corner is 0 —
         !! pair with `vertex_geomean_kdmin` (see below).
      real(wp) :: vertex_geomean_kdmin = 0.0_wp
         !! Floor (m^2/s) applied to each corner Kd BEFORE the geometric
         !! mean (MOM6 `VERTEX_SHEAR_GEOMETRIC_MEAN_KDMIN`; inert unless
         !! `vertex_geometric_mean`).  With 0 the geometric mean hard-
         !! zeros Kd at every shear-zone edge; OM5 configs use 1e-9.

      ! ---- EOS hookup (shared handle from the eos slot) ----
      ! Value copy of the flat-POD `eos_t` set at configure
      ! from `ocean_state%eos` — the buoyancy derivatives use the
      ! SAME EOS the dyn-core runs.  Maps onto the device with the
      ! parent for free; one source of truth (no drifting copies).
      type(eos_t) :: eos
      real(wp) :: rho0 = 1035.0_wp
         !! Boussinesq reference density (kg/m^3).

      ! ---- Persistent fields ----
      real(wp), allocatable :: f_centre(:, :)
         !! |f| at cell centres (1/s); filled by `set_f_centre`.
      real(wp), allocatable :: kd_int(:, :, :)
         !! Kappa-shear diffusivity at interfaces (m^2/s), (nx,ny,nz+1),
         !! global bottom-up: zero at bed (K=1) and surface (K=nz+1).
      real(wp), allocatable :: tke_int(:, :, :)
         !! Time-mean TKE at interfaces (m^2/s^2), same shape/convention
         !! — diagnostic (currently filled to 0; reserved for the TKE
         !! budget diag).  Vertex mode zeroes it (corner TKE not carried).

      ! ---- Vertex-mode persistent fields (allocated by `init_vertex`,
      ! only when `at_vertex` — the corner carrier is nz+1 full planes) ----
      real(wp), allocatable :: f_corner(:, :)
         !! SIGNED Coriolis f at C-grid corners (1/s), `(nx+1,ny+1)`;
         !! corner (i,j) is the SW corner of cell (i,j).  The kernel
         !! squares it (MOM6 vertex form takes f^2 straight at the
         !! corner, no 4-point average).  Filled by `set_f_corner`
         !! (beta-plane) or `fill_coriolis_corner` at configure.
      real(wp), allocatable :: kd_corner(:, :, :)
         !! Corner diffusivity at interfaces (m^2/s), `(nx+1,ny+1,nz+1)`,
         !! global bottom-up.  The Pass-B -> Pass-C carrier (MOM6
         !! `kappa_vertex`): the corner solve writes it, the scatter
         !! kernel averages it to tracer points, and the momentum vdiff
         !! reads it as the corner Kv source (`kv_corner_source`,
         !! scaled by `prandtl_turb` — the corner->face viscosity seam,
         !! MOM6 `Kv_shear_Bu`).  A REQUIRED snapshot —
         !! fusing solve+scatter would be a read-neighbour/write-own
         !! `do concurrent` race.  Ring corners (ic=1, ic=nx+1, jc=1,
         !! jc=ny+1) are never solved and stay exactly 0.
   contains
      procedure, non_overridable :: init => ocean_kappa_shear_init
      procedure, non_overridable :: init_vertex => ocean_kappa_shear_init_vertex
      procedure, non_overridable :: destroy => ocean_kappa_shear_destroy
      procedure, non_overridable :: enter_data => ocean_kappa_shear_enter_data
      procedure, non_overridable :: exit_data => ocean_kappa_shear_exit_data
      procedure, non_overridable :: set_f_centre => ocean_kappa_shear_set_f_centre
      procedure, non_overridable :: set_f_corner => ocean_kappa_shear_set_f_corner
      procedure, non_overridable :: bytes => ocean_kappa_shear_bytes
   end type ocean_kappa_shear_t