source: LMDZ4/trunk/libf/phylmd/stdlevvar.F90 @ 557

Last change on this file since 557 was 539, checked in by lmdzadmin, 20 years ago

Petit bug sur le signe de la valeur seuil du gradient tair-tsol
IM

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 8.9 KB
Line 
1!
2! $Header$
3!
4      SUBROUTINE stdlevvar(klon, knon, nsrf, zxli, &
5 &                         u1, v1, t1, q1, z1, &
6 &                         ts1, qsurf, rugos, psol, pat1, &
7 &                         t_2m, q_2m, u_10m)
8      IMPLICIT NONE
9!-------------------------------------------------------------------------
10!
11! Objet : calcul de la temperature et l'humidite relative a 2m et du
12!         module du vent a 10m a partir des relations de Dyer-Businger et
13!         des equations de Louis.
14!
15! Reference : Hess, Colman et McAvaney (1995)       
16!
17! I. Musat, 01.07.2002
18!-------------------------------------------------------------------------
19!
20! klon----input-I- dimension de la grille physique (= nb_pts_latitude X nb_pts_longitude)
21! knon----input-I- nombre de points pour un type de surface
22! nsrf----input-I- indice pour le type de surface; voir indicesol.inc
23! zxli----input-L- TRUE si calcul des cdrags selon Laurent Li
24! u1------input-R- vent zonal au 1er niveau du modele
25! v1------input-R- vent meridien au 1er niveau du modele
26! t1------input-R- temperature de l'air au 1er niveau du modele
27! q1------input-R- humidite relative au 1er niveau du modele
28! z1------input-R- geopotentiel au 1er niveau du modele
29! ts1-----input-R- temperature de l'air a la surface
30! qsurf---input-R- humidite relative a la surface
31! rugos---input-R- rugosite
32! psol----input-R- pression au sol
33! pat1----input-R- pression au 1er niveau du modele
34!
35! t_2m---output-R- temperature de l'air a 2m
36! q_2m---output-R- humidite relative a 2m
37! u_10m--output-R- vitesse du vent a 10m
38!
39      INTEGER, intent(in) :: klon, knon, nsrf
40      LOGICAL, intent(in) :: zxli
41      REAL, dimension(klon), intent(in) :: u1, v1, t1, q1, z1, ts1
42      REAL, dimension(klon), intent(in) :: qsurf, rugos
43      REAL, dimension(klon), intent(in) :: psol, pat1
44!
45      REAL, dimension(klon), intent(out) :: t_2m, q_2m, u_10m
46!-------------------------------------------------------------------------
47#include "YOMCST.inc"
48!IM PLUS
49#include "YOETHF.inc"
50!
51! Quelques constantes et options:
52!
53! RKAR : constante de von Karman
54      REAL, PARAMETER :: RKAR=0.40
55! niter : nombre iterations calcul "corrector"
56!     INTEGER, parameter :: niter=6, ncon=niter-1
57      INTEGER, parameter :: niter=2, ncon=niter-1
58!
59! Variables locales
60      INTEGER :: i, n
61      REAL :: zref
62      REAL, dimension(klon) :: speed
63! tpot : temperature potentielle
64      REAL, dimension(klon) :: tpot
65      REAL, dimension(klon) :: zri1, cdran
66      REAL, dimension(klon) :: cdram, cdrah
67! ri1 : nb. de Richardson entre la surface --> la 1ere couche
68      REAL, dimension(klon) :: ri1
69      REAL, dimension(klon) :: ustar, testar, qstar
70      REAL, dimension(klon) :: zdte, zdq   
71! lmon : longueur de Monin-Obukhov selon Hess, Colman and McAvaney
72      DOUBLE PRECISION, dimension(klon) :: lmon
73      DOUBLE PRECISION, parameter :: eps=1.0D-20
74      REAL, dimension(klon) :: delu, delte, delq
75      REAL, dimension(klon) :: u_zref, te_zref, q_zref 
76      REAL, dimension(klon) :: temp, pref
77      LOGICAL :: okri
78      REAL, dimension(klon) :: u_zref_p, te_zref_p, temp_p, q_zref_p
79!convertgence
80      REAL, dimension(klon) :: te_zref_con, q_zref_con
81      REAL, dimension(klon) :: u_zref_c, te_zref_c, temp_c, q_zref_c
82      REAL, dimension(klon) :: ok_pred, ok_corr
83!     REAL, dimension(klon) :: conv_te, conv_q
84!-------------------------------------------------------------------------
85      DO i=1, knon
86       speed(i)=SQRT(u1(i)**2+v1(i)**2)
87       ri1(i) = 0.0
88      ENDDO
89!
90      okri=.FALSE.
91      CALL coefcdrag(klon, knon, nsrf, zxli, &
92 &                   speed, t1, q1, z1, psol, &
93 &                   ts1, qsurf, rugos, okri, ri1,  &         
94 &                   cdram, cdrah, cdran, zri1, pref)           
95!
96!---------Star variables----------------------------------------------------
97!
98      DO i = 1, knon
99        ri1(i) = zri1(i)
100        tpot(i) = t1(i)* (psol(i)/pat1(i))**RKAPPA
101        ustar(i) = sqrt(cdram(i) * speed(i) * speed(i))
102        zdte(i) = tpot(i) - ts1(i)
103!       PRINT*,'AVANT i,zdte',i,zdte(i)
104!IM cf FH : on prend le max : pour eviter le plantage sur SUN
105!IM BUG BUG BUG       zdte(i) = max(zdte(i),1.e-10)
106        zdte(i) = sign(max(abs(zdte(i)),1.e-10),zdte(i))
107!       PRINT*,'APRES i,zdte',i,zdte(i)
108!
109        zdq(i) = max(q1(i),0.0) - max(qsurf(i),0.0)
110!
111        testar(i) = (cdrah(i) * zdte(i) * speed(i))/ustar(i)
112        qstar(i) = (cdrah(i) * zdq(i) * speed(i))/ustar(i)
113        lmon(i) = (ustar(i) * ustar(i) * tpot(i))/ &
114 &                (RKAR * RG * testar(i))
115      ENDDO
116!
117!----------First aproximation of variables at zref --------------------------
118      zref = 2.0
119      CALL screenp(klon, knon, nsrf, speed, tpot, q1, &
120 &                 ts1, qsurf, rugos, lmon, &
121 &                 ustar, testar, qstar, zref, &
122 &                 delu, delte, delq)
123!
124      DO i = 1, knon
125        u_zref(i) = delu(i)
126        q_zref(i) = max(qsurf(i),0.0) + delq(i)
127        te_zref(i) = ts1(i) + delte(i)
128        temp(i) = te_zref(i) * (psol(i)/pat1(i))**(-RKAPPA)
129        q_zref_p(i) = q_zref(i)
130!       te_zref_p(i) = te_zref(i)
131        temp_p(i) = temp(i)
132      ENDDO
133!
134! Iteration of the variables at the reference level zref : corrector calculation ; see Hess & McAvaney, 1995
135!
136      DO n = 1, niter
137!
138        okri=.TRUE.
139        CALL screenc(klon, knon, nsrf, zxli, &
140 &                   u_zref, temp, q_zref, zref, &
141 &                   ts1, qsurf, rugos, psol, &           
142 &                   ustar, testar, qstar, okri, ri1, &
143 &                   pref, delu, delte, delq)
144!
145        DO i = 1, knon
146          u_zref(i) = delu(i)
147          q_zref(i) = delq(i) + max(qsurf(i),0.0)
148          te_zref(i) = delte(i) + ts1(i)
149!
150! return to normal temperature
151!
152          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
153!         temp(i) = te_zref(i) - (zref* RG)/RCPD/ &
154!                 (1 + RVTMP2 * max(q_zref(i),0.0))
155!
156!IM +++
157!         IF(temp(i).GT.350.) THEN
158!           WRITE(*,*) 'temp(i) GT 350 K !!',i,nsrf,temp(i)
159!         ENDIF
160!IM ---
161!
162        IF(n.EQ.ncon) THEN
163          te_zref_con(i) = te_zref(i)
164          q_zref_con(i) = q_zref(i)
165        ENDIF
166!
167        ENDDO
168!
169      ENDDO
170!
171! verifier le critere de convergence : 0.25% pour te_zref et 5% pour qe_zref
172!
173!       DO i = 1, knon
174!         conv_te(i) = (te_zref(i) - te_zref_con(i))/te_zref_con(i)
175!         conv_q(i) = (q_zref(i) - q_zref_con(i))/q_zref_con(i)
176!IM +++
177!         IF(abs(conv_te(i)).GE.0.0025.AND.abs(conv_q(i)).GE.0.05) THEN
178!           PRINT*,'DIV','i=',i,te_zref_con(i),te_zref(i),conv_te(i), &
179!           q_zref_con(i),q_zref(i),conv_q(i)
180!         ENDIF
181!IM ---
182!       ENDDO
183!
184      DO i = 1, knon
185        q_zref_c(i) = q_zref(i)
186        temp_c(i) = temp(i)
187!
188!       IF(zri1(i).LT.0.) THEN
189!         IF(nsrf.EQ.1) THEN
190!           ok_pred(i)=1.
191!           ok_corr(i)=0.
192!         ELSE
193!           ok_pred(i)=0.
194!           ok_corr(i)=1.
195!         ENDIF
196!       ELSE
197!         ok_pred(i)=0.
198!         ok_corr(i)=1.
199!       ENDIF
200!
201        ok_pred(i)=0.
202        ok_corr(i)=1.
203!
204        t_2m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
205        q_2m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
206!IM +++
207!       IF(n.EQ.niter) THEN
208!       IF(t_2m(i).LT.t1(i).AND.t_2m(i).LT.ts1(i)) THEN
209!         PRINT*,' BAD t2m LT ',i,nsrf,t_2m(i),t1(i),ts1(i)
210!       ELSEIF(t_2m(i).GT.t1(i).AND.t_2m(i).GT.ts1(i)) THEN
211!         PRINT*,' BAD t2m GT ',i,nsrf,t_2m(i),t1(i),ts1(i)
212!       ENDIF
213!       ENDIF
214!IM ---
215      ENDDO
216!
217!
218!----------First aproximation of variables at zref --------------------------
219!
220      zref = 10.0
221      CALL screenp(klon, knon, nsrf, speed, tpot, q1, &
222 &                 ts1, qsurf, rugos, lmon, &
223 &                 ustar, testar, qstar, zref, &
224 &                 delu, delte, delq)
225!
226      DO i = 1, knon
227        u_zref(i) = delu(i)
228        q_zref(i) = max(qsurf(i),0.0) + delq(i)
229        te_zref(i) = ts1(i) + delte(i)
230        temp(i) = te_zref(i) * (psol(i)/pat1(i))**(-RKAPPA)
231!       temp(i) = te_zref(i) - (zref* RG)/RCPD/ &
232!                 (1 + RVTMP2 * max(q_zref(i),0.0))
233        u_zref_p(i) = u_zref(i)
234      ENDDO
235!
236! Iteration of the variables at the reference level zref : corrector ; see Hess & McAvaney, 1995
237!
238      DO n = 1, niter
239!
240        okri=.TRUE.
241        CALL screenc(klon, knon, nsrf, zxli, &
242 &                   u_zref, temp, q_zref, zref, &
243 &                   ts1, qsurf, rugos, psol, &
244 &                   ustar, testar, qstar, okri, ri1, &
245 &                   pref, delu, delte, delq)
246!
247        DO i = 1, knon
248          u_zref(i) = delu(i)
249          q_zref(i) = delq(i) + max(qsurf(i),0.0)
250          te_zref(i) = delte(i) + ts1(i)
251          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
252!         temp(i) = te_zref(i) - (zref* RG)/RCPD/ &
253!                   (1 + RVTMP2 * max(q_zref(i),0.0))
254        ENDDO
255!
256      ENDDO
257!
258      DO i = 1, knon
259        u_zref_c(i) = u_zref(i)
260!
261        u_10m(i) = u_zref_p(i) * ok_pred(i) + u_zref_c(i) * ok_corr(i)
262      ENDDO
263!
264      RETURN
265      END subroutine stdlevvar
Note: See TracBrowser for help on using the repository browser.