Ignore:
Timestamp:
Nov 28, 2014, 4:36:29 PM (10 years ago)
Author:
Laurent Fairhead
Message:

Merged trunk changes -r2070:2158 into testing branch. Compilation problems introduced by revision r2155 have been corrected by hand

Location:
LMDZ5/branches/testing
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • LMDZ5/branches/testing

  • LMDZ5/branches/testing/libf/phylmd/phytrac_mod.F90

    r2056 r2160  
    5454CONTAINS
    5555
    56 SUBROUTINE phytrac(                                 &
    57      nstep,     julien,   gmtime,   debutphy,       &
    58      lafin,     pdtphys,  u, v,     t_seri,         &
    59      paprs,     pplay,    pmfu,     pmfd,           &
    60      pen_u,     pde_u,    pen_d,    pde_d,          &
    61      cdragh,    coefh,    fm_therm, entr_therm,     &
    62      yu1,       yv1,      ftsol,    pctsrf,         &
    63      ustar,     u10m,      v10m,                    &
    64      wstar,     ale_bl,      ale_wake,              &
    65      xlat,      xlon,                               &
    66      frac_impa,frac_nucl,beta_fisrt,beta_v1,        &
    67      presnivs,  pphis,    pphi,     albsol,         &
    68      sh,        rh,       cldfra,   rneb,           &
    69      diafra,    cldliq,   itop_con, ibas_con,       &
    70      pmflxr,    pmflxs,   prfl,     psfl,           &
    71      da,        phi,      mp,       upwd,           &
    72      phi2,      d1a,      dam,      sij, wght_cvfd, &   ! RomP +RL
    73      wdtrainA,  wdtrainM, sigd,     clw, elij,      &   ! RomP
    74      evap,      ep,       epmlmMm,  eplaMm,         &   ! RomP
    75      dnwd,      aerosol_couple,     flxmass_w,      &
    76      tau_aero,  piz_aero,  cg_aero, ccm,            &
    77      rfname,                                        &
    78      d_tr_dyn,                                      &   ! RomP
    79      tr_seri)         
    80 !
    81 !======================================================================
    82 ! Auteur(s) FH
    83 ! Objet: Moniteur general des tendances traceurs
    84 ! Modification R. Pilon 01 janvier 2012 transport+scavenging KE scheme : cvltr
    85 ! Modification R. Pilon 10 octobre 2012 large scale scavenging incloud_scav + bc_scav
    86 !======================================================================
    87 
    88   USE ioipsl
    89   USE phys_cal_mod, only : hour
    90   USE dimphy
    91   USE infotrac
    92   USE mod_grid_phy_lmdz
    93   USE mod_phys_lmdz_para
    94   USE comgeomphy
    95   USE iophy
    96   USE traclmdz_mod
    97   USE tracinca_mod
    98   USE tracreprobus_mod
    99   USE control_mod
    100   USE indice_sol_mod
    101 
    102   IMPLICIT NONE
    103 
    104   INCLUDE "YOMCST.h"
    105   INCLUDE "dimensions.h"
    106   INCLUDE "clesphys.h"
    107   INCLUDE "temps.h"
    108   INCLUDE "paramet.h"
    109   INCLUDE "thermcell.h"
    110   INCLUDE "iniprint.h"
    111 !==========================================================================
    112 !                   -- ARGUMENT DESCRIPTION --
    113 !==========================================================================
    114 
    115 ! Input arguments
    116 !----------------
    117 !Configuration grille,temps:
    118   INTEGER,INTENT(IN) :: nstep      ! Appel physique
    119   INTEGER,INTENT(IN) :: julien     ! Jour julien
    120   REAL,INTENT(IN)    :: gmtime     ! Heure courante
    121   REAL,INTENT(IN)    :: pdtphys    ! Pas d'integration pour la physique (seconde)
    122   LOGICAL,INTENT(IN) :: debutphy   ! le flag de l'initialisation de la physique
    123   LOGICAL,INTENT(IN) :: lafin      ! le flag de la fin de la physique
    124  
    125   REAL,DIMENSION(klon),INTENT(IN) :: xlat    ! latitudes pour chaque point
    126   REAL,DIMENSION(klon),INTENT(IN) :: xlon    ! longitudes pour chaque point
    127 !
    128 !Physique:
    129 !--------
    130   REAL,DIMENSION(klon,klev),INTENT(IN)   :: t_seri  ! Temperature
    131   REAL,DIMENSION(klon,klev),INTENT(IN)   :: u       ! variable not used
    132   REAL,DIMENSION(klon,klev),INTENT(IN)   :: v       ! variable not used
    133   REAL,DIMENSION(klon,klev),INTENT(IN)   :: sh      ! humidite specifique
    134   REAL,DIMENSION(klon,klev),INTENT(IN)   :: rh      ! humidite relative
    135   REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs   ! pression pour chaque inter-couche (en Pa)
    136   REAL,DIMENSION(klon,klev),INTENT(IN)   :: pplay   ! pression pour le mileu de chaque couche (en Pa)
    137   REAL,DIMENSION(klon,klev),INTENT(IN)   :: pphi    ! geopotentiel
    138   REAL,DIMENSION(klon),INTENT(IN)        :: pphis
    139   REAL,DIMENSION(klev),INTENT(IN)        :: presnivs
    140   REAL,DIMENSION(klon,klev),INTENT(IN)   :: cldliq  ! eau liquide nuageuse
    141   REAL,DIMENSION(klon,klev),INTENT(IN)   :: cldfra  ! fraction nuageuse (tous les nuages)
    142   REAL,DIMENSION(klon,klev),INTENT(IN)   :: diafra  ! fraction nuageuse (convection ou stratus artificiels)
    143   REAL,DIMENSION(klon,klev),INTENT(IN)   :: rneb    ! fraction nuageuse (grande echelle)
    144 !
    145   REAL                                   :: ql_incl ! contenu en eau liquide nuageuse dans le nuage ! ql_incl=oliq/rneb
    146   REAL,DIMENSION(klon,klev),INTENT(IN)   :: beta_fisrt ! taux de conversion de l'eau cond (de fisrtilp)
    147   REAL,DIMENSION(klon,klev),INTENT(out)  :: beta_v1    ! -- (originale version)
    148 
    149 !
    150   INTEGER,DIMENSION(klon),INTENT(IN)     :: itop_con
    151   INTEGER,DIMENSION(klon),INTENT(IN)     :: ibas_con
    152   REAL,DIMENSION(klon),INTENT(IN)        :: albsol  ! albedo surface
    153 !
    154 !Dynamique
    155 !--------
    156   REAL,DIMENSION(klon,klev,nbtr),INTENT(IN)    :: d_tr_dyn
    157 !
    158 !Convection:
    159 !----------
    160   REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfu  ! flux de masse dans le panache montant
    161   REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfd  ! flux de masse dans le panache descendant
    162   REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_u ! flux entraine dans le panache montant
    163   REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_u ! flux detraine dans le panache montant
    164   REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_d ! flux entraine dans le panache descendant
    165   REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_d ! flux detraine dans le panache descendant
    166 
    167 !...Tiedke     
    168   REAL,DIMENSION(klon,klev+1),INTENT(IN)   :: pmflxr, pmflxs ! Flux precipitant de pluie, neige aux interfaces [convection]
    169   REAL,DIMENSION(klon,klev+1),INTENT(IN)   :: prfl, psfl ! Flux precipitant de pluie, neige aux interfaces [large-scale]
    170 
    171   LOGICAL,INTENT(IN)                       :: aerosol_couple
    172   REAL,DIMENSION(klon,klev),INTENT(IN)     :: flxmass_w
    173   REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: tau_aero
    174   REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: piz_aero
    175   REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: cg_aero
    176   CHARACTER(len=4),DIMENSION(9),INTENT(IN) :: rfname
    177   REAL,DIMENSION(klon,klev,2),INTENT(IN)   :: ccm
    178 !... K.Emanuel
    179   REAL,DIMENSION(klon,klev),INTENT(IN)     :: da
    180   REAL,DIMENSION(klon,klev,klev),INTENT(IN):: phi
    181 ! RomP >>>
    182   REAL,DIMENSION(klon,klev),INTENT(IN)      :: d1a,dam
    183   REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: phi2
    184 !
    185   REAL,DIMENSION(klon,klev),INTENT(IN)      :: wdtrainA
    186   REAL,DIMENSION(klon,klev),INTENT(IN)      :: wdtrainM
    187   REAL,DIMENSION(klon),INTENT(IN)           :: sigd
    188 ! ---- RomP flux entraine, detraine et precipitant kerry Emanuel
    189   REAL,DIMENSION(klon,klev),INTENT(IN)      :: evap
    190   REAL,DIMENSION(klon,klev),INTENT(IN)      :: ep
    191   REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: sij
    192   REAL,DIMENSION(klon,klev),INTENT(IN)      :: wght_cvfd          !RL
    193   REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: elij
    194   REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: epmlmMm
    195   REAL,DIMENSION(klon,klev),INTENT(IN)      :: eplaMm
    196   REAL,DIMENSION(klon,klev),INTENT(IN)      :: clw
    197 ! RomP <<<
    198 
    199 !
    200   REAL,DIMENSION(klon,klev),INTENT(IN)     :: mp
    201   REAL,DIMENSION(klon,klev),INTENT(IN)     :: upwd      ! saturated updraft mass flux
    202   REAL,DIMENSION(klon,klev),INTENT(IN)     :: dnwd      ! saturated downdraft mass flux
    203 !
    204 !Thermiques:
    205 !----------
    206   REAL,DIMENSION(klon,klev+1),INTENT(IN)   :: fm_therm
    207   REAL,DIMENSION(klon,klev),INTENT(IN)     :: entr_therm
    208 !
    209 !Couche limite:
    210 !--------------
    211 !
    212 !
    213   REAL,DIMENSION(:),INTENT(IN)   :: cdragh          ! (klon) coeff drag pour T et Q
    214   REAL,DIMENSION(:,:),INTENT(IN) :: coefh           ! (klon,klev) coeff melange CL (m**2/s)
    215   REAL,DIMENSION(:),INTENT(IN)   :: ustar,u10m,v10m ! (klon) u* & vent a 10m (m/s)
    216   REAL,DIMENSION(:),INTENT(IN)   :: wstar,ale_bl,ale_wake ! (klon) w* and Avail. Lifting Ener.
    217   REAL,DIMENSION(:),INTENT(IN)   :: yu1             ! (klon) vents au premier niveau
    218   REAL,DIMENSION(:),INTENT(IN)   :: yv1             ! (klon) vents au premier niveau
    219 
    220 !
    221 !Lessivage:
    222 !----------
    223 !
    224 ! pour le ON-LINE
    225 !
    226   REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_impa ! fraction d'aerosols non impactes
    227   REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_nucl ! fraction d'aerosols non nuclees
    228 
    229 ! Arguments necessaires pour les sources et puits de traceur:
    230   REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: ftsol  ! Temperature du sol (surf)(Kelvin)
    231   REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: pctsrf ! Pourcentage de sol (nature du sol)
    232 
    233 
    234 ! Output argument
    235 !----------------
    236   REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA]
    237   REAL,DIMENSION(klon,klev)                    :: sourceBE
    238 !=======================================================================================
    239 !                        -- LOCAL VARIABLES --
    240 !=======================================================================================
    241 
    242   INTEGER :: i, k, it
    243   INTEGER :: nsplit
    244 
    245 !Sources et Reservoirs de traceurs (ex:Radon):
    246 !--------------------------------------------
    247 !
    248   REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: source  ! a voir lorsque le flux de surface est prescrit
     56  SUBROUTINE phytrac(                                 &
     57       nstep,     julien,   gmtime,   debutphy,       &
     58       lafin,     pdtphys,  u, v,     t_seri,         &
     59       paprs,     pplay,    pmfu,     pmfd,           &
     60       pen_u,     pde_u,    pen_d,    pde_d,          &
     61       cdragh,    coefh,    fm_therm, entr_therm,     &
     62       yu1,       yv1,      ftsol,    pctsrf,         &
     63       ustar,     u10m,      v10m,                    &
     64       wstar,     ale_bl,      ale_wake,              &
     65       xlat,      xlon,                               &
     66       frac_impa,frac_nucl,beta_fisrt,beta_v1,        &
     67       presnivs,  pphis,    pphi,     albsol,         &
     68       sh,        rh,       cldfra,   rneb,           &
     69       diafra,    cldliq,   itop_con, ibas_con,       &
     70       pmflxr,    pmflxs,   prfl,     psfl,           &
     71       da,        phi,      mp,       upwd,           &
     72       phi2,      d1a,      dam,      sij, wght_cvfd, &   ! RomP +RL
     73       wdtrainA,  wdtrainM, sigd,     clw, elij,      &   ! RomP
     74       evap,      ep,       epmlmMm,  eplaMm,         &   ! RomP
     75       dnwd,      aerosol_couple,     flxmass_w,      &
     76       tau_aero,  piz_aero,  cg_aero, ccm,            &
     77       rfname,                                        &
     78       d_tr_dyn,                                      &   ! RomP
     79       tr_seri)         
     80    !
     81    !======================================================================
     82    ! Auteur(s) FH
     83    ! Objet: Moniteur general des tendances traceurs
     84    ! Modification R. Pilon 01 janvier 2012 transport+scavenging KE scheme : cvltr
     85    ! Modification R. Pilon 10 octobre 2012 large scale scavenging incloud_scav + bc_scav
     86    !======================================================================
     87
     88    USE ioipsl
     89    USE phys_cal_mod, only : hour
     90    USE dimphy
     91    USE infotrac
     92    USE mod_grid_phy_lmdz
     93    USE mod_phys_lmdz_para
     94    USE comgeomphy
     95    USE iophy
     96    USE traclmdz_mod
     97    USE tracinca_mod
     98    USE tracreprobus_mod
     99    USE control_mod
     100    USE indice_sol_mod
     101
     102    USE mod_phys_lmdz_mpi_data, ONLY :  is_mpi_root
     103
     104    IMPLICIT NONE
     105
     106    INCLUDE "YOMCST.h"
     107    INCLUDE "dimensions.h"
     108    INCLUDE "clesphys.h"
     109    INCLUDE "temps.h"
     110    INCLUDE "paramet.h"
     111    INCLUDE "thermcell.h"
     112    INCLUDE "iniprint.h"
     113    !==========================================================================
     114    !                   -- ARGUMENT DESCRIPTION --
     115    !==========================================================================
     116
     117    ! Input arguments
     118    !----------------
     119    !Configuration grille,temps:
     120    INTEGER,INTENT(IN) :: nstep      ! Appel physique
     121    INTEGER,INTENT(IN) :: julien     ! Jour julien
     122    REAL,INTENT(IN)    :: gmtime     ! Heure courante
     123    REAL,INTENT(IN)    :: pdtphys    ! Pas d'integration pour la physique (seconde)
     124    LOGICAL,INTENT(IN) :: debutphy   ! le flag de l'initialisation de la physique
     125    LOGICAL,INTENT(IN) :: lafin      ! le flag de la fin de la physique
     126
     127    REAL,DIMENSION(klon),INTENT(IN) :: xlat    ! latitudes pour chaque point
     128    REAL,DIMENSION(klon),INTENT(IN) :: xlon    ! longitudes pour chaque point
     129    !
     130    !Physique:
     131    !--------
     132    REAL,DIMENSION(klon,klev),INTENT(IN)   :: t_seri  ! Temperature
     133    REAL,DIMENSION(klon,klev),INTENT(IN)   :: u       ! variable not used
     134    REAL,DIMENSION(klon,klev),INTENT(IN)   :: v       ! variable not used
     135    REAL,DIMENSION(klon,klev),INTENT(IN)   :: sh      ! humidite specifique
     136    REAL,DIMENSION(klon,klev),INTENT(IN)   :: rh      ! humidite relative
     137    REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs   ! pression pour chaque inter-couche (en Pa)
     138    REAL,DIMENSION(klon,klev),INTENT(IN)   :: pplay   ! pression pour le mileu de chaque couche (en Pa)
     139    REAL,DIMENSION(klon,klev),INTENT(IN)   :: pphi    ! geopotentiel
     140    REAL,DIMENSION(klon),INTENT(IN)        :: pphis
     141    REAL,DIMENSION(klev),INTENT(IN)        :: presnivs
     142    REAL,DIMENSION(klon,klev),INTENT(IN)   :: cldliq  ! eau liquide nuageuse
     143    REAL,DIMENSION(klon,klev),INTENT(IN)   :: cldfra  ! fraction nuageuse (tous les nuages)
     144    REAL,DIMENSION(klon,klev),INTENT(IN)   :: diafra  ! fraction nuageuse (convection ou stratus artificiels)
     145    REAL,DIMENSION(klon,klev),INTENT(IN)   :: rneb    ! fraction nuageuse (grande echelle)
     146    !
     147    REAL                                   :: ql_incl ! contenu en eau liquide nuageuse dans le nuage ! ql_incl=oliq/rneb
     148    REAL,DIMENSION(klon,klev),INTENT(IN)   :: beta_fisrt ! taux de conversion de l'eau cond (de fisrtilp)
     149    REAL,DIMENSION(klon,klev),INTENT(out)  :: beta_v1    ! -- (originale version)
     150
     151    !
     152    INTEGER,DIMENSION(klon),INTENT(IN)     :: itop_con
     153    INTEGER,DIMENSION(klon),INTENT(IN)     :: ibas_con
     154    REAL,DIMENSION(klon),INTENT(IN)        :: albsol  ! albedo surface
     155    !
     156    !Dynamique
     157    !--------
     158    REAL,DIMENSION(klon,klev,nbtr),INTENT(IN)    :: d_tr_dyn
     159    !
     160    !Convection:
     161    !----------
     162    REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfu  ! flux de masse dans le panache montant
     163    REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfd  ! flux de masse dans le panache descendant
     164    REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_u ! flux entraine dans le panache montant
     165    REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_u ! flux detraine dans le panache montant
     166    REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_d ! flux entraine dans le panache descendant
     167    REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_d ! flux detraine dans le panache descendant
     168
     169    !...Tiedke     
     170    REAL,DIMENSION(klon,klev+1),INTENT(IN)   :: pmflxr, pmflxs ! Flux precipitant de pluie, neige aux interfaces [convection]
     171    REAL,DIMENSION(klon,klev+1),INTENT(IN)   :: prfl, psfl ! Flux precipitant de pluie, neige aux interfaces [large-scale]
     172
     173    LOGICAL,INTENT(IN)                       :: aerosol_couple
     174    REAL,DIMENSION(klon,klev),INTENT(IN)     :: flxmass_w
     175    REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: tau_aero
     176    REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: piz_aero
     177    REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: cg_aero
     178    CHARACTER(len=4),DIMENSION(9),INTENT(IN) :: rfname
     179    REAL,DIMENSION(klon,klev,2),INTENT(IN)   :: ccm
     180    !... K.Emanuel
     181    REAL,DIMENSION(klon,klev),INTENT(IN)     :: da
     182    REAL,DIMENSION(klon,klev,klev),INTENT(IN):: phi
     183    ! RomP >>>
     184    REAL,DIMENSION(klon,klev),INTENT(IN)      :: d1a,dam
     185    REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: phi2
     186    !
     187    REAL,DIMENSION(klon,klev),INTENT(IN)      :: wdtrainA
     188    REAL,DIMENSION(klon,klev),INTENT(IN)      :: wdtrainM
     189    REAL,DIMENSION(klon),INTENT(IN)           :: sigd
     190    ! ---- RomP flux entraine, detraine et precipitant kerry Emanuel
     191    REAL,DIMENSION(klon,klev),INTENT(IN)      :: evap
     192    REAL,DIMENSION(klon,klev),INTENT(IN)      :: ep
     193    REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: sij
     194    REAL,DIMENSION(klon,klev),INTENT(IN)      :: wght_cvfd          !RL
     195    REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: elij
     196    REAL,DIMENSION(klon,klev,klev),INTENT(IN) :: epmlmMm
     197    REAL,DIMENSION(klon,klev),INTENT(IN)      :: eplaMm
     198    REAL,DIMENSION(klon,klev),INTENT(IN)      :: clw
     199    ! RomP <<<
     200
     201    !
     202    REAL,DIMENSION(klon,klev),INTENT(IN)     :: mp
     203    REAL,DIMENSION(klon,klev),INTENT(IN)     :: upwd      ! saturated updraft mass flux
     204    REAL,DIMENSION(klon,klev),INTENT(IN)     :: dnwd      ! saturated downdraft mass flux
     205    !
     206    !Thermiques:
     207    !----------
     208    REAL,DIMENSION(klon,klev+1),INTENT(IN)   :: fm_therm
     209    REAL,DIMENSION(klon,klev),INTENT(IN)     :: entr_therm
     210    !
     211    !Couche limite:
     212    !--------------
     213    !
     214    !
     215    REAL,DIMENSION(:),INTENT(IN)   :: cdragh          ! (klon) coeff drag pour T et Q
     216    REAL,DIMENSION(:,:),INTENT(IN) :: coefh           ! (klon,klev) coeff melange CL (m**2/s)
     217    REAL,DIMENSION(:),INTENT(IN)   :: ustar,u10m,v10m ! (klon) u* & vent a 10m (m/s)
     218    REAL,DIMENSION(:),INTENT(IN)   :: wstar,ale_bl,ale_wake ! (klon) w* and Avail. Lifting Ener.
     219    REAL,DIMENSION(:),INTENT(IN)   :: yu1             ! (klon) vents au premier niveau
     220    REAL,DIMENSION(:),INTENT(IN)   :: yv1             ! (klon) vents au premier niveau
     221
     222    !
     223    !Lessivage:
     224    !----------
     225    !
     226    REAL, DIMENSION(:), ALLOCATABLE, SAVE :: ccntrAA
     227    REAL, DIMENSION(:), ALLOCATABLE, SAVE :: ccntrENV
     228    REAL, DIMENSION(:), ALLOCATABLE, SAVE :: coefcoli
     229    LOGICAL, DIMENSION(:), ALLOCATABLE, SAVE :: flag_cvltr
     230!$OMP THREADPRIVATE(ccntrAA,ccntrENV,coefcoli,flag_cvltr)
     231    REAL, DIMENSION(klon,klev) :: ccntrAA_3d
     232    REAL, DIMENSION(klon,klev) :: ccntrENV_3d
     233    REAL, DIMENSION(klon,klev) :: coefcoli_3d
     234    !
     235    ! pour le ON-LINE
     236    !
     237    REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_impa ! fraction d'aerosols non impactes
     238    REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_nucl ! fraction d'aerosols non nuclees
     239
     240    ! Arguments necessaires pour les sources et puits de traceur:
     241    REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: ftsol  ! Temperature du sol (surf)(Kelvin)
     242    REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: pctsrf ! Pourcentage de sol (nature du sol)
     243
     244
     245    ! Output argument
     246    !----------------
     247    REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA]
     248    REAL,DIMENSION(klon,klev)                    :: sourceBE
     249    !=======================================================================================
     250    !                        -- LOCAL VARIABLES --
     251    !=======================================================================================
     252
     253    INTEGER :: i, k, it
     254    INTEGER :: nsplit
     255
     256    !Sources et Reservoirs de traceurs (ex:Radon):
     257    !--------------------------------------------
     258    !
     259    REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: source  ! a voir lorsque le flux de surface est prescrit
    249260!$OMP THREADPRIVATE(source)
    250261
    251 !
    252 !Entrees/Sorties: (cf ini_histrac.h et write_histrac.h) 
    253 !---------------
    254   INTEGER                   :: iiq, ierr
    255   INTEGER                   :: nhori, nvert
    256   REAL                      :: zsto, zout, zjulian
    257   INTEGER,SAVE              :: nid_tra     ! pointe vers le fichier histrac.nc         
     262    !
     263    !Entrees/Sorties: (cf ini_histrac.h et write_histrac.h) 
     264    !---------------
     265    INTEGER                   :: iiq, ierr
     266    INTEGER                   :: nhori, nvert
     267    REAL                      :: zsto, zout, zjulian
     268    INTEGER,SAVE              :: nid_tra     ! pointe vers le fichier histrac.nc         
    258269!$OMP THREADPRIVATE(nid_tra)
    259   REAL,DIMENSION(klon)      :: zx_tmp_fi2d ! variable temporaire grille physique
    260   INTEGER                   :: itau_w      ! pas de temps ecriture = nstep + itau_phy
    261   LOGICAL,PARAMETER         :: ok_sync=.TRUE.
    262 
    263 !
    264 ! Nature du traceur
    265 !------------------
    266   LOGICAL,DIMENSION(:),ALLOCATABLE,SAVE :: aerosol  ! aerosol(it) = true  => aerosol => lessivage
     270    REAL,DIMENSION(klon)      :: zx_tmp_fi2d ! variable temporaire grille physique
     271    INTEGER                   :: itau_w      ! pas de temps ecriture = nstep + itau_phy
     272    LOGICAL,PARAMETER         :: ok_sync=.TRUE.
     273
     274    !
     275    ! Nature du traceur
     276    !------------------
     277    LOGICAL,DIMENSION(:),ALLOCATABLE,SAVE :: aerosol  ! aerosol(it) = true  => aerosol => lessivage
    267278!$OMP THREADPRIVATE(aerosol)                        ! aerosol(it) = false => gaz
    268   REAL,DIMENSION(klon,klev)             :: delp     ! epaisseur de couche (Pa)
    269 !
    270 ! Tendances de traceurs (Td) et flux de traceurs:
    271 !------------------------
    272   REAL,DIMENSION(klon,klev)      :: d_tr     ! Td dans l'atmosphere
    273   REAL,DIMENSION(klon,klev)      :: Mint
    274   REAL,DIMENSION(klon,klev,nbtr) :: zmfd1a
    275   REAL,DIMENSION(klon,klev,nbtr) :: zmfdam
    276   REAL,DIMENSION(klon,klev,nbtr) :: zmfphi2
    277 ! Physique
    278 !----------
    279   REAL,DIMENSION(klon,klev,nbtr) :: flestottr ! flux de lessivage dans chaque couche
    280   REAL,DIMENSION(klon,klev)      :: zmasse    ! densité atmosphérique Kg/m2
    281   REAL,DIMENSION(klon,klev)      :: ztra_th
    282 !PhH
    283   REAL,DIMENSION(klon,klev)      :: zrho
    284   REAL,DIMENSION(klon,klev)      :: zdz
    285   REAL                           :: evaplsc,dx,beta ! variable pour lessivage Genthon
    286   REAL,DIMENSION(klon)           :: his_dh          ! ---
    287 ! in-cloud scav variables
    288   REAL           :: ql_incloud_ref     ! ref value of in-cloud condensed water content
    289  
    290 !Controles:
    291 !---------
    292   INTEGER,SAVE :: iflag_vdf_trac,iflag_con_trac,iflag_the_trac
    293   INTEGER,SAVE  :: iflag_con_trac_omp, iflag_vdf_trac_omp,iflag_the_trac_omp
     279    REAL,DIMENSION(klon,klev)             :: delp     ! epaisseur de couche (Pa)
     280    !
     281    ! Tendances de traceurs (Td) et flux de traceurs:
     282    !------------------------
     283    REAL,DIMENSION(klon,klev)      :: d_tr     ! Td dans l'atmosphere
     284    REAL,DIMENSION(klon,klev)      :: Mint
     285    REAL,DIMENSION(klon,klev,nbtr) :: zmfd1a
     286    REAL,DIMENSION(klon,klev,nbtr) :: zmfdam
     287    REAL,DIMENSION(klon,klev,nbtr) :: zmfphi2
     288    ! Physique
     289    !----------
     290    REAL,DIMENSION(klon,klev,nbtr) :: flestottr ! flux de lessivage dans chaque couche
     291    REAL,DIMENSION(klon,klev)      :: zmasse    ! densité atmosphérique Kg/m2
     292    REAL,DIMENSION(klon,klev)      :: ztra_th
     293    !PhH
     294    REAL,DIMENSION(klon,klev)      :: zrho
     295    REAL,DIMENSION(klon,klev)      :: zdz
     296    REAL                           :: evaplsc,dx,beta ! variable pour lessivage Genthon
     297    REAL,DIMENSION(klon)           :: his_dh          ! ---
     298    ! in-cloud scav variables
     299    REAL           :: ql_incloud_ref     ! ref value of in-cloud condensed water content
     300
     301    !Controles:
     302    !---------
     303    INTEGER,SAVE :: iflag_vdf_trac,iflag_con_trac,iflag_the_trac
     304    INTEGER,SAVE  :: iflag_con_trac_omp, iflag_vdf_trac_omp,iflag_the_trac_omp
    294305!$OMP THREADPRIVATE(iflag_vdf_trac,iflag_con_trac,iflag_the_trac)
    295306
    296   LOGICAL,SAVE :: lessivage
     307    LOGICAL,SAVE :: lessivage
    297308!$OMP THREADPRIVATE(lessivage)
    298309
    299   CHARACTER(len=8),DIMENSION(nbtr) :: solsym
    300 !RomP >>>
    301   INTEGER,SAVE  :: iflag_lscav_omp,iflag_lscav
    302   LOGICAL,SAVE  :: convscav_omp,convscav
     310    CHARACTER(len=8),DIMENSION(nbtr) :: solsym
     311    !RomP >>>
     312    INTEGER,SAVE  :: iflag_lscav_omp,iflag_lscav
     313    LOGICAL,SAVE  :: convscav_omp,convscav
    303314!$OMP THREADPRIVATE(iflag_lscav)
    304315!$OMP THREADPRIVATE(convscav)
    305 !RomP <<<
    306 !######################################################################
    307 !                    -- INITIALIZATION --
    308 !######################################################################
    309 IF (debutphy) THEN
    310 ALLOCATE(d_tr_cl(klon,klev,nbtr),d_tr_dry(klon,nbtr))
    311 ALLOCATE(flux_tr_dry(klon,nbtr),d_tr_dec(klon,klev,nbtr),d_tr_cv(klon,klev,nbtr))
    312 ALLOCATE(d_tr_insc(klon,klev,nbtr),d_tr_bcscav(klon,klev,nbtr))
    313 ALLOCATE(d_tr_evapls(klon,klev,nbtr),d_tr_ls(klon,klev,nbtr))
    314 ALLOCATE(qPrls(klon,nbtr),d_tr_trsp(klon,klev,nbtr))
    315 ALLOCATE(d_tr_sscav(klon,klev,nbtr),d_tr_sat(klon,klev,nbtr))
    316 ALLOCATE(d_tr_uscav(klon,klev,nbtr),qPr(klon,klev,nbtr),qDi(klon,klev,nbtr))
    317 ALLOCATE(qPa(klon,klev,nbtr),qMel(klon,klev,nbtr))
    318 ALLOCATE(qTrdi(klon,klev,nbtr),dtrcvMA(klon,klev,nbtr))
    319 ALLOCATE(d_tr_th(klon,klev,nbtr))
    320 ALLOCATE(d_tr_lessi_impa(klon,klev,nbtr),d_tr_lessi_nucl(klon,klev,nbtr))
    321 ENDIF
    322 
    323   DO k=1,klev
    324      DO i=1,klon
    325       sourceBE(i,k)=0.
    326       Mint(i,k)=0.
    327       zrho(i,k)=0.
    328       zdz(i,k)=0.
    329      END DO
    330   END DO
    331 
    332   DO it=1, nbtr
    333    DO k=1,klev
    334     DO i=1,klon
    335     d_tr_insc(i,k,it)=0.
    336     d_tr_bcscav(i,k,it)=0.
    337     d_tr_evapls(i,k,it)=0.
    338     d_tr_ls(i,k,it)=0.
    339     d_tr_cv(i,k,it)=0.
    340     d_tr_cl(i,k,it)=0.
    341     d_tr_trsp(i,k,it)=0.
    342     d_tr_sscav(i,k,it)=0.
    343     d_tr_sat(i,k,it)=0.
    344     d_tr_uscav(i,k,it)=0.
    345     d_tr_lessi_impa(i,k,it)=0.
    346     d_tr_lessi_nucl(i,k,it)=0.
    347     qDi(i,k,it)=0.
    348     qPr(i,k,it)=0.
    349     qPa(i,k,it)=0.
    350     qMel(i,k,it)=0.
    351     qTrdi(i,k,it)=0.
    352     dtrcvMA(i,k,it)=0.
    353     zmfd1a(i,k,it)=0.
    354     zmfdam(i,k,it)=0.
    355     zmfphi2(i,k,it)=0.
     316    !RomP <<<
     317    !######################################################################
     318    !                    -- INITIALIZATION --
     319    !######################################################################
     320    IF (debutphy) THEN
     321       ALLOCATE(d_tr_cl(klon,klev,nbtr),d_tr_dry(klon,nbtr))
     322       ALLOCATE(flux_tr_dry(klon,nbtr),d_tr_dec(klon,klev,nbtr),d_tr_cv(klon,klev,nbtr))
     323       ALLOCATE(d_tr_insc(klon,klev,nbtr),d_tr_bcscav(klon,klev,nbtr))
     324       ALLOCATE(d_tr_evapls(klon,klev,nbtr),d_tr_ls(klon,klev,nbtr))
     325       ALLOCATE(qPrls(klon,nbtr),d_tr_trsp(klon,klev,nbtr))
     326       ALLOCATE(d_tr_sscav(klon,klev,nbtr),d_tr_sat(klon,klev,nbtr))
     327       ALLOCATE(d_tr_uscav(klon,klev,nbtr),qPr(klon,klev,nbtr),qDi(klon,klev,nbtr))
     328       ALLOCATE(qPa(klon,klev,nbtr),qMel(klon,klev,nbtr))
     329       ALLOCATE(qTrdi(klon,klev,nbtr),dtrcvMA(klon,klev,nbtr))
     330       ALLOCATE(d_tr_th(klon,klev,nbtr))
     331       ALLOCATE(d_tr_lessi_impa(klon,klev,nbtr),d_tr_lessi_nucl(klon,klev,nbtr))
     332    ENDIF
     333
     334    DO k=1,klev
     335       DO i=1,klon
     336          sourceBE(i,k)=0.
     337          Mint(i,k)=0.
     338          zrho(i,k)=0.
     339          zdz(i,k)=0.
     340       END DO
    356341    END DO
    357    END DO
    358   END DO
    359   IF (debutphy) THEN
    360 !!jyg
     342
     343    DO it=1, nbtr
     344       DO k=1,klev
     345          DO i=1,klon
     346             d_tr_insc(i,k,it)=0.
     347             d_tr_bcscav(i,k,it)=0.
     348             d_tr_evapls(i,k,it)=0.
     349             d_tr_ls(i,k,it)=0.
     350             d_tr_cv(i,k,it)=0.
     351             d_tr_cl(i,k,it)=0.
     352             d_tr_trsp(i,k,it)=0.
     353             d_tr_sscav(i,k,it)=0.
     354             d_tr_sat(i,k,it)=0.
     355             d_tr_uscav(i,k,it)=0.
     356             d_tr_lessi_impa(i,k,it)=0.
     357             d_tr_lessi_nucl(i,k,it)=0.
     358             qDi(i,k,it)=0.
     359             qPr(i,k,it)=0.
     360             qPa(i,k,it)=0.
     361             qMel(i,k,it)=0.
     362             qTrdi(i,k,it)=0.
     363             dtrcvMA(i,k,it)=0.
     364             zmfd1a(i,k,it)=0.
     365             zmfdam(i,k,it)=0.
     366             zmfphi2(i,k,it)=0.
     367          END DO
     368       END DO
     369    END DO
     370
     371    DO k = 1, klev
     372       DO i = 1, klon
     373          delp(i,k) = paprs(i,k)-paprs(i,k+1)
     374       END DO
     375    END DO
     376
     377    IF (debutphy) THEN
     378       !!jyg
    361379!$OMP BARRIER
    362    ecrit_tra=86400. ! frequence de stokage en dur
    363                     ! obsolete car remplace par des ecritures dans phys_output_write
    364 !RomP >>>
    365 !
    366 !Config Key  = convscav
    367 !Config Desc = Convective scavenging switch: 0=off, 1=on.
    368 !Config Def  = .false.
    369 !Config Help =
    370 !
     380       ecrit_tra=86400. ! frequence de stokage en dur
     381       ! obsolete car remplace par des ecritures dans phys_output_write
     382       !RomP >>>
     383       !
     384       !Config Key  = convscav
     385       !Config Desc = Convective scavenging switch: 0=off, 1=on.
     386       !Config Def  = .false.
     387       !Config Help =
     388       !
    371389!$OMP MASTER
    372   convscav_omp=.false.
    373   call getin('convscav', convscav_omp)
    374   iflag_vdf_trac_omp=1
    375   call getin('iflag_vdf_trac', iflag_vdf_trac_omp)
    376   iflag_con_trac_omp=1
    377   call getin('iflag_con_trac', iflag_con_trac_omp)
    378   iflag_the_trac_omp=1
    379   call getin('iflag_the_trac', iflag_the_trac_omp)
     390       convscav_omp=.false.
     391       call getin('convscav', convscav_omp)
     392       iflag_vdf_trac_omp=1
     393       call getin('iflag_vdf_trac', iflag_vdf_trac_omp)
     394       iflag_con_trac_omp=1
     395       call getin('iflag_con_trac', iflag_con_trac_omp)
     396       iflag_the_trac_omp=1
     397       call getin('iflag_the_trac', iflag_the_trac_omp)
    380398!$OMP END MASTER
    381399!$OMP BARRIER
    382   convscav=convscav_omp
    383   iflag_vdf_trac=iflag_vdf_trac_omp
    384   iflag_con_trac=iflag_con_trac_omp
    385   iflag_the_trac=iflag_the_trac_omp
    386   print*,'phytrac passage dans routine conv avec lessivage', convscav
    387 !
    388 !Config Key  = iflag_lscav
    389 !Config Desc = Large scale scavenging parametrization: 0=none, 1=old(Genthon92),
    390 !              2=1+PHeinrich, 3=Reddy_Boucher2004, 4=3+RPilon.
    391 !Config Def  = 1
    392 !Config Help =
    393 !
     400       convscav=convscav_omp
     401       iflag_vdf_trac=iflag_vdf_trac_omp
     402       iflag_con_trac=iflag_con_trac_omp
     403       iflag_the_trac=iflag_the_trac_omp
     404       write(lunout,*) 'phytrac passage dans routine conv avec lessivage', convscav
     405       !
     406       !Config Key  = iflag_lscav
     407       !Config Desc = Large scale scavenging parametrization: 0=none, 1=old(Genthon92),
     408       !              2=1+PHeinrich, 3=Reddy_Boucher2004, 4=3+RPilon.
     409       !Config Def  = 1
     410       !Config Help =
     411       !
    394412!$OMP MASTER
    395   iflag_lscav_omp=1
    396   call getin('iflag_lscav', iflag_lscav_omp)
     413       iflag_lscav_omp=1
     414       call getin('iflag_lscav', iflag_lscav_omp)
    397415!$OMP END MASTER
    398416!$OMP BARRIER
    399   iflag_lscav=iflag_lscav_omp
     417       iflag_lscav=iflag_lscav_omp
     418       !
     419       SELECT CASE(iflag_lscav)
     420       CASE(0)
     421          WRITE(lunout,*)  'Large scale scavenging: none'
     422       CASE(1)
     423          WRITE(lunout,*)  'Large scale scavenging: C. Genthon, Tellus(1992), 44B, 371-389'
     424       CASE(2)
     425          WRITE(lunout,*)  'Large scale scavenging: C. Genthon, modified P. Heinrich'
     426       CASE(3)
     427          WRITE(lunout,*)  'Large scale scavenging: M. Shekkar Reddy and O. Boucher, JGR(2004), 109, D14202'
     428       CASE(4)
     429          WRITE(lunout,*)  'Large scale scavenging: Reddy and Boucher, modified R. Pilon'
     430       END SELECT
     431       !RomP <<<
     432       WRITE(*,*) 'FIRST TIME IN PHYTRAC : pdtphys(sec) = ',pdtphys,'ecrit_tra (sec) = ',ecrit_tra
     433       ALLOCATE( source(klon,nbtr), stat=ierr)
     434       IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 1',1)
     435
     436       ALLOCATE( aerosol(nbtr), stat=ierr)
     437       IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 2',1)
     438
     439
     440       ! Initialize module for specific tracers
     441       SELECT CASE(type_trac)
     442       CASE('lmdz')
     443          CALL traclmdz_init(pctsrf,xlat,xlon,ftsol,tr_seri,t_seri,pplay,sh,pdtphys,aerosol,lessivage)
     444       CASE('inca')
     445          source(:,:)=0.
     446          CALL tracinca_init(aerosol,lessivage)
     447       CASE('repr')
     448          source(:,:)=0.
     449       END SELECT
     450
     451       !
     452       !--initialising coefficients for scavenging in the case of NP
     453       !
     454       ALLOCATE(flag_cvltr(nbtr))
     455       IF (iflag_con.EQ.3) THEN
     456          !
     457          ALLOCATE(ccntrAA(nbtr))
     458          ALLOCATE(ccntrENV(nbtr))
     459          ALLOCATE(coefcoli(nbtr))
     460          !
     461          DO it=1, nbtr
     462             SELECT CASE(type_trac)
     463             CASE('lmdz')
     464                IF (convscav.and.aerosol(it)) THEN
     465                   flag_cvltr(it)=.true.
     466                   ccntrAA(it) =1.0         !--a modifier par JYG a lire depuis fichier
     467                   ccntrENV(it)=1.0
     468                   coefcoli(it)=0.001
     469                ELSE
     470                   flag_cvltr(it)=.false.
     471                ENDIF
     472
     473             CASE('inca')
     474!                IF ((it.EQ.id_Rn222) .OR. ((it.GE.id_SO2) .AND. (it.LE.id_NH3)) ) THEN
     475!                   !--gas-phase species
     476!                   flag_cvltr(it)=.false.
    400477!
    401   SELECT CASE(iflag_lscav)
    402   CASE(0)
    403    PRINT*, 'Large scale scavenging: none'
    404   CASE(1)
    405    PRINT*, 'Large scale scavenging: C. Genthon, Tellus(1992), 44B, 371-389'
    406   CASE(2)
    407    PRINT*, 'Large scale scavenging: C. Genthon, modified P. Heinrich'
    408   CASE(3)
    409    PRINT*, 'Large scale scavenging: M. Shekkar Reddy and O. Boucher, JGR(2004), 109, D14202'
    410   CASE(4)
    411    PRINT*, 'Large scale scavenging: Reddy and Boucher, modified R. Pilon'
    412   END SELECT
    413 !RomP <<<
    414      WRITE(*,*) 'FIRST TIME IN PHYTRAC : pdtphys(sec) = ',pdtphys,'ecrit_tra (sec) = ',ecrit_tra
    415      ALLOCATE( source(klon,nbtr), stat=ierr)
    416      IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 1',1)
    417      
    418      ALLOCATE( aerosol(nbtr), stat=ierr)
    419      IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 2',1)
    420      
    421 
    422      ! Initialize module for specific tracers
    423      SELECT CASE(type_trac)
    424      CASE('lmdz')
    425         CALL traclmdz_init(pctsrf, xlat, xlon, ftsol, tr_seri, t_seri, pplay, sh, pdtphys, aerosol, lessivage)
    426      CASE('inca')
    427         source(:,:)=0.
    428         CALL tracinca_init(aerosol,lessivage)
    429      CASE('repr')
    430         source(:,:)=0.
    431      END SELECT
    432 !
    433 ! Initialize diagnostic output
    434 ! ----------------------------
     478!                ELSEIF ( (it.GE.id_CIDUSTM) .AND. (it.LE.id_AIN) ) THEN
     479!                   !--insoluble aerosol species
     480!                   flag_cvltr(it)=.true.
     481!                   ccntrAA(it)=0.7
     482!                   ccntrENV(it)=0.7
     483!                   coefcoli(it)=0.001
     484!                ELSEIF ( (it.EQ.id_Pb210) .OR. ((it.GE.id_CSSSM) .AND. (it.LE.id_SSN))) THEN
     485!                   !--soluble aerosol species
     486!                   flag_cvltr(it)=.true.
     487!                   ccntrAA(it)=0.9
     488!                   ccntrENV(it)=0.9
     489!                   coefcoli(it)=0.001
     490!                ELSE
     491!                   WRITE(lunout,*) 'pb it=', it
     492!                   CALL abort_gcm('phytrac','pb it scavenging',1)
     493!                ENDIF
     494                !--test OB
     495                !--for now we do not scavenge in cvltr
     496                flag_cvltr(it)=.false.
     497             END SELECT
     498          ENDDO
     499          !
     500       ELSE ! iflag_con .ne. 3
     501          flag_cvltr(:) = .false.
     502       ENDIF
     503       !
     504       ! Initialize diagnostic output
     505       ! ----------------------------
    435506#ifdef CPP_IOIPSL
    436 !     INCLUDE "ini_histrac.h"
     507       !     INCLUDE "ini_histrac.h"
    437508#endif
    438   END IF ! of IF (debutphy)
    439 !############################################ END INITIALIZATION #######
    440 
    441   DO k=1,klev
    442      DO i=1,klon
    443         zmasse(i,k)=(paprs(i,k)-paprs(i,k+1))/rg
    444      END DO
    445   END DO
    446 !
    447   IF (id_be .GT. 0) THEN
    448   DO k=1,klev
    449      DO i=1,klon
    450         sourceBE(i,k)=srcbe(i,k)       !RomP  -> pour sortie histrac
    451      END DO
    452   END DO
    453   ENDIF
    454 
    455 !===============================================================================
    456 !    -- Do specific treatment according to chemestry model or local LMDZ tracers
    457 !     
    458 !===============================================================================
    459   SELECT CASE(type_trac)
    460   CASE('lmdz')
    461      !    -- Traitement des traceurs avec traclmdz
    462      CALL traclmdz(nstep, julien, gmtime, pdtphys, t_seri, paprs, pplay, &
    463           cdragh,  coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon,iflag_vdf_trac>=0,sh, &
    464            rh, pphi, ustar, wstar, ale_bl, ale_wake,  u10m, v10m, &
    465            tr_seri, source, solsym, d_tr_cl,d_tr_dec, zmasse)               !RomP
    466 
    467   CASE('inca')
    468      !    -- CHIMIE INCA  config_inca = aero or chem --
    469 
    470      CALL tracinca(&
    471           nstep,    julien,   gmtime,         lafin,     &
    472           pdtphys,  t_seri,   paprs,          pplay,     &
    473           pmfu,     ftsol,    pctsrf,         pphis,     &
    474           pphi,     albsol,   sh,             rh,        &
    475           cldfra,   rneb,     diafra,         cldliq,    &
    476           itop_con, ibas_con, pmflxr,         pmflxs,    &
    477           prfl,     psfl,     aerosol_couple, flxmass_w, &
    478           tau_aero, piz_aero, cg_aero,        ccm,       &
    479           rfname,                                        &
    480           tr_seri,  source,   solsym)     
    481 
    482   CASE('repr')
    483      !   -- CHIMIE REPROBUS --
    484 
    485      CALL tracreprobus(pdtphys, gmtime, debutphy, julien, &
    486           presnivs, xlat, xlon, pphis, pphi, &
    487           t_seri, pplay, paprs, sh , &
    488           tr_seri, solsym)
    489      
    490   END SELECT
    491 !======================================================================
    492 !       -- Calcul de l'effet de la convection --
    493 !======================================================================
    494 
    495   IF (iflag_con_trac==1) THEN
    496      DO it=1, nbtr
    497         IF ( conv_flg(it) == 0 ) CYCLE
    498         IF (iflag_con.LT.2) THEN
    499            d_tr_cv(:,:,it)=0.
    500         ELSE IF (iflag_con.EQ.2) THEN
    501 !..Tiedke
    502            CALL nflxtr(pdtphys, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &
    503                 pplay, paprs, tr_seri(:,:,it), d_tr_cv(:,:,it))
    504 ! RomP >>>               
    505         ELSE   
    506 !..K.Emanuel                  !RomP modif arg
    507         if (convscav.and.aerosol(it)) then    ! lessivage convectif pour aerosol
    508 !
    509           CALL cvltr(pdtphys, da, phi,phi2,d1a,dam, mp,ep,              &
    510 !!               sigd,sij,clw,elij,epmlmMm,eplaMm,                        &   !RL
    511                sigd,sij,wght_cvfd,clw,elij,epmlmMm,eplaMm,              &     !RL
    512                pmflxr,pmflxs,evap,t_seri,wdtrainA,wdtrainM,             &   
    513                paprs,it,tr_seri,upwd,dnwd,itop_con,ibas_con,            &
    514                d_tr_cv,d_tr_trsp,d_tr_sscav,d_tr_sat,d_tr_uscav,qDi,qPr,&
    515                qPa,qMel,qTrdi,dtrcvMA,Mint,                             &
    516                zmfd1a,zmfphi2,zmfdam)
    517         else  !pas de lessivage convectif ou n'est pas un aerosol
    518 !!           CALL cvltrorig(it,pdtphys, da, phi,mp,paprs,pplay,tr_seri,&      !jyg
    519 !!                    upwd,dnwd,d_tr_cv)                                      !jyg
    520            CALL cvltr_noscav(it,pdtphys, da, phi,mp,wght_cvfd,paprs,pplay, &  !jyg
    521                     tr_seri,upwd,dnwd,d_tr_cv)                                !jyg
    522         endif
    523         END IF
    524 ! RomP <<<
    525 
    526         DO k = 1, klev
    527            DO i = 1, klon       
    528               tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cv(i,k,it)
    529            END DO
    530         END DO
    531 
    532         CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'convection it = '//solsym(it))
    533              
    534      END DO ! nbtr
    535   END IF ! convection
    536 
    537 !======================================================================
    538 !    -- Calcul de l'effet des thermiques --
    539 !======================================================================
    540 
    541   DO it=1,nbtr
    542      DO k=1,klev
    543         DO i=1,klon
    544            d_tr_th(i,k,it)=0.
    545            tr_seri(i,k,it)=MAX(tr_seri(i,k,it),0.)
    546            tr_seri(i,k,it)=MIN(tr_seri(i,k,it),1.e10)
    547         END DO
    548      END DO
    549   END DO
    550  
    551   IF (iflag_thermals.GT.0.AND.iflag_the_trac>0) THEN   
    552      nsplit=10
    553      DO it=1, nbtr
    554         DO isplit=1,nsplit
    555 
    556            CALL dqthermcell(klon,klev,pdtphys/nsplit, &
    557                 fm_therm,entr_therm,zmasse, &
    558                 tr_seri(1:klon,1:klev,it),d_tr,ztra_th)
    559 
    560            DO k=1,klev
    561               DO i=1,klon
    562                  d_tr(i,k)=pdtphys*d_tr(i,k)/nsplit
    563                  d_tr_th(i,k,it)=d_tr_th(i,k,it)+d_tr(i,k)
    564                  tr_seri(i,k,it)=MAX(tr_seri(i,k,it)+d_tr(i,k),0.)
    565               END DO
    566            END DO
    567         END DO ! nsplit
    568      END DO ! it
    569   END IF ! Thermiques
    570 
    571 !======================================================================
    572 !     -- Calcul de l'effet de la couche limite --
    573 !======================================================================
    574 
    575   DO k = 1, klev
    576      DO i = 1, klon
    577         delp(i,k) = paprs(i,k)-paprs(i,k+1)
    578      END DO
    579   END DO
    580 
    581   IF (iflag_vdf_trac==1) THEN
    582      DO it=1, nbtr
    583         if (prt_level > 20) write(lunout,*)'trac pbl ',it,pbl_flg(it)
    584         IF( pbl_flg(it) /= 0 ) THEN
    585            CALL cltrac(pdtphys, coefh,t_seri,       &
    586                 tr_seri(:,:,it), source(:,it),      &
    587                 paprs, pplay, delp,                 &
    588                 d_tr_cl(:,:,it),d_tr_dry(:,it),flux_tr_dry(:,it))
    589            tr_seri(:,:,it)=tr_seri(:,:,it)+d_tr_cl(:,:,it)
    590         END IF
    591      END DO
    592   ELSE IF (iflag_vdf_trac==0) THEN
    593 !   Injection of source in the first model layer
     509       !
     510       ! print out all tracer flags
     511       !
     512       WRITE(lunout,*) 'print out all tracer flags'
     513       WRITE(lunout,*) 'type_trac      =', type_trac
     514       WRITE(lunout,*) 'config_inca    =', config_inca
     515       WRITE(lunout,*) 'iflag_con_trac =', iflag_con_trac
     516       WRITE(lunout,*) 'iflag_con      =', iflag_con
     517       WRITE(lunout,*) 'convscav       =', convscav
     518       WRITE(lunout,*) 'iflag_lscav    =', iflag_lscav
     519       WRITE(lunout,*) 'aerosol        =', aerosol
     520       WRITE(lunout,*) 'iflag_the_trac =', iflag_the_trac
     521       WRITE(lunout,*) 'iflag_thermals =', iflag_thermals
     522       WRITE(lunout,*) 'iflag_vdf_trac =', iflag_vdf_trac
     523       WRITE(lunout,*) 'pbl_flg        =', pbl_flg
     524       WRITE(lunout,*) 'lessivage      =', lessivage
     525       write(lunout,*)  'flag_cvltr    = ', flag_cvltr
     526
     527       IF (lessivage.AND.config_inca.EQ.'inca') THEN
     528          CALL abort_gcm('phytrac', 'lessivage=T config_inca=inca impossible',1)
     529          STOP
     530       ENDIF
     531       !
     532    END IF ! of IF (debutphy)
     533    !############################################ END INITIALIZATION #######
     534
     535    DO k=1,klev
     536       DO i=1,klon
     537          zmasse(i,k)=(paprs(i,k)-paprs(i,k+1))/rg
     538       END DO
     539    END DO
     540    !
     541    IF (id_be .GT. 0) THEN
     542       DO k=1,klev
     543          DO i=1,klon
     544             sourceBE(i,k)=srcbe(i,k)       !RomP  -> pour sortie histrac
     545          END DO
     546       END DO
     547    ENDIF
     548
     549    !===============================================================================
     550    !    -- Do specific treatment according to chemestry model or local LMDZ tracers
     551    !     
     552    !===============================================================================
     553    SELECT CASE(type_trac)
     554    CASE('lmdz')
     555       !    -- Traitement des traceurs avec traclmdz
     556       CALL traclmdz(nstep, julien, gmtime, pdtphys, t_seri, paprs, pplay, &
     557            cdragh,  coefh, yu1, yv1, ftsol, pctsrf, xlat, xlon,iflag_vdf_trac>=0,sh, &
     558            rh, pphi, ustar, wstar, ale_bl, ale_wake,  u10m, v10m, &
     559            tr_seri, source, solsym, d_tr_cl,d_tr_dec, zmasse)               !RomP
     560
     561    CASE('inca')
     562       !    -- CHIMIE INCA  config_inca = aero or chem --
     563
     564       CALL tracinca(&
     565            nstep,    julien,   gmtime,         lafin,     &
     566            pdtphys,  t_seri,   paprs,          pplay,     &
     567            pmfu,     upwd,     ftsol,  pctsrf, pphis,     &
     568            pphi,     albsol,   sh,             rh,        &
     569            cldfra,   rneb,     diafra,         cldliq,    &
     570            itop_con, ibas_con, pmflxr,         pmflxs,    &
     571            prfl,     psfl,     aerosol_couple, flxmass_w, &
     572            tau_aero, piz_aero, cg_aero,        ccm,       &
     573            rfname,                                        &
     574            tr_seri,  source,   solsym)     
     575
     576    CASE('repr')
     577       !   -- CHIMIE REPROBUS --
     578
     579       CALL tracreprobus(pdtphys, gmtime, debutphy, julien, &
     580            presnivs, xlat, xlon, pphis, pphi, &
     581            t_seri, pplay, paprs, sh , &
     582            tr_seri, solsym)
     583
     584    END SELECT
     585    !======================================================================
     586    !       -- Calcul de l'effet de la convection --
     587    !======================================================================
     588
     589    IF (iflag_con_trac==1) THEN
     590
     591       DO it=1, nbtr
     592          IF ( conv_flg(it) == 0 ) CYCLE
     593          IF (iflag_con.LT.2) THEN
     594             !--pas de transport convectif
     595
     596             d_tr_cv(:,:,it)=0.
     597          ELSE IF (iflag_con.EQ.2) THEN
     598             !--ancien transport convectif de Tiedtke
     599
     600             CALL nflxtr(pdtphys, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &
     601                  pplay, paprs, tr_seri(:,:,it), d_tr_cv(:,:,it))
     602          ELSE   
     603             !--nouveau transport convectif de Emanuel
     604
     605             IF (flag_cvltr(it)) THEN
     606                !--nouveau transport convectif de Emanuel avec lessivage convectif
     607                !
     608                !
     609                ccntrAA_3d(:,:) =ccntrAA(it)
     610                ccntrENV_3d(:,:)=ccntrENV(it)
     611                coefcoli_3d(:,:)=coefcoli(it)
     612
     613                !--beware this interface is a bit weird because it is called for each tracer
     614                !--with the full array tr_seri even if only item it is processed
     615
     616                CALL cvltr_scav(pdtphys, da, phi,phi2,d1a,dam, mp,ep,         &
     617                     sigd,sij,wght_cvfd,clw,elij,epmlmMm,eplaMm,              &     
     618                     pmflxr,pmflxs,evap,t_seri,wdtrainA,wdtrainM,             &   
     619                     paprs,it,tr_seri,upwd,dnwd,itop_con,ibas_con,            &
     620                     ccntrAA_3d,ccntrENV_3d,coefcoli_3d,                      &
     621                     d_tr_cv,d_tr_trsp,d_tr_sscav,d_tr_sat,d_tr_uscav,qDi,qPr,&
     622                     qPa,qMel,qTrdi,dtrcvMA,Mint,                             &
     623                     zmfd1a,zmfphi2,zmfdam)
     624
     625
     626             ELSE  !---flag_cvltr(it).EQ.FALSE
     627                !--nouveau transport convectif de Emanuel mais pas de lessivage convectif
     628
     629                !--beware this interface is a bit weird because it is called for each tracer
     630                !--with the full array tr_seri even if only item it is processed
     631                !
     632                CALL cvltr_noscav(it,pdtphys, da, phi,mp,wght_cvfd,paprs,pplay, &  !jyg
     633                     tr_seri,upwd,dnwd,d_tr_cv)                                !jyg
     634
     635             ENDIF
     636
     637          ENDIF !--iflag
     638
     639          !--on ajoute les tendances
     640
     641          DO k = 1, klev
     642             DO i = 1, klon       
     643                tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cv(i,k,it)
     644             END DO
     645          END DO
     646
     647          CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'convection it = '//solsym(it))
     648
     649       END DO ! nbtr
     650
     651    END IF ! convection
     652
     653    !======================================================================
     654    !    -- Calcul de l'effet des thermiques --
     655    !======================================================================
     656
    594657    DO it=1,nbtr
    595        d_tr_cl(:,1,it)=source(:,it)*rg/delp(:,1)*pdtphys
    596        tr_seri(:,1,it)=tr_seri(:,1,it)+d_tr_cl(:,1,it)
    597     ENDDO
    598     d_tr_cl(:,2:klev,1:nbtr)=0.
    599   ELSE IF (iflag_vdf_trac==-1) THEN
    600     d_tr_cl=0.
    601   ELSE
    602     CALL abort_gcm('iflag_vdf_trac', 'cas non prevu',1)
    603   END IF ! couche limite
    604 
    605 
    606 
    607 !======================================================================
    608 !   Calcul de l'effet de la precipitation grande echelle
    609 !======================================================================
    610   IF (lessivage) THEN
    611 
    612    ql_incloud_ref = 10.e-4
    613    ql_incloud_ref =  5.e-4
    614 
    615 
    616 ! calcul du contenu en eau liquide au sein du nuage
    617      ql_incl = ql_incloud_ref
    618 ! choix du lessivage
    619 !
    620   IF (iflag_lscav .EQ. 3 .OR. iflag_lscav .EQ. 4) THEN
    621 ! ********  Olivier Boucher version (3) possibly with modified ql_incl (4)
    622 !
    623    DO it = 1, nbtr
    624 !  incloud scavenging and removal by large scale rain ! orig : ql_incl was replaced by 0.5e-3 kg/kg
    625 ! the value 0.5e-3 kg/kg is from Giorgi and Chameides (1986), JGR
    626 ! Liu (2001) proposed to use 1.5e-3 kg/kg
    627 
    628     CALL lsc_scav(pdtphys,it,iflag_lscav,ql_incl,prfl,psfl,rneb,beta_fisrt,  &
     658       DO k=1,klev
     659          DO i=1,klon
     660             d_tr_th(i,k,it)=0.
     661             tr_seri(i,k,it)=MAX(tr_seri(i,k,it),0.)
     662             tr_seri(i,k,it)=MIN(tr_seri(i,k,it),1.e10)
     663          END DO
     664       END DO
     665    END DO
     666
     667    IF (iflag_thermals.GT.0.AND.iflag_the_trac>0) THEN   
     668
     669       DO it=1, nbtr
     670
     671          CALL thermcell_dq(klon,klev,1,pdtphys,fm_therm,entr_therm, &
     672               zmasse,tr_seri(1:klon,1:klev,it),        &
     673               d_tr_th(1:klon,1:klev,it),ztra_th,0 )
     674
     675          DO k=1,klev
     676             DO i=1,klon
     677                d_tr_th(i,k,it)=pdtphys*d_tr_th(i,k,it)
     678                tr_seri(i,k,it)=MAX(tr_seri(i,k,it)+d_tr_th(i,k,it),0.)
     679             END DO
     680          END DO
     681
     682       END DO ! it
     683
     684    END IF ! Thermiques
     685
     686    !======================================================================
     687    !     -- Calcul de l'effet de la couche limite --
     688    !======================================================================
     689
     690    IF (iflag_vdf_trac==1) THEN
     691
     692       !  Injection during BL mixing
     693       !
     694       DO it=1, nbtr
     695          !
     696          IF( pbl_flg(it) /= 0 ) THEN
     697             !
     698             CALL cltrac(pdtphys, coefh,t_seri,       &
     699                  tr_seri(:,:,it), source(:,it),      &
     700                  paprs, pplay, delp,                 &
     701                  d_tr_cl(:,:,it),d_tr_dry(:,it),flux_tr_dry(:,it))
     702             !
     703             tr_seri(:,:,it)=tr_seri(:,:,it)+d_tr_cl(:,:,it)
     704             !
     705          END IF
     706          !
     707       END DO
     708       !
     709    ELSE IF (iflag_vdf_trac==0) THEN
     710       !
     711       !   Injection of source in the first model layer
     712       !
     713       DO it=1,nbtr
     714          d_tr_cl(:,1,it)=source(:,it)*RG/delp(:,1)*pdtphys
     715          tr_seri(:,1,it)=tr_seri(:,1,it)+d_tr_cl(:,1,it)
     716       ENDDO
     717       d_tr_cl(:,2:klev,1:nbtr)=0.
     718       !
     719    ELSE IF (iflag_vdf_trac==-1) THEN
     720       !
     721       ! Nothing happens
     722       !
     723       d_tr_cl=0.
     724       !
     725    ELSE
     726       !
     727       CALL abort_gcm('iflag_vdf_trac', 'cas non prevu',1)
     728       !
     729    END IF ! couche limite
     730
     731    !======================================================================
     732    !   Calcul de l'effet de la precipitation grande echelle
     733    !   POUR INCA le lessivage est fait directement dans INCA
     734    !======================================================================
     735
     736    IF (lessivage) THEN
     737
     738       ql_incloud_ref = 10.e-4
     739       ql_incloud_ref =  5.e-4
     740
     741
     742       ! calcul du contenu en eau liquide au sein du nuage
     743       ql_incl = ql_incloud_ref
     744       ! choix du lessivage
     745       !
     746       IF (iflag_lscav .EQ. 3 .OR. iflag_lscav .EQ. 4) THEN
     747          ! ********  Olivier Boucher version (3) possibly with modified ql_incl (4)
     748          !
     749          DO it = 1, nbtr
     750             !  incloud scavenging and removal by large scale rain ! orig : ql_incl was replaced by 0.5e-3 kg/kg
     751             ! the value 0.5e-3 kg/kg is from Giorgi and Chameides (1986), JGR
     752             ! Liu (2001) proposed to use 1.5e-3 kg/kg
     753
     754             CALL lsc_scav(pdtphys,it,iflag_lscav,ql_incl,prfl,psfl,rneb,beta_fisrt,  &
    629755                  beta_v1,pplay,paprs,t_seri,tr_seri,d_tr_insc,   &
    630756                  d_tr_bcscav,d_tr_evapls,qPrls)
    631757
    632 !large scale scavenging tendency
    633    DO k = 1, klev
    634     DO i = 1, klon
    635     d_tr_ls(i,k,it)=d_tr_insc(i,k,it)+d_tr_bcscav(i,k,it)+d_tr_evapls(i,k,it)
    636     tr_seri(i,k,it)=tr_seri(i,k,it)+d_tr_ls(i,k,it)
    637     ENDDO
    638    ENDDO
    639      CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'lsc scav it = '//solsym(it))
    640    END DO  !tr
    641 
    642  ELSE IF (iflag_lscav .EQ. 2) THEN ! frac_impa, frac_nucl
    643 ! *********   modified  old version
    644 
    645      d_tr_lessi_nucl(:,:,:) = 0.
    646      d_tr_lessi_impa(:,:,:) = 0.
    647      flestottr(:,:,:) = 0.
    648 ! Tendance des aerosols nuclees et impactes
    649      DO it = 1, nbtr
    650         IF (aerosol(it)) THEN
    651         his_dh(:)=0.
    652            DO k = 1, klev
    653               DO i = 1, klon
    654 !PhH
    655               zrho(i,k)=pplay(i,k)/t_seri(i,k)/RD
    656               zdz(i,k)=(paprs(i,k)-paprs(i,k+1))/zrho(i,k)/RG
    657 !
    658               END DO
    659            END DO
    660 
    661           DO k=klev-1, 1, -1
    662             DO i=1, klon
    663 !             d_tr_ls(i,k,it)=tr_seri(i,k,it)*(frac_impa(i,k)*frac_nucl(i,k)-1.)
    664              dx=d_tr_ls(i,k,it)
    665              his_dh(i)=his_dh(i)-dx*zrho(i,k)*zdz(i,k)/pdtphys !  kg/m2/s
    666              evaplsc = prfl(i,k) - prfl(i,k+1) + psfl(i,k) - psfl(i,k+1)
    667 ! Evaporation Partielle -> Liberation Partielle 0.5*evap
    668             IF ( evaplsc .LT.0..and.abs(prfl(i,k+1)+psfl(i,k+1)).gt.1.e-10) THEN
    669                 evaplsc = (-evaplsc)/(prfl(i,k+1)+psfl(i,k+1))
    670 ! evaplsc est donc positif, his_dh(i) est positif
    671 !-------------- 
    672              d_tr_evapls(i,k,it)=0.5*evaplsc*(d_tr_lessi_nucl(i,k+1,it) &
    673                                   +d_tr_lessi_impa(i,k+1,it))
    674 !-------------   d_tr_evapls(i,k,it)=-0.5*evaplsc*(d_tr_lsc(i,k+1,it))
    675              beta=0.5*evaplsc
    676               if ((prfl(i,k)+psfl(i,k)).lt.1.e-10) THEN
    677                beta=1.0*evaplsc
    678               endif
    679             dx=beta*his_dh(i)/zrho(i,k)/zdz(i,k)*pdtphys
    680             his_dh(i)=(1.-beta)*his_dh(i)   ! tracer from
    681             d_tr_evapls(i,k,it)=dx
    682             ENDIF
    683             d_tr_ls(i,k,it)=tr_seri(i,k,it)*(frac_impa(i,k)*frac_nucl(i,k)-1.) &
    684                             +d_tr_evapls(i,k,it)
    685 
    686 !--------------
    687                  d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) +    &
    688                       ( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
    689                  d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) +    &
    690                       ( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
    691 !
    692 ! Flux lessivage total
    693                  flestottr(i,k,it) = flestottr(i,k,it) -   &
    694                       ( d_tr_lessi_nucl(i,k,it)   +        &
    695                       d_tr_lessi_impa(i,k,it) ) *          &
    696                       ( paprs(i,k)-paprs(i,k+1) ) /        &
    697                       (RG * pdtphys)
    698 !! Mise a jour des traceurs due a l'impaction,nucleation
    699 !                 tr_seri(i,k,it)=tr_seri(i,k,it)*frac_impa(i,k)*frac_nucl(i,k)
    700 !!  calcul de la tendance liee au lessivage stratiforme
    701 !                 d_tr_ls(i,k,it)=tr_seri(i,k,it)*&
    702 !                                (1.-1./(frac_impa(i,k)*frac_nucl(i,k)))
    703 !--------------
    704               END DO
    705            END DO
    706         END IF
    707      END DO
    708 ! *********   end modified old version
    709 
    710  ELSE IF (iflag_lscav .EQ. 1) THEN ! frac_impa, frac_nucl
    711 ! *********    old version
    712 
    713      d_tr_lessi_nucl(:,:,:) = 0.
    714      d_tr_lessi_impa(:,:,:) = 0.
    715      flestottr(:,:,:) = 0.
    716 !=========================
    717 ! LESSIVAGE LARGE SCALE :
    718 !=========================
    719 
    720 ! Tendance des aerosols nuclees et impactes
    721 ! -----------------------------------------
    722      DO it = 1, nbtr
    723         IF (aerosol(it)) THEN
    724            DO k = 1, klev
    725               DO i = 1, klon
    726                  d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) +    &
    727                       ( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
    728                  d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) +    &
    729                       ( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
    730 
    731 !
    732 ! Flux lessivage total
    733 ! ------------------------------------------------------------
    734                  flestottr(i,k,it) = flestottr(i,k,it) -   &
    735                       ( d_tr_lessi_nucl(i,k,it)   +        &
    736                       d_tr_lessi_impa(i,k,it) ) *          &
    737                       ( paprs(i,k)-paprs(i,k+1) ) /        &
    738                       (RG * pdtphys)
    739 !
    740 ! Mise a jour des traceurs due a l'impaction,nucleation
    741 ! ----------------------------------------------------------------------
    742                  tr_seri(i,k,it)=tr_seri(i,k,it)*frac_impa(i,k)*frac_nucl(i,k)
    743               END DO
    744            END DO
    745         END IF
    746      END DO
    747      
    748 ! *********   end old version
    749   ENDIF  !  iflag_lscav . EQ. 1, 2, 3 or 4
    750 !
    751   END IF !  lessivage
    752 
    753 !=============================================================
    754 !   Ecriture des sorties
    755 !=============================================================
     758             !large scale scavenging tendency
     759             DO k = 1, klev
     760                DO i = 1, klon
     761                   d_tr_ls(i,k,it)=d_tr_insc(i,k,it)+d_tr_bcscav(i,k,it)+d_tr_evapls(i,k,it)
     762                   tr_seri(i,k,it)=tr_seri(i,k,it)+d_tr_ls(i,k,it)
     763                ENDDO
     764             ENDDO
     765             CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'lsc scav it = '//solsym(it))
     766          END DO  !tr
     767
     768       ELSE IF (iflag_lscav .EQ. 2) THEN ! frac_impa, frac_nucl
     769          ! *********   modified  old version
     770
     771          d_tr_lessi_nucl(:,:,:) = 0.
     772          d_tr_lessi_impa(:,:,:) = 0.
     773          flestottr(:,:,:) = 0.
     774          ! Tendance des aerosols nuclees et impactes
     775          DO it = 1, nbtr
     776             IF (aerosol(it)) THEN
     777                his_dh(:)=0.
     778                DO k = 1, klev
     779                   DO i = 1, klon
     780                      !PhH
     781                      zrho(i,k)=pplay(i,k)/t_seri(i,k)/RD
     782                      zdz(i,k)=(paprs(i,k)-paprs(i,k+1))/zrho(i,k)/RG
     783                      !
     784                   END DO
     785                END DO
     786
     787                DO k=klev-1, 1, -1
     788                   DO i=1, klon
     789                      !             d_tr_ls(i,k,it)=tr_seri(i,k,it)*(frac_impa(i,k)*frac_nucl(i,k)-1.)
     790                      dx=d_tr_ls(i,k,it)
     791                      his_dh(i)=his_dh(i)-dx*zrho(i,k)*zdz(i,k)/pdtphys !  kg/m2/s
     792                      evaplsc = prfl(i,k) - prfl(i,k+1) + psfl(i,k) - psfl(i,k+1)
     793                      ! Evaporation Partielle -> Liberation Partielle 0.5*evap
     794                      IF ( evaplsc .LT.0..and.abs(prfl(i,k+1)+psfl(i,k+1)).gt.1.e-10) THEN
     795                         evaplsc = (-evaplsc)/(prfl(i,k+1)+psfl(i,k+1))
     796                         ! evaplsc est donc positif, his_dh(i) est positif
     797                         !-------------- 
     798                         d_tr_evapls(i,k,it)=0.5*evaplsc*(d_tr_lessi_nucl(i,k+1,it) &
     799                              +d_tr_lessi_impa(i,k+1,it))
     800                         !-------------   d_tr_evapls(i,k,it)=-0.5*evaplsc*(d_tr_lsc(i,k+1,it))
     801                         beta=0.5*evaplsc
     802                         if ((prfl(i,k)+psfl(i,k)).lt.1.e-10) THEN
     803                            beta=1.0*evaplsc
     804                         endif
     805                         dx=beta*his_dh(i)/zrho(i,k)/zdz(i,k)*pdtphys
     806                         his_dh(i)=(1.-beta)*his_dh(i)   ! tracer from
     807                         d_tr_evapls(i,k,it)=dx
     808                      ENDIF
     809                      d_tr_ls(i,k,it)=tr_seri(i,k,it)*(frac_impa(i,k)*frac_nucl(i,k)-1.) &
     810                           +d_tr_evapls(i,k,it)
     811
     812                      !--------------
     813                      d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) +    &
     814                           ( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
     815                      d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) +    &
     816                           ( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
     817                      !
     818                      ! Flux lessivage total
     819                      flestottr(i,k,it) = flestottr(i,k,it) -   &
     820                           ( d_tr_lessi_nucl(i,k,it)   +        &
     821                           d_tr_lessi_impa(i,k,it) ) *          &
     822                           ( paprs(i,k)-paprs(i,k+1) ) /        &
     823                           (RG * pdtphys)
     824                      !! Mise a jour des traceurs due a l'impaction,nucleation
     825                      !                 tr_seri(i,k,it)=tr_seri(i,k,it)*frac_impa(i,k)*frac_nucl(i,k)
     826                      !!  calcul de la tendance liee au lessivage stratiforme
     827                      !                 d_tr_ls(i,k,it)=tr_seri(i,k,it)*&
     828                      !                                (1.-1./(frac_impa(i,k)*frac_nucl(i,k)))
     829                      !--------------
     830                   END DO
     831                END DO
     832             END IF
     833          END DO
     834          ! *********   end modified old version
     835
     836       ELSE IF (iflag_lscav .EQ. 1) THEN ! frac_impa, frac_nucl
     837          ! *********    old version
     838
     839          d_tr_lessi_nucl(:,:,:) = 0.
     840          d_tr_lessi_impa(:,:,:) = 0.
     841          flestottr(:,:,:) = 0.
     842          !=========================
     843          ! LESSIVAGE LARGE SCALE :
     844          !=========================
     845
     846          ! Tendance des aerosols nuclees et impactes
     847          ! -----------------------------------------
     848          DO it = 1, nbtr
     849             IF (aerosol(it)) THEN
     850                DO k = 1, klev
     851                   DO i = 1, klon
     852                      d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) +    &
     853                           ( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
     854                      d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) +    &
     855                           ( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
     856
     857                      !
     858                      ! Flux lessivage total
     859                      ! ------------------------------------------------------------
     860                      flestottr(i,k,it) = flestottr(i,k,it) -   &
     861                           ( d_tr_lessi_nucl(i,k,it)   +        &
     862                           d_tr_lessi_impa(i,k,it) ) *          &
     863                           ( paprs(i,k)-paprs(i,k+1) ) /        &
     864                           (RG * pdtphys)
     865                      !
     866                      ! Mise a jour des traceurs due a l'impaction,nucleation
     867                      ! ----------------------------------------------------------------------
     868                      tr_seri(i,k,it)=tr_seri(i,k,it)*frac_impa(i,k)*frac_nucl(i,k)
     869                   END DO
     870                END DO
     871             END IF
     872          END DO
     873
     874          ! *********   end old version
     875       ENDIF  !  iflag_lscav . EQ. 1, 2, 3 or 4
     876       !
     877    END IF !  lessivage
     878
     879    !=============================================================
     880    !   Ecriture des sorties
     881    !=============================================================
    756882#ifdef CPP_IOIPSL
    757 ! INCLUDE "write_histrac.h"
     883    ! INCLUDE "write_histrac.h"
    758884#endif
    759885
    760 END SUBROUTINE phytrac
     886  END SUBROUTINE phytrac
    761887
    762888END MODULE
Note: See TracChangeset for help on using the changeset viewer.