From 3c0e7ee9b0862c6ebd735cd77ae409d652e5b59b Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Beno=C3=AEt=20Vi=C3=A9?= <benoit.vie@meteo.fr>
Date: Wed, 30 Mar 2022 17:21:15 +0200
Subject: [PATCH] Add LSNOW_T switch in ICE3
---
src/MNH/ice4_fast_rg.f90 | 10 +-
src/MNH/ice4_fast_rh.f90 | 11 +-
src/MNH/ice4_fast_rs.f90 | 38 ++-----
src/MNH/ice4_rsrimcg_old.f90 | 5 +-
src/MNH/ice4_sedimentation_split.f90 | 54 +++++-----
src/MNH/ice4_sedimentation_stat.f90 | 34 +++----
src/MNH/ice4_slow.f90 | 16 +--
src/MNH/ice4_tendencies.f90 | 70 +++++++------
src/MNH/ini_rain_ice.f90 | 146 ++++++++++++++-------------
src/MNH/ini_rain_ice_elec.f90 | 140 +++++++++++++++----------
src/MNH/modd_rain_ice_descr.f90 | 6 +-
src/MNH/modd_rain_ice_param.f90 | 3 -
src/MNH/rain_ice_red.f90 | 56 +++++-----
13 files changed, 297 insertions(+), 292 deletions(-)
diff --git a/src/MNH/ice4_fast_rg.f90 b/src/MNH/ice4_fast_rg.f90
index c78a154bb..0f7dea217 100644
--- a/src/MNH/ice4_fast_rg.f90
+++ b/src/MNH/ice4_fast_rg.f90
@@ -94,6 +94,7 @@ SUBROUTINE ICE4_FAST_RG(KSIZE, LDSOFT, PCOMPUTE, KRR, &
!!
! P. Wautelet 26/04/2019: replace non-standard FLOAT function by REAL function
! P. Wautelet 29/05/2019: remove PACK/UNPACK intrinsics (to get more performance and better OpenACC support)
+! J. Wurtz 03/2022: New snow characteristics with LSNOW_T
!
!
!* 0. DECLARATIONS
@@ -335,14 +336,7 @@ ELSE
END DO
!
WHERE(GDRY(:))
-! PRG_TEND(:, IRSWETG)=XFSDRYG*ZZW(:) & ! RSDRYG
-! / XCOLSG &
-! *(PLBDAS(:)**(XCXS-XBS))*( PLBDAG(:)**XCXG ) &
-! *(PRHODREF(:)**(-XCEXVT-1.)) &
-! *( XLBSDRYG1/( PLBDAG(:)**2 ) + &
-! XLBSDRYG2/( PLBDAG(:) * PLBDAS(:) ) + &
-! XLBSDRYG3/( PLBDAS(:)**2))
- PRG_TEND(:, IRSWETG)=XFSDRYG*ZZW(:) & ! RSDRYG ! Modif Wurtz snow diag
+ PRG_TEND(:, IRSWETG)=XFSDRYG*ZZW(:) & ! RSDRYG
/ XCOLSG &
*(PRST(:))*( PLBDAG(:)**XCXG ) &
*(PRHODREF(:)**(-XCEXVT)) &
diff --git a/src/MNH/ice4_fast_rh.f90 b/src/MNH/ice4_fast_rh.f90
index 98106308d..a21c85f19 100644
--- a/src/MNH/ice4_fast_rh.f90
+++ b/src/MNH/ice4_fast_rh.f90
@@ -84,6 +84,7 @@ SUBROUTINE ICE4_FAST_RH(KSIZE, LDSOFT, PCOMPUTE, PWETG, &
!!
! P. Wautelet 26/04/2019: replace non-standard FLOAT function by REAL function
! P. Wautelet 29/05/2019: remove PACK/UNPACK intrinsics (to get more performance and better OpenACC support)
+! J. Wurtz 03/2022: New snow characteristics with LSNOW_T
!
!
!* 0. DECLARATIONS
@@ -269,16 +270,10 @@ ELSE
END DO
!
WHERE(GWET(:))
-! PRH_TEND(:, IRSWETH)=XFSWETH*ZZW(:) & ! RSWETH
-! *( PLBDAS(:)**(XCXS-XBS) )*( PLBDAH(:)**XCXH ) &
-! *( PRHODREF(:)**(-XCEXVT-1.) ) &
-! *( XLBSWETH1/( PLBDAH(:)**2 ) + &
-! XLBSWETH2/( PLBDAH(:) * PLBDAS(:) ) + &
-! XLBSWETH3/( PLBDAS(:)**2) )
- PRH_TEND(:, IRSWETH)=XFSWETH*ZZW(:) & ! RSWETH Modif Wurtz conc snow
+ PRH_TEND(:, IRSWETH)=XFSWETH*ZZW(:) &
*( PRST(:))*( PLBDAH(:)**XCXH ) &
*( PRHODREF(:)**(-XCEXVT) ) &
- *( XLBSWETH1/( PLBDAH(:)**2 ) + & ! Il s'agit de moment (?)
+ *( XLBSWETH1/( PLBDAH(:)**2 ) + &
XLBSWETH2/( PLBDAH(:) * PLBDAS(:) ) + &
XLBSWETH3/( PLBDAS(:)**2) )
PRH_TEND(:, IRSDRYH)=PRH_TEND(:, IRSWETH)*(XCOLSH*EXP(XCOLEXSH*(PT(:)-XTT)))
diff --git a/src/MNH/ice4_fast_rs.f90 b/src/MNH/ice4_fast_rs.f90
index b3df38c8c..50d0cae7d 100644
--- a/src/MNH/ice4_fast_rs.f90
+++ b/src/MNH/ice4_fast_rs.f90
@@ -77,6 +77,7 @@ SUBROUTINE ICE4_FAST_RS(KSIZE, LDSOFT, PCOMPUTE, &
!!
! P. Wautelet 26/04/2019: replace non-standard FLOAT function by REAL function
! P. Wautelet 29/05/2019: remove PACK/UNPACK intrinsics (to get more performance and better OpenACC support)
+! J. Wurtz 03/2022: New snow characteristics with LSNOW_T
!
!
!* 0. DECLARATIONS
@@ -244,13 +245,10 @@ ELSE
! 5.1.4 riming of the small sized aggregates
!
WHERE (GRIM(:))
-! PRS_TEND(:, IRCRIMSS) = XCRIMSS * ZZW(:) * PRCT(:) & ! RCRIMSS
-! * PLBDAS(:)**XEXCRIMSS &
-! * PRHODREF(:)**(-XCEXVT)
- PRS_TEND(:, IRCRIMSS) = XCRIMSS * ZZW(:) * PRCT(:) & ! RCRIMSS !Wurtz ! Thompson
+ PRS_TEND(:, IRCRIMSS) = XCRIMSS * ZZW(:) * PRCT(:) & ! RCRIMSS
* PRST(:)*(1+(XFVELOS/PLBDAS(:))**XALPHAS)**(-XNUS+XEXCRIMSS/XALPHAS) &
* PRHODREF(:)**(-XCEXVT+1.) &
- * (PLBDAS(:)) ** (XEXCRIMSS+XBS) ! Thompson
+ * (PLBDAS(:)) ** (XEXCRIMSS+XBS)
END WHERE
!
! 5.1.5 perform the linear interpolation of the normalized
@@ -275,21 +273,17 @@ ELSE
!
!
WHERE(GRIM(:))
-! PRS_TEND(:, IRCRIMS)=XCRIMSG * PRCT(:) & ! RCRIMS
-! * PLBDAS(:)**XEXCRIMSG &
-! * PRHODREF(:)**(-XCEXVT)
PRS_TEND(:, IRCRIMS) = XCRIMSG * PRCT(:) & ! RCRIMS
* PRST(:)*(1+(XFVELOS/PLBDAS(:))**(XALPHAS))**(-XNUS+XEXCRIMSG/XALPHAS) &
* PRHODREF(:)**(-XCEXVT+1.) &
- * PLBDAS(:)**(XBS+XEXCRIMSG) ! GAMMAGEN LH_EXTENDED
+ * PLBDAS(:)**(XBS+XEXCRIMSG)
ZZW6(:) = PRS_TEND(:, IRCRIMS) - PRS_TEND(:, IRCRIMSS) ! RCRIMSG
END WHERE
IF(CSNOWRIMING=='M90 ')THEN
!Murakami 1990
WHERE(GRIM(:))
- ! PRS_TEND(:, IRSRIMCG)=XSRIMCG * PLBDAS(:)**XEXSRIMCG*(1.0-ZZW(:))
- PRS_TEND(:, IRSRIMCG)=XSRIMCG * PRST(:)*PRHODREF(:)*PLBDAS(:)**(XEXSRIMCG+XBS)*(1.0-ZZW(:)) !Wurtz
+ PRS_TEND(:, IRSRIMCG)=XSRIMCG * PRST(:)*PRHODREF(:)*PLBDAS(:)**(XEXSRIMCG+XBS)*(1.0-ZZW(:))
PRS_TEND(:, IRSRIMCG)=ZZW6(:)*PRS_TEND(:, IRSRIMCG)/ &
MAX(1.E-20, &
@@ -394,13 +388,8 @@ ELSE
! 5.2.4 raindrop accretion on the small sized aggregates
!
WHERE(GACC(:))
-! ZZW6(:) = & !! coef of RRACCS
-! XFRACCSS*( PLBDAS(:)**XCXS )*( PRHODREF(:)**(-XCEXVT-1.) ) &
-! *( XLBRACCS1/((PLBDAS(:)**2) ) + &
-! XLBRACCS2/( PLBDAS(:) * PLBDAR(:) ) + &
-! XLBRACCS3/( (PLBDAR(:)**2)) )/PLBDAR(:)**4
ZZW6(:) = & !! coef of RRACCS
- XFRACCSS*( PRST(:)*PLBDAS(:)**XBS )*( PRHODREF(:)**(-XCEXVT) ) & ! Wurtz
+ XFRACCSS*( PRST(:)*PLBDAS(:)**XBS )*( PRHODREF(:)**(-XCEXVT) ) &
*( XLBRACCS1/((PLBDAS(:)**2) ) + &
XLBRACCS2/( PLBDAS(:) * PLBDAR(:) ) + &
XLBRACCS3/( (PLBDAR(:)**2)) )/PLBDAR(:)**4
@@ -445,12 +434,7 @@ ELSE
! into graupeln
!
WHERE(GACC(:))
-! PRS_TEND(:, IRSACCRG) = XFSACCRG*ZZW(:)* & ! RSACCRG
-! ( PLBDAS(:)**(XCXS-XBS) )*( PRHODREF(:)**(-XCEXVT-1.) ) &
-! *( XLBSACCR1/((PLBDAR(:)**2) ) + &
-! XLBSACCR2/( PLBDAR(:) * PLBDAS(:) ) + &
-! XLBSACCR3/( (PLBDAS(:)**2)) )/PLBDAR(:)
- PRS_TEND(:, IRSACCRG) = XFSACCRG*ZZW(:)* & ! RSACCRG ! Modif Wurtz
+ PRS_TEND(:, IRSACCRG) = XFSACCRG*ZZW(:)* & ! RSACCRG
( PRST(:))*( PRHODREF(:)**(-XCEXVT) ) &
*( XLBSACCR1/((PLBDAR(:)**2) ) + &
XLBSACCR2/( PLBDAR(:) * PLBDAS(:) ) + &
@@ -516,13 +500,7 @@ ELSE
!
! compute RSMLT
!
-! PRSMLTG(:) = XFSCVMG*MAX( 0.0,( -PRSMLTG(:) * &
-! ( X0DEPS* PLBDAS(:)**XEX0DEPS + &
-! X1DEPS*PCJ(:)*PLBDAS(:)**XEX1DEPS ) - &
-! ( PRS_TEND(:, IRCRIMS) + PRS_TEND(:, IRRACCS) ) * &
-! ( PRHODREF(:)*XCL*(XTT-PT(:))) ) / &
-! ( PRHODREF(:)*XLMTT ) )
- PRSMLTG(:) = XFSCVMG*MAX( 0.0,( -PRSMLTG(:) * & ! Modif GAMMAGEN LH_EXTENDED
+ PRSMLTG(:) = XFSCVMG*MAX( 0.0,( -PRSMLTG(:) * &
PRST(:)*PRHODREF(:) * &
( X0DEPS *PLBDAS(:)**XEX0DEPS + &
X1DEPS*PCJ(:)*(1+(XFVELOS/(2*PLBDAS(:))**XALPHAS))**(XNUS+XEX1DEPS/XALPHAS)*((PLBDAS(:))**(XBS+XEX1DEPS))) - &
diff --git a/src/MNH/ice4_rsrimcg_old.f90 b/src/MNH/ice4_rsrimcg_old.f90
index 5fc932536..d250ddc0e 100644
--- a/src/MNH/ice4_rsrimcg_old.f90
+++ b/src/MNH/ice4_rsrimcg_old.f90
@@ -44,6 +44,7 @@ SUBROUTINE ICE4_RSRIMCG_OLD(KSIZE, ODSOFT, ODCOMPUTE, &
!!
! P. Wautelet 26/04/2019: replace non-standard FLOAT function by REAL function
! P. Wautelet 29/05/2019: remove PACK/UNPACK intrinsics (to get more performance and better OpenACC support)
+! J. Wurtz 03/2022: New snow characteristics with LSNOW_T
!
!
!* 0. DECLARATIONS
@@ -131,9 +132,7 @@ IF(.NOT. ODSOFT) THEN
!
!
WHERE(GRIM(:))
-! PRSRIMCG_MR(:) = XSRIMCG * PLBDAS(:)**XEXSRIMCG & ! RSRIMCG
-! * (1.0 - ZZW(:) )/PRHODREF(:)
- PRSRIMCG_MR(:) = XSRIMCG * PLBDAS(:)**XEXSRIMCG & ! RSRIMCG !Modif Wurtz concentration snow
+ PRSRIMCG_MR(:) = XSRIMCG * PLBDAS(:)**XEXSRIMCG & ! RSRIMCG
* (1.0 - ZZW(:) )*PRST(:)
PRSRIMCG_MR(:)=MIN(PRST(:), PRSRIMCG_MR(:))
END WHERE
diff --git a/src/MNH/ice4_sedimentation_split.f90 b/src/MNH/ice4_sedimentation_split.f90
index 9e93d59e5..f9369a878 100644
--- a/src/MNH/ice4_sedimentation_split.f90
+++ b/src/MNH/ice4_sedimentation_split.f90
@@ -6,12 +6,12 @@
MODULE MODI_ICE4_SEDIMENTATION_SPLIT
INTERFACE
SUBROUTINE ICE4_SEDIMENTATION_SPLIT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKTB, KKTE, KKT, KKL, &
- &PTSTEP, KRR, OSEDIC, ODEPOSC, PVDEPOSC, PDZZ, &
- &PRHODREF, PPABST, PTHT, PRHODJ, &
- & PLBDAS, & ! Modif Wurtz
- &PRCS, PRCT, PRRS, PRRT, PRIS, PRIT, PRSS, PRST, PRGS, PRGT,&
- &PINPRC, PINDEP, PINPRR, PINPRI, PINPRS, PINPRG, &
- &PSEA, PTOWN, &
+ &PTSTEP, KRR, OSEDIC, ODEPOSC, PVDEPOSC, PDZZ, &
+ &PRHODREF, PPABST, PTHT, PRHODJ, &
+ & PLBDAS, &
+ &PRCS, PRCT, PRRS, PRRT, PRIS, PRIT, PRSS, PRST, PRGS, PRGT, &
+ &PINPRC, PINDEP, PINPRR, PINPRI, PINPRS, PINPRG, &
+ &PSEA, PTOWN, &
&PINPRH, PRHT, PRHS, PFPR)
IMPLICIT NONE
INTEGER, INTENT(IN) :: KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKTB, KKTE, KKT
@@ -26,7 +26,7 @@ REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRHODREF! Reference den
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PPABST ! absolute pressure at t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PTHT ! Theta at time t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRHODJ ! Dry density * Jacobian
-REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PLBDAS ! lambda parameter for snow ! Modif Wurtz
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PLBDAS ! lambda parameter for snow
REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PRCS ! Cloud water m.r. source
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRCT ! Cloud water m.r. at t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PRRS ! Rain water m.r. source
@@ -53,12 +53,12 @@ END SUBROUTINE ICE4_SEDIMENTATION_SPLIT
END INTERFACE
END MODULE MODI_ICE4_SEDIMENTATION_SPLIT
SUBROUTINE ICE4_SEDIMENTATION_SPLIT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKTB, KKTE, KKT, KKL, &
- &PTSTEP, KRR, OSEDIC, ODEPOSC, PVDEPOSC, PDZZ, &
- &PRHODREF, PPABST, PTHT, PRHODJ, &
- & PLBDAS, & ! Modif Wurtz
- &PRCS, PRCT, PRRS, PRRT, PRIS, PRIT, PRSS, PRST, PRGS, PRGT,&
- &PINPRC, PINDEP, PINPRR, PINPRI, PINPRS, PINPRG, &
- &PSEA, PTOWN, &
+ &PTSTEP, KRR, OSEDIC, ODEPOSC, PVDEPOSC, PDZZ, &
+ &PRHODREF, PPABST, PTHT, PRHODJ, &
+ & PLBDAS, &
+ &PRCS, PRCT, PRRS, PRRT, PRIS, PRIT, PRSS, PRST, PRGS, PRGT, &
+ &PINPRC, PINDEP, PINPRR, PINPRI, PINPRS, PINPRG, &
+ &PSEA, PTOWN, &
&PINPRH, PRHT, PRHS, PFPR)
!!
!!** PURPOSE
@@ -75,6 +75,7 @@ SUBROUTINE ICE4_SEDIMENTATION_SPLIT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKTB
!!
! P. Wautelet 11/02/2019: dimensions of PINPRC and PINDEP not necessarily KIT,KJT
! P. Wautelet 10/04/2019: replace ABORT and STOP calls by Print_msg
+! J. Wurtz 03/2022: New snow characteristics with LSNOW_T
!
!
!* 0. DECLARATIONS
@@ -105,7 +106,7 @@ REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRHODREF! Reference den
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PPABST ! absolute pressure at t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PTHT ! Theta at time t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRHODJ ! Dry density * Jacobian
-REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PLBDAS ! lambda parameter for snow ! Modif Wurtz
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PLBDAS ! lambda parameter for snow
REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PRCS ! Cloud water m.r. source
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRCT ! Cloud water m.r. at t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PRRS ! Rain water m.r. source
@@ -243,7 +244,7 @@ IF (GSEDIC) THEN
&XSPLIT_MAXCFL, &
&PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
&2, &
- &PLBDAS, & ! Modif Wurtz
+ &PLBDAS, &
&ZRCT, PRCS, PINPRC, ZPRCS, &
&ZRAY, ZLBC, ZFSEDC, ZCONC3D, PFPR=PFPR)
ENDIF
@@ -270,7 +271,7 @@ END IF
&XSPLIT_MAXCFL, &
&PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
&3, &
- &PLBDAS, & ! Modif Wurtz
+ &PLBDAS, &
&ZRRT, PRRS, PINPRR, ZPRRS, &
&PFPR=PFPR)
!
@@ -280,7 +281,7 @@ END IF
&XSPLIT_MAXCFL, &
&PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
&4, &
- &PLBDAS, & ! Modif Wurtz
+ &PLBDAS, &
&ZRIT, PRIS, PINPRI, ZPRIS, &
PFPR=PFPR)
!
@@ -290,7 +291,7 @@ END IF
&XSPLIT_MAXCFL, &
&PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
&5, &
- &PLBDAS, & ! Modif Wurtz
+ &PLBDAS, &
&ZRST, PRSS, PINPRS, ZPRSS, &
PFPR=PFPR)
!
@@ -300,7 +301,7 @@ END IF
&XSPLIT_MAXCFL, &
&PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
&6, &
- &PLBDAS, & ! Modif Wurtz
+ &PLBDAS, &
&ZRGT, PRGS, PINPRG, ZPRGS, &
PFPR=PFPR)
!
@@ -311,7 +312,7 @@ IF (IRR==7) THEN
&XSPLIT_MAXCFL, &
&PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
&7, &
- &PLBDAS, & ! Modif Wurtz
+ &PLBDAS, &
&ZRHT, PRHS, PINPRH, ZPRHS, &
PFPR=PFPR)
ENDIF
@@ -332,9 +333,9 @@ SUBROUTINE INTERNAL_SEDIM_SPLI(KIB,KIE,KIT,KJB,KJE,KJT,KKB,KKTB,KKTE,KKT,KKL,KRR
! ------------
!
USE MODD_CST, ONLY: XCPD,XP00,XRD
-USE MODD_RAIN_ICE_DESCR, ONLY: XCC,XCEXVT,XDC,XLBEXC,XRTMIN,XALPHAS,XNUS,XBS
-USE MODD_RAIN_ICE_PARAM, ONLY: XEXCSEDI,XEXSEDG,XEXSEDH,XEXSEDR,XEXSEDS,XFSEDG,XFSEDH,XFSEDI,XFSEDR,XFSEDS,&
- XFVELOS
+USE MODD_RAIN_ICE_DESCR, ONLY: XCC,XCEXVT,XDC,XLBEXC,XRTMIN,XALPHAS,XNUS,XBS,XFVELOS
+USE MODD_RAIN_ICE_PARAM, ONLY: XEXCSEDI,XEXSEDG,XEXSEDH,XEXSEDR,XEXSEDS,XFSEDG,XFSEDH,XFSEDI,XFSEDR,XFSEDS
+
!
IMPLICIT NONE
!
@@ -456,8 +457,8 @@ DO WHILE (ANY(ZREMAINT>0.))
ZWSED(JI, JJ, JK) = XFSEDS * &
& PRXT(JI,JJ,JK)* &
- & PRHODREF(JI,JJ,JK)**(1-XCEXVT) * & ! Modif Wurtz snow
- & (1 + (XFVELOS/PLBDAS(JI, JJ, JK))**XALPHAS)** (-XNUS+XEXSEDS/XALPHAS) * & ! GAMMAGEN LH_EXTENDED
+ & PRHODREF(JI,JJ,JK)**(1-XCEXVT) * &
+ & (1 + (XFVELOS/PLBDAS(JI, JJ, JK))**XALPHAS)** (-XNUS+XEXSEDS/XALPHAS) * &
& PLBDAS(JI, JJ, JK) ** (XBS+XEXSEDS)
! GAMMAGEN_LH_EXTENDED
@@ -469,9 +470,6 @@ DO WHILE (ANY(ZREMAINT>0.))
CASE(3)
ZFSED=XFSEDR
ZEXSED=XEXSEDR
-! CASE(5)
-! ZFSED=XFSEDS
-! ZEXSED=XEXSEDS
CASE(6)
ZFSED=XFSEDG
ZEXSED=XEXSEDG
diff --git a/src/MNH/ice4_sedimentation_stat.f90 b/src/MNH/ice4_sedimentation_stat.f90
index 4daa043da..891223920 100644
--- a/src/MNH/ice4_sedimentation_stat.f90
+++ b/src/MNH/ice4_sedimentation_stat.f90
@@ -8,7 +8,7 @@ INTERFACE
SUBROUTINE ICE4_SEDIMENTATION_STAT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKTB, KKTE, KKT, KKL, &
&PTSTEP, KRR, OSEDIC, ODEPOSC, PVDEPOSC, PDZZ, &
&PRHODREF, PPABST, PTHT, PRHODJ, &
- & PLBDAS, & ! Modif Wurtz
+ & PLBDAS, &
&PRCS, PRCT, PRRS, PRRT, PRIS, PRIT,&
&PRSS, PRST, PRGS, PRGT,&
&PINPRC, PINDEP, PINPRR, PINPRI, PINPRS, PINPRG, &
@@ -27,7 +27,7 @@ REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRHODREF! Reference den
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PPABST ! absolute pressure at t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PTHT ! Theta at time t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRHODJ ! Dry density * Jacobian
-REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PLBDAS ! lambda parameter for snow ! Modif Wurtz
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PLBDAS ! lambda parameter for snow
REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PRCS ! Cloud water m.r. source
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRCT ! Cloud water m.r. at t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PRRS ! Rain water m.r. source
@@ -56,7 +56,7 @@ END MODULE MODI_ICE4_SEDIMENTATION_STAT
SUBROUTINE ICE4_SEDIMENTATION_STAT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKTB, KKTE, KKT, KKL, &
&PTSTEP, KRR, OSEDIC, ODEPOSC, PVDEPOSC, PDZZ, &
&PRHODREF, PPABST, PTHT, PRHODJ, &
- & PLBDAS, & ! Modif Wurtz
+ & PLBDAS, &
&PRCS, PRCT, PRRS, PRRT, PRIS, PRIT, &
&PRSS, PRST, PRGS, PRGT,&
&PINPRC, PINDEP, PINPRR, PINPRI, PINPRS, PINPRG, &
@@ -78,6 +78,7 @@ SUBROUTINE ICE4_SEDIMENTATION_STAT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKTB,
! P. Wautelet 10/04/2019: replace ABORT and STOP calls by Print_msg
! P. Wautelet 28/05/2019: move COUNTJV function to tools.f90
! P. Wautelet 21/01/2021: initialize untouched part of PFPR
+! J. Wurtz 03/2022: New snow characteristics with LSNOW_T
!
!
!* 0. DECLARATIONS
@@ -104,7 +105,7 @@ REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRHODREF! Reference den
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PPABST ! absolute pressure at t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PTHT ! Theta at time t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRHODJ ! Dry density * Jacobian
-REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PLBDAS ! lambda parameter for snow ! Modif Wurtz
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PLBDAS ! lambda parameter for snow
REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PRCS ! Cloud water m.r. source
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRCT ! Cloud water m.r. at t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PRRS ! Rain water m.r. source
@@ -193,7 +194,7 @@ END IF
CALL INTERNAL_SEDIM_STAT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKTB, KKTE, KKL, &
&PRHODREF, PDZZ, ZW, PPABST, PTHT, PTSTEP, &
&3, &
- &PLBDAS, &!Modif Wurtz
+ &PLBDAS, &
&PRRT, PRRS, ZWSED)
IF (PRESENT(PFPR)) THEN
DO JK = KKTB , KKTE
@@ -207,7 +208,7 @@ PINPRR(:,:) = ZWSED(:,:,KKB)/XRHOLW ! in m/s
CALL INTERNAL_SEDIM_STAT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKTB, KKTE, KKL, &
&PRHODREF, PDZZ, ZW, PPABST, PTHT, PTSTEP, &
&4, &
- &PLBDAS, &!Modif Wurtz
+ &PLBDAS, &
&PRIT, PRIS, ZWSED)
IF (PRESENT(PFPR)) THEN
DO JK = KKTB , KKTE
@@ -221,7 +222,7 @@ PINPRI(:,:) = ZWSED(:,:,KKB)/XRHOLW ! in m/s
CALL INTERNAL_SEDIM_STAT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKTB, KKTE, KKL, &
&PRHODREF, PDZZ, ZW, PPABST, PTHT, PTSTEP, &
&5, &
- &PLBDAS, &!Modif Wurtz
+ &PLBDAS, &
&PRST, PRSS, ZWSED)
IF (PRESENT(PFPR)) THEN
DO JK = KKTB , KKTE
@@ -235,7 +236,7 @@ PINPRS(:,:) = ZWSED(:,:,KKB)/XRHOLW ! in m/s
CALL INTERNAL_SEDIM_STAT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKTB, KKTE, KKL, &
&PRHODREF, PDZZ, ZW, PPABST, PTHT, PTSTEP, &
&6, &
- &PLBDAS, &!Modif Wurtz
+ &PLBDAS, &
&PRGT, PRGS, ZWSED)
IF (PRESENT(PFPR)) THEN
DO JK = KKTB , KKTE
@@ -250,7 +251,7 @@ IF ( KRR == 7 ) THEN
CALL INTERNAL_SEDIM_STAT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKTB, KKTE, KKL, &
&PRHODREF, PDZZ, ZW, PPABST, PTHT, PTSTEP, &
&7, &
- &PLBDAS, &!Modif Wurtz
+ &PLBDAS, &
&PRHT, PRHS, ZWSED)
IF (PRESENT(PFPR)) THEN
DO JK = KKTB , KKTE
@@ -265,7 +266,7 @@ CONTAINS
SUBROUTINE INTERNAL_SEDIM_STAT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKTB, KKTE, KKL, &
&PRHODREF, PDZZ, PTSORHODZ, PPABST, PTHT, PTSTEP, &
&KSPE, &
- &PLBDAS, & !Modif Wurtz
+ &PLBDAS, &
&PRXT, PRXS, PWSED, PSEA, PTOWN)
!
!* 0. DECLARATIONS
@@ -284,7 +285,7 @@ CONTAINS
!
INTEGER, INTENT(IN) :: KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKT, KKE, KKTB, KKTE, KKL
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRHODREF ! Reference density
- REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PLBDAS ! lambda parameter for snow ! Modif Wurtz
+ REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PLBDAS ! lambda parameter for snow
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PDZZ ! Layer thikness (m)
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PTSORHODZ ! TimeStep Over (Rhodref times delta Z)
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PPABST
@@ -409,14 +410,14 @@ CONTAINS
IF(PRXT(JI,JJ,JK) > XRTMIN(KSPE)) THEN
ZWSEDW1(JI,JJ,JK)= XFSEDS * &
- & PRHODREF(JI,JJ,JK)**(-XCEXVT) * & ! Modif Wurtz snow
- & (1+(XFVELOS/PLBDAS(JI,JJ,JK))**XALPHAS)**(-XNUS+XEXSEDS/XALPHAS) * & ! GAMMAGEN LH_EXTENDED
+ & PRHODREF(JI,JJ,JK)**(-XCEXVT) * &
+ & (1+(XFVELOS/PLBDAS(JI,JJ,JK))**XALPHAS)**(-XNUS+XEXSEDS/XALPHAS) * &
& PLBDAS(JI,JJ,JK)**(XBS+XEXSEDS)
ENDIF
IF ( ZQP(JI,JJ) > XRTMIN(KSPE)) THEN
ZWSEDW2(JI,JJ,JK)= XFSEDS * &
- & PRHODREF(JI,JJ,JK)**(-XCEXVT) * & ! Modif Wurtz snow Thompson
- & (1+(XFVELOS/PLBDAS(JI,JJ,JK))**XALPHAS)**(-XNUS+XEXSEDS/XALPHAS) * & ! GAMMAGEN LH_EXTENDED
+ & PRHODREF(JI,JJ,JK)**(-XCEXVT) * &
+ & (1+(XFVELOS/PLBDAS(JI,JJ,JK))**XALPHAS)**(-XNUS+XEXSEDS/XALPHAS) * &
& PLBDAS(JI,JJ,JK)**(XBS+XEXSEDS)
ENDIF
@@ -427,9 +428,6 @@ CONTAINS
IF(KSPE==3) THEN
ZFSED=XFSEDR
ZEXSED=XEXSEDR
- ! ELSEIF(KSPE==5) THEN
- ! ZFSED=XFSEDS
- ! ZEXSED=XEXSEDS
ELSEIF(KSPE==6) THEN
ZFSED=XFSEDG
ZEXSED=XEXSEDG
diff --git a/src/MNH/ice4_slow.f90 b/src/MNH/ice4_slow.f90
index df434ef76..521096775 100644
--- a/src/MNH/ice4_slow.f90
+++ b/src/MNH/ice4_slow.f90
@@ -65,6 +65,7 @@ SUBROUTINE ICE4_SLOW(KSIZE, LDSOFT, PCOMPUTE, PRHODREF, PT, &
!! MODIFICATIONS
!! -------------
!!
+! J. Wurtz 03/2022: New snow characteristics with LSNOW_T
!
!
!* 0. DECLARATIONS
@@ -173,12 +174,9 @@ IF(LDSOFT) THEN
ELSE
PRVDEPS(:) = 0.
WHERE(ZMASK(:)==1.)
-! PRVDEPS(:) = ( PSSI(:)/(PRHODREF(:)*PAI(:)) ) * &
-! ( X0DEPS*PLBDAS(:)**XEX0DEPS + X1DEPS*PCJ(:)*PLBDAS(:)**XEX1DEPS )
- PRVDEPS(:) = ( PRST(:)*PSSI(:)/PAI(:)) * & !Modif Wurtz
- ! ( X0DEPS*PLBDAS(:)**XEX0DEPS + (X1DEPS*PCJ(:)*(PLBDAS(:)+XFVELOS/2.)**(XEX1DEPS)*(PLBDAS(:))**(XNUS+XBS))) ! Thompson
+ PRVDEPS(:) = ( PRST(:)*PSSI(:)/PAI(:)) * &
( X0DEPS*PLBDAS(:)**XEX0DEPS + (X1DEPS*PCJ(:)*(1+(PLBDAS(:)/(2*XFVELOS)**XALPHAS))**(-XNUS+XEX1DEPS) &
- *(PLBDAS(:))**(XBS+XEX1DEPS))) ! GAMMAGEN LH_EXTENDED
+ *(PLBDAS(:))**(XBS+XEX1DEPS)))
END WHERE
ENDIF
DO JL=1, KSIZE
@@ -201,15 +199,11 @@ IF(LDSOFT) THEN
ELSE
PRIAGGS(:) = 0.
WHERE(ZMASK(:)==1)
-! PRIAGGS(:) = XFIAGGS * EXP( XCOLEXIS*(PT(:)-XTT) ) &
-! * PRIT(:) &
-! * PLBDAS(:)**XEXIAGGS &
-! * PRHODREF(:)**(-XCEXVT)
- PRIAGGS(:) = XFIAGGS * EXP( XCOLEXIS*(PT(:)-XTT) ) & !Modif Wurtz GAMMAGEN LH_EXTENDED
+ PRIAGGS(:) = XFIAGGS * EXP( XCOLEXIS*(PT(:)-XTT) ) &
* PRIT(:) &
* PRST(:) * (1+(XFVELOS/PLBDAS(:))**XALPHAS)**(-XNUS+XEXIAGGS/XALPHAS) &
* PRHODREF(:)**(-XCEXVT+1.) &
- * ((PLBDAS(:))**(XBS+XEXIAGGS)) ! Thompson
+ * ((PLBDAS(:))**(XBS+XEXIAGGS))
END WHERE
ENDIF
DO JL=1, KSIZE
diff --git a/src/MNH/ice4_tendencies.f90 b/src/MNH/ice4_tendencies.f90
index 5e5f90006..f012b2184 100644
--- a/src/MNH/ice4_tendencies.f90
+++ b/src/MNH/ice4_tendencies.f90
@@ -12,7 +12,6 @@ SUBROUTINE ICE4_TENDENCIES(KSIZE, KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, K
&PEXN, PRHODREF, PLVFACT, PLSFACT, K1, K2, K3, &
&PPRES, PCF, PSIGMA_RC, &
&PCIT, &
- &PLBDAS ,& ! Wurtz
&PT, PTHT, &
&PRVT, PRCT, PRRT, PRIT, PRST, PRGT, PRHT, &
&PRVHENI_MR, PRRHONG_MR, PRIMLTC_MR, PRSRIMCG_MR, &
@@ -136,7 +135,6 @@ REAL, DIMENSION(KSIZE), INTENT(OUT) :: PHLI_LCF
REAL, DIMENSION(KSIZE), INTENT(OUT) :: PHLI_HRI
REAL, DIMENSION(KSIZE), INTENT(OUT) :: PHLI_LRI
REAL, DIMENSION(KIT,KJT,KKT), INTENT(OUT) :: PRAINFR ! Rain fraction
-REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PLBDAS ! Wurtz
END SUBROUTINE ICE4_TENDENCIES
END INTERFACE
END MODULE MODI_ICE4_TENDENCIES
@@ -147,7 +145,6 @@ SUBROUTINE ICE4_TENDENCIES(KSIZE, KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, K
&PEXN, PRHODREF, PLVFACT, PLSFACT, K1, K2, K3, &
&PPRES, PCF, PSIGMA_RC, &
&PCIT, &
- &PLBDAS ,& ! Wurtz
&PT, PTHT, &
&PRVT, PRCT, PRRT, PRIT, PRST, PRGT, PRHT, &
&PRVHENI_MR, PRRHONG_MR, PRIMLTC_MR, PRSRIMCG_MR, &
@@ -177,6 +174,7 @@ SUBROUTINE ICE4_TENDENCIES(KSIZE, KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, K
!! -------------
!
! P. Wautelet 29/05/2019: remove PACK/UNPACK intrinsics (to get more performance and better OpenACC support)
+! J. Wurtz 03/2022: New snow characteristics with LSNOW_T
!
!
!* 0. DECLARATIONS
@@ -310,7 +308,6 @@ REAL, DIMENSION(KSIZE), INTENT(OUT) :: PHLI_LCF
REAL, DIMENSION(KSIZE), INTENT(OUT) :: PHLI_HRI
REAL, DIMENSION(KSIZE), INTENT(OUT) :: PHLI_LRI
REAL, DIMENSION(KIT,KJT,KKT), INTENT(OUT) :: PRAINFR ! Rain fraction
-REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PLBDAS ! Wurtz
!
!* 0.2 declaration of local variables
!
@@ -321,7 +318,7 @@ REAL, DIMENSION(KSIZE) :: ZRVT, ZRCT, ZRRT, ZRIT, ZRST, ZRGT, &
& ZRF, &
& ZLBDAR, ZLBDAS, ZLBDAG, ZLBDAH, ZLBDAR_RF, &
& ZRGSI, ZRGSI_MR
-REAL, DIMENSION(KIT,KJT,KKT) :: ZRRT3D, ZRST3D, ZRGT3D, ZRHT3D
+REAL, DIMENSION(KIT,KJT,KKT) :: ZLBDAS3D, ZRRT3D, ZRST3D, ZRGT3D, ZRHT3D
INTEGER :: JL
REAL, DIMENSION(KSIZE) :: ZWETG ! 1. if graupel growths in wet mode, 0. otherwise
@@ -404,22 +401,28 @@ ELSE
! 5.1.6 riming-conversion of the large sized aggregates into graupel (old parametrisation)
!
IF(CSNOWRIMING=='OLD ') THEN
-! ZLBDAS(:)=0.
-! WHERE(ZRST(:)>0.)
-! ZLBDAS(:) = MIN(XLBDAS_MAX, XLBS*(PRHODREF(:)*MAX(ZRST(:), XRTMIN(5)))**XLBEXS)
-! END WHERE
- CALL ICE4_RSRIMCG_OLD(KSIZE, ODSOFT, PCOMPUTE==1., &
- &PRHODREF, &
-! &ZLBDAS, &
- &PLBDAS, & ! Wurtz
- &ZT, ZRCT, ZRST, &
- &PRSRIMCG_MR, PB_RS, PB_RG)
- DO JL=1, KSIZE
- ZRST(JL) = ZRST(JL) - PRSRIMCG_MR(JL)
- ZRGT(JL) = ZRGT(JL) + PRSRIMCG_MR(JL)
- ENDDO
+ ZLBDAS(:)=0.
+ IF (LSNOW_T) THEN
+ WHERE (ZRST(:)>XRTMIN(5) .AND. ZT(:)>263.15)
+ ZLBDAS(:) = MAX(MIN(XLBDAS_MAX, 10**(14.554-0.0423*ZT(:))),XLBDAS_MIN)*XTRANS_MP_GAMMAS
+ END WHERE
+ WHERE (ZRST(:)>XRTMIN(5) .AND. ZT(:)<=263.15)
+ ZLBDAS(:) = MAX(MIN(XLBDAS_MAX, 10**(6.226-0.0106*ZT(:))),XLBDAS_MIN)*XTRANS_MP_GAMMAS
+ END WHERE
+ ELSE
+ ZLBDAS(:) = MAX(MIN(XLBDAS_MAX,XLBS*(PRHODREF(:)*MAX(ZRST(:), XRTMIN(5)))**XLBEXS,XLBDAS_MIN)
+ END IF
+ CALL ICE4_RSRIMCG_OLD(KSIZE, ODSOFT, PCOMPUTE==1., &
+ &PRHODREF, &
+ &ZLBDAS, &
+ &ZT, ZRCT, ZRST, &
+ &PRSRIMCG_MR, PB_RS, PB_RG)
+ DO JL=1, KSIZE
+ ZRST(JL) = ZRST(JL) - PRSRIMCG_MR(JL)
+ ZRGT(JL) = ZRGT(JL) + PRSRIMCG_MR(JL)
+ ENDDO
ELSE
- PRSRIMCG_MR(:) = 0.
+ PRSRIMCG_MR(:) = 0.
ENDIF
ENDIF
!
@@ -474,10 +477,17 @@ IF(KSIZE>0) THEN
!
!* compute the slope parameters
!
-! ZLBDAS(:)=0.
-! WHERE(ZRST(:)>0.)
-! ZLBDAS(:) = MIN(XLBDAS_MAX, XLBS*(PRHODREF(:)*MAX(ZRST(:), XRTMIN(5)))**XLBEXS)
-! END WHERE
+ ZLBDAS(:)=0.
+ IF (LSNOW_T) THEN
+ WHERE (ZRST(:)>XRTMIN(5) .AND. ZT(:)>263.15)
+ ZLBDAS(:) = MAX(MIN(XLBDAS_MAX, 10**(14.554-0.0423*ZT(:))),XLBDAS_MIN)*XTRANS_MP_GAMMAS
+ END WHERE
+ WHERE (ZRST(:)>XRTMIN(5) .AND. ZT(:)<=263.15)
+ ZLBDAS(:) = MAX(MIN(XLBDAS_MAX, 10**(6.226-0.0106*ZT(:))),XLBDAS_MIN)*XTRANS_MP_GAMMAS
+ END WHERE
+ ELSE
+ ZLBDAS(:) = MAX(MIN(XLBDAS_MAX,XLBS*(PRHODREF(:)*MAX(ZRST(:), XRTMIN(5)))**XLBEXS,XLBDAS_MIN)
+ END IF
ZLBDAG(:)=0.
WHERE(ZRGT(:)>0.)
ZLBDAG(:) = XLBG*(PRHODREF(:)*MAX(ZRGT(:), XRTMIN(6)))**XLBEXG
@@ -508,8 +518,7 @@ ENDIF
CALL ICE4_SLOW(KSIZE, ODSOFT, PCOMPUTE, PRHODREF, ZT, &
&PSSI, PLVFACT, PLSFACT, &
&ZRVT, ZRCT, ZRIT, ZRST, ZRGT, &
- ! &ZLBDAS, ZLBDAG, &
- &PLBDAS, ZLBDAG, & ! Wurtz
+ &ZLBDAS, ZLBDAG, &
&ZAI, ZCJ, PHLI_HCF, PHLI_HRI, &
&PRCHONI, PRVDEPS, PRIAGGS, PRIAUTS, PRVDEPG, &
&PA_TH, PA_RV, PA_RC, PA_RI, PA_RS, PA_RG)
@@ -546,8 +555,7 @@ END IF
CALL ICE4_FAST_RS(KSIZE, ODSOFT, PCOMPUTE, &
&PRHODREF, PLVFACT, PLSFACT, PPRES, &
&ZDV, ZKA, ZCJ, &
- ! &ZLBDAR, ZLBDAS, &
- &ZLBDAR, PLBDAS, & ! Wurtz
+ &ZLBDAR, ZLBDAS, &
&ZT, ZRVT, ZRCT, ZRRT, ZRST, &
&PRIAGGS, &
&PRCRIMSS, PRCRIMSG, PRSRIMCG, &
@@ -570,8 +578,7 @@ ENDDO
CALL ICE4_FAST_RG(KSIZE, ODSOFT, PCOMPUTE, KRR, &
&PRHODREF, PLVFACT, PLSFACT, PPRES, &
&ZDV, ZKA, ZCJ, PCIT, &
- ! &ZLBDAR, ZLBDAS, ZLBDAG, &
- &ZLBDAR, PLBDAS, ZLBDAG, & ! Wurtz
+ &ZLBDAR, ZLBDAS, ZLBDAG, &
&ZT, ZRVT, ZRCT, ZRRT, ZRIT, ZRST, ZRGT, &
&ZRGSI, ZRGSI_MR(:), &
&ZWETG, &
@@ -590,8 +597,7 @@ IF (KRR==7) THEN
CALL ICE4_FAST_RH(KSIZE, ODSOFT, PCOMPUTE, ZWETG, &
&PRHODREF, PLVFACT, PLSFACT, PPRES, &
&ZDV, ZKA, ZCJ, &
- ! &ZLBDAS, ZLBDAG, ZLBDAR, ZLBDAH, &
- &PLBDAS, ZLBDAG, ZLBDAR, ZLBDAH, & ! Wurtz
+ &ZLBDAS, ZLBDAG, ZLBDAR, ZLBDAH, &
&ZT, ZRVT, ZRCT, ZRRT, ZRIT, ZRST, ZRGT, PRHT, &
&PRCWETH, PRIWETH, PRSWETH, PRGWETH, PRRWETH, &
&PRCDRYH, PRIDRYH, PRSDRYH, PRRDRYH, PRGDRYH, PRDRYHG, PRHMLTR, &
diff --git a/src/MNH/ini_rain_ice.f90 b/src/MNH/ini_rain_ice.f90
index b2d117da3..5b0b0fc66 100644
--- a/src/MNH/ini_rain_ice.f90
+++ b/src/MNH/ini_rain_ice.f90
@@ -102,6 +102,7 @@ END MODULE MODI_INI_RAIN_ICE
!! S. Riette 2016-11: new ICE3/ICE4 options
!! P. Wautelet 22/01/2019 bug correction: incorrect write
! P. Wautelet 26/04/2019: replace non-standard FLOAT function by REAL function
+! J. Wurtz 03/2022: New snow characteristics with LSNOW_T
!
!-------------------------------------------------------------------------------
!
@@ -297,13 +298,25 @@ XF2I = 0.14
!
XAS = 0.02
XBS = 1.9
-!Values of Thompson (WRF)
-XCS = 12.5 ! Wurtz GAMMAGEN LH_EXTENDED
-XDS = 0.4 ! Wurtz GAMMAGEN LH_EXTENDED
-XFVELOS = .5 ! Wurtz GAMMAGEN LH_EXTENDED
+IF (LSNOW_T) THEN
+!Cas Gamma generalisee
+XCS = 11.52
+XDS = 0.39
+XFVELOS =0.097
+!Cas MP
+!XCS = 13.2
+!XDS = 0.423
+!XFVELOS = 25.14
+ELSE
+XCS = 5.
+XDS = 0.27
+XFVELOS = 0.
+END IF
!
+IF (.NOT. LSNOW_T) THEN
XCCS = 5.0
XCXS = 1.0
+END IF
!
XF0S = 0.86
XF1S = 0.28
@@ -379,11 +392,17 @@ XNUR = 1.0 ! Exponential law
XALPHAI = 3.0 ! Gamma law for the ice crystal volume
XNUI = 3.0 ! Gamma law with little dispersion
!
-XALPHAS = .214 ! Exponential law
-XNUS = 43.7 ! Exponential law
-XTRANS_MP_GAMMAS = SQRT( ( GAMMA(XNUS + 2./XALPHAS)*GAMMA(XNUS + 4./XALPHAS) ) / &
- ( 8.* GAMMA(XNUS + 1./XALPHAS)*GAMMA(XNUS + 3./XALPHAS) ) ) ! Wurtz
-!0.5 * GAMMA(XNUS + 2./XALPHAS) / GAMMA(XNUS+1./XALPHAS) ! Wurtz permet de convertir gamma MP en gamma generalisee
+IF (LSNOW_T) THEN
+!Cas GAMMAGEN
+ XALPHAS = .214 ! Generalized gamma law
+ XNUS = 43.7 ! Generalized gamma law
+ XTRANS_MP_GAMMAS = SQRT( ( GAMMA(XNUS + 2./XALPHAS)*GAMMA(XNUS + 4./XALPHAS) ) / &
+ ( 8.* GAMMA(XNUS + 1./XALPHAS)*GAMMA(XNUS + 3./XALPHAS) ) )
+ELSE
+ XALPHAS = 1.0 ! Exponential law
+ XNUS = 1.0 ! Exponential law
+ XTRANS_MP_GAMMAS = 1.
+END IF
!
XALPHAG = 1.0 ! Exponential law
XNUG = 1.0 ! Exponential law
@@ -405,9 +424,13 @@ XLBR = ( XAR*XCCR*MOMG(XALPHAR,XNUR,XBR) )**(-XLBEXR)
XLBEXI = 1.0/(-XBI)
XLBI = ( XAI*MOMG(XALPHAI,XNUI,XBI) )**(-XLBEXI)
!
-!XLBEXS = 1.0/(XCXS-XBS) ! Wurtz Not used new snow diagnostic
-!XLBS = ( XAS*XCCS*MOMG(XALPHAS,XNUS,XBS) )**(-XLBEXS) !Old constant ! OPER
-XLBS = 1.0/(XAS*MOMG(XALPHAS,XNUS,XBS)) !Wurtz for new diag
+IF (LSNOW_T) THEN
+ XLBEXS = 0. ! Not used
+ XLBS = 1.0/(XAS*MOMG(XALPHAS,XNUS,XBS))
+ELSE
+ XLBEXS = 1.0/(XCXS-XBS)
+ XLBS = ( XAS*XCCS*MOMG(XALPHAS,XNUS,XBS) )**(-XLBEXS)
+END IF
!
XLBEXG = 1.0/(XCXG-XBG)
XLBG = ( XAG*XCCG*MOMG(XALPHAG,XNUG,XBG) )**(-XLBEXG)
@@ -419,10 +442,9 @@ XLBH = ( XAH*XCCH*MOMG(XALPHAH,XNUH,XBH) )**(-XLBEXH)
!
XLBDAS_MAX = 100000.0
!
-ZCONC_MAX = 1.E6 ! Maximal concentration for falling particules set to 1 per cc ! Wurtz verifier par rapport a 1E7
-!XLBDAS_MAX = ( ZCONC_MAX/XCCS )**(1./XCXS) ! OLD CANIAUX
-XLBDAS_MAX = 1.E6 ! Wurtz
-XLBDAS_MIN = 1000. ! Valeurs litterature ! Modif Wurtz ! inutile si on se debrouille bien ?
+ZCONC_MAX = 1.E6 ! Maximal concentration for falling particules set to 1 per cc
+XLBDAS_MAX = 1.E6
+XLBDAS_MIN = 1000.
!
IF (HCLOUD == 'ICE4') THEN
ALLOCATE( XRTMIN(7) )
@@ -489,15 +511,20 @@ XEXCSEDI =-0.9324*3.0
WRITE (KLUOUT,FMT=*)' PRISTINE ICE SEDIMENTATION for columns XFSEDI =',XFSEDI
!
!
-!XEXSEDS = (XBS+XDS-XCXS)/(XBS-XCXS) ! Wurtz
-!XFSEDS = XCS*XAS*XCCS*MOMG(XALPHAS,XNUS,XBS+XDS)* & ! Wurtz
-! (XAS*XCCS*MOMG(XALPHAS,XNUS,XBS))**(-XEXSEDS)*(ZRHO00)**XCEXVT ! Wurtz
-
-!XEXSEDS = -XDS !(2*XBS+XDS) ! Modif Wurtz
-!XEXSEDS = -XDS-XBS-XNUS !(2*XBS+XDS) ! Modif Wurtz Thompson
-XEXSEDS = -XDS-XBS !(2*XBS+XDS) ! Modif Wurtz GAMMAGEN LH_EXTENDED
-XFSEDS = XCS*MOMG(XALPHAS,XNUS,XBS+XDS)/(MOMG(XALPHAS,XNUS,XBS)) & ! Modif Wurtz
+IF (LSNOW_T) THEN
+!HOUZE/HAIC
+ !XEXSEDS = -XDS !(2*XBS+XDS)
+ !XFSEDS = XCS*MOMG(XALPHAS,XNUS,XBS+XDS)/(MOMG(XALPHAS,XNUS,XBS)) &
+ ! *(ZRHO00)**XCEXVT
+!LH_EXTENDED
+ XEXSEDS = -XDS-XBS
+ XFSEDS = XCS*MOMG(XALPHAS,XNUS,XBS+XDS)/(MOMG(XALPHAS,XNUS,XBS)) &
*(ZRHO00)**XCEXVT
+ELSE
+ XEXSEDS = (XBS+XDS-XCXS)/(XBS-XCXS)
+ XFSEDS = XCS*XAS*XCCS*MOMG(XALPHAS,XNUS,XBS+XDS)* &
+ (XAS*XCCS*MOMG(XALPHAS,XNUS,XBS))**(-XEXSEDS)*(ZRHO00)**XCEXVT
+END IF
!
XEXSEDG = (XBG+XDG-XCXG)/(XBG-XCXG)
XFSEDG = XCG*XAG*XCCG*MOMG(XALPHAG,XNUG,XBG+XDG)* &
@@ -570,17 +597,10 @@ XSCFAC = (0.63**(1./3.))*SQRT((ZRHO00)**XCEXVT) ! One assumes Sc=0.63
X0DEPI = (4.0*XPI)*XC1I*XF0I*MOMG(XALPHAI,XNUI,1.)
X2DEPI = (4.0*XPI)*XC1I*XF2I*XC_I*MOMG(XALPHAI,XNUI,XDI+2.0)
!
-!X0DEPS = (4.0*XPI)*XCCS*XC1S*XF0S*MOMG(XALPHAS,XNUS,1.)
-!X1DEPS = (4.0*XPI)*XCCS*XC1S*XF1S*SQRT(XCS)*MOMG(XALPHAS,XNUS,0.5*XDS+1.5)
-
-X0DEPS = XLBS*(4.0*XPI)*XC1S*XF0S*MOMG(XALPHAS,XNUS,1.) ! Wurtz
-X1DEPS = XLBS*(4.0*XPI)*XC1S*XF1S*SQRT(XCS)*MOMG(XALPHAS,XNUS,0.5*XDS+1.5) ! Wurtz
-!XEX0DEPS = XCXS-1.0
-!XEX1DEPS = XCXS-0.5*(XDS+3.0)
-XEX0DEPS = XBS-1.0 ! Wurtz
-!XEX1DEPS = XBS-0.5*(XDS+3.0) ! Wurtz Houze
-!XEX1DEPS = -0.5*(XDS+XNUS+3.0) ! Wurtz Thompson
-XEX1DEPS = -0.5*(XDS+3.0) ! Wurtz GAMMGEN LH_EXTENDED
+X0DEPS = XLBS*(4.0*XPI)*XC1S*XF0S*MOMG(XALPHAS,XNUS,1.)
+X1DEPS = XLBS*(4.0*XPI)*XC1S*XF1S*SQRT(XCS)*MOMG(XALPHAS,XNUS,0.5*XDS+1.5)
+XEX0DEPS = XBS-1.0
+XEX1DEPS = -0.5*(XDS+3.0)
!
X0DEPG = (4.0*XPI)*XCCG*XC1G*XF0G*MOMG(XALPHAG,XNUG,1.)
X1DEPG = (4.0*XPI)*XCCG*XC1G*XF1G*SQRT(XCG)*MOMG(XALPHAG,XNUG,0.5*XDG+1.5)
@@ -620,11 +640,7 @@ END IF
!
XCOLIS = 0.25 ! Collection efficiency of I+S
XCOLEXIS = 0.05 ! Temperature factor of the I+S collection efficiency
-!XFIAGGS = (XPI/4.0)*XCOLIS*XCCS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0)
-XFIAGGS = XLBS*(XPI/4.0)*XCOLIS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0) ! Wurtz
-!XEXIAGGS = XCXS-XDS-2.0
-!XEXIAGGS = XBS-XDS-2.0 ! Wurtz Houze
-!XEXIAGGS = -XDS-XNUS-2.0 ! TO MODIFY THOMPSON WURTZ
+XFIAGGS = XLBS*(XPI/4.0)*XCOLIS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0)
XEXIAGGS = -XDS - 2.0 ! GAMMGEN LH_EXTENDED
!
GFLAG = .TRUE.
@@ -676,24 +692,15 @@ XEX1EVAR = -1.0-0.5*(XDR+3.0)
!
XDCSLIM = 0.007 ! D_cs^lim = 7 mm as suggested by Farley et al. (1989)
XCOLCS = 1.0
-!XEXCRIMSS= XCXS-XDS-2.0
-!XEXCRIMSS= XBS-XDS-2.0 ! Wurtz Houze
-!XEXCRIMSS= -XDS-XNUS-2.0 ! Modif Wurtz THOMPSON
-XEXCRIMSS= -XDS-2.0 ! Modif Wurtz GAMMAGEN LH_EXTENDED
-!XCRIMSS = (XPI/4.0)*XCOLCS*XCCS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0)
-XCRIMSS = XLBS * (XPI/4.0)*XCOLCS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0) ! Wurtz
+XEXCRIMSS= -XDS-2.0
+XCRIMSS = XLBS * (XPI/4.0)*XCOLCS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0)
XEXCRIMSG= XEXCRIMSS
XCRIMSG = XCRIMSS
-!XSRIMCG = XCCS*XAS*MOMG(XALPHAS,XNUS,XBS)
-XSRIMCG = XLBS*XAS*MOMG(XALPHAS,XNUS,XBS) ! Wurtz
-!XEXSRIMCG= XCXS-XBS
-!XEXSRIMCG= 0. !XBS-XBS Wurtz
-XEXSRIMCG = -XBS ! Wurtz GAMMAGEN LH_EXTENDED
-!XSRIMCG2 = XCCS*XAG*MOMG(XALPHAS,XNUS,XBG)
-XSRIMCG2 = XLBS*XAG*MOMG(XALPHAS,XNUS,XBG) ! Wurtz
+XSRIMCG = XLBS*XAS*MOMG(XALPHAS,XNUS,XBS)
+XEXSRIMCG = -XBS
+XSRIMCG2 = XLBS*XAG*MOMG(XALPHAS,XNUS,XBG)
XSRIMCG3 = XFRACM90
-!XEXSRIMCG2=XCXS-XBG
-XEXSRIMCG2=XBS-XBG ! Wurtz
+XEXSRIMCG2=XBS-XBG
!
GFLAG = .TRUE.
IF (GFLAG) THEN
@@ -723,15 +730,13 @@ XRIMINTP2 = 1.0 + XRIMINTP1*LOG( XDCSLIM/(XGAMINC_BOUND_MIN)**(1.0/XALPHAS) )
!
!* 7.2 Constants for the accretion of raindrops onto aggregates
!
-!XFRACCSS = ((XPI**2)/24.0)*XCCS*XCCR*XRHOLW*(ZRHO00**XCEXVT)
-XFRACCSS = XLBS*((XPI**2)/24.0)*XCCR*XRHOLW*(ZRHO00**XCEXVT) ! Wurtz
+XFRACCSS = XLBS*((XPI**2)/24.0)*XCCR*XRHOLW*(ZRHO00**XCEXVT)
!
XLBRACCS1 = MOMG(XALPHAS,XNUS,2.)*MOMG(XALPHAR,XNUR,3.)
XLBRACCS2 = 2.*MOMG(XALPHAS,XNUS,1.)*MOMG(XALPHAR,XNUR,4.)
XLBRACCS3 = MOMG(XALPHAR,XNUR,5.)
!
-!XFSACCRG = (XPI/4.0)*XAS*XCCS*XCCR*(ZRHO00**XCEXVT)
-XFSACCRG = XLBS*(XPI/4.0)*XAS*XCCR*(ZRHO00**XCEXVT) ! Wurtz
+XFSACCRG = XLBS*(XPI/4.0)*XAS*XCCR*(ZRHO00**XCEXVT)
!
XLBSACCR1 = MOMG(XALPHAR,XNUR,2.)*MOMG(XALPHAS,XNUS,XBS)
XLBSACCR2 = 2.*MOMG(XALPHAR,XNUR,1.)*MOMG(XALPHAS,XNUS,XBS+1.)
@@ -778,15 +783,15 @@ IF( (KACCLBDAS/=NACCLBDAS) .OR. (KACCLBDAR/=NACCLBDAR) .OR. (KND/=IND) .OR. &
(PACCLBDAS_MIN/=XACCLBDAS_MIN) .OR. (PACCLBDAR_MIN/=XACCLBDAR_MIN) .OR. &
(PFDINFTY/=ZFDINFTY) ) THEN
CALL RRCOLSS ( IND, XALPHAS, XNUS, XALPHAR, XNUR, &
- ZESR, XBR, XCS, XDS, XFVELOS, XCR, XDR, & ! Wurtz Thompson
+ ZESR, XBR, XCS, XDS, XFVELOS, XCR, XDR, &
XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
ZFDINFTY, XKER_RACCSS, XAG, XBS, XAS )
CALL RZCOLX ( IND, XALPHAS, XNUS, XALPHAR, XNUR, &
- ZESR, XBR, XCS, XDS, XFVELOS, XCR, XDR, 0., & ! Wurtz Thompson
+ ZESR, XBR, XCS, XDS, XFVELOS, XCR, XDR, 0., &
XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
ZFDINFTY, XKER_RACCS )
CALL RSCOLRG ( IND, XALPHAS, XNUS, XALPHAR, XNUR, &
- ZESR, XBS, XCS, XDS, XFVELOS, XCR, XDR, & ! Wurtz Thompson
+ ZESR, XBS, XCS, XDS, XFVELOS, XCR, XDR, &
XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
ZFDINFTY, XKER_SACCRG, XAG, XBS, XAS )
WRITE(UNIT=KLUOUT,FMT='("*****************************************")')
@@ -923,8 +928,7 @@ XCOLSG = 0.01 ! Collection efficiency of S+G
XCOLEXSG = 0.1 ! Temperature factor of the S+G collection efficiency
WRITE (KLUOUT, FMT=*) ' NEW Constants for the aggregate collection by the graupeln'
WRITE (KLUOUT, FMT=*) ' XCOLSG, XCOLEXSG = ',XCOLSG,XCOLEXSG
-!XFSDRYG = (XPI/4.0)*XCOLSG*XCCG*XCCS*XAS*(ZRHO00**XCEXVT)
-XFSDRYG = XLBS*(XPI/4.0)*XCOLSG*XCCG*XAS*(ZRHO00**XCEXVT) ! Wurtz
+XFSDRYG = XLBS*(XPI/4.0)*XCOLSG*XCCG*XAS*(ZRHO00**XCEXVT)
!
XLBSDRYG1 = MOMG(XALPHAG,XNUG,2.)*MOMG(XALPHAS,XNUS,XBS)
XLBSDRYG2 = 2.*MOMG(XALPHAG,XNUG,1.)*MOMG(XALPHAS,XNUS,XBS+1.)
@@ -987,7 +991,7 @@ IF( (KDRYLBDAG/=NDRYLBDAG) .OR. (KDRYLBDAS/=NDRYLBDAS) .OR. (KND/=IND) .OR. &
(PDRYLBDAG_MIN/=XDRYLBDAG_MIN) .OR. (PDRYLBDAS_MIN/=XDRYLBDAS_MIN) .OR. &
(PFDINFTY/=ZFDINFTY) ) THEN
CALL RZCOLX ( IND, XALPHAG, XNUG, XALPHAS, XNUS, &
- ZEGS, XBS, XCG, XDG,0., XCS, XDS, XFVELOS, & ! Thompson Wurtz
+ ZEGS, XBS, XCG, XDG,0., XCS, XDS, XFVELOS, &
XDRYLBDAG_MAX, XDRYLBDAS_MAX, XDRYLBDAG_MIN, XDRYLBDAS_MIN, &
ZFDINFTY, XKER_SDRYG )
WRITE(UNIT=KLUOUT,FMT='("*****************************************")')
@@ -1053,7 +1057,7 @@ IF( (KDRYLBDAG/=NDRYLBDAG) .OR. (KDRYLBDAR/=NDRYLBDAR) .OR. (KND/=IND) .OR. &
(PDRYLBDAG_MIN/=XDRYLBDAG_MIN) .OR. (PDRYLBDAR_MIN/=XDRYLBDAR_MIN) .OR. &
(PFDINFTY/=ZFDINFTY) ) THEN
CALL RZCOLX ( IND, XALPHAG, XNUG, XALPHAR, XNUR, &
- ZEGR, XBR, XCG, XDG,0., XCR, XDR, 0., & ! Thompson
+ ZEGR, XBR, XCG, XDG,0., XCR, XDR, 0., &
XDRYLBDAG_MAX, XDRYLBDAR_MAX, XDRYLBDAG_MIN, XDRYLBDAR_MIN, &
ZFDINFTY, XKER_RDRYG )
WRITE(UNIT=KLUOUT,FMT='("*****************************************")')
@@ -1192,7 +1196,7 @@ IF( (KWETLBDAH/=NWETLBDAH) .OR. (KWETLBDAS/=NWETLBDAS) .OR. (KND/=IND) .OR. &
(PWETLBDAH_MIN/=XWETLBDAH_MIN) .OR. (PWETLBDAS_MIN/=XWETLBDAS_MIN) .OR. &
(PFDINFTY/=ZFDINFTY) ) THEN
CALL RZCOLX ( IND, XALPHAH, XNUH, XALPHAS, XNUS, &
- ZEHS, XBS, XCH, XDH,0., XCS, XDS, XFVELOS, & ! Thompson Wurtz
+ ZEHS, XBS, XCH, XDH,0., XCS, XDS, XFVELOS, &
XWETLBDAH_MAX, XWETLBDAS_MAX, XWETLBDAH_MIN, XWETLBDAS_MIN, &
ZFDINFTY, XKER_SWETH )
WRITE(UNIT=KLUOUT,FMT='("*****************************************")')
@@ -1258,7 +1262,7 @@ IF( (KWETLBDAH/=NWETLBDAH) .OR. (KWETLBDAG/=NWETLBDAG) .OR. (KND/=IND) .OR. &
(PWETLBDAH_MIN/=XWETLBDAH_MIN) .OR. (PWETLBDAG_MIN/=XWETLBDAG_MIN) .OR. &
(PFDINFTY/=ZFDINFTY) ) THEN
CALL RZCOLX ( IND, XALPHAH, XNUH, XALPHAG, XNUG, &
- ZEHG, XBG, XCH, XDH, 0., XCG, XDG, 0., & ! Thompson
+ ZEHG, XBG, XCH, XDH, 0., XCG, XDG, 0., &
XWETLBDAH_MAX, XWETLBDAG_MAX, XWETLBDAH_MIN, XWETLBDAG_MIN, &
ZFDINFTY, XKER_GWETH )
WRITE(UNIT=KLUOUT,FMT='("*****************************************")')
@@ -1324,7 +1328,7 @@ IF( (KWETLBDAH/=NWETLBDAH) .OR. (KWETLBDAR/=NWETLBDAR) .OR. (KND/=IND) .OR. &
(PWETLBDAH_MIN/=XWETLBDAH_MIN) .OR. (PWETLBDAR_MIN/=XWETLBDAR_MIN) .OR. &
(PFDINFTY/=ZFDINFTY) ) THEN
CALL RZCOLX ( IND, XALPHAH, XNUH, XALPHAR, XNUR, &
- ZEHR, XBR, XCH, XDH,0., XCR, XDR,0., & ! Thompson
+ ZEHR, XBR, XCH, XDH,0., XCR, XDR,0., &
XWETLBDAH_MAX, XWETLBDAR_MAX, XWETLBDAH_MIN, XWETLBDAR_MIN, &
ZFDINFTY, XKER_RWETH )
WRITE(UNIT=KLUOUT,FMT='("*****************************************")')
@@ -1392,8 +1396,8 @@ IF (GFLAG) THEN
XAS,XBS
WRITE(UNIT=KLUOUT,FMT='(" vitesse: C=",E13.6," D=",E13.6)') &
XCS,XDS
- !WRITE(UNIT=KLUOUT,FMT='(" concentration:CC=",E13.6," x=",E13.6)') & ! Modif Wurtz unused
- ! XCCS,XCXS
+ WRITE(UNIT=KLUOUT,FMT='(" concentration:CC=",E13.6," x=",E13.6)') &
+ XCCS,XCXS
WRITE(UNIT=KLUOUT,FMT='(" distribution:AL=",E13.6,"NU=",E13.6)') &
XALPHAS,XNUS
WRITE(UNIT=KLUOUT,FMT='(" GRAUPEL")')
diff --git a/src/MNH/ini_rain_ice_elec.f90 b/src/MNH/ini_rain_ice_elec.f90
index b99ef0f4c..105c64c68 100644
--- a/src/MNH/ini_rain_ice_elec.f90
+++ b/src/MNH/ini_rain_ice_elec.f90
@@ -87,6 +87,7 @@ END MODULE MODI_INI_RAIN_ICE_ELEC
!! Modifications:
!! C. Barthe 20/11/09 update to version 4.8.1
! P. Wautelet 26/04/2019: replace non-standard FLOAT function by REAL function
+! J. Wurtz 03/2022: New snow characteristics with LSNOW_T
!
!-------------------------------------------------------------------------------
!
@@ -271,13 +272,25 @@ XF2I = 0.14
!
XAS = 0.02
XBS = 1.9
-!Values of Thompson (WRF)
-XCS = 12.5 ! Wurtz GAMMAGEN LH_EXTENDED
-XDS = 0.4 ! Wurtz GAMMAGEN LH_EXTENDED
-XFVELOS = .5 ! Wurtz GAMMAGEN LH_EXTENDED
+IF (LSNOW_T) THEN
+!Cas Gamma generalisee
+XCS = 11.52
+XDS = 0.39
+XFVELOS =0.097
+!Cas MP
+!XCS = 13.2
+!XDS = 0.423
+!XFVELOS = 25.14
+ELSE
+XCS = 5.
+XDS = 0.27
+XFVELOS = 0.
+END IF
!
+IF (.NOT. LSNOW_T) THEN
XCCS = 5.0
XCXS = 1.0
+END IF
!
XF0S = 0.86
XF1S = 0.28
@@ -344,11 +357,17 @@ XNUR = 1.0 ! Exponential law
XALPHAI = 3.0 ! Gamma law for the ice crystal volume
XNUI = 3.0 ! Gamma law with little dispersion
!
-XALPHAS = .214 ! Exponential law
-XNUS = 43.7 ! Exponential law
-XTRANS_MP_GAMMAS = SQRT( ( GAMMA(XNUS + 2./XALPHAS)*GAMMA(XNUS + 4./XALPHAS) ) / &
- ( 8.* GAMMA(XNUS + 1./XALPHAS)*GAMMA(XNUS + 3./XALPHAS) ) ) ! Wurtz
-!0.5 * GAMMA(XNUS + 2./XALPHAS) / GAMMA(XNUS+1./XALPHAS) ! Wurtz permet de convertir gamma MP en gamma generalisee
+IF (LSNOW_T) THEN
+!Cas GAMMAGEN
+ XALPHAS = .214 ! Generalized gamma law
+ XNUS = 43.7 ! Generalized gamma law
+ XTRANS_MP_GAMMAS = SQRT( ( GAMMA(XNUS + 2./XALPHAS)*GAMMA(XNUS + 4./XALPHAS) ) / &
+ ( 8.* GAMMA(XNUS + 1./XALPHAS)*GAMMA(XNUS + 3./XALPHAS) ) )
+ELSE
+ XALPHAS = 1.0 ! Exponential law
+ XNUS = 1.0 ! Exponential law
+ XTRANS_MP_GAMMAS = 1.
+END IF
!
XALPHAG = 1.0 ! Exponential law
XNUG = 1.0 ! Exponential law
@@ -370,9 +389,13 @@ XLBR = (XAR * XCCR * MOMG(XALPHAR,XNUR,XBR))**(-XLBEXR)
XLBEXI = 1.0 / (-XBI)
XLBI = (XAI * MOMG(XALPHAI,XNUI,XBI))**(-XLBEXI)
!
-!XLBEXS = 1.0/(XCXS-XBS) ! Wurtz Not used new snow diagnostic
-!XLBS = ( XAS*XCCS*MOMG(XALPHAS,XNUS,XBS) )**(-XLBEXS) !Old constant ! OPER
-XLBS = 1.0/(XAS*MOMG(XALPHAS,XNUS,XBS)) !Wurtz for new diag
+IF (LSNOW_T) THEN
+ XLBEXS = 0. ! Not used
+ XLBS = 1.0/(XAS*MOMG(XALPHAS,XNUS,XBS))
+ELSE
+ XLBEXS = 1.0/(XCXS-XBS)
+ XLBS = ( XAS*XCCS*MOMG(XALPHAS,XNUS,XBS) )**(-XLBEXS)
+END IF
!
XLBEXG = 1.0 / (XCXG - XBG)
XLBG = (XAG * XCCG * MOMG(XALPHAG,XNUG,XBG))**(-XLBEXG)
@@ -386,10 +409,9 @@ XLBDAR_MAX = 100000.0
XLBDAS_MAX = 100000.0
XLBDAG_MAX = 100000.0
!
-ZCONC_MAX = 1.E6 ! Maximal concentration for falling particules set to 1 per cc ! Wurtz verifier par rapport a 1E7
-!XLBDAS_MAX = ( ZCONC_MAX/XCCS )**(1./XCXS) ! OLD CANIAUX
-XLBDAS_MAX = 1.E6 ! Wurtz
-XLBDAS_MIN = 1000. ! Valeurs litterature ! Modif Wurtz ! inutile si on se debrouille bien ?
+ZCONC_MAX = 1.E6 ! Maximal concentration for falling particules set to 1 per cc
+XLBDAS_MAX = 1.E6
+XLBDAS_MIN = 1000.
!
IF (HCLOUD == 'ICE4') THEN
ALLOCATE( XRTMIN(7) )
@@ -449,17 +471,29 @@ XFSEDI = 3.89745E11 * MOMG(XALPHAI,XNUI,3.285) * &
XEXCSEDI =-0.9324 * 3.0
WRITE (KLUOUT,FMT=*)' PRISTINE ICE SEDIMENTATION for columns XFSEDI =',XFSEDI
!
-XEXSEDS = (XBS + XDS - XCXS) / (XBS - XCXS)
-XFSEDS = XCS * XAS * XCCS * MOMG(XALPHAS,XNUS,XBS+XDS) * &
- (XAS * XCCS * MOMG(XALPHAS,XNUS,XBS))**(-XEXSEDS) * (ZRHO00)**XCEXVT
!
-XEXSEDG = (XBG + XDG - XCXG) / (XBG - XCXG)
-XFSEDG = XCG * XAG * XCCG * MOMG(XALPHAG,XNUG,XBG+XDG) * &
- (XAG * XCCG * MOMG(XALPHAG,XNUG,XBG))**(-XEXSEDG) * (ZRHO00)**XCEXVT
+IF (LSNOW_T) THEN
+!HOUZE/HAIC
+ !XEXSEDS = -XDS !(2*XBS+XDS)
+ !XFSEDS = XCS*MOMG(XALPHAS,XNUS,XBS+XDS)/(MOMG(XALPHAS,XNUS,XBS)) &
+ ! *(ZRHO00)**XCEXVT
+!LH_EXTENDED
+ XEXSEDS = -XDS-XBS
+ XFSEDS = XCS*MOMG(XALPHAS,XNUS,XBS+XDS)/(MOMG(XALPHAS,XNUS,XBS)) &
+ *(ZRHO00)**XCEXVT
+ELSE
+ XEXSEDS = (XBS+XDS-XCXS)/(XBS-XCXS)
+ XFSEDS = XCS*XAS*XCCS*MOMG(XALPHAS,XNUS,XBS+XDS)* &
+ (XAS*XCCS*MOMG(XALPHAS,XNUS,XBS))**(-XEXSEDS)*(ZRHO00)**XCEXVT
+END IF
+!
+XEXSEDG = (XBG+XDG-XCXG)/(XBG-XCXG)
+XFSEDG = XCG*XAG*XCCG*MOMG(XALPHAG,XNUG,XBG+XDG)* &
+ (XAG*XCCG*MOMG(XALPHAG,XNUG,XBG))**(-XEXSEDG)*(ZRHO00)**XCEXVT
!
-XEXSEDH = (XBH + XDH - XCXH) / (XBH - XCXH)
-XFSEDH = XCH * XAH * XCCH * MOMG(XALPHAH,XNUH,XBH+XDH) * &
- (XAH * XCCH * MOMG(XALPHAH,XNUH,XBH))**(-XEXSEDH) * (ZRHO00)**XCEXVT
+XEXSEDH = (XBH+XDH-XCXH)/(XBH-XCXH)
+XFSEDH = XCH*XAH*XCCH*MOMG(XALPHAH,XNUH,XBH+XDH)* &
+ (XAH*XCCH*MOMG(XALPHAH,XNUH,XBH))**(-XEXSEDH)*(ZRHO00)**XCEXVT
!
!
!-------------------------------------------------------------------------------
@@ -524,10 +558,10 @@ XSCFAC = (0.63**(1./3.)) * SQRT((ZRHO00)**XCEXVT) ! One assumes Sc=0.63
X0DEPI = (4.0 * XPI) * XC1I * XF0I * MOMG(XALPHAI,XNUI,1.)
X2DEPI = (4.0 * XPI) * XC1I * XF2I * XC_I * MOMG(XALPHAI,XNUI,XDI+2.0)
!
-X0DEPS = (4.0 * XPI) * XCCS * XC1S * XF0S * MOMG(XALPHAS,XNUS,1.)
-X1DEPS = (4.0 * XPI) * XCCS * XC1S * XF1S * SQRT(XCS) * MOMG(XALPHAS,XNUS,0.5*XDS+1.5)
-XEX0DEPS = XCXS - 1.0
-XEX1DEPS = XCXS - 0.5 * (XDS + 3.0)
+X0DEPS = XLBS*(4.0*XPI)*XC1S*XF0S*MOMG(XALPHAS,XNUS,1.)
+X1DEPS = XLBS*(4.0*XPI)*XC1S*XF1S*SQRT(XCS)*MOMG(XALPHAS,XNUS,0.5*XDS+1.5)
+XEX0DEPS = XBS-1.0
+XEX1DEPS = -0.5*(XDS+3.0)
!
X0DEPG = (4.0 * XPI) * XCCG * XC1G * XF0G * MOMG(XALPHAG,XNUG,1.)
X1DEPG = (4.0 * XPI) * XCCG * XC1G * XF1G * SQRT(XCG) * MOMG(XALPHAG,XNUG,0.5*XDG+1.5)
@@ -560,9 +594,8 @@ END IF
!
XCOLIS = 0.25 ! Collection efficiency of I+S
XCOLEXIS = 0.05 ! Temperature factor of the I+S collection efficiency
-XFIAGGS = (XPI / 4.0) * XCOLIS * XCCS * XCS * (ZRHO00**XCEXVT) * &
- MOMG(XALPHAS,XNUS,XDS+2.0)
-XEXIAGGS = XCXS - XDS - 2.0
+XFIAGGS = XLBS*(XPI/4.0)*XCOLIS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0)
+XEXIAGGS = -XDS - 2.0 ! GAMMGEN LH_EXTENDED
!
GFLAG = .TRUE.
IF (GFLAG) THEN
@@ -612,12 +645,15 @@ XEX1EVAR = -1.0 - 0.5 * (XDR + 3.0)
!
XDCSLIM = 0.007 ! D_cs^lim = 7 mm as suggested by Farley et al. (1989)
XCOLCS = 1.0
-XEXCRIMSS= XCXS - XDS - 2.0
-XCRIMSS = (XPI / 4.0) * XCOLCS * XCCS * XCS * (ZRHO00**XCEXVT) * MOMG(XALPHAS,XNUS,XDS+2.0)
+XEXCRIMSS= -XDS-2.0
+XCRIMSS = XLBS * (XPI/4.0)*XCOLCS*XCS*(ZRHO00**XCEXVT)*MOMG(XALPHAS,XNUS,XDS+2.0)
XEXCRIMSG= XEXCRIMSS
XCRIMSG = XCRIMSS
-XSRIMCG = XCCS * XAS * MOMG(XALPHAS,XNUS,XBS)
-XEXSRIMCG= XCXS - XBS
+XSRIMCG = XLBS*XAS*MOMG(XALPHAS,XNUS,XBS)
+XEXSRIMCG = -XBS
+XSRIMCG2 = XLBS*XAG*MOMG(XALPHAS,XNUS,XBG)
+XSRIMCG3 = XFRACM90
+XEXSRIMCG2=XBS-XBG
!
GFLAG = .TRUE.
IF (GFLAG) THEN
@@ -648,13 +684,13 @@ XRIMINTP2 = 1.0 + XRIMINTP1 * LOG(XDCSLIM/(XGAMINC_BOUND_MIN)**(1.0/XALPHAS))
!
!* 7.2 Constants for the accretion of raindrops onto aggregates
!
-XFRACCSS = ((XPI**2) / 24.0) * XCCS * XCCR * XRHOLW * (ZRHO00**XCEXVT)
+XFRACCSS = XLBS*((XPI**2)/24.0)*XCCR*XRHOLW*(ZRHO00**XCEXVT)
!
XLBRACCS1 = MOMG(XALPHAS,XNUS,2.) * MOMG(XALPHAR,XNUR,3.)
XLBRACCS2 = 2. * MOMG(XALPHAS,XNUS,1.) * MOMG(XALPHAR,XNUR,4.)
XLBRACCS3 = MOMG(XALPHAR,XNUR,5.)
!
-XFSACCRG = (XPI / 4.0) * XAS * XCCS * XCCR * (ZRHO00**XCEXVT)
+XFSACCRG = XLBS*(XPI/4.0)*XAS*XCCR*(ZRHO00**XCEXVT)
!
XLBSACCR1 = MOMG(XALPHAR,XNUR,2.) * MOMG(XALPHAS,XNUS,XBS)
XLBSACCR2 = 2. * MOMG(XALPHAR,XNUR,1.) * MOMG(XALPHAS,XNUS,XBS+1.)
@@ -702,15 +738,15 @@ IF( (KACCLBDAS/=NACCLBDAS) .OR. (KACCLBDAR/=NACCLBDAR) .OR. (KND/=IND) .OR. &
(PACCLBDAS_MIN/=XACCLBDAS_MIN) .OR. (PACCLBDAR_MIN/=XACCLBDAR_MIN) .OR. &
(PFDINFTY/=ZFDINFTY) ) THEN
CALL RRCOLSS ( IND, XALPHAS, XNUS, XALPHAR, XNUR, &
- ZESR, XBR, XCS, XDS, XFVELOS, XCR, XDR, & ! Wurtz Thompson
+ ZESR, XBR, XCS, XDS, XFVELOS, XCR, XDR, &
XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
ZFDINFTY, XKER_RACCSS, XAG, XBS, XAS )
CALL RZCOLX ( IND, XALPHAS, XNUS, XALPHAR, XNUR, &
- ZESR, XBR, XCS, XDS, XFVELOS, XCR, XDR, 0., & ! Wurtz Thompson
+ ZESR, XBR, XCS, XDS, XFVELOS, XCR, XDR, 0., &
XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
ZFDINFTY, XKER_RACCS )
CALL RSCOLRG ( IND, XALPHAS, XNUS, XALPHAR, XNUR, &
- ZESR, XBS, XCS, XDS, XFVELOS, XCR, XDR, & ! Wurtz Thompson
+ ZESR, XBS, XCS, XDS, XFVELOS, XCR, XDR, &
XACCLBDAS_MAX, XACCLBDAR_MAX, XACCLBDAS_MIN, XACCLBDAR_MIN, &
ZFDINFTY, XKER_SACCRG, XAG, XBS, XAS )
WRITE(UNIT=KLUOUT,FMT='("*****************************************")')
@@ -828,9 +864,9 @@ XCOLIG = 0.01 ! Collection efficiency of I+G
XCOLEXIG = 0.1 ! Temperature factor of the I+G collection efficiency
WRITE (KLUOUT, FMT=*) ' NEW Constants for the cloud ice collection by the graupeln'
WRITE (KLUOUT, FMT=*) ' XCOLIG, XCOLEXIG = ',XCOLIG,XCOLEXIG
-!
-XFIDRYG = (XPI / 4.0) * XCOLIG * XCCG * XCG * (ZRHO00**XCEXVT) * &
- MOMG(XALPHAG,XNUG,XDG+2.0)
+XFIDRYG = (XPI/4.0)*XCOLIG*XCCG*XCG*(ZRHO00**XCEXVT)*MOMG(XALPHAG,XNUG,XDG+2.0)
+XEXFIDRYG=(XCXG-XDG-2.)/(XCXG-XBG)
+XFIDRYG2=XFIDRYG/XCOLIG*(XAG*XCCG*MOMG(XALPHAG,XNUG,XBG))**(-XEXFIDRYG)
!
GFLAG = .TRUE.
IF (GFLAG) THEN
@@ -847,8 +883,7 @@ XCOLSG = 0.01 ! Collection efficiency of S+G
XCOLEXSG = 0.1 ! Temperature factor of the S+G collection efficiency
WRITE (KLUOUT, FMT=*) ' NEW Constants for the aggregate collection by the graupeln'
WRITE (KLUOUT, FMT=*) ' XCOLSG, XCOLEXSG = ',XCOLSG,XCOLEXSG
-!
-XFSDRYG = (XPI / 4.0) * XCOLSG * XCCG * XCCS * XAS * (ZRHO00**XCEXVT)
+XFSDRYG = XLBS*(XPI/4.0)*XCOLSG*XCCG*XAS*(ZRHO00**XCEXVT)
!
XLBSDRYG1 = MOMG(XALPHAG,XNUG,2.) * MOMG(XALPHAS,XNUS,XBS)
XLBSDRYG2 = 2. * MOMG(XALPHAG,XNUG,1.) * MOMG(XALPHAS,XNUS,XBS+1.)
@@ -913,7 +948,7 @@ IF( (KDRYLBDAG/=NDRYLBDAG) .OR. (KDRYLBDAS/=NDRYLBDAS) .OR. (KND/=IND) .OR. &
(PDRYLBDAG_MIN/=XDRYLBDAG_MIN) .OR. (PDRYLBDAS_MIN/=XDRYLBDAS_MIN) .OR. &
(PFDINFTY/=ZFDINFTY) ) THEN
CALL RZCOLX ( IND, XALPHAG, XNUG, XALPHAS, XNUS, &
- ZEGS, XBS, XCG, XDG,0., XCS, XDS, XFVELOS, & ! Thompson Wurtz
+ ZEGS, XBS, XCG, XDG,0., XCS, XDS, XFVELOS, &
XDRYLBDAG_MAX, XDRYLBDAS_MAX, XDRYLBDAG_MIN, XDRYLBDAS_MIN, &
ZFDINFTY, XKER_SDRYG )
WRITE(UNIT=KLUOUT,FMT='("*****************************************")')
@@ -979,7 +1014,7 @@ IF( (KDRYLBDAG/=NDRYLBDAG) .OR. (KDRYLBDAR/=NDRYLBDAR) .OR. (KND/=IND) .OR. &
(PDRYLBDAG_MIN/=XDRYLBDAG_MIN) .OR. (PDRYLBDAR_MIN/=XDRYLBDAR_MIN) .OR. &
(PFDINFTY/=ZFDINFTY) ) THEN
CALL RZCOLX ( IND, XALPHAG, XNUG, XALPHAR, XNUR, &
- ZEGR, XBR, XCG, XDG,0., XCR, XDR, 0., & ! Thompson
+ ZEGR, XBR, XCG, XDG,0., XCR, XDR, 0., &
XDRYLBDAG_MAX, XDRYLBDAR_MAX, XDRYLBDAG_MIN, XDRYLBDAR_MIN, &
ZFDINFTY, XKER_RDRYG )
WRITE(UNIT=KLUOUT,FMT='("*****************************************")')
@@ -1041,7 +1076,10 @@ XFWETH = (XPI / 4.0) * XCCH * XCH * (ZRHO00**XCEXVT) * MOMG(XALPHAH,XNUH,XDH+2.0
!
!* 9.2.2 Constants for the aggregate collection by the hailstones
!
-XFSWETH = (XPI/4.0) * XCCH * XCCS * XAS * (ZRHO00**XCEXVT)
+XCOLSH = 0.01 ! Collection efficiency of S+H
+XCOLEXSH = 0.1 ! Temperature factor of the S+H collection efficiency
+!XFSWETH = (XPI/4.0)*XCCH*XCCS*XAS*(ZRHO00**XCEXVT)
+XFSWETH = XLBS*(XPI/4.0)*XCCH*XAS*(ZRHO00**XCEXVT) ! Wurtz
!
XLBSWETH1 = MOMG(XALPHAH,XNUH,2.) * MOMG(XALPHAS,XNUS,XBS)
XLBSWETH2 = 2. * MOMG(XALPHAH,XNUH,1.) * MOMG(XALPHAS,XNUS,XBS+1.)
@@ -1097,7 +1135,7 @@ IF( (KWETLBDAH/=NWETLBDAH) .OR. (KWETLBDAS/=NWETLBDAS) .OR. (KND/=IND) .OR. &
(PWETLBDAH_MIN/=XWETLBDAH_MIN) .OR. (PWETLBDAS_MIN/=XWETLBDAS_MIN) .OR. &
(PFDINFTY/=ZFDINFTY) ) THEN
CALL RZCOLX ( IND, XALPHAH, XNUH, XALPHAS, XNUS, &
- ZEHS, XBS, XCH, XDH,0., XCS, XDS, XFVELOS, & ! Thompson Wurtz
+ ZEHS, XBS, XCH, XDH,0., XCS, XDS, XFVELOS, &
XWETLBDAH_MAX, XWETLBDAS_MAX, XWETLBDAH_MIN, XWETLBDAS_MIN, &
ZFDINFTY, XKER_SWETH )
WRITE(UNIT=KLUOUT,FMT='("*****************************************")')
@@ -1163,7 +1201,7 @@ IF( (KWETLBDAH/=NWETLBDAH) .OR. (KWETLBDAG/=NWETLBDAG) .OR. (KND/=IND) .OR. &
(PWETLBDAH_MIN/=XWETLBDAH_MIN) .OR. (PWETLBDAG_MIN/=XWETLBDAG_MIN) .OR. &
(PFDINFTY/=ZFDINFTY) ) THEN
CALL RZCOLX ( IND, XALPHAH, XNUH, XALPHAG, XNUG, &
- ZEHG, XBG, XCH, XDH, 0., XCG, XDG, 0., & ! Thompson
+ ZEHG, XBG, XCH, XDH, 0., XCG, XDG, 0., &
XWETLBDAH_MAX, XWETLBDAG_MAX, XWETLBDAH_MIN, XWETLBDAG_MIN, &
ZFDINFTY, XKER_GWETH )
WRITE(UNIT=KLUOUT,FMT='("*****************************************")')
diff --git a/src/MNH/modd_rain_ice_descr.f90 b/src/MNH/modd_rain_ice_descr.f90
index 96295e4e7..e754d5a8f 100644
--- a/src/MNH/modd_rain_ice_descr.f90
+++ b/src/MNH/modd_rain_ice_descr.f90
@@ -76,7 +76,11 @@ REAL,SAVE :: XALPHAS,XNUS,XLBEXS,XLBS ! Snow/agg. distribution parameters
REAL,SAVE :: XALPHAG,XNUG,XLBEXG,XLBG ! Graupel distribution parameters
REAL,SAVE :: XALPHAH,XNUH,XLBEXH,XLBH ! Hail distribution parameters
!
-REAL,SAVE :: XLBDAR_MAX,XLBDAS_MAX,XLBDAG_MAX ! Max values allowed for the shape
+LOGICAL,SAVE :: LSNOW_T ! Snow parameterization from Wurtz (2021)
+REAL,SAVE :: XFVELOS ! factor for snow fall speed after Thompson (2008)
+REAL,SAVE :: XTRANS_MP_GAMMAS ! coefficient to convert lambdas for gamma function
+!
+REAL,SAVE :: XLBDAR_MAX,XLBDAS_MIN,XLBDAS_MAX,XLBDAG_MAX ! Max values allowed for the shape
! parameters (rain,snow,graupeln)
!
REAL,DIMENSION(:),SAVE,ALLOCATABLE :: XRTMIN ! Min values allowed for the mixing ratios
diff --git a/src/MNH/modd_rain_ice_param.f90 b/src/MNH/modd_rain_ice_param.f90
index 4394deea2..434c4bc76 100644
--- a/src/MNH/modd_rain_ice_param.f90
+++ b/src/MNH/modd_rain_ice_param.f90
@@ -39,9 +39,6 @@
! ------------
!
IMPLICIT NONE
-!certaines constantes a deplacer dans modd_rain_ice_descr.F90 ... peut etre ...
-REAL,SAVE :: XFVELOS,XLBDAS_MIN ! Wurtz
-REAL,SAVE :: XTRANS_MP_GAMMAS ! coefficient to convert lambdas for gamma function Wurtz
REAL,DIMENSION(2),SAVE :: XFSEDC ! Constants for sedimentation fluxes of C
REAL,SAVE :: XFSEDR,XEXSEDR, & ! Constants for sedimentation
XFSEDI,XEXCSEDI,XEXRSEDI, & ! fluxes of R, I, S and G
diff --git a/src/MNH/rain_ice_red.f90 b/src/MNH/rain_ice_red.f90
index 153617ef6..a44c3af43 100644
--- a/src/MNH/rain_ice_red.f90
+++ b/src/MNH/rain_ice_red.f90
@@ -257,6 +257,7 @@ END MODULE MODI_RAIN_ICE_RED
! P. Wautelet 17/01/2020: move Quicksort to tools.f90
! P. Wautelet 02/2020: use the new data structures and subroutines for budgets
! P. Wautelet 25/02/2020: bugfix: add missing budget: WETH_BU_RRG
+! J. Wurtz 03/2022: New snow characteristics with LSNOW_T
!-----------------------------------------------------------------
!
!* 0. DECLARATIONS
@@ -273,8 +274,7 @@ USE MODD_CST, ONLY: XCI,XCL,XCPD,XCPV,XLSTT,XLVTT,XTT
USE MODD_PARAMETERS, ONLY: JPVEXT,XUNDEF
USE MODD_PARAM_ICE, ONLY: CSUBG_PR_PDF,CSUBG_RC_RR_ACCR,CSUBG_RR_EVAP,LDEPOSC,LFEEDBACKT,LSEDIM_AFTER, &
NMAXITER,XMRSTEP,XTSTEP_TS,XVDEPOSC
-USE MODD_RAIN_ICE_DESCR, ONLY: XRTMIN,XLBDAS_MAX
-USE MODD_RAIN_ICE_PARAM, ONLY: XLBDAS_MIN,XTRANS_MP_GAMMAS
+USE MODD_RAIN_ICE_DESCR, ONLY: XRTMIN,XLBDAS_MIN,XLBDAS_MAX,XTRANS_MP_GAMMAS,LSNOW_T
USE MODD_VAR_ll, ONLY: IP
use mode_budget, only: Budget_store_add, Budget_store_init, Budget_store_end
@@ -357,7 +357,7 @@ REAL, DIMENSION(KIT,KJT,KKT,KRR), OPTIONAL, INTENT(OUT) :: PFPR ! upper-air pre
!
!* 0.2 Declarations of local variables :
!
-REAL, DIMENSION(SIZE(PRST,1),SIZE(PRST,2),SIZE(PRST,3)) :: PLBDAS ! Modif !lbda parameter snow
+REAL, DIMENSION(SIZE(PRST,1),SIZE(PRST,2),SIZE(PRST,3)) :: ZLBDAS ! Modif !lbda parameter snow
INTEGER :: IIB ! Define the domain where is
INTEGER :: IIE ! the microphysical sources have to be computed
@@ -402,7 +402,6 @@ REAL, DIMENSION(KSIZE) :: ZRVT, & ! Water vapor m.r. at t
& ZRRT, & ! Rain water m.r. at t
& ZRIT, & ! Pristine ice m.r. at t
& ZRST, & ! Snow/aggregate m.r. at t
- & ZLBDAS, & ! Snow/aggregate lbdas ! Modif Wurtz
& ZRGT, & ! Graupel m.r. at t
& ZRHT, & ! Hail m.r. at t
& ZCIT, & ! Pristine ice conc. at t
@@ -576,22 +575,25 @@ ELSE
ENDDO
ENDIF
-!Compute lambda_snow parameter ! Modif Wurtz
+!Compute lambda_snow parameter
!ZT en KELVIN
-PLBDAS(:,:,:)=1000. !Wurtz
+ZLBDAS(:,:,:)=1000.
DO JK = 1, KKT
- DO JJ = 1, KJT
- DO JI = 1, KIT
- IF (PRST(JI,JJ,JK)>XRTMIN(5)) THEN
-
- IF(ZT(JI,JJ,JK)>263.15) THEN
- PLBDAS(JI,JJ,JK) = MAX(MIN(XLBDAS_MAX, 10**(14.554-0.0423*ZT(JI,JJ,JK))),XLBDAS_MIN)*XTRANS_MP_GAMMAS ! On le calcul ici car il est utilise dans la sedim en 3D et dans tendencies
- ELSE
- PLBDAS(JI,JJ,JK) = MAX(MIN(XLBDAS_MAX, 10**(6.226-0.0106*ZT(JI,JJ,JK))),XLBDAS_MIN)*XTRANS_MP_GAMMAS
- END IF
- END IF
- END DO
- END DO
+ DO JJ = 1, KJT
+ DO JI = 1, KIT
+ IF (LSNOW_T) THEN
+ IF (PRST(JI,JJ,JK)>XRTMIN(5)) THEN
+ IF(ZT(JI,JJ,JK)>263.15) THEN
+ ZLBDAS(JI,JJ,JK) = MAX(MIN(XLBDAS_MAX, 10**(14.554-0.0423*ZT(JI,JJ,JK))),XLBDAS_MIN)*XTRANS_MP_GAMMAS
+ ELSE
+ ZLBDAS(JI,JJ,JK) = MAX(MIN(XLBDAS_MAX, 10**(6.226-0.0106*ZT(JI,JJ,JK))),XLBDAS_MIN)*XTRANS_MP_GAMMAS
+ END IF
+ END IF
+ ELSE
+ ZLBDAS(JI,JJ,JK) = MAX(MIN(XLBDAS_MAX,XLBS*(PRHODREF(JI,JK,JL)*MAX(PRST(JI,JK,JL), XRTMIN(5)))**XLBEXS,XLBDAS_MIN)
+ END IF
+ END DO
+ END DO
END DO
!
!-------------------------------------------------------------------------------
@@ -622,7 +624,7 @@ IF(.NOT. LSEDIM_AFTER) THEN
CALL ICE4_SEDIMENTATION_STAT(IIB, IIE, KIT, IJB, IJE, KJT, IKB, IKE, IKTB, IKTE, KKT, KKL, &
&PTSTEP, KRR, OSEDIC, LDEPOSC, XVDEPOSC, PDZZ, &
&PRHODREF, PPABST, PTHT, PRHODJ, &
- &PLBDAS, & ! Modif Wurtz
+ &ZLBDAS, &
&PRCS, PRCS*PTSTEP, PRRS, PRRS*PTSTEP, PRIS, PRIS*PTSTEP,&
&PRSS, PRSS*PTSTEP, PRGS, PRGS*PTSTEP,&
&PINPRC, PINDEP, PINPRR, ZINPRI, PINPRS, PINPRG, &
@@ -632,7 +634,7 @@ IF(.NOT. LSEDIM_AFTER) THEN
CALL ICE4_SEDIMENTATION_STAT(IIB, IIE, KIT, IJB, IJE, KJT, IKB, IKE, IKTB, IKTE, KKT, KKL, &
&PTSTEP, KRR, OSEDIC, LDEPOSC, XVDEPOSC, PDZZ, &
&PRHODREF, PPABST, PTHT, PRHODJ, &
- &PLBDAS, & ! Modif Wurtz
+ &ZLBDAS, &
&PRCS, PRCS*PTSTEP, PRRS, PRRS*PTSTEP, PRIS, PRIS*PTSTEP,&
&PRSS, PRSS*PTSTEP, PRGS, PRGS*PTSTEP,&
&PINPRC, PINDEP, PINPRR, ZINPRI, PINPRS, PINPRG, &
@@ -647,7 +649,7 @@ IF(.NOT. LSEDIM_AFTER) THEN
CALL ICE4_SEDIMENTATION_SPLIT(IIB, IIE, KIT, IJB, IJE, KJT, IKB, IKE, IKTB, IKTE, KKT, KKL, &
&PTSTEP, KRR, OSEDIC, LDEPOSC, XVDEPOSC, PDZZ, &
&PRHODREF, PPABST, PTHT, PRHODJ, &
- &PLBDAS, & ! Modif Wurtz
+ &ZLBDAS, &
&PRCS, PRCT, PRRS, PRRT, PRIS, PRIT, PRSS, PRST, PRGS, PRGT,&
&PINPRC, PINDEP, PINPRR, ZINPRI, PINPRS, PINPRG, &
&PSEA=PSEA, PTOWN=PTOWN, &
@@ -656,7 +658,7 @@ IF(.NOT. LSEDIM_AFTER) THEN
CALL ICE4_SEDIMENTATION_SPLIT(IIB, IIE, KIT, IJB, IJE, KJT, IKB, IKE, IKTB, IKTE, KKT, KKL, &
&PTSTEP, KRR, OSEDIC, LDEPOSC, XVDEPOSC, PDZZ, &
&PRHODREF, PPABST, PTHT, PRHODJ, &
- &PLBDAS, & ! Modif Wurtz
+ &ZLBDAS, &
&PRCS, PRCT, PRRS, PRRT, PRIS, PRIT, PRSS, PRST, PRGS, PRGT,&
&PINPRC, PINDEP, PINPRR, ZINPRI, PINPRS, PINPRG, &
&PSEA=PSEA, PTOWN=PTOWN, &
@@ -736,7 +738,6 @@ IF(IMICRO>0) THEN
ELSE
ZHLI_LCF(JL)=0.
ENDIF
- ZLBDAS(JL) = PLBDAS(I1(JL),I2(JL),I3(JL))
ENDDO
IF(GEXT_TEND) THEN
DO JL=1, IMICRO
@@ -886,7 +887,6 @@ DO WHILE(ANY(ZTIME(:)<PTSTEP)) ! Loop to *really* compute tendencies
&ZEXN, ZRHODREF, ZLVFACT, ZLSFACT, I1, I2, I3, &
&ZPRES, ZCF, ZSIGMA_RC,&
&ZCIT, &
- &ZLBDAS, & ! Wurtz
&ZZT, ZTHT, &
&ZRVT, ZRCT, ZRRT, ZRIT, ZRST, ZRGT, ZRHT, &
&ZRVHENI_MR, ZRRHONG_MR, ZRIMLTC_MR, ZRSRIMCG_MR, &
@@ -1712,7 +1712,7 @@ IF(LSEDIM_AFTER) THEN
CALL ICE4_SEDIMENTATION_STAT(IIB, IIE, KIT, IJB, IJE, KJT, IKB, IKE, IKTB, IKTE, KKT, KKL, &
&PTSTEP, KRR, OSEDIC, LDEPOSC, XVDEPOSC, PDZZ, &
&PRHODREF, PPABST, PTHT, PRHODJ, &
- &PLBDAS, & ! Modif Wurtz
+ &ZLBDAS, &
&PRCS, PRCS*PTSTEP, PRRS, PRRS*PTSTEP, PRIS, PRIS*PTSTEP,&
&PRSS, PRSS*PTSTEP, PRGS, PRGS*PTSTEP,&
&PINPRC, PINDEP, PINPRR, ZINPRI, PINPRS, PINPRG, &
@@ -1722,7 +1722,7 @@ IF(LSEDIM_AFTER) THEN
CALL ICE4_SEDIMENTATION_STAT(IIB, IIE, KIT, IJB, IJE, KJT, IKB, IKE, IKTB, IKTE, KKT, KKL, &
&PTSTEP, KRR, OSEDIC, LDEPOSC, XVDEPOSC, PDZZ,&
&PRHODREF, PPABST, PTHT, PRHODJ, &
- &PLBDAS, & ! Modif Wurtz
+ &ZLBDAS, &
&PRCS, PRCS*PTSTEP, PRRS, PRRS*PTSTEP, PRIS, PRIS*PTSTEP,&
&PRSS, PRSS*PTSTEP, PRGS, PRGS*PTSTEP,&
&PINPRC, PINDEP, PINPRR, ZINPRI, PINPRS, PINPRG, &
@@ -1737,7 +1737,7 @@ IF(LSEDIM_AFTER) THEN
CALL ICE4_SEDIMENTATION_SPLIT(IIB, IIE, KIT, IJB, IJE, KJT, IKB, IKE, IKTB, IKTE, KKT, KKL, &
&PTSTEP, KRR, OSEDIC, LDEPOSC, XVDEPOSC, PDZZ, &
&PRHODREF, PPABST, PTHT, PRHODJ, &
- &PLBDAS, & ! Modif Wurtz
+ &ZLBDAS, &
&PRCS, PRCT, PRRS, PRRT, PRIS, PRIT, PRSS, PRST, PRGS, PRGT,&
&PINPRC, PINDEP, PINPRR, ZINPRI, PINPRS, PINPRG, &
&PSEA=PSEA, PTOWN=PTOWN, &
@@ -1746,7 +1746,7 @@ IF(LSEDIM_AFTER) THEN
CALL ICE4_SEDIMENTATION_SPLIT(IIB, IIE, KIT, IJB, IJE, KJT, IKB, IKE, IKTB, IKTE, KKT, KKL, &
&PTSTEP, KRR, OSEDIC, LDEPOSC, XVDEPOSC, PDZZ, &
&PRHODREF, PPABST, PTHT, PRHODJ, &
- &PLBDAS, & ! Modif Wurtz
+ &ZLBDAS, &
&PRCS, PRCT, PRRS, PRRT, PRIS, PRIT, PRSS, PRST, PRGS, PRGT,&
&PINPRC, PINDEP, PINPRR, ZINPRI, PINPRS, PINPRG, &
&PSEA=PSEA, PTOWN=PTOWN, &
--
GitLab