source: LMDZ6/trunk/libf/phylmdiso/reevap.F90 @ 5282

  SUBROUTINE reevap (klon,klev,iflag_ice_thermo,t_seri,qx, &
   &         d_t_eva,d_qx_eva)

    ! flag to include modifications to ensure energy conservation (if flag >0)
    USE add_phys_tend_mod, only : fl_cor_ebil 
#ifdef ISO
    USE infotrac_phy, ONLY: ntiso,nqtot,ivap,iliq,isol,iqWIsoPha
#ifdef ISOVERIF
    USE isotopes_verif_mod
!, ONLY: errmax,errmaxrel, iso_verif_o18_aberrant_nostop,deltaD,deltaO
    USE isotopes_mod, ONLY: iso_eau,iso_hdo,iso_o18,ridicule
#ifdef ISOTRAC
    USE isotrac_routines_mod, ONLY: iso_verif_traceur_pbidouille    
#endif
#endif 
#endif
    USE yomcst_mod_h, ONLY: RPI, RCLUM, RHPLA, RKBOL, RNAVO                   &
          , RDAY, REA, REPSM, RSIYEA, RSIDAY, ROMEGA                  &
          , R_ecc, R_peri, R_incl                                      &
          , RA, RG, R1SA                                         &
          , RSIGMA                                                     &
          , R, RMD, RMV, RD, RV, RCPD                    &
          , RMO3, RMCO2, RMC, RMCH4, RMN2O, RMCFC11, RMCFC12        &
          , RCPV, RCVD, RCVV, RKAPPA, RETV, eps_w                    &
          , RCW, RCS                                                 &
          , RLVTT, RLSTT, RLMLT, RTT, RATM                           &
          , RESTT, RALPW, RBETW, RGAMW, RALPS, RBETS, RGAMS            &
          , RALPD, RBETD, RGAMD
IMPLICIT none
    !>======================================================================

    INTEGER klon,klev,iflag_ice_thermo
    REAL, DIMENSION(klon,klev), INTENT(in) :: t_seri
    REAL, DIMENSION(klon,klev,nqtot), INTENT(in) ::     qx
    REAL, DIMENSION(klon,klev), INTENT(out) :: d_t_eva
    REAL, DIMENSION(klon,klev,nqtot), INTENT(out) ::    d_qx_eva

    REAL za,zb,zdelta,zlvdcp,zlsdcp
    INTEGER i,k,ixt,ivapcur,iliqcur,isolcur


    !--------Stochastic Boundary Layer Triggering: ALE_BL--------
    !---Propri\'et\'es du thermiques au LCL

    include "YOETHF.h"
    include "FCTTRE.h"
    !IM 100106 BEG : pouvoir sortir les ctes de la physique
    !
    DO ixt = 1, 1+ntiso
    ! Re-evaporer l'eau liquide nuageuse
    !
    iliqcur= iqWIsoPha(ixt,iliq)    
    ivapcur= iqWIsoPha(ixt,ivap)    
    isolcur= iqWIsoPha(ixt,isol)   
!print *,'rrevap ; fl_cor_ebil:',fl_cor_ebil,' iflag_ice_thermo:',iflag_ice_thermo,' RVTMP2',RVTMP2
    DO k = 1, klev  ! re-evaporation de l'eau liquide nuageuse
      DO i = 1, klon

        IF (ixt == 1) THEN ! water
          IF (fl_cor_ebil > 0) THEN
            !zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*(q_seri(i,k)+ql_seri(i,k)+qs_seri(i,k)))
            !zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*(q_seri(i,k)+ql_seri(i,k)+qs_seri(i,k)))
            zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*(qx(i,k,ivapcur)+qx(i,k,iliqcur)+qx(i,k,isolcur)))
            zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*(qx(i,k,ivapcur)+qx(i,k,iliqcur)+qx(i,k,isolcur)))
          ELSE
            zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*qx(i,k,ivapcur))
            !jyg<
            !  Attention : Arnaud a propose des formules completement differentes
            !                  A verifier !!!
            zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*qx(i,k,ivapcur))
         ENDIF
          IF (iflag_ice_thermo == 0) THEN
            zlsdcp=zlvdcp
          ENDIF
          !>jyg
        ENDIF
        IF (iflag_ice_thermo == 0) THEN   
           !pas necessaire a priori

            zdelta = MAX(0.,SIGN(1.,RTT-t_seri(i,k)))
            zdelta = 0.
            zb = MAX(0.0,qx(i,k,iliqcur))
            IF (ixt == 1) THEN
              za = - MAX(0.0,qx(i,k,iliqcur)) &
                   * (zlvdcp*(1.-zdelta)+zlsdcp*zdelta)
              d_t_eva(i,k) = za
            ENDIF
            d_qx_eva(i,k,ivapcur)  = zb
            d_qx_eva(i,k,iliqcur) = -qx(i,k,iliqcur)
            d_qx_eva(i,k,isolcur) = 0.

        ELSE
              
             !CR: on r\'e-\'evapore eau liquide et glace

             !        zdelta = MAX(0.,SIGN(1.,RTT-t_seri(i,k)))
             !        zb = MAX(0.0,ql_seri(i,k))
             !        za = - MAX(0.0,ql_seri(i,k)) &
             !             * (zlvdcp*(1.-zdelta)+zlsdcp*zdelta)
            IF (ixt == 1) THEN
              za = - MAX(0.0,qx(i,k,iliqcur))*zlvdcp & 
                   - MAX(0.0,qx(i,k,iliqcur))*zlsdcp
              d_t_eva(i,k) = za
            ENDIF
            !zb = MAX(0.0,ql_seri(i,k)+qs_seri(i,k))
            !d_q_eva(i,k) = zb
            !d_ql_eva(i,k) = -ql_seri(i,k)
            !d_qs_eva(i,k) = -qs_seri(i,k)

            zb = MAX(0.0,qx(i,k,iliqcur)+qx(i,k,isolcur))
            d_qx_eva(i,k,ivapcur) = zb
            d_qx_eva(i,k,iliqcur) = -qx(i,k,iliqcur)
            d_qx_eva(i,k,isolcur) = -qx(i,k,isolcur)
        ENDIF


      ENDDO
    ENDDO

    ENDDO ! DO ixt = 1, 1+niso*(nzone +1)

RETURN

END SUBROUTINE reevap