source: LMDZ6/branches/LMDZ_cdrag_LSCE/libf/phylmd/stdlevvar_mod.F90 @ 5304

Last change on this file since 5304 was 4661, checked in by Laurent Fairhead, 14 months ago

Inclusion and merge of Olivier's modifications

  • Property svn:executable set to *
File size: 26.2 KB
Line 
1!
2MODULE stdlevvar_mod
3!
4! This module contains main procedures for calculation
5! of temperature, specific humidity and wind at a reference level
6!
7  USE cdrag_mod
8  USE screenp_mod
9  USE screenc_mod
10  IMPLICIT NONE
11
12CONTAINS
13!
14!****************************************************************************************
15!
16!r original routine svn3623
17!
18      SUBROUTINE stdlevvar(klon, knon, nsrf, zxli, &
19                           u1, v1, t1, q1, z1, &
20                           ts1, qsurf, z0m, z0h, psol, pat1, &
21                           t_2m, q_2m, t_10m, q_10m, u_10m, ustar, s_pblh, prain, tsol)
22      IMPLICIT NONE
23!-------------------------------------------------------------------------
24!
25! Objet : calcul de la temperature et l'humidite relative a 2m et du
26!         module du vent a 10m a partir des relations de Dyer-Businger et
27!         des equations de Louis.
28!
29! Reference : Hess, Colman et McAvaney (1995)       
30!
31! I. Musat, 01.07.2002
32!
33!AM On rajoute en sortie t et q a 10m pr le calcule d'hbtm2 dans clmain
34!
35!-------------------------------------------------------------------------
36!
37! klon----input-I- dimension de la grille physique (= nb_pts_latitude X nb_pts_longitude)
38! knon----input-I- nombre de points pour un type de surface
39! nsrf----input-I- indice pour le type de surface; voir indice_sol_mod.F90
40! zxli----input-L- TRUE si calcul des cdrags selon Laurent Li
41! u1------input-R- vent zonal au 1er niveau du modele
42! v1------input-R- vent meridien au 1er niveau du modele
43! t1------input-R- temperature de l'air au 1er niveau du modele
44! q1------input-R- humidite relative au 1er niveau du modele
45! z1------input-R- geopotentiel au 1er niveau du modele
46! ts1-----input-R- temperature de l'air a la surface
47! qsurf---input-R- humidite relative a la surface
48! z0m, z0h---input-R- rugosite
49! psol----input-R- pression au sol
50! pat1----input-R- pression au 1er niveau du modele
51!
52! t_2m---output-R- temperature de l'air a 2m
53! q_2m---output-R- humidite relative a 2m
54! u_10m--output-R- vitesse du vent a 10m
55!AM
56! t_10m--output-R- temperature de l'air a 10m
57! q_10m--output-R- humidite specifique a 10m
58! ustar--output-R- u*
59!
60      INTEGER, intent(in) :: klon, knon, nsrf
61      LOGICAL, intent(in) :: zxli
62      REAL, dimension(klon), intent(in) :: u1, v1, t1, q1, z1, ts1
63      REAL, dimension(klon), intent(in) :: qsurf
64      REAL, dimension(klon), intent(inout) :: z0m, z0h
65      REAL, dimension(klon), intent(in) :: psol, pat1
66!
67      REAL, dimension(klon), intent(out) :: t_2m, q_2m, ustar
68      REAL, dimension(klon), intent(out) :: u_10m, t_10m, q_10m
69      REAL, DIMENSION(klon), INTENT(INOUT) :: s_pblh
70      REAL, DIMENSION(klon), INTENT(IN) :: prain
71      REAL, DIMENSION(klon), INTENT(IN) :: tsol
72!-------------------------------------------------------------------------
73      include "flux_arp.h"
74      include "YOMCST.h"
75!IM PLUS
76      include "YOETHF.h"
77!
78! Quelques constantes et options:
79!
80! RKAR : constante de von Karman
81      REAL, PARAMETER :: RKAR=0.40
82! niter : nombre iterations calcul "corrector"
83!     INTEGER, parameter :: niter=6, ncon=niter-1
84      INTEGER, parameter :: niter=2, ncon=niter-1
85!
86! Variables locales
87      INTEGER :: i, n
88      REAL :: zref
89      REAL, dimension(klon) :: speed
90! tpot : temperature potentielle
91      REAL, dimension(klon) :: tpot
92      REAL, dimension(klon) :: zri1, cdran
93      REAL, dimension(klon) :: cdram, cdrah
94! ri1 : nb. de Richardson entre la surface --> la 1ere couche
95      REAL, dimension(klon) :: ri1
96      REAL, dimension(klon) :: testar, qstar
97      REAL, dimension(klon) :: zdte, zdq   
98! lmon : longueur de Monin-Obukhov selon Hess, Colman and McAvaney
99      DOUBLE PRECISION, dimension(klon) :: lmon
100      DOUBLE PRECISION, parameter :: eps=1.0D-20
101      REAL, dimension(klon) :: delu, delte, delq
102      REAL, dimension(klon) :: u_zref, te_zref, q_zref 
103      REAL, dimension(klon) :: temp, pref
104      LOGICAL :: okri
105      REAL, dimension(klon) :: u_zref_p, te_zref_p, temp_p, q_zref_p
106!convertgence
107      REAL, dimension(klon) :: te_zref_con, q_zref_con
108      REAL, dimension(klon) :: u_zref_c, te_zref_c, temp_c, q_zref_c
109      REAL, dimension(klon) :: ok_pred, ok_corr, zri_zero
110!     REAL, dimension(klon) :: conv_te, conv_q
111!-------------------------------------------------------------------------
112      DO i=1, knon
113       speed(i)=SQRT(u1(i)**2+v1(i)**2)
114       ri1(i) = 0.0
115      ENDDO
116!
117      okri=.FALSE.
118!      CALL coefcdrag(klon, knon, nsrf, zxli, &
119! &                   speed, t1, q1, z1, psol, &
120! &                   ts1, qsurf, rugos, okri, ri1,  &         
121! &                   cdram, cdrah, cdran, zri1, pref)           
122! Fuxing WANG, 04/03/2015, replace the coefcdrag by the merged version: cdrag
123
124      CALL cdrag(knon, nsrf, &
125 &                   speed, t1, q1, z1, &
126 &                   psol, s_pblh, ts1, qsurf, z0m, z0h, &
127 &                   zri_zero, 0, &
128 &                   cdram, cdrah, zri1, pref, prain, tsol, pat1)
129
130! --- special Dice: on force cdragm ( a defaut de forcer ustar) MPL 05082013
131     IF (ok_prescr_ust) then
132      DO i = 1, knon
133       print *,'cdram avant=',cdram(i)
134       cdram(i) = ust*ust/speed(i)/speed(i)
135       print *,'cdram ust speed apres=',cdram(i),ust,speed
136      ENDDO
137     ENDIF
138!
139!---------Star variables----------------------------------------------------
140!
141      DO i = 1, knon
142        ri1(i) = zri1(i)
143        tpot(i) = t1(i)* (psol(i)/pat1(i))**RKAPPA
144        ustar(i) = sqrt(cdram(i) * speed(i) * speed(i))
145        zdte(i) = tpot(i) - ts1(i)
146        zdq(i) = max(q1(i),0.0) - max(qsurf(i),0.0)
147!
148!
149!IM BUG BUG BUG       zdte(i) = max(zdte(i),1.e-10)
150        zdte(i) = sign(max(abs(zdte(i)),1.e-10),zdte(i))
151!
152        testar(i) = (cdrah(i) * zdte(i) * speed(i))/ustar(i)
153        qstar(i) = (cdrah(i) * zdq(i) * speed(i))/ustar(i)
154        lmon(i) = (ustar(i) * ustar(i) * tpot(i))/ &
155 &                (RKAR * RG * testar(i))
156      ENDDO
157!
158!----------First aproximation of variables at zref --------------------------
159      zref = 2.0
160      CALL screenp(klon, knon, nsrf, speed, tpot, q1, &
161 &                 ts1, qsurf, z0m, lmon, &
162 &                 ustar, testar, qstar, zref, &
163 &                 delu, delte, delq)
164!
165      DO i = 1, knon
166        u_zref(i) = delu(i)
167        q_zref(i) = max(qsurf(i),0.0) + delq(i)
168        te_zref(i) = ts1(i) + delte(i)
169        temp(i) = te_zref(i) * (psol(i)/pat1(i))**(-RKAPPA)
170        q_zref_p(i) = q_zref(i)
171!       te_zref_p(i) = te_zref(i)
172        temp_p(i) = temp(i)
173      ENDDO
174!
175! Iteration of the variables at the reference level zref : corrector calculation ; see Hess & McAvaney, 1995
176!
177      DO n = 1, niter
178!
179        okri=.TRUE.
180        CALL screenc(klon, knon, nsrf, zxli, &
181 &                   u_zref, temp, q_zref, zref, &
182 &                   ts1, qsurf, z0m, z0h, psol, &           
183 &                   ustar, testar, qstar, okri, ri1, &
184 &                   pref, delu, delte, delq, s_pblh ,prain, tsol, pat1)
185!
186        DO i = 1, knon
187          u_zref(i) = delu(i)
188          q_zref(i) = delq(i) + max(qsurf(i),0.0)
189          te_zref(i) = delte(i) + ts1(i)
190!
191! return to normal temperature
192!
193          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
194!         temp(i) = te_zref(i) - (zref* RG)/RCPD/ &
195!                 (1 + RVTMP2 * max(q_zref(i),0.0))
196!
197!IM +++
198!         IF(temp(i).GT.350.) THEN
199!           WRITE(*,*) 'temp(i) GT 350 K !!',i,nsrf,temp(i)
200!         ENDIF
201!IM ---
202!
203        IF(n.EQ.ncon) THEN
204          te_zref_con(i) = te_zref(i)
205          q_zref_con(i) = q_zref(i)
206        ENDIF
207!
208        ENDDO
209!
210      ENDDO
211!
212! verifier le critere de convergence : 0.25% pour te_zref et 5% pour qe_zref
213!
214!       DO i = 1, knon
215!         conv_te(i) = (te_zref(i) - te_zref_con(i))/te_zref_con(i)
216!         conv_q(i) = (q_zref(i) - q_zref_con(i))/q_zref_con(i)
217!IM +++
218!         IF(abs(conv_te(i)).GE.0.0025.AND.abs(conv_q(i)).GE.0.05) THEN
219!           PRINT*,'DIV','i=',i,te_zref_con(i),te_zref(i),conv_te(i), &
220!           q_zref_con(i),q_zref(i),conv_q(i)
221!         ENDIF
222!IM ---
223!       ENDDO
224!
225      DO i = 1, knon
226        q_zref_c(i) = q_zref(i)
227        temp_c(i) = temp(i)
228!
229!       IF(zri1(i).LT.0.) THEN
230!         IF(nsrf.EQ.1) THEN
231!           ok_pred(i)=1.
232!           ok_corr(i)=0.
233!         ELSE
234!           ok_pred(i)=0.
235!           ok_corr(i)=1.
236!         ENDIF
237!       ELSE
238!         ok_pred(i)=0.
239!         ok_corr(i)=1.
240!       ENDIF
241!
242        ok_pred(i)=0.
243        ok_corr(i)=1.
244!
245        t_2m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
246        q_2m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
247!IM +++
248!       IF(n.EQ.niter) THEN
249!       IF(t_2m(i).LT.t1(i).AND.t_2m(i).LT.ts1(i)) THEN
250!         PRINT*,' BAD t2m LT ',i,nsrf,t_2m(i),t1(i),ts1(i)
251!       ELSEIF(t_2m(i).GT.t1(i).AND.t_2m(i).GT.ts1(i)) THEN
252!         PRINT*,' BAD t2m GT ',i,nsrf,t_2m(i),t1(i),ts1(i)
253!       ENDIF
254!       ENDIF
255!IM ---
256      ENDDO
257!
258!
259!----------First aproximation of variables at zref --------------------------
260!
261      zref = 10.0
262      CALL screenp(klon, knon, nsrf, speed, tpot, q1, &
263 &                 ts1, qsurf, z0m, lmon, &
264 &                 ustar, testar, qstar, zref, &
265 &                 delu, delte, delq)
266!
267      DO i = 1, knon
268        u_zref(i) = delu(i)
269        q_zref(i) = max(qsurf(i),0.0) + delq(i)
270        te_zref(i) = ts1(i) + delte(i)
271        temp(i) = te_zref(i) * (psol(i)/pat1(i))**(-RKAPPA)
272!       temp(i) = te_zref(i) - (zref* RG)/RCPD/ &
273!                 (1 + RVTMP2 * max(q_zref(i),0.0))
274        u_zref_p(i) = u_zref(i)
275      ENDDO
276!
277! Iteration of the variables at the reference level zref : corrector ; see Hess & McAvaney, 1995
278!
279      DO n = 1, niter
280!
281        okri=.TRUE.
282        CALL screenc(klon, knon, nsrf, zxli, &
283 &                   u_zref, temp, q_zref, zref, &
284 &                   ts1, qsurf, z0m, z0h, psol, &
285 &                   ustar, testar, qstar, okri, ri1, &
286 &                   pref, delu, delte, delq, s_pblh ,prain, tsol, pat1)
287!
288        DO i = 1, knon
289          u_zref(i) = delu(i)
290          q_zref(i) = delq(i) + max(qsurf(i),0.0)
291          te_zref(i) = delte(i) + ts1(i)
292          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
293!         temp(i) = te_zref(i) - (zref* RG)/RCPD/ &
294!                   (1 + RVTMP2 * max(q_zref(i),0.0))
295        ENDDO
296!
297      ENDDO
298!
299      DO i = 1, knon
300        u_zref_c(i) = u_zref(i)
301!
302        u_10m(i) = u_zref_p(i) * ok_pred(i) + u_zref_c(i) * ok_corr(i)
303!
304!AM
305        q_zref_c(i) = q_zref(i)
306        temp_c(i) = temp(i)
307        t_10m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
308        q_10m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
309!MA
310      ENDDO
311!
312      RETURN
313      END subroutine stdlevvar
314!
315      SUBROUTINE stdlevvarn(klon, knon, nsrf, zxli, &
316                           u1, v1, t1, q1, z1, &
317                           ts1, qsurf, z0m, z0h, psol, pat1, &
318                           t_2m, q_2m, t_10m, q_10m, u_10m, ustar, &
319                           n2mout)
320!
321      USE ioipsl_getin_p_mod, ONLY : getin_p
322      IMPLICIT NONE
323!-------------------------------------------------------------------------
324!
325! Objet : calcul de la temperature et l'humidite relative a 2m et du
326!         module du vent a 10m a partir des relations de Dyer-Businger et
327!         des equations de Louis.
328!
329! Reference : Hess, Colman et McAvaney (1995)       
330!
331! I. Musat, 01.07.2002
332!
333!AM On rajoute en sortie t et q a 10m pr le calcule d'hbtm2 dans clmain
334!
335!-------------------------------------------------------------------------
336!
337! klon----input-I- dimension de la grille physique (= nb_pts_latitude X nb_pts_longitude)
338! knon----input-I- nombre de points pour un type de surface
339! nsrf----input-I- indice pour le type de surface; voir indice_sol_mod.F90
340! zxli----input-L- TRUE si calcul des cdrags selon Laurent Li
341! u1------input-R- vent zonal au 1er niveau du modele
342! v1------input-R- vent meridien au 1er niveau du modele
343! t1------input-R- temperature de l'air au 1er niveau du modele
344! q1------input-R- humidite relative au 1er niveau du modele
345! z1------input-R- geopotentiel au 1er niveau du modele
346! ts1-----input-R- temperature de l'air a la surface
347! qsurf---input-R- humidite relative a la surface
348! z0m, z0h---input-R- rugosite
349! psol----input-R- pression au sol
350! pat1----input-R- pression au 1er niveau du modele
351!
352! t_2m---output-R- temperature de l'air a 2m
353! q_2m---output-R- humidite relative a 2m
354! u_2m--output-R- vitesse du vent a 2m
355! u_10m--output-R- vitesse du vent a 10m
356! ustar--output-R- u*
357!AM
358! t_10m--output-R- temperature de l'air a 10m
359! q_10m--output-R- humidite specifique a 10m
360!
361      INTEGER, intent(in) :: klon, knon, nsrf
362      LOGICAL, intent(in) :: zxli
363      REAL, dimension(klon), intent(in) :: u1, v1, t1, q1, ts1, z1
364      REAL, dimension(klon), intent(inout) :: z0m, z0h
365      REAL, dimension(klon), intent(in) :: qsurf
366      REAL, dimension(klon), intent(in) :: psol, pat1
367!
368      REAL, dimension(klon), intent(out) :: t_2m, q_2m, ustar
369      REAL, dimension(klon), intent(out) :: u_10m, t_10m, q_10m
370      INTEGER, dimension(klon, 6), intent(out) :: n2mout
371!
372      REAL, dimension(klon) :: u_2m
373      REAL, dimension(klon) :: cdrm2m, cdrh2m, ri2m
374      REAL, dimension(klon) :: cdram, cdrah, zri1
375      REAL, dimension(klon) :: cdmn1, cdhn1, fm1, fh1
376      REAL, dimension(klon) :: cdmn2m, cdhn2m, fm2m, fh2m
377      REAL, dimension(klon) :: ri2m_new
378      REAL, DIMENSION(klon) :: s_pblh
379      REAL, DIMENSION(klon) :: prain
380      REAL, DIMENSION(klon) :: tsol
381!-------------------------------------------------------------------------
382      include "flux_arp.h"
383      include "YOMCST.h"
384!IM PLUS
385      include "YOETHF.h"
386!
387! Quelques constantes et options:
388!
389! RKAR : constante de von Karman
390      REAL, PARAMETER :: RKAR=0.40
391! niter : nombre iterations calcul "corrector"
392!     INTEGER, parameter :: niter=6, ncon=niter-1
393!IM 071020     INTEGER, parameter :: niter=2, ncon=niter-1
394      INTEGER, parameter :: niter=2, ncon=niter
395!     INTEGER, parameter :: niter=6, ncon=niter
396!
397! Variables locales
398      INTEGER :: i, n
399      REAL :: zref
400      REAL, dimension(klon) :: speed
401! tpot : temperature potentielle
402      REAL, dimension(klon) :: tpot
403      REAL, dimension(klon) :: cdran
404! ri1 : nb. de Richardson entre la surface --> la 1ere couche
405      REAL, dimension(klon) :: ri1
406      DOUBLE PRECISION, parameter :: eps=1.0D-20
407      REAL, dimension(klon) :: delu, delte, delq
408      REAL, dimension(klon) :: delh, delm
409      REAL, dimension(klon) :: delh_new, delm_new
410      REAL, dimension(klon) :: u_zref, te_zref, q_zref 
411      REAL, dimension(klon) :: u_zref_pnew, te_zref_pnew, q_zref_pnew
412      REAL, dimension(klon) :: temp, pref
413      REAL, dimension(klon) :: temp_new, pref_new
414      LOGICAL :: okri
415      REAL, dimension(klon) :: u_zref_p, te_zref_p, temp_p, q_zref_p
416      REAL, dimension(klon) :: u_zref_p_new, te_zref_p_new, temp_p_new, q_zref_p_new
417!convergence
418      REAL, dimension(klon) :: te_zref_con, q_zref_con
419      REAL, dimension(klon) :: u_zref_c, te_zref_c, temp_c, q_zref_c
420      REAL, dimension(klon) :: ok_pred, ok_corr
421!
422      REAL, dimension(klon) :: cdrm10m, cdrh10m, ri10m
423      REAL, dimension(klon) :: cdmn10m, cdhn10m, fm10m, fh10m
424      REAL, dimension(klon) :: cdn2m, cdn1, zri_zero
425      REAL :: CEPDUE,zdu2
426      INTEGER :: nzref, nz1
427      LOGICAL, dimension(klon) :: ok_t2m_toosmall, ok_t2m_toobig
428      LOGICAL, dimension(klon) :: ok_q2m_toosmall, ok_q2m_toobig
429      LOGICAL, dimension(klon) :: ok_u2m_toobig
430      LOGICAL, dimension(klon) :: ok_t10m_toosmall, ok_t10m_toobig
431      LOGICAL, dimension(klon) :: ok_q10m_toosmall, ok_q10m_toobig
432      LOGICAL, dimension(klon) :: ok_u10m_toobig
433      INTEGER, dimension(klon, 6) :: n10mout
434
435!-------------------------------------------------------------------------
436      CEPDUE=0.1
437!
438! n2mout : compteur des pas de temps ou t2m,q2m ou u2m sont en dehors des intervalles
439!          [tsurf, temp], [qsurf, q1] ou [0, speed]
440! n10mout : compteur des pas de temps ou t10m,q10m ou u10m sont en dehors des intervalles
441!          [tsurf, temp], [qsurf, q1] ou [0, speed]
442!
443      n2mout(:,:)=0
444      n10mout(:,:)=0
445     
446      DO i=1, knon
447       speed(i)=MAX(SQRT(u1(i)**2+v1(i)**2),CEPDUE)
448       ri1(i) = 0.0
449      ENDDO
450!
451      okri=.FALSE.
452      CALL cdrag(knon, nsrf, &
453 &                   speed, t1, q1, z1, &
454 &                   psol, s_pblh, ts1, qsurf, z0m, z0h, &
455 &                   zri_zero, 0, &
456 &                   cdram, cdrah, zri1, pref, prain, tsol, pat1)
457
458!
459      DO i = 1, knon
460        ri1(i) = zri1(i)
461        tpot(i) = t1(i)* (psol(i)/pat1(i))**RKAPPA
462        zdu2 = MAX(CEPDUE*CEPDUE, speed(i)**2)
463        ustar(i) = sqrt(cdram(i) * zdu2)
464!
465      ENDDO
466!
467!----------First aproximation of variables at zref --------------------------
468      zref = 2.0
469!
470! calcul first-guess en utilisant dans les calculs à 2m
471! le Richardson de la premiere couche atmospherique
472!
473       CALL screencn(klon, knon, nsrf, zxli, &
474 &                   speed, tpot, q1, zref, &
475 &                   ts1, qsurf, z0m, z0h, psol, &           
476 &                   cdram, cdrah,  okri, &
477 &                   ri1, 1, &
478 &                   pref_new, delm_new, delh_new, ri2m, &
479 &                   s_pblh, prain, tsol, pat1      )
480!
481       DO i = 1, knon
482         u_zref(i) = delm_new(i)*speed(i)
483         u_zref_p(i) = u_zref(i)
484         q_zref(i) = delh_new(i)*max(q1(i),0.0) + &
485         &           max(qsurf(i),0.0)*(1-delh_new(i))
486         q_zref_p(i) = q_zref(i)
487         te_zref(i) = delh_new(i)*tpot(i) + ts1(i)*(1-delh_new(i))
488         te_zref_p(i) = te_zref(i)
489!
490! return to normal temperature
491         temp(i) = te_zref(i) * (psol(i)/pref_new(i))**(-RKAPPA)
492         temp_p(i) = temp(i)
493!
494! compteurs ici
495!
496         ok_t2m_toosmall(i)=te_zref(i).LT.tpot(i).AND. &
497         & te_zref(i).LT.ts1(i)
498         ok_t2m_toobig(i)=te_zref(i).GT.tpot(i).AND. &
499         & te_zref(i).GT.ts1(i)
500         ok_q2m_toosmall(i)=q_zref(i).LT.q1(i).AND. &
501         & q_zref(i).LT.qsurf(i)
502         ok_q2m_toobig(i)=q_zref(i).GT.q1(i).AND. &
503         & q_zref(i).GT.qsurf(i)
504         ok_u2m_toobig(i)=u_zref(i).GT.speed(i)
505!
506         IF(ok_t2m_toosmall(i).OR.ok_t2m_toobig(i)) THEN
507             n2mout(i,1)=n2mout(i,1)+1
508         ENDIF
509         IF(ok_q2m_toosmall(i).OR.ok_q2m_toobig(i)) THEN
510             n2mout(i,3)=n2mout(i,3)+1
511         ENDIF
512         IF(ok_u2m_toobig(i)) THEN
513             n2mout(i,5)=n2mout(i,5)+1
514         ENDIF
515!
516         IF(ok_t2m_toosmall(i).OR.ok_t2m_toobig(i).OR. &
517          & ok_q2m_toosmall(i).OR.ok_q2m_toobig(i).OR. &
518          & ok_u2m_toobig(i)) THEN
519             delm_new(i)=min(max(delm_new(i),0.),1.)
520             delh_new(i)=min(max(delh_new(i),0.),1.)
521             u_zref(i) = delm_new(i)*speed(i)
522             u_zref_p(i) = u_zref(i)
523             q_zref(i) = delh_new(i)*max(q1(i),0.0) + &
524         &               max(qsurf(i),0.0)*(1-delh_new(i))
525             q_zref_p(i) = q_zref(i)
526             te_zref(i) = delh_new(i)*tpot(i) + ts1(i)*(1-delh_new(i))
527             te_zref_p(i) = te_zref(i)
528!
529! return to normal temperature
530             temp(i) = te_zref(i) * (psol(i)/pref_new(i))**(-RKAPPA)
531             temp_p(i) = temp(i)
532         ENDIF
533!
534       ENDDO
535!
536! Iteration of the variables at the reference level zref : corrector calculation ; see Hess & McAvaney, 1995
537!
538      DO n = 1, niter
539!
540        okri=.TRUE.
541        CALL screencn(klon, knon, nsrf, zxli, &
542 &                   u_zref, temp, q_zref, zref, &
543 &                   ts1, qsurf, z0m, z0h, psol, &
544 &                   cdram, cdrah,  okri, &
545 &                   ri1, 0, &
546 &                   pref, delm, delh, ri2m, &
547 &                   s_pblh, prain, tsol, pat1      )
548!
549        DO i = 1, knon
550          u_zref(i) = delm(i)*speed(i)
551          q_zref(i) = delh(i)*max(q1(i),0.0) + &
552          &           max(qsurf(i),0.0)*(1-delh(i))
553          te_zref(i) = delh(i)*tpot(i) + ts1(i)*(1-delh(i))
554!
555! return to normal temperature
556          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
557!
558! compteurs ici
559!
560          ok_t2m_toosmall(i)=te_zref(i).LT.tpot(i).AND. &
561          & te_zref(i).LT.ts1(i)
562          ok_t2m_toobig(i)=te_zref(i).GT.tpot(i).AND. &
563          & te_zref(i).GT.ts1(i)
564          ok_q2m_toosmall(i)=q_zref(i).LT.q1(i).AND. &
565          & q_zref(i).LT.qsurf(i)
566          ok_q2m_toobig(i)=q_zref(i).GT.q1(i).AND. &
567          & q_zref(i).GT.qsurf(i)
568          ok_u2m_toobig(i)=u_zref(i).GT.speed(i)
569!
570          IF(ok_t2m_toosmall(i).OR.ok_t2m_toobig(i)) THEN
571              n2mout(i,2)=n2mout(i,2)+1
572          ENDIF
573          IF(ok_q2m_toosmall(i).OR.ok_q2m_toobig(i)) THEN
574              n2mout(i,4)=n2mout(i,4)+1
575          ENDIF
576          IF(ok_u2m_toobig(i)) THEN
577              n2mout(i,6)=n2mout(i,6)+1
578          ENDIF
579!
580          IF(ok_t2m_toosmall(i).OR.ok_t2m_toobig(i).OR. &
581           & ok_q2m_toosmall(i).OR.ok_q2m_toobig(i).OR. &
582           & ok_u2m_toobig(i)) THEN
583              delm(i)=min(max(delm(i),0.),1.)
584              delh(i)=min(max(delh(i),0.),1.)
585              u_zref(i) = delm(i)*speed(i)
586              q_zref(i) = delh(i)*max(q1(i),0.0) + &
587          &           max(qsurf(i),0.0)*(1-delh(i))
588              te_zref(i) = delh(i)*tpot(i) + ts1(i)*(1-delh(i))
589              temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
590          ENDIF
591!
592!
593          IF(n.EQ.ncon) THEN
594            te_zref_con(i) = te_zref(i)
595            q_zref_con(i) = q_zref(i)
596          ENDIF
597!
598        ENDDO
599!
600      ENDDO
601!
602      DO i = 1, knon
603        q_zref_c(i) = q_zref(i)
604        temp_c(i) = temp(i)
605!
606        ok_pred(i)=0.
607        ok_corr(i)=1.
608!
609        t_2m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
610        q_2m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
611!
612        u_zref_c(i) = u_zref(i)
613        u_2m(i) = u_zref_p(i) * ok_pred(i) + u_zref_c(i) * ok_corr(i)
614      ENDDO
615!
616!
617!----------First aproximation of variables at zref --------------------------
618!
619      zref = 10.0
620!
621       CALL screencn(klon, knon, nsrf, zxli, &
622 &                   speed, tpot, q1, zref, &
623 &                   ts1, qsurf, z0m, z0h, psol, &           
624 &                   cdram, cdrah,  okri, &
625 &                   ri1, 1, &
626 &                   pref_new, delm_new, delh_new, ri10m, &
627 &                   s_pblh, prain, tsol, pat1      )
628!
629       DO i = 1, knon
630         u_zref(i) = delm_new(i)*speed(i)
631         q_zref(i) = delh_new(i)*max(q1(i),0.0) + &
632         &           max(qsurf(i),0.0)*(1-delh_new(i))
633         te_zref(i) = delh_new(i)*tpot(i) + ts1(i)*(1-delh_new(i))
634         temp(i) = te_zref(i) * (psol(i)/pref_new(i))**(-RKAPPA)
635         u_zref_p(i) = u_zref(i)
636!
637! compteurs ici
638!
639         ok_t10m_toosmall(i)=te_zref(i).LT.tpot(i).AND. &
640         & te_zref(i).LT.ts1(i)
641         ok_t10m_toobig(i)=te_zref(i).GT.tpot(i).AND. &
642         & te_zref(i).GT.ts1(i)
643         ok_q10m_toosmall(i)=q_zref(i).LT.q1(i).AND. &
644         & q_zref(i).LT.qsurf(i)
645         ok_q10m_toobig(i)=q_zref(i).GT.q1(i).AND. &
646         & q_zref(i).GT.qsurf(i)
647         ok_u10m_toobig(i)=u_zref(i).GT.speed(i)
648!
649         IF(ok_t10m_toosmall(i).OR.ok_t10m_toobig(i)) THEN
650             n10mout(i,1)=n10mout(i,1)+1
651         ENDIF
652         IF(ok_q10m_toosmall(i).OR.ok_q10m_toobig(i)) THEN
653             n10mout(i,3)=n10mout(i,3)+1
654         ENDIF
655         IF(ok_u10m_toobig(i)) THEN
656             n10mout(i,5)=n10mout(i,5)+1
657         ENDIF
658!
659         IF(ok_t10m_toosmall(i).OR.ok_t10m_toobig(i).OR. &
660          & ok_q10m_toosmall(i).OR.ok_q10m_toobig(i).OR. &
661          & ok_u10m_toobig(i)) THEN
662             delm_new(i)=min(max(delm_new(i),0.),1.)
663             delh_new(i)=min(max(delh_new(i),0.),1.)
664             u_zref(i) = delm_new(i)*speed(i)
665             u_zref_p(i) = u_zref(i)
666             q_zref(i) = delh_new(i)*max(q1(i),0.0) + &
667         &               max(qsurf(i),0.0)*(1-delh_new(i))
668             te_zref(i) = delh_new(i)*tpot(i) + ts1(i)*(1-delh_new(i))
669             temp(i) = te_zref(i) * (psol(i)/pref_new(i))**(-RKAPPA)
670         ENDIF
671!
672       ENDDO
673!
674! Iteration of the variables at the reference level zref : corrector calculation ; see Hess & McAvaney, 1995
675!
676      DO n = 1, niter
677!
678        okri=.TRUE.
679        CALL screencn(klon, knon, nsrf, zxli, &
680 &                   u_zref, temp, q_zref, zref, &
681 &                   ts1, qsurf, z0m, z0h, psol, &
682 &                   cdram, cdrah,  okri, &
683 &                   ri1, 0, &
684 &                   pref, delm, delh, ri10m, &
685 &                   s_pblh, prain, tsol, pat1      )
686!
687        DO i = 1, knon
688          u_zref(i) = delm(i)*speed(i)
689          q_zref(i) = delh(i)*max(q1(i),0.0) + &
690          &           max(qsurf(i),0.0)*(1-delh(i))
691          te_zref(i) = delh(i)*tpot(i) + ts1(i)*(1-delh(i))
692!
693! return to normal temperature
694          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
695!
696! compteurs ici
697!
698          ok_t10m_toosmall(i)=te_zref(i).LT.tpot(i).AND. &
699          & te_zref(i).LT.ts1(i)
700          ok_t10m_toobig(i)=te_zref(i).GT.tpot(i).AND. &
701          & te_zref(i).GT.ts1(i)
702          ok_q10m_toosmall(i)=q_zref(i).LT.q1(i).AND. &
703          & q_zref(i).LT.qsurf(i)
704          ok_q10m_toobig(i)=q_zref(i).GT.q1(i).AND. &
705          & q_zref(i).GT.qsurf(i)
706          ok_u10m_toobig(i)=u_zref(i).GT.speed(i)
707!
708          IF(ok_t10m_toosmall(i).OR.ok_t10m_toobig(i)) THEN
709              n10mout(i,2)=n10mout(i,2)+1
710          ENDIF
711          IF(ok_q10m_toosmall(i).OR.ok_q10m_toobig(i)) THEN
712              n10mout(i,4)=n10mout(i,4)+1
713          ENDIF
714          IF(ok_u10m_toobig(i)) THEN
715              n10mout(i,6)=n10mout(i,6)+1
716          ENDIF
717!
718          IF(ok_t10m_toosmall(i).OR.ok_t10m_toobig(i).OR. &
719           & ok_q10m_toosmall(i).OR.ok_q10m_toobig(i).OR. &
720           & ok_u10m_toobig(i)) THEN
721              delm(i)=min(max(delm(i),0.),1.)
722              delh(i)=min(max(delh(i),0.),1.)
723              u_zref(i) = delm(i)*speed(i)
724              q_zref(i) = delh(i)*max(q1(i),0.0) + &
725          &           max(qsurf(i),0.0)*(1-delh(i))
726              te_zref(i) = delh(i)*tpot(i) + ts1(i)*(1-delh(i))
727              temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
728          ENDIF
729!
730!
731          IF(n.EQ.ncon) THEN
732            te_zref_con(i) = te_zref(i)
733            q_zref_con(i) = q_zref(i)
734          ENDIF
735!
736        ENDDO
737!
738      ENDDO
739!
740      DO i = 1, knon
741        q_zref_c(i) = q_zref(i)
742        temp_c(i) = temp(i)
743!
744        ok_pred(i)=0.
745        ok_corr(i)=1.
746!
747        t_10m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
748        q_10m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
749!
750        u_zref_c(i) = u_zref(i)
751        u_10m(i) = u_zref_p(i) * ok_pred(i) + u_zref_c(i) * ok_corr(i)
752      ENDDO
753!
754      RETURN
755      END subroutine stdlevvarn
756
757END MODULE stdlevvar_mod
Note: See TracBrowser for help on using the repository browser.