source: LMDZ6/trunk/libf/dyn3d_common/advyp.f90 @ 5272

Last change on this file since 5272 was 5272, checked in by abarral, 23 hours ago

Turn paramet.h into a module

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