/[lmdze]/trunk/Sources/phylmd/stdlevvar.f
ViewVC logotype

Diff of /trunk/Sources/phylmd/stdlevvar.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

trunk/libf/phylmd/stdlevvar.f90 revision 38 by guez, Thu Jan 6 17:52:19 2011 UTC trunk/phylmd/stdlevvar.f revision 108 by guez, Tue Sep 16 14:00:41 2014 UTC
# Line 1  Line 1 
1  !  module stdlevvar_m
2  ! $Header: /home/cvsroot/LMDZ4/libf/phylmd/stdlevvar.F90,v 1.3 2005/05/25 13:10:09 fairhead Exp $  
3  !   IMPLICIT NONE
4        SUBROUTINE stdlevvar(klon, knon, nsrf, zxli, &  
5                             u1, v1, t1, q1, z1, &  contains
6                             ts1, qsurf, rugos, psol, pat1, &  
7                             t_2m, q_2m, t_10m, q_10m, u_10m, ustar)    SUBROUTINE stdlevvar(klon, knon, nsrf, zxli, u1, v1, t1, q1, z1, ts1, &
8        use SUPHEC_M         qsurf, rugos, psol, pat1, t_2m, q_2m, t_10m, q_10m, u_10m, ustar)
9              use yoethf_m  
10        IMPLICIT NONE      ! From LMDZ4/libf/phylmd/stdlevvar.F90, version 1.3 2005/05/25 13:10:09
11  !-------------------------------------------------------------------------  
12  !      use coefcdrag_m, only: coefcdrag
13  ! Objet : calcul de la temperature et l'humidite relative a 2m et du      USE suphec_m, ONLY: rg, rkappa
14  !         module du vent a 10m a partir des relations de Dyer-Businger et  
15  !         des equations de Louis.      ! Objet : calcul de la température et de l'humidité relative à 2 m
16  !      ! et du module du vent à 10 m à partir des relations de
17  ! Reference : Hess, Colman et McAvaney (1995)              ! Dyer-Businger et des équations de Louis.
18  !  
19  ! I. Musat, 01.07.2002      ! Reference: Hess, Colman and McAvaney (1995)
20  !  
21  !AM On rajoute en sortie t et q a 10m pr le calcule d'hbtm2 dans clmain      ! Author: I. Musat, 01.07.2002
22  !  
23  !-------------------------------------------------------------------------      INTEGER, intent(in):: klon
24  !      ! dimension de la grille physique (= nb_pts_latitude X nb_pts_longitude)
25  ! klon----input-I- dimension de la grille physique (= nb_pts_latitude X nb_pts_longitude)  
26  ! knon----input-I- nombre de points pour un type de surface      INTEGER, intent(in):: knon
27  ! nsrf----input-I- indice pour le type de surface; voir indicesol.inc      ! knon----input-I- nombre de points pour un type de surface
28  ! zxli----input-L- TRUE si calcul des cdrags selon Laurent Li      INTEGER, intent(in):: nsrf
29  ! u1------input-R- vent zonal au 1er niveau du modele      ! nsrf----input-I- indice pour le type de surface; voir indicesol.inc
30  ! v1------input-R- vent meridien au 1er niveau du modele      LOGICAL, intent(in):: zxli
31  ! t1------input-R- temperature de l'air au 1er niveau du modele      ! zxli----input-L- TRUE si calcul des cdrags selon Laurent Li
32  ! q1------input-R- humidite relative au 1er niveau du modele      REAL, dimension(klon), intent(in):: u1
33  ! z1------input-R- geopotentiel au 1er niveau du modele      ! u1------input-R- vent zonal au 1er niveau du modele
34  ! ts1-----input-R- temperature de l'air a la surface      REAL, dimension(klon), intent(in):: v1
35  ! qsurf---input-R- humidite relative a la surface      ! v1------input-R- vent meridien au 1er niveau du modele
36  ! rugos---input-R- rugosite      REAL, dimension(klon), intent(in):: t1
37  ! psol----input-R- pression au sol      ! t1------input-R- temperature de l'air au 1er niveau du modele
38  ! pat1----input-R- pression au 1er niveau du modele      REAL, dimension(klon), intent(in):: q1
39  !      ! q1------input-R- humidite relative au 1er niveau du modele
40  ! t_2m---output-R- temperature de l'air a 2m      REAL, dimension(klon), intent(in):: z1
41  ! q_2m---output-R- humidite relative a 2m      ! z1------input-R- geopotentiel au 1er niveau du modele
42  ! u_10m--output-R- vitesse du vent a 10m      REAL, dimension(klon), intent(in):: ts1
43  !AM      ! ts1-----input-R- temperature de l'air a la surface
44  ! t_10m--output-R- temperature de l'air a 10m      REAL, dimension(klon), intent(in):: qsurf
45  ! q_10m--output-R- humidite specifique a 10m       ! qsurf---input-R- humidite relative a la surface
46  ! ustar--output-R- u*      REAL, dimension(klon), intent(in):: rugos
47  !      ! rugos---input-R- rugosite
48        INTEGER, intent(in) :: klon, knon, nsrf     REAL, dimension(klon), intent(in):: psol
49        LOGICAL, intent(in) :: zxli      ! psol----input-R- pression au sol
50        REAL, dimension(klon), intent(in) :: u1, v1, t1, q1, z1, ts1     REAL, dimension(klon), intent(in):: pat1
51        REAL, dimension(klon), intent(in) :: qsurf, rugos      ! pat1----input-R- pression au 1er niveau du modele
52        REAL, dimension(klon), intent(in) :: psol, pat1  
53  !      REAL, dimension(klon), intent(out):: t_2m
54        REAL, dimension(klon), intent(out) :: t_2m, q_2m, ustar      ! t_2m---output-R- temperature de l'air a 2m
55        REAL, dimension(klon), intent(out) :: u_10m, t_10m, q_10m      REAL, dimension(klon), intent(out):: q_2m
56  !-------------------------------------------------------------------------      ! q_2m---output-R- humidite relative a 2m
57  !IM PLUS      REAL, dimension(klon), intent(out):: t_10m
58  !      ! t_10m--output-R- temperature de l'air a 10m
59  ! Quelques constantes et options:      REAL, dimension(klon), intent(out):: q_10m
60  !      ! q_10m--output-R- humidite specifique a 10m
61  ! RKAR : constante de von Karman      REAL, dimension(klon), intent(out):: u_10m
62        REAL, PARAMETER :: RKAR=0.40      ! u_10m--output-R- vitesse du vent a 10m
63  ! niter : nombre iterations calcul "corrector"      REAL, intent(out):: ustar(klon) ! u*
64  !     INTEGER, parameter :: niter=6, ncon=niter-1  
65        INTEGER, parameter :: niter=2, ncon=niter-1      ! Local:
66  !  
67  ! Variables locales      ! RKAR : constante de von Karman
68        INTEGER :: i, n      REAL, PARAMETER:: RKAR=0.40
69        REAL :: zref      ! niter : nombre iterations calcul "corrector"
70        REAL, dimension(klon) :: speed      INTEGER, parameter:: niter=2, ncon=niter-1
71  ! tpot : temperature potentielle  
72        REAL, dimension(klon) :: tpot      ! Variables locales
73        REAL, dimension(klon) :: zri1, cdran      INTEGER i, n
74        REAL, dimension(klon) :: cdram, cdrah      REAL zref
75  ! ri1 : nb. de Richardson entre la surface --> la 1ere couche      REAL, dimension(klon):: speed
76        REAL, dimension(klon) :: ri1      ! tpot : temperature potentielle
77        REAL, dimension(klon) :: testar, qstar      REAL, dimension(klon):: tpot
78        REAL, dimension(klon) :: zdte, zdq        REAL, dimension(klon):: zri1, cdran
79  ! lmon : longueur de Monin-Obukhov selon Hess, Colman and McAvaney      REAL cdram(klon), cdrah(klon)
80        DOUBLE PRECISION, dimension(klon) :: lmon      ! ri1 : nb. de Richardson entre la surface --> la 1ere couche
81        DOUBLE PRECISION, parameter :: eps=1.0D-20      REAL, dimension(klon):: ri1
82        REAL, dimension(klon) :: delu, delte, delq      REAL, dimension(klon):: testar, qstar
83        REAL, dimension(klon) :: u_zref, te_zref, q_zref        REAL, dimension(klon):: zdte, zdq
84        REAL, dimension(klon) :: temp, pref      ! lmon : longueur de Monin-Obukhov selon Hess, Colman and McAvaney
85        LOGICAL :: okri      DOUBLE PRECISION, dimension(klon):: lmon
86        REAL, dimension(klon) :: u_zref_p, te_zref_p, temp_p, q_zref_p      DOUBLE PRECISION, parameter:: eps=1.0D-20
87  !convertgence      REAL, dimension(klon):: delu, delte, delq
88        REAL, dimension(klon) :: te_zref_con, q_zref_con      REAL, dimension(klon):: u_zref, te_zref, q_zref
89        REAL, dimension(klon) :: u_zref_c, te_zref_c, temp_c, q_zref_c      REAL, dimension(klon):: temp, pref
90        REAL, dimension(klon) :: ok_pred, ok_corr      LOGICAL okri
91  !     REAL, dimension(klon) :: conv_te, conv_q      REAL, dimension(klon):: u_zref_p, temp_p, q_zref_p
92  !-------------------------------------------------------------------------      !convertgence
93        DO i=1, knon      REAL, dimension(klon):: te_zref_con, q_zref_con
94        REAL, dimension(klon):: u_zref_c, temp_c, q_zref_c
95        REAL, dimension(klon):: ok_pred, ok_corr
96    
97        !-------------------------------------------------------------------------
98    
99        DO i=1, knon
100         speed(i)=SQRT(u1(i)**2+v1(i)**2)         speed(i)=SQRT(u1(i)**2+v1(i)**2)
101         ri1(i) = 0.0         ri1(i) = 0.0
102        ENDDO      ENDDO
103  !  
104        okri=.FALSE.      okri=.FALSE.
105        CALL coefcdrag(klon, knon, nsrf, zxli, &      CALL coefcdrag(klon, knon, nsrf, zxli, speed, t1, q1, z1, psol, ts1, &
106   &                   speed, t1, q1, z1, psol, &           qsurf, rugos, okri, ri1, cdram, cdrah, cdran, zri1, pref)
107   &                   ts1, qsurf, rugos, okri, ri1,  &          
108   &                   cdram, cdrah, cdran, zri1, pref)                  ! Star variables
109  !  
110  !---------Star variables----------------------------------------------------      DO i = 1, knon
111  !         ri1(i) = zri1(i)
112        DO i = 1, knon         tpot(i) = t1(i)* (psol(i)/pat1(i))**RKAPPA
113          ri1(i) = zri1(i)         ustar(i) = sqrt(cdram(i) * speed(i) * speed(i))
114          tpot(i) = t1(i)* (psol(i)/pat1(i))**RKAPPA         zdte(i) = tpot(i) - ts1(i)
115          ustar(i) = sqrt(cdram(i) * speed(i) * speed(i))         zdq(i) = max(q1(i), 0.0) - max(qsurf(i), 0.0)
116          zdte(i) = tpot(i) - ts1(i)  
117          zdq(i) = max(q1(i),0.0) - max(qsurf(i),0.0)         zdte(i) = sign(max(abs(zdte(i)), 1.e-10), zdte(i))
118  !  
119  !         testar(i) = (cdrah(i) * zdte(i) * speed(i))/ustar(i)
120  !IM BUG BUG BUG       zdte(i) = max(zdte(i),1.e-10)         qstar(i) = (cdrah(i) * zdq(i) * speed(i))/ustar(i)
121          zdte(i) = sign(max(abs(zdte(i)),1.e-10),zdte(i))         lmon(i) = (ustar(i) * ustar(i) * tpot(i))/ &
122  !              (RKAR * RG * testar(i))
123          testar(i) = (cdrah(i) * zdte(i) * speed(i))/ustar(i)      ENDDO
124          qstar(i) = (cdrah(i) * zdq(i) * speed(i))/ustar(i)  
125          lmon(i) = (ustar(i) * ustar(i) * tpot(i))/ &      ! First aproximation of variables at zref  
126   &                (RKAR * RG * testar(i))      zref = 2.0
127        ENDDO      CALL screenp(klon, knon, nsrf, speed, tpot, q1, &
128  !           ts1, qsurf, rugos, lmon, &
129  !----------First aproximation of variables at zref --------------------------           ustar, testar, qstar, zref, &
130        zref = 2.0           delu, delte, delq)
131        CALL screenp(klon, knon, nsrf, speed, tpot, q1, &  
132   &                 ts1, qsurf, rugos, lmon, &      DO i = 1, knon
133   &                 ustar, testar, qstar, zref, &         u_zref(i) = delu(i)
134   &                 delu, delte, delq)         q_zref(i) = max(qsurf(i), 0.0) + delq(i)
135  !         te_zref(i) = ts1(i) + delte(i)
136        DO i = 1, knon         temp(i) = te_zref(i) * (psol(i)/pat1(i))**(-RKAPPA)
137          u_zref(i) = delu(i)         q_zref_p(i) = q_zref(i)
138          q_zref(i) = max(qsurf(i),0.0) + delq(i)         temp_p(i) = temp(i)
139          te_zref(i) = ts1(i) + delte(i)      ENDDO
140          temp(i) = te_zref(i) * (psol(i)/pat1(i))**(-RKAPPA)  
141          q_zref_p(i) = q_zref(i)      ! Iteration of the variables at the reference level zref :
142  !       te_zref_p(i) = te_zref(i)      ! corrector calculation ; see Hess & McAvaney, 1995
143          temp_p(i) = temp(i)  
144        ENDDO      DO n = 1, niter
145  !         okri=.TRUE.
146  ! Iteration of the variables at the reference level zref : corrector calculation ; see Hess & McAvaney, 1995         CALL screenc(klon, knon, nsrf, zxli, &
147  !              u_zref, temp, q_zref, zref, &
148        DO n = 1, niter              ts1, qsurf, rugos, psol, &
149  !              ustar, testar, qstar, okri, ri1, &
150          okri=.TRUE.              pref, delu, delte, delq)
151          CALL screenc(klon, knon, nsrf, zxli, &  
152   &                   u_zref, temp, q_zref, zref, &         DO i = 1, knon
  &                   ts1, qsurf, rugos, psol, &            
  &                   ustar, testar, qstar, okri, ri1, &  
  &                   pref, delu, delte, delq)  
 !  
         DO i = 1, knon  
153            u_zref(i) = delu(i)            u_zref(i) = delu(i)
154            q_zref(i) = delq(i) + max(qsurf(i),0.0)            q_zref(i) = delq(i) + max(qsurf(i), 0.0)
155            te_zref(i) = delte(i) + ts1(i)            te_zref(i) = delte(i) + ts1(i)
156  !  
157  ! return to normal temperature            ! return to normal temperature
158  !  
159            temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)            temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
160  !         temp(i) = te_zref(i) - (zref* RG)/RCPD/ &  
161  !                 (1 + RVTMP2 * max(q_zref(i),0.0))            IF(n == ncon) THEN
162  !               te_zref_con(i) = te_zref(i)
163  !IM +++               q_zref_con(i) = q_zref(i)
164  !         IF(temp(i).GT.350.) THEN            ENDIF
165  !           WRITE(*,*) 'temp(i) GT 350 K !!',i,nsrf,temp(i)         ENDDO
166  !         ENDIF      ENDDO
167  !IM ---  
168  !      ! verifier le critere de convergence : 0.25% pour te_zref et 5% pour qe_zref
169          IF(n.EQ.ncon) THEN  
170            te_zref_con(i) = te_zref(i)      DO i = 1, knon
171            q_zref_con(i) = q_zref(i)         q_zref_c(i) = q_zref(i)
172          ENDIF         temp_c(i) = temp(i)
173  !  
174          ENDDO         ok_pred(i)=0.
175  !         ok_corr(i)=1.
176        ENDDO  
177  !         t_2m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
178  ! verifier le critere de convergence : 0.25% pour te_zref et 5% pour qe_zref         q_2m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
179  !      ENDDO
180  !       DO i = 1, knon  
181  !         conv_te(i) = (te_zref(i) - te_zref_con(i))/te_zref_con(i)      ! First aproximation of variables at zref  
182  !         conv_q(i) = (q_zref(i) - q_zref_con(i))/q_zref_con(i)  
183  !IM +++      zref = 10.0
184  !         IF(abs(conv_te(i)).GE.0.0025.AND.abs(conv_q(i)).GE.0.05) THEN      CALL screenp(klon, knon, nsrf, speed, tpot, q1, &
185  !           PRINT*,'DIV','i=',i,te_zref_con(i),te_zref(i),conv_te(i), &           ts1, qsurf, rugos, lmon, &
186  !           q_zref_con(i),q_zref(i),conv_q(i)           ustar, testar, qstar, zref, &
187  !         ENDIF           delu, delte, delq)
188  !IM ---  
189  !       ENDDO      DO i = 1, knon
190  !         u_zref(i) = delu(i)
191        DO i = 1, knon         q_zref(i) = max(qsurf(i), 0.0) + delq(i)
192          q_zref_c(i) = q_zref(i)         te_zref(i) = ts1(i) + delte(i)
193          temp_c(i) = temp(i)         temp(i) = te_zref(i) * (psol(i)/pat1(i))**(-RKAPPA)
194  !         u_zref_p(i) = u_zref(i)
195  !       IF(zri1(i).LT.0.) THEN      ENDDO
196  !         IF(nsrf.EQ.1) THEN  
197  !           ok_pred(i)=1.      ! Iteration of the variables at the reference level zref:
198  !           ok_corr(i)=0.      ! corrector ; see Hess & McAvaney, 1995
199  !         ELSE  
200  !           ok_pred(i)=0.      DO n = 1, niter
201  !           ok_corr(i)=1.         okri=.TRUE.
202  !         ENDIF         CALL screenc(klon, knon, nsrf, zxli, &
203  !       ELSE              u_zref, temp, q_zref, zref, &
204  !         ok_pred(i)=0.              ts1, qsurf, rugos, psol, &
205  !         ok_corr(i)=1.              ustar, testar, qstar, okri, ri1, &
206  !       ENDIF              pref, delu, delte, delq)
207  !  
208          ok_pred(i)=0.         DO i = 1, knon
         ok_corr(i)=1.  
 !  
         t_2m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)  
         q_2m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)  
 !IM +++  
 !       IF(n.EQ.niter) THEN  
 !       IF(t_2m(i).LT.t1(i).AND.t_2m(i).LT.ts1(i)) THEN  
 !         PRINT*,' BAD t2m LT ',i,nsrf,t_2m(i),t1(i),ts1(i)  
 !       ELSEIF(t_2m(i).GT.t1(i).AND.t_2m(i).GT.ts1(i)) THEN  
 !         PRINT*,' BAD t2m GT ',i,nsrf,t_2m(i),t1(i),ts1(i)  
 !       ENDIF  
 !       ENDIF  
 !IM ---  
       ENDDO  
 !  
 !  
 !----------First aproximation of variables at zref --------------------------  
 !  
       zref = 10.0  
       CALL screenp(klon, knon, nsrf, speed, tpot, q1, &  
  &                 ts1, qsurf, rugos, lmon, &  
  &                 ustar, testar, qstar, zref, &  
  &                 delu, delte, delq)  
 !  
       DO i = 1, knon  
         u_zref(i) = delu(i)  
         q_zref(i) = max(qsurf(i),0.0) + delq(i)  
         te_zref(i) = ts1(i) + delte(i)  
         temp(i) = te_zref(i) * (psol(i)/pat1(i))**(-RKAPPA)  
 !       temp(i) = te_zref(i) - (zref* RG)/RCPD/ &  
 !                 (1 + RVTMP2 * max(q_zref(i),0.0))  
         u_zref_p(i) = u_zref(i)  
       ENDDO  
 !  
 ! Iteration of the variables at the reference level zref : corrector ; see Hess & McAvaney, 1995  
 !  
       DO n = 1, niter  
 !  
         okri=.TRUE.  
         CALL screenc(klon, knon, nsrf, zxli, &  
  &                   u_zref, temp, q_zref, zref, &  
  &                   ts1, qsurf, rugos, psol, &  
  &                   ustar, testar, qstar, okri, ri1, &  
  &                   pref, delu, delte, delq)  
 !  
         DO i = 1, knon  
209            u_zref(i) = delu(i)            u_zref(i) = delu(i)
210            q_zref(i) = delq(i) + max(qsurf(i),0.0)            q_zref(i) = delq(i) + max(qsurf(i), 0.0)
211            te_zref(i) = delte(i) + ts1(i)            te_zref(i) = delte(i) + ts1(i)
212            temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)            temp(i) = te_zref(i) * (psol(i)/pref(i))**(-RKAPPA)
213  !         temp(i) = te_zref(i) - (zref* RG)/RCPD/ &         ENDDO
214  !                   (1 + RVTMP2 * max(q_zref(i),0.0))      ENDDO
215          ENDDO  
216  !      DO i = 1, knon
217        ENDDO         u_zref_c(i) = u_zref(i)
218  !  
219        DO i = 1, knon         u_10m(i) = u_zref_p(i) * ok_pred(i) + u_zref_c(i) * ok_corr(i)
220          u_zref_c(i) = u_zref(i)  
221  !         q_zref_c(i) = q_zref(i)
222          u_10m(i) = u_zref_p(i) * ok_pred(i) + u_zref_c(i) * ok_corr(i)         temp_c(i) = temp(i)
223  !         t_10m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)
224  !AM         q_10m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)
225          q_zref_c(i) = q_zref(i)      ENDDO
226          temp_c(i) = temp(i)  
227          t_10m(i) = temp_p(i) * ok_pred(i) + temp_c(i) * ok_corr(i)    END subroutine stdlevvar
228          q_10m(i) = q_zref_p(i) * ok_pred(i) + q_zref_c(i) * ok_corr(i)  
229  !MA  end module stdlevvar_m
       ENDDO  
 !  
       RETURN  
       END subroutine stdlevvar  
Legend:
Removed from v.38  
changed lines
  Added in v.108
  ViewVC Help
Powered by ViewVC 1.1.21