rdb_ocean_pressure_force Module

Carries scheme-variant flags + reusable workspace for the layered hydrostatic pressure-force kernels on C-grid face arrays: Montgomery potential, finite-volume (lite / Wright / MOM6) and the NK=2 reduced-gravity form.

The FV variants integrate pressure per column from the surface (p=0) down to the bed, p_edge(k) = p_edge(k+1) + g * rho_layer(k) * h_layer(k) p_centre(k) = 0.5 * (p_edge(k) + p_edge(k+1)) and difference cell-centred pressure across each face at constant z (pressure difference + g*rho_face*dz_centre correction) for the acceleration du/dt = -(1/rho0) * dp/dx.

The Montgomery variant instead builds the Boussinesq Montgomery potential M = p/rho0 + (g*rho/rho0)*z — which is CONSTANT within a layer of uniform density — by a vertical recursion, and takes ONE horizontal difference of M. See the OPGF_VARIANT_* constants below.

MONT, FV_LITE and FV_WRIGHT are measured from the FREE SURFACE (z = 0 at the surface, negative below), so they carry no barotropic -g*grad(eta) term: the split-explicit barotropic substep owns that. FV_MOM6 (pa(nz+1) = rho_ref*g*eta_geo) and GPRIME (-g_FS*grad(eta) in the top layer) DO carry it; the split sheds exactly that term from the barotropic forcing (pgf_free_surface_gravity, set_fast_forcing_eta_pf) and keeps the rest of the depth mean.


Uses

  • module~~rdb_ocean_pressure_force~~UsesGraph module~rdb_ocean_pressure_force rdb_ocean_pressure_force iso_fortran_env iso_fortran_env module~rdb_ocean_pressure_force->iso_fortran_env module~rdb_constants rdb_constants module~rdb_ocean_pressure_force->module~rdb_constants module~rdb_eos rdb_eos module~rdb_ocean_pressure_force->module~rdb_eos module~rdb_grid rdb_grid module~rdb_ocean_pressure_force->module~rdb_grid module~rdb_mem_report rdb_mem_report module~rdb_ocean_pressure_force->module~rdb_mem_report module~rdb_multilayer_state rdb_multilayer_state module~rdb_ocean_pressure_force->module~rdb_multilayer_state module~rdb_ocean_metrics rdb_ocean_metrics module~rdb_ocean_pressure_force->module~rdb_ocean_metrics module~rdb_ocean_pgf_reconstruct rdb_ocean_pgf_reconstruct module~rdb_ocean_pressure_force->module~rdb_ocean_pgf_reconstruct module~rdb_scratch_3d rdb_scratch_3d module~rdb_ocean_pressure_force->module~rdb_scratch_3d 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_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_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_status rdb_ocean_status 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_pgf_reconstruct->module~rdb_constants module~rdb_ocean_pgf_reconstruct->module~rdb_eos module~rdb_scratch_3d->iso_fortran_env module~rdb_scratch_3d->module~rdb_constants module~rdb_scratch_3d->module~rdb_mem_report module~rdb_efp->iso_fortran_env ieee_arithmetic ieee_arithmetic module~rdb_efp->ieee_arithmetic 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_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_pressure_force~~UsedByGraph module~rdb_ocean_pressure_force rdb_ocean_pressure_force module~rdb_barotropic_coupling rdb_barotropic_coupling module~rdb_barotropic_coupling->module~rdb_ocean_pressure_force module~rdb_ocean_bt_budget_probe rdb_ocean_bt_budget_probe module~rdb_ocean_bt_budget_probe->module~rdb_ocean_pressure_force module~rdb_ocean_dyn rdb_ocean_dyn module~rdb_ocean_dyn->module~rdb_ocean_pressure_force module~rdb_ocean_dyn->module~rdb_barotropic_coupling module~rdb_ocean_dyn->module~rdb_ocean_bt_budget_probe module~rdb_ocean_setup rdb_ocean_setup module~rdb_ocean_setup->module~rdb_ocean_pressure_force module~rdb_ocean_setup->module~rdb_ocean_dyn module~rdb_ocean_state rdb_ocean_state module~rdb_ocean_setup->module~rdb_ocean_state module~rdb_ocean_state->module~rdb_ocean_pressure_force module~rdb_ocean_state->module~rdb_ocean_dyn proc~validate_config validate_config proc~validate_config->module~rdb_ocean_pressure_force 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_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_setup module~rdb_ocean_engine->module~rdb_ocean_state module~rdb_ocean_engine->module~rdb_ocean_diag_derived module~rdb_ocean_engine->module~rdb_ocean_diag_fills

Variables

Type Visibility Attributes Name Initial
integer, public, parameter :: OPGF_VARIANT_FV_LITE = 2

Finite-volume PGF with the z-position correction PGF_x = -(1/rho_0) * [(p_R - p_L)/dx + grho_face(z_R - z_L)/dx] that differences pressure at constant z rather than at constant layer index. Reduces to MONT for aligned columns; cancels the spurious bottom-current PGF when h_layer varies across columns. Reuses the layer-mean rho_layer from the EOS (no in-layer quadrature — use FV_WRIGHT for that).
integer, public, parameter :: OPGF_VARIANT_FV_MOM6 = 5

Faithful layer-integrated FV-Bouss PGF (Boussinesq, per-layer Rlay density). Unlike FV_LITE (layer-centre pressure with a z-correction), uses layer-integrated pressure differences divided by face-averaged thickness:

PFu(I,j,k) = [ (pa·h + intz_dpa)_L − (pa·h + intz_dpa)_R + (h_R − h_L) · intx_pa − (e_bot_R − e_bot_L) · intx_dpa ] · (2 · I_Rho0 · IdxCu) / (h_L + h_R + h_neglect)

pa is the pressure anomaly stack relative to rho_ref·g·z, intz_dpa the in-layer vertical integral, intx_pa/intx_dpa the horizontal face integrals. The face-thickness divisor self- regulates bed-layer force at thin shelf-break cells (the mechanism FV_LITE+rho_init lacked, driving a 33× WBC bed overshoot).

Density per layer from ms%rho_layer. rho_ref defaults to rho0 so the surface layer’s anomaly contribution vanishes. per-cell convention: our k=1 is bed, MOM6’s is surf.
integer, public, parameter :: OPGF_VARIANT_FV_WRIGHT = 3

FV_LITE z-correction plus in-situ density along the integration path: each layer’s centre density is re-evaluated by Wright at the pressure from a single Picard step (ms%rho_layer seeds the half-layer pressure). The pressure stack and rho_face both use the in-situ value, adding the compressibility piece FV_LITE misses; reduces to FV_LITE for incompressible water.

Caller must have run ocean_eos_compute with EOS_VARIANT_WRIGHT_97 first (the Picard seed is the nonlinear surface ρ, not a Boussinesq constant).
integer, public, parameter :: OPGF_VARIANT_GPRIME = 4

Reduced-gravity / gprime PGF (reduced-gravity analogue). Per-layer acceleration a_k = -g_FS · ∇η - Σ_{j>k} g’_j · ∇η_j_interface. NK = 2 only (asserted at init): a_top = -g_FS · ∂(h_1 + h_2)/∂x a_bot = -g_FS · ∂(h_1 + h_2)/∂x - g’_int · ∂h_1/∂x pgf%gprime_gfs is g_FS (m/s²), pgf%gprime_gint the single internal g’. No EOS, no Picard — densities fixed by GFS/GINT.
integer, public, parameter :: OPGF_VARIANT_MONT = 1

Boussinesq Montgomery-potential form (cheapest, layered, no EOS calls along the integration path, and no pressure stack at all).

For a layer of horizontally uniform density the Boussinesq Montgomery potential M = p/rho0 + rho_starz, rho_star = grho_layer/rho0 is CONSTANT through the layer (moving up by dz costs -g*rho*dz/rho0 of p/rho0 and gains exactly rho_star*dz), and -(1/rho0) * dp/dx|_z == -dM/dx along the layer. That is what makes ONE horizontal difference of M legitimate: M already carries the geopotential.

Built by a vertical recursion, BOTTOM-UP here (k=1 bed, k=nz surface — MOM6 runs the same recursion top-down). Seeded at the free surface, where the surface-relative interface height and the pressure are both zero:

e_edge(k) = z_centre(k) + 0.5*h_layer(k) ! TOP of layer k M(nz) = 0 M(k) = M(k+1) + (rho_star(k) - rho_star(k+1)) * e_edge(k+1)

(e_edge(k+1) is the interface SHARED by layers k and k+1, where p and z agree, so the whole jump in M is the jump in rho_star.) Then, per face,

PGF_x = -(M_R - M_L)idxCu + (rho_star_R - rho_star_L)z_effidxCu z_eff = (e_Lh_R + e_Rh_L - h_Lh_R) / (h_L + h_R)

There is NO 1/rho0 multiplying dM: M already has units of geopotential (m^2 s^-2), so its horizontal gradient IS an acceleration. The 1/rho0 lives inside rho_star.

The second term is load-bearing, not a refinement. -dM/dx is the PGF only where rho_layer is horizontally uniform within the layer; where it is not, the exact relation picks up + z * d(rho_star)/dx, and z_eff is the thickness-weighted height at which to evaluate it. Drop it and the scheme gets the SIGN of a horizontal density contrast wrong. On aligned columns (h_L == h_R) z_eff collapses to the arithmetic mean layer centre and the whole expression reduces ALGEBRAICALLY to FV_LITE.

Exact at rest in an isopycnal (VCOORD_LAGRANGIAN) column over ANY bathymetry as long as the layer is present on both sides of the face: flat isopycnals make every e_edge horizontally uniform and every rho_star difference zero, so M is uniform and the PGF is zero to round-off. Where an isopycnal layer has GROUNDED on one side the two columns share no common depth for that layer and the residual returns — skip_nonoverlap gates exactly that face, and honours mont for the same geometric reason it honours the FV variants.
integer, private, parameter :: N_BOOLE = 5

Sub-points of the in-layer Boole (5-point closed Newton-Cotes) rule; the same value as rdb_ocean_pgf_reconstruct’s.

Derived Types

type, public ::  ocean_pressure_force_t

Components

Type Visibility Attributes Name Initial
real(kind=wp), public, allocatable :: b(:,:)
type(scratch_3d_buffer_t), public :: conc_S

Layer-mean salinity, as conc_T.

Read more…
type(scratch_3d_buffer_t), public :: conc_T

Layer-mean temperature the in-situ PCM and reconstruction kernels read, shape (nx, ny, nz): hT/h on a live layer, the I1′ donor’s on a vanished one (Pass C of both kernels). Filled once per call by a per-column pass, so no later pass divides by a thickness or walks a column. Unused by the rho_layer path.
type(scratch_3d_buffer_t), public :: dpdx_face

East-face PGF acceleration: -(1/rho0) * dp/dx at u-face. Shape (nx+1, ny, nz_ml). Apply adds dt*dpdx to u_face_x_layer.
type(scratch_3d_buffer_t), public :: dpdy_face

North-face PGF acceleration. Shape (nx, ny+1, nz_ml).

Read more…
type(scratch_3d_buffer_t), public :: e_face

Interface heights at cell centres, shape (nx, ny, nz+1). Pressure anomaly stack at interfaces, units Pa. Relative to the rho_ref · g · z baseline so the surface-pressure contribution is just pa(top) = rho_ref · g · η. Built by marching down from the surface; pa(k) − pa(k+1) = dpa(k) where dpa(k) = (Rlay(k) − rho_ref) · g · h(k).
real(kind=wp), public :: gfs_scale = 1.0_wp

Free-surface gravity scaling (= GFS / G_EARTH). Default 1.0 ⇒ pure FV_MOM6, bit-identical. When < 1, Pass 5 applies the Montgomery dM correction subtracting (1 - gfs_scale)·(g/ρ₀)·ρ_surf·∇η from every layer’s PGF (depth-independent, so the BT mass-flux invariant survives); the driver also drops bt_work%g_bt to gfs_scale·GRAVITY. Combined ⇒ wave speed sqrt(gfs_scale · g · H).
real(kind=wp), public :: gprime_gfs = 9.81_wp

Free-surface reduced gravity (m/s²). Full physical g reduces to a standard free-surface model; reducing it slows the BT mode (wave speed sqrt(gfs · H)) for a slower BT CFL.
real(kind=wp), public :: gprime_gint = 0.0098_wp

Internal reduced gravity (m/s²) at the single layer-1/layer-2 interface. Only used by OPGF_VARIANT_GPRIME.
real(kind=wp), public :: h_neglect = 1.0e-10_wp

Face-thickness floor in the FV_MOM6 denominator (h_L + h_R + h_neglect). Prevents division by zero when both adjacent cells have vanishing bed layers.
logical, public :: insitu_density = .true.

FV_MOM6 constant-by-layer (PCM) density at its IN-SITU pressure (&ocean_pgf_nml insitu_density, MOM6 parity). .true. (default): each layer’s dpa/intz_dpa and the cross-face intx_dpa/inty_dpa integrate EOS(T, S, p = -g·rho0·z) with the layer-mean T/S (compute_fv_mom6_insitu_pcm_impl). Under Wright the vertical integral is ANALYTIC (MOM6 int_density_dz_wright: one polynomial evaluation per layer and per lateral sub-column); under Roquet the 5-point Boole rule of MOM6 int_density_dz_generic_pcm with the (T, S) part of the EOS hoisted out of the pressure points (roquet_pcm_dpa_intz: 1 T/S polynomial + 5 pressure Horners per layer, 3 + 15 per face). .false.: the legacy PCM integral of ms%rho_layer, a POTENTIAL density at the single &ocean_eos_nml p_ref, which drops the pressure dependence of the horizontal density gradient below the reference level. Consulted only by FV_MOM6 with reconstruct_for_pressure = .false., an EOS handle and T/S, and only for a PRESSURE-DEPENDENT EOS (Wright, Roquet): for the linear EOS in-situ and potential density coincide, so the legacy path runs and answers are bit-identical.
type(scratch_3d_buffer_t), public :: intx_dpa

u-face ∫ dpa, shape (nx+1, ny, nz), units Pa·m.
type(scratch_3d_buffer_t), public :: intx_pa

u-face ∫ pa across x, shape (nx+1, ny, nz+1), units Pa·m.
type(scratch_3d_buffer_t), public :: inty_dpa

v-face ∫ dpa, shape (nx, ny+1, nz), units Pa·m.
type(scratch_3d_buffer_t), public :: inty_pa

v-face ∫ pa across y, shape (nx, ny+1, nz+1), units Pa·m.
type(scratch_3d_buffer_t), public :: intz_dpa

Per-layer ∫ dpa dz, shape (nx, ny, nz), units Pa·m.
logical, public :: is_init = .false.

True between init and destroy. Prefer this to allocated(...) — tracks GPU device attachment too.
logical, public :: mass_weight = .false.

FV_MOM6 shelf-break hWght mass-weighting toggle. When .true. Pass-2’s horizontal pressure integral biases the face density toward the thinner column at unequal-depth faces, cancelling the spurious bottom-layer shelf-break PGF. Equal-depth columns ⇒ reduces exactly to the midpoint average ⇒ bit-identical.
type(scratch_3d_buffer_t), public :: mont_M

Boussinesq Montgomery potential M at layer centres (m² s⁻²). Shape (nx, ny, nz_ml). Written by the MONT column recursion and read by its two face passes; no other variant touches it.
real(kind=wp), public :: nonoverlap_vanish_tol = 2.0_wp*H_VANISHED

Thickness (m) at or below which a layer counts as GROUNDED for the skip_nonoverlap gate: a face is zeroed only where the layer’s z-extents do not overlap AND the layer is this thin on at least one side. Non-overlap alone is not grounding — a layer that is MASSIVE on both sides but sits at different depths (a sigma-seeded stack over a bathymetric step, where a 10 m layer at 90-100 m faces a 25 m layer at 225-250 m) is a steep coordinate surface the FV forms are built for, exactly as under vcoord_type="sigma". Zeroing its PGF while continuity keeps moving its mass across the face breaks the PGF-work / PE exchange and grows energy without bound (the compat-matrix staircase, MaxCFL panic at step 178). Set by the driver from nonoverlap_vanish_tol_for(angstrom_h); the default matches angstrom_h = 0.
type(scratch_3d_buffer_t), public :: p_edge

Layer-edge pressure stack at cell centres. Shape (nx, ny, nz_ml+1). k=1 bed, k=nz_ml+1 surface. Reused across the two horizontal-gradient passes.
logical, public :: p_top_in_bc = .false.

FV_MOM6 top-of-column pressure in the surface boundary condition (&ocean_pgf_nml p_top_in_bc). .false. (default): pa(nz+1) = rho_ref·g·eta_geo, bit-identical. .true.: the load multilayer_state_t%p_top (Pa) is ADDED there, so the pressure stack measures down from the loaded surface — pa(nz+1) = rho_ref·g·eta_geo + p_top. Only consulted by FV_MOM6 (both the PCM and the reconstruct_for_pressure branch); fail-loud at configure for any other variant, which carries no injectable pa stack.

Read more…
type(scratch_3d_buffer_t), public :: pa

Pressure anomaly at interfaces, shape (nx, ny, nz+1). Per-layer vertical integral of dpa from layer top inward. For Boussinesq Rlay path: intz_dpa(k) = 0.5 · (Rlay(k) − rho_ref) · g · h(k)² (mid-point rule).
type(scratch_3d_buffer_t), public :: recon_S_b
type(scratch_3d_buffer_t), public :: recon_S_t
type(scratch_3d_buffer_t), public :: recon_T_b
type(scratch_3d_buffer_t), public :: recon_T_t
integer, public :: recon_scheme = PGF_RECON_PLM

In-layer reconstruction scheme: 1 = PLM, 2 = PPM. Only consulted when reconstruct_for_pressure = .true..
logical, public :: reconstruct_for_pressure = .false.

FV_MOM6 in-layer T/S reconstruction toggle. .false. (default): layer-mean (PCM) density ⇒ bit-identical. .true.: per-layer dpa / intz_dpa from a 5-point Boole quadrature of a monotone PLM/PPM sub-layer T/S profile (Adcroft, Hallberg & Harrison 2008; White, Adcroft & Hallberg 2009), removing the spurious-PGF error on thick/sloped layers. Only consulted by FV_MOM6.
real(kind=wp), public :: rho0 = 1035.0_wp

BOUSSINESQ reference density (kg/m³) — the ρ₀ that divides the pressure gradient into an acceleration, du/dt = −(1/ρ₀)·∂p/∂x. Read by EVERY variant (it is inv_rho0 in the face passes and g_over_rho0 in the Montgomery recursion), and by compute_pbce in the barotropic coupling.

Read more…
type(scratch_3d_buffer_t), public :: rho_insitu

In-situ layer-centre density from the Wright column-sweep Picard step. Shape (nx, ny, nz_ml). Populated only by the FV_WRIGHT branch; the other variants source ρ from ms%rho_layer.
real(kind=wp), public :: rho_ref = 1035.0_wp

ANOMALY reference density (kg/m³) subtracted from layer densities when building the pa pressure-anomaly stack (FV_MOM6), and the surface-layer g·ρ_ref/ρ₀ in compute_pbce. rho_ref = rho0 = ρ_surf makes the surface layer’s anomaly vanish; only denser bed layers contribute.

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

ALLOCATION GATE for the variant-specific scratch buffers.

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

Grounded-layer PGF gate (&ocean_isopycnal_nml pgf_skip_nonoverlap, set by the driver ONLY under VCOORD_LAGRANGIAN). When .true. the face PGF is zeroed wherever the layer’s z-extents in the two abutting columns do NOT overlap — a layer that has wedged out against the bed on one side, where the two-point Jacobian Δp_centre + g·ρ_layer·Δz_centre has no common depth to difference across and leaves g·(ρ_layer − ρ̄_ambient)·∂z/∂x of acceleration on a RESTING ocean. .false. ⇒ the gate branch is never taken ⇒ bit-identical.

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

Convenience flag — re-evaluate density at each integration sample rather than at layer centres. Phase 5d.
integer, public :: variant = OPGF_VARIANT_MONT

Active pressure-force scheme variant. MONT is the default here AND the &ocean_pgf_nml form default, so a state built straight from this type and one built through the config agree. MONT is general-purpose: valid over variable bathymetry and every vcoord.
type(scratch_3d_buffer_t), public :: z_centre

Layer-centre z at cell centres, from the free surface downward (z=0 surface, negative below). Shape (nx, ny, nz_ml). Filled for MONT and the FV_LITE-family branches. Surface-relative (NOT bed-relative) so the σ-coord Jacobian cancels the spurious cross-bathymetry pressure gradient — and so that neither form double-counts the barotropic -g·∇η the BT substep already carries. MONT reads it as the interface height e_edge(k) = z_centre(k) + 0.5·h_layer(k).

Type-Bound Procedures

procedure, public, non_overridable :: bytes => ocean_pressure_force_bytes
procedure, public, non_overridable :: destroy => ocean_pressure_force_destroy
procedure, public, non_overridable :: enter_data => ocean_pressure_force_enter_data
procedure, public, non_overridable :: exit_data => ocean_pressure_force_exit_data
procedure, public, non_overridable :: init => ocean_pressure_force_init
procedure, public, non_overridable :: set_bathymetry => ocean_pressure_force_set_bathymetry

Functions

public pure function gprime_nz_is_supported(variant, nz) result(ok)

.true. unless variant == OPGF_VARIANT_GPRIME with nz /= 2 (PR-6 fail-loud). The reduced-gravity gprime PGF hard-writes ONLY k=1 (bottom) and k=2 (top) — with nz > 2 layers k=3..nz carry zero pressure gradient AND the “top” branch lands on layer 2 of nz (the abyss under the bottom-up convention); with nz < 2 the kernel early-returns leaving the whole PGF zero. The guard is gprime-specific — every other variant supports general nz, so this returns .true. for them regardless of nz. Wired into validate_config against cfg%nz_layers (config-time).

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: variant
integer, intent(in) :: nz

Return Value logical

public pure function nonoverlap_vanish_tol_for(angstrom_h) result(tol)

The grounded-layer gate’s “vanished on one side” threshold, 2·max(angstrom_h, H_VANISHED): a layer within a factor two of the floor it can rest on. The factor is the margin that keeps every floor the Lagrangian path parks a grounded layer on inside the band — the uniform_z seed’s max(angstrom_h, 2·H_VANISHED) collapse, the conservative-floor borrow’s angstrom_h (to round-off), and the positive-definite continuity’s [angstrom_h, angstrom_h + H_VANISHED] at-floor band — while staying orders of magnitude below any layer that carries real mass.

Arguments

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

&ocean_isopycnal_nml angstrom_h (m); 0 when the floor is off.

Return Value real(kind=wp)

public pure function parse_opgf_variant(name) result(code)

Translate a namelist string into an OPGF_VARIANT_* code. Unrecognised values fall back to OPGF_VARIANT_FV_LITE (the production default).

Arguments

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

Return Value integer

private pure function ocean_pressure_force_bytes(this) result(nbytes)

Counted allocatable footprint of the pressure force 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_pressure_force_t), intent(in) :: this

Return Value integer(kind=int64)

private pure function recon_rho_surf(pa_k, pa_kp1, h_surf, rho_ref) result(rho_surf)

Recover the layer-mean surface density from the reconstructed pressure-anomaly stack: dpa(nz) = pa(nz) - pa(nz+1) = (rho_surf - rho_ref)gh_surf, so rho_surf = rho_ref + dpa/(g*h). Used only by the gfs_scale Montgomery correction (Pass 5).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: pa_k
real(kind=wp), intent(in) :: pa_kp1
real(kind=wp), intent(in) :: h_surf
real(kind=wp), intent(in) :: rho_ref

Return Value real(kind=wp)

private pure function use_insitu_pcm(pgf, ms, eos) result(yes)

Does the FV_MOM6 constant-by-layer branch take the IN-SITU density path (compute_fv_mom6_insitu_pcm_impl)? Only when it can change the answer: the knob is on, there is an EOS handle and T/S to evaluate it on, and the EOS depends on pressure. For the linear EOS ms%rho_layer already IS the in-situ density, so the legacy path runs, bit-identical.

Arguments

Type IntentOptional Attributes Name
type(ocean_pressure_force_t), intent(in) :: pgf
type(multilayer_state_t), intent(in) :: ms
type(eos_t), intent(in), optional :: eos

Return Value logical

private pure function wright_rho(t, s, p) result(rho)

Wright (1997) in-situ density (kg/m^3) – term for term the Wright branch of eos_density_point, without the handle.

Arguments

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

Temperature (degC), salinity (PSU), pressure (Pa).
real(kind=wp), intent(in) :: s

Temperature (degC), salinity (PSU), pressure (Pa).
real(kind=wp), intent(in) :: p

Temperature (degC), salinity (PSU), pressure (Pa).

Return Value real(kind=wp)


Subroutines

public pure subroutine boole_dpa_face_wright(rho0, rho_ref, e_top_l, e_top_r, dz_l, dz_r, t_t_l, t_b_l, t_m_l, t_t_r, t_b_r, t_m_r, s_t_l, s_b_l, s_m_l, s_t_r, s_b_r, s_m_r, dpa_l, dpa_r, parabolic, dpa_face)

boole_dpa_face with the Wright (1997) vertical rule boole_dpa_intz_layer_wright at each of the five sub-columns – identical sub-columns, weights and summation order, no eos_t handle. 15 inline Wright density evaluations per face per layer (the end points are the columns’ own dpa).

Arguments

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

Pressure-estimate density / anomaly reference (kg/m^3).
real(kind=wp), intent(in) :: rho_ref

Pressure-estimate density / anomaly reference (kg/m^3).
real(kind=wp), intent(in) :: e_top_l

Left / right shallower-interface heights and thicknesses (m).
real(kind=wp), intent(in) :: e_top_r

Left / right shallower-interface heights and thicknesses (m).
real(kind=wp), intent(in) :: dz_l

Left / right shallower-interface heights and thicknesses (m).
real(kind=wp), intent(in) :: dz_r

Left / right shallower-interface heights and thicknesses (m).
real(kind=wp), intent(in) :: t_t_l

Left / right temperature triples (top, bottom, mean).
real(kind=wp), intent(in) :: t_b_l

Left / right temperature triples (top, bottom, mean).
real(kind=wp), intent(in) :: t_m_l

Left / right temperature triples (top, bottom, mean).
real(kind=wp), intent(in) :: t_t_r

Left / right temperature triples (top, bottom, mean).
real(kind=wp), intent(in) :: t_b_r

Left / right temperature triples (top, bottom, mean).
real(kind=wp), intent(in) :: t_m_r

Left / right temperature triples (top, bottom, mean).
real(kind=wp), intent(in) :: s_t_l

Left / right salinity triples.
real(kind=wp), intent(in) :: s_b_l

Left / right salinity triples.
real(kind=wp), intent(in) :: s_m_l

Left / right salinity triples.
real(kind=wp), intent(in) :: s_t_r

Left / right salinity triples.
real(kind=wp), intent(in) :: s_b_r

Left / right salinity triples.
real(kind=wp), intent(in) :: s_m_r

Left / right salinity triples.
real(kind=wp), intent(in) :: dpa_l

The left / right columns’ own layer dpa (Pa), the end points.
real(kind=wp), intent(in) :: dpa_r

The left / right columns’ own layer dpa (Pa), the end points.
logical, intent(in) :: parabolic

.true. -> PPM curvature.
real(kind=wp), intent(out) :: dpa_face

Along-face mean of g * int rho' dz over the layer (Pa).

public pure subroutine boole_dpa_intz_layer_wright(rho0, rho_ref, e_top, dz, t_t, t_b, t_mean, s_t, s_b, s_mean, parabolic, dpa, intz_dpa)

boole_dpa_intz_layer specialised to Wright (1997): the same five sub-points and weights, with the density written inline (wright_rho, the expression eos_density_point evaluates) – no eos_t handle, no per-point variant dispatch. The caller selects this twin ONCE, outside its loops.

Arguments

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

Pressure-estimate density / anomaly reference (kg/m^3).
real(kind=wp), intent(in) :: rho_ref

Pressure-estimate density / anomaly reference (kg/m^3).
real(kind=wp), intent(in) :: e_top

Shallower-interface height / layer thickness (m).
real(kind=wp), intent(in) :: dz

Shallower-interface height / layer thickness (m).
real(kind=wp), intent(in) :: t_t

T / S top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: t_b

T / S top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: t_mean

T / S top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: s_t

T / S top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: s_b

T / S top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: s_mean

T / S top edge, bottom edge, layer mean.
logical, intent(in) :: parabolic

.true. -> PPM curvature.
real(kind=wp), intent(out) :: dpa

As boole_dpa_intz_layer.
real(kind=wp), intent(out) :: intz_dpa

As boole_dpa_intz_layer.

public subroutine ocean_pressure_force_apply(pgf, ms, dt, no_wait)

Forward-Euler accumulation of the PGF acceleration onto the face velocities. Additive (not overwriting), so apply ordering vs the Coriolis apply doesn’t matter before the next tendency-compute. no_wait (optional, default .false. ⇒ blocking): when .true. the apply loops run on OpenACC queue 1 and return WITHOUT syncing, so the batched velocity-apply chain can !$acc wait(1) ONCE. Not pure (async/wait directives); still functionally pure.

Arguments

Type IntentOptional Attributes Name
type(ocean_pressure_force_t), intent(in) :: pgf
type(multilayer_state_t), intent(inout) :: ms
real(kind=wp), intent(in) :: dt
logical, intent(in), optional :: no_wait

public pure subroutine ocean_pressure_force_compute(grid, metrics, pgf, ms, eos)

Compute the hydrostatic pressure-gradient acceleration at every C-grid face. Variants (see the OPGF_VARIANT_* / GPRIME / FV_MOM6 constants for the per-variant formulas): MONT (layer-mean ρgh), FV_LITE (+ z-correction), FV_WRIGHT (+ in-situ Wright density), GPRIME, FV_MOM6.

Read more…

Arguments

Type IntentOptional Attributes Name
type(hgrid_t), intent(in) :: grid
type(ocean_metrics_t), intent(in) :: metrics

Curvilinear horizontal metrics. The u-face gradient divides by idxCu(i,j), the v-face by idyCv(i,j). On uniform Cartesian idxCu == 1/dx bitwise (byte-identical to scalar inv_dx/inv_dy).
type(ocean_pressure_force_t), intent(inout) :: pgf
type(multilayer_state_t), intent(in) :: ms
type(eos_t), intent(in), optional :: eos

EOS handle — REQUIRED when pgf%reconstruct_for_pressure is on (the in-layer Boole quadrature evaluates the EOS at each sub-point). Optional so the legacy PCM call sites (and the non-reconstruct variant tests) need not thread it through.

public pure subroutine plm_edges_layer(k, nz, h_dn, h_c, h_up, q_dn, q_c, q_up, q_t, q_b)

PLM top/bottom edge values of ONE layer k of a layer-mean field q, via a two-stage h-weighted van-Leer slope (White, Adcroft & Hallberg 2009 §2). Returns the SHALLOWER edge in q_t (toward k+1) and the DEEPER edge in q_b (toward k-1), bottom-up. Boundary layers (k=1, k=nz) -> boundary_edges_linear, the linear-exact one-sided pair.

Read more…

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: k

Layer index (1 = bed, nz = surface).
integer, intent(in) :: nz

Number of layers in the column.
real(kind=wp), intent(in) :: h_dn

Thicknesses (m) of layers k-1 (deeper), k, k+1 (shallower).
real(kind=wp), intent(in) :: h_c

Thicknesses (m) of layers k-1 (deeper), k, k+1 (shallower).
real(kind=wp), intent(in) :: h_up

Thicknesses (m) of layers k-1 (deeper), k, k+1 (shallower).
real(kind=wp), intent(in) :: q_dn

Layer means of layers k-1, k, k+1.
real(kind=wp), intent(in) :: q_c

Layer means of layers k-1, k, k+1.
real(kind=wp), intent(in) :: q_up

Layer means of layers k-1, k, k+1.
real(kind=wp), intent(out) :: q_t

Top (shallower) edge value of layer k.
real(kind=wp), intent(out) :: q_b

Bottom (deeper) edge value of layer k.

public pure subroutine ppm_edges_layer(k, nz, h_m2, h_m1, h_c, h_p1, h_p2, s_m2, s_m1, s_c, s_p1, s_p2, t_m2, t_m1, t_c, t_p1, t_p2, s_top, s_bot, t_top, t_bot)

PPM top/bottom edge values of salinity and temperature in ONE layer k: the implicit-h4 interface estimates (ppm_interface_values) + the Colella & Woodward (1984) limiter (ppm_limit_edges). Returns the SHALLOWER edges in *_top, the DEEPER in *_bot (bottom-up).

Read more…

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: k

Layer index (1 = bed, nz = surface).
integer, intent(in) :: nz

Number of layers in the column.
real(kind=wp), intent(in) :: h_m2

Thicknesses (m) of layers k-2 .. k+2.
real(kind=wp), intent(in) :: h_m1

Thicknesses (m) of layers k-2 .. k+2.
real(kind=wp), intent(in) :: h_c

Thicknesses (m) of layers k-2 .. k+2.
real(kind=wp), intent(in) :: h_p1

Thicknesses (m) of layers k-2 .. k+2.
real(kind=wp), intent(in) :: h_p2

Thicknesses (m) of layers k-2 .. k+2.
real(kind=wp), intent(in) :: s_m2

Salinity layer means of layers k-2 .. k+2.
real(kind=wp), intent(in) :: s_m1

Salinity layer means of layers k-2 .. k+2.
real(kind=wp), intent(in) :: s_c

Salinity layer means of layers k-2 .. k+2.
real(kind=wp), intent(in) :: s_p1

Salinity layer means of layers k-2 .. k+2.
real(kind=wp), intent(in) :: s_p2

Salinity layer means of layers k-2 .. k+2.
real(kind=wp), intent(in) :: t_m2

Temperature layer means of layers k-2 .. k+2.
real(kind=wp), intent(in) :: t_m1

Temperature layer means of layers k-2 .. k+2.
real(kind=wp), intent(in) :: t_c

Temperature layer means of layers k-2 .. k+2.
real(kind=wp), intent(in) :: t_p1

Temperature layer means of layers k-2 .. k+2.
real(kind=wp), intent(in) :: t_p2

Temperature layer means of layers k-2 .. k+2.
real(kind=wp), intent(out) :: s_top

Top (shallower) / bottom (deeper) edge values of layer k.
real(kind=wp), intent(out) :: s_bot

Top (shallower) / bottom (deeper) edge values of layer k.
real(kind=wp), intent(out) :: t_top

Top (shallower) / bottom (deeper) edge values of layer k.
real(kind=wp), intent(out) :: t_bot

Top (shallower) / bottom (deeper) edge values of layer k.

public pure subroutine roquet_pcm_dpa_face(e_top_l, e_top_r, dz_l, dz_r, t_l, t_r, s_l, s_r, dpa_l, dpa_r, hwt_ll, hwt_lr, hwt_rr, hwt_rl, rho0, rho_ref, dpa_face)

Cross-face 5-point Boole quadrature of the layer dpa for a PCM column pair under Roquet SpV (intx_dpa / inty_dpa): the Roquet twin of wright_pcm_dpa_face, and boole_dpa_face_pcm with the factored vertical rule roquet_pcm_dpa_intz at the three interior sub-columns (same sub-column T/S / height / thickness interpolation, same weights and summation order). 3 T/S polynomials + 15 pressure Horners per face per layer, against 15 full generic-EOS calls.

Arguments

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

Shallower-interface heights in the left / right column (m).
real(kind=wp), intent(in) :: e_top_r

Shallower-interface heights in the left / right column (m).
real(kind=wp), intent(in) :: dz_l

Layer thicknesses in the left / right column (m, >= 0).
real(kind=wp), intent(in) :: dz_r

Layer thicknesses in the left / right column (m, >= 0).
real(kind=wp), intent(in) :: t_l

Layer-mean temperature / salinity in the left / right column.
real(kind=wp), intent(in) :: t_r

Layer-mean temperature / salinity in the left / right column.
real(kind=wp), intent(in) :: s_l

Layer-mean temperature / salinity in the left / right column.
real(kind=wp), intent(in) :: s_r

Layer-mean temperature / salinity in the left / right column.
real(kind=wp), intent(in) :: dpa_l

The columns’ own g * int rho' dz over the layer (Pa).
real(kind=wp), intent(in) :: dpa_r

The columns’ own g * int rho' dz over the layer (Pa).
real(kind=wp), intent(in) :: hwt_ll

MOM6 hWt_LL/LR/RR/RL mass-weighting fractions.
real(kind=wp), intent(in) :: hwt_lr

MOM6 hWt_LL/LR/RR/RL mass-weighting fractions.
real(kind=wp), intent(in) :: hwt_rr

MOM6 hWt_LL/LR/RR/RL mass-weighting fractions.
real(kind=wp), intent(in) :: hwt_rl

MOM6 hWt_LL/LR/RR/RL mass-weighting fractions.
real(kind=wp), intent(in) :: rho0

Boussinesq pressure-estimate density / anomaly reference.
real(kind=wp), intent(in) :: rho_ref

Boussinesq pressure-estimate density / anomaly reference.
real(kind=wp), intent(out) :: dpa_face

Along-face mean of g * int rho' dz over the layer (Pa).

public pure subroutine roquet_pcm_dpa_intz(t, s, e_top, dz, rho0, rho_ref, dpa, intz_dpa)

Vertical integral of the Roquet et al. (2015) SpV in-situ density anomaly over one constant-T/S (PCM) layer: the 5-point Boole rule of boole_dpa_intz_layer (MOM6 int_density_dz_generic_pcm, which is also what MOM6 runs for ROQUET_SPV under a Boussinesq PGF), with the EOS FACTORED.

Read more…

Arguments

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

Layer potential temperature (degC), constant through the layer.
real(kind=wp), intent(in) :: s

Layer practical salinity (PSU), constant through the layer.
real(kind=wp), intent(in) :: e_top

Height of the SHALLOWER interface (m, geopotential).
real(kind=wp), intent(in) :: dz

Layer thickness (m, >= 0); the layer spans [e_top - dz, e_top].
real(kind=wp), intent(in) :: rho0

Boussinesq reference density of the pressure estimate (kg/m^3).
real(kind=wp), intent(in) :: rho_ref

Anomaly reference subtracted from the in-situ density (kg/m^3).
real(kind=wp), intent(out) :: dpa

g * int (rho - rho_ref) dz over the layer (Pa).
real(kind=wp), intent(out) :: intz_dpa

First moment from the top (Pa*m), as boole_dpa_intz_layer.

public pure subroutine roquet_recon_dpa_face(rho0, rho_ref, e_top_l, e_top_r, dz_l, dz_r, t_t_l, t_b_l, t_m_l, t_t_r, t_b_r, t_m_r, s_t_l, s_b_l, s_m_l, s_t_r, s_b_r, s_m_r, dpa_l, dpa_r, parabolic, dpa_face)

boole_dpa_face (the reconstruct-for-pressure cross-face 5-point Boole rule) with the Roquet vertical rule roquet_recon_dpa_intz at the three interior sub-columns – identical sub-columns, weights and summation order; the end points are the columns’ own dpa. 15 inlined Roquet evaluations per face per layer.

Arguments

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

Pressure-estimate density / anomaly reference (kg/m^3).
real(kind=wp), intent(in) :: rho_ref

Pressure-estimate density / anomaly reference (kg/m^3).
real(kind=wp), intent(in) :: e_top_l

Left / right shallower-interface heights and thicknesses (m).
real(kind=wp), intent(in) :: e_top_r

Left / right shallower-interface heights and thicknesses (m).
real(kind=wp), intent(in) :: dz_l

Left / right shallower-interface heights and thicknesses (m).
real(kind=wp), intent(in) :: dz_r

Left / right shallower-interface heights and thicknesses (m).
real(kind=wp), intent(in) :: t_t_l

Left / right temperature triples (top, bottom, mean).
real(kind=wp), intent(in) :: t_b_l

Left / right temperature triples (top, bottom, mean).
real(kind=wp), intent(in) :: t_m_l

Left / right temperature triples (top, bottom, mean).
real(kind=wp), intent(in) :: t_t_r

Left / right temperature triples (top, bottom, mean).
real(kind=wp), intent(in) :: t_b_r

Left / right temperature triples (top, bottom, mean).
real(kind=wp), intent(in) :: t_m_r

Left / right temperature triples (top, bottom, mean).
real(kind=wp), intent(in) :: s_t_l

Left / right salinity triples.
real(kind=wp), intent(in) :: s_b_l

Left / right salinity triples.
real(kind=wp), intent(in) :: s_m_l

Left / right salinity triples.
real(kind=wp), intent(in) :: s_t_r

Left / right salinity triples.
real(kind=wp), intent(in) :: s_b_r

Left / right salinity triples.
real(kind=wp), intent(in) :: s_m_r

Left / right salinity triples.
real(kind=wp), intent(in) :: dpa_l

The left / right columns’ own layer dpa (Pa), the end points.
real(kind=wp), intent(in) :: dpa_r

The left / right columns’ own layer dpa (Pa), the end points.
logical, intent(in) :: parabolic

.true. -> PPM curvature.
real(kind=wp), intent(out) :: dpa_face

Along-face mean of g * int rho' dz over the layer (Pa).

public pure subroutine roquet_recon_dpa_intz(rho0, rho_ref, e_top, dz, t_t, t_b, t_mean, s_t, s_b, s_mean, parabolic, dpa, intz_dpa)

boole_dpa_intz_layer (the reconstruct-for-pressure in-layer 5-point Boole rule over a PLM / PPM T/S profile) specialised to Roquet SpV: the same five sub-points, weights and summation order, with the density 1/SV from this module’s copy of the SpV value (rdb_roq_ts_coeffs + rdb_roq_spv_p, rdb_roquet_spv.inc) – no eos_t handle, no per-point variant dispatch, and a body the compiler inlines into the kernel. T and S vary through the layer, so unlike roquet_pcm_dpa_intz there is no (T, S) hoist: each point is a full EOS evaluation.

Arguments

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

Pressure-estimate density / anomaly reference (kg/m^3).
real(kind=wp), intent(in) :: rho_ref

Pressure-estimate density / anomaly reference (kg/m^3).
real(kind=wp), intent(in) :: e_top

Shallower-interface height / layer thickness (m).
real(kind=wp), intent(in) :: dz

Shallower-interface height / layer thickness (m).
real(kind=wp), intent(in) :: t_t

T / S top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: t_b

T / S top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: t_mean

T / S top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: s_t

T / S top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: s_b

T / S top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: s_mean

T / S top edge, bottom edge, layer mean.
logical, intent(in) :: parabolic

.true. -> PPM curvature.
real(kind=wp), intent(out) :: dpa

As boole_dpa_intz_layer.
real(kind=wp), intent(out) :: intz_dpa

As boole_dpa_intz_layer.

public pure subroutine wright_pcm_dpa_face(e_top_l, e_top_r, dz_l, dz_r, t_l, t_r, s_l, s_r, dpa_l, dpa_r, hwt_ll, hwt_lr, hwt_rr, hwt_rl, rho0, rho_ref, dpa_face)

Cross-face 5-point Boole quadrature of the layer dpa for a PCM column pair with the ANALYTIC Wright vertical integral at every lateral sub-column (MOM6 int_density_dz_wright, intx_dpa / inty_dpa). The Wright twin of boole_dpa_face_pcm, with the same sub-column definition: the end points are the columns’ own dpa_l / dpa_r; the three interior sub-columns interpolate the interface height and thickness linearly in the cross-face fraction and T/S with the near-bottom mass-weighted fractions wtT_L = wl*hwt_ll + wr*hwt_rl, wtT_R = wl*hwt_lr + wr*hwt_rr.

Read more…

Arguments

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

Shallower-interface heights in the left / right column (m).
real(kind=wp), intent(in) :: e_top_r

Shallower-interface heights in the left / right column (m).
real(kind=wp), intent(in) :: dz_l

Layer thicknesses in the left / right column (m, >= 0).
real(kind=wp), intent(in) :: dz_r

Layer thicknesses in the left / right column (m, >= 0).
real(kind=wp), intent(in) :: t_l

Layer-mean temperature / salinity in the left / right column.
real(kind=wp), intent(in) :: t_r

Layer-mean temperature / salinity in the left / right column.
real(kind=wp), intent(in) :: s_l

Layer-mean temperature / salinity in the left / right column.
real(kind=wp), intent(in) :: s_r

Layer-mean temperature / salinity in the left / right column.
real(kind=wp), intent(in) :: dpa_l

The columns’ own g * int rho' dz over the layer (Pa).
real(kind=wp), intent(in) :: dpa_r

The columns’ own g * int rho' dz over the layer (Pa).
real(kind=wp), intent(in) :: hwt_ll

MOM6 hWt_LL/LR/RR/RL mass-weighting fractions.
real(kind=wp), intent(in) :: hwt_lr

MOM6 hWt_LL/LR/RR/RL mass-weighting fractions.
real(kind=wp), intent(in) :: hwt_rr

MOM6 hWt_LL/LR/RR/RL mass-weighting fractions.
real(kind=wp), intent(in) :: hwt_rl

MOM6 hWt_LL/LR/RR/RL mass-weighting fractions.
real(kind=wp), intent(in) :: rho0

Boussinesq pressure-estimate density / anomaly reference.
real(kind=wp), intent(in) :: rho_ref

Boussinesq pressure-estimate density / anomaly reference.
real(kind=wp), intent(out) :: dpa_face

Along-face mean of g * int rho' dz over the layer (Pa).

public pure subroutine wright_pcm_dpa_intz(t, s, e_top, dz, rho0, rho_ref, dpa, intz_dpa)

ANALYTIC vertical integral of the Wright (1997) in-situ density anomaly over one constant-T/S (PCM) layer — MOM6 int_density_dz_wright (reduced-range coefficients, MOM6 EQN_OF_STATE = "WRIGHT" / "WRIGHT_RED", the set rdb_eos carries). Replaces the 5-point Boole quadrature of boole_dpa_intz_layer (5 generic-EOS evaluations) with one polynomial evaluation, one division pair and a short series.

Read more…

Arguments

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

Layer temperature (degC), constant through the layer.
real(kind=wp), intent(in) :: s

Layer salinity (PSU), constant through the layer.
real(kind=wp), intent(in) :: e_top

Height of the SHALLOWER interface (m, geopotential, negative below the datum).
real(kind=wp), intent(in) :: dz

Layer thickness (m, >= 0); the layer spans [e_top - dz, e_top].
real(kind=wp), intent(in) :: rho0

Boussinesq reference density of the pressure estimate (kg/m^3).
real(kind=wp), intent(in) :: rho_ref

Anomaly reference subtracted from the in-situ density (kg/m^3).
real(kind=wp), intent(out) :: dpa

g * int (rho - rho_ref) dz over the layer (Pa).
real(kind=wp), intent(out) :: intz_dpa

First moment from the top (Pa*m).

private pure subroutine boole_layer_combine(r5, dz, dpa, intz_dpa)

Boole weights of the five sub-point density anomalies r5 (top to bottom): dpa = g*dz*<rho'> and the first moment from the top.

Arguments

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

Density anomaly at the five sub-points (kg/m^3).
real(kind=wp), intent(in) :: dz

Layer thickness (m).
real(kind=wp), intent(out) :: dpa

g * int rho’ dz over the layer (Pa).
real(kind=wp), intent(out) :: intz_dpa

0.5 * g * dz^2 * bracket (Pa*m), first moment from the top.

private pure subroutine boole_layer_points(rho0, e_top, dz, t_t, t_b, t_mean, s_t, s_b, s_mean, parabolic, t5, s5, p5)

The five sub-point (T, S, p) triples of the in-layer Boole rule, top (n = 1) to bottom (n = 5), for the per-EOS twins – the same points boole_dpa_intz_layer writes out in place.

Arguments

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

Boussinesq reference density used in the pressure estimate.
real(kind=wp), intent(in) :: e_top

Surface-relative height of the SHALLOWER interface.
real(kind=wp), intent(in) :: dz

Layer thickness (m), dz >= 0.
real(kind=wp), intent(in) :: t_t

Temperature: top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: t_b

Temperature: top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: t_mean

Temperature: top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: s_t

Salinity: top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: s_b

Salinity: top edge, bottom edge, layer mean.
real(kind=wp), intent(in) :: s_mean

Salinity: top edge, bottom edge, layer mean.
logical, intent(in) :: parabolic

.true. -> add the PPM curvature (s6/t6) term.
real(kind=wp), intent(out) :: t5(N_BOOLE)

Sub-point temperature, salinity and Boussinesq pressure (Pa).
real(kind=wp), intent(out) :: s5(N_BOOLE)

Sub-point temperature, salinity and Boussinesq pressure (Pa).
real(kind=wp), intent(out) :: p5(N_BOOLE)

Sub-point temperature, salinity and Boussinesq pressure (Pa).

private pure subroutine boundary_edges_linear(h_self, h_nbr, q_self, dq_up, q_t, q_b)

Linear-exact one-sided edge pair for a BOUNDARY layer (k=1 or k=nz), where a centred slope has no second neighbour.

Read more…

Arguments

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

Thickness of the boundary layer itself (m).
real(kind=wp), intent(in) :: h_nbr

Thickness of its single interior neighbour (m).
real(kind=wp), intent(in) :: q_self

Layer mean of the boundary layer.
real(kind=wp), intent(in) :: dq_up

Layer-mean increment toward the surface, neighbour -> self at the surface layer, self -> neighbour at the bed layer.
real(kind=wp), intent(out) :: q_t

Top (shallower) edge value.
real(kind=wp), intent(out) :: q_b

Bottom (deeper) edge value.

private pure subroutine compute_fv_mom6_impl(h_layer, rho_layer, b, e_face, pa, intz_dpa, intx_pa, inty_pa, intx_dpa, inty_dpa, dpdx_face, dpdy_face, rho0, rho_ref, h_neglect, gfs_scale, mass_weight, p_top, p_top_in_bc, idxCu, idyCv, nx, ny, nz)

Faithful port of MOM6’s PressureForce_FV_Bouss per-layer PGF for the Boussinesq + per-layer Rlay path.

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: h_layer(nx,ny,nz)
real(kind=wp), intent(in) :: rho_layer(nx,ny,nz)
real(kind=wp), intent(in) :: b(nx,ny)
real(kind=wp), intent(inout) :: e_face(nx,ny,nz+1)
real(kind=wp), intent(inout) :: pa(nx,ny,nz+1)
real(kind=wp), intent(inout) :: intz_dpa(nx,ny,nz)
real(kind=wp), intent(inout) :: intx_pa(nx+1,ny,nz+1)
real(kind=wp), intent(inout) :: inty_pa(nx,ny+1,nz+1)
real(kind=wp), intent(inout) :: intx_dpa(nx+1,ny,nz)
real(kind=wp), intent(inout) :: inty_dpa(nx,ny+1,nz)
real(kind=wp), intent(inout) :: dpdx_face(nx+1,ny,nz)
real(kind=wp), intent(inout) :: dpdy_face(nx,ny+1,nz)
real(kind=wp), intent(in) :: rho0
real(kind=wp), intent(in) :: rho_ref
real(kind=wp), intent(in) :: h_neglect
real(kind=wp), intent(in) :: gfs_scale
logical, intent(in) :: mass_weight
real(kind=wp), intent(in) :: p_top(nx,ny)

Top-of-column pressure (Pa, >= 0), multilayer_state_t%p_top. Consulted only when p_top_in_bc; the zero array otherwise.
logical, intent(in) :: p_top_in_bc

Add p_top to the Pass-1 surface BC. .false. ⇒ the assignment is character-for-character the pre-knob one ⇒ bit-identical (same branch-on-a-scalar-knob shape as mass_weight in Pass 2).
real(kind=wp), intent(in) :: idxCu(nx+1,ny)
real(kind=wp), intent(in) :: idyCv(nx,ny+1)
integer, intent(in) :: nx
integer, intent(in) :: ny
integer, intent(in) :: nz

private pure subroutine compute_fv_mom6_insitu_pcm_impl(h_layer, hS, hT, b, conc_T, conc_S, e_face, pa, intz_dpa, intx_pa, inty_pa, intx_dpa, inty_dpa, dpdx_face, dpdy_face, rho0, rho_ref, h_neglect, gfs_scale, mass_weight, eos_variant, p_top, p_top_in_bc, idxCu, idyCv, nx, ny, nz)

FV_MOM6 pressure gradient, constant-by-layer (PCM) T/S, density at the IN-SITU pressure — MOM6 PressureForce_FV_Bouss with RECONSTRUCT_FOR_PRESSURE = False (int_density_dz_generic_pcm).

Read more…

Arguments

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

Salinity * thickness (PSU*m) — layer-mean S = hS / h.
real(kind=wp), intent(in) :: hT(nx,ny,nz)

Temperature * thickness (degC*m) — layer-mean T = hT / h.
real(kind=wp), intent(in) :: b(nx,ny)
real(kind=wp), intent(inout) :: conc_T(nx,ny,nz)

Layer-mean T / S as every pass below reads them (Pass C).
real(kind=wp), intent(inout) :: conc_S(nx,ny,nz)

Layer-mean T / S as every pass below reads them (Pass C).
real(kind=wp), intent(inout) :: e_face(nx,ny,nz+1)
real(kind=wp), intent(inout) :: pa(nx,ny,nz+1)
real(kind=wp), intent(inout) :: intz_dpa(nx,ny,nz)
real(kind=wp), intent(inout) :: intx_pa(nx+1,ny,nz+1)
real(kind=wp), intent(inout) :: inty_pa(nx,ny+1,nz+1)
real(kind=wp), intent(inout) :: intx_dpa(nx+1,ny,nz)
real(kind=wp), intent(inout) :: inty_dpa(nx,ny+1,nz)
real(kind=wp), intent(inout) :: dpdx_face(nx+1,ny,nz)
real(kind=wp), intent(inout) :: dpdy_face(nx,ny+1,nz)
real(kind=wp), intent(in) :: rho0
real(kind=wp), intent(in) :: rho_ref
real(kind=wp), intent(in) :: h_neglect
real(kind=wp), intent(in) :: gfs_scale
logical, intent(in) :: mass_weight

MOM6 MASS_WEIGHT_IN_PRESSURE_GRADIENT (near-bottom hWght).
integer, intent(in) :: eos_variant

eos%variant, selected ONCE here, outside the loops – each variant has its own loop copies with its integrals inlined, no generic eos_t dispatch in the hot loop. use_insitu_pcm admits exactly two: * EOS_VARIANT_WRIGHT_97: the ANALYTIC Wright layer integral (wright_pcm_dpa_intz / wright_pcm_dpa_face, MOM6 int_density_dz_wright). * EOS_VARIANT_ROQUET_SPV (the else copies): the 5-point Boole rule of MOM6 int_density_dz_generic_pcm, with the (T, S) part of the EOS evaluated once per sub-column (roquet_pcm_dpa_intz / roquet_pcm_dpa_face).
real(kind=wp), intent(in) :: p_top(nx,ny)

Top-of-column pressure (Pa, >= 0), multilayer_state_t%p_top.
logical, intent(in) :: p_top_in_bc

Add p_top to the Pass-1 surface BC (.false. ⇒ the plain rho_ref·g·eta seed).
real(kind=wp), intent(in) :: idxCu(nx+1,ny)
real(kind=wp), intent(in) :: idyCv(nx,ny+1)
integer, intent(in) :: nx
integer, intent(in) :: ny
integer, intent(in) :: nz

private pure subroutine compute_fv_mom6_reconstruct_impl(h_layer, hS, hT, b, eos, S_t, S_b, T_t, T_b, conc_T, conc_S, e_face, pa, intz_dpa, intx_pa, inty_pa, intx_dpa, inty_dpa, dpdx_face, dpdy_face, rho0, rho_ref, h_neglect, gfs_scale, recon_scheme, p_top, p_top_in_bc, idxCu, idyCv, nx, ny, nz)

FV_MOM6 pressure-gradient with in-layer T/S reconstruction.

Read more…

Arguments

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

Salinity * thickness (PSU*m) — layer-mean S = hS / h.
real(kind=wp), intent(in) :: hT(nx,ny,nz)

Temperature * thickness (degC*m) — layer-mean T = hT / h.
real(kind=wp), intent(in) :: b(nx,ny)
type(eos_t), intent(in) :: eos
real(kind=wp), intent(inout) :: S_t(nx,ny,nz)
real(kind=wp), intent(inout) :: S_b(nx,ny,nz)
real(kind=wp), intent(inout) :: T_t(nx,ny,nz)
real(kind=wp), intent(inout) :: T_b(nx,ny,nz)
real(kind=wp), intent(inout) :: conc_T(nx,ny,nz)

Layer-mean T / S as every pass below reads them (Pass C).
real(kind=wp), intent(inout) :: conc_S(nx,ny,nz)

Layer-mean T / S as every pass below reads them (Pass C).
real(kind=wp), intent(inout) :: e_face(nx,ny,nz+1)
real(kind=wp), intent(inout) :: pa(nx,ny,nz+1)
real(kind=wp), intent(inout) :: intz_dpa(nx,ny,nz)
real(kind=wp), intent(inout) :: intx_pa(nx+1,ny,nz+1)
real(kind=wp), intent(inout) :: inty_pa(nx,ny+1,nz+1)
real(kind=wp), intent(inout) :: intx_dpa(nx+1,ny,nz)
real(kind=wp), intent(inout) :: inty_dpa(nx,ny+1,nz)
real(kind=wp), intent(inout) :: dpdx_face(nx+1,ny,nz)
real(kind=wp), intent(inout) :: dpdy_face(nx,ny+1,nz)
real(kind=wp), intent(in) :: rho0
real(kind=wp), intent(in) :: rho_ref
real(kind=wp), intent(in) :: h_neglect
real(kind=wp), intent(in) :: gfs_scale
integer, intent(in) :: recon_scheme
real(kind=wp), intent(in) :: p_top(nx,ny)

Top-of-column pressure (Pa, >= 0), multilayer_state_t%p_top.
logical, intent(in) :: p_top_in_bc

Add p_top to the Pass-1 surface BC (.false. ⇒ bit-identical).
real(kind=wp), intent(in) :: idxCu(nx+1,ny)
real(kind=wp), intent(in) :: idyCv(nx,ny+1)
integer, intent(in) :: nx
integer, intent(in) :: ny
integer, intent(in) :: nz

private pure subroutine compute_gprime_impl(h_layer, b, dpdx_face, dpdy_face, gfs, gint, idxCu, idyCv, nx, ny, nz)

Reduced-gravity / gprime PGF for NK = 2.

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: h_layer(nx,ny,nz)
real(kind=wp), intent(in) :: b(nx,ny)
real(kind=wp), intent(inout) :: dpdx_face(nx+1,ny,nz)
real(kind=wp), intent(inout) :: dpdy_face(nx,ny+1,nz)
real(kind=wp), intent(in) :: gfs
real(kind=wp), intent(in) :: gint
real(kind=wp), intent(in) :: idxCu(nx+1,ny)
real(kind=wp), intent(in) :: idyCv(nx,ny+1)
integer, intent(in) :: nx
integer, intent(in) :: ny
integer, intent(in) :: nz

private pure subroutine fv_mom6_mass_weights(mass_weight, e_bed_l, e_bed_r, e_top_l, e_top_r, h_l, h_r, h_neglect, hwt_ll, hwt_lr, hwt_rr, hwt_rl)

MOM6 near-bottom hWght mass-weighting fractions of one face (MASS_WEIGHT_IN_PRESSURE_GRADIENT): the distance by which the layer top of one column sits below the other column’s bed, scaled by the squared relative thickness contrast. Zero ⇒ plain linear interpolation (hwt_ll = hwt_rr = 1, hwt_lr = hwt_rl = 0).

Arguments

Type IntentOptional Attributes Name
logical, intent(in) :: mass_weight

Mass weighting on?
real(kind=wp), intent(in) :: e_bed_l

Bed heights of the left / right column (m).
real(kind=wp), intent(in) :: e_bed_r

Bed heights of the left / right column (m).
real(kind=wp), intent(in) :: e_top_l

Layer top heights of the left / right column (m).
real(kind=wp), intent(in) :: e_top_r

Layer top heights of the left / right column (m).
real(kind=wp), intent(in) :: h_l

Layer thicknesses of the left / right column (m).
real(kind=wp), intent(in) :: h_r

Layer thicknesses of the left / right column (m).
real(kind=wp), intent(in) :: h_neglect

Negligible thickness (m).
real(kind=wp), intent(out) :: hwt_ll

MOM6 hWt_LL/LR/RR/RL.
real(kind=wp), intent(out) :: hwt_lr

MOM6 hWt_LL/LR/RR/RL.
real(kind=wp), intent(out) :: hwt_rr

MOM6 hWt_LL/LR/RR/RL.
real(kind=wp), intent(out) :: hwt_rl

MOM6 hWt_LL/LR/RR/RL.

private subroutine ocean_pressure_force_destroy(this)

Arguments

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

private subroutine ocean_pressure_force_enter_data(this)

Arguments

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

private subroutine ocean_pressure_force_enter_data_impl(this)

Arguments

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

private subroutine ocean_pressure_force_exit_data(this)

Arguments

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

private subroutine ocean_pressure_force_exit_data_impl(this)

Arguments

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

private subroutine ocean_pressure_force_init(this, grid, nz_ml)

Allocate the scratch buffers. Default nz_ml=1 keeps the barotropic-only path constructible; passing nz_ml sizes them for the multilayer kernel.

Read more…

Arguments

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

private subroutine ocean_pressure_force_set_bathymetry(this, b)

Copy b(:, :) into this%b on the host. Must be called BEFORE enter_data (the device copy is taken from the host values). Issues no update device; to refresh post enter_data the caller must issue !$acc update device(this%b) itself.

Arguments

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

private pure subroutine ppm_interface_values(m, nz, h0, h1, h2, h3, s0, s1, s2, s3, t0, t1, t2, t3, edge_s, edge_t)

The PPM estimates of salinity and temperature at the interface between layer m (deeper) and m+1 (shallower), 1 <= m <= nz-1, from the four-layer stencil m-1 .. m+2 (thicknesses h0..h3, layer means s0..s3, t0..t3). Interior interfaces (2 <= m <= nz-2) take the implicit-h4 estimate in its explicit form (White & Adcroft 2008 non-uniform stencil, exactly 4th-order on non-uniform layers; matches remap_column_ppm_h4 Step 1 to round-off). The two near-boundary interfaces (m = 1, m = nz-1) take the thickness-weighted (h2) estimate — the linear-exact value at the shared face of two piecewise-linear cells, (q_m h_{m+1} + q_{m+1} h_m)/(h_m + h_{m+1}) (a plain mean biases the thick interior layer’s edge on non-uniform thicknesses); layer m-1 (resp. m+2) is then not referenced.

Read more…

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: m

Interface index: between layers m and m+1.
integer, intent(in) :: nz

Number of layers in the column.
real(kind=wp), intent(in) :: h0

Thicknesses of layers m-1, m, m+1, m+2 (m).
real(kind=wp), intent(in) :: h1

Thicknesses of layers m-1, m, m+1, m+2 (m).
real(kind=wp), intent(in) :: h2

Thicknesses of layers m-1, m, m+1, m+2 (m).
real(kind=wp), intent(in) :: h3

Thicknesses of layers m-1, m, m+1, m+2 (m).
real(kind=wp), intent(in) :: s0

Salinity layer means of layers m-1 .. m+2.
real(kind=wp), intent(in) :: s1

Salinity layer means of layers m-1 .. m+2.
real(kind=wp), intent(in) :: s2

Salinity layer means of layers m-1 .. m+2.
real(kind=wp), intent(in) :: s3

Salinity layer means of layers m-1 .. m+2.
real(kind=wp), intent(in) :: t0

Temperature layer means of layers m-1 .. m+2.
real(kind=wp), intent(in) :: t1

Temperature layer means of layers m-1 .. m+2.
real(kind=wp), intent(in) :: t2

Temperature layer means of layers m-1 .. m+2.
real(kind=wp), intent(in) :: t3

Temperature layer means of layers m-1 .. m+2.
real(kind=wp), intent(out) :: edge_s

Interface salinity / temperature.
real(kind=wp), intent(out) :: edge_t

Interface salinity / temperature.

private pure subroutine ppm_limit_edges(q_m1, q_c, q_p1, ql_raw, qr_raw, q_t, q_b)

The PPM edge limiter of one interior layer: clip both interface estimates into the monotone bounds of the three adjacent means, flatten a local extremum to PCM, and apply the Colella & Woodward (1984) parabola limiter.

Arguments

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

Layer means of layers k-1, k, k+1.
real(kind=wp), intent(in) :: q_c

Layer means of layers k-1, k, k+1.
real(kind=wp), intent(in) :: q_p1

Layer means of layers k-1, k, k+1.
real(kind=wp), intent(in) :: ql_raw

Interface estimates at the layer’s deeper / shallower interface.
real(kind=wp), intent(in) :: qr_raw

Interface estimates at the layer’s deeper / shallower interface.
real(kind=wp), intent(out) :: q_t

Limited top (shallower) / bottom (deeper) edge values.
real(kind=wp), intent(out) :: q_b

Limited top (shallower) / bottom (deeper) edge values.