CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
dcyc2t3.f
1
5
6!! This module contains the CCPP-compliant dcyc2t3 codes that fits
7!! radiative fluxes and heating rates from a coarse radiation
8!! calculation time interval into model's more frequent time steps.
9 module dcyc2t3
10
11 implicit none
12
13 private
14
15 public :: dcyc2t3_run
16
17 contains
18
19! ===================================================================== !
20! description: !
21! !
22! dcyc2t3 fits radiative fluxes and heating rates from a coarse !
23! radiation calc time interval into model's more frequent time steps.!
24! solar heating rates and fluxes are scaled by the ratio of cosine !
25! of zenith angle at the current time to the mean value used in !
26! radiation calc. surface downward lw flux is scaled by the ratio !
27! of current surface air temperature (temp**4) to the corresponding !
28! temperature saved during lw radiation calculation. upward lw flux !
29! at the surface is computed by current ground surface temperature. !
30! surface emissivity effect will be taken in other part of the model.!
31! !
32! usage: !
33! !
34! call dcyc2t3 !
35! inputs: !
36! ( solhr,slag,sdec,cdec,sinlat,coslat, !
37! xlon,coszen,tsfc_lnd,tsfc_ice,tsfc_wat, !
38! tf,tsflw,sfcemis_lnd,sfcemis_ice,sfcemis_wat, !
39! sfcdsw,sfcnsw,sfcdlw,sfculw,swh,swhc,hlw,hlwc, !
40! sfcnirbmu,sfcnirdfu,sfcvisbmu,sfcvisdfu, !
41! sfcnirbmd,sfcnirdfd,sfcvisbmd,sfcvisdfd, !
42! im, levs, deltim, fhswr, !
43! dry, icy, wet !
44! input/output: !
45! dtdt,dtdtnp, !
46! outputs: !
47! adjsfcdsw,adjsfcnsw,adjsfcdlw, !
48! adjsfculw_lnd,adjsfculw_ice,adjsfculw_wat,xmu,xcosz, !
49! adjnirbmu,adjnirdfu,adjvisbmu,adjvisdfu, !
50! adjdnnbmd,adjdnndfd,adjdnvbmd,adjdnvdfd) !
51! !
52! !
53!
54! inputs:
55! solhr - real, forecast time in 24-hour form (hr)
56! slag - real, equation of time in radians
57! sdec, cdec - real, sin and cos of the solar declination angle
58! sinlat(im), coslat(im): !
59! - real, sin and cos of latitude !
60! xlon (im) - real, longitude in radians !
61! coszen (im) - real, avg of cosz over daytime sw call interval !
62! tsfc_lnd (im) - real, bottom surface temperature over land (k) !
63! tsfc_ice (im) - real, bottom surface temperature over ice (k) !
64! tsfc_wat (im) - real, bottom surface temperature over ocean (k) !
65! tf (im) - real, surface air (layer 1) temperature (k) !
66! sfcemis_lnd(im) - real, surface emissivity (fraction) o. land (k) !
67! sfcemis_ice(im) - real, surface emissivity (fraction) o. ice (k) !
68! sfcemis_wat(im) - real, surface emissivity (fraction) o. ocean (k)!
69! tsflw (im) - real, sfc air (layer 1) temp in k saved in lw call !
70! sfcdsw (im) - real, total sky sfc downward sw flux ( w/m**2 ) !
71! sfcnsw (im) - real, total sky sfc net sw into ground (w/m**2) !
72! sfcdlw (im) - real, total sky sfc downward lw flux ( w/m**2 ) !
73! sfculw (im) - real, total sky sfc upward lw flux ( w/m**2 ) !
74! swh(im,levs) - real, total sky sw heating rates ( k/s ) !
75! swhc(im,levs) - real, clear sky sw heating rates ( k/s ) !
76! hlw(im,levs) - real, total sky lw heating rates ( k/s ) !
77! hlwc(im,levs) - real, clear sky lw heating rates ( k/s ) !
78! sfcnirbmu(im)- real, tot sky sfc nir-beam sw upward flux (w/m2) !
79! sfcnirdfu(im)- real, tot sky sfc nir-diff sw upward flux (w/m2) !
80! sfcvisbmu(im)- real, tot sky sfc uv+vis-beam sw upward flux (w/m2)!
81! sfcvisdfu(im)- real, tot sky sfc uv+vis-diff sw upward flux (w/m2)!
82! sfcnirbmd(im)- real, tot sky sfc nir-beam sw downward flux (w/m2) !
83! sfcnirdfd(im)- real, tot sky sfc nir-diff sw downward flux (w/m2) !
84! sfcvisbmd(im)- real, tot sky sfc uv+vis-beam sw dnward flux (w/m2)!
85! sfcvisdfd(im)- real, tot sky sfc uv+vis-diff sw dnward flux (w/m2)!
86! im - integer, horizontal dimension !
87! levs - integer, vertical layer dimension !
88! deltim - real, physics time step in seconds !
89! fhswr - real, Short wave radiation time step in seconds !
90! dry - logical, true over land !
91! icy - logical, true over ice !
92! wet - logical, true over water !
93! !
94! input/output: !
95! dtdt(im,levs)- real, model time step adjusted total radiation !
96! heating rates ( k/s ) !
97! dtdtnp(im,levs)- real, heating rate adjustment for SPPT !
98! !
99! outputs: !
100! adjsfcdsw(im)- real, time step adjusted sfc dn sw flux (w/m**2) !
101! adjsfcnsw(im)- real, time step adj sfc net sw into ground (w/m**2)!
102! adjsfcdlw(im)- real, time step adjusted sfc dn lw flux (w/m**2) !
103! adjsfculw_lnd(im)- real, sfc upw. lw flux at current time (w/m**2)!
104! adjsfculw_ice(im)- real, sfc upw. lw flux at current time (w/m**2)!
105! adjsfculw_wat(im)- real, sfc upw. lw flux at current time (w/m**2)!
106! adjnirbmu(im)- real, t adj sfc nir-beam sw upward flux (w/m2) !
107! adjnirdfu(im)- real, t adj sfc nir-diff sw upward flux (w/m2) !
108! adjvisbmu(im)- real, t adj sfc uv+vis-beam sw upward flux (w/m2) !
109! adjvisdfu(im)- real, t adj sfc uv+vis-diff sw upward flux (w/m2) !
110! adjnirbmd(im)- real, t adj sfc nir-beam sw downward flux (w/m2) !
111! adjnirdfd(im)- real, t adj sfc nir-diff sw downward flux (w/m2) !
112! adjvisbmd(im)- real, t adj sfc uv+vis-beam sw dnward flux (w/m2) !
113! adjvisdfd(im)- real, t adj sfc uv+vis-diff sw dnward flux (w/m2) !
114! xmu (im) - real, time step zenith angle adjust factor for sw !
115! xcosz (im) - real, cosine of zenith angle at current time step !
116! !
117! ==================== end of description ===================== !
118
166 subroutine dcyc2t3_run &
167! --- inputs:
168 & ( solhr,slag,sdec,cdec,sinlat,coslat, &
169 & con_g, con_cp, con_pi, con_sbc, &
170 & xlon,coszen,tsfc_lnd,tsfc_ice,tsfc_wat,tf,tsflw,tsfc, &
171 & sfcemis_lnd, sfcemis_ice, sfcemis_wat, &
172 & sfcdsw,sfcnsw,sfcdlw,swh,swhc,hlw,hlwc, &
173 & sfcnirbmu,sfcnirdfu,sfcvisbmu,sfcvisdfu, &
174 & sfcnirbmd,sfcnirdfd,sfcvisbmd,sfcvisdfd, &
175 & im, levs, deltim, fhswr, &
176 & dry, icy, wet, damp_lw_fluxadj, lfnc_k, lfnc_p0, &
177 & use_lw_jacobian, sfculw, use_med_flux, sfculw_med, &
178 & fluxlwup_jac, t_lay, p_lay, p_lev, flux2d_lwup, &
179 & flux2d_lwdown,pert_radtend,do_sppt,ca_global,tsfc_radtime, &
180! & dry, icy, wet, lprnt, ipr, &
181! --- input/output:
182 & dtdt,dtdtnp,htrlw, &
183! --- outputs:
184 & adjsfcdsw,adjsfcnsw,adjsfcdlw, &
185 & adjsfculw_lnd,adjsfculw_ice,adjsfculw_wat,xmu,xcosz, &
186 & adjnirbmu,adjnirdfu,adjvisbmu,adjvisdfu, &
187 & adjnirbmd,adjnirdfd,adjvisbmd,adjvisdfd, &
188 & errmsg,errflg &
189 & )
190!
191 use machine, only : kind_phys
192
193 implicit none
194!
195! --- constant parameters:
196 real(kind=kind_phys), parameter :: f_eps = 0.0001_kind_phys, &
197 & zero = 0.0d0, one = 1.0d0, &
198 & hour12 = 12.0_kind_phys, &
199 & f3600 = one/3600.0_kind_phys, &
200 & f7200 = one/7200.0_kind_phys, &
201 & czlimt = 0.0001_kind_phys ! ~ cos(89.99427)
202
203! --- inputs:
204 integer, intent(in) :: im, levs
205
206! integer, intent(in) :: ipr
207! logical lprnt
208 logical, dimension(:), intent(in) :: dry, icy, wet
209 logical, intent(in) :: use_lw_jacobian, damp_lw_fluxadj, &
210 & pert_radtend, use_med_flux
211 logical, intent(in) :: do_sppt,ca_global
212 real(kind=kind_phys), intent(in) :: solhr, slag, cdec, sdec, &
213 & deltim, fhswr, lfnc_k, lfnc_p0
214
215 real(kind=kind_phys), dimension(:), intent(in) :: &
216 & sinlat, coslat, xlon, coszen, tf, tsflw, sfcdlw, &
217 & sfcdsw, sfcnsw, sfculw, tsfc
218 real(kind=kind_phys), dimension(:), intent(in), optional :: &
219 & sfculw_med, tsfc_radtime
220 real(kind=kind_phys), dimension(:), intent(in) :: &
221 & tsfc_lnd, tsfc_ice, tsfc_wat, &
222 & sfcemis_lnd, sfcemis_ice, sfcemis_wat
223
224 real(kind=kind_phys), dimension(:), intent(in) :: &
225 & sfcnirbmu, sfcnirdfu, sfcvisbmu, sfcvisdfu, &
226 & sfcnirbmd, sfcnirdfd, sfcvisbmd, sfcvisdfd
227
228 real(kind=kind_phys), dimension(:,:), intent(in) :: swh, hlw, &
229 & swhc, hlwc, p_lay, t_lay
230
231 real(kind=kind_phys), dimension(:,:), intent(in) :: p_lev
232 real(kind=kind_phys), dimension(:,:), intent(in), optional :: &
233 & flux2d_lwup, flux2d_lwdown, fluxlwup_jac
234
235 real(kind_phys), intent(in ) :: con_g, con_cp, &
236 & con_pi, con_sbc
237
238 real(kind_phys) :: pid12
239
240
241! --- input/output:
242 real(kind=kind_phys), dimension(:,:), intent(inout) :: dtdt
243 real(kind=kind_phys), dimension(:,:), intent(inout), optional :: &
244 & dtdtnp, htrlw
245
246! --- outputs:
247 real(kind=kind_phys), dimension(:), intent(out) :: &
248 & adjsfcdsw, adjsfcnsw, adjsfcdlw, xmu, xcosz, &
249 & adjnirbmu, adjnirdfu, adjvisbmu, adjvisdfu, &
250 & adjnirbmd, adjnirdfd, adjvisbmd, adjvisdfd
251
252 real(kind=kind_phys), dimension(:), intent(out) :: &
253 & adjsfculw_lnd, adjsfculw_ice, adjsfculw_wat
254
255 character(len=*), intent(out) :: errmsg
256 integer, intent(out) :: errflg
257
258! --- locals:
259 integer :: i, k, nstp, nstl, it, istsun(im),isfc,itoa
260 real(kind=kind_phys) :: cns, coszn, tem1, tem2, anginc, &
261 & rstl, solang, dt
262 real(kind=kind_phys), dimension(im,levs+1) :: flxlwup_adj, &
263 & flxlwdn_adj
264 real(kind=kind_phys) :: fluxlwnet_adj,fluxlwnet,dt_sfc, &
265 &fluxlwdown_jac,lfnc,c1
266 ! Length scale for flux-adjustment scaling
267 real(kind=kind_phys), parameter :: &
268 & l = 1.
269 ! Scaling factor for downwelling LW Jacobian profile.
270 real(kind=kind_phys), parameter :: &
271 & gamma = 0.2
272!
273!===> ... begin here
274!
275 ! Initialize CCPP error handling variables
276 errmsg = ''
277 errflg = 0
278
279! Vertical ordering?
280 if (p_lev(1,1) .lt. p_lev(1, levs)) then
281 isfc = levs + 1
282 itoa = 1
283 else
284 isfc = 1
285 itoa = levs + 1
286 endif
287
288 tem1 = fhswr / deltim
289 nstp = max(6, nint(tem1))
290 nstl = max(1, nint(nstp/tem1))
291 pid12 = con_pi / hour12
292!
293! --- ... sw time-step adjustment for current cosine of zenith angle
294! ----------------------------------------------------------
295 if (nstl == 1) then
296 cns = pid12 * (solhr + deltim*f7200 - hour12) + slag
297 do i = 1, im
298 xcosz(i) = sdec*sinlat(i) + cdec*coslat(i)*cos(cns+xlon(i))
299 enddo
300 elseif (nstl == nstp) then
301 do i = 1, im
302 xcosz(i) = coszen(i)
303 enddo
304 else
305 rstl = one / float(nstl)
306 solang = pid12 * (solhr - hour12)
307 anginc = pid12 * deltim * f3600 * rstl
308 do i = 1, im
309 xcosz(i) = zero
310 istsun(i) = zero
311 enddo
312 do it=1,nstl
313 cns = solang + (float(it)-0.5_kind_phys)*anginc + slag
314 do i = 1, im
315 coszn = sdec*sinlat(i) + cdec*coslat(i)*cos(cns+xlon(i))
316 xcosz(i) = xcosz(i) + max(zero, coszn)
317 if (coszn > czlimt) istsun(i) = istsun(i) + 1
318 enddo
319 enddo
320 do i = 1, im
321 if (istsun(i) > 0) xcosz(i) = xcosz(i) / istsun(i) ! mean cosine of solar zenith angle at current time
322 enddo
323 endif
324!
325
326 do i = 1, im
328 tem1 = tf(i) / tsflw(i)
329 tem2 = tem1 * tem1
330 adjsfcdlw(i) = sfcdlw(i) * tem2 * tem2
331!! - adjust \a sfc downward LW flux to account for t changes in the lowest model layer.
332!! compute 4th power of the ratio of \c tf in the lowest model layer over the mean value \c tsflw.
333 if (dry(i)) then
334 tem2 = tsfc_lnd(i) * tsfc_lnd(i)
335 adjsfculw_lnd(i) = sfcemis_lnd(i) * con_sbc * tem2 * tem2
336 & + (one - sfcemis_lnd(i)) * adjsfcdlw(i)
337 endif
338 if (icy(i)) then
339 tem2 = tsfc_ice(i) * tsfc_ice(i)
340 adjsfculw_ice(i) = sfcemis_ice(i) * con_sbc * tem2 * tem2
341 & + (one - sfcemis_ice(i)) * adjsfcdlw(i)
342 endif
343 if (wet(i)) then
344 tem2 = tsfc_wat(i) * tsfc_wat(i)
345 adjsfculw_wat(i) = sfcemis_wat(i) * con_sbc *
346 & tem2 * tem2
347 & + (one - sfcemis_wat(i)) * adjsfcdlw(i)
349 if (use_med_flux) then
350 if (sfculw_med(i) > f_eps) then
351 adjsfculw_wat(i) = sfculw_med(i)
352 end if
353 end if
354 endif
355
356! if (lprnt .and. i == ipr) write(0,*)' in dcyc3: dry==',dry(i)
357! &,' wet=',wet(i),' icy=',icy(i),' tsfc3=',tsfc3(i,:)
358! &,' sfcemis=',sfcemis(i,:)
359!
360
362 if ( xcosz(i) > f_eps .and. coszen(i) > f_eps ) then
363 xmu(i) = xcosz(i) / coszen(i)
364 else
365 xmu(i) = zero
366 endif
367
369! note: sfc emiss effect will not be appied here
370
371 adjsfcnsw(i) = sfcnsw(i) * xmu(i)
372 adjsfcdsw(i) = sfcdsw(i) * xmu(i)
373
374 adjnirbmu(i) = sfcnirbmu(i) * xmu(i)
375 adjnirdfu(i) = sfcnirdfu(i) * xmu(i)
376 adjvisbmu(i) = sfcvisbmu(i) * xmu(i)
377 adjvisdfu(i) = sfcvisdfu(i) * xmu(i)
378
379 adjnirbmd(i) = sfcnirbmd(i) * xmu(i)
380 adjnirdfd(i) = sfcnirdfd(i) * xmu(i)
381 adjvisbmd(i) = sfcvisbmd(i) * xmu(i)
382 adjvisdfd(i) = sfcvisdfd(i) * xmu(i)
383 enddo
384
385 ! Adjust the LW and SW heating-rates.
386 ! For LW, optionally scale using the Jacobian of the upward LW flux. *RRTMGP ONLY*
387 ! For SW, adjust heating rates with zenith angle change.
388 if (use_lw_jacobian) then
389 ! Compute adjusted net LW flux foillowing Hogan and Bozzo 2015 (10.1002/2015MS000455)
390 ! Here we assume that the profile of the downwelling LW Jaconiam has the same shape
391 ! as the upwelling, but scaled and offset.
392 ! The scaling factor is 0.2
393 ! The profile of the downwelling Jacobian (J) is offset so that
394 ! J_dn_sfc / J_up_sfc = scaling_factor
395 ! J_dn_toa / J_up_sfc = 0
396 !
397 ! Optionally, the flux adjustment can be damped with height using a logistic function
398 ! fx ~ L / (1 + exp(-k*dp)), where dp = p - p0
399 ! L = 1, fix scale between 0-1. - Fixed
400 ! k = 1 / pressure decay length (Pa) - Controlled by namelist
401 ! p0 = Transition pressure (Pa) - Controlled by namelsit
402 do i = 1, im
403 c1 = fluxlwup_jac(i,itoa) / fluxlwup_jac(i,isfc)
404 !dT_sfc = t_lev2(i,iSFC) - t_lev(i,iSFC)
405 dt_sfc = tsfc(i) - tsfc_radtime(i)
406 do k = 1, levs
407 ! LW net flux
408 fluxlwnet = (flux2d_lwup(i, k+1) - flux2d_lwup(i, k) - &
409 & flux2d_lwdown(i,k+1) + flux2d_lwdown(i,k))
410 ! Downward LW Jacobian (Eq. 9)
411 fluxlwdown_jac = gamma * &
412 & (fluxlwup_jac(i,k)/fluxlwup_jac(i,isfc) - c1) / &
413 & (1 - c1)
414 ! Adjusted LW net flux(Eq. 10)
415 fluxlwnet_adj = fluxlwnet + dt_sfc* &
416 & (fluxlwup_jac(i,k)/fluxlwup_jac(i,isfc) - &
417 & fluxlwdown_jac)
418 ! Adjusted LW heating rate
419 htrlw(i,k) = fluxlwnet_adj * con_g / &
420 & (con_cp * (p_lev(i,k+1) - p_lev(i,k)))
421
422 ! Add radiative heating rates to physics heating rate. Optionally, scaled w/ height
423 ! using a logistic function
424 if (damp_lw_fluxadj) then
425 lfnc = l / (1+exp(-(p_lev(i,k) - lfnc_p0)/lfnc_k))
426 else
427 lfnc = 1.
428 endif
429 dtdt(i,k) = dtdt(i,k) + swh(i,k)*xmu(i) + &
430 & htrlw(i,k)*lfnc + (1.-lfnc)*hlw(i,k)
431 enddo
432 enddo
433 else
434 do k = 1, levs
435 do i = 1, im
436 dtdt(i,k) = dtdt(i,k) + swh(i,k)*xmu(i) + hlw(i,k)
437 enddo
438 enddo
439 endif
440
441 if (do_sppt .or. ca_global) then
442 if (pert_radtend) then
443! clear sky
444 do k = 1, levs
445 do i = 1, im
446 dtdtnp(i,k) = dtdtnp(i,k) + swhc(i,k)*xmu(i) + hlwc(i,k)
447 enddo
448 enddo
449 else
450! all sky
451 do k = 1, levs
452 do i = 1, im
453 dtdtnp(i,k) = dtdtnp(i,k) + swh(i,k)*xmu(i) + hlw(i,k)
454 enddo
455 enddo
456 endif
457 endif
458!
459 return
460!...................................
461 end subroutine dcyc2t3_run
463!-----------------------------------
464 end module dcyc2t3
subroutine, public dcyc2t3_run(solhr, slag, sdec, cdec, sinlat, coslat, con_g, con_cp, con_pi, con_sbc, xlon, coszen, tsfc_lnd, tsfc_ice, tsfc_wat, tf, tsflw, tsfc, sfcemis_lnd, sfcemis_ice, sfcemis_wat, sfcdsw, sfcnsw, sfcdlw, swh, swhc, hlw, hlwc, sfcnirbmu, sfcnirdfu, sfcvisbmu, sfcvisdfu, sfcnirbmd, sfcnirdfd, sfcvisbmd, sfcvisdfd, im, levs, deltim, fhswr, dry, icy, wet, damp_lw_fluxadj, lfnc_k, lfnc_p0, use_lw_jacobian, sfculw, use_med_flux, sfculw_med, fluxlwup_jac, t_lay, p_lay, p_lev, flux2d_lwup, flux2d_lwdown, pert_radtend, do_sppt, ca_global, tsfc_radtime, dtdt, dtdtnp, htrlw, adjsfcdsw, adjsfcnsw, adjsfcdlw, adjsfculw_lnd, adjsfculw_ice, adjsfculw_wat, xmu, xcosz, adjnirbmu, adjnirdfu, adjvisbmu, adjvisdfu, adjnirbmd, adjnirdfd, adjvisbmd, adjvisdfd, errmsg, errflg)
Definition dcyc2t3.f:190