source: LMDZ4/branches/LMDZ4-dev/libf/phylmd/phytrac.F90 @ 1214

Last change on this file since 1214 was 1212, checked in by jghattas, 15 years ago
  • Correction dans la recherche du traceur Berelium
  • Correction des messages d'erreurs
  • Enleve allocate pas necessaire
  • Ajoute argument d'entree pour traclmdz bientot necessaire (pour parametrization Cariolle)
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 15.0 KB
RevLine 
[1191]1!$Id $
2
3SUBROUTINE phytrac(                            &
4     nstep,     julien,   gmtime,   debutphy,  &
5     lafin,     pdtphys,  u, v,     t_seri,     &
6     paprs,     pplay,    pmfu,     pmfd,      &
7     pen_u,     pde_u,    pen_d,    pde_d,     &
8     cdragh,    coefh,    fm_therm, entr_therm,&
9     yu1,       yv1,      ftsol,    pctsrf,    &
10     xlat,      frac_impa,frac_nucl,xlon,      &
11     presnivs,  pphis,    pphi,     albsol,    &
12     sh,        rh,       cldfra,   rneb,      &
13     diafra,    cldliq,   itop_con, ibas_con,  &
14     pmflxr,    pmflxs,   prfl,     psfl,      &
15     da,        phi,      mp,       upwd,      &
16     dnwd,      aerosol_couple,     flxmass_w, &
17     tau_aero,  piz_aero,  cg_aero, ccm,       &
18     rfname,                                   &
19     tr_seri)         
[524]20!
[1191]21!======================================================================
22! Auteur(s) FH
23! Objet: Moniteur general des tendances traceurs
24!======================================================================
[524]25
[1191]26  USE ioipsl
27  USE dimphy
28  USE infotrac
29  USE mod_grid_phy_lmdz
30  USE mod_phys_lmdz_para
31  USE comgeomphy
32  USE iophy
33  USE traclmdz_mod
34  USE tracinca_mod
[524]35
36
[1191]37  IMPLICIT NONE
[959]38
[1191]39  INCLUDE "YOMCST.h"
40  INCLUDE "dimensions.h"
41  INCLUDE "indicesol.h"
42  INCLUDE "clesphys.h"
43  INCLUDE "temps.h"
44  INCLUDE "paramet.h"
45  INCLUDE "control.h"
46  INCLUDE "thermcell.h"
47!==========================================================================
48!                   -- ARGUMENT DESCRIPTION --
49!==========================================================================
[541]50
[1191]51! Input arguments
52!----------------
53!Configuration grille,temps:
54  INTEGER,INTENT(IN) :: nstep      ! Appel physique
55  INTEGER,INTENT(IN) :: julien     ! Jour julien
56  REAL,INTENT(IN)    :: gmtime
57  REAL,INTENT(IN)    :: pdtphys    ! Pas d'integration pour la physique (seconde)
58  LOGICAL,INTENT(IN) :: debutphy   ! le flag de l'initialisation de la physique
59  LOGICAL,INTENT(IN) :: lafin      ! le flag de la fin de la physique
60 
61  REAL,DIMENSION(klon),INTENT(IN) :: xlat    ! latitudes pour chaque point
62  REAL,DIMENSION(klon),INTENT(IN) :: xlon    ! longitudes pour chaque point
63!
64!Physique:
65!--------
66  REAL,DIMENSION(klon,klev),INTENT(IN)   :: t_seri  ! Temperature
67  REAL,DIMENSION(klon,klev),INTENT(IN)   :: u       !
68  REAL,DIMENSION(klon,klev),INTENT(IN)   :: v       !
69  REAL,DIMENSION(klon,klev),INTENT(IN)   :: sh      ! humidite specifique
70  REAL,DIMENSION(klon,klev),INTENT(IN)   :: rh      ! humidite relative
71  REAL,DIMENSION(klon,klev+1),INTENT(IN) :: paprs   ! pression pour chaque inter-couche (en Pa)
72  REAL,DIMENSION(klon,klev),INTENT(IN)   :: pplay   ! pression pour le mileu de chaque couche (en Pa)
73  REAL,DIMENSION(klon,klev),INTENT(IN)   :: pphi    ! geopotentiel
74  REAL,DIMENSION(klon),INTENT(IN)        :: pphis
75  REAL,DIMENSION(klev),INTENT(IN)        :: presnivs
76  REAL,DIMENSION(klon,klev),INTENT(IN)   :: cldliq  ! eau liquide nuageuse
77  REAL,DIMENSION(klon,klev),INTENT(IN)   :: cldfra  ! fraction nuageuse (tous les nuages)
78  REAL,DIMENSION(klon,klev),INTENT(IN)   :: diafra  ! fraction nuageuse (convection ou stratus artificiels)
79  REAL,DIMENSION(klon,klev),INTENT(IN)   :: rneb    ! fraction nuageuse (grande echelle)
80  INTEGER,DIMENSION(klon),INTENT(IN)     :: itop_con
81  INTEGER,DIMENSION(klon),INTENT(IN)     :: ibas_con
82  REAL,DIMENSION(klon),INTENT(IN)        :: albsol  ! albedo surface
83!
84!Convection:
85!----------
86  REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfu  ! flux de masse dans le panache montant
87  REAL,DIMENSION(klon,klev),INTENT(IN) :: pmfd  ! flux de masse dans le panache descendant
88  REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_u ! flux entraine dans le panache montant
89  REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_u ! flux detraine dans le panache montant
90  REAL,DIMENSION(klon,klev),INTENT(IN) :: pen_d ! flux entraine dans le panache descendant
91  REAL,DIMENSION(klon,klev),INTENT(IN) :: pde_d ! flux detraine dans le panache descendant
[619]92
[1191]93!...Tiedke     
94  REAL,DIMENSION(klon,klev+1),INTENT(IN)   :: pmflxr, pmflxs ! Flux precipitant de pluie, neige aux interfaces [convection]
95  REAL,DIMENSION(klon,klev+1),INTENT(IN)   :: prfl, psfl ! Flux precipitant de pluie, neige aux interfaces [large-scale]
[524]96
[1191]97  LOGICAL,INTENT(IN)                       :: aerosol_couple
98  REAL,DIMENSION(klon,klev),INTENT(IN)     :: flxmass_w
99  REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: tau_aero
100  REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: piz_aero
101  REAL,DIMENSION(klon,klev,9,2),INTENT(IN) :: cg_aero
102  CHARACTER(len=4),DIMENSION(9),INTENT(IN) :: rfname
103  REAL,DIMENSION(klon,klev,2),INTENT(IN)   :: ccm
104!... K.Emanuel
105  REAL,DIMENSION(klon,klev),INTENT(IN)     :: da
106  REAL,DIMENSION(klon,klev,klev),INTENT(IN):: phi
107  REAL,DIMENSION(klon,klev),INTENT(IN)     :: mp
108  REAL,DIMENSION(klon,klev),INTENT(IN)     :: upwd      ! saturated updraft mass flux
109  REAL,DIMENSION(klon,klev),INTENT(IN)     :: dnwd      ! saturated downdraft mass flux
110!
111!Thermiques:
112!----------
113  REAL,DIMENSION(klon,klev+1),INTENT(IN)   :: fm_therm
114  REAL,DIMENSION(klon,klev),INTENT(IN)     :: entr_therm
115!
116!Couche limite:
117!--------------
118!
119  REAL,DIMENSION(klon,klev),INTENT(IN) :: cdragh ! coeff drag pour T et Q
120  REAL,DIMENSION(klon,klev),INTENT(IN) :: coefh  ! coeff melange CL (m**2/s)
121  REAL,DIMENSION(klon),INTENT(IN)      :: yu1    ! vents au premier niveau
122  REAL,DIMENSION(klon),INTENT(IN)      :: yv1    ! vents au premier niveau
123!
124!Lessivage:
125!----------
126!
127! pour le ON-LINE
128!
129  REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_impa ! fraction d'aerosols non impactes
130  REAL,DIMENSION(klon,klev),INTENT(IN) :: frac_nucl ! fraction d'aerosols non nuclees
[524]131
[1191]132! Arguments necessaires pour les sources et puits de traceur:
133  REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: ftsol  ! Temperature du sol (surf)(Kelvin)
134  REAL,DIMENSION(klon,nbsrf),INTENT(IN) :: pctsrf ! Pourcentage de sol (nature du sol)
[524]135
136
[1191]137! Output argument
138!----------------
139  REAL,DIMENSION(klon,klev,nbtr),INTENT(INOUT) :: tr_seri ! Concentration Traceur [U/KgA] 
[524]140
[1191]141!=======================================================================================
142!                        -- LOCAL VARIABLES --
143!=======================================================================================
[766]144
[1191]145  INTEGER :: i, k, it
146  INTEGER :: nsplit
[524]147
[1191]148!Sources et Reservoirs de traceurs (ex:Radon):
149!--------------------------------------------
150!
151  REAL,DIMENSION(:,:),ALLOCATABLE,SAVE :: source  ! a voir lorsque le flux de surface est prescrit
152!$OMP THREADPRIVATE(source)
[524]153
[1191]154!
155!Entrees/Sorties: (cf ini_histrac.h et write_histrac.h) 
156!---------------
157  INTEGER                   :: iiq, ierr
158  INTEGER                   :: nhori, nvert
159  REAL                      :: zsto, zout, zjulian
160  INTEGER,SAVE              :: nid_tra     ! pointe vers le fichier histrac.nc         
161!$OMP THREADPRIVATE(nid_tra)
162  REAL,DIMENSION(klon)      :: zx_tmp_fi2d ! variable temporaire grille physique
163  INTEGER                   :: itau_w      ! pas de temps ecriture = nstep + itau_phy
164  LOGICAL,PARAMETER :: ok_sync=.TRUE.
[524]165
[1191]166!
167! Nature du traceur
168!------------------
169  LOGICAL,DIMENSION(:),ALLOCATABLE,SAVE :: aerosol  ! aerosol(it) = true  => aerosol => lessivage
170!$OMP THREADPRIVATE(aerosol)                        ! aerosol(it) = false => gaz
171  REAL,DIMENSION(klon,klev)             :: delp     ! epaisseur de couche (Pa)
172!
173! Tendances de traceurs (Td):
174!------------------------
175!
176  REAL,DIMENSION(klon,klev)      :: d_tr     ! Td dans l'atmosphere
177  REAL,DIMENSION(klon,klev,nbtr) :: d_tr_cl  ! Td couche limite/traceur
178  REAL,DIMENSION(klon,klev,nbtr) :: d_tr_cv  ! Td convection/traceur
179  REAL,DIMENSION(klon,klev,nbtr) :: d_tr_th  ! Td thermique
180  REAL,DIMENSION(klon,klev,nbtr) :: d_tr_lessi_impa ! Td du lessivage par impaction
181  REAL,DIMENSION(klon,klev,nbtr) :: d_tr_lessi_nucl ! Td du lessivage par nucleation
182!
183! Physique
184!----------   
185  REAL,DIMENSION(klon,klev,nbtr) :: flestottr ! flux de lessivage dans chaque couche
186  REAL,DIMENSION(klon,klev)      :: zmasse    ! densité atmosphérique Kg/m2
187  REAL,DIMENSION(klon,klev)      :: ztra_th
188 
189!Controles:
190!---------
191  LOGICAL,SAVE :: couchelimite=.TRUE.
192  LOGICAL,SAVE :: convection=.TRUE.
193  LOGICAL,SAVE :: lessivage
194!$OMP THREADPRIVATE(couchelimite,convection,lessivage)
[766]195
[1191]196  CHARACTER(len=8),DIMENSION(nbtr) :: solsym
[524]197
198
[1191]199!######################################################################
200!                    -- INITIALIZATION --
201!######################################################################
202  IF (debutphy) THEN
203     WRITE(*,*) 'FIRST TIME IN PHYTRAC : pdtphys(sec) = ',pdtphys,'ecrit_tra (sec) = ',ecrit_tra
204     ALLOCATE( source(klon,nbtr), stat=ierr)
205     IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 1',1)
206     
207     ALLOCATE( aerosol(nbtr), stat=ierr)
208     IF (ierr /= 0) CALL abort_gcm('phytrac', 'pb in allocation 2',1)
209     
[524]210
[1191]211     ! Initialize module for specific tracers
212     SELECT CASE(type_trac)
213     CASE('lmdz')
214        CALL traclmdz_init(pctsrf, ftsol, aerosol, lessivage)
215     CASE('inca')
216        source(:,:)=0.
217        CALL tracinca_init(aerosol,lessivage)
218     END SELECT
219!
220! Initialize diagnostic output
221! ----------------------------
[1192]222#ifdef CPP_IOIPSL
[1191]223     INCLUDE "ini_histrac.h"
[1192]224#endif
[1191]225  END IF
226!############################################ END INITIALIZATION #######
[524]227
[1191]228!===============================================================================
229!    -- Do specific treatment according to chemestry model or local LMDZ tracers
230!     
231!===============================================================================
232  SELECT CASE(type_trac)
233  CASE('lmdz')
234     !    -- Traitement des traceurs avec traclmdz
235     
236     CALL traclmdz(&
[1212]237          nstep,    pdtphys,      t_seri,           &
[1191]238          paprs,    pplay,        cdragh,  coefh,   &
239          yu1,      yv1,          ftsol,   pctsrf,  &
240          xlat,     couchelimite,                   &
241          tr_seri,  source,       solsym,  d_tr_cl)
242     
243  CASE('inca')
244     !    -- CHIMIE INCA  config_inca = aero or chem --
[776]245
[1191]246     CALL tracinca(&
247          nstep,    julien,   gmtime,         lafin,     &
248          pdtphys,  t_seri,   paprs,          pplay,     &
249          pmfu,     ftsol,    pctsrf,         pphis,     &
250          pphi,     albsol,   sh,             rh,        &
251          cldfra,   rneb,     diafra,         cldliq,    &
252          itop_con, ibas_con, pmflxr,         pmflxs,    &
253          prfl,     psfl,     aerosol_couple, flxmass_w, &
254          tau_aero, piz_aero, cg_aero,        ccm,       &
255          rfname,                                        &
256          tr_seri,  source,   solsym)     
257  END SELECT
[524]258
[1191]259!======================================================================
260!       -- Calcul de l'effet de la convection --
261!======================================================================
262  IF (convection) THEN
263     DO it=1, nbtr
264        IF ( conv_flg(it) == 0 ) CYCLE
265       
266        IF (iflag_con.LT.2) THEN
267           d_tr_cv(:,:,:)=0.
268        ELSE IF (iflag_con.EQ.2) THEN
269!..Tiedke
270           CALL nflxtr(pdtphys, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &
271                pplay, paprs, tr_seri(:,:,it), d_tr_cv(:,:,it))
272        ELSE
273!..K.Emanuel
274           CALL cvltr(pdtphys, da, phi, mp, paprs,pplay, tr_seri(:,:,it),&
275                upwd,dnwd,d_tr_cv(:,:,it))
276        END IF
[959]277
[1191]278        DO k = 1, klev
279           DO i = 1, klon       
280              tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cv(i,k,it)
281           END DO
282        END DO
[959]283
[1191]284        CALL minmaxqfi(tr_seri(:,:,it),0.,1.e33,'convection it = '//solsym(it))
285             
286     END DO ! nbtr
287  END IF ! convection
[959]288
[524]289!======================================================================
[1191]290!    -- Calcul de l'effet des thermiques --
[524]291!======================================================================
292
[1191]293  DO k=1,klev
294     DO i=1,klon
295        zmasse(i,k)=(paprs(i,k)-paprs(i,k+1))/rg
296     END DO
297  END DO
[524]298
[1191]299  DO it=1,nbtr
300     DO k=1,klev
301        DO i=1,klon
302           d_tr_th(i,k,it)=0.
303           tr_seri(i,k,it)=MAX(tr_seri(i,k,it),0.)
304           tr_seri(i,k,it)=MIN(tr_seri(i,k,it),1.e10)
[524]305        END DO
[1191]306     END DO
307  END DO
308 
309  IF (iflag_thermals.GT.0) THEN   
310     nsplit=10
311     DO it=1, nbtr
312        DO isplit=1,nsplit
[524]313
[1191]314           CALL dqthermcell(klon,klev,pdtphys/nsplit, &
315                fm_therm,entr_therm,zmasse, &
316                tr_seri(1:klon,1:klev,it),d_tr,ztra_th)
[524]317
[1191]318           DO k=1,klev
319              DO i=1,klon
320                 d_tr(i,k)=pdtphys*d_tr(i,k)/nsplit
321                 d_tr_th(i,k,it)=d_tr_th(i,k,it)+d_tr(i,k)
322                 tr_seri(i,k,it)=MAX(tr_seri(i,k,it)+d_tr(i,k),0.)
323              END DO
324           END DO
325        END DO ! nsplit
326     END DO ! it
327  END IF ! Thermiques
[682]328
[1191]329!======================================================================
330!     -- Calcul de l'effet de la couche limite --
331!======================================================================
[682]332
[1191]333  IF (couchelimite) THEN
[524]334
[1191]335     DO k = 1, klev
336        DO i = 1, klon
337           delp(i,k) = paprs(i,k)-paprs(i,k+1)
338        END DO
339     END DO
[524]340
[1191]341     DO it=1, nbtr
342       
343        IF( pbl_flg(it) /= 0 ) THEN
344       
345           CALL cltrac(pdtphys, coefh,t_seri,       &
346                tr_seri(:,:,it), source(:,it),      &
347                paprs, pplay, delp,                 &
348                d_tr_cl(:,:,it))
349           
350           DO k = 1, klev
351              DO i = 1, klon
352                 tr_seri(i,k,it) = tr_seri(i,k,it) + d_tr_cl(i,k,it)
353              END DO
354           END DO
355        END IF
[959]356
[1191]357     END DO
358     
359  END IF ! couche limite
[959]360
[619]361
[1191]362!======================================================================
363!   Calcul de l'effet de la precipitation
364!======================================================================
[959]365
[1191]366  IF (lessivage) THEN
367     
368     d_tr_lessi_nucl(:,:,:) = 0.
369     d_tr_lessi_impa(:,:,:) = 0.
370     flestottr(:,:,:) = 0.
371!=========================
372! LESSIVAGE LARGE SCALE :
373!=========================
[524]374
[1191]375! Tendance des aerosols nuclees et impactes
376! -----------------------------------------
377     DO it = 1, nbtr
378        IF (aerosol(it)) THEN
379           DO k = 1, klev
[524]380              DO i = 1, klon
[1191]381                 d_tr_lessi_nucl(i,k,it) = d_tr_lessi_nucl(i,k,it) +    &
382                      ( 1 - frac_nucl(i,k) )*tr_seri(i,k,it)
383                 d_tr_lessi_impa(i,k,it) = d_tr_lessi_impa(i,k,it) +    &
384                      ( 1 - frac_impa(i,k) )*tr_seri(i,k,it)
[524]385
[1191]386!
387! Flux lessivage total
388! ------------------------------------------------------------
389                 flestottr(i,k,it) = flestottr(i,k,it) -   &
390                      ( d_tr_lessi_nucl(i,k,it)   +        &
391                      d_tr_lessi_impa(i,k,it) ) *          &
392                      ( paprs(i,k)-paprs(i,k+1) ) /        &
393                      (RG * pdtphys)
394!
395! Mise a jour des traceurs due a l'impaction,nucleation
396! ----------------------------------------------------------------------
397                 tr_seri(i,k,it)=tr_seri(i,k,it)*frac_impa(i,k)*frac_nucl(i,k)
398              END DO
399           END DO
400        END IF
401     END DO
402     
403  END IF ! lessivage
[524]404
[1191]405!=============================================================
406!   Ecriture des sorties
407!=============================================================
[1192]408#ifdef CPP_IOIPSL
[1191]409  INCLUDE "write_histrac.h"
[1192]410#endif
[524]411
[1191]412END SUBROUTINE phytrac
Note: See TracBrowser for help on using the repository browser.