8 INTEGER,
PARAMETER ::
nz = 30
14 L_MOD, sfr, FAI, StabilityMethod, &
15 avcp, lv_J_kg, avdens, &
16 avU1, Temp_C, avRH, Press_hPa, zMeas, qh, qe, & ! input
17 T2_C, q2_gkg, U10_ms, RH2, &!output
32 REAL(KIND(1d0)),
DIMENSION(nsurf),
INTENT(in) ::sfr
33 REAL(KIND(1d0)),
INTENT(in):: zMeas
34 REAL(KIND(1d0)),
INTENT(in):: avU1
35 REAL(KIND(1d0)),
INTENT(in):: Temp_C
36 REAL(KIND(1d0)),
INTENT(in):: avRH
37 REAL(KIND(1d0)),
INTENT(in):: Press_hPa
38 REAL(KIND(1d0)),
INTENT(in):: L_MOD
39 REAL(KIND(1d0)),
INTENT(in):: avcp
40 REAL(KIND(1d0)),
INTENT(in):: lv_J_kg
41 REAL(KIND(1d0)),
INTENT(in):: avdens
42 REAL(KIND(1d0)),
INTENT(in):: qh
43 REAL(KIND(1d0)),
INTENT(in):: qe
44 REAL(KIND(1d0)),
INTENT(in):: Zh
45 REAL(KIND(1d0)),
INTENT(in):: z0m
46 REAL(KIND(1d0)),
INTENT(in):: zdm
47 REAL(KIND(1d0)),
INTENT(in):: FAI
49 INTEGER,
INTENT(in)::StabilityMethod
51 REAL(KIND(1d0)),
INTENT(out):: T2_C
52 REAL(KIND(1d0)),
INTENT(out):: q2_gkg
53 REAL(KIND(1d0)),
INTENT(out):: U10_ms
54 REAL(KIND(1d0)),
INTENT(out):: RH2
56 INTEGER,
PARAMETER :: nz = 30
58 REAL(KIND(1d0)),
PARAMETER:: cd_tree = 1.2, &
63 pi = 4.*atan(1.0), r = 0.1, &
64 a1 = 4., a2 = -0.1, a3 = 1.5, a4 = -1.
67 REAL(KIND(1d0)),
INTENT(out),
DIMENSION(nz*4):: dataoutLineRSL
68 REAL(KIND(1d0)),
DIMENSION(nz):: psihatm_z
69 REAL(KIND(1d0)),
DIMENSION(nz):: psihath_z
70 REAL(KIND(1d0)),
DIMENSION(nz):: dif
72 REAL(KIND(1d0)),
DIMENSION(nz):: zarray
73 REAL(KIND(1d0)),
DIMENSION(nz):: dataoutLineURSL
74 REAL(KIND(1d0)),
DIMENSION(nz):: dataoutLineTRSL
75 REAL(KIND(1d0)),
DIMENSION(nz):: dataoutLineqRSL
80 REAL(KIND(1d0))::Lc_build, Lc_tree, Lc
82 REAL(KIND(1d0))::phim, psimz, psimZh, psimz0, psimza, phi_hatmZh, phi_hathZh, phimzp, phimz, phihzp, phihz, psihz, psihza
83 REAL(KIND(1d0))::betaHF, betaNL, beta, betaN2
85 REAL(KIND(1d0))::xxm1, xxm1_2, xxh1, xxh1_2, err, z01, dphi, dphih
86 REAL(KIND(1d0))::f, cm, c2, ch, c2h
87 REAL(KIND(1d0))::t_h, q_h
88 REAL(KIND(1d0))::TStar_RSL
89 REAL(KIND(1d0))::UStar_RSL
90 REAL(KIND(1d0))::sfr_zh
91 REAL(KIND(1d0))::sfr_tr
92 REAL(KIND(1d0))::L_MOD_RSL
93 REAL(KIND(1d0))::zH_RSL
95 REAL(KIND(1d0)),
parameter::Zh_min = 0.15
96 REAL(KIND(1d0)),
parameter::ratio_dz = 1.618
98 REAL(KIND(1d0))::qa_gkg, qStar_RSL
99 INTEGER :: I, z, it, idx_can, idx_za, idx_2m, idx_10m
116 stabilitymethod, zh, l_mod, sfr, fai, &
117 l_mod_rsl, zh_rsl, lc, beta, zd, z0, elm, scc, f)
121 IF (zh_rsl <= 2)
THEN 123 ELSE IF (zh_rsl <= 10)
THEN 134 dz = (zmeas - zh_rsl)/(nz - nz_can)
135 do i = nz_can + 1, nz
136 zarray(i) = zh_rsl + (i - nz_can)*dz
142 dif(z) = abs(zarray(z) - 2)
144 idx_2m = minloc(dif, dim=1)
148 dif(z) = abs(zarray(z) - 10)
150 idx_10m = minloc(dif, dim=1)
155 dif(z) = abs(zarray(z) - zh_rsl)
157 idx_can = minloc(dif, dim=1)
158 zarray(idx_can) = zh_rsl
163 dif(z) = abs(zarray(z) - zmeas)
165 idx_za = minloc(dif, dim=1)
166 zarray(idx_za) = zmeas
168 if (zh_rsl - zd < z0 .or. zh < zh_min)
then 186 xxm1 =
stab_phi_mom(stabilitymethod, (zh_rsl - zd)/l_mod_rsl)
187 xxm1_2 =
stab_phi_mom(stabilitymethod, (zh_rsl - zd + 1.)/l_mod_rsl)
189 xxh1 =
stab_phi_heat(stabilitymethod, (zh_rsl - zd)/l_mod_rsl)
190 xxh1_2 =
stab_phi_heat(stabilitymethod, (zh_rsl - zd + 1.)/l_mod_rsl)
192 phi_hatmzh = kappa/(2.*beta*xxm1)
193 phi_hathzh = kappa*scc/(2.*beta*xxh1)
195 dphih = xxh1_2 - xxh1
196 IF (phi_hatmzh > 1.)
THEN 200 c2 = (kappa*(3.-(2.*beta**2.*lc/xxm1*dphi)))/(2.*beta*xxm1 - kappa)
201 c2h = (kappa*scc*(2.+f - (dphih*2.*beta**2.*lc/xxh1)))/(2.*beta*xxh1 - kappa*scc)
203 cm = (1.-phi_hatmzh)*exp(c2/2.)
204 ch = (1.-phi_hathzh)*exp(c2h/2.)
209 psihatm_z = 0.*zarray
210 psihath_z = 0.*zarray
213 DO z = nz - 1, idx_can - 1, -1
214 phimz =
stab_phi_mom(stabilitymethod, (zarray(z) - zd)/l_mod_rsl)
215 phimzp =
stab_phi_mom(stabilitymethod, (zarray(z + 1) - zd)/l_mod_rsl)
216 phihz =
stab_phi_heat(stabilitymethod, (zarray(z) - zd)/l_mod_rsl)
217 phihzp =
stab_phi_heat(stabilitymethod, (zarray(z + 1) - zd)/l_mod_rsl)
219 psihatm_z(z) = psihatm_z(z + 1) + dz/2.*phimzp*(cm*exp(-1.*c2*beta*(zarray(z + 1) - zd)/elm)) &
220 /(zarray(z + 1) - zd)
221 psihatm_z(z) = psihatm_z(z) + dz/2.*phimz*(cm*exp(-1.*c2*beta*(zarray(z) - zd)/elm)) &
223 psihath_z(z) = psihath_z(z + 1) + dz/2.*phihzp*(ch*exp(-1.*c2h*beta*(zarray(z + 1) - zd)/elm)) &
224 /(zarray(z + 1) - zd)
225 psihath_z(z) = psihath_z(z) + dz/2.*phihz*(ch*exp(-1.*c2h*beta*(zarray(z) - zd)/elm)) &
236 psimza =
stab_psi_mom(stabilitymethod, (zmeas - zd)/l_mod_rsl)
237 psihza =
stab_psi_heat(stabilitymethod, (zmeas - zd)/l_mod_rsl)
238 ustar_rsl = avu1*kappa/(log((zmeas - zd)/z0) - psimza + psimz0 + psihatm_z(nz))
239 tstar_rsl = -1.*(qh/(avcp*avdens))/ustar_rsl
240 qstar_rsl = -1.*(qe/lv_j_kg*avdens)/ustar_rsl
241 qa_gkg =
rh2qa(avrh/100, press_hpa, temp_c)
244 psimz =
stab_psi_mom(stabilitymethod, (zarray(z) - zd)/l_mod_rsl)
245 psihz =
stab_psi_heat(stabilitymethod, (zarray(z) - zd)/l_mod_rsl)
247 dataoutlineursl(z) = (log((zarray(z) - zd)/z0) - psimz + psimz0 + psihatm_z(z))/kappa
249 dataoutlinetrsl(z) = (log((zarray(z) - zd)/(zmeas - zd)) - psihz + psihza + psihath_z(z) - psihath_z(idx_za))/kappa
250 dataoutlineqrsl(z) = (log((zarray(z) - zd)/(zmeas - zd)) - psihz + psihza + psihath_z(z) - psihath_z(idx_za))/kappa
256 t_h = scc*tstar_rsl/(beta*f)
257 q_h = scc*qstar_rsl/(beta*f)
259 dataoutlineursl(z) = dataoutlineursl(idx_can)*exp(beta*(zarray(z) - zh_rsl)/elm)
260 dataoutlinetrsl(z) = dataoutlinetrsl(idx_can) + (t_h*exp(beta*f*(zarray(z) - zh_rsl)/elm) - t_h)/tstar_rsl
261 dataoutlineqrsl(z) = dataoutlineqrsl(idx_can) + (q_h*exp(beta*f*(zarray(z) - zh_rsl)/elm) - q_h)/qstar_rsl
264 dataoutlineursl = dataoutlineursl*ustar_rsl
265 dataoutlinetrsl = dataoutlinetrsl*tstar_rsl + temp_c
266 dataoutlineqrsl = (dataoutlineqrsl*qstar_rsl + qa_gkg/1000.)*1000.
268 dataoutlinersl = (/zarray, dataoutlineursl, dataoutlinetrsl, dataoutlineqrsl/)
274 t2_c = dataoutlinetrsl(idx_2m)
275 q2_gkg = dataoutlineqrsl(idx_2m)
276 u10_ms = dataoutlineursl(idx_10m)
278 rh2 =
qa2rh(q2_gkg, press_hpa, t2_c)
292 recursive function cal_psim_hat(StabilityMethod, z, zh_RSL, L_MOD_RSL, beta, Lc)
result(psim_hat_z)
296 integer,
intent(in) :: StabilityMethod
297 real(KIND(1D0)),
intent(in) :: z
298 real(KIND(1D0)),
intent(in) :: zh_RSL
299 real(KIND(1D0)),
intent(in) :: Lc
300 real(KIND(1D0)),
intent(in) :: beta
301 real(KIND(1D0)),
intent(in) :: L_MOD_RSL
304 real(KIND(1D0)) ::psim_hat_z
309 real(KIND(1D0)) ::phim_lc
310 real(KIND(1D0)) ::phim_z
311 real(KIND(1D0)) ::phim_zp
312 real(KIND(1D0)) ::psim_hat_zp
313 real(KIND(1D0)) ::elm
314 real(KIND(1D0)) ::xxm1
315 real(KIND(1D0)) ::xxm1_2
316 real(KIND(1D0)) ::dphi
317 real(KIND(1D0)) ::phi_hatmZh
321 REAL(KIND(1d0)),
PARAMETER::kappa = 0.40
322 REAL(KIND(1d0)),
PARAMETER::dz = 0.1
331 zd = zh_rsl - (beta**2.)*lc
337 phim_z =
stab_phi_mom(stabilitymethod, (z - zd)/l_mod_rsl)
338 phim_zp =
stab_phi_mom(stabilitymethod, (zp - zd)/l_mod_rsl)
340 psim_hat_zp =
cal_psim_hat(stabilitymethod, zp, zh_rsl, l_mod_rsl, beta, lc)
342 xxm1 =
stab_phi_mom(stabilitymethod, (zh_rsl - zd)/l_mod_rsl)
343 xxm1_2 =
stab_phi_mom(stabilitymethod, (zh_rsl - zd + dz)/l_mod_rsl)
345 dphi = (xxm1_2 - xxm1)/dz
347 phi_hatmzh = kappa/(2.*beta*xxm1)
349 IF (phi_hatmzh > 1.)
THEN 354 c2 = (kappa*(3.-(2.*beta**2.*lc/xxm1*dphi)))/(2.*beta*xxm1 - kappa)
357 cm = (1.-phi_hatmzh)*exp(c2/2.)
360 psim_hat_z = psim_hat_zp + dz/2.*phim_zp*(cm*exp(-1.*c2*beta*(zp - zd)/elm))/(zp - zd)
361 psim_hat_z = psim_hat_z + dz/2.*phim_z*(cm*exp(-1.*c2*beta*(z - zd)/elm))/(z - zd)
365 function cal_psihatm_z(StabilityMethod, nz, zarray, L_MOD_RSL, zH_RSL, Lc, beta, zd, elm)
result(psihatm_z)
370 integer,
intent(in) :: StabilityMethod
371 integer,
intent(in) :: nz
372 real(KIND(1D0)),
DIMENSION(nz),
intent(in) :: zarray
373 real(KIND(1D0)),
intent(in) :: zh_RSL
374 real(KIND(1D0)),
intent(in) :: Lc
375 real(KIND(1D0)),
intent(in) :: beta
376 real(KIND(1D0)),
intent(in) :: L_MOD_RSL
377 real(KIND(1D0)),
intent(in) ::zd
378 real(KIND(1D0)),
intent(in) ::elm
381 real(KIND(1D0)),
DIMENSION(nz) ::psihatm_z
386 real(KIND(1D0)) ::phimz
387 real(KIND(1D0)) ::phimzp
389 real(KIND(1D0)) ::xxm1
390 real(KIND(1D0)) ::xxm1_2
391 real(KIND(1D0)) ::dphi
392 real(KIND(1D0)) ::phi_hatmZh
395 REAL(KIND(1d0)),
DIMENSION(nz):: dif
397 REAL(KIND(1d0)),
PARAMETER::kappa = 0.40
398 REAL(KIND(1d0)),
PARAMETER::dz = 0.1
402 psihatm_z = 0.*zarray
406 dif(z) = abs(zarray(z) - zh_rsl)
408 idx_can = minloc(dif, dim=1)
412 xxm1 =
stab_phi_mom(stabilitymethod, (zh_rsl - zd)/l_mod_rsl)
413 xxm1_2 =
stab_phi_mom(stabilitymethod, (zh_rsl - zd + 1.)/l_mod_rsl)
415 phi_hatmzh = kappa/(2.*beta*xxm1)
418 IF (phi_hatmzh > 1.)
THEN 422 c2 = (kappa*(3.-(2.*beta**2.*lc/xxm1*dphi)))/(2.*beta*xxm1 - kappa)
425 cm = (1.-phi_hatmzh)*exp(c2/2.)
428 DO z = nz - 1, idx_can - 1, -1
429 phimz =
stab_phi_mom(stabilitymethod, (zarray(z) - zd)/l_mod_rsl)
430 phimzp =
stab_phi_mom(stabilitymethod, (zarray(z + 1) - zd)/l_mod_rsl)
432 psihatm_z(z) = psihatm_z(z + 1) + dz/2.*phimzp*(cm*exp(-1.*c2*beta*(zarray(z + 1) - zd)/elm)) &
433 /(zarray(z + 1) - zd)
434 psihatm_z(z) = psihatm_z(z) + dz/2.*phimz*(cm*exp(-1.*c2*beta*(zarray(z) - zd)/elm)) &
441 function cal_psihath_z(StabilityMethod, nz, zarray, L_MOD_RSL, zH_RSL, Lc, beta, zd, elm, Scc, f)
result(psihath_z)
446 integer,
intent(in) :: StabilityMethod
447 integer,
intent(in) :: nz
449 real(KIND(1D0)),
DIMENSION(nz),
intent(in) :: zarray
450 real(KIND(1D0)),
intent(in) :: zh_RSL
451 real(KIND(1D0)),
intent(in) :: Lc
452 real(KIND(1D0)),
intent(in) :: beta
453 real(KIND(1D0)),
intent(in) :: Scc
454 real(KIND(1D0)),
intent(in) :: f
455 real(KIND(1D0)),
intent(in) :: L_MOD_RSL
456 real(KIND(1D0)),
intent(in) :: elm
457 real(KIND(1D0)),
intent(in) ::zd
460 real(KIND(1D0)),
DIMENSION(nz) ::psihath_z
465 real(KIND(1D0)) ::phihz
466 real(KIND(1D0)) ::phihzp
468 real(KIND(1D0)) ::xxh1
469 real(KIND(1D0)) ::xxh1_2
470 real(KIND(1D0)) ::dphih
471 real(KIND(1D0)) ::phi_hathZh
473 real(KIND(1D0)) ::c2h
474 REAL(KIND(1d0)),
DIMENSION(nz):: dif
476 REAL(KIND(1d0)),
PARAMETER::kappa = 0.40
477 REAL(KIND(1d0)),
PARAMETER::dz = 0.1
481 psihath_z = 0.*zarray
485 dif(z) = abs(zarray(z) - zh_rsl)
487 idx_can = minloc(dif, dim=1)
491 xxh1 =
stab_phi_heat(stabilitymethod, (zh_rsl - zd)/l_mod_rsl)
492 xxh1_2 =
stab_phi_heat(stabilitymethod, (zh_rsl - zd + 1.)/l_mod_rsl)
494 phi_hathzh = kappa*scc/(2.*beta*xxh1)
495 dphih = xxh1_2 - xxh1
497 IF (phi_hathzh > 1.)
THEN 502 c2h = (kappa*scc*(2.+f - (dphih*2.*beta**2.*lc/xxh1)))/(2.*beta*xxh1 - kappa*scc)
505 ch = (1.-phi_hathzh)*exp(c2h/2.)
507 DO z = nz - 1, idx_can - 1, -1
508 phihz =
stab_phi_heat(stabilitymethod, (zarray(z) - zd)/l_mod_rsl)
509 phihzp =
stab_phi_heat(stabilitymethod, (zarray(z + 1) - zd)/l_mod_rsl)
511 psihath_z(z) = psihath_z(z + 1) + dz/2.*phihzp*(ch*exp(-1.*c2h*beta*(zarray(z + 1) - zd)/elm)) &
512 /(zarray(z + 1) - zd)
513 psihath_z(z) = psihath_z(z) + dz/2.*phihz*(ch*exp(-1.*c2h*beta*(zarray(z) - zd)/elm)) &
520 function rsl_cal_z0(StabilityMethod, zH_RSL, zd, beta, L_MOD_RSL, Lc)
result(z0)
524 integer,
intent(in) ::StabilityMethod
525 real(KIND(1D0)),
intent(in) :: zH_RSL
526 real(KIND(1D0)),
intent(in) :: zd
527 real(KIND(1D0)),
intent(in) :: L_MOD_RSL
528 real(KIND(1D0)),
intent(in) :: Lc
529 real(KIND(1D0)),
intent(in) :: beta
535 real(KIND(1D0)) ::psimZh, psimz0, z01, psihatm_Zh
536 real(KIND(1D0)) ::err
539 REAL(KIND(1d0)),
PARAMETER::kappa = 0.40
543 psimzh =
stab_psi_mom(stabilitymethod, (zh_rsl - zd)/l_mod_rsl)
544 psihatm_zh =
cal_psim_hat(stabilitymethod, zh_rsl, zh_rsl, l_mod_rsl, beta, lc)
553 DO WHILE ((err > 0.001) .AND. (it < 10))
556 z0 = (zh_rsl - zd)*exp(-1.*kappa/beta)*exp(-1.*psimzh + psimz0)*exp(psihatm_zh)
565 StabilityMethod, zh, L_MOD, sfr, planF, &!input
566 L_MOD_RSL, zH_RSL, Lc, beta, zd, z0, elm, Scc, f)
570 integer,
intent(in) :: StabilityMethod
571 real(KIND(1D0)),
intent(in) :: zh_min
572 real(KIND(1D0)),
intent(in) :: z0m
573 real(KIND(1D0)),
intent(in) :: zdm
574 real(KIND(1D0)),
intent(in) :: zh
575 real(KIND(1D0)),
intent(in) :: planF
576 real(KIND(1D0)),
intent(in) :: L_MOD
577 real(KIND(1D0)),
DIMENSION(nsurf),
intent(in) :: sfr
580 real(KIND(1D0)),
intent(out) ::L_MOD_RSL
581 real(KIND(1D0)),
intent(out) ::zH_RSL
582 real(KIND(1D0)),
intent(out) ::Lc
583 real(KIND(1D0)),
intent(out) ::beta
584 real(KIND(1D0)),
intent(out) ::zd
585 real(KIND(1D0)),
intent(out) ::z0
586 real(KIND(1D0)),
intent(out) ::elm
587 real(KIND(1D0)),
intent(out) ::Scc
588 real(KIND(1D0)),
intent(out) ::f
591 real(KIND(1D0)) ::sfr_zh
592 real(KIND(1D0)) ::sfr_tr
593 real(KIND(1D0)) ::phim
594 real(KIND(1D0)) ::betaN2
595 real(KIND(1D0)) ::betaHF
596 real(KIND(1D0)) ::betaNL
598 REAL(KIND(1d0)),
PARAMETER::planF_low = 1e-6
599 REAL(KIND(1d0)),
PARAMETER::kappa = 0.40
601 REAL(KIND(1d0)),
PARAMETER::r = 0.1
602 REAL(KIND(1d0)),
PARAMETER::a1 = 4., a2 = -0.1, a3 = 1.5, a4 = -1.
606 l_mod_rsl = merge(l_mod, max(300., l_mod), l_mod < 0)
609 zh_rsl = max(zh, zh_min)
615 sfr_zh = min(sfr_zh, 0.8)
627 lc = (1.-sfr_zh)/planf*zh_rsl
632 scc = 0.5 + 0.3*tanh(2.*lc/l_mod_rsl)
633 f = 0.5*((1.+4.*r*scc)**0.5) - 0.5
641 betan2 = 0.30*sfr_tr/sfr_zh + (sfr_zh - sfr_tr)/sfr_zh*0.4
647 betanl = (kappa/2.)/phim
649 IF (lc/l_mod_rsl > a2)
THEN 652 beta = betanl + ((betahf - betanl)/(1.+a1*abs(lc/l_mod_rsl - a2)**a3))
658 zd = zh_rsl - (beta**2.)*lc
662 z0 =
rsl_cal_z0(stabilitymethod, zh_rsl, zd, beta, l_mod_rsl, lc)
real(kind(1d0)) function stab_phi_mom(StabilityMethod, ZL)
real(kind(1d0)) function, dimension(nz) cal_psihatm_z(StabilityMethod, nz, zarray, L_MOD_RSL, zH_RSL, Lc, beta, zd, elm)
real(kind(1d0)) function qa2rh(qa_gkg, pres_hPa, Ta_degC)
subroutine rsl_cal_prms(zh_min, z0m, zdm, StabilityMethod, zh, L_MOD, sfr, planF, L_MOD_RSL, zH_RSL, Lc, beta, zd, z0, elm, Scc, f)
real(kind(1d0)) function, dimension(nz) cal_psihath_z(StabilityMethod, nz, zarray, L_MOD_RSL, zH_RSL, Lc, beta, zd, elm, Scc, f)
subroutine cal_stab(StabilityMethod, zzd, z0m, zdm, avU1, Temp_C, QH_init, avdens, avcp, L_MOD, TStar, UStar, zL)
integer, parameter conifsurf
real(kind(1d0)) function stab_psi_heat(StabilityMethod, ZL)
real(kind(1d0)) function rsl_cal_z0(StabilityMethod, zH_RSL, zd, beta, L_MOD_RSL, Lc)
recursive real(kind(1d0)) function cal_psim_hat(StabilityMethod, z, zh_RSL, L_MOD_RSL, beta, Lc)
real(kind(1d0)) function rh2qa(RH_dec, pres_hPa, Ta_degC)
real(kind(1d0)) function stab_phi_heat(StabilityMethod, ZL)
subroutine rslprofile(Zh, z0m, zdm, L_MOD, sfr, FAI, StabilityMethod, avcp, lv_J_kg, avdens, avU1, Temp_C, avRH, Press_hPa, zMeas, qh, qe, T2_C, q2_gkg, U10_ms, RH2, dataoutLineRSL)
real(kind(1d0)) function stab_psi_mom(StabilityMethod, ZL)
integer, parameter decidsurf
integer, parameter bldgsurf