/[lmdze]/trunk/libf/phylmd/Conflx/flxflux.f90
ViewVC logotype

Contents of /trunk/libf/phylmd/Conflx/flxflux.f90

Parent Directory Parent Directory | Revision Log Revision Log


Revision 71 - (show annotations)
Mon Jul 8 18:12:18 2013 UTC (10 years, 9 months ago) by guez
File size: 8821 byte(s)
No reason to call inidissip in ce0l.

In inidissip, set random seed to 1 beacuse PGI compiler does not
accept all zeros.

dq was computed needlessly in caladvtrac. Arguments masse and dq of
calfis not used.

Replaced real*8 by double precision.

Pass arrays with inverted order of vertical levels to conflx instead
of creating local variables for this inside conflx.

1 module flxflux_m
2
3 IMPLICIT none
4
5 contains
6
7 SUBROUTINE flxflux(pdtime, pqen, pqsen, ptenh, pqenh, pap, paph, ldland, &
8 pgeoh, kcbot, kctop, lddraf, kdtop, ktype, ldcum, mfu, pmfd, mfus, &
9 pmfds, mfuq, pmfdq, mful, plude, pdmfup, pdmfdp, pten, prfl, psfl, &
10 pdpmel, ktopm2, pmflxr, pmflxs)
11
12 ! This routine does the final calculation of convective fluxes in
13 ! the cloud layer and in the subcloud layer.
14
15 USE dimphy, ONLY: klev, klon
16 USE suphec_m, ONLY: rcpd, retv, rg, rlmlt, rtt
17 USE yoethf_m, ONLY: r2es
18 USE fcttre, ONLY: foeew
19
20 REAL, intent(in):: pdtime
21 REAL, intent(in):: pqen(klon, klev)
22 real, intent(inout):: pqsen(klon, klev)
23 REAL, intent(in):: ptenh(klon, klev)
24 REAL, intent(in):: pqenh(klon, klev)
25 REAL pap(klon, klev)
26 REAL paph(klon, klev+1)
27 LOGICAL ldland(klon)
28 real pgeoh(klon, klev)
29 INTEGER, intent(in):: kcbot(klon), kctop(klon)
30 LOGICAL lddraf(klon)
31 INTEGER kdtop(klon)
32 INTEGER, intent(inout):: ktype(klon)
33 LOGICAL, intent(in):: ldcum(klon)
34 REAL mfu(klon, klev)
35 REAL pmfd(klon, klev)
36 real, intent(inout):: mfus(klon, klev)
37 real pmfds(klon, klev)
38 REAL mfuq(klon, klev)
39 REAL pmfdq(klon, klev)
40 real mful(klon, klev)
41 REAL plude(klon, klev)
42 REAL pdmfup(klon, klev)
43 REAL pdmfdp(klon, klev)
44 REAL, intent(in):: pten(klon, klev)
45 REAL prfl(klon), psfl(klon)
46 real pdpmel(klon, klev)
47 INTEGER, intent(out):: ktopm2
48 REAL pmflxr(klon, klev+1), pmflxs(klon, klev+1)
49
50 ! Local:
51 REAL cevapcu(klev)
52 REAL ztmsmlt, zdelta, zqsat
53
54 !jq The variable maxpdmfdp(klon) has been introduced by Olivier Boucher
55 !jq 14/11/00 to fix the problem with the negative precipitation.
56 real maxpdmfdp(klon, klev)
57
58 INTEGER k, kp, i
59 REAL zcons1, zcons2, zcucov, ztmelp2
60 real zdp, zzp, zfac, zsnmlt, zrfl, zrnew
61 REAL zrmin, zrfln, zdrfl
62 REAL zpds, zpdr, zdenom
63 INTEGER ikb
64
65 !---------------------------------------------------------------
66
67 DO k=1, klev
68 CEVAPCU(k)=1.93E-6*261.*SQRT(1.E3/(38.3*0.293) &
69 *SQRT(0.5*(paph(1, k)+paph(1, k+1))/paph(1, klev+1))) * 0.5/RG
70 end DO
71
72 ! SPECIFY CONSTANTS
73
74 zcons1 = RCPD/(RLMLT*RG*pdtime)
75 zcons2 = 1./(RG*pdtime)
76 zcucov = 0.05
77 ztmelp2 = RTT + 2.
78
79 ! DETERMINE FINAL CONVECTIVE FLUXES
80
81 DO i = 1, klon
82 IF (.NOT. ldcum(i) .OR. kdtop(i) < kctop(i)) lddraf(i)=.FALSE.
83 IF (.NOT. ldcum(i)) ktype(i) = 0
84 end DO
85
86 ktopm2 = min(klev, minval(kctop)) - 2
87
88 DO k = ktopm2, klev
89 DO i = 1, klon
90 IF (ldcum(i) .AND. k >= kctop(i) - 1) THEN
91 mfus(i, k) = mfus(i, k) &
92 - mfu(i, k) * (RCPD * ptenh(i, k) + pgeoh(i, k))
93 mfuq(i, k)=mfuq(i, k)-mfu(i, k)*pqenh(i, k)
94 zdp = 1.5E4
95 IF (ldland(i)) zdp = 3.E4
96
97 ! L'eau liquide détrainée est precipitée quand certaines
98 ! conditions sont réunies (sinon, elle est considérée
99 ! évaporée dans l'environnement)
100
101 IF (paph(i, kcbot(i)) - paph(i, kctop(i)) >= zdp &
102 .AND. pqen(i, k - 1) > 0.8 * pqsen(i, k - 1)) &
103 pdmfup(i, k - 1) = pdmfup(i, k - 1) + plude(i, k - 1)
104
105 IF (lddraf(i).AND.k >= kdtop(i)) THEN
106 pmfds(i, k)=pmfds(i, k)-pmfd(i, k)* &
107 (RCPD*ptenh(i, k)+pgeoh(i, k))
108 pmfdq(i, k)=pmfdq(i, k)-pmfd(i, k)*pqenh(i, k)
109 ELSE
110 pmfd(i, k)=0.
111 pmfds(i, k)=0.
112 pmfdq(i, k)=0.
113 pdmfdp(i, k-1)=0.
114 END IF
115 ELSE
116 mfu(i, k)=0.
117 mfus(i, k)=0.
118 mfuq(i, k)=0.
119 mful(i, k)=0.
120 pdmfup(i, k-1)=0.
121 plude(i, k-1)=0.
122 pmfd(i, k)=0.
123 pmfds(i, k)=0.
124 pmfdq(i, k)=0.
125 pdmfdp(i, k-1)=0.
126 ENDIF
127 end DO
128 end DO
129
130 DO k=ktopm2, klev
131 DO i = 1, klon
132 IF (ldcum(i) .AND. k > kcbot(i)) THEN
133 ikb=kcbot(i)
134 zzp = ((paph(i, klev + 1) - paph(i, k)) &
135 / (paph(i, klev + 1) - paph(i, ikb)))
136 IF (ktype(i) == 3) zzp = zzp**2
137 mfu(i, k)=mfu(i, ikb)*zzp
138 mfus(i, k)=mfus(i, ikb)*zzp
139 mfuq(i, k)=mfuq(i, ikb)*zzp
140 mful(i, k)=mful(i, ikb)*zzp
141 ENDIF
142 end DO
143 end DO
144
145 ! CALCULATE RAIN/SNOW FALL RATES
146 ! CALCULATE MELTING OF SNOW
147 ! CALCULATE EVAPORATION OF PRECIP
148
149 DO k = 1, klev+1
150 DO i = 1, klon
151 pmflxr(i, k) = 0.0
152 pmflxs(i, k) = 0.0
153 ENDDO
154 ENDDO
155 DO k = ktopm2, klev
156 DO i = 1, klon
157 IF (ldcum(i)) THEN
158 IF (pmflxs(i, k) > 0.0 .AND. pten(i, k) > ztmelp2) THEN
159 zfac=zcons1*(paph(i, k+1)-paph(i, k))
160 zsnmlt=MIN(pmflxs(i, k), zfac*(pten(i, k)-ztmelp2))
161 pdpmel(i, k)=zsnmlt
162 ztmsmlt=pten(i, k)-zsnmlt/zfac
163 zdelta=MAX(0., SIGN(1., RTT-ztmsmlt))
164 zqsat = R2ES * FOEEW(ztmsmlt, zdelta) / pap(i, k)
165 zqsat = MIN(0.5, zqsat)
166 zqsat = zqsat / (1. - RETV * zqsat)
167 pqsen(i, k) = zqsat
168 ENDIF
169 IF (pten(i, k) > RTT) THEN
170 pmflxr(i, k + 1) = pmflxr(i, k) + pdmfup(i, k) + pdmfdp(i, k) &
171 + pdpmel(i, k)
172 pmflxs(i, k+1)=pmflxs(i, k)-pdpmel(i, k)
173 ELSE
174 pmflxs(i, k+1)=pmflxs(i, k)+pdmfup(i, k)+pdmfdp(i, k)
175 pmflxr(i, k+1)=pmflxr(i, k)
176 ENDIF
177 ! si la precipitation est negative, on ajuste le plux du
178 ! panache descendant pour eliminer la negativite
179 IF ((pmflxr(i, k+1)+pmflxs(i, k+1)) < 0.0) THEN
180 pdmfdp(i, k) = -pmflxr(i, k)-pmflxs(i, k)-pdmfup(i, k)
181 pmflxr(i, k+1) = 0.0
182 pmflxs(i, k+1) = 0.0
183 pdpmel(i, k) = 0.0
184 ENDIF
185 ENDIF
186 ENDDO
187 ENDDO
188
189 ! The new variable is initialized here. It contains the
190 ! humidity which is fed to the downdraft by evaporation of
191 ! precipitation in the column below the base of convection.
192
193 ! In the former version, this term has been subtracted from precip
194 ! as well as the evaporation.
195
196 DO k = 1, klev
197 DO i = 1, klon
198 maxpdmfdp(i, k)=0.0
199 ENDDO
200 ENDDO
201 DO k = 1, klev
202 DO kp = k, klev
203 DO i = 1, klon
204 maxpdmfdp(i, k)=maxpdmfdp(i, k)+pdmfdp(i, kp)
205 ENDDO
206 ENDDO
207 ENDDO
208 ! End of initialization
209
210 DO k = ktopm2, klev
211 DO i = 1, klon
212 IF (ldcum(i) .AND. k >= kcbot(i)) THEN
213 zrfl = pmflxr(i, k) + pmflxs(i, k)
214 IF (zrfl > 1.0E-20) THEN
215 zrnew = (MAX(0., SQRT(zrfl / zcucov) - CEVAPCU(k) &
216 * (paph(i, k + 1) - paph(i, k)) &
217 * MAX(0., pqsen(i, k) - pqen(i, k))))**2 * zcucov
218 zrmin = zrfl - zcucov &
219 * MAX(0., 0.8 * pqsen(i, k) - pqen(i, k)) &
220 * zcons2 * (paph(i, k + 1) - paph(i, k))
221 zrnew=MAX(zrnew, zrmin)
222 zrfln=MAX(zrnew, 0.)
223 zdrfl=MIN(0., zrfln-zrfl)
224 ! At least the amount of precipiation needed to feed
225 ! the downdraft with humidity below the base of
226 ! convection has to be left and can't be evaporated
227 ! (surely the evaporation can't be positive):
228 zdrfl=MAX(zdrfl, &
229 MIN(-pmflxr(i, k)-pmflxs(i, k)-maxpdmfdp(i, k), 0.0))
230
231 zdenom=1.0/MAX(1.0E-20, pmflxr(i, k)+pmflxs(i, k))
232 IF (pten(i, k) > RTT) THEN
233 zpdr = pdmfdp(i, k)
234 zpds = 0.0
235 ELSE
236 zpdr = 0.0
237 zpds = pdmfdp(i, k)
238 ENDIF
239 pmflxr(i, k+1) = pmflxr(i, k) + zpdr + pdpmel(i, k) &
240 + zdrfl*pmflxr(i, k)*zdenom
241 pmflxs(i, k+1) = pmflxs(i, k) + zpds - pdpmel(i, k) &
242 + zdrfl*pmflxs(i, k)*zdenom
243 pdmfup(i, k) = pdmfup(i, k) + zdrfl
244 ELSE
245 pmflxr(i, k+1) = 0.0
246 pmflxs(i, k+1) = 0.0
247 pdmfdp(i, k) = 0.0
248 pdpmel(i, k) = 0.0
249 ENDIF
250 if (pmflxr(i, k) + pmflxs(i, k) < -1.e-26) &
251 write(*, *) 'precip. < 1e-16 ', pmflxr(i, k) + pmflxs(i, k)
252 ENDIF
253 ENDDO
254 ENDDO
255
256 DO i = 1, klon
257 prfl(i) = pmflxr(i, klev+1)
258 psfl(i) = pmflxs(i, klev+1)
259 end DO
260
261 RETURN
262 END SUBROUTINE flxflux
263
264 end module flxflux_m

  ViewVC Help
Powered by ViewVC 1.1.21