Context Navigation

advyp.f90 @ 5201

Last change on this file since 5201 was 5159, checked in by abarral, 3 months ago
Put dimensions.h and paramet.h into modules
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: 20.8 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format