/[lmdze]/trunk/Sources/phylmd/yamada.f
ViewVC logotype

Contents of /trunk/Sources/phylmd/yamada.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 134 - (show annotations)
Wed Apr 29 15:47:56 2015 UTC (9 years, 1 month ago) by guez
File size: 5205 byte(s)
Sources inside, compilation outside.
1 module yamada_m
2
3 IMPLICIT NONE
4
5 contains
6
7 SUBROUTINE yamada(ngrid, g, rconst, plev, temp, zlev, zlay, u, v, teta, q2, &
8 km, kn, ustar, l_mix)
9
10 ! From LMDZ4/libf/phylmd/yamada.F,v 1.1 2004/06/22 11:45:36
11
12 USE dimens_m
13 USE dimphy
14 ! .......................................................................
15 ! .......................................................................
16
17 ! g : g
18 ! zlev : altitude a chaque niveau (interface inferieure de la couche
19 ! de meme indice)
20 ! zlay : altitude au centre de chaque couche
21 ! u,v : vitesse au centre de chaque couche
22 ! (en entree : la valeur au debut du pas de temps)
23 ! teta : temperature potentielle au centre de chaque couche
24 ! (en entree : la valeur au debut du pas de temps)
25 ! q2 : $q^2$ au bas de chaque couche
26 ! (en entree : la valeur au debut du pas de temps)
27 ! (en sortie : la valeur a la fin du pas de temps)
28 ! km : diffusivite turbulente de quantite de mouvement (au bas de chaque
29 ! couche)
30 ! (en sortie : la valeur a la fin du pas de temps)
31 ! kn : diffusivite turbulente des scalaires (au bas de chaque couche)
32 ! (en sortie : la valeur a la fin du pas de temps)
33
34 ! .......................................................................
35 REAL, INTENT (IN) :: g
36 REAL rconst
37 REAL plev(klon, klev+1), temp(klon, klev)
38 REAL ustar(klon), snstable
39 REAL zlev(klon, klev+1)
40 REAL zlay(klon, klev)
41 REAL u(klon, klev)
42 REAL v(klon, klev)
43 REAL teta(klon, klev)
44 REAL q2(klon, klev+1)
45 REAL km(klon, klev+1)
46 REAL kn(klon, klev+1)
47 INTEGER l_mix, ngrid
48
49
50 INTEGER nlay, nlev
51 PARAMETER (nlay=klev)
52 PARAMETER (nlev=klev+1)
53
54 LOGICAL first
55 SAVE first
56 DATA first/.TRUE./
57
58
59 INTEGER ig, k
60
61 REAL ri, zrif, zalpha, zsm
62 REAL rif(klon, klev+1), sm(klon, klev+1), alpha(klon, klev)
63
64 REAL m2(klon, klev+1), dz(klon, klev+1), zq, n2(klon, klev+1)
65 REAL l(klon, klev+1), l0(klon)
66
67 REAL sq(klon), sqz(klon), zz(klon, klev+1)
68 INTEGER iter
69
70 REAL ric, rifc, b1, kap
71 SAVE ric, rifc, b1, kap
72 DATA ric, rifc, b1, kap/0.195, 0.191, 16.6, 0.3/
73
74 IF (0==1 .AND. first) THEN
75 DO ig = 1, 1000
76 ri = (ig-800.)/500.
77 IF (ri<ric) THEN
78 zrif = frif(ri)
79 ELSE
80 zrif = rifc
81 END IF
82 IF (zrif<0.16) THEN
83 zalpha = falpha(zrif)
84 zsm = fsm(zrif)
85 ELSE
86 zalpha = 1.12
87 zsm = 0.085
88 END IF
89 PRINT *, ri, rif, zalpha, zsm
90 END DO
91 first = .FALSE.
92 END IF
93
94 ! Correction d'un bug sauvage a verifier.
95 ! do k=2,nlev
96 DO k = 2, nlay
97 DO ig = 1, ngrid
98 dz(ig, k) = zlay(ig, k) - zlay(ig, k-1)
99 m2(ig, k) = ((u(ig,k)-u(ig,k-1))**2+(v(ig,k)-v(ig, &
100 k-1))**2)/(dz(ig,k)*dz(ig,k))
101 n2(ig, k) = g*2.*(teta(ig,k)-teta(ig,k-1))/(teta(ig,k-1)+teta(ig,k))/ &
102 dz(ig, k)
103 ri = n2(ig, k)/max(m2(ig,k), 1.E-10)
104 IF (ri<ric) THEN
105 rif(ig, k) = frif(ri)
106 ELSE
107 rif(ig, k) = rifc
108 END IF
109 IF (rif(ig,k)<0.16) THEN
110 alpha(ig, k) = falpha(rif(ig,k))
111 sm(ig, k) = fsm(rif(ig,k))
112 ELSE
113 alpha(ig, k) = 1.12
114 sm(ig, k) = 0.085
115 END IF
116 zz(ig, k) = b1*m2(ig, k)*(1.-rif(ig,k))*sm(ig, k)
117 END DO
118 END DO
119
120 ! iterration pour determiner la longueur de melange
121
122 DO ig = 1, ngrid
123 l0(ig) = 100.
124 END DO
125 DO k = 2, klev - 1
126 DO ig = 1, ngrid
127 l(ig, k) = l0(ig)*kap*zlev(ig, k)/(kap*zlev(ig,k)+l0(ig))
128 END DO
129 END DO
130
131 DO iter = 1, 10
132 DO ig = 1, ngrid
133 sq(ig) = 1.E-10
134 sqz(ig) = 1.E-10
135 END DO
136 DO k = 2, klev - 1
137 DO ig = 1, ngrid
138 q2(ig, k) = l(ig, k)**2*zz(ig, k)
139 l(ig, k) = min(l0(ig)*kap*zlev(ig,k)/(kap*zlev(ig, &
140 k)+l0(ig)), 0.5*sqrt(q2(ig,k))/sqrt(max(n2(ig,k),1.E-10)))
141 zq = sqrt(q2(ig,k))
142 sqz(ig) = sqz(ig) + zq*zlev(ig, k)*(zlay(ig,k)-zlay(ig,k-1))
143 sq(ig) = sq(ig) + zq*(zlay(ig,k)-zlay(ig,k-1))
144 END DO
145 END DO
146 DO ig = 1, ngrid
147 l0(ig) = 0.2*sqz(ig)/sq(ig)
148 END DO
149 ! (abd 3 5 2) print*,'ITER=',iter,' L0=',l0
150
151 END DO
152
153 DO k = 2, klev
154 DO ig = 1, ngrid
155 l(ig, k) = min(l0(ig)*kap*zlev(ig,k)/(kap*zlev(ig, &
156 k)+l0(ig)), 0.5*sqrt(q2(ig,k))/sqrt(max(n2(ig,k),1.E-10)))
157 q2(ig, k) = l(ig, k)**2*zz(ig, k)
158 km(ig, k) = l(ig, k)*sqrt(q2(ig,k))*sm(ig, k)
159 kn(ig, k) = km(ig, k)*alpha(ig, k)
160 END DO
161 END DO
162
163 contains
164
165 REAL function frif(ri)
166 real ri
167 frif = 0.6588*(ri+0.1776-sqrt(ri*ri-0.3221*ri+0.03156))
168 end function frif
169
170 REAL function falpha(ri)
171 real ri
172 falpha = 1.318*(0.2231-ri)/(0.2341-ri)
173 end function falpha
174
175 REAL function fsm(ri)
176 real ri
177 fsm = 1.96*(0.1912-ri)*(0.2341-ri)/((1.-ri)*(0.2231-ri))
178 end function fsm
179
180 END SUBROUTINE yamada
181
182 end module yamada_m

  ViewVC Help
Powered by ViewVC 1.1.21