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