Context Navigation

calfis_p.F @ 1114

Last change on this file since 1114 was 1000, checked in by Laurent Fairhead, 16 years ago

Modifs sur le parallelisme: masquage dans la physique
Inclusion strato
mise en coherence etat0
le mode offline fonctionne maintenant en parallele,
les fichiers de la dynamiques sont correctement sortis et peuvent etre reconstruit avec rebuild
la version parallele de la dynamique peut s'executer sans MPI (sur 1 proc)
L'OPENMP fonctionne maintenant sans la parallelisation MPI.

YM
LF

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

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

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

Download in other formats:

Original Format