/[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 66 by guez, Thu Sep 20 13:00:41 2012 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 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 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   ------  
39        LOGICAL limit    ! Local:
40        real zq(iip1,jjp1,llm)    ! ------
41        REAL sm ( iip1,jjp1, llm )    LOGICAL limit
42        REAL s0( iip1,jjp1,llm ),  sx( iip1,jjp1,llm )    REAL sm(iip1, jjp1, llm)
43        REAL sy( iip1,jjp1,llm ),  sz( iip1,jjp1,llm )    REAL s0(iip1, jjp1, llm), sx(iip1, jjp1, llm)
44        REAL sxx( iip1,jjp1,llm)    REAL sy(iip1, jjp1, llm), sz(iip1, jjp1, llm)
45        REAL sxy( iip1,jjp1,llm)    REAL sxx(iip1, jjp1, llm)
46        REAL sxz( iip1,jjp1,llm)    REAL sxy(iip1, jjp1, llm)
47        REAL syy( iip1,jjp1,llm )    REAL sxz(iip1, jjp1, llm)
48        REAL syz( iip1,jjp1,llm )    REAL syy(iip1, jjp1, llm)
49        REAL szz( iip1,jjp1,llm ),zz    REAL syz(iip1, jjp1, llm)
50        INTEGER i,j,l,indice    REAL szz(iip1, jjp1, llm), zz
51        real sxn(iip1),sxs(iip1)    INTEGER i, j, l, indice
52      REAL sxn(iip1), sxs(iip1)
53        real sinlon(iip1),sinlondlon(iip1)  
54        real coslon(iip1),coslondlon(iip1)    REAL sinlon(iip1), sinlondlon(iip1)
55        real qmin,qmax    REAL coslon(iip1), coslondlon(iip1)
56        save qmin,qmax    SAVE sinlon, coslon, sinlondlon, coslondlon
57        save sinlon,coslon,sinlondlon,coslondlon    REAL dyn1, dyn2, dys1, dys2, qpn, qps, dqzpn, dqzps
58        real dyn1,dyn2,dys1,dys2,qpn,qps,dqzpn,dqzps    REAL masn, mass
59        real masn,mass  
60  c    REAL ssum
61        REAL      SSUM    INTEGER ismax, ismin
62        integer ismax,ismin    EXTERNAL ssum, convflu, ismin, ismax
63        EXTERNAL  SSUM, convflu,ismin,ismax    LOGICAL first
64        logical first    SAVE first
65        save first    EXTERNAL advxp, advyp, advzp
66        EXTERNAL advxp,advyp,advzp  
67    
68      DATA first/.TRUE./
69        data first/.true./  
70        data qmin,qmax/-1.e33,1.e33/    ! ==========================================================================
71      ! ==========================================================================
72      ! MODIFICATION POUR PAS DE TEMPS ADAPTATIF, dtvr remplace par dt
73  c==========================================================================    ! ==========================================================================
74  c==========================================================================    ! ==========================================================================
75  c     MODIFICATION POUR PAS DE TEMPS ADAPTATIF, dtvr remplace par dt    REAL dt
76  c==========================================================================    ! ==========================================================================
77  c==========================================================================    limit = .TRUE.
78        REAL dt  
79  c==========================================================================    IF (first) THEN
80        limit = .TRUE.      PRINT *, 'SCHEMA PRATHER'
81        first = .FALSE.
82        if(first) then      DO i = 2, iip1
83           print*,'SCHEMA PRATHER'        coslon(i) = cos(rlonv(i))
84           first=.false.        sinlon(i) = sin(rlonv(i))
85           do i=2,iip1        coslondlon(i) = coslon(i)*(rlonu(i)-rlonu(i-1))/pi
86              coslon(i)=cos(rlonv(i))        sinlondlon(i) = sinlon(i)*(rlonu(i)-rlonu(i-1))/pi
87              sinlon(i)=sin(rlonv(i))      END DO
88              coslondlon(i)=coslon(i)*(rlonu(i)-rlonu(i-1))/pi      coslon(1) = coslon(iip1)
89              sinlondlon(i)=sinlon(i)*(rlonu(i)-rlonu(i-1))/pi      coslondlon(1) = coslondlon(iip1)
90           enddo      sinlon(1) = sinlon(iip1)
91           coslon(1)=coslon(iip1)      sinlondlon(1) = sinlondlon(iip1)
92           coslondlon(1)=coslondlon(iip1)  
93           sinlon(1)=sinlon(iip1)      DO l = 1, llm
94           sinlondlon(1)=sinlondlon(iip1)        DO j = 1, jjp1
95            DO i = 1, iip1
96          DO l = 1,llm            q(i, j, l, 1) = 0.
97          DO j = 1,jjp1            q(i, j, l, 2) = 0.
98          DO i = 1,iip1            q(i, j, l, 3) = 0.
99          q( i,j,l,1 )=0.            q(i, j, l, 4) = 0.
100          q( i,j,l,2)=0.            q(i, j, l, 5) = 0.
101          q( i,j,l,3)=0.            q(i, j, l, 6) = 0.
102          q( i,j,l,4)=0.            q(i, j, l, 7) = 0.
103          q( i,j,l,5)=0.            q(i, j, l, 8) = 0.
104          q( i,j,l,6)=0.            q(i, j, l, 9) = 0.
105          q( i,j,l,7)=0.          END DO
106          q( i,j,l,8)=0.        END DO
107          q( i,j,l,9)=0.      END DO
108          ENDDO    END IF
109          ENDDO    ! Fin modif Fred
110          ENDDO  
111        endif    ! *** On calcule la masse d'air en kg
112  c   Fin modif Fred  
113      DO l = 1, llm
114  c *** On calcule la masse d'air en kg      DO j = 1, jjp1
115          DO i = 1, iip1
116         DO l = 1,llm          sm(i, j, llm+1-l) = masse(i, j, l)
117          DO j = 1,jjp1        END DO
118           DO i = 1,iip1      END DO
119           sm( i,j,llm+1-l ) =masse(i,j,l)    END DO
120           ENDDO  
121          ENDDO    ! *** q contient les qqtes de traceur avant l'advection
122         ENDDO  
123      ! *** Affectation des tableaux S a partir de Q
124  c *** q contient les qqtes de traceur avant l'advection  
125      DO l = 1, llm
126  c *** Affectation des tableaux S a partir de Q      DO j = 1, jjp1
127          DO i = 1, iip1
128         DO l = 1,llm          s0(i, j, l) = q(i, j, llm+1-l, 0)*sm(i, j, l)
129          DO j = 1,jjp1          sx(i, j, l) = q(i, j, llm+1-l, 1)*sm(i, j, l)
130           DO i = 1,iip1          sy(i, j, l) = q(i, j, llm+1-l, 2)*sm(i, j, l)
131         s0( i,j,l) = q ( i,j,llm+1-l,0 )*sm(i,j,l)          sz(i, j, l) = q(i, j, llm+1-l, 3)*sm(i, j, l)
132         sx( i,j,l) = q( i,j,llm+1-l,1 )*sm(i,j,l)          sxx(i, j, l) = q(i, j, llm+1-l, 4)*sm(i, j, l)
133         sy( i,j,l) = q( i,j,llm+1-l,2)*sm(i,j,l)          sxy(i, j, l) = q(i, j, llm+1-l, 5)*sm(i, j, l)
134         sz( i,j,l) = q( i,j,llm+1-l,3)*sm(i,j,l)          sxz(i, j, l) = q(i, j, llm+1-l, 6)*sm(i, j, l)
135         sxx( i,j,l) = q( i,j,llm+1-l,4)*sm(i,j,l)          syy(i, j, l) = q(i, j, llm+1-l, 7)*sm(i, j, l)
136         sxy( i,j,l) = q( i,j,llm+1-l,5)*sm(i,j,l)          syz(i, j, l) = q(i, j, llm+1-l, 8)*sm(i, j, l)
137         sxz( i,j,l) = q( i,j,llm+1-l,6)*sm(i,j,l)          szz(i, j, l) = q(i, j, llm+1-l, 9)*sm(i, j, l)
138         syy( i,j,l) = q( i,j,llm+1-l,7)*sm(i,j,l)        END DO
139         syz( i,j,l) = q( i,j,llm+1-l,8)*sm(i,j,l)      END DO
140         szz( i,j,l) = q( i,j,llm+1-l,9)*sm(i,j,l)    END DO
141           ENDDO    ! *** Appel des subroutines d'advection en X, en Y et en Z
142          ENDDO    ! *** Advection avec "time-splitting"
143         ENDDO  
144  c *** Appel des subroutines d'advection en X, en Y et en Z    ! -----------------------------------------------------------
145  c *** Advection avec "time-splitting"    DO indice = 1, nt
146              CALL advxp(limit, 0.5*dt, pbaru, sm, s0, sx, sy, sz, sxx, sxy, sxz, syy, &
147  c-----------------------------------------------------------        syz, szz, 1)
148         do indice =1,nt    END DO
149         call advxp( limit,0.5*dt,pbaru,sm,s0,sx,sy,sz    DO l = 1, llm
150       .             ,sxx,sxy,sxz,syy,syz,szz,1 )      DO i = 1, iip1
151          end do        sy(i, 1, l) = 0.
152          do l=1,llm        sy(i, jjp1, l) = 0.
153          do i=1,iip1      END DO
154          sy(i,1,l)=0.    END DO
155          sy(i,jjp1,l)=0.    ! ---------------------------------------------------------
156          enddo    CALL advyp(limit, .5*dt*nt, pbarv, sm, s0, sx, sy, sz, sxx, sxy, sxz, syy, &
157          enddo      syz, szz, 1)
158  c---------------------------------------------------------    ! ---------------------------------------------------------
159         call advyp( limit,.5*dt*nt,pbarv,sm,s0,sx,sy,sz  
160       .             ,sxx,sxy,sxz,syy,syz,szz,1 )    ! ---------------------------------------------------------
161  c---------------------------------------------------------    DO j = 1, jjp1
162        DO i = 1, iip1
163  c---------------------------------------------------------        sz(i, j, 1) = 0.
164         do j=1,jjp1        sz(i, j, llm) = 0.
165            do i=1,iip1        sxz(i, j, 1) = 0.
166               sz(i,j,1)=0.        sxz(i, j, llm) = 0.
167               sz(i,j,llm)=0.        syz(i, j, 1) = 0.
168               sxz(i,j,1)=0.        syz(i, j, llm) = 0.
169               sxz(i,j,llm)=0.        szz(i, j, 1) = 0.
170               syz(i,j,1)=0.        szz(i, j, llm) = 0.
171               syz(i,j,llm)=0.      END DO
172               szz(i,j,1)=0.    END DO
173               szz(i,j,llm)=0.    CALL advzp(limit, dt*nt, w, sm, s0, sx, sy, sz, sxx, sxy, sxz, syy, syz, &
174            enddo      szz, 1)
175         enddo    DO l = 1, llm
176         call advzp( limit,dt*nt,w,sm,s0,sx,sy,sz      DO i = 1, iip1
177       .             ,sxx,sxy,sxz,syy,syz,szz,1 )        sy(i, 1, l) = 0.
178          do l=1,llm        sy(i, jjp1, l) = 0.
179          do i=1,iip1      END DO
180          sy(i,1,l)=0.    END DO
181          sy(i,jjp1,l)=0.  
182          enddo    ! ---------------------------------------------------------
183          enddo  
184      ! ---------------------------------------------------------
185  c---------------------------------------------------------    CALL advyp(limit, .5*dt*nt, pbarv, sm, s0, sx, sy, sz, sxx, sxy, sxz, syy, &
186        syz, szz, 1)
187  c---------------------------------------------------------    ! ---------------------------------------------------------
188         call advyp( limit,.5*dt*nt,pbarv,sm,s0,sx,sy,sz    DO l = 1, llm
189       .             ,sxx,sxy,sxz,syy,syz,szz,1 )      DO j = 1, jjp1
190  c---------------------------------------------------------        s0(iip1, j, l) = s0(1, j, l)
191         DO l = 1,llm        sx(iip1, j, l) = sx(1, j, l)
192          DO j = 1,jjp1        sy(iip1, j, l) = sy(1, j, l)
193               s0( iip1,j,l)=s0( 1,j,l )        sz(iip1, j, l) = sz(1, j, l)
194               sx( iip1,j,l)=sx( 1,j,l )        sxx(iip1, j, l) = sxx(1, j, l)
195               sy( iip1,j,l)=sy( 1,j,l )        sxy(iip1, j, l) = sxy(1, j, l)
196               sz( iip1,j,l)=sz( 1,j,l )        sxz(iip1, j, l) = sxz(1, j, l)
197               sxx( iip1,j,l)=sxx( 1,j,l )        syy(iip1, j, l) = syy(1, j, l)
198               sxy( iip1,j,l)=sxy( 1,j,l)        syz(iip1, j, l) = syz(1, j, l)
199               sxz( iip1,j,l)=sxz( 1,j,l )        szz(iip1, j, l) = szz(1, j, l)
200               syy( iip1,j,l)=syy( 1,j,l )      END DO
201               syz( iip1,j,l)=syz( 1,j,l)    END DO
202               szz( iip1,j,l)=szz( 1,j,l )    DO indice = 1, nt
203          ENDDO      CALL advxp(limit, 0.5*dt, pbaru, sm, s0, sx, sy, sz, sxx, sxy, sxz, syy, &
204         ENDDO        syz, szz, 1)
205         do indice=1,nt    END DO
206         call advxp( limit,0.5*dt,pbaru,sm,s0,sx,sy,sz    ! ---------------------------------------------------------
207       .             ,sxx,sxy,sxz,syy,syz,szz,1 )    ! ---------------------------------------------------------
208          end do    ! ***   On repasse les S dans la variable qpr
209  c---------------------------------------------------------    ! ***   On repasse les S dans la variable q directement 14/10/94
210  c---------------------------------------------------------  
211  c ***   On repasse les S dans la variable qpr    DO l = 1, llm
212  c ***   On repasse les S dans la variable q directement 14/10/94      DO j = 1, jjp1
213          DO i = 1, iip1
214         DO  l = 1,llm          q(i, j, llm+1-l, 0) = s0(i, j, l)/sm(i, j, l)
215          DO  j = 1,jjp1          q(i, j, llm+1-l, 1) = sx(i, j, l)/sm(i, j, l)
216           DO  i = 1,iip1          q(i, j, llm+1-l, 2) = sy(i, j, l)/sm(i, j, l)
217        q( i,j,llm+1-l,0 )=s0( 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,1 ) = sx( 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,2 ) = sy( 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,3 ) = sz( 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,4 ) = sxx( 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,5 ) = sxy( 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,6 ) = sxz( 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,7 ) = syy( i,j,l )/sm(i,j,l)        END DO
225        q( i,j,llm+1-l,8 ) = syz( i,j,l )/sm(i,j,l)      END DO
226        q( i,j,llm+1-l,9 ) = szz( i,j,l )/sm(i,j,l)    END DO
227        ENDDO  
228        ENDDO    ! ---------------------------------------------------------
229        ENDDO    ! go to  777
230      ! filtrages aux poles
231  c---------------------------------------------------------  
232  c      go to  777    ! Traitements specifiques au pole
233  c   filtrages aux poles  
234      ! filtrages aux poles
235  c Traitements specifiques au pole    DO l = 1, llm
236        ! filtrages aux poles
237  c   filtrages aux poles      masn = ssum(iim, sm(1,1,l), 1)
238           DO l=1,llm      mass = ssum(iim, sm(1,jjp1,l), 1)
239  c   filtrages aux poles      qpn = ssum(iim, s0(1,1,l), 1)/masn
240           masn=ssum(iim,sm(1,1,l),1)      qps = ssum(iim, s0(1,jjp1,l), 1)/mass
241           mass=ssum(iim,sm(1,jjp1,l),1)      dqzpn = ssum(iim, sz(1,1,l), 1)/masn
242           qpn=ssum(iim,s0(1,1,l),1)/masn      dqzps = ssum(iim, sz(1,jjp1,l), 1)/mass
243           qps=ssum(iim,s0(1,jjp1,l),1)/mass      DO i = 1, iip1
244           dqzpn=ssum(iim,sz(1,1,l),1)/masn        q(i, 1, llm+1-l, 3) = dqzpn
245           dqzps=ssum(iim,sz(1,jjp1,l),1)/mass        q(i, jjp1, llm+1-l, 3) = dqzps
246           do i=1,iip1        q(i, 1, llm+1-l, 0) = qpn
247            q( i,1,llm+1-l,3)=dqzpn        q(i, jjp1, llm+1-l, 0) = qps
248            q( i,jjp1,llm+1-l,3)=dqzps      END DO
249            q( i,1,llm+1-l,0)=qpn      dyn1 = 0.
250            q( i,jjp1,llm+1-l,0)=qps      dys1 = 0.
251           enddo      dyn2 = 0.
252             dyn1=0.      dys2 = 0.
253             dys1=0.      DO i = 1, iim
254             dyn2=0.        zz = s0(i, 2, l)/sm(i, 2, l) - q(i, 1, llm+1-l, 0)
255             dys2=0.        dyn1 = dyn1 + sinlondlon(i)*zz
256          do i=1,iim        dyn2 = dyn2 + coslondlon(i)*zz
257          zz=s0(i,2,l)/sm(i,2,l)-q(i,1,llm+1-l,0)        zz = q(i, jjp1, llm+1-l, 0) - s0(i, jjm, l)/sm(i, jjm, l)
258          dyn1=dyn1+sinlondlon(i)*zz        dys1 = dys1 + sinlondlon(i)*zz
259          dyn2=dyn2+coslondlon(i)*zz        dys2 = dys2 + coslondlon(i)*zz
260          zz=q(i,jjp1,llm+1-l,0)-s0(i,jjm,l)/sm(i,jjm,l)      END DO
261          dys1=dys1+sinlondlon(i)*zz      DO i = 1, iim
262          dys2=dys2+coslondlon(i)*zz        q(i, 1, llm+1-l, 2) = (sinlon(i)*dyn1+coslon(i)*dyn2)/2.
263          enddo        q(i, 1, llm+1-l, 0) = q(i, 1, llm+1-l, 0) + q(i, 1, llm+1-l, 2)
264           do i=1,iim        q(i, jjp1, llm+1-l, 2) = (sinlon(i)*dys1+coslon(i)*dys2)/2.
265           q(i,1,llm+1-l,2)=        q(i, jjp1, llm+1-l, 0) = q(i, jjp1, llm+1-l, 0) - &
266       $   (sinlon(i)*dyn1+coslon(i)*dyn2)/2.          q(i, jjp1, llm+1-l, 2)
267           q(i,1,llm+1-l,0)=q(i,1,llm+1-l,0)      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           q(i,jjp1,llm+1-l,2)=      q(iip1, jjp1, llm+1-l, 0) = q(1, jjp1, llm+1-l, 0)
270       $   (sinlon(i)*dys1+coslon(i)*dys2)/2.      DO i = 1, iim
271           q(i,jjp1,llm+1-l,0)=q(i,jjp1,llm+1-l,0)        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           enddo      END DO
274        q(iip1,1,llm+1-l,0)=q(1,1,llm+1-l,0)      sxn(iip1) = sxn(1)
275        q(iip1,jjp1,llm+1-l,0)=q(1,jjp1,llm+1-l,0)      sxs(iip1) = sxs(1)
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)        q(i+1, 1, llm+1-l, 1) = 0.25*(sxn(i)+sxn(i+1))
278        sxs(i)=q(i+1,jjp1,llm+1-l,0)-q(i,jjp1,llm+1-l,0)        q(i+1, jjp1, llm+1-l, 1) = 0.25*(sxs(i)+sxs(i+1))
279        enddo      END DO
280        sxn(iip1)=sxn(1)      q(1, 1, llm+1-l, 1) = q(iip1, 1, llm+1-l, 1)
281        sxs(iip1)=sxs(1)      q(1, jjp1, llm+1-l, 1) = q(iip1, jjp1, llm+1-l, 1)
282        do i=1,iim    END DO
283        q(i+1,1,llm+1-l,1)=0.25*(sxn(i)+sxn(i+1))    DO l = 1, llm
284        q(i+1,jjp1,llm+1-l,1)=0.25*(sxs(i)+sxs(i+1))      DO i = 1, iim
285          q(i, 1, llm+1-l, 4) = 0.
286          q(i, jjp1, llm+1-l, 4) = 0.
287          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)  
                      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.66  
changed lines
  Added in v.157
  ViewVC Help
Powered by ViewVC 1.1.21