Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

calfis_loc.F90 @ 5248

Last change on this file since 5248 was 5246, checked in by abarral, 21 hours ago
Convert fixed-form to free-form sources .F -> .{f,F}90 (WIP: some .F remain, will be handled in subsequent commits)
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.1 KB

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

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

Download in other formats:

Original Format