Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Options to enforce KHTR_MIN and KHTH_MIN in the whole water column and fix naming inconsistency #294

Merged
merged 3 commits into from
Aug 9, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
52 changes: 43 additions & 9 deletions src/parameterizations/lateral/MOM_thickness_diffuse.F90
Original file line number Diff line number Diff line change
Expand Up @@ -49,6 +49,9 @@ module MOM_thickness_diffuse
real :: kappa_smooth !< Vertical diffusivity used to interpolate more sensible values
!! of T & S into thin layers [H Z T-1 ~> m2 s-1 or kg m-1 s-1]
logical :: thickness_diffuse !< If true, interfaces heights are diffused.
logical :: full_depth_khth_min !< If true, KHTH_MIN is enforced throughout the whole water column.
!! Otherwise, KHTH_MIN is only enforced at the surface. This parameter
!! is only available when KHTH_USE_EBT_STRUCT=True and KHTH_MIN>0.
logical :: use_FGNV_streamfn !< If true, use the streamfunction formulation of
!! Ferrari et al., 2010, which effectively emphasizes
!! graver vertical modes by smoothing in the vertical.
Expand Down Expand Up @@ -294,10 +297,18 @@ subroutine thickness_diffuse(h, uhtr, vhtr, tv, dt, G, GV, US, MEKE, VarMix, CDp
enddo ; enddo

if (khth_use_ebt_struct) then
!$OMP do
do K=2,nz+1 ; do j=js,je ; do I=is-1,ie
KH_u(I,j,K) = KH_u(I,j,1) * 0.5 * ( VarMix%ebt_struct(i,j,k-1) + VarMix%ebt_struct(i+1,j,k-1) )
enddo ; enddo ; enddo
if (CS%full_depth_khth_min) then
!$OMP do
do K=2,nz+1 ; do j=js,je ; do I=is-1,ie
KH_u(I,j,K) = KH_u(I,j,1) * 0.5 * ( VarMix%ebt_struct(i,j,k-1) + VarMix%ebt_struct(i+1,j,k-1) )
KH_u(I,j,K) = max(KH_u(I,j,K), CS%Khth_Min)
enddo ; enddo ; enddo
else
!$OMP do
do K=2,nz+1 ; do j=js,je ; do I=is-1,ie
KH_u(I,j,K) = KH_u(I,j,1) * 0.5 * ( VarMix%ebt_struct(i,j,k-1) + VarMix%ebt_struct(i+1,j,k-1) )
enddo ; enddo ; enddo
endif
else
!$OMP do
do K=2,nz+1 ; do j=js,je ; do I=is-1,ie
Expand Down Expand Up @@ -390,10 +401,18 @@ subroutine thickness_diffuse(h, uhtr, vhtr, tv, dt, G, GV, US, MEKE, VarMix, CDp
endif

if (khth_use_ebt_struct) then
!$OMP do
do K=2,nz+1 ; do J=js-1,je ; do i=is,ie
KH_v(i,J,K) = KH_v(i,J,1) * 0.5 * ( VarMix%ebt_struct(i,j,k-1) + VarMix%ebt_struct(i,j+1,k-1) )
enddo ; enddo ; enddo
if (CS%full_depth_khth_min) then
!$OMP do
do K=2,nz+1 ; do J=js-1,je ; do i=is,ie
KH_v(i,J,K) = KH_v(i,J,1) * 0.5 * ( VarMix%ebt_struct(i,j,k-1) + VarMix%ebt_struct(i,j+1,k-1) )
KH_v(i,J,K) = max(KH_v(i,J,K), CS%Khth_Min)
enddo ; enddo ; enddo
else
!$OMP do
do K=2,nz+1 ; do J=js-1,je ; do i=is,ie
KH_v(i,J,K) = KH_v(i,J,1) * 0.5 * ( VarMix%ebt_struct(i,j,k-1) + VarMix%ebt_struct(i,j+1,k-1) )
enddo ; enddo ; enddo
endif
else
!$OMP do
do K=2,nz+1 ; do J=js-1,je ; do i=is,ie
Expand Down Expand Up @@ -2132,7 +2151,11 @@ subroutine thickness_diffuse_init(Time, G, GV, US, param_file, diag, CDp, CS)
! rotation [nondim].
real :: Stanley_coeff ! Coefficient relating the temperature gradient and sub-gridscale
! temperature variance [nondim]
integer :: default_answer_date ! The default setting for the various ANSWER_DATE flags.
logical :: khth_use_ebt_struct ! If true, uses the equivalent barotropic structure
! as the vertical structure of thickness diffusivity.
! Used to determine if FULL_DEPTH_KHTH_MIN should be
! available.
integer :: default_answer_date ! The default setting for the various ANSWER_DATE flags.
integer :: i, j

CS%initialized = .true.
Expand Down Expand Up @@ -2178,6 +2201,17 @@ subroutine thickness_diffuse_init(Time, G, GV, US, param_file, diag, CDp, CS)
call get_param(param_file, mdl, "KHTH_MIN", CS%KHTH_Min, &
"The minimum horizontal thickness diffusivity.", &
default=0.0, units="m2 s-1", scale=US%m_to_L**2*US%T_to_s)
call get_param(param_file, mdl, "KHTH_USE_EBT_STRUCT", khth_use_ebt_struct, &
"If true, uses the equivalent barotropic structure "//&
"as the vertical structure of thickness diffusivity.",&
default=.false., do_not_log=.true.)
if (khth_use_ebt_struct .and. CS%KHTH_Min>0.0) then
call get_param(param_file, mdl, "FULL_DEPTH_KHTH_MIN", CS%full_depth_khth_min, &
"If true, KHTH_MIN is enforced throughout the whole water column. "//&
"Otherwise, KHTH_MIN is only enforced at the surface. This parameter "//&
"is only available when KHTH_USE_EBT_STRUCT=True and KHTH_MIN>0.", &
default=.false.)
endif
call get_param(param_file, mdl, "KHTH_MAX", CS%KHTH_Max, &
"The maximum horizontal thickness diffusivity.", &
default=0.0, units="m2 s-1", scale=US%m_to_L**2*US%T_to_s)
Expand Down
63 changes: 46 additions & 17 deletions src/tracer/MOM_tracer_hor_diff.F90
Original file line number Diff line number Diff line change
Expand Up @@ -52,8 +52,11 @@ module MOM_tracer_hor_diff
real :: max_diff_CFL !< If positive, locally limit the along-isopycnal
!! tracer diffusivity to keep the diffusive CFL
!! locally at or below this value [nondim].
logical :: KhTh_use_ebt_struct !< If true, uses the equivalent barotropic structure
logical :: KhTr_use_ebt_struct !< If true, uses the equivalent barotropic structure
!! as the vertical structure of tracer diffusivity.
logical :: full_depth_khtr_min !< If true, KHTR_MIN is enforced throughout the whole water column.
!! Otherwise, KHTR_MIN is only enforced at the surface. This parameter
!! is only available when KHTR_USE_EBT_STRUCT=True and KHTR_MIN>0.
logical :: Diffuse_ML_interior !< If true, diffuse along isopycnals between
!! the mixed layer and the interior.
logical :: check_diffusive_CFL !< If true, automatically iterate the diffusion
Expand Down Expand Up @@ -412,21 +415,40 @@ subroutine tracer_hordiff(h, dt, MEKE, VarMix, G, GV, US, CS, Reg, tv, do_online
enddo
enddo
enddo
if (CS%KhTh_use_ebt_struct) then
do K=2,nz+1
do J=js-1,je
do i=is,ie
Coef_y(i,J,K) = Coef_y(i,J,1) * 0.5 * ( VarMix%ebt_struct(i,j,k-1) + VarMix%ebt_struct(i,j+1,k-1) )
if (CS%KhTr_use_ebt_struct) then
if (CS%full_depth_khtr_min) then
do K=2,nz+1
do J=js-1,je
do i=is,ie
Coef_y(i,J,K) = Coef_y(i,J,1) * 0.5 * ( VarMix%ebt_struct(i,j,k-1) + VarMix%ebt_struct(i,j+1,k-1) )
Coef_y(i,J,K) = max(Coef_y(i,J,K), CS%KhTr_min)
enddo
enddo
enddo
enddo
do k=2,nz+1
do j=js,je
do I=is-1,ie
Coef_x(I,j,K) = Coef_x(I,j,1) * 0.5 * ( VarMix%ebt_struct(i,j,k-1) + VarMix%ebt_struct(i+1,j,k-1) )
do k=2,nz+1
do j=js,je
do I=is-1,ie
Coef_x(I,j,K) = Coef_x(I,j,1) * 0.5 * ( VarMix%ebt_struct(i,j,k-1) + VarMix%ebt_struct(i+1,j,k-1) )
Coef_x(I,j,K) = max(Coef_x(I,j,K), CS%KhTr_min)
enddo
enddo
enddo
enddo
else
do K=2,nz+1
do J=js-1,je
do i=is,ie
Coef_y(i,J,K) = Coef_y(i,J,1) * 0.5 * ( VarMix%ebt_struct(i,j,k-1) + VarMix%ebt_struct(i,j+1,k-1) )
enddo
enddo
enddo
do k=2,nz+1
do j=js,je
do I=is-1,ie
Coef_x(I,j,K) = Coef_x(I,j,1) * 0.5 * ( VarMix%ebt_struct(i,j,k-1) + VarMix%ebt_struct(i+1,j,k-1) )
enddo
enddo
enddo
endif
endif

do itt=1,num_itts
Expand Down Expand Up @@ -468,7 +490,7 @@ subroutine tracer_hordiff(h, dt, MEKE, VarMix, G, GV, US, CS, Reg, tv, do_online
enddo
enddo
enddo
if (CS%KhTh_use_ebt_struct) then
if (CS%KhTr_use_ebt_struct) then
do K=2,nz+1
do J=js-1,je
do i=is,ie
Expand Down Expand Up @@ -594,7 +616,7 @@ subroutine tracer_hordiff(h, dt, MEKE, VarMix, G, GV, US, CS, Reg, tv, do_online
do j=js,je ; do I=is-1,ie
Kh_u(I,j,:) = G%mask2dCu(I,j)*Kh_u(I,j,1)
enddo ; enddo
if (CS%KhTh_use_ebt_struct) then
if (CS%KhTr_use_ebt_struct) then
do K=2,nz+1
do j=js,je
do I=is-1,ie
Expand All @@ -610,7 +632,7 @@ subroutine tracer_hordiff(h, dt, MEKE, VarMix, G, GV, US, CS, Reg, tv, do_online
do J=js-1,je ; do i=is,ie
Kh_v(i,J,:) = G%mask2dCv(i,J)*Kh_v(i,J,1)
enddo ; enddo
if (CS%KhTh_use_ebt_struct) then
if (CS%KhTr_use_ebt_struct) then
do K=2,nz+1
do J=js-1,je
do i=is,ie
Expand All @@ -636,7 +658,7 @@ subroutine tracer_hordiff(h, dt, MEKE, VarMix, G, GV, US, CS, Reg, tv, do_online
(G%mask2dCv(i,J-1)+G%mask2dCv(i,J)) + 1.0e-37)
Kh_h(i,j,:) = normalize*G%mask2dT(i,j)*((Kh_u(I-1,j,1)+Kh_u(I,j,1)) + &
(Kh_v(i,J-1,1)+Kh_v(i,J,1)))
if (CS%KhTh_use_ebt_struct) then
if (CS%KhTr_use_ebt_struct) then
do K=2,nz+1
Kh_h(i,j,K) = normalize*G%mask2dT(i,j)*VarMix%ebt_struct(i,j,k-1)*((Kh_u(I-1,j,1)+Kh_u(I,j,1)) + &
(Kh_v(i,J-1,1)+Kh_v(i,J,1)))
Expand Down Expand Up @@ -1561,7 +1583,7 @@ subroutine tracer_hor_diff_init(Time, G, GV, US, param_file, diag, EOS, diabatic
call get_param(param_file, mdl, "KHTR", CS%KhTr, &
"The background along-isopycnal tracer diffusivity.", &
units="m2 s-1", default=0.0, scale=US%m_to_L**2*US%T_to_s)
call get_param(param_file, mdl, "KHTR_USE_EBT_STRUCT", CS%KhTh_use_ebt_struct, &
call get_param(param_file, mdl, "KHTR_USE_EBT_STRUCT", CS%KhTr_use_ebt_struct, &
"If true, uses the equivalent barotropic structure "//&
"as the vertical structure of the tracer diffusivity.",&
default=.false.)
Expand All @@ -1573,6 +1595,13 @@ subroutine tracer_hor_diff_init(Time, G, GV, US, param_file, diag, EOS, diabatic
call get_param(param_file, mdl, "KHTR_MIN", CS%KhTr_Min, &
"The minimum along-isopycnal tracer diffusivity.", &
units="m2 s-1", default=0.0, scale=US%m_to_L**2*US%T_to_s)
if (CS%KhTr_use_ebt_struct .and. CS%KhTr_Min > 0.0) then
call get_param(param_file, mdl, "FULL_DEPTH_KHTR_MIN", CS%full_depth_khtr_min, &
"If true, KHTR_MIN is enforced throughout the whole water column. "//&
"Otherwise, KHTR_MIN is only enforced at the surface. This parameter "//&
"is only available when KHTR_USE_EBT_STRUCT=True and KHTR_MIN>0.", &
default=.false.)
endif
call get_param(param_file, mdl, "KHTR_MAX", CS%KhTr_Max, &
"The maximum along-isopycnal tracer diffusivity.", &
units="m2 s-1", default=0.0, scale=US%m_to_L**2*US%T_to_s)
Expand Down
Loading