Context Navigation

calfis_p.F @ 977

Last change on this file since 977 was 961, checked in by Laurent Fairhead, 17 years ago
Un peu de menage pour avoir moins de sorties a l'execution LF
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 27.5 KB

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

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

Download in other formats:

Original Format