/[lmdze]/trunk/Sources/dyn3d/prather.f
ViewVC logotype

Diff of /trunk/Sources/dyn3d/prather.f

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

trunk/libf/dyn3d/prather.f revision 3 by guez, Wed Feb 27 13:16:39 2008 UTC trunk/Sources/dyn3d/prather.f revision 157 by guez, Mon Jul 20 16:01:49 2015 UTC
# Line 1  Line 1 
1  !  
2  ! $Header: /home/cvsroot/LMDZ4/libf/dyn3d/prather.F,v 1.1.1.1 2004/05/19 12:53:07 lmdzadmin Exp $  ! $Header: /home/cvsroot/LMDZ4/libf/dyn3d/prather.F,v 1.1.1.1 2004/05/19
3  !  ! 12:53:07 lmdzadmin Exp $
4        SUBROUTINE prather (q,w,masse,pbaru,pbarv,nt,dt)  
5        use dimens_m  SUBROUTINE prather(q, w, masse, pbaru, pbarv, nt, dt)
6        use paramet_m  
7        use comconst    USE comconst
8        use comvert    USE dimens_m
9        use comgeom    USE disvert_m
10        IMPLICIT NONE    USE dynetat0_m, only: rlonv, rlonu
11      USE nr_util, ONLY: pi
12  c=======================================================================    USE paramet_m
13  c   Adaptation LMDZ:  A.Armengaud (LGGE)  
14  c   ----------------    IMPLICIT NONE
15  c  
16  c   ************************************************    ! =======================================================================
17  c   Transport des traceurs par la methode de prather    ! Adaptation LMDZ:  A.Armengaud (LGGE)
18  c   Ref :    ! ----------------
19  c  
20  c   ************************************************    ! ************************************************
21  c   q,w,pext,pbaru et pbarv : arguments d'entree  pour le s-pg    ! Transport des traceurs par la methode de prather
22  c    ! Ref :
23  c=======================================================================  
24      ! ************************************************
25      ! q,w,pext,pbaru et pbarv : arguments d'entree  pour le s-pg
26    
27  c   Arguments:    ! =======================================================================
28  c   ----------  
29        INTEGER iq,nt  
30        REAL pbaru( ip1jmp1,llm ),pbarv( ip1jm,llm )  
31        REAL masse(iip1,jjp1,llm)    ! Arguments:
32        REAL q( iip1,jjp1,llm,0:9)    ! ----------
33        REAL w( ip1jmp1,llm )    INTEGER nt
34        integer ordre,ilim    REAL, INTENT (IN) :: pbaru(ip1jmp1, llm), pbarv(ip1jm, llm)
35      REAL masse(iip1, jjp1, llm)
36  c   Local:    REAL q(iip1, jjp1, llm, 0:9)
37  c   ------    REAL w(ip1jmp1, llm)
38        LOGICAL limit  
39        real zq(iip1,jjp1,llm)    ! Local:
40        REAL sm ( iip1,jjp1, llm )    ! ------
41        REAL s0( iip1,jjp1,llm ),  sx( iip1,jjp1,llm )    LOGICAL limit
42        REAL sy( iip1,jjp1,llm ),  sz( iip1,jjp1,llm )    REAL sm(iip1, jjp1, llm)
43        REAL sxx( iip1,jjp1,llm)    REAL s0(iip1, jjp1, llm), sx(iip1, jjp1, llm)
44        REAL sxy( iip1,jjp1,llm)    REAL sy(iip1, jjp1, llm), sz(iip1, jjp1, llm)
45        REAL sxz( iip1,jjp1,llm)    REAL sxx(iip1, jjp1, llm)
46        REAL syy( iip1,jjp1,llm )    REAL sxy(iip1, jjp1, llm)
47        REAL syz( iip1,jjp1,llm )    REAL sxz(iip1, jjp1, llm)
48        REAL szz( iip1,jjp1,llm ),zz    REAL syy(iip1, jjp1, llm)
49        INTEGER i,j,l,indice    REAL syz(iip1, jjp1, llm)
50        real sxn(iip1),sxs(iip1)    REAL szz(iip1, jjp1, llm), zz
51      INTEGER i, j, l, indice
52        real sinlon(iip1),sinlondlon(iip1)    REAL sxn(iip1), sxs(iip1)
53        real coslon(iip1),coslondlon(iip1)  
54        real qmin,qmax    REAL sinlon(iip1), sinlondlon(iip1)
55        save qmin,qmax    REAL coslon(iip1), coslondlon(iip1)
56        save sinlon,coslon,sinlondlon,coslondlon    SAVE sinlon, coslon, sinlondlon, coslondlon
57        real dyn1,dyn2,dys1,dys2,qpn,qps,dqzpn,dqzps    REAL dyn1, dyn2, dys1, dys2, qpn, qps, dqzpn, dqzps
58        real masn,mass    REAL masn, mass
59  c  
60        REAL      SSUM    REAL ssum
61        integer ismax,ismin    INTEGER ismax, ismin
62        EXTERNAL  SSUM, convflu,ismin,ismax    EXTERNAL ssum, convflu, ismin, ismax
63        logical first    LOGICAL first
64        save first    SAVE first
65        EXTERNAL advxp,advyp,advzp    EXTERNAL advxp, advyp, advzp
66    
67    
68        data first/.true./    DATA first/.TRUE./
69        data qmin,qmax/-1.e33,1.e33/  
70      ! ==========================================================================
71      ! ==========================================================================
72  c==========================================================================    ! MODIFICATION POUR PAS DE TEMPS ADAPTATIF, dtvr remplace par dt
73  c==========================================================================    ! ==========================================================================
74  c     MODIFICATION POUR PAS DE TEMPS ADAPTATIF, dtvr remplace par dt    ! ==========================================================================
75  c==========================================================================    REAL dt
76  c==========================================================================    ! ==========================================================================
77        REAL dt    limit = .TRUE.
78  c==========================================================================  
79        limit = .TRUE.    IF (first) THEN
80        PRINT *, 'SCHEMA PRATHER'
81        if(first) then      first = .FALSE.
82           print*,'SCHEMA PRATHER'      DO i = 2, iip1
83           first=.false.        coslon(i) = cos(rlonv(i))
84           do i=2,iip1        sinlon(i) = sin(rlonv(i))
85              coslon(i)=cos(rlonv(i))        coslondlon(i) = coslon(i)*(rlonu(i)-rlonu(i-1))/pi
86              sinlon(i)=sin(rlonv(i))        sinlondlon(i) = sinlon(i)*(rlonu(i)-rlonu(i-1))/pi
87              coslondlon(i)=coslon(i)*(rlonu(i)-rlonu(i-1))/pi      END DO
88              sinlondlon(i)=sinlon(i)*(rlonu(i)-rlonu(i-1))/pi      coslon(1) = coslon(iip1)
89           enddo      coslondlon(1) = coslondlon(iip1)
90           coslon(1)=coslon(iip1)      sinlon(1) = sinlon(iip1)
91           coslondlon(1)=coslondlon(iip1)      sinlondlon(1) = sinlondlon(iip1)
92           sinlon(1)=sinlon(iip1)  
93           sinlondlon(1)=sinlondlon(iip1)      DO l = 1, llm
94          DO j = 1, jjp1
95          DO l = 1,llm          DO i = 1, iip1
96          DO j = 1,jjp1            q(i, j, l, 1) = 0.
97          DO i = 1,iip1            q(i, j, l, 2) = 0.
98          q( i,j,l,1 )=0.            q(i, j, l, 3) = 0.
99          q( i,j,l,2)=0.            q(i, j, l, 4) = 0.
100          q( i,j,l,3)=0.            q(i, j, l, 5) = 0.
101          q( i,j,l,4)=0.            q(i, j, l, 6) = 0.
102          q( i,j,l,5)=0.            q(i, j, l, 7) = 0.
103          q( i,j,l,6)=0.            q(i, j, l, 8) = 0.
104          q( i,j,l,7)=0.            q(i, j, l, 9) = 0.
105          q( i,j,l,8)=0.          END DO
106          q( i,j,l,9)=0.        END DO
107          ENDDO      END DO
108          ENDDO    END IF
109          ENDDO    ! Fin modif Fred
110        endif  
111  c   Fin modif Fred    ! *** On calcule la masse d'air en kg
112    
113  c *** On calcule la masse d'air en kg    DO l = 1, llm
114        DO j = 1, jjp1
115         DO l = 1,llm        DO i = 1, iip1
116          DO j = 1,jjp1          sm(i, j, llm+1-l) = masse(i, j, l)
117           DO i = 1,iip1        END DO
118           sm( i,j,llm+1-l ) =masse(i,j,l)      END DO
119           ENDDO    END DO
120          ENDDO  
121         ENDDO    ! *** q contient les qqtes de traceur avant l'advection
122    
123  c *** q contient les qqtes de traceur avant l'advection    ! *** Affectation des tableaux S a partir de Q
124    
125  c *** Affectation des tableaux S a partir de Q    DO l = 1, llm
126        DO j = 1, jjp1
127         DO l = 1,llm        DO i = 1, iip1
128          DO j = 1,jjp1          s0(i, j, l) = q(i, j, llm+1-l, 0)*sm(i, j, l)
129           DO i = 1,iip1          sx(i, j, l) = q(i, j, llm+1-l, 1)*sm(i, j, l)
130         s0( i,j,l) = q ( i,j,llm+1-l,0 )*sm(i,j,l)          sy(i, j, l) = q(i, j, llm+1-l, 2)*sm(i, j, l)
131         sx( i,j,l) = q( i,j,llm+1-l,1 )*sm(i,j,l)          sz(i, j, l) = q(i, j, llm+1-l, 3)*sm(i, j, l)
132         sy( i,j,l) = q( i,j,llm+1-l,2)*sm(i,j,l)          sxx(i, j, l) = q(i, j, llm+1-l, 4)*sm(i, j, l)
133         sz( i,j,l) = q( i,j,llm+1-l,3)*sm(i,j,l)          sxy(i, j, l) = q(i, j, llm+1-l, 5)*sm(i, j, l)
134         sxx( i,j,l) = q( i,j,llm+1-l,4)*sm(i,j,l)          sxz(i, j, l) = q(i, j, llm+1-l, 6)*sm(i, j, l)
135         sxy( i,j,l) = q( i,j,llm+1-l,5)*sm(i,j,l)          syy(i, j, l) = q(i, j, llm+1-l, 7)*sm(i, j, l)
136         sxz( i,j,l) = q( i,j,llm+1-l,6)*sm(i,j,l)          syz(i, j, l) = q(i, j, llm+1-l, 8)*sm(i, j, l)
137         syy( i,j,l) = q( i,j,llm+1-l,7)*sm(i,j,l)          szz(i, j, l) = q(i, j, llm+1-l, 9)*sm(i, j, l)
138         syz( i,j,l) = q( i,j,llm+1-l,8)*sm(i,j,l)        END DO
139         szz( i,j,l) = q( i,j,llm+1-l,9)*sm(i,j,l)      END DO
140           ENDDO    END DO
141          ENDDO    ! *** Appel des subroutines d'advection en X, en Y et en Z
142         ENDDO    ! *** Advection avec "time-splitting"
143  c *** Appel des subroutines d'advection en X, en Y et en Z  
144  c *** Advection avec "time-splitting"    ! -----------------------------------------------------------
145            DO indice = 1, nt
146  c-----------------------------------------------------------      CALL advxp(limit, 0.5*dt, pbaru, sm, s0, sx, sy, sz, sxx, sxy, sxz, syy, &
147         do indice =1,nt        syz, szz, 1)
148         call advxp( limit,0.5*dt,pbaru,sm,s0,sx,sy,sz    END DO
149       .             ,sxx,sxy,sxz,syy,syz,szz,1 )    DO l = 1, llm
150          end do      DO i = 1, iip1
151          do l=1,llm        sy(i, 1, l) = 0.
152          do i=1,iip1        sy(i, jjp1, l) = 0.
153          sy(i,1,l)=0.      END DO
154          sy(i,jjp1,l)=0.    END DO
155          enddo    ! ---------------------------------------------------------
156          enddo    CALL advyp(limit, .5*dt*nt, pbarv, sm, s0, sx, sy, sz, sxx, sxy, sxz, syy, &
157  c---------------------------------------------------------      syz, szz, 1)
158         call advyp( limit,.5*dt*nt,pbarv,sm,s0,sx,sy,sz    ! ---------------------------------------------------------
159       .             ,sxx,sxy,sxz,syy,syz,szz,1 )  
160  c---------------------------------------------------------    ! ---------------------------------------------------------
161      DO j = 1, jjp1
162  c---------------------------------------------------------      DO i = 1, iip1
163         do j=1,jjp1        sz(i, j, 1) = 0.
164            do i=1,iip1        sz(i, j, llm) = 0.
165               sz(i,j,1)=0.        sxz(i, j, 1) = 0.
166               sz(i,j,llm)=0.        sxz(i, j, llm) = 0.
167               sxz(i,j,1)=0.        syz(i, j, 1) = 0.
168               sxz(i,j,llm)=0.        syz(i, j, llm) = 0.
169               syz(i,j,1)=0.        szz(i, j, 1) = 0.
170               syz(i,j,llm)=0.        szz(i, j, llm) = 0.
171               szz(i,j,1)=0.      END DO
172               szz(i,j,llm)=0.    END DO
173            enddo    CALL advzp(limit, dt*nt, w, sm, s0, sx, sy, sz, sxx, sxy, sxz, syy, syz, &
174         enddo      szz, 1)
175         call advzp( limit,dt*nt,w,sm,s0,sx,sy,sz    DO l = 1, llm
176       .             ,sxx,sxy,sxz,syy,syz,szz,1 )      DO i = 1, iip1
177          do l=1,llm        sy(i, 1, l) = 0.
178          do i=1,iip1        sy(i, jjp1, l) = 0.
179          sy(i,1,l)=0.      END DO
180          sy(i,jjp1,l)=0.    END DO
181          enddo  
182          enddo    ! ---------------------------------------------------------
183    
184  c---------------------------------------------------------    ! ---------------------------------------------------------
185      CALL advyp(limit, .5*dt*nt, pbarv, sm, s0, sx, sy, sz, sxx, sxy, sxz, syy, &
186  c---------------------------------------------------------      syz, szz, 1)
187         call advyp( limit,.5*dt*nt,pbarv,sm,s0,sx,sy,sz    ! ---------------------------------------------------------
188       .             ,sxx,sxy,sxz,syy,syz,szz,1 )    DO l = 1, llm
189  c---------------------------------------------------------      DO j = 1, jjp1
190         DO l = 1,llm        s0(iip1, j, l) = s0(1, j, l)
191          DO j = 1,jjp1        sx(iip1, j, l) = sx(1, j, l)
192               s0( iip1,j,l)=s0( 1,j,l )        sy(iip1, j, l) = sy(1, j, l)
193               sx( iip1,j,l)=sx( 1,j,l )        sz(iip1, j, l) = sz(1, j, l)
194               sy( iip1,j,l)=sy( 1,j,l )        sxx(iip1, j, l) = sxx(1, j, l)
195               sz( iip1,j,l)=sz( 1,j,l )        sxy(iip1, j, l) = sxy(1, j, l)
196               sxx( iip1,j,l)=sxx( 1,j,l )        sxz(iip1, j, l) = sxz(1, j, l)
197               sxy( iip1,j,l)=sxy( 1,j,l)        syy(iip1, j, l) = syy(1, j, l)
198               sxz( iip1,j,l)=sxz( 1,j,l )        syz(iip1, j, l) = syz(1, j, l)
199               syy( iip1,j,l)=syy( 1,j,l )        szz(iip1, j, l) = szz(1, j, l)
200               syz( iip1,j,l)=syz( 1,j,l)      END DO
201               szz( iip1,j,l)=szz( 1,j,l )    END DO
202          ENDDO    DO indice = 1, nt
203         ENDDO      CALL advxp(limit, 0.5*dt, pbaru, sm, s0, sx, sy, sz, sxx, sxy, sxz, syy, &
204         do indice=1,nt        syz, szz, 1)
205         call advxp( limit,0.5*dt,pbaru,sm,s0,sx,sy,sz    END DO
206       .             ,sxx,sxy,sxz,syy,syz,szz,1 )    ! ---------------------------------------------------------
207          end do    ! ---------------------------------------------------------
208  c---------------------------------------------------------    ! ***   On repasse les S dans la variable qpr
209  c---------------------------------------------------------    ! ***   On repasse les S dans la variable q directement 14/10/94
210  c ***   On repasse les S dans la variable qpr  
211  c ***   On repasse les S dans la variable q directement 14/10/94    DO l = 1, llm
212        DO j = 1, jjp1
213         DO  l = 1,llm        DO i = 1, iip1
214          DO  j = 1,jjp1          q(i, j, llm+1-l, 0) = s0(i, j, l)/sm(i, j, l)
215           DO  i = 1,iip1          q(i, j, llm+1-l, 1) = sx(i, j, l)/sm(i, j, l)
216        q( i,j,llm+1-l,0 )=s0( i,j,l )/sm(i,j,l)          q(i, j, llm+1-l, 2) = sy(i, j, l)/sm(i, j, l)
217        q( i,j,llm+1-l,1 ) = sx( i,j,l )/sm(i,j,l)          q(i, j, llm+1-l, 3) = sz(i, j, l)/sm(i, j, l)
218        q( i,j,llm+1-l,2 ) = sy( i,j,l )/sm(i,j,l)          q(i, j, llm+1-l, 4) = sxx(i, j, l)/sm(i, j, l)
219        q( i,j,llm+1-l,3 ) = sz( i,j,l )/sm(i,j,l)          q(i, j, llm+1-l, 5) = sxy(i, j, l)/sm(i, j, l)
220        q( i,j,llm+1-l,4 ) = sxx( i,j,l )/sm(i,j,l)          q(i, j, llm+1-l, 6) = sxz(i, j, l)/sm(i, j, l)
221        q( i,j,llm+1-l,5 ) = sxy( i,j,l )/sm(i,j,l)          q(i, j, llm+1-l, 7) = syy(i, j, l)/sm(i, j, l)
222        q( i,j,llm+1-l,6 ) = sxz( i,j,l )/sm(i,j,l)          q(i, j, llm+1-l, 8) = syz(i, j, l)/sm(i, j, l)
223        q( i,j,llm+1-l,7 ) = syy( i,j,l )/sm(i,j,l)          q(i, j, llm+1-l, 9) = szz(i, j, l)/sm(i, j, l)
224        q( i,j,llm+1-l,8 ) = syz( i,j,l )/sm(i,j,l)        END DO
225        q( i,j,llm+1-l,9 ) = szz( i,j,l )/sm(i,j,l)      END DO
226        ENDDO    END DO
227        ENDDO  
228        ENDDO    ! ---------------------------------------------------------
229      ! go to  777
230  c---------------------------------------------------------    ! filtrages aux poles
231  c      go to  777  
232  c   filtrages aux poles    ! Traitements specifiques au pole
233    
234  c Traitements specifiques au pole    ! filtrages aux poles
235      DO l = 1, llm
236  c   filtrages aux poles      ! filtrages aux poles
237           DO l=1,llm      masn = ssum(iim, sm(1,1,l), 1)
238  c   filtrages aux poles      mass = ssum(iim, sm(1,jjp1,l), 1)
239           masn=ssum(iim,sm(1,1,l),1)      qpn = ssum(iim, s0(1,1,l), 1)/masn
240           mass=ssum(iim,sm(1,jjp1,l),1)      qps = ssum(iim, s0(1,jjp1,l), 1)/mass
241           qpn=ssum(iim,s0(1,1,l),1)/masn      dqzpn = ssum(iim, sz(1,1,l), 1)/masn
242           qps=ssum(iim,s0(1,jjp1,l),1)/mass      dqzps = ssum(iim, sz(1,jjp1,l), 1)/mass
243           dqzpn=ssum(iim,sz(1,1,l),1)/masn      DO i = 1, iip1
244           dqzps=ssum(iim,sz(1,jjp1,l),1)/mass        q(i, 1, llm+1-l, 3) = dqzpn
245           do i=1,iip1        q(i, jjp1, llm+1-l, 3) = dqzps
246            q( i,1,llm+1-l,3)=dqzpn        q(i, 1, llm+1-l, 0) = qpn
247            q( i,jjp1,llm+1-l,3)=dqzps        q(i, jjp1, llm+1-l, 0) = qps
248            q( i,1,llm+1-l,0)=qpn      END DO
249            q( i,jjp1,llm+1-l,0)=qps      dyn1 = 0.
250           enddo      dys1 = 0.
251  c       enddo      dyn2 = 0.
252  c         print*,'qpn',qpn,'qps',qps      dys2 = 0.
253  c          print*,'dqzpn',dqzpn,'dqzps',dqzps      DO i = 1, iim
254  c       enddo        zz = s0(i, 2, l)/sm(i, 2, l) - q(i, 1, llm+1-l, 0)
255             dyn1=0.        dyn1 = dyn1 + sinlondlon(i)*zz
256             dys1=0.        dyn2 = dyn2 + coslondlon(i)*zz
257             dyn2=0.        zz = q(i, jjp1, llm+1-l, 0) - s0(i, jjm, l)/sm(i, jjm, l)
258             dys2=0.        dys1 = dys1 + sinlondlon(i)*zz
259          do i=1,iim        dys2 = dys2 + coslondlon(i)*zz
260          zz=s0(i,2,l)/sm(i,2,l)-q(i,1,llm+1-l,0)      END DO
261          dyn1=dyn1+sinlondlon(i)*zz      DO i = 1, iim
262          dyn2=dyn2+coslondlon(i)*zz        q(i, 1, llm+1-l, 2) = (sinlon(i)*dyn1+coslon(i)*dyn2)/2.
263          zz=q(i,jjp1,llm+1-l,0)-s0(i,jjm,l)/sm(i,jjm,l)        q(i, 1, llm+1-l, 0) = q(i, 1, llm+1-l, 0) + q(i, 1, llm+1-l, 2)
264          dys1=dys1+sinlondlon(i)*zz        q(i, jjp1, llm+1-l, 2) = (sinlon(i)*dys1+coslon(i)*dys2)/2.
265          dys2=dys2+coslondlon(i)*zz        q(i, jjp1, llm+1-l, 0) = q(i, jjp1, llm+1-l, 0) - &
266          enddo          q(i, jjp1, llm+1-l, 2)
267           do i=1,iim      END DO
268           q(i,1,llm+1-l,2)=      q(iip1, 1, llm+1-l, 0) = q(1, 1, llm+1-l, 0)
269       $   (sinlon(i)*dyn1+coslon(i)*dyn2)/2.      q(iip1, jjp1, llm+1-l, 0) = q(1, jjp1, llm+1-l, 0)
270           q(i,1,llm+1-l,0)=q(i,1,llm+1-l,0)      DO i = 1, iim
271       $          +q(i,1,llm+1-l,2)        sxn(i) = q(i+1, 1, llm+1-l, 0) - q(i, 1, llm+1-l, 0)
272           q(i,jjp1,llm+1-l,2)=        sxs(i) = q(i+1, jjp1, llm+1-l, 0) - q(i, jjp1, llm+1-l, 0)
273       $   (sinlon(i)*dys1+coslon(i)*dys2)/2.      END DO
274           q(i,jjp1,llm+1-l,0)=q(i,jjp1,llm+1-l,0)      sxn(iip1) = sxn(1)
275       $      -q(i,jjp1,llm+1-l,2)      sxs(iip1) = sxs(1)
276           enddo      DO i = 1, iim
277        q(iip1,1,llm+1-l,0)=q(1,1,llm+1-l,0)        q(i+1, 1, llm+1-l, 1) = 0.25*(sxn(i)+sxn(i+1))
278        q(iip1,jjp1,llm+1-l,0)=q(1,jjp1,llm+1-l,0)        q(i+1, jjp1, llm+1-l, 1) = 0.25*(sxs(i)+sxs(i+1))
279        do i=1,iim      END DO
280        sxn(i)=q(i+1,1,llm+1-l,0)-q(i,1,llm+1-l,0)      q(1, 1, llm+1-l, 1) = q(iip1, 1, llm+1-l, 1)
281        sxs(i)=q(i+1,jjp1,llm+1-l,0)-q(i,jjp1,llm+1-l,0)      q(1, jjp1, llm+1-l, 1) = q(iip1, jjp1, llm+1-l, 1)
282        enddo    END DO
283        sxn(iip1)=sxn(1)    DO l = 1, llm
284        sxs(iip1)=sxs(1)      DO i = 1, iim
285        do i=1,iim        q(i, 1, llm+1-l, 4) = 0.
286        q(i+1,1,llm+1-l,1)=0.25*(sxn(i)+sxn(i+1))        q(i, jjp1, llm+1-l, 4) = 0.
287        q(i+1,jjp1,llm+1-l,1)=0.25*(sxs(i)+sxs(i+1))        q(i, 1, llm+1-l, 5) = 0.
288          q(i, jjp1, llm+1-l, 5) = 0.
289          q(i, 1, llm+1-l, 6) = 0.
290          q(i, jjp1, llm+1-l, 6) = 0.
291          q(i, 1, llm+1-l, 7) = 0.
292          q(i, jjp1, llm+1-l, 7) = 0.
293          q(i, 1, llm+1-l, 8) = 0.
294          q(i, jjp1, llm+1-l, 8) = 0.
295          q(i, 1, llm+1-l, 9) = 0.
296          q(i, jjp1, llm+1-l, 9) = 0.
297        END DO
298      END DO
299    
300      ! bouclage en longitude
301      DO l = 1, llm
302        DO j = 1, jjp1
303          q(iip1, j, l, 0) = q(1, j, l, 0)
304          q(iip1, j, llm+1-l, 0) = q(1, j, llm+1-l, 0)
305          q(iip1, j, llm+1-l, 1) = q(1, j, llm+1-l, 1)
306          q(iip1, j, llm+1-l, 2) = q(1, j, llm+1-l, 2)
307          q(iip1, j, llm+1-l, 3) = q(1, j, llm+1-l, 3)
308          q(iip1, j, llm+1-l, 4) = q(1, j, llm+1-l, 4)
309          q(iip1, j, llm+1-l, 5) = q(1, j, llm+1-l, 5)
310          q(iip1, j, llm+1-l, 6) = q(1, j, llm+1-l, 6)
311          q(iip1, j, llm+1-l, 7) = q(1, j, llm+1-l, 7)
312          q(iip1, j, llm+1-l, 8) = q(1, j, llm+1-l, 8)
313          q(iip1, j, llm+1-l, 9) = q(1, j, llm+1-l, 9)
314        END DO
315      END DO
316      DO l = 1, llm
317        DO j = 2, jjm
318          DO i = 1, iip1
319            IF (q(i,j,l,0)<0.) THEN
320              PRINT *, '------------ BIP-----------'
321              PRINT *, 'S0(', i, j, l, ')=', q(i, j, l, 0), q(i, j-1, l, 0)
322              PRINT *, 'SX(', i, j, l, ')=', q(i, j, l, 1)
323              PRINT *, 'SY(', i, j, l, ')=', q(i, j, l, 2), q(i, j-1, l, 2)
324              PRINT *, 'SZ(', i, j, l, ')=', q(i, j, l, 3)
325              q(i, j, l, 0) = 0.
326            END IF
327          END DO
328        END DO
329        DO j = 1, jjp1, jjm
330          DO i = 1, iip1
331            IF (q(i,j,l,0)<0.) THEN
332              PRINT *, '------------ BIP 2-----------'
333              PRINT *, 'S0(', i, j, l, ')=', q(i, j, l, 0)
334              PRINT *, 'SX(', i, j, l, ')=', q(i, j, l, 1)
335              PRINT *, 'SY(', i, j, l, ')=', q(i, j, l, 2)
336              PRINT *, 'SZ(', i, j, l, ')=', q(i, j, l, 3)
337    
338              q(i, j, l, 0) = 0.
339              ! STOP
340            END IF
341        END DO        END DO
342        q(1,1,llm+1-l,1)=q(iip1,1,llm+1-l,1)      END DO
343        q(1,jjp1,llm+1-l,1)=    END DO
344       $   q(iip1,jjp1,llm+1-l,1)    RETURN
345          enddo  END SUBROUTINE prather
          do l=1,llm  
            do i=1,iim  
             q( i,1,llm+1-l,4)=0.  
             q( i,jjp1,llm+1-l,4)=0.  
             q( i,1,llm+1-l,5)=0.  
             q( i,jjp1,llm+1-l,5)=0.  
             q( i,1,llm+1-l,6)=0.  
             q( i,jjp1,llm+1-l,6)=0.  
             q( i,1,llm+1-l,7)=0.  
             q( i,jjp1,llm+1-l,7)=0.  
             q( i,1,llm+1-l,8)=0.  
             q( i,jjp1,llm+1-l,8)=0.  
             q( i,1,llm+1-l,9)=0.  
             q( i,jjp1,llm+1-l,9)=0.  
           enddo  
          ENDDO  
   
 777      continue  
 c  
 c   bouclage en longitude  
       do l=1,llm  
       do j=1,jjp1  
       q(iip1,j,l,0)=q(1,j,l,0)  
       q(iip1,j,llm+1-l,0)=q(1,j,llm+1-l,0)  
       q(iip1,j,llm+1-l,1)=q(1,j,llm+1-l,1)  
       q(iip1,j,llm+1-l,2)=q(1,j,llm+1-l,2)  
       q(iip1,j,llm+1-l,3)=q(1,j,llm+1-l,3)  
       q(iip1,j,llm+1-l,4)=q(1,j,llm+1-l,4)  
       q(iip1,j,llm+1-l,5)=q(1,j,llm+1-l,5)  
       q(iip1,j,llm+1-l,6)=q(1,j,llm+1-l,6)  
       q(iip1,j,llm+1-l,7)=q(1,j,llm+1-l,7)  
       q(iip1,j,llm+1-l,8)=q(1,j,llm+1-l,8)  
       q(iip1,j,llm+1-l,9)=q(1,j,llm+1-l,9)  
       enddo  
       enddo  
         DO l = 1,llm  
          DO j = 2,jjm  
            DO i = 1,iip1  
          IF (q(i,j,l,0).lt.0.)  THEN  
          PRINT*,'------------ BIP-----------'  
          PRINT*,'S0(',i,j,l,')=',q(i,j,l,0),  
      $          q(i,j-1,l,0)  
          PRINT*,'SX(',i,j,l,')=',q(i,j,l,1)  
          PRINT*,'SY(',i,j,l,')=',q(i,j,l,2),  
      $   q(i,j-1,l,2)    
          PRINT*,'SZ(',i,j,l,')=',q(i,j,l,3)  
 c                    PRINT*,' PBL EN SORTIE D'' ADVZP'  
                      q(i,j,l,0)=0.  
 c                  STOP  
                ENDIF  
            ENDDO  
          ENDDO  
          do j=1,jjp1,jjm  
          do i=1,iip1  
                IF (q(i,j,l,0).lt.0.)  THEN  
                PRINT*,'------------ BIP 2-----------'  
          PRINT*,'S0(',i,j,l,')=',q(i,j,l,0)  
          PRINT*,'SX(',i,j,l,')=',q(i,j,l,1)  
          PRINT*,'SY(',i,j,l,')=',q(i,j,l,2)  
          PRINT*,'SZ(',i,j,l,')=',q(i,j,l,3)  
   
                      q(i,j,l,0)=0.  
 c                  STOP  
                ENDIF  
          enddo  
          enddo  
         ENDDO  
       RETURN  
       END  

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

  ViewVC Help
Powered by ViewVC 1.1.21