Context Navigation

calfis_p.F @ 960

Last change on this file since 960 was 960, checked in by lsce, 16 years ago

Ajoute du parametre config_inca dans conf_gcm.F config_inca='none'(sans INCA, par defaut) config_inca='chem'(avec INCA config chemie) config_inca='aero'(avec INCA config aerosol)
Menage parmis les cles CPP INCA
Enleve le calcul d'omega dans calfis.F et active le calcul correspondant dans advtrac.F(avant uniquement pour INCA).

Property svn:eol-style set to native
Property svn:keywords set to Author Date Id Revision

File size: 27.5 KB

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

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

Download in other formats:

Original Format