Context Navigation

bilan_dyn_loc.f90 @ 5411

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

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

Download in other formats:

Original Format