diff --git a/src/core/MOM_variables.F90 b/src/core/MOM_variables.F90 index efd558cdb3..cabda1cf58 100644 --- a/src/core/MOM_variables.F90 +++ b/src/core/MOM_variables.F90 @@ -239,6 +239,9 @@ module MOM_variables !! h points [H T-1 ~> m s-1 or kg m-2 s-1]. BBL_meanKE_loss, & !< The viscous loss of mean kinetic energy in the bottom boundary layer !! [H L2 T-3 ~> m3 s-3 or W m-2]. + BBL_meanKE_loss_sqrtCd, & !< The viscous loss of mean kinetic energy in the bottom boundary layer + !! divided by the square root of the drag coefficient [H L2 T-3 ~> m3 s-3 or W m-2]. + !! This is being set only to retain old answers, and should be phased out. taux_shelf, & !< The zonal stresses on the ocean under shelves [R Z L T-2 ~> Pa]. tauy_shelf !< The meridional stresses on the ocean under shelves [R Z L T-2 ~> Pa]. real, allocatable, dimension(:,:) :: tbl_thick_shelf_u diff --git a/src/parameterizations/vertical/MOM_energetic_PBL.F90 b/src/parameterizations/vertical/MOM_energetic_PBL.F90 index 6e3cdfa1a5..6caaace9e9 100644 --- a/src/parameterizations/vertical/MOM_energetic_PBL.F90 +++ b/src/parameterizations/vertical/MOM_energetic_PBL.F90 @@ -198,6 +198,9 @@ module MOM_energetic_PBL !! energetics that accounts for an exponential decay of TKE from a !! near-bottom source and an assumed piecewise linear linear profile !! of the buoyancy flux response to a change in a diffusivity. + logical :: BBL_effic_bug !< If true, overestimate the efficiency of the non-tidal ePBL bottom boundary + !! layer diffusivity by a factor of 1/sqrt(CDRAG), which is often a factor of + !! about 18.3. !/ Options for documenting differences from parameter choices integer :: options_diff !< If positive, this is a coded integer indicating a pair of @@ -662,7 +665,11 @@ subroutine energetic_PBL(h_3d, u_3d, v_3d, tv, fluxes, visc, dt, Kd_int, G, GV, if (BBL_mixing) then if (CS%MLD_iteration_guess .and. (CS%BBL_depth(i,j) > 0.0)) BBLD_io = CS%BBL_depth(i,j) BBLD_in = BBLD_io - BBL_TKE = CS%ePBL_BBL_effic * GV%H_to_RZ * dt * visc%BBL_meanKE_loss(i,j) + if (CS%BBL_effic_bug) then + BBL_TKE = CS%ePBL_BBL_effic * GV%H_to_RZ * dt * visc%BBL_meanKE_loss_sqrtCd(i,j) + else + BBL_TKE = CS%ePBL_BBL_effic * GV%H_to_RZ * dt * visc%BBL_meanKE_loss(i,j) + endif u_star_BBL = max(visc%ustar_BBL(i,j), CS%ustar_min*GV%Z_to_H) ! Add in tidal dissipation energy at the bottom, noting that fluxes%BBL_tidal_dis is @@ -3752,6 +3759,10 @@ subroutine energetic_PBL_init(Time, G, GV, US, param_file, diag, CS) "length scale by rotation in the bottom boundary layer. Making this larger "//& "decreases the bottom boundary layer diffusivity.", & units="nondim", default=CS%Ekman_scale_coef, do_not_log=no_BBL) + call get_param(param_file, mdl, "EPBL_BBL_EFFIC_BUG", CS%BBL_effic_bug, & + "If true, overestimate the efficiency of the non-tidal ePBL bottom boundary "//& + "layer diffusivity by a factor of 1/sqrt(CDRAG), which is often a factor of "//& + "about 18.3.", default=.false., do_not_log=(CS%ePBL_BBL_effic<=0.0)) call get_param(param_file, mdl, "DECAY_ADJUSTED_BBL_TKE", CS%decay_adjusted_BBL_TKE, & "If true, include an adjustment factor in the bottom boundary layer energetics "//& diff --git a/src/parameterizations/vertical/MOM_set_diffusivity.F90 b/src/parameterizations/vertical/MOM_set_diffusivity.F90 index 056a89a190..8035b88b96 100644 --- a/src/parameterizations/vertical/MOM_set_diffusivity.F90 +++ b/src/parameterizations/vertical/MOM_set_diffusivity.F90 @@ -163,11 +163,17 @@ module MOM_set_diffusivity !! calculations. Values below 20190101 recover the answers from the !! end of 2018, while higher values use updated and more robust forms !! of the same expressions. Values above 20240630 use more accurate - !! expressions for cases where USE_LOTW_BBL_DIFFUSIVITY is true. + !! expressions for cases where USE_LOTW_BBL_DIFFUSIVITY is true. Values + !! above 20250301 use less confusing expressions to set the bottom-drag + !! generated diffusivity when USE_LOTW_BBL_DIFFUSIVITY is false. integer :: LOTW_BBL_answer_date !< The vintage of the order of arithmetic and expressions !! in the LOTW_BBL calculations. Values below 20240630 recover the !! original answers, while higher values use more accurate expressions. !! This only applies when USE_LOTW_BBL_DIFFUSIVITY is true. + integer :: drag_diff_answer_date !< The vintage of the order of arithmetic in the drag diffusivity + !! calculations. Values above 20250301 use less confusing expressions + !! to set the bottom-drag generated diffusivity when + !! USE_LOTW_BBL_DIFFUSIVITY is false. character(len=200) :: inputdir !< The directory in which input files are found type(user_change_diff_CS), pointer :: user_change_diff_CSp => NULL() !< Control structure for a child module @@ -1414,7 +1420,11 @@ subroutine add_drag_diffusivity(h, u, v, tv, fluxes, visc, j, TKE_to_Kd, maxTKE, ! If ustar_h = 0, this is land so this value doesn't matter. I2decay(i) = 0.5*CS%IMax_decay endif - TKE(i) = ((CS%BBL_effic * cdrag_sqrt) * exp(-I2decay(i)*h(i,j,nz)) ) * visc%BBL_meanKE_loss(i,j) + if (CS%drag_diff_answer_date <= 20250301) then + TKE(i) = ((CS%BBL_effic * cdrag_sqrt) * exp(-I2decay(i)*h(i,j,nz)) ) * visc%BBL_meanKE_loss_sqrtCd(i,j) + else + TKE(i) = (CS%BBL_effic * exp(-I2decay(i)*h(i,j,nz)) ) * visc%BBL_meanKE_loss(i,j) + endif if (associated(fluxes%BBL_tidal_dis)) & TKE(i) = TKE(i) + fluxes%BBL_tidal_dis(i,j) * GV%RZ_to_H * & @@ -1644,8 +1654,7 @@ subroutine add_LOTW_BBL_diffusivity(h, u, v, tv, fluxes, visc, j, N2_int, Rho_bo ! Energy input at the bottom [H Z2 T-3 ~> m3 s-3 or W m-2]. ! (Note that visc%BBL_meanKE_loss is in [H L2 T-3 ~> m3 s-3 or W m-2], set in set_BBL_TKE().) - !### I am still unsure about sqrt(cdrag) in this expressions - AJA - BBL_meanKE_dis = cdrag_sqrt * visc%BBL_meanKE_loss(i,j) + BBL_meanKE_dis = visc%BBL_meanKE_loss(i,j) ! Add in tidal dissipation energy at the bottom [H Z2 T-3 ~> m3 s-3 or W m-2]. ! Note that BBL_tidal_dis is in [R Z L2 T-3 ~> W m-2]. if (associated(fluxes%BBL_tidal_dis)) & @@ -1956,6 +1965,9 @@ subroutine set_BBL_TKE(u, v, h, tv, fluxes, visc, G, GV, US, CS, OBC) if (allocated(visc%BBL_meanKE_loss)) then do j=js,je ; do i=is,ie ; visc%BBL_meanKE_loss(i,j) = 0.0 ; enddo ; enddo endif + if (allocated(visc%BBL_meanKE_loss_sqrtCd)) then + do j=js,je ; do i=is,ie ; visc%BBL_meanKE_loss_sqrtCd(i,j) = 0.0 ; enddo ; enddo + endif return endif @@ -2079,7 +2091,13 @@ subroutine set_BBL_TKE(u, v, h, tv, fluxes, visc, G, GV, US, CS, OBC) (G%areaCu(I,j)*(ustar(I)*ustar(I)))) + & ((G%areaCv(i,J-1)*(vstar(i,J-1)*vstar(i,J-1))) + & (G%areaCv(i,J)*(vstar(i,J)*vstar(i,J)))) ) ) - visc%BBL_meanKE_loss(i,j) = & + visc%BBL_meanKE_loss(i,j) = cdrag_sqrt * & + ((((G%areaCu(I-1,j)*(ustar(I-1)*u2_bbl(I-1))) + & + (G%areaCu(I,j) * (ustar(I)*u2_bbl(I)))) + & + ((G%areaCv(i,J-1)*(vstar(i,J-1)*v2_bbl(i,J-1))) + & + (G%areaCv(i,J) * (vstar(i,J)*v2_bbl(i,J)))) )*G%IareaT(i,j)) + ! The following line could be omitted if SET_DIFF_ANSWER_DATE > 20250301 and EPBL_BBL_EFFIC_BUG is false. + visc%BBL_meanKE_loss_sqrtCd(i,j) = & ((((G%areaCu(I-1,j)*(ustar(I-1)*u2_bbl(I-1))) + & (G%areaCu(I,j) * (ustar(I)*u2_bbl(I)))) + & ((G%areaCv(i,J-1)*(vstar(i,J-1)*v2_bbl(i,J-1))) + & @@ -2289,7 +2307,9 @@ subroutine set_diffusivity_init(Time, G, GV, US, param_file, diag, CS, int_tide_ call get_param(param_file, mdl, "SET_DIFF_ANSWER_DATE", CS%answer_date, & "The vintage of the order of arithmetic and expressions in the set diffusivity "//& "calculations. Values below 20190101 recover the answers from the end of 2018, "//& - "while higher values use updated and more robust forms of the same expressions.", & + "while higher values use updated and more robust forms of the same expressions. "//& + "Values above 20250301 also use less confusing expressions to set the bottom-drag "//& + "generated diffusivity when USE_LOTW_BBL_DIFFUSIVITY is false.", & default=default_answer_date, do_not_log=.not.GV%Boussinesq) if (.not.GV%Boussinesq) CS%answer_date = max(CS%answer_date, 20230701) @@ -2405,6 +2425,12 @@ subroutine set_diffusivity_init(Time, G, GV, US, param_file, diag, CS, int_tide_ "USE_LOTW_BBL_DIFFUSIVITY is true.", & default=20190101, do_not_log=.not.CS%use_LOTW_BBL_diffusivity) !### Set default as default=default_answer_date, or use SET_DIFF_ANSWER_DATE. + call get_param(param_file, mdl, "DRAG_DIFFUSIVITY_ANSWER_DATE", CS%drag_diff_answer_date, & + "The vintage of the order of arithmetic in the drag diffusivity calculations. "//& + "Values above 20250301 use less confusing expressions to set the bottom-drag "//& + "generated diffusivity when USE_LOTW_BBL_DIFFUSIVITY is false. ", & + default=20250101, do_not_log=CS%use_LOTW_BBL_diffusivity.or.(CS%BBL_effic<=0.0)) + !### Set default as default=default_answer_date, or use SET_DIFF_ANSWER_DATE. CS%id_Kd_BBL = register_diag_field('ocean_model', 'Kd_BBL', diag%axesTi, Time, & 'Bottom Boundary Layer Diffusivity', 'm2 s-1', conversion=GV%HZ_T_to_m2_s) diff --git a/src/parameterizations/vertical/MOM_set_viscosity.F90 b/src/parameterizations/vertical/MOM_set_viscosity.F90 index 963c8c4525..c6a389e158 100644 --- a/src/parameterizations/vertical/MOM_set_viscosity.F90 +++ b/src/parameterizations/vertical/MOM_set_viscosity.F90 @@ -3140,6 +3140,7 @@ subroutine set_visc_init(Time, G, GV, US, param_file, diag, visc, CS, restart_CS allocate(visc%kv_bbl_v(isd:ied,JsdB:JedB), source=0.0) allocate(visc%ustar_bbl(isd:ied,jsd:jed), source=0.0) allocate(visc%BBL_meanKE_loss(isd:ied,jsd:jed), source=0.0) + allocate(visc%BBL_meanKE_loss_sqrtCd(isd:ied,jsd:jed), source=0.0) CS%id_bbl_thick_u = register_diag_field('ocean_model', 'bbl_thick_u', & diag%axesCu1, Time, 'BBL thickness at u points', 'm', conversion=US%Z_to_m) @@ -3215,6 +3216,7 @@ subroutine set_visc_end(visc, CS) if (associated(visc%Kv_shear_Bu)) deallocate(visc%Kv_shear_Bu) if (allocated(visc%ustar_bbl)) deallocate(visc%ustar_bbl) if (allocated(visc%BBL_meanKE_loss)) deallocate(visc%BBL_meanKE_loss) + if (allocated(visc%BBL_meanKE_loss_sqrtCd)) deallocate(visc%BBL_meanKE_loss_sqrtCd) if (allocated(visc%taux_shelf)) deallocate(visc%taux_shelf) if (allocated(visc%tauy_shelf)) deallocate(visc%tauy_shelf) if (allocated(visc%tbl_thick_shelf_u)) deallocate(visc%tbl_thick_shelf_u)