25 subroutine gfs_surface_composites_post_run ( &
26 im, kice, km, rd, rvrdm1, cplflx, cplwav2atm, frac_grid, flag_cice, thsfc_loc, islmsk, dry, wet, icy, wind, t1, q1, prsl1, &
27 landfrac, lakefrac, oceanfrac, zorl, zorlo, zorll, zorli, garea, frac_ice, &
28 cd, cd_wat, cd_lnd, cd_ice, cdq, cdq_wat, cdq_lnd, cdq_ice, rb, rb_wat, rb_lnd, rb_ice, stress, stress_wat, stress_lnd, &
29 stress_ice, ffmm, ffmm_wat, ffmm_lnd, ffmm_ice, ffhh, ffhh_wat, ffhh_lnd, ffhh_ice, uustar, uustar_wat, uustar_lnd, &
30 uustar_ice, fm10, fm10_wat, fm10_lnd, fm10_ice, fh2, fh2_wat, fh2_lnd, fh2_ice, tsurf_wat, tsurf_lnd, tsurf_ice, &
31 cmm, cmm_wat, cmm_lnd, cmm_ice, chh, chh_wat, chh_lnd, chh_ice, gflx, gflx_wat, gflx_lnd, gflx_ice, ep1d, ep1d_wat, &
32 ep1d_lnd, ep1d_ice, weasd, weasd_lnd, weasd_ice, snowd, snowd_lnd, snowd_ice, tprcp, tprcp_wat, &
33 tprcp_lnd, tprcp_ice, evap, evap_wat, evap_lnd, evap_ice, hflx, hflx_wat, hflx_lnd, hflx_ice, qss, qss_wat, qss_lnd, &
34 qss_ice, tsfc, tsfco, tsfcl, tsfc_wat, tisfc, hice, cice, tiice, &
35 sigmaf, zvfun, lheatstrg, h0facu, h0facs, hflxq, hffac, stc, lkm, iopt_lake, iopt_lake_clm, use_lake_model, &
36 grav, prsik1, prslk1, prslki, z1, ztmax_wat, ztmax_lnd, ztmax_ice, huge, errmsg, errflg)
40 integer,
intent(in) :: im, kice, km, lkm, iopt_lake, iopt_lake_clm
41 logical,
intent(in) :: cplflx, frac_grid, cplwav2atm, frac_ice
42 logical,
intent(in) :: lheatstrg
43 logical,
dimension(:),
intent(in) :: flag_cice, dry, icy
44 logical,
dimension(:),
intent(in) :: wet
45 integer,
dimension(:),
intent(in) :: islmsk, use_lake_model
46 real(kind=kind_phys),
dimension(:),
intent(in) :: wind, t1, q1, prsl1, landfrac, lakefrac, oceanfrac, &
47 cd_wat, cd_lnd, cd_ice, cdq_wat, cdq_lnd, cdq_ice, rb_wat, rb_lnd, rb_ice, stress_wat, &
48 stress_lnd, stress_ice, ffmm_wat, ffmm_lnd, ffmm_ice, ffhh_wat, ffhh_lnd, ffhh_ice, uustar_wat, uustar_lnd, uustar_ice, &
49 fm10_wat, fm10_lnd, fm10_ice, fh2_wat, fh2_lnd, fh2_ice, tsurf_wat, tsurf_lnd, tsurf_ice, cmm_wat, cmm_lnd, cmm_ice, &
50 chh_wat, chh_lnd, chh_ice, gflx_wat, gflx_lnd, gflx_ice, ep1d_wat, ep1d_lnd, ep1d_ice, weasd_lnd, weasd_ice, &
51 snowd_lnd, snowd_ice, tprcp_wat, tprcp_lnd, tprcp_ice, evap_wat, evap_lnd, evap_ice, hflx_wat, hflx_lnd, &
52 hflx_ice, qss_wat, qss_lnd, qss_ice, tsfc_wat, zorlo, zorll, zorli, garea
54 real(kind=kind_phys),
dimension(:),
intent(inout) :: zorl, cd, cdq, rb, stress, ffmm, ffhh, uustar, fm10, &
55 fh2, cmm, chh, gflx, ep1d, weasd, snowd, tprcp, evap, hflx, qss, tsfc, tsfco, tsfcl, tisfc
57 real(kind=kind_phys),
dimension(:),
intent(inout) :: hice, cice
58 real(kind=kind_phys),
dimension(:),
intent(inout) :: sigmaf, zvfun, hflxq, hffac
59 real(kind=kind_phys),
intent(in ) :: h0facu, h0facs
61 real(kind=kind_phys),
intent(in ) :: rd, rvrdm1, huge
63 real(kind=kind_phys),
dimension(:,:),
intent(in ) :: tiice
64 real(kind=kind_phys),
dimension(:,:),
intent(inout) :: stc
67 logical,
intent(in ) :: thsfc_loc
68 real(kind=kind_phys),
intent(in ) :: grav
69 real(kind=kind_phys),
dimension(:),
intent(in ) :: prsik1, prslk1, prslki, z1
70 real(kind=kind_phys),
dimension(:),
intent(in ) :: ztmax_wat, ztmax_lnd, ztmax_ice
72 character(len=*),
intent(out) :: errmsg
73 integer,
intent(out) :: errflg
77 real(kind=kind_phys) :: txl, txi, txo, wfrac, q0, rho
79 real(kind=kind_phys) :: tsurf, virtfac, tv1, thv1, tvs, z0max, ztmax
80 real(kind=kind_phys) :: lnzorll, lnzorli, lnzorlo
82 real(kind=kind_phys) :: tem1, tem2, gdx
83 real(kind=kind_phys),
parameter :: z0lo=0.1, z0up=1.0
92 fractional_grid:
if (frac_grid)
then
100 txo = max(zero, wfrac-txi)
103 ep1d(i) = txl*ep1d_lnd(i) + txi*ep1d_ice(i) + txo*ep1d_wat(i)
104 weasd(i) = txl*weasd_lnd(i) + txi*weasd_ice(i)
105 snowd(i) = txl*snowd_lnd(i) + txi*snowd_ice(i)
108 sigmaf(i) = txl*sigmaf(i)
110 evap(i) = txl*evap_lnd(i) + txi*evap_ice(i) + txo*evap_wat(i)
111 hflx(i) = txl*hflx_lnd(i) + txi*hflx_ice(i) + txo*hflx_wat(i)
112 qss(i) = txl*qss_lnd(i) + txi*qss_ice(i) + txo*qss_wat(i)
113 gflx(i) = txl*gflx_lnd(i) + txi*gflx_ice(i) + txo*gflx_wat(i)
121 tsfc(i) = ( txl * cdq_lnd(i) * tsfcl(i) &
122 + txi * cdq_ice(i) * tisfc(i) &
123 + txo * cdq_wat(i) * tsfc_wat(i)) &
124 / (txl * cdq_lnd(i) + txi * cdq_ice(i) + txo * cdq_wat(i) )
125 tsurf = ( txl * cdq_lnd(i) * tsurf_lnd(i) &
126 + txi * cdq_ice(i) * tsurf_ice(i) &
127 + txo * cdq_wat(i) * tsurf_wat(i)) &
128 / (txl * cdq_lnd(i) + txi * cdq_ice(i) + txo * cdq_wat(i) )
130 q0 = max( q1(i), qmin )
131 virtfac = one + rvrdm1 * q0
133 tv1 = t1(i) * virtfac
135 thv1 = t1(i) * prslki(i) * virtfac
136 tvs = half * (tsfc(i)+tsurf) * virtfac
138 thv1 = t1(i) / prslk1(i) * virtfac
139 tvs = half * (tsfc(i)+tsurf)/prsik1(i) * virtfac
142 lnzorll = zero ; lnzorli = zero ; lnzorlo = zero
143 if (zorll(i) /= huge)
then
144 lnzorll = log(zorll(i))
146 if (zorli(i) /= huge)
then
147 lnzorli = log(zorli(i))
149 if (zorlo(i) /= huge)
then
150 lnzorlo = log(zorlo(i))
152 zorl(i) = exp(txl*lnzorll + txi*lnzorli + txo*lnzorlo)
154 z0max = 0.01_kind_phys * zorl(i)
155 ztmax = exp(txl*log(ztmax_lnd(i)) + txi*log(ztmax_ice(i)) + txo*log(ztmax_wat(i)))
160 ffmm(i) = ffmm_lnd(i)
161 ffhh(i) = ffhh_lnd(i)
162 fm10(i) = fm10_lnd(i)
166 stress(i) = stress_lnd(i)
167 uustar(i) = uustar_lnd(i)
168 elseif(txo == one)
then
170 ffmm(i) = ffmm_wat(i)
171 ffhh(i) = ffhh_wat(i)
172 fm10(i) = fm10_wat(i)
176 stress(i) = stress_wat(i)
177 uustar(i) = uustar_wat(i)
178 elseif(txi == one)
then
180 ffmm(i) = ffmm_ice(i)
181 ffhh(i) = ffhh_ice(i)
182 fm10(i) = fm10_ice(i)
186 stress(i) = stress_ice(i)
187 uustar(i) = uustar_ice(i)
192 tem1 = (z0max - z0lo) / (z0up - z0lo)
193 tem1 = min(max(tem1, zero), one)
194 tem2 = max(sigmaf(i), 0.1)
195 zvfun(i) = sqrt(tem1 * tem2)
204 if(hflx(i) > 0.)
then
205 hffac(i) = h0facu * zvfun(i)
207 hffac(i) = h0facs * zvfun(i)
209 hffac(i) = 1. + hffac(i)
210 hflxq(i) = hflx(i) / hffac(i)
213 call stability(z1(i), zvfun(i), gdx, tv1, thv1, wind(i), &
214 z0max, ztmax, tvs, grav, thsfc_loc, &
215 rb(i), ffmm(i), ffhh(i), fm10(i), fh2(i), cd(i), cdq(i), &
216 stress(i), uustar(i))
220 rho = prsl1(i) / (rd*tv1)
221 cmm(i) = cd(i)*wind(i)
222 chh(i) = rho*cdq(i)*wind(i)
228 tsfcl(i) = tsfc_wat(i)
233 tsfco(i) = tsfc_wat(i)
242 tisfc(i) = tsfc_wat(i)
258 if (.not. icy(i))
then
268 if (use_lake_model(i)>0)
then
269 if(frac_ice .and. icy(i))
then
270 call composite_icy(iopt_lake==iopt_lake_clm)
271 call composite_wet_and_icy
275 else if (islmsk(i) == 1)
then
278 elseif (islmsk(i) == 0)
then
283 call composite_icy(.false.)
284 call composite_wet_and_icy
288 endif fractional_grid
294 subroutine composite_land
300 stress(i) = stress_lnd(i)
301 ffmm(i) = ffmm_lnd(i)
302 ffhh(i) = ffhh_lnd(i)
303 uustar(i) = uustar_lnd(i)
304 fm10(i) = fm10_lnd(i)
311 gflx(i) = gflx_lnd(i)
312 ep1d(i) = ep1d_lnd(i)
313 weasd(i) = weasd_lnd(i)
314 snowd(i) = snowd_lnd(i)
315 evap(i) = evap_lnd(i)
316 hflx(i) = hflx_lnd(i)
320 end subroutine composite_land
322 subroutine composite_wet
328 stress(i) = stress_wat(i)
329 ffmm(i) = ffmm_wat(i)
330 ffhh(i) = ffhh_wat(i)
331 uustar(i) = uustar_wat(i)
332 fm10(i) = fm10_wat(i)
334 tsfco(i) = tsfc_wat(i)
340 gflx(i) = gflx_wat(i)
341 ep1d(i) = ep1d_wat(i)
344 evap(i) = evap_wat(i)
345 hflx(i) = hflx_wat(i)
349 end subroutine composite_wet
351 subroutine composite_icy(is_clm)
353 logical,
intent(in) :: is_clm
358 ffmm(i) = ffmm_ice(i)
359 ffhh(i) = ffhh_ice(i)
360 uustar(i) = uustar_ice(i)
361 fm10(i) = fm10_ice(i)
363 stress(i) = stress_ice(i)
366 gflx(i) = gflx_ice(i)
367 ep1d(i) = ep1d_ice(i)
369 weasd(i) = weasd_ice(i)
370 snowd(i) = snowd_ice(i)
372 weasd(i) = weasd_ice(i) * cice(i)
373 snowd(i) = snowd_ice(i) * cice(i)
376 evap(i) = evap_ice(i)
377 hflx(i) = hflx_ice(i)
378 end subroutine composite_icy
380 subroutine composite_wet_and_icy
384 evap(i) = txi * evap_ice(i) + txo * evap_wat(i)
385 hflx(i) = txi * hflx_ice(i) + txo * hflx_wat(i)
386 tsfc(i) = txi * tisfc(i) + txo * tsfc_wat(i)
387 stress(i) = txi * stress_ice(i) + txo * stress_wat(i)
388 qss(i) = txi * qss_ice(i) + txo * qss_wat(i)
389 ep1d(i) = txi * ep1d_ice(i) + txo * ep1d_wat(i)
391 lnzorli = zero ; lnzorlo = zero
392 if (zorli(i) /= huge)
then
393 lnzorli = log(zorli(i))
395 if (zorlo(i) /= huge)
then
396 lnzorlo = log(zorlo(i))
398 zorl(i) = exp(txi*lnzorli + txo*lnzorlo)
402 tsfco(i) = tsfc_wat(i)
408 stc(i,k) = tiice(i,k)
410 end subroutine composite_wet_and_icy