Changeset 5408 for LMDZ6/trunk/libf
- Timestamp:
- Dec 15, 2024, 3:48:56 PM (6 hours ago)
- Location:
- LMDZ6/trunk/libf/phylmd
- Files:
-
- 1 added
- 1 deleted
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
LMDZ6/trunk/libf/phylmd/lmdz_lscp.f90
r5406 r5408 108 108 USE lmdz_lscp_tools, ONLY : fallice_velocity, distance_to_cloud_top 109 109 USE lmdz_lscp_condensation, ONLY : condensation_lognormal, condensation_ice_supersat 110 USE lmdz_lscp_p oprecip, ONLY :poprecip_precld, poprecip_postcld110 USE lmdz_lscp_precip, ONLY : histprecip_precld, poprecip_precld, poprecip_postcld 111 111 112 112 ! Use du module lmdz_lscp_ini contenant les constantes 113 113 USE lmdz_lscp_ini, ONLY : prt_level, lunout, eps 114 USE lmdz_lscp_ini, ONLY : seuil_neb, niter_lscp, iflag_evap_prec, t_coup, DDT0, ztfondue, rain_int_min114 USE lmdz_lscp_ini, ONLY : seuil_neb, niter_lscp, iflag_evap_prec, DDT0, rain_int_min 115 115 USE lmdz_lscp_ini, ONLY : ok_radocond_snow, a_tr_sca, cld_expo_con, cld_expo_lsc 116 116 USE lmdz_lscp_ini, ONLY : iflag_cloudth_vert, iflag_rain_incloud_vol, iflag_t_glace, t_glace_min 117 USE lmdz_lscp_ini, ONLY : c oef_eva, coef_sub,cld_tau_lsc, cld_tau_con, cld_lc_lsc, cld_lc_con117 USE lmdz_lscp_ini, ONLY : cld_tau_lsc, cld_tau_con, cld_lc_lsc, cld_lc_con 118 118 USE lmdz_lscp_ini, ONLY : iflag_bergeron, iflag_fisrtilp_qsat, iflag_vice, cice_velo, dice_velo 119 119 USE lmdz_lscp_ini, ONLY : iflag_autoconversion, ffallv_con, ffallv_lsc, min_frac_th_cld … … 267 267 ! LOCAL VARIABLES: 268 268 !---------------- 269 REAL,DIMENSION(klon) :: qsl, qsi ! saturation threshold at current vertical level270 269 REAL,DIMENSION(klon) :: qice_ini 271 270 REAL :: zct, zcl,zexpo … … 284 283 REAL, DIMENSION(klon) :: zfice_th 285 284 REAL, DIMENSION(klon) :: qcloud, qincloud_mpc 286 REAL, DIMENSION(klon) :: zrfl, zrfln 287 REAL :: zqev, zqevt 288 REAL, DIMENSION(klon) :: zifl, zifln, ziflprev 289 REAL :: zqev0,zqevi, zqevti 285 REAL, DIMENSION(klon) :: zrfl 286 REAL, DIMENSION(klon) :: zifl, ziflprev 290 287 REAL, DIMENSION(klon) :: zoliq, zcond, zq, zqn 291 288 REAL, DIMENSION(klon) :: zoliql, zoliqi … … 294 291 REAL, DIMENSION(klon) :: zdz,iwc 295 292 REAL :: zchau,zfroi 296 REAL, DIMENSION(klon) :: zfice,zneb ,znebprecip297 REAL :: z melt,zrain,zsnow,zprecip293 REAL, DIMENSION(klon) :: zfice,zneb 294 REAL :: zrain,zsnow,zprecip 298 295 REAL, DIMENSION(klon) :: dzfice 299 296 REAL, DIMENSION(klon) :: zfice_turb, dzfice_turb 300 297 REAL :: zsolid 301 298 REAL, DIMENSION(klon) :: qtot, qzero 302 REAL, DIMENSION(klon) :: dqsl,dqsi303 299 REAL :: smallestreal 304 300 ! Variables for Bergeron process … … 310 306 REAL :: zfrac_lessi 311 307 REAL, DIMENSION(klon) :: zprec_cond 312 REAL :: zmair313 REAL :: zcpair, zcpeau314 308 REAL, DIMENSION(klon) :: zlh_solid 315 309 REAL, DIMENSION(klon) :: ztupnew … … 318 312 REAL, DIMENSION(klon) :: zrflclr, zrflcld 319 313 REAL, DIMENSION(klon) :: d_zrfl_clr_cld, d_zifl_clr_cld 320 REAL, DIMENSION(klon) :: d_zrfl_cld_clr, d_zifl_cld_clr 314 REAL, DIMENSION(klon) :: d_zrfl_cld_clr, d_zifl_cld_clr 321 315 REAL, DIMENSION(klon) :: ziflclr, ziflcld 322 316 REAL, DIMENSION(klon) :: znebprecipclr, znebprecipcld … … 369 363 ! Few initialisations 370 364 371 znebprecip(:)=0.0372 365 ctot_vol(1:klon,1:klev)=0.0 373 366 rneblsvol(1:klon,1:klev)=0.0 … … 488 481 DO k = klev, 1, -1 489 482 483 qice_ini = qi_seri(:,k) 484 490 485 IF (k.LE.klev-1) THEN 491 486 iftop=.false. … … 511 506 !c_iso init of iso 512 507 ENDDO 508 509 ! -------------------------------------------------------------------- 510 ! P1> Precipitation processes, before cloud formation: 511 ! Thermalization of precipitation falling from the overlying layer AND 512 ! Precipitation evaporation/sublimation/melting 513 !--------------------------------------------------------------------- 513 514 514 515 !================================================================ 515 516 ! Flag for the new and more microphysical treatment of precipitation from Atelier Nuage (R) 516 IF ( ok_poprecip) THEN517 518 519 520 521 522 523 dqreva(:,k),dqssub(:,k) &524 517 IF ( ok_poprecip .OR. ( iflag_evap_prec .EQ. 6 ) ) THEN 518 519 CALL poprecip_precld(klon, dtime, iftop, paprs(:,k), paprs(:,k+1), pplay(:,k), & 520 zt, ztupnew, zq, zmqc, znebprecipclr, znebprecipcld, & 521 zqvapclr, zqupnew, & 522 zrfl, zrflclr, zrflcld, & 523 zifl, ziflclr, ziflcld, & 524 dqreva(:,k), dqssub(:,k) & 525 ) 525 526 526 527 !================================================================ 527 528 ELSE 528 529 529 ! -------------------------------------------------------------------- 530 ! P1> Thermalization of precipitation falling from the overlying layer 531 ! -------------------------------------------------------------------- 532 ! Computes air temperature variation due to enthalpy transported by 533 ! precipitation. Precipitation is then thermalized with the air in the 534 ! layer. 535 ! The precipitation should remain thermalized throughout the different 536 ! thermodynamical transformations. 537 ! The corresponding water mass should 538 ! be added when calculating the layer's enthalpy change with 539 ! temperature 540 ! See lmdzpedia page todoan 541 ! todoan: check consistency with ice phase 542 ! todoan: understand why several steps 543 ! --------------------------------------------------------------------- 544 545 IF (iftop) THEN 546 547 DO i = 1, klon 548 zmqc(i) = 0. 549 ENDDO 550 551 ELSE 552 553 DO i = 1, klon 554 555 zmair=(paprs(i,k)-paprs(i,k+1))/RG 556 ! no condensed water so cp=cp(vapor+dry air) 557 ! RVTMP2=rcpv/rcpd-1 558 zcpair=RCPD*(1.0+RVTMP2*zq(i)) 559 zcpeau=RCPD*RVTMP2 560 561 ! zmqc: precipitation mass that has to be thermalized with 562 ! layer's air so that precipitation at the ground has the 563 ! same temperature as the lowermost layer 564 zmqc(i) = (zrfl(i)+zifl(i))*dtime/zmair 565 ! t(i,k+1)+d_t(i,k+1): new temperature of the overlying layer 566 zt(i) = ( ztupnew(i)*zmqc(i)*zcpeau + zcpair*zt(i) ) & 567 / (zcpair + zmqc(i)*zcpeau) 568 569 ENDDO 570 571 ENDIF 572 573 qice_ini = qi_seri(:,k) 574 ! -------------------------------------------------------------------- 575 ! P2> Precipitation evaporation/sublimation/melting 576 ! -------------------------------------------------------------------- 577 ! A part of the precipitation coming from above is evaporated/sublimated/melted. 578 ! -------------------------------------------------------------------- 579 580 IF (iflag_evap_prec.GE.1) THEN 581 582 ! Calculation of saturation specific humidity 583 ! depending on temperature: 584 CALL calc_qsat_ecmwf(klon,zt,qzero,pplay(:,k),RTT,0,.false.,zqs,zdqs) 585 ! wrt liquid water 586 CALL calc_qsat_ecmwf(klon,zt,qzero,pplay(:,k),RTT,1,.false.,qsl,dqsl) 587 ! wrt ice 588 CALL calc_qsat_ecmwf(klon,zt,qzero,pplay(:,k),RTT,2,.false.,qsi,dqsi) 589 590 DO i = 1, klon 591 592 ! if precipitation 593 IF (zrfl(i)+zifl(i).GT.0.) THEN 594 595 ! LudoTP: we only account for precipitation evaporation in the clear-sky (iflag_evap_prec>=4). 596 ! c_iso: likely important to distinguish cs from neb isotope precipitation 597 598 IF (iflag_evap_prec.GE.4) THEN 599 zrfl(i) = zrflclr(i) 600 zifl(i) = ziflclr(i) 601 ENDIF 602 603 IF (iflag_evap_prec.EQ.1) THEN 604 znebprecip(i)=zneb(i) 605 ELSE 606 znebprecip(i)=MAX(zneb(i),znebprecip(i)) 607 ENDIF 608 609 IF (iflag_evap_prec.GT.4) THEN 610 ! Max evaporation not to saturate the clear sky precip fraction 611 ! i.e. the fraction where evaporation occurs 612 zqev0 = MAX(0.0, (zqs(i)-zq(i))*znebprecipclr(i)) 613 ELSEIF (iflag_evap_prec .EQ. 4) THEN 614 ! Max evaporation not to saturate the whole mesh 615 ! Pay attention -> lead to unrealistic and excessive evaporation 616 zqev0 = MAX(0.0, zqs(i)-zq(i)) 617 ELSE 618 ! Max evap not to saturate the fraction below the cloud 619 zqev0 = MAX(0.0, (zqs(i)-zq(i))*znebprecip(i)) 620 ENDIF 621 622 ! Evaporation of liquid precipitation coming from above 623 ! dP/dz=beta*(1-q/qsat)*sqrt(P) 624 ! formula from Sundquist 1988, Klemp & Wilhemson 1978 625 ! LTP: evaporation only in the clear sky part (iflag_evap_prec>=4) 626 627 IF (iflag_evap_prec.EQ.3) THEN 628 zqevt = znebprecip(i)*coef_eva*(1.0-zq(i)/qsl(i)) & 629 *SQRT(zrfl(i)/max(1.e-4,znebprecip(i))) & 630 *(paprs(i,k)-paprs(i,k+1))/pplay(i,k)*zt(i)*RD/RG 631 ELSE IF (iflag_evap_prec.GE.4) THEN 632 zqevt = znebprecipclr(i)*coef_eva*(1.0-zq(i)/qsl(i)) & 633 *SQRT(zrfl(i)/max(1.e-8,znebprecipclr(i))) & 634 *(paprs(i,k)-paprs(i,k+1))/pplay(i,k)*zt(i)*RD/RG 635 ELSE 636 zqevt = 1.*coef_eva*(1.0-zq(i)/qsl(i))*SQRT(zrfl(i)) & 637 *(paprs(i,k)-paprs(i,k+1))/pplay(i,k)*zt(i)*RD/RG 638 ENDIF 639 640 zqevt = MAX(0.0,MIN(zqevt,zrfl(i))) & 641 *RG*dtime/(paprs(i,k)-paprs(i,k+1)) 642 643 ! sublimation of the solid precipitation coming from above 644 IF (iflag_evap_prec.EQ.3) THEN 645 zqevti = znebprecip(i)*coef_sub*(1.0-zq(i)/qsi(i)) & 646 *SQRT(zifl(i)/max(1.e-4,znebprecip(i))) & 647 *(paprs(i,k)-paprs(i,k+1))/pplay(i,k)*zt(i)*RD/RG 648 ELSE IF (iflag_evap_prec.GE.4) THEN 649 zqevti = znebprecipclr(i)*coef_sub*(1.0-zq(i)/qsi(i)) & 650 *SQRT(zifl(i)/max(1.e-8,znebprecipclr(i))) & 651 *(paprs(i,k)-paprs(i,k+1))/pplay(i,k)*zt(i)*RD/RG 652 ELSE 653 zqevti = 1.*coef_sub*(1.0-zq(i)/qsi(i))*SQRT(zifl(i)) & 654 *(paprs(i,k)-paprs(i,k+1))/pplay(i,k)*zt(i)*RD/RG 655 ENDIF 656 657 zqevti = MAX(0.0,MIN(zqevti,zifl(i))) & 658 *RG*dtime/(paprs(i,k)-paprs(i,k+1)) 659 660 ! A. JAM 661 ! Evaporation limit: we ensure that the layer's fraction below 662 ! the cloud or the whole mesh (depending on iflag_evap_prec) 663 ! does not reach saturation. In this case, we 664 ! redistribute zqev0 conserving the ratio liquid/ice 665 666 IF (zqevt+zqevti.GT.zqev0) THEN 667 zqev=zqev0*zqevt/(zqevt+zqevti) 668 zqevi=zqev0*zqevti/(zqevt+zqevti) 669 ELSE 670 zqev=zqevt 671 zqevi=zqevti 672 ENDIF 673 674 !--Diagnostics 675 dqreva(i,k) = - zqev / dtime 676 dqssub(i,k) = - zqevti / dtime 677 678 ! New solid and liquid precipitation fluxes after evap and sublimation 679 zrfln(i) = Max(0.,zrfl(i) - zqev*(paprs(i,k)-paprs(i,k+1)) & 680 /RG/dtime) 681 zifln(i) = Max(0.,zifl(i) - zqevi*(paprs(i,k)-paprs(i,k+1)) & 682 /RG/dtime) 683 684 685 ! vapor, temperature, precip fluxes update 686 ! vapor is updated after evaporation/sublimation (it is increased) 687 zq(i) = zq(i) - (zrfln(i)+zifln(i)-zrfl(i)-zifl(i)) & 688 * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime 689 ! zmqc is the total condensed water in the precip flux (it is decreased) 690 zmqc(i) = zmqc(i) + (zrfln(i)+zifln(i)-zrfl(i)-zifl(i)) & 691 * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime 692 ! air and precip temperature (i.e., gridbox temperature) 693 ! is updated due to latent heat cooling 694 zt(i) = zt(i) + (zrfln(i)-zrfl(i)) & 695 * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime & 696 * RLVTT/RCPD/(1.0+RVTMP2*(zq(i)+zmqc(i))) & 697 + (zifln(i)-zifl(i)) & 698 * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime & 699 * RLSTT/RCPD/(1.0+RVTMP2*(zq(i)+zmqc(i))) 700 701 ! New values of liquid and solid precipitation 702 zrfl(i) = zrfln(i) 703 zifl(i) = zifln(i) 704 705 ! c_iso here call_reevap that updates isotopic zrfl, zifl (in inout) 706 ! due to evap + sublim 707 708 709 IF (iflag_evap_prec.GE.4) THEN 710 zrflclr(i) = zrfl(i) 711 ziflclr(i) = zifl(i) 712 IF(zrflclr(i) + ziflclr(i).LE.0) THEN 713 znebprecipclr(i) = 0.0 714 ENDIF 715 zrfl(i) = zrflclr(i) + zrflcld(i) 716 zifl(i) = ziflclr(i) + ziflcld(i) 717 ENDIF 718 719 ! c_iso duplicate for isotopes or loop on isotopes 720 721 ! Melting: 722 zmelt = ((zt(i)-RTT)/(ztfondue-RTT)) ! JYG 723 ! precip fraction that is melted 724 zmelt = MIN(MAX(zmelt,0.),1.) 725 726 ! update of rainfall and snowfall due to melting 727 IF (iflag_evap_prec.GE.4) THEN 728 zrflclr(i)=zrflclr(i)+zmelt*ziflclr(i) 729 zrflcld(i)=zrflcld(i)+zmelt*ziflcld(i) 730 zrfl(i)=zrflclr(i)+zrflcld(i) 731 ELSE 732 zrfl(i)=zrfl(i)+zmelt*zifl(i) 733 ENDIF 734 735 736 ! c_iso: melting of isotopic precipi with zmelt (no fractionation) 737 738 ! Latent heat of melting because of precipitation melting 739 ! NB: the air + precip temperature is simultaneously updated 740 zt(i)=zt(i)-zifl(i)*zmelt*(RG*dtime)/(paprs(i,k)-paprs(i,k+1)) & 741 *RLMLT/RCPD/(1.0+RVTMP2*(zq(i)+zmqc(i))) 742 743 IF (iflag_evap_prec.GE.4) THEN 744 ziflclr(i)=ziflclr(i)*(1.-zmelt) 745 ziflcld(i)=ziflcld(i)*(1.-zmelt) 746 zifl(i)=ziflclr(i)+ziflcld(i) 747 ELSE 748 zifl(i)=zifl(i)*(1.-zmelt) 749 ENDIF 750 751 ELSE 752 ! if no precip, we reinitialize the cloud fraction used for the precip to 0 753 znebprecip(i)=0. 754 755 ENDIF ! (zrfl(i)+zifl(i).GT.0.) 756 757 ENDDO ! loop on klon 758 759 ENDIF ! (iflag_evap_prec>=1) 530 CALL histprecip_precld(klon, dtime, iftop, paprs(:,k), paprs(:,k+1), pplay(:,k), & 531 zt, ztupnew, zq, zmqc, zneb, znebprecipclr, & 532 zrfl, zrflclr, zrflcld, & 533 zifl, ziflclr, ziflcld, & 534 dqreva(:,k), dqssub(:,k) & 535 ) 760 536 761 537 ENDIF ! (ok_poprecip) 762 763 ! --------------------------------------------------------------------764 ! End precip evaporation765 ! --------------------------------------------------------------------766 538 767 539 ! Calculation of qsat, L/Cp*dqsat/dT and ncoreczq counter … … 782 554 783 555 ! -------------------------------------------------------------------- 784 ! P 3> Cloud formation556 ! P2> Cloud formation 785 557 !--------------------------------------------------------------------- 786 558 ! … … 793 565 ! ------------------------------------------------------------------- 794 566 795 ! P 3.1> With the PDFs (log-normal, bigaussian)567 ! P2.1> With the PDFs (log-normal, bigaussian) 796 568 ! cloud properties calculation with the initial values of t and q 797 569 ! ---------------------------------------------------------------- … … 1038 810 ENDDO ! iter=1,iflag_fisrtilp_qsat+1 1039 811 1040 ! P 3.2> Final quantities calculation812 ! P2.2> Final quantities calculation 1041 813 ! Calculated variables: 1042 814 ! rneb : cloud fraction … … 1136 908 1137 909 ! ---------------------------------------------------------------- 1138 ! P4> Precipitation formation 910 ! P3> Precipitation processes, after cloud formation 911 ! - precipitation formation, melting/freezing 1139 912 ! ---------------------------------------------------------------- 1140 913 1141 914 !================================================================ 915 ! Flag for the new and more microphysical treatment of precipitation from Atelier Nuage (R) 1142 916 IF (ok_poprecip) THEN 1143 917 … … 1468 1242 ENDIF ! ok_poprecip 1469 1243 1470 ! End of precipitation formation1471 ! -------------------------------- 1244 ! End of precipitation processes after cloud formation 1245 ! ---------------------------------------------------- 1472 1246 1473 1247 1474 1248 !---------------------------------------------------------------------- 1475 ! P 5> Calculation of cloud condensates amount seen by radiative scheme1249 ! P4> Calculation of cloud condensates amount seen by radiative scheme 1476 1250 !---------------------------------------------------------------------- 1477 1251 … … 1506 1280 1507 1281 ! ---------------------------------------------------------------- 1508 ! P 6> Wet scavenging1282 ! P5> Wet scavenging 1509 1283 ! ---------------------------------------------------------------- 1510 1284 … … 1563 1337 1564 1338 !------------------------------------------------------------ 1565 ! P 7> write diagnostics and outputs1339 ! P6 > write diagnostics and outputs 1566 1340 !------------------------------------------------------------ 1567 1341
Note: See TracChangeset
for help on using the changeset viewer.