Context Navigation

advy.f90 @ 5319

Last change on this file since 5319 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: 10.9 KB

Line
1	! $Header$
2
3	SUBROUTINE advy(limit, dty, pbarv, sm, s0, sx, sy, sz)
4	USE lmdz_comgeom2
5
6	USE lmdz_dimensions, ONLY: iim, jjm, llm, ndm
7	USE lmdz_paramet
8	IMPLICIT NONE
9
10	!CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
11	! C
12	! first-order moments (SOM) advection of tracer in Y direction C
13	! C
14	! Source : Pascal Simon ( Meteo, CNRM ) C
15	! Adaptation : A.A. (LGGE) C
16	! Derniere Modif : 15/12/94 LAST
17	! C
18	! sont les arguments d'entree pour le s-pg C
19	! C
20	! argument de sortie du s-pg C
21	! C
22	!CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
23	!CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
24
25	! Rem : Probleme aux poles il faut reecrire ce cas specifique
26	! Attention au sens de l'indexation
27
28	! parametres principaux du modele
29
30	!
31
32
33
34	! Arguments :
35	! ----------
36	! dty : frequence fictive d'appel du transport
37	! parbu,pbarv : flux de masse en x et y en Pa.m2.s-1
38
39	INTEGER :: lon, lat, niv
40	INTEGER :: i, j, jv, k, kp, l
41	INTEGER :: ntra
42	PARAMETER (ntra = 1)
43
44	REAL :: dty
45	REAL :: pbarv (iip1, jjm, llm)
46
47	! moments: SM total mass in each grid box
48	! S0 mass of tracer in each grid box
49	! Si 1rst order moment in i direction
50
51	REAL :: SM(iip1, jjp1, llm) &
52	, S0(iip1, jjp1, llm, ntra)
53	REAL :: sx(iip1, jjp1, llm, ntra) &
54	, sy(iip1, jjp1, llm, ntra) &
55	, sz(iip1, jjp1, llm, ntra)
56
57
58	! Local :
59	! -------
60
61	! mass fluxes across the boundaries (UGRI,VGRI,WGRI)
62	! mass fluxes in kg
63	! declaration :
64
65	REAL :: VGRI(iip1, 0:jjp1, llm)
66
67	! Rem : UGRI et WGRI ne sont pas utilises dans
68	! cette SUBROUTINE ( advection en y uniquement )
69	! Rem 2 :le dimensionnement de VGRI depend de celui de pbarv
70
71	! the moments F are similarly defined and used as temporary
72	! storage for portions of the grid boxes in transit
73
74	REAL :: F0(iim, 0:jjp1, ntra), FM(iim, 0:jjp1)
75	REAL :: FX(iim, jjm, ntra), FY(iim, jjm, ntra)
76	REAL :: FZ(iim, jjm, ntra)
77	REAL :: S00(ntra)
78	REAL :: SM0 ! Just temporal variable
79
80	! work arrays
81
82	REAL :: ALF(iim, 0:jjp1), ALF1(iim, 0:jjp1)
83	REAL :: ALFQ(iim, 0:jjp1), ALF1Q(iim, 0:jjp1)
84	REAL :: TEMPTM ! Just temporal variable
85
86	! Special pour poles
87
88	REAL :: sbms, sfms, sfzs, sbmn, sfmn, sfzn
89	REAL :: sns0(ntra), snsz(ntra), snsm
90	REAL :: s1v(llm), slatv(llm)
91	REAL :: qy1(iim, llm, ntra), qylat(iim, llm, ntra)
92	REAL :: cx1(llm, ntra), cxLAT(llm, ntra)
93	REAL :: cy1(llm, ntra), cyLAT(llm, ntra)
94	REAL :: z1(iim), zcos(iim), zsin(iim)
95	REAL :: smpn, smps, s0pn, s0ps
96
97	REAL :: sqi, sqf
98	LOGICAL :: LIMIT
99
100	lon = iim ! rem : Il est possible qu'un pbl. arrive ici
101	lat = jjp1 ! a cause des dim. differentes entre les
102	niv = llm
103
104
105	! the moments Fi are used as temporary storage for
106	! portions of the grid boxes in transit at the current level
107
108	! work arrays
109	!
110
111	DO l = 1, llm
112	DO j = 1, jjm
113	DO i = 1, iip1
114	vgri (i, j, llm + 1 - l) = -1. * pbarv(i, j, l)
115	enddo
116	enddo
117	DO i = 1, iip1
118	vgri(i, 0, l) = 0.
119	vgri(i, jjp1, l) = 0.
120	enddo
121	enddo
122
123	DO L = 1, NIV
124
125	! place limits on appropriate moments before transport
126	! (if flux-limiting is to be applied)
127
128	IF(.NOT.LIMIT) GO TO 11
129
130	DO JV = 1, NTRA
131	DO K = 1, LAT
132	DO I = 1, LON
133	sy(I, K, L, JV) = SIGN(AMIN1(AMAX1(S0(I, K, L, JV), 0.), &
134	ABS(sy(I, K, L, JV))), sy(I, K, L, JV))
135	END DO
136	END DO
137	END DO
138
139	11 CONTINUE
140
141	! le flux a travers le pole Nord est traite separement
142
143	SM0 = 0.
144	DO JV = 1, NTRA
145	S00(JV) = 0.
146	END DO
147
148	DO I = 1, LON
149
150	IF(VGRI(I, 0, L)<=0.) THEN
151	FM(I, 0) = -VGRI(I, 0, L) * DTY
152	ALF(I, 0) = FM(I, 0) / SM(I, 1, L)
153	SM(I, 1, L) = SM(I, 1, L) - FM(I, 0)
154	SM0 = SM0 + FM(I, 0)
155	ENDIF
156
157	ALFQ(I, 0) = ALF(I, 0) * ALF(I, 0)
158	ALF1(I, 0) = 1. - ALF(I, 0)
159	ALF1Q(I, 0) = ALF1(I, 0) * ALF1(I, 0)
160
161	END DO
162
163	DO JV = 1, NTRA
164	DO I = 1, LON
165
166	IF(VGRI(I, 0, L)<=0.) THEN
167
168	F0(I, 0, JV) = ALF(I, 0) * &
169	(S0(I, 1, L, JV) - ALF1(I, 0) * sy(I, 1, L, JV))
170
171	S00(JV) = S00(JV) + F0(I, 0, JV)
172	S0(I, 1, L, JV) = S0(I, 1, L, JV) - F0(I, 0, JV)
173	sy(I, 1, L, JV) = ALF1Q(I, 0) * sy(I, 1, L, JV)
174	sx(I, 1, L, JV) = ALF1 (I, 0) * sx(I, 1, L, JV)
175	sz(I, 1, L, JV) = ALF1 (I, 0) * sz(I, 1, L, JV)
176
177	ENDIF
178
179	END DO
180	END DO
181
182	DO I = 1, LON
183	IF(VGRI(I, 0, L)>0.) THEN
184	FM(I, 0) = VGRI(I, 0, L) * DTY
185	ALF(I, 0) = FM(I, 0) / SM0
186	ENDIF
187	END DO
188
189	DO JV = 1, NTRA
190	DO I = 1, LON
191	IF(VGRI(I, 0, L)>0.) THEN
192	F0(I, 0, JV) = ALF(I, 0) * S00(JV)
193	ENDIF
194	END DO
195	END DO
196
197	! puts the temporary moments Fi into appropriate neighboring boxes
198
199	DO I = 1, LON
200
201	IF(VGRI(I, 0, L)>0.) THEN
202	SM(I, 1, L) = SM(I, 1, L) + FM(I, 0)
203	ALF(I, 0) = FM(I, 0) / SM(I, 1, L)
204	ENDIF
205
206	ALF1(I, 0) = 1. - ALF(I, 0)
207
208	END DO
209
210	DO JV = 1, NTRA
211	DO I = 1, LON
212
213	IF(VGRI(I, 0, L)>0.) THEN
214
215	TEMPTM = ALF(I, 0) * S0(I, 1, L, JV) - ALF1(I, 0) * F0(I, 0, JV)
216	S0(I, 1, L, JV) = S0(I, 1, L, JV) + F0(I, 0, JV)
217	sy(I, 1, L, JV) = ALF1(I, 0) * sy(I, 1, L, JV) + 3. * TEMPTM
218
219	ENDIF
220
221	END DO
222	END DO
223
224	! calculate flux and moments between adjacent boxes
225	! 1- create temporary moments/masses for partial boxes in transit
226	! 2- reajusts moments remaining in the box
227
228	! flux from KP to K if V(K).lt.0 and from K to KP if V(K).gt.0
229
230	DO K = 1, LAT - 1
231	KP = K + 1
232	DO I = 1, LON
233
234	IF(VGRI(I, K, L)<0.) THEN
235	FM(I, K) = -VGRI(I, K, L) * DTY
236	ALF(I, K) = FM(I, K) / SM(I, KP, L)
237	SM(I, KP, L) = SM(I, KP, L) - FM(I, K)
238	ELSE
239	FM(I, K) = VGRI(I, K, L) * DTY
240	ALF(I, K) = FM(I, K) / SM(I, K, L)
241	SM(I, K, L) = SM(I, K, L) - FM(I, K)
242	ENDIF
243
244	ALFQ(I, K) = ALF(I, K) * ALF(I, K)
245	ALF1(I, K) = 1. - ALF(I, K)
246	ALF1Q(I, K) = ALF1(I, K) * ALF1(I, K)
247
248	END DO
249	END DO
250
251	DO JV = 1, NTRA
252	DO K = 1, LAT - 1
253	KP = K + 1
254	DO I = 1, LON
255
256	IF(VGRI(I, K, L)<0.) THEN
257
258	F0(I, K, JV) = ALF (I, K) * &
259	(S0(I, KP, L, JV) - ALF1(I, K) * sy(I, KP, L, JV))
260	FY(I, K, JV) = ALFQ(I, K) * sy(I, KP, L, JV)
261	FX(I, K, JV) = ALF (I, K) * sx(I, KP, L, JV)
262	FZ(I, K, JV) = ALF (I, K) * sz(I, KP, L, JV)
263
264	S0(I, KP, L, JV) = S0(I, KP, L, JV) - F0(I, K, JV)
265	sy(I, KP, L, JV) = ALF1Q(I, K) * sy(I, KP, L, JV)
266	sx(I, KP, L, JV) = sx(I, KP, L, JV) - FX(I, K, JV)
267	sz(I, KP, L, JV) = sz(I, KP, L, JV) - FZ(I, K, JV)
268
269	ELSE
270
271	F0(I, K, JV) = ALF (I, K) * &
272	(S0(I, K, L, JV) + ALF1(I, K) * sy(I, K, L, JV))
273	FY(I, K, JV) = ALFQ(I, K) * sy(I, K, L, JV)
274	FX(I, K, JV) = ALF(I, K) * sx(I, K, L, JV)
275	FZ(I, K, JV) = ALF(I, K) * sz(I, K, L, JV)
276
277	S0(I, K, L, JV) = S0(I, K, L, JV) - F0(I, K, JV)
278	sy(I, K, L, JV) = ALF1Q(I, K) * sy(I, K, L, JV)
279	sx(I, K, L, JV) = sx(I, K, L, JV) - FX(I, K, JV)
280	sz(I, K, L, JV) = sz(I, K, L, JV) - FZ(I, K, JV)
281
282	ENDIF
283
284	END DO
285	END DO
286	END DO
287
288	! puts the temporary moments Fi into appropriate neighboring boxes
289
290	DO K = 1, LAT - 1
291	KP = K + 1
292	DO I = 1, LON
293
294	IF(VGRI(I, K, L)<0.) THEN
295	SM(I, K, L) = SM(I, K, L) + FM(I, K)
296	ALF(I, K) = FM(I, K) / SM(I, K, L)
297	ELSE
298	SM(I, KP, L) = SM(I, KP, L) + FM(I, K)
299	ALF(I, K) = FM(I, K) / SM(I, KP, L)
300	ENDIF
301
302	ALF1(I, K) = 1. - ALF(I, K)
303
304	END DO
305	END DO
306
307	DO JV = 1, NTRA
308	DO K = 1, LAT - 1
309	KP = K + 1
310	DO I = 1, LON
311
312	IF(VGRI(I, K, L)<0.) THEN
313
314	TEMPTM = -ALF(I, K) * S0(I, K, L, JV) + ALF1(I, K) * F0(I, K, JV)
315	S0(I, K, L, JV) = S0(I, K, L, JV) + F0(I, K, JV)
316	sy(I, K, L, JV) = ALF(I, K) * FY(I, K, JV) + ALF1(I, K) * sy(I, K, L, JV) &
317	+ 3. * TEMPTM
318	sx(I, K, L, JV) = sx(I, K, L, JV) + FX(I, K, JV)
319	sz(I, K, L, JV) = sz(I, K, L, JV) + FZ(I, K, JV)
320
321	ELSE
322
323	TEMPTM = ALF(I, K) * S0(I, KP, L, JV) - ALF1(I, K) * F0(I, K, JV)
324	S0(I, KP, L, JV) = S0(I, KP, L, JV) + F0(I, K, JV)
325	sy(I, KP, L, JV) = ALF(I, K) * FY(I, K, JV) + ALF1(I, K) * sy(I, KP, L, JV) &
326	+ 3. * TEMPTM
327	sx(I, KP, L, JV) = sx(I, KP, L, JV) + FX(I, K, JV)
328	sz(I, KP, L, JV) = sz(I, KP, L, JV) + FZ(I, K, JV)
329
330	ENDIF
331
332	END DO
333	END DO
334	END DO
335
336	! traitement special pour le pole Sud (idem pole Nord)
337
338	K = LAT
339
340	SM0 = 0.
341	DO JV = 1, NTRA
342	S00(JV) = 0.
343	END DO
344
345	DO I = 1, LON
346
347	IF(VGRI(I, K, L)>=0.) THEN
348	FM(I, K) = VGRI(I, K, L) * DTY
349	ALF(I, K) = FM(I, K) / SM(I, K, L)
350	SM(I, K, L) = SM(I, K, L) - FM(I, K)
351	SM0 = SM0 + FM(I, K)
352	ENDIF
353
354	ALFQ(I, K) = ALF(I, K) * ALF(I, K)
355	ALF1(I, K) = 1. - ALF(I, K)
356	ALF1Q(I, K) = ALF1(I, K) * ALF1(I, K)
357
358	END DO
359
360	DO JV = 1, NTRA
361	DO I = 1, LON
362
363	IF(VGRI(I, K, L)>=0.) THEN
364	F0 (I, K, JV) = ALF(I, K) * &
365	(S0(I, K, L, JV) + ALF1(I, K) * sy(I, K, L, JV))
366	S00(JV) = S00(JV) + F0(I, K, JV)
367
368	S0(I, K, L, JV) = S0 (I, K, L, JV) - F0 (I, K, JV)
369	sy(I, K, L, JV) = ALF1Q(I, K) * sy(I, K, L, JV)
370	sx(I, K, L, JV) = ALF1(I, K) * sx(I, K, L, JV)
371	sz(I, K, L, JV) = ALF1(I, K) * sz(I, K, L, JV)
372	ENDIF
373
374	END DO
375	END DO
376
377	DO I = 1, LON
378	IF(VGRI(I, K, L)<0.) THEN
379	FM(I, K) = -VGRI(I, K, L) * DTY
380	ALF(I, K) = FM(I, K) / SM0
381	ENDIF
382	END DO
383
384	DO JV = 1, NTRA
385	DO I = 1, LON
386	IF(VGRI(I, K, L)<0.) THEN
387	F0(I, K, JV) = ALF(I, K) * S00(JV)
388	ENDIF
389	END DO
390	END DO
391
392	! puts the temporary moments Fi into appropriate neighboring boxes
393
394	DO I = 1, LON
395
396	IF(VGRI(I, K, L)<0.) THEN
397	SM(I, K, L) = SM(I, K, L) + FM(I, K)
398	ALF(I, K) = FM(I, K) / SM(I, K, L)
399	ENDIF
400
401	ALF1(I, K) = 1. - ALF(I, K)
402
403	END DO
404
405	DO JV = 1, NTRA
406	DO I = 1, LON
407
408	IF(VGRI(I, K, L)<0.) THEN
409
410	TEMPTM = -ALF(I, K) * S0(I, K, L, JV) + ALF1(I, K) * F0(I, K, JV)
411	S0(I, K, L, JV) = S0(I, K, L, JV) + F0(I, K, JV)
412	sy(I, K, L, JV) = ALF1(I, K) * sy(I, K, L, JV) + 3. * TEMPTM
413
414	ENDIF
415
416	END DO
417	END DO
418
419	END DO
420
421	RETURN
422	END SUBROUTINE advy
423

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

Download in other formats:

Original Format