Context Navigation

calfis_p.F @ 1137

Last change on this file since 1137 was 1114, checked in by jghattas, 15 years ago

Creation du module infotrac:

contient les variables de advtrac.h
contient la subroutine iniadvtrac renommer en infotrac_init
le nombre des traceurs est lu dans tracer.def en dynamique (ou par default ou recu par INCA)
ce module est utilise dans la dynamique et la physique
contient aussi la variable nbtr qui avant etait stockee dans dimphy

Le fichier advtrac.h n'existe plus.
La compilation ne prend plus en compte le nombre de traceur.

/JG

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 28.1 KB

Rev	Line
[630]	1	!
	2	! $Header$
	3	!
	4	C
	5	C
[1114]	6	SUBROUTINE calfis_p(lafin,
[630]	7	$ rdayvrai,
	8	$ heure,
	9	$ pucov,
	10	$ pvcov,
	11	$ pteta,
	12	$ pq,
	13	$ pmasse,
	14	$ pps,
	15	$ pp,
	16	$ ppk,
	17	$ pphis,
	18	$ pphi,
	19	$ pducov,
	20	$ pdvcov,
	21	$ pdteta,
	22	$ pdq,
	23	$ flxw,
	24	$ clesphy0,
	25	$ pdufi,
	26	$ pdvfi,
	27	$ pdhfi,
	28	$ pdqfi,
	29	$ pdpsfi)
	30	c
	31	c Auteur : P. Le Van, F. Hourdin
	32	c .........
	33	USE dimphy
[774]	34	USE mod_phys_lmdz_para, mpi_root_xx=>mpi_root
[1000]	35	USE parallel, ONLY : omp_chunk, using_mpi
[774]	36	USE mod_interface_dyn_phys
[630]	37	USE Write_Field
	38	Use Write_field_p
	39	USE Times
[764]	40	USE IOPHY
[1114]	41	USE infotrac
	42
[630]	43	IMPLICIT NONE
	44	c=======================================================================
	45	c
	46	c 1. rearrangement des tableaux et transformation
	47	c variables dynamiques > variables physiques
	48	c 2. calcul des termes physiques
	49	c 3. retransformation des tendances physiques en tendances dynamiques
	50	c
	51	c remarques:
	52	c ----------
	53	c
	54	c - les vents sont donnes dans la physique par leurs composantes
	55	c naturelles.
	56	c - la variable thermodynamique de la physique est une variable
	57	c intensive : T
	58	c pour la dynamique on prend T * ( preff / p(l) ) **kappa
	59	c - les deux seules variables dependant de la geometrie necessaires
	60	c pour la physique sont la latitude pour le rayonnement et
	61	c l'aire de la maille quand on veut integrer une grandeur
	62	c horizontalement.
	63	c - les points de la physique sont les points scalaires de la
	64	c la dynamique; numerotation:
	65	c 1 pour le pole nord
	66	c (jjm-1)*iim pour l'interieur du domaine
	67	c ngridmx pour le pole sud
	68	c ---> ngridmx=2+(jjm-1)*iim
	69	c
	70	c Input :
	71	c -------
	72	c ecritphy frequence d'ecriture (en jours)de histphy
	73	c pucov covariant zonal velocity
	74	c pvcov covariant meridional velocity
	75	c pteta potential temperature
	76	c pps surface pressure
	77	c pmasse masse d'air dans chaque maille
	78	c pts surface temperature (K)
	79	c callrad clef d'appel au rayonnement
	80	c
	81	c Output :
	82	c --------
	83	c pdufi tendency for the natural zonal velocity (ms-1)
	84	c pdvfi tendency for the natural meridional velocity
	85	c pdhfi tendency for the potential temperature
	86	c pdtsfi tendency for the surface temperature
	87	c
	88	c pdtrad radiative tendencies \ both input
	89	c pfluxrad radiative fluxes / and output
	90	c
	91	c=======================================================================
	92	c
	93	c-----------------------------------------------------------------------
	94	c
	95	c 0. Declarations :
	96	c ------------------
	97
	98	#include "dimensions.h"
	99	#include "paramet.h"
	100	#include "temps.h"
	101
[1114]	102	INTEGER ngridmx
[630]	103	PARAMETER( ngridmx = 2+(jjm-1)*iim - 1/jjm )
	104
	105	#include "comconst.h"
	106	#include "comvert.h"
	107	#include "comgeom2.h"
	108	#include "control.h"
[1000]	109	#ifdef CPP_MPI
[630]	110	include 'mpif.h'
[1000]	111	#endif
[630]	112	c Arguments :
	113	c -----------
	114	LOGICAL lafin
	115	REAL heure
	116
	117	REAL pvcov(iip1,jjm,llm)
	118	REAL pucov(iip1,jjp1,llm)
	119	REAL pteta(iip1,jjp1,llm)
	120	REAL pmasse(iip1,jjp1,llm)
[1114]	121	REAL pq(iip1,jjp1,llm,nqtot)
[630]	122	REAL pphis(iip1,jjp1)
	123	REAL pphi(iip1,jjp1,llm)
	124	c
	125	REAL pdvcov(iip1,jjm,llm)
	126	REAL pducov(iip1,jjp1,llm)
	127	REAL pdteta(iip1,jjp1,llm)
[1114]	128	REAL pdq(iip1,jjp1,llm,nqtot)
[630]	129	c
	130	REAL pps(iip1,jjp1)
	131	REAL pp(iip1,jjp1,llmp1)
	132	REAL ppk(iip1,jjp1,llm)
	133	c
	134	REAL pdvfi(iip1,jjm,llm)
	135	REAL pdufi(iip1,jjp1,llm)
	136	REAL pdhfi(iip1,jjp1,llm)
[1114]	137	REAL pdqfi(iip1,jjp1,llm,nqtot)
[630]	138	REAL pdpsfi(iip1,jjp1)
	139
	140	INTEGER longcles
	141	PARAMETER ( longcles = 20 )
	142	REAL clesphy0( longcles )
	143
	144
	145	c Local variables :
	146	c -----------------
	147
	148	INTEGER i,j,l,ig0,ig,iq,iiq
[764]	149	REAL,ALLOCATABLE,SAVE :: zpsrf(:)
	150	REAL,ALLOCATABLE,SAVE :: zplev(:,:),zplay(:,:)
	151	REAL,ALLOCATABLE,SAVE :: zphi(:,:),zphis(:)
[630]	152	c
[764]	153	REAL,ALLOCATABLE,SAVE :: zufi(:,:), zvfi(:,:)
	154	REAL,ALLOCATABLE,SAVE :: ztfi(:,:),zqfi(:,:,:)
[630]	155	c
[764]	156	REAL,ALLOCATABLE,SAVE :: pcvgu(:,:), pcvgv(:,:)
	157	REAL,ALLOCATABLE,SAVE :: pcvgt(:,:), pcvgq(:,:,:)
[630]	158	c
[960]	159	c REAL,ALLOCATABLE,SAVE :: pvervel(:,:)
[630]	160	c
[764]	161	REAL,ALLOCATABLE,SAVE :: zdufi(:,:),zdvfi(:,:)
	162	REAL,ALLOCATABLE,SAVE :: zdtfi(:,:),zdqfi(:,:,:)
	163	REAL,ALLOCATABLE,SAVE :: zdpsrf(:)
[985]	164	REAL,SAVE,ALLOCATABLE :: flxwfi(:,:) ! Flux de masse verticale sur la grille physiq
	165
[630]	166	c
[764]	167	REAL,ALLOCATABLE,SAVE :: zplev_omp(:,:)
	168	REAL,ALLOCATABLE,SAVE :: zplay_omp(:,:)
	169	REAL,ALLOCATABLE,SAVE :: zphi_omp(:,:)
	170	REAL,ALLOCATABLE,SAVE :: zphis_omp(:)
	171	REAL,ALLOCATABLE,SAVE :: presnivs_omp(:)
	172	REAL,ALLOCATABLE,SAVE :: zufi_omp(:,:)
	173	REAL,ALLOCATABLE,SAVE :: zvfi_omp(:,:)
	174	REAL,ALLOCATABLE,SAVE :: ztfi_omp(:,:)
	175	REAL,ALLOCATABLE,SAVE :: zqfi_omp(:,:,:)
[960]	176	c REAL,ALLOCATABLE,SAVE :: pvervel_omp(:,:)
[764]	177	REAL,ALLOCATABLE,SAVE :: zdufi_omp(:,:)
	178	REAL,ALLOCATABLE,SAVE :: zdvfi_omp(:,:)
	179	REAL,ALLOCATABLE,SAVE :: zdtfi_omp(:,:)
	180	REAL,ALLOCATABLE,SAVE :: zdqfi_omp(:,:,:)
	181	REAL,ALLOCATABLE,SAVE :: zdpsrf_omp(:)
[985]	182	REAL,SAVE,ALLOCATABLE :: flxwfi_omp(:,:) ! Flux de masse verticale sur la grille physiq
[764]	183
	184	c$OMP THREADPRIVATE(zplev_omp,zplay_omp,zphi_omp,zphis_omp,
	185	c$OMP+ presnivs_omp,zufi_omp,zvfi_omp,ztfi_omp,
[985]	186	c$OMP+ zqfi_omp,zdufi_omp,zdvfi_omp,
	187	c$OMP+ zdtfi_omp,zdqfi_omp,zdpsrf_omp,flxwfi_omp)
[764]	188
	189	LOGICAL,SAVE :: first_omp=.true.
	190	c$OMP THREADPRIVATE(first_omp)
	191
[630]	192	REAL zsin(iim),zcos(iim),z1(iim)
	193	REAL zsinbis(iim),zcosbis(iim),z1bis(iim)
	194	REAL unskap, pksurcp
[764]	195	c
	196	cIM diagnostique PVteta, Amip2
	197	INTEGER ntetaSTD
	198	PARAMETER(ntetaSTD=3)
	199	REAL rtetaSTD(ntetaSTD)
	200	DATA rtetaSTD/350., 380., 405./
	201	REAL PVteta(klon,ntetaSTD)
	202
[960]	203	REAL flxw(iip1,jjp1,llm) ! Flux de masse verticale sur la grille dynamique
[630]	204
	205	REAL SSUM
	206
	207	LOGICAL firstcal, debut
	208	DATA firstcal/.true./
	209	SAVE firstcal,debut
[764]	210	c$OMP THREADPRIVATE(firstcal,debut)
[630]	211	REAL rdayvrai
	212
[764]	213	REAL,SAVE,dimension(1:iim,1:llm):: du_send,du_recv,dv_send,dv_recv
[630]	214	INTEGER :: ierr
[1000]	215	#ifdef CPP_MPI
[630]	216	INTEGER,dimension(MPI_STATUS_SIZE,4) :: Status
[1000]	217	#else
	218	INTEGER,dimension(1,4) :: Status
	219	#endif
[630]	220	INTEGER, dimension(4) :: Req
[764]	221	REAL,ALLOCATABLE,SAVE:: zdufi2(:,:),zdvfi2(:,:)
	222	integer :: k,kstart,kend
	223	INTEGER :: offset
[630]	224	c
	225	c-----------------------------------------------------------------------
	226	c
	227	c 1. Initialisations :
	228	c --------------------
	229	c
	230
[764]	231	klon=klon_mpi
	232
	233	PVteta(:,:)=0.
	234
[630]	235	IF (ngridmx.NE.2+(jjm-1)*iim) THEN
	236	PRINT*,'STOP dans calfis'
	237	PRINT,'La dimension ngridmx doit etre egale a 2 + (jjm-1)iim'
	238	PRINT*,' ngridmx jjm iim '
	239	PRINT*,ngridmx,jjm,iim
	240	STOP
	241	ENDIF
	242
	243	c-----------------------------------------------------------------------
	244	c latitude, longitude et aires des mailles pour la physique:
	245	c ----------------------------------------------------------
	246
	247	c
	248	IF ( firstcal ) THEN
	249	debut = .TRUE.
[764]	250	c$OMP MASTER
	251	ALLOCATE(zpsrf(klon))
	252	ALLOCATE(zplev(klon,llm+1),zplay(klon,llm))
	253	ALLOCATE(zphi(klon,llm),zphis(klon))
	254	ALLOCATE(zufi(klon,llm), zvfi(klon,llm))
[1114]	255	ALLOCATE(ztfi(klon,llm),zqfi(klon,llm,nqtot))
[764]	256	ALLOCATE(pcvgu(klon,llm), pcvgv(klon,llm))
	257	ALLOCATE(pcvgt(klon,llm), pcvgq(klon,llm,2))
[960]	258	c ALLOCATE(pvervel(klon,llm))
[764]	259	ALLOCATE(zdufi(klon,llm),zdvfi(klon,llm))
[1114]	260	ALLOCATE(zdtfi(klon,llm),zdqfi(klon,llm,nqtot))
[764]	261	ALLOCATE(zdpsrf(klon))
	262	ALLOCATE(zdufi2(klon+iim,llm),zdvfi2(klon+iim,llm))
[985]	263	ALLOCATE(flxwfi(klon,llm))
[764]	264	c$OMP END MASTER
	265	c$OMP BARRIER
[630]	266	ELSE
	267	debut = .FALSE.
	268	ENDIF
	269
	270	c
	271	c
	272	c-----------------------------------------------------------------------
	273	c 40. transformation des variables dynamiques en variables physiques:
	274	c ---------------------------------------------------------------
	275
	276	c 41. pressions au sol (en Pascals)
	277	c ----------------------------------
	278
[764]	279	c$OMP MASTER
[630]	280	call start_timer(timer_physic)
[764]	281	c$OMP END MASTER
	282
	283	c$OMP MASTER
[630]	284	do ig0=1,klon
[774]	285	i=index_i(ig0)
	286	j=index_j(ig0)
[630]	287	zpsrf(ig0)=pps(i,j)
	288	enddo
[764]	289	c$OMP END MASTER
[630]	290
	291
	292	c 42. pression intercouches :
	293	c
	294	c -----------------------------------------------------------------
	295	c .... zplev definis aux (llm +1) interfaces des couches ....
	296	c .... zplay definis aux ( llm ) milieux des couches ....
	297	c -----------------------------------------------------------------
	298
	299	c ... Exner = cp * ( p(l) / preff ) ** kappa ....
	300	c
	301	unskap = 1./ kappa
	302	c
[961]	303	c print *,omp_rank,'klon--->',klon
[764]	304	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	305	DO l = 1, llmp1
	306	do ig0=1,klon
[774]	307	i=index_i(ig0)
	308	j=index_j(ig0)
[630]	309	zplev( ig0,l ) = pp(i,j,l)
	310	enddo
	311	ENDDO
[764]	312	c$OMP END DO NOWAIT
[630]	313	c
	314	c
	315
	316	c 43. temperature naturelle (en K) et pressions milieux couches .
	317	c ---------------------------------------------------------------
[764]	318	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	319	DO l=1,llm
	320
	321	do ig0=1,klon
[774]	322	i=index_i(ig0)
	323	j=index_j(ig0)
[630]	324	pksurcp = ppk(i,j,l) / cpp
	325	zplay(ig0,l) = preff * pksurcp ** unskap
	326	ztfi(ig0,l) = pteta(i,j,l) * pksurcp
	327	c pcvgt(ig0,l) = pdteta(i,j,l) * pksurcp / pmasse(i,j,l)
	328	enddo
	329
	330	ENDDO
[764]	331	c$OMP END DO NOWAIT
[630]	332
	333	c 43.bis traceurs
	334	c ---------------
	335	c
	336
[1114]	337	DO iq=1,nqtot
[630]	338	iiq=niadv(iq)
[764]	339	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	340	DO l=1,llm
	341	do ig0=1,klon
[774]	342	i=index_i(ig0)
	343	j=index_j(ig0)
[630]	344	zqfi(ig0,l,iq) = pq(i,j,l,iiq)
	345	enddo
	346	ENDDO
[764]	347	c$OMP END DO NOWAIT
[630]	348	ENDDO
	349
	350	c convergence dynamique pour les traceurs "EAU"
	351
	352	DO iq=1,2
[764]	353	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	354	DO l=1,llm
	355	do ig0=1,klon
[774]	356	i=index_i(ig0)
	357	j=index_j(ig0)
[630]	358	c pcvgq(ig0,l,iq) = pdq(i,j,l,iq) / pmasse(i,j,l)
	359	enddo
	360	ENDDO
[764]	361	c$OMP END DO NOWAIT
[630]	362	ENDDO
	363
	364
	365
	366	c Geopotentiel calcule par rapport a la surface locale:
	367	c -----------------------------------------------------
	368
	369	CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,pphi,zphi)
[764]	370
	371	c$OMP MASTER
[630]	372	CALL gr_dyn_fi_p(1,iip1,jjp1,klon,pphis,zphis)
[764]	373	c$OMP END MASTER
	374	c$OMP BARRIER
[630]	375
[764]	376	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	377	DO l=1,llm
	378	DO ig=1,klon
	379	zphi(ig,l)=zphi(ig,l)-zphis(ig)
	380	ENDDO
	381	ENDDO
[960]	382	c$OMP END DO NOWAIT
	383
[630]	384	c .... Calcul de la vitesse verticale ( en Pams ou Kg/s ) ....
[960]	385	c JG : ancien calcule de omega utilise dans physiq.F. Maintenant le flux
	386	c de masse est calclue dans advtrac_p.F
[630]	387	c
[960]	388	cc$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	389	c DO l=1,llm
	390	c do ig0=1,klon
	391	c i=index_i(ig0)
	392	c j=index_j(ig0)
	393	c pvervel(ig0,l) = pw(i,j,l)g unsaire(i,j)
	394	c enddo
	395	c if (is_north_pole) pvervel(1,l)=pw(1,1,l)*g /apoln
	396	c if (is_south_pole) pvervel(klon,l)=pw(1,jjp1,l)*g/apols
	397	c ENDDO
	398	cc$OMP END DO NOWAIT
[630]	399
	400	c
	401	c 45. champ u:
	402	c ------------
	403
	404	kstart=1
	405	kend=klon
	406
[774]	407	if (is_north_pole) kstart=2
	408	if (is_south_pole) kend=klon-1
[630]	409
[764]	410	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	411	DO l=1,llm
	412	do ig0=kstart,kend
[774]	413	i=index_i(ig0)
	414	j=index_j(ig0)
[630]	415	if (i==1) then
	416	zufi(ig0,l)= 0.5 *( pucov(iim,j,l)/cu(iim,j)
	417	$ + pucov(1,j,l)/cu(1,j) )
	418	c pcvgu(ig0,l)= 0.5*( pducov(iim,j,l)/cu(iim,j)
	419	c $ + pducov(1,j,l)/cu(1,j) )
	420	else
	421	zufi(ig0,l)= 0.5*( pucov(i-1,j,l)/cu(i-1,j)
	422	$ + pucov(i,j,l)/cu(i,j) )
	423	c pcvgu(ig0,l)= 0.5*( pducov(i-1,j,l)/cu(i-1,j)
	424	c $ + pducov(i,j,l)/cu(i,j) )
	425	endif
	426	enddo
	427	ENDDO
[764]	428	c$OMP END DO NOWAIT
[630]	429	c 46.champ v:
	430	c -----------
[764]	431	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	432	DO l=1,llm
	433	DO ig0=kstart,kend
[774]	434	i=index_i(ig0)
	435	j=index_j(ig0)
[630]	436	zvfi(ig0,l)= 0.5 *( pvcov(i,j-1,l)/cv(i,j-1)
	437	$ + pvcov(i,j,l)/cv(i,j) )
	438
	439	c pcvgv(ig0+i,l)= 0.5 * ( pdvcov(i,j-1,l)/cv(i,j-1)
	440	c $ + pdvcov(i,j,l)/cv(i,j) )
	441	ENDDO
	442	ENDDO
[764]	443	c$OMP END DO NOWAIT
[630]	444
	445	c 47. champs de vents aux pole nord
	446	c ------------------------------
	447	c U = 1 / pi * integrale [ v * cos(long) * d long ]
	448	c V = 1 / pi * integrale [ v * sin(long) * d long ]
	449
[774]	450	if (is_north_pole) then
[764]	451	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	452	DO l=1,llm
	453
	454	z1(1) =(rlonu(1)-rlonu(iim)+2.pi)pvcov(1,1,l)/cv(1,1)
	455	c z1bis(1)=(rlonu(1)-rlonu(iim)+2.pi)pdvcov(1,1,l)/cv(1,1)
	456	DO i=2,iim
	457	z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,1,l)/cv(i,1)
	458	c z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,1,l)/cv(i,1)
	459	ENDDO
	460
	461	DO i=1,iim
	462	zcos(i) = COS(rlonv(i))*z1(i)
	463	c zcosbis(i)= COS(rlonv(i))*z1bis(i)
	464	zsin(i) = SIN(rlonv(i))*z1(i)
	465	c zsinbis(i)= SIN(rlonv(i))*z1bis(i)
	466	ENDDO
	467
	468	zufi(1,l) = SSUM(iim,zcos,1)/pi
	469	c pcvgu(1,l) = SSUM(iim,zcosbis,1)/pi
	470	zvfi(1,l) = SSUM(iim,zsin,1)/pi
	471	c pcvgv(1,l) = SSUM(iim,zsinbis,1)/pi
	472
	473	ENDDO
[764]	474	c$OMP END DO NOWAIT
[630]	475	endif
	476
	477
	478	c 48. champs de vents aux pole sud:
	479	c ---------------------------------
	480	c U = 1 / pi * integrale [ v * cos(long) * d long ]
	481	c V = 1 / pi * integrale [ v * sin(long) * d long ]
	482
[774]	483	if (is_south_pole) then
[764]	484	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	485	DO l=1,llm
	486
	487	z1(1) =(rlonu(1)-rlonu(iim)+2.pi)pvcov(1,jjm,l)/cv(1,jjm)
	488	c z1bis(1)=(rlonu(1)-rlonu(iim)+2.pi)pdvcov(1,jjm,l)/cv(1,jjm)
	489	DO i=2,iim
	490	z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,jjm,l)/cv(i,jjm)
	491	c z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,jjm,l)/cv(i,jjm)
	492	ENDDO
	493
	494	DO i=1,iim
	495	zcos(i) = COS(rlonv(i))*z1(i)
	496	c zcosbis(i) = COS(rlonv(i))*z1bis(i)
	497	zsin(i) = SIN(rlonv(i))*z1(i)
	498	c zsinbis(i) = SIN(rlonv(i))*z1bis(i)
	499	ENDDO
	500
	501	zufi(klon,l) = SSUM(iim,zcos,1)/pi
	502	c pcvgu(klon,l) = SSUM(iim,zcosbis,1)/pi
	503	zvfi(klon,l) = SSUM(iim,zsin,1)/pi
	504	c pcvgv(klon,l) = SSUM(iim,zsinbis,1)/pi
	505
	506	ENDDO
[764]	507	c$OMP END DO NOWAIT
[630]	508	endif
	509
	510
[774]	511	IF (is_sequential) THEN
[764]	512	c
	513	cIM calcul PV a teta=350, 380, 405K
	514	CALL PVtheta(ngridmx,llm,pucov,pvcov,pteta,
	515	$ ztfi,zplay,zplev,
	516	$ ntetaSTD,rtetaSTD,PVteta)
	517	c
	518	ENDIF
[960]	519
	520	c On change de grille, dynamique vers physiq, pour le flux de masse verticale
[630]	521	CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,flxw,flxwfi)
	522
	523	c-----------------------------------------------------------------------
	524	c Appel de la physique:
	525	c ---------------------
	526
[764]	527	cc$OMP PARALLEL DEFAULT(NONE)
	528	cc$OMP+ PRIVATE(i,l,offset,iq)
[1114]	529	cc$OMP+ SHARED(klon_omp_nb,nqtot,klon_omp_begin,
[764]	530	cc$OMP+ debut,lafin,rdayvrai,heure,dtphys,zplev,zplay,
	531	cc$OMP+ zphi,zphis,presnivs,clesphy0,zufi,zvfi,ztfi,
	532	cc$OMP+ zqfi,pvervel,zdufi,zdvfi,zdtfi,zdqfi,zdpsrf)
[630]	533
[764]	534	c PRIVATE(zplev_omp,zplay_omp,zphi_omp,zphis_omp,
	535	c c$OMP+ presnivs_omp,zufi_omp,zvfi_omp,ztfi_omp,
	536	c c$OMP+ zqfi_omp,pvervel_omp,zdufi_omp,zdvfi_omp,
	537	c c$OMP+ zdtfi_omp,zdqfi_omp,zdpsrf_omp)
	538
	539	c$OMP BARRIER
	540	if (first_omp) then
[774]	541	klon=klon_omp
[764]	542
	543	allocate(zplev_omp(klon,llm+1))
	544	allocate(zplay_omp(klon,llm))
	545	allocate(zphi_omp(klon,llm))
	546	allocate(zphis_omp(klon))
	547	allocate(presnivs_omp(llm))
	548	allocate(zufi_omp(klon,llm))
	549	allocate(zvfi_omp(klon,llm))
	550	allocate(ztfi_omp(klon,llm))
[1114]	551	allocate(zqfi_omp(klon,llm,nqtot))
[960]	552	c allocate(pvervel_omp(klon,llm))
[764]	553	allocate(zdufi_omp(klon,llm))
	554	allocate(zdvfi_omp(klon,llm))
	555	allocate(zdtfi_omp(klon,llm))
[1114]	556	allocate(zdqfi_omp(klon,llm,nqtot))
[764]	557	allocate(zdpsrf_omp(klon))
[985]	558	allocate(flxwfi_omp(klon,llm))
[764]	559	first_omp=.false.
	560	endif
	561
	562
[774]	563	klon=klon_omp
	564	offset=klon_omp_begin-1
[764]	565
	566	do l=1,llm+1
	567	do i=1,klon
	568	zplev_omp(i,l)=zplev(offset+i,l)
	569	enddo
	570	enddo
	571
	572	do l=1,llm
	573	do i=1,klon
	574	zplay_omp(i,l)=zplay(offset+i,l)
	575	enddo
	576	enddo
	577
	578	do l=1,llm
	579	do i=1,klon
	580	zphi_omp(i,l)=zphi(offset+i,l)
	581	enddo
	582	enddo
	583
	584	do i=1,klon
	585	zphis_omp(i)=zphis(offset+i)
	586	enddo
	587
	588
	589	do l=1,llm
	590	presnivs_omp(l)=presnivs(l)
	591	enddo
	592
	593	do l=1,llm
	594	do i=1,klon
	595	zufi_omp(i,l)=zufi(offset+i,l)
	596	enddo
	597	enddo
	598
	599	do l=1,llm
	600	do i=1,klon
	601	zvfi_omp(i,l)=zvfi(offset+i,l)
	602	enddo
	603	enddo
	604
	605	do l=1,llm
	606	do i=1,klon
	607	ztfi_omp(i,l)=ztfi(offset+i,l)
	608	enddo
	609	enddo
	610
[1114]	611	do iq=1,nqtot
[764]	612	do l=1,llm
	613	do i=1,klon
	614	zqfi_omp(i,l,iq)=zqfi(offset+i,l,iq)
	615	enddo
	616	enddo
	617	enddo
	618
[960]	619	c do l=1,llm
	620	c do i=1,klon
	621	c pvervel_omp(i,l)=pvervel(offset+i,l)
	622	c enddo
	623	c enddo
[764]	624
	625	do l=1,llm
	626	do i=1,klon
	627	zdufi_omp(i,l)=zdufi(offset+i,l)
	628	enddo
	629	enddo
	630
	631	do l=1,llm
	632	do i=1,klon
	633	zdvfi_omp(i,l)=zdvfi(offset+i,l)
	634	enddo
	635	enddo
	636
	637	do l=1,llm
	638	do i=1,klon
	639	zdtfi_omp(i,l)=zdtfi(offset+i,l)
	640	enddo
	641	enddo
	642
[1114]	643	do iq=1,nqtot
[764]	644	do l=1,llm
	645	do i=1,klon
	646	zdqfi_omp(i,l,iq)=zdqfi(offset+i,l,iq)
	647	enddo
	648	enddo
	649	enddo
	650
	651	do i=1,klon
	652	zdpsrf_omp(i)=zdpsrf(offset+i)
	653	enddo
[985]	654
	655	do l=1,llm
	656	do i=1,klon
	657	flxwfi_omp(i,l)=flxwfi(offset+i,l)
	658	enddo
	659	enddo
[764]	660
	661	c$OMP BARRIER
	662	cym call WriteField_phy_p('zdtfi_omp',zdtfi_omp(:,:),llm)
	663
[630]	664	CALL physiq (klon,
	665	. llm,
	666	. debut,
	667	. lafin,
	668	. rdayvrai,
	669	. heure,
	670	. dtphys,
[764]	671	. zplev_omp,
	672	. zplay_omp,
	673	. zphi_omp,
	674	. zphis_omp,
	675	. presnivs_omp,
[630]	676	. clesphy0,
[764]	677	. zufi_omp,
	678	. zvfi_omp,
	679	. ztfi_omp,
	680	. zqfi_omp,
[960]	681	c . pvervel_omp,
	682	c#ifdef INCA
[985]	683	. flxwfi_omp,
[960]	684	c#endif
[764]	685	. zdufi_omp,
	686	. zdvfi_omp,
	687	. zdtfi_omp,
	688	. zdqfi_omp,
	689	. zdpsrf_omp,
	690	cIM diagnostique PVteta, Amip2
	691	. pducov,
	692	. PVteta)
[630]	693
[764]	694	cym call WriteField_phy_p('zdtfi_omp',zdtfi_omp(:,:),llm)
	695
	696	c$OMP BARRIER
	697
	698	do l=1,llm+1
	699	do i=1,klon
	700	zplev(offset+i,l)=zplev_omp(i,l)
	701	enddo
	702	enddo
	703
	704	do l=1,llm
	705	do i=1,klon
	706	zplay(offset+i,l)=zplay_omp(i,l)
	707	enddo
	708	enddo
	709
	710	do l=1,llm
	711	do i=1,klon
	712	zphi(offset+i,l)=zphi_omp(i,l)
	713	enddo
	714	enddo
	715
	716
	717	do i=1,klon
	718	zphis(offset+i)=zphis_omp(i)
	719	enddo
	720
	721
	722	do l=1,llm
	723	presnivs(l)=presnivs_omp(l)
	724	enddo
	725
	726	do l=1,llm
	727	do i=1,klon
	728	zufi(offset+i,l)=zufi_omp(i,l)
	729	enddo
	730	enddo
	731
	732	do l=1,llm
	733	do i=1,klon
	734	zvfi(offset+i,l)=zvfi_omp(i,l)
	735	enddo
	736	enddo
	737
	738	do l=1,llm
	739	do i=1,klon
	740	ztfi(offset+i,l)=ztfi_omp(i,l)
	741	enddo
	742	enddo
	743
[1114]	744	do iq=1,nqtot
[764]	745	do l=1,llm
	746	do i=1,klon
	747	zqfi(offset+i,l,iq)=zqfi_omp(i,l,iq)
	748	enddo
	749	enddo
	750	enddo
	751
[960]	752	c do l=1,llm
	753	c do i=1,klon
	754	c pvervel(offset+i,l)=pvervel_omp(i,l)
	755	c enddo
	756	c enddo
[764]	757
	758	do l=1,llm
	759	do i=1,klon
	760	zdufi(offset+i,l)=zdufi_omp(i,l)
	761	enddo
	762	enddo
	763
	764	do l=1,llm
	765	do i=1,klon
	766	zdvfi(offset+i,l)=zdvfi_omp(i,l)
	767	enddo
	768	enddo
	769
	770	do l=1,llm
	771	do i=1,klon
	772	zdtfi(offset+i,l)=zdtfi_omp(i,l)
	773	enddo
	774	enddo
	775
[1114]	776	do iq=1,nqtot
[764]	777	do l=1,llm
	778	do i=1,klon
	779	zdqfi(offset+i,l,iq)=zdqfi_omp(i,l,iq)
	780	enddo
	781	enddo
	782	enddo
	783
	784	do i=1,klon
	785	zdpsrf(offset+i)=zdpsrf_omp(i)
	786	enddo
	787
	788
	789	cc$OMP END PARALLEL
	790	klon=klon_mpi
[630]	791	500 CONTINUE
[764]	792	c$OMP BARRIER
[630]	793
[764]	794	c$OMP MASTER
	795	cym call WriteField_phy('zdtfi',zdtfi(:,:),llm)
[630]	796	call stop_timer(timer_physic)
[764]	797	c$OMP END MASTER
[1000]	798
	799	IF (using_mpi) THEN
	800
[630]	801	if (MPI_rank>0) then
[764]	802
	803	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	804	DO l=1,llm
	805	du_send(1:iim,l)=zdufi(1:iim,l)
	806	dv_send(1:iim,l)=zdvfi(1:iim,l)
	807	ENDDO
	808	c$OMP END DO NOWAIT
	809
	810	c$OMP BARRIER
[1000]	811	#ifdef CPP_MPI
[764]	812	c$OMP MASTER
[985]	813	!$OMP CRITICAL (MPI)
[630]	814	call MPI_ISSEND(du_send,iim*llm,MPI_REAL8,MPI_Rank-1,401,
[764]	815	& COMM_LMDZ,Req(1),ierr)
[630]	816	call MPI_ISSEND(dv_send,iim*llm,MPI_REAL8,MPI_Rank-1,402,
[764]	817	& COMM_LMDZ,Req(2),ierr)
[985]	818	!$OMP END CRITICAL (MPI)
[764]	819	c$OMP END MASTER
[1000]	820	#endif
[764]	821	c$OMP BARRIER
[630]	822
	823	endif
	824
	825	if (MPI_rank<MPI_Size-1) then
[764]	826	c$OMP BARRIER
[1000]	827	#ifdef CPP_MPI
[764]	828	c$OMP MASTER
[985]	829	!$OMP CRITICAL (MPI)
[630]	830	call MPI_IRECV(du_recv,iim*llm,MPI_REAL8,MPI_Rank+1,401,
[764]	831	& COMM_LMDZ,Req(3),ierr)
[630]	832	call MPI_IRECV(dv_recv,iim*llm,MPI_REAL8,MPI_Rank+1,402,
[764]	833	& COMM_LMDZ,Req(4),ierr)
[985]	834	!$OMP END CRITICAL (MPI)
[764]	835	c$OMP END MASTER
[1000]	836	#endif
[630]	837	endif
[764]	838
	839	c$OMP BARRIER
[1000]	840
	841
	842	#ifdef CPP_MPI
[764]	843	c$OMP MASTER
[985]	844	!$OMP CRITICAL (MPI)
[630]	845	if (MPI_rank>0 .and. MPI_rank< MPI_Size-1) then
	846	call MPI_WAITALL(4,Req(1),Status,ierr)
	847	else if (MPI_rank>0) then
	848	call MPI_WAITALL(2,Req(1),Status,ierr)
	849	else if (MPI_rank <MPI_Size-1) then
	850	call MPI_WAITALL(2,Req(3),Status,ierr)
	851	endif
[985]	852	!$OMP END CRITICAL (MPI)
[764]	853	c$OMP END MASTER
[1000]	854	#endif
	855
[764]	856	c$OMP BARRIER
	857
[1000]	858	ENDIF ! using_mpi
	859
	860
[764]	861	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	862	DO l=1,llm
	863
	864	zdufi2(1:klon,l)=zdufi(1:klon,l)
	865	zdufi2(klon+1:klon+iim,l)=du_recv(1:iim,l)
	866
	867	zdvfi2(1:klon,l)=zdvfi(1:klon,l)
	868	zdvfi2(klon+1:klon+iim,l)=dv_recv(1:iim,l)
	869
[774]	870	pdhfi(:,jj_begin,l)=0
	871	pdqfi(:,jj_begin,l,:)=0
	872	pdufi(:,jj_begin,l)=0
	873	pdvfi(:,jj_begin,l)=0
[764]	874
[774]	875	if (.not. is_south_pole) then
	876	pdhfi(:,jj_end,l)=0
	877	pdqfi(:,jj_end,l,:)=0
	878	pdufi(:,jj_end,l)=0
	879	pdvfi(:,jj_end,l)=0
[764]	880	endif
[630]	881
[764]	882	ENDDO
	883	c$OMP END DO NOWAIT
[630]	884
[764]	885	c$OMP MASTER
[774]	886	pdpsfi(:,jj_begin)=0
	887	if (.not. is_south_pole) then
	888	pdpsfi(:,jj_end)=0
[630]	889	endif
[764]	890	c$OMP END MASTER
[630]	891	c-----------------------------------------------------------------------
	892	c transformation des tendances physiques en tendances dynamiques:
	893	c ---------------------------------------------------------------
	894
	895	c tendance sur la pression :
	896	c -----------------------------------
[764]	897	c$OMP MASTER
[630]	898	CALL gr_fi_dyn_p(1,klon,iip1,jjp1,zdpsrf,pdpsfi)
[764]	899	c$OMP END MASTER
[630]	900	c
	901	c 62. enthalpie potentielle
	902	c ---------------------
	903
	904	kstart=1
	905	kend=klon
	906
[774]	907	if (is_north_pole) kstart=2
	908	if (is_south_pole) kend=klon-1
[630]	909
[764]	910	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	911	DO l=1,llm
	912
[764]	913	!!cdir NODEP
[630]	914	do ig0=kstart,kend
[774]	915	i=index_i(ig0)
	916	j=index_j(ig0)
[630]	917	pdhfi(i,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l)
	918	if (i==1) pdhfi(iip1,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l)
	919	enddo
	920
[774]	921	if (is_north_pole) then
[630]	922	DO i=1,iip1
	923	pdhfi(i,1,l) = cpp * zdtfi(1,l) / ppk(i, 1 ,l)
	924	enddo
	925	endif
	926
[774]	927	if (is_south_pole) then
[630]	928	DO i=1,iip1
	929	pdhfi(i,jjp1,l) = cpp * zdtfi(klon,l)/ ppk(i,jjp1,l)
	930	ENDDO
	931	endif
	932	ENDDO
[764]	933	c$OMP END DO NOWAIT
[630]	934
	935	c 62. humidite specifique
	936	c ---------------------
	937
[1114]	938	DO iq=1,nqtot
[764]	939	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	940	DO l=1,llm
[764]	941	!!cdir NODEP
[630]	942	do ig0=kstart,kend
[774]	943	i=index_i(ig0)
	944	j=index_j(ig0)
[630]	945	pdqfi(i,j,l,iq) = zdqfi(ig0,l,iq)
	946	if (i==1) pdqfi(iip1,j,l,iq) = zdqfi(ig0,l,iq)
	947	enddo
	948
[774]	949	if (is_north_pole) then
[630]	950	do i=1,iip1
	951	pdqfi(i,1,l,iq) = zdqfi(1,l,iq)
	952	enddo
	953	endif
	954
[774]	955	if (is_south_pole) then
[630]	956	do i=1,iip1
	957	pdqfi(i,jjp1,l,iq) = zdqfi(klon,l,iq)
	958	enddo
	959	endif
	960
	961	ENDDO
[764]	962	c$OMP END DO NOWAIT
[630]	963	ENDDO
	964
	965	c 63. traceurs
	966	c ------------
	967	C initialisation des tendances
[764]	968
	969	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	970	DO l=1,llm
	971	pdqfi(:,:,l,:)=0.
	972	ENDDO
	973	c$OMP END DO NOWAIT
	974
[630]	975	C
	976
[1114]	977	DO iq=1,nqtot
[630]	978	iiq=niadv(iq)
[764]	979	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	980	DO l=1,llm
	981
[764]	982	!!cdir NODEP
[630]	983	DO ig0=kstart,kend
[774]	984	i=index_i(ig0)
	985	j=index_j(ig0)
[630]	986	pdqfi(i,j,l,iiq) = zdqfi(ig0,l,iq)
	987	if (i==1) pdqfi(iip1,j,l,iiq) = zdqfi(ig0,l,iq)
	988	ENDDO
	989
[774]	990	IF (is_north_pole) then
[630]	991	DO i=1,iip1
	992	pdqfi(i,1,l,iiq) = zdqfi(1,l,iq)
	993	ENDDO
	994	ENDIF
	995
[774]	996	IF (is_south_pole) then
[630]	997	DO i=1,iip1
	998	pdqfi(i,jjp1,l,iiq) = zdqfi(klon,l,iq)
	999	ENDDO
	1000	ENDIF
	1001
	1002	ENDDO
[764]	1003	c$OMP END DO NOWAIT
[630]	1004	ENDDO
	1005
	1006	c 65. champ u:
	1007	c ------------
[764]	1008	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	1009	DO l=1,llm
[764]	1010	!!cdir NODEP
[630]	1011	do ig0=kstart,kend
[774]	1012	i=index_i(ig0)
	1013	j=index_j(ig0)
[630]	1014
	1015	if (i/=iim) then
	1016	pdufi(i,j,l)=0.5(zdufi2(ig0,l)+zdufi2(ig0+1,l))cu(i,j)
	1017	endif
	1018
	1019	if (i==1) then
	1020	pdufi(iim,j,l)=0.5*( zdufi2(ig0,l)
	1021	$ + zdufi2(ig0+iim-1,l))*cu(iim,j)
[764]	1022	pdufi(iip1,j,l)=0.5(zdufi2(ig0,l)+zdufi2(ig0+1,l))cu(i,j)
[630]	1023	endif
	1024
	1025	enddo
	1026
[774]	1027	if (is_north_pole) then
[630]	1028	DO i=1,iip1
	1029	pdufi(i,1,l) = 0.
	1030	ENDDO
	1031	endif
	1032
[774]	1033	if (is_south_pole) then
[630]	1034	DO i=1,iip1
	1035	pdufi(i,jjp1,l) = 0.
	1036	ENDDO
	1037	endif
	1038
	1039	ENDDO
[764]	1040	c$OMP END DO NOWAIT
[630]	1041
	1042	c 67. champ v:
	1043	c ------------
	1044
	1045	kstart=1
	1046	kend=klon
	1047
[774]	1048	if (is_north_pole) kstart=2
	1049	if (is_south_pole) kend=klon-1-iim
[630]	1050
[764]	1051	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	1052	DO l=1,llm
[764]	1053	!!cdir NODEP
[630]	1054	do ig0=kstart,kend
[774]	1055	i=index_i(ig0)
	1056	j=index_j(ig0)
[630]	1057	pdvfi(i,j,l)=0.5(zdvfi2(ig0,l)+zdvfi2(ig0+iim,l))cv(i,j)
	1058	if (i==1) pdvfi(iip1,j,l) = 0.5*(zdvfi2(ig0,l)+
	1059	$ zdvfi2(ig0+iim,l))
	1060	$ *cv(i,j)
	1061	enddo
	1062
	1063	ENDDO
[764]	1064	c$OMP END DO NOWAIT
[630]	1065
	1066
	1067	c 68. champ v pres des poles:
	1068	c ---------------------------
	1069	c v = U * cos(long) + V * SIN(long)
	1070
[774]	1071	if (is_north_pole) then
[764]	1072
	1073	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	1074	DO l=1,llm
	1075
	1076	DO i=1,iim
	1077	pdvfi(i,1,l)=
	1078	$ zdufi(1,l)COS(rlonv(i))+zdvfi(1,l)SIN(rlonv(i))
	1079
	1080	pdvfi(i,1,l)=
	1081	$ 0.5(pdvfi(i,1,l)+zdvfi(i+1,l))cv(i,1)
	1082	ENDDO
	1083
	1084	pdvfi(iip1,1,l) = pdvfi(1,1,l)
	1085
	1086	ENDDO
[764]	1087	c$OMP END DO NOWAIT
[630]	1088
	1089	endif
	1090
[774]	1091	if (is_south_pole) then
[764]	1092
	1093	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	1094	DO l=1,llm
[630]	1095
	1096	DO i=1,iim
	1097	pdvfi(i,jjm,l)=zdufi(klon,l)*COS(rlonv(i))
	1098	$ +zdvfi(klon,l)*SIN(rlonv(i))
	1099
	1100	pdvfi(i,jjm,l)=
	1101	$ 0.5(pdvfi(i,jjm,l)+zdvfi(klon-iip1+i,l))cv(i,jjm)
	1102	ENDDO
	1103
	1104	pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l)
	1105
	1106	ENDDO
[764]	1107	c$OMP END DO NOWAIT
[630]	1108
	1109	endif
	1110	c-----------------------------------------------------------------------
	1111
	1112	700 CONTINUE
	1113
	1114	firstcal = .FALSE.
	1115
	1116	RETURN
	1117	END

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

Download in other formats:

Original Format