Changeset 5641


Ignore:
Timestamp:
May 1, 2025, 6:00:03 PM (5 weeks ago)
Author:
aborella
Message:

Changed the contrails classes that are advected

Location:
LMDZ6/branches/contrails
Files:
17 edited

Legend:

Unmodified
Added
Removed

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

    r5618 r5641  
    912912        <field id="qsati"      long_name="Saturation with respect to ice"    unit="kg/kg" />
    913913
    914         <field id="cfaseri"    long_name="Linear contrail fraction" unit="-" />
    915         <field id="dcfadyn"    long_name="Dynamics linear contrail fraction tendency" unit="s-1" />
    916         <field id="pcfseri"    long_name="Contrail induced cirrus fraction" unit="-" />
    917         <field id="dpcfdyn"    long_name="Dynamics contrail induced cirrus fraction tendency" unit="s-1" />
    918         <field id="qvaseri"    long_name="Linear contrail vapor specific humidity" unit="kg/kg" />
    919         <field id="dqvadyn"    long_name="Dynamics linear contrail vapor specific humidity tendency" unit="kg/kg/s" />
    920         <field id="qiaseri"    long_name="Linear contrail ice specific humidity" unit="kg/kg" />
    921         <field id="dqiadyn"    long_name="Dynamics linear contrail ice specific humidity tendency" unit="kg/kg/s" />
     914        <field id="cflseri"    long_name="Linear contrail fraction" unit="-" />
     915        <field id="dcfldyn"    long_name="Dynamics linear contrail fraction tendency" unit="s-1" />
     916        <field id="cfcseri"    long_name="Contrail induced cirrus fraction" unit="-" />
     917        <field id="dcfcdyn"    long_name="Dynamics contrail induced cirrus fraction tendency" unit="s-1" />
     918        <field id="qtlseri"    long_name="Linear contrail total specific humidity" unit="kg/kg" />
     919        <field id="dqtldyn"    long_name="Dynamics linear contrail total specific humidity tendency" unit="kg/kg/s" />
     920        <field id="qtcseri"    long_name="Contrail cirrus total specific humidity" unit="kg/kg" />
     921        <field id="dqtcdyn"    long_name="Dynamics contrail cirrus total specific humidity tendency" unit="kg/kg/s" />
    922922        <field id="Tcritcont"  long_name="Temperature threshold for contrail formation"    unit="K" />
    923923        <field id="qcritcont"  long_name="Specific humidity threshold for contrail formation"    unit="kg/kg" />
    924924        <field id="potcontfraP"  long_name="Fraction with pontential persistent contrail"    unit="-" />
    925925        <field id="potcontfraNP" long_name="Fraction with potential non-persistent contrail"    unit="-" />
    926         <field id="qice_perscont" long_name="Contrail induced cirrus ice specific humidity"    unit="kg/kg" />
    927         <field id="dcfaini"    long_name="Initial formation linear contrail fraction tendency"    unit="s-1" />
    928         <field id="dqiaini"    long_name="Initial formation linear contrail ice specific humidity tendency"    unit="kg/kg/s" />
    929         <field id="dqtaini"    long_name="Initial formation linear contrail total specific humidity tendency"    unit="kg/kg/s" />
    930         <field id="dcfacir"    long_name="Conversion to cirrus linear contrail fraction tendency"    unit="s-1" />
    931         <field id="dqtacir"    long_name="Conversion to cirrus linear contrail total specific humidity tendency"    unit="kg/kg/s" />
    932         <field id="dcfasub"    long_name="Sublimation linear contrail fraction tendency"    unit="s-1" />
    933         <field id="dqiasub"    long_name="Sublimation linear contrail ice specific humidity tendency"    unit="kg/kg/s" />
    934         <field id="dqtasub"    long_name="Sublimation linear contrail total specific humidity tendency"    unit="kg/kg/s" />
    935         <field id="dcfamix"    long_name="Mixing linear contrail fraction tendency"    unit="s-1" />
    936         <field id="dqiamix"    long_name="Mixing linear contrail ice specific humidity tendency"    unit="kg/kg/s" />
    937         <field id="dqtamix"    long_name="Mixing linear contrail total specific humidity tendency"    unit="kg/kg/s" />
     926        <field id="qice_lincont" long_name="Linear contrails ice specific humidity"    unit="kg/kg" />
     927        <field id="qice_circont" long_name="Contrail induced cirrus ice specific humidity"    unit="kg/kg" />
     928        <field id="dcflini"    long_name="Initial formation linear contrail fraction tendency"    unit="s-1" />
     929        <field id="dqilini"    long_name="Initial formation linear contrail ice specific humidity tendency"    unit="kg/kg/s" />
     930        <field id="dqtlini"    long_name="Initial formation linear contrail total specific humidity tendency"    unit="kg/kg/s" />
     931        <field id="dcflcir"    long_name="Conversion to cirrus linear contrail fraction tendency"    unit="s-1" />
     932        <field id="dqtlcir"    long_name="Conversion to cirrus linear contrail total specific humidity tendency"    unit="kg/kg/s" />
     933        <field id="dcflsub"    long_name="Sublimation linear contrail fraction tendency"    unit="s-1" />
     934        <field id="dqilsub"    long_name="Sublimation linear contrail ice specific humidity tendency"    unit="kg/kg/s" />
     935        <field id="dqtlsub"    long_name="Sublimation linear contrail total specific humidity tendency"    unit="kg/kg/s" />
     936        <field id="dcflmix"    long_name="Mixing linear contrail fraction tendency"    unit="s-1" />
     937        <field id="dqilmix"    long_name="Mixing linear contrail ice specific humidity tendency"    unit="kg/kg/s" />
     938        <field id="dqtlmix"    long_name="Mixing linear contrail total specific humidity tendency"    unit="kg/kg/s" />
     939        <field id="dcfcsub"    long_name="Sublimation contrail cirrus fraction tendency"    unit="s-1" />
     940        <field id="dqicsub"    long_name="Sublimation contrail cirrus ice specific humidity tendency"    unit="kg/kg/s" />
     941        <field id="dqtcsub"    long_name="Sublimation contrail cirrus total specific humidity tendency"    unit="kg/kg/s" />
     942        <field id="dcfcmix"    long_name="Mixing contrail cirrus fraction tendency"    unit="s-1" />
     943        <field id="dqicmix"    long_name="Mixing contrail cirrus ice specific humidity tendency"    unit="kg/kg/s" />
     944        <field id="dqtcmix"    long_name="Mixing contrail cirrus total specific humidity tendency"    unit="kg/kg/s" />
    938945        <field id="flightdist" long_name="Aviation flown distance concentration"    unit="m/s/m3" />
    939946        <field id="flighth2o"  long_name="Aviation emitted H2O concentration"    unit="kg H2O/s/m3" />

  • LMDZ6/branches/contrails/libf/dynphy_lonlat/phylmd/etat0phys_netcdf.f90

    r5626 r5641  
    5050    ale_bl, ale_bl_trig, alp_bl, &
    5151    ale_wake, ale_bl_stat, AWAKE_S, &
    52     cf_ancien, qvc_ancien, qvcon, qccon, cfa_ancien, pcf_ancien, qva_ancien, qia_ancien
     52    cf_ancien, qvc_ancien, qvcon, qccon, cfl_ancien, cfc_ancien, qtl_ancien, qtc_ancien
    5353
    5454  USE comconst_mod, ONLY: pi, dtvr
    … …  
    245245  qvcon = 0.
    246246  qccon = 0.
    247   cfa_ancien = 0.
    248   pcf_ancien = 0.
    249   qva_ancien = 0.
    250   qia_ancien = 0.
     247  cfl_ancien = 0.
     248  cfc_ancien = 0.
     249  qtl_ancien = 0.
     250  qtc_ancien = 0.
    251251 
    252252  z0m(:,:)=0 ! ym missing 5th subsurface initialization

  • LMDZ6/branches/contrails/libf/misc/readTracFiles_mod.f90

    r5623 r5641  
    121121  !--- SOME PARAMETERS THAT ARE NOT LIKELY TO CHANGE OFTEN
    122122  CHARACTER(LEN=maxlen), SAVE      :: tran0        = 'air'      !--- Default transporting fluid
    123   CHARACTER(LEN=maxlen), PARAMETER :: old_phases   = 'vlibfcapqt'     !--- Old phases for water (no separator)
    124   CHARACTER(LEN=maxlen), PARAMETER :: known_phases = 'glsbfcapqt'     !--- Known phases initials
     123  CHARACTER(LEN=maxlen), PARAMETER :: old_phases   = 'vlibfczyxw'     !--- Old phases for water (no separator)
     124  CHARACTER(LEN=maxlen), PARAMETER :: known_phases = 'glsbfczyxw'     !--- Known phases initials
    125125  INTEGER, PARAMETER :: nphases = LEN_TRIM(known_phases)        !--- Number of phases
    126126  CHARACTER(LEN=maxlen), SAVE      :: phases_names(nphases) &   !--- Known phases names
    127                                 = ['gaseous  ', 'liquid   ', 'solid    ','blownSnow', 'fraccld  ', 'cldvap   ', 'avifrac  ', 'pavifra  ', 'qvapavi  ', 'ticeavi  ']
     127                                = ['gaseous  ', 'liquid   ', 'solid    ','blownSnow', 'fraccld  ', 'cldvap   ', 'zlincont ', 'ycircont ', 'xlincontw', 'wcircontw']
    128128  CHARACTER(LEN=1),      SAVE :: phases_sep  =  '_'             !--- Phase separator
    129129  CHARACTER(LEN=maxlen), SAVE :: isoFile = 'isotopes_params.def'!--- Name of the isotopes parameters file

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

    r5623 r5641  
    2727   !=== FOR ISOTOPES: Specific to water
    2828   PUBLIC :: iH2O                                          !--- Value of "ixIso" for "H2O" isotopes class
    29    PUBLIC :: ivap, iliq, isol, ibs, icf, iqvc, icfa, ipcf, iqia, iqva
     29   PUBLIC :: ivap, iliq, isol, ibs, icf, iqvc, icfl, icfc, iqtl, iqtc
    3030   !=== FOR ISOTOPES: Depending on the selected isotopes family
    3131   PUBLIC :: isotope                                       !--- Selected isotopes database (argument of getKey)
    … …  
    103103
    104104   !=== INDICES FOR WATER
    105    INTEGER, SAVE :: ivap, iliq, isol, ibs, icf, iqvc, icfa, ipcf, iqia, iqva
    106 !$OMP THREADPRIVATE(ivap, iliq, isol, ibs, icf, iqvc, icfa, ipcf, iqia, iqva)
     105   INTEGER, SAVE :: ivap, iliq, isol, ibs, icf, iqvc, icfl, icfc, iqtl, iqtc
     106!$OMP THREADPRIVATE(ivap, iliq, isol, ibs, icf, iqvc, icfl, icfc, iqtl, iqtc)
    107107
    108108   !=== DIMENSIONS OF THE TRACERS TABLES AND OTHER SCALAR VARIABLES
    … …  
    364364   ibs  = strIdx(tracers(:)%name, addPhase('H2O', 'b'))
    365365   icf  = strIdx(tracers(:)%name, 'CLDFRA')
    366    iqvc = strIdx(tracers(:)%name, 'CLDVAP_g')
    367    icfa = strIdx(tracers(:)%name, 'CONTFRA')
    368    ipcf = strIdx(tracers(:)%name, 'PERSCONTFRA')
    369    iqva = strIdx(tracers(:)%name, 'CONTWATER_g')
    370    iqia = strIdx(tracers(:)%name, 'CONTWATER_s')
     366   iqvc = strIdx(tracers(:)%name, 'CLDVAP')
     367   icfl = strIdx(tracers(:)%name, 'LINCONTFRA')
     368   icfc = strIdx(tracers(:)%name, 'CIRCONTFRA')
     369   iqtl = strIdx(tracers(:)%name, 'LINCONTWATER')
     370   iqtc = strIdx(tracers(:)%name, 'CIRCONTWATER')
    371371   !--Two ways of declaring tracers - the way below should be deleted in the future
    372372   IF (icf.EQ.0) icf  = strIdx(tracers(:)%name, addPhase('H2O', 'f'))
    373373   IF (iqvc.EQ.0) iqvc = strIdx(tracers(:)%name, addPhase('H2O', 'c'))
    374    IF (icfa.EQ.0) icfa = strIdx(tracers(:)%name, addPhase('H2O', 'a'))
    375    IF (ipcf.EQ.0) ipcf = strIdx(tracers(:)%name, addPhase('H2O', 'p'))
    376    IF (iqva.EQ.0) iqva = strIdx(tracers(:)%name, addPhase('H2O', 'q'))
    377    IF (iqia.EQ.0) iqia = strIdx(tracers(:)%name, addPhase('H2O', 't'))
     374   IF (icfl.EQ.0) icfl = strIdx(tracers(:)%name, addPhase('H2O', 'z'))
     375   IF (icfc.EQ.0) icfc = strIdx(tracers(:)%name, addPhase('H2O', 'y'))
     376   IF (iqtl.EQ.0) iqtl = strIdx(tracers(:)%name, addPhase('H2O', 'x'))
     377   IF (iqtc.EQ.0) iqtc = strIdx(tracers(:)%name, addPhase('H2O', 'w'))
    378378
    379379   IF(CPPKEY_STRATAER .AND. type_trac == 'coag') THEN

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

    r5637 r5641  
    7575      clrfra, qclr, pdf_scale, pdf_alpha, pdf_gamma, keepgoing, pt_pron_clds, &
    7676      Tcritcont, qcritcont, potcontfraP, potcontfraNP, &
    77       dcfa_ini, dqia_ini, dqta_ini)
     77      dcfl_ini, dqil_ini, dqtl_ini)
    7878
    7979USE lmdz_lscp_ini, ONLY: RPI, RCPD, EPS_W, RTT
    … …  
    110110REAL,     INTENT(INOUT), DIMENSION(klon) :: potcontfraP  ! potential persistent contrail fraction [-]
    111111REAL,     INTENT(INOUT), DIMENSION(klon) :: potcontfraNP ! potential non-persistent contrail fraction [-]
    112 REAL,     INTENT(INOUT), DIMENSION(klon) :: dcfa_ini     ! contrails cloud fraction tendency bec ause of initial formation [s-1]
    113 REAL,     INTENT(INOUT), DIMENSION(klon) :: dqia_ini     ! contrails ice specific humidity tende ncy because of initial formation [kg/kg/s]
    114 REAL,     INTENT(INOUT), DIMENSION(klon) :: dqta_ini     ! contrails total specific humidity ten dency because of initial formation [kg/kg/s]
     112REAL,     INTENT(INOUT), DIMENSION(klon) :: dcfl_ini     ! linear contrails cloud fraction tendency bec ause of initial formation [s-1]
     113REAL,     INTENT(INOUT), DIMENSION(klon) :: dqil_ini     ! linear contrails ice specific humidity tende ncy because of initial formation [kg/kg/s]
     114REAL,     INTENT(INOUT), DIMENSION(klon) :: dqtl_ini     ! linear contrails total specific humidity ten dency because of initial formation [kg/kg/s]
    115115!
    116116! Local
    … …  
    217217      contrail_cross_section = CONTRAIL_CROSS_SECTION_ONERA()
    218218      contfra_new = MIN(1., flight_dist(i) * dtime * contrail_cross_section)
    219       dcfa_ini(i) = potcontfraP(i) * contfra_new
    220       dqta_ini(i) = qpotcontP * contfra_new
    221       dqia_ini(i) = dqta_ini(i) - qsat(i) * dcfa_ini(i)
     219      dcfl_ini(i) = potcontfraP(i) * contfra_new
     220      dqtl_ini(i) = qpotcontP * contfra_new
     221      dqil_ini(i) = dqtl_ini(i) - qsat(i) * dcfl_ini(i)
    222222    ENDIF
    223223

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

    r5609 r5641  
    1212    icefrac_optics, dNovrN, ptconv,rnebcon, ccwcon, &
    1313    !--AB contrails
    14     contfra, perscontfra, qice_cont, qice_perscont, pclc_nocont, &
     14    lincontfra, circontfra, qice_lincont, qice_circont, pclc_nocont, &
    1515    pcltau_nocont, pclemi_nocont, pcltau_cont, pclemi_cont, pch_nocont, pct_cont, &
    1616    xfiwp_nocont, xfiwc_nocont, reice_nocont)
    … …  
    9595
    9696  !--AB for contrails. All these are used / outputed only if ok_plane_contrail=y
    97   REAL, INTENT(IN)  :: contfra(klon, klev)       ! linear contrails fraction [-]
    98   REAL, INTENT(IN)  :: perscontfra(klon, klev)   ! contrail induced cirrus fraction [-]
    99   REAL, INTENT(IN)  :: qice_cont(klon, klev)     ! linear contrails condensed water [kg/kg]
    100   REAL, INTENT(IN)  :: qice_perscont(klon, klev) ! contrail induced cirrus condensed water [kg/kg]
     97  REAL, INTENT(IN)  :: lincontfra(klon, klev)    ! linear contrails fraction [-]
     98  REAL, INTENT(IN)  :: circontfra(klon, klev)    ! contrail induced cirrus fraction [-]
     99  REAL, INTENT(IN)  :: qice_lincont(klon, klev)  ! linear contrails condensed water [kg/kg]
     100  REAL, INTENT(IN)  :: qice_circont(klon, klev) ! contrail induced cirrus condensed water [kg/kg]
    101101  REAL, INTENT(OUT) :: pch_nocont(klon)          ! 2D high cloud cover without contrails[-]
    102102  REAL, INTENT(OUT) :: pct_cont(klon)            ! 2D total contrails cover[-]
    … …  
    139139    icefrac_optics, dNovrN, ptconv,rnebcon, ccwcon, &
    140140    !--AB for contrails
    141     contfra, perscontfra, qice_cont, qice_perscont, pclc_nocont, pcltau_nocont, &
     141    lincontfra, circontfra, qice_lincont, qice_circont, pclc_nocont, pcltau_nocont, &
    142142    pclemi_nocont, pcltau_cont, pclemi_cont, pch_nocont, pct_cont, &
    143143    xfiwp_nocont, xfiwc_nocont, reice_nocont)

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

    r5639 r5641  
    1111    icefrac_optics, dNovrN, ptconv, rnebcon, ccwcon, &
    1212    !--AB contrails
    13     contfra, perscontfra, qice_cont, qice_perscont, pclc_nocont, pcltau_nocont, &
     13    lincontfra, circontfra, qice_lincont, qice_circont, pclc_nocont, pcltau_nocont, &
    1414    pclemi_nocont, pcltau_cont, pclemi_cont, pch_nocont, pct_cont, &
    1515    xfiwp_nocont, xfiwc_nocont, reice_nocont)
    … …  
    119119
    120120  !--AB for contrails. All these are used / outputed only if ok_plane_contrail=y
    121   REAL, INTENT(IN)  :: contfra(klon, klev)       ! linear contrails fraction [-]
    122   REAL, INTENT(IN)  :: perscontfra(klon, klev)   ! contrail induced cirrus fraction [-]
    123   REAL, INTENT(IN)  :: qice_cont(klon, klev)     ! linear contrails condensed water [kg/kg]
    124   REAL, INTENT(IN)  :: qice_perscont(klon, klev) ! contrail induced cirrus condensed water [kg/kg]
     121  REAL, INTENT(IN)  :: lincontfra(klon, klev)    ! linear contrails fraction [-]
     122  REAL, INTENT(IN)  :: circontfra(klon, klev)    ! contrail induced cirrus fraction [-]
     123  REAL, INTENT(IN)  :: qice_lincont(klon, klev)  ! linear contrails condensed water [kg/kg]
     124  REAL, INTENT(IN)  :: qice_circont(klon, klev) ! contrail induced cirrus condensed water [kg/kg]
    125125  REAL, INTENT(OUT) :: pch_nocont(klon)          ! 2D high cloud cover without contrails[-]
    126126  REAL, INTENT(OUT) :: pct_cont(klon)            ! 2D total contrails cover[-]
    … …  
    173173  REAL zflwp_var, zfiwp_var
    174174  REAL d_rei_dt
    175   REAL pclc_cont(klon,klev), pclc_perscont(klon,klev)
    176   REAL rei_cont, rei_perscont
     175  REAL pclc_lincont(klon,klev), pclc_circont(klon,klev)
     176  REAL rei_lincont, rei_circont
    177177
    178178
    … …  
    255255      DO k = 1, klev
    256256        DO i = 1, klon
    257           pclc_nocont(i,k) = pclc(i, k) - contfra(i, k) - perscontfra(i, k)
    258           xfiwc_nocont(i, k) = xfiwc(i, k) - qice_cont(i, k) - qice_perscont(i, k)
     257          pclc_nocont(i,k) = pclc(i, k) - lincontfra(i, k) - circontfra(i, k)
     258          xfiwc_nocont(i, k) = xfiwc(i, k) - qice_lincont(i, k) - qice_circont(i, k)
    259259        ENDDO
    260260      ENDDO
    … …  
    574574        reice_nocont(i,k) = 0.
    575575        pclc_nocont(i,k) = 0.
    576         pclc_cont(i,k) = 0.
    577         pclc_perscont(i,k) = 0.
     576        pclc_lincont(i,k) = 0.
     577        pclc_circont(i,k) = 0.
    578578        pcltau_cont(i,k) = 0.
    579579        pclemi_cont(i,k) = 0.
    … …  
    583583      ELSE
    584584
    585         IF ( pclc_nocont(i,k) .GT. 0.01 * seuil_neb ) THEN
     585        IF ( pclc_nocont(i,k) .GT. zepsec ) THEN
    586586          ! -- liquid/ice cloud water paths:
    587587
    … …  
    682682          reice_nocont(i,k) = 1.
    683683          pclc_nocont(i,k) = 0.
    684           pclc(i,k) = contfra(i,k) + perscontfra(i,k)
     684          pclc(i,k) = lincontfra(i,k) + circontfra(i,k)
    685685          pcltau_nocont(i, k) = 0.
    686686          pclemi_nocont(i, k) = 0.
    … …  
    689689
    690690
    691         IF ( contfra(i,k) .GT. 0.01 * seuil_neb ) THEN
    692           pclc_cont(i,k) = contfra(i,k)
    693           rei_cont = re_ice_crystals_contrails * 1.E6
     691        IF ( lincontfra(i,k) .GT. zepsec ) THEN
     692          pclc_lincont(i,k) = lincontfra(i,k)
     693          rei_lincont = re_ice_crystals_contrails * 1.E6
    694694        ELSE
    695           pclc_cont(i,k) = 0.
    696           rei_cont = 1.
    697         ENDIF
    698 
    699 
    700         IF ( perscontfra(i,k) .GT. 0.01 * seuil_neb ) THEN
     695          pclc_lincont(i,k) = 0.
     696          rei_lincont = 1.
     697        ENDIF
     698
     699
     700        IF ( circontfra(i,k) .GT. zepsec ) THEN
    701701          !--Everything is the same but with contrails
    702           pclc_perscont(i,k) = perscontfra(i,k)
     702          pclc_circont(i,k) = circontfra(i,k)
    703703
    704704          tc = temp(i, k) - 273.15
    705           rei_perscont = d_rei_dt*tc + rei_max
    706           IF (tc<=-81.4) rei_perscont = rei_min
     705          rei_circont = d_rei_dt*tc + rei_max
     706          IF (tc<=-81.4) rei_circont = rei_min
    707707
    708708        ELSE
    709           pclc_perscont(i,k) = 0.
    710           rei_perscont = 1.
    711         ENDIF
    712 
    713         IF ( MAX(contfra(i,k), perscontfra(i,k)) .GT. 0.01 * seuil_neb ) THEN
    714 
    715           rei = ( rei_cont * pclc_cont(i,k) + rei_perscont * pclc_perscont(i,k) ) &
    716               / ( pclc_cont(i,k) + pclc_perscont(i,k) )
     709          pclc_circont(i,k) = 0.
     710          rei_circont = 1.
     711        ENDIF
     712
     713        IF ( MAX(lincontfra(i,k), circontfra(i,k)) .GT. zepsec ) THEN
     714
     715          rei = ( rei_lincont * pclc_lincont(i,k) + rei_circont * pclc_circont(i,k) ) &
     716              / ( pclc_lincont(i,k) + pclc_circont(i,k) )
    717717          zfiwp_var = 1000.*(xfiwc(i, k)-xfiwc_nocont(i, k))&
    718718                    / (pclc(i, k)-pclc_nocont(i, k))*rhodz(i, k)
    … …  
    732732
    733733
    734         rei = ( rei_cont * pclc_cont(i,k) + rei_perscont * pclc_perscont(i,k) &
     734        rei = ( rei_lincont * pclc_lincont(i,k) + rei_circont * pclc_circont(i,k) &
    735735            + reice_nocont(i, k) * pclc_nocont(i, k) ) &
    736             / ( pclc_cont(i,k) + pclc_perscont(i,k) + pclc_nocont(i,k) )
     736            / ( pclc_lincont(i,k) + pclc_circont(i,k) + pclc_nocont(i,k) )
    737737
    738738        zflwp_var = 1000.*xflwc(i, k)/pclc(i, k)*rhodz(i, k)
    … …  
    881881    DO k = klev, 1, -1
    882882      DO i = 1, klon
    883         zclear(i) = zclear(i)*(1.-max(pclc_cont(i,k)+pclc_perscont(i,k),zcloud(i)))/&
     883        zclear(i) = zclear(i)*(1.-max(pclc_lincont(i,k)+pclc_circont(i,k),zcloud(i)))/&
    884884          (1.-min(real(zcloud(i),kind=8),1.-zepsec))
    885885        pct_cont(i) = 1. - zclear(i)
    886         zcloud(i) = pclc_cont(i,k) + pclc_perscont(i,k)
     886        zcloud(i) = pclc_lincont(i,k) + pclc_circont(i,k)
    887887        IF (paprs(i,k)<prmhc) THEN
    888888          pch_nocont(i) = pch_nocont(i)*(1.-max(pclc_nocont(i,k),zcloudh(i)))/(1.-min(real( &

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

    r5631 r5641  
    2525     dqi_adj, dqi_sub, dqi_con, dqi_mix, dqvc_adj,      &
    2626     dqvc_sub, dqvc_con, dqvc_mix, qsatl, qsati,        &
    27      cfa_seri, pcf_seri, qva_seri, qia_seri,flight_dist,&
    28      flight_h2o, qice_perscont, Tcritcont,              &
     27     cfl_seri, cfc_seri, qtl_seri, qtc_seri,flight_dist,&
     28     flight_h2o, qice_lincont, qice_circont, Tcritcont, &
    2929     qcritcont, potcontfraP, potcontfraNP,              &
    3030     cloudth_sth,                                       &
    … …  
    128128USE phys_local_var_mod, ONLY : issrfra100to150, issrfra150to200, issrfra200to250
    129129USE phys_local_var_mod, ONLY : issrfra250to300, issrfra300to400, issrfra400to500
    130 USE phys_local_var_mod, ONLY : dcfa_ini, dqia_ini, dqta_ini, dcfa_sub, dqia_sub, dqta_sub
    131 USE phys_local_var_mod, ONLY : dcfa_cir, dqta_cir, dcfa_mix, dqia_mix, dqta_mix
     130USE phys_local_var_mod, ONLY : dcfl_ini, dqil_ini, dqtl_ini, dcfl_sub, dqil_sub, dqtl_sub
     131USE phys_local_var_mod, ONLY : dcfl_cir, dqtl_cir, dcfl_mix, dqil_mix, dqtl_mix
     132USE phys_local_var_mod, ONLY : dcfc_sub, dqic_sub, dqtc_sub, dcfc_mix, dqic_mix, dqtc_mix
    132133
    133134IMPLICIT NONE
    … …  
    190191  ! INPUT/OUTPUT aviation
    191192  !--------------------------------------------------
    192   REAL, DIMENSION(klon,klev),      INTENT(INOUT):: cfa_seri         ! linear contrails fraction [-]
    193   REAL, DIMENSION(klon,klev),      INTENT(INOUT):: pcf_seri         ! contrails induced cirrus fraction [-]
    194   REAL, DIMENSION(klon,klev),      INTENT(INOUT):: qva_seri         ! linear contrails total specific humidity [kg/kg]
    195   REAL, DIMENSION(klon,klev),      INTENT(INOUT):: qia_seri         ! linear contrails total specific humidity [kg/kg]
     193  REAL, DIMENSION(klon,klev),      INTENT(INOUT):: cfl_seri         ! linear contrails fraction [-]
     194  REAL, DIMENSION(klon,klev),      INTENT(INOUT):: cfc_seri         ! contrail cirrus fraction [-]
     195  REAL, DIMENSION(klon,klev),      INTENT(INOUT):: qtl_seri         ! linear contrails total specific humidity [kg/kg]
     196  REAL, DIMENSION(klon,klev),      INTENT(INOUT):: qtc_seri         ! contrail cirrus total specific humidity [kg/kg]
    196197  REAL, DIMENSION(klon,klev),      INTENT(IN)   :: flight_dist      ! aviation distance flown within the mesh [m/s/mesh]
    197198  REAL, DIMENSION(klon,klev),      INTENT(IN)   :: flight_h2o       ! aviation H2O emitted within the mesh [kgH2O/s/mesh]
    … …  
    251252  ! for contrails and aviation
    252253
    253   REAL, DIMENSION(klon,klev),      INTENT(OUT)  :: qice_perscont  !--condensed water in contrail induced cirrus used in the radiation scheme [kg/kg]
     254  REAL, DIMENSION(klon,klev),      INTENT(OUT)  :: qice_lincont   !--condensed water in linear contrails used in the radiation scheme [kg/kg]
     255  REAL, DIMENSION(klon,klev),      INTENT(OUT)  :: qice_circont   !--condensed water in contrail cirrus used in the radiation scheme [kg/kg]
    254256  REAL, DIMENSION(klon,klev),      INTENT(OUT)  :: Tcritcont      !--critical temperature for contrail formation [K]
    255257  REAL, DIMENSION(klon,klev),      INTENT(OUT)  :: qcritcont      !--critical specific humidity for contrail formation [kg/kg]
    … …  
    332334  REAL :: delta_z
    333335  ! for contrails
    334   REAL, DIMENSION(klon) :: contfra, perscontfra, qcont
     336  REAL, DIMENSION(klon) :: lincontfra, circontfra, qlincont, qcircont
    335337  REAL, DIMENSION(klon) :: zq_nodeep
    336338  LOGICAL, DIMENSION(klon) :: pt_pron_clds
    … …  
    662664          !--Initialisation
    663665          IF ( ok_plane_contrail ) THEN
    664             contfra(:)     = 0.
    665             qcont(:)       = 0.
    666             perscontfra(:) = 0.
     666            lincontfra(:) = 0.
     667            circontfra(:) = 0.
     668            qlincont(:)   = 0.
     669            qcircont(:)   = 0.
    667670          ENDIF
    668671
    … …  
    816819                        dqi_adj(:,k), dqi_sub(:,k), dqi_con(:,k), dqi_mix(:,k), dqised(:,k), &
    817820                        dqvc_adj(:,k), dqvc_sub(:,k), dqvc_con(:,k), dqvc_mix(:,k), dqvcsed(:,k), &
    818                         cfa_seri(:,k), pcf_seri(:,k), qva_seri(:,k), qia_seri(:,k), &
    819                         flight_dist(:,k), flight_h2o(:,k), contfra, perscontfra, qcont, &
     821                        cfl_seri(:,k), cfc_seri(:,k), qtl_seri(:,k), qtc_seri(:,k), &
     822                        flight_dist(:,k), flight_h2o(:,k), &
     823                        lincontfra, circontfra, qlincont, qcircont, &
    820824                        Tcritcont(:,k), qcritcont(:,k), potcontfraP(:,k), potcontfraNP(:,k), &
    821                         dcfa_ini(:,k), dqia_ini(:,k), dqta_ini(:,k), &
    822                         dcfa_sub(:,k), dqia_sub(:,k), dqta_sub(:,k), &
    823                         dcfa_cir(:,k), dqta_cir(:,k), &
    824                         dcfa_mix(:,k), dqia_mix(:,k), dqta_mix(:,k))
     825                        dcfl_ini(:,k), dqil_ini(:,k), dqtl_ini(:,k), &
     826                        dcfl_sub(:,k), dqil_sub(:,k), dqtl_sub(:,k), &
     827                        dcfl_cir(:,k), dqtl_cir(:,k), &
     828                        dcfl_mix(:,k), dqil_mix(:,k), dqtl_mix(:,k), &
     829                        dcfc_sub(:,k), dqic_sub(:,k), dqtc_sub(:,k), &
     830                        dcfc_mix(:,k), dqic_mix(:,k), dqtc_mix(:,k))
    825831
    826832                    DO i = 1, klon
    … …  
    10491055      !--Contrails precipitate as natural clouds. We save the partition of ice
    10501056      !--between natural clouds and contrails
    1051       !--NB. we use qcont as a temporary variable to save this partition
     1057      !--NB. we use qlincont / qcircont as a temporary variable to save this partition
    10521058      DO i = 1, klon
     1059        dcf_sub(i,k) = lincontfra(i)
     1060        dqi_sub(i,k) = qlincont(i)
    10531061        IF ( zoliqi(i) .GT. 0. ) THEN
    1054           qcont(i) = ( qcont(i) - zqs(i) * contfra(i) ) / zoliqi(i)
     1062          qlincont(i) = ( qlincont(i) - zqs(i) * lincontfra(i) ) / zoliqi(i)
     1063          qcircont(i) = ( qcircont(i) - zqs(i) * circontfra(i) ) / zoliqi(i)
    10551064        ELSE
    1056           qcont(i) = 0.
     1065          qlincont(i) = 0.
     1066          qcircont(i) = 0.
    10571067        ENDIF
    10581068      ENDDO
    … …  
    10951105      !--Contrails fraction is left unchanged, but contrails water has changed
    10961106      DO i = 1, klon
    1097         IF ( zoliqi(i) .LE. 0. ) THEN
    1098           contfra(i) = 0.
    1099           qcont(i) = 0.
     1107        IF ( zoliqi(i) .GT. 0. ) THEN
     1108          qlincont(i) = zqs(i) * lincontfra(i) + zoliqi(i) * qlincont(i)
     1109          qcircont(i) = zqs(i) * circontfra(i) + zoliqi(i) * qcircont(i)
    11001110        ELSE
    1101           qcont(i) = zqs(i) * contfra(i) + zoliqi(i) * qcont(i)
     1111          lincontfra(i) = 0.
     1112          circontfra(i) = 0.
     1113          qlincont(i) = 0.
     1114          qcircont(i) = 0.
    11021115        ENDIF
    11031116      ENDDO
    … …  
    11981211    !--AB Write diagnostics and tracers for ice supersaturation
    11991212    IF ( ok_plane_contrail ) THEN
    1200       cfa_seri(:,k) = contfra(:)
    1201       pcf_seri(:,k) = perscontfra(:)
    1202       qva_seri(:,k) = zqs(:) * contfra(:)
     1213      cfl_seri(:,k) = lincontfra(:)
     1214      cfc_seri(:,k) = circontfra(:)
     1215      qtl_seri(:,k) = qlincont(:)
     1216      qtc_seri(:,k) = qcircont(:)
    12031217      DO i = 1, klon
    12041218        IF ( zoliqi(i) .GT. 0. ) THEN
    1205           !--The quantity of ice in linear contrails has been reduced by precipitation
    1206           !--with the same factor as the rest of the clouds
    1207           qia_seri(i,k) = ( qcont(i) - zqs(i) * contfra(i) ) / zoliqi(i) * radocond(i,k)
     1219          !--The ice water content seen by radiation is higher than the actual ice
     1220          !--water content. We take into account this difference
     1221          qice_lincont(i,k) = ( qlincont(i) - zqs(i) * lincontfra(i) ) &
     1222              / zoliqi(i) * radocond(i,k)
     1223          qice_circont(i,k) = ( qcircont(i) - zqs(i) * circontfra(i) ) &
     1224              / zoliqi(i) * radocond(i,k)
    12081225        ELSE
    1209           qia_seri(i,k) = 0.
    1210         ENDIF
    1211         IF ( ( rneb(i,k) - cfa_seri(i,k) ) .GT. eps ) THEN
    1212           !--The in-cloud quantity of ice in persistent contrail cirrus is the same as
    1213           !--the one from all natural cirrus clouds
    1214           qice_perscont(i,k) = ( radocond(i,k) - qia_seri(i,k) ) &
    1215               * perscontfra(i) / ( rneb(i,k) - cfa_seri(i,k) )
    1216         ELSE
    1217           qice_perscont(i,k) = 0.
     1226          qice_lincont(i,k) = 0.
     1227          qice_circont(i,k) = 0.
    12181228        ENDIF
    12191229      ENDDO

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

    r5639 r5641  
    101101      dqi_adj, dqi_sub, dqi_con, dqi_mix, dqi_sed, &
    102102      dqvc_adj, dqvc_sub, dqvc_con, dqvc_mix, dqvc_sed, &
    103       contfra_in, perscontfra_in, qva_in, qia_in, flight_dist, flight_h2o, &
    104       contfra, perscontfra, qcont, Tcritcont, qcritcont, potcontfraP, potcontfraNP, &
    105       dcfa_ini, dqia_ini, dqta_ini, dcfa_sub, dqia_sub, dqta_sub, &
    106       dcfa_cir, dqta_cir, dcfa_mix, dqia_mix, dqta_mix)
     103      lincontfra_in, circontfra_in, qtl_in, qtc_in, flight_dist, flight_h2o, &
     104      lincontfra, circontfra, qlincont, qcircont, &
     105      Tcritcont, qcritcont, potcontfraP, potcontfraNP, &
     106      dcfl_ini, dqil_ini, dqtl_ini, dcfl_sub, dqil_sub, dqtl_sub, &
     107      dcfl_cir, dqtl_cir, dcfl_mix, dqil_mix, dqtl_mix, &
     108      dcfc_sub, dqic_sub, dqtc_sub, dcfc_mix, dqic_mix, dqtc_mix)
    107109
    108110!----------------------------------------------------------------------
    … …  
    127129USE lmdz_lscp_ini,   ONLY: beta_pdf_lscp, temp_thresh_pdf_lscp
    128130USE lmdz_lscp_ini,   ONLY: std100_pdf_lscp, k0_pdf_lscp, kappa_pdf_lscp
    129 USE lmdz_lscp_ini,   ONLY: coef_mixing_lscp, coef_shear_lscp, cooling_rate_ice_thresh
    130 
    131 USE lmdz_lscp_ini,   ONLY: ok_plane_contrail, aspect_ratio_contrails
    132 USE lmdz_lscp_ini,   ONLY: coef_mixing_contrails, coef_shear_contrails
    133 USE lmdz_lscp_ini,   ONLY: chi_mixing_contrails, linear_contrails_lifetime
     131USE lmdz_lscp_ini,   ONLY: coef_mixing_lscp, coef_shear_lscp
     132USE lmdz_lscp_ini,   ONLY: aspect_ratio_cirrus, cooling_rate_ice_thresh
     133
     134USE lmdz_lscp_ini,   ONLY: ok_plane_contrail, aspect_ratio_lincontrails
     135USE lmdz_lscp_ini,   ONLY: coef_mixing_lincontrails, coef_shear_lincontrails
     136USE lmdz_lscp_ini,   ONLY: chi_mixing_lincontrails, linear_contrails_lifetime
    134137USE lmdz_aviation,   ONLY: contrails_formation
    135138
    … …  
    165168! Input for aviation
    166169!
    167 REAL,     INTENT(IN)   , DIMENSION(klon) :: contfra_in    ! input linear contrails fraction [-]
    168 REAL,     INTENT(IN)   , DIMENSION(klon) :: perscontfra_in! input contrail induced cirrus fraction [-]
    169 REAL,     INTENT(IN)   , DIMENSION(klon) :: qva_in        ! input linear contrails total specific humidity [kg/kg]
    170 REAL,     INTENT(IN)   , DIMENSION(klon) :: qia_in        ! input linear contrails total specific humidity [kg/kg]
     170REAL,     INTENT(IN)   , DIMENSION(klon) :: lincontfra_in ! input linear contrails fraction [-]
     171REAL,     INTENT(IN)   , DIMENSION(klon) :: circontfra_in ! input contrail cirrus fraction [-]
     172REAL,     INTENT(IN)   , DIMENSION(klon) :: qtl_in        ! input linear contrails total specific humidity [kg/kg]
     173REAL,     INTENT(IN)   , DIMENSION(klon) :: qtc_in        ! input contrail cirrus total specific humidity [kg/kg]
    171174REAL,     INTENT(IN)   , DIMENSION(klon) :: flight_dist   ! aviation distance flown concentration [m/s/m3]
    172175REAL,     INTENT(IN)   , DIMENSION(klon) :: flight_h2o    ! aviation emitted H2O concentration [kgH2O/s/m3]
    … …  
    206209!  NB. idem for the INOUT
    207210!
    208 REAL,     INTENT(INOUT), DIMENSION(klon) :: contfra      ! linear contrail fraction [-]
    209 REAL,     INTENT(INOUT), DIMENSION(klon) :: perscontfra  ! linear contrail induced cirrus fraction [-]
    210 REAL,     INTENT(INOUT), DIMENSION(klon) :: qcont        ! linear contrail specific humidity [kg/kg]
     211REAL,     INTENT(INOUT), DIMENSION(klon) :: lincontfra   ! linear contrail fraction [-]
     212REAL,     INTENT(INOUT), DIMENSION(klon) :: circontfra   ! contrail cirrus fraction [-]
     213REAL,     INTENT(INOUT), DIMENSION(klon) :: qlincont     ! linear contrail specific humidity [kg/kg]
     214REAL,     INTENT(INOUT), DIMENSION(klon) :: qcircont     ! contrail cirrus specific humidity [kg/kg]
    211215REAL,     INTENT(INOUT), DIMENSION(klon) :: Tcritcont    ! critical temperature for contrail formation [K]
    212216REAL,     INTENT(INOUT), DIMENSION(klon) :: qcritcont    ! critical specific humidity for contrail formation [kg/kg]
    213217REAL,     INTENT(INOUT), DIMENSION(klon) :: potcontfraP  ! potential persistent contrail fraction [-]
    214218REAL,     INTENT(INOUT), DIMENSION(klon) :: potcontfraNP ! potential non-persistent contrail fraction [-]
    215 REAL,     INTENT(INOUT), DIMENSION(klon) :: dcfa_ini     ! contrails cloud fraction tendency because of initial formation [s-1]
    216 REAL,     INTENT(INOUT), DIMENSION(klon) :: dqia_ini     ! contrails ice specific humidity tendency because of initial formation [kg/kg/s]
    217 REAL,     INTENT(INOUT), DIMENSION(klon) :: dqta_ini     ! contrails total specific humidity tendency because of initial formation [kg/kg/s]
    218 REAL,     INTENT(INOUT), DIMENSION(klon) :: dcfa_sub     ! contrails cloud fraction tendency because of sublimation [s-1]
    219 REAL,     INTENT(INOUT), DIMENSION(klon) :: dqia_sub     ! contrails ice specific humidity tendency because of sublimation [kg/kg/s]
    220 REAL,     INTENT(INOUT), DIMENSION(klon) :: dqta_sub     ! contrails total specific humidity tendency because of sublimation [kg/kg/s]
    221 REAL,     INTENT(INOUT), DIMENSION(klon) :: dcfa_cir     ! contrails cloud fraction tendency because of conversion in cirrus [s-1]
    222 REAL,     INTENT(INOUT), DIMENSION(klon) :: dqta_cir     ! contrails total specific humidity tendency because of conversion in cirrus [kg/kg/s]
    223 REAL,     INTENT(INOUT), DIMENSION(klon) :: dcfa_mix     ! contrails cloud fraction tendency because of mixing [s-1]
    224 REAL,     INTENT(INOUT), DIMENSION(klon) :: dqia_mix     ! contrails ice specific humidity tendency because of mixing [kg/kg/s]
    225 REAL,     INTENT(INOUT), DIMENSION(klon) :: dqta_mix     ! contrails total specific humidity tendency because of mixing [kg/kg/s]
     219REAL,     INTENT(INOUT), DIMENSION(klon) :: dcfl_ini     ! linear contrails cloud fraction tendency because of initial formation [s-1]
     220REAL,     INTENT(INOUT), DIMENSION(klon) :: dqil_ini     ! linear contrails ice specific humidity tendency because of initial formation [kg/kg/s]
     221REAL,     INTENT(INOUT), DIMENSION(klon) :: dqtl_ini     ! linear contrails total specific humidity tendency because of initial formation [kg/kg/s]
     222REAL,     INTENT(INOUT), DIMENSION(klon) :: dcfl_sub     ! linear contrails cloud fraction tendency because of sublimation [s-1]
     223REAL,     INTENT(INOUT), DIMENSION(klon) :: dqil_sub     ! linear contrails ice specific humidity tendency because of sublimation [kg/kg/s]
     224REAL,     INTENT(INOUT), DIMENSION(klon) :: dqtl_sub     ! linear contrails total specific humidity tendency because of sublimation [kg/kg/s]
     225REAL,     INTENT(INOUT), DIMENSION(klon) :: dcfl_cir     ! linear contrails cloud fraction tendency because of conversion in cirrus [s-1]
     226REAL,     INTENT(INOUT), DIMENSION(klon) :: dqtl_cir     ! linear contrails total specific humidity tendency because of conversion in cirrus [kg/kg/s]
     227REAL,     INTENT(INOUT), DIMENSION(klon) :: dcfl_mix     ! linear contrails cloud fraction tendency because of mixing [s-1]
     228REAL,     INTENT(INOUT), DIMENSION(klon) :: dqil_mix     ! linear contrails ice specific humidity tendency because of mixing [kg/kg/s]
     229REAL,     INTENT(INOUT), DIMENSION(klon) :: dqtl_mix     ! linear contrails total specific humidity tendency because of mixing [kg/kg/s]
     230REAL,     INTENT(INOUT), DIMENSION(klon) :: dcfc_sub     ! contrail cirrus cloud fraction tendency because of sublimation [s-1]
     231REAL,     INTENT(INOUT), DIMENSION(klon) :: dqic_sub     ! contrail cirrus ice specific humidity tendency because of sublimation [kg/kg/s]
     232REAL,     INTENT(INOUT), DIMENSION(klon) :: dqtc_sub     ! contrail cirrus total specific humidity tendency because of sublimation [kg/kg/s]
     233REAL,     INTENT(INOUT), DIMENSION(klon) :: dcfc_mix     ! contrail cirrus cloud fraction tendency because of mixing [s-1]
     234REAL,     INTENT(INOUT), DIMENSION(klon) :: dqic_mix     ! contrail cirrus ice specific humidity tendency because of mixing [kg/kg/s]
     235REAL,     INTENT(INOUT), DIMENSION(klon) :: dqtc_mix     ! contrail cirrus total specific humidity tendency because of mixing [kg/kg/s]
    226236!
    227237! Local
    … …  
    273283!
    274284! for contrails
    275 REAL :: perscontfra_ratio
    276285REAL :: contrails_conversion_factor
    277286
    … …  
    455464      !--If contrails are activated
    456465      IF ( ok_plane_contrail ) THEN
    457         contfra(i) = MAX(0., MIN(cldfra(i), contfra_in(i)))
    458         perscontfra(i) = MAX(0., MIN(cldfra(i) - contfra(i), perscontfra_in(i)))
    459         qcont(i) = MAX(0., MIN(qcld(i), qva_in(i) + qia_in(i)))
    460 
    461         dcfa_ini(i) = 0.
    462         dqia_ini(i) = 0.
    463         dqta_ini(i) = 0.
    464         dcfa_sub(i) = 0.
    465         dqia_sub(i) = 0.
    466         dqta_sub(i) = 0.
    467         dcfa_cir(i) = 0.
    468         dqta_cir(i) = 0.
    469         dcfa_mix(i) = 0.
    470         dqia_mix(i) = 0.
    471         dqta_mix(i) = 0.
     466        lincontfra(i) = MAX(0., lincontfra_in(i))
     467        circontfra(i) = MAX(0., circontfra_in(i))
     468        qlincont(i) = MAX(0., qtl_in(i))
     469        qcircont(i) = MAX(0., qtc_in(i))
     470        !--The following barriers are needed since the advection scheme does not
     471        !--conserve order relations
     472        mixed_fraction = lincontfra(i) + circontfra(i)
     473        IF ( mixed_fraction .GT. cldfra(i) ) THEN
     474          mixed_fraction = cldfra(i) / mixed_fraction
     475          lincontfra(i) = lincontfra(i) * mixed_fraction
     476          circontfra(i) = circontfra(i) * mixed_fraction
     477          qlincont(i)   = qlincont(i)   * mixed_fraction
     478          qcircont(i)   = qcircont(i)   * mixed_fraction
     479        ENDIF
     480        mixed_fraction = qlincont(i) + qcircont(i)
     481        IF ( mixed_fraction .GT. qcld(i) ) THEN
     482          mixed_fraction = qcld(i) / mixed_fraction
     483          lincontfra(i) = lincontfra(i) * mixed_fraction
     484          circontfra(i) = circontfra(i) * mixed_fraction
     485          qlincont(i)   = qlincont(i)   * mixed_fraction
     486          qcircont(i)   = qcircont(i)   * mixed_fraction
     487        ENDIF
     488
     489
     490        dcfl_ini(i) = 0.
     491        dqil_ini(i) = 0.
     492        dqtl_ini(i) = 0.
     493        dcfl_sub(i) = 0.
     494        dqil_sub(i) = 0.
     495        dqtl_sub(i) = 0.
     496        dcfl_cir(i) = 0.
     497        dqtl_cir(i) = 0.
     498        dcfl_mix(i) = 0.
     499        dqil_mix(i) = 0.
     500        dqtl_mix(i) = 0.
     501        dcfc_sub(i) = 0.
     502        dqic_sub(i) = 0.
     503        dqtc_sub(i) = 0.
     504        dcfc_mix(i) = 0.
     505        dqic_mix(i) = 0.
     506        dqtc_mix(i) = 0.
    472507      ELSE
    473         contfra(i) = 0.
    474         perscontfra(i) = 0.
    475         qcont(i) = 0.
     508        lincontfra(i) = 0.
     509        circontfra(i) = 0.
     510        qlincont(i) = 0.
     511        qcircont(i) = 0.
    476512      ENDIF
    477513
    … …  
    481517      !----------------------------------------------------------------------
    482518
    483       !--If there is a contrail
    484       IF ( contfra(i) .GT. eps ) THEN
     519      !--If there is a linear contrail
     520      IF ( lincontfra(i) .GT. eps ) THEN
    485521        !--We remove contrails from the main class
    486         cldfra(i) = cldfra(i) - contfra(i)
    487         qcld(i) = qcld(i) - qcont(i)
    488         qvc(i) = qvc(i) - qsat(i) * contfra(i)
     522        cldfra(i) = cldfra(i) - lincontfra(i)
     523        qcld(i) = qcld(i) - qlincont(i)
     524        qvc(i) = qvc(i) - qsat(i) * lincontfra(i)
    489525
    490526        !--The contrail is always adjusted to saturation
    491         qiceincld = ( qcont(i) / contfra(i) - qsat(i) )
     527        qiceincld = ( qlincont(i) / lincontfra(i) - qsat(i) )
    492528       
    493529        !--If the ice water content is too low, the cloud is purely sublimated
    494530        IF ( qiceincld .LT. eps ) THEN
    495           dcfa_sub(i) = - contfra(i)
    496           dqia_sub(i) = - qiceincld * contfra(i)
    497           dqta_sub(i) = - qcont(i)
    498           contfra(i) = 0.
    499           qcont(i) = 0.
    500           clrfra(i) = MIN(1., clrfra(i) - dcfa_sub(i))
    501           qclr(i) = qclr(i) - dqta_sub(i)
     531          dcfl_sub(i) = - lincontfra(i)
     532          dqil_sub(i) = - qiceincld * lincontfra(i)
     533          dqtl_sub(i) = - qlincont(i)
     534          lincontfra(i) = 0.
     535          qlincont(i) = 0.
     536          clrfra(i) = MIN(1., clrfra(i) - dcfl_sub(i))
     537          qclr(i) = qclr(i) - dqtl_sub(i)
    502538        ENDIF   ! qiceincld .LT. eps
    503       ENDIF ! contfra(i) .GT. eps
    504 
    505       !--If there is a contrail induced cirrus, we save the ratio
    506       perscontfra_ratio = perscontfra(i) / MAX(eps, cldfra(i))
     539      ENDIF ! lincontfra(i) .GT. eps
     540
     541      !--If there is a contrail cirrus
     542      IF ( circontfra(i) .GT. eps ) THEN
     543        !--We remove contrails from the main class
     544        cldfra(i) = cldfra(i) - circontfra(i)
     545        qcld(i) = qcld(i) - qcircont(i)
     546        qvc(i) = qvc(i) - qsat(i) * circontfra(i)
     547
     548        !--The contrail is always adjusted to saturation
     549        qiceincld = ( qcircont(i) / circontfra(i) - qsat(i) )
     550       
     551        !--If the ice water content is too low, the cloud is purely sublimated
     552        IF ( qiceincld .LT. eps ) THEN
     553          dcfc_sub(i) = - circontfra(i)
     554          dqic_sub(i) = - qiceincld * circontfra(i)
     555          dqtc_sub(i) = - qcircont(i)
     556          circontfra(i) = 0.
     557          qcircont(i) = 0.
     558          clrfra(i) = MIN(1., clrfra(i) - dcfc_sub(i))
     559          qclr(i) = qclr(i) - dqtc_sub(i)
     560        ENDIF   ! qiceincld .LT. eps
     561      ENDIF ! circontfra(i) .GT. eps
    507562
    508563
    … …  
    640695      ENDIF ! cldfra(i) .GT. eps
    641696
    642       !--If there is a contrail induced cirrus, we restore it
    643       perscontfra(i) = perscontfra_ratio * cldfra(i)
    644 
    645697
    646698      !--------------------------------------------------------------------------
    … …  
    750802      ENDIF ! clrfra(i) .GT. eps
    751803
    752 
    753       !--If there is a contrail induced cirrus, we save the ratio
    754       perscontfra_ratio = perscontfra(i) / MAX(eps, cldfra(i))
    755804
    756805      !--------------------------------------
    … …  
    787836        !-- b / a = bovera = MAX(0.1, cldfra)
    788837        !bovera = MAX(0.1, cldfra(i))
    789         bovera = 0.111111111
     838        bovera = aspect_ratio_cirrus
    790839        !--P / a is a function of b / a only, that we can calculate
    791840        !-- P / a = RPI * ( 3. * ( 1. + b / a ) - SQRT( (3. + b / a) * (1. + 3. * b / a) ) )
    … …  
    925974      !--------------------------------------
    926975
    927       IF ( ( contfra(i) .GT. eps ) .AND. ( clrfra(i) .GT.  eps ) ) THEN
     976      IF ( ( lincontfra(i) .GT. eps ) .AND. ( clrfra(i) .GT.  eps ) ) THEN
    928977
    929978        !-- PART 1 - TURBULENT DIFFUSION
    … …  
    950999        !--are spherical.
    9511000        !-- b / a = bovera = MAX(0.1, cldfra)
    952         bovera = aspect_ratio_contrails
     1001        bovera = aspect_ratio_lincontrails
    9531002        !--P / a is a function of b / a only, that we can calculate
    9541003        !-- P / a = RPI * ( 3. * ( 1. + b / a ) - SQRT( (3. + b / a) * (1. + 3. * b / a) ) )
    … …  
    9611010        !-- cldfra = a * b / a * RPI / (P / a) * coef_mix_lscp * cldfra * ( 1. - cldfra )
    9621011        !-- a = (P / a) / ( coef_mix_lscp * RPI * ( 1. - cldfra ) * (b / a) )
    963         a_mix = Povera / coef_mixing_contrails / RPI / ( 1. - contfra(i) ) / bovera
     1012        a_mix = Povera / coef_mixing_lincontrails / RPI / ( 1. - lincontfra(i) ) / bovera
    9641013        !--and finally,
    9651014        !-- N_cld_mix = coef_mix_lscp * cldfra * ( 1. - cldfra ) * cell_area / ( P / a * a )
    966         N_cld_mix = coef_mixing_contrails * contfra(i) * ( 1. - contfra(i) ) * cell_area(i) &
    967                   / Povera / a_mix
     1015        N_cld_mix = coef_mixing_lincontrails * lincontfra(i) * ( 1. - lincontfra(i) ) &
     1016                  * cell_area(i) / Povera / a_mix
    9681017       
    9691018        !--The time required for turbulent diffusion to homogenize a region of size
    … …  
    9941043        !--With this, the clouds increase in size along b only, by a factor
    9951044        !--L_shear * SQRT(2.) / 2. (to account for the 45 degrees orientation of the wind)
    996         L_shear = coef_shear_contrails * shear(i) * dz * dtime
     1045        L_shear = coef_shear_lincontrails * shear(i) * dz * dtime
    9971046        !--therefore, the fraction of clear sky mixed is
    9981047        !-- N_cld_mix * ( a * (b + L_shear * SQRT(2.) / 2.) - a * b ) * RPI / 2. / cell_area
    … …  
    10211070        !--cloud's vapor without condensing or sublimating ice crystals
    10221071        qvapinclr_lim = qsat(i) * ( 1. + cldfra_mix / clrfra_mix ) &
    1023                       - qcont(i) / contfra(i) * cldfra_mix / clrfra_mix
     1072                      - qlincont(i) / lincontfra(i) * cldfra_mix / clrfra_mix
    10241073
    10251074        IF ( qvapinclr_lim .LT. 0. ) THEN
    10261075          !--Whatever we do, the cloud will increase in size
    1027           dcfa_mix(i) = clrfra_mix
    1028           dqta_mix(i) = clrfra_mix * qclr(i) / clrfra(i)
     1076          dcfl_mix(i) = clrfra_mix
     1077          dqtl_mix(i) = clrfra_mix * qclr(i) / clrfra(i)
    10291078        ELSE
    10301079          !--We then calculate the clear sky part where the humidity is lower than
    … …  
    10621111          !--decreases the size of the clouds
    10631112          !--For aviation, we increase the chance that the air surrounding contrails
    1064           !--is moist. This is quantified with chi_mixing_contrails
    1065           sigma_mix = chi_mixing_contrails * pdf_fra_above_lim &
    1066                     / ( pdf_fra_below_lim + chi_mixing_contrails * pdf_fra_above_lim )
     1113          !--is moist. This is quantified with chi_mixing_lincontrails
     1114          sigma_mix = chi_mixing_lincontrails * pdf_fra_above_lim &
     1115                    / ( pdf_fra_below_lim + chi_mixing_lincontrails * pdf_fra_above_lim )
    10671116
    10681117          IF ( pdf_fra_above_lim .GT. eps ) THEN
    1069             dcfa_mix(i) = clrfra_mix * sigma_mix
    1070             dqta_mix(i) = clrfra_mix * sigma_mix * pdf_q_above_lim / pdf_fra_above_lim
     1118            dcfl_mix(i) = clrfra_mix * sigma_mix
     1119            dqtl_mix(i) = clrfra_mix * sigma_mix * pdf_q_above_lim / pdf_fra_above_lim
    10711120          ENDIF
    10721121
    10731122          IF ( pdf_fra_below_lim .GT. eps ) THEN
    1074             qvapincld = qcont(i) / contfra(i)
     1123            qvapincld = qlincont(i) / lincontfra(i)
    10751124            qiceincld = qvapincld - qsat(i)
    1076             dcfa_mix(i) = dcfa_mix(i) - cldfra_mix * ( 1. - sigma_mix )
    1077             dqta_mix(i) = dqta_mix(i) - cldfra_mix * ( 1. - sigma_mix ) * qvapincld
    1078             dqia_mix(i) = dqia_mix(i) - cldfra_mix * ( 1. - sigma_mix ) * qiceincld
     1125            dcfl_mix(i) = dcfl_mix(i) - cldfra_mix * ( 1. - sigma_mix )
     1126            dqtl_mix(i) = dqtl_mix(i) - cldfra_mix * ( 1. - sigma_mix ) * qvapincld
     1127            dqil_mix(i) = dqil_mix(i) - cldfra_mix * ( 1. - sigma_mix ) * qiceincld
    10791128          ENDIF
    10801129
    10811130        ENDIF
    1082       ENDIF ! ( contfra(i) .GT. eps ) .AND. ( clrfra(i) .GT.  eps )
    1083 
    1084       IF ( contfra(i) .GT. eps ) THEN
    1085         !--We balance the mixing tendencies between the different cloud classes
    1086         mixed_fraction = dcf_mix(i) + dcfa_mix(i)
    1087         IF ( mixed_fraction .GT. clrfra(i) ) THEN
    1088           mixed_fraction = clrfra(i) / mixed_fraction
    1089           dcf_mix(i)  = dcf_mix(i)  * mixed_fraction
    1090           dqvc_mix(i) = dqvc_mix(i) * mixed_fraction
    1091           dqi_mix(i)  = dqi_mix(i)  * mixed_fraction
    1092           dcfa_mix(i) = dcfa_mix(i) * mixed_fraction
    1093           dqta_mix(i) = dqta_mix(i) * mixed_fraction
    1094         ENDIF
    1095 
    1096         !--Add tendencies
    1097         contfra(i) = contfra(i) + dcfa_mix(i)
    1098         qcont(i)   = qcont(i) + dqta_mix(i)
    1099         clrfra(i)  = clrfra(i) - dcfa_mix(i)
    1100         qclr(i)    = qclr(i) - dqta_mix(i)
    1101       ENDIF
    1102 
    1103       !--Add tendencies
    1104       cldfra(i)  = cldfra(i) + dcf_mix(i)
    1105       qcld(i)    = qcld(i) + dqvc_mix(i) + dqi_mix(i)
    1106       qvc(i)     = qvc(i) + dqvc_mix(i)
    1107       clrfra(i)  = clrfra(i) - dcf_mix(i)
    1108       qclr(i)    = qclr(i) - dqvc_mix(i) - dqi_mix(i)
    1109 
    1110       !--If there is a contrail induced cirrus, we restore it
    1111       perscontfra(i) = perscontfra_ratio * cldfra(i)
    1112 
    1113 
    1114       !---------------------------------------
    1115       !--         ICE SEDIMENTATION         --
    1116       !---------------------------------------
    1117       !
    1118       !--If ice supersaturation is activated, the cloud properties are prognostic.
    1119       !--The falling ice is then partly considered a new cloud in the gridbox.
    1120       !--BEWARE with this parameterisation, we can create a new cloud with a much
    1121       !--different ice water content and water vapor content than the existing cloud
    1122       !--(if it exists). This implies than unphysical fluxes of ice and vapor
    1123       !--occur between the existing cloud and the newly formed cloud (from sedimentation).
    1124       !--Note also that currently, the precipitation scheme assume a maximum
    1125       !--random overlap, meaning that the new formed clouds will not be affected
    1126       !--by vertical wind shear.
    1127       !
    1128       IF ( icesed_flux(i) .GT. 0. ) THEN
    1129 
    1130         clrfra_sed = MIN(clrfra(i), cldfra_above(i) - cldfra(i))
    1131 
    1132         IF ( ( clrfra_sed .GT. eps ) .AND. ( clrfra(i) .GT. eps ) ) THEN
    1133 
    1134           qiceincld = qice_sedim(i) / cldfra_above(i)
    1135           IF ( ok_unadjusted_clouds .AND. .NOT. ok_warm_cloud ) THEN
    1136             tauinv_depsub = qiceincld**(1./3.) * kappa_depsub
    1137             qvapinclr_lim = qsat(i) - qiceincld / ( 1. - EXP( - dtime * tauinv_depsub ) )
    1138           ELSE
    1139             qvapinclr_lim = qsat(i) - qiceincld
    1140           ENDIF
     1131      ENDIF ! ( lincontfra(i) .GT. eps ) .AND. ( clrfra(i) .GT.  eps )
     1132
     1133      IF ( ( circontfra(i) .GT. eps ) .AND. ( clrfra(i) .GT.  eps ) ) THEN
     1134
     1135        !-- PART 1 - TURBULENT DIFFUSION
     1136       
     1137        !--Clouds within the mesh are assumed to be ellipses. The length of the
     1138        !--semi-major axis is a and the length of the semi-minor axis is b.
     1139        !--N_cld_mix is the number of clouds in contact with clear sky, and can be non-integer.
     1140        !--In particular, it is 0 if cldfra = 1.
     1141        !--clouds_perim is the total perimeter of the clouds within the mesh,
     1142        !--not considering interfaces with other meshes (only the interfaces with clear
     1143        !--sky are taken into account).
     1144        !--
     1145        !--The area of each cloud is A = a * b * RPI,
     1146        !--and the perimeter of each cloud is
     1147        !-- P ~= RPI * ( 3 * (a + b) - SQRT( (3 * a + b) * (a + 3 * b) ) )
     1148        !--
     1149        !--With cell_area the area of the cell, we have:
     1150        !-- cldfra = A * N_cld_mix / cell_area
     1151        !-- clouds_perim = P * N_cld_mix
     1152        !--
     1153        !--We assume that the ratio between b and a is a function of
     1154        !--cldfra such that it is 1 for cldfra = 1 and it is low for little cldfra, because
     1155        !--if cldfra is low the clouds are linear, and if cldfra is high, the clouds
     1156        !--are spherical.
     1157        !-- b / a = bovera = MAX(0.1, cldfra)
     1158        bovera = aspect_ratio_cirrus
     1159        !--P / a is a function of b / a only, that we can calculate
     1160        !-- P / a = RPI * ( 3. * ( 1. + b / a ) - SQRT( (3. + b / a) * (1. + 3. * b / a) ) )
     1161        Povera = RPI * ( 3. * (1. + bovera) - SQRT( (3. + bovera) * (1. + 3. * bovera) ) )
     1162     
     1163        !--The clouds perimeter is imposed using the formula from Morcrette 2012,
     1164        !--based on observations.
     1165        !-- clouds_perim / cell_area = N_cld_mix * ( P / a * a ) / cell_area = coef_mix_lscp * cldfra * ( 1. - cldfra )
     1166        !--With cldfra = a * ( b / a * a ) * RPI * N_cld_mix / cell_area, we have:
     1167        !-- cldfra = a * b / a * RPI / (P / a) * coef_mix_lscp * cldfra * ( 1. - cldfra )
     1168        !-- a = (P / a) / ( coef_mix_lscp * RPI * ( 1. - cldfra ) * (b / a) )
     1169        a_mix = Povera / coef_mixing_lscp / RPI / ( 1. - circontfra(i) ) / bovera
     1170        !--and finally,
     1171        !-- N_cld_mix = coef_mix_lscp * cldfra * ( 1. - cldfra ) * cell_area / ( P / a * a )
     1172        N_cld_mix = coef_mixing_lscp * circontfra(i) * ( 1. - circontfra(i) ) &
     1173                  * cell_area(i) / Povera / a_mix
     1174       
     1175        !--The time required for turbulent diffusion to homogenize a region of size
     1176        !--L_mix is defined as (L_mix**2/tke_dissip)**(1./3.) (Pope, 2000; Field et al., 2014)
     1177        !--We compute L_mix and assume that the cloud is mixed over this length
     1178        L_mix = SQRT( dtime**3 * pbl_eps(i) )
     1179        !--The mixing length cannot be greater than the semi-minor axis. In this case,
     1180        !--the entire cloud is mixed.
     1181        L_mix = MIN(L_mix, a_mix * bovera)
     1182       
     1183        !--The fraction of clear sky mixed is
     1184        !-- N_cld_mix * ( (a + L_mix) * (b + L_mix) - a * b ) * RPI / cell_area
     1185        clrfra_mix = N_cld_mix * RPI / cell_area(i) &
     1186                   * ( a_mix * ( 1. + bovera ) * L_mix + L_mix**2 )
     1187        !--The fraction of clear sky mixed is
     1188        !-- N_cld_mix * ( a * b - (a - L_mix) * (b - L_mix) ) * RPI / cell_area
     1189        cldfra_mix = N_cld_mix * RPI / cell_area(i) &
     1190                   * ( a_mix * ( 1. + bovera ) * L_mix - L_mix**2 )
     1191       
     1192       
     1193        !-- PART 2 - SHEARING
     1194       
     1195        !--The clouds are then sheared. We keep the shape and number
     1196        !--assumptions from before. The clouds are sheared along their
     1197        !--semi-major axis (a_mix), on the entire cell heigh dz.
     1198        !--The increase in size is
     1199        L_shear = coef_shear_lscp * shear(i) * dz * dtime
     1200        !--therefore, the fraction of clear sky mixed is
     1201        !-- N_cld_mix * ( (a + L_shear) * b - a * b ) * RPI / 2. / cell_area
     1202        !--and the fraction of cloud mixed is
     1203        !-- N_cld_mix * ( (a * b) - (a - L_shear) * b ) * RPI / 2. / cell_area
     1204        !--(note that they are equal)
     1205        shear_fra = RPI * L_shear * a_mix * bovera / 2. * N_cld_mix / cell_area(i)
     1206        !--and the environment and cloud mixed fractions are the same,
     1207        !--which we add to the previous calculated mixed fractions.
     1208        !--We therefore assume that the sheared clouds and the turbulent
     1209        !--mixed clouds are different.
     1210        clrfra_mix = clrfra_mix + shear_fra
     1211        cldfra_mix = cldfra_mix + shear_fra
     1212       
     1213       
     1214        !-- PART 3 - CALCULATION OF THE MIXING PROPERTIES
     1215       
     1216        clrfra_mix = MAX(eps, MIN(clrfra(i), clrfra_mix))
     1217        cldfra_mix = MAX(eps, MIN(cldfra(i), cldfra_mix))
     1218       
     1219        !--We compute the limit vapor in clear sky where the mixed cloud could not
     1220        !--survive if all the ice crystals were sublimated. Note that here we assume,
     1221        !--for growth or reduction of the cloud, that the mixed cloud is adjusted
     1222        !--to saturation, ie the vapor in the mixed cloud is qsat. This is only a
     1223        !--diagnostic, and if the cloud size is increased, we add the new vapor to the
     1224        !--cloud's vapor without condensing or sublimating ice crystals
     1225        qvapinclr_lim = qsat(i) * ( 1. + cldfra_mix / clrfra_mix ) &
     1226                      - qcircont(i) / circontfra(i) * cldfra_mix / clrfra_mix
     1227
     1228        IF ( qvapinclr_lim .LT. 0. ) THEN
     1229          !--Whatever we do, the cloud will increase in size
     1230          dcfl_mix(i) = clrfra_mix
     1231          dqtl_mix(i) = clrfra_mix * qclr(i) / clrfra(i)
     1232        ELSE
     1233          !--We then calculate the clear sky part where the humidity is lower than
     1234          !--qvapinclr_lim, and the part where it is higher than qvapinclr_lim
     1235          !--This is the clear-sky PDF calculated in the condensation section. Note
     1236          !--that if we are here, we necessarily went through the condensation part
     1237          !--because the clear sky fraction can only be reduced by condensation.
     1238          !--Thus the `pdf_xxx` variables are well defined.
    11411239
    11421240          rhl_clr = qclr(i) / clrfra(i) / qsatl(i) * 100.
    … …  
    11601258          ENDIF
    11611259
     1260          pdf_fra_below_lim = clrfra(i) - pdf_fra_above_lim
     1261
     1262          !--sigma_mix is the ratio of the surroundings of the clouds where mixing
     1263          !--increases the size of the cloud, to the total surroundings of the clouds.
     1264          !--This implies that ( 1. - sigma_mix ) quantifies the ratio where mixing
     1265          !--decreases the size of the clouds
     1266          sigma_mix = pdf_fra_above_lim / ( pdf_fra_below_lim + pdf_fra_above_lim )
     1267
     1268          IF ( pdf_fra_above_lim .GT. eps ) THEN
     1269            dcfc_mix(i) = clrfra_mix * sigma_mix
     1270            dqtc_mix(i) = clrfra_mix * sigma_mix * pdf_q_above_lim / pdf_fra_above_lim
     1271          ENDIF
     1272
     1273          IF ( pdf_fra_below_lim .GT. eps ) THEN
     1274            qvapincld = qcircont(i) / circontfra(i)
     1275            qiceincld = qvapincld - qsat(i)
     1276            dcfc_mix(i) = dcfc_mix(i) - cldfra_mix * ( 1. - sigma_mix )
     1277            dqtc_mix(i) = dqtc_mix(i) - cldfra_mix * ( 1. - sigma_mix ) * qvapincld
     1278            dqic_mix(i) = dqic_mix(i) - cldfra_mix * ( 1. - sigma_mix ) * qiceincld
     1279          ENDIF
     1280
     1281        ENDIF
     1282      ENDIF ! ( circontfra(i) .GT. eps ) .AND. ( clrfra(i) .GT.  eps )
     1283
     1284      !--We balance the mixing tendencies between the different cloud classes
     1285      mixed_fraction = dcf_mix(i) + dcfl_mix(i) + dcfc_mix(i)
     1286      IF ( mixed_fraction .GT. clrfra(i) ) THEN
     1287        mixed_fraction = clrfra(i) / mixed_fraction
     1288        dcf_mix(i)  = dcf_mix(i)  * mixed_fraction
     1289        dqvc_mix(i) = dqvc_mix(i) * mixed_fraction
     1290        dqi_mix(i)  = dqi_mix(i)  * mixed_fraction
     1291        dcfl_mix(i) = dcfl_mix(i) * mixed_fraction
     1292        dqtl_mix(i) = dqtl_mix(i) * mixed_fraction
     1293        dqil_mix(i) = dqil_mix(i) * mixed_fraction
     1294        dcfc_mix(i) = dcfc_mix(i) * mixed_fraction
     1295        dqtc_mix(i) = dqtc_mix(i) * mixed_fraction
     1296        dqic_mix(i) = dqic_mix(i) * mixed_fraction
     1297      ENDIF
     1298
     1299      IF ( lincontfra(i) .GT. eps ) THEN
     1300        !--Add tendencies
     1301        lincontfra(i) = lincontfra(i) + dcfl_mix(i)
     1302        qlincont(i)   = qlincont(i) + dqtl_mix(i)
     1303        clrfra(i)  = clrfra(i) - dcfl_mix(i)
     1304        qclr(i)    = qclr(i) - dqtl_mix(i)
     1305      ENDIF
     1306
     1307      IF ( circontfra(i) .GT. eps ) THEN
     1308        !--Add tendencies
     1309        circontfra(i) = circontfra(i) + dcfc_mix(i)
     1310        qcircont(i)   = qcircont(i) + dqtc_mix(i)
     1311        clrfra(i)  = clrfra(i) - dcfc_mix(i)
     1312        qclr(i)    = qclr(i) - dqtc_mix(i)
     1313      ENDIF
     1314
     1315      !--Add tendencies
     1316      cldfra(i)  = cldfra(i) + dcf_mix(i)
     1317      qcld(i)    = qcld(i) + dqvc_mix(i) + dqi_mix(i)
     1318      qvc(i)     = qvc(i) + dqvc_mix(i)
     1319      clrfra(i)  = clrfra(i) - dcf_mix(i)
     1320      qclr(i)    = qclr(i) - dqvc_mix(i) - dqi_mix(i)
     1321
     1322
     1323      !---------------------------------------
     1324      !--         ICE SEDIMENTATION         --
     1325      !---------------------------------------
     1326      !
     1327      !--If ice supersaturation is activated, the cloud properties are prognostic.
     1328      !--The falling ice is then partly considered a new cloud in the gridbox.
     1329      !--BEWARE with this parameterisation, we can create a new cloud with a much
     1330      !--different ice water content and water vapor content than the existing cloud
     1331      !--(if it exists). This implies than unphysical fluxes of ice and vapor
     1332      !--occur between the existing cloud and the newly formed cloud (from sedimentation).
     1333      !--Note also that currently, the precipitation scheme assume a maximum
     1334      !--random overlap, meaning that the new formed clouds will not be affected
     1335      !--by vertical wind shear.
     1336      !
     1337      IF ( icesed_flux(i) .GT. 0. ) THEN
     1338
     1339        clrfra_sed = MIN(clrfra(i), cldfra_above(i) - cldfra(i))
     1340
     1341        IF ( ( clrfra_sed .GT. eps ) .AND. ( clrfra(i) .GT. eps ) ) THEN
     1342
     1343          qiceincld = qice_sedim(i) / cldfra_above(i)
     1344          IF ( ok_unadjusted_clouds .AND. .NOT. ok_warm_cloud ) THEN
     1345            tauinv_depsub = qiceincld**(1./3.) * kappa_depsub
     1346            qvapinclr_lim = qsat(i) - qiceincld / ( 1. - EXP( - dtime * tauinv_depsub ) )
     1347          ELSE
     1348            qvapinclr_lim = qsat(i) - qiceincld
     1349          ENDIF
     1350
     1351          rhl_clr = qclr(i) / clrfra(i) / qsatl(i) * 100.
     1352          pdf_x = qvapinclr_lim / qsatl(i) * 100.
     1353          pdf_loc = rhl_clr - pdf_scale(i) * pdf_gamma(i)
     1354          pdf_x = qsat(i) / qsatl(i) * 100.
     1355          pdf_y = LOG( MAX( ( pdf_x - pdf_loc ) / pdf_scale(i), eps) ) * pdf_alpha(i)
     1356          IF ( pdf_y .GT. 10. ) THEN !--Avoid overflows
     1357            pdf_fra_above_lim = 0.
     1358            pdf_q_above_lim = 0.
     1359          ELSEIF ( pdf_y .LT. -10. ) THEN
     1360            pdf_fra_above_lim = clrfra(i)
     1361            pdf_q_above_lim = qclr(i)
     1362          ELSE
     1363            pdf_y = EXP( pdf_y )
     1364            pdf_e3 = GAMMAINC ( 1. + 1. / pdf_alpha(i) , pdf_y )
     1365            pdf_e3 = pdf_scale(i) * ( 1. - pdf_e3 ) * pdf_gamma(i)
     1366            pdf_fra_above_lim = EXP( - pdf_y ) * clrfra(i)
     1367            pdf_q_above_lim = ( pdf_e3 * clrfra(i) &
     1368                              + pdf_loc * pdf_fra_above_lim ) * qsatl(i) / 100.
     1369          ENDIF
     1370
    11621371          IF ( pdf_fra_above_lim .GT. eps ) THEN
    11631372            dcf_sed(i)  = clrfra_sed * pdf_fra_above_lim / clrfra(i)
    … …  
    11861395
    11871396
    1188       !--We put back contrails in the clouds class
    1189       cldfra(i) = cldfra(i) + contfra(i)
    1190       qcld(i) = qcld(i) + qcont(i)
    1191       qvc(i) = qvc(i) + qsat(i) * contfra(i)
     1397      IF ( ( lincontfra(i) + circontfra(i) ) .GT. 0. ) THEN
     1398        !--We put back contrails in the clouds class
     1399        cldfra(i) = cldfra(i) + lincontfra(i) + circontfra(i)
     1400        qcld(i) = qcld(i) + qlincont(i) + qcircont(i)
     1401        qvc(i) = qvc(i) + qsat(i) * ( lincontfra(i) + circontfra(i) )
     1402      ENDIF
    11921403
    11931404
    … …  
    12671478      keepgoing, pt_pron_clds, &
    12681479      Tcritcont, qcritcont, potcontfraP, potcontfraNP, &
    1269       dcfa_ini, dqia_ini, dqta_ini)
     1480      dcfl_ini, dqil_ini, dqtl_ini)
    12701481
    12711482  DO i = 1, klon
    … …  
    12751486      IF ( cldfra(i) .GE. ( 1. - cfcon(i) - eps ) ) THEN
    12761487        contrails_conversion_factor = 1.
    1277       ELSEIF ( contfra(i) .LT. 1.e-6 ) THEN
     1488      ELSEIF ( lincontfra(i) .LT. 1.e-6 ) THEN
    12781489        contrails_conversion_factor = 1.
    12791490      ELSE
    … …  
    12861497            * ( 1. - cldfra(i) / ( 1. - cfcon(i) ) )**0.1 )
    12871498      ENDIF
    1288       dcfa_cir(i) = - contrails_conversion_factor * contfra(i)
    1289       dqta_cir(i) = - contrails_conversion_factor * qcont(i)
     1499      dcfl_cir(i) = - contrails_conversion_factor * lincontfra(i)
     1500      dqtl_cir(i) = - contrails_conversion_factor * qlincont(i)
     1501
     1502      dcfl_ini(i) = MIN(MIN(dcfl_ini(i), issrfra(i)), 1. - cfcon(i) - cldfra(i))
     1503      dqtl_ini(i) = MIN(MIN(dqtl_ini(i), qissr(i)), qtot_in(i) - qcld(i))
    12901504
    12911505      !--Add tendencies
    1292       issrfra(i) = MAX(0., issrfra(i) - dcfa_ini(i))
    1293       qissr(i)   = MAX(0., qissr(i) - dqta_ini(i))
    1294       clrfra(i)  = MAX(0., clrfra(i) - dcfa_ini(i))
    1295       qclr(i)    = MAX(0., qclr(i) - dqta_ini(i))
    1296 
    1297       cldfra(i)  = MAX(0., MIN(1. - cfcon(i), cldfra(i) + dcfa_ini(i)))
    1298       qcld(i)    = MAX(0., MIN(qtot_in(i), qcld(i) + dqta_ini(i)))
    1299       qvc(i)     = MAX(0., MIN(qcld(i), qvc(i) + dcfa_ini(i) * qsat(i)))
    1300       contfra(i) = MAX(0., MIN(cldfra(i), contfra(i) + dcfa_cir(i) + dcfa_ini(i)))
    1301       qcont(i)   = MAX(0., MIN(qcld(i), qcont(i) + dqta_cir(i) + dqta_ini(i)))
    1302       perscontfra(i) = perscontfra(i) - dcfa_cir(i)
     1506      issrfra(i) = issrfra(i) - dcfl_ini(i)
     1507      qissr(i)   = qissr(i) - dqtl_ini(i)
     1508      clrfra(i)  = clrfra(i) - dcfl_ini(i)
     1509      qclr(i)    = qclr(i) - dqtl_ini(i)
     1510
     1511      cldfra(i)  = cldfra(i) + dcfl_ini(i)
     1512      qcld(i)    = qcld(i) + dqtl_ini(i)
     1513      qvc(i)     = MIN(qcld(i), qvc(i) + dcfl_ini(i) * qsat(i))
     1514      lincontfra(i) = lincontfra(i) + dcfl_cir(i) + dcfl_ini(i)
     1515      qlincont(i)   = qlincont(i) + dqtl_cir(i) + dqtl_ini(i)
     1516      circontfra(i) = circontfra(i) - dcfl_cir(i)
     1517      qcircont(i)   = qcircont(i) - dqtl_cir(i)
    13031518
    13041519      !--Diagnostics
    1305       dcfa_ini(i) = dcfa_ini(i) / dtime
    1306       dqia_ini(i) = dqia_ini(i) / dtime
    1307       dqta_ini(i) = dqta_ini(i) / dtime
    1308       dcfa_sub(i) = dcfa_sub(i) / dtime
    1309       dqia_sub(i) = dqta_sub(i) / dtime
    1310       dqta_sub(i) = dqta_sub(i) / dtime
    1311       dcfa_cir(i) = dcfa_cir(i) / dtime
    1312       dqta_cir(i) = dqta_cir(i) / dtime
    1313       dcfa_mix(i) = dcfa_mix(i) / dtime
    1314       dqia_mix(i) = dqia_mix(i) / dtime
    1315       dqta_mix(i) = dqta_mix(i) / dtime
     1520      dcfl_ini(i) = dcfl_ini(i) / dtime
     1521      dqil_ini(i) = dqil_ini(i) / dtime
     1522      dqtl_ini(i) = dqtl_ini(i) / dtime
     1523      dcfl_sub(i) = dcfl_sub(i) / dtime
     1524      dqil_sub(i) = dqil_sub(i) / dtime
     1525      dqtl_sub(i) = dqtl_sub(i) / dtime
     1526      dcfl_cir(i) = dcfl_cir(i) / dtime
     1527      dqtl_cir(i) = dqtl_cir(i) / dtime
     1528      dcfl_mix(i) = dcfl_mix(i) / dtime
     1529      dqil_mix(i) = dqil_mix(i) / dtime
     1530      dqtl_mix(i) = dqtl_mix(i) / dtime
     1531      dcfc_sub(i) = dcfc_sub(i) / dtime
     1532      dqic_sub(i) = dqic_sub(i) / dtime
     1533      dqtc_sub(i) = dqtc_sub(i) / dtime
     1534      dcfc_mix(i) = dcfc_mix(i) / dtime
     1535      dqic_mix(i) = dqic_mix(i) / dtime
     1536      dqtc_mix(i) = dqtc_mix(i) / dtime
    13161537
    13171538      !-------------------------------------------
    … …  
    13221543        !--If the cloud is too small, it is sublimated.
    13231544        cldfra(i) = 0.
    1324         contfra(i)= 0.
    1325         perscontfra(i) = 0.
    13261545        qcld(i)   = 0.
    13271546        qvc(i)    = 0.
    1328         qcont(i)  = 0.
     1547        lincontfra(i) = 0.
     1548        qlincont(i)   = 0.
     1549        circontfra(i) = 0.
     1550        qcircont(i)   = 0.
    13291551        qincld(i) = qsat(i)
    13301552      ELSE
    … …  
    13321554      ENDIF ! cldfra .LT. eps
    13331555
    1334       IF ( contfra(i) .LT. eps ) THEN
    1335         contfra(i) = 0.
    1336         qcont(i) = 0.
     1556      IF ( lincontfra(i) .LT. eps ) THEN
     1557        lincontfra(i) = 0.
     1558        qlincont(i) = 0.
    13371559      ENDIF
    13381560
    1339       IF ( perscontfra(i) .LT. eps ) perscontfra(i) = 0.
     1561      IF ( circontfra(i) .LT. eps ) THEN
     1562        circontfra(i) = 0.
     1563        qcircont(i) = 0.
     1564      ENDIF
    13401565
    13411566    ENDIF ! keepgoing

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

    r5631 r5641  
    206206  REAL, SAVE, PROTECTED :: delta_hetfreez=0.85               ! [-] value between 0 and 1 to simulate for heterogeneous freezing.
    207207  !$OMP THREADPRIVATE(delta_hetfreez)
     208
     209  REAL, SAVE, PROTECTED :: aspect_ratio_cirrus=1./9.         ! [-] aspect ratio of natural cirrus clouds
     210  !$OMP THREADPRIVATE(aspect_ratio_cirrus)
    208211 
    209212  REAL, SAVE, PROTECTED :: coef_mixing_lscp=1.E-3            ! [-] tuning coefficient for the mixing process
    … …  
    219222  !$OMP THREADPRIVATE(ok_plane_contrail)
    220223
    221   REAL, SAVE, PROTECTED :: aspect_ratio_contrails=.1         ! [-] aspect ratio of the contrails clouds
    222   !$OMP THREADPRIVATE(aspect_ratio_contrails)
    223 
    224   REAL, SAVE, PROTECTED :: coef_mixing_contrails             ! [-] tuning coefficient for the contrails mixing process
    225   !$OMP THREADPRIVATE(coef_mixing_contrails)
    226  
    227   REAL, SAVE, PROTECTED :: coef_shear_contrails              ! [-] additional coefficient for the contrails shearing process (subprocess of the contrails mixing process)
    228   !$OMP THREADPRIVATE(coef_shear_contrails)
    229  
    230   REAL, SAVE, PROTECTED :: chi_mixing_contrails=1.           ! [-] factor for increasing the chance that moist air is surrounding contrails
    231   !$OMP THREADPRIVATE(chi_mixing_contrails)
    232 
    233   REAL, SAVE, PROTECTED :: rm_ice_crystals_contrails=7.5E-6  ! [m] geometric radius of ice crystals in contrails
    234   !$OMP THREADPRIVATE(rm_ice_crystals_contrails)
     224  REAL, SAVE, PROTECTED :: aspect_ratio_lincontrails=.1      ! [-] aspect ratio of linear contrails
     225  !$OMP THREADPRIVATE(aspect_ratio_lincontrails)
     226
     227  REAL, SAVE, PROTECTED :: coef_mixing_lincontrails          ! [-] tuning coefficient for the linear contrails mixing process
     228  !$OMP THREADPRIVATE(coef_mixing_lincontrails)
     229 
     230  REAL, SAVE, PROTECTED :: coef_shear_lincontrails           ! [-] additional coefficient for the linear contrails shearing process (subprocess of the contrails mixing process)
     231  !$OMP THREADPRIVATE(coef_shear_lincontrails)
     232 
     233  REAL, SAVE, PROTECTED :: chi_mixing_lincontrails=3.        ! [-] factor for increasing the chance that moist air is surrounding linear contrails
     234  !$OMP THREADPRIVATE(chi_mixing_lincontrails)
    235235
    236236  REAL, SAVE, PROTECTED :: EI_H2O_aviation=1.25              ! [kgH2O/kg] emission index of water vapor for a given fuel type
    … …  
    497497    CALL getin_p('b_homofreez',b_homofreez)
    498498    CALL getin_p('delta_hetfreez',delta_hetfreez)
     499    CALL getin_p('aspect_ratio_cirrus',aspect_ratio_cirrus)
    499500    CALL getin_p('coef_mixing_lscp',coef_mixing_lscp)
    500501    CALL getin_p('coef_shear_lscp',coef_shear_lscp)
    501502    ! for aviation
    502     CALL getin_p('aspect_ratio_contrails',aspect_ratio_contrails)
    503     coef_mixing_contrails=coef_mixing_lscp
    504     CALL getin_p('coef_mixing_contrails',coef_mixing_contrails)
    505     coef_shear_contrails=coef_shear_lscp
    506     CALL getin_p('coef_shear_contrails',coef_shear_contrails)
    507     CALL getin_p('chi_mixing_contrails',chi_mixing_contrails)
    508     CALL getin_p('rm_ice_crystals_contrails',rm_ice_crystals_contrails)
     503    CALL getin_p('aspect_ratio_lincontrails',aspect_ratio_lincontrails)
     504    coef_mixing_lincontrails=coef_mixing_lscp
     505    CALL getin_p('coef_mixing_lincontrails',coef_mixing_lincontrails)
     506    coef_shear_lincontrails=coef_shear_lscp
     507    CALL getin_p('coef_shear_lincontrails',coef_shear_lincontrails)
     508    CALL getin_p('chi_mixing_lincontrails',chi_mixing_lincontrails)
    509509    CALL getin_p('EI_H2O_aviation',EI_H2O_aviation)
    510510    CALL getin_p('qheat_fuel_aviation',qheat_fuel_aviation)
    … …  
    597597    WRITE(lunout,*) 'lscp_ini, b_homofreez:', b_homofreez
    598598    WRITE(lunout,*) 'lscp_ini, delta_hetfreez', delta_hetfreez
     599    WRITE(lunout,*) 'lscp_ini, aspect_ratio_cirrus:', aspect_ratio_cirrus
    599600    WRITE(lunout,*) 'lscp_ini, coef_mixing_lscp:', coef_mixing_lscp
    600601    WRITE(lunout,*) 'lscp_ini, coef_shear_lscp:', coef_shear_lscp
    601602    ! for aviation
    602     WRITE(lunout,*) 'lscp_ini, aspect_ratio_contrails:', aspect_ratio_contrails
    603     WRITE(lunout,*) 'lscp_ini, coef_mixing_contrails:', coef_mixing_contrails
    604     WRITE(lunout,*) 'lscp_ini, coef_shear_contrails:', coef_shear_contrails
    605     WRITE(lunout,*) 'lscp_ini, chi_mixing_contrails:', chi_mixing_contrails
    606     WRITE(lunout,*) 'lscp_ini, rm_ice_crystals_contrails:', rm_ice_crystals_contrails
     603    WRITE(lunout,*) 'lscp_ini, aspect_ratio_lincontrails:', aspect_ratio_lincontrails
     604    WRITE(lunout,*) 'lscp_ini, coef_mixing_lincontrails:', coef_mixing_lincontrails
     605    WRITE(lunout,*) 'lscp_ini, coef_shear_lincontrails:', coef_shear_lincontrails
     606    WRITE(lunout,*) 'lscp_ini, chi_mixing_lincontrails:', chi_mixing_lincontrails
    607607    WRITE(lunout,*) 'lscp_ini, EI_H2O_aviation:', EI_H2O_aviation
    608608    WRITE(lunout,*) 'lscp_ini, qheat_fuel_aviation:', qheat_fuel_aviation

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

    r5626 r5641  
    2323       falb_dir, falb_dif, prw_ancien, prlw_ancien, prsw_ancien, prbsw_ancien, &
    2424       ftsol, pbl_tke, pctsrf, q_ancien, ql_ancien, qs_ancien, qbs_ancien, &
    25        cf_ancien, qvc_ancien, qvcon, qccon, cfa_ancien, pcf_ancien, qva_ancien, qia_ancien, &
     25       cf_ancien, qvc_ancien, qvcon, qccon, cfl_ancien, cfc_ancien, qtl_ancien, qtc_ancien, &
    2626       radpas, radsol, rain_fall, &
    2727       ratqs, rnebcon, rugoro, sig1, snow_fall, bs_fall, solaire_etat0, sollw, sollwdown, &
    … …  
    425425  ! cas specifique des variables de l'aviation
    426426  IF ( ok_plane_contrail ) THEN
    427     ancien_ok=ancien_ok.AND.phyetat0_get(cfa_ancien,"CFAANCIEN","CFAANCIEN",0.)
    428     ancien_ok=ancien_ok.AND.phyetat0_get(pcf_ancien,"PCFANCIEN","PCFANCIEN",0.)
    429     ancien_ok=ancien_ok.AND.phyetat0_get(qva_ancien,"QVAANCIEN","QVAANCIEN",0.)
    430     ancien_ok=ancien_ok.AND.phyetat0_get(qia_ancien,"QIAANCIEN","QIAANCIEN",0.)
    431   ELSE
    432     cfa_ancien(:,:)=0.
    433     pcf_ancien(:,:)=0.
    434     qva_ancien(:,:)=0.
    435     qia_ancien(:,:)=0.
     427    ancien_ok=ancien_ok.AND.phyetat0_get(cfl_ancien,"CFLANCIEN","CFLANCIEN",0.)
     428    ancien_ok=ancien_ok.AND.phyetat0_get(cfc_ancien,"CFCANCIEN","CFCANCIEN",0.)
     429    ancien_ok=ancien_ok.AND.phyetat0_get(qtl_ancien,"QTLANCIEN","QTLANCIEN",0.)
     430    ancien_ok=ancien_ok.AND.phyetat0_get(qtc_ancien,"QTCANCIEN","QTCANCIEN",0.)
     431  ELSE
     432    cfl_ancien(:,:)=0.
     433    cfc_ancien(:,:)=0.
     434    qtl_ancien(:,:)=0.
     435    qtc_ancien(:,:)=0.
    436436  ENDIF
    437437
    … …  
    464464
    465465  IF ( ok_plane_contrail ) THEN
    466     IF ( ( maxval(cfa_ancien).EQ.minval(cfa_ancien) ) .OR. &
    467          ( maxval(pcf_ancien).EQ.minval(pcf_ancien) ) .OR. &
    468          ( maxval(qva_ancien).EQ.minval(qva_ancien) ) .OR. &
    469          ( maxval(qia_ancien).EQ.minval(qia_ancien) ) ) THEN
     466    IF ( ( maxval(cfl_ancien).EQ.minval(cfl_ancien) ) .OR. &
     467         ( maxval(cfc_ancien).EQ.minval(cfc_ancien) ) .OR. &
     468         ( maxval(qtl_ancien).EQ.minval(qtl_ancien) ) .OR. &
     469         ( maxval(qtc_ancien).EQ.minval(qtc_ancien) ) ) THEN
    470470       ancien_ok=.false.
    471471     ENDIF

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

    r5626 r5641  
    2424                                ql_ancien, qs_ancien, qbs_ancien, cf_ancien, &
    2525                                qvcon, qccon,                                &
    26                                 qvc_ancien, cfa_ancien, pcf_ancien,          &
    27                                 qva_ancien, qia_ancien, u_ancien, v_ancien,  &
     26                                qvc_ancien, cfl_ancien, cfc_ancien,          &
     27                                qtl_ancien, qtc_ancien, u_ancien, v_ancien,  &
    2828                                clwcon, rnebcon, ratqs, pbl_tke,             &
    2929                                wake_delta_pbl_tke, zmax0, f0, sig1, w01,    &
    … …  
    261261
    262262    IF ( ok_plane_contrail ) THEN
    263       CALL put_field(pass,"CFAANCIEN", "CFAANCIEN", cfa_ancien)
    264       CALL put_field(pass,"PCFANCIEN", "PCFANCIEN", pcf_ancien)
    265       CALL put_field(pass,"QVAANCIEN", "QVAANCIEN", qva_ancien)
    266       CALL put_field(pass,"QIAANCIEN", "QIAANCIEN", qia_ancien)
     263      CALL put_field(pass,"CFLANCIEN", "CFLANCIEN", cfl_ancien)
     264      CALL put_field(pass,"CFCANCIEN", "CFCANCIEN", cfc_ancien)
     265      CALL put_field(pass,"QTLANCIEN", "QTLANCIEN", qtl_ancien)
     266      CALL put_field(pass,"QTCANCIEN", "QTCANCIEN", qtc_ancien)
    267267    ENDIF
    268268

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

    r5618 r5641  
    671671      REAL, SAVE, ALLOCATABLE :: d_q_avi(:,:)
    672672      !$OMP THREADPRIVATE(d_q_avi)
    673       REAL, SAVE, ALLOCATABLE :: cfa_seri(:,:), d_cfa_dyn(:,:)
    674       !$OMP THREADPRIVATE(cfa_seri, d_cfa_dyn)
    675       REAL, SAVE, ALLOCATABLE :: pcf_seri(:,:), d_pcf_dyn(:,:)
    676       !$OMP THREADPRIVATE(pcf_seri, d_pcf_dyn)
    677       REAL, SAVE, ALLOCATABLE :: qva_seri(:,:), d_qva_dyn(:,:)
    678       !$OMP THREADPRIVATE(qva_seri, d_qva_dyn)
    679       REAL, SAVE, ALLOCATABLE :: qia_seri(:,:), d_qia_dyn(:,:)
    680       !$OMP THREADPRIVATE(qia_seri, d_qia_dyn)
     673      REAL, SAVE, ALLOCATABLE :: cfl_seri(:,:), d_cfl_dyn(:,:)
     674      !$OMP THREADPRIVATE(cfl_seri, d_cfl_dyn)
     675      REAL, SAVE, ALLOCATABLE :: cfc_seri(:,:), d_cfc_dyn(:,:)
     676      !$OMP THREADPRIVATE(cfc_seri, d_cfc_dyn)
     677      REAL, SAVE, ALLOCATABLE :: qtl_seri(:,:), d_qtl_dyn(:,:)
     678      !$OMP THREADPRIVATE(qtl_seri, d_qtl_dyn)
     679      REAL, SAVE, ALLOCATABLE :: qtc_seri(:,:), d_qtc_dyn(:,:)
     680      !$OMP THREADPRIVATE(qtc_seri, d_qtc_dyn)
    681681      REAL, SAVE, ALLOCATABLE :: flight_dist(:,:), flight_h2o(:,:)
    682682      !$OMP THREADPRIVATE(flight_dist, flight_h2o)
    … …  
    685685      REAL, SAVE, ALLOCATABLE :: potcontfraP(:,:), potcontfraNP(:,:)
    686686      !$OMP THREADPRIVATE(potcontfraP, potcontfraNP)
    687       REAL, SAVE, ALLOCATABLE :: qice_perscont(:,:)
    688       !$OMP THREADPRIVATE(qice_perscont)
    689       REAL, SAVE, ALLOCATABLE :: dcfa_cir(:,:), dqta_cir(:,:)
    690       !$OMP THREADPRIVATE(dcfa_cir, dqta_cir)
    691       REAL, SAVE, ALLOCATABLE :: dcfa_ini(:,:), dqia_ini(:,:), dqta_ini(:,:)
    692       !$OMP THREADPRIVATE(dcfa_ini, dqia_ini, dqta_ini)
    693       REAL, SAVE, ALLOCATABLE :: dcfa_sub(:,:), dqia_sub(:,:), dqta_sub(:,:)
    694       !$OMP THREADPRIVATE(dcfa_sub, dqia_sub, dqta_sub)
    695       REAL, SAVE, ALLOCATABLE :: dcfa_mix(:,:), dqia_mix(:,:), dqta_mix(:,:)
    696       !$OMP THREADPRIVATE(dcfa_mix, dqia_mix, dqta_mix)
     687      REAL, SAVE, ALLOCATABLE :: qice_lincont(:,:), qice_circont(:,:)
     688      !$OMP THREADPRIVATE(qice_lincont, qice_circont)
     689      REAL, SAVE, ALLOCATABLE :: dcfl_cir(:,:), dqtl_cir(:,:)
     690      !$OMP THREADPRIVATE(dcfl_cir, dqtl_cir)
     691      REAL, SAVE, ALLOCATABLE :: dcfl_ini(:,:), dqil_ini(:,:), dqtl_ini(:,:)
     692      !$OMP THREADPRIVATE(dcfl_ini, dqil_ini, dqtl_ini)
     693      REAL, SAVE, ALLOCATABLE :: dcfl_sub(:,:), dqil_sub(:,:), dqtl_sub(:,:)
     694      !$OMP THREADPRIVATE(dcfl_sub, dqil_sub, dqtl_sub)
     695      REAL, SAVE, ALLOCATABLE :: dcfl_mix(:,:), dqil_mix(:,:), dqtl_mix(:,:)
     696      !$OMP THREADPRIVATE(dcfl_mix, dqil_mix, dqtl_mix)
     697      REAL, SAVE, ALLOCATABLE :: dcfc_sub(:,:), dqic_sub(:,:), dqtc_sub(:,:)
     698      !$OMP THREADPRIVATE(dcfc_sub, dqic_sub, dqtc_sub)
     699      REAL, SAVE, ALLOCATABLE :: dcfc_mix(:,:), dqic_mix(:,:), dqtc_mix(:,:)
     700      !$OMP THREADPRIVATE(dcfc_mix, dqic_mix, dqtc_mix)
    697701      REAL, SAVE, ALLOCATABLE :: cldfra_nocont(:,:), cldtau_nocont(:,:), cldemi_nocont(:,:)
    698702      !$OMP THREADPRIVATE(cldfra_nocont, cldtau_nocont, cldemi_nocont)
    … …  
    12561260!-- LSCP - aviation and contrails variables
    12571261      ALLOCATE(d_q_avi(klon,klev))
    1258       ALLOCATE(cfa_seri(klon,klev), d_cfa_dyn(klon,klev))
    1259       ALLOCATE(pcf_seri(klon,klev), d_pcf_dyn(klon,klev))
    1260       ALLOCATE(qva_seri(klon,klev), d_qva_dyn(klon,klev))
    1261       ALLOCATE(qia_seri(klon,klev), d_qia_dyn(klon,klev))
     1262      ALLOCATE(cfl_seri(klon,klev), d_cfl_dyn(klon,klev))
     1263      ALLOCATE(cfc_seri(klon,klev), d_cfc_dyn(klon,klev))
     1264      ALLOCATE(qtl_seri(klon,klev), d_qtl_dyn(klon,klev))
     1265      ALLOCATE(qtc_seri(klon,klev), d_qtc_dyn(klon,klev))
    12621266      ALLOCATE(flight_dist(klon,klev), flight_h2o(klon,klev))
    12631267      ALLOCATE(Tcritcont(klon,klev), qcritcont(klon,klev))
    12641268      ALLOCATE(potcontfraP(klon,klev), potcontfraNP(klon,klev))
    1265       ALLOCATE(qice_perscont(klon,klev))
    1266       ALLOCATE(dcfa_cir(klon,klev), dqta_cir(klon,klev))
    1267       ALLOCATE(dcfa_ini(klon,klev), dqia_ini(klon,klev), dqta_ini(klon,klev))
    1268       ALLOCATE(dcfa_sub(klon,klev), dqia_sub(klon,klev), dqta_sub(klon,klev))
    1269       ALLOCATE(dcfa_mix(klon,klev), dqia_mix(klon,klev), dqta_mix(klon,klev))
     1269      ALLOCATE(qice_lincont(klon,klev), qice_circont(klon,klev))
     1270      ALLOCATE(dcfl_cir(klon,klev), dqtl_cir(klon,klev))
     1271      ALLOCATE(dcfl_ini(klon,klev), dqil_ini(klon,klev), dqtl_ini(klon,klev))
     1272      ALLOCATE(dcfl_sub(klon,klev), dqil_sub(klon,klev), dqtl_sub(klon,klev))
     1273      ALLOCATE(dcfl_mix(klon,klev), dqil_mix(klon,klev), dqtl_mix(klon,klev))
     1274      ALLOCATE(dcfc_sub(klon,klev), dqic_sub(klon,klev), dqtc_sub(klon,klev))
     1275      ALLOCATE(dcfc_mix(klon,klev), dqic_mix(klon,klev), dqtc_mix(klon,klev))
    12701276      ALLOCATE(cldfra_nocont(klon,klev), cldtau_nocont(klon,klev), cldemi_nocont(klon,klev))
    12711277      ALLOCATE(cldh_nocont(klon), contcov(klon), conttau(klon,klev), contemi(klon,klev))
    … …  
    16781684
    16791685!-- LSCP - aviation and contrails variables
    1680       DEALLOCATE(cfa_seri, d_cfa_dyn, pcf_seri, d_pcf_dyn)
    1681       DEALLOCATE(qva_seri, d_qva_dyn, qia_seri, d_qia_dyn)
     1686      DEALLOCATE(cfl_seri, d_cfl_dyn, cfc_seri, d_cfc_dyn)
     1687      DEALLOCATE(qtl_seri, d_qtl_dyn, qtc_seri, d_qtc_dyn)
    16821688      DEALLOCATE(d_q_avi, flight_dist, flight_h2o)
    16831689      DEALLOCATE(Tcritcont, qcritcont, potcontfraP, potcontfraNP)
    1684       DEALLOCATE(qice_perscont, dcfa_cir, dqta_cir, dcfa_ini, dqia_ini, dqta_ini)
    1685       DEALLOCATE(dcfa_sub, dqia_sub, dqta_sub, dcfa_mix, dqia_mix, dqta_mix)
     1690      DEALLOCATE(qice_lincont, qice_circont, dcfl_cir, dqtl_cir, dcfl_ini, dqil_ini)
     1691      DEALLOCATE(dqtl_ini, dcfl_sub, dqil_sub, dqtl_sub, dcfl_mix, dqil_mix, dqtl_mix)
     1692      DEALLOCATE(dcfc_sub, dqic_sub, dqtc_sub, dcfc_mix, dqic_mix, dqtc_mix)
    16861693      DEALLOCATE(cldfra_nocont, cldtau_nocont, cldemi_nocont, conttau, contemi)
    16871694      DEALLOCATE(cldh_nocont, contcov, fiwp_nocont, fiwc_nocont, ref_ice_nocont)

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

    r5618 r5641  
    21972197  TYPE(ctrl_out), SAVE :: o_dqavi = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
    21982198    'dqavi', 'Water vapor emissions from aviation tendency', 'kg/kg/s', (/ ('',i=1,10) /))
    2199   TYPE(ctrl_out), SAVE :: o_cfaseri = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
    2200     'cfaseri', 'Linear contrail fraction', '-', (/ ('',i=1,10) /))
    2201   TYPE(ctrl_out), SAVE :: o_dcfadyn = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
    2202     'dcfadyn', 'Dynamics linear contrail fraction tendency', 's-1', (/ ('',i=1,10) /))
    2203   TYPE(ctrl_out), SAVE :: o_pcfseri = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
    2204     'pcfseri', 'Contrail induced cirrus fraction', '-', (/ ('',i=1,10) /))
    2205   TYPE(ctrl_out), SAVE :: o_dpcfdyn = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
    2206     'dpcfdyn', 'Dynamics contrail induced cirrus fraction tendency', 's-1', (/ ('',i=1,10) /))
    2207   TYPE(ctrl_out), SAVE :: o_qvaseri = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
    2208     'qvaseri', 'Linear contrail total specific humidity', 'kg/kg', (/ ('',i=1,10) /))
    2209   TYPE(ctrl_out), SAVE :: o_dqvadyn = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
    2210     'dqvadyn', 'Dynamics linear contrail total specific humidity tendency', 'kg/kg/s', (/ ('',i=1,10) /))
    2211   TYPE(ctrl_out), SAVE :: o_qiaseri = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
    2212     'qiaseri', 'Linear contrail total specific humidity', 'kg/kg', (/ ('',i=1,10) /))
    2213   TYPE(ctrl_out), SAVE :: o_dqiadyn = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
    2214     'dqiadyn', 'Dynamics linear contrail total specific humidity tendency', 'kg/kg/s', (/ ('',i=1,10) /))
     2199  TYPE(ctrl_out), SAVE :: o_cflseri = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
     2200    'cflseri', 'Linear contrail fraction', '-', (/ ('',i=1,10) /))
     2201  TYPE(ctrl_out), SAVE :: o_dcfldyn = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
     2202    'dcfldyn', 'Dynamics linear contrail fraction tendency', 's-1', (/ ('',i=1,10) /))
     2203  TYPE(ctrl_out), SAVE :: o_cfcseri = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
     2204    'cfcseri', 'Contrail cirrus fraction', '-', (/ ('',i=1,10) /))
     2205  TYPE(ctrl_out), SAVE :: o_dcfcdyn = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
     2206    'dcfcdyn', 'Dynamics contrail cirrus fraction tendency', 's-1', (/ ('',i=1,10) /))
     2207  TYPE(ctrl_out), SAVE :: o_qtlseri = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
     2208    'qtlseri', 'Linear contrail total specific humidity', 'kg/kg', (/ ('',i=1,10) /))
     2209  TYPE(ctrl_out), SAVE :: o_dqtldyn = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
     2210    'dqtldyn', 'Dynamics linear contrail total specific humidity tendency', 'kg/kg/s', (/ ('',i=1,10) /))
     2211  TYPE(ctrl_out), SAVE :: o_qtcseri = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
     2212    'qtcseri', 'Contrail cirrus total specific humidity', 'kg/kg', (/ ('',i=1,10) /))
     2213  TYPE(ctrl_out), SAVE :: o_dqtcdyn = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
     2214    'dqtcdyn', 'Dynamics contrail cirrus total specific humidity tendency', 'kg/kg/s', (/ ('',i=1,10) /))
    22152215  TYPE(ctrl_out), SAVE :: o_Tcritcont = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
    22162216    'Tcritcont', 'Temperature threshold for contrail formation', 'K', (/ ('',i=1,10) /))
    … …  
    22212221  TYPE(ctrl_out), SAVE :: o_potcontfraNP = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/),&
    22222222    'potcontfraNP', 'Potential non-persistent contrail fraction', '-', (/ ('', i=1,10)/))
    2223   TYPE(ctrl_out), SAVE :: o_qice_perscont = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    2224     'qice_perscont', 'Contrail induced cirrus ice specific humidity', 'kg/kg', (/ ('', i=1, 10) /))
    2225   TYPE(ctrl_out), SAVE :: o_dcfaini = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    2226     'dcfaini', 'Initial formation linear contrail fraction tendency', 's-1', (/ ('', i=1, 10) /))
    2227   TYPE(ctrl_out), SAVE :: o_dqiaini = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    2228     'dqiaini', 'Initial formation linear contrail ice specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
    2229   TYPE(ctrl_out), SAVE :: o_dqtaini = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    2230     'dqtaini', 'Initial formation linear contrail total specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
    2231   TYPE(ctrl_out), SAVE :: o_dcfacir = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    2232     'dcfacir', 'Conversion to cirrus linear contrail fraction tendency', 's-1', (/ ('', i=1, 10) /))
    2233   TYPE(ctrl_out), SAVE :: o_dqtacir = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    2234     'dqtacir', 'Conversion to cirrus linear contrail total specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
    2235   TYPE(ctrl_out), SAVE :: o_dcfasub = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    2236     'dcfasub', 'Sublimation linear contrail fraction tendency', 's-1', (/ ('', i=1, 10) /))
    2237   TYPE(ctrl_out), SAVE :: o_dqiasub = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    2238     'dqiasub', 'Sublimation linear contrail ice specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
    2239   TYPE(ctrl_out), SAVE :: o_dqtasub = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    2240     'dqtasub', 'Sublimation linear contrail total specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
    2241   TYPE(ctrl_out), SAVE :: o_dcfamix = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    2242     'dcfamix', 'Mixing linear contrail fraction tendency', 's-1', (/ ('', i=1, 10) /))
    2243   TYPE(ctrl_out), SAVE :: o_dqiamix = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    2244     'dqiamix', 'Mixing linear contrail ice specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
    2245   TYPE(ctrl_out), SAVE :: o_dqtamix = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    2246     'dqtamix', 'Mixing linear contrail total specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
     2223  TYPE(ctrl_out), SAVE :: o_qice_lincont = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2224    'qice_lincont', 'Linear contrails ice specific humidity', 'kg/kg', (/ ('', i=1, 10) /))
     2225  TYPE(ctrl_out), SAVE :: o_qice_circont = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2226    'qice_circont', 'Contrail cirrus ice specific humidity', 'kg/kg', (/ ('', i=1, 10) /))
     2227  TYPE(ctrl_out), SAVE :: o_dcflini = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2228    'dcflini', 'Initial formation linear contrail fraction tendency', 's-1', (/ ('', i=1, 10) /))
     2229  TYPE(ctrl_out), SAVE :: o_dqilini = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2230    'dqilini', 'Initial formation linear contrail ice specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
     2231  TYPE(ctrl_out), SAVE :: o_dqtlini = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2232    'dqtlini', 'Initial formation linear contrail total specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
     2233  TYPE(ctrl_out), SAVE :: o_dcflcir = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2234    'dcflcir', 'Conversion to cirrus linear contrail fraction tendency', 's-1', (/ ('', i=1, 10) /))
     2235  TYPE(ctrl_out), SAVE :: o_dqtlcir = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2236    'dqtlcir', 'Conversion to cirrus linear contrail total specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
     2237  TYPE(ctrl_out), SAVE :: o_dcflsub = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2238    'dcflsub', 'Sublimation linear contrail fraction tendency', 's-1', (/ ('', i=1, 10) /))
     2239  TYPE(ctrl_out), SAVE :: o_dqilsub = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2240    'dqilsub', 'Sublimation linear contrail ice specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
     2241  TYPE(ctrl_out), SAVE :: o_dqtlsub = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2242    'dqtlsub', 'Sublimation linear contrail total specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
     2243  TYPE(ctrl_out), SAVE :: o_dcflmix = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2244    'dcflmix', 'Mixing linear contrail fraction tendency', 's-1', (/ ('', i=1, 10) /))
     2245  TYPE(ctrl_out), SAVE :: o_dqilmix = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2246    'dqilmix', 'Mixing linear contrail ice specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
     2247  TYPE(ctrl_out), SAVE :: o_dqtlmix = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2248    'dqtlmix', 'Mixing linear contrail total specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
     2249  TYPE(ctrl_out), SAVE :: o_dcfcsub = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2250    'dcfcsub', 'Sublimation contrail cirrus fraction tendency', 's-1', (/ ('', i=1, 10) /))
     2251  TYPE(ctrl_out), SAVE :: o_dqicsub = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2252    'dqicsub', 'Sublimation contrail cirrus ice specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
     2253  TYPE(ctrl_out), SAVE :: o_dqtcsub = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2254    'dqtcsub', 'Sublimation contrail cirrus total specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
     2255  TYPE(ctrl_out), SAVE :: o_dcfcmix = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2256    'dcfcmix', 'Mixing contrail cirrus fraction tendency', 's-1', (/ ('', i=1, 10) /))
     2257  TYPE(ctrl_out), SAVE :: o_dqicmix = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2258    'dqicmix', 'Mixing contrail cirrus ice specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
     2259  TYPE(ctrl_out), SAVE :: o_dqtcmix = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     2260    'dqtcmix', 'Mixing contrail cirrus total specific humidity tendency', 'kg/kg/s', (/ ('', i=1, 10) /))
    22472261  TYPE(ctrl_out), SAVE :: o_flight_dist = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
    22482262    'flightdist', 'Aviation flown distance', 'm/s/m^3', (/ ('', i=1, 10) /))

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

    r5639 r5641  
    232232         o_issrfra250to300, o_issrfra300to400, o_issrfra400to500, &
    233233!-- LSCP - aviation variables
    234          o_cfaseri, o_dcfadyn, o_pcfseri, o_dpcfdyn, &
    235          o_qvaseri, o_dqvadyn, o_qiaseri, o_dqiadyn, o_dqavi, &
     234         o_cflseri, o_dcfldyn, o_cfcseri, o_dcfcdyn, &
     235         o_qtlseri, o_dqtldyn, o_qtcseri, o_dqtcdyn, o_dqavi, &
    236236         o_Tcritcont, o_qcritcont, o_potcontfraP, o_potcontfraNP, &
    237          o_flight_dist, o_flight_h2o, o_qice_perscont, o_dcfacir, o_dqtacir, &
    238          o_dcfaini, o_dqiaini, o_dqtaini, o_dcfasub, o_dqiasub, o_dqtasub, &
    239          o_dcfamix, o_dqiamix, o_dqtamix, &
     237         o_flight_dist, o_flight_h2o, o_qice_lincont, o_qice_circont, o_dcflcir, o_dqtlcir, &
     238         o_dcflini, o_dqilini, o_dqtlini, o_dcflsub, o_dqilsub, o_dqtlsub, &
     239         o_dcflmix, o_dqilmix, o_dqtlmix, &
     240         o_dcfcsub, o_dqicsub, o_dqtcsub, o_dcfcmix, o_dqicmix, o_dqtcmix, &
    240241         o_cldfra_nocont, o_cldtau_nocont, o_cldemi_nocont, o_cldh_nocont, &
    241242         o_contcov, o_conttau, o_contemi, o_iwp_nocont, o_iwc_nocont, o_ref_ice_nocont, &
    … …  
    365366         issrfra100to150, issrfra150to200, issrfra200to250, &
    366367         issrfra250to300, issrfra300to400, issrfra400to500, &
    367          cfa_seri, d_cfa_dyn, pcf_seri, d_pcf_dyn, &
    368          qva_seri, d_qva_dyn, qia_seri, d_qia_dyn, d_q_avi, &
     368         cfl_seri, d_cfl_dyn, cfc_seri, d_cfc_dyn, &
     369         qtl_seri, d_qtl_dyn, qtc_seri, d_qtc_dyn, d_q_avi, &
    369370         Tcritcont, qcritcont, potcontfraP, potcontfraNP, &
    370          flight_dist, flight_h2o, qice_perscont, dcfa_cir, dqta_cir, &
    371          dcfa_ini, dqia_ini, dqta_ini, dcfa_sub, dqia_sub, dqta_sub, &
    372          dcfa_mix, dqia_mix, dqta_mix, &
     371         flight_dist, flight_h2o, qice_lincont, qice_circont, dcfl_cir, dqtl_cir, &
     372         dcfl_ini, dqil_ini, dqtl_ini, dcfl_sub, dqil_sub, dqtl_sub, &
     373         dcfl_mix, dqil_mix, dqtl_mix, &
     374         dcfc_sub, dqic_sub, dqtc_sub, dcfc_mix, dqic_mix, dqtc_mix, &
    373375         cldfra_nocont, cldtau_nocont, cldemi_nocont, cldh_nocont, &
    374376         contcov, conttau, contemi, fiwp_nocont, fiwc_nocont, ref_ice_nocont, &
    … …  
    21612163!-- LSCP - aviation variables
    21622164       IF (ok_plane_contrail) THEN
    2163          CALL histwrite_phy(o_cfaseri, cfa_seri)
    2164          CALL histwrite_phy(o_dcfadyn, d_cfa_dyn)
    2165          CALL histwrite_phy(o_pcfseri, pcf_seri)
    2166          CALL histwrite_phy(o_dpcfdyn, d_pcf_dyn)
    2167          CALL histwrite_phy(o_qvaseri, qva_seri)
    2168          CALL histwrite_phy(o_dqvadyn, d_qva_dyn)
    2169          CALL histwrite_phy(o_qiaseri, qia_seri)
    2170          CALL histwrite_phy(o_dqiadyn, d_qia_dyn)
     2165         CALL histwrite_phy(o_cflseri, cfl_seri)
     2166         CALL histwrite_phy(o_dcfldyn, d_cfl_dyn)
     2167         CALL histwrite_phy(o_cfcseri, cfc_seri)
     2168         CALL histwrite_phy(o_dcfcdyn, d_cfc_dyn)
     2169         CALL histwrite_phy(o_qtlseri, qtl_seri)
     2170         CALL histwrite_phy(o_dqtldyn, d_qtl_dyn)
     2171         CALL histwrite_phy(o_qtcseri, qtc_seri)
     2172         CALL histwrite_phy(o_dqtcdyn, d_qtc_dyn)
    21712173         CALL histwrite_phy(o_flight_dist, flight_dist)
    21722174         CALL histwrite_phy(o_Tcritcont, Tcritcont)
    … …  
    21742176         CALL histwrite_phy(o_potcontfraP, potcontfraP)
    21752177         CALL histwrite_phy(o_potcontfraNP, potcontfraNP)
    2176          CALL histwrite_phy(o_qice_perscont, qice_perscont)
    2177          CALL histwrite_phy(o_dcfacir, dcfa_cir)
    2178          CALL histwrite_phy(o_dqtacir, dqta_cir)
    2179          CALL histwrite_phy(o_dcfaini, dcfa_ini)
    2180          CALL histwrite_phy(o_dqiaini, dqia_ini)
    2181          CALL histwrite_phy(o_dqtaini, dqta_ini)
    2182          CALL histwrite_phy(o_dcfasub, dcfa_sub)
    2183          CALL histwrite_phy(o_dqiasub, dqia_sub)
    2184          CALL histwrite_phy(o_dqtasub, dqta_sub)
    2185          CALL histwrite_phy(o_dcfamix, dcfa_mix)
    2186          CALL histwrite_phy(o_dqiamix, dqia_mix)
    2187          CALL histwrite_phy(o_dqtamix, dqta_mix)
     2178         CALL histwrite_phy(o_qice_lincont, qice_lincont)
     2179         CALL histwrite_phy(o_qice_circont, qice_circont)
     2180         CALL histwrite_phy(o_dcflcir, dcfl_cir)
     2181         CALL histwrite_phy(o_dqtlcir, dqtl_cir)
     2182         CALL histwrite_phy(o_dcflini, dcfl_ini)
     2183         CALL histwrite_phy(o_dqilini, dqil_ini)
     2184         CALL histwrite_phy(o_dqtlini, dqtl_ini)
     2185         CALL histwrite_phy(o_dcflsub, dcfl_sub)
     2186         CALL histwrite_phy(o_dqilsub, dqil_sub)
     2187         CALL histwrite_phy(o_dqtlsub, dqtl_sub)
     2188         CALL histwrite_phy(o_dcflmix, dcfl_mix)
     2189         CALL histwrite_phy(o_dqilmix, dqil_mix)
     2190         CALL histwrite_phy(o_dqtlmix, dqtl_mix)
     2191         CALL histwrite_phy(o_dcfcsub, dcfc_sub)
     2192         CALL histwrite_phy(o_dqicsub, dqic_sub)
     2193         CALL histwrite_phy(o_dqtcsub, dqtc_sub)
     2194         CALL histwrite_phy(o_dcfcmix, dcfc_mix)
     2195         CALL histwrite_phy(o_dqicmix, dqic_mix)
     2196         CALL histwrite_phy(o_dqtcmix, dqtc_mix)
    21882197         CALL histwrite_phy(o_cldfra_nocont, cldfra_nocont)
    21892198         CALL histwrite_phy(o_cldtau_nocont, cldtau_nocont)

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

    r5626 r5641  
    9696      REAL, ALLOCATABLE, SAVE :: qvcon(:,:), qccon(:,:)
    9797!$OMP THREADPRIVATE(qvcon, qccon)
    98       REAL, ALLOCATABLE, SAVE :: cfa_ancien(:,:), pcf_ancien(:,:)
    99 !$OMP THREADPRIVATE(cfa_ancien, pcf_ancien)
    100       REAL, ALLOCATABLE, SAVE :: qva_ancien(:,:), qia_ancien(:,:)
    101 !$OMP THREADPRIVATE(qva_ancien, qia_ancien)
     98      REAL, ALLOCATABLE, SAVE :: cfl_ancien(:,:), cfc_ancien(:,:)
     99!$OMP THREADPRIVATE(cfl_ancien, cfc_ancien)
     100      REAL, ALLOCATABLE, SAVE :: qtl_ancien(:,:), qtc_ancien(:,:)
     101!$OMP THREADPRIVATE(qtl_ancien, qtc_ancien)
    102102!!! RomP >>>
    103103      REAL, ALLOCATABLE, SAVE :: tr_ancien(:,:,:)
    … …  
    594594      ALLOCATE(u_ancien(klon,klev), v_ancien(klon,klev))
    595595      ALLOCATE(cf_ancien(klon,klev), qvc_ancien(klon,klev))
    596       ALLOCATE(cfa_ancien(klon,klev), pcf_ancien(klon,klev))
    597       ALLOCATE(qva_ancien(klon,klev), qia_ancien(klon,klev))
     596      ALLOCATE(cfl_ancien(klon,klev), cfc_ancien(klon,klev))
     597      ALLOCATE(qtl_ancien(klon,klev), qtc_ancien(klon,klev))
    598598      ALLOCATE(qvcon(klon,klev), qccon(klon,klev))
    599599!!! Rom P >>>
    … …  
    824824      DEALLOCATE(qtc_cv,sigt_cv,detrain_cv,fm_cv)
    825825      DEALLOCATE(u_ancien, v_ancien)
    826       DEALLOCATE(cf_ancien, qvc_ancien, cfa_ancien, pcf_ancien)
    827       DEALLOCATE(qva_ancien, qia_ancien)
     826      DEALLOCATE(cf_ancien, qvc_ancien)
     827      DEALLOCATE(cfl_ancien, cfc_ancien, qtl_ancien, qtc_ancien)
    828828      DEALLOCATE(qvcon, qccon)
    829829      DEALLOCATE(tr_ancien)                           !RomP

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

    r5637 r5641  
    3939    USE ioipsl_getin_p_mod, ONLY : getin_p
    4040    USE indice_sol_mod
    41     USE infotrac_phy, ONLY: nqtot, nbtr, nqo, tracers, type_trac, addPhase, ivap, iliq, isol, ibs, icf, iqvc, icfa, ipcf, iqva, iqia
     41    USE infotrac_phy, ONLY: nqtot, nbtr, nqo, tracers, type_trac, addPhase, ivap, iliq, isol, ibs, icf, iqvc, icfl, icfc, iqtl, iqtc
    4242    USE strings_mod,  ONLY: strIdx
    4343    USE iophy
    … …  
    333333       dqvc_adj, dqvc_sub, dqvc_con, dqvc_mix, qsatliq, qsatice, &
    334334       !-- LSCP - aviation and contrails variables
    335        cfa_seri, pcf_seri, qva_seri, qia_seri, d_cfa_dyn, d_pcf_dyn, d_qva_dyn, d_qia_dyn, &
    336        d_q_avi, flight_dist, flight_h2o, qice_perscont, &
     335       cfl_seri, cfc_seri, qtl_seri, qtc_seri, d_cfl_dyn, d_cfc_dyn, d_qtl_dyn, d_qtc_dyn, &
     336       d_q_avi, flight_dist, flight_h2o, qice_lincont, qice_circont, &
    337337       Tcritcont, qcritcont, potcontfraP, potcontfraNP, &
    338338       cldfra_nocont, cldtau_nocont, cldemi_nocont, cldh_nocont, &
    … …  
    14321432       ENDIF
    14331433
    1434        IF (ok_plane_contrail.AND.((icfa.EQ.0).OR.(iqva.EQ.0).OR.(iqia.EQ.0).OR.(ipcf.EQ.0))) THEN
    1435           WRITE (lunout, *) ' ok_plane_contrail=y requires 8 tracers ', &
    1436               '(H2O_g, H2O_l, H2O_s, CLDFRA, CLDVAP_g, CONTFRA, PERSCONTFRA, CONTWATER_s) ', &
     1434       IF (ok_plane_contrail.AND.((icfl.EQ.0).OR.(icfc.EQ.0).OR.(iqtl.EQ.0).OR.(iqtc.EQ.0))) THEN
     1435          WRITE (lunout, *) ' ok_plane_contrail=y requires 9 tracers ', &
     1436              '(H2O_g, H2O_l, H2O_s, CLDFRA, CLDVAP, LINCONTFRA, ', &
     1437              'CIRCONTFRA, LINCONTWATER, CIRCONTWATER) ', &
    14371438              'but not all are provided. Might as well stop here.'
    14381439          abort_message='see above'
    … …  
    24612462          cf_seri(i,k) = 0.
    24622463          qvc_seri(i,k)= 0.
    2463           cfa_seri(i,k)= 0.
    2464           pcf_seri(i,k)= 0.
    2465           qva_seri(i,k)= 0.
    2466           qia_seri(i,k)= 0.
     2464          cfl_seri(i,k)= 0.
     2465          cfc_seri(i,k)= 0.
     2466          qtl_seri(i,k)= 0.
     2467          qtc_seri(i,k)= 0.
    24672468          !CR: ATTENTION, on rajoute la variable glace
    24682469          IF (nqo .GE. 3) THEN        !--vapour, liquid and ice
    … …  
    24792480          ENDIF
    24802481          IF (ok_plane_contrail) THEN
    2481              cfa_seri(i,k) = qx(i,k,icfa)
    2482              pcf_seri(i,k) = qx(i,k,ipcf)
    2483              qva_seri(i,k) = qx(i,k,iqva)
    2484              qia_seri(i,k) = qx(i,k,iqia)
     2482             cfl_seri(i,k) = qx(i,k,icfl)
     2483             cfc_seri(i,k) = qx(i,k,icfc)
     2484             qtl_seri(i,k) = qx(i,k,iqtl)
     2485             qtc_seri(i,k) = qx(i,k,iqtc)
    24852486             !--DYNAMICO can return NaNs for children tracers
    2486              IF (ISNAN(cfa_seri(i,k))) cfa_seri(i,k) = 0.
    2487              IF (ISNAN(pcf_seri(i,k))) pcf_seri(i,k) = 0.
    2488              IF (ISNAN(qva_seri(i,k))) qva_seri(i,k) = 0.
    2489              IF (ISNAN(qia_seri(i,k))) qia_seri(i,k) = 0.
     2487             IF (ISNAN(cfl_seri(i,k))) cfl_seri(i,k) = 0.
     2488             IF (ISNAN(cfc_seri(i,k))) cfc_seri(i,k) = 0.
     2489             IF (ISNAN(qtl_seri(i,k))) qtl_seri(i,k) = 0.
     2490             IF (ISNAN(qtc_seri(i,k))) qtc_seri(i,k) = 0.
    24902491          ENDIF
    24912492       ENDDO
    … …  
    25622563       d_cf_dyn(:,:) = (cf_seri(:,:)-cf_ancien(:,:))/phys_tstep
    25632564       d_qvc_dyn(:,:)= (qvc_seri(:,:)-qvc_ancien(:,:))/phys_tstep
    2564        d_cfa_dyn(:,:)= (cfa_seri(:,:)-cfa_ancien(:,:))/phys_tstep
    2565        d_pcf_dyn(:,:)= (pcf_seri(:,:)-pcf_ancien(:,:))/phys_tstep
    2566        d_qva_dyn(:,:)= (qva_seri(:,:)-qva_ancien(:,:))/phys_tstep
    2567        d_qia_dyn(:,:)= (qia_seri(:,:)-qia_ancien(:,:))/phys_tstep
     2565       d_cfl_dyn(:,:)= (cfl_seri(:,:)-cfl_ancien(:,:))/phys_tstep
     2566       d_cfc_dyn(:,:)= (cfc_seri(:,:)-cfc_ancien(:,:))/phys_tstep
     2567       d_qtl_dyn(:,:)= (qtl_seri(:,:)-qtl_ancien(:,:))/phys_tstep
     2568       d_qtc_dyn(:,:)= (qtc_seri(:,:)-qtc_ancien(:,:))/phys_tstep
    25682569       CALL water_int(klon,klev,q_seri,zmasse,zx_tmp_fi2d)
    25692570       d_q_dyn2d(:)=(zx_tmp_fi2d(:)-prw_ancien(:))/phys_tstep
    … …  
    25872588       d_cf_dyn(:,:) = 0.0
    25882589       d_qvc_dyn(:,:)= 0.0
    2589        d_cfa_dyn(:,:)= 0.0
    2590        d_pcf_dyn(:,:)= 0.0
    2591        d_qva_dyn(:,:)= 0.0
    2592        d_qia_dyn(:,:)= 0.0
     2590       d_cfl_dyn(:,:)= 0.0
     2591       d_cfc_dyn(:,:)= 0.0
     2592       d_qtl_dyn(:,:)= 0.0
     2593       d_qtc_dyn(:,:)= 0.0
    25932594       d_q_dyn2d(:)  = 0.0
    25942595       d_ql_dyn2d(:) = 0.0
    … …  
    27152716            qi_seri_lscp(i,k) = qs_seri(i,k) / ( q_seri(i,k) + ql_seri(i,k) + qs_seri(i,k) )
    27162717            qvc_seri(i,k) = qvc_seri(i,k) / ( q_seri(i,k) + ql_seri(i,k) + qs_seri(i,k) )
    2717             qva_seri(i,k) = qva_seri(i,k) / ( q_seri(i,k) + ql_seri(i,k) + qs_seri(i,k) )
    2718             qia_seri(i,k) = qia_seri(i,k) / ( q_seri(i,k) + ql_seri(i,k) + qs_seri(i,k) )
     2718            qtl_seri(i,k) = qtl_seri(i,k) / ( q_seri(i,k) + ql_seri(i,k) + qs_seri(i,k) )
     2719            qtc_seri(i,k) = qtc_seri(i,k) / ( q_seri(i,k) + ql_seri(i,k) + qs_seri(i,k) )
    27192720          ELSE
    27202721            ql_seri_lscp(i,k) = 0.
    27212722            qi_seri_lscp(i,k) = 0.
    27222723            qvc_seri(i,k) = 0.
    2723             qva_seri(i,k) = 0.
    2724             qia_seri(i,k) = 0.
     2724            qtl_seri(i,k) = 0.
     2725            qtc_seri(i,k) = 0.
    27252726          ENDIF
    27262727        ENDDO
    … …  
    39583959      DO k = 1, klev
    39593960        DO i = 1, klon
    3960           qva_seri(i,k) = qva_seri(i,k) * q_seri(i,k)
    3961           qia_seri(i,k) = qia_seri(i,k) * q_seri(i,k)
     3961          qtl_seri(i,k) = qtl_seri(i,k) * q_seri(i,k)
     3962          qtc_seri(i,k) = qtc_seri(i,k) * q_seri(i,k)
    39623963        ENDDO
    39633964      ENDDO
    … …  
    39863987         dcf_sub, dcf_con, dcf_mix, dqi_adj, dqi_sub, dqi_con, dqi_mix, &
    39873988         dqvc_adj, dqvc_sub, dqvc_con, dqvc_mix, qsatliq, qsatice, &
    3988          cfa_seri, pcf_seri, qva_seri, qia_seri, flight_dist, flight_h2o, &
    3989          qice_perscont, Tcritcont, qcritcont, potcontfraP, potcontfraNP, &
     3989         cfl_seri, cfc_seri, qtl_seri, qtc_seri, flight_dist, flight_h2o, &
     3990         qice_lincont, qice_circont, Tcritcont, qcritcont, potcontfraP, potcontfraNP, &
    39903991         cloudth_sth,cloudth_senv,cloudth_sigmath,cloudth_sigmaenv, &
    39913992         qraindiag, qsnowdiag, dqreva, dqssub, dqrauto, dqrcol, dqrmelt, &
    … …  
    45394540               zfice, dNovrN, ptconv, rnebcon, clwcon, &
    45404541               !--AB contrails
    4541                cfa_seri, pcf_seri, qia_seri, qice_perscont, cldfra_nocont, &
     4542               cfl_seri, cfc_seri, qice_lincont, qice_circont, cldfra_nocont, &
    45424543               cldtau_nocont, cldemi_nocont, conttau, contemi, cldh_nocont, contcov, &
    45434544               fiwp_nocont, fiwc_nocont, ref_ice_nocont)
    … …  
    57465747          ENDIF
    57475748          IF (ok_plane_contrail) THEN
    5748              d_qx(i,k,icfa) = ( cfa_seri(i,k) - qx(i,k,icfa) ) / phys_tstep
    5749              d_qx(i,k,ipcf) = ( pcf_seri(i,k) - qx(i,k,ipcf) ) / phys_tstep
    5750              d_qx(i,k,iqva) = ( qva_seri(i,k) - qx(i,k,iqva) ) / phys_tstep
    5751              d_qx(i,k,iqia) = ( qia_seri(i,k) - qx(i,k,iqia) ) / phys_tstep
     5749             d_qx(i,k,icfl) = ( cfl_seri(i,k) - qx(i,k,icfl) ) / phys_tstep
     5750             d_qx(i,k,icfc) = ( cfc_seri(i,k) - qx(i,k,icfc) ) / phys_tstep
     5751             d_qx(i,k,iqtl) = ( qtl_seri(i,k) - qx(i,k,iqtl) ) / phys_tstep
     5752             d_qx(i,k,iqtc) = ( qtc_seri(i,k) - qx(i,k,iqtc) ) / phys_tstep
    57525753          ENDIF
    57535754       ENDDO
    … …  
    57815782    cf_ancien(:,:) = cf_seri(:,:)
    57825783    qvc_ancien(:,:)= qvc_seri(:,:)
    5783     cfa_ancien(:,:)= cfa_seri(:,:)
    5784     pcf_ancien(:,:)= pcf_seri(:,:)
    5785     qva_ancien(:,:)= qva_seri(:,:)
    5786     qia_ancien(:,:)= qia_seri(:,:)
     5784    cfl_ancien(:,:)= cfl_seri(:,:)
     5785    cfc_ancien(:,:)= cfc_seri(:,:)
     5786    qtl_ancien(:,:)= qtl_seri(:,:)
     5787    qtc_ancien(:,:)= qtc_seri(:,:)
    57875788    CALL water_int(klon,klev,q_ancien,zmasse,prw_ancien)
    57885789    CALL water_int(klon,klev,ql_ancien,zmasse,prlw_ancien)
Note: See TracChangeset for help on using the changeset viewer.