Sources/phylmd/coef_diff_turb.f

module coef_diff_turb_m

  implicit none

contains

  subroutine coef_diff_turb(dtime, nsrf, ni, paprs, pplay, u, v, q, t, ts, &
       cdragm, zgeop, coefm, coefh, q2)

    ! Computes coefficients for turbulent diffusion in the atmosphere.

    USE clesphys, ONLY: ok_kzmin
    use coefkz_m, only: coefkz
    use coefkzmin_m, only: coefkzmin
    use coefkz2_m, only: coefkz2
    USE conf_phys_m, ONLY: iflag_pbl
    USE dimphy, ONLY: klev
    use nr_util, only: assert
    USE suphec_m, ONLY: rd, rg, rkappa
    use ustarhb_m, only: ustarhb
    use yamada4_m, only: yamada4

    REAL, INTENT(IN):: dtime ! interval du temps (secondes)
    INTEGER, INTENT(IN):: nsrf
    INTEGER, INTENT(IN):: ni(:) ! (knon)
    REAL, INTENT(IN):: paprs(:, :) ! (knon, klev + 1)
    REAL, INTENT(IN):: pplay(:, :) ! (knon, klev)
    REAL, INTENT(IN), dimension(:, :):: u, v, q, t ! (knon, klev)
    REAL, INTENT(IN):: ts(:), cdragm(:) ! (knon)
    REAL, INTENT(IN):: zgeop(:, :) ! (knon, klev)
    REAL, intent(out):: coefm(:, 2:) ! (knon, 2:klev) coefficient, vitesse

    real, intent(out):: coefh(:, 2:) ! (knon, 2:klev) 
    ! coefficient, chaleur et humidité

    real, intent(inout):: q2(:, :) ! (knon, klev + 1)

    ! Local:
    REAL coefm0(size(ni), 2:klev), coefh0(size(ni), 2:klev) ! (knon, 2:klev)
    REAL zlay(size(ni), klev), teta(size(ni), klev) ! (knon, klev)
    real zlev(size(ni), klev + 1) ! (knon, klev + 1)
    REAL ustar(size(ni)) ! (knon)
    integer k

    !-------------------------------------------------------------------------

    call assert(size(ni) == [size(paprs, 1), size(pplay, 1), size(u, 1), &
         size(v, 1), size(q, 1), size(t, 1), size(ts), size(cdragm), &
         size(zgeop, 1), size(coefm, 1), size(coefh, 1), size(q2, 1)], &
         "coef_diff_turb knon")
    
    CALL coefkz(nsrf, paprs, pplay, ts, u, v, t, q, zgeop, coefm, coefh)

    IF (iflag_pbl == 1) THEN
       CALL coefkz2(nsrf, paprs, pplay, t, coefm0, coefh0)
       coefm = max(coefm, coefm0)
       coefh = max(coefh, coefh0)
    END IF

    IF (ok_kzmin) THEN
       ! Calcul d'une diffusion minimale pour les conditions tres stables
       CALL coefkzmin(paprs, pplay, u, v, t, q, cdragm, coefh0)
       coefm0 = coefh0
       coefm = max(coefm, coefm0)
       coefh = max(coefh, coefh0)
    END IF

    IF (iflag_pbl >= 6) THEN
       ! Mellor et Yamada adapt\'e \`a Mars, Richard Fournier et
       ! Fr\'ed\'eric Hourdin
       zlay(:, 1) = rd * t(:, 1) / (0.5 * (paprs(:, 1) &
            + pplay(:, 1))) * (paprs(:, 1) - pplay(:, 1)) / rg

       DO k = 2, klev
          zlay(:, k) = zlay(:, k-1) + rd * 0.5 * (t(:, k-1) + t(:, k)) &
               / paprs(:, k) * (pplay(:, k-1) - pplay(:, k)) / rg
       END DO

       DO k = 1, klev
          teta(:, k) = t(:, k) * (paprs(:, 1) / pplay(:, k))**rkappa &
               * (1. + 0.61 * q(:, k))
       END DO

       zlev(:, 1) = 0.
       zlev(:, klev + 1) = 2. * zlay(:, klev) - zlay(:, klev - 1)

       DO k = 2, klev
          zlev(:, k) = 0.5 * (zlay(:, k) + zlay(:, k-1))
       END DO

       ustar = ustarhb(u(:, 1), v(:, 1), cdragm)
       CALL yamada4(dtime, zlev, zlay, u, v, teta, q2, coefm, coefh, ustar)
    END IF

  end subroutine coef_diff_turb

end module coef_diff_turb_m
1	guez	250	module coef_diff_turb_m
2
3			implicit none
4
5			contains
6
7	guez	251	subroutine coef_diff_turb(dtime, nsrf, ni, paprs, pplay, u, v, q, t, ts, &
8			cdragm, zgeop, coefm, coefh, q2)
9	guez	250
10	guez	251	! Computes coefficients for turbulent diffusion in the atmosphere.
11
12	guez	250	USE clesphys, ONLY: ok_kzmin
13			use coefkz_m, only: coefkz
14			use coefkzmin_m, only: coefkzmin
15			use coefkz2_m, only: coefkz2
16			USE conf_phys_m, ONLY: iflag_pbl
17	guez	251	USE dimphy, ONLY: klev
18			use nr_util, only: assert
19	guez	250	USE suphec_m, ONLY: rd, rg, rkappa
20			use ustarhb_m, only: ustarhb
21			use yamada4_m, only: yamada4
22
23			REAL, INTENT(IN):: dtime ! interval du temps (secondes)
24			INTEGER, INTENT(IN):: nsrf
25			INTEGER, INTENT(IN):: ni(:) ! (knon)
26	guez	251	REAL, INTENT(IN):: paprs(:, :) ! (knon, klev + 1)
27			REAL, INTENT(IN):: pplay(:, :) ! (knon, klev)
28			REAL, INTENT(IN), dimension(:, :):: u, v, q, t ! (knon, klev)
29			REAL, INTENT(IN):: ts(:), cdragm(:) ! (knon)
30	guez	250	REAL, INTENT(IN):: zgeop(:, :) ! (knon, klev)
31	guez	251	REAL, intent(out):: coefm(:, 2:) ! (knon, 2:klev) coefficient, vitesse
32	guez	250
33	guez	251	real, intent(out):: coefh(:, 2:) ! (knon, 2:klev)
34	guez	250	! coefficient, chaleur et humidité
35
36	guez	251	real, intent(inout):: q2(:, :) ! (knon, klev + 1)
37	guez	250
38			! Local:
39	guez	251	REAL coefm0(size(ni), 2:klev), coefh0(size(ni), 2:klev) ! (knon, 2:klev)
40			REAL zlay(size(ni), klev), teta(size(ni), klev) ! (knon, klev)
41			real zlev(size(ni), klev + 1) ! (knon, klev + 1)
42			REAL ustar(size(ni)) ! (knon)
43			integer k
44	guez	250
45			!-------------------------------------------------------------------------
46
47	guez	251	call assert(size(ni) == [size(paprs, 1), size(pplay, 1), size(u, 1), &
48			size(v, 1), size(q, 1), size(t, 1), size(ts), size(cdragm), &
49			size(zgeop, 1), size(coefm, 1), size(coefh, 1), size(q2, 1)], &
50			"coef_diff_turb knon")
51
52			CALL coefkz(nsrf, paprs, pplay, ts, u, v, t, q, zgeop, coefm, coefh)
53	guez	250
54			IF (iflag_pbl == 1) THEN
55	guez	251	CALL coefkz2(nsrf, paprs, pplay, t, coefm0, coefh0)
56			coefm = max(coefm, coefm0)
57			coefh = max(coefh, coefh0)
58	guez	250	END IF
59
60			IF (ok_kzmin) THEN
61			! Calcul d'une diffusion minimale pour les conditions tres stables
62	guez	251	CALL coefkzmin(paprs, pplay, u, v, t, q, cdragm, coefh0)
63			coefm0 = coefh0
64			coefm = max(coefm, coefm0)
65			coefh = max(coefh, coefh0)
66	guez	250	END IF
67
68			IF (iflag_pbl >= 6) THEN
69			! Mellor et Yamada adapt\'e \`a Mars, Richard Fournier et
70			! Fr\'ed\'eric Hourdin
71	guez	251	zlay(:, 1) = rd * t(:, 1) / (0.5 * (paprs(:, 1) &
72			+ pplay(:, 1))) * (paprs(:, 1) - pplay(:, 1)) / rg
73	guez	250
74			DO k = 2, klev
75	guez	251	zlay(:, k) = zlay(:, k-1) + rd * 0.5 * (t(:, k-1) + t(:, k)) &
76			/ paprs(:, k) * (pplay(:, k-1) - pplay(:, k)) / rg
77	guez	250	END DO
78
79			DO k = 1, klev
80	guez	251	teta(:, k) = t(:, k) * (paprs(:, 1) / pplay(:, k))**rkappa &
81			* (1. + 0.61 * q(:, k))
82	guez	250	END DO
83
84	guez	251	zlev(:, 1) = 0.
85			zlev(:, klev + 1) = 2. * zlay(:, klev) - zlay(:, klev - 1)
86	guez	250
87			DO k = 2, klev
88	guez	251	zlev(:, k) = 0.5 * (zlay(:, k) + zlay(:, k-1))
89	guez	250	END DO
90
91	guez	251	ustar = ustarhb(u(:, 1), v(:, 1), cdragm)
92			CALL yamada4(dtime, zlev, zlay, u, v, teta, q2, coefm, coefh, ustar)
93	guez	250	END IF
94
95			end subroutine coef_diff_turb
96
97			end module coef_diff_turb_m