From 9edba42112880661934599f311baccc2042f142c Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Beno=C3=AEt=20Vi=C3=A9?= <benoit.vie@meteo.fr>
Date: Thu, 10 Feb 2022 11:25:25 +0100
Subject: [PATCH] initial commit for LIMA_JW

---
 src/MNH/ini_ice_c1r3.f90                 |  40 +++----
 src/MNH/ini_lima_cold_mixed.f90          | 144 +++++++++++++++--------
 src/MNH/init_aerosol_concentration.f90   |   2 +-
 src/MNH/lima_cold_slow_processes.f90     |  46 ++++++--
 src/MNH/lima_conversion_melting_snow.f90 |  14 ++-
 src/MNH/lima_droplets_riming_snow.f90    |  14 ++-
 src/MNH/lima_graupel.f90                 |  11 +-
 src/MNH/lima_ice_aggregation_snow.f90    |   8 +-
 src/MNH/lima_ice_snow_deposition.f90     |  13 +-
 src/MNH/lima_mixed.f90                   |  21 +++-
 src/MNH/lima_mixed_fast_processes.f90    | 100 ++++++++++++----
 src/MNH/lima_rain_accr_snow.f90          |  22 +++-
 src/MNH/lima_sedimentation.f90           |  20 +++-
 src/MNH/lima_snow_deposition.f90         |  13 +-
 src/MNH/lima_tendencies.f90              |  17 ++-
 src/MNH/modd_param_lima.f90              |   4 +
 src/MNH/modd_param_lima_mixed.f90        |   1 +
 src/MNH/rrcolss.f90                      |  10 +-
 src/MNH/rscolrg.f90                      |  10 +-
 src/MNH/rzcolx.f90                       |  12 +-
 20 files changed, 362 insertions(+), 160 deletions(-)

diff --git a/src/MNH/ini_ice_c1r3.f90 b/src/MNH/ini_ice_c1r3.f90
index b0a35554a..771069f6b 100644
--- a/src/MNH/ini_ice_c1r3.f90
+++ b/src/MNH/ini_ice_c1r3.f90
@@ -724,18 +724,18 @@ IF( (KACCLBDAS/=NACCLBDAS) .OR. (KACCLBDAR/=NACCLBDAR) .OR. (KND/=IND) .OR. &
     (PACCLBDAS_MAX/=XACCLBDAS_MAX) .OR. (PACCLBDAR_MAX/=XACCLBDAR_MAX) .OR. &
     (PACCLBDAS_MIN/=XACCLBDAS_MIN) .OR. (PACCLBDAR_MIN/=XACCLBDAR_MIN) .OR. &
     (PFDINFTY/=ZFDINFTY)                                               ) THEN
-  CALL RRCOLSS ( IND, XALPHAS, XNUS, XALPHAR, XNUR,                          &
-                 ZESR, XBR, XCS, XDS, XCR, XDR,                              &
-                 XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
-                 ZFDINFTY, XKER_RACCSS, XAG, XBS, XAS                        )
-  CALL RZCOLX  ( IND, XALPHAS, XNUS, XALPHAR, XNUR,                          &
-                 ZESR, XBR, XCS, XDS, XCR, XDR,                              &
-                 XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
-                 ZFDINFTY, XKER_RACCS                                        )
-  CALL RSCOLRG ( IND, XALPHAS, XNUS, XALPHAR, XNUR,                          &
-                 ZESR, XBS, XCS, XDS, XCR, XDR,                              &
-                 XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
-                 ZFDINFTY, XKER_SACCRG,XAG, XBS, XAS                         )
+!  CALL RRCOLSS ( IND, XALPHAS, XNUS, XALPHAR, XNUR,                          &
+!                 ZESR, XBR, XCS, XDS, XCR, XDR,                              &
+!                 XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
+!                 ZFDINFTY, XKER_RACCSS, XAG, XBS, XAS                        )
+!  CALL RZCOLX  ( IND, XALPHAS, XNUS, XALPHAR, XNUR,                          &
+!                 ZESR, XBR, XCS, XDS, XCR, XDR,                              &
+!                 XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
+!                 ZFDINFTY, XKER_RACCS                                        )
+!  CALL RSCOLRG ( IND, XALPHAS, XNUS, XALPHAR, XNUR,                          &
+!                 ZESR, XBS, XCS, XDS, XCR, XDR,                              &
+!                 XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
+!                 ZFDINFTY, XKER_SACCRG,XAG, XBS, XAS                         )
   WRITE(UNIT=ILUOUT0,FMT='("*****************************************")')
   WRITE(UNIT=ILUOUT0,FMT='("**** UPDATE NEW SET OF RACSS KERNELS ****")')
   WRITE(UNIT=ILUOUT0,FMT='("**** UPDATE NEW SET OF RACS  KERNELS ****")')
@@ -937,10 +937,10 @@ IF( (KDRYLBDAG/=NDRYLBDAG) .OR. (KDRYLBDAS/=NDRYLBDAS) .OR. (KND/=IND) .OR. &
     (PDRYLBDAG_MAX/=XDRYLBDAG_MAX) .OR. (PDRYLBDAS_MAX/=XDRYLBDAS_MAX) .OR. &
     (PDRYLBDAG_MIN/=XDRYLBDAG_MIN) .OR. (PDRYLBDAS_MIN/=XDRYLBDAS_MIN) .OR. &
     (PFDINFTY/=ZFDINFTY)                                               ) THEN
-  CALL RZCOLX ( IND, XALPHAG, XNUG, XALPHAS, XNUS,                          &
-                ZEGS, XBS, XCG, XDG, XCS, XDS,                              &
-                XDRYLBDAG_MAX, XDRYLBDAS_MAX, XDRYLBDAG_MIN, XDRYLBDAS_MIN, &
-                ZFDINFTY, XKER_SDRYG                                        )
+!  CALL RZCOLX ( IND, XALPHAG, XNUG, XALPHAS, XNUS,                          &
+!                ZEGS, XBS, XCG, XDG, XCS, XDS,                              &
+!                XDRYLBDAG_MAX, XDRYLBDAS_MAX, XDRYLBDAG_MIN, XDRYLBDAS_MIN, &
+!                ZFDINFTY, XKER_SDRYG                                        )
   WRITE(UNIT=ILUOUT0,FMT='("*****************************************")')
   WRITE(UNIT=ILUOUT0,FMT='("**** UPDATE NEW SET OF SDRYG KERNELS ****")')
   WRITE(UNIT=ILUOUT0,FMT='("*****************************************")')
@@ -1003,10 +1003,10 @@ IF( (KDRYLBDAG/=NDRYLBDAG) .OR. (KDRYLBDAR/=NDRYLBDAR) .OR. (KND/=IND) .OR. &
     (PDRYLBDAG_MAX/=XDRYLBDAG_MAX) .OR. (PDRYLBDAR_MAX/=XDRYLBDAR_MAX) .OR. &
     (PDRYLBDAG_MIN/=XDRYLBDAG_MIN) .OR. (PDRYLBDAR_MIN/=XDRYLBDAR_MIN) .OR. &
     (PFDINFTY/=ZFDINFTY)                                               ) THEN
-  CALL RZCOLX ( IND, XALPHAG, XNUG, XALPHAR, XNUR,                          &
-                ZEGR, XBR, XCG, XDG, XCR, XDR,                              &
-                XDRYLBDAG_MAX, XDRYLBDAR_MAX, XDRYLBDAG_MIN, XDRYLBDAR_MIN, &
-                ZFDINFTY, XKER_RDRYG                                        )
+!  CALL RZCOLX ( IND, XALPHAG, XNUG, XALPHAR, XNUR,                          &
+!                ZEGR, XBR, XCG, XDG, XCR, XDR,                              &
+!                XDRYLBDAG_MAX, XDRYLBDAR_MAX, XDRYLBDAG_MIN, XDRYLBDAR_MIN, &
+!                ZFDINFTY, XKER_RDRYG                                        )
   WRITE(UNIT=ILUOUT0,FMT='("*****************************************")')
   WRITE(UNIT=ILUOUT0,FMT='("**** UPDATE NEW SET OF RDRYG KERNELS ****")')
   WRITE(UNIT=ILUOUT0,FMT='("*****************************************")')
diff --git a/src/MNH/ini_lima_cold_mixed.f90 b/src/MNH/ini_lima_cold_mixed.f90
index cb427cdb4..defb832c1 100644
--- a/src/MNH/ini_lima_cold_mixed.f90
+++ b/src/MNH/ini_lima_cold_mixed.f90
@@ -171,12 +171,27 @@ XF1IS = 0.28
 !
 XAS = 0.02
 XBS = 1.9
-XCS = 5.
-XDS = 0.27
-!
-XCCS = 5.0
-XCXS = 1.0
-!
+
+!OPER
+!XCS = 5.      ! Wurtz
+!XDS = 0.27    ! Wurtz
+!XFVELOS = 0.  ! Wurtz
+
+!Cas GAMMAGEN
+XCS = 11.52     ! Wurtz
+XDS = 0.39      ! Wurtz
+XFVELOS =0.097  ! Wurtz
+
+!Cas MP
+!XCS = 13.2        ! Wurtz
+!XDS = 0.423       ! Wurtz        
+!XFVELOS = 25.14   ! Wurtz
+
+!
+!XCCS = 5.0 ! Wurtz unused
+!XCXS = 1.0 ! Wurtz unused
+!
+
 XF0S = 0.86
 XF1S = 0.28
 !
@@ -234,8 +249,17 @@ XC1H = 1./2.
 XALPHAI = 3.0  ! Gamma law for the ice crystal volume
 XNUI    = 3.0  ! Gamma law with little dispersion
 !
-XALPHAS = 1.0  ! Exponential law
-XNUS    = 1.0  ! Exponential law
+!Wurtz
+!Cas MP
+!XALPHAS = 1.0  ! Exponential law
+!XNUS    = 1.0  ! Exponential law
+
+!Wurtz
+!Cas GAMMAGEN
+XALPHAS = .214   ! Generalized gamma law    !Wurtz
+XNUS    = 43.7   ! Generalized gamma law    !Wurtz
+XTRANS_MP_GAMMAS = SQRT( ( GAMMA(XNUS + 2./XALPHAS)*GAMMA(XNUS + 4./XALPHAS) ) / ( 8.* GAMMA(XNUS + 1./XALPHAS)*GAMMA(XNUS + 3./XALPHAS) ) )    ! Wurtz
+
 !
 XALPHAG = 1.0  ! Exponential law
 XNUG    = 1.0  ! Exponential law
@@ -248,8 +272,10 @@ XNUH    = 8.0  ! Gamma law with little dispersion
 XLBEXI = 1.0/XBI
 XLBI   = XAI*MOMG(XALPHAI,XNUI,XBI)
 !
-XLBEXS = 1.0/(XCXS-XBS)
-XLBS   = ( XAS*XCCS*MOMG(XALPHAS,XNUS,XBS) )**(-XLBEXS)
+!!XLBEXS = 1.0/(XCXS-XBS) ! Wurtz unused
+XLBEXS = 0. !1.0/(XCXS-XBS) ! Wurtz unused but used for XLB
+!XLBS   = ( XAS*XCCS*MOMG(XALPHAS,XNUS,XBS) )**(-XLBEXS)
+XLBS   = 1.0/(XAS*MOMG(XALPHAS,XNUS,XBS)) ! Wurtz
 !
 XLBEXG = 1.0/(XCXG-XBG)
 XLBG   = ( XAG*XCCG*MOMG(XALPHAG,XNUG,XBG))**(-XLBEXG)
@@ -266,7 +292,8 @@ IF (GFLAG) THEN
   WRITE(UNIT=ILUOUT0,FMT='(" XLBEXH =",E13.6," XLBH =",E13.6)') XLBEXH,XLBH
 END IF
 !
-XLBDAS_MAX = 500000
+XLBDAS_MAX = 500000. ! Attention il faut veuiller que LBDAS_MAX soit bien comparer avec un LBDAS avec une forme de Marshall-Palmer
+XLBDAS_MIN = 1000.
 XLBDAG_MAX = 100000.0
 !
 ZCONC_MAX  = 1.E6 ! Maximal concentration for falling particules set to 1 per cc
@@ -312,9 +339,21 @@ XFSEDRI  = XC_I*GAMMA_X0D(XNUI+(XDI+XBI)/XALPHAI)/GAMMA_X0D(XNUI+XBI/XALPHAI)*
 XFSEDCI  = XC_I*GAMMA_X0D(XNUI+XDI/XALPHAI)/GAMMA_X0D(XNUI)*     &
             (ZRHO00)**XCEXVT
 !
-XEXSEDS = (XBS+XDS-XCXS)/(XBS-XCXS)
-XFSEDS  = XCS*XAS*XCCS*MOMG(XALPHAS,XNUS,XBS+XDS)*                         &
-            (XAS*XCCS*MOMG(XALPHAS,XNUS,XBS))**(-XEXSEDS)*(ZRHO00)**XCEXVT
+!OPER
+!XEXSEDS = (XBS+XDS-XCXS)/(XBS-XCXS)    ! Wurtz
+!XFSEDS  = XCS*XAS*XCCS*MOMG(XALPHAS,XNUS,XBS+XDS)*                         & ! Wurtz
+!            (XAS*XCCS*MOMG(XALPHAS,XNUS,XBS))**(-XEXSEDS)*(ZRHO00)**XCEXVT  ! Wurtz
+
+!HOUZE/HAIC
+!XEXSEDS = -XDS !(2*XBS+XDS) ! Modif Wurtz
+!XFSEDS  = XCS*MOMG(XALPHAS,XNUS,XBS+XDS)/(MOMG(XALPHAS,XNUS,XBS))    & ! Modif Wurtz
+!            *(ZRHO00)**XCEXVT
+
+!LH_EXTENDED
+XEXSEDS = -XDS-XBS !-XNUS GAMMAGEN !(2*XBS+XDS) ! Modif Wurtz 
+XFSEDS  = XCS*MOMG(XALPHAS,XNUS,XBS+XDS)/(MOMG(XALPHAS,XNUS,XBS))    & ! Modif Wurtz
+            *(ZRHO00)**XCEXVT
+!
 !
 XEXSEDG = (XBG+XDG-XCXG)/(XBG-XCXG)
 XFSEDG  = XCG*XAG*XCCG*MOMG(XALPHAG,XNUG,XBG+XDG)*                         &
@@ -561,7 +600,7 @@ XR1DEPIS = XC1DEPIS *(XAI*XDICNVS_LIM**XBI)
 !
 ! Harrington parameterization for snow to ice conversion
 !
-XLBDASCNVI_MAX = 6000. ! lbdas max after Field (1999)
+XLBDASCNVI_MAX = 6000.*XTRANS_MP_GAMMAS ! lbdas max after Field (1999) ! Wurtz
 !
 XDSCNVI_LIM    = 125.E-6  ! size in microns
 XLBDASCNVI_LIM = (50.0**(1.0/(XALPHAS)))/XDSCNVI_LIM  ! ZLBDAS Limitation
@@ -574,10 +613,15 @@ XR1DEPSI = XC1DEPSI *(XAS*XDSCNVI_LIM**XBS)
 !
 ! Vapor deposition on snow and graupel and hail
 !
-X0DEPS = (4.0*XPI)*XCCS*XC1S*XF0S*MOMG(XALPHAS,XNUS,1.)
-X1DEPS = (4.0*XPI)*XCCS*XC1S*XF1S*SQRT(XCS)*MOMG(XALPHAS,XNUS,0.5*XDS+1.5)
-XEX0DEPS = XCXS-1.0
-XEX1DEPS = XCXS-0.5*(XDS+3.0)
+!X0DEPS = (4.0*XPI)*XCCS*XC1S*XF0S*MOMG(XALPHAS,XNUS,1.)
+!X1DEPS = (4.0*XPI)*XCCS*XC1S*XF1S*SQRT(XCS)*MOMG(XALPHAS,XNUS,0.5*XDS+1.5)
+X0DEPS = XLBS*(4.0*XPI)*XC1S*XF0S*MOMG(XALPHAS,XNUS,1.) ! Wurtz 
+X1DEPS = XLBS*(4.0*XPI)*XC1S*XF1S*SQRT(XCS)*MOMG(XALPHAS,XNUS,0.5*XDS+1.5) ! Wurtz 
+!XEX0DEPS = XCXS-1.0
+XEX0DEPS = XBS-1.0 ! Wurtz 
+!XEX1DEPS = XCXS-0.5*(XDS+3.0)
+XEX1DEPS = -0.5*(XDS+3.0) !GAMMAGEN !-0.5*(XDS+3.0)-XNUS ! Wurtz 
+
 !
 X0DEPG = (4.0*XPI)*XCCG*XC1G*XF0G*MOMG(XALPHAG,XNUG,1.)
 X1DEPG = (4.0*XPI)*XCCG*XC1G*XF1G*SQRT(XCG)*MOMG(XALPHAG,XNUG,0.5*XDG+1.5)
@@ -669,12 +713,16 @@ END IF
 !
 XDCSLIM  = 0.007 ! D_cs^lim = 7 mm as suggested by Farley et al. (1989)
 XCOLCS   = 1.0
-XEXCRIMSS= XCXS-XDS-2.0
-XCRIMSS  = (XPI/4.0)*XCOLCS*XCCS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0)
-XEXCRIMSG= XEXCRIMSS
-XCRIMSG  = XCRIMSS
-XSRIMCG  = XCCS*XAS*MOMG(XALPHAS,XNUS,XBS)
-XEXSRIMCG= XCXS-XBS
+!XEXCRIMSS= XBS-XDS-2.0 ! Wurtz 
+XEXCRIMSS = -XDS-2.0 ! Wurtz GAMMAGEN
+!XCRIMSS  = (XPI/4.0)*XCOLCS*XCCS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0)
+XCRIMSS  = XLBS * (XPI/4.0)*XCOLCS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0) ! Wurtz
+
+XEXCRIMSG= XEXCRIMSS ! Done RHO
+XCRIMSG  = XCRIMSS   ! Done RHO
+!XSRIMCG  = XCCS*XAS*MOMG(XALPHAS,XNUS,XBS) 
+XSRIMCG  = XLBS*XAS*MOMG(XALPHAS,XNUS,XBS) ! Wurtz
+XEXSRIMCG=-XBS ! Wurtz GAMMAGEN
 !
 GFLAG = .TRUE.
 IF (GFLAG) THEN
@@ -684,9 +732,9 @@ IF (GFLAG) THEN
 END IF
 !
 NGAMINC = 80
-XGAMINC_BOUND_MIN = 1.0E-1 ! Minimal value of (Lbda * D_cs^lim)**alpha
-XGAMINC_BOUND_MAX = 1.0E7  ! Maximal value of (Lbda * D_cs^lim)**alpha
-ZRATE = EXP(LOG(XGAMINC_BOUND_MAX/XGAMINC_BOUND_MIN)/REAL(NGAMINC-1))
+XGAMINC_BOUND_MIN = (1000.*XTRANS_MP_GAMMAS*XDCSLIM)**XALPHAS !1.0E-1 !Correction Wurtz ! Minimal value of (Lbda * D_cs^lim)**alpha
+XGAMINC_BOUND_MAX = (50000.*XTRANS_MP_GAMMAS*XDCSLIM)**XALPHAS !1.0E7 !Correction Wurtz ! Maximal value of (Lbda * D_cs^lim)**alpha
+ZRATE = EXP(LOG(XGAMINC_BOUND_MAX/XGAMINC_BOUND_MIN)/FLOAT(NGAMINC-1))
 !
 ALLOCATE( XGAMINC_RIM1(NGAMINC) )
 ALLOCATE( XGAMINC_RIM2(NGAMINC) )
@@ -732,13 +780,15 @@ XHMLINTP2 = 1.0 + XHMLINTP1*LOG( 25.E-6/(XGAMINC_HMC_BOUND_MIN)**(1.0/XALPHAC) )
 !
 !*       7.2    Constants for the accretion of raindrops onto aggregates
 !
-XFRACCSS = ((XPI**2)/24.0)*XCCS*XRHOLW*(ZRHO00**XCEXVT)
+!XFRACCSS = ((XPI**2)/24.0)*XCCS*XCCR*XRHOLW*(ZRHO00**XCEXVT)
+XFRACCSS = XLBS*((XPI**2)/24.0)*XCCR*XRHOLW*(ZRHO00**XCEXVT) ! Wurtz 
 !
 XLBRACCS1   =    MOMG(XALPHAS,XNUS,2.)*MOMG(XALPHAR,XNUR,3.)
 XLBRACCS2   = 2.*MOMG(XALPHAS,XNUS,1.)*MOMG(XALPHAR,XNUR,4.)
 XLBRACCS3   =                          MOMG(XALPHAR,XNUR,5.)
 !
-XFSACCRG = (XPI/4.0)*XAS*XCCS*(ZRHO00**XCEXVT)
+!XFSACCRG = (XPI/4.0)*XAS*XCCS*XCCR*(ZRHO00**XCEXVT)
+XFSACCRG = XLBS*(XPI/4.0)*XAS*XCCR*(ZRHO00**XCEXVT) ! Wurtz 
 !
 XLBSACCR1   =    MOMG(XALPHAR,XNUR,2.)*MOMG(XALPHAS,XNUS,XBS)
 XLBSACCR2   = 2.*MOMG(XALPHAR,XNUR,1.)*MOMG(XALPHAS,XNUS,XBS+1.)
@@ -750,9 +800,9 @@ XLBSACCR3   =                          MOMG(XALPHAS,XNUS,XBS+2.)
 ! Notice: One magnitude of lambda discretized over 10 points for snow
 !
 NACCLBDAS = 40
-XACCLBDAS_MIN = 5.0E1 ! Minimal value of Lbda_s to tabulate XKER_RACCS
-XACCLBDAS_MAX = 5.0E5 ! Maximal value of Lbda_s to tabulate XKER_RACCS
-ZRATE = LOG(XACCLBDAS_MAX/XACCLBDAS_MIN)/REAL(NACCLBDAS-1)
+XACCLBDAS_MIN = 5.0E2*XTRANS_MP_GAMMAS !5.0E1*XTRANS_MP_GAMMAS ! Minimal value of Lbda_s to tabulate XKER_RACCS
+XACCLBDAS_MAX = 1.0E5*XTRANS_MP_GAMMAS !5.0E5*XTRANS_MP_GAMMAS ! Maximal value of Lbda_s to tabulate XKER_RACCS
+ZRATE = LOG(XACCLBDAS_MAX/XACCLBDAS_MIN)/FLOAT(NACCLBDAS-1)
 XACCINTP1S = 1.0 / ZRATE
 XACCINTP2S = 1.0 - LOG( XACCLBDAS_MIN ) / ZRATE
 NACCLBDAR = 40
@@ -785,15 +835,15 @@ IF( (KACCLBDAS/=NACCLBDAS) .OR. (KACCLBDAR/=NACCLBDAR) .OR. (KND/=IND) .OR. &
     (PACCLBDAS_MIN/=XACCLBDAS_MIN) .OR. (PACCLBDAR_MIN/=XACCLBDAR_MIN) .OR. &
     (PFDINFTY/=ZFDINFTY)                                               ) THEN
   CALL RRCOLSS ( IND, XALPHAS, XNUS, XALPHAR, XNUR,                          &
-                 ZESR, XBR, XCS, XDS, XCR, XDR,                              &
+                 ZESR, XBR, XCS, XDS,XFVELOS, XCR, XDR,                      & ! Wurtz
                  XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
                  ZFDINFTY, XKER_RACCSS, XAG, XBS, XAS                        )
   CALL RZCOLX  ( IND, XALPHAS, XNUS, XALPHAR, XNUR,                          &
-                 ZESR, XBR, XCS, XDS, XCR, XDR,                              &
+                 ZESR, XBR, XCS, XDS, XFVELOS, XCR, XDR, 0.,                 & ! Wurtz
                  XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
                  ZFDINFTY, XKER_RACCS                                        )
   CALL RSCOLRG ( IND, XALPHAS, XNUS, XALPHAR, XNUR,                          &
-                 ZESR, XBS, XCS, XDS, XCR, XDR,                              &
+                 ZESR, XBS, XCS, XDS, XFVELOS, XCR, XDR,                     & ! Wurtz
                  XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
                  ZFDINFTY, XKER_SACCRG,XAG, XBS, XAS                         )
   WRITE(UNIT=ILUOUT0,FMT='("*****************************************")')
@@ -936,7 +986,8 @@ XCOLSG   = 0.01 ! Collection efficiency of S+G
 XCOLEXSG = 0.1  ! Temperature factor of the S+G collection efficiency
 WRITE (ILUOUT0, FMT=*) ' NEW Constants for the aggregate collection by the graupeln'
 WRITE (ILUOUT0, FMT=*) ' XCOLSG, XCOLEXSG  = ',XCOLSG,XCOLEXSG
-XFSDRYG = (XPI/4.0)*XCOLSG*XCCG*XCCS*XAS*(ZRHO00**XCEXVT)
+!XFSDRYG = (XPI/4.0)*XCOLSG*XCCG*XCCS*XAS*(ZRHO00**XCEXVT)
+XFSDRYG = XLBS*(XPI/4.0)*XCOLSG*XCCG*XAS*(ZRHO00**XCEXVT) ! Wurtz 
 !
 XLBSDRYG1   =    MOMG(XALPHAG,XNUG,2.)*MOMG(XALPHAS,XNUS,XBS)
 XLBSDRYG2   = 2.*MOMG(XALPHAG,XNUG,1.)*MOMG(XALPHAS,XNUS,XBS+1.)
@@ -966,8 +1017,8 @@ ZRATE = LOG(XDRYLBDAR_MAX/XDRYLBDAR_MIN)/REAL(NDRYLBDAR-1)
 XDRYINTP1R = 1.0 / ZRATE
 XDRYINTP2R = 1.0 - LOG( XDRYLBDAR_MIN ) / ZRATE
 NDRYLBDAS = 80
-XDRYLBDAS_MIN = 2.5E1 ! Minimal value of Lbda_s to tabulate XKER_SDRYG
-XDRYLBDAS_MAX = 2.5E9 ! Maximal value of Lbda_s to tabulate XKER_SDRYG
+XDRYLBDAS_MIN = 5.0E2*XTRANS_MP_GAMMAS ! Minimal value of Lbda_s to tabulate XKER_SDRYG
+XDRYLBDAS_MAX = 1.0E5*XTRANS_MP_GAMMAS ! Maximal value of Lbda_s to tabulate XKER_SDRYG
 ZRATE = LOG(XDRYLBDAS_MAX/XDRYLBDAS_MIN)/REAL(NDRYLBDAS-1)
 XDRYINTP1S = 1.0 / ZRATE
 XDRYINTP2S = 1.0 - LOG( XDRYLBDAS_MIN ) / ZRATE
@@ -999,7 +1050,7 @@ IF( (KDRYLBDAG/=NDRYLBDAG) .OR. (KDRYLBDAS/=NDRYLBDAS) .OR. (KND/=IND) .OR. &
     (PDRYLBDAG_MIN/=XDRYLBDAG_MIN) .OR. (PDRYLBDAS_MIN/=XDRYLBDAS_MIN) .OR. &
     (PFDINFTY/=ZFDINFTY)                                               ) THEN
   CALL RZCOLX ( IND, XALPHAG, XNUG, XALPHAS, XNUS,                          &
-                ZEGS, XBS, XCG, XDG, XCS, XDS,                              &
+                ZEGS, XBS, XCG, XDG,0., XCS, XDS,XFVELOS,                   & !Wurtz
                 XDRYLBDAG_MAX, XDRYLBDAS_MAX, XDRYLBDAG_MIN, XDRYLBDAS_MIN, &
                 ZFDINFTY, XKER_SDRYG                                        )
   WRITE(UNIT=ILUOUT0,FMT='("*****************************************")')
@@ -1065,7 +1116,7 @@ IF( (KDRYLBDAG/=NDRYLBDAG) .OR. (KDRYLBDAR/=NDRYLBDAR) .OR. (KND/=IND) .OR. &
     (PDRYLBDAG_MIN/=XDRYLBDAG_MIN) .OR. (PDRYLBDAR_MIN/=XDRYLBDAR_MIN) .OR. &
     (PFDINFTY/=ZFDINFTY)                                               ) THEN
   CALL RZCOLX ( IND, XALPHAG, XNUG, XALPHAR, XNUR,                          &
-                ZEGR, XBR, XCG, XDG, XCR, XDR,                              &
+                ZEGR, XBR, XCG, XDG,0., XCR, XDR, 0.,                             &
                 XDRYLBDAG_MAX, XDRYLBDAR_MAX, XDRYLBDAG_MIN, XDRYLBDAR_MIN, &
                 ZFDINFTY, XKER_RDRYG                                        )
   WRITE(UNIT=ILUOUT0,FMT='("*****************************************")')
@@ -1126,7 +1177,8 @@ XFWETH = (XPI/4.0)*XCCH*XCH*(ZRHO00**XCEXVT)*MOMG(XALPHAH,XNUH,XDH+2.0)
 !
 !*       9.2.2  Constants for the aggregate collection by the hailstones
 !
-XFSWETH = (XPI/4.0)*XCCH*XCCS*XAS*(ZRHO00**XCEXVT)
+!XFSWETH = (XPI/4.0)*XCCH*XCCS*XAS*(ZRHO00**XCEXVT)
+XFSWETH = XLBS*(XPI/4.0)*XCCH*XAS*(ZRHO00**XCEXVT) ! Wurtz  
 !
 XLBSWETH1   =    MOMG(XALPHAH,XNUH,2.)*MOMG(XALPHAS,XNUS,XBS)
 XLBSWETH2   = 2.*MOMG(XALPHAH,XNUH,1.)*MOMG(XALPHAS,XNUS,XBS+1.)
@@ -1143,8 +1195,8 @@ XLBGWETH3   =                          MOMG(XALPHAG,XNUG,XBG+2.)
 ! Notice: One magnitude of lambda discretized over 10 points
 !
 NWETLBDAS = 80
-XWETLBDAS_MIN = 2.5E1 ! Minimal value of Lbda_s to tabulate XKER_SWETH
-XWETLBDAS_MAX = 2.5E9 ! Maximal value of Lbda_s to tabulate XKER_SWETH
+XWETLBDAS_MIN = 5.0E2*XTRANS_MP_GAMMAS ! Minimal value of Lbda_s to tabulate XKER_SWETH
+XWETLBDAS_MAX = 1.0E5*XTRANS_MP_GAMMAS ! Maximal value of Lbda_s to tabulate XKER_SWETH
 ZRATE = LOG(XWETLBDAS_MAX/XWETLBDAS_MIN)/REAL(NWETLBDAS-1)
 XWETINTP1S = 1.0 / ZRATE
 XWETINTP2S = 1.0 - LOG( XWETLBDAS_MIN ) / ZRATE
@@ -1182,7 +1234,7 @@ IF( (KWETLBDAH/=NWETLBDAH) .OR. (KWETLBDAS/=NWETLBDAS) .OR. (KND/=IND) .OR. &
     (PWETLBDAH_MIN/=XWETLBDAH_MIN) .OR. (PWETLBDAS_MIN/=XWETLBDAS_MIN) .OR. &
     (PFDINFTY/=ZFDINFTY)                                               ) THEN
   CALL RZCOLX ( IND, XALPHAH, XNUH, XALPHAS, XNUS,                          &
-                ZEHS, XBS, XCH, XDH, XCS, XDS,                              &
+                ZEHS, XBS, XCH, XDH,0., XCS, XDS, XFVELOS,                  &
                 XWETLBDAH_MAX, XWETLBDAS_MAX, XWETLBDAH_MIN, XWETLBDAS_MIN, &
                 ZFDINFTY, XKER_SWETH                                        )
   WRITE(UNIT=ILUOUT0,FMT='("*****************************************")')
@@ -1248,7 +1300,7 @@ IF( (KWETLBDAH/=NWETLBDAH) .OR. (KWETLBDAG/=NWETLBDAG) .OR. (KND/=IND) .OR. &
     (PWETLBDAH_MIN/=XWETLBDAH_MIN) .OR. (PWETLBDAG_MIN/=XWETLBDAG_MIN) .OR. &
     (PFDINFTY/=ZFDINFTY)                                               ) THEN
   CALL RZCOLX ( IND, XALPHAH, XNUH, XALPHAG, XNUG,                          &
-                ZEHG, XBG, XCH, XDH, XCG, XDG,                              &
+                ZEHG, XBG, XCH, XDH,0., XCG, XDG,   0.,                     &
                 XWETLBDAH_MAX, XWETLBDAG_MAX, XWETLBDAH_MIN, XWETLBDAG_MIN, &
                 ZFDINFTY, XKER_GWETH                                        )
   WRITE(UNIT=ILUOUT0,FMT='("*****************************************")')
diff --git a/src/MNH/init_aerosol_concentration.f90 b/src/MNH/init_aerosol_concentration.f90
index e86998c4b..90d16461f 100644
--- a/src/MNH/init_aerosol_concentration.f90
+++ b/src/MNH/init_aerosol_concentration.f90
@@ -79,7 +79,7 @@ INTEGER                                 :: IKB, IKE
 !
 !
 IF ( LWARM .AND. LACTI ) THEN
-  DO JSV = NSV_LIMA_CCN_FREE, NSV_LIMA_CCN_ACTI+NMOD_CCN-1         
+  DO JSV = NSV_LIMA_CCN_FREE, NSV_LIMA_CCN_ACTI+NMOD_CCN-1
     PSVT(:,:,:,JSV) = 0.0
   ENDDO
   IKB = 1+JPVEXT
diff --git a/src/MNH/lima_cold_slow_processes.f90 b/src/MNH/lima_cold_slow_processes.f90
index 9fcacdd5a..098d3a406 100644
--- a/src/MNH/lima_cold_slow_processes.f90
+++ b/src/MNH/lima_cold_slow_processes.f90
@@ -95,7 +95,7 @@ USE MODD_CST,             ONLY: XP00, XRD, XRV, XMV, XMD, XCPD, XCPV,        &
 USE MODD_NSV,             ONLY: NSV_LIMA_NI
 USE MODD_PARAMETERS,      ONLY: JPHEXT, JPVEXT
 USE MODD_PARAM_LIMA,      ONLY: LSNOW, XRTMIN, XCTMIN, XALPHAI, XALPHAS,     &
-                                XNUI
+                                XNUI, XNUS, XLBDAS_MIN,XFVELOS,XTRANS_MP_GAMMAS
 USE MODD_PARAM_LIMA_COLD, ONLY: XLBI, XLBEXI, XLBS, XLBEXS, XBI, XCXS, XCCS, &
                                 XLBDAS_MAX, XDSCNVI_LIM, XLBDASCNVI_MAX,     &
                                 XC0DEPSI, XC1DEPSI, XR0DEPSI, XR1DEPSI,      &
@@ -103,7 +103,7 @@ USE MODD_PARAM_LIMA_COLD, ONLY: XLBI, XLBEXI, XLBS, XLBEXS, XBI, XCXS, XCCS, &
                                 XDICNVS_LIM, XLBDAICNVS_LIM,                 &
                                 XC0DEPIS, XC1DEPIS, XR0DEPIS, XR1DEPIS,      &
                                 XCOLEXIS, XAGGS_CLARGE1, XAGGS_CLARGE2,      &
-                                XAGGS_RLARGE1, XAGGS_RLARGE2
+                                XAGGS_RLARGE1, XAGGS_RLARGE2, XBS
 
 use mode_budget,          only: Budget_store_init, Budget_store_end
 use mode_tools,           only: Countjv
@@ -317,9 +317,15 @@ IF( IMICRO >= 1 ) THEN
    END WHERE
    ZLBDAS(:)  = 1.E10
    WHERE (ZRST(:)>XRTMIN(5) )
-      ZLBDAS(:) = XLBS*( ZRHODREF(:)*ZRST(:) )**XLBEXS
+WHERE(ZZT(:)>263.15) 
+ZLBDAS(:) = MAX(MIN(XLBDAS_MAX, 10**(14.554-0.0423*ZZT(:))),XLBDAS_MIN)
+END WHERE
+WHERE(ZZT(:)<=263.15) 
+ZLBDAS(:) = MAX(MIN(XLBDAS_MAX, 10**(6.226-0.0106*ZZT(:))),XLBDAS_MIN)
+END WHERE
    END WHERE
-!
+ZLBDAS(:) = ZLBDAS(:) * XTRANS_MP_GAMMAS
+
    ZKA(:) = 2.38E-2 + 0.0071E-2 * ( ZZT(:) - XTT )          ! k_a
    ZDV(:) = 0.211E-4 * (ZZT(:)/XTT)**1.94 * (XP00/ZPRES(:)) ! D_v
 !
@@ -342,16 +348,21 @@ IF( IMICRO >= 1 ) THEN
           call Budget_store_init( tbudgets(NBUDGET_SV1 - 1 + nsv_lima_ni), 'CNVI', pcis(:, :, :) * prhodj(:, :, :) )
       end if
 
-      WHERE ( ZRST(:)>XRTMIN(5) )
-         ZLBDAS(:)  = MIN( XLBDAS_MAX,                                           &
-                           XLBS*( ZRHODREF(:)*MAX( ZRST(:),XRTMIN(5) ) )**XLBEXS )
-      END WHERE
+!      WHERE ( ZRST(:)>XRTMIN(5) )
+!         ZLBDAS(:)  = MIN( XLBDAS_MAX,                                           &
+!                           XLBS*( ZRHODREF(:)*MAX( ZRST(:),XRTMIN(5) ) )**XLBEXS )
+!      END WHERE
       ZZW(:) = 0.0
       WHERE ( ZLBDAS(:)<XLBDASCNVI_MAX .AND. (ZRST(:)>XRTMIN(5)) &
                                        .AND. (ZSSI(:)<0.0)       )
          ZZW(:) = (ZLBDAS(:)*XDSCNVI_LIM)**(XALPHAS)
-         ZZX(:) = ( -ZSSI(:)/ZAI(:) ) * (XCCS*ZLBDAS(:)**XCXS)/ZRHODREF(:) * (ZZW(:)**XNUI) &
-                                                               * EXP(-ZZW(:))
+       !  ZZX(:) = ( -ZSSI(:)/ZAI(:) ) * (XCCS*ZLBDAS(:)**XCXS) * (ZZW(:)**XNUI) &
+                                                              ! * EXP(-ZZW(:))
+         ZZX(:) = ( -ZSSI(:)/ZAI(:) ) * (XLBS*ZRST(:)*ZLBDAS(:)**XBS) * (ZZW(:)**XNUI)  &
+
+                                                              * EXP(-ZZW(:))
+!MODIF CNVI -> concentration de la neige en kg-1 et non en m-3
+
 !
          ZZW(:) = MIN( ( XR0DEPSI+XR1DEPSI*ZCJ(:) )*ZZX(:),ZRSS(:) )
          ZRIS(:) = ZRIS(:) + ZZW(:)
@@ -384,8 +395,13 @@ IF( IMICRO >= 1 ) THEN
 
       ZZW(:) = 0.0
       WHERE ( (ZRST(:)>XRTMIN(5)) .AND. (ZRSS(:)>ZRTMIN(5)) )
-         ZZW(:) = ( ZSSI(:)/(ZRHODREF(:)*ZAI(:)) ) *                               &
-                  ( X0DEPS*ZLBDAS(:)**XEX0DEPS + X1DEPS*ZCJ(:)*ZLBDAS(:)**XEX1DEPS )
+!Correction BVIE rhodref
+!         ZZW(:) = ( ZSSI(:)/(ZRHODREF(:)*ZAI(:)) ) *                               &
+     !    ZZW(:) = ( ZSSI(:)/(ZAI(:)) ) *                               &
+      !            ( X0DEPS*ZLBDAS(:)**XEX0DEPS + X1DEPS*ZCJ(:)*ZLBDAS(:)**XEX1DEPS )
+         ZZW(:) = ( ZRHODREF(:)*ZRST(:)*ZSSI(:)/(ZAI(:)) ) *                               & ! Wurtz
+                 ( X0DEPS*ZLBDAS(:)**XEX0DEPS + (X1DEPS*ZCJ(:)*(1+(XFVELOS/(2.*ZLBDAS))**XALPHAS)**(-XNUS+XEX1DEPS/XALPHAS)*(ZLBDAS(:))**(XEX1DEPS+XBS))) ! Wurtz
+
          ZZW(:) =    MIN( ZRVS(:),ZZW(:)      )*(0.5+SIGN(0.5,ZZW(:))) &
                    - MIN( ZRSS(:),ABS(ZZW(:)) )*(0.5-SIGN(0.5,ZZW(:)))
          ZRSS(:) = ZRSS(:) + ZZW(:)
@@ -458,7 +474,11 @@ IF( IMICRO >= 1 ) THEN
       WHERE ( (ZRIT(:)>XRTMIN(4)) .AND. (ZRST(:)>XRTMIN(5)) .AND. (ZRIS(:)>ZRTMIN(4)) &
                                                             .AND. (ZCIS(:)>ZCTMIN(4)) )
          ZZW1(:,3) = (ZLBDAI(:) / ZLBDAS(:))**3
-         ZZW1(:,1) = (ZCIT(:)*(XCCS*ZLBDAS(:)**XCXS)/ZRHODREF(:)*EXP( XCOLEXIS*(ZZT(:)-XTT) )) &
+        ! ZZW1(:,1) = (ZCIT(:)*(XCCS*ZLBDAS(:)**XCXS)*EXP( XCOLEXIS*(ZZT(:)-XTT) )) &
+         !                                           / (ZLBDAI(:)**3)
+        ! ZZW1(:,1) = (ZCIT(:)*(XLBS*ZRHODREF(:)*ZRST(:)*ZLBDAS(:)**XBS)*EXP(XCOLEXIS*(ZZT(:)-XTT) )) &
+         !                                           / (ZLBDAI(:)**3)
+         ZZW1(:,1) = (ZCIT(:)*(XLBS*ZRST(:)*ZLBDAS(:)**XBS)*EXP(XCOLEXIS*(ZZT(:)-XTT) )) &
                                                     / (ZLBDAI(:)**3)
          ZZW1(:,2) = MIN( ZZW1(:,1)*(XAGGS_CLARGE1+XAGGS_CLARGE2*ZZW1(:,3)),ZCIS(:) )
          ZCIS(:) = ZCIS(:) - ZZW1(:,2)
diff --git a/src/MNH/lima_conversion_melting_snow.f90 b/src/MNH/lima_conversion_melting_snow.f90
index ff5a69146..a078380b3 100644
--- a/src/MNH/lima_conversion_melting_snow.f90
+++ b/src/MNH/lima_conversion_melting_snow.f90
@@ -61,9 +61,9 @@ END MODULE MODI_LIMA_CONVERSION_MELTING_SNOW
 !              ------------
 !
 USE MODD_CST,              ONLY : XTT, XMV, XMD, XLVTT, XCPV, XCL, XESTT, XRV
-USE MODD_PARAM_LIMA,       ONLY : XRTMIN
+USE MODD_PARAM_LIMA,       ONLY : XRTMIN, XFVELOS, XNUS, XALPHAS
 USE MODD_PARAM_LIMA_MIXED, ONLY : XFSCVMG
-USE MODD_PARAM_LIMA_COLD,  ONLY : X0DEPS, XEX0DEPS, X1DEPS, XEX1DEPS
+USE MODD_PARAM_LIMA_COLD,  ONLY : X0DEPS, XEX0DEPS, X1DEPS, XEX1DEPS, XBS
 !
 IMPLICIT NONE
 !
@@ -106,9 +106,13 @@ WHERE( (PRST(:)>XRTMIN(5)) .AND. (PT(:)>XTT) .AND. LDCOMPUTE(:) )
 !
 ! compute RSMLT
 !
-   ZW(:)  = XFSCVMG*MAX( 0.0,( -ZW(:) *             &
-                          ( X0DEPS*           PLBDS(:)**XEX0DEPS +     &
-                            X1DEPS*PCJ(:)*PLBDS(:)**XEX1DEPS ) ))!-    &
+!   ZW(:)  = XFSCVMG*MAX( 0.0,( -ZW(:) *             &
+!                          ( X0DEPS*           PLBDS(:)**XEX0DEPS +     &
+!                            X1DEPS*PCJ(:)*PLBDS(:)**XEX1DEPS ) ))!-    &
+
+   ZW(:)  = XFSCVMG*MAX( 0.0,( -ZW(:) *             &   ! Wurtz
+                          PRHODREF(:)*PRST(:)*( X0DEPS*            PLBDS(:)**XEX0DEPS +     &
+                            X1DEPS*PCJ(:)*(1+(XFVELOS/(2.*PLBDS(:)))**XALPHAS)**(-XNUS+XEX1DEPS/XALPHAS)*(PLBDS(:))**(XEX1DEPS+XBS))))
 ! On ne tient pas compte de la collection de pluie et gouttelettes par la neige si T>0 !!!! 
 ! Note that no heat is exchanged because the graupeln produced are still icy!!!
    P_RS_CMEL(:) = - ZW(:)
diff --git a/src/MNH/lima_droplets_riming_snow.f90 b/src/MNH/lima_droplets_riming_snow.f90
index 6bef29df3..d8c3a4bb6 100644
--- a/src/MNH/lima_droplets_riming_snow.f90
+++ b/src/MNH/lima_droplets_riming_snow.f90
@@ -73,11 +73,11 @@ END MODULE MODI_LIMA_DROPLETS_RIMING_SNOW
 !              ------------
 !
 USE MODD_CST,              ONLY : XTT
-USE MODD_PARAM_LIMA,       ONLY : XRTMIN, XCEXVT
+USE MODD_PARAM_LIMA,       ONLY : XRTMIN, XCEXVT, XFVELOS, XNUS,XALPHAS
 USE MODD_PARAM_LIMA_MIXED, ONLY : NGAMINC, XRIMINTP1, XRIMINTP2, XGAMINC_RIM1, XGAMINC_RIM2, &
-                                  XCRIMSS, XEXCRIMSS, XSRIMCG, XEXSRIMCG, &
+                                  XCRIMSS, XEXCRIMSS, XSRIMCG, XEXSRIMCG, XSRIMCG2, XSRIMCG3, XEXSRIMCG2, &
                                   XHMLINTP1, XHMLINTP2, XGAMINC_HMC, XHM_FACTS, XHMTMIN, XHMTMAX
-USE MODD_PARAM_LIMA_COLD,  ONLY : XMNU0
+USE MODD_PARAM_LIMA_COLD,  ONLY : XMNU0, XBS
 !
 IMPLICIT NONE
 !
@@ -171,7 +171,10 @@ WHERE( GRIM )
 !        4.     riming
 !
    ! Cloud droplets collected
-   P_RC_RIM(:) = - XCRIMSS * PRCT(:) * PLBDS(:)**XEXCRIMSS * PRHODREF(:)**(-XCEXVT)
+!   P_RC_RIM(:) = - XCRIMSS * PRCT(:) * PLBDS(:)**XEXCRIMSS * PRHODREF(:)**(-XCEXVT)
+   P_RC_RIM(:) = - XCRIMSS  * PRCT(:) * PRST(:)*(1+(XFVELOS/PLBDS(:))**XALPHAS)**(-XNUS+XEXCRIMSS/XALPHAS) &
+                                      * PRHODREF(:)**(-XCEXVT+1) &
+                                      * (PLBDS(:)) ** (XEXCRIMSS+XBS) ! Wurtz
    P_CC_RIM(:) = P_RC_RIM(:) *(PCCT(:)/PRCT(:)) ! Lambda_c**3
    !
    ! Cloud droplets collected on small aggregates add to snow
@@ -181,7 +184,8 @@ WHERE( GRIM )
    P_RG_RIM(:) = - P_RC_RIM(:) - P_RS_RIM(:) 
    !
    ! Large aggregates collecting droplets add to graupel (instant process ???)
-   ZZW3(:) = XSRIMCG * PLBDS(:)**XEXSRIMCG * (1.0 - ZZW2(:))/(PTSTEP*PRHODREF(:))
+!   ZZW3(:) = XSRIMCG * PLBDS(:)**XEXSRIMCG * (1.0 - ZZW2(:))/(PTSTEP*PRHODREF(:))
+   ZZW3(:) = XSRIMCG *PRHODREF(:)*PRST(:)* PLBDS(:)**(XEXSRIMCG+XBS) * (1.0 - ZZW2(:))!/(PTSTEP*PRHODREF(:)) ! Wurtz
    P_RS_RIM(:) = P_RS_RIM(:) - ZZW3(:) 
    P_RG_RIM(:) = P_RG_RIM(:) + ZZW3(:) 
    !
diff --git a/src/MNH/lima_graupel.f90 b/src/MNH/lima_graupel.f90
index ad114da36..6806322cd 100644
--- a/src/MNH/lima_graupel.f90
+++ b/src/MNH/lima_graupel.f90
@@ -330,9 +330,16 @@ WHERE( GDRY )
        			                                    * (ZVEC1(:) - 1.0)
    ZZW(:) = ZVEC3(:)
 !
+!   ZZW3(:) = XFSDRYG * ZZW(:) * EXP( XCOLEXSG*(PT(:)-XTT) )  & ! RSDRYG - rs collected by graupel in dry mode
+!                    *( PLBDS(:)**(XCXS-XBS) )*( PLBDG(:)**XCXG )    &
+!                    *( PRHODREF(:)**(-XCEXVT-1.) )                      &
+!                    *( XLBSDRYG1/( PLBDG(:)**2                 ) + &
+!                       XLBSDRYG2/( PLBDG(:)   * PLBDS(:)   ) + &
+!                       XLBSDRYG3/(                  PLBDS(:)**2) )
+! Wurtz 
    ZZW3(:) = XFSDRYG * ZZW(:) * EXP( XCOLEXSG*(PT(:)-XTT) )  & ! RSDRYG - rs collected by graupel in dry mode
-                    *( PLBDS(:)**(XCXS-XBS) )*( PLBDG(:)**XCXG )    &
-                    *( PRHODREF(:)**(-XCEXVT-1.) )                      &
+                    *( PRST(:))*( PLBDG(:)**XCXG )    &
+                    *( PRHODREF(:)**(-XCEXVT) )                      &
                     *( XLBSDRYG1/( PLBDG(:)**2                 ) + &
                        XLBSDRYG2/( PLBDG(:)   * PLBDS(:)   ) + &
                        XLBSDRYG3/(                  PLBDS(:)**2) )
diff --git a/src/MNH/lima_ice_aggregation_snow.f90 b/src/MNH/lima_ice_aggregation_snow.f90
index 15e01ec84..6c7be6eee 100644
--- a/src/MNH/lima_ice_aggregation_snow.f90
+++ b/src/MNH/lima_ice_aggregation_snow.f90
@@ -61,7 +61,7 @@ END MODULE MODI_LIMA_ICE_AGGREGATION_SNOW
 USE MODD_CST,             ONLY : XTT
 USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN
 USE MODD_PARAM_LIMA_COLD, ONLY : XBI, XCCS, XCXS, XCOLEXIS, XAGGS_CLARGE1, XAGGS_CLARGE2, &
-                                 XAGGS_RLARGE1, XAGGS_RLARGE2
+                                 XAGGS_RLARGE1, XAGGS_RLARGE2, XBS, XLBS
 !
 IMPLICIT NONE
 !
@@ -101,8 +101,10 @@ P_CI_AGGS(:) = 0.
 !
 WHERE ( (PRIT(:)>XRTMIN(4)) .AND. (PRST(:)>XRTMIN(5)) .AND. LDCOMPUTE(:) )
    ZZW1(:) = (PLBDI(:) / PLBDS(:))**3
-   ZZW2(:) = (PCIT(:)*(XCCS*PLBDS(:)**XCXS)/PRHODREF(:)*EXP( XCOLEXIS*(PT(:)-XTT) )) &
-        / (PLBDI(:)**3)
+!   ZZW2(:) = (PCIT(:)*(XCCS*PLBDS(:)**XCXS)/PRHODREF(:)*EXP( XCOLEXIS*(PT(:)-XTT) )) &
+!        / (PLBDI(:)**3)
+   ZZW2(:) = (PCIT(:)*(XLBS*PRST(:)*PLBDS(:)**XBS)*EXP(XCOLEXIS*(PT(:)-XTT) )) & ! Wurtz
+        / (PLBDI(:)**3)                                                          ! Wurtz
    ZZW3(:) = ZZW2(:)*(XAGGS_CLARGE1+XAGGS_CLARGE2*ZZW1(:))
 !
    P_CI_AGGS(:) = - ZZW3(:)
diff --git a/src/MNH/lima_ice_snow_deposition.f90 b/src/MNH/lima_ice_snow_deposition.f90
index 4d92b528a..8da16334f 100644
--- a/src/MNH/lima_ice_snow_deposition.f90
+++ b/src/MNH/lima_ice_snow_deposition.f90
@@ -84,8 +84,8 @@ SUBROUTINE LIMA_ICE_SNOW_DEPOSITION (PTSTEP, LDCOMPUTE,                        &
 !*       0.    DECLARATIONS
 !              ------------
 !
-USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN, XALPHAI, XALPHAS, XNUI, XNUS 
-USE MODD_PARAM_LIMA_COLD, ONLY : XCXS, XCCS, &
+USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN, XALPHAI, XALPHAS, XNUI, XNUS, XFVELOS
+USE MODD_PARAM_LIMA_COLD, ONLY : XCXS, XCCS, XLBS, XBS, &
                                  XLBDAS_MAX, XDSCNVI_LIM, XLBDASCNVI_MAX,     &
                                  XC0DEPSI, XC1DEPSI, XR0DEPSI, XR1DEPSI,      &
                                  XSCFAC, X1DEPS, X0DEPS, XEX1DEPS, XEX0DEPS,  &
@@ -166,7 +166,8 @@ WHERE( GMICRO )
    WHERE ( PLBDS(:)<XLBDASCNVI_MAX .AND. (PRST(:)>XRTMIN(5)) &
                                    .AND. (PSSI(:)<0.0)       )
       ZZW(:) = (PLBDS(:)*XDSCNVI_LIM)**(XALPHAS)
-      ZZX(:) = ( -PSSI(:)/PAI(:) ) * (XCCS*PLBDS(:)**XCXS)/PRHODREF(:) * (ZZW(:)**XNUS) * EXP(-ZZW(:))
+!      ZZX(:) = ( -PSSI(:)/PAI(:) ) * (XCCS*PLBDS(:)**XCXS)/PRHODREF(:) * (ZZW(:)**XNUS) * EXP(-ZZW(:))
+      ZZX(:) = ( -PSSI(:)/PAI(:) ) * (XLBS*PRST(:)*PLBDS(:)**XBS) * (ZZW(:)**XNUS) * EXP(-ZZW(:))
 !
       ZZW(:) = ( XR0DEPSI+XR1DEPSI*PCJ(:) )*ZZX(:)
 !
@@ -187,8 +188,10 @@ WHERE( GMICRO )
 !
    ZZW(:) = 0.0
    WHERE ( (PRST(:)>XRTMIN(5)) )
-      ZZW(:) = ( PSSI(:)/(PAI(:))/PRHODREF(:) ) * &
-           ( X0DEPS*PLBDS(:)**XEX0DEPS + X1DEPS*PCJ(:)*PLBDS(:)**XEX1DEPS )
+!      ZZW(:) = ( PSSI(:)/(PAI(:))/PRHODREF(:) ) * &
+!           ( X0DEPS*PLBDS(:)**XEX0DEPS + X1DEPS*PCJ(:)*PLBDS(:)**XEX1DEPS )
+      ZZW(:) =10* ( PRHODREF(:) * PRST(:)*PSSI(:)/PAI(:)) *                               &       !Modif Wurtz
+                 ( X0DEPS*PLBDS(:)**XEX0DEPS + (X1DEPS*PCJ(:)*(1+(XFVELOS/(2.*PLBDS(:)))**XALPHAS)**(-XNUS+XEX1DEPS/XALPHAS)*(PLBDS(:))**(XBS+XEX1DEPS))) ! Wurtz
       ZZW(:) =    ZZW(:)*(0.5+SIGN(0.5,ZZW(:))) - ABS(ZZW(:))*(0.5-SIGN(0.5,ZZW(:)))
    END WHERE
 !
diff --git a/src/MNH/lima_mixed.f90 b/src/MNH/lima_mixed.f90
index 49024b7b5..7ddedc901 100644
--- a/src/MNH/lima_mixed.f90
+++ b/src/MNH/lima_mixed.f90
@@ -108,9 +108,9 @@ USE MODD_CST,              ONLY: XP00, XRD, XRV, XMV, XMD, XCPD, XCPV,       &
 USE MODD_NSV
 USE MODD_PARAMETERS,       ONLY: JPHEXT, JPVEXT
 USE MODD_PARAM_LIMA,       ONLY: NMOD_IFN, XRTMIN, XCTMIN, LWARM, LCOLD,     &
-                                 NMOD_CCN, NMOD_IMM, LRAIN, LSNOW, LHAIL
+                                 NMOD_CCN, NMOD_IMM, LRAIN, LSNOW, LHAIL, XLBDAS_MIN, XTRANS_MP_GAMMAS
 USE MODD_PARAM_LIMA_WARM,  ONLY: XLBC, XLBEXC, XLBR, XLBEXR
-USE MODD_PARAM_LIMA_COLD,  ONLY: XLBI, XLBEXI, XLBS, XLBEXS, XSCFAC
+USE MODD_PARAM_LIMA_COLD,  ONLY: XLBI, XLBEXI, XLBS, XLBEXS, XSCFAC, XLBDAS_MAX
 USE MODD_PARAM_LIMA_MIXED, ONLY: XLBG, XLBEXG, XLBH, XLBEXH
 
 use mode_tools,            only: Countjv
@@ -467,9 +467,20 @@ IF( IMICRO >= 1 ) THEN
       ZLBDAI(:) = ( XLBI*ZCIT(:) / ZRIT(:) )**XLBEXI
    END WHERE
    ZLBDAS(:)  = 1.E10
-   WHERE (ZRST(:)>XRTMIN(5) )
-      ZLBDAS(:) = XLBS*( ZRHODREF(:)*ZRST(:) )**XLBEXS
-   END WHERE
+ !  WHERE (ZRST(:)>XRTMIN(5) )
+  !    ZLBDAS(:) = XLBS*( ZRHODREF(:)*ZRST(:) )**XLBEXS
+   !END WHERE
+!Wurtz 
+      WHERE (ZRST(:)>XRTMIN(5) )
+                        WHERE(ZZT(:)>263.15) 
+                                ZLBDAS(:) = MAX(MIN(XLBDAS_MAX, 10**(14.554-0.0423*ZZT(:))),XLBDAS_MIN)
+                        END WHERE
+
+                       WHERE(ZZT(:)<=263.15) 
+                                ZLBDAS(:) = MAX(MIN(XLBDAS_MAX, 10**(6.226-0.0106*ZZT(:))),XLBDAS_MIN)
+                       END WHERE
+      END WHERE
+   ZLBDAS(:) = ZLBDAS(:)*XTRANS_MP_GAMMAS
    ZLBDAG(:)  = 1.E10
    WHERE (ZRGT(:)>XRTMIN(6) )
       ZLBDAG(:) = XLBG*( ZRHODREF(:)*ZRGT(:) )**XLBEXG
diff --git a/src/MNH/lima_mixed_fast_processes.f90 b/src/MNH/lima_mixed_fast_processes.f90
index 09c86c8a2..e12dfcad5 100644
--- a/src/MNH/lima_mixed_fast_processes.f90
+++ b/src/MNH/lima_mixed_fast_processes.f90
@@ -296,10 +296,17 @@ IF( IGRIM>0 ) THEN
 !        1.1.4  riming of the small sized aggregates
 !
    WHERE ( GRIM(:) )
-      ZZW1(:,1) = MIN( ZRCS(:),                              &
-  	           XCRIMSS * ZZW(:) * ZRCT(:)                & ! RCRIMSS
-  	                            *   ZLBDAS(:)**XEXCRIMSS &
-    			            * ZRHODREF(:)**(-XCEXVT) )
+    !  ZZW1(:,1) = MIN( ZRCS(:),                              &
+      !           XCRIMSS * ZZW(:) * ZRCT(:)                & ! RCRIMSS
+         !                           *   ZLBDAS(:)**XEXCRIMSS &
+    !		            * ZRHODREF(:)**(-XCEXVT) )
+!
+      ZZW1(:,1) = MIN( ZRCS(:),                              & ! Wurtz
+             XCRIMSS * ZZW(:) * ZRCT(:)                & ! RCRIMSS
+                                      *  ZRST(:)*(1+(XFVELOS/ZLBDAS(:))**XALPHAS)**(-XNUS+XEXCRIMSS/XALPHAS) &
+                                      * ZRHODREF(:)**(-XCEXVT+1) &
+                                      * (ZLBDAS(:)) ** (XEXCRIMSS+XBS) )! Wurtz
+
       ZRCS(:) = ZRCS(:) - ZZW1(:,1)
       ZRSS(:) = ZRSS(:) + ZZW1(:,1)
       ZTHS(:) = ZTHS(:) + ZZW1(:,1)*(ZLSFACT(:)-ZLVFACT(:)) ! f(L_f*(RCRIMSS))
@@ -318,14 +325,24 @@ IF( IGRIM>0 ) THEN
 !
 !
    WHERE ( GRIM(:) .AND. (ZRSS(:)>XRTMIN(5)/PTSTEP) )
-      ZZW1(:,2) = MIN( ZRCS(:),                     &
-    	           XCRIMSG * ZRCT(:)                & ! RCRIMSG
-    	                   *  ZLBDAS(:)**XEXCRIMSG  &
-  	                   * ZRHODREF(:)**(-XCEXVT) &
-    		           - ZZW1(:,1)              )
+     ! ZZW1(:,2) = MIN( ZRCS(:),                     &
+     !           XCRIMSG * ZRCT(:)                & ! RCRIMSG
+     !                  *  ZLBDAS(:)**XEXCRIMSG  &
+  	  !                 * ZRHODREF(:)**(-XCEXVT) &
+    !            - ZZW1(:,1)              )
+
+      ZZW1(:,2) = MIN( ZRCS(:),                     & ! Modif Wurtz
+                   XCRIMSG * ZRCT(:)*  ZRST(:)                & ! RCRIMSG
+                           *(1+(XFVELOS/ZLBDAS(:))**XALPHAS)**(-XNUS+XEXCRIMSG/XALPHAS)*ZLBDAS(:)**(XBS+XEXCRIMSG)  &
+                           * ZRHODREF(:)**(-XCEXVT+1) &
+                           - ZZW1(:,1)              )
+
       ZZW1(:,3) = MIN( ZRSS(:),                         &
-                       XSRIMCG * ZLBDAS(:)**XEXSRIMCG   & ! RSRIMCG
-   	                       * (1.0 - ZZW(:) )/(PTSTEP*ZRHODREF(:)))
+                     !  XSRIMCG * ZLBDAS(:)**XEXSRIMCG   & ! RSRIMCG
+                 !         * (1.0 - ZZW(:) )/(PTSTEP*ZRHODREF(:)))
+                       XSRIMCG * ZRHODREF(:) *ZRST(:)* ZLBDAS(:)**(XEXSRIMCG+XBS)   & ! RSRIMCG ! Wurtz
+                               * (1.0 - ZZW(:) )/(PTSTEP*ZRHODREF(:)))
+
       ZRCS(:) = ZRCS(:) - ZZW1(:,2)
       ZRSS(:) = ZRSS(:) - ZZW1(:,3)
       ZRGS(:) = ZRGS(:) + ZZW1(:,2) + ZZW1(:,3)
@@ -476,8 +493,13 @@ IF( IGACC>0 .AND. LRAIN) THEN
 !        1.3.4  raindrop accretion on the small sized aggregates
 !
    WHERE ( GACC(:) )
-      ZZW1(:,2) = ZCRT(:) *                                           & !! coef of RRACCS
-              XFRACCSS*( ZLBDAS(:)**XCXS )*( ZRHODREF(:)**(-XCEXVT-1.) ) &
+ !     ZZW1(:,2) = ZCRT(:) *                                           & !! coef of RRACCS
+  !            XFRACCSS*( ZLBDAS(:)**XCXS )*( ZRHODREF(:)**(-XCEXVT-1.) ) &
+   !      *( XLBRACCS1/((ZLBDAS(:)**2)               ) +                  &
+    !        XLBRACCS2/( ZLBDAS(:)    * ZLBDAR(:)    ) +                  &
+     !       XLBRACCS3/(               (ZLBDAR(:)**2)) )/ZLBDAR(:)**3
+      ZZW1(:,2) = ZCRT(:) *                                           & !! coef of RRACCS ! Wurtz
+              XFRACCSS*( ZRST(:)*ZLBDAS(:)**XBS )*( ZRHODREF(:)**(-XCEXVT) ) &
          *( XLBRACCS1/((ZLBDAS(:)**2)               ) +                  &
             XLBRACCS2/( ZLBDAS(:)    * ZLBDAR(:)    ) +                  &
             XLBRACCS3/(               (ZLBDAR(:)**2)) )/ZLBDAR(:)**3
@@ -520,8 +542,14 @@ IF( IGACC>0 .AND. LRAIN) THEN
 !
    WHERE ( GACC(:) .AND. (ZRSS(:)>XRTMIN(5)/PTSTEP) )
       ZZW1(:,2) = MAX( MIN( ZRRS(:),ZZW1(:,2)-ZZW1(:,4) ) , 0. )      ! RRACCSG
-      ZZW1(:,3) = MIN( ZRSS(:),XFSACCRG*ZZW(:)*                     & ! RSACCRG
-            ( ZLBDAS(:)**(XCXS-XBS) )*( ZRHODREF(:)**(-XCEXVT-1.) ) &
+   !   ZZW1(:,3) = MIN( ZRSS(:),XFSACCRG*ZZW(:)*                     & ! RSACCRG
+    !        ( ZLBDAS(:)**(XCXS-XBS) )*( ZRHODREF(:)**(-XCEXVT-1.) ) &
+     !      *( XLBSACCR1/((ZLBDAR(:)**2)               ) +           &
+      !        XLBSACCR2/( ZLBDAR(:)    * ZLBDAS(:)    ) +           &
+       !       XLBSACCR3/(               (ZLBDAS(:)**2)) ) )
+
+      ZZW1(:,3) = MIN( ZRSS(:),XFSACCRG*ZZW(:)*                     & ! RSACCRG ! Wurtz
+            ( ZRST(:) )*( ZRHODREF(:)**(-XCEXVT) ) &
            *( XLBSACCR1/((ZLBDAR(:)**2)               ) +           &
               XLBSACCR2/( ZLBDAR(:)    * ZLBDAS(:)    ) +           &
               XLBSACCR3/(               (ZLBDAS(:)**2)) ) )
@@ -572,9 +600,16 @@ WHERE( (ZRST(:)>XRTMIN(5)) .AND. (ZRSS(:)>XRTMIN(5)/PTSTEP) .AND. (ZZT(:)>XTT) )
 !
 ! compute RSMLT
 !
-   ZZW(:)  = MIN( ZRSS(:), XFSCVMG*MAX( 0.0,( -ZZW(:) *             &
-                          ( X0DEPS*       ZLBDAS(:)**XEX0DEPS +     &
-                            X1DEPS*ZCJ(:)*ZLBDAS(:)**XEX1DEPS ) -   &
+  ! 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 ) ) )
+
+   ZZW(:)  = MIN( ZRSS(:), XFSCVMG*MAX( 0.0,( -ZZW(:) *             & ! Wurtz
+                          ZRHODREF(:) * ZRST(:)*( X0DEPS*            ZLBDAS(:)**XEX0DEPS +     &
+                            X1DEPS*ZCJ(:)*(1+(XFVELOS/(2.*ZLBDAS(:)))**XALPHAS)**(-XNUS+XEX1DEPS/XALPHAS)*(ZLBDAS(:))**(XEX1DEPS+XBS))-    &
                                     ( ZZW1(:,1)+ZZW1(:,4) ) *       &
                              ( ZRHODREF(:)*XCL*(XTT-ZZT(:))) ) /    &
                                             ( ZRHODREF(:)*XLMTT ) ) )
@@ -738,14 +773,23 @@ IF( IGDRY>0 ) THEN
    END DO
    ZZW(:) = UNPACK( VECTOR=ZVEC3(:),MASK=GDRY,FIELD=0.0 )
 !
-   WHERE( GDRY(:) )
+   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) ) )
+
+
       ZZW1(:,3) = MIN( ZRSS(:),XFSDRYG*ZZW(:)                         & ! RSDRYG
                                       * EXP( XCOLEXSG*(ZZT(:)-XTT) )  &
-                    *( ZLBDAS(:)**(XCXS-XBS) )*( ZLBDAG(:)**XCXG )    &
-                    *( ZRHODREF(:)**(-XCEXVT-1.) )                    &
+                    *( ZRST(:)) )*( ZLBDAG(:)**XCXG )    &
+                    *( ZRHODREF(:)**(-XCEXVT) )                    &
                          *( XLBSDRYG1/( ZLBDAG(:)**2              ) + &
                             XLBSDRYG2/( ZLBDAG(:)   * ZLBDAS(:)   ) + &
-                            XLBSDRYG3/(               ZLBDAS(:)**2) ) )
+                            XLBSDRYG3/(               ZLBDAS(:)**2) ) 
    END WHERE
    DEALLOCATE(IVEC2)
    DEALLOCATE(IVEC1)
@@ -1184,9 +1228,17 @@ IF( IHAIL>0 ) THEN
       ZZW(:) = UNPACK( VECTOR=ZVEC3(:),MASK=GWET,FIELD=0.0 )
 !
       WHERE( GWET(:) )
+!        ZZW1(:,3) = MIN( ZRSS(:),XFSWETH*ZZW(:)                       & ! RSWETH
+ !                     *( ZLBDAS(:)**(XCXS-XBS) )*( ZLBDAH(:)**XCXH )  &
+  !     	                 *( ZRHODREF(:)**(-XCEXVT-1.) )               &
+   !                      *( XLBSWETH1/( ZLBDAH(:)**2              ) + &
+     !                       XLBSWETH2/( ZLBDAH(:)   * ZLBDAS(:)   ) + &
+    !                        XLBSWETH3/(               ZLBDAS(:)**2) ) )
+
+! Wurtz 
         ZZW1(:,3) = MIN( ZRSS(:),XFSWETH*ZZW(:)                       & ! RSWETH
-                      *( ZLBDAS(:)**(XCXS-XBS) )*( ZLBDAH(:)**XCXH )  &
-       	                 *( ZRHODREF(:)**(-XCEXVT-1.) )               &
+                      *( ZRST(:))*( ZLBDAH(:)**XCXH )  &
+                        *( ZRHODREF(:)**(-XCEXVT) )               &
                          *( XLBSWETH1/( ZLBDAH(:)**2              ) + &
                             XLBSWETH2/( ZLBDAH(:)   * ZLBDAS(:)   ) + &
                             XLBSWETH3/(               ZLBDAS(:)**2) ) )
diff --git a/src/MNH/lima_rain_accr_snow.f90 b/src/MNH/lima_rain_accr_snow.f90
index 01c31afbe..eed58c987 100644
--- a/src/MNH/lima_rain_accr_snow.f90
+++ b/src/MNH/lima_rain_accr_snow.f90
@@ -67,7 +67,7 @@ END MODULE MODI_LIMA_RAIN_ACCR_SNOW
 !              ------------
 !
 USE MODD_CST,              ONLY : XTT
-USE MODD_PARAM_LIMA,       ONLY : XRTMIN, XCEXVT
+USE MODD_PARAM_LIMA,       ONLY : XRTMIN, XCEXVT, XTRANS_MP_GAMMAS
 USE MODD_PARAM_LIMA_COLD,  ONLY : XBS, XCXS
 USE MODD_PARAM_LIMA_MIXED, ONLY : NACCLBDAS, XACCINTP1S, XACCINTP2S,         &
                                   NACCLBDAR, XACCINTP1R, XACCINTP2R,         &
@@ -141,7 +141,7 @@ WHERE( GACC )
 !
 !        1.3.1  select the (ZLBDAS,ZLBDAR) couplet
    !
-   ZVEC1(:) = MAX(MIN(PLBDS(:),5.E5),5.E1)
+   ZVEC1(:) = MAX(MIN(PLBDS(:),5.E5*XTRANS_MP_GAMMAS),5.E1*XTRANS_MP_GAMMAS)
    ZVEC2(:) = PLBDR(:)
 !
 !        1.3.2  find the next lower indice for the ZLBDAS and for the ZLBDAR
@@ -212,8 +212,13 @@ WHERE( GACC )
 !      
 ! BVIE manque PCRT ???????????????????????????????????
 !      ZZW4(:) =                                            & !! coef of RRACCS and RRACCS
-   ZZW4(:) = PCRT(:)                                     & !! coef of RRACCS and RRACCS
-         *  XFRACCSS *( PLBDS(:)**XCXS )*( PRHODREF(:)**(-XCEXVT-1.) ) &
+!   ZZW4(:) = PCRT(:)                                     & !! coef of RRACCS and RRACCS
+!         *  XFRACCSS *( PLBDS(:)**XCXS )*( PRHODREF(:)**(-XCEXVT-1.) ) &
+!         *( XLBRACCS1/( PLBDS(:)**2 )                +                     &
+!            XLBRACCS2/( PLBDS(:) * PLBDR(:)    ) +                     &
+!            XLBRACCS3/(                PLBDR(:)**2 ) ) / PLBDR(:)**3
+   ZZW4(:) = PCRT(:)                                     & !! coef of RRACCS and RRACCS ! Wurtz
+         *  XFRACCSS *( PRST(:)*PLBDS(:)**XBS )*( PRHODREF(:)**(-XCEXVT) ) &
          *( XLBRACCS1/( PLBDS(:)**2 )                +                     &
             XLBRACCS2/( PLBDS(:) * PLBDR(:)    ) +                     &
             XLBRACCS3/(                PLBDR(:)**2 ) ) / PLBDR(:)**3
@@ -227,8 +232,13 @@ WHERE( GACC )
 !        1.3.6  raindrop accretion-conversion of the large sized aggregates
 !               into graupeln
 !
-   ZZW5(:) = XFSACCRG*ZZW3(:) *                             & ! RSACCRG
-            ( PLBDS(:)**(XCXS-XBS) )*( PRHODREF(:)**(-XCEXVT-1.) ) &
+!   ZZW5(:) = XFSACCRG*ZZW3(:) *                             & ! RSACCRG
+!            ( PLBDS(:)**(XCXS-XBS) )*( PRHODREF(:)**(-XCEXVT-1.) ) &
+!           *( XLBSACCR1/((PLBDR(:)**2)               ) +           &
+!              XLBSACCR2/( PLBDR(:)    * PLBDS(:) ) +           &
+!              XLBSACCR3/(                  (PLBDS(:)**2)) )
+   ZZW5(:) = XFSACCRG*ZZW3(:) *                             & ! RSACCRG ! Wurtz
+            ( PRST(:) )*( PRHODREF(:)**(-XCEXVT) ) &
            *( XLBSACCR1/((PLBDR(:)**2)               ) +           &
               XLBSACCR2/( PLBDR(:)    * PLBDS(:) ) +           &
               XLBSACCR3/(                  (PLBDS(:)**2)) )
diff --git a/src/MNH/lima_sedimentation.f90 b/src/MNH/lima_sedimentation.f90
index 365ae0f23..d60655b96 100644
--- a/src/MNH/lima_sedimentation.f90
+++ b/src/MNH/lima_sedimentation.f90
@@ -183,11 +183,23 @@ DO JN = 1 ,  NSPLITSED(KID)
          IF (KMOMENTS==2) ZCS(JL) = PCS(I1(JL),I2(JL),I3(JL))
       END DO
 !
-      IF (KMOMENTS==1) ZLBDA(:) = XLB(KID) * ( ZRHODREF(:) * ZRS(:) )**XLBEX(KID)
-      IF (KMOMENTS==2) ZLBDA(:) = ( XLB(KID)*ZCS(:) / ZRS(:) )**XLBEX(KID)
+      IF (KID == 5) THEN
+         ZLBDA(:) = 1.E10
+         WHERE(ZT(:)>263.15 .AND. ZRS(:)>XRTMIN(5))
+            ZLBDA(:) = MAX(MIN(XLBDAS_MAX, 10**(14.554-0.0423*ZT(:))),XLBDAS_MIN)
+         END WHERE
+         WHERE(ZT(:)<=263.15 .AND. ZRS(:)>XRTMIN(5))
+            ZLBDA(:) = MAX(MIN(XLBDAS_MAX, 10**(6.226-0.0106*ZT(:))),XLBDAS_MIN)
+         END WHERE
+         ZLBDA(:) = ZLBDA(:)*XTRANS_MP_GAMMAS
+         ZZW(:) = XFSEDR(KID) * ZRHODREF(:)**(1.-XCEXVT)*ZRS(:)*(1 + (XFVELOS/ZLBDA(:))**XALPHAS)**(-XNUS+XEXSEDS/XALPHAS) * ZLBDA(:)**(XBS+XEXSEDS)
+      ELSE
+         IF (KMOMENTS==1) ZLBDA(:) = XLB(KID) * ( ZRHODREF(:) * ZRS(:) )**XLBEX(KID)
+         IF (KMOMENTS==2) ZLBDA(:) = ( XLB(KID)*ZCS(:) / ZRS(:) )**XLBEX(KID)
+         ZZY(:) = ZRHODREF(:)**(-XCEXVT) * ZLBDA(:)**(-XD(KID))
+         ZZW(:) = XFSEDR(KID) * ZRS(:) * ZZY(:) * ZRHODREF(:)
+      END IF ! Wurtz
 !
-      ZZY(:) = ZRHODREF(:)**(-XCEXVT) * ZLBDA(:)**(-XD(KID))
-      ZZW(:) = XFSEDR(KID) * ZRS(:) * ZZY(:) * ZRHODREF(:)
       IF (KMOMENTS==2) ZZX(:) = XFSEDC(KID) * ZCS(:) * ZZY(:) * ZRHODREF(:)
 
       IF (KID==2) THEN
diff --git a/src/MNH/lima_snow_deposition.f90 b/src/MNH/lima_snow_deposition.f90
index d39f714d8..e29b71173 100644
--- a/src/MNH/lima_snow_deposition.f90
+++ b/src/MNH/lima_snow_deposition.f90
@@ -68,8 +68,8 @@ SUBROUTINE LIMA_SNOW_DEPOSITION (LDCOMPUTE,                                &
 !*       0.    DECLARATIONS
 !              ------------
 !
-USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN, XALPHAI, XALPHAS, XNUI, XNUS 
-USE MODD_PARAM_LIMA_COLD, ONLY : XCXS, XCCS, &
+USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN, XALPHAI, XALPHAS, XNUI, XNUS, XFVELOS 
+USE MODD_PARAM_LIMA_COLD, ONLY : XCXS, XCCS,XLBS,XBS, &
                                  XLBDAS_MAX, XDSCNVI_LIM, XLBDASCNVI_MAX,     &
                                  XC0DEPSI, XC1DEPSI, XR0DEPSI, XR1DEPSI,      &
                                  XSCFAC, X1DEPS, X0DEPS, XEX1DEPS, XEX0DEPS,  &
@@ -132,7 +132,8 @@ WHERE( GMICRO )
    WHERE ( PLBDS(:)<XLBDASCNVI_MAX .AND. (PRST(:)>XRTMIN(5)) &
                                    .AND. (PSSI(:)<0.0)       )
       ZZW(:) = (PLBDS(:)*XDSCNVI_LIM)**(XALPHAS)
-      ZZX(:) = ( -PSSI(:)/PAI(:) ) * (XCCS*PLBDS(:)**XCXS) * (ZZW(:)**XNUS) * EXP(-ZZW(:))
+!      ZZX(:) = ( -PSSI(:)/PAI(:) ) * (XCCS*PLBDS(:)**XCXS) * (ZZW(:)**XNUS) * EXP(-ZZW(:))
+      ZZX(:) = ( -PSSI(:)/PAI(:) ) * (XLBS*PRST(:)*PLBDS(:)**XBS) * (ZZW(:)**XNUS) * EXP(-ZZW(:))
 !
       ZZW(:) = ( XR0DEPSI+XR1DEPSI*PCJ(:) )*ZZX(:)
 !
@@ -149,8 +150,10 @@ WHERE( GMICRO )
 !
    ZZW(:) = 0.0
    WHERE ( (PRST(:)>XRTMIN(5)) )
-      ZZW(:) = ( PSSI(:)/(PAI(:)*PRHODREF(:)) ) * &
-           ( X0DEPS*PLBDS(:)**XEX0DEPS + X1DEPS*PCJ(:)*PLBDS(:)**XEX1DEPS )
+!      ZZW(:) = ( PSSI(:)/(PAI(:)*PRHODREF(:)) ) * &
+!           ( X0DEPS*PLBDS(:)**XEX0DEPS + X1DEPS*PCJ(:)*PLBDS(:)**XEX1DEPS )
+      ZZW(:) = ( PRHODREF(:) * PRST(:)*PSSI(:)/(PAI(:)) ) * & ! Wurtz
+                ( X0DEPS*PLBDS(:)**XEX0DEPS + (X1DEPS*PCJ(:)*(1+(XFVELOS/(2.*PLBDS))**XALPHAS)**(-XNUS+XEX1DEPS/XALPHAS)*(PLBDS(:))**(XEX1DEPS+XBS)) ) ! Wurtz
       ZZW(:) =    ZZW(:)*(0.5+SIGN(0.5,ZZW(:))) - ABS(ZZW(:))*(0.5-SIGN(0.5,ZZW(:)))
    END WHERE
 !
diff --git a/src/MNH/lima_tendencies.f90 b/src/MNH/lima_tendencies.f90
index bd98d503c..0aac7dec6 100644
--- a/src/MNH/lima_tendencies.f90
+++ b/src/MNH/lima_tendencies.f90
@@ -239,7 +239,7 @@ SUBROUTINE LIMA_TENDENCIES (PTSTEP, LDCOMPUTE,
 USE MODD_CST,              ONLY : XP00, XRD, XRV, XMD, XMV, XCPD, XCPV, XCL, XCI, XLVTT, XLSTT, XTT, &
                                   XALPW, XBETAW, XGAMW, XALPI, XBETAI, XGAMI
 USE MODD_PARAM_LIMA,       ONLY : XRTMIN, XCTMIN,                                                    &
-                                  LCOLD, LNUCL, LSNOW, LHAIL, LWARM, LACTI, LRAIN
+                                  LCOLD, LNUCL, LSNOW, LHAIL, LWARM, LACTI, LRAIN,XNUS, XTRANS_MP_GAMMAS,XFVELOS,XLBDAS_MIN
 USE MODD_PARAM_LIMA_WARM,  ONLY : XLBC, XLBEXC, XLBR, XLBEXR
 USE MODD_PARAM_LIMA_MIXED, ONLY : XLBG, XLBEXG, XLBH, XLBEXH, XLBDAG_MAX
 USE MODD_PARAM_LIMA_COLD,  ONLY : XSCFAC, XLBI, XLBEXI, XLBS, XLBEXS, XLBDAS_MAX
@@ -502,9 +502,18 @@ WHERE (PRIT(:)>XRTMIN(4) .AND. PCIT(:)>XCTMIN(4) .AND. LDCOMPUTE(:))
    ZLBDI(:) = ( XLBI*PCIT(:) / PRIT(:) )**XLBEXI
 END WHERE
 ZLBDS(:)  = 1.E10
+!WHERE (PRST(:)>XRTMIN(5) .AND. LDCOMPUTE(:) )
+!   ZLBDS(:) = XLBS*( PRHODREF(:)*PRST(:) )**XLBEXS
+!END WHERE
 WHERE (PRST(:)>XRTMIN(5) .AND. LDCOMPUTE(:) )
-   ZLBDS(:) = XLBS*( PRHODREF(:)*PRST(:) )**XLBEXS
-END WHERE
+   WHERE(ZT(:)>263.15)
+      ZLBDS(:) = MAX(MIN(XLBDAS_MAX, 10**(14.554-0.0423*ZT(:))),XLBDAS_MIN)
+   END WHERE
+   WHERE(ZT(:)<=263.15)
+      ZLBDS(:) = MAX(MIN(XLBDAS_MAX, 10**(6.226-0.0106*ZT(:))),XLBDAS_MIN)
+   END WHERE
+END WHERE                               ! Wurtz
+ZLBDS(:) =  ZLBDS(:) * XTRANS_MP_GAMMAS ! Wurtz
 ZLBDG(:)  = 1.E10
 WHERE (PRGT(:)>XRTMIN(6) .AND. LDCOMPUTE(:) )
    ZLBDG(:) = XLBG*( PRHODREF(:)*PRGT(:) )**XLBEXG
@@ -639,7 +648,7 @@ END IF
 ! Lambda_s limited for collection processes to prevent too high concentrations
 ! must be changed or removed if C and x modified
 !
-ZLBDS(:) = MIN( XLBDAS_MAX, ZLBDS(:))
+!ZLBDS(:) = MIN( XLBDAS_MAX, ZLBDS(:))
 !
 !
 IF (LCOLD .AND. LSNOW) THEN
diff --git a/src/MNH/modd_param_lima.f90 b/src/MNH/modd_param_lima.f90
index 66156a056..5ef7eaaca 100644
--- a/src/MNH/modd_param_lima.f90
+++ b/src/MNH/modd_param_lima.f90
@@ -80,6 +80,10 @@ REAL, DIMENSION(:),    SAVE, ALLOCATABLE :: XFRAC_REF       ! AP compostion in P
 !
 ! 1.3 Ice characteristics
 !
+REAL,SAVE              :: XFVELOS !Wurtz
+REAL,SAVE              :: XLBDAS_MIN
+
+REAL,SAVE              :: XTRANS_MP_GAMMAS
 CHARACTER(LEN=4), SAVE :: CPRISTINE_ICE_LIMA ! Pristine type PLAT, COLU or BURO
 CHARACTER(LEN=4), SAVE :: CHEVRIMED_ICE_LIMA ! Heavily rimed type GRAU or HAIL
 REAL,SAVE              :: XALPHAI,XNUI,    & ! Pristine ice   distribution parameters
diff --git a/src/MNH/modd_param_lima_mixed.f90 b/src/MNH/modd_param_lima_mixed.f90
index f13accfc6..2a2b75491 100644
--- a/src/MNH/modd_param_lima_mixed.f90
+++ b/src/MNH/modd_param_lima_mixed.f90
@@ -67,6 +67,7 @@ REAL,SAVE :: XDCSLIM,XCOLCS,                   & ! Constants for the riming of
     	     XEXCRIMSS,XCRIMSS,                & ! the aggregates : RIM
     	     XEXCRIMSG,XCRIMSG,                & !
     	     XEXSRIMCG,XSRIMCG,                & !
+             XSRIMCG2, XSRIMCG3, XEXSRIMCG2,   & ! Murakami 1990
     	     XGAMINC_BOUND_MIN,                & ! Min val. of Lbda_s for RIM
     	     XGAMINC_BOUND_MAX,                & ! Max val. of Lbda_s for RIM
     	     XRIMINTP1,XRIMINTP2                 ! Csts for lin. interpol. of 
diff --git a/src/MNH/rrcolss.f90 b/src/MNH/rrcolss.f90
index 527165111..ccb2db545 100644
--- a/src/MNH/rrcolss.f90
+++ b/src/MNH/rrcolss.f90
@@ -10,7 +10,7 @@
 INTERFACE
 !
       SUBROUTINE RRCOLSS( KND, PALPHAS, PNUS, PALPHAR, PNUR,                 &
-                         PESR, PEXMASSR, PFALLS, PEXFALLS, PFALLR, PEXFALLR, &
+                         PESR, PEXMASSR, PFALLS, PEXFALLS, PFALLEXPS, PFALLR, PEXFALLR, & ! Wurtz
                          PLBDASMAX, PLBDARMAX, PLBDASMIN, PLBDARMIN,         &
                          PDINFTY, PRRCOLSS, PAG, PBS, PAS                    )
 !
@@ -28,6 +28,7 @@ REAL, INTENT(IN) :: PESR      ! Efficiency of aggregates collecting rain
 REAL, INTENT(IN) :: PEXMASSR  ! Mass exponent of rain 
 REAL, INTENT(IN) :: PFALLS    ! Fall speed constant of aggregates
 REAL, INTENT(IN) :: PEXFALLS  ! Fall speed exponent of aggregates
+REAL, INTENT(IN) :: PFALLEXPS ! Fall speed exponential of aggregates ! Wurtz
 REAL, INTENT(IN) :: PFALLR    ! Fall speed constant of rain 
 REAL, INTENT(IN) :: PEXFALLR  ! Fall speed exponent of rain 
 REAL, INTENT(IN) :: PLBDASMAX ! Maximun slope of size distribution of aggregates
@@ -49,7 +50,7 @@ END INTERFACE
       END MODULE MODI_RRCOLSS
 !     ########################################################################
       SUBROUTINE RRCOLSS( KND, PALPHAS, PNUS, PALPHAR, PNUR,                 &
-                         PESR, PEXMASSR, PFALLS, PEXFALLS, PFALLR, PEXFALLR, &
+                         PESR, PEXMASSR, PFALLS, PEXFALLS, PFALLEXPS, PFALLR, PEXFALLR, &
                          PLBDASMAX, PLBDARMAX, PLBDASMIN, PLBDARMIN,         &
                          PDINFTY, PRRCOLSS, PAG, PBS, PAS                    )
 !     ########################################################################
@@ -151,6 +152,7 @@ REAL, INTENT(IN) :: PESR      ! Efficiency of aggregates collecting rain
 REAL, INTENT(IN) :: PEXMASSR  ! Mass exponent of rain 
 REAL, INTENT(IN) :: PFALLS    ! Fall speed constant of aggregates
 REAL, INTENT(IN) :: PEXFALLS  ! Fall speed exponent of aggregates
+REAL, INTENT(IN) :: PFALLEXPS ! Fall speed exponential of aggregates ! Wurtz
 REAL, INTENT(IN) :: PFALLR    ! Fall speed constant of rain 
 REAL, INTENT(IN) :: PEXFALLR  ! Fall speed exponent of rain 
 REAL, INTENT(IN) :: PLBDASMAX ! Maximun slope of size distribution of aggregates
@@ -277,11 +279,11 @@ DO JLBDAS = 1,SIZE(PRRCOLSS(:,:),1)
           DO JDR = 1,INR-1
             ZDR = ZDDCOLLR * REAL(JDR)
             ZCOLLR = ZCOLLR + (ZDS+ZDR)**2 * ZDR**PEXMASSR                     &
-                       * PESR * ABS(PFALLS*ZDS**PEXFALLS-PFALLR*ZDR**PEXFALLR) &
+                       * PESR * ABS(PFALLS*ZDS**PEXFALLS * EXP(-(PFALLEXPS*XDS)**PALPHAS)-PFALLR*ZDR**PEXFALLR) & ! Wurtz
                                       * GENERAL_GAMMA(PALPHAR,PNUR,ZLBDAR,ZDR)
           END DO
           ZCOLLDRMAX = (ZDS+ZDRMAX)**2 * ZDRMAX**PEXMASSR                      &
-                    * PESR * ABS(PFALLS*ZDS**PEXFALLS-PFALLR*ZDRMAX**PEXFALLR) &
+                    * PESR * ABS(PFALLS*ZDS**PEXFALLS* EXP(-(PFALLEXPS*XDS)**PALPHAS)-PFALLR*ZDRMAX**PEXFALLR) & ! Wurtz
                                    * GENERAL_GAMMA(PALPHAR,PNUR,ZLBDAR,ZDRMAX)
           ZCOLLR = (ZCOLLR + 0.5*ZCOLLDRMAX)*(ZDDCOLLR/ZDDSCALR)
 !
diff --git a/src/MNH/rscolrg.f90 b/src/MNH/rscolrg.f90
index caa868e91..a1508649b 100644
--- a/src/MNH/rscolrg.f90
+++ b/src/MNH/rscolrg.f90
@@ -10,7 +10,7 @@
 INTERFACE
 !
       SUBROUTINE RSCOLRG( KND, PALPHAS, PZNUS, PALPHAR, PNUR,                &
-                         PESR, PEXMASSS, PFALLS, PEXFALLS, PFALLR, PEXFALLR, &
+                         PESR, PEXMASSS, PFALLS, PEXFALLS, PFALLEXPS, PFALLR, PEXFALLR, & ! Wurtz
                          PLBDASMAX, PLBDARMAX, PLBDASMIN, PLBDARMIN,         &
                          PDINFTY, PRSCOLRG,PAG, PBS, PAS                     )
 !
@@ -28,6 +28,7 @@ REAL, INTENT(IN) :: PESR      ! Efficiency of the aggregates collecting rain
 REAL, INTENT(IN) :: PEXMASSS  ! Mass exponent of the aggregates
 REAL, INTENT(IN) :: PFALLS    ! Fall speed constant of the aggregates
 REAL, INTENT(IN) :: PEXFALLS  ! Fall speed exponent of the aggregates
+REAL, INTENT(IN) :: PFALLEXPS  ! Fall speed exponential constant of the aggregates ! Wurtz
 REAL, INTENT(IN) :: PFALLR    ! Fall speed constant of rain 
 REAL, INTENT(IN) :: PEXFALLR  ! Fall speed exponent of rain 
 REAL, INTENT(IN) :: PLBDASMAX ! Maximun slope of size distribution of the aggregates
@@ -49,7 +50,7 @@ END INTERFACE
       END MODULE MODI_RSCOLRG
 !     ########################################################################
       SUBROUTINE RSCOLRG( KND, PALPHAS, PZNUS, PALPHAR, PNUR,                &
-                         PESR, PEXMASSS, PFALLS, PEXFALLS, PFALLR, PEXFALLR, &
+                         PESR, PEXMASSS, PFALLS, PEXFALLS, PFALLEXPS, PFALLR, PEXFALLR, &
                          PLBDASMAX, PLBDARMAX, PLBDASMIN, PLBDARMIN,         &
                          PDINFTY, PRSCOLRG,PAG, PBS, PAS                     )
 !     ########################################################################
@@ -149,6 +150,7 @@ REAL, INTENT(IN) :: PESR      ! Efficiency of the aggregates collecting rain
 REAL, INTENT(IN) :: PEXMASSS  ! Mass exponent of the aggregates
 REAL, INTENT(IN) :: PFALLS    ! Fall speed constant of the aggregates
 REAL, INTENT(IN) :: PEXFALLS  ! Fall speed exponent of the aggregates
+REAL, INTENT(IN) :: PFALLEXPS ! Fall speed exponential constant of the aggregates ! Wurtz
 REAL, INTENT(IN) :: PFALLR    ! Fall speed constant of rain 
 REAL, INTENT(IN) :: PEXFALLR  ! Fall speed exponent of rain 
 REAL, INTENT(IN) :: PLBDASMAX ! Maximun slope of size distribution of the aggregates
@@ -271,12 +273,12 @@ DO JLBDAR = 1,SIZE(PRSCOLRG(:,:),1)
             ZDR = ZDDCOLLR * REAL(JDR) + ZDRMIN
             ZCOLLR = ZCOLLR + (ZDS+ZDR)**2                                     &
                        * GENERAL_GAMMA(PALPHAR,PNUR,ZLBDAR,ZDR)                &
-                         * PESR * ABS(PFALLS*ZDS**PEXFALLS-PFALLR*ZDR**PEXFALLR)
+                         * PESR * ABS(PFALLS*ZDS**PEXFALLS*EXP(-(ZDS*PFALLEXPS)**PALPHAS)-PFALLR*ZDR**PEXFALLR) ! Wurtz
           END DO
           IF( ZDRMIN>0.0 ) THEN
             ZCOLLDRMIN = (ZDS+ZDRMIN)**2                                       &
                       * GENERAL_GAMMA(PALPHAR,PNUR,ZLBDAR,ZDRMIN)              &
-                      * PESR * ABS(PFALLS*ZDS**PEXFALLS-PFALLR*ZDRMIN**PEXFALLR)
+                      * PESR * ABS(PFALLS*ZDS**PEXFALLS*EXP(-(ZDS*PFALLEXPS)**PALPHAS)-PFALLR*ZDRMIN**PEXFALLR) ! Wurtz
             ELSE
             ZCOLLDRMIN = 0.0
           END IF 
diff --git a/src/MNH/rzcolx.f90 b/src/MNH/rzcolx.f90
index 28658241c..7aa3d97ed 100644
--- a/src/MNH/rzcolx.f90
+++ b/src/MNH/rzcolx.f90
@@ -10,7 +10,7 @@
 INTERFACE
 !
       SUBROUTINE RZCOLX( KND, PALPHAX, PNUX, PALPHAZ, PNUZ,                  &
-			 PEXZ, PEXMASSZ, PFALLX, PEXFALLX, PFALLZ, PEXFALLZ, &
+			 PEXZ, PEXMASSZ, PFALLX, PEXFALLX, PFALLEXPX, PFALLZ, PEXFALLZ, PFALLEXPZ, & ! Wurtz
 		         PLBDAXMAX, PLBDAZMAX, PLBDAXMIN, PLBDAZMIN,         &
 		         PDINFTY, PRZCOLX                                    )
 !
@@ -29,8 +29,10 @@ REAL, INTENT(IN) :: PEXZ      ! Efficiency of specy X collecting specy Z
 REAL, INTENT(IN) :: PEXMASSZ  ! Mass exponent of specy Z
 REAL, INTENT(IN) :: PFALLX    ! Fall speed constant of specy X
 REAL, INTENT(IN) :: PEXFALLX  ! Fall speed exponent of specy X
+REAL, INTENT(IN) :: PFALLEXPX ! Fall speed exponential constant of specy X ! Wurtz
 REAL, INTENT(IN) :: PFALLZ    ! Fall speed constant of specy Z
 REAL, INTENT(IN) :: PEXFALLZ  ! Fall speed exponent of specy Z
+REAL, INTENT(IN) :: PFALLEXPZ ! Fall speed exponent of specy Z  ! Wurtz
 REAL, INTENT(IN) :: PLBDAXMAX ! Maximun slope of size distribution of specy X
 REAL, INTENT(IN) :: PLBDAZMAX ! Maximun slope of size distribution of specy Z
 REAL, INTENT(IN) :: PLBDAXMIN ! Minimun slope of size distribution of specy X
@@ -49,7 +51,7 @@ END INTERFACE
       END MODULE MODI_RZCOLX
 !     ########################################################################
       SUBROUTINE RZCOLX( KND, PALPHAX, PNUX, PALPHAZ, PNUZ,                  &
-			 PEXZ, PEXMASSZ, PFALLX, PEXFALLX, PFALLZ, PEXFALLZ, &
+			 PEXZ, PEXMASSZ, PFALLX, PEXFALLX, PFALLEXPX, PFALLZ, PEXFALLZ, PFALLEXPZ, & ! Wurtz
 		         PLBDAXMAX, PLBDAZMAX, PLBDAXMIN, PLBDAZMIN,         &
 		         PDINFTY, PRZCOLX                                    )
 !     ########################################################################
@@ -152,8 +154,10 @@ REAL, INTENT(IN) :: PEXZ      ! Efficiency of specy X collecting specy Z
 REAL, INTENT(IN) :: PEXMASSZ  ! Mass exponent of specy Z
 REAL, INTENT(IN) :: PFALLX    ! Fall speed constant of specy X
 REAL, INTENT(IN) :: PEXFALLX  ! Fall speed exponent of specy X
+REAL, INTENT(IN) :: PFALLEXPX ! Fall speed exponential constant of specy X ! Wurtz
 REAL, INTENT(IN) :: PFALLZ    ! Fall speed constant of specy Z
 REAL, INTENT(IN) :: PEXFALLZ  ! Fall speed exponent of specy Z
+REAL, INTENT(IN) :: PFALLEXPZ ! Fall speed exponent of specy Z  ! Wurtz
 REAL, INTENT(IN) :: PLBDAXMAX ! Maximun slope of size distribution of specy X
 REAL, INTENT(IN) :: PLBDAZMAX ! Maximun slope of size distribution of specy Z
 REAL, INTENT(IN) :: PLBDAXMIN ! Minimun slope of size distribution of specy X
@@ -246,8 +250,8 @@ DO JLBDAX = 1,SIZE(PRZCOLX(:,:),1)
 !*       1.7     Compute the scaled fall speed difference by integration over
 !                the dimensional spectrum of specy Z
 !
-	ZCOLLZ = ZCOLLZ + ZFUNC                                               &
-		        * PEXZ * ABS(PFALLX*ZDX**PEXFALLX-PFALLZ*ZDZ**PEXFALLZ)
+        ZCOLLZ = ZCOLLZ + ZFUNC   &
+                      &  * PEXZ * ABS(PFALLX*ZDX**PEXFALLX *EXP(-(ZDX*PFALLEXPX)**PALPHAX)-PFALLZ*ZDZ**PEXFALLZ * EXP(-(ZDZ*PFALLEXPZ)**PALPHAZ)) ! Wurtz
       END DO
 !
 !*       1.8     Compute the normalization factor by integration over the
-- 
GitLab