source: LMDZ4/branches/unlabeled-1.1.1/libf/phylmd/stdlevvar.F90 @ 525

Last change on this file since 525 was 524, checked in by lmdzadmin, 21 years ago

Initial revision

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 8.7 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!IM cf FH : on prend le max : pour eviter le plantage sur SUN
104        zdte(i) = max(zdte(i),1.e-10)
105        zdq(i) = max(q1(i),0.0) - max(qsurf(i),0.0)
106!
107        testar(i) = (cdrah(i) * zdte(i) * speed(i))/ustar(i)
108        qstar(i) = (cdrah(i) * zdq(i) * speed(i))/ustar(i)
109        lmon(i) = (ustar(i) * ustar(i) * tpot(i))/ &
110 &                (RKAR * RG * testar(i))
111      ENDDO
112!
113!----------First aproximation of variables at zref --------------------------
114      zref = 2.0
115      CALL screenp(klon, knon, nsrf, speed, tpot, q1, &
116 &                 ts1, qsurf, rugos, lmon, &
117 &                 ustar, testar, qstar, zref, &
118 &                 delu, delte, delq)
119!
120      DO i = 1, knon
121        u_zref(i) = delu(i)
122        q_zref(i) = max(qsurf(i),0.0) + delq(i)
123        te_zref(i) = ts1(i) + delte(i)
124        temp(i) = te_zref(i) * (psol(i)/pat1(i))**(-RKAPPA)
125        q_zref_p(i) = q_zref(i)
126!       te_zref_p(i) = te_zref(i)
127        temp_p(i) = temp(i)
128      ENDDO
129!
130! Iteration of the variables at the reference level zref : corrector calculation ; see Hess & McAvaney, 1995
131!
132      DO n = 1, niter
133!
134        okri=.TRUE.
135        CALL screenc(klon, knon, nsrf, zxli, &
136 &                   u_zref, temp, q_zref, zref, &
137 &                   ts1, qsurf, rugos, psol, &           
138 &                   ustar, testar, qstar, okri, ri1, &
139 &                   pref, delu, delte, delq)
140!
141        DO i = 1, knon
142          u_zref(i) = delu(i)
143          q_zref(i) = delq(i) + max(qsurf(i),0.0)
144          te_zref(i) = delte(i) + ts1(i)
145!
146! return to normal temperature
147!
148          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
149!         temp(i) = te_zref(i) - (zref* RG)/RCPD/ &
150!                 (1 + RVTMP2 * max(q_zref(i),0.0))
151!
152!IM +++
153!         IF(temp(i).GT.350.) THEN
154!           WRITE(*,*) 'temp(i) GT 350 K !!',i,nsrf,temp(i)
155!         ENDIF
156!IM ---
157!
158        IF(n.EQ.ncon) THEN
159          te_zref_con(i) = te_zref(i)
160          q_zref_con(i) = q_zref(i)
161        ENDIF
162!
163        ENDDO
164!
165      ENDDO
166!
167! verifier le critere de convergence : 0.25% pour te_zref et 5% pour qe_zref
168!
169!       DO i = 1, knon
170!         conv_te(i) = (te_zref(i) - te_zref_con(i))/te_zref_con(i)
171!         conv_q(i) = (q_zref(i) - q_zref_con(i))/q_zref_con(i)
172!IM +++
173!         IF(abs(conv_te(i)).GE.0.0025.AND.abs(conv_q(i)).GE.0.05) THEN
174!           PRINT*,'DIV','i=',i,te_zref_con(i),te_zref(i),conv_te(i), &
175!           q_zref_con(i),q_zref(i),conv_q(i)
176!         ENDIF
177!IM ---
178!       ENDDO
179!
180      DO i = 1, knon
181        q_zref_c(i) = q_zref(i)
182        temp_c(i) = temp(i)
183!
184!       IF(zri1(i).LT.0.) THEN
185!         IF(nsrf.EQ.1) THEN
186!           ok_pred(i)=1.
187!           ok_corr(i)=0.
188!         ELSE
189!           ok_pred(i)=0.
190!           ok_corr(i)=1.
191!         ENDIF
192!       ELSE
193!         ok_pred(i)=0.
194!         ok_corr(i)=1.
195!       ENDIF
196!
197        ok_pred(i)=0.
198        ok_corr(i)=1.
199!
200        t_2m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
201        q_2m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
202!IM +++
203!       IF(n.EQ.niter) THEN
204!       IF(t_2m(i).LT.t1(i).AND.t_2m(i).LT.ts1(i)) THEN
205!         PRINT*,' BAD t2m LT ',i,nsrf,t_2m(i),t1(i),ts1(i)
206!       ELSEIF(t_2m(i).GT.t1(i).AND.t_2m(i).GT.ts1(i)) THEN
207!         PRINT*,' BAD t2m GT ',i,nsrf,t_2m(i),t1(i),ts1(i)
208!       ENDIF
209!       ENDIF
210!IM ---
211      ENDDO
212!
213!
214!----------First aproximation of variables at zref --------------------------
215!
216      zref = 10.0
217      CALL screenp(klon, knon, nsrf, speed, tpot, q1, &
218 &                 ts1, qsurf, rugos, lmon, &
219 &                 ustar, testar, qstar, zref, &
220 &                 delu, delte, delq)
221!
222      DO i = 1, knon
223        u_zref(i) = delu(i)
224        q_zref(i) = max(qsurf(i),0.0) + delq(i)
225        te_zref(i) = ts1(i) + delte(i)
226        temp(i) = te_zref(i) * (psol(i)/pat1(i))**(-RKAPPA)
227!       temp(i) = te_zref(i) - (zref* RG)/RCPD/ &
228!                 (1 + RVTMP2 * max(q_zref(i),0.0))
229        u_zref_p(i) = u_zref(i)
230      ENDDO
231!
232! Iteration of the variables at the reference level zref : corrector ; see Hess & McAvaney, 1995
233!
234      DO n = 1, niter
235!
236        okri=.TRUE.
237        CALL screenc(klon, knon, nsrf, zxli, &
238 &                   u_zref, temp, q_zref, zref, &
239 &                   ts1, qsurf, rugos, psol, &
240 &                   ustar, testar, qstar, okri, ri1, &
241 &                   pref, delu, delte, delq)
242!
243        DO i = 1, knon
244          u_zref(i) = delu(i)
245          q_zref(i) = delq(i) + max(qsurf(i),0.0)
246          te_zref(i) = delte(i) + ts1(i)
247          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
248!         temp(i) = te_zref(i) - (zref* RG)/RCPD/ &
249!                   (1 + RVTMP2 * max(q_zref(i),0.0))
250        ENDDO
251!
252      ENDDO
253!
254      DO i = 1, knon
255        u_zref_c(i) = u_zref(i)
256!
257        u_10m(i) = u_zref_p(i) * ok_pred(i) + u_zref_c(i) * ok_corr(i)
258      ENDDO
259!
260      RETURN
261      END subroutine stdlevvar
Note: See TracBrowser for help on using the repository browser.