CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_phys_time_vary.scm.F90
1
3
9
10 use machine, only : kind_phys, kind_dbl_prec, kind_sngl_prec
11
13
14 use module_ozphys, only: ty_ozphys
15
16 use h2o_def, only : levh2o, h2o_coeff, h2o_lat, h2o_pres, h2o_time, h2oplin
17 use h2ointerp, only : read_h2odata, setindxh2o, h2ointerpol
18
19 use aerclm_def, only : aerin, aer_pres, ntrcaer, ntrcaerm, iamin, iamax, jamin, jamax
20 use aerinterp, only : read_aerdata, setindxaer, aerinterpol, read_aerdataf
21
22 use iccn_def, only : ciplin, ccnin, ci_pres
23 use iccninterp, only : read_cidata, setindxci, ciinterpol
24
25 use cires_tauamf_data, only: cires_indx_ugwp, read_tau_amf, tau_amf_interp
26 use cires_tauamf_data, only: tau_limb, days_limb, ugwp_taulat
27
28 !--- variables needed for calculating 'sncovr'
29 use namelist_soilveg, only: salp_data, snupx
31
32 ! --- needed for Noah MP init
33 use noahmp_tables, only: read_mp_table_parameters, &
34 laim_table,saim_table,sla_table, &
35 bexp_table,smcmax_table,smcwlt_table, &
36 dwsat_table,dksat_table,psisat_table, &
37 isurban_table,isbarren_table, &
38 isice_table,iswater_table
39 implicit none
40
41 private
42
44
45 logical :: is_initialized = .false.
46
47 real(kind=kind_phys), parameter :: con_hr = 3600.0_kind_phys
48 real(kind=kind_phys), parameter :: con_99 = 99.0_kind_phys
49 real(kind=kind_phys), parameter :: con_100 = 100.0_kind_phys
50 real(kind=kind_phys), parameter :: missing_value = 9.99e20_kind_phys
51 real(kind=kind_phys), parameter :: drythresh = 1.e-4_kind_phys
52 real(kind=kind_phys), parameter :: zero = 0.0_kind_phys
53 real(kind=kind_phys), parameter :: one = 1.0_kind_phys
54
55 contains
56
63 me, master, ntoz, h2o_phys, iaerclm, iccn, iflip, im, nx, ny, idate, xlat_d, xlon_d, &
64 jindx1_o3, jindx2_o3, ddy_o3, jindx1_h, jindx2_h, ddy_h, h2opl,fhour, &
65 jindx1_aer, jindx2_aer, ddy_aer, iindx1_aer, iindx2_aer, ddx_aer, aer_nm, &
66 jindx1_ci, jindx2_ci, ddy_ci, iindx1_ci, iindx2_ci, ddx_ci, imap, jmap, &
67 do_ugwp_v1, jindx1_tau, jindx2_tau, ddy_j1tau, ddy_j2tau, &
68 isot, ivegsrc, nlunit, sncovr, sncovr_ice, lsm, lsm_noahmp, lsm_ruc, min_seaice, &
69 fice, landfrac, vtype, weasd, lsoil, zs, dzs, lsnow_lsm_lbound, lsnow_lsm_ubound, &
70 tvxy, tgxy, tahxy, canicexy, canliqxy, eahxy, cmxy, chxy, fwetxy, sneqvoxy, alboldxy,&
71 qsnowxy, wslakexy, albdvis_lnd, albdnir_lnd, albivis_lnd, albinir_lnd, albdvis_ice, &
72 albdnir_ice, albivis_ice, albinir_ice, emiss_lnd, emiss_ice, taussxy, waxy, wtxy, &
73 zwtxy, xlaixy, xsaixy, lfmassxy, stmassxy, rtmassxy, woodxy, stblcpxy, fastcpxy, &
74 smcwtdxy, deeprechxy, rechxy, snowxy, snicexy, snliqxy, tsnoxy , smoiseq, zsnsoxy, &
75 slc, smc, stc, tsfcl, snowd, canopy, tg3, stype, con_t0c, lsm_cold_start, nthrds, &
76 ozphys, errmsg, errflg)
77
78 implicit none
79
80 ! Interface variables
81 integer, intent(in) :: me, master, ntoz, iccn, iflip, im, nx, ny
82 logical, intent(in) :: h2o_phys, iaerclm, lsm_cold_start
83 integer, intent(in) :: idate(:)
84 real(kind_phys), intent(in) :: fhour
85 real(kind_phys), intent(in) :: xlat_d(:), xlon_d(:)
86
87 integer, intent(inout), optional :: jindx1_o3(:), jindx2_o3(:), jindx1_h(:), jindx2_h(:)
88 real(kind_phys), intent(inout), optional :: ddy_o3(:), ddy_h(:)
89 real(kind_phys), intent(in) :: h2opl(:,:,:)
90 integer, intent(inout), optional :: jindx1_aer(:), jindx2_aer(:), iindx1_aer(:), iindx2_aer(:)
91 real(kind_phys), intent(inout), optional :: ddy_aer(:), ddx_aer(:)
92 real(kind_phys), intent(in) :: aer_nm(:,:,:)
93 integer, intent(inout), optional :: jindx1_ci(:), jindx2_ci(:), iindx1_ci(:), iindx2_ci(:)
94 real(kind_phys), intent(inout), optional :: ddy_ci(:), ddx_ci(:)
95 integer, intent(inout) :: imap(:), jmap(:)
96 logical, intent(in) :: do_ugwp_v1
97 real(kind_phys), intent(inout), optional :: ddy_j1tau(:), ddy_j2tau(:)
98 integer, intent(inout), optional :: jindx1_tau(:), jindx2_tau(:)
99
100 integer, intent(in) :: isot, ivegsrc, nlunit
101 real(kind_phys), intent(inout) :: sncovr(:), sncovr_ice(:)
102 integer, intent(in) :: lsm, lsm_noahmp, lsm_ruc, vtype(:)
103 real(kind_phys), intent(in) :: min_seaice, fice(:)
104 real(kind_phys), intent(in) :: landfrac(:)
105 real(kind_phys), intent(inout) :: weasd(:)
106 type(ty_ozphys), intent(in) :: ozphys
107
108 ! NoahMP - only allocated when NoahMP is used
109 integer, intent(in) :: lsoil, lsnow_lsm_lbound, lsnow_lsm_ubound
110 real(kind_phys), intent(in) :: zs(:)
111 real(kind_phys), intent(in) :: dzs(:)
112 real(kind_phys), intent(inout), optional :: tvxy(:)
113 real(kind_phys), intent(inout), optional :: tgxy(:)
114 real(kind_phys), intent(inout), optional :: tahxy(:)
115 real(kind_phys), intent(inout), optional :: canicexy(:)
116 real(kind_phys), intent(inout), optional :: canliqxy(:)
117 real(kind_phys), intent(inout), optional :: eahxy(:)
118 real(kind_phys), intent(inout), optional :: cmxy(:)
119 real(kind_phys), intent(inout), optional :: chxy(:)
120 real(kind_phys), intent(inout), optional :: fwetxy(:)
121 real(kind_phys), intent(inout), optional :: sneqvoxy(:)
122 real(kind_phys), intent(inout), optional :: alboldxy(:)
123 real(kind_phys), intent(inout), optional :: qsnowxy(:)
124 real(kind_phys), intent(inout), optional :: wslakexy(:)
125 real(kind_phys), intent(inout) :: albdvis_lnd(:)
126 real(kind_phys), intent(inout) :: albdnir_lnd(:)
127 real(kind_phys), intent(inout) :: albivis_lnd(:)
128 real(kind_phys), intent(inout) :: albinir_lnd(:)
129 real(kind_phys), intent(inout), optional :: albdvis_ice(:)
130 real(kind_phys), intent(inout), optional :: albdnir_ice(:)
131 real(kind_phys), intent(inout), optional :: albivis_ice(:)
132 real(kind_phys), intent(inout), optional :: albinir_ice(:)
133 real(kind_phys), intent(inout) :: emiss_lnd(:)
134 real(kind_phys), intent(inout) :: emiss_ice(:)
135 real(kind_phys), intent(inout), optional :: taussxy(:)
136 real(kind_phys), intent(inout), optional :: waxy(:)
137 real(kind_phys), intent(inout), optional :: wtxy(:)
138 real(kind_phys), intent(inout), optional :: zwtxy(:)
139 real(kind_phys), intent(inout), optional :: xlaixy(:)
140 real(kind_phys), intent(inout), optional :: xsaixy(:)
141 real(kind_phys), intent(inout), optional :: lfmassxy(:)
142 real(kind_phys), intent(inout), optional :: stmassxy(:)
143 real(kind_phys), intent(inout), optional :: rtmassxy(:)
144 real(kind_phys), intent(inout), optional :: woodxy(:)
145 real(kind_phys), intent(inout), optional :: stblcpxy(:)
146 real(kind_phys), intent(inout), optional :: fastcpxy(:)
147 real(kind_phys), intent(inout), optional :: smcwtdxy(:)
148 real(kind_phys), intent(inout), optional :: deeprechxy(:)
149 real(kind_phys), intent(inout), optional :: rechxy(:)
150 real(kind_phys), intent(inout), optional :: snowxy(:)
151 real(kind_phys), intent(inout), optional :: snicexy(:,lsnow_lsm_lbound:)
152 real(kind_phys), intent(inout), optional :: snliqxy(:,lsnow_lsm_lbound:)
153 real(kind_phys), intent(inout), optional :: tsnoxy (:,lsnow_lsm_lbound:)
154 real(kind_phys), intent(inout), optional :: smoiseq(:,:)
155 real(kind_phys), intent(inout), optional :: zsnsoxy(:,lsnow_lsm_lbound:)
156 real(kind_phys), intent(inout) :: slc(:,:)
157 real(kind_phys), intent(inout) :: smc(:,:)
158 real(kind_phys), intent(inout) :: stc(:,:)
159 real(kind_phys), intent(in) :: tsfcl(:)
160 real(kind_phys), intent(in) :: snowd(:)
161 real(kind_phys), intent(in) :: canopy(:)
162 real(kind_phys), intent(in) :: tg3(:)
163 integer, intent(in) :: stype(:)
164 real(kind_phys), intent(in) :: con_t0c
165
166 integer, intent(in) :: nthrds
167 character(len=*), intent(out) :: errmsg
168 integer, intent(out) :: errflg
169
170 ! Local variables
171 integer :: i, j, ix, vegtyp
172 real(kind_phys) :: rsnow
173
174 !--- Noah MP
175 integer :: soiltyp, isnow, is, imn
176 real(kind=kind_phys) :: masslai, masssai, snd
177 real(kind=kind_phys) :: bexp, ddz, smcmax, smcwlt, dwsat, dksat, psisat
178
179 real(kind=kind_phys), dimension(:), allocatable :: dzsno
180 real(kind=kind_phys), dimension(:), allocatable :: dzsnso
181
182 ! Initialize CCPP error handling variables
183 errmsg = ''
184 errflg = 0
185
186 if (is_initialized) return
187 iamin=999
188 iamax=-999
189 jamin=999
190 jamax=-999
191
193 call read_h2odata (h2o_phys, me, master)
194
195 ! Consistency check that the hardcoded values for levh2o and
196 ! h2o_coeff in GFS_typedefs.F90 match what is set by read_h2odata
197 ! in GFS_typedefs.F90: allocate (Tbd%h2opl (IM,levh2o,h2o_coeff))
198 if (size(h2opl, dim=2).ne.levh2o) then
199 write(errmsg,'(2a,i0,a,i0)') "Value error in GFS_phys_time_vary_init: ", &
200 "levh2o from read_h2odata does not match value in GFS_typedefs.F90: ", &
201 levh2o, " /= ", size(h2opl, dim=2)
202 errflg = 1
203 end if
204 if (size(h2opl, dim=3).ne.h2o_coeff) then
205 write(errmsg,'(2a,i0,a,i0)') "Value error in GFS_phys_time_vary_init: ", &
206 "h2o_coeff from read_h2odata does not match value in GFS_typedefs.F90: ", &
207 h2o_coeff, " /= ", size(h2opl, dim=3)
208 errflg = 1
209 end if
210
212 if (iaerclm) then
213 ! Consistency check that the value for ntrcaerm set in GFS_typedefs.F90
214 ! and used to allocate aer_nm matches the value defined in aerclm_def
215 if (size(aer_nm, dim=3).ne.ntrcaerm) then
216 write(errmsg,'(2a,i0,a,i0)') "Value error in GFS_phys_time_vary_init: ", &
217 "ntrcaerm from aerclm_def does not match value in GFS_typedefs.F90: ", &
218 ntrcaerm, " /= ", size(aer_nm, dim=3)
219 errflg = 1
220 else
221 ! Update the value of ntrcaer in aerclm_def with the value defined
222 ! in GFS_typedefs.F90 that is used to allocate the Tbd DDT.
223 ! If iaerclm is .true., then ntrcaer == ntrcaerm
224 ntrcaer = size(aer_nm, dim=3)
225 ! Read aerosol climatology
226 call read_aerdata (me,master,iflip,idate,errmsg,errflg)
227 endif
228 else
229 ! Update the value of ntrcaer in aerclm_def with the value defined
230 ! in GFS_typedefs.F90 that is used to allocate the Tbd DDT.
231 ! If iaerclm is .false., then ntrcaer == 1
232 ntrcaer = size(aer_nm, dim=3)
233 endif
234
236 if (iccn == 1) then
237 call read_cidata (me,master)
238 ! No consistency check needed for in/ccn data, all values are
239 ! hardcoded in module iccn_def.F and GFS_typedefs.F90
240 endif
241
243 if (do_ugwp_v1) then
244 call read_tau_amf(me, master, errmsg, errflg)
245 endif
246
248 call set_soilveg(me, isot, ivegsrc, nlunit, errmsg, errflg)
249
251 if(lsm == lsm_noahmp) then
252 call read_mp_table_parameters(errmsg, errflg)
253 endif
254
256 if (ntoz > 0) then
257 call ozphys%setup_o3prog(xlat_d, jindx1_o3, jindx2_o3, ddy_o3)
258 endif
259
261 if (h2o_phys) then
262 call setindxh2o (im, xlat_d, jindx1_h, jindx2_h, ddy_h)
263 endif
264
266 if (iaerclm) then
267 call setindxaer (im, xlat_d, jindx1_aer, &
268 jindx2_aer, ddy_aer, xlon_d, &
269 iindx1_aer, iindx2_aer, ddx_aer, &
270 me, master)
271 iamin = min(minval(iindx1_aer), iamin)
272 iamax = max(maxval(iindx2_aer), iamax)
273 jamin = min(minval(jindx1_aer), jamin)
274 jamax = max(maxval(jindx2_aer), jamax)
275 endif
276
278 if (iccn == 1) then
279 call setindxci (im, xlat_d, jindx1_ci, &
280 jindx2_ci, ddy_ci, xlon_d, &
281 iindx1_ci, iindx2_ci, ddx_ci)
282 endif
283
285 if (do_ugwp_v1) then
286 call cires_indx_ugwp (im, me, master, xlat_d, jindx1_tau, jindx2_tau, &
287 ddy_j1tau, ddy_j2tau)
288 endif
289
290 !--- initial calculation of maps local ix -> global i and j
291 ix = 0
292 do j = 1,ny
293 do i = 1,nx
294 ix = ix + 1
295 jmap(ix) = j
296 imap(ix) = i
297 enddo
298 enddo
299
300 !--- if sncovr does not exist in the restart, need to create it
301 if (all(sncovr < zero)) then
302 if (me == master ) write(*,'(a)') 'GFS_phys_time_vary_init: compute sncovr from weasd and soil vegetation parameters'
303 !--- compute sncovr from existing variables
304 !--- code taken directly from read_fix.f
305 sncovr(:) = zero
306 do ix=1,im
307 if (landfrac(ix) >= drythresh .or. fice(ix) >= min_seaice) then
308 vegtyp = vtype(ix)
309 if (vegtyp == 0) vegtyp = 7
310 rsnow = 0.001_kind_phys*weasd(ix)/snupx(vegtyp)
311 if (0.001_kind_phys*weasd(ix) < snupx(vegtyp)) then
312 sncovr(ix) = one - (exp(-salp_data*rsnow) - rsnow*exp(-salp_data))
313 else
314 sncovr(ix) = one
315 endif
316 endif
317 enddo
318 endif
319
320 !--- For RUC LSM: create sncovr_ice from sncovr
321 if (lsm == lsm_ruc) then
322 if (all(sncovr_ice < zero)) then
323 if (me == master ) write(*,'(a)') 'GFS_phys_time_vary_init: fill sncovr_ice with sncovr for RUC LSM'
324 sncovr_ice(:) = sncovr(:)
325 endif
326 endif
327
328 if (errflg/=0) return
329
330 if (iaerclm) then
331 call read_aerdataf (me, master, iflip, idate, fhour, errmsg, errflg)
332 if (errflg/=0) return
333 end if
334
335 !--- For Noah MP or RUC LSMs: initialize four components of albedo for
336 !--- land and ice - not for restart runs
337 lsm_init: if (lsm_cold_start) then
338 if (lsm == lsm_noahmp .or. lsm == lsm_ruc) then
339 if (me == master ) write(*,'(a)') 'GFS_phys_time_vary_init: initialize albedo for land and ice'
340 do ix=1,im
341 albdvis_lnd(ix) = 0.2_kind_phys
342 albdnir_lnd(ix) = 0.2_kind_phys
343 albivis_lnd(ix) = 0.2_kind_phys
344 albinir_lnd(ix) = 0.2_kind_phys
345 emiss_lnd(ix) = 0.95_kind_phys
346 enddo
347 endif
348 if (lsm == lsm_ruc) then
349 do ix=1,im
350 albdvis_ice(ix) = 0.6_kind_phys
351 albdnir_ice(ix) = 0.6_kind_phys
352 albivis_ice(ix) = 0.6_kind_phys
353 albinir_ice(ix) = 0.6_kind_phys
354 emiss_ice(ix) = 0.97_kind_phys
355 enddo
356 endif
357
358 noahmp_init: if (lsm == lsm_noahmp) then
359 allocate(dzsno(lsnow_lsm_lbound:lsnow_lsm_ubound))
360 allocate(dzsnso(lsnow_lsm_lbound:lsoil) )
361 dzsno(:) = missing_value
362 dzsnso(:) = missing_value
363
364 tvxy(:) = missing_value
365 tgxy(:) = missing_value
366 tahxy(:) = missing_value
367 canicexy(:) = missing_value
368 canliqxy(:) = missing_value
369 eahxy(:) = missing_value
370 cmxy(:) = missing_value
371 chxy(:) = missing_value
372 fwetxy(:) = missing_value
373 sneqvoxy(:) = missing_value
374 alboldxy(:) = missing_value
375 qsnowxy(:) = missing_value
376 wslakexy(:) = missing_value
377 taussxy(:) = missing_value
378 waxy(:) = missing_value
379 wtxy(:) = missing_value
380 zwtxy(:) = missing_value
381 xlaixy(:) = missing_value
382 xsaixy(:) = missing_value
383
384 lfmassxy(:) = missing_value
385 stmassxy(:) = missing_value
386 rtmassxy(:) = missing_value
387 woodxy(:) = missing_value
388 stblcpxy(:) = missing_value
389 fastcpxy(:) = missing_value
390 smcwtdxy(:) = missing_value
391 deeprechxy(:) = missing_value
392 rechxy(:) = missing_value
393
394 snowxy(:) = missing_value
395 snicexy(:,:) = missing_value
396 snliqxy(:,:) = missing_value
397 tsnoxy(:,:) = missing_value
398 smoiseq(:,:) = missing_value
399 zsnsoxy(:,:) = missing_value
400
401 imn = idate(2)
402
403 do ix=1,im
404 if (landfrac(ix) >= drythresh) then
405 tvxy(ix) = tsfcl(ix)
406 tgxy(ix) = tsfcl(ix)
407 tahxy(ix) = tsfcl(ix)
408
409 if (snowd(ix) > 0.01_kind_phys .and. tsfcl(ix) > con_t0c ) then
410 tvxy(ix) = con_t0c
411 tgxy(ix) = con_t0c
412 tahxy(ix) = con_t0c
413 end if
414
415 canicexy(ix) = 0.0_kind_phys
416 canliqxy(ix) = canopy(ix)
417
418 eahxy(ix) = 2000.0_kind_phys
419
420 cmxy(ix) = zero
421 chxy(ix) = zero
422 fwetxy(ix) = zero
423 sneqvoxy(ix) = weasd(ix) ! mm
424 alboldxy(ix) = 0.65_kind_phys
425 qsnowxy(ix) = zero
426
427! if (srflag(ix) > 0.001) qsnowxy(ix) = tprcp(ix)/dtp
428 ! already set to 0.0
429 wslakexy(ix) = zero
430 taussxy(ix) = zero
431
432 waxy(ix) = 4900.0_kind_phys
433 wtxy(ix) = waxy(ix)
434 zwtxy(ix) = (25.0_kind_phys + 2.0_kind_phys) - waxy(ix) / 1000.0_kind_phys / 0.2_kind_phys
435
436 vegtyp = vtype(ix)
437 if (vegtyp == 0) vegtyp = 7
438
439 if ((vegtyp == isbarren_table) .or. (vegtyp == isice_table) .or. (vegtyp == isurban_table) .or. (vegtyp == iswater_table)) then
440
441 xlaixy(ix) = zero
442 xsaixy(ix) = zero
443
444 lfmassxy(ix) = zero
445 stmassxy(ix) = zero
446 rtmassxy(ix) = zero
447
448 woodxy(ix) = zero
449 stblcpxy(ix) = zero
450 fastcpxy(ix) = zero
451
452 else
453
454 xlaixy(ix) = max(laim_table(vegtyp, imn),0.05_kind_phys)
455! xsaixy(ix) = max(saim_table(vegtyp, imn),0.05)
456 xsaixy(ix) = max(xlaixy(ix)*0.1_kind_phys,0.05_kind_phys)
457
458 masslai = 1000.0_kind_phys / max(sla_table(vegtyp),one)
459 lfmassxy(ix) = xlaixy(ix)*masslai
460 masssai = 1000.0_kind_phys / 3.0_kind_phys
461 stmassxy(ix) = xsaixy(ix)* masssai
462
463 rtmassxy(ix) = 500.0_kind_phys
464
465 woodxy(ix) = 500.0_kind_phys
466 stblcpxy(ix) = 1000.0_kind_phys
467 fastcpxy(ix) = 1000.0_kind_phys
468
469 endif ! non urban ...
470
471 if (vegtyp == isice_table) then
472 do is = 1,lsoil
473 stc(ix,is) = min(stc(ix,is),min(tg3(ix),263.15_kind_phys))
474 smc(ix,is) = one
475 slc(ix,is) = zero
476 enddo
477 endif
478
479 snd = snowd(ix)/1000.0_kind_phys ! go to m from snwdph
480
481 if (weasd(ix) /= zero .and. snd == zero ) then
482 snd = weasd(ix)/1000.0
483 endif
484
485 if (vegtyp == 15) then ! land ice in MODIS/IGBP
486 if (weasd(ix) < 0.1_kind_phys) then
487 weasd(ix) = 0.1_kind_phys
488 snd = 0.01_kind_phys
489 endif
490 endif
491
492 if (snd < 0.025_kind_phys ) then
493 snowxy(ix) = zero
494 dzsno(-2:0) = zero
495 elseif (snd >= 0.025_kind_phys .and. snd <= 0.05_kind_phys ) then
496 snowxy(ix) = -1.0_kind_phys
497 dzsno(0) = snd
498 elseif (snd > 0.05_kind_phys .and. snd <= 0.10_kind_phys ) then
499 snowxy(ix) = -2.0_kind_phys
500 dzsno(-1) = 0.5_kind_phys*snd
501 dzsno(0) = 0.5_kind_phys*snd
502 elseif (snd > 0.10_kind_phys .and. snd <= 0.25_kind_phys ) then
503 snowxy(ix) = -2.0_kind_phys
504 dzsno(-1) = 0.05_kind_phys
505 dzsno(0) = snd - 0.05_kind_phys
506 elseif (snd > 0.25_kind_phys .and. snd <= 0.45_kind_phys ) then
507 snowxy(ix) = -3.0_kind_phys
508 dzsno(-2) = 0.05_kind_phys
509 dzsno(-1) = 0.5_kind_phys*(snd-0.05_kind_phys)
510 dzsno(0) = 0.5_kind_phys*(snd-0.05_kind_phys)
511 elseif (snd > 0.45_kind_phys) then
512 snowxy(ix) = -3.0_kind_phys
513 dzsno(-2) = 0.05_kind_phys
514 dzsno(-1) = 0.20_kind_phys
515 dzsno(0) = snd - 0.05_kind_phys - 0.20_kind_phys
516 else
517 errmsg = 'Error in GFS_phys_time_vary.scm.F90: Problem with the logic assigning snow layers in Noah MP initialization'
518 errflg = 1
519 endif
520
521! Now we have the snowxy field
522! snice + snliq + tsno allocation and compute them from what we have
523
524 tsnoxy(ix,:) = zero
525 snicexy(ix,:) = zero
526 snliqxy(ix,:) = zero
527 zsnsoxy(ix,:) = zero
528
529 isnow = nint(snowxy(ix))+1 ! snowxy <=0.0, dzsno >= 0.0
530
531 do is = isnow,0
532 tsnoxy(ix,is) = tgxy(ix)
533 snliqxy(ix,is) = zero
534 snicexy(ix,is) = one * dzsno(is) * weasd(ix)/snd
535 enddo
536!
537!zsnsoxy, all negative ?
538!
539 do is = isnow,0
540 dzsnso(is) = -dzsno(is)
541 enddo
542
543 do is = 1,4
544 dzsnso(is) = -dzs(is)
545 enddo
546!
547! Assign to zsnsoxy
548!
549 zsnsoxy(ix,isnow) = dzsnso(isnow)
550 do is = isnow+1,4
551 zsnsoxy(ix,is) = zsnsoxy(ix,is-1) + dzsnso(is)
552 enddo
553!
554! smoiseq
555! Init water table related quantities here
556!
557 soiltyp = stype(ix)
558 if (soiltyp /= 0) then
559 bexp = bexp_table(soiltyp)
560 smcmax = smcmax_table(soiltyp)
561 smcwlt = smcwlt_table(soiltyp)
562 dwsat = dwsat_table(soiltyp)
563 dksat = dksat_table(soiltyp)
564 psisat = -psisat_table(soiltyp)
565 endif
566
567 if (vegtyp == isurban_table) then
568 smcmax = 0.45_kind_phys
569 smcwlt = 0.40_kind_phys
570 endif
571
572 if ((bexp > zero) .and. (smcmax > zero) .and. (-psisat > zero)) then
573 do is = 1, lsoil
574 if ( is == 1 )then
575 ddz = -zs(is+1) * 0.5_kind_phys
576 elseif ( is < lsoil ) then
577 ddz = ( zs(is-1) - zs(is+1) ) * 0.5_kind_phys
578 else
579 ddz = zs(is-1) - zs(is)
580 endif
581 smoiseq(ix,is) = min(max(find_eq_smc(bexp, dwsat, dksat, ddz, smcmax),1.e-4_kind_phys),smcmax*0.99_kind_phys)
582 enddo
583 else ! bexp <= 0.0
584 smoiseq(ix,1:4) = smcmax
585 endif ! end the bexp condition
586
587 smcwtdxy(ix) = smcmax
588 deeprechxy(ix) = zero
589 rechxy(ix) = zero
590
591 endif
592
593 enddo ! ix
594
595 if (errflg/=0) return
596
597 deallocate(dzsno)
598 deallocate(dzsnso)
599
600 endif noahmp_init
601 endif lsm_init
602
603 is_initialized = .true.
604
605 contains
606
607!
608! Use newton-raphson method to find eq soil moisture
609!
610 function find_eq_smc(bexp, dwsat, dksat, ddz, smcmax) result(smc)
611 implicit none
612 real(kind=kind_phys), intent(in) :: bexp, dwsat, dksat, ddz, smcmax
613 real(kind=kind_phys) :: smc
614 real(kind=kind_phys) :: expon, aa, bb, func, dfunc, dx
615 integer :: iter
616 !
617 expon = bexp + 1.
618 aa = dwsat / ddz
619 bb = dksat / smcmax ** expon
620 smc = 0.5 * smcmax
621 !
622 do iter = 1,100
623 func = (smc - smcmax) * aa + bb * smc ** expon
624 dfunc = aa + bb * expon * smc ** bexp
625 dx = func / dfunc
626 smc = smc - dx
627 if ( abs(dx) < 1.e-6_kind_phys) return
628 enddo
629 end function find_eq_smc
630
631 end subroutine gfs_phys_time_vary_init
632!! @}
633
640 me, master, cnx, cny, isc, jsc, nrcm, im, levs, kdt, idate, nsswr, fhswr, lsswr, fhour, &
641 imfdeepcnv, cal_pre, random_clds, ozphys, ntoz, h2o_phys, iaerclm, iccn, clstp, &
642 jindx1_o3, jindx2_o3, ddy_o3, ozpl, jindx1_h, jindx2_h, ddy_h, h2opl, iflip, &
643 jindx1_aer, jindx2_aer, ddy_aer, iindx1_aer, iindx2_aer, ddx_aer, aer_nm, &
644 jindx1_ci, jindx2_ci, ddy_ci, iindx1_ci, iindx2_ci, ddx_ci, in_nm, ccn_nm, &
645 imap, jmap, prsl, seed0, rann, do_ugwp_v1, jindx1_tau, jindx2_tau, ddy_j1tau, ddy_j2tau,&
646 tau_amf, nthrds, errmsg, errflg)
647
648 implicit none
649
650 ! Interface variables
651 integer, intent(in) :: me, master, cnx, cny, isc, jsc, nrcm, im, levs, kdt, &
652 nsswr, imfdeepcnv, iccn, ntoz, iflip
653 integer, intent(in) :: idate(:)
654 real(kind_phys), intent(in) :: fhswr, fhour
655 logical, intent(in) :: lsswr, cal_pre, random_clds, h2o_phys, iaerclm
656 real(kind_phys), intent(out) :: clstp
657 integer, intent(in), optional :: jindx1_o3(:), jindx2_o3(:), jindx1_h(:), jindx2_h(:)
658 real(kind_phys), intent(in), optional :: ddy_o3(:), ddy_h(:)
659 real(kind_phys), intent(inout) :: ozpl(:,:,:), h2opl(:,:,:)
660 integer, intent(in), optional :: jindx1_aer(:), jindx2_aer(:), iindx1_aer(:), iindx2_aer(:)
661 real(kind_phys), intent(in), optional :: ddy_aer(:), ddx_aer(:)
662 real(kind_phys), intent(inout) :: aer_nm(:,:,:)
663 integer, intent(in), optional :: jindx1_ci(:), jindx2_ci(:), iindx1_ci(:), iindx2_ci(:)
664 real(kind_phys), intent(in), optional :: ddy_ci(:), ddx_ci(:)
665 real(kind_phys), intent(inout) :: in_nm(:,:), ccn_nm(:,:)
666 integer, intent(in) :: imap(:), jmap(:)
667 real(kind_phys), intent(in) :: prsl(:,:)
668 integer, intent(in) :: seed0
669 real(kind_phys), intent(inout) :: rann(:,:)
670
671 logical, intent(in) :: do_ugwp_v1
672 integer, intent(in), optional :: jindx1_tau(:), jindx2_tau(:)
673 real(kind_phys), intent(in), optional :: ddy_j1tau(:), ddy_j2tau(:)
674 real(kind_phys), intent(inout) :: tau_amf(:)
675 type(ty_ozphys), intent(in) :: ozphys
676 integer, intent(in) :: nthrds
677 character(len=*), intent(out) :: errmsg
678 integer, intent(out) :: errflg
679
680 ! Local variables
681 integer :: i, j, k, iseed, iskip, ix, idat(8), jdat(8), iday, j1, j2, nc, n1, n2, jdow, &
682 jdoy, jday, w3kindreal, w3kindint
683 real(kind_phys) :: wrk(1), tem, tx1, tx2, rjday
684 real(kind_phys) :: rannie(cny)
685 real(kind_phys) :: rndval(cnx*cny*nrcm)
686 real(kind_dbl_prec) :: rinc(5)
687 real(kind_sngl_prec) :: rinc4(5)
688
689 ! Initialize CCPP error handling variables
690 errmsg = ''
691 errflg = 0
692
693 ! Check initialization status
694 if (.not.is_initialized) then
695 write(errmsg,'(*(a))') "Logic error: GFS_phys_time_vary_timestep_init called before GFS_phys_time_vary_init"
696 errflg = 1
697 return
698 end if
699
700 !--- switch for saving convective clouds - cnvc90.f
701 !--- aka Ken Campana/Yu-Tai Hou legacy
702 if ((mod(kdt,nsswr) == 0) .and. (lsswr)) then
703 !--- initialize,accumulate,convert
704 clstp = 1100 + min(fhswr/con_hr,fhour,con_99)
705 elseif (mod(kdt,nsswr) == 0) then
706 !--- accumulate,convert
707 clstp = 0100 + min(fhswr/con_hr,fhour,con_99)
708 elseif (lsswr) then
709 !--- initialize,accumulate
710 clstp = 1100
711 else
712 !--- accumulate
713 clstp = 0100
714 endif
715
716 !--- random number needed for RAS and old SAS and when cal_pre=.true.
717 ! imfdeepcnv < 0 when ras = .true.
718 if ( (imfdeepcnv <= 0 .or. cal_pre) .and. random_clds ) then
719
720 iseed = mod(con_100*sqrt(fhour*con_hr),1.0d9) + seed0
721 call random_setseed(iseed)
722 call random_number(wrk)
723 do i = 1,cnx*nrcm
724 iseed = iseed + nint(wrk(1)*1000.0) * i
725 call random_setseed(iseed)
726 call random_number(rannie)
727 rndval(1+(i-1)*cny:i*cny) = rannie(1:cny)
728 enddo
729
730 do k = 1,nrcm
731 iskip = (k-1)*cnx*cny
732 do ix=1,im
733 j = jmap(ix)
734 i = imap(ix)
735 rann(ix,k) = rndval(i+isc-1 + (j+jsc-2)*cnx + iskip)
736 enddo
737 enddo
738
739 endif ! imfdeepcnv, cal_re, random_clds
740
742 idat=0
743 idat(1)=idate(4)
744 idat(2)=idate(2)
745 idat(3)=idate(3)
746 idat(5)=idate(1)
747 rinc=0.
748 rinc(2)=fhour
749 call w3kind(w3kindreal,w3kindint)
750 if(w3kindreal==4) then
751 rinc4=rinc
752 CALL w3movdat(rinc4,idat,jdat)
753 else
754 CALL w3movdat(rinc,idat,jdat)
755 endif
756 jdow = 0
757 jdoy = 0
758 jday = 0
759 call w3doxdat(jdat,jdow,jdoy,jday)
760 rjday = jdoy + jdat(5) / 24.
761 if (rjday < ozphys%time(1)) rjday = rjday + 365.
762
763 n2 = ozphys%ntime + 1
764 do j=2,ozphys%ntime
765 if (rjday < ozphys%time(j)) then
766 n2 = j
767 exit
768 endif
769 enddo
770 n1 = n2 - 1
771 if (n2 > ozphys%ntime) n2 = n2 - ozphys%ntime
772
774 if (ntoz > 0) then
775 call ozphys%update_o3prog(jindx1_o3, jindx2_o3, ddy_o3, rjday, n1, n2, ozpl)
776 endif
777
779 if (h2o_phys) then
780 call h2ointerpol (me, im, idate, fhour, &
781 jindx1_h, jindx2_h, &
782 h2opl, ddy_h)
783 endif
784
786 if (iccn == 1) then
787 call ciinterpol (me, im, idate, fhour, &
788 jindx1_ci, jindx2_ci, &
789 ddy_ci, iindx1_ci, &
790 iindx2_ci, ddx_ci, &
791 levs, prsl, in_nm, ccn_nm)
792 endif
793
795 if (do_ugwp_v1) then
796 call tau_amf_interp(me, master, im, idate, fhour, &
797 jindx1_tau, jindx2_tau, &
798 ddy_j1tau, ddy_j2tau, tau_amf)
799 endif
800
802 if (iaerclm) then
803 ! aerinterpol is using threading inside, don't
804 ! move into OpenMP parallel section above
805 call aerinterpol (me, master, nthrds, im, idate, &
806 fhour, iflip, jindx1_aer, jindx2_aer, &
807 ddy_aer, iindx1_aer, &
808 iindx2_aer, ddx_aer, &
809 levs, prsl, aer_nm, errmsg, errflg)
810 if(errflg /= 0) then
811 return
812 endif
813 endif
814
815! Not needed for SCM:
817 ! if (nscyc > 0) then
818 ! if (mod(kdt,nscyc) == 1) THEN
819 ! call gcycle (me, nthrds, nx, ny, isc, jsc, nsst, tile_num, nlunit, &
820 ! input_nml_file, lsoil, lsoil_lsm, kice, idate, ialb, isot, ivegsrc, &
821 ! use_ufo, nst_anl, fhcyc, phour, landfrac, lakefrac, min_seaice, min_lakeice,&
822 ! frac_grid, smc, slc, stc, smois, sh2o, tslb, tiice, tg3, tref, tsfc, &
823 ! tsfco, tisfc, hice, fice, facsf, facwf, alvsf, alvwf, alnsf, alnwf, &
824 ! zorli, zorll, zorlo, weasd, slope, snoalb, canopy, vfrac, vtype, &
825 ! stype, shdmin, shdmax, snowd, cv, cvb, cvt, oro, oro_uf, &
826 ! xlat_d, xlon_d, slmsk, imap, jmap)
827 ! endif
828 ! endif
829
831!! @}
832
838 subroutine gfs_phys_time_vary_timestep_finalize (errmsg, errflg)
839
840 implicit none
841
842 ! Interface variables
843 character(len=*), intent(out) :: errmsg
844 integer, intent(out) :: errflg
845
846 ! Initialize CCPP error handling variables
847 errmsg = ''
848 errflg = 0
849
851!! @}
852
856 subroutine gfs_phys_time_vary_finalize(errmsg, errflg)
857
858 implicit none
859
860 ! Interface variables
861 character(len=*), intent(out) :: errmsg
862 integer, intent(out) :: errflg
863
864 ! Initialize CCPP error handling variables
865 errmsg = ''
866 errflg = 0
867
868 if (.not.is_initialized) return
869
870 ! Deallocate h2o arrays
871 if (allocated(h2o_lat) ) deallocate(h2o_lat)
872 if (allocated(h2o_pres)) deallocate(h2o_pres)
873 if (allocated(h2o_time)) deallocate(h2o_time)
874 if (allocated(h2oplin) ) deallocate(h2oplin)
875
876 ! Deallocate aerosol arrays
877 if (allocated(aerin) ) deallocate(aerin)
878 if (allocated(aer_pres)) deallocate(aer_pres)
879
880 ! Deallocate IN and CCN arrays
881 if (allocated(ciplin) ) deallocate(ciplin)
882 if (allocated(ccnin) ) deallocate(ccnin)
883 if (allocated(ci_pres) ) deallocate(ci_pres)
884
885 ! Deallocate UGWP-input arrays
886 if (allocated(ugwp_taulat)) deallocate(ugwp_taulat)
887 if (allocated(tau_limb )) deallocate(tau_limb)
888 if (allocated(days_limb )) deallocate(days_limb)
889
890 is_initialized = .false.
891
892 end subroutine gfs_phys_time_vary_finalize
893
894 end module gfs_phys_time_vary
subroutine, public set_soilveg(me, isot, ivet, nlunit, errmsg, errflg)
This subroutine initializes soil and vegetation.
Definition set_soilveg.f:17
real(kind=kind_phys), parameter drythresh
real(kind=kind_phys), parameter one
real(kind=kind_phys), parameter con_99
subroutine, public gfs_phys_time_vary_init(me, master, ntoz, h2o_phys, iaerclm, iccn, iflip, im, nx, ny, idate, xlat_d, xlon_d, jindx1_o3, jindx2_o3, ddy_o3, jindx1_h, jindx2_h, ddy_h, h2opl, fhour, jindx1_aer, jindx2_aer, ddy_aer, iindx1_aer, iindx2_aer, ddx_aer, aer_nm, jindx1_ci, jindx2_ci, ddy_ci, iindx1_ci, iindx2_ci, ddx_ci, imap, jmap, do_ugwp_v1, jindx1_tau, jindx2_tau, ddy_j1tau, ddy_j2tau, isot, ivegsrc, nlunit, sncovr, sncovr_ice, lsm, lsm_noahmp, lsm_ruc, min_seaice, fice, landfrac, vtype, weasd, lsoil, zs, dzs, lsnow_lsm_lbound, lsnow_lsm_ubound, tvxy, tgxy, tahxy, canicexy, canliqxy, eahxy, cmxy, chxy, fwetxy, sneqvoxy, alboldxy, qsnowxy, wslakexy, albdvis_lnd, albdnir_lnd, albivis_lnd, albinir_lnd, albdvis_ice, albdnir_ice, albivis_ice, albinir_ice, emiss_lnd, emiss_ice, taussxy, waxy, wtxy, zwtxy, xlaixy, xsaixy, lfmassxy, stmassxy, rtmassxy, woodxy, stblcpxy, fastcpxy, smcwtdxy, deeprechxy, rechxy, snowxy, snicexy, snliqxy, tsnoxy, smoiseq, zsnsoxy, slc, smc, stc, tsfcl, snowd, canopy, tg3, stype, con_t0c, lsm_cold_start, nthrds, ozphys, errmsg, errflg)
subroutine, public gfs_phys_time_vary_finalize(errmsg, errflg)
real(kind=kind_phys), parameter missing_value
real(kind=kind_phys), parameter con_hr
real(kind=kind_phys), parameter con_100
subroutine, public gfs_phys_time_vary_timestep_finalize(errmsg, errflg)
real(kind=kind_phys), parameter zero
real(kind=kind_phys) function find_eq_smc(bexp, dwsat, dksat, ddz, smcmax)
subroutine, public gfs_phys_time_vary_timestep_init(me, master, cnx, cny, isc, jsc, nrcm, im, levs, kdt, idate, nsswr, fhswr, lsswr, fhour, imfdeepcnv, cal_pre, random_clds, ozphys, ntoz, h2o_phys, iaerclm, iccn, clstp, jindx1_o3, jindx2_o3, ddy_o3, ozpl, jindx1_h, jindx2_h, ddy_h, h2opl, iflip, jindx1_aer, jindx2_aer, ddy_aer, iindx1_aer, iindx2_aer, ddx_aer, aer_nm, jindx1_ci, jindx2_ci, ddy_ci, iindx1_ci, iindx2_ci, ddx_ci, in_nm, ccn_nm, imap, jmap, prsl, seed0, rann, do_ugwp_v1, jindx1_tau, jindx2_tau, ddy_j1tau, ddy_j2tau, tau_amf, nthrds, errmsg, errflg)
This module defines arrays in H2O scheme.
Definition h2o_def.f:6
This module contains subroutines of reading and interpolating h2o coefficients.
Definition h2ointerp.f90:8
This module defines IN and CCN arrays.
Definition iccn_def.F:6
This module contains subroutines of reading and interplating IN and CCN data.
Definition iccninterp.F90:8
The operational GFS currently parameterizes ozone production and destruction based on monthly mean co...
Data from MPTABLE.TBL, SOILPARM.TBL, GENPARM.TBL for NoahMP.
This module contains set_soilveg subroutine.
Definition set_soilveg.f:4
Derived type containing data and procedures needed by ozone photochemistry parameterization Note All ...