source: LMDZ6/trunk/libf/phylmd/stdlevvar_mod.F90 @ 4842

Last change on this file since 4842 was 4722, checked in by Laurent Fairhead, 8 months ago

Modification by O. Torres to the cdrag routines to include different bulk formulae
to calculate cdrag coefficients over ocean as well as an iteration of that
calculation.
The iteration is controlled by flag ok_cdrag_iter which if set to FALSE by default
to converge with previous results.
The choice of bulk formulae is set with the choix_bulk parameter
The number of iterations to run is set with nit_bulk
OT, PB, CD, LF

  • 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.