Changeset 5408 for LMDZ6


Ignore:
Timestamp:
Dec 15, 2024, 3:48:56 PM (6 hours ago)
Author:
evignon
Message:

debut de l'externalisation des processus de precipitation de lmdz_lscp avec ici, les processus precedents
la formation des nuages.
A Borella, E Vignon

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  
    108108USE lmdz_lscp_tools, ONLY : fallice_velocity, distance_to_cloud_top
    109109USE lmdz_lscp_condensation, ONLY : condensation_lognormal, condensation_ice_supersat
    110 USE lmdz_lscp_poprecip, ONLY : poprecip_precld, poprecip_postcld
     110USE lmdz_lscp_precip, ONLY : histprecip_precld, poprecip_precld, poprecip_postcld
    111111
    112112! Use du module lmdz_lscp_ini contenant les constantes
    113113USE 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_min
     114USE lmdz_lscp_ini, ONLY : seuil_neb, niter_lscp, iflag_evap_prec, DDT0, rain_int_min
    115115USE lmdz_lscp_ini, ONLY : ok_radocond_snow, a_tr_sca, cld_expo_con, cld_expo_lsc
    116116USE lmdz_lscp_ini, ONLY : iflag_cloudth_vert, iflag_rain_incloud_vol, iflag_t_glace, t_glace_min
    117 USE lmdz_lscp_ini, ONLY : coef_eva, coef_sub,cld_tau_lsc, cld_tau_con, cld_lc_lsc, cld_lc_con
     117USE lmdz_lscp_ini, ONLY : cld_tau_lsc, cld_tau_con, cld_lc_lsc, cld_lc_con
    118118USE lmdz_lscp_ini, ONLY : iflag_bergeron, iflag_fisrtilp_qsat, iflag_vice, cice_velo, dice_velo
    119119USE lmdz_lscp_ini, ONLY : iflag_autoconversion, ffallv_con, ffallv_lsc, min_frac_th_cld
     
    267267  ! LOCAL VARIABLES:
    268268  !----------------
    269   REAL,DIMENSION(klon) :: qsl, qsi                                ! saturation threshold at current vertical level
    270269  REAL,DIMENSION(klon) :: qice_ini
    271270  REAL :: zct, zcl,zexpo
     
    284283  REAL, DIMENSION(klon) :: zfice_th
    285284  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
    290287  REAL, DIMENSION(klon) :: zoliq, zcond, zq, zqn
    291288  REAL, DIMENSION(klon) :: zoliql, zoliqi
     
    294291  REAL, DIMENSION(klon) :: zdz,iwc
    295292  REAL :: zchau,zfroi
    296   REAL, DIMENSION(klon) :: zfice,zneb,znebprecip
    297   REAL :: zmelt,zrain,zsnow,zprecip
     293  REAL, DIMENSION(klon) :: zfice,zneb
     294  REAL :: zrain,zsnow,zprecip
    298295  REAL, DIMENSION(klon) :: dzfice
    299296  REAL, DIMENSION(klon) :: zfice_turb, dzfice_turb
    300297  REAL :: zsolid
    301298  REAL, DIMENSION(klon) :: qtot, qzero
    302   REAL, DIMENSION(klon) :: dqsl,dqsi
    303299  REAL :: smallestreal
    304300  !  Variables for Bergeron process
     
    310306  REAL :: zfrac_lessi
    311307  REAL, DIMENSION(klon) :: zprec_cond
    312   REAL :: zmair
    313   REAL :: zcpair, zcpeau
    314308  REAL, DIMENSION(klon) :: zlh_solid
    315309  REAL, DIMENSION(klon) :: ztupnew
     
    318312  REAL, DIMENSION(klon) :: zrflclr, zrflcld
    319313  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
    321315  REAL, DIMENSION(klon) :: ziflclr, ziflcld
    322316  REAL, DIMENSION(klon) :: znebprecipclr, znebprecipcld
     
    369363! Few initialisations
    370364
    371 znebprecip(:)=0.0
    372365ctot_vol(1:klon,1:klev)=0.0
    373366rneblsvol(1:klon,1:klev)=0.0
     
    488481DO k = klev, 1, -1
    489482
     483    qice_ini = qi_seri(:,k)
     484
    490485    IF (k.LE.klev-1) THEN
    491486       iftop=.false.
     
    511506        !c_iso init of iso
    512507    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    !---------------------------------------------------------------------
    513514   
    514515    !================================================================
    515516    ! Flag for the new and more microphysical treatment of precipitation from Atelier Nuage (R)
    516     IF (ok_poprecip) THEN
    517 
    518             CALL poprecip_precld(klon, dtime, iftop, paprs(:,k), paprs(:,k+1), pplay(:,k), &
    519                               zt, ztupnew, zq, zmqc, znebprecipclr, znebprecipcld, &
    520                               zqvapclr, zqupnew, &
    521                               zrfl, zrflclr, zrflcld, &
    522                               zifl, ziflclr, ziflcld, &
    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                        )
    525526
    526527    !================================================================
    527528    ELSE
    528529
    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                        )
    760536
    761537    ENDIF ! (ok_poprecip)
    762 
    763     ! --------------------------------------------------------------------
    764     ! End precip evaporation
    765     ! --------------------------------------------------------------------
    766538   
    767539    ! Calculation of qsat, L/Cp*dqsat/dT and ncoreczq counter
     
    782554   
    783555    ! --------------------------------------------------------------------
    784     ! P3> Cloud formation
     556    ! P2> Cloud formation
    785557    !---------------------------------------------------------------------
    786558    !
     
    793565    ! -------------------------------------------------------------------
    794566
    795         ! P3.1> With the PDFs (log-normal, bigaussian)
     567        ! P2.1> With the PDFs (log-normal, bigaussian)
    796568        ! cloud properties calculation with the initial values of t and q
    797569        ! ----------------------------------------------------------------
     
    1038810        ENDDO         ! iter=1,iflag_fisrtilp_qsat+1
    1039811
    1040         ! P3.2> Final quantities calculation
     812        ! P2.2> Final quantities calculation
    1041813        ! Calculated variables:
    1042814        !   rneb : cloud fraction
     
    1136908
    1137909    ! ----------------------------------------------------------------
    1138     ! P4> Precipitation formation
     910    ! P3> Precipitation processes, after cloud formation
     911    !     - precipitation formation, melting/freezing
    1139912    ! ----------------------------------------------------------------
    1140913
    1141914    !================================================================
     915    ! Flag for the new and more microphysical treatment of precipitation from Atelier Nuage (R)
    1142916    IF (ok_poprecip) THEN
    1143917
     
    14681242    ENDIF ! ok_poprecip
    14691243
    1470     ! End of precipitation formation
    1471     ! --------------------------------
     1244    ! End of precipitation processes after cloud formation
     1245    ! ----------------------------------------------------
    14721246
    14731247
    14741248    !----------------------------------------------------------------------
    1475     ! P5> Calculation of cloud condensates amount seen by radiative scheme
     1249    ! P4> Calculation of cloud condensates amount seen by radiative scheme
    14761250    !----------------------------------------------------------------------
    14771251
     
    15061280
    15071281    ! ----------------------------------------------------------------
    1508     ! P6> Wet scavenging
     1282    ! P5> Wet scavenging
    15091283    ! ----------------------------------------------------------------
    15101284
     
    15631337   
    15641338    !------------------------------------------------------------
    1565     ! P7 > write diagnostics and outputs
     1339    ! P6 > write diagnostics and outputs
    15661340    !------------------------------------------------------------
    15671341   
Note: See TracChangeset for help on using the changeset viewer.