Changeset 4809 for LMDZ6/trunk

Timestamp:

Feb 8, 2024, 4:49:33 PM (10 months ago)

Author:

evignon

Message:

mise a jour des routines proprecip suite à l'atelier nuage

Location:

LMDZ6/trunk/libf/phylmd

Files:

• lmdz_lscp_ini.F90 (modified) (1 diff)
• lmdz_lscp_poprecip.F90 (modified) (11 diffs)

LMDZ6/trunk/libf/phylmd/lmdz_lscp_ini.F90

@@ -158 +158 @@
   !$OMP THREADPRIVATE(rho_snow)
 
-  REAL, SAVE, PROTECTED :: r_rain=100.E-6                   ! A COMMENTER TODO [m]
+  REAL, SAVE, PROTECTED :: r_rain=500.E-6                   ! A COMMENTER TODO [m]
   !$OMP THREADPRIVATE(r_rain)
 
-  REAL, SAVE, PROTECTED :: r_snow=100.E-6                    ! A COMMENTER TODO [m]
+  REAL, SAVE, PROTECTED :: r_snow=1.E-3                    ! A COMMENTER TODO [m]
   !$OMP THREADPRIVATE(r_snow)
 

LMDZ6/trunk/libf/phylmd/lmdz_lscp_poprecip.F90

@@ -48 +48 @@
 REAL,    INTENT(INOUT), DIMENSION(klon) :: precipfraccld  !--fraction of precipitation in the cloudy air IN THE LAYER ABOVE [-]
 
-REAL,    INTENT(INOUT), DIMENSION(klon) :: rain           !--flux of rain coming from the layer above [kg/s/m2]
-REAL,    INTENT(INOUT), DIMENSION(klon) :: rainclr        !--flux of rain in clear sky coming from the layer above [kg/s/m2]
-REAL,    INTENT(IN),    DIMENSION(klon) :: raincld        !--flux of rain in cloudy air coming from the layer above [kg/s/m2]
-REAL,    INTENT(INOUT), DIMENSION(klon) :: snow           !--flux of snow coming from the layer above [kg/s/m2]
-REAL,    INTENT(INOUT), DIMENSION(klon) :: snowclr        !--flux of snow in clear sky coming from the layer above [kg/s/m2]
-REAL,    INTENT(IN),    DIMENSION(klon) :: snowcld        !--flux of snow in cloudy air coming from the layer above [kg/s/m2]
+REAL,    INTENT(INOUT), DIMENSION(klon) :: rain           !--flux of rain gridbox-mean coming from the layer above [kg/s/m2]
+REAL,    INTENT(INOUT), DIMENSION(klon) :: rainclr        !--flux of rain gridbox-mean in clear sky coming from the layer above [kg/s/m2]
+REAL,    INTENT(IN),    DIMENSION(klon) :: raincld        !--flux of rain gridbox-mean in cloudy air coming from the layer above [kg/s/m2]
+REAL,    INTENT(INOUT), DIMENSION(klon) :: snow           !--flux of snow gridbox-mean coming from the layer above [kg/s/m2]
+REAL,    INTENT(INOUT), DIMENSION(klon) :: snowclr        !--flux of snow gridbox-mean in clear sky coming from the layer above [kg/s/m2]
+REAL,    INTENT(IN),    DIMENSION(klon) :: snowcld        !--flux of snow gridbox-mean in cloudy air coming from the layer above [kg/s/m2]
 
 REAL,    INTENT(OUT),   DIMENSION(klon) :: dqreva         !--rain tendency due to evaporation [kg/kg/s]
@@ -110 +110 @@
        ! layer's air so that precipitation at the ground has the
        ! same temperature as the lowermost layer
+       ! we convert the flux into a specific quantity qprecip
        qprecip(i) = (rain(i)+snow(i))*dtime/((paprsdn(i)-paprsup(i))/RG)
        ! t(i,k+1)+d_t(i,k+1): new temperature of the overlying layer
@@ -269 +270 @@
                                                           !--in the current layer
 
-REAL,    INTENT(INOUT), DIMENSION(klon) :: rain           !--flux of rain coming from the layer above [kg/s/m2]
-REAL,    INTENT(INOUT), DIMENSION(klon) :: rainclr        !--flux of rain in clear sky coming from the layer above [kg/s/m2]
-REAL,    INTENT(INOUT), DIMENSION(klon) :: raincld        !--flux of rain in cloudy air coming from the layer above [kg/s/m2]
-REAL,    INTENT(INOUT), DIMENSION(klon) :: snow           !--flux of snow coming from the layer above [kg/s/m2]
-REAL,    INTENT(INOUT), DIMENSION(klon) :: snowclr        !--flux of snow in clear sky coming from the layer above [kg/s/m2]
-REAL,    INTENT(INOUT), DIMENSION(klon) :: snowcld        !--flux of snow in cloudy air coming from the layer above [kg/s/m2]
+REAL,    INTENT(INOUT), DIMENSION(klon) :: rain           !--flux of rain gridbox-mean coming from the layer above [kg/s/m2]
+REAL,    INTENT(INOUT), DIMENSION(klon) :: rainclr        !--flux of rain gridbox-mean in clear sky coming from the layer above [kg/s/m2]
+REAL,    INTENT(INOUT), DIMENSION(klon) :: raincld        !--flux of rain gridbox-mean in cloudy air coming from the layer above [kg/s/m2]
+REAL,    INTENT(INOUT), DIMENSION(klon) :: snow           !--flux of snow gridbox-mean coming from the layer above [kg/s/m2]
+REAL,    INTENT(INOUT), DIMENSION(klon) :: snowclr        !--flux of snow gridbox-mean in clear sky coming from the layer above [kg/s/m2]
+REAL,    INTENT(INOUT), DIMENSION(klon) :: snowcld        !--flux of snow gridbox-mean in cloudy air coming from the layer above [kg/s/m2]
 
 REAL,    INTENT(OUT),   DIMENSION(klon) :: dqrcol         !-- rain tendendy due to collection by rain of liquid cloud droplets [kg/kg/s]
@@ -299 +300 @@
 REAL :: draincld, dsnowcld
 REAL :: eff_cldfra
+REAL :: coef_col, coef_agg, coef_tmp, qrain
 REAL :: qthresh_auto_rain, tau_auto_rain, expo_auto_rain
 REAL :: qthresh_auto_snow, tau_auto_snow, expo_auto_snow
@@ -309 +311 @@
 
 !--Initialisation of variables
-
-
-dqlrim=0.0
-dqlcol=0.0
-dqiagg=0.0
-dqiauto=0.0
-dqlauto=0.0
-
 
 dqrcol(:) = 0.
@@ -329 +323 @@
 
 
-
-!Partitionning between precipitation coming from clouds and that coming from CS
-
-!0) Calculate tot_zneb, total cloud fraction above the cloud
-!assuming a maximum-random overlap (voir Jakob and Klein, 2000)
-
 DO i = 1, klon
 
+
+  ! variables initialisation
+  dqlrim = 0.0
+  dqlcol = 0.0
+  dqiagg = 0.0
+  dqiauto = 0.0
+  dqlauto = 0.0
+
+  !------------------------------------------------------------
+  !--     PRECIPITATION FRACTIONS UPDATE
+  !------------------------------------------------------------
+  !--The goal of this routine is to reattribute precipitation fractions
+  !--and fluxes to clear or cloudy air, depending on the variation of
+  !--the cloud fraction on the vertical dimension. We assume a
+  !--maximum-random overlap of the cloud cover (see Jakob and Klein, 2000,
+  !--and LTP thesis, 2021)
+  !--NB. in fact, we assume a maximum-random overlap of the total precip. frac
+
   !--Initialisation
-  drainclr = 0.
-  dsnowclr = 0.
+  !--hum_to_flux: coef to convert a specific quantity to a flux 
+  !-- hum_to_flux = rho * dz/dt = 1 / g * dP/dt
   hum_to_flux = ( paprsdn(i) - paprsup(i) ) / RG / dtime
   precipfractot = precipfracclr(i) + precipfraccld(i)
 
-  ! formule maximum random
+  !--Instead of using the cloud cover which was use in LTP thesis, we use the
+  !--total precip. fraction to compute the maximum-random overlap. This is
+  !--because all the information of the cloud cover is embedded into
+  !--precipfractot, and this allows for taking into account the potential
+  !--reduction of the precipitation fraction because either the flux is too
+  !--small (see barrier at the end of poprecip_postcld) or the flux is completely
+  !--evaporated (see barrier at the end of poprecip_precld)
+  !--NB. precipfraccld(i) is here the cloud fraction of the layer above
   precipfractot = 1. - ( 1. - precipfractot ) * &
                  ( 1. - MAX( cldfra(i), precipfraccld(i) ) ) &
@@ -349 +362 @@
 
 
-  ! precipfraccld est égal a rneb au niveau d'au dessus
+  !--precipfraccld(i) is here the cloud fraction of the layer above
   dcldfra = cldfra(i) - precipfraccld(i)
-  ! add max to guaranttee precip fraction >=0
-  dprecipfracclr = MAX(0., ( precipfractot - cldfra(i) ) ) - precipfracclr(i)
-  dprecipfraccld = MAX( dcldfra , -precipfraccld(i) ) 
-
-
+  !--Tendency of the clear-sky precipitation fraction. We add a MAX on the
+  !--calculation of the current CS precip. frac.
+  dprecipfracclr = MAX( 0., ( precipfractot - cldfra(i) ) ) - precipfracclr(i)
+  !--Tendency of the cloudy precipitation fraction. We add a MAX on the
+  !--calculation of the current CS precip. frac.
+  !dprecipfraccld = MAX( dcldfra , - precipfraccld(i) ) 
+  !--We remove it, because cldfra is guaranteed to be > O (the MAX is activated
+  !--if cldfra < 0)
+  dprecipfraccld = dcldfra
+
+
+  !--If the cloud extends
   IF ( dprecipfraccld .GT. 0. ) THEN
-    IF ( precipfracclr(i) .GT. 0. ) THEN
+    !--If there is no CS precip, nothing happens.
+    !--If there is, we reattribute some of the CS precip flux
+    !--to the cloud precip flux, proportionnally to the
+    !--decrease of the CS precip fraction
+    IF ( precipfracclr(i) .LE. 0. ) THEN
+      drainclr = 0.
+      dsnowclr = 0.
+    ELSE
       drainclr = dprecipfracclr / precipfracclr(i) * rainclr(i) 
       dsnowclr = dprecipfracclr / precipfracclr(i) * snowclr(i) 
     ENDIF
+  !--If the cloud narrows
+  ELSEIF ( dprecipfraccld .LT. 0. ) THEN
+    !--We reattribute some of the cloudy precip flux
+    !--to the CS precip flux, proportionnally to the
+    !--decrease of the cloud precip fraction
+    draincld = dprecipfraccld / precipfraccld(i) * raincld(i) 
+    dsnowcld = dprecipfraccld / precipfraccld(i) * snowcld(i) 
+    drainclr = - draincld
+    dsnowclr = - dsnowcld
+  !--If the cloud stays the same or if there is no cloud above and
+  !--in the current layer, nothing happens
   ELSE
-    IF  ( precipfraccld(i) .GT. 0. ) THEN
-      draincld = dprecipfraccld / precipfraccld(i) * raincld(i) 
-      dsnowcld = dprecipfraccld / precipfraccld(i) * snowcld(i) 
-      drainclr = - draincld
-      dsnowclr = - dsnowcld
-    ENDIF
+    drainclr = 0.
+    dsnowclr = 0.
   ENDIF
 
-
+  !--We add the tendencies
   precipfraccld(i) = precipfraccld(i) + dprecipfraccld
   precipfracclr(i) = precipfracclr(i) + dprecipfracclr
@@ -378 +412 @@
   snowcld(i) = snowcld(i) - dsnowclr
   
-  
+
+  ! if vertical heterogeneity is taken into account, we use
+  ! the "true" volume fraction instead of a modified 
+  ! surface fraction (which is larger and artificially
+  ! reduces the in-cloud water).
   IF ( ( iflag_cloudth_vert .GE. 3 ) .AND. ( iflag_rain_incloud_vol .EQ. 1 ) ) THEN
     eff_cldfra = ctot_vol(i)
@@ -385 +423 @@
   ENDIF 
 
+
+  ! Start precipitation growth processes
+
+  !--If the cloud is big enough, the precipitation processes activate
   IF ( cldfra(i) .GE. seuil_neb ) THEN
 
-    ! if vertical heterogeneity is taken into account, we use
-    ! the "true" volume fraction instead of a modified 
-    ! surface fraction (which is larger and artificially
-    ! reduces the in-cloud water).
-
-    ! agg et col avant autoconv car on ne veut pas agg et col la nouvelle precip
     !---------------------------------------------------------
     !--           COLLECTION AND AGGREGATION
     !---------------------------------------------------------
-
-    ! remplacer la premiere ligne par "coef_col" ?
-    dqlcol = - gamma_col * 3. / 4. / rho_rain / r_rain * Eff_rain_liq &
-              * qliq(i) / eff_cldfra * raincld(i) * dtime
-
-    ! remplacer la premiere ligne par "coef_agg" ?
-    dqiagg = - gamma_agg * 3. / 4. / rho_snow / r_snow * Eff_snow_ice &
-              * qice(i) / eff_cldfra * snowcld(i) * dtime
-
+    !--Collection: processus through which rain collects small liquid droplets
+    !--in suspension, and add it to the rain flux
+    !--Aggregation: same for snow (precip flux) and ice crystals (in suspension)
+    !--Those processes are treated before autoconversion because we do not
+    !--want to collect/aggregate the newly formed fluxes, which already
+    !--"saw" the cloud as they come from it
+    !--gamma_col: tuning coefficient [-]
+    !--rho_rain: volumic mass of rain [kg/m3]
+    !--r_rain: size of the rain droplets [m]
+    !--Eff_rain_liq: efficiency of the collection process [-] (between 0 and 1)
+    !--dqlcol is a gridbox-mean quantity, as is qliq and raincld. They are
+    !--divided by respectively eff_cldfra, eff_cldfra and precipfraccld to
+    !--get in-cloud mean quantities. The two divisions by eff_cldfra are
+    !--then simplified.
+
+    coef_col = gamma_col * 3. / 4. / rho_rain / r_rain * Eff_rain_liq
+    IF ((raincld(i) .GT. 0.) .AND. (coef_col .GT. 0.)) THEN
+       !--Explicit version
+       !dqlcol = - coef_col * qliq(i) * raincld(i) / precipfraccld(i) *dtime
+       !--Semi-implicit version
+       !dqlcol = qliq(i) * ( 1. / ( 1. + coef_col * raincld(i) / precipfraccld(i)*dtime ) - 1. )
+       !--Implicit version
+       !qrain = raincld(i) / hum_to_flux
+       !coef_tmp = coef_col * dtime *( qrain / precipfraccld(i) + qliq(i) / eff_cldfra )
+       !dqlcol = qliq(i) * ( 1. / ( 1. + 0.5 * ( coef_tmp - 1. + SQRT( &
+       !    ( 1. - coef_tmp )**2. + 4. * coef_col * dtime *qrain / precipfraccld(i) ) &
+       !                   ) ) - 1. )
+       !dqlcol=max(dqlcol,-qliq(i))
+       ! Exact version
+       dqlcol=qliq(i)*(exp(-dtime * coef_col * raincld(i) / precipfraccld(i))-1.)
+    ENDIF
+
+    !--Same as for aggregation
+    coef_agg=gamma_agg * 3. / 4. / rho_snow / r_snow * Eff_snow_ice
+    IF ((snowcld(i) .GT. 0.) .AND. (coef_agg .GT. 0.)) THEN
+        !--Explicit version     
+        !dqiagg = - coef_agg &
+        !      * qice(i) * snowcld(i) / precipfraccld(i) * dtime
+        ! Exact version
+        dqiagg=qice(i)*(exp(-dtime * coef_agg * snowcld(i) / precipfraccld(i))-1.)
+    ENDIF
+    !--Barriers so that the processes do not consume more liquid/ice than
+    !--available.
+    dqlcol = MAX( - qliq(i), dqlcol )
     dqiagg = MAX( - qice(i), dqiagg )
-    dqlcol = MAX( - qliq(i), dqlcol )
-
+
+    !--Add tendencies
     qliq(i) = qliq(i) + dqlcol
     qice(i) = qice(i) + dqiagg
-
     raincld(i) = raincld(i) - dqlcol * hum_to_flux
     snowcld(i) = snowcld(i) - dqiagg * hum_to_flux
@@ -443 +513 @@
 
 
-    ! Liquid water quantity to remove: zchau (Sundqvist, 1978)
+    ! Liquid water quantity to remove according to (Sundqvist, 1978)
     ! dqliq/dt=-qliq/tau*(1-exp(-qcin/clw)**2)
     !.........................................................
@@ -473 +543 @@
     !---------------------------------------------------------
 
+    dqlrim=0.0
+
     ! remplacer la premiere ligne par "coef_rim" ?
-    dqlrim = - gamma_rim * 3. / 4. / rho_snow / r_snow * Eff_snow_liq &
-              * qliq(i) / eff_cldfra * snowcld(i) * dtime
-
+    IF (snowcld(i) .GT. 0.) THEN
+       dqlrim = - gamma_rim * 3. / 4. / rho_snow / r_snow * Eff_snow_liq &
+              * qliq(i) * snowcld(i) /  precipfraccld(i) * dtime
+    ENDIF
     dqlrim = MAX( - qliq(i), dqlrim )