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