14 use mom_time_manager,
only : time_type,
operator(+),
operator(/), time_type_to_real
18 implicit none ;
private
20 #include <MOM_memory.h>
35 logical :: usewindstress
36 logical :: useheatflux
37 logical :: useevaporation
38 logical :: usediurnalsw
48 #include "version_variable.h"
50 character(len=40) ::
mdl =
"SCM_CVMix_tests"
56 real,
dimension(NIMEM_,NJMEM_, NKMEM_),
intent(out) :: t
57 real,
dimension(NIMEM_,NJMEM_, NKMEM_),
intent(out) :: s
58 real,
dimension(NIMEM_,NJMEM_, NKMEM_),
intent(in) :: h
63 logical,
optional,
intent(in) :: just_read_params
66 real :: upperlayertempmld
67 real :: upperlayersaltmld
68 real :: upperlayertemp
69 real :: upperlayersalt
70 real :: lowerlayertemp
71 real :: lowerlayersalt
72 real :: lowerlayerdtdz
73 real :: lowerlayerdsdz
74 real :: lowerlayermintemp
75 real :: zc, dz, top, bottom
77 integer :: i, j, k, is, ie, js, je, isd, ied, jsd, jed, nz
79 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec ; nz = g%ke
80 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
83 just_read = .false. ;
if (
present(just_read_params)) just_read = just_read_params
86 call get_param(param_file,
mdl,
"SCM_TEMP_MLD", upperlayertempmld, &
87 'Initial temp mixed layer depth', &
88 units=
'm', default=0.0, scale=us%m_to_Z, do_not_log=just_read)
89 call get_param(param_file,
mdl,
"SCM_SALT_MLD", upperlayersaltmld, &
90 'Initial salt mixed layer depth', &
91 units=
'm', default=0.0, scale=us%m_to_Z, do_not_log=just_read)
92 call get_param(param_file,
mdl,
"SCM_L1_SALT", upperlayersalt, &
93 'Layer 2 surface salinity', units=
'1e-3', default=35.0, do_not_log=just_read)
94 call get_param(param_file,
mdl,
"SCM_L1_TEMP", upperlayertemp, &
95 'Layer 1 surface temperature', units=
'C', default=20.0, do_not_log=just_read)
96 call get_param(param_file,
mdl,
"SCM_L2_SALT", lowerlayersalt, &
97 'Layer 2 surface salinity', units=
'1e-3', default=35.0, do_not_log=just_read)
98 call get_param(param_file,
mdl,
"SCM_L2_TEMP", lowerlayertemp, &
99 'Layer 2 surface temperature', units=
'C', default=20.0, do_not_log=just_read)
100 call get_param(param_file,
mdl,
"SCM_L2_DTDZ", lowerlayerdtdz, &
101 'Initial temperature stratification in layer 2', &
102 units=
'C/m', default=0.0, scale=us%Z_to_m, do_not_log=just_read)
103 call get_param(param_file,
mdl,
"SCM_L2_DSDZ", lowerlayerdsdz, &
104 'Initial salinity stratification in layer 2', &
105 units=
'PPT/m', default=0.0, scale=us%Z_to_m, do_not_log=just_read)
106 call get_param(param_file,
mdl,
"SCM_L2_MINTEMP",lowerlayermintemp, &
107 'Layer 2 minimum temperature', units=
'C', default=4.0, do_not_log=just_read)
109 if (just_read)
return
111 do j=js,je ;
do i=is,ie
115 bottom = bottom - h(i,j,k)*gv%H_to_Z
116 zc = 0.5*( top + bottom )
117 dz = min(0., zc + upperlayertempmld)
118 t(i,j,k) = max(lowerlayermintemp,lowerlayertemp + lowerlayerdtdz * dz)
119 dz = min(0., zc + upperlayersaltmld)
120 s(i,j,k) = lowerlayersalt + lowerlayerdsdz * dz
129 type(time_type),
intent(in) :: time
136 # include "version_variable.h"
141 if (
associated(cs))
then
142 call mom_error(fatal,
"SCM_CVMix_tests_surface_forcing_init called with an associated "// &
143 "control structure.")
150 call get_param(param_file,
mdl,
"SCM_USE_WIND_STRESS", &
151 cs%UseWindStress,
"Wind Stress switch "// &
152 "used in the SCM CVMix surface forcing.", &
153 units=
'', default=.false.)
155 cs%UseHeatFlux,
"Heat flux switch "// &
156 "used in the SCM CVMix test surface forcing.", &
157 units=
'', default=.false.)
158 call get_param(param_file,
mdl,
"SCM_USE_EVAPORATION", &
159 cs%UseEvaporation,
"Evaporation switch "// &
160 "used in the SCM CVMix test surface forcing.", &
161 units=
'', default=.false.)
163 cs%UseDiurnalSW,
"Diurnal sw radation switch "// &
164 "used in the SCM CVMix test surface forcing.", &
165 units=
'', default=.false.)
166 if (cs%UseWindStress)
then
168 cs%tau_x,
"Constant X-dir wind stress "// &
169 "used in the SCM CVMix test surface forcing.", &
170 units=
'N/m2', scale=us%kg_m3_to_R*us%m_s_to_L_T**2*us%L_to_Z, fail_if_missing=.true.)
172 cs%tau_y,
"Constant y-dir wind stress "// &
173 "used in the SCM CVMix test surface forcing.", &
174 units=
'N/m2', scale=us%kg_m3_to_R*us%m_s_to_L_T**2*us%L_to_Z, fail_if_missing=.true.)
176 if (cs%UseHeatFlux)
then
178 cs%surf_HF,
"Constant surface heat flux "// &
179 "used in the SCM CVMix test surface forcing.", &
180 units=
'm K/s', fail_if_missing=.true.)
182 if (cs%UseEvaporation)
then
184 cs%surf_evap,
"Constant surface evaporation "// &
185 "used in the SCM CVMix test surface forcing.", &
186 units=
'm/s', fail_if_missing=.true.)
188 if (cs%UseDiurnalSW)
then
190 cs%Max_sw,
"Maximum diurnal sw radiation "// &
191 "used in the SCM CVMix test surface forcing.", &
192 units=
'm K/s', fail_if_missing=.true.)
200 type(time_type),
intent(in) :: day
205 integer :: i, j, is, ie, js, je, isq, ieq, jsq, jeq
206 integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
209 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
210 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
211 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
212 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
214 do j=js,je ;
do i=isq,ieq
215 forces%taux(i,j) = cs%tau_x
217 do j=jsq,jeq ;
do i=is,ie
218 forces%tauy(i,j) = cs%tau_y
222 mag_tau = sqrt(cs%tau_x*cs%tau_x + cs%tau_y*cs%tau_y)
223 if (
associated(forces%ustar))
then ;
do j=js,je ;
do i=is,ie
224 forces%ustar(i,j) = sqrt( us%L_to_Z * mag_tau / (us%kg_m3_to_R*cs%Rho0) )
225 enddo ;
enddo ;
endif
232 type(
forcing),
intent(inout) :: fluxes
233 type(time_type),
intent(in) :: day
239 integer :: i, j, is, ie, js, je, isq, ieq, jsq, jeq
240 integer :: isd, ied, jsd, jed, isdb, iedb, jsdb, jedb
246 is = g%isc ; ie = g%iec ; js = g%jsc ; je = g%jec
247 isq = g%IscB ; ieq = g%IecB ; jsq = g%JscB ; jeq = g%JecB
248 isd = g%isd ; ied = g%ied ; jsd = g%jsd ; jed = g%jed
249 isdb = g%IsdB ; iedb = g%IedB ; jsdb = g%JsdB ; jedb = g%JedB
251 if (cs%UseHeatFlux)
then
255 do j=jsq,jeq ;
do i=is,ie
256 fluxes%sens(i,j) = cs%surf_HF * cs%Rho0 * fluxes%C_p
260 if (cs%UseEvaporation)
then
261 do j=jsq,jeq ;
do i=is,ie
265 fluxes%evap(i,j) = cs%surf_evap * us%kg_m3_to_R*us%m_to_Z*us%T_to_s * cs%Rho0
269 if (cs%UseDiurnalSW)
then
270 do j=jsq,jeq ;
do i=is,ie
275 fluxes%sw(i,j) = cs%Max_sw * max(0.0,cos(2*pi* &
276 (time_type_to_real(day)/86400.-0.5))) * cs%RHO0 * fluxes%C_p