Context Navigation

calfis_loc.F90 @ 5272

Last change on this file since 5272 was 5272, checked in by abarral, 8 weeks ago
Turn paramet.h into a module
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.2 KB

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

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

Download in other formats:

Original Format