Context Navigation

advtrac.F @ 1171

Last change on this file since 1171 was 1114, checked in by jghattas, 16 years ago

Creation du module infotrac:

contient les variables de advtrac.h
contient la subroutine iniadvtrac renommer en infotrac_init
le nombre des traceurs est lu dans tracer.def en dynamique (ou par default ou recu par INCA)
ce module est utilise dans la dynamique et la physique
contient aussi la variable nbtr qui avant etait stockee dans dimphy

Le fichier advtrac.h n'existe plus.
La compilation ne prend plus en compte le nombre de traceur.

/JG

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 11.6 KB

Line
1	!
2	! $Header$
3	!
4	c
5	c
6	SUBROUTINE advtrac(pbaru,pbarv ,
7	* p, masse,q,iapptrac,teta,
8	* flxw,
9	* pk)
10	c Auteur : F. Hourdin
11	c
12	c Modif. P. Le Van (20/12/97)
13	c F. Codron (10/99)
14	c D. Le Croller (07/2001)
15	c M.A Filiberti (04/2002)
16	c
17	USE infotrac
18
19	IMPLICIT NONE
20	c
21	#include "dimensions.h"
22	#include "paramet.h"
23	#include "comconst.h"
24	#include "comvert.h"
25	#include "comdissip.h"
26	#include "comgeom2.h"
27	#include "logic.h"
28	#include "temps.h"
29	#include "control.h"
30	#include "ener.h"
31	#include "description.h"
32
33	c-------------------------------------------------------------------
34	c Arguments
35	c-------------------------------------------------------------------
36	c Ajout PPM
37	c--------------------------------------------------------
38	REAL massebx(ip1jmp1,llm),masseby(ip1jm,llm)
39	c--------------------------------------------------------
40	INTEGER iapptrac
41	REAL pbaru(ip1jmp1,llm),pbarv(ip1jm,llm)
42	REAL q(ip1jmp1,llm,nqtot),masse(ip1jmp1,llm)
43	REAL p( ip1jmp1,llmp1 ),teta(ip1jmp1,llm)
44	REAL pk(ip1jmp1,llm)
45	REAL flxw(ip1jmp1,llm)
46
47	c-------------------------------------------------------------
48	c Variables locales
49	c-------------------------------------------------------------
50
51	REAL pbaruc(ip1jmp1,llm),pbarvc(ip1jm,llm)
52	REAL massem(ip1jmp1,llm),zdp(ip1jmp1)
53	REAL pbarug(ip1jmp1,llm),pbarvg(ip1jm,llm),wg(ip1jmp1,llm)
54	REAL (kind=kind(1.d0)) :: t_initial, t_final, tps_cpu
55	INTEGER iadvtr
56	INTEGER ij,l,iq,iiq
57	REAL zdpmin, zdpmax
58	EXTERNAL minmax
59	SAVE iadvtr, massem, pbaruc, pbarvc
60	DATA iadvtr/0/
61	c----------------------------------------------------------
62	c Rajouts pour PPM
63	c----------------------------------------------------------
64	INTEGER indice,n
65	REAL dtbon ! Pas de temps adaptatif pour que CFL<1
66	REAL CFLmaxz,aaa,bbb ! CFL maximum
67	REAL psppm(iim,jjp1) ! pression au sol
68	REAL unatppm(iim,jjp1,llm),vnatppm(iim,jjp1,llm)
69	REAL qppm(iim*jjp1,llm,nqtot)
70	REAL fluxwppm(iim,jjp1,llm)
71	REAL apppm(llmp1), bpppm(llmp1)
72	LOGICAL dum,fill
73	DATA fill/.true./
74	DATA dum/.true./
75
76
77	IF(iadvtr.EQ.0) THEN
78	CALL initial0(ijp1llm,pbaruc)
79	CALL initial0(ijmllm,pbarvc)
80	ENDIF
81
82	c accumulation des flux de masse horizontaux
83	DO l=1,llm
84	DO ij = 1,ip1jmp1
85	pbaruc(ij,l) = pbaruc(ij,l) + pbaru(ij,l)
86	ENDDO
87	DO ij = 1,ip1jm
88	pbarvc(ij,l) = pbarvc(ij,l) + pbarv(ij,l)
89	ENDDO
90	ENDDO
91
92	c selection de la masse instantannee des mailles avant le transport.
93	IF(iadvtr.EQ.0) THEN
94
95	CALL SCOPY(ip1jmp1*llm,masse,1,massem,1)
96	ccc CALL filtreg ( massem ,jjp1, llm,-2, 2, .TRUE., 1 )
97	c
98	ENDIF
99
100	iadvtr = iadvtr+1
101	iapptrac = iadvtr
102
103
104	c Test pour savoir si on advecte a ce pas de temps
105	IF ( iadvtr.EQ.iapp_tracvl ) THEN
106
107	cc .. Modif P.Le Van ( 20/12/97 ) ....
108	cc
109
110	c traitement des flux de masse avant advection.
111	c 1. calcul de w
112	c 2. groupement des mailles pres du pole.
113
114	CALL groupe( massem, pbaruc,pbarvc, pbarug,pbarvg,wg )
115
116	! ... Flux de masse diaganostiques traceurs
117	flxw = wg / FLOAT(iapp_tracvl)
118
119	c test sur l'eventuelle creation de valeurs negatives de la masse
120	DO l=1,llm-1
121	DO ij = iip2+1,ip1jm
122	zdp(ij) = pbarug(ij-1,l) - pbarug(ij,l)
123	s - pbarvg(ij-iip1,l) + pbarvg(ij,l)
124	s + wg(ij,l+1) - wg(ij,l)
125	ENDDO
126	CALL SCOPY( jjm -1 ,zdp(iip1+iip1),iip1,zdp(iip2),iip1 )
127	DO ij = iip2,ip1jm
128	zdp(ij)= zdp(ij)*dtvr/ massem(ij,l)
129	ENDDO
130
131
132	CALL minmax ( ip1jm-iip1, zdp(iip2), zdpmin,zdpmax )
133
134	IF(MAX(ABS(zdpmin),ABS(zdpmax)).GT.0.5) THEN
135	PRINT*,'WARNING DP/P l=',l,' MIN:',zdpmin,
136	s ' MAX:', zdpmax
137	ENDIF
138
139	ENDDO
140
141	c-------------------------------------------------------------------
142	c Advection proprement dite (Modification Le Croller (07/2001)
143	c-------------------------------------------------------------------
144
145	c----------------------------------------------------
146	c Calcul des moyennes basées sur la masse
147	c----------------------------------------------------
148	call massbar(massem,massebx,masseby)
149
150	c-----------------------------------------------------------
151	c Appel des sous programmes d'advection
152	c-----------------------------------------------------------
153	do iq=1,nqtot
154	c call clock(t_initial)
155	if(iadv(iq) == 0) cycle
156	c ----------------------------------------------------------------
157	c Schema de Van Leer I MUSCL
158	c ----------------------------------------------------------------
159	if(iadv(iq).eq.10) THEN
160	call vlsplt(q(1,1,iq),2.,massem,wg,pbarug,pbarvg,dtvr)
161
162
163	c ----------------------------------------------------------------
164	c Schema "pseudo amont" + test sur humidite specifique
165	C pour la vapeur d'eau. F. Codron
166	c ----------------------------------------------------------------
167	else if(iadv(iq).eq.14) then
168	c
169	CALL vlspltqs( q(1,1,1), 2., massem, wg ,
170	* pbarug,pbarvg,dtvr,p,pk,teta )
171	c ----------------------------------------------------------------
172	c Schema de Frederic Hourdin
173	c ----------------------------------------------------------------
174	else if(iadv(iq).eq.12) then
175	c Pas de temps adaptatif
176	call adaptdt(iadv(iq),dtbon,n,pbarug,massem)
177	if (n.GT.1) then
178	write(,) 'WARNING horizontal dt=',dtbon,'dtvr=',
179	s dtvr,'n=',n
180	endif
181	do indice=1,n
182	call advn(q(1,1,iq),massem,wg,pbarug,pbarvg,dtbon,1)
183	end do
184	else if(iadv(iq).eq.13) then
185	c Pas de temps adaptatif
186	call adaptdt(iadv(iq),dtbon,n,pbarug,massem)
187	if (n.GT.1) then
188	write(,) 'WARNING horizontal dt=',dtbon,'dtvr=',
189	s dtvr,'n=',n
190	endif
191	do indice=1,n
192	call advn(q(1,1,iq),massem,wg,pbarug,pbarvg,dtbon,2)
193	end do
194	c ----------------------------------------------------------------
195	c Schema de pente SLOPES
196	c ----------------------------------------------------------------
197	else if (iadv(iq).eq.20) then
198	call pentes_ini (q(1,1,iq),wg,massem,pbarug,pbarvg,0)
199
200	c ----------------------------------------------------------------
201	c Schema de Prather
202	c ----------------------------------------------------------------
203	else if (iadv(iq).eq.30) then
204	c Pas de temps adaptatif
205	call adaptdt(iadv(iq),dtbon,n,pbarug,massem)
206	if (n.GT.1) then
207	write(,) 'WARNING horizontal dt=',dtbon,'dtvr=',
208	s dtvr,'n=',n
209	endif
210	call prather(q(1,1,iq),wg,massem,pbarug,pbarvg,
211	s n,dtbon)
212
213	c ----------------------------------------------------------------
214	c Schemas PPM Lin et Rood
215	c ----------------------------------------------------------------
216	else if (iadv(iq).eq.11.OR.(iadv(iq).GE.16.AND.
217	s iadv(iq).LE.18)) then
218
219	c Test sur le flux horizontal
220	c Pas de temps adaptatif
221	call adaptdt(iadv(iq),dtbon,n,pbarug,massem)
222	if (n.GT.1) then
223	write(,) 'WARNING horizontal dt=',dtbon,'dtvr=',
224	s dtvr,'n=',n
225	endif
226	c Test sur le flux vertical
227	CFLmaxz=0.
228	do l=2,llm
229	do ij=iip2,ip1jm
230	aaa=wg(ij,l)*dtvr/massem(ij,l)
231	CFLmaxz=max(CFLmaxz,aaa)
232	bbb=-wg(ij,l)*dtvr/massem(ij,l-1)
233	CFLmaxz=max(CFLmaxz,bbb)
234	enddo
235	enddo
236	if (CFLmaxz.GE.1) then
237	write(,) 'WARNING vertical','CFLmaxz=', CFLmaxz
238	endif
239
240	c-----------------------------------------------------------
241	c Ss-prg interface LMDZ.4->PPM3d
242	c-----------------------------------------------------------
243
244	call interpre(q(1,1,iq),qppm(1,1,iq),wg,fluxwppm,massem,
245	s apppm,bpppm,massebx,masseby,pbarug,pbarvg,
246	s unatppm,vnatppm,psppm)
247
248	do indice=1,n
249	c---------------------------------------------------------------------
250	c VL (version PPM) horiz. et PPM vert.
251	c---------------------------------------------------------------------
252	if (iadv(iq).eq.11) then
253	c Ss-prg PPM3d de Lin
254	call ppm3d(1,qppm(1,1,iq),
255	s psppm,psppm,
256	s unatppm,vnatppm,fluxwppm,dtbon,2,2,2,1,
257	s iim,jjp1,2,llm,apppm,bpppm,0.01,6400000,
258	s fill,dum,220.)
259
260	c----------------------------------------------------------------------
261	c Monotonic PPM
262	c----------------------------------------------------------------------
263	else if (iadv(iq).eq.16) then
264	c Ss-prg PPM3d de Lin
265	call ppm3d(1,qppm(1,1,iq),
266	s psppm,psppm,
267	s unatppm,vnatppm,fluxwppm,dtbon,3,3,3,1,
268	s iim,jjp1,2,llm,apppm,bpppm,0.01,6400000,
269	s fill,dum,220.)
270	c---------------------------------------------------------------------
271
272	c---------------------------------------------------------------------
273	c Semi Monotonic PPM
274	c---------------------------------------------------------------------
275	else if (iadv(iq).eq.17) then
276	c Ss-prg PPM3d de Lin
277	call ppm3d(1,qppm(1,1,iq),
278	s psppm,psppm,
279	s unatppm,vnatppm,fluxwppm,dtbon,4,4,4,1,
280	s iim,jjp1,2,llm,apppm,bpppm,0.01,6400000,
281	s fill,dum,220.)
282	c---------------------------------------------------------------------
283
284	c---------------------------------------------------------------------
285	c Positive Definite PPM
286	c---------------------------------------------------------------------
287	else if (iadv(iq).eq.18) then
288	c Ss-prg PPM3d de Lin
289	call ppm3d(1,qppm(1,1,iq),
290	s psppm,psppm,
291	s unatppm,vnatppm,fluxwppm,dtbon,5,5,5,1,
292	s iim,jjp1,2,llm,apppm,bpppm,0.01,6400000,
293	s fill,dum,220.)
294	c---------------------------------------------------------------------
295	endif
296	enddo
297	c-----------------------------------------------------------------
298	c Ss-prg interface PPM3d-LMDZ.4
299	c-----------------------------------------------------------------
300	call interpost(q(1,1,iq),qppm(1,1,iq))
301	endif
302	c----------------------------------------------------------------------
303
304	c-----------------------------------------------------------------
305	c On impose une seule valeur du traceur au pôle Sud j=jjm+1=jjp1
306	c et Nord j=1
307	c-----------------------------------------------------------------
308
309	c call traceurpole(q(1,1,iq),massem)
310
311	c calcul du temps cpu pour un schema donne
312
313	c call clock(t_final)
314	cym tps_cpu=t_final-t_initial
315	cym cpuadv(iq)=cpuadv(iq)+tps_cpu
316
317	end DO
318
319
320	c------------------------------------------------------------------
321	c on reinitialise a zero les flux de masse cumules
322	c---------------------------------------------------
323	iadvtr=0
324
325	ENDIF ! if iadvtr.EQ.iapp_tracvl
326
327	RETURN
328	END
329

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

Download in other formats:

Original Format