Context Navigation

bilan_dyn_loc.f90 @ 5115

Last change on this file since 5115 was 5106, checked in by abarral, 4 months ago
Turn coefils.h into lmdz_coefils.f90 Put filtreg.F90 inside lmdz_filtreg.F90 Turn mod_filtreg_p.F90 into lmdz_filtreg_p.F90 Delete obsolete parafilt.h* (lint) remove spaces between routine name and args
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
File size: 21.9 KB

Line
1
2	! $Id: bilan_dyn_p.F 1299 2010-01-20 14:27:21Z fairhead $
3
4	SUBROUTINE bilan_dyn_loc(ntrac,dt_app,dt_cum, &
5	ps,masse,pk,flux_u,flux_v,teta,phi,ucov,vcov,trac)
6
7	! AFAIRE
8	! Prevoir en champ nq+1 le diagnostique de l'energie
9	! en faisant Qzon=Cv T + L * ...
10	! vQ..A=Cp T + L * ...
11
12	USE IOIPSL
13	USE parallel_lmdz
14	USE mod_hallo
15	use misc_mod
16	USE write_field_loc
17	USE comconst_mod, ONLY: cpp, pi
18	USE comvert_mod, ONLY: presnivs
19	USE temps_mod, ONLY: annee_ref, day_ref, itau_dyn
20
21	IMPLICIT NONE
22
23	include "dimensions.h"
24	include "paramet.h"
25	include "comgeom2.h"
26	include "iniprint.h"
27
28	!====================================================================
29	!
30	! Sous-programme consacre à des diagnostics dynamiques de base
31	!
32	!
33	! De facon generale, les moyennes des scalaires Q sont ponderees par
34	! la masse.
35	!
36	! Les flux de masse sont eux simplement moyennes.
37	!
38	!====================================================================
39
40	! Arguments :
41	! ===========
42
43	integer :: ntrac
44	real :: dt_app,dt_cum
45	real :: ps(iip1,jjb_u:jje_u)
46	real :: masse(iip1,jjb_u:jje_u,llm),pk(iip1,jjb_u:jje_u,llm)
47	real :: flux_u(iip1,jjb_u:jje_u,llm)
48	real :: flux_v(iip1,jjb_v:jje_v,llm)
49	real :: teta(iip1,jjb_u:jje_u,llm)
50	real :: phi(iip1,jjb_u:jje_u,llm)
51	real :: ucov(iip1,jjb_u:jje_u,llm)
52	real :: vcov(iip1,jjb_v:jje_v,llm)
53	real :: trac(iip1,jjb_u:jje_u,llm,ntrac)
54
55	! Local :
56	! =======
57
58	integer,SAVE :: icum,ncum
59	!$OMP THREADPRIVATE(icum,ncum)
60	LOGICAL,SAVE :: first=.TRUE.
61	!$OMP THREADPRIVATE(first)
62
63	real :: zz,zqy
64	REAl,SAVE,ALLOCATABLE :: zfactv(:,:)
65
66	INTEGER,PARAMETER :: nQ=7
67
68
69	!ym character*6 nom(nQ)
70	!ym character*6 unites(nQ)
71	character(len=6),save :: nom(nQ)
72	character(len=6),save :: unites(nQ)
73
74	character(len=10) file
75	integer :: ifile
76	parameter (ifile=4)
77
78	integer,PARAMETER :: itemp=1,igeop=2,iecin=3,iang=4,iu=5
79	INTEGER,PARAMETER :: iovap=6,iun=7
80	integer,PARAMETER :: i_sortie=1
81
82	real,SAVE :: time=0.
83	integer,SAVE :: itau=0.
84	!$OMP THREADPRIVATE(time,itau)
85
86	real :: ww
87
88	! variables dynamiques intermédiaires
89	REAL,SAVE,ALLOCATABLE :: vcont(:,:,:),ucont(:,:,:)
90	REAL,SAVE,ALLOCATABLE :: ang(:,:,:),unat(:,:,:)
91	REAL,SAVE,ALLOCATABLE :: massebx(:,:,:),masseby(:,:,:)
92	REAL,SAVE,ALLOCATABLE :: vorpot(:,:,:)
93	REAL,SAVE,ALLOCATABLE :: w(:,:,:),ecin(:,:,:),convm(:,:,:)
94	REAL,SAVE,ALLOCATABLE :: bern(:,:,:)
95
96	! champ contenant les scalaires advectés.
97	real,SAVE,ALLOCATABLE :: Q(:,:,:,:)
98
99	! champs cumulés
100	real,SAVE,ALLOCATABLE :: ps_cum(:,:)
101	real,SAVE,ALLOCATABLE :: masse_cum(:,:,:)
102	real,SAVE,ALLOCATABLE :: flux_u_cum(:,:,:)
103	real,SAVE,ALLOCATABLE :: flux_v_cum(:,:,:)
104	real,SAVE,ALLOCATABLE :: Q_cum(:,:,:,:)
105	real,SAVE,ALLOCATABLE :: flux_uQ_cum(:,:,:,:)
106	real,SAVE,ALLOCATABLE :: flux_vQ_cum(:,:,:,:)
107	real,SAVE,ALLOCATABLE :: flux_wQ_cum(:,:,:,:)
108	real,SAVE,ALLOCATABLE :: dQ(:,:,:,:)
109
110
111	! champs de tansport en moyenne zonale
112	integer :: ntr,itr
113	parameter (ntr=5)
114
115	!ym character*10 znom(ntr,nQ)
116	!ym character*20 znoml(ntr,nQ)
117	!ym character*10 zunites(ntr,nQ)
118	character*10,save :: znom(ntr,nQ)
119	character*20,save :: znoml(ntr,nQ)
120	character*10,save :: zunites(ntr,nQ)
121
122	INTEGER,PARAMETER :: iave=1,itot=2,immc=3,itrs=4,istn=5
123
124	character(len=3) :: ctrs(ntr)
125	data ctrs/' ','TOT','MMC','TRS','STN'/
126
127	real,SAVE,ALLOCATABLE :: zvQ(:,:,:,:),zvQtmp(:,:)
128	real,SAVE,ALLOCATABLE :: zavQ(:,:,:),psiQ(:,:,:)
129	real,SAVE,ALLOCATABLE :: zmasse(:,:),zamasse(:)
130
131	real,SAVE,ALLOCATABLE :: zv(:,:),psi(:,:)
132
133	integer :: i,j,l,iQ
134
135
136	! Initialisation du fichier contenant les moyennes zonales.
137	! ---------------------------------------------------------
138
139	character(len=10) :: infile
140
141	integer, save :: fileid
142	integer :: thoriid, zvertiid
143
144	INTEGER,SAVE,ALLOCATABLE :: ndex3d(:)
145
146	! Variables locales
147	!
148	integer :: tau0
149	real :: zjulian
150	character(len=3) :: str
151	character(len=10) :: ctrac
152	integer :: ii,jj
153	integer :: zan, dayref
154	!
155	real,SAVE,ALLOCATABLE :: rlong(:),rlatg(:)
156	integer :: jjb,jje,jjn,ijb,ije
157	type(Request),SAVE :: Req
158	!$OMP THREADPRIVATE(Req)
159
160	! definition du domaine d'ecriture pour le rebuild
161
162	INTEGER,DIMENSION(1) :: ddid
163	INTEGER,DIMENSION(1) :: dsg
164	INTEGER,DIMENSION(1) :: dsl
165	INTEGER,DIMENSION(1) :: dpf
166	INTEGER,DIMENSION(1) :: dpl
167	INTEGER,DIMENSION(1) :: dhs
168	INTEGER,DIMENSION(1) :: dhe
169
170	INTEGER :: bilan_dyn_domain_id
171
172	!=====================================================================
173	! Initialisation
174	!=====================================================================
175	if (adjust) return
176
177	time=time+dt_app
178	itau=itau+1
179
180	if (first) then
181	!$OMP BARRIER
182	!$OMP MASTER
183	ALLOCATE(zfactv(jjb_v:jje_v,llm))
184	ALLOCATE(vcont(iip1,jjb_v:jje_v,llm))
185	ALLOCATE(ucont(iip1,jjb_u:jje_u,llm))
186	ALLOCATE(ang(iip1,jjb_u:jje_u,llm))
187	ALLOCATE(unat(iip1,jjb_u:jje_u,llm))
188	ALLOCATE(massebx(iip1,jjb_u:jje_u,llm))
189	ALLOCATE(masseby(iip1,jjb_v:jje_v,llm))
190	ALLOCATE(vorpot(iip1,jjb_v:jje_v,llm))
191	ALLOCATE(w(iip1,jjb_u:jje_u,llm))
192	ALLOCATE(ecin(iip1,jjb_u:jje_u,llm))
193	ALLOCATE(convm(iip1,jjb_u:jje_u,llm))
194	ALLOCATE(bern(iip1,jjb_u:jje_u,llm))
195	ALLOCATE(Q(iip1,jjb_u:jje_u,llm,nQ))
196	ALLOCATE(ps_cum(iip1,jjb_u:jje_u))
197	ALLOCATE(masse_cum(iip1,jjb_u:jje_u,llm))
198	ALLOCATE(flux_u_cum(iip1,jjb_u:jje_u,llm))
199	ALLOCATE(flux_v_cum(iip1,jjb_v:jje_v,llm))
200	ALLOCATE(Q_cum(iip1,jjb_u:jje_u,llm,nQ))
201	ALLOCATE(flux_uQ_cum(iip1,jjb_u:jje_u,llm,nQ))
202	ALLOCATE(flux_vQ_cum(iip1,jjb_v:jje_v,llm,nQ))
203	ALLOCATE(flux_wQ_cum(iip1,jjb_u:jje_u,llm,nQ))
204	ALLOCATE(dQ(iip1,jjb_u:jje_u,llm,nQ))
205	ALLOCATE(zvQ(jjb_v:jje_v,llm,ntr,nQ))
206	ALLOCATE(zvQtmp(jjb_v:jje_v,llm))
207	ALLOCATE(zavQ(jjb_v:jje_v,ntr,nQ))
208	ALLOCATE(psiQ(jjb_v:jje_v,llm+1,nQ))
209	ALLOCATE(zmasse(jjb_v:jje_v,llm))
210	ALLOCATE(zamasse(jjb_v:jje_v))
211	ALLOCATE(zv(jjb_v:jje_v,llm))
212	ALLOCATE(psi(jjb_v:jje_v,llm+1))
213	ALLOCATE(ndex3d(jjb_v:jje_v*llm))
214	ndex3d=0
215	ALLOCATE(rlong(1))
216	ALLOCATE(rlatg(jjm))
217
218	!$OMP END MASTER
219	!$OMP BARRIER
220	icum=0
221	! initialisation des fichiers
222	first=.FALSE.
223	! ncum est la frequence de stokage en pas de temps
224	ncum=dt_cum/dt_app
225	if (abs(ncumdt_app-dt_cum)>1.e-5dt_app) then
226	WRITE(lunout,*) &
227	'Pb : le pas de cumule doit etre multiple du pas'
228	WRITE(lunout,*)'dt_app=',dt_app
229	WRITE(lunout,*)'dt_cum=',dt_cum
230	CALL abort_gcm("conf_gcmbilan_dyn_loc","stopped",1)
231	endif
232
233	!$OMP MASTER
234	nom(itemp)='T'
235	nom(igeop)='gz'
236	nom(iecin)='K'
237	nom(iang)='ang'
238	nom(iu)='u'
239	nom(iovap)='ovap'
240	nom(iun)='un'
241
242	unites(itemp)='K'
243	unites(igeop)='m2/s2'
244	unites(iecin)='m2/s2'
245	unites(iang)='ang'
246	unites(iu)='m/s'
247	unites(iovap)='kg/kg'
248	unites(iun)='un'
249
250
251	! Initialisation du fichier contenant les moyennes zonales.
252	! ---------------------------------------------------------
253
254	infile='dynzon'
255
256	zan = annee_ref
257	dayref = day_ref
258	CALL ymds2ju(zan, 1, dayref, 0.0, zjulian)
259	tau0 = itau_dyn
260
261	rlong=0.
262	rlatg=rlatv*180./pi
263
264	jjb=jj_begin
265	jje=jj_end
266	jjn=jj_nb
267	IF (pole_sud) THEN
268	jjn=jj_nb-1
269	jje=jj_end-1
270	ENDIF
271
272	ddid=(/ 2 /)
273	dsg=(/ jjm /)
274	dsl=(/ jjn /)
275	dpf=(/ jjb /)
276	dpl=(/ jje /)
277	dhs=(/ 0 /)
278	dhe=(/ 0 /)
279
280	CALL flio_dom_set(mpi_size,mpi_rank,ddid,dsg,dsl,dpf,dpl,dhs,dhe, &
281	'box',bilan_dyn_domain_id)
282
283	CALL histbeg(trim(infile), &
284	1, rlong, jjn, rlatg(jjb:jje), &
285	1, 1, 1, jjn, &
286	tau0, zjulian, dt_cum, thoriid, fileid, &
287	bilan_dyn_domain_id)
288
289	!
290	! Appel a histvert pour la grille verticale
291	!
292	CALL histvert(fileid, 'presnivs', 'Niveaux sigma','mb', &
293	llm, presnivs, zvertiid)
294	!
295	! Appels a histdef pour la definition des variables a sauvegarder
296	do iQ=1,nQ
297	do itr=1,ntr
298	if(itr==1) then
299	znom(itr,iQ)=nom(iQ)
300	znoml(itr,iQ)=nom(iQ)
301	zunites(itr,iQ)=unites(iQ)
302	else
303	znom(itr,iQ)=ctrs(itr)//'v'//nom(iQ)
304	znoml(itr,iQ)='transport : v * '//nom(iQ)//' '//ctrs(itr)
305	zunites(itr,iQ)='m/s * '//unites(iQ)
306	endif
307	enddo
308	enddo
309
310	! Declarations des champs avec dimension verticale
311	! PRINT*,'1HISTDEF'
312	do iQ=1,nQ
313	do itr=1,ntr
314	IF (prt_level > 5) &
315	WRITE(lunout,*)'var ',itr,iQ &
316	,znom(itr,iQ),znoml(itr,iQ),zunites(itr,iQ)
317	CALL histdef(fileid,znom(itr,iQ),znoml(itr,iQ), &
318	zunites(itr,iQ),1,jjn,thoriid,llm,1,llm,zvertiid, &
319	32,'ave(X)',dt_cum,dt_cum)
320	enddo
321	! Declarations pour les fonctions de courant
322	! PRINT*,'2HISTDEF'
323	CALL histdef(fileid,'psi'//nom(iQ) &
324	,'stream fn. '//znoml(itot,iQ), &
325	zunites(itot,iQ),1,jjn,thoriid,llm,1,llm,zvertiid, &
326	32,'ave(X)',dt_cum,dt_cum)
327	enddo
328
329
330	! Declarations pour les champs de transport d'air
331	! PRINT*,'3HISTDEF'
332	CALL histdef(fileid, 'masse', 'masse', &
333	'kg', 1, jjn, thoriid, llm, 1, llm, zvertiid, &
334	32, 'ave(X)', dt_cum, dt_cum)
335	CALL histdef(fileid, 'v', 'v', &
336	'm/s', 1, jjn, thoriid, llm, 1, llm, zvertiid, &
337	32, 'ave(X)', dt_cum, dt_cum)
338	! Declarations pour les fonctions de courant
339	! PRINT*,'4HISTDEF'
340	CALL histdef(fileid,'psi','stream fn. MMC ','mega t/s', &
341	1,jjn,thoriid,llm,1,llm,zvertiid, &
342	32,'ave(X)',dt_cum,dt_cum)
343
344
345	! Declaration des champs 1D de transport en latitude
346	! PRINT*,'5HISTDEF'
347	do iQ=1,nQ
348	do itr=2,ntr
349	CALL histdef(fileid,'a'//znom(itr,iQ),znoml(itr,iQ), &
350	zunites(itr,iQ),1,jjn,thoriid,1,1,1,-99, &
351	32,'ave(X)',dt_cum,dt_cum)
352	enddo
353	enddo
354
355
356	! PRINT*,'8HISTDEF'
357	CALL histend(fileid)
358
359	!$OMP END MASTER
360	endif
361
362
363	!=====================================================================
364	! Calcul des champs dynamiques
365	! ----------------------------
366
367	jjb=jj_begin
368	jje=jj_end
369
370	! énergie cinétique
371	! ucont(:,jjb:jje,:)=0
372
373	CALL Register_Hallo_u(ucov,llm,1,1,1,1,Req)
374	CALL Register_Hallo_v(vcov,llm,1,1,1,1,Req)
375	CALL SendRequest(Req)
376	!$OMP BARRIER
377	CALL WaitRequest(Req)
378
379	CALL covcont_loc(llm,ucov,vcov,ucont,vcont)
380	CALL enercin_loc(vcov,ucov,vcont,ucont,ecin)
381
382	! moment cinétique
383	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
384	do l=1,llm
385	ang(:,jjb:jje,l)=ucov(:,jjb:jje,l)+constang(:,jjb:jje)
386	unat(:,jjb:jje,l)=ucont(:,jjb:jje,l)*cu(:,jjb:jje)
387	enddo
388	!$OMP END DO NOWAIT
389
390	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
391	DO l=1,llm
392	Q(:,jjb:jje,l,itemp)=teta(:,jjb:jje,l)*pk(:,jjb:jje,l)/cpp
393	Q(:,jjb:jje,l,igeop)=phi(:,jjb:jje,l)
394	Q(:,jjb:jje,l,iecin)=ecin(:,jjb:jje,l)
395	Q(:,jjb:jje,l,iang)=ang(:,jjb:jje,l)
396	Q(:,jjb:jje,l,iu)=unat(:,jjb:jje,l)
397	Q(:,jjb:jje,l,iovap)=trac(:,jjb:jje,l,1)
398	Q(:,jjb:jje,l,iun)=1.
399	ENDDO
400	!$OMP END DO NOWAIT
401
402	!=====================================================================
403	! Cumul
404	!=====================================================================
405	!
406	if(icum==0) then
407	jjb=jj_begin
408	jje=jj_end
409
410	!$OMP MASTER
411	ps_cum(:,jjb:jje)=0.
412	!$OMP END MASTER
413
414
415	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
416	DO l=1,llm
417	masse_cum(:,jjb:jje,l)=0.
418	flux_u_cum(:,jjb:jje,l)=0.
419	Q_cum(:,jjb:jje,l,:)=0.
420	flux_uQ_cum(:,jjb:jje,l,:)=0.
421	if (pole_sud) jje=jj_end-1
422	flux_v_cum(:,jjb:jje,l)=0.
423	flux_vQ_cum(:,jjb:jje,l,:)=0.
424	ENDDO
425	!$OMP END DO NOWAIT
426	endif
427
428	IF (prt_level > 5) &
429	WRITE(lunout,*)'dans bilan_dyn ',icum,'->',icum+1
430	icum=icum+1
431
432	! accumulation des flux de masse horizontaux
433	jjb=jj_begin
434	jje=jj_end
435
436	!$OMP MASTER
437	ps_cum(:,jjb:jje)=ps_cum(:,jjb:jje)+ps(:,jjb:jje)
438	!$OMP END MASTER
439
440
441	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
442	DO l=1,llm
443	masse_cum(:,jjb:jje,l)=masse_cum(:,jjb:jje,l)+masse(:,jjb:jje,l)
444	flux_u_cum(:,jjb:jje,l)=flux_u_cum(:,jjb:jje,l) &
445	+flux_u(:,jjb:jje,l)
446	ENDDO
447	!$OMP END DO NOWAIT
448
449	if (pole_sud) jje=jj_end-1
450
451	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
452	DO l=1,llm
453	flux_v_cum(:,jjb:jje,l)=flux_v_cum(:,jjb:jje,l) &
454	+flux_v(:,jjb:jje,l)
455	ENDDO
456	!$OMP END DO NOWAIT
457
458	jjb=jj_begin
459	jje=jj_end
460
461	do iQ=1,nQ
462	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
463	DO l=1,llm
464	Q_cum(:,jjb:jje,l,iQ)=Q_cum(:,jjb:jje,l,iQ) &
465	+Q(:,jjb:jje,l,iQ)*masse(:,jjb:jje,l)
466	ENDDO
467	!$OMP END DO NOWAIT
468	enddo
469
470	!=====================================================================
471	! FLUX ET TENDANCES
472	!=====================================================================
473
474	! Flux longitudinal
475	! -----------------
476	do iQ=1,nQ
477	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
478	do l=1,llm
479	do j=jjb,jje
480	do i=1,iim
481	flux_uQ_cum(i,j,l,iQ)=flux_uQ_cum(i,j,l,iQ) &
482	+flux_u(i,j,l)0.5(Q(i,j,l,iQ)+Q(i+1,j,l,iQ))
483	enddo
484	flux_uQ_cum(iip1,j,l,iQ)=flux_uQ_cum(1,j,l,iQ)
485	enddo
486	enddo
487	!$OMP END DO NOWAIT
488	enddo
489
490	! flux méridien
491	! -------------
492	do iQ=1,nQ
493	CALL Register_Hallo_u(Q(1,jjb_u,1,iQ),llm,0,1,1,0,Req)
494	enddo
495	CALL SendRequest(Req)
496	!$OMP BARRIER
497	CALL WaitRequest(Req)
498
499	jjb=jj_begin
500	jje=jj_end
501	if (pole_sud) jje=jj_end-1
502
503	do iQ=1,nQ
504	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
505	do l=1,llm
506	do j=jjb,jje
507	do i=1,iip1
508	flux_vQ_cum(i,j,l,iQ)=flux_vQ_cum(i,j,l,iQ) &
509	+flux_v(i,j,l)0.5(Q(i,j,l,iQ)+Q(i,j+1,l,iQ))
510	enddo
511	enddo
512	enddo
513	!$OMP ENDDO NOWAIT
514	!$OMP BARRIER
515	enddo
516
517	! tendances
518	! ---------
519
520	! convergence horizontale
521	CALL Register_Hallo_u(flux_uQ_cum,llm,2,2,2,2,Req)
522	CALL Register_Hallo_v(flux_vQ_cum,llm,2,2,2,2,Req)
523	CALL SendRequest(Req)
524	!$OMP BARRIER
525	CALL WaitRequest(Req)
526
527	CALL convflu_loc(flux_uQ_cum,flux_vQ_cum,llm*nQ,dQ)
528
529	! calcul de la vitesse verticale
530	CALL Register_Hallo_u(flux_u_cum,llm,2,2,2,2,Req)
531	CALL Register_Hallo_v(flux_v_cum,llm,2,2,2,2,Req)
532	CALL SendRequest(Req)
533	!$OMP BARRIER
534	CALL WaitRequest(Req)
535
536	CALL convmas_loc(flux_u_cum,flux_v_cum,convm)
537	CALL vitvert_loc(convm,w)
538	!$OMP BARRIER
539
540
541	jjb=jj_begin
542	jje=jj_end
543
544	! do iQ=1,nQ
545	! do l=1,llm-1
546	! do j=jjb,jje
547	! do i=1,iip1
548	! ww=-0.5w(i,j,l+1)(Q(i,j,l,iQ)+Q(i,j,l+1,iQ))
549	! dQ(i,j,l ,iQ)=dQ(i,j,l ,iQ)-ww
550	! dQ(i,j,l+1,iQ)=dQ(i,j,l+1,iQ)+ww
551	! enddo
552	! enddo
553	! enddo
554	! enddo
555
556	do iQ=1,nQ
557	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
558	do l=1,llm
559	IF (l<llm) THEN
560	do j=jjb,jje
561	do i=1,iip1
562	ww=-0.5w(i,j,l+1)(Q(i,j,l,iQ)+Q(i,j,l+1,iQ))
563	dQ(i,j,l ,iQ)=dQ(i,j,l ,iQ)-ww
564	dQ(i,j,l+1,iQ)=dQ(i,j,l+1,iQ)+ww
565	enddo
566	enddo
567	ENDIF
568	IF (l>2) THEN
569	do j=jjb,jje
570	do i=1,iip1
571	ww=-0.5w(i,j,l)(Q(i,j,l-1,iQ)+Q(i,j,l,iQ))
572	dQ(i,j,l,iQ)=dQ(i,j,l,iQ)+ww
573	enddo
574	enddo
575	ENDIF
576	enddo
577	!$OMP ENDDO NOWAIT
578	enddo
579	IF (prt_level > 5) &
580	WRITE(lunout,*)'Apres les calculs fait a chaque pas'
581	!=====================================================================
582	! PAS DE TEMPS D'ECRITURE
583	!=====================================================================
584	if (icum==ncum) then
585	!=====================================================================
586
587	IF (prt_level > 5) &
588	WRITE(lunout,*)'Pas d ecriture'
589
590	jjb=jj_begin
591	jje=jj_end
592
593	! Normalisation
594	do iQ=1,nQ
595	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
596	do l=1,llm
597	Q_cum(:,jjb:jje,l,iQ)=Q_cum(:,jjb:jje,l,iQ) &
598	/masse_cum(:,jjb:jje,l)
599	enddo
600	!$OMP ENDDO NOWAIT
601	enddo
602
603	zz=1./REAL(ncum)
604
605	!$OMP MASTER
606	ps_cum(:,jjb:jje)=ps_cum(:,jjb:jje)*zz
607	!$OMP END MASTER
608
609	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
610	DO l=1,llm
611	masse_cum(:,jjb:jje,l)=masse_cum(:,jjb:jje,l)*zz
612	flux_u_cum(:,jjb:jje,l)=flux_u_cum(:,jjb:jje,l)*zz
613	flux_uQ_cum(:,jjb:jje,l,:)=flux_uQ_cum(:,jjb:jje,l,:)*zz
614	dQ(:,jjb:jje,l,:)=dQ(:,jjb:jje,l,:)*zz
615	ENDDO
616	!$OMP ENDDO NOWAIT
617
618	IF (pole_sud) jje=jj_end-1
619	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
620	DO l=1,llm
621	flux_v_cum(:,jjb:jje,l)=flux_v_cum(:,jjb:jje,l)*zz
622	flux_vQ_cum(:,jjb:jje,l,:)=flux_vQ_cum(:,jjb:jje,l,:)*zz
623	ENDDO
624	!$OMP ENDDO NOWAIT
625	!$OMP BARRIER
626
627	jjb=jj_begin
628	jje=jj_end
629
630
631	! A retravailler eventuellement
632	! division de dQ par la masse pour revenir aux bonnes grandeurs
633	do iQ=1,nQ
634	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
635	DO l=1,llm
636	dQ(:,jjb:jje,l,iQ)=dQ(:,jjb:jje,l,iQ)/masse_cum(:,jjb:jje,l)
637	ENDDO
638	!$OMP ENDDO NOWAIT
639	enddo
640
641	!=====================================================================
642	! Transport méridien
643	!=====================================================================
644
645	! cumul zonal des masses des mailles
646	! ----------------------------------
647	jjb=jj_begin
648	jje=jj_end
649	if (pole_sud) jje=jj_end-1
650
651	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
652	DO l=1,llm
653	zv(jjb:jje,l)=0.
654	zmasse(jjb:jje,l)=0.
655	ENDDO
656	!$OMP ENDDO NOWAIT
657	!$OMP BARRIER
658
659	CALL Register_Hallo_u(masse_cum,llm,1,1,1,1,Req)
660	do iQ=1,nQ
661	CALL Register_Hallo_u(Q_cum(1,jjb_u,1,iQ),llm,0,1,1,0,Req)
662	enddo
663
664	CALL SendRequest(Req)
665	!$OMP BARRIER
666	CALL WaitRequest(Req)
667
668	CALL massbar_loc(masse_cum,massebx,masseby)
669
670	jjb=jj_begin
671	jje=jj_end
672	if (pole_sud) jje=jj_end-1
673
674	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
675	do l=1,llm
676	do j=jjb,jje
677	do i=1,iim
678	zmasse(j,l)=zmasse(j,l)+masseby(i,j,l)
679	zv(j,l)=zv(j,l)+flux_v_cum(i,j,l)
680	enddo
681	zfactv(j,l)=cv(1,j)/zmasse(j,l)
682	enddo
683	enddo
684	!$OMP ENDDO NOWAIT
685	!$OMP BARRIER
686
687	! PRINT*,'3OK'
688	! --------------------------------------------------------------
689	! calcul de la moyenne zonale du transport :
690	! ------------------------------------------
691	!
692	! --
693	! TOT : la circulation totale [ vq ]
694	!
695	! - -
696	! MMC : mean meridional circulation [ v ] [ q ]
697	!
698	! ---- -- - -
699	! TRS : transitoires [ v'q'] = [ vq ] - [ v q ]
700	!
701	! - * - * - - - -
702	! STT : stationaires [ v q ] = [ v q ] - [ v ] [ q ]
703	!
704	! - -
705	! on utilise aussi l'intermediaire TMP : [ v q ]
706	!
707	! la variable zfactv transforme un transport meridien cumule
708	! en kg/s * unte-du-champ-transporte en m/s * unite-du-champ-transporte
709	!
710	! --------------------------------------------------------------
711
712
713	! ----------------------------------------
714	! Transport dans le plan latitude-altitude
715	! ----------------------------------------
716
717	jjb=jj_begin
718	jje=jj_end
719	if (pole_sud) jje=jj_end-1
720
721	zvQ=0.
722	psiQ=0.
723	do iQ=1,nQ
724	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
725	do l=1,llm
726	zvQtmp(:,l)=0.
727	do j=jjb,jje
728	! PRINT*,'j,l,iQ=',j,l,iQ
729	! Calcul des moyennes zonales du transort total et de zvQtmp
730	do i=1,iim
731	zvQ(j,l,itot,iQ)=zvQ(j,l,itot,iQ) &
732	+flux_vQ_cum(i,j,l,iQ)
733	zqy= 0.5(Q_cum(i,j,l,iQ)masse_cum(i,j,l)+ &
734	Q_cum(i,j+1,l,iQ)*masse_cum(i,j+1,l))
735	zvQtmp(j,l)=zvQtmp(j,l)+flux_v_cum(i,j,l)*zqy &
736	/(0.5*(masse_cum(i,j,l)+masse_cum(i,j+1,l)))
737	zvQ(j,l,iave,iQ)=zvQ(j,l,iave,iQ)+zqy
738	enddo
739	! PRINT*,'aOK'
740	! Decomposition
741	zvQ(j,l,iave,iQ)=zvQ(j,l,iave,iQ)/zmasse(j,l)
742	zvQ(j,l,itot,iQ)=zvQ(j,l,itot,iQ)*zfactv(j,l)
743	zvQtmp(j,l)=zvQtmp(j,l)*zfactv(j,l)
744	zvQ(j,l,immc,iQ)=zv(j,l)zvQ(j,l,iave,iQ)zfactv(j,l)
745	zvQ(j,l,itrs,iQ)=zvQ(j,l,itot,iQ)-zvQtmp(j,l)
746	zvQ(j,l,istn,iQ)=zvQtmp(j,l)-zvQ(j,l,immc,iQ)
747	enddo
748	enddo
749	!$OMP ENDDO NOWAIT
750	! fonction de courant meridienne pour la quantite Q
751	!$OMP BARRIER
752	!$OMP MASTER
753	do l=llm,1,-1
754	do j=jjb,jje
755	psiQ(j,l,iQ)=psiQ(j,l+1,iQ)+zvQ(j,l,itot,iQ)
756	enddo
757	enddo
758	!$OMP END MASTER
759	!$OMP BARRIER
760	enddo
761
762	! fonction de courant pour la circulation meridienne moyenne
763	!$OMP BARRIER
764	!$OMP MASTER
765	psi(jjb:jje,:)=0.
766	do l=llm,1,-1
767	do j=jjb,jje
768	psi(j,l)=psi(j,l+1)+zv(j,l)
769	zv(j,l)=zv(j,l)*zfactv(j,l)
770	enddo
771	enddo
772	!$OMP END MASTER
773	!$OMP BARRIER
774
775	! PRINT*,'4OK'
776	! sorties proprement dites
777	!$OMP MASTER
778	if (i_sortie==1) then
779	jjb=jj_begin
780	jje=jj_end
781	jjn=jj_nb
782	if (pole_sud) jje=jj_end-1
783	if (pole_sud) jjn=jj_nb-1
784	do iQ=1,nQ
785	do itr=1,ntr
786	CALL histwrite(fileid,znom(itr,iQ),itau, &
787	zvQ(jjb:jje,:,itr,iQ) &
788	,jjn*llm,ndex3d)
789	enddo
790	CALL histwrite(fileid,'psi'//nom(iQ), &
791	itau,psiQ(jjb:jje,1:llm,iQ) &
792	,jjn*llm,ndex3d)
793	enddo
794
795	CALL histwrite(fileid,'masse',itau,zmasse(jjb:jje,1:llm) &
796	,jjn*llm,ndex3d)
797	CALL histwrite(fileid,'v',itau,zv(jjb:jje,1:llm) &
798	,jjn*llm,ndex3d)
799	psi(jjb:jje,:)=psi(jjb:jje,:)*1.e-9
800	CALL histwrite(fileid,'psi',itau,psi(jjb:jje,1:llm), &
801	jjn*llm,ndex3d)
802
803	endif
804
805
806	! -----------------
807	! Moyenne verticale
808	! -----------------
809
810	zamasse(jjb:jje)=0.
811	do l=1,llm
812	zamasse(jjb:jje)=zamasse(jjb:jje)+zmasse(jjb:jje,l)
813	enddo
814
815	zavQ(jjb:jje,:,:)=0.
816	do iQ=1,nQ
817	do itr=2,ntr
818	do l=1,llm
819	zavQ(jjb:jje,itr,iQ)=zavQ(jjb:jje,itr,iQ) &
820	+zvQ(jjb:jje,l,itr,iQ) &
821	*zmasse(jjb:jje,l)
822	enddo
823	zavQ(jjb:jje,itr,iQ)=zavQ(jjb:jje,itr,iQ)/zamasse(jjb:jje)
824	CALL histwrite(fileid,'a'//znom(itr,iQ),itau, &
825	zavQ(jjb:jje,itr,iQ),jjn*llm,ndex3d)
826	enddo
827	enddo
828	!$OMP END MASTER
829	! on doit pouvoir tracer systematiquement la fonction de courant.
830
831	!=====================================================================
832	!/////////////////////////////////////////////////////////////////////
833	icum=0 !///////////////////////////////////////
834	endif ! icum.eq.ncum !///////////////////////////////////////
835	!/////////////////////////////////////////////////////////////////////
836	!=====================================================================
837	!$OMP MASTER
838	CALL histsync(fileid)
839	!$OMP END MASTER
840
841
842
843	end subroutine bilan_dyn_loc

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

Download in other formats:

Original Format