[13] | 1 | !********************************************************************** |
---|
| 2 | ! Copyright 1998,1999,2000,2001,2002,2005,2007,2008,2009,2010 * |
---|
| 3 | ! Andreas Stohl, Petra Seibert, A. Frank, Gerhard Wotawa, * |
---|
| 4 | ! Caroline Forster, Sabine Eckhardt, John Burkhart, Harald Sodemann * |
---|
| 5 | ! * |
---|
| 6 | ! This file is part of FLEXPART. * |
---|
| 7 | ! * |
---|
| 8 | ! FLEXPART is free software: you can redistribute it and/or modify * |
---|
| 9 | ! it under the terms of the GNU General Public License as published by* |
---|
| 10 | ! the Free Software Foundation, either version 3 of the License, or * |
---|
| 11 | ! (at your option) any later version. * |
---|
| 12 | ! * |
---|
| 13 | ! FLEXPART is distributed in the hope that it will be useful, * |
---|
| 14 | ! but WITHOUT ANY WARRANTY; without even the implied warranty of * |
---|
| 15 | ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * |
---|
| 16 | ! GNU General Public License for more details. * |
---|
| 17 | ! * |
---|
| 18 | ! You should have received a copy of the GNU General Public License * |
---|
| 19 | ! along with FLEXPART. If not, see <http://www.gnu.org/licenses/>. * |
---|
| 20 | !********************************************************************** |
---|
| 21 | |
---|
| 22 | !************************************************************************** |
---|
| 23 | !**** SUBROUTINE CONVECT ***** |
---|
| 24 | !**** VERSION 4.3c ***** |
---|
| 25 | !**** 20 May, 2002 ***** |
---|
| 26 | !**** Kerry Emanuel ***** |
---|
| 27 | !************************************************************************** |
---|
| 28 | ! |
---|
| 29 | SUBROUTINE CONVECT & |
---|
| 30 | (ND, NL, DELT, IFLAG, & |
---|
| 31 | PRECIP, WD, TPRIME, QPRIME, CBMF ) |
---|
| 32 | ! |
---|
| 33 | !-cv ************************************************************************* |
---|
| 34 | !-cv C. Forster, November 2003 - May 2004: |
---|
| 35 | !-cv |
---|
| 36 | !-cv The subroutine has been downloaded from Kerry Emanuel's homepage, |
---|
| 37 | !-cv where further infos on the convection scheme can be found |
---|
| 38 | !-cv http://www-paoc.mit.edu/~emanuel/home.html |
---|
| 39 | !-cv |
---|
| 40 | !-cv The following changes have been made to integrate this subroutine |
---|
| 41 | !-cv into FLEXPART |
---|
| 42 | !-cv |
---|
| 43 | !-cv Putting most of the variables in a new common block |
---|
| 44 | !-cv renaming eps to eps0 because there is some eps already in includepar |
---|
| 45 | !-cv |
---|
| 46 | !-cv removing the arrays U,V,TRA and related arrays |
---|
| 47 | !-cv |
---|
| 48 | !-cv renaming the original arrays T,Q,QS,P,PH to |
---|
| 49 | !-cv TCONV,QCONV,QSCONV,PCONV_HPA,PHCONV_HPA |
---|
| 50 | !-cv |
---|
| 51 | !-cv Initialization of variables has been put into parameter statements |
---|
| 52 | !-cv instead of assignment of values at each call, in order to save |
---|
| 53 | !-cv computation time. |
---|
| 54 | !*************************************************************************** |
---|
| 55 | ! |
---|
| 56 | !----------------------------------------------------------------------------- |
---|
| 57 | ! *** On input: *** |
---|
| 58 | ! |
---|
| 59 | !T: Array of absolute temperature (K) of dimension ND, with first |
---|
| 60 | ! index corresponding to lowest model level. Note that this array |
---|
| 61 | ! will be altered by the subroutine if dry convective adjustment |
---|
| 62 | ! occurs and if IPBL is not equal to 0. |
---|
| 63 | ! |
---|
| 64 | !Q: Array of specific humidity (gm/gm) of dimension ND, with first |
---|
| 65 | ! index corresponding to lowest model level. Must be defined |
---|
| 66 | ! at same grid levels as T. Note that this array will be altered |
---|
| 67 | ! if dry convective adjustment occurs and if IPBL is not equal to 0. |
---|
| 68 | ! |
---|
| 69 | !QS: Array of saturation specific humidity of dimension ND, with first |
---|
| 70 | ! index corresponding to lowest model level. Must be defined |
---|
| 71 | ! at same grid levels as T. Note that this array will be altered |
---|
| 72 | ! if dry convective adjustment occurs and if IPBL is not equal to 0. |
---|
| 73 | ! |
---|
| 74 | !U: Array of zonal wind velocity (m/s) of dimension ND, witth first |
---|
| 75 | ! index corresponding with the lowest model level. Defined at |
---|
| 76 | ! same levels as T. Note that this array will be altered if |
---|
| 77 | ! dry convective adjustment occurs and if IPBL is not equal to 0. |
---|
| 78 | ! |
---|
| 79 | !V: Same as U but for meridional velocity. |
---|
| 80 | ! |
---|
| 81 | !TRA: Array of passive tracer mixing ratio, of dimensions (ND,NTRA), |
---|
| 82 | ! where NTRA is the number of different tracers. If no |
---|
| 83 | ! convective tracer transport is needed, define a dummy |
---|
| 84 | ! input array of dimension (ND,1). Tracers are defined at |
---|
| 85 | ! same vertical levels as T. Note that this array will be altered |
---|
| 86 | ! if dry convective adjustment occurs and if IPBL is not equal to 0. |
---|
| 87 | ! |
---|
| 88 | !P: Array of pressure (mb) of dimension ND, with first |
---|
| 89 | ! index corresponding to lowest model level. Must be defined |
---|
| 90 | ! at same grid levels as T. |
---|
| 91 | ! |
---|
| 92 | !PH: Array of pressure (mb) of dimension ND+1, with first index |
---|
| 93 | ! corresponding to lowest level. These pressures are defined at |
---|
| 94 | ! levels intermediate between those of P, T, Q and QS. The first |
---|
| 95 | ! value of PH should be greater than (i.e. at a lower level than) |
---|
| 96 | ! the first value of the array P. |
---|
| 97 | ! |
---|
| 98 | !ND: The dimension of the arrays T,Q,QS,P,PH,FT and FQ |
---|
| 99 | ! |
---|
| 100 | !NL: The maximum number of levels to which convection can |
---|
| 101 | ! penetrate, plus 1. |
---|
| 102 | ! NL MUST be less than or equal to ND-1. |
---|
| 103 | ! |
---|
| 104 | !NTRA:The number of different tracers. If no tracer transport |
---|
| 105 | ! is needed, set this equal to 1. (On most compilers, setting |
---|
| 106 | ! NTRA to 0 will bypass tracer calculation, saving some CPU.) |
---|
| 107 | ! |
---|
| 108 | !DELT: The model time step (sec) between calls to CONVECT |
---|
| 109 | ! |
---|
| 110 | !---------------------------------------------------------------------------- |
---|
| 111 | ! *** On Output: *** |
---|
| 112 | ! |
---|
| 113 | !IFLAG: An output integer whose value denotes the following: |
---|
| 114 | ! |
---|
| 115 | ! VALUE INTERPRETATION |
---|
| 116 | ! ----- -------------- |
---|
| 117 | ! 0 No moist convection; atmosphere is not |
---|
| 118 | ! unstable, or surface temperature is less |
---|
| 119 | ! than 250 K or surface specific humidity |
---|
| 120 | ! is non-positive. |
---|
| 121 | ! |
---|
| 122 | ! 1 Moist convection occurs. |
---|
| 123 | ! |
---|
| 124 | ! 2 No moist convection: lifted condensation |
---|
| 125 | ! level is above the 200 mb level. |
---|
| 126 | ! |
---|
| 127 | ! 3 No moist convection: cloud base is higher |
---|
| 128 | ! then the level NL-1. |
---|
| 129 | ! |
---|
| 130 | ! 4 Moist convection occurs, but a CFL condition |
---|
| 131 | ! on the subsidence warming is violated. This |
---|
| 132 | ! does not cause the scheme to terminate. |
---|
| 133 | ! |
---|
| 134 | !FT: Array of temperature tendency (K/s) of dimension ND, defined at same |
---|
| 135 | ! grid levels as T, Q, QS and P. |
---|
| 136 | ! |
---|
| 137 | !FQ: Array of specific humidity tendencies ((gm/gm)/s) of dimension ND, |
---|
| 138 | ! defined at same grid levels as T, Q, QS and P. |
---|
| 139 | ! |
---|
| 140 | !FU: Array of forcing of zonal velocity (m/s^2) of dimension ND, |
---|
| 141 | ! defined at same grid levels as T. |
---|
| 142 | ! |
---|
| 143 | !FV: Same as FU, but for forcing of meridional velocity. |
---|
| 144 | ! |
---|
| 145 | !FTRA: Array of forcing of tracer content, in tracer mixing ratio per |
---|
| 146 | ! second, defined at same levels as T. Dimensioned (ND,NTRA). |
---|
| 147 | ! |
---|
| 148 | !PRECIP: Scalar convective precipitation rate (mm/day). |
---|
| 149 | ! |
---|
| 150 | !WD: A convective downdraft velocity scale. For use in surface |
---|
| 151 | ! flux parameterizations. See convect.ps file for details. |
---|
| 152 | ! |
---|
| 153 | !TPRIME: A convective downdraft temperature perturbation scale (K). |
---|
| 154 | ! For use in surface flux parameterizations. See convect.ps |
---|
| 155 | ! file for details. |
---|
| 156 | ! |
---|
| 157 | !QPRIME: A convective downdraft specific humidity |
---|
| 158 | ! perturbation scale (gm/gm). |
---|
| 159 | ! For use in surface flux parameterizations. See convect.ps |
---|
| 160 | ! file for details. |
---|
| 161 | ! |
---|
| 162 | !CBMF: The cloud base mass flux ((kg/m**2)/s). THIS SCALAR VALUE MUST |
---|
| 163 | ! BE STORED BY THE CALLING PROGRAM AND RETURNED TO CONVECT AT |
---|
| 164 | ! ITS NEXT CALL. That is, the value of CBMF must be "remembered" |
---|
| 165 | ! by the calling program between calls to CONVECT. |
---|
| 166 | ! |
---|
| 167 | !----------------------------------------------------------------------------- |
---|
| 168 | ! |
---|
| 169 | ! *** THE PARAMETER NA SHOULD IN GENERAL BE GREATER THAN *** |
---|
| 170 | ! *** OR EQUAL TO ND + 1 *** |
---|
| 171 | ! |
---|
| 172 | ! |
---|
| 173 | use par_mod |
---|
| 174 | use conv_mod |
---|
| 175 | |
---|
| 176 | implicit none |
---|
| 177 | ! |
---|
| 178 | !-cv====>Begin Module CONVECT File convect.f Undeclared variables |
---|
| 179 | ! |
---|
| 180 | !Argument variables |
---|
| 181 | ! |
---|
| 182 | integer :: iflag, nd, nl |
---|
| 183 | ! |
---|
| 184 | real :: cbmf, delt, precip, qprime, tprime, wd |
---|
| 185 | ! |
---|
| 186 | !Local variables |
---|
| 187 | ! |
---|
| 188 | integer :: i, icb, ihmin, inb, inb1, j, jtt, k |
---|
| 189 | integer :: nk |
---|
| 190 | ! |
---|
| 191 | real :: ad, afac, ahmax, ahmin, alt, altem |
---|
| 192 | real :: am, amp1, anum, asij, awat, b6, bf2, bsum, by |
---|
| 193 | real :: byp, c6, cape, capem, cbmfold, chi, coeff |
---|
| 194 | real :: cpinv, cwat, damps, dbo, dbosum |
---|
| 195 | real :: defrac, dei, delm, delp, delt0, delti, denom, dhdp |
---|
| 196 | real :: dpinv, dtma, dtmin, dtpbl, elacrit, ents |
---|
| 197 | real :: epmax, fac, fqold, frac, ftold |
---|
| 198 | real :: plcl, qp1, qsm, qstm, qti, rat |
---|
| 199 | real :: rdcp, revap, rh, scrit, sigt, sjmax |
---|
| 200 | real :: sjmin, smid, smin, stemp, tca |
---|
| 201 | real :: tvaplcl, tvpplcl, tvx, tvy, wdtrain |
---|
| 202 | |
---|
| 203 | !integer jc,jn |
---|
| 204 | !real alvnew,a2,ahm,alv,rm,sum,qnew,dphinv,tc,thbar,tnew,x |
---|
| 205 | |
---|
| 206 | real :: FUP(NA),FDOWN(NA) |
---|
| 207 | ! |
---|
| 208 | !-cv====>End Module CONVECT File convect.f |
---|
| 209 | |
---|
| 210 | INTEGER :: NENT(NA) |
---|
| 211 | REAL :: M(NA),MP(NA),MENT(NA,NA),QENT(NA,NA),ELIJ(NA,NA) |
---|
| 212 | REAL :: SIJ(NA,NA),TVP(NA),TV(NA),WATER(NA) |
---|
| 213 | REAL :: QP(NA),EP(NA),TH(NA),WT(NA),EVAP(NA),CLW(NA) |
---|
| 214 | REAL :: SIGP(NA),TP(NA),CPN(NA) |
---|
| 215 | REAL :: LV(NA),LVCP(NA),H(NA),HP(NA),GZ(NA),HM(NA) |
---|
| 216 | !REAL TOLD(NA) |
---|
| 217 | ! |
---|
| 218 | ! ----------------------------------------------------------------------- |
---|
| 219 | ! |
---|
| 220 | ! *** Specify Switches *** |
---|
| 221 | ! |
---|
| 222 | ! *** IPBL: Set to zero to bypass dry adiabatic adjustment *** |
---|
| 223 | ! *** Any other value results in dry adiabatic adjustment *** |
---|
| 224 | ! *** (Zero value recommended for use in models with *** |
---|
| 225 | ! *** boundary layer schemes) *** |
---|
| 226 | ! |
---|
| 227 | ! *** MINORIG: Lowest level from which convection may originate *** |
---|
| 228 | ! *** (Should be first model level at which T is defined *** |
---|
| 229 | ! *** for models using bulk PBL schemes; otherwise, it should *** |
---|
| 230 | ! *** be the first model level at which T is defined above *** |
---|
| 231 | ! *** the surface layer) *** |
---|
| 232 | ! |
---|
| 233 | INTEGER,PARAMETER :: IPBL=0 |
---|
| 234 | INTEGER,PARAMETER :: MINORIG=1 |
---|
| 235 | ! |
---|
| 236 | !------------------------------------------------------------------------------ |
---|
| 237 | ! |
---|
| 238 | ! *** SPECIFY PARAMETERS *** |
---|
| 239 | ! |
---|
| 240 | ! *** ELCRIT IS THE AUTOCONVERSION THERSHOLD WATER CONTENT (gm/gm) *** |
---|
| 241 | ! *** TLCRIT IS CRITICAL TEMPERATURE BELOW WHICH THE AUTO- *** |
---|
| 242 | ! *** CONVERSION THRESHOLD IS ASSUMED TO BE ZERO *** |
---|
| 243 | ! *** (THE AUTOCONVERSION THRESHOLD VARIES LINEARLY *** |
---|
| 244 | ! *** BETWEEN 0 C AND TLCRIT) *** |
---|
| 245 | ! *** ENTP IS THE COEFFICIENT OF MIXING IN THE ENTRAINMENT *** |
---|
| 246 | ! *** FORMULATION *** |
---|
| 247 | ! *** SIGD IS THE FRACTIONAL AREA COVERED BY UNSATURATED DNDRAFT *** |
---|
| 248 | ! *** SIGS IS THE FRACTION OF PRECIPITATION FALLING OUTSIDE *** |
---|
| 249 | ! *** OF CLOUD *** |
---|
| 250 | ! *** OMTRAIN IS THE ASSUMED FALL SPEED (P/s) OF RAIN *** |
---|
| 251 | ! *** OMTSNOW IS THE ASSUMED FALL SPEED (P/s) OF SNOW *** |
---|
| 252 | ! *** COEFFR IS A COEFFICIENT GOVERNING THE RATE OF EVAPORATION *** |
---|
| 253 | ! *** OF RAIN *** |
---|
| 254 | ! *** COEFFS IS A COEFFICIENT GOVERNING THE RATE OF EVAPORATION *** |
---|
| 255 | ! *** OF SNOW *** |
---|
| 256 | ! *** CU IS THE COEFFICIENT GOVERNING CONVECTIVE MOMENTUM *** |
---|
| 257 | ! *** TRANSPORT *** |
---|
| 258 | ! *** DTMAX IS THE MAXIMUM NEGATIVE TEMPERATURE PERTURBATION *** |
---|
| 259 | ! *** A LIFTED PARCEL IS ALLOWED TO HAVE BELOW ITS LFC *** |
---|
| 260 | ! *** ALPHA AND DAMP ARE PARAMETERS THAT CONTROL THE RATE OF *** |
---|
| 261 | ! *** APPROACH TO QUASI-EQUILIBRIUM *** |
---|
| 262 | ! *** (THEIR STANDARD VALUES ARE 0.20 AND 0.1, RESPECTIVELY) *** |
---|
| 263 | ! *** (DAMP MUST BE LESS THAN 1) *** |
---|
| 264 | ! |
---|
| 265 | REAL,PARAMETER :: ELCRIT=.0011 |
---|
| 266 | REAL,PARAMETER :: TLCRIT=-55.0 |
---|
| 267 | REAL,PARAMETER :: ENTP=1.5 |
---|
| 268 | REAL,PARAMETER :: SIGD=0.05 |
---|
| 269 | REAL,PARAMETER :: SIGS=0.12 |
---|
| 270 | REAL,PARAMETER :: OMTRAIN=50.0 |
---|
| 271 | REAL,PARAMETER :: OMTSNOW=5.5 |
---|
| 272 | REAL,PARAMETER :: COEFFR=1.0 |
---|
| 273 | REAL,PARAMETER :: COEFFS=0.8 |
---|
| 274 | REAL,PARAMETER :: CU=0.7 |
---|
| 275 | REAL,PARAMETER :: BETA=10.0 |
---|
| 276 | REAL,PARAMETER :: DTMAX=0.9 |
---|
| 277 | REAL,PARAMETER :: ALPHA=0.025 !original 0.2 |
---|
| 278 | REAL,PARAMETER :: DAMP=0.1 |
---|
| 279 | ! |
---|
| 280 | ! *** ASSIGN VALUES OF THERMODYNAMIC CONSTANTS, *** |
---|
| 281 | ! *** GRAVITY, AND LIQUID WATER DENSITY. *** |
---|
| 282 | ! *** THESE SHOULD BE CONSISTENT WITH *** |
---|
| 283 | ! *** THOSE USED IN CALLING PROGRAM *** |
---|
| 284 | ! *** NOTE: THESE ARE ALSO SPECIFIED IN SUBROUTINE TLIFT *** |
---|
| 285 | ! |
---|
| 286 | REAL,PARAMETER :: CPD=1005.7 |
---|
| 287 | REAL,PARAMETER :: CPV=1870.0 |
---|
| 288 | REAL,PARAMETER :: CL=2500.0 |
---|
| 289 | REAL,PARAMETER :: RV=461.5 |
---|
| 290 | REAL,PARAMETER :: RD=287.04 |
---|
| 291 | REAL,PARAMETER :: LV0=2.501E6 |
---|
| 292 | REAL,PARAMETER :: G=9.81 |
---|
| 293 | REAL,PARAMETER :: ROWL=1000.0 |
---|
| 294 | ! |
---|
| 295 | REAL,PARAMETER :: CPVMCL=CL-CPV |
---|
| 296 | REAL,PARAMETER :: EPS0=RD/RV |
---|
| 297 | REAL,PARAMETER :: EPSI=1./EPS0 |
---|
| 298 | REAL,PARAMETER :: GINV=1.0/G |
---|
| 299 | REAL,PARAMETER :: EPSILON=1.e-20 |
---|
| 300 | |
---|
| 301 | ! EPSILON IS A SMALL NUMBER USED TO EXCLUDE MASS FLUXES OF ZERO |
---|
| 302 | ! |
---|
| 303 | DELTI=1.0/DELT |
---|
| 304 | ! |
---|
| 305 | ! *** INITIALIZE OUTPUT ARRAYS AND PARAMETERS *** |
---|
| 306 | ! |
---|
| 307 | |
---|
| 308 | DO I=1,NL+1 |
---|
| 309 | FT(I)=0.0 |
---|
| 310 | FQ(I)=0.0 |
---|
| 311 | FDOWN(I)=0.0 |
---|
| 312 | SUB(I)=0.0 |
---|
| 313 | FUP(I)=0.0 |
---|
| 314 | M(I)=0.0 |
---|
| 315 | MP(I)=0.0 |
---|
| 316 | DO J=1,NL+1 |
---|
| 317 | FMASS(I,J)=0.0 |
---|
| 318 | MENT(I,J)=0.0 |
---|
| 319 | END DO |
---|
| 320 | END DO |
---|
| 321 | DO I=1,NL+1 |
---|
| 322 | RDCP=(RD*(1.-QCONV(I))+QCONV(I)*RV)/ & |
---|
| 323 | (CPD*(1.-QCONV(I))+QCONV(I)*CPV) |
---|
| 324 | TH(I)=TCONV(I)*(1000.0/PCONV_HPA(I))**RDCP |
---|
| 325 | END DO |
---|
| 326 | PRECIP=0.0 |
---|
| 327 | WD=0.0 |
---|
| 328 | TPRIME=0.0 |
---|
| 329 | QPRIME=0.0 |
---|
| 330 | IFLAG=0 |
---|
| 331 | ! |
---|
| 332 | ! IF(IPBL.NE.0)THEN |
---|
| 333 | ! |
---|
| 334 | !*** PERFORM DRY ADIABATIC ADJUSTMENT *** |
---|
| 335 | ! |
---|
| 336 | ! JC=0 |
---|
| 337 | ! DO 30 I=NL-1,1,-1 |
---|
| 338 | ! JN=0 |
---|
| 339 | ! SUM=TH(I)*(1.+QCONV(I)*EPSI-QCONV(I)) |
---|
| 340 | ! DO 10 J=I+1,NL |
---|
| 341 | ! SUM=SUM+TH(J)*(1.+QCONV(J)*EPSI-QCONV(J)) |
---|
| 342 | ! THBAR=SUM/REAL(J+1-I) |
---|
| 343 | ! IF((TH(J)*(1.+QCONV(J)*EPSI-QCONV(J))).LT.THBAR)JN=J |
---|
| 344 | ! 10 CONTINUE |
---|
| 345 | ! IF(I.EQ.1)JN=MAX(JN,2) |
---|
| 346 | ! IF(JN.EQ.0)GOTO 30 |
---|
| 347 | ! 12 CONTINUE |
---|
| 348 | ! AHM=0.0 |
---|
| 349 | ! RM=0.0 |
---|
| 350 | ! DO 15 J=I,JN |
---|
| 351 | ! AHM=AHM+(CPD*(1.-QCONV(J))+QCONV(J)*CPV)*TCONV(J)* |
---|
| 352 | ! + (PHCONV_HPA(J)-PHCONV_HPA(J+1)) |
---|
| 353 | ! RM=RM+QCONV(J)*(PHCONV_HPA(J)-PHCONV_HPA(J+1)) |
---|
| 354 | ! 15 CONTINUE |
---|
| 355 | ! DPHINV=1./(PHCONV_HPA(I)-PHCONV_HPA(JN+1)) |
---|
| 356 | ! RM=RM*DPHINV |
---|
| 357 | ! A2=0.0 |
---|
| 358 | ! DO 20 J=I,JN |
---|
| 359 | ! QCONV(J)=RM |
---|
| 360 | ! RDCP=(RD*(1.-QCONV(J))+QCONV(J)*RV)/ |
---|
| 361 | ! 1 (CPD*(1.-QCONV(J))+QCONV(J)*CPV) |
---|
| 362 | ! X=(0.001*PCONV_HPA(J))**RDCP |
---|
| 363 | ! TOLD(J)=TCONV(J) |
---|
| 364 | ! TCONV(J)=X |
---|
| 365 | ! A2=A2+(CPD*(1.-QCONV(J))+QCONV(J)*CPV)*X* |
---|
| 366 | ! 1 (PHCONV_HPA(J)-PHCONV_HPA(J+1)) |
---|
| 367 | ! 20 CONTINUE |
---|
| 368 | ! DO 25 J=I,JN |
---|
| 369 | ! TH(J)=AHM/A2 |
---|
| 370 | ! TCONV(J)=TCONV(J)*TH(J) |
---|
| 371 | ! TC=TOLD(J)-273.15 |
---|
| 372 | ! ALV=LV0-CPVMCL*TC |
---|
| 373 | ! QSCONV(J)=QSCONV(J)+QSCONV(J)*(1.+QSCONV(J)*(EPSI-1.))*ALV* |
---|
| 374 | ! 1 (TCONV(J)- TOLD(J))/(RV*TOLD(J)*TOLD(J)) |
---|
| 375 | ! if (qslev(j) .lt. 0.) then |
---|
| 376 | ! write(*,*) 'qslev.lt.0 ',j,qslev |
---|
| 377 | ! endif |
---|
| 378 | ! 25 CONTINUE |
---|
| 379 | ! IF((TH(JN+1)*(1.+QCONV(JN+1)*EPSI-QCONV(JN+1))).LT. |
---|
| 380 | ! 1 (TH(JN)*(1.+QCONV(JN)*EPSI-QCONV(JN))))THEN |
---|
| 381 | ! JN=JN+1 |
---|
| 382 | ! GOTO 12 |
---|
| 383 | ! END IF |
---|
| 384 | ! IF(I.EQ.1)JC=JN |
---|
| 385 | ! 30 CONTINUE |
---|
| 386 | ! |
---|
| 387 | ! *** Remove any supersaturation that results from adjustment *** |
---|
| 388 | ! |
---|
| 389 | !IF(JC.GT.1)THEN |
---|
| 390 | ! DO 38 J=1,JC |
---|
| 391 | ! IF(QSCONV(J).LT.QCONV(J))THEN |
---|
| 392 | ! ALV=LV0-CPVMCL*(TCONV(J)-273.15) |
---|
| 393 | ! TNEW=TCONV(J)+ALV*(QCONV(J)-QSCONV(J))/(CPD*(1.-QCONV(J))+ |
---|
| 394 | ! 1 CL*QCONV(J)+QSCONV(J)*(CPV-CL+ALV*ALV/(RV*TCONV(J)*TCONV(J)))) |
---|
| 395 | ! ALVNEW=LV0-CPVMCL*(TNEW-273.15) |
---|
| 396 | ! QNEW=(ALV*QCONV(J)-(TNEW-TCONV(J))*(CPD*(1.-QCONV(J)) |
---|
| 397 | ! 1 +CL*QCONV(J)))/ALVNEW |
---|
| 398 | ! PRECIP=PRECIP+24.*3600.*1.0E5*(PHCONV_HPA(J)-PHCONV_HPA(J+1))* |
---|
| 399 | ! 1 (QCONV(J)-QNEW)/(G*DELT*ROWL) |
---|
| 400 | ! TCONV(J)=TNEW |
---|
| 401 | ! QCONV(J)=QNEW |
---|
| 402 | ! QSCONV(J)=QNEW |
---|
| 403 | ! END IF |
---|
| 404 | ! 38 CONTINUE |
---|
| 405 | !END IF |
---|
| 406 | ! |
---|
| 407 | !END IF |
---|
| 408 | ! |
---|
| 409 | ! *** CALCULATE ARRAYS OF GEOPOTENTIAL, HEAT CAPACITY AND STATIC ENERGY |
---|
| 410 | ! |
---|
| 411 | GZ(1)=0.0 |
---|
| 412 | CPN(1)=CPD*(1.-QCONV(1))+QCONV(1)*CPV |
---|
| 413 | H(1)=TCONV(1)*CPN(1) |
---|
| 414 | LV(1)=LV0-CPVMCL*(TCONV(1)-273.15) |
---|
| 415 | HM(1)=LV(1)*QCONV(1) |
---|
| 416 | TV(1)=TCONV(1)*(1.+QCONV(1)*EPSI-QCONV(1)) |
---|
| 417 | AHMIN=1.0E12 |
---|
| 418 | IHMIN=NL |
---|
| 419 | DO I=2,NL+1 |
---|
| 420 | TVX=TCONV(I)*(1.+QCONV(I)*EPSI-QCONV(I)) |
---|
| 421 | TVY=TCONV(I-1)*(1.+QCONV(I-1)*EPSI-QCONV(I-1)) |
---|
| 422 | GZ(I)=GZ(I-1)+0.5*RD*(TVX+TVY)*(PCONV_HPA(I-1)-PCONV_HPA(I))/ & |
---|
| 423 | PHCONV_HPA(I) |
---|
| 424 | CPN(I)=CPD*(1.-QCONV(I))+CPV*QCONV(I) |
---|
| 425 | H(I)=TCONV(I)*CPN(I)+GZ(I) |
---|
| 426 | LV(I)=LV0-CPVMCL*(TCONV(I)-273.15) |
---|
| 427 | HM(I)=(CPD*(1.-QCONV(I))+CL*QCONV(I))*(TCONV(I)-TCONV(1))+ & |
---|
| 428 | LV(I)*QCONV(I)+GZ(I) |
---|
| 429 | TV(I)=TCONV(I)*(1.+QCONV(I)*EPSI-QCONV(I)) |
---|
| 430 | ! |
---|
| 431 | ! *** Find level of minimum moist static energy *** |
---|
| 432 | ! |
---|
| 433 | IF(I.GE.MINORIG.AND.HM(I).LT.AHMIN.AND.HM(I).LT.HM(I-1))THEN |
---|
| 434 | AHMIN=HM(I) |
---|
| 435 | IHMIN=I |
---|
| 436 | END IF |
---|
| 437 | END DO |
---|
| 438 | IHMIN=MIN(IHMIN, NL-1) |
---|
| 439 | ! |
---|
| 440 | ! *** Find that model level below the level of minimum moist *** |
---|
| 441 | ! *** static energy that has the maximum value of moist static energy *** |
---|
| 442 | ! |
---|
| 443 | AHMAX=0.0 |
---|
| 444 | ! *** bug fixed: need to assign an initial value to NK |
---|
| 445 | ! HSO, 05.08.2009 |
---|
| 446 | NK=MINORIG |
---|
| 447 | DO I=MINORIG,IHMIN |
---|
| 448 | IF(HM(I).GT.AHMAX)THEN |
---|
| 449 | NK=I |
---|
| 450 | AHMAX=HM(I) |
---|
| 451 | END IF |
---|
| 452 | END DO |
---|
| 453 | ! |
---|
| 454 | ! *** CHECK WHETHER PARCEL LEVEL TEMPERATURE AND SPECIFIC HUMIDITY *** |
---|
| 455 | ! *** ARE REASONABLE *** |
---|
| 456 | ! *** Skip convection if HM increases monotonically upward *** |
---|
| 457 | ! |
---|
| 458 | IF(TCONV(NK).LT.250.0.OR.QCONV(NK).LE.0.0.OR.IHMIN.EQ.(NL-1)) & |
---|
| 459 | THEN |
---|
| 460 | IFLAG=0 |
---|
| 461 | CBMF=0.0 |
---|
| 462 | RETURN |
---|
| 463 | END IF |
---|
| 464 | ! |
---|
| 465 | ! *** CALCULATE LIFTED CONDENSATION LEVEL OF AIR AT PARCEL ORIGIN LEVEL *** |
---|
| 466 | ! *** (WITHIN 0.2% OF FORMULA OF BOLTON, MON. WEA. REV.,1980) *** |
---|
| 467 | ! |
---|
| 468 | RH=QCONV(NK)/QSCONV(NK) |
---|
| 469 | CHI=TCONV(NK)/(1669.0-122.0*RH-TCONV(NK)) |
---|
| 470 | PLCL=PCONV_HPA(NK)*(RH**CHI) |
---|
| 471 | IF(PLCL.LT.200.0.OR.PLCL.GE.2000.0)THEN |
---|
| 472 | IFLAG=2 |
---|
| 473 | CBMF=0.0 |
---|
| 474 | RETURN |
---|
| 475 | END IF |
---|
| 476 | ! |
---|
| 477 | ! *** CALCULATE FIRST LEVEL ABOVE LCL (=ICB) *** |
---|
| 478 | ! |
---|
| 479 | ICB=NL-1 |
---|
| 480 | DO I=NK+1,NL |
---|
| 481 | IF(PCONV_HPA(I).LT.PLCL)THEN |
---|
| 482 | ICB=MIN(ICB,I) |
---|
| 483 | END IF |
---|
| 484 | END DO |
---|
| 485 | IF(ICB.GE.(NL-1))THEN |
---|
| 486 | IFLAG=3 |
---|
| 487 | CBMF=0.0 |
---|
| 488 | RETURN |
---|
| 489 | END IF |
---|
| 490 | ! |
---|
| 491 | ! *** FIND TEMPERATURE UP THROUGH ICB AND TEST FOR INSTABILITY *** |
---|
| 492 | ! |
---|
| 493 | ! *** SUBROUTINE TLIFT CALCULATES PART OF THE LIFTED PARCEL VIRTUAL *** |
---|
| 494 | ! *** TEMPERATURE, THE ACTUAL TEMPERATURE AND THE ADIABATIC *** |
---|
| 495 | ! *** LIQUID WATER CONTENT *** |
---|
| 496 | ! |
---|
| 497 | CALL TLIFT(GZ,ICB,NK,TVP,TP,CLW,ND,NL,1) |
---|
| 498 | DO I=NK,ICB |
---|
| 499 | TVP(I)=TVP(I)-TP(I)*QCONV(NK) |
---|
| 500 | END DO |
---|
| 501 | ! |
---|
| 502 | ! *** If there was no convection at last time step and parcel *** |
---|
| 503 | ! *** is stable at ICB then skip rest of calculation *** |
---|
| 504 | ! |
---|
| 505 | IF(CBMF.EQ.0.0.AND.TVP(ICB).LE.(TV(ICB)-DTMAX))THEN |
---|
| 506 | IFLAG=0 |
---|
| 507 | RETURN |
---|
| 508 | END IF |
---|
| 509 | ! |
---|
| 510 | ! *** IF THIS POINT IS REACHED, MOIST CONVECTIVE ADJUSTMENT IS NECESSARY *** |
---|
| 511 | ! |
---|
| 512 | IF(IFLAG.NE.4)IFLAG=1 |
---|
| 513 | ! |
---|
| 514 | ! *** FIND THE REST OF THE LIFTED PARCEL TEMPERATURES *** |
---|
| 515 | ! |
---|
| 516 | CALL TLIFT(GZ,ICB,NK,TVP,TP,CLW,ND,NL,2) |
---|
| 517 | ! |
---|
| 518 | ! *** SET THE PRECIPITATION EFFICIENCIES AND THE FRACTION OF *** |
---|
| 519 | ! *** PRECIPITATION FALLING OUTSIDE OF CLOUD *** |
---|
| 520 | ! *** THESE MAY BE FUNCTIONS OF TP(I), PCONV_HPA(I) AND CLW(I) *** |
---|
| 521 | ! |
---|
| 522 | DO I=1,NK |
---|
| 523 | EP(I)=0.0 |
---|
| 524 | SIGP(I)=SIGS |
---|
| 525 | END DO |
---|
| 526 | DO I=NK+1,NL |
---|
| 527 | TCA=TP(I)-273.15 |
---|
| 528 | IF(TCA.GE.0.0)THEN |
---|
| 529 | ELACRIT=ELCRIT |
---|
| 530 | ELSE |
---|
| 531 | ELACRIT=ELCRIT*(1.0-TCA/TLCRIT) |
---|
| 532 | END IF |
---|
| 533 | ELACRIT=MAX(ELACRIT,0.0) |
---|
| 534 | EPMAX=0.999 |
---|
| 535 | EP(I)=EPMAX*(1.0-ELACRIT/MAX(CLW(I),1.0E-8)) |
---|
| 536 | EP(I)=MAX(EP(I),0.0) |
---|
| 537 | EP(I)=MIN(EP(I),EPMAX) |
---|
| 538 | SIGP(I)=SIGS |
---|
| 539 | END DO |
---|
| 540 | ! |
---|
| 541 | ! *** CALCULATE VIRTUAL TEMPERATURE AND LIFTED PARCEL *** |
---|
| 542 | ! *** VIRTUAL TEMPERATURE *** |
---|
| 543 | ! |
---|
| 544 | DO I=ICB+1,NL |
---|
| 545 | TVP(I)=TVP(I)-TP(I)*QCONV(NK) |
---|
| 546 | END DO |
---|
| 547 | TVP(NL+1)=TVP(NL)-(GZ(NL+1)-GZ(NL))/CPD |
---|
| 548 | ! |
---|
| 549 | ! *** NOW INITIALIZE VARIOUS ARRAYS USED IN THE COMPUTATIONS *** |
---|
| 550 | ! |
---|
| 551 | DO I=1,NL+1 |
---|
| 552 | HP(I)=H(I) |
---|
| 553 | NENT(I)=0 |
---|
| 554 | WATER(I)=0.0 |
---|
| 555 | EVAP(I)=0.0 |
---|
| 556 | WT(I)=OMTSNOW |
---|
| 557 | LVCP(I)=LV(I)/CPN(I) |
---|
| 558 | DO J=1,NL+1 |
---|
| 559 | QENT(I,J)=QCONV(J) |
---|
| 560 | ELIJ(I,J)=0.0 |
---|
| 561 | SIJ(I,J)=0.0 |
---|
| 562 | END DO |
---|
| 563 | END DO |
---|
| 564 | QP(1)=QCONV(1) |
---|
| 565 | DO I=2,NL+1 |
---|
| 566 | QP(I)=QCONV(I-1) |
---|
| 567 | END DO |
---|
| 568 | ! |
---|
| 569 | ! *** FIND THE FIRST MODEL LEVEL (INB1) ABOVE THE PARCEL'S *** |
---|
| 570 | ! *** HIGHEST LEVEL OF NEUTRAL BUOYANCY *** |
---|
| 571 | ! *** AND THE HIGHEST LEVEL OF POSITIVE CAPE (INB) *** |
---|
| 572 | ! |
---|
| 573 | CAPE=0.0 |
---|
| 574 | CAPEM=0.0 |
---|
| 575 | INB=ICB+1 |
---|
| 576 | INB1=INB |
---|
| 577 | BYP=0.0 |
---|
| 578 | DO I=ICB+1,NL-1 |
---|
| 579 | BY=(TVP(I)-TV(I))*(PHCONV_HPA(I)-PHCONV_HPA(I+1))/PCONV_HPA(I) |
---|
| 580 | CAPE=CAPE+BY |
---|
| 581 | IF(BY.GE.0.0)INB1=I+1 |
---|
| 582 | IF(CAPE.GT.0.0)THEN |
---|
| 583 | INB=I+1 |
---|
| 584 | BYP=(TVP(I+1)-TV(I+1))*(PHCONV_HPA(I+1)-PHCONV_HPA(I+2))/ & |
---|
| 585 | PCONV_HPA(I+1) |
---|
| 586 | CAPEM=CAPE |
---|
| 587 | END IF |
---|
| 588 | END DO |
---|
| 589 | INB=MAX(INB,INB1) |
---|
| 590 | CAPE=CAPEM+BYP |
---|
| 591 | DEFRAC=CAPEM-CAPE |
---|
| 592 | DEFRAC=MAX(DEFRAC,0.001) |
---|
| 593 | FRAC=-CAPE/DEFRAC |
---|
| 594 | FRAC=MIN(FRAC,1.0) |
---|
| 595 | FRAC=MAX(FRAC,0.0) |
---|
| 596 | ! |
---|
| 597 | ! *** CALCULATE LIQUID WATER STATIC ENERGY OF LIFTED PARCEL *** |
---|
| 598 | ! |
---|
| 599 | DO I=ICB,INB |
---|
| 600 | HP(I)=H(NK)+(LV(I)+(CPD-CPV)*TCONV(I))*EP(I)*CLW(I) |
---|
| 601 | END DO |
---|
| 602 | ! |
---|
| 603 | ! *** CALCULATE CLOUD BASE MASS FLUX AND RATES OF MIXING, M(I), *** |
---|
| 604 | ! *** AT EACH MODEL LEVEL *** |
---|
| 605 | ! |
---|
| 606 | DBOSUM=0.0 |
---|
| 607 | ! |
---|
| 608 | ! *** INTERPOLATE DIFFERENCE BETWEEN LIFTED PARCEL AND *** |
---|
| 609 | ! *** ENVIRONMENTAL TEMPERATURES TO LIFTED CONDENSATION LEVEL *** |
---|
| 610 | ! |
---|
| 611 | TVPPLCL=TVP(ICB-1)-RD*TVP(ICB-1)*(PCONV_HPA(ICB-1)-PLCL)/ & |
---|
| 612 | (CPN(ICB-1)*PCONV_HPA(ICB-1)) |
---|
| 613 | TVAPLCL=TV(ICB)+(TVP(ICB)-TVP(ICB+1))*(PLCL-PCONV_HPA(ICB))/ & |
---|
| 614 | (PCONV_HPA(ICB)-PCONV_HPA(ICB+1)) |
---|
| 615 | DTPBL=0.0 |
---|
| 616 | DO I=NK,ICB-1 |
---|
| 617 | DTPBL=DTPBL+(TVP(I)-TV(I))*(PHCONV_HPA(I)-PHCONV_HPA(I+1)) |
---|
| 618 | END DO |
---|
| 619 | DTPBL=DTPBL/(PHCONV_HPA(NK)-PHCONV_HPA(ICB)) |
---|
| 620 | DTMIN=TVPPLCL-TVAPLCL+DTMAX+DTPBL |
---|
| 621 | DTMA=DTMIN |
---|
| 622 | ! |
---|
| 623 | ! *** ADJUST CLOUD BASE MASS FLUX *** |
---|
| 624 | ! |
---|
| 625 | CBMFOLD=CBMF |
---|
| 626 | ! *** C. Forster: adjustment of CBMF is not allowed to depend on FLEXPART timestep |
---|
| 627 | DELT0=DELT/3. |
---|
| 628 | DAMPS=DAMP*DELT/DELT0 |
---|
| 629 | CBMF=(1.-DAMPS)*CBMF+0.1*ALPHA*DTMA |
---|
| 630 | CBMF=MAX(CBMF,0.0) |
---|
| 631 | ! |
---|
| 632 | ! *** If cloud base mass flux is zero, skip rest of calculation *** |
---|
| 633 | ! |
---|
| 634 | IF(CBMF.EQ.0.0.AND.CBMFOLD.EQ.0.0)THEN |
---|
| 635 | RETURN |
---|
| 636 | END IF |
---|
| 637 | |
---|
| 638 | ! |
---|
| 639 | ! *** CALCULATE RATES OF MIXING, M(I) *** |
---|
| 640 | ! |
---|
| 641 | M(ICB)=0.0 |
---|
| 642 | DO I=ICB+1,INB |
---|
| 643 | K=MIN(I,INB1) |
---|
| 644 | DBO=ABS(TV(K)-TVP(K))+ & |
---|
| 645 | ENTP*0.02*(PHCONV_HPA(K)-PHCONV_HPA(K+1)) |
---|
| 646 | DBOSUM=DBOSUM+DBO |
---|
| 647 | M(I)=CBMF*DBO |
---|
| 648 | END DO |
---|
| 649 | DO I=ICB+1,INB |
---|
| 650 | M(I)=M(I)/DBOSUM |
---|
| 651 | END DO |
---|
| 652 | ! |
---|
| 653 | ! *** CALCULATE ENTRAINED AIR MASS FLUX (MENT), TOTAL WATER MIXING *** |
---|
| 654 | ! *** RATIO (QENT), TOTAL CONDENSED WATER (ELIJ), AND MIXING *** |
---|
| 655 | ! *** FRACTION (SIJ) *** |
---|
| 656 | ! |
---|
| 657 | DO I=ICB+1,INB |
---|
| 658 | QTI=QCONV(NK)-EP(I)*CLW(I) |
---|
| 659 | DO J=ICB,INB |
---|
| 660 | BF2=1.+LV(J)*LV(J)*QSCONV(J)/(RV*TCONV(J)*TCONV(J)*CPD) |
---|
| 661 | ANUM=H(J)-HP(I)+(CPV-CPD)*TCONV(J)*(QTI-QCONV(J)) |
---|
| 662 | DENOM=H(I)-HP(I)+(CPD-CPV)*(QCONV(I)-QTI)*TCONV(J) |
---|
| 663 | DEI=DENOM |
---|
| 664 | IF(ABS(DEI).LT.0.01)DEI=0.01 |
---|
| 665 | SIJ(I,J)=ANUM/DEI |
---|
| 666 | SIJ(I,I)=1.0 |
---|
| 667 | ALTEM=SIJ(I,J)*QCONV(I)+(1.-SIJ(I,J))*QTI-QSCONV(J) |
---|
| 668 | ALTEM=ALTEM/BF2 |
---|
| 669 | CWAT=CLW(J)*(1.-EP(J)) |
---|
| 670 | STEMP=SIJ(I,J) |
---|
| 671 | IF((STEMP.LT.0.0.OR.STEMP.GT.1.0.OR. & |
---|
| 672 | ALTEM.GT.CWAT).AND.J.GT.I)THEN |
---|
| 673 | ANUM=ANUM-LV(J)*(QTI-QSCONV(J)-CWAT*BF2) |
---|
| 674 | DENOM=DENOM+LV(J)*(QCONV(I)-QTI) |
---|
| 675 | IF(ABS(DENOM).LT.0.01)DENOM=0.01 |
---|
| 676 | SIJ(I,J)=ANUM/DENOM |
---|
| 677 | ALTEM=SIJ(I,J)*QCONV(I)+(1.-SIJ(I,J))*QTI-QSCONV(J) |
---|
| 678 | ALTEM=ALTEM-(BF2-1.)*CWAT |
---|
| 679 | END IF |
---|
| 680 | IF(SIJ(I,J).GT.0.0.AND.SIJ(I,J).LT.0.9)THEN |
---|
| 681 | QENT(I,J)=SIJ(I,J)*QCONV(I)+(1.-SIJ(I,J))*QTI |
---|
| 682 | ELIJ(I,J)=ALTEM |
---|
| 683 | ELIJ(I,J)=MAX(0.0,ELIJ(I,J)) |
---|
| 684 | MENT(I,J)=M(I)/(1.-SIJ(I,J)) |
---|
| 685 | NENT(I)=NENT(I)+1 |
---|
| 686 | END IF |
---|
| 687 | SIJ(I,J)=MAX(0.0,SIJ(I,J)) |
---|
| 688 | SIJ(I,J)=MIN(1.0,SIJ(I,J)) |
---|
| 689 | END DO |
---|
| 690 | ! |
---|
| 691 | ! *** IF NO AIR CAN ENTRAIN AT LEVEL I ASSUME THAT UPDRAFT DETRAINS *** |
---|
| 692 | ! *** AT THAT LEVEL AND CALCULATE DETRAINED AIR FLUX AND PROPERTIES *** |
---|
| 693 | ! |
---|
| 694 | IF(NENT(I).EQ.0)THEN |
---|
| 695 | MENT(I,I)=M(I) |
---|
| 696 | QENT(I,I)=QCONV(NK)-EP(I)*CLW(I) |
---|
| 697 | ELIJ(I,I)=CLW(I) |
---|
| 698 | SIJ(I,I)=1.0 |
---|
| 699 | END IF |
---|
| 700 | END DO |
---|
| 701 | SIJ(INB,INB)=1.0 |
---|
| 702 | ! |
---|
| 703 | ! *** NORMALIZE ENTRAINED AIR MASS FLUXES TO REPRESENT EQUAL *** |
---|
| 704 | ! *** PROBABILITIES OF MIXING *** |
---|
| 705 | ! |
---|
| 706 | DO I=ICB+1,INB |
---|
| 707 | IF(NENT(I).NE.0)THEN |
---|
| 708 | QP1=QCONV(NK)-EP(I)*CLW(I) |
---|
| 709 | ANUM=H(I)-HP(I)-LV(I)*(QP1-QSCONV(I)) |
---|
| 710 | DENOM=H(I)-HP(I)+LV(I)*(QCONV(I)-QP1) |
---|
| 711 | IF(ABS(DENOM).LT.0.01)DENOM=0.01 |
---|
| 712 | SCRIT=ANUM/DENOM |
---|
| 713 | ALT=QP1-QSCONV(I)+SCRIT*(QCONV(I)-QP1) |
---|
| 714 | IF(ALT.LT.0.0)SCRIT=1.0 |
---|
| 715 | SCRIT=MAX(SCRIT,0.0) |
---|
| 716 | ASIJ=0.0 |
---|
| 717 | SMIN=1.0 |
---|
| 718 | DO J=ICB,INB |
---|
| 719 | IF(SIJ(I,J).GT.0.0.AND.SIJ(I,J).LT.0.9)THEN |
---|
| 720 | IF(J.GT.I)THEN |
---|
| 721 | SMID=MIN(SIJ(I,J),SCRIT) |
---|
| 722 | SJMAX=SMID |
---|
| 723 | SJMIN=SMID |
---|
| 724 | IF(SMID.LT.SMIN.AND.SIJ(I,J+1).LT.SMID)THEN |
---|
| 725 | SMIN=SMID |
---|
| 726 | SJMAX=MIN(SIJ(I,J+1),SIJ(I,J),SCRIT) |
---|
| 727 | SJMIN=MAX(SIJ(I,J-1),SIJ(I,J)) |
---|
| 728 | SJMIN=MIN(SJMIN,SCRIT) |
---|
| 729 | END IF |
---|
| 730 | ELSE |
---|
| 731 | SJMAX=MAX(SIJ(I,J+1),SCRIT) |
---|
| 732 | SMID=MAX(SIJ(I,J),SCRIT) |
---|
| 733 | SJMIN=0.0 |
---|
| 734 | IF(J.GT.1)SJMIN=SIJ(I,J-1) |
---|
| 735 | SJMIN=MAX(SJMIN,SCRIT) |
---|
| 736 | END IF |
---|
| 737 | DELP=ABS(SJMAX-SMID) |
---|
| 738 | DELM=ABS(SJMIN-SMID) |
---|
| 739 | ASIJ=ASIJ+(DELP+DELM)*(PHCONV_HPA(J)-PHCONV_HPA(J+1)) |
---|
| 740 | MENT(I,J)=MENT(I,J)*(DELP+DELM)* & |
---|
| 741 | (PHCONV_HPA(J)-PHCONV_HPA(J+1)) |
---|
| 742 | END IF |
---|
| 743 | END DO |
---|
| 744 | ASIJ=MAX(1.0E-21,ASIJ) |
---|
| 745 | ASIJ=1.0/ASIJ |
---|
| 746 | DO J=ICB,INB |
---|
| 747 | MENT(I,J)=MENT(I,J)*ASIJ |
---|
| 748 | END DO |
---|
| 749 | BSUM=0.0 |
---|
| 750 | DO J=ICB,INB |
---|
| 751 | BSUM=BSUM+MENT(I,J) |
---|
| 752 | END DO |
---|
| 753 | IF(BSUM.LT.1.0E-18)THEN |
---|
| 754 | NENT(I)=0 |
---|
| 755 | MENT(I,I)=M(I) |
---|
| 756 | QENT(I,I)=QCONV(NK)-EP(I)*CLW(I) |
---|
| 757 | ELIJ(I,I)=CLW(I) |
---|
| 758 | SIJ(I,I)=1.0 |
---|
| 759 | END IF |
---|
| 760 | END IF |
---|
| 761 | END DO |
---|
| 762 | ! |
---|
| 763 | ! *** CHECK WHETHER EP(INB)=0, IF SO, SKIP PRECIPITATING *** |
---|
| 764 | ! *** DOWNDRAFT CALCULATION *** |
---|
| 765 | ! |
---|
| 766 | IF(EP(INB).LT.0.0001)GOTO 405 |
---|
| 767 | ! |
---|
| 768 | ! *** INTEGRATE LIQUID WATER EQUATION TO FIND CONDENSED WATER *** |
---|
| 769 | ! *** AND CONDENSED WATER FLUX *** |
---|
| 770 | ! |
---|
| 771 | JTT=2 |
---|
| 772 | ! |
---|
| 773 | ! *** BEGIN DOWNDRAFT LOOP *** |
---|
| 774 | ! |
---|
| 775 | DO I=INB,1,-1 |
---|
| 776 | ! |
---|
| 777 | ! *** CALCULATE DETRAINED PRECIPITATION *** |
---|
| 778 | ! |
---|
| 779 | WDTRAIN=G*EP(I)*M(I)*CLW(I) |
---|
| 780 | IF(I.GT.1)THEN |
---|
| 781 | DO J=1,I-1 |
---|
| 782 | AWAT=ELIJ(J,I)-(1.-EP(I))*CLW(I) |
---|
| 783 | AWAT=MAX(0.0,AWAT) |
---|
| 784 | WDTRAIN=WDTRAIN+G*AWAT*MENT(J,I) |
---|
| 785 | END DO |
---|
| 786 | END IF |
---|
| 787 | ! |
---|
| 788 | ! *** FIND RAIN WATER AND EVAPORATION USING PROVISIONAL *** |
---|
| 789 | ! *** ESTIMATES OF QP(I)AND QP(I-1) *** |
---|
| 790 | ! |
---|
| 791 | ! |
---|
| 792 | ! *** Value of terminal velocity and coefficient of evaporation for snow *** |
---|
| 793 | ! |
---|
| 794 | COEFF=COEFFS |
---|
| 795 | WT(I)=OMTSNOW |
---|
| 796 | ! |
---|
| 797 | ! *** Value of terminal velocity and coefficient of evaporation for rain *** |
---|
| 798 | ! |
---|
| 799 | IF(TCONV(I).GT.273.0)THEN |
---|
| 800 | COEFF=COEFFR |
---|
| 801 | WT(I)=OMTRAIN |
---|
| 802 | END IF |
---|
| 803 | QSM=0.5*(QCONV(I)+QP(I+1)) |
---|
| 804 | AFAC=COEFF*PHCONV_HPA(I)*(QSCONV(I)-QSM)/ & |
---|
| 805 | (1.0E4+2.0E3*PHCONV_HPA(I)*QSCONV(I)) |
---|
| 806 | AFAC=MAX(AFAC,0.0) |
---|
| 807 | SIGT=SIGP(I) |
---|
| 808 | SIGT=MAX(0.0,SIGT) |
---|
| 809 | SIGT=MIN(1.0,SIGT) |
---|
| 810 | B6=100.*(PHCONV_HPA(I)-PHCONV_HPA(I+1))*SIGT*AFAC/WT(I) |
---|
| 811 | C6=(WATER(I+1)*WT(I+1)+WDTRAIN/SIGD)/WT(I) |
---|
| 812 | REVAP=0.5*(-B6+SQRT(B6*B6+4.*C6)) |
---|
| 813 | EVAP(I)=SIGT*AFAC*REVAP |
---|
| 814 | WATER(I)=REVAP*REVAP |
---|
| 815 | ! |
---|
| 816 | ! *** CALCULATE PRECIPITATING DOWNDRAFT MASS FLUX UNDER *** |
---|
| 817 | ! *** HYDROSTATIC APPROXIMATION *** |
---|
| 818 | ! |
---|
| 819 | IF(I.EQ.1)GOTO 360 |
---|
| 820 | DHDP=(H(I)-H(I-1))/(PCONV_HPA(I-1)-PCONV_HPA(I)) |
---|
| 821 | DHDP=MAX(DHDP,10.0) |
---|
| 822 | MP(I)=100.*GINV*LV(I)*SIGD*EVAP(I)/DHDP |
---|
| 823 | MP(I)=MAX(MP(I),0.0) |
---|
| 824 | ! |
---|
| 825 | ! *** ADD SMALL AMOUNT OF INERTIA TO DOWNDRAFT *** |
---|
| 826 | ! |
---|
| 827 | FAC=20.0/(PHCONV_HPA(I-1)-PHCONV_HPA(I)) |
---|
| 828 | MP(I)=(FAC*MP(I+1)+MP(I))/(1.+FAC) |
---|
| 829 | ! |
---|
| 830 | ! *** FORCE MP TO DECREASE LINEARLY TO ZERO *** |
---|
| 831 | ! *** BETWEEN ABOUT 950 MB AND THE SURFACE *** |
---|
| 832 | ! |
---|
| 833 | IF(PCONV_HPA(I).GT.(0.949*PCONV_HPA(1)))THEN |
---|
| 834 | JTT=MAX(JTT,I) |
---|
| 835 | MP(I)=MP(JTT)*(PCONV_HPA(1)-PCONV_HPA(I))/(PCONV_HPA(1)- & |
---|
| 836 | PCONV_HPA(JTT)) |
---|
| 837 | END IF |
---|
| 838 | 360 CONTINUE |
---|
| 839 | ! |
---|
| 840 | ! *** FIND MIXING RATIO OF PRECIPITATING DOWNDRAFT *** |
---|
| 841 | ! |
---|
| 842 | IF(I.EQ.INB)GOTO 400 |
---|
| 843 | IF(I.EQ.1)THEN |
---|
| 844 | QSTM=QSCONV(1) |
---|
| 845 | ELSE |
---|
| 846 | QSTM=QSCONV(I-1) |
---|
| 847 | END IF |
---|
| 848 | IF(MP(I).GT.MP(I+1))THEN |
---|
| 849 | RAT=MP(I+1)/MP(I) |
---|
| 850 | QP(I)=QP(I+1)*RAT+QCONV(I)*(1.0-RAT)+100.*GINV* & |
---|
| 851 | SIGD*(PHCONV_HPA(I)-PHCONV_HPA(I+1))*(EVAP(I)/MP(I)) |
---|
| 852 | ELSE |
---|
| 853 | IF(MP(I+1).GT.0.0)THEN |
---|
| 854 | QP(I)=(GZ(I+1)-GZ(I)+QP(I+1)*(LV(I+1)+TCONV(I+1)*(CL-CPD))+ & |
---|
| 855 | CPD*(TCONV(I+1)-TCONV(I)))/(LV(I)+TCONV(I)*(CL-CPD)) |
---|
| 856 | END IF |
---|
| 857 | END IF |
---|
| 858 | QP(I)=MIN(QP(I),QSTM) |
---|
| 859 | QP(I)=MAX(QP(I),0.0) |
---|
| 860 | 400 CONTINUE |
---|
| 861 | END DO |
---|
| 862 | ! |
---|
| 863 | ! *** CALCULATE SURFACE PRECIPITATION IN MM/DAY *** |
---|
| 864 | ! |
---|
| 865 | PRECIP=PRECIP+WT(1)*SIGD*WATER(1)*3600.*24000./(ROWL*G) |
---|
| 866 | ! |
---|
| 867 | 405 CONTINUE |
---|
| 868 | ! |
---|
| 869 | ! *** CALCULATE DOWNDRAFT VELOCITY SCALE AND SURFACE TEMPERATURE AND *** |
---|
| 870 | ! *** WATER VAPOR FLUCTUATIONS *** |
---|
| 871 | ! |
---|
| 872 | WD=BETA*ABS(MP(ICB))*0.01*RD*TCONV(ICB)/(SIGD*PCONV_HPA(ICB)) |
---|
| 873 | QPRIME=0.5*(QP(1)-QCONV(1)) |
---|
| 874 | TPRIME=LV0*QPRIME/CPD |
---|
| 875 | ! |
---|
| 876 | ! *** CALCULATE TENDENCIES OF LOWEST LEVEL POTENTIAL TEMPERATURE *** |
---|
| 877 | ! *** AND MIXING RATIO *** |
---|
| 878 | ! |
---|
| 879 | DPINV=0.01/(PHCONV_HPA(1)-PHCONV_HPA(2)) |
---|
| 880 | AM=0.0 |
---|
| 881 | IF(NK.EQ.1)THEN |
---|
| 882 | DO K=2,INB |
---|
| 883 | AM=AM+M(K) |
---|
| 884 | END DO |
---|
| 885 | END IF |
---|
| 886 | ! save saturated upward mass flux for first level |
---|
| 887 | FUP(1)=AM |
---|
| 888 | IF((2.*G*DPINV*AM).GE.DELTI)IFLAG=4 |
---|
| 889 | FT(1)=FT(1)+G*DPINV*AM*(TCONV(2)-TCONV(1)+(GZ(2)-GZ(1))/CPN(1)) |
---|
| 890 | FT(1)=FT(1)-LVCP(1)*SIGD*EVAP(1) |
---|
| 891 | FT(1)=FT(1)+SIGD*WT(2)*(CL-CPD)*WATER(2)*(TCONV(2)- & |
---|
| 892 | TCONV(1))*DPINV/CPN(1) |
---|
| 893 | FQ(1)=FQ(1)+G*MP(2)*(QP(2)-QCONV(1))* & |
---|
| 894 | DPINV+SIGD*EVAP(1) |
---|
| 895 | FQ(1)=FQ(1)+G*AM*(QCONV(2)-QCONV(1))*DPINV |
---|
| 896 | DO J=2,INB |
---|
| 897 | FQ(1)=FQ(1)+G*DPINV*MENT(J,1)*(QENT(J,1)-QCONV(1)) |
---|
| 898 | END DO |
---|
| 899 | ! |
---|
| 900 | ! *** CALCULATE TENDENCIES OF POTENTIAL TEMPERATURE AND MIXING RATIO *** |
---|
| 901 | ! *** AT LEVELS ABOVE THE LOWEST LEVEL *** |
---|
| 902 | ! |
---|
| 903 | ! *** FIRST FIND THE NET SATURATED UPDRAFT AND DOWNDRAFT MASS FLUXES *** |
---|
| 904 | ! *** THROUGH EACH LEVEL *** |
---|
| 905 | ! |
---|
| 906 | DO I=2,INB |
---|
| 907 | DPINV=0.01/(PHCONV_HPA(I)-PHCONV_HPA(I+1)) |
---|
| 908 | CPINV=1.0/CPN(I) |
---|
| 909 | AMP1=0.0 |
---|
| 910 | AD=0.0 |
---|
| 911 | IF(I.GE.NK)THEN |
---|
| 912 | DO K=I+1,INB+1 |
---|
| 913 | AMP1=AMP1+M(K) |
---|
| 914 | END DO |
---|
| 915 | END IF |
---|
| 916 | DO K=1,I |
---|
| 917 | DO J=I+1,INB+1 |
---|
| 918 | AMP1=AMP1+MENT(K,J) |
---|
| 919 | END DO |
---|
| 920 | END DO |
---|
| 921 | ! save saturated upward mass flux |
---|
| 922 | FUP(I)=AMP1 |
---|
| 923 | IF((2.*G*DPINV*AMP1).GE.DELTI)IFLAG=4 |
---|
| 924 | DO K=1,I-1 |
---|
| 925 | DO J=I,INB |
---|
| 926 | AD=AD+MENT(J,K) |
---|
| 927 | END DO |
---|
| 928 | END DO |
---|
| 929 | ! save saturated downward mass flux |
---|
| 930 | FDOWN(I)=AD |
---|
| 931 | FT(I)=FT(I)+G*DPINV*(AMP1*(TCONV(I+1)-TCONV(I)+(GZ(I+1)-GZ(I))* & |
---|
| 932 | CPINV)-AD*(TCONV(I)-TCONV(I-1)+(GZ(I)-GZ(I-1))*CPINV)) & |
---|
| 933 | -SIGD*LVCP(I)*EVAP(I) |
---|
| 934 | FT(I)=FT(I)+G*DPINV*MENT(I,I)*(HP(I)-H(I)+ & |
---|
| 935 | TCONV(I)*(CPV-CPD)*(QCONV(I)-QENT(I,I)))*CPINV |
---|
| 936 | FT(I)=FT(I)+SIGD*WT(I+1)*(CL-CPD)*WATER(I+1)* & |
---|
| 937 | (TCONV(I+1)-TCONV(I))*DPINV*CPINV |
---|
| 938 | FQ(I)=FQ(I)+G*DPINV*(AMP1*(QCONV(I+1)-QCONV(I))- & |
---|
| 939 | AD*(QCONV(I)-QCONV(I-1))) |
---|
| 940 | DO K=1,I-1 |
---|
| 941 | AWAT=ELIJ(K,I)-(1.-EP(I))*CLW(I) |
---|
| 942 | AWAT=MAX(AWAT,0.0) |
---|
| 943 | FQ(I)=FQ(I)+G*DPINV*MENT(K,I)*(QENT(K,I)-AWAT-QCONV(I)) |
---|
| 944 | END DO |
---|
| 945 | DO K=I,INB |
---|
| 946 | FQ(I)=FQ(I)+G*DPINV*MENT(K,I)*(QENT(K,I)-QCONV(I)) |
---|
| 947 | END DO |
---|
| 948 | FQ(I)=FQ(I)+SIGD*EVAP(I)+G*(MP(I+1)* & |
---|
| 949 | (QP(I+1)-QCONV(I))-MP(I)*(QP(I)-QCONV(I-1)))*DPINV |
---|
| 950 | END DO |
---|
| 951 | ! |
---|
| 952 | ! *** Adjust tendencies at top of convection layer to reflect *** |
---|
| 953 | ! *** actual position of the level zero CAPE *** |
---|
| 954 | ! |
---|
| 955 | FQOLD=FQ(INB) |
---|
| 956 | FQ(INB)=FQ(INB)*(1.-FRAC) |
---|
| 957 | FQ(INB-1)=FQ(INB-1)+FRAC*FQOLD*((PHCONV_HPA(INB)- & |
---|
| 958 | PHCONV_HPA(INB+1))/ & |
---|
| 959 | (PHCONV_HPA(INB-1)-PHCONV_HPA(INB)))*LV(INB)/LV(INB-1) |
---|
| 960 | FTOLD=FT(INB) |
---|
| 961 | FT(INB)=FT(INB)*(1.-FRAC) |
---|
| 962 | FT(INB-1)=FT(INB-1)+FRAC*FTOLD*((PHCONV_HPA(INB)- & |
---|
| 963 | PHCONV_HPA(INB+1))/ & |
---|
| 964 | (PHCONV_HPA(INB-1)-PHCONV_HPA(INB)))*CPN(INB)/CPN(INB-1) |
---|
| 965 | ! |
---|
| 966 | ! *** Very slightly adjust tendencies to force exact *** |
---|
| 967 | ! *** enthalpy, momentum and tracer conservation *** |
---|
| 968 | ! |
---|
| 969 | ENTS=0.0 |
---|
| 970 | DO I=1,INB |
---|
| 971 | ENTS=ENTS+(CPN(I)*FT(I)+LV(I)*FQ(I))* & |
---|
| 972 | (PHCONV_HPA(I)-PHCONV_HPA(I+1)) |
---|
| 973 | END DO |
---|
| 974 | ENTS=ENTS/(PHCONV_HPA(1)-PHCONV_HPA(INB+1)) |
---|
| 975 | DO I=1,INB |
---|
| 976 | FT(I)=FT(I)-ENTS/CPN(I) |
---|
| 977 | END DO |
---|
| 978 | |
---|
| 979 | ! ************************************************ |
---|
| 980 | ! **** DETERMINE MASS DISPLACEMENT MATRIX |
---|
| 981 | ! ***** AND COMPENSATING SUBSIDENCE |
---|
| 982 | ! ************************************************ |
---|
| 983 | |
---|
| 984 | ! mass displacement matrix due to saturated up-and downdrafts |
---|
| 985 | ! inside the cloud and determine compensating subsidence |
---|
| 986 | ! FUP(I) (saturated updrafts), FDOWN(I) (saturated downdrafts) are assumed to be |
---|
| 987 | ! balanced by compensating subsidence (SUB(I)) |
---|
| 988 | ! FDOWN(I) and SUB(I) defined positive downwards |
---|
| 989 | |
---|
| 990 | ! NCONVTOP IS THE TOP LEVEL AT WHICH CONVECTIVE MASS FLUXES ARE DIAGNOSED |
---|
| 991 | ! EPSILON IS A SMALL NUMBER |
---|
| 992 | |
---|
| 993 | SUB(1)=0. |
---|
| 994 | NCONVTOP=1 |
---|
| 995 | do i=1,INB+1 |
---|
| 996 | do j=1,INB+1 |
---|
| 997 | if (j.eq.NK) then |
---|
| 998 | FMASS(j,i)=FMASS(j,i)+M(i) |
---|
| 999 | endif |
---|
| 1000 | FMASS(j,i)=FMASS(j,i)+MENT(j,i) |
---|
| 1001 | IF (FMASS(J,I).GT.EPSILON) NCONVTOP=MAX(NCONVTOP,I,J) |
---|
| 1002 | end do |
---|
| 1003 | if (i.gt.1) then |
---|
| 1004 | SUB(i)=FUP(i-1)-FDOWN(i) |
---|
| 1005 | endif |
---|
| 1006 | end do |
---|
| 1007 | NCONVTOP=NCONVTOP+1 |
---|
| 1008 | |
---|
| 1009 | RETURN |
---|
| 1010 | ! |
---|
| 1011 | END SUBROUTINE CONVECT |
---|
| 1012 | ! |
---|
| 1013 | ! --------------------------------------------------------------------------- |
---|
| 1014 | ! |
---|
| 1015 | SUBROUTINE TLIFT(GZ,ICB,NK,TVP,TPK,CLW,ND,NL,KK) |
---|
| 1016 | ! |
---|
| 1017 | !-cv |
---|
| 1018 | use par_mod |
---|
| 1019 | use conv_mod |
---|
| 1020 | |
---|
| 1021 | implicit none |
---|
| 1022 | !-cv |
---|
| 1023 | !====>Begin Module TLIFT File convect.f Undeclared variables |
---|
| 1024 | ! |
---|
| 1025 | !Argument variables |
---|
| 1026 | ! |
---|
| 1027 | integer :: icb, kk, nd, nk, nl |
---|
| 1028 | ! |
---|
| 1029 | !Local variables |
---|
| 1030 | ! |
---|
| 1031 | integer :: i, j, nsb, nst |
---|
| 1032 | ! |
---|
| 1033 | real :: ah0, ahg, alv, cpinv, cpp, denom |
---|
| 1034 | real :: es, qg, rg, s, tc, tg |
---|
| 1035 | ! |
---|
| 1036 | !====>End Module TLIFT File convect.f |
---|
| 1037 | |
---|
| 1038 | REAL :: GZ(ND),TPK(ND),CLW(ND) |
---|
| 1039 | REAL :: TVP(ND) |
---|
| 1040 | ! |
---|
| 1041 | ! *** ASSIGN VALUES OF THERMODYNAMIC CONSTANTS *** |
---|
| 1042 | ! |
---|
| 1043 | REAL,PARAMETER :: CPD=1005.7 |
---|
| 1044 | REAL,PARAMETER :: CPV=1870.0 |
---|
| 1045 | REAL,PARAMETER :: CL=2500.0 |
---|
| 1046 | REAL,PARAMETER :: RV=461.5 |
---|
| 1047 | REAL,PARAMETER :: RD=287.04 |
---|
| 1048 | REAL,PARAMETER :: LV0=2.501E6 |
---|
| 1049 | ! |
---|
| 1050 | REAL,PARAMETER :: CPVMCL=CL-CPV |
---|
| 1051 | REAL,PARAMETER :: EPS0=RD/RV |
---|
| 1052 | REAL,PARAMETER :: EPSI=1./EPS0 |
---|
| 1053 | ! |
---|
| 1054 | ! *** CALCULATE CERTAIN PARCEL QUANTITIES, INCLUDING STATIC ENERGY *** |
---|
| 1055 | ! |
---|
| 1056 | AH0=(CPD*(1.-QCONV(NK))+CL*QCONV(NK))*TCONV(NK)+QCONV(NK)* & |
---|
| 1057 | (LV0-CPVMCL*( & |
---|
| 1058 | TCONV(NK)-273.15))+GZ(NK) |
---|
| 1059 | CPP=CPD*(1.-QCONV(NK))+QCONV(NK)*CPV |
---|
| 1060 | CPINV=1./CPP |
---|
| 1061 | ! |
---|
| 1062 | IF(KK.EQ.1)THEN |
---|
| 1063 | ! |
---|
| 1064 | ! *** CALCULATE LIFTED PARCEL QUANTITIES BELOW CLOUD BASE *** |
---|
| 1065 | ! |
---|
| 1066 | DO I=1,ICB-1 |
---|
| 1067 | CLW(I)=0.0 |
---|
| 1068 | END DO |
---|
| 1069 | DO I=NK,ICB-1 |
---|
| 1070 | TPK(I)=TCONV(NK)-(GZ(I)-GZ(NK))*CPINV |
---|
| 1071 | TVP(I)=TPK(I)*(1.+QCONV(NK)*EPSI) |
---|
| 1072 | END DO |
---|
| 1073 | END IF |
---|
| 1074 | ! |
---|
| 1075 | ! *** FIND LIFTED PARCEL QUANTITIES ABOVE CLOUD BASE *** |
---|
| 1076 | ! |
---|
| 1077 | NST=ICB |
---|
| 1078 | NSB=ICB |
---|
| 1079 | IF(KK.EQ.2)THEN |
---|
| 1080 | NST=NL |
---|
| 1081 | NSB=ICB+1 |
---|
| 1082 | END IF |
---|
| 1083 | DO I=NSB,NST |
---|
| 1084 | TG=TCONV(I) |
---|
| 1085 | QG=QSCONV(I) |
---|
| 1086 | ALV=LV0-CPVMCL*(TCONV(I)-273.15) |
---|
| 1087 | DO J=1,2 |
---|
| 1088 | S=CPD+ALV*ALV*QG/(RV*TCONV(I)*TCONV(I)) |
---|
| 1089 | S=1./S |
---|
| 1090 | AHG=CPD*TG+(CL-CPD)*QCONV(NK)*TCONV(I)+ALV*QG+GZ(I) |
---|
| 1091 | TG=TG+S*(AH0-AHG) |
---|
| 1092 | TG=MAX(TG,35.0) |
---|
| 1093 | TC=TG-273.15 |
---|
| 1094 | DENOM=243.5+TC |
---|
| 1095 | IF(TC.GE.0.0)THEN |
---|
| 1096 | ES=6.112*EXP(17.67*TC/DENOM) |
---|
| 1097 | ELSE |
---|
| 1098 | ES=EXP(23.33086-6111.72784/TG+0.15215*LOG(TG)) |
---|
| 1099 | END IF |
---|
| 1100 | QG=EPS0*ES/(PCONV_HPA(I)-ES*(1.-EPS0)) |
---|
| 1101 | END DO |
---|
| 1102 | ALV=LV0-CPVMCL*(TCONV(I)-273.15) |
---|
| 1103 | TPK(I)=(AH0-(CL-CPD)*QCONV(NK)*TCONV(I)-GZ(I)-ALV*QG)/CPD |
---|
| 1104 | CLW(I)=QCONV(NK)-QG |
---|
| 1105 | CLW(I)=MAX(0.0,CLW(I)) |
---|
| 1106 | RG=QG/(1.-QCONV(NK)) |
---|
| 1107 | TVP(I)=TPK(I)*(1.+RG*EPSI) |
---|
| 1108 | END DO |
---|
| 1109 | RETURN |
---|
| 1110 | END SUBROUTINE TLIFT |
---|