147 subroutine cires_ugwpv1_init (me, master, nlunit, logunit, jdat_gfs, con_pi, &
148 con_rerth, fn_nml2, input_nml_file, lonr, latr, levs, ak, bk, &
149 pref, dtp, errmsg, errflg)
166 integer,
intent (in) :: me
167 integer,
intent (in) :: master
168 integer,
intent (in) :: nlunit
169 integer,
intent (in) :: logunit
170 integer,
intent (in) :: lonr
171 integer,
intent (in) :: levs
172 integer,
intent (in) :: latr
173 integer,
intent (in) :: jdat_gfs(8)
174 real(kind=kind_phys),
intent (in) :: ak(levs+1), bk(levs+1), pref
175 real(kind=kind_phys),
intent (in) :: dtp
179 real(kind=kind_phys),
intent (in) :: con_pi, con_rerth
181 character(len=64),
intent (in) :: fn_nml2
182 character(len=*),
intent (in) :: input_nml_file(:)
184 character(len=*),
intent(out) :: errmsg
185 integer,
intent(out) :: errflg
192 integer :: ncid, iernc, vid, dimid, status
194 integer :: ddd_ugwp, curday_ugwp
201#ifdef INTERNAL_FILE_NML
202 read (input_nml_file, nml = cires_ugwp_nml)
204 if (me == master) print *, trim(fn_nml2),
' GW-namelist file '
205 inquire (file =trim(fn_nml2) , exist = exists)
207 if (.not. exists)
then
208 write(errmsg,
'(3a)')
'cires_ugwpv1_init: namelist file: ', trim(fn_nml2),
' does not exist'
212 open (unit = nlunit, file = trim(fn_nml2), action =
'read', status =
'old', iostat = ios)
215 read (nlunit, nml = cires_ugwp_nml)
219 strsolver= knob_ugwp_orosolv
221 curday_ugwp = jdat_gfs(1)*10000 + jdat_gfs(2)*100 +jdat_gfs(3)
222 call calendar_ugwp(jdat_gfs(1), jdat_gfs(2), jdat_gfs(3), ddd_ugwp)
225 if (me == master)
then
226 write (logunit, *)
" ================================================================== "
227 write (logunit, *)
"CCPP cires_ugwp_namelist_extended_v1"
228 write (logunit, nml = cires_ugwp_nml)
229 write (logunit, *)
" ================================================================== "
231 write (6, *)
" ================================================================== "
232 write (6, *)
"CCPP cires_ugwp_namelist_extended_v1"
233 write (6, nml = cires_ugwp_nml)
234 write (6, *)
" ================================================================== "
235 write (6, *)
"calendar_ugwp ddd_ugwp=", ddd_ugwp
236 write (6, *)
"calendar_ugwp curday_ugwp=", curday_ugwp
237 write (6, *)
" ================================================================== "
238 write (6, *) ddd_ugwp,
' jdat_gfs ddd of year '
244 kxw = pi2/knob_ugwp_lhmet
252 allocate( kvg(levs+1), ktg(levs+1) )
253 allocate( krad(levs+1), kion(levs+1) )
254 allocate( zkm(levs), pmb(levs) )
261 pmb(k) = ak(k) + pref*bk(k)
262 zkm(k) = -hpskm*alog(pmb(k)/pref)
268 if (knob_ugwp_palaunch > 900.e2)
then
269 write(errmsg,
'(a,e16.7)')
'cires_ugwpv1_init: unrealistic value of knob_ugwp_palaunch', knob_ugwp_palaunch
275 if (pmb(k) .gt. knob_ugwp_palaunch )
exit
278 launch_level = max(k-1, 5)
279 if (me == master)
then
280 print *,
'cires_ugwpv1 klev_ngw =', launch_level, nint(pmb(launch_level))
285 call init_global_gwdis(levs, zkm, pmb, kvg, ktg, krad, kion, me, master)
294 call init_oro_gws( knob_ugwp_wvspec(1), knob_ugwp_azdir(1), &
295 knob_ugwp_stoch(1), knob_ugwp_effac(1), lonr, kxw )
299 do_physb_gwsrcs=.true.
301 IF (do_physb_gwsrcs)
THEN
304 if (knob_ugwp_wvspec(4) > 0)
then
306 call init_okw_gws(knob_ugwp_wvspec(4), knob_ugwp_azdir(4), &
307 knob_ugwp_stoch(4), knob_ugwp_effac(4), &
312 if (knob_ugwp_wvspec(3) > 0)
then
314 call init_fjet_gws(knob_ugwp_wvspec(3), knob_ugwp_azdir(3), &
315 knob_ugwp_stoch(3), knob_ugwp_effac(3), &
320 if (knob_ugwp_wvspec(2) > 0)
then
322 call init_conv_gws(knob_ugwp_wvspec(2), knob_ugwp_azdir(2), &
323 knob_ugwp_stoch(2), knob_ugwp_effac(2), &
338 if (knob_ugwp_solver==1)
then
341 call initsolv_lsatdis(me, master, knob_ugwp_wvspec(2), knob_ugwp_azdir(2), &
342 knob_ugwp_stoch(2), knob_ugwp_effac(2), do_physb_gwsrcs, kxw )
344 if (knob_ugwp_solver == 2)
then
348 nslope = knob_ugwp_nslope
349 lzstar = knob_ugwp_lzstar
350 lzmax = knob_ugwp_lzmax
351 lzmin = knob_ugwp_lzmin
352 lhmet = knob_ugwp_lhmet
353 tamp_mpa =knob_ugwp_tauamp
354 tau_min =knob_ugwp_taumin
355 ilaunch = launch_level
359 call initsolv_wmsdis(me, master, knob_ugwp_wvspec(2), knob_ugwp_azdir(2), &
360 knob_ugwp_stoch(2), knob_ugwp_effac(2), do_physb_gwsrcs, kxw, knob_ugwp_version)
366 module_is_initialized = .true.
367 if (me == master) print *,
' CIRES_ugwpV1 is initialized ', module_is_initialized
448 subroutine ngwflux_update(me, master, im, levs, kdt, ddd, curdate, &
449 tau_ddd, xlatd, sinlat,coslat, rain, tau_ngw)
451 use machine,
only: kind_phys
455 integer,
intent(in) :: me, master
456 integer,
intent(in) :: im, levs, kdt
457 integer,
intent(in) :: ddd, curdate
462 real(kind=kind_phys),
intent(in),
dimension(im) :: xlatd, sinlat,coslat
463 real(kind=kind_phys),
intent(in),
dimension(im) :: rain, tau_ddd
465 real(kind=kind_phys),
intent(inout),
dimension(im) :: tau_ngw
470 integer :: i, j1, j2, k, it1, it2, iday
471 real(kind=kind_phys) :: tem, tx1, tx2, w1, w2, wlat, rw1, rw2
472 real(kind=kind_phys) :: tau_rain, flat_rain, tau_3dt
502 tau_3dt = tau_ngw(i) * w_merra + w_nomerra *tau_ddd(i)
504 if (w_rain > 0. .and. rain(i) > 0.)
then
508 if (wlat <= rain_lat .and. rain(i) > rain_lim)
then
509 flat_rain = wlat/rain_lat
510 rw1 = 0.75 * flat_rain ; rw2 = 1.-rw1
512 tau_rain = tau_3dt * rw1 + rw2 * mtau_rain*min(rain_max, rain(i))/rain_lim
513 tau_rain = tau_3dt*(1.-w_rain) + w_rain* tau_rain
519 if (tau_rain < ft_min *tau_3dt) tau_rain = ft_min *tau_3dt
520 if (tau_rain > ft_max *tau_3dt) tau_rain = ft_max *tau_3dt
525 if (metoum_rain == 1)
then
526 tau_rain = min( sqrt(rain(i))*pbase_um, tau_ngw_max)
527 tau_3dt = max(tau_ngw_min, tau_rain)