Changeset 4803 for LMDZ6/trunk/libf/phylmd/lmdz_lscp.F90
- Timestamp:
- Feb 5, 2024, 10:16:07 PM (9 months ago)
- File:
-
- 1 edited
-
LMDZ6/trunk/libf/phylmd/lmdz_lscp.F90 (modified) (23 diffs)
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/trunk/libf/phylmd/lmdz_lscp.F90
r4686 r4803 18 18 qclr, qcld, qss, qvc, rnebclr, rnebss, gamma_ss, & 19 19 Tcontr, qcontr, qcontr2, fcontrN, fcontrP, & 20 cloudth_sth,cloudth_senv,cloudth_sigmath,cloudth_sigmaenv) 20 cloudth_sth,cloudth_senv,cloudth_sigmath,cloudth_sigmaenv, & 21 dqreva,dqssub,dqrauto,dqrcol,dqrmelt,dqrfreez,dqsauto, & 22 dqsagg,dqsrim,dqsmelt,dqsfreez) 21 23 22 24 !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ … … 96 98 USE lmdz_lscp_tools, ONLY : fallice_velocity, distance_to_cloud_top 97 99 USE ice_sursat_mod, ONLY : ice_sursat 100 USE lmdz_lscp_poprecip, ONLY : poprecip_evapsub, poprecip_postcld 98 101 99 102 ! Use du module lmdz_lscp_ini contenant les constantes … … 106 109 USE lmdz_lscp_ini, ONLY : iflag_autoconversion, ffallv_con, ffallv_lsc 107 110 USE lmdz_lscp_ini, ONLY : RCPD, RLSTT, RLVTT, RLMLT, RVTMP2, RTT, RD, RG 111 USE lmdz_lscp_ini, ONLY : ok_poprecip 108 112 109 113 IMPLICIT NONE … … 197 201 REAL, DIMENSION(klon,klev), INTENT(OUT) :: cloudth_sigmaenv ! std of saturation deficit in environment 198 202 203 ! for POPRECIP 204 205 REAL, DIMENSION(klon,klev), INTENT(OUT) :: dqreva !-- rain tendendy due to evaporation [kg/kg/s] 206 REAL, DIMENSION(klon,klev), INTENT(OUT) :: dqssub !-- snow tendency due to sublimation [kg/kg/s] 207 REAL, DIMENSION(klon,klev), INTENT(OUT) :: dqrcol !-- rain tendendy due to collection by rain of liquid cloud droplets [kg/kg/s] 208 REAL, DIMENSION(klon,klev), INTENT(OUT) :: dqsagg !-- snow tendency due to collection of lcoud ice by aggregation [kg/kg/s] 209 REAL, DIMENSION(klon,klev), INTENT(OUT) :: dqrauto !-- rain tendency due to autoconversion of cloud liquid [kg/kg/s] 210 REAL, DIMENSION(klon,klev), INTENT(OUT) :: dqsauto !-- snow tendency due to autoconversion of cloud ice [kg/kg/s] 211 REAL, DIMENSION(klon,klev), INTENT(OUT) :: dqsrim !-- snow tendency due to riming [kg/kg/s] 212 REAL, DIMENSION(klon,klev), INTENT(OUT) :: dqsmelt !-- snow tendency due to melting [kg/kg/s] 213 REAL, DIMENSION(klon,klev), INTENT(OUT) :: dqrmelt !-- rain tendency due to melting [kg/kg/s] 214 REAL, DIMENSION(klon,klev), INTENT(OUT) :: dqsfreez !-- snow tendency due to freezing [kg/kg/s] 215 REAL, DIMENSION(klon,klev), INTENT(OUT) :: dqrfreez !-- rain tendency due to freezing [kg/kg/s] 199 216 200 217 … … 243 260 REAL :: zfrac_lessi 244 261 REAL, DIMENSION(klon) :: zprec_cond 245 REAL , DIMENSION(klon):: zmair262 REAL :: zmair 246 263 REAL :: zcpair, zcpeau 247 264 REAL, DIMENSION(klon) :: zlh_solid 265 REAL, DIMENSION(klon) :: ztupnew 248 266 REAL :: zm_solid ! for liquid -> solid conversion 249 267 REAL, DIMENSION(klon) :: zrflclr, zrflcld … … 266 284 INTEGER ncoreczq 267 285 INTEGER, DIMENSION(klon,klev) :: mpc_bl_points 286 LOGICAL iftop 268 287 269 288 LOGICAL, DIMENSION(klon) :: lognormale … … 339 358 distcltop1D(:)=0.0 340 359 temp_cltop1D(:) = 0.0 360 ztupnew(:)=0.0 341 361 !--ice supersaturation 342 362 gamma_ss(:,:) = 1.0 … … 350 370 distcltop(:,:)=0. 351 371 temp_cltop(:,:)=0. 372 !-- poprecip 373 dqreva(:,:)=0.0 374 dqrauto(:,:)=0.0 375 dqrmelt(:,:)=0.0 376 dqrfreez(:,:)=0.0 377 dqrcol(:,:)=0.0 378 dqssub(:,:)=0.0 379 dqsauto(:,:)=0.0 380 dqsrim(:,:)=0.0 381 dqsagg(:,:)=0.0 382 dqsfreez(:,:)=0.0 383 dqsmelt(:,:)=0.0 384 385 386 352 387 353 388 !c_iso: variable initialisation for iso … … 361 396 362 397 DO k = klev, 1, -1 398 399 IF (k.LE.klev-1) THEN 400 iftop=.false. 401 ELSE 402 iftop=.true. 403 ENDIF 363 404 364 405 ! Initialisation temperature and specific humidity … … 366 407 zt(i)=temp(i,k) 367 408 zq(i)=qt(i,k) 409 IF (.not. iftop) THEN 410 ztupnew(i)=temp(i,k+1)+d_t(i,k+1) 411 ENDIF 368 412 !c_iso init of iso 369 413 ENDDO 370 414 415 !================================================================ 416 ! Flag for the new and more microphysical treatment of precipitation from Atelier Nuage (R) 417 IF (ok_poprecip) THEN 418 419 CALL poprecip_evapsub(klon, dtime, iftop, paprs(:,k), paprs(:,k+1), pplay(:,k), & 420 zt, ztupnew, zq, zmqc, znebprecipclr, znebprecipcld, & 421 zrfl, zrflclr, zrflcld, & 422 zifl, ziflclr, ziflcld, & 423 dqreva(:,k),dqssub(:,k) & 424 ) 425 426 !================================================================ 427 ELSE 428 371 429 ! -------------------------------------------------------------------- 372 ! P 0> Thermalization of precipitation falling from the overlying layer430 ! P1> Thermalization of precipitation falling from the overlying layer 373 431 ! -------------------------------------------------------------------- 374 432 ! Computes air temperature variation due to enthalpy transported by … … 385 443 ! --------------------------------------------------------------------- 386 444 387 IF (k.LE.klev-1) THEN 445 IF (iftop) THEN 446 447 DO i = 1, klon 448 zmqc(i) = 0. 449 ENDDO 450 451 ELSE 388 452 389 453 DO i = 1, klon 390 454 391 zmair (i)=(paprs(i,k)-paprs(i,k+1))/RG455 zmair=(paprs(i,k)-paprs(i,k+1))/RG 392 456 ! no condensed water so cp=cp(vapor+dry air) 393 457 ! RVTMP2=rcpv/rcpd-1 … … 398 462 ! layer's air so that precipitation at the ground has the 399 463 ! same temperature as the lowermost layer 400 zmqc(i) = (zrfl(i)+zifl(i))*dtime/zmair (i)464 zmqc(i) = (zrfl(i)+zifl(i))*dtime/zmair 401 465 ! t(i,k+1)+d_t(i,k+1): new temperature of the overlying layer 402 zt(i) = ( (temp(i,k+1)+d_t(i,k+1))*zmqc(i)*zcpeau + zcpair*zt(i) ) &466 zt(i) = ( ztupnew(i)*zmqc(i)*zcpeau + zcpair*zt(i) ) & 403 467 / (zcpair + zmqc(i)*zcpeau) 404 468 405 ENDDO406 407 ELSE408 409 DO i = 1, klon410 zmair(i)=(paprs(i,k)-paprs(i,k+1))/RG411 zmqc(i) = 0.412 469 ENDDO 413 470 … … 415 472 416 473 ! -------------------------------------------------------------------- 417 ! P 1> Precipitation evaporation/sublimation474 ! P2> Precipitation evaporation/sublimation/melting 418 475 ! -------------------------------------------------------------------- 419 ! A part of the precipitation coming from above is evaporated/sublimated .476 ! A part of the precipitation coming from above is evaporated/sublimated/melted. 420 477 ! -------------------------------------------------------------------- 421 478 422 479 IF (iflag_evap_prec.GE.1) THEN 423 480 … … 506 563 ! redistribute zqev0 conserving the ratio liquid/ice 507 564 508 ! todoan: check the consistency of this formula509 565 IF (zqevt+zqevti.GT.zqev0) THEN 510 566 zqev=zqev0*zqevt/(zqevt+zqevti) … … 524 580 525 581 ! vapor, temperature, precip fluxes update 582 ! vapor is updated after evaporation/sublimation (it is increased) 526 583 zq(i) = zq(i) - (zrfln(i)+zifln(i)-zrfl(i)-zifl(i)) & 527 584 * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime 585 ! zmqc is the total condensed water in the precip flux (it is decreased) 528 586 zmqc(i) = zmqc(i) + (zrfln(i)+zifln(i)-zrfl(i)-zifl(i)) & 529 587 * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime 588 ! air and precip temperature (i.e., gridbox temperature) 589 ! is updated due to latent heat cooling 530 590 zt(i) = zt(i) + (zrfln(i)-zrfl(i)) & 531 591 * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime & … … 579 639 ! c_iso: melting of isotopic precipi with zmelt (no fractionation) 580 640 581 ! Latent heat of melting with precipitation thermalization 641 ! Latent heat of melting because of precipitation melting 642 ! NB: the air + precip temperature is simultaneously updated 582 643 zt(i)=zt(i)-zifl(i)*zmelt*(RG*dtime)/(paprs(i,k)-paprs(i,k+1)) & 583 644 *RLMLT/RCPD/(1.0+RVTMP2*(zq(i)+zmqc(i))) … … 590 651 ENDIF ! (zrfl(i)+zifl(i).GT.0.) 591 652 592 ENDDO 653 ENDDO ! loop on klon 593 654 594 655 ENDIF ! (iflag_evap_prec>=1) 656 657 ENDIF ! (ok_poprecip) 595 658 596 659 ! -------------------------------------------------------------------- … … 647 710 648 711 649 650 651 652 712 ELSEIF (iflag_cloudth_vert.GE.3 .AND. iflag_cloudth_vert.LE.5) THEN 653 713 … … 861 921 862 922 ! Partition function in stratiform clouds (will be overwritten in boundary-layer MPCs) 863 CALL icefrac_lscp(klon,zt (:),iflag_ice_thermo,distcltop1D(:),temp_cltop1D(:),zfice(:), dzfice(:))923 CALL icefrac_lscp(klon,zt,iflag_ice_thermo,distcltop1D,temp_cltop1D,zfice, dzfice) 864 924 865 925 … … 918 978 ! ---------------------------------------------------------------- 919 979 980 !================================================================ 981 IF (ok_poprecip) THEN 982 983 DO i = 1, klon 984 zoliql(i) = zcond(i) * ( 1. - zfice(i) ) 985 zoliqi(i) = zcond(i) * zfice(i) 986 ENDDO 987 988 CALL poprecip_postcld(klon, dtime, paprs(:,k), paprs(:,k+1), pplay(:,k), & 989 ctot_vol(:,k), ptconv(:,k), & 990 zt, zq, zoliql, zoliqi, zfice, & 991 rneb(:,k), znebprecipclr, znebprecipcld, & 992 zrfl, zrflclr, zrflcld, & 993 zifl, ziflclr, ziflcld, & 994 dqrauto(:,k),dqrcol(:,k),dqrmelt(:,k),dqrfreez(:,k), & 995 dqsauto(:,k),dqsagg(:,k),dqsrim(:,k),dqsmelt(:,k),dqsfreez(:,k) & 996 ) 997 998 DO i = 1, klon 999 zoliq(i) = zoliql(i) + zoliqi(i) 1000 IF ( zoliq(i) .GT. 0. ) THEN 1001 zfice(i) = zoliqi(i)/zoliq(i) 1002 ELSE 1003 zfice(i) = 0.0 1004 ENDIF 1005 ! when poprecip activated, radiation does not see any precipitation content 1006 radocond(i,k) = zoliq(i) 1007 radocondl(i,k)= radocond(i,k)*(1.-zfice(i)) 1008 radocondi(i,k)= radocond(i,k)*zfice(i) 1009 ENDDO 1010 1011 !================================================================ 1012 ELSE 1013 920 1014 ! LTP: 921 1015 IF (iflag_evap_prec .GE. 4) THEN … … 982 1076 zneb(i) = MAX(rneb(i,k),seuil_neb) 983 1077 radocond(i,k) = zoliq(i)/REAL(niter_lscp+1) 984 radicefrac(i,k) = zfice(i)985 1078 radocondi(i,k)=zoliq(i)*zfice(i)/REAL(niter_lscp+1) 986 1079 radocondl(i,k)=zoliq(i)*(1.-zfice(i))/REAL(niter_lscp+1) … … 1195 1288 ENDDO 1196 1289 1197 1198 1290 ENDIF 1291 1292 ENDIF ! ok_poprecip 1199 1293 1200 1294 ! End of precipitation formation 1201 1295 ! -------------------------------- 1202 1296 1203 ! Outputs: 1204 ! Precipitation fluxes at layer interfaces 1205 ! + precipitation fractions + 1206 ! temperature and water species tendencies 1207 DO i = 1, klon 1208 psfl(i,k)=zifl(i) 1209 prfl(i,k)=zrfl(i) 1210 pfraclr(i,k)=znebprecipclr(i) 1211 pfracld(i,k)=znebprecipcld(i) 1212 d_ql(i,k) = (1-zfice(i))*zoliq(i) 1213 d_qi(i,k) = zfice(i)*zoliq(i) 1214 d_q(i,k) = zq(i) - qt(i,k) 1215 ! c_iso: same for isotopes 1216 d_t(i,k) = zt(i) - temp(i,k) 1217 ENDDO 1297 1298 ! Calculation of cloud condensates amount seen by radiative scheme 1299 !----------------------------------------------------------------- 1218 1300 1219 1301 ! Calculation of the concentration of condensates seen by the radiation scheme 1220 1302 ! for liquid phase, we take radocondl 1221 1303 ! for ice phase, we take radocondi if we neglect snowfall, otherwise (ok_radocond_snow=true) 1222 ! we recaulate radocondi to account for contributions from the precipitation flux 1223 1224 IF ((ok_radocond_snow) .AND. (k .LT. klev)) THEN 1304 ! we recalculate radocondi to account for contributions from the precipitation flux 1305 ! TODO: for the moment, we deactivate this option when ok_poprecip=.true. 1306 1307 IF ((ok_radocond_snow) .AND. (k .LT. klev) .AND. (.NOT. ok_poprecip)) THEN 1225 1308 ! for the solid phase (crystals + snowflakes) 1226 1309 ! we recalculate radocondi by summing … … 1319 1402 CALL calc_qsat_ecmwf(klon,Tbef(:),qzero(:),pplay(:,k),RTT,2,.false.,qsats(:,k),zdqs(:)) 1320 1403 1404 1405 1406 ! Outputs: 1407 !------------------------------- 1408 ! Precipitation fluxes at layer interfaces 1409 ! + precipitation fractions + 1410 ! temperature and water species tendencies 1411 DO i = 1, klon 1412 psfl(i,k)=zifl(i) 1413 prfl(i,k)=zrfl(i) 1414 pfraclr(i,k)=znebprecipclr(i) 1415 pfracld(i,k)=znebprecipcld(i) 1416 d_ql(i,k) = (1-zfice(i))*zoliq(i) 1417 d_qi(i,k) = zfice(i)*zoliq(i) 1418 d_q(i,k) = zq(i) - qt(i,k) 1419 ! c_iso: same for isotopes 1420 d_t(i,k) = zt(i) - temp(i,k) 1421 ENDDO 1422 1423 1321 1424 ENDDO 1322 1425 1323 !======================================================================1324 ! END OF VERTICAL LOOP1325 !======================================================================1326 1426 1327 1427 ! Rain or snow at the surface (depending on the first layer temperature)
Note: See TracChangeset
for help on using the changeset viewer.
