From 70bfe8f0c48b22857db358f3033215ce5f0603a2 Mon Sep 17 00:00:00 2001 From: Kate Hedstrom Date: Sun, 2 Jun 2024 01:55:12 -0800 Subject: [PATCH] Spatially variable fields for MLE%Bodner (#617) * Spatially variable fields for MLE%Bodner - Allows reading in 2d fields for Cr and for MLD_decaying_Tfilt. * Finish the job? * Renamed one variable, fixed some units --- .../lateral/MOM_mixed_layer_restrat.F90 | 61 ++++++++++++++++--- 1 file changed, 52 insertions(+), 9 deletions(-) diff --git a/src/parameterizations/lateral/MOM_mixed_layer_restrat.F90 b/src/parameterizations/lateral/MOM_mixed_layer_restrat.F90 index e7ada31430..9f0061f104 100644 --- a/src/parameterizations/lateral/MOM_mixed_layer_restrat.F90 +++ b/src/parameterizations/lateral/MOM_mixed_layer_restrat.F90 @@ -15,6 +15,7 @@ module MOM_mixed_layer_restrat use MOM_grid, only : ocean_grid_type use MOM_hor_index, only : hor_index_type use MOM_intrinsic_functions, only : cuberoot +use MOM_io, only : MOM_read_data use MOM_lateral_mixing_coeffs, only : VarMix_CS use MOM_restart, only : register_restart_field, query_initialized, MOM_restart_CS use MOM_unit_scaling, only : unit_scale_type @@ -98,11 +99,15 @@ module MOM_mixed_layer_restrat real :: Kv_restrat !< A viscosity that sets a floor on the momentum mixing rate !! during restratification, rescaled into thickness-based !! units [H2 T-1 ~> m2 s-1 or kg2 m-4 s-1] + logical :: MLD_grid !< If true, read a spacially varying field for MLD_decaying_Tfilt + logical :: Cr_grid !< If true, read a spacially varying field for Cr real, dimension(:,:), allocatable :: & MLD_filtered, & !< Time-filtered MLD [H ~> m or kg m-2] MLD_filtered_slow, & !< Slower time-filtered MLD [H ~> m or kg m-2] - wpup_filtered !< Time-filtered vertical momentum flux [H L T-2 ~> m2 s-2 or kg m-1 s-2] + wpup_filtered, & !< Time-filtered vertical momentum flux [H L T-2 ~> m2 s-2 or kg m-1 s-2] + MLD_Tfilt_space, & !< Spatially varying time scale for MLD filter [T ~> s] + Cr_space !< Spatially varying Cr coefficient [nondim] !>@{ !! Diagnostic identifier @@ -887,11 +892,19 @@ subroutine mixedlayer_restrat_Bodner(CS, G, GV, US, h, uhtr, vhtr, tv, forces, d enddo ; enddo ! Calculate "big H", representative of the mixed layer depth, used in B22 formula (eq 27). - do j=js-1,je+1 ; do i=is-1,ie+1 - big_H(i,j) = rmean2ts(little_h(i,j), CS%MLD_filtered_slow(i,j), & - CS%MLD_growing_Tfilt, CS%MLD_decaying_Tfilt, dt) - CS%MLD_filtered_slow(i,j) = big_H(i,j) - enddo ; enddo + if (CS%MLD_grid) then + do j=js-1,je+1 ; do i=is-1,ie+1 + big_H(i,j) = rmean2ts(little_h(i,j), CS%MLD_filtered_slow(i,j), & + CS%MLD_growing_Tfilt, CS%MLD_Tfilt_space(i,j), dt) + CS%MLD_filtered_slow(i,j) = big_H(i,j) + enddo ; enddo + else + do j=js-1,je+1 ; do i=is-1,ie+1 + big_H(i,j) = rmean2ts(little_h(i,j), CS%MLD_filtered_slow(i,j), & + CS%MLD_growing_Tfilt, CS%MLD_decaying_Tfilt, dt) + CS%MLD_filtered_slow(i,j) = big_H(i,j) + enddo ; enddo + endif ! Estimate w'u' at h-points, with a floor to avoid division by zero later. if (allocated(tv%SpV_avg) .and. .not.(GV%Boussinesq .or. GV%semi_Boussinesq)) then @@ -1050,7 +1063,7 @@ subroutine mixedlayer_restrat_Bodner(CS, G, GV, US, h, uhtr, vhtr, tv, forces, d h_big = 0.5*( big_H(i,j) + big_H(i+1,j) ) ! H ~> m or kg m-3 grd_b = ( buoy_av(i+1,j) - buoy_av(i,j) ) * G%IdxCu(I,j) ! L H-1 T-2 ~> s-2 or m3 kg-1 s-2 r_wpup = 2. / ( wpup(i,j) + wpup(i+1,j) ) ! T2 L-1 H-1 ~> s2 m-2 or m s2 kg-1 - psi_mag = ( ( ( CS%Cr * grid_dsd ) * ( absf * h_sml ) ) & ! L2 H T-1 ~> m3 s-1 or kg s-1 + psi_mag = ( ( ( CS%Cr_space(i,j) * grid_dsd ) * ( absf * h_sml ) ) & ! L2 H T-1 ~> m3 s-1 or kg s-1 * ( ( h_big**2 ) * grd_b ) ) * r_wpup else ! There is no flux on land and no gradient at open boundary points. psi_mag = 0.0 @@ -1091,7 +1104,7 @@ subroutine mixedlayer_restrat_Bodner(CS, G, GV, US, h, uhtr, vhtr, tv, forces, d h_big = 0.5*( big_H(i,j) + big_H(i,j+1) ) ! H ~> m or kg m-3 grd_b = ( buoy_av(i,j+1) - buoy_av(i,j) ) * G%IdyCv(I,j) ! L H-1 T-2 ~> s-2 or m3 kg-1 s-2 r_wpup = 2. / ( wpup(i,j) + wpup(i,j+1) ) ! T2 L-1 H-1 ~> s2 m-2 or m s2 kg-1 - psi_mag = ( ( ( CS%Cr * grid_dsd ) * ( absf * h_sml ) ) & ! L2 H T-1 ~> m3 s-1 or kg s-1 + psi_mag = ( ( ( CS%Cr_space(i,j) * grid_dsd ) * ( absf * h_sml ) ) & ! L2 H T-1 ~> m3 s-1 or kg s-1 * ( ( h_big**2 ) * grd_b ) ) * r_wpup else ! There is no flux on land and no gradient at open boundary points. psi_mag = 0.0 @@ -1549,6 +1562,7 @@ logical function mixedlayer_restrat_init(Time, G, GV, US, param_file, diag, CS, ! This include declares and sets the variable "version". # include "version_variable.h" integer :: i, j + character(len=200) :: filename, inputdir, varname ! Read all relevant parameters and write them to the model log. call get_param(param_file, mdl, "MIXEDLAYER_RESTRAT", mixedlayer_restrat_init, & @@ -1571,11 +1585,14 @@ logical function mixedlayer_restrat_init(Time, G, GV, US, param_file, diag, CS, CS%MLE_MLD_stretch = -9.e9 CS%use_Stanley_ML = .false. CS%use_Bodner = .false. + CS%MLD_grid = .false. + CS%Cr_grid = .false. call get_param(param_file, mdl, "DEBUG", CS%debug, default=.false., do_not_log=.true.) call get_param(param_file, mdl, "DEFAULT_ANSWER_DATE", default_answer_date, & "This sets the default value for the various _ANSWER_DATE parameters.", & default=99991231, do_not_log=.true.) + call get_param(param_file, mdl, "INPUTDIR", inputdir, default=".") call openParameterBlock(param_file,'MLE') ! Prepend MLE% to all parameters if (GV%nkml==0) then call get_param(param_file, mdl, "USE_BODNER23", CS%use_Bodner, & @@ -1584,7 +1601,7 @@ logical function mixedlayer_restrat_init(Time, G, GV, US, param_file, diag, CS, default=.false.) endif if (CS%use_Bodner) then - call get_param(param_file, mdl, "CR", CS%CR, & + call get_param(param_file, mdl, "CR", CS%Cr, & "The efficiency coefficient in eq 27 of Bodner et al., 2023.", & units="nondim", default=0.0) call get_param(param_file, mdl, "BODNER_NSTAR", CS%Nstar, & @@ -1638,6 +1655,32 @@ logical function mixedlayer_restrat_init(Time, G, GV, US, param_file, diag, CS, call get_param(param_file, mdl, "USE_STANLEY_TVAR", CS%use_Stanley_ML, & "If true, turn on Stanley SGS T variance parameterization "// & "in ML restrat code.", default=.false.) + call get_param(param_file, mdl, "USE_CR_GRID", CS%Cr_grid, & + "If true, read in a spatially varying Cr field.", default=.false.) + call get_param(param_file, mdl, "USE_MLD_GRID", CS%MLD_grid, & + "If true, read in a spatially varying MLD_decaying_Tfilt field.", default=.false.) + if (CS%MLD_grid) then + call get_param(param_file, mdl, "MLD_TFILT_FILE", filename, & + "The path to the file containing the MLD_decaying_Tfilt fields.", & + default="") + call get_param(param_file, mdl, "MLD_TFILT_VAR", varname, & + "The variable name for MLD_decaying_Tfilt field.", & + default="MLD_tfilt") + filename = trim(inputdir) // "/" // trim(filename) + allocate(CS%MLD_Tfilt_space(G%isd:G%ied,G%jsd:G%jed), source=0.0) + call MOM_read_data(filename, varname, CS%MLD_Tfilt_space, G%domain, scale=US%s_to_T) + endif + allocate(CS%Cr_space(G%isd:G%ied,G%jsd:G%jed), source=CS%Cr) + if (CS%Cr_grid) then + call get_param(param_file, mdl, "CR_FILE", filename, & + "The path to the file containing the Cr fields.", & + default="") + call get_param(param_file, mdl, "CR_VAR", varname, & + "The variable name for Cr field.", & + default="Cr") + filename = trim(inputdir) // "/" // trim(filename) + call MOM_read_data(filename, varname, CS%Cr_space, G%domain) + endif call closeParameterBlock(param_file) ! The remaining parameters do not have MLE% prepended call get_param(param_file, mdl, "MLE_USE_PBL_MLD", CS%MLE_use_PBL_MLD, & "If true, the MLE parameterization will use the mixed-layer "//&