rdb_ocean_vcoord Module

Holds the vertical-coordinate configuration and the per-step target grid that the ALE remap step relamps multilayer.h_layer + every multilayer.tracers(t)%hTr onto. The coastal path has VCOORD_SIGMA, VCOORD_ZSIGMA, VCOORD_ZSTAR, VCOORD_ZSTAR_SIGMA, VCOORD_ZSTAR_FULL plus rdb_remap_column; this module is the C-grid counterpart for the ocean dynamical core.

Phase 5g build-out status (this commit, Layer 2):

  • target_h(nx, ny, nz_ml) is allocated on init + bound to the device via enter_data / exit_data.
  • compute_target_h(this, total_h, eta) populates target_h per the current coord_type. Working bodies: VCOORD_EULERIAN_Z — H · dsig(k) (η ignored) VCOORD_SIGMA — (H + η) · dsig(k) VCOORD_ZSTAR — MOM6 z: the z_fixed nominal profile dilated by (H + η)/H over bed fillers VCOORD_ZSIGMA — smoothstep blend(sigma, fixed z-levels) VCOORD_ZSTAR_SIGMA — smoothstep blend(sigma, z-lite) VCOORD_ZSTAR_FULL — per-column z_ref + vanishing-layer floors
  • build_zref_full(this, h_bed) populates z_ref(:, :, 0:nz_ml) per column from local bathymetry. Call once at init (or any time the bathymetry changes); the per-step compute_target_h then walks the cached z_ref table.
  • Driver call site (remap kernel invocation between slow stages in ocean_dyn_step_split) is NOT wired yet — simulation state is untouched. External callers can invoke compute_target_h / build_zref_full to inspect the target grid without advancing dynamics.

Collaborator hand-off — remaining Phase 5g work (Layer 3): 1. Driver wiring in ocean_dyn_step_split: call compute_target_h(total_h_2d, dyn%bt_work%bt_eta) then invoke the remap kernel on multilayer.h_layer + each tracers(t)%hTr + face velocities. Recompute bt_eta from the new sum of h_layer - bt_H_ref. 2. Face-velocity remap adapter — MOM6 reconstructs u at centres, remaps, then projects back to faces with a divergence-free correction. Donor-cell on (h·u)_face is the simpler fallback. 3. VANISHING_LAYER_TOL-gated CWC for VCOORD_ZSTAR_FULL, lifted from the coastal multilayer kernels.

See docs/ROADMAP_OCEAN.md Phase 5g for the full scope.


Uses

  • module~~rdb_ocean_vcoord~~UsesGraph module~rdb_ocean_vcoord rdb_ocean_vcoord iso_fortran_env iso_fortran_env module~rdb_ocean_vcoord->iso_fortran_env module~rdb_constants rdb_constants module~rdb_ocean_vcoord->module~rdb_constants module~rdb_eos rdb_eos module~rdb_ocean_vcoord->module~rdb_eos module~rdb_grid rdb_grid module~rdb_ocean_vcoord->module~rdb_grid module~rdb_mem_report rdb_mem_report module~rdb_ocean_vcoord->module~rdb_mem_report module~rdb_vcoord rdb_vcoord module~rdb_ocean_vcoord->module~rdb_vcoord pic_types pic_types module~rdb_constants->pic_types module~rdb_eos->module~rdb_constants module~rdb_eos->module~rdb_grid 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_vcoord->module~rdb_constants module~rdb_vcoord->pic_logger module~rdb_vcoord->pic_strings

Used by

  • module~~rdb_ocean_vcoord~~UsedByGraph module~rdb_ocean_vcoord rdb_ocean_vcoord module~rdb_ocean_diag_fills rdb_ocean_diag_fills module~rdb_ocean_diag_fills->module~rdb_ocean_vcoord module~rdb_ocean_state rdb_ocean_state module~rdb_ocean_diag_fills->module~rdb_ocean_state module~rdb_ocean_dyn rdb_ocean_dyn module~rdb_ocean_dyn->module~rdb_ocean_vcoord module~rdb_ocean_remap rdb_ocean_remap module~rdb_ocean_dyn->module~rdb_ocean_remap module~rdb_ocean_engine rdb_ocean_engine module~rdb_ocean_engine->module~rdb_ocean_vcoord module~rdb_ocean_engine->module~rdb_ocean_diag_fills module~rdb_ocean_engine->module~rdb_ocean_dyn module~rdb_ocean_setup rdb_ocean_setup module~rdb_ocean_engine->module~rdb_ocean_setup module~rdb_ocean_engine->module~rdb_ocean_state module~rdb_ocean_diag_derived rdb_ocean_diag_derived module~rdb_ocean_engine->module~rdb_ocean_diag_derived module~rdb_ocean_remap->module~rdb_ocean_vcoord module~rdb_ocean_setup->module~rdb_ocean_vcoord module~rdb_ocean_setup->module~rdb_ocean_dyn module~rdb_ocean_setup->module~rdb_ocean_state module~rdb_ocean_state->module~rdb_ocean_vcoord 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_engine module~rdb_driver->module~rdb_ocean_state module~rdb_handle rdb_handle module~rdb_handle->module~rdb_ocean_engine module~rdb_handle->module~rdb_ocean_state module~rdb_ocean_api rdb_ocean_api module~rdb_ocean_api->module~rdb_ocean_diag_fills module~rdb_ocean_api->module~rdb_ocean_dyn module~rdb_ocean_api->module~rdb_ocean_engine module~rdb_ocean_api->module~rdb_handle module~rdb_ocean_api->module~rdb_ocean_diag_derived module~rdb_ocean_diag_derived->module~rdb_ocean_diag_fills module~rdb_ocean_diag_derived->module~rdb_ocean_state

Variables

Type Visibility Attributes Name Initial
integer, private, parameter :: HYCOM_FLOOR_PROFILE = 2

Floor at Σ z_fixed_dz(k)·(H+η)/H — a stretched z* resolution in metres (z_fixed_profile = "list" | "tanh").
integer, private, parameter :: HYCOM_FLOOR_SIGMA = 0

Floor at Σ dsig·(H+η) — a column FRACTION, i.e. a sigma floor. Only the fallback for a slot no setup path configured (z_fixed_h_ref <= 0); it was the only floor before 2026-10-02, when it set 95 % of the 1-degree Southern Ocean’s interfaces and made hycom terrain-following there (audit finding H1).
integer, private, parameter :: HYCOM_FLOOR_UNIFORM = 1

Floor at Σ (z_fixed_h_ref/nz)·(H+η)/H — uniform z* in metres (&vcoord_nml z_fixed_profile = "uniform", max_depth/nz).
integer, private, parameter :: NR_ITERS = 8

Fixed (GPU-uniform) Newton iteration budget for the density→depth inversion. Unrolled, no data-dependent while.
real(kind=wp), private, parameter :: NR_OFFSET = 1.0e-6_wp

Out-of-range nudge applied only when the boundary gradient ≈ 0.
real(kind=wp), private, parameter :: NR_TOL = 1.0e-12_wp

Newton convergence tolerance — tested on |delta| AFTER xi += delta.
real(kind=wp), private, parameter :: Z_FIXED_BED_PARTIAL_MIN = H_VANISHED

Minimum LIVE partial BOTTOM cell thickness for VCOORD_Z_FIXED — the bed mirror of Z_FIXED_TOP_PARTIAL_FRAC, expressed as an absolute thickness rather than a fraction of the spacing.

What it closes. The bed branch used to floor its partial cell at h_min (zstar_h_min) alone, and the contract for this family is zstar_h_min <= H_VANISHED. That leaves the half-open band (zstar_h_min, H_VANISHED] = (1e-4, 1.5e-4] with the default h_min, in which a layer is live to the coordinate and vanished to every consumer — the EOS substitutes reference T/S, the remap reads its concentration as zero, k_top-style scans walk past it. A millimetre of η is enough to flip a bed remainder across it, and on validation_examples/ocean/isomip_plus/ocean0_idealised_zfixed.nml five columns did so 165 times in five days. With this floor the bed partial cell is either a live layer strictly ABOVE the marker or an inert filler at exactly h_min; nothing lands in between.

Why not Z_FIXED_TOP_PARTIAL_FRAC*h_nominal at the bed too. That is the full mirror, and it is a different, larger change: with a 20 m spacing it would merge away every bed cell thinner than 2 m, moving the level structure — and the effective bathymetry — of every shipped z_fixed configuration. The thin-layer CFL argument that motivates the top’s 10 % applies at the bed as well and that mirror may still be worth taking, but it is a bathymetry change and belongs in its own slice with its own validation. This constant fixes the defect that is actually diagnosed: a 0.05 mm window, so in practice only a column whose bed remainder lands inside it moves at all.
real(kind=wp), private, parameter :: Z_FIXED_TOP_PARTIAL_FRAC = 0.1_wp

Minimum PARTIAL TOP CELL thickness, as a fraction of the nominal layer spacing h_nominal, for VCOORD_Z_FIXED under a rigid top (z_top > 0, i.e. an ice-shelf cavity). The ice base cuts the first live layer wherever the draft crosses a nominal level; left unguarded that cut can leave a sliver of arbitrarily small thickness (0 < h << h_nominal) sitting against the ice. A layer that thin is live (it clears H_VANISHED) but violates every thin-layer CFL the column has, so when the cut would leave less than Z_FIXED_TOP_PARTIAL_FRAC*h_nominal the sliver is merged into the layer BELOW and the vacated index becomes an inert filler — exactly mirroring how the BED side treats its own slivers today (a would-be sub-h_min bed cell collapses to h_min and hands its water to the layer above).

The two thresholds differ on purpose and the asymmetry is the honest one: the bed’s threshold is Z_FIXED_BED_PARTIAL_MIN (below), the top’s is a fraction of the spacing. 0.1 is MITgcm’s hFacMin default, which is the minimum partial-cell fraction Losch (2008, JGR 113 C08043, §2.1) used for exactly this ice-shelf partial-top-cell problem. Not a namelist knob: it is inert unless z_top > 0, and a cavity on a z-like coordinate is itself a new, fenced configuration.

Derived Types

type, public ::  ocean_vcoord_t

Components

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

&vcoord_nml check_vanished_content — the I1′ tripwire. Carried on this slot (rather than on ocean_dyn_t) because the vertical coordinate is what MAKES vanished layers, so the knob that polices them belongs beside zstar_h_min and the filler contract. A plain scalar on the type: it rides the existing copyin(this) and adds no device array. Read by check_vanished_invariant_or_die in rdb_ocean_dyn. Default .false. ⇒ no scan, no cost.
integer, public :: coord_type = VCOORD_EULERIAN_Z

Selected vertical-coordinate variant.
real(kind=wp), public, allocatable :: dsig(:)

Per-layer σ-fraction. Sums to 1.0; size nz_ml.
logical, public :: is_init = .false.

True between init and destroy. Prefer this to allocated(...) — tracks GPU device attachment too.
integer, public :: nx_total = 0

Total i-extent of target_h (incl. halos).
integer, public :: ny_total = 0

Total j-extent of target_h (incl. halos).
integer, public :: nz_ml = 0

Number of active layers.
real(kind=wp), public :: regrid_time_scale = 0.0_wp

Grid time-filter timescale τ (s) for the ALE regrid. After compute_target_h builds the new target grid, the remap step relaxes the coordinate a fraction dt/(τ+dt) toward that target each outer step rather than jumping to it — damping the per-step grid-motion shock that drives the σ/z* PGE (White & Adcroft 2008, the grid time-filter). Scalar on the type, reaches the device through the existing copyin(this); no new device array. Default 0.0 ⇒ wtd = 1 ⇒ jump to target ⇒ bit-identical to the no-filter remap.
logical, public :: remap_boundary_extrap = .false.

Close the ALE remap’s reconstruction at the two boundary cells (k=1, k=nz) with the linear-exact one-sided edge pair instead of the PCM flatten (MOM6 BOUNDARY_EXTRAPOLATION).

Read more…
logical, public :: remap_check_preconditions = .false.

Assert the ALE remap’s column preconditions once per remap and fail loud on a violation (audit findings V5, V6).

Read more…
real(kind=wp), public, allocatable :: remap_conc_s(:,:,:)

Layer-mean S concentration scratch for VCOORD_RHO, shape (nx, ny, nz).
real(kind=wp), public, allocatable :: remap_conc_t(:,:,:)

Layer-mean T concentration scratch for the VCOORD_RHO density inversion, shape (nx, ny, nz). Built from hTr / remap_h_old (vanishing-layer-guarded) once per remap.
real(kind=wp), public, allocatable :: remap_h_old(:,:,:)

Snapshot of h_layer before the remap, shape (nx, ny, nz).
real(kind=wp), public, allocatable :: remap_h_ref(:,:)

H reference (total_h − bt_eta) scratch, shape (nx, ny).
integer, public :: remap_method = REMAP_PPM

ALE remap reconstruction order (REMAP_PCM/PLM/PPM/PPM_H4/PQM). PQM falls back to PPM for nz < 5 (see remap_column_pqm).
logical, public :: remap_nonuniform_weights = .false.

Use the non-uniform-grid reconstruction weights in the ALE remap’s PLM slope and PPM edge estimate — Colella & Woodward (1984) eqs (1.6)-(1.8) — instead of their equal-thickness specialisations (0.5·minmod and (7/12, -1/12)).

Read more…
real(kind=wp), public, allocatable :: remap_total_h(:,:)

Column-total h_layer scratch, shape (nx, ny).
logical, public :: remap_vel_conserve_ke = .false.

Enable the KE-conserving rescale of the remapped layer velocities. After the per-face column remap (which already conserves u·h, i.e. momentum), rescale the BAROCLINIC velocity anomaly per column so column KE Σ ½ h·u² is preserved (Adcroft & Hallberg 2006 layer-velocity remap), capped at a 1.25× rescale factor. The barotropic/depth-mean component is never touched (mode-split consistency). Default .false. ⇒ velocities unchanged ⇒ bit-identical.
real(kind=wp), public :: rho_ref_pressure = 2.0e7_wp

Reference pressure (Pa, default 2e7 = 2000 dbar) for the potential density that defines the VCOORD_RHO coordinate. A rdb convention (not MOM6-inherited).
real(kind=wp), public, allocatable :: rho_target(:)

Target interface potential densities (kg/m³), shape 0:nz_ml.
real(kind=wp), public, allocatable :: target_h(:,:,:)

Target layer thickness (m), shape (nx, ny, nz_ml).
real(kind=wp), public, allocatable :: z_fixed_dz(:)

VCOORD_Z_FIXED nominal layer thicknesses (m), shape nz_ml, bottom-up: z_fixed_zi(k-1) - z_fixed_zi(k), stored separately so the partial-top-cell threshold uses the exact namelist value. Read only when z_fixed_use_profile.
real(kind=wp), public :: z_fixed_h_ref = 0.0_wp

Total reference depth (m) for VCOORD_Z_FIXED and VCOORD_ZSTAR (the MOM6 z* nominal profile is the z_fixed one, dilated per column — ocean_vcoord_zstar_target). Layer interfaces sit at z = k · h_ref / nz_ml from the surface, same as MOM6’s COORD_CONFIG = "gprime" with MAXIMUM_DEPTH = h_ref. Driver writes from cfg%ocean%topo%max_depth at init. When 0 (default) the compute_target_h Z_FIXED branch falls back to a uniform H · dsig(k) target so the path stays sane in tests that don’t explicitly set this knob. Under a stretched profile (z_fixed_use_profile) it is the profile’s total depth, z_fixed_zi(0), and the nominal interfaces come from z_fixed_zi instead of h_ref/nz.
logical, public :: z_fixed_use_profile = .false.

&vcoord_nml z_fixed_profile /= "uniform": the VCOORD_Z_FIXED (and VCOORD_ZSTAR) nominal interfaces come from z_fixed_zi / z_fixed_dz (set by ocean_vcoord_set_z_fixed_profile) rather than from the uniform z_fixed_h_ref/nz_ml. Scalar, rides copyin(this). Default .false. ⇒ the uniform arithmetic, byte-identical.
real(kind=wp), public, allocatable :: z_fixed_zi(:)

VCOORD_Z_FIXED nominal interface depths (m, positive down, below z = 0), shape 0:nz_ml, BOTTOM-UP like the state: z_fixed_zi(k) is the TOP interface of layer k, so z_fixed_zi(nz_ml) = 0 (the surface) and z_fixed_zi(0) is the profile’s total depth. Allocated at init (zeros), so it is never a placeholder; read only when z_fixed_use_profile.
real(kind=wp), public, allocatable :: z_ref(:,:,:)

Per-column z* reference (m), shape (nx, ny, 0:nz_ml).
real(kind=wp), public, allocatable :: z_ref_global(:)

Global reference z-interfaces (m, positive-down), shape 0:nz_ml. Used by ZSIGMA / ZSTAR_SIGMA (NOT by ZSTAR, whose nominal profile is the z_fixed one, z_fixed_zi).
real(kind=wp), public, allocatable :: z_top(:,:)

Geopotential depth of the column top (m, positive down), shape (nx, ny). 0 = the free surface at z = 0.
logical, public :: zfixed_closed_faces = .false.

&vcoord_nml zfixed_closed_faces — partial-step z-level face closure. Meaningful on the three GEOMETRIC families that vanish bed-side layers, VCOORD_Z_FIXED, VCOORD_ZSTAR (MOM6 z*) and VCOORD_ZSTAR_FULL, where a layer whose reference range lies inside the bed (or, under z_fixed, the ice draft) is an inert FILLER; a velocity face at which that layer is a filler on EITHER side is a z-level WALL, not a thin passage (Adcroft, Hill & Marshall 1997; Losch 2008).

Read more…
real(kind=wp), public :: zsigma_blend_width = 100.0_wp

Smoothstep blend width (m) above the transition depth.
real(kind=wp), public :: zsigma_depth_transition = 200.0_wp

Sigma → z* transition depth (m) for VCOORD_ZSTAR_SIGMA.
real(kind=wp), public :: zstar_h_min = 1.0e-4_wp

Bed-side vanishing-layer floor (m). Two contracts, picked by the coordinate family — see rdb_vcoord :: vcoord_h_min_role.

Read more…
real(kind=wp), public :: zstar_h_surf_target = 5.0_wp

Surface-layer thickness anchor for VCOORD_ZSTAR_FULL (m).
integer, public :: zstar_n_surf = 0

Number of fine near-surface layers for ZSTAR_FULL. ≤ 0 = auto-pick (max(1, nz_ml/3)).
integer, public :: zstar_stretching = STRETCH_UNIFORM

Stretching mode for ZSTAR_FULL. STRETCH_UNIFORM (default) or STRETCH_LOG for a geometric near-surface fine zone.

Type-Bound Procedures

procedure, public, non_overridable :: build_zref_full => ocean_vcoord_build_zref_full
procedure, public, non_overridable :: bytes => ocean_vcoord_bytes
procedure, public, non_overridable :: compute_target_h => ocean_vcoord_compute_target_h
procedure, public, non_overridable :: compute_target_h_rho => ocean_vcoord_compute_target_h_rho
procedure, public, non_overridable :: destroy => ocean_vcoord_destroy
procedure, public, non_overridable :: enter_data => ocean_vcoord_enter_data
procedure, public, non_overridable :: exit_data => ocean_vcoord_exit_data
procedure, public, non_overridable :: init => ocean_vcoord_init

Functions

public pure function ocean_vcoord_count_bed_steps(target_h, total_h, dy_cu, dx_cv, nx, ny, nz, i0, i1, j0, j1, h_vanished) result(n_step)

Count the wet velocity faces at which the two columns’ BED falls in different nominal z_fixed layers — a bed staircase step that crosses a nominal interface, so that an OPEN face pairs a live layer on one side with a bed FILLER on the other. That is the face &vcoord_nml zfixed_closed_faces closes; left open, the FV pressure gradient across the step drives the flow from rest (configure_ocean_closed_faces refuses the configuration).

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: target_h(nx,ny,nz)

The z_fixed target thickness at η = 0.
real(kind=wp), intent(in) :: total_h(nx,ny)

Column reference thickness (bt_H_ref, ghost-filled).
real(kind=wp), intent(in) :: dy_cu(nx+1,ny)

Land-masked u-face width.
real(kind=wp), intent(in) :: dx_cv(nx,ny+1)

Land-masked v-face width.
integer, intent(in) :: nx

i-extent of the CENTRE arrays (total, incl. halos).
integer, intent(in) :: ny

j-extent of the CENTRE arrays (total, incl. halos).
integer, intent(in) :: nz

Number of layers; k = 1 is the bed, k = nz the top.
integer, intent(in) :: i0

Owned cell range (nghost+1 : nghost+n_phys).
integer, intent(in) :: i1

Owned cell range (nghost+1 : nghost+n_phys).
integer, intent(in) :: j0

Owned cell range (nghost+1 : nghost+n_phys).
integer, intent(in) :: j1

Owned cell range (nghost+1 : nghost+n_phys).
real(kind=wp), intent(in) :: h_vanished

Inert-filler marker (H_VANISHED).

Return Value integer

public pure function ocean_vcoord_count_ledges(open_u, open_v, target_h, nx, ny, nz, h_vanished) result(n_ledge)

Count LEDGE cells: a cell that is LIVE at layer k but all four of whose own-layer faces are closed, i.e. water the mask has isolated. A ledge needs a one-cell-wide spike in the bed or the draft; it is inert by construction (no flux in or out, and its velocity is zeroed every stage), but a non-zero count is worth saying out loud once at configure, because it means the mask is walling off real water.

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: open_u(nx+1,ny,nz)
real(kind=wp), intent(in) :: open_v(nx,ny+1,nz)
real(kind=wp), intent(in) :: target_h(nx,ny,nz)
integer, intent(in) :: nx
integer, intent(in) :: ny
integer, intent(in) :: nz
real(kind=wp), intent(in) :: h_vanished

Return Value integer

public pure function parse_ocean_vcoord_type(name) result(code)

Ocean-path wrapper around the canonical parse_vcoord_type in rdb_vcoord. Pins the unrecognised-string fallback to VCOORD_EULERIAN_Z — the ocean path’s “do nothing” default, distinct from the coastal path’s VCOORD_SIGMA fallback. Kept as a thin name-preserving wrapper so the ocean-only semantic (fallback choice) is visible at the call site.

Arguments

Type IntentOptional Attributes Name
character(len=*), intent(in) :: name

Return Value integer

private pure function ocean_vcoord_bytes(this) result(nbytes)

Counted allocatable footprint of the vertical coordinate 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_vcoord_t), intent(in) :: this

Return Value integer(kind=int64)


Subroutines

public pure subroutine invert_density_targets(nk, hc, rhoc, n_int, rho_tgt, z_new)

Density-space interface inversion — the single source of truth for the RHO vcoord regrid (ocean_vcoord_compute_target_h_rho_impl) AND the DENSITY diagnostic remap (rdb_ocean_diag_fills).

Read more…

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: nk
real(kind=wp), intent(in) :: hc(nk)
real(kind=wp), intent(in) :: rhoc(nk)
integer, intent(in) :: n_int
real(kind=wp), intent(in) :: rho_tgt(n_int)
real(kind=wp), intent(out) :: z_new(n_int+2)

public pure subroutine ocean_vcoord_closed_face_masks(open_u, open_v, target_h, nx, ny, nz, h_vanished)

Partial-step z-level FACE CLOSURE mask for VCOORD_Z_FIXED (&vcoord_nml zfixed_closed_faces; Adcroft, Hill & Marshall 1997; Losch 2008 §2.1 for the ice-shelf cavity).

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(out) :: open_u(nx+1,ny,nz)

u-face 0/1 open mask.
real(kind=wp), intent(out) :: open_v(nx,ny+1,nz)

v-face 0/1 open mask.
real(kind=wp), intent(in) :: target_h(nx,ny,nz)

The z_fixed target thickness at η = 0.
integer, intent(in) :: nx

i-extent of the CENTRE arrays (total, incl. halos).
integer, intent(in) :: ny

j-extent of the CENTRE arrays (total, incl. halos).
integer, intent(in) :: nz

Number of layers; k = 1 is the bed, k = nz the top.
real(kind=wp), intent(in) :: h_vanished

Inert-filler marker (H_VANISHED). A layer is LIVE iff its target thickness is strictly greater than this.

public pure subroutine ocean_vcoord_eta0_target(vc, target_h, total_h, nx, ny, nz)

The target layer thickness at eta = 0 of a GEOMETRIC family that vanishes layers — the ONE definition of “live” that the partial-step face mask (configure_ocean_closed_faces) and the on-target initial seed share with the running ALE regrid.

Read more…

Arguments

Type IntentOptional Attributes Name
type(ocean_vcoord_t), intent(in) :: vc

The vertical-coordinate slot (coord_type + every table).
real(kind=wp), intent(inout) :: target_h(nx,ny,nz)

Target thickness at eta = 0 (m), bottom-up.
real(kind=wp), intent(in) :: total_h(nx,ny)

Reference column thickness (m) — bt_H_ref at configure, the seed’s water column at IC time.
integer, intent(in) :: nx

i-extent (total, incl. halos).
integer, intent(in) :: ny

j-extent (total, incl. halos).
integer, intent(in) :: nz

Number of layers; k = 1 is the bed.

public pure subroutine ocean_vcoord_k_bot_from_target(k_bot, k_bot_u, k_bot_v, target_h, nx, ny, nz, h_vanished)

The shared FIRST-LIVE-LAYER index counting UP from the bed — multilayer_state_t%k_bot and its two face twins. The bed-side mirror of ocean_vcoord_k_top_from_target, built from the same layer-thickness field by the same strict > h_vanished test.

Read more…

Arguments

Type IntentOptional Attributes Name
integer, intent(out) :: k_bot(nx,ny)

Cell-centred first live layer counting up from the bed.
integer, intent(out) :: k_bot_u(nx+1,ny)

u-face twin.
integer, intent(out) :: k_bot_v(nx,ny+1)

v-face twin.
real(kind=wp), intent(in) :: target_h(nx,ny,nz)

Layer thickness (m) the live/filler pattern is read from — the z_fixed target at eta = 0 at configure time.
integer, intent(in) :: nx

i-extent of the CENTRE arrays (total, incl. halos).
integer, intent(in) :: ny

j-extent of the CENTRE arrays (total, incl. halos).
integer, intent(in) :: nz

Number of layers; k = 1 is the bed, k = nz the top.
real(kind=wp), intent(in) :: h_vanished

Inert-filler marker (H_VANISHED). LIVE iff strictly greater.

public pure subroutine ocean_vcoord_k_top_from_target(k_top, k_top_u, k_top_v, target_h, nx, ny, nz, h_vanished)

The shared FIRST-LIVE-LAYER index, counting down from the top — multilayer_state_t%k_top and its two face twins — built from a layer-thickness field.

Read more…

Arguments

Type IntentOptional Attributes Name
integer, intent(out) :: k_top(nx,ny)

Cell-centred first live layer.
integer, intent(out) :: k_top_u(nx+1,ny)

u-face twin.
integer, intent(out) :: k_top_v(nx,ny+1)

v-face twin.
real(kind=wp), intent(in) :: target_h(nx,ny,nz)

Layer thickness (m) the live/filler pattern is read from — the z_fixed target at eta = 0 at configure time.
integer, intent(in) :: nx

i-extent of the CENTRE arrays (total, incl. halos).
integer, intent(in) :: ny

j-extent of the CENTRE arrays (total, incl. halos).
integer, intent(in) :: nz

Number of layers; k = 1 is the bed, k = nz the top.
real(kind=wp), intent(in) :: h_vanished

Inert-filler marker (H_VANISHED). A layer is LIVE iff its thickness is strictly greater than this.

public pure subroutine ocean_vcoord_set_z_fixed_profile(this, dz_surface_first)

Install a stretched VCOORD_Z_FIXED nominal profile: flip the surface-first thicknesses into the bottom-up z_fixed_dz, build the interface table z_fixed_zi by accumulating from the surface (z_fixed_zi(nz) = 0 exactly), set z_fixed_h_ref to the total and raise z_fixed_use_profile. Setup-time host code; must run BEFORE enter_data (the copyin captures the tables).

Arguments

Type IntentOptional Attributes Name
type(ocean_vcoord_t), intent(inout) :: this
real(kind=wp), intent(in) :: dz_surface_first(:)

Nominal thicknesses (m), dz_surface_first(1) = top layer; size must be nz_ml.

public pure subroutine ocean_vcoord_z_fixed_target(target_h, total_h, eta, z_top, nx, ny, nz, h_nominal, use_profile, zi, dz_nom, h_min)

VCOORD_Z_FIXED target grid — quasi-geopotential interfaces under a rigid top, with inert fillers and a partial cell at BOTH ends (Yung, Hallberg, Adcroft & Morrison 2026, JAMES, Fig. 1b: quasi-z layers are geopotential and VANISH where they outcrop into the ice base; Asay-Davis et al. 2016 §3.1.5: z-level models use both partial top and bottom cells).

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(out) :: target_h(nx,ny,nz)

Target layer thickness (m).
real(kind=wp), intent(in) :: total_h(nx,ny)

Column reference thickness H (m) — Σ h_layer − η.
real(kind=wp), intent(in) :: eta(nx,ny)

Free-surface anomaly η (m).
real(kind=wp), intent(in) :: z_top(nx,ny)

Geopotential depth of the column top (m, positive down, >= 0). 0 ⇒ the pre-cavity arithmetic, bit-for-bit.
integer, intent(in) :: nx

i-extent of every array (total, incl. halos).
integer, intent(in) :: ny

j-extent of every array (total, incl. halos).
integer, intent(in) :: nz

Number of layers; k = 1 is the bed, k = nz the top.
real(kind=wp), intent(in) :: h_nominal

Nominal layer spacing z_fixed_h_ref/nz (m), > 0. Unused when use_profile.
logical, intent(in) :: use_profile

Take the nominal interfaces from zi / dz_nom.
real(kind=wp), intent(in) :: zi(0:nz)

Nominal interface depths (m), bottom-up, zi(k) = top of layer k, zi(nz) = 0. Read only when use_profile.
real(kind=wp), intent(in) :: dz_nom(nz)

Nominal layer thicknesses (m), bottom-up. Read only when use_profile.
real(kind=wp), intent(in) :: h_min

Inert-filler thickness (zstar_h_min, <= H_VANISHED).

public pure subroutine ocean_vcoord_z_fixed_target_uniform(target_h, total_h, eta, z_top, nx, ny, nz, h_nominal, h_min)

ocean_vcoord_z_fixed_target on the UNIFORM nominal spacing h_nominal — the historical signature, for callers (tests, setup code without a vcoord slot) that have no profile tables. Same kernel, use_profile = .false., so the arithmetic is the uniform branch’s exactly; the two tables are never read.

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(out) :: target_h(nx,ny,nz)

Target layer thickness (m).
real(kind=wp), intent(in) :: total_h(nx,ny)

Column reference thickness H (m).
real(kind=wp), intent(in) :: eta(nx,ny)

Free-surface anomaly η (m).
real(kind=wp), intent(in) :: z_top(nx,ny)

Geopotential depth of the column top (m, positive down).
integer, intent(in) :: nx

i-extent of every array (total, incl. halos).
integer, intent(in) :: ny

j-extent of every array (total, incl. halos).
integer, intent(in) :: nz

Number of layers; k = 1 is the bed, k = nz the top.
real(kind=wp), intent(in) :: h_nominal

Nominal layer spacing (m), > 0.
real(kind=wp), intent(in) :: h_min

Inert-filler thickness.

public pure subroutine ocean_vcoord_zstar_target(target_h, total_h, eta, nx, ny, nz, h_nominal, use_profile, zi, h_min)

VCOORD_ZSTAR target grid — MOM6 z* (REGRIDDING_COORDINATE_MODE = "Z*", build_zstar_column, MOM6 src/ALE/coord_zlike.F90 lines 65-146): the FIXED nominal z profile of z_fixed (&vcoord_nml z_fixed_profile uniform / list / tanh — the same z_fixed_zi table, not a parallel one) DILATED per column by the free-surface stretching, over a partial bed cell and inert bed fillers.

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(out) :: target_h(nx,ny,nz)

Target layer thickness (m), bottom-up.
real(kind=wp), intent(in) :: total_h(nx,ny)

Column reference thickness H (m) — Σ h_layer − η.
real(kind=wp), intent(in) :: eta(nx,ny)

Free-surface anomaly η (m).
integer, intent(in) :: nx

i-extent of every array (total, incl. halos).
integer, intent(in) :: ny

j-extent of every array (total, incl. halos).
integer, intent(in) :: nz

Number of layers; k = 1 is the bed, k = nz the surface.
real(kind=wp), intent(in) :: h_nominal

Uniform nominal spacing z_fixed_h_ref/nz (m). Unused when use_profile.
logical, intent(in) :: use_profile

Take the nominal interfaces from zi.
real(kind=wp), intent(in) :: zi(0:nz)

Nominal interface depths (m), bottom-up, zi(k) = top of layer k, zi(nz) = 0. Read only when use_profile.
real(kind=wp), intent(in) :: h_min

Inert-filler thickness (zstar_h_min, <= H_VANISHED).

private pure subroutine ocean_vcoord_build_zref_full(this, h_bed)

Populate z_ref(i, j, 0:nz_ml) per column from the local bathymetry h_bed(i, j). Mirrors zstar_full_build_column from src/ALE/rdb_vcoord.F90 but as a 2D loop owned by this slot — keeps the coastal helper untouched while letting the ocean path own its z_ref lifecycle.

Read more…

Arguments

Type IntentOptional Attributes Name
class(ocean_vcoord_t), intent(inout) :: this
real(kind=wp), intent(in) :: h_bed(:,:)

Bed depth at cell centres (m, positive-down).

private pure subroutine ocean_vcoord_compute_target_h(this, total_h, eta)

Thin polymorphic wrapper. A type-bound procedure’s passed object must be class(...), but mapping a polymorphic list item into a target / offload region is unspecified behaviour (gfortran -Wopenmp; see FORTRAN_STYLE.md) — so the do concurrent kernels live in the type(ocean_vcoord_t) _impl and this wrapper only resolves the concrete type. ocean_vcoord_t is never extended, so the dynamic type is always the declared type. Mirrors the enter_data/exit_data split.

Arguments

Type IntentOptional Attributes Name
class(ocean_vcoord_t), intent(inout) :: this
real(kind=wp), intent(in) :: total_h(:,:)
real(kind=wp), intent(in) :: eta(:,:)

private pure subroutine ocean_vcoord_compute_target_h_impl(this, total_h, eta)

Populate target_h(i,j,k) from the column-total depth (H, constant per column for the ocean path; coastal would pass the bathymetry) and the free-surface anomaly η.

Read more…

Arguments

Type IntentOptional Attributes Name
type(ocean_vcoord_t), intent(inout) :: this
real(kind=wp), intent(in) :: total_h(:,:)

Column-total depth H(i, j) (m).
real(kind=wp), intent(in) :: eta(:,:)

Free-surface anomaly η(i, j) (m).

private pure subroutine ocean_vcoord_compute_target_h_rho(this, total_h, eta, T, S, eos, hybrid)

Thin polymorphic wrapper for the isopycnal (VCOORD_RHO) and hybrid z*/isopycnal (VCOORD_HYCOM) target-grid build. Mirrors ocean_vcoord_compute_target_h (resolves the concrete type so the do concurrent kernel runs on type(ocean_vcoord_t), never a polymorphic list item). Separate from compute_target_h because the RHO/HYCOM branch needs the per-layer T/S concentrations and the device-resident EOS coefficients, which the shared pure (total_h, eta) TBP cannot carry.

Read more…

Arguments

Type IntentOptional Attributes Name
class(ocean_vcoord_t), intent(inout) :: this
real(kind=wp), intent(in) :: total_h(:,:)
real(kind=wp), intent(in) :: eta(:,:)
real(kind=wp), intent(in) :: T(:,:,:)
real(kind=wp), intent(in) :: S(:,:,:)
type(eos_t), intent(in) :: eos
logical, intent(in), optional :: hybrid

Enable the HYCOM hybrid deltas (default .false. = pure RHO).

private pure subroutine ocean_vcoord_compute_target_h_rho_impl(this, total_h, eta, T, S, eos, hybrid)

Isopycnal regrid: place layer interfaces on the prescribed rho_target(0:nz) potential-density surfaces. ONE do concurrent(j,i) over columns, each running the full density-space inversion on NZ_STACK_MAX fixed-size locals — no host loop, no per-call allocate.

Read more…

Arguments

Type IntentOptional Attributes Name
type(ocean_vcoord_t), intent(inout) :: this
real(kind=wp), intent(in) :: total_h(:,:)

Column reference depth H(i, j) (m). Caller passes the live column total (sum of h_layer) so the new grid spans it exactly.
real(kind=wp), intent(in) :: eta(:,:)

Free-surface anomaly η(i, j) (m). Added to total_h to form the column extent the new interfaces span.
real(kind=wp), intent(in) :: T(:,:,:)

Layer-mean potential temperature concentration (°C), (nx,ny,nz).
real(kind=wp), intent(in) :: S(:,:,:)

Layer-mean salinity concentration (PSU), (nx,ny,nz).
type(eos_t), intent(in) :: eos

Shared device-resident EOS handle (flat POD, by value).
logical, intent(in) :: hybrid

.true. = HYCOM (apply the monotonize + z*-floor deltas); .false. = pure RHO (bit-identical with the P2 kernel).

private subroutine ocean_vcoord_destroy(this)

Arguments

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

private subroutine ocean_vcoord_enter_data(this)

Map every host allocatable onto the device. Idempotent guard via is_init.

Arguments

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

private subroutine ocean_vcoord_enter_data_impl(this)

Arguments

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

private subroutine ocean_vcoord_exit_data(this)

Arguments

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

private subroutine ocean_vcoord_exit_data_impl(this)

Arguments

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

private pure subroutine ocean_vcoord_geometric_target(coord_type, nx, ny, nz, target_h, total_h, eta, dsig, z_ref_global, z_ref, zsigma_depth_transition, zsigma_blend_width, zstar_h_min)

Geometric target-grid kernels (EULERIAN_Z, SIGMA, ZSIGMA, ZSTAR_SIGMA, ZSTAR_FULL; formulae documented on ocean_vcoord_compute_target_h_impl). Flat on purpose — every array an explicit-shape dummy, every knob a scalar dummy — so no derived-type component and no associate-name reaches a do concurrent (see ocean_vcoord_rho_target for the GPU fault that shape caused there).

Arguments

Type IntentOptional Attributes Name
integer, intent(in), value :: coord_type

VCOORD_* family (not LAGRANGIAN / Z_FIXED / ZSTAR: the dispatcher owns those).
integer, intent(in), value :: nx

i-extent of every horizontal array (total, incl. halos).
integer, intent(in), value :: ny

j-extent of every horizontal array (total, incl. halos).
integer, intent(in), value :: nz

Number of layers; k = 1 is the bed, k = nz the surface.
real(kind=wp), intent(inout) :: target_h(nx,ny,nz)

Target layer thickness (m), bottom-up.
real(kind=wp), intent(in) :: total_h(nx,ny)

Column-total depth H (m).
real(kind=wp), intent(in) :: eta(nx,ny)

Free-surface anomaly η (m).
real(kind=wp), intent(in) :: dsig(nz)

Nominal layer fractions, bottom-up.
real(kind=wp), intent(in) :: z_ref_global(0:nz)

Global reference interface depths (ZSIGMA / ZSTAR_SIGMA).
real(kind=wp), intent(in) :: z_ref(nx,ny,0:nz)

Per-column reference interface depths (ZSTAR_FULL).
real(kind=wp), intent(in), value :: zsigma_depth_transition

Sigma → z* transition depth (m).
real(kind=wp), intent(in), value :: zsigma_blend_width

Smoothstep blend width (m).
real(kind=wp), intent(in), value :: zstar_h_min

Vanished-layer thickness (m).

private subroutine ocean_vcoord_init(this, grid, nz_ml)

Allocate every per-column array the slot owns: dsig, z_ref_global, target_h, z_ref. All sized once at init — grids don’t resize. Host allocations only; enter_data ships them to the device.

Arguments

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

private pure subroutine ocean_vcoord_rho_target(nx, ny, nz, target_h, remap_h_old, total_h, eta, t_conc, s_conc, dsig, floor_dz, rho_target, eos, p_ref, h_min, h_floor_eff, hybrid, floor_mode, h_nominal)

Column kernel of the RHO / HYCOM regrid (algorithm: see ocean_vcoord_compute_target_h_rho_impl). Flat on purpose: every array is an explicit-shape dummy and every knob a scalar dummy — no derived-type component and no associate reaches the do concurrent.

Read more…

Arguments

Type IntentOptional Attributes Name
integer, intent(in), value :: nx

i-extent of every horizontal array (total, incl. halos).
integer, intent(in), value :: ny

j-extent of every horizontal array (total, incl. halos).
integer, intent(in), value :: nz

Number of layers; k = 1 is the bed, k = nz the surface.
real(kind=wp), intent(inout) :: target_h(nx,ny,nz)

Target layer thickness (m), bottom-up.
real(kind=wp), intent(in) :: remap_h_old(nx,ny,nz)

Pre-remap layer thickness snapshot (m), bottom-up.
real(kind=wp), intent(in) :: total_h(nx,ny)

Column reference depth H (m).
real(kind=wp), intent(in) :: eta(nx,ny)

Free-surface anomaly η (m).
real(kind=wp), intent(in) :: t_conc(nx,ny,nz)

Layer-mean potential temperature (°C), bottom-up.
real(kind=wp), intent(in) :: s_conc(nx,ny,nz)

Layer-mean salinity (PSU), bottom-up.
real(kind=wp), intent(in) :: dsig(nz)

Nominal layer fractions, bottom-up (dsig(nz) = surface) — the HYCOM floor only under HYCOM_FLOOR_SIGMA.
real(kind=wp), intent(in) :: floor_dz(nz)

HYCOM z* nominal layer thicknesses (m), bottom-up (floor_dz(nz) = surface) — read only under HYCOM_FLOOR_PROFILE.
real(kind=wp), intent(in) :: rho_target(0:nz)

Target potential densities (kg/m³), 0 = lightest = surface.
type(eos_t), intent(in) :: eos

Shared EOS handle (flat POD).
real(kind=wp), intent(in), value :: p_ref

Coordinate reference pressure (Pa) — rho_ref_pressure.
real(kind=wp), intent(in), value :: h_min

Pre-compaction strip threshold (m) — zstar_h_min.
real(kind=wp), intent(in), value :: h_floor_eff

Min-thickness inflation floor (m), > H_VANISHED.
logical, intent(in), value :: hybrid

.true. = HYCOM deltas; .false. = pure RHO.
integer, intent(in), value :: floor_mode

HYCOM z* floor source: HYCOM_FLOOR_PROFILE (floor_dz), HYCOM_FLOOR_UNIFORM (h_nominal) or HYCOM_FLOOR_SIGMA (dsig·H, the unconfigured-slot fallback).
real(kind=wp), intent(in), value :: h_nominal

Uniform z* nominal thickness (m) for HYCOM_FLOOR_UNIFORM.

private pure subroutine ocean_vcoord_rho_target_column(i, j, nx, ny, nz, target_h, remap_h_old, total_h, eta, t_conc, s_conc, dsig, floor_dz, rho_target, eos, p_ref, h_min, h_floor_eff, hybrid, floor_mode, h_nominal)

One column of the RHO / HYCOM regrid (steps 0-5 of ocean_vcoord_compute_target_h_rho_impl) — the per-thread body of ocean_vcoord_rho_target. Same module as its caller (the project rule for !$acc routine seq callees of a do concurrent).

Arguments

Type IntentOptional Attributes Name
integer, intent(in), value :: i

Column i-index.
integer, intent(in), value :: j

Column j-index.
integer, intent(in), value :: nx

i-extent of every horizontal array (total, incl. halos).
integer, intent(in), value :: ny

j-extent of every horizontal array (total, incl. halos).
integer, intent(in), value :: nz

Number of layers; k = 1 is the bed, k = nz the surface.
real(kind=wp), intent(inout) :: target_h(nx,ny,nz)

Target layer thickness (m), bottom-up.
real(kind=wp), intent(in) :: remap_h_old(nx,ny,nz)

Pre-remap layer thickness snapshot (m), bottom-up.
real(kind=wp), intent(in) :: total_h(nx,ny)

Column reference depth H (m).
real(kind=wp), intent(in) :: eta(nx,ny)

Free-surface anomaly η (m).
real(kind=wp), intent(in) :: t_conc(nx,ny,nz)

Layer-mean potential temperature (°C), bottom-up.
real(kind=wp), intent(in) :: s_conc(nx,ny,nz)

Layer-mean salinity (PSU), bottom-up.
real(kind=wp), intent(in) :: dsig(nz)

Nominal layer fractions, bottom-up (dsig(nz) = surface) — the HYCOM floor only under HYCOM_FLOOR_SIGMA.
real(kind=wp), intent(in) :: floor_dz(nz)

HYCOM z* nominal layer thicknesses (m), bottom-up (floor_dz(nz) = surface) — read only under HYCOM_FLOOR_PROFILE.
real(kind=wp), intent(in) :: rho_target(0:nz)

Target potential densities (kg/m³), 0 = lightest = surface.
type(eos_t), intent(in) :: eos

Shared EOS handle (flat POD).
real(kind=wp), intent(in), value :: p_ref

Coordinate reference pressure (Pa) — rho_ref_pressure.
real(kind=wp), intent(in), value :: h_min

Pre-compaction strip threshold (m) — zstar_h_min.
real(kind=wp), intent(in), value :: h_floor_eff

Min-thickness inflation floor (m), > H_VANISHED.
logical, intent(in), value :: hybrid

.true. = HYCOM deltas; .false. = pure RHO.
integer, intent(in), value :: floor_mode

HYCOM z* floor source: HYCOM_FLOOR_PROFILE (floor_dz), HYCOM_FLOOR_UNIFORM (h_nominal) or HYCOM_FLOOR_SIGMA (dsig·H, the unconfigured-slot fallback).
real(kind=wp), intent(in), value :: h_nominal

Uniform z* nominal thickness (m) for HYCOM_FLOOR_UNIFORM.