From 6e1a2ebb03b760f86d806ce58eb6e8bdcfa77dd6 Mon Sep 17 00:00:00 2001 From: Robert Hallberg Date: Wed, 23 Nov 2022 08:22:48 -0500 Subject: [PATCH] Document units in 9 vertical param modules Documented numerous internal variables and their units in 9 vertical parameterization modules (MOM_bkgnd_mixing, MOM_bulk_mixed_layer, MOM_diabatic_aux, MOM_diabatic_driver, MOM_entrain_diffusive, MOM_kappa_shear, MOM_regularize_layers, MOM_set_diffusivity and MOM_set_viscosity). This commit includes the addition of units arguments in 9 unlogged get_param calls while 4 of these calls are now being scaled into the units used in MOM6 and then unscaled when used as the default for other parameters. A number of spelling errors were also corrected in comments. All answers and output are bitwise identical. --- .../vertical/MOM_bkgnd_mixing.F90 | 20 +- .../vertical/MOM_bulk_mixed_layer.F90 | 200 ++++++++++-------- .../vertical/MOM_diabatic_aux.F90 | 68 +++--- .../vertical/MOM_diabatic_driver.F90 | 67 +++--- .../vertical/MOM_entrain_diffusive.F90 | 143 +++++++------ .../vertical/MOM_kappa_shear.F90 | 85 ++++---- .../vertical/MOM_regularize_layers.F90 | 24 ++- .../vertical/MOM_set_diffusivity.F90 | 71 +++---- .../vertical/MOM_set_viscosity.F90 | 47 ++-- 9 files changed, 400 insertions(+), 325 deletions(-) diff --git a/src/parameterizations/vertical/MOM_bkgnd_mixing.F90 b/src/parameterizations/vertical/MOM_bkgnd_mixing.F90 index ba47f281e8..5a39f83c5d 100644 --- a/src/parameterizations/vertical/MOM_bkgnd_mixing.F90 +++ b/src/parameterizations/vertical/MOM_bkgnd_mixing.F90 @@ -57,11 +57,11 @@ module MOM_bkgnd_mixing real :: omega !< The Earth's rotation rate [T-1 ~> s-1]. real :: N0_2Omega !< ratio of the typical Buoyancy frequency to !! twice the Earth's rotation period, used with the - !! Henyey scaling from the mixing + !! Henyey scaling from the mixing [nondim] real :: prandtl_bkgnd !< Turbulent Prandtl number used to convert - !! vertical background diffusivity into viscosity + !! vertical background diffusivity into viscosity [nondim] real :: Kd_tanh_lat_scale !< A nondimensional scaling for the range of - !! diffusivities with Kd_tanh_lat_fn. Valid values + !! diffusivities with Kd_tanh_lat_fn [nondim]. Valid values !! are in the range of -2 to 2; 0.4 reproduces CM2M. real :: Kd_tot_ml !< The mixed layer diapycnal diffusivity [Z2 T-1 ~> m2 s-1] !! when no other physically based mixed layer turbulence @@ -151,10 +151,12 @@ subroutine bkgnd_mixing_init(Time, G, GV, US, param_file, diag, CS, physical_OBL CS%physical_OBL_scheme = physical_OBL_scheme if (CS%physical_OBL_scheme) then ! Check that Kdml is not set when using bulk mixed layer - call get_param(param_file, mdl, "KDML", CS%Kd_tot_ml, default=-1., do_not_log=.true.) + call get_param(param_file, mdl, "KDML", CS%Kd_tot_ml, & + units="m2 s-1", default=-1., scale=US%m2_s_to_Z2_T, do_not_log=.true.) if (CS%Kd_tot_ml>0.) call MOM_error(FATAL, & "bkgnd_mixing_init: KDML is a depricated parameter that should not be used.") - call get_param(param_file, mdl, "KD_ML_TOT", CS%Kd_tot_ml, default=-1., do_not_log=.true.) + call get_param(param_file, mdl, "KD_ML_TOT", CS%Kd_tot_ml, & + units="m2 s-1", default=-1., scale=US%m2_s_to_Z2_T, do_not_log=.true.) if (CS%Kd_tot_ml>0.) call MOM_error(FATAL, & "bkgnd_mixing_init: KD_ML_TOT cannot be set when using a physically based ocean "//& "boundary layer mixing parameterization.") @@ -338,8 +340,8 @@ subroutine calculate_bkgnd_mixing(h, tv, N2_lay, Kd_lay, Kd_int, Kv_bkgnd, j, G, real :: I_2Omega !< 1/(2 Omega) [T ~> s] real :: N_2Omega ! The ratio of the stratification to the Earth's rotation rate [nondim] real :: N02_N2 ! The ratio a reference stratification to the actual stratification [nondim] - real :: I_x30 !< 2/acos(2) = 1/(sin(30 deg) * acosh(1/sin(30 deg))) - real :: deg_to_rad !< factor converting degrees to radians, pi/180. + real :: I_x30 !< 2/acos(2) = 1/(sin(30 deg) * acosh(1/sin(30 deg))) [nondim] + real :: deg_to_rad !< factor converting degrees to radians [radians degree-1], pi/180. real :: abs_sinlat !< absolute value of sine of latitude [nondim] real :: min_sinlat ! The minimum value of the sine of latitude [nondim] real :: bckgrnd_vdc_psin !< PSI diffusivity in northern hemisphere [Z2 T-1 ~> m2 s-1] @@ -455,7 +457,7 @@ subroutine calculate_bkgnd_mixing(h, tv, N2_lay, Kd_lay, Kd_int, Kv_bkgnd, j, G, enddo endif - ! Now set background diffusivies based on these surface values, possibly with vertical structure. + ! Now set background diffusivities based on these surface values, possibly with vertical structure. if ((.not.CS%physical_OBL_scheme) .and. (CS%Kd /= CS%Kd_tot_ml)) then ! This is a crude way to put in a diffusive boundary layer without an explicit boundary ! layer turbulence scheme. It should not be used for any realistic ocean models. @@ -527,7 +529,7 @@ subroutine check_bkgnd_scheme(CS, str) end subroutine -!> Clear pointers and dealocate memory +!> Clear pointers and deallocate memory subroutine bkgnd_mixing_end(CS) type(bkgnd_mixing_cs), pointer :: CS !< Control structure for this module that !! will be deallocated in this subroutine diff --git a/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 b/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 index 49d62bbde4..2f8e03480e 100644 --- a/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 +++ b/src/parameterizations/vertical/MOM_bulk_mixed_layer.F90 @@ -35,7 +35,7 @@ module MOM_bulk_mixed_layer integer :: nkbl !< The number of buffer layers. integer :: nsw !< The number of bands of penetrating shortwave radiation. real :: mstar !< The ratio of the friction velocity cubed to the - !! TKE input to the mixed layer, nondimensional. + !! TKE input to the mixed layer [nondim]. real :: nstar !< The fraction of the TKE input to the mixed layer !! available to drive entrainment [nondim]. real :: nstar2 !< The fraction of potential energy released by @@ -43,7 +43,7 @@ module MOM_bulk_mixed_layer logical :: absorb_all_SW !< If true, all shortwave radiation is absorbed by the !! ocean, instead of passing through to the bottom mud. real :: TKE_decay !< The ratio of the natural Ekman depth to the TKE - !! decay scale, nondimensional. + !! decay scale [nondim]. real :: bulk_Ri_ML !< The efficiency with which mean kinetic energy !! released by mechanically forced entrainment of !! the mixed layer is converted to TKE [nondim]. @@ -84,9 +84,9 @@ module MOM_bulk_mixed_layer integer :: ML_presort_nz_conv_adj !< If ML_resort is true, do convective !! adjustment on this many layers (starting from the !! top) before sorting the remaining layers. - real :: omega_frac !< When setting the decay scale for turbulence, use - !! this fraction of the absolute rotation rate blended - !! with the local value of f, as sqrt((1-of)*f^2 + of*4*omega^2). + real :: omega_frac !< When setting the decay scale for turbulence, use this fraction + !! of the absolute rotation rate blended with the local value of f, + !! as sqrt((1-of)*f^2 + of*4*omega^2) [nondim]. logical :: correct_absorption !< If true, the depth at which penetrating !! shortwave radiation is absorbed is corrected by !! moving some of the heating upward in the water @@ -105,9 +105,8 @@ module MOM_bulk_mixed_layer !! points of the surface region (mixed & buffer !! layer) thickness [nondim]. 0.5 by default. real :: lim_det_dH_bathy !< The fraction of the total depth by which the - !! thickness of the surface region (mixed & buffer - !! layer) is allowed to change between grid points. - !! Nondimensional, 0.2 by default. + !! thickness of the surface region (mixed & buffer layers) is allowed + !! to change between grid points [nondim]. 0.2 by default. logical :: use_river_heat_content !< If true, use the fluxes%runoff_Hflx field !! to set the heat carried by runoff, instead of !! using SST for temperature of liq_runoff @@ -118,21 +117,21 @@ module MOM_bulk_mixed_layer type(diag_ctrl), pointer :: diag => NULL() !< A structure that is used to regulate the !! timing of diagnostic output. real :: Allowed_T_chg !< The amount by which temperature is allowed - !! to exceed previous values during detrainment, K. + !! to exceed previous values during detrainment [C ~> degC] real :: Allowed_S_chg !< The amount by which salinity is allowed !! to exceed previous values during detrainment [S ~> ppt] ! These are terms in the mixed layer TKE budget, all in [Z L2 T-3 ~> m3 s-3] except as noted. real, allocatable, dimension(:,:) :: & ML_depth, & !< The mixed layer depth [H ~> m or kg m-2]. - diag_TKE_wind, & !< The wind source of TKE. - diag_TKE_RiBulk, & !< The resolved KE source of TKE. - diag_TKE_conv, & !< The convective source of TKE. - diag_TKE_pen_SW, & !< The TKE sink required to mix penetrating shortwave heating. - diag_TKE_mech_decay, & !< The decay of mechanical TKE. - diag_TKE_conv_decay, & !< The decay of convective TKE. - diag_TKE_mixing, & !< The work done by TKE to deepen the mixed layer. - diag_TKE_conv_s2, & !< The convective source of TKE due to to mixing in sigma2. + diag_TKE_wind, & !< The wind source of TKE [Z L2 T-3 ~> m3 s-3]. + diag_TKE_RiBulk, & !< The resolved KE source of TKE [Z L2 T-3 ~> m3 s-3]. + diag_TKE_conv, & !< The convective source of TKE [Z L2 T-3 ~> m3 s-3]. + diag_TKE_pen_SW, & !< The TKE sink required to mix penetrating shortwave heating [Z L2 T-3 ~> m3 s-3]. + diag_TKE_mech_decay, & !< The decay of mechanical TKE [Z L2 T-3 ~> m3 s-3]. + diag_TKE_conv_decay, & !< The decay of convective TKE [Z L2 T-3 ~> m3 s-3]. + diag_TKE_mixing, & !< The work done by TKE to deepen the mixed layer [Z L2 T-3 ~> m3 s-3]. + diag_TKE_conv_s2, & !< The convective source of TKE due to to mixing in sigma2 [Z L2 T-3 ~> m3 s-3]. diag_PE_detrain, & !< The spurious source of potential energy due to mixed layer !! detrainment [R Z L2 T-3 ~> W m-2]. diag_PE_detrain2 !< The spurious source of potential energy due to mixed layer only @@ -171,10 +170,10 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt, ea, eb, G, GV, US, C !! [L T-1 ~> m s-1]. type(thermo_var_ptrs), intent(inout) :: tv !< A structure containing pointers to any !! available thermodynamic fields. Absent - !! fields have NULL ptrs. + !! fields have NULL pointers. type(forcing), intent(inout) :: fluxes !< A structure containing pointers to any !! possible forcing fields. Unused fields - !! have NULL ptrs. + !! have NULL pointers. real, intent(in) :: dt !< Time increment [T ~> s]. real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), & intent(inout) :: ea !< The amount of fluid moved downward into a @@ -184,7 +183,7 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt, ea, eb, G, GV, US, C intent(inout) :: eb !< The amount of fluid moved upward into a !! layer; this should be increased due to !! mixed layer entrainment [H ~> m or kg m-2]. - type(bulkmixedlayer_CS), intent(inout) :: CS !< Bulk mixed layer control struct + type(bulkmixedlayer_CS), intent(inout) :: CS !< Bulk mixed layer control structure type(optics_type), pointer :: optics !< The structure that can be queried for the !! inverse of the vertical absorption decay !! scale for penetrating shortwave radiation. @@ -195,7 +194,7 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt, ea, eb, G, GV, US, C !! being applied separately. real, optional, intent(in) :: dt_diag !< The diagnostic time step, !! which may be less than dt if there are - !! two callse to mixedlayer [T ~> s]. + !! two calls to mixedlayer [T ~> s]. logical, optional, intent(in) :: last_call !< if true, this is the last call !! to mixedlayer in the current time step, so !! diagnostics will be written. The default is @@ -247,8 +246,8 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt, ea, eb, G, GV, US, C ! entrained [C H ~> degC m or degC kg m-2]. Stot, & ! The integrated salt of layers which are fully entrained ! [H S ~> m ppt or ppt kg m-2]. - uhtot, & ! The depth integrated zonal and meridional velocities in the - vhtot, & ! mixed layer [H L T-1 ~> m2 s-1 or kg m-1 s-1]. + uhtot, & ! The depth integrated zonal velocity in the mixed layer [H L T-1 ~> m2 s-1 or kg m-1 s-1] + vhtot, & ! The depth integrated meridional velocity in the mixed layer [H L T-1 ~> m2 s-1 or kg m-1 s-1] netMassInOut, & ! The net mass flux (if non-Boussinsq) or volume flux (if ! Boussinesq - i.e. the fresh water flux (P+R-E)) into the @@ -278,7 +277,7 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt, ea, eb, G, GV, US, C Pen_SW_bnd ! The penetrating fraction of the shortwave heating integrated ! over a time step in each band [C H ~> degC m or degC kg m-2]. real, dimension(max(CS%nsw,1),SZI_(G),SZK_(GV)) :: & - opacity_band ! The opacity in each band [H-1 ~> m-1 or m2 kg-1]. The indicies are band, i, k. + opacity_band ! The opacity in each band [H-1 ~> m-1 or m2 kg-1]. The indices are band, i, k. real :: cMKE(2,SZI_(G)) ! Coefficients of HpE and HpE^2 used in calculating the ! denominator of MKE_rate; the two elements have differing @@ -318,7 +317,7 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt, ea, eb, G, GV, US, C real :: absf_x_H ! The absolute value of f times the mixed layer thickness [Z T-1 ~> m s-1]. real :: kU_star ! Ustar times the Von Karman constant [Z T-1 ~> m s-1]. - real :: dt__diag ! A recaled copy of dt_diag (if present) or dt [T ~> s]. + real :: dt__diag ! A rescaled copy of dt_diag (if present) or dt [T ~> s]. logical :: write_diags ! If true, write out diagnostics with this step. logical :: reset_diags ! If true, zero out the accumulated diagnostics. integer, dimension(2) :: EOSdom ! The i-computational domain for the equation of state @@ -585,9 +584,9 @@ subroutine bulkmixedlayer(h_3d, u_3d, v_3d, tv, fluxes, dt, ea, eb, G, GV, US, C enddo ; endif endif -! Move water left in the former mixed layer into the buffer layer and -! from the buffer layer into the interior. These steps might best be -! treated in conjuction. + ! Move water left in the former mixed layer into the buffer layer and + ! from the buffer layer into the interior. These steps might best be + ! treated in conjunction. if (CS%nkbl == 1) then call mixedlayer_detrain_1(h(:,0:), T(:,0:), S(:,0:), R0(:,0:), Rcv(:,0:), & GV%Rlay(:), dt, dt__diag, d_ea, d_eb, j, G, GV, US, CS, & @@ -777,7 +776,7 @@ subroutine convective_adjustment(h, u, v, R0, Rcv, T, S, eps, d_eb, & !! [Z L2 T-2 ~> m3 s-2]. integer, intent(in) :: j !< The j-index to work on. type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type - type(bulkmixedlayer_CS), intent(in) :: CS !< Bulk mixed layer control struct + type(bulkmixedlayer_CS), intent(in) :: CS !< Bulk mixed layer control structure integer, optional, intent(in) :: nz_conv !< If present, the number of layers !! over which to do convective adjustment !! (perhaps CS%nkml). @@ -952,13 +951,13 @@ subroutine mixedlayer_convection(h, d_eb, htot, Ttot, Stot, uhtot, vhtot, & integer, dimension(SZI_(G),SZK_(GV)), & intent(in) :: ksort !< The density-sorted k-indices. type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type - type(bulkmixedlayer_CS), intent(in) :: CS !< Bulk mixed layer control struct + type(bulkmixedlayer_CS), intent(in) :: CS !< Bulk mixed layer control structure type(thermo_var_ptrs), intent(inout) :: tv !< A structure containing pointers to any !! available thermodynamic fields. Absent - !! fields have NULL ptrs. + !! fields have NULL pointers. type(forcing), intent(inout) :: fluxes !< A structure containing pointers to any !! possible forcing fields. Unused fields - !! have NULL ptrs. + !! have NULL pointers. real, intent(in) :: dt !< Time increment [T ~> s]. logical, intent(in) :: aggregate_FW_forcing !< If true, the net incoming and !! outgoing surface freshwater fluxes are @@ -1261,7 +1260,7 @@ subroutine find_starting_TKE(htot, h_CA, fluxes, Conv_En, cTKE, dKE_FC, dKE_CA, !! adjustment [H ~> m or kg m-2]. type(forcing), intent(in) :: fluxes !< A structure containing pointers to any !! possible forcing fields. Unused fields - !! have NULL ptrs. + !! have NULL pointers. real, dimension(SZI_(G)), intent(inout) :: Conv_En !< The buoyant turbulent kinetic energy source !! due to free convection [Z L2 T-2 ~> m3 s-2]. real, dimension(SZI_(G)), intent(in) :: dKE_FC !< The vertically integrated change in @@ -1290,8 +1289,8 @@ subroutine find_starting_TKE(htot, h_CA, fluxes, Conv_En, cTKE, dKE_FC, dKE_CA, !! time interval [T-1 ~> s-1]. integer, intent(in) :: j !< The j-index to work on. integer, dimension(SZI_(G),SZK_(GV)), & - intent(in) :: ksort !< The density-sorted k-indicies. - type(bulkmixedlayer_CS), intent(inout) :: CS !< Bulk mixed layer control struct + intent(in) :: ksort !< The density-sorted k-indices. + type(bulkmixedlayer_CS), intent(inout) :: CS !< Bulk mixed layer control structure ! This subroutine determines the TKE available at the depth of free ! convection to drive mechanical entrainment. @@ -1500,14 +1499,14 @@ subroutine mechanical_entrainment(h, d_eb, htot, Ttot, Stot, uhtot, vhtot, & real, dimension(SZI_(G)), intent(inout) :: Idecay_len_TKE !< The vertical TKE decay rate [H-1 ~> m-1 or m2 kg-1]. integer, intent(in) :: j !< The j-index to work on. integer, dimension(SZI_(G),SZK_(GV)), & - intent(in) :: ksort !< The density-sorted k-indicies. - type(bulkmixedlayer_CS), intent(inout) :: CS !< Bulk mixed layer control struct + intent(in) :: ksort !< The density-sorted k-indices. + type(bulkmixedlayer_CS), intent(inout) :: CS !< Bulk mixed layer control structure ! This subroutine calculates mechanically driven entrainment. ! Local variables real :: SW_trans ! The fraction of shortwave radiation that is not - ! absorbed in a layer, nondimensional. + ! absorbed in a layer [nondim]. real :: Pen_absorbed ! The amount of penetrative shortwave radiation ! that is absorbed in a layer [C H ~> degC m or degC kg m-2]. real :: h_avail ! The thickness in a layer available for entrainment [H ~> m or kg m-2]. @@ -1517,7 +1516,7 @@ subroutine mechanical_entrainment(h, d_eb, htot, Ttot, Stot, uhtot, vhtot, & ! h_ent between iterations [H ~> m or kg m-2]. real :: MKE_rate ! The fraction of the energy in resolved shears ! within the mixed layer that will be eliminated - ! within a timestep, nondim, 0 to 1. + ! within a timestep [nondim], 0 to 1. real :: HpE ! The current thickness plus entrainment [H ~> m or kg m-2]. real :: g_H_2Rho0 ! Half the gravitational acceleration times the ! conversion from H to m divided by the mean density, @@ -1541,17 +1540,17 @@ subroutine mechanical_entrainment(h, d_eb, htot, Ttot, Stot, uhtot, vhtot, & 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 :: dEF4_dh ! The partial derivative of EF4 with h [H-2 ~> m-2 or m4 kg-2]. - real :: Pen_En1 ! A nondimensional temporary variable. - real :: kh, exp_kh ! Nondimensional temporary variables related to the - real :: f1_kh ! fractional decay of TKE across a layer. - real :: x1, e_x1 ! Nondimensional temporary variables related to - real :: f1_x1, f2_x1 ! the relative decay of TKE and SW radiation across - real :: f3_x1 ! a layer, and exponential-related functions of x1. + real :: Pen_En1 ! A nondimensional temporary variable [nondim]. + real :: kh, exp_kh, f1_kh ! Nondimensional temporary variables related to the + ! fractional decay of TKE across a layer [nondim]. + real :: x1, e_x1 ! Nondimensional temporary variables related to the relative decay + ! of TKE and SW radiation across a layer [nondim] + real :: f1_x1, f2_x1, f3_x1 ! Exponential-related functions of x1 [nondim]. real :: E_HxHpE ! Entrainment divided by the product of the new and old ! thicknesses [H-1 ~> m-1 or m2 kg-1]. real :: Hmix_min ! The minimum mixed layer depth [H ~> m or kg m-2]. - real :: opacity - real :: C1_3, C1_6, C1_24 ! 1/3, 1/6, and 1/24. + real :: opacity ! The opacity of a layer in a band of shortwave radiation [H-1 ~> m-1 or m2 kg-1] + real :: C1_3, C1_6, C1_24 ! 1/3, 1/6, and 1/24. [nondim] integer :: is, ie, nz, i, k, ks, itt, n C1_3 = 1.0/3.0 ; C1_6 = 1.0/6.0 ; C1_24 = 1.0/24.0 @@ -1784,12 +1783,12 @@ subroutine sort_ML(h, R0, eps, G, GV, CS, ksort) !! the layers [R ~> kg m-3]. real, dimension(SZI_(G),SZK_(GV)), intent(in) :: eps !< The (small) thickness that must !! remain in each layer [H ~> m or kg m-2]. - type(bulkmixedlayer_CS), intent(in) :: CS !< Bulk mixed layer control struct + type(bulkmixedlayer_CS), intent(in) :: CS !< Bulk mixed layer control structure integer, dimension(SZI_(G),SZK_(GV)), intent(out) :: ksort !< The k-index to use in the sort. ! Local variables - real :: R0sort(SZI_(G),SZK_(GV)) - integer :: nsort(SZI_(G)) + real :: R0sort(SZI_(G),SZK_(GV)) ! The sorted potential density [R ~> kg m-3] + integer :: nsort(SZI_(G)) ! The number of layers left to sort logical :: done_sorting(SZI_(G)) integer :: i, k, ks, is, ie, nz, nkmb @@ -1852,14 +1851,14 @@ subroutine resort_ML(h, T, S, R0, Rcv, RcvTgt, eps, d_ea, d_eb, ksort, G, GV, CS !! layer in the entrainment from !! below [H ~> m or kg m-2]. Positive values go !! with mass gain by a layer. - integer, dimension(SZI_(G),SZK_(GV)), intent(in) :: ksort !< The density-sorted k-indicies. - type(bulkmixedlayer_CS), intent(in) :: CS !< Bulk mixed layer control struct + integer, dimension(SZI_(G),SZK_(GV)), intent(in) :: ksort !< The density-sorted k-indices. + type(bulkmixedlayer_CS), intent(in) :: CS !< Bulk mixed layer control structure real, dimension(SZI_(G)), intent(in) :: dR0_dT !< The partial derivative of !! potential density referenced !! to the surface with potential !! temperature [R C-1 ~> kg m-3 degC-1]. real, dimension(SZI_(G)), intent(in) :: dR0_dS !< The partial derivative of - !! cpotential density referenced + !! potential density referenced !! to the surface with salinity, !! [R S-1 ~> kg m-3 ppt-1]. real, dimension(SZI_(G)), intent(in) :: dRcv_dT !< The partial derivative of @@ -1880,21 +1879,38 @@ subroutine resort_ML(h, T, S, R0, Rcv, RcvTgt, eps, d_ea, d_eb, ksort, G, GV, CS ! and the coordinate density (sigma-2)) between the newly forming mixed layer ! and a residual buffer- or mixed layer, and the number of massive layers above ! the deepest massive buffer or mixed layer is greater than nkbl, then split -! those buffer layers into peices that match the target density of the two +! those buffer layers into pieces that match the target density of the two ! nearest interior layers. ! Otherwise, if there are more than nkbl+1 remaining massive layers ! Local variables - real :: h_move, h_tgt_old, I_hnew - real :: dT_dS_wt2, dT_dR, dS_dR, I_denom - real :: Rcv_int - real :: T_up, S_up, R0_up, I_hup, h_to_up - real :: T_dn, S_dn, R0_dn, I_hdn, h_to_dn - real :: wt_dn - real :: dR1, dR2 - real :: dPE, hmin, min_dPE, min_hmin - real, dimension(SZK_(GV)) :: & - h_tmp, R0_tmp, T_tmp, S_tmp, Rcv_tmp + real :: h_move ! The thickness of water being moved between layers [H ~> m or kg m-2] + real :: h_tgt_old ! The previous thickness of the recipient layer [H ~> m or kg m-2] + real :: I_hnew ! The inverse of a new layer thickness [H-1 ~> m-1 or m3 kg-1] + real :: dT_dS_wt2 ! The square of the relative weighting of temperature and salinity changes + ! when extraploating to match a target density [C2 S-2 ~> degC2 ppt-2] + real :: dT_dR ! The ratio of temperature changes to density changes when + ! extrapolating [C R-1 ~> degC m3 kg-1] + real :: dS_dR ! The ratio of salinity changes to density changes when + ! extrapolating [S R-1 ~> ppt m3 kg-1] + real :: I_denom ! A work variable with units of [S2 R-2 ~> ppt2 m6 kg-2]. + real :: Rcv_int ! The target coordinate density of an interior layer [R ~> kg m-3] + real :: T_up, T_dn ! Temperatures projected to match the target densities of two layers [C ~> degC] + real :: S_up, S_dn ! Salinities projected to match the target densities of two layers [S ~> ppt] + real :: R0_up, R0_dn ! Potential densities projected to match the target coordinate + ! densities of two layers [R ~> kg m-3] + real :: I_hup, I_hdn ! Inverse of the new thicknesses of the two layers [H-1 ~> m-1 or m2 kg-1] + real :: h_to_up, h_to_dn ! Thickness transferred to two layers [H ~> m or kg m-2] + real :: wt_dn ! Fraction of the thickness transferred to the deeper layer [nondim] + real :: dR1, dR2 ! Density difference with the target densities of two layers [R ~> kg m-3] + real :: dPE, min_dPE ! Values proportional to the potential energy change due to the merging + ! of a pair of layers [R H2 ~> kg m-1 or kg3 m-6] + real :: hmin, min_hmin ! The thickness of the thinnest layer [H ~> m or kg m-2] + real :: h_tmp(SZK_(GV)) ! A copy of the original layer thicknesses [H ~> m or kg m-2] + real :: R0_tmp(SZK_(GV)) ! A copy of the original layer potential densities [R ~> kg m-3] + real :: T_tmp(SZK_(GV)) ! A copy of the original layer temperatures [C ~> degC] + real :: S_tmp(SZK_(GV)) ! A copy of the original layer salinities [S ~> ppt] + real :: Rcv_tmp(SZK_(GV)) ! A copy of the original layer coordinate densities [R ~> kg m-3] integer :: ks_min logical :: sorted, leave_in_layer integer :: ks_deep(SZI_(G)), k_count(SZI_(G)), ks2_reverse(SZI_(G), SZK_(GV)) @@ -2168,13 +2184,13 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt, dt_diag, d_ea, j, !! goes with layer thickness increases. integer, intent(in) :: j !< The meridional row to work on. type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type - type(bulkmixedlayer_CS), intent(inout) :: CS !< Bulk mixed layer control struct + type(bulkmixedlayer_CS), intent(inout) :: CS !< Bulk mixed layer control structure real, dimension(SZI_(G)), intent(in) :: dR0_dT !< The partial derivative of !! potential density referenced to the !! surface with potential temperature, !! [R C-1 ~> kg m-3 degC-1]. real, dimension(SZI_(G)), intent(in) :: dR0_dS !< The partial derivative of - !! cpotential density referenced to the + !! potential density referenced to the !! surface with salinity !! [R S-1 ~> kg m-3 ppt-1]. real, dimension(SZI_(G)), intent(in) :: dRcv_dT !< The partial derivative of @@ -2224,10 +2240,11 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt, dt_diag, d_ea, j, real :: stays_min_merge ! The minimum allowed value of stays_merge [H ~> m or kg m-2]. real :: dR0_2dz, dRcv_2dz ! Half the vertical gradients of R0 and Rcv [R H-1 ~> kg m-4 or m-1] -! real :: dT_2dz, dS_2dz ! Half the vertical gradients of T and S, in degC H-1, and ppt H-1. +! real :: dT_2dz ! Half the vertical gradient of T [C H-1 ~> degC m-1 or degC m2 kg-1] +! real :: dS_2dz ! Half the vertical gradient of S [S H-1 ~> ppt m-1 or ppt m2 kg-1] real :: scale_slope ! A nondimensional number < 1 used to scale down ! the slope within the upper buffer layer when - ! water MUST be detrained to the lower layer. + ! water MUST be detrained to the lower layer [nondim]. real :: dPE_extrap ! The potential energy change due to dispersive ! advection or mixing layers, divided by @@ -2264,9 +2281,9 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt, dt_diag, d_ea, j, real :: dPE_ratio ! Multiplier of dPE_det at which merging is ! permitted - here (detrainment_per_day/dt)*30 - ! days? + ! days? [nondim] real :: num_events ! The number of detrainment events over which - ! to prefer merging the buffer layers. + ! to prefer merging the buffer layers [nondim]. real :: dPE_time_ratio ! Larger of 1 and the detrainment timescale over dt [nondim]. real :: dT_dS_gauge, dS_dT_gauge ! The relative scales of temperature and ! salinity changes in defining spiciness, in @@ -2287,14 +2304,16 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt, dt_diag, d_ea, j, real :: s1, s2, bh0 ! Work variables [H ~> m or kg m-2]. real :: s3sq ! A work variable [H2 ~> m2 or kg2 m-4]. real :: I_ya, b1 ! Nondimensional work variables [nondim] - real :: Ih, Ihdet, Ih1f, Ih2f ! Assorted inverse thickness work variables, - real :: Ihk0, Ihk1, Ih12 ! all in [H-1 ~> m-1 or m2 kg-1]. - real :: dR1, dR2, dR2b, dRk1 ! Assorted density difference work variables, - real :: dR0, dR21, dRcv ! all in [R ~> kg m-3]. + real :: Ih, Ihdet, Ih1f, Ih2f ! Assorted inverse thickness work variables [H-1 ~> m-1 or m2 kg-1] + real :: Ihk0, Ihk1, Ih12 ! Assorted inverse thickness work variables [H-1 ~> m-1 or m2 kg-1] + real :: dR1, dR2, dR2b, dRk1 ! Assorted density difference work variables [R ~> kg m-3] + real :: dR0, dR21, dRcv ! Assorted density difference work variables [R ~> kg m-3] real :: dRcv_stays, dRcv_det, dRcv_lim ! Assorted densities [R ~> kg m-3] - real :: Angstrom ! The minumum layer thickness [H ~> m or kg m-2]. + real :: Angstrom ! The minimum layer thickness [H ~> m or kg m-2]. - real :: h2_to_k1_lim, T_new, S_new, T_max, T_min, S_max, S_min + real :: h2_to_k1_lim ! A limit on the thickness that can be detrained to layer k1 [H ~> m or kg m-2] + real :: T_new, T_max, T_min ! Temperature of the detrained water and limits on it [C ~> degC] + real :: S_new, S_max, S_min ! Salinity of the detrained water and limits on it [S ~> ppt] integer :: i, k, k0, k1, is, ie, nz, kb1, kb2, nkmb is = G%isc ; ie = G%iec ; nz = GV%ke @@ -2352,7 +2371,7 @@ subroutine mixedlayer_detrain_2(h, T, S, R0, Rcv, RcvTgt, dt, dt_diag, d_ea, j, ! (3) The lower buffer layer density extrapolated to its base with a ! linear fit between the two layers must exceed the density of the ! next denser interior layer. - ! (4) The average extroplated coordinate density that is moved into the + ! (4) The average extrapolated coordinate density that is moved into the ! isopycnal interior matches the target value for that layer. ! (5) The potential energy change is calculated and might be used later ! to allow the upper buffer layer to mix more into the lower buffer @@ -3062,7 +3081,7 @@ subroutine mixedlayer_detrain_1(h, T, S, R0, Rcv, RcvTgt, dt, dt_diag, d_ea, d_e !! a layer. integer, intent(in) :: j !< The meridional row to work on. type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type - type(bulkmixedlayer_CS), intent(inout) :: CS !< Bulk mixed layer control struct + type(bulkmixedlayer_CS), intent(inout) :: CS !< Bulk mixed layer control structure real, dimension(SZI_(G)), intent(in) :: dRcv_dT !< The partial derivative of !! coordinate defining potential density !! with potential temperature @@ -3081,9 +3100,17 @@ subroutine mixedlayer_detrain_1(h, T, S, R0, Rcv, RcvTgt, dt, dt_diag, d_ea, d_e real :: max_det_rem(SZI_(G)) ! Remaining permitted detrainment [H ~> m or kg m-2]. real :: detrain(SZI_(G)) ! The thickness of fluid to detrain ! from the mixed layer [H ~> m or kg m-2]. - real :: dT_dR, dS_dR, dRml, dR0_dRcv, dT_dS_wt2 + real :: dT_dS_wt2 ! The square of the relative weighting of temperature and salinity changes + ! when extraploating to match a target density [C2 S-2 ~> degC2 ppt-2] + real :: dT_dR ! The ratio of temperature changes to density changes when + ! extrapolating [C R-1 ~> degC m3 kg-1] + real :: dS_dR ! The ratio of salinity changes to density changes when + ! extrapolating [S R-1 ~> ppt m3 kg-1] + real :: dRml ! The density range within the extent of the mixed layers [R ~> kg m-3] + real :: dR0_dRcv ! The relative changes in the potential density and the coordinate density [nondim] real :: I_denom ! A work variable [S2 R-2 ~> ppt2 m6 kg-2]. - real :: Sdown, Tdown ! A salinity [S ~> ppt] and a temperature [C ~> degC] + real :: Sdown ! The salinity of the detrained water [S ~> ppt] + real :: Tdown ! The temperature of the detrained water [C ~> degC] real :: dt_Time ! The timestep divided by the detrainment timescale [nondim]. real :: g_H2_2Rho0dt ! Half the gravitational acceleration times the square of the ! conversion from H to m divided by the mean density times the time @@ -3091,11 +3118,10 @@ subroutine mixedlayer_detrain_1(h, T, S, R0, Rcv, RcvTgt, dt, dt_diag, d_ea, d_e real :: g_H2_2dt ! Half the gravitational acceleration times the square of the ! conversion from H to Z divided by the diagnostic time step ! [L2 Z H-2 T-3 ~> m s-3 or m7 kg-2 s-3]. - + real :: x1 ! A temporary work variable [various] logical :: splittable_BL(SZI_(G)), orthogonal_extrap - real :: x1 - integer :: i, is, ie, k, k1, nkmb, nz + is = G%isc ; ie = G%iec ; nz = GV%ke nkmb = CS%nkml+CS%nkbl if (CS%nkbl /= 1) call MOM_error(FATAL,"MOM_mixed_layer: "// & @@ -3329,13 +3355,15 @@ subroutine bulkmixedlayer_init(Time, G, GV, US, param_file, diag, CS) !! parameters. type(diag_ctrl), target, intent(inout) :: diag !< A structure that is used to regulate diagnostic !! output. - type(bulkmixedlayer_CS), intent(inout) :: CS !< Bulk mixed layer control struct + type(bulkmixedlayer_CS), intent(inout) :: CS !< Bulk mixed layer control structure ! This include declares and sets the variable "version". #include "version_variable.h" character(len=40) :: mdl = "MOM_mixed_layer" ! This module's name. real :: BL_detrain_time_dflt ! The default value for BUFFER_LAY_DETRAIN_TIME [s] - real :: omega_frac_dflt, ustar_min_dflt, Hmix_min_m + real :: omega_frac_dflt ! The default value for ML_OMEGA_FRAC [nondim] + real :: ustar_min_dflt ! The default value for BML_USTAR_MIN [m s-1] + real :: Hmix_min_m ! The unscaled value of HMIX_MIN [m] integer :: isd, ied, jsd, jed logical :: use_temperature, use_omega isd = G%isd ; ied = G%ied ; jsd = G%jsd ; jed = G%jed @@ -3600,7 +3628,7 @@ function EF4(Ht, En, I_L, dR_de) real :: EF4 !< The integral [H-1 ~> m-1 or m2 kg-1]. ! Local variables - real :: exp_LHpE ! A nondimensional exponential decay. + real :: exp_LHpE ! A nondimensional exponential decay [nondim]. real :: I_HpE ! An inverse thickness plus entrainment [H-1 ~> m-1 or m2 kg-1]. real :: Res ! The result of the integral above [H-1 ~> m-1 or m2 kg-1]. diff --git a/src/parameterizations/vertical/MOM_diabatic_aux.F90 b/src/parameterizations/vertical/MOM_diabatic_aux.F90 index a3450bd6e4..d51f796df1 100644 --- a/src/parameterizations/vertical/MOM_diabatic_aux.F90 +++ b/src/parameterizations/vertical/MOM_diabatic_aux.F90 @@ -242,20 +242,19 @@ subroutine differential_diffuse_T_S(h, T, S, Kd_T, Kd_S, dt, G, GV) ! local variables real, dimension(SZI_(G)) :: & - b1_T, b1_S, & ! Variables used by the tridiagonal solvers of T & S [H ~> m or kg m-2]. - d1_T, d1_S ! Variables used by the tridiagonal solvers [nondim]. + b1_T, b1_S, & ! Variables used by the tridiagonal solvers of T & S [H ~> m or kg m-2]. + d1_T, d1_S ! Variables used by the tridiagonal solvers [nondim]. real, dimension(SZI_(G),SZK_(GV)) :: & - c1_T, c1_S ! Variables used by the tridiagonal solvers [H ~> m or kg m-2]. + c1_T, c1_S ! Variables used by the tridiagonal solvers [H ~> m or kg m-2]. real, dimension(SZI_(G),SZK_(GV)+1) :: & - mix_T, mix_S ! Mixing distances in both directions across each interface [H ~> m or kg m-2]. - real :: h_tr ! h_tr is h at tracer points with a tiny thickness - ! added to ensure positive definiteness [H ~> m or kg m-2]. - 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 :: I_h_int ! The inverse of the thickness associated with an - ! interface [H-1 ~> m-1 or m2 kg-1]. - real :: b_denom_T ! The first term in the denominators for the expressions - real :: b_denom_S ! for b1_T and b1_S, both [H ~> m or kg m-2]. + mix_T, mix_S ! Mixing distances in both directions across each interface [H ~> m or kg m-2]. + real :: h_tr ! h_tr is h at tracer points with a tiny thickness + ! added to ensure positive definiteness [H ~> m or kg m-2]. + 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 :: I_h_int ! The inverse of the thickness associated with an interface [H-1 ~> m-1 or m2 kg-1]. + real :: b_denom_T ! The first term in the denominator for the expression for b1_T [H ~> m or kg m-2]. + real :: b_denom_S ! The first term in the denominator for the expression for b1_S [H ~> m or kg m-2]. integer :: i, j, k, is, ie, js, je, nz is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke @@ -497,17 +496,18 @@ subroutine find_uv_at_h(u, v, h, u_h, v_h, G, GV, US, ea, eb, zero_mix) !! v_h as though ea and eb were being supplied with !! uniformly zero values. - ! local variables + ! Local variables real :: b_denom_1 ! The first term in the denominator of b1 [H ~> m or kg m-2]. 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 :: b1(SZI_(G)) ! A thickness used in the tridiagonal solver [H ~> m or kg m-2] real :: c1(SZI_(G),SZK_(GV)) ! A variable used in the tridiagonal solver [nondim] real :: d1(SZI_(G)) ! The complement of c1 [nondim] - real :: a_n(SZI_(G)), a_s(SZI_(G)) ! Fractional weights of the neighboring - real :: a_e(SZI_(G)), a_w(SZI_(G)) ! velocity points, ~1/2 in the open - ! ocean, nondimensional. - real :: sum_area, Idenom + ! Fractional weights of the neighboring velocity points, ~1/2 in the open ocean. + real :: a_n(SZI_(G)), a_s(SZI_(G)) ! Fractional weights of the neighboring velocity points [nondim] + real :: a_e(SZI_(G)), a_w(SZI_(G)) ! Fractional weights of the neighboring velocity points [nondim] + real :: sum_area ! A sum of adjacent areas [L2 ~> m2] + real :: Idenom ! The inverse of the denomninator in a weighted average [L-2 ~> m-2] logical :: mix_vertically, zero_mixing integer :: i, j, k, is, ie, js, je, nz is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec ; nz = GV%ke @@ -527,6 +527,12 @@ subroutine find_uv_at_h(u, v, h, u_h, v_h, G, GV, US, ea, eb, zero_mix) do i=is,ie sum_area = G%areaCu(I-1,j) + G%areaCu(I,j) if (sum_area > 0.0) then + ! If this were a simple area weighted average, this would just be I_denom = 1.0 / sum_area. + ! The other factor of sqrt(0.5*sum_area*G%IareaT(i,j)) is 1 for open ocean points on a + ! Cartesian grid. This construct predates the initial commit of the MOM6 code, and was + ! present in the GOLD code before February, 2010. I do not recall why this was added, and + ! the GOLD CVS server that contained the relevant history and logs appears to have been + ! decommissioned. Idenom = sqrt(0.5*G%IareaT(i,j) / sum_area) a_w(i) = G%areaCu(I-1,j) * Idenom a_e(i) = G%areaCu(I,j) * Idenom @@ -803,7 +809,7 @@ subroutine diagnoseMLDbyEnergy(id_MLD, h, tv, G, GV, US, Mixing_Energy, diagPtr) ! converges extremely quickly (usually 1 guess) since this equation turns out to be rather ! linear for PE change with increasing X. ! Input parameters: - integer, dimension(3), intent(in) :: id_MLD !< Energy output diag IDs + integer, dimension(3), intent(in) :: id_MLD !< Energy output diagnostic IDs type(ocean_grid_type), intent(in) :: G !< Grid type type(verticalGrid_type), intent(in) :: GV !< ocean vertical grid structure type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type @@ -1045,13 +1051,19 @@ subroutine applyBoundaryFluxesInOut(CS, G, GV, US, dt, fluxes, optics, nsw, h, t ! Local variables integer, parameter :: maxGroundings = 5 integer :: numberOfGroundings, iGround(maxGroundings), jGround(maxGroundings) - real :: H_limit_fluxes - real :: IforcingDepthScale + real :: H_limit_fluxes ! Surface fluxes are scaled down fluxes when the total depth of the ocean + ! drops below this value [H ~> m or kg m-2] + real :: IforcingDepthScale ! The inverse of the layer thickness below which mass losses are + ! shifted to the next deeper layer [H ~> m or kg m-2] real :: Idt ! The inverse of the timestep [T-1 ~> s-1] - real :: dThickness, dTemp, dSalt - real :: fractionOfForcing, hOld, Ithickness + real :: dThickness ! The change in layer thickness [H ~> m or kg m-2] + real :: dTemp ! The integrated change in layer temperature [C H ~> degC m or degC kg m-2] + real :: dSalt ! The integrated change in layer salinity [S H ~> ppt m or ppt kg m-2] + real :: fractionOfForcing ! THe fraction of the remaining forcing applied to a layer [nondim] + real :: hOld ! The original thickness of a layer [H ~> m or kg m-2] + real :: Ithickness ! The inverse of the new layer thickness [H-1 ~> m-1 or m2 kg-1] real :: RivermixConst ! A constant used in implementing river mixing [R Z2 T-1 ~> Pa s]. - real :: EnthalpyConst ! A constant used to control the enthalpy calculation + real :: EnthalpyConst ! A constant used to control the enthalpy calculation [nondim] ! By default EnthalpyConst = 1.0. If fluxes%heat_content_evap ! is associated enthalpy is provided via coupler and EnthalpyConst = 0.0. real, dimension(SZI_(G)) :: & @@ -1092,13 +1104,17 @@ subroutine applyBoundaryFluxesInOut(CS, G, GV, US, dt, fluxes, optics, nsw, h, t opacityBand ! The opacity (inverse of the exponential absorption length) of each frequency ! band of shortwave radiation in each layer [H-1 ~> m-1 or m2 kg-1] real, dimension(maxGroundings) :: hGrounding ! Thickness added by each grounding event [H ~> m or kg m-2] - real :: Temp_in, Salin_in + real :: Temp_in ! The initial temperature of a layer [C ~> degC] + real :: Salin_in ! The initial salinity of a layer [S ~> ppt] real :: g_Hconv2 ! A conversion factor for use in the TKE calculation ! in units of [Z3 R2 T-2 H-2 ~> kg2 m-5 s-2 or m s-2]. real :: GoRho ! g_Earth times a unit conversion factor divided by density ! [Z T-2 R-1 ~> m4 s-2 kg-1] - logical :: calculate_energetics - logical :: calculate_buoyancy + logical :: calculate_energetics ! If true, calculate the energy required to mix the newly added + ! water over the topmost grid cell, assuming that the fluxes of heat and salt + ! and rejected brine are initially applied in vanishingly thin layers at the + ! top of the layer before being mixed throughout the layer. + logical :: calculate_buoyancy ! If true, calculate the surface buoyancy flux. integer, dimension(2) :: EOSdom ! The i-computational domain for the equation of state integer :: i, j, is, ie, js, je, k, nz, nb character(len=45) :: mesg diff --git a/src/parameterizations/vertical/MOM_diabatic_driver.F90 b/src/parameterizations/vertical/MOM_diabatic_driver.F90 index ddafbc3274..7cfcbfab07 100644 --- a/src/parameterizations/vertical/MOM_diabatic_driver.F90 +++ b/src/parameterizations/vertical/MOM_diabatic_driver.F90 @@ -122,7 +122,7 @@ module MOM_diabatic_driver !! other diffusivities. Otherwise, the larger of kappa- !! shear and ePBL diffusivities are used. real :: ePBL_Prandtl !< The Prandtl number used by ePBL to convert vertical - !! diffusivities into viscosities. + !! diffusivities into viscosities [nondim]. integer :: nMode = 1 !< Number of baroclinic modes to consider real :: uniform_test_cg !< Uniform group velocity of internal tide !! for testing internal tides [L T-1 ~> m s-1] @@ -133,7 +133,7 @@ module MOM_diabatic_driver !! FW fluxes are applied separately or combined before !! being applied. real :: ML_mix_first !< The nondimensional fraction of the mixed layer - !! algorithm that is applied before diffusive mixing. + !! algorithm that is applied before diffusive mixing [nondim]. !! The default is 0, while 0.5 gives Strang splitting !! and 1 is a sensible value too. Note that if there !! are convective instabilities in the initial state, @@ -174,8 +174,8 @@ module MOM_diabatic_driver real :: MLD_EN_VALS(3) !< Energy values for energy mixed layer diagnostics [R Z L2 T-2 ~> J m-2] !>@{ Diagnostic IDs - integer :: id_cg1 = -1 ! diag handle for mode-1 speed - integer, allocatable, dimension(:) :: id_cn ! diag handle for all mode speeds + integer :: id_cg1 = -1 ! diagnostic handle for mode-1 speed + integer, allocatable, dimension(:) :: id_cn ! diagnostic handle for all mode speeds integer :: id_ea = -1, id_eb = -1 ! used by layer diabatic integer :: id_ea_t = -1, id_eb_t = -1, id_ea_s = -1, id_eb_s = -1 integer :: id_Kd_heat = -1, id_Kd_salt = -1, id_Kd_int = -1, id_Kd_ePBL = -1 @@ -231,14 +231,14 @@ module MOM_diabatic_driver type(KPP_CS), pointer :: KPP_CSp => NULL() !< Control structure for a child module type(diapyc_energy_req_CS), pointer :: diapyc_en_rec_CSp => NULL() !< Control structure for a child module type(oda_incupd_CS), pointer :: oda_incupd_CSp => NULL() !< Control structure for a child module - type(bulkmixedlayer_CS) :: bulkmixedlayer !< Bulk mixed layer control struct - type(CVMix_conv_CS) :: CVMix_conv !< CVMix convection control struct - type(energetic_PBL_CS) :: ePBL !< Energetic PBL control struct - type(entrain_diffusive_CS) :: entrain_diffusive !< Diffusive entrainment control struct - type(geothermal_CS) :: geothermal !< Geothermal control struct - type(int_tide_CS) :: int_tide !< Internal tide control struct - type(opacity_CS) :: opacity !< Opacity control struct - type(regularize_layers_CS) :: regularize_layers !< Regularize layer control struct + type(bulkmixedlayer_CS) :: bulkmixedlayer !< Bulk mixed layer control structure + type(CVMix_conv_CS) :: CVMix_conv !< CVMix convection control structure + type(energetic_PBL_CS) :: ePBL !< Energetic PBL control structure + type(entrain_diffusive_CS) :: entrain_diffusive !< Diffusive entrainment control structure + type(geothermal_CS) :: geothermal !< Geothermal control structure + type(int_tide_CS) :: int_tide !< Internal tide control structure + type(opacity_CS) :: opacity !< Opacity control structure + type(regularize_layers_CS) :: regularize_layers !< Regularize layer control structure type(group_pass_type) :: pass_hold_eb_ea !< For group halo pass type(group_pass_type) :: pass_Kv !< For group halo pass @@ -1659,9 +1659,10 @@ subroutine layered_diabatic(u, v, h, tv, Hml, fluxes, visc, ADp, CDp, dt, Time_e real, dimension(SZI_(G),SZJ_(G),SZK_(GV)), target :: & ! These are targets so that the space can be shared with eaml & ebml. - eatr, & ! The equivalent of ea and eb for tracers, which differ from ea and - ebtr ! eb in that they tend to homogenize tracers in massless layers - ! near the boundaries [H ~> m or kg m-2] (for Bous or non-Bouss) + eatr, & ! The equivalent of ea for tracers, which differs from ea in that it tends to + ! homogenize tracers in massless layers near the boundaries [H ~> m or kg m-2] + ebtr ! The equivalent of eb for tracers, which differs from eb in that it tends to + ! homogenize tracers in massless layers near the boundaries [H ~> m or kg m-2] real, dimension(SZI_(G),SZJ_(G),SZK_(GV)+1) :: & Kd_int, & ! diapycnal diffusivity of interfaces [Z2 T-1 ~> m2 s-1] @@ -2620,7 +2621,7 @@ subroutine adiabatic(h, tv, fluxes, dt, G, GV, US, CS) type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(diabatic_CS), pointer :: CS !< module control structure - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: zeros ! An array of zeros. + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: zeros ! An array of zeros with units of [H ~> m or kg m-2] zeros(:,:,:) = 0.0 @@ -2646,8 +2647,8 @@ subroutine diagnose_diabatic_diff_tendency(tv, h, temp_old, saln_old, dt, G, GV, type(diabatic_CS), pointer :: CS !< module control structure ! Local variables - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: work_3d - real, dimension(SZI_(G),SZJ_(G)) :: work_2d + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: work_3d ! A 3-d work array for diagnostics [various] + real, dimension(SZI_(G),SZJ_(G)) :: work_2d ! A 2-d work array for diagnostics [various] real :: Idt ! The inverse of the timestep [T-1 ~> s-1] real :: ppt2mks ! Conversion factor from S to kg/kg [S-1 ~> ppt-1]. integer :: i, j, k, is, ie, js, je, nz @@ -2741,8 +2742,8 @@ subroutine diagnose_boundary_forcing_tendency(tv, h, temp_old, saln_old, h_old, type(diabatic_CS), pointer :: CS !< module control structure ! Local variables - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: work_3d - real, dimension(SZI_(G),SZJ_(G)) :: work_2d + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: work_3d ! A 3-d work array for diagnostics [various] + real, dimension(SZI_(G),SZJ_(G)) :: work_2d ! A 2-d work array for diagnostics [various] real :: Idt ! The inverse of the timestep [T-1 ~> s-1] real :: ppt2mks ! Conversion factor from S to kg/kg [S-1 ~> ppt-1]. integer :: i, j, k, is, ie, js, je, nz @@ -2828,8 +2829,8 @@ subroutine diagnose_frazil_tendency(tv, h, temp_old, dt, G, GV, US, CS) type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type type(diabatic_CS), pointer :: CS !< module control structure - real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: work_3d - real, dimension(SZI_(G),SZJ_(G)) :: work_2d + real, dimension(SZI_(G),SZJ_(G),SZK_(GV)) :: work_3d ! A 3-d work array for diagnostics [various] + real, dimension(SZI_(G),SZJ_(G)) :: work_2d ! A 2-d work array for diagnostics [various] real :: Idt ! The inverse of the timestep [T-1 ~> s-1] integer :: i, j, k, is, ie, js, je, nz @@ -2942,10 +2943,11 @@ subroutine diabatic_driver_init(Time, G, GV, US, param_file, useALEalgorithm, di !! tracer flow control module type(sponge_CS), pointer :: sponge_CSp !< pointer to the sponge module control structure type(ALE_sponge_CS), pointer :: ALE_sponge_CSp !< pointer to the ALE sponge module control structure - type(oda_incupd_CS), pointer :: oda_incupd_CSp !< pointer to the oda incupd module control structure + type(oda_incupd_CS), pointer :: oda_incupd_CSp !< pointer to the ocean data assimilation incremental + !! update module control structure ! Local variables - real :: Kd ! A diffusivity used in the default for other tracer diffusivities, in MKS units [m2 s-1] + real :: Kd ! A diffusivity used in the default for other tracer diffusivities [Z2 T-1 ~> m2 s-1] logical :: use_temperature character(len=20) :: EN1, EN2, EN3 @@ -3082,11 +3084,12 @@ subroutine diabatic_driver_init(Time, G, GV, US, param_file, useALEalgorithm, di "KD_MIN_TR were operating.", default=.false., do_not_log=.not.CS%useALEalgorithm) if (CS%mix_boundary_tracers .or. CS%mix_boundary_tracer_ALE) then - call get_param(param_file, mdl, "KD", Kd, default=0.0) + call get_param(param_file, mdl, "KD", Kd, units="m2 s-1", default=0.0, scale=US%m2_s_to_Z2_T) call get_param(param_file, mdl, "KD_MIN_TR", CS%Kd_min_tr, & "A minimal diffusivity that should always be applied to "//& "tracers, especially in massless layers near the bottom. "//& - "The default is 0.1*KD.", units="m2 s-1", default=0.1*Kd, scale=US%m2_s_to_Z2_T) + "The default is 0.1*KD.", & + units="m2 s-1", default=0.1*Kd*US%Z2_T_to_m2_s, scale=US%m2_s_to_Z2_T) call get_param(param_file, mdl, "KD_BBL_TR", CS%Kd_BBL_tr, & "A bottom boundary layer tracer diffusivity that will "//& "allow for explicitly specified bottom fluxes. The "//& @@ -3280,9 +3283,9 @@ subroutine diabatic_driver_init(Time, G, GV, US, param_file, useALEalgorithm, di endif - ! diagnostics for tendencies of temp and saln due to diabatic processes + ! Diagnostics for tendencies of temperature and salinity due to diabatic processes, ! available only for ALE algorithm. - ! diagnostics for tendencies of temp and heat due to frazil + ! Diagnostics for tendencies of temperature and heat due to frazil CS%id_diabatic_diff_h = register_diag_field('ocean_model', 'diabatic_diff_h', diag%axesTL, Time, & 'Cell thickness used during diabatic diffusion', & thickness_units, conversion=GV%H_to_MKS, v_extensive=.true.) @@ -3354,9 +3357,9 @@ subroutine diabatic_driver_init(Time, G, GV, US, param_file, useALEalgorithm, di CS%diabatic_diff_tendency_diag = .true. endif - ! diagnostics for tendencies of thickness temp and saln due to boundary forcing + ! Diagnostics for tendencies of thickness temperature and salinity due to boundary forcing, ! available only for ALE algorithm. - ! diagnostics for tendencies of temp and heat due to frazil + ! Diagnostics for tendencies of temperature and heat due to frazil CS%id_boundary_forcing_h = register_diag_field('ocean_model', 'boundary_forcing_h', diag%axesTL, Time, & 'Cell thickness after applying boundary forcing', & thickness_units, conversion=GV%H_to_MKS, v_extensive=.true.) @@ -3593,8 +3596,8 @@ end subroutine diabatic_driver_end !! calculated flux of the layer above and an estimated flux in the !! layer below. This flux is subject to the following conditions: !! (1) the flux in the top and bottom layers are set by the boundary -!! conditions, and (2) no layer may be driven below an Angstrom thick- -!! ness. If there is a bulk mixed layer, the buffer layer is treated +!! conditions, and (2) no layer may be driven below a minimal thickness. +!! If there is a bulk mixed layer, the buffer layer is treated !! as a fixed density layer with vanishingly small diffusivity. !! !! diabatic takes 5 arguments: the two velocities (u and v), the diff --git a/src/parameterizations/vertical/MOM_entrain_diffusive.F90 b/src/parameterizations/vertical/MOM_entrain_diffusive.F90 index beb207624a..792c30cc98 100644 --- a/src/parameterizations/vertical/MOM_entrain_diffusive.F90 +++ b/src/parameterizations/vertical/MOM_entrain_diffusive.F90 @@ -113,7 +113,8 @@ subroutine entrainment_diffusive(h, tv, fluxes, dt, G, GV, US, CS, ea, eb, & diff_work ! The work actually done by diffusion across each ! interface [R Z3 T-3 ~> W m-2]. Sum vertically for the total work. - real :: hm, fm, fr, fk ! Work variables with units of H, H, H, and H2. + real :: hm, fm, fr ! Work variables with units of [H ~> m or kg m-2]. + real :: fk ! A Work variable with units of [H2 ~> m2 or kg2 m-4] real :: b1(SZI_(G)) ! A variable used by the tridiagonal solver [H ~> m or kg m-2] real :: c1(SZI_(G),SZK_(GV)) ! A variable used by the tridiagonal solver [nondim] @@ -140,9 +141,11 @@ subroutine entrainment_diffusive(h, tv, fluxes, dt, G, GV, US, CS, ea, eb, & zeros, & ! An array of all zeros. (Usually used with [H ~> m or kg m-2].) max_eakb, & ! The maximum value of eakb that might be realized [H ~> m or kg m-2]. min_eakb, & ! The minimum value of eakb that might be realized [H ~> m or kg m-2]. - err_max_eakb0, & ! The value of error returned by determine_Ea_kb - err_min_eakb0, & ! when eakb = min_eakb and max_eakb and ea_kbp1 = 0. - err_eakb0, & ! A value of error returned by determine_Ea_kb. + err_max_eakb0, & ! The value of error returned by determine_Ea_kb when eakb = max_eakb + ! and ea_kbp1 = 0 [H2 ~> m2 or kg2 m-4]. + err_min_eakb0, & ! The value of error returned by determine_Ea_kb when eakb = min_eakb + ! and ea_kbp1 = 0 [H2 ~> m2 or kg2 m-4]. + err_eakb0, & ! A value of error returned by determine_Ea_kb [H2 ~> m2 or kg2 m-4]. F_kb, & ! The value of F in layer kb, or equivalently the entrainment ! from below by layer kb [H ~> m or kg m-2]. dFdfm_kb, & ! The partial derivative of F with fm [nondim]. See dFdfm. @@ -187,7 +190,7 @@ subroutine entrainment_diffusive(h, tv, fluxes, dt, G, GV, US, CS, ea, eb, & ! entrain from the layer above [H ~> m or kg m-2]. real :: Kd_here ! The effective diapycnal diffusivity times the timestep [H2 ~> m2 or kg2 m-4]. real :: h_avail ! The thickness that is available for entrainment [H ~> m or kg m-2]. - real :: dS_kb_eff ! The value of dS_kb after limiting is taken into account. + real :: dS_kb_eff ! The value of dS_kb after limiting is taken into account [R ~> kg m-3]. real :: Rho_cor ! The depth-integrated potential density anomaly that ! needs to be corrected for [H R ~> kg m-2 or kg2 m-5]. real :: ea_cor ! The corrective adjustment to eakb [H ~> m or kg m-2]. @@ -752,7 +755,7 @@ subroutine entrainment_diffusive(h, tv, fluxes, dt, G, GV, US, CS, ea, eb, & ea(i,j,k) = ea(i,j,k) + ea_cor eb(i,j,k) = eb(i,j,k) - (dS_kb(i) * I_dSkbp1(i)) * ea_cor elseif (k < kb(i)) then - ! Repetative, unless ea(kb) has been corrected. + ! Repetitive, unless ea(kb) has been corrected. ea(i,j,k) = ea(i,j,k+1) endif enddo @@ -761,7 +764,7 @@ subroutine entrainment_diffusive(h, tv, fluxes, dt, G, GV, US, CS, ea, eb, & ea(i,j,k) = ea(i,j,k+1) enddo ; enddo - ! Repetative, unless ea(kb) has been corrected. + ! Repetitive, unless ea(kb) has been corrected. k=kmb do i=is,ie ! Do not adjust eb through the base of the buffer layers, but it @@ -909,7 +912,7 @@ subroutine F_to_ent(F, h, kb, kmb, j, G, GV, CS, dsp1_ds, eakb, Ent_bl, ea, eb) real, dimension(SZI_(G),SZK_(GV)), intent(in) :: dsp1_ds !< The ratio of coordinate variable !! differences across the interfaces below !! a layer over the difference across the - !! interface above the layer. + !! interface above the layer [nondim]. real, dimension(SZI_(G)), intent(in) :: eakb !< The entrainment from above by the layer !! below the buffer layer [H ~> m or kg m-2]. real, dimension(SZI_(G),SZK_(GV)), intent(in) :: Ent_bl !< The average entrainment upward and @@ -1232,13 +1235,14 @@ subroutine determine_dSkb(h_bl, Sref, Ent_bl, E_kb, is, ie, kmb, G, GV, limit, & ! Local variables real, dimension(SZI_(G),SZK_(GV)) :: & - b1, c1, & ! b1 and c1 are variables used by the tridiagonal solver. - S, dS_dE, & ! The coordinate density [R ~> kg m-3] and its derivative with E. - ea, dea_dE, & ! The entrainment from above and its derivative with E. - eb, deb_dE ! The entrainment from below and its derivative with E. - real :: deriv_dSkb(SZI_(G)) - real :: d1(SZI_(G)) ! d1 = 1.0-c1 is also used by the tridiagonal solver. - real :: src ! A source term for dS_dR. + b1, c1, & ! b1 [H-1 ~> m-1 or m2 kg-1] and c1 [nondim] are variables used by the tridiagonal solver. + S, dS_dE, & ! The coordinate density [R ~> kg m-3] and its derivative with E [R H-1 ~> kg m-4 or m-1]. + ea, dea_dE, & ! The entrainment from above [H ~> m or kg m-2] and its derivative with E [nondim]. + eb, deb_dE ! The entrainment from below [H ~> m or kg m-2] and its derivative with E [nondim]. + real :: deriv_dSkb(SZI_(G)) ! The limited derivative of the new density difference across the base of + ! the buffer layers with the new density of the bottommost buffer layer [nondim] + real :: d1(SZI_(G)) ! d1 = 1.0-c1 is also used by the tridiagonal solver [nondim]. + real :: src ! A source term for dS_dR [R ~> kg m-3]. real :: h1 ! The thickness in excess of the minimum that will remain ! after exchange with the layer below [H ~> m or kg m-2]. logical, dimension(SZI_(G)) :: do_i @@ -1247,13 +1251,15 @@ subroutine determine_dSkb(h_bl, Sref, Ent_bl, E_kb, is, ie, kmb, G, GV, limit, & real :: h_tr ! h_tr is h at tracer points with a tiny thickness ! added to ensure positive definiteness [H ~> m or kg m-2]. real :: b_denom_1 ! The first term in the denominator of b1 [H ~> m or kg m-2]. - real :: rat - real :: dS_kbp1, IdS_kbp1 - real :: deriv_dSLay - real :: Inv_term ! [nondim] + real :: rat ! A ratio of density differences [nondim] + real :: dS_kbp1 ! The density difference between the top two interior layers [R ~> kg m-3]. + real :: IdS_kbp1 ! The inverse of dS_kbp1 [R-1 ~> m3 kg-1] + real :: deriv_dSLay ! The derivative of the projected density difference across the topmost interior + ! layer with the density difference across the interface above it [nondim] + real :: Inv_term ! The inverse of a nondimensional expression [nondim] real :: f1, df1_drat ! Temporary variables [nondim]. real :: z, dz_drat, f2, df2_dz, expz ! Temporary variables [nondim]. - real :: eps_dSLay, eps_dSkb ! Small nondimensional constants. + real :: eps_dSLay, eps_dSkb ! Small nondimensional constants [nondim]. integer :: i, k if (present(ddSlay_dE) .and. .not.present(dSlay)) call MOM_error(FATAL, & @@ -1447,16 +1453,21 @@ subroutine F_kb_to_ea_kb(h_bl, Sref, Ent_bl, I_dSkbp1, F_kb, kmb, i, & real, optional, intent(in) :: tol_in !< A tolerance for the iterative determination !! of the entrainment [H ~> m or kg m-2]. - real :: max_ea, min_ea - real :: err, err_min, err_max - real :: derr_dea - real :: val, tolerance, tol1 - real :: ea_prev - real :: dS_kbp1 - logical :: bisect_next, Newton - real, dimension(SZI_(G)) :: dS_kb - real, dimension(SZI_(G)) :: maxF, ent_maxF, zeros - real, dimension(SZI_(G)) :: ddSkb_dE + real :: max_ea, min_ea ! Bounds on the estimated entraiment [H ~> m or kg m-2] + real :: err, err_min, err_max ! Errors in the mass flux balance [H R ~> kg m-2 or kg2 m-5] + real :: derr_dea ! The change in error with the change in ea [R ~> kg m-3] + real :: val ! An estimate mass flux [H R ~> kg m-2 or kg2 m-5] + real :: tolerance, tol1 ! Tolerances for the determination of the entrainment [H ~> m or kg m-2] + real :: ea_prev ! A previous estimate of ea_kb [H ~> m or kg m-2] + real :: dS_kbp1 ! The density difference between two interior layers [R ~> kg m-3] + real :: dS_kb(SZI_(G)) ! The limited potential density difference across the interface + ! between the bottommost buffer layer and the topmost interior layer [R ~> kg m-3] + real :: maxF(SZI_(G)) ! The maximum value of F (the density flux divided by density + ! differences) found in the range min_ent < ent < max_ent [H ~> m or kg m-2]. + real :: ent_maxF(SZI_(G)) ! The value of entrainment that gives maxF [H ~> m or kg m-2] + real :: zeros(SZI_(G)) ! An array of zero entrainments [H ~> m or kg m-2] + real :: ddSkb_dE(SZI_(G)) ! The partial derivative of dS_kb with ea_kb [R H-1 ~> kg m-4 or m-1] + logical :: bisect_next, Newton ! These indicate what method the next iteration should use integer :: it integer, parameter :: MAXIT = 30 @@ -1589,13 +1600,15 @@ subroutine determine_Ea_kb(h_bl, dtKd_kb, Sref, I_dSkbp1, Ent_bl, ea_kbp1, & !! The input value is the first guess. real, dimension(SZI_(G)), optional, intent(out) :: error !< The error (locally defined in this !! routine) associated with the returned - !! solution. + !! solution [H2 ~> m2 or kg2 m-4] real, dimension(SZI_(G)), optional, intent(in) :: err_min_eakb0 !< The errors (locally defined) !! associated with min_eakb when ea_kbp1 = 0, - !! returned from a previous call to this fn. + !! returned from a previous call to this + !! subroutine [H2 ~> m2 or kg2 m-4]. real, dimension(SZI_(G)), optional, intent(in) :: err_max_eakb0 !< The errors (locally defined) !! associated with min_eakb when ea_kbp1 = 0, - !! returned from a previous call to this fn. + !! returned from a previous call to this + !! subroutine [H2 ~> m2 or kg2 m-4]. real, dimension(SZI_(G)), optional, intent(out) :: F_kb !< The entrainment from below by the !! uppermost interior layer !! corresponding to the returned @@ -1719,7 +1732,7 @@ subroutine determine_Ea_kb(h_bl, dtKd_kb, Sref, I_dSkbp1, Ent_bl, ea_kbp1, & Ent(i) = Ent(i) - err(i) / derror_dE(i) elseif (false_position(i) .and. & (error_maxE(i) - error_minE(i) < 0.9*large_err)) then - ! Use the false postion method if there are decent error estimates. + ! Use the false position method if there are decent error estimates. Ent(i) = E_min(i) + (E_max(i)-E_min(i)) * & (-error_minE(i)/(error_maxE(i) - error_minE(i))) false_position(i) = .false. @@ -1813,17 +1826,21 @@ subroutine find_maxF_kb(h_bl, Sref, Ent_bl, I_dSkbp1, min_ent_in, max_ent_in, & ! negative) value. It is faster to find the true maximum by first finding the ! unlimited maximum and comparing it to the limited value at max_ent_in. real, dimension(SZI_(G)) :: & - ent, & - minent, maxent, ent_best, & - F_max_ent_in, & - F_maxent, F_minent, F, F_best, & - dF_dent, dF_dE_max, dF_dE_min, dF_dE_best, & - dS_kb, dS_kb_lim, ddSkb_dE, dS_anom_lim, & - chg_prev, chg_pre_prev - real :: dF_dE_mean, maxslope, minslope - real :: tolerance - real :: ratio_select_end - real :: rat, max_chg, min_chg, chg1, chg2, chg + ent, & ! The updated estimate of the entrainment [H ~> m or kg m-2] + minent, maxent, ent_best, & ! Various previous estimates of the entrainment [H ~> m or kg m-2] + F_max_ent_in, & ! The value of F that gives the input maximum value of ent [H ~> m or kg m-2] + F_maxent, F_minent, F, F_best, & ! Various estimates of F [H ~> m or kg m-2] + dF_dent, dF_dE_max, dF_dE_min, dF_dE_best, & ! Various derivatives of F with ent [nondim] + dS_kb, & ! The density difference across the interface between the bottommost + ! buffer layer and the topmost interior layer [R ~> kg m-3] + dS_kb_lim, dS_anom_lim, & ! Various limits on dS_kb [R ~> kg m-3] + ddSkb_dE, & ! The partial derivative of dS_kb with ent [R H-1 ~> kg m-4 or m-1]. + chg_prev, chg_pre_prev ! Changes in estimates of the entrainment from previous iterations [H ~> m or kg m-2] + real :: dF_dE_mean, maxslope, minslope ! Various derivatives of F with ent [nondim] + real :: tolerance ! The tolerance within which ent must be converged [H ~> m or kg m-2] + real :: ratio_select_end, rat ! Fractional changes in the value of ent to use for the next iteration + ! relative to its bounded range [nondim] + real :: max_chg, min_chg, chg1, chg2, chg ! Changes in entrainment estimates [H ~> m or kg m-2] logical, dimension(SZI_(G)) :: do_i, last_it, need_bracket, may_use_best logical :: doany, OK1, OK2, bisect, new_min_bound integer :: i, it, is1, ie1 @@ -1876,14 +1893,14 @@ subroutine find_maxF_kb(h_bl, Sref, Ent_bl, I_dSkbp1, min_ent_in, max_ent_in, & maxslope = MAX(dF_dE_mean, dF_dE_min(i), dF_dE_max(i)) minslope = MIN(dF_dE_mean, dF_dE_min(i), dF_dE_max(i)) if (F_minent(i) >= F_maxent(i)) then - if (dF_dE_min(i) > 0.0) then ; rat = 0.02 ! A small step should bracket the soln. + if (dF_dE_min(i) > 0.0) then ; rat = 0.02 ! A small step should bracket the solution. elseif (maxslope < ratio_select_end*minslope) then ! The maximum of F is at minent. F_best(i) = F_minent(i) ; ent_best(i) = minent(i) ; rat = 0.0 do_i(i) = .false. else ; rat = 0.382 ; endif ! Use the golden ratio else - if (dF_dE_max(i) < 0.0) then ; rat = 0.98 ! A small step should bracket the soln. + if (dF_dE_max(i) < 0.0) then ; rat = 0.98 ! A small step should bracket the solution. elseif (minslope > ratio_select_end*maxslope) then ! The maximum of F is at maxent. F_best(i) = F_maxent(i) ; ent_best(i) = maxent(i) ; rat = 1.0 @@ -1979,7 +1996,7 @@ subroutine find_maxF_kb(h_bl, Sref, Ent_bl, I_dSkbp1, min_ent_in, max_ent_in, & new_min_bound = .true. ! We have a new minimum bound. elseif ((F(i) <= F_maxent(i)) .and. (F(i) > F_minent(i))) then new_min_bound = .false. ! We have a new maximum bound. - else ! This case would bracket a minimum. Wierd. + else ! This case would bracket a minimum. Weird. ! Unless the derivative indicates that there is a maximum near the ! lower bound, try keeping the end with the larger value of F ! in a tie keep the minimum as the answer here will be compared @@ -2068,14 +2085,14 @@ subroutine entrain_diffusive_init(Time, G, GV, US, param_file, diag, CS, just_re !! parameters. type(diag_ctrl), target, intent(inout) :: diag !< A structure that is used to regulate diagnostic !! output. - type(entrain_diffusive_CS), intent(inout) :: CS !< Entrainment diffusion control struct + type(entrain_diffusive_CS), intent(inout) :: CS !< Entrainment diffusion control structure logical, intent(in) :: just_read_params !< If true, this call will only read !! and log parameters without registering !! any diagnostics ! Local variables - real :: dt ! The dynamics timestep, used here in the default for TOLERANCE_ENT, in MKS units [s] - real :: Kd ! A diffusivity used in the default for TOLERANCE_ENT, in MKS units [m2 s-1] + real :: dt ! The dynamics timestep, used here in the default for TOLERANCE_ENT [T ~> s] + real :: Kd ! A diffusivity used in the default for TOLERANCE_ENT [Z2 T-1 ~> m2 s-1] ! This include declares and sets the variable "version". # include "version_variable.h" character(len=40) :: mdl = "MOM_entrain_diffusive" ! This module's name. @@ -2090,14 +2107,14 @@ subroutine entrain_diffusive_init(Time, G, GV, US, param_file, diag, CS, just_re call get_param(param_file, mdl, "MAX_ENT_IT", CS%max_ent_it, & "The maximum number of iterations that may be used to "//& "calculate the interior diapycnal entrainment.", default=5, do_not_log=just_read_params) - ! In this module, KD is only used to set the default for TOLERANCE_ENT. [m2 s-1] - call get_param(param_file, mdl, "KD", Kd, default=0.0) + ! In this module, KD is only used to set the default for TOLERANCE_ENT. [Z2 T-1 ~> m2 s-1] + call get_param(param_file, mdl, "KD", Kd, units="m2 s-1", default=0.0, scale=US%m2_s_to_Z2_T) call get_param(param_file, mdl, "DT", dt, & - "The (baroclinic) dynamics time step.", units = "s", & - fail_if_missing=.true., do_not_log=just_read_params) + "The (baroclinic) dynamics time step.", & + units="s", scale=US%s_to_T, fail_if_missing=.true., do_not_log=just_read_params) call get_param(param_file, mdl, "TOLERANCE_ENT", CS%Tolerance_Ent, & "The tolerance with which to solve for entrainment values.", & - units="m", default=MAX(100.0*GV%Angstrom_m,1.0e-4*sqrt(dt*Kd)), scale=GV%m_to_H, & + units="m", default=MAX(100.0*GV%Angstrom_m,1.0e-4*sqrt(dt*Kd)*US%Z_to_m), scale=GV%m_to_H, & do_not_log=just_read_params) CS%Rho_sig_off = 1000.0*US%kg_m3_to_R @@ -2119,10 +2136,10 @@ end subroutine entrain_diffusive_init !! mixing and advection in isopycnal layers. The main subroutine, !! calculate_entrainment, returns the entrainment by each layer !! across the interfaces above and below it. These are calculated -!! subject to the constraints that no layers can be driven to neg- -!! ative thickness and that the each layer maintains its target -!! density, using the scheme described in Hallberg (MWR 2000). There -!! may or may not be a bulk mixed layer above the isopycnal layers. +!! subject to the constraints that no layers can be driven to negative +!! thickness and that the each layer maintains its target density, +!! using the scheme described in Hallberg (MWR 2000). There may or +!! may not be a bulk mixed layer above the isopycnal layers. !! The solution is iterated until the change in the entrainment !! between successive iterations is less than some small tolerance. !! @@ -2134,9 +2151,9 @@ end subroutine entrain_diffusive_init !! diffusion, so the fully implicit upwind differencing scheme that !! is used is entirely appropriate. The downward buoyancy flux in !! each layer is determined from an implicit calculation based on -!! the previously calculated flux of the layer above and an estim- -!! ated flux in the layer below. This flux is subject to the foll- -!! owing conditions: (1) the flux in the top and bottom layers are +!! the previously calculated flux of the layer above and an estimated +!! flux in the layer below. This flux is subject to the following +!! conditions: (1) the flux in the top and bottom layers are !! set by the boundary conditions, and (2) no layer may be driven !! below an Angstrom thickness. If there is a bulk mixed layer, the !! mixed and buffer layers are treated as Eulerian layers, whose diff --git a/src/parameterizations/vertical/MOM_kappa_shear.F90 b/src/parameterizations/vertical/MOM_kappa_shear.F90 index c088eea5bb..1e8015eacd 100644 --- a/src/parameterizations/vertical/MOM_kappa_shear.F90 +++ b/src/parameterizations/vertical/MOM_kappa_shear.F90 @@ -356,8 +356,8 @@ subroutine Calc_kappa_shear_vertex(u_in, v_in, h, T_in, S_in, tv, p_surf, kappa_ real, dimension(SZIB_(G),SZJB_(G),SZK_(GV)+1), & intent(inout) :: kv_io !< The vertical viscosity at each interface [Z2 T-1 ~> m2 s-1]. !! The previous value is used to initialize kappa - !! in the vertex columes as Kappa = Kv/Prandtl - !! to accelerate the iteration toward covergence. + !! in the vertex columns as Kappa = Kv/Prandtl + !! to accelerate the iteration toward convergence. real, intent(in) :: dt !< Time increment [T ~> s]. type(Kappa_shear_CS), pointer :: CS !< The control structure returned by a previous !! call to kappa_shear_init. @@ -650,7 +650,7 @@ subroutine kappa_shear_column(kappa, tke, dt, nzc, f2, surface_pres, dz, & S2, & ! The squared shear at an interface [T-2 ~> s-2]. a1, & ! a1 is the coupling between adjacent interfaces in the TKE, ! velocity, and density equations [Z s-1 ~> m s-1] or [Z ~> m] - c1, & ! c1 is used in the tridiagonal (and similar) solvers. + c1, & ! c1 is used in the tridiagonal (and similar) solvers [nondim]. k_src, & ! The shear-dependent source term in the kappa equation [T-1 ~> s-1]. kappa_src, & ! The shear-dependent source term in the kappa equation [T-1 ~> s-1]. kappa_out, & ! The kappa that results from the kappa equation [Z2 T-1 ~> m2 s-1]. @@ -675,9 +675,9 @@ subroutine kappa_shear_column(kappa, tke, dt, nzc, f2, surface_pres, dz, & ! sources from the elliptic term [T-1 ~> s-1]. real :: dist_from_bot ! The distance from the bottom surface [Z ~> m]. - real :: b1 ! The inverse of the pivot in the tridiagonal equations. - real :: bd1 ! A term in the denominator of b1. - real :: d1 ! 1 - c1 in the tridiagonal equations. + real :: b1 ! The inverse of the pivot in the tridiagonal equations [Z-1 ~> m-1]. + real :: bd1 ! A term in the denominator of b1 [Z ~> m]. + real :: d1 ! 1 - c1 in the tridiagonal equations [nondim] real :: gR0 ! A conversion factor from Z to pressure, given by Rho_0 times g ! [R L2 T-2 Z-1 ~> kg m-2 s-2]. real :: g_R0 ! g_R0 is a rescaled version of g/Rho [Z R-1 T-2 ~> m4 kg-1 s-2]. @@ -1060,10 +1060,13 @@ subroutine calculate_projected_state(kappa, u0, v0, T0, S0, dt, nz, dz, I_dz_int !! diffusivity. ! Local variables - real, dimension(nz+1) :: c1 - real :: L2_to_Z2 ! A conversion factor from horizontal length units to vertical depth - ! units squared [Z2 s2 T-2 m-2 ~> 1]. - real :: a_a, a_b, b1, d1, bd1, b1nz_0 + real, dimension(nz+1) :: c1 ! A tridiagonal variable [nondim] + real :: L2_to_Z2 ! A conversion factor from horizontal length units to vertical depth + ! units squared [Z2 s2 T-2 m-2 ~> 1]. + real :: a_a, a_b ! Tridiagonal coupling coefficients [Z ~> m] + real :: b1, b1nz_0 ! Tridiagonal variables [Z-1 ~> m-1] + real :: bd1 ! A term in the denominator of b1 [Z ~> m] + real :: d1 ! A tridiagonal variable [nondim] integer :: k, ks, ke ks = 1 ; ke = nz @@ -1166,7 +1169,7 @@ subroutine find_kappa_tke(N2, S2, kappa_in, Idz, dz_Int, I_L2_bdry, f2, & real, dimension(nz+1), intent(in) :: kappa_in !< The initial guess at the diffusivity !! [Z2 T-1 ~> m2 s-1]. real, dimension(nz+1), intent(in) :: dz_Int !< The thicknesses associated with interfaces - !! [Z-1 ~> m-1]. + !! [Z ~> m]. real, dimension(nz+1), intent(in) :: I_L2_bdry !< The inverse of the squared distance to !! boundaries [Z-2 ~> m-2]. real, dimension(nz), intent(in) :: Idz !< The inverse grid spacing of layers [Z-1 ~> m-1]. @@ -1203,7 +1206,7 @@ subroutine find_kappa_tke(N2, S2, kappa_in, Idz, dz_Int, I_L2_bdry, f2, & dQmdK, & ! With Newton's method the change in dQ(k-1) due to dK(k) [T ~> s]. dKdQ, & ! With Newton's method the change in dK(k) due to dQ(k) [T-1 ~> s-1]. e1 ! The fractional change in a layer TKE due to a change in the - ! TKE of the layer above when all the kappas below are 0. + ! TKE of the layer above when all the kappas below are 0 [nondim]. ! e1 is nondimensional, and 0 < e1 < 1. real :: tke_src ! The net source of TKE due to mixing against the shear ! and stratification [Z2 T-3 ~> m2 s-3]. (For convenience, @@ -1213,13 +1216,13 @@ subroutine find_kappa_tke(N2, S2, kappa_in, Idz, dz_Int, I_L2_bdry, f2, & real :: bK ! The inverse of the pivot in the tridiagonal equations [Z-1 ~> m-1]. real :: bQd1 ! A term in the denominator of bQ [Z T-1 ~> m s-1]. real :: bKd1 ! A term in the denominator of bK [Z ~> m]. - real :: cQcomp, cKcomp ! 1 - cQ or 1 - cK in the tridiagonal equations. + real :: cQcomp, cKcomp ! 1 - cQ or 1 - cK in the tridiagonal equations [nondim]. real :: c_s2 ! The coefficient for the decay of TKE due to - ! shear (i.e. proportional to |S|*tke), nondimensional. + ! shear (i.e. proportional to |S|*tke) [nondim]. real :: c_n2 ! The coefficient for the decay of TKE due to ! stratification (i.e. proportional to N*tke) [nondim]. real :: Ri_crit ! The critical shear Richardson number for shear- - ! driven mixing. The theoretical value is 0.25. + ! driven mixing [nondim]. The theoretical value is 0.25. real :: q0 ! The background level of TKE [Z2 T-2 ~> m2 s-2]. real :: Ilambda2 ! 1.0 / CS%lambda**2 [nondim] real :: TKE_min ! The minimum value of shear-driven TKE that can be @@ -1227,31 +1230,33 @@ subroutine find_kappa_tke(N2, S2, kappa_in, Idz, dz_Int, I_L2_bdry, f2, & real :: kappa0 ! The background diapycnal diffusivity [Z2 T-1 ~> m2 s-1]. real :: kappa_trunc ! Diffusivities smaller than this are rounded to 0 [Z2 T-1 ~> m2 s-1]. - real :: eden1, eden2, I_eden, ome ! Variables used in calculating e1. + real :: eden1, eden2 ! Variables used in calculating e1 [Z-1 ~> m-1] + real :: I_eden ! The inverse of the denominator in e1 [Z ~> m] + real :: ome ! Variables used in calculating e1 [nondim] real :: diffusive_src ! The diffusive source in the kappa equation [Z T-1 ~> m s-1]. real :: chg_by_k0 ! The value of k_src that leads to an increase of ! kappa_0 if only the diffusive term is a sink [T-1 ~> s-1]. real :: kappa_mean ! A mean value of kappa [Z2 T-1 ~> m2 s-1]. real :: Newton_test ! The value of relative error that will cause the next - ! iteration to use Newton's method. + ! iteration to use Newton's method [nondim]. ! Temporary variables used in the Newton's method iterations. - real :: decay_term_k ! The decay term in the diffusivity equation + real :: decay_term_k ! The decay term in the diffusivity equation [Z-1 ~> m-1] real :: decay_term_Q ! The decay term in the TKE equation - proportional to [T-1 ~> s-1] real :: I_Q ! The inverse of TKE [T2 Z-2 ~> s2 m-2] - real :: kap_src + real :: kap_src ! A source term in the kappa equation [Z T-1 ~> m s-1] real :: v1 ! A temporary variable proportional to [T-1 ~> s-1] - real :: v2 - real :: tol_err ! The tolerance for max_err that determines when to - ! stop iterating. - real :: Newton_err ! The tolerance for max_err that determines when to - ! start using Newton's method. Empirically, an initial - ! value of about 0.2 seems to be most efficient. - real, parameter :: roundoff = 1.0e-16 ! A negligible fractional change in TKE. - ! This could be larger but performance gains are small. + real :: v2 ! A temporary variable in [Z T-2 ~> m s-2] + real :: tol_err ! The tolerance for max_err that determines when to + ! stop iterating [nondim]. + real :: Newton_err ! The tolerance for max_err that determines when to + ! start using Newton's method [nondim]. Empirically, an initial + ! value of about 0.2 seems to be most efficient. + real, parameter :: roundoff = 1.0e-16 ! A negligible fractional change in TKE [nondim]. + ! This could be larger but performance gains are small. logical :: tke_noflux_bottom_BC = .false. ! Specify the boundary conditions - logical :: tke_noflux_top_BC = .false. ! that are applied to the TKE eqns. + logical :: tke_noflux_top_BC = .false. ! that are applied to the TKE equations. logical :: do_Newton ! If .true., use Newton's method for the next iteration. logical :: abort_Newton ! If .true., an Newton's method has encountered a 0 ! pivot, and should not have been used. @@ -1265,7 +1270,8 @@ subroutine find_kappa_tke(N2, S2, kappa_in, Idz, dz_Int, I_L2_bdry, f2, & ! These variables are used only for debugging. logical, parameter :: debug_soln = .false. - real :: K_err_lin, Q_err_lin + real :: K_err_lin ! The imbalance in the K equation [Z T-1 ~> m s-1] + real :: Q_err_lin ! The imbalance in the Q equation [Z2 T-3 ~> m2 s-3] real, dimension(nz+1) :: & I_Ld2_debug, & ! A separate version of I_Ld2 for debugging [Z-2 ~> m-2]. kappa_prev, & ! The value of kappa at the start of the current iteration [Z2 T-1 ~> m2 s-1]. @@ -1726,15 +1732,15 @@ function kappa_shear_init(Time, G, GV, US, param_file, diag, CS) logical :: kappa_shear_init !< True if module is to be used, False otherwise ! Local variables + real :: KD_normal ! The KD of the main model, read here only as a parameter + ! for setting the default of KD_SMOOTH [Z2 T-1 ~> m2 s-1] + real :: kappa_0_default ! The default value for KD_KAPPA_SHEAR_0 [Z2 T-1 ~> m2 s-1] logical :: merge_mixedlayer logical :: debug_shear logical :: just_read ! If true, this module is not used, so only read the parameters. ! This include declares and sets the variable "version". # include "version_variable.h" character(len=40) :: mdl = "MOM_kappa_shear" ! This module's name. - real :: kappa_0_unscaled ! The value of kappa_0 in MKS units [m2 s-1] - real :: KD_normal ! The KD of the main model, read here only as a parameter - ! for setting the default of KD_SMOOTH in MKS units [m2 s-1] if (associated(CS)) then call MOM_error(WARNING, "kappa_shear_init called with an associated "// & @@ -1775,18 +1781,21 @@ function kappa_shear_init(Time, G, GV, US, param_file, diag, CS) call get_param(param_file, mdl, "MAX_RINO_IT", CS%max_RiNo_it, & "The maximum number of iterations that may be used to "//& "estimate the Richardson number driven mixing.", & - default=50, do_not_log=just_read) - call get_param(param_file, mdl, "KD", KD_normal, default=0.0, do_not_log=.true.) + units="nondim", default=50, do_not_log=just_read) + call get_param(param_file, mdl, "KD", KD_normal, & + units="m2 s-1", scale=US%m2_s_to_Z2_T, default=0.0, do_not_log=.true.) + kappa_0_default = max(Kd_normal, 1.0e-7*US%m2_s_to_Z2_T) call get_param(param_file, mdl, "KD_KAPPA_SHEAR_0", CS%kappa_0, & "The background diffusivity that is used to smooth the "//& "density and shear profiles before solving for the "//& "diffusivities. The default is the greater of KD and 1e-7 m2 s-1.", & - units="m2 s-1", default=max(KD_normal, 1.0e-7), scale=US%m2_s_to_Z2_T, & - unscaled=kappa_0_unscaled, do_not_log=just_read) + units="m2 s-1", default=kappa_0_default*US%Z2_T_to_m2_s, scale=US%m2_s_to_Z2_T, & + do_not_log=just_read) call get_param(param_file, mdl, "KD_TRUNC_KAPPA_SHEAR", CS%kappa_trunc, & "The value of shear-driven diffusivity that is considered negligible "//& "and is rounded down to 0. The default is 1% of KD_KAPPA_SHEAR_0.", & - units="m2 s-1", default=0.01*kappa_0_unscaled, scale=US%m2_s_to_Z2_T, do_not_log=just_read) + units="m2 s-1", default=0.01*CS%kappa_0*US%Z2_T_to_m2_s, scale=US%m2_s_to_Z2_T, & + do_not_log=just_read) call get_param(param_file, mdl, "FRI_CURVATURE", CS%FRi_curvature, & "The nondimensional curvature of the function of the "//& "Richardson number in the kappa source term in the "//& @@ -1950,7 +1959,7 @@ end function kappa_shear_at_vertex !! TKE with shear and stratification fixed, then marches the density !! and velocities forward with an adaptive (and aggressive) time step !! in a predictor-corrector-corrector emulation of a trapezoidal -!! scheme. Run-time-settable parameters determine the tolerence to +!! scheme. Run-time-settable parameters determine the tolerance to !! which the kappa and TKE equations are solved and the minimum time !! step that can be taken. diff --git a/src/parameterizations/vertical/MOM_regularize_layers.F90 b/src/parameterizations/vertical/MOM_regularize_layers.F90 index deb1c90ca9..8966c12b79 100644 --- a/src/parameterizations/vertical/MOM_regularize_layers.F90 +++ b/src/parameterizations/vertical/MOM_regularize_layers.F90 @@ -151,10 +151,11 @@ subroutine regularize_surface(h, tv, dt, ea, eb, G, GV, US, CS) ! the coordinate variable, set to P_Ref [R L2 T-2 ~> Pa]. Rcv_tol, & ! A tolerence, relative to the target density differences ! between layers, for detraining into the interior [nondim]. - h_add_tgt, h_add_tot, & - h_tot1, Th_tot1, Sh_tot1, & - h_tot3, Th_tot3, Sh_tot3, & - h_tot2, Th_tot2, Sh_tot2 + h_add_tgt, & ! The target for the thickness to add to the mixed layers [H ~> m or kg m-2] + h_add_tot, & ! The net thickness added to the mixed layers [H ~> m or kg m-2] + h_tot1, h_tot2, h_tot3, & ! Debugging diagnostics of total thicknesses [H ~> m or kg m-2] + Th_tot1, Th_tot2, Th_tot3, & ! Debugging diagnostics of integrated temperatures [C H ~> degC m or degC kg m-2] + Sh_tot1, Sh_tot2, Sh_tot3 ! Debugging diagnostics of integrated salinities [S H ~> ppt m or ppt kg m-2] real, dimension(SZK_(GV)) :: & h_prev_1d ! The previous thicknesses [H ~> m or kg m-2]. real :: I_dtol ! The inverse of the tolerance changes [nondim]. @@ -168,16 +169,17 @@ subroutine regularize_surface(h, tv, dt, ea, eb, G, GV, US, CS) int_flux, & ! Mass flux across the interfaces [H ~> m or kg m-2] int_Tflux, & ! Temperature flux across the interfaces [C H ~> degC m or degC kg m-2] int_Sflux ! Salinity flux across the interfaces [S H ~> ppt m or ppt kg m-2] - real :: h_add - real :: h_det_tot - real :: max_def_rat + real :: h_add ! The thickness to add to the layers above an interface [H ~> m or kg m-2] + real :: h_det_tot ! The total thickness detrained by the mixed layers [H ~> m or kg m-2] + real :: max_def_rat ! The maximum value of the ratio of the thickness deficit to the minimum depth [nondim] real :: Rcv_min_det ! The lightest (min) and densest (max) coordinate density real :: Rcv_max_det ! that can detrain into a layer [R ~> kg m-3]. - real :: int_top, int_bot - real :: h_predicted - real :: h_prev - real :: h_deficit + real :: int_top, int_bot ! The interface depths above and below a layer [H ~> m or kg m-2], positive upward. + real :: h_predicted ! An updated thickness [H ~> m or kg m-2] + real :: h_prev ! The previous thickness [H ~> m or kg m-2] + real :: h_deficit ! The difference between the layer thickness and the value estimated from the + ! filtered interface depths [H ~> m or kg m-2] logical :: cols_left, ent_any, more_ent_i(SZI_(G)), ent_i(SZI_(G)) logical :: det_any, det_i(SZI_(G)) diff --git a/src/parameterizations/vertical/MOM_set_diffusivity.F90 b/src/parameterizations/vertical/MOM_set_diffusivity.F90 index 6d35616b3a..8257d19bd3 100644 --- a/src/parameterizations/vertical/MOM_set_diffusivity.F90 +++ b/src/parameterizations/vertical/MOM_set_diffusivity.F90 @@ -69,8 +69,7 @@ module MOM_set_diffusivity !! drag law c_drag*|u|*u. logical :: BBL_mixing_as_max !< If true, take the maximum of the diffusivity !! from the BBL mixing and the other diffusivities. - !! Otherwise, diffusivities from the BBL_mixing is - !! added. + !! Otherwise, diffusivities from the BBL_mixing is added. logical :: use_LOTW_BBL_diffusivity !< If true, use simpler/less precise, BBL diffusivity. logical :: LOTW_BBL_use_omega !< If true, use simpler/less precise, BBL diffusivity. real :: Von_Karm !< The von Karman constant as used in the BBL diffusivity calculation @@ -115,10 +114,9 @@ module MOM_set_diffusivity !! is the rotation rate of the earth squared. real :: ML_rad_kd_max !< Maximum diapycnal diffusivity due to turbulence !! radiated from the base of the mixed layer [Z2 T-1 ~> m2 s-1]. - real :: ML_rad_efold_coeff !< non-dim coefficient to scale penetration depth - real :: ML_rad_coeff !< coefficient, which scales MSTAR*USTAR^3 to - !! obtain energy available for mixing below - !! mixed layer base [nondim] + real :: ML_rad_efold_coeff !< Coefficient to scale penetration depth [nondim] + real :: ML_rad_coeff !< Coefficient which scales MSTAR*USTAR^3 to obtain energy + !! available for mixing below mixed layer base [nondim] logical :: ML_rad_bug !< If true use code with a bug that reduces the energy available !! in the transition layer by a factor of the inverse of the energy !! deposition lenthscale (in m). @@ -135,7 +133,7 @@ module MOM_set_diffusivity !! of the vertical component of rotation when !! setting the decay scale for mixed layer turbulence. real :: ML_omega_frac !< When setting the decay scale for turbulence, use - !! this fraction of the absolute rotation rate blended + !! this fraction [nondim] of the absolute rotation rate blended !! with the local value of f, as f^2 ~= (1-of)*f^2 + of*4*omega^2. logical :: user_change_diff !< If true, call user-defined code to change diffusivity. logical :: useKappaShear !< If true, use the kappa_shear module to find the @@ -149,9 +147,9 @@ module MOM_set_diffusivity logical :: use_tidal_mixing !< If true, activate tidal mixing diffusivity. logical :: simple_TKE_to_Kd !< If true, uses a simple estimate of Kd/TKE that !! does not rely on a layer-formulation. - real :: Max_Rrho_salt_fingers !< max density ratio for salt fingering + real :: Max_Rrho_salt_fingers !< max density ratio for salt fingering [nondim] real :: Max_salt_diff_salt_fingers !< max salt diffusivity for salt fingers [Z2 T-1 ~> m2 s-1] - real :: Kv_molecular !< molecular visc for double diff convect [Z2 T-1 ~> m2 s-1] + real :: Kv_molecular !< Molecular viscosity for double diffusive convection [Z2 T-1 ~> m2 s-1] integer :: answer_date !< The vintage of the order of arithmetic and expressions in this module's !! calculations. Values below 20190101 recover the answers from the @@ -185,9 +183,9 @@ module MOM_set_diffusivity Kd_work => NULL(), & !< layer integrated work by diapycnal mixing [R Z3 T-3 ~> W m-2] maxTKE => NULL(), & !< energy required to entrain to h_max [Z3 T-3 ~> m3 s-3] Kd_bkgnd => NULL(), & !< Background diffusivity at interfaces [Z2 T-1 ~> m2 s-1] - Kv_bkgnd => NULL(), & !< Viscosity from ackground diffusivity at interfaces [Z2 T-1 ~> m2 s-1] - KT_extra => NULL(), & !< double diffusion diffusivity for temp [Z2 T-1 ~> m2 s-1]. - KS_extra => NULL(), & !< double diffusion diffusivity for saln [Z2 T-1 ~> m2 s-1]. + Kv_bkgnd => NULL(), & !< Viscosity from background diffusivity at interfaces [Z2 T-1 ~> m2 s-1] + KT_extra => NULL(), & !< Double diffusion diffusivity for temperature [Z2 T-1 ~> m2 s-1]. + KS_extra => NULL(), & !< Double diffusion diffusivity for salinity [Z2 T-1 ~> m2 s-1]. drho_rat => NULL() !< The density difference ratio used in double diffusion [nondim]. real, pointer, dimension(:,:,:) :: TKE_to_Kd => NULL() !< conversion rate (~1.0 / (G_Earth + dRho_lay)) between TKE @@ -262,8 +260,8 @@ subroutine set_diffusivity(u, v, h, u_h, v_h, tv, fluxes, optics, visc, dt, Kd_i Kd_int_2d, & !< The interface diffusivities [Z2 T-1 ~> m2 s-1] Kv_bkgnd, & !< The background diffusion related interface viscosities [Z2 T-1 ~> m2 s-1] dRho_int, & !< Locally referenced potential density difference across interfaces [R ~> kg m-3] - KT_extra, & !< Double difusion diffusivity of temperature [Z2 T-1 ~> m2 s-1] - KS_extra !< Double difusion diffusivity of salinity [Z2 T-1 ~> m2 s-1] + KT_extra, & !< Double diffusion diffusivity of temperature [Z2 T-1 ~> m2 s-1] + KS_extra !< Double diffusion diffusivity of salinity [Z2 T-1 ~> m2 s-1] real :: dissip ! local variable for dissipation calculations [Z2 R T-3 ~> W m-3] real :: Omega2 ! squared absolute rotation rate [T-2 ~> s-2] @@ -673,7 +671,7 @@ subroutine find_TKE_to_Kd(h, tv, dRho_int, N2_lay, j, dt, G, GV, US, CS, & type(set_diffusivity_CS), pointer :: CS !< Diffusivity control structure real, dimension(SZI_(G),SZK_(GV)), intent(out) :: TKE_to_Kd !< The conversion rate between the !! TKE dissipated within a layer and the - !! diapycnal diffusivity witin that layer, + !! diapycnal diffusivity within that layer, !! usually (~Rho_0 / (G_Earth * dRho_lay)) !! [Z2 T-1 / Z3 T-3 = T2 Z-1 ~> s2 m-1] real, dimension(SZI_(G),SZK_(GV)), intent(out) :: maxTKE !< The energy required to for a layer to entrain @@ -701,12 +699,11 @@ subroutine find_TKE_to_Kd(h, tv, dRho_int, N2_lay, j, dt, G, GV, US, CS, & Rcv_kmb, & ! coordinate density in the lowest buffer layer [R ~> kg m-3] p_0 ! An array of 0 pressures [R L2 T-2 ~> Pa] - real :: dh_max ! maximum amount of entrainment a layer could - ! undergo before entraining all fluid in the layers - ! above or below [Z ~> m]. + real :: dh_max ! maximum amount of entrainment a layer could undergo before + ! entraining all fluid in the layers above or below [Z ~> m]. real :: dRho_lay ! density change across a layer [R ~> kg m-3] real :: Omega2 ! rotation rate squared [T-2 ~> s-2] - real :: G_Rho0 ! gravitation accel divided by Bouss ref density [Z T-2 R-1 ~> m4 s-2 kg-1] + real :: G_Rho0 ! Gravitational acceleration divided by Boussinesq reference density [Z T-2 R-1 ~> m4 s-2 kg-1] real :: G_IRho0 ! Alternate calculation of G_Rho0 for reproducibility [Z T-2 R-1 ~> m4 s-2 kg-1] real :: I_Rho0 ! inverse of Boussinesq reference density [R-1 ~> m3 kg-1] real :: I_dt ! 1/dt [T-1 ~> s-1] @@ -911,16 +908,16 @@ subroutine find_N2(h, tv, T_f, S_f, fluxes, j, G, GV, US, CS, dRho_int, & z_from_bot ! The hieght above the bottom [Z ~> m]. real :: dz_int ! thickness associated with an interface [Z ~> m]. - real :: G_Rho0 ! gravitation acceleration divided by Bouss reference density + real :: G_Rho0 ! Gravitational acceleration divided by Boussinesq reference density ! times some unit conversion factors [Z T-2 R-1 ~> m4 s-2 kg-1]. - real :: H_neglect ! negligibly small thickness, in the same units as h. + real :: H_neglect ! A negligibly small thickness [H ~> m or kg m-2] logical :: do_i(SZI_(G)), do_any integer, dimension(2) :: EOSdom ! The i-computational domain for the equation of state integer :: i, k, is, ie, nz is = G%isc ; ie = G%iec ; nz = GV%ke - G_Rho0 = (US%L_to_Z**2 * GV%g_Earth) / (GV%Rho0) + G_Rho0 = (US%L_to_Z**2 * GV%g_Earth) / GV%Rho0 H_neglect = GV%H_subroundoff ! Find the (limited) density jump across each interface. @@ -1064,8 +1061,8 @@ subroutine double_diffusion(tv, h, T_f, S_f, j, G, GV, US, CS, Kd_T_dd, Kd_S_dd) !! diffusivity for saln [Z2 T-1 ~> m2 s-1]. real, dimension(SZI_(G)) :: & - dRho_dT, & ! partial derivatives of density wrt temp [R C-1 ~> kg m-3 degC-1] - dRho_dS, & ! partial derivatives of density wrt saln [R S-1 ~> kg m-3 ppt-1] + dRho_dT, & ! partial derivatives of density with respect to temperature [R C-1 ~> kg m-3 degC-1] + dRho_dS, & ! partial derivatives of density with respect to salinity [R S-1 ~> kg m-3 ppt-1] pres, & ! pressure at each interface [R L2 T-2 ~> Pa] Temp_int, & ! temperature at interfaces [C ~> degC] Salin_int ! Salinity at interfaces [S ~> ppt] @@ -1076,7 +1073,7 @@ subroutine double_diffusion(tv, h, T_f, S_f, j, G, GV, US, CS, Kd_T_dd, Kd_S_dd) real :: Rrho ! vertical density ratio [nondim] real :: diff_dd ! factor for double-diffusion [nondim] real :: Kd_dd ! The dominant double diffusive diffusivity [Z2 T-1 ~> m2 s-1] - real :: prandtl ! flux ratio for diffusive convection regime + real :: prandtl ! flux ratio for diffusive convection regime [nondim] real, parameter :: Rrho0 = 1.9 ! limit for double-diffusive density ratio [nondim] @@ -1146,7 +1143,7 @@ subroutine add_drag_diffusivity(h, u, v, tv, fluxes, visc, j, TKE_to_Kd, & integer, intent(in) :: j !< j-index of row to work on real, dimension(SZI_(G),SZK_(GV)), intent(in) :: TKE_to_Kd !< The conversion rate between the TKE !! TKE dissipated within a layer and the - !! diapycnal diffusivity witin that layer, + !! diapycnal diffusivity within that layer, !! usually (~Rho_0 / (G_Earth * dRho_lay)) !! [Z2 T-1 / Z3 T-3 = T2 Z-1 ~> s2 m-1] real, dimension(SZI_(G),SZK_(GV)), intent(in) :: maxTKE !< The energy required to for a layer to entrain @@ -1274,7 +1271,7 @@ subroutine add_drag_diffusivity(h, u, v, tv, fluxes, visc, j, TKE_to_Kd, & ! if (maxEnt(i,k) <= 0.0) cycle if (maxTKE(i,k) <= 0.0) cycle - ! This is an analytic integral where diffusity is a quadratic function of + ! This is an analytic integral where diffusivity is a quadratic function of ! rho that goes asymptotically to 0 at Rho_top (vaguely following KPP?). if (TKE(i) > 0.0) then if (Rint(K) <= Rho_top(i)) then @@ -1395,7 +1392,7 @@ subroutine add_LOTW_BBL_diffusivity(h, u, v, tv, fluxes, visc, j, N2_int, Kd_int real :: ustar ! value of ustar at a thickness point [Z T-1 ~> m s-1]. real :: ustar2 ! square of ustar, for convenience [Z2 T-2 ~> m2 s-2] real :: absf ! average absolute value of Coriolis parameter around a thickness point [T-1 ~> s-1] - real :: dh, dhm1 ! thickness of layers k and k-1, respecitvely [Z ~> m]. + real :: dh, dhm1 ! thickness of layers k and k-1, respectively [Z ~> m]. real :: z_bot ! distance to interface k from bottom [Z ~> m]. real :: D_minus_z ! distance to interface k from surface [Z ~> m]. real :: total_thickness ! total thickness of water column [Z ~> m]. @@ -1550,7 +1547,7 @@ subroutine add_MLrad_diffusivity(h, fluxes, j, Kd_int, G, GV, US, CS, TKE_to_Kd, real :: h_ml_sq ! The square of the mixed layer thickness [Z2 ~> m2]. real :: ustar_sq ! ustar squared [Z2 T-2 ~> m2 s-2] real :: Kd_mlr ! A diffusivity associated with mixed layer turbulence radiation [Z2 T-1 ~> m2 s-1]. - real :: C1_6 ! 1/6 + real :: C1_6 ! 1/6 [nondim] real :: Omega2 ! rotation rate squared [T-2 ~> s-2]. real :: z1 ! layer thickness times I_decay [nondim] real :: dzL ! thickness converted to heights [Z ~> m]. @@ -1623,7 +1620,7 @@ subroutine add_MLrad_diffusivity(h, fluxes, j, Kd_int, G, GV, US, CS, TKE_to_Kd, do i=is,ie ; if (do_i(i)) then dzL = GV%H_to_Z*h(i,j,k) ; z1 = dzL*I_decay(i) if (CS%ML_Rad_bug) then - ! These expresssions are dimensionally inconsistent. -RWH + ! These expressions are dimensionally inconsistent. -RWH ! This is supposed to be the integrated energy deposited in the layer, ! not the average over the layer as in these expressions. if (z1 > 1e-5) then @@ -1881,8 +1878,8 @@ subroutine set_density_ratios(h, tv, kb, G, GV, US, CS, j, ds_dsp1, rho_0) ! Local variables real :: g_R0 ! g_R0 is a rescaled version of g/Rho [L2 Z-1 R-1 T-2 ~> m4 kg-1 s-2] - real :: eps, tmp ! nondimensional temporary variables - real :: a(SZK_(GV)), a_0(SZK_(GV)) ! nondimensional temporary variables + real :: eps, tmp ! nondimensional temporary variables [nondim] + real :: a(SZK_(GV)), a_0(SZK_(GV)) ! nondimensional temporary variables [nondim] real :: p_ref(SZI_(G)) ! an array of tv%P_Ref pressures [R L2 T-2 ~> Pa] real :: Rcv(SZI_(G),SZK_(GV)) ! coordinate density in the mixed and buffer layers [R ~> kg m-3] real :: I_Drho ! temporary variable [R-1 ~> m3 kg-1] @@ -1950,7 +1947,7 @@ subroutine set_density_ratios(h, tv, kb, G, GV, US, CS, j, ds_dsp1, rho_0) do k3=2,kmb ! ds_dsp1(i,k3) = MAX(a(k3),1e-5) - ! Deliberately treat convective instabilies of the upper mixed + ! Deliberately treat convective instabilities of the upper mixed ! and buffer layers with respect to the deepest buffer layer as ! though they don't exist. They will be eliminated by the upcoming ! call to the mixedlayer code anyway. @@ -1974,7 +1971,7 @@ subroutine set_diffusivity_init(Time, G, GV, US, param_file, diag, CS, int_tide_ type(diag_ctrl), target, intent(inout) :: diag !< A structure used to regulate diagnostic output. type(set_diffusivity_CS), pointer :: CS !< pointer set to point to the module control !! structure. - type(int_tide_CS), intent(in), target :: int_tide_CSp !< Internal tide control struct + type(int_tide_CS), intent(in), target :: int_tide_CSp !< Internal tide control structure integer, intent(out) :: halo_TS !< The halo size of tracer points that must be !! valid for the calculations in set_diffusivity. logical, intent(out) :: double_diffuse !< This indicates whether some version @@ -1986,7 +1983,7 @@ subroutine set_diffusivity_init(Time, G, GV, US, param_file, diag, CS, int_tide_ !! surface boundary layer. ! Local variables - real :: decay_length + real :: decay_length ! The maximum decay scale for the BBL diffusion [Z ~> m] logical :: ML_use_omega integer :: default_answer_date ! The default setting for the various ANSWER_DATE flags. logical :: default_2018_answers ! The default setting for the various 2018_ANSWERS flags. @@ -2176,7 +2173,7 @@ subroutine set_diffusivity_init(Time, G, GV, US, param_file, diag, CS, int_tide_ "for an isopycnal layer-formulation.", & default=.false., do_not_log=.not.TKE_to_Kd_used) - ! set params related to the background mixing + ! set parameters related to the background mixing call bkgnd_mixing_init(Time, G, GV, US, param_file, CS%diag, CS%bkgnd_mixing_csp, physical_OBL_scheme) call get_param(param_file, mdl, "KV", CS%Kv, & @@ -2340,7 +2337,7 @@ subroutine set_diffusivity_init(Time, G, GV, US, param_file, diag, CS, int_tide_ end subroutine set_diffusivity_init -!> Clear pointers and dealocate memory +!> Clear pointers and deallocate memory subroutine set_diffusivity_end(CS) type(set_diffusivity_CS), intent(inout) :: CS !< Control structure for this module diff --git a/src/parameterizations/vertical/MOM_set_viscosity.F90 b/src/parameterizations/vertical/MOM_set_viscosity.F90 index e9497d0e92..249cbb5777 100644 --- a/src/parameterizations/vertical/MOM_set_viscosity.F90 +++ b/src/parameterizations/vertical/MOM_set_viscosity.F90 @@ -48,10 +48,10 @@ module MOM_set_visc logical :: initialized = .false. !< True if this control structure has been initialized. real :: Hbbl !< The static bottom boundary layer thickness [H ~> m or kg m-2]. !! Runtime parameter `HBBL`. - real :: cdrag !< The quadratic drag coefficient. + real :: cdrag !< The quadratic drag coefficient [nondim]. !! Runtime parameter `CDRAG`. real :: c_Smag !< The Laplacian Smagorinsky coefficient for - !! calculating the drag in channels. + !! calculating the drag in channels [nondim]. real :: drag_bg_vel !< An assumed unresolved background velocity for !! calculating the bottom drag [L T-1 ~> m s-1]. !! Runtime parameter `DRAG_BG_VEL`. @@ -233,7 +233,7 @@ subroutine set_viscous_BBL(u, v, h, tv, visc, G, GV, US, CS, pbv) real :: C24_a ! 24/a [H-1 ~> m-1 or m2 kg-1]. real :: slope ! The absolute value of the bottom depth slope across ! a cell times the cell width [H ~> m or kg m-2]. - real :: apb_4a, ax2_3apb ! Various nondimensional ratios of a and slope. + real :: apb_4a, ax2_3apb ! Various nondimensional ratios of a and slope [nondim]. real :: a2x48_apb3, Iapb, Ibma_2 ! Combinations of a and slope [H-1 ~> m-1 or m2 kg-1]. ! All of the following "volumes" have units of thickness because they are normalized ! by the full horizontal area of a velocity cell. @@ -282,12 +282,12 @@ subroutine set_viscous_BBL(u, v, h, tv, visc, G, GV, US, CS, pbv) real :: h_bbl_fr ! The fraction of the bottom boundary layer in a layer [nondim]. real :: h_sum ! The sum of the thicknesses of the layers below the one being ! worked on [H ~> m or kg m-2]. - real, parameter :: C1_3 = 1.0/3.0, C1_6 = 1.0/6.0, C1_12 = 1.0/12.0 - real :: C2pi_3 ! An irrational constant, 2/3 pi. - real :: tmp ! A temporary variable. - real :: tmp_val_m1_to_p1 + real, parameter :: C1_3 = 1.0/3.0, C1_6 = 1.0/6.0, C1_12 = 1.0/12.0 ! Rational constants [nondim] + real :: C2pi_3 ! An irrational constant, 2/3 pi. [nondim] + real :: tmp ! A temporary variable, sometimes in [Z ~> m] + real :: tmp_val_m1_to_p1 ! A temporary variable [nondim] real :: curv_tol ! Numerator of curvature cubed, used to estimate - ! accuracy of a single L(:) Newton iteration + ! accuracy of a single L(:) Newton iteration [H5 ~> m3 or kg5 m-10] logical :: use_L0, do_one_L_iter ! Control flags for L(:) Newton iteration logical :: use_BBL_EOS, do_i(SZIB_(G)) integer, dimension(2) :: EOSdom ! The computational domain for the equation of state @@ -1099,7 +1099,7 @@ function set_v_at_u(v, h, G, GV, i, j, k, mask2dCv, OBC) integer, intent(in) :: j !< The j-index of the u-location to work on. integer, intent(in) :: k !< The k-index of the u-location to work on. real, dimension(SZI_(G),SZJB_(G)),& - intent(in) :: mask2dCv !< A multiplicative mask of the v-points + intent(in) :: mask2dCv !< A multiplicative mask of the v-points [nondim] type(ocean_OBC_type), pointer :: OBC !< A pointer to an open boundary condition structure real :: set_v_at_u !< The return value of v at u points points in the !! same units as u, i.e. [L T-1 ~> m s-1] or other units. @@ -1144,7 +1144,7 @@ function set_u_at_v(u, h, G, GV, i, j, k, mask2dCu, OBC) integer, intent(in) :: j !< The j-index of the u-location to work on. integer, intent(in) :: k !< The k-index of the u-location to work on. real, dimension(SZIB_(G),SZJ_(G)), & - intent(in) :: mask2dCu !< A multiplicative mask of the u-points + intent(in) :: mask2dCu !< A multiplicative mask of the u-points [nondim] type(ocean_OBC_type), pointer :: OBC !< A pointer to an open boundary condition structure real :: set_u_at_v !< The return value of u at v points in the !! same units as u, i.e. [L T-1 ~> m s-1] or other units. @@ -1213,8 +1213,10 @@ subroutine set_viscous_ML(u, v, h, tv, forces, visc, dt, G, GV, US, CS) Rhtot, & ! The integrated density of layers that are within the surface mixed layer ! [H R ~> kg m-2 or kg2 m-5]. Rhtot is only used if no ! equation of state is used. - uhtot, & ! The depth integrated zonal and meridional velocities within - vhtot, & ! the surface mixed layer [H L T-1 ~> m2 s-1 or kg m-1 s-1]. + uhtot, & ! The depth integrated zonal velocity within the surface + ! mixed layer [H L T-1 ~> m2 s-1 or kg m-1 s-1]. + vhtot, & ! The depth integrated meridional velocity within the surface + ! mixed layer [H L T-1 ~> m2 s-1 or kg m-1 s-1]. Idecay_len_TKE, & ! The inverse of a turbulence decay length scale [H-1 ~> m-1 or m2 kg-1]. dR_dT, & ! Partial derivative of the density at the base of layer nkml ! (roughly the base of the mixed layer) with temperature [R C-1 ~> kg m-3 degC-1]. @@ -1261,8 +1263,7 @@ subroutine set_viscous_ML(u, v, h, tv, forces, visc, dt, G, GV, US, CS) ! on the mixed layer thickness and density difference across ! the base of the mixed layer [L2 T-2 ~> m2 s-2]. real :: RiBulk ! The bulk Richardson number below which water is in the - ! viscous mixed layer, including reduction for turbulent - ! decay. Nondimensional. + ! viscous mixed layer, including reduction for turbulent decay [nondim] real :: dt_Rho0 ! The time step divided by the conversion from the layer ! thickness to layer mass [T H Z-1 R-1 ~> s m3 kg-1 or s]. real :: g_H_Rho0 ! The gravitational acceleration times the conversion from H to m divided @@ -1871,7 +1872,7 @@ subroutine set_visc_register_restarts(HI, GV, US, param_file, visc, restart_CS) logical :: adiabatic, useKPP, useEPBL logical :: use_CVMix_shear, MLE_use_PBL_MLD, use_CVMix_conv integer :: isd, ied, jsd, jed, nz - real :: hfreeze !< If hfreeze > 0 [m], melt potential will be computed. + real :: hfreeze !< If hfreeze > 0 [Z ~> m], melt potential will be computed. character(len=40) :: mdl = "MOM_set_visc" ! This module's name. isd = HI%isd ; ied = HI%ied ; jsd = HI%jsd ; jed = HI%jed ; nz = GV%ke @@ -1929,7 +1930,7 @@ subroutine set_visc_register_restarts(HI, GV, US, param_file, visc, restart_CS) default=.false., do_not_log=.true.) ! visc%MLD needs to be allocated when melt potential is computed (HFREEZE>0) call get_param(param_file, mdl, "HFREEZE", hfreeze, & - default=-1.0, do_not_log=.true.) + units="m", default=-1.0, scale=US%m_to_Z, do_not_log=.true.) if (MLE_use_PBL_MLD) then call safe_alloc_ptr(visc%MLD, isd, ied, jsd, jed) @@ -2001,11 +2002,11 @@ subroutine set_visc_init(Time, G, GV, US, param_file, diag, visc, CS, restart_CS real :: Chan_max_thick_dflt ! The default value for CHANNEL_DRAG_MAX_THICK [m] real :: Z_rescale ! A rescaling factor for heights from the representation in - ! a restart file to the internal representation in this run. + ! a restart file to the internal representation in this run [nondim]? real :: I_T_rescale ! A rescaling factor for time from the internal representation in this run - ! to the representation in a restart file. + ! to the representation in a restart file [nondim]? real :: Z2_T_rescale ! A rescaling factor for vertical diffusivities and viscosities from the - ! representation in a restart file to the internal representation in this run. + ! representation in a restart file to the internal representation in this run [nondim]? integer :: i, j, k, is, ie, js, je integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB, nz integer :: default_answer_date ! The default setting for the various ANSWER_DATE flags. @@ -2100,14 +2101,14 @@ subroutine set_visc_init(Time, G, GV, US, param_file, diag, visc, CS, restart_CS "determine the mixed layer thickness for viscosity.", & default=.false.) if (CS%dynamic_viscous_ML) then - call get_param(param_file, mdl, "BULK_RI_ML", bulk_Ri_ML_dflt, default=0.0) + call get_param(param_file, mdl, "BULK_RI_ML", bulk_Ri_ML_dflt, units="nondim", default=0.0) call get_param(param_file, mdl, "BULK_RI_ML_VISC", CS%bulk_Ri_ML, & "The efficiency with which mean kinetic energy released "//& "by mechanically forced entrainment of the mixed layer "//& "is converted to turbulent kinetic energy. By default, "//& "BULK_RI_ML_VISC = BULK_RI_ML or 0.", units="nondim", & default=bulk_Ri_ML_dflt) - call get_param(param_file, mdl, "TKE_DECAY", TKE_decay_dflt, default=0.0) + call get_param(param_file, mdl, "TKE_DECAY", TKE_decay_dflt, units="nondim", default=0.0) call get_param(param_file, mdl, "TKE_DECAY_VISC", CS%TKE_decay, & "TKE_DECAY_VISC relates the vertical rate of decay of "//& "the TKE available for mechanical entrainment to the "//& @@ -2118,7 +2119,7 @@ subroutine set_visc_init(Time, G, GV, US, param_file, diag, visc, CS, restart_CS call get_param(param_file, mdl, "ML_USE_OMEGA", use_omega, & "If true, use the absolute rotation rate instead of the "//& "vertical component of rotation when setting the decay "//& - "scale for turbulence.", default=.false., do_not_log=.true.) + "scale for turbulence.", default=.false., do_not_log=.true.) omega_frac_dflt = 0.0 if (use_omega) then call MOM_error(WARNING, "ML_USE_OMEGA is deprecated; use ML_OMEGA_FRAC=1.0 instead.") @@ -2223,7 +2224,7 @@ subroutine set_visc_init(Time, G, GV, US, param_file, diag, visc, CS, restart_CS default=.false.) if (CS%Channel_drag) then - call get_param(param_file, mdl, "SMAG_LAP_CONST", smag_const1, default=-1.0) + call get_param(param_file, mdl, "SMAG_LAP_CONST", smag_const1, units="nondim", default=-1.0) cSmag_chan_dflt = 0.15 if (smag_const1 >= 0.0) cSmag_chan_dflt = smag_const1