Context Navigation

calfis.F @ 5097

Last change on this file since 5097 was 5090, checked in by abarral, 2 months ago
Move lmdz_netcdf_format.F90 -> lmdz_cppkeys_wrapper.F90 to handle other CPP keys Replace all (except wrapper) use of CPP_PHYS by fortran logical
Property copyright set to Name of program: LMDZ Creation date: 1984 Version: LMDZ5 License: CeCILL version 2 Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539 See the license file in the root directory Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 21.2 KB

Rev	Line
[524]	1	!
[1279]	2	! $Id: calfis.F 5090 2024-07-20 16:08:57Z snguyen $
[524]	3	!
	4	C
	5	C
[1146]	6	SUBROUTINE calfis(lafin,
[1279]	7	$ jD_cur, jH_cur,
[524]	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	$ pdufi,
	24	$ pdvfi,
	25	$ pdhfi,
	26	$ pdqfi,
	27	$ pdpsfi)
	28	c
[5090]	29	c Auteur : P. Le Van, F. Hourdin
[524]	30	c .........
[4056]	31	USE infotrac, ONLY: nqtot, tracers
[1987]	32	USE control_mod, ONLY: planet_type, nsplit_phys
[2418]	33	USE callphysiq_mod, ONLY: call_physiq
[5090]	34	USE lmdz_cppkeys_wrapper, ONLY: CPPKEY_PHYS
[2597]	35	USE comconst_mod, ONLY: cpp, daysec, dtphys, dtvr, kappa, pi
[2600]	36	USE comvert_mod, ONLY: preff, presnivs
[5090]	37
[524]	38	IMPLICIT NONE
	39	c=======================================================================
	40	c
	41	c 1. rearrangement des tableaux et transformation
	42	c variables dynamiques > variables physiques
	43	c 2. calcul des termes physiques
	44	c 3. retransformation des tendances physiques en tendances dynamiques
	45	c
	46	c remarques:
	47	c ----------
	48	c
[5090]	49	c - les vents sont donnes dans la physique par leurs composantes
[524]	50	c naturelles.
	51	c - la variable thermodynamique de la physique est une variable
[5090]	52	c intensive : T
[524]	53	c pour la dynamique on prend T * ( preff / p(l) ) **kappa
	54	c - les deux seules variables dependant de la geometrie necessaires
[5090]	55	c pour la physique sont la latitude pour le rayonnement et
	56	c l'aire de la maille quand on veut integrer une grandeur
[524]	57	c horizontalement.
[5090]	58	c - les points de la physique sont les points scalaires de la
[524]	59	c la dynamique; numerotation:
	60	c 1 pour le pole nord
	61	c (jjm-1)*iim pour l'interieur du domaine
	62	c ngridmx pour le pole sud
	63	c ---> ngridmx=2+(jjm-1)*iim
	64	c
	65	c Input :
	66	c -------
	67	c pucov covariant zonal velocity
[5090]	68	c pvcov covariant meridional velocity
[524]	69	c pteta potential temperature
	70	c pps surface pressure
	71	c pmasse masse d'air dans chaque maille
	72	c pts surface temperature (K)
	73	c callrad clef d'appel au rayonnement
	74	c
	75	c Output :
	76	c --------
	77	c pdufi tendency for the natural zonal velocity (ms-1)
[5090]	78	c pdvfi tendency for the natural meridional velocity
[524]	79	c pdhfi tendency for the potential temperature
	80	c pdtsfi tendency for the surface temperature
	81	c
	82	c pdtrad radiative tendencies \ both input
	83	c pfluxrad radiative fluxes / and output
	84	c
	85	c=======================================================================
	86	c
	87	c-----------------------------------------------------------------------
	88	c
	89	c 0. Declarations :
	90	c ------------------
	91
[2600]	92	include "dimensions.h"
	93	include "paramet.h"
[524]	94
[1146]	95	INTEGER ngridmx
[524]	96	PARAMETER( ngridmx = 2+(jjm-1)*iim - 1/jjm )
	97
[2600]	98	include "comgeom2.h"
	99	include "iniprint.h"
[524]	100
	101	c Arguments :
	102	c -----------
[1987]	103	LOGICAL,INTENT(IN) :: lafin ! .true. for the very last call to physics
	104	REAL,INTENT(IN):: jD_cur, jH_cur
	105	REAL,INTENT(IN) :: pvcov(iip1,jjm,llm) ! covariant meridional velocity
	106	REAL,INTENT(IN) :: pucov(iip1,jjp1,llm) ! covariant zonal velocity
	107	REAL,INTENT(IN) :: pteta(iip1,jjp1,llm) ! potential temperature
	108	REAL,INTENT(IN) :: pmasse(iip1,jjp1,llm) ! mass in each cell ! not used
	109	REAL,INTENT(IN) :: pq(iip1,jjp1,llm,nqtot) ! tracers
	110	REAL,INTENT(IN) :: pphis(iip1,jjp1) ! surface geopotential
	111	REAL,INTENT(IN) :: pphi(iip1,jjp1,llm) ! geopotential
[524]	112
[1987]	113	REAL,INTENT(IN) :: pdvcov(iip1,jjm,llm) ! dynamical tendency on vcov
	114	REAL,INTENT(IN) :: pducov(iip1,jjp1,llm) ! dynamical tendency on ucov
	115	REAL,INTENT(IN) :: pdteta(iip1,jjp1,llm) ! dynamical tendency on teta
	116	! NB: pdteta is used only to compute pcvgt which is in fact not used...
	117	REAL,INTENT(IN) :: pdq(iip1,jjp1,llm,nqtot) ! dynamical tendency on tracers
	118	! NB: pdq is only used to compute pcvgq which is in fact not used...
[1279]	119
[1987]	120	REAL,INTENT(IN) :: pps(iip1,jjp1) ! surface pressure (Pa)
	121	REAL,INTENT(IN) :: pp(iip1,jjp1,llmp1) ! pressure at mesh interfaces (Pa)
	122	REAL,INTENT(IN) :: ppk(iip1,jjp1,llm) ! Exner at mid-layer
[2418]	123	REAL,INTENT(IN) :: flxw(iip1,jjp1,llm) ! Vertical mass flux on lower mesh interfaces (kg/s) (on llm because flxw(:,:,llm+1)=0)
[524]	124
[1987]	125	! tendencies (in */s) from the physics
	126	REAL,INTENT(OUT) :: pdvfi(iip1,jjm,llm) ! tendency on covariant meridional wind
	127	REAL,INTENT(OUT) :: pdufi(iip1,jjp1,llm) ! tendency on covariant zonal wind
	128	REAL,INTENT(OUT) :: pdhfi(iip1,jjp1,llm) ! tendency on potential temperature (K/s)
	129	REAL,INTENT(OUT) :: pdqfi(iip1,jjp1,llm,nqtot) ! tendency on tracers
	130	REAL,INTENT(OUT) :: pdpsfi(iip1,jjp1) ! tendency on surface pressure (Pa/s)
[524]	131
	132
	133	c Local variables :
	134	c -----------------
	135
[4056]	136	INTEGER i,j,l,ig0,ig,iq,itr
[524]	137	REAL zpsrf(ngridmx)
	138	REAL zplev(ngridmx,llm+1),zplay(ngridmx,llm)
	139	REAL zphi(ngridmx,llm),zphis(ngridmx)
	140	c
[2333]	141	REAL zrot(iip1,jjm,llm) ! AdlC May 2014
[2604]	142	REAL zufi(ngridmx,llm), zvfi(ngridmx,llm)
	143	REAL zrfi(ngridmx,llm) ! relative wind vorticity
[1146]	144	REAL ztfi(ngridmx,llm),zqfi(ngridmx,llm,nqtot)
[2604]	145	REAL zpk(ngridmx,llm)
[524]	146	c
	147	REAL pcvgu(ngridmx,llm), pcvgv(ngridmx,llm)
	148	REAL pcvgt(ngridmx,llm), pcvgq(ngridmx,llm,2)
	149	c
	150	REAL zdufi(ngridmx,llm),zdvfi(ngridmx,llm)
[1146]	151	REAL zdtfi(ngridmx,llm),zdqfi(ngridmx,llm,nqtot)
[524]	152	REAL zdpsrf(ngridmx)
	153	c
[1403]	154	REAL zdufic(ngridmx,llm),zdvfic(ngridmx,llm)
	155	REAL zdtfic(ngridmx,llm),zdqfic(ngridmx,llm,nqtot)
	156	REAL jH_cur_split,zdt_split
	157	LOGICAL debut_split,lafin_split
	158	INTEGER isplit
	159
[524]	160	REAL zsin(iim),zcos(iim),z1(iim)
	161	REAL zsinbis(iim),zcosbis(iim),z1bis(iim)
	162	REAL unskap, pksurcp
[644]	163	c
[960]	164	REAL flxwfi(ngridmx,llm) ! Flux de masse verticale sur la grille physiq
[524]	165	c
[5090]	166
[524]	167	REAL SSUM
	168
[1987]	169	LOGICAL,SAVE :: firstcal=.true., debut=.true.
[1279]	170	! REAL rdayvrai
[1654]	171
[524]	172	c
	173	c-----------------------------------------------------------------------
	174	c
	175	c 1. Initialisations :
	176	c --------------------
	177	c
[1279]	178	c
	179	IF ( firstcal ) THEN
	180	debut = .TRUE.
[5082]	181	IF (ngridmx/=2+(jjm-1)*iim) THEN
[1403]	182	write(lunout,*) 'STOP dans calfis'
	183	write(lunout,*)
	184	& 'La dimension ngridmx doit etre egale a 2 + (jjm-1)*iim'
	185	write(lunout,*) ' ngridmx jjm iim '
	186	write(lunout,*) ngridmx,jjm,iim
[4464]	187	call abort_gcm("calfis", "", 1)
[1279]	188	ENDIF
[524]	189	ELSE
[1279]	190	debut = .FALSE.
	191	ENDIF ! of IF (firstcal)
[524]	192
	193	c
	194	c
	195	c-----------------------------------------------------------------------
	196	c 40. transformation des variables dynamiques en variables physiques:
	197	c ---------------------------------------------------------------
	198
	199	c 41. pressions au sol (en Pascals)
	200	c ----------------------------------
	201
[5090]	202
[524]	203	zpsrf(1) = pps(1,1)
	204
	205	ig0 = 2
	206	DO j = 2,jjm
	207	CALL SCOPY( iim,pps(1,j),1,zpsrf(ig0), 1 )
	208	ig0 = ig0+iim
	209	ENDDO
	210
	211	zpsrf(ngridmx) = pps(1,jjp1)
	212
	213
[2604]	214	c 42. pression intercouches et fonction d'Exner:
[524]	215	c
	216	c -----------------------------------------------------------------
	217	c .... zplev definis aux (llm +1) interfaces des couches ....
[5090]	218	c .... zplay definis aux ( llm ) milieux des couches ....
[524]	219	c -----------------------------------------------------------------
	220
	221	c ... Exner = cp * ( p(l) / preff ) ** kappa ....
	222	c
	223	unskap = 1./ kappa
	224	c
[2604]	225	DO l = 1, llm
	226	zpk( 1,l ) = ppk(1,1,l)
[524]	227	zplev( 1,l ) = pp(1,1,l)
	228	ig0 = 2
	229	DO j = 2, jjm
	230	DO i =1, iim
[2604]	231	zpk( ig0,l ) = ppk(i,j,l)
[524]	232	zplev( ig0,l ) = pp(i,j,l)
	233	ig0 = ig0 +1
	234	ENDDO
	235	ENDDO
[2604]	236	zpk( ngridmx,l ) = ppk(1,jjp1,l)
[524]	237	zplev( ngridmx,l ) = pp(1,jjp1,l)
	238	ENDDO
[2604]	239	zplev( 1,llmp1 ) = pp(1,1,llmp1)
	240	ig0 = 2
	241	DO j = 2, jjm
	242	DO i =1, iim
	243	zplev( ig0,llmp1 ) = pp(i,j,llmp1)
	244	ig0 = ig0 +1
	245	ENDDO
	246	ENDDO
	247	zplev( ngridmx,llmp1 ) = pp(1,jjp1,llmp1)
[524]	248	c
	249	c
	250
	251	c 43. temperature naturelle (en K) et pressions milieux couches .
	252	c ---------------------------------------------------------------
	253
	254	DO l=1,llm
	255
	256	pksurcp = ppk(1,1,l) / cpp
	257	zplay(1,l) = preff * pksurcp ** unskap
	258	ztfi(1,l) = pteta(1,1,l) * pksurcp
	259	pcvgt(1,l) = pdteta(1,1,l) * pksurcp / pmasse(1,1,l)
	260	ig0 = 2
	261
	262	DO j = 2, jjm
	263	DO i = 1, iim
	264	pksurcp = ppk(i,j,l) / cpp
	265	zplay(ig0,l) = preff * pksurcp ** unskap
	266	ztfi(ig0,l) = pteta(i,j,l) * pksurcp
	267	pcvgt(ig0,l) = pdteta(i,j,l) * pksurcp / pmasse(i,j,l)
	268	ig0 = ig0 + 1
	269	ENDDO
	270	ENDDO
	271
	272	pksurcp = ppk(1,jjp1,l) / cpp
	273	zplay(ig0,l) = preff * pksurcp ** unskap
	274	ztfi (ig0,l) = pteta(1,jjp1,l) * pksurcp
	275	pcvgt(ig0,l) = pdteta(1,jjp1,l) * pksurcp/ pmasse(1,jjp1,l)
	276
	277	ENDDO
	278
	279	c 43.bis traceurs
	280	c ---------------
	281	c
[4056]	282	itr=0
[1146]	283	DO iq=1,nqtot
[4056]	284	IF(.NOT.tracers(iq)%isAdvected) CYCLE
	285	itr = itr + 1
[524]	286	DO l=1,llm
[4056]	287	zqfi(1,l,itr) = pq(1,1,l,iq)
	288	ig0 = 2
[524]	289	DO j=2,jjm
	290	DO i = 1, iim
[4056]	291	zqfi(ig0,l,itr) = pq(i,j,l,iq)
[524]	292	ig0 = ig0 + 1
	293	ENDDO
	294	ENDDO
[4056]	295	zqfi(ig0,l,itr) = pq(1,jjp1,l,iq)
[524]	296	ENDDO
	297	ENDDO
	298
	299	c convergence dynamique pour les traceurs "EAU"
[1279]	300	! Earth-specific treatment of first 2 tracers (water)
	301	if (planet_type=="earth") then
	302	DO iq=1,2
[524]	303	DO l=1,llm
	304	pcvgq(1,l,iq)= pdq(1,1,l,iq) / pmasse(1,1,l)
	305	ig0 = 2
	306	DO j=2,jjm
	307	DO i = 1, iim
	308	pcvgq(ig0,l,iq) = pdq(i,j,l,iq) / pmasse(i,j,l)
	309	ig0 = ig0 + 1
	310	ENDDO
	311	ENDDO
	312	pcvgq(ig0,l,iq)= pdq(1,jjp1,l,iq) / pmasse(1,jjp1,l)
	313	ENDDO
[1279]	314	ENDDO
	315	endif ! of if (planet_type=="earth")
[524]	316
	317
	318	c Geopotentiel calcule par rapport a la surface locale:
	319	c -----------------------------------------------------
	320
	321	CALL gr_dyn_fi(llm,iip1,jjp1,ngridmx,pphi,zphi)
	322	CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,pphis,zphis)
	323	DO l=1,llm
[1403]	324	DO ig=1,ngridmx
	325	zphi(ig,l)=zphi(ig,l)-zphis(ig)
	326	ENDDO
[524]	327	ENDDO
	328
	329	c .... Calcul de la vitesse verticale ( en Pams ou Kg/s ) ....
[5090]	330	c JG : ancien calcule de omega utilise dans physiq.F. Maintenant le flux
	331	c de masse est calclue dans advtrac.F
[960]	332	c DO l=1,llm
	333	c pvervel(1,l)=pw(1,1,l) * g /apoln
	334	c ig0=2
	335	c DO j=2,jjm
	336	c DO i = 1, iim
	337	c pvervel(ig0,l) = pw(i,j,l) * g * unsaire(i,j)
	338	c ig0 = ig0 + 1
	339	c ENDDO
	340	c ENDDO
	341	c pvervel(ig0,l)=pw(1,jjp1,l) * g /apols
	342	c ENDDO
[524]	343
	344	c
	345	c 45. champ u:
	346	c ------------
	347
[5086]	348	DO l=1,llm
[524]	349
[5086]	350	DO j=2,jjm
[524]	351	ig0 = 1+(j-2)*iim
[5090]	352	zufi(ig0+1,l)= 0.5 *
[524]	353	$ ( pucov(iim,j,l)/cu(iim,j) + pucov(1,j,l)/cu(1,j) )
[5090]	354	pcvgu(ig0+1,l)= 0.5 *
[524]	355	$ ( pducov(iim,j,l)/cu(iim,j) + pducov(1,j,l)/cu(1,j) )
[5086]	356	DO i=2,iim
[524]	357	zufi(ig0+i,l)= 0.5 *
	358	$ ( pucov(i-1,j,l)/cu(i-1,j) + pucov(i,j,l)/cu(i,j) )
	359	pcvgu(ig0+i,l)= 0.5 *
	360	$ ( pducov(i-1,j,l)/cu(i-1,j) + pducov(i,j,l)/cu(i,j) )
[5086]	361	END DO
	362	END DO
[524]	363
[5086]	364	END DO
[524]	365
	366
[2333]	367	C Alvaro de la Camara (May 2014)
	368	C 46.1 Calcul de la vorticite et passage sur la grille physique
	369	C --------------------------------------------------------------
	370	DO l=1,llm
	371	do i=1,iim
	372	do j=1,jjm
	373	zrot(i,j,l) = (pvcov(i+1,j,l) - pvcov(i,j,l)
[5090]	374	$ + pucov(i,j+1,l) - pucov(i,j,l))
	375	$ / (cu(i,j)+cu(i,j+1))
[2333]	376	$ / (cv(i+1,j)+cv(i,j)) *4
	377	enddo
	378	enddo
	379	ENDDO
	380
[524]	381	c 46.champ v:
	382	c -----------
	383
	384	DO l=1,llm
	385	DO j=2,jjm
	386	ig0=1+(j-2)*iim
	387	DO i=1,iim
	388	zvfi(ig0+i,l)= 0.5 *
	389	$ ( pvcov(i,j-1,l)/cv(i,j-1) + pvcov(i,j,l)/cv(i,j) )
	390	pcvgv(ig0+i,l)= 0.5 *
	391	$ ( pdvcov(i,j-1,l)/cv(i,j-1) + pdvcov(i,j,l)/cv(i,j) )
	392	ENDDO
[2333]	393	zrfi(ig0 + 1,l)= 0.25 *(zrot(iim,j-1,l)+zrot(iim,j,l)
	394	& +zrot(1,j-1,l)+zrot(1,j,l))
	395	DO i=2,iim
	396	zrfi(ig0 + i,l)= 0.25 *(zrot(i-1,j-1,l)+zrot(i-1,j,l)
	397	$ +zrot(i,j-1,l)+zrot(i,j,l)) ! AdlC MAY 2014
	398	ENDDO
[524]	399	ENDDO
	400	ENDDO
	401
	402
[5090]	403	c 47. champs de vents aux pole nord
[524]	404	c ------------------------------
	405	c U = 1 / pi * integrale [ v * cos(long) * d long ]
	406	c V = 1 / pi * integrale [ v * sin(long) * d long ]
	407
	408	DO l=1,llm
	409
	410	z1(1) =(rlonu(1)-rlonu(iim)+2.pi)pvcov(1,1,l)/cv(1,1)
	411	z1bis(1)=(rlonu(1)-rlonu(iim)+2.pi)pdvcov(1,1,l)/cv(1,1)
	412	DO i=2,iim
	413	z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,1,l)/cv(i,1)
	414	z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,1,l)/cv(i,1)
	415	ENDDO
	416
	417	DO i=1,iim
	418	zcos(i) = COS(rlonv(i))*z1(i)
	419	zcosbis(i)= COS(rlonv(i))*z1bis(i)
	420	zsin(i) = SIN(rlonv(i))*z1(i)
	421	zsinbis(i)= SIN(rlonv(i))*z1bis(i)
	422	ENDDO
	423
	424	zufi(1,l) = SSUM(iim,zcos,1)/pi
	425	pcvgu(1,l) = SSUM(iim,zcosbis,1)/pi
	426	zvfi(1,l) = SSUM(iim,zsin,1)/pi
	427	pcvgv(1,l) = SSUM(iim,zsinbis,1)/pi
[2333]	428	zrfi(1, l) = 0.
[524]	429	ENDDO
	430
	431
	432	c 48. champs de vents aux pole sud:
	433	c ---------------------------------
	434	c U = 1 / pi * integrale [ v * cos(long) * d long ]
	435	c V = 1 / pi * integrale [ v * sin(long) * d long ]
	436
	437	DO l=1,llm
	438
	439	z1(1) =(rlonu(1)-rlonu(iim)+2.pi)pvcov(1,jjm,l)/cv(1,jjm)
	440	z1bis(1)=(rlonu(1)-rlonu(iim)+2.pi)pdvcov(1,jjm,l)/cv(1,jjm)
	441	DO i=2,iim
	442	z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,jjm,l)/cv(i,jjm)
	443	z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,jjm,l)/cv(i,jjm)
[1403]	444	ENDDO
[524]	445
	446	DO i=1,iim
	447	zcos(i) = COS(rlonv(i))*z1(i)
	448	zcosbis(i) = COS(rlonv(i))*z1bis(i)
	449	zsin(i) = SIN(rlonv(i))*z1(i)
	450	zsinbis(i) = SIN(rlonv(i))*z1bis(i)
[1403]	451	ENDDO
[524]	452
	453	zufi(ngridmx,l) = SSUM(iim,zcos,1)/pi
	454	pcvgu(ngridmx,l) = SSUM(iim,zcosbis,1)/pi
	455	zvfi(ngridmx,l) = SSUM(iim,zsin,1)/pi
	456	pcvgv(ngridmx,l) = SSUM(iim,zsinbis,1)/pi
[2333]	457	zrfi(ngridmx, l) = 0.
[524]	458	ENDDO
[644]	459	c
[960]	460	c On change de grille, dynamique vers physiq, pour le flux de masse verticale
[524]	461	CALL gr_dyn_fi(llm,iip1,jjp1,ngridmx,flxw,flxwfi)
	462
	463	c-----------------------------------------------------------------------
	464	c Appel de la physique:
	465	c ---------------------
	466
	467
[1403]	468
[1407]	469	! write(lunout,*) 'PHYSIQUE AVEC NSPLIT_PHYS=',nsplit_phys
[1403]	470	zdt_split=dtphys/nsplit_phys
	471	zdufic(:,:)=0.
	472	zdvfic(:,:)=0.
	473	zdtfic(:,:)=0.
	474	zdqfic(:,:,:)=0.
	475
[5090]	476	IF (CPPKEY_PHYS) THEN
[1403]	477
[1615]	478	do isplit=1,nsplit_phys
	479
[1403]	480	jH_cur_split=jH_cur+(isplit-1) * dtvr / (daysec *nsplit_phys)
	481	debut_split=debut.and.isplit==1
	482	lafin_split=lafin.and.isplit==nsplit_phys
	483
[4056]	484	! if (planet_type=="earth") then
[4046]	485	CALL call_physiq(ngridmx,llm,nqtot,tracers(:)%name,
[2418]	486	& debut_split,lafin_split,
	487	& jD_cur,jH_cur_split,zdt_split,
	488	& zplev,zplay,
[2604]	489	& zpk,zphi,zphis,
[2418]	490	& presnivs,
	491	& zufi,zvfi,zrfi,ztfi,zqfi,
	492	& flxwfi,pducov,
	493	& zdufi,zdvfi,zdtfi,zdqfi,zdpsrf)
	494	!
	495	! else if ( planet_type=="generic" ) then
	496	!
	497	! CALL physiq (ngridmx, !! ngrid
	498	! . llm, !! nlayer
	499	! . nqtot, !! nq
[4046]	500	! . tracers(:)%name,!! tracer names from dynamical core (given in infotrac)
[5090]	501	! . debut_split, !! firstcall
[2418]	502	! . lafin_split, !! lastcall
	503	! . jD_cur, !! pday. see leapfrog
	504	! . jH_cur_split, !! ptime "fraction of day"
	505	! . zdt_split, !! ptimestep
	506	! . zplev, !! pplev
	507	! . zplay, !! pplay
	508	! . zphi, !! pphi
	509	! . zufi, !! pu
	510	! . zvfi, !! pv
	511	! . ztfi, !! pt
	512	! . zqfi, !! pq
	513	! . flxwfi, !! pw !! or 0. anyway this is for diagnostic. not used in physiq.
	514	! . zdufi, !! pdu
	515	! . zdvfi, !! pdv
	516	! . zdtfi, !! pdt
	517	! . zdqfi, !! pdq
	518	! . zdpsrf, !! pdpsrf
	519	! . tracerdyn) !! tracerdyn <-- utilite ???
	520	!
	521	! endif ! of if (planet_type=="earth")
[1654]	522
[1403]	523	zufi(:,:)=zufi(:,:)+zdufi(:,:)*zdt_split
	524	zvfi(:,:)=zvfi(:,:)+zdvfi(:,:)*zdt_split
	525	ztfi(:,:)=ztfi(:,:)+zdtfi(:,:)*zdt_split
	526	zqfi(:,:,:)=zqfi(:,:,:)+zdqfi(:,:,:)*zdt_split
	527
	528	zdufic(:,:)=zdufic(:,:)+zdufi(:,:)
	529	zdvfic(:,:)=zdvfic(:,:)+zdvfi(:,:)
	530	zdtfic(:,:)=zdtfic(:,:)+zdtfi(:,:)
	531	zdqfic(:,:,:)=zdqfic(:,:,:)+zdqfi(:,:,:)
	532
[1615]	533	enddo ! of do isplit=1,nsplit_phys
	534
[5090]	535	END IF
[1615]	536
[1403]	537	zdufi(:,:)=zdufic(:,:)/nsplit_phys
	538	zdvfi(:,:)=zdvfic(:,:)/nsplit_phys
	539	zdtfi(:,:)=zdtfic(:,:)/nsplit_phys
	540	zdqfi(:,:,:)=zdqfic(:,:,:)/nsplit_phys
	541
[524]	542	c-----------------------------------------------------------------------
	543	c transformation des tendances physiques en tendances dynamiques:
	544	c ---------------------------------------------------------------
	545
	546	c tendance sur la pression :
	547	c -----------------------------------
	548
	549	CALL gr_fi_dyn(1,ngridmx,iip1,jjp1,zdpsrf,pdpsfi)
	550	c
	551	c 62. enthalpie potentielle
	552	c ---------------------
	553
	554	DO l=1,llm
	555
	556	DO i=1,iip1
	557	pdhfi(i,1,l) = cpp * zdtfi(1,l) / ppk(i, 1 ,l)
	558	pdhfi(i,jjp1,l) = cpp * zdtfi(ngridmx,l)/ ppk(i,jjp1,l)
	559	ENDDO
	560
	561	DO j=2,jjm
	562	ig0=1+(j-2)*iim
	563	DO i=1,iim
	564	pdhfi(i,j,l) = cpp * zdtfi(ig0+i,l) / ppk(i,j,l)
	565	ENDDO
	566	pdhfi(iip1,j,l) = pdhfi(1,j,l)
	567	ENDDO
	568
	569	ENDDO
	570
	571
	572	c 62. humidite specifique
	573	c ---------------------
[1279]	574	! Ehouarn: removed this useless bit: was overwritten at step 63 anyways
	575	! DO iq=1,nqtot
	576	! DO l=1,llm
	577	! DO i=1,iip1
	578	! pdqfi(i,1,l,iq) = zdqfi(1,l,iq)
	579	! pdqfi(i,jjp1,l,iq) = zdqfi(ngridmx,l,iq)
	580	! ENDDO
	581	! DO j=2,jjm
	582	! ig0=1+(j-2)*iim
	583	! DO i=1,iim
	584	! pdqfi(i,j,l,iq) = zdqfi(ig0+i,l,iq)
	585	! ENDDO
	586	! pdqfi(iip1,j,l,iq) = pdqfi(1,j,l,iq)
	587	! ENDDO
	588	! ENDDO
	589	! ENDDO
[524]	590
	591	c 63. traceurs
	592	c ------------
	593	C initialisation des tendances
[1279]	594	pdqfi(:,:,:,:)=0.
[524]	595	C
[4056]	596	itr = 0
[1146]	597	DO iq=1,nqtot
[4056]	598	IF(.NOT.tracers(iq)%isAdvected) CYCLE
	599	itr = itr + 1
[524]	600	DO l=1,llm
	601	DO i=1,iip1
[4056]	602	pdqfi(i,1,l,iq) = zdqfi(1,l,itr)
	603	pdqfi(i,jjp1,l,iq) = zdqfi(ngridmx,l,itr)
[524]	604	ENDDO
	605	DO j=2,jjm
	606	ig0=1+(j-2)*iim
	607	DO i=1,iim
[4056]	608	pdqfi(i,j,l,iq) = zdqfi(ig0+i,l,itr)
[524]	609	ENDDO
[4056]	610	pdqfi(iip1,j,l,iq) = pdqfi(1,j,l,itr)
[524]	611	ENDDO
	612	ENDDO
	613	ENDDO
	614
	615	c 65. champ u:
	616	c ------------
	617
	618	DO l=1,llm
	619
	620	DO i=1,iip1
	621	pdufi(i,1,l) = 0.
	622	pdufi(i,jjp1,l) = 0.
	623	ENDDO
	624
	625	DO j=2,jjm
	626	ig0=1+(j-2)*iim
	627	DO i=1,iim-1
	628	pdufi(i,j,l)=
	629	$ 0.5(zdufi(ig0+i,l)+zdufi(ig0+i+1,l))cu(i,j)
	630	ENDDO
	631	pdufi(iim,j,l)=
	632	$ 0.5(zdufi(ig0+1,l)+zdufi(ig0+iim,l))cu(iim,j)
	633	pdufi(iip1,j,l)=pdufi(1,j,l)
	634	ENDDO
	635
	636	ENDDO
	637
	638
	639	c 67. champ v:
	640	c ------------
	641
	642	DO l=1,llm
	643
	644	DO j=2,jjm-1
	645	ig0=1+(j-2)*iim
	646	DO i=1,iim
	647	pdvfi(i,j,l)=
	648	$ 0.5(zdvfi(ig0+i,l)+zdvfi(ig0+i+iim,l))cv(i,j)
	649	ENDDO
	650	pdvfi(iip1,j,l) = pdvfi(1,j,l)
	651	ENDDO
	652	ENDDO
	653
	654
	655	c 68. champ v pres des poles:
	656	c ---------------------------
	657	c v = U * cos(long) + V * SIN(long)
	658
	659	DO l=1,llm
	660
	661	DO i=1,iim
	662	pdvfi(i,1,l)=
	663	$ zdufi(1,l)COS(rlonv(i))+zdvfi(1,l)SIN(rlonv(i))
	664	pdvfi(i,jjm,l)=zdufi(ngridmx,l)*COS(rlonv(i))
	665	$ +zdvfi(ngridmx,l)*SIN(rlonv(i))
	666	pdvfi(i,1,l)=
	667	$ 0.5(pdvfi(i,1,l)+zdvfi(i+1,l))cv(i,1)
	668	pdvfi(i,jjm,l)=
	669	$ 0.5(pdvfi(i,jjm,l)+zdvfi(ngridmx-iip1+i,l))cv(i,jjm)
	670	ENDDO
	671
	672	pdvfi(iip1,1,l) = pdvfi(1,1,l)
	673	pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l)
	674
	675	ENDDO
	676
	677	c-----------------------------------------------------------------------
	678	firstcal = .FALSE.
	679
	680	RETURN
	681	END

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

Download in other formats:

Original Format