Context Navigation

calfis_p.F @ 1114

Last change on this file since 1114 was 1114, checked in by jghattas, 15 years ago

Creation du module infotrac:

contient les variables de advtrac.h
contient la subroutine iniadvtrac renommer en infotrac_init
le nombre des traceurs est lu dans tracer.def en dynamique (ou par default ou recu par INCA)
ce module est utilise dans la dynamique et la physique
contient aussi la variable nbtr qui avant etait stockee dans dimphy

Le fichier advtrac.h n'existe plus.
La compilation ne prend plus en compte le nombre de traceur.

/JG

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(lafin,
7	$ rdayvrai,
8	$ heure,
9	$ pucov,
10	$ pvcov,
11	$ pteta,
12	$ pq,
13	$ pmasse,
14	$ pps,
15	$ pp,
16	$ ppk,
17	$ pphis,
18	$ pphi,
19	$ pducov,
20	$ pdvcov,
21	$ pdteta,
22	$ pdq,
23	$ flxw,
24	$ clesphy0,
25	$ pdufi,
26	$ pdvfi,
27	$ pdhfi,
28	$ pdqfi,
29	$ pdpsfi)
30	c
31	c Auteur : P. Le Van, F. Hourdin
32	c .........
33	USE dimphy
34	USE mod_phys_lmdz_para, mpi_root_xx=>mpi_root
35	USE parallel, ONLY : omp_chunk, using_mpi
36	USE mod_interface_dyn_phys
37	USE Write_Field
38	Use Write_field_p
39	USE Times
40	USE IOPHY
41	USE infotrac
42
43	IMPLICIT NONE
44	c=======================================================================
45	c
46	c 1. rearrangement des tableaux et transformation
47	c variables dynamiques > variables physiques
48	c 2. calcul des termes physiques
49	c 3. retransformation des tendances physiques en tendances dynamiques
50	c
51	c remarques:
52	c ----------
53	c
54	c - les vents sont donnes dans la physique par leurs composantes
55	c naturelles.
56	c - la variable thermodynamique de la physique est une variable
57	c intensive : T
58	c pour la dynamique on prend T * ( preff / p(l) ) **kappa
59	c - les deux seules variables dependant de la geometrie necessaires
60	c pour la physique sont la latitude pour le rayonnement et
61	c l'aire de la maille quand on veut integrer une grandeur
62	c horizontalement.
63	c - les points de la physique sont les points scalaires de la
64	c la dynamique; numerotation:
65	c 1 pour le pole nord
66	c (jjm-1)*iim pour l'interieur du domaine
67	c ngridmx pour le pole sud
68	c ---> ngridmx=2+(jjm-1)*iim
69	c
70	c Input :
71	c -------
72	c ecritphy frequence d'ecriture (en jours)de histphy
73	c pucov covariant zonal velocity
74	c pvcov covariant meridional velocity
75	c pteta potential temperature
76	c pps surface pressure
77	c pmasse masse d'air dans chaque maille
78	c pts surface temperature (K)
79	c callrad clef d'appel au rayonnement
80	c
81	c Output :
82	c --------
83	c pdufi tendency for the natural zonal velocity (ms-1)
84	c pdvfi tendency for the natural meridional velocity
85	c pdhfi tendency for the potential temperature
86	c pdtsfi tendency for the surface temperature
87	c
88	c pdtrad radiative tendencies \ both input
89	c pfluxrad radiative fluxes / and output
90	c
91	c=======================================================================
92	c
93	c-----------------------------------------------------------------------
94	c
95	c 0. Declarations :
96	c ------------------
97
98	#include "dimensions.h"
99	#include "paramet.h"
100	#include "temps.h"
101
102	INTEGER ngridmx
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,nqtot)
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,nqtot)
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,nqtot)
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,nqtot))
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,nqtot))
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,nqtot
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,nqtot,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,nqtot))
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,nqtot))
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,nqtot
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,nqtot
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	. debut,
667	. lafin,
668	. rdayvrai,
669	. heure,
670	. dtphys,
671	. zplev_omp,
672	. zplay_omp,
673	. zphi_omp,
674	. zphis_omp,
675	. presnivs_omp,
676	. clesphy0,
677	. zufi_omp,
678	. zvfi_omp,
679	. ztfi_omp,
680	. zqfi_omp,
681	c . pvervel_omp,
682	c#ifdef INCA
683	. flxwfi_omp,
684	c#endif
685	. zdufi_omp,
686	. zdvfi_omp,
687	. zdtfi_omp,
688	. zdqfi_omp,
689	. zdpsrf_omp,
690	cIM diagnostique PVteta, Amip2
691	. pducov,
692	. PVteta)
693
694	cym call WriteField_phy_p('zdtfi_omp',zdtfi_omp(:,:),llm)
695
696	c$OMP BARRIER
697
698	do l=1,llm+1
699	do i=1,klon
700	zplev(offset+i,l)=zplev_omp(i,l)
701	enddo
702	enddo
703
704	do l=1,llm
705	do i=1,klon
706	zplay(offset+i,l)=zplay_omp(i,l)
707	enddo
708	enddo
709
710	do l=1,llm
711	do i=1,klon
712	zphi(offset+i,l)=zphi_omp(i,l)
713	enddo
714	enddo
715
716
717	do i=1,klon
718	zphis(offset+i)=zphis_omp(i)
719	enddo
720
721
722	do l=1,llm
723	presnivs(l)=presnivs_omp(l)
724	enddo
725
726	do l=1,llm
727	do i=1,klon
728	zufi(offset+i,l)=zufi_omp(i,l)
729	enddo
730	enddo
731
732	do l=1,llm
733	do i=1,klon
734	zvfi(offset+i,l)=zvfi_omp(i,l)
735	enddo
736	enddo
737
738	do l=1,llm
739	do i=1,klon
740	ztfi(offset+i,l)=ztfi_omp(i,l)
741	enddo
742	enddo
743
744	do iq=1,nqtot
745	do l=1,llm
746	do i=1,klon
747	zqfi(offset+i,l,iq)=zqfi_omp(i,l,iq)
748	enddo
749	enddo
750	enddo
751
752	c do l=1,llm
753	c do i=1,klon
754	c pvervel(offset+i,l)=pvervel_omp(i,l)
755	c enddo
756	c enddo
757
758	do l=1,llm
759	do i=1,klon
760	zdufi(offset+i,l)=zdufi_omp(i,l)
761	enddo
762	enddo
763
764	do l=1,llm
765	do i=1,klon
766	zdvfi(offset+i,l)=zdvfi_omp(i,l)
767	enddo
768	enddo
769
770	do l=1,llm
771	do i=1,klon
772	zdtfi(offset+i,l)=zdtfi_omp(i,l)
773	enddo
774	enddo
775
776	do iq=1,nqtot
777	do l=1,llm
778	do i=1,klon
779	zdqfi(offset+i,l,iq)=zdqfi_omp(i,l,iq)
780	enddo
781	enddo
782	enddo
783
784	do i=1,klon
785	zdpsrf(offset+i)=zdpsrf_omp(i)
786	enddo
787
788
789	cc$OMP END PARALLEL
790	klon=klon_mpi
791	500 CONTINUE
792	c$OMP BARRIER
793
794	c$OMP MASTER
795	cym call WriteField_phy('zdtfi',zdtfi(:,:),llm)
796	call stop_timer(timer_physic)
797	c$OMP END MASTER
798
799	IF (using_mpi) THEN
800
801	if (MPI_rank>0) then
802
803	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
804	DO l=1,llm
805	du_send(1:iim,l)=zdufi(1:iim,l)
806	dv_send(1:iim,l)=zdvfi(1:iim,l)
807	ENDDO
808	c$OMP END DO NOWAIT
809
810	c$OMP BARRIER
811	#ifdef CPP_MPI
812	c$OMP MASTER
813	!$OMP CRITICAL (MPI)
814	call MPI_ISSEND(du_send,iim*llm,MPI_REAL8,MPI_Rank-1,401,
815	& COMM_LMDZ,Req(1),ierr)
816	call MPI_ISSEND(dv_send,iim*llm,MPI_REAL8,MPI_Rank-1,402,
817	& COMM_LMDZ,Req(2),ierr)
818	!$OMP END CRITICAL (MPI)
819	c$OMP END MASTER
820	#endif
821	c$OMP BARRIER
822
823	endif
824
825	if (MPI_rank<MPI_Size-1) then
826	c$OMP BARRIER
827	#ifdef CPP_MPI
828	c$OMP MASTER
829	!$OMP CRITICAL (MPI)
830	call MPI_IRECV(du_recv,iim*llm,MPI_REAL8,MPI_Rank+1,401,
831	& COMM_LMDZ,Req(3),ierr)
832	call MPI_IRECV(dv_recv,iim*llm,MPI_REAL8,MPI_Rank+1,402,
833	& COMM_LMDZ,Req(4),ierr)
834	!$OMP END CRITICAL (MPI)
835	c$OMP END MASTER
836	#endif
837	endif
838
839	c$OMP BARRIER
840
841
842	#ifdef CPP_MPI
843	c$OMP MASTER
844	!$OMP CRITICAL (MPI)
845	if (MPI_rank>0 .and. MPI_rank< MPI_Size-1) then
846	call MPI_WAITALL(4,Req(1),Status,ierr)
847	else if (MPI_rank>0) then
848	call MPI_WAITALL(2,Req(1),Status,ierr)
849	else if (MPI_rank <MPI_Size-1) then
850	call MPI_WAITALL(2,Req(3),Status,ierr)
851	endif
852	!$OMP END CRITICAL (MPI)
853	c$OMP END MASTER
854	#endif
855
856	c$OMP BARRIER
857
858	ENDIF ! using_mpi
859
860
861	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
862	DO l=1,llm
863
864	zdufi2(1:klon,l)=zdufi(1:klon,l)
865	zdufi2(klon+1:klon+iim,l)=du_recv(1:iim,l)
866
867	zdvfi2(1:klon,l)=zdvfi(1:klon,l)
868	zdvfi2(klon+1:klon+iim,l)=dv_recv(1:iim,l)
869
870	pdhfi(:,jj_begin,l)=0
871	pdqfi(:,jj_begin,l,:)=0
872	pdufi(:,jj_begin,l)=0
873	pdvfi(:,jj_begin,l)=0
874
875	if (.not. is_south_pole) then
876	pdhfi(:,jj_end,l)=0
877	pdqfi(:,jj_end,l,:)=0
878	pdufi(:,jj_end,l)=0
879	pdvfi(:,jj_end,l)=0
880	endif
881
882	ENDDO
883	c$OMP END DO NOWAIT
884
885	c$OMP MASTER
886	pdpsfi(:,jj_begin)=0
887	if (.not. is_south_pole) then
888	pdpsfi(:,jj_end)=0
889	endif
890	c$OMP END MASTER
891	c-----------------------------------------------------------------------
892	c transformation des tendances physiques en tendances dynamiques:
893	c ---------------------------------------------------------------
894
895	c tendance sur la pression :
896	c -----------------------------------
897	c$OMP MASTER
898	CALL gr_fi_dyn_p(1,klon,iip1,jjp1,zdpsrf,pdpsfi)
899	c$OMP END MASTER
900	c
901	c 62. enthalpie potentielle
902	c ---------------------
903
904	kstart=1
905	kend=klon
906
907	if (is_north_pole) kstart=2
908	if (is_south_pole) kend=klon-1
909
910	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
911	DO l=1,llm
912
913	!!cdir NODEP
914	do ig0=kstart,kend
915	i=index_i(ig0)
916	j=index_j(ig0)
917	pdhfi(i,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l)
918	if (i==1) pdhfi(iip1,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l)
919	enddo
920
921	if (is_north_pole) then
922	DO i=1,iip1
923	pdhfi(i,1,l) = cpp * zdtfi(1,l) / ppk(i, 1 ,l)
924	enddo
925	endif
926
927	if (is_south_pole) then
928	DO i=1,iip1
929	pdhfi(i,jjp1,l) = cpp * zdtfi(klon,l)/ ppk(i,jjp1,l)
930	ENDDO
931	endif
932	ENDDO
933	c$OMP END DO NOWAIT
934
935	c 62. humidite specifique
936	c ---------------------
937
938	DO iq=1,nqtot
939	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
940	DO l=1,llm
941	!!cdir NODEP
942	do ig0=kstart,kend
943	i=index_i(ig0)
944	j=index_j(ig0)
945	pdqfi(i,j,l,iq) = zdqfi(ig0,l,iq)
946	if (i==1) pdqfi(iip1,j,l,iq) = zdqfi(ig0,l,iq)
947	enddo
948
949	if (is_north_pole) then
950	do i=1,iip1
951	pdqfi(i,1,l,iq) = zdqfi(1,l,iq)
952	enddo
953	endif
954
955	if (is_south_pole) then
956	do i=1,iip1
957	pdqfi(i,jjp1,l,iq) = zdqfi(klon,l,iq)
958	enddo
959	endif
960
961	ENDDO
962	c$OMP END DO NOWAIT
963	ENDDO
964
965	c 63. traceurs
966	c ------------
967	C initialisation des tendances
968
969	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
970	DO l=1,llm
971	pdqfi(:,:,l,:)=0.
972	ENDDO
973	c$OMP END DO NOWAIT
974
975	C
976
977	DO iq=1,nqtot
978	iiq=niadv(iq)
979	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
980	DO l=1,llm
981
982	!!cdir NODEP
983	DO ig0=kstart,kend
984	i=index_i(ig0)
985	j=index_j(ig0)
986	pdqfi(i,j,l,iiq) = zdqfi(ig0,l,iq)
987	if (i==1) pdqfi(iip1,j,l,iiq) = zdqfi(ig0,l,iq)
988	ENDDO
989
990	IF (is_north_pole) then
991	DO i=1,iip1
992	pdqfi(i,1,l,iiq) = zdqfi(1,l,iq)
993	ENDDO
994	ENDIF
995
996	IF (is_south_pole) then
997	DO i=1,iip1
998	pdqfi(i,jjp1,l,iiq) = zdqfi(klon,l,iq)
999	ENDDO
1000	ENDIF
1001
1002	ENDDO
1003	c$OMP END DO NOWAIT
1004	ENDDO
1005
1006	c 65. champ u:
1007	c ------------
1008	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1009	DO l=1,llm
1010	!!cdir NODEP
1011	do ig0=kstart,kend
1012	i=index_i(ig0)
1013	j=index_j(ig0)
1014
1015	if (i/=iim) then
1016	pdufi(i,j,l)=0.5(zdufi2(ig0,l)+zdufi2(ig0+1,l))cu(i,j)
1017	endif
1018
1019	if (i==1) then
1020	pdufi(iim,j,l)=0.5*( zdufi2(ig0,l)
1021	$ + zdufi2(ig0+iim-1,l))*cu(iim,j)
1022	pdufi(iip1,j,l)=0.5(zdufi2(ig0,l)+zdufi2(ig0+1,l))cu(i,j)
1023	endif
1024
1025	enddo
1026
1027	if (is_north_pole) then
1028	DO i=1,iip1
1029	pdufi(i,1,l) = 0.
1030	ENDDO
1031	endif
1032
1033	if (is_south_pole) then
1034	DO i=1,iip1
1035	pdufi(i,jjp1,l) = 0.
1036	ENDDO
1037	endif
1038
1039	ENDDO
1040	c$OMP END DO NOWAIT
1041
1042	c 67. champ v:
1043	c ------------
1044
1045	kstart=1
1046	kend=klon
1047
1048	if (is_north_pole) kstart=2
1049	if (is_south_pole) kend=klon-1-iim
1050
1051	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1052	DO l=1,llm
1053	!!cdir NODEP
1054	do ig0=kstart,kend
1055	i=index_i(ig0)
1056	j=index_j(ig0)
1057	pdvfi(i,j,l)=0.5(zdvfi2(ig0,l)+zdvfi2(ig0+iim,l))cv(i,j)
1058	if (i==1) pdvfi(iip1,j,l) = 0.5*(zdvfi2(ig0,l)+
1059	$ zdvfi2(ig0+iim,l))
1060	$ *cv(i,j)
1061	enddo
1062
1063	ENDDO
1064	c$OMP END DO NOWAIT
1065
1066
1067	c 68. champ v pres des poles:
1068	c ---------------------------
1069	c v = U * cos(long) + V * SIN(long)
1070
1071	if (is_north_pole) then
1072
1073	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1074	DO l=1,llm
1075
1076	DO i=1,iim
1077	pdvfi(i,1,l)=
1078	$ zdufi(1,l)COS(rlonv(i))+zdvfi(1,l)SIN(rlonv(i))
1079
1080	pdvfi(i,1,l)=
1081	$ 0.5(pdvfi(i,1,l)+zdvfi(i+1,l))cv(i,1)
1082	ENDDO
1083
1084	pdvfi(iip1,1,l) = pdvfi(1,1,l)
1085
1086	ENDDO
1087	c$OMP END DO NOWAIT
1088
1089	endif
1090
1091	if (is_south_pole) then
1092
1093	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1094	DO l=1,llm
1095
1096	DO i=1,iim
1097	pdvfi(i,jjm,l)=zdufi(klon,l)*COS(rlonv(i))
1098	$ +zdvfi(klon,l)*SIN(rlonv(i))
1099
1100	pdvfi(i,jjm,l)=
1101	$ 0.5(pdvfi(i,jjm,l)+zdvfi(klon-iip1+i,l))cv(i,jjm)
1102	ENDDO
1103
1104	pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l)
1105
1106	ENDDO
1107	c$OMP END DO NOWAIT
1108
1109	endif
1110	c-----------------------------------------------------------------------
1111
1112	700 CONTINUE
1113
1114	firstcal = .FALSE.
1115
1116	RETURN
1117	END

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

Download in other formats:

Original Format