17 use cvmix_convection,
only : cvmix_init_conv, cvmix_coeffs_conv
18 use cvmix_kpp,
only : cvmix_kpp_compute_kobl_depth
20 implicit none ;
private
22 #include <MOM_memory.h>
40 integer :: id_n2 = -1, id_kd_conv = -1, id_kv_conv = -1
44 real,
allocatable,
dimension(:,:,:) :: n2
45 real,
allocatable,
dimension(:,:,:) :: kd_conv
46 real,
allocatable,
dimension(:,:,:) :: kv_conv
50 character(len=40) ::
mdl =
"MOM_CVMix_conv"
57 type(time_type),
intent(in) :: time
62 type(
diag_ctrl),
target,
intent(inout) :: diag
69 #include "version_variable.h"
71 if (
associated(cs))
then
72 call mom_error(warning,
"CVMix_conv_init called with an associated "// &
80 "Parameterization of enhanced mixing due to convection via CVMix")
82 "If true, turns on the enhanced mixing due to convection "//&
83 "via CVMix. This scheme increases diapycnal diffs./viscs. "//&
84 "at statically unstable interfaces. Relevant parameters are "//&
85 "contained in the CVMix_CONVECTION% parameter block.", &
90 call get_param(param_file,
mdl,
"ENERGETICS_SFC_PBL", useepbl, default=.false., &
96 call mom_error(warning,
'MOM_CVMix_conv_init: '// &
97 'CVMix convection may not be properly applied when ENERGETICS_SFC_PBL = True'//&
98 'as convective mixing might occur in the boundary layer.')
101 call get_param(param_file,
mdl,
'DEBUG', cs%debug, default=.false., do_not_log=.true.)
103 call get_param(param_file,
mdl,
'MIN_THICKNESS', cs%min_thickness, default=0.001, do_not_log=.true.)
105 call openparameterblock(param_file,
'CVMix_CONVECTION')
107 call get_param(param_file,
mdl,
"PRANDTL_CONV", prandtl_conv, &
108 "The turbulent Prandtl number applied to convective "//&
109 "instabilities (i.e., used to convert KD_CONV into KV_CONV)", &
110 units=
"nondim", default=1.0)
112 call get_param(param_file,
mdl,
'KD_CONV', cs%kd_conv_const, &
113 "Diffusivity used in convective regime. Corresponding viscosity "//&
114 "(KV_CONV) will be set to KD_CONV * PRANDTL_TURB.", &
115 units=
'm2/s', default=1.00)
117 call get_param(param_file,
mdl,
'BV_SQR_CONV', cs%bv_sqr_conv, &
118 "Threshold for squared buoyancy frequency needed to trigger "//&
119 "Brunt-Vaisala parameterization.", &
120 units=
'1/s^2', default=0.0)
122 call closeparameterblock(param_file)
125 cs%kv_conv_const = cs%kd_conv_const * prandtl_conv
128 allocate(cs%N2(szi_(g), szj_(g), szk_(g)+1)); cs%N2(:,:,:) = 0.
129 allocate(cs%kd_conv(szi_(g), szj_(g), szk_(g)+1)); cs%kd_conv(:,:,:) = 0.
130 allocate(cs%kv_conv(szi_(g), szj_(g), szk_(g)+1)); cs%kv_conv(:,:,:) = 0.
134 cs%id_N2 = register_diag_field(
'ocean_model',
'N2_conv', diag%axesTi, time, &
135 'Square of Brunt-Vaisala frequency used by MOM_CVMix_conv module',
'1/s2')
136 cs%id_kd_conv = register_diag_field(
'ocean_model',
'kd_conv', diag%axesTi, time, &
137 'Additional diffusivity added by MOM_CVMix_conv module',
'm2/s', conversion=us%Z2_T_to_m2_s)
138 cs%id_kv_conv = register_diag_field(
'ocean_model',
'kv_conv', diag%axesTi, time, &
139 'Additional viscosity added by MOM_CVMix_conv module',
'm2/s', conversion=us%Z2_T_to_m2_s)
141 call cvmix_init_conv(convect_diff=cs%kd_conv_const, &
142 convect_visc=cs%kv_conv_const, &
143 lbruntvaisala=.true., &
144 bvsqr_convect=cs%bv_sqr_conv)
155 real,
dimension(SZI_(G),SZJ_(G),SZK_(G)),
intent(in) :: h
159 real,
dimension(:,:),
optional,
pointer :: hbl
161 real,
dimension(SZK_(G)) :: rho_lwr
164 real,
dimension(SZK_(G)) :: rho_1d
166 real,
dimension(SZK_(G)+1) :: kv_col
167 real,
dimension(SZK_(G)+1) :: kd_col
168 real,
dimension(SZK_(G)+1) :: ifaceheight
169 real,
dimension(SZK_(G)) :: cellheight
172 real :: pref, rhok, rhokm1, dz, dh, hcorr
175 g_o_rho0 = gv%mks_g_Earth / (us%R_to_kg_m3*gv%Rho0)
178 rho_lwr(:) = 0.0; rho_1d(:) = 0.0
180 if (.not.
associated(hbl))
then
181 allocate(hbl(szi_(g), szj_(g)))
190 cs%N2(i,j,g%ke+1) = 0.
200 pref = pref + gv%H_to_Pa * h(i,j,k)
202 call calculate_density(tv%t(i,j,k-1), tv%s(i,j,k-1), pref, rhokm1, tv%eqn_of_state)
204 dz = ((0.5*(h(i,j,k-1) + h(i,j,k))+gv%H_subroundoff)*gv%H_to_m)
205 cs%N2(i,j,k) = g_o_rho0 * (rhok - rhokm1) / dz
213 dh = h(i,j,k) * gv%H_to_m
215 hcorr = min( dh - cs%min_thickness, 0. )
216 dh = max( dh, cs%min_thickness )
217 cellheight(k) = ifaceheight(k) - 0.5 * dh
218 ifaceheight(k+1) = ifaceheight(k) - dh
222 kobl = cvmix_kpp_compute_kobl_depth(ifaceheight, cellheight,hbl(i,j))
224 kv_col(:) = 0.0 ; kd_col(:) = 0.0
225 call cvmix_coeffs_conv(mdiff_out=kv_col(:), &
226 tdiff_out=kd_col(:), &
229 dens_lwr=rho_lwr(:), &
235 cs%kv_conv(i,j,k) = us%m2_s_to_Z2_T * kv_col(k)
236 cs%Kd_conv(i,j,k) = us%m2_s_to_Z2_T * kd_col(k)
240 cs%kv_conv(i,j,k) = 0.0
241 cs%kd_conv(i,j,k) = 0.0
248 call hchksum(cs%N2,
"MOM_CVMix_conv: N2",g%HI,haloshift=0)
249 call hchksum(cs%kd_conv,
"MOM_CVMix_conv: kd_conv",g%HI,haloshift=0,scale=us%Z2_T_to_m2_s)
250 call hchksum(cs%kv_conv,
"MOM_CVMix_conv: kv_conv",g%HI,haloshift=0,scale=us%m2_s_to_Z2_T)
254 if (cs%id_N2 > 0)
call post_data(cs%id_N2, cs%N2, cs%diag)
255 if (cs%id_kd_conv > 0)
call post_data(cs%id_kd_conv, cs%kd_conv, cs%diag)
256 if (cs%id_kv_conv > 0)
call post_data(cs%id_kv_conv, cs%kv_conv, cs%diag)
266 default=.false., do_not_log = .true.)
275 if (.not.
associated(cs))
return
278 deallocate(cs%kd_conv)
279 deallocate(cs%kv_conv)