From ddabe7cb6a2eaa2134506b57b0537e2cbb3f9d7b Mon Sep 17 00:00:00 2001
From: Quentin Rodier <quentin.rodier@meteo.fr>
Date: Thu, 21 Dec 2023 14:01:41 +0100
Subject: [PATCH] M. Taufour 18/21/2023: Estimation of effective radius of ice
and liquid from microphysics hypotheses options NRADLP 4 and NRADIP 4
---
.../ifs/ice_effective_radius.F90 | 29 ++++++++---
.../ifs/liquid_effective_radius.F90 | 52 +++++++++++++++++++
.../ecrad-1.4.0_mnh/ifs/radiation_scheme.F90 | 11 +++-
src/MNH/ecrad_interface.f90 | 1 +
4 files changed, 85 insertions(+), 8 deletions(-)
diff --git a/src/LIB/RAD/ecrad-1.4.0_mnh/ifs/ice_effective_radius.F90 b/src/LIB/RAD/ecrad-1.4.0_mnh/ifs/ice_effective_radius.F90
index bf4906471..6e1958014 100644
--- a/src/LIB/RAD/ecrad-1.4.0_mnh/ifs/ice_effective_radius.F90
+++ b/src/LIB/RAD/ecrad-1.4.0_mnh/ifs/ice_effective_radius.F90
@@ -5,8 +5,8 @@ INTERFACE
SUBROUTINE ICE_EFFECTIVE_RADIUS &
& (KIDIA, KFDIA, KLON, KLEV, &
- & PPRESSURE, PTEMPERATURE, PCLOUD_FRAC, PQ_ICE, PQ_SNOW, PGEMU, &
- & PRE_UM)
+ & PPRESSURE, PTEMPERATURE, PCLOUD_FRAC, PQ_ICE, PQ_SNOW, PCIT_LIMA, &
+ & PGEMU, PRE_UM)
USE PARKIND1 , ONLY : JPIM, JPRB
! INPUT ARGUMENTS
@@ -24,6 +24,7 @@ REAL(KIND=JPRB), INTENT(IN) :: PCLOUD_FRAC(KLON,KLEV) ! (kg/kg)
REAL(KIND=JPRB), INTENT(IN) :: PQ_ICE(KLON,KLEV) ! (kg/kg)
REAL(KIND=JPRB), INTENT(IN) :: PQ_SNOW(KLON,KLEV) ! (kg/kg)
+REAL(KIND=JPRB), INTENT(IN) :: PCIT_LIMA(KLON,KLEV) ! (/kg)
! *** Single level variable
REAL(KIND=JPRB), INTENT(IN) :: PGEMU(KLON) ! Sine of latitude
@@ -37,8 +38,8 @@ END MODULE MODI_ICE_EFFECTIVE_RADIUS
SUBROUTINE ICE_EFFECTIVE_RADIUS &
& (KIDIA, KFDIA, KLON, KLEV, &
- & PPRESSURE, PTEMPERATURE, PCLOUD_FRAC, PQ_ICE, PQ_SNOW, PGEMU, &
- & PRE_UM)
+ & PPRESSURE, PTEMPERATURE, PCLOUD_FRAC, PQ_ICE, PQ_SNOW, PCIT_LIMA, &
+ & PGEMU, PRE_UM)
! ICE_EFFECTIVE_RADIUS
!
@@ -76,7 +77,8 @@ USE YOERDU , ONLY : REPLOG, REPSCW
USE YOMLUN , ONLY : NULERR
USE YOMCST , ONLY : RD, RPI, RTT
USE MODD_PARAM_ECRAD_n , ONLY : NRADIP, NMINICE, XRE2DE, XRMINICE ! ice optical properties model
-
+USE MODD_PARAM_LIMA , ONLY : XALPHAI,XNUI ! Ice crystals distribution parameters
+USE MODD_PARAM_LIMA_COLD , ONLY : XDELTAI,XBI,XREFFI,XLBI,XLBEXI ! Ice crystals distribution parameters
! -------------------------------------------------------------------
IMPLICIT NONE
@@ -96,6 +98,8 @@ REAL(KIND=JPRB), INTENT(IN) :: PCLOUD_FRAC(KLON,KLEV) ! (kg/kg)
REAL(KIND=JPRB), INTENT(IN) :: PQ_ICE(KLON,KLEV) ! (kg/kg)
REAL(KIND=JPRB), INTENT(IN) :: PQ_SNOW(KLON,KLEV) ! (kg/kg)
+REAL(KIND=JPRB), INTENT(IN) :: PCIT_LIMA(KLON,KLEV) ! (/kg)
+
! *** Single level variable
REAL(KIND=JPRB), INTENT(IN) :: PGEMU(KLON) ! Sine of latitude
@@ -116,6 +120,8 @@ REAL(KIND=JPRB) :: ZDIAMETER_UM ! Effective diameter in microns
! Min effective diameter in microns; may vary with latitude
REAL(KIND=JPRB) :: ZMIN_DIAMETER_UM(KLON)
+REAL(KIND=JPRB) :: ZLBDAI(KLON,KLEV) ! (/kg)
+
INTEGER :: JL, JK
@@ -199,7 +205,18 @@ CASE(3)
ENDIF
ENDDO
ENDDO
-
+
+CASE(4)
+ !!!! Ice effective radius from Yang et al (2013) reff = V/A
+ ! Distribution parameters for ice and snow
+ DO JK = 1,KLEV
+ DO JL = KIDIA,KFDIA
+ IF (PCLOUD_FRAC(JL,JK) > 0.001_JPRB &
+ & .AND. (PQ_ICE(JL,JK)+PQ_SNOW(JL,JK)) > 0.0_JPRB) THEN
+ PRE_UM(JL,JK) = XREFFI*( XLBI*PCIT_LIMA(JL,JK) / PQ_ICE(JL,JK) )**(XLBEXI*(XDELTAI-XBI))
+ ENDIF
+ ENDDO
+ ENDDO
CASE DEFAULT
WRITE(NULERR,'(A,I0,A)') 'ICE EFFECTIVE RADIUS OPTION NRADLP=',NRADIP,' NOT AVAILABLE'
CALL ABOR1('ERROR IN ICE_EFFECTIVE_RADIUS')
diff --git a/src/LIB/RAD/ecrad-1.4.0_mnh/ifs/liquid_effective_radius.F90 b/src/LIB/RAD/ecrad-1.4.0_mnh/ifs/liquid_effective_radius.F90
index 0478665d5..d3e1a888f 100644
--- a/src/LIB/RAD/ecrad-1.4.0_mnh/ifs/liquid_effective_radius.F90
+++ b/src/LIB/RAD/ecrad-1.4.0_mnh/ifs/liquid_effective_radius.F90
@@ -5,6 +5,7 @@ INTERFACE
SUBROUTINE LIQUID_EFFECTIVE_RADIUS &
& (KIDIA, KFDIA, KLON, KLEV, &
& PPRESSURE, PTEMPERATURE, PCLOUD_FRAC, PQ_LIQ, PQ_RAIN, &
+ & PCCT_LIMA, PCRT_LIMA, &
& PLAND_FRAC, PCCN_LAND, PCCN_SEA, &
& PRE_UM)
@@ -27,6 +28,9 @@ REAL(KIND=JPRB), INTENT(IN) :: PCLOUD_FRAC(KLON,KLEV)
REAL(KIND=JPRB), INTENT(IN) :: PQ_LIQ(KLON,KLEV) ! (kg/kg)
REAL(KIND=JPRB), INTENT(IN) :: PQ_RAIN(KLON,KLEV) ! (kg/kg)
+REAL(KIND=JPRB), INTENT (IN) :: PCCT_LIMA(KLON,KLEV) ! cloud water concentration (LIMA)
+REAL(KIND=JPRB), INTENT (IN) :: PCRT_LIMA(KLON,KLEV) ! rain water concentration (LIMA)
+
! *** Single-level variables
REAL(KIND=JPRB), INTENT(IN) :: PLAND_FRAC(KLON) ! 1=land, 0=sea
REAL(KIND=JPRB), INTENT(IN) :: PCCN_LAND(KLON)
@@ -44,6 +48,7 @@ END MODULE MODI_LIQUID_EFFECTIVE_RADIUS
SUBROUTINE LIQUID_EFFECTIVE_RADIUS &
& (KIDIA, KFDIA, KLON, KLEV, &
& PPRESSURE, PTEMPERATURE, PCLOUD_FRAC, PQ_LIQ, PQ_RAIN, &
+ & PCCT_LIMA, PCRT_LIMA, &
& PLAND_FRAC, PCCN_LAND, PCCN_SEA, &
& PRE_UM)
@@ -84,6 +89,8 @@ USE YOM_YGFL , ONLY : YGFL
USE YOERDU , ONLY : REPLOG, REPSCW
USE YOMLUN , ONLY : NULERR
USE YOMCST , ONLY : RD, RPI
+USE MODD_PARAM_LIMA , ONLY : XALPHAR,XNUR, & ! Raindrop distribution parameters
+ & XALPHAC,XNUC ! Cloud droplet distribution parameters
! -------------------------------------------------------------------
@@ -104,6 +111,9 @@ REAL(KIND=JPRB), INTENT(IN) :: PCLOUD_FRAC(KLON,KLEV)
REAL(KIND=JPRB), INTENT(IN) :: PQ_LIQ(KLON,KLEV) ! (kg/kg)
REAL(KIND=JPRB), INTENT(IN) :: PQ_RAIN(KLON,KLEV) ! (kg/kg)
+REAL(KIND=JPRB), INTENT (IN) :: PCCT_LIMA(KLON,KLEV) ! cloud water concentration (LIMA)
+REAL(KIND=JPRB), INTENT (IN) :: PCRT_LIMA(KLON,KLEV) ! rain water concentration (LIMA)
+
! *** Single-level variables
REAL(KIND=JPRB), INTENT(IN) :: PLAND_FRAC(KLON) ! 1=land, 0=sea
REAL(KIND=JPRB), INTENT(IN) :: PCCN_LAND(KLON)
@@ -129,6 +139,7 @@ REAL(KIND=JPRB) :: ZSPECTRAL_DISPERSION
REAL(KIND=JPRB) :: ZNTOT_CM3 ! Number conc in cm-3
REAL(KIND=JPRB) :: ZRE_CUBED
REAL(KIND=JPRB) :: ZLWC_GM3, ZRWC_GM3 ! In-cloud liquid, rain content in g m-3
+REAL(KIND=JPRB) :: PCCT_CM3 ! Droplet concentration(from lima) in cm-3
REAL(KIND=JPRB) :: ZAIR_DENSITY_GM3 ! Air density in g m-3
REAL(KIND=JPRB) :: ZRAIN_RATIO ! Ratio of rain to liquid water content
REAL(KIND=JPRB) :: ZWOOD_FACTOR, ZRATIO
@@ -242,12 +253,53 @@ CASE(2)
! Cloud fraction or liquid+rain water content too low to
! consider this a cloud
PRE_UM(JL,JK) = PP_MIN_RE_UM
+
ENDIF
ENDDO
ENDDO
+
+CASE(4)
+
+ ! Calculates effective radius from LIMA
+ DO JL = KIDIA,KFDIA
+ ! First compute the cloud droplet concentration
+
+ DO JK = 1,KLEV
+
+ ! Only consider cloudy regions
+ IF (PCLOUD_FRAC(JL,JK) >= 0.001_JPRB &
+ & .AND. (PQ_LIQ(JL,JK)+PQ_RAIN(JL,JK)) > 0.0_JPRB) THEN
+
+ ! Compute liquid and rain water contents
+ ZAIR_DENSITY_GM3 = 1000.0_JPRB * PPRESSURE(JL,JK) &
+ & / (RD*PTEMPERATURE(JL,JK))
+ ! In-cloud mean water contents found by dividing by cloud
+ ! fraction
+ PCCT_CM3 = PCCT_LIMA (JL, JK) * ZAIR_DENSITY_GM3 * 1.E-9_JPRB
+ ZLWC_GM3 = ZAIR_DENSITY_GM3 * PQ_LIQ(JL,JK) / PCLOUD_FRAC(JL,JK)
+ PRE_UM(JL,JK)=100*(4.*ZLWC_GM3*(XNUC+2)**2/(PCCT_CM3*3.*RPI*(XNUC**2+XNUC)))**(1./3.)
+! print*,"r_eff=",PRE_UM(JL,JK)
+ if (PRE_UM(JL,JK)>50) then
+ PRE_UM(JL,JK)=50.
+ else if (PRE_UM(JL,JK)<1.) then
+ PRE_UM(JL,JK)=1.
+ end if
+
+
+ ELSE
+ ! Cloud fraction or liquid+rain water content too low to
+ ! consider this a cloud
+ PRE_UM(JL,JK) = 1.
+
+ ENDIF
+
+ ENDDO
+
+ ENDDO
+
CASE DEFAULT
WRITE(NULERR,'(A,I0,A)') 'LIQUID EFFECTIVE RADIUS OPTION IRADLP=',IRADLP,' NOT AVAILABLE'
diff --git a/src/LIB/RAD/ecrad-1.4.0_mnh/ifs/radiation_scheme.F90 b/src/LIB/RAD/ecrad-1.4.0_mnh/ifs/radiation_scheme.F90
index 06cd9e271..270f990cf 100644
--- a/src/LIB/RAD/ecrad-1.4.0_mnh/ifs/radiation_scheme.F90
+++ b/src/LIB/RAD/ecrad-1.4.0_mnh/ifs/radiation_scheme.F90
@@ -25,6 +25,7 @@ SUBROUTINE RADIATION_SCHEME &
& PPRESSURE_H, PTEMPERATURE_H, &
& PQ, PCO2, PCH4, PN2O, PNO2, PCFC11, PCFC12, PHCFC22, PCCL4, PO3_DP, &
& PCLOUD_FRAC, PQ_LIQUID, PQ_ICE, PQ_RAIN, PQ_SNOW, &
+ & PCCT_LIMA, PCRT_LIMA, PCIT_LIMA, &
& PAEROSOL_OLD, PAEROSOL, &
& PFLUX_SW, PFLUX_LW, PFLUX_SW_CLEAR, PFLUX_LW_CLEAR, &
& PFLUX_SW_SURF, PFLUX_LW_SURF, PFLUX_SW_SURF_CLEAR, PFLUX_LW_SURF_CLEAR, &
@@ -186,6 +187,11 @@ REAL(KIND=JPRB), INTENT(IN) :: PQ_LIQUID(KLON,KLEV)
REAL(KIND=JPRB), INTENT(IN) :: PQ_ICE(KLON,KLEV)
REAL(KIND=JPRB), INTENT(IN) :: PQ_RAIN(KLON,KLEV)
REAL(KIND=JPRB), INTENT(IN) :: PQ_SNOW(KLON,KLEV)
+! *** Particle concentrations
+REAL(KIND=JPRB), INTENT (IN) :: PCCT_LIMA(KLON,KLEV) ! cloud water concentration (LIMA)
+REAL(KIND=JPRB), INTENT (IN) :: PCRT_LIMA(KLON,KLEV) ! rain water concentration (LIMA)
+REAL(KIND=JPRB), INTENT (IN) :: PCIT_LIMA(KLON,KLEV) ! ice crystal concentration (LIMA)
+
! *** Aerosol mass mixing ratios
REAL(KIND=JPRB), INTENT(IN) :: PAEROSOL_OLD(KLON,6,KLEV)
@@ -436,13 +442,14 @@ cloud%fraction(KIDIA:KFDIA,:) = PCLOUD_FRAC(KIDIA:KFDIA,:)
! Compute effective radii and convert to metres
CALL LIQUID_EFFECTIVE_RADIUS(KIDIA, KFDIA, KLON, KLEV, &
& PPRESSURE, PTEMPERATURE, PCLOUD_FRAC, PQ_LIQUID, PQ_RAIN, &
+ & PCCT_LIMA, PCRT_LIMA, &
& PLAND_SEA_MASK, PCCN_LAND, PCCN_SEA, &
& PRE_LIQUID_UM)
cloud%re_liq(KIDIA:KFDIA,:) = PRE_LIQUID_UM(KIDIA:KFDIA,:) * 1.0e-6_JPRB
CALL ICE_EFFECTIVE_RADIUS(KIDIA, KFDIA, KLON, KLEV, &
- & PPRESSURE, PTEMPERATURE, PCLOUD_FRAC, PQ_ICE, PQ_SNOW, PGEMU, &
- & PRE_ICE_UM)
+ & PPRESSURE, PTEMPERATURE, PCLOUD_FRAC, PQ_ICE, PQ_SNOW, PCIT_LIMA, &
+ & PGEMU, PRE_ICE_UM)
cloud%re_ice(KIDIA:KFDIA,:) = PRE_ICE_UM(KIDIA:KFDIA,:) * 1.0e-6_JPRB
! Get the cloud overlap decorrelation length (for cloud boundaries),
diff --git a/src/MNH/ecrad_interface.f90 b/src/MNH/ecrad_interface.f90
index 19830f989..96a2821ab 100644
--- a/src/MNH/ecrad_interface.f90
+++ b/src/MNH/ecrad_interface.f90
@@ -441,6 +441,7 @@ CALL RADIATION_SCHEME (IKIDIA, IKFDIA,KLON, KLEV, IAER, &
PAP, PT, PAPH, PTH, &
PQ, ZCCO2, ZCCH4, ZCN2O, ZCNO2, ZCCFC11, ZCCFC12, ZCCFC22, ZCCCL4, ZOZON, &
PCLFR, PQLWC, PQIWC, PQRWC, ZQSWC, &
+ PCCT_C2R2, PCRT_C2R2, PCIT_C1R3, &
ZAEROSOL_OLD, ZAEROSOL, &
ZFLUX_SW, ZFLUX_LW, ZFLUX_SW_CLEAR, ZFLUX_LW_CLEAR, &
ZFLUX_SW_SURF, ZFLUX_LW_SURF, ZFLUX_SW_SURF_CLEAR, ZFLUX_LW_SURF_CLEAR, &
--
GitLab