9 INTEGER,
PARAMETER ::
nz = 30
16 L_MOD, sfr_surf, FAI, PAI, &
17 StabilityMethod, RA_h, &
18 avcp, lv_J_kg, avdens, &
19 avU1, Temp_C, avRH, Press_hPa, zMeas, qh, qe, & ! input
20 T2_C, q2_gkg, U10_ms, RH2, & !output
35 REAL(KIND(1D0)),
DIMENSION(nsurf),
INTENT(in) :: sfr_surf
36 REAL(KIND(1D0)),
INTENT(in) :: zMeas
37 REAL(KIND(1D0)),
INTENT(in) :: avU1
38 REAL(KIND(1D0)),
INTENT(in) :: Temp_C
39 REAL(KIND(1D0)),
INTENT(in) :: avRH
40 REAL(KIND(1D0)),
INTENT(in) :: Press_hPa
41 REAL(KIND(1D0)),
INTENT(in) :: L_MOD
42 REAL(KIND(1D0)),
INTENT(in) :: RA_h
43 REAL(KIND(1D0)),
INTENT(in) :: avcp
44 REAL(KIND(1D0)),
INTENT(in) :: lv_J_kg
45 REAL(KIND(1D0)),
INTENT(in) :: avdens
46 REAL(KIND(1D0)),
INTENT(in) :: qh
47 REAL(KIND(1D0)),
INTENT(in) :: qe
48 REAL(KIND(1D0)),
INTENT(in) :: Zh
49 REAL(KIND(1D0)),
INTENT(in) :: z0m
50 REAL(KIND(1D0)),
INTENT(in) :: z0v
51 REAL(KIND(1D0)),
INTENT(in) :: zdm
52 REAL(KIND(1D0)),
INTENT(in) :: FAI
53 REAL(KIND(1D0)),
INTENT(in) :: PAI
57 INTEGER,
INTENT(in) :: StabilityMethod
58 INTEGER,
INTENT(in) :: DiagMethod
60 REAL(KIND(1D0)),
INTENT(out) :: T2_C
61 REAL(KIND(1D0)),
INTENT(out) :: q2_gkg
62 REAL(KIND(1D0)),
INTENT(out) :: U10_ms
63 REAL(KIND(1D0)),
INTENT(out) :: RH2
65 INTEGER,
PARAMETER :: nz = 30
67 REAL(KIND(1D0)),
PARAMETER :: cd_tree = 1.2
68 REAL(KIND(1D0)),
PARAMETER :: a_tree = 0.05
69 REAL(KIND(1D0)),
PARAMETER :: kappa = 0.40
71 REAL(KIND(1D0)),
PARAMETER :: beta_N = 0.40
72 REAL(KIND(1D0)),
PARAMETER :: pi = 4.*atan(1.0), r = 0.1
73 REAL(KIND(1D0)),
PARAMETER :: a1 = 4.,
a2 = -0.1,
a3 = 1.5, a4 = -1.
75 REAL(KIND(1D0)),
PARAMETER :: porosity_evetr = 0.32
81 REAL(KIND(1D0)),
INTENT(out),
DIMENSION(ncolumnsDataOutRSL - 5) :: dataoutLineRSL
82 REAL(KIND(1D0)),
DIMENSION(nz) :: psihatm_z
83 REAL(KIND(1D0)),
DIMENSION(nz) :: psihath_z
86 REAL(KIND(1D0)),
DIMENSION(nz) :: zarray
87 REAL(KIND(1D0)),
DIMENSION(nz) :: dataoutLineURSL
88 REAL(KIND(1D0)),
DIMENSION(nz) :: dataoutLineTRSL
89 REAL(KIND(1D0)),
DIMENSION(nz) :: dataoutLineqRSL
91 REAL(KIND(1D0)) :: z0_RSL
92 REAL(KIND(1D0)) :: zd_RSL
98 REAL(KIND(1D0)) :: Scc
99 REAL(KIND(1D0)) :: psimz, psimz0, psimza, psihz, psihz0, psihza
101 REAL(KIND(1D0)) :: beta
102 REAL(KIND(1D0)) :: elm
104 REAL(KIND(1D0)) :: fx
105 REAL(KIND(1D0)) :: t_h, q_h
106 REAL(KIND(1D0)) :: TStar_RSL
107 REAL(KIND(1D0)) :: UStar_RSL
108 REAL(KIND(1D0)) :: UStar_heat
111 REAL(KIND(1D0)) :: L_MOD_RSL
115 REAL(KIND(1D0)) :: zH_RSL
116 REAL(KIND(1D0)) :: dz_above
117 REAL(KIND(1D0)) :: dz_can
122 REAL(KIND(1D0)) :: qa_gkg, qStar_RSL
148 IF (diagmethod == 0)
THEN
151 ELSEIF (diagmethod == 1)
THEN
154 ELSEIF (diagmethod == 2)
THEN
159 fai > 0.25*(1 - pai) .AND. fai < 0.45 .AND. &
161 pai > 0.1 .AND. pai < 0.68 .AND. &
177 ELSE IF (zh <= 10)
THEN
189 IF (dz_can > 2) zarray(1) = 1.999
192 nz_above = nz - nz_can
193 dz_above = (zmeas - zh)/nz_above
194 DO i = nz_can + 1, nz
195 zarray(i) = zh + (i - nz_can)*dz_above
215 psihatm_z = 0.*zarray
216 psihath_z = 0.*zarray
220 nz_above, zarray(nz_can + 1:nz), &
221 zh, l_mod, sfr_surf, fai, pai, &
222 psihatm_z(nz_can + 1:nz), psihath_z(nz_can + 1:nz), &
224 lc, beta, zd_rsl, z0_rsl, elm, scc, fx)
275 psimz0 =
stab_psi_mom(stabilitymethod, z0_rsl/l_mod_rsl)
276 psimza =
stab_psi_mom(stabilitymethod, (zmeas - zd_rsl)/l_mod_rsl)
277 psihza =
stab_psi_heat(stabilitymethod, (zmeas - zd_rsl)/l_mod_rsl)
279 ustar_rsl = avu1*kappa/(log((zmeas - zd_rsl)/z0_rsl) - psimza + psimz0 + psihatm_z(nz))
283 ustar_rsl = max(0.001, ustar_rsl)
285 IF ((zmeas - zd_rsl)/l_mod_rsl < -
neut_limit) ustar_rsl = max(0.15, ustar_rsl)
290 ustar_heat = max(0.15, ustar_rsl)
294 ustar_heat = 1/(kappa*ra_h)*(log((zmeas - zd_rsl)/z0v) - psihza + psihz0)
296 tstar_rsl = -1.*(qh/(avcp*avdens))/ustar_heat
298 qstar_rsl = 10.**(-10)
300 qstar_rsl = -1.*(qe/lv_j_kg*avdens)/ustar_heat
302 qa_gkg =
rh2qa(avrh/100, press_hpa, temp_c)
305 psimz =
stab_psi_mom(stabilitymethod, (zarray(z) - zd_rsl)/l_mod_rsl)
306 psihz =
stab_psi_heat(stabilitymethod, (zarray(z) - zd_rsl)/l_mod_rsl)
307 dataoutlineursl(z) = (log((zarray(z) - zd_rsl)/z0_rsl) - psimz + psimz0 + psihatm_z(z))/kappa
308 dataoutlinetrsl(z) = (log((zarray(z) - zd_rsl)/(zmeas - zd_rsl)) - psihz + psihza + psihath_z(z) - psihath_z(nz))/kappa
309 dataoutlineqrsl(z) = dataoutlinetrsl(z)
317 t_h = scc*tstar_rsl/(beta*fx)
318 q_h = scc*qstar_rsl/(beta*fx)
320 dataoutlineursl(z) = dataoutlineursl(nz_can)*exp(beta*(zarray(z) - zh_rsl)/elm)
321 dataoutlinetrsl(z) = dataoutlinetrsl(nz_can) + (t_h*exp(beta*fx*(zarray(z) - zh_rsl)/elm) - t_h)/tstar_rsl
322 dataoutlineqrsl(z) = dataoutlineqrsl(nz_can) + (q_h*exp(beta*fx*(zarray(z) - zh_rsl)/elm) - q_h)/qstar_rsl
329 IF (zarray(z) <= zd_rsl) zarray(z) = 1.01*(zd_rsl + z0_rsl)
330 psimz =
stab_psi_mom(stabilitymethod, (zarray(z) - zd_rsl)/l_mod_rsl)
331 psihz =
stab_psi_heat(stabilitymethod, (zarray(z) - zd_rsl)/l_mod_rsl)
332 dataoutlineursl(z) = (log((zarray(z) - zd_rsl)/z0_rsl) - psimz + psimz0)/kappa
333 dataoutlinetrsl(z) = (log((zarray(z) - zd_rsl)/(zmeas - zd_rsl)) - psihz + psihza)/kappa
334 dataoutlineqrsl(z) = dataoutlinetrsl(z)
339 dataoutlineursl = dataoutlineursl*ustar_rsl
340 dataoutlinetrsl = dataoutlinetrsl*tstar_rsl + temp_c
341 dataoutlineqrsl = (dataoutlineqrsl*qstar_rsl + qa_gkg/1000.)*1000.
344 dataoutlinersl = [zarray, dataoutlineursl, dataoutlinetrsl, dataoutlineqrsl, &
351 beta, zd_rsl, z0_rsl, elm, scc, fx, &
352 ustar_rsl, ustar_heat, tstar_rsl, fai, pai, merge(1.d0, 0.d0, flag_rsl) &
361 t2_c =
interp_z(2d0, zarray, dataoutlinetrsl)
362 q2_gkg =
interp_z(2d0, zarray, dataoutlineqrsl)
363 u10_ms =
interp_z(10d0, zarray, dataoutlineursl)
366 t2_c =
interp_z(2d0 + zd_rsl + z0_rsl, zarray, dataoutlinetrsl)
367 q2_gkg =
interp_z(2d0 + zd_rsl + z0_rsl, zarray, dataoutlineqrsl)
368 u10_ms =
interp_z(10d0 + zd_rsl + z0_rsl, zarray, dataoutlineursl)
371 rh2 =
qa2rh(q2_gkg, press_hpa, t2_c)
377 REAL(kind(1d0)),
INTENT(in) :: z_x
378 REAL(kind(1d0)),
DIMENSION(nz),
INTENT(in) :: z
379 REAL(kind(1d0)),
DIMENSION(nz),
INTENT(in) :: v
382 REAL(kind(1d0)) :: v_x
385 REAL(kind(1d0)) :: slope
386 REAL(kind(1d0)) :: dz
387 REAL(kind(1d0)),
DIMENSION(nz) :: dif
392 INTEGER,
PARAMETER ::
nz = 30
398 idx_x = maxloc(dif, 1, abs(dif) < 1.d-6)
399 idx_low = maxloc(dif, 1, dif < 0.)
400 idx_high = minloc(dif, 1, dif > 0.)
407 dz = z(idx_high) - z(idx_low)
408 slope = (v(idx_high) - v(idx_low))/dz
409 v_x = v(idx_low) + (z_x - z(idx_low))*slope
416 REAL(kind(1d0)),
INTENT(in) :: lc
417 REAL(kind(1d0)),
INTENT(in) :: beta
420 REAL(kind(1d0)) :: elm
427 psihatm_top, psihatm_mid, &
428 z_top, z_mid, z_btm, &
430 zh_RSL, zd_RSL, L_MOD, beta, elm, Lc) &
435 INTEGER,
INTENT(in) :: stabilitymethod
436 REAL(kind(1d0)),
INTENT(in) :: psihatm_top
437 REAL(kind(1d0)),
INTENT(in) :: z_top
438 REAL(kind(1d0)),
INTENT(in) :: psihatm_mid
439 REAL(kind(1d0)),
INTENT(in) :: z_mid
440 REAL(kind(1d0)),
INTENT(in) :: z_btm
441 REAL(kind(1d0)),
INTENT(in) :: zh_rsl
442 REAL(kind(1d0)),
INTENT(in) :: zd_rsl
443 REAL(kind(1d0)),
INTENT(in) :: lc
444 REAL(kind(1d0)),
INTENT(in) :: beta
445 REAL(kind(1d0)),
INTENT(in) :: elm
446 REAL(kind(1d0)),
INTENT(in) :: l_mod
447 REAL(kind(1d0)),
INTENT(in) :: cm
448 REAL(kind(1d0)),
INTENT(in) :: c2
451 REAL(kind(1d0)) :: psihatm_btm
454 REAL(kind(1d0)) :: phim_btm
455 REAL(kind(1d0)) :: phim_mid
456 REAL(kind(1d0)) :: slope_top
457 REAL(kind(1d0)) :: slope_btm
458 REAL(kind(1d0)) :: z_plus
459 REAL(kind(1d0)) :: psihatm_plus
461 REAL(kind(1d0)),
PARAMETER :: tol = 0.5
462 REAL(kind(1d0)),
PARAMETER :: kappa = 0.40
463 REAL(kind(1d0)) :: dz_above
465 dz_above = z_mid - z_btm
467 phim_mid =
stab_phi_mom(stabilitymethod, (z_mid - zd_rsl)/l_mod)
468 phim_btm =
stab_phi_mom(stabilitymethod, (z_btm - zd_rsl)/l_mod)
470 psihatm_btm = psihatm_mid + dz_above/2.*phim_mid*(cm*exp(-1.*c2*beta*(z_mid - zd_rsl)/elm)) &
472 psihatm_btm = psihatm_btm + dz_above/2.*phim_btm*(cm*exp(-1.*c2*beta*(z_btm - zd_rsl)/elm)) &
474 IF (dz_above/zd_rsl < 1e-2)
THEN
478 IF (abs(psihatm_btm) > 1e-3)
THEN
480 slope_top = (psihatm_top - psihatm_mid)/(z_top - z_mid)
481 slope_btm = (psihatm_mid - psihatm_btm)/(z_mid - z_btm)
482 IF (abs(slope_btm) > (1 + tol)*abs(slope_top))
THEN
483 z_plus = (z_mid + z_btm)/2
486 psihatm_top, psihatm_mid, &
487 z_top, z_mid, z_plus, &
489 zh_rsl, zd_rsl, l_mod, beta, elm, lc)
492 psihatm_mid, psihatm_plus, &
493 z_mid, z_plus, z_btm, &
495 zh_rsl, zd_rsl, l_mod, beta, elm, lc)
502 psihath_top, psihath_mid, &
503 z_top, z_mid, z_btm, &
505 zh_RSL, zd_RSL, L_MOD, beta, elm, Lc) &
510 INTEGER,
INTENT(in) :: stabilitymethod
511 REAL(kind(1d0)),
INTENT(in) :: psihath_top
512 REAL(kind(1d0)),
INTENT(in) :: z_top
513 REAL(kind(1d0)),
INTENT(in) :: psihath_mid
514 REAL(kind(1d0)),
INTENT(in) :: z_mid
515 REAL(kind(1d0)),
INTENT(in) :: z_btm
516 REAL(kind(1d0)),
INTENT(in) :: zh_rsl
517 REAL(kind(1d0)),
INTENT(in) :: zd_rsl
518 REAL(kind(1d0)),
INTENT(in) :: lc
519 REAL(kind(1d0)),
INTENT(in) :: beta
520 REAL(kind(1d0)),
INTENT(in) :: elm
521 REAL(kind(1d0)),
INTENT(in) :: l_mod
522 REAL(kind(1d0)),
INTENT(in) :: ch
523 REAL(kind(1d0)),
INTENT(in) :: c2h
526 REAL(kind(1d0)) :: psihath_btm
529 REAL(kind(1d0)) :: phih_btm
530 REAL(kind(1d0)) :: phih_mid
532 REAL(kind(1d0)) :: slope_top
533 REAL(kind(1d0)) :: slope_btm
534 REAL(kind(1d0)) :: z_plus
535 REAL(kind(1d0)) :: psihath_plus
537 REAL(kind(1d0)),
PARAMETER :: tol = 0.1
538 REAL(kind(1d0)),
PARAMETER :: kappa = 0.40
539 REAL(kind(1d0)) :: dz_above
541 dz_above = z_mid - z_btm
543 phih_mid =
stab_phi_heat(stabilitymethod, (z_mid - zd_rsl)/l_mod)
544 phih_btm =
stab_phi_heat(stabilitymethod, (z_btm - zd_rsl)/l_mod)
546 psihath_btm = psihath_mid + dz_above/2.*phih_mid*(ch*exp(-1.*c2h*beta*(z_mid - zd_rsl)/elm)) &
548 psihath_btm = psihath_btm + dz_above/2.*phih_btm*(ch*exp(-1.*c2h*beta*(z_btm - zd_rsl)/elm)) &
550 IF (dz_above/zd_rsl < 1e-2)
THEN
554 IF (abs(psihath_btm) > 1e-3)
THEN
556 slope_top = (psihath_top - psihath_mid)/(z_top - z_mid)
557 slope_btm = (psihath_mid - psihath_btm)/(z_mid - z_btm)
558 IF (abs(slope_btm) > (1 + tol)*abs(slope_top))
THEN
559 z_plus = (z_mid + z_btm)/2
562 psihath_top, psihath_mid, &
563 z_top, z_mid, z_plus, &
565 zh_rsl, zd_rsl, l_mod, beta, elm, lc)
568 psihath_mid, psihath_plus, &
569 z_mid, z_plus, z_btm, &
571 zh_rsl, zd_rsl, l_mod, beta, elm, lc)
577 FUNCTION cal_phim_hat(StabilityMethod, z, zh_RSL, L_MOD, beta, lc)
RESULT(phim_hat)
579 INTEGER,
INTENT(in) :: stabilitymethod
580 REAL(kind(1d0)),
INTENT(in) :: z
581 REAL(kind(1d0)),
INTENT(in) :: zh_rsl
582 REAL(kind(1d0)),
INTENT(in) :: l_mod
583 REAL(kind(1d0)),
INTENT(in) :: beta
584 REAL(kind(1d0)),
INTENT(in) :: lc
585 REAL(kind(1d0)) :: phim_hat
586 REAL(kind(1d0)) :: zd_rsl
588 REAL(kind(1d0)) :: elm
589 REAL(kind(1d0)) :: c2
590 REAL(kind(1d0)) :: cm
596 CALL cal_cm(stabilitymethod, zh_rsl, zd_rsl, lc, beta, l_mod, c2, cm)
598 phim_hat = 1 - cm*exp(-1.*c2*beta*(z - zd_rsl)/elm)
602 SUBROUTINE cal_cm(StabilityMethod, zh_RSL, zd_RSL, Lc, beta, L_MOD, c2, cm)
605 INTEGER,
INTENT(in) :: StabilityMethod
607 REAL(KIND(1D0)),
INTENT(in) :: zh_RSL
608 REAL(KIND(1D0)),
INTENT(in) :: zd_RSL
609 REAL(KIND(1D0)),
INTENT(in) :: Lc
610 REAL(KIND(1D0)),
INTENT(in) :: beta
611 REAL(KIND(1D0)),
INTENT(in) :: L_MOD
614 REAL(KIND(1D0)),
INTENT(out) :: c2
615 REAL(KIND(1D0)),
INTENT(out) :: cm
619 REAL(KIND(1D0)) :: phi_hatmZh
620 REAL(KIND(1D0)) :: phim_zh
621 REAL(KIND(1D0)) :: phim_zhdz
622 REAL(KIND(1D0)) :: dphi
625 REAL(KIND(1D0)),
PARAMETER :: kappa = 0.40
626 REAL(KIND(1D0)),
PARAMETER :: dz = 0.1
628 phim_zh =
stab_phi_mom(stabilitymethod, (zh_rsl - zd_rsl)/l_mod)
629 phim_zhdz =
stab_phi_mom(stabilitymethod, (zh_rsl - zd_rsl + dz)/l_mod)
631 dphi = (phim_zhdz - phim_zh)/dz
632 IF (phim_zh /= 0.)
THEN
633 phi_hatmzh = kappa/(2.*beta*phim_zh)
639 IF (phi_hatmzh >= 1.)
THEN
645 c2 = (kappa*(3.-(2.*beta**2.*lc/phim_zh*dphi)))/(2.*beta*phim_zh - kappa)
651 cm = (1.-phi_hatmzh)*exp(c2/2.)
655 SUBROUTINE cal_ch(StabilityMethod, zh_RSL, zd_RSL, Lc, beta, L_MOD, Scc, f, c2h, ch)
658 INTEGER,
INTENT(in) :: StabilityMethod
660 REAL(KIND(1D0)),
INTENT(in) :: zh_RSL
661 REAL(KIND(1D0)),
INTENT(in) :: zd_RSL
662 REAL(KIND(1D0)),
INTENT(in) :: Scc
663 REAL(KIND(1D0)),
INTENT(in) :: f
664 REAL(KIND(1D0)),
INTENT(in) :: Lc
665 REAL(KIND(1D0)),
INTENT(in) :: beta
666 REAL(KIND(1D0)),
INTENT(in) :: L_MOD
669 REAL(KIND(1D0)),
INTENT(out) :: ch
670 REAL(KIND(1D0)),
INTENT(out) :: c2h
673 REAL(KIND(1D0)) :: phih_zh
674 REAL(KIND(1D0)) :: phih_zhdz
675 REAL(KIND(1D0)) :: dphih
676 REAL(KIND(1D0)) :: phi_hathZh
678 REAL(KIND(1D0)),
PARAMETER :: kappa = 0.40
679 REAL(KIND(1D0)),
PARAMETER :: dz = 0.1
681 phih_zh =
stab_phi_heat(stabilitymethod, (zh_rsl - zd_rsl)/l_mod)
682 phih_zhdz =
stab_phi_heat(stabilitymethod, (zh_rsl - zd_rsl + 1.)/l_mod)
684 dphih = phih_zhdz - phih_zh
685 IF (phih_zh /= 0.)
THEN
686 phi_hathzh = kappa*scc/(2.*beta*phih_zh)
691 IF (phi_hathzh >= 1.)
THEN
697 c2h = (kappa*scc*(2.+f - (dphih*2.*beta**2.*lc/phih_zh)))/(2.*beta*phih_zh - kappa*scc)
703 ch = (1.-phi_hathzh)*exp(c2h/2.)
864 REAL(kind(1d0)),
INTENT(in) :: zh_rsl
865 REAL(kind(1d0)),
INTENT(in) :: lc
866 REAL(kind(1d0)),
INTENT(in) :: beta
869 REAL(kind(1d0)) :: zd_rsl
871 zd_rsl = zh_rsl - (beta**2.)*lc
877 FUNCTION cal_z0_rsl(StabilityMethod, zH_RSL, zd_RSL, beta, L_MOD_RSL, psihatm_Zh)
RESULT(z0_RSL)
881 INTEGER,
INTENT(in) :: stabilitymethod
882 REAL(kind(1d0)),
INTENT(in) :: zh_rsl
883 REAL(kind(1d0)),
INTENT(in) :: zd_rsl
884 REAL(kind(1d0)),
INTENT(in) :: l_mod_rsl
886 REAL(kind(1d0)),
INTENT(in) :: beta
887 REAL(kind(1d0)),
INTENT(in) :: psihatm_zh
890 REAL(kind(1d0)) :: z0_rsl
893 REAL(kind(1d0)) :: psimzh, psimz0, z0_rsl_x
894 REAL(kind(1d0)) :: err
897 REAL(kind(1d0)),
PARAMETER :: kappa = 0.40
901 psimzh =
stab_psi_mom(stabilitymethod, (zh_rsl - zd_rsl)/l_mod_rsl)
909 DO WHILE ((err > 0.001) .AND. (it < 10))
911 psimz0 =
stab_psi_mom(stabilitymethod, z0_rsl_x/l_mod_rsl)
912 z0_rsl = (zh_rsl - zd_rsl)*exp(-1.*kappa/beta)*exp(-1.*psimzh + psimz0)*exp(psihatm_zh)
913 err = abs(z0_rsl_x - z0_rsl)
918 z0_rsl = merge(z0_rsl, 1d-2, z0_rsl > 1d-2)
923 StabilityMethod, nz_above, z_array, zh, L_MOD, sfr_surf, FAI, PAI, & !input
924 psihatm_array, psihath_array, zH_RSL, L_MOD_RSL, Lc, beta, zd_RSL, z0_RSL, elm, Scc, fx)
927 INTEGER,
INTENT(in) :: StabilityMethod
928 INTEGER,
INTENT(IN) :: nz_above
929 REAL(KIND(1D0)),
DIMENSION(nz_above),
INTENT(in) :: z_array
930 REAL(KIND(1D0)),
DIMENSION(nz_above),
INTENT(out) :: psihatm_array
931 REAL(KIND(1D0)),
DIMENSION(nz_above),
INTENT(out) :: psihath_array
932 REAL(KIND(1D0)),
INTENT(in) :: zh
933 REAL(KIND(1D0)),
INTENT(in) :: FAI
935 REAL(KIND(1D0)),
INTENT(in) :: PAI
937 REAL(KIND(1D0)),
INTENT(in) :: L_MOD
938 REAL(KIND(1D0)),
DIMENSION(nsurf),
INTENT(in) :: sfr_surf
943 REAL(KIND(1D0)),
INTENT(out) :: L_MOD_RSL
944 REAL(KIND(1D0)),
INTENT(out) :: zH_RSL
947 REAL(KIND(1D0)),
INTENT(out) :: Lc
948 REAL(KIND(1D0)),
INTENT(out) :: beta
949 REAL(KIND(1D0)),
INTENT(out) :: zd_RSL
950 REAL(KIND(1D0)),
INTENT(out) :: z0_RSL
951 REAL(KIND(1D0)),
INTENT(out) :: elm
952 REAL(KIND(1D0)),
INTENT(out) :: Scc
953 REAL(KIND(1D0)),
INTENT(out) :: fx
957 REAL(KIND(1D0)) :: sfr_tr
958 REAL(KIND(1D0)) :: psihatm_Zh
961 REAL(KIND(1D0)) :: lc_over_L
962 REAL(KIND(1D0)) :: z_top, z_mid, z_btm
963 REAL(KIND(1D0)) :: psihatm_top, psihatm_mid, psihatm_btm
964 REAL(KIND(1D0)) :: psihath_top, psihath_mid, psihath_btm
965 REAL(KIND(1D0)) :: cm, ch, c2m, c2h
968 REAL(KIND(1D0)) :: Lc_min
969 REAL(KIND(1D0)) :: dim_bluffbody
971 REAL(KIND(1D0)),
PARAMETER :: cd_tree = 1.2
972 REAL(KIND(1D0)),
PARAMETER :: a_tree = 0.05
974 REAL(KIND(1D0)),
PARAMETER :: beta_N = 0.40
975 REAL(KIND(1D0)),
PARAMETER :: pi = 4.*atan(1.0)
977 REAL(KIND(1D0)),
PARAMETER :: planF_low = 1e-6
978 REAL(KIND(1D0)),
PARAMETER :: kappa = 0.40
980 REAL(KIND(1D0)),
PARAMETER :: r = 0.1
981 REAL(KIND(1D0)),
PARAMETER :: a1 = 4.,
a2 = -0.1,
a3 = 1.5, a4 = -1.
982 REAL(KIND(1D0)),
PARAMETER :: Zh_min = 0.4
983 REAL(KIND(1D0)),
PARAMETER :: porosity_evetr = 0.32
989 zh_rsl = max(zh, zh_min)
1008 lc = (1.-pai)/fai*zh_rsl
1014 dim_bluffbody = zh_rsl*pai/fai
1015 lc_min = dim_bluffbody*pai**(-0.5)/3.
1016 lc = merge(lc, lc_min, lc > lc_min)
1019 lc_over_l = lc/l_mod
1021 IF (lc_over_l > 0)
THEN
1022 lc_over_l = min(2., lc_over_l)
1024 lc_over_l = max(-2., lc_over_l)
1027 l_mod_rsl = lc/lc_over_l
1031 beta =
cal_beta_rsl(stabilitymethod, pai, sfr_tr, lc_over_l)
1034 scc = 0.5 + 0.3*tanh(2.*lc_over_l)
1035 fx = 0.5*((1.+4.*r*scc)**0.5) - 0.5
1043 CALL cal_ch(stabilitymethod, zh_rsl, zd_rsl, lc, beta, l_mod_rsl, scc, fx, c2h, ch)
1044 CALL cal_cm(stabilitymethod, zh_rsl, zd_rsl, lc, beta, l_mod_rsl, c2m, cm)
1049 psihatm_array(nz_above) = 0
1050 psihatm_array(nz_above - 1) = 0
1054 psihath_array(nz_above) = 0
1055 psihath_array(nz_above - 1) = 0
1057 DO iz = nz_above, 3, -1
1059 z_mid = z_array(iz - 1)
1060 z_btm = z_array(iz - 2)
1064 psihatm_top, psihatm_mid, &
1065 z_top, z_mid, z_btm, &
1067 zh_rsl, zd_rsl, l_mod, beta, elm, lc)
1068 psihatm_array(iz - 2) = psihatm_btm
1069 psihatm_top = psihatm_mid
1070 psihatm_mid = psihatm_btm
1074 psihath_top, psihath_mid, &
1075 z_top, z_mid, z_btm, &
1077 zh_rsl, zd_rsl, l_mod, beta, elm, lc)
1078 psihath_top = psihath_mid
1079 psihath_mid = psihath_btm
1080 psihath_array(iz - 2) = psihath_btm
1085 psihatm_zh = psihatm_array(1)
1086 z0_rsl =
cal_z0_rsl(stabilitymethod, zh_rsl, zd_rsl, beta, l_mod_rsl, psihatm_zh)
1090 FUNCTION cal_beta_rsl(StabilityMethod, PAI, sfr_tr, lc_over_L)
RESULT(beta)
1095 INTEGER,
INTENT(in) :: stabilitymethod
1096 REAL(kind(1d0)),
INTENT(in) :: pai
1097 REAL(kind(1d0)),
INTENT(in) :: sfr_tr
1098 REAL(kind(1d0)),
INTENT(in) :: lc_over_l
1101 REAL(kind(1d0)) :: beta
1104 REAL(kind(1d0)) :: betahf
1105 REAL(kind(1d0)) :: betanl
1107 REAL(kind(1d0)),
PARAMETER :: kappa = 0.4
1108 REAL(kind(1d0)),
PARAMETER ::
a1 = 4.,
a2 = -0.1,
a3 = 1.5, a4 = -1.
1112 REAL(kind(1d0)) :: betan2
1120 betan2 = 0.30*sfr_tr/pai + (pai - sfr_tr)/pai*0.4
1125 betahf =
cal_beta_lc(stabilitymethod, betan2, lc_over_l)
1127 betanl =
cal_beta_lc(stabilitymethod, kappa/2., lc_over_l)
1129 IF (lc_over_l >
a2)
THEN
1132 beta = betanl + ((betahf - betanl)/(1.+
a1*abs(lc_over_l -
a2)**
a3))
1135 IF (beta > 0.5)
THEN
1141 FUNCTION cal_beta_lc(stabilityMethod, beta0, lc_over_l)
RESULT(beta_x)
1146 INTEGER,
INTENT(in) :: stabilitymethod
1147 REAL(kind(1d0)),
INTENT(in) :: beta0
1148 REAL(kind(1d0)),
INTENT(in) :: lc_over_l
1149 REAL(kind(1d0)) :: beta_x
1151 REAL(kind(1d0)) :: phim, err, beta_x0
1161 DO WHILE ((err > 0.001) .AND. (it < 20))
1162 beta_x0 = beta0/phim
1163 phim =
stab_phi_heat(stabilitymethod, (beta_x0**2)*lc_over_l)
1167 err = abs(beta_x - beta_x0)
integer, parameter bldgsurf
integer, parameter conifsurf
integer, parameter ncolumnsdataoutrsl
integer, parameter decidsurf
real(kind(1d0)) function stab_phi_heat(StabilityMethod, ZL)
real(kind(1d0)) function stab_phi_mom(StabilityMethod, ZL)
subroutine cal_stab(StabilityMethod, zzd, z0m, zdm, avU1, Temp_C, QH_init, avdens, avcp, L_MOD, TStar, UStar, zL)
real(kind(1d0)), parameter neut_limit
real(kind(1d0)) function stab_psi_heat(StabilityMethod, ZL)
real(kind(1d0)) function stab_psi_mom(StabilityMethod, ZL)
real(kind(1d0)) function rh2qa(RH_dec, pres_hPa, Ta_degC)
real(kind(1d0)) function qa2rh(qa_gkg, pres_hPa, Ta_degC)
subroutine suews_cal_roughnessparameters(RoughLenMomMethod, sfr_surf, bldgH, EveTreeH, DecTreeH, porosity_dectr, FAIBldg, FAIEveTree, FAIDecTree, z0m_in, zdm_in, Z, FAI, PAI, Zh, z0m, zdm, ZZD)
real(kind(1d0)) function cal_beta_rsl(StabilityMethod, PAI, sfr_tr, lc_over_L)
real(kind(1d0)) function interp_z(z_x, z, v)
recursive real(kind(1d0)) function cal_psim_hat(StabilityMethod, psihatm_top, psihatm_mid, z_top, z_mid, z_btm, cm, c2, zh_RSL, zd_RSL, L_MOD, beta, elm, Lc)
subroutine cal_ch(StabilityMethod, zh_RSL, zd_RSL, Lc, beta, L_MOD, Scc, f, c2h, ch)
real(kind(1d0)) function cal_zd_rsl(zh_RSL, beta, Lc)
subroutine rsl_cal_prms(StabilityMethod, nz_above, z_array, zh, L_MOD, sfr_surf, FAI, PAI, psihatm_array, psihath_array, zH_RSL, L_MOD_RSL, Lc, beta, zd_RSL, z0_RSL, elm, Scc, fx)
subroutine cal_cm(StabilityMethod, zh_RSL, zd_RSL, Lc, beta, L_MOD, c2, cm)
recursive real(kind(1d0)) function cal_psih_hat(StabilityMethod, psihath_top, psihath_mid, z_top, z_mid, z_btm, ch, c2h, zh_RSL, zd_RSL, L_MOD, beta, elm, Lc)
real(kind(1d0)) function cal_z0_rsl(StabilityMethod, zH_RSL, zd_RSL, beta, L_MOD_RSL, psihatm_Zh)
subroutine rslprofile(DiagMethod, Zh, z0m, zdm, z0v, L_MOD, sfr_surf, FAI, PAI, StabilityMethod, RA_h, 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 cal_beta_lc(stabilityMethod, beta0, lc_over_l)
real(kind(1d0)) function cal_elm_rsl(beta, Lc)
real(kind(1d0)) function cal_phim_hat(StabilityMethod, z, zh_RSL, L_MOD, beta, lc)