source: LMDZ4/branches/LMDZ4V5.0-LF/libf/phylmd/phytrac.F90 @ 3603

Last change on this file since 3603 was 1299, checked in by Laurent Fairhead, 15 years ago

Nettoyage general pour se rapprocher des normes et éviter des erreurs a la
compilation:

  • tous les FLOAT() sont remplacés par des REAL()
  • tous les STOP dans phylmd sont remplacés par des appels à abort_gcm
  • le common control défini dans le fichier control.h est remplacé par le module control_mod pour éviter des messages sur l'alignement des variables dans les déclarations
  • des $Header$ remplacés par des $Id$ pour svn

Quelques remplacements à faire ont pu m'échapper


General cleanup of the code to try and adhere to norms and to prevent some
compilation errors:

  • all FLOAT() instructions have been replaced by REAL() instructions
  • all STOP instructions in phylmd have been replaced by calls to abort_gcm
  • the common block control defined in the control.h file has been replaced by the control_mod to prevent compilation warnings on the alignement of declared variables
  • $Header$ replaced by $Id$ for svn

Some changes which should have been made might have escaped me

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