Context Navigation

calfis_p.F @ 1146

Last change on this file since 1146 was 1146, checked in by Laurent Fairhead, 15 years ago
Réintegration dans le tronc des modifications issues de la branche LMDZ-dev comprises entre la révision 1074 et 1145 Validation: une simulation de 1 jour en séquentiel sur PC donne les mêmes résultats entre la trunk et la dev LF
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
[1146]	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
[1146]	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
[1146]	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)
[1146]	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)
[1146]	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)
[1146]	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))
[1146]	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))
[1146]	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
[1146]	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
[630]	371	CALL gr_dyn_fi_p(1,iip1,jjp1,klon,pphis,zphis)
[1146]	372
[764]	373	c$OMP BARRIER
[630]	374
[764]	375	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	376	DO l=1,llm
	377	DO ig=1,klon
	378	zphi(ig,l)=zphi(ig,l)-zphis(ig)
	379	ENDDO
	380	ENDDO
[960]	381	c$OMP END DO NOWAIT
	382
[630]	383	c .... Calcul de la vitesse verticale ( en Pams ou Kg/s ) ....
[960]	384	c JG : ancien calcule de omega utilise dans physiq.F. Maintenant le flux
	385	c de masse est calclue dans advtrac_p.F
[630]	386	c
[960]	387	cc$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	388	c DO l=1,llm
	389	c do ig0=1,klon
	390	c i=index_i(ig0)
	391	c j=index_j(ig0)
	392	c pvervel(ig0,l) = pw(i,j,l)g unsaire(i,j)
	393	c enddo
	394	c if (is_north_pole) pvervel(1,l)=pw(1,1,l)*g /apoln
	395	c if (is_south_pole) pvervel(klon,l)=pw(1,jjp1,l)*g/apols
	396	c ENDDO
	397	cc$OMP END DO NOWAIT
[630]	398
	399	c
	400	c 45. champ u:
	401	c ------------
	402
	403	kstart=1
	404	kend=klon
	405
[774]	406	if (is_north_pole) kstart=2
	407	if (is_south_pole) kend=klon-1
[630]	408
[764]	409	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	410	DO l=1,llm
	411	do ig0=kstart,kend
[774]	412	i=index_i(ig0)
	413	j=index_j(ig0)
[630]	414	if (i==1) then
	415	zufi(ig0,l)= 0.5 *( pucov(iim,j,l)/cu(iim,j)
	416	$ + pucov(1,j,l)/cu(1,j) )
	417	c pcvgu(ig0,l)= 0.5*( pducov(iim,j,l)/cu(iim,j)
	418	c $ + pducov(1,j,l)/cu(1,j) )
	419	else
	420	zufi(ig0,l)= 0.5*( pucov(i-1,j,l)/cu(i-1,j)
	421	$ + pucov(i,j,l)/cu(i,j) )
	422	c pcvgu(ig0,l)= 0.5*( pducov(i-1,j,l)/cu(i-1,j)
	423	c $ + pducov(i,j,l)/cu(i,j) )
	424	endif
	425	enddo
	426	ENDDO
[764]	427	c$OMP END DO NOWAIT
[630]	428	c 46.champ v:
	429	c -----------
[764]	430	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	431	DO l=1,llm
	432	DO ig0=kstart,kend
[774]	433	i=index_i(ig0)
	434	j=index_j(ig0)
[630]	435	zvfi(ig0,l)= 0.5 *( pvcov(i,j-1,l)/cv(i,j-1)
	436	$ + pvcov(i,j,l)/cv(i,j) )
	437
	438	c pcvgv(ig0+i,l)= 0.5 * ( pdvcov(i,j-1,l)/cv(i,j-1)
	439	c $ + pdvcov(i,j,l)/cv(i,j) )
	440	ENDDO
	441	ENDDO
[764]	442	c$OMP END DO NOWAIT
[630]	443
	444	c 47. champs de vents aux pole nord
	445	c ------------------------------
	446	c U = 1 / pi * integrale [ v * cos(long) * d long ]
	447	c V = 1 / pi * integrale [ v * sin(long) * d long ]
	448
[774]	449	if (is_north_pole) then
[764]	450	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	451	DO l=1,llm
	452
	453	z1(1) =(rlonu(1)-rlonu(iim)+2.pi)pvcov(1,1,l)/cv(1,1)
	454	c z1bis(1)=(rlonu(1)-rlonu(iim)+2.pi)pdvcov(1,1,l)/cv(1,1)
	455	DO i=2,iim
	456	z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,1,l)/cv(i,1)
	457	c z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,1,l)/cv(i,1)
	458	ENDDO
	459
	460	DO i=1,iim
	461	zcos(i) = COS(rlonv(i))*z1(i)
	462	c zcosbis(i)= COS(rlonv(i))*z1bis(i)
	463	zsin(i) = SIN(rlonv(i))*z1(i)
	464	c zsinbis(i)= SIN(rlonv(i))*z1bis(i)
	465	ENDDO
	466
	467	zufi(1,l) = SSUM(iim,zcos,1)/pi
	468	c pcvgu(1,l) = SSUM(iim,zcosbis,1)/pi
	469	zvfi(1,l) = SSUM(iim,zsin,1)/pi
	470	c pcvgv(1,l) = SSUM(iim,zsinbis,1)/pi
	471
	472	ENDDO
[764]	473	c$OMP END DO NOWAIT
[630]	474	endif
	475
	476
	477	c 48. champs de vents aux pole sud:
	478	c ---------------------------------
	479	c U = 1 / pi * integrale [ v * cos(long) * d long ]
	480	c V = 1 / pi * integrale [ v * sin(long) * d long ]
	481
[774]	482	if (is_south_pole) then
[764]	483	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	484	DO l=1,llm
	485
	486	z1(1) =(rlonu(1)-rlonu(iim)+2.pi)pvcov(1,jjm,l)/cv(1,jjm)
	487	c z1bis(1)=(rlonu(1)-rlonu(iim)+2.pi)pdvcov(1,jjm,l)/cv(1,jjm)
	488	DO i=2,iim
	489	z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,jjm,l)/cv(i,jjm)
	490	c z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,jjm,l)/cv(i,jjm)
	491	ENDDO
	492
	493	DO i=1,iim
	494	zcos(i) = COS(rlonv(i))*z1(i)
	495	c zcosbis(i) = COS(rlonv(i))*z1bis(i)
	496	zsin(i) = SIN(rlonv(i))*z1(i)
	497	c zsinbis(i) = SIN(rlonv(i))*z1bis(i)
	498	ENDDO
	499
	500	zufi(klon,l) = SSUM(iim,zcos,1)/pi
	501	c pcvgu(klon,l) = SSUM(iim,zcosbis,1)/pi
	502	zvfi(klon,l) = SSUM(iim,zsin,1)/pi
	503	c pcvgv(klon,l) = SSUM(iim,zsinbis,1)/pi
	504
	505	ENDDO
[764]	506	c$OMP END DO NOWAIT
[630]	507	endif
	508
	509
[774]	510	IF (is_sequential) THEN
[764]	511	c
	512	cIM calcul PV a teta=350, 380, 405K
	513	CALL PVtheta(ngridmx,llm,pucov,pvcov,pteta,
	514	$ ztfi,zplay,zplev,
	515	$ ntetaSTD,rtetaSTD,PVteta)
	516	c
	517	ENDIF
[960]	518
	519	c On change de grille, dynamique vers physiq, pour le flux de masse verticale
[630]	520	CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,flxw,flxwfi)
	521
	522	c-----------------------------------------------------------------------
	523	c Appel de la physique:
	524	c ---------------------
	525
[764]	526	cc$OMP PARALLEL DEFAULT(NONE)
	527	cc$OMP+ PRIVATE(i,l,offset,iq)
[1146]	528	cc$OMP+ SHARED(klon_omp_nb,nqtot,klon_omp_begin,
[764]	529	cc$OMP+ debut,lafin,rdayvrai,heure,dtphys,zplev,zplay,
	530	cc$OMP+ zphi,zphis,presnivs,clesphy0,zufi,zvfi,ztfi,
	531	cc$OMP+ zqfi,pvervel,zdufi,zdvfi,zdtfi,zdqfi,zdpsrf)
[630]	532
[764]	533	c PRIVATE(zplev_omp,zplay_omp,zphi_omp,zphis_omp,
	534	c c$OMP+ presnivs_omp,zufi_omp,zvfi_omp,ztfi_omp,
	535	c c$OMP+ zqfi_omp,pvervel_omp,zdufi_omp,zdvfi_omp,
	536	c c$OMP+ zdtfi_omp,zdqfi_omp,zdpsrf_omp)
	537
	538	c$OMP BARRIER
	539	if (first_omp) then
[774]	540	klon=klon_omp
[764]	541
	542	allocate(zplev_omp(klon,llm+1))
	543	allocate(zplay_omp(klon,llm))
	544	allocate(zphi_omp(klon,llm))
	545	allocate(zphis_omp(klon))
	546	allocate(presnivs_omp(llm))
	547	allocate(zufi_omp(klon,llm))
	548	allocate(zvfi_omp(klon,llm))
	549	allocate(ztfi_omp(klon,llm))
[1146]	550	allocate(zqfi_omp(klon,llm,nqtot))
[960]	551	c allocate(pvervel_omp(klon,llm))
[764]	552	allocate(zdufi_omp(klon,llm))
	553	allocate(zdvfi_omp(klon,llm))
	554	allocate(zdtfi_omp(klon,llm))
[1146]	555	allocate(zdqfi_omp(klon,llm,nqtot))
[764]	556	allocate(zdpsrf_omp(klon))
[985]	557	allocate(flxwfi_omp(klon,llm))
[764]	558	first_omp=.false.
	559	endif
	560
	561
[774]	562	klon=klon_omp
	563	offset=klon_omp_begin-1
[764]	564
	565	do l=1,llm+1
	566	do i=1,klon
	567	zplev_omp(i,l)=zplev(offset+i,l)
	568	enddo
	569	enddo
	570
	571	do l=1,llm
	572	do i=1,klon
	573	zplay_omp(i,l)=zplay(offset+i,l)
	574	enddo
	575	enddo
	576
	577	do l=1,llm
	578	do i=1,klon
	579	zphi_omp(i,l)=zphi(offset+i,l)
	580	enddo
	581	enddo
	582
	583	do i=1,klon
	584	zphis_omp(i)=zphis(offset+i)
	585	enddo
	586
	587
	588	do l=1,llm
	589	presnivs_omp(l)=presnivs(l)
	590	enddo
	591
	592	do l=1,llm
	593	do i=1,klon
	594	zufi_omp(i,l)=zufi(offset+i,l)
	595	enddo
	596	enddo
	597
	598	do l=1,llm
	599	do i=1,klon
	600	zvfi_omp(i,l)=zvfi(offset+i,l)
	601	enddo
	602	enddo
	603
	604	do l=1,llm
	605	do i=1,klon
	606	ztfi_omp(i,l)=ztfi(offset+i,l)
	607	enddo
	608	enddo
	609
[1146]	610	do iq=1,nqtot
[764]	611	do l=1,llm
	612	do i=1,klon
	613	zqfi_omp(i,l,iq)=zqfi(offset+i,l,iq)
	614	enddo
	615	enddo
	616	enddo
	617
[960]	618	c do l=1,llm
	619	c do i=1,klon
	620	c pvervel_omp(i,l)=pvervel(offset+i,l)
	621	c enddo
	622	c enddo
[764]	623
	624	do l=1,llm
	625	do i=1,klon
	626	zdufi_omp(i,l)=zdufi(offset+i,l)
	627	enddo
	628	enddo
	629
	630	do l=1,llm
	631	do i=1,klon
	632	zdvfi_omp(i,l)=zdvfi(offset+i,l)
	633	enddo
	634	enddo
	635
	636	do l=1,llm
	637	do i=1,klon
	638	zdtfi_omp(i,l)=zdtfi(offset+i,l)
	639	enddo
	640	enddo
	641
[1146]	642	do iq=1,nqtot
[764]	643	do l=1,llm
	644	do i=1,klon
	645	zdqfi_omp(i,l,iq)=zdqfi(offset+i,l,iq)
	646	enddo
	647	enddo
	648	enddo
	649
	650	do i=1,klon
	651	zdpsrf_omp(i)=zdpsrf(offset+i)
	652	enddo
[985]	653
	654	do l=1,llm
	655	do i=1,klon
	656	flxwfi_omp(i,l)=flxwfi(offset+i,l)
	657	enddo
	658	enddo
[764]	659
	660	c$OMP BARRIER
	661	cym call WriteField_phy_p('zdtfi_omp',zdtfi_omp(:,:),llm)
	662
[630]	663	CALL physiq (klon,
	664	. llm,
	665	. debut,
	666	. lafin,
	667	. rdayvrai,
	668	. heure,
	669	. dtphys,
[764]	670	. zplev_omp,
	671	. zplay_omp,
	672	. zphi_omp,
	673	. zphis_omp,
	674	. presnivs_omp,
[630]	675	. clesphy0,
[764]	676	. zufi_omp,
	677	. zvfi_omp,
	678	. ztfi_omp,
	679	. zqfi_omp,
[960]	680	c . pvervel_omp,
	681	c#ifdef INCA
[985]	682	. flxwfi_omp,
[960]	683	c#endif
[764]	684	. zdufi_omp,
	685	. zdvfi_omp,
	686	. zdtfi_omp,
	687	. zdqfi_omp,
	688	. zdpsrf_omp,
	689	cIM diagnostique PVteta, Amip2
	690	. pducov,
	691	. PVteta)
[630]	692
[764]	693	cym call WriteField_phy_p('zdtfi_omp',zdtfi_omp(:,:),llm)
	694
	695	c$OMP BARRIER
	696
	697	do l=1,llm+1
	698	do i=1,klon
	699	zplev(offset+i,l)=zplev_omp(i,l)
	700	enddo
	701	enddo
	702
	703	do l=1,llm
	704	do i=1,klon
	705	zplay(offset+i,l)=zplay_omp(i,l)
	706	enddo
	707	enddo
	708
	709	do l=1,llm
	710	do i=1,klon
	711	zphi(offset+i,l)=zphi_omp(i,l)
	712	enddo
	713	enddo
	714
	715
	716	do i=1,klon
	717	zphis(offset+i)=zphis_omp(i)
	718	enddo
	719
	720
	721	do l=1,llm
	722	presnivs(l)=presnivs_omp(l)
	723	enddo
	724
	725	do l=1,llm
	726	do i=1,klon
	727	zufi(offset+i,l)=zufi_omp(i,l)
	728	enddo
	729	enddo
	730
	731	do l=1,llm
	732	do i=1,klon
	733	zvfi(offset+i,l)=zvfi_omp(i,l)
	734	enddo
	735	enddo
	736
	737	do l=1,llm
	738	do i=1,klon
	739	ztfi(offset+i,l)=ztfi_omp(i,l)
	740	enddo
	741	enddo
	742
[1146]	743	do iq=1,nqtot
[764]	744	do l=1,llm
	745	do i=1,klon
	746	zqfi(offset+i,l,iq)=zqfi_omp(i,l,iq)
	747	enddo
	748	enddo
	749	enddo
	750
[960]	751	c do l=1,llm
	752	c do i=1,klon
	753	c pvervel(offset+i,l)=pvervel_omp(i,l)
	754	c enddo
	755	c enddo
[764]	756
	757	do l=1,llm
	758	do i=1,klon
	759	zdufi(offset+i,l)=zdufi_omp(i,l)
	760	enddo
	761	enddo
	762
	763	do l=1,llm
	764	do i=1,klon
	765	zdvfi(offset+i,l)=zdvfi_omp(i,l)
	766	enddo
	767	enddo
	768
	769	do l=1,llm
	770	do i=1,klon
	771	zdtfi(offset+i,l)=zdtfi_omp(i,l)
	772	enddo
	773	enddo
	774
[1146]	775	do iq=1,nqtot
[764]	776	do l=1,llm
	777	do i=1,klon
	778	zdqfi(offset+i,l,iq)=zdqfi_omp(i,l,iq)
	779	enddo
	780	enddo
	781	enddo
	782
	783	do i=1,klon
	784	zdpsrf(offset+i)=zdpsrf_omp(i)
	785	enddo
	786
	787
	788	cc$OMP END PARALLEL
	789	klon=klon_mpi
[630]	790	500 CONTINUE
[764]	791	c$OMP BARRIER
[630]	792
[764]	793	c$OMP MASTER
	794	cym call WriteField_phy('zdtfi',zdtfi(:,:),llm)
[630]	795	call stop_timer(timer_physic)
[764]	796	c$OMP END MASTER
[1000]	797
	798	IF (using_mpi) THEN
	799
[630]	800	if (MPI_rank>0) then
[764]	801
	802	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	803	DO l=1,llm
	804	du_send(1:iim,l)=zdufi(1:iim,l)
	805	dv_send(1:iim,l)=zdvfi(1:iim,l)
	806	ENDDO
	807	c$OMP END DO NOWAIT
	808
	809	c$OMP BARRIER
[1000]	810	#ifdef CPP_MPI
[764]	811	c$OMP MASTER
[985]	812	!$OMP CRITICAL (MPI)
[630]	813	call MPI_ISSEND(du_send,iim*llm,MPI_REAL8,MPI_Rank-1,401,
[764]	814	& COMM_LMDZ,Req(1),ierr)
[630]	815	call MPI_ISSEND(dv_send,iim*llm,MPI_REAL8,MPI_Rank-1,402,
[764]	816	& COMM_LMDZ,Req(2),ierr)
[985]	817	!$OMP END CRITICAL (MPI)
[764]	818	c$OMP END MASTER
[1000]	819	#endif
[764]	820	c$OMP BARRIER
[630]	821
	822	endif
	823
	824	if (MPI_rank<MPI_Size-1) then
[764]	825	c$OMP BARRIER
[1000]	826	#ifdef CPP_MPI
[764]	827	c$OMP MASTER
[985]	828	!$OMP CRITICAL (MPI)
[630]	829	call MPI_IRECV(du_recv,iim*llm,MPI_REAL8,MPI_Rank+1,401,
[764]	830	& COMM_LMDZ,Req(3),ierr)
[630]	831	call MPI_IRECV(dv_recv,iim*llm,MPI_REAL8,MPI_Rank+1,402,
[764]	832	& COMM_LMDZ,Req(4),ierr)
[985]	833	!$OMP END CRITICAL (MPI)
[764]	834	c$OMP END MASTER
[1000]	835	#endif
[630]	836	endif
[764]	837
	838	c$OMP BARRIER
[1000]	839
	840
	841	#ifdef CPP_MPI
[764]	842	c$OMP MASTER
[985]	843	!$OMP CRITICAL (MPI)
[630]	844	if (MPI_rank>0 .and. MPI_rank< MPI_Size-1) then
	845	call MPI_WAITALL(4,Req(1),Status,ierr)
	846	else if (MPI_rank>0) then
	847	call MPI_WAITALL(2,Req(1),Status,ierr)
	848	else if (MPI_rank <MPI_Size-1) then
	849	call MPI_WAITALL(2,Req(3),Status,ierr)
	850	endif
[985]	851	!$OMP END CRITICAL (MPI)
[764]	852	c$OMP END MASTER
[1000]	853	#endif
	854
[764]	855	c$OMP BARRIER
	856
[1000]	857	ENDIF ! using_mpi
	858
	859
[764]	860	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	861	DO l=1,llm
	862
	863	zdufi2(1:klon,l)=zdufi(1:klon,l)
	864	zdufi2(klon+1:klon+iim,l)=du_recv(1:iim,l)
	865
	866	zdvfi2(1:klon,l)=zdvfi(1:klon,l)
	867	zdvfi2(klon+1:klon+iim,l)=dv_recv(1:iim,l)
	868
[774]	869	pdhfi(:,jj_begin,l)=0
	870	pdqfi(:,jj_begin,l,:)=0
	871	pdufi(:,jj_begin,l)=0
	872	pdvfi(:,jj_begin,l)=0
[764]	873
[774]	874	if (.not. is_south_pole) then
	875	pdhfi(:,jj_end,l)=0
	876	pdqfi(:,jj_end,l,:)=0
	877	pdufi(:,jj_end,l)=0
	878	pdvfi(:,jj_end,l)=0
[764]	879	endif
[630]	880
[764]	881	ENDDO
	882	c$OMP END DO NOWAIT
[630]	883
[764]	884	c$OMP MASTER
[774]	885	pdpsfi(:,jj_begin)=0
	886	if (.not. is_south_pole) then
	887	pdpsfi(:,jj_end)=0
[630]	888	endif
[764]	889	c$OMP END MASTER
[630]	890	c-----------------------------------------------------------------------
	891	c transformation des tendances physiques en tendances dynamiques:
	892	c ---------------------------------------------------------------
	893
	894	c tendance sur la pression :
	895	c -----------------------------------
	896	CALL gr_fi_dyn_p(1,klon,iip1,jjp1,zdpsrf,pdpsfi)
	897	c
	898	c 62. enthalpie potentielle
	899	c ---------------------
	900
	901	kstart=1
	902	kend=klon
	903
[774]	904	if (is_north_pole) kstart=2
	905	if (is_south_pole) kend=klon-1
[630]	906
[764]	907	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	908	DO l=1,llm
	909
[764]	910	!!cdir NODEP
[630]	911	do ig0=kstart,kend
[774]	912	i=index_i(ig0)
	913	j=index_j(ig0)
[630]	914	pdhfi(i,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l)
	915	if (i==1) pdhfi(iip1,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l)
	916	enddo
	917
[774]	918	if (is_north_pole) then
[630]	919	DO i=1,iip1
	920	pdhfi(i,1,l) = cpp * zdtfi(1,l) / ppk(i, 1 ,l)
	921	enddo
	922	endif
	923
[774]	924	if (is_south_pole) then
[630]	925	DO i=1,iip1
	926	pdhfi(i,jjp1,l) = cpp * zdtfi(klon,l)/ ppk(i,jjp1,l)
	927	ENDDO
	928	endif
	929	ENDDO
[764]	930	c$OMP END DO NOWAIT
[630]	931
	932	c 62. humidite specifique
	933	c ---------------------
	934
[1146]	935	DO iq=1,nqtot
[764]	936	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	937	DO l=1,llm
[764]	938	!!cdir NODEP
[630]	939	do ig0=kstart,kend
[774]	940	i=index_i(ig0)
	941	j=index_j(ig0)
[630]	942	pdqfi(i,j,l,iq) = zdqfi(ig0,l,iq)
	943	if (i==1) pdqfi(iip1,j,l,iq) = zdqfi(ig0,l,iq)
	944	enddo
	945
[774]	946	if (is_north_pole) then
[630]	947	do i=1,iip1
	948	pdqfi(i,1,l,iq) = zdqfi(1,l,iq)
	949	enddo
	950	endif
	951
[774]	952	if (is_south_pole) then
[630]	953	do i=1,iip1
	954	pdqfi(i,jjp1,l,iq) = zdqfi(klon,l,iq)
	955	enddo
	956	endif
	957
	958	ENDDO
[764]	959	c$OMP END DO NOWAIT
[630]	960	ENDDO
	961
	962	c 63. traceurs
	963	c ------------
	964	C initialisation des tendances
[764]	965
	966	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	967	DO l=1,llm
	968	pdqfi(:,:,l,:)=0.
	969	ENDDO
	970	c$OMP END DO NOWAIT
	971
[630]	972	C
	973
[1146]	974	DO iq=1,nqtot
[630]	975	iiq=niadv(iq)
[764]	976	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	977	DO l=1,llm
	978
[764]	979	!!cdir NODEP
[630]	980	DO ig0=kstart,kend
[774]	981	i=index_i(ig0)
	982	j=index_j(ig0)
[630]	983	pdqfi(i,j,l,iiq) = zdqfi(ig0,l,iq)
	984	if (i==1) pdqfi(iip1,j,l,iiq) = zdqfi(ig0,l,iq)
	985	ENDDO
	986
[774]	987	IF (is_north_pole) then
[630]	988	DO i=1,iip1
	989	pdqfi(i,1,l,iiq) = zdqfi(1,l,iq)
	990	ENDDO
	991	ENDIF
	992
[774]	993	IF (is_south_pole) then
[630]	994	DO i=1,iip1
	995	pdqfi(i,jjp1,l,iiq) = zdqfi(klon,l,iq)
	996	ENDDO
	997	ENDIF
	998
	999	ENDDO
[764]	1000	c$OMP END DO NOWAIT
[630]	1001	ENDDO
	1002
	1003	c 65. champ u:
	1004	c ------------
[764]	1005	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	1006	DO l=1,llm
[764]	1007	!!cdir NODEP
[630]	1008	do ig0=kstart,kend
[774]	1009	i=index_i(ig0)
	1010	j=index_j(ig0)
[630]	1011
	1012	if (i/=iim) then
	1013	pdufi(i,j,l)=0.5(zdufi2(ig0,l)+zdufi2(ig0+1,l))cu(i,j)
	1014	endif
	1015
	1016	if (i==1) then
	1017	pdufi(iim,j,l)=0.5*( zdufi2(ig0,l)
	1018	$ + zdufi2(ig0+iim-1,l))*cu(iim,j)
[764]	1019	pdufi(iip1,j,l)=0.5(zdufi2(ig0,l)+zdufi2(ig0+1,l))cu(i,j)
[630]	1020	endif
	1021
	1022	enddo
	1023
[774]	1024	if (is_north_pole) then
[630]	1025	DO i=1,iip1
	1026	pdufi(i,1,l) = 0.
	1027	ENDDO
	1028	endif
	1029
[774]	1030	if (is_south_pole) then
[630]	1031	DO i=1,iip1
	1032	pdufi(i,jjp1,l) = 0.
	1033	ENDDO
	1034	endif
	1035
	1036	ENDDO
[764]	1037	c$OMP END DO NOWAIT
[630]	1038
	1039	c 67. champ v:
	1040	c ------------
	1041
	1042	kstart=1
	1043	kend=klon
	1044
[774]	1045	if (is_north_pole) kstart=2
	1046	if (is_south_pole) kend=klon-1-iim
[630]	1047
[764]	1048	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	1049	DO l=1,llm
[764]	1050	!!cdir NODEP
[630]	1051	do ig0=kstart,kend
[774]	1052	i=index_i(ig0)
	1053	j=index_j(ig0)
[630]	1054	pdvfi(i,j,l)=0.5(zdvfi2(ig0,l)+zdvfi2(ig0+iim,l))cv(i,j)
	1055	if (i==1) pdvfi(iip1,j,l) = 0.5*(zdvfi2(ig0,l)+
	1056	$ zdvfi2(ig0+iim,l))
	1057	$ *cv(i,j)
	1058	enddo
	1059
	1060	ENDDO
[764]	1061	c$OMP END DO NOWAIT
[630]	1062
	1063
	1064	c 68. champ v pres des poles:
	1065	c ---------------------------
	1066	c v = U * cos(long) + V * SIN(long)
	1067
[774]	1068	if (is_north_pole) then
[764]	1069
	1070	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
[630]	1071	DO l=1,llm
	1072
	1073	DO i=1,iim
	1074	pdvfi(i,1,l)=
	1075	$ zdufi(1,l)COS(rlonv(i))+zdvfi(1,l)SIN(rlonv(i))
	1076
	1077	pdvfi(i,1,l)=
	1078	$ 0.5(pdvfi(i,1,l)+zdvfi(i+1,l))cv(i,1)
	1079	ENDDO
	1080
	1081	pdvfi(iip1,1,l) = pdvfi(1,1,l)
	1082
	1083	ENDDO
[764]	1084	c$OMP END DO NOWAIT
[630]	1085
	1086	endif
	1087
[774]	1088	if (is_south_pole) then
[764]	1089
	1090	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
	1091	DO l=1,llm
[630]	1092
	1093	DO i=1,iim
	1094	pdvfi(i,jjm,l)=zdufi(klon,l)*COS(rlonv(i))
	1095	$ +zdvfi(klon,l)*SIN(rlonv(i))
	1096
	1097	pdvfi(i,jjm,l)=
	1098	$ 0.5(pdvfi(i,jjm,l)+zdvfi(klon-iip1+i,l))cv(i,jjm)
	1099	ENDDO
	1100
	1101	pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l)
	1102
	1103	ENDDO
[764]	1104	c$OMP END DO NOWAIT
[630]	1105
	1106	endif
	1107	c-----------------------------------------------------------------------
	1108
	1109	700 CONTINUE
	1110
	1111	firstcal = .FALSE.
	1112
	1113	RETURN
	1114	END

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

Download in other formats:

Original Format