/[lmdze]/trunk/dyn3d/advyp.f
ViewVC logotype

Diff of /trunk/dyn3d/advyp.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

trunk/libf/dyn3d/advyp.f revision 3 by guez, Wed Feb 27 13:16:39 2008 UTC trunk/dyn3d/advyp.f revision 113 by guez, Thu Sep 18 19:56:46 2014 UTC
# Line 1  Line 1 
1  !  
2  ! $Header: /home/cvsroot/LMDZ4/libf/dyn3d/advyp.F,v 1.1.1.1 2004/05/19 12:53:06 lmdzadmin Exp $  ! $Header: /home/cvsroot/LMDZ4/libf/dyn3d/advyp.F,v 1.1.1.1 2004/05/19
3  !  ! 12:53:06 lmdzadmin Exp $
4        SUBROUTINE ADVYP(LIMIT,DTY,PBARV,SM,S0,SSX,SY,SZ  
5       .                 ,SSXX,SSXY,SSXZ,SYY,SYZ,SZZ,ntra )  SUBROUTINE advyp(limit, dty, pbarv, sm, s0, ssx, sy, sz, ssxx, ssxy, ssxz, &
6        use dimens_m      syy, syz, szz, ntra)
7        use comconst    USE dimens_m
8        use paramet_m    USE comconst
9        use comvert    USE paramet_m
10        use comgeom    USE disvert_m
11        IMPLICIT NONE    USE comgeom
12  CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC    IMPLICIT NONE
13  C                                                                 C    ! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
14  C  second-order moments (SOM) advection of tracer in Y direction  C    ! C
15  C                                                                 C    ! second-order moments (SOM) advection of tracer in Y direction  C
16  CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC    ! C
17  C                                                                C    ! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
18  C  Source : Pascal Simon ( Meteo, CNRM )                         C    ! C
19  C  Adaptation : A.A. (LGGE)                                      C    ! Source : Pascal Simon ( Meteo, CNRM )                        C
20  C  Derniere Modif : 19/10/95 LAST    ! Adaptation : A.A. (LGGE)                                     C
21  C                                                                C    ! Derniere Modif : 19/10/95 LAST
22  C  sont les arguments d'entree pour le s-pg                      C    ! C
23  C                                                                C    ! sont les arguments d'entree pour le s-pg                     C
24  C  argument de sortie du s-pg                                    C    ! C
25  C                                                                C    ! argument de sortie du s-pg                                   C
26  CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC    ! C
27  CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC    ! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
28  C    ! CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
29  C  Rem : Probleme aux poles il faut reecrire ce cas specifique  
30  C        Attention au sens de l'indexation    ! Rem : Probleme aux poles il faut reecrire ce cas specifique
31  C    ! Attention au sens de l'indexation
32  C  parametres principaux du modele  
33  C    ! parametres principaux du modele
34  C  
35    
36  C  Arguments :  
37  C  ----------    ! Arguments :
38  C  dty : frequence fictive d'appel du transport    ! ----------
39  C  parbu,pbarv : flux de masse en x et y en Pa.m2.s-1    ! dty : frequence fictive d'appel du transport
40      ! parbu,pbarv : flux de masse en x et y en Pa.m2.s-1
41        INTEGER lon,lat,niv  
42        INTEGER i,j,jv,k,kp,l    INTEGER lon, lat, niv
43        INTEGER ntra    INTEGER i, j, jv, k, kp, l
44  C      PARAMETER (ntra = 1)    INTEGER ntra
45      ! PARAMETER (ntra = 1)
46        REAL dty  
47        REAL pbarv ( iip1,jjm, llm )    REAL dty
48      REAL, INTENT (IN) :: pbarv(iip1, jjm, llm)
49  C  moments: SM  total mass in each grid box  
50  C           S0  mass of tracer in each grid box    ! moments: SM  total mass in each grid box
51  C           Si  1rst order moment in i direction    ! S0  mass of tracer in each grid box
52  C    ! Si  1rst order moment in i direction
53        REAL SM(iip1,jjp1,llm)  
54       +    ,S0(iip1,jjp1,llm,ntra)    REAL sm(iip1, jjp1, llm), s0(iip1, jjp1, llm, ntra)
55        REAL SSX(iip1,jjp1,llm,ntra)    REAL ssx(iip1, jjp1, llm, ntra), sy(iip1, jjp1, llm, ntra), &
56       +    ,SY(iip1,jjp1,llm,ntra)      sz(iip1, jjp1, llm, ntra), ssxx(iip1, jjp1, llm, ntra), &
57       +    ,SZ(iip1,jjp1,llm,ntra)      ssxy(iip1, jjp1, llm, ntra), ssxz(iip1, jjp1, llm, ntra), &
58       +    ,SSXX(iip1,jjp1,llm,ntra)      syy(iip1, jjp1, llm, ntra), syz(iip1, jjp1, llm, ntra), &
59       +    ,SSXY(iip1,jjp1,llm,ntra)      szz(iip1, jjp1, llm, ntra)
60       +    ,SSXZ(iip1,jjp1,llm,ntra)  
61       +    ,SYY(iip1,jjp1,llm,ntra)    ! Local :
62       +    ,SYZ(iip1,jjp1,llm,ntra)    ! -------
63       +    ,SZZ(iip1,jjp1,llm,ntra)  
64  C    ! mass fluxes across the boundaries (UGRI,VGRI,WGRI)
65  C  Local :    ! mass fluxes in kg
66  C  -------    ! declaration :
67    
68  C  mass fluxes across the boundaries (UGRI,VGRI,WGRI)    REAL vgri(iip1, 0:jjp1, llm)
69  C  mass fluxes in kg  
70  C  declaration :    ! Rem : UGRI et WGRI ne sont pas utilises dans
71      ! cette subroutine ( advection en y uniquement )
72        REAL VGRI(iip1,0:jjp1,llm)    ! Rem 2 :le dimensionnement de VGRI depend de celui de pbarv
73    
74  C  Rem : UGRI et WGRI ne sont pas utilises dans    ! the moments F are similarly defined and used as temporary
75  C  cette subroutine ( advection en y uniquement )    ! storage for portions of the grid boxes in transit
76  C  Rem 2 :le dimensionnement de VGRI depend de celui de pbarv  
77  C    ! the moments Fij are used as temporary storage for
78  C  the moments F are similarly defined and used as temporary    ! portions of the grid boxes in transit at the current level
79  C  storage for portions of the grid boxes in transit  
80  C    ! work arrays
81  C  the moments Fij are used as temporary storage for  
82  C  portions of the grid boxes in transit at the current level  
83  C    REAL f0(iim, 0:jjp1, ntra), fm(iim, 0:jjp1)
84  C  work arrays    REAL fx(iim, jjm, ntra), fy(iim, jjm, ntra)
85  C    REAL fz(iim, jjm, ntra)
86  C    REAL fxx(iim, jjm, ntra), fxy(iim, jjm, ntra)
87        REAL F0(iim,0:jjp1,ntra),FM(iim,0:jjp1)    REAL fxz(iim, jjm, ntra), fyy(iim, jjm, ntra)
88        REAL FX(iim,jjm,ntra),FY(iim,jjm,ntra)    REAL fyz(iim, jjm, ntra), fzz(iim, jjm, ntra)
89        REAL FZ(iim,jjm,ntra)    REAL s00(ntra)
90        REAL FXX(iim,jjm,ntra),FXY(iim,jjm,ntra)    REAL sm0 ! Just temporal variable
91        REAL FXZ(iim,jjm,ntra),FYY(iim,jjm,ntra)  
92        REAL FYZ(iim,jjm,ntra),FZZ(iim,jjm,ntra)    ! work arrays
93        REAL S00(ntra)  
94        REAL SM0             ! Just temporal variable    REAL alf(iim, 0:jjp1), alf1(iim, 0:jjp1)
95  C    REAL alfq(iim, 0:jjp1), alf1q(iim, 0:jjp1)
96  C  work arrays    REAL alf2(iim, 0:jjp1), alf3(iim, 0:jjp1)
97  C    REAL alf4(iim, 0:jjp1)
98        REAL ALF(iim,0:jjp1),ALF1(iim,0:jjp1)    REAL temptm ! Just temporal variable
99        REAL ALFQ(iim,0:jjp1),ALF1Q(iim,0:jjp1)    REAL slpmax, s1max, s1new, s2new
100        REAL ALF2(iim,0:jjp1),ALF3(iim,0:jjp1)  
101        REAL ALF4(iim,0:jjp1)    ! Special pour poles
102        REAL TEMPTM          ! Just temporal variable  
103        REAL SLPMAX,S1MAX,S1NEW,S2NEW    REAL sbms, sfms, sfzs, sbmn, sfmn, sfzn
104  c    REAL ssum
105  C  Special pour poles    EXTERNAL ssum
106  c  
107        REAL sbms,sfms,sfzs,sbmn,sfmn,sfzn    REAL sqi, sqf
108        REAL sns0(ntra),snsz(ntra),snsm    LOGICAL limit
109        REAL qy1(iim,llm,ntra),qylat(iim,llm,ntra)  
110        REAL cx1(llm,ntra), cxLAT(llm,ntra)    lon = iim ! rem : Il est possible qu'un pbl. arrive ici
111        REAL cy1(llm,ntra), cyLAT(llm,ntra)    lat = jjp1 ! a cause des dim. differentes entre les
112        REAL z1(iim), zcos(iim), zsin(iim)    niv = llm !       tab. S et VGRI
113        REAL SSUM  
114        EXTERNAL SSUM    ! -----------------------------------------------------------------
115  C    ! initialisations
116        REAL sqi,sqf  
117        LOGICAL LIMIT    sbms = 0.
118      sfms = 0.
119        lon = iim         ! rem : Il est possible qu'un pbl. arrive ici    sfzs = 0.
120        lat = jjp1        ! a cause des dim. differentes entre les    sbmn = 0.
121        niv = llm         !       tab. S et VGRI    sfmn = 0.
122                          sfzn = 0.
123  c-----------------------------------------------------------------  
124  C initialisations    ! -----------------------------------------------------------------
125      ! *** Test : diag de la qtite totale de traceur dans
126        sbms = 0.    ! l'atmosphere avant l'advection en Y
127        sfms = 0.  
128        sfzs = 0.    sqi = 0.
129        sbmn = 0.    sqf = 0.
130        sfmn = 0.  
131        sfzn = 0.    DO l = 1, llm
132        DO j = 1, jjp1
133  c-----------------------------------------------------------------        DO i = 1, iim
134  C *** Test : diag de la qtite totale de traceur dans          sqi = sqi + s0(i, j, l, ntra)
135  C            l'atmosphere avant l'advection en Y        END DO
136  c      END DO
137        sqi = 0.    END DO
138        sqf = 0.    PRINT *, '---------- DIAG DANS ADVY - ENTREE --------'
139      PRINT *, 'sqi=', sqi
140        DO l = 1,llm  
141           DO j = 1,jjp1    ! -----------------------------------------------------------------
142             DO i = 1,iim    ! Interface : adaptation nouveau modele
143                sqi = sqi + S0(i,j,l,ntra)    ! -------------------------------------
144             END DO  
145           END DO    ! Conversion des flux de masses en kg
146        END DO    ! -AA 20/10/94  le signe -1 est necessaire car indexation opposee
147        PRINT*,'---------- DIAG DANS ADVY - ENTREE --------'  
148        PRINT*,'sqi=',sqi    DO l = 1, llm
149        DO j = 1, jjm
150  c-----------------------------------------------------------------        DO i = 1, iip1
151  C  Interface : adaptation nouveau modele          vgri(i, j, llm+1-l) = -1.*pbarv(i, j, l)
152  C  -------------------------------------        END DO
153  C      END DO
154  C  Conversion des flux de masses en kg    END DO
155  C-AA 20/10/94  le signe -1 est necessaire car indexation opposee  
156      ! AA Initialisation de flux fictifs aux bords sup. des boites pol.
157        DO 500 l = 1,llm  
158           DO 500 j = 1,jjm    DO l = 1, llm
159              DO 500 i = 1,iip1        DO i = 1, iip1
160              vgri (i,j,llm+1-l)=-1.*pbarv (i,j,l)        vgri(i, 0, l) = 0.
161    500 CONTINUE        vgri(i, jjp1, l) = 0.
162        END DO
163  CAA Initialisation de flux fictifs aux bords sup. des boites pol.    END DO
164    
165        DO l = 1,llm    ! ----------------- START HERE -----------------------
166           DO i = 1,iip1      ! boucle sur les niveaux
167               vgri(i,0,l) = 0.  
168               vgri(i,jjp1,l) = 0.    DO l = 1, niv
169           ENDDO  
170        ENDDO      ! place limits on appropriate moments before transport
171  c      ! (if flux-limiting is to be applied)
172  c----------------- START HERE -----------------------  
173  C  boucle sur les niveaux      IF (.NOT. limit) GO TO 11
174  C  
175        DO 1 L=1,NIV      DO jv = 1, ntra
176  C        DO k = 1, lat
177  C  place limits on appropriate moments before transport          DO i = 1, lon
178  C      (if flux-limiting is to be applied)            IF (s0(i,k,l,jv)>0.) THEN
179  C              slpmax = amax1(s0(i,k,l,jv), 0.)
180        IF(.NOT.LIMIT) GO TO 11              s1max = 1.5*slpmax
181  C              s1new = amin1(s1max, amax1(-s1max,sy(i,k,l,jv)))
182        DO 10 JV=1,NTRA              s2new = amin1(2.*slpmax-abs(s1new)/3., amax1(abs( &
183        DO 10 K=1,LAT                s1new)-slpmax,syy(i,k,l,jv)))
184        DO 100 I=1,LON              sy(i, k, l, jv) = s1new
185           IF(S0(I,K,L,JV).GT.0.) THEN              syy(i, k, l, jv) = s2new
186             SLPMAX=AMAX1(S0(I,K,L,JV),0.)              ssxy(i, k, l, jv) = amin1(slpmax, amax1(-slpmax,ssxy(i,k,l,jv)))
187             S1MAX=1.5*SLPMAX              syz(i, k, l, jv) = amin1(slpmax, amax1(-slpmax,syz(i,k,l,jv)))
188             S1NEW=AMIN1(S1MAX,AMAX1(-S1MAX,SY(I,K,L,JV)))            ELSE
189             S2NEW=AMIN1( 2.*SLPMAX-ABS(S1NEW)/3. ,              sy(i, k, l, jv) = 0.
190       +                  AMAX1(ABS(S1NEW)-SLPMAX,SYY(I,K,L,JV)) )              syy(i, k, l, jv) = 0.
191             SY (I,K,L,JV)=S1NEW              ssxy(i, k, l, jv) = 0.
192             SYY(I,K,L,JV)=S2NEW              syz(i, k, l, jv) = 0.
193         SSXY(I,K,L,JV)=AMIN1(SLPMAX,AMAX1(-SLPMAX,SSXY(I,K,L,JV)))            END IF
194         SYZ(I,K,L,JV)=AMIN1(SLPMAX,AMAX1(-SLPMAX,SYZ(I,K,L,JV)))          END DO
195           ELSE        END DO
196             SY (I,K,L,JV)=0.      END DO
197             SYY(I,K,L,JV)=0.  
198             SSXY(I,K,L,JV)=0.  11  CONTINUE
199             SYZ(I,K,L,JV)=0.  
200           ENDIF      ! le flux a travers le pole Nord est traite separement
201   100  CONTINUE  
202   10   CONTINUE      sm0 = 0.
203  C      DO jv = 1, ntra
204   11   CONTINUE        s00(jv) = 0.
205  C      END DO
206  C  le flux a travers le pole Nord est traite separement  
207  C      DO i = 1, lon
208        SM0=0.  
209        DO 20 JV=1,NTRA        IF (vgri(i,0,l)<=0.) THEN
210           S00(JV)=0.          fm(i, 0) = -vgri(i, 0, l)*dty
211   20   CONTINUE          alf(i, 0) = fm(i, 0)/sm(i, 1, l)
212  C          sm(i, 1, l) = sm(i, 1, l) - fm(i, 0)
213        DO 21 I=1,LON          sm0 = sm0 + fm(i, 0)
214  C        END IF
215           IF(VGRI(I,0,L).LE.0.) THEN  
216             FM(I,0)=-VGRI(I,0,L)*DTY        alfq(i, 0) = alf(i, 0)*alf(i, 0)
217             ALF(I,0)=FM(I,0)/SM(I,1,L)        alf1(i, 0) = 1. - alf(i, 0)
218             SM(I,1,L)=SM(I,1,L)-FM(I,0)        alf1q(i, 0) = alf1(i, 0)*alf1(i, 0)
219             SM0=SM0+FM(I,0)        alf2(i, 0) = alf1(i, 0) - alf(i, 0)
220           ENDIF        alf3(i, 0) = alf(i, 0)*alfq(i, 0)
221  C        alf4(i, 0) = alf1(i, 0)*alf1q(i, 0)
222           ALFQ(I,0)=ALF(I,0)*ALF(I,0)  
223           ALF1(I,0)=1.-ALF(I,0)      END DO
224           ALF1Q(I,0)=ALF1(I,0)*ALF1(I,0)      ! print*,'ADVYP 21'
225           ALF2(I,0)=ALF1(I,0)-ALF(I,0)  
226           ALF3(I,0)=ALF(I,0)*ALFQ(I,0)      DO jv = 1, ntra
227           ALF4(I,0)=ALF1(I,0)*ALF1Q(I,0)        DO i = 1, lon
228  C  
229   21   CONTINUE          IF (vgri(i,0,l)<=0.) THEN
230  c     print*,'ADVYP 21'  
231  C            f0(i, 0, jv) = alf(i, 0)*(s0(i,1,l,jv)-alf1(i,0)*(sy(i,1,l, &
232        DO 22 JV=1,NTRA              jv)-alf2(i,0)*syy(i,1,l,jv)))
233        DO 220 I=1,LON  
234  C            s00(jv) = s00(jv) + f0(i, 0, jv)
235           IF(VGRI(I,0,L).LE.0.) THEN            s0(i, 1, l, jv) = s0(i, 1, l, jv) - f0(i, 0, jv)
236  C            sy(i, 1, l, jv) = alf1q(i, 0)*(sy(i,1,l,jv)+3.*alf(i,0)*syy(i,1,l, &
237             F0(I,0,JV)=ALF(I,0)* ( S0(I,1,L,JV)-ALF1(I,0)*              jv))
238       +        ( SY(I,1,L,JV)-ALF2(I,0)*SYY(I,1,L,JV) ) )            syy(i, 1, l, jv) = alf4(i, 0)*syy(i, 1, l, jv)
239  C            ssx(i, 1, l, jv) = alf1(i, 0)*(ssx(i,1,l,jv)+alf(i,0)*ssxy(i,1,l,jv &
240             S00(JV)=S00(JV)+F0(I,0,JV)              ))
241             S0 (I,1,L,JV)=S0(I,1,L,JV)-F0(I,0,JV)            sz(i, 1, l, jv) = alf1(i, 0)*(sz(i,1,l,jv)+alf(i,0)*ssxz(i,1,l,jv))
242             SY (I,1,L,JV)=ALF1Q(I,0)*            ssxx(i, 1, l, jv) = alf1(i, 0)*ssxx(i, 1, l, jv)
243       +            (SY(I,1,L,JV)+3.*ALF(I,0)*SYY(I,1,L,JV))            ssxz(i, 1, l, jv) = alf1(i, 0)*ssxz(i, 1, l, jv)
244             SYY(I,1,L,JV)=ALF4 (I,0)*SYY(I,1,L,JV)            szz(i, 1, l, jv) = alf1(i, 0)*szz(i, 1, l, jv)
245             SSX (I,1,L,JV)=ALF1 (I,0)*            ssxy(i, 1, l, jv) = alf1q(i, 0)*ssxy(i, 1, l, jv)
246       +            (SSX(I,1,L,JV)+ALF(I,0)*SSXY(I,1,L,JV) )            syz(i, 1, l, jv) = alf1q(i, 0)*syz(i, 1, l, jv)
247             SZ (I,1,L,JV)=ALF1 (I,0)*  
248       +            (SZ(I,1,L,JV)+ALF(I,0)*SSXZ(I,1,L,JV) )          END IF
249             SSXX(I,1,L,JV)=ALF1 (I,0)*SSXX(I,1,L,JV)  
250             SSXZ(I,1,L,JV)=ALF1 (I,0)*SSXZ(I,1,L,JV)        END DO
251             SZZ(I,1,L,JV)=ALF1 (I,0)*SZZ(I,1,L,JV)      END DO
252             SSXY(I,1,L,JV)=ALF1Q(I,0)*SSXY(I,1,L,JV)  
253             SYZ(I,1,L,JV)=ALF1Q(I,0)*SYZ(I,1,L,JV)      DO i = 1, lon
254  C        IF (vgri(i,0,l)>0.) THEN
255           ENDIF          fm(i, 0) = vgri(i, 0, l)*dty
256  C          alf(i, 0) = fm(i, 0)/sm0
257   220  CONTINUE        END IF
258   22   CONTINUE      END DO
259  C  
260        DO 23 I=1,LON      DO jv = 1, ntra
261           IF(VGRI(I,0,L).GT.0.) THEN        DO i = 1, lon
262             FM(I,0)=VGRI(I,0,L)*DTY          IF (vgri(i,0,l)>0.) THEN
263             ALF(I,0)=FM(I,0)/SM0            f0(i, 0, jv) = alf(i, 0)*s00(jv)
264           ENDIF          END IF
265   23   CONTINUE        END DO
266  C      END DO
267        DO 24 JV=1,NTRA  
268        DO 240 I=1,LON      ! puts the temporary moments Fi into appropriate neighboring boxes
269           IF(VGRI(I,0,L).GT.0.) THEN  
270             F0(I,0,JV)=ALF(I,0)*S00(JV)      ! print*,'av ADVYP 25'
271           ENDIF      DO i = 1, lon
272   240  CONTINUE  
273   24   CONTINUE        IF (vgri(i,0,l)>0.) THEN
274  C          sm(i, 1, l) = sm(i, 1, l) + fm(i, 0)
275  C  puts the temporary moments Fi into appropriate neighboring boxes          alf(i, 0) = fm(i, 0)/sm(i, 1, l)
276  C        END IF
277  c     print*,'av ADVYP 25'  
278        DO 25 I=1,LON        alfq(i, 0) = alf(i, 0)*alf(i, 0)
279  C        alf1(i, 0) = 1. - alf(i, 0)
280           IF(VGRI(I,0,L).GT.0.) THEN        alf1q(i, 0) = alf1(i, 0)*alf1(i, 0)
281             SM(I,1,L)=SM(I,1,L)+FM(I,0)        alf2(i, 0) = alf1(i, 0) - alf(i, 0)
282             ALF(I,0)=FM(I,0)/SM(I,1,L)        alf3(i, 0) = alf1(i, 0)*alf(i, 0)
283           ENDIF  
284  C      END DO
285           ALFQ(I,0)=ALF(I,0)*ALF(I,0)      ! print*,'av ADVYP 25'
286           ALF1(I,0)=1.-ALF(I,0)  
287           ALF1Q(I,0)=ALF1(I,0)*ALF1(I,0)      DO jv = 1, ntra
288           ALF2(I,0)=ALF1(I,0)-ALF(I,0)        DO i = 1, lon
289           ALF3(I,0)=ALF1(I,0)*ALF(I,0)  
290  C          IF (vgri(i,0,l)>0.) THEN
291   25   CONTINUE  
292  c     print*,'av ADVYP 25'            temptm = alf(i, 0)*s0(i, 1, l, jv) - alf1(i, 0)*f0(i, 0, jv)
293  C            s0(i, 1, l, jv) = s0(i, 1, l, jv) + f0(i, 0, jv)
294        DO 26 JV=1,NTRA            syy(i, 1, l, jv) = alf1q(i, 0)*syy(i, 1, l, jv) + &
295        DO 260 I=1,LON              5.*(alf3(i,0)*sy(i,1,l,jv)-alf2(i,0)*temptm)
296  C            sy(i, 1, l, jv) = alf1(i, 0)*sy(i, 1, l, jv) + 3.*temptm
297           IF(VGRI(I,0,L).GT.0.) THEN            ssxy(i, 1, l, jv) = alf1(i, 0)*ssxy(i, 1, l, jv) + &
298  C              3.*alf(i, 0)*ssx(i, 1, l, jv)
299           TEMPTM=ALF(I,0)*S0(I,1,L,JV)-ALF1(I,0)*F0(I,0,JV)            syz(i, 1, l, jv) = alf1(i, 0)*syz(i, 1, l, jv) + &
300           S0 (I,1,L,JV)=S0(I,1,L,JV)+F0(I,0,JV)              3.*alf(i, 0)*sz(i, 1, l, jv)
301           SYY(I,1,L,JV)=ALF1Q(I,0)*SYY(I,1,L,JV)  
302       +        +5.*( ALF3 (I,0)*SY (I,1,L,JV)-ALF2(I,0)*TEMPTM )          END IF
303           SY (I,1,L,JV)=ALF1 (I,0)*SY (I,1,L,JV)+3.*TEMPTM  
304        SSXY(I,1,L,JV)=ALF1 (I,0)*SSXY(I,1,L,JV)+3.*ALF(I,0)*SSX(I,1,L,JV)        END DO
305        SYZ(I,1,L,JV)=ALF1 (I,0)*SYZ(I,1,L,JV)+3.*ALF(I,0)*SZ(I,1,L,JV)      END DO
306  C  
307           ENDIF      ! calculate flux and moments between adjacent boxes
308  C      ! 1- create temporary moments/masses for partial boxes in transit
309   260  CONTINUE      ! 2- reajusts moments remaining in the box
310   26   CONTINUE  
311  C      ! flux from KP to K if V(K).lt.0 and from K to KP if V(K).gt.0
312  C  calculate flux and moments between adjacent boxes  
313  C  1- create temporary moments/masses for partial boxes in transit      ! print*,'av ADVYP 30'
314  C  2- reajusts moments remaining in the box      DO k = 1, lat - 1
315  C        kp = k + 1
316  C  flux from KP to K if V(K).lt.0 and from K to KP if V(K).gt.0        DO i = 1, lon
317  C  
318  c     print*,'av ADVYP 30'          IF (vgri(i,k,l)<0.) THEN
319        DO 30 K=1,LAT-1            fm(i, k) = -vgri(i, k, l)*dty
320        KP=K+1            alf(i, k) = fm(i, k)/sm(i, kp, l)
321        DO 300 I=1,LON            sm(i, kp, l) = sm(i, kp, l) - fm(i, k)
322  C          ELSE
323           IF(VGRI(I,K,L).LT.0.) THEN            fm(i, k) = vgri(i, k, l)*dty
324             FM(I,K)=-VGRI(I,K,L)*DTY            alf(i, k) = fm(i, k)/sm(i, k, l)
325             ALF(I,K)=FM(I,K)/SM(I,KP,L)            sm(i, k, l) = sm(i, k, l) - fm(i, k)
326             SM(I,KP,L)=SM(I,KP,L)-FM(I,K)          END IF
327           ELSE  
328             FM(I,K)=VGRI(I,K,L)*DTY          alfq(i, k) = alf(i, k)*alf(i, k)
329             ALF(I,K)=FM(I,K)/SM(I,K,L)          alf1(i, k) = 1. - alf(i, k)
330             SM(I,K,L)=SM(I,K,L)-FM(I,K)          alf1q(i, k) = alf1(i, k)*alf1(i, k)
331           ENDIF          alf2(i, k) = alf1(i, k) - alf(i, k)
332  C          alf3(i, k) = alf(i, k)*alfq(i, k)
333           ALFQ(I,K)=ALF(I,K)*ALF(I,K)          alf4(i, k) = alf1(i, k)*alf1q(i, k)
334           ALF1(I,K)=1.-ALF(I,K)  
335           ALF1Q(I,K)=ALF1(I,K)*ALF1(I,K)        END DO
336           ALF2(I,K)=ALF1(I,K)-ALF(I,K)      END DO
337           ALF3(I,K)=ALF(I,K)*ALFQ(I,K)      ! print*,'ap ADVYP 30'
338           ALF4(I,K)=ALF1(I,K)*ALF1Q(I,K)  
339  C      DO jv = 1, ntra
340   300  CONTINUE        DO k = 1, lat - 1
341   30   CONTINUE          kp = k + 1
342  c     print*,'ap ADVYP 30'          DO i = 1, lon
343  C  
344        DO 31 JV=1,NTRA            IF (vgri(i,k,l)<0.) THEN
345        DO 31 K=1,LAT-1  
346        KP=K+1              f0(i, k, jv) = alf(i, k)*(s0(i,kp,l,jv)-alf1(i,k)*(sy(i,kp,l, &
347        DO 310 I=1,LON                jv)-alf2(i,k)*syy(i,kp,l,jv)))
348  C              fy(i, k, jv) = alfq(i, k)*(sy(i,kp,l,jv)-3.*alf1(i,k)*syy(i,kp,l, &
349           IF(VGRI(I,K,L).LT.0.) THEN                jv))
350  C              fyy(i, k, jv) = alf3(i, k)*syy(i, kp, l, jv)
351             F0 (I,K,JV)=ALF (I,K)* ( S0(I,KP,L,JV)-ALF1(I,K)*              fx(i, k, jv) = alf(i, k)*(ssx(i,kp,l,jv)-alf1(i,k)*ssxy(i,kp,l,jv &
352       +        ( SY(I,KP,L,JV)-ALF2(I,K)*SYY(I,KP,L,JV) ) )                ))
353             FY (I,K,JV)=ALFQ(I,K)*              fz(i, k, jv) = alf(i, k)*(sz(i,kp,l,jv)-alf1(i,k)*syz(i,kp,l,jv))
354       +                 (SY(I,KP,L,JV)-3.*ALF1(I,K)*SYY(I,KP,L,JV))              fxy(i, k, jv) = alfq(i, k)*ssxy(i, kp, l, jv)
355             FYY(I,K,JV)=ALF3(I,K)*SYY(I,KP,L,JV)              fyz(i, k, jv) = alfq(i, k)*syz(i, kp, l, jv)
356             FX (I,K,JV)=ALF (I,K)*              fxx(i, k, jv) = alf(i, k)*ssxx(i, kp, l, jv)
357       +                 (SSX(I,KP,L,JV)-ALF1(I,K)*SSXY(I,KP,L,JV))              fxz(i, k, jv) = alf(i, k)*ssxz(i, kp, l, jv)
358             FZ (I,K,JV)=ALF (I,K)*              fzz(i, k, jv) = alf(i, k)*szz(i, kp, l, jv)
359       +                 (SZ(I,KP,L,JV)-ALF1(I,K)*SYZ(I,KP,L,JV))  
360             FXY(I,K,JV)=ALFQ(I,K)*SSXY(I,KP,L,JV)              s0(i, kp, l, jv) = s0(i, kp, l, jv) - f0(i, k, jv)
361             FYZ(I,K,JV)=ALFQ(I,K)*SYZ(I,KP,L,JV)              sy(i, kp, l, jv) = alf1q(i, k)*(sy(i,kp,l,jv)+3.*alf(i,k)*syy(i, &
362             FXX(I,K,JV)=ALF (I,K)*SSXX(I,KP,L,JV)                kp,l,jv))
363             FXZ(I,K,JV)=ALF (I,K)*SSXZ(I,KP,L,JV)              syy(i, kp, l, jv) = alf4(i, k)*syy(i, kp, l, jv)
364             FZZ(I,K,JV)=ALF (I,K)*SZZ(I,KP,L,JV)              ssx(i, kp, l, jv) = ssx(i, kp, l, jv) - fx(i, k, jv)
365  C              sz(i, kp, l, jv) = sz(i, kp, l, jv) - fz(i, k, jv)
366             S0 (I,KP,L,JV)=S0(I,KP,L,JV)-F0(I,K,JV)              ssxx(i, kp, l, jv) = ssxx(i, kp, l, jv) - fxx(i, k, jv)
367             SY (I,KP,L,JV)=ALF1Q(I,K)*              ssxz(i, kp, l, jv) = ssxz(i, kp, l, jv) - fxz(i, k, jv)
368       +                 (SY(I,KP,L,JV)+3.*ALF(I,K)*SYY(I,KP,L,JV))              szz(i, kp, l, jv) = szz(i, kp, l, jv) - fzz(i, k, jv)
369             SYY(I,KP,L,JV)=ALF4(I,K)*SYY(I,KP,L,JV)              ssxy(i, kp, l, jv) = alf1q(i, k)*ssxy(i, kp, l, jv)
370             SSX (I,KP,L,JV)=SSX (I,KP,L,JV)-FX (I,K,JV)              syz(i, kp, l, jv) = alf1q(i, k)*syz(i, kp, l, jv)
371             SZ (I,KP,L,JV)=SZ (I,KP,L,JV)-FZ (I,K,JV)  
372             SSXX(I,KP,L,JV)=SSXX(I,KP,L,JV)-FXX(I,K,JV)            ELSE
373             SSXZ(I,KP,L,JV)=SSXZ(I,KP,L,JV)-FXZ(I,K,JV)  
374             SZZ(I,KP,L,JV)=SZZ(I,KP,L,JV)-FZZ(I,K,JV)              f0(i, k, jv) = alf(i, k)*(s0(i,k,l,jv)+alf1(i,k)*(sy(i,k,l, &
375             SSXY(I,KP,L,JV)=ALF1Q(I,K)*SSXY(I,KP,L,JV)                jv)+alf2(i,k)*syy(i,k,l,jv)))
376             SYZ(I,KP,L,JV)=ALF1Q(I,K)*SYZ(I,KP,L,JV)              fy(i, k, jv) = alfq(i, k)*(sy(i,k,l,jv)+3.*alf1(i,k)*syy(i,k,l,jv &
377  C                ))
378           ELSE              fyy(i, k, jv) = alf3(i, k)*syy(i, k, l, jv)
379  C              fx(i, k, jv) = alf(i, k)*(ssx(i,k,l,jv)+alf1(i,k)*ssxy(i,k,l,jv))
380             F0 (I,K,JV)=ALF (I,K)* ( S0(I,K,L,JV)+ALF1(I,K)*              fz(i, k, jv) = alf(i, k)*(sz(i,k,l,jv)+alf1(i,k)*syz(i,k,l,jv))
381       +             ( SY(I,K,L,JV)+ALF2(I,K)*SYY(I,K,L,JV) ) )              fxy(i, k, jv) = alfq(i, k)*ssxy(i, k, l, jv)
382             FY (I,K,JV)=ALFQ(I,K)*              fyz(i, k, jv) = alfq(i, k)*syz(i, k, l, jv)
383       +                 (SY(I,K,L,JV)+3.*ALF1(I,K)*SYY(I,K,L,JV))              fxx(i, k, jv) = alf(i, k)*ssxx(i, k, l, jv)
384             FYY(I,K,JV)=ALF3(I,K)*SYY(I,K,L,JV)              fxz(i, k, jv) = alf(i, k)*ssxz(i, k, l, jv)
385        FX (I,K,JV)=ALF (I,K)*(SSX(I,K,L,JV)+ALF1(I,K)*SSXY(I,K,L,JV))              fzz(i, k, jv) = alf(i, k)*szz(i, k, l, jv)
386        FZ (I,K,JV)=ALF (I,K)*(SZ(I,K,L,JV)+ALF1(I,K)*SYZ(I,K,L,JV))  
387             FXY(I,K,JV)=ALFQ(I,K)*SSXY(I,K,L,JV)              s0(i, k, l, jv) = s0(i, k, l, jv) - f0(i, k, jv)
388             FYZ(I,K,JV)=ALFQ(I,K)*SYZ(I,K,L,JV)              sy(i, k, l, jv) = alf1q(i, k)*(sy(i,k,l,jv)-3.*alf(i,k)*syy(i,k,l &
389             FXX(I,K,JV)=ALF (I,K)*SSXX(I,K,L,JV)                ,jv))
390             FXZ(I,K,JV)=ALF (I,K)*SSXZ(I,K,L,JV)              syy(i, k, l, jv) = alf4(i, k)*syy(i, k, l, jv)
391             FZZ(I,K,JV)=ALF (I,K)*SZZ(I,K,L,JV)              ssx(i, k, l, jv) = ssx(i, k, l, jv) - fx(i, k, jv)
392  C              sz(i, k, l, jv) = sz(i, k, l, jv) - fz(i, k, jv)
393             S0 (I,K,L,JV)=S0 (I,K,L,JV)-F0 (I,K,JV)              ssxx(i, k, l, jv) = ssxx(i, k, l, jv) - fxx(i, k, jv)
394             SY (I,K,L,JV)=ALF1Q(I,K)*              ssxz(i, k, l, jv) = ssxz(i, k, l, jv) - fxz(i, k, jv)
395       +                  (SY(I,K,L,JV)-3.*ALF(I,K)*SYY(I,K,L,JV))              szz(i, k, l, jv) = szz(i, k, l, jv) - fzz(i, k, jv)
396             SYY(I,K,L,JV)=ALF4(I,K)*SYY(I,K,L,JV)              ssxy(i, k, l, jv) = alf1q(i, k)*ssxy(i, k, l, jv)
397             SSX (I,K,L,JV)=SSX (I,K,L,JV)-FX (I,K,JV)              syz(i, k, l, jv) = alf1q(i, k)*syz(i, k, l, jv)
398             SZ (I,K,L,JV)=SZ (I,K,L,JV)-FZ (I,K,JV)  
399             SSXX(I,K,L,JV)=SSXX(I,K,L,JV)-FXX(I,K,JV)            END IF
400             SSXZ(I,K,L,JV)=SSXZ(I,K,L,JV)-FXZ(I,K,JV)  
401             SZZ(I,K,L,JV)=SZZ(I,K,L,JV)-FZZ(I,K,JV)          END DO
402             SSXY(I,K,L,JV)=ALF1Q(I,K)*SSXY(I,K,L,JV)        END DO
403             SYZ(I,K,L,JV)=ALF1Q(I,K)*SYZ(I,K,L,JV)      END DO
404  C      ! print*,'ap ADVYP 31'
405           ENDIF  
406  C      ! puts the temporary moments Fi into appropriate neighboring boxes
407   310  CONTINUE  
408   31   CONTINUE      DO k = 1, lat - 1
409  c     print*,'ap ADVYP 31'        kp = k + 1
410  C        DO i = 1, lon
411  C  puts the temporary moments Fi into appropriate neighboring boxes  
412  C          IF (vgri(i,k,l)<0.) THEN
413        DO 32 K=1,LAT-1            sm(i, k, l) = sm(i, k, l) + fm(i, k)
414        KP=K+1            alf(i, k) = fm(i, k)/sm(i, k, l)
415        DO 320 I=1,LON          ELSE
416  C            sm(i, kp, l) = sm(i, kp, l) + fm(i, k)
417           IF(VGRI(I,K,L).LT.0.) THEN            alf(i, k) = fm(i, k)/sm(i, kp, l)
418             SM(I,K,L)=SM(I,K,L)+FM(I,K)          END IF
419             ALF(I,K)=FM(I,K)/SM(I,K,L)  
420           ELSE          alfq(i, k) = alf(i, k)*alf(i, k)
421             SM(I,KP,L)=SM(I,KP,L)+FM(I,K)          alf1(i, k) = 1. - alf(i, k)
422             ALF(I,K)=FM(I,K)/SM(I,KP,L)          alf1q(i, k) = alf1(i, k)*alf1(i, k)
423           ENDIF          alf2(i, k) = alf1(i, k) - alf(i, k)
424  C          alf3(i, k) = alf1(i, k)*alf(i, k)
425           ALFQ(I,K)=ALF(I,K)*ALF(I,K)  
426           ALF1(I,K)=1.-ALF(I,K)        END DO
427           ALF1Q(I,K)=ALF1(I,K)*ALF1(I,K)      END DO
428           ALF2(I,K)=ALF1(I,K)-ALF(I,K)      ! print*,'ap ADVYP 32'
429           ALF3(I,K)=ALF1(I,K)*ALF(I,K)  
430  C      DO jv = 1, ntra
431   320  CONTINUE        DO k = 1, lat - 1
432   32   CONTINUE          kp = k + 1
433  c     print*,'ap ADVYP 32'          DO i = 1, lon
434  C  
435        DO 33 JV=1,NTRA            IF (vgri(i,k,l)<0.) THEN
436        DO 33 K=1,LAT-1  
437        KP=K+1              temptm = -alf(i, k)*s0(i, k, l, jv) + alf1(i, k)*f0(i, k, jv)
438        DO 330 I=1,LON              s0(i, k, l, jv) = s0(i, k, l, jv) + f0(i, k, jv)
439  C              syy(i, k, l, jv) = alfq(i, k)*fyy(i, k, jv) + &
440           IF(VGRI(I,K,L).LT.0.) THEN                alf1q(i, k)*syy(i, k, l, jv) + 5.*(alf3(i,k)*(fy(i,k,jv)-sy(i, &
441  C                k,l,jv))+alf2(i,k)*temptm)
442           TEMPTM=-ALF(I,K)*S0(I,K,L,JV)+ALF1(I,K)*F0(I,K,JV)              sy(i, k, l, jv) = alf(i, k)*fy(i, k, jv) + &
443           S0 (I,K,L,JV)=S0(I,K,L,JV)+F0(I,K,JV)                alf1(i, k)*sy(i, k, l, jv) + 3.*temptm
444         SYY(I,K,L,JV)=ALFQ(I,K)*FYY(I,K,JV)+ALF1Q(I,K)*SYY(I,K,L,JV)              ssxy(i, k, l, jv) = alf(i, k)*fxy(i, k, jv) + &
445       +  +5.*( ALF3(I,K)*(FY(I,K,JV)-SY(I,K,L,JV))+ALF2(I,K)*TEMPTM )                alf1(i, k)*ssxy(i, k, l, jv) + 3.*(alf1(i,k)*fx(i,k,jv)-alf(i,k &
446           SY (I,K,L,JV)=ALF(I,K)*FY(I,K,JV)+ALF1(I,K)*SY(I,K,L,JV)                )*ssx(i,k,l,jv))
447       +            +3.*TEMPTM              syz(i, k, l, jv) = alf(i, k)*fyz(i, k, jv) + &
448         SSXY(I,K,L,JV)=ALF (I,K)*FXY(I,K,JV)+ALF1(I,K)*SSXY(I,K,L,JV)                alf1(i, k)*syz(i, k, l, jv) + 3.*(alf1(i,k)*fz(i,k,jv)-alf(i,k) &
449       +         +3.*(ALF1(I,K)*FX (I,K,JV)-ALF (I,K)*SSX (I,K,L,JV))                *sz(i,k,l,jv))
450         SYZ(I,K,L,JV)=ALF (I,K)*FYZ(I,K,JV)+ALF1(I,K)*SYZ(I,K,L,JV)              ssx(i, k, l, jv) = ssx(i, k, l, jv) + fx(i, k, jv)
451       +         +3.*(ALF1(I,K)*FZ (I,K,JV)-ALF (I,K)*SZ (I,K,L,JV))              sz(i, k, l, jv) = sz(i, k, l, jv) + fz(i, k, jv)
452           SSX (I,K,L,JV)=SSX (I,K,L,JV)+FX (I,K,JV)              ssxx(i, k, l, jv) = ssxx(i, k, l, jv) + fxx(i, k, jv)
453           SZ (I,K,L,JV)=SZ (I,K,L,JV)+FZ (I,K,JV)              ssxz(i, k, l, jv) = ssxz(i, k, l, jv) + fxz(i, k, jv)
454           SSXX(I,K,L,JV)=SSXX(I,K,L,JV)+FXX(I,K,JV)              szz(i, k, l, jv) = szz(i, k, l, jv) + fzz(i, k, jv)
455           SSXZ(I,K,L,JV)=SSXZ(I,K,L,JV)+FXZ(I,K,JV)  
456           SZZ(I,K,L,JV)=SZZ(I,K,L,JV)+FZZ(I,K,JV)            ELSE
457  C  
458           ELSE              temptm = alf(i, k)*s0(i, kp, l, jv) - alf1(i, k)*f0(i, k, jv)
459  C              s0(i, kp, l, jv) = s0(i, kp, l, jv) + f0(i, k, jv)
460           TEMPTM=ALF(I,K)*S0(I,KP,L,JV)-ALF1(I,K)*F0(I,K,JV)              syy(i, kp, l, jv) = alfq(i, k)*fyy(i, k, jv) + &
461           S0 (I,KP,L,JV)=S0(I,KP,L,JV)+F0(I,K,JV)                alf1q(i, k)*syy(i, kp, l, jv) + 5.*(alf3(i,k)*(sy(i,kp,l, &
462         SYY(I,KP,L,JV)=ALFQ(I,K)*FYY(I,K,JV)+ALF1Q(I,K)*SYY(I,KP,L,JV)                jv)-fy(i,k,jv))-alf2(i,k)*temptm)
463       +  +5.*( ALF3(I,K)*(SY(I,KP,L,JV)-FY(I,K,JV))-ALF2(I,K)*TEMPTM )              sy(i, kp, l, jv) = alf(i, k)*fy(i, k, jv) + &
464           SY (I,KP,L,JV)=ALF(I,K)*FY(I,K,JV)+ALF1(I,K)*SY(I,KP,L,JV)                alf1(i, k)*sy(i, kp, l, jv) + 3.*temptm
465       +                 +3.*TEMPTM              ssxy(i, kp, l, jv) = alf(i, k)*fxy(i, k, jv) + &
466         SSXY(I,KP,L,JV)=ALF(I,K)*FXY(I,K,JV)+ALF1(I,K)*SSXY(I,KP,L,JV)                alf1(i, k)*ssxy(i, kp, l, jv) + 3.*(alf(i,k)*ssx(i,kp,l,jv)- &
467       +             +3.*(ALF(I,K)*SSX(I,KP,L,JV)-ALF1(I,K)*FX(I,K,JV))                alf1(i,k)*fx(i,k,jv))
468           SYZ(I,KP,L,JV)=ALF(I,K)*FYZ(I,K,JV)+ALF1(I,K)*SYZ(I,KP,L,JV)              syz(i, kp, l, jv) = alf(i, k)*fyz(i, k, jv) + &
469       +             +3.*(ALF(I,K)*SZ(I,KP,L,JV)-ALF1(I,K)*FZ(I,K,JV))                alf1(i, k)*syz(i, kp, l, jv) + 3.*(alf(i,k)*sz(i,kp,l,jv)-alf1( &
470           SSX (I,KP,L,JV)=SSX (I,KP,L,JV)+FX (I,K,JV)                i,k)*fz(i,k,jv))
471           SZ (I,KP,L,JV)=SZ (I,KP,L,JV)+FZ (I,K,JV)              ssx(i, kp, l, jv) = ssx(i, kp, l, jv) + fx(i, k, jv)
472           SSXX(I,KP,L,JV)=SSXX(I,KP,L,JV)+FXX(I,K,JV)              sz(i, kp, l, jv) = sz(i, kp, l, jv) + fz(i, k, jv)
473           SSXZ(I,KP,L,JV)=SSXZ(I,KP,L,JV)+FXZ(I,K,JV)              ssxx(i, kp, l, jv) = ssxx(i, kp, l, jv) + fxx(i, k, jv)
474           SZZ(I,KP,L,JV)=SZZ(I,KP,L,JV)+FZZ(I,K,JV)              ssxz(i, kp, l, jv) = ssxz(i, kp, l, jv) + fxz(i, k, jv)
475  C              szz(i, kp, l, jv) = szz(i, kp, l, jv) + fzz(i, k, jv)
476           ENDIF  
477  C            END IF
478   330  CONTINUE  
479   33   CONTINUE          END DO
480  c     print*,'ap ADVYP 33'        END DO
481  C      END DO
482  C  traitement special pour le pole Sud (idem pole Nord)      ! print*,'ap ADVYP 33'
483  C  
484        K=LAT      ! traitement special pour le pole Sud (idem pole Nord)
485  C  
486        SM0=0.      k = lat
487        DO 40 JV=1,NTRA  
488           S00(JV)=0.      sm0 = 0.
489   40   CONTINUE      DO jv = 1, ntra
490  C        s00(jv) = 0.
491        DO 41 I=1,LON      END DO
492  C  
493           IF(VGRI(I,K,L).GE.0.) THEN      DO i = 1, lon
494             FM(I,K)=VGRI(I,K,L)*DTY  
495             ALF(I,K)=FM(I,K)/SM(I,K,L)        IF (vgri(i,k,l)>=0.) THEN
496             SM(I,K,L)=SM(I,K,L)-FM(I,K)          fm(i, k) = vgri(i, k, l)*dty
497             SM0=SM0+FM(I,K)          alf(i, k) = fm(i, k)/sm(i, k, l)
498           ENDIF          sm(i, k, l) = sm(i, k, l) - fm(i, k)
499  C          sm0 = sm0 + fm(i, k)
500           ALFQ(I,K)=ALF(I,K)*ALF(I,K)        END IF
501           ALF1(I,K)=1.-ALF(I,K)  
502           ALF1Q(I,K)=ALF1(I,K)*ALF1(I,K)        alfq(i, k) = alf(i, k)*alf(i, k)
503           ALF2(I,K)=ALF1(I,K)-ALF(I,K)        alf1(i, k) = 1. - alf(i, k)
504           ALF3(I,K)=ALF(I,K)*ALFQ(I,K)        alf1q(i, k) = alf1(i, k)*alf1(i, k)
505           ALF4(I,K)=ALF1(I,K)*ALF1Q(I,K)        alf2(i, k) = alf1(i, k) - alf(i, k)
506  C        alf3(i, k) = alf(i, k)*alfq(i, k)
507   41   CONTINUE        alf4(i, k) = alf1(i, k)*alf1q(i, k)
508  c     print*,'ap ADVYP 41'  
509  C      END DO
510        DO 42 JV=1,NTRA      ! print*,'ap ADVYP 41'
511        DO 420 I=1,LON  
512  C      DO jv = 1, ntra
513           IF(VGRI(I,K,L).GE.0.) THEN        DO i = 1, lon
514             F0 (I,K,JV)=ALF(I,K)* ( S0(I,K,L,JV)+ALF1(I,K)*  
515       +             ( SY(I,K,L,JV)+ALF2(I,K)*SYY(I,K,L,JV) ) )          IF (vgri(i,k,l)>=0.) THEN
516             S00(JV)=S00(JV)+F0(I,K,JV)            f0(i, k, jv) = alf(i, k)*(s0(i,k,l,jv)+alf1(i,k)*(sy(i,k,l, &
517  C              jv)+alf2(i,k)*syy(i,k,l,jv)))
518             S0 (I,K,L,JV)=S0 (I,K,L,JV)-F0 (I,K,JV)            s00(jv) = s00(jv) + f0(i, k, jv)
519             SY (I,K,L,JV)=ALF1Q(I,K)*  
520       +                  (SY(I,K,L,JV)-3.*ALF(I,K)*SYY(I,K,L,JV))            s0(i, k, l, jv) = s0(i, k, l, jv) - f0(i, k, jv)
521             SYY(I,K,L,JV)=ALF4 (I,K)*SYY(I,K,L,JV)            sy(i, k, l, jv) = alf1q(i, k)*(sy(i,k,l,jv)-3.*alf(i,k)*syy(i,k,l, &
522        SSX (I,K,L,JV)=ALF1(I,K)*(SSX(I,K,L,JV)-ALF(I,K)*SSXY(I,K,L,JV))              jv))
523        SZ (I,K,L,JV)=ALF1(I,K)*(SZ(I,K,L,JV)-ALF(I,K)*SYZ(I,K,L,JV))            syy(i, k, l, jv) = alf4(i, k)*syy(i, k, l, jv)
524             SSXX(I,K,L,JV)=ALF1 (I,K)*SSXX(I,K,L,JV)            ssx(i, k, l, jv) = alf1(i, k)*(ssx(i,k,l,jv)-alf(i,k)*ssxy(i,k,l,jv &
525             SSXZ(I,K,L,JV)=ALF1 (I,K)*SSXZ(I,K,L,JV)              ))
526             SZZ(I,K,L,JV)=ALF1 (I,K)*SZZ(I,K,L,JV)            sz(i, k, l, jv) = alf1(i, k)*(sz(i,k,l,jv)-alf(i,k)*syz(i,k,l,jv))
527             SSXY(I,K,L,JV)=ALF1Q(I,K)*SSXY(I,K,L,JV)            ssxx(i, k, l, jv) = alf1(i, k)*ssxx(i, k, l, jv)
528             SYZ(I,K,L,JV)=ALF1Q(I,K)*SYZ(I,K,L,JV)            ssxz(i, k, l, jv) = alf1(i, k)*ssxz(i, k, l, jv)
529           ENDIF            szz(i, k, l, jv) = alf1(i, k)*szz(i, k, l, jv)
530  C            ssxy(i, k, l, jv) = alf1q(i, k)*ssxy(i, k, l, jv)
531   420  CONTINUE            syz(i, k, l, jv) = alf1q(i, k)*syz(i, k, l, jv)
532   42   CONTINUE          END IF
533  c     print*,'ap ADVYP 42'  
534  C        END DO
535        DO 43 I=1,LON      END DO
536           IF(VGRI(I,K,L).LT.0.) THEN      ! print*,'ap ADVYP 42'
537             FM(I,K)=-VGRI(I,K,L)*DTY  
538             ALF(I,K)=FM(I,K)/SM0      DO i = 1, lon
539           ENDIF        IF (vgri(i,k,l)<0.) THEN
540   43   CONTINUE          fm(i, k) = -vgri(i, k, l)*dty
541  c     print*,'ap ADVYP 43'          alf(i, k) = fm(i, k)/sm0
542  C        END IF
543        DO 44 JV=1,NTRA      END DO
544        DO 440 I=1,LON      ! print*,'ap ADVYP 43'
545           IF(VGRI(I,K,L).LT.0.) THEN  
546             F0(I,K,JV)=ALF(I,K)*S00(JV)      DO jv = 1, ntra
547           ENDIF        DO i = 1, lon
548   440  CONTINUE          IF (vgri(i,k,l)<0.) THEN
549   44   CONTINUE            f0(i, k, jv) = alf(i, k)*s00(jv)
550  C          END IF
551  C  puts the temporary moments Fi into appropriate neighboring boxes        END DO
552  C      END DO
553        DO 45 I=1,LON  
554  C      ! puts the temporary moments Fi into appropriate neighboring boxes
555           IF(VGRI(I,K,L).LT.0.) THEN  
556             SM(I,K,L)=SM(I,K,L)+FM(I,K)      DO i = 1, lon
557             ALF(I,K)=FM(I,K)/SM(I,K,L)  
558           ENDIF        IF (vgri(i,k,l)<0.) THEN
559  C          sm(i, k, l) = sm(i, k, l) + fm(i, k)
560           ALFQ(I,K)=ALF(I,K)*ALF(I,K)          alf(i, k) = fm(i, k)/sm(i, k, l)
561           ALF1(I,K)=1.-ALF(I,K)        END IF
562           ALF1Q(I,K)=ALF1(I,K)*ALF1(I,K)  
563           ALF2(I,K)=ALF1(I,K)-ALF(I,K)        alfq(i, k) = alf(i, k)*alf(i, k)
564           ALF3(I,K)=ALF1(I,K)*ALF(I,K)        alf1(i, k) = 1. - alf(i, k)
565  C        alf1q(i, k) = alf1(i, k)*alf1(i, k)
566   45   CONTINUE        alf2(i, k) = alf1(i, k) - alf(i, k)
567  c     print*,'ap ADVYP 45'        alf3(i, k) = alf1(i, k)*alf(i, k)
568  C  
569        DO 46 JV=1,NTRA      END DO
570        DO 460 I=1,LON      ! print*,'ap ADVYP 45'
571  C  
572           IF(VGRI(I,K,L).LT.0.) THEN      DO jv = 1, ntra
573  C        DO i = 1, lon
574           TEMPTM=-ALF(I,K)*S0(I,K,L,JV)+ALF1(I,K)*F0(I,K,JV)  
575           S0 (I,K,L,JV)=S0(I,K,L,JV)+F0(I,K,JV)          IF (vgri(i,k,l)<0.) THEN
576           SYY(I,K,L,JV)=ALF1Q(I,K)*SYY(I,K,L,JV)  
577       +           +5.*(-ALF3 (I,K)*SY (I,K,L,JV)+ALF2(I,K)*TEMPTM )            temptm = -alf(i, k)*s0(i, k, l, jv) + alf1(i, k)*f0(i, k, jv)
578           SY (I,K,L,JV)=ALF1(I,K)*SY (I,K,L,JV)+3.*TEMPTM            s0(i, k, l, jv) = s0(i, k, l, jv) + f0(i, k, jv)
579        SSXY(I,K,L,JV)=ALF1(I,K)*SSXY(I,K,L,JV)-3.*ALF(I,K)*SSX(I,K,L,JV)            syy(i, k, l, jv) = alf1q(i, k)*syy(i, k, l, jv) + &
580        SYZ(I,K,L,JV)=ALF1(I,K)*SYZ(I,K,L,JV)-3.*ALF(I,K)*SZ(I,K,L,JV)              5.*(-alf3(i,k)*sy(i,k,l,jv)+alf2(i,k)*temptm)
581  C            sy(i, k, l, jv) = alf1(i, k)*sy(i, k, l, jv) + 3.*temptm
582           ENDIF            ssxy(i, k, l, jv) = alf1(i, k)*ssxy(i, k, l, jv) - &
583  C              3.*alf(i, k)*ssx(i, k, l, jv)
584   460  CONTINUE            syz(i, k, l, jv) = alf1(i, k)*syz(i, k, l, jv) - &
585   46   CONTINUE              3.*alf(i, k)*sz(i, k, l, jv)
586  c     print*,'ap ADVYP 46'  
587  C          END IF
588   1    CONTINUE  
589          END DO
590  c--------------------------------------------------      END DO
591  C     bouclage cyclique horizontal .      ! print*,'ap ADVYP 46'
592        
593        DO l = 1,llm    END DO
594           DO jv = 1,ntra  
595              DO j = 1,jjp1    ! --------------------------------------------------
596                 SM(iip1,j,l) = SM(1,j,l)    ! bouclage cyclique horizontal .
597                 S0(iip1,j,l,jv) = S0(1,j,l,jv)  
598                 SSX(iip1,j,l,jv) = SSX(1,j,l,jv)      DO l = 1, llm
599                 SY(iip1,j,l,jv) = SY(1,j,l,jv)      DO jv = 1, ntra
600                 SZ(iip1,j,l,jv) = SZ(1,j,l,jv)        DO j = 1, jjp1
601              END DO          sm(iip1, j, l) = sm(1, j, l)
602           END DO          s0(iip1, j, l, jv) = s0(1, j, l, jv)
603        END DO          ssx(iip1, j, l, jv) = ssx(1, j, l, jv)
604            sy(iip1, j, l, jv) = sy(1, j, l, jv)
605  c -------------------------------------------------------------------          sz(iip1, j, l, jv) = sz(1, j, l, jv)
606  C *** Test  negativite:        END DO
607        END DO
608  c      DO jv = 1,ntra    END DO
609  c       DO l = 1,llm  
610  c         DO j = 1,jjp1    ! -------------------------------------------------------------------
611  c           DO i = 1,iip1    ! *** Test  negativite:
612  c              IF (s0( i,j,l,jv ).lt.0.) THEN  
613  c                 PRINT*, '------ S0 < 0 en FIN ADVYP ---'    ! DO jv = 1,ntra
614  c                 PRINT*, 'S0(',i,j,l,jv,')=', S0(i,j,l,jv)    ! DO l = 1,llm
615  cc                 STOP    ! DO j = 1,jjp1
616  c              ENDIF    ! DO i = 1,iip1
617  c           ENDDO    ! IF (s0( i,j,l,jv ).lt.0.) THEN
618  c         ENDDO    ! PRINT*, '------ S0 < 0 en FIN ADVYP ---'
619  c       ENDDO    ! PRINT*, 'S0(',i,j,l,jv,')=', S0(i,j,l,jv)
620  c      ENDDO    ! c                 STOP
621      ! ENDIF
622        ! ENDDO
623  c -------------------------------------------------------------------    ! ENDDO
624  C *** Test : diag de la qtite totale de traceur dans    ! ENDDO
625  C            l'atmosphere avant l'advection en Y    ! ENDDO
626    
627         DO l = 1,llm  
628           DO j = 1,jjp1    ! -------------------------------------------------------------------
629             DO i = 1,iim    ! *** Test : diag de la qtite totale de traceur dans
630                sqf = sqf + S0(i,j,l,ntra)    ! l'atmosphere avant l'advection en Y
631             END DO  
632           END DO    DO l = 1, llm
633         END DO      DO j = 1, jjp1
634        PRINT*,'---------- DIAG DANS ADVY - SORTIE --------'        DO i = 1, iim
635        PRINT*,'sqf=',sqf          sqf = sqf + s0(i, j, l, ntra)
636  c     print*,'ap ADVYP fin'        END DO
637        END DO
638  c-----------------------------------------------------------------    END DO
639  C    PRINT *, '---------- DIAG DANS ADVY - SORTIE --------'
640        RETURN    PRINT *, 'sqf=', sqf
641        END    ! print*,'ap ADVYP fin'
642    
643      ! -----------------------------------------------------------------
644    
645      RETURN
646    END SUBROUTINE advyp
647    
648    
649    

Legend:
Removed from v.3  
changed lines
  Added in v.113

  ViewVC Help
Powered by ViewVC 1.1.21