Context Navigation

calfis_p.F @ 988

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

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

Download in other formats:

Original Format