Context Navigation

bilan_dyn_loc.F @ 5210

Last change on this file since 5210 was 4605, checked in by fhourdin, 17 months ago
Synchronisation des dynzon.nc a chaque ecriture
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: 23.9 KB

Rev	Line
[1632]	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	s ps,masse,pk,flux_u,flux_v,teta,phi,ucov,vcov,trac)
	6
	7	c AFAIRE
	8	c Prevoir en champ nq+1 le diagnostique de l'energie
	9	c en faisant Qzon=Cv T + L * ...
	10	c vQ..A=Cp T + L * ...
	11
	12	#ifdef CPP_IOIPSL
	13	USE IOIPSL
	14	#endif
[1823]	15	USE parallel_lmdz
[1632]	16	USE mod_hallo
	17	use misc_mod
[2475]	18	USE write_field_loc
[2597]	19	USE comconst_mod, ONLY: cpp, pi
[2600]	20	USE comvert_mod, ONLY: presnivs
[2601]	21	USE temps_mod, ONLY: annee_ref, day_ref, itau_dyn
[2600]	22
[1632]	23	IMPLICIT NONE
	24
[2597]	25	include "dimensions.h"
	26	include "paramet.h"
	27	include "comgeom2.h"
	28	include "iniprint.h"
[1632]	29
	30	c====================================================================
	31	c
	32	c Sous-programme consacre à des diagnostics dynamiques de base
	33	c
	34	c
	35	c De facon generale, les moyennes des scalaires Q sont ponderees par
	36	c la masse.
	37	c
	38	c Les flux de masse sont eux simplement moyennes.
	39	c
	40	c====================================================================
	41
	42	c Arguments :
	43	c ===========
	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	c Local :
	58	c =======
	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	cym character*6 nom(nQ)
	72	cym 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	c 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	c champ contenant les scalaires advectés.
	99	real,SAVE,ALLOCATABLE :: Q(:,:,:,:)
	100
	101	c 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	c champs de tansport en moyenne zonale
	114	integer ntr,itr
	115	parameter (ntr=5)
	116
	117	cym character*10 znom(ntr,nQ)
	118	cym character*20 znoml(ntr,nQ)
	119	cym 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*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	c Initialisation du fichier contenant les moyennes zonales.
	139	c ---------------------------------------------------------
	140
	141	character*10 infile
	142
[4605]	143	integer, save :: fileid
[1632]	144	integer thoriid, zvertiid
	145
	146	INTEGER,SAVE,ALLOCATABLE :: ndex3d(:)
	147
	148	C Variables locales
	149	C
	150	integer tau0
	151	real zjulian
	152	character*3 str
	153	character*10 ctrac
	154	integer ii,jj
	155	integer zan, dayref
	156	C
	157	real,SAVE,ALLOCATABLE :: rlong(:),rlatg(:)
	158	integer :: jjb,jje,jjn,ijb,ije
[1848]	159	type(Request),SAVE :: Req
	160	!$OMP THREADPRIVATE(Req)
[1632]	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	c=====================================================================
	175	c Initialisation
	176	c=====================================================================
	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
[2475]	217	ALLOCATE(rlong(1))
	218	ALLOCATE(rlatg(jjm))
[1632]	219
	220	!$OMP END MASTER
	221	!$OMP BARRIER
	222	icum=0
	223	c initialisation des fichiers
	224	first=.false.
	225	c ncum est la frequence de stokage en pas de temps
	226	ncum=dt_cum/dt_app
	227	if (abs(ncumdt_app-dt_cum).gt.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
[4469]	232	CALL abort_gcm("conf_gcmbilan_dyn_loc","stopped",1)
[1632]	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	c Initialisation du fichier contenant les moyennes zonales.
	254	c ---------------------------------------------------------
	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),
[2475]	286	. 1, rlong, jjn, rlatg(jjb:jje),
[1632]	287	. 1, 1, 1, jjn,
	288	. tau0, zjulian, dt_cum, thoriid, fileid,
	289	. bilan_dyn_domain_id)
	290
	291	C
	292	C Appel a histvert pour la grille verticale
	293	C
	294	call histvert(fileid, 'presnivs', 'Niveaux sigma','mb',
	295	. llm, presnivs, zvertiid)
	296	C
	297	C Appels a histdef pour la definition des variables a sauvegarder
	298	do iQ=1,nQ
	299	do itr=1,ntr
	300	if(itr.eq.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	c Declarations des champs avec dimension verticale
	313	c 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	c Declarations pour les fonctions de courant
	324	c 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	c Declarations pour les champs de transport d'air
	333	c 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	c Declarations pour les fonctions de courant
	341	c 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	c Declaration des champs 1D de transport en latitude
	348	c 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	c print*,'8HISTDEF'
	359	CALL histend(fileid)
	360
	361	!$OMP END MASTER
	362	endif
	363
	364
	365	c=====================================================================
	366	c Calcul des champs dynamiques
	367	c ----------------------------
	368
	369	jjb=jj_begin
	370	jje=jj_end
	371
	372	c é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	c$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	c 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	c=====================================================================
	405	c Cumul
	406	c=====================================================================
	407	c
	408	if(icum.EQ.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.
[1874]	421	Q_cum(:,jjb:jje,l,:)=0.
[1632]	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	c 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
[1885]	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)
[1632]	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
[1885]	455	flux_v_cum(:,jjb:jje,l)=flux_v_cum(:,jjb:jje,l)
	456	. +flux_v(:,jjb:jje,l)
[1632]	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
[1885]	466	Q_cum(:,jjb:jje,l,iQ)=Q_cum(:,jjb:jje,l,iQ)
	467	. +Q(:,jjb:jje,l,iQ)*masse(:,jjb:jje,l)
[1632]	468	ENDDO
	469	!$OMP END DO NOWAIT
	470	enddo
	471
	472	c=====================================================================
	473	c FLUX ET TENDANCES
	474	c=====================================================================
	475
	476	c Flux longitudinal
	477	c -----------------
	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	s +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	c flux méridien
	493	c -------------
	494	do iQ=1,nQ
[1874]	495	call Register_Hallo_u(Q(1,jjb_u,1,iQ),llm,0,1,1,0,Req)
[1632]	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	s +flux_v(i,j,l)0.5(Q(i,j,l,iQ)+Q(i,j+1,l,iQ))
	512	enddo
	513	enddo
	514	enddo
[2475]	515	!$OMP ENDDO NOWAIT
	516	!$OMP BARRIER
[1632]	517	enddo
	518
	519	c tendances
	520	c ---------
	521
	522	c 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	c 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
[2475]	542
[1632]	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	c=====================================================================
	584	c PAS DE TEMPS D'ECRITURE
	585	c=====================================================================
	586	if (icum.eq.ncum) then
	587	c=====================================================================
	588
	589	IF (prt_level > 5)
	590	. WRITE(lunout,*)'Pas d ecriture'
	591
	592	jjb=jj_begin
	593	jje=jj_end
	594
	595	c 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)
[2600]	600	. /masse_cum(:,jjb:jje,l)
[1632]	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
[2475]	627	!$OMP BARRIER
[1632]	628
	629	jjb=jj_begin
	630	jje=jj_end
	631
	632
	633	c A retravailler eventuellement
	634	c 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
[2475]	642
[1632]	643	c=====================================================================
	644	c Transport méridien
	645	c=====================================================================
	646
	647	c cumul zonal des masses des mailles
	648	c ----------------------------------
	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
[2475]	659	!$OMP BARRIER
[1632]	660
	661	call Register_Hallo_u(masse_cum,llm,1,1,1,1,Req)
	662	do iQ=1,nQ
[1874]	663	call Register_Hallo_u(Q_cum(1,jjb_u,1,iQ),llm,0,1,1,0,Req)
[1632]	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
[2475]	687	!$OMP BARRIER
[1632]	688
	689	c print*,'3OK'
	690	c --------------------------------------------------------------
	691	c calcul de la moyenne zonale du transport :
	692	c ------------------------------------------
	693	c
	694	c --
	695	c TOT : la circulation totale [ vq ]
	696	c
	697	c - -
	698	c MMC : mean meridional circulation [ v ] [ q ]
	699	c
	700	c ---- -- - -
	701	c TRS : transitoires [ v'q'] = [ vq ] - [ v q ]
	702	c
	703	c - * - * - - - -
	704	c STT : stationaires [ v q ] = [ v q ] - [ v ] [ q ]
	705	c
	706	c - -
	707	c on utilise aussi l'intermediaire TMP : [ v q ]
	708	c
	709	c la variable zfactv transforme un transport meridien cumule
	710	c en kg/s * unte-du-champ-transporte en m/s * unite-du-champ-transporte
	711	c
	712	c --------------------------------------------------------------
	713
	714
	715	c ----------------------------------------
	716	c Transport dans le plan latitude-altitude
	717	c ----------------------------------------
	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	c print*,'j,l,iQ=',j,l,iQ
	731	c 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	s +flux_vQ_cum(i,j,l,iQ)
	735	zqy= 0.5(Q_cum(i,j,l,iQ)masse_cum(i,j,l)+
	736	s 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	s /(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	c print*,'aOK'
	742	c 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	c 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	c 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	c print*,'4OK'
	778	c sorties proprement dites
	779	!$OMP MASTER
	780	if (i_sortie.eq.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	s zvQ(jjb:jje,:,itr,iQ)
	790	s ,jjn*llm,ndex3d)
	791	enddo
	792	call histwrite(fileid,'psi'//nom(iQ),
	793	s itau,psiQ(jjb:jje,1:llm,iQ)
	794	s ,jjn*llm,ndex3d)
	795	enddo
	796
	797	call histwrite(fileid,'masse',itau,zmasse(jjb:jje,1:llm)
	798	s ,jjn*llm,ndex3d)
	799	call histwrite(fileid,'v',itau,zv(jjb:jje,1:llm)
	800	s ,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	s jjn*llm,ndex3d)
	804
	805	endif
	806
	807
	808	c -----------------
	809	c Moyenne verticale
	810	c -----------------
	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	s +zvQ(jjb:jje,l,itr,iQ)
	823	s *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	s zavQ(jjb:jje,itr,iQ),jjn*llm,ndex3d)
	828	enddo
	829	enddo
	830	!$OMP END MASTER
	831	c on doit pouvoir tracer systematiquement la fonction de courant.
	832
	833	c=====================================================================
	834	c/////////////////////////////////////////////////////////////////////
	835	icum=0 !///////////////////////////////////////
	836	endif ! icum.eq.ncum !///////////////////////////////////////
	837	c/////////////////////////////////////////////////////////////////////
	838	c=====================================================================
[4605]	839	!$OMP MASTER
	840	call histsync(fileid)
	841	!$OMP END MASTER
[1632]	842
[4605]	843
[1632]	844	return
	845	end

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

Download in other formats:

Original Format