diff --git a/src/MNH/aer_monitorn.f90 b/src/MNH/aer_monitorn.f90
index 22d5f69a9b8e052d1f28a46fe8a2d524a9e5913c..ec800306eab3ece233568cc174bb489d21e5799f 100644
--- a/src/MNH/aer_monitorn.f90
+++ b/src/MNH/aer_monitorn.f90
@@ -347,8 +347,6 @@ SELECT CASE (CCLOUD)
                               XRHODREF(IIB:IIE,IJB:IJE,IKB:IKE),   &
                               XRT(IIB:IIE,IJB:IJE,IKB:IKE,2),      &
                               XRT(IIB:IIE,IJB:IJE,IKB:IKE,3),      &
-                              ZRCS(IIB:IIE,IJB:IJE,IKB:IKE),       &
-                              ZRRS(IIB:IIE,IJB:IJE,IKB:IKE),       &
                               ZSVDST(IIB:IIE,IJB:IJE,IKB:IKE,:),   &
                               XTHT(IIB:IIE,IJB:IJE,IKB:IKE),       &
                               XPABST(IIB:IIE,IJB:IJE,IKB:IKE),     &
@@ -366,8 +364,6 @@ SELECT CASE (CCLOUD)
                               XRHODREF(IIB:IIE,IJB:IJE,IKB:IKE),   &
                               XRT(IIB:IIE,IJB:IJE,IKB:IKE,2),      &
                               XRT(IIB:IIE,IJB:IJE,IKB:IKE,3),      &
-                              ZRCS(IIB:IIE,IJB:IJE,IKB:IKE),       &
-                              ZRRS(IIB:IIE,IJB:IJE,IKB:IKE),       &
                               ZSVDST(IIB:IIE,IJB:IJE,IKB:IKE,:),   &
                               XTHT(IIB:IIE,IJB:IJE,IKB:IKE),       &
                               XPABST(IIB:IIE,IJB:IJE,IKB:IKE),     &
@@ -385,8 +381,6 @@ CASE ('LIMA')
                               XRHODREF(IIB:IIE,IJB:IJE,IKB:IKE),   &
                               XRT(IIB:IIE,IJB:IJE,IKB:IKE,2),      &
                               XRT(IIB:IIE,IJB:IJE,IKB:IKE,3),      &
-                              ZRCS(IIB:IIE,IJB:IJE,IKB:IKE),       &
-                              ZRRS(IIB:IIE,IJB:IJE,IKB:IKE),       &
                               ZSVDST(IIB:IIE,IJB:IJE,IKB:IKE,:),   &
                               XTHT(IIB:IIE,IJB:IJE,IKB:IKE),       &
                               XPABST(IIB:IIE,IJB:IJE,IKB:IKE),     &
@@ -537,8 +531,6 @@ SELECT CASE (CCLOUD)
                               XRHODREF(IIB:IIE,IJB:IJE,IKB:IKE),  &
                               XRT(IIB:IIE,IJB:IJE,IKB:IKE,2),     &
                               XRT(IIB:IIE,IJB:IJE,IKB:IKE,3),     &
-                              ZRCS(IIB:IIE,IJB:IJE,IKB:IKE),      &
-                              ZRRS(IIB:IIE,IJB:IJE,IKB:IKE),      &
                               ZSVSLT(IIB:IIE,IJB:IJE,IKB:IKE,:),  &
                               XTHT(IIB:IIE,IJB:IJE,IKB:IKE),      &
                               XPABST(IIB:IIE,IJB:IJE,IKB:IKE),    &
@@ -556,8 +548,6 @@ SELECT CASE (CCLOUD)
                               XRHODREF(IIB:IIE,IJB:IJE,IKB:IKE),   &
                               XRT(IIB:IIE,IJB:IJE,IKB:IKE,2),      &
                               XRT(IIB:IIE,IJB:IJE,IKB:IKE,3),      &
-                              ZRCS(IIB:IIE,IJB:IJE,IKB:IKE),       &
-                              ZRRS(IIB:IIE,IJB:IJE,IKB:IKE),       &
                               ZSVSLT(IIB:IIE,IJB:IJE,IKB:IKE,:),   &
                               XTHT(IIB:IIE,IJB:IJE,IKB:IKE),       &
                               XPABST(IIB:IIE,IJB:IJE,IKB:IKE),     &
@@ -574,8 +564,6 @@ SELECT CASE (CCLOUD)
                               XRHODREF(IIB:IIE,IJB:IJE,IKB:IKE),   &
                               XRT(IIB:IIE,IJB:IJE,IKB:IKE,2),      &
                               XRT(IIB:IIE,IJB:IJE,IKB:IKE,3),      &
-                              ZRCS(IIB:IIE,IJB:IJE,IKB:IKE),       &
-                              ZRRS(IIB:IIE,IJB:IJE,IKB:IKE),       &
                               ZSVSLT(IIB:IIE,IJB:IJE,IKB:IKE,:),   &
                               XTHT(IIB:IIE,IJB:IJE,IKB:IKE),       &
                               XPABST(IIB:IIE,IJB:IJE,IKB:IKE),     &
diff --git a/src/MNH/aer_wet_dep_kmt_warm.f90 b/src/MNH/aer_wet_dep_kmt_warm.f90
index 40083c0c1e3c9336156114b3148db8d5f8f99354..d0ebbdfff6edd447323c9d85a65c73b6aad43aef 100644
--- a/src/MNH/aer_wet_dep_kmt_warm.f90
+++ b/src/MNH/aer_wet_dep_kmt_warm.f90
@@ -12,7 +12,7 @@ INTERFACE
 !!
 SUBROUTINE AER_WET_DEP_KMT_WARM(KSPLITR, PTSTEP, PZZ, PRHODREF,       &
                                 PRCT, PRRT,                           &
-                                PRCS, PRRS,  PSVT, PTHT,              &
+                                PSVT, PTHT,                           &
                                 PPABST, PRGAER, PEVAP3D, KMODE,       &
                                 PDENSITY_AER, PMASSMIN, PSEA, PTOWN,  &
                                 PCCT, PCRT )
@@ -30,8 +30,6 @@ 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(INOUT) :: PSVT    ! Tracer m.r. at t
 !
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRCS    ! Cloud water conc derived from source term
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRRS    ! Rain water conc derifed from source term
 REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PEVAP3D ! Instantaneous 3D Rain Evaporation flux (KG/KG/S)
 !
 REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PTHT       !Potential temp
@@ -53,7 +51,7 @@ END MODULE MODI_AER_WET_DEP_KMT_WARM
 !     ###############################################################
       SUBROUTINE AER_WET_DEP_KMT_WARM (KSPLITR, PTSTEP, PZZ,      &
                             PRHODREF, PRCT, PRRT,                 &
-                            PRCS, PRRS,  PSVT, PTHT,              &
+                            PSVT, PTHT,                           &
                             PPABST, PRGAER, PEVAP3D, KMODE,       &
                             PDENSITY_AER, PMASSMIN, PSEA, PTOWN,  &
                             PCCT, PCRT )
@@ -123,14 +121,16 @@ END MODULE MODI_AER_WET_DEP_KMT_WARM
 !              ------------
 !
 USE MODD_CST
-USE MODD_RAIN_ICE_PARAM
+USE MODD_RAIN_ICE_PARAM, ONLY : YEXCACCR=>XEXCACCR, XFSEDC, XFCACCR,&
+                                XEXSEDR, XCRIAUTC, XFSEDR, XTIMAUTC,&
+                                YFCACCR => XFCACCR
 !++th++ 10/05/17
 USE MODD_RAIN_ICE_DESCR, ONLY : YRTMIN => XRTMIN, YCEXVT => XCEXVT, &
                                 XCONC_LAND, XCONC_SEA, XCONC_URBAN, &
                                 XNUC2, XALPHAC2, XNUC, XALPHAC,     &
                                 YLBC => XLBC, XLBEXC,               &
                                 XCCR,                               &
-                                YLBR => XLBR, YLBEXR => XLBEXR 
+                                YLBR => XLBR, YLBEXR => XLBEXR
 !--th--
 USE MODD_PRECIP_n
 USE MODI_AER_VELGRAV
@@ -145,7 +145,8 @@ USE MODD_PARAM_LIMA_WARM, ONLY : WLBR => XLBR, WLBEXR => XLBEXR,          & ! fo
                                  WLBC => XLBC,                            &
                                  XACCR1, XACCR2, XACCR3, XACCR4, XACCR5,  & ! for
                                  XACCR_RLARGE1, XACCR_RLARGE2,            & ! accr.
-                                 XACCR_RSMALL1, XACCR_RSMALL2
+                                 XACCR_RSMALL1, XACCR_RSMALL2,            &
+                                 WEXCACCR=>XEXCACCR, WFCACCR=>XFCACCR
 USE MODD_PARAM_n,         ONLY: CCLOUD
 !--th--
 
@@ -165,8 +166,6 @@ 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(INOUT) :: PSVT    ! Tracer m.r. at t
 !
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRCS    ! Cloud water m.r. from source term
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRRS    ! Rain water m.r. from source term 
 REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PEVAP3D ! Instantaneous 3D Rain Evaporation flux (KG/KG/S) 
 !
 REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PTHT     ! Potential temp
@@ -263,7 +262,7 @@ INTEGER :: IKE
 REAL, DIMENSION(:), ALLOCATABLE :: KRTMIN 
 REAL :: KCEXVT, KLBR, KLBEXR, KLBC, ZLBEXC
 REAL, DIMENSION(2) :: ZXLBC
-REAL :: ZEXSEDR, ZDR
+REAL :: ZEXSEDR, ZDR, ZEXCACCR, ZFCACCR
 !  
 !-------------------------------------------------------------------------------
 !
@@ -282,6 +281,8 @@ CASE('ICE3')
   KLBEXR = YLBEXR
   ZXLBC(:) = YLBC(:)
   ZLBEXC = XLBEXC
+  ZEXCACCR = YEXCACCR
+  ZFCACCR = YFCACCR
 CASE('LIMA')
   ALLOCATE(KRTMIN(SIZE(WRTMIN)))
   KRTMIN = WRTMIN
@@ -291,6 +292,8 @@ CASE('LIMA')
   KLBC   = WLBC
   ZLBEXC = 1.0 / 3.0
   ZDR = 0.8
+  ZEXCACCR = WEXCACCR
+  ZFCACCR = WFCACCR
 END SELECT
 !--cb--
 !
@@ -361,9 +364,7 @@ CALL AER_WET_DEP_KMT_ICE_WARM
 !
 CALL AER_WET_DEP_KMT_EVAP   
 !
-!++cb++
 DEALLOCATE(KRTMIN)
-!--cb--
 !
 !-------------------------------------------------------------------------------
 !
@@ -396,12 +397,10 @@ INTEGER                           :: JKAQ     ! counter for chemistry
 GCLOUD(:,:,:) = .FALSE.
 !
 IF (PRESENT(PCCT)) THEN ! case KHKO, C2R2, C3R5, LIMA (2-moment schemes)
-  GCLOUD(:,:,:) = PRCS(:,:,:) > KRTMIN(2) .AND. PCCT(:,:,:) > XCTMIN(2)
+  GCLOUD(:,:,:) = PRCT(:,:,:) > KRTMIN(2) .AND. PCCT(:,:,:) > XCTMIN(2)
 ELSE                    ! Case ICE3, REVE, KESS, ... (1-moment schemes)
-  GCLOUD(:,:,:) = PRCS(:,:,:) > KRTMIN(2)
+  GCLOUD(:,:,:) = PRCT(:,:,:) > KRTMIN(2)
 END IF
-!--cb--
-!--th--
 
 ICLOUD = COUNTJV( GCLOUD(:,:,:),I1C(:),I2C(:),I3C(:))
 IF( ICLOUD >= 1 ) THEN   
@@ -437,7 +436,7 @@ IF( ICLOUD >= 1 ) THEN
     ZTHT(JL)           = PTHT(I1C(JL),I2C(JL),I3C(JL))
     ZRC(JL)            = ZRAY(I1C(JL),I2C(JL),I3C(JL))
     ZPABST(JL)         = PPABST(I1C(JL),I2C(JL),I3C(JL))
-    ZRCT(JL)           = PRCS(I1C(JL),I2C(JL),I3C(JL)) 
+    ZRCT(JL)           = PRCT(I1C(JL),I2C(JL),I3C(JL)) 
     ZRHODREF(JL)       = PRHODREF(I1C(JL),I2C(JL),I3C(JL))   
     ZMASSMIN(JL,:)     = PMASSMIN(I1C(JL),I2C(JL),I3C(JL),:)   
     ZWLBDC(JL)         = ZLBC(I1C(JL),I2C(JL),I3C(JL))
@@ -711,7 +710,7 @@ DO JN = 1 , KSPLITR
         ZSVT(JL,KMODE*2+JKAQ) = PSVT(IR1(JL),IR2(JL),IR3(JL),KMODE*2+JKAQ)
       END DO
 !
-      IF (PRESENT(PCRT)) ZCRT(JL) = PCRT(IR1(JL),IR2(JL),IR2(JL))
+      IF (PRESENT(PCRT)) ZCRT(JL) = PCRT(IR1(JL),IR2(JL),IR3(JL))
       ZRRT(JL)           = PRRT(IR1(JL),IR2(JL),IR3(JL))
       ZRHODREF(JL)       = PRHODREF(IR1(JL),IR2(JL),IR3(JL))
      ENDDO
@@ -803,10 +802,12 @@ ZZRCT(:,:,:) = MAX(ZZRCT(:,:,:), KRTMIN(2)/2.)
 IF (PRESENT(PCRT)) THEN  ! 2-moment schemes
 !
 ! from lima_warm_coal.f90 (AUTO)
-  ZLBDC3(:,:,:) = XMNH_HUGE
+  ZLBDC3(:,:,:) = 1E40
+  ! ZLBDC3(:,:,:) = XMNH_HUGE
   ZLBDC(:,:,:)  = 1.E15
   WHERE (ZZRCT(:,:,:) > KRTMIN(2) .AND. PCCT(:,:,:) > XCTMIN(2))
-    ZLBDC3(:,:,:) = KLBC * PCCT(:,:,:) / PRCT(:,:,:)
+    ZLBDC3(:,:,:) = KLBC * PCCT(:,:,:) / ZZRCT(:,:,:)
+    ! ZLBDC3(:,:,:) = KLBC * PCCT(:,:,:) / PRCT(:,:,:)
     ZLBDC(:,:,:)  = ZLBDC3(:,:,:)**ZLBEXC
   END WHERE
 !
@@ -814,14 +815,14 @@ IF (PRESENT(PCRT)) THEN  ! 2-moment schemes
   WHERE (ZZRCT(:,:,:) > KRTMIN(2))
     ZZW3(:,:,:) = MAX(0.0, XLAUTR*PRHODREF(:,:,:)*ZZRCT(:,:,:)*             &
                            (XAUTO1/ZLBDC3(:,:,:)**4-XLAUTR_THRESHOLD)) ! L 
-    ZZW4(:,:,:) = MIN(PRCS(:,:,:), MAX(0.0, XITAUTR*ZZW3(:,:,:)*ZZRCT(:,:,:)*  &
+    ZZW4(:,:,:) = MIN(PRCT(:,:,:), MAX(0.0, XITAUTR*ZZW3(:,:,:)*ZZRCT(:,:,:)*  &
                            (XAUTO2/ZLBDC3(:,:,:)-XITAUTR_THRESHOLD))) ! L/tau
   END WHERE
 !
 ELSE                     ! 1-moment scheme
 !
-  WHERE ((ZZRCT(:,:,:) > KRTMIN(2)) .AND. (PRCS(:,:,:) > 0.0)) 
-    ZZW4(:,:,:) = MIN(PRCS(:,:,:), XTIMAUTC* &
+  WHERE ((ZZRCT(:,:,:) > KRTMIN(2)) .AND. (PRCT(:,:,:) > 0.0)) 
+    ZZW4(:,:,:) = MIN(PRCT(:,:,:), XTIMAUTC* &
                     MAX((ZZRCT(:,:,:)-XCRIAUTC/PRHODREF(:,:,:)), 0.0))
   END WHERE
 !
@@ -853,12 +854,14 @@ IF (PRESENT(PCRT)) THEN ! 2-moment schemes
 ! from lima_warm_coal.f90 (ACCR)
   ZLBDR3(:,:,:) = 1.E30
   ZLBDR(:,:,:)  = 1.E10
+
+
   WHERE (PRRT(:,:,:) > KRTMIN(3) .AND. PCRT(:,:,:) > XCTMIN(3))
     ZLBDAR(:,:,:) = KLBR * (PRHODREF(:,:,:) * PRRT(:,:,:))**KLBEXR
     ZLBDR3(:,:,:) = KLBR * PCRT(:,:,:) / PRRT(:,:,:)
     ZLBDR(:,:,:)  = ZLBDR3(:,:,:)**KLBEXR
-    ZZW4(:,:,:) = MIN(PRCS(:,:,:), XFCACCR * ZZRCT(:,:,:)            &
-                                         * ZLBDAR(:,:,:)**XEXCACCR &
+    ZZW4(:,:,:) = MIN(PRCT(:,:,:), ZFCACCR * ZZRCT(:,:,:)            &
+                                         * ZLBDAR(:,:,:)**ZEXCACCR &
                                          * PRHODREF(:,:,:)**(-KCEXVT) )
     ZDIM(:,:,:) = XACCR1 / ZLBDAR(:,:,:)
   END WHERE
@@ -871,7 +874,7 @@ IF (PRESENT(PCRT)) THEN ! 2-moment schemes
     ZZW5(:,:,:) = ZLBDC3(:,:,:) / ZLBDR3(:,:,:)
     ZZW1(:,:,:) = (PCCT(:,:,:) * PCRT(:,:,:) / ZLBDC3(:,:,:)**2) * PRHODREF(:,:,:)
     ZZW4(:,:,:) = MIN(ZZW1(:,:,:)*(XACCR_RLARGE1+XACCR_RLARGE2*ZZW5(:,:,:)), &
-                      PRCS(:,:,:))
+                      PRCT(:,:,:))
   END WHERE
 ! Accretion for D < 100 10-6 m
   WHERE (PRRT(:,:,:) > KRTMIN(3)  .AND. PCRT(:,:,:) > XCTMIN(3) .AND. &
@@ -881,17 +884,17 @@ IF (PRESENT(PCRT)) THEN ! 2-moment schemes
     ZZW5(:,:,:) = (ZLBDC3(:,:,:) / ZLBDR3(:,:,:))**2
     ZZW1(:,:,:) = (PCCT(:,:,:) * PCRT(:,:,:) / ZLBDC3(:,:,:)**3) * PRHODREF(:,:,:)
     ZZW4(:,:,:) = MIN(ZZW1(:,:,:)*(XACCR_RSMALL1+XACCR_RSMALL2*ZZW5(:,:,:)), &
-                      PRCS(:,:,:))
+                      PRCT(:,:,:))
   END WHERE
 !
 ELSE                    ! 1-moment schemes
 !
   ZLBDR(:,:,:) = 0.0
   WHERE ((ZZRCT(:,:,:) > KRTMIN(2)) .AND. (PRRT(:,:,:) > KRTMIN(3)) &
-                                    .AND. (PRCS(:,:,:) > 0.0))      
+                                    .AND. (PRCT(:,:,:) > 0.0))      
     ZLBDR(:,:,:) = KLBR * (PRHODREF(:,:,:) * PRRT(:,:,:))**KLBEXR
-    ZZW4(:,:,:) = MIN(PRCS(:,:,:), XFCACCR * ZZRCT(:,:,:)            &
-                                           * ZLBDR(:,:,:)**XEXCACCR &
+    ZZW4(:,:,:) = MIN(PRCT(:,:,:), ZFCACCR * ZZRCT(:,:,:)            &
+                                           * ZLBDR(:,:,:)**ZEXCACCR &
                                            * PRHODREF(:,:,:)**(-KCEXVT) )
   END WHERE
 END IF
@@ -960,7 +963,7 @@ ZWEVAP(:,:,:) = MAX(ZWEVAP(:,:,:), 0.0)
 !              no partial cloud evaporation at this stage
 !
 ZMASK(:,:,:) = 0.
-WHERE(PRCS(:,:,:) .LT. KRTMIN(2))
+WHERE(PRCT(:,:,:) .LT. KRTMIN(2))
   ZMASK(:,:,:) = 1.
 END WHERE
 !
diff --git a/src/MNH/aerocamsn.f90 b/src/MNH/aerocamsn.f90
index b3ceb1d48989bf77722b0996b3845d3735e14c80..3b7dc1cd7e87b856e1ff4401412e68b43ee00074 100644
--- a/src/MNH/aerocamsn.f90
+++ b/src/MNH/aerocamsn.f90
@@ -52,6 +52,7 @@ END MODULE MODI_AEROCAMS_n
 !!
 
 USE MODE_AERO_PSD
+USE MODD_CH_AEROSOL
 !!
 IMPLICIT NONE
 !!
@@ -72,7 +73,14 @@ INTEGER  :: JN
 
 DO JN =1,SIZE(PSV, 4)
  PSV(:,:,:,JN) = PSV(:,:,:,JN) * 1E9 / PRHODREF(:,:,:)
+ !print*, CAERONAMES(JN),' =',MINVAL(PSV(:,:,:,JN)), MAXVAL(PSV(:,:,:,JN))
 ENDDO
+PSV(:,:,:,JP_CH_BCi) = MAX(PSV(:,:,:,JP_CH_BCi), 1E-4)
+PSV(:,:,:,JP_CH_BCj) = MAX(PSV(:,:,:,JP_CH_BCj), 1E-3)
+PSV(:,:,:,JP_CH_OCi) = MAX(PSV(:,:,:,JP_CH_OCi), 1E-4)
+PSV(:,:,:,JP_CH_OCj) = MAX(PSV(:,:,:,JP_CH_OCj), 1E-3)
+PSV(:,:,:,JP_CH_DSTi) = MAX(PSV(:,:,:,JP_CH_DSTi), 1E-4)
+PSV(:,:,:,JP_CH_DSTj) = MAX(PSV(:,:,:,JP_CH_DSTj), 1E-3)
 
 ! Compute moment from aerosol mass and conversion  SV aerosols variables into ppv
 
diff --git a/src/MNH/ch_aer_coag.f90 b/src/MNH/ch_aer_coag.f90
index 87a97af5cb13d02cfe303f8e431a97330568cfb0..4e50174f9f637fed1dac8a6c648e773e8cfb21b4 100644
--- a/src/MNH/ch_aer_coag.f90
+++ b/src/MNH/ch_aer_coag.f90
@@ -14,15 +14,15 @@
 !!
 INTERFACE
 !!
-SUBROUTINE CH_AER_COAG(PM, PSIG0, PRG0, PN0,PDMINTRA,PDMINTER,PTGAS,PMU,&
-                         PLAMBDA,PRHOP0)
+SUBROUTINE CH_AER_COAG(PM,PLNSIG,PRG,PN,PDMINTRA,PDMINTER,&
+                       PTEMP,PMU,PLAMBDA,PRHOP            )
 IMPLICIT NONE
-REAL, DIMENSION(:,:), INTENT(INOUT) :: PM,PRHOP0 
+REAL, DIMENSION(:,:), INTENT(INOUT) :: PM,PRHOP 
 REAL, DIMENSION(:),   INTENT(INOUT) :: PLAMBDA, PMU
 REAL, DIMENSION(:,:), INTENT(INOUT) :: PDMINTRA
 REAL, DIMENSION(:,:), INTENT(INOUT) :: PDMINTER
-REAL, DIMENSION(:),INTENT(IN) :: PTGAS
-REAL,   DIMENSION(:,:), INTENT(IN) :: PSIG0, PRG0, PN0
+REAL, DIMENSION(:),   INTENT(IN)    :: PTEMP
+REAL, DIMENSION(:,:), INTENT(IN)    :: PLNSIG, PRG, PN
 END SUBROUTINE CH_AER_COAG
 !!
 END INTERFACE
@@ -30,59 +30,43 @@ END INTERFACE
 END MODULE MODI_CH_AER_COAG
 !!
 !!   #############################################
-     SUBROUTINE CH_AER_COAG(PM, PSIG0, PRG0, PN0,PDMINTRA,PDMINTER,PTGAS,PMU,&
-                              PLAMBDA,PRHOP0)
+     SUBROUTINE CH_AER_COAG(PM,PLNSIG,PRG,PN,PDMINTRA,PDMINTER,&
+                            PTEMP,PMU,PLAMBDA,PRHOP            )
 !!   #############################################
 !!
 !!   PURPOSE
 !!   -------
-!!
-!!   compute the terms due to Brownian, turbulent and Gravitational
-!!   coagulation:
+!!   Compute the terms due to Brownian, turbulent and gravitational
+!!   coagulation: 
 !!   a set of arrays are used to evaluate the double integral
-!!    REFERENCE
-!!    ---------
-!!    none
+!!   Based on Whitby et al. 1991 : Appendix H
 !!
-!!    AUTHOR
-!!    ------
-!!    Vincent Crassier (LA)
+!!   METHOD
+!!   ------
+!!   * Arrays of numerical evaluation of coagulation terms
+!!     in the free-molecule regime (computed from the ESMAP code)
+!!     ZINTRA     - Intamodal coagulation, mode i,j 0th and 6th Moment
+!!     ZINTER0I   - Intermodal coagulation, mode i, 0th Moment
+!!     ZINTER3I   - Intermodal coagulation, mode i, 3rd Moment
+!!     ZINTER6I   - Intermodal coagulation, mode i, 6th Moment
+!!     ZINTER6J   - Intermodal coagulation, mode j, 6th Moment
 !!
-!!    MODIFICATIONS
-!!    -------------
-!*****************************************************************
-! * Arrays of numerical evaluation of coagulation terms
-!   in the free-molecule regime (computed from the ESMAP code)
-!
-! ZINTRA     - Intamodal coagulation, mode i,j 0th and 6th Moment
-!
-! ZINTER0I   - Intermodal coagulation, mode i, 0th Moment
-! ZINTER3I   - Intermodal coagulation, mode i, 3rd Moment
-! ZINTER6I   - Intermodal coagulation, mode i, 6th Moment
-! ZINTER6J   - Intermodal coagulation, mode j, 6th Moment
-!
-! * Variables used during the coefficients evaluation
-! ZXI(i)     - Variables values at the array nodes
-! ZXINT(i)   - Variables values where the interpolation
-!             is to be made
-!
-! intramodal coagulation
-!
-! ZXINTRAMIN     - Minimal value of ln(sigma)
-! ZXINTRAMAX     - Maximal value of ln(sigma)
-! ZDXINTRA       - Step of ln(sigma) in the array
-!
-! intermodal coagulation
-!
-! ZXINTERMIN(i)  - Minimal value of the variable i
-! ZXINTERMAX(i)  - Maximal value of the variable i
-! ZDXINTER(i)    - Step of the variable i in the arrays
-!
-! i=1           - ln(sigmaj)
-! i=2           - ln(sigmai)
-! i=3           - ln((ZR=Rgj/Rgi)**2)
-!
-!***************************************************************
+!!   * Variables used during the coefficients evaluation
+!!     ZXINT(i)   - Variables values where the interpolation
+!!                  is to be made
+!!
+!!   * intramodal coagulation terms
+!!     ZXINTRAMIN     - Minimal value of ln(sigma)
+!!     ZXINTRAMAX     - Maximal value of ln(sigma)
+!!     ZDXINTRA       - Step of ln(sigma) in the array
+!!
+!!   * intermodal coagulation terms:
+!!     ZXINTERMIN(i)  - Minimal value of the variable i
+!!     ZXINTERMAX(i)  - Maximal value of the variable i
+!!     ZDXINTER(i)    - Step of the variable i in the arrays
+!!     i=1           - ln(sigmaj)
+!!     i=2           - ln(sigmai)
+!!     i=3           - ln((ZR=Rgj/Rgi)**2)
 !!
 !!   EXTERNAL
 !!   -------
@@ -90,336 +74,261 @@ END MODULE MODI_CH_AER_COAG
 !!   IMPLICIT ARGUMENTS
 !!   ------------------
 !!
-USE MODD_CH_AEROSOL
-USE MODD_CST, ONLY :    &
-       XPI              & !Definition of pi
-      ,XBOLTZ            ! Boltzman constant
+!!   REFERENCE
+!!   ---------
+!!   none
 !!
+!!   AUTHOR
+!!   ------
+!!   Vincent Crassier (LA)
+!!
+!!   MODIFICATIONS
+!!   -------------
+!!   October 2018 J. Pianezze - Add comments, cleaning and debug
+!!                              + move mode merging into ch_aer_driver
+!!
+!!
+!-------------------------------------------------------------------------------
+!
+!*       0.     DECLARATIONS
+!               ------------
+! 
+USE MODD_CH_AEROSOL
+USE MODD_CST,       ONLY : XPI, XBOLTZ
+USE MODD_CONF,      ONLY : NVERB
+!
 IMPLICIT NONE
 !
 !*       0.1   Declarations of arguments
 !
-REAL, DIMENSION(:,:), INTENT(INOUT) :: PM,PRHOP0 
+REAL, DIMENSION(:,:), INTENT(INOUT) :: PM,PRHOP 
 REAL, DIMENSION(:),   INTENT(INOUT) :: PLAMBDA, PMU
 REAL, DIMENSION(:,:), INTENT(INOUT) :: PDMINTRA
 REAL, DIMENSION(:,:), INTENT(INOUT) :: PDMINTER
-REAL, DIMENSION(:),INTENT(IN) :: PTGAS
-REAL,   DIMENSION(:,:), INTENT(IN) :: PSIG0, PRG0, PN0
+REAL, DIMENSION(:),   INTENT(IN)    :: PTEMP
+REAL, DIMENSION(:,:), INTENT(IN)    :: PLNSIG, PRG, PN
 !
 !*       0.2   Declarations of local variables
 !
-INTEGER :: JI,JJ
+INTEGER                             :: JI,JJ
 !
-REAL :: ZTURBDS ! Rate of dissipation of kinetic energy per unit mass (m2/s3)
+REAL, DIMENSION(SIZE(PM,1))         :: ZKFM,ZKNC
+REAL, DIMENSION(SIZE(PM,1))         :: ZR,ZR2,ZR4
+REAL, DIMENSION(SIZE(PM,1))         :: ZRM,ZRM2,ZRM3
+REAL, DIMENSION(SIZE(PM,1))         :: ZKNG
+REAL, DIMENSION(SIZE(PM,1))         :: ZAI,ZKNGI,ZAJ,ZKNGJ
+REAL, DIMENSION(SIZE(PM,1))         :: ZINTRA0NC,ZINTRA0FM,ZINTRA0
+REAL, DIMENSION(SIZE(PM,1))         :: ZINTRA6NC,ZINTRA6FM,ZINTRA6
+REAL, DIMENSION(SIZE(PM,1))         :: ZINTERNC,ZINTERFM,ZINTER
+REAL, DIMENSION(SIZE(PM,1))         :: ZAPPROX
 !
-REAL, DIMENSION(SIZE(PM,1)) :: ZKFM,ZKNC
-!REAL, DIMENSION(SIZE(PM,1)) :: ZKTURB,ZKGRAV,ZR3,ZRM4
-REAL, DIMENSION(SIZE(PM,1)) :: ZR,ZR2,ZR4
-REAL, DIMENSION(SIZE(PM,1)) :: ZRM,ZRM2,ZRM3
-REAL, DIMENSION(SIZE(PM,1)) :: ZKNG
-REAL, DIMENSION(SIZE(PM,1)) :: ZAI,ZKNGI,ZAJ,ZKNGJ
-REAL, DIMENSION(SIZE(PM,1)) :: ZINTRA0NC,ZINTRA0FM,ZINTRA0
-REAL, DIMENSION(SIZE(PM,1)) :: ZINTRA3NC,ZINTRA3FM,ZINTRA3
-REAL, DIMENSION(SIZE(PM,1)) :: ZINTRA6NC,ZINTRA6FM,ZINTRA6
-REAL, DIMENSION(SIZE(PM,1)) :: ZINTERNC,ZINTERFM,ZINTER
-REAL, DIMENSION(SIZE(PM,1)) :: ZAPPROX
+REAL, DIMENSION(SIZE(PM,1))         :: ZRGJ, ZRGI, ZRG
 !
-REAL, DIMENSION(SIZE(PM,1)) :: ZA,ZB,ZC,ZD
-REAL, DIMENSION(SIZE(PM,1)) :: ZRGJ, ZRGI, ZRG
+REAL, DIMENSION(SIZE(PM,1),JPMODE)  :: ZESG01,ZESG04,ZESG05,ZESG08,ZESG09
+REAL, DIMENSION(SIZE(PM,1),JPMODE)  :: ZESG12,ZESG16
+REAL, DIMENSION(SIZE(PM,1),JPMODE)  :: ZESG20,ZESG24,ZESG25,ZESG28
+REAL, DIMENSION(SIZE(PM,1),JPMODE)  :: ZESG32,ZESG36
+REAL, DIMENSION(SIZE(PM,1),JPMODE)  :: ZESG49
+REAL, DIMENSION(SIZE(PM,1),JPMODE)  :: ZESG52
+REAL, DIMENSION(SIZE(PM,1),JPMODE)  :: ZESG64
+REAL, DIMENSION(SIZE(PM,1),JPMODE)  :: ZESG81,ZESG85
+REAL, DIMENSION(SIZE(PM,1),JPMODE)  :: ZESG100,ZESG121,ZESG144,ZESG169,ZESG196
+REAL, DIMENSION(SIZE(PM,1),JPMODE)  :: ZESG256
+REAL, DIMENSION(SIZE(PM,1))         :: ZRB0,ZRB6
+REAL, DIMENSION(SIZE(PM,1))         :: ZRES
 !
-REAL, DIMENSION(SIZE(PM,1)) :: ZERF0,ZPHI0,ZXi,ZSOL
-REAL, DIMENSION(SIZE(PM,1)) :: ZERF3,ZPHI3
-REAL, DIMENSION(SIZE(PM,1)) :: ZERF6,ZPHI6
-REAL, DIMENSION(SIZE(PM,1),JPMODE) :: ZINVSIG,ZLNDG
+!-------------------------------------------------------------------------------
 !
-REAL, DIMENSION(SIZE(PM,1),JPMODE) :: ZESG01,ZESG04,ZESG05,ZESG08,ZESG09
-REAL, DIMENSION(SIZE(PM,1),JPMODE) :: ZESG12,ZESG16
-REAL, DIMENSION(SIZE(PM,1),JPMODE) :: ZESG20,ZESG24,ZESG25,ZESG28
-REAL, DIMENSION(SIZE(PM,1),JPMODE) :: ZESG32,ZESG36
-REAL, DIMENSION(SIZE(PM,1),JPMODE) :: ZESG49
-REAL, DIMENSION(SIZE(PM,1),JPMODE) :: ZESG52
-REAL, DIMENSION(SIZE(PM,1),JPMODE) :: ZESG64
-REAL, DIMENSION(SIZE(PM,1),JPMODE) :: ZESG81,ZESG85
-REAL, DIMENSION(SIZE(PM,1),JPMODE) :: ZESG100,ZESG121,ZESG144,ZESG169,ZESG196
-REAL, DIMENSION(SIZE(PM,1),JPMODE) :: ZESG256
-REAL, DIMENSION(SIZE(PM,1)) :: ZRB0,ZRB6
-REAL, DIMENSION(SIZE(PM,1)) :: ZRES
-!-------------------------------------------------------------------------------!
-ZTURBDS=0.001
-ZKNC(:)=2.*XBOLTZ*PTGAS(:)/(3.*PMU(:))
+!*       1.    INITIALIZATION
+!              --------------
 !
-PDMINTRA(:,:)=0.
-PDMINTER(:,:)=0.
+PDMINTRA(:,:) = 0.0
+PDMINTER(:,:) = 0.0
 !
-!****************************************************************
-! Initialisation des variables utilisees dans le calcul des
-! coefficients de coagulation
-!****************************************************************
-
-ZESG01(:,:) = exp(0.125*PSIG0(:,1:JPMODE)**2)            
-ZESG04(:,:)  = ZESG01(:,:) ** 4            
-ZESG05(:,:)  = ZESG04(:,:) * ZESG01(:,:)                        
-ZESG08(:,:)  = ZESG04(:,:) * ZESG04(:,:)                        
-ZESG09(:,:)  = ZESG04(:,:) * ZESG05(:,:)            
-ZESG12(:,:)  = ZESG04(:,:) * ZESG04(:,:) * ZESG04(:,:)
-ZESG16(:,:)  = ZESG08(:,:) * ZESG08(:,:)
-ZESG20(:,:)  = ZESG16(:,:) * ZESG04(:,:)
-ZESG24(:,:)  = ZESG12(:,:) * ZESG12(:,:)
-ZESG25(:,:)  = ZESG16(:,:) * ZESG09(:,:)
-ZESG28(:,:)  = ZESG20(:,:) * ZESG08(:,:)
-ZESG32(:,:)  = ZESG16(:,:) * ZESG16(:,:)
-ZESG36(:,:)  = ZESG16(:,:) * ZESG20(:,:)
-ZESG49(:,:)  = ZESG25(:,:) * ZESG20(:,:) * ZESG04(:,:)
-ZESG52(:,:)  = ZESG16(:,:) * ZESG36(:,:)
-ZESG64(:,:)  = ZESG32(:,:) * ZESG32(:,:)
-ZESG81(:,:)  = ZESG49(:,:) * ZESG32(:,:)
-ZESG85(:,:)  = ZESG64(:,:) * ZESG20(:,:) * ZESG01(:,:)
-ZESG100(:,:) = ZESG36(:,:) * ZESG64(:,:) 
-ZESG121(:,:) = ZESG85(:,:) * ZESG36(:,:)
-ZESG144(:,:) = ZESG100(:,:) * ZESG36(:,:) * ZESG08(:,:)
-ZESG169(:,:) = ZESG144(:,:) * ZESG25(:,:)
-ZESG196(:,:) = ZESG144(:,:) * ZESG52(:,:)
+ZKNC(:)       = 2.*XBOLTZ*PTEMP(:)/(3.*PMU(:))
+!
+! Compute coagulation coefficients : Whitby et al. 1991 : Appendix H
+!
+ZESG01 (:,:) = EXP(0.125*PLNSIG(:,:)**2)            
+ZESG04 (:,:) = ZESG01 (:,:) ** 4            
+ZESG05 (:,:) = ZESG04 (:,:) * ZESG01 (:,:)                        
+ZESG08 (:,:) = ZESG04 (:,:) * ZESG04 (:,:)                        
+ZESG09 (:,:) = ZESG04 (:,:) * ZESG05 (:,:)            
+ZESG12 (:,:) = ZESG04 (:,:) * ZESG04 (:,:) * ZESG04(:,:)
+ZESG16 (:,:) = ZESG08 (:,:) * ZESG08 (:,:)
+ZESG20 (:,:) = ZESG16 (:,:) * ZESG04 (:,:)
+ZESG24 (:,:) = ZESG12 (:,:) * ZESG12 (:,:)
+ZESG25 (:,:) = ZESG16 (:,:) * ZESG09 (:,:)
+ZESG28 (:,:) = ZESG20 (:,:) * ZESG08 (:,:)
+ZESG32 (:,:) = ZESG16 (:,:) * ZESG16 (:,:)
+ZESG36 (:,:) = ZESG16 (:,:) * ZESG20 (:,:)
+ZESG49 (:,:) = ZESG25 (:,:) * ZESG20 (:,:) * ZESG04(:,:)
+ZESG52 (:,:) = ZESG16 (:,:) * ZESG36 (:,:)
+ZESG64 (:,:) = ZESG32 (:,:) * ZESG32 (:,:)
+ZESG81 (:,:) = ZESG49 (:,:) * ZESG32 (:,:)
+ZESG85 (:,:) = ZESG64 (:,:) * ZESG20 (:,:) * ZESG01(:,:)
+ZESG100(:,:) = ZESG36 (:,:) * ZESG64 (:,:) 
+ZESG121(:,:) = ZESG85 (:,:) * ZESG36 (:,:)
+ZESG144(:,:) = ZESG100(:,:) * ZESG36 (:,:) * ZESG08(:,:)
+ZESG169(:,:) = ZESG144(:,:) * ZESG25 (:,:)
+ZESG196(:,:) = ZESG144(:,:) * ZESG52 (:,:)
 ZESG256(:,:) = ZESG144(:,:) * ZESG100(:,:) * ZESG12(:,:)
-
-!***************************************************************
-! Transfert de moments entre les modes i et j
-!***************************************************************
-
-ZINVSIG(:,:)=1./PSIG0(:,1:JPMODE)**2
-ZLNDG(:,:)=log(2.*PRG0(:,1:JPMODE))
-
-ZA(:)=0.5*(ZINVSIG(:,1)-ZINVSIG(:,2))
-ZD(:) = 0.
-ZXi(:)= 0.
-
-WHERE (ABS(ZA(:)) > 1E-4)
-  ZB(:)=ZINVSIG(:,2)*ZLNDG(:,2)-ZINVSIG(:,1)*ZLNDG(:,1)
-  ZC(:)=0.5*(ZINVSIG(:,1)*ZLNDG(:,1)**2-ZINVSIG(:,2)*ZLNDG(:,2)**2) - &
-       &log((PN0(:,1)*PSIG0(:,2))/(PN0(:,2)*PSIG0(:,1)))
-
-  ZD(:)=ZB(:)**2-4.*ZA(:)*ZC(:)
-
-  ZSOL(:)=(-ZB(:)+sqrt(ABS(ZD(:))))/(2.*ZA(:))
-  WHERE (ZSOL(:) < 5.E+2)
-    ZSOL(:)=exp(ZSOL(:))/2.
-    ZXi(:)=log(ZSOL(:)/PRG0(:,1))/(sqrt(2.)*PSIG0(:,1))
-  ENDWHERE
-ENDWHERE
-
-!*********************************************************************
-!      calculate the intramodal moment coefficients (log-normal model)
-!*********************************************************************
-       
-do JI=1,JPMODE
-
-  ZKFM(:)=sqrt(3.*XBOLTZ*PTGAS(:)/PRHOP0(:,JI))*1.e-3
-  !ZKTURB(:)=sqrt(XPI*ZTURBDS*PMU(:)/(120.*PRHOP0(:,JI)))*1.e-18
-  !ZKGRAV(:)=1.5/4.*0.544*XPI*PRHOP0(:,JI)/PMU(:)*1.e-24
-!*************************************************************
-!      calculate ZVG,ln2(sigma) and sigma
-!      (log-normal model)
-!*************************************************************
-
-  ZRG(:)=PRG0(:,JI)
-  ZKNG(:)=PLAMBDA(:)/ZRG(:)
-  ZAI(:)=1.392*ZKNG(:)**0.0783
-
-!***********************
-! Brownian Coagulation  
-!***********************
-       
-  ZRB0(:)=0.8
-  ZRB6(:)=ZRB0
-       
-  ZINTRA0FM(:)=ZKFM(:)*ZRB0(:)*sqrt(2.*ZRG(:))*(ZESG01(:,JI)+ZESG25(:,JI)+2.*ZESG05(:,JI))
-  ZINTRA3FM(:)=ZKFM(:)*ZRB0(:)*sqrt(ZRG(:))**7*sqrt(2.)*(ZESG49(:,JI)+ZESG36(:,JI)*ZESG01(:,JI)+&
-              &2.*ZESG25(:,JI)*ZESG04(:,JI)+ZESG09(:,JI)*ZESG16(:,JI)+ZESG100(:,JI)*ZESG09(:,JI)+&
-              &2.*ZESG64(:,JI)*ZESG01(:,JI))
-  ZINTRA6FM(:)=ZKFM(:)*ZRB6(:)*sqrt(ZRG(:))**13*sqrt(2.)*ZESG85(:,JI)*&
-              (1.+2.*ZESG04(:,JI)+ZESG24(:,JI))
-  ZINTRA0NC(:)=ZKNC(:)*(1.+ZESG08(:,JI)+ZAI(:)*ZKNG(:)*(ZESG20(:,JI)+ZESG04(:,JI)))
-  ZINTRA3NC(:)=ZKNC(:)*ZRG(:)**3*(2.*ZESG36(:,JI)+ZAI(:)*ZKNG(:)*(ZESG16(:,JI)+ZESG04(:,JI)*ZESG04(:,JI)+&
-              &ZESG36(:,JI)*ZESG04(:,JI)+ZESG64(:,JI)*ZESG16(:,JI))+ZESG16(:,JI)*ZESG04(:,JI)+&
-              &ZESG64(:,JI)*ZESG04(:,JI))
-  ZINTRA6NC(:)=2.*ZKNC(:)*(ZRG(:))**6*ZESG52(:,JI)*(ZESG20(:,JI)+ZESG28(:,JI)+ZAI(:)*ZKNG(:)*(1.+ZESG16(:,JI)))
-  ZINTRA0(:)=ZINTRA0FM(:)*(ZINTRA0NC(:)/(ZINTRA0FM(:)+ZINTRA0NC(:)))
-  ZINTRA3(:)=ZINTRA3FM(:)*(ZINTRA3NC(:)/(ZINTRA3FM(:)+ZINTRA3NC(:)))
-  ZINTRA6(:)=ZINTRA6FM(:)*(ZINTRA6NC(:)/(ZINTRA6FM(:)+ZINTRA6NC(:)))
-   
-  PDMINTRA(:,NM0(JI))=ZINTRA0(:)
-  PDMINTRA(:,NM3(JI))=ZINTRA3(:)
-  PDMINTRA(:,NM6(JI))=ZINTRA6(:)
-  !print*,'PDMINTRA(:,NM0(',JI,') =',MINVAL(PDMINTRA(:,NM0(JI))), MAXVAL(PDMINTRA(:,NM0(JI)))
-  !print*,'PDMINTRA(:,NM3(',JI,') =',MINVAL(PDMINTRA(:,NM3(JI))), MAXVAL(PDMINTRA(:,NM3(JI)))
-  !print*,'PDMINTRA(:,NM6(',JI,') =',MINVAL(PDMINTRA(:,NM6(JI))), MAXVAL(PDMINTRA(:,NM6(JI)))
-
-enddo
-!print*,'=============================='
-!print*,'=============================='
-
-WHERE (ZD(:) > 0. .AND. ZXi(:) > (6.*PSIG0(:,1)/sqrt(2.)))
-
-! transfert du moment d'ordre 0 (nombre)
-!**************************************
-
-  ZERF0(:)=sqrt(1.-exp(-4.*(ZXi(:))**2/XPI))
-  ZPHI0(:)=0.5*(1.+ZERF0(:))
-
-! transfert du moment d'ordre 3 (masse)
-!**************************************
-
-  ZERF3(:)=sqrt(1.-exp(-4.*(ZXi(:)-3.*PSIG0(:,1)/sqrt(2.))**2/XPI))
-  ZPHI3(:)=0.5*(1.+ZERF3(:))
-  
-! transfert du moment d'ordre 6 (dispersion)
-!**************************************
-
-  ZERF6(:)=sqrt(1.-exp(-4.*(ZXi(:)-6.*PSIG0(:,1)/sqrt(2.))**2/XPI))
-  ZPHI6(:)=0.5*(1.+ZERF6(:))
-  
-  PDMINTRA(:,NM0(2))=PDMINTRA(:,NM0(2))-(1.-ZPHI0(:)**2)*PDMINTRA(:,NM0(1))*(PM(:,NM0(1))/PM(:,NM0(2)))**2
-  PDMINTRA(:,NM0(1))=(2.-ZPHI0(:)**2)*PDMINTRA(:,NM0(1))
-
-  PDMINTRA(:,NM3(2))=PDMINTRA(:,NM3(1))*(1.-ZPHI0(:)*ZPHI3(:))*PM(:,NM0(1))**2
-  PDMINTRA(:,NM3(1))=PDMINTRA(:,NM3(1))*(ZPHI0(:)*ZPHI3(:)-1.)*PM(:,NM0(1))**2
-  
-  ZKFM(:)=sqrt(3.*XBOLTZ*PTGAS(:)/PRHOP0(:,1))*1.e-3
-  ZRG(:)=PRG0(:,1)
-  ZKNG(:)=PLAMBDA(:)/ZRG(:)
-  ZAI(:)=1.392*ZKNG(:)**0.0783
-  
-  ZINTRA6FM(:)=ZKFM(:)*sqrt(2.)*sqrt(ZRG(:))**13*(ZESG169(:,1)+ZESG144(:,1)*ZESG01(:,1)+&
-               2.*ZESG121(:,1)*ZESG04(:,1)+ZESG81(:,1)*ZESG16(:,1)+&
-               ZESG256(:,1)*ZESG09(:,1)+ZESG196(:,1)*ZESG01(:,1))
-
-  ZINTRA6NC(:)=ZKNC(:)*(ZRG(:))**6*(2.*ZESG144(:,1)+ZAI(:)*ZKNG(:)*(ZESG100(:,1)+&
-             ZESG64(:,1)*ZESG04(:,1))+ZAI(:)*ZKNG(:)*(ZESG144(:,1)*ZESG04(:,1)+&
-             ZESG196(:,1)*ZESG16(:,1))+ZESG100(:,1)*ZESG04(:,1)+&
-             ZESG196(:,1)*ZESG04(:,1))
-             
-  ZINTRA6(:)=ZINTRA6FM(:)*(ZINTRA6NC(:)/(ZINTRA6FM(:)+ZINTRA6NC(:)))
-
-  PDMINTRA(:,NM6(2))=PDMINTRA(:,NM6(2))+(PDMINTRA(:,NM6(1))*(1.-ZPHI3(:)**2)+ZINTRA6(:)*(1.-ZPHI0(:)*ZPHI6(:)))&
-                    &*(PM(:,NM0(1))/PM(:,NM0(2)))**2
-                   
-  PDMINTRA(:,NM6(1))=PDMINTRA(:,NM6(1))*(ZPHI3(:)**2)+ZINTRA6(:)*(ZPHI0(:)*ZPHI6(:)-1.)
-  
-ELSEWHERE
-
-  PDMINTRA(:,NM3(1))=0.
-  PDMINTRA(:,NM3(2))=0.
-
-
-ENDWHERE
-
-do JI=1,JPMODE
-!print*,'2.-ZPHI0(:)**2 =',MINVAL(2.-ZPHI0(:)**2), MAXVAL(2.-ZPHI0(:)**2)
-!  print*,'apres corr PDMINTRA(:,NM0(',JI,') =',MINVAL(PDMINTRA(:,NM0(JI))), MAXVAL(PDMINTRA(:,NM0(JI)))
-!  print*,'apres corr PDMINTRA(:,NM3(',JI,') =',MINVAL(PDMINTRA(:,NM3(JI))), MAXVAL(PDMINTRA(:,NM3(JI)))
-!  print*,'apres corr PDMINTRA(:,NM6(',JI,') =',MINVAL(PDMINTRA(:,NM6(JI))), MAXVAL(PDMINTRA(:,NM6(JI)))
- enddo
-
-!*********************************************************************
-!   calculate the intermodal moment coefficients (log-normal model)
-!*********************************************************************
-
-do JI=1,(JPMODE-1)
-  do JJ=(JI+1),JPMODE
-
-    ZRGI(:)=PRG0(:,JI)
-    ZKNGI(:)=PLAMBDA(:)/ZRGI(:)
-    ZAI(:)=1.392*ZKNGI(:)**0.0783
-
-    ZRGJ(:)=PRG0(:,JJ)
-    ZKNGJ(:)=PLAMBDA(:)/ZRGJ(:)
-    ZAJ(:)=1.392*ZKNGJ(:)**0.0783
-          
-    ZR(:)=sqrt(ZRGJ(:)/ZRGI(:))
-    ZR2(:)=ZR(:)*ZR(:)
-    !ZR3(:)=ZR(:)*ZR2(:)
-    ZR4(:)=ZR2(:)*ZR2(:)
-    ZRM(:)=1./ZR(:)
-    ZRM2(:)=ZRM(:)*ZRM(:)
-    ZRM3(:)=ZRM(:)*ZRM2(:)
-    !ZRM4(:)=ZRM2(:)*ZRM2(:)
-
-!**********************
-! Brownian Coagulation
-!**********************
-
-      ZRES(:)=0.9
-
-      ZAPPROX(:)=sqrt(2.*ZRGI(:))*(ZESG01(:,JI)+ZR(:)*ZESG01(:,JJ)+2.*ZR2(:)*ZESG01(:,JI)*ZESG04(:,JJ)&
-                 +ZR4(:)*ZESG09(:,JI)*ZESG16(:,JJ)+ZRM3(:)*ZESG16(:,JI)*ZESG09(:,JJ)+&
-                 2.*ZRM(:)*ZESG04(:,JI)*ZESG01(:,JJ))
+!
+!-------------------------------------------------------------------------------
+!
+!*       2.    COMPUTE INTRA-MODAL COAGULATION TERMS
+!              -------------------------------------
+!
+DO JI=1,JPMODE
+  !
+  !*          2.0 INITIALIZATION
+  !               --------------
+  !
+  ZRG(:)  = PRG(:,JI)
+  ZKNG(:) = PLAMBDA(:)/ZRG(:)
+  ZAI(:)  = 1.392*ZKNG(:)**0.0783
+  ZKFM(:) = SQRT(3.*XBOLTZ*PTEMP(:)/PRHOP(:,JI))*1.E-3
+  !   
+  ZRB0(:) = 0.8
+  ZRB6(:) = ZRB0
+  !    
+  !*          2.1 FREE-MOLECULE REGIME (KN > 10) 
+  !               ------------------------------
+  !
+  ZINTRA0FM(:) = ZKFM(:) * ZRB0(:) * SQRT(2.*ZRG(:)) *     &
+               & (ZESG01(:,JI)+ZESG25(:,JI)+2.*ZESG05(:,JI))
+  ZINTRA6FM(:) = ZKFM(:)*ZRB6(:)*sqrt(ZRG(:))**13*sqrt(2.)*ZESG85(:,JI) * &
+               & (1.+2.*ZESG04(:,JI)+ZESG24(:,JI))
+  !
+  !*          2.2 NEAR-CONTINUUM (KN < 0.1)
+  !               -------------------------
+  !
+  ZINTRA0NC(:) = ZKNC(:)*(1.+ZESG08(:,JI)+ZAI(:)*ZKNG(:)*(ZESG20(:,JI)+ZESG04(:,JI)))
+  ZINTRA6NC(:) = 2.*ZKNC(:)*(ZRG(:))**6*ZESG52(:,JI)*(ZESG20(:,JI)+ZESG28(:,JI)+&
+                 ZAI(:)*ZKNG(:)*(1.+ZESG16(:,JI)))
+  !
+  !*         2.3 HARMONIC MEAN
+  !              -------------
+  !
+  ZINTRA0(:) = ZINTRA0FM(:)*(ZINTRA0NC(:) / (ZINTRA0FM(:)+ZINTRA0NC(:)))
+  ZINTRA6(:) = ZINTRA6FM(:)*(ZINTRA6NC(:) / (ZINTRA6FM(:)+ZINTRA6NC(:)))
+  !
+  PDMINTRA(:,NM0(JI)) = - ZINTRA0(:) * PM(:,NM0(JI))**2.0
+  PDMINTRA(:,NM3(JI)) =   0.0
+  PDMINTRA(:,NM6(JI)) =   ZINTRA6(:) * PM(:,NM0(JI))**2.0
+  !
+ENDDO
+!
+!-------------------------------------------------------------------------------
+!
+!*       3.    COMPUTE INTER-MODAL COAGULATION TERMS
+!              -------------------------------------
+!
+! JPMODE = 2
+! donc :
+!   - JI = 1
+!   - JJ = 2
+!
+DO JI=1,(JPMODE-1)
+  DO JJ=(JI+1),JPMODE
+    !
+    ZRGI (:) = PRG(:,JI)
+    ZKNGI(:) = PLAMBDA(:)/ZRGI(:)
+    ZAI  (:) = 1.392*ZKNGI(:)**0.0783
+    !
+    ZRGJ (:) = PRG(:,JJ)
+    ZKNGJ(:) = PLAMBDA(:)/ZRGJ(:)
+    ZAJ  (:) = 1.392*ZKNGJ(:)**0.0783
+    !      
+    ZR  (:)  = SQRT(ZRGJ(:)/ZRGI(:))
+    ZR2 (:)  = ZR(:)*ZR(:)
+    ZR4 (:)  = ZR2(:)*ZR2(:)
+    ZRM (:)  = 1./ZR(:)
+    ZRM2(:)  = ZRM(:)*ZRM(:)
+    ZRM3(:)  = ZRM(:)*ZRM2(:)
+    !
+    !    
+    ! * 3.1 Free-Molecule Regime (Kn > 10) : appendix H.2.2.1 - Whitby et al. 1991
+    !       ----------------------------------------------------------------------
+    !
+    ZRES(:)=0.9
+    !
+    !       moment 0
+    !
+    ZAPPROX(:)=sqrt(2.*ZRGI(:))*(ZESG01(:,JI)+ZR(:)*ZESG01(:,JJ)+2.*ZR2(:)*ZESG01(:,JI)*ZESG04(:,JJ)&
+               +ZR4(:)*ZESG09(:,JI)*ZESG16(:,JJ)+ZRM3(:)*ZESG16(:,JI)*ZESG09(:,JJ)+&
+               2.*ZRM(:)*ZESG04(:,JI)*ZESG01(:,JJ))
        
-      ZINTERFM(:)=ZKFM(:)*ZRES(:)*ZAPPROX(:)
+    ZINTERFM(:)=ZKFM(:)*ZRES(:)*ZAPPROX(:)
 
-      ZAPPROX(:)=2.+ZAI(:)*ZKNGI(:)*(ZESG04(:,JI)+ZR2(:)*ZESG16(:,JI)*ZESG04(:,JJ))+&
+    ZAPPROX(:)=2.+ZAI(:)*ZKNGI(:)*(ZESG04(:,JI)+ZR2(:)*ZESG16(:,JI)*ZESG04(:,JJ))+&
                 ZAJ(:)*ZKNGJ(:)*(ZESG04(:,JJ)+ZRM2(:)*ZESG16(:,JJ)*ZESG04(:,JI))+&
                 (ZR2(:)+ZRM2(:))*(ZESG04(:,JI)*ZESG04(:,JJ))
 
-      ZINTERNC(:)=ZKNC(:)*ZAPPROX(:)
+    ZINTERNC(:)=ZKNC(:)*ZAPPROX(:)
 
-      ZINTER(:)=ZINTERNC(:)*(ZINTERFM(:)/(ZINTERNC(:)+ZINTERFM(:)))
+    ZINTER(:)=ZINTERNC(:)*(ZINTERFM(:)/(ZINTERNC(:)+ZINTERFM(:)))
 
-      PDMINTER(:,NM0(JI))=PM(:,NM0(JJ))*ZINTER(:)
-      PDMINTER(:,NM0(JJ))=-PM(:,NM0(JJ))*ZINTER(:)
-
-      ZAPPROX(:)=sqrt(2.)*sqrt(ZRGI(:))**7*(ZESG49(:,JI)+ZR(:)*ZESG36(:,JI)*ZESG01(:,JJ)+2.*ZR2(:)*&
+    PDMINTER(:,NM0(JI))=-PM(:,NM0(JI))*PM(:,NM0(JJ))*ZINTER(:)
+    PDMINTER(:,NM0(JJ))= 0.0
+    !
+    !       moment 3
+    !
+    ZAPPROX(:)=sqrt(2.)*sqrt(ZRGI(:))**7*(ZESG49(:,JI)+ZR(:)*ZESG36(:,JI)*ZESG01(:,JJ)+2.*ZR2(:)*&
                  ZESG25(:,JI)*ZESG04(:,JJ)+ZR4(:)*ZESG09(:,JI)*ZESG16(:,JJ)+ZRM3(:)*&
                  ZESG100(:,JI)*ZESG09(:,JJ)+2.*ZRM(:)*ZESG64(:,JI)*ZESG01(:,JJ))
 
-      ZINTERFM(:)=ZKFM(:)*ZRES(:)*ZAPPROX(:)
+    ZINTERFM(:)=ZKFM(:)*ZRES(:)*ZAPPROX(:)
 
-      ZAPPROX(:)=(2.*ZESG36(:,JI)+ZAI(:)*ZKNGI(:)*(ZESG16(:,JI)+ZR2(:)*ZESG04(:,JI)*ZESG04(:,JJ))+&
-      ZAJ(:)*ZKNGJ(:)*(ZESG36(:,JI)*ZESG04(:,JJ)+ZRM2(:)*ZESG16(:,JJ)*ZESG64(:,JI))+&
-      ZR2(:)*ZESG16(:,JI)*ZESG04(:,JJ)+ZRM2(:)*ZESG64(:,JI)*ZESG04(:,JJ))*(ZRGI(:))**3      
+    ZAPPROX(:)=(2.*ZESG36(:,JI)+ZAI(:)*ZKNGI(:)*(ZESG16(:,JI)+ZR2(:)*ZESG04(:,JI)*ZESG04(:,JJ))+&
+    ZAJ(:)*ZKNGJ(:)*(ZESG36(:,JI)*ZESG04(:,JJ)+ZRM2(:)*ZESG16(:,JJ)*ZESG64(:,JI))+&
+    ZR2(:)*ZESG16(:,JI)*ZESG04(:,JJ)+ZRM2(:)*ZESG64(:,JI)*ZESG04(:,JJ))*(ZRGI(:))**3      
 
-      ZINTERNC(:)=ZKNC(:)*ZAPPROX(:)
+    ZINTERNC(:)=ZKNC(:)*ZAPPROX(:)
 
-      ZINTER(:)=ZINTERNC(:)*(ZINTERFM(:)/(ZINTERNC(:)+ZINTERFM(:)))
+    ZINTER(:)=ZINTERNC(:)*(ZINTERFM(:)/(ZINTERNC(:)+ZINTERFM(:)))
        
-      PDMINTER(:,NM3(JI))=-PM(:,NM0(JI))*PM(:,NM0(JJ))*ZINTER(:)
-      PDMINTER(:,NM3(JJ))=PM(:,NM0(JI))*PM(:,NM0(JJ))*ZINTER(:)
+    PDMINTER(:,NM3(JI))=-PM(:,NM0(JI))*PM(:,NM0(JJ))*ZINTER(:)
+    PDMINTER(:,NM3(JJ))=-PDMINTER(:,NM3(JI))
+
+    !       moment 6
 
-      ZAPPROX(:)=sqrt(2.)*sqrt(ZRGI(:))**13*(ZESG169(:,JI)+ZR(:)*ZESG144(:,JI)*ZESG01(:,JJ)+&
-             2.*ZR2(:)*ZESG121(:,JI)*ZESG04(:,JJ)+ZR4(:)*ZESG81(:,JI)*ZESG16(:,JJ)+&
-             ZRM3(:)*ZESG256(:,JI)*ZESG09(:,JJ)+2*ZRM(:)*ZESG196(:,JI)*ZESG01(:,JJ))
+    ZAPPROX(:)=sqrt(2.)*sqrt(ZRGI(:))**13*(ZESG169(:,JI)+ZR(:)*ZESG144(:,JI)*ZESG01(:,JJ)+&
+           2.*ZR2(:)*ZESG121(:,JI)*ZESG04(:,JJ)+ZR4(:)*ZESG81(:,JI)*ZESG16(:,JJ)+&
+           ZRM3(:)*ZESG256(:,JI)*ZESG09(:,JJ)+2*ZRM(:)*ZESG196(:,JI)*ZESG01(:,JJ))
       
-      ZINTERFM(:)=ZKFM(:)*ZRES(:)*ZAPPROX(:)
+    ZINTERFM(:)=ZKFM(:)*ZRES(:)*ZAPPROX(:)
 
-      ZAPPROX(:)=(ZRGI(:))**6*(2.*ZESG144(:,JI)+ZAI(:)*ZKNGI(:)*(ZESG100(:,JI)+&
+    ZAPPROX(:)=(ZRGI(:))**6*(2.*ZESG144(:,JI)+ZAI(:)*ZKNGI(:)*(ZESG100(:,JI)+&
        ZR2(:)*ZESG64(:,JI)*ZESG04(:,JJ))+ZAJ(:)*ZKNGJ(:)*(ZESG144(:,JI)*ZESG04(:,JJ)+&
        ZRM2(:)*ZESG196(:,JI)*ZESG16(:,JJ))+ZR2(:)*ZESG100(:,JI)*ZESG04(:,JJ)+&
        ZRM2(:)*ZESG196(:,JI)*ZESG04(:,JJ))
 
-      ZINTERNC(:)=ZKNC(:)*ZAPPROX(:)
+    ZINTERNC(:)=ZKNC(:)*ZAPPROX(:)
 
-      ZINTER(:)=ZINTERNC(:)*(ZINTERFM(:)/(ZINTERNC(:)+ZINTERFM(:)))
+    ZINTER(:)=ZINTERNC(:)*(ZINTERFM(:)/(ZINTERNC(:)+ZINTERFM(:)))
 
-      PDMINTER(:,NM6(JI))=-PM(:,NM0(JI))*PM(:,NM0(JJ))*ZINTER(:)
-      PDMINTER(:,NM6(JJ))=PM(:,NM0(JI))*PM(:,NM0(JJ))*ZINTER(:)
+    PDMINTER(:,NM6(JI))=-PM(:,NM0(JI))*PM(:,NM0(JJ))*ZINTER(:)
 
-      ZAPPROX(:)=sqrt(2.)*sqrt(ZRGI(:))**7*sqrt(ZRGJ(:))**6*(ZESG49(:,JI)*&
+    ZAPPROX(:)=sqrt(2.)*sqrt(ZRGI(:))**7*sqrt(ZRGJ(:))**6*(ZESG49(:,JI)*&
            ZESG36(:,JJ)+ZR(:)*ZESG36(:,JI)*ZESG49(:,JJ)+2.*ZR2(:)*ZESG25(:,JI)*&
            ZESG64(:,JJ)+ZR4(:)*ZESG09(:,JI)*ZESG100(:,JJ)+ZRM3(:)*ZESG100(:,JI)*&
            ZESG09(:,JJ)+2.*ZRM(:)*ZESG64(:,JI)*ZESG25(:,JJ))
        
-      ZINTERFM(:)=ZKFM(:)*ZRES(:)*ZAPPROX(:)
+    ZINTERFM(:)=ZKFM(:)*ZRES(:)*ZAPPROX(:)
 
-      ZAPPROX(:)=(ZRGI(:))**3*(ZRGJ(:))**3*(2.*ZESG36(:,JI)*ZESG36(:,JJ)+&
+    ZAPPROX(:)=(ZRGI(:))**3*(ZRGJ(:))**3*(2.*ZESG36(:,JI)*ZESG36(:,JJ)+&
        ZAI(:)*ZKNGI(:)*(ZESG16(:,JI)*ZESG16(:,JJ)+ZR2(:)*ZESG04(:,JI)*ZESG64(:,JJ))+&
        ZAJ(:)*ZKNGJ(:)*(ZESG36(:,JI)*ZESG16(:,JJ)+ZRM2(:)*ZESG64(:,JI)*ZESG04(:,JJ))+&
        ZR2(:)*ZESG16(:,JI)*ZESG64(:,JJ)+ZRM2(:)*ZESG64(:,JI)*ZESG16(:,JJ))
 
-      ZINTERNC(:)=ZKNC(:)*ZAPPROX(:)
+    ZINTERNC(:)=ZKNC(:)*ZAPPROX(:)
 
-      ZINTER(:)=ZINTERNC(:)*(ZINTERFM(:)/(ZINTERNC(:)+ZINTERFM(:)))
+    ZINTER(:)=ZINTERNC(:)*(ZINTERFM(:)/(ZINTERNC(:)+ZINTERFM(:)))
        
-      PDMINTER(:,NM6(JJ))=PDMINTER(:,NM6(JJ))+2.*PM(:,NM0(JI))*PM(:,NM0(JJ))*ZINTER(:)
+    PDMINTER(:,NM6(JJ))=-PDMINTER(:,NM6(JI))+2.*PM(:,NM0(JI))*PM(:,NM0(JJ))*ZINTER(:)
  
-  enddo
-enddo
+  ENDDO
+ENDDO
 
 END SUBROUTINE CH_AER_COAG 
diff --git a/src/MNH/ch_aer_depos.f90 b/src/MNH/ch_aer_depos.f90
index cf5659aa332d24dec88b471aae1fd4ae678c0ac7..37ca0645fe6e357302967318c408f59e0c2469a3 100644
--- a/src/MNH/ch_aer_depos.f90
+++ b/src/MNH/ch_aer_depos.f90
@@ -116,12 +116,12 @@ END IF
 !Get minimum values possible
 ZPMIN(1) = XN0IMIN
 ZRGMIN = ZINIRADIUSI
-ZPMIN(2) = ZPMIN(1) * (ZRGMIN**3)*EXP(4.5 * LOG(XSIGIMIN)**2) 
-ZPMIN(3) = ZPMIN(1) * (ZRGMIN**6)*EXP(18. * LOG(XSIGIMIN)**2)
+ZPMIN(2) = ZPMIN(1) * (ZRGMIN**3)*EXP(4.5 * LOG(XINISIGI)**2) 
+ZPMIN(3) = ZPMIN(1) * (ZRGMIN**6)*EXP(18. * LOG(XINISIGI)**2)
 ZPMIN(4) = XN0JMIN
 ZRGMIN = ZINIRADIUSJ
-ZPMIN(5) = ZPMIN(4) * (ZRGMIN**3)*EXP(4.5 * LOG(XSIGJMIN)**2) 
-ZPMIN(6) = ZPMIN(4) * (ZRGMIN**6)*EXP(18. * LOG(XSIGJMIN)**2)
+ZPMIN(5) = ZPMIN(4) * (ZRGMIN**3)*EXP(4.5 * LOG(XINISIGJ)**2) 
+ZPMIN(6) = ZPMIN(4) * (ZRGMIN**6)*EXP(18. * LOG(XINISIGJ)**2)
 !
 !
 CALL  PPP2AERO(ZSVT, PRHODREF, &
@@ -167,15 +167,15 @@ SELECT CASE (CCLOUD)
   CASE ('KESS','REVE','ICE3','ICE4')
 ! One moment cloud scheme
   CALL AER_WET_DEP_KMT_WARM  (NSPLITR, PTSTEP, PZZ, PRHODREF,    &
-                              PRT(:,:,:,2), PRT(:,:,:,3), ZRCS,  &
-                              ZRRS, ZSVAER, PTHT, PPABST,  ZRG,  &
+                              PRT(:,:,:,2), PRT(:,:,:,3),        &
+                              ZSVAER, PTHT, PPABST,  ZRG,        &
                               PEVAP3D, JPMODE,ZRHOP, ZMASSMIN,   &
                               PSEA=ZSEA, PTOWN=ZTOWN)
   CASE ('KHKO','C2R2','C3R5')
 ! Two moment cloud scheme
   CALL AER_WET_DEP_KMT_WARM  (NSPLITR, PTSTEP, PZZ, PRHODREF,    &
-                              PRT(:,:,:,2), PRT(:,:,:,3), ZRCS,  &
-                              ZRRS, ZSVAER, PTHT, PPABST,  ZRG,  &
+                              PRT(:,:,:,2), PRT(:,:,:,3),        &
+                              ZSVAER, PTHT, PPABST,  ZRG,        &
                               PEVAP3D, JPMODE,ZRHOP, ZMASSMIN,   &
                               PSEA=ZSEA, PTOWN=ZTOWN,            &
                               PCCT=PSVT(:,:,:,NSV_C2R2BEG+1),    &
@@ -183,8 +183,8 @@ SELECT CASE (CCLOUD)
   CASE ('LIMA')
 ! Two moment cloud scheme
   CALL AER_WET_DEP_KMT_WARM  (NSPLITR, PTSTEP, PZZ, PRHODREF,    &
-                              PRT(:,:,:,2), PRT(:,:,:,3), ZRCS,  &
-                              ZRRS, ZSVAER, PTHT, PPABST,  ZRG,  &
+                              PRT(:,:,:,2), PRT(:,:,:,3),        &
+                              ZSVAER, PTHT, PPABST,  ZRG,        &
                               PEVAP3D, JPMODE,ZRHOP, ZMASSMIN,   &
                               PSEA=ZSEA, PTOWN=ZTOWN,            &
                               PCCT=PSVT(:,:,:,NSV_LIMA_NC),      &
diff --git a/src/MNH/ch_aer_eqm_initn.f90 b/src/MNH/ch_aer_eqm_initn.f90
index 0759d63cebf4fb9ecacb825551b4dc67bb1bec1e..d9b1a77621153cfe89c0d5365b9c2b267761147a 100644
--- a/src/MNH/ch_aer_eqm_initn.f90
+++ b/src/MNH/ch_aer_eqm_initn.f90
@@ -403,16 +403,16 @@ XSVMIN(NSV_AERBEG:NSV_AEREND) = XMNH_TINY
 XSVMIN(NSV_CHEMBEG-1+JP_CH_CO) =  1E-10
 ! For i mode
 ZRHODREFMIN = MAX_ll( PRHODREF(:,:,:), IINFO_ll)
-ZMASS  = XN0IMIN *  ((ZMINRGI**3)*EXP(4.5 * (LOG(XSIGIMIN))**2))
-ZM6MIN = XN0IMIN *  ((ZMINRGI**6)*EXP(18. * (LOG(XSIGIMIN))**2))
+ZMASS  = XN0IMIN *  ((ZMINRGI**3)*EXP(4.5 * (LOG(XINISIGI))**2))
+ZM6MIN = XN0IMIN *  ((ZMINRGI**6)*EXP(18. * (LOG(XINISIGI))**2))
 XSVMIN(NSV_AERBEG-1+JP_CH_BCi) = 0.5*ZMASS * XFAC(JP_AER_BC) * 6.0221367E+11/(ZDEN2MOL*12.*ZRHODREFMIN)
 XSVMIN(NSV_AERBEG-1+JP_CH_OCi) = 0.5*ZMASS * XFAC(JP_AER_OC) * 6.0221367E+11/(ZDEN2MOL*12.*ZRHODREFMIN)
 XSVMIN(NSV_AERBEG-1+JP_CH_M0i) = XN0IMIN * 1E-6 / (ZDEN2MOL*ZRHODREFMIN)
 IF (LVARSIGI) XSVMIN(NSV_AERBEG-1+JP_CH_M6i) = ZM6MIN  / (ZDEN2MOL*ZRHODREFMIN)
 !
 ! For j mode
-ZMASS  = XN0JMIN *  ((ZMINRGJ**3)*EXP(4.5 * (LOG(XSIGJMIN))**2))
-ZM6MIN = XN0JMIN *  ((ZMINRGJ**6)*EXP(18. * (LOG(XSIGJMIN))**2))
+ZMASS  = XN0JMIN *  ((ZMINRGJ**3)*EXP(4.5 * (LOG(XINISIGJ))**2))
+ZM6MIN = XN0JMIN *  ((ZMINRGJ**6)*EXP(18. * (LOG(XINISIGJ))**2))
 XSVMIN(NSV_AERBEG-1+JP_CH_BCj) = 0.5*ZMASS * XFAC(JP_AER_BC) * 6.0221367E+11/(ZDEN2MOL*12.*ZRHODREFMIN)
 XSVMIN(NSV_AERBEG-1+JP_CH_OCj) = 0.5*ZMASS * XFAC(JP_AER_OC) * 6.0221367E+11/(ZDEN2MOL*12.*ZRHODREFMIN)
 XSVMIN(NSV_AERBEG-1+JP_CH_M0j) = XN0JMIN * 1E-6 / (ZDEN2MOL*ZRHODREFMIN)
diff --git a/src/MNH/ch_aer_reallfin.f90 b/src/MNH/ch_aer_reallfin.f90
index 505f3acfae89b18fac882c45bc4a2fe83b953ac0..101835c6c8a9ce8afeef7ffa24f52c4a1aaa3b8b 100644
--- a/src/MNH/ch_aer_reallfin.f90
+++ b/src/MNH/ch_aer_reallfin.f90
@@ -134,7 +134,7 @@ ZRHOI(:) = 1.8e3
 ZRHOI(JP_AER_H2O) = 1.0e3   ! water
 ZRHOI(JP_AER_DST) = XDENSITY_DUST
 !
-PSV(:,:,:,:)  = MAX(PSV(:,:,:,:), 0.)
+PSV(:,:,:,:)  = MAX(PSV(:,:,:,:), 1E-20)
 PCO(:,:,:)    = MAX(PCO(:,:,:), 10.E-9*PRHODREF(:,:,:))
 !
 ! Special treatment for BC and OC (link to CO gaseous concentration)
@@ -310,6 +310,7 @@ DO JN=1,JPMODE
 !
   ZRG(:,:,:,JN)=(ZM(:,:,:,NM3(JN))/ZN(:,:,:,JN))**(1./3.)*EXP(-1.5*ZSIGMA(:,:,:))
 !
+
 ENDDO
 !
 !conversion into ppv
diff --git a/src/MNH/ch_aer_trans.f90 b/src/MNH/ch_aer_trans.f90
index d8844d117baa87b4a65456a4af08737169242bb2..400ccec1f6e4b79c458fe1a55b8502cb8d3180f7 100644
--- a/src/MNH/ch_aer_trans.f90
+++ b/src/MNH/ch_aer_trans.f90
@@ -95,16 +95,16 @@ CHARACTER(LEN=10),        INTENT(IN) :: HSCHEME
 !
 !*      0.2    declarations local variables
 !
-INTEGER :: JJ, JN  ! loop counter
+INTEGER :: JJ, JN, II  ! loop counter
 !   variables for the aerosol module
 !
 REAL, DIMENSION(SIZE(PM,1)) :: ZSIGMA
 REAL, DIMENSION(SIZE(PM,1)) :: ZSUM
 REAL, SAVE, DIMENSION(JPIN) :: ZPMIN
 LOGICAL, SAVE               :: GPHYSLIM = .TRUE. ! flag
-REAL :: ZRGMIN
+REAL, SAVE    :: ZRGMIN
 REAL, PARAMETER :: ZCSTAVOG=6.0221367E+11
-REAL    :: ZINIRADIUSI, ZINIRADIUSJ
+REAL, SAVE    :: ZINIRADIUSI, ZINIRADIUSJ
 
 
 !
@@ -128,13 +128,13 @@ END IF
 
 ZPMIN(1) = XN0IMIN
 ZRGMIN = ZINIRADIUSI
-ZPMIN(2) = ZPMIN(1) * (ZRGMIN**3)*EXP(4.5 * LOG(XSIGIMIN)**2) 
-ZPMIN(3) = ZPMIN(1) * (ZRGMIN**6)*EXP(18. * LOG(XSIGIMIN)**2)
+ZPMIN(2) = ZPMIN(1) * (ZRGMIN**3)*EXP(4.5 * LOG(XINISIGI)**2) 
+ZPMIN(3) = ZPMIN(1) * (ZRGMIN**6)*EXP(18. * LOG(XINISIGI)**2)
 
 ZPMIN(4) = XN0JMIN
 ZRGMIN = ZINIRADIUSJ
-ZPMIN(5) = ZPMIN(4) * (ZRGMIN**3)*EXP(4.5 * LOG(XSIGJMIN)**2) 
-ZPMIN(6) = ZPMIN(4) * (ZRGMIN**6)*EXP(18. * LOG(XSIGJMIN)**2)
+ZPMIN(5) = ZPMIN(4) * (ZRGMIN**3)*EXP(4.5 * LOG(XINISIGJ)**2) 
+ZPMIN(6) = ZPMIN(4) * (ZRGMIN**6)*EXP(18. * LOG(XINISIGJ)**2)
 
 END IF
 
@@ -768,6 +768,8 @@ END IF
 !
 !*       1.n    transfer moment 0 from gas to aerosol variable
 !
+!print*,'aer_trans N0i =',MINVAL(PAERO(:,JP_CH_M0i)), MAXVAL(PAERO(:,JP_CH_M0i))
+!print*,'aer_trans N0j =',MINVAL(PAERO(:,JP_CH_M0j)), MAXVAL(PAERO(:,JP_CH_M0j))
   PM(:,1) = MAX(PAERO(:,JP_CH_M0i) * 1E+6, XMNH_TINY) ! convert from 1/cc to 1/m3
   PM(:,4) = MAX(PAERO(:,JP_CH_M0j) * 1E+6, XMNH_TINY) ! convert from 1/cc to 1/m3
 
@@ -892,7 +894,7 @@ PMASK(:,JN) = 1.
 END WHERE
 
 ENDDO
-
+  
 ELSE
 !
 !*       2.     TRANSFER FROM AEROSOL TO GAS  MODULE
diff --git a/src/MNH/ch_aer_wetdepn.f90 b/src/MNH/ch_aer_wetdepn.f90
index b1041d55d5171b6a39b85289f99ee1b09a47536b..dd051acd8e62e7409f874285fe408355fa6c0793 100644
--- a/src/MNH/ch_aer_wetdepn.f90
+++ b/src/MNH/ch_aer_wetdepn.f90
@@ -99,12 +99,12 @@ END IF
 !Get minimum values possible
 ZPMIN(1) = XN0IMIN
 ZRGMIN =  ZINIRADIUSI
-ZPMIN(2) = ZPMIN(1) * (ZRGMIN**3)*EXP(4.5 * LOG(XSIGIMIN)**2)
-ZPMIN(3) = ZPMIN(1) * (ZRGMIN**6)*EXP(18. * LOG(XSIGIMIN)**2)
+ZPMIN(2) = ZPMIN(1) * (ZRGMIN**3)*EXP(4.5 * LOG(XINISIGI)**2)
+ZPMIN(3) = ZPMIN(1) * (ZRGMIN**6)*EXP(18. * LOG(XINISIGI)**2)
 ZPMIN(4) = XN0JMIN
 ZRGMIN =  ZINIRADIUSJ
-ZPMIN(5) = ZPMIN(4) * (ZRGMIN**3)*EXP(4.5 * LOG(XSIGJMIN)**2)
-ZPMIN(6) = ZPMIN(4) * (ZRGMIN**6)*EXP(18. * LOG(XSIGJMIN)**2)
+ZPMIN(5) = ZPMIN(4) * (ZRGMIN**3)*EXP(4.5 * LOG(XINISIGJ)**2)
+ZPMIN(6) = ZPMIN(4) * (ZRGMIN**6)*EXP(18. * LOG(XINISIGJ)**2)
 !
 CALL  PPP2AERO(PSVT, PRHODREF, PSIG3D=ZSIG, PRG3D=ZRG, PM3D=ZPMOLD)
 CALL  PPP2AERO(PCWETDEP, PRHODREF, PSIG3D=ZSIGN, PRG3D=ZRGN, PM3D=ZPM)
diff --git a/src/MNH/ch_monitorn.f90 b/src/MNH/ch_monitorn.f90
index b61f19288e3472fd3542f6188ba7cd5b26472f24..e4eb5b0d7c4c3ccc9d4630bd3af39a20d7a5e612 100644
--- a/src/MNH/ch_monitorn.f90
+++ b/src/MNH/ch_monitorn.f90
@@ -392,7 +392,6 @@ REAL, DIMENSION(:,:), ALLOCATABLE :: ZRHOP, ZSOLORG
 REAL, DIMENSION(:),   ALLOCATABLE :: ZSO4RAT
 REAL, DIMENSION(:),   ALLOCATABLE :: ZJNUC, ZJ2RAT
 
-REAL,DIMENSION(SIZE(XSVT,1),SIZE(XSVT,2),SIZE(XSVT,3),SIZE(XSVT,4)) :: ZSVT
 REAL,DIMENSION(SIZE(XSVT,1),SIZE(XSVT,2),SIZE(XSVT,3),NSV_AER) :: ZCWETAERO
 !
 !-------------------------------------------------------------------------------
@@ -609,25 +608,6 @@ ENDIF  ! first time step
 !
 ZDTSOLVER = PTSTEP / NCH_SUBSTEPS
 !
-!*       1.3   give minimum value and conserve mass for aerosols
-!
-!
-IF (LORILAM) THEN
-
-  IF (CPROGRAM/='DIAG  ') THEN
-   DO JSV = 1, SIZE(XSVT,4)
-    ZSVT(:,:,:,JSV) =  XRSVS(:,:,:,JSV) *PTSTEP / XRHODJ(:,:,:) 
-   END DO
-  ELSE
-   DO JSV = 1, SIZE(XSVT,4)
-    ZSVT(:,:,:,JSV) =  XSVT(:,:,:,JSV)
-   END DO
-  END IF
-  ZSVT(:,:,:,NSV_CHEMBEG:NSV_CHEMEND) = MAX(ZSVT(:,:,:,NSV_CHEMBEG:NSV_CHEMEND), XMNH_TINY)
-  ZSVT(:,:,:,NSV_AERBEG:NSV_AEREND)   = MAX(ZSVT(:,:,:,NSV_AERBEG:NSV_AEREND), XMNH_TINY)
-!
-END IF
-!
 !*       1.4 compute conversion factor ppv/m3 --> molec/cm3
 !
 ZDEN2MOL = 1E-6 * XAVOGADRO / XMD
@@ -774,10 +754,10 @@ IF (LORILAM) THEN
   !
   IF ((LSEDIMAERO).AND.(CPROGRAM/='DIAG  ')) THEN
     CALL CH_AER_SEDIM_n(PTSTEP,                                                &
-            ZSVT(IIB:IIE,IJB:IJE,IKB:IKE,NSV_AERBEG:NSV_AEREND),               &
+            XSVT(IIB:IIE,IJB:IJE,IKB:IKE,NSV_AERBEG:NSV_AEREND),               &
             XTHT(IIB:IIE,IJB:IJE,IKB:IKE),  XRHODREF(IIB:IIE,IJB:IJE,IKB:IKE), &
-          XPABST(IIB:IIE,IJB:IJE,IKB:IKE), XVDEPAERO(IIB:IIE,IJB:IJE,:),       &
-             XZZ(IIB:IIE,IJB:IJE,IKB:IKE),     XSEDA(IIB:IIE,IJB:IJE,IKB:IKE,:))
+            XPABST(IIB:IIE,IJB:IJE,IKB:IKE), XVDEPAERO(IIB:IIE,IJB:IJE,:),     &
+            XZZ(IIB:IIE,IJB:IJE,IKB:IKE),     XSEDA(IIB:IIE,IJB:IJE,IKB:IKE,:))
   ENDIF
 ! implicit  wet deposition
   IF ((LCH_CONV_SCAV).AND.(CPROGRAM/='DIAG  ')) THEN
@@ -786,7 +766,7 @@ IF (LORILAM) THEN
     END DO
     ZCWETAERO(:,:,:,:)= MAX(ZCWETAERO(:,:,:,:), XMNH_TINY)
 
-    CALL CH_AER_WETDEP_n(PTSTEP, ZSVT(IIB:IIE,IJB:IJE,IKB:IKE,NSV_AERBEG:NSV_AEREND),             &
+    CALL CH_AER_WETDEP_n(PTSTEP, XSVT(IIB:IIE,IJB:IJE,IKB:IKE,NSV_AERBEG:NSV_AEREND),             &
                          ZCWETAERO(IIB:IIE,IJB:IJE,IKB:IKE,:), XRHODREF(IIB:IIE,IJB:IJE,IKB:IKE), &
                          XSEDA(IIB:IIE,IJB:IJE,IKB:IKE,:))
   ENDIF
@@ -798,7 +778,7 @@ IF (LORILAM) THEN
                         XRT(IIB:IIE,IJB:IJE,IKB:IKE,:),     &
                         XRRS(IIB:IIE,IJB:IJE,IKB:IKE,:),    &
                         XRHODJ(IIB:IIE,IJB:IJE,IKB:IKE),    &
-                        ZSVT(IIB:IIE,IJB:IJE,IKB:IKE,:),    &
+                        XSVT(IIB:IIE,IJB:IJE,IKB:IKE,:),    &
                         XMI(IIB:IIE,IJB:IJE,IKB:IKE,:),     &
                         XTHT(IIB:IIE,IJB:IJE,IKB:IKE),      &
                         XPABST(IIB:IIE,IJB:IJE,IKB:IKE),    &
@@ -808,7 +788,7 @@ IF (LORILAM) THEN
   ENDIF
 ! Update aerosol tendency before aerosol solver
   DO JSV = 1, SIZE(XSVT,4)
-    XRSVS(:,:,:,JSV) = ZSVT(:,:,:,JSV) * XRHODJ(:,:,:) / PTSTEP
+    XRSVS(:,:,:,JSV) = XSVT(:,:,:,JSV) * XRHODJ(:,:,:) / PTSTEP
   END DO
 ENDIF
 !
diff --git a/src/MNH/default_desfmn.f90 b/src/MNH/default_desfmn.f90
index 6678d6e268ccc85266d79418972f7814fb54923f..22911ad2b5bee50a891e215ac87a3d700dcfd4d4 100644
--- a/src/MNH/default_desfmn.f90
+++ b/src/MNH/default_desfmn.f90
@@ -613,7 +613,7 @@ XLAT_PROF(:)  = XUNDEF
 XLON_PROF(:)  = XUNDEF
 CNAME_PROF(:) = ''
 CFILE_PROF    = 'NO_INPUT_CSV'
-LDIAG_SURFRAD_PROF = .TRUE.
+! LDIAG_SURFRAD = .TRUE.
 !------------------------------------------------------------------------------
 !*      10.f   SET DEFAULT VALUES FOR MODD_ALLSTATION_n :
 !             ----------------------------------
@@ -627,7 +627,7 @@ XLAT_STAT(:)  = XUNDEF
 XLON_STAT(:)  = XUNDEF
 CNAME_STAT(:) = ''
 CFILE_STAT    = 'NO_INPUT_CSV'
-LDIAG_SURFRAD_STAT = .TRUE.
+LDIAG_SURFRAD = .TRUE.
 !
 !-------------------------------------------------------------------------------
 !
@@ -1088,7 +1088,7 @@ LCH_PH              = .FALSE.
 LCH_RET_ICE         = .FALSE.
 XCH_PHINIT          = 5.2
 XRTMIN_AQ           = 5.e-8
-CCHEM_INPUT_FILE    = 'EXSEG1.nam'
+CCHEM_INPUT_FILE    = 'MNHC.input'
 CCH_TDISCRETIZATION = 'SPLIT'
 NCH_SUBSTEPS        = 1
 LCH_TUV_ONLINE      = .FALSE.
@@ -1168,8 +1168,8 @@ LHETEROSO4 = .FALSE.  ! switch to active sulfates heteronegeous
                       ! production
 LSEDIMAERO = .FALSE.  ! switch to active aerosol sedimentation
 LAERINIT   = .FALSE.  ! switch to initialize aerosol in arome
-CMINERAL      = "NONE"   ! mineral equilibrium scheme
-CORGANIC      = "NONE"   ! mineral equilibrium scheme
+CMINERAL      = "EQSAM"   ! mineral equilibrium scheme
+CORGANIC      = "MPMPO"   ! mineral equilibrium scheme
 CNUCLEATION   = "NONE" ! sulfates nucleation scheme
 LDEPOS_AER(:) = .FALSE.
 
diff --git a/src/MNH/dustcamsn.f90 b/src/MNH/dustcamsn.f90
index 33966adfd96f0d7b3e6d5684f374c77cc7721e38..94732c05c0006a92918710c24f4e1a329ebb660d 100644
--- a/src/MNH/dustcamsn.f90
+++ b/src/MNH/dustcamsn.f90
@@ -125,13 +125,12 @@ DO JN = 1, NMODE_DST
   ZMMIN(IM3(JN)) = XN0MIN(IMODEIDX) * (ZRGMIN**3)*EXP(4.5 * LOG(ZINISIGMA(JN))**2) 
   ZMMIN(IM6(JN)) = XN0MIN(IMODEIDX) * (ZRGMIN**6)*EXP(18. * LOG(ZINISIGMA(JN))**2)
 
-  IF (JPDUSTORDER(JN) == 1) ZMASS(:,:,:,JN) = PMASSCAMS(:,:,:,1) ! fin mode 
-  IF (JPDUSTORDER(JN) == 2) ZMASS(:,:,:,JN) = PMASSCAMS(:,:,:,2) ! median mode 
-  IF (JPDUSTORDER(JN) == 3) ZMASS(:,:,:,JN) = PMASSCAMS(:,:,:,3) ! large mode
+  IF (JPDUSTORDER(JN) == 1) ZMASS(:,:,:,JN) = MAX(PMASSCAMS(:,:,:,1), 1E-16) ! fin mode 
+  IF (JPDUSTORDER(JN) == 2) ZMASS(:,:,:,JN) = MAX(PMASSCAMS(:,:,:,2), 1E-15) ! median mode 
+  IF (JPDUSTORDER(JN) == 3) ZMASS(:,:,:,JN) = MAX(PMASSCAMS(:,:,:,3), 1E-15) ! large mode
 
 ENDDO
 
-ZMASS(:,:,:,:) = MAX(ZMASS(:,:,:,:), 1E-40)
 !
 !
 ZRHOI = XDENSITY_DUST !1.8e3 !++changed alfgr
@@ -158,14 +157,14 @@ DO JN = 1, NMODE_DST
                       (ZINIRADIUS(JN)**3) * &
                       EXP(4.5*LOG(ZINISIGMA(JN))**2) 
 
-  ZM(:,:,:,IM3(JN)) = MAX(ZMMIN(IM3(JN)), ZM(:,:,:,IM3(JN)))
+!  ZM(:,:,:,IM3(JN)) = MAX(ZMMIN(IM3(JN)), ZM(:,:,:,IM3(JN)))
 !
 !*       1.3    calculate moment 6 from m0,  RG and SIG 
 !
   ZM(:,:,:,IM6(JN))= ZM(:,:,:,IM0(JN)) * ((ZINIRADIUS(JN)**6) * &
                      EXP(18.*(LOG(ZINISIGMA(JN)))**2))
 !
-  ZM(:,:,:,IM6(JN)) = MAX(ZMMIN(IM6(JN)), ZM(:,:,:,IM6(JN)))
+!  ZM(:,:,:,IM6(JN)) = MAX(ZMMIN(IM6(JN)), ZM(:,:,:,IM6(JN)))
 !
 !*       1.4    output concentration
 !
diff --git a/src/MNH/init_salt.f90 b/src/MNH/init_salt.f90
index ab14998b35e89dc56ce1e4b224b74877ce471401..bb99ab50f8e781aa38ef3af7ddf7cb790406ad9b 100644
--- a/src/MNH/init_salt.f90
+++ b/src/MNH/init_salt.f90
@@ -45,7 +45,7 @@ XINIRADIUS_SLT=  (/0.009, 0.021, 0.045, 0.115, 0.415,2.5, 7.0, 25.0/)
 !Initial, standard deviation from  Ova et al., 2014
 XINISIG_SLT =  (/ 1.37, 1.5, 1.42, 1.53, 1.85,1.7, 1.8, 2.1 /)
 !Minimum allowed number concentration for any mode (#/m3)
-XN0MIN_SLT  = (/1.e1 , 1.e1, 1.e1, 1., 1.e-4,1.e-20 , 1.e-20, 1.e-20 /)
+XN0MIN_SLT  = (/1.e1 , 1.e1, 1.e1, 1., 1.e-4, 1.e-5, 1.e-6, 1.e-7 /)
 
 ELSE IF ( NMODE_SLT == 3) THEN
 
diff --git a/src/MNH/modd_dust.f90 b/src/MNH/modd_dust.f90
index 540de108e225683892178339a0ca054f7064f0c2..2f2a7dca2de875bddd3dd308f7aa875cc26c8d21 100644
--- a/src/MNH/modd_dust.f90
+++ b/src/MNH/modd_dust.f90
@@ -100,7 +100,7 @@ REAL, DIMENSION(3)          :: XINIRADIUS= 0.5*(/0.078, 0.641, 5.00 /)
 !Initial, standard deviation from Alfaro et al 1998
 REAL, DIMENSION(3)          :: XINISIG =  (/1.75, 1.76, 1.70/)
 !Minimum allowed number concentration for any mode (#/m3)
-REAL, DIMENSION(3)          :: XN0MIN  = (/1.e3 , 1.e1 , 1.e-2 /)
+REAL, DIMENSION(3)          :: XN0MIN  = (/1.e1 , 1.e-1 , 1.e-4 /)
 CHARACTER(LEN=9),DIMENSION(:),ALLOCATABLE :: CDEDSTNAMES
 CHARACTER(LEN=9),DIMENSION(6), PARAMETER  :: YPDEDST_INI = &
      (/'DEDSTM31C','DEDSTM32C','DEDSTM33C' &
diff --git a/src/MNH/modd_salt.f90 b/src/MNH/modd_salt.f90
index e111b15db085287e012c5afaf137f9be3bd03185..c8ebbebd09a3525e588a3ce3b226f22ac60aaa44 100644
--- a/src/MNH/modd_salt.f90
+++ b/src/MNH/modd_salt.f90
@@ -83,7 +83,7 @@ REAL,DIMENSION(8)    :: XINIRADIUS_SLT=  (/0.009, 0.021, 0.045, 0.115,0.415,2.5,
 REAL,DIMENSION(8)      :: XINISIG_SLT =  (/ 1.37, 1.5, 1.42, 1.53, 1.85,1.7,1.8, 2.9 /)
 
 !Minimum allowed number concentration for any mode (#/m3)
-REAL,DIMENSION(8)  :: XN0MIN_SLT  = (/1.e1 , 1.e1, 1.e1, 1., 1.e-4,1.e-20, 1.e-20,1.e-20 /)
+REAL,DIMENSION(8)  :: XN0MIN_SLT  = (/1.e1 , 1.e1, 1.e1, 1., 1.e-4,1.e-5, 1.e-6,1.e-7 /)
 !Test Thomas
 REAL, DIMENSION(:,:,:), ALLOCATABLE :: XSLTMSS   ! [kg/m3] total mass concentration of sea salt
 !
diff --git a/src/MNH/mode_aero_psd.f90 b/src/MNH/mode_aero_psd.f90
index d5d72fd29799702f2504a7badf56c22b208b1800..f03bdfcee1b9b84623a350d908857a96eff8ec4d 100644
--- a/src/MNH/mode_aero_psd.f90
+++ b/src/MNH/mode_aero_psd.f90
@@ -101,6 +101,7 @@ REAL                                 :: ZDEN2MOL
 REAL,DIMENSION(JPMODE*3)             :: ZPMIN               ! [aerosol units] minimum values for N, sigma, M
 INTEGER                              :: JI,JJ,JK,JSV, JN              ! [idx] loop counters
 REAL    :: ZINIRADIUSI, ZINIRADIUSJ
+INTEGER :: II,IJ,IK
 !
 !-------------------------------------------------------------------------------
 !
@@ -253,6 +254,7 @@ END IF
 !
 !*       5    set moment 6  ==> um6_{aer}/m3_{air}
 !
+
 IF (LVARSIGI) THEN ! set M6 variable standard deviation
   IF ((CPROGRAM=="REAL  ").OR.(CPROGRAM=="IDEAL ")) THEN
   ZM(:,:,:,3)= ZM(:,:,:,1) * (ZINIRADIUSJ**6)*EXP(18. * LOG(XINISIGJ)**2)
@@ -328,6 +330,7 @@ DO JN=1,JPMODE
 
 ENDDO
 !
+
 IF(PRESENT(PCTOTA)) PCTOTA(:,:,:,:,:) = ZCTOTA(:,:,:,:,:)
 IF(PRESENT(PM3D)) PM3D(:,:,:,:) = ZM(:,:,:,:)
 !
@@ -463,7 +466,9 @@ ENDDO
   ZCTOTA(:,:,:,JP_AER_SOA10,1) = PSVT(:,:,:,JP_CH_SOA10i)
   ZCTOTA(:,:,:,JP_AER_SOA10,2) = PSVT(:,:,:,JP_CH_SOA10j)
   END IF
-  ZCTOTA(:,:,:,:,:) = MAX(ZCTOTA(:,:,:,:,:),XMNH_TINY)
+  ZCTOTA(:,:,:,:,1) = MAX(ZCTOTA(:,:,:,:,1),1E-6)
+  ZCTOTA(:,:,:,:,2) = MAX(ZCTOTA(:,:,:,:,2),1E-4)
+
 
 !
 !*       3    calculate moment 3 from total aerosol mass
@@ -547,8 +552,6 @@ IF (LVARSIGJ) PSVT(:,:,:,JP_CH_M6j) = ZM(:,:,:,6)
 DO JJ=1,SIZE(PSVT,4)
   PSVT(:,:,:,JJ) =  PSVT(:,:,:,JJ) / (ZDEN2MOL * PRHODREF(:,:,:))
 ENDDO
-
-
 END SUBROUTINE CON2MIX
 
 !   ############################################################
diff --git a/src/MNH/read_chem_data_netcdf_case.f90 b/src/MNH/read_chem_data_netcdf_case.f90
index 2648709e0c93aed6ecd7ca7c21a4727a94d9881b..f55fe1f28bfeeac1d67c04d6975e47be7ccce3df 100644
--- a/src/MNH/read_chem_data_netcdf_case.f90
+++ b/src/MNH/read_chem_data_netcdf_case.f90
@@ -175,7 +175,7 @@ INTEGER(kind=CDFINT)               :: ind_netcdf    ! Indice for netcdf var.
 INTEGER                                       :: ICHANNEL
 CHARACTER(LEN=8)                              :: YMOZ="MOZ1.nam"
 integer                                       :: IMOZ
-CHARACTER(LEN=68)                             :: YFORMAT
+CHARACTER(LEN=100)                            :: YFORMAT
 CHARACTER(LEN=40), DIMENSION(:), ALLOCATABLE  :: YSPCMNH
 integer, dimension(:), ALLOCATABLE            :: ISPCMOZ
 CHARACTER(LEN=9)                              :: YA
diff --git a/src/MNH/saltcamsn.f90 b/src/MNH/saltcamsn.f90
index 1747bddc7934c3ada4267f7df353af05de3ab6e1..96dcde95d634c1a05a2ce9a284648bfd704ecc4c 100644
--- a/src/MNH/saltcamsn.f90
+++ b/src/MNH/saltcamsn.f90
@@ -135,6 +135,8 @@ DELTA_CAMS_3 = RAY_CAMS_4 - RAY_CAMS_3
 ! puis calcul de la masse correspondante avec facteur correctif pour eviter
 ! la surestimation des concentrations en aerosols
 
+ZMASS(:,:,:,1) = PMASSCAMS(:,:,:,1) * 1E-3
+
 DELTA_1 = RAY_2 - RAY_CAMS_1
 RATIO_1 = DELTA_1 / DELTA_CAMS_1
 ZMASS(:,:,:,2) = PMASSCAMS(:,:,:,1) * RATIO_1 ! * 1E-2 ! Attribution Mode 2 ORILAM
@@ -159,18 +161,16 @@ DELTA_6 = 10 - RAY_CAMS_3
 RATIO_6 = DELTA_3 / DELTA_CAMS_1
 ZMASS(:,:,:,5) = (PMASSCAMS(:,:,:,3) * RATIO_6) + ZMASS(:,:,:,5) ! Attribution Mode 5 ter ORILAM
 
-ZMASS(:,:,:,5) = ZMASS(:,:,:,5) * 1E-1
-
-! Hyp : the ultrafine mode is neglected for orilam-lima coupling
-ZMASS(:,:,:,1) = PMASSCAMS(:,:,:,1) * 1E-5 ! ultrafin mode
-!
 !========================================================
 ! Adjust the mass / SSA emissions after a few hours
-ZMASS(:,:,:,1) = ZMASS(:,:,:,1) * 1.
-ZMASS(:,:,:,2) = ZMASS(:,:,:,2) * 1.
-ZMASS(:,:,:,3) = ZMASS(:,:,:,3) * 1.
-ZMASS(:,:,:,4) = ZMASS(:,:,:,4) * 1.
-ZMASS(:,:,:,5) = ZMASS(:,:,:,5) * 1.
+ZMASS(:,:,:,1) = MAX(ZMASS(:,:,:,1) * 0.16, 1E-18)
+ZMASS(:,:,:,2) = MAX(ZMASS(:,:,:,2) * 0.1, 1E-17)
+ZMASS(:,:,:,3) = MAX(ZMASS(:,:,:,3) * 0.5, 1E-16)
+ZMASS(:,:,:,4) = MAX(ZMASS(:,:,:,4) * 0.1, 1E-16)
+ZMASS(:,:,:,5) = MAX(ZMASS(:,:,:,5), 1E-17) 
+IF (NMODE_SLT >= 6) ZMASS(:,:,:,6) = MAX(ZMASS(:,:,:,5) * 0.01, 1E-16)
+IF (NMODE_SLT >= 7) ZMASS(:,:,:,7) = MAX(ZMASS(:,:,:,5) * 0.001, 1E-16)
+IF (NMODE_SLT >= 8) ZMASS(:,:,:,8) = MAX(ZMASS(:,:,:,5) * 0.0001, 1E-16)
 !========================================================
 
 DO JN = 1, NMODE_SLT
@@ -197,15 +197,11 @@ DO JN = 1, NMODE_SLT
   ZMMIN(IM6(JN)) = XN0MIN_SLT(IMODEIDX) * (ZRGMIN**6)*EXP(18. * LOG(ZINISIGMA(JN))**2)
 
 END DO
-
-ZMASS(:,:,:,:) = MAX(ZMASS(:,:,:,:), 1E-40)
-!
 !
 ZRHOI = XDENSITY_SALT
 ZMI   = XMOLARWEIGHT_SALT
 ZDEN2MOL = 1E-6 * XAVOGADRO / XMD
 ZFAC = (4. / 3.) * XPI * ZRHOI * 1.e-9
-
 !
 DO JN = 1, NMODE_SLT
 
@@ -225,13 +221,13 @@ DO JN = 1, NMODE_SLT
                       (ZINIRADIUS(JN)**3) * & 
                       EXP(4.5*LOG(ZINISIGMA(JN))**2) 
 
-  ZM(:,:,:,IM3(JN)) = MAX(ZMMIN(IM3(JN)), ZM(:,:,:,IM3(JN)))
+!  ZM(:,:,:,IM3(JN)) = MAX(ZMMIN(IM3(JN)), ZM(:,:,:,IM3(JN)))
 !
 !*       1.3    calculate moment 6 from m0,  RG and SIG 
 !
   ZM(:,:,:,IM6(JN))= ZM(:,:,:,IM0(JN)) * ((ZINIRADIUS(JN)**6)*&
                         EXP(18. * (LOG(ZINISIGMA(JN)))**2))
-  ZM(:,:,:,IM6(JN)) = MAX(ZMMIN(IM6(JN)), ZM(:,:,:,IM6(JN)))
+!  ZM(:,:,:,IM6(JN)) = MAX(ZMMIN(IM6(JN)), ZM(:,:,:,IM6(JN)))
 !
 !*       1.4    output concentration (in ppv)
 !
diff --git a/src/MNH/write_lfifm1_for_diag.f90 b/src/MNH/write_lfifm1_for_diag.f90
index a6099e6a0f4eb779347699adbcf1e6f85fc896ca..01c1c560e957fed8f28c7445895ae399ea2487da 100644
--- a/src/MNH/write_lfifm1_for_diag.f90
+++ b/src/MNH/write_lfifm1_for_diag.f90
@@ -1348,8 +1348,26 @@ IF ((LCHEMDIAG).AND.(LORILAM).AND.(LUSECHEM)) THEN
   IF (.NOT.(ASSOCIATED(XSIG3D))) &
     ALLOCATE(XSIG3D(SIZE(XSVT,1),SIZE(XSVT,2),SIZE(XSVT,3),JPMODE))
   !
-  CALL  PPP2AERO(XSVT(:,:,:,NSV_AERBEG:NSV_AEREND), XRHODREF, &
-                 PSIG3D=XSIG3D, PRG3D=XRG3D, PN3D=XN3D, PCTOTA=ZPTOTA)
+  IF (CRGUNIT=="MASS") THEN
+  XRG3D(:,:,:,1) = XINIRADIUSI * EXP(-3.*(LOG(XINISIGI))**2)
+  XRG3D(:,:,:,2) = XINIRADIUSJ * EXP(-3.*(LOG(XINISIGJ))**2)
+  ELSE
+  XRG3D(:,:,:,1) = XINIRADIUSI
+  XRG3D(:,:,:,2) = XINIRADIUSJ
+  END IF
+  XSIG3D(:,:,:,1) = XINISIGI
+  XSIG3D(:,:,:,2) = XINISIGJ
+  XN3D(:,:,:,1) = XN0IMIN
+  XN3D(:,:,:,2) = XN0JMIN
+  
+  ZPTOTA(:,:,:,:,:) = 0.
+
+  CALL  PPP2AERO(XSVT(IIB:IIE,IJB:IJE,IKB:IKE,NSV_AERBEG:NSV_AEREND),&
+                 XRHODREF(IIB:IIE,IJB:IJE,IKB:IKE), &
+                 PSIG3D=XSIG3D(IIB:IIE,IJB:IJE,IKB:IKE,:),&
+                 PRG3D=XRG3D(IIB:IIE,IJB:IJE,IKB:IKE,:),&
+                 PN3D=XN3D(IIB:IIE,IJB:IJE,IKB:IKE,:),& 
+                 PCTOTA=ZPTOTA(IIB:IIE,IJB:IJE,IKB:IKE,:,:))
 
   TZFIELD = TFIELDMETADATA(                   &
     CMNHNAME   = 'generic for aerosol modes', &
diff --git a/src/SURFEX/coupling_sltn.F90 b/src/SURFEX/coupling_sltn.F90
index 3db0c7b1a987464590ce42376196f115d2188399..cdb6c10ee14606387b1c2ef4ade7427feb3e2155 100644
--- a/src/SURFEX/coupling_sltn.F90
+++ b/src/SURFEX/coupling_sltn.F90
@@ -192,7 +192,7 @@ IF ((CEMISPARAM_SLT .eq. "Ova14").OR.(CEMISPARAM_SLT .eq. "OvB21a").OR.(CEMISPAR
 
   ! Tableau d'interpolation pour calculer ZNUWATER en fonction de la SST
   ! Cas ou 0 < SST < 10 C
-  WHERE ((PSST(:) >= 273.15).AND.(PSST(:) < 283.15))
+  WHERE (PSST(:) < 283.15)
     ZVISCO(:) = ZNUWATER(1) + (PSST(:) - ZWT(1)) * (ZNUWATER(2)-ZNUWATER(1)) / &
                 (ZWT(2) - ZWT(1))
   ENDWHERE
@@ -210,7 +210,7 @@ IF ((CEMISPARAM_SLT .eq. "Ova14").OR.(CEMISPARAM_SLT .eq. "OvB21a").OR.(CEMISPAR
   ENDWHERE
 
   ! Cas ou 30 < SST < 40 C
-  WHERE ((PSST(:) >= 303.15).AND.(PSST(:) < 313.15))
+  WHERE (PSST(:) >= 303.15)
     ZVISCO(:) = ZNUWATER(4) + (PSST(:) - ZWT(4)) * (ZNUWATER(5)-ZNUWATER(4)) / &
                 (ZWT(5) - ZWT(4))
   ENDWHERE
diff --git a/src/SURFEX/default_slt.F90 b/src/SURFEX/default_slt.F90
index ca96f03f98d0b5c56c5cbd4b31838ac4d8c27045..1ac9363d66c36f06f83a25449d1b2d1467410d60 100644
--- a/src/SURFEX/default_slt.F90
+++ b/src/SURFEX/default_slt.F90
@@ -53,10 +53,8 @@ IMPLICIT NONE
 REAL(KIND=JPRB) :: ZHOOK_HANDLE
 !
 IF (LHOOK) CALL DR_HOOK('DEFAULT_SLT',0,ZHOOK_HANDLE)
-! ++ PIERRE / MARINE SSA - MODIF ++
-CEMISPARAM_SLT = 'Ova14'
-JPMODE_SLT     = 5
-! -- PIERRE / MARINE SSA - MODIF --
+CEMISPARAM_SLT = 'OvB21b'
+JPMODE_SLT     = 8
 LVARSIG_SLT    = .FALSE.
 LRGFIX_SLT     = .TRUE.
 IF (LHOOK) CALL DR_HOOK('DEFAULT_SLT',1,ZHOOK_HANDLE)
diff --git a/src/SURFEX/init_slt.F90 b/src/SURFEX/init_slt.F90
index 2182cbeb851c052aa64199836e92cc7eba49711f..51f0eeb3f1b6c55b71106f7c8a7bcd2302a84cf9 100644
--- a/src/SURFEX/init_slt.F90
+++ b/src/SURFEX/init_slt.F90
@@ -37,12 +37,11 @@ IF (LHOOK) CALL DR_HOOK('INIT_SLT',0,ZHOOK_HANDLE)
 !Get initial size distributions. This is cut and pasted
 !from dead routine dstpsd.F90
 !Check for different source parameterizations
-! Default : CEMISPARAM_SLT.eq."Ova14"
-
-  NSLTMDE = 5
+! Default : CEMISPARAM_SLT.eq."OvB21b"
+  NSLTMDE = 8
   CRGUNITS   = 'NUMB'
-  XEMISRADIUS_INI_SLT = (/0.009, 0.021, 0.045, 0.115, 0.415, 0.0, 0.0, 0.0/)
-  XEMISSIG_INI_SLT = (/1.37, 1.5, 1.42, 1.53, 1.85,0.,0.,0./)
+  XEMISRADIUS_INI_SLT = (/0.009, 0.021, 0.045, 0.115, 0.415, 2.5, 7.0, 25.0/)
+  XEMISSIG_INI_SLT = (/1.37, 1.5, 1.42, 1.53,1.70,1.80, 1.85, 2.1/)
 
 IF ((CEMISPARAM_SLT.eq."OvB21a").OR.(CEMISPARAM_SLT.eq."OvB21b")) THEN
   NSLTMDE = 8
@@ -50,6 +49,11 @@ IF ((CEMISPARAM_SLT.eq."OvB21a").OR.(CEMISPARAM_SLT.eq."OvB21b")) THEN
   XEMISRADIUS_INI_SLT = (/0.009, 0.021, 0.045, 0.115, 0.415, 2.5, 7.0, 25.0/)
   XEMISSIG_INI_SLT = (/1.37, 1.5, 1.42, 1.53,1.70,1.80, 1.85, 2.1/)
 
+ELSE IF (CEMISPARAM_SLT.eq."Ova14") THEN
+  NSLTMDE = 5
+  CRGUNITS   = 'NUMB'
+  XEMISRADIUS_INI_SLT = (/0.009, 0.021, 0.045, 0.115, 0.415, 0.0, 0.0, 0.0/)
+  XEMISSIG_INI_SLT = (/1.37, 1.5, 1.42, 1.53, 1.85,0.,0.,0./)
 
 ELSE IF (CEMISPARAM_SLT.eq."Vig01") THEN
    NSLTMDE = 5
diff --git a/src/SURFEX/modd_ch_isban.F90 b/src/SURFEX/modd_ch_isban.F90
index 6f97ac6cb325b72e30ec9459622f0a8365fb8786..b817f5810da674a597ab9be8eb737ff68bd2e07d 100644
--- a/src/SURFEX/modd_ch_isban.F90
+++ b/src/SURFEX/modd_ch_isban.F90
@@ -83,7 +83,7 @@ NULLIFY(YCH_ISBA%CAER_NAMES)
 NULLIFY(YCH_ISBA%CDSTNAMES)
 NULLIFY(YCH_ISBA%CSLTNAMES)
 NULLIFY(YCH_ISBA%CSNWNAMES)
-YCH_ISBA%CCHEM_SURF_FILE=' '
+YCH_ISBA%CCHEM_SURF_FILE='MNHC.input'
 YCH_ISBA%CCH_DRY_DEP=' '
 YCH_ISBA%CPARAMBVOC=' '
 YCH_ISBA%LCH_BIO_FLUX=.FALSE.
diff --git a/src/SURFEX/modd_ch_surfn.F90 b/src/SURFEX/modd_ch_surfn.F90
index 6114c9132946c60dc2147da017b12579bbafe55b..3108882b02d72437eaa285a5ebbe5cfa6e76d8fa 100644
--- a/src/SURFEX/modd_ch_surfn.F90
+++ b/src/SURFEX/modd_ch_surfn.F90
@@ -85,7 +85,7 @@ IF (LHOOK) CALL DR_HOOK("MODD_CH_SURF_N:CH_SURF_INIT",0,ZHOOK_HANDLE)
 YCH_SURF%CCH_EMIS='NONE'
 YCH_SURF%CCH_BIOEMIS='NONE'
 YCH_SURF%CCH_DMSEMIS='NONE'
-YCH_SURF%CCHEM_SURF_FILE='EXSEG1.nam'
+YCH_SURF%CCHEM_SURF_FILE='MNHC.input'
 YCH_SURF%LCH_SURF_EMIS=.FALSE.
 YCH_SURF%LCH_EMIS=.FALSE.
 YCH_SURF%LCH_BIOEMIS=.FALSE.
diff --git a/src/SURFEX/modd_ch_tebn.F90 b/src/SURFEX/modd_ch_tebn.F90
index 511d142c50b1756d0d906b28650751676e6662fb..3e43a15ad3d37f04480de0c9ee6a12a8edaafe9f 100644
--- a/src/SURFEX/modd_ch_tebn.F90
+++ b/src/SURFEX/modd_ch_tebn.F90
@@ -81,7 +81,7 @@ NULLIFY(YCH_TEB%CCH_NAMES)
 NULLIFY(YCH_TEB%CAER_NAMES)
 NULLIFY(YCH_TEB%CDSTNAMES)
 NULLIFY(YCH_TEB%CSLTNAMES)
-YCH_TEB%CCHEM_SURF_FILE=' '
+YCH_TEB%CCHEM_SURF_FILE='MNHC.input'
 YCH_TEB%CCH_DRY_DEP=' '
 YCH_TEB%CPARAMBVOC=' '
 YCH_TEB%LCH_BIO_FLUX=.FALSE.