Changeset 5790 for LMDZ6/branches/contrails/libf/phylmd/lmdz_lscp_ini.f90

Timestamp:

Jul 28, 2025, 6:44:28 PM (4 months ago)

Author:

aborella

Message:

Major modifs to treatment of contrails (from 2 classes to 2 moments) + diagnostics. Increased numerical efficiency

File:

• LMDZ6/branches/contrails/libf/phylmd/lmdz_lscp_ini.f90 (modified) (7 diffs)

LMDZ6/branches/contrails/libf/phylmd/lmdz_lscp_ini.f90

@@ -232 +232 @@
   !$OMP THREADPRIVATE(ok_plane_contrail)
 
-  LOGICAL, SAVE, PROTECTED :: ok_precip_lincontrails=.TRUE.  ! if True, linear contrails can be autoconverted to snow
-  !$OMP THREADPRIVATE(ok_precip_lincontrails)
-
-  LOGICAL, SAVE, PROTECTED :: ok_precip_circontrails=.FALSE. ! if True, cirrus contrails can be autoconverted to snow
-  !$OMP THREADPRIVATE(ok_precip_circontrails)
-
-  REAL, SAVE, PROTECTED :: aspect_ratio_lincontrails=.1      ! [-] aspect ratio of linear contrails
-  !$OMP THREADPRIVATE(aspect_ratio_lincontrails)
-
-  REAL, SAVE, PROTECTED :: coef_mixing_lincontrails          ! [-] tuning coefficient for the linear contrails mixing process
-  !$OMP THREADPRIVATE(coef_mixing_lincontrails)
-  
-  REAL, SAVE, PROTECTED :: coef_shear_lincontrails           ! [-] additional coefficient for the linear contrails shearing process (subprocess of the contrails mixing process)
-  !$OMP THREADPRIVATE(coef_shear_lincontrails)
-  
-  REAL, SAVE, PROTECTED :: chi_mixing_lincontrails=3.        ! [-] factor for increasing the chance that moist air is surrounding linear contrails
-  !$OMP THREADPRIVATE(chi_mixing_lincontrails)
+  REAL, SAVE, PROTECTED :: aspect_ratio_contrails            ! [-] aspect ratio of contrails
+  !$OMP THREADPRIVATE(aspect_ratio_contrails)
+
+  REAL, SAVE, PROTECTED :: coef_mixing_contrails             ! [-] tuning coefficient for the contrails mixing process
+  !$OMP THREADPRIVATE(coef_mixing_contrails)
+  
+  REAL, SAVE, PROTECTED :: coef_shear_contrails              ! [-] additional coefficient for the contrails shearing process (subprocess of the contrails mixing process)
+  !$OMP THREADPRIVATE(coef_shear_contrails)
+  
+  REAL, SAVE, PROTECTED :: chi_mixing_contrails              ! [-] factor for increasing the chance that moist air is surrounding contrails
+  !$OMP THREADPRIVATE(chi_mixing_contrails)
+
+  REAL, SAVE, PROTECTED :: Nice_init_contrails=100.          ! [#/cm3] initial ice crystals concentration in contrails
+  !$OMP THREADPRIVATE(Nice_init_contrails)
+
+  REAL, SAVE, PROTECTED :: N_Brunt_Vaisala_aviation=0.01     ! [s-1] average Brunt Vaisala frequency, for contrail formation
+  !$OMP THREADPRIVATE(N_Brunt_Vaisala_aviation)
+
+  REAL, SAVE, PROTECTED :: circ_0_loss                       ! [m/s] circulation for contrail formation, will be set automatically
+  !$OMP THREADPRIVATE(circ_0_loss)
+
+  REAL, SAVE, PROTECTED :: plume_area_loss                   ! [m2] plume area for contrail formation, will be set automatically
+  !$OMP THREADPRIVATE(plume_area_loss)
+
+  REAL, SAVE, PROTECTED :: nice_init_ref_loss                ! [#/m3] reference ice crystals concentration for contrail formation, will be set automatically
+  !$OMP THREADPRIVATE(nice_init_ref_loss)
+
+  REAL, SAVE, PROTECTED :: Naer_amb=600.                     ! [#/cm3] ambiant background aerosol number concentration, for contrail formation
+  !$OMP THREADPRIVATE(Naer_amb)
+
+  REAL, SAVE, PROTECTED :: raer_amb_mean=15.                 ! [nm] average radius of the background aerosols, for contrail formation
+  !$OMP THREADPRIVATE(raer_amb_mean)
+
+  REAL, SAVE, PROTECTED :: raer_amb_std=2.2                  ! [?] standard deviation of the radius of the background aerosols, for contrail formation
+  !$OMP THREADPRIVATE(raer_amb_std)
+
+  REAL, SAVE, PROTECTED :: r_soot_mean=15.                   ! [nm] average radius of the background aerosols, for contrail formation
+  !$OMP THREADPRIVATE(r_soot_mean)
+
+  REAL, SAVE, PROTECTED :: r_soot_std=1.6                    ! [?] standard deviation of the radius of the background aerosols, for contrail formation
+  !$OMP THREADPRIVATE(r_soot_std)
+
+  REAL, SAVE, PROTECTED :: air_to_fuel_ratio_engine=70.      ! [-] air to fuel ratio engine, for contrail formation
+  !$OMP THREADPRIVATE(air_to_fuel_ratio_engine)
+
+  REAL, SAVE, PROTECTED :: wingspan=60.                      ! [m] average aircraft wingspan, for contrail formation
+  !$OMP THREADPRIVATE(wingspan)
+
+  REAL, SAVE, PROTECTED :: EI_soot_aviation=1.5e15           ! [#/kg] emission index of soot number for a given fuel type
+  !$OMP THREADPRIVATE(EI_soot_aviation)
 
   REAL, SAVE, PROTECTED :: EI_H2O_aviation=1.25              ! [kgH2O/kg] emission index of water vapor for a given fuel type
@@ -273 +306 @@
   REAL, SAVE, PROTECTED :: fallice_cirrus_contrails=1.       ! [m/s] Ice fallspeed velocity in cirrus contrails
   !$OMP THREADPRIVATE(fallice_cirrus_contrails)
+
+  REAL, SAVE, PROTECTED :: eff2vol_radius_contrails=0.7      ! [-] 
+  !$OMP THREADPRIVATE(eff2vol_radius_contrails)
 
   REAL, SAVE, PROTECTED :: aviation_coef=1.                  ! [-] scaling factor for aviation emissions and flown distance
@@ -568 +604 @@
     CALL getin_p('chi_mixing',chi_mixing)
     ! for aviation
-    CALL getin_p('ok_precip_lincontrails',ok_precip_lincontrails)
-    CALL getin_p('ok_precip_circontrails',ok_precip_circontrails)
-    CALL getin_p('aspect_ratio_lincontrails',aspect_ratio_lincontrails)
-    coef_mixing_lincontrails=coef_mixing_lscp
-    CALL getin_p('coef_mixing_lincontrails',coef_mixing_lincontrails)
-    coef_shear_lincontrails=coef_shear_lscp
-    CALL getin_p('coef_shear_lincontrails',coef_shear_lincontrails)
-    CALL getin_p('chi_mixing_lincontrails',chi_mixing_lincontrails)
+    aspect_ratio_contrails=aspect_ratio_cirrus
+    CALL getin_p('aspect_ratio_contrails',aspect_ratio_contrails)
+    coef_mixing_contrails=coef_mixing_lscp
+    CALL getin_p('coef_mixing_contrails',coef_mixing_contrails)
+    coef_shear_contrails=coef_shear_lscp
+    CALL getin_p('coef_shear_contrails',coef_shear_contrails)
+    chi_mixing_contrails=chi_mixing
+    CALL getin_p('chi_mixing_contrails',chi_mixing_contrails)
+    CALL getin_p('Nice_init_contrails',Nice_init_contrails)
     CALL getin_p('EI_H2O_aviation',EI_H2O_aviation)
+    CALL getin_p('EI_soot_aviation',EI_soot_aviation)
+    CALL getin_p('air_to_fuel_ratio_engine',air_to_fuel_ratio_engine)
+    CALL getin_p('wingspan',wingspan)
+    CALL getin_p('Naer_amb',Naer_amb)
+    CALL getin_p('raer_amb_mean',raer_amb_mean)
+    CALL getin_p('raer_amb_std',raer_amb_std)
+    CALL getin_p('r_soot_mean',r_soot_mean)
+    CALL getin_p('r_soot_std',r_soot_std)
+    CALL getin_p('N_Brunt_Vaisala_aviation',N_Brunt_Vaisala_aviation)
     CALL getin_p('qheat_fuel_aviation',qheat_fuel_aviation)
     CALL getin_p('prop_efficiency_aviation',prop_efficiency_aviation)
@@ -584 +630 @@
     CALL getin_p('fallice_linear_contrails',fallice_linear_contrails)
     CALL getin_p('fallice_cirrus_contrails',fallice_cirrus_contrails)
+    CALL getin_p('eff2vol_radius_contrails',eff2vol_radius_contrails)
     CALL getin_p('aviation_coef',aviation_coef)
     ! for cloudth routines
@@ -689 +736 @@
     WRITE(lunout,*) 'lscp_ini, chi_mixing:', chi_mixing
     ! for aviation
-    WRITE(lunout,*) 'lscp_ini, ok_precip_lincontrails:', ok_precip_lincontrails
-    WRITE(lunout,*) 'lscp_ini, ok_precip_circontrails:', ok_precip_circontrails
-    WRITE(lunout,*) 'lscp_ini, aspect_ratio_lincontrails:', aspect_ratio_lincontrails
-    WRITE(lunout,*) 'lscp_ini, coef_mixing_lincontrails:', coef_mixing_lincontrails
-    WRITE(lunout,*) 'lscp_ini, coef_shear_lincontrails:', coef_shear_lincontrails
-    WRITE(lunout,*) 'lscp_ini, chi_mixing_lincontrails:', chi_mixing_lincontrails
+    WRITE(lunout,*) 'lscp_ini, aspect_ratio_contrails:', aspect_ratio_contrails
+    WRITE(lunout,*) 'lscp_ini, coef_mixing_contrails:', coef_mixing_contrails
+    WRITE(lunout,*) 'lscp_ini, coef_shear_contrails:', coef_shear_contrails
+    WRITE(lunout,*) 'lscp_ini, chi_mixing_contrails:', chi_mixing_contrails
+    WRITE(lunout,*) 'lscp_ini, Nice_init_contrails:', Nice_init_contrails
     WRITE(lunout,*) 'lscp_ini, EI_H2O_aviation:', EI_H2O_aviation
+    WRITE(lunout,*) 'lscp_ini, EI_soot_aviation:', EI_soot_aviation
+    WRITE(lunout,*) 'lscp_ini, air_to_fuel_ratio_engine:', air_to_fuel_ratio_engine
+    WRITE(lunout,*) 'lscp_ini, wingspan:', wingspan
+    WRITE(lunout,*) 'lscp_ini, Naer_amb:', Naer_amb
+    WRITE(lunout,*) 'lscp_ini, raer_amb_mean:', raer_amb_mean
+    WRITE(lunout,*) 'lscp_ini, raer_amb_std:', raer_amb_std
+    WRITE(lunout,*) 'lscp_ini, r_soot_mean:', r_soot_mean
+    WRITE(lunout,*) 'lscp_ini, r_soot_std:', r_soot_std
+    WRITE(lunout,*) 'lscp_ini, N_Brunt_Vaisala_aviation:', N_Brunt_Vaisala_aviation
     WRITE(lunout,*) 'lscp_ini, qheat_fuel_aviation:', qheat_fuel_aviation
     WRITE(lunout,*) 'lscp_ini, prop_efficiency_aviation:', prop_efficiency_aviation
@@ -703 +758 @@
     WRITE(lunout,*) 'lscp_ini, fallice_linear_contrails:', fallice_linear_contrails
     WRITE(lunout,*) 'lscp_ini, fallice_cirrus_contrails:', fallice_cirrus_contrails
+    WRITE(lunout,*) 'lscp_ini, eff2vol_radius_contrails:', eff2vol_radius_contrails
     WRITE(lunout,*) 'lscp_ini, aviation_coef:', aviation_coef
     ! for cloudth routines
@@ -765 +821 @@
                       / nu_iwc_pdf_lscp**(1./3.)
 
+    IF ( ok_plane_contrail ) THEN
+      !--Calculated here to lighten calculations
+      !--See Unterstrasser (2016), Lottermoser and Unterstrasser (2025), and lmdz_aviation routine
+      ! U16, Eq. A16 and LU25, section A2
+      plume_area_loss = 2. * RPI * (1.5 + 0.314 * wingspan)**2
+      ! initial circulation [m/s], U16, Eq. A5
+      circ_0_loss = 10. * wingspan - 70.
+      ! reference ice crystals number concentration [#/m3], LU25, section A1
+      nice_init_ref_loss = 3.38e12 / (2. * RPI * (1.5 + 0.314 * 60.3)**2)
+    ENDIF
+
     !AA Temporary initialisation
     a_tr_sca(1) = -0.5