Context Navigation

calfis.f90 @ 5875

Last change on this file since 5875 was 5481, checked in by dcugnet, 10 months ago
Remove tracers attributes "isAdvected" and "isInPhysics" from infotrac (iadv is enough). Remove tracers attribute "isAdvected" from infotrac_phy (isInPhysics is now equivalent to former isInPhysics .AND. iadv > 0
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: 19.7 KB

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

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

Download in other formats:

Original Format