Context Navigation

calfis_p.F @ 2422

Last change on this file since 2422 was 2418, checked in by Ehouarn Millour, 9 years ago

Improving the physics/dynamics interface:

added module callphysiq_mod.F90 in dynphy_lonlat/phy* which contains the routine "call_physiq" which is called by calfis* and calls the physics. This way different "physiq" routine from different physics packages may be called: The calfis* routines now exposes all available fields that might be transmitted to physiq but which is actually send (ie: expected/needed by physiq) is decided in call_physiq.
turned "physiq.F90" into module "physiq_mod.F90" for better control of "physiq" arguments. Extracted embeded "gr_fi_ecrit" as self-standing routine (but note that this routine actually only works in serial mode).

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: 29.9 KB

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

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

Download in other formats:

Original Format