Context Navigation

← Previous Revision
Next Revision →
Blame
Revision Log

calfis_p.F

Last change on this file was 2604, checked in by Ehouarn Millour, 8 years ago
Add Exner function to the call_physiq arguments (not used by the Earth physics) to harmonize physics/dynamics interface. EM
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: 31.2 KB

Line
1	!
2	! $Id: calfis_p.F 2604 2016-07-26 15:37:18Z aclsce $
3	!
4	C
5	C
6	SUBROUTINE calfis_p(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	#ifdef CPP_PHYS
29	! If using physics
30	USE dimphy
31	USE mod_phys_lmdz_mpi_data, mpi_root_xx=>mpi_master
32	USE mod_phys_lmdz_omp_data, ONLY: klon_omp, klon_omp_begin
33	USE mod_const_mpi, ONLY: COMM_LMDZ
34	USE mod_interface_dyn_phys
35	USE IOPHY
36	#endif
37	#ifdef CPP_PARA
38	USE parallel_lmdz,ONLY:omp_chunk,using_mpi,jjb_u,jje_u,jjb_v,jje_v
39	$ ,jj_begin_dyn=>jj_begin,jj_end_dyn=>jj_end
40	USE Write_Field
41	Use Write_field_p
42	USE Times
43	#endif
44	USE infotrac, ONLY: nqtot, niadv, tname
45	USE control_mod, ONLY: planet_type, nsplit_phys
46	#ifdef CPP_PHYS
47	USE callphysiq_mod, ONLY: call_physiq
48	#endif
49	USE comvert_mod, ONLY: preff, presnivs
50	USE comconst_mod, ONLY: cpp, daysec, dtphys, dtvr, kappa, pi
51
52	IMPLICIT NONE
53	c=======================================================================
54	c
55	c 1. rearrangement des tableaux et transformation
56	c variables dynamiques > variables physiques
57	c 2. calcul des termes physiques
58	c 3. retransformation des tendances physiques en tendances dynamiques
59	c
60	c remarques:
61	c ----------
62	c
63	c - les vents sont donnes dans la physique par leurs composantes
64	c naturelles.
65	c - la variable thermodynamique de la physique est une variable
66	c intensive : T
67	c pour la dynamique on prend T * ( preff / p(l) ) **kappa
68	c - les deux seules variables dependant de la geometrie necessaires
69	c pour la physique sont la latitude pour le rayonnement et
70	c l'aire de la maille quand on veut integrer une grandeur
71	c horizontalement.
72	c - les points de la physique sont les points scalaires de la
73	c la dynamique; numerotation:
74	c 1 pour le pole nord
75	c (jjm-1)*iim pour l'interieur du domaine
76	c ngridmx pour le pole sud
77	c ---> ngridmx=2+(jjm-1)*iim
78	c
79	c Input :
80	c -------
81	c ecritphy frequence d'ecriture (en jours)de histphy
82	c pucov covariant zonal velocity
83	c pvcov covariant meridional velocity
84	c pteta potential temperature
85	c pps surface pressure
86	c pmasse masse d'air dans chaque maille
87	c pts surface temperature (K)
88	c callrad clef d'appel au rayonnement
89	c
90	c Output :
91	c --------
92	c pdufi tendency for the natural zonal velocity (ms-1)
93	c pdvfi tendency for the natural meridional velocity
94	c pdhfi tendency for the potential temperature
95	c pdtsfi tendency for the surface temperature
96	c
97	c pdtrad radiative tendencies \ both input
98	c pfluxrad radiative fluxes / and output
99	c
100	c=======================================================================
101	c
102	c-----------------------------------------------------------------------
103	c
104	c 0. Declarations :
105	c ------------------
106
107	#include "dimensions.h"
108	#include "paramet.h"
109
110	INTEGER ngridmx
111	PARAMETER( ngridmx = 2+(jjm-1)*iim - 1/jjm )
112
113	#include "comgeom2.h"
114	#include "iniprint.h"
115	#ifdef CPP_MPI
116	include 'mpif.h'
117	#endif
118	c Arguments :
119	c -----------
120	LOGICAL,INTENT(IN) :: lafin ! .true. for the very last call to physics
121	REAL,INTENT(IN) :: jD_cur, jH_cur
122	REAL,INTENT(IN) :: pvcov(iip1,jjm,llm) ! covariant meridional velocity
123	REAL,INTENT(IN) :: pucov(iip1,jjp1,llm) ! covariant zonal velocity
124	REAL,INTENT(IN) :: pteta(iip1,jjp1,llm) ! potential temperature
125	REAL,INTENT(IN) :: pmasse(iip1,jjp1,llm) ! mass in each cell ! not used
126	REAL,INTENT(IN) :: pq(iip1,jjp1,llm,nqtot) ! tracers
127	REAL,INTENT(IN) :: pphis(iip1,jjp1) ! surface geopotential
128	REAL,INTENT(IN) :: pphi(iip1,jjp1,llm) ! geopotential
129
130	REAL,INTENT(IN) :: pdvcov(iip1,jjm,llm) ! dynamical tendency on vcov ! not used
131	REAL,INTENT(IN) :: pducov(iip1,jjp1,llm) ! dynamical tendency on ucov
132	REAL,INTENT(IN) :: pdteta(iip1,jjp1,llm) ! dynamical tendency on teta
133	! NB: pdteta is used only to compute pcvgt which is in fact not used...
134	REAL,INTENT(IN) :: pdq(iip1,jjp1,llm,nqtot) ! dynamical tendency on tracers
135	! NB: pdq is only used to compute pcvgq which is in fact not used...
136
137	REAL,INTENT(IN) :: pps(iip1,jjp1) ! surface pressure (Pa)
138	REAL,INTENT(IN) :: pp(iip1,jjp1,llmp1) ! pressure at mesh interfaces (Pa)
139	REAL,INTENT(IN) :: ppk(iip1,jjp1,llm) ! Exner at mid-layer
140	REAL,INTENT(IN) :: flxw(iip1,jjp1,llm) ! Vertical mass flux on lower mesh interfaces (kg/s) (on llm because flxw(:,:,llm+1)=0)
141
142	! tendencies (in */s) from the physics
143	REAL,INTENT(OUT) :: pdvfi(iip1,jjm,llm) ! tendency on covariant meridional wind
144	REAL,INTENT(OUT) :: pdufi(iip1,jjp1,llm) ! tendency on covariant zonal wind
145	REAL,INTENT(OUT) :: pdhfi(iip1,jjp1,llm) ! tendency on potential temperature (K/s)
146	REAL,INTENT(OUT) :: pdqfi(iip1,jjp1,llm,nqtot) ! tendency on tracers
147	REAL,INTENT(OUT) :: pdpsfi(iip1,jjp1) ! tendency on surface pressure (Pa/s)
148
149	#ifdef CPP_PARA
150	#ifdef CPP_PHYS
151	! Ehouarn: for now calfis_p needs some informations from physics to compile
152	c Local variables :
153	c -----------------
154
155	INTEGER i,j,l,ig0,ig,iq,iiq
156	REAL,ALLOCATABLE,SAVE :: zpsrf(:)
157	REAL,ALLOCATABLE,SAVE :: zplev(:,:),zplay(:,:)
158	REAL,ALLOCATABLE,SAVE :: zphi(:,:),zphis(:)
159	c
160	! REAL zrot(iip1,jjb_v:jje_v,llm) ! AdlC May 2014
161	REAL :: zrot(iip1,jjm,llm)
162	REAL,ALLOCATABLE,SAVE :: zufi(:,:), zvfi(:,:), zrfi(:,:)
163	REAL,ALLOCATABLE,SAVE :: ztfi(:,:),zqfi(:,:,:)
164	REAL,ALLOCATABLE,SAVE :: zpk(:,:)
165	c
166	REAL,ALLOCATABLE,SAVE :: pcvgu(:,:), pcvgv(:,:)
167	REAL,ALLOCATABLE,SAVE :: pcvgt(:,:), pcvgq(:,:,:)
168	c
169	REAL,ALLOCATABLE,SAVE :: zdufi(:,:),zdvfi(:,:)
170	REAL,ALLOCATABLE,SAVE :: zdtfi(:,:),zdqfi(:,:,:)
171	REAL,ALLOCATABLE,SAVE :: zdpsrf(:)
172	REAL,SAVE,ALLOCATABLE :: flxwfi(:,:) ! Flux de masse verticale sur la grille physiq
173
174	c
175	REAL,ALLOCATABLE,SAVE :: zplev_omp(:,:)
176	REAL,ALLOCATABLE,SAVE :: zplay_omp(:,:)
177	REAL,ALLOCATABLE,SAVE :: zpk_omp(:,:)
178	REAL,ALLOCATABLE,SAVE :: zphi_omp(:,:)
179	REAL,ALLOCATABLE,SAVE :: zphis_omp(:)
180	REAL,ALLOCATABLE,SAVE :: presnivs_omp(:)
181	REAL,ALLOCATABLE,SAVE :: zufi_omp(:,:)
182	REAL,ALLOCATABLE,SAVE :: zvfi_omp(:,:)
183	REAL,ALLOCATABLE,SAVE :: zrfi_omp(:,:)
184	REAL,ALLOCATABLE,SAVE :: ztfi_omp(:,:)
185	REAL,ALLOCATABLE,SAVE :: zqfi_omp(:,:,:)
186	REAL,ALLOCATABLE,SAVE :: zdufi_omp(:,:)
187	REAL,ALLOCATABLE,SAVE :: zdvfi_omp(:,:)
188	REAL,ALLOCATABLE,SAVE :: zdtfi_omp(:,:)
189	REAL,ALLOCATABLE,SAVE :: zdqfi_omp(:,:,:)
190	REAL,ALLOCATABLE,SAVE :: zdpsrf_omp(:)
191	REAL,SAVE,ALLOCATABLE :: flxwfi_omp(:,:) ! Flux de masse verticale sur la grille physiq
192
193	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
194	! Introduction du splitting (FH)
195	! Question pour Yann :
196	! J'ai été surpris au début que les tableaux zufi_omp, zdufi_omp n'co soitent
197	! en SAVE. Je crois comprendre que c'est parce que tu voulais qu'il
198	! soit allocatable (plutot par exemple que de passer une dimension
199	! dépendant du process en argument des routines) et que, du coup,
200	! le SAVE évite d'avoir à refaire l'allocation à chaque appel.
201	! Tu confirmes ?
202	! J'ai suivi le même principe pour les zdufic_omp
203	! Mais c'est surement bien que tu controles.
204	!
205
206	REAL,ALLOCATABLE,SAVE :: zdufic_omp(:,:)
207	REAL,ALLOCATABLE,SAVE :: zdvfic_omp(:,:)
208	REAL,ALLOCATABLE,SAVE :: zdtfic_omp(:,:)
209	REAL,ALLOCATABLE,SAVE :: zdqfic_omp(:,:,:)
210	REAL jH_cur_split,zdt_split
211	LOGICAL debut_split,lafin_split
212	INTEGER isplit
213	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
214
215	c$OMP THREADPRIVATE(zplev_omp,zplay_omp,zpk_omp,zphi_omp,zphis_omp,
216	c$OMP+ presnivs_omp,zufi_omp,zvfi_omp,ztfi_omp,
217	c$OMP+ zrfi_omp,zqfi_omp,zdufi_omp,zdvfi_omp,
218	c$OMP+ zdtfi_omp,zdqfi_omp,zdpsrf_omp,flxwfi_omp,
219	c$OMP+ zdufic_omp,zdvfic_omp,zdtfic_omp,zdqfic_omp)
220
221	LOGICAL,SAVE :: first_omp=.true.
222	c$OMP THREADPRIVATE(first_omp)
223
224	REAL zsin(iim),zcos(iim),z1(iim)
225	REAL zsinbis(iim),zcosbis(iim),z1bis(iim)
226	REAL unskap, pksurcp
227	c
228	REAL SSUM
229
230	LOGICAL,SAVE :: firstcal=.true., debut=.true.
231	c$OMP THREADPRIVATE(firstcal,debut)
232
233	REAL,SAVE,dimension(1:iim,1:llm):: du_send,du_recv,dv_send,dv_recv
234	INTEGER :: ierr
235	#ifdef CPP_MPI
236	INTEGER,dimension(MPI_STATUS_SIZE,4) :: Status
237	#else
238	INTEGER,dimension(1,4) :: Status
239	#endif
240	INTEGER, dimension(4) :: Req
241	REAL,ALLOCATABLE,SAVE:: zdufi2(:,:),zdvfi2(:,:)
242	integer :: k,kstart,kend
243	INTEGER :: offset
244	INTEGER :: jjb,jje
245
246	c
247	c-----------------------------------------------------------------------
248	c
249	c 1. Initialisations :
250	c --------------------
251	c
252
253	klon=klon_mpi
254
255	c
256	IF ( firstcal ) THEN
257	debut = .TRUE.
258	IF (ngridmx.NE.2+(jjm-1)*iim) THEN
259	write(lunout,*) 'STOP dans calfis'
260	write(lunout,*)
261	& 'La dimension ngridmx doit etre egale a 2 + (jjm-1)*iim'
262	write(lunout,*) ' ngridmx jjm iim '
263	write(lunout,*) ngridmx,jjm,iim
264	STOP
265	ENDIF
266	c$OMP MASTER
267	ALLOCATE(zpsrf(klon))
268	ALLOCATE(zplev(klon,llm+1),zplay(klon,llm))
269	ALLOCATE(zphi(klon,llm),zphis(klon))
270	ALLOCATE(zufi(klon,llm), zvfi(klon,llm),zrfi(klon,llm))
271	ALLOCATE(ztfi(klon,llm),zqfi(klon,llm,nqtot))
272	ALLOCATE(pcvgu(klon,llm), pcvgv(klon,llm))
273	ALLOCATE(pcvgt(klon,llm), pcvgq(klon,llm,2))
274	ALLOCATE(zdufi(klon,llm),zdvfi(klon,llm))
275	ALLOCATE(zdtfi(klon,llm),zdqfi(klon,llm,nqtot))
276	ALLOCATE(zdpsrf(klon))
277	ALLOCATE(zdufi2(klon+iim,llm),zdvfi2(klon+iim,llm))
278	ALLOCATE(flxwfi(klon,llm))
279	ALLOCATE(zpk(klon,llm))
280	c$OMP END MASTER
281	c$OMP BARRIER
282	ELSE
283	debut = .FALSE.
284	ENDIF
285
286	c
287	c
288	c-----------------------------------------------------------------------
289	c 40. transformation des variables dynamiques en variables physiques:
290	c ---------------------------------------------------------------
291
292	c 41. pressions au sol (en Pascals)
293	c ----------------------------------
294
295	c$OMP MASTER
296	call start_timer(timer_physic)
297	c$OMP END MASTER
298
299	c$OMP MASTER
300	!CDIR ON_ADB(index_i)
301	!CDIR ON_ADB(index_j)
302	do ig0=1,klon
303	i=index_i(ig0)
304	j=index_j(ig0)
305	zpsrf(ig0)=pps(i,j)
306	enddo
307	c$OMP END MASTER
308
309
310	c 42. pression intercouches :
311	c
312	c -----------------------------------------------------------------
313	c .... zplev definis aux (llm +1) interfaces des couches ....
314	c .... zplay definis aux ( llm ) milieux des couches ....
315	c -----------------------------------------------------------------
316
317	c ... Exner = cp * ( p(l) / preff ) ** kappa ....
318	c
319	unskap = 1./ kappa
320	c
321	c print *,omp_rank,'klon--->',klon
322	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
323	DO l = 1, llmp1
324	!CDIR ON_ADB(index_i)
325	!CDIR ON_ADB(index_j)
326	do ig0=1,klon
327	i=index_i(ig0)
328	j=index_j(ig0)
329	zplev( ig0,l ) = pp(i,j,l)
330	enddo
331	ENDDO
332	c$OMP END DO NOWAIT
333
334	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
335	DO l=1,llm
336	do ig0=1,klon
337	i=index_i(ig0)
338	j=index_j(ig0)
339	zpk(ig0,l)=ppk(i,j,l)
340	enddo
341	ENDDO
342	c$OMP END DO NOWAIT
343
344	c
345	c
346
347	c 43. temperature naturelle (en K) et pressions milieux couches .
348	c ---------------------------------------------------------------
349	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
350	DO l=1,llm
351	!CDIR ON_ADB(index_i)
352	!CDIR ON_ADB(index_j)
353	do ig0=1,klon
354	i=index_i(ig0)
355	j=index_j(ig0)
356	pksurcp = ppk(i,j,l) / cpp
357	zplay(ig0,l) = preff * pksurcp ** unskap
358	ztfi(ig0,l) = pteta(i,j,l) * pksurcp
359	enddo
360
361	ENDDO
362	c$OMP END DO NOWAIT
363
364	c 43.bis traceurs
365	c ---------------
366	c
367
368	DO iq=1,nqtot
369	iiq=niadv(iq)
370	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
371	DO l=1,llm
372	!CDIR ON_ADB(index_i)
373	!CDIR ON_ADB(index_j)
374	do ig0=1,klon
375	i=index_i(ig0)
376	j=index_j(ig0)
377	zqfi(ig0,l,iq) = pq(i,j,l,iiq)
378	enddo
379	ENDDO
380	c$OMP END DO NOWAIT
381	ENDDO
382
383
384	c Geopotentiel calcule par rapport a la surface locale:
385	c -----------------------------------------------------
386
387	CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,pphi,zphi)
388
389	CALL gr_dyn_fi_p(1,iip1,jjp1,klon,pphis,zphis)
390
391	c$OMP BARRIER
392
393	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
394	DO l=1,llm
395	DO ig=1,klon
396	zphi(ig,l)=zphi(ig,l)-zphis(ig)
397	ENDDO
398	ENDDO
399	c$OMP END DO NOWAIT
400
401
402	c
403	c 45. champ u:
404	c ------------
405
406	kstart=1
407	kend=klon
408
409	if (is_north_pole_dyn) kstart=2
410	if (is_south_pole_dyn) kend=klon-1
411
412	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
413	DO l=1,llm
414	!CDIR ON_ADB(index_i)
415	!CDIR ON_ADB(index_j)
416	!CDIR SPARSE
417	do ig0=kstart,kend
418	i=index_i(ig0)
419	j=index_j(ig0)
420	if (i==1) then
421	zufi(ig0,l)= 0.5 *( pucov(iim,j,l)/cu(iim,j)
422	$ + pucov(1,j,l)/cu(1,j) )
423	else
424	zufi(ig0,l)= 0.5*( pucov(i-1,j,l)/cu(i-1,j)
425	$ + pucov(i,j,l)/cu(i,j) )
426	endif
427	enddo
428	ENDDO
429	c$OMP END DO NOWAIT
430
431	c
432	C Alvaro de la Camara (May 2014)
433	C 46.1 Calcul de la vorticite et passage sur la grille physique
434	C --------------------------------------------------------------
435
436	jjb=jj_begin_dyn-1
437	jje=jj_end_dyn+1
438	if (is_north_pole_dyn) jjb=1
439	if (is_south_pole_dyn) jje=jjm
440
441	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
442
443	DO l=1,llm
444	do i=1,iim
445	do j=jjb,jje
446	zrot(i,j,l) = (pvcov(i+1,j,l) - pvcov(i,j,l)
447	$ + pucov(i,j+1,l) - pucov(i,j,l))
448	$ / (cu(i,j)+cu(i,j+1))
449	$ / (cv(i+1,j)+cv(i,j)) *4
450	enddo
451	enddo
452	ENDDO
453
454
455	c 46.2champ v:
456	c -----------
457
458	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
459	DO l=1,llm
460	!CDIR ON_ADB(index_i)
461	!CDIR ON_ADB(index_j)
462	DO ig0=kstart,kend
463	i=index_i(ig0)
464	j=index_j(ig0)
465	zvfi(ig0,l)= 0.5 *( pvcov(i,j-1,l)/cv(i,j-1)
466	$ + pvcov(i,j,l)/cv(i,j) )
467
468	if (j==1 .OR. j==jjp1) then ! AdlC MAY 2014
469	zrfi(ig0,l) = 0 ! AdlC MAY 2014
470	else
471	if(i==1)then
472	zrfi(ig0,l)= 0.25 *(zrot(iim,j-1,l)+zrot(iim,j,l)
473	$ +zrot(1,j-1,l)+zrot(1,j,l)) ! AdlC MAY 2014
474	else
475	zrfi(ig0,l)= 0.25 *(zrot(i-1,j-1,l)+zrot(i-1,j,l)
476	$ +zrot(i,j-1,l)+zrot(i,j,l)) ! AdlC MAY 2014
477	endif
478	endif
479
480
481	ENDDO
482	ENDDO
483	c$OMP END DO NOWAIT
484
485	c 47. champs de vents aux pole nord
486	c ------------------------------
487	c U = 1 / pi * integrale [ v * cos(long) * d long ]
488	c V = 1 / pi * integrale [ v * sin(long) * d long ]
489
490	if (is_north_pole_dyn) then
491	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
492	DO l=1,llm
493
494	z1(1) =(rlonu(1)-rlonu(iim)+2.pi)pvcov(1,1,l)/cv(1,1)
495	DO i=2,iim
496	z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,1,l)/cv(i,1)
497	ENDDO
498
499	DO i=1,iim
500	zcos(i) = COS(rlonv(i))*z1(i)
501	zsin(i) = SIN(rlonv(i))*z1(i)
502	ENDDO
503
504	zufi(1,l) = SSUM(iim,zcos,1)/pi
505	zvfi(1,l) = SSUM(iim,zsin,1)/pi
506	zrfi(1,l) = 0.
507
508	ENDDO
509	c$OMP END DO NOWAIT
510	endif
511
512
513	c 48. champs de vents aux pole sud:
514	c ---------------------------------
515	c U = 1 / pi * integrale [ v * cos(long) * d long ]
516	c V = 1 / pi * integrale [ v * sin(long) * d long ]
517
518	if (is_south_pole_dyn) then
519	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
520	DO l=1,llm
521
522	z1(1) =(rlonu(1)-rlonu(iim)+2.pi)pvcov(1,jjm,l)/cv(1,jjm)
523	DO i=2,iim
524	z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,jjm,l)/cv(i,jjm)
525	ENDDO
526
527	DO i=1,iim
528	zcos(i) = COS(rlonv(i))*z1(i)
529	zsin(i) = SIN(rlonv(i))*z1(i)
530	ENDDO
531
532	zufi(klon,l) = SSUM(iim,zcos,1)/pi
533	zvfi(klon,l) = SSUM(iim,zsin,1)/pi
534	zrfi(klon,l) = 0.
535	ENDDO
536	c$OMP END DO NOWAIT
537	endif
538
539	c On change de grille, dynamique vers physiq, pour le flux de masse verticale
540	CALL gr_dyn_fi_p(llm,iip1,jjp1,klon,flxw,flxwfi)
541
542	c-----------------------------------------------------------------------
543	c Appel de la physique:
544	c ---------------------
545
546
547	c$OMP BARRIER
548	if (first_omp) then
549	klon=klon_omp
550
551	allocate(zplev_omp(klon,llm+1))
552	allocate(zplay_omp(klon,llm))
553	allocate(zpk_omp(klon,llm))
554	allocate(zphi_omp(klon,llm))
555	allocate(zphis_omp(klon))
556	allocate(presnivs_omp(llm))
557	allocate(zufi_omp(klon,llm))
558	allocate(zvfi_omp(klon,llm))
559	allocate(zrfi_omp(klon,llm)) ! LG Ari 2014
560	allocate(ztfi_omp(klon,llm))
561	allocate(zqfi_omp(klon,llm,nqtot))
562	allocate(zdufi_omp(klon,llm))
563	allocate(zdvfi_omp(klon,llm))
564	allocate(zdtfi_omp(klon,llm))
565	allocate(zdqfi_omp(klon,llm,nqtot))
566	allocate(zdufic_omp(klon,llm))
567	allocate(zdvfic_omp(klon,llm))
568	allocate(zdtfic_omp(klon,llm))
569	allocate(zdqfic_omp(klon,llm,nqtot))
570	allocate(zdpsrf_omp(klon))
571	allocate(flxwfi_omp(klon,llm))
572	first_omp=.false.
573	endif
574
575
576	klon=klon_omp
577	offset=klon_omp_begin-1
578
579	do l=1,llm+1
580	do i=1,klon
581	zplev_omp(i,l)=zplev(offset+i,l)
582	enddo
583	enddo
584
585	do l=1,llm
586	do i=1,klon
587	zplay_omp(i,l)=zplay(offset+i,l)
588	enddo
589	enddo
590
591	do l=1,llm
592	do i=1,klon
593	zpk_omp(i,l)=zpk(offset+i,l)
594	enddo
595	enddo
596
597	do l=1,llm
598	do i=1,klon
599	zphi_omp(i,l)=zphi(offset+i,l)
600	enddo
601	enddo
602
603	do i=1,klon
604	zphis_omp(i)=zphis(offset+i)
605	enddo
606
607
608	do l=1,llm
609	presnivs_omp(l)=presnivs(l)
610	enddo
611
612	do l=1,llm
613	do i=1,klon
614	zufi_omp(i,l)=zufi(offset+i,l)
615	enddo
616	enddo
617
618	do l=1,llm
619	do i=1,klon
620	zvfi_omp(i,l)=zvfi(offset+i,l)
621	enddo
622	enddo
623
624	do l=1,llm
625	do i=1,klon
626	zrfi_omp(i,l)=zrfi(offset+i,l)
627	enddo
628	enddo
629
630
631	do l=1,llm
632	do i=1,klon
633	ztfi_omp(i,l)=ztfi(offset+i,l)
634	enddo
635	enddo
636
637	do iq=1,nqtot
638	do l=1,llm
639	do i=1,klon
640	zqfi_omp(i,l,iq)=zqfi(offset+i,l,iq)
641	enddo
642	enddo
643	enddo
644
645	do l=1,llm
646	do i=1,klon
647	zdufi_omp(i,l)=zdufi(offset+i,l)
648	enddo
649	enddo
650
651	do l=1,llm
652	do i=1,klon
653	zdvfi_omp(i,l)=zdvfi(offset+i,l)
654	enddo
655	enddo
656
657	do l=1,llm
658	do i=1,klon
659	zdtfi_omp(i,l)=zdtfi(offset+i,l)
660	enddo
661	enddo
662
663	do iq=1,nqtot
664	do l=1,llm
665	do i=1,klon
666	zdqfi_omp(i,l,iq)=zdqfi(offset+i,l,iq)
667	enddo
668	enddo
669	enddo
670
671	do i=1,klon
672	zdpsrf_omp(i)=zdpsrf(offset+i)
673	enddo
674
675	do l=1,llm
676	do i=1,klon
677	flxwfi_omp(i,l)=flxwfi(offset+i,l)
678	enddo
679	enddo
680
681	c$OMP BARRIER
682
683	!$OMP MASTER
684	! write(lunout,*) 'PHYSIQUE AVEC NSPLIT_PHYS=',nsplit_phys
685	!$OMP END MASTER
686	zdt_split=dtphys/nsplit_phys
687	zdufic_omp(:,:)=0.
688	zdvfic_omp(:,:)=0.
689	zdtfic_omp(:,:)=0.
690	zdqfic_omp(:,:,:)=0.
691
692	#ifdef CPP_PHYS
693	do isplit=1,nsplit_phys
694
695	jH_cur_split=jH_cur+(isplit-1) * dtvr / (daysec *nsplit_phys)
696	debut_split=debut.and.isplit==1
697	lafin_split=lafin.and.isplit==nsplit_phys
698
699	CALL call_physiq(klon,llm,nqtot,tname,
700	& debut_split,lafin_split,
701	& jD_cur,jH_cur_split,zdt_split,
702	& zplev_omp,zplay_omp,
703	& zpk_omp,zphi_omp,zphis_omp,
704	& presnivs_omp,
705	& zufi_omp,zvfi_omp,zrfi_omp,ztfi_omp,zqfi_omp,
706	& flxwfi_omp,pducov,
707	& zdufi_omp,zdvfi_omp,zdtfi_omp,zdqfi_omp,
708	& zdpsrf_omp)
709
710
711	zufi_omp(:,:)=zufi_omp(:,:)+zdufi_omp(:,:)*zdt_split
712	zvfi_omp(:,:)=zvfi_omp(:,:)+zdvfi_omp(:,:)*zdt_split
713	ztfi_omp(:,:)=ztfi_omp(:,:)+zdtfi_omp(:,:)*zdt_split
714	zqfi_omp(:,:,:)=zqfi_omp(:,:,:)+zdqfi_omp(:,:,:)*zdt_split
715
716	zdufic_omp(:,:)=zdufic_omp(:,:)+zdufi_omp(:,:)
717	zdvfic_omp(:,:)=zdvfic_omp(:,:)+zdvfi_omp(:,:)
718	zdtfic_omp(:,:)=zdtfic_omp(:,:)+zdtfi_omp(:,:)
719	zdqfic_omp(:,:,:)=zdqfic_omp(:,:,:)+zdqfi_omp(:,:,:)
720
721	enddo
722
723	#endif
724	! of #ifdef CPP_PHYS
725
726	zdufi_omp(:,:)=zdufic_omp(:,:)/nsplit_phys
727	zdvfi_omp(:,:)=zdvfic_omp(:,:)/nsplit_phys
728	zdtfi_omp(:,:)=zdtfic_omp(:,:)/nsplit_phys
729	zdqfi_omp(:,:,:)=zdqfic_omp(:,:,:)/nsplit_phys
730	c$OMP BARRIER
731
732	do l=1,llm+1
733	do i=1,klon
734	zplev(offset+i,l)=zplev_omp(i,l)
735	enddo
736	enddo
737
738	do l=1,llm
739	do i=1,klon
740	zplay(offset+i,l)=zplay_omp(i,l)
741	enddo
742	enddo
743
744	do l=1,llm
745	do i=1,klon
746	zphi(offset+i,l)=zphi_omp(i,l)
747	enddo
748	enddo
749
750
751	do i=1,klon
752	zphis(offset+i)=zphis_omp(i)
753	enddo
754
755
756	do l=1,llm
757	presnivs(l)=presnivs_omp(l)
758	enddo
759
760	do l=1,llm
761	do i=1,klon
762	zufi(offset+i,l)=zufi_omp(i,l)
763	enddo
764	enddo
765
766	do l=1,llm
767	do i=1,klon
768	zvfi(offset+i,l)=zvfi_omp(i,l)
769	enddo
770	enddo
771
772	do l=1,llm
773	do i=1,klon
774	ztfi(offset+i,l)=ztfi_omp(i,l)
775	enddo
776	enddo
777
778	do iq=1,nqtot
779	do l=1,llm
780	do i=1,klon
781	zqfi(offset+i,l,iq)=zqfi_omp(i,l,iq)
782	enddo
783	enddo
784	enddo
785
786	do l=1,llm
787	do i=1,klon
788	zdufi(offset+i,l)=zdufi_omp(i,l)
789	enddo
790	enddo
791
792	do l=1,llm
793	do i=1,klon
794	zdvfi(offset+i,l)=zdvfi_omp(i,l)
795	enddo
796	enddo
797
798	do l=1,llm
799	do i=1,klon
800	zdtfi(offset+i,l)=zdtfi_omp(i,l)
801	enddo
802	enddo
803
804	do iq=1,nqtot
805	do l=1,llm
806	do i=1,klon
807	zdqfi(offset+i,l,iq)=zdqfi_omp(i,l,iq)
808	enddo
809	enddo
810	enddo
811
812	do i=1,klon
813	zdpsrf(offset+i)=zdpsrf_omp(i)
814	enddo
815
816
817	klon=klon_mpi
818	500 CONTINUE
819	c$OMP BARRIER
820
821	c$OMP MASTER
822	call stop_timer(timer_physic)
823	c$OMP END MASTER
824
825	IF (using_mpi) THEN
826
827	if (MPI_rank>0) then
828
829	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
830	DO l=1,llm
831	du_send(1:iim,l)=zdufi(1:iim,l)
832	dv_send(1:iim,l)=zdvfi(1:iim,l)
833	ENDDO
834	c$OMP END DO NOWAIT
835
836	c$OMP BARRIER
837	#ifdef CPP_MPI
838	c$OMP MASTER
839	!$OMP CRITICAL (MPI)
840	call MPI_ISSEND(du_send,iim*llm,MPI_REAL8,MPI_Rank-1,401,
841	& COMM_LMDZ,Req(1),ierr)
842	call MPI_ISSEND(dv_send,iim*llm,MPI_REAL8,MPI_Rank-1,402,
843	& COMM_LMDZ,Req(2),ierr)
844	!$OMP END CRITICAL (MPI)
845	c$OMP END MASTER
846	#endif
847	c$OMP BARRIER
848
849	endif
850
851	if (MPI_rank<MPI_Size-1) then
852	c$OMP BARRIER
853	#ifdef CPP_MPI
854	c$OMP MASTER
855	!$OMP CRITICAL (MPI)
856	call MPI_IRECV(du_recv,iim*llm,MPI_REAL8,MPI_Rank+1,401,
857	& COMM_LMDZ,Req(3),ierr)
858	call MPI_IRECV(dv_recv,iim*llm,MPI_REAL8,MPI_Rank+1,402,
859	& COMM_LMDZ,Req(4),ierr)
860	!$OMP END CRITICAL (MPI)
861	c$OMP END MASTER
862	#endif
863	endif
864
865	c$OMP BARRIER
866
867
868	#ifdef CPP_MPI
869	c$OMP MASTER
870	!$OMP CRITICAL (MPI)
871	if (MPI_rank>0 .and. MPI_rank< MPI_Size-1) then
872	call MPI_WAITALL(4,Req(1),Status,ierr)
873	else if (MPI_rank>0) then
874	call MPI_WAITALL(2,Req(1),Status,ierr)
875	else if (MPI_rank <MPI_Size-1) then
876	call MPI_WAITALL(2,Req(3),Status,ierr)
877	endif
878	!$OMP END CRITICAL (MPI)
879	c$OMP END MASTER
880	#endif
881
882	c$OMP BARRIER
883
884	ENDIF ! using_mpi
885
886
887	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
888	DO l=1,llm
889
890	zdufi2(1:klon,l)=zdufi(1:klon,l)
891	zdufi2(klon+1:klon+iim,l)=du_recv(1:iim,l)
892
893	zdvfi2(1:klon,l)=zdvfi(1:klon,l)
894	zdvfi2(klon+1:klon+iim,l)=dv_recv(1:iim,l)
895
896	pdhfi(:,jj_begin,l)=0
897	pdqfi(:,jj_begin,l,:)=0
898	pdufi(:,jj_begin,l)=0
899	pdvfi(:,jj_begin,l)=0
900
901	if (.not. is_south_pole_dyn) then
902	pdhfi(:,jj_end,l)=0
903	pdqfi(:,jj_end,l,:)=0
904	pdufi(:,jj_end,l)=0
905	pdvfi(:,jj_end,l)=0
906	endif
907
908	ENDDO
909	c$OMP END DO NOWAIT
910
911	c$OMP MASTER
912	pdpsfi(:,jj_begin)=0
913	if (.not. is_south_pole_dyn) then
914	pdpsfi(:,jj_end)=0
915	endif
916	c$OMP END MASTER
917	c-----------------------------------------------------------------------
918	c transformation des tendances physiques en tendances dynamiques:
919	c ---------------------------------------------------------------
920
921	c tendance sur la pression :
922	c -----------------------------------
923	CALL gr_fi_dyn_p(1,klon,iip1,jjp1,zdpsrf,pdpsfi)
924	c
925	c 62. enthalpie potentielle
926	c ---------------------
927
928	kstart=1
929	kend=klon
930
931	if (is_north_pole_dyn) kstart=2
932	if (is_south_pole_dyn) kend=klon-1
933
934	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
935	DO l=1,llm
936
937	!CDIR ON_ADB(index_i)
938	!CDIR ON_ADB(index_j)
939	!cdir NODEP
940	do ig0=kstart,kend
941	i=index_i(ig0)
942	j=index_j(ig0)
943	pdhfi(i,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l)
944	if (i==1) pdhfi(iip1,j,l) = cpp * zdtfi(ig0,l) / ppk(i,j,l)
945	enddo
946
947	if (is_north_pole_dyn) then
948	DO i=1,iip1
949	pdhfi(i,1,l) = cpp * zdtfi(1,l) / ppk(i, 1 ,l)
950	enddo
951	endif
952
953	if (is_south_pole_dyn) then
954	DO i=1,iip1
955	pdhfi(i,jjp1,l) = cpp * zdtfi(klon,l)/ ppk(i,jjp1,l)
956	ENDDO
957	endif
958	ENDDO
959	c$OMP END DO NOWAIT
960
961	c 62. humidite specifique
962	c ---------------------
963	! Ehouarn: removed this useless bit: was overwritten at step 63 anyways
964	! DO iq=1,nqtot
965	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
966	! DO l=1,llm
967	!!!cdir NODEP
968	! do ig0=kstart,kend
969	! i=index_i(ig0)
970	! j=index_j(ig0)
971	! pdqfi(i,j,l,iq) = zdqfi(ig0,l,iq)
972	! if (i==1) pdqfi(iip1,j,l,iq) = zdqfi(ig0,l,iq)
973	! enddo
974	!
975	! if (is_north_pole_dyn) then
976	! do i=1,iip1
977	! pdqfi(i,1,l,iq) = zdqfi(1,l,iq)
978	! enddo
979	! endif
980	!
981	! if (is_south_pole_dyn) then
982	! do i=1,iip1
983	! pdqfi(i,jjp1,l,iq) = zdqfi(klon,l,iq)
984	! enddo
985	! endif
986	! ENDDO
987	!c$OMP END DO NOWAIT
988	! ENDDO
989
990	c 63. traceurs
991	c ------------
992	C initialisation des tendances
993
994	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
995	DO l=1,llm
996	pdqfi(:,:,l,:)=0.
997	ENDDO
998	c$OMP END DO NOWAIT
999
1000	C
1001	!cdir NODEP
1002	DO iq=1,nqtot
1003	iiq=niadv(iq)
1004	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1005	DO l=1,llm
1006	!CDIR ON_ADB(index_i)
1007	!CDIR ON_ADB(index_j)
1008	!cdir NODEP
1009	DO ig0=kstart,kend
1010	i=index_i(ig0)
1011	j=index_j(ig0)
1012	pdqfi(i,j,l,iiq) = zdqfi(ig0,l,iq)
1013	if (i==1) pdqfi(iip1,j,l,iiq) = zdqfi(ig0,l,iq)
1014	ENDDO
1015
1016	IF (is_north_pole_dyn) then
1017	DO i=1,iip1
1018	pdqfi(i,1,l,iiq) = zdqfi(1,l,iq)
1019	ENDDO
1020	ENDIF
1021
1022	IF (is_south_pole_dyn) then
1023	DO i=1,iip1
1024	pdqfi(i,jjp1,l,iiq) = zdqfi(klon,l,iq)
1025	ENDDO
1026	ENDIF
1027
1028	ENDDO
1029	c$OMP END DO NOWAIT
1030	ENDDO
1031
1032	c 65. champ u:
1033	c ------------
1034	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1035	DO l=1,llm
1036	!CDIR ON_ADB(index_i)
1037	!CDIR ON_ADB(index_j)
1038	!cdir NODEP
1039	do ig0=kstart,kend
1040	i=index_i(ig0)
1041	j=index_j(ig0)
1042
1043	if (i/=iim) then
1044	pdufi(i,j,l)=0.5(zdufi2(ig0,l)+zdufi2(ig0+1,l))cu(i,j)
1045	endif
1046
1047	if (i==1) then
1048	pdufi(iim,j,l)=0.5*( zdufi2(ig0,l)
1049	$ + zdufi2(ig0+iim-1,l))*cu(iim,j)
1050	pdufi(iip1,j,l)=0.5(zdufi2(ig0,l)+zdufi2(ig0+1,l))cu(i,j)
1051	endif
1052
1053	enddo
1054
1055	if (is_north_pole_dyn) then
1056	DO i=1,iip1
1057	pdufi(i,1,l) = 0.
1058	ENDDO
1059	endif
1060
1061	if (is_south_pole_dyn) then
1062	DO i=1,iip1
1063	pdufi(i,jjp1,l) = 0.
1064	ENDDO
1065	endif
1066
1067	ENDDO
1068	c$OMP END DO NOWAIT
1069
1070	c 67. champ v:
1071	c ------------
1072
1073	kstart=1
1074	kend=klon
1075
1076	if (is_north_pole_dyn) kstart=2
1077	if (is_south_pole_dyn) kend=klon-1-iim
1078
1079	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1080	DO l=1,llm
1081	!CDIR ON_ADB(index_i)
1082	!CDIR ON_ADB(index_j)
1083	!cdir NODEP
1084	do ig0=kstart,kend
1085	i=index_i(ig0)
1086	j=index_j(ig0)
1087	pdvfi(i,j,l)=0.5(zdvfi2(ig0,l)+zdvfi2(ig0+iim,l))cv(i,j)
1088	if (i==1) pdvfi(iip1,j,l) = 0.5*(zdvfi2(ig0,l)+
1089	$ zdvfi2(ig0+iim,l))
1090	$ *cv(i,j)
1091	enddo
1092
1093	ENDDO
1094	c$OMP END DO NOWAIT
1095
1096
1097	c 68. champ v pres des poles:
1098	c ---------------------------
1099	c v = U * cos(long) + V * SIN(long)
1100
1101	if (is_north_pole_dyn) then
1102
1103	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1104	DO l=1,llm
1105
1106	DO i=1,iim
1107	pdvfi(i,1,l)=
1108	$ zdufi(1,l)COS(rlonv(i))+zdvfi(1,l)SIN(rlonv(i))
1109
1110	pdvfi(i,1,l)=
1111	$ 0.5(pdvfi(i,1,l)+zdvfi(i+1,l))cv(i,1)
1112	ENDDO
1113
1114	pdvfi(iip1,1,l) = pdvfi(1,1,l)
1115
1116	ENDDO
1117	c$OMP END DO NOWAIT
1118
1119	endif
1120
1121	if (is_south_pole_dyn) then
1122
1123	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1124	DO l=1,llm
1125
1126	DO i=1,iim
1127	pdvfi(i,jjm,l)=zdufi(klon,l)*COS(rlonv(i))
1128	$ +zdvfi(klon,l)*SIN(rlonv(i))
1129
1130	pdvfi(i,jjm,l)=
1131	$ 0.5(pdvfi(i,jjm,l)+zdvfi(klon-iip1+i,l))cv(i,jjm)
1132	ENDDO
1133
1134	pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l)
1135
1136	ENDDO
1137	c$OMP END DO NOWAIT
1138
1139	endif
1140	c-----------------------------------------------------------------------
1141
1142	700 CONTINUE
1143
1144	firstcal = .FALSE.
1145
1146	#else
1147	write(lunout,*)
1148	& "calfis_p: for now can only work with parallel physics"
1149	stop
1150	#endif
1151	! of #ifdef CPP_PHYS
1152	#endif
1153	! of #ifdef CPP_PARA
1154	END

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

Download in other formats:

Original Format