/[lmdze]/trunk/dyn3d/advyp.f
ViewVC logotype

Annotation of /trunk/dyn3d/advyp.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 113 - (hide annotations)
Thu Sep 18 19:56:46 2014 UTC (9 years, 8 months ago) by guez
File size: 19750 byte(s)
Moved the call to read_serre out of conf_gcm so that it can be called
only in the program ce0l, not in gcm. In gcm, variables of module
serre are read from start file. Added reading of dzoomx, dzoomy, taux,
tauy from start file, in dynetat0. Those variables were written by
dynredem0 but not read.

Removed possibility fxyhypb = false, because the geometric part of the
program is such a mess. Could then remove variables transx, transy,
alphax, alphay, pxo, pyo of module serre.

Bug fix in tau2alpha: missing save attributes. The first call to
tau2alpha needs to compute dxdyu and dxdyv regardless of value of
argument type, because they will be needed for subsequent calls to
tau2alpha with various values of argument type.

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

  ViewVC Help
Powered by ViewVC 1.1.21