Context Navigation

advx.F @ 1907

Last change on this file since 1907 was 1907, checked in by lguez, 10 years ago

Added a copyright property to every file of the distribution, except
for the fcm files (which have their own copyright). Use svn propget on
a file to see the copyright. For instance:

$ svn propget copyright libf/phylmd/physiq.F90
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

Also added the files defining the CeCILL version 2 license, in French
and English, at the top of the LMDZ tree.

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: 12.1 KB

Line
1	!
2	! $Header$
3	!
4	SUBROUTINE advx(limit,dtx,pbaru,sm,s0,
5	$ sx,sy,sz,lati,latf)
6	IMPLICIT NONE
7
8	CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
9	C C
10	C first-order moments (FOM) advection of tracer in X direction C
11	C C
12	C Source : Pascal Simon (Meteo,CNRM) C
13	C Adaptation : A.Armengaud (LGGE) juin 94 C
14	C C
15	C limit,dtx,pbaru,pbarv,sm,s0,sx,sy,sz C
16	C sont des arguments d'entree pour le s-pg... C
17	C C
18	C sm,s0,sx,sy,sz C
19	C sont les arguments de sortie pour le s-pg C
20	C C
21	CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
22	C
23	C parametres principaux du modele
24	C
25	#include "dimensions.h"
26	#include "paramet.h"
27	#include "comconst.h"
28	#include "comvert.h"
29
30	C Arguments :
31	C -----------
32	C dtx : frequence fictive d'appel du transport
33	C pbaru, pbarv : flux de masse en x et y en Pa.m2.s-1
34
35	INTEGER ntra
36	PARAMETER (ntra = 1)
37
38	C ATTENTION partout ou on trouve ntra, insertion de boucle
39	C possible dans l'avenir.
40
41	REAL dtx
42	REAL pbaru ( iip1,jjp1,llm )
43
44	C moments: SM total mass in each grid box
45	C S0 mass of tracer in each grid box
46	C Si 1rst order moment in i direction
47	C
48	REAL SM(iip1,jjp1,llm),S0(iip1,jjp1,llm,ntra)
49	REAL sx(iip1,jjp1,llm,ntra)
50	$ ,sy(iip1,jjp1,llm,ntra)
51	REAL sz(iip1,jjp1,llm,ntra)
52
53	C Local :
54	C -------
55
56	C mass fluxes across the boundaries (UGRI,VGRI,WGRI)
57	C mass fluxes in kg
58	C declaration :
59
60	REAL UGRI(iip1,jjp1,llm)
61
62	C Rem : VGRI et WGRI ne sont pas utilises dans
63	C cette subroutine ( advection en x uniquement )
64	C
65	C Ti are the moments for the current latitude and level
66	C
67	REAL TM(iim)
68	REAL T0(iim,ntra),TX(iim,ntra)
69	REAL TY(iim,ntra),TZ(iim,ntra)
70	REAL TEMPTM ! just a temporary variable
71	C
72	C the moments F are similarly defined and used as temporary
73	C storage for portions of the grid boxes in transit
74	C
75	REAL FM(iim)
76	REAL F0(iim,ntra),FX(iim,ntra)
77	REAL FY(iim,ntra),FZ(iim,ntra)
78	C
79	C work arrays
80	C
81	REAL ALF(iim),ALF1(iim),ALFQ(iim),ALF1Q(iim)
82	C
83	REAL SMNEW(iim),UEXT(iim)
84	C
85	REAL sqi,sqf
86
87	LOGICAL LIMIT
88	INTEGER NUM(jjp1),LONK,NUMK
89	INTEGER lon,lati,latf,niv
90	INTEGER i,i2,i3,j,jv,l,k,itrac
91
92	lon = iim
93	niv = llm
94
95	C * Test de passage d'arguments ****
96
97
98	C -------------------------------------
99	DO 300 j = 1,jjp1
100	NUM(j) = 1
101	300 CONTINUE
102	sqi = 0.
103	sqf = 0.
104
105	DO l = 1,llm
106	DO j = 1,jjp1
107	DO i = 1,iim
108	cIM 240305 sqi = sqi + S0(i,j,l,9)
109	sqi = sqi + S0(i,j,l,ntra)
110	ENDDO
111	ENDDO
112	ENDDO
113	PRINT*,'-------- DIAG DANS ADVX - ENTREE ---------'
114	PRINT*,'sqi=',sqi
115
116
117	C Interface : adaptation nouveau modele
118	C -------------------------------------
119	C
120	C ---------------------------------------------------------
121	C Conversion des flux de masses en kg/s
122	C pbaru est en N/s d'ou :
123	C ugri est en kg/s
124
125	DO 500 l = 1,llm
126	DO 500 j = 1,jjm+1
127	DO 500 i = 1,iip1
128	C ugri (i,j,llm+1-l) = pbaru (i,j,l) * ( dsig(l) / g )
129	ugri (i,j,llm+1-l) = pbaru (i,j,l)
130	500 CONTINUE
131
132
133	C ---------------------------------------------------------
134	C ---------------------------------------------------------
135	C ---------------------------------------------------------
136
137	C start here
138	C
139	C boucle principale sur les niveaux et les latitudes
140	C
141	DO 1 L=1,NIV
142	DO 1 K=lati,latf
143	C
144	C initialisation
145	C
146	C program assumes periodic boundaries in X
147	C
148	DO 10 I=2,LON
149	SMNEW(I)=SM(I,K,L)+(UGRI(I-1,K,L)-UGRI(I,K,L))*DTX
150	10 CONTINUE
151	SMNEW(1)=SM(1,K,L)+(UGRI(LON,K,L)-UGRI(1,K,L))*DTX
152	C
153	C modifications for extended polar zones
154	C
155	NUMK=NUM(K)
156	LONK=LON/NUMK
157	C
158	IF(NUMK.GT.1) THEN
159	C
160	DO 111 I=1,LON
161	TM(I)=0.
162	111 CONTINUE
163	DO 112 JV=1,NTRA
164	DO 1120 I=1,LON
165	T0(I,JV)=0.
166	TX(I,JV)=0.
167	TY(I,JV)=0.
168	TZ(I,JV)=0.
169	1120 CONTINUE
170	112 CONTINUE
171	C
172	DO 11 I2=1,NUMK
173	C
174	DO 113 I=1,LONK
175	I3=(I-1)*NUMK+I2
176	TM(I)=TM(I)+SM(I3,K,L)
177	ALF(I)=SM(I3,K,L)/TM(I)
178	ALF1(I)=1.-ALF(I)
179	113 CONTINUE
180	C
181	DO JV=1,NTRA
182	DO I=1,LONK
183	I3=(I-1)*NUMK+I2
184	TEMPTM=-ALF(I)*T0(I,JV)+ALF1(I)
185	$ *S0(I3,K,L,JV)
186	T0(I,JV)=T0(I,JV)+S0(I3,K,L,JV)
187	TX(I,JV)=ALF(I) *sx(I3,K,L,JV)+
188	$ ALF1(I)TX(I,JV) +3.TEMPTM
189	TY(I,JV)=TY(I,JV)+sy(I3,K,L,JV)
190	TZ(I,JV)=TZ(I,JV)+sz(I3,K,L,JV)
191	ENDDO
192	ENDDO
193	C
194	11 CONTINUE
195	C
196	ELSE
197	C
198	DO 115 I=1,LON
199	TM(I)=SM(I,K,L)
200	115 CONTINUE
201	DO 116 JV=1,NTRA
202	DO 1160 I=1,LON
203	T0(I,JV)=S0(I,K,L,JV)
204	TX(I,JV)=sx(I,K,L,JV)
205	TY(I,JV)=sy(I,K,L,JV)
206	TZ(I,JV)=sz(I,K,L,JV)
207	1160 CONTINUE
208	116 CONTINUE
209	C
210	ENDIF
211	C
212	DO 117 I=1,LONK
213	UEXT(I)=UGRI(I*NUMK,K,L)
214	117 CONTINUE
215	C
216	C place limits on appropriate moments before transport
217	C (if flux-limiting is to be applied)
218	C
219	IF(.NOT.LIMIT) GO TO 13
220	C
221	DO 12 JV=1,NTRA
222	DO 120 I=1,LONK
223	TX(I,JV)=SIGN(AMIN1(AMAX1(T0(I,JV),0.),ABS(TX(I,JV))),TX(I,JV))
224	120 CONTINUE
225	12 CONTINUE
226	C
227	13 CONTINUE
228	C
229	C calculate flux and moments between adjacent boxes
230	C 1- create temporary moments/masses for partial boxes in transit
231	C 2- reajusts moments remaining in the box
232	C
233	C flux from IP to I if U(I).lt.0
234	C
235	DO 140 I=1,LONK-1
236	IF(UEXT(I).LT.0.) THEN
237	FM(I)=-UEXT(I)*DTX
238	ALF(I)=FM(I)/TM(I+1)
239	TM(I+1)=TM(I+1)-FM(I)
240	ENDIF
241	140 CONTINUE
242	C
243	I=LONK
244	IF(UEXT(I).LT.0.) THEN
245	FM(I)=-UEXT(I)*DTX
246	ALF(I)=FM(I)/TM(1)
247	TM(1)=TM(1)-FM(I)
248	ENDIF
249	C
250	C flux from I to IP if U(I).gt.0
251	C
252	DO 141 I=1,LONK
253	IF(UEXT(I).GE.0.) THEN
254	FM(I)=UEXT(I)*DTX
255	ALF(I)=FM(I)/TM(I)
256	TM(I)=TM(I)-FM(I)
257	ENDIF
258	141 CONTINUE
259	C
260	DO 142 I=1,LONK
261	ALFQ(I)=ALF(I)*ALF(I)
262	ALF1(I)=1.-ALF(I)
263	ALF1Q(I)=ALF1(I)*ALF1(I)
264	142 CONTINUE
265	C
266	DO 150 JV=1,NTRA
267	DO 1500 I=1,LONK-1
268	C
269	IF(UEXT(I).LT.0.) THEN
270	C
271	F0(I,JV)=ALF (I)* ( T0(I+1,JV)-ALF1(I)*TX(I+1,JV) )
272	FX(I,JV)=ALFQ(I)*TX(I+1,JV)
273	FY(I,JV)=ALF (I)*TY(I+1,JV)
274	FZ(I,JV)=ALF (I)*TZ(I+1,JV)
275	C
276	T0(I+1,JV)=T0(I+1,JV)-F0(I,JV)
277	TX(I+1,JV)=ALF1Q(I)*TX(I+1,JV)
278	TY(I+1,JV)=TY(I+1,JV)-FY(I,JV)
279	TZ(I+1,JV)=TZ(I+1,JV)-FZ(I,JV)
280	C
281	ENDIF
282	C
283	1500 CONTINUE
284	150 CONTINUE
285	C
286	I=LONK
287	IF(UEXT(I).LT.0.) THEN
288	C
289	DO 151 JV=1,NTRA
290	C
291	F0 (I,JV)=ALF (I)* ( T0(1,JV)-ALF1(I)*TX(1,JV) )
292	FX (I,JV)=ALFQ(I)*TX(1,JV)
293	FY (I,JV)=ALF (I)*TY(1,JV)
294	FZ (I,JV)=ALF (I)*TZ(1,JV)
295	C
296	T0(1,JV)=T0(1,JV)-F0(I,JV)
297	TX(1,JV)=ALF1Q(I)*TX(1,JV)
298	TY(1,JV)=TY(1,JV)-FY(I,JV)
299	TZ(1,JV)=TZ(1,JV)-FZ(I,JV)
300	C
301	151 CONTINUE
302	C
303	ENDIF
304	C
305	DO 152 JV=1,NTRA
306	DO 1520 I=1,LONK
307	C
308	IF(UEXT(I).GE.0.) THEN
309	C
310	F0(I,JV)=ALF (I)* ( T0(I,JV)+ALF1(I)*TX(I,JV) )
311	FX(I,JV)=ALFQ(I)*TX(I,JV)
312	FY(I,JV)=ALF (I)*TY(I,JV)
313	FZ(I,JV)=ALF (I)*TZ(I,JV)
314	C
315	T0(I,JV)=T0(I,JV)-F0(I,JV)
316	TX(I,JV)=ALF1Q(I)*TX(I,JV)
317	TY(I,JV)=TY(I,JV)-FY(I,JV)
318	TZ(I,JV)=TZ(I,JV)-FZ(I,JV)
319	C
320	ENDIF
321	C
322	1520 CONTINUE
323	152 CONTINUE
324	C
325	C puts the temporary moments Fi into appropriate neighboring boxes
326	C
327	DO 160 I=1,LONK
328	IF(UEXT(I).LT.0.) THEN
329	TM(I)=TM(I)+FM(I)
330	ALF(I)=FM(I)/TM(I)
331	ENDIF
332	160 CONTINUE
333	C
334	DO 161 I=1,LONK-1
335	IF(UEXT(I).GE.0.) THEN
336	TM(I+1)=TM(I+1)+FM(I)
337	ALF(I)=FM(I)/TM(I+1)
338	ENDIF
339	161 CONTINUE
340	C
341	I=LONK
342	IF(UEXT(I).GE.0.) THEN
343	TM(1)=TM(1)+FM(I)
344	ALF(I)=FM(I)/TM(1)
345	ENDIF
346	C
347	DO 162 I=1,LONK
348	ALF1(I)=1.-ALF(I)
349	162 CONTINUE
350	C
351	DO 170 JV=1,NTRA
352	DO 1700 I=1,LONK
353	C
354	IF(UEXT(I).LT.0.) THEN
355	C
356	TEMPTM=-ALF(I)T0(I,JV)+ALF1(I)F0(I,JV)
357	T0(I,JV)=T0(I,JV)+F0(I,JV)
358	TX(I,JV)=ALF(I)FX(I,JV)+ALF1(I)TX(I,JV)+3.*TEMPTM
359	TY(I,JV)=TY(I,JV)+FY(I,JV)
360	TZ(I,JV)=TZ(I,JV)+FZ(I,JV)
361	C
362	ENDIF
363	C
364	1700 CONTINUE
365	170 CONTINUE
366	C
367	DO 171 JV=1,NTRA
368	DO 1710 I=1,LONK-1
369	C
370	IF(UEXT(I).GE.0.) THEN
371	C
372	TEMPTM=ALF(I)T0(I+1,JV)-ALF1(I)F0(I,JV)
373	T0(I+1,JV)=T0(I+1,JV)+F0(I,JV)
374	TX(I+1,JV)=ALF(I)FX(I,JV)+ALF1(I)TX(I+1,JV)+3.*TEMPTM
375	TY(I+1,JV)=TY(I+1,JV)+FY(I,JV)
376	TZ(I+1,JV)=TZ(I+1,JV)+FZ(I,JV)
377	C
378	ENDIF
379	C
380	1710 CONTINUE
381	171 CONTINUE
382	C
383	I=LONK
384	IF(UEXT(I).GE.0.) THEN
385	DO 172 JV=1,NTRA
386	TEMPTM=ALF(I)T0(1,JV)-ALF1(I)F0(I,JV)
387	T0(1,JV)=T0(1,JV)+F0(I,JV)
388	TX(1,JV)=ALF(I)FX(I,JV)+ALF1(I)TX(1,JV)+3.*TEMPTM
389	TY(1,JV)=TY(1,JV)+FY(I,JV)
390	TZ(1,JV)=TZ(1,JV)+FZ(I,JV)
391	172 CONTINUE
392	ENDIF
393	C
394	C retour aux mailles d'origine (passage des Tij aux Sij)
395	C
396	IF(NUMK.GT.1) THEN
397	C
398	DO 180 I2=1,NUMK
399	C
400	DO 180 I=1,LONK
401	C
402	I3=I2+(I-1)*NUMK
403	SM(I3,K,L)=SMNEW(I3)
404	ALF(I)=SMNEW(I3)/TM(I)
405	TM(I)=TM(I)-SMNEW(I3)
406	C
407	ALFQ(I)=ALF(I)*ALF(I)
408	ALF1(I)=1.-ALF(I)
409	ALF1Q(I)=ALF1(I)*ALF1(I)
410	C
411	180 CONTINUE
412	C
413	DO JV=1,NTRA
414	DO I=1,LONK
415	C
416	I3=I2+(I-1)*NUMK
417	S0(I3,K,L,JV)=ALF (I)
418	$ * (T0(I,JV)-ALF1(I)*TX(I,JV))
419	sx(I3,K,L,JV)=ALFQ(I)*TX(I,JV)
420	sy(I3,K,L,JV)=ALF (I)*TY(I,JV)
421	sz(I3,K,L,JV)=ALF (I)*TZ(I,JV)
422	C
423	C reajusts moments remaining in the box
424	C
425	T0(I,JV)=T0(I,JV)-S0(I3,K,L,JV)
426	TX(I,JV)=ALF1Q(I)*TX(I,JV)
427	TY(I,JV)=TY(I,JV)-sy(I3,K,L,JV)
428	TZ(I,JV)=TZ(I,JV)-sz(I3,K,L,JV)
429	ENDDO
430	ENDDO
431	C
432	C
433	ELSE
434	C
435	DO 190 I=1,LON
436	SM(I,K,L)=TM(I)
437	190 CONTINUE
438	DO 191 JV=1,NTRA
439	DO 1910 I=1,LON
440	S0(I,K,L,JV)=T0(I,JV)
441	sx(I,K,L,JV)=TX(I,JV)
442	sy(I,K,L,JV)=TY(I,JV)
443	sz(I,K,L,JV)=TZ(I,JV)
444	1910 CONTINUE
445	191 CONTINUE
446	C
447	ENDIF
448	C
449	1 CONTINUE
450	C
451	C ----------- AA Test en fin de ADVX ------ Controle des S*
452	c OK
453	c DO 9998 l = 1, llm
454	c DO 9998 j = 1, jjp1
455	c DO 9998 i = 1, iip1
456	c IF (S0(i,j,l,ntra).lt.0..and.LIMIT) THEN
457	c PRINT*, '-------------------'
458	c PRINT*, 'En fin de ADVX'
459	c PRINT*,'SM(',i,j,l,')=',SM(i,j,l)
460	c PRINT*,'S0(',i,j,l,')=',S0(i,j,l,ntra)
461	c print*, 'sx(',i,j,l,')=',sx(i,j,l,ntra)
462	c print*, 'sy(',i,j,l,')=',sy(i,j,l,ntra)
463	c print*, 'sz(',i,j,l,')=',sz(i,j,l,ntra)
464	c WRITE (,) 'On arrete !! - pbl en fin de ADVX1'
465	cc STOP
466	c ENDIF
467	c 9998 CONTINUE
468	c
469	C ---------- bouclage cyclique
470	DO itrac=1,ntra
471	DO l = 1,llm
472	DO j = lati,latf
473	SM(iip1,j,l) = SM(1,j,l)
474	S0(iip1,j,l,itrac) = S0(1,j,l,itrac)
475	sx(iip1,j,l,itrac) = sx(1,j,l,itrac)
476	sy(iip1,j,l,itrac) = sy(1,j,l,itrac)
477	sz(iip1,j,l,itrac) = sz(1,j,l,itrac)
478	END DO
479	END DO
480	ENDDO
481
482	c ----------- qqtite totale de traceur dans tte l'atmosphere
483	DO l = 1, llm
484	DO j = 1, jjp1
485	DO i = 1, iim
486	cIM 240405 sqf = sqf + S0(i,j,l,9)
487	sqf = sqf + S0(i,j,l,ntra)
488	END DO
489	END DO
490	END DO
491	c
492	PRINT*,'------ DIAG DANS ADVX - SORTIE -----'
493	PRINT*,'sqf=',sqf
494	c-------------
495
496	RETURN
497	END
498	C_________________________________________________________________
499	C_________________________________________________________________

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

Download in other formats:

Original Format