Context Navigation

lmdz_calfis_loc.F90 @ 5072

Last change on this file since 5072 was 5066, checked in by abarral, 2 months ago
Transform gr_dyn_fi_p.F, gr_fi_dyn_p.F, calfis_loc.F into free-form modules. Reorder CPP_PARA keys in lmdz_call_calfis.F90, lmdz_calfis_loc.F90, lmdz_gr_dyn_fi_p.F90, lmdz_gr_fi_dyn_p.F90 to avoid implicit declarations. Remove redundant `-cpp -D..` on arch. Correct "!OMP" -> "!$OMP" Correct typo in lmdz_xios.F90, wstats.F90
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:keywords set to `Id`
File size: 29.4 KB

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

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

Download in other formats:

Original Format