Ignore:
Timestamp:
Sep 20, 2024, 12:32:04 PM (11 hours ago)
Author:
Laurent Fairhead
Message:

Updating cirrus branch to trunk revision 5171

Location:
LMDZ6/branches/cirrus
Files:
3 deleted
12 edited
14 copied

Legend:

Unmodified
Added
Removed
  • LMDZ6/branches/cirrus

  • LMDZ6/branches/cirrus/libf/phylmdiso/add_phys_tend_mod.F90

    r4523 r5202  
    957957      bilh_bnd = (-(rcw-rcpd)*t_seri(1,1) + rlvtt) * rain_lsc(1) &
    958958    &         + (-(rcs-rcpd)*t_seri(1,1) + rlstt) * snow_lsc(1)
    959   CASE("bs") param
     959  CASE("bsss") param
    960960      bilq_bnd = - bs_fall(1)
    961961      bilh_bnd = (-(rcs-rcpd)*t_seri(1,1) + rlstt) * bs_fall(1)
  • LMDZ6/branches/cirrus/libf/phylmdiso/add_wake_tend.F90

    r4143 r5202  
    1 SUBROUTINE add_wake_tend(zddeltat, zddeltaq, zds, zddensaw, zddensw, zoccur, text, abortphy &
     1SUBROUTINE add_wake_tend(zddeltat, zddeltaq, zds, zdas, zddensw, zddensaw, zoccur, text, abortphy &
    22#ifdef ISO
    33        , zddeltaxt &
     
    1313
    1414USE dimphy, ONLY: klon, klev
    15 USE phys_state_var_mod, ONLY: wake_deltat, wake_deltaq, wake_s, &
    16                               awake_dens, wake_dens
     15USE phys_state_var_mod, ONLY: wake_deltat, wake_deltaq, wake_s, awake_s, &
     16                              wake_dens, awake_dens
    1717
    1818USE print_control_mod, ONLY: prt_level
     
    2626!------------
    2727  REAL, DIMENSION(klon, klev),   INTENT (IN)         :: zddeltat, zddeltaq
    28   REAL, DIMENSION(klon),         INTENT (IN)         :: zds, zddensaw, zddensw
     28  REAL, DIMENSION(klon),         INTENT (IN)         :: zds, zdas, zddensw, zddensaw
    2929  INTEGER, DIMENSION(klon),      INTENT (IN)         :: zoccur
    3030  CHARACTER*(*),                 INTENT (IN)         :: text
     
    7979           IF (zoccur(i) .GE. 1) THEN
    8080             wake_s(i)     = wake_s(i)    + zds(i)
     81             awake_s(i)    = awake_s(i)    + zdas(i)
     82             wake_dens(i)  = wake_dens(i) + zddensw(i)
    8183             awake_dens(i) = awake_dens(i) + zddensaw(i)
    82              wake_dens(i)  = wake_dens(i) + zddensw(i)
    8384           ELSE
    8485             wake_s(i)     = 0.
     86             awake_s(i)    = 0.
     87             wake_dens(i)  = 0.
    8588             awake_dens(i) = 0.
    86              wake_dens(i)  = 0.
    8789           ENDIF   ! (zoccur(i) .GE. 1)
    8890         END DO
  • LMDZ6/branches/cirrus/libf/phylmdiso/isotopes_mod.F90

    r4491 r5202  
    2020
    2121   !--- Variables not depending on isotopes
    22    REAL,    SAVE :: pxtmelt, pxtice, pxtmin, pxtmax
    23 !$OMP THREADPRIVATE(pxtmelt, pxtice, pxtmin, pxtmax)
    24    REAL,    SAVE :: tdifexp, tv0cin, thumxt1
    25 !$OMP THREADPRIVATE(tdifexp, tv0cin, thumxt1)
     22   REAL,    SAVE :: thumxt1
     23   !$OMP THREADPRIVATE(thumxt1)
    2624   INTEGER, SAVE :: ntot
    2725!$OMP THREADPRIVATE(ntot)
     
    3028   REAL,    SAVE :: P_veg
    3129!$OMP THREADPRIVATE(P_veg)
    32    REAL,    SAVE :: musi, lambda_sursat
    33 !$OMP THREADPRIVATE(musi, lambda_sursat)
    34    REAL,    SAVE :: Kd
    35 !$OMP THREADPRIVATE(Kd)
    36    REAL,    SAVE :: rh_cste_surf_cond, T_cste_surf_cond
    37 !$OMP THREADPRIVATE(rh_cste_surf_cond, T_cste_surf_cond)
     30   REAL,    SAVE :: lambda_sursat
     31!$OMP THREADPRIVATE(lambda_sursat)
    3832   LOGICAL, SAVE :: bidouille_anti_divergence    ! T: regularly, xteau <- q to avoid slow drifts
    3933!$OMP THREADPRIVATE(bidouille_anti_divergence)
     
    5448   REAL,    SAVE :: sstlatcrit, dsstlatcrit
    5549!$OMP THREADPRIVATE(sstlatcrit, dsstlatcrit)
    56    REAL,    SAVE :: deltaO18_oce
    57 !$OMP THREADPRIVATE(deltaO18_oce)
    5850   INTEGER, SAVE :: albedo_prescrit              ! 0: default ; 1: constant albedo
    5951!$OMP THREADPRIVATE(albedo_prescrit)
     
    8880   REAL,    SAVE :: fac_modif_evaoce
    8981!$OMP THREADPRIVATE(fac_modif_evaoce)
     82   REAL,    SAVE :: deltaO18_oce
     83!$OMP THREADPRIVATE(deltaO18_oce)
    9084   INTEGER, SAVE :: ok_bidouille_wake
    9185!$OMP THREADPRIVATE(ok_bidouille_wake)
     
    106100                    alpha_liq_sol, Rdefault, Rmethox
    107101!$OMP THREADPRIVATE(alpha_liq_sol, Rdefault, Rmethox)
    108    REAL, SAVE ::    fac_coeff_eq17_liq, fac_coeff_eq17_ice
    109 !$OMP THREADPRIVATE(fac_coeff_eq17_liq, fac_coeff_eq17_ice)
     102!   REAL, SAVE ::    fac_coeff_eq17_liq, fac_coeff_eq17_ice
     103!!$OMP THREADPRIVATE(fac_coeff_eq17_liq, fac_coeff_eq17_ice)
     104
     105   !--- H2[18]O reference
     106   REAL, PARAMETER :: fac_enrichoce18=0.0005
     107   REAL, PARAMETER :: alpha_liq_sol_O18=1.00291
     108   REAL, PARAMETER :: talph1_O18=1137.
     109   REAL, PARAMETER :: talph2_O18=-0.4156
     110   REAL, PARAMETER :: talph3_O18=-2.0667E-3
     111   REAL, PARAMETER :: talps1_O18=11.839
     112   REAL, PARAMETER :: talps2_O18=-0.028244
     113   REAL, PARAMETER :: tdifrel_O18=1./0.9723
     114   REAL, PARAMETER :: tkcin0_O18=0.006
     115   REAL, PARAMETER :: tkcin1_O18=0.000285
     116   REAL, PARAMETER :: tkcin2_O18=0.00082
     117   REAL, PARAMETER :: fac_coeff_eq17_liq=0.529
     118   REAL, PARAMETER :: fac_coeff_eq17_ice=0.529
     119
     120   !---- Parameters that do not depend on the nature of water isotopes:
     121   REAL, PARAMETER :: pxtmelt = 273.15 ! temperature at which ice formation starts
     122   REAL, PARAMETER :: pxtice  = 273.15-10.0 ! -- temperature at which all condensate is ice:
     123   REAL, PARAMETER :: pxtmin = 273.15 - 120.0   ! On ne calcule qu'au dessus de -120°C
     124   REAL, PARAMETER :: pxtmax = 273.15 +  60.0   ! On ne calcule qu'au dessus de +60°C
     125   REAL, PARAMETER :: tdifexp = 0.58 ! -- a constant for alpha_eff for equilibrium below cloud base:
     126   REAL, PARAMETER :: tv0cin  = 7.0 ! wind threshold (m/s) for smooth/rough regime (Merlivat and Jouzel 1979)
     127   REAL, PARAMETER :: musi=1.0  ! facteurs lambda et mu dans Si=musi-lambda*T
     128   REAL, PARAMETER :: Kd=2.5e-9 ! m2/s ! diffusion dans le sol
     129   REAL, PARAMETER :: rh_cste_surf_cond = 0.6 ! cas où cste_surf_cond: on met rhs ou/et Ts cste pour voir
     130   REAL, PARAMETER :: T_cste_surf_cond = 288.0
     131
    110132
    111133   !--- Negligible lower thresholds: no need to check for absurd values under these lower limits
     
    140162   INTEGER :: ixt
    141163
    142    !--- H2[18]O reference
    143    REAL :: fac_enrichoce18, alpha_liq_sol_O18, &
    144            talph1_O18, talph2_O18, talph3_O18, talps1_O18, talps2_O18, &
    145            tkcin0_O18, tkcin1_O18, tkcin2_O18, tdifrel_O18 
    146 
     164 
    147165   !--- For H2[17]O
    148166   REAL    :: fac_kcin, pente_MWL
     
    152170   LOGICAL, PARAMETER ::   ok_nocinsat = .FALSE. ! if T: no sursaturation effect for ice
    153171   LOGICAL, PARAMETER :: Rdefault_smow = .FALSE. ! if T: Rdefault=smow; if F: nul
     172   LOGICAL, PARAMETER :: tnat1 = .TRUE. ! If T: all tnats are 1.
    154173
    155174   !--- For [3]H
     
    157176
    158177   CHARACTER(LEN=maxlen) :: modname, sxt
    159    REAL, ALLOCATABLE :: tmp(:)
    160178
    161179   modname = 'iso_init'
     
    214232      CALL get_in('lat_max_albedo', lat_max_albedo,  100.)
    215233   END IF
    216    deltaO18_oce=0.0
    217234   CALL get_in('deltaP_BL',           deltaP_BL,     10.0)
    218235   CALL get_in('ruissellement_pluie', ruissellement_pluie, 0)
     
    249266   CALL get_in('ok_prod_nucl_tritium', ok_prod_nucl_tritium, .FALSE., .FALSE.)
    250267
    251    !--------------------------------------------------------------
    252    ! Parameters that do not depend on the nature of water isotopes:
    253    !--------------------------------------------------------------
    254    ! -- temperature at which ice condensate starts to form (valeur ECHAM?):
    255    pxtmelt = 273.15
    256 
    257    ! -- temperature at which all condensate is ice:
    258    pxtice  = 273.15-10.0
    259 
    260    !- -- test PHASE
    261 !   pxtmelt = 273.15 - 10.0
    262 !   pxtice  = 273.15 - 30.0
    263 
    264    ! -- minimum temperature to calculate fractionation coeff
    265    pxtmin = 273.15 - 120.0   ! On ne calcule qu'au dessus de -120°C
    266    pxtmax = 273.15 +  60.0   ! On ne calcule qu'au dessus de +60°C
    267    !    Remarque: les coeffs ont ete mesures seulement jusq'à -40!
    268 
    269    ! -- a constant for alpha_eff for equilibrium below cloud base:
    270    tdifexp = 0.58
    271    tv0cin  = 7.0
    272 
    273    ! facteurs lambda et mu dans Si=musi-lambda*T
    274    musi=1.0
    275    if (ok_nocinsat) lambda_sursat = 0.0          ! no sursaturation effect
    276 
    277    ! diffusion dans le sol
    278    Kd=2.5e-9 ! m2/s   
    279 
    280    ! cas où cste_surf_cond: on met rhs ou/et Ts cste pour voir
    281    rh_cste_surf_cond = 0.6
    282     T_cste_surf_cond = 288.0
     268   ! Ocean composition
     269   CALL get_in('deltaO18_oce',  deltaO18_oce, 0.0)
    283270   
    284271   CALL msg('iso_O18, iso_HDO, iso_eau = '//TRIM(strStack(int2str([iso_O18, iso_HDO, iso_eau]))), modname)
    285272
    286273   !--------------------------------------------------------------
    287    ! Parameters that depend on the nature of water isotopes:
     274   ! Isotope fractionation factors and a few isotopic constants
    288275   !--------------------------------------------------------------
    289    IF(getKey('tnat',    tnat,    isoName)) CALL abort_physic(modname, 'can''t get tnat',    1)
    290    IF(getKey('toce',    toce,    isoName)) CALL abort_physic(modname, 'can''t get toce',    1)
    291    IF(getKey('tcorr',   tcorr,   isoName)) CALL abort_physic(modname, 'can''t get tcorr',   1)
    292    IF(getKey('talph1',  talph1,  isoName)) CALL abort_physic(modname, 'can''t get talph1',  1)
    293    IF(getKey('talph2',  talph2,  isoName)) CALL abort_physic(modname, 'can''t get talph2',  1)
    294    IF(getKey('talph3',  talph3,  isoName)) CALL abort_physic(modname, 'can''t get talph3',  1)
    295    IF(getKey('talps1',  talps1,  isoName)) CALL abort_physic(modname, 'can''t get talps1',  1)
    296    IF(getKey('talps2',  talps2,  isoName)) CALL abort_physic(modname, 'can''t get talps2',  1)
    297    IF(getKey('tkcin0',  tkcin0,  isoName)) CALL abort_physic(modname, 'can''t get tkcin0',  1)
    298    IF(getKey('tkcin1',  tkcin1,  isoName)) CALL abort_physic(modname, 'can''t get tkcin1',  1)
    299    IF(getKey('tkcin2',  tkcin2,  isoName)) CALL abort_physic(modname, 'can''t get tkcin2',  1)
    300    IF(getKey('tdifrel', tdifrel, isoName)) CALL abort_physic(modname, 'can''t get tdifrel', 1)
    301    IF(getKey('alpha_liq_sol', alpha_liq_sol, isoName)) CALL abort_physic(modname, 'can''t get alpha_liq_sol',  1)
    302    IF(getKey('Rdefault',Rdefault,isoName)) CALL abort_physic(modname, 'can''t get Rdefault',1)
    303    IF(getKey('Rmethox', Rmethox, isoName)) CALL abort_physic(modname, 'can''t get Rmethox', 1)
     276   ALLOCATE(tkcin0(niso))
     277   ALLOCATE(tkcin1(niso))
     278   ALLOCATE(tkcin2(niso))
     279   ALLOCATE(tnat(niso))
     280   ALLOCATE(tdifrel(niso))
     281   ALLOCATE(toce(niso))
     282   ALLOCATE(tcorr(niso))
     283   ALLOCATE(talph1(niso))
     284   ALLOCATE(talph2(niso))
     285   ALLOCATE(talph3(niso))
     286   ALLOCATE(talps1(niso))
     287   ALLOCATE(talps2(niso))
     288   ALLOCATE(alpha_liq_sol(niso))
     289   ALLOCATE(Rdefault(niso))
     290   ALLOCATE(Rmethox(niso))
     291
     292   do ixt=1,niso
     293     if (ixt.eq.iso_HTO) then  ! Tritium
     294       tkcin0(ixt) = 0.01056
     295       tkcin1(ixt) = 0.0005016
     296       tkcin2(ixt) = 0.0014432
     297       if (tnat1) then
     298               tnat(ixt)=1
     299       else
     300               tnat(ixt)=0.
     301       endif
     302       toce(ixt)=4.0E-19 ! rapport T/H = 0.2 TU Dreisigacker and Roether 1978
     303       tcorr(ixt)=1.
     304       tdifrel(ixt)=1./0.968
     305       talph1(ixt)=46480.
     306       talph2(ixt)=-103.87
     307       talph3(ixt)=0.
     308       talps1(ixt)=46480.
     309       talps2(ixt)=-103.87
     310       alpha_liq_sol(ixt)=1.
     311       Rmethox(ixt)=0.0
     312     endif
     313     if (ixt.eq.iso_O17) then  ! O17
     314       pente_MWL=0.528
     315       tdifrel(ixt)=1./0.98555 ! valeur utilisée en 1D et dans modèle de LdG ! tdifrel(ixt)=1./0.985452 ! donné par Amaelle
     316       fac_kcin= (tdifrel(ixt)-1.0)/(tdifrel_O18-1.0) ! fac_kcin=0.5145 ! donné par Amaelle
     317       tkcin0(ixt) = tkcin0_O18*fac_kcin
     318       tkcin1(ixt) = tkcin1_O18*fac_kcin
     319       tkcin2(ixt) = tkcin2_O18*fac_kcin
     320       if (tnat1) then
     321               tnat(ixt)=1
     322       else
     323               tnat(ixt)=0.004/100. ! O17 représente 0.004% de l'oxygène
     324       endif
     325       toce(ixt)=tnat(ixt)*(1.0+deltaO18_oce/1000.0)**pente_MWL
     326       tcorr(ixt)=1.0+fac_enrichoce18*pente_MWL ! donné par Amaelle           
     327       talph1(ixt)=talph1_O18
     328       talph2(ixt)=talph2_O18
     329       talph3(ixt)=talph3_O18
     330       talps1(ixt)=talps1_O18
     331       talps2(ixt)=talps2_O18     
     332       alpha_liq_sol(ixt)=(alpha_liq_sol_O18)**fac_coeff_eq17_liq
     333       Rdefault(ixt)=tnat(ixt)*(-3.15/1000.0+1.0)
     334       Rmethox(ixt)=(230./1000.+1.)*tnat(ixt) !Zahn et al 2006
     335     endif
     336     if (ixt.eq.iso_O18) then  ! Oxygene18
     337       tkcin0(ixt) = tkcin0_O18
     338       tkcin1(ixt) = tkcin1_O18
     339       tkcin2(ixt) = tkcin2_O18
     340       if (tnat1) then
     341               tnat(ixt)=1
     342       else
     343               tnat(ixt)=2005.2E-6
     344       endif
     345       toce(ixt)=tnat(ixt)*(1.0+deltaO18_oce/1000.0)
     346       tcorr(ixt)=1.0+fac_enrichoce18
     347       tdifrel(ixt)=tdifrel_O18
     348       talph1(ixt)=talph1_O18
     349       talph2(ixt)=talph2_O18
     350       talph3(ixt)=talph3_O18
     351       talps1(ixt)=talps1_O18
     352       talps2(ixt)=talps2_O18
     353       alpha_liq_sol(ixt)=alpha_liq_sol_O18   
     354       Rdefault(ixt)=tnat(ixt)*(-6.0/1000.0+1.0)
     355       Rmethox(ixt)=(130./1000.+1.)*tnat(ixt) !Zahn et al 2006 
     356     endif
     357     if (ixt.eq.iso_HDO) then ! Deuterium
     358       pente_MWL=8.0
     359       tdifrel(ixt)=1./0.9755 !          fac_kcin=0.88
     360       fac_kcin= (tdifrel(ixt)-1)/(tdifrel_O18-1)
     361       tkcin0(ixt) = tkcin0_O18*fac_kcin
     362       tkcin1(ixt) = tkcin1_O18*fac_kcin
     363       tkcin2(ixt) = tkcin2_O18*fac_kcin
     364       if (tnat1) then
     365               tnat(ixt)=1
     366       else
     367               tnat(ixt)=155.76E-6
     368       endif
     369       toce(ixt)=tnat(ixt)*(1.0+pente_MWL*deltaO18_oce/1000.0)
     370       tcorr(ixt)=1.0+fac_enrichoce18*pente_MWL         
     371       talph1(ixt)=24844.
     372       talph2(ixt)=-76.248
     373       talph3(ixt)=52.612E-3
     374       talps1(ixt)=16288.
     375       talps2(ixt)=-0.0934
     376       !ZXalpha_liq_sol=1.0192 ! Weston, Ralph, 1955
     377       alpha_liq_sol(ixt)=1.0212
     378       ! valeur de Lehmann & Siegenthaler, 1991, Journal of
     379       ! Glaciology, vol 37, p 23
     380       Rdefault(ixt)=tnat(ixt)*((-6.0*pente_MWL+10.0)/1000.0+1.0)
     381       Rmethox(ixt)=tnat(ixt)*(-25.0/1000.+1.) ! Zahn et al 2006
     382     endif
     383     if (ixt.eq.iso_eau) then ! Oxygene16
     384       tkcin0(ixt) = 0.0
     385       tkcin1(ixt) = 0.0
     386       tkcin2(ixt) = 0.0
     387       tnat(ixt)=1.
     388       toce(ixt)=tnat(ixt)
     389       tcorr(ixt)=1.0
     390       tdifrel(ixt)=1.
     391       talph1(ixt)=0.
     392       talph2(ixt)=0.
     393       talph3(ixt)=0.
     394       talps1(ixt)=0.
     395       talph3(ixt)=0.
     396       alpha_liq_sol(ixt)=1.
     397       Rdefault(ixt)=tnat(ixt)*1.0
     398       Rmethox(ixt)=1.0
     399     endif
     400   enddo ! ixt=1,niso
    304401
    305402   IF(.NOT.Rdefault_smow) then
     
    308405   ENDIF
    309406   write(*,*) 'Rdefault=',Rdefault
     407   write(*,*) 'toce=',toce
    310408
    311409   !--- Sensitivity test: no kinetic effect in sfc evaporation
  • LMDZ6/branches/cirrus/libf/phylmdiso/isotopes_routines_mod.F90

    r4491 r5202  
    15251525#endif       
    15261526       pxtfra=max(min(pxtfra,alpha_max),0.0)
    1527 
    15281527
    15291528      end subroutine fractcalk_liq
     
    1592215921
    1592315922      ! verif
    15924 !      text="phyisoetat0 67"
    15925 !      write(*,*) 'snow(8,1)=',snow(8,1)
    15926 !      write(*,*) 'xtsnow(4,8,1)=',xtsnow(4,8,1)
    1592715923#ifdef ISOVERIF
    1592815924      do i=1,klon
     
    1593415930         enddo !do ixt=1,niso
    1593515931      enddo !do i=1,klon
    15936 #endif     
    15937 #ifdef ISOVERIF
    1593815932      do i=1,klon
    1593915933         if (iso_eau.gt.0) then
     
    1602116015         endif
    1602216016       enddo !do i=1,klon
    16023 
    1602416017#endif
    1602516018      !end verif
     
    1612816121          deltaD_run_off_lic_0(ixt)=deltaD_sol(ixt)
    1612916122          deltaD_land_ice(ixt)=deltaD_snow(ixt)
    16130           call fractcalk_liq(ixt, 283.0, alpha(ixt))           
     16123          call fractcalk_liq(ixt, 283.0, alpha(ixt))   
    1613116124        enddo !do ixt=1,niso
    1613216125        call calcul_kcin(2.0,kcin)
     
    1883018823        if ((iso_HDO.gt.0).and.(ixt.eq.iso_HDO)) then
    1883118824            if (q.gt.ridicule) then
     18825                    write(*,*) 'xt,q=',xt,q
     18826                    write(*,*) 'alpha=',alpha
     18827                    write(*,*) 'toce,kcin,h0=',toce,kcin,h0
     18828                    write(*,*) 'RMerlivat=',RMerlivat
    1883218829                call iso_verif_aberrant_encadre( xt/q, 'isotopes_routines_mod 18930b: iso_init_ideal')
    1883318830            endif
     
    1890218899end subroutine appel_stewart_debug
    1890318900
     18901
     18902subroutine dispatch(klon,klev,qx,q_seri,xt_seri,ql_seri,xtl_seri,qs_seri,xts_seri)
     18903
     18904use infotrac_phy, ONLY: nqtot,nqo,ivap,iliq,isol,iqIsoPha,ntraciso=>ntiso
     18905implicit none
     18906
     18907! inputs
     18908integer, intent(in) :: klon,klev
     18909real,dimension(klon,klev,nqtot), intent(in) ::qx
     18910
     18911! outputs
     18912real,dimension(klon,klev), intent(out) ::q_seri,ql_seri,qs_seri
     18913real,dimension(ntraciso,klon,klev), intent(out) :: xt_seri,xtl_seri,xts_seri
     18914
     18915! locals
     18916integer :: i,k,ixt
     18917
     18918do k=1,klev
     18919do i=1,klon
     18920    q_seri(i,k)  = qx(i,k,ivap)
     18921    ql_seri(i,k) = qx(i,k,iliq)
     18922    IF (nqo.EQ.2) THEN             !--vapour and liquid only
     18923             qs_seri(i,k) = 0.
     18924    ELSE IF (nqo.ge.3) THEN        !--vapour, liquid and ice
     18925             qs_seri(i,k) = qx(i,k,isol)
     18926    ENDIF
     18927    do ixt=1,ntraciso
     18928          xt_seri(ixt,i,k)  = qx(i,k,iqIsoPha(ixt,ivap))
     18929          xtl_seri(ixt,i,k) = qx(i,k,iqIsoPha(ixt,iliq))
     18930          if (nqo.eq.2) then
     18931             xts_seri(ixt,i,k) = 0.
     18932          else if (nqo.eq.3) then
     18933             xts_seri(ixt,i,k) = qx(i,k,iqIsoPha(ixt,isol))
     18934          endif
     18935    enddo !do ixt=1,niso
     18936
     18937enddo
     18938enddo
     18939
     18940end subroutine dispatch
     18941
     18942subroutine together(klon,klev,qx,q_seri,xt_seri,ql_seri,xtl_seri,qs_seri,xts_seri)
     18943
     18944use infotrac_phy, ONLY: nqtot,nqo,ivap,iliq,isol,iqIsoPha,ntraciso=>ntiso
     18945implicit none
     18946
     18947! inputs
     18948integer, intent(in) :: klon,klev
     18949real,dimension(klon,klev), intent(in) ::q_seri,ql_seri,qs_seri
     18950real,dimension(ntraciso,klon,klev), intent(in) :: xt_seri,xtl_seri,xts_seri
     18951
     18952! inputs
     18953real,dimension(klon,klev,nqtot), intent(out) ::qx
     18954
     18955! locals
     18956integer :: i,k,ixt
     18957
     18958do k=1,klev
     18959do i=1,klon 
     18960    qx(i,k,ivap)  = q_seri(i,k)
     18961    qx(i,k,iliq) = ql_seri(i,k)
     18962    IF (nqo.ge.3) THEN        !--vapour, liquid and ice
     18963             qx(i,k,isol) = qs_seri(i,k)
     18964    ENDIF
     18965    do ixt=1,ntraciso
     18966          qx(i,k,iqIsoPha(ixt,ivap)) = xt_seri(ixt,i,k) 
     18967          qx(i,k,iqIsoPha(ixt,iliq)) = xtl_seri(ixt,i,k)
     18968          if (nqo.ge.3) then
     18969             qx(i,k,iqIsoPha(ixt,isol)) = xts_seri(ixt,i,k)
     18970          endif
     18971    enddo !do ixt=1,niso
     18972
     18973enddo
     18974enddo
     18975
     18976end subroutine together
     18977
     18978
    1890418979END MODULE isotopes_routines_mod
    1890518980#endif
  • LMDZ6/branches/cirrus/libf/phylmdiso/isotopes_verif_mod.F90

    r4491 r5202  
    10421042            write(*,*) 'deltaD=',deltaD
    10431043            write(*,*) 'Dexcess=',dexcess
     1044            write(*,*) 'tnat=',tnat
    10441045!            stop
    10451046            iso_verif_O18_aberrant_nostop=1
  • LMDZ6/branches/cirrus/libf/phylmdiso/isotrac_routines_mod.F90

    r4491 r5202  
    681681              Eqi_prime_cas(il)=Eqi_prime(cas(il)) &
    682682     &           *(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
    683               Eqi_cas(il)=Eqi(il) &
     683              Eqi_cas(il)=Eqi(cas(il)) & ! corr bug Camille 15 juin 2024
    684684     &           *(Pxtisup(ieau,cas(il))/Pqisup(cas(il)))
    685685            else
  • LMDZ6/branches/cirrus/libf/phylmdiso/phyetat0_mod.F90

    r5055 r5202  
    3030       solsw, solswfdiff, t_ancien, u_ancien, v_ancien, w01, wake_cstar, wake_deltaq, &
    3131       wake_deltat, wake_delta_pbl_TKE, delta_tsurf, beta_aridity, wake_fip, wake_pe, &
    32        wake_s, wake_dens, awake_dens, cv_gen, zgam, zmax0, zmea, zpic, zsig, &
     32       wake_s, awake_s, wake_dens, awake_dens, cv_gen, zgam, zmax0, zmea, zpic, zsig, &
    3333#ifdef ISO
    3434       fxtevap, xtsol, xt_ancien, xtl_ancien, xts_ancien, wake_deltaxt, &
     
    4949  USE time_phylmdz_mod, ONLY: init_iteration, pdtphys, itau_phy
    5050  USE wxios, ONLY: missing_val_xios => missing_val, using_xios
    51   use netcdf, only: nf90_fill_real
     51  use netcdf, only: missing_val_netcdf => nf90_fill_real
    5252  use config_ocean_skin_m, only: activate_ocean_skin
    5353#ifdef ISO
     
    112112  REAL Rland_ice(niso,klon)
    113113#endif
     114
     115  IF (using_xios) THEN
     116    missing_val=missing_val_xios
     117  ELSE
     118    missing_val=missing_val_netcdf
     119  ENDIF
     120
    114121  ! FH1D
    115122  !     real iolat(jjm+1)
     
    117124
    118125  ! Ouvrir le fichier contenant l'etat initial:
    119   IF (using_xios) THEN
    120     missing_val = missing_val_xios
    121   ELSE
    122     missing_val =  nf90_fill_real
    123   ENDIF
    124126
    125127  CALL open_startphy(fichnom)
     
    324326
    325327!===================================================================
     328! Lecture dans le cas iflag_pbl_surface =1
     329!===================================================================
     330
     331   if ( iflag_physiq <= 1 ) then
     332!===================================================================
    326333  ! Lecture des temperatures du sol profond:
    327334!===================================================================
     
    351358  found=phyetat0_get(snow_fall,"snow_f","snow fall",0.)
    352359  found=phyetat0_get(rain_fall,"rain_f","rain fall",0.)
    353 
    354360  IF (ok_bs) THEN
    355361     found=phyetat0_get(bs_fall,"bs_f","blowing snow fall",0.)
     
    405411  ENDIF
    406412
     413  endif ! iflag_physiq <= 1
     414
    407415  ! Lecture de l'age de la neige:
    408416  found=phyetat0_srf(agesno,"AGESNO","SNOW AGE",0.001)
     
    498506  found=phyetat0_get(wake_deltaq,"WAKE_DELTAQ","Delta hum. wake/env",0.)
    499507  found=phyetat0_get(wake_s,"WAKE_S","Wake frac. area",0.)
     508  found=phyetat0_get(awake_s,"AWAKE_S","Active Wake frac. area",0.)
    500509!jyg<
    501510!  Set wake_dens to -1000. when there is no restart so that the actual
     
    677686!        write(*,*) 'xtsnow(:,994,2)=',xtsnow(:,994,2)
    678687!#endif
    679 
     688  if ( iflag_physiq <= 1 ) then
    680689  CALL pbl_surface_init(fder, snow, qsurf, tsoil)
    681690#ifdef ISO
    682691  CALL pbl_surface_init_iso(xtsnow,Rland_ice)
    683692#endif
     693  endif
    684694
    685695  ! Initialize module ocean_cpl_mod for the case of coupled ocean
  • LMDZ6/branches/cirrus/libf/phylmdiso/phys_local_var_mod.F90

    r5055 r5202  
    1414      REAL, SAVE, ALLOCATABLE :: ql_seri(:,:),qs_seri(:,:)
    1515      !$OMP THREADPRIVATE(ql_seri,qs_seri)
     16! SN 15/07/2024 ISO 4D
     17      REAL, SAVE, ALLOCATABLE :: qx_seri(:,:,:)
     18      !$OMP THREADPRIVATE(qx_seri)
     19! SN
    1620      REAL, SAVE, ALLOCATABLE :: qbs_seri(:,:)
    1721      !$OMP THREADPRIVATE(qbs_seri)
     
    2428      REAL, SAVE, ALLOCATABLE :: pbl_eps(:,:,:)
    2529      !$OMP THREADPRIVATE(pbl_eps)
     30      REAL, SAVE, ALLOCATABLE :: tke_shear(:,:,:), tke_buoy(:,:,:), tke_trans(:,:,:)
     31      !$OMP THREADPRIVATE(tke_shear,tke_buoy,tke_trans)
    2632      REAL, SAVE, ALLOCATABLE :: tr_seri(:,:,:)
    2733      !$OMP THREADPRIVATE(tr_seri)
     
    6470      REAL, SAVE, ALLOCATABLE :: d_t_eva(:,:),d_q_eva(:,:),d_ql_eva(:,:),d_qi_eva(:,:)
    6571      !$OMP THREADPRIVATE(d_t_eva,d_q_eva,d_ql_eva,d_qi_eva)
     72! SN 15/07/2024 ISO 4D
     73      REAL, SAVE, ALLOCATABLE :: d_qx_eva(:,:,:)
     74      !$OMP THREADPRIVATE(d_qx_eva)
     75! SN
    6676      REAL, SAVE, ALLOCATABLE :: d_t_lscst(:,:),d_q_lscst(:,:)
    6777      !$OMP THREADPRIVATE(d_t_lscst,d_q_lscst)
     
    8494      REAL, SAVE, ALLOCATABLE :: d_t_vdf_x(:,:), d_q_vdf_x(:,:)
    8595      !$OMP THREADPRIVATE( d_t_vdf_x, d_q_vdf_x)
    86       REAL, SAVE, ALLOCATABLE :: d_t_bs(:,:), d_q_bs(:,:), d_qbs_bs(:,:)
    87       !$OMP THREADPRIVATE( d_t_bs,d_q_bs, d_qbs_bs)
     96      REAL, SAVE, ALLOCATABLE :: d_t_bsss(:,:), d_q_bsss(:,:), d_qbs_bsss(:,:)
     97      !$OMP THREADPRIVATE( d_t_bsss,d_q_bsss, d_qbs_bsss)
    8898!>nrlmd+jyg
    8999      REAL, SAVE, ALLOCATABLE :: d_t_oro(:,:)
     
    124134      REAL, SAVE, ALLOCATABLE :: xts_seri(:,:,:)
    125135      !$OMP THREADPRIVATE( xts_seri)
     136      REAL, SAVE, ALLOCATABLE :: xtbs_seri(:,:,:)
     137      !$OMP THREADPRIVATE( xtbs_seri)
    126138      REAL, SAVE, ALLOCATABLE :: d_xt_eva(:,:,:)
    127139      !$OMP THREADPRIVATE( d_xt_eva)
     
    134146      REAL, SAVE, ALLOCATABLE :: d_xt_dyn(:,:,:)
    135147      !$OMP THREADPRIVATE( d_xt_dyn)
    136       REAL, SAVE, ALLOCATABLE :: d_xtl_dyn(:,:,:), d_xts_dyn(:,:,:)
    137       !$OMP THREADPRIVATE(d_xtl_dyn, d_xts_dyn)
     148      REAL, SAVE, ALLOCATABLE :: d_xtl_dyn(:,:,:), d_xts_dyn(:,:,:), d_xtbs_dyn(:,:,:)
     149      !$OMP THREADPRIVATE(d_xtl_dyn, d_xts_dyn, d_xtbs_dyn)
    138150      REAL, SAVE, ALLOCATABLE :: d_xt_con(:,:,:)
    139151      !$OMP THREADPRIVATE( d_xt_con)
     
    166178      !$OMP THREADPRIVATE(d_ts, d_tr)
    167179
    168 ! aerosols
    169       REAL, SAVE, ALLOCATABLE :: m_allaer (:,:,:)
    170       !$OMP THREADPRIVATE(m_allaer)
    171180! diagnostique pour le rayonnement
    172181      REAL, SAVE, ALLOCATABLE :: topswad_aero(:),  solswad_aero(:)      ! diag
     
    292301!$OMP THREADPRIVATE(toplwad0_aerop, sollwad0_aerop)
    293302
     303!AI 08 2023 ajout pour Ecrad
     304      REAL,ALLOCATABLE,SAVE :: topswad_aero_s2(:), solswad_aero_s2(:)
     305!$OMP THREADPRIVATE(topswad_aero_s2, solswad_aero_s2)
     306      REAL,ALLOCATABLE,SAVE :: topswai_aero_s2(:), solswai_aero_s2(:)
     307!$OMP THREADPRIVATE(topswai_aero_s2, solswai_aero_s2)
     308      REAL,ALLOCATABLE,SAVE :: topswad0_aero_s2(:), solswad0_aero_s2(:)
     309!$OMP THREADPRIVATE(topswad0_aero_s2, solswad0_aero_s2)
     310      REAL,ALLOCATABLE,SAVE :: topsw_aero_s2(:,:), topsw0_aero_s2(:,:)
     311!$OMP THREADPRIVATE(topsw_aero_s2, topsw0_aero_s2)
     312      REAL,ALLOCATABLE,SAVE :: solsw_aero_s2(:,:), solsw0_aero_s2(:,:)
     313!$OMP THREADPRIVATE(solsw_aero_s2, solsw0_aero_s2)
     314      REAL,ALLOCATABLE,SAVE :: topswcf_aero_s2(:,:), solswcf_aero_s2(:,:)
     315!$OMP THREADPRIVATE(topswcf_aero_s2, solswcf_aero_s2)
     316! additional LW variables CK
     317      REAL,ALLOCATABLE,SAVE :: toplwad_aero_s2(:), sollwad_aero_s2(:)
     318!$OMP THREADPRIVATE(toplwad_aero_s2, sollwad_aero_s2)
     319      REAL,ALLOCATABLE,SAVE :: toplwai_aero_s2(:), sollwai_aero_s2(:)
     320!$OMP THREADPRIVATE(toplwai_aero_s2, sollwai_aero_s2)
     321      REAL,ALLOCATABLE,SAVE :: toplwad0_aero_s2(:), sollwad0_aero_s2(:)
     322!$OMP THREADPRIVATE(toplwad0_aero_s2, sollwad0_aero_s2)
     323
    294324!Ajout de celles n??cessaires au phys_output_write_mod
    295325      REAL, SAVE, ALLOCATABLE :: tal1(:), pal1(:), pab1(:), pab2(:)
     
    300330!$OMP THREADPRIVATE(sens, flwp, fiwp)
    301331!!
    302 !FC
     332!FC 
    303333      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: zxfluxt, zxfluxq
    304334!$OMP THREADPRIVATE(zxfluxt, zxfluxq)
     
    315345    REAL, SAVE, ALLOCATABLE,DIMENSION(:,:)          :: d_deltat_wk, d_deltaq_wk
    316346!$OMP THREADPRIVATE(d_deltat_wk, d_deltaq_wk)
    317       REAL,ALLOCATABLE,SAVE,DIMENSION(:)            :: d_s_wk, d_dens_a_wk, d_dens_wk
    318 !$OMP THREADPRIVATE(d_s_wk, d_dens_a_wk, d_dens_wk)
     347      REAL,ALLOCATABLE,SAVE,DIMENSION(:)            :: d_s_wk, d_s_a_wk, d_dens_wk, d_dens_a_wk
     348!$OMP THREADPRIVATE(d_s_wk, d_s_a_wk, d_dens_wk, d_dens_a_wk)
    319349    REAL, SAVE, ALLOCATABLE,DIMENSION(:,:)          :: d_deltat_wk_gw, d_deltaq_wk_gw
    320350!$OMP THREADPRIVATE(d_deltat_wk_gw, d_deltaq_wk_gw)
     
    328358!!!OMP THREADPRIVATE(d_s_the, d_dens_the)
    329359      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:)           :: d_deltat_ajs_cv, d_deltaq_ajs_cv
    330 !$OMP THREADPRIVATE(d_deltat_ajs_cv, d_deltaq_ajs_cv)                       
     360!$OMP THREADPRIVATE(d_deltat_ajs_cv, d_deltaq_ajs_cv)
    331361#ifdef ISO
    332362    REAL, SAVE, ALLOCATABLE,DIMENSION(:,:,:)          :: d_deltaxt_wk
     
    376406      REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: zxfluxlat, zxtsol, snow_lsc, zxfqfonte
    377407!$OMP THREADPRIVATE(zxfluxlat, zxtsol, snow_lsc, zxfqfonte)
    378       REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: zxrunofflic
    379 !$OMP THREADPRIVATE(zxrunofflic)
     408!SN runoffdiag
     409      REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: zxrunofflic, runoff_diag
     410!$OMP THREADPRIVATE(zxrunofflic, runoff_diag)
    380411      REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: zxqsurf, rain_lsc, rain_num
    381412!$OMP THREADPRIVATE(zxqsurf, rain_lsc, rain_num)
     
    383414      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: xtevap,xtprw
    384415!$OMP THREADPRIVATE(xtevap,xtprw)
    385       REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: h1_diag,runoff_diag
     416      REAL,ALLOCATABLE,SAVE,DIMENSION(:) :: h1_diag
    386417      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: xtrunoff_diag
    387 !$OMP THREADPRIVATE(h1_diag,runoff_diag,xtrunoff_diag)
     418!$OMP THREADPRIVATE(h1_diagv,xtrunoff_diag)
    388419      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: zxfxtcalving
    389420!$OMP THREADPRIVATE(zxfxtcalving)
     
    581612      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: pfraclr,pfracld
    582613!$OMP THREADPRIVATE(pfraclr,pfracld)
     614      REAL, SAVE, ALLOCATABLE :: cldfraliq(:,:)
     615!$OMP THREADPRIVATE(cldfraliq)
     616      REAL, SAVE, ALLOCATABLE ::mean_icefracturb(:,:)
     617!$OMP THREADPRIVATE(mean_icefracturb)
     618      REAL, SAVE, ALLOCATABLE :: sigma2_icefracturb(:,:)
     619!$OMP THREADPRIVATE(sigma2_icefracturb)
    583620
    584621! variables de sorties MM
     
    671708!
    672709! variables for stratospheric aerosol
     710      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: d_q_emiss
     711!$OMP THREADPRIVATE(d_q_emiss)
    673712      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: R2SO4
    674713!$OMP THREADPRIVATE(R2SO4)
     714      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: R2SO4B
     715!$OMP THREADPRIVATE(R2SO4B)
    675716      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: DENSO4
    676717!$OMP THREADPRIVATE(DENSO4)
     718      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:, :) :: DENSO4B
     719!$OMP THREADPRIVATE(DENSO4B)     
    677720      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: f_r_wet
    678721!$OMP THREADPRIVATE(f_r_wet)
     722      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:, :) :: f_r_wetB
     723!$OMP THREADPRIVATE(f_r_wetB)
    679724      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: decfluxaer
    680725!$OMP THREADPRIVATE(decfluxaer)
     
    685730      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: SO2_lifetime
    686731!$OMP THREADPRIVATE(SO2_lifetime)
     732      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: H2SO4_lifetime
     733!$OMP THREADPRIVATE(H2SO4_lifetime)
     734      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: O3_clim
     735!$OMP THREADPRIVATE(O3_clim)
    687736      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: alpha_bin
    688737!$OMP THREADPRIVATE(alpha_bin)
     
    701750      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: vsed_aer
    702751!$OMP THREADPRIVATE(vsed_aer)
     752!     Sulfate aerosol concentration (dry mixing ratio) (condensed H2SO4 mmr)
     753      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sulfmmr
     754!$OMP THREADPRIVATE(sulfmmr)
     755!     SAD all aerosols (cm2/cm3)
     756      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: SAD_sulfate
     757!$OMP THREADPRIVATE(SAD_sulfate)
     758!     Effective radius of wet surface aerosols (cm)
     759      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:) :: reff_sulfate
     760!$OMP THREADPRIVATE(reff_sulfate)
     761!     sulfate MMR in different modes (based on sulfmmr, it must be dry mmr)
     762      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: sulfmmr_mode
     763!$OMP THREADPRIVATE(sulfmmr_mode)
     764!     particle concentration in different modes (part/m3)
     765      REAL, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: nd_mode
     766!$OMP THREADPRIVATE(nd_mode)
    703767!
    704768!---3D budget variables
     
    749813SUBROUTINE phys_local_var_init
    750814USE dimphy
    751 USE infotrac_phy, ONLY : nbtr
     815USE infotrac_phy, ONLY : nbtr,nqtot
    752816#ifdef ISO
    753817USE infotrac_phy, ONLY : ntraciso=>ntiso,niso
     
    757821USE phys_output_var_mod
    758822USE phys_state_var_mod
     823#ifdef CPP_StratAer
     824USE infotrac_phy, ONLY : nbtr_bin
     825#endif
    759826
    760827IMPLICIT NONE
    761828      ALLOCATE(t_seri(klon,klev),q_seri(klon,klev),ql_seri(klon,klev),qs_seri(klon,klev), qbs_seri(klon,klev))
     829! SN 4D ISO
     830      ALLOCATE(qx_seri(klon,klev,nqtot))
     831! SN
    762832      ALLOCATE(u_seri(klon,klev),v_seri(klon,klev))
    763833      ALLOCATE(cf_seri(klon,klev),rvc_seri(klon,klev))
    764834      ALLOCATE(l_mixmin(klon,klev+1,nbsrf),l_mix(klon,klev+1,nbsrf),wprime(klon,klev+1,nbsrf))
    765835      ALLOCATE(pbl_eps(klon,klev+1,nbsrf+1))
     836      ALLOCATE(tke_shear(klon,klev+1,nbsrf), tke_buoy(klon,klev+1,nbsrf), tke_trans(klon,klev+1,nbsrf))
    766837      pbl_eps(:,:,:)=0.
     838      tke_shear(:,:,:)=0.; tke_buoy(:,:,:)=0.; tke_trans(:,:,:)=0.
    767839      l_mix(:,:,:)=0.;l_mixmin(:,:,:)=0.;wprime(:,:,:)=0. ! doit etre initialse car pas toujours remplis
    768840      ALLOCATE(rhcl(klon,klev))
     
    789861      ALLOCATE(d_u_ajs(klon,klev),d_v_ajs(klon,klev))
    790862      ALLOCATE(d_t_eva(klon,klev),d_q_eva(klon,klev))
     863! SN 4D ISO
     864      ALLOCATE(d_qx_eva(klon,klev,nqtot))
     865! SN
    791866      ALLOCATE(d_ql_eva(klon,klev),d_qi_eva(klon,klev))
    792867      ALLOCATE(d_t_lscst(klon,klev),d_q_lscst(klon,klev))
     
    795870      ALLOCATE(d_t_vdf(klon,klev),d_q_vdf(klon,klev),d_t_diss(klon,klev))
    796871      ALLOCATE (d_qbs_vdf(klon,klev))
    797       ALLOCATE(d_t_bs(klon,klev),d_q_bs(klon,klev),d_qbs_bs(klon,klev))
     872      ALLOCATE(d_t_bsss(klon,klev),d_q_bsss(klon,klev),d_qbs_bsss(klon,klev))
    798873      ALLOCATE(d_t_vdf_w(klon,klev),d_q_vdf_w(klon,klev))
    799874      ALLOCATE(d_t_vdf_x(klon,klev),d_q_vdf_x(klon,klev))
     
    802877      allocate(xtl_seri(ntraciso,klon,klev))
    803878      allocate(xts_seri(ntraciso,klon,klev))
     879      allocate(xtbs_seri(ntraciso,klon,klev))
    804880      allocate(d_xt_dyn(ntraciso,klon,klev))
    805881      allocate(d_xtl_dyn(ntraciso,klon,klev))
    806882      allocate(d_xts_dyn(ntraciso,klon,klev))
     883      allocate(d_xtbs_dyn(ntraciso,klon,klev))
    807884      allocate(d_xt_con(ntraciso,klon,klev))
    808885      allocate(d_xt_wake(ntraciso,klon,klev))
     
    835912      ALLOCATE(d_u_lif(klon,klev),d_v_lif(klon,klev))
    836913      ALLOCATE(d_ts(klon,nbsrf), d_tr(klon,klev,nbtr))
     914
    837915! Special RRTM
    838916      ALLOCATE(ZLWFT0_i(klon,klev+1),ZSWFT0_i(klon,klev+1),ZFLDN0(klon,klev+1))
     
    913991      ALLOCATE(toplwad0_aerop(klon), sollwad0_aerop(klon))
    914992
     993!AI Ajout Ecrad (3Deffect)
     994      ALLOCATE(topswad_aero_s2(klon), solswad_aero_s2(klon))
     995      ALLOCATE(topswai_aero_s2(klon), solswai_aero_s2(klon))
     996      ALLOCATE(topswad0_aero_s2(klon), solswad0_aero_s2(klon))
     997      ALLOCATE(topsw_aero_s2(klon,naero_grp), topsw0_aero_s2(klon,naero_grp))
     998      ALLOCATE(solsw_aero_s2(klon,naero_grp), solsw0_aero_s2(klon,naero_grp))
     999      ALLOCATE(topswcf_aero_s2(klon,naero_grp), solswcf_aero_s2(klon,naero_grp))
     1000! additional LW variables CK
     1001      ALLOCATE(toplwad_aero_s2(klon), sollwad_aero_s2(klon))
     1002      ALLOCATE(toplwai_aero_s2(klon), sollwai_aero_s2(klon))
     1003      ALLOCATE(toplwad0_aero_s2(klon), sollwad0_aero_s2(klon))
     1004
    9151005! FH Ajout de celles necessaires au phys_output_write_mod
    9161006
     
    9231013      ALLOCATE(wake_omg(klon, klev))
    9241014      ALLOCATE(d_deltat_wk(klon, klev), d_deltaq_wk(klon, klev))
    925       ALLOCATE(d_s_wk(klon), d_dens_a_wk(klon), d_dens_wk(klon))
     1015      ALLOCATE(d_s_wk(klon), d_s_a_wk(klon), d_dens_wk(klon), d_dens_a_wk(klon))
    9261016      ALLOCATE(d_deltat_wk_gw(klon, klev), d_deltaq_wk_gw(klon, klev))
    9271017      ALLOCATE(d_deltat_vdf(klon, klev), d_deltaq_vdf(klon, klev))
     
    9581048      ALLOCATE(zxfqcalving(klon), zxfluxlat(klon))
    9591049      ALLOCATE(zxtsol(klon), snow_lsc(klon), zxfqfonte(klon), zxqsurf(klon))
    960       ALLOCATE(zxrunofflic(klon))
     1050! SN add runoff_diag
     1051      ALLOCATE(zxrunofflic(klon), runoff_diag(klon))
     1052      runoff_diag(:)=0.
    9611053      ALLOCATE(zxustartlic(klon), zxrhoslic(klon), zxqsaltlic(klon))
    9621054      zxustartlic(:)=0. ; zxrhoslic(:)=0. ; zxqsaltlic(:)=0.
     
    9731065      ALLOCATE(xtrain_lsc(ntraciso,klon))
    9741066      ALLOCATE(xtrunoff_diag(niso,klon))
    975       ALLOCATE(h1_diag(klon),runoff_diag(klon))
     1067      ALLOCATE(h1_diag(klon))
     1068!SN
     1069      xtrunoff_diag(:,:)=0. ! because variables are only given values on knon grid points
    9761070#endif
    9771071!
     
    10321126      ALLOCATE(wfevap(klon, nbsrf))
    10331127      ALLOCATE(evap_pot(klon, nbsrf))
    1034 ! FC
     1128! FC 
    10351129      ALLOCATE(zxfluxq(klon,klev),zxfluxt(klon,klev))
    1036 !
    10371130!
    10381131!  Deep convective variables used in phytrac
    10391132      ALLOCATE(pmflxr(klon, klev+1), pmflxs(klon, klev+1))
    10401133      ALLOCATE(wdtrainA(klon,klev),wdtrainS(klon,klev),wdtrainM(klon,klev))
    1041       ALLOCATE(dnwd(klon, klev), upwd(klon, klev) )
     1134      ALLOCATE(dnwd(klon, klev), upwd(klon, klev))
    10421135      ALLOCATE(ep(klon,klev))                          ! epmax_cape
    1043       ALLOCATE(da(klon,klev), mp(klon,klev) )
    1044       ALLOCATE(phi(klon,klev,klev) )
    1045       ALLOCATE(wght_cvfd(klon,klev) )
    1046       ALLOCATE(phi2(klon,klev,klev) )
     1136      ALLOCATE(da(klon,klev), mp(klon,klev))
     1137      ALLOCATE(phi(klon,klev,klev))
     1138      ALLOCATE(wght_cvfd(klon,klev))
     1139      ALLOCATE(phi2(klon,klev,klev))
    10471140      ALLOCATE(d1a(klon,klev), dam(klon,klev))
    1048       ALLOCATE(ev(klon,klev) )
    1049       ALLOCATE(elij(klon,klev,klev) )
    1050       ALLOCATE(qtaa(klon,klev) )
    1051       ALLOCATE(clw(klon,klev) )
    1052       ALLOCATE(epmlmMm(klon,klev,klev), eplaMm(klon,klev) )
    1053       ALLOCATE(sij(klon,klev,klev) )
     1141      ALLOCATE(ev(klon,klev))
     1142      ALLOCATE(elij(klon,klev,klev))
     1143      ALLOCATE(qtaa(klon,klev))
     1144      ALLOCATE(clw(klon,klev))
     1145      ALLOCATE(epmlmMm(klon,klev,klev), eplaMm(klon,klev))
     1146      ALLOCATE(sij(klon,klev,klev))
    10541147#ifdef ISO
    10551148      ALLOCATE(xtwdtrainA(ntraciso,klon,klev))
     
    10941187      ALLOCATE(pfraclr(klon,klev),pfracld(klon,klev))
    10951188      pfraclr(:,:)=0. ; pfracld(:,:)=0. ! because not always defined
     1189      ALLOCATE(cldfraliq(klon,klev))
     1190      ALLOCATE(sigma2_icefracturb(klon,klev))
     1191      ALLOCATE(mean_icefracturb(klon,klev))
    10961192      ALLOCATE(distcltop(klon,klev))
    10971193      ALLOCATE(temp_cltop(klon,klev))
     
    11341230
    11351231#ifdef CPP_StratAer
     1232      ALLOCATE (d_q_emiss(klon,klev))
    11361233      ALLOCATE (R2SO4(klon,klev))
     1234      ALLOCATE (R2SO4B(klon,klev,nbtr_bin))
    11371235      ALLOCATE (DENSO4(klon,klev))
     1236      ALLOCATE (DENSO4B(klon,klev,nbtr_bin))
    11381237      ALLOCATE (f_r_wet(klon,klev))
     1238      ALLOCATE (f_r_wetB(klon,klev,nbtr_bin))
    11391239      ALLOCATE (decfluxaer(klon,nbtr))
    11401240      ALLOCATE (mdw(nbtr))
     
    11471247      ALLOCATE (OCS_lifetime(klon,klev))
    11481248      ALLOCATE (SO2_lifetime(klon,klev))
     1249      ALLOCATE (H2SO4_lifetime(klon,klev))
     1250      ALLOCATE (O3_clim(klon,klev))
    11491251      ALLOCATE (alpha_bin(nbands_sw_rrtm+nbands_lw_rrtm+nwave,nbtr))
    11501252      ALLOCATE (piz_bin(nbands_sw_rrtm+nbands_lw_rrtm+nwave,nbtr))
     
    11711273      ALLOCATE (surf_PM25_sulf(klon))
    11721274      ALLOCATE (vsed_aer(klon,klev))
     1275      ALLOCATE (sulfmmr(klon,klev))
     1276      ALLOCATE (SAD_sulfate(klon,klev))
     1277      ALLOCATE (reff_sulfate(klon,klev))
     1278      ALLOCATE (sulfmmr_mode(klon,klev,nbtr_bin))
     1279      ALLOCATE (nd_mode(klon,klev,nbtr_bin))
    11731280#endif
    11741281
     
    11811288IMPLICIT NONE
    11821289      DEALLOCATE(t_seri,q_seri,ql_seri,qs_seri, qbs_seri)
     1290! SN 4D ISO
     1291      DEALLOCATE(qx_seri)
     1292! SN
    11831293      DEALLOCATE(u_seri,v_seri)
    11841294      DEALLOCATE(cf_seri,rvc_seri)
    11851295      DEALLOCATE(l_mixmin,l_mix,wprime)
     1296      DEALLOCATE(tke_shear,tke_buoy,tke_trans)
    11861297      DEALLOCATE(pbl_eps)
    11871298      DEALLOCATE(rhcl)
     
    12081319      DEALLOCATE(d_u_ajs,d_v_ajs)
    12091320      DEALLOCATE(d_t_eva,d_q_eva)
     1321! SN 4D ISO
     1322      DEALLOCATE(d_qx_eva)
     1323! SN
    12101324      DEALLOCATE(d_ql_eva,d_qi_eva)
    12111325      DEALLOCATE(d_t_lscst,d_q_lscst)
     
    12141328      DEALLOCATE(d_t_vdf,d_q_vdf,d_t_diss)
    12151329      DEALLOCATE(d_qbs_vdf)
    1216       DEALLOCATE(d_t_bs,d_q_bs,d_qbs_bs)
    1217 #ifdef ISO
    1218       deallocate(xt_seri,xtl_seri,xts_seri)
     1330      DEALLOCATE(d_t_bsss,d_q_bsss,d_qbs_bsss)
     1331#ifdef ISO
     1332      deallocate(xt_seri,xtl_seri,xts_seri,xtbs_seri)
    12191333      DEALLOCATE(d_xtl_eva,d_xti_eva)
    1220       deallocate(d_xt_dyn,d_xtl_dyn,d_xts_dyn)
     1334      deallocate(d_xt_dyn,d_xtl_dyn,d_xts_dyn,d_xtbs_dyn)
    12211335      deallocate(d_xt_con)
    12221336      deallocate(d_xt_wake)
     
    13081422      DEALLOCATE(solsw_aerop, solsw0_aerop)
    13091423      DEALLOCATE(topswcf_aerop, solswcf_aerop)
    1310 
    13111424!CK LW diagnostics
    13121425      DEALLOCATE(toplwad_aerop, sollwad_aerop)
     
    13141427      DEALLOCATE(toplwad0_aerop, sollwad0_aerop)
    13151428
     1429!AI Ajout pour Ecrad (3Deffect)
     1430      DEALLOCATE(topswad_aero_s2, solswad_aero_s2)
     1431      DEALLOCATE(topswai_aero_s2, solswai_aero_s2)
     1432      DEALLOCATE(topswad0_aero_s2, solswad0_aero_s2)
     1433      DEALLOCATE(topsw_aero_s2, topsw0_aero_s2)
     1434      DEALLOCATE(solsw_aero_s2, solsw0_aero_s2)
     1435      DEALLOCATE(topswcf_aero_s2, solswcf_aero_s2)
     1436!CK LW diagnostics
     1437      DEALLOCATE(toplwad_aero_s2, sollwad_aero_s2)
     1438      DEALLOCATE(toplwai_aero_s2, sollwai_aero_s2)
     1439      DEALLOCATE(toplwad0_aero_s2, sollwad0_aero_s2)     
     1440
    13161441! FH Ajout de celles necessaires au phys_output_write_mod
    13171442      DEALLOCATE(tal1, pal1, pab1, pab2)
     
    13221447      DEALLOCATE(wake_omg)
    13231448      DEALLOCATE(d_deltat_wk, d_deltaq_wk)
    1324       DEALLOCATE(d_s_wk, d_dens_a_wk, d_dens_wk)
     1449      DEALLOCATE(d_s_wk, d_s_a_wk, d_dens_wk, d_dens_a_wk)
    13251450      DEALLOCATE(d_deltat_wk_gw, d_deltaq_wk_gw)
    13261451      DEALLOCATE(d_deltat_vdf, d_deltaq_vdf)
     
    13531478      DEALLOCATE(uwat, vwat)
    13541479      DEALLOCATE(zxfqcalving, zxfluxlat)
    1355       DEALLOCATE(zxrunofflic)
     1480! SN runoff_diag
     1481      DEALLOCATE(zxrunofflic, runoff_diag)
    13561482      DEALLOCATE(zxustartlic, zxrhoslic, zxqsaltlic)
    13571483      DEALLOCATE(zxtsol, snow_lsc, zxfqfonte, zxqsurf)
     
    13821508      DEALLOCATE(dxtvdf_x, dxtvdf_w)
    13831509      DEALLOCATE(xt_therm)
    1384       DEALLOCATE(h1_diag,runoff_diag,xtrunoff_diag)
     1510      DEALLOCATE(h1_diag,xtrunoff_diag)
    13851511#endif
    13861512!
     
    14221548      DEALLOCATE(upwd, dnwd)
    14231549      DEALLOCATE(ep)
    1424       DEALLOCATE(da, mp )
    1425       DEALLOCATE(phi )
    1426       DEALLOCATE(wght_cvfd )
    1427       DEALLOCATE(phi2 )
     1550      DEALLOCATE(da, mp)
     1551      DEALLOCATE(phi)
     1552      DEALLOCATE(wght_cvfd)
     1553      DEALLOCATE(phi2)
    14281554      DEALLOCATE(d1a, dam)
    1429       DEALLOCATE(ev )
    1430       DEALLOCATE(elij )
    1431       DEALLOCATE(qtaa )
    1432       DEALLOCATE(clw )
    1433       DEALLOCATE(epmlmMm, eplaMm )
    1434       DEALLOCATE(sij )
     1555      DEALLOCATE(ev)
     1556      DEALLOCATE(elij)
     1557      DEALLOCATE(qtaa)
     1558      DEALLOCATE(clw)
     1559      DEALLOCATE(epmlmMm, eplaMm)
     1560      DEALLOCATE(sij)
    14351561#ifdef ISO
    14361562      DEALLOCATE(xtwdtrainA)
     
    14721598      DEALLOCATE(rneb)
    14731599      DEALLOCATE(pfraclr,pfracld)
     1600      DEALLOCATE(cldfraliq)
     1601      DEALLOCATE(sigma2_icefracturb)
     1602      DEALLOCATE(mean_icefracturb)
     1603      DEALLOCATE (zxsnow,snowhgt,qsnow,to_ice,sissnow,runoff,albsol3_lic)
    14741604      DEALLOCATE(distcltop)
    14751605      DEALLOCATE(temp_cltop)
    1476       DEALLOCATE (zxsnow,snowhgt,qsnow,to_ice,sissnow,runoff,albsol3_lic)
    14771606#ifdef ISO
    14781607      DEALLOCATE (zxxtsnow,xtVprecip,xtVprecipi,pxtrfl,pxtsfl)
     
    15071636#ifdef CPP_StratAer
    15081637! variables for strat. aerosol CK
    1509       DEALLOCATE (R2SO4)
    1510       DEALLOCATE (DENSO4)
    1511       DEALLOCATE (f_r_wet)
     1638      DEALLOCATE (d_q_emiss)
     1639      DEALLOCATE (R2SO4, R2SO4B)
     1640      DEALLOCATE (DENSO4, DENSO4B)
     1641      DEALLOCATE (f_r_wet, f_r_wetB)
    15121642      DEALLOCATE (decfluxaer)
    15131643      DEALLOCATE (mdw)
    15141644      DEALLOCATE (SO2_lifetime)
    15151645      DEALLOCATE (OCS_lifetime)
     1646      DEALLOCATE (H2SO4_lifetime)
     1647      DEALLOCATE (O3_clim)
    15161648      DEALLOCATE (alpha_bin)
    15171649      DEALLOCATE (piz_bin)
     
    15221654      DEALLOCATE (surf_PM25_sulf)
    15231655      DEALLOCATE (vsed_aer)
     1656      DEALLOCATE (sulfmmr)
     1657      DEALLOCATE (SAD_sulfate)
     1658      DEALLOCATE (reff_sulfate)
     1659      DEALLOCATE (sulfmmr_mode)
     1660      DEALLOCATE (nd_mode)
    15241661      DEALLOCATE (budg_3D_ocs_to_so2)
    15251662      DEALLOCATE (budg_3D_so2_to_h2so4)
  • LMDZ6/branches/cirrus/libf/phylmdiso/phys_output_ctrlout_mod.F90

    r5055 r5202  
    11121112  TYPE(ctrl_out), SAVE :: o_tke = ctrl_out((/ 4, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    11131113    'tke ', 'TKE', 'm2/s2', (/ ('', i=1, 10) /))
     1114  TYPE(ctrl_out), SAVE :: o_tke_shear = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     1115    'tke_shear ', 'TKE shear term', 'm2/s3', (/ ('', i=1, 10) /)) 
     1116  TYPE(ctrl_out), SAVE :: o_tke_buoy = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     1117    'tke_buoy ', 'TKE buoyancy term', 'm2/s3', (/ ('', i=1, 10) /))
     1118  TYPE(ctrl_out), SAVE :: o_tke_trans = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     1119    'tke_trans ', 'TKE transport term', 'm2/s3', (/ ('', i=1, 10) /))
    11141120  TYPE(ctrl_out), SAVE :: o_tke_dissip = ctrl_out((/ 10, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    1115     'tke_dissip ', 'TKE DISSIPATION', 'm2/s3', (/ ('', i=1, 10) /))   
     1121    'tke_dissip ', 'TKE dissipation term', 'm2/s3', (/ ('', i=1, 10) /))
     1122
    11161123  TYPE(ctrl_out), SAVE :: o_tke_max = ctrl_out((/ 4, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    11171124    'tke_max', 'TKE max', 'm2/s2',                                  &
     
    14421449  TYPE(ctrl_out), SAVE :: o_tau_strat_1020 = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 1/), &
    14431450    'OD1020_strat_only', 'Stratospheric Aerosol Optical depth at 1020 nm ', '1', (/ ('', i=1, 10) /))
     1451  TYPE(ctrl_out), SAVE :: o_SAD_sulfate = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 1/), &
     1452    'SAD_sulfate', 'SAD WET sulfate aerosols', 'cm2/cm3', (/ ('', i=1, 10) /))
     1453  TYPE(ctrl_out), SAVE :: o_reff_sulfate = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 1/), &
     1454    'reff_sulfate', 'Effective radius of WET sulfate aerosols', 'cm', (/ ('', i=1, 10) /))
     1455  TYPE(ctrl_out), SAVE :: o_sulfmmr = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 1/), &
     1456    'sulfMMR', 'Sulfate aerosol concentration (dry mass mixing ratio)', 'kg(H2SO4)/kg(air)', (/ ('', i=1, 10) /))
     1457  TYPE(ctrl_out), SAVE, ALLOCATABLE :: o_nd_mode(:)
     1458  TYPE(ctrl_out), SAVE, ALLOCATABLE :: o_sulfmmr_mode(:)
    14441459!--chemistry
    14451460  TYPE(ctrl_out), SAVE :: o_R2SO4 = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 1/), &
     
    15511566  TYPE(ctrl_out), SAVE :: o_rneb = ctrl_out((/ 2, 5, 10, 10, 10, 10, 11, 11, 11, 11/), &
    15521567    'rneb', 'Cloud fraction', '-', (/ ('', i=1, 10) /))
     1568  TYPE(ctrl_out), SAVE :: o_cldfraliq = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     1569    'cldfraliq', 'Liquid fraction of the cloud', '-', (/ ('', i=1, 10) /))
     1570  TYPE(ctrl_out), SAVE :: o_sigma2_icefracturb = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     1571    'sigma2_icefracturb', 'Variance of the diagnostic supersaturation distribution (icefrac_turb) [-]', '-', (/ ('', i=1, 10) /))
     1572  TYPE(ctrl_out), SAVE :: o_mean_icefracturb = ctrl_out((/ 11, 11, 11, 11, 11, 11, 11, 11, 11, 11/), &
     1573    'mean_icefracturb', 'Mean of the diagnostic supersaturation distribution (icefrac_turb) [-]', '-', (/ ('', i=1, 10) /))
     1574 
    15531575  TYPE(ctrl_out), SAVE :: o_rnebjn = ctrl_out((/ 2, 5, 10, 10, 10, 10, 11, 11,11, 11/), &     
    15541576    'rnebjn', 'Cloud fraction in day', '-', (/ ('', i=1, 10) /))
     
    19812003  TYPE(ctrl_out), SAVE, ALLOCATABLE :: o_dtr_dry(:)
    19822004
    1983 
    19842005#ifdef ISO
    19852006  TYPE(ctrl_out), SAVE, ALLOCATABLE :: o_xtprecip(:)
     
    19912012  TYPE(ctrl_out), SAVE, ALLOCATABLE :: o_xtoliq(:)
    19922013  TYPE(ctrl_out), SAVE, ALLOCATABLE :: o_xtcond(:)
     2014  TYPE(ctrl_out), SAVE, ALLOCATABLE :: o_xtrunoff_diag(:)
    19932015  TYPE(ctrl_out), SAVE, ALLOCATABLE :: o_dxtdyn(:)
    19942016  TYPE(ctrl_out), SAVE, ALLOCATABLE :: o_dxtldyn(:)
     
    20882110  TYPE(ctrl_out), SAVE :: o_runoff = ctrl_out((/ 1, 1, 10, 1, 10, 10, 11, 11, 11, 11/), &
    20892111    'runoff', 'Run-off rate land ice', 'kg/m2/s', (/ ('', i=1, 10) /))
     2112! SN add runoff_diag
     2113!#ifdef ISO
     2114  TYPE(ctrl_out), SAVE :: o_runoff_diag = ctrl_out((/ 1, 1, 10, 1, 10, 10, 11, 11, 11, 11/), &
     2115    'runoffland', 'Run-off rate land for bucket', 'kg/m2/s', (/ ('', i=1, 10) /))
     2116!#endif
    20902117  TYPE(ctrl_out), SAVE :: o_albslw3 = ctrl_out((/ 1, 1, 1, 1, 10, 10, 11, 11, 11, 11/), &
    20912118    'albslw3', 'Surface albedo LW3', '-', (/ ('', i=1, 10) /))
  • LMDZ6/branches/cirrus/libf/phylmdiso/physiq_mod.F90

    r5055 r5202  
    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
     41    USE infotrac_phy, ONLY: nqtot, nbtr, nqo, tracers, type_trac,ivap,iliq,isol
    4242    USE readTracFiles_mod, ONLY: addPhase
    4343    USE strings_mod,  ONLY: strIdx
     
    9393    USE phys_output_var_mod, ONLY : cloud_cover_sw, cloud_cover_sw_s2
    9494
     95
    9596    !USE cmp_seri_mod
    9697!    USE add_phys_tend_mod, only : add_pbl_tend, add_phys_tend, diag_phys_tend, prt_enerbil, &
     
    117118    USE chem_rep, ONLY: Init_chem_rep_xjour, d_q_rep, d_ql_rep, d_qi_rep, &
    118119                        ptrop, ttrop, ztrop, gravit, itroprep, Z1, Z2, fac, B
     120    USE strataer_local_var_mod
     121    USE strataer_emiss_mod, ONLY: strataer_emiss_init
    119122#endif
    120123#if defined INCA || defined REPROBUS
     
    131134
    132135#ifdef CPP_StratAer
     136    USE phys_local_var_mod, ONLY: d_q_emiss
    133137    USE strataer_local_var_mod
    134138    USE strataer_nuc_mod, ONLY: strataer_nuc_init
    135139    USE strataer_emiss_mod, ONLY: strataer_emiss_init
    136140#endif
    137 
    138141
    139142    USE lmdz_xios, ONLY: xios_update_calendar, xios_context_finalize
     
    153156        & modif_ratqs,essai_convergence,iso_init,ridicule_rain,tnat, &
    154157        & ridicule,ridicule_snow
    155     USE isotopes_routines_mod, ONLY: iso_tritium
     158    USE isotopes_routines_mod, ONLY: iso_tritium,dispatch,together
    156159#ifdef ISOVERIF
    157160    USE isotopes_verif_mod, ONLY: errmax,errmaxrel, &
     
    188191!!!!!!!!!!!!!!!!!!  END "USE" for CPP keys !!!!!!!!!!!!!!!!!!!!!!
    189192
     193USE physiqex_mod, ONLY : physiqex
    190194USE phys_local_var_mod, ONLY: phys_local_var_init, phys_local_var_end, &
    191195       ! [Variables internes non sauvegardees de la physique]
     
    193197       t_seri,q_seri,ql_seri,qs_seri,qbs_seri, &
    194198       u_seri,v_seri,cf_seri,rvc_seri,tr_seri, &
     199       rhcl, &       
     200       qx_seri, & ! CR
    195201       rhcl, &       
    196202       ! Dynamic tendencies (diagnostics)
     
    209215       !
    210216       d_t_eva,d_q_eva,d_ql_eva,d_qi_eva, &
     217       d_qx_eva, &
    211218       d_t_lsc,d_q_lsc,d_ql_lsc,d_qi_lsc, &
    212219       d_t_lscst,d_q_lscst, &
     
    219226       d_ts, &
    220227       !
    221        d_t_bs,d_q_bs,d_qbs_bs, &
     228       d_t_bsss,d_q_bsss,d_qbs_bsss, &
    222229       !
    223230!       d_t_oli,d_u_oli,d_v_oli, &
     
    247254       toplwai_aero,sollwai_aero,   &
    248255       toplwad0_aero,sollwad0_aero, &
     256       !pour Ecrad
     257       topswad_aero_s2, solswad_aero_s2,   &
     258       topswai_aero_s2, solswai_aero_s2,   &
     259       topswad0_aero_s2, solswad0_aero_s2, &
     260       topsw_aero_s2, topsw0_aero_s2,      &
     261       solsw_aero_s2, solsw0_aero_s2,      &
     262       topswcf_aero_s2, solswcf_aero_s2,   &
     263       !LW diagnostics
     264       toplwad_aero_s2, sollwad_aero_s2,   &
     265       toplwai_aero_s2, sollwai_aero_s2,   &
     266       toplwad0_aero_s2, sollwad0_aero_s2, &
    249267       !
    250268       topsw_aero,solsw_aero,       &
     
    265283       toplwai_aerop, sollwai_aerop,   &
    266284       toplwad0_aerop, sollwad0_aerop, &
     285       !pour Ecrad
     286       topswad_aero_s2, solswad_aero_s2,   &
     287       topswai_aero_s2, solswai_aero_s2,   &
     288       topswad0_aero_s2, solswad0_aero_s2, &
     289       topsw_aero_s2, topsw0_aero_s2,      &
     290       solsw_aero_s2, solsw0_aero_s2,      &
     291       topswcf_aero_s2, solswcf_aero_s2,   &
     292       !LW diagnostics
     293       toplwad_aero_s2, sollwad_aero_s2,   &
     294       toplwai_aero_s2, sollwai_aero_s2,   &
     295       toplwad0_aero_s2, sollwad0_aero_s2, &
    267296       !
    268297       ptstar, pt0, slp, &
     
    346375       !
    347376       rneblsvol, &
    348        pfraclr,pfracld, &
    349        distcltop,temp_cltop, &
     377       pfraclr, pfracld, cldfraliq, sigma2_icefracturb, mean_icefracturb, &
     378       distcltop, temp_cltop, &
    350379       !-- LSCP - condensation and ice supersaturation variables
    351380       qsub, qissr, qcld, subfra, issrfra, gamma_cond, ratio_qi_qtot, &
     
    384413
    385414#ifdef ISO
    386        USE phys_local_var_mod, ONLY: xt_seri,xtl_seri,xts_seri, &
     415       USE phys_local_var_mod, ONLY: xt_seri,xtl_seri,xts_seri,xtbs_seri, &
    387416       d_xt_eva,d_xtl_eva,d_xti_eva,           &
    388        d_xt_dyn,d_xtl_dyn,d_xts_dyn,          &
     417       d_xt_dyn,d_xtl_dyn,d_xts_dyn,d_xtbs_dyn, &
    389418       d_xt_con, d_xt_wake,                    &
    390419       d_xt_ajsb, d_xt_ajs,                    &
     
    412441       USE phys_output_var_mod, ONLY: scdnc, cldncl, reffclwtop, lcc, reffclws, &
    413442       reffclwc, cldnvi, lcc3d, lcc3dcon, lcc3dstra, icc3dcon, icc3dstra
     443       USE output_physiqex_mod, ONLY: output_physiqex
    414444
    415445
     
    556586    !
    557587    ! indices de traceurs eau vapeur, liquide, glace, fraction nuageuse LS (optional), blowing snow (optional)
    558     INTEGER,SAVE :: ivap, iliq, isol, ibs, icf, irvc
    559 !$OMP THREADPRIVATE(ivap, iliq, isol, ibs, icf, irvc)
    560     !
     588!    INTEGER,SAVE :: ivap, iliq, isol, irneb, ibs
     589!!$OMP THREADPRIVATE(ivap, iliq, isol, irneb, ibs)
     590! Camille Risi 25 juillet 2023: ivap,iliq,isol deja definis dans infotrac_phy.
     591! Et ils sont utiles ailleurs que dans physiq_mod (ex:
     592! reevap -> je commente les 2 lignes au dessus et je laisse la definition
     593! plutot dans infotrac_phy
     594    INTEGER,SAVE :: irneb, ibs, icf,irvc
     595!$OMP THREADPRIVATE(irneb, ibs, icf,irvc)
     596!
    561597    !
    562598    ! Variables argument:
     
    812848    real therm_tke_max0(klon)   ! TKE dans les thermiques au LCL
    813849    real env_tke_max0(klon)     ! TKE dans l'environnement au LCL
     850    INTEGER, SAVE :: iflag_thermcell_tke ! transtport TKE by thermals
     851    !$OMP THREADPRIVATE(iflag_thermcell_tke)
    814852
    815853!JLD    !---D\'eclenchement stochastique
     
    900938    EXTERNAL ajsec     ! ajustement sec
    901939    EXTERNAL conlmd    ! convection (schema LMD)
    902     !KE43
    903940    EXTERNAL conema3  ! convect4.3
    904     !AA
    905     ! JBM (3/14) fisrtilp_tr not loaded
    906     ! EXTERNAL fisrtilp_tr ! schema de condensation a grande echelle (pluie)
    907     !                          ! stockage des coefficients necessaires au
    908     !                          ! lessivage OFF-LINE et ON-LINE
    909941    EXTERNAL hgardfou  ! verifier les temperatures
    910942    EXTERNAL nuage     ! calculer les proprietes radiatives
     
    960992    REAL conv_t(klon,klev) ! convergence de la temperature(K/s)
    961993    !
    962 #ifdef INCA
    963     REAL zxsnow_dummy(klon)
    964 #endif
    965994    REAL zsav_tsol(klon)
    966995    !
     
    10681097    !$OMP THREADPRIVATE(ok_bug_split_th)
    10691098
     1099    ! Logical switch to a bug : modifying directly wake_deltat  by adding
     1100    ! the (w) dry adjustment tendency to wake_deltat
     1101    LOGICAL, SAVE :: ok_bug_ajs_cv = .TRUE.
     1102    !$OMP THREADPRIVATE(ok_bug_ajs_cv)
    10701103    !
    10711104    !********************************************************
     
    12051238    REAL, DIMENSION(klon,klev)     :: mass_solu_aero_pi
    12061239    ! - " - (pre-industrial value)
     1240    REAL, DIMENSION(klon,klev,naero_tot) :: m_allaer
    12071241
    12081242    ! Parameters
     
    12711305    ! Declarations pour Simulateur COSP
    12721306    !============================================================
     1307    ! AI 10-22
     1308#ifdef CPP_COSP
     1309    include "ini_COSP.h"
     1310#endif
     1311#ifdef CPP_COSPV2
     1312    include "ini_COSP.h"
     1313#endif
    12731314    real :: mr_ozone(klon,klev), phicosp(klon,klev)
    12741315
     
    13461387
    13471388    REAL, dimension(klon,klev) :: t_env,q_env
     1389#ifdef ISO
     1390    real, dimension(ntraciso,klon,klev) ::  xt_env
     1391#endif
    13481392
    13491393    REAL, dimension(klon) :: pr_et
     
    13561400    !AI namelist pour gerer le double appel de Ecrad
    13571401    CHARACTER(len=512) :: namelist_ecrad_file
     1402
     1403    !======================================================================!
     1404    ! Bifurcation vers un nouveau moniteur physique pour experimenter      !
     1405    ! des solutions et préparer le couplage avec la physique de MesoNH     !
     1406    ! 14 mai 2023                                                          !
     1407    !======================================================================!
     1408    if (debut) then                                                        !
     1409       iflag_physiq=0
     1410       call getin_p('iflag_physiq', iflag_physiq)                          !
     1411    endif                                                                  !
     1412    if ( iflag_physiq == 2 ) then
     1413#ifdef ISO
     1414       abort_message='isotopes pas encore dans physiqex'
     1415       CALL abort_physic(modname,abort_message,1)
     1416#endif
     1417       call physiqex (nlon,nlev, &                                         !
     1418       debut,lafin,pdtphys_, &                                             !
     1419       paprs,pplay,pphi,pphis,presnivs, &                                  !
     1420       u,v,rot,t,qx, &                                                     !
     1421       flxmass_w, &                                                        !
     1422       d_u, d_v, d_t, d_qx, d_ps)                                          !
     1423       return                                                              !
     1424    endif                                                                  !
     1425    !======================================================================!
     1426
    13581427
    13591428    pi = 4. * ATAN(1.)
     
    13721441    phys_tstep=NINT(pdtphys)
    13731442    IF (.NOT. using_xios) missing_val=nf90_fill_real
    1374 #ifdef CPP_XIOS
    1375 ! switch to XIOS LMDZ physics context
    1376     IF (.NOT. debut .AND. is_omp_master) THEN
    1377        CALL wxios_set_context()
    1378        CALL xios_update_calendar(itap+1)
     1443
     1444    IF (using_xios) THEN
     1445      ! switch to XIOS LMDZ physics context
     1446      IF (.NOT. debut .AND. is_omp_master) THEN
     1447        CALL wxios_set_context()
     1448        CALL xios_update_calendar(itap+1)
     1449      ENDIF
    13791450    ENDIF
    1380 #endif
    13811451
    13821452    !======================================================================
     
    13841454    ! Utilise notamment en 1D mais peut etre active egalement en 3D
    13851455    ! en imposant la valeur de igout.
    1386     !======================================================================d
     1456    !======================================================================
    13871457    IF (prt_level.ge.1) THEN
    13881458       igout=klon/2+1/klon
     
    14411511            read_climoz, &
    14421512            alp_offset)
     1513       CALL init_etat0_limit_unstruct
     1514       IF (.NOT. create_etat0_limit) CALL init_limit_read(days_elapsed)
    14431515       CALL phys_state_var_init(read_climoz)
    14441516       CALL phys_output_var_init
    14451517       IF (read_climoz>=1 .AND. create_etat0_limit .AND. grid_type==unstructured) &
    14461518          CALL regr_horiz_time_climoz(read_climoz,ok_daily_climoz)
     1519
     1520#ifdef REPROBUS
     1521       CALL strataer_init
     1522       CALL strataer_emiss_init
     1523#endif
    14471524
    14481525#ifdef CPP_StratAer
     
    14881565
    14891566        IF (ok_bs) THEN
     1567#ifdef ISO
    14901568          abort_message='blowing snow cannot be activated with water isotopes yet'
    14911569          CALL abort_physic(modname,abort_message, 1)
     
    14971575         ENDIF
    14981576        ENDIF
     1577
    14991578       Ncvpaseq1 = 0
    15001579       dnwd0=0.0
     
    15381617       tau_gl=86400.*tau_gl
    15391618       WRITE(lunout,*) 'debut physiq_mod tau_gl=',tau_gl
     1619       iflag_thermcell_tke=0
     1620       call getin_p('iflag_thermcell_tke', iflag_thermcell_tke)                          !
    15401621
    15411622       CALL getin_p('iflag_alp_wk_cond', iflag_alp_wk_cond)
     
    15601641       CALL getin_p('ok_bug_cv_trac',ok_bug_cv_trac)
    15611642       CALL getin_p('ok_bug_split_th',ok_bug_split_th)
     1643       CALL getin_p('ok_bug_ajs_cv',ok_bug_ajs_cv)
    15621644       fl_ebil = 0 ! by default, conservation diagnostics are desactivated
    15631645       CALL getin_p('fl_ebil',fl_ebil)
     
    15961678       CALL infocfields_init
    15971679
     1680       !AI 08 2023
    15981681#ifdef CPP_ECRAD
    15991682       ok_3Deffect=.false.
     
    18751958      IF (.NOT. create_etat0_limit) CALL init_readaerosolstrato(flag_aerosol_strat)  !! initialise aero strato from file for XIOS interpolation (unstructured_grid)
    18761959
     1960      ! A.I : Initialisations pour le 1er passage a Cosp
     1961      if (ok_cosp) then
     1962
    18771963#ifdef CPP_COSP
    1878       IF (ok_cosp) THEN
    1879 !           DO k = 1, klev
    1880 !             DO i = 1, klon
    1881 !               phicosp(i,k) = pphi(i,k) + pphis(i)
    1882 !             ENDDO
    1883 !           ENDDO
     1964        CALL ini_COSP(ref_liq_cosp0,ref_ice_cosp0,pctsrf_cosp0,zu10m_cosp0,zv10m_cosp0, &
     1965               zxtsol_cosp0,zx_rh_cosp0,cldfra_cosp0,rnebcon_cosp0,flwc_cosp0, &
     1966               fiwc_cosp0,prfl_cosp0,psfl_cosp0,pmflxr_cosp0,pmflxs_cosp0, &
     1967               mr_ozone_cosp0,cldtau_cosp0,cldemi_cosp0,JrNt_cosp0)
     1968
    18841969        CALL phys_cosp(itap,phys_tstep,freq_cosp, &
    18851970               ok_mensuelCOSP,ok_journeCOSP,ok_hfCOSP, &
    18861971               ecrit_mth,ecrit_day,ecrit_hf, ok_all_xml, missing_val, &
    18871972               klon,klev,longitude_deg,latitude_deg,presnivs,overlap, &
    1888                JrNt,ref_liq,ref_ice, &
    1889                pctsrf(:,is_ter)+pctsrf(:,is_lic), &
    1890                zu10m,zv10m,pphis, &
    1891                zphi,paprs(:,1:klev),pplay,zxtsol,t_seri, &
    1892                qx(:,:,ivap),zx_rh,cldfra,rnebcon,flwc,fiwc, &
    1893                prfl(:,1:klev),psfl(:,1:klev), &
    1894                pmflxr(:,1:klev),pmflxs(:,1:klev), &
    1895                mr_ozone,cldtau, cldemi)
    1896       ENDIF
    1897 #endif
    1898 
    1899 #ifdef CPP_COSP2
    1900         IF (ok_cosp) THEN
    1901 !           DO k = 1, klev
    1902 !             DO i = 1, klon
    1903 !               phicosp(i,k) = pphi(i,k) + pphis(i)
    1904 !             ENDDO
    1905 !           ENDDO
    1906           CALL phys_cosp2(itap,phys_tstep,freq_cosp, &
    1907                ok_mensuelCOSP,ok_journeCOSP,ok_hfCOSP, &
    1908                ecrit_mth,ecrit_day,ecrit_hf, ok_all_xml, missing_val, &
    1909                klon,klev,longitude_deg,latitude_deg,presnivs,overlap, &
    1910                JrNt,ref_liq,ref_ice, &
    1911                pctsrf(:,is_ter)+pctsrf(:,is_lic), &
    1912                zu10m,zv10m,pphis, &
    1913                zphi,paprs(:,1:klev),pplay,zxtsol,t_seri, &
    1914                qx(:,:,ivap),zx_rh,cldfra,rnebcon,flwc,fiwc, &
    1915                prfl(:,1:klev),psfl(:,1:klev), &
    1916                pmflxr(:,1:klev),pmflxs(:,1:klev), &
    1917                mr_ozone,cldtau, cldemi)
    1918        ENDIF
     1973               JrNt_cosp0,ref_liq_cosp0,ref_ice_cosp0, &
     1974               pctsrf_cosp0, &
     1975               zu10m_cosp0,zv10m_cosp0,pphis, &
     1976               pphi,paprs(:,1:klev),pplay,zxtsol_cosp0,t, &
     1977               qx(:,:,ivap),zx_rh_cosp0,cldfra_cosp0,rnebcon_cosp0,flwc_cosp0,fiwc_cosp0, &
     1978               prfl_cosp0(:,1:klev),psfl_cosp0(:,1:klev), &
     1979               pmflxr_cosp0(:,1:klev),pmflxs_cosp0(:,1:klev), &
     1980               mr_ozone_cosp0,cldtau_cosp0, cldemi_cosp0)
    19191981#endif
    19201982
    19211983#ifdef CPP_COSPV2
    1922         IF (ok_cosp) THEN
    1923            DO k = 1, klev
    1924              DO i = 1, klon
    1925                phicosp(i,k) = pphi(i,k) + pphis(i)
    1926              ENDDO
    1927            ENDDO
     1984          CALL ini_COSP(ref_liq_cosp0,ref_ice_cosp0,pctsrf_cosp0,zu10m_cosp0,zv10m_cosp0, &
     1985               zxtsol_cosp0,zx_rh_cosp0,cldfra_cosp0,rnebcon_cosp0,flwc_cosp0, &
     1986               fiwc_cosp0,prfl_cosp0,psfl_cosp0,pmflxr_cosp0,pmflxs_cosp0, &
     1987               mr_ozone_cosp0,cldtau_cosp0,cldemi_cosp0,JrNt_cosp0)
     1988
    19281989          CALL lmdz_cosp_interface(itap,phys_tstep,freq_cosp, &
    19291990               ok_mensuelCOSP,ok_journeCOSP,ok_hfCOSP, &
    19301991               ecrit_mth,ecrit_day,ecrit_hf, ok_all_xml, missing_val, &
    19311992               klon,klev,longitude_deg,latitude_deg,presnivs,overlap, &
    1932                JrNt,ref_liq,ref_ice, &
    1933                pctsrf(:,is_ter)+pctsrf(:,is_lic), &
    1934                zu10m,zv10m,pphis, &
    1935                phicosp,paprs(:,1:klev),pplay,zxtsol,t_seri, &
    1936                qx(:,:,ivap),zx_rh,cldfra,rnebcon,flwc,fiwc, &
    1937                prfl(:,1:klev),psfl(:,1:klev), &
    1938                pmflxr(:,1:klev),pmflxs(:,1:klev), &
    1939                mr_ozone,cldtau, cldemi)
    1940        ENDIF
    1941 #endif
     1993               JrNt_cosp0,ref_liq_cosp0,ref_ice_cosp0, &
     1994               pctsrf_cosp0, &
     1995               zu10m_cosp0,zv10m_cosp0,pphis, &
     1996               pphi,paprs(:,1:klev),pplay,zxtsol_cosp0,t, &
     1997               qx(:,:,ivap),zx_rh_cosp0,cldfra_cosp0,rnebcon_cosp0,flwc_cosp0,fiwc_cosp0, &
     1998               prfl_cosp0(:,1:klev),psfl_cosp0(:,1:klev), &
     1999               pmflxr_cosp0(:,1:klev),pmflxs_cosp0(:,1:klev), &
     2000               mr_ozone_cosp0,cldtau_cosp0, cldemi_cosp0)
     2001#endif
     2002
     2003      endif !ok_cosp
    19422004
    19432005       !
     
    20142076
    20152077
    2016 #ifdef CPP_XIOS
    2017        IF (is_omp_master) CALL xios_update_calendar(1)
    2018 #endif
     2078       IF (using_xios) THEN
     2079         IF (is_omp_master) CALL xios_update_calendar(1)
     2080       ENDIF
     2081       
    20192082       IF(read_climoz>=1 .AND. create_etat0_limit) CALL regr_horiz_time_climoz(read_climoz,ok_daily_climoz)
    20202083       CALL create_etat0_limit_unstruct
     
    22132276       ENDIF
    22142277
    2215       IF (using_xios) THEN
    2216 ! Need to put this initialisation after phyetat0 as in the coupled model the XIOS context is only
    2217 ! initialised at that moment
    2218        ! Get "missing_val" value from XML files (from temperature variable)
    2219         IF (is_omp_master) CALL xios_get_field_attr("temp",default_value=missing_val)
    2220         CALL bcast_omp(missing_val)
    2221        
     2278       IF (using_xios) THEN   
     2279         ! Need to put this initialisation after phyetat0 as in the coupled model the XIOS context is only
     2280         ! initialised at that moment
     2281         ! Get "missing_val" value from XML files (from temperature variable)
     2282         IF (is_omp_master) CALL xios_get_field_attr("temp",default_value=missing_val)
     2283         CALL bcast_omp(missing_val)
    22222284       !
    22232285       ! Now we activate some double radiation call flags only if some
    22242286       ! diagnostics are requested, otherwise there is no point in doing this
    2225        IF (is_master) THEN
    2226          !--setting up swaero_diag to TRUE in XIOS case
    2227          IF (xios_field_is_active("topswad").OR.xios_field_is_active("topswad0").OR. &
    2228             xios_field_is_active("solswad").OR.xios_field_is_active("solswad0").OR. &
    2229             xios_field_is_active("topswai").OR.xios_field_is_active("solswai").OR.  &
    2230               (iflag_rrtm==1.AND.(xios_field_is_active("toplwad").OR.xios_field_is_active("toplwad0").OR. &
    2231                                   xios_field_is_active("sollwad").OR.xios_field_is_active("sollwad0"))))  &
    2232             !!!--for now these fields are not in the XML files so they are omitted
    2233             !!!  xios_field_is_active("toplwai").OR.xios_field_is_active("sollwai") !))) &
    2234             swaero_diag=.TRUE.
     2287         IF (is_master) THEN
     2288           !--setting up swaero_diag to TRUE in XIOS case
     2289           IF (xios_field_is_active("topswad").OR.xios_field_is_active("topswad0").OR. &
     2290              xios_field_is_active("solswad").OR.xios_field_is_active("solswad0").OR. &
     2291              xios_field_is_active("topswai").OR.xios_field_is_active("solswai").OR.  &
     2292                (iflag_rrtm==1.AND.(xios_field_is_active("toplwad").OR.xios_field_is_active("toplwad0").OR. &
     2293                                    xios_field_is_active("sollwad").OR.xios_field_is_active("sollwad0"))))  &
     2294              !!!--for now these fields are not in the XML files so they are omitted
     2295              !!!  xios_field_is_active("toplwai").OR.xios_field_is_active("sollwai") !))) &
     2296              swaero_diag=.TRUE.
    22352297 
    2236          !--setting up swaerofree_diag to TRUE in XIOS case
    2237          IF (xios_field_is_active("SWdnSFCcleanclr").OR.xios_field_is_active("SWupSFCcleanclr").OR. &
    2238             xios_field_is_active("SWupTOAcleanclr").OR.xios_field_is_active("rsucsaf").OR.   &
    2239             xios_field_is_active("rsdcsaf") .OR. xios_field_is_active("LWdnSFCcleanclr").OR. &
    2240             xios_field_is_active("LWupTOAcleanclr")) &
    2241             swaerofree_diag=.TRUE.
     2298           !--setting up swaerofree_diag to TRUE in XIOS case
     2299           IF (xios_field_is_active("SWdnSFCcleanclr").OR.xios_field_is_active("SWupSFCcleanclr").OR. &
     2300              xios_field_is_active("SWupTOAcleanclr").OR.xios_field_is_active("rsucsaf").OR.   &
     2301              xios_field_is_active("rsdcsaf") .OR. xios_field_is_active("LWdnSFCcleanclr").OR. &
     2302              xios_field_is_active("LWupTOAcleanclr")) &
     2303              swaerofree_diag=.TRUE.
    22422304 
    2243          !--setting up dryaod_diag to TRUE in XIOS case
    2244          DO naero = 1, naero_tot-1
    2245           IF (xios_field_is_active("dryod550_"//name_aero_tau(naero))) dryaod_diag=.TRUE.
    2246          ENDDO
    2247          !
    2248          !--setting up ok_4xCO2atm to TRUE in XIOS case
    2249          IF (xios_field_is_active("rsut4co2").OR.xios_field_is_active("rlut4co2").OR. &
    2250             xios_field_is_active("rsutcs4co2").OR.xios_field_is_active("rlutcs4co2").OR. &
    2251             xios_field_is_active("rsu4co2").OR.xios_field_is_active("rsucs4co2").OR. &
    2252             xios_field_is_active("rsd4co2").OR.xios_field_is_active("rsdcs4co2").OR. &
    2253             xios_field_is_active("rlu4co2").OR.xios_field_is_active("rlucs4co2").OR. &
    2254             xios_field_is_active("rld4co2").OR.xios_field_is_active("rldcs4co2")) &
    2255             ok_4xCO2atm=.TRUE.
    2256        ENDIF
    2257        !$OMP BARRIER
    2258        CALL bcast(swaero_diag)
    2259        CALL bcast(swaerofree_diag)
    2260        CALL bcast(dryaod_diag)
    2261        CALL bcast(ok_4xCO2atm)
    2262 
    2263      ENDIF !using_xios
    2264 
     2305           !--setting up dryaod_diag to TRUE in XIOS case
     2306           DO naero = 1, naero_tot-1
     2307             IF (xios_field_is_active("dryod550_"//name_aero_tau(naero))) dryaod_diag=.TRUE.
     2308           ENDDO
     2309           !
     2310          !--setting up ok_4xCO2atm to TRUE in XIOS case
     2311           IF (xios_field_is_active("rsut4co2").OR.xios_field_is_active("rlut4co2").OR. &
     2312              xios_field_is_active("rsutcs4co2").OR.xios_field_is_active("rlutcs4co2").OR. &
     2313              xios_field_is_active("rsu4co2").OR.xios_field_is_active("rsucs4co2").OR. &
     2314              xios_field_is_active("rsd4co2").OR.xios_field_is_active("rsdcs4co2").OR. &
     2315              xios_field_is_active("rlu4co2").OR.xios_field_is_active("rlucs4co2").OR. &
     2316              xios_field_is_active("rld4co2").OR.xios_field_is_active("rldcs4co2")) &
     2317              ok_4xCO2atm=.TRUE.
     2318           ENDIF
     2319           !$OMP BARRIER
     2320           CALL bcast(swaero_diag)
     2321           CALL bcast(swaerofree_diag)
     2322           CALL bcast(dryaod_diag)
     2323           CALL bcast(ok_4xCO2atm)
     2324         ENDIF !using_xios
    22652325       !
    22662326       CALL printflag( tabcntr0,radpas,ok_journe, &
     
    25492609          u_seri(i,k)  = u(i,k)
    25502610          v_seri(i,k)  = v(i,k)
     2611          qx_seri(i,k,:)  = qx(i,k,:)
    25512612          q_seri(i,k)  = qx(i,k,ivap)
    25522613          ql_seri(i,k) = qx(i,k,iliq)
     
    25902651      DO k = 1, klev
    25912652       DO i = 1, klon
     2653          xtbs_seri(ixt,i,k) = 0.
    25922654          xt_seri(ixt,i,k)  = qx(i,k,iqIsoPha(ixt,ivap))
    25932655          xtl_seri(ixt,i,k) = qx(i,k,iqIsoPha(ixt,iliq))
     
    26102672    qql1(:)=0.0
    26112673    DO k = 1, klev
    2612       qql1(:)=qql1(:)+(q_seri(:,k)+ql_seri(:,k)+qs_seri(:,k)+qbs_seri(:,k))*zmasse(:,k)
     2674      qql1(:)=qql1(:)+(q_seri(:,k)+ql_seri(:,k))*zmasse(:,k)
     2675      IF (nqo >= 3) THEN
     2676        qql1(:)=qql1(:)+qs_seri(:,k)*zmasse(:,k)
     2677      ENDIF
     2678      IF (ok_bs) THEN
     2679        qql1(:)=qql1(:)+qbs_seri(:,k)*zmasse(:,k)
     2680      ENDIF
    26132681    ENDDO
    26142682#ifdef ISO
    2615 #ifdef ISOVERIF
    2616         write(*,*) 'physiq tmp 1913'
    2617 #endif
    2618     do ixt=1,ntraciso
     2683    DO ixt=1,ntraciso
    26192684    xtql1(ixt,:)=0.0
    26202685    DO k = 1, klev
    2621       xtql1(ixt,:)=xtql1(ixt,:)+(xt_seri(ixt,:,k)+xtl_seri(ixt,:,k)+xts_seri(ixt,:,k))*zmasse(:,k)
    2622     ENDDO
    2623     enddo !do ixt=1,ntraciso
     2686      xtql1(ixt,:)=xtql1(ixt,:)+(xt_seri(ixt,:,k)+xtl_seri(ixt,:,k))*zmasse(:,k)
     2687      IF (nqo >= 3) THEN
     2688        xtql1(ixt,:)=xtql1(ixt,:)+xts_seri(ixt,:,k)*zmasse(:,k)
     2689      ENDIF
     2690      IF (ok_bs) THEN
     2691        xtql1(ixt,:)=xtql1(ixt,:)+xtbs_seri(ixt,:,k)*zmasse(:,k)
     2692      ENDIF
     2693    ENDDO !DO k = 1, klev
     2694    ENDDO !DO ixt=1,ntraciso
    26242695#endif
    26252696    ENDIF
     
    26332704         IF(.NOT.tracers(iq)%isInPhysics) CYCLE
    26342705         itr = itr+1
    2635 !#ifdef ISOVERIF
    2636 !         write(*,*) 'physiq 1973: itr,iq=',itr,iq
    2637 !         write(*,*) 'qx(1,1,iq)=',qx(1,1,iq)
    2638 !#endif
    2639          DO  k = 1, klev
     2706          DO  k = 1, klev
    26402707             DO  i = 1, klon
    26412708                tr_seri(i,k,itr) = qx(i,k,iq)
     
    27522819              d_xts_dyn(ixt,i,k) =  &
    27532820     &           (xts_seri(ixt,i,k)-xts_ancien(ixt,i,k))/phys_tstep
     2821              d_xtbs_dyn(ixt,i,k) =  &
     2822     &           (xtbs_seri(ixt,i,k)-xtbs_ancien(ixt,i,k))/phys_tstep
    27542823            enddo ! do ixt=1,ntraciso
    27552824         ENDDO
     
    27652834           call iso_verif_noNaN(d_xtl_dyn(ixt,i,k),'physiq 2220d')
    27662835           call iso_verif_noNaN(d_xts_dyn(ixt,i,k),'physiq 2220e')
     2836           call iso_verif_noNaN(d_xtbs_dyn(ixt,i,k),'physiq 2220f')
    27672837           enddo ! do ixt=1,ntraciso
    27682838         enddo
     
    28482918       ! !! RomP <<<
    28492919       ancien_ok = .TRUE.
     2920#ifdef ISO
     2921       d_xtbs_dyn(:,:,:)=0.0
     2922#endif
    28502923    ENDIF
    28512924    !
     
    29863059    ! Re-evaporer l'eau liquide nuageuse
    29873060    !
    2988      CALL reevap (klon,klev,iflag_ice_thermo,t_seri,q_seri,ql_seri,qs_seri, &
    2989    &         d_t_eva,d_q_eva,d_ql_eva,d_qi_eva &
    2990 #ifdef ISO
    2991              ,xt_seri,xtl_seri,xts_seri,d_xt_eva,d_xtl_eva,d_xti_eva &
    2992 #endif
    2993    &     )
     3061     CALL reevap (klon,klev,iflag_ice_thermo,t_seri,qx_seri, &
     3062   &         d_t_eva,d_qx_eva)
     3063
     3064     call dispatch(klon,klev,qx_seri,q_seri,xt_seri,ql_seri,xtl_seri,qs_seri,xts_seri)
     3065     call dispatch(klon,klev,d_qx_eva,d_q_eva,d_xt_eva,d_ql_eva,d_xtl_eva,d_qi_eva,d_xti_eva)
     3066
     3067
     3068#ifdef ISO
     3069#ifdef ISOVERIF
     3070 DO k = 1, klev
     3071     DO i = 1, klon
     3072      do ixt=1,ntraciso
     3073      call iso_verif_noNaN(xt_seri(ixt,i,k), &
     3074     &     'reevap 2417: apres evap tot')
     3075      enddo
     3076      if (iso_eau.gt.0) then
     3077              call iso_verif_egalite_choix( &
     3078     &           xt_seri(iso_eau,i,k),q_seri(i,k), &
     3079     &          'reevap 1891, après réévap totale',errmax,errmaxrel)
     3080              call iso_verif_egalite_choix( &
     3081     &           xtl_seri(iso_eau,i,k),ql_seri(i,k), &
     3082     &          'reevap 2209, après réévap totale',errmax,errmaxrel)
     3083              call iso_verif_egalite_choix( &
     3084     &           xts_seri(iso_eau,i,k),qs_seri(i,k), &
     3085     &          'reevap 2209b, après réévap totale',errmax,errmaxrel)
     3086       endif !if (iso_eau.gt.0) then
     3087     
     3088      if ((iso_HDO.gt.0).and.(iso_O18.gt.0)) then 
     3089            if (q_seri(i,k).gt.ridicule) then 
     3090               if (iso_verif_o18_aberrant_nostop( &
     3091     &              xt_seri(iso_HDO,i,k)/q_seri(i,k), &
     3092     &              xt_seri(iso_O18,i,k)/q_seri(i,k), &
     3093     &              'reevap 2408: apres reevap totale').eq.1) then
     3094                  write(*,*) 'i,k,q_seri(i,k)=',i,k,q_seri(i,k)                       
     3095                  stop
     3096              endif !  if (iso_verif_o18_aberrant_nostop
     3097            endif !if (q_seri(i,k).gt.errmax) then   
     3098        endif !if ((iso_HDO.gt.0).and.(iso_O18.gt.0)) then       
     3099#ifdef ISOTRAC     
     3100             call iso_verif_traceur(xt_seri(1,i,k), &
     3101     &           'reevap 2165c')
     3102             call iso_verif_traceur_pbidouille(xt_seri(1,i,k), &
     3103     &           'reevap 2165d')
     3104#endif
     3105       ENDDO
     3106    ENDDO               
     3107#endif                 
     3108#endif
     3109
    29943110
    29953111     CALL add_phys_tend &
     
    31233239    ! Calcul de l'humidite de saturation au niveau du sol
    31243240
     3241! Tests Fredho, instensibilite au pas de temps -------------------------------
     3242! A detruire en 2024 une fois les tests documentes et les choix faits        !
     3243! Conservation des variables avant l'appel à l a diffusion pour les tehrmic  !
     3244    if (iflag_thermals_tenv / 10 == 1 ) then                                 !
     3245        do k=1,klev                                                          !
     3246           do i=1,klon                                                       !
     3247              t_env(i,k)=t_seri(i,k)                                         !
     3248              q_env(i,k)=q_seri(i,k)   
     3249#ifdef ISO
     3250              do ixt=1,ntraciso
     3251                xt_env(ixt,i,k)=xt_seri(ixt,i,k)
     3252              enddo
     3253#endif
     3254           enddo                                                             !
     3255        enddo                                                                !
     3256    else if (iflag_thermals_tenv / 10 == 2 ) then                            !
     3257        do k=1,klev                                                          !
     3258           do i=1,klon                                                       !
     3259              t_env(i,k)=t_seri(i,k)                                         !
     3260           enddo                                                             !
     3261        enddo                                                                !
     3262    endif                                                                    !
     3263! Tests Fredho, instensibilite au pas de temps -------------------------------
    31253264
    31263265
     
    33113450          d_deltaq_vdf(:,:) = d_q_vdf_w(:,:)-d_q_vdf_x(:,:)
    33123451          CALL add_wake_tend &
    3313              (d_deltat_vdf, d_deltaq_vdf, dsig0, ddens0, ddens0, wkoccur1, 'vdf', abortphy &
     3452             (d_deltat_vdf, d_deltaq_vdf, dsig0, dsig0, ddens0, ddens0, wkoccur1, 'vdf', abortphy &
    33143453#ifdef ISO
    33153454               ,d_deltaxt_vdf &
     
    33443483     &   )
    33453484       ENDIF
    3346 #ifdef ISOVERIF
    3347         write(*,*) 'physiq tmp 2736'
    3348 #endif
    3349 
    33503485       CALL prt_enerbil('vdf',itap)
     3486
    33513487       !--------------------------------------------------------------------
    33523488
     
    34033539    ! Blowing snow sublimation and sedimentation
    34043540
    3405     d_t_bs(:,:)=0.
    3406     d_q_bs(:,:)=0.
    3407     d_qbs_bs(:,:)=0.
     3541    d_t_bsss(:,:)=0.
     3542    d_q_bsss(:,:)=0.
     3543    d_qbs_bsss(:,:)=0.
    34083544    bsfl(:,:)=0.
    34093545    bs_fall(:)=0.
     
    34113547
    34123548     CALL call_blowing_snow_sublim_sedim(klon,klev,phys_tstep,t_seri,q_seri,qbs_seri,pplay,paprs, &
    3413                                         d_t_bs,d_q_bs,d_qbs_bs,bsfl,bs_fall)
     3549                                        d_t_bsss,d_q_bsss,d_qbs_bsss,bsfl,bs_fall)
    34143550
    34153551     CALL add_phys_tend &
    3416                (du0,dv0,d_t_bs,d_q_bs,dql0,dqi0,d_qbs_bs,paprs,&
     3552               (du0,dv0,d_t_bsss,d_q_bsss,dql0,dqi0,d_qbs_bsss,paprs,&
    34173553               'bs',abortphy,flag_inhib_tend,itap,0  &
    34183554#ifdef ISO                                       
     
    37133849                ENDDO
    37143850             ENDDO
    3715              IF (iflag_adjwk == 2) THEN
     3851             IF (iflag_adjwk == 2 .AND. OK_bug_ajs_cv) THEN
    37163852               CALL add_wake_tend &
    3717                  (d_deltat_ajs_cv, d_deltaq_ajs_cv, dsig0, ddens0, ddens0, wkoccur1, 'ajs_cv', abortphy &
     3853                 (d_deltat_ajs_cv, d_deltaq_ajs_cv, dsig0, dsig0, ddens0, ddens0, wkoccur1, 'ajs_cv', abortphy &
    37183854#ifdef ISO
    37193855                      ,d_deltaxt_ajs_cv  &
    37203856#endif
    37213857                 )
    3722              ENDIF  ! (iflag_adjwk == 2)
     3858             ENDIF  ! (iflag_adjwk == 2 .AND. OK_bug_ajs_cv)
    37233859          ENDIF  ! (iflag_adjwk >= 1)
    37243860       ENDIF ! (iflag_wake>=1)
     
    44244560       !  ====
    44254561       IF (prt_level>9) WRITE(lunout,*)'pas de convection seche'
     4562       WRITE(lunout,*) 'WARNING : running without dry convection. Somme intermediate variables are not properly defined in physiq_mod.F90'
     4563       ! Reprendre proprement les initialisation ci dessouds si on veut vraiment utiliser l'option (FH)
     4564          fraca(:,:)=0.
     4565          fm_therm(:,:)=0.
     4566          ztv(:,:)=t_seri(:,:)
     4567          zqasc(:,:)=q_seri(:,:)
     4568          ztla(:,:)=0.
     4569          zthl(:,:)=0.
     4570          zpspsk(:,:)=(pplay(:,:)/100000.)**RKAPPA
    44264571
    44274572
     
    45154660
    45164661       IF (iflag_thermals>=1) THEN
     4662
     4663! Tests Fredho, instensibilite au pas de temps -------------------------------
     4664! A detruire en 2024 une fois les tests documentes et les choix faits        !
     4665          if (iflag_thermals_tenv /10 == 0 ) then                            !
     4666            do k=1,klev                                                      !
     4667               do i=1,klon                                                   !
     4668                  t_env(i,k)=t_seri(i,k)                                     !
     4669                  q_env(i,k)=q_seri(i,k)                                     !
     4670#ifdef ISO
     4671                  do ixt=1,ntraciso
     4672                        xt_env(ixt,i,k)=xt_seri(ixt,i,k)
     4673                  enddo
     4674#endif
     4675               enddo                                                         !
     4676            enddo                                                            !
     4677          else if (iflag_thermals_tenv / 10 == 2 ) then                      !
     4678            do k=1,klev                                                      !
     4679               do i=1,klon                                                   !
     4680                  q_env(i,k)=q_seri(i,k)                                     !
     4681#ifdef ISO
     4682                  do ixt=1,ntraciso
     4683                        xt_env(ixt,i,k)=xt_seri(ixt,i,k)
     4684                  enddo
     4685#endif
     4686               enddo                                                         !
     4687            enddo                                                            !
     4688          else if (iflag_thermals_tenv / 10 == 3 ) then                      !
     4689            do k=1,klev                                                      !
     4690               do i=1,klon                                                   !
     4691                  t_env(i,k)=t(i,k)                                          !
     4692                  q_env(i,k)=qx(i,k,1)                                       !
     4693#ifdef ISO
     4694                  do ixt=1,ntraciso
     4695                        xt_env(ixt,i,k)=xt_seri(ixt,i,k)
     4696                  enddo
     4697#endif
     4698               enddo                                                         !
     4699            enddo                                                            !
     4700          endif                                                              !
     4701! Tests Fredho, instensibilite au pas de temps ------------------------------
     4702
    45174703          !jyg<
    45184704!!       IF (mod(iflag_pbl_split/2,2) .EQ. 1) THEN
     
    45234709                   t_therm(i,k) = t_seri(i,k) - wake_s(i)*wake_deltat(i,k)
    45244710                   q_therm(i,k) = q_seri(i,k) - wake_s(i)*wake_deltaq(i,k)
     4711                   t_env(i,k)   = t_env(i,k) - wake_s(i)*wake_deltat(i,k)
     4712                   q_env(i,k)   = q_env(i,k) - wake_s(i)*wake_deltaq(i,k)
    45254713                   u_therm(i,k) = u_seri(i,k)
    45264714                   v_therm(i,k) = v_seri(i,k)
     
    45284716                   do ixt=1,ntraciso
    45294717                     xt_therm(ixt,i,k) = xt_seri(ixt,i,k) - wake_s(i)*wake_deltaxt(ixt,i,k)
     4718                     xt_env(ixt,i,k) = xt_env(ixt,i,k) - wake_s(i)*wake_deltaxt(ixt,i,k)
    45304719                   enddo !do ixt=1,ntraciso
    45314720#endif
     
    45704759               ,pplay,paprs,pphi,weak_inversion &
    45714760                        ! ,u_seri,v_seri,t_seri,q_seri,zqsat,debut & !jyg
    4572                ,u_therm,v_therm,t_therm,q_therm,t_therm,q_therm,zqsat,debut &  !jyg
     4761               ,u_therm,v_therm,t_therm,q_therm,t_env,q_env,zqsat,debut &  !jyg
    45734762               ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs &
    45744763               ,fm_therm,entr_therm,detr_therm &
     
    45894778               ,zqla,ztva &
    45904779#ifdef ISO         
    4591      &      ,xt_therm,d_xt_ajs &
     4780     &      ,xt_env,d_xt_ajs &
    45924781#ifdef DIAGISO         
    45934782     &      ,q_the,xt_the &
     
    46244813             IF (ok_bug_split_th) THEN
    46254814               CALL add_wake_tend &
    4626                    (d_deltat_the, d_deltaq_the, dsig0, ddens0, ddens0, wkoccur1, 'the', abortphy &
     4815                   (d_deltat_the, d_deltaq_the, dsig0, dsig0, ddens0, ddens0, wkoccur1, 'the', abortphy &
    46274816#ifdef ISO
    46284817                   ,d_deltaxt_the &
     
    46314820             ELSE
    46324821               CALL add_wake_tend &
    4633                    (d_deltat_the, d_deltaq_the, dsig0, ddens0, ddens0, wake_k, 'the', abortphy &
     4822                   (d_deltat_the, d_deltaq_the, dsig0, dsig0, ddens0, ddens0, wake_k, 'the', abortphy &
    46344823#ifdef ISO
    46354824                   ,d_deltaxt_the &
     
    46604849          ! Transport de la TKE par les panaches thermiques.
    46614850          ! FH : 2010/02/01
    4662           !     if (iflag_pbl.eq.10) then
    4663           !     call thermcell_dtke(klon,klev,nbsrf,pdtphys,fm_therm,entr_therm,
    4664           !    s           rg,paprs,pbl_tke)
    4665           !     endif
     4851               if (iflag_thermcell_tke==1) then
     4852               call thermcell_dtke(klon,klev,nbsrf,pdtphys,fm_therm,entr_therm,rg,paprs,pbl_tke)
     4853               endif
    46664854          ! -------------------------------------------------------------------
    46674855
     
    49025090
    49035091    CALL lscp(klon,klev,phys_tstep,missing_val,paprs,pplay, &
    4904          t_seri, q_seri,ptconv,ratqs, &
     5092         t_seri, q_seri,qs_ancien,ptconv,ratqs, &
    49055093         d_t_lsc, d_q_lsc, d_ql_lsc, d_qi_lsc, rneb, rneblsvol, &
    4906          pfraclr,pfracld, &
     5094         pfraclr, pfracld, cldfraliq, sigma2_icefracturb, mean_icefracturb, &
    49075095         radocond, picefra, rain_lsc, snow_lsc, &
    49085096         frac_impa, frac_nucl, beta_prec_fisrt, &
    49095097         prfl, psfl, rhcl,  &
    49105098         zqasc, fraca,ztv,zpspsk,ztla,zthl,iflag_cld_th, &
    4911          iflag_ice_thermo, distcltop, temp_cltop, cell_area, &
    4912          cf_seri, rvc_seri, u_seri, v_seri, pbl_eps(:,:,is_ave), &
     5099         iflag_ice_thermo, distcltop, temp_cltop,
     5100         pbl_tke(:,:,is_ave), pbl_eps(:,:,is_ave), &
     5101         cell_area, &
     5102         cf_seri, rvc_seri, u_seri, v_seri, &
    49135103         qsub, qissr, qcld, subfra, issrfra, gamma_cond, ratio_qi_qtot, &
    49145104         dcf_sub, dcf_con, dcf_mix, dqi_adj, dqi_sub, dqi_con, dqi_mix, &
     
    49235113    ELSE
    49245114
     5115    ! Camille Risi mai 2024: on ne met pas à jour ici pour ne pas s'mbêter à modifier fisrtilp
    49255116    CALL fisrtilp(phys_tstep,paprs,pplay, &
    49265117         t_seri, q_seri,ptconv,ratqs, &
     
    55225713                     tausum_aero, tau3d_aero)
    55235714             ENDIF
    5524           ELSE                       ! RRTM radiation
     5715          ELSE IF (iflag_rrtm .EQ.1) THEN  ! RRTM radiation
    55255716             IF (aerosol_couple .AND. config_inca == 'aero' ) THEN
    55265717                abort_message='config_inca=aero et rrtm=1 impossible'
     
    55885779                !
    55895780             ENDIF
     5781          ELSE IF (iflag_rrtm .EQ.2) THEN    ! ecrad RADIATION
     5782#ifdef CPP_ECRAD
     5783             !--climatologies or INCA aerosols
     5784             CALL readaerosol_optic_ecrad( debut, aerosol_couple, ok_alw, ok_volcan, &
     5785                  flag_aerosol, flag_bc_internal_mixture, itap, jD_cur-jD_ref, &
     5786                  pdtphys, pplay, paprs, t_seri, rhcl, presnivs,  &
     5787                  tr_seri, mass_solu_aero, mass_solu_aero_pi, m_allaer)
     5788#else
     5789                abort_message='You should compile with -rad ecrad if running with iflag_rrtm=2'
     5790                CALL abort_physic(modname,abort_message,1)
     5791#endif
    55905792          ENDIF
    55915793       ELSE   !--flag_aerosol = 0
     
    58286030               ! Rajoute par OB pour RRTM
    58296031               tau_aero_lw_rrtm, &
    5830                cldtaupirad, &
     6032               cldtaupirad, m_allaer, &
    58316033!              zqsat, flwcrad, fiwcrad, &
    58326034               zqsat, flwc, fiwc, &
     
    59076109                                ! Rajoute par OB pour RRTM
    59086110                     tau_aero_lw_rrtm, &
    5909                      cldtaupi, &
     6111                     cldtaupi, m_allaer, &
    59106112!                    zqsat, flwcrad, fiwcrad, &
    59116113                     zqsat, flwc, fiwc, &
     
    59346136                     cloud_cover_sw)
    59356137          ENDIF !ok_4xCO2atm
     6138
     6139! A.I aout 2023
     6140! Effet 3D des nuages Ecrad
     6141! a passer : nom du ficher namelist et cles ok_3Deffect
     6142! a declarer comme iflag_rrtm et a lire dans physiq.def
     6143#ifdef CPP_ECRAD
     6144          IF (ok_3Deffect) then
     6145!                print*,'ok_3Deffect = ',ok_3Deffect 
     6146                namelist_ecrad_file='namelist_ecrad_s2'
     6147                CALL radlwsw &
     6148                     (debut, dist, rmu0, fract,  &
     6149                     paprs, pplay,zxtsol,SFRWL,albsol_dir, albsol_dif, &
     6150                     t_seri,q_seri,wo, &
     6151                     cldfrarad, cldemirad, cldtaurad, &
     6152                     ok_ade.OR.flag_aerosol_strat.GT.0, ok_aie,  ok_volcan, flag_volc_surfstrat, &
     6153                     flag_aerosol, flag_aerosol_strat, flag_aer_feedback, &
     6154                     tau_aero, piz_aero, cg_aero, &
     6155                     tau_aero_sw_rrtm, piz_aero_sw_rrtm, cg_aero_sw_rrtm, &
     6156                     tau_aero_lw_rrtm, &
     6157                     cldtaupi, &
     6158                     zqsat, flwc, fiwc, &
     6159                     ref_liq, ref_ice, ref_liq_pi, ref_ice_pi, &
     6160                     namelist_ecrad_file, &
     6161! A modifier             
     6162                     heat_s2,heat0_s2,cool_s2,cool0_s2,albpla_s2, &
     6163                     heat_volc,cool_volc, &
     6164                     topsw_s2,toplw_s2,solsw_s2,solswfdiff_s2,sollw_s2, &
     6165                     sollwdown_s2, &
     6166                     topsw0_s2,toplw0_s2,solsw0_s2,sollw0_s2, &
     6167                     lwdnc0_s2, lwdn0_s2, lwdn_s2, lwupc0_s2, lwup0_s2, lwup_s2,  &
     6168                     swdnc0_s2, swdn0_s2, swdn_s2, swupc0_s2, swup0_s2, swup_s2, &
     6169                     topswad_aero_s2, solswad_aero_s2, &
     6170                     topswai_aero_s2, solswai_aero_s2, &
     6171                     topswad0_aero_s2, solswad0_aero_s2, &
     6172                     topsw_aero_s2, topsw0_aero_s2, &
     6173                     solsw_aero_s2, solsw0_aero_s2, &
     6174                     topswcf_aero_s2, solswcf_aero_s2, &
     6175                                !-C. Kleinschmitt for LW diagnostics
     6176                     toplwad_aero_s2, sollwad_aero_s2,&
     6177                     toplwai_aero_s2, sollwai_aero_s2, &
     6178                     toplwad0_aero_s2, sollwad0_aero_s2,&
     6179                                !-end
     6180                     ZLWFT0_i, ZFLDN0, ZFLUP0, &
     6181                     ZSWFT0_i, ZFSDN0, ZFSUP0, &
     6182                     cloud_cover_sw_s2)
     6183          ENDIF ! ok_3Deffect
     6184#endif
     6185
    59366186       ENDIF ! aerosol_couple
    59376187       itaprad = 0
     
    61576407       d_t_hin(:, :)=0.
    61586408       CALL add_phys_tend(du_gwd_hines, dv_gwd_hines, d_t_hin, dq0, dql0, &
    6159             dqi0, dqbs0,paprs, 'hin', abortphy,flag_inhib_tend,itap,0 &
     6409            dqi0, dqbs0, paprs, 'hin', abortphy,flag_inhib_tend,itap,0 &
    61606410#ifdef ISO
    61616411     &    ,dxt0,dxtl0,dxti0 &
     
    62806530       
    62816531       SELECT CASE(flag_emit)
    6282        CASE(1) ! emission volc H2O dans LMDZ
     6532       CASE(1) ! emission volc H2O in LMDZ
    62836533          DO ieru=1, nErupt
    62846534             IF (year_cur==year_emit_vol(ieru).AND.&
     
    62886538               
    62896539                IF(flag_verbose_strataer) print *,'IN physiq_mod: date=',year_cur,mth_cur,day_cur
    6290                 ! initialisation tendance q emission
     6540                ! initialisation of q tendency emission
    62916541                d_q_emiss(:,:)=0.
    62926542                ! daily injection mass emission - NL
     
    62956545                !
    62966546                CALL STRATEMIT(pdtphys,pdtphys,latitude_deg,longitude_deg,t_seri,&
    6297                      pplay,paprs,tr_seri,&
    6298                      m_H2O_emiss_vol_daily,&
    6299                      xlat_min_vol(ieru),xlat_max_vol(ieru),&
    6300                      xlon_min_vol(ieru),xlon_max_vol(ieru),&
    6301                      altemiss_vol(ieru),&
    6302                      sigma_alt_vol(ieru),1,&
    6303                      1,nAerErupt+1,0)
     6547                    pplay,paprs,tr_seri,&
     6548                    m_H2O_emiss_vol_daily,&
     6549                    xlat_min_vol(ieru),xlat_max_vol(ieru),&
     6550                    xlon_min_vol(ieru),xlon_max_vol(ieru),&
     6551                    altemiss_vol(ieru),sigma_alt_vol(ieru),1,1.,&
     6552                    nAerErupt+1,0)
    63046553               
    63056554                IF(flag_verbose_strataer) print *,'IN physiq_mod: min max d_q_emiss=',&
     
    63156564    ENDIF
    63166565#endif
    6317 
    63186566
    63196567!===============================================================
     
    67547002    t_seri(:,:)=t_seri(:,:)+d_t_ec(:,:)
    67557003
    6756     !=======================================================================
    6757     !   SORTIES
    6758     !=======================================================================
    6759     !
    6760     !IM initialisation + calculs divers diag AMIP2
    6761     !
    6762     include "calcul_divers.h"
    6763     !
    6764     !IM Interpolation sur les niveaux de pression du NMC
    6765     !   -------------------------------------------------
    6766     !
    6767     include "calcul_STDlev.h"
    6768     !
    6769     ! slp sea level pressure derived from Arpege-IFS : CALL ctstar + CALL pppmer
    6770     CALL diag_slp(klon,t_seri,paprs,pplay,pphis,ptstar,pt0,slp)
    6771     !
     7004    !==================================================================
     7005    !--OB water mass fixer for the physics
     7006    !--water profiles are corrected to force mass conservation of water
     7007    !--currently flag is turned off
     7008    !==================================================================
     7009    IF (mass_fixer) THEN
     7010#ifdef ISO
     7011      CALL abort_gcm('physiq 6936','isos pas prevus dans le mass fixer',1)
     7012      ! Camille Risi mai 2024: on attend d'avoir la 4e dimension qui rendra tout plus simple.
     7013#endif
     7014    qql2(:)=0.0
     7015    DO k = 1, klev
     7016      qql2(:)=qql2(:)+(q_seri(:,k)+ql_seri(:,k))*zmasse(:,k)
     7017      IF (nqo >= 3) THEN
     7018        qql2(:)=qql2(:)+qs_seri(:,k)*zmasse(:,k)
     7019      ENDIF
     7020      IF (ok_bs) THEN
     7021        qql2(:)=qql2(:)+qbs_seri(:,k)*zmasse(:,k)
     7022      ENDIF
     7023    ENDDO
     7024
     7025#ifdef CPP_StratAer
     7026    IF (ok_qemiss) THEN
     7027       DO k = 1, klev
     7028          qql1(:) = qql1(:)+d_q_emiss(:,k)*zmasse(:,k)
     7029       ENDDO
     7030    ENDIF
     7031#endif
     7032    IF (ok_qch4) THEN
     7033       DO k = 1, klev
     7034          qql1(:) = qql1(:)+d_q_ch4_dtime(:,k)*zmasse(:,k)
     7035       ENDDO
     7036    ENDIF
     7037   
     7038    DO i = 1, klon
     7039      !--compute ratio of what q+ql should be with conservation to what it is
     7040      IF (ok_bs) THEN
     7041        corrqql=(qql1(i)+(evap(i)-rain_fall(i)-snow_fall(i)-bs_fall(i))*pdtphys)/qql2(i)
     7042      ELSE
     7043        corrqql=(qql1(i)+(evap(i)-rain_fall(i)-snow_fall(i))*pdtphys)/qql2(i)
     7044      ENDIF
     7045      DO k = 1, klev
     7046        q_seri(i,k) =q_seri(i,k)*corrqql
     7047        ql_seri(i,k)=ql_seri(i,k)*corrqql
     7048        IF (nqo >= 3) THEN
     7049          qs_seri(i,k)=qs_seri(i,k)*corrqql
     7050        ENDIF
     7051        IF (ok_bs) THEN
     7052          qbs_seri(i,k)=qbs_seri(i,k)*corrqql
     7053        ENDIF
     7054      ENDDO
     7055    ENDDO
     7056    ENDIF
     7057    !--fin mass fixer
     7058
    67727059    !cc prw  = eau precipitable
    67737060    !   prlw = colonne eau liquide
    67747061    !   prlw = colonne eau solide
    67757062    !   prbsw = colonne neige soufflee
     7063    !   water_budget = non-conservation residual from the LMDZ physics
     7064    !                  (should be equal to machine precision if mass fixer is activated)
    67767065    prw(:) = 0.
    67777066    prlw(:) = 0.
    67787067    prsw(:) = 0.
    67797068    prbsw(:) = 0.
     7069    water_budget(:) = 0.0
    67807070    DO k = 1, klev
    67817071       prw(:)  = prw(:)  + q_seri(:,k)*zmasse(:,k)
    67827072       prlw(:) = prlw(:) + ql_seri(:,k)*zmasse(:,k)
    6783        prsw(:) = prsw(:) + qs_seri(:,k)*zmasse(:,k)
    6784        prbsw(:)= prbsw(:) + qbs_seri(:,k)*zmasse(:,k)
     7073       water_budget(:) = water_budget(:) + (q_seri(:,k)-qx(:,k,ivap)+ql_seri(:,k)-qx(:,k,iliq))*zmasse(:,k)
     7074       IF (nqo >= 3) THEN
     7075         prsw(:) = prsw(:) + qs_seri(:,k)*zmasse(:,k)
     7076         water_budget(:) = water_budget(:) + (qs_seri(:,k)-qx(:,k,isol))*zmasse(:,k)
     7077       ENDIF
     7078       IF (nqo >= 4 .AND. ok_bs) THEN
     7079         prbsw(:)= prbsw(:) + qbs_seri(:,k)*zmasse(:,k)
     7080         water_budget(:) = water_budget(:) + (qbs_seri(:,k)-qx(:,k,ibs))*zmasse(:,k)
     7081       ENDIF
    67857082    ENDDO
    6786 
    6787 #ifdef ISO
    6788       DO i = 1, klon
    6789       do ixt=1,ntraciso
    6790        xtprw(ixt,i) = 0.
    6791        DO k = 1, klev
    6792         xtprw(ixt,i) = xtprw(ixt,i) + &
    6793      &           xt_seri(ixt,i,k)*(paprs(i,k)-paprs(i,k+1))/RG
    6794        ENDDO !DO k = 1, klev
    6795       enddo !do ixt=1,ntraciso
    6796       enddo !DO i = 1, klon
    6797 #endif
     7083    water_budget(:)=water_budget(:)+(rain_fall(:)+snow_fall(:)-evap(:))*pdtphys
     7084    IF (ok_bs) THEN
     7085      water_budget(:)=water_budget(:)+bs_fall(:)*pdtphys
     7086    ENDIF
     7087    ! Camille Risi mai 2024: pour les isotopes, on attend d'avoir la 4e dimension, ça rendra tout plus facile
     7088    ! ces variables sont diagnostiques, donc pas indispensables
     7089
     7090    !=======================================================================
     7091    !   SORTIES
     7092    !=======================================================================
     7093    !
     7094    !IM initialisation + calculs divers diag AMIP2
     7095    !
     7096    include "calcul_divers.h"
     7097    !
     7098    !IM Interpolation sur les niveaux de pression du NMC
     7099    !   -------------------------------------------------
     7100    !
     7101    include "calcul_STDlev.h"
     7102    !
     7103    ! slp sea level pressure derived from Arpege-IFS : CALL ctstar + CALL pppmer
     7104    CALL diag_slp(klon,t_seri,paprs,pplay,pphis,ptstar,pt0,slp)
     7105    !
    67987106    !
    67997107    IF (ANY(type_trac == ['inca','inco'])) THEN
     
    68987206    !IM global posePB      include "write_bilKP_ave.h"
    68997207    !
    6900 
    6901     !--OB mass fixer
    6902     !--profile is corrected to force mass conservation of water
    6903     IF (mass_fixer) THEN
    6904     qql2(:)=0.0
    6905     DO k = 1, klev
    6906       qql2(:)=qql2(:)+(q_seri(:,k)+ql_seri(:,k)+qs_seri(:,k))*zmasse(:,k)
    6907     ENDDO
    6908    
    6909 #ifdef CPP_StratAer
    6910     IF (ok_qemiss) THEN
    6911        DO k = 1, klev
    6912           qql1(:) = qql1(:)+d_q_emiss(:,k)*zmasse(:,k)
    6913        ENDDO
    6914     ENDIF
    6915 #endif
    6916     IF (ok_qch4) THEN
    6917        DO k = 1, klev
    6918           qql1(:) = qql1(:)+d_q_ch4_dtime(:,k)*zmasse(:,k)
    6919        ENDDO
    6920     ENDIF
    6921    
    6922     DO i = 1, klon
    6923       !--compute ratio of what q+ql should be with conservation to what it is
    6924       corrqql=(qql1(i)+(evap(i)-rain_fall(i)-snow_fall(i))*pdtphys)/qql2(i)
    6925       DO k = 1, klev
    6926         q_seri(i,k) =q_seri(i,k)*corrqql
    6927         ql_seri(i,k)=ql_seri(i,k)*corrqql
    6928       ENDDO
    6929     ENDDO
    6930 #ifdef ISO
    6931     do ixt=1,ntraciso
    6932     xtql2(ixt,:)=0.0
    6933     DO k = 1, klev
    6934       xtql2(ixt,:)=xtql2(ixt,:)+(xt_seri(ixt,:,k)+xtl_seri(ixt,:,k)+xts_seri(ixt,:,k))*zmasse(:,k)
    6935     ENDDO
    6936     DO i = 1, klon
    6937       !--compute ratio of what q+ql should be with conservation to what it is
    6938       corrxtql(ixt)=(xtql1(ixt,i)+(xtevap(ixt,i)-xtrain_fall(ixt,i)-xtsnow_fall(ixt,i))*pdtphys)/xtql2(ixt,i)
    6939       DO k = 1, klev
    6940         xt_seri(ixt,i,k) =xt_seri(ixt,i,k)*corrxtql(ixt)
    6941         xtl_seri(ixt,i,k)=xtl_seri(ixt,i,k)*corrxtql(ixt)
    6942       ENDDO
    6943     ENDDO   
    6944     enddo !do ixt=1,ntraciso
    6945 #endif
    6946     ENDIF
    6947     !--fin mass fixer
    6948 
    69497208    ! Sauvegarder les valeurs de t et q a la fin de la physique:
    69507209    !
     
    69627221    xtl_ancien(:,:,:)=xtl_seri(:,:,:)
    69637222    xts_ancien(:,:,:)=xts_seri(:,:,:)
     7223    xtbs_ancien(:,:,:)=xtbs_seri(:,:,:)
    69647224#endif
    69657225    CALL water_int(klon,klev,q_ancien,zmasse,prw_ancien)
     
    70987358         ok_sync, ptconv, read_climoz, clevSTD,          &
    70997359         ptconvth, d_u, d_t, qx, d_qx, zmasse,           &
    7100          flag_aerosol, flag_aerosol_strat, ok_cdnc,t,u1,v1)
     7360         flag_aerosol, flag_aerosol_strat, ok_cdnc,t, u1, v1)
    71017361#endif
    71027362
    71037363#ifndef CPP_XIOS
    7104     CALL write_paramLMDZ_phy(itap,nid_ctesGCM,ok_sync)
    7105 #endif
    7106 
    7107 #endif
    7108 
    7109 ! Pour XIOS : On remet des variables a .false. apres un premier appel
    7110     IF (debut) THEN
    7111 
    7112       IF (using_xios) THEN
    7113         swaero_diag=.FALSE.
    7114         swaerofree_diag=.FALSE.
    7115         dryaod_diag=.FALSE.
    7116         ok_4xCO2atm= .FALSE.
    7117 !       write (lunout,*)'ok_4xCO2atm= ',swaero_diag, swaerofree_diag, dryaod_diag, ok_4xCO2atm
    7118 
    7119         IF (is_master) THEN
    7120           !--setting up swaero_diag to TRUE in XIOS case
    7121           IF (xios_field_is_active("topswad").OR.xios_field_is_active("topswad0").OR. &
    7122              xios_field_is_active("solswad").OR.xios_field_is_active("solswad0").OR. &
    7123              xios_field_is_active("topswai").OR.xios_field_is_active("solswai").OR.  &
    7124                (iflag_rrtm==1.AND.(xios_field_is_active("toplwad").OR.xios_field_is_active("toplwad0").OR. &
    7125                                    xios_field_is_active("sollwad").OR.xios_field_is_active("sollwad0"))))  &
    7126              !!!--for now these fields are not in the XML files so they are omitted
    7127              !!!  xios_field_is_active("toplwai").OR.xios_field_is_active("sollwai") !))) &
    7128              swaero_diag=.TRUE.
    7129 
    7130           !--setting up swaerofree_diag to TRUE in XIOS case
    7131           IF (xios_field_is_active("SWdnSFCcleanclr").OR.xios_field_is_active("SWupSFCcleanclr").OR. &
    7132              xios_field_is_active("SWupTOAcleanclr").OR.xios_field_is_active("rsucsaf").OR.   &
    7133              xios_field_is_active("rsdcsaf") .OR. xios_field_is_active("LWdnSFCcleanclr").OR. &
    7134              xios_field_is_active("LWupTOAcleanclr")) &
    7135              swaerofree_diag=.TRUE.
    7136 
    7137           !--setting up dryaod_diag to TRUE in XIOS case
    7138           DO naero = 1, naero_tot-1
    7139            IF (xios_field_is_active("dryod550_"//name_aero_tau(naero))) dryaod_diag=.TRUE.
    7140           ENDDO
    7141           !
    7142           !--setting up ok_4xCO2atm to TRUE in XIOS case
    7143           IF (xios_field_is_active("rsut4co2").OR.xios_field_is_active("rlut4co2").OR. &
    7144              xios_field_is_active("rsutcs4co2").OR.xios_field_is_active("rlutcs4co2").OR. &
    7145              xios_field_is_active("rsu4co2").OR.xios_field_is_active("rsucs4co2").OR. &
    7146              xios_field_is_active("rsd4co2").OR.xios_field_is_active("rsdcs4co2").OR. &
    7147              xios_field_is_active("rlu4co2").OR.xios_field_is_active("rlucs4co2").OR. &
    7148              xios_field_is_active("rld4co2").OR.xios_field_is_active("rldcs4co2")) &
    7149              ok_4xCO2atm=.TRUE.
    7150         ENDIF
    7151         !$OMP BARRIER
    7152         CALL bcast(swaero_diag)
    7153         CALL bcast(swaerofree_diag)
    7154         CALL bcast(dryaod_diag)
    7155         CALL bcast(ok_4xCO2atm)
    7156 !        write (lunout,*)'ok_4xCO2atm= ',swaero_diag, swaerofree_diag, dryaod_diag, ok_4xCO2atm
    7157       ENDIF !using_xios
    7158     ENDIF
     7364      CALL write_paramLMDZ_phy(itap,nid_ctesGCM,ok_sync)
     7365#endif
     7366
     7367#endif
     7368    ! Petit appelle de sorties pour accompagner le travail sur phyex
     7369    if ( iflag_physiq == 1 ) then
     7370        call output_physiqex(debut,jD_eq,pdtphys,presnivs,paprs,u,v,t,qx,cldfra,0.*t,0.*t,0.*t,pbl_tke,theta)
     7371    endif
    71597372
    71607373    !====================================================================
     
    72007413    ! Disabling calls to the prt_alerte function
    72017414    alert_first_call = .FALSE.
     7415
    72027416   
    72037417    IF (lafin) THEN
     
    72127426         IF (read_climoz >= 1) THEN
    72137427           IF (is_mpi_root) CALL nf95_close(ncid_climoz)
    7214             DEALLOCATE(press_edg_climoz) ! pointer
    7215             DEALLOCATE(press_cen_climoz) ! pointer
     7428            DEALLOCATE(press_edg_climoz)
     7429            DEALLOCATE(press_cen_climoz)
    72167430         ENDIF
    72177431       
    72187432       ENDIF
     7433
     7434       IF (using_xios) THEN
     7435
     7436#ifdef INCA
     7437          IF (type_trac == 'inca') THEN
     7438             IF (is_omp_master .AND. grid_type==unstructured) THEN
     7439                CALL finalize_inca
     7440             ENDIF
     7441          ENDIF
     7442#endif
     7443
     7444          IF (is_omp_master .and. grid_type==unstructured) CALL xios_context_finalize
     7445       ENDIF
     7446
     7447       WRITE(lunout,*) ' physiq fin, nombre de steps ou cvpas = 1 : ', Ncvpaseq1
    72197448       
    7220        IF (using_xios) THEN
    7221          IF (is_omp_master) CALL xios_context_finalize
    7222 
    7223 #ifdef INCA
    7224          if (type_trac == 'inca') then
    7225             IF (is_omp_master .and. grid_type==unstructured) THEN
    7226                CALL finalize_inca
    7227             ENDIF
    7228          endif
    7229 #endif
    7230        ENDIF !using_xios
    7231        WRITE(lunout,*) ' physiq fin, nombre de steps ou cvpas = 1 : ', Ncvpaseq1
    72327449    ENDIF
    72337450
  • LMDZ6/branches/cirrus/libf/phylmdiso/reevap.F90

    r4491 r5202  
    1   SUBROUTINE reevap (klon,klev,iflag_ice_thermo,t_seri,q_seri,ql_seri,qs_seri, &
    2    &         d_t_eva,d_q_eva,d_ql_eva,d_qs_eva &
    3 #ifdef ISO
    4              ,xt_seri,xtl_seri,xts_seri,d_xt_eva,d_xtl_eva,d_xts_eva &
    5 #endif
    6    &     )
     1  SUBROUTINE reevap (klon,klev,iflag_ice_thermo,t_seri,qx, &
     2   &         d_t_eva,d_qx_eva)
    73
    84    ! flag to include modifications to ensure energy conservation (if flag >0)
    95    USE add_phys_tend_mod, only : fl_cor_ebil
    106#ifdef ISO
    11     USE infotrac_phy, ONLY: ntiso   
     7    USE infotrac_phy, ONLY: ntiso,nqtot,ivap,iliq,isol,iqWIsoPha
    128#ifdef ISOVERIF
    139    USE isotopes_verif_mod
     
    2319
    2420    INTEGER klon,klev,iflag_ice_thermo
    25     REAL, DIMENSION(klon,klev), INTENT(in) :: t_seri,q_seri,ql_seri,qs_seri
    26     REAL, DIMENSION(klon,klev), INTENT(out) :: d_t_eva,d_q_eva,d_ql_eva,d_qs_eva
     21    REAL, DIMENSION(klon,klev), INTENT(in) :: t_seri
     22    REAL, DIMENSION(klon,klev,nqtot), INTENT(in) ::     qx
     23    REAL, DIMENSION(klon,klev), INTENT(out) :: d_t_eva
     24    REAL, DIMENSION(klon,klev,nqtot), INTENT(out) ::    d_qx_eva
    2725
    2826    REAL za,zb,zdelta,zlvdcp,zlsdcp
    29     INTEGER i,k
     27    INTEGER i,k,ixt,ivapcur,iliqcur,isolcur   
    3028
    31 #ifdef ISO
    32     REAL, DIMENSION(ntiso,klon,klev), INTENT(in) :: xt_seri,xtl_seri,xts_seri
    33     REAL, DIMENSION(ntiso,klon,klev), INTENT(out) :: d_xt_eva,d_xtl_eva,d_xts_eva
    34     integer ixt
    35 #endif
    3629
    3730    !--------Stochastic Boundary Layer Triggering: ALE_BL--------
     
    4235    !IM 100106 BEG : pouvoir sortir les ctes de la physique
    4336    !
     37    DO ixt = 1, 1+ntiso
    4438    ! Re-evaporer l'eau liquide nuageuse
    4539    !
     40    iliqcur= iqWIsoPha(ixt,iliq)   
     41    ivapcur= iqWIsoPha(ixt,ivap)   
     42    isolcur= iqWIsoPha(ixt,isol)   
    4643!print *,'rrevap ; fl_cor_ebil:',fl_cor_ebil,' iflag_ice_thermo:',iflag_ice_thermo,' RVTMP2',RVTMP2
    4744    DO k = 1, klev  ! re-evaporation de l'eau liquide nuageuse
    48        DO i = 1, klon
    49         if (fl_cor_ebil .GT. 0) then
    50           zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*(q_seri(i,k)+ql_seri(i,k)+qs_seri(i,k)))
    51           zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*(q_seri(i,k)+ql_seri(i,k)+qs_seri(i,k)))
    52         else
    53           zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
    54           !jyg<
    55           !  Attention : Arnaud a propose des formules completement differentes
    56           !                  A verifier !!!
    57           zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*q_seri(i,k))
    58         end if
    59           IF (iflag_ice_thermo .EQ. 0) THEN
    60              zlsdcp=zlvdcp
     45      DO i = 1, klon
     46
     47        IF (ixt == 1) THEN ! water
     48          IF (fl_cor_ebil > 0) THEN
     49            !zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*(q_seri(i,k)+ql_seri(i,k)+qs_seri(i,k)))
     50            !zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*(q_seri(i,k)+ql_seri(i,k)+qs_seri(i,k)))
     51            zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*(qx(i,k,ivapcur)+qx(i,k,iliqcur)+qx(i,k,isolcur)))
     52            zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*(qx(i,k,ivapcur)+qx(i,k,iliqcur)+qx(i,k,isolcur)))
     53          ELSE
     54            zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*qx(i,k,ivapcur))
     55            !jyg<
     56            !  Attention : Arnaud a propose des formules completement differentes
     57            !                  A verifier !!!
     58            zlsdcp=RLSTT/RCPD/(1.0+RVTMP2*qx(i,k,ivapcur))
     59         ENDIF
     60          IF (iflag_ice_thermo == 0) THEN
     61            zlsdcp=zlvdcp
    6162          ENDIF
    6263          !>jyg
     64        ENDIF
     65        IF (iflag_ice_thermo == 0) THEN   
     66           !pas necessaire a priori
    6367
    64           IF (iflag_ice_thermo.eq.0) THEN   
    65              !pas necessaire a priori
     68            zdelta = MAX(0.,SIGN(1.,RTT-t_seri(i,k)))
     69            zdelta = 0.
     70            zb = MAX(0.0,qx(i,k,iliqcur))
     71            IF (ixt == 1) THEN
     72              za = - MAX(0.0,qx(i,k,iliqcur)) &
     73                   * (zlvdcp*(1.-zdelta)+zlsdcp*zdelta)
     74              d_t_eva(i,k) = za
     75            ENDIF
     76            d_qx_eva(i,k,ivapcur)  = zb
     77            d_qx_eva(i,k,iliqcur) = -qx(i,k,iliqcur)
     78            d_qx_eva(i,k,isolcur) = 0.
    6679
    67              zdelta = MAX(0.,SIGN(1.,RTT-t_seri(i,k)))
    68   zdelta = 0.
    69              zb = MAX(0.0,ql_seri(i,k))
    70              za = - MAX(0.0,ql_seri(i,k)) &
    71                   * (zlvdcp*(1.-zdelta)+zlsdcp*zdelta)
    72              d_t_eva(i,k) = za
    73              d_q_eva(i,k) = zb
    74              d_ql_eva(i,k) = -ql_seri(i,k)
    75              d_qs_eva(i,k) = 0.
    76 
    77 #ifdef ISO
    78          do ixt=1,ntiso
    79             zb = MAX(0.0,xtl_seri(ixt,i,k))
    80             d_xt_eva(ixt,i,k) = zb
    81             d_xtl_eva(ixt,i,k) = -xtl_seri(ixt,i,k)
    82             d_xts_eva(ixt,i,k) = 0.
    83          enddo
    84 #ifdef ISOVERIF
    85       do ixt=1,ntiso
    86         call iso_verif_noNaN(xt_seri(ixt,i,k), &
    87      &     'reevap 2417: apres evap tot')
    88       enddo
    89       if (iso_eau.gt.0) then
    90               call iso_verif_egalite_choix( &
    91      &           xt_seri(iso_eau,i,k),q_seri(i,k), &
    92      &          'reevap 1891+, après reevap totale',errmax,errmaxrel)
    93               call iso_verif_egalite_choix( &
    94      &           xtl_seri(iso_eau,i,k),ql_seri(i,k), &
    95      &          'reevap 2209+, après reevap totale',errmax,errmaxrel)
    96        endif !if (iso_eau.gt.0) then       
    97       if ((iso_HDO.gt.0).and.(iso_O18.gt.0)) then 
    98             if (q_seri(i,k).gt.ridicule) then 
    99                if (iso_verif_o18_aberrant_nostop( &
    100      &              xt_seri(iso_HDO,i,k)/q_seri(i,k), &
    101      &              xt_seri(iso_O18,i,k)/q_seri(i,k), &
    102      &              'reevap 2315: apres reevap totale').eq.1) then
    103                   write(*,*) 'i,k,q_seri(i,k)=',i,k,q_seri(i,k)
    104                   write(*,*) 'd_q_eva(i,k)=',d_q_eva(i,k)
    105                   write(*,*) 'deltaD(d_q_eva(i,k))=',deltaD(d_xt_eva(iso_HDO,i,k)/d_q_eva(i,k))
    106                   write(*,*) 'deltaO18(d_q_eva(i,k))=',deltaO(d_xt_eva(iso_O18,i,k)/d_q_eva(i,k))
    107                   stop
    108               endif !  if (iso_verif_o18_aberrant_nostop
    109             endif !if (q_seri(i,k).gt.errmax) then   
    110         endif !if ((iso_HDO.gt.0).and.(iso_O18.gt.0)) then
    111 #ifdef ISOTRAC     
    112              call iso_verif_traceur(xt_seri(1,i,k), &
    113      &           'reevap 2165a')
    114              call iso_verif_traceur_pbidouille(xt_seri(1,i,k), &
    115      &           'reevap 2165b')
    116 #endif               
    117          
    118 #endif           
    119 #endif
    120 
    121           ELSE
    122 
     80        ELSE
     81             
    12382             !CR: on r\'e-\'evapore eau liquide et glace
    12483
     
    12786             !        za = - MAX(0.0,ql_seri(i,k)) &
    12887             !             * (zlvdcp*(1.-zdelta)+zlsdcp*zdelta)
    129              zb = MAX(0.0,ql_seri(i,k)+qs_seri(i,k))
    130              za = - MAX(0.0,ql_seri(i,k))*zlvdcp &
    131                   - MAX(0.0,qs_seri(i,k))*zlsdcp
    132              d_t_eva(i,k) = za
    133              d_q_eva(i,k) = zb
    134              d_ql_eva(i,k) = -ql_seri(i,k)
    135              d_qs_eva(i,k) = -qs_seri(i,k)
     88            IF (ixt == 1) THEN
     89              za = - MAX(0.0,qx(i,k,iliqcur))*zlvdcp &
     90                   - MAX(0.0,qx(i,k,iliqcur))*zlsdcp
     91              d_t_eva(i,k) = za
     92            ENDIF
     93            !zb = MAX(0.0,ql_seri(i,k)+qs_seri(i,k))
     94            !d_q_eva(i,k) = zb
     95            !d_ql_eva(i,k) = -ql_seri(i,k)
     96            !d_qs_eva(i,k) = -qs_seri(i,k)
    13697
    137 #ifdef ISO
    138          do ixt=1,ntiso
    139             zb = MAX(0.0,xtl_seri(ixt,i,k)+xts_seri(ixt,i,k))
    140             d_xt_eva(ixt,i,k) = zb
    141             d_xtl_eva(ixt,i,k) = -xtl_seri(ixt,i,k)
    142             d_xts_eva(ixt,i,k) = -xts_seri(ixt,i,k)
    143          enddo
     98            zb = MAX(0.0,qx(i,k,iliqcur)+qx(i,k,isolcur))
     99            d_qx_eva(i,k,ivapcur) = zb
     100            d_qx_eva(i,k,iliqcur) = -qx(i,k,iliqcur)
     101            d_qx_eva(i,k,isolcur) = -qx(i,k,isolcur)
     102        ENDIF
    144103
    145 #ifdef ISOVERIF
    146       do ixt=1,ntiso
    147       call iso_verif_noNaN(xt_seri(ixt,i,k), &
    148      &     'reevap 2417: apres evap tot')
    149       enddo
    150       if (iso_eau.gt.0) then
    151               call iso_verif_egalite_choix( &
    152      &           xt_seri(iso_eau,i,k),q_seri(i,k), &
    153      &          'reevap 1891, après réévap totale',errmax,errmaxrel)
    154               call iso_verif_egalite_choix( &
    155      &           xtl_seri(iso_eau,i,k),ql_seri(i,k), &
    156      &          'reevap 2209, après réévap totale',errmax,errmaxrel)
    157               call iso_verif_egalite_choix( &
    158      &           xts_seri(iso_eau,i,k),qs_seri(i,k), &
    159      &          'reevap 2209b, après réévap totale',errmax,errmaxrel)
    160        endif !if (iso_eau.gt.0) then
    161      
    162       if ((iso_HDO.gt.0).and.(iso_O18.gt.0)) then 
    163             if (q_seri(i,k).gt.ridicule) then 
    164                if (iso_verif_o18_aberrant_nostop( &
    165      &              xt_seri(iso_HDO,i,k)/q_seri(i,k), &
    166      &              xt_seri(iso_O18,i,k)/q_seri(i,k), &
    167      &              'reevap 2408: apres reevap totale').eq.1) then
    168                   write(*,*) 'i,k,q_seri(i,k)=',i,k,q_seri(i,k)                       
    169                   stop
    170               endif !  if (iso_verif_o18_aberrant_nostop
    171             endif !if (q_seri(i,k).gt.errmax) then   
    172         endif !if ((iso_HDO.gt.0).and.(iso_O18.gt.0)) then       
    173 #ifdef ISOTRAC     
    174              call iso_verif_traceur(xt_seri(1,i,k), &
    175      &           'reevap 2165c')
    176              call iso_verif_traceur_pbidouille(xt_seri(1,i,k), &
    177      &           'reevap 2165d')
    178 #endif                 
    179 #endif                 
    180 #endif
    181104
    182           ENDIF
     105      ENDDO
     106    ENDDO
    183107
    184        ENDDO
    185     ENDDO
     108    ENDDO ! DO ixt = 1, 1+niso*(nzone +1)
    186109
    187110RETURN
Note: See TracChangeset for help on using the changeset viewer.