GCC Code Coverage Report
Directory: ./ Exec Total Coverage
File: phylmd/stdlevvar_mod.F90 Lines: 175 267 65.5 %
Date: 2023-06-30 12:56:34 Branches: 159 202 78.7 %

Line Branch Exec Source
1
!
2
MODULE 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
12
CONTAINS
13
!
14
!****************************************************************************************
15
!
16
!r original routine svn3623
17
!
18
4349870
      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)
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, z0m, z0h
64
      REAL, dimension(klon), intent(in) :: psol, pat1
65
!
66
      REAL, dimension(klon), intent(out) :: t_2m, q_2m, ustar
67
      REAL, dimension(klon), intent(out) :: u_10m, t_10m, q_10m
68
!-------------------------------------------------------------------------
69
      include "flux_arp.h"
70
      include "YOMCST.h"
71
!IM PLUS
72
      include "YOETHF.h"
73
!
74
! Quelques constantes et options:
75
!
76
! RKAR : constante de von Karman
77
      REAL, PARAMETER :: RKAR=0.40
78
! niter : nombre iterations calcul "corrector"
79
!     INTEGER, parameter :: niter=6, ncon=niter-1
80
      INTEGER, parameter :: niter=2, ncon=niter-1
81
!
82
! Variables locales
83
      INTEGER :: i, n
84
      REAL :: zref
85
      REAL, dimension(klon) :: speed
86
! tpot : temperature potentielle
87
      REAL, dimension(klon) :: tpot
88
      REAL, dimension(klon) :: zri1, cdran
89
      REAL, dimension(klon) :: cdram, cdrah
90
! ri1 : nb. de Richardson entre la surface --> la 1ere couche
91
      REAL, dimension(klon) :: ri1
92
      REAL, dimension(klon) :: testar, qstar
93
      REAL, dimension(klon) :: zdte, zdq
94
! lmon : longueur de Monin-Obukhov selon Hess, Colman and McAvaney
95
      DOUBLE PRECISION, dimension(klon) :: lmon
96
      DOUBLE PRECISION, parameter :: eps=1.0D-20
97
      REAL, dimension(klon) :: delu, delte, delq
98
      REAL, dimension(klon) :: u_zref, te_zref, q_zref
99
      REAL, dimension(klon) :: temp, pref
100
      LOGICAL :: okri
101
      REAL, dimension(klon) :: u_zref_p, te_zref_p, temp_p, q_zref_p
102
!convertgence
103
      REAL, dimension(klon) :: te_zref_con, q_zref_con
104
      REAL, dimension(klon) :: u_zref_c, te_zref_c, temp_c, q_zref_c
105
      REAL, dimension(klon) :: ok_pred, ok_corr, zri_zero
106
!     REAL, dimension(klon) :: conv_te, conv_q
107
!-------------------------------------------------------------------------
108
      DO i=1, knon
109
       speed(i)=SQRT(u1(i)**2+v1(i)**2)
110
       ri1(i) = 0.0
111
      ENDDO
112
!
113
      okri=.FALSE.
114
!      CALL coefcdrag(klon, knon, nsrf, zxli, &
115
! &                   speed, t1, q1, z1, psol, &
116
! &                   ts1, qsurf, rugos, okri, ri1,  &
117
! &                   cdram, cdrah, cdran, zri1, pref)
118
! Fuxing WANG, 04/03/2015, replace the coefcdrag by the merged version: cdrag
119
120
      CALL cdrag(knon, nsrf, &
121
 &                   speed, t1, q1, z1, &
122
 &                   psol, ts1, qsurf, z0m, z0h, &
123
 &                   zri_zero, 0, &
124
 &                   cdram, cdrah, zri1, pref)
125
126
! --- special Dice: on force cdragm ( a defaut de forcer ustar) MPL 05082013
127
     IF (ok_prescr_ust) then
128
      DO i = 1, knon
129
       print *,'cdram avant=',cdram(i)
130
       cdram(i) = ust*ust/speed(i)/speed(i)
131
       print *,'cdram ust speed apres=',cdram(i),ust,speed
132
      ENDDO
133
     ENDIF
134
!
135
!---------Star variables----------------------------------------------------
136
!
137
      DO i = 1, knon
138
        ri1(i) = zri1(i)
139
        tpot(i) = t1(i)* (psol(i)/pat1(i))**RKAPPA
140
        ustar(i) = sqrt(cdram(i) * speed(i) * speed(i))
141
        zdte(i) = tpot(i) - ts1(i)
142
        zdq(i) = max(q1(i),0.0) - max(qsurf(i),0.0)
143
!
144
!
145
!IM BUG BUG BUG       zdte(i) = max(zdte(i),1.e-10)
146
        zdte(i) = sign(max(abs(zdte(i)),1.e-10),zdte(i))
147
!
148
        testar(i) = (cdrah(i) * zdte(i) * speed(i))/ustar(i)
149
        qstar(i) = (cdrah(i) * zdq(i) * speed(i))/ustar(i)
150
        lmon(i) = (ustar(i) * ustar(i) * tpot(i))/ &
151
 &                (RKAR * RG * testar(i))
152
      ENDDO
153
!
154
!----------First aproximation of variables at zref --------------------------
155
      zref = 2.0
156
      CALL screenp(klon, knon, nsrf, speed, tpot, q1, &
157
 &                 ts1, qsurf, z0m, lmon, &
158
 &                 ustar, testar, qstar, zref, &
159
 &                 delu, delte, delq)
160
!
161
      DO i = 1, knon
162
        u_zref(i) = delu(i)
163
        q_zref(i) = max(qsurf(i),0.0) + delq(i)
164
        te_zref(i) = ts1(i) + delte(i)
165
        temp(i) = te_zref(i) * (psol(i)/pat1(i))**(-RKAPPA)
166
        q_zref_p(i) = q_zref(i)
167
!       te_zref_p(i) = te_zref(i)
168
        temp_p(i) = temp(i)
169
      ENDDO
170
!
171
! Iteration of the variables at the reference level zref : corrector calculation ; see Hess & McAvaney, 1995
172
!
173
      DO n = 1, niter
174
!
175
        okri=.TRUE.
176
        CALL screenc(klon, knon, nsrf, zxli, &
177
 &                   u_zref, temp, q_zref, zref, &
178
 &                   ts1, qsurf, z0m, z0h, psol, &
179
 &                   ustar, testar, qstar, okri, ri1, &
180
 &                   pref, delu, delte, delq)
181
!
182
        DO i = 1, knon
183
          u_zref(i) = delu(i)
184
          q_zref(i) = delq(i) + max(qsurf(i),0.0)
185
          te_zref(i) = delte(i) + ts1(i)
186
!
187
! return to normal temperature
188
!
189
          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
190
!         temp(i) = te_zref(i) - (zref* RG)/RCPD/ &
191
!                 (1 + RVTMP2 * max(q_zref(i),0.0))
192
!
193
!IM +++
194
!         IF(temp(i).GT.350.) THEN
195
!           WRITE(*,*) 'temp(i) GT 350 K !!',i,nsrf,temp(i)
196
!         ENDIF
197
!IM ---
198
!
199
        IF(n.EQ.ncon) THEN
200
          te_zref_con(i) = te_zref(i)
201
          q_zref_con(i) = q_zref(i)
202
        ENDIF
203
!
204
        ENDDO
205
!
206
      ENDDO
207
!
208
! verifier le critere de convergence : 0.25% pour te_zref et 5% pour qe_zref
209
!
210
!       DO i = 1, knon
211
!         conv_te(i) = (te_zref(i) - te_zref_con(i))/te_zref_con(i)
212
!         conv_q(i) = (q_zref(i) - q_zref_con(i))/q_zref_con(i)
213
!IM +++
214
!         IF(abs(conv_te(i)).GE.0.0025.AND.abs(conv_q(i)).GE.0.05) THEN
215
!           PRINT*,'DIV','i=',i,te_zref_con(i),te_zref(i),conv_te(i), &
216
!           q_zref_con(i),q_zref(i),conv_q(i)
217
!         ENDIF
218
!IM ---
219
!       ENDDO
220
!
221
      DO i = 1, knon
222
        q_zref_c(i) = q_zref(i)
223
        temp_c(i) = temp(i)
224
!
225
!       IF(zri1(i).LT.0.) THEN
226
!         IF(nsrf.EQ.1) THEN
227
!           ok_pred(i)=1.
228
!           ok_corr(i)=0.
229
!         ELSE
230
!           ok_pred(i)=0.
231
!           ok_corr(i)=1.
232
!         ENDIF
233
!       ELSE
234
!         ok_pred(i)=0.
235
!         ok_corr(i)=1.
236
!       ENDIF
237
!
238
        ok_pred(i)=0.
239
        ok_corr(i)=1.
240
!
241
        t_2m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
242
        q_2m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
243
!IM +++
244
!       IF(n.EQ.niter) THEN
245
!       IF(t_2m(i).LT.t1(i).AND.t_2m(i).LT.ts1(i)) THEN
246
!         PRINT*,' BAD t2m LT ',i,nsrf,t_2m(i),t1(i),ts1(i)
247
!       ELSEIF(t_2m(i).GT.t1(i).AND.t_2m(i).GT.ts1(i)) THEN
248
!         PRINT*,' BAD t2m GT ',i,nsrf,t_2m(i),t1(i),ts1(i)
249
!       ENDIF
250
!       ENDIF
251
!IM ---
252
      ENDDO
253
!
254
!
255
!----------First aproximation of variables at zref --------------------------
256
!
257
      zref = 10.0
258
      CALL screenp(klon, knon, nsrf, speed, tpot, q1, &
259
 &                 ts1, qsurf, z0m, lmon, &
260
 &                 ustar, testar, qstar, zref, &
261
 &                 delu, delte, delq)
262
!
263
      DO i = 1, knon
264
        u_zref(i) = delu(i)
265
        q_zref(i) = max(qsurf(i),0.0) + delq(i)
266
        te_zref(i) = ts1(i) + delte(i)
267
        temp(i) = te_zref(i) * (psol(i)/pat1(i))**(-RKAPPA)
268
!       temp(i) = te_zref(i) - (zref* RG)/RCPD/ &
269
!                 (1 + RVTMP2 * max(q_zref(i),0.0))
270
        u_zref_p(i) = u_zref(i)
271
      ENDDO
272
!
273
! Iteration of the variables at the reference level zref : corrector ; see Hess & McAvaney, 1995
274
!
275
      DO n = 1, niter
276
!
277
        okri=.TRUE.
278
        CALL screenc(klon, knon, nsrf, zxli, &
279
 &                   u_zref, temp, q_zref, zref, &
280
 &                   ts1, qsurf, z0m, z0h, psol, &
281
 &                   ustar, testar, qstar, okri, ri1, &
282
 &                   pref, delu, delte, delq)
283
!
284
        DO i = 1, knon
285
          u_zref(i) = delu(i)
286
          q_zref(i) = delq(i) + max(qsurf(i),0.0)
287
          te_zref(i) = delte(i) + ts1(i)
288
          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
289
!         temp(i) = te_zref(i) - (zref* RG)/RCPD/ &
290
!                   (1 + RVTMP2 * max(q_zref(i),0.0))
291
        ENDDO
292
!
293
      ENDDO
294
!
295
      DO i = 1, knon
296
        u_zref_c(i) = u_zref(i)
297
!
298
        u_10m(i) = u_zref_p(i) * ok_pred(i) + u_zref_c(i) * ok_corr(i)
299
!
300
!AM
301
        q_zref_c(i) = q_zref(i)
302
        temp_c(i) = temp(i)
303
        t_10m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
304
        q_10m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
305
!MA
306
      ENDDO
307
!
308
      RETURN
309
      END subroutine stdlevvar
310
!
311
1440
      SUBROUTINE stdlevvarn(klon, knon, nsrf, zxli, &
312
1440
                           u1, v1, t1, q1, z1, &
313
                           ts1, qsurf, z0m, z0h, psol, pat1, &
314
                           t_2m, q_2m, t_10m, q_10m, u_10m, ustar, &
315
1440
                           n2mout)
316
!
317
      USE ioipsl_getin_p_mod, ONLY : getin_p
318
      IMPLICIT NONE
319
!-------------------------------------------------------------------------
320
!
321
! Objet : calcul de la temperature et l'humidite relative a 2m et du
322
!         module du vent a 10m a partir des relations de Dyer-Businger et
323
!         des equations de Louis.
324
!
325
! Reference : Hess, Colman et McAvaney (1995)
326
!
327
! I. Musat, 01.07.2002
328
!
329
!AM On rajoute en sortie t et q a 10m pr le calcule d'hbtm2 dans clmain
330
!
331
!-------------------------------------------------------------------------
332
!
333
! klon----input-I- dimension de la grille physique (= nb_pts_latitude X nb_pts_longitude)
334
! knon----input-I- nombre de points pour un type de surface
335
! nsrf----input-I- indice pour le type de surface; voir indice_sol_mod.F90
336
! zxli----input-L- TRUE si calcul des cdrags selon Laurent Li
337
! u1------input-R- vent zonal au 1er niveau du modele
338
! v1------input-R- vent meridien au 1er niveau du modele
339
! t1------input-R- temperature de l'air au 1er niveau du modele
340
! q1------input-R- humidite relative au 1er niveau du modele
341
! z1------input-R- geopotentiel au 1er niveau du modele
342
! ts1-----input-R- temperature de l'air a la surface
343
! qsurf---input-R- humidite relative a la surface
344
! z0m, z0h---input-R- rugosite
345
! psol----input-R- pression au sol
346
! pat1----input-R- pression au 1er niveau du modele
347
!
348
! t_2m---output-R- temperature de l'air a 2m
349
! q_2m---output-R- humidite relative a 2m
350
! u_2m--output-R- vitesse du vent a 2m
351
! u_10m--output-R- vitesse du vent a 10m
352
! ustar--output-R- u*
353
!AM
354
! t_10m--output-R- temperature de l'air a 10m
355
! q_10m--output-R- humidite specifique a 10m
356
!
357
      INTEGER, intent(in) :: klon, knon, nsrf
358
      LOGICAL, intent(in) :: zxli
359
      REAL, dimension(klon), intent(in) :: u1, v1, t1, q1, z1, ts1
360
      REAL, dimension(klon), intent(in) :: qsurf, z0m, z0h
361
      REAL, dimension(klon), intent(in) :: psol, pat1
362
!
363
      REAL, dimension(klon), intent(out) :: t_2m, q_2m, ustar
364
      REAL, dimension(klon), intent(out) :: u_10m, t_10m, q_10m
365
      INTEGER, dimension(klon, 6), intent(out) :: n2mout
366
!
367
2880
      REAL, dimension(klon) :: u_2m
368
2880
      REAL, dimension(klon) :: cdrm2m, cdrh2m, ri2m
369
2880
      REAL, dimension(klon) :: cdram, cdrah, zri1
370
      REAL, dimension(klon) :: cdmn1, cdhn1, fm1, fh1
371
      REAL, dimension(klon) :: cdmn2m, cdhn2m, fm2m, fh2m
372
      REAL, dimension(klon) :: ri2m_new
373
!-------------------------------------------------------------------------
374
      include "flux_arp.h"
375
      include "YOMCST.h"
376
!IM PLUS
377
      include "YOETHF.h"
378
!
379
! Quelques constantes et options:
380
!
381
! RKAR : constante de von Karman
382
      REAL, PARAMETER :: RKAR=0.40
383
! niter : nombre iterations calcul "corrector"
384
!     INTEGER, parameter :: niter=6, ncon=niter-1
385
!IM 071020     INTEGER, parameter :: niter=2, ncon=niter-1
386
      INTEGER, parameter :: niter=2, ncon=niter
387
!     INTEGER, parameter :: niter=6, ncon=niter
388
!
389
! Variables locales
390
      INTEGER :: i, n
391
      REAL :: zref
392
2880
      REAL, dimension(klon) :: speed
393
! tpot : temperature potentielle
394
2880
      REAL, dimension(klon) :: tpot
395
      REAL, dimension(klon) :: cdran
396
! ri1 : nb. de Richardson entre la surface --> la 1ere couche
397
2880
      REAL, dimension(klon) :: ri1
398
      DOUBLE PRECISION, parameter :: eps=1.0D-20
399
      REAL, dimension(klon) :: delu, delte, delq
400
2880
      REAL, dimension(klon) :: delh, delm
401
2880
      REAL, dimension(klon) :: delh_new, delm_new
402
2880
      REAL, dimension(klon) :: u_zref, te_zref, q_zref
403
      REAL, dimension(klon) :: u_zref_pnew, te_zref_pnew, q_zref_pnew
404
2880
      REAL, dimension(klon) :: temp, pref
405
2880
      REAL, dimension(klon) :: temp_new, pref_new
406
      LOGICAL :: okri
407
2880
      REAL, dimension(klon) :: u_zref_p, te_zref_p, temp_p, q_zref_p
408
      REAL, dimension(klon) :: u_zref_p_new, te_zref_p_new, temp_p_new, q_zref_p_new
409
!convergence
410
2880
      REAL, dimension(klon) :: te_zref_con, q_zref_con
411
2880
      REAL, dimension(klon) :: u_zref_c, te_zref_c, temp_c, q_zref_c
412
2880
      REAL, dimension(klon) :: ok_pred, ok_corr
413
!
414
2880
      REAL, dimension(klon) :: cdrm10m, cdrh10m, ri10m
415
      REAL, dimension(klon) :: cdmn10m, cdhn10m, fm10m, fh10m
416
2880
      REAL, dimension(klon) :: cdn2m, cdn1, zri_zero
417
      REAL :: CEPDUE,zdu2
418
      INTEGER :: nzref, nz1
419
2880
      LOGICAL, dimension(klon) :: ok_t2m_toosmall, ok_t2m_toobig
420
2880
      LOGICAL, dimension(klon) :: ok_q2m_toosmall, ok_q2m_toobig
421
2880
      LOGICAL, dimension(klon) :: ok_u2m_toobig
422
2880
      LOGICAL, dimension(klon) :: ok_t10m_toosmall, ok_t10m_toobig
423
2880
      LOGICAL, dimension(klon) :: ok_q10m_toosmall, ok_q10m_toobig
424
2880
      LOGICAL, dimension(klon) :: ok_u10m_toobig
425
1440
      INTEGER, dimension(klon, 6) :: n10mout
426
427
!-------------------------------------------------------------------------
428
      CEPDUE=0.1
429
!
430
! n2mout : compteur des pas de temps ou t2m,q2m ou u2m sont en dehors des intervalles
431
!          [tsurf, temp], [qsurf, q1] ou [0, speed]
432
! n10mout : compteur des pas de temps ou t10m,q10m ou u10m sont en dehors des intervalles
433
!          [tsurf, temp], [qsurf, q1] ou [0, speed]
434
!
435

8598240
      n2mout(:,:)=0
436

8598240
      n10mout(:,:)=0
437
438
622850
      DO i=1, knon
439
621410
       speed(i)=MAX(SQRT(u1(i)**2+v1(i)**2),CEPDUE)
440
622850
       ri1(i) = 0.0
441
      ENDDO
442
!
443
1440
      okri=.FALSE.
444
      CALL cdrag(knon, nsrf, &
445
 &                   speed, t1, q1, z1, &
446
 &                   psol, ts1, qsurf, z0m, z0h, &
447
 &                   zri_zero, 0, &
448
1440
 &                   cdram, cdrah, zri1, pref)
449
450
!
451
622850
      DO i = 1, knon
452
621410
        ri1(i) = zri1(i)
453
621410
        tpot(i) = t1(i)* (psol(i)/pat1(i))**RKAPPA
454
621410
        zdu2 = MAX(CEPDUE*CEPDUE, speed(i)**2)
455
622850
        ustar(i) = sqrt(cdram(i) * zdu2)
456
!
457
      ENDDO
458
!
459
!----------First aproximation of variables at zref --------------------------
460
1440
      zref = 2.0
461
!
462
! calcul first-guess en utilisant dans les calculs à 2m
463
! le Richardson de la premiere couche atmospherique
464
!
465
       CALL screencn(klon, knon, nsrf, zxli, &
466
 &                   speed, tpot, q1, zref, &
467
 &                   ts1, qsurf, z0m, z0h, psol, &
468
 &                   cdram, cdrah,  okri, &
469
 &                   ri1, 1, &
470
1440
 &                   pref_new, delm_new, delh_new, ri2m)
471
!
472
622850
       DO i = 1, knon
473
621410
         u_zref(i) = delm_new(i)*speed(i)
474
621410
         u_zref_p(i) = u_zref(i)
475
         q_zref(i) = delh_new(i)*max(q1(i),0.0) + &
476
621410
         &           max(qsurf(i),0.0)*(1-delh_new(i))
477
621410
         q_zref_p(i) = q_zref(i)
478
621410
         te_zref(i) = delh_new(i)*tpot(i) + ts1(i)*(1-delh_new(i))
479
621410
         te_zref_p(i) = te_zref(i)
480
!
481
! return to normal temperature
482
621410
         temp(i) = te_zref(i) * (psol(i)/pref_new(i))**(-RKAPPA)
483
621410
         temp_p(i) = temp(i)
484
!
485
! compteurs ici
486
!
487
         ok_t2m_toosmall(i)=te_zref(i).LT.tpot(i).AND. &
488

621410
         & te_zref(i).LT.ts1(i)
489
         ok_t2m_toobig(i)=te_zref(i).GT.tpot(i).AND. &
490

621410
         & te_zref(i).GT.ts1(i)
491
         ok_q2m_toosmall(i)=q_zref(i).LT.q1(i).AND. &
492

621410
         & q_zref(i).LT.qsurf(i)
493
         ok_q2m_toobig(i)=q_zref(i).GT.q1(i).AND. &
494

621410
         & q_zref(i).GT.qsurf(i)
495
621410
         ok_u2m_toobig(i)=u_zref(i).GT.speed(i)
496
!
497

621410
         IF(ok_t2m_toosmall(i).OR.ok_t2m_toobig(i)) THEN
498
             n2mout(i,1)=n2mout(i,1)+1
499
         ENDIF
500

621410
         IF(ok_q2m_toosmall(i).OR.ok_q2m_toobig(i)) THEN
501
1611
             n2mout(i,3)=n2mout(i,3)+1
502
         ENDIF
503
621410
         IF(ok_u2m_toobig(i)) THEN
504
             n2mout(i,5)=n2mout(i,5)+1
505
         ENDIF
506
!
507
         IF(ok_t2m_toosmall(i).OR.ok_t2m_toobig(i).OR. &
508


621410
          & ok_q2m_toosmall(i).OR.ok_q2m_toobig(i).OR. &
509
1440
          & ok_u2m_toobig(i)) THEN
510
1611
             delm_new(i)=min(max(delm_new(i),0.),1.)
511
1611
             delh_new(i)=min(max(delh_new(i),0.),1.)
512
1611
             u_zref(i) = delm_new(i)*speed(i)
513
1611
             u_zref_p(i) = u_zref(i)
514
             q_zref(i) = delh_new(i)*max(q1(i),0.0) + &
515
1611
         &               max(qsurf(i),0.0)*(1-delh_new(i))
516
1611
             q_zref_p(i) = q_zref(i)
517
1611
             te_zref(i) = delh_new(i)*tpot(i) + ts1(i)*(1-delh_new(i))
518
1611
             te_zref_p(i) = te_zref(i)
519
!
520
! return to normal temperature
521
1611
             temp(i) = te_zref(i) * (psol(i)/pref_new(i))**(-RKAPPA)
522
1611
             temp_p(i) = temp(i)
523
         ENDIF
524
!
525
       ENDDO
526
!
527
! Iteration of the variables at the reference level zref : corrector calculation ; see Hess & McAvaney, 1995
528
!
529
4320
      DO n = 1, niter
530
!
531
2880
        okri=.TRUE.
532
        CALL screencn(klon, knon, nsrf, zxli, &
533
 &                   u_zref, temp, q_zref, zref, &
534
 &                   ts1, qsurf, z0m, z0h, psol, &
535
 &                   cdram, cdrah,  okri, &
536
 &                   ri1, 0, &
537
2880
 &                   pref, delm, delh, ri2m)
538
!
539
1247140
        DO i = 1, knon
540
1242820
          u_zref(i) = delm(i)*speed(i)
541
          q_zref(i) = delh(i)*max(q1(i),0.0) + &
542
1242820
          &           max(qsurf(i),0.0)*(1-delh(i))
543
1242820
          te_zref(i) = delh(i)*tpot(i) + ts1(i)*(1-delh(i))
544
!
545
! return to normal temperature
546
1242820
          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
547
!
548
! compteurs ici
549
!
550
          ok_t2m_toosmall(i)=te_zref(i).LT.tpot(i).AND. &
551

1242820
          & te_zref(i).LT.ts1(i)
552
          ok_t2m_toobig(i)=te_zref(i).GT.tpot(i).AND. &
553

1242820
          & te_zref(i).GT.ts1(i)
554
          ok_q2m_toosmall(i)=q_zref(i).LT.q1(i).AND. &
555

1242820
          & q_zref(i).LT.qsurf(i)
556
          ok_q2m_toobig(i)=q_zref(i).GT.q1(i).AND. &
557

1242820
          & q_zref(i).GT.qsurf(i)
558
1242820
          ok_u2m_toobig(i)=u_zref(i).GT.speed(i)
559
!
560

1242820
          IF(ok_t2m_toosmall(i).OR.ok_t2m_toobig(i)) THEN
561
407
              n2mout(i,2)=n2mout(i,2)+1
562
          ENDIF
563

1242820
          IF(ok_q2m_toosmall(i).OR.ok_q2m_toobig(i)) THEN
564
3615
              n2mout(i,4)=n2mout(i,4)+1
565
          ENDIF
566
1242820
          IF(ok_u2m_toobig(i)) THEN
567
274
              n2mout(i,6)=n2mout(i,6)+1
568
          ENDIF
569
!
570
          IF(ok_t2m_toosmall(i).OR.ok_t2m_toobig(i).OR. &
571


1242820
           & ok_q2m_toosmall(i).OR.ok_q2m_toobig(i).OR. &
572
           & ok_u2m_toobig(i)) THEN
573
3618
              delm(i)=min(max(delm(i),0.),1.)
574
3618
              delh(i)=min(max(delh(i),0.),1.)
575
3618
              u_zref(i) = delm(i)*speed(i)
576
              q_zref(i) = delh(i)*max(q1(i),0.0) + &
577
3618
          &           max(qsurf(i),0.0)*(1-delh(i))
578
3618
              te_zref(i) = delh(i)*tpot(i) + ts1(i)*(1-delh(i))
579
3618
              temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
580
          ENDIF
581
!
582
!
583
1245700
          IF(n.EQ.ncon) THEN
584
621410
            te_zref_con(i) = te_zref(i)
585
621410
            q_zref_con(i) = q_zref(i)
586
          ENDIF
587
!
588
        ENDDO
589
!
590
      ENDDO
591
!
592
622850
      DO i = 1, knon
593
621410
        q_zref_c(i) = q_zref(i)
594
621410
        temp_c(i) = temp(i)
595
!
596
621410
        ok_pred(i)=0.
597
621410
        ok_corr(i)=1.
598
!
599
621410
        t_2m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
600
621410
        q_2m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
601
!
602
621410
        u_zref_c(i) = u_zref(i)
603
622850
        u_2m(i) = u_zref_p(i) * ok_pred(i) + u_zref_c(i) * ok_corr(i)
604
      ENDDO
605
!
606
!
607
!----------First aproximation of variables at zref --------------------------
608
!
609
1440
      zref = 10.0
610
!
611
       CALL screencn(klon, knon, nsrf, zxli, &
612
 &                   speed, tpot, q1, zref, &
613
 &                   ts1, qsurf, z0m, z0h, psol, &
614
 &                   cdram, cdrah,  okri, &
615
 &                   ri1, 1, &
616
1440
 &                   pref_new, delm_new, delh_new, ri10m)
617
!
618
622850
       DO i = 1, knon
619
621410
         u_zref(i) = delm_new(i)*speed(i)
620
         q_zref(i) = delh_new(i)*max(q1(i),0.0) + &
621
621410
         &           max(qsurf(i),0.0)*(1-delh_new(i))
622
621410
         te_zref(i) = delh_new(i)*tpot(i) + ts1(i)*(1-delh_new(i))
623
621410
         temp(i) = te_zref(i) * (psol(i)/pref_new(i))**(-RKAPPA)
624
621410
         u_zref_p(i) = u_zref(i)
625
!
626
! compteurs ici
627
!
628
         ok_t10m_toosmall(i)=te_zref(i).LT.tpot(i).AND. &
629

621410
         & te_zref(i).LT.ts1(i)
630
         ok_t10m_toobig(i)=te_zref(i).GT.tpot(i).AND. &
631

621410
         & te_zref(i).GT.ts1(i)
632
         ok_q10m_toosmall(i)=q_zref(i).LT.q1(i).AND. &
633

621410
         & q_zref(i).LT.qsurf(i)
634
         ok_q10m_toobig(i)=q_zref(i).GT.q1(i).AND. &
635

621410
         & q_zref(i).GT.qsurf(i)
636
621410
         ok_u10m_toobig(i)=u_zref(i).GT.speed(i)
637
!
638

621410
         IF(ok_t10m_toosmall(i).OR.ok_t10m_toobig(i)) THEN
639
             n10mout(i,1)=n10mout(i,1)+1
640
         ENDIF
641

621410
         IF(ok_q10m_toosmall(i).OR.ok_q10m_toobig(i)) THEN
642
775
             n10mout(i,3)=n10mout(i,3)+1
643
         ENDIF
644
621410
         IF(ok_u10m_toobig(i)) THEN
645
             n10mout(i,5)=n10mout(i,5)+1
646
         ENDIF
647
!
648
         IF(ok_t10m_toosmall(i).OR.ok_t10m_toobig(i).OR. &
649


621410
          & ok_q10m_toosmall(i).OR.ok_q10m_toobig(i).OR. &
650
1440
          & ok_u10m_toobig(i)) THEN
651
775
             delm_new(i)=min(max(delm_new(i),0.),1.)
652
775
             delh_new(i)=min(max(delh_new(i),0.),1.)
653
775
             u_zref(i) = delm_new(i)*speed(i)
654
775
             u_zref_p(i) = u_zref(i)
655
             q_zref(i) = delh_new(i)*max(q1(i),0.0) + &
656
775
         &               max(qsurf(i),0.0)*(1-delh_new(i))
657
775
             te_zref(i) = delh_new(i)*tpot(i) + ts1(i)*(1-delh_new(i))
658
775
             temp(i) = te_zref(i) * (psol(i)/pref_new(i))**(-RKAPPA)
659
         ENDIF
660
!
661
       ENDDO
662
!
663
! Iteration of the variables at the reference level zref : corrector calculation ; see Hess & McAvaney, 1995
664
!
665
4320
      DO n = 1, niter
666
!
667
2880
        okri=.TRUE.
668
        CALL screencn(klon, knon, nsrf, zxli, &
669
 &                   u_zref, temp, q_zref, zref, &
670
 &                   ts1, qsurf, z0m, z0h, psol, &
671
 &                   cdram, cdrah,  okri, &
672
 &                   ri1, 0, &
673
2880
 &                   pref, delm, delh, ri10m)
674
!
675
1247140
        DO i = 1, knon
676
1242820
          u_zref(i) = delm(i)*speed(i)
677
          q_zref(i) = delh(i)*max(q1(i),0.0) + &
678
1242820
          &           max(qsurf(i),0.0)*(1-delh(i))
679
1242820
          te_zref(i) = delh(i)*tpot(i) + ts1(i)*(1-delh(i))
680
!
681
! return to normal temperature
682
1242820
          temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
683
!
684
! compteurs ici
685
!
686
          ok_t10m_toosmall(i)=te_zref(i).LT.tpot(i).AND. &
687

1242820
          & te_zref(i).LT.ts1(i)
688
          ok_t10m_toobig(i)=te_zref(i).GT.tpot(i).AND. &
689

1242820
          & te_zref(i).GT.ts1(i)
690
          ok_q10m_toosmall(i)=q_zref(i).LT.q1(i).AND. &
691

1242820
          & q_zref(i).LT.qsurf(i)
692
          ok_q10m_toobig(i)=q_zref(i).GT.q1(i).AND. &
693

1242820
          & q_zref(i).GT.qsurf(i)
694
1242820
          ok_u10m_toobig(i)=u_zref(i).GT.speed(i)
695
!
696

1242820
          IF(ok_t10m_toosmall(i).OR.ok_t10m_toobig(i)) THEN
697
511
              n10mout(i,2)=n10mout(i,2)+1
698
          ENDIF
699

1242820
          IF(ok_q10m_toosmall(i).OR.ok_q10m_toobig(i)) THEN
700
2278
              n10mout(i,4)=n10mout(i,4)+1
701
          ENDIF
702
1242820
          IF(ok_u10m_toobig(i)) THEN
703
292
              n10mout(i,6)=n10mout(i,6)+1
704
          ENDIF
705
!
706
          IF(ok_t10m_toosmall(i).OR.ok_t10m_toobig(i).OR. &
707


1242820
           & ok_q10m_toosmall(i).OR.ok_q10m_toobig(i).OR. &
708
           & ok_u10m_toobig(i)) THEN
709
2290
              delm(i)=min(max(delm(i),0.),1.)
710
2290
              delh(i)=min(max(delh(i),0.),1.)
711
2290
              u_zref(i) = delm(i)*speed(i)
712
              q_zref(i) = delh(i)*max(q1(i),0.0) + &
713
2290
          &           max(qsurf(i),0.0)*(1-delh(i))
714
2290
              te_zref(i) = delh(i)*tpot(i) + ts1(i)*(1-delh(i))
715
2290
              temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
716
          ENDIF
717
!
718
!
719
1245700
          IF(n.EQ.ncon) THEN
720
621410
            te_zref_con(i) = te_zref(i)
721
621410
            q_zref_con(i) = q_zref(i)
722
          ENDIF
723
!
724
        ENDDO
725
!
726
      ENDDO
727
!
728
622850
      DO i = 1, knon
729
621410
        q_zref_c(i) = q_zref(i)
730
621410
        temp_c(i) = temp(i)
731
!
732
621410
        ok_pred(i)=0.
733
621410
        ok_corr(i)=1.
734
!
735
621410
        t_10m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
736
621410
        q_10m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
737
!
738
621410
        u_zref_c(i) = u_zref(i)
739
622850
        u_10m(i) = u_zref_p(i) * ok_pred(i) + u_zref_c(i) * ok_corr(i)
740
      ENDDO
741
!
742
1440
      RETURN
743
      END subroutine stdlevvarn
744
745
END MODULE stdlevvar_mod