Context Navigation

calfis_p.F @ 1211

Last change on this file since 1211 was 1207, checked in by Laurent Fairhead, 16 years ago
Probleme sur la forme des commentaires laisses par emacs, openmp pourrait mal les comprendre L'argument heure disparait de la liste des arguments d'appel a calfis
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 28.0 KB

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

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

Download in other formats:

Original Format