Context Navigation

leapfrog_loc.F @ 1763

Last change on this file since 1763 was 1705, checked in by Ehouarn Millour, 12 years ago

Added arch files for ADA (IDRIS IBMx3750) and made the following code modifications:
phylmd/printflag.F : removed "print" of unset variable (radpas0)
dyn3dmem/integrd_loc.F : removed unecessary "include mpif.h"
dyn3dmem/leapfrog_loc.F : removed unecessary "include mpif.h" and allocate saved variables at first call
dyn3dmem/mod_filtreg_p.F : added matmul() alternatives to call to BLAS routine SGEMM (which was incorectly set as DGEMM; which would fail if running with -r4)
filtrez/filtreg.F: changed calls to DGEMM into calls to SGEMM, so that code works with either -r4 or -r8 (the later being used in conjunction with "BLAS SGEMV=DGEMV SGEMM=DGEMM" preprocessing statements)
EM

File size: 51.7 KB

Line
1	!
2	! $Id$
3	!
4	c
5	c
6	#define DEBUG_IO
7	#undef DEBUG_IO
8
9
10	SUBROUTINE leapfrog_loc(ucov0,vcov0,teta0,ps0,
11	& masse0,phis0,q0,clesphy0,
12	& time_0)
13
14	USE misc_mod
15	USE parallel
16	USE times
17	USE mod_hallo
18	USE Bands
19	USE Write_Field
20	USE Write_Field_p
21	USE vampir
22	USE timer_filtre, ONLY : print_filtre_timer
23	USE infotrac
24	USE guide_loc_mod, ONLY : guide_main
25	USE getparam
26	USE control_mod
27	USE mod_filtreg_p
28	USE write_field_loc
29	USE allocate_field
30	USE call_dissip_mod, ONLY : call_dissip
31	USE call_calfis_mod, ONLY : call_calfis
32	USE leapfrog_mod
33	IMPLICIT NONE
34
35	c ...... Version du 10/01/98 ..........
36
37	c avec coordonnees verticales hybrides
38	c avec nouveaux operat. dissipation * ( gradiv2,divgrad2,nxgraro2 )
39
40	c=======================================================================
41	c
42	c Auteur: P. Le Van /L. Fairhead/F.Hourdin
43	c -------
44	c
45	c Objet:
46	c ------
47	c
48	c GCM LMD nouvelle grille
49	c
50	c=======================================================================
51	c
52	c ... Dans inigeom , nouveaux calculs pour les elongations cu , cv
53	c et possibilite d'appeler une fonction f(y) a derivee tangente
54	c hyperbolique a la place de la fonction a derivee sinusoidale.
55
56	c ... Possibilite de choisir le shema pour l'advection de
57	c q , en modifiant iadv dans traceur.def (10/02) .
58	c
59	c Pour Van-Leer + Vapeur d'eau saturee, iadv(1)=4. (F.Codron,10/99)
60	c Pour Van-Leer iadv=10
61	c
62	c-----------------------------------------------------------------------
63	c Declarations:
64	c -------------
65
66	#include "dimensions.h"
67	#include "paramet.h"
68	#include "comconst.h"
69	#include "comdissnew.h"
70	#include "comvert.h"
71	#include "comgeom.h"
72	#include "logic.h"
73	#include "temps.h"
74	#include "ener.h"
75	#include "description.h"
76	#include "serre.h"
77	!#include "com_io_dyn.h"
78	#include "iniprint.h"
79	#include "academic.h"
80	! include "mpif.h"
81
82	INTEGER longcles
83	PARAMETER ( longcles = 20 )
84	REAL clesphy0( longcles )
85
86	real zqmin,zqmax
87
88	c variables dynamiques
89	REAL :: vcov0(ijb_v:ije_v,llm),ucov0(ijb_u:ije_u,llm) ! vents covariants
90	REAL :: teta0(ijb_u:ije_u,llm) ! temperature potentielle
91	REAL :: q0(ijb_u:ije_u,llm,nqtot) ! champs advectes
92	REAL :: ps0(ijb_u:ije_u) ! pression au sol
93	REAL :: masse0(ijb_u:ije_u,llm) ! masse d'air
94	REAL :: phis0(ijb_u:ije_u) ! geopotentiel au sol
95
96	! REAL,SAVE,ALLOCATABLE :: p (:,: ) ! pression aux interfac.des couches
97	! REAL,SAVE,ALLOCATABLE :: pks(:) ! exner au sol
98	! REAL,SAVE,ALLOCATABLE :: pk(:,:) ! exner au milieu des couches
99	! REAL,SAVE,ALLOCATABLE :: pkf(:,:) ! exner filt.au milieu des couches
100	! REAL,SAVE,ALLOCATABLE :: phi(:,:) ! geopotentiel
101	! REAL,SAVE,ALLOCATABLE :: w(:,:) ! vitesse verticale
102
103	c variables dynamiques intermediaire pour le transport
104	! REAL,SAVE,ALLOCATABLE :: pbaru(:,:),pbarv(:,:) !flux de masse
105
106	c variables dynamiques au pas -1
107	! REAL,SAVE,ALLOCATABLE :: vcovm1(:,:),ucovm1(:,:)
108	! REAL,SAVE,ALLOCATABLE :: tetam1(:,:),psm1(:)
109	! REAL,SAVE,ALLOCATABLE :: massem1(:,:)
110
111	c tendances dynamiques
112	! REAL,SAVE,ALLOCATABLE :: dv(:,:),du(:,:)
113	! REAL,SAVE,ALLOCATABLE :: dteta(:,:),dp(:)
114	! REAL,DIMENSION(:,:,:), ALLOCATABLE, SAVE :: dq
115
116	c tendances de la dissipation
117	! REAL,SAVE,ALLOCATABLE :: dvdis(:,:),dudis(:,:)
118	! REAL,SAVE,ALLOCATABLE :: dtetadis(:,:)
119
120	c tendances physiques
121	REAL,SAVE,ALLOCATABLE :: dvfi(:,:),dufi(:,:)
122	REAL,SAVE,ALLOCATABLE :: dtetafi(:,:)
123	REAL,SAVE,ALLOCATABLE :: dpfi(:)
124	REAL,DIMENSION(:,:,:),ALLOCATABLE,SAVE :: dqfi
125
126	c variables pour le fichier histoire
127	REAL dtav ! intervalle de temps elementaire
128
129	REAL tppn(iim),tpps(iim),tpn,tps
130	c
131	INTEGER itau,itaufinp1,iav
132	! INTEGER iday ! jour julien
133	REAL time
134
135	REAL SSUM
136	REAL time_0
137	! REAL,SAVE,ALLOCATABLE :: finvmaold(:,:)
138
139	cym LOGICAL lafin
140	LOGICAL :: lafin
141	INTEGER ij,iq,l
142	INTEGER ik
143
144	real time_step, t_wrt, t_ops
145
146	! jD_cur: jour julien courant
147	! jH_cur: heure julienne courante
148	REAL :: jD_cur, jH_cur
149	INTEGER :: an, mois, jour
150	REAL :: secondes
151
152	logical :: physic
153	LOGICAL first,callinigrads
154
155	data callinigrads/.true./
156	character*10 string10
157
158	! REAL,SAVE,ALLOCATABLE :: alpha(:,:),beta(:,:)
159	! REAL,SAVE,ALLOCATABLE :: flxw(:,:) ! flux de masse verticale
160
161	c+jld variables test conservation energie
162	! REAL,SAVE,ALLOCATABLE :: ecin(:,:),ecin0(:,:)
163	C Tendance de la temp. potentiel d (theta)/ d t due a la
164	C tansformation d'energie cinetique en energie thermique
165	C cree par la dissipation
166	! REAL,SAVE,ALLOCATABLE :: dtetaecdt(:,:)
167	! REAL,SAVE,ALLOCATABLE :: vcont(:,:),ucont(:,:)
168	! REAL,SAVE,ALLOCATABLE :: vnat(:,:),unat(:,:)
169	REAL d_h_vcol, d_qt, d_qw, d_ql, d_ec
170	CHARACTER*15 ztit
171	!! INTEGER ip_ebil_dyn ! PRINT level for energy conserv. diag.
172	! SAVE ip_ebil_dyn
173	! DATA ip_ebil_dyn/0/
174	c-jld
175
176	character*80 dynhist_file, dynhistave_file
177	character(len=*),parameter :: modname="leapfrog_loc"
178	character*80 abort_message
179
180
181	logical,PARAMETER :: dissip_conservative=.TRUE.
182
183	INTEGER testita
184	PARAMETER (testita = 9)
185
186	logical , parameter :: flag_verif = .false.
187
188	c declaration liees au parallelisme
189	INTEGER :: ierr
190	LOGICAL :: FirstCaldyn
191	LOGICAL :: FirstPhysic
192	INTEGER :: ijb,ije,j,i
193	type(Request) :: TestRequest
194	type(Request) :: Request_Dissip
195	type(Request) :: Request_physic
196
197	INTEGER :: true_itau
198	INTEGER :: iapptrac
199	INTEGER :: AdjustCount
200	! INTEGER :: var_time
201	LOGICAL :: ok_start_timer=.FALSE.
202	LOGICAL, SAVE :: firstcall=.TRUE.
203	TYPE(distrib),SAVE :: new_dist
204
205	c$OMP MASTER
206	ItCount=0
207	c$OMP END MASTER
208	true_itau=0
209	FirstCaldyn=.TRUE.
210	FirstPhysic=.TRUE.
211	iapptrac=0
212	AdjustCount = 0
213	lafin=.false.
214
215	itaufin = nday*day_step
216	itaufinp1 = itaufin +1
217
218	itau = 0
219	physic=.true.
220	if (iflag_phys==0.or.iflag_phys==2) physic=.false.
221	CALL init_nan
222	CALL leapfrog_allocate
223	ucov=ucov0
224	vcov=vcov0
225	teta=teta0
226	ps=ps0
227	masse=masse0
228	phis=phis0
229	q=q0
230
231	! iday = day_ini+itau/day_step
232	! time = REAL(itau-(iday-day_ini)*day_step)/day_step+time_0
233	! IF(time.GT.1.) THEN
234	! time = time-1.
235	! iday = iday+1
236	! ENDIF
237
238	c Allocate variables depending on dynamic variable nqtot
239	!$OMP MASTER
240	if (firstcall) then
241	!
242	! ALLOCATE(p(ijb_u:ije_u,llmp1))
243	! ALLOCATE(pks(ijb_u:ije_u))
244	! ALLOCATE(pk(ijb_u:ije_u,llm))
245	! ALLOCATE(pkf(ijb_u:ije_u,llm))
246	! ALLOCATE(phi(ijb_u:ije_u,llm))
247	! ALLOCATE(w(ijb_u:ije_u,llm))
248	! ALLOCATE(pbaru(ip1jmp1,llm),pbarv(ip1jm,llm))
249	! ALLOCATE(vcovm1(ijb_v:ije_v,llm),ucovm1(ijb_u:ije_u,llm))
250	! ALLOCATE(tetam1(ijb_u:ije_u,llm),psm1(ijb_u:ije_u))
251	! ALLOCATE(massem1(ijb_u:ije_u,llm))
252	! ALLOCATE(dv(ijb_v:ije_v,llm),du(ijb_u:ije_u,llm))
253	! ALLOCATE(dteta(ijb_u:ije_u,llm),dp(ijb_u:ije_u))
254	! ALLOCATE(dvdis(ijb_v:ije_v,llm),dudis(ijb_u:ije_u,llm))
255	! ALLOCATE(dtetadis(ijb_u:ije_u,llm))
256	ALLOCATE(dvfi(ijb_v:ije_v,llm),dufi(ijb_u:ije_u,llm))
257	ALLOCATE(dtetafi(ijb_u:ije_u,llm))
258	ALLOCATE(dpfi(ijb_u:ije_u))
259	! ALLOCATE(dq(ijb_u:ije_u,llm,nqtot))
260	ALLOCATE(dqfi(ijb_u:ije_u,llm,nqtot))
261	! ALLOCATE(dqfi_tmp(iip1,llm,nqtot))
262	! ALLOCATE(finvmaold(ijb_u:ije_u,llm))
263	! ALLOCATE(alpha(ijb_u:ije_u,llm),beta(ijb_u:ije_u,llm))
264	! ALLOCATE(flxw(ijb_u:ije_u,llm))
265	! ALLOCATE(ecin(ijb_u:ije_u,llm),ecin0(ijb_u:ije_u,llm))
266	! ALLOCATE(dtetaecdt(ijb_u:ije_u,llm))
267	! ALLOCATE(vcont(ijb_v:ije_v,llm),ucont(ijb_u:ije_u,llm))
268	! ALLOCATE(vnat(ijb_v:ije_v,llm),unat(ijb_u:ije_u,llm))
269	endif
270	!$OMP END MASTER
271	!$OMP BARRIER
272
273	! CALL dynredem1_loc("restart.nc",0.0,
274	! & vcov,ucov,teta,q,masse,ps)
275
276
277	c-----------------------------------------------------------------------
278	c On initialise la pression et la fonction d'Exner :
279	c --------------------------------------------------
280
281	c$OMP MASTER
282	dq(:,:,:)=0.
283	CALL pression ( ijnb_u, ap, bp, ps, p )
284	c$OMP END MASTER
285	if (pressure_exner) then
286	CALL exner_hyb_loc( ijnb_u, ps, p,alpha,beta, pks, pk, pkf)
287	else
288	CALL exner_milieu_loc( ijnb_u, ps, p, beta, pks, pk, pkf )
289	endif
290	c-----------------------------------------------------------------------
291	c Debut de l'integration temporelle:
292	c ----------------------------------
293	c et du parallelisme !!
294
295	1 CONTINUE ! Matsuno Forward step begins here
296
297	jD_cur = jD_ref + day_ini - day_ref + &
298	& itau/day_step
299	jH_cur = jH_ref + start_time + &
300	& mod(itau,day_step)/float(day_step)
301	if (jH_cur > 1.0 ) then
302	jD_cur = jD_cur +1.
303	jH_cur = jH_cur -1.
304	endif
305
306
307	#ifdef CPP_IOIPSL
308	if (ok_guide) then
309	!$OMP MASTER
310	call guide_main(itau,ucov,vcov,teta,q,masse,ps)
311	!$OMP END MASTER
312	!$OMP BARRIER
313	endif
314	#endif
315
316
317	c
318	c IF( MOD( itau, 10* day_step ).EQ.0 ) THEN
319	c CALL test_period ( ucov,vcov,teta,q,p,phis )
320	c PRINT *,' ---- Test_period apres continue OK ! -----', itau
321	c ENDIF
322	c
323	cym CALL SCOPY( ijmllm ,vcov , 1, vcovm1 , 1 )
324	cym CALL SCOPY( ijp1llm,ucov , 1, ucovm1 , 1 )
325	cym CALL SCOPY( ijp1llm,teta , 1, tetam1 , 1 )
326	cym CALL SCOPY( ijp1llm,masse, 1, massem1, 1 )
327	cym CALL SCOPY( ip1jmp1, ps , 1, psm1 , 1 )
328
329	if (FirstCaldyn) then
330	c$OMP MASTER
331	ucovm1=ucov
332	vcovm1=vcov
333	tetam1= teta
334	massem1= masse
335	psm1= ps
336
337	! Ehouarn: finvmaold is actually not used
338	! finvmaold = masse
339	c$OMP END MASTER
340	c$OMP BARRIER
341	! CALL filtreg_p ( finvmaold ,jjb_u,jje_u,jjb_u,jje_u,jjp1, llm,
342	! & -2,2, .TRUE., 1 )
343	else
344	! Save fields obtained at previous time step as '...m1'
345	ijb=ij_begin
346	ije=ij_end
347
348	c$OMP MASTER
349	psm1 (ijb:ije) = ps (ijb:ije)
350	c$OMP END MASTER
351
352	c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
353	DO l=1,llm
354	ije=ij_end
355	ucovm1 (ijb:ije,l) = ucov (ijb:ije,l)
356	tetam1 (ijb:ije,l) = teta (ijb:ije,l)
357	massem1 (ijb:ije,l) = masse (ijb:ije,l)
358	! finvmaold(ijb:ije,l)=masse(ijb:ije,l)
359
360	if (pole_sud) ije=ij_end-iip1
361	vcovm1(ijb:ije,l) = vcov (ijb:ije,l)
362
363
364	ENDDO
365	c$OMP ENDDO
366
367
368	! Ehouarn: finvmaold not used
369	! CALL filtreg_p(finvmaold ,jjb_u,jje_u,jj_begin,jj_end,jjp1,
370	! . llm, -2,2, .TRUE., 1 )
371
372	endif ! of if (FirstCaldyn)
373
374	forward = .TRUE.
375	leapf = .FALSE.
376	dt = dtvr
377
378	c ... P.Le Van .26/04/94 ....
379
380	cym CALL SCOPY ( ijp1llm, masse, 1, finvmaold, 1 )
381	cym CALL filtreg ( finvmaold ,jjp1, llm, -2,2, .TRUE., 1 )
382
383	cym ne sert a rien
384	cym call minmax(ijp1llm,q(:,:,3),zqmin,zqmax)
385
386	2 CONTINUE ! Matsuno backward or leapfrog step begins here
387
388	c$OMP MASTER
389	ItCount=ItCount+1
390	if (MOD(ItCount,1)==1) then
391	debug=.true.
392	else
393	debug=.false.
394	endif
395	c$OMP END MASTER
396	c-----------------------------------------------------------------------
397
398	c date:
399	c -----
400
401
402	c gestion des appels de la physique et des dissipations:
403	c ------------------------------------------------------
404	c
405	c ... P.Le Van ( 6/02/95 ) ....
406
407	apphys = .FALSE.
408	statcl = .FALSE.
409	conser = .FALSE.
410	apdiss = .FALSE.
411
412	IF( purmats ) THEN
413	! Purely Matsuno time stepping
414	IF( MOD(itau,iconser) .EQ.0.AND. forward ) conser = .TRUE.
415	IF( MOD(itau,dissip_period ).EQ.0.AND..NOT.forward )
416	s apdiss = .TRUE.
417	IF( MOD(itau,iphysiq ).EQ.0.AND..NOT.forward
418	s .and. physic ) apphys = .TRUE.
419	ELSE
420	! Leapfrog/Matsuno time stepping
421	IF( MOD(itau ,iconser) .EQ. 0 ) conser = .TRUE.
422	IF( MOD(itau+1,dissip_period).EQ.0 .AND. .NOT. forward )
423	s apdiss = .TRUE.
424	IF( MOD(itau+1,iphysiq).EQ.0.AND.physic) apphys=.TRUE.
425	END IF
426
427	! Ehouarn: for Shallow Water case (ie: 1 vertical layer),
428	! supress dissipation step
429	if (llm.eq.1) then
430	apdiss=.false.
431	endif
432
433	cym ---> Pour le moment
434	cym apphys = .FALSE.
435	statcl = .FALSE.
436	conser = .FALSE. ! ie: no output of control variables to stdout in //
437
438	if (firstCaldyn) then
439	c$OMP MASTER
440	call Set_Distrib(distrib_caldyn)
441	c$OMP END MASTER
442	c$OMP BARRIER
443	firstCaldyn=.FALSE.
444	cym call InitTime
445	c$OMP MASTER
446	call Init_timer
447	c$OMP END MASTER
448	endif
449
450	c$OMP MASTER
451	IF (ok_start_timer) THEN
452	CALL InitTime
453	! ok_start_timer=.FALSE.
454	ok_start_timer=.TRUE.
455	ENDIF
456	c$OMP END MASTER
457
458
459	!ym PAS D'AJUSTEMENT POUR LE MOMENT
460	if (Adjust) then
461	AdjustCount=AdjustCount+1
462	! if (iapptrac==iapp_tracvl .and. (forward. OR . leapf)
463	! & .and. itau/iphysiq>2 .and. Adjustcount>30) then
464	if (Adjustcount>1) then
465	AdjustCount=0
466	c$OMP MASTER
467	call allgather_timer_average
468
469	if (prt_level > 9) then
470
471	print ,'********************************'
472	print ,'*** TIMER CALDYN ****'
473	do i=0,mpi_size-1
474	print *,'proc',i,' : Nb Bandes :',jj_nb_caldyn(i),
475	& ' : temps moyen :',
476	& timer_average(jj_nb_caldyn(i),timer_caldyn,i),
477	& '+-',timer_delta(jj_nb_caldyn(i),timer_caldyn,i)
478	enddo
479
480	print ,'********************************'
481	print ,'*** TIMER VANLEER ****'
482	do i=0,mpi_size-1
483	print *,'proc',i,' : Nb Bandes :',jj_nb_vanleer(i),
484	& ' : temps moyen :',
485	& timer_average(jj_nb_vanleer(i),timer_vanleer,i),
486	& '+-',timer_delta(jj_nb_vanleer(i),timer_vanleer,i)
487	enddo
488
489	print ,'********************************'
490	print ,'*** TIMER DISSIP ****'
491	do i=0,mpi_size-1
492	print *,'proc',i,' : Nb Bandes :',jj_nb_dissip(i),
493	& ' : temps moyen :',
494	& timer_average(jj_nb_dissip(i),timer_dissip,i),
495	& '+-',timer_delta(jj_nb_dissip(i),timer_dissip,i)
496	enddo
497
498	! if (mpi_rank==0) call WriteBands
499
500	endif
501
502	call AdjustBands_caldyn(new_dist)
503	!$OMP END MASTER
504	!$OMP BARRIER
505	CALL leapfrog_switch_caldyn(new_dist)
506	!$OMP BARRIER
507
508
509	!$OMP MASTER
510	distrib_caldyn=new_dist
511	CALL set_distrib(distrib_caldyn)
512	!$OMP END MASTER
513	!$OMP BARRIER
514	! call Register_SwapFieldHallo(ucov,ucov,ip1jmp1,llm,
515	! & jj_Nb_caldyn,0,0,TestRequest)
516	! call Register_SwapFieldHallo(ucovm1,ucovm1,ip1jmp1,llm,
517	! & jj_Nb_caldyn,0,0,TestRequest)
518	! call Register_SwapFieldHallo(vcov,vcov,ip1jm,llm,
519	! & jj_Nb_caldyn,0,0,TestRequest)
520	! call Register_SwapFieldHallo(vcovm1,vcovm1,ip1jm,llm,
521	! & jj_Nb_caldyn,0,0,TestRequest)
522	! call Register_SwapFieldHallo(teta,teta,ip1jmp1,llm,
523	! & jj_Nb_caldyn,0,0,TestRequest)
524	! call Register_SwapFieldHallo(tetam1,tetam1,ip1jmp1,llm,
525	! & jj_Nb_caldyn,0,0,TestRequest)
526	! call Register_SwapFieldHallo(masse,masse,ip1jmp1,llm,
527	! & jj_Nb_caldyn,0,0,TestRequest)
528	! call Register_SwapFieldHallo(massem1,massem1,ip1jmp1,llm,
529	! & jj_Nb_caldyn,0,0,TestRequest)
530	! call Register_SwapFieldHallo(ps,ps,ip1jmp1,1,
531	! & jj_Nb_caldyn,0,0,TestRequest)
532	! call Register_SwapFieldHallo(psm1,psm1,ip1jmp1,1,
533	! & jj_Nb_caldyn,0,0,TestRequest)
534	! call Register_SwapFieldHallo(pkf,pkf,ip1jmp1,llm,
535	! & jj_Nb_caldyn,0,0,TestRequest)
536	! call Register_SwapFieldHallo(pk,pk,ip1jmp1,llm,
537	! & jj_Nb_caldyn,0,0,TestRequest)
538	! call Register_SwapFieldHallo(pks,pks,ip1jmp1,1,
539	! & jj_Nb_caldyn,0,0,TestRequest)
540	! call Register_SwapFieldHallo(phis,phis,ip1jmp1,1,
541	! & jj_Nb_caldyn,0,0,TestRequest)
542	! call Register_SwapFieldHallo(phi,phi,ip1jmp1,llm,
543	! & jj_Nb_caldyn,0,0,TestRequest)
544	! call Register_SwapFieldHallo(finvmaold,finvmaold,ip1jmp1,llm,
545	! & jj_Nb_caldyn,0,0,TestRequest)
546	!
547	! do j=1,nqtot
548	! call Register_SwapFieldHallo(q(:,:,j),q(:,:,j),ip1jmp1,llm,
549	! & jj_nb_caldyn,0,0,TestRequest)
550	! enddo
551	!
552	! call Set_Distrib(distrib_caldyn)
553	! call SendRequest(TestRequest)
554	! call WaitRequest(TestRequest)
555
556	!$OMP MASTER
557	call AdjustBands_dissip(new_dist)
558	!$OMP END MASTER
559	!$OMP BARRIER
560	CALL leapfrog_switch_dissip(new_dist)
561	!$OMP BARRIER
562	!$OMP MASTER
563	distrib_dissip=new_dist
564	!$OMP END MASTER
565	!$OMP BARRIER
566	! call AdjustBands_physic
567
568	c$OMP MASTER
569	if (mpi_rank==0) call WriteBands
570	c$OMP END MASTER
571
572
573	endif
574	endif
575
576
577
578	c-----------------------------------------------------------------------
579	c calcul des tendances dynamiques:
580	c --------------------------------
581	c$OMP BARRIER
582	c$OMP MASTER
583	call VTb(VThallo)
584	c$OMP END MASTER
585
586	call Register_Hallo_u(ucov,llm,1,1,1,1,TestRequest)
587	call Register_Hallo_v(vcov,llm,1,1,1,1,TestRequest)
588	call Register_Hallo_u(teta,llm,1,1,1,1,TestRequest)
589	call Register_Hallo_u(ps,1,1,2,2,1,TestRequest)
590	call Register_Hallo_u(pkf,llm,1,1,1,1,TestRequest)
591	call Register_Hallo_u(pk,llm,1,1,1,1,TestRequest)
592	call Register_Hallo_u(pks,1,1,1,1,1,TestRequest)
593	call Register_Hallo_u(p,llmp1,1,1,1,1,TestRequest)
594
595	c do j=1,nqtot
596	c call Register_Hallo(q(1,1,j),ip1jmp1,llm,1,1,1,1,
597	c * TestRequest)
598	c enddo
599
600	call SendRequest(TestRequest)
601	c$OMP BARRIER
602	call WaitRequest(TestRequest)
603
604	c$OMP MASTER
605	call VTe(VThallo)
606	c$OMP END MASTER
607	c$OMP BARRIER
608
609	if (debug) then
610	call WriteField_u('ucov',ucov)
611	call WriteField_v('vcov',vcov)
612	call WriteField_u('teta',teta)
613	call WriteField_u('ps',ps)
614	call WriteField_u('masse',masse)
615	call WriteField_u('pk',pk)
616	call WriteField_u('pks',pks)
617	call WriteField_u('pkf',pkf)
618	call WriteField_u('phis',phis)
619	do j=1,nqtot
620	call WriteField_u('q'//trim(int2str(j)),
621	. q(:,:,j))
622	enddo
623	endif
624
625
626	True_itau=True_itau+1
627
628	c$OMP MASTER
629	PRINT *,"---> itau=",itau," True_itau=",True_itau
630	c$OMP END MASTER
631
632	c$OMP MASTER
633	IF (prt_level>9) THEN
634	WRITE(lunout,*)"leapfrog_p: Iteration No",True_itau
635	ENDIF
636
637
638	call start_timer(timer_caldyn)
639
640	! compute geopotential phi()
641	CALL geopot_loc ( ip1jmp1, teta , pk , pks, phis , phi )
642
643
644	call VTb(VTcaldyn)
645	c$OMP END MASTER
646	! var_time=time+iday-day_ini
647
648	c$OMP BARRIER
649	! CALL FTRACE_REGION_BEGIN("caldyn")
650	time = jD_cur + jH_cur
651	CALL caldyn_loc
652	$ ( itau,ucov,vcov,teta,ps,masse,pk,pkf,phis ,
653	$ phi,conser,du,dv,dteta,dp,w, pbaru,pbarv, time )
654
655	! CALL FTRACE_REGION_END("caldyn")
656
657	c$OMP MASTER
658	call VTe(VTcaldyn)
659	c$OMP END MASTER
660
661	#ifdef DEBUG_IO
662	call WriteField_u('du',du)
663	call WriteField_v('dv',dv)
664	call WriteField_u('dteta',dteta)
665	call WriteField_u('dp',dp)
666	call WriteField_u('w',w)
667	call WriteField_u('pbaru',pbaru)
668	call WriteField_v('pbarv',pbarv)
669	call WriteField_u('p',p)
670	call WriteField_u('masse',masse)
671	call WriteField_u('pk',pk)
672	#endif
673	c-----------------------------------------------------------------------
674	c calcul des tendances advection des traceurs (dont l'humidite)
675	c -------------------------------------------------------------
676
677
678	IF( forward. OR . leapf ) THEN
679
680
681	CALL caladvtrac_loc(q,pbaru,pbarv,
682	* p, masse, dq, teta,
683	. flxw,pk, iapptrac)
684
685	! do j=1,nqtot
686	! call WriteField_u('qadv'//trim(int2str(j)),q(:,:,j))
687	! enddo
688
689	c
690	ENDIF ! of IF( forward. OR . leapf )
691
692
693	c-----------------------------------------------------------------------
694	c integrations dynamique et traceurs:
695	c ----------------------------------
696
697	c$OMP MASTER
698	call VTb(VTintegre)
699	c$OMP END MASTER
700	#ifdef DEBUG_IO
701	if (true_itau>20) then
702	call WriteField_u('ucovm1',ucovm1)
703	call WriteField_v('vcovm1',vcovm1)
704	call WriteField_u('tetam1',tetam1)
705	call WriteField_u('psm1',psm1)
706	call WriteField_u('ucov_int',ucov)
707	call WriteField_v('vcov_int',vcov)
708	call WriteField_u('teta_int',teta)
709	call WriteField_u('ps_int',ps)
710	endif
711	#endif
712	c$OMP BARRIER
713	! CALL FTRACE_REGION_BEGIN("integrd")
714
715	CALL integrd_loc ( 2,vcovm1,ucovm1,tetam1,psm1,massem1 ,
716	$ dv,du,dteta,dq,dp,vcov,ucov,teta,q,ps,masse,phis)
717	! $ finvmaold )
718
719	! CALL FTRACE_REGION_END("integrd")
720	c$OMP BARRIER
721	#ifdef DEBUG_IO
722	call WriteField_u('ucovm1',ucovm1)
723	call WriteField_v('vcovm1',vcovm1)
724	call WriteField_u('tetam1',tetam1)
725	call WriteField_u('psm1',psm1)
726	call WriteField_u('ucov_int',ucov)
727	call WriteField_v('vcov_int',vcov)
728	call WriteField_u('teta_int',teta)
729	call WriteField_u('ps_int',ps)
730	#endif
731	c do j=1,nqtot
732	c call WriteField_p('q'//trim(int2str(j)),
733	c . reshape(q(:,:,j),(/iip1,jmp1,llm/)))
734	c call WriteField_p('dq'//trim(int2str(j)),
735	c . reshape(dq(:,:,j),(/iip1,jmp1,llm/)))
736	c enddo
737
738
739	c$OMP MASTER
740	call VTe(VTintegre)
741	c$OMP END MASTER
742	c .P.Le Van (26/04/94 ajout de finvpold dans l'appel d'integrd)
743	c
744	c-----------------------------------------------------------------------
745	c calcul des tendances physiques:
746	c -------------------------------
747	c ######## P.Le Van ( Modif le 6/02/95 ) ###########
748	c
749	IF( purmats ) THEN
750	IF( itau.EQ.itaufin.AND..NOT.forward ) lafin = .TRUE.
751	ELSE
752	IF( itau+1. EQ. itaufin ) lafin = .TRUE.
753	ENDIF
754
755	cc$OMP END PARALLEL
756
757	c
758	c
759	IF( apphys ) THEN
760
761	CALL call_calfis(itau,lafin,clesphy0,ucov,vcov,teta,masse,ps,
762	& phis,q,flxw)
763	! #ifdef DEBUG_IO
764	! call WriteField_u('ucovfi',ucov)
765	! call WriteField_v('vcovfi',vcov)
766	! call WriteField_u('tetafi',teta)
767	! call WriteField_u('pfi',p)
768	! call WriteField_u('pkfi',pk)
769	! do j=1,nqtot
770	! call WriteField_u('qfi'//trim(int2str(j)),q(:,:,j))
771	! enddo
772	! #endif
773	! c
774	! c ....... Ajout P.Le Van ( 17/04/96 ) ...........
775	! c
776	! cc$OMP PARALLEL DEFAULT(SHARED)
777	! cc$OMP+ PRIVATE(rdaym_ini,rdayvrai,ijb,ije)
778
779	! c$OMP MASTER
780	! call suspend_timer(timer_caldyn)
781
782	! write(lunout,*)
783	! & 'leapfrog_p: Entree dans la physique : Iteration No ',true_itau
784	! c$OMP END MASTER
785
786	! CALL pression_loc ( ip1jmp1, ap, bp, ps, p )
787
788	! c$OMP BARRIER
789	! CALL exner_hyb_loc( ip1jmp1, ps, p,alpha,beta,pks, pk, pkf )
790	! c$OMP BARRIER
791	! jD_cur = jD_ref + day_ini - day_ref
792	! $ + int (itau * dtvr / daysec)
793	! jH_cur = jH_ref + &
794	! & (itau * dtvr / daysec - int(itau * dtvr / daysec))
795	! ! call ju2ymds(jD_cur+jH_cur, an, mois, jour, secondes)
796
797	! c rajout debug
798	! c lafin = .true.
799
800
801	! c Inbterface avec les routines de phylmd (phymars ... )
802	! c -----------------------------------------------------
803
804	! c+jld
805
806	! c Diagnostique de conservation de l'energie : initialisation
807	!
808	! c-jld
809	! c$OMP BARRIER
810	! c$OMP MASTER
811	! call VTb(VThallo)
812	! c$OMP END MASTER
813
814	! #ifdef DEBUG_IO
815	! call WriteField_u('ucovfi',ucov)
816	! call WriteField_v('vcovfi',vcov)
817	! call WriteField_u('tetafi',teta)
818	! call WriteField_u('pfi',p)
819	! call WriteField_u('pkfi',pk)
820	! #endif
821	! call SetTag(Request_physic,800)
822	!
823	! call Register_SwapField_u(ucov,ucov,distrib_physic,
824	! * Request_physic,up=2,down=2)
825	!
826	! call Register_SwapField_v(vcov,vcov,distrib_physic,
827	! * Request_physic,up=2,down=2)
828
829	! call Register_SwapField_u(teta,teta,distrib_physic,
830	! * Request_physic,up=2,down=2)
831	!
832	! call Register_SwapField_u(masse,masse,distrib_physic,
833	! * Request_physic,up=1,down=2)
834
835	! call Register_SwapField_u(p,p,distrib_physic,
836	! * Request_physic,up=2,down=2)
837	!
838	! call Register_SwapField_u(pk,pk,distrib_physic,
839	! * Request_physic,up=2,down=2)
840	!
841	! call Register_SwapField_u(phis,phis,distrib_physic,
842	! * Request_physic,up=2,down=2)
843	!
844	! call Register_SwapField_u(phi,phi,distrib_physic,
845	! * Request_physic,up=2,down=2)
846	!
847	! call Register_SwapField_u(w,w,distrib_physic,
848	! * Request_physic,up=2,down=2)
849	!
850	! call Register_SwapField_u(q,q,distrib_physic,
851	! * Request_physic,up=2,down=2)
852
853	! call Register_SwapField_u(flxw,flxw,distrib_physic,
854	! * Request_physic,up=2,down=2)
855	!
856	! call SendRequest(Request_Physic)
857	! c$OMP BARRIER
858	! call WaitRequest(Request_Physic)
859
860	! c$OMP BARRIER
861	! c$OMP MASTER
862	! call Set_Distrib(distrib_Physic)
863	! call VTe(VThallo)
864	!
865	! call VTb(VTphysiq)
866	! c$OMP END MASTER
867	! c$OMP BARRIER
868
869	! #ifdef DEBUG_IO
870	! call WriteField_u('ucovfi',ucov)
871	! call WriteField_v('vcovfi',vcov)
872	! call WriteField_u('tetafi',teta)
873	! call WriteField_u('pfi',p)
874	! call WriteField_u('pkfi',pk)
875	! do j=1,nqtot
876	! call WriteField_u('qfi'//trim(int2str(j)),q(:,:,j))
877	! enddo
878	! #endif
879	! STOP
880	! c$OMP BARRIER
881	! ! CALL FTRACE_REGION_BEGIN("calfis")
882	! CALL calfis_loc(lafin ,jD_cur, jH_cur,
883	! $ ucov,vcov,teta,q,masse,ps,p,pk,phis,phi ,
884	! $ du,dv,dteta,dq,
885	! $ flxw,
886	! $ clesphy0, dufi,dvfi,dtetafi,dqfi,dpfi )
887	! ! CALL FTRACE_REGION_END("calfis")
888	! ! ijb=ij_begin
889	! ! ije=ij_end
890	! ! if ( .not. pole_nord) then
891	! !c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
892	! ! DO l=1,llm
893	! ! dufi_tmp(1:iip1,l) = dufi(ijb:ijb+iim,l)
894	! ! dvfi_tmp(1:iip1,l) = dvfi(ijb:ijb+iim,l)
895	! ! dtetafi_tmp(1:iip1,l)= dtetafi(ijb:ijb+iim,l)
896	! ! dqfi_tmp(1:iip1,l,:) = dqfi(ijb:ijb+iim,l,:)
897	! ! ENDDO
898	! !c$OMP END DO NOWAIT
899	! !
900	! !c$OMP MASTER
901	! ! dpfi_tmp(1:iip1) = dpfi(ijb:ijb+iim)
902	! !c$OMP END MASTER
903	! ! endif ! of if ( .not. pole_nord)
904
905	! !c$OMP BARRIER
906	! !c$OMP MASTER
907	! ! call Set_Distrib(distrib_physic_bis)
908
909	! ! call VTb(VThallo)
910	! !c$OMP END MASTER
911	! !c$OMP BARRIER
912	! !
913	! ! call Register_Hallo_u(dufi,llm,
914	! ! * 1,0,0,1,Request_physic)
915	! !
916	! ! call Register_Hallo_v(dvfi,llm,
917	! ! * 1,0,0,1,Request_physic)
918	! !
919	! ! call Register_Hallo_u(dtetafi,llm,
920	! ! * 1,0,0,1,Request_physic)
921	! !
922	! ! call Register_Hallo_u(dpfi,1,
923	! ! * 1,0,0,1,Request_physic)
924	! !
925	! ! do j=1,nqtot
926	! ! call Register_Hallo_u(dqfi(ijb_u,1,j),llm,
927	! ! * 1,0,0,1,Request_physic)
928	! ! enddo
929	! !
930	! ! call SendRequest(Request_Physic)
931	! !c$OMP BARRIER
932	! ! call WaitRequest(Request_Physic)
933	! !
934	! !c$OMP BARRIER
935	! !c$OMP MASTER
936	! ! call VTe(VThallo)
937	! !
938	! ! call set_Distrib(distrib_Physic)
939	! !c$OMP END MASTER
940	! !c$OMP BARRIER
941	! ! ijb=ij_begin
942	! ! if (.not. pole_nord) then
943	! !
944	! !c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
945	! ! DO l=1,llm
946	! ! dufi(ijb:ijb+iim,l) = dufi(ijb:ijb+iim,l)+dufi_tmp(1:iip1,l)
947	! ! dvfi(ijb:ijb+iim,l) = dvfi(ijb:ijb+iim,l)+dvfi_tmp(1:iip1,l)
948	! ! dtetafi(ijb:ijb+iim,l) = dtetafi(ijb:ijb+iim,l)
949	! ! & +dtetafi_tmp(1:iip1,l)
950	! ! dqfi(ijb:ijb+iim,l,:) = dqfi(ijb:ijb+iim,l,:)
951	! ! & + dqfi_tmp(1:iip1,l,:)
952	! ! ENDDO
953	! !c$OMP END DO NOWAIT
954	! !
955	! !c$OMP MASTER
956	! ! dpfi(ijb:ijb+iim) = dpfi(ijb:ijb+iim)+ dpfi_tmp(1:iip1)
957	! !c$OMP END MASTER
958	! !
959	! ! endif ! of if (.not. pole_nord)
960
961	! #ifdef DEBUG_IO
962	! call WriteField_u('dufi',dufi)
963	! call WriteField_v('dvfi',dvfi)
964	! call WriteField_u('dtetafi',dtetafi)
965	! call WriteField_u('dpfi',dpfi)
966	! do j=1,nqtot
967	! call WriteField_u('dqfi'//trim(int2str(j)),dqfi(:,:,j))
968	! enddo
969	! #endif
970
971	! c$OMP BARRIER
972
973	! c ajout des tendances physiques:
974	! c ------------------------------
975	! #ifdef DEBUG_IO
976	! call WriteField_u('ucovfi',ucov)
977	! call WriteField_v('vcovfi',vcov)
978	! call WriteField_u('tetafi',teta)
979	! call WriteField_u('psfi',ps)
980	! do j=1,nqtot
981	! call WriteField_u('qfi'//trim(int2str(j)),q(:,:,j))
982	! enddo
983	! #endif
984
985	! IF (ok_strato) THEN
986	! CALL top_bound_loc( vcov,ucov,teta,masse,dufi,dvfi,dtetafi)
987	! ENDIF
988
989	! #ifdef DEBUG_IO
990	! call WriteField_u('ucovfi',ucov)
991	! call WriteField_v('vcovfi',vcov)
992	! call WriteField_u('tetafi',teta)
993	! call WriteField_u('psfi',ps)
994	! do j=1,nqtot
995	! call WriteField_u('qfi'//trim(int2str(j)),q(:,:,j))
996	! enddo
997	! #endif
998
999	! CALL addfi_loc( dtphys, leapf, forward ,
1000	! $ ucov, vcov, teta , q ,ps ,
1001	! $ dufi, dvfi, dtetafi , dqfi ,dpfi )
1002
1003	! #ifdef DEBUG_IO
1004	! call WriteField_u('ucovfi',ucov)
1005	! call WriteField_v('vcovfi',vcov)
1006	! call WriteField_u('tetafi',teta)
1007	! call WriteField_u('psfi',ps)
1008	! do j=1,nqtot
1009	! call WriteField_u('qfi'//trim(int2str(j)),q(:,:,j))
1010	! enddo
1011	! #endif
1012
1013	! c$OMP BARRIER
1014	! c$OMP MASTER
1015	! call VTe(VTphysiq)
1016
1017	! call VTb(VThallo)
1018	! c$OMP END MASTER
1019
1020	! call SetTag(Request_physic,800)
1021	! call Register_SwapField_u(ucov,ucov,
1022	! * distrib_caldyn,Request_physic)
1023	!
1024	! call Register_SwapField_v(vcov,vcov,
1025	! * distrib_caldyn,Request_physic)
1026	!
1027	! call Register_SwapField_u(teta,teta,
1028	! * distrib_caldyn,Request_physic)
1029	!
1030	! call Register_SwapField_u(masse,masse,
1031	! * distrib_caldyn,Request_physic)
1032
1033	! call Register_SwapField_u(p,p,
1034	! * distrib_caldyn,Request_physic)
1035	!
1036	! call Register_SwapField_u(pk,pk,
1037	! * distrib_caldyn,Request_physic)
1038	!
1039	! call Register_SwapField_u(phis,phis,
1040	! * distrib_caldyn,Request_physic)
1041	!
1042	! call Register_SwapField_u(phi,phi,
1043	! * distrib_caldyn,Request_physic)
1044	!
1045	! call Register_SwapField_u(w,w,
1046	! * distrib_caldyn,Request_physic)
1047
1048	! call Register_SwapField_u(q,q,
1049	! * distrib_caldyn,Request_physic)
1050	!
1051	! call SendRequest(Request_Physic)
1052	! c$OMP BARRIER
1053	! call WaitRequest(Request_Physic)
1054
1055	! c$OMP BARRIER
1056	! c$OMP MASTER
1057	! call VTe(VThallo)
1058	! call set_distrib(distrib_caldyn)
1059	! c$OMP END MASTER
1060	! c$OMP BARRIER
1061	! c
1062	! c Diagnostique de conservation de l'energie : difference
1063	! IF (ip_ebil_dyn.ge.1 ) THEN
1064	! ztit='bil phys'
1065	! CALL diagedyn(ztit,2,1,1,dtphys
1066	! e , ucov , vcov , ps, p ,pk , teta , q(:,:,1), q(:,:,2))
1067	! ENDIF
1068
1069	! #ifdef DEBUG_IO
1070	! call WriteField_u('ucovfi',ucov)
1071	! call WriteField_v('vcovfi',vcov)
1072	! call WriteField_u('tetafi',teta)
1073	! call WriteField_u('psfi',ps)
1074	! do j=1,nqtot
1075	! call WriteField_u('qfi'//trim(int2str(j)),q(:,:,j))
1076	! enddo
1077	! #endif
1078
1079
1080	! c-jld
1081	c$OMP MASTER
1082	if (FirstPhysic) then
1083	ok_start_timer=.TRUE.
1084	FirstPhysic=.false.
1085	endif
1086	c$OMP END MASTER
1087	ENDIF ! of IF( apphys )
1088
1089	IF(iflag_phys.EQ.2) THEN ! "Newtonian" case
1090	c Calcul academique de la physique = Rappel Newtonien + fritcion
1091	c --------------------------------------------------------------
1092	cym teta(:,:)=teta(:,:)
1093	cym s -iphysiqdtvr(teta(:,:)-tetarappel(:,:))/taurappel
1094	ijb=ij_begin
1095	ije=ij_end
1096	!LF teta(ijb:ije,:)=teta(ijb:ije,:)
1097	!LF s -iphysiqdtvr(teta(ijb:ije,:)-tetarappel(ijb:ije,:))/taurappel
1098	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1099	do l=1,llm
1100	teta(ijb:ije,l)=teta(ijb:ije,l) -dtvr*
1101	& (teta(ijb:ije,l)-tetarappel(ijb:ije,l))*
1102	& (knewt_g+knewt_t(l)*clat4(ijb:ije))
1103	enddo
1104	!$OMP END DO
1105
1106	!$OMP MASTER
1107	if (planet_type.eq."giant") then
1108	! add an intrinsic heat flux at the base of the atmosphere
1109	teta(ijb:ije,1) = teta(ijb:ije,1)
1110	& + dtvr * aire(ijb:ije) * ihf / cpp / masse(ijb:ije,1)
1111	endif
1112	!$OMP END MASTER
1113	!$OMP BARRIER
1114
1115
1116	call Register_Hallo_u(ucov,llm,0,1,1,0,Request_Physic)
1117	call Register_Hallo_v(vcov,llm,1,1,1,1,Request_Physic)
1118	call SendRequest(Request_Physic)
1119	c$OMP BARRIER
1120	call WaitRequest(Request_Physic)
1121	c$OMP BARRIER
1122	call friction_loc(ucov,vcov,dtvr)
1123	!$OMP BARRIER
1124
1125	! Sponge layer (if any)
1126	IF (ok_strato) THEN
1127	! set dufi,dvfi,... to zero
1128	ijb=ij_begin
1129	ije=ij_end
1130	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1131	do l=1,llm
1132	dufi(ijb:ije,l)=0
1133	dtetafi(ijb:ije,l)=0
1134	dqfi(ijb:ije,l,1:nqtot)=0
1135	enddo
1136	!$OMP END DO
1137	!$OMP MASTER
1138	dpfi(ijb:ije)=0
1139	!$OMP END MASTER
1140	ijb=ij_begin
1141	ije=ij_end
1142	if (pole_sud) ije=ije-iip1
1143	!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1144	do l=1,llm
1145	dvfi(ijb:ije,l)=0
1146	enddo
1147	!$OMP END DO
1148
1149	CALL top_bound_loc(vcov,ucov,teta,masse,dufi,dvfi,dtetafi)
1150	CALL addfi_loc( dtvr, leapf, forward ,
1151	$ ucov, vcov, teta , q ,ps ,
1152	$ dufi, dvfi, dtetafi , dqfi ,dpfi )
1153	!$OMP BARRIER
1154	ENDIF ! of IF (ok_strato)
1155	ENDIF ! of IF(iflag_phys.EQ.2)
1156
1157
1158	CALL pression_loc ( ip1jmp1, ap, bp, ps, p )
1159	c$OMP BARRIER
1160	if (pressure_exner) then
1161	CALL exner_hyb_loc( ijnb_u, ps, p,alpha,beta, pks, pk, pkf )
1162	else
1163	CALL exner_milieu_loc( ijnb_u, ps, p, beta, pks, pk, pkf )
1164	endif
1165	c$OMP BARRIER
1166
1167	cc$OMP END PARALLEL
1168
1169	c-----------------------------------------------------------------------
1170	c dissipation horizontale et verticale des petites echelles:
1171	c ----------------------------------------------------------
1172
1173	IF(apdiss) THEN
1174
1175	CALL call_dissip(ucov,vcov,teta,p,pk,ps)
1176	!cc$OMP PARALLEL DEFAULT(SHARED)
1177	!cc$OMP+ PRIVATE(ijb,ije,tppn,tpn,tpps,tps)
1178	!c$OMP MASTER
1179	! call suspend_timer(timer_caldyn)
1180	!
1181	!c print*,'Entree dans la dissipation : Iteration No ',true_itau
1182	!c calcul de l'energie cinetique avant dissipation
1183	!c print *,'Passage dans la dissipation'
1184
1185	! call VTb(VThallo)
1186	!c$OMP END MASTER
1187
1188	!c$OMP BARRIER
1189
1190	! call Register_SwapField_u(ucov,ucov,distrib_dissip,
1191	! * Request_dissip,up=1,down=1)
1192
1193	! call Register_SwapField_v(vcov,vcov,distrib_dissip,
1194	! * Request_dissip,up=1,down=1)
1195
1196	! call Register_SwapField_u(teta,teta,distrib_dissip,
1197	! * Request_dissip)
1198
1199	! call Register_SwapField_u(p,p,distrib_dissip,
1200	! * Request_dissip)
1201
1202	! call Register_SwapField_u(pk,pk,distrib_dissip,
1203	! * Request_dissip)
1204
1205	! call SendRequest(Request_dissip)
1206	!c$OMP BARRIER
1207	! call WaitRequest(Request_dissip)
1208
1209	!c$OMP BARRIER
1210	!c$OMP MASTER
1211	! call set_distrib(distrib_dissip)
1212	! call VTe(VThallo)
1213	! call VTb(VTdissipation)
1214	! call start_timer(timer_dissip)
1215	!c$OMP END MASTER
1216	!c$OMP BARRIER
1217
1218	! call covcont_loc(llm,ucov,vcov,ucont,vcont)
1219	! call enercin_loc(vcov,ucov,vcont,ucont,ecin0)
1220
1221	!c dissipation
1222
1223	!! CALL FTRACE_REGION_BEGIN("dissip")
1224	! CALL dissip_loc(vcov,ucov,teta,p,dvdis,dudis,dtetadis)
1225
1226	!#ifdef DEBUG_IO
1227	! call WriteField_u('dudis',dudis)
1228	! call WriteField_v('dvdis',dvdis)
1229	! call WriteField_u('dtetadis',dtetadis)
1230	!#endif
1231	!
1232	!! CALL FTRACE_REGION_END("dissip")
1233	!
1234	! ijb=ij_begin
1235	! ije=ij_end
1236	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1237	! DO l=1,llm
1238	! ucov(ijb:ije,l)=ucov(ijb:ije,l)+dudis(ijb:ije,l)
1239	! ENDDO
1240	!c$OMP END DO NOWAIT
1241	! if (pole_sud) ije=ije-iip1
1242	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1243	! DO l=1,llm
1244	! vcov(ijb:ije,l)=vcov(ijb:ije,l)+dvdis(ijb:ije,l)
1245	! ENDDO
1246	!c$OMP END DO NOWAIT
1247
1248	!c teta=teta+dtetadis
1249
1250
1251	!c------------------------------------------------------------------------
1252	! if (dissip_conservative) then
1253	!C On rajoute la tendance due a la transform. Ec -> E therm. cree
1254	!C lors de la dissipation
1255	!c$OMP BARRIER
1256	!c$OMP MASTER
1257	! call suspend_timer(timer_dissip)
1258	! call VTb(VThallo)
1259	!c$OMP END MASTER
1260	! call Register_Hallo_u(ucov,llm,1,1,1,1,Request_Dissip)
1261	! call Register_Hallo_v(vcov,llm,1,1,1,1,Request_Dissip)
1262	! call SendRequest(Request_Dissip)
1263	!c$OMP BARRIER
1264	! call WaitRequest(Request_Dissip)
1265	!c$OMP MASTER
1266	! call VTe(VThallo)
1267	! call resume_timer(timer_dissip)
1268	!c$OMP END MASTER
1269	!c$OMP BARRIER
1270	! call covcont_loc(llm,ucov,vcov,ucont,vcont)
1271	! call enercin_loc(vcov,ucov,vcont,ucont,ecin)
1272	!
1273	! ijb=ij_begin
1274	! ije=ij_end
1275	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1276	! do l=1,llm
1277	! do ij=ijb,ije
1278	! dtetaecdt(ij,l)= (ecin0(ij,l)-ecin(ij,l))/ pk(ij,l)
1279	! dtetadis(ij,l)=dtetadis(ij,l)+dtetaecdt(ij,l)
1280	! enddo
1281	! enddo
1282	!c$OMP END DO NOWAIT
1283	! endif
1284
1285	! ijb=ij_begin
1286	! ije=ij_end
1287	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1288	! do l=1,llm
1289	! do ij=ijb,ije
1290	! teta(ij,l)=teta(ij,l)+dtetadis(ij,l)
1291	! enddo
1292	! enddo
1293	!c$OMP END DO NOWAIT
1294	!c------------------------------------------------------------------------
1295
1296
1297	!c ....... P. Le Van ( ajout le 17/04/96 ) ...........
1298	!c ... Calcul de la valeur moyenne, unique de h aux poles .....
1299	!c
1300
1301	! ijb=ij_begin
1302	! ije=ij_end
1303	!
1304	! if (pole_nord) then
1305	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1306	! DO l = 1, llm
1307	! DO ij = 1,iim
1308	! tppn(ij) = aire( ij ) * teta( ij ,l)
1309	! ENDDO
1310	! tpn = SSUM(iim,tppn,1)/apoln
1311
1312	! DO ij = 1, iip1
1313	! teta( ij ,l) = tpn
1314	! ENDDO
1315	! ENDDO
1316	!c$OMP END DO NOWAIT
1317
1318	!c$OMP MASTER
1319	! DO ij = 1,iim
1320	! tppn(ij) = aire( ij ) * ps ( ij )
1321	! ENDDO
1322	! tpn = SSUM(iim,tppn,1)/apoln
1323	!
1324	! DO ij = 1, iip1
1325	! ps( ij ) = tpn
1326	! ENDDO
1327	!c$OMP END MASTER
1328	! endif
1329	!
1330	! if (pole_sud) then
1331	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1332	! DO l = 1, llm
1333	! DO ij = 1,iim
1334	! tpps(ij) = aire(ij+ip1jm) * teta(ij+ip1jm,l)
1335	! ENDDO
1336	! tps = SSUM(iim,tpps,1)/apols
1337
1338	! DO ij = 1, iip1
1339	! teta(ij+ip1jm,l) = tps
1340	! ENDDO
1341	! ENDDO
1342	!c$OMP END DO NOWAIT
1343
1344	!c$OMP MASTER
1345	! DO ij = 1,iim
1346	! tpps(ij) = aire(ij+ip1jm) * ps (ij+ip1jm)
1347	! ENDDO
1348	! tps = SSUM(iim,tpps,1)/apols
1349	!
1350	! DO ij = 1, iip1
1351	! ps(ij+ip1jm) = tps
1352	! ENDDO
1353	!c$OMP END MASTER
1354	! endif
1355
1356
1357	!c$OMP BARRIER
1358	!c$OMP MASTER
1359	! call VTe(VTdissipation)
1360
1361	! call stop_timer(timer_dissip)
1362	!
1363	! call VTb(VThallo)
1364	!c$OMP END MASTER
1365	! call Register_SwapField_u(ucov,ucov,distrib_caldyn,
1366	! * Request_dissip)
1367
1368	! call Register_SwapField_v(vcov,vcov,distrib_caldyn,
1369	! * Request_dissip)
1370
1371	! call Register_SwapField_u(teta,teta,distrib_caldyn,
1372	! * Request_dissip)
1373
1374	! call Register_SwapField_u(p,p,distrib_caldyn,
1375	! * Request_dissip)
1376
1377	! call Register_SwapField_u(pk,pk,distrib_caldyn,
1378	! * Request_dissip)
1379
1380	! call SendRequest(Request_dissip)
1381	!c$OMP BARRIER
1382	! call WaitRequest(Request_dissip)
1383
1384	!c$OMP BARRIER
1385	!c$OMP MASTER
1386	! call set_distrib(distrib_caldyn)
1387	! call VTe(VThallo)
1388	! call resume_timer(timer_caldyn)
1389	!c print *,'fin dissipation'
1390	!c$OMP END MASTER
1391	!c$OMP BARRIER
1392	END IF ! of IF(apdiss)
1393
1394	cc$OMP END PARALLEL
1395
1396	c ajout debug
1397	c IF( lafin ) then
1398	c abort_message = 'Simulation finished'
1399	c call abort_gcm(modname,abort_message,0)
1400	c ENDIF
1401
1402	c ********************************************************************
1403	c ********************************************************************
1404	c .... fin de l'integration dynamique et physique pour le pas itau ..
1405	c ********************************************************************
1406	c ********************************************************************
1407
1408	c preparation du pas d'integration suivant ......
1409	cym call WriteField('ucov',reshape(ucov,(/iip1,jmp1,llm/)))
1410	cym call WriteField('vcov',reshape(vcov,(/iip1,jjm,llm/)))
1411	c$OMP MASTER
1412	call stop_timer(timer_caldyn)
1413	c$OMP END MASTER
1414	IF (itau==itaumax) then
1415	c$OMP MASTER
1416	call allgather_timer_average
1417
1418	if (mpi_rank==0) then
1419
1420	print ,'********************************'
1421	print ,'*** TIMER CALDYN ****'
1422	do i=0,mpi_size-1
1423	print *,'proc',i,' : Nb Bandes :',jj_nb_caldyn(i),
1424	& ' : temps moyen :',
1425	& timer_average(jj_nb_caldyn(i),timer_caldyn,i)
1426	enddo
1427
1428	print ,'********************************'
1429	print ,'*** TIMER VANLEER ****'
1430	do i=0,mpi_size-1
1431	print *,'proc',i,' : Nb Bandes :',jj_nb_vanleer(i),
1432	& ' : temps moyen :',
1433	& timer_average(jj_nb_vanleer(i),timer_vanleer,i)
1434	enddo
1435
1436	print ,'********************************'
1437	print ,'*** TIMER DISSIP ****'
1438	do i=0,mpi_size-1
1439	print *,'proc',i,' : Nb Bandes :',jj_nb_dissip(i),
1440	& ' : temps moyen :',
1441	& timer_average(jj_nb_dissip(i),timer_dissip,i)
1442	enddo
1443
1444	print ,'********************************'
1445	print ,'*** TIMER PHYSIC ****'
1446	do i=0,mpi_size-1
1447	print *,'proc',i,' : Nb Bandes :',jj_nb_physic(i),
1448	& ' : temps moyen :',
1449	& timer_average(jj_nb_physic(i),timer_physic,i)
1450	enddo
1451
1452	endif
1453
1454	print ,'Taille du Buffer MPI (REAL8)',MaxBufferSize
1455	print ,'Taille du Buffer MPI utilise (REAL8)',MaxBufferSize_Used
1456	print *, 'Temps total ecoule sur la parallelisation :',DiffTime()
1457	print *, 'Temps CPU ecoule sur la parallelisation :',DiffCpuTime()
1458	CALL print_filtre_timer
1459	c$OMP END MASTER
1460	CALL dynredem1_loc("restart.nc",0.0,
1461	. vcov,ucov,teta,q,masse,ps)
1462	c$OMP MASTER
1463	call fin_getparam
1464	call finalize_parallel
1465	c$OMP END MASTER
1466	c$OMP BARRIER
1467	RETURN
1468	ENDIF
1469
1470	IF ( .NOT.purmats ) THEN
1471	c ........................................................
1472	c .............. schema matsuno + leapfrog ..............
1473	c ........................................................
1474
1475	IF(forward. OR. leapf) THEN
1476	itau= itau + 1
1477	! iday= day_ini+itau/day_step
1478	! time= REAL(itau-(iday-day_ini)*day_step)/day_step+time_0
1479	! IF(time.GT.1.) THEN
1480	! time = time-1.
1481	! iday = iday+1
1482	! ENDIF
1483	ENDIF
1484
1485
1486	IF( itau. EQ. itaufinp1 ) then
1487
1488	if (flag_verif) then
1489	write(79,*) 'ucov',ucov
1490	write(80,*) 'vcov',vcov
1491	write(81,*) 'teta',teta
1492	write(82,*) 'ps',ps
1493	write(83,*) 'q',q
1494	WRITE(85,*) 'q1 = ',q(:,:,1)
1495	WRITE(86,*) 'q3 = ',q(:,:,3)
1496	endif
1497
1498
1499	c$OMP MASTER
1500	call fin_getparam
1501	call finalize_parallel
1502	c$OMP END MASTER
1503	abort_message = 'Simulation finished'
1504	call abort_gcm(modname,abort_message,0)
1505	RETURN
1506	ENDIF
1507	c-----------------------------------------------------------------------
1508	c ecriture du fichier histoire moyenne:
1509	c -------------------------------------
1510
1511	IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) THEN
1512	c$OMP BARRIER
1513	IF(itau.EQ.itaufin) THEN
1514	iav=1
1515	ELSE
1516	iav=0
1517	ENDIF
1518
1519	#ifdef CPP_IOIPSL
1520	IF (ok_dynzon) THEN
1521
1522	CALL bilan_dyn_loc(2,dtvriperiod,dtvrday_step*periodav,
1523	, ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,q)
1524
1525	ENDIF !ok_dynzon
1526
1527	IF (ok_dyn_ave) THEN
1528	CALL writedynav_loc(itau,vcov,
1529	& ucov,teta,pk,phi,q,masse,ps,phis)
1530	ENDIF
1531	#endif
1532	ENDIF
1533
1534	c-----------------------------------------------------------------------
1535	c ecriture de la bande histoire:
1536	c ------------------------------
1537
1538	IF( MOD(itau,iecri).EQ.0) THEN
1539	! Ehouarn: output only during LF or Backward Matsuno
1540	if (leapf.or.(.not.leapf.and.(.not.forward))) then
1541
1542	c$OMP BARRIER
1543	c$OMP MASTER
1544	CALL geopot_loc(ip1jmp1,teta,pk,pks,phis,phi)
1545	c$OMP END MASTER
1546	c$OMP BARRIER
1547
1548	#ifdef CPP_IOIPSL
1549	if (ok_dyn_ins) then
1550	CALL writehist_loc(itau,vcov,ucov,teta,phi,q,
1551	& masse,ps,phis)
1552	endif
1553	#endif
1554	endif ! of if (leapf.or.(.not.leapf.and.(.not.forward)))
1555	ENDIF ! of IF(MOD(itau,iecri).EQ.0)
1556
1557	IF(itau.EQ.itaufin) THEN
1558
1559	c$OMP BARRIER
1560
1561	! if (planet_type.eq."earth") then
1562	! Write an Earth-format restart file
1563	CALL dynredem1_loc("restart.nc",0.0,
1564	& vcov,ucov,teta,q,masse,ps)
1565	! endif ! of if (planet_type.eq."earth")
1566
1567	! CLOSE(99)
1568	ENDIF ! of IF (itau.EQ.itaufin)
1569
1570	c-----------------------------------------------------------------------
1571	c gestion de l'integration temporelle:
1572	c ------------------------------------
1573
1574	IF( MOD(itau,iperiod).EQ.0 ) THEN
1575	GO TO 1
1576	ELSE IF ( MOD(itau-1,iperiod). EQ. 0 ) THEN
1577
1578	IF( forward ) THEN
1579	c fin du pas forward et debut du pas backward
1580
1581	forward = .FALSE.
1582	leapf = .FALSE.
1583	GO TO 2
1584
1585	ELSE
1586	c fin du pas backward et debut du premier pas leapfrog
1587
1588	leapf = .TRUE.
1589	dt = 2.*dtvr
1590	GO TO 2
1591	END IF
1592	ELSE
1593
1594	c ...... pas leapfrog .....
1595
1596	leapf = .TRUE.
1597	dt = 2.*dtvr
1598	GO TO 2
1599	END IF ! of IF (MOD(itau,iperiod).EQ.0)
1600	! ELSEIF (MOD(itau-1,iperiod).EQ.0)
1601
1602
1603	ELSE ! of IF (.not.purmats)
1604
1605	c ........................................................
1606	c .............. schema matsuno ...............
1607	c ........................................................
1608	IF( forward ) THEN
1609
1610	itau = itau + 1
1611	! iday = day_ini+itau/day_step
1612	! time = REAL(itau-(iday-day_ini)*day_step)/day_step+time_0
1613	!
1614	! IF(time.GT.1.) THEN
1615	! time = time-1.
1616	! iday = iday+1
1617	! ENDIF
1618
1619	forward = .FALSE.
1620	IF( itau. EQ. itaufinp1 ) then
1621	c$OMP MASTER
1622	call fin_getparam
1623	call finalize_parallel
1624	c$OMP END MASTER
1625	abort_message = 'Simulation finished'
1626	call abort_gcm(modname,abort_message,0)
1627	RETURN
1628	ENDIF
1629	GO TO 2
1630
1631	ELSE ! of IF(forward) i.e. backward step
1632
1633	IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) THEN
1634	IF(itau.EQ.itaufin) THEN
1635	iav=1
1636	ELSE
1637	iav=0
1638	ENDIF
1639
1640	#ifdef CPP_IOIPSL
1641	IF (ok_dynzon) THEN
1642	CALL bilan_dyn_loc(2,dtvriperiod,dtvrday_step*periodav,
1643	, ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,q)
1644	ENDIF
1645
1646	IF (ok_dyn_ave) THEN
1647	CALL writedynav_loc(itau,vcov,
1648	& ucov,teta,pk,phi,q,masse,ps,phis)
1649	ENDIF
1650	#endif
1651	ENDIF ! of IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin)
1652
1653
1654	IF(MOD(itau,iecri ).EQ.0) THEN
1655
1656	c$OMP BARRIER
1657	c$OMP MASTER
1658	CALL geopot_loc(ip1jmp1,teta,pk,pks,phis,phi)
1659	c$OMP END MASTER
1660	c$OMP BARRIER
1661
1662
1663	#ifdef CPP_IOIPSL
1664	if (ok_dyn_ins) then
1665	CALL writehist_loc(itau,vcov,ucov,teta,phi,q,
1666	& masse,ps,phis)
1667	endif ! of if (ok_dyn_ins)
1668	#endif
1669	ENDIF ! of IF(MOD(itau,iecri).EQ.0)
1670
1671
1672	IF(itau.EQ.itaufin) THEN
1673	! if (planet_type.eq."earth") then
1674	CALL dynredem1_loc("restart.nc",0.0,
1675	. vcov,ucov,teta,q,masse,ps)
1676	! endif ! of if (planet_type.eq."earth")
1677	ENDIF ! of IF(itau.EQ.itaufin)
1678
1679	forward = .TRUE.
1680	GO TO 1
1681
1682	ENDIF ! of IF (forward)
1683
1684	END IF ! of IF(.not.purmats)
1685	c$OMP MASTER
1686	call fin_getparam
1687	call finalize_parallel
1688	c$OMP END MASTER
1689	abort_message = 'Simulation finished'
1690	call abort_gcm(modname,abort_message,0)
1691	RETURN
1692	END

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

Download in other formats:

Original Format