| 1 | ! $Header$ |
|---|
| 2 | |
|---|
| 3 | SUBROUTINE convect3(dtime, epmax, ok_adj, t1, r1, rs, u, v, tra, p, ph, nd, & |
|---|
| 4 | ndp1, nl, ntra, delt, iflag, ft, fr, fu, fv, ftra, precip, icb, inb, & |
|---|
| 5 | upwd, dnwd, dnwd0, sig, w0, mike, mke, ma, ments, qents, tps, tls, sigij, & |
|---|
| 6 | cape, tvp, pbase, buoybase, & ! ccc * |
|---|
| 7 | ! DTVPDT1,DTVPDQ1,DPLCLDT,DPLCLDR) |
|---|
| 8 | dtvpdt1, dtvpdq1, dplcldt, dplcldr, & ! sbl |
|---|
| 9 | ft2, fr2, fu2, fv2, wd, qcond, qcondc) ! sbl |
|---|
| 10 | |
|---|
| 11 | ! *** THE PARAMETER NA SHOULD IN GENERAL EQUAL ND *** |
|---|
| 12 | |
|---|
| 13 | ! ################################################################# |
|---|
| 14 | ! Fleur Introduction des traceurs dans convect3 le 6 juin 200 |
|---|
| 15 | ! ################################################################# |
|---|
| 16 | USE dimphy |
|---|
| 17 | USE infotrac_phy, ONLY: nbtr |
|---|
| 18 | USE lmdz_yomcst |
|---|
| 19 | |
|---|
| 20 | IMPLICIT NONE |
|---|
| 21 | INTEGER na |
|---|
| 22 | PARAMETER (na = 60) |
|---|
| 23 | |
|---|
| 24 | REAL deltac ! cld |
|---|
| 25 | PARAMETER (deltac = 0.01) ! cld |
|---|
| 26 | |
|---|
| 27 | INTEGER nent(na) |
|---|
| 28 | INTEGER nd, ndp1, nl, ntra, iflag, icb, inb |
|---|
| 29 | REAL dtime, epmax, delt, precip, cape |
|---|
| 30 | REAL dplcldt, dplcldr |
|---|
| 31 | REAL t1(nd), r1(nd), rs(nd), u(nd), v(nd), tra(nd, ntra) |
|---|
| 32 | REAL p(nd), ph(ndp1) |
|---|
| 33 | REAL ft(nd), fr(nd), fu(nd), fv(nd), ftra(nd, ntra) |
|---|
| 34 | REAL sig(nd), w0(nd) |
|---|
| 35 | REAL uent(na, na), vent(na, na), traent(na, na, nbtr), tratm(na) |
|---|
| 36 | REAL up(na), vp(na), trap(na, nbtr) |
|---|
| 37 | REAL m(na), mp(na), ment(na, na), qent(na, na), elij(na, na) |
|---|
| 38 | REAL sij(na, na), tvp(na), tv(na), water(na) |
|---|
| 39 | REAL rp(na), ep(na), th(na), wt(na), evap(na), clw(na) |
|---|
| 40 | REAL sigp(na), b(na), tp(na), cpn(na) |
|---|
| 41 | REAL lv(na), lvcp(na), h(na), hp(na), gz(na) |
|---|
| 42 | REAL t(na), rr(na) |
|---|
| 43 | |
|---|
| 44 | REAL ft2(nd), fr2(nd), fu2(nd), fv2(nd) ! added sbl |
|---|
| 45 | REAL u1(nd), v1(nd) ! added sbl |
|---|
| 46 | |
|---|
| 47 | REAL buoy(na) ! Lifted parcel buoyancy |
|---|
| 48 | REAL dtvpdt1(nd), dtvpdq1(nd) ! Derivatives of parcel virtual |
|---|
| 49 | ! temperature wrt T1 and Q1 |
|---|
| 50 | REAL clw_new(na), qi(na) |
|---|
| 51 | |
|---|
| 52 | REAL wd, betad ! for gust factor (sb) |
|---|
| 53 | REAL qcondc(nd) ! interface cld param (sb) |
|---|
| 54 | REAL qcond(nd), nqcond(na), wa(na), maa(na), siga(na), axc(na) ! cld |
|---|
| 55 | |
|---|
| 56 | LOGICAL ice_conv, ok_adj |
|---|
| 57 | PARAMETER (ice_conv = .TRUE.) |
|---|
| 58 | |
|---|
| 59 | ! ccccccccccccccccccccccccccccccccccccccccccccc |
|---|
| 60 | ! declaration des variables a sortir |
|---|
| 61 | ! cccccccccccccccccccccccccccccccccccccccccccc |
|---|
| 62 | REAL mke(nd) |
|---|
| 63 | REAL mike(nd) |
|---|
| 64 | REAL ma(nd) |
|---|
| 65 | REAL tps(nd) !temperature dans les ascendances non diluees |
|---|
| 66 | REAL tls(nd) !temperature potentielle |
|---|
| 67 | REAL ments(nd, nd) |
|---|
| 68 | REAL qents(nd, nd) |
|---|
| 69 | REAL sigij(klev, klev) |
|---|
| 70 | REAL pbase ! pressure at the cloud base level |
|---|
| 71 | REAL buoybase ! buoyancy at the cloud base level |
|---|
| 72 | ! ccccccccccccccccccccccccccccccccccccccccccccc |
|---|
| 73 | |
|---|
| 74 | REAL :: cpv, cl, cpvmcl, eps, alv0, rdcp, pbcrit, ptcrit, sigd, spfac |
|---|
| 75 | REAL :: tau, beta, alpha, dtcrit, dtovsh, ahm, rm, um, vm, dphinv |
|---|
| 76 | REAL :: a2, x, tvx, tvy, plcl, pden, dpbase, tvpbase, tvbase, tdif |
|---|
| 77 | REAL :: ath1, ath, delti, deltap, dcape, dlnp, sigold, dtmin, fac, w |
|---|
| 78 | REAL :: amu, rti, cpd, bf2, anum, denom, dei, altem, cwat, stemp, qp |
|---|
| 79 | REAL :: scrit, alt, smax, asij, wgh, sjmax, sjmin, smid, delp, delm |
|---|
| 80 | REAL :: asum, bsum, csum, wflux, tinv, wdtrain, awat, afac, afac1, afac2 |
|---|
| 81 | REAL :: bfac, pr1, pr2, sigt, b6, c6, revap, tevap, delth, amfac, amp2 |
|---|
| 82 | REAL :: xf, tf, af, bf, fac2, ur, sru, d, ampmax, dpinv, am, amde, cpinv |
|---|
| 83 | REAL :: amp1, ad, rat, ax, bx, cx, dx, ex, dsum |
|---|
| 84 | INTEGER :: nk, i, j, nopt, jn, k, im, jm, n |
|---|
| 85 | |
|---|
| 86 | REAL dnwd0(nd) ! precipitation driven unsaturated downdraft flux |
|---|
| 87 | REAL dnwd(nd), dn1 ! in-cloud saturated downdraft mass flux |
|---|
| 88 | REAL upwd(nd), up1 ! in-cloud saturated updraft mass flux |
|---|
| 89 | |
|---|
| 90 | ! *** ASSIGN VALUES OF THERMODYNAMIC CONSTANTS *** |
|---|
| 91 | ! *** THESE SHOULD BE CONSISTENT WITH *** |
|---|
| 92 | ! *** THOSE USED IN CALLING PROGRAM *** |
|---|
| 93 | ! *** NOTE: THESE ARE ALSO SPECIFIED IN SUBROUTINE TLIFT *** |
|---|
| 94 | |
|---|
| 95 | ! sb CPD=1005.7 |
|---|
| 96 | ! sb CPV=1870.0 |
|---|
| 97 | ! sb CL=4190.0 |
|---|
| 98 | ! sb CPVMCL=CL-CPV |
|---|
| 99 | ! sb RV=461.5 |
|---|
| 100 | ! sb RD=287.04 |
|---|
| 101 | ! sb EPS=RD/RV |
|---|
| 102 | ! sb ALV0=2.501E6 |
|---|
| 103 | ! cccccccccccccccccccccc |
|---|
| 104 | ! constantes coherentes avec le modele du Centre Europeen |
|---|
| 105 | ! sb RD = 1000.0 * 1.380658E-23 * 6.0221367E+23 / 28.9644 |
|---|
| 106 | ! sb RV = 1000.0 * 1.380658E-23 * 6.0221367E+23 / 18.0153 |
|---|
| 107 | ! sb CPD = 3.5 * RD |
|---|
| 108 | ! sb CPV = 4.0 * RV |
|---|
| 109 | ! sb CL = 4218.0 |
|---|
| 110 | ! sb CPVMCL=CL-CPV |
|---|
| 111 | ! sb EPS=RD/RV |
|---|
| 112 | ! sb ALV0=2.5008E+06 |
|---|
| 113 | ! ccccccccccccccccccccc |
|---|
| 114 | ! on utilise les constantes thermo du Centre Europeen: (SB) |
|---|
| 115 | |
|---|
| 116 | cpd = rcpd |
|---|
| 117 | cpv = rcpv |
|---|
| 118 | cl = rcw |
|---|
| 119 | cpvmcl = cl - cpv |
|---|
| 120 | eps = rd / rv |
|---|
| 121 | alv0 = rlvtt |
|---|
| 122 | |
|---|
| 123 | nk = 1 ! origin level of the lifted parcel |
|---|
| 124 | |
|---|
| 125 | ! ccccccccccccccccccccc |
|---|
| 126 | |
|---|
| 127 | ! *** INITIALIZE OUTPUT ARRAYS AND PARAMETERS *** |
|---|
| 128 | |
|---|
| 129 | DO i = 1, nd |
|---|
| 130 | ft(i) = 0.0 |
|---|
| 131 | fr(i) = 0.0 |
|---|
| 132 | fu(i) = 0.0 |
|---|
| 133 | fv(i) = 0.0 |
|---|
| 134 | |
|---|
| 135 | ft2(i) = 0.0 |
|---|
| 136 | fr2(i) = 0.0 |
|---|
| 137 | fu2(i) = 0.0 |
|---|
| 138 | fv2(i) = 0.0 |
|---|
| 139 | |
|---|
| 140 | DO j = 1, ntra |
|---|
| 141 | ftra(i, j) = 0.0 |
|---|
| 142 | END DO |
|---|
| 143 | |
|---|
| 144 | qcondc(i) = 0.0 ! cld |
|---|
| 145 | qcond(i) = 0.0 ! cld |
|---|
| 146 | nqcond(i) = 0.0 ! cld |
|---|
| 147 | |
|---|
| 148 | t(i) = t1(i) |
|---|
| 149 | rr(i) = r1(i) |
|---|
| 150 | u1(i) = u(i) ! added sbl |
|---|
| 151 | v1(i) = v(i) ! added sbl |
|---|
| 152 | END DO |
|---|
| 153 | DO i = 1, nl |
|---|
| 154 | rdcp = (rd * (1. - rr(i)) + rr(i) * rv) / (cpd * (1. - rr(i)) + rr(i) * cpv) |
|---|
| 155 | th(i) = t(i) * (1000.0 / p(i))**rdcp |
|---|
| 156 | END DO |
|---|
| 157 | |
|---|
| 158 | ! ************************************************************ |
|---|
| 159 | ! * CALCUL DES TEMPERATURES POTENTIELLES A SORTIR |
|---|
| 160 | ! ************************************************************ |
|---|
| 161 | DO i = 1, nd |
|---|
| 162 | rdcp = (rd * (1. - rr(i)) + rr(i) * rv) / (cpd * (1. - rr(i)) + rr(i) * cpv) |
|---|
| 163 | |
|---|
| 164 | tls(i) = t(i) * (1000.0 / p(i))**rdcp |
|---|
| 165 | END DO |
|---|
| 166 | |
|---|
| 167 | |
|---|
| 168 | |
|---|
| 169 | |
|---|
| 170 | ! *********************************************************** |
|---|
| 171 | |
|---|
| 172 | precip = 0.0 |
|---|
| 173 | wd = 0.0 ! sb |
|---|
| 174 | iflag = 1 |
|---|
| 175 | |
|---|
| 176 | ! *** SPECIFY PARAMETERS *** |
|---|
| 177 | ! *** PBCRIT IS THE CRITICAL CLOUD DEPTH (MB) BENEATH WHICH THE *** |
|---|
| 178 | ! *** PRECIPITATION EFFICIENCY IS ASSUMED TO BE ZERO *** |
|---|
| 179 | ! *** PTCRIT IS THE CLOUD DEPTH (MB) ABOVE WHICH THE PRECIP. *** |
|---|
| 180 | ! *** EFFICIENCY IS ASSUMED TO BE UNITY *** |
|---|
| 181 | ! *** SIGD IS THE FRACTIONAL AREA COVERED BY UNSATURATED DNDRAFT *** |
|---|
| 182 | ! *** SPFAC IS THE FRACTION OF PRECIPITATION FALLING OUTSIDE *** |
|---|
| 183 | ! *** OF CLOUD *** |
|---|
| 184 | ! *** ALPHA AND BETA ARE PARAMETERS THAT CONTROL THE RATE OF *** |
|---|
| 185 | ! *** APPROACH TO QUASI-EQUILIBRIUM *** |
|---|
| 186 | ! *** (THEIR STANDARD VALUES ARE 1.0 AND 0.96, RESPECTIVELY) *** |
|---|
| 187 | ! *** (BETA MUST BE LESS THAN OR EQUAL TO 1) *** |
|---|
| 188 | ! *** DTCRIT IS THE CRITICAL BUOYANCY (K) USED TO ADJUST THE *** |
|---|
| 189 | ! *** APPROACH TO QUASI-EQUILIBRIUM *** |
|---|
| 190 | ! *** IT MUST BE LESS THAN 0 *** |
|---|
| 191 | |
|---|
| 192 | pbcrit = 150.0 |
|---|
| 193 | ptcrit = 500.0 |
|---|
| 194 | sigd = 0.01 |
|---|
| 195 | spfac = 0.15 |
|---|
| 196 | ! sb: |
|---|
| 197 | ! EPMAX=0.993 ! precip efficiency less than unity |
|---|
| 198 | ! EPMAX=1. ! precip efficiency less than unity |
|---|
| 199 | |
|---|
| 200 | ! jyg |
|---|
| 201 | ! CC BETA=0.96 |
|---|
| 202 | ! Beta is now expressed as a function of the characteristic time |
|---|
| 203 | ! of the convective process. |
|---|
| 204 | ! CC Old value : TAU = 15000. !(for dtime = 600.s) |
|---|
| 205 | ! CC Other value (inducing little change) :TAU = 8000. |
|---|
| 206 | tau = 8000. |
|---|
| 207 | beta = 1. - dtime / tau |
|---|
| 208 | ! jyg |
|---|
| 209 | ! CC ALPHA=1.0 |
|---|
| 210 | alpha = 1.5E-3 * dtime / tau |
|---|
| 211 | ! Increase alpha in order to compensate W decrease |
|---|
| 212 | alpha = alpha * 1.5 |
|---|
| 213 | |
|---|
| 214 | ! jyg (voir CONVECT 3) |
|---|
| 215 | ! CC DTCRIT=-0.2 |
|---|
| 216 | dtcrit = -2. |
|---|
| 217 | ! gf&jyg |
|---|
| 218 | ! CC DT pour l'overshoot. |
|---|
| 219 | dtovsh = -0.2 |
|---|
| 220 | |
|---|
| 221 | |
|---|
| 222 | ! *** INCREMENT THE COUNTER *** |
|---|
| 223 | |
|---|
| 224 | sig(nd) = sig(nd) + 1.0 |
|---|
| 225 | sig(nd) = amin1(sig(nd), 12.1) |
|---|
| 226 | |
|---|
| 227 | ! *** IF NOPT IS AN INTEGER OTHER THAN 0, CONVECT *** |
|---|
| 228 | ! *** RETURNS ARRAYS T AND R THAT MAY HAVE BEEN *** |
|---|
| 229 | ! *** ALTERED BY DRY ADIABATIC ADJUSTMENT; OTHERWISE *** |
|---|
| 230 | ! *** THE RETURNED ARRAYS ARE UNALTERED. *** |
|---|
| 231 | |
|---|
| 232 | nopt = 0 |
|---|
| 233 | ! NOPT=1 ! sbl |
|---|
| 234 | |
|---|
| 235 | ! *** PERFORM DRY ADIABATIC ADJUSTMENT *** |
|---|
| 236 | |
|---|
| 237 | ! *** DO NOT BYPASS THIS EVEN IF THE CALLING PROGRAM HAS A *** |
|---|
| 238 | ! *** BOUNDARY LAYER SCHEME !!! *** |
|---|
| 239 | |
|---|
| 240 | IF (ok_adj) THEN ! added sbl |
|---|
| 241 | |
|---|
| 242 | DO i = nl - 1, 1, -1 |
|---|
| 243 | jn = 0 |
|---|
| 244 | DO j = i + 1, nl |
|---|
| 245 | IF (th(j)<th(i)) jn = j |
|---|
| 246 | END DO |
|---|
| 247 | IF (jn==0) GO TO 30 |
|---|
| 248 | ahm = 0.0 |
|---|
| 249 | rm = 0.0 |
|---|
| 250 | um = 0.0 |
|---|
| 251 | vm = 0.0 |
|---|
| 252 | DO k = 1, ntra |
|---|
| 253 | tratm(k) = 0.0 |
|---|
| 254 | END DO |
|---|
| 255 | DO j = i, jn |
|---|
| 256 | ahm = ahm + (cpd * (1. - rr(j)) + rr(j) * cpv) * t(j) * (ph(j) - ph(j + 1)) |
|---|
| 257 | rm = rm + rr(j) * (ph(j) - ph(j + 1)) |
|---|
| 258 | um = um + u(j) * (ph(j) - ph(j + 1)) |
|---|
| 259 | vm = vm + v(j) * (ph(j) - ph(j + 1)) |
|---|
| 260 | DO k = 1, ntra |
|---|
| 261 | tratm(k) = tratm(k) + tra(j, k) * (ph(j) - ph(j + 1)) |
|---|
| 262 | END DO |
|---|
| 263 | END DO |
|---|
| 264 | dphinv = 1. / (ph(i) - ph(jn + 1)) |
|---|
| 265 | rm = rm * dphinv |
|---|
| 266 | um = um * dphinv |
|---|
| 267 | vm = vm * dphinv |
|---|
| 268 | DO k = 1, ntra |
|---|
| 269 | tratm(k) = tratm(k) * dphinv |
|---|
| 270 | END DO |
|---|
| 271 | a2 = 0.0 |
|---|
| 272 | DO j = i, jn |
|---|
| 273 | rr(j) = rm |
|---|
| 274 | u(j) = um |
|---|
| 275 | v(j) = vm |
|---|
| 276 | DO k = 1, ntra |
|---|
| 277 | tra(j, k) = tratm(k) |
|---|
| 278 | END DO |
|---|
| 279 | rdcp = (rd * (1. - rr(j)) + rr(j) * rv) / (cpd * (1. - rr(j)) + rr(j) * cpv) |
|---|
| 280 | x = (0.001 * p(j))**rdcp |
|---|
| 281 | t(j) = x |
|---|
| 282 | a2 = a2 + (cpd * (1. - rr(j)) + rr(j) * cpv) * x * (ph(j) - ph(j + 1)) |
|---|
| 283 | END DO |
|---|
| 284 | DO j = i, jn |
|---|
| 285 | th(j) = ahm / a2 |
|---|
| 286 | t(j) = t(j) * th(j) |
|---|
| 287 | END DO |
|---|
| 288 | 30 END DO |
|---|
| 289 | |
|---|
| 290 | END IF ! added sbl |
|---|
| 291 | |
|---|
| 292 | ! *** RESET INPUT ARRAYS IF ok_adj 0 *** |
|---|
| 293 | |
|---|
| 294 | IF (ok_adj) THEN |
|---|
| 295 | DO i = 1, nd |
|---|
| 296 | |
|---|
| 297 | ft2(i) = (t(i) - t1(i)) / delt ! sbl |
|---|
| 298 | fr2(i) = (rr(i) - r1(i)) / delt ! sbl |
|---|
| 299 | fu2(i) = (u(i) - u1(i)) / delt ! sbl |
|---|
| 300 | fv2(i) = (v(i) - v1(i)) / delt ! sbl |
|---|
| 301 | |
|---|
| 302 | ! T1(I)=T(I) ! commente sbl |
|---|
| 303 | ! R1(I)=RR(I) ! commente sbl |
|---|
| 304 | END DO |
|---|
| 305 | END IF |
|---|
| 306 | |
|---|
| 307 | ! *** CALCULATE ARRAYS OF GEOPOTENTIAL, HEAT CAPACITY AND STATIC ENERGY |
|---|
| 308 | |
|---|
| 309 | gz(1) = 0.0 |
|---|
| 310 | cpn(1) = cpd * (1. - rr(1)) + rr(1) * cpv |
|---|
| 311 | h(1) = t(1) * cpn(1) |
|---|
| 312 | DO i = 2, nl |
|---|
| 313 | tvx = t(i) * (1. + rr(i) / eps - rr(i)) |
|---|
| 314 | tvy = t(i - 1) * (1. + rr(i - 1) / eps - rr(i - 1)) |
|---|
| 315 | gz(i) = gz(i - 1) + 0.5 * rd * (tvx + tvy) * (p(i - 1) - p(i)) / ph(i) |
|---|
| 316 | cpn(i) = cpd * (1. - rr(i)) + cpv * rr(i) |
|---|
| 317 | h(i) = t(i) * cpn(i) + gz(i) |
|---|
| 318 | END DO |
|---|
| 319 | |
|---|
| 320 | ! *** CALCULATE LIFTED CONDENSATION LEVEL OF AIR AT LOWEST MODEL LEVEL *** |
|---|
| 321 | ! *** (WITHIN 0.2% OF FORMULA OF BOLTON, MON. WEA. REV.,1980) *** |
|---|
| 322 | |
|---|
| 323 | IF (t(1)<250.0 .OR. rr(1)<=0.0) THEN |
|---|
| 324 | iflag = 0 |
|---|
| 325 | ! sb3d PRINT*,'je suis passe par 366' |
|---|
| 326 | RETURN |
|---|
| 327 | END IF |
|---|
| 328 | |
|---|
| 329 | ! jyg1 Utilisation de la SUBROUTINE CLIFT |
|---|
| 330 | ! C RH=RR(1)/RS(1) |
|---|
| 331 | ! C CHI=T(1)/(1669.0-122.0*RH-T(1)) |
|---|
| 332 | ! C PLCL=P(1)*(RH**CHI) |
|---|
| 333 | CALL clift(p(1), t(1), rr(1), rs(1), plcl, dplcldt, dplcldr) |
|---|
| 334 | ! jyg2 |
|---|
| 335 | ! sb3d PRINT *,' em_plcl,p1,t1,r1,rs1,rh ' |
|---|
| 336 | ! sb3d $ ,PLCL,P(1),T(1),RR(1),RS(1),RH |
|---|
| 337 | |
|---|
| 338 | IF (plcl<200.0 .OR. plcl>=2000.0) THEN |
|---|
| 339 | iflag = 2 |
|---|
| 340 | RETURN |
|---|
| 341 | END IF |
|---|
| 342 | ! jyg1 |
|---|
| 343 | ! Essais de modification de ICB |
|---|
| 344 | |
|---|
| 345 | ! *** CALCULATE FIRST LEVEL ABOVE LCL (=ICB) *** |
|---|
| 346 | |
|---|
| 347 | ! C ICB=NL-1 |
|---|
| 348 | ! C DO 50 I=2,NL-1 |
|---|
| 349 | ! C IF(P(I).LT.PLCL)THEN |
|---|
| 350 | ! C ICB=MIN(ICB,I) ! ICB sup ou egal a 2 |
|---|
| 351 | ! C END IF |
|---|
| 352 | ! C 50 CONTINUE |
|---|
| 353 | ! C IF(ICB.EQ.(NL-1))THEN |
|---|
| 354 | ! C IFLAG=3 |
|---|
| 355 | ! C RETURN |
|---|
| 356 | ! C END IF |
|---|
| 357 | |
|---|
| 358 | ! *** CALCULATE LAYER CONTAINING LCL (=ICB) *** |
|---|
| 359 | |
|---|
| 360 | icb = nl - 1 |
|---|
| 361 | ! sb DO 50 I=2,NL-1 |
|---|
| 362 | DO i = 3, nl - 1 ! modif sb pour que ICB soit sup/egal a 2 |
|---|
| 363 | ! la modification consiste a comparer PLCL a PH et non a P: |
|---|
| 364 | ! ICB est defini par : PH(ICB)<PLCL<PH(ICB-!) |
|---|
| 365 | IF (ph(i)<plcl) THEN |
|---|
| 366 | icb = min(icb, i) |
|---|
| 367 | END IF |
|---|
| 368 | END DO |
|---|
| 369 | IF (icb==(nl - 1)) THEN |
|---|
| 370 | iflag = 3 |
|---|
| 371 | RETURN |
|---|
| 372 | END IF |
|---|
| 373 | icb = icb - 1 ! ICB sup ou egal a 2 |
|---|
| 374 | ! jyg2 |
|---|
| 375 | |
|---|
| 376 | |
|---|
| 377 | |
|---|
| 378 | ! *** SUBROUTINE TLIFT CALCULATES PART OF THE LIFTED PARCEL VIRTUAL |
|---|
| 379 | ! *** |
|---|
| 380 | ! *** TEMPERATURE, THE ACTUAL TEMPERATURE AND THE ADIABATIC |
|---|
| 381 | ! *** |
|---|
| 382 | ! *** LIQUID WATER CONTENT |
|---|
| 383 | ! *** |
|---|
| 384 | |
|---|
| 385 | |
|---|
| 386 | ! jyg1 |
|---|
| 387 | ! make sure that "Cloud base" seen by TLIFT is actually the |
|---|
| 388 | ! fisrt level where adiabatic ascent is saturated |
|---|
| 389 | IF (plcl>p(icb)) THEN |
|---|
| 390 | ! sb CALL TLIFT(P,T,RR,RS,GZ,PLCL,ICB,TVP,TP,CLW,ND,NL) |
|---|
| 391 | CALL tlift(p, t, rr, rs, gz, plcl, icb, nk, tvp, tp, clw, nd, nl, & |
|---|
| 392 | dtvpdt1, dtvpdq1) |
|---|
| 393 | ELSE |
|---|
| 394 | ! sb CALL TLIFT(P,T,RR,RS,GZ,PLCL,ICB+1,TVP,TP,CLW,ND,NL) |
|---|
| 395 | CALL tlift(p, t, rr, rs, gz, plcl, icb + 1, nk, tvp, tp, clw, nd, nl, & |
|---|
| 396 | dtvpdt1, dtvpdq1) |
|---|
| 397 | END IF |
|---|
| 398 | ! jyg2 |
|---|
| 399 | |
|---|
| 400 | ! ***************************************************************************** |
|---|
| 401 | ! *** SORTIE DE LA TEMPERATURE DE L ASCENDANCE NON DILUE |
|---|
| 402 | ! ***************************************************************************** |
|---|
| 403 | DO i = 1, nd |
|---|
| 404 | tps(i) = tp(i) |
|---|
| 405 | END DO |
|---|
| 406 | |
|---|
| 407 | |
|---|
| 408 | ! ***************************************************************************** |
|---|
| 409 | |
|---|
| 410 | |
|---|
| 411 | ! *** SET THE PRECIPITATION EFFICIENCIES AND THE FRACTION OF *** |
|---|
| 412 | ! *** PRECIPITATION FALLING OUTSIDE OF CLOUD *** |
|---|
| 413 | ! *** THESE MAY BE FUNCTIONS OF TP(I), P(I) AND CLW(I) *** |
|---|
| 414 | |
|---|
| 415 | DO i = 1, nl |
|---|
| 416 | pden = ptcrit - pbcrit |
|---|
| 417 | |
|---|
| 418 | ! jyg |
|---|
| 419 | ! cc EP(I)=(P(ICB)-P(I)-PBCRIT)/PDEN |
|---|
| 420 | ! sb EP(I)=(PLCL-P(I)-PBCRIT)/PDEN |
|---|
| 421 | ep(i) = (plcl - p(i) - pbcrit) / pden * epmax ! sb |
|---|
| 422 | |
|---|
| 423 | ep(i) = amax1(ep(i), 0.0) |
|---|
| 424 | ! sb EP(I)=AMIN1(EP(I),1.0) |
|---|
| 425 | ep(i) = amin1(ep(i), epmax) ! sb |
|---|
| 426 | sigp(i) = spfac |
|---|
| 427 | |
|---|
| 428 | ! *** CALCULATE VIRTUAL TEMPERATURE AND LIFTED PARCEL *** |
|---|
| 429 | ! *** VIRTUAL TEMPERATURE *** |
|---|
| 430 | |
|---|
| 431 | tv(i) = t(i) * (1. + rr(i) / eps - rr(i)) |
|---|
| 432 | ! cd1 |
|---|
| 433 | ! . Keep all liquid water in lifted parcel (-> adiabatic CAPE) |
|---|
| 434 | |
|---|
| 435 | ! cc TVP(I)=TVP(I)-TP(I)*(RR(1)-EP(I)*CLW(I)) |
|---|
| 436 | ! !!! sb TVP(I)=TVP(I)-TP(I)*RR(1) ! calcule dans tlift |
|---|
| 437 | ! cd2 |
|---|
| 438 | |
|---|
| 439 | ! *** Calculate first estimate of buoyancy |
|---|
| 440 | |
|---|
| 441 | buoy(i) = tvp(i) - tv(i) |
|---|
| 442 | END DO |
|---|
| 443 | |
|---|
| 444 | ! *** Set Cloud Base Buoyancy at (Plcl+DPbase) level buoyancy |
|---|
| 445 | |
|---|
| 446 | dpbase = -40. !That is 400m above LCL |
|---|
| 447 | pbase = plcl + dpbase |
|---|
| 448 | tvpbase = tvp(icb) * (pbase - p(icb + 1)) / (p(icb) - p(icb + 1)) + & |
|---|
| 449 | tvp(icb + 1) * (p(icb) - pbase) / (p(icb) - p(icb + 1)) |
|---|
| 450 | tvbase = tv(icb) * (pbase - p(icb + 1)) / (p(icb) - p(icb + 1)) + & |
|---|
| 451 | tv(icb + 1) * (p(icb) - pbase) / (p(icb) - p(icb + 1)) |
|---|
| 452 | |
|---|
| 453 | ! test sb: |
|---|
| 454 | ! @ WRITE(*,*) '++++++++++++++++++++++++++++++++++++++++' |
|---|
| 455 | ! @ WRITE(*,*) 'plcl,dpbas,tvpbas,tvbas,tvp(icb),tvp(icb1) |
|---|
| 456 | ! @ : ,tv(icb),tv(icb1)' |
|---|
| 457 | ! @ WRITE(*,*) plcl,dpbase,tvpbase,tvbase,tvp(icb) |
|---|
| 458 | ! @ L ,tvp(icb+1),tv(icb),tv(icb+1) |
|---|
| 459 | ! @ WRITE(*,*) '++++++++++++++++++++++++++++++++++++++++' |
|---|
| 460 | ! fin test sb |
|---|
| 461 | buoybase = tvpbase - tvbase |
|---|
| 462 | |
|---|
| 463 | ! C Set buoyancy = BUOYBASE for all levels below BASE. |
|---|
| 464 | ! C For safety, set : BUOY(ICB) = BUOYBASE |
|---|
| 465 | DO i = icb, nl |
|---|
| 466 | IF (p(i)>=pbase) THEN |
|---|
| 467 | buoy(i) = buoybase |
|---|
| 468 | END IF |
|---|
| 469 | END DO |
|---|
| 470 | buoy(icb) = buoybase |
|---|
| 471 | |
|---|
| 472 | ! sb3d print *,'buoybase,tvp_tv,tvpbase,tvbase,pbase,plcl' |
|---|
| 473 | ! sb3d $, buoybase,tvp(icb)-tv(icb),tvpbase,tvbase,pbase,plcl |
|---|
| 474 | ! sb3d print *,'TVP ',(tvp(i),i=1,nl) |
|---|
| 475 | ! sb3d print *,'TV ',(tv(i),i=1,nl) |
|---|
| 476 | ! sb3d print *, 'P ',(p(i),i=1,nl) |
|---|
| 477 | ! sb3d print *,'ICB ',icb |
|---|
| 478 | ! test sb: |
|---|
| 479 | ! @ WRITE(*,*) '@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@' |
|---|
| 480 | ! @ WRITE(*,*) 'icb,icbs,inb,buoybase:' |
|---|
| 481 | ! @ WRITE(*,*) icb,icb+1,inb,buoybase |
|---|
| 482 | ! @ WRITE(*,*) 'k,tvp,tv,tp,buoy,ep: ' |
|---|
| 483 | ! @ do k=1,nl |
|---|
| 484 | ! @ WRITE(*,*) k,tvp(k),tv(k),tp(k),buoy(k),ep(k) |
|---|
| 485 | ! @ enddo |
|---|
| 486 | ! @ WRITE(*,*) '@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@' |
|---|
| 487 | ! fin test sb |
|---|
| 488 | |
|---|
| 489 | |
|---|
| 490 | |
|---|
| 491 | ! *** MAKE SURE THAT COLUMN IS DRY ADIABATIC BETWEEN THE SURFACE *** |
|---|
| 492 | ! *** AND CLOUD BASE, AND THAT LIFTED AIR IS POSITIVELY BUOYANT *** |
|---|
| 493 | ! *** AT CLOUD BASE *** |
|---|
| 494 | ! *** IF NOT, RETURN TO CALLING PROGRAM AFTER RESETTING *** |
|---|
| 495 | ! *** SIG(I) AND W0(I) *** |
|---|
| 496 | |
|---|
| 497 | ! jyg |
|---|
| 498 | ! CC TDIF=TVP(ICB)-TV(ICB) |
|---|
| 499 | tdif = buoy(icb) |
|---|
| 500 | ath1 = th(1) |
|---|
| 501 | ! jyg |
|---|
| 502 | ! CC ATH=TH(ICB-1)-1.0 |
|---|
| 503 | ath = th(icb - 1) - 5.0 |
|---|
| 504 | ! ATH=0. ! ajout sb |
|---|
| 505 | ! IF (ICB.GT.1) ATH=TH(ICB-1)-5.0 ! modif sb |
|---|
| 506 | IF (tdif<dtcrit .OR. ath>ath1) THEN |
|---|
| 507 | DO i = 1, nl |
|---|
| 508 | sig(i) = beta * sig(i) - 2. * alpha * tdif * tdif |
|---|
| 509 | sig(i) = amax1(sig(i), 0.0) |
|---|
| 510 | w0(i) = beta * w0(i) |
|---|
| 511 | END DO |
|---|
| 512 | iflag = 0 |
|---|
| 513 | RETURN |
|---|
| 514 | END IF |
|---|
| 515 | |
|---|
| 516 | |
|---|
| 517 | |
|---|
| 518 | ! *** IF THIS POINT IS REACHED, MOIST CONVECTIVE ADJUSTMENT IS NECESSARY |
|---|
| 519 | ! *** |
|---|
| 520 | ! *** NOW INITIALIZE VARIOUS ARRAYS USED IN THE COMPUTATIONS |
|---|
| 521 | ! *** |
|---|
| 522 | |
|---|
| 523 | DO i = 1, nl |
|---|
| 524 | hp(i) = h(i) |
|---|
| 525 | wt(i) = 0.001 |
|---|
| 526 | rp(i) = rr(i) |
|---|
| 527 | up(i) = u(i) |
|---|
| 528 | vp(i) = v(i) |
|---|
| 529 | DO j = 1, ntra |
|---|
| 530 | trap(i, j) = tra(i, j) |
|---|
| 531 | END DO |
|---|
| 532 | nent(i) = 0 |
|---|
| 533 | water(i) = 0.0 |
|---|
| 534 | evap(i) = 0.0 |
|---|
| 535 | b(i) = 0.0 |
|---|
| 536 | mp(i) = 0.0 |
|---|
| 537 | m(i) = 0.0 |
|---|
| 538 | lv(i) = alv0 - cpvmcl * (t(i) - 273.15) |
|---|
| 539 | lvcp(i) = lv(i) / cpn(i) |
|---|
| 540 | DO j = 1, nl |
|---|
| 541 | qent(i, j) = rr(j) |
|---|
| 542 | elij(i, j) = 0.0 |
|---|
| 543 | ment(i, j) = 0.0 |
|---|
| 544 | sij(i, j) = 0.0 |
|---|
| 545 | uent(i, j) = u(j) |
|---|
| 546 | vent(i, j) = v(j) |
|---|
| 547 | DO k = 1, ntra |
|---|
| 548 | traent(i, j, k) = tra(j, k) |
|---|
| 549 | END DO |
|---|
| 550 | END DO |
|---|
| 551 | END DO |
|---|
| 552 | |
|---|
| 553 | delti = 1.0 / delt |
|---|
| 554 | |
|---|
| 555 | ! *** FIND THE FIRST MODEL LEVEL (INB) ABOVE THE PARCEL'S *** |
|---|
| 556 | ! *** LEVEL OF NEUTRAL BUOYANCY *** |
|---|
| 557 | |
|---|
| 558 | inb = nl - 1 |
|---|
| 559 | DO i = icb, nl - 1 |
|---|
| 560 | ! jyg |
|---|
| 561 | ! CC IF((TVP(I)-TV(I)).LT.DTCRIT)THEN |
|---|
| 562 | IF (buoy(i)<dtovsh) THEN |
|---|
| 563 | inb = min(inb, i) |
|---|
| 564 | END IF |
|---|
| 565 | END DO |
|---|
| 566 | |
|---|
| 567 | |
|---|
| 568 | |
|---|
| 569 | |
|---|
| 570 | |
|---|
| 571 | ! *** RESET SIG(I) AND W0(I) FOR I>INB AND I<ICB *** |
|---|
| 572 | |
|---|
| 573 | IF (inb<(nl - 1)) THEN |
|---|
| 574 | DO i = inb + 1, nl - 1 |
|---|
| 575 | ! jyg |
|---|
| 576 | ! CC SIG(I)=BETA*SIG(I)-2.0E-4*ALPHA*(TV(INB)-TVP(INB))* |
|---|
| 577 | ! CC 1 ABS(TV(INB)-TVP(INB)) |
|---|
| 578 | sig(i) = beta * sig(i) + 2. * alpha * buoy(inb) * abs(buoy(inb)) |
|---|
| 579 | sig(i) = amax1(sig(i), 0.0) |
|---|
| 580 | w0(i) = beta * w0(i) |
|---|
| 581 | END DO |
|---|
| 582 | END IF |
|---|
| 583 | DO i = 1, icb |
|---|
| 584 | ! jyg |
|---|
| 585 | ! CC SIG(I)=BETA*SIG(I)-2.0E-4*ALPHA*(TV(ICB)-TVP(ICB))* |
|---|
| 586 | ! CC 1 (TV(ICB)-TVP(ICB)) |
|---|
| 587 | sig(i) = beta * sig(i) - 2. * alpha * buoy(icb) * buoy(icb) |
|---|
| 588 | sig(i) = amax1(sig(i), 0.0) |
|---|
| 589 | w0(i) = beta * w0(i) |
|---|
| 590 | END DO |
|---|
| 591 | |
|---|
| 592 | ! *** RESET FRACTIONAL AREAS OF UPDRAFTS AND W0 AT INITIAL TIME *** |
|---|
| 593 | ! *** AND AFTER 10 TIME STEPS OF NO CONVECTION *** |
|---|
| 594 | |
|---|
| 595 | IF (sig(nd)<1.5 .OR. sig(nd)>12.0) THEN |
|---|
| 596 | DO i = 1, nl - 1 |
|---|
| 597 | sig(i) = 0.0 |
|---|
| 598 | w0(i) = 0.0 |
|---|
| 599 | END DO |
|---|
| 600 | END IF |
|---|
| 601 | |
|---|
| 602 | ! *** CALCULATE LIQUID WATER STATIC ENERGY OF LIFTED PARCEL *** |
|---|
| 603 | |
|---|
| 604 | DO i = icb, inb |
|---|
| 605 | hp(i) = h(1) + (lv(i) + (cpd - cpv) * t(i)) * ep(i) * clw(i) |
|---|
| 606 | END DO |
|---|
| 607 | |
|---|
| 608 | ! *** CALCULATE CONVECTIVE AVAILABLE POTENTIAL ENERGY (CAPE), *** |
|---|
| 609 | ! *** VERTICAL VELOCITY (W), FRACTIONAL AREA COVERED BY *** |
|---|
| 610 | ! *** UNDILUTE UPDRAFT (SIG), AND UPDRAFT MASS FLUX (M) *** |
|---|
| 611 | |
|---|
| 612 | cape = 0.0 |
|---|
| 613 | |
|---|
| 614 | DO i = icb + 1, inb |
|---|
| 615 | ! jyg1 |
|---|
| 616 | ! CC CAPE=CAPE+RD*(TVP(I-1)-TV(I-1))*(PH(I-1)-PH(I))/P(I-1) |
|---|
| 617 | ! CC DCAPE=RD*BUOY(I-1)*(PH(I-1)-PH(I))/P(I-1) |
|---|
| 618 | ! CC DLNP=(PH(I-1)-PH(I))/P(I-1) |
|---|
| 619 | ! The interval on which CAPE is computed starts at PBASE : |
|---|
| 620 | deltap = min(pbase, ph(i - 1)) - min(pbase, ph(i)) |
|---|
| 621 | cape = cape + rd * buoy(i - 1) * deltap / p(i - 1) |
|---|
| 622 | dcape = rd * buoy(i - 1) * deltap / p(i - 1) |
|---|
| 623 | dlnp = deltap / p(i - 1) |
|---|
| 624 | ! jyg2 |
|---|
| 625 | ! sb3d print *,'buoy,dlnp,dcape,cape',buoy(i-1),dlnp,dcape,cape |
|---|
| 626 | ! test sb: |
|---|
| 627 | ! @ WRITE(*,*) '############################################' |
|---|
| 628 | ! @ WRITE(*,*) 'cape,rrd,buoy,deltap,p,pbase,ph:' |
|---|
| 629 | ! @ : ,cape,rd,buoy(i-1),deltap,p(i-1),pbase,ph(i) |
|---|
| 630 | ! @ WRITE(*,*) '############################################' |
|---|
| 631 | |
|---|
| 632 | ! fin test sb |
|---|
| 633 | cape = amax1(0.0, cape) |
|---|
| 634 | |
|---|
| 635 | sigold = sig(i) |
|---|
| 636 | dtmin = 100.0 |
|---|
| 637 | DO j = icb, i - 1 |
|---|
| 638 | ! jyg |
|---|
| 639 | ! CC DTMIN=AMIN1(DTMIN,(TVP(J)-TV(J))) |
|---|
| 640 | dtmin = amin1(dtmin, buoy(j)) |
|---|
| 641 | END DO |
|---|
| 642 | ! sb3d print *, 'DTMIN, BETA, ALPHA, SIG = ',DTMIN,BETA,ALPHA,SIG(I) |
|---|
| 643 | sig(i) = beta * sig(i) + alpha * dtmin * abs(dtmin) |
|---|
| 644 | sig(i) = amax1(sig(i), 0.0) |
|---|
| 645 | sig(i) = amin1(sig(i), 0.01) |
|---|
| 646 | fac = amin1(((dtcrit - dtmin) / dtcrit), 1.0) |
|---|
| 647 | ! jyg |
|---|
| 648 | ! C Essais de reduction de la vitesse |
|---|
| 649 | ! C FAC = FAC*.5 |
|---|
| 650 | |
|---|
| 651 | w = (1. - beta) * fac * sqrt(cape) + beta * w0(i) |
|---|
| 652 | amu = 0.5 * (sig(i) + sigold) * w |
|---|
| 653 | m(i) = amu * 0.007 * p(i) * (ph(i) - ph(i + 1)) / tv(i) |
|---|
| 654 | |
|---|
| 655 | ! --------- test sb: |
|---|
| 656 | ! WRITE(*,*) '############################################' |
|---|
| 657 | ! WRITE(*,*) 'k,amu,buoy(k-1),deltap,w,beta,fac,cape,w0(k)' |
|---|
| 658 | ! WRITE(*,*) i,amu,buoy(i-1),deltap |
|---|
| 659 | ! : ,w,beta,fac,cape,w0(i) |
|---|
| 660 | ! WRITE(*,*) '############################################' |
|---|
| 661 | ! --------- |
|---|
| 662 | |
|---|
| 663 | w0(i) = w |
|---|
| 664 | END DO |
|---|
| 665 | w0(icb) = 0.5 * w0(icb + 1) |
|---|
| 666 | m(icb) = 0.5 * m(icb + 1) * (ph(icb) - ph(icb + 1)) / (ph(icb + 1) - ph(icb + 2)) |
|---|
| 667 | sig(icb) = sig(icb + 1) |
|---|
| 668 | sig(icb - 1) = sig(icb) |
|---|
| 669 | ! jyg1 |
|---|
| 670 | ! sb3d print *, 'Cloud base, c. top, CAPE',ICB,INB,cape |
|---|
| 671 | ! sb3d print *, 'SIG ',(sig(i),i=1,inb) |
|---|
| 672 | ! sb3d print *, 'W ',(w0(i),i=1,inb) |
|---|
| 673 | ! sb3d print *, 'M ',(m(i), i=1,inb) |
|---|
| 674 | ! sb3d print *, 'Dt1 ',(tvp(i)-tv(i),i=1,inb) |
|---|
| 675 | ! sb3d print *, 'Dt_vrai ',(buoy(i),i=1,inb) |
|---|
| 676 | ! jyg2 |
|---|
| 677 | |
|---|
| 678 | ! *** CALCULATE ENTRAINED AIR MASS FLUX (MENT), TOTAL WATER MIXING *** |
|---|
| 679 | ! *** RATIO (QENT), TOTAL CONDENSED WATER (ELIJ), AND MIXING *** |
|---|
| 680 | ! *** FRACTION (SIJ) *** |
|---|
| 681 | |
|---|
| 682 | DO i = icb, inb |
|---|
| 683 | rti = rr(1) - ep(i) * clw(i) |
|---|
| 684 | DO j = icb - 1, inb |
|---|
| 685 | bf2 = 1. + lv(j) * lv(j) * rs(j) / (rv * t(j) * t(j) * cpd) |
|---|
| 686 | anum = h(j) - hp(i) + (cpv - cpd) * t(j) * (rti - rr(j)) |
|---|
| 687 | denom = h(i) - hp(i) + (cpd - cpv) * (rr(i) - rti) * t(j) |
|---|
| 688 | dei = denom |
|---|
| 689 | IF (abs(dei)<0.01) dei = 0.01 |
|---|
| 690 | sij(i, j) = anum / dei |
|---|
| 691 | sij(i, i) = 1.0 |
|---|
| 692 | altem = sij(i, j) * rr(i) + (1. - sij(i, j)) * rti - rs(j) |
|---|
| 693 | altem = altem / bf2 |
|---|
| 694 | cwat = clw(j) * (1. - ep(j)) |
|---|
| 695 | stemp = sij(i, j) |
|---|
| 696 | IF ((stemp<0.0 .OR. stemp>1.0 .OR. altem>cwat) .AND. j>i) THEN |
|---|
| 697 | anum = anum - lv(j) * (rti - rs(j) - cwat * bf2) |
|---|
| 698 | denom = denom + lv(j) * (rr(i) - rti) |
|---|
| 699 | IF (abs(denom)<0.01) denom = 0.01 |
|---|
| 700 | sij(i, j) = anum / denom |
|---|
| 701 | altem = sij(i, j) * rr(i) + (1. - sij(i, j)) * rti - rs(j) |
|---|
| 702 | altem = altem - (bf2 - 1.) * cwat |
|---|
| 703 | END IF |
|---|
| 704 | |
|---|
| 705 | IF (sij(i, j)>0.0 .AND. sij(i, j)<0.95) THEN |
|---|
| 706 | qent(i, j) = sij(i, j) * rr(i) + (1. - sij(i, j)) * rti |
|---|
| 707 | uent(i, j) = sij(i, j) * u(i) + (1. - sij(i, j)) * u(nk) |
|---|
| 708 | vent(i, j) = sij(i, j) * v(i) + (1. - sij(i, j)) * v(nk) |
|---|
| 709 | DO k = 1, ntra |
|---|
| 710 | traent(i, j, k) = sij(i, j) * tra(i, k) + (1. - sij(i, j)) * tra(nk, k) |
|---|
| 711 | END DO |
|---|
| 712 | elij(i, j) = altem |
|---|
| 713 | elij(i, j) = amax1(0.0, elij(i, j)) |
|---|
| 714 | ment(i, j) = m(i) / (1. - sij(i, j)) |
|---|
| 715 | nent(i) = nent(i) + 1 |
|---|
| 716 | END IF |
|---|
| 717 | sij(i, j) = amax1(0.0, sij(i, j)) |
|---|
| 718 | sij(i, j) = amin1(1.0, sij(i, j)) |
|---|
| 719 | END DO |
|---|
| 720 | |
|---|
| 721 | ! *** IF NO AIR CAN ENTRAIN AT LEVEL I ASSUME THAT UPDRAFT DETRAINS |
|---|
| 722 | ! *** |
|---|
| 723 | ! *** AT THAT LEVEL AND CALCULATE DETRAINED AIR FLUX AND PROPERTIES |
|---|
| 724 | ! *** |
|---|
| 725 | |
|---|
| 726 | IF (nent(i)==0) THEN |
|---|
| 727 | ment(i, i) = m(i) |
|---|
| 728 | qent(i, i) = rr(nk) - ep(i) * clw(i) |
|---|
| 729 | uent(i, i) = u(nk) |
|---|
| 730 | vent(i, i) = v(nk) |
|---|
| 731 | DO j = 1, ntra |
|---|
| 732 | traent(i, i, j) = tra(nk, j) |
|---|
| 733 | END DO |
|---|
| 734 | elij(i, i) = clw(i) |
|---|
| 735 | sij(i, i) = 1.0 |
|---|
| 736 | END IF |
|---|
| 737 | |
|---|
| 738 | DO j = icb - 1, inb |
|---|
| 739 | sigij(i, j) = sij(i, j) |
|---|
| 740 | END DO |
|---|
| 741 | |
|---|
| 742 | END DO |
|---|
| 743 | |
|---|
| 744 | ! *** NORMALIZE ENTRAINED AIR MASS FLUXES TO REPRESENT EQUAL *** |
|---|
| 745 | ! *** PROBABILITIES OF MIXING *** |
|---|
| 746 | |
|---|
| 747 | DO i = icb, inb |
|---|
| 748 | IF (nent(i)/=0) THEN |
|---|
| 749 | qp = rr(1) - ep(i) * clw(i) |
|---|
| 750 | anum = h(i) - hp(i) - lv(i) * (qp - rs(i)) + (cpv - cpd) * t(i) * (qp - rr(i)) |
|---|
| 751 | denom = h(i) - hp(i) + lv(i) * (rr(i) - qp) + (cpd - cpv) * t(i) * (rr(i) - qp) |
|---|
| 752 | IF (abs(denom)<0.01) denom = 0.01 |
|---|
| 753 | scrit = anum / denom |
|---|
| 754 | alt = qp - rs(i) + scrit * (rr(i) - qp) |
|---|
| 755 | IF (scrit<=0.0 .OR. alt<=0.0) scrit = 1.0 |
|---|
| 756 | smax = 0.0 |
|---|
| 757 | asij = 0.0 |
|---|
| 758 | DO j = inb, icb - 1, -1 |
|---|
| 759 | IF (sij(i, j)>1.0E-16 .AND. sij(i, j)<0.95) THEN |
|---|
| 760 | wgh = 1.0 |
|---|
| 761 | IF (j>i) THEN |
|---|
| 762 | sjmax = amax1(sij(i, j + 1), smax) |
|---|
| 763 | sjmax = amin1(sjmax, scrit) |
|---|
| 764 | smax = amax1(sij(i, j), smax) |
|---|
| 765 | sjmin = amax1(sij(i, j - 1), smax) |
|---|
| 766 | sjmin = amin1(sjmin, scrit) |
|---|
| 767 | IF (sij(i, j)<(smax - 1.0E-16)) wgh = 0.0 |
|---|
| 768 | smid = amin1(sij(i, j), scrit) |
|---|
| 769 | ELSE |
|---|
| 770 | sjmax = amax1(sij(i, j + 1), scrit) |
|---|
| 771 | smid = amax1(sij(i, j), scrit) |
|---|
| 772 | sjmin = 0.0 |
|---|
| 773 | IF (j>1) sjmin = sij(i, j - 1) |
|---|
| 774 | sjmin = amax1(sjmin, scrit) |
|---|
| 775 | END IF |
|---|
| 776 | delp = abs(sjmax - smid) |
|---|
| 777 | delm = abs(sjmin - smid) |
|---|
| 778 | asij = asij + wgh * (delp + delm) |
|---|
| 779 | ment(i, j) = ment(i, j) * (delp + delm) * wgh |
|---|
| 780 | END IF |
|---|
| 781 | END DO |
|---|
| 782 | asij = amax1(1.0E-16, asij) |
|---|
| 783 | asij = 1.0 / asij |
|---|
| 784 | DO j = icb - 1, inb |
|---|
| 785 | ment(i, j) = ment(i, j) * asij |
|---|
| 786 | END DO |
|---|
| 787 | asum = 0.0 |
|---|
| 788 | bsum = 0.0 |
|---|
| 789 | DO j = icb - 1, inb |
|---|
| 790 | asum = asum + ment(i, j) |
|---|
| 791 | ment(i, j) = ment(i, j) * sig(j) |
|---|
| 792 | bsum = bsum + ment(i, j) |
|---|
| 793 | END DO |
|---|
| 794 | bsum = amax1(bsum, 1.0E-16) |
|---|
| 795 | bsum = 1.0 / bsum |
|---|
| 796 | DO j = icb - 1, inb |
|---|
| 797 | ment(i, j) = ment(i, j) * asum * bsum |
|---|
| 798 | END DO |
|---|
| 799 | csum = 0.0 |
|---|
| 800 | DO j = icb - 1, inb |
|---|
| 801 | csum = csum + ment(i, j) |
|---|
| 802 | END DO |
|---|
| 803 | |
|---|
| 804 | IF (csum<m(i)) THEN |
|---|
| 805 | nent(i) = 0 |
|---|
| 806 | ment(i, i) = m(i) |
|---|
| 807 | qent(i, i) = rr(1) - ep(i) * clw(i) |
|---|
| 808 | uent(i, i) = u(nk) |
|---|
| 809 | vent(i, i) = v(nk) |
|---|
| 810 | DO j = 1, ntra |
|---|
| 811 | traent(i, i, j) = tra(nk, j) |
|---|
| 812 | END DO |
|---|
| 813 | elij(i, i) = clw(i) |
|---|
| 814 | sij(i, i) = 1.0 |
|---|
| 815 | END IF |
|---|
| 816 | END IF |
|---|
| 817 | END DO |
|---|
| 818 | |
|---|
| 819 | |
|---|
| 820 | |
|---|
| 821 | ! ************************************************************** |
|---|
| 822 | ! * CALCUL DES MENTS(I,J) ET DES QENTS(I,J) |
|---|
| 823 | ! ************************************************************* |
|---|
| 824 | |
|---|
| 825 | DO im = 1, nd |
|---|
| 826 | DO jm = 1, nd |
|---|
| 827 | |
|---|
| 828 | qents(im, jm) = qent(im, jm) |
|---|
| 829 | ments(im, jm) = ment(im, jm) |
|---|
| 830 | END DO |
|---|
| 831 | END DO |
|---|
| 832 | |
|---|
| 833 | ! ********************************************************** |
|---|
| 834 | ! --- test sb: |
|---|
| 835 | ! @ WRITE(*,*) '^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^' |
|---|
| 836 | ! @ WRITE(*,*) 'inb,m(inb),ment(inb,inb),sigij(inb,inb):' |
|---|
| 837 | ! @ WRITE(*,*) inb,m(inb),ment(inb,inb),sigij(inb,inb) |
|---|
| 838 | ! @ WRITE(*,*) '^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^' |
|---|
| 839 | ! --- |
|---|
| 840 | |
|---|
| 841 | |
|---|
| 842 | |
|---|
| 843 | |
|---|
| 844 | |
|---|
| 845 | ! *** CHECK WHETHER EP(INB)=0, IF SO, SKIP PRECIPITATING *** |
|---|
| 846 | ! *** DOWNDRAFT CALCULATION *** |
|---|
| 847 | |
|---|
| 848 | IF (ep(inb)<0.0001) GO TO 405 |
|---|
| 849 | |
|---|
| 850 | ! *** INTEGRATE LIQUID WATER EQUATION TO FIND CONDENSED WATER *** |
|---|
| 851 | ! *** AND CONDENSED WATER FLUX *** |
|---|
| 852 | |
|---|
| 853 | wflux = 0.0 |
|---|
| 854 | tinv = 1. / 3. |
|---|
| 855 | |
|---|
| 856 | ! *** BEGIN DOWNDRAFT LOOP *** |
|---|
| 857 | |
|---|
| 858 | DO i = inb, 1, -1 |
|---|
| 859 | |
|---|
| 860 | ! *** CALCULATE DETRAINED PRECIPITATION *** |
|---|
| 861 | |
|---|
| 862 | wdtrain = 10.0 * ep(i) * m(i) * clw(i) |
|---|
| 863 | IF (i>1) THEN |
|---|
| 864 | DO j = 1, i - 1 |
|---|
| 865 | awat = elij(j, i) - (1. - ep(i)) * clw(i) |
|---|
| 866 | awat = amax1(awat, 0.0) |
|---|
| 867 | wdtrain = wdtrain + 10.0 * awat * ment(j, i) |
|---|
| 868 | END DO |
|---|
| 869 | END IF |
|---|
| 870 | |
|---|
| 871 | ! *** FIND RAIN WATER AND EVAPORATION USING PROVISIONAL *** |
|---|
| 872 | ! *** ESTIMATES OF RP(I)AND RP(I-1) *** |
|---|
| 873 | |
|---|
| 874 | wt(i) = 45.0 |
|---|
| 875 | IF (i<inb) THEN |
|---|
| 876 | rp(i) = rp(i + 1) + (cpd * (t(i + 1) - t(i)) + gz(i + 1) - gz(i)) / lv(i) |
|---|
| 877 | rp(i) = 0.5 * (rp(i) + rr(i)) |
|---|
| 878 | END IF |
|---|
| 879 | rp(i) = amax1(rp(i), 0.0) |
|---|
| 880 | rp(i) = amin1(rp(i), rs(i)) |
|---|
| 881 | rp(inb) = rr(inb) |
|---|
| 882 | IF (i==1) THEN |
|---|
| 883 | afac = p(1) * (rs(1) - rp(1)) / (1.0E4 + 2000.0 * p(1) * rs(1)) |
|---|
| 884 | ELSE |
|---|
| 885 | rp(i - 1) = rp(i) + (cpd * (t(i) - t(i - 1)) + gz(i) - gz(i - 1)) / lv(i) |
|---|
| 886 | rp(i - 1) = 0.5 * (rp(i - 1) + rr(i - 1)) |
|---|
| 887 | rp(i - 1) = amin1(rp(i - 1), rs(i - 1)) |
|---|
| 888 | rp(i - 1) = amax1(rp(i - 1), 0.0) |
|---|
| 889 | afac1 = p(i) * (rs(i) - rp(i)) / (1.0E4 + 2000.0 * p(i) * rs(i)) |
|---|
| 890 | afac2 = p(i - 1) * (rs(i - 1) - rp(i - 1)) / (1.0E4 + 2000.0 * p(i - 1) * rs(i - 1)) |
|---|
| 891 | afac = 0.5 * (afac1 + afac2) |
|---|
| 892 | END IF |
|---|
| 893 | IF (i==inb) afac = 0.0 |
|---|
| 894 | afac = amax1(afac, 0.0) |
|---|
| 895 | bfac = 1. / (sigd * wt(i)) |
|---|
| 896 | |
|---|
| 897 | ! jyg1 |
|---|
| 898 | ! CC SIGT=1.0 |
|---|
| 899 | ! CC IF(I.GE.ICB)SIGT=SIGP(I) |
|---|
| 900 | ! Prise en compte de la variation progressive de SIGT dans |
|---|
| 901 | ! les couches ICB et ICB-1: |
|---|
| 902 | ! Pour PLCL<PH(I+1), PR1=0 & PR2=1 |
|---|
| 903 | ! Pour PLCL>PH(I), PR1=1 & PR2=0 |
|---|
| 904 | ! Pour PH(I+1)<PLCL<PH(I), PR1 est la proportion a cheval |
|---|
| 905 | ! sur le nuage, et PR2 est la proportion sous la base du |
|---|
| 906 | ! nuage. |
|---|
| 907 | pr1 = (plcl - ph(i + 1)) / (ph(i) - ph(i + 1)) |
|---|
| 908 | pr1 = max(0., min(1., pr1)) |
|---|
| 909 | pr2 = (ph(i) - plcl) / (ph(i) - ph(i + 1)) |
|---|
| 910 | pr2 = max(0., min(1., pr2)) |
|---|
| 911 | sigt = sigp(i) * pr1 + pr2 |
|---|
| 912 | ! sb3d print *,'i,sigt,pr1,pr2', i,sigt,pr1,pr2 |
|---|
| 913 | ! jyg2 |
|---|
| 914 | |
|---|
| 915 | b6 = bfac * 50. * sigd * (ph(i) - ph(i + 1)) * sigt * afac |
|---|
| 916 | c6 = water(i + 1) + bfac * wdtrain - 50. * sigd * bfac * (ph(i) - ph(i + 1)) * evap(i + 1) |
|---|
| 917 | IF (c6>0.0) THEN |
|---|
| 918 | revap = 0.5 * (-b6 + sqrt(b6 * b6 + 4. * c6)) |
|---|
| 919 | evap(i) = sigt * afac * revap |
|---|
| 920 | water(i) = revap * revap |
|---|
| 921 | ELSE |
|---|
| 922 | evap(i) = -evap(i + 1) + 0.02 * (wdtrain + sigd * wt(i) * water(i + 1)) / (sigd * (ph(i & |
|---|
| 923 | ) - ph(i + 1))) |
|---|
| 924 | END IF |
|---|
| 925 | |
|---|
| 926 | |
|---|
| 927 | |
|---|
| 928 | ! *** CALCULATE PRECIPITATING DOWNDRAFT MASS FLUX UNDER *** |
|---|
| 929 | ! *** HYDROSTATIC APPROXIMATION *** |
|---|
| 930 | |
|---|
| 931 | IF (i==1) GO TO 360 |
|---|
| 932 | tevap = amax1(0.0, evap(i)) |
|---|
| 933 | delth = amax1(0.001, (th(i) - th(i - 1))) |
|---|
| 934 | mp(i) = 10. * lvcp(i) * sigd * tevap * (p(i - 1) - p(i)) / delth |
|---|
| 935 | |
|---|
| 936 | ! *** IF HYDROSTATIC ASSUMPTION FAILS, *** |
|---|
| 937 | ! *** SOLVE CUBIC DIFFERENCE EQUATION FOR DOWNDRAFT THETA *** |
|---|
| 938 | ! *** AND MASS FLUX FROM TWO SIMULTANEOUS DIFFERENTIAL EQNS *** |
|---|
| 939 | |
|---|
| 940 | amfac = sigd * sigd * 70.0 * ph(i) * (p(i - 1) - p(i)) * (th(i) - th(i - 1)) / (tv(i) * th(i)) |
|---|
| 941 | amp2 = abs(mp(i + 1) * mp(i + 1) - mp(i) * mp(i)) |
|---|
| 942 | IF (amp2>(0.1 * amfac)) THEN |
|---|
| 943 | xf = 100.0 * sigd * sigd * sigd * (ph(i) - ph(i + 1)) |
|---|
| 944 | tf = b(i) - 5.0 * (th(i) - th(i - 1)) * t(i) / (lvcp(i) * sigd * th(i)) |
|---|
| 945 | af = xf * tf + mp(i + 1) * mp(i + 1) * tinv |
|---|
| 946 | bf = 2. * (tinv * mp(i + 1))**3 + tinv * mp(i + 1) * xf * tf + & |
|---|
| 947 | 50. * (p(i - 1) - p(i)) * xf * tevap |
|---|
| 948 | fac2 = 1.0 |
|---|
| 949 | IF (bf<0.0) fac2 = -1.0 |
|---|
| 950 | bf = abs(bf) |
|---|
| 951 | ur = 0.25 * bf * bf - af * af * af * tinv * tinv * tinv |
|---|
| 952 | IF (ur>=0.0) THEN |
|---|
| 953 | sru = sqrt(ur) |
|---|
| 954 | fac = 1.0 |
|---|
| 955 | IF ((0.5 * bf - sru)<0.0) fac = -1.0 |
|---|
| 956 | mp(i) = mp(i + 1) * tinv + (0.5 * bf + sru)**tinv + & |
|---|
| 957 | fac * (abs(0.5 * bf - sru))**tinv |
|---|
| 958 | ELSE |
|---|
| 959 | d = atan(2. * sqrt(-ur) / (bf + 1.0E-28)) |
|---|
| 960 | IF (fac2<0.0) d = 3.14159 - d |
|---|
| 961 | mp(i) = mp(i + 1) * tinv + 2. * sqrt(af * tinv) * cos(d * tinv) |
|---|
| 962 | END IF |
|---|
| 963 | mp(i) = amax1(0.0, mp(i)) |
|---|
| 964 | b(i - 1) = b(i) + 100.0 * (p(i - 1) - p(i)) * tevap / (mp(i) + sigd * 0.1) - & |
|---|
| 965 | 10.0 * (th(i) - th(i - 1)) * t(i) / (lvcp(i) * sigd * th(i)) |
|---|
| 966 | b(i - 1) = amax1(b(i - 1), 0.0) |
|---|
| 967 | END IF |
|---|
| 968 | |
|---|
| 969 | |
|---|
| 970 | |
|---|
| 971 | ! *** LIMIT MAGNITUDE OF MP(I) TO MEET CFL CONDITION *** |
|---|
| 972 | |
|---|
| 973 | ampmax = 2.0 * (ph(i) - ph(i + 1)) * delti |
|---|
| 974 | amp2 = 2.0 * (ph(i - 1) - ph(i)) * delti |
|---|
| 975 | ampmax = amin1(ampmax, amp2) |
|---|
| 976 | mp(i) = amin1(mp(i), ampmax) |
|---|
| 977 | |
|---|
| 978 | ! *** FORCE MP TO DECREASE LINEARLY TO ZERO *** |
|---|
| 979 | ! *** BETWEEN CLOUD BASE AND THE SURFACE *** |
|---|
| 980 | |
|---|
| 981 | IF (p(i)>p(icb)) THEN |
|---|
| 982 | mp(i) = mp(icb) * (p(1) - p(i)) / (p(1) - p(icb)) |
|---|
| 983 | END IF |
|---|
| 984 | 360 CONTINUE |
|---|
| 985 | |
|---|
| 986 | ! *** FIND MIXING RATIO OF PRECIPITATING DOWNDRAFT *** |
|---|
| 987 | |
|---|
| 988 | IF (i==inb) GO TO 400 |
|---|
| 989 | rp(i) = rr(i) |
|---|
| 990 | IF (mp(i)>mp(i + 1)) THEN |
|---|
| 991 | rp(i) = rp(i + 1) * mp(i + 1) + rr(i) * (mp(i) - mp(i + 1)) + & |
|---|
| 992 | 5. * sigd * (ph(i) - ph(i + 1)) * (evap(i + 1) + evap(i)) |
|---|
| 993 | rp(i) = rp(i) / mp(i) |
|---|
| 994 | up(i) = up(i + 1) * mp(i + 1) + u(i) * (mp(i) - mp(i + 1)) |
|---|
| 995 | up(i) = up(i) / mp(i) |
|---|
| 996 | vp(i) = vp(i + 1) * mp(i + 1) + v(i) * (mp(i) - mp(i + 1)) |
|---|
| 997 | vp(i) = vp(i) / mp(i) |
|---|
| 998 | DO j = 1, ntra |
|---|
| 999 | trap(i, j) = trap(i + 1, j) * mp(i + 1) + trap(i, j) * (mp(i) - mp(i + 1)) |
|---|
| 1000 | trap(i, j) = trap(i, j) / mp(i) |
|---|
| 1001 | END DO |
|---|
| 1002 | ELSE |
|---|
| 1003 | IF (mp(i + 1)>1.0E-16) THEN |
|---|
| 1004 | rp(i) = rp(i + 1) + 5.0 * sigd * (ph(i) - ph(i + 1)) * (evap(i + 1) + evap(i)) / mp(i + 1 & |
|---|
| 1005 | ) |
|---|
| 1006 | up(i) = up(i + 1) |
|---|
| 1007 | vp(i) = vp(i + 1) |
|---|
| 1008 | DO j = 1, ntra |
|---|
| 1009 | trap(i, j) = trap(i + 1, j) |
|---|
| 1010 | END DO |
|---|
| 1011 | END IF |
|---|
| 1012 | END IF |
|---|
| 1013 | rp(i) = amin1(rp(i), rs(i)) |
|---|
| 1014 | rp(i) = amax1(rp(i), 0.0) |
|---|
| 1015 | 400 END DO |
|---|
| 1016 | |
|---|
| 1017 | ! *** CALCULATE SURFACE PRECIPITATION IN MM/DAY *** |
|---|
| 1018 | |
|---|
| 1019 | precip = wt(1) * sigd * water(1) * 8640.0 |
|---|
| 1020 | |
|---|
| 1021 | ! sb *** Calculate downdraft velocity scale and surface temperature and |
|---|
| 1022 | ! *** |
|---|
| 1023 | ! sb *** water vapor fluctuations |
|---|
| 1024 | ! *** |
|---|
| 1025 | ! sb (inspire de convect 4.3) |
|---|
| 1026 | |
|---|
| 1027 | ! BETAD=10.0 |
|---|
| 1028 | betad = 5.0 |
|---|
| 1029 | wd = betad * abs(mp(icb)) * 0.01 * rd * t(icb) / (sigd * p(icb)) |
|---|
| 1030 | |
|---|
| 1031 | 405 CONTINUE |
|---|
| 1032 | |
|---|
| 1033 | ! *** CALCULATE TENDENCIES OF LOWEST LEVEL POTENTIAL TEMPERATURE *** |
|---|
| 1034 | ! *** AND MIXING RATIO *** |
|---|
| 1035 | |
|---|
| 1036 | dpinv = 1.0 / (ph(1) - ph(2)) |
|---|
| 1037 | am = 0.0 |
|---|
| 1038 | DO k = 2, inb |
|---|
| 1039 | am = am + m(k) |
|---|
| 1040 | END DO |
|---|
| 1041 | IF ((0.1 * dpinv * am)>=delti) iflag = 4 |
|---|
| 1042 | ft(1) = 0.1 * dpinv * am * (t(2) - t(1) + (gz(2) - gz(1)) / cpn(1)) |
|---|
| 1043 | ft(1) = ft(1) - 0.5 * lvcp(1) * sigd * (evap(1) + evap(2)) |
|---|
| 1044 | ft(1) = ft(1) - 0.09 * sigd * mp(2) * t(1) * b(1) * dpinv |
|---|
| 1045 | ft(1) = ft(1) + 0.01 * sigd * wt(1) * (cl - cpd) * water(2) * (t(2) - t(1)) * dpinv / cpn(1) |
|---|
| 1046 | fr(1) = 0.1 * mp(2) * (rp(2) - rr(1)) * & ! correction bug conservation eau |
|---|
| 1047 | ! 1 DPINV+SIGD*0.5*(EVAP(1)+EVAP(2)) |
|---|
| 1048 | dpinv + sigd * 0.5 * (evap(1) + evap(2)) |
|---|
| 1049 | ! IM cf. SBL |
|---|
| 1050 | ! 1 DPINV+SIGD*EVAP(1) |
|---|
| 1051 | fr(1) = fr(1) + 0.1 * am * (rr(2) - rr(1)) * dpinv |
|---|
| 1052 | fu(1) = fu(1) + 0.1 * dpinv * (mp(2) * (up(2) - u(1)) + am * (u(2) - u(1))) |
|---|
| 1053 | fv(1) = fv(1) + 0.1 * dpinv * (mp(2) * (vp(2) - v(1)) + am * (v(2) - v(1))) |
|---|
| 1054 | DO j = 1, ntra |
|---|
| 1055 | ftra(1, j) = ftra(1, j) + 0.1 * dpinv * (mp(2) * (trap(2, j) - tra(1, & |
|---|
| 1056 | j)) + am * (tra(2, j) - tra(1, j))) |
|---|
| 1057 | END DO |
|---|
| 1058 | amde = 0.0 |
|---|
| 1059 | DO j = 2, inb |
|---|
| 1060 | fr(1) = fr(1) + 0.1 * dpinv * ment(j, 1) * (qent(j, 1) - rr(1)) |
|---|
| 1061 | fu(1) = fu(1) + 0.1 * dpinv * ment(j, 1) * (uent(j, 1) - u(1)) |
|---|
| 1062 | fv(1) = fv(1) + 0.1 * dpinv * ment(j, 1) * (vent(j, 1) - v(1)) |
|---|
| 1063 | DO k = 1, ntra |
|---|
| 1064 | ftra(1, k) = ftra(1, k) + 0.1 * dpinv * ment(j, 1) * (traent(j, 1, k) - tra(1, k)) |
|---|
| 1065 | END DO |
|---|
| 1066 | END DO |
|---|
| 1067 | |
|---|
| 1068 | ! *** CALCULATE TENDENCIES OF POTENTIAL TEMPERATURE AND MIXING RATIO *** |
|---|
| 1069 | ! *** AT LEVELS ABOVE THE LOWEST LEVEL *** |
|---|
| 1070 | |
|---|
| 1071 | ! *** FIRST FIND THE NET SATURATED UPDRAFT AND DOWNDRAFT MASS FLUXES *** |
|---|
| 1072 | ! *** THROUGH EACH LEVEL *** |
|---|
| 1073 | |
|---|
| 1074 | DO i = 2, inb |
|---|
| 1075 | dpinv = 1.0 / (ph(i) - ph(i + 1)) |
|---|
| 1076 | cpinv = 1.0 / cpn(i) |
|---|
| 1077 | amp1 = 0.0 |
|---|
| 1078 | DO k = i + 1, inb + 1 |
|---|
| 1079 | amp1 = amp1 + m(k) |
|---|
| 1080 | END DO |
|---|
| 1081 | DO k = 1, i |
|---|
| 1082 | DO j = i + 1, inb + 1 |
|---|
| 1083 | amp1 = amp1 + ment(k, j) |
|---|
| 1084 | END DO |
|---|
| 1085 | END DO |
|---|
| 1086 | IF ((0.1 * dpinv * amp1)>=delti) iflag = 4 |
|---|
| 1087 | ad = 0.0 |
|---|
| 1088 | DO k = 1, i - 1 |
|---|
| 1089 | DO j = i, inb |
|---|
| 1090 | ad = ad + ment(j, k) |
|---|
| 1091 | END DO |
|---|
| 1092 | END DO |
|---|
| 1093 | ft(i) = 0.1 * dpinv * (amp1 * (t(i + 1) - t(i) + (gz(i + 1) - gz(i)) * cpinv) - ad * (t(i) - t(i - & |
|---|
| 1094 | 1) + (gz(i) - gz(i - 1)) * cpinv)) - 0.5 * sigd * lvcp(i) * (evap(i) + evap(i + 1)) |
|---|
| 1095 | rat = cpn(i - 1) * cpinv |
|---|
| 1096 | ft(i) = ft(i) - 0.09 * sigd * (mp(i + 1) * t(i) * b(i) - mp(i) * t(i - 1) * rat * b(i - 1)) * & |
|---|
| 1097 | dpinv |
|---|
| 1098 | ft(i) = ft(i) + 0.1 * dpinv * ment(i, i) * (hp(i) - h(i) + t(i) * (cpv - cpd) * (rr(i) - & |
|---|
| 1099 | qent(i, i))) * cpinv |
|---|
| 1100 | ft(i) = ft(i) + 0.01 * sigd * wt(i) * (cl - cpd) * water(i + 1) * (t(i + 1) - t(i)) * dpinv * & |
|---|
| 1101 | cpinv |
|---|
| 1102 | fr(i) = 0.1 * dpinv * (amp1 * (rr(i + 1) - rr(i)) - ad * (rr(i) - rr(i - 1))) |
|---|
| 1103 | fu(i) = fu(i) + 0.1 * dpinv * (amp1 * (u(i + 1) - u(i)) - ad * (u(i) - u(i - 1))) |
|---|
| 1104 | fv(i) = fv(i) + 0.1 * dpinv * (amp1 * (v(i + 1) - v(i)) - ad * (v(i) - v(i - 1))) |
|---|
| 1105 | DO k = 1, ntra |
|---|
| 1106 | ftra(i, k) = ftra(i, k) + 0.1 * dpinv * (amp1 * (tra(i + 1, k) - tra(i, & |
|---|
| 1107 | k)) - ad * (tra(i, k) - tra(i - 1, k))) |
|---|
| 1108 | END DO |
|---|
| 1109 | DO k = 1, i - 1 |
|---|
| 1110 | awat = elij(k, i) - (1. - ep(i)) * clw(i) |
|---|
| 1111 | awat = amax1(awat, 0.0) |
|---|
| 1112 | fr(i) = fr(i) + 0.1 * dpinv * ment(k, i) * (qent(k, i) - awat - rr(i)) |
|---|
| 1113 | fu(i) = fu(i) + 0.1 * dpinv * ment(k, i) * (uent(k, i) - u(i)) |
|---|
| 1114 | fv(i) = fv(i) + 0.1 * dpinv * ment(k, i) * (vent(k, i) - v(i)) |
|---|
| 1115 | ! (saturated updrafts resulting from mixing) ! cld |
|---|
| 1116 | qcond(i) = qcond(i) + (elij(k, i) - awat) ! cld |
|---|
| 1117 | nqcond(i) = nqcond(i) + 1. ! cld |
|---|
| 1118 | DO j = 1, ntra |
|---|
| 1119 | ftra(i, j) = ftra(i, j) + 0.1 * dpinv * ment(k, i) * (traent(k, i, j) - tra(i, j & |
|---|
| 1120 | )) |
|---|
| 1121 | END DO |
|---|
| 1122 | END DO |
|---|
| 1123 | DO k = i, inb |
|---|
| 1124 | fr(i) = fr(i) + 0.1 * dpinv * ment(k, i) * (qent(k, i) - rr(i)) |
|---|
| 1125 | fu(i) = fu(i) + 0.1 * dpinv * ment(k, i) * (uent(k, i) - u(i)) |
|---|
| 1126 | fv(i) = fv(i) + 0.1 * dpinv * ment(k, i) * (vent(k, i) - v(i)) |
|---|
| 1127 | DO j = 1, ntra |
|---|
| 1128 | ftra(i, j) = ftra(i, j) + 0.1 * dpinv * ment(k, i) * (traent(k, i, j) - tra(i, j & |
|---|
| 1129 | )) |
|---|
| 1130 | END DO |
|---|
| 1131 | END DO |
|---|
| 1132 | fr(i) = fr(i) + 0.5 * sigd * (evap(i) + evap(i + 1)) + 0.1 * (mp(i + 1) * (rp(i + & |
|---|
| 1133 | 1) - rr(i)) - mp(i) * (rp(i) - rr(i - 1))) * dpinv |
|---|
| 1134 | fu(i) = fu(i) + 0.1 * (mp(i + 1) * (up(i + 1) - u(i)) - mp(i) * (up(i) - u(i - 1))) * dpinv |
|---|
| 1135 | fv(i) = fv(i) + 0.1 * (mp(i + 1) * (vp(i + 1) - v(i)) - mp(i) * (vp(i) - v(i - 1))) * dpinv |
|---|
| 1136 | DO j = 1, ntra |
|---|
| 1137 | ftra(i, j) = ftra(i, j) + 0.1 * dpinv * (mp(i + 1) * (trap(i + 1, j) - tra(i, & |
|---|
| 1138 | j)) - mp(i) * (trap(i, j) - trap(i - 1, j))) |
|---|
| 1139 | END DO |
|---|
| 1140 | ! (saturated downdrafts resulting from mixing) ! cld |
|---|
| 1141 | DO k = i + 1, inb ! cld |
|---|
| 1142 | qcond(i) = qcond(i) + elij(k, i) ! cld |
|---|
| 1143 | nqcond(i) = nqcond(i) + 1. ! cld |
|---|
| 1144 | END DO ! cld |
|---|
| 1145 | ! (particular case: no detraining level is found) ! cld |
|---|
| 1146 | IF (nent(i)==0) THEN ! cld |
|---|
| 1147 | qcond(i) = qcond(i) + (1 - ep(i)) * clw(i) ! cld |
|---|
| 1148 | nqcond(i) = nqcond(i) + 1. ! cld |
|---|
| 1149 | END IF ! cld |
|---|
| 1150 | IF (nqcond(i)/=0.) THEN ! cld |
|---|
| 1151 | qcond(i) = qcond(i) / nqcond(i) ! cld |
|---|
| 1152 | END IF ! cld |
|---|
| 1153 | END DO |
|---|
| 1154 | |
|---|
| 1155 | |
|---|
| 1156 | |
|---|
| 1157 | |
|---|
| 1158 | ! *** MOVE THE DETRAINMENT AT LEVEL INB DOWN TO LEVEL INB-1 *** |
|---|
| 1159 | ! *** IN SUCH A WAY AS TO PRESERVE THE VERTICALLY *** |
|---|
| 1160 | ! *** INTEGRATED ENTHALPY AND WATER TENDENCIES *** |
|---|
| 1161 | |
|---|
| 1162 | ! test sb: |
|---|
| 1163 | ! @ WRITE(*,*) '--------------------------------------------' |
|---|
| 1164 | ! @ WRITE(*,*) 'inb,ft,hp,h,t,rr,qent,ment,water,waterp,wt,mp,b' |
|---|
| 1165 | ! @ WRITE(*,*) inb,ft(inb),hp(inb),h(inb) |
|---|
| 1166 | ! @ : ,t(inb),rr(inb),qent(inb,inb) |
|---|
| 1167 | ! @ : ,ment(inb,inb),water(inb) |
|---|
| 1168 | ! @ : ,water(inb+1),wt(inb),mp(inb),b(inb) |
|---|
| 1169 | ! @ WRITE(*,*) '--------------------------------------------' |
|---|
| 1170 | ! fin test sb: |
|---|
| 1171 | |
|---|
| 1172 | ax = 0.1 * ment(inb, inb) * (hp(inb) - h(inb) + t(inb) * (cpv - cpd) * (rr(inb) - qent(inb, & |
|---|
| 1173 | inb))) / (cpn(inb) * (ph(inb) - ph(inb + 1))) |
|---|
| 1174 | ft(inb) = ft(inb) - ax |
|---|
| 1175 | ft(inb - 1) = ft(inb - 1) + ax * cpn(inb) * (ph(inb) - ph(inb + 1)) / (cpn(inb - 1) * (ph(inb & |
|---|
| 1176 | - 1) - ph(inb))) |
|---|
| 1177 | bx = 0.1 * ment(inb, inb) * (qent(inb, inb) - rr(inb)) / (ph(inb) - ph(inb + 1)) |
|---|
| 1178 | fr(inb) = fr(inb) - bx |
|---|
| 1179 | fr(inb - 1) = fr(inb - 1) + bx * (ph(inb) - ph(inb + 1)) / (ph(inb - 1) - ph(inb)) |
|---|
| 1180 | cx = 0.1 * ment(inb, inb) * (uent(inb, inb) - u(inb)) / (ph(inb) - ph(inb + 1)) |
|---|
| 1181 | fu(inb) = fu(inb) - cx |
|---|
| 1182 | fu(inb - 1) = fu(inb - 1) + cx * (ph(inb) - ph(inb + 1)) / (ph(inb - 1) - ph(inb)) |
|---|
| 1183 | dx = 0.1 * ment(inb, inb) * (vent(inb, inb) - v(inb)) / (ph(inb) - ph(inb + 1)) |
|---|
| 1184 | fv(inb) = fv(inb) - dx |
|---|
| 1185 | fv(inb - 1) = fv(inb - 1) + dx * (ph(inb) - ph(inb + 1)) / (ph(inb - 1) - ph(inb)) |
|---|
| 1186 | DO j = 1, ntra |
|---|
| 1187 | ex = 0.1 * ment(inb, inb) * (traent(inb, inb, j) - tra(inb, j)) / & |
|---|
| 1188 | (ph(inb) - ph(inb + 1)) |
|---|
| 1189 | ftra(inb, j) = ftra(inb, j) - ex |
|---|
| 1190 | ftra(inb - 1, j) = ftra(inb - 1, j) + ex * (ph(inb) - ph(inb + 1)) / (ph(inb - 1) - ph(& |
|---|
| 1191 | inb)) |
|---|
| 1192 | END DO |
|---|
| 1193 | |
|---|
| 1194 | ! *** HOMOGINIZE TENDENCIES BELOW CLOUD BASE *** |
|---|
| 1195 | |
|---|
| 1196 | asum = 0.0 |
|---|
| 1197 | bsum = 0.0 |
|---|
| 1198 | csum = 0.0 |
|---|
| 1199 | dsum = 0.0 |
|---|
| 1200 | DO i = 1, icb - 1 |
|---|
| 1201 | asum = asum + ft(i) * (ph(i) - ph(i + 1)) |
|---|
| 1202 | bsum = bsum + fr(i) * (lv(i) + (cl - cpd) * (t(i) - t(1))) * (ph(i) - ph(i + 1)) |
|---|
| 1203 | csum = csum + (lv(i) + (cl - cpd) * (t(i) - t(1))) * (ph(i) - ph(i + 1)) |
|---|
| 1204 | dsum = dsum + t(i) * (ph(i) - ph(i + 1)) / th(i) |
|---|
| 1205 | END DO |
|---|
| 1206 | DO i = 1, icb - 1 |
|---|
| 1207 | ft(i) = asum * t(i) / (th(i) * dsum) |
|---|
| 1208 | fr(i) = bsum / csum |
|---|
| 1209 | END DO |
|---|
| 1210 | |
|---|
| 1211 | ! *** RESET COUNTER AND RETURN *** |
|---|
| 1212 | |
|---|
| 1213 | sig(nd) = 2.0 |
|---|
| 1214 | |
|---|
| 1215 | DO i = 1, nd |
|---|
| 1216 | upwd(i) = 0.0 |
|---|
| 1217 | dnwd(i) = 0.0 |
|---|
| 1218 | ! sb dnwd0(i) = - mp(i) |
|---|
| 1219 | END DO |
|---|
| 1220 | |
|---|
| 1221 | DO i = 1, nl |
|---|
| 1222 | dnwd0(i) = -mp(i) |
|---|
| 1223 | END DO |
|---|
| 1224 | DO i = nl + 1, nd |
|---|
| 1225 | dnwd0(i) = 0. |
|---|
| 1226 | END DO |
|---|
| 1227 | |
|---|
| 1228 | DO i = icb, inb |
|---|
| 1229 | upwd(i) = 0.0 |
|---|
| 1230 | dnwd(i) = 0.0 |
|---|
| 1231 | |
|---|
| 1232 | DO k = i, inb |
|---|
| 1233 | up1 = 0.0 |
|---|
| 1234 | dn1 = 0.0 |
|---|
| 1235 | DO n = 1, i - 1 |
|---|
| 1236 | up1 = up1 + ment(n, k) |
|---|
| 1237 | dn1 = dn1 - ment(k, n) |
|---|
| 1238 | END DO |
|---|
| 1239 | upwd(i) = upwd(i) + m(k) + up1 |
|---|
| 1240 | dnwd(i) = dnwd(i) + dn1 |
|---|
| 1241 | END DO |
|---|
| 1242 | END DO |
|---|
| 1243 | |
|---|
| 1244 | ! ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc |
|---|
| 1245 | ! DETERMINATION DE LA VARIATION DE FLUX ASCENDANT ENTRE |
|---|
| 1246 | ! DEUX NIVEAU NON DILUE Mike |
|---|
| 1247 | ! ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc |
|---|
| 1248 | |
|---|
| 1249 | |
|---|
| 1250 | ! sb do i=1,ND |
|---|
| 1251 | ! sb Mike(i)=M(i) |
|---|
| 1252 | ! sb enddo |
|---|
| 1253 | |
|---|
| 1254 | DO i = 1, nl |
|---|
| 1255 | mike(i) = m(i) |
|---|
| 1256 | END DO |
|---|
| 1257 | DO i = nl + 1, nd |
|---|
| 1258 | mike(i) = 0. |
|---|
| 1259 | END DO |
|---|
| 1260 | |
|---|
| 1261 | DO i = 1, nd |
|---|
| 1262 | ma(i) = 0 |
|---|
| 1263 | END DO |
|---|
| 1264 | |
|---|
| 1265 | ! sb do i=1,nd |
|---|
| 1266 | ! sb do j=i,nd |
|---|
| 1267 | ! sb Ma(i)=Ma(i)+M(j) |
|---|
| 1268 | ! sb enddo |
|---|
| 1269 | ! sb enddo |
|---|
| 1270 | |
|---|
| 1271 | DO i = 1, nl |
|---|
| 1272 | DO j = i, nl |
|---|
| 1273 | ma(i) = ma(i) + m(j) |
|---|
| 1274 | END DO |
|---|
| 1275 | END DO |
|---|
| 1276 | |
|---|
| 1277 | DO i = nl + 1, nd |
|---|
| 1278 | ma(i) = 0. |
|---|
| 1279 | END DO |
|---|
| 1280 | |
|---|
| 1281 | DO i = 1, icb - 1 |
|---|
| 1282 | ma(i) = 0 |
|---|
| 1283 | END DO |
|---|
| 1284 | |
|---|
| 1285 | |
|---|
| 1286 | |
|---|
| 1287 | ! cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc |
|---|
| 1288 | ! ICB REPRESENTE DE NIVEAU OU SE TROUVE LA |
|---|
| 1289 | ! BASE DU NUAGE , ET INB LE TOP DU NUAGE |
|---|
| 1290 | ! ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc |
|---|
| 1291 | |
|---|
| 1292 | DO i = 1, nd |
|---|
| 1293 | mke(i) = upwd(i) + dnwd(i) |
|---|
| 1294 | END DO |
|---|
| 1295 | |
|---|
| 1296 | |
|---|
| 1297 | ! *** Diagnose the in-cloud mixing ratio *** ! cld |
|---|
| 1298 | ! *** of condensed water *** ! cld |
|---|
| 1299 | ! cld |
|---|
| 1300 | DO i = 1, nd ! cld |
|---|
| 1301 | maa(i) = 0.0 ! cld |
|---|
| 1302 | wa(i) = 0.0 ! cld |
|---|
| 1303 | siga(i) = 0.0 ! cld |
|---|
| 1304 | END DO ! cld |
|---|
| 1305 | DO i = nk, inb ! cld |
|---|
| 1306 | DO k = i + 1, inb + 1 ! cld |
|---|
| 1307 | maa(i) = maa(i) + m(k) ! cld |
|---|
| 1308 | END DO ! cld |
|---|
| 1309 | END DO ! cld |
|---|
| 1310 | DO i = icb, inb - 1 ! cld |
|---|
| 1311 | axc(i) = 0. ! cld |
|---|
| 1312 | DO j = icb, i ! cld |
|---|
| 1313 | axc(i) = axc(i) + rd * (tvp(j) - tv(j)) * (ph(j) - ph(j + 1)) / p(j) ! cld |
|---|
| 1314 | END DO ! cld |
|---|
| 1315 | IF (axc(i)>0.0) THEN ! cld |
|---|
| 1316 | wa(i) = sqrt(2. * axc(i)) ! cld |
|---|
| 1317 | END IF ! cld |
|---|
| 1318 | END DO ! cld |
|---|
| 1319 | DO i = 1, nl ! cld |
|---|
| 1320 | IF (wa(i)>0.0) & ! cld |
|---|
| 1321 | siga(i) = maa(i) / wa(i) * rd * tvp(i) / p(i) / 100. / deltac ! cld |
|---|
| 1322 | siga(i) = min(siga(i), 1.0) ! cld |
|---|
| 1323 | qcondc(i) = siga(i) * clw(i) * (1. - ep(i)) & ! cld |
|---|
| 1324 | + (1. - siga(i)) * qcond(i) ! cld |
|---|
| 1325 | END DO ! cld |
|---|
| 1326 | |
|---|
| 1327 | |
|---|
| 1328 | ! @$$cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc |
|---|
| 1329 | ! @$$ CALL writeg1d(1,klev,ma,'ma ','ma ') |
|---|
| 1330 | ! @$$ CALL writeg1d(1,klev,upwd,'upwd ','upwd ') |
|---|
| 1331 | ! @$$ CALL writeg1d(1,klev,dnwd,'dnwd ','dnwd ') |
|---|
| 1332 | ! @$$ CALL writeg1d(1,klev,dnwd0,'dnwd0 ','dnwd0 ') |
|---|
| 1333 | ! @$$ CALL writeg1d(1,klev,tvp,'tvp ','tvp ') |
|---|
| 1334 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,3),'tra3 ','tra3 ') |
|---|
| 1335 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,4),'tra4 ','tra4 ') |
|---|
| 1336 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,5),'tra5 ','tra5 ') |
|---|
| 1337 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,6),'tra6 ','tra6 ') |
|---|
| 1338 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,7),'tra7 ','tra7 ') |
|---|
| 1339 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,8),'tra8 ','tra8 ') |
|---|
| 1340 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,9),'tra9 ','tra9 ') |
|---|
| 1341 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,10),'tra10','tra10') |
|---|
| 1342 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,11),'tra11','tra11') |
|---|
| 1343 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,12),'tra12','tra12') |
|---|
| 1344 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,13),'tra13','tra13') |
|---|
| 1345 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,14),'tra14','tra14') |
|---|
| 1346 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,15),'tra15','tra15') |
|---|
| 1347 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,16),'tra16','tra16') |
|---|
| 1348 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,17),'tra17','tra17') |
|---|
| 1349 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,18),'tra18','tra18') |
|---|
| 1350 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,19),'tra19','tra19') |
|---|
| 1351 | ! @$$ CALL writeg1d(1,klev,tra(1:klev,20),'tra20','tra20 ') |
|---|
| 1352 | ! @$$ CALL writeg1d(1,klev,trap(1:klev,1),'trp1','trp1') |
|---|
| 1353 | ! @$$ CALL writeg1d(1,klev,trap(1:klev,2),'trp2','trp2') |
|---|
| 1354 | ! @$$ CALL writeg1d(1,klev,trap(1:klev,3),'trp3','trp3') |
|---|
| 1355 | ! @$$ CALL writeg1d(1,klev,trap(1:klev,4),'trp4','trp4') |
|---|
| 1356 | ! @$$ CALL writeg1d(1,klev,trap(1:klev,5),'trp5','trp5') |
|---|
| 1357 | ! @$$ CALL writeg1d(1,klev,trap(1:klev,10),'trp10','trp10') |
|---|
| 1358 | ! @$$ CALL writeg1d(1,klev,trap(1:klev,12),'trp12','trp12') |
|---|
| 1359 | ! @$$ CALL writeg1d(1,klev,trap(1:klev,15),'trp15','trp15') |
|---|
| 1360 | ! @$$ CALL writeg1d(1,klev,trap(1:klev,20),'trp20','trp20') |
|---|
| 1361 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,1),'ftr1 ','ftr1 ') |
|---|
| 1362 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,2),'ftr2 ','ftr2 ') |
|---|
| 1363 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,3),'ftr3 ','ftr3 ') |
|---|
| 1364 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,4),'ftr4 ','ftr4 ') |
|---|
| 1365 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,5),'ftr5 ','ftr5 ') |
|---|
| 1366 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,6),'ftr6 ','ftr6 ') |
|---|
| 1367 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,7),'ftr7 ','ftr7 ') |
|---|
| 1368 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,8),'ftr8 ','ftr8 ') |
|---|
| 1369 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,9),'ftr9 ','ftr9 ') |
|---|
| 1370 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,10),'ftr10','ftr10') |
|---|
| 1371 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,11),'ftr11','ftr11') |
|---|
| 1372 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,12),'ftr12','ftr12') |
|---|
| 1373 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,13),'ftr13','ftr13') |
|---|
| 1374 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,14),'ftr14','ftr14') |
|---|
| 1375 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,15),'ftr15','ftr15') |
|---|
| 1376 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,16),'ftr16','ftr16') |
|---|
| 1377 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,17),'ftr17','ftr17') |
|---|
| 1378 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,18),'ftr18','ftr18') |
|---|
| 1379 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,19),'ftr19','ftr19') |
|---|
| 1380 | ! @$$ CALL writeg1d(1,klev,ftra(1:klev,20),'ftr20','ftr20 ') |
|---|
| 1381 | ! @$$ CALL writeg1d(1,klev,mp,'mp ','mp ') |
|---|
| 1382 | ! @$$ CALL writeg1d(1,klev,Mke,'Mke ','Mke ') |
|---|
| 1383 | |
|---|
| 1384 | |
|---|
| 1385 | |
|---|
| 1386 | ! cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc |
|---|
| 1387 | |
|---|
| 1388 | END SUBROUTINE convect3 |
|---|
| 1389 | ! --------------------------------------------------------------------------- |
|---|
