mom_neutral_diffusion module reference¶

A column-wise toolbox for implementing neutral diffusion.

Data Types¶

neutral_diffusion_cs

The control structure for the MOM_neutral_diffusion module.

Functions/Subroutines¶

`neutral_diffusion_init()`	Read parameters and allocate control structure for neutral_diffusion module.
`neutral_diffusion_calc_coeffs()`	Calculate remapping factors for u/v columns used to map adjoining columns to a shared coordinate space.
`neutral_diffusion()`	Update tracer concentration due to neutral diffusion; layer thickness unchanged by this update.
`compute_tapering_coeffs()`	Computes linear tapering coefficients at interfaces of the left and right columns within a region defined by the boundary layer depths in the two columns.
`interface_scalar()`	Returns interface scalar, Si, for a column of layer values, S.
`ppm_edge()`	Returns the PPM quasi-fourth order edge value at k+1/2 following equation 1.6 in Colella & Woodward, 1984: JCP 54, 174-201.
`ppm_ave()`	Returns the average of a PPM reconstruction between two fractional positions in the same arbitrary concentration units as aMean (e.g.
`signum()`	A true signum function that returns either -abs(a), when x<0; or abs(a) when x>0; or 0 when x=0.
`plm_diff()`	Returns PLM slopes for a column where the slopes are the difference in value across each cell.
`fv_diff()`	Returns the cell-centered second-order finite volume (unlimited PLM) slope using three consecutive cell widths and average values.
`fvlsq_slope()`	Returns the cell-centered second-order weighted least squares slope using three consecutive cell widths and average values.
`find_neutral_surface_positions_continuous()`	Returns positions within left/right columns of combined interfaces using continuous reconstructions of T/S.
`interpolate_for_nondim_position()`	Returns the non-dimensional position between Pneg and Ppos where the interpolated density difference equals zero.
`find_neutral_surface_positions_discontinuous()`	Higher order version of find_neutral_surface_positions.
`mark_unstable_cells()`	Sweep down through the column and mark as stable if the bottom interface of a cell is denser than the top.
`search_other_column()`	Searches the “other” (searched) column for the position of the neutral surface.
`increment_interface()`	Increments the interface which was just connected and also set flags if the bottom is reached.
`find_neutral_pos_linear()`	Search a layer to find where delta_rho = 0 based on a linear interpolation of alpha and beta of the top and bottom being searched and polynomial reconstructions of T and S.
`find_neutral_pos_full()`	Use the full equation of state to calculate the difference in locally referenced potential density.
`calc_delta_rho_and_derivs()`	Calculate the difference in density between two points in a variety of ways.
`delta_rho_from_derivs()`	Calculate delta rho from derivatives and gradients of properties \(\Delta \rho = \frac{1}{2}\left[ (\alpha_1 + \alpha_2)(T_1-T_2) + (\beta_1 + \beta_2)(S_1-S_2) + (\gamma^{-1}_1 + \gamma^{-1}_2)*(P_1-P_2) \right]\).
`absolute_position()`	Converts non-dimensional position within a layer to absolute position (for debugging)
`absolute_positions()`	Converts non-dimensional positions within layers to absolute positions (for debugging)
`neutral_surface_flux()`	Returns a single column of neutral diffusion fluxes of a tracer.
`neutral_surface_t_eval()`	Evaluate various parts of the reconstructions to calculate gradient-based flux limiter.
`ppm_left_right_edge_values()`	Discontinuous PPM reconstructions of the left/right edge values within a cell.
`neutral_diffusion_unit_tests()`	Returns true if unit tests of neutral_diffusion functions fail.
`ndiff_unit_tests_continuous()`	Returns true if unit tests of neutral_diffusion functions fail.
`ndiff_unit_tests_discontinuous()`
`test_fv_diff()`	Returns true if a test of `fv_diff()` fails, and conditionally writes results to stream.
`test_fvlsq_slope()`	Returns true if a test of `fvlsq_slope()` fails, and conditionally writes results to stream.
`test_ifndp()`	Returns true if a test of `interpolate_for_nondim_position()` fails, and conditionally writes results to stream.
`test_data1d()`	Returns true if comparison of Po and Ptrue fails, and conditionally writes results to stream.
`test_data1di()`	Returns true if comparison of Po and Ptrue fails, and conditionally writes results to stream.
`test_nsp()`	Returns true if output of find_neutral_surface_positions() does not match correct values, and conditionally writes results to stream.
`compare_nsp_row()`	Compares a single row, k, of output from find_neutral_surface_positions()
`test_rnp()`	Compares output position from refine_nondim_position with an expected value.
`neutral_diffusion_end()`	Deallocates neutral_diffusion control structure.

Detailed Description¶

A column-wise toolbox for implementing neutral diffusion.

Type Documentation¶

type mom_neutral_diffusion/neutral_diffusion_cs¶

The control structure for the MOM_neutral_diffusion module.

Type fields:

% nkp1 [integer] :: Number of interfaces for a column = nk + 1.
% nsurf [integer] :: Number of neutral surfaces.
% deg [integer] :: Degree of polynomial used for reconstructions.
% continuous_reconstruction [logical] :: True if using continuous PPM reconstruction at interfaces.
% debug [logical] :: If true, write verbose debugging messages.
% hard_fail_heff [logical] :: Bring down the model if a problem with heff is detected.
% max_iter [integer] :: Maximum number of iterations if refine_position is defined.
% drho_tol [real] :: Convergence criterion representing density difference from true neutrality [R ~> kg m-3].
% x_tol [real] :: Convergence criterion for how small an update of the position can be [nondim].
% ref_pres [real] :: Reference pressure, negative if using locally referenced neutral density [R L2 T-2 ~> Pa].
% interior_only [logical] :: If true, only applies neutral diffusion in the ocean interior. That is, the algorithm will exclude the surface and bottom boundary layers.
% tapering [logical] :: If true, neutral diffusion linearly decays towards zero within a transition zone defined using boundary layer depths. Only available when interior_only=true.
% khth_use_ebt_struct [logical] :: If true, uses the equivalent barotropic structure as the vertical structure of tracer diffusivity.
% use_unmasked_transport_bug [logical] :: If true, use an older form for the accumulation of neutral-diffusion transports that were unmasked, as used prior to Jan 2018.
% hbl [real(:,:),allocatable] :: Boundary layer depth [H ~> m or kg m-2].
% coeff_l [real(:),allocatable] :: Non-dimensional coefficient in the left column, at cell interfaces [nondim].
% coeff_r [real(:),allocatable] :: Non-dimensional coefficient in the right column, at cell interfaces [nondim].
% coef_h [real(:,:,:),allocatable] :: Coef_x and Coef_y averaged at t-points [L2 ~> m2].
% upol [real(:,:,:),allocatable] :: Non-dimensional position with left layer uKoL-1, u-point [nondim].
% upor [real(:,:,:),allocatable] :: Non-dimensional position with right layer uKoR-1, u-point [nondim].
% ukol [integer(:,:,:),allocatable] :: Index of left interface corresponding to neutral surface, at a u-point.
% ukor [integer(:,:,:),allocatable] :: Index of right interface corresponding to neutral surface, at a u-point.
% uheff [real(:,:,:),allocatable] :: Effective thickness at u-point [H ~> m or kg m-2].
% vpol [real(:,:,:),allocatable] :: Non-dimensional position with left layer uKoL-1, v-point [nondim].
% vpor [real(:,:,:),allocatable] :: Non-dimensional position with right layer uKoR-1, v-point [nondim].
% vkol [integer(:,:,:),allocatable] :: Index of left interface corresponding to neutral surface, at a v-point.
% vkor [integer(:,:,:),allocatable] :: Index of right interface corresponding to neutral surface, at a v-point.
% vheff [real(:,:,:),allocatable] :: Effective thickness at v-point [H ~> m or kg m-2].
% ppoly_coeffs_t [real(:,:,:,:),allocatable] :: Polynomial coefficients of the sub-gridscale temperatures [C ~> degC].
% ppoly_coeffs_s [real(:,:,:,:),allocatable] :: Polynomial coefficients of the sub-gridscale salinity [S ~> ppt].
% drdt [real(:,:,:),allocatable] :: dRho/dT [R C-1 ~> kg m-3 degC-1] at interfaces
% drds [real(:,:,:),allocatable] :: dRho/dS [R S-1 ~> kg m-3 ppt-1] at interfaces
% tint [real(:,:,:),allocatable] :: Interface T [C ~> degC].
% sint [real(:,:,:),allocatable] :: Interface S [S ~> ppt].
% pint [real(:,:,:),allocatable] :: Interface pressure [R L2 T-2 ~> Pa].
% t_i [real(:,:,:,:),allocatable] :: Top edge reconstruction of temperature [C ~> degC].
% s_i [real(:,:,:,:),allocatable] :: Top edge reconstruction of salinity [S ~> ppt].
% p_i [real(:,:,:,:),allocatable] :: Interface pressures [R L2 T-2 ~> Pa].
% drdt_i [real(:,:,:,:),allocatable] :: dRho/dT [R C-1 ~> kg m-3 degC-1] at top edge
% drds_i [real(:,:,:,:),allocatable] :: dRho/dS [R S-1 ~> kg m-3 ppt-1] at top edge
% ns [integer(:,:),allocatable] :: Number of interfaces in a column.
% stable_cell [logical(:,:,:),allocatable] :: True if the cell is stably stratified wrt to the next cell.
% r_to_kg_m3 [real] :: A rescaling factor translating density to kg m-3 for use in diagnostic messages [kg m-3 R-1 ~> 1].
% diag [type( diag_ctrl ),pointer] :: A structure that is used to regulate the timing of diagnostic output.
% neutral_pos_method [integer] :: Method to find the position of a neutral surface within the layer.
% delta_rho_form [character (len=40)] :: Determine which (if any) approximation is made to the equation describing the difference in density.
% id_uheff_2d [integer] :: Diagnostic IDs.
% id_vheff_2d [integer] :: Diagnostic IDs.
% eos [type( eos_type ),pointer] :: Equation of state parameters.
% remap_cs [type( remapping_cs )] :: Remapping control structure used to create sublayers.
% remap_answer_date [integer] :: The vintage of the order of arithmetic and expressions to use for remapping. Values below 20190101 recover the remapping answers from 2018, while higher values use more robust forms of the same remapping expressions.
% ndiff_answer_date [integer] :: The vintage of the order of arithmetic to use for the neutral diffusion. Values of 20240330 or below recover the answers from the original form of this code, while higher values use mathematically equivalent expressions that recover rotational symmetry.
% kpp_csp [type( kpp_cs ),pointer] :: KPP control structure needed to get BLD.
% energetic_pbl_csp [type( energetic_pbl_cs ),pointer] :: ePBL control structure needed to get MLD

Function/Subroutine Documentation¶

function mom_neutral_diffusion/neutral_diffusion_init(Time, G, GV, US, param_file, diag, EOS, diabatic_CSp, CS) [logical]¶

Read parameters and allocate control structure for neutral_diffusion module.

Parameters:

time :: [in] Time structure
g :: [in] Grid structure
gv :: [in] The ocean’s vertical grid structure
us :: [in] A dimensional unit scaling type
diag :: [inout] Diagnostics control structure
param_file :: [in] Parameter file structure
eos :: [in] Equation of state
diabatic_csp :: KPP control structure needed to get BLD
cs :: Neutral diffusion control structure

Call to:

mom_diabatic_driver::extract_diabatic_member mdl mom_error_handler::mom_error mom_remapping::remappingdefaultscheme mom_remapping::remappingschemesdoc

subroutine mom_neutral_diffusion/neutral_diffusion_calc_coeffs(G, GV, US, h, T, S, CS, p_surf)¶

Calculate remapping factors for u/v columns used to map adjoining columns to a shared coordinate space.

Parameters:

g :: [in] Ocean grid structure
gv :: [in] ocean vertical grid structure
us :: [in] A dimensional unit scaling type
h :: [in] Layer thickness [H ~> m or kg m-2]
t :: [in] Potential temperature [C ~> degC]
s :: [in] Salinity [S ~> ppt]
cs :: Neutral diffusion control structure
p_surf :: [in] Surface pressure to include in pressures used for equation of state calculations [R L2 T-2 ~> Pa]

Call to:

mom_hor_bnd_diffusion::boundary_k_range mom_energetic_pbl::energetic_pbl_get_mld polynomial_functions::evaluation_polynomial find_neutral_surface_positions_continuous find_neutral_surface_positions_discontinuous interface_scalar mom_cvmix_kpp::kpp_get_bld mark_unstable_cells mom_hor_bnd_diffusion::surface

Called from:

mom_tracer_hor_diff::tracer_hordiff

subroutine mom_neutral_diffusion/neutral_diffusion(G, GV, h, Coef_x, Coef_y, dt, Reg, US, CS)¶

Update tracer concentration due to neutral diffusion; layer thickness unchanged by this update.

Parameters:

g :: [in] Ocean grid structure
gv :: [in] ocean vertical grid structure
h :: [in] Layer thickness [H ~> m or kg m-2]
coef_x :: [in] dt * Kh * dy / dx at u-points [L2 ~> m2]
coef_y :: [in] dt * Kh * dx / dy at v-points [L2 ~> m2]
dt :: [in] Tracer time step * I_numitts [T ~> s] (I_numitts is in tracer_hordiff)
reg :: Tracer registry
us :: [in] A dimensional unit scaling type
cs :: Neutral diffusion control structure

Call to:

compute_tapering_coeffs neutral_surface_flux

Called from:

mom_tracer_hor_diff::tracer_hordiff

subroutine mom_neutral_diffusion/compute_tapering_coeffs(ne, bld_l, bld_r, coeff_l, coeff_r, h_l, h_r)¶

Computes linear tapering coefficients at interfaces of the left and right columns within a region defined by the boundary layer depths in the two columns.

Parameters:

ne :: [in] Number of interfaces
bld_l :: [in] Boundary layer depth, left column [H ~> m or kg m-2]
bld_r :: [in] Boundary layer depth, right column [H ~> m or kg m-2]
h_l :: [in] Layer thickness, left column [H ~> m or kg m-2]
h_r :: [in] Layer thickness, right column [H ~> m or kg m-2]
coeff_l :: [inout] Tapering coefficient, left column [nondim]
coeff_r :: [inout] Tapering coefficient, right column [nondim]

Call to:

mom_hor_bnd_diffusion::boundary_k_range mom_hor_bnd_diffusion::surface

Called from:

neutral_diffusion

subroutine mom_neutral_diffusion/interface_scalar(nk, h, S, Si, i_method, h_neglect)¶

Returns interface scalar, Si, for a column of layer values, S.

Parameters:

nk :: [in] Number of levels
h :: [in] Layer thickness [H ~> m or kg m-2]
s :: [in] Layer scalar (or concentrations) in arbitrary concentration units (e.g. [C ~> degC] for temperature)
si :: [inout] Interface scalar (or concentrations) in arbitrary concentration units (e.g. [C ~> degC] for temperature)
i_method :: [in] =1 use average of PLM edges =2 use continuous PPM edge interpolation
h_neglect :: [in] A negligibly small thickness [H ~> m or kg m-2]

Call to:

plm_diff ppm_edge

Called from:

ndiff_unit_tests_continuous neutral_diffusion_calc_coeffs neutral_surface_flux

function mom_neutral_diffusion/ppm_edge(hkm1, hk, hkp1, hkp2, Ak, Akp1, Pk, Pkp1, h_neglect) [real]¶

Returns the PPM quasi-fourth order edge value at k+1/2 following equation 1.6 in Colella & Woodward, 1984: JCP 54, 174-201.

Parameters:

hkm1 :: [in] Width of cell k-1 in [H ~> m or kg m-2] or other units
hk :: [in] Width of cell k in [H ~> m or kg m-2] or other units
hkp1 :: [in] Width of cell k+1 in [H ~> m or kg m-2] or other units
hkp2 :: [in] Width of cell k+2 in [H ~> m or kg m-2] or other units
ak :: [in] Average scalar value of cell k in arbitrary concentration units (e.g. [C ~> degC] for temperature)
akp1 :: [in] Average scalar value of cell k+1 in arbitrary concentration units (e.g. [C ~> degC] for temperature)
pk :: [in] PLM slope for cell k in arbitrary concentration units (e.g. [C ~> degC] for temperature)
pkp1 :: [in] PLM slope for cell k+1 in arbitrary concentration units (e.g. [C ~> degC] for temperature)
h_neglect :: [in] A negligibly small thickness [H ~> m or kg m-2]

Called from:

interface_scalar

function mom_neutral_diffusion/ppm_ave(xL, xR, aL, aR, aMean) [real]¶

Returns the average of a PPM reconstruction between two fractional positions in the same arbitrary concentration units as aMean (e.g. usually [C ~> degC] for temperature)

Parameters:

xl :: [in] Fraction position of left bound (0,1) [nondim]
xr :: [in] Fraction position of right bound (0,1) [nondim]
al :: [in] Left edge scalar value, at x=0, in arbitrary concentration units (e.g. usually [C ~> degC] for temperature)
ar :: [in] Right edge scalar value, at x=1 in arbitrary concentration units (e.g. usually [C ~> degC] for temperature)
amean :: [in] Average scalar value of cell in arbitrary concentration units (e.g. usually [C ~> degC] for temperature)

Called from:

neutral_surface_flux

function mom_neutral_diffusion/signum(a, x) [real]¶

A true signum function that returns either -abs(a), when x<0; or abs(a) when x>0; or 0 when x=0. The returned units are the same as those of a [arbitrary].

Parameters:

a :: [in] The magnitude argument in arbitrary units [arbitrary]
x :: [in] The sign (or zero) argument [arbitrary]

Called from:

neutral_surface_flux plm_diff ppm_left_right_edge_values

subroutine mom_neutral_diffusion/plm_diff(nk, h, S, c_method, b_method, diff)¶

Returns PLM slopes for a column where the slopes are the difference in value across each cell. The limiting follows equation 1.8 in Colella & Woodward, 1984: JCP 54, 174-201.

Parameters:

nk :: [in] Number of levels
h :: [in] Layer thickness [H ~> m or kg m-2] or other units
s :: [in] Layer salinity (conc, e.g. ppt) or other tracer concentration in arbitrary units [A ~> a]
c_method :: [in] Method to use for the centered difference
b_method :: [in] =1, use PCM in first/last cell, =2 uses linear extrapolation
diff :: [inout] Scalar difference across layer (conc, e.g. ppt) in the same arbitrary units as S [A ~> a], determined by the following values for c_method:

Second order finite difference (not recommended)
Second order finite volume (used in original PPM)
Finite-volume weighted least squares linear fit

Todo

The use of c_method to choose a scheme is inefficient and should eventually be moved up the call tree.

Call to:: fv_diff fvlsq_slope signum
Called from:: interface_scalar

function mom_neutral_diffusion/fv_diff(hkm1, hk, hkp1, Skm1, Sk, Skp1) [real]¶

Returns the cell-centered second-order finite volume (unlimited PLM) slope using three consecutive cell widths and average values. Slope is returned as a difference across the central cell (i.e. units of scalar S). Discretization follows equation 1.7 in Colella & Woodward, 1984: JCP 54, 174-201.

Parameters:

hkm1 :: [in] Left cell width [H ~> m or kg m-2] or other arbitrary units
hk :: [in] Center cell width [H ~> m or kg m-2] or other arbitrary units
hkp1 :: [in] Right cell width [H ~> m or kg m-2] or other arbitrary units
skm1 :: [in] Left cell average value in arbitrary concentration units (e.g. [C ~> degC] for temperature)
sk :: [in] Center cell average value in arbitrary concentration units (e.g. [C ~> degC] for temperature)
skp1 :: [in] Right cell average value in arbitrary concentration units (e.g. [C ~> degC] for temperature)

Called from:

plm_diff test_fv_diff

function mom_neutral_diffusion/fvlsq_slope(hkm1, hk, hkp1, Skm1, Sk, Skp1) [real]¶

Returns the cell-centered second-order weighted least squares slope using three consecutive cell widths and average values. Slope is returned as a gradient (i.e. units of scalar S over width units). For example, for temperature fvlsq_slope would usually be returned in units of [C H-1 ~> degC m-1 or degC m2 kg-1].

Parameters:

hkm1 :: [in] Left cell width [H ~> m or kg m-2] or other arbitrary units
hk :: [in] Center cell width [H ~> m or kg m-2] or other arbitrary units
hkp1 :: [in] Right cell width [H ~> m or kg m-2] or other arbitrary units
skm1 :: [in] Left cell average value in arbitrary concentration units (e.g. [C ~> degC] for temperature)
sk :: [in] Center cell average value often in arbitrary concentration units (e.g. [C ~> degC] for temperature)
skp1 :: [in] Right cell average value often in arbitrary concentration units (e.g. [C ~> degC] for temperature)

Called from:

plm_diff test_fvlsq_slope

subroutine mom_neutral_diffusion/find_neutral_surface_positions_continuous(nk, Pl, Tl, Sl, dRdTl, dRdSl, Pr, Tr, Sr, dRdTr, dRdSr, PoL, PoR, KoL, KoR, hEff, bl_kl, bl_kr, bl_zl, bl_zr)¶

Returns positions within left/right columns of combined interfaces using continuous reconstructions of T/S.

Parameters:

nk :: [in] Number of levels
pl :: [in] Left-column interface pressure [R L2 T-2 ~> Pa] or other units
tl :: [in] Left-column interface potential temperature [C ~> degC]
sl :: [in] Left-column interface salinity [S ~> ppt]
drdtl :: [in] Left-column dRho/dT [R C-1 ~> kg m-3 degC-1]
drdsl :: [in] Left-column dRho/dS [R S-1 ~> kg m-3 ppt-1]
pr :: [in] Right-column interface pressure [R L2 T-2 ~> Pa] or other units
tr :: [in] Right-column interface potential temperature [C ~> degC]
sr :: [in] Right-column interface salinity [S ~> ppt]
drdtr :: [in] Left-column dRho/dT [R C-1 ~> kg m-3 degC-1]
drdsr :: [in] Left-column dRho/dS [R S-1 ~> kg m-3 ppt-1]
pol :: [inout] Fractional position of neutral surface within layer KoL of left column [nondim]
por :: [inout] Fractional position of neutral surface within layer KoR of right column [nondim]
kol :: [inout] Index of first left interface above neutral surface
kor :: [inout] Index of first right interface above neutral surface
heff :: [inout] Effective thickness between two neutral surfaces [R L2 T-2 ~> Pa] or other units following Pl and Pr.
bl_kl :: [in] Layer index of the boundary layer (left)
bl_kr :: [in] Layer index of the boundary layer (right)
bl_zl :: [in] Fractional position of the boundary layer (left) [nondim]
bl_zr :: [in] Fractional position of the boundary layer (right) [nondim]

Call to:

absolute_position interpolate_for_nondim_position

Called from:

ndiff_unit_tests_continuous neutral_diffusion_calc_coeffs

function mom_neutral_diffusion/interpolate_for_nondim_position(dRhoNeg, Pneg, dRhoPos, Ppos) [real]¶

Returns the non-dimensional position between Pneg and Ppos where the interpolated density difference equals zero. The result is always bounded to be between 0 and 1.

Parameters:

drhoneg :: [in] Negative density difference [R ~> kg m-3]
pneg :: [in] Position of negative density difference [R L2 T-2 ~> Pa] or [nondim]
drhopos :: [in] Positive density difference [R ~> kg m-3]
ppos :: [in] Position of positive density difference [R L2 T-2 ~> Pa] or [nondim]

Call to:

mom_error_handler::mom_error mom_io::stderr

Called from:

find_neutral_surface_positions_continuous search_other_column test_ifndp

subroutine mom_neutral_diffusion/find_neutral_surface_positions_discontinuous(CS, nk, Pres_l, hcol_l, Tl, Sl, ppoly_T_l, ppoly_S_l, stable_l, Pres_r, hcol_r, Tr, Sr, ppoly_T_r, ppoly_S_r, stable_r, PoL, PoR, KoL, KoR, hEff, zeta_bot_L, zeta_bot_R, k_bot_L, k_bot_R, hard_fail_heff)¶

Higher order version of find_neutral_surface_positions. Returns positions within left/right columns of combined interfaces using intracell reconstructions of T/S. Note that the polynomial reconstructions of T and S are optional to aid with unit testing, but will always be passed otherwise.

Parameters:

cs :: [inout] Neutral diffusion control structure
nk :: [in] Number of levels
pres_l :: [in] Left-column interface pressure [R L2 T-2 ~> Pa]
hcol_l :: [in] Left-column layer thicknesses [H ~> m or kg m-2] or other units
tl :: [in] Left-column top interface potential temperature [C ~> degC]
sl :: [in] Left-column top interface salinity [S ~> ppt]
ppoly_t_l :: [in] Left-column coefficients of T reconstruction [C ~> degC]
ppoly_s_l :: [in] Left-column coefficients of S reconstruction [S ~> ppt]
stable_l :: [in] True where the left-column is stable
pres_r :: [in] Right-column interface pressure [R L2 T-2 ~> Pa]
hcol_r :: [in] Left-column layer thicknesses [H ~> m or kg m-2] or other units
tr :: [in] Right-column top interface potential temperature [C ~> degC]
sr :: [in] Right-column top interface salinity [S ~> ppt]
ppoly_t_r :: [in] Right-column coefficients of T reconstruction [C ~> degC]
ppoly_s_r :: [in] Right-column coefficients of S reconstruction [S ~> ppt]
stable_r :: [in] True where the right-column is stable
pol :: [inout] Fractional position of neutral surface within layer KoL of left column [nondim]
por :: [inout] Fractional position of neutral surface within layer KoR of right column [nondim]
kol :: [inout] Index of first left interface above neutral surface
kor :: [inout] Index of first right interface above neutral surface
heff :: [inout] Effective thickness between two neutral surfaces [H ~> m or kg m-2] or other units taken from hcol_l
zeta_bot_l :: [in] Non-dimensional distance to where the boundary layer intersects the cell (left) [nondim]
zeta_bot_r :: [in] Non-dimensional distance to where the boundary layer intersects the cell (right) [nondim]
k_bot_l :: [in] k-index for the boundary layer (left) [nondim]
k_bot_r :: [in] k-index for the boundary layer (right) [nondim]
hard_fail_heff :: [in] If true (default) bring down the model if the neutral surfaces ever cross

Call to:

calc_delta_rho_and_derivs increment_interface mom_error_handler::mom_error search_other_column mom_io::stdout

Called from:

ndiff_unit_tests_discontinuous neutral_diffusion_calc_coeffs

subroutine mom_neutral_diffusion/mark_unstable_cells(CS, nk, T, S, P, stable_cell)¶

Sweep down through the column and mark as stable if the bottom interface of a cell is denser than the top.

Parameters:

cs :: [inout] Neutral diffusion control structure
nk :: [in] Number of levels in a column
t :: [in] Temperature at interfaces [C ~> degC]
s :: [in] Salinity at interfaces [S ~> ppt]
p :: [in] Pressure at interfaces [R L2 T-2 ~> Pa]
stable_cell :: [out] True if this cell is unstably stratified

Call to:

calc_delta_rho_and_derivs

Called from:

ndiff_unit_tests_discontinuous neutral_diffusion_calc_coeffs

function mom_neutral_diffusion/search_other_column(CS, ksurf, pos_last, T_from, S_from, P_from, T_top, S_top, P_top, T_bot, S_bot, P_bot, T_poly, S_poly) [real]¶

Searches the “other” (searched) column for the position of the neutral surface.

Parameters:

cs :: [in] Neutral diffusion control structure
ksurf :: [in] Current index of neutral surface
pos_last :: [in] Last position within the current layer, used as the lower bound in the root finding algorithm [nondim]
t_from :: [in] Temperature at the searched from interface [C ~> degC]
s_from :: [in] Salinity at the searched from interface [S ~> ppt]
p_from :: [in] Pressure at the searched from interface [R L2 T-2 ~> Pa]
t_top :: [in] Temperature at the searched to top interface [C ~> degC]
s_top :: [in] Salinity at the searched to top interface [S ~> ppt]
p_top :: [in] Pressure at the searched to top interface [R L2 T-2 ~> Pa] interface [R L2 T-2 ~> Pa]
t_bot :: [in] Temperature at the searched to bottom interface [C ~> degC]
s_bot :: [in] Salinity at the searched to bottom interface [S ~> ppt]
p_bot :: [in] Pressure at the searched to bottom interface [R L2 T-2 ~> Pa]
t_poly :: [in] Temperature polynomial reconstruction coefficients [C ~> degC]
s_poly :: [in] Salinity polynomial reconstruction coefficients [S ~> ppt]

Call to:

calc_delta_rho_and_derivs find_neutral_pos_full find_neutral_pos_linear interpolate_for_nondim_position

Called from:

find_neutral_surface_positions_discontinuous

subroutine mom_neutral_diffusion/increment_interface(nk, kl, ki, reached_bottom, searching_this_column, searching_other_column)¶

Increments the interface which was just connected and also set flags if the bottom is reached.

Parameters:

nk :: [in] Number of vertical levels
kl :: [inout] Current layer (potentially updated)
ki :: [inout] Current interface
reached_bottom :: [inout] Updated when kl == nk and ki == 2
searching_this_column :: [inout] Updated when kl == nk and ki == 2
searching_other_column :: [inout] Updated when kl == nk and ki == 2

Call to:

mom_error_handler::mom_error

Called from:

find_neutral_surface_positions_discontinuous

function mom_neutral_diffusion/find_neutral_pos_linear(CS, z0, T_ref, S_ref, dRdT_ref, dRdS_ref, dRdT_top, dRdS_top, dRdT_bot, dRdS_bot, ppoly_T, ppoly_S) [real]¶

Search a layer to find where delta_rho = 0 based on a linear interpolation of alpha and beta of the top and bottom being searched and polynomial reconstructions of T and S. Compressibility is not needed because either, we are assuming incompressibility in the equation of state for this module or alpha and beta are calculated having been displaced to the average pressures of the two pressures We need Newton’s method because the T and S reconstructions make delta_rho a polynomial function of z if using PPM or higher. If Newton’s method would search fall out of the interval [0,1], a bisection step would be taken instead. Also this linearization of alpha, beta means that second derivatives of the EOS are not needed. Note that delta in variable names below refers to horizontal differences and ‘d’ refers to vertical differences.

Parameters:

cs :: [in] Control structure with parameters for this module
z0 :: [in] Lower bound of position, also serves as the initial guess [nondim]
t_ref :: [in] Temperature at the searched from interface [C ~> degC]
s_ref :: [in] Salinity at the searched from interface [S ~> ppt]
drdt_ref :: [in] dRho/dT at the searched from interface [R C-1 ~> kg m-3 degC-1]
drds_ref :: [in] dRho/dS at the searched from interface [R S-1 ~> kg m-3 ppt-1]
drdt_top :: [in] dRho/dT at top of layer being searched [R C-1 ~> kg m-3 degC-1]
drds_top :: [in] dRho/dS at top of layer being searched [R S-1 ~> kg m-3 ppt-1]
drdt_bot :: [in] dRho/dT at bottom of layer being searched [R C-1 ~> kg m-3 degC-1]
drds_bot :: [in] dRho/dS at bottom of layer being searched [R S-1 ~> kg m-3 ppt-1]
ppoly_t :: [in] Coefficients of the polynomial reconstruction of T within the layer to be searched [C ~> degC].
ppoly_s :: [in] Coefficients of the polynomial reconstruction of S within the layer to be searched [S ~> ppt].

Return:

undefined :: Position where drho = 0 [nondim]

Call to:

polynomial_functions::evaluation_polynomial polynomial_functions::first_derivative_polynomial mom_error_handler::mom_error

Called from:

ndiff_unit_tests_discontinuous search_other_column

function mom_neutral_diffusion/find_neutral_pos_full(CS, z0, T_ref, S_ref, P_ref, P_top, P_bot, ppoly_T, ppoly_S) [real]¶

Use the full equation of state to calculate the difference in locally referenced potential density. The derivatives in this case are not trivial to calculate, so instead we use a regula falsi method.

Parameters:

cs :: [in] Control structure with parameters for this module
z0 :: [in] Lower bound of position, also serves as the initial guess [nondim]
t_ref :: [in] Temperature at the searched from interface [C ~> degC]
s_ref :: [in] Salinity at the searched from interface [S ~> ppt]
p_ref :: [in] Pressure at the searched from interface [R L2 T-2 ~> Pa]
p_top :: [in] Pressure at top of layer being searched [R L2 T-2 ~> Pa]
p_bot :: [in] Pressure at bottom of layer being searched [R L2 T-2 ~> Pa]
ppoly_t :: [in] Coefficients of the polynomial reconstruction of T within the layer to be searched [C ~> degC]
ppoly_s :: [in] Coefficients of the polynomial reconstruction of T within the layer to be searched [S ~> ppt]

Return:

undefined :: Position where drho = 0 [nondim]

Call to:

calc_delta_rho_and_derivs polynomial_functions::evaluation_polynomial

Called from:

search_other_column

subroutine mom_neutral_diffusion/calc_delta_rho_and_derivs(CS, T1, S1, p1_in, T2, S2, p2_in, drho, drdt1_out, drds1_out, drdt2_out, drds2_out)¶

Calculate the difference in density between two points in a variety of ways.

Parameters:

cs :: Neutral diffusion control structure
t1 :: [in] Temperature at point 1 [C ~> degC]
s1 :: [in] Salinity at point 1 [S ~> ppt]
p1_in :: [in] Pressure at point 1 [R L2 T-2 ~> Pa]
t2 :: [in] Temperature at point 2 [C ~> degC]
s2 :: [in] Salinity at point 2 [S ~> ppt]
p2_in :: [in] Pressure at point 2 [R L2 T-2 ~> Pa]
drho :: [out] Difference in density between the two points [R ~> kg m-3]
drdt1_out :: [out] drho_dt at point 1 [R C-1 ~> kg m-3 degC-1]
drds1_out :: [out] drho_ds at point 1 [R S-1 ~> kg m-3 ppt-1]
drdt2_out :: [out] drho_dt at point 2 [R C-1 ~> kg m-3 degC-1]
drds2_out :: [out] drho_ds at point 2 [R S-1 ~> kg m-3 ppt-1]

Call to:

delta_rho_from_derivs mom_error_handler::mom_error

Called from:

find_neutral_pos_full find_neutral_surface_positions_discontinuous mark_unstable_cells search_other_column

function mom_neutral_diffusion/delta_rho_from_derivs(T1, S1, P1, dRdT1, dRdS1, T2, S2, P2, dRdT2, dRdS2) [real]¶

Calculate delta rho from derivatives and gradients of properties \(\Delta \rho = \frac{1}{2}\left[ (\alpha_1 + \alpha_2)*(T_1-T_2) + (\beta_1 + \beta_2)*(S_1-S_2) + (\gamma^{-1}_1 + \gamma^{-1}_2)*(P_1-P_2) \right]\).

Parameters:

t1 :: Temperature at point 1 [C ~> degC]
s1 :: Salinity at point 1 [S ~> ppt]
p1 :: Pressure at point 1 [R L2 T-2 ~> Pa]
drdt1 :: The partial derivative of density with temperature at point 1 [R C-1 ~> kg m-3 degC-1]
drds1 :: The partial derivative of density with salinity at point 1 [R S-1 ~> kg m-3 ppt-1]
t2 :: Temperature at point 2 [C ~> degC]
s2 :: Salinity at point 2 [S ~> ppt]
p2 :: Pressure at point 2 [R L2 T-2 ~> Pa]
drdt2 :: The partial derivative of density with temperature at point 2 [R C-1 ~> kg m-3 degC-1]
drds2 :: The partial derivative of density with salinity at point 2 [R S-1 ~> kg m-3 ppt-1]

Called from:

calc_delta_rho_and_derivs

function mom_neutral_diffusion/absolute_position(n, ns, Pint, Karr, NParr, k_surface) [real]¶

Converts non-dimensional position within a layer to absolute position (for debugging)

Parameters:

n :: [in] Number of levels
ns :: [in] Number of neutral surfaces
pint :: [in] Position of interfaces [R L2 T-2 ~> Pa] or other units
karr :: [in] Index of interface above position
nparr :: [in] Non-dimensional position within layer Karr(:) [nondim]
k_surface :: [in] k-interface to query

Return:

undefined :: The absolute position of a location [R L2 T-2 ~> Pa] or other units following Pint

Called from:

absolute_positions find_neutral_surface_positions_continuous

function mom_neutral_diffusion/absolute_positions(n, ns, Pint, Karr, NParr) [real]¶

Converts non-dimensional positions within layers to absolute positions (for debugging)

Parameters:

n :: [in] Number of levels
ns :: [in] Number of neutral surfaces
pint :: [in] Position of interface [R L2 T-2 ~> Pa] or other units
karr :: [in] Indexes of interfaces about positions
nparr :: [in] Non-dimensional positions within layers Karr(:) [nondim]

Return:

undefined :: Absolute positions [R L2 T-2 ~> Pa] or other units following Pint

Call to:

absolute_position

Called from:

ndiff_unit_tests_continuous

subroutine mom_neutral_diffusion/neutral_surface_flux(nk, nsurf, deg, hl, hr, Tl, Tr, PiL, PiR, KoL, KoR, hEff, Flx, continuous, h_neglect, remap_CS, h_neglect_edge, coeff_l, coeff_r)¶

Returns a single column of neutral diffusion fluxes of a tracer.

Parameters:

nk :: [in] Number of levels
nsurf :: [in] Number of neutral surfaces
deg :: [in] Degree of polynomial reconstructions
hl :: [in] Left-column layer thickness [H ~> m or kg m-2]
hr :: [in] Right-column layer thickness [H ~> m or kg m-2]
tl :: [in] Left-column layer tracer in arbitrary concentration units (e.g. [C ~> degC] for temperature)
tr :: [in] Right-column layer tracer in arbitrary concentration units (e.g. [C ~> degC] for temperature)
pil :: [in] Fractional position of neutral surface within layer KoL of left column [nondim]
pir :: [in] Fractional position of neutral surface within layer KoR of right column [nondim]
kol :: [in] Index of first left interface above neutral surface
kor :: [in] Index of first right interface above neutral surface
heff :: [in] Effective thickness between two neutral surfaces [H ~> m or kg m-2]
flx :: [inout] Flux of tracer between pairs of neutral layers in units (conc H or conc H L2) that depend on the presence and units of coeff_l and coeff_r. If the tracer is temperature, this could have units of [C H ~> degC m or degC kg m-2] or [C H L2 ~> degC m3 or degC kg] if coeff_l has units of [L2 ~> m2]
continuous :: [in] True if using continuous reconstruction
h_neglect :: [in] A negligibly small width for the purpose of cell reconstructions [H ~> m or kg m-2]
remap_cs :: [in] Remapping control structure used to create sublayers
h_neglect_edge :: [in] A negligibly small width used for edge value calculations if continuous is false [H ~> m or kg m-2]
coeff_l :: [in] Left-column diffusivity [L2 ~> m2] or [nondim]
coeff_r :: [in] Right-column diffusivity [L2 ~> m2] or [nondim]

Call to:

interface_scalar neutral_surface_t_eval ppm_ave ppm_left_right_edge_values signum

Called from:

ndiff_unit_tests_continuous neutral_diffusion

subroutine mom_neutral_diffusion/neutral_surface_t_eval(nk, ns, k_sub, Ks, Ps, T_mean, T_int, deg, iMethod, T_poly, T_top, T_bot, T_sub, T_top_int, T_bot_int, T_layer)¶

Evaluate various parts of the reconstructions to calculate gradient-based flux limiter.

Parameters:

nk :: [in] Number of cell averages
ns :: [in] Number of neutral surfaces
k_sub :: [in] Index of current neutral layer
ks :: [in] List of the layers associated with each neutral surface
ps :: [in] List of the positions within a layer of each surface [nondim]
t_mean :: [in] Layer average of tracer in arbitrary concentration units (e.g. [C ~> degC] for temperature)
t_int :: [in] Layer interface values of tracer from reconstruction in concentration units (e.g. [C ~> degC] for temperature)
deg :: [in] Degree of reconstruction polynomial (e.g. 1 is linear)
imethod :: [in] Method of integration to use
t_poly :: [in] Coefficients of polynomial reconstructions in arbitrary concentration units (e.g. [C ~> degC] for temperature)
t_top :: [out] Tracer value at top (across discontinuity if necessary) in concentration units (e.g. [C ~> degC] for temperature)
t_bot :: [out] Tracer value at bottom (across discontinuity if necessary) in concentration units (e.g. [C ~> degC] for temperature)
t_sub :: [out] Average of the tracer value over the sublayer in arbitrary concentration units (e.g. [C ~> degC] for temperature)
t_top_int :: [out] Tracer value at the top interface of a neutral layer in concentration units (e.g. [C ~> degC] for temperature)
t_bot_int :: [out] Tracer value at the bottom interface of a neutral layer in concentration units (e.g. [C ~> degC] for temperature)
t_layer :: [out] Cell-average tracer concentration in a layer that the reconstruction belongs to in concentration units (e.g. [C ~> degC] for temperature)

Call to:

mom_remapping::average_value_ppoly polynomial_functions::evaluation_polynomial mom_error_handler::mom_error

Called from:

neutral_surface_flux

subroutine mom_neutral_diffusion/ppm_left_right_edge_values(nk, Tl, Ti, aL, aR)¶

Discontinuous PPM reconstructions of the left/right edge values within a cell.

Parameters:

nk :: [in] Number of levels
tl :: [in] Layer tracer (conc, e.g. degC) in arbitrary units [A ~> a]
ti :: [in] Interface tracer (conc, e.g. degC) in arbitrary units [A ~> a]
al :: [inout] Left edge value of tracer (conc, e.g. degC) in the same arbitrary units as Tl and Ti [A ~> a]
ar :: [inout] Right edge value of tracer (conc, e.g. degC) in the same arbitrary units as Tl and Ti [A ~> a]

Call to:

signum

Called from:

neutral_surface_flux

function mom_neutral_diffusion/neutral_diffusion_unit_tests(verbose) [logical]¶

Returns true if unit tests of neutral_diffusion functions fail. Otherwise returns false.

Parameters:: verbose :: [in] If true, write results to stdout
Call to:: ndiff_unit_tests_continuous ndiff_unit_tests_discontinuous
Called from:: mom_unit_tests::unit_tests

function mom_neutral_diffusion/ndiff_unit_tests_continuous(verbose) [logical]¶

Returns true if unit tests of neutral_diffusion functions fail. Otherwise returns false.

Parameters:: verbose :: [in] If true, write results to stdout
Call to:: absolute_positions find_neutral_surface_positions_continuous interface_scalar neutral_surface_flux mom_io::stdout test_data1d test_fv_diff test_fvlsq_slope test_ifndp test_nsp
Called from:: neutral_diffusion_unit_tests

function mom_neutral_diffusion/ndiff_unit_tests_discontinuous(verbose) [logical]¶

Parameters:: verbose :: [in] It true, write results to stdout
Call to:: mom_eos::eos_linear find_neutral_pos_linear find_neutral_surface_positions_discontinuous mark_unstable_cells mom_io::stdout test_nsp test_rnp
Called from:: neutral_diffusion_unit_tests

function mom_neutral_diffusion/test_fv_diff(verbose, hkm1, hk, hkp1, Skm1, Sk, Skp1, Ptrue, title) [logical]¶

Returns true if a test of fv_diff() fails, and conditionally writes results to stream. fails, and conditionally writes results to stream.

Parameters:

verbose :: [in] If true, write results to stdout
hkm1 :: [in] Left cell width [nondim]
hk :: [in] Center cell width [nondim]
hkp1 :: [in] Right cell width [nondim]
skm1 :: [in] Left cell average value in arbitrary units [arbitrary]
sk :: [in] Center cell average value in arbitrary units [arbitrary]
skp1 :: [in] Right cell average value in arbitrary units [arbitrary]
ptrue :: [in] True answer in arbitrary units [arbitrary]
title :: [in] Title for messages

Call to:

fv_diff mom_io::stderr mom_io::stdout

Called from:

ndiff_unit_tests_continuous

function mom_neutral_diffusion/test_fvlsq_slope(verbose, hkm1, hk, hkp1, Skm1, Sk, Skp1, Ptrue, title) [logical]¶

Returns true if a test of fvlsq_slope() fails, and conditionally writes results to stream. fails, and conditionally writes results to stream.

Parameters:

verbose :: [in] If true, write results to stdout
hkm1 :: [in] Left cell width in arbitrary units [B ~> b]
hk :: [in] Center cell width in arbitrary units [B ~> b]
hkp1 :: [in] Right cell width in arbitrary units [B ~> b]
skm1 :: [in] Left cell average value in arbitrary units [A ~> a]
sk :: [in] Center cell average value in arbitrary units [A ~> a]
skp1 :: [in] Right cell average value in arbitrary units [A ~> a]
ptrue :: [in] True answer in arbitrary units [A B-1 ~> a b-1]
title :: [in] Title for messages

Call to:

fvlsq_slope mom_io::stderr mom_io::stdout

Called from:

ndiff_unit_tests_continuous

function mom_neutral_diffusion/test_ifndp(verbose, rhoNeg, Pneg, rhoPos, Ppos, Ptrue, title) [logical]¶

Returns true if a test of interpolate_for_nondim_position() fails, and conditionally writes results to stream. fails, and conditionally writes results to stream.

Parameters:

verbose :: [in] If true, write results to stdout
rhoneg :: [in] Lighter density [R ~> kg m-3]
pneg :: [in] Interface position of lighter density [nondim]
rhopos :: [in] Heavier density [R ~> kg m-3]
ppos :: [in] Interface position of heavier density [nondim]
ptrue :: [in] True answer [nondim]
title :: [in] Title for messages

Call to:

interpolate_for_nondim_position mom_io::stderr mom_io::stdout

Called from:

ndiff_unit_tests_continuous

function mom_neutral_diffusion/test_data1d(verbose, nk, Po, Ptrue, title) [logical]¶

Returns true if comparison of Po and Ptrue fails, and conditionally writes results to stream.

Parameters:

verbose :: [in] If true, write results to stdout
nk :: [in] Number of layers
po :: [in] Calculated answer [arbitrary]
ptrue :: [in] True answer [arbitrary]
title :: [in] Title for messages

Call to:

mom_io::stderr mom_io::stdout

Called from:

ndiff_unit_tests_continuous

function mom_neutral_diffusion/test_data1di(verbose, nk, Po, Ptrue, title) [logical]¶

Returns true if comparison of Po and Ptrue fails, and conditionally writes results to stream.

Parameters:

verbose :: [in] If true, write results to stdout
nk :: [in] Number of layers
po :: [in] Calculated answer [arbitrary]
ptrue :: [in] True answer [arbitrary]
title :: [in] Title for messages

Call to:

mom_io::stderr mom_io::stdout

function mom_neutral_diffusion/test_nsp(verbose, ns, KoL, KoR, pL, pR, hEff, KoL0, KoR0, pL0, pR0, hEff0, title) [logical]¶

Returns true if output of find_neutral_surface_positions() does not match correct values, and conditionally writes results to stream.

Parameters:

verbose :: [in] If true, write results to stdout
ns :: [in] Number of surfaces
kol :: [in] Index of first left interface above neutral surface
kor :: [in] Index of first right interface above neutral surface
pl :: [in] Fractional position of neutral surface within layer KoL of left column [nondim]
pr :: [in] Fractional position of neutral surface within layer KoR of right column [nondim]
heff :: [in] Effective thickness between two neutral surfaces [R L2 T-2 ~> Pa]
kol0 :: [in] Correct value for KoL
kor0 :: [in] Correct value for KoR
pl0 :: [in] Correct value for pL [nondim]
pr0 :: [in] Correct value for pR [nondim]
heff0 :: [in] Correct value for hEff [R L2 T-2 ~> Pa]
title :: [in] Title for messages

Call to:

compare_nsp_row mom_error_handler::mom_error mom_io::stderr mom_io::stdout

Called from:

ndiff_unit_tests_continuous ndiff_unit_tests_discontinuous

function mom_neutral_diffusion/compare_nsp_row(KoL, KoR, pL, pR, KoL0, KoR0, pL0, pR0) [logical]¶

Compares a single row, k, of output from find_neutral_surface_positions()

Parameters:

kol :: [in] Index of first left interface above neutral surface
kor :: [in] Index of first right interface above neutral surface
pl :: [in] Fractional position of neutral surface within layer KoL of left column [nondim]
pr :: [in] Fractional position of neutral surface within layer KoR of right column [nondim]
kol0 :: [in] Correct value for KoL
kor0 :: [in] Correct value for KoR
pl0 :: [in] Correct value for pL [nondim]
pr0 :: [in] Correct value for pR [nondim]

Called from:

test_nsp

function mom_neutral_diffusion/test_rnp(expected_pos, test_pos, title) [logical]¶

Compares output position from refine_nondim_position with an expected value.

Parameters:

expected_pos :: [in] The expected position [arbitrary]
test_pos :: [in] The position returned by the code [arbitrary]
title :: [in] A label for this test

Call to:

mom_io::stdout

Called from:

ndiff_unit_tests_discontinuous

subroutine mom_neutral_diffusion/neutral_diffusion_end(CS)¶

Deallocates neutral_diffusion control structure.

Parameters:: cs :: Neutral diffusion control structure