Context Navigation

calfis_p.F @ 1422

Last change on this file since 1422 was 1422, checked in by milmd, 10 years ago
In GENERIC, MARS and COMMON models replace some include files by modules (usefull for decoupling physics with dynamics).
File size: 40.0 KB

Rev	Line
[1]	1	!
	2	! $Id: calfis_p.F 1407 2010-07-07 10:31:52Z fairhead $
	3	!
	4	C
	5	C
	6	SUBROUTINE calfis_p(lafin,
	7	$ jD_cur, jH_cur,
	8	$ pucov,
	9	$ pvcov,
	10	$ pteta,
	11	$ pq,
	12	$ pmasse,
	13	$ pps,
	14	$ pp,
	15	$ ppk,
	16	$ pphis,
	17	$ pphi,
	18	$ pducov,
	19	$ pdvcov,
	20	$ pdteta,
	21	$ pdq,
	22	$ flxw,
	23	$ pdufi,
	24	$ pdvfi,
	25	$ pdhfi,
	26	$ pdqfi,
	27	$ pdpsfi)
[8]	28	#ifdef CPP_PHYS
[776]	29	! Ehouarn: if using (parallelized) physics
[1]	30	USE dimphy
	31	USE mod_phys_lmdz_para, mpi_root_xx=>mpi_root
	32	USE mod_interface_dyn_phys
[1107]	33	! USE IOPHY
[1]	34	#endif
[1403]	35	#ifdef CPP_PARA
[1056]	36	USE parallel_lmdz, ONLY : omp_chunk, using_mpi, AllGather_Field
[1]	37	USE Write_Field
	38	Use Write_field_p
	39	USE Times
[1086]	40	USE cpdet_mod, only: tpot2t_p, t2tpot_p
[1056]	41	! used only for zonal averages
	42	USE moyzon_mod
[1403]	43	#endif
	44	USE infotrac, ONLY: nqtot, niadv, tname
	45	USE control_mod, ONLY: planet_type, nsplit_phys
[1422]	46	USE comvert_mod, ONLY: preff,presnivs
	47	USE comconst_mod, ONLY: daysec,dtvr,dtphys,kappa,cpp,g,rad,pi
	48	USE logic_mod, ONLY: moyzon_ch,moyzon_mu
[1056]	49
[1]	50	IMPLICIT NONE
	51	c=======================================================================
	52	c
	53	c 1. rearrangement des tableaux et transformation
	54	c variables dynamiques > variables physiques
	55	c 2. calcul des termes physiques
	56	c 3. retransformation des tendances physiques en tendances dynamiques
	57	c
	58	c remarques:
	59	c ----------
	60	c
	61	c - les vents sont donnes dans la physique par leurs composantes
	62	c naturelles.
	63	c - la variable thermodynamique de la physique est une variable
	64	c intensive : T
	65	c pour la dynamique on prend T * ( preff / p(l) ) **kappa
	66	c - les deux seules variables dependant de la geometrie necessaires
	67	c pour la physique sont la latitude pour le rayonnement et
	68	c l'aire de la maille quand on veut integrer une grandeur
	69	c horizontalement.
	70	c - les points de la physique sont les points scalaires de la
	71	c la dynamique; numerotation:
	72	c 1 pour le pole nord
	73	c (jjm-1)*iim pour l'interieur du domaine
	74	c ngridmx pour le pole sud
	75	c ---> ngridmx=2+(jjm-1)*iim
	76	c
	77	c Input :
	78	c -------
	79	c ecritphy frequence d'ecriture (en jours)de histphy
	80	c pucov covariant zonal velocity
	81	c pvcov covariant meridional velocity
	82	c pteta potential temperature
	83	c pps surface pressure
	84	c pmasse masse d'air dans chaque maille
	85	c pts surface temperature (K)
	86	c callrad clef d'appel au rayonnement
	87	c
	88	c Output :
	89	c --------
	90	c pdufi tendency for the natural zonal velocity (ms-1)
	91	c pdvfi tendency for the natural meridional velocity
[8]	92	c pdhfi tendency for the potential temperature (K/s)
[1]	93	c pdtsfi tendency for the surface temperature
	94	c
	95	c pdtrad radiative tendencies \ both input
	96	c pfluxrad radiative fluxes / and output
	97	c
	98	c=======================================================================
	99	c
	100	c-----------------------------------------------------------------------
	101	c
	102	c 0. Declarations :
	103	c ------------------
	104
	105	#include "dimensions.h"
	106	#include "paramet.h"
	107
	108	INTEGER ngridmx
	109	PARAMETER( ngridmx = 2+(jjm-1)*iim - 1/jjm )
	110
	111	#include "comgeom2.h"
	112	#include "iniprint.h"
	113	#ifdef CPP_MPI
	114	include 'mpif.h'
	115	#endif
	116	c Arguments :
	117	c -----------
[1189]	118	LOGICAL,INTENT(IN) :: lafin ! .true. for the very last call to physics
	119	REAL,INTENT(IN) :: jD_cur, jH_cur
	120	REAL,INTENT(IN) :: pvcov(iip1,jjm,llm) ! covariant meridional velocity
	121	REAL,INTENT(IN) :: pucov(iip1,jjp1,llm) ! covariant zonal velocity
	122	REAL,INTENT(IN) :: pteta(iip1,jjp1,llm) ! potential temperature
	123	REAL,INTENT(IN) :: pmasse(iip1,jjp1,llm) ! mass in each cell ! not used
	124	REAL,INTENT(IN) :: pq(iip1,jjp1,llm,nqtot) ! tracers
	125	REAL,INTENT(IN) :: pphis(iip1,jjp1) ! surface geopotential
	126	REAL,INTENT(IN) :: pphi(iip1,jjp1,llm) ! geopotential
[1056]	127
[1189]	128	REAL,INTENT(IN) :: pdvcov(iip1,jjm,llm) ! dynamical tendency on vcov
	129	REAL,INTENT(IN) :: pducov(iip1,jjp1,llm) ! dynamical tendency on ucov
	130	REAL,INTENT(IN) :: pdteta(iip1,jjp1,llm) ! dynamical tendency on teta
[8]	131	! commentaire SL: pdq ne sert que pour le calcul de pcvgq,
	132	! qui lui meme ne sert a rien dans la routine telle qu'elle est
	133	! ecrite, et que j'ai donc commente....
[1189]	134	REAL,INTENT(IN) :: pdq(iip1,jjp1,llm,nqtot) ! dynamical tendency on tracers
	135	! NB: pdq is only used to compute pcvgq which is in fact not used...
[1056]	136
[1189]	137	REAL,INTENT(IN) :: pps(iip1,jjp1) ! surface pressure (Pa)
	138	REAL,INTENT(IN) :: pp(iip1,jjp1,llmp1) ! pressure at mesh interfaces (Pa)
	139	REAL,INTENT(IN) :: ppk(iip1,jjp1,llm) ! Exner at mid-layer
[1312]	140	REAL,INTENT(IN) :: flxw(iip1,jjp1,llm) ! Vertical mass flux on lower mesh interfaces (kg/s) (on llm because flxw(:,:,llm+1)=0)
[1056]	141
[1189]	142	! tendencies (in */s) from the physics
	143	REAL,INTENT(OUT) :: pdvfi(iip1,jjm,llm) ! tendency on covariant meridional wind
	144	REAL,INTENT(OUT) :: pdufi(iip1,jjp1,llm) ! tendency on covariant zonal wind
	145	REAL,INTENT(OUT) :: pdhfi(iip1,jjp1,llm) ! tendency on potential temperature (K/s)
	146	REAL,INTENT(OUT) :: pdqfi(iip1,jjp1,llm,nqtot) ! tendency on tracers
	147	REAL,INTENT(OUT) :: pdpsfi(iip1,jjp1) ! tendency on surface pressure (Pa/s)
[1]	148
[1403]	149	#ifdef CPP_PARA
[8]	150	#ifdef CPP_PHYS
[1]	151	c Local variables :
	152	c -----------------
	153
	154	INTEGER i,j,l,ig0,ig,iq,iiq
	155	REAL,ALLOCATABLE,SAVE :: zpsrf(:)
	156	REAL,ALLOCATABLE,SAVE :: zplev(:,:),zplay(:,:)
	157	REAL,ALLOCATABLE,SAVE :: zphi(:,:),zphis(:)
[1056]	158
[1]	159	REAL,ALLOCATABLE,SAVE :: zufi(:,:), zvfi(:,:)
	160	REAL,ALLOCATABLE,SAVE :: ztfi(:,:),zqfi(:,:,:)
[8]	161	! ADAPTATION GCM POUR CP(T)
	162	REAL,ALLOCATABLE,SAVE :: zteta(:,:)
	163	REAL,ALLOCATABLE,SAVE :: zpk(:,:)
[1056]	164
[8]	165	! Ces calculs ne servent pas.
	166	! Si necessaire, decommenter ces variables et les calculs...
	167	! REAL,ALLOCATABLE,SAVE :: pcvgu(:,:), pcvgv(:,:)
	168	! REAL,ALLOCATABLE,SAVE :: pcvgt(:,:), pcvgq(:,:,:)
[1]	169	c
	170	REAL,ALLOCATABLE,SAVE :: zdufi(:,:),zdvfi(:,:)
	171	REAL,ALLOCATABLE,SAVE :: zdtfi(:,:),zdqfi(:,:,:)
	172	REAL,ALLOCATABLE,SAVE :: zdpsrf(:)
	173	REAL,SAVE,ALLOCATABLE :: flxwfi(:,:) ! Flux de masse verticale sur la grille physiq
	174
	175	REAL,ALLOCATABLE,SAVE :: zplev_omp(:,:)
	176	REAL,ALLOCATABLE,SAVE :: zplay_omp(:,:)
[892]	177	REAL,ALLOCATABLE,SAVE :: zpk_omp(:,:)
[1]	178	REAL,ALLOCATABLE,SAVE :: zphi_omp(:,:)
	179	REAL,ALLOCATABLE,SAVE :: zphis_omp(:)
	180	REAL,ALLOCATABLE,SAVE :: presnivs_omp(:)
	181	REAL,ALLOCATABLE,SAVE :: zufi_omp(:,:)
	182	REAL,ALLOCATABLE,SAVE :: zvfi_omp(:,:)
	183	REAL,ALLOCATABLE,SAVE :: ztfi_omp(:,:)
	184	REAL,ALLOCATABLE,SAVE :: zqfi_omp(:,:,:)
	185	REAL,ALLOCATABLE,SAVE :: zdufi_omp(:,:)
	186	REAL,ALLOCATABLE,SAVE :: zdvfi_omp(:,:)
	187	REAL,ALLOCATABLE,SAVE :: zdtfi_omp(:,:)
	188	REAL,ALLOCATABLE,SAVE :: zdqfi_omp(:,:,:)
	189	REAL,ALLOCATABLE,SAVE :: zdpsrf_omp(:)
	190	REAL,SAVE,ALLOCATABLE :: flxwfi_omp(:,:) ! Flux de masse verticale sur la grille physiq
	191
	192	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	193	! Introduction du splitting (FH)
	194	! Question pour Yann :
	195	! J'ai été surpris au début que les tableaux zufi_omp, zdufi_omp n'co soitent
	196	! en SAVE. Je crois comprendre que c'est parce que tu voulais qu'il
	197	! soit allocatable (plutot par exemple que de passer une dimension
	198	! dépendant du process en argument des routines) et que, du coup,
	199	! le SAVE évite d'avoir à refaire l'allocation à chaque appel.
	200	! Tu confirmes ?
	201	! J'ai suivi le même principe pour les zdufic_omp
	202	! Mais c'est surement bien que tu controles.
	203	!
	204
	205	REAL,ALLOCATABLE,SAVE :: zdufic_omp(:,:)
	206	REAL,ALLOCATABLE,SAVE :: zdvfic_omp(:,:)
	207	REAL,ALLOCATABLE,SAVE :: zdtfic_omp(:,:)
	208	REAL,ALLOCATABLE,SAVE :: zdqfic_omp(:,:,:)
	209	REAL jH_cur_split,zdt_split
	210	LOGICAL debut_split,lafin_split
	211	INTEGER isplit
	212	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
	213
[892]	214	c$OMP THREADPRIVATE(zplev_omp,zplay_omp,zpk_omp,zphi_omp,zphis_omp,
[1]	215	c$OMP+ presnivs_omp,zufi_omp,zvfi_omp,ztfi_omp,
	216	c$OMP+ zqfi_omp,zdufi_omp,zdvfi_omp,
	217	c$OMP+ zdtfi_omp,zdqfi_omp,zdpsrf_omp,flxwfi_omp,
	218	c$OMP+ zdufic_omp,zdvfic_omp,zdtfic_omp,zdqfic_omp)
	219
	220	LOGICAL,SAVE :: first_omp=.true.
	221	c$OMP THREADPRIVATE(first_omp)
	222
	223	REAL zsin(iim),zcos(iim),z1(iim)
	224	REAL zsinbis(iim),zcosbis(iim),z1bis(iim)
	225	REAL unskap, pksurcp
[1056]	226	save unskap
	227
[1]	228	REAL SSUM
	229
[1189]	230	LOGICAL,SAVE :: firstcal=.true., debut=.true.
[1]	231	c$OMP THREADPRIVATE(firstcal,debut)
	232
	233	REAL,SAVE,dimension(1:iim,1:llm):: du_send,du_recv,dv_send,dv_recv
	234	INTEGER :: ierr
	235	#ifdef CPP_MPI
	236	INTEGER,dimension(MPI_STATUS_SIZE,4) :: Status
	237	#else
	238	INTEGER,dimension(1,4) :: Status
	239	#endif
	240	INTEGER, dimension(4) :: Req
	241	REAL,ALLOCATABLE,SAVE:: zdufi2(:,:),zdvfi2(:,:)
	242	integer :: k,kstart,kend
	243	INTEGER :: offset
[841]	244
	245	LOGICAL tracerdyn ! for generic/mars physics call ; possibly to get rid of
[1056]	246
	247	! For Titan only right now:
	248	! to allow for 2D computation of microphys and chemistry
	249	LOGICAL,save :: flag_moyzon
[1076]	250	REAL,allocatable,save :: tmpvar(:,:)
	251	REAL,allocatable,save :: tmpvarp1(:,:)
	252	REAL,allocatable,save :: tmpvarbar(:)
	253	REAL,allocatable,save :: tmpvarbarp1(:)
[1126]	254	real :: zz1,zz2
[1056]	255
[1]	256	c-----------------------------------------------------------------------
[1056]	257
[1]	258	c 1. Initialisations :
	259	c --------------------
	260
[1056]	261
[1]	262	klon=klon_mpi
	263
	264	IF ( firstcal ) THEN
	265	debut = .TRUE.
	266	IF (ngridmx.NE.2+(jjm-1)*iim) THEN
	267	write(lunout,*) 'STOP dans calfis'
	268	write(lunout,*)
	269	& 'La dimension ngridmx doit etre egale a 2 + (jjm-1)*iim'
	270	write(lunout,*) ' ngridmx jjm iim '
	271	write(lunout,*) ngridmx,jjm,iim
	272	STOP
	273	ENDIF
[1056]	274
	275	unskap = 1./ kappa
	276
[1076]	277	c$OMP MASTER
[1056]	278	flag_moyzon = .false.
	279	if(moyzon_ch.or.moyzon_mu) then
	280	flag_moyzon = .true.
[1071]	281	allocate(tmpvar(iip1,llm))
	282	allocate(tmpvarp1(iip1,llmp1))
	283	allocate(tmpvarbar(llm))
	284	allocate(tmpvarbarp1(llmp1))
[1056]	285	endif
	286
[1]	287	ALLOCATE(zpsrf(klon))
	288	ALLOCATE(zplev(klon,llm+1),zplay(klon,llm))
	289	ALLOCATE(zphi(klon,llm),zphis(klon))
	290	ALLOCATE(zufi(klon,llm), zvfi(klon,llm))
	291	ALLOCATE(ztfi(klon,llm),zqfi(klon,llm,nqtot))
[8]	292	! ALLOCATE(pcvgu(klon,llm), pcvgv(klon,llm))
	293	! ALLOCATE(pcvgt(klon,llm), pcvgq(klon,llm,2))
[1]	294	ALLOCATE(zdufi(klon,llm),zdvfi(klon,llm))
	295	ALLOCATE(zdtfi(klon,llm),zdqfi(klon,llm,nqtot))
	296	ALLOCATE(zdpsrf(klon))
	297	ALLOCATE(zdufi2(klon+iim,llm),zdvfi2(klon+iim,llm))
	298	ALLOCATE(flxwfi(klon,llm))
[8]	299	! ADAPTATION GCM POUR CP(T)
[52]	300	ALLOCATE(zteta(klon,llm))
[8]	301	ALLOCATE(zpk(klon,llm))
[1056]	302
[1251]	303	! zonal means. horizontal dimension should be iip1
[1256]	304	if (flag_moyzon) call moyzon_init(klon,llm,nqtot)
[1056]	305
	306	c------------------------------------------------------------------
	307	c moyennes globales pour les profils de pression et de temperature
[1126]	308	if(planet_type.eq."titan".or.planet_type.eq."venus") then
[1056]	309	call AllGather_Field(pp,iip1*jjp1,llmp1)
	310	call AllGather_Field(pteta,iip1*jjp1,llm)
	311	call AllGather_Field(ppk,iip1*jjp1,llm)
[1126]	312	call AllGather_Field(pphi,iip1*jjp1,llm)
	313	call AllGather_Field(pphis,iip1*jjp1,1)
[1056]	314	ALLOCATE(plevmoy(llm+1))
	315	ALLOCATE(playmoy(llm))
	316	ALLOCATE(tmoy(llm))
	317	ALLOCATE(tetamoy(llm))
	318	ALLOCATE(pkmoy(llm))
[1126]	319	ALLOCATE(phimoy(0:llm))
	320	ALLOCATE(zlevmoy(llm+1))
	321	ALLOCATE(zlaymoy(llm))
[1056]	322	plevmoy=0.
	323	do l=1,llmp1
	324	do i=1,iip1
	325	do j=1,jjp1
	326	plevmoy(l)=plevmoy(l)+pp(i,j,l)/(iip1*jjp1)
	327	enddo
	328	enddo
	329	enddo
	330	tetamoy=0.
	331	pkmoy=0.
[1126]	332	phimoy=0.
	333	do i=1,iip1
	334	do j=1,jjp1
	335	phimoy(0)=phimoy(0)+pphis(i,j)/(iip1*jjp1)
	336	enddo
	337	enddo
[1056]	338	do l=1,llm
	339	do i=1,iip1
	340	do j=1,jjp1
	341	tetamoy(l)=tetamoy(l)+pteta(i,j,l)/(iip1*jjp1)
	342	pkmoy(l)=pkmoy(l)+ppk(i,j,l)/(iip1*jjp1)
[1126]	343	phimoy(l)=phimoy(l)+pphi(i,j,l)/(iip1*jjp1)
[1056]	344	enddo
	345	enddo
	346	enddo
[1126]	347	playmoy(:) = preff * (pkmoy(:)/cpp) ** unskap
[1056]	348	call tpot2t_p(1,llm,tetamoy,tmoy,pkmoy)
[1126]	349	c SI ON TIENT COMPTE DE LA VARIATION DE G AVEC L'ALTITUDE:
	350	zlaymoy(:) = gradrad/(g*rad-phimoy(:))-rad
	351	zlevmoy(1) = phimoy(0)/g
	352	DO l=2,llm
	353	zz1=(playmoy(l-1)+plevmoy(l))/(playmoy(l-1)-plevmoy(l))
	354	zz2=(plevmoy(l) +playmoy(l))/(plevmoy(l) -playmoy(l))
	355	zlevmoy(l)=(zz1zlaymoy(l-1)+zz2zlaymoy(l))/(zz1+zz2)
	356	ENDDO
	357	zlevmoy(llmp1)=zlaymoy(llm)+(zlaymoy(llm)-zlevmoy(llm))
[1071]	358	c-------------------
[1056]	359	c + lat index
[1071]	360	allocate(klat(klon))
	361	do ig0=1,klon
[1126]	362	j=index_j(ig0)
	363	klat(ig0)=j
[1071]	364	enddo
	365	endif ! planet_type=titan
[1056]	366	c------------------------------------------------------------------
	367
[1]	368	c$OMP END MASTER
	369	c$OMP BARRIER
	370	ELSE
	371	debut = .FALSE.
	372	ENDIF
	373
[1056]	374
[1]	375	c-----------------------------------------------------------------------
	376	c 40. transformation des variables dynamiques en variables physiques:
	377	c ---------------------------------------------------------------
	378
	379	c 41. pressions au sol (en Pascals)
	380	c ----------------------------------
	381
	382	c$OMP MASTER
	383	call start_timer(timer_physic)
	384	c$OMP END MASTER
	385
	386	c$OMP MASTER
	387	!CDIR ON_ADB(index_i)
	388	!CDIR ON_ADB(index_j)
	389	do ig0=1,klon
	390	i=index_i(ig0)
	391	j=index_j(ig0)
	392	zpsrf(ig0)=pps(i,j)
	393	enddo
	394	c$OMP END MASTER
	395
	396
[8]	397	c 42. pression intercouches et fonction d'Exner:
[1056]	398
[1]	399	c -----------------------------------------------------------------
	400	c .... zplev definis aux (llm +1) interfaces des couches ....
	401	c .... zplay definis aux ( llm ) milieux des couches ....
	402	c -----------------------------------------------------------------
	403
	404	c ... Exner = cp * ( p(l) / preff ) ** kappa ....
[1056]	405
[1]	406	c print *,omp_rank,'klon--->',klon
	407	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	408	DO l = 1, llmp1
	409	!CDIR ON_ADB(index_i)
	410	!CDIR ON_ADB(index_j)
	411	do ig0=1,klon
	412	i=index_i(ig0)
	413	j=index_j(ig0)
	414	zplev( ig0,l ) = pp(i,j,l)
	415	enddo
	416	ENDDO
	417	c$OMP END DO NOWAIT
[1056]	418	if (flag_moyzon) then
	419	call AllGather_Field(pp,iip1*jjp1,llmp1)
	420	j=index_j(1)
	421	tmpvarp1(:,:) = pp(:,j,:)
	422	call moyzon(llmp1,tmpvarp1,tmpvarbarp1)
	423	zplevbar_mpi(1,:) = tmpvarbarp1
	424	do ig0=2,klon
	425	j=index_j(ig0)
	426	if (j.ne.index_j(ig0-1)) then
	427	tmpvarp1(:,:) = pp(:,j,:)
	428	call moyzon(llmp1,tmpvarp1,tmpvarbarp1)
	429	zplevbar_mpi(ig0,:) = tmpvarbarp1
	430	else
	431	zplevbar_mpi(ig0,:) = zplevbar_mpi(ig0-1,:)
	432	endif
	433	enddo
	434	endif
	435
[8]	436	! ADAPTATION GCM POUR CP(T)
[1056]	437	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[8]	438	DO l=1,llm
	439	!CDIR ON_ADB(index_i)
	440	!CDIR ON_ADB(index_j)
	441	do ig0=1,klon
	442	i=index_i(ig0)
	443	j=index_j(ig0)
	444	zpk(ig0,l)=ppk(i,j,l)
	445	zteta(ig0,l)=pteta(i,j,l)
	446	enddo
	447	ENDDO
	448	c$OMP END DO NOWAIT
[1056]	449	if (flag_moyzon) then
	450	call AllGather_Field(pteta,iip1*jjp1,llm)
	451	call AllGather_Field(ppk,iip1*jjp1,llm)
	452	j=index_j(1)
	453	tmpvar(:,:) = pteta(:,j,:)
	454	call moyzon(llm,tmpvar,tmpvarbar)
	455	ztetabar_mpi(1,:) = tmpvarbar
	456	tmpvar(:,:) = ppk(:,j,:)
	457	call moyzon(llm,tmpvar,tmpvarbar)
	458	zpkbar_mpi(1,:) = tmpvarbar
	459	call tpot2t_p(1,llm,ztetabar_mpi(1,:),ztfibar_mpi(1,:),
	460	& zpkbar_mpi(1,:))
	461	do ig0=2,klon
	462	j=index_j(ig0)
	463	if (j.ne.index_j(ig0-1)) then
	464	tmpvar(:,:) = pteta(:,j,:)
	465	call moyzon(llm,tmpvar,tmpvarbar)
	466	ztetabar_mpi(ig0,:) = tmpvarbar
	467	tmpvar(:,:) = ppk(:,j,:)
	468	call moyzon(llm,tmpvar,tmpvarbar)
	469	zpkbar_mpi(ig0,:) = tmpvarbar
	470	call tpot2t_p(1,llm,ztetabar_mpi(ig0,:),ztfibar_mpi(ig0,:),
	471	& zpkbar_mpi(ig0,:))
	472	else
	473	zpkbar_mpi(ig0,:) = zpkbar_mpi(ig0-1,:)
	474	ztetabar_mpi(ig0,:) = ztetabar_mpi(ig0-1,:)
	475	ztfibar_mpi(ig0,:) = ztfibar_mpi(ig0-1,:)
	476	endif
	477	enddo
	478	endif
[8]	479
[1]	480	c 43. temperature naturelle (en K) et pressions milieux couches .
	481	c ---------------------------------------------------------------
[8]	482
	483	! ADAPTATION GCM POUR CP(T)
[847]	484	call tpot2t_p(klon,llm,zteta,ztfi,zpk)
[8]	485
[1]	486	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	487	DO l=1,llm
	488	!CDIR ON_ADB(index_i)
	489	!CDIR ON_ADB(index_j)
	490	do ig0=1,klon
	491	i=index_i(ig0)
	492	j=index_j(ig0)
	493	pksurcp = ppk(i,j,l) / cpp
	494	zplay(ig0,l) = preff * pksurcp ** unskap
	495	enddo
	496	ENDDO
	497	c$OMP END DO NOWAIT
[1056]	498	if (flag_moyzon) then
	499	zplaybar_mpi(:,:) = preff * (zpkbar_mpi(:,:)/cpp)**unskap
	500	endif
[1]	501
[8]	502	c 43.bis traceurs (tous intensifs)
[1]	503	c ---------------
	504
	505	DO iq=1,nqtot
	506	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	507	DO l=1,llm
	508	!CDIR ON_ADB(index_i)
	509	!CDIR ON_ADB(index_j)
	510	do ig0=1,klon
	511	i=index_i(ig0)
	512	j=index_j(ig0)
[1056]	513	zqfi(ig0,l,iq) = pq(i,j,l,iq)
[1]	514	enddo
	515	ENDDO
	516	c$OMP END DO NOWAIT
[8]	517	ENDDO ! of DO iq=1,nqtot
[1056]	518	if (flag_moyzon) then
	519	DO iq=1,nqtot
	520	call AllGather_Field(pq(:,:,:,iq),iip1*jjp1,llm)
	521	j=index_j(1)
	522	tmpvar(:,:) = pq(:,j,:,iq)
	523	call moyzon(llm,tmpvar,tmpvarbar)
	524	zqfibar_mpi(1,:,iq) = tmpvarbar
	525	do ig0=2,klon
	526	j=index_j(ig0)
	527	if (j.ne.index_j(ig0-1)) then
	528	tmpvar(:,:) = pq(:,j,:,iq)
	529	call moyzon(llm,tmpvar,tmpvarbar)
	530	zqfibar_mpi(ig0,:,iq) = tmpvarbar
	531	else
	532	zqfibar_mpi(ig0,:,iq) = zqfibar_mpi(ig0-1,:,iq)
	533	endif
	534	enddo
	535	ENDDO ! of DO iq=1,nqtot
	536	endif
[1]	537
	538
	539	c Geopotentiel calcule par rapport a la surface locale:
	540	c -----------------------------------------------------
	541
	542	CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,pphi,zphi)
	543
	544	CALL gr_dyn_fi_p(1,iip1,jjp1,klon,pphis,zphis)
	545
	546	c$OMP BARRIER
	547	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	548	DO l=1,llm
	549	DO ig=1,klon
	550	zphi(ig,l)=zphi(ig,l)-zphis(ig)
	551	ENDDO
	552	ENDDO
	553	c$OMP END DO NOWAIT
	554
[1056]	555	if (flag_moyzon) then
	556	call AllGather_Field(pphis,iip1*jjp1,1)
	557	call AllGather_Field(pphi,iip1*jjp1,llm)
	558	j=index_j(1)
	559	tmpvar(:,1) = pphis(:,j)
	560	call moyzon(1,tmpvar(:,1),tmpvarbar(1))
	561	zphisbar_mpi(1) = tmpvarbar(1)
	562	tmpvar(:,:) = pphi(:,j,:)
	563	call moyzon(llm,tmpvar,tmpvarbar)
	564	zphibar_mpi(1,:) = tmpvarbar-zphisbar_mpi(1)
	565	do ig0=2,klon
	566	j=index_j(ig0)
	567	if (j.ne.index_j(ig0-1)) then
	568	tmpvar(:,1) = pphis(:,j)
	569	call moyzon(1,tmpvar(:,1),tmpvarbar(1))
	570	zphisbar_mpi(ig0) = tmpvarbar(1)
	571	tmpvar(:,:) = pphi(:,j,:)
	572	call moyzon(llm,tmpvar,tmpvarbar)
	573	zphibar_mpi(ig0,:) = tmpvarbar-zphisbar_mpi(ig0)
	574	else
	575	zphisbar_mpi(ig0) = zphisbar_mpi(ig0-1)
	576	zphibar_mpi(ig0,:) = zphibar_mpi(ig0-1,:)
	577	endif
	578	enddo
	579	endif
[1]	580
	581	c
	582	c 45. champ u:
	583	c ------------
	584
	585	kstart=1
	586	kend=klon
	587
	588	if (is_north_pole) kstart=2
	589	if (is_south_pole) kend=klon-1
	590
	591	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	592	DO l=1,llm
	593	!CDIR ON_ADB(index_i)
	594	!CDIR ON_ADB(index_j)
	595	!CDIR SPARSE
	596	do ig0=kstart,kend
	597	i=index_i(ig0)
	598	j=index_j(ig0)
	599	if (i==1) then
	600	zufi(ig0,l)= 0.5 *( pucov(iim,j,l)/cu(iim,j)
	601	$ + pucov(1,j,l)/cu(1,j) )
	602	else
	603	zufi(ig0,l)= 0.5*( pucov(i-1,j,l)/cu(i-1,j)
	604	$ + pucov(i,j,l)/cu(i,j) )
	605	endif
	606	enddo
	607	ENDDO
	608	c$OMP END DO NOWAIT
	609
	610	c 46.champ v:
	611	c -----------
	612
	613	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	614	DO l=1,llm
	615	!CDIR ON_ADB(index_i)
	616	!CDIR ON_ADB(index_j)
	617	DO ig0=kstart,kend
	618	i=index_i(ig0)
	619	j=index_j(ig0)
	620	zvfi(ig0,l)= 0.5 *( pvcov(i,j-1,l)/cv(i,j-1)
	621	$ + pvcov(i,j,l)/cv(i,j) )
	622
	623	ENDDO
	624	ENDDO
	625	c$OMP END DO NOWAIT
	626
	627	c 47. champs de vents aux pole nord
	628	c ------------------------------
	629	c U = 1 / pi * integrale [ v * cos(long) * d long ]
	630	c V = 1 / pi * integrale [ v * sin(long) * d long ]
	631
	632	if (is_north_pole) then
	633	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	634	DO l=1,llm
	635
	636	z1(1) =(rlonu(1)-rlonu(iim)+2.pi)pvcov(1,1,l)/cv(1,1)
	637	DO i=2,iim
	638	z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,1,l)/cv(i,1)
	639	ENDDO
	640
	641	DO i=1,iim
	642	zcos(i) = COS(rlonv(i))*z1(i)
	643	zsin(i) = SIN(rlonv(i))*z1(i)
	644	ENDDO
	645
	646	zufi(1,l) = SSUM(iim,zcos,1)/pi
	647	zvfi(1,l) = SSUM(iim,zsin,1)/pi
	648
	649	ENDDO
	650	c$OMP END DO NOWAIT
	651	endif
	652
	653
	654	c 48. champs de vents aux pole sud:
	655	c ---------------------------------
	656	c U = 1 / pi * integrale [ v * cos(long) * d long ]
	657	c V = 1 / pi * integrale [ v * sin(long) * d long ]
	658
	659	if (is_south_pole) then
	660	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	661	DO l=1,llm
	662
	663	z1(1) =(rlonu(1)-rlonu(iim)+2.pi)pvcov(1,jjm,l)/cv(1,jjm)
	664	DO i=2,iim
	665	z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,jjm,l)/cv(i,jjm)
	666	ENDDO
	667
	668	DO i=1,iim
	669	zcos(i) = COS(rlonv(i))*z1(i)
	670	zsin(i) = SIN(rlonv(i))*z1(i)
	671	ENDDO
	672
	673	zufi(klon,l) = SSUM(iim,zcos,1)/pi
	674	zvfi(klon,l) = SSUM(iim,zsin,1)/pi
	675	ENDDO
	676	c$OMP END DO NOWAIT
	677	endif
	678
	679	c On change de grille, dynamique vers physiq, pour le flux de masse verticale
	680	CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,flxw,flxwfi)
	681
	682	c-----------------------------------------------------------------------
	683	c Appel de la physique:
	684	c ---------------------
	685
[8]	686	! Appel de la physique: pose probleme quand on tourne
	687	! SANS physique, car physiq.F est dans le repertoire phy[]...
	688	! Il faut une cle CPP_PHYS
[1]	689
[8]	690	! Le fait que les arguments de physiq soient differents selon les planetes
	691	! ne pose pas de probleme a priori.
	692
	693
[1]	694	c$OMP BARRIER
	695	if (first_omp) then
	696	klon=klon_omp
	697
	698	allocate(zplev_omp(klon,llm+1))
	699	allocate(zplay_omp(klon,llm))
[892]	700	allocate(zpk_omp(klon,llm))
[1]	701	allocate(zphi_omp(klon,llm))
	702	allocate(zphis_omp(klon))
	703	allocate(presnivs_omp(llm))
	704	allocate(zufi_omp(klon,llm))
	705	allocate(zvfi_omp(klon,llm))
	706	allocate(ztfi_omp(klon,llm))
	707	allocate(zqfi_omp(klon,llm,nqtot))
	708	allocate(zdufi_omp(klon,llm))
	709	allocate(zdvfi_omp(klon,llm))
	710	allocate(zdtfi_omp(klon,llm))
	711	allocate(zdqfi_omp(klon,llm,nqtot))
	712	allocate(zdufic_omp(klon,llm))
	713	allocate(zdvfic_omp(klon,llm))
	714	allocate(zdtfic_omp(klon,llm))
	715	allocate(zdqfic_omp(klon,llm,nqtot))
	716	allocate(zdpsrf_omp(klon))
	717	allocate(flxwfi_omp(klon,llm))
[1056]	718
[1251]	719	if (flag_moyzon) call moyzon_init_omp(klon,llm,nqtot)
[1056]	720
[1]	721	first_omp=.false.
	722	endif
	723
	724
	725	klon=klon_omp
	726	offset=klon_omp_begin-1
	727
	728	do l=1,llm+1
	729	do i=1,klon
	730	zplev_omp(i,l)=zplev(offset+i,l)
	731	enddo
	732	enddo
	733
	734	do l=1,llm
	735	do i=1,klon
	736	zplay_omp(i,l)=zplay(offset+i,l)
	737	enddo
	738	enddo
	739
[892]	740	do l=1,llm
	741	do i=1,klon
	742	zpk_omp(i,l)=zpk(offset+i,l)
	743	enddo
	744	enddo
	745
[1]	746	do l=1,llm
	747	do i=1,klon
	748	zphi_omp(i,l)=zphi(offset+i,l)
	749	enddo
	750	enddo
	751
	752	do i=1,klon
	753	zphis_omp(i)=zphis(offset+i)
	754	enddo
	755
	756
	757	do l=1,llm
	758	presnivs_omp(l)=presnivs(l)
	759	enddo
	760
	761	do l=1,llm
	762	do i=1,klon
	763	zufi_omp(i,l)=zufi(offset+i,l)
	764	enddo
	765	enddo
	766
	767	do l=1,llm
	768	do i=1,klon
	769	zvfi_omp(i,l)=zvfi(offset+i,l)
	770	enddo
	771	enddo
	772
	773	do l=1,llm
	774	do i=1,klon
	775	ztfi_omp(i,l)=ztfi(offset+i,l)
	776	enddo
	777	enddo
	778
	779	do iq=1,nqtot
	780	do l=1,llm
	781	do i=1,klon
	782	zqfi_omp(i,l,iq)=zqfi(offset+i,l,iq)
	783	enddo
	784	enddo
	785	enddo
	786
	787	do l=1,llm
	788	do i=1,klon
	789	zdufi_omp(i,l)=zdufi(offset+i,l)
	790	enddo
	791	enddo
	792
	793	do l=1,llm
	794	do i=1,klon
	795	zdvfi_omp(i,l)=zdvfi(offset+i,l)
	796	enddo
	797	enddo
	798
	799	do l=1,llm
	800	do i=1,klon
	801	zdtfi_omp(i,l)=zdtfi(offset+i,l)
	802	enddo
	803	enddo
	804
	805	do iq=1,nqtot
	806	do l=1,llm
	807	do i=1,klon
	808	zdqfi_omp(i,l,iq)=zdqfi(offset+i,l,iq)
	809	enddo
	810	enddo
	811	enddo
	812
	813	do i=1,klon
	814	zdpsrf_omp(i)=zdpsrf(offset+i)
	815	enddo
	816
	817	do l=1,llm
	818	do i=1,klon
	819	flxwfi_omp(i,l)=flxwfi(offset+i,l)
	820	enddo
	821	enddo
[1056]	822
	823	if (flag_moyzon) then
	824	do l=1,llm+1
	825	do i=1,klon
	826	zplevbar(i,l)=zplevbar_mpi(offset+i,l)
	827	enddo
	828	enddo
	829
	830	do l=1,llm
	831	do i=1,klon
	832	zplaybar(i,l)=zplaybar_mpi(offset+i,l)
	833	enddo
	834	enddo
	835
	836	do l=1,llm
	837	do i=1,klon
	838	zphibar(i,l)=zphibar_mpi(offset+i,l)
	839	enddo
	840	enddo
	841
	842	do i=1,klon
	843	zphisbar(i)=zphisbar_mpi(offset+i)
	844	enddo
	845
	846	do l=1,llm
	847	do i=1,klon
	848	ztfibar(i,l)=ztfibar_mpi(offset+i,l)
	849	enddo
	850	enddo
	851
	852	do iq=1,nqtot
	853	do l=1,llm
	854	do i=1,klon
	855	zqfibar(i,l,iq)=zqfibar_mpi(offset+i,l,iq)
	856	enddo
	857	enddo
	858	enddo
	859	endif
	860
[1]	861
	862	c$OMP BARRIER
	863
	864	!$OMP MASTER
	865	! write(lunout,*) 'PHYSIQUE AVEC NSPLIT_PHYS=',nsplit_phys
	866	!$OMP END MASTER
	867	zdt_split=dtphys/nsplit_phys
	868	zdufic_omp(:,:)=0.
	869	zdvfic_omp(:,:)=0.
	870	zdtfic_omp(:,:)=0.
	871	zdqfic_omp(:,:,:)=0.
	872
[776]	873	#ifdef CPP_PHYS
[1]	874	do isplit=1,nsplit_phys
	875
	876	jH_cur_split=jH_cur+(isplit-1) * dtvr / (daysec *nsplit_phys)
	877	debut_split=debut.and.isplit==1
	878	lafin_split=lafin.and.isplit==nsplit_phys
	879
	880
[8]	881	if (planet_type=="earth") then
	882	CALL physiq (klon,
[1]	883	. llm,
	884	. debut_split,
	885	. lafin_split,
	886	. jD_cur,
	887	. jH_cur_split,
	888	. zdt_split,
	889	. zplev_omp,
	890	. zplay_omp,
	891	. zphi_omp,
	892	. zphis_omp,
	893	. presnivs_omp,
	894	. zufi_omp,
	895	. zvfi_omp,
	896	. ztfi_omp,
	897	. zqfi_omp,
	898	. flxwfi_omp,
	899	. zdufi_omp,
	900	. zdvfi_omp,
	901	. zdtfi_omp,
	902	. zdqfi_omp,
	903	. zdpsrf_omp,
[1302]	904	. pducov)
[841]	905
	906	else if ( planet_type=="generic" ) then
	907
	908	CALL physiq (klon, !! ngrid
	909	. llm, !! nlayer
	910	. nqtot, !! nq
	911	. tname, !! tracer names from dynamical core (given in infotrac)
	912	. debut_split, !! firstcall
	913	. lafin_split, !! lastcall
	914	. jD_cur, !! pday. see leapfrog_p
	915	. jH_cur_split, !! ptime "fraction of day"
	916	. zdt_split, !! ptimestep
	917	. zplev_omp, !! pplev
	918	. zplay_omp, !! pplay
	919	. zphi_omp, !! pphi
	920	. zufi_omp, !! pu
	921	. zvfi_omp, !! pv
	922	. ztfi_omp, !! pt
	923	. zqfi_omp, !! pq
	924	. flxwfi_omp, !! pw !! or 0. anyway this is for diagnostic. not used in physiq.
	925	. zdufi_omp, !! pdu
	926	. zdvfi_omp, !! pdv
	927	. zdtfi_omp, !! pdt
	928	. zdqfi_omp, !! pdq
	929	. zdpsrf_omp, !! pdpsrf
	930	. tracerdyn) !! tracerdyn <-- utilite ???
	931
[1107]	932	else if ( planet_type=="mars" ) then
	933
	934	CALL physiq (klon, ! ngrid
	935	. llm, ! nlayer
	936	. nqtot, ! nq
	937	. debut_split, ! firstcall
	938	. lafin_split, ! lastcall
	939	. jD_cur, ! pday
	940	. jH_cur_split, ! ptime
	941	. zdt_split, ! ptimestep
	942	. zplev_omp, ! pplev
	943	. zplay_omp, ! pplay
	944	. zphi_omp, ! pphi
	945	. zufi_omp, ! pu
	946	. zvfi_omp, ! pv
	947	. ztfi_omp, ! pt
	948	. zqfi_omp, ! pq
	949	. flxwfi_omp, ! pw
	950	. zdufi_omp, ! pdu
	951	. zdvfi_omp, ! pdv
	952	. zdtfi_omp, ! pdt
	953	. zdqfi_omp, ! pdq
	954	. zdpsrf_omp, ! pdpsrf
	955	. tracerdyn) ! tracerdyn (somewhat obsolete)
	956
[1129]	957	else if ((planet_type=="titan").or.(planet_type=="venus")) then
	958
	959	CALL physiq (klon,
	960	. llm,
	961	. nqtot,
	962	. debut_split,
	963	. lafin_split,
	964	. jD_cur,
	965	. jH_cur_split,
	966	. zdt_split,
	967	. zplev_omp,
	968	. zplay_omp,
	969	. zpk_omp,
	970	. zphi_omp,
	971	. zphis_omp,
	972	. presnivs_omp,
	973	. zufi_omp,
	974	. zvfi_omp,
	975	. ztfi_omp,
	976	. zqfi_omp,
	977	. flxwfi_omp,
	978	. zdufi_omp,
	979	. zdvfi_omp,
	980	. zdtfi_omp,
	981	. zdqfi_omp,
	982	. zdpsrf_omp)
	983
[1107]	984	else ! unknown "planet_type"
	985
	986	write(lunout,*) "calfis_p: error, unknown planet_type: ",
	987	& trim(planet_type)
	988	stop
	989
[8]	990	endif ! planet_type
[1]	991	zufi_omp(:,:)=zufi_omp(:,:)+zdufi_omp(:,:)*zdt_split
	992	zvfi_omp(:,:)=zvfi_omp(:,:)+zdvfi_omp(:,:)*zdt_split
	993	ztfi_omp(:,:)=ztfi_omp(:,:)+zdtfi_omp(:,:)*zdt_split
	994	zqfi_omp(:,:,:)=zqfi_omp(:,:,:)+zdqfi_omp(:,:,:)*zdt_split
	995
	996	zdufic_omp(:,:)=zdufic_omp(:,:)+zdufi_omp(:,:)
	997	zdvfic_omp(:,:)=zdvfic_omp(:,:)+zdvfi_omp(:,:)
	998	zdtfic_omp(:,:)=zdtfic_omp(:,:)+zdtfi_omp(:,:)
	999	zdqfic_omp(:,:,:)=zdqfic_omp(:,:,:)+zdqfi_omp(:,:,:)
	1000
	1001	enddo
	1002
[776]	1003	#endif
	1004	! of #ifdef CPP_PHYS
	1005
[1]	1006	zdufi_omp(:,:)=zdufic_omp(:,:)/nsplit_phys
	1007	zdvfi_omp(:,:)=zdvfic_omp(:,:)/nsplit_phys
	1008	zdtfi_omp(:,:)=zdtfic_omp(:,:)/nsplit_phys
	1009	zdqfi_omp(:,:,:)=zdqfic_omp(:,:,:)/nsplit_phys
	1010
	1011	c$OMP BARRIER
	1012
	1013	do l=1,llm+1
	1014	do i=1,klon
	1015	zplev(offset+i,l)=zplev_omp(i,l)
	1016	enddo
	1017	enddo
	1018
	1019	do l=1,llm
	1020	do i=1,klon
	1021	zplay(offset+i,l)=zplay_omp(i,l)
	1022	enddo
	1023	enddo
	1024
	1025	do l=1,llm
	1026	do i=1,klon
	1027	zphi(offset+i,l)=zphi_omp(i,l)
	1028	enddo
	1029	enddo
	1030
	1031
	1032	do i=1,klon
	1033	zphis(offset+i)=zphis_omp(i)
	1034	enddo
	1035
	1036
	1037	do l=1,llm
	1038	presnivs(l)=presnivs_omp(l)
	1039	enddo
	1040
	1041	do l=1,llm
	1042	do i=1,klon
	1043	zufi(offset+i,l)=zufi_omp(i,l)
	1044	enddo
	1045	enddo
	1046
	1047	do l=1,llm
	1048	do i=1,klon
	1049	zvfi(offset+i,l)=zvfi_omp(i,l)
	1050	enddo
	1051	enddo
	1052
	1053	do l=1,llm
	1054	do i=1,klon
	1055	ztfi(offset+i,l)=ztfi_omp(i,l)
	1056	enddo
	1057	enddo
	1058
	1059	do iq=1,nqtot
	1060	do l=1,llm
	1061	do i=1,klon
	1062	zqfi(offset+i,l,iq)=zqfi_omp(i,l,iq)
	1063	enddo
	1064	enddo
	1065	enddo
	1066
	1067	do l=1,llm
	1068	do i=1,klon
	1069	zdufi(offset+i,l)=zdufi_omp(i,l)
	1070	enddo
	1071	enddo
	1072
	1073	do l=1,llm
	1074	do i=1,klon
	1075	zdvfi(offset+i,l)=zdvfi_omp(i,l)
	1076	enddo
	1077	enddo
	1078
	1079	do l=1,llm
	1080	do i=1,klon
	1081	zdtfi(offset+i,l)=zdtfi_omp(i,l)
	1082	enddo
	1083	enddo
	1084
	1085	do iq=1,nqtot
	1086	do l=1,llm
	1087	do i=1,klon
	1088	zdqfi(offset+i,l,iq)=zdqfi_omp(i,l,iq)
	1089	enddo
	1090	enddo
	1091	enddo
	1092
	1093	do i=1,klon
	1094	zdpsrf(offset+i)=zdpsrf_omp(i)
	1095	enddo
	1096
	1097
	1098	klon=klon_mpi
	1099	500 CONTINUE
	1100	c$OMP BARRIER
	1101
	1102	c$OMP MASTER
	1103	call stop_timer(timer_physic)
	1104	c$OMP END MASTER
	1105
	1106	IF (using_mpi) THEN
	1107
	1108	if (MPI_rank>0) then
	1109
	1110	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	1111	DO l=1,llm
	1112	du_send(1:iim,l)=zdufi(1:iim,l)
	1113	dv_send(1:iim,l)=zdvfi(1:iim,l)
	1114	ENDDO
	1115	c$OMP END DO NOWAIT
	1116
	1117	c$OMP BARRIER
	1118	#ifdef CPP_MPI
	1119	c$OMP MASTER
	1120	!$OMP CRITICAL (MPI)
	1121	call MPI_ISSEND(du_send,iim*llm,MPI_REAL8,MPI_Rank-1,401,
	1122	& COMM_LMDZ,Req(1),ierr)
	1123	call MPI_ISSEND(dv_send,iim*llm,MPI_REAL8,MPI_Rank-1,402,
	1124	& COMM_LMDZ,Req(2),ierr)
	1125	!$OMP END CRITICAL (MPI)
	1126	c$OMP END MASTER
	1127	#endif
	1128	c$OMP BARRIER
	1129
	1130	endif
	1131
	1132	if (MPI_rank<MPI_Size-1) then
	1133	c$OMP BARRIER
	1134	#ifdef CPP_MPI
	1135	c$OMP MASTER
	1136	!$OMP CRITICAL (MPI)
	1137	call MPI_IRECV(du_recv,iim*llm,MPI_REAL8,MPI_Rank+1,401,
	1138	& COMM_LMDZ,Req(3),ierr)
	1139	call MPI_IRECV(dv_recv,iim*llm,MPI_REAL8,MPI_Rank+1,402,
	1140	& COMM_LMDZ,Req(4),ierr)
	1141	!$OMP END CRITICAL (MPI)
	1142	c$OMP END MASTER
	1143	#endif
	1144	endif
	1145
	1146	c$OMP BARRIER
	1147
	1148
	1149	#ifdef CPP_MPI
	1150	c$OMP MASTER
	1151	!$OMP CRITICAL (MPI)
	1152	if (MPI_rank>0 .and. MPI_rank< MPI_Size-1) then
	1153	call MPI_WAITALL(4,Req(1),Status,ierr)
	1154	else if (MPI_rank>0) then
	1155	call MPI_WAITALL(2,Req(1),Status,ierr)
	1156	else if (MPI_rank <MPI_Size-1) then
	1157	call MPI_WAITALL(2,Req(3),Status,ierr)
	1158	endif
	1159	!$OMP END CRITICAL (MPI)
	1160	c$OMP END MASTER
	1161	#endif
	1162
	1163	c$OMP BARRIER
	1164
	1165	ENDIF ! using_mpi
	1166
	1167
	1168	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	1169	DO l=1,llm
	1170
	1171	zdufi2(1:klon,l)=zdufi(1:klon,l)
	1172	zdufi2(klon+1:klon+iim,l)=du_recv(1:iim,l)
	1173
	1174	zdvfi2(1:klon,l)=zdvfi(1:klon,l)
	1175	zdvfi2(klon+1:klon+iim,l)=dv_recv(1:iim,l)
	1176
	1177	pdhfi(:,jj_begin,l)=0
	1178	pdqfi(:,jj_begin,l,:)=0
	1179	pdufi(:,jj_begin,l)=0
	1180	pdvfi(:,jj_begin,l)=0
	1181
	1182	if (.not. is_south_pole) then
	1183	pdhfi(:,jj_end,l)=0
	1184	pdqfi(:,jj_end,l,:)=0
	1185	pdufi(:,jj_end,l)=0
	1186	pdvfi(:,jj_end,l)=0
	1187	endif
	1188
	1189	ENDDO
	1190	c$OMP END DO NOWAIT
	1191
	1192	c$OMP MASTER
	1193	pdpsfi(:,jj_begin)=0
	1194	if (.not. is_south_pole) then
	1195	pdpsfi(:,jj_end)=0
	1196	endif
	1197	c$OMP END MASTER
	1198	c-----------------------------------------------------------------------
	1199	c transformation des tendances physiques en tendances dynamiques:
	1200	c ---------------------------------------------------------------
	1201
	1202	c tendance sur la pression :
	1203	c -----------------------------------
	1204	CALL gr_fi_dyn_p(1,klon,iip1,jjp1,zdpsrf,pdpsfi)
	1205	c
	1206	c 62. enthalpie potentielle
	1207	c ---------------------
[8]	1208
[1]	1209	kstart=1
	1210	kend=klon
	1211
	1212	if (is_north_pole) kstart=2
	1213	if (is_south_pole) kend=klon-1
[8]	1214
	1215	! ADAPTATION GCM POUR CP(T)
[847]	1216	call t2tpot_p(klon,llm,ztfi,zteta,zpk)
[1]	1217
[8]	1218
[1]	1219	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	1220	DO l=1,llm
	1221	!CDIR ON_ADB(index_i)
	1222	!CDIR ON_ADB(index_j)
	1223	!cdir NODEP
	1224	do ig0=kstart,kend
	1225	i=index_i(ig0)
	1226	j=index_j(ig0)
[8]	1227	! pdhfi(i,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l)
	1228	pdhfi(i,j,l) = (zteta(ig0,l) - pteta(i,j,l))/dtphys
	1229	! if (i==1) pdhfi(iip1,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l)
	1230	if (i==1) then
	1231	pdhfi(iip1,j,l) = (zteta(ig0,l) - pteta(i,j,l))/dtphys
	1232	endif
[1]	1233	enddo
	1234
	1235	if (is_north_pole) then
	1236	DO i=1,iip1
[8]	1237	! pdhfi(i,1,l) = cpp * zdtfi(1,l) / ppk(i, 1 ,l)
	1238	pdhfi(i,1,l) = (zteta(1,l) - pteta(i,1,l))/dtphys
[1]	1239	enddo
	1240	endif
	1241
	1242	if (is_south_pole) then
	1243	DO i=1,iip1
[8]	1244	! pdhfi(i,jjp1,l) = cpp * zdtfi(klon,l)/ ppk(i,jjp1,l)
	1245	pdhfi(i,jjp1,l) = (zteta(klon,l) - pteta(i,jjp1,l))/dtphys
[1]	1246	ENDDO
	1247	endif
	1248	ENDDO
	1249	c$OMP END DO NOWAIT
	1250
	1251	c 62. humidite specifique
	1252	c ---------------------
	1253	! Ehouarn: removed this useless bit: was overwritten at step 63 anyways
	1254	! DO iq=1,nqtot
	1255	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	1256	! DO l=1,llm
	1257	!!!cdir NODEP
	1258	! do ig0=kstart,kend
	1259	! i=index_i(ig0)
	1260	! j=index_j(ig0)
	1261	! pdqfi(i,j,l,iq) = zdqfi(ig0,l,iq)
	1262	! if (i==1) pdqfi(iip1,j,l,iq) = zdqfi(ig0,l,iq)
	1263	! enddo
	1264	!
	1265	! if (is_north_pole) then
	1266	! do i=1,iip1
	1267	! pdqfi(i,1,l,iq) = zdqfi(1,l,iq)
	1268	! enddo
	1269	! endif
	1270	!
	1271	! if (is_south_pole) then
	1272	! do i=1,iip1
	1273	! pdqfi(i,jjp1,l,iq) = zdqfi(klon,l,iq)
	1274	! enddo
	1275	! endif
	1276	! ENDDO
	1277	!c$OMP END DO NOWAIT
	1278	! ENDDO
	1279
[8]	1280	c 63. traceurs (tous en intensifs)
[1]	1281	c ------------
	1282	C initialisation des tendances
	1283
	1284	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	1285	DO l=1,llm
	1286	pdqfi(:,:,l,:)=0.
	1287	ENDDO
	1288	c$OMP END DO NOWAIT
	1289
	1290	C
	1291	!cdir NODEP
	1292	DO iq=1,nqtot
	1293	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	1294	DO l=1,llm
	1295	!CDIR ON_ADB(index_i)
	1296	!CDIR ON_ADB(index_j)
	1297	!cdir NODEP
	1298	DO ig0=kstart,kend
	1299	i=index_i(ig0)
	1300	j=index_j(ig0)
[1056]	1301	pdqfi(i,j,l,iq) = zdqfi(ig0,l,iq)
	1302	if (i==1) pdqfi(iip1,j,l,iq) = zdqfi(ig0,l,iq)
[1]	1303	ENDDO
	1304
	1305	IF (is_north_pole) then
	1306	DO i=1,iip1
[1056]	1307	pdqfi(i,1,l,iq) = zdqfi(1,l,iq)
[1]	1308	ENDDO
	1309	ENDIF
	1310
	1311	IF (is_south_pole) then
	1312	DO i=1,iip1
[1056]	1313	pdqfi(i,jjp1,l,iq) = zdqfi(klon,l,iq)
[1]	1314	ENDDO
	1315	ENDIF
	1316
	1317	ENDDO
	1318	c$OMP END DO NOWAIT
	1319	ENDDO
	1320
	1321	c 65. champ u:
	1322	c ------------
	1323	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	1324	DO l=1,llm
	1325	!CDIR ON_ADB(index_i)
	1326	!CDIR ON_ADB(index_j)
	1327	!cdir NODEP
	1328	do ig0=kstart,kend
	1329	i=index_i(ig0)
	1330	j=index_j(ig0)
	1331
	1332	if (i/=iim) then
	1333	pdufi(i,j,l)=0.5(zdufi2(ig0,l)+zdufi2(ig0+1,l))cu(i,j)
	1334	endif
	1335
	1336	if (i==1) then
	1337	pdufi(iim,j,l)=0.5*( zdufi2(ig0,l)
	1338	$ + zdufi2(ig0+iim-1,l))*cu(iim,j)
	1339	pdufi(iip1,j,l)=0.5(zdufi2(ig0,l)+zdufi2(ig0+1,l))cu(i,j)
	1340	endif
	1341
	1342	enddo
	1343
	1344	if (is_north_pole) then
	1345	DO i=1,iip1
	1346	pdufi(i,1,l) = 0.
	1347	ENDDO
	1348	endif
	1349
	1350	if (is_south_pole) then
	1351	DO i=1,iip1
	1352	pdufi(i,jjp1,l) = 0.
	1353	ENDDO
	1354	endif
	1355
	1356	ENDDO
	1357	c$OMP END DO NOWAIT
	1358
	1359	c 67. champ v:
	1360	c ------------
	1361
	1362	kstart=1
	1363	kend=klon
	1364
	1365	if (is_north_pole) kstart=2
	1366	if (is_south_pole) kend=klon-1-iim
	1367
	1368	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	1369	DO l=1,llm
	1370	!CDIR ON_ADB(index_i)
	1371	!CDIR ON_ADB(index_j)
	1372	!cdir NODEP
	1373	do ig0=kstart,kend
	1374	i=index_i(ig0)
	1375	j=index_j(ig0)
	1376	pdvfi(i,j,l)=0.5(zdvfi2(ig0,l)+zdvfi2(ig0+iim,l))cv(i,j)
	1377	if (i==1) pdvfi(iip1,j,l) = 0.5*(zdvfi2(ig0,l)+
	1378	$ zdvfi2(ig0+iim,l))
	1379	$ *cv(i,j)
	1380	enddo
	1381
	1382	ENDDO
	1383	c$OMP END DO NOWAIT
	1384
	1385
	1386	c 68. champ v pres des poles:
	1387	c ---------------------------
	1388	c v = U * cos(long) + V * SIN(long)
	1389
	1390	if (is_north_pole) then
	1391
	1392	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	1393	DO l=1,llm
	1394
	1395	DO i=1,iim
	1396	pdvfi(i,1,l)=
	1397	$ zdufi(1,l)COS(rlonv(i))+zdvfi(1,l)SIN(rlonv(i))
	1398
	1399	pdvfi(i,1,l)=
	1400	$ 0.5(pdvfi(i,1,l)+zdvfi(i+1,l))cv(i,1)
	1401	ENDDO
	1402
	1403	pdvfi(iip1,1,l) = pdvfi(1,1,l)
	1404
	1405	ENDDO
	1406	c$OMP END DO NOWAIT
	1407
	1408	endif
	1409
	1410	if (is_south_pole) then
	1411
	1412	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	1413	DO l=1,llm
	1414
	1415	DO i=1,iim
	1416	pdvfi(i,jjm,l)=zdufi(klon,l)*COS(rlonv(i))
	1417	$ +zdvfi(klon,l)*SIN(rlonv(i))
	1418
	1419	pdvfi(i,jjm,l)=
	1420	$ 0.5(pdvfi(i,jjm,l)+zdvfi(klon-iip1+i,l))cv(i,jjm)
	1421	ENDDO
	1422
	1423	pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l)
	1424
	1425	ENDDO
	1426	c$OMP END DO NOWAIT
	1427
	1428	endif
	1429	c-----------------------------------------------------------------------
	1430
	1431	700 CONTINUE
	1432
	1433	firstcal = .FALSE.
	1434
	1435	#else
	1436	write(lunout,*)
	1437	& "calfis_p: for now can only work with parallel physics"
	1438	stop
	1439	#endif
[8]	1440	! of #ifdef CPP_PHYS
[1403]	1441	#endif
	1442	! of #ifdef CPP_PARA
[1]	1443	END

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

Download in other formats:

Original Format