From 819a47fe3926d2b4826a0bdd0ccb43b57b9354c4 Mon Sep 17 00:00:00 2001
From: VIE Benoit <vie@sxphynh>
Date: Thu, 19 May 2022 16:59:20 +0200
Subject: [PATCH] LIMA1M : first run passed

---
 src/MNH/default_desfmn.f90               |  3 ++-
 src/MNH/ini_lima_warm.f90                |  4 ++--
 src/MNH/lima.f90                         |  9 ++++++--
 src/MNH/lima_adjust_split.f90            | 16 +++++++--------
 src/MNH/lima_droplets_autoconversion.f90 |  6 ++++--
 src/MNH/lima_droplets_hom_freezing.f90   | 26 +++---------------------
 src/MNH/lima_tendencies.f90              | 13 ++++++++----
 src/MNH/modd_param_lima.f90              |  1 +
 src/MNH/modn_param_lima.f90              |  2 +-
 src/MNH/set_conc_lima.f90                |  9 ++++----
 src/MNH/sources_neg_correct.f90          |  2 +-
 11 files changed, 43 insertions(+), 48 deletions(-)

diff --git a/src/MNH/default_desfmn.f90 b/src/MNH/default_desfmn.f90
index 1925ef8ba..58ac77649 100644
--- a/src/MNH/default_desfmn.f90
+++ b/src/MNH/default_desfmn.f90
@@ -283,7 +283,7 @@ USE MODD_PARAM_LIMA, ONLY : LCOLD, LNUCL, LSEDI, LHHONI, LSNOW, LHAIL, LMEYERS,
                             XFACTNUC_DEP, XFACTNUC_CON,                                 &
                             OWARM=>LWARM, LACTI, ORAIN=>LRAIN, OSEDC=>LSEDC,            &
                             OACTIT=>LACTIT, LBOUND, LSPRO, LADJ, LKHKO, NMOM_C, NMOM_R, &
-                            NMOD_CCN, XCCN_CONC,                                        &
+                            NMOD_CCN, XCCN_CONC, LKESSLERAC,                            &
                             LCCN_HOM, CCCN_MODES,                                       &
                             YALPHAR=>XALPHAR, YNUR=>XNUR,                               &
                             YALPHAC=>XALPHAC, YNUC=>XNUC, CINI_CCN=>HINI_CCN,           &
@@ -1000,6 +1000,7 @@ IF (KMI == 1) THEN
    ODEPOC = .FALSE.
    LBOUND = .FALSE.
    OACTTKE = .TRUE.
+   LKESSLERAC = .FALSE.
 !
    NMOM_C = 2
    NMOM_R = 2
diff --git a/src/MNH/ini_lima_warm.f90 b/src/MNH/ini_lima_warm.f90
index a8e8cdcce..08d1e85a5 100644
--- a/src/MNH/ini_lima_warm.f90
+++ b/src/MNH/ini_lima_warm.f90
@@ -164,10 +164,10 @@ XLBC   = XAR*ZGAMC(2)
 XLBEXC = 1.0/XBC
 !
 IF (NMOM_R.EQ.1) THEN
-   XLBEXR = 1.0/XBR
-   XLBR   = ( XAR*XCCR*MOMG(XALPHAR,XNUR,XBR) )**(-XLBEXR)
    XCCR   = 8.E6
    XCXR   = -1.
+   XLBEXR = 1.0/XBR
+   XLBR   = ( XAR*XCCR*MOMG(XALPHAR,XNUR,XBR) )**(-XLBEXR)
 ELSE
    XLBR   = XAR*ZGAMR(2)
    XLBEXR = 1.0/XBR
diff --git a/src/MNH/lima.f90 b/src/MNH/lima.f90
index 66e16a895..5dd354385 100644
--- a/src/MNH/lima.f90
+++ b/src/MNH/lima.f90
@@ -548,8 +548,13 @@ IF ( KRR .GE. 7 ) ZRHS(:,:,:) = PRS(:,:,:,7)
 !
 ! Concentrations
 !
-IF ( LWARM .AND. NMOM_C.GE.2)             ZCCT(:,:,:)   = PSVS(:,:,:,NSV_LIMA_NC) * PTSTEP
-IF ( LWARM .AND. NMOM_C.GE.2)             ZCCS(:,:,:)   = PSVS(:,:,:,NSV_LIMA_NC)
+IF ( NMOM_C.GE.2) THEN
+   ZCCT(:,:,:)   = PSVS(:,:,:,NSV_LIMA_NC) * PTSTEP
+   ZCCS(:,:,:)   = PSVS(:,:,:,NSV_LIMA_NC)
+ELSE
+   ZCCT(:,:,:)   = 300.E6 / PRHODREF(:,:,:)
+   ZCCS(:,:,:)   = ZCCT(:,:,:) / PTSTEP
+END IF
 IF ( LWARM .AND. LRAIN .AND. NMOM_R.GE.2) ZCRT(:,:,:)   = PSVS(:,:,:,NSV_LIMA_NR) * PTSTEP
 IF ( LWARM .AND. LRAIN .AND. NMOM_R.GE.2) ZCRS(:,:,:)   = PSVS(:,:,:,NSV_LIMA_NR)
 IF ( LCOLD .AND. NMOM_I.GE.2)             ZCIT(:,:,:)   = PSVS(:,:,:,NSV_LIMA_NI) * PTSTEP
diff --git a/src/MNH/lima_adjust_split.f90 b/src/MNH/lima_adjust_split.f90
index 19b6b3778..a1afd1350 100644
--- a/src/MNH/lima_adjust_split.f90
+++ b/src/MNH/lima_adjust_split.f90
@@ -374,11 +374,11 @@ PCIT(:,:,:) = 0.
 PCCS(:,:,:) = 0.
 ! PCIS(:,:,:) = 0.
 !
-IF ( LWARM ) PCCT(:,:,:) = PSVS(:,:,:,NSV_LIMA_NC)*PTSTEP
-IF ( LCOLD ) PCIT(:,:,:) = PSVT(:,:,:,NSV_LIMA_NI)
+IF ( LWARM .AND. NMOM_C.GE.2 ) PCCT(:,:,:) = PSVS(:,:,:,NSV_LIMA_NC)*PTSTEP
+IF ( LCOLD .AND. NMOM_I.GE.2 ) PCIT(:,:,:) = PSVT(:,:,:,NSV_LIMA_NI)
 !
-IF ( LWARM ) PCCS(:,:,:) = PSVS(:,:,:,NSV_LIMA_NC)
-! IF ( LCOLD ) PCIS(:,:,:) = PSVS(:,:,:,NSV_LIMA_NI)
+IF ( LWARM .AND. NMOM_C.GE.2 ) PCCS(:,:,:) = PSVS(:,:,:,NSV_LIMA_NC)
+! IF ( LCOLD .AND. NMOM_I.GE.2 ) PCIS(:,:,:) = PSVS(:,:,:,NSV_LIMA_NI)
 !
 IF ( LSCAV .AND. LAERO_MASS ) PMAS(:,:,:) = PSVS(:,:,:,NSV_LIMA_SCAVMASS)
 ! 
@@ -519,7 +519,7 @@ DO JITER =1,ITERMAX
            Z_SIGS, PMFCONV, PCLDFR, Z_SRCS, GUSERI, G_SIGMAS,                  &
            Z_SIGQSAT, PLV=ZLV, PLS=ZLS, PCPH=ZCPH )
    END IF
-   IF (OSUBG_COND) THEN
+   IF (OSUBG_COND .AND. NMOM_C.GE.2) THEN
       PSRCS=Z_SRCS
       ZW_MF=0.
       CALL LIMA_CCN_ACTIVATION (TPFILE,                          &
@@ -630,7 +630,7 @@ END IF
 !
 ZMASK(:,:,:) = 0.0
 ZW(:,:,:) = 0.
-WHERE (PRCS(:,:,:) <= ZRTMIN(2) .OR. PCCS(:,:,:) <= ZCTMIN(2)) 
+WHERE (PRCS(:,:,:) <= ZRTMIN(2) .OR. PCCS(:,:,:) <0.) 
    PRVS(:,:,:) = PRVS(:,:,:) + PRCS(:,:,:) 
    PTHS(:,:,:) = PTHS(:,:,:) - PRCS(:,:,:)*ZLV(:,:,:)/(ZCPH(:,:,:)*ZEXNS(:,:,:))
    PRCS(:,:,:) = 0.0
@@ -696,8 +696,8 @@ IF ( KRR .GE. 6 ) PRS(:,:,:,6) = PRGS(:,:,:)
 !
 ! Prepare 3D number concentrations
 !
-IF ( LWARM ) PSVS(:,:,:,NSV_LIMA_NC) = PCCS(:,:,:)
-! IF ( LCOLD ) PSVS(:,:,:,NSV_LIMA_NI) = PCIS(:,:,:)
+IF ( LWARM .AND. NMOM_C.GE.2 ) PSVS(:,:,:,NSV_LIMA_NC) = PCCS(:,:,:)
+! IF ( LCOLD .AND. NMOM_I.GE.2 ) PSVS(:,:,:,NSV_LIMA_NI) = PCIS(:,:,:)
 !
 IF ( LSCAV .AND. LAERO_MASS ) PSVS(:,:,:,NSV_LIMA_SCAVMASS) = PMAS(:,:,:)
 ! 
diff --git a/src/MNH/lima_droplets_autoconversion.f90 b/src/MNH/lima_droplets_autoconversion.f90
index 270906662..3fa32e7a6 100644
--- a/src/MNH/lima_droplets_autoconversion.f90
+++ b/src/MNH/lima_droplets_autoconversion.f90
@@ -59,7 +59,7 @@ END MODULE MODI_LIMA_DROPLETS_AUTOCONVERSION
 !*       0.    DECLARATIONS
 !              ------------
 !
-USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN, LKHKO
+USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN, LKHKO, LKESSLERAC, NMOM_C
 USE MODD_PARAM_LIMA_WARM, ONLY : XLAUTR, XAUTO1, XLAUTR_THRESHOLD, &
                                  XITAUTR, XAUTO2, XITAUTR_THRESHOLD, &
                                  XACCR4, XACCR5, XACCR3, XACCR1, XAC, XR0
@@ -96,7 +96,9 @@ ZW3(:) = 0.0
 ZW2(:) = 0.0
 ZW1(:) = 0.0
 !
-IF (LKHKO) THEN
+IF (NMOM_C.EQ.1 .AND. LKESSLERAC) THEN
+   P_RC_AUTO(:) = - 1.E-3 * MAX ( PRCT(:) - 0.5E-3 / PRHODREF(:), 0. )
+ELSE IF (LKHKO) THEN
 !
 !        1. Autoconversion of cloud droplets (Berry-Reinhardt parameterization)
 !   	 ----------------------------------------------------------------------
diff --git a/src/MNH/lima_droplets_hom_freezing.f90 b/src/MNH/lima_droplets_hom_freezing.f90
index db27f466e..58e12b5f1 100644
--- a/src/MNH/lima_droplets_hom_freezing.f90
+++ b/src/MNH/lima_droplets_hom_freezing.f90
@@ -10,8 +10,7 @@ INTERFACE
    SUBROUTINE LIMA_DROPLETS_HOM_FREEZING (PTSTEP, LDCOMPUTE,                &
                                           PT, PLVFACT, PLSFACT,             &
                                           PRCT, PCCT, PLBDC,                &
-                                          P_TH_HONC, P_RC_HONC, P_CC_HONC,  &
-                                          PA_TH, PA_RC, PA_CC, PA_RI, PA_CI )
+                                          P_TH_HONC, P_RC_HONC, P_CC_HONC   )
 !
 REAL,                 INTENT(IN)    :: PTSTEP 
 LOGICAL, DIMENSION(:),INTENT(IN)    :: LDCOMPUTE
@@ -28,12 +27,6 @@ REAL, DIMENSION(:),   INTENT(INOUT) :: P_TH_HONC
 REAL, DIMENSION(:),   INTENT(INOUT) :: P_RC_HONC
 REAL, DIMENSION(:),   INTENT(INOUT) :: P_CC_HONC
 !
-REAL, DIMENSION(:),   INTENT(INOUT) :: PA_TH
-REAL, DIMENSION(:),   INTENT(INOUT) :: PA_RC
-REAL, DIMENSION(:),   INTENT(INOUT) :: PA_CC
-REAL, DIMENSION(:),   INTENT(INOUT) :: PA_RI
-REAL, DIMENSION(:),   INTENT(INOUT) :: PA_CI
-!
 END SUBROUTINE LIMA_DROPLETS_HOM_FREEZING
 END INTERFACE
 END MODULE MODI_LIMA_DROPLETS_HOM_FREEZING
@@ -42,8 +35,7 @@ END MODULE MODI_LIMA_DROPLETS_HOM_FREEZING
       SUBROUTINE LIMA_DROPLETS_HOM_FREEZING (PTSTEP,  LDCOMPUTE,               &
                                              PT, PLVFACT, PLSFACT,             &
                                              PRCT, PCCT, PLBDC,                &
-                                             P_TH_HONC, P_RC_HONC, P_CC_HONC,  &
-                                             PA_TH, PA_RC, PA_CC, PA_RI, PA_CI )
+                                             P_TH_HONC, P_RC_HONC, P_CC_HONC   )
 !     ##########################################################################
 !
 !!    PURPOSE
@@ -91,12 +83,6 @@ REAL, DIMENSION(:),   INTENT(INOUT) :: P_TH_HONC
 REAL, DIMENSION(:),   INTENT(INOUT) :: P_RC_HONC
 REAL, DIMENSION(:),   INTENT(INOUT) :: P_CC_HONC
 !
-REAL, DIMENSION(:),   INTENT(INOUT) :: PA_TH
-REAL, DIMENSION(:),   INTENT(INOUT) :: PA_RC
-REAL, DIMENSION(:),   INTENT(INOUT) :: PA_CC
-REAL, DIMENSION(:),   INTENT(INOUT) :: PA_RI
-REAL, DIMENSION(:),   INTENT(INOUT) :: PA_CI
-!
 !*       0.2   Declarations of local variables :
 !
 REAL, DIMENSION(SIZE(PT)) ::  ZZW, ZZX, ZZY, ZTCELSIUS
@@ -130,13 +116,7 @@ WHERE ( (PT(:)<XTT-35.0) .AND. (PCCT(:)>XCTMIN(2)) .AND. (PRCT(:)>XRTMIN(2)) )
 !
    P_RC_HONC(:) = - ZZY(:)/PTSTEP
    P_CC_HONC(:) = - ZZW(:)/PTSTEP
-   P_TH_HONC(:) = P_RC_HONC(:) * (PLSFACT(:)-PLVFACT(:))
-!
-   PA_TH(:) = PA_TH(:) + P_TH_HONC(:)
-   PA_RC(:) = PA_RC(:) + P_RC_HONC(:)
-   PA_CC(:) = PA_CC(:) + P_CC_HONC(:)
-   PA_RI(:) = PA_RI(:) - P_RC_HONC(:)
-   PA_CI(:) = PA_CI(:) - P_CC_HONC(:)
+!   P_TH_HONC(:) = P_RC_HONC(:) * (PLSFACT(:)-PLVFACT(:))
 !
 END WHERE
 !
diff --git a/src/MNH/lima_tendencies.f90 b/src/MNH/lima_tendencies.f90
index 4977e2109..ea1068aaa 100644
--- a/src/MNH/lima_tendencies.f90
+++ b/src/MNH/lima_tendencies.f90
@@ -606,11 +606,16 @@ IF (LCOLD .AND. LWARM) THEN
    CALL LIMA_DROPLETS_HOM_FREEZING (PTSTEP, LDCOMPUTE,                 & ! independent from CF,IF,PF
                                     ZT, ZLVFACT, ZLSFACT,              &
                                     PRCT, PCCT, ZLBDC,                 &
-                                    P_TH_HONC, P_RC_HONC, P_CC_HONC,   &
-                                    PA_TH, PA_RC, PA_CC, PA_RI, PA_CI  )
+                                    P_TH_HONC, P_RC_HONC, P_CC_HONC    )
+   PA_RC(:) = PA_RC(:) + P_RC_HONC(:)
+   PA_CC(:) = PA_CC(:) + P_CC_HONC(:)
+   PA_RI(:) = PA_RI(:) - P_RC_HONC(:)
+   PA_CI(:) = PA_CI(:) - P_CC_HONC(:) 
+   P_TH_ACC(:) = - P_CC_HONC(:) * (ZLSFACT(:)-ZLVFACT(:))
+   PA_TH(:) = PA_TH(:) + P_TH_HONC(:)
 END IF
 !
-IF (LWARM .AND. LRAIN .AND. (.NOT. LKHKO)) THEN
+IF (LWARM .AND. LRAIN .AND. (.NOT. LKHKO) .AND. NMOM_C.GE.2) THEN
    CALL LIMA_DROPLETS_SELF_COLLECTION (LDCOMPUTE,          & ! depends on CF
                                        PRHODREF,           &
                                        PCCT/ZCF1D, ZLBDC3, &
@@ -649,7 +654,7 @@ IF (LWARM .AND. LRAIN) THEN
    PA_RR(:) = PA_RR(:) - P_RC_ACCR(:)
 END IF
 !
-IF (LWARM .AND. LRAIN .AND. (.NOT. LKHKO)) THEN 
+IF (LWARM .AND. LRAIN .AND. (.NOT. LKHKO) .AND. NMOM_R.GE.2) THEN 
    CALL LIMA_DROPS_SELF_COLLECTION (LDCOMPUTE,           & ! depends on PF
                                     PRHODREF,            &
                                     PCRT/ZPF1D(:), ZLBDR, ZLBDR3, &
diff --git a/src/MNH/modd_param_lima.f90 b/src/MNH/modd_param_lima.f90
index e3c53475e..d199c2a66 100644
--- a/src/MNH/modd_param_lima.f90
+++ b/src/MNH/modd_param_lima.f90
@@ -148,6 +148,7 @@ LOGICAL, SAVE :: LACTTKE       ! TRUE to take into account TKE in W for activati
 LOGICAL, SAVE :: LADJ          ! TRUE for adjustment procedure + Smax (false for diagnostic supersaturation)
 LOGICAL, SAVE :: LSPRO         ! TRUE for prognostic supersaturation                     
 LOGICAL, SAVE :: LKHKO         ! TRUE for Scu simulation (replicates the previous KHKO scheme)                     
+LOGICAL, SAVE :: LKESSLERAC    ! TRUE for Kessler autoconversion (if NMOM_C=1)
 !
 INTEGER, SAVE :: NMOM_C        ! Number of moments for cloud droplets
 INTEGER, SAVE :: NMOM_R        ! Number of moments for rain drops
diff --git a/src/MNH/modn_param_lima.f90 b/src/MNH/modn_param_lima.f90
index 6edddc884..44d467145 100644
--- a/src/MNH/modn_param_lima.f90
+++ b/src/MNH/modn_param_lima.f90
@@ -26,7 +26,7 @@ NAMELIST/NAM_PARAM_LIMA/LCOLD, LNUCL, LSEDI, LSNOW, LHAIL, LHHONI, LMEYERS,&
                         XFACTNUC_DEP, XFACTNUC_CON, NPHILLIPS,             &
                         LCIBU, XNDEBRIS_CIBU, LRDSF, LMURAKAMI,            &
                         LWARM, LACTI, LRAIN, LSEDC, LACTIT, LBOUND, LSPRO, &
-                        LADJ, LKHKO, NMOM_C, NMOM_R,                       &
+                        LADJ, LKHKO, LKESSLERAC, NMOM_C, NMOM_R,           &
                         NMOD_CCN, XCCN_CONC,                               &
                         LCCN_HOM, CCCN_MODES, HINI_CCN, HTYPE_CCN,         &
                         XALPHAC, XNUC, XALPHAR, XNUR,                      &
diff --git a/src/MNH/set_conc_lima.f90 b/src/MNH/set_conc_lima.f90
index 3eeda7f5a..c7dce39bd 100644
--- a/src/MNH/set_conc_lima.f90
+++ b/src/MNH/set_conc_lima.f90
@@ -73,7 +73,8 @@ contains
 !*       0.    DECLARATIONS
 !              ------------
 !
-USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN, LCOLD, LWARM, LRAIN, NMOD_CCN, NMOD_IFN
+USE MODD_PARAM_LIMA,      ONLY : XRTMIN, XCTMIN, LCOLD, LWARM, LRAIN, NMOD_CCN, NMOD_IFN, &
+                                 NMOM_C, NMOM_R, NMOM_I
 USE MODD_PARAM_LIMA_COLD, ONLY : XAI, XBI
 USE MODD_NSV,             ONLY : NSV_LIMA_BEG_A, NSV_LIMA_NC_A, NSV_LIMA_NR_A, NSV_LIMA_CCN_ACTI_A, &
                                  NSV_LIMA_NI_A, NSV_LIMA_IFN_NUCL_A
@@ -110,7 +111,7 @@ ILUOUT = TLUOUT%NLU
 !*       2.    INITIALIZATION
 !              --------------
 !
-IF (LWARM .AND. NRR.GE.2) THEN
+IF (LWARM .AND. NRR.GE.2 .AND. NMOM_C.GE.2) THEN
 !
 !  droplets
 !
@@ -138,7 +139,7 @@ IF (LWARM .AND. NRR.GE.2) THEN
    END IF
 END IF
 !
-IF (LWARM .AND. LRAIN .AND. NRR.GE.3) THEN
+IF (LWARM .AND. LRAIN .AND. NRR.GE.3 .AND. NMOM_R.GE.2) THEN
 !
 !  drops
 !
@@ -161,7 +162,7 @@ IF (LWARM .AND. LRAIN .AND. NRR.GE.3) THEN
    END IF
 END IF
 !
-IF (LCOLD .AND. NRR.GE.4) THEN
+IF (LCOLD .AND. NRR.GE.4 .AND. NMOM_I.GE.2) THEN
 !
 ! ice crystals
 !
diff --git a/src/MNH/sources_neg_correct.f90 b/src/MNH/sources_neg_correct.f90
index 81d49856c..c17d1dbc1 100644
--- a/src/MNH/sources_neg_correct.f90
+++ b/src/MNH/sources_neg_correct.f90
@@ -235,7 +235,7 @@ CLOUD: select case ( hcloud )
     end where
 !
 !
-  case( 'LIMA' )
+  case( 'LIMA_OFF' )
 ! Correction where rc<0 or Nc<0
     if ( lwarm_lima ) then
       where ( prrs(:, :, :, 2) < xrtmin_lima(2) / ptstep .or. prsvs(:, :, :, nsv_lima_nc) < xctmin_lima(2) / ptstep )
-- 
GitLab