Context Navigation

leapfrog_loc.F @ 2270

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

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

Download in other formats:

Original Format