source: codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/config/ppsrc/dyn3d_common/advx.f @ 224

Last change on this file since 224 was 224, checked in by ymipsl, 10 years ago
File size: 17.0 KB
Line 
1!
2! $Header$
3!
4      SUBROUTINE  advx(limit,dtx,pbaru,sm,s0,
5     $     sx,sy,sz,lati,latf)
6      IMPLICIT NONE
7
8CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
9C                                                                C
10C  first-order moments (FOM) advection of tracer in X direction  C
11C                                                                C
12C  Source : Pascal Simon (Meteo,CNRM)                            C
13C  Adaptation : A.Armengaud (LGGE) juin 94                       C
14C                                                                C
15C  limit,dtx,pbaru,pbarv,sm,s0,sx,sy,sz                       C
16C  sont des arguments d'entree pour le s-pg...                   C
17C                                                                C
18C  sm,s0,sx,sy,sz                                                C
19C  sont les arguments de sortie pour le s-pg                     C
20C                                                                C
21CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
22C
23C  parametres principaux du modele
24C
25!-----------------------------------------------------------------------
26!   INCLUDE 'dimensions.h'
27!
28!   dimensions.h contient les dimensions du modele
29!   ndm est tel que iim=2**ndm
30!-----------------------------------------------------------------------
31
32      INTEGER iim,jjm,llm,ndm
33
34      PARAMETER (iim= 128,jjm=96,llm=64,ndm=1)
35
36!-----------------------------------------------------------------------
37!
38! $Header$
39!
40!
41!  ATTENTION!!!!: ce fichier include est compatible format fixe/format libre
42!                 veillez  n'utiliser que des ! pour les commentaires
43!                 et  bien positionner les & des lignes de continuation
44!                 (les placer en colonne 6 et en colonne 73)
45!
46!
47!-----------------------------------------------------------------------
48!   INCLUDE 'paramet.h'
49
50      INTEGER  iip1,iip2,iip3,jjp1,llmp1,llmp2,llmm1
51      INTEGER  kftd,ip1jm,ip1jmp1,ip1jmi1,ijp1llm
52      INTEGER  ijmllm,mvar
53      INTEGER jcfil,jcfllm
54
55      PARAMETER( iip1= iim+1,iip2=iim+2,iip3=iim+3                       &
56     &    ,jjp1=jjm+1-1/jjm)
57      PARAMETER( llmp1 = llm+1,  llmp2 = llm+2, llmm1 = llm-1 )
58      PARAMETER( kftd  = iim/2 -ndm )
59      PARAMETER( ip1jm  = iip1*jjm,  ip1jmp1= iip1*jjp1 )
60      PARAMETER( ip1jmi1= ip1jm - iip1 )
61      PARAMETER( ijp1llm= ip1jmp1 * llm, ijmllm= ip1jm * llm )
62      PARAMETER( mvar= ip1jmp1*( 2*llm+1) + ijmllm )
63      PARAMETER( jcfil=jjm/2+5, jcfllm=jcfil*llm )
64
65!-----------------------------------------------------------------------
66!
67! $Id: comconst.h 1437 2010-09-30 08:29:10Z emillour $
68!
69!-----------------------------------------------------------------------
70! INCLUDE comconst.h
71
72      COMMON/comconsti/im,jm,lllm,imp1,jmp1,lllmm1,lllmp1,lcl,          &
73     &                 iflag_top_bound,mode_top_bound
74      COMMON/comconstr/dtvr,daysec,                                     &
75     & pi,dtphys,dtdiss,rad,r,kappa,cotot,unsim,g,omeg                  &
76     & ,dissip_fac_mid,dissip_fac_up,dissip_deltaz,dissip_hdelta        &
77     & ,dissip_pupstart  ,tau_top_bound,                                &
78     & daylen,molmass, ihf
79      COMMON/cpdetvenus/cpp,nu_venus,t0_venus
80
81      INTEGER im,jm,lllm,imp1,jmp1,lllmm1,lllmp1,lcl
82      REAL dtvr ! dynamical time step (in s)
83      REAL daysec !length (in s) of a standard day
84      REAL pi    ! something like 3.14159....
85      REAL dtphys ! (s) time step for the physics
86      REAL dtdiss ! (s) time step for the dissipation
87      REAL rad ! (m) radius of the planet
88      REAL r ! Reduced Gas constant r=R/mu
89             ! with R=8.31.. J.K-1.mol-1, mu: mol mass of atmosphere (kg/mol)
90      REAL cpp   ! Cp
91      REAL kappa ! kappa=R/Cp
92      REAL cotot
93      REAL unsim ! = 1./iim
94      REAL g ! (m/s2) gravity
95      REAL omeg ! (rad/s) rotation rate of the planet
96! Dissipation factors, for Earth model:
97      REAL dissip_factz,dissip_zref !dissip_deltaz
98! Dissipation factors, for other planets:
99      REAL dissip_fac_mid,dissip_fac_up,dissip_deltaz,dissip_hdelta
100      REAL dissip_pupstart
101      INTEGER iflag_top_bound,mode_top_bound
102      REAL tau_top_bound
103      REAL daylen ! length of solar day, in 'standard' day length
104      REAL molmass ! (g/mol) molar mass of the atmosphere
105
106      REAL nu_venus,t0_venus ! coeffs needed for Cp(T), Venus atmosphere
107      REAL ihf  ! (W/m2) intrinsic heat flux for giant planets
108
109
110!-----------------------------------------------------------------------
111!
112! $Id: comvert.h 1654 2012-09-24 15:07:18Z aslmd $
113!
114!-----------------------------------------------------------------------
115!   INCLUDE 'comvert.h'
116
117      COMMON/comvertr/ap(llm+1),bp(llm+1),presnivs(llm),dpres(llm),     &
118     &               pa,preff,nivsigs(llm),nivsig(llm+1),               &
119     &               aps(llm),bps(llm),scaleheight,pseudoalt(llm)
120
121      common/comverti/disvert_type, pressure_exner
122
123      real ap     ! hybrid pressure contribution at interlayers
124      real bp     ! hybrid sigma contribution at interlayer
125      real presnivs ! (reference) pressure at mid-layers
126      real dpres
127      real pa     ! reference pressure (Pa) at which hybrid coordinates
128                  ! become purely pressure
129      real preff  ! reference surface pressure (Pa)
130      real nivsigs
131      real nivsig
132      real aps    ! hybrid pressure contribution at mid-layers
133      real bps    ! hybrid sigma contribution at mid-layers
134      real scaleheight ! atmospheric (reference) scale height (km)
135      real pseudoalt ! pseudo-altitude of model levels (km), based on presnivs(),
136                     ! preff and scaleheight
137
138      integer disvert_type ! type of vertical discretization:
139                           ! 1: Earth (default for planet_type==earth),
140                           !     automatic generation
141                           ! 2: Planets (default for planet_type!=earth),
142                           !     using 'z2sig.def' (or 'esasig.def) file
143
144      logical pressure_exner
145!     compute pressure inside layers using Exner function, else use mean
146!     of pressure values at interfaces
147
148 !-----------------------------------------------------------------------
149
150C  Arguments :
151C  -----------
152C  dtx : frequence fictive d'appel du transport
153C  pbaru, pbarv : flux de masse en x et y en Pa.m2.s-1
154
155       INTEGER ntra
156       PARAMETER (ntra = 1)
157
158C ATTENTION partout ou on trouve ntra, insertion de boucle
159C           possible dans l'avenir.
160
161      REAL dtx
162      REAL pbaru ( iip1,jjp1,llm )
163
164C  moments: SM  total mass in each grid box
165C           S0  mass of tracer in each grid box
166C           Si  1rst order moment in i direction
167C
168      REAL SM(iip1,jjp1,llm),S0(iip1,jjp1,llm,ntra)
169      REAL sx(iip1,jjp1,llm,ntra)
170     $    ,sy(iip1,jjp1,llm,ntra)
171      REAL sz(iip1,jjp1,llm,ntra)
172
173C  Local :
174C  ------- 
175
176C  mass fluxes across the boundaries (UGRI,VGRI,WGRI)
177C  mass fluxes in kg
178C  declaration :
179
180      REAL UGRI(iip1,jjp1,llm)
181
182C  Rem : VGRI et WGRI ne sont pas utilises dans 
183C  cette subroutine ( advection en x uniquement )
184C
185C  Ti are the moments for the current latitude and level
186C
187      REAL TM(iim)
188      REAL T0(iim,ntra),TX(iim,ntra)
189      REAL TY(iim,ntra),TZ(iim,ntra)
190      REAL TEMPTM                ! just a temporary variable
191C
192C  the moments F are similarly defined and used as temporary
193C  storage for portions of the grid boxes in transit
194C
195      REAL FM(iim)
196      REAL F0(iim,ntra),FX(iim,ntra)
197      REAL FY(iim,ntra),FZ(iim,ntra)
198C
199C  work arrays
200C
201      REAL ALF(iim),ALF1(iim),ALFQ(iim),ALF1Q(iim)
202C
203      REAL SMNEW(iim),UEXT(iim)
204C
205      REAL sqi,sqf
206
207      LOGICAL LIMIT
208      INTEGER NUM(jjp1),LONK,NUMK
209      INTEGER lon,lati,latf,niv
210      INTEGER i,i2,i3,j,jv,l,k,itrac 
211
212      lon = iim 
213      niv = llm 
214
215C *** Test de passage d'arguments ******
216
217
218C  -------------------------------------
219      DO 300 j = 1,jjp1
220         NUM(j) = 1
221  300 CONTINUE
222      sqi = 0.
223      sqf = 0.
224
225      DO l = 1,llm
226         DO j = 1,jjp1
227            DO i = 1,iim
228cIM 240305            sqi = sqi + S0(i,j,l,9)
229               sqi = sqi + S0(i,j,l,ntra)
230            ENDDO
231         ENDDO
232      ENDDO
233      PRINT*,'-------- DIAG DANS ADVX - ENTREE ---------'
234      PRINT*,'sqi=',sqi
235
236
237C  Interface : adaptation nouveau modele
238C  -------------------------------------
239C
240C  ---------------------------------------------------------
241C  Conversion des flux de masses en kg/s
242C  pbaru est en N/s d'ou :
243C  ugri est en kg/s
244
245      DO 500 l = 1,llm
246         DO 500 j = 1,jjm+1
247            DO 500 i = 1,iip1 
248C            ugri (i,j,llm+1-l) = pbaru (i,j,l) * ( dsig(l) / g )
249             ugri (i,j,llm+1-l) = pbaru (i,j,l)
250  500 CONTINUE
251
252
253C  ---------------------------------------------------------
254C  ---------------------------------------------------------
255C  ---------------------------------------------------------
256 
257C  start here         
258C
259C  boucle principale sur les niveaux et les latitudes
260C
261      DO 1 L=1,NIV
262      DO 1 K=lati,latf
263C
264C  initialisation
265C
266C  program assumes periodic boundaries in X
267C
268      DO 10 I=2,LON
269         SMNEW(I)=SM(I,K,L)+(UGRI(I-1,K,L)-UGRI(I,K,L))*DTX
270 10   CONTINUE
271      SMNEW(1)=SM(1,K,L)+(UGRI(LON,K,L)-UGRI(1,K,L))*DTX
272C
273C  modifications for extended polar zones
274C
275      NUMK=NUM(K)
276      LONK=LON/NUMK
277C
278      IF(NUMK.GT.1) THEN
279C
280      DO 111 I=1,LON
281         TM(I)=0.
282 111  CONTINUE
283      DO 112 JV=1,NTRA
284      DO 1120 I=1,LON
285         T0(I,JV)=0.
286         TX(I,JV)=0.
287         TY(I,JV)=0.
288         TZ(I,JV)=0.
289 1120 CONTINUE
290 112  CONTINUE
291C
292      DO 11 I2=1,NUMK
293C
294         DO 113 I=1,LONK
295            I3=(I-1)*NUMK+I2
296            TM(I)=TM(I)+SM(I3,K,L)
297            ALF(I)=SM(I3,K,L)/TM(I)
298            ALF1(I)=1.-ALF(I)
299 113     CONTINUE
300C
301         DO  JV=1,NTRA
302         DO  I=1,LONK
303            I3=(I-1)*NUMK+I2
304            TEMPTM=-ALF(I)*T0(I,JV)+ALF1(I)
305     $          *S0(I3,K,L,JV)
306            T0(I,JV)=T0(I,JV)+S0(I3,K,L,JV)
307            TX(I,JV)=ALF(I)  *sx(I3,K,L,JV)+
308     $       ALF1(I)*TX(I,JV) +3.*TEMPTM
309            TY(I,JV)=TY(I,JV)+sy(I3,K,L,JV)
310            TZ(I,JV)=TZ(I,JV)+sz(I3,K,L,JV)
311         ENDDO 
312         ENDDO
313C
314 11   CONTINUE
315C
316      ELSE
317C
318      DO 115 I=1,LON
319         TM(I)=SM(I,K,L)
320 115  CONTINUE
321      DO 116 JV=1,NTRA
322      DO 1160 I=1,LON
323         T0(I,JV)=S0(I,K,L,JV)
324         TX(I,JV)=sx(I,K,L,JV)
325         TY(I,JV)=sy(I,K,L,JV)
326         TZ(I,JV)=sz(I,K,L,JV)
327 1160 CONTINUE
328 116  CONTINUE
329C
330      ENDIF
331C
332      DO 117 I=1,LONK
333         UEXT(I)=UGRI(I*NUMK,K,L)
334 117  CONTINUE
335C
336C  place limits on appropriate moments before transport
337C      (if flux-limiting is to be applied)
338C
339      IF(.NOT.LIMIT) GO TO 13
340C
341      DO 12 JV=1,NTRA
342      DO 120 I=1,LONK
343        TX(I,JV)=SIGN(AMIN1(AMAX1(T0(I,JV),0.),ABS(TX(I,JV))),TX(I,JV))
344 120  CONTINUE
345 12   CONTINUE
346C
347 13   CONTINUE
348C
349C  calculate flux and moments between adjacent boxes
350C  1- create temporary moments/masses for partial boxes in transit
351C  2- reajusts moments remaining in the box
352C
353C  flux from IP to I if U(I).lt.0
354C
355      DO 140 I=1,LONK-1
356         IF(UEXT(I).LT.0.) THEN
357           FM(I)=-UEXT(I)*DTX
358           ALF(I)=FM(I)/TM(I+1)
359           TM(I+1)=TM(I+1)-FM(I)
360         ENDIF
361 140  CONTINUE
362C
363      I=LONK
364      IF(UEXT(I).LT.0.) THEN
365        FM(I)=-UEXT(I)*DTX
366        ALF(I)=FM(I)/TM(1)
367        TM(1)=TM(1)-FM(I)
368      ENDIF
369C
370C  flux from I to IP if U(I).gt.0
371C
372      DO 141 I=1,LONK
373         IF(UEXT(I).GE.0.) THEN
374           FM(I)=UEXT(I)*DTX
375           ALF(I)=FM(I)/TM(I)
376           TM(I)=TM(I)-FM(I)
377         ENDIF
378 141  CONTINUE
379C
380      DO 142 I=1,LONK
381         ALFQ(I)=ALF(I)*ALF(I)
382         ALF1(I)=1.-ALF(I)
383         ALF1Q(I)=ALF1(I)*ALF1(I)
384 142  CONTINUE
385C
386      DO 150 JV=1,NTRA
387      DO 1500 I=1,LONK-1
388C
389         IF(UEXT(I).LT.0.) THEN
390C
391           F0(I,JV)=ALF (I)* ( T0(I+1,JV)-ALF1(I)*TX(I+1,JV) )
392           FX(I,JV)=ALFQ(I)*TX(I+1,JV)
393           FY(I,JV)=ALF (I)*TY(I+1,JV)
394           FZ(I,JV)=ALF (I)*TZ(I+1,JV)
395C
396           T0(I+1,JV)=T0(I+1,JV)-F0(I,JV)
397           TX(I+1,JV)=ALF1Q(I)*TX(I+1,JV)
398           TY(I+1,JV)=TY(I+1,JV)-FY(I,JV)
399           TZ(I+1,JV)=TZ(I+1,JV)-FZ(I,JV)
400C
401         ENDIF
402C
403 1500 CONTINUE
404 150  CONTINUE
405C
406      I=LONK
407      IF(UEXT(I).LT.0.) THEN
408C
409        DO 151 JV=1,NTRA
410C
411           F0 (I,JV)=ALF (I)* ( T0(1,JV)-ALF1(I)*TX(1,JV) )
412           FX (I,JV)=ALFQ(I)*TX(1,JV)
413           FY (I,JV)=ALF (I)*TY(1,JV)
414           FZ (I,JV)=ALF (I)*TZ(1,JV)
415C
416           T0(1,JV)=T0(1,JV)-F0(I,JV)
417           TX(1,JV)=ALF1Q(I)*TX(1,JV)
418           TY(1,JV)=TY(1,JV)-FY(I,JV)
419           TZ(1,JV)=TZ(1,JV)-FZ(I,JV)
420C
421 151    CONTINUE
422C
423      ENDIF
424C
425      DO 152 JV=1,NTRA
426      DO 1520 I=1,LONK
427C
428         IF(UEXT(I).GE.0.) THEN
429C
430           F0(I,JV)=ALF (I)* ( T0(I,JV)+ALF1(I)*TX(I,JV) )
431           FX(I,JV)=ALFQ(I)*TX(I,JV)
432           FY(I,JV)=ALF (I)*TY(I,JV)
433           FZ(I,JV)=ALF (I)*TZ(I,JV)
434C
435           T0(I,JV)=T0(I,JV)-F0(I,JV)
436           TX(I,JV)=ALF1Q(I)*TX(I,JV)
437           TY(I,JV)=TY(I,JV)-FY(I,JV)
438           TZ(I,JV)=TZ(I,JV)-FZ(I,JV)
439C
440         ENDIF
441C
442 1520 CONTINUE
443 152  CONTINUE
444C
445C  puts the temporary moments Fi into appropriate neighboring boxes
446C
447      DO 160 I=1,LONK
448         IF(UEXT(I).LT.0.) THEN
449           TM(I)=TM(I)+FM(I)
450           ALF(I)=FM(I)/TM(I)
451         ENDIF
452 160  CONTINUE
453C
454      DO 161 I=1,LONK-1
455         IF(UEXT(I).GE.0.) THEN
456           TM(I+1)=TM(I+1)+FM(I)
457           ALF(I)=FM(I)/TM(I+1)
458         ENDIF
459 161  CONTINUE
460C
461      I=LONK
462      IF(UEXT(I).GE.0.) THEN
463        TM(1)=TM(1)+FM(I)
464        ALF(I)=FM(I)/TM(1)
465      ENDIF
466C
467      DO 162 I=1,LONK
468         ALF1(I)=1.-ALF(I)
469 162  CONTINUE
470C
471      DO 170 JV=1,NTRA
472      DO 1700 I=1,LONK
473C
474         IF(UEXT(I).LT.0.) THEN
475C
476           TEMPTM=-ALF(I)*T0(I,JV)+ALF1(I)*F0(I,JV)
477           T0(I,JV)=T0(I,JV)+F0(I,JV)
478           TX(I,JV)=ALF(I)*FX(I,JV)+ALF1(I)*TX(I,JV)+3.*TEMPTM
479           TY(I,JV)=TY(I,JV)+FY(I,JV)
480           TZ(I,JV)=TZ(I,JV)+FZ(I,JV)
481C
482         ENDIF
483C
484 1700 CONTINUE
485 170  CONTINUE
486C
487      DO 171 JV=1,NTRA
488      DO 1710 I=1,LONK-1
489C
490         IF(UEXT(I).GE.0.) THEN
491C
492           TEMPTM=ALF(I)*T0(I+1,JV)-ALF1(I)*F0(I,JV)
493           T0(I+1,JV)=T0(I+1,JV)+F0(I,JV)
494           TX(I+1,JV)=ALF(I)*FX(I,JV)+ALF1(I)*TX(I+1,JV)+3.*TEMPTM
495           TY(I+1,JV)=TY(I+1,JV)+FY(I,JV)
496           TZ(I+1,JV)=TZ(I+1,JV)+FZ(I,JV)
497C
498         ENDIF
499C
500 1710 CONTINUE
501 171  CONTINUE
502C
503      I=LONK
504      IF(UEXT(I).GE.0.) THEN
505        DO 172 JV=1,NTRA
506           TEMPTM=ALF(I)*T0(1,JV)-ALF1(I)*F0(I,JV)
507           T0(1,JV)=T0(1,JV)+F0(I,JV)
508           TX(1,JV)=ALF(I)*FX(I,JV)+ALF1(I)*TX(1,JV)+3.*TEMPTM
509           TY(1,JV)=TY(1,JV)+FY(I,JV)
510           TZ(1,JV)=TZ(1,JV)+FZ(I,JV)
511 172    CONTINUE
512      ENDIF
513C
514C  retour aux mailles d'origine (passage des Tij aux Sij)
515C
516      IF(NUMK.GT.1) THEN
517C
518      DO 180 I2=1,NUMK
519C
520         DO 180 I=1,LONK
521C
522            I3=I2+(I-1)*NUMK
523            SM(I3,K,L)=SMNEW(I3)
524            ALF(I)=SMNEW(I3)/TM(I)
525            TM(I)=TM(I)-SMNEW(I3)
526C
527            ALFQ(I)=ALF(I)*ALF(I)
528            ALF1(I)=1.-ALF(I)
529            ALF1Q(I)=ALF1(I)*ALF1(I)
530C
531 180     CONTINUE
532C
533         DO  JV=1,NTRA
534         DO  I=1,LONK
535C
536            I3=I2+(I-1)*NUMK
537            S0(I3,K,L,JV)=ALF (I)
538     $       * (T0(I,JV)-ALF1(I)*TX(I,JV))
539            sx(I3,K,L,JV)=ALFQ(I)*TX(I,JV)
540            sy(I3,K,L,JV)=ALF (I)*TY(I,JV)
541            sz(I3,K,L,JV)=ALF (I)*TZ(I,JV)
542C
543C   reajusts moments remaining in the box
544C
545            T0(I,JV)=T0(I,JV)-S0(I3,K,L,JV)
546            TX(I,JV)=ALF1Q(I)*TX(I,JV)
547            TY(I,JV)=TY(I,JV)-sy(I3,K,L,JV)
548            TZ(I,JV)=TZ(I,JV)-sz(I3,K,L,JV)
549          ENDDO
550          ENDDO
551C
552C
553      ELSE
554C
555      DO 190 I=1,LON
556         SM(I,K,L)=TM(I)
557 190  CONTINUE
558      DO 191 JV=1,NTRA
559      DO 1910 I=1,LON
560         S0(I,K,L,JV)=T0(I,JV)
561         sx(I,K,L,JV)=TX(I,JV)
562         sy(I,K,L,JV)=TY(I,JV)
563         sz(I,K,L,JV)=TZ(I,JV)
564 1910 CONTINUE
565 191  CONTINUE
566C
567      ENDIF
568C
569 1    CONTINUE
570C
571C ----------- AA Test en fin de ADVX ------ Controle des S*
572c OK
573c      DO 9998 l = 1, llm
574c      DO 9998 j = 1, jjp1
575c      DO 9998 i = 1, iip1
576c         IF (S0(i,j,l,ntra).lt.0..and.LIMIT) THEN
577c            PRINT*, '-------------------'
578c            PRINT*, 'En fin de ADVX'
579c            PRINT*,'SM(',i,j,l,')=',SM(i,j,l)
580c            PRINT*,'S0(',i,j,l,')=',S0(i,j,l,ntra)
581c            print*, 'sx(',i,j,l,')=',sx(i,j,l,ntra)
582c            print*, 'sy(',i,j,l,')=',sy(i,j,l,ntra)
583c            print*, 'sz(',i,j,l,')=',sz(i,j,l,ntra)
584c            WRITE (*,*) 'On arrete !! - pbl en fin de ADVX1'
585cc            STOP
586c         ENDIF
587c 9998 CONTINUE
588c
589C ---------- bouclage cyclique 
590      DO itrac=1,ntra
591      DO l = 1,llm
592        DO j = lati,latf
593           SM(iip1,j,l) = SM(1,j,l)
594           S0(iip1,j,l,itrac) = S0(1,j,l,itrac)
595           sx(iip1,j,l,itrac) = sx(1,j,l,itrac)
596           sy(iip1,j,l,itrac) = sy(1,j,l,itrac)
597           sz(iip1,j,l,itrac) = sz(1,j,l,itrac)
598        END DO
599      END DO
600      ENDDO 
601
602c ----------- qqtite totale de traceur dans tte l'atmosphere
603      DO l = 1, llm
604        DO j = 1, jjp1
605          DO i = 1, iim
606cIM 240405          sqf = sqf + S0(i,j,l,9)
607             sqf = sqf + S0(i,j,l,ntra)
608          END DO 
609        END DO
610      END DO
611c
612      PRINT*,'------ DIAG DANS ADVX - SORTIE -----'
613      PRINT*,'sqf=',sqf
614c-------------
615
616      RETURN
617      END
618C_________________________________________________________________
619C_________________________________________________________________
Note: See TracBrowser for help on using the repository browser.