Skip to content
Snippets Groups Projects
rain_ice.f90 145 KiB
Newer Older
  • Learn to ignore specific revisions
  • !
     IF ( KRR == 7 ) THEN
    
         PRHS(:,:,:) = PRHS(:,:,:) * PTSTEP
    
         ZWSED(:,:,:) = 0.
         ZWSEDW1(:,:,:) = 0.
         ZWSEDW2(:,:,:) = 0.
    ! calculation of ZP1, ZP2 and sedimentation flux
         DO JK = IKE , IKB, -1*KKL
         !estimation of q' taking into account incomming ZWSED
         ZQP(:,:)=ZWSED(:,:,JK+KKL)*ZW(:,:,JK)
    
         JCOUNT=COUNTJV2((PRHS(:,:,JK)+ZQP(JI,JJ) > ZRTMIN(7)) .OR. &
                         (ZQP(:,:) > ZRTMIN(7)),I1(:),I2(:))
         DO JL=1, JCOUNT
           JI=I1(JL)
           JJ=I2(JL)
             !calculation of w
             IF ((PRHS(JI,JJ,JK)+ZQP(JI,JJ)) > ZRTMIN(7) ) THEN
               ZWSEDW1 (JI,JJ,JK)= XFSEDH  * (PRHS(JI,JJ,JK))**(XEXSEDH-1) *   &
                                        PRHODREF(JI,JJ,JK)**(XEXSEDH-XCEXVT-1)
             ENDIF
             IF ( ZQP(JI,JJ) > ZRTMIN(7) ) THEN
               ZWSEDW2 (JI,JJ,JK)= XFSEDH  * ZQP(JI,JJ)**(XEXSEDH-1) *   &
                                        PRHODREF(JI,JJ,JK)**(XEXSEDH-XCEXVT-1)
             ENDIF
           ENDDO
           DO JJ = IJB, IJE
             DO JI = IIB, IIE
               ZH=PDZZ(JI,JJ,JK)
               ZP1 = MIN(1., ZWSEDW1(JI,JJ,JK) * PTSTEP / ZH)
               IF (ZWSEDW2(JI,JJ,JK) /= 0.) THEN
                 ZP2 = MAX(0.,1 - ZH &
               &  / (PTSTEP*ZWSEDW2(JI,JJ,JK)) )
               ELSE
                 ZP2 = 0.
               ENDIF
               ZWSED (JI,JJ,JK)=ZP1*PRHODREF(JI,JJ,JK)*&
               &ZH*PRHS(JI,JJ,JK)&
               &* ZINVTSTEP+ ZP2 * ZWSED (JI,JJ,JK+KKL)
             ENDDO
           ENDDO
         ENDDO
    
         DO JK = IKTB , IKTE
           PRHS(:,:,JK) = PRHS(:,:,JK) + ZW(:,:,JK)*(ZWSED(:,:,JK+KKL)-ZWSED(:,:,JK))
         ENDDO
    
         IF (PRESENT(PFPR)) THEN
           DO JK = IKTB , IKTE
             PFPR(:,:,JK,7)=ZWSED(:,:,JK)
           ENDDO
         ENDIF
    
    
         PINPRH(:,:) = ZWSED(:,:,IKB)/XRHOLW                        ! in m/s
    
         PRHS(:,:,:) = PRHS(:,:,:) * ZINVTSTEP
    
     ENDIF
    !
    
    !
    !*       2.3     budget storage
    !
    IF (LBUDGET_RC .AND. OSEDIC) &
                    CALL BUDGET (PRCS(:,:,:)*PRHODJ(:,:,:),7 ,'SEDI_BU_RRC')
    IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:)*PRHODJ(:,:,:),8 ,'SEDI_BU_RRR')
    IF (LBUDGET_RI) CALL BUDGET (PRIS(:,:,:)*PRHODJ(:,:,:),9 ,'SEDI_BU_RRI')
    IF (LBUDGET_RS) CALL BUDGET (PRSS(:,:,:)*PRHODJ(:,:,:),10,'SEDI_BU_RRS')
    IF (LBUDGET_RG) CALL BUDGET (PRGS(:,:,:)*PRHODJ(:,:,:),11,'SEDI_BU_RRG')
    
    IF ( KRR == 7 .AND. LBUDGET_RH) &
    
                    CALL BUDGET (PRHS(:,:,:)*PRHODJ(:,:,:),12,'SEDI_BU_RRH')
    !
    
    !
    !*       2.4  DROPLET DEPOSITION AT THE 1ST LEVEL ABOVE GROUND
    !
    IF (LDEPOSC) THEN
      GDEP(:,:) = .FALSE.
    
      GDEP(IIB:IIE,IJB:IJE) =    PRCS(IIB:IIE,IJB:IJE,IKB) >0 
    
         PRCS(:,:,IKB) = PRCS(:,:,IKB) - XVDEPOSC * PRCT(:,:,IKB) / PDZZ(:,:,IKB)
         PINPRC(:,:) = PINPRC(:,:) + XVDEPOSC * PRCT(:,:,IKB) * PRHODREF(:,:,IKB) /XRHOLW 
         PINDEP(:,:) = XVDEPOSC * PRCT(:,:,IKB) * PRHODREF(:,:,IKB) /XRHOLW 
    
      END WHERE
    END IF
    !
    !*       2.5     budget storage
    !
    IF ( LBUDGET_RC .AND. LDEPOSC ) &
       CALL BUDGET (PRCS(:,:,:)*PRHODJ(:,:,:),7 ,'DEPO_BU_RRC')
    
    
      END SUBROUTINE RAIN_ICE_SEDIMENTATION_STAT
    !
    !-------------------------------------------------------------------------------
    !
    
    !
      SUBROUTINE RAIN_ICE_NUCLEATION
    !
    !*      0. DECLARATIONS
    !          ------------
    !
    IMPLICIT NONE
    !
    !*       0.2  declaration of local variables
    !
    INTEGER , DIMENSION(SIZE(GNEGT))  :: I1,I2,I3 ! Used to replace the COUNT
    INTEGER                           :: JL       ! and PACK intrinsics
    !
    !-------------------------------------------------------------------------------
    !
    !
    !  compute the temperature and the pressure
    !
    ZT(:,:,:) = PTHT(:,:,:) * ( PPABST(:,:,:) / XP00 ) ** (XRD/XCPD)
    !
    !  optimization by looking for locations where
    !  the temperature is negative only !!!
    !
    GNEGT(:,:,:) = .FALSE.
    GNEGT(IIB:IIE,IJB:IJE,IKTB:IKTE) = ZT(IIB:IIE,IJB:IJE,IKTB:IKTE)<XTT
    INEGT = COUNTJV( GNEGT(:,:,:),I1(:),I2(:),I3(:))
    IF( INEGT >= 1 ) THEN
      ALLOCATE(ZRVT(INEGT)) ;
      ALLOCATE(ZCIT(INEGT)) ;
      ALLOCATE(ZZT(INEGT))  ;
      ALLOCATE(ZPRES(INEGT));
      DO JL=1,INEGT
        ZRVT(JL) = PRVT(I1(JL),I2(JL),I3(JL))
        ZCIT(JL) = PCIT(I1(JL),I2(JL),I3(JL))
        ZZT(JL) = ZT(I1(JL),I2(JL),I3(JL))
        ZPRES(JL) = PPABST(I1(JL),I2(JL),I3(JL))
      ENDDO
      ALLOCATE(ZZW(INEGT))
      ALLOCATE(ZUSW(INEGT))
      ALLOCATE(ZSSI(INEGT))
        ZZW(:) = EXP( XALPI - XBETAI/ZZT(:) - XGAMI*ALOG(ZZT(:) ) )           ! es_i
        ZZW(:) = MIN(ZPRES(:)/2., ZZW(:))             ! safety limitation
        ZSSI(:) = ZRVT(:)*( ZPRES(:)-ZZW(:) ) / ( (XMV/XMD) * ZZW(:) ) - 1.0
                                                          ! Supersaturation over ice
        ZUSW(:) = EXP( XALPW - XBETAW/ZZT(:) - XGAMW*ALOG(ZZT(:) ) )          ! es_w
        ZUSW(:) = MIN(ZPRES(:)/2.,ZUSW(:))            ! safety limitation
        ZUSW(:) = ( ZUSW(:)/ZZW(:) )*( (ZPRES(:)-ZZW(:))/(ZPRES(:)-ZUSW(:)) ) - 1.0
                                 ! Supersaturation of saturated water vapor over ice
    !
    !*       3.1     compute the heterogeneous nucleation source: RVHENI
    !
    !*       3.1.1   compute the cloud ice concentration
    !
      ZZW(:) = 0.0
      ZSSI(:) = MIN( ZSSI(:), ZUSW(:) ) ! limitation of SSi according to SSw=0
      WHERE( (ZZT(:)<XTT-5.0) .AND. (ZSSI(:)>0.0) )
        ZZW(:) = XNU20 * EXP( XALPHA2*ZSSI(:)-XBETA2 )
      END WHERE
      WHERE( (ZZT(:)<=XTT-2.0) .AND. (ZZT(:)>=XTT-5.0) .AND. (ZSSI(:)>0.0) )
        ZZW(:) = MAX( XNU20 * EXP( -XBETA2 ),XNU10 * EXP( -XBETA1*(ZZT(:)-XTT) ) * &
                                   ( ZSSI(:)/ZUSW(:) )**XALPHA1 )
      END WHERE
      ZZW(:) = ZZW(:) - ZCIT(:)
      IF( MAXVAL(ZZW(:)) > 0.0 ) THEN
    !
    !*       3.1.2   update the r_i and r_v mixing ratios
    !
        ZZW(:) = MIN( ZZW(:),50.E3 ) ! limitation provisoire a 50 l^-1
        ZW(:,:,:) = UNPACK( ZZW(:),MASK=GNEGT(:,:,:),FIELD=0.0 )
        ZW(:,:,:) = MAX( ZW(:,:,:) ,0.0 ) *XMNU0/(PRHODREF(:,:,:)*PTSTEP)
        PRIS(:,:,:) = PRIS(:,:,:) + ZW(:,:,:)
        PRVS(:,:,:) = PRVS(:,:,:) - ZW(:,:,:)
        IF ( KRR == 7 ) THEN
            PTHS(:,:,:) = PTHS(:,:,:) + ZW(:,:,:)*(XLSTT+(XCPV-XCI)*(ZT(:,:,:)-XTT))   &
                     /( (XCPD + XCPV*PRVT(:,:,:) + XCL*(PRCT(:,:,:)+PRRT(:,:,:))   &
           + XCI*(PRIT(:,:,:)+PRST(:,:,:)+PRGT(:,:,:)+PRHT(:,:,:)))*PEXNREF(:,:,:) )
          ELSE IF( KRR == 6 ) THEN
            PTHS(:,:,:) = PTHS(:,:,:) + ZW(:,:,:)*(XLSTT+(XCPV-XCI)*(ZT(:,:,:)-XTT))   &
                     /( (XCPD + XCPV*PRVT(:,:,:) + XCL*(PRCT(:,:,:)+PRRT(:,:,:))   &
                       + XCI*(PRIT(:,:,:)+PRST(:,:,:)+PRGT(:,:,:)))*PEXNREF(:,:,:) )
        END IF
    
        ZZW(:) = MAX( ZZW(:)+ZCIT(:),ZCIT(:) )
        PCIT(:,:,:) = MAX( UNPACK( ZZW(:),MASK=GNEGT(:,:,:),FIELD=0.0 ) , &
                           PCIT(:,:,:) )
      END IF
      DEALLOCATE(ZSSI)
      DEALLOCATE(ZUSW)
      DEALLOCATE(ZZW)
      DEALLOCATE(ZPRES)
      DEALLOCATE(ZZT)
      DEALLOCATE(ZCIT)
      DEALLOCATE(ZRVT)
    END IF
    !
    !*       3.1.3   budget storage
    !
    IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'HENU_BU_RTH')
    IF (LBUDGET_RV) CALL BUDGET (PRVS(:,:,:)*PRHODJ(:,:,:),6,'HENU_BU_RRV')
    IF (LBUDGET_RI) CALL BUDGET (PRIS(:,:,:)*PRHODJ(:,:,:),9,'HENU_BU_RRI')
    !
      END SUBROUTINE RAIN_ICE_NUCLEATION
    !
    !-------------------------------------------------------------------------------
    !
    !
      SUBROUTINE RAIN_ICE_SLOW
    !
    !*      0. DECLARATIONS
    !          ------------
    !
    IMPLICIT NONE
    !
    !-------------------------------------------------------------------------------
    !
    !
    !*       3.2     compute the homogeneous nucleation source: RCHONI
    !
      ZZW(:) = 0.0
      WHERE( (ZZT(:)<XTT-35.0) .AND. (ZRCT(:)>XRTMIN(2)) .AND. (ZRCS(:)>0.) )
        ZZW(:) = MIN( ZRCS(:),XHON*ZRHODREF(:)*ZRCT(:)       &
    
                                     *EXP( MIN(XMNH_HUGE_12_LOG,XALPHA3*(ZZT(:)-XTT)-XBETA3) ) )
                                     ! *EXP( XALPHA3*(ZZT(:)-XTT)-XBETA3 ) )
    
        ZRIS(:) = ZRIS(:) + ZZW(:)
        ZRCS(:) = ZRCS(:) - ZZW(:)
        ZTHS(:) = ZTHS(:) + ZZW(:)*(ZLSFACT(:)-ZLVFACT(:)) ! f(L_f*(RCHONI))
      ENDWHERE
    !
      IF (LBUDGET_TH) CALL BUDGET (                                                &
                     UNPACK(ZTHS(:),MASK=GMICRO(:,:,:),FIELD=PTHS)*PRHODJ(:,:,:),  &
                                                                  4,'HON_BU_RTH')
    
      IF (LBUDGET_RC) CALL BUDGET (                                                   &
                         UNPACK(ZRCS(:),MASK=GMICRO(:,:,:),FIELD=PRCS)*PRHODJ(:,:,:), &
    
                                                                  7,'HON_BU_RRC')
    
      IF (LBUDGET_RI) CALL BUDGET (                                                   &
                         UNPACK(ZRIS(:),MASK=GMICRO(:,:,:),FIELD=PRIS)*PRHODJ(:,:,:), &
    
                                                                  9,'HON_BU_RRI')
    !
    !*       3.3     compute the spontaneous freezing source: RRHONG
    !
      ZZW(:) = 0.0
      WHERE( (ZZT(:)<XTT-35.0) .AND. (ZRRT(:)>XRTMIN(3)) .AND. (ZRRS(:)>0.) )
        ZZW(:) = MIN( ZRRS(:),ZRRT(:)* ZINVTSTEP )
        ZRGS(:) = ZRGS(:) + ZZW(:)
        ZRRS(:) = ZRRS(:) - ZZW(:)
        ZTHS(:) = ZTHS(:) + ZZW(:)*(ZLSFACT(:)-ZLVFACT(:)) ! f(L_f*(RRHONG))
      ENDWHERE
    !
      IF (LBUDGET_TH) CALL BUDGET (                                                &
                     UNPACK(ZTHS(:),MASK=GMICRO(:,:,:),FIELD=PTHS)*PRHODJ(:,:,:),  &
                                                                  4,'SFR_BU_RTH')
    
      IF (LBUDGET_RR) CALL BUDGET (                                                   &
                         UNPACK(ZRRS(:),MASK=GMICRO(:,:,:),FIELD=PRRS)*PRHODJ(:,:,:), &
    
                                                                  8,'SFR_BU_RRR')
    
      IF (LBUDGET_RG) CALL BUDGET (                                                   &
                         UNPACK(ZRGS(:),MASK=GMICRO(:,:,:),FIELD=PRGS)*PRHODJ(:,:,:), &
    
                                                                 11,'SFR_BU_RRG')
    !
    !*       3.4    compute the deposition, aggregation and autoconversion sources
    !
      ZKA(:) = 2.38E-2 + 0.0071E-2 * ( ZZT(:) - XTT )          ! k_a
      ZDV(:) = 0.211E-4 * (ZZT(:)/XTT)**1.94 * (XP00/ZPRES(:)) ! D_v
    !
    !*       3.4.1  compute the thermodynamical function A_i(T,P)
    !*              and the c^prime_j (in the ventilation factor)
    !
    
      ZAI(:) = EXP( XALPI - XBETAI/ZZT(:) - XGAMI*ALOG(ZZT(:) ) ) ! es_i
      ZAI(:) = ( XLSTT + (XCPV-XCI)*(ZZT(:)-XTT) )**2 / (ZKA(:)*XRV*ZZT(:)**2) &
                                     + ( XRV*ZZT(:) ) / (ZDV(:)*ZAI(:))
      ZCJ(:) = XSCFAC * ZRHODREF(:)**0.3 / SQRT( 1.718E-5+0.0049E-5*(ZZT(:)-XTT) )
    !
    !*       3.4.2  compute the riming-conversion of r_c for r_i production: RCAUTI
    !
    !  ZZW(:) = 0.0
    !  ZTIMAUTIC = SQRT( XTIMAUTI*XTIMAUTC )
    !  WHERE ( (ZRCT(:)>0.0) .AND. (ZRIT(:)>0.0) .AND. (ZRCS(:)>0.0) )
    !    ZZW(:) = MIN( ZRCS(:),ZTIMAUTIC * MAX( SQRT( ZRIT(:)*ZRCT(:) ),0.0 ) )
    !    ZRIS(:) = ZRIS(:) + ZZW(:)
    !    ZRCS(:) = ZRCS(:) - ZZW(:)
    !    ZTHS(:) = ZTHS(:) + ZZW(:)*(ZLSFACT(:)-ZLVFACT(:)) ! f(L_f*(RCAUTI))
    !  END WHERE
    !
    !*       3.4.3  compute the deposition on r_s: RVDEPS
    !
      WHERE ( ZRST(:)>0.0 )
        ZLBDAS(:)  = MIN( XLBDAS_MAX,                                           &
                          XLBS*( ZRHODREF(:)*MAX( ZRST(:),XRTMIN(5) ) )**XLBEXS )
      END WHERE
      ZZW(:) = 0.0
      WHERE ( (ZRST(:)>XRTMIN(5)) .AND. (ZRSS(:)>0.0) )
        ZZW(:) = ( ZSSI(:)/(ZRHODREF(:)*ZAI(:)) ) *                               &
                 ( X0DEPS*ZLBDAS(:)**XEX0DEPS + X1DEPS*ZCJ(:)*ZLBDAS(:)**XEX1DEPS )
        ZZW(:) =         MIN( ZRVS(:),ZZW(:)      )*(0.5+SIGN(0.5,ZZW(:))) &
                       - MIN( ZRSS(:),ABS(ZZW(:)) )*(0.5-SIGN(0.5,ZZW(:)))
        ZRSS(:) = ZRSS(:) + ZZW(:)
        ZRVS(:) = ZRVS(:) - ZZW(:)
        ZTHS(:) = ZTHS(:) + ZZW(:)*ZLSFACT(:)
      END WHERE
      IF (LBUDGET_TH) CALL BUDGET (                                                &
                     UNPACK(ZTHS(:),MASK=GMICRO(:,:,:),FIELD=PTHS)*PRHODJ(:,:,:),  &
                                                                  4,'DEPS_BU_RTH')
      IF (LBUDGET_RV) CALL BUDGET (                                                &
                     UNPACK(ZRVS(:),MASK=GMICRO(:,:,:),FIELD=PRVS)*PRHODJ(:,:,:),  &
                                                                  6,'DEPS_BU_RRV')
    
      IF (LBUDGET_RS) CALL BUDGET (                                                   &
                         UNPACK(ZRSS(:),MASK=GMICRO(:,:,:),FIELD=PRSS)*PRHODJ(:,:,:), &
    
                                                                 10,'DEPS_BU_RRS')
    !
    !*       3.4.4  compute the aggregation on r_s: RIAGGS
    !
      ZZW(:) = 0.0
      WHERE ( (ZRIT(:)>XRTMIN(4)) .AND. (ZRST(:)>XRTMIN(5)) .AND. (ZRIS(:)>0.0) )
        ZZW(:) = MIN( ZRIS(:),XFIAGGS * EXP( XCOLEXIS*(ZZT(:)-XTT) ) &
    
                                      * ZRIT(:)                      &
                                      * ZLBDAS(:)**XEXIAGGS          &
                                      * ZRHODREF(:)**(-XCEXVT)       )
    
        ZRSS(:)  = ZRSS(:)  + ZZW(:)
        ZRIS(:)  = ZRIS(:)  - ZZW(:)
      END WHERE
    
      IF (LBUDGET_RI) CALL BUDGET (                                                   &
                         UNPACK(ZRIS(:),MASK=GMICRO(:,:,:),FIELD=PRIS)*PRHODJ(:,:,:), &
    
                                                                  9,'AGGS_BU_RRI')
    
      IF (LBUDGET_RS) CALL BUDGET (                                                   &
                         UNPACK(ZRSS(:),MASK=GMICRO(:,:,:),FIELD=PRSS)*PRHODJ(:,:,:), &
    
                                                                 10,'AGGS_BU_RRS')
    !
    !*       3.4.5  compute the autoconversion of r_i for r_s production: RIAUTS
    !
      ALLOCATE(ZCRIAUTI(IMICRO))
    
    !  ZCRIAUTI(:)=MIN(XCRIAUTI,10**(0.06*(ZZT(:)-XTT)-3.5))
    
      ZCRIAUTI(:)=MIN(XCRIAUTI,10**(XACRIAUTI*(ZZT(:)-XTT)+XBCRIAUTI))
    
      ZZW(:) = 0.0
      WHERE ( (ZRIT(:)>XRTMIN(4)) .AND. (ZRIS(:)>0.0) )
        ZZW(:) = MIN( ZRIS(:),XTIMAUTI * EXP( XTEXAUTI*(ZZT(:)-XTT) ) &
                                * MAX( ZRIT(:)-ZCRIAUTI(:),0.0 ) )
        ZRSS(:)  = ZRSS(:)  + ZZW(:)
        ZRIS(:)  = ZRIS(:)  - ZZW(:)
      END WHERE
      DEALLOCATE(ZCRIAUTI)
    
      IF (LBUDGET_RI) CALL BUDGET (                                                   &
                         UNPACK(ZRIS(:),MASK=GMICRO(:,:,:),FIELD=PRIS)*PRHODJ(:,:,:), &
    
                                                                  9,'AUTS_BU_RRI')
    
      IF (LBUDGET_RS) CALL BUDGET (                                                   &
                         UNPACK(ZRSS(:),MASK=GMICRO(:,:,:),FIELD=PRSS)*PRHODJ(:,:,:), &
    
                                                                 10,'AUTS_BU_RRS')
    !
    !*       3.4.6  compute the deposition on r_g: RVDEPG
    !
    !
      WHERE ( ZRGT(:)>0.0 )
        ZLBDAG(:)  = XLBG*( ZRHODREF(:)*MAX( ZRGT(:),XRTMIN(6) ) )**XLBEXG
      END WHERE
      ZZW(:) = 0.0
      WHERE ( (ZRGT(:)>XRTMIN(6)) .AND. (ZRGS(:)>0.0) )
        ZZW(:) = ( ZSSI(:)/(ZRHODREF(:)*ZAI(:)) ) *                               &
                 ( X0DEPG*ZLBDAG(:)**XEX0DEPG + X1DEPG*ZCJ(:)*ZLBDAG(:)**XEX1DEPG )
        ZZW(:) =         MIN( ZRVS(:),ZZW(:)      )*(0.5+SIGN(0.5,ZZW(:))) &
                       - MIN( ZRGS(:),ABS(ZZW(:)) )*(0.5-SIGN(0.5,ZZW(:)))
        ZRGS(:) = ZRGS(:) + ZZW(:)
        ZRVS(:) = ZRVS(:) - ZZW(:)
        ZTHS(:) = ZTHS(:) + ZZW(:)*ZLSFACT(:)
      END WHERE
      IF (LBUDGET_TH) CALL BUDGET (                                                 &
                     UNPACK(ZTHS(:),MASK=GMICRO(:,:,:),FIELD=PTHS)*PRHODJ(:,:,:),   &
                                                                  4,'DEPG_BU_RTH')
      IF (LBUDGET_RV) CALL BUDGET (                                                 &
                     UNPACK(ZRVS(:),MASK=GMICRO(:,:,:),FIELD=PRVS)*PRHODJ(:,:,:),   &
                                                                  6,'DEPG_BU_RRV')
    
      IF (LBUDGET_RG) CALL BUDGET (                                                   &
                         UNPACK(ZRGS(:),MASK=GMICRO(:,:,:),FIELD=PRGS)*PRHODJ(:,:,:), &
    
                                                                 11,'DEPG_BU_RRG')
    !
      END SUBROUTINE RAIN_ICE_SLOW
    !
    !-------------------------------------------------------------------------------
    !
    !
      SUBROUTINE RAIN_ICE_WARM
    !
    !*      0. DECLARATIONS
    !          ------------
    !
    IMPLICIT NONE
    !
    !
    !-------------------------------------------------------------------------------
    !
    !*       4.2    compute the autoconversion of r_c for r_r production: RCAUTR
    !
    
    
        WHERE( ZRCS(:)>0.0 .AND. ZHLC_HCF(:).GT.0.0 )
          ZZW(:) = XTIMAUTC*MAX( ZHLC_HRC(:)/ZHLC_HCF(:)  - XCRIAUTC/ZRHODREF(:),0.0)
          ZZW(:) = MIN( ZRCS(:),ZHLC_HCF(:)*ZZW(:))
          ZRCS(:) = ZRCS(:) - ZZW(:)
          ZRRS(:) = ZRRS(:) + ZZW(:)
        END WHERE
    
        IF (LBUDGET_RC) CALL BUDGET (                                                   &
                           UNPACK(ZRCS(:),MASK=GMICRO(:,:,:),FIELD=PRCS)*PRHODJ(:,:,:), &
                                                                  7,'AUTO_BU_RRC')
        IF (LBUDGET_RR) CALL BUDGET (                                                   &
                           UNPACK(ZRRS(:),MASK=GMICRO(:,:,:),FIELD=PRRS)*PRHODJ(:,:,:), &
    
                                                                    8,'AUTO_BU_RRR')
    !
    !*       4.3    compute the accretion of r_c for r_r production: RCACCR
    !
    
        IF (CSUBG_RC_RR_ACCR=='NONE') THEN
          !CLoud water and rain are diluted over the grid box
          WHERE( ZRCT(:)>XRTMIN(2) .AND. ZRRT(:)>XRTMIN(3) .AND. ZRCS(:)>0.0 )
            ZZW(:) = MIN( ZRCS(:), XFCACCR * ZRCT(:)                &
                     * ZLBDAR(:)**XEXCACCR    &
                     * ZRHODREF(:)**(-XCEXVT) )
            ZRCS(:) = ZRCS(:) - ZZW(:)
            ZRRS(:) = ZRRS(:) + ZZW(:)
          END WHERE
    
        ELSEIF (CSUBG_RC_RR_ACCR=='PRFR') THEN
          !Cloud water is concentrated over its fraction with possibly to parts with high and low content as set for autoconversion
          !Rain is concnetrated over its fraction
          !Rain in high content area fraction: ZHLC_HCF
          !Rain in low content area fraction:
          ! if ZRF<ZCF (rain is entirely falling in cloud): ZRF-ZHLC_HCF
          ! if ZRF>ZCF (rain is falling in cloud and in clear sky): ZCF-ZHLC_HCF
          ! => min(ZCF, ZRF)-ZHLC_HCF
          ZZW(:) = 0.
          WHERE( ZHLC_HRC(:)>XRTMIN(2) .AND. ZRRT(:)>XRTMIN(3) .AND. ZRCS(:)>0.0 &
                .AND. ZHLC_HCF(:)>0 )
            !Accretion due to rain falling in high cloud content
            ZZW(:) = XFCACCR * ( ZHLC_HRC(:)/ZHLC_HCF(:) )     &
                   * ZLBDAR_RF(:)**XEXCACCR &
                   * ZRHODREF(:)**(-XCEXVT) &
                   * ZHLC_HCF
          END WHERE
          WHERE( ZHLC_LRC(:)>XRTMIN(2) .AND. ZRRT(:)>XRTMIN(3) .AND. ZRCS(:)>0.0 &
                .AND. ZHLC_LCF(:)>0 )
            !We add acrretion due to rain falling in low cloud content
            ZZW(:) = ZZW(:) + XFCACCR * ( ZHLC_LRC(:)/ZHLC_LCF(:) )     &
                            * ZLBDAR_RF(:)**XEXCACCR &
                            * ZRHODREF(:)**(-XCEXVT) &
                            * (MIN(ZCF(:), ZRF(:))-ZHLC_HCF(:))
          END WHERE
          ZZW(:)=MIN(ZRCS(:), ZZW(:))
    
          ZRCS(:) = ZRCS(:) - ZZW(:)
          ZRRS(:) = ZRRS(:) + ZZW(:)
    
    
        ELSE
          !wrong CSUBG_RC_RR_ACCR case
          WRITE(*,*) 'wrong CSUBG_RC_RR_ACCR case'
                CALL PRINT_MSG(NVERB_FATAL,'GEN','RAIN_ICE_WARM','')  
        ENDIF
    
    
        IF (LBUDGET_RC) CALL BUDGET (                                                   &
                           UNPACK(ZRCS(:),MASK=GMICRO(:,:,:),FIELD=PRCS)*PRHODJ(:,:,:), &
                                                                    7,'ACCR_BU_RRC')
        IF (LBUDGET_RR) CALL BUDGET (                                                   &
                           UNPACK(ZRRS(:),MASK=GMICRO(:,:,:),FIELD=PRRS)*PRHODJ(:,:,:), &
                                                                    8,'ACCR_BU_RRR')
    
    !
    !*       4.4    compute the evaporation of r_r: RREVAV
    !
        ZZW(:) = 0.0
    
    
        IF (CSUBG_RR_EVAP=='NONE') THEN
          !Evaporation only when there's no cloud (RC must be 0)
           WHERE( (ZRRT(:)>XRTMIN(3)) .AND. (ZRCT(:)<=XRTMIN(2)) )
              ZZW(:)  = EXP( XALPW - XBETAW/ZZT(:) - XGAMW*ALOG(ZZT(:) ) ) ! es_w
              ZUSW(:) = 1.0 - ZRVT(:)*( ZPRES(:)-ZZW(:) ) / ( (XMV/XMD) * ZZW(:) )
                                                            ! Undersaturation over water
              ZZW(:) = ( XLVTT+(XCPV-XCL)*(ZZT(:)-XTT) )**2 / ( ZKA(:)*XRV*ZZT(:)**2 ) &
                   + ( XRV*ZZT(:) ) / ( ZDV(:)*ZZW(:) )
              ZZW(:) = MIN( ZRRS(:),( MAX( 0.0,ZUSW(:) )/(ZRHODREF(:)*ZZW(:)) ) *      &
                ( X0EVAR*ZLBDAR(:)**XEX0EVAR+X1EVAR*ZCJ(:)*ZLBDAR(:)**XEX1EVAR ) )
              ZRRS(:) = ZRRS(:) - ZZW(:)
              ZRVS(:) = ZRVS(:) + ZZW(:)
              ZTHS(:) = ZTHS(:) - ZZW(:)*ZLVFACT(:)
           END WHERE
    
        ELSEIF (CSUBG_RR_EVAP=='CLFR' .OR. CSUBG_RR_EVAP=='PRFR') THEN
          !Evaporation in clear sky part
          !With CLFR, rain is diluted over the grid box
          !With PRFR, rain is concentrated in its fraction
          !Use temperature and humidity in clear sky part like Bechtold et al. (1993)
          IF (CSUBG_RR_EVAP=='CLFR') THEN
            ZZW4(:)=1. !Precipitation fraction
            ZZW3(:)=ZLBDAR(:)
          ELSE
            ZZW4(:)=ZRF(:) !Precipitation fraction
            ZZW3(:)=ZLBDAR_RF(:)
          ENDIF
    
          !ATTENTION
          !Il faudrait recalculer les variables ZKA, ZDV, ZCJ en tenant compte de la température T^u
          !Ces variables devraient être sorties de rain_ice_slow et on mettrait le calcul de T^u, T^s
          !et plusieurs versions (comme actuellement, en ciel clair, en ciel nuageux) de ZKA, ZDV, ZCJ dans rain_ice
          !On utiliserait la bonne version suivant l'option NONE, CLFR... dans l'évaporation et ailleurs
    
          WHERE(  (ZRRT(:)>XRTMIN(3)) .AND. ( ZZW4(:) > ZCF(:) ) )
            ! outside the cloud (environment) the use of T^u (unsaturated) instead of T
            ! Bechtold et al. 1993
            !
            ! T^u = T_l = theta_l * (T/theta)
            ZZW2(:) =  ZTHLT(:) * ZZT(:) / ZTHT(:)
            !
            ! es_w with new T^u
            ZZW(:)  = EXP( XALPW - XBETAW/ZZW2(:) - XGAMW*ALOG(ZZW2(:) ) )
            !
            ! S, Undersaturation over water (with new theta^u)
            ZUSW(:) = 1.0 - ZRVT(:)*( ZPRES(:)-ZZW(:) ) / ( (XMV/XMD) * ZZW(:) )
            !
            ZZW(:) = ( XLVTT+(XCPV-XCL)*(ZZW2(:)-XTT) )**2 / ( ZKA(:)*XRV*ZZW2(:)**2 ) &
                   + ( XRV*ZZW2(:) ) / ( ZDV(:)*ZZW(:) )
            !
            ZZW(:) = MAX( 0.0,ZUSW(:) )/(ZRHODREF(:)*ZZW(:))  *      &
                   ( X0EVAR*ZZW3(:)**XEX0EVAR+X1EVAR*ZCJ(:)*ZZW3(:)**XEX1EVAR )
            !
            ZZW(:) = MIN( ZRRS(:),  ZZW(:) *( ZZW4(:) - ZCF(:) ) )
            !
            ZRRS(:) = ZRRS(:) - ZZW(:)
            ZRVS(:) = ZRVS(:) + ZZW(:)
            ZTHS(:) = ZTHS(:) - ZZW(:)*ZLVFACT(:)
          END WHERE
    
        ELSE
          !wrong CSUBG_RR_EVAP case
          WRITE(*,*) 'wrong CSUBG_RR_EVAP case'
          CALL PRINT_MSG(NVERB_FATAL,'GEN','RAIN_ICE_WARM','')  
        END IF
    
    
        IF (LBUDGET_TH) CALL BUDGET (                                               &
                     UNPACK(ZTHS(:),MASK=GMICRO(:,:,:),FIELD=PTHS)*PRHODJ(:,:,:),   &
                                                                  4,'REVA_BU_RTH')
        IF (LBUDGET_RV) CALL BUDGET (                                               &
                     UNPACK(ZRVS(:),MASK=GMICRO(:,:,:),FIELD=PRVS)*PRHODJ(:,:,:),   &
                                                                  6,'REVA_BU_RRV')
    
        IF (LBUDGET_RR) CALL BUDGET (                                                   &
                           UNPACK(ZRRS(:),MASK=GMICRO(:,:,:),FIELD=PRRS)*PRHODJ(:,:,:), &
    
                                                                  8,'REVA_BU_RRR')
        ZW(:,:,:)=PEVAP3D(:,:,:)
        PEVAP3D(:,:,:)=UNPACK(ZZW(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:))
    !
      END SUBROUTINE RAIN_ICE_WARM
    !
    !-------------------------------------------------------------------------------
    !
    !
      SUBROUTINE RAIN_ICE_FAST_RS
    !
    !*      0. DECLARATIONS
    !          ------------
    !
    IMPLICIT NONE
    !
    !-------------------------------------------------------------------------------
    !
    !*       5.1    cloud droplet riming of the aggregates
    !
      ZZW1(:,:) = 0.0
    !
      ALLOCATE(GRIM(IMICRO))
    ! GRIM(:) = (ZRCT(:)>0.0) .AND. (ZRST(:)>0.0) .AND.            &
      GRIM(:) = (ZRCT(:)>XRTMIN(2)) .AND. (ZRST(:)>XRTMIN(5)) .AND.            &
                                    (ZRCS(:)>0.0) .AND. (ZZT(:)<XTT)
      IGRIM = COUNT( GRIM(:) )
    !
      IF( IGRIM>0 ) THEN
    !
    !        5.1.0  allocations
    !
        ALLOCATE(ZVEC1(IGRIM))
        ALLOCATE(ZVEC2(IGRIM))
        ALLOCATE(IVEC1(IGRIM))
        ALLOCATE(IVEC2(IGRIM))
    !
    !        5.1.1  select the ZLBDAS
    !
        ZVEC1(:) = PACK( ZLBDAS(:),MASK=GRIM(:) )
    !
    !        5.1.2  find the next lower indice for the ZLBDAS in the geometrical
    
    !               set of Lbda_s used to tabulate some moments of the incomplete
    
    !               gamma function
    !
        ZVEC2(1:IGRIM) = MAX( 1.00001, MIN( FLOAT(NGAMINC)-0.00001,           &
                              XRIMINTP1 * LOG( ZVEC1(1:IGRIM) ) + XRIMINTP2 ) )
        IVEC2(1:IGRIM) = INT( ZVEC2(1:IGRIM) )
        ZVEC2(1:IGRIM) = ZVEC2(1:IGRIM) - FLOAT( IVEC2(1:IGRIM) )
    !
    !        5.1.3  perform the linear interpolation of the normalized
    !               "2+XDS"-moment of the incomplete gamma function
    !
        ZVEC1(1:IGRIM) =   XGAMINC_RIM1( IVEC2(1:IGRIM)+1 )* ZVEC2(1:IGRIM)      &
                         - XGAMINC_RIM1( IVEC2(1:IGRIM)   )*(ZVEC2(1:IGRIM) - 1.0)
        ZZW(:) = UNPACK( VECTOR=ZVEC1(:),MASK=GRIM,FIELD=0.0 )
    !
    !        5.1.4  riming of the small sized aggregates
    !
        WHERE ( GRIM(:) )
          ZZW1(:,1) = MIN( ZRCS(:),                                &
    
                         XCRIMSS * ZZW(:) * ZRCT(:)                & ! RCRIMSS
                                          *   ZLBDAS(:)**XEXCRIMSS &
                                          * ZRHODREF(:)**(-XCEXVT) )
    
          ZRCS(:) = ZRCS(:) - ZZW1(:,1)
          ZRSS(:) = ZRSS(:) + ZZW1(:,1)
          ZTHS(:) = ZTHS(:) + ZZW1(:,1)*(ZLSFACT(:)-ZLVFACT(:)) ! f(L_f*(RCRIMSS))
        END WHERE
    !
    !        5.1.5  perform the linear interpolation of the normalized
    !               "XBS"-moment of the incomplete gamma function
    !
        ZVEC1(1:IGRIM) =  XGAMINC_RIM2( IVEC2(1:IGRIM)+1 )* ZVEC2(1:IGRIM)      &
    
                        - XGAMINC_RIM2( IVEC2(1:IGRIM)   )*(ZVEC2(1:IGRIM) - 1.0)
    
        ZZW(:) = UNPACK( VECTOR=ZVEC1(:),MASK=GRIM,FIELD=0.0 )
    !
    !        5.1.6  riming-conversion of the large sized aggregates into graupeln
    !
    !
        WHERE ( GRIM(:) .AND. (ZRSS(:)>0.0) )
          ZZW1(:,2) = MIN( ZRCS(:),                     &
    
                       XCRIMSG * ZRCT(:)                & ! RCRIMSG
                               *  ZLBDAS(:)**XEXCRIMSG  &
                               * ZRHODREF(:)**(-XCEXVT) &
                               - ZZW1(:,1)              )
    
          ZZW1(:,3) = MIN( ZRSS(:),                         &
                           XSRIMCG * ZLBDAS(:)**XEXSRIMCG   & ! RSRIMCG
    
                                   * (1.0 - ZZW(:) )/(PTSTEP*ZRHODREF(:)) )
    
          ZRCS(:) = ZRCS(:) - ZZW1(:,2)
          ZRSS(:) = ZRSS(:) - ZZW1(:,3)
          ZRGS(:) = ZRGS(:) + ZZW1(:,2)+ZZW1(:,3)
          ZTHS(:) = ZTHS(:) + ZZW1(:,2)*(ZLSFACT(:)-ZLVFACT(:)) ! f(L_f*(RCRIMSG))
        END WHERE
        DEALLOCATE(IVEC2)
        DEALLOCATE(IVEC1)
        DEALLOCATE(ZVEC2)
        DEALLOCATE(ZVEC1)
      END IF
      IF (LBUDGET_TH) CALL BUDGET (                                               &
                   UNPACK(ZTHS(:),MASK=GMICRO(:,:,:),FIELD=PTHS)*PRHODJ(:,:,:),   &
                                                                 4,'RIM_BU_RTH')
    
      IF (LBUDGET_RC) CALL BUDGET (                                                   &
                         UNPACK(ZRCS(:),MASK=GMICRO(:,:,:),FIELD=PRCS)*PRHODJ(:,:,:), &
                                                                  7,'RIM_BU_RRC')
      IF (LBUDGET_RS) CALL BUDGET (                                                   &
                         UNPACK(ZRSS(:),MASK=GMICRO(:,:,:),FIELD=PRSS)*PRHODJ(:,:,:), &
                                                                 10,'RIM_BU_RRS')
      IF (LBUDGET_RG) CALL BUDGET (                                                   &
                         UNPACK(ZRGS(:),MASK=GMICRO(:,:,:),FIELD=PRGS)*PRHODJ(:,:,:), &
                                                                 11,'RIM_BU_RRG')
    
      DEALLOCATE(GRIM)
    !
    !*       5.2    rain accretion onto the aggregates
    !
      ZZW1(:,2:3) = 0.0
      ALLOCATE(GACC(IMICRO))
       GACC(:) = (ZRRT(:)>XRTMIN(3)) .AND. (ZRST(:)>XRTMIN(5)) .AND.            &
    
                                (ZRRS(:)>0.0) .AND. (ZZT(:)<XTT)
    
      IGACC = COUNT( GACC(:) )
    !
      IF( IGACC>0 ) THEN
    !
    !        5.2.0  allocations
    !
        ALLOCATE(ZVEC1(IGACC))
        ALLOCATE(ZVEC2(IGACC))
        ALLOCATE(ZVEC3(IGACC))
        ALLOCATE(IVEC1(IGACC))
        ALLOCATE(IVEC2(IGACC))
    !
    !        5.2.1  select the (ZLBDAS,ZLBDAR) couplet
    !
        ZVEC1(:) = PACK( ZLBDAS(:),MASK=GACC(:) )
        ZVEC2(:) = PACK( ZLBDAR(:),MASK=GACC(:) )
    !
    !        5.2.2  find the next lower indice for the ZLBDAS and for the ZLBDAR
    !               in the geometrical set of (Lbda_s,Lbda_r) couplet use to
    !               tabulate the RACCSS-kernel
    !
        ZVEC1(1:IGACC) = MAX( 1.00001, MIN( FLOAT(NACCLBDAS)-0.00001,           &
                              XACCINTP1S * LOG( ZVEC1(1:IGACC) ) + XACCINTP2S ) )
        IVEC1(1:IGACC) = INT( ZVEC1(1:IGACC) )
        ZVEC1(1:IGACC) = ZVEC1(1:IGACC) - FLOAT( IVEC1(1:IGACC) )
    !
        ZVEC2(1:IGACC) = MAX( 1.00001, MIN( FLOAT(NACCLBDAR)-0.00001,           &
                              XACCINTP1R * LOG( ZVEC2(1:IGACC) ) + XACCINTP2R ) )
        IVEC2(1:IGACC) = INT( ZVEC2(1:IGACC) )
        ZVEC2(1:IGACC) = ZVEC2(1:IGACC) - FLOAT( IVEC2(1:IGACC) )
    !
    !        5.2.3  perform the bilinear interpolation of the normalized
    !               RACCSS-kernel
    !
        DO JJ = 1,IGACC
          ZVEC3(JJ) =  (  XKER_RACCSS(IVEC1(JJ)+1,IVEC2(JJ)+1)* ZVEC2(JJ)          &
                        - XKER_RACCSS(IVEC1(JJ)+1,IVEC2(JJ)  )*(ZVEC2(JJ) - 1.0) ) &
    
                     - (  XKER_RACCSS(IVEC1(JJ)  ,IVEC2(JJ)+1)* ZVEC2(JJ)          &
                        - XKER_RACCSS(IVEC1(JJ)  ,IVEC2(JJ)  )*(ZVEC2(JJ) - 1.0) ) &
    
        END DO
        ZZW(:) = UNPACK( VECTOR=ZVEC3(:),MASK=GACC,FIELD=0.0 )
    !
    !        5.2.4  raindrop accretion on the small sized aggregates
    !
        WHERE ( GACC(:) )
          ZZW1(:,2) =                                            & !! coef of RRACCS
                  XFRACCSS*( ZLBDAS(:)**XCXS )*( ZRHODREF(:)**(-XCEXVT-1.) ) &
             *( XLBRACCS1/((ZLBDAS(:)**2)               ) +                  &
                XLBRACCS2/( ZLBDAS(:)    * ZLBDAR(:)    ) +                  &
                XLBRACCS3/(               (ZLBDAR(:)**2)) )/ZLBDAR(:)**4
          ZZW1(:,4) = MIN( ZRRS(:),ZZW1(:,2)*ZZW(:) )           ! RRACCSS
          ZRRS(:) = ZRRS(:) - ZZW1(:,4)
          ZRSS(:) = ZRSS(:) + ZZW1(:,4)
          ZTHS(:) = ZTHS(:) + ZZW1(:,4)*(ZLSFACT(:)-ZLVFACT(:)) ! f(L_f*(RRACCSS))
        END WHERE
    !
    !        5.2.4b perform the bilinear interpolation of the normalized
    !               RACCS-kernel
    !
        DO JJ = 1,IGACC
          ZVEC3(JJ) =  (   XKER_RACCS(IVEC2(JJ)+1,IVEC1(JJ)+1)* ZVEC1(JJ)          &
                        -  XKER_RACCS(IVEC2(JJ)+1,IVEC1(JJ)  )*(ZVEC1(JJ) - 1.0) ) &
                                                                       * ZVEC2(JJ) &
                     - (   XKER_RACCS(IVEC2(JJ)  ,IVEC1(JJ)+1)* ZVEC1(JJ)          &
                        -  XKER_RACCS(IVEC2(JJ)  ,IVEC1(JJ)  )*(ZVEC1(JJ) - 1.0) ) &
                                                               * (ZVEC2(JJ) - 1.0)
        END DO
    
        ZZW1(:,2) = ZZW1(:,2)*UNPACK( VECTOR=ZVEC3(:),MASK=GACC(:),FIELD=0.0 )
    
                                                                           !! RRACCS!
    !        5.2.5  perform the bilinear interpolation of the normalized
    !               SACCRG-kernel
    !
        DO JJ = 1,IGACC
          ZVEC3(JJ) =  (  XKER_SACCRG(IVEC2(JJ)+1,IVEC1(JJ)+1)* ZVEC1(JJ)          &
                        - XKER_SACCRG(IVEC2(JJ)+1,IVEC1(JJ)  )*(ZVEC1(JJ) - 1.0) ) &
    
                     - (  XKER_SACCRG(IVEC2(JJ)  ,IVEC1(JJ)+1)* ZVEC1(JJ)          &
                        - XKER_SACCRG(IVEC2(JJ)  ,IVEC1(JJ)  )*(ZVEC1(JJ) - 1.0) ) &
    
        END DO
        ZZW(:) = UNPACK( VECTOR=ZVEC3(:),MASK=GACC,FIELD=0.0 )
    !
    !        5.2.6  raindrop accretion-conversion of the large sized aggregates
    !               into graupeln
    !
        WHERE ( GACC(:) .AND. (ZRSS(:)>0.0) )
          ZZW1(:,2) = MAX( MIN( ZRRS(:),ZZW1(:,2)-ZZW1(:,4) ),0.0 )       ! RRACCSG
        END WHERE
        WHERE ( GACC(:) .AND. (ZRSS(:)>0.0) .AND. ZZW1(:,2)>0.0 )
          ZZW1(:,3) = MIN( ZRSS(:),XFSACCRG*ZZW(:)*                     & ! RSACCRG
                ( ZLBDAS(:)**(XCXS-XBS) )*( ZRHODREF(:)**(-XCEXVT-1.) ) &
               *( XLBSACCR1/((ZLBDAR(:)**2)               ) +           &
                  XLBSACCR2/( ZLBDAR(:)    * ZLBDAS(:)    ) +           &
                  XLBSACCR3/(               (ZLBDAS(:)**2)) )/ZLBDAR(:) )
          ZRRS(:) = ZRRS(:) - ZZW1(:,2)
          ZRSS(:) = ZRSS(:) - ZZW1(:,3)
          ZRGS(:) = ZRGS(:) + ZZW1(:,2)+ZZW1(:,3)
          ZTHS(:) = ZTHS(:) + ZZW1(:,2)*(ZLSFACT(:)-ZLVFACT(:)) !
    
        END WHERE
        DEALLOCATE(IVEC2)
        DEALLOCATE(IVEC1)
        DEALLOCATE(ZVEC3)
        DEALLOCATE(ZVEC2)
        DEALLOCATE(ZVEC1)
      END IF
      DEALLOCATE(GACC)
      IF (LBUDGET_TH) CALL BUDGET (                                               &
                   UNPACK(ZTHS(:),MASK=GMICRO(:,:,:),FIELD=PTHS)*PRHODJ(:,:,:),   &
                                                                 4,'ACC_BU_RTH')
    
      IF (LBUDGET_RR) CALL BUDGET (                                                   &
                         UNPACK(ZRRS(:),MASK=GMICRO(:,:,:),FIELD=PRRS)*PRHODJ(:,:,:), &
                                                                  8,'ACC_BU_RRR')
      IF (LBUDGET_RS) CALL BUDGET (                                                   &
                         UNPACK(ZRSS(:),MASK=GMICRO(:,:,:),FIELD=PRSS)*PRHODJ(:,:,:), &
                                                                 10,'ACC_BU_RRS')
      IF (LBUDGET_RG) CALL BUDGET (                                                   &
                         UNPACK(ZRGS(:),MASK=GMICRO(:,:,:),FIELD=PRGS)*PRHODJ(:,:,:), &
                                                                 11,'ACC_BU_RRG')
    
    !
    !*       5.3    Conversion-Melting of the aggregates
    !
      ZZW(:) = 0.0
      WHERE( (ZRST(:)>XRTMIN(5)) .AND. (ZRSS(:)>0.0) .AND. (ZZT(:)>XTT) )
        ZZW(:) = ZRVT(:)*ZPRES(:)/((XMV/XMD)+ZRVT(:)) ! Vapor pressure
        ZZW(:) =  ZKA(:)*(XTT-ZZT(:)) +                                 &
    
                   ( ZDV(:)*(XLVTT + ( XCPV - XCL ) * ( ZZT(:) - XTT )) &
    
                               *(XESTT-ZZW(:))/(XRV*ZZT(:))             )
    !
    ! compute RSMLT
    !
        ZZW(:)  = MIN( ZRSS(:), XFSCVMG*MAX( 0.0,( -ZZW(:) *             &
                               ( X0DEPS*       ZLBDAS(:)**XEX0DEPS +     &
                                 X1DEPS*ZCJ(:)*ZLBDAS(:)**XEX1DEPS ) -   &
                                         ( ZZW1(:,1)+ZZW1(:,4) ) *       &
    
                                  ( ZRHODREF(:)*XCL*(XTT-ZZT(:))) ) /    &
    
                                                 ( ZRHODREF(:)*XLMTT ) ) )
    !
    ! note that RSCVMG = RSMLT*XFSCVMG but no heat is exchanged (at the rate RSMLT)
    ! because the graupeln produced by this process are still icy!!!
    !
        ZRSS(:) = ZRSS(:) - ZZW(:)
        ZRGS(:) = ZRGS(:) + ZZW(:)
      END WHERE
    
      IF (LBUDGET_RS) CALL BUDGET (                                                   &
                         UNPACK(ZRSS(:),MASK=GMICRO(:,:,:),FIELD=PRSS)*PRHODJ(:,:,:), &
    
                                                                 10,'CMEL_BU_RRS')
    
      IF (LBUDGET_RG) CALL BUDGET (                                                   &
                         UNPACK(ZRGS(:),MASK=GMICRO(:,:,:),FIELD=PRGS)*PRHODJ(:,:,:), &
    
                                                                 11,'CMEL_BU_RRG')
    !
      END SUBROUTINE RAIN_ICE_FAST_RS
    !
    !-------------------------------------------------------------------------------
    !
    !
      SUBROUTINE RAIN_ICE_FAST_RG
    !
    !*      0. DECLARATIONS
    !          ------------
    !
    IMPLICIT NONE
    !
    !-------------------------------------------------------------------------------
    !
    !*       6.1    rain contact freezing
    !
      ZZW1(:,3:4) = 0.0
      WHERE( (ZRIT(:)>XRTMIN(4)) .AND. (ZRRT(:)>XRTMIN(3)) .AND.  &
                                 (ZRIS(:)>0.0) .AND. (ZRRS(:)>0.0) )
        ZZW1(:,3) = MIN( ZRIS(:),XICFRR * ZRIT(:)                & ! RICFRRG
    
                                        * ZLBDAR(:)**XEXICFRR    &
                                        * ZRHODREF(:)**(-XCEXVT) )
    
        ZZW1(:,4) = MIN( ZRRS(:),XRCFRI * ZCIT(:)                & ! RRCFRIG
    
                                        * ZLBDAR(:)**XEXRCFRI    &
                                        * ZRHODREF(:)**(-XCEXVT-1.) )
    
        ZRIS(:) = ZRIS(:) - ZZW1(:,3)
        ZRRS(:) = ZRRS(:) - ZZW1(:,4)
        ZRGS(:) = ZRGS(:) + ZZW1(:,3)+ZZW1(:,4)
        ZTHS(:) = ZTHS(:) + ZZW1(:,4)*(ZLSFACT(:)-ZLVFACT(:)) ! f(L_f*RRCFRIG)
      END WHERE
      IF (LBUDGET_TH) CALL BUDGET (                                                 &
                     UNPACK(ZTHS(:),MASK=GMICRO(:,:,:),FIELD=PTHS)*PRHODJ(:,:,:),   &
                                                                  4,'CFRZ_BU_RTH')
    
      IF (LBUDGET_RR) CALL BUDGET (                                                   &
                         UNPACK(ZRRS(:),MASK=GMICRO(:,:,:),FIELD=PRRS)*PRHODJ(:,:,:), &
    
                                                                  8,'CFRZ_BU_RRR')
    
      IF (LBUDGET_RI) CALL BUDGET (                                                   &
                         UNPACK(ZRIS(:),MASK=GMICRO(:,:,:),FIELD=PRIS)*PRHODJ(:,:,:), &
    
                                                                  9,'CFRZ_BU_RRI')
    
      IF (LBUDGET_RG) CALL BUDGET (                                                   &
                         UNPACK(ZRGS(:),MASK=GMICRO(:,:,:),FIELD=PRGS)*PRHODJ(:,:,:), &
    
                                                                 11,'CFRZ_BU_RRG')
    !
    !*       6.2    compute the Dry growth case
    !
      ZZW1(:,:) = 0.0
      WHERE( (ZRGT(:)>XRTMIN(6)) .AND. ((ZRCT(:)>XRTMIN(2) .AND. ZRCS(:)>0.0)) )
        ZZW(:) = ZLBDAG(:)**(XCXG-XDG-2.0) * ZRHODREF(:)**(-XCEXVT)
        ZZW1(:,1) = MIN( ZRCS(:),XFCDRYG * ZRCT(:) * ZZW(:) )             ! RCDRYG
      END WHERE
      WHERE( (ZRGT(:)>XRTMIN(6)) .AND. ((ZRIT(:)>XRTMIN(4) .AND. ZRIS(:)>0.0)) )
        ZZW(:) = ZLBDAG(:)**(XCXG-XDG-2.0) * ZRHODREF(:)**(-XCEXVT)
        ZZW1(:,2) = MIN( ZRIS(:),XFIDRYG * EXP( XCOLEXIG*(ZZT(:)-XTT) ) &
                                         * ZRIT(:) * ZZW(:) )             ! RIDRYG
      END WHERE
    !
    !*       6.2.1  accretion of aggregates on the graupeln
    !
      ALLOCATE(GDRY(IMICRO))
      GDRY(:) = (ZRST(:)>XRTMIN(5)) .AND. (ZRGT(:)>XRTMIN(6)) .AND. (ZRSS(:)>0.0)
      IGDRY = COUNT( GDRY(:) )
    !
      IF( IGDRY>0 ) THEN
    !
    !*       6.2.2  allocations
    !
        ALLOCATE(ZVEC1(IGDRY))
        ALLOCATE(ZVEC2(IGDRY))
        ALLOCATE(ZVEC3(IGDRY))
        ALLOCATE(IVEC1(IGDRY))
        ALLOCATE(IVEC2(IGDRY))
    !
    !*       6.2.3  select the (ZLBDAG,ZLBDAS) couplet
    !
        ZVEC1(:) = PACK( ZLBDAG(:),MASK=GDRY(:) )
        ZVEC2(:) = PACK( ZLBDAS(:),MASK=GDRY(:) )
    !
    !*       6.2.4  find the next lower indice for the ZLBDAG and for the ZLBDAS
    !               in the geometrical set of (Lbda_g,Lbda_s) couplet use to
    !               tabulate the SDRYG-kernel
    !
        ZVEC1(1:IGDRY) = MAX( 1.00001, MIN( FLOAT(NDRYLBDAG)-0.00001,           &
                              XDRYINTP1G * LOG( ZVEC1(1:IGDRY) ) + XDRYINTP2G ) )
        IVEC1(1:IGDRY) = INT( ZVEC1(1:IGDRY) )
        ZVEC1(1:IGDRY) = ZVEC1(1:IGDRY) - FLOAT( IVEC1(1:IGDRY) )
    !
        ZVEC2(1:IGDRY) = MAX( 1.00001, MIN( FLOAT(NDRYLBDAS)-0.00001,           &
                              XDRYINTP1S * LOG( ZVEC2(1:IGDRY) ) + XDRYINTP2S ) )
        IVEC2(1:IGDRY) = INT( ZVEC2(1:IGDRY) )
        ZVEC2(1:IGDRY) = ZVEC2(1:IGDRY) - FLOAT( IVEC2(1:IGDRY) )
    !
    !*       6.2.5  perform the bilinear interpolation of the normalized
    !               SDRYG-kernel
    !
        DO JJ = 1,IGDRY
          ZVEC3(JJ) =  (  XKER_SDRYG(IVEC1(JJ)+1,IVEC2(JJ)+1)* ZVEC2(JJ)          &
                        - XKER_SDRYG(IVEC1(JJ)+1,IVEC2(JJ)  )*(ZVEC2(JJ) - 1.0) ) &
    
                     - (  XKER_SDRYG(IVEC1(JJ)  ,IVEC2(JJ)+1)* ZVEC2(JJ)          &
                        - XKER_SDRYG(IVEC1(JJ)  ,IVEC2(JJ)  )*(ZVEC2(JJ) - 1.0) ) &
    
        END DO
        ZZW(:) = UNPACK( VECTOR=ZVEC3(:),MASK=GDRY,FIELD=0.0 )
    !
        WHERE( GDRY(:) )
          ZZW1(:,3) = MIN( ZRSS(:),XFSDRYG*ZZW(:)                         & ! RSDRYG
                                          * EXP( XCOLEXSG*(ZZT(:)-XTT) )  &
                        *( ZLBDAS(:)**(XCXS-XBS) )*( ZLBDAG(:)**XCXG )    &
    
                        *( ZRHODREF(:)**(-XCEXVT-1.) )                    &
    
                             *( XLBSDRYG1/( ZLBDAG(:)**2              ) + &
                                XLBSDRYG2/( ZLBDAG(:)   * ZLBDAS(:)   ) + &
                                XLBSDRYG3/(               ZLBDAS(:)**2) ) )
        END WHERE
        DEALLOCATE(IVEC2)
        DEALLOCATE(IVEC1)
        DEALLOCATE(ZVEC3)
        DEALLOCATE(ZVEC2)
        DEALLOCATE(ZVEC1)
      END IF
    !
    !*       6.2.6  accretion of raindrops on the graupeln
    !
      GDRY(:) = (ZRRT(:)>XRTMIN(3)) .AND. (ZRGT(:)>XRTMIN(6)) .AND. (ZRRS(:)>0.0)
      IGDRY = COUNT( GDRY(:) )
    !
      IF( IGDRY>0 ) THEN
    !
    !*       6.2.7  allocations
    !
        ALLOCATE(ZVEC1(IGDRY))
        ALLOCATE(ZVEC2(IGDRY))
        ALLOCATE(ZVEC3(IGDRY))
        ALLOCATE(IVEC1(IGDRY))
        ALLOCATE(IVEC2(IGDRY))
    !
    !*       6.2.8  select the (ZLBDAG,ZLBDAR) couplet
    !
        ZVEC1(:) = PACK( ZLBDAG(:),MASK=GDRY(:) )
        ZVEC2(:) = PACK( ZLBDAR(:),MASK=GDRY(:) )
    !
    !*       6.2.9  find the next lower indice for the ZLBDAG and for the ZLBDAR
    !               in the geometrical set of (Lbda_g,Lbda_r) couplet use to
    !               tabulate the RDRYG-kernel
    !
        ZVEC1(1:IGDRY) = MAX( 1.00001, MIN( FLOAT(NDRYLBDAG)-0.00001,           &
                              XDRYINTP1G * LOG( ZVEC1(1:IGDRY) ) + XDRYINTP2G ) )
        IVEC1(1:IGDRY) = INT( ZVEC1(1:IGDRY) )
        ZVEC1(1:IGDRY) = ZVEC1(1:IGDRY) - FLOAT( IVEC1(1:IGDRY) )
    !
        ZVEC2(1:IGDRY) = MAX( 1.00001, MIN( FLOAT(NDRYLBDAR)-0.00001,           &
                              XDRYINTP1R * LOG( ZVEC2(1:IGDRY) ) + XDRYINTP2R ) )
        IVEC2(1:IGDRY) = INT( ZVEC2(1:IGDRY) )
        ZVEC2(1:IGDRY) = ZVEC2(1:IGDRY) - FLOAT( IVEC2(1:IGDRY) )
    !
    !*       6.2.10 perform the bilinear interpolation of the normalized
    !               RDRYG-kernel
    !
        DO JJ = 1,IGDRY
          ZVEC3(JJ) =  (  XKER_RDRYG(IVEC1(JJ)+1,IVEC2(JJ)+1)* ZVEC2(JJ)          &
                        - XKER_RDRYG(IVEC1(JJ)+1,IVEC2(JJ)  )*(ZVEC2(JJ) - 1.0) ) &
    
                     - (  XKER_RDRYG(IVEC1(JJ)  ,IVEC2(JJ)+1)* ZVEC2(JJ)          &
                        - XKER_RDRYG(IVEC1(JJ)  ,IVEC2(JJ)  )*(ZVEC2(JJ) - 1.0) ) &
    
        END DO
        ZZW(:) = UNPACK( VECTOR=ZVEC3(:),MASK=GDRY,FIELD=0.0 )
    !
        WHERE( GDRY(:) )
          ZZW1(:,4) = MIN( ZRRS(:),XFRDRYG*ZZW(:)                    & ! RRDRYG
                            *( ZLBDAR(:)**(-4) )*( ZLBDAG(:)**XCXG ) &
    
                                   *( ZRHODREF(:)**(-XCEXVT-1.) )   &
    
                        *( XLBRDRYG1/( ZLBDAG(:)**2              ) + &
                           XLBRDRYG2/( ZLBDAG(:)   * ZLBDAR(:)   ) + &
                           XLBRDRYG3/(               ZLBDAR(:)**2) ) )
        END WHERE
        DEALLOCATE(IVEC2)
        DEALLOCATE(IVEC1)
        DEALLOCATE(ZVEC3)
        DEALLOCATE(ZVEC2)
        DEALLOCATE(ZVEC1)
      END IF
    !
      ZRDRYG(:) = ZZW1(:,1) + ZZW1(:,2) + ZZW1(:,3) + ZZW1(:,4)
      DEALLOCATE(GDRY)
    !
    !*       6.3    compute the Wet growth case
    !
      ZZW(:) = 0.0
      ZRWETG(:) = 0.0
      WHERE( ZRGT(:)>XRTMIN(6) )
        ZZW1(:,5) = MIN( ZRIS(:),                                    &
                    ZZW1(:,2) / (XCOLIG*EXP(XCOLEXIG*(ZZT(:)-XTT)) ) ) ! RIWETG
        ZZW1(:,6) = MIN( ZRSS(:),                                    &
                    ZZW1(:,3) / (XCOLSG*EXP(XCOLEXSG*(ZZT(:)-XTT)) ) ) ! RSWETG
    !
        ZZW(:) = ZRVT(:)*ZPRES(:)/((XMV/XMD)+ZRVT(:)) ! Vapor pressure
        ZZW(:) =   ZKA(:)*(XTT-ZZT(:)) +                              &