Skip to content
Snippets Groups Projects
rain_ice.f90 131 KiB
Newer Older
  • Learn to ignore specific revisions
  • !MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
    
    !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
    
    !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
    
    !MNH_LIC for details. version 1.
    
    !-----------------------------------------------------------------
    !--------------- special set of characters for RCS information
    !-----------------------------------------------------------------
    ! $Source$ $Revision$
    ! masdev4_7 BUG1 2007/06/15 17:47:18
    !-----------------------------------------------------------------
    !      ####################
           MODULE MODI_RAIN_ICE
    !      ####################
    !
    INTERFACE
    
          SUBROUTINE RAIN_ICE ( OSEDIC,HSEDIM, HSUBG_AUCV, OWARM, KKA, KKU, KKL,      &
                                KSPLITR, PTSTEP, KMI, KRR,                            &
                                PDZZ, PRHODJ, PRHODREF, PEXNREF, PPABST, PCIT, PCLDFR,&
                                PTHT, PRVT, PRCT, PRRT, PRIT, PRST,                   &
                                PRGT, PTHS, PRVS, PRCS, PRRS, PRIS, PRSS, PRGS,       &
                                PINPRC,PINPRC3D, PINPRR, PINPRR3D, PEVAP3D,           &
                                PINPRS,PINPRS3D, PINPRG,PINPRG3D, PSIGS,              &
                                PSPEEDC, PSPEEDR, PSPEEDS, PSPEEDG, PSPEEDH,          &
                                PSEA, PTOWN,                                          &
                                PRHT,  PRHS, PINPRH,PINPRH3D,OCONVHG             )
    
    !
    !
    LOGICAL,                  INTENT(IN)    :: OSEDIC ! Switch for droplet sedim.
    CHARACTER(LEN=4),         INTENT(IN)    :: HSEDIM ! Sedimentation scheme
    CHARACTER(LEN=4),         INTENT(IN)    :: HSUBG_AUCV ! Switch for Subgrid autoconversion
                                            ! Kind of Subgrid autoconversion method
    LOGICAL,                  INTENT(IN)    :: OWARM   ! .TRUE. allows raindrops to
                                                       !   form by warm processes
                                                       !      (Kessler scheme)
    !
    INTEGER,                  INTENT(IN)    :: KKA   !near ground array index  
    INTEGER,                  INTENT(IN)    :: KKU   !uppest atmosphere array index
    INTEGER,                  INTENT(IN)    :: KKL   !vert. levels type 1=MNH -1=ARO
    INTEGER,                  INTENT(IN)    :: KSPLITR ! Number of small time step 
                                          ! integration for  rain sedimendation
    REAL,                     INTENT(IN)    :: PTSTEP  ! Double Time step
                                                       ! (single if cold start)
    INTEGER,                  INTENT(IN)    :: KMI     ! Model index 
    INTEGER,                  INTENT(IN)    :: KRR     ! Number of moist variable
    !
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PDZZ     ! Layer thikness (m)
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRHODJ  ! Dry density * Jacobian
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRHODREF! Reference density
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PEXNREF ! Reference Exner function
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PPABST  ! absolute pressure at t
    !
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PCIT    ! Pristine ice n.c. at t
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PCLDFR  ! Cloud fraction
    !
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PTHT    ! Theta at time t
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRVT    ! Water vapor m.r. at t 
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRCT    ! Cloud water m.r. at t 
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRRT    ! Rain water m.r. at t 
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRIT    ! Pristine ice m.r. at t
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRST    ! Snow/aggregate m.r. at t
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRGT    ! Graupel/hail m.r. at t
    !
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PSIGS   ! Sigma_s at t
    !
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PTHS    ! Theta source
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRVS    ! Water vapor m.r. source
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRCS    ! Cloud water m.r. source
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRRS    ! Rain water m.r. source
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRIS    ! Pristine ice m.r. source
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRSS    ! Snow/aggregate m.r. source
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRGS    ! Graupel m.r. source
    
    
    REAL, DIMENSION(:,:,:),   INTENT(OUT)   :: PSPEEDC ! Cloud sedimentation speed
    REAL, DIMENSION(:,:,:),   INTENT(OUT)   :: PSPEEDR ! Rain sedimentation speed
    REAL, DIMENSION(:,:,:),   INTENT(OUT)   :: PSPEEDS ! Snow sedimentation speed
    REAL, DIMENSION(:,:,:),   INTENT(OUT)   :: PSPEEDG ! Graupel sedimentation speed
    REAL, DIMENSION(:,:,:),   INTENT(OUT)   :: PSPEEDH ! Hail sedimentation speed
    
    !
    REAL, DIMENSION(:,:), INTENT(INOUT)     :: PINPRC! Cloud instant precip
    
    REAL, DIMENSION(:,:,:), INTENT(INOUT)   :: PINPRC3D! Cloud inst precip 3D
    
    REAL, DIMENSION(:,:), INTENT(INOUT)     :: PINPRR! Rain instant precip
    REAL, DIMENSION(:,:,:), INTENT(INOUT)   :: PINPRR3D! Rain inst precip 3D
    REAL, DIMENSION(:,:,:), INTENT(INOUT)   :: PEVAP3D! Rain evap profile
    REAL, DIMENSION(:,:), INTENT(INOUT)     :: PINPRS! Snow instant precip
    
    REAL, DIMENSION(:,:,:), INTENT(INOUT)   :: PINPRS3D! Snow inst precip 3D
    
    REAL, DIMENSION(:,:), INTENT(INOUT)     :: PINPRG! Graupel instant precip
    
    REAL, DIMENSION(:,:,:), INTENT(INOUT)   :: PINPRG3D! Graupel inst precip 3D
    
    REAL, DIMENSION(:,:),OPTIONAL,INTENT(IN)         :: PSEA
    REAL, DIMENSION(:,:),OPTIONAL,INTENT(IN)         :: PTOWN
    REAL, DIMENSION(:,:,:), OPTIONAL,  INTENT(IN)    :: PRHT    ! Hail m.r. at t
    REAL, DIMENSION(:,:,:), OPTIONAL,  INTENT(INOUT) :: PRHS    ! Hail m.r. source
    
    REAL, DIMENSION(:,:), OPTIONAL, INTENT(INOUT)    :: PINPRH! Hail instant precip
    
    REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(INOUT)     :: PINPRH3D! Hail instant precip 3D
    
    LOGICAL, OPTIONAL,                 INTENT(IN)    :: OCONVHG! Switch for conversion from
                                                      ! hail to graupel
    
    
    !
    END SUBROUTINE RAIN_ICE
    END INTERFACE
    END MODULE MODI_RAIN_ICE
    !     ######spl
    
          SUBROUTINE RAIN_ICE ( OSEDIC,HSEDIM, HSUBG_AUCV, OWARM, KKA, KKU, KKL,          &
                                KSPLITR, PTSTEP, KMI, KRR,                                &
                                PDZZ, PRHODJ, PRHODREF, PEXNREF, PPABST, PCIT, PCLDFR,    &
                                PTHT, PRVT, PRCT, PRRT, PRIT, PRST,                       &
                                PRGT, PTHS, PRVS, PRCS, PRRS, PRIS, PRSS, PRGS,           &
    
                                PINPRC,PINPRC3D, PINPRR, PINPRR3D, PEVAP3D,               &
                                PINPRS,PINPRS3D, PINPRG,PINPRG3D, PSIGS,                  &
                                PSPEEDC, PSPEEDR, PSPEEDS, PSPEEDG, PSPEEDH,              &
                                PSEA, PTOWN,                                              &
                                PRHT,  PRHS, PINPRH,PINPRH3D,OCONVHG                 )
    !     #####################################################################
    
    !
    !!****  * -  compute the explicit microphysical sources
    !!
    !!    PURPOSE
    !!    -------
    !!      The purpose of this routine is to compute the slow microphysical sources
    !!    which can be computed explicitly
    !!
    !!
    !!**  METHOD
    !!    ------
    !!      The autoconversion computation follows Kessler (1969).
    !!      The sedimentation rate is computed with a time spliting technique and 
    !!    an upstream scheme, written as a difference of non-advective fluxes. This
    !!    source term is added to the future instant ( split-implicit process ).
    !!      The others microphysical processes are evaluated at the central instant 
    !!    (split-explicit process ): autoconversion, accretion and rain evaporation.
    !!      These last 3 terms are bounded in order not to create negative values 
    !!    for the water species at the future instant.
    !!
    !!    EXTERNAL
    !!    --------
    !!      None
    !!
    !!
    !!    IMPLICIT ARGUMENTS
    !!    ------------------
    !!      Module MODD_PARAMETERS
    !!          JPHEXT       : Horizontal external points number
    !!          JPVEXT       : Vertical external points number
    !!      Module MODD_CONF :
    !!          CCONF configuration of the model for the first time step
    !!      Module MODD_CST
    !!          XP00               ! Reference pressure
    !!          XRD,XRV            ! Gaz  constant for dry air, vapor
    !!          XMD,XMV            ! Molecular weight for dry air, vapor
    !!          XCPD               ! Cpd (dry air)
    !!          XCL                ! Cl (liquid)
    !!          XCI                ! Ci (solid)
    !!          XTT                ! Triple point temperature
    !!          XLVTT              ! Vaporization heat constant
    !!          XALPW,XBETAW,XGAMW ! Constants for saturation vapor pressure
    !!                               function over liquid water
    !!          XALPI,XBETAI,XGAMI ! Constants for saturation vapor pressure
    !!                               function over solid ice
    !!      Module MODD_BUDGET:
    !!         NBUMOD       : model in which budget is calculated
    !!         CBUTYPE      : type of desired budget
    !!                          'CART' for cartesian box configuration
    !!                          'MASK' for budget zone defined by a mask 
    !!                          'NONE'  ' for no budget
    !!         NBUPROCCTR   : process counter used for each budget variable
    !!         LBU_RTH      : logical for budget of RTH (potential temperature)
    !!                        .TRUE. = budget of RTH
    !!                        .FALSE. = no budget of RTH
    !!         LBU_RRV      : logical for budget of RRV (water vapor)
    !!                        .TRUE. = budget of RRV
    !!                        .FALSE. = no budget of RRV
    !!         LBU_RRC      : logical for budget of RRC (cloud water)
    !!                        .TRUE. = budget of RRC
    !!                        .FALSE. = no budget of RRC
    !!         LBU_RRI      : logical for budget of RRI (cloud ice)
    !!                        .TRUE. = budget of RRI
    !!                        .FALSE. = no budget of RRI
    !!         LBU_RRR      : logical for budget of RRR (rain water)
    !!                        .TRUE. = budget of RRR
    !!                        .FALSE. = no budget of RRR
    !!         LBU_RRS      : logical for budget of RRS (aggregates)
    !!                        .TRUE. = budget of RRS
    !!                        .FALSE. = no budget of RRS
    !!         LBU_RRG      : logical for budget of RRG (graupeln)
    !!                        .TRUE. = budget of RRG
    !!                        .FALSE. = no budget of RRG
    !!
    !!    REFERENCE
    !!    ---------
    !!
    !!      Book1 and Book2 of documentation ( routine RAIN_ICE )
    !!
    !!    AUTHOR
    !!    ------
    !!      J.-P. Pinty      * Laboratoire d'Aerologie*
    !!
    !!    MODIFICATIONS
    !!    -------------
    !!      Original    02/11/95
    !!      (J.Viviand) 04/02/97  debug accumulated prcipitation & convert
    !!                            precipitation rate in m/s
    !!      (J.-P. Pinty) 17/02/97  add budget calls
    !!      (J.-P. Pinty) 17/11/97  set ice sedim. for cirrus ice, reset RCHONI
    !!                              and RRHONG, reverse order for DEALLOCATE
    !!      (J.-P. Pinty) 11/02/98  correction of the air dynamical viscosity and
    !!                              add advance of the budget calls
    !!      (J.-P. Pinty) 18/05/98  correction of the air density in the RIAUTS
    !!                              process
    !!      (J.-P. Pinty) 18/11/98  split the main routine
    !!      (V. Masson)   18/11/98  bug in IVEC1 and IVEC2 upper limits
    !!      (J. Escobar & J.-P. Pinty)
    !!                    11/12/98  contains and rewrite count+pack
    !!      (J. Stein & J.-P. Pinty)
    !!                    14/10/99  correction for very small RIT
    !!      (J. Escobar & J.-P. Pinty)
    !!                    24/07/00  correction for very samll m.r. in
    !!                              the sedimentation subroutine
    !!      (M. Tomasini) 11/05/01  Autoconversion of rc into rr modification to take 
    !!                              into account the subgrid variance 
    !!                              (cf Redelsperger & Sommeria JAS 86)
    !!      (G. Molinie)  21/05/99  bug in RRCFRIG process, RHODREF**(-1) missing
    !!                              in RSRIMCG
    !!      (G. Molinie & J.-P. Pinty)
    !!                    21/06/99  bug in RACCS process
    !!      (P. Jabouille) 27/05/04 safety test for case where esw/i(T)> pabs (~Z>40km)
    !!      (J-.P. Chaboureau) 12/02/05  temperature depending ice-to-snow autocon-
    !                              version threshold (Chaboureau and Pinty GRL 2006)
    !!      (J.-P. Pinty) 01/01/O1  add the hail category and correction of the
    !!                              wet growth rate of the graupeln
    !!      (S.Remy & C.Lac) 06/06 Add the cloud sedimentation
    !!      (S.Remy & C.Lac) 06/06 Sedimentation becoming the last process
    !!      to settle the precipitating species created during the current time step
    !!      (S.Remy & C.Lac) 06/06 Modification of the algorithm of sedimentation
    !!      to settle n times the precipitating species created during Dt/n instead
    !!      of Dt
    !!      (C.Lac) 11/06 Optimization of the sedimentation loop for NEC
    !!      (J.Escobar) 18/01/2008 Parallel Bug in Budget when IMICRO >= 1
    !!                  --> Path inhibit this test by IMICRO >= 0 allway true
    !!      (Y.Seity) 03/2008 Add Statistic sedimentation
    !!      (Y.Seity) 10/2009 Added condition for the raindrop accretion of the aggregates
    !!         into graupeln process (5.2.6) to avoid negative graupel mixing ratio
    !!      (V.Masson, C.Lac) 09/2010 Correction in split sedimentation for
    !!                                reproducibility
    !!      (S. Riette) Oct 2010 Better vectorisation of RAIN_ICE_SEDIMENTATION_STAT
    !!      (Y. Seity), 02-2012  add possibility to run with reversed vertical levels
    !!      Juan 24/09/2012: for BUG Pgi rewrite PACK function on mode_pack_pgi
    
    !!      (C. Lac) FIT temporal scheme : instant M removed
    
    !!      (JP Pinty), 01-2014 : ICE4 : partial reconversion of hail to graupel
    
    !!              July, 2015 (O.Nuissier/F.Duffourg) Add microphysics diagnostic for
    !!                                      aircraft, ballon and profiler
    
    !-------------------------------------------------------------------------------
    !
    !*       0.    DECLARATIONS
    !              ------------
    !
    USE MODD_PARAMETERS
    USE MODD_CST
    USE MODD_CONF
    USE MODD_RAIN_ICE_DESCR
    USE MODD_RAIN_ICE_PARAM
    USE MODD_PARAM_ICE
    USE MODD_BUDGET
    USE MODD_LES
    USE MODI_BUDGET
    USE MODI_GAMMA
    
    USE MODE_FMWRIT
    
    !
    #ifdef MNH_PGI
    USE MODE_PACK_PGI
    #endif
    !
    IMPLICIT NONE
    !
    !*       0.1   Declarations of dummy arguments :
    !
    !
    !
    LOGICAL,                  INTENT(IN)    :: OSEDIC ! Switch for droplet sedim.
    CHARACTER(LEN=4),         INTENT(IN)    :: HSEDIM ! Sedimentation scheme
    CHARACTER(LEN=4),         INTENT(IN)    :: HSUBG_AUCV
                                            ! Kind of Subgrid autoconversion method
    LOGICAL,                  INTENT(IN)    :: OWARM   ! .TRUE. allows raindrops to
                                                       !   form by warm processes
                                                       !      (Kessler scheme)
    !
    INTEGER,                  INTENT(IN)    :: KKA   !near ground array index  
    INTEGER,                  INTENT(IN)    :: KKU   !uppest atmosphere array index
    INTEGER,                  INTENT(IN)    :: KKL   !vert. levels type 1=MNH -1=ARO
    INTEGER,                  INTENT(IN)    :: KSPLITR ! Number of small time step 
                                          ! integration for  rain sedimendation
    REAL,                     INTENT(IN)    :: PTSTEP  ! Double Time step
                                                       ! (single if cold start)
    INTEGER,                  INTENT(IN)    :: KMI     ! Model index 
    INTEGER,                  INTENT(IN)    :: KRR     ! Number of moist variable
    !
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PDZZ    ! Layer thikness (m)
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRHODJ  ! Dry density * Jacobian
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRHODREF! Reference density
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PEXNREF ! Reference Exner function
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PPABST  ! absolute pressure at t
    !
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PCIT    ! Pristine ice n.c. at t
    REAL, DIMENSION(:,:,:),     INTENT(IN)  :: PCLDFR! Convective Mass Flux Cloud fraction
    !
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PTHT    ! Theta at time t
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRVT    ! Water vapor m.r. at t 
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRCT    ! Cloud water m.r. at t 
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRRT    ! Rain water m.r. at t 
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRIT    ! Pristine ice m.r. at t
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRST    ! Snow/aggregate m.r. at t
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRGT    ! Graupel/hail m.r. at t
    REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PSIGS   ! Sigma_s at t
    !
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PTHS    ! Theta source
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRVS    ! Water vapor m.r. source
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRCS    ! Cloud water m.r. source
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRRS    ! Rain water m.r. source
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRIS    ! Pristine ice m.r. source
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRSS    ! Snow/aggregate m.r. source
    REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRGS    ! Graupel m.r. source
    !
    
    REAL, DIMENSION(:,:,:),   INTENT(OUT)   :: PSPEEDC ! Cloud sedimentation speed
    REAL, DIMENSION(:,:,:),   INTENT(OUT)   :: PSPEEDR ! Rain sedimentation speed
    REAL, DIMENSION(:,:,:),   INTENT(OUT)   :: PSPEEDS ! Snow sedimentation speed
    REAL, DIMENSION(:,:,:),   INTENT(OUT)   :: PSPEEDG ! Graupel sedimentation speed
    REAL, DIMENSION(:,:,:),   INTENT(OUT)   :: PSPEEDH ! Hail sedimentation speed
    !
    
    REAL, DIMENSION(:,:), INTENT(INOUT)     :: PINPRC! Cloud instant precip
    
    REAL, DIMENSION(:,:,:), INTENT(INOUT)   :: PINPRC3D! Cloud inst precip 3D
    
    REAL, DIMENSION(:,:), INTENT(INOUT)     :: PINPRR! Rain instant precip
    REAL, DIMENSION(:,:,:), INTENT(INOUT)   :: PINPRR3D! Rain inst precip 3D
    REAL, DIMENSION(:,:,:), INTENT(INOUT)   :: PEVAP3D! Rain evap profile
    REAL, DIMENSION(:,:), INTENT(INOUT)     :: PINPRS! Snow instant precip
    
    REAL, DIMENSION(:,:,:), INTENT(INOUT)   :: PINPRS3D! Snow inst precip 3D
    
    REAL, DIMENSION(:,:), INTENT(INOUT)     :: PINPRG! Graupel instant precip
    
    REAL, DIMENSION(:,:,:), INTENT(INOUT)   :: PINPRG3D! Graupel inst precip 3D
    
    REAL, DIMENSION(:,:),OPTIONAL,INTENT(IN)         :: PSEA
    REAL, DIMENSION(:,:),OPTIONAL,INTENT(IN)         :: PTOWN
    REAL, DIMENSION(:,:,:), OPTIONAL,  INTENT(IN)    :: PRHT    ! Hail m.r. at t
    REAL, DIMENSION(:,:,:), OPTIONAL,  INTENT(INOUT) :: PRHS    ! Hail m.r. source
    
    REAL, DIMENSION(:,:), OPTIONAL, INTENT(INOUT)    :: PINPRH! Hail instant precip
    
    REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(INOUT)     :: PINPRH3D! Hail instant precip 3D
    
    LOGICAL, OPTIONAL,                 INTENT(IN)    :: OCONVHG! Switch for conversion from
                                                      ! hail to graupel
    
    !
    !*       0.2   Declarations of local variables :
    !
    INTEGER :: JK            ! Vertical loop index for the rain sedimentation 
    INTEGER :: JN            ! Temporal loop index for the rain sedimentation
    INTEGER :: JJ            ! Loop index for the interpolation
    INTEGER :: JI            ! Loop index for the interpolation
    INTEGER :: IIB           !  Define the domain where is 
    INTEGER :: IIE           !  the microphysical sources have to be computed
    INTEGER :: IJB           !
    INTEGER :: IJE           !
    INTEGER :: IKB,IKTB,IKT  !
    INTEGER :: IKE,IKTE      !
    !
    REAL    :: ZTSPLITR      ! Small time step for rain sedimentation
    !
    !
    INTEGER :: ISEDIMR,ISEDIMC, ISEDIMI, ISEDIMS, ISEDIMG, ISEDIMH, &
      INEGT, IMICRO ! Case number of sedimentation, T>0 (for HEN)
    
    INTEGER :: IGRIM, IGACC, IGDRY ! Case number of riming, accretion and dry growth
                                   ! locations
    INTEGER :: IGWET, IHAIL   ! wet growth locations and case number
    LOGICAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3)) &
        :: GSEDIMR,GSEDIMC, GSEDIMI, GSEDIMS, GSEDIMG, GSEDIMH ! Test where to compute the SED processes
    LOGICAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3)) &
    
      :: GNEGT  ! Test where to compute the HEN process
    
    LOGICAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3)) &
    
      :: GMICRO ! Test where to compute all processes
    
    LOGICAL, DIMENSION(:), ALLOCATABLE :: GRIM ! Test where to compute riming
    LOGICAL, DIMENSION(:), ALLOCATABLE :: GACC ! Test where to compute accretion
    LOGICAL, DIMENSION(:), ALLOCATABLE :: GDRY ! Test where to compute dry growth
    LOGICAL, DIMENSION(:), ALLOCATABLE :: GWET  ! Test where to compute wet growth
    LOGICAL, DIMENSION(:), ALLOCATABLE :: GHAIL ! Test where to compute hail growth
    INTEGER, DIMENSION(:), ALLOCATABLE :: IVEC1,IVEC2       ! Vectors of indices for
    
    REAL,    DIMENSION(:), ALLOCATABLE :: ZVEC1,ZVEC2,ZVEC3 ! Work vectors for 
    
    REAL,    DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3))   &
                                      :: ZW ! work array
    REAL,    DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3))   &
          :: ZPRCS,ZPRRS,ZPRSS,ZPRGS,ZPRHS   ! Mixing ratios created during the time step
    REAL,    DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),0:SIZE(PEXNREF,3)+1)   &
                                      :: ZWSED        ! sedimentation fluxes
    REAL,    DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),0:SIZE(PEXNREF,3)+1)   &
                                      :: ZWSEDW1       ! sedimentation speed
    REAL,    DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),0:SIZE(PEXNREF,3)+1)   &
                                      :: ZWSEDW2       ! sedimentation speed
    REAL,    DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2))                   &
                                      :: ZCONC_TMP    ! Weighted concentration
    REAL,    DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3))   &
                                      :: ZT ! Temperature
    REAL,    DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2),SIZE(PRHODREF,3)) ::  &
                                         ZRAY,   & ! Cloud Mean radius
                                         ZLBC,   & ! XLBC weighted by sea fraction
                                         ZFSEDC
    REAL, DIMENSION(:), ALLOCATABLE :: ZRVT    ! Water vapor m.r. at t 
    REAL, DIMENSION(:), ALLOCATABLE :: ZRCT    ! Cloud water m.r. at t 
    REAL, DIMENSION(:), ALLOCATABLE :: ZRRT    ! Rain water m.r. at t 
    REAL, DIMENSION(:), ALLOCATABLE :: ZRIT    ! Pristine ice m.r. at t
    REAL, DIMENSION(:), ALLOCATABLE :: ZRST    ! Snow/aggregate m.r. at t
    REAL, DIMENSION(:), ALLOCATABLE :: ZRGT    ! Graupel m.r. at t
    REAL, DIMENSION(:), ALLOCATABLE :: ZRHT    ! Hail m.r. at t
    REAL, DIMENSION(:), ALLOCATABLE :: ZCIT    ! Pristine ice conc. at t
    !
    REAL, DIMENSION(:), ALLOCATABLE :: ZRVS    ! Water vapor m.r. source
    REAL, DIMENSION(:), ALLOCATABLE :: ZRCS    ! Cloud water m.r. source
    REAL, DIMENSION(:), ALLOCATABLE :: ZRRS    ! Rain water m.r. source
    REAL, DIMENSION(:), ALLOCATABLE :: ZRIS    ! Pristine ice m.r. source
    REAL, DIMENSION(:), ALLOCATABLE :: ZRSS    ! Snow/aggregate m.r. source
    REAL, DIMENSION(:), ALLOCATABLE :: ZRGS    ! Graupel m.r. source
    REAL, DIMENSION(:), ALLOCATABLE :: ZRHS    ! Hail m.r. source
    REAL, DIMENSION(:), ALLOCATABLE :: ZTHS    ! Theta source
    REAL, DIMENSION(:), ALLOCATABLE :: ZCRIAUTI ! Snow-to-ice autoconversion thres.
    !
    REAL, DIMENSION(:), ALLOCATABLE &
                   :: ZRHODREF, & ! RHO Dry REFerence
                      ZRHODREFC,& ! RHO Dry REFerence
                      ZRHODREFR,& ! RHO Dry REFerence
                      ZRHODREFI,& ! RHO Dry REFerence
                      ZRHODREFS,& ! RHO Dry REFerence
                      ZRHODREFG,& ! RHO Dry REFerence
                      ZRHODREFH,& ! RHO Dry REFerence
                      ZRHODJ,   & ! RHO times Jacobian
                      ZZT,      & ! Temperature
                      ZPRES,    & ! Pressure
                      ZEXNREF,  & ! EXNer Pressure REFerence
    
                      ZZW,      & ! Work array
                      ZLSFACT,  & ! L_s/(Pi_ref*C_ph)
                      ZLVFACT,  & ! L_v/(Pi_ref*C_ph)
    
                      ZUSW,     & ! Undersaturation over water
                      ZSSI,     & ! Supersaturation over ice
                      ZLBDAR,   & ! Slope parameter of the raindrop  distribution
                      ZLBDAS,   & ! Slope parameter of the aggregate distribution
                      ZLBDAG,   & ! Slope parameter of the graupel   distribution
                      ZLBDAH,   & ! Slope parameter of the hail      distribution
                      ZRDRYG,   & ! Dry growth rate of the graupeln
                      ZRWETG,   & ! Wet growth rate of the graupeln
                      ZAI,      & ! Thermodynamical function
                      ZCJ,      & ! Function to compute the ventilation coefficient
                      ZKA,      & ! Thermal conductivity of the air
                      ZDV,      & ! Diffusivity of water vapor in the air
                      ZSIGMA_RC,& ! Standard deviation of rc at time t  
                      ZCF,      & ! Cloud fraction
                      ZCC,      & ! terminal velocity
                      ZFSEDC1D, & ! For cloud sedimentation
                      ZWLBDC,   & ! Slope parameter of the droplet  distribution
                      ZCONC,    & ! Concentration des aerosols
                      ZRAY1D,   & ! Mean radius
                      ZWLBDA      ! Libre parcours moyen
    REAL, DIMENSION(:,:), ALLOCATABLE :: ZZW1 ! Work arrays
    
    REAL            :: ZINVTSTEP
    REAL, DIMENSION(SIZE(XRTMIN))     :: ZRTMIN
    ! XRTMIN = Minimum value for the mixing ratio
    ! ZRTMIN = Minimum value for the source (tendency)
    !
    INTEGER , DIMENSION(SIZE(GMICRO)) :: I1,I2,I3 ! Used to replace the COUNT
    INTEGER                           :: JL       ! and PACK intrinsics
    !
    !-------------------------------------------------------------------------------
    !
    !*       1.     COMPUTE THE LOOP BOUNDS
    !   	        -----------------------
    !
    IIB=1+JPHEXT
    IIE=SIZE(PDZZ,1) - JPHEXT
    IJB=1+JPHEXT
    IJE=SIZE(PDZZ,2) - JPHEXT
    IKB=KKA+JPVEXT*KKL
    IKE=KKU-JPVEXT*KKL
    IKT=SIZE(PDZZ,3)          
    IKTB=1+JPVEXT              
    IKTE=IKT-JPVEXT
    
    !
    !
    ZINVTSTEP=1./PTSTEP
    !
    !
    !*       2.     COMPUTES THE SLOW COLD PROCESS SOURCES
    !   	        --------------------------------------
    !
    CALL RAIN_ICE_NUCLEATION
    !
    !
    !  optimization by looking for locations where
    !  the microphysical fields are larger than a minimal value only !!!
    !
    GMICRO(:,:,:) = .FALSE.
    
     IF ( KRR == 7 ) THEN
      GMICRO(IIB:IIE,IJB:IJE,IKTB:IKTE) =                          &
                    PRCT(IIB:IIE,IJB:IJE,IKTB:IKTE)>XRTMIN(2) .OR. &
                    PRRT(IIB:IIE,IJB:IJE,IKTB:IKTE)>XRTMIN(3) .OR. &
                    PRIT(IIB:IIE,IJB:IJE,IKTB:IKTE)>XRTMIN(4) .OR. &
                    PRST(IIB:IIE,IJB:IJE,IKTB:IKTE)>XRTMIN(5) .OR. &
                    PRGT(IIB:IIE,IJB:IJE,IKTB:IKTE)>XRTMIN(6) .OR. &
                    PRHT(IIB:IIE,IJB:IJE,IKTB:IKTE)>XRTMIN(7)
     ELSE IF( KRR == 6 ) THEN
      GMICRO(IIB:IIE,IJB:IJE,IKTB:IKTE) =                          &
                    PRCT(IIB:IIE,IJB:IJE,IKTB:IKTE)>XRTMIN(2) .OR. &
                    PRRT(IIB:IIE,IJB:IJE,IKTB:IKTE)>XRTMIN(3) .OR. &
                    PRIT(IIB:IIE,IJB:IJE,IKTB:IKTE)>XRTMIN(4) .OR. &
                    PRST(IIB:IIE,IJB:IJE,IKTB:IKTE)>XRTMIN(5) .OR. &
                    PRGT(IIB:IIE,IJB:IJE,IKTB:IKTE)>XRTMIN(6)
     END IF
    
    IMICRO = COUNTJV( GMICRO(:,:,:),I1(:),I2(:),I3(:))
    IF( IMICRO >= 0 ) THEN
      ALLOCATE(ZRVT(IMICRO)) 
      ALLOCATE(ZRCT(IMICRO))
      ALLOCATE(ZRRT(IMICRO))
      ALLOCATE(ZRIT(IMICRO))
      ALLOCATE(ZRST(IMICRO))
      ALLOCATE(ZRGT(IMICRO))
     IF ( KRR == 7 ) ALLOCATE(ZRHT(IMICRO))
      ALLOCATE(ZCIT(IMICRO))
      ALLOCATE(ZRVS(IMICRO))
      ALLOCATE(ZRCS(IMICRO))
      ALLOCATE(ZRRS(IMICRO))
      ALLOCATE(ZRIS(IMICRO))
      ALLOCATE(ZRSS(IMICRO))
      ALLOCATE(ZRGS(IMICRO))
     IF ( KRR == 7 ) ALLOCATE(ZRHS(IMICRO))
      ALLOCATE(ZTHS(IMICRO))
      ALLOCATE(ZRHODREF(IMICRO))
      ALLOCATE(ZZT(IMICRO))
      ALLOCATE(ZPRES(IMICRO))
      ALLOCATE(ZEXNREF(IMICRO))
      ALLOCATE(ZSIGMA_RC(IMICRO))
      ALLOCATE(ZCF(IMICRO))
      DO JL=1,IMICRO   
        ZRVT(JL) = PRVT(I1(JL),I2(JL),I3(JL))
        ZRCT(JL) = PRCT(I1(JL),I2(JL),I3(JL))
        ZRRT(JL) = PRRT(I1(JL),I2(JL),I3(JL))
        ZRIT(JL) = PRIT(I1(JL),I2(JL),I3(JL))
        ZRST(JL) = PRST(I1(JL),I2(JL),I3(JL))
        ZRGT(JL) = PRGT(I1(JL),I2(JL),I3(JL))
        IF ( KRR == 7 ) ZRHT(JL) = PRHT(I1(JL),I2(JL),I3(JL))
        ZCIT(JL) = PCIT(I1(JL),I2(JL),I3(JL))
        IF ( HSUBG_AUCV == 'SIGM') THEN
             ZSIGMA_RC(JL) = PSIGS(I1(JL),I2(JL),I3(JL)) * 2.
          ELSE IF ( HSUBG_AUCV == 'CLFR') THEN
             ZCF(JL) = PCLDFR(I1(JL),I2(JL),I3(JL))
        END IF
    !
        ZRVS(JL) = PRVS(I1(JL),I2(JL),I3(JL))
        ZRCS(JL) = PRCS(I1(JL),I2(JL),I3(JL))
        ZRRS(JL) = PRRS(I1(JL),I2(JL),I3(JL))
        ZRIS(JL) = PRIS(I1(JL),I2(JL),I3(JL))
        ZRSS(JL) = PRSS(I1(JL),I2(JL),I3(JL))
        ZRGS(JL) = PRGS(I1(JL),I2(JL),I3(JL))
        IF ( KRR == 7 ) ZRHS(JL) = PRHS(I1(JL),I2(JL),I3(JL))
        ZTHS(JL) = PTHS(I1(JL),I2(JL),I3(JL))
    !
        ZRHODREF(JL) = PRHODREF(I1(JL),I2(JL),I3(JL))
        ZZT(JL) = ZT(I1(JL),I2(JL),I3(JL))
        ZPRES(JL) = PPABST(I1(JL),I2(JL),I3(JL))
        ZEXNREF(JL) = PEXNREF(I1(JL),I2(JL),I3(JL))
      ENDDO
      ALLOCATE(ZZW(IMICRO))
      ALLOCATE(ZLSFACT(IMICRO))
      ALLOCATE(ZLVFACT(IMICRO))
        ZZW(:)  = ZEXNREF(:)*( XCPD+XCPV*ZRVT(:)+XCL*(ZRCT(:)+ZRRT(:)) &
                                        +XCI*(ZRIT(:)+ZRST(:)+ZRGT(:)) )
        ZLSFACT(:) = (XLSTT+(XCPV-XCI)*(ZZT(:)-XTT))/ZZW(:) ! L_s/(Pi_ref*C_ph)
        ZLVFACT(:) = (XLVTT+(XCPV-XCL)*(ZZT(:)-XTT))/ZZW(:) ! L_v/(Pi_ref*C_ph)
      ALLOCATE(ZUSW(IMICRO))
      ALLOCATE(ZSSI(IMICRO))
        ZZW(:) = EXP( XALPI - XBETAI/ZZT(:) - XGAMI*ALOG(ZZT(:) ) )
        ZSSI(:) = ZRVT(:)*( ZPRES(:)-ZZW(:) ) / ( (XMV/XMD) * ZZW(:) ) - 1.0
                                                          ! Supersaturation over ice
    !
      ALLOCATE(ZLBDAR(IMICRO))
      ALLOCATE(ZLBDAS(IMICRO))
      ALLOCATE(ZLBDAG(IMICRO))
     IF ( KRR == 7 ) ALLOCATE(ZLBDAH(IMICRO))
      ALLOCATE(ZRDRYG(IMICRO))
      ALLOCATE(ZRWETG(IMICRO))
      ALLOCATE(ZAI(IMICRO))
      ALLOCATE(ZCJ(IMICRO))
      ALLOCATE(ZKA(IMICRO))
      ALLOCATE(ZDV(IMICRO))
    !
     IF ( KRR == 7 ) THEN
      ALLOCATE(ZZW1(IMICRO,7))
        ELSE IF( KRR == 6 ) THEN
            ALLOCATE(ZZW1(IMICRO,6))
     ENDIF
    !
    
      IF (LBU_ENABLE .OR. LLES_CALL .OR. LCHECK) THEN
    
        ALLOCATE(ZRHODJ(IMICRO))
        ZRHODJ(:) = PACK( PRHODJ(:,:,:),MASK=GMICRO(:,:,:) )
      END IF
    !
      CALL RAIN_ICE_SLOW
    !
    !-------------------------------------------------------------------------------
    !
    !
    !*       3.     COMPUTES THE SLOW WARM PROCESS SOURCES
    !   	        --------------------------------------
    !
    !*       3.1    compute the slope parameter Lbda_r
    !
      WHERE( ZRRT(:)>0.0 )
        ZLBDAR(:)  = XLBR*( ZRHODREF(:)*MAX( ZRRT(:),XRTMIN(3) ) )**XLBEXR
      END WHERE
    !
      IF( OWARM ) THEN    !  Check if the formation of the raindrops by the slow
                          !  warm processes is allowed
        PEVAP3D(:,:,:)= 0.
        CALL RAIN_ICE_WARM
      END IF
    !
    !-------------------------------------------------------------------------------
    !
    !
    !*       4.     COMPUTES THE FAST COLD PROCESS SOURCES FOR r_s
    !               ----------------------------------------------
    !
      CALL RAIN_ICE_FAST_RS
    !
    !-------------------------------------------------------------------------------
    !
    !
    !*       5.     COMPUTES THE FAST COLD PROCESS SOURCES FOR r_g
    !               ----------------------------------------------
    !
      CALL RAIN_ICE_FAST_RG
    !
    !-------------------------------------------------------------------------------
    !
    !
    !*       6.     COMPUTES THE FAST COLD PROCESS SOURCES FOR r_h
    !               ----------------------------------------------
    !
     IF ( KRR == 7 ) THEN
      CALL RAIN_ICE_FAST_RH
     END IF
    !
    !-------------------------------------------------------------------------------
    !
    !
    !*       7.     COMPUTES SPECIFIC SOURCES OF THE WARM AND COLD CLOUDY SPECIES
    !               -------------------------------------------------------------
    !
      CALL RAIN_ICE_FAST_RI
    !
    !
    !-------------------------------------------------------------------------------
    !
    !
    !
      ZW(:,:,:) = PRVS(:,:,:)
      PRVS(:,:,:) = UNPACK( ZRVS(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:) )
      ZW(:,:,:) = PRCS(:,:,:)
      PRCS(:,:,:) = UNPACK( ZRCS(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:) )
      ZW(:,:,:) = PRRS(:,:,:)
      PRRS(:,:,:) = UNPACK( ZRRS(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:) )
      ZW(:,:,:) = PRIS(:,:,:)
      PRIS(:,:,:) = UNPACK( ZRIS(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:) )
      ZW(:,:,:) = PRSS(:,:,:)
      PRSS(:,:,:) = UNPACK( ZRSS(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:) )
      ZW(:,:,:) = PRGS(:,:,:)
      PRGS(:,:,:) = UNPACK( ZRGS(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:) )
      IF ( KRR == 7 ) THEN
        ZW(:,:,:) = PRHS(:,:,:)
        PRHS(:,:,:) = UNPACK( ZRHS(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:) )
      END IF
      ZW(:,:,:) = PTHS(:,:,:)
      PTHS(:,:,:) = UNPACK( ZTHS(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:) )
      ZW(:,:,:) = PCIT(:,:,:)
      PCIT(:,:,:) = UNPACK( ZCIT(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:) )
    !
    !
    !
      DEALLOCATE(ZZW1)
      DEALLOCATE(ZDV)
      DEALLOCATE(ZCJ)
      DEALLOCATE(ZRDRYG)
      DEALLOCATE(ZRWETG)
      DEALLOCATE(ZLBDAG)
     IF ( KRR == 7 ) DEALLOCATE(ZLBDAH)
      DEALLOCATE(ZLBDAS)
      DEALLOCATE(ZLBDAR)
      DEALLOCATE(ZSSI)
      DEALLOCATE(ZUSW)
      DEALLOCATE(ZLVFACT)
      DEALLOCATE(ZLSFACT)
      DEALLOCATE(ZZW)
      DEALLOCATE(ZEXNREF)
      DEALLOCATE(ZPRES)
      DEALLOCATE(ZRHODREF)
      DEALLOCATE(ZZT)
      IF(LBU_ENABLE .OR. LLES_CALL) DEALLOCATE(ZRHODJ)
      DEALLOCATE(ZTHS)
      IF ( KRR == 7 ) DEALLOCATE(ZRHS)
      DEALLOCATE(ZRGS)
      DEALLOCATE(ZRSS)
      DEALLOCATE(ZRIS)
      DEALLOCATE(ZRRS)
      DEALLOCATE(ZRCS)
      DEALLOCATE(ZRVS)
      DEALLOCATE(ZCIT)
      DEALLOCATE(ZRGT)
      IF ( KRR == 7 ) DEALLOCATE(ZRHT)
      DEALLOCATE(ZRST)
      DEALLOCATE(ZRIT)
      DEALLOCATE(ZRRT)
      DEALLOCATE(ZAI)
      DEALLOCATE(ZRCT)
      DEALLOCATE(ZKA)
      DEALLOCATE(ZRVT)
      DEALLOCATE(ZSIGMA_RC)
      DEALLOCATE(ZCF)
    !
      ELSE
    !
    ! Advance the budget calls
    !
    ! Reordered for compability with flexible structures like in AROME
    
     ! rain_ice_slow
     IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'HON_BU_RTH')
     IF (LBUDGET_RC) CALL BUDGET (PRCS(:,:,:)*PRHODJ(:,:,:),7,'HON_BU_RRC')
     IF (LBUDGET_RI) CALL BUDGET (PRIS(:,:,:)*PRHODJ(:,:,:),9,'HON_BU_RRI')
     IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'SFR_BU_RTH')
     IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:)*PRHODJ(:,:,:),8,'SFR_BU_RRR')
     IF (LBUDGET_RG) CALL BUDGET (PRGS(:,:,:)*PRHODJ(:,:,:),11,'SFR_BU_RRG')
     IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'DEPS_BU_RTH')
     IF (LBUDGET_RV) CALL BUDGET (PRVS(:,:,:)*PRHODJ(:,:,:),6,'DEPS_BU_RRV')
     IF (LBUDGET_RS) CALL BUDGET (PRSS(:,:,:)*PRHODJ(:,:,:),10,'DEPS_BU_RRS')
     IF (LBUDGET_RI) CALL BUDGET (PRIS(:,:,:)*PRHODJ(:,:,:),9,'AGGS_BU_RRI')
     IF (LBUDGET_RS) CALL BUDGET (PRSS(:,:,:)*PRHODJ(:,:,:),10,'AGGS_BU_RRS')
     IF (LBUDGET_RI) CALL BUDGET (PRIS(:,:,:)*PRHODJ(:,:,:),9,'AUTS_BU_RRI')
     IF (LBUDGET_RS) CALL BUDGET (PRSS(:,:,:)*PRHODJ(:,:,:),10,'AUTS_BU_RRS')
     IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'DEPG_BU_RTH')
     IF (LBUDGET_RV) CALL BUDGET (PRVS(:,:,:)*PRHODJ(:,:,:),6,'DEPG_BU_RRV')
     IF (LBUDGET_RG) CALL BUDGET (PRGS(:,:,:)*PRHODJ(:,:,:),11,'DEPG_BU_RRG')
    
     IF (OWARM) THEN ! rain_ice_warm
       IF (LBUDGET_RC) CALL BUDGET (PRCS(:,:,:)*PRHODJ(:,:,:),7,'AUTO_BU_RRC')
       IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:)*PRHODJ(:,:,:),8,'AUTO_BU_RRR')
       IF (LBUDGET_RC) CALL BUDGET (PRCS(:,:,:)*PRHODJ(:,:,:),7,'ACCR_BU_RRC')
       IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:)*PRHODJ(:,:,:),8,'ACCR_BU_RRR')
       IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'REVA_BU_RTH')
       IF (LBUDGET_RV) CALL BUDGET (PRVS(:,:,:)*PRHODJ(:,:,:),6,'REVA_BU_RRV')
       IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:)*PRHODJ(:,:,:),8,'REVA_BU_RRR')
     ENDIF
    
     !rain_ice_fast_rs
     IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'RIM_BU_RTH')
     IF (LBUDGET_RC) CALL BUDGET (PRCS(:,:,:)*PRHODJ(:,:,:),7,'RIM_BU_RRC')
     IF (LBUDGET_RS) CALL BUDGET (PRSS(:,:,:)*PRHODJ(:,:,:),10,'RIM_BU_RRS')
     IF (LBUDGET_RG) CALL BUDGET (PRGS(:,:,:)*PRHODJ(:,:,:),11,'RIM_BU_RRG')
     IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'ACC_BU_RTH')
     IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:)*PRHODJ(:,:,:),8,'ACC_BU_RRR')
     IF (LBUDGET_RS) CALL BUDGET (PRSS(:,:,:)*PRHODJ(:,:,:),10,'ACC_BU_RRS')
     IF (LBUDGET_RG) CALL BUDGET (PRGS(:,:,:)*PRHODJ(:,:,:),11,'ACC_BU_RRG')
     IF (LBUDGET_RS) CALL BUDGET (PRSS(:,:,:)*PRHODJ(:,:,:),10,'CMEL_BU_RRS')
     IF (LBUDGET_RG) CALL BUDGET (PRGS(:,:,:)*PRHODJ(:,:,:),11,'CMEL_BU_RRG')
    
     !rain_ice_fast_rg
     IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'CFRZ_BU_RTH')
     IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:)*PRHODJ(:,:,:),8,'CFRZ_BU_RRR')
     IF (LBUDGET_RI) CALL BUDGET (PRIS(:,:,:)*PRHODJ(:,:,:),9,'CFRZ_BU_RRI')
     IF (LBUDGET_RG) CALL BUDGET (PRGS(:,:,:)*PRHODJ(:,:,:),11,'CFRZ_BU_RRG')
     IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'WETG_BU_RTH')
     IF (LBUDGET_RC) CALL BUDGET (PRCS(:,:,:)*PRHODJ(:,:,:),7,'WETG_BU_RRC')
     IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:)*PRHODJ(:,:,:),8,'WETG_BU_RRR')
     IF (LBUDGET_RI) CALL BUDGET (PRIS(:,:,:)*PRHODJ(:,:,:),9,'WETG_BU_RRI')
     IF (LBUDGET_RS) CALL BUDGET (PRSS(:,:,:)*PRHODJ(:,:,:),10,'WETG_BU_RRS')
     IF (LBUDGET_RG) CALL BUDGET (PRGS(:,:,:)*PRHODJ(:,:,:),11,'WETG_BU_RRG')
     IF (LBUDGET_RH) CALL BUDGET (PRHS(:,:,:)*PRHODJ(:,:,:),12,'WETG_BU_RRH')
     IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'DRYG_BU_RTH')
     IF (LBUDGET_RC) CALL BUDGET (PRCS(:,:,:)*PRHODJ(:,:,:),7,'DRYG_BU_RRC')
     IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:)*PRHODJ(:,:,:),8,'DRYG_BU_RRR')
     IF (LBUDGET_RI) CALL BUDGET (PRIS(:,:,:)*PRHODJ(:,:,:),9,'DRYG_BU_RRI')
     IF (LBUDGET_RS) CALL BUDGET (PRSS(:,:,:)*PRHODJ(:,:,:),10,'DRYG_BU_RRS')
     IF (LBUDGET_RG) CALL BUDGET (PRGS(:,:,:)*PRHODJ(:,:,:),11,'DRYG_BU_RRG')
     IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'GMLT_BU_RTH')
     IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:)*PRHODJ(:,:,:),8,'GMLT_BU_RRR')
     IF (LBUDGET_RG) CALL BUDGET (PRGS(:,:,:)*PRHODJ(:,:,:),11,'GMLT_BU_RRG')
    
     IF(KRR==7) THEN ! rain_ice_fast_rh
       IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'WETH_BU_RTH')
       IF (LBUDGET_RC) CALL BUDGET (PRCS(:,:,:)*PRHODJ(:,:,:),7,'WETH_BU_RRC')
       IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:)*PRHODJ(:,:,:),8,'WETH_BU_RRR')
       IF (LBUDGET_RI) CALL BUDGET (PRIS(:,:,:)*PRHODJ(:,:,:),9,'WETH_BU_RRI')
       IF (LBUDGET_RS) CALL BUDGET (PRSS(:,:,:)*PRHODJ(:,:,:),10,'WETH_BU_RRS')
       IF (LBUDGET_RG) CALL BUDGET (PRGS(:,:,:)*PRHODJ(:,:,:),11,'WETH_BU_RRG')
       IF (LBUDGET_RH) CALL BUDGET (PRHS(:,:,:)*PRHODJ(:,:,:),12,'WETH_BU_RRH')
    
       IF (LBUDGET_RH) CALL BUDGET (PRGS(:,:,:)*PRHODJ(:,:,:),11,'COHG_BU_RRG')
       IF (LBUDGET_RH) CALL BUDGET (PRHS(:,:,:)*PRHODJ(:,:,:),12,'COHG_BU_RRH')
    
       IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'HMLT_BU_RTH')
       IF (LBUDGET_RR) CALL BUDGET (PRRS(:,:,:)*PRHODJ(:,:,:),8,'HMLT_BU_RRR')
       IF (LBUDGET_RH) CALL BUDGET (PRHS(:,:,:)*PRHODJ(:,:,:),12,'HMLT_BU_RRH')
     ENDIF
    
     !rain_ice_fast_ri
     IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'IMLT_BU_RTH')
     IF (LBUDGET_RC) CALL BUDGET (PRCS(:,:,:)*PRHODJ(:,:,:),7,'IMLT_BU_RRC')
     IF (LBUDGET_RI) CALL BUDGET (PRIS(:,:,:)*PRHODJ(:,:,:),9,'IMLT_BU_RRI')
     IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:)*PRHODJ(:,:,:),4,'BERFI_BU_RTH')
     IF (LBUDGET_RC) CALL BUDGET (PRCS(:,:,:)*PRHODJ(:,:,:),7,'BERFI_BU_RRC')
     IF (LBUDGET_RI) CALL BUDGET (PRIS(:,:,:)*PRHODJ(:,:,:),9,'BERFI_BU_RRI')
    !
    END IF
    !
    !-------------------------------------------------------------------------------
    !
    !*       8.     COMPUTE THE SEDIMENTATION (RS) SOURCE
    !	        -------------------------------------
    !
    !*       8.1    time splitting loop initialization
    !
    ZTSPLITR= PTSTEP / FLOAT(KSPLITR)
    !
    !
    IF (HSEDIM == 'STAT') THEN
      CALL RAIN_ICE_SEDIMENTATION_STAT
    ELSEIF (HSEDIM == 'SPLI') THEN
      CALL RAIN_ICE_SEDIMENTATION_SPLIT
    ELSE
      WRITE(*,*) ' STOP'                                                     
      WRITE(*,*) ' NO SEDIMENTATION SCHEME FOR HSEDIM=',HSEDIM 
      CALL ABORT
      STOP
    END IF
    
    !
    !
    !-------------------------------------------------------------------------------
    !
    !-------------------------------------------------------------------------------
    !
    !
    CONTAINS
    !
    !
    !-------------------------------------------------------------------------------
    !
    !
      SUBROUTINE RAIN_ICE_SEDIMENTATION_SPLIT
    !
    !*      0. DECLARATIONS
    !          ------------
    !
    IMPLICIT NONE
    !
    !*       0.2  declaration of local variables
    !
    !
    INTEGER , DIMENSION(SIZE(GSEDIMC)) :: IC1,IC2,IC3 ! Used to replace the COUNT
    INTEGER , DIMENSION(SIZE(GSEDIMR)) :: IR1,IR2,IR3 ! Used to replace the COUNT
    INTEGER , DIMENSION(SIZE(GSEDIMS)) :: IS1,IS2,IS3 ! Used to replace the COUNT
    INTEGER , DIMENSION(SIZE(GSEDIMI)) :: II1,II2,II3 ! Used to replace the COUNT
    INTEGER , DIMENSION(SIZE(GSEDIMG)) :: IG1,IG2,IG3 ! Used to replace the COUNT
    INTEGER , DIMENSION(SIZE(GSEDIMH)) :: IH1,IH2,IH3 ! Used to replace the COUNT
    INTEGER   :: ILENALLOCC,ILENALLOCR,ILENALLOCI,ILENALLOCS,ILENALLOCG,ILENALLOCH
    INTEGER   :: ILISTLENC,ILISTLENR,ILISTLENI,ILISTLENS,ILISTLENG,ILISTLENH
    INTEGER, ALLOCATABLE :: ILISTR(:),ILISTC(:),ILISTI(:),ILISTS(:),ILISTG(:),ILISTH(:)
    ! Optimization for NEC
    !INTEGER, SAVE :: IOLDALLOCC = SIZE(PEXNREF,1)*SIZE(PEXNREF,2)*SIZE(PEXNREF,3)/10
    !INTEGER, SAVE :: IOLDALLOCR = SIZE(PEXNREF,1)*SIZE(PEXNREF,2)*SIZE(PEXNREF,3)/10
    !INTEGER, SAVE :: IOLDALLOCI = SIZE(PEXNREF,1)*SIZE(PEXNREF,2)*SIZE(PEXNREF,3)/10
    !INTEGER, SAVE :: IOLDALLOCS = SIZE(PEXNREF,1)*SIZE(PEXNREF,2)*SIZE(PEXNREF,3)/10
    !INTEGER, SAVE :: IOLDALLOCG = SIZE(PEXNREF,1)*SIZE(PEXNREF,2)*SIZE(PEXNREF,3)/10
    !INTEGER, SAVE :: IOLDALLOCH = SIZE(PEXNREF,1)*SIZE(PEXNREF,2)*SIZE(PEXNREF,3)/10
    INTEGER, SAVE :: IOLDALLOCC = 6000
    INTEGER, SAVE :: IOLDALLOCR = 6000
    INTEGER, SAVE :: IOLDALLOCI = 6000
    INTEGER, SAVE :: IOLDALLOCS = 6000
    INTEGER, SAVE :: IOLDALLOCG = 6000
    INTEGER, SAVE :: IOLDALLOCH = 6000
    !
    REAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2),SIZE(PRHODREF,3)) :: ZCONC3D !  droplet condensation
    !-------------------------------------------------------------------------------
    !
    !
    !        O. Initialization of for sedimentation                  
    !
    IF (OSEDIC) PINPRC (:,:) = 0.
    
    PINPRR (:,:) = 0.
    PINPRR3D (:,:,:) = 0.
    PINPRS (:,:) = 0.
    
    PINPRG (:,:) = 0.
    
    PINPRG3D (:,:,:) = 0.
    IF ( KRR == 7 ) PINPRH (:,:) = 0.   
    IF ( KRR == 7 ) PINPRH3D (:,:,:) = 0.
    PSPEEDC(:,:,:) = 0.
    PSPEEDR(:,:,:) = 0.
    PSPEEDS(:,:,:) = 0.
    PSPEEDG(:,:,:) = 0.
    PSPEEDH(:,:,:) = 0.
    
    !
    !*       1. Parameters for cloud sedimentation
    !
       IF (OSEDIC) THEN
        ZRAY(:,:,:)   = 0.
        ZLBC(:,:,:)   = XLBC(1)
        ZFSEDC(:,:,:) = XFSEDC(1)
        ZCONC3D(:,:,:)= XCONC_LAND
        ZCONC_TMP(:,:)= XCONC_LAND
        IF (PRESENT(PSEA)) THEN
          ZCONC_TMP(:,:)=PSEA(:,:)*XCONC_SEA+(1.-PSEA(:,:))*XCONC_LAND
          DO JK=IKTB,IKTE
            ZLBC(:,:,JK)   = PSEA(:,:)*XLBC(2)+(1.-PSEA(:,:))*XLBC(1)
            ZFSEDC(:,:,JK) = (PSEA(:,:)*XFSEDC(2)+(1.-PSEA(:,:))*XFSEDC(1))
            ZFSEDC(:,:,JK) = MAX(MIN(XFSEDC(1),XFSEDC(2)),ZFSEDC(:,:,JK))
            ZCONC3D(:,:,JK)= (1.-PTOWN(:,:))*ZCONC_TMP(:,:)+PTOWN(:,:)*XCONC_URBAN
            ZRAY(:,:,JK)   = 0.5*((1.-PSEA(:,:))*GAMMA(XNUC+1.0/XALPHAC)/(GAMMA(XNUC)) + &
                    PSEA(:,:)*GAMMA(XNUC2+1.0/XALPHAC2)/(GAMMA(XNUC2)))
          END DO
        ELSE
            ZCONC3D(:,:,:) = XCONC_LAND
            ZRAY(:,:,:)  = 0.5*(GAMMA(XNUC+1.0/XALPHAC)/(GAMMA(XNUC)))
        END IF
        ZRAY(:,:,:)      = MAX(1.,ZRAY(:,:,:))
        ZLBC(:,:,:)      = MAX(MIN(XLBC(1),XLBC(2)),ZLBC(:,:,:))
       ENDIF
    !
    !*       2.    compute the fluxes
    !
    !  optimization by looking for locations where
    !  the precipitating fields are larger than a minimal value only !!!
    !  For optimization we consider each variable separately
    
    ZRTMIN(:)    = XRTMIN(:) * ZINVTSTEP
    IF (OSEDIC) GSEDIMC(:,:,:) = .FALSE.
    GSEDIMR(:,:,:) = .FALSE.
    GSEDIMI(:,:,:) = .FALSE.
    GSEDIMS(:,:,:) = .FALSE.
    GSEDIMG(:,:,:) = .FALSE.
    IF ( KRR == 7 ) GSEDIMH(:,:,:) = .FALSE.
    !
    ILENALLOCR = 0
    IF (OSEDIC) ILENALLOCC = 0
    ILENALLOCI = 0
    ILENALLOCS = 0
    ILENALLOCG = 0
    IF ( KRR == 7 ) ILENALLOCH = 0
    !
    ! ZPiS = Specie i source creating during the current time step
    ! PRiS = Source of the previous time step
    !
    IF (OSEDIC) THEN
      ZPRCS(:,:,:) = 0.0
    
      ZPRCS(:,:,:) = PRCS(:,:,:)-PRCT(:,:,:)* ZINVTSTEP
      PRCS(:,:,:)  = PRCT(:,:,:)* ZINVTSTEP
    
    END IF
    ZPRRS(:,:,:) = 0.0
    ZPRSS(:,:,:) = 0.0
    ZPRGS(:,:,:) = 0.0
    IF ( KRR == 7 ) ZPRHS(:,:,:) = 0.0
    !
    
    ZPRRS(:,:,:) = PRRS(:,:,:)-PRRT(:,:,:)* ZINVTSTEP
    ZPRSS(:,:,:) = PRSS(:,:,:)-PRST(:,:,:)* ZINVTSTEP
    ZPRGS(:,:,:) = PRGS(:,:,:)-PRGT(:,:,:)* ZINVTSTEP
    IF ( KRR == 7 ) ZPRHS(:,:,:) = PRHS(:,:,:)-PRHT(:,:,:)* ZINVTSTEP
    PRRS(:,:,:)  = PRRT(:,:,:)* ZINVTSTEP
    PRSS(:,:,:)  = PRST(:,:,:)* ZINVTSTEP
    PRGS(:,:,:)  = PRGT(:,:,:)* ZINVTSTEP
    IF ( KRR == 7 ) PRHS(:,:,:)  = PRHT(:,:,:)* ZINVTSTEP
    
    !
    ! PRiS = Source of the previous time step + source created during the subtime
    ! step
    !
    DO JN = 1 , KSPLITR
      IF( JN==1 ) THEN
       IF (OSEDIC) PRCS(:,:,:) = PRCS(:,:,:) + ZPRCS(:,:,:)/KSPLITR
       PRRS(:,:,:) = PRRS(:,:,:) + ZPRRS(:,:,:)/KSPLITR
       PRSS(:,:,:) = PRSS(:,:,:) + ZPRSS(:,:,:)/KSPLITR
       PRGS(:,:,:) = PRGS(:,:,:) + ZPRGS(:,:,:)/KSPLITR
       IF ( KRR == 7 ) PRHS(:,:,:) = PRHS(:,:,:) + ZPRHS(:,:,:)/KSPLITR
       DO JK = IKTB , IKTE
         ZW(:,:,JK) =ZTSPLITR/(PRHODREF(:,:,JK)* PDZZ(:,:,JK))
       END DO
     ELSE
       IF (OSEDIC) PRCS(:,:,:) = PRCS(:,:,:) + ZPRCS(:,:,:)*ZTSPLITR
       PRRS(:,:,:) = PRRS(:,:,:) + ZPRRS(:,:,:)*ZTSPLITR
       PRSS(:,:,:) = PRSS(:,:,:) + ZPRSS(:,:,:)*ZTSPLITR
       PRGS(:,:,:) = PRGS(:,:,:) + ZPRGS(:,:,:)*ZTSPLITR
       IF ( KRR == 7 ) PRHS(:,:,:) = PRHS(:,:,:) + ZPRHS(:,:,:)*ZTSPLITR
     END IF
     !
     IF (OSEDIC) GSEDIMC(IIB:IIE,IJB:IJE,IKTB:IKTE) =                &
                      PRCS(IIB:IIE,IJB:IJE,IKTB:IKTE)>ZRTMIN(2)
     GSEDIMR(IIB:IIE,IJB:IJE,IKTB:IKTE) =                            &
                      PRRS(IIB:IIE,IJB:IJE,IKTB:IKTE)>ZRTMIN(3)
     GSEDIMI(IIB:IIE,IJB:IJE,IKTB:IKTE) =                            &
                      PRIS(IIB:IIE,IJB:IJE,IKTB:IKTE)>ZRTMIN(4)
     GSEDIMS(IIB:IIE,IJB:IJE,IKTB:IKTE) =                            &