From b97034eee98b6bd1a3d2dac0a3b6ebe56cd02443 Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Thu, 25 Apr 2024 16:49:44 -0400 Subject: [PATCH] +*Pass visc%h_ML to call_tracer_column_fns Added the optional argument h_BL to call_tracer_column_fns, and then pass visc%h_ML to this routine as a part of the calls from the diabatic routines. This argument is not provided in adiabatic mode or in offline tracer mode, which is why it has been made optional. All answers are bitwise identical in Boussinesq mode but ideal age tracers can change in non-Boussinesq mode due to the updates to the layer specific volumes between the calculation of the boundary layer properties and the calls to call_tracer_column_fns. --- .../vertical/MOM_diabatic_driver.F90 | 13 +++++-------- .../vertical/MOM_set_viscosity.F90 | 6 ++++-- src/tracer/MOM_tracer_flow_control.F90 | 18 +++++++++++++++--- 3 files changed, 24 insertions(+), 13 deletions(-) diff --git a/src/parameterizations/vertical/MOM_diabatic_driver.F90 b/src/parameterizations/vertical/MOM_diabatic_driver.F90 index f967c005e6..a384f3b327 100644 --- a/src/parameterizations/vertical/MOM_diabatic_driver.F90 +++ b/src/parameterizations/vertical/MOM_diabatic_driver.F90 @@ -1078,7 +1078,7 @@ subroutine diabatic_ALE_legacy(u, v, h, tv, BLD, fluxes, visc, ADp, CDp, dt, Tim KPP_CSp=CS%KPP_CSp, & nonLocalTrans=CS%KPP_NLTscalar, & evap_CFL_limit=CS%evap_CFL_limit, & - minimum_forcing_depth=CS%minimum_forcing_depth) + minimum_forcing_depth=CS%minimum_forcing_depth, h_BL=visc%h_ML) call cpu_clock_end(id_clock_tracers) @@ -1588,7 +1588,7 @@ subroutine diabatic_ALE(u, v, h, tv, BLD, fluxes, visc, ADp, CDp, dt, Time_end, KPP_CSp=CS%KPP_CSp, & nonLocalTrans=CS%KPP_NLTscalar, & evap_CFL_limit=CS%evap_CFL_limit, & - minimum_forcing_depth=CS%minimum_forcing_depth) + minimum_forcing_depth=CS%minimum_forcing_depth, h_BL=visc%h_ML) call cpu_clock_end(id_clock_tracers) @@ -2378,8 +2378,7 @@ subroutine layered_diabatic(u, v, h, tv, BLD, fluxes, visc, ADp, CDp, dt, Time_e call call_tracer_column_fns(hold, h, eatr, ebtr, fluxes, BLD, dt, G, GV, US, tv, & CS%optics, CS%tracer_flow_CSp, CS%debug, & - KPP_CSp=CS%KPP_CSp, & - nonLocalTrans=CS%KPP_NLTscalar) + KPP_CSp=CS%KPP_CSp, nonLocalTrans=CS%KPP_NLTscalar, h_BL=visc%h_ML) elseif (CS%double_diffuse) then ! extra diffusivity for passive tracers @@ -2400,14 +2399,12 @@ subroutine layered_diabatic(u, v, h, tv, BLD, fluxes, visc, ADp, CDp, dt, Time_e call call_tracer_column_fns(hold, h, eatr, ebtr, fluxes, BLD, dt, G, GV, US, tv, & CS%optics, CS%tracer_flow_CSp, CS%debug, & - KPP_CSp=CS%KPP_CSp, & - nonLocalTrans=CS%KPP_NLTscalar) + KPP_CSp=CS%KPP_CSp, nonLocalTrans=CS%KPP_NLTscalar, h_BL=visc%h_ML) else call call_tracer_column_fns(hold, h, ea, eb, fluxes, BLD, dt, G, GV, US, tv, & CS%optics, CS%tracer_flow_CSp, CS%debug, & - KPP_CSp=CS%KPP_CSp, & - nonLocalTrans=CS%KPP_NLTscalar) + KPP_CSp=CS%KPP_CSp, nonLocalTrans=CS%KPP_NLTscalar, h_BL=visc%h_ML) endif ! (CS%mix_boundary_tracers) diff --git a/src/parameterizations/vertical/MOM_set_viscosity.F90 b/src/parameterizations/vertical/MOM_set_viscosity.F90 index b2821a4ea8..59f869f458 100644 --- a/src/parameterizations/vertical/MOM_set_viscosity.F90 +++ b/src/parameterizations/vertical/MOM_set_viscosity.F90 @@ -2679,7 +2679,7 @@ subroutine set_visc_register_restarts(HI, G, GV, US, param_file, visc, restart_C logical, intent(in) :: use_ice_shelf !< if true, register tau_shelf restarts ! Local variables logical :: use_kappa_shear, KS_at_vertex - logical :: adiabatic, useKPP, useEPBL + logical :: adiabatic, useKPP, useEPBL, use_ideal_age logical :: do_brine_plume, use_hor_bnd_diff, use_neutral_diffusion, use_fpmix logical :: use_CVMix_shear, MLE_use_PBL_MLD, MLE_use_Bodner, use_CVMix_conv integer :: isd, ied, jsd, jed, nz @@ -2763,12 +2763,14 @@ subroutine set_visc_register_restarts(HI, G, GV, US, param_file, visc, restart_C default=.false., do_not_log=.true.) call get_param(param_file, mdl, "FPMIX", use_fpmix, & default=.false., do_not_log=.true.) + call get_param(param_file, mdl, "USE_IDEAL_AGE_TRACER", use_ideal_age, & + default=.false., do_not_log=.true.) if (MLE_use_PBL_MLD) then call safe_alloc_ptr(visc%MLD, isd, ied, jsd, jed) endif if ((hfreeze >= 0.0) .or. MLE_use_PBL_MLD .or. do_brine_plume .or. use_fpmix .or. & - use_neutral_diffusion .or. use_hor_bnd_diff) then + use_neutral_diffusion .or. use_hor_bnd_diff .or. use_ideal_age) then call safe_alloc_ptr(visc%h_ML, isd, ied, jsd, jed) endif diff --git a/src/tracer/MOM_tracer_flow_control.F90 b/src/tracer/MOM_tracer_flow_control.F90 index 10aba675da..ca85fc234f 100644 --- a/src/tracer/MOM_tracer_flow_control.F90 +++ b/src/tracer/MOM_tracer_flow_control.F90 @@ -429,7 +429,7 @@ end subroutine call_tracer_set_forcing !> This subroutine calls all registered tracer column physics subroutines. subroutine call_tracer_column_fns(h_old, h_new, ea, eb, fluxes, mld, dt, G, GV, US, tv, optics, CS, & - debug, KPP_CSp, nonLocalTrans, evap_CFL_limit, minimum_forcing_depth) + debug, KPP_CSp, nonLocalTrans, evap_CFL_limit, minimum_forcing_depth, h_BL) type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure. type(verticalGrid_type), intent(in) :: GV !< The ocean's vertical grid structure. real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), intent(in) :: h_old !< Layer thickness before entrainment @@ -464,9 +464,11 @@ subroutine call_tracer_column_fns(h_old, h_new, ea, eb, fluxes, mld, dt, G, GV, !! of the top layer in a timestep [nondim] real, optional, intent(in) :: minimum_forcing_depth !< The smallest depth over !! which fluxes can be applied [H ~> m or kg m-2] + real, dimension(:,:), optional, pointer :: h_BL !< Thickness of active mixing layer [H ~> m or kg m-2] ! Local variables real :: Hbl(SZI_(G),SZJ_(G)) !< Boundary layer thickness [H ~> m or kg m-2] + logical :: use_h_BL if (.not. associated(CS)) call MOM_error(FATAL, "call_tracer_column_fns: "// & "Module must be initialized via call_tracer_register before it is used.") @@ -493,7 +495,12 @@ subroutine call_tracer_column_fns(h_old, h_new, ea, eb, fluxes, mld, dt, G, GV, evap_CFL_limit=evap_CFL_limit, & minimum_forcing_depth=minimum_forcing_depth) if (CS%use_ideal_age) then - call convert_MLD_to_ML_thickness(mld, h_new, Hbl, tv, G, GV) + use_h_BL = .false. ; if (present(h_BL)) use_h_BL = associated(h_BL) + if (present(h_BL)) then + Hbl(:,:) = h_BL(:,:) + else ! This option is here mostly to support the offline tracers. + call convert_MLD_to_ML_thickness(mld, h_new, Hbl, tv, G, GV) + endif call ideal_age_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, & G, GV, US, CS%ideal_age_tracer_CSp, & evap_CFL_limit=evap_CFL_limit, & @@ -573,7 +580,12 @@ subroutine call_tracer_column_fns(h_old, h_new, ea, eb, fluxes, mld, dt, G, GV, call RGC_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, & G, GV, US, CS%RGC_tracer_CSp) if (CS%use_ideal_age) then - call convert_MLD_to_ML_thickness(mld, h_new, Hbl, tv, G, GV) + use_h_BL = .false. ; if (present(h_BL)) use_h_BL = associated(h_BL) + if (present(h_BL)) then + Hbl(:,:) = h_BL(:,:) + else ! This option is here mostly to support the offline tracers. + call convert_MLD_to_ML_thickness(mld, h_new, Hbl, tv, G, GV) + endif call ideal_age_tracer_column_physics(h_old, h_new, ea, eb, fluxes, dt, & G, GV, US, CS%ideal_age_tracer_CSp, Hbl=Hbl) endif