Changeset 5551 for LMDZ6/branches/contrails/libf/phylmd/lmdz_aviation.f90

Timestamp:

Feb 19, 2025, 5:40:06 PM (3 days ago)

Author:

aborella

Message:

Changed input: from absolute to concentration (per mesh to per m3)

File:

• LMDZ6/branches/contrails/libf/phylmd/lmdz_aviation.f90 (modified) (9 diffs)

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

@@ -8 +8 @@
 
 SUBROUTINE aviation_water_emissions( &
-      klon, klev, dtime, paprs, pplay, temp, qtot, cell_area, &
+      klon, klev, dtime, pplay, temp, qtot, &
       flight_h2o, d_q_avi &
       )
 
-USE lmdz_lscp_ini, ONLY: RD, RG
+USE lmdz_lscp_ini, ONLY: RD
 
 IMPLICIT NONE
@@ -18 +18 @@
 INTEGER,                      INTENT(IN)  :: klon, klev ! number of horizontal grid points and vertical levels
 REAL,                         INTENT(IN)  :: dtime      ! time step [s]
-REAL, DIMENSION(klon,klev+1), INTENT(IN)  :: paprs      ! inter-layer pressure [Pa]
 REAL, DIMENSION(klon,klev),   INTENT(IN)  :: pplay      ! mid-layer pressure [Pa]
 REAL, DIMENSION(klon,klev),   INTENT(IN)  :: temp       ! temperature (K)
 REAL, DIMENSION(klon,klev),   INTENT(IN)  :: qtot       ! total specific humidity (in vapor phase) [kg/kg] 
-REAL, DIMENSION(klon),        INTENT(IN)  :: cell_area  ! area of each cell [m2]
-REAL, DIMENSION(klon,klev),   INTENT(IN)  :: flight_h2o ! aviation H2O emitted within the mesh [kgH2O/s/mesh]
+REAL, DIMENSION(klon,klev),   INTENT(IN)  :: flight_h2o ! aviation emitted H2O concentration [kgH2O/s/m3]
 REAL, DIMENSION(klon,klev),   INTENT(OUT) :: d_q_avi    ! water vapor tendency from aviation [kg/kg]
 ! Local
 INTEGER :: i, k
-REAL :: rho, rhodz, dz, M_cell
+REAL :: rho
 
 DO i=1, klon
@@ -34 +32 @@
     rho = pplay(i,k) / temp(i,k) / RD
     !--Dry air mass [kg/m2]
-    rhodz = ( paprs(i,k) - paprs(i,k+1) ) / RG
+    !rhodz = ( paprs(i,k) - paprs(i,k+1) ) / RG
     !--Cell thickness [m]
-    dz = rhodz / rho
+    !dz = rhodz / rho
     !--Cell dry air mass [kg]
-    M_cell = rhodz * cell_area(i)
+    !M_cell = rhodz * cell_area(i)
 
     !--q is the specific humidity (kg/kg humid air) hence the complicated equation to update q
@@ -44 +42 @@
     !      = ( m_dry_air * qold + dm_h2O * (1-qold) ) / (m_dry_air + dm_H2O * (1-qold) ) 
     !--The equation is derived by writing m_humid_air = m_dry_air + m_H2O = m_dry_air / (1-q) 
-    !--flight_h2O is in kg H2O / s / mesh
+    !--flight_h2O is in kg H2O / s / m3
     
-    !d_q_avi(i,k) = ( M_cell * qtot(i,k) + flight_h2o(i,k) * dtime * ( 1. - qtot(i,k) ) ) &
-    !             / ( M_cell             + flight_h2o(i,k) * dtime * ( 1. - qtot(i,k) ) ) &
+    !d_q_avi(i,k) = ( M_cell * qtot(i,k) + V_cell * flight_h2o(i,k) * dtime * ( 1. - qtot(i,k) ) ) &
+    !             / ( M_cell             + V_cell * flight_h2o(i,k) * dtime * ( 1. - qtot(i,k) ) ) &
+    !             - qtot(i,k)
+    !--NB., M_cell = V_cell * rho
+    !d_q_avi(i,k) = ( rho * qtot(i,k) + flight_h2o(i,k) * dtime * ( 1. - qtot(i,k) ) ) &
+    !             / ( rho             + flight_h2o(i,k) * dtime * ( 1. - qtot(i,k) ) ) &
     !             - qtot(i,k)
     !--Same formula, more computationally effective but less readable
     d_q_avi(i,k) = flight_h2o(i,k) * ( 1. - qtot(i,k) ) &
-                 / ( M_cell / dtime / ( 1. - qtot(i,k) ) + flight_h2o(i,k) )
+                 / ( rho / dtime / ( 1. - qtot(i,k) ) + flight_h2o(i,k) )
   ENDDO
 ENDDO
@@ -61 +63 @@
 SUBROUTINE contrails_formation_evolution( &
       dtime, pplay, temp, qsat, qsatl, gamma_cond, rcont_seri, flight_dist, &
-      cldfra, qvc, dz, V_cell, pdf_loc, pdf_scale, pdf_alpha, &
+      cldfra, qvc, dz, pdf_loc, pdf_scale, pdf_alpha, &
       Tcritcont, qcritcont, potcontfraP, potcontfraNP, contfra, &
       dcontfra_cir, dcf_avi, dqvc_avi, dqi_avi &
@@ -85 +87 @@
 REAL, INTENT(IN)  :: gamma_cond   ! condensation threshold w.r.t. qsat [-]
 REAL, INTENT(IN)  :: rcont_seri   ! ratio of contrails fraction to total cloud fraction [-]
-REAL, INTENT(IN)  :: flight_dist  ! aviation distance flown within the mesh [m/s/mesh]
+REAL, INTENT(IN)  :: flight_dist  ! aviation distance flown concentration [m/s/m3]
 REAL, INTENT(IN)  :: cldfra       ! cloud fraction [-]
 REAL, INTENT(IN)  :: qvc          ! gridbox-mean vapor in the cloud [kg/kg]
 REAL, INTENT(IN)  :: dz           ! cell width [m]
-REAL, INTENT(IN)  :: V_cell       ! cell volume [m3]
 REAL, INTENT(IN)  :: pdf_loc      ! location parameter of the clear sky PDF [%]
 REAL, INTENT(IN)  :: pdf_scale    ! scale parameter of the clear sky PDF [%]
@@ -197 +198 @@
 IF ( potcontfraP .GT. eps ) THEN
   contrail_cross_section = CONTRAIL_CROSS_SECTION_ONERA(dz)
-  contfra_new = MIN(1., flight_dist * dtime * contrail_cross_section / V_cell)
+  contfra_new = MIN(1., flight_dist * dtime * contrail_cross_section)
   dcf_avi = potcontfraP * contfra_new
   IF ( cldfra .GT. eps ) THEN
@@ -634 +635 @@
 
     INTEGER,                    INTENT(IN)  :: klon, klev  ! number of horizontal grid points and vertical levels
-    REAL, DIMENSION(klon,klev), INTENT(OUT) :: flight_dist ! Aviation distance flown within the mesh [m/s/mesh]
-    REAL, DIMENSION(klon,klev), INTENT(OUT) :: flight_h2o  ! Aviation H2O emitted within the mesh [kgH2O/s/mesh]
+    REAL, DIMENSION(klon,klev), INTENT(OUT) :: flight_dist ! Aviation distance flown concentration [m/s/m3]
+    REAL, DIMENSION(klon,klev), INTENT(OUT) :: flight_h2o  ! Aviation emitted H2O [kgH2O/s/m3]
 
     !----------------------------------------------------
@@ -650 +651 @@
     IF (is_omp_master) CALL xios_recv_field("KMFLOWN_interp", flight_dist_mpi(:,:,1))
 
-   ! Propagate to other OMP threads: flight_dist_mpi(klon_mpi,klev) to flight_dist(klon,klev)
-   ! (klon_mpi,klon) = (200,50) avec 80 MPI, 4 OMP, nbp40
-   CALL scatter_omp(flight_dist_mpi(:,:,1), flight_dist)
+    ! Propagate to other OMP threads: flight_dist_mpi(klon_mpi,klev) to flight_dist(klon,klev)
+    ! (klon_mpi,klon) = (200,50) avec 80 MPI, 4 OMP, nbp40
+    CALL scatter_omp(flight_dist_mpi(:,:,1), flight_dist)
 
 END SUBROUTINE read_aviation_emissions