source: LMDZ6/branches/PBLSURF_GPUPORT/libf/phylmd/stdlevvar_mod.f90 @ 5914

Last change on this file since 5914 was 5882, checked in by yann meurdesoif, 5 weeks ago

GPU port of stdlevvar & stdlevvarn
YM

  • Property svn:executable set to *
File size: 27.7 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 yomcst_mod_h
8      USE cdrag_mod
9  USE screenp_mod
10  USE screenc_mod
11  IMPLICIT NONE
12
13CONTAINS
14!
15!****************************************************************************************
16!
17!r original routine svn3623
18!
19      SUBROUTINE stdlevvar(klon, knon, nsrf, zxli, &
20                           u1, v1, t1, q1, z1, &
21                           ts1, qsurf, z0m, z0h, psol, pat1, &
22                           t_2m, q_2m, t_10m, q_10m, u_10m, ustar, s_pblh, prain, tsol)
23!$gpum horizontal knon klon
24        USE yoethf_mod_h
25        USE flux_arp_mod_h
26      IMPLICIT NONE
27!-------------------------------------------------------------------------
28!
29! Objet : calcul de la temperature et l'humidite relative a 2m et du
30!         module du vent a 10m a partir des relations de Dyer-Businger et
31!         des equations de Louis.
32!
33! Reference : Hess, Colman et McAvaney (1995)
34!
35! I. Musat, 01.07.2002
36!
37!AM On rajoute en sortie t et q a 10m pr le calcule d'hbtm2 dans clmain
38!
39!-------------------------------------------------------------------------
40!
41! klon----input-I- dimension de la grille physique (= nb_pts_latitude X nb_pts_longitude)
42! knon----input-I- nombre de points pour un type de surface
43! nsrf----input-I- indice pour le type de surface; voir indice_sol_mod.F90
44! zxli----input-L- TRUE si calcul des cdrags selon Laurent Li
45! u1------input-R- vent zonal au 1er niveau du modele
46! v1------input-R- vent meridien au 1er niveau du modele
47! t1------input-R- temperature de l'air au 1er niveau du modele
48! q1------input-R- humidite relative au 1er niveau du modele
49! z1------input-R- geopotentiel au 1er niveau du modele
50! ts1-----input-R- temperature de l'air a la surface
51! qsurf---input-R- humidite relative a la surface
52! z0m, z0h---input-R- rugosite
53! psol----input-R- pression au sol
54! pat1----input-R- pression au 1er niveau du modele
55!
56! t_2m---output-R- temperature de l'air a 2m
57! q_2m---output-R- humidite relative a 2m
58! u_10m--output-R- vitesse du vent a 10m
59!AM
60! t_10m--output-R- temperature de l'air a 10m
61! q_10m--output-R- humidite specifique a 10m
62! ustar--output-R- u*
63!
64      INTEGER, intent(in) :: klon, knon, nsrf
65      LOGICAL, intent(in) :: zxli
66      REAL, dimension(klon), intent(in) :: u1, v1, t1, q1, z1, ts1
67      REAL, dimension(klon), intent(in) :: qsurf
68      REAL, dimension(klon), intent(inout) :: z0m, z0h
69      REAL, dimension(klon), intent(in) :: psol, pat1
70!
71      REAL, dimension(klon), intent(out) :: t_2m, q_2m, ustar
72      REAL, dimension(klon), intent(out) :: u_10m, t_10m, q_10m
73      REAL, DIMENSION(klon), INTENT(INOUT) :: s_pblh
74      REAL, DIMENSION(klon), INTENT(IN) :: prain
75      REAL, DIMENSION(klon), INTENT(IN) :: tsol
76!-------------------------------------------------------------------------
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!$gpum horizontal knon klon
321
322      USE indice_sol_mod, ONLY : nbsrf
323      USE yomcst_mod_h
324      USE ioipsl_getin_p_mod, ONLY : getin_p
325      USE yoethf_mod_h
326      USE flux_arp_mod_h
327      USE lmdz_checksum
328      IMPLICIT NONE
329!-------------------------------------------------------------------------
330!
331! Objet : calcul de la temperature et l'humidite relative a 2m et du
332!         module du vent a 10m a partir des relations de Dyer-Businger et
333!         des equations de Louis.
334!
335! Reference : Hess, Colman et McAvaney (1995)
336!
337! I. Musat, 01.07.2002
338!
339!AM On rajoute en sortie t et q a 10m pr le calcule d'hbtm2 dans clmain
340!
341!-------------------------------------------------------------------------
342!
343! klon----input-I- dimension de la grille physique (= nb_pts_latitude X nb_pts_longitude)
344! knon----input-I- nombre de points pour un type de surface
345! nsrf----input-I- indice pour le type de surface; voir indice_sol_mod.F90
346! zxli----input-L- TRUE si calcul des cdrags selon Laurent Li
347! u1------input-R- vent zonal au 1er niveau du modele
348! v1------input-R- vent meridien au 1er niveau du modele
349! t1------input-R- temperature de l'air au 1er niveau du modele
350! q1------input-R- humidite relative au 1er niveau du modele
351! z1------input-R- geopotentiel au 1er niveau du modele
352! ts1-----input-R- temperature de l'air a la surface
353! qsurf---input-R- humidite relative a la surface
354! z0m, z0h---input-R- rugosite
355! psol----input-R- pression au sol
356! pat1----input-R- pression au 1er niveau du modele
357!
358! t_2m---output-R- temperature de l'air a 2m
359! q_2m---output-R- humidite relative a 2m
360! u_2m--output-R- vitesse du vent a 2m
361! u_10m--output-R- vitesse du vent a 10m
362! ustar--output-R- u*
363!AM
364! t_10m--output-R- temperature de l'air a 10m
365! q_10m--output-R- humidite specifique a 10m
366!
367      INTEGER, intent(in) :: klon, knon, nsrf
368      LOGICAL, intent(in) :: zxli
369      REAL, dimension(klon), intent(in) :: u1, v1, t1, q1, ts1, z1
370      REAL, dimension(klon), intent(inout) :: z0m, z0h
371      REAL, dimension(klon), intent(in) :: qsurf
372      REAL, dimension(klon), intent(in) :: psol, pat1
373!
374      REAL, dimension(klon), intent(out) :: t_2m, q_2m, ustar
375      REAL, dimension(klon), intent(out) :: u_10m, t_10m, q_10m
376      INTEGER, dimension(klon, nbsrf, 6), intent(out) :: n2mout
377!
378      REAL, dimension(klon) :: u_2m
379      REAL, dimension(klon) :: cdrm2m, cdrh2m, ri2m
380      REAL, dimension(klon) :: cdram, cdrah, zri1
381      REAL, dimension(klon) :: cdmn1, cdhn1, fm1, fh1
382      REAL, dimension(klon) :: cdmn2m, cdhn2m, fm2m, fh2m
383      REAL, dimension(klon) :: ri2m_new
384      REAL, DIMENSION(klon) :: s_pblh
385      REAL, DIMENSION(klon) :: prain
386      REAL, DIMENSION(klon) :: tsol
387!-------------------------------------------------------------------------
388!
389! Quelques constantes et options:
390!
391! RKAR : constante de von Karman
392      REAL, PARAMETER :: RKAR=0.40
393! niter : nombre iterations calcul "corrector"
394!     INTEGER, parameter :: niter=6, ncon=niter-1
395!IM 071020     INTEGER, parameter :: niter=2, ncon=niter-1
396      INTEGER, parameter :: niter=2, ncon=niter
397!     INTEGER, parameter :: niter=6, ncon=niter
398!
399! Variables locales
400      INTEGER :: i, n
401      REAL :: zref
402      REAL, dimension(klon) :: speed
403! tpot : temperature potentielle
404      REAL, dimension(klon) :: tpot
405      REAL, dimension(klon) :: cdran
406! ri1 : nb. de Richardson entre la surface --> la 1ere couche
407      REAL, dimension(klon) :: ri1
408      DOUBLE PRECISION, parameter :: eps=1.0D-20
409      REAL, dimension(klon) :: delu, delte, delq
410      REAL, dimension(klon) :: delh, delm
411      REAL, dimension(klon) :: delh_new, delm_new
412      REAL, dimension(klon) :: u_zref, te_zref, q_zref 
413      REAL, dimension(klon) :: u_zref_pnew, te_zref_pnew, q_zref_pnew
414      REAL, dimension(klon) :: temp, pref
415      REAL, dimension(klon) :: temp_new, pref_new
416      LOGICAL :: okri
417      REAL, dimension(klon) :: u_zref_p, te_zref_p, temp_p, q_zref_p
418      REAL, dimension(klon) :: u_zref_p_new, te_zref_p_new, temp_p_new, q_zref_p_new
419!convergence
420      REAL, dimension(klon) :: te_zref_con, q_zref_con
421      REAL, dimension(klon) :: u_zref_c, te_zref_c, temp_c, q_zref_c
422      REAL, dimension(klon) :: ok_pred, ok_corr
423!
424      REAL, dimension(klon) :: cdrm10m, cdrh10m, ri10m
425      REAL, dimension(klon) :: cdmn10m, cdhn10m, fm10m, fh10m
426      REAL, dimension(klon) :: cdn2m, cdn1, zri_zero
427      REAL :: CEPDUE,zdu2
428      INTEGER :: nzref, nz1
429      LOGICAL, dimension(klon) :: ok_t2m_toosmall, ok_t2m_toobig
430      LOGICAL, dimension(klon) :: ok_q2m_toosmall, ok_q2m_toobig
431      LOGICAL, dimension(klon) :: ok_u2m_toobig
432      LOGICAL, dimension(klon) :: ok_t10m_toosmall, ok_t10m_toobig
433      LOGICAL, dimension(klon) :: ok_q10m_toosmall, ok_q10m_toobig
434      LOGICAL, dimension(klon) :: ok_u10m_toobig
435      INTEGER, dimension(klon, 6) :: n10mout
436
437!$gpum nocall checksum
438      CALL checksum("u1", u1)
439      CALL checksum("v1", v1)
440      CALL checksum("t1", t1)
441      CALL checksum("q1", q1)
442      CALL checksum("z1", z1)
443      CALL checksum("ts1", ts1)
444      CALL checksum("qsurf", qsurf)
445      CALL checksum("z0m", z0m)
446      CALL checksum("z0h", z0h)
447      CALL checksum("psol", psol)
448      CALL checksum("pat1", pat1)
449      CALL checksum("t_2m", t_2m)
450      CALL checksum("q_2m", q_2m)
451      CALL checksum("t_10m", t_10m)
452      CALL checksum("q_10m", q_10m)
453      CALL checksum("u_10m", u_10m)
454      CALL checksum("ustar", ustar)
455      CALL checksum("s_pblh", s_pblh)
456      CALL checksum("prain", prain)
457      CALL checksum("tsol", tsol)
458
459!-------------------------------------------------------------------------
460      CEPDUE=0.1
461!
462! n2mout : compteur des pas de temps ou t2m,q2m ou u2m sont en dehors des intervalles
463!          [tsurf, temp], [qsurf, q1] ou [0, speed]
464! n10mout : compteur des pas de temps ou t10m,q10m ou u10m sont en dehors des intervalles
465!          [tsurf, temp], [qsurf, q1] ou [0, speed]
466!
467      n2mout(:, nsrf, :)=0
468      n10mout(:, :)=0
469     
470      DO i=1, knon
471       speed(i)=MAX(SQRT(u1(i)**2+v1(i)**2),CEPDUE)
472       ri1(i) = 0.0
473      ENDDO
474!
475      okri=.FALSE.
476      CALL cdrag(knon, nsrf, &
477 &                   speed, t1, q1, z1, &
478 &                   psol, s_pblh, ts1, qsurf, z0m, z0h, &
479 &                   zri_zero, 0, &
480 &                   cdram, cdrah, zri1, pref, prain, tsol, pat1)
481
482!
483      DO i = 1, knon
484        ri1(i) = zri1(i)
485        tpot(i) = t1(i)* (psol(i)/pat1(i))**RKAPPA
486        zdu2 = MAX(CEPDUE*CEPDUE, speed(i)**2)
487        ustar(i) = sqrt(cdram(i) * zdu2)
488!
489      ENDDO
490!
491!----------First aproximation of variables at zref --------------------------
492      zref = 2.0
493!
494! calcul first-guess en utilisant dans les calculs à 2m
495! le Richardson de la premiere couche atmospherique
496!
497       CALL screencn(klon, knon, nsrf, zxli, &
498 &                   speed, tpot, q1, zref, &
499 &                   ts1, qsurf, z0m, z0h, psol, &           
500 &                   cdram, cdrah,  okri, &
501 &                   ri1, 1, &
502 &                   pref_new, delm_new, delh_new, ri2m, &
503 &                   s_pblh, prain, tsol, pat1      )
504!
505       DO i = 1, knon
506         u_zref(i) = delm_new(i)*speed(i)
507         u_zref_p(i) = u_zref(i)
508         q_zref(i) = delh_new(i)*max(q1(i),0.0) + &
509         &           max(qsurf(i),0.0)*(1-delh_new(i))
510         q_zref_p(i) = q_zref(i)
511         te_zref(i) = delh_new(i)*tpot(i) + ts1(i)*(1-delh_new(i))
512         te_zref_p(i) = te_zref(i)
513!
514! return to normal temperature
515         temp(i) = te_zref(i) * (psol(i)/pref_new(i))**(-RKAPPA)
516         temp_p(i) = temp(i)
517!
518! compteurs ici
519!
520         ok_t2m_toosmall(i)=te_zref(i).LT.tpot(i).AND. &
521         & te_zref(i).LT.ts1(i)
522         ok_t2m_toobig(i)=te_zref(i).GT.tpot(i).AND. &
523         & te_zref(i).GT.ts1(i)
524         ok_q2m_toosmall(i)=q_zref(i).LT.q1(i).AND. &
525         & q_zref(i).LT.qsurf(i)
526         ok_q2m_toobig(i)=q_zref(i).GT.q1(i).AND. &
527         & q_zref(i).GT.qsurf(i)
528         ok_u2m_toobig(i)=u_zref(i).GT.speed(i)
529!
530         IF(ok_t2m_toosmall(i).OR.ok_t2m_toobig(i)) THEN
531             n2mout(i, nsrf, 1)=n2mout(i, nsrf, 1)+1
532         ENDIF
533         IF(ok_q2m_toosmall(i).OR.ok_q2m_toobig(i)) THEN
534             n2mout(i, nsrf, 3)=n2mout(i, nsrf, 3)+1
535         ENDIF
536         IF(ok_u2m_toobig(i)) THEN
537             n2mout(i, nsrf, 5)=n2mout(i, nsrf, 5)+1
538         ENDIF
539!
540         IF(ok_t2m_toosmall(i).OR.ok_t2m_toobig(i).OR. &
541          & ok_q2m_toosmall(i).OR.ok_q2m_toobig(i).OR. &
542          & ok_u2m_toobig(i)) THEN
543             delm_new(i)=min(max(delm_new(i),0.),1.)
544             delh_new(i)=min(max(delh_new(i),0.),1.)
545             u_zref(i) = delm_new(i)*speed(i)
546             u_zref_p(i) = u_zref(i)
547             q_zref(i) = delh_new(i)*max(q1(i),0.0) + &
548         &               max(qsurf(i),0.0)*(1-delh_new(i))
549             q_zref_p(i) = q_zref(i)
550             te_zref(i) = delh_new(i)*tpot(i) + ts1(i)*(1-delh_new(i))
551             te_zref_p(i) = te_zref(i)
552!
553! return to normal temperature
554             temp(i) = te_zref(i) * (psol(i)/pref_new(i))**(-RKAPPA)
555             temp_p(i) = temp(i)
556         ENDIF
557!
558       ENDDO
559!
560! Iteration of the variables at the reference level zref : corrector calculation ; see Hess & McAvaney, 1995
561!
562      DO n = 1, niter
563!
564        okri=.TRUE.
565        CALL screencn(klon, knon, nsrf, zxli, &
566 &                   u_zref, temp, q_zref, zref, &
567 &                   ts1, qsurf, z0m, z0h, psol, &
568 &                   cdram, cdrah,  okri, &
569 &                   ri1, 0, &
570 &                   pref, delm, delh, ri2m, &
571 &                   s_pblh, prain, tsol, pat1      )
572!
573        DO i = 1, knon
574          u_zref(i) = delm(i)*speed(i)
575          q_zref(i) = delh(i)*max(q1(i),0.0) + &
576          &           max(qsurf(i),0.0)*(1-delh(i))
577          te_zref(i) = delh(i)*tpot(i) + ts1(i)*(1-delh(i))
578!
579! return to normal temperature
580          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
581!
582! compteurs ici
583!
584          ok_t2m_toosmall(i)=te_zref(i).LT.tpot(i).AND. &
585          & te_zref(i).LT.ts1(i)
586          ok_t2m_toobig(i)=te_zref(i).GT.tpot(i).AND. &
587          & te_zref(i).GT.ts1(i)
588          ok_q2m_toosmall(i)=q_zref(i).LT.q1(i).AND. &
589          & q_zref(i).LT.qsurf(i)
590          ok_q2m_toobig(i)=q_zref(i).GT.q1(i).AND. &
591          & q_zref(i).GT.qsurf(i)
592          ok_u2m_toobig(i)=u_zref(i).GT.speed(i)
593!
594          IF(ok_t2m_toosmall(i).OR.ok_t2m_toobig(i)) THEN
595              n2mout(i, nsrf, 2)=n2mout(i, nsrf, 2)+1
596          ENDIF
597          IF(ok_q2m_toosmall(i).OR.ok_q2m_toobig(i)) THEN
598              n2mout(i, nsrf, 4)=n2mout(i, nsrf, 4)+1
599          ENDIF
600          IF(ok_u2m_toobig(i)) THEN
601              n2mout(i, nsrf, 6)=n2mout(i, nsrf, 6)+1
602          ENDIF
603!
604          IF(ok_t2m_toosmall(i).OR.ok_t2m_toobig(i).OR. &
605           & ok_q2m_toosmall(i).OR.ok_q2m_toobig(i).OR. &
606           & ok_u2m_toobig(i)) THEN
607              delm(i)=min(max(delm(i),0.),1.)
608              delh(i)=min(max(delh(i),0.),1.)
609              u_zref(i) = delm(i)*speed(i)
610              q_zref(i) = delh(i)*max(q1(i),0.0) + &
611          &           max(qsurf(i),0.0)*(1-delh(i))
612              te_zref(i) = delh(i)*tpot(i) + ts1(i)*(1-delh(i))
613              temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
614          ENDIF
615!
616!
617          IF(n.EQ.ncon) THEN
618            te_zref_con(i) = te_zref(i)
619            q_zref_con(i) = q_zref(i)
620          ENDIF
621!
622        ENDDO
623!
624      ENDDO
625!
626      DO i = 1, knon
627        q_zref_c(i) = q_zref(i)
628        temp_c(i) = temp(i)
629!
630        ok_pred(i)=0.
631        ok_corr(i)=1.
632!
633        t_2m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
634        q_2m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
635!
636        u_zref_c(i) = u_zref(i)
637        u_2m(i) = u_zref_p(i) * ok_pred(i) + u_zref_c(i) * ok_corr(i)
638      ENDDO
639!
640!
641!----------First aproximation of variables at zref --------------------------
642!
643      zref = 10.0
644!
645       CALL screencn(klon, knon, nsrf, zxli, &
646 &                   speed, tpot, q1, zref, &
647 &                   ts1, qsurf, z0m, z0h, psol, &           
648 &                   cdram, cdrah,  okri, &
649 &                   ri1, 1, &
650 &                   pref_new, delm_new, delh_new, ri10m, &
651 &                   s_pblh, prain, tsol, pat1      )
652!
653       DO i = 1, knon
654         u_zref(i) = delm_new(i)*speed(i)
655         q_zref(i) = delh_new(i)*max(q1(i),0.0) + &
656         &           max(qsurf(i),0.0)*(1-delh_new(i))
657         te_zref(i) = delh_new(i)*tpot(i) + ts1(i)*(1-delh_new(i))
658         temp(i) = te_zref(i) * (psol(i)/pref_new(i))**(-RKAPPA)
659         u_zref_p(i) = u_zref(i)
660!
661! compteurs ici
662!
663         ok_t10m_toosmall(i)=te_zref(i).LT.tpot(i).AND. &
664         & te_zref(i).LT.ts1(i)
665         ok_t10m_toobig(i)=te_zref(i).GT.tpot(i).AND. &
666         & te_zref(i).GT.ts1(i)
667         ok_q10m_toosmall(i)=q_zref(i).LT.q1(i).AND. &
668         & q_zref(i).LT.qsurf(i)
669         ok_q10m_toobig(i)=q_zref(i).GT.q1(i).AND. &
670         & q_zref(i).GT.qsurf(i)
671         ok_u10m_toobig(i)=u_zref(i).GT.speed(i)
672!
673         IF(ok_t10m_toosmall(i).OR.ok_t10m_toobig(i)) THEN
674             n10mout(i,1)=n10mout(i,1)+1
675         ENDIF
676         IF(ok_q10m_toosmall(i).OR.ok_q10m_toobig(i)) THEN
677             n10mout(i,3)=n10mout(i,3)+1
678         ENDIF
679         IF(ok_u10m_toobig(i)) THEN
680             n10mout(i,5)=n10mout(i,5)+1
681         ENDIF
682!
683         IF(ok_t10m_toosmall(i).OR.ok_t10m_toobig(i).OR. &
684          & ok_q10m_toosmall(i).OR.ok_q10m_toobig(i).OR. &
685          & ok_u10m_toobig(i)) THEN
686             delm_new(i)=min(max(delm_new(i),0.),1.)
687             delh_new(i)=min(max(delh_new(i),0.),1.)
688             u_zref(i) = delm_new(i)*speed(i)
689             u_zref_p(i) = u_zref(i)
690             q_zref(i) = delh_new(i)*max(q1(i),0.0) + &
691         &               max(qsurf(i),0.0)*(1-delh_new(i))
692             te_zref(i) = delh_new(i)*tpot(i) + ts1(i)*(1-delh_new(i))
693             temp(i) = te_zref(i) * (psol(i)/pref_new(i))**(-RKAPPA)
694         ENDIF
695!
696       ENDDO
697!
698! Iteration of the variables at the reference level zref : corrector calculation ; see Hess & McAvaney, 1995
699!
700      DO n = 1, niter
701!
702        okri=.TRUE.
703        CALL screencn(klon, knon, nsrf, zxli, &
704 &                   u_zref, temp, q_zref, zref, &
705 &                   ts1, qsurf, z0m, z0h, psol, &
706 &                   cdram, cdrah,  okri, &
707 &                   ri1, 0, &
708 &                   pref, delm, delh, ri10m, &
709 &                   s_pblh, prain, tsol, pat1      )
710!
711        DO i = 1, knon
712          u_zref(i) = delm(i)*speed(i)
713          q_zref(i) = delh(i)*max(q1(i),0.0) + &
714          &           max(qsurf(i),0.0)*(1-delh(i))
715          te_zref(i) = delh(i)*tpot(i) + ts1(i)*(1-delh(i))
716!
717! return to normal temperature
718          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
719!
720! compteurs ici
721!
722          ok_t10m_toosmall(i)=te_zref(i).LT.tpot(i).AND. &
723          & te_zref(i).LT.ts1(i)
724          ok_t10m_toobig(i)=te_zref(i).GT.tpot(i).AND. &
725          & te_zref(i).GT.ts1(i)
726          ok_q10m_toosmall(i)=q_zref(i).LT.q1(i).AND. &
727          & q_zref(i).LT.qsurf(i)
728          ok_q10m_toobig(i)=q_zref(i).GT.q1(i).AND. &
729          & q_zref(i).GT.qsurf(i)
730          ok_u10m_toobig(i)=u_zref(i).GT.speed(i)
731!
732          IF(ok_t10m_toosmall(i).OR.ok_t10m_toobig(i)) THEN
733              n10mout(i,2)=n10mout(i,2)+1
734          ENDIF
735          IF(ok_q10m_toosmall(i).OR.ok_q10m_toobig(i)) THEN
736              n10mout(i,4)=n10mout(i,4)+1
737          ENDIF
738          IF(ok_u10m_toobig(i)) THEN
739              n10mout(i,6)=n10mout(i,6)+1
740          ENDIF
741!
742          IF(ok_t10m_toosmall(i).OR.ok_t10m_toobig(i).OR. &
743           & ok_q10m_toosmall(i).OR.ok_q10m_toobig(i).OR. &
744           & ok_u10m_toobig(i)) THEN
745              delm(i)=min(max(delm(i),0.),1.)
746              delh(i)=min(max(delh(i),0.),1.)
747              u_zref(i) = delm(i)*speed(i)
748              q_zref(i) = delh(i)*max(q1(i),0.0) + &
749          &           max(qsurf(i),0.0)*(1-delh(i))
750              te_zref(i) = delh(i)*tpot(i) + ts1(i)*(1-delh(i))
751              temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
752          ENDIF
753!
754!
755          IF(n.EQ.ncon) THEN
756            te_zref_con(i) = te_zref(i)
757            q_zref_con(i) = q_zref(i)
758          ENDIF
759!
760        ENDDO
761!
762      ENDDO
763!
764      DO i = 1, knon
765        q_zref_c(i) = q_zref(i)
766        temp_c(i) = temp(i)
767!
768        ok_pred(i)=0.
769        ok_corr(i)=1.
770!
771        t_10m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
772        q_10m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
773!
774        u_zref_c(i) = u_zref(i)
775        u_10m(i) = u_zref_p(i) * ok_pred(i) + u_zref_c(i) * ok_corr(i)
776      ENDDO
777!
778      CALL checksum("u1", u1)
779      CALL checksum("v1", v1)
780      CALL checksum("t1", t1)
781      CALL checksum("q1", q1)
782      CALL checksum("z1", z1)
783      CALL checksum("ts1", ts1)
784      CALL checksum("qsurf", qsurf)
785      CALL checksum("z0m", z0m)
786      CALL checksum("z0h", z0h)
787      CALL checksum("psol", psol)
788      CALL checksum("pat1", pat1)
789      CALL checksum("t_2m", t_2m)
790      CALL checksum("q_2m", q_2m)
791      CALL checksum("t_10m", t_10m)
792      CALL checksum("q_10m", q_10m)
793      CALL checksum("u_10m", u_10m)
794      CALL checksum("ustar", ustar)
795      CALL checksum("s_pblh", s_pblh)
796      CALL checksum("prain", prain)
797      CALL checksum("tsol", tsol)
798      RETURN
799      END subroutine stdlevvarn
800
801END MODULE stdlevvar_mod
Note: See TracBrowser for help on using the repository browser.