From 3476a7df87848827b898152368eefb46ad91ca50 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Beno=C3=AEt=20Vi=C3=A9?= <benoit.vie@meteo.fr>
Date: Tue, 28 Sep 2021 11:59:30 +0200
Subject: [PATCH] Initial version

---
 src/MNH/ini_lima_cold_mixed.f90       |   1 +
 src/MNH/lima_adjust_split.f90         | 154 ++++++--------------------
 src/MNH/lima_bergeron.f90             |   6 +-
 src/MNH/lima_ice_aggregation_snow.f90 |  31 +++---
 src/MNH/lima_ice_deposition.f90       |  58 +++++-----
 src/MNH/lima_nucleation_procs.f90     |  39 +++++--
 src/MNH/lima_snow_deposition.f90      |  28 ++---
 src/MNH/lima_tendencies.f90           |  31 +++---
 src/MNH/modd_param_lima_cold.f90      |   3 +-
 9 files changed, 155 insertions(+), 196 deletions(-)

diff --git a/src/MNH/ini_lima_cold_mixed.f90 b/src/MNH/ini_lima_cold_mixed.f90
index cb427cdb4..9bcc16a27 100644
--- a/src/MNH/ini_lima_cold_mixed.f90
+++ b/src/MNH/ini_lima_cold_mixed.f90
@@ -647,6 +647,7 @@ XITAUTS_THRESHOLD = 7.5
 !*       6.4    Constants for snow aggregation
 !
 XCOLEXIS = 0.05    ! Temperature factor of the I+S collection efficiency
+XFIAGGS  = (XPI/4.0)*0.25*XCCS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0)
 XAGGS_CLARGE1 = XKER_ZRNIC_A2*ZGAMI(2)
 XAGGS_CLARGE2 = XKER_ZRNIC_A2*ZGAMS(2)
 XAGGS_RLARGE1 = XKER_ZRNIC_A2*ZGAMI(6)*XAI
diff --git a/src/MNH/lima_adjust_split.f90 b/src/MNH/lima_adjust_split.f90
index 6f7dbb738..6d8bd3dae 100644
--- a/src/MNH/lima_adjust_split.f90
+++ b/src/MNH/lima_adjust_split.f90
@@ -273,7 +273,7 @@ REAL, DIMENSION(SIZE(PRHODJ,1),SIZE(PRHODJ,2),SIZE(PRHODJ,3)) &
                             ZRV, ZRV2,  &
                             ZRC, ZRC2,  &
                             ZRI,  &
-                            ZSIGS, &
+                            ZSIGS, ZSRCS, &
                             ZW_MF
 LOGICAL, DIMENSION(SIZE(PRHODJ,1),SIZE(PRHODJ,2),SIZE(PRHODJ,3)) &
                          :: GMICRO ! Test where to compute cond/dep proc.
@@ -523,124 +523,42 @@ DO JITER =1,ITERMAX
 !
 !-------------------------------------------------------------------------------
 !
-!
-!
-!*              FULLY IMPLICIT CONDENSATION SCHEME
-!               ---------------------------------
-! 
-!*              select cases where r_c>0
-! 
-!
-      GMICRO(:,:,:) = .FALSE.
-      GMICRO(IIB:IIE,IJB:IJE,IKB:IKE) =( PRCS(IIB:IIE,IJB:IJE,IKB:IKE)>0. .AND.        &
-                                         PCCS(IIB:IIE,IJB:IJE,IKB:IKE)>0.      )
-      IMICRO = COUNTJV( GMICRO(:,:,:),I1(:),I2(:),I3(:))
-      IF( IMICRO >= 1 ) THEN
-         ALLOCATE(ZRVT(IMICRO))
-         ALLOCATE(ZRCT(IMICRO))
-!
-         ALLOCATE(ZRVS(IMICRO))
-         ALLOCATE(ZRCS(IMICRO))
-         ALLOCATE(ZCCS(IMICRO))
-         ALLOCATE(ZTHS(IMICRO))
-!
-         ALLOCATE(ZRHODREF(IMICRO))
-         ALLOCATE(ZZT(IMICRO))
-         ALLOCATE(ZPRES(IMICRO))
-         ALLOCATE(ZEXNREF(IMICRO))
-         ALLOCATE(ZZCPH(IMICRO))
-         DO JL=1,IMICRO
-            ZRVT(JL) = PRVT(I1(JL),I2(JL),I3(JL))
-            ZRCT(JL) = PRCT(I1(JL),I2(JL),I3(JL))
-            !
-            ZRVS(JL) = PRVS(I1(JL),I2(JL),I3(JL))
-            ZRCS(JL) = PRCS(I1(JL),I2(JL),I3(JL))
-            ZCCS(JL) = PCCS(I1(JL),I2(JL),I3(JL))
-            ZTHS(JL) = PTHS(I1(JL),I2(JL),I3(JL))
-            !
-            ZRHODREF(JL) = PRHODREF(I1(JL),I2(JL),I3(JL))
-            ZZT(JL) = ZT(I1(JL),I2(JL),I3(JL))
-            ZPRES(JL) = 2.0*PPABST(I1(JL),I2(JL),I3(JL))-PPABSM(I1(JL),I2(JL),I3(JL))
-            ZEXNREF(JL) = PEXNREF(I1(JL),I2(JL),I3(JL))
-            ZZCPH(JL) = ZCPH(I1(JL),I2(JL),I3(JL))
-         ENDDO
-         ALLOCATE(ZZW(IMICRO))
-         ALLOCATE(ZLVFACT(IMICRO))
-         ALLOCATE(ZRVSATW(IMICRO))
-         ALLOCATE(ZCND(IMICRO))
-         ZLVFACT(:) = (XLVTT+(XCPV-XCL)*(ZZT(:)-XTT))/ZZCPH(:) ! L_v/C_ph
-         ZZW(:) = EXP( XALPW - XBETAW/ZZT(:) - XGAMW*ALOG(ZZT(:) ) ) ! es_w
-         ZRVSATW(:) = ZEPS*ZZW(:) / ( ZPRES(:)-ZZW(:) )              ! r_sw
+      ZRV=PRVS*PTSTEP
+      ZRC=PRCS*PTSTEP
+      ZRI=PRIS*PTSTEP
+      ZSIGS=0
+      CALL CONDENSATION(IIU, IJU, IKU, IIB, IIE, IJB, IJE, IKB, IKE, 1, 'S',   &
+           HCONDENS, HLAMBDA3, &
+           PPABST, PZZ, PRHODREF, ZT, ZRV, ZRC, ZRI, PRSS*PTSTEP, PRGS*PTSTEP, &
+           ZSIGS, PMFCONV, PCLDFR, ZSRCS, .TRUE., OSIGMAS=.TRUE., &
+           PSIGQSAT=0., PLV=ZLV, PLS=ZLS, PCPH=ZCPH )
+
+      ZW1(:,:,:) = (ZRC(:,:,:) - PRCS(:,:,:)*PTSTEP) / PTSTEP       ! Pcon = ----------
+                                                      !         2 Delta t
 
-         IF (LADJ) THEN
-            ALLOCATE(ZRVSATW_PRIME(IMICRO))
-            ALLOCATE(ZAWW(IMICRO))
-            ALLOCATE(ZDELT1(IMICRO))
-            ALLOCATE(ZDELT2(IMICRO))
-            ZRVSATW_PRIME(:) = (( XBETAW/ZZT(:) - XGAMW ) / ZZT(:))  &  ! r'_sw
-                               * ZRVSATW(:) * ( 1. + ZRVSATW(:)/ZEPS )
-            ZAWW(:) = 1.0 + ZRVSATW_PRIME(:)*ZLVFACT(:)
-            ZDELT2(:) = (ZRVSATW_PRIME(:)*ZLVFACT(:)/ZAWW(:)) *                     &
-                        ( ((-2.*XBETAW+XGAMW*ZZT(:))/(XBETAW-XGAMW*ZZT(:))          &
-                        + (XBETAW/ZZT(:)-XGAMW)*(1.0+2.0*ZRVSATW(:)/ZEPS))/ZZT(:) )
-            ZDELT1(:) = (ZLVFACT(:)/ZAWW(:)) * ( ZRVSATW(:) - ZRVS(:)*ZDT )
-            ZCND(:) = - ZDELT1(:)*( 1.0 + 0.5*ZDELT1(:)*ZDELT2(:) ) / (ZLVFACT(:)*ZDT)
-            DEALLOCATE(ZRVSATW_PRIME)
-            DEALLOCATE(ZAWW)
-            DEALLOCATE(ZDELT1)
-            DEALLOCATE(ZDELT2)
-         ELSE
-            ALLOCATE(ZS(IMICRO))
-            ALLOCATE(ZZW2(IMICRO))
-            ALLOCATE(ZVEC1(IMICRO))
-            ALLOCATE(IVEC1(IMICRO))
-            ZVEC1(:) = MAX( 1.0001, MIN( REAL(NAHEN)-0.0001, XAHENINTP1 * ZZT(:) + XAHENINTP2 ) )
-            IVEC1(:) = INT( ZVEC1(:) )
-            ZVEC1(:) = ZVEC1(:) - REAL( IVEC1(:) )
-            ZS(:) = ZRVS(:)*PTSTEP / ZRVSATW(:) - 1.
-            ZZW(:) = ZCCS(:)*PTSTEP/(XLBC*ZCCS(:)/ZRCS(:))**XLBEXC
-            ZZW2(:) = XAHENG3(IVEC1(:)+1)*ZVEC1(:)-XAHENG3(IVEC1(:))*(ZVEC1(:)-1.)
-            ZCND(:) = 2.*3.14*1000.*ZZW2(:)*ZS(:)*ZZW(:)
-            DEALLOCATE(ZS)
-            DEALLOCATE(ZZW2)
-            DEALLOCATE(ZVEC1)
-            DEALLOCATE(IVEC1)
-         END IF
-!
-!
-! Integration
-!
-         WHERE( ZCND(:) < 0.0 )
-            ZCND(:) = MAX ( ZCND(:), -ZRCS(:) )
-         ELSEWHERE
-            ZCND(:) = MIN ( ZCND(:),  ZRVS(:) )
-         END WHERE
-         ZRVS(:) = ZRVS(:) - ZCND(:)
-         ZRCS(:) = ZRCS(:) + ZCND(:)
-         ZTHS(:) = ZTHS(:) + ZCND(:) * ZLVFACT(:) / ZEXNREF(:)
-!
-         ZW(:,:,:) = PRVS(:,:,:)
-         PRVS(:,:,:) = UNPACK( ZRVS(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:) )
-         ZW(:,:,:) = PRCS(:,:,:)
-         PRCS(:,:,:) = UNPACK( ZRCS(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:) )
-         ZW(:,:,:) = PTHS(:,:,:)
-         PTHS(:,:,:) = UNPACK( ZTHS(:),MASK=GMICRO(:,:,:),FIELD=ZW(:,:,:) )
-!
-         DEALLOCATE(ZRVT)
-         DEALLOCATE(ZRCT)
-         DEALLOCATE(ZRVS)
-         DEALLOCATE(ZRCS)
-         DEALLOCATE(ZTHS)
-         DEALLOCATE(ZRHODREF)
-         DEALLOCATE(ZZT)
-         DEALLOCATE(ZPRES)
-         DEALLOCATE(ZEXNREF)
-         DEALLOCATE(ZZCPH)
-         DEALLOCATE(ZZW)
-         DEALLOCATE(ZLVFACT)
-         DEALLOCATE(ZRVSATW)
-         DEALLOCATE(ZCND)
-      END IF ! IMICRO
+      ZW2(:,:,:) = (ZRI(:,:,:) - PRIS(:,:,:)*PTSTEP) / PTSTEP       ! idem ZW1 but for Ri
+!
+!*       5.1    compute the sources
+!
+      WHERE( ZW1(:,:,:) < 0.0 )
+         ZW1(:,:,:) = MAX ( ZW1(:,:,:), -PRCS(:,:,:) )
+      ELSEWHERE
+         ZW1(:,:,:) = MIN ( ZW1(:,:,:),  PRVS(:,:,:) )
+      END WHERE
+      PRVS(:,:,:) = PRVS(:,:,:) - ZW1(:,:,:)
+      PRCS(:,:,:) = PRCS(:,:,:) + ZW1(:,:,:)
+      PTHS(:,:,:) = PTHS(:,:,:) +        &
+           ZW1(:,:,:) * ZLV(:,:,:) / (ZCPH(:,:,:) * PEXNREF(:,:,:))
+!
+      WHERE( ZW2(:,:,:) < 0.0 )
+         ZW2(:,:,:) = MAX ( ZW2(:,:,:), -PRIS(:,:,:) )
+      ELSEWHERE
+         ZW2(:,:,:) = MIN ( ZW2(:,:,:),  PRVS(:,:,:) )
+      END WHERE
+      PRVS(:,:,:) = PRVS(:,:,:) - ZW2(:,:,:)
+      PRIS(:,:,:) = PRIS(:,:,:) + ZW2(:,:,:)
+      PTHS(:,:,:) = PTHS(:,:,:) +        &
+           ZW2(:,:,:) * ZLS(:,:,:) / (ZCPH(:,:,:) * PEXNREF(:,:,:))
 !
    END IF ! end of adjustment procedure (test on OSUBG_COND)
 !
diff --git a/src/MNH/lima_bergeron.f90 b/src/MNH/lima_bergeron.f90
index 9105c78d6..9afd9ac36 100644
--- a/src/MNH/lima_bergeron.f90
+++ b/src/MNH/lima_bergeron.f90
@@ -111,9 +111,9 @@ WHERE( (PRCT(:)>XRTMIN(2)) .AND. (PRIT(:)>XRTMIN(4)) .AND. (PCIT(:)>XCTMIN(4)) .
    P_TH_BERFI(:) = - P_RC_BERFI(:)*(PLSFACT(:)-PLVFACT(:))
 END WHERE
 !
-PA_RC(:) = PA_RC(:) + P_RC_BERFI(:)
-PA_RI(:) = PA_RI(:) - P_RC_BERFI(:)
-PA_TH(:) = PA_TH(:) + P_TH_BERFI(:)
+!PA_RC(:) = PA_RC(:) + P_RC_BERFI(:)
+!PA_RI(:) = PA_RI(:) - P_RC_BERFI(:)
+!PA_TH(:) = PA_TH(:) + P_TH_BERFI(:)
 !
 !
 !-------------------------------------------------------------------------------
diff --git a/src/MNH/lima_ice_aggregation_snow.f90 b/src/MNH/lima_ice_aggregation_snow.f90
index 15e01ec84..698697a2c 100644
--- a/src/MNH/lima_ice_aggregation_snow.f90
+++ b/src/MNH/lima_ice_aggregation_snow.f90
@@ -59,9 +59,9 @@ END MODULE MODI_LIMA_ICE_AGGREGATION_SNOW
 !              ------------
 !
 USE MODD_CST,             ONLY : XTT
-USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN
+USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN, XCEXVT
 USE MODD_PARAM_LIMA_COLD, ONLY : XBI, XCCS, XCXS, XCOLEXIS, XAGGS_CLARGE1, XAGGS_CLARGE2, &
-                                 XAGGS_RLARGE1, XAGGS_RLARGE2
+                                 XAGGS_RLARGE1, XAGGS_RLARGE2, XFIAGGS
 !
 IMPLICIT NONE
 !
@@ -100,17 +100,22 @@ 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)
-   ZZW3(:) = ZZW2(:)*(XAGGS_CLARGE1+XAGGS_CLARGE2*ZZW1(:))
-!
-   P_CI_AGGS(:) = - ZZW3(:)
-!
-   ZZW2(:) = ZZW2(:) / PLBDI(:)**XBI
-   ZZW2(:) = ZZW2(:)*(XAGGS_RLARGE1+XAGGS_RLARGE2*ZZW1(:))
-!
-   P_RI_AGGS(:) = - ZZW2(:)
+!!$   ZZW1(:) = (PLBDI(:) / PLBDS(:))**3
+!!$   ZZW2(:) = (PCIT(:)*(XCCS*PLBDS(:)**XCXS)/PRHODREF(:)*EXP( XCOLEXIS*(PT(:)-XTT) )) &
+!!$        / (PLBDI(:)**3)
+!!$   ZZW3(:) = ZZW2(:)*(XAGGS_CLARGE1+XAGGS_CLARGE2*ZZW1(:))
+!!$!
+!!$   P_CI_AGGS(:) = - ZZW3(:)
+!!$!
+!!$   ZZW2(:) = ZZW2(:) / PLBDI(:)**XBI
+!!$   ZZW2(:) = ZZW2(:)*(XAGGS_RLARGE1+XAGGS_RLARGE2*ZZW1(:))
+!
+   ZZW1(:) = XFIAGGS * EXP( XCOLEXIS*(PT(:)-XTT) ) &
+                     * PRIT(:)                     &
+                     * PLBDS(:)**(1.-0.27-2.)      &
+                     * PRHODREF(:)**(-XCEXVT)
+!
+   P_RI_AGGS(:) = - ZZW1(:)
 END WHERE
 !
 !
diff --git a/src/MNH/lima_ice_deposition.f90 b/src/MNH/lima_ice_deposition.f90
index 752feb93b..34dd0c2c1 100644
--- a/src/MNH/lima_ice_deposition.f90
+++ b/src/MNH/lima_ice_deposition.f90
@@ -8,16 +8,17 @@
 !      #####################
 !
 INTERFACE
-      SUBROUTINE LIMA_ICE_DEPOSITION (PTSTEP, LDCOMPUTE,                 &
-                                      PRHODREF, PSSI, PAI, PCJ, PLSFACT, &
-                                      PRIT, PCIT, PLBDI,                 &
-                                      P_TH_DEPI, P_RI_DEPI,              &
-                                      P_RI_CNVS, P_CI_CNVS               )
+      SUBROUTINE LIMA_ICE_DEPOSITION (PTSTEP, LDCOMPUTE,                     &
+                                      PRHODREF, PT, PSSI, PAI, PCJ, PLSFACT, &
+                                      PRIT, PCIT, PLBDI,                     &
+                                      P_TH_DEPI, P_RI_DEPI,                  &
+                                      P_RI_CNVS, P_CI_CNVS                   )
 !
 REAL,                 INTENT(IN)    :: PTSTEP
 LOGICAL, DIMENSION(:),INTENT(IN)    :: LDCOMPUTE
 !
 REAL, DIMENSION(:),   INTENT(IN)    :: PRHODREF! Reference density
+REAL, DIMENSION(:),   INTENT(IN)    :: PT  ! abs. pressure at time t
 REAL, DIMENSION(:),   INTENT(IN)    :: PSSI  ! abs. pressure at time t
 REAL, DIMENSION(:),   INTENT(IN)    :: PAI  ! abs. pressure at time t
 REAL, DIMENSION(:),   INTENT(IN)    :: PCJ  ! abs. pressure at time t
@@ -40,7 +41,7 @@ END MODULE MODI_LIMA_ICE_DEPOSITION
 !
 !     ##########################################################################
 SUBROUTINE LIMA_ICE_DEPOSITION (PTSTEP, LDCOMPUTE,                        &
-                                PRHODREF,  PSSI, PAI, PCJ, PLSFACT,       &
+                                PRHODREF, PT,  PSSI, PAI, PCJ, PLSFACT,   &
                                 PRIT, PCIT, PLBDI,                        &
                                 P_TH_DEPI, P_RI_DEPI,                     &
                                 P_RI_CNVS, P_CI_CNVS                      )
@@ -72,7 +73,7 @@ SUBROUTINE LIMA_ICE_DEPOSITION (PTSTEP, LDCOMPUTE,                        &
 !*       0.    DECLARATIONS
 !              ------------
 !
-USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN, XALPHAI, XALPHAS, XNUI, XNUS 
+USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN, XALPHAI, XALPHAS, XNUI, XNUS, LSNOW 
 USE MODD_PARAM_LIMA_COLD, ONLY : XCXS, XCCS, &
                                  XLBDAS_MAX, XDSCNVI_LIM, XLBDASCNVI_MAX,     &
                                  XC0DEPSI, XC1DEPSI, XR0DEPSI, XR1DEPSI,      &
@@ -82,7 +83,7 @@ USE MODD_PARAM_LIMA_COLD, ONLY : XCXS, XCCS, &
                                  XCOLEXIS, XAGGS_CLARGE1, XAGGS_CLARGE2,      &
                                  XAGGS_RLARGE1, XAGGS_RLARGE2,                &
                                  XDI, X0DEPI, X2DEPI
-
+USE MODD_CST,             ONLY : XTT
 !
 IMPLICIT NONE
 !
@@ -92,6 +93,7 @@ REAL,                 INTENT(IN)    :: PTSTEP
 LOGICAL, DIMENSION(:),INTENT(IN)    :: LDCOMPUTE
 !
 REAL, DIMENSION(:),   INTENT(IN)    :: PRHODREF! Reference density
+REAL, DIMENSION(:),   INTENT(IN)    :: PT  ! abs. pressure at time t
 REAL, DIMENSION(:),   INTENT(IN)    :: PSSI  ! abs. pressure at time t
 REAL, DIMENSION(:),   INTENT(IN)    :: PAI  ! abs. pressure at time t
 REAL, DIMENSION(:),   INTENT(IN)    :: PCJ  ! abs. pressure at time t
@@ -111,7 +113,7 @@ REAL, DIMENSION(:),   INTENT(OUT)   :: P_CI_CNVS
 !*       0.2   Declarations of local variables :
 !
 LOGICAL, DIMENSION(SIZE(PRHODREF)) :: GMICRO ! Computations only where necessary
-REAL,    DIMENSION(SIZE(PRHODREF)) :: ZZW, ZZW2, ZZX ! Work array
+REAL,    DIMENSION(SIZE(PRHODREF)) :: ZZW, ZZW2, ZZX, ZCRIAUTI ! Work array
 !
 !
 !-------------------------------------------------------------------------------
@@ -136,13 +138,13 @@ WHERE( GMICRO )
 !        -----------------------------------------------
 !
 !
-   ZZW(:) = 0.0
-   WHERE ( (PRIT(:)>XRTMIN(4)) .AND. (PCIT(:)>XCTMIN(4)) )
-      ZZW(:) = ( PSSI(:) / PAI(:) ) * PCIT(:) *        &
-        ( X0DEPI/PLBDI(:)+X2DEPI*PCJ(:)*PCJ(:)/PLBDI(:)**(XDI+2.0) )
-   END WHERE
+!   ZZW(:) = 0.0
+!   WHERE ( (PRIT(:)>XRTMIN(4)) .AND. (PCIT(:)>XCTMIN(4)) )
+!      ZZW(:) = ( PSSI(:) / PAI(:) ) * PCIT(:) *        &
+!        ( X0DEPI/PLBDI(:)+X2DEPI*PCJ(:)*PCJ(:)/PLBDI(:)**(XDI+2.0) )
+!   END WHERE
 !
-   P_RI_DEPI(:) = ZZW(:)
+!   P_RI_DEPI(:) = ZZW(:)
 !!$   P_TH_DEPI(:) = P_RI_DEPI(:) * PLSFACT(:)
 !
 !!$   PA_TH(:) = PA_TH(:) + P_TH_DEPI(:)
@@ -154,20 +156,26 @@ WHERE( GMICRO )
 !        ------------------------------------------------
 !
 !
-   ZZW(:) = 0.0
-   ZZW2(:) = 0.0
-   WHERE ( (PLBDI(:)<XLBDAICNVS_LIM) .AND. (PCIT(:)>XCTMIN(4)) &
-                                     .AND. (PSSI(:)>0.0)       )
-      ZZW(:) = (PLBDI(:)*XDICNVS_LIM)**(XALPHAI)
-      ZZX(:) = ( PSSI(:)/PAI(:) )*PCIT(:) * (ZZW(:)**XNUI) *EXP(-ZZW(:))
+!!$   ZZW(:) = 0.0
+!!$   ZZW2(:) = 0.0
+!!$   WHERE ( (PLBDI(:)<XLBDAICNVS_LIM) .AND. (PCIT(:)>XCTMIN(4)) &
+!!$                                     .AND. (PSSI(:)>0.0)       )
+!!$      ZZW(:) = (PLBDI(:)*XDICNVS_LIM)**(XALPHAI)
+!!$      ZZX(:) = ( PSSI(:)/PAI(:) )*PCIT(:) * (ZZW(:)**XNUI) *EXP(-ZZW(:))
+!!$!
+!!$      ZZW(:) = ( XR0DEPIS + XR1DEPIS*PCJ(:) )*ZZX(:)                             
+!!$!
+!!$      ZZW2(:) = ZZW(:) * (XC0DEPIS+XC1DEPIS*PCJ(:)) / (XR0DEPIS+XR1DEPIS*PCJ(:))
+!!$   END WHERE
 !
-      ZZW(:) = ( XR0DEPIS + XR1DEPIS*PCJ(:) )*ZZX(:)                             
-!
-      ZZW2(:) = ZZW(:) * (XC0DEPIS+XC1DEPIS*PCJ(:)) / (XR0DEPIS+XR1DEPIS*PCJ(:))
+   ZCRIAUTI(:)=MIN(0.2E-4,10**(0.06*(PT(:)-XTT)-3.5))
+   ZZW(:) = 0.0
+   WHERE ( (PRIT(:)>XRTMIN(4)))
+      ZZW(:)   = 1.E-3 * EXP( 0.015*(PT(:)-XTT) ) * MAX( PRIT(:)-ZCRIAUTI(:),0.0 )
    END WHERE
 !
 P_RI_CNVS(:) = - ZZW(:)
-P_CI_CNVS(:) = - ZZW2(:)
+!P_CI_CNVS(:) = - ZZW2(:)
 !
 !
 END WHERE
diff --git a/src/MNH/lima_nucleation_procs.f90 b/src/MNH/lima_nucleation_procs.f90
index 122d4b3c8..430099268 100644
--- a/src/MNH/lima_nucleation_procs.f90
+++ b/src/MNH/lima_nucleation_procs.f90
@@ -86,6 +86,7 @@ use modd_budget,     only: lbu_enable, lbudget_th, lbudget_rv, lbudget_rc, lbudg
                            NBUDGET_TH, NBUDGET_RV, NBUDGET_RC, NBUDGET_RI, NBUDGET_SV1, &
                            tbudgets
 USE MODD_IO,         ONLY: TFILEDATA
+USE MODD_PARAMETERS, ONLY : JPHEXT, JPVEXT
 USE MODD_NSV,        ONLY : NSV_LIMA_NC, NSV_LIMA_NR, NSV_LIMA_CCN_FREE, NSV_LIMA_CCN_ACTI, &
                             NSV_LIMA_NI, NSV_LIMA_IFN_FREE, NSV_LIMA_IFN_NUCL, NSV_LIMA_IMM_NUCL, NSV_LIMA_HOM_HAZE
 USE MODD_PARAM_LIMA, ONLY : LCOLD, LNUCL, LMEYERS, LSNOW, LWARM, LACTI, LRAIN, LHHONI,  &
@@ -98,6 +99,7 @@ USE MODI_LIMA_CCN_ACTIVATION
 USE MODI_LIMA_CCN_HOM_FREEZING
 USE MODI_LIMA_MEYERS_NUCLEATION
 USE MODI_LIMA_PHILLIPS_IFN_NUCLEATION
+USE MODE_RAIN_ICE_NUCLEATION
 !
 !-------------------------------------------------------------------------------
 !
@@ -142,6 +144,7 @@ REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PPRCFR     ! Precipitation fraction
 !-------------------------------------------------------------------------------
 !
 REAL, DIMENSION(SIZE(PT,1),SIZE(PT,2),SIZE(PT,3))          :: Z_TH_HIND, Z_RI_HIND, Z_CI_HIND, Z_TH_HINC, Z_RC_HINC, Z_CC_HINC
+REAL, DIMENSION(SIZE(PT,1),SIZE(PT,2),SIZE(PT,3))          :: ZTHS, ZRIS, ZRVS, ZRHT, ZCIT, ZT
 !
 integer :: idx
 INTEGER :: JL
@@ -258,15 +261,35 @@ END IF
 !-------------------------------------------------------------------------------
 !
 IF (LCOLD .AND. LNUCL .AND. LMEYERS) THEN
-   CALL LIMA_MEYERS_NUCLEATION (PTSTEP,                                     &
-                                PRHODREF, PEXNREF, PPABST,                  &
-                                PTHT, PRVT, PRCT, PRRT, PRIT, PRST, PRGT,   &
-                                PCCT, PCIT, PINT,                           &
-                                Z_TH_HIND, Z_RI_HIND, Z_CI_HIND,            &
-                                Z_TH_HINC, Z_RC_HINC, Z_CC_HINC,            &
-                                PICEFR                                      )
+!   CALL LIMA_MEYERS_NUCLEATION (PTSTEP,                                     &
+!                                PRHODREF, PEXNREF, PPABST,                  &
+!                                PTHT, PRVT, PRCT, PRRT, PRIT, PRST, PRGT,   &
+!                                PCCT, PCIT, PINT,                           &
+!                                Z_TH_HIND, Z_RI_HIND, Z_CI_HIND,            &
+!                                Z_TH_HINC, Z_RC_HINC, Z_CC_HINC,            &
+!                                PICEFR                                      )
   WHERE(PICEFR(:,:,:)<1.E-10 .AND. PRIT(:,:,:)>XRTMIN(4) .AND. PCIT(:,:,:)>XCTMIN(4)) PICEFR(:,:,:)=1.
-!
+  !
+  ZTHS=PTHT/PTSTEP
+  ZRVS=PRVT/PTSTEP
+  ZRIS=PRIT/PTSTEP
+  ZRHT=0.
+  ZCIT=PCIT
+  ZT=PT
+  CALL RAIN_ICE_NUCLEATION(1+JPHEXT, SIZE(PT,1)-JPHEXT, 1+JPHEXT, SIZE(PT,2)-JPHEXT, 1+JPVEXT, SIZE(PT,3)-JPVEXT, 6, &
+                           PTSTEP, PTHT, PPABST, PRHODJ, PRHODREF, PRVT, PRCT, PRRT, PRIT, PRST, PRGT,               &
+                           ZCIT, PEXNREF, ZTHS, ZRVS, ZRIS, ZT, ZRHT)
+  !
+  Z_TH_HIND=ZTHS*PTSTEP-PTHT
+  Z_RI_HIND=ZRIS*PTSTEP-PRIT
+  Z_CI_HIND=ZCIT-PCIT
+  PCIT=ZCIT
+  PRIT=ZRIS*PTSTEP
+  PTHT=ZTHS*PTSTEP
+  Z_TH_HINC=0.
+  Z_RC_HINC=0.
+  Z_CC_HINC=0.
+  !
   if ( lbu_enable ) then
     if ( lbudget_th ) call Budget_store_add( tbudgets(NBUDGET_TH), 'HIND',  z_th_hind(:, :, :) * prhodj(:, :, :) / ptstep )
     if ( lbudget_rv ) call Budget_store_add( tbudgets(NBUDGET_RV), 'HIND', -z_ri_hind(:, :, :) * prhodj(:, :, :) / ptstep )
diff --git a/src/MNH/lima_snow_deposition.f90 b/src/MNH/lima_snow_deposition.f90
index d39f714d8..7d923f1c9 100644
--- a/src/MNH/lima_snow_deposition.f90
+++ b/src/MNH/lima_snow_deposition.f90
@@ -127,20 +127,20 @@ WHERE( GMICRO )
 !        ----------------------------------------
 !
 !
-   ZZW2(:) = 0.0
-   ZZW(:) = 0.0
-   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(:))
-!
-      ZZW(:) = ( XR0DEPSI+XR1DEPSI*PCJ(:) )*ZZX(:)
-!
-      ZZW2(:) = ZZW(:)*( XC0DEPSI+XC1DEPSI*PCJ(:) )/( XR0DEPSI+XR1DEPSI*PCJ(:) )
-   END WHERE
-!
-   P_RI_CNVI(:) = ZZW(:)
-   P_CI_CNVI(:) = ZZW2(:)
+!!$   ZZW2(:) = 0.0
+!!$   ZZW(:) = 0.0
+!!$   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(:))
+!!$!
+!!$      ZZW(:) = ( XR0DEPSI+XR1DEPSI*PCJ(:) )*ZZX(:)
+!!$!
+!!$      ZZW2(:) = ZZW(:)*( XC0DEPSI+XC1DEPSI*PCJ(:) )/( XR0DEPSI+XR1DEPSI*PCJ(:) )
+!!$   END WHERE
+!!$!
+!!$   P_RI_CNVI(:) = ZZW(:)
+!!$   P_CI_CNVI(:) = ZZW2(:)
 !
 !
 !*       2.2    Deposition of water vapor on r_s: RVDEPS
diff --git a/src/MNH/lima_tendencies.f90 b/src/MNH/lima_tendencies.f90
index bd98d503c..bcf8e41ef 100644
--- a/src/MNH/lima_tendencies.f90
+++ b/src/MNH/lima_tendencies.f90
@@ -595,7 +595,7 @@ IF (LCOLD) THEN
    ! Includes vapour deposition on ice, ice -> snow conversion
    !
    CALL LIMA_ICE_DEPOSITION (PTSTEP, LDCOMPUTE,                 & ! depends on IF, PF
-                             PRHODREF, ZSSI, ZAI, ZCJ, ZLSFACT, &
+                             PRHODREF, ZT, ZSSI, ZAI, ZCJ, ZLSFACT, &
                              PRIT/ZIF1D, PCIT/ZIF1D, ZLBDI,     &
                              P_TH_DEPI, P_RI_DEPI,              &
                              P_RI_CNVS, P_CI_CNVS               )
@@ -617,11 +617,11 @@ IF (LCOLD .AND. LSNOW) THEN
    !
    ! Includes vapour deposition on snow, snow -> ice conversion
    !
-   CALL LIMA_SNOW_DEPOSITION (LDCOMPUTE,                         & ! depends on IF, PF
+   CALL LIMA_SNOW_DEPOSITION (LDCOMPUTE,                             & ! depends on IF, PF
                               PRHODREF, ZSSI, ZAI, ZCJ, ZLSFACT, &
-                              PRST/ZPF1D, ZLBDS,                 &
-                              P_RI_CNVI, P_CI_CNVI,              &
-                              P_TH_DEPS, P_RS_DEPS               )
+                              PRST/ZPF1D, ZLBDS,                     &
+                              P_RI_CNVI, P_CI_CNVI,                  &
+                              P_TH_DEPS, P_RS_DEPS                   )
    !
    P_RI_CNVI(:) = P_RI_CNVI(:) * ZPF1D(:)
    P_CI_CNVI(:) = P_CI_CNVI(:) * ZPF1D(:)
@@ -667,16 +667,19 @@ IF (LWARM .AND. LCOLD) THEN
    PA_TH(:) = PA_TH(:) + P_TH_DEPG(:)
 END IF
 !
-!!$IF (LWARM .AND. LCOLD) THEN
-!!$   CALL LIMA_BERGERON (LDCOMPUTE,                         & ! depends on CF, IF
-!!$                       PRCT, PRIT, PCIT, ZLBDI,           &
-!!$                       ZSSIW, ZAI, ZCJ, ZLVFACT, ZLSFACT, &
-!!$                       P_TH_BERFI, P_RC_BERFI,            &
-!!$                       PA_TH, PA_RC, PA_RI                )
-!!$END IF
-P_TH_BERFI(:) = 0.
-P_RC_BERFI(:) = 0.
+IF (LWARM .AND. LCOLD) THEN
+   CALL LIMA_BERGERON (LDCOMPUTE,                         & ! depends on CF, IF
+                       PRCT, PRIT, PCIT, ZLBDI,           &
+                       ZSSIW, ZAI, ZCJ, ZLVFACT, ZLSFACT, &
+                       P_TH_BERFI, P_RC_BERFI,            &
+                       PA_TH, PA_RC, PA_RI                )
+END IF
+!P_TH_BERFI(:) = 0.
+!P_RC_BERFI(:) = 0.
 !
+PA_RC(:) = PA_RC(:) + P_RC_BERFI(:)
+PA_RI(:) = PA_RI(:) - P_RC_BERFI(:)
+PA_TH(:) = PA_TH(:) + P_TH_BERFI(:)
 !
 IF (LWARM .AND. LCOLD .AND. LSNOW) THEN
      !
diff --git a/src/MNH/modd_param_lima_cold.f90 b/src/MNH/modd_param_lima_cold.f90
index 64494219e..0618a4ebc 100644
--- a/src/MNH/modd_param_lima_cold.f90
+++ b/src/MNH/modd_param_lima_cold.f90
@@ -106,7 +106,8 @@ REAL,SAVE :: XDICNVS_LIM, XLBDAICNVS_LIM,      & ! Constants for pristine ice
 !
 REAL,SAVE :: XCOLEXIS,                         & ! Constants for snow 
     	     XAGGS_CLARGE1,XAGGS_CLARGE2,      & ! aggregation : AGG
-             XAGGS_RLARGE1,XAGGS_RLARGE2
+             XAGGS_RLARGE1,XAGGS_RLARGE2,      &
+             XFIAGGS
 !
 !??????????????????
 REAL,SAVE :: XKER_ZRNIC_A1,XKER_ZRNIC_A2         ! Long-Zrnic Kernels (ini_ice_coma)
-- 
GitLab