Changeset 5575 for LMDZ6


Ignore:
Timestamp:
Mar 14, 2025, 9:20:21 PM (3 months ago)
Author:
aborella
Message:

Fixed the interpolation of aviation data + added the ok_no_issr_strato option + added additional outputs to compare ISSRs to Lamquin et al (2012)

Location:
LMDZ6/branches/contrails
Files:
12 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/branches/contrails/DefLists/field_def_lmdz.xml

    r5551 r5575  
    588588        <field id="taur" long_name="momentum flux due to rain"         unit="Pa"               detect_missing_value=".true." />
    589589        <field id="SSS" long_name="bulk sea-surface salinity" unit="ppt" detect_missing_value=".true." />
     590        <field id="issrfra100to150" long_name="Supersaturated fraction in the 100to150 hPa layer"    unit="kg/kg" />
     591        <field id="issrfra150to200" long_name="Supersaturated fraction in the 150to200 hPa layer"    unit="kg/kg" />
     592        <field id="issrfra200to250" long_name="Supersaturated fraction in the 200to250 hPa layer"    unit="kg/kg" />
     593        <field id="issrfra250to300" long_name="Supersaturated fraction in the 250to300 hPa layer"    unit="kg/kg" />
     594        <field id="issrfra300to400" long_name="Supersaturated fraction in the 300to400 hPa layer"    unit="kg/kg" />
     595        <field id="issrfra400to500" long_name="Supersaturated fraction in the 400to500 hPa layer"    unit="kg/kg" />
    590596            <!-- Begin Added SN isotopes 2D fields 07 2023 -->
    591597            <!-- water oxygen H216O H217O H218O -->

  • LMDZ6/branches/contrails/libf/phylmd/clesphys_mod_h.f90

    r5536 r5575  
    4444          , ok_lic_melt, ok_lic_cond, aer_type                         &
    4545          , iflag_rrtm, ok_strato, ok_hines, ok_qch4                    &
    46           , iflag_ice_thermo, ok_ice_supersat                            &
     46          , iflag_ice_thermo, ok_ice_supersat, ok_no_issr_strato       &
    4747          , ok_plane_h2o, ok_plane_contrail                            &
    4848          , ok_gwd_rando, NSW, iflag_albedo                            &
    … …  
    141141  LOGICAL :: ok_airs
    142142  INTEGER :: ip_ebil_phy, iflag_rrtm, iflag_ice_thermo, NSW, iflag_albedo
    143   LOGICAL :: ok_ice_supersat, ok_plane_h2o, ok_plane_contrail
     143  LOGICAL :: ok_ice_supersat, ok_no_issr_strato
     144  LOGICAL :: ok_plane_h2o, ok_plane_contrail
    144145  LOGICAL :: ok_chlorophyll
    145146  LOGICAL :: ok_strato
    … …  
    201202  !$OMP      , ok_lic_melt, ok_lic_cond, aer_type                         &
    202203  !$OMP      , iflag_rrtm, ok_strato, ok_hines, ok_qch4                    &
    203   !$OMP      , iflag_ice_thermo, ok_ice_supersat                            &
     204  !$OMP      , iflag_ice_thermo, ok_ice_supersat, ok_no_issr_strato       &
    204205  !$OMP      , ok_plane_h2o, ok_plane_contrail                            &
    205206  !$OMP      , ok_gwd_rando, NSW, iflag_albedo                            &

  • LMDZ6/branches/contrails/libf/phylmd/conf_phys_m.f90

    r5536 r5575  
    172172    INTEGER,SAVE :: iflag_clw_omp
    173173    INTEGER,SAVE :: iflag_ice_thermo_omp
    174     LOGICAL,SAVE :: ok_ice_supersat_omp
     174    LOGICAL,SAVE :: ok_ice_supersat_omp, ok_no_issr_strato_omp
    175175    LOGICAL,SAVE :: ok_plane_h2o_omp, ok_plane_contrail_omp
    176176    REAL,SAVE :: rad_froid_omp, rad_chau1_omp, rad_chau2_omp
    … …  
    13541354    CALL getin('ok_ice_supersat',ok_ice_supersat_omp)
    13551355
     1356    !
     1357    !Config Key  = ok_no_issr_strato
     1358    !Config Desc = deactivates ice supersaturation in the stratosphere
     1359    !Config Def  = .FALSE.
     1360    !Config Help =
     1361    !
     1362    ok_no_issr_strato_omp = .FALSE.
     1363    CALL getin('ok_no_issr_strato',ok_no_issr_strato_omp)
     1364
    13561365    !Config Key  = ok_plane_h2o
    13571366    !Config Desc = include H2O emissions from aviation
    … …  
    23452354    iflag_ice_thermo = iflag_ice_thermo_omp
    23462355    ok_ice_supersat = ok_ice_supersat_omp
     2356    ok_no_issr_strato = ok_no_issr_strato_omp
    23472357    ok_plane_h2o = ok_plane_h2o_omp
    23482358    ok_plane_contrail = ok_plane_contrail_omp
    … …  
    27702780    WRITE(lunout,*) ' iflag_ice_thermo = ',iflag_ice_thermo
    27712781    WRITE(lunout,*) ' ok_ice_supersat = ',ok_ice_supersat
     2782    WRITE(lunout,*) ' ok_no_issr_strato = ',ok_no_issr_strato
    27722783    WRITE(lunout,*) ' ok_plane_h2o = ',ok_plane_h2o
    27732784    WRITE(lunout,*) ' ok_plane_contrail = ',ok_plane_contrail

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

    r5573 r5575  
    44
    55IMPLICIT NONE
     6
     7! Arrays for the lecture of aviation files
     8! The allocation is done in the read_aviation module
     9! The size is (klon, nleva, 1) where
     10! nleva            is the size of the vertical axis (read from file)
     11! flight_dist_read is the number of km per second
     12! flight_h2o_read  is the water content added to the air
     13! aviation_lev     is the value of the levels
     14INTEGER, SAVE :: nleva  ! Size of the vertical axis in the file
     15!$OMP THREADPRIVATE(nleva)
     16REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE, PRIVATE :: flight_dist_read ! Aviation distance flown within the mesh [m/s/mesh]
     17!$OMP THREADPRIVATE(flight_dist_read)
     18REAL, ALLOCATABLE, DIMENSION(:,:,:), SAVE, PRIVATE :: flight_h2o_read  ! Aviation H2O emitted within the mesh [kgH2O/s/mesh]
     19!$OMP THREADPRIVATE(flight_h2o_read)
     20REAL, ALLOCATABLE, DIMENSION(:),     SAVE, PRIVATE :: aviation_lev     ! Pressure in the middle of the layers [Pa]
     21!$OMP THREADPRIVATE(aviation_lev)
    622
    723CONTAINS
    … …  
    3147    !--Dry density [kg/m3]
    3248    rho = pplay(i,k) / temp(i,k) / RD
    33     !--Dry air mass [kg/m2]
    34     !rhodz = ( paprs(i,k) - paprs(i,k+1) ) / RG
    35     !--Cell thickness [m]
    36     !dz = rhodz / rho
    37     !--Cell dry air mass [kg]
    38     !M_cell = rhodz * cell_area(i)
    3949
    4050    !--q is the specific humidity (kg/kg humid air) hence the complicated equation to update q
    … …  
    622632END FUNCTION contrail_cross_section_onera
    623633
    624 SUBROUTINE read_aviation_emissions(klon, klev, flight_dist_read, flight_h2o_read, aviation_lev, nleva)
     634SUBROUTINE read_aviation_emissions(klon, klev)
    625635    ! This subroutine allows to read the traffic density data read in the file aviation.nc
    626636    ! This file is defined in ./COMP/lmdz.card
    … …  
    634644    IMPLICIT NONE
    635645
    636     INTEGER,                    INTENT(IN)  :: klon, klev  ! number of horizontal grid points and vertical levels
    637     INTEGER, INTENT(out) :: nleva  ! Size of the vertical axis in the file
    638     !REAL, DIMENSION(klon,klev), INTENT(OUT) :: flight_dist ! Aviation distance flown within the mesh [m/s/mesh]
    639     !REAL, DIMENSION(klon,klev), INTENT(OUT) :: flight_h2o  ! Aviation H2O emitted within the mesh [kgH2O/s/mesh]
    640     REAL, ALLOCATABLE, INTENT(OUT) :: flight_dist_read(:,:,:) ! Aviation distance flown within the mesh [m/s/mesh]
    641     REAL, ALLOCATABLE, INTENT(OUT) :: flight_h2o_read(:,:,:)  ! Aviation H2O emitted within the mesh [kgH2O/s/mesh]
    642     REAL, ALLOCATABLE, INTENT(OUT) :: aviation_lev(:)  ! Pressure in the middle of the layers [Pa]
     646    INTEGER, INTENT(IN) :: klon, klev  ! number of horizontal grid points and vertical levels
    643647
    644648    !----------------------------------------------------
    645649    ! Local variable
    646650    !----------------------------------------------------
    647     !REAL, DIMENSION(klon_mpi,klev,1) :: flight_dist_mpi
    648651    REAL, ALLOCATABLE :: flight_dist_mpi(:,:,:)
    649652    INTEGER :: ierr
    650653
    651654    ! Get number of vertical levels and level values
    652     CALL xios_get_axis_attr( "aviation_lev", n_glo=nleva )
     655    IF (is_omp_master) CALL xios_get_axis_attr( "aviation_lev", n_glo=nleva )
     656    CALL bcast_omp(nleva)
    653657
    654658    ! Allocation of arrays
    655     !$OMP MASTER
    656     ALLOCATE(aviation_lev(nleva), STAT=ierr)
    657     IF (ierr /= 0) CALL abort_physic('read_aviation_emissions', 'problem to allocate aviation_lev',1)
    658659    ALLOCATE(flight_dist_read(klon, nleva,1), STAT=ierr)
    659660    IF (ierr /= 0) CALL abort_physic('read_aviation_emissions', 'problem to allocate flight_dist',1)
    660661    ALLOCATE(flight_h2o_read(klon, nleva,1), STAT=ierr)
    661662    IF (ierr /= 0) CALL abort_physic('read_aviation_emissions', 'problem to allocate flight_h2o',1)
    662     ALLOCATE(flight_dist_mpi(klon_mpi, nleva,1), STAT=ierr)
    663     IF (ierr /= 0) CALL abort_physic('read_aviation_emissions', 'problem to allocate flight_dist_mpi',1)
    664     !$OMP END MASTER
    665 
    666     !$OMP BARRIER  ! Ensure all threads wait until the arrays are allocated
    667 
    668     !--Initialisation
    669     aviation_lev(:) = 0.
    670     flight_dist_read(:,:,1) = 0.
    671     flight_h2o_read(:,:,1) = 0.
    672 
    673     ! Get number of vertical levels and level values
    674     CALL xios_get_axis_attr( "aviation_lev", value=aviation_lev(:))
     663    ALLOCATE(aviation_lev(nleva), STAT=ierr)
     664    IF (ierr /= 0) CALL abort_physic('read_aviation_emissions', 'problem to allocate aviation_lev',1)
    675665
    676666    ! Read the data from the file
    677667    ! is_omp_master is necessary to make XIOS works
    678     IF (is_omp_master) CALL xios_recv_field("KMFLOWN_interp", flight_dist_mpi(:,:,1))
    679 
    680    ! Propagate to other OMP threads: flight_dist_mpi(klon_mpi,klev) to flight_dist(klon,klev)
    681    ! (klon_mpi,klon) = (200,50) avec 80 MPI, 4 OMP, nbp40
    682    CALL scatter_omp(flight_dist_mpi(:,:,1), flight_dist_read(:,:,1))
     668    IF (is_omp_master) THEN
     669        ALLOCATE(flight_dist_mpi(klon_mpi, nleva,1), STAT=ierr)
     670        IF (ierr /= 0) CALL abort_physic('read_aviation_emissions', 'problem to allocate flight_dist_mpi',1)
     671        CALL xios_recv_field("KMFLOWN_interp", flight_dist_mpi(:,:,1))
     672        ! Get number of vertical levels and level values
     673        CALL xios_get_axis_attr( "aviation_lev", value=aviation_lev(:))
     674    ENDIF
     675
     676    ! Propagate to other OMP threads: flight_dist_mpi(klon_mpi,klev) to flight_dist(klon,klev)
     677    ! (klon_mpi,klon) = (200,50) avec 80 MPI, 4 OMP, nbp40
     678    CALL scatter_omp(flight_dist_mpi, flight_dist_read)
     679    CALL bcast_omp(aviation_lev)
    683680
    684681END SUBROUTINE read_aviation_emissions
    685682
    686 SUBROUTINE vertical_interpolation_aviation(klon, klev, paprs, pplay, flight_dist_read, &
    687                                      flight_h2o_read, aviation_lev, nleva, flight_dist, &
    688                                      flight_h2o)
     683SUBROUTINE vertical_interpolation_aviation(klon, klev, paprs, pplay, temp, flight_dist, flight_h2o)
    689684    ! This subroutine performs the vertical interpolation from the read data in aviation.nc
    690685    ! where there are nleva vertical levels described in aviation_lev to the klev levels or
    … …  
    693688    ! flight_h2o_read(klon,nleva) -> flight_h2o(klon, klev)
    694689    USE print_control_mod, ONLY: lunout
     690    USE lmdz_lscp_ini, ONLY: RD, RG
     691
    695692    IMPLICIT NONE
    696693
    … …  
    698695    REAL, INTENT(IN)    :: paprs(klon, klev+1) ! inter-layer pressure [Pa]
    699696    REAL, INTENT(IN)    :: pplay(klon, klev) ! mid-layer pressure [Pa]
    700     INTEGER, INTENT(IN) :: nleva  ! Size of the vertical axis in the file
    701     REAL, INTENT(OUT) :: flight_dist(klon,klev,1) ! Aviation distance flown within the mesh [m/s/mesh]
    702     REAL, INTENT(OUT) :: flight_h2o(klon,klev,1)  ! Aviation H2O emitted within the mesh [kgH2O/s/mesh]
    703     REAL, INTENT(IN) :: flight_dist_read(klon,nleva,1) ! Aviation distance flown within the mesh in file [m/s/mesh]
    704     REAL, INTENT(IN) :: flight_h2o_read(klon,nleva,1)  ! Aviation H2O emitted within the mesh in file [kgH2O/s/mesh]
    705     REAL, INTENT(IN) :: aviation_lev(nleva)  !  Pressure in the middle of the layers [Pa]
     697    REAL, INTENT(IN)    :: temp(klon, klev) ! temperature [K]
     698    REAL, INTENT(OUT)   :: flight_dist(klon,klev) ! Aviation distance flown within the mesh [m/s/mesh]
     699    REAL, INTENT(OUT)   :: flight_h2o(klon,klev)  ! Aviation H2O emitted within the mesh [kgH2O/s/mesh]
    706700
    707701    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    … …  
    713707    REAL :: zfrac ! Fraction of layer kori in layer k
    714708    REAL :: width_read_layer(1:nleva) ! width of a given layer [ Pa ]
     709    REAL :: rho, rhodz, dz
    715710
    716711    ! Initialisation
    717     flight_dist(:,:,1) = 0.
    718     flight_h2o(:,:,1) = 0.
     712    flight_dist(:,:) = 0.
     713    flight_h2o(:,:) = 0.
    719714
    720715    ! Compute the array with the vertical interface
    … …  
    751746                 
    752747                 ! Vertical reprojection for each desired array
    753                  flight_dist(i,k,1) = flight_dist(i,k,1) + zfrac * flight_dist_read(i,kori,1)
    754                  flight_h2o(i,k,1)  = flight_h2o(i,k,1) + zfrac * flight_h2o(i,kori,1)
     748                 flight_dist(i,k) = flight_dist(i,k) + zfrac * flight_dist_read(i,kori,1)
     749                 flight_h2o(i,k)  = flight_h2o(i,k) + zfrac * flight_h2o_read(i,kori,1)
    755750            ENDDO
     751
     752            !--Dry density [kg/m3]
     753            rho = pplay(i,k) / temp(i,k) / RD
     754            !--Dry air mass [kg/m2]
     755            rhodz = ( paprs(i,k) - paprs(i,k+1) ) / RG
     756            !--Cell thickness [m]
     757            dz = rhodz / rho
     758
     759            !--Normalisation with the cell thickness
     760            flight_dist(i,k) = flight_dist(i,k) / dz
     761            flight_h2o(i,k) = flight_h2o(i,k) / dz
    756762        ENDDO
    757763    ENDDO

  • LMDZ6/branches/contrails/libf/phylmd/lmdz_lscp.f90

    r5573 r5575  
    1818     iflag_ice_thermo, distcltop, temp_cltop,           &
    1919     tke, tke_dissip,                                   &
    20      cell_area,                                         &
     20     cell_area, stratomask,                             &
    2121     cf_seri, rvc_seri, u_seri, v_seri,                 &
    2222     qsub, qissr, qcld, subfra, issrfra, gamma_cond,    &
    … …  
    2525     dqvc_sub, dqvc_con, dqvc_mix, qsatl, qsati,        &
    2626     rcont_seri, flight_dist, flight_h2o,               &
    27      flight_dist_read, flight_h2o_read,                 &
    28      aviation_lev, nleva, contfra, Tcritcont, qcritcont,&
     27     contfra, Tcritcont, qcritcont,                     &
    2928     potcontfraP, potcontfraNP, dcontfra_cir, dcf_avi,  &
    3029     dqi_avi, dqvc_avi, cloudth_sth,cloudth_senv,       &
    … …  
    124123USE lmdz_lscp_ini, ONLY : ok_plane_contrail
    125124
    126 ! aviation module
    127 USE mod_phys_lmdz_para, ONLY : is_omp_master
    128 USE lmdz_aviation, ONLY : vertical_interpolation_aviation
     125! Temporary call for Lamquin et al (2012) diagnostics
     126USE phys_local_var_mod, ONLY : issrfra100to150, issrfra150to200, issrfra200to250
     127USE phys_local_var_mod, ONLY : issrfra250to300, issrfra300to400, issrfra400to500
    129128
    130129IMPLICIT NONE
    … …  
    177176  REAL, DIMENSION(klon,klev),      INTENT(IN)   :: v_seri           ! northward wind [m/s]
    178177  REAL, DIMENSION(klon),           INTENT(IN)   :: cell_area        ! area of each cell [m2]
     178  REAL, DIMENSION(klon,klev),      INTENT(IN)   :: stratomask       ! fraction of stratosphere (0 or 1)
    179179
    180180  ! INPUT/OUTPUT aviation
    181181  !--------------------------------------------------
    182182  REAL, DIMENSION(klon,klev),      INTENT(INOUT):: rcont_seri       ! ratio of contrails fraction to total cloud fraction [-]
    183   REAL, DIMENSION(klon,klev,1),      INTENT(OUT)   :: flight_dist      ! aviation distance flown within the mesh [m/s/mesh]
    184   REAL, DIMENSION(klon,klev,1),      INTENT(OUT)   :: flight_h2o       ! aviation H2O emitted within the mesh [kgH2O/s/mesh]
    185   ! Read from file - nleva vertical levels - define in state_var
    186   INTEGER,                         INTENT(IN)   :: nleva         
    187   REAL, DIMENSION(nleva),          INTENT(IN)   :: aviation_lev     ! vertical levels [km]
    188   REAL, DIMENSION(klon,nleva,1),     INTENT(IN)   :: flight_dist_read  ! aviation distance flown within the mesh [m/s/mesh]
    189   REAL, DIMENSION(klon,nleva,1),     INTENT(IN)   :: flight_h2o_read       ! aviation H2O emitted within the mesh [kgH2O/s/mesh]
     183  REAL, DIMENSION(klon,klev),      INTENT(IN)   :: flight_dist      ! aviation distance flown within the mesh [m/s/mesh]
     184  REAL, DIMENSION(klon,klev),      INTENT(IN)   :: flight_h2o       ! aviation H2O emitted within the mesh [kgH2O/s/mesh]
    190185 
    191186  ! OUTPUT variables
    … …  
    331326  ! dyn3d_common/infotrac.F90
    332327  REAL :: min_qParent, min_ratio
     328  !--for Lamquin et al 2012 diagnostics
     329  REAL, DIMENSION(klon) :: issrfra100to150UP, issrfra150to200UP, issrfra200to250UP
     330  REAL, DIMENSION(klon) :: issrfra250to300UP, issrfra300to400UP, issrfra400to500UP
    333331
    334332  INTEGER i, k, kk, iter
    … …  
    442440min_ratio       = 1.e-16
    443441
     442!--for Lamquin et al (2012) diagnostics
     443issrfra100to150(:)   = 0.
     444issrfra100to150UP(:) = 0.
     445issrfra150to200(:)   = 0.
     446issrfra150to200UP(:) = 0.
     447issrfra200to250(:)   = 0.
     448issrfra200to250UP(:) = 0.
     449issrfra250to300(:)   = 0.
     450issrfra250to300UP(:) = 0.
     451issrfra300to400(:)   = 0.
     452issrfra300to400UP(:) = 0.
     453issrfra400to500(:)   = 0.
     454issrfra400to500UP(:) = 0.
     455
    444456!-- poprecip
    445457qraindiag(:,:)= 0.
    … …  
    462474
    463475!c_iso: variable initialisation for iso
    464 
    465 IF ( ok_plane_contrail ) THEN
    466     ! Vertical interpolation is done at each physical timestep
    467     !
    468     IF (is_omp_master) CALL vertical_interpolation_aviation(klon, klev, paprs, pplay, flight_dist_read, &
    469                    flight_h2o_read, aviation_lev, nleva, flight_dist, flight_h2o)
    470 ENDIF
    471476
    472477!===============================================================================
    … …  
    742747                        pplay(:,k), paprs(:,k), paprs(:,k+1), &
    743748                        cf_seri(:,k), rvc_seri(:,k), ql_seri(:,k), qi_seri(:,k), &
    744                         shear, tke_dissip(:,k), cell_area, &
     749                        shear, tke_dissip(:,k), cell_area, stratomask(:,k), &
    745750                        Tbef, zq, zqs, gammasat, ratqs(:,k), keepgoing, &
    746751                        rneb(:,k), zqn, qvc, issrfra(:,k), qissr(:,k), &
    … …  
    748753                        dqi_adj(:,k), dqi_sub(:,k), dqi_con(:,k), dqi_mix(:,k), &
    749754                        dqvc_adj(:,k), dqvc_sub(:,k), dqvc_con(:,k), dqvc_mix(:,k), &
    750                         rcont_seri(:,k), flight_dist(:,k,1), flight_h2o(:,k,1), contfra(:,k), &
     755                        rcont_seri(:,k), flight_dist(:,k), flight_h2o(:,k), contfra(:,k), &
    751756                        Tcritcont(:,k), qcritcont(:,k), potcontfraP(:,k), potcontfraNP(:,k), &
    752757                        dcontfra_cir(:,k), dcf_avi(:,k), dqi_avi(:,k), dqvc_avi(:,k))
    … …  
    11091114        qsub(i,k) = zq(i) - qvc(i) - qissr(i,k)
    11101115        qcld(i,k) = qvc(i) + zoliq(i)
     1116
     1117        !--Calculation of the ice supersaturated fraction following Lamquin et al (2012)
     1118        !--methodology: in each layer, we make a maximum random overlap assumption for
     1119        !--ice supersaturation
     1120        IF ( ( paprs(i,k) .GT. 10000. ) .AND. ( paprs(i,k) .LE. 15000. ) ) THEN
     1121                IF ( issrfra100to150UP(i) .GT. ( 1. - eps ) ) THEN
     1122                        issrfra100to150(i) = 1.
     1123                ELSE
     1124                        issrfra100to150(i) = 1. - ( 1. - issrfra100to150(i) ) * &
     1125                                ( 1. - MAX( issrfra(i,k), issrfra100to150UP(i) ) ) &
     1126                              / ( 1. - issrfra100to150UP(i) )
     1127                        issrfra100to150UP(i) = issrfra(i,k)
     1128                ENDIF
     1129        ELSEIF ( ( paprs(i,k) .GT. 15000. ) .AND. ( paprs(i,k) .LE. 20000. ) ) THEN
     1130                IF ( issrfra150to200UP(i) .GT. ( 1. - eps ) ) THEN
     1131                        issrfra150to200(i) = 1.
     1132                ELSE
     1133                        issrfra150to200(i) = 1. - ( 1. - issrfra150to200(i) ) * &
     1134                                ( 1. - MAX( issrfra(i,k), issrfra150to200UP(i) ) ) &
     1135                              / ( 1. - issrfra150to200UP(i) )
     1136                        issrfra150to200UP(i) = issrfra(i,k)
     1137                ENDIF
     1138        ELSEIF ( ( paprs(i,k) .GT. 20000. ) .AND. ( paprs(i,k) .LE. 25000. ) ) THEN
     1139                IF ( issrfra200to250UP(i) .GT. ( 1. - eps ) ) THEN
     1140                        issrfra200to250(i) = 1.
     1141                ELSE
     1142                        issrfra200to250(i) = 1. - ( 1. - issrfra200to250(i) ) * &
     1143                                ( 1. - MAX( issrfra(i,k), issrfra200to250UP(i) ) ) &
     1144                              / ( 1. - issrfra200to250UP(i) )
     1145                        issrfra200to250UP(i) = issrfra(i,k)
     1146                ENDIF
     1147        ELSEIF ( ( paprs(i,k) .GT. 25000. ) .AND. ( paprs(i,k) .LE. 30000. ) ) THEN
     1148                IF ( issrfra250to300UP(i) .GT. ( 1. - eps ) ) THEN
     1149                        issrfra250to300(i) = 1.
     1150                ELSE
     1151                        issrfra250to300(i) = 1. - ( 1. - issrfra250to300(i) ) * &
     1152                                ( 1. - MAX( issrfra(i,k), issrfra250to300UP(i) ) ) &
     1153                              / ( 1. - issrfra250to300UP(i) )
     1154                        issrfra250to300UP(i) = issrfra(i,k)
     1155                ENDIF
     1156        ELSEIF ( ( paprs(i,k) .GT. 30000. ) .AND. ( paprs(i,k) .LE. 40000. ) ) THEN
     1157                IF ( issrfra300to400UP(i) .GT. ( 1. - eps ) ) THEN
     1158                        issrfra300to400(i) = 1.
     1159                ELSE
     1160                        issrfra300to400(i) = 1. - ( 1. - issrfra300to400(i) ) * &
     1161                                ( 1. - MAX( issrfra(i,k), issrfra300to400UP(i) ) ) &
     1162                              / ( 1. - issrfra300to400UP(i) )
     1163                        issrfra300to400UP(i) = issrfra(i,k)
     1164                ENDIF
     1165        ELSEIF ( ( paprs(i,k) .GT. 40000. ) .AND. ( paprs(i,k) .LE. 50000. ) ) THEN
     1166                IF ( issrfra400to500UP(i) .GT. ( 1. - eps ) ) THEN
     1167                        issrfra400to500(i) = 1.
     1168                ELSE
     1169                        issrfra400to500(i) = 1. - ( 1. - issrfra400to500(i) ) * &
     1170                                ( 1. - MAX( issrfra(i,k), issrfra400to500UP(i) ) ) &
     1171                              / ( 1. - issrfra400to500UP(i) )
     1172                        issrfra400to500UP(i) = issrfra(i,k)
     1173                ENDIF
     1174        ENDIF
     1175
    11111176      ENDDO
    11121177    ENDIF

  • LMDZ6/branches/contrails/libf/phylmd/lmdz_lscp_condensation.f90

    r5551 r5575  
    9595SUBROUTINE condensation_ice_supersat( &
    9696      klon, dtime, missing_val, pplay, paprsdn, paprsup, &
    97       cf_seri, rvc_seri, ql_seri, qi_seri, shear, pbl_eps, cell_area, &
     97      cf_seri, rvc_seri, ql_seri, qi_seri, shear, pbl_eps, cell_area, stratomask, &
    9898      temp, qtot, qsat, gamma_cond, ratqs, keepgoing, &
    9999      cldfra, qincld, qvc, issrfra, qissr, dcf_sub, dcf_con, dcf_mix, &
    … …  
    120120USE lmdz_lscp_ini,        ONLY: eps, temp_nowater, ok_weibull_warm_clouds
    121121USE lmdz_lscp_ini,        ONLY: ok_unadjusted_clouds, iflag_cloud_sublim_pdf
    122 USE lmdz_lscp_ini,        ONLY: ok_plane_contrail
     122USE lmdz_lscp_ini,        ONLY: ok_plane_contrail, ok_no_issr_strato
    123123
    124124USE lmdz_lscp_ini,        ONLY: depo_coef_cirrus, capa_cond_cirrus, std_subl_pdf_lscp, &
    … …  
    148148REAL,     INTENT(IN)   , DIMENSION(klon) :: pbl_eps       ! TKE dissipation [m2/s3]
    149149REAL,     INTENT(IN)   , DIMENSION(klon) :: cell_area     ! cell area [m2]
     150REAL,     INTENT(IN)   , DIMENSION(klon) :: stratomask    ! fraction of stratosphere in the mesh (1 or 0)
    150151REAL,     INTENT(IN)   , DIMENSION(klon) :: temp          ! temperature [K]
    151152REAL,     INTENT(IN)   , DIMENSION(klon) :: qtot          ! total specific humidity (without precip) [kg/kg]
    … …  
    266267    !--If ok_weibull_warm_clouds = .TRUE., the Weibull law is used for
    267268    !--all clouds, and the lognormal scheme is not activated
    268     IF ( ( temp(i) .GT. temp_nowater ) .AND. .NOT. ok_weibull_warm_clouds ) THEN
     269    IF ( ( ( temp(i) .GT. temp_nowater ) .AND. .NOT. ok_weibull_warm_clouds ) .OR. &
     270         ( ok_no_issr_strato .AND. ( stratomask(i) .EQ. 1. ) ) ) THEN
    269271
    270272      pdf_std   = ratqs(i) * qtot(i)

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

    r5453 r5575  
    150150  !$OMP THREADPRIVATE(ok_ice_supersat)
    151151
     152  LOGICAL, SAVE, PROTECTED :: ok_no_issr_strato=.FALSE.      ! deactivates the ice supersaturation scheme in the stratosphere
     153  !$OMP THREADPRIVATE(ok_no_issr_strato)
     154
    152155  LOGICAL, SAVE, PROTECTED :: ok_unadjusted_clouds=.FALSE.   ! if True, relax the saturation adjustment assumption for ice clouds
    153156  !$OMP THREADPRIVATE(ok_unadjusted_clouds)
    … …  
    339342CONTAINS
    340343
    341 SUBROUTINE lscp_ini(dtime, lunout_in, prt_level_in, ok_ice_supersat_in, ok_plane_contrail_in, &
     344SUBROUTINE lscp_ini(dtime, lunout_in, prt_level_in, ok_ice_supersat_in, &
     345                    ok_no_issr_strato_in, ok_plane_contrail_in, &
    342346                    iflag_ratqs, fl_cor_ebil_in, &
    343347                    RCPD_in, RLSTT_in, RLVTT_in, RLMLT_in, RVTMP2_in, &
    … …  
    350354   REAL, INTENT(IN)      :: dtime
    351355   INTEGER, INTENT(IN)   :: lunout_in,prt_level_in,iflag_ratqs,fl_cor_ebil_in
    352    LOGICAL, INTENT(IN)   :: ok_ice_supersat_in, ok_plane_contrail_in
     356   LOGICAL, INTENT(IN)   :: ok_ice_supersat_in, ok_no_issr_strato_in, ok_plane_contrail_in
    353357
    354358   REAL, INTENT(IN)      :: RCPD_in, RLSTT_in, RLVTT_in, RLMLT_in
    … …  
    363367
    364368    ok_ice_supersat=ok_ice_supersat_in
     369    ok_no_issr_strato=ok_no_issr_strato_in
    365370    ok_plane_contrail=ok_plane_contrail_in
    366371

  • LMDZ6/branches/contrails/libf/phylmd/phys_local_var_mod.F90

    r5536 r5575  
    663663      REAL, SAVE, ALLOCATABLE :: qsatliq(:,:), qsatice(:,:)
    664664      !$OMP THREADPRIVATE(qsatliq, qsatice)
     665      REAL, SAVE, ALLOCATABLE :: issrfra100to150(:), issrfra150to200(:), issrfra200to250(:)
     666      !$OMP THREADPRIVATE(issrfra100to150, issrfra150to200, issrfra200to250)
     667      REAL, SAVE, ALLOCATABLE :: issrfra250to300(:), issrfra300to400(:), issrfra400to500(:)
     668      !$OMP THREADPRIVATE(issrfra250to300, issrfra300to400, issrfra400to500)
    665669
    666670!-- LSCP - aviation and contrails variables
    … …  
    12171221      ALLOCATE(dqvc_adj(klon,klev), dqvc_sub(klon,klev), dqvc_con(klon,klev), dqvc_mix(klon,klev))
    12181222      ALLOCATE(qsatliq(klon,klev), qsatice(klon,klev))
     1223      ALLOCATE(issrfra100to150(klon), issrfra150to200(klon), issrfra200to250(klon))
     1224      ALLOCATE(issrfra250to300(klon), issrfra300to400(klon), issrfra400to500(klon))
    12191225
    12201226!-- LSCP - aviation and contrails variables
    … …  
    16251631      DEALLOCATE(dqvc_adj, dqvc_sub, dqvc_con, dqvc_mix)
    16261632      DEALLOCATE(qsatliq, qsatice)
     1633      DEALLOCATE(issrfra100to150, issrfra150to200, issrfra200to250)
     1634      DEALLOCATE(issrfra250to300, issrfra300to400, issrfra400to500)
    16271635
    16281636!-- LSCP - aviation and contrails variables

  • LMDZ6/branches/contrails/libf/phylmd/phys_output_ctrlout_mod.F90

    r5573 r5575  
    21662166  TYPE(ctrl_out), SAVE :: o_qsati = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    21672167    'qsati', 'Saturation with respect to ice', 'kg/kg', (/ ('', i=1, 10) /))
     2168  TYPE(ctrl_out), SAVE :: o_issrfra100to150 = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2169    'issrfra100to150', 'Supersaturated fraction in the 100to150 hPa layer', '-', (/ ('', i=1, 10) /))
     2170  TYPE(ctrl_out), SAVE :: o_issrfra150to200 = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2171    'issrfra150to200', 'Supersaturated fraction in the 150to200 hPa layer', '-', (/ ('', i=1, 10) /))
     2172  TYPE(ctrl_out), SAVE :: o_issrfra200to250 = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2173    'issrfra200to250', 'Supersaturated fraction in the 200to250 hPa layer', '-', (/ ('', i=1, 10) /))
     2174  TYPE(ctrl_out), SAVE :: o_issrfra250to300 = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2175    'issrfra250to300', 'Supersaturated fraction in the 250to300 hPa layer', '-', (/ ('', i=1, 10) /))
     2176  TYPE(ctrl_out), SAVE :: o_issrfra300to400 = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2177    'issrfra300to400', 'Supersaturated fraction in the 300to400 hPa layer', '-', (/ ('', i=1, 10) /))
     2178  TYPE(ctrl_out), SAVE :: o_issrfra400to500 = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2179    'issrfra400to500', 'Supersaturated fraction in the 400to500 hPa layer', '-', (/ ('', i=1, 10) /))
    21682180
    21692181!-- LSCP - aviation variables

  • LMDZ6/branches/contrails/libf/phylmd/phys_output_write_mod.F90

    r5536 r5575  
    225225         o_dcfsub, o_dcfcon, o_dcfmix, o_dqiadj, o_dqisub, o_dqicon, o_dqimix, &
    226226         o_dqvcadj, o_dqvcsub, o_dqvccon, o_dqvcmix, o_qsatl, o_qsati, &
     227         o_issrfra100to150, o_issrfra150to200, o_issrfra200to250, &
     228         o_issrfra250to300, o_issrfra300to400, o_issrfra400to500, &
    227229!-- LSCP - aviation variables
    228230         o_rcontseri, o_drcontdyn, o_dqavi, o_contfra, &
    … …  
    350352         dqvc_adj, dqvc_sub, dqvc_con, dqvc_mix, &
    351353         qsatliq, qsatice, &
     354         issrfra100to150, issrfra150to200, issrfra200to250, &
     355         issrfra250to300, issrfra300to400, issrfra400to500, &
    352356         rcont_seri, d_rcont_dyn, d_q_avi, contfra, &
    353357         Tcritcont, qcritcont, potcontfraP, potcontfraNP, &
    … …  
    21222126         CALL histwrite_phy(o_qsatl, qsatliq)
    21232127         CALL histwrite_phy(o_qsati, qsatice)
     2128         CALL histwrite_phy(o_issrfra100to150, issrfra100to150)
     2129         CALL histwrite_phy(o_issrfra150to200, issrfra150to200)
     2130         CALL histwrite_phy(o_issrfra200to250, issrfra200to250)
     2131         CALL histwrite_phy(o_issrfra250to300, issrfra250to300)
     2132         CALL histwrite_phy(o_issrfra300to400, issrfra300to400)
     2133         CALL histwrite_phy(o_issrfra400to500, issrfra400to500)
    21242134       ENDIF
    21252135!-- LSCP - aviation variables

  • LMDZ6/branches/contrails/libf/phylmd/phys_state_var_mod.F90

    r5573 r5575  
    528528      !$OMP THREADPRIVATE(delta_sal, ds_ns, dt_ns, delta_sst, dter, dser, dt_ds)
    529529
    530       ! Arrays for the lecture of aviation files
    531       ! The allocation is done in the read_aviation module
    532       ! The size is (klon, nleva, 1) where
    533       ! nleva            is the size of the vertical axis (read from file)
    534       ! flight_dist_read is the number of km per second
    535       ! flight_h2o_read  is the water content added to the air
    536       ! aviation_lev     is the value of the levels
    537       REAL, SAVE, ALLOCATABLE :: flight_dist_read(:,:,:), flight_h2o_read(:,:,:)
    538       REAL, SAVE, ALLOCATABLE :: aviation_lev(:)
    539       !$OMP THREADPRIVATE(flight_dist_read, flight_h2o_read, aviation_lev)
    540       INTEGER, SAVE :: nleva
    541       !$OMP THREADPRIVATE(nleva)
    542 
    543530
    544531    CONTAINS
    … …  
    956943      DEALLOCATE(ratqs_inter_,sigma_qtherm)
    957944
    958       ! DEALLOCATE aviation arrays
    959       DEALLOCATE(flight_dist_read, flight_h2o_read, aviation_lev)
    960 
    961945      if (activate_ocean_skin >= 1) then
    962946         deALLOCATE(delta_sal, ds_ns, dt_ns, delta_sst, dter, dser)

  • LMDZ6/branches/contrails/libf/phylmd/physiq_mod.F90

    r5573 r5575  
    7575    USE write_field_phy
    7676    use wxios_mod, ONLY: g_ctx, wxios_set_context
    77     USE lmdz_aviation, ONLY : read_aviation_emissions, aviation_water_emissions
     77    USE lmdz_aviation, ONLY : read_aviation_emissions, aviation_water_emissions, vertical_interpolation_aviation
    7878    USE lmdz_lscp, ONLY : lscp
    7979    USE lmdz_call_cloud_optics_prop, ONLY : call_cloud_optics_prop
    … …  
    331331       contfra, Tcritcont, qcritcont, potcontfraP, potcontfraNP, &
    332332       dcontfra_cir, dcf_avi, dqi_avi, dqvc_avi, &
     333       !
     334       stratomask, &
    333335       !
    334336       cldemi,  &
    … …  
    18531855   &    RG,RD,RCPD,RKAPPA,RLVTT,RETV)
    18541856       CALL ratqs_ini(klon,klev,iflag_thermals,lunout,nbsrf,is_lic,is_ter,RG,RV,RD,RCPD,RLSTT,RLVTT,RTT)
    1855        CALL lscp_ini(pdtphys,lunout,prt_level,ok_ice_supersat,&
     1857       CALL lscp_ini(pdtphys,lunout,prt_level,ok_ice_supersat,ok_no_issr_strato,&
    18561858                     ok_plane_contrail,iflag_ratqs,fl_cor_ebil, &
    18571859                     RCPD,RLSTT,RLVTT,RLMLT,RVTMP2,RTT,RD,RV,RG,RPI,EPS_W)
    … …  
    21232125
    21242126         !--Read the aviation emissions
    2125          IF ( ok_plane_h2o .OR. ok_plane_contrail ) THEN
    2126            CALL read_aviation_emissions(klon, klev, flight_dist_read, flight_h2o_read, &
    2127                                         aviation_lev, nleva)
    2128          ENDIF
     2127         IF ( ok_plane_h2o .OR. ok_plane_contrail ) CALL read_aviation_emissions(klon, klev)
    21292128       !
    21302129       ! Now we activate some double radiation call flags only if some
    … …  
    38923891
    38933892    IF (ok_new_lscp) THEN
    3894  
     3893
     3894    IF (ok_no_issr_strato) CALL stratosphere_mask(missing_val, pphis, t_seri, pplay, latitude_deg)
     3895    ! Vertical interpolation is done at each physical timestep
     3896    IF (ok_plane_contrail) CALL vertical_interpolation_aviation(klon, klev, paprs, pplay, t_seri, flight_dist, flight_h2o)
     3897
    38953898    DO k = 1, klev
    38963899      DO i = 1, klon
    … …  
    39103913         iflag_ice_thermo, distcltop, temp_cltop,   &
    39113914         pbl_tke(:,:,is_ave), pbl_eps(:,:,is_ave), &
    3912          cell_area, &
     3915         cell_area, stratomask, &
    39133916         cf_seri, rvc_seri, u_seri, v_seri, &
    39143917         qsub, qissr, qcld, subfra, issrfra, gamma_cond,  &
    39153918         dcf_sub, dcf_con, dcf_mix, dqi_adj, dqi_sub, dqi_con, dqi_mix, &
    39163919         dqvc_adj, dqvc_sub, dqvc_con, dqvc_mix, qsatliq, qsatice, &
    3917          rcont_seri, flight_dist, flight_h2o, flight_dist_read, flight_h2o_read, &
    3918          aviation_lev, nleva, contfra, Tcritcont, qcritcont, potcontfraP, &
     3920         rcont_seri, flight_dist, flight_h2o, &
     3921         contfra, Tcritcont, qcritcont, potcontfraP, &
    39193922         potcontfraNP, dcontfra_cir, dcf_avi, dqi_avi, dqvc_avi, &
    39203923         cloudth_sth,cloudth_senv,cloudth_sigmath,cloudth_sigmaenv, &
Note: See TracChangeset for help on using the changeset viewer.