Skip to content

Commit

Permalink
Merge pull request #1351 from adcroft/units-for-slope
Browse files Browse the repository at this point in the history 
Change units of slope returned from calc_isoneutral_slopes() to "Z L-1"
  • Loading branch information
@marshallward
marshallward authored Mar 12, 2021
2 parents 4255ada + ab241b6 commit a6f680c
Show file tree
Hide file tree
Showing 3 changed files with 38 additions and 36 deletions.
26 changes: 13 additions & 13 deletions src/core/MOM_isopycnal_slopes.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -37,14 +37,14 @@ subroutine calc_isoneutral_slopes(G, GV, US, h, e, tv, dt_kappa_smooth, &
!! thermodynamic variables
real, intent(in) :: dt_kappa_smooth !< A smoothing vertical diffusivity
!! times a smoothing timescale [Z2 ~> m2].
real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)+1), intent(inout) :: slope_x !< Isopycnal slope in i-direction [nondim]
real, dimension(SZI_(G),SZJB_(G),SZK_(GV)+1), intent(inout) :: slope_y !< Isopycnal slope in j-direction [nondim]
real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)+1), intent(inout) :: slope_x !< Isopycnal slope in i-dir [Z L-1 ~> nondim]
real, dimension(SZI_(G),SZJB_(G),SZK_(GV)+1), intent(inout) :: slope_y !< Isopycnal slope in j-dir [Z L-1 ~> nondim]
real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)+1), &
optional, intent(inout) :: N2_u !< Brunt-Vaisala frequency squared at
!! interfaces between u-points [T-2 ~> s-2]
!! interfaces between u-points [L2 Z-2 T-2 ~> s-2]
real, dimension(SZI_(G),SZJB_(G),SZK_(GV)+1), &
optional, intent(inout) :: N2_v !< Brunt-Vaisala frequency squared at
!! interfaces between u-points [T-2 ~> s-2]
!! interfaces between v-points [L2 Z-2 T-2 ~> s-2]
integer, optional, intent(in) :: halo !< Halo width over which to compute
type(ocean_OBC_type), optional, pointer :: OBC !< Open boundaries control structure.

Expand Down Expand Up @@ -86,15 +86,15 @@ subroutine calc_isoneutral_slopes(G, GV, US, h, e, tv, dt_kappa_smooth, &
real :: drdz ! Vertical density gradient [R Z-1 ~> kg m-4].
real :: Slope ! The slope of density surfaces, calculated in a way
! that is always between -1 and 1.
real :: mag_grad2 ! The squared magnitude of the 3-d density gradient [R2 L-2 ~> kg2 m-8].
real :: mag_grad2 ! The squared magnitude of the 3-d density gradient [R2 Z-2 ~> kg2 m-8].
real :: slope2_Ratio ! The ratio of the slope squared to slope_max squared.
real :: h_neglect ! A thickness that is so small it is usually lost
! in roundoff and can be neglected [H ~> m or kg m-2].
real :: h_neglect2 ! h_neglect^2 [H2 ~> m2 or kg2 m-4].
real :: dz_neglect ! A change in interface heighs that is so small it is usually lost
! in roundoff and can be neglected [Z ~> m].
logical :: use_EOS ! If true, density is calculated from T & S using an equation of state.
real :: G_Rho0 ! The gravitational acceleration divided by density [Z2 T-2 R-1 ~> m5 kg-2 s-2]
real :: G_Rho0 ! The gravitational acceleration divided by density [L2 Z-1 T-2 R-1 ~> m4 s-2 kg-1]
real :: Z_to_L ! A conversion factor between from units for e to the
! units for lateral distances.
real :: L_to_Z ! A conversion factor between from units for lateral distances
Expand Down Expand Up @@ -134,7 +134,7 @@ subroutine calc_isoneutral_slopes(G, GV, US, h, e, tv, dt_kappa_smooth, &

present_N2_u = PRESENT(N2_u)
present_N2_v = PRESENT(N2_v)
G_Rho0 = (US%L_to_Z*L_to_Z*GV%g_Earth) / GV%Rho0
G_Rho0 = GV%g_Earth / GV%Rho0
if (present_N2_u) then
do j=js,je ; do I=is-1,ie
N2_u(I,j,1) = 0.
Expand Down Expand Up @@ -248,17 +248,17 @@ subroutine calc_isoneutral_slopes(G, GV, US, h, e, tv, dt_kappa_smooth, &

! This estimate of slope is accurate for small slopes, but bounded
! to be between -1 and 1.
mag_grad2 = drdx**2 + (L_to_Z*drdz)**2
mag_grad2 = (Z_to_L*drdx)**2 + drdz**2
if (mag_grad2 > 0.0) then
slope_x(I,j,K) = drdx / sqrt(mag_grad2)
else ! Just in case mag_grad2 = 0 ever.
slope_x(I,j,K) = 0.0
endif

if (present_N2_u) N2_u(I,j,k) = G_Rho0 * drdz * G%mask2dCu(I,j) ! Square of buoyancy frequency [T-2 ~> s-2]
if (present_N2_u) N2_u(I,j,k) = G_Rho0 * drdz * G%mask2dCu(I,j) ! Square of buoyancy freq. [L2 Z-2 T-2 ~> s-2]

else ! With .not.use_EOS, the layers are constant density.
slope_x(I,j,K) = (Z_to_L*(e(i,j,K)-e(i+1,j,K))) * G%IdxCu(I,j)
slope_x(I,j,K) = (e(i,j,K)-e(i+1,j,K)) * G%IdxCu(I,j)
endif
if (local_open_u_BC) then
l_seg = OBC%segnum_u(I,j)
Expand Down Expand Up @@ -351,17 +351,17 @@ subroutine calc_isoneutral_slopes(G, GV, US, h, e, tv, dt_kappa_smooth, &

! This estimate of slope is accurate for small slopes, but bounded
! to be between -1 and 1.
mag_grad2 = drdy**2 + (L_to_Z*drdz)**2
mag_grad2 = (Z_to_L*drdy)**2 + drdz**2
if (mag_grad2 > 0.0) then
slope_y(i,J,K) = drdy / sqrt(mag_grad2)
else ! Just in case mag_grad2 = 0 ever.
slope_y(i,J,K) = 0.0
endif

if (present_N2_v) N2_v(i,J,k) = G_Rho0 * drdz * G%mask2dCv(i,J) ! Square of buoyancy frequency [T-2 ~> s-2]
if (present_N2_v) N2_v(i,J,k) = G_Rho0 * drdz * G%mask2dCv(i,J) ! Square of buoyancy freq. [L2 Z-2 T-2 ~> s-2]

else ! With .not.use_EOS, the layers are constant density.
slope_y(i,J,K) = (Z_to_L*(e(i,j,K)-e(i,j+1,K))) * G%IdyCv(i,J)
slope_y(i,J,K) = (e(i,j,K)-e(i,j+1,K)) * G%IdyCv(i,J)
endif
if (local_open_v_BC) then
l_seg = OBC%segnum_v(i,J)
Expand Down
10 changes: 6 additions & 4 deletions src/parameterizations/lateral/MOM_lateral_mixing_coeffs.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1125,16 +1125,18 @@ subroutine VarMix_init(Time, G, GV, US, param_file, diag, CS)
if (CS%calculate_Eady_growth_rate .and. CS%use_stored_slopes) then
CS%id_N2_u = register_diag_field('ocean_model', 'N2_u', diag%axesCui, Time, &
'Square of Brunt-Vaisala frequency, N^2, at u-points, as used in Visbeck et al.', &
's-2', conversion=US%s_to_T**2)
's-2', conversion=(US%L_to_Z*US%s_to_T)**2)
CS%id_N2_v = register_diag_field('ocean_model', 'N2_v', diag%axesCvi, Time, &
'Square of Brunt-Vaisala frequency, N^2, at v-points, as used in Visbeck et al.', &
's-2', conversion=US%s_to_T**2)
's-2', conversion=(US%L_to_Z*US%s_to_T)**2)
endif
if (CS%use_stored_slopes) then
CS%id_S2_u = register_diag_field('ocean_model', 'S2_u', diag%axesCu1, Time, &
'Depth average square of slope magnitude, S^2, at u-points, as used in Visbeck et al.', 'nondim')
'Depth average square of slope magnitude, S^2, at u-points, as used in Visbeck et al.', &
'nondim', conversion=US%Z_to_L**2)
CS%id_S2_v = register_diag_field('ocean_model', 'S2_v', diag%axesCv1, Time, &
'Depth average square of slope magnitude, S^2, at v-points, as used in Visbeck et al.', 'nondim')
'Depth average square of slope magnitude, S^2, at v-points, as used in Visbeck et al.', &
'nondim', conversion=US%Z_to_L**2)
endif

oneOrTwo = 1.0
Expand Down
38 changes: 19 additions & 19 deletions src/parameterizations/lateral/MOM_thickness_diffuse.F90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -40,7 +40,7 @@ module MOM_thickness_diffuse
real :: max_Khth_CFL !< Maximum value of the diffusive CFL for thickness diffusion
real :: Khth_Min !< Minimum value of Khth [L2 T-1 ~> m2 s-1]
real :: Khth_Max !< Maximum value of Khth [L2 T-1 ~> m2 s-1], or 0 for no max
real :: slope_max !< Slopes steeper than slope_max are limited in some way [nondim].
real :: slope_max !< Slopes steeper than slope_max are limited in some way [Z L-1 ~> nondim].
real :: kappa_smooth !< Vertical diffusivity used to interpolate more
!! sensible values of T & S into thin layers [Z2 T-1 ~> m2 s-1].
logical :: thickness_diffuse !< If true, interfaces heights are diffused.
Expand Down Expand Up @@ -83,8 +83,8 @@ module MOM_thickness_diffuse

type(diag_ctrl), pointer :: diag => NULL() !< structure used to regulate timing of diagnostics
real, pointer :: GMwork(:,:) => NULL() !< Work by thickness diffusivity [R Z L2 T-3 ~> W m-2]
real, pointer :: diagSlopeX(:,:,:) => NULL() !< Diagnostic: zonal neutral slope [nondim]
real, pointer :: diagSlopeY(:,:,:) => NULL() !< Diagnostic: zonal neutral slope [nondim]
real, pointer :: diagSlopeX(:,:,:) => NULL() !< Diagnostic: zonal neutral slope [Z L-1 ~> nondim]
real, pointer :: diagSlopeY(:,:,:) => NULL() !< Diagnostic: zonal neutral slope [Z L-1 ~> nondim]

real, dimension(:,:,:), pointer :: &
KH_u_GME => NULL(), & !< interface height diffusivities in u-columns [L2 T-1 ~> m2 s-1]
Expand Down Expand Up @@ -578,8 +578,8 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV
!! the isopycnal slopes are taken directly from
!! the interface slopes without consideration of
!! density gradients [nondim].
real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)+1), optional, intent(in) :: slope_x !< Isopycnal slope at u-points
real, dimension(SZI_(G),SZJB_(G),SZK_(GV)+1), optional, intent(in) :: slope_y !< Isopycnal slope at v-points
real, dimension(SZIB_(G),SZJ_(G),SZK_(GV)+1), optional, intent(in) :: slope_x !< Isopyc. slope at u [Z L-1 ~> nondim]
real, dimension(SZI_(G),SZJB_(G),SZK_(GV)+1), optional, intent(in) :: slope_y !< Isopyc. slope at v [Z L-1 ~> nondim]
! Local variables
real, dimension(SZI_(G), SZJ_(G), SZK_(GV)) :: &
T, & ! The temperature (or density) [degC], with the values in
Expand Down Expand Up @@ -660,7 +660,7 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV
real :: Slope ! The slope of density surfaces, calculated in a way
! that is always between -1 and 1, nondimensional.
real :: mag_grad2 ! The squared magnitude of the 3-d density gradient [R2 L-2 ~> kg2 m-8].
real :: I_slope_max2 ! The inverse of slope_max squared, nondimensional.
real :: I_slope_max2 ! The inverse of slope_max squared [L2 Z-2 ~> nondim].
real :: h_neglect ! A thickness that is so small it is usually lost
! in roundoff and can be neglected [H ~> m or kg m-2].
real :: h_neglect2 ! h_neglect^2 [H2 ~> m2 or kg2 m-4].
Expand Down Expand Up @@ -919,7 +919,7 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV

! This estimate of slope is accurate for small slopes, but bounded
! to be between -1 and 1.
mag_grad2 = drdx**2 + (US%L_to_Z*drdz)**2
mag_grad2 = (US%Z_to_L*drdx)**2 + drdz**2
if (mag_grad2 > 0.0) then
Slope = drdx / sqrt(mag_grad2)
slope2_Ratio_u(I,K) = Slope**2 * I_slope_max2
Expand All @@ -933,7 +933,7 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV
! that ignore density gradients along layers.
if (present_int_slope_u) then
Slope = (1.0 - int_slope_u(I,j,K)) * Slope + &
int_slope_u(I,j,K) * US%Z_to_L*((e(i+1,j,K)-e(i,j,K)) * G%IdxCu(I,j))
int_slope_u(I,j,K) * ((e(i+1,j,K)-e(i,j,K)) * G%IdxCu(I,j))
slope2_Ratio_u(I,K) = (1.0 - int_slope_u(I,j,K)) * slope2_Ratio_u(I,K)
endif

Expand All @@ -942,7 +942,7 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV
if (CS%id_slope_x > 0) CS%diagSlopeX(I,j,k) = Slope

! Estimate the streamfunction at each interface [Z L2 T-1 ~> m3 s-1].
Sfn_unlim_u(I,K) = -((KH_u(I,j,K)*G%dy_Cu(I,j))*US%L_to_Z*Slope)
Sfn_unlim_u(I,K) = -((KH_u(I,j,K)*G%dy_Cu(I,j))*Slope)

! Avoid moving dense water upslope from below the level of
! the bottom on the receiving side.
Expand All @@ -968,10 +968,10 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV
if (present_slope_x) then
Slope = slope_x(I,j,k)
else
Slope = US%Z_to_L*((e(i,j,K)-e(i+1,j,K))*G%IdxCu(I,j)) * G%mask2dCu(I,j)
Slope = ((e(i,j,K)-e(i+1,j,K))*G%IdxCu(I,j)) * G%mask2dCu(I,j)
endif
if (CS%id_slope_x > 0) CS%diagSlopeX(I,j,k) = Slope
Sfn_unlim_u(I,K) = ((KH_u(I,j,K)*G%dy_Cu(I,j))*US%L_to_Z*Slope)
Sfn_unlim_u(I,K) = ((KH_u(I,j,K)*G%dy_Cu(I,j))*Slope)
hN2_u(I,K) = GV%g_prime(K)
endif ! if (use_EOS)
else ! if (k > nk_linear)
Expand Down Expand Up @@ -1185,7 +1185,7 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV

! This estimate of slope is accurate for small slopes, but bounded
! to be between -1 and 1.
mag_grad2 = drdy**2 + (US%L_to_Z*drdz)**2
mag_grad2 = (US%Z_to_L*drdy)**2 + drdz**2
if (mag_grad2 > 0.0) then
Slope = drdy / sqrt(mag_grad2)
slope2_Ratio_v(i,K) = Slope**2 * I_slope_max2
Expand All @@ -1199,7 +1199,7 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV
! that ignore density gradients along layers.
if (present_int_slope_v) then
Slope = (1.0 - int_slope_v(i,J,K)) * Slope + &
int_slope_v(i,J,K) * US%Z_to_L*((e(i,j+1,K)-e(i,j,K)) * G%IdyCv(i,J))
int_slope_v(i,J,K) * ((e(i,j+1,K)-e(i,j,K)) * G%IdyCv(i,J))
slope2_Ratio_v(i,K) = (1.0 - int_slope_v(i,J,K)) * slope2_Ratio_v(i,K)
endif

Expand All @@ -1208,7 +1208,7 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV
if (CS%id_slope_y > 0) CS%diagSlopeY(I,j,k) = Slope

! Estimate the streamfunction at each interface [Z L2 T-1 ~> m3 s-1].
Sfn_unlim_v(i,K) = -((KH_v(i,J,K)*G%dx_Cv(i,J))*US%L_to_Z*Slope)
Sfn_unlim_v(i,K) = -((KH_v(i,J,K)*G%dx_Cv(i,J))*Slope)

! Avoid moving dense water upslope from below the level of
! the bottom on the receiving side.
Expand All @@ -1234,10 +1234,10 @@ subroutine thickness_diffuse_full(h, e, Kh_u, Kh_v, tv, uhD, vhD, cg1, dt, G, GV
if (present_slope_y) then
Slope = slope_y(i,J,k)
else
Slope = US%Z_to_L*((e(i,j,K)-e(i,j+1,K))*G%IdyCv(i,J)) * G%mask2dCv(i,J)
Slope = ((e(i,j,K)-e(i,j+1,K))*G%IdyCv(i,J)) * G%mask2dCv(i,J)
endif
if (CS%id_slope_y > 0) CS%diagSlopeY(I,j,k) = Slope
Sfn_unlim_v(i,K) = ((KH_v(i,J,K)*G%dx_Cv(i,J))*US%L_to_Z*Slope)
Sfn_unlim_v(i,K) = ((KH_v(i,J,K)*G%dx_Cv(i,J))*Slope)
hN2_v(i,K) = GV%g_prime(K)
endif ! if (use_EOS)
else ! if (k > nk_linear)
Expand Down Expand Up @@ -1947,7 +1947,7 @@ subroutine thickness_diffuse_init(Time, G, GV, US, param_file, diag, CDp, CS)
"longer than DT, or 0 to use DT.", units="s", default=0.0, scale=US%s_to_T)
call get_param(param_file, mdl, "KHTH_SLOPE_MAX", CS%slope_max, &
"A slope beyond which the calculated isopycnal slope is "//&
"not reliable and is scaled away.", units="nondim", default=0.01)
"not reliable and is scaled away.", units="nondim", default=0.01, scale=US%L_to_Z)
call get_param(param_file, mdl, "KD_SMOOTH", CS%kappa_smooth, &
"A diapycnal diffusivity that is used to interpolate "//&
"more sensible values of T & S into thin layers.", &
Expand Down Expand Up @@ -2065,10 +2065,10 @@ subroutine thickness_diffuse_init(Time, G, GV, US, param_file, diag, CDp, CS)
'm2 s-1', conversion=US%L_to_m**2*US%s_to_T)

CS%id_slope_x = register_diag_field('ocean_model', 'neutral_slope_x', diag%axesCui, Time, &
'Zonal slope of neutral surface', 'nondim')
'Zonal slope of neutral surface', 'nondim', conversion=US%Z_to_L)
if (CS%id_slope_x > 0) call safe_alloc_ptr(CS%diagSlopeX,G%IsdB,G%IedB,G%jsd,G%jed,GV%ke+1)
CS%id_slope_y = register_diag_field('ocean_model', 'neutral_slope_y', diag%axesCvi, Time, &
'Meridional slope of neutral surface', 'nondim')
'Meridional slope of neutral surface', 'nondim', conversion=US%Z_to_L)
if (CS%id_slope_y > 0) call safe_alloc_ptr(CS%diagSlopeY,G%isd,G%ied,G%JsdB,G%JedB,GV%ke+1)
CS%id_sfn_x = register_diag_field('ocean_model', 'GM_sfn_x', diag%axesCui, Time, &
'Parameterized Zonal Overturning Streamfunction', &
Expand Down

0 comments on commit a6f680c

Please sign in to comment.