trunk/phylmd/coefkzmin.f

module coefkzmin_m

  IMPLICIT NONE

contains

  SUBROUTINE coefkzmin(ypaprs, ypplay, yu, yv, yt, yq, ycoefm, kn)

    ! From LMDZ4/libf/phylmd/coefkzmin.F, version 1.1.1.1, 2004/05/19 12:53:08

    ! Entrées modifiées en attendant une version où les zlev et zlay
    ! soient disponibles.

    USE dimphy, ONLY: klev
    USE suphec_m, ONLY: rd, rg, rkappa

    REAL, intent(in):: ypaprs(:, :) ! (knon, klev+1)
    REAL, intent(in):: ypplay(:, :) ! (knon, klev)
    REAL, intent(in):: yu(:, :), yv(:, :) ! (knon, klev) wind, in m s-1
    REAL, intent(in):: yt(:, :) ! (knon, klev) temperature, in K
    REAL, intent(in):: yq(:, :) ! (knon, klev)
    REAL, intent(in):: ycoefm(:) ! (knon) drag coefficient

    REAL, intent(out):: kn(:, 2:) ! (knon, 2:klev) coefficient de
    ! diffusion turbulente de la quantité de mouvement et des
    ! scalaires (au bas de chaque couche) (en sortie : la valeur à la
    ! fin du pas de temps), m2 s-1

    ! Local:

    integer knon
    real ustar(size(ypaprs, 1)) ! (knon) u*
    real zlay(size(ypaprs, 1), klev) ! (knon, klev) in m
    integer i, k
    real pblhmin(size(ypaprs, 1)) ! (knon)
    real, parameter:: coriol = 1e-4

    REAL zlev(size(ypaprs, 1), 2: klev) ! (knon, 2: klev)
    ! altitude at level (interface between layer with same index), in m

    REAL teta(size(ypaprs, 1), klev) ! (knon, klev)
    ! température potentielle au centre de chaque couche (la valeur au
    ! debut du pas de temps)

    real, PARAMETER:: kap = 0.4

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

    knon = size(ypaprs, 1)
    
    ! Debut de la partie qui doit etre incluse a terme dans pbl_surface.

    do i = 1, knon
       zlay(i, 1) = RD * yt(i, 1) * 2 / (ypaprs(i, 1) + ypplay(i, 1)) &
            * (ypaprs(i, 1) - ypplay(i, 1)) / RG
    enddo

    do k = 2, klev
       do i = 1, knon
          zlay(i, k) = zlay(i, k-1) + RD * 0.5 * (yt(i, k - 1) + yt(i, k)) &
               / ypaprs(i, k) * (ypplay(i, k - 1) - ypplay(i, k)) / RG
       enddo
    enddo

    do k=1, klev
       do i = 1, knon
          ! Attention : on passe la temperature potentielle virtuelle
          ! pour le calcul de K.
          teta(i, k) = yt(i, k) * (ypaprs(i, 1) / ypplay(i, k))**rkappa &
               * (1. + 0.61 * yq(i, k))
       enddo
    enddo

    forall (k = 2: klev) zlev(:, k) = 0.5 * (zlay(:, k) + zlay(:, k-1))
    ustar = SQRT(ycoefm * (yu(:, 1)**2 + yv(:, 1)**2))

    ! Fin de la partie qui doit être incluse à terme dans pbl_surface

    ! Cette routine est ecrite pour avoir en entree ustar, teta et zlev
    ! Ici, on a inclus le calcul de ces trois variables dans la routine
    ! coefkzmin en attendant une nouvelle version de la couche limite
    ! ou ces variables seront disponibles.

    ! Debut de la routine coefkzmin proprement dite

    pblhmin = 0.07 * ustar / coriol

    do k = 2, klev
       do i = 1, knon
          if (teta(i, 2) > teta(i, 1)) then
             kn(i, k) = kap * zlev(i, k) * ustar(i) &
                  * (max(1. - zlev(i, k) / pblhmin(i), 0.))**2
          else
             kn(i, k) = 0. ! min n'est utilisé que pour les SL stables
          endif
       enddo
    enddo

  end SUBROUTINE coefkzmin

end module coefkzmin_m
1	guez	47	module coefkzmin_m
2	guez	3
3	guez	47	IMPLICIT NONE
4	guez	3
5	guez	47	contains
6	guez	3
7	guez	251	SUBROUTINE coefkzmin(ypaprs, ypplay, yu, yv, yt, yq, ycoefm, kn)
8	guez	3
9	guez	251	! From LMDZ4/libf/phylmd/coefkzmin.F, version 1.1.1.1, 2004/05/19 12:53:08
10	guez	3
11	guez	47	! Entrées modifiées en attendant une version où les zlev et zlay
12			! soient disponibles.
13	guez	3
14	guez	251	USE dimphy, ONLY: klev
15	guez	47	USE suphec_m, ONLY: rd, rg, rkappa
16	guez	3
17	guez	251	REAL, intent(in):: ypaprs(:, :) ! (knon, klev+1)
18			REAL, intent(in):: ypplay(:, :) ! (knon, klev)
19			REAL, intent(in):: yu(:, :), yv(:, :) ! (knon, klev) wind, in m s-1
20			REAL, intent(in):: yt(:, :) ! (knon, klev) temperature, in K
21			REAL, intent(in):: yq(:, :) ! (knon, klev)
22	guez	239	REAL, intent(in):: ycoefm(:) ! (knon) drag coefficient
23	guez	3
24	guez	239	REAL, intent(out):: kn(:, 2:) ! (knon, 2:klev) coefficient de
25			! diffusion turbulente de la quantité de mouvement et des
26			! scalaires (au bas de chaque couche) (en sortie : la valeur à la
27			! fin du pas de temps), m2 s-1
28	guez	3
29	guez	47	! Local:
30	guez	3
31	guez	251	integer knon
32			real ustar(size(ypaprs, 1)) ! (knon) u*
33			real zlay(size(ypaprs, 1), klev) ! (knon, klev) in m
34	guez	47	integer i, k
35	guez	251	real pblhmin(size(ypaprs, 1)) ! (knon)
36	guez	47	real, parameter:: coriol = 1e-4
37	guez	3
38	guez	251	REAL zlev(size(ypaprs, 1), 2: klev) ! (knon, 2: klev)
39	guez	47	! altitude at level (interface between layer with same index), in m
40	guez	3
41	guez	251	REAL teta(size(ypaprs, 1), klev) ! (knon, klev)
42	guez	47	! température potentielle au centre de chaque couche (la valeur au
43			! debut du pas de temps)
44	guez	3
45	guez	47	real, PARAMETER:: kap = 0.4
46	guez	3
47	guez	47	!---------------------------------------------------------------------
48	guez	3
49	guez	251	knon = size(ypaprs, 1)
50
51	guez	267	! Debut de la partie qui doit etre incluse a terme dans pbl_surface.
52	guez	3
53	guez	239	do i = 1, knon
54	guez	47	zlay(i, 1) = RD * yt(i, 1) * 2 / (ypaprs(i, 1) + ypplay(i, 1)) &
55			* (ypaprs(i, 1) - ypplay(i, 1)) / RG
56			enddo
57	guez	3
58	guez	47	do k = 2, klev
59	guez	239	do i = 1, knon
60	guez	47	zlay(i, k) = zlay(i, k-1) + RD * 0.5 * (yt(i, k - 1) + yt(i, k)) &
61			/ ypaprs(i, k) * (ypplay(i, k - 1) - ypplay(i, k)) / RG
62			enddo
63			enddo
64	guez	3
65	guez	47	do k=1, klev
66	guez	239	do i = 1, knon
67	guez	47	! Attention : on passe la temperature potentielle virtuelle
68			! pour le calcul de K.
69			teta(i, k) = yt(i, k) * (ypaprs(i, 1) / ypplay(i, k))**rkappa &
70			* (1. + 0.61 * yq(i, k))
71			enddo
72			enddo
73	guez	3
74	guez	47	forall (k = 2: klev) zlev(:, k) = 0.5 * (zlay(:, k) + zlay(:, k-1))
75	guez	251	ustar = SQRT(ycoefm * (yu(:, 1)2 + yv(:, 1)2))
76	guez	3
77	guez	267	! Fin de la partie qui doit être incluse à terme dans pbl_surface
78	guez	3
79	guez	47	! Cette routine est ecrite pour avoir en entree ustar, teta et zlev
80			! Ici, on a inclus le calcul de ces trois variables dans la routine
81			! coefkzmin en attendant une nouvelle version de la couche limite
82			! ou ces variables seront disponibles.
83	guez	3
84	guez	47	! Debut de la routine coefkzmin proprement dite
85
86			pblhmin = 0.07 * ustar / coriol
87
88			do k = 2, klev
89	guez	239	do i = 1, knon
90	guez	47	if (teta(i, 2) > teta(i, 1)) then
91			kn(i, k) = kap * zlev(i, k) * ustar(i) &
92			* (max(1. - zlev(i, k) / pblhmin(i), 0.))**2
93			else
94			kn(i, k) = 0. ! min n'est utilisé que pour les SL stables
95			endif
96			enddo
97			enddo
98
99			end SUBROUTINE coefkzmin
100
101			end module coefkzmin_m