Context Navigation

vlsplt_p.F @ 5166

Last change on this file since 5166 was 2603, checked in by Ehouarn Millour, 8 years ago
Cleanup in the dynamics: turn logic.h into module logic_mod.F90 EM
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: 29.6 KB

Rev	Line
[1520]	1	c
	2	c $Id: vlsplt_p.F 2603 2016-07-25 09:31:56Z idelkadi $
	3	c
	4
[630]	5	SUBROUTINE vlsplt_p(q,pente_max,masse,w,pbaru,pbarv,pdt)
	6	c
	7	c Auteurs: P.Le Van, F.Hourdin, F.Forget
	8	c
	9	c ********************************************************************
	10	c Shema d'advection " pseudo amont " .
	11	c ********************************************************************
	12	c q,pbaru,pbarv,w sont des arguments d'entree pour le s-pg ....
	13	c
	14	c pente_max facteur de limitation des pentes: 2 en general
	15	c 0 pour un schema amont
	16	c pbaru,pbarv,w flux de masse en u ,v ,w
	17	c pdt pas de temps
	18	c
	19	c --------------------------------------------------------------------
[1823]	20	USE parallel_lmdz
[630]	21	USE mod_hallo
	22	USE Vampir
	23	IMPLICIT NONE
	24	c
	25	#include "dimensions.h"
	26	#include "paramet.h"
	27
	28	c
	29	c Arguments:
	30	c ----------
	31	REAL masse(ip1jmp1,llm),pente_max
	32	c REAL masse(iip1,jjp1,llm),pente_max
	33	REAL pbaru( ip1jmp1,llm ),pbarv( ip1jm,llm)
	34	REAL q(ip1jmp1,llm)
	35	c REAL q(iip1,jjp1,llm)
	36	REAL w(ip1jmp1,llm),pdt
	37	c
	38	c Local
	39	c ---------
	40	c
	41	INTEGER i,ij,l,j,ii
	42	INTEGER ijlqmin,iqmin,jqmin,lqmin
	43	c
	44	REAL zm(ip1jmp1,llm),newmasse
	45	REAL mu(ip1jmp1,llm)
	46	REAL mv(ip1jm,llm)
	47	REAL mw(ip1jmp1,llm+1)
	48	REAL zq(ip1jmp1,llm),zz
	49	REAL dqx(ip1jmp1,llm),dqy(ip1jmp1,llm),dqz(ip1jmp1,llm)
	50	REAL second,temps0,temps1,temps2,temps3
	51	REAL ztemps1,ztemps2,ztemps3
	52	REAL zzpbar, zzw
	53	LOGICAL testcpu
	54	SAVE testcpu
	55	SAVE temps1,temps2,temps3
	56	INTEGER iminn,imaxx
	57
	58	REAL qmin,qmax
	59	DATA qmin,qmax/0.,1.e33/
	60	DATA testcpu/.false./
	61	DATA temps1,temps2,temps3/0.,0.,0./
	62	INTEGER ijb,ije
	63	type(request) :: MyRequest1
	64	type(request) :: MyRequest2
	65
	66	call SetTag(MyRequest1,100)
	67	call SetTag(MyRequest2,101)
	68
	69	zzpbar = 0.5 * pdt
	70	zzw = pdt
	71
	72	ijb=ij_begin
	73	ije=ij_end
	74	if (pole_nord) ijb=ijb+iip1
	75	if (pole_sud) ije=ije-iip1
	76
	77	DO l=1,llm
	78	DO ij = ijb,ije
	79	mu(ij,l)=pbaru(ij,l) * zzpbar
	80	ENDDO
	81	ENDDO
	82
	83	ijb=ij_begin-iip1
	84	ije=ij_end
	85	if (pole_nord) ijb=ij_begin
	86	if (pole_sud) ije=ij_end-iip1
	87
	88	DO l=1,llm
	89	DO ij=ijb,ije
	90	mv(ij,l)=pbarv(ij,l) * zzpbar
	91	ENDDO
	92	ENDDO
	93
	94	ijb=ij_begin
	95	ije=ij_end
	96
	97	DO l=1,llm
	98	DO ij=ijb,ije
	99	mw(ij,l)=w(ij,l) * zzw
	100	ENDDO
	101	ENDDO
	102
	103	DO ij=ijb,ije
	104	mw(ij,llm+1)=0.
	105	ENDDO
	106
	107	c CALL SCOPY(ijp1llm,q,1,zq,1)
	108	c CALL SCOPY(ijp1llm,masse,1,zm,1)
	109
	110	ijb=ij_begin
	111	ije=ij_end
	112	zq(ijb:ije,:)=q(ijb:ije,:)
	113	zm(ijb:ije,:)=masse(ijb:ije,:)
	114
	115
	116	c print*,'Entree vlx1'
	117	c call minmaxq(zq,qmin,qmax,'avant vlx ')
	118	call vlx_p(zq,pente_max,zm,mu,ij_begin,ij_begin+2*iip1-1)
	119	call vlx_p(zq,pente_max,zm,mu,ij_end-2*iip1+1,ij_end)
	120	call VTb(VTHallo)
	121	call Register_Hallo(zq,ip1jmp1,llm,2,2,2,2,MyRequest1)
	122	call Register_Hallo(zm,ip1jmp1,llm,1,1,1,1,MyRequest1)
	123	call SendRequest(MyRequest1)
	124	call VTe(VTHallo)
	125	call vlx_p(zq,pente_max,zm,mu,ij_begin+2iip1,ij_end-2iip1)
	126	c call vlx_p(zq,pente_max,zm,mu,ij_begin,ij_end)
	127	call VTb(VTHallo)
	128	call WaitRecvRequest(MyRequest1)
	129	call VTe(VTHallo)
	130
	131
	132	c print*,'Sortie vlx1'
	133	c call minmaxq(zq,qmin,qmax,'apres vlx1 ')
	134
	135	c print*,'Entree vly1'
	136	c call exchange_hallo(zq,ip1jmp1,llm,2,2)
	137	c call exchange_hallo(zm,ip1jmp1,llm,1,1)
	138
	139	call vly_p(zq,pente_max,zm,mv)
	140	c call minmaxq(zq,qmin,qmax,'apres vly1 ')
	141	c print*,'Sortie vly1'
	142	call vlz_p(zq,pente_max,zm,mw,ij_begin,ij_begin+2*iip1-1)
	143	call vlz_p(zq,pente_max,zm,mw,ij_end-2*iip1+1,ij_end)
	144	call VTb(VTHallo)
	145	call Register_Hallo(zq,ip1jmp1,llm,2,2,2,2,MyRequest2)
	146	call Register_Hallo(zm,ip1jmp1,llm,1,1,1,1,MyRequest2)
	147	call SendRequest(MyRequest2)
	148	call VTe(VTHallo)
	149	call vlz_p(zq,pente_max,zm,mw,ij_begin+2iip1,ij_end-2iip1)
	150	call VTb(VTHallo)
	151	call WaitRecvRequest(MyRequest2)
	152
	153	call VTe(VTHallo)
	154
	155	c call minmaxq(zq,qmin,qmax,'apres vlz ')
	156
	157
	158
	159
	160	call vly_p(zq,pente_max,zm,mv)
	161	c call minmaxq(zq,qmin,qmax,'apres vly ')
	162
	163
	164	call vlx_p(zq,pente_max,zm,mu,ij_begin,ij_end)
	165	c call minmaxq(zq,qmin,qmax,'apres vlx2 ')
	166
	167
	168	ijb=ij_begin
	169	ije=ij_end
	170
	171	DO l=1,llm
	172	DO ij=ijb,ije
	173	q(ij,l)=zq(ij,l)
	174	ENDDO
	175	ENDDO
	176
	177
	178	DO l=1,llm
	179	DO ij=ijb,ije-iip1+1,iip1
	180	q(ij+iim,l)=q(ij,l)
	181	ENDDO
	182	ENDDO
	183
	184	call WaitSendRequest(MyRequest1)
	185	call WaitSendRequest(MyRequest2)
	186
	187	RETURN
	188	END
	189
	190
	191	SUBROUTINE vlx_p(q,pente_max,masse,u_m,ijb_x,ije_x)
	192
	193	c Auteurs: P.Le Van, F.Hourdin, F.Forget
	194	c
	195	c ********************************************************************
	196	c Shema d'advection " pseudo amont " .
	197	c ********************************************************************
	198	c nq,iq,q,pbaru,pbarv,w sont des arguments d'entree pour le s-pg ....
	199	c
	200	c
	201	c --------------------------------------------------------------------
[1823]	202	USE parallel_lmdz
[630]	203	IMPLICIT NONE
	204	c
	205	#include "dimensions.h"
	206	#include "paramet.h"
	207	c
	208	c
	209	c Arguments:
	210	c ----------
	211	REAL masse(ip1jmp1,llm),pente_max
	212	REAL u_m( ip1jmp1,llm ),pbarv( iip1,jjm,llm)
	213	REAL q(ip1jmp1,llm)
	214	REAL w(ip1jmp1,llm)
	215	c
	216	c Local
	217	c ---------
	218	c
	219	INTEGER ij,l,j,i,iju,ijq,indu(ip1jmp1),niju
	220	INTEGER n0,iadvplus(ip1jmp1,llm),nl(llm)
	221	c
	222	REAL new_m,zu_m,zdum(ip1jmp1,llm)
	223	REAL sigu(ip1jmp1),dxq(ip1jmp1,llm),dxqu(ip1jmp1)
	224	REAL zz(ip1jmp1)
	225	REAL adxqu(ip1jmp1),dxqmax(ip1jmp1,llm)
	226	REAL u_mq(ip1jmp1,llm)
	227
[764]	228	Logical extremum
[630]	229
	230	REAL SSUM
	231	EXTERNAL SSUM
	232
	233	REAL z1,z2,z3
	234
	235	INTEGER ijb,ije,ijb_x,ije_x
	236
	237	c calcul de la pente a droite et a gauche de la maille
	238
	239	ijb=ijb_x
	240	ije=ije_x
	241
	242	if (pole_nord.and.ijb==1) ijb=ijb+iip1
	243	if (pole_sud.and.ije==ip1jmp1) ije=ije-iip1
	244
	245	IF (pente_max.gt.-1.e-5) THEN
	246	c IF (pente_max.gt.10) THEN
	247
	248	c calcul des pentes avec limitation, Van Leer scheme I:
	249	c -----------------------------------------------------
	250
	251	c calcul de la pente aux points u
[764]	252	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	253	DO l = 1, llm
	254
	255	DO ij=ijb,ije-1
	256	dxqu(ij)=q(ij+1,l)-q(ij,l)
	257	c IF(u_m(ij,l).lt.0.) stop'limx n admet pas les U<0'
	258	c sigu(ij)=u_m(ij,l)/masse(ij,l)
	259	ENDDO
	260	DO ij=ijb+iip1-1,ije,iip1
	261	dxqu(ij)=dxqu(ij-iim)
	262	c sigu(ij)=sigu(ij-iim)
	263	ENDDO
	264
	265	DO ij=ijb,ije
	266	adxqu(ij)=abs(dxqu(ij))
	267	ENDDO
	268
	269	c calcul de la pente maximum dans la maille en valeur absolue
	270
	271	DO ij=ijb+1,ije
	272	dxqmax(ij,l)=pente_max*
	273	, min(adxqu(ij-1),adxqu(ij))
	274	c limitation subtile
	275	c , min(adxqu(ij-1)/sigu(ij-1),adxqu(ij)/(1.-sigu(ij)))
	276
	277
	278	ENDDO
	279
	280	DO ij=ijb+iip1-1,ije,iip1
	281	dxqmax(ij-iim,l)=dxqmax(ij,l)
	282	ENDDO
	283
	284	DO ij=ijb+1,ije
	285	#ifdef CRAY
	286	dxq(ij,l)=
	287	, cvmgp(dxqu(ij-1)+dxqu(ij),0.,dxqu(ij-1)*dxqu(ij))
	288	#else
	289	IF(dxqu(ij-1)*dxqu(ij).gt.0) THEN
	290	dxq(ij,l)=dxqu(ij-1)+dxqu(ij)
	291	ELSE
	292	c extremum local
	293	dxq(ij,l)=0.
	294	ENDIF
	295	#endif
	296	dxq(ij,l)=0.5*dxq(ij,l)
	297	dxq(ij,l)=
	298	, sign(min(abs(dxq(ij,l)),dxqmax(ij,l)),dxq(ij,l))
	299	ENDDO
	300
	301	ENDDO ! l=1,llm
[764]	302	c$OMP END DO NOWAIT
[630]	303	c print*,'Ok calcul des pentes'
	304
	305	ELSE ! (pente_max.lt.-1.e-5)
	306
	307	c Pentes produits:
	308	c ----------------
[764]	309	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	310	DO l = 1, llm
	311	DO ij=ijb,ije-1
	312	dxqu(ij)=q(ij+1,l)-q(ij,l)
	313	ENDDO
	314	DO ij=ijb+iip1-1,ije,iip1
	315	dxqu(ij)=dxqu(ij-iim)
	316	ENDDO
	317
	318	DO ij=ijb+1,ije
	319	zz(ij)=dxqu(ij-1)*dxqu(ij)
	320	zz(ij)=zz(ij)+zz(ij)
	321	IF(zz(ij).gt.0) THEN
	322	dxq(ij,l)=zz(ij)/(dxqu(ij-1)+dxqu(ij))
	323	ELSE
	324	c extremum local
	325	dxq(ij,l)=0.
	326	ENDIF
	327	ENDDO
	328
	329	ENDDO
[764]	330	c$OMP END DO NOWAIT
[630]	331	ENDIF ! (pente_max.lt.-1.e-5)
	332
	333	c bouclage de la pente en iip1:
	334	c -----------------------------
[764]	335	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	336	DO l=1,llm
	337	DO ij=ijb+iip1-1,ije,iip1
	338	dxq(ij-iim,l)=dxq(ij,l)
	339	ENDDO
[764]	340	DO ij=ijb,ije
[630]	341	iadvplus(ij,l)=0
	342	ENDDO
	343
	344	ENDDO
[764]	345	c$OMP END DO NOWAIT
[630]	346	c print*,'Bouclage en iip1'
	347
	348	c calcul des flux a gauche et a droite
	349
	350	#ifdef CRAY
[764]	351	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	352	DO l=1,llm
	353	DO ij=ijb,ije-1
	354	zdum(ij,l)=cvmgp(1.-u_m(ij,l)/masse(ij,l),
	355	, 1.+u_m(ij,l)/masse(ij+1,l),
	356	, u_m(ij,l))
	357	zdum(ij,l)=0.5*zdum(ij,l)
	358	u_mq(ij,l)=cvmgp(
	359	, q(ij,l)+zdum(ij,l)*dxq(ij,l),
	360	, q(ij+1,l)-zdum(ij,l)*dxq(ij+1,l),
	361	, u_m(ij,l))
	362	u_mq(ij,l)=u_m(ij,l)*u_mq(ij,l)
	363	ENDDO
	364	ENDDO
[764]	365	c$OMP END DO NOWAIT
[630]	366	#else
	367	c on cumule le flux correspondant a toutes les mailles dont la masse
	368	c au travers de la paroi pENDant le pas de temps.
	369	c print*,'Cumule ....'
[764]	370	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	371	DO l=1,llm
	372	DO ij=ijb,ije-1
	373	c print*,'masse(',ij,')=',masse(ij,l)
	374	IF (u_m(ij,l).gt.0.) THEN
	375	zdum(ij,l)=1.-u_m(ij,l)/masse(ij,l)
	376	u_mq(ij,l)=u_m(ij,l)(q(ij,l)+0.5zdum(ij,l)*dxq(ij,l))
	377	ELSE
	378	zdum(ij,l)=1.+u_m(ij,l)/masse(ij+1,l)
	379	u_mq(ij,l)=u_m(ij,l)(q(ij+1,l)-0.5zdum(ij,l)*dxq(ij+1,l))
	380	ENDIF
	381	ENDDO
	382	ENDDO
[764]	383	c$OMP END DO NOWAIT
[630]	384	#endif
	385	c stop
	386
	387	c go to 9999
	388	c detection des points ou on advecte plus que la masse de la
	389	c maille
[764]	390	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	391	DO l=1,llm
	392	DO ij=ijb,ije-1
	393	IF(zdum(ij,l).lt.0) THEN
	394	iadvplus(ij,l)=1
	395	u_mq(ij,l)=0.
	396	ENDIF
	397	ENDDO
	398	ENDDO
[764]	399	c$OMP END DO NOWAIT
[630]	400	c print*,'Ok test 1'
[764]	401	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	402	DO l=1,llm
	403	DO ij=ijb+iip1-1,ije,iip1
	404	iadvplus(ij,l)=iadvplus(ij-iim,l)
	405	ENDDO
	406	ENDDO
[764]	407	c$OMP END DO NOWAIT
[630]	408	c print*,'Ok test 2'
	409
	410
	411	c traitement special pour le cas ou on advecte en longitude plus que le
	412	c contenu de la maille.
	413	c cette partie est mal vectorisee.
	414
	415	c calcul du nombre de maille sur lequel on advecte plus que la maille.
	416
	417	n0=0
[764]	418	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	419	DO l=1,llm
	420	nl(l)=0
	421	DO ij=ijb,ije
	422	nl(l)=nl(l)+iadvplus(ij,l)
	423	ENDDO
	424	n0=n0+nl(l)
	425	ENDDO
[764]	426	c$OMP END DO NOWAIT
[630]	427	cym IF(n0.gt.1) THEN
[764]	428	cym IF(n0.gt.0) THEN
[630]	429
	430	c PRINT*,'Nombre de points pour lesquels on advect plus que le'
	431	c & ,'contenu de la maille : ',n0
[764]	432	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	433	DO l=1,llm
	434	IF(nl(l).gt.0) THEN
	435	iju=0
	436	c indicage des mailles concernees par le traitement special
	437	DO ij=ijb,ije
	438	IF(iadvplus(ij,l).eq.1.and.mod(ij,iip1).ne.0) THEN
	439	iju=iju+1
	440	indu(iju)=ij
	441	ENDIF
	442	ENDDO
	443	niju=iju
	444	c PRINT*,'niju,nl',niju,nl(l)
	445
	446	c traitement des mailles
	447	DO iju=1,niju
	448	ij=indu(iju)
	449	j=(ij-1)/iip1+1
	450	zu_m=u_m(ij,l)
	451	u_mq(ij,l)=0.
	452	IF(zu_m.gt.0.) THEN
	453	ijq=ij
	454	i=ijq-(j-1)*iip1
	455	c accumulation pour les mailles completements advectees
	456	do while(zu_m.gt.masse(ijq,l))
	457	u_mq(ij,l)=u_mq(ij,l)+q(ijq,l)*masse(ijq,l)
	458	zu_m=zu_m-masse(ijq,l)
	459	i=mod(i-2+iim,iim)+1
	460	ijq=(j-1)*iip1+i
	461	ENDDO
	462	c ajout de la maille non completement advectee
	463	u_mq(ij,l)=u_mq(ij,l)+zu_m*
	464	& (q(ijq,l)+0.5(1.-zu_m/masse(ijq,l))dxq(ijq,l))
	465	ELSE
	466	ijq=ij+1
	467	i=ijq-(j-1)*iip1
	468	c accumulation pour les mailles completements advectees
	469	do while(-zu_m.gt.masse(ijq,l))
	470	u_mq(ij,l)=u_mq(ij,l)-q(ijq,l)*masse(ijq,l)
	471	zu_m=zu_m+masse(ijq,l)
	472	i=mod(i,iim)+1
	473	ijq=(j-1)*iip1+i
	474	ENDDO
	475	c ajout de la maille non completement advectee
	476	u_mq(ij,l)=u_mq(ij,l)+zu_m*(q(ijq,l)-
	477	& 0.5(1.+zu_m/masse(ijq,l))dxq(ijq,l))
	478	ENDIF
	479	ENDDO
	480	ENDIF
	481	ENDDO
[764]	482	c$OMP END DO NOWAIT
	483	cym ENDIF ! n0.gt.0
[630]	484	9999 continue
	485
	486
	487	c bouclage en latitude
	488	c print*,'Avant bouclage en latitude'
[764]	489	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	490	DO l=1,llm
	491	DO ij=ijb+iip1-1,ije,iip1
	492	u_mq(ij,l)=u_mq(ij-iim,l)
	493	ENDDO
	494	ENDDO
[764]	495	c$OMP END DO NOWAIT
[630]	496
	497	c calcul des tENDances
[764]	498	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	499	DO l=1,llm
	500	DO ij=ijb+1,ije
	501	new_m=masse(ij,l)+u_m(ij-1,l)-u_m(ij,l)
	502	q(ij,l)=(q(ij,l)*masse(ij,l)+
	503	& u_mq(ij-1,l)-u_mq(ij,l))
	504	& /new_m
	505	masse(ij,l)=new_m
	506	ENDDO
	507	c ModIF Fred 22 03 96 correction d'un bug (les scopy ci-dessous)
	508	DO ij=ijb+iip1-1,ije,iip1
	509	q(ij-iim,l)=q(ij,l)
	510	masse(ij-iim,l)=masse(ij,l)
	511	ENDDO
	512	ENDDO
[764]	513	c$OMP END DO NOWAIT
[630]	514	c CALL SCOPY((jjm-1)*llm,q(iip1+iip1,1),iip1,q(iip2,1),iip1)
	515	c CALL SCOPY((jjm-1)*llm,masse(iip1+iip1,1),iip1,masse(iip2,1),iip1)
	516
	517
	518	RETURN
	519	END
	520
	521
	522	SUBROUTINE vly_p(q,pente_max,masse,masse_adv_v)
	523	c
	524	c Auteurs: P.Le Van, F.Hourdin, F.Forget
	525	c
	526	c ********************************************************************
	527	c Shema d'advection " pseudo amont " .
	528	c ********************************************************************
	529	c q,masse_adv_v,w sont des arguments d'entree pour le s-pg ....
	530	c dq sont des arguments de sortie pour le s-pg ....
	531	c
	532	c
	533	c --------------------------------------------------------------------
[1823]	534	USE parallel_lmdz
[2597]	535	USE comconst_mod, ONLY: pi
[630]	536	IMPLICIT NONE
	537	c
	538	#include "dimensions.h"
	539	#include "paramet.h"
	540	#include "comgeom.h"
	541	c
	542	c
	543	c Arguments:
	544	c ----------
	545	REAL masse(ip1jmp1,llm),pente_max
	546	REAL masse_adv_v( ip1jm,llm)
	547	REAL q(ip1jmp1,llm), dq( ip1jmp1,llm)
	548	c
	549	c Local
	550	c ---------
	551	c
	552	INTEGER i,ij,l
	553	c
	554	REAL airej2,airejjm,airescb(iim),airesch(iim)
	555	REAL dyq(ip1jmp1,llm),dyqv(ip1jm),zdvm(ip1jmp1,llm)
	556	REAL adyqv(ip1jm),dyqmax(ip1jmp1)
	557	REAL qbyv(ip1jm,llm)
	558
[1520]	559	REAL qpns,qpsn,appn,apps,dyn1,dys1,dyn2,dys2,newmasse,fn,fs
[630]	560	c REAL newq,oldmasse
	561	Logical extremum,first,testcpu
	562	REAL temps0,temps1,temps2,temps3,temps4,temps5,second
	563	SAVE temps0,temps1,temps2,temps3,temps4,temps5
[764]	564	c$OMP THREADPRIVATE(temps0,temps1,temps2,temps3,temps4,temps5)
[630]	565	SAVE first,testcpu
[764]	566	c$OMP THREADPRIVATE(first,testcpu)
[630]	567
	568	REAL convpn,convps,convmpn,convmps
	569	real massepn,masseps,qpn,qps
	570	REAL sinlon(iip1),sinlondlon(iip1)
	571	REAL coslon(iip1),coslondlon(iip1)
	572	SAVE sinlon,coslon,sinlondlon,coslondlon
[764]	573	c$OMP THREADPRIVATE(sinlon,coslon,sinlondlon,coslondlon)
[630]	574	SAVE airej2,airejjm
[764]	575	c$OMP THREADPRIVATE(airej2,airejjm)
[630]	576	c
	577	c
	578	REAL SSUM
	579	EXTERNAL SSUM
	580
	581	DATA first,testcpu/.true.,.false./
	582	DATA temps0,temps1,temps2,temps3,temps4,temps5/0.,0.,0.,0.,0.,0./
	583	INTEGER ijb,ije
	584
	585	IF(first) THEN
	586	c PRINT*,'Shema Amont nouveau appele dans Vanleer '
	587	first=.false.
	588	do i=2,iip1
	589	coslon(i)=cos(rlonv(i))
	590	sinlon(i)=sin(rlonv(i))
	591	coslondlon(i)=coslon(i)*(rlonu(i)-rlonu(i-1))/pi
	592	sinlondlon(i)=sinlon(i)*(rlonu(i)-rlonu(i-1))/pi
	593	ENDDO
	594	coslon(1)=coslon(iip1)
	595	coslondlon(1)=coslondlon(iip1)
	596	sinlon(1)=sinlon(iip1)
	597	sinlondlon(1)=sinlondlon(iip1)
	598	airej2 = SSUM( iim, aire(iip2), 1 )
	599	airejjm= SSUM( iim, aire(ip1jm -iim), 1 )
	600	ENDIF
	601
	602	c
	603	c PRINT*,'CALCUL EN LATITUDE'
	604
[764]	605	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	606	DO l = 1, llm
	607	c
	608	c --------------------------------
	609	c CALCUL EN LATITUDE
	610	c --------------------------------
	611
	612	c On commence par calculer la valeur du traceur moyenne sur le premier cercle
	613	c de latitude autour du pole (qpns pour le pole nord et qpsn pour
	614	c le pole nord) qui sera utilisee pour evaluer les pentes au pole.
	615
	616	if (pole_nord) then
	617	DO i = 1, iim
	618	airescb(i) = aire(i+ iip1) * q(i+ iip1,l)
	619	ENDDO
	620	qpns = SSUM( iim, airescb ,1 ) / airej2
	621	endif
	622
	623	if (pole_sud) then
	624	DO i = 1, iim
	625	airesch(i) = aire(i+ ip1jm- iip1) * q(i+ ip1jm- iip1,l)
	626	ENDDO
	627	qpsn = SSUM( iim, airesch ,1 ) / airejjm
	628	endif
	629
	630
	631
	632	c calcul des pentes aux points v
	633
	634	ijb=ij_begin-2*iip1
	635	ije=ij_end+iip1
	636	if (pole_nord) ijb=ij_begin
	637	if (pole_sud) ije=ij_end-iip1
	638
	639	DO ij=ijb,ije
	640	dyqv(ij)=q(ij,l)-q(ij+iip1,l)
	641	adyqv(ij)=abs(dyqv(ij))
	642	ENDDO
	643
	644	c calcul des pentes aux points scalaires
	645	ijb=ij_begin-iip1
	646	ije=ij_end+iip1
	647	if (pole_nord) ijb=ij_begin+iip1
	648	if (pole_sud) ije=ij_end-iip1
	649
	650	DO ij=ijb,ije
	651	dyq(ij,l)=.5*(dyqv(ij-iip1)+dyqv(ij))
	652	dyqmax(ij)=min(adyqv(ij-iip1),adyqv(ij))
	653	dyqmax(ij)=pente_max*dyqmax(ij)
	654	ENDDO
	655
	656	c calcul des pentes aux poles
	657	IF (pole_nord) THEN
	658	DO ij=1,iip1
	659	dyq(ij,l)=qpns-q(ij+iip1,l)
	660	ENDDO
	661
	662	dyn1=0.
	663	dyn2=0.
	664	DO ij=1,iim
	665	dyn1=dyn1+sinlondlon(ij)*dyq(ij,l)
	666	dyn2=dyn2+coslondlon(ij)*dyq(ij,l)
	667	ENDDO
	668	DO ij=1,iip1
	669	dyq(ij,l)=dyn1sinlon(ij)+dyn2coslon(ij)
	670	ENDDO
	671
	672	DO ij=1,iip1
	673	dyq(ij,l)=0.
	674	ENDDO
	675	c ym tout cela ne sert pas a grand chose
	676	ENDIF
	677
	678	IF (pole_sud) THEN
	679
	680	DO ij=1,iip1
	681	dyq(ip1jm+ij,l)=q(ip1jm+ij-iip1,l)-qpsn
	682	ENDDO
	683
	684	dys1=0.
	685	dys2=0.
	686
	687	DO ij=1,iim
	688	dys1=dys1+sinlondlon(ij)*dyq(ip1jm+ij,l)
	689	dys2=dys2+coslondlon(ij)*dyq(ip1jm+ij,l)
	690	ENDDO
	691
	692	DO ij=1,iip1
	693	dyq(ip1jm+ij,l)=dys1sinlon(ij)+dys2coslon(ij)
	694	ENDDO
	695
	696	DO ij=1,iip1
	697	dyq(ip1jm+ij,l)=0.
	698	ENDDO
	699	c ym tout cela ne sert pas a grand chose
	700	ENDIF
	701
	702	c filtrage de la derivee
	703
	704	c calcul des pentes limites aux poles
	705	c ym partie inutile
	706	c goto 8888
	707	c fn=1.
	708	c fs=1.
	709	c DO ij=1,iim
	710	c IF(pente_max*adyqv(ij).lt.abs(dyq(ij,l))) THEN
	711	c fn=min(pente_max*adyqv(ij)/abs(dyq(ij,l)),fn)
	712	c ENDIF
	713	c IF(pente_max*adyqv(ij+ip1jm-iip1).lt.abs(dyq(ij+ip1jm,l))) THEN
	714	c fs=min(pente_max*adyqv(ij+ip1jm-iip1)/abs(dyq(ij+ip1jm,l)),fs)
	715	c ENDIF
	716	c ENDDO
	717	c DO ij=1,iip1
	718	c dyq(ij,l)=fn*dyq(ij,l)
	719	c dyq(ip1jm+ij,l)=fs*dyq(ip1jm+ij,l)
	720	c ENDDO
	721	c 8888 continue
	722
	723
	724	CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
	725	C En memoire de dIFferents tests sur la
	726	C limitation des pentes aux poles.
	727	CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
	728	C PRINT*,dyq(1)
	729	C PRINT*,dyqv(iip1+1)
[1520]	730	C appn=abs(dyq(1)/dyqv(iip1+1))
[630]	731	C PRINT*,dyq(ip1jm+1)
	732	C PRINT*,dyqv(ip1jm-iip1+1)
[1520]	733	C apps=abs(dyq(ip1jm+1)/dyqv(ip1jm-iip1+1))
[630]	734	C DO ij=2,iim
[1520]	735	C appn=amax1(abs(dyq(ij)/dyqv(ij)),appn)
	736	C apps=amax1(abs(dyq(ip1jm+ij)/dyqv(ip1jm-iip1+ij)),apps)
[630]	737	C ENDDO
[1520]	738	C appn=min(pente_max/appn,1.)
	739	C apps=min(pente_max/apps,1.)
[630]	740	C
	741	C
	742	C cas ou on a un extremum au pole
	743	C
	744	C IF(dyqv(ismin(iim,dyqv,1))*dyqv(ismax(iim,dyqv,1)).le.0.)
[1520]	745	C & appn=0.
[630]	746	C IF(dyqv(ismax(iim,dyqv(ip1jm-iip1+1),1)+ip1jm-iip1+1)*
	747	C & dyqv(ismin(iim,dyqv(ip1jm-iip1+1),1)+ip1jm-iip1+1).le.0.)
[1520]	748	C & apps=0.
[630]	749	C
	750	C limitation des pentes aux poles
	751	C DO ij=1,iip1
[1520]	752	C dyq(ij)=appn*dyq(ij)
	753	C dyq(ip1jm+ij)=apps*dyq(ip1jm+ij)
[630]	754	C ENDDO
	755	C
	756	C test
	757	C DO ij=1,iip1
	758	C dyq(iip1+ij)=0.
	759	C dyq(ip1jm+ij-iip1)=0.
	760	C ENDDO
	761	C DO ij=1,ip1jmp1
	762	C dyq(ij)=dyq(ij)*cos(rlatu((ij-1)/iip1+1))
	763	C ENDDO
	764	C
	765	C changement 10 07 96
	766	C IF(dyqv(ismin(iim,dyqv,1))*dyqv(ismax(iim,dyqv,1)).le.0.)
	767	C & THEN
	768	C DO ij=1,iip1
	769	C dyqmax(ij)=0.
	770	C ENDDO
	771	C ELSE
	772	C DO ij=1,iip1
	773	C dyqmax(ij)=pente_max*abs(dyqv(ij))
	774	C ENDDO
	775	C ENDIF
	776	C
	777	C IF(dyqv(ismax(iim,dyqv(ip1jm-iip1+1),1)+ip1jm-iip1+1)*
	778	C & dyqv(ismin(iim,dyqv(ip1jm-iip1+1),1)+ip1jm-iip1+1).le.0.)
	779	C &THEN
	780	C DO ij=ip1jm+1,ip1jmp1
	781	C dyqmax(ij)=0.
	782	C ENDDO
	783	C ELSE
	784	C DO ij=ip1jm+1,ip1jmp1
	785	C dyqmax(ij)=pente_max*abs(dyqv(ij-iip1))
	786	C ENDDO
	787	C ENDIF
	788	C fin changement 10 07 96
	789	CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
	790
	791	c calcul des pentes limitees
	792	ijb=ij_begin-iip1
	793	ije=ij_end+iip1
	794	if (pole_nord) ijb=ij_begin+iip1
	795	if (pole_sud) ije=ij_end-iip1
	796
	797	DO ij=ijb,ije
	798	IF(dyqv(ij)*dyqv(ij-iip1).gt.0.) THEN
	799	dyq(ij,l)=sign(min(abs(dyq(ij,l)),dyqmax(ij)),dyq(ij,l))
	800	ELSE
	801	dyq(ij,l)=0.
	802	ENDIF
	803	ENDDO
	804
	805	ENDDO
[764]	806	c$OMP END DO NOWAIT
[630]	807
	808	ijb=ij_begin-iip1
	809	ije=ij_end
	810	if (pole_nord) ijb=ij_begin
	811	if (pole_sud) ije=ij_end-iip1
	812
[764]	813	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	814	DO l=1,llm
	815	DO ij=ijb,ije
	816	IF(masse_adv_v(ij,l).gt.0) THEN
	817	qbyv(ij,l)=q(ij+iip1,l)+dyq(ij+iip1,l)*
	818	, 0.5*(1.-masse_adv_v(ij,l)/masse(ij+iip1,l))
	819	ELSE
	820	qbyv(ij,l)=q(ij,l)-dyq(ij,l)*
	821	, 0.5*(1.+masse_adv_v(ij,l)/masse(ij,l))
	822	ENDIF
	823	qbyv(ij,l)=masse_adv_v(ij,l)*qbyv(ij,l)
	824	ENDDO
	825	ENDDO
[764]	826	c$OMP END DO NOWAIT
[630]	827
	828	ijb=ij_begin
	829	ije=ij_end
	830	if (pole_nord) ijb=ij_begin+iip1
	831	if (pole_sud) ije=ij_end-iip1
	832
[764]	833	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	834	DO l=1,llm
	835	DO ij=ijb,ije
	836	newmasse=masse(ij,l)
	837	& +masse_adv_v(ij,l)-masse_adv_v(ij-iip1,l)
	838
	839	q(ij,l)=(q(ij,l)*masse(ij,l)+qbyv(ij,l)-qbyv(ij-iip1,l))
	840	& /newmasse
	841	masse(ij,l)=newmasse
	842	ENDDO
	843	c.-. ancienne version
	844	c convpn=SSUM(iim,qbyv(1,l),1)/apoln
	845	c convmpn=ssum(iim,masse_adv_v(1,l),1)/apoln
	846	if (pole_nord) then
	847	convpn=SSUM(iim,qbyv(1,l),1)
	848	convmpn=ssum(iim,masse_adv_v(1,l),1)
	849	massepn=ssum(iim,masse(1,l),1)
	850	qpn=0.
	851	do ij=1,iim
	852	qpn=qpn+masse(ij,l)*q(ij,l)
	853	enddo
	854	qpn=(qpn+convpn)/(massepn+convmpn)
	855	do ij=1,iip1
	856	q(ij,l)=qpn
	857	enddo
	858	endif
	859
	860	c convps=-SSUM(iim,qbyv(ip1jm-iim,l),1)/apols
	861	c convmps=-ssum(iim,masse_adv_v(ip1jm-iim,l),1)/apols
	862
	863	if (pole_sud) then
	864
	865	convps=-SSUM(iim,qbyv(ip1jm-iim,l),1)
	866	convmps=-ssum(iim,masse_adv_v(ip1jm-iim,l),1)
	867	masseps=ssum(iim, masse(ip1jm+1,l),1)
	868	qps=0.
	869	do ij = ip1jm+1,ip1jmp1-1
	870	qps=qps+masse(ij,l)*q(ij,l)
	871	enddo
	872	qps=(qps+convps)/(masseps+convmps)
	873	do ij=ip1jm+1,ip1jmp1
	874	q(ij,l)=qps
	875	enddo
	876	endif
	877	c.-. fin ancienne version
	878
	879	c._. nouvelle version
	880	c convpn=SSUM(iim,qbyv(1,l),1)
	881	c convmpn=ssum(iim,masse_adv_v(1,l),1)
	882	c oldmasse=ssum(iim,masse(1,l),1)
	883	c newmasse=oldmasse+convmpn
	884	c newq=(q(1,l)*oldmasse+convpn)/newmasse
	885	c newmasse=newmasse/apoln
	886	c DO ij = 1,iip1
	887	c q(ij,l)=newq
	888	c masse(ij,l)=newmasse*aire(ij)
	889	c ENDDO
	890	c convps=-SSUM(iim,qbyv(ip1jm-iim,l),1)
	891	c convmps=-ssum(iim,masse_adv_v(ip1jm-iim,l),1)
	892	c oldmasse=ssum(iim,masse(ip1jm-iim,l),1)
	893	c newmasse=oldmasse+convmps
	894	c newq=(q(ip1jmp1,l)*oldmasse+convps)/newmasse
	895	c newmasse=newmasse/apols
	896	c DO ij = ip1jm+1,ip1jmp1
	897	c q(ij,l)=newq
	898	c masse(ij,l)=newmasse*aire(ij)
	899	c ENDDO
	900	c._. fin nouvelle version
	901	ENDDO
[764]	902	c$OMP END DO NOWAIT
[630]	903
	904	RETURN
	905	END
	906
	907
	908
	909	SUBROUTINE vlz_p(q,pente_max,masse,w,ijb_x,ije_x)
	910	c
	911	c Auteurs: P.Le Van, F.Hourdin, F.Forget
	912	c
	913	c ********************************************************************
	914	c Shema d'advection " pseudo amont " .
	915	c ********************************************************************
	916	c q,pbaru,pbarv,w sont des arguments d'entree pour le s-pg ....
	917	c dq sont des arguments de sortie pour le s-pg ....
	918	c
	919	c
	920	c --------------------------------------------------------------------
[1823]	921	USE parallel_lmdz
[630]	922	IMPLICIT NONE
	923	c
	924	#include "dimensions.h"
	925	#include "paramet.h"
	926	c
	927	c
	928	c Arguments:
	929	c ----------
	930	REAL masse(ip1jmp1,llm),pente_max
	931	REAL q(ip1jmp1,llm)
	932	REAL w(ip1jmp1,llm+1)
	933	c
	934	c Local
	935	c ---------
	936	c
	937	INTEGER i,ij,l,j,ii
	938	c
[764]	939	REAL,SAVE :: wq(ip1jmp1,llm+1)
	940	REAL newmasse
[630]	941
[764]	942	REAL,SAVE :: dzq(ip1jmp1,llm),dzqw(ip1jmp1,llm),adzqw(ip1jmp1,llm)
	943	REAL dzqmax
[630]	944	REAL sigw
	945
	946	LOGICAL testcpu
	947	SAVE testcpu
[764]	948	c$OMP THREADPRIVATE(testcpu)
[630]	949	REAL temps0,temps1,temps2,temps3,temps4,temps5,second
	950	SAVE temps0,temps1,temps2,temps3,temps4,temps5
[764]	951	c$OMP THREADPRIVATE(temps0,temps1,temps2,temps3,temps4,temps5)
	952
[630]	953	REAL SSUM
[764]	954	EXTERNAL SSUM
[630]	955
	956	DATA testcpu/.false./
	957	DATA temps0,temps1,temps2,temps3,temps4,temps5/0.,0.,0.,0.,0.,0./
	958	INTEGER ijb,ije,ijb_x,ije_x
	959	c On oriente tout dans le sens de la pression c'est a dire dans le
	960	c sens de W
	961
	962	#ifdef BIDON
	963	IF(testcpu) THEN
	964	temps0=second(0.)
	965	ENDIF
	966	#endif
	967
	968	ijb=ijb_x
	969	ije=ije_x
[764]	970
	971	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	972	DO l=2,llm
	973	DO ij=ijb,ije
	974	dzqw(ij,l)=q(ij,l-1)-q(ij,l)
	975	adzqw(ij,l)=abs(dzqw(ij,l))
	976	ENDDO
	977	ENDDO
[764]	978	c$OMP END DO
[630]	979
[764]	980	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	981	DO l=2,llm-1
	982	DO ij=ijb,ije
	983	#ifdef CRAY
	984	dzq(ij,l)=0.5*
	985	, cvmgp(dzqw(ij,l)+dzqw(ij,l+1),0.,dzqw(ij,l)*dzqw(ij,l+1))
	986	#else
	987	IF(dzqw(ij,l)*dzqw(ij,l+1).gt.0.) THEN
	988	dzq(ij,l)=0.5*(dzqw(ij,l)+dzqw(ij,l+1))
	989	ELSE
	990	dzq(ij,l)=0.
	991	ENDIF
	992	#endif
	993	dzqmax=pente_max*min(adzqw(ij,l),adzqw(ij,l+1))
	994	dzq(ij,l)=sign(min(abs(dzq(ij,l)),dzqmax),dzq(ij,l))
	995	ENDDO
	996	ENDDO
[764]	997	c$OMP END DO NOWAIT
[630]	998
[764]	999	c$OMP MASTER
[630]	1000	DO ij=ijb,ije
	1001	dzq(ij,1)=0.
	1002	dzq(ij,llm)=0.
	1003	ENDDO
[764]	1004	c$OMP END MASTER
	1005	c$OMP BARRIER
[630]	1006	#ifdef BIDON
	1007	IF(testcpu) THEN
	1008	temps1=temps1+second(0.)-temps0
	1009	ENDIF
	1010	#endif
	1011	c ---------------------------------------------------------------
	1012	c .... calcul des termes d'advection verticale .......
	1013	c ---------------------------------------------------------------
	1014
	1015	c calcul de - d( q * w )/ d(sigma) qu'on ajoute a dq pour calculer dq
	1016
[764]	1017	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	1018	DO l = 1,llm-1
	1019	do ij = ijb,ije
	1020	IF(w(ij,l+1).gt.0.) THEN
	1021	sigw=w(ij,l+1)/masse(ij,l+1)
	1022	wq(ij,l+1)=w(ij,l+1)(q(ij,l+1)+0.5(1.-sigw)*dzq(ij,l+1))
	1023	ELSE
	1024	sigw=w(ij,l+1)/masse(ij,l)
	1025	wq(ij,l+1)=w(ij,l+1)(q(ij,l)-0.5(1.+sigw)*dzq(ij,l))
	1026	ENDIF
	1027	ENDDO
	1028	ENDDO
[764]	1029	c$OMP END DO NOWAIT
[630]	1030
[764]	1031	c$OMP MASTER
[630]	1032	DO ij=ijb,ije
	1033	wq(ij,llm+1)=0.
	1034	wq(ij,1)=0.
	1035	ENDDO
[764]	1036	c$OMP END MASTER
	1037	c$OMP BARRIER
[630]	1038
[764]	1039	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	1040	DO l=1,llm
	1041	DO ij=ijb,ije
	1042	newmasse=masse(ij,l)+w(ij,l+1)-w(ij,l)
	1043	q(ij,l)=(q(ij,l)*masse(ij,l)+wq(ij,l+1)-wq(ij,l))
	1044	& /newmasse
	1045	masse(ij,l)=newmasse
	1046	ENDDO
	1047	ENDDO
[764]	1048	c$OMP END DO NOWAIT
[630]	1049
	1050
	1051	RETURN
	1052	END
	1053	c SUBROUTINE minmaxq(zq,qmin,qmax,comment)
	1054	c
	1055	c#include "dimensions.h"
	1056	c#include "paramet.h"
	1057
	1058	c CHARACTER() comment
	1059	c real qmin,qmax
	1060	c real zq(ip1jmp1,llm)
	1061
	1062	c INTEGER jadrs(ip1jmp1), jbad, k, i
	1063
	1064
	1065	c DO k = 1, llm
	1066	c jbad = 0
	1067	c DO i = 1, ip1jmp1
	1068	c IF (zq(i,k).GT.qmax .OR. zq(i,k).LT.qmin) THEN
	1069	c jbad = jbad + 1
	1070	c jadrs(jbad) = i
	1071	c ENDIF
	1072	c ENDDO
	1073	c IF (jbad.GT.0) THEN
	1074	c PRINT*, comment
	1075	c DO i = 1, jbad
	1076	cc PRINT*, "i,k,zq=", jadrs(i),k,zq(jadrs(i),k)
	1077	c ENDDO
	1078	c ENDIF
	1079	c ENDDO
	1080
	1081	c return
	1082	c end
	1083
	1084
	1085	subroutine minmaxq_p(zq,qmin,qmax,comment)
	1086
	1087	#include "dimensions.h"
	1088	#include "paramet.h"
	1089
	1090	character*20 comment
	1091	real qmin,qmax
	1092	real zq(ip1jmp1,llm)
	1093	real zzq(iip1,jjp1,llm)
	1094
	1095	integer imin,jmin,lmin,ijlmin
	1096	integer imax,jmax,lmax,ijlmax
	1097
	1098	integer ismin,ismax
	1099
	1100	#ifdef isminismax
	1101	call scopy (ip1jmp1*llm,zq,1,zzq,1)
	1102
	1103	ijlmin=ismin(ijp1llm,zq,1)
	1104	lmin=(ijlmin-1)/ip1jmp1+1
	1105	ijlmin=ijlmin-(lmin-1.)*ip1jmp1
	1106	jmin=(ijlmin-1)/iip1+1
	1107	imin=ijlmin-(jmin-1.)*iip1
	1108	zqmin=zq(ijlmin,lmin)
	1109
	1110	ijlmax=ismax(ijp1llm,zq,1)
	1111	lmax=(ijlmax-1)/ip1jmp1+1
	1112	ijlmax=ijlmax-(lmax-1.)*ip1jmp1
	1113	jmax=(ijlmax-1)/iip1+1
	1114	imax=ijlmax-(jmax-1.)*iip1
	1115	zqmax=zq(ijlmax,lmax)
	1116
	1117	if(zqmin.lt.qmin)
	1118	c s write(*,9999) comment,
	1119	s write(,) comment,
	1120	s imin,jmin,lmin,zqmin,zzq(imin,jmin,lmin)
	1121	if(zqmax.gt.qmax)
	1122	c s write(*,9999) comment,
	1123	s write(,) comment,
	1124	s imax,jmax,lmax,zqmax,zzq(imax,jmax,lmax)
	1125	#endif
	1126	return
	1127	9999 format(a20,' q(',i3,',',i2,',',i2,')=',e12.5,e12.5)
	1128	end
	1129
	1130
	1131

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

Download in other formats:

Original Format