Context Navigation

calfis_p.F @ 1

Last change on this file since 1 was 1, checked in by emillour, 14 years ago
Import initial LMDZ5
File size: 28.0 KB

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

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

Download in other formats:

Original Format