mom_density_integrals module reference¶
Provides integrals of density.
Functions/Subroutines¶
Calls the appropriate subroutine to calculate analytical and nearly-analytical integrals in z across layers of pressure anomalies, which are required for calculating the finite-volume form pressure accelerations in a Boussinesq model. |
|
Calculates (by numerical quadrature) integrals of pressure anomalies across layers, which are required for calculating the finite-volume form pressure accelerations in a Boussinesq model. |
|
Compute pressure gradient force integrals by quadrature for the case where T and S are linear profiles. |
|
Compute pressure gradient force integrals for layer “k” and the case where T and S are parabolic profiles. |
|
Calls the appropriate subroutine to calculate analytical and nearly-analytical integrals in pressure across layers of geopotential anomalies, which are required for calculating the finite-volume form pressure accelerations in a non-Boussinesq model. |
|
This subroutine calculates integrals of specific volume anomalies in pressure across layers, which are required for calculating the finite-volume form pressure accelerations in a non-Boussinesq model. |
|
This subroutine calculates integrals of specific volume anomalies in pressure across layers, which are required for calculating the finite-volume form pressure accelerations in a non-Boussinesq model. |
|
Find the depth at which the reconstructed pressure matches P_tgt. |
|
Calculate the average in situ specific volume across layers. |
|
Returns change in anomalous pressure change from top to non-dimensional position pos between z_t and z_b [R L2 T-2 ~> Pa]. |
Detailed Description¶
Provides integrals of density.
Function/Subroutine Documentation¶
-
subroutine
mom_density_integrals/
int_density_dz
(T, S, z_t, z_b, rho_ref, rho_0, G_e, HI, EOS, US, dpa, intz_dpa, intx_dpa, inty_dpa, bathyT, dz_neglect, useMassWghtInterp, Z_0p)¶ Calls the appropriate subroutine to calculate analytical and nearly-analytical integrals in z across layers of pressure anomalies, which are required for calculating the finite-volume form pressure accelerations in a Boussinesq model.
- Parameters:
hi :: [in] Ocean horizontal index structures for the arrays
t :: [in] Potential temperature referenced to the surface [C ~> degC]
s :: [in] Salinity [S ~> ppt]
z_t :: [in] Height at the top of the layer in depth units [Z ~> m]
z_b :: [in] Height at the bottom of the layer [Z ~> m]
rho_ref :: [in] A mean density [R ~> kg m-3], that is subtracted out to reduce the magnitude of each of the integrals.
rho_0 :: [in] A density [R ~> kg m-3], that is used to calculate the pressure (as p~=-z*rho_0*G_e) used in the equation of state.
g_e :: [in] The Earth’s gravitational acceleration [L2 Z-1 T-2 ~> m s-2]
eos :: [in] Equation of state structure
us :: [in] A dimensional unit scaling type
dpa :: [inout] The change in the pressure anomaly
intz_dpa :: [inout] The integral through the thickness of the
intx_dpa :: [inout] The integral in x of the difference between
inty_dpa :: [inout] The integral in y of the difference between
bathyt :: [in] The depth of the bathymetry [Z ~> m]
dz_neglect :: [in] A minuscule thickness change [Z ~> m]
usemasswghtinterp :: [in] If true, uses mass weighting to interpolate T/S for top and bottom integrals.
z_0p :: [in] The height at which the pressure is 0 [Z ~> m]
- Call to:
- Called from:
-
subroutine
mom_density_integrals/
int_density_dz_generic_pcm
(T, S, z_t, z_b, rho_ref, rho_0, G_e, HI, EOS, US, dpa, intz_dpa, intx_dpa, inty_dpa, bathyT, dz_neglect, useMassWghtInterp, use_inaccurate_form, Z_0p)¶ Calculates (by numerical quadrature) integrals of pressure anomalies across layers, which are required for calculating the finite-volume form pressure accelerations in a Boussinesq model.
- Parameters:
hi :: [in] Horizontal index type for input variables.
t :: [in] Potential temperature of the layer [C ~> degC]
s :: [in] Salinity of the layer [S ~> ppt]
z_t :: [in] Height at the top of the layer in depth units [Z ~> m]
z_b :: [in] Height at the bottom of the layer [Z ~> m]
rho_ref :: [in] A mean density [R ~> kg m-3], that is subtracted out to reduce the magnitude of each of the integrals.
rho_0 :: [in] A density [R ~> kg m-3], that is used to calculate the pressure (as p~=-z*rho_0*G_e) used in the equation of state.
g_e :: [in] The Earth’s gravitational acceleration [L2 Z-1 T-2 ~> m s-2]
eos :: [in] Equation of state structure
us :: [in] A dimensional unit scaling type
dpa :: [inout] The change in the pressure anomaly
intz_dpa :: [inout] The integral through the thickness of the
intx_dpa :: [inout] The integral in x of the difference between
inty_dpa :: [inout] The integral in y of the difference between
bathyt :: [in] The depth of the bathymetry [Z ~> m]
dz_neglect :: [in] A minuscule thickness change [Z ~> m]
usemasswghtinterp :: [in] If true, uses mass weighting to interpolate T/S for top and bottom integrals.
use_inaccurate_form :: [in] If true, uses an inaccurate form of density anomalies, as was used prior to March 2018.
z_0p :: [in] The height at which the pressure is 0 [Z ~> m]
- Call to:
- Called from:
-
subroutine
mom_density_integrals/
int_density_dz_generic_plm
(k, tv, T_t, T_b, S_t, S_b, e, rho_ref, rho_0, G_e, dz_subroundoff, bathyT, HI, GV, EOS, US, use_stanley_eos, dpa, intz_dpa, intx_dpa, inty_dpa, useMassWghtInterp, use_inaccurate_form, Z_0p)¶ Compute pressure gradient force integrals by quadrature for the case where T and S are linear profiles.
- Parameters:
k :: [in] Layer index to calculate integrals for
hi :: [in] Ocean horizontal index structures for the input arrays
gv :: [in] Vertical grid structure
tv :: [in] Thermodynamic variables
t_t :: [in] Potential temperature at the cell top [C ~> degC]
t_b :: [in] Potential temperature at the cell bottom [C ~> degC]
s_t :: [in] Salinity at the cell top [S ~> ppt]
s_b :: [in] Salinity at the cell bottom [S ~> ppt]
e :: [in] Height of interfaces [Z ~> m]
rho_ref :: [in] A mean density [R ~> kg m-3], that is subtracted out to reduce the magnitude of each of the integrals.
rho_0 :: [in] A density [R ~> kg m-3], that is used to calculate the pressure (as p~=-z*rho_0*G_e) used in the equation of state.
g_e :: [in] The Earth’s gravitational acceleration [L2 Z-1 T-2 ~> m s-2]
dz_subroundoff :: [in] A minuscule thickness change [Z ~> m]
bathyt :: [in] The depth of the bathymetry [Z ~> m]
eos :: [in] Equation of state structure
us :: [in] A dimensional unit scaling type
use_stanley_eos :: [in] If true, turn on Stanley SGS T variance parameterization
dpa :: [inout] The change in the pressure anomaly across the layer [R L2 T-2 ~> Pa]
intz_dpa :: [inout] The integral through the thickness of the layer of
intx_dpa :: [inout] The integral in x of the difference between the
inty_dpa :: [inout] The integral in y of the difference between the
usemasswghtinterp :: [in] If true, uses mass weighting to interpolate T/S for top and bottom integrals.
use_inaccurate_form :: [in] If true, uses an inaccurate form of density anomalies, as was used prior to March 2018.
z_0p :: [in] The height at which the pressure is 0 [Z ~> m]
-
subroutine
mom_density_integrals/
int_density_dz_generic_ppm
(k, tv, T_t, T_b, S_t, S_b, e, rho_ref, rho_0, G_e, dz_subroundoff, bathyT, HI, GV, EOS, US, use_stanley_eos, dpa, intz_dpa, intx_dpa, inty_dpa, useMassWghtInterp, Z_0p)¶ Compute pressure gradient force integrals for layer “k” and the case where T and S are parabolic profiles.
- Parameters:
k :: [in] Layer index to calculate integrals for
hi :: [in] Ocean horizontal index structures for the input arrays
gv :: [in] Vertical grid structure
tv :: [in] Thermodynamic variables
t_t :: [in] Potential temperature at the cell top [C ~> degC]
t_b :: [in] Potential temperature at the cell bottom [C ~> degC]
s_t :: [in] Salinity at the cell top [S ~> ppt]
s_b :: [in] Salinity at the cell bottom [S ~> ppt]
e :: [in] Height of interfaces [Z ~> m]
rho_ref :: [in] A mean density [R ~> kg m-3], that is subtracted out to reduce the magnitude of each of the integrals.
rho_0 :: [in] A density [R ~> kg m-3], that is used to calculate the pressure (as p~=-z*rho_0*G_e) used in the equation of state.
g_e :: [in] The Earth’s gravitational acceleration [L2 Z-1 T-2 ~> m s-2]
dz_subroundoff :: [in] A minuscule thickness change [Z ~> m]
bathyt :: [in] The depth of the bathymetry [Z ~> m]
eos :: [in] Equation of state structure
us :: [in] A dimensional unit scaling type
use_stanley_eos :: [in] If true, turn on Stanley SGS T variance parameterization
dpa :: [inout] The change in the pressure anomaly across the layer [R L2 T-2 ~> Pa]
intz_dpa :: [inout] The integral through the thickness of the layer of
intx_dpa :: [inout] The integral in x of the difference between the
inty_dpa :: [inout] The integral in y of the difference between the
usemasswghtinterp :: [in] If true, uses mass weighting to interpolate T/S for top and bottom integrals.
z_0p :: [in] The height at which the pressure is 0 [Z ~> m]
-
subroutine
mom_density_integrals/
int_specific_vol_dp
(T, S, p_t, p_b, alpha_ref, HI, EOS, US, dza, intp_dza, intx_dza, inty_dza, halo_size, bathyP, dP_tiny, useMassWghtInterp)¶ Calls the appropriate subroutine to calculate analytical and nearly-analytical integrals in pressure across layers of geopotential anomalies, which are required for calculating the finite-volume form pressure accelerations in a non-Boussinesq model. There are essentially no free assumptions, apart from the use of Boole’s rule to do the horizontal integrals, and from a truncation in the series for log(1-eps/1+eps) that assumes that |eps| < 0.34.
- Parameters:
hi :: [in] The horizontal index structure
t :: [in] Potential temperature referenced to the surface [C ~> degC]
s :: [in] Salinity [S ~> ppt]
p_t :: [in] Pressure at the top of the layer [R L2 T-2 ~> Pa]
p_b :: [in] Pressure at the bottom of the layer [R L2 T-2 ~> Pa]
alpha_ref :: [in] A mean specific volume that is subtracted out to reduce the magnitude of each of the integrals [R-1 ~> m3 kg-1] The calculation is mathematically identical with different values of alpha_ref, but this reduces the effects of roundoff.
eos :: [in] Equation of state structure
us :: [in] A dimensional unit scaling type
dza :: [inout] The change in the geopotential anomaly across
intp_dza :: [inout] The integral in pressure through the layer of the
intx_dza :: [inout] The integral in x of the difference between the
inty_dza :: [inout] The integral in y of the difference between the
halo_size :: [in] The width of halo points on which to calculate dza.
bathyp :: [in] The pressure at the bathymetry [R L2 T-2 ~> Pa]
dp_tiny :: [in] A minuscule pressure change with the same units as p_t [R L2 T-2 ~> Pa]
usemasswghtinterp :: [in] If true, uses mass weighting to interpolate T/S for top and bottom integrals.
- Call to:
- Called from:
mom_interface_heights::dz_to_thickness::dz_to_thickness_eos
mom_interface_heights::find_eta::find_eta_2d
mom_interface_heights::find_eta::find_eta_3d
mom_pressureforce_mont::pressureforce_mont_nonbouss
-
subroutine
mom_density_integrals/
int_spec_vol_dp_generic_pcm
(T, S, p_t, p_b, alpha_ref, HI, EOS, US, dza, intp_dza, intx_dza, inty_dza, halo_size, bathyP, dP_neglect, useMassWghtInterp)¶ This subroutine calculates integrals of specific volume anomalies in pressure across layers, which are required for calculating the finite-volume form pressure accelerations in a non-Boussinesq model. There are essentially no free assumptions, apart from the use of Boole’s rule quadrature to do the integrals.
- Parameters:
hi :: [in] A horizontal index type structure.
t :: [in] Potential temperature of the layer [C ~> degC]
s :: [in] Salinity of the layer [S ~> ppt]
p_t :: [in] Pressure atop the layer [R L2 T-2 ~> Pa]
p_b :: [in] Pressure below the layer [R L2 T-2 ~> Pa]
alpha_ref :: [in] A mean specific volume that is subtracted out to reduce the magnitude of each of the integrals [R-1 ~> m3 kg-1] The calculation is mathematically identical with different values of alpha_ref, but alpha_ref alters the effects of roundoff, and answers do change.
eos :: [in] Equation of state structure
us :: [in] A dimensional unit scaling type
dza :: [inout] The change in the geopotential anomaly
intp_dza :: [inout] The integral in pressure through the layer of
intx_dza :: [inout] The integral in x of the difference between
inty_dza :: [inout] The integral in y of the difference between
halo_size :: [in] The width of halo points on which to calculate dza.
bathyp :: [in] The pressure at the bathymetry [R L2 T-2 ~> Pa]
dp_neglect :: [in] A minuscule pressure change with the same units as p_t [R L2 T-2 ~> Pa]
usemasswghtinterp :: [in] If true, uses mass weighting to interpolate T/S for top and bottom integrals.
- Call to:
- Called from:
-
subroutine
mom_density_integrals/
int_spec_vol_dp_generic_plm
(T_t, T_b, S_t, S_b, p_t, p_b, alpha_ref, dP_neglect, bathyP, HI, EOS, US, dza, intp_dza, intx_dza, inty_dza, useMassWghtInterp)¶ This subroutine calculates integrals of specific volume anomalies in pressure across layers, which are required for calculating the finite-volume form pressure accelerations in a non-Boussinesq model. There are essentially no free assumptions, apart from the use of Boole’s rule quadrature to do the integrals.
- Parameters:
hi :: [in] A horizontal index type structure.
t_t :: [in] Potential temperature at the top of the layer [C ~> degC]
t_b :: [in] Potential temperature at the bottom of the layer [C ~> degC]
s_t :: [in] Salinity at the top the layer [S ~> ppt]
s_b :: [in] Salinity at the bottom the layer [S ~> ppt]
p_t :: [in] Pressure atop the layer [R L2 T-2 ~> Pa]
p_b :: [in] Pressure below the layer [R L2 T-2 ~> Pa]
alpha_ref :: [in] A mean specific volume that is subtracted out to reduce the magnitude of each of the integrals [R-1 ~> m3 kg-1] The calculation is mathematically identical with different values of alpha_ref, but alpha_ref alters the effects of roundoff, and answers do change.
dp_neglect :: [in] !< A miniscule pressure change with the same units as p_t [R L2 T-2 ~> Pa]
bathyp :: [in] The pressure at the bathymetry [R L2 T-2 ~> Pa]
eos :: [in] Equation of state structure
us :: [in] A dimensional unit scaling type
dza :: [inout] The change in the geopotential anomaly
intp_dza :: [inout] The integral in pressure through the layer of
intx_dza :: [inout] The integral in x of the difference between
inty_dza :: [inout] The integral in y of the difference between
usemasswghtinterp :: [in] If true, uses mass weighting to interpolate T/S for top and bottom integrals.
-
subroutine
mom_density_integrals/
find_depth_of_pressure_in_cell
(T_t, T_b, S_t, S_b, z_t, z_b, P_t, P_tgt, rho_ref, G_e, EOS, US, P_b, z_out, z_tol)¶ Find the depth at which the reconstructed pressure matches P_tgt.
- Parameters:
t_t :: [in] Potential temperature at the cell top [C ~> degC]
t_b :: [in] Potential temperature at the cell bottom [C ~> degC]
s_t :: [in] Salinity at the cell top [S ~> ppt]
s_b :: [in] Salinity at the cell bottom [S ~> ppt]
z_t :: [in] Absolute height of top of cell [Z ~> m] (Boussinesq ??]
z_b :: [in] Absolute height of bottom of cell [Z ~> m]
p_t :: [in] Anomalous pressure of top of cell, relative to g*rho_ref*z_t [R L2 T-2 ~> Pa]
p_tgt :: [in] Target pressure at height z_out, relative to g*rho_ref*z_out [R L2 T-2 ~> Pa]
rho_ref :: [in] Reference density with which calculation are anomalous to [R ~> kg m-3]
g_e :: [in] Gravitational acceleration [L2 Z-1 T-2 ~> m s-2]
eos :: [in] Equation of state structure
us :: [in] A dimensional unit scaling type
p_b :: [out] Pressure at the bottom of the cell [R L2 T-2 ~> Pa]
z_out :: [out] Absolute depth at which anomalous pressure = p_tgt [Z ~> m]
z_tol :: [in] The tolerance in finding z_out [Z ~> m]
- Call to:
- Called from:
-
subroutine
mom_density_integrals/
avg_specific_vol
(T, S, p_t, dp, HI, EOS, SpV_avg, halo_size)¶ Calculate the average in situ specific volume across layers.
- Parameters:
hi :: [in] The horizontal index structure
t :: [in] Potential temperature of the layer [C ~> degC]
s :: [in] Salinity of the layer [S ~> ppt]
p_t :: [in] Pressure at the top of the layer [R L2 T-2 ~> Pa]
dp :: [in] Pressure change in the layer [R L2 T-2 ~> Pa]
eos :: [in] Equation of state structure
spv_avg :: [inout] The vertical average specific volume
halo_size :: [in] The number of halo points in which to work.
- Call to:
- Called from:
-
function
mom_density_integrals/
frac_dp_at_pos
(T_t, T_b, S_t, S_b, z_t, z_b, rho_ref, G_e, pos, EOS) [real]¶ Returns change in anomalous pressure change from top to non-dimensional position pos between z_t and z_b [R L2 T-2 ~> Pa].
- Parameters:
t_t :: [in] Potential temperature at the cell top [C ~> degC]
t_b :: [in] Potential temperature at the cell bottom [C ~> degC]
s_t :: [in] Salinity at the cell top [S ~> ppt]
s_b :: [in] Salinity at the cell bottom [S ~> ppt]
z_t :: [in] The geometric height at the top of the layer [Z ~> m]
z_b :: [in] The geometric height at the bottom of the layer [Z ~> m]
rho_ref :: [in] A mean density [R ~> kg m-3], that is subtracted out to reduce the magnitude of each of the integrals.
g_e :: [in] The Earth’s gravitational acceleration [L2 Z-1 T-2 ~> m s-2]
pos :: [in] The fractional vertical position, 0 to 1 [nondim]
eos :: [in] Equation of state structure
- Called from: