Context Navigation

leapfrog_loc.F @ 1795

Last change on this file since 1795 was 1795, checked in by Ehouarn Millour, 11 years ago

Version testing basee sur la r1794

Testing release based on r1794

File size: 51.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,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	CALL top_bound_loc(vcov,ucov,teta,masse,dtvr)
1128	!$OMP BARRIER
1129	ENDIF ! of IF (ok_strato)
1130	ENDIF ! of IF(iflag_phys.EQ.2)
1131
1132
1133	CALL pression_loc ( ip1jmp1, ap, bp, ps, p )
1134	c$OMP BARRIER
1135	if (pressure_exner) then
1136	CALL exner_hyb_loc( ijnb_u, ps, p,alpha,beta, pks, pk, pkf )
1137	else
1138	CALL exner_milieu_loc( ijnb_u, ps, p, beta, pks, pk, pkf )
1139	endif
1140	c$OMP BARRIER
1141
1142	cc$OMP END PARALLEL
1143
1144	c-----------------------------------------------------------------------
1145	c dissipation horizontale et verticale des petites echelles:
1146	c ----------------------------------------------------------
1147
1148	IF(apdiss) THEN
1149
1150	CALL call_dissip(ucov,vcov,teta,p,pk,ps)
1151	!cc$OMP PARALLEL DEFAULT(SHARED)
1152	!cc$OMP+ PRIVATE(ijb,ije,tppn,tpn,tpps,tps)
1153	!c$OMP MASTER
1154	! call suspend_timer(timer_caldyn)
1155	!
1156	!c print*,'Entree dans la dissipation : Iteration No ',true_itau
1157	!c calcul de l'energie cinetique avant dissipation
1158	!c print *,'Passage dans la dissipation'
1159
1160	! call VTb(VThallo)
1161	!c$OMP END MASTER
1162
1163	!c$OMP BARRIER
1164
1165	! call Register_SwapField_u(ucov,ucov,distrib_dissip,
1166	! * Request_dissip,up=1,down=1)
1167
1168	! call Register_SwapField_v(vcov,vcov,distrib_dissip,
1169	! * Request_dissip,up=1,down=1)
1170
1171	! call Register_SwapField_u(teta,teta,distrib_dissip,
1172	! * Request_dissip)
1173
1174	! call Register_SwapField_u(p,p,distrib_dissip,
1175	! * Request_dissip)
1176
1177	! call Register_SwapField_u(pk,pk,distrib_dissip,
1178	! * Request_dissip)
1179
1180	! call SendRequest(Request_dissip)
1181	!c$OMP BARRIER
1182	! call WaitRequest(Request_dissip)
1183
1184	!c$OMP BARRIER
1185	!c$OMP MASTER
1186	! call set_distrib(distrib_dissip)
1187	! call VTe(VThallo)
1188	! call VTb(VTdissipation)
1189	! call start_timer(timer_dissip)
1190	!c$OMP END MASTER
1191	!c$OMP BARRIER
1192
1193	! call covcont_loc(llm,ucov,vcov,ucont,vcont)
1194	! call enercin_loc(vcov,ucov,vcont,ucont,ecin0)
1195
1196	!c dissipation
1197
1198	!! CALL FTRACE_REGION_BEGIN("dissip")
1199	! CALL dissip_loc(vcov,ucov,teta,p,dvdis,dudis,dtetadis)
1200
1201	!#ifdef DEBUG_IO
1202	! call WriteField_u('dudis',dudis)
1203	! call WriteField_v('dvdis',dvdis)
1204	! call WriteField_u('dtetadis',dtetadis)
1205	!#endif
1206	!
1207	!! CALL FTRACE_REGION_END("dissip")
1208	!
1209	! ijb=ij_begin
1210	! ije=ij_end
1211	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1212	! DO l=1,llm
1213	! ucov(ijb:ije,l)=ucov(ijb:ije,l)+dudis(ijb:ije,l)
1214	! ENDDO
1215	!c$OMP END DO NOWAIT
1216	! if (pole_sud) ije=ije-iip1
1217	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1218	! DO l=1,llm
1219	! vcov(ijb:ije,l)=vcov(ijb:ije,l)+dvdis(ijb:ije,l)
1220	! ENDDO
1221	!c$OMP END DO NOWAIT
1222
1223	!c teta=teta+dtetadis
1224
1225
1226	!c------------------------------------------------------------------------
1227	! if (dissip_conservative) then
1228	!C On rajoute la tendance due a la transform. Ec -> E therm. cree
1229	!C lors de la dissipation
1230	!c$OMP BARRIER
1231	!c$OMP MASTER
1232	! call suspend_timer(timer_dissip)
1233	! call VTb(VThallo)
1234	!c$OMP END MASTER
1235	! call Register_Hallo_u(ucov,llm,1,1,1,1,Request_Dissip)
1236	! call Register_Hallo_v(vcov,llm,1,1,1,1,Request_Dissip)
1237	! call SendRequest(Request_Dissip)
1238	!c$OMP BARRIER
1239	! call WaitRequest(Request_Dissip)
1240	!c$OMP MASTER
1241	! call VTe(VThallo)
1242	! call resume_timer(timer_dissip)
1243	!c$OMP END MASTER
1244	!c$OMP BARRIER
1245	! call covcont_loc(llm,ucov,vcov,ucont,vcont)
1246	! call enercin_loc(vcov,ucov,vcont,ucont,ecin)
1247	!
1248	! ijb=ij_begin
1249	! ije=ij_end
1250	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1251	! do l=1,llm
1252	! do ij=ijb,ije
1253	! dtetaecdt(ij,l)= (ecin0(ij,l)-ecin(ij,l))/ pk(ij,l)
1254	! dtetadis(ij,l)=dtetadis(ij,l)+dtetaecdt(ij,l)
1255	! enddo
1256	! enddo
1257	!c$OMP END DO NOWAIT
1258	! endif
1259
1260	! ijb=ij_begin
1261	! ije=ij_end
1262	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1263	! do l=1,llm
1264	! do ij=ijb,ije
1265	! teta(ij,l)=teta(ij,l)+dtetadis(ij,l)
1266	! enddo
1267	! enddo
1268	!c$OMP END DO NOWAIT
1269	!c------------------------------------------------------------------------
1270
1271
1272	!c ....... P. Le Van ( ajout le 17/04/96 ) ...........
1273	!c ... Calcul de la valeur moyenne, unique de h aux poles .....
1274	!c
1275
1276	! ijb=ij_begin
1277	! ije=ij_end
1278	!
1279	! if (pole_nord) then
1280	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1281	! DO l = 1, llm
1282	! DO ij = 1,iim
1283	! tppn(ij) = aire( ij ) * teta( ij ,l)
1284	! ENDDO
1285	! tpn = SSUM(iim,tppn,1)/apoln
1286
1287	! DO ij = 1, iip1
1288	! teta( ij ,l) = tpn
1289	! ENDDO
1290	! ENDDO
1291	!c$OMP END DO NOWAIT
1292
1293	!c$OMP MASTER
1294	! DO ij = 1,iim
1295	! tppn(ij) = aire( ij ) * ps ( ij )
1296	! ENDDO
1297	! tpn = SSUM(iim,tppn,1)/apoln
1298	!
1299	! DO ij = 1, iip1
1300	! ps( ij ) = tpn
1301	! ENDDO
1302	!c$OMP END MASTER
1303	! endif
1304	!
1305	! if (pole_sud) then
1306	!c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
1307	! DO l = 1, llm
1308	! DO ij = 1,iim
1309	! tpps(ij) = aire(ij+ip1jm) * teta(ij+ip1jm,l)
1310	! ENDDO
1311	! tps = SSUM(iim,tpps,1)/apols
1312
1313	! DO ij = 1, iip1
1314	! teta(ij+ip1jm,l) = tps
1315	! ENDDO
1316	! ENDDO
1317	!c$OMP END DO NOWAIT
1318
1319	!c$OMP MASTER
1320	! DO ij = 1,iim
1321	! tpps(ij) = aire(ij+ip1jm) * ps (ij+ip1jm)
1322	! ENDDO
1323	! tps = SSUM(iim,tpps,1)/apols
1324	!
1325	! DO ij = 1, iip1
1326	! ps(ij+ip1jm) = tps
1327	! ENDDO
1328	!c$OMP END MASTER
1329	! endif
1330
1331
1332	!c$OMP BARRIER
1333	!c$OMP MASTER
1334	! call VTe(VTdissipation)
1335
1336	! call stop_timer(timer_dissip)
1337	!
1338	! call VTb(VThallo)
1339	!c$OMP END MASTER
1340	! call Register_SwapField_u(ucov,ucov,distrib_caldyn,
1341	! * Request_dissip)
1342
1343	! call Register_SwapField_v(vcov,vcov,distrib_caldyn,
1344	! * Request_dissip)
1345
1346	! call Register_SwapField_u(teta,teta,distrib_caldyn,
1347	! * Request_dissip)
1348
1349	! call Register_SwapField_u(p,p,distrib_caldyn,
1350	! * Request_dissip)
1351
1352	! call Register_SwapField_u(pk,pk,distrib_caldyn,
1353	! * Request_dissip)
1354
1355	! call SendRequest(Request_dissip)
1356	!c$OMP BARRIER
1357	! call WaitRequest(Request_dissip)
1358
1359	!c$OMP BARRIER
1360	!c$OMP MASTER
1361	! call set_distrib(distrib_caldyn)
1362	! call VTe(VThallo)
1363	! call resume_timer(timer_caldyn)
1364	!c print *,'fin dissipation'
1365	!c$OMP END MASTER
1366	!c$OMP BARRIER
1367	END IF ! of IF(apdiss)
1368
1369	cc$OMP END PARALLEL
1370
1371	c ajout debug
1372	c IF( lafin ) then
1373	c abort_message = 'Simulation finished'
1374	c call abort_gcm(modname,abort_message,0)
1375	c ENDIF
1376
1377	c ********************************************************************
1378	c ********************************************************************
1379	c .... fin de l'integration dynamique et physique pour le pas itau ..
1380	c ********************************************************************
1381	c ********************************************************************
1382
1383	c preparation du pas d'integration suivant ......
1384	cym call WriteField('ucov',reshape(ucov,(/iip1,jmp1,llm/)))
1385	cym call WriteField('vcov',reshape(vcov,(/iip1,jjm,llm/)))
1386	c$OMP MASTER
1387	call stop_timer(timer_caldyn)
1388	c$OMP END MASTER
1389	IF (itau==itaumax) then
1390	c$OMP MASTER
1391	call allgather_timer_average
1392
1393	if (mpi_rank==0) then
1394
1395	print ,'********************************'
1396	print ,'*** TIMER CALDYN ****'
1397	do i=0,mpi_size-1
1398	print *,'proc',i,' : Nb Bandes :',jj_nb_caldyn(i),
1399	& ' : temps moyen :',
1400	& timer_average(jj_nb_caldyn(i),timer_caldyn,i)
1401	enddo
1402
1403	print ,'********************************'
1404	print ,'*** TIMER VANLEER ****'
1405	do i=0,mpi_size-1
1406	print *,'proc',i,' : Nb Bandes :',jj_nb_vanleer(i),
1407	& ' : temps moyen :',
1408	& timer_average(jj_nb_vanleer(i),timer_vanleer,i)
1409	enddo
1410
1411	print ,'********************************'
1412	print ,'*** TIMER DISSIP ****'
1413	do i=0,mpi_size-1
1414	print *,'proc',i,' : Nb Bandes :',jj_nb_dissip(i),
1415	& ' : temps moyen :',
1416	& timer_average(jj_nb_dissip(i),timer_dissip,i)
1417	enddo
1418
1419	print ,'********************************'
1420	print ,'*** TIMER PHYSIC ****'
1421	do i=0,mpi_size-1
1422	print *,'proc',i,' : Nb Bandes :',jj_nb_physic(i),
1423	& ' : temps moyen :',
1424	& timer_average(jj_nb_physic(i),timer_physic,i)
1425	enddo
1426
1427	endif
1428
1429	print ,'Taille du Buffer MPI (REAL8)',MaxBufferSize
1430	print ,'Taille du Buffer MPI utilise (REAL8)',MaxBufferSize_Used
1431	print *, 'Temps total ecoule sur la parallelisation :',DiffTime()
1432	print *, 'Temps CPU ecoule sur la parallelisation :',DiffCpuTime()
1433	CALL print_filtre_timer
1434	c$OMP END MASTER
1435	CALL dynredem1_loc("restart.nc",0.0,
1436	. vcov,ucov,teta,q,masse,ps)
1437	c$OMP MASTER
1438	call fin_getparam
1439	call finalize_parallel
1440	c$OMP END MASTER
1441	c$OMP BARRIER
1442	RETURN
1443	ENDIF
1444
1445	IF ( .NOT.purmats ) THEN
1446	c ........................................................
1447	c .............. schema matsuno + leapfrog ..............
1448	c ........................................................
1449
1450	IF(forward. OR. leapf) THEN
1451	itau= itau + 1
1452	! iday= day_ini+itau/day_step
1453	! time= REAL(itau-(iday-day_ini)*day_step)/day_step+time_0
1454	! IF(time.GT.1.) THEN
1455	! time = time-1.
1456	! iday = iday+1
1457	! ENDIF
1458	ENDIF
1459
1460
1461	IF( itau. EQ. itaufinp1 ) then
1462
1463	if (flag_verif) then
1464	write(79,*) 'ucov',ucov
1465	write(80,*) 'vcov',vcov
1466	write(81,*) 'teta',teta
1467	write(82,*) 'ps',ps
1468	write(83,*) 'q',q
1469	WRITE(85,*) 'q1 = ',q(:,:,1)
1470	WRITE(86,*) 'q3 = ',q(:,:,3)
1471	endif
1472
1473
1474	c$OMP MASTER
1475	call fin_getparam
1476	call finalize_parallel
1477	c$OMP END MASTER
1478	abort_message = 'Simulation finished'
1479	call abort_gcm(modname,abort_message,0)
1480	RETURN
1481	ENDIF
1482	c-----------------------------------------------------------------------
1483	c ecriture du fichier histoire moyenne:
1484	c -------------------------------------
1485
1486	IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) THEN
1487	c$OMP BARRIER
1488	IF(itau.EQ.itaufin) THEN
1489	iav=1
1490	ELSE
1491	iav=0
1492	ENDIF
1493
1494	#ifdef CPP_IOIPSL
1495	IF (ok_dynzon) THEN
1496
1497	CALL bilan_dyn_loc(2,dtvriperiod,dtvrday_step*periodav,
1498	, ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,q)
1499
1500	ENDIF !ok_dynzon
1501
1502	IF (ok_dyn_ave) THEN
1503	CALL writedynav_loc(itau,vcov,
1504	& ucov,teta,pk,phi,q,masse,ps,phis)
1505	ENDIF
1506	#endif
1507	ENDIF
1508
1509	c-----------------------------------------------------------------------
1510	c ecriture de la bande histoire:
1511	c ------------------------------
1512
1513	IF( MOD(itau,iecri).EQ.0) THEN
1514	! Ehouarn: output only during LF or Backward Matsuno
1515	if (leapf.or.(.not.leapf.and.(.not.forward))) then
1516
1517	c$OMP BARRIER
1518	c$OMP MASTER
1519	CALL geopot_loc(ip1jmp1,teta,pk,pks,phis,phi)
1520	c$OMP END MASTER
1521	c$OMP BARRIER
1522
1523	#ifdef CPP_IOIPSL
1524	if (ok_dyn_ins) then
1525	CALL writehist_loc(itau,vcov,ucov,teta,phi,q,
1526	& masse,ps,phis)
1527	endif
1528	#endif
1529	endif ! of if (leapf.or.(.not.leapf.and.(.not.forward)))
1530	ENDIF ! of IF(MOD(itau,iecri).EQ.0)
1531
1532	IF(itau.EQ.itaufin) THEN
1533
1534	c$OMP BARRIER
1535
1536	! if (planet_type.eq."earth") then
1537	! Write an Earth-format restart file
1538	CALL dynredem1_loc("restart.nc",0.0,
1539	& vcov,ucov,teta,q,masse,ps)
1540	! endif ! of if (planet_type.eq."earth")
1541
1542	! CLOSE(99)
1543	ENDIF ! of IF (itau.EQ.itaufin)
1544
1545	c-----------------------------------------------------------------------
1546	c gestion de l'integration temporelle:
1547	c ------------------------------------
1548
1549	IF( MOD(itau,iperiod).EQ.0 ) THEN
1550	GO TO 1
1551	ELSE IF ( MOD(itau-1,iperiod). EQ. 0 ) THEN
1552
1553	IF( forward ) THEN
1554	c fin du pas forward et debut du pas backward
1555
1556	forward = .FALSE.
1557	leapf = .FALSE.
1558	GO TO 2
1559
1560	ELSE
1561	c fin du pas backward et debut du premier pas leapfrog
1562
1563	leapf = .TRUE.
1564	dt = 2.*dtvr
1565	GO TO 2
1566	END IF
1567	ELSE
1568
1569	c ...... pas leapfrog .....
1570
1571	leapf = .TRUE.
1572	dt = 2.*dtvr
1573	GO TO 2
1574	END IF ! of IF (MOD(itau,iperiod).EQ.0)
1575	! ELSEIF (MOD(itau-1,iperiod).EQ.0)
1576
1577
1578	ELSE ! of IF (.not.purmats)
1579
1580	c ........................................................
1581	c .............. schema matsuno ...............
1582	c ........................................................
1583	IF( forward ) THEN
1584
1585	itau = itau + 1
1586	! iday = day_ini+itau/day_step
1587	! time = REAL(itau-(iday-day_ini)*day_step)/day_step+time_0
1588	!
1589	! IF(time.GT.1.) THEN
1590	! time = time-1.
1591	! iday = iday+1
1592	! ENDIF
1593
1594	forward = .FALSE.
1595	IF( itau. EQ. itaufinp1 ) then
1596	c$OMP MASTER
1597	call fin_getparam
1598	call finalize_parallel
1599	c$OMP END MASTER
1600	abort_message = 'Simulation finished'
1601	call abort_gcm(modname,abort_message,0)
1602	RETURN
1603	ENDIF
1604	GO TO 2
1605
1606	ELSE ! of IF(forward) i.e. backward step
1607
1608	IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin) THEN
1609	IF(itau.EQ.itaufin) THEN
1610	iav=1
1611	ELSE
1612	iav=0
1613	ENDIF
1614
1615	#ifdef CPP_IOIPSL
1616	IF (ok_dynzon) THEN
1617	CALL bilan_dyn_loc(2,dtvriperiod,dtvrday_step*periodav,
1618	, ps,masse,pk,pbaru,pbarv,teta,phi,ucov,vcov,q)
1619	ENDIF
1620
1621	IF (ok_dyn_ave) THEN
1622	CALL writedynav_loc(itau,vcov,
1623	& ucov,teta,pk,phi,q,masse,ps,phis)
1624	ENDIF
1625	#endif
1626	ENDIF ! of IF(MOD(itau,iperiod).EQ.0 .OR. itau.EQ.itaufin)
1627
1628
1629	IF(MOD(itau,iecri ).EQ.0) THEN
1630
1631	c$OMP BARRIER
1632	c$OMP MASTER
1633	CALL geopot_loc(ip1jmp1,teta,pk,pks,phis,phi)
1634	c$OMP END MASTER
1635	c$OMP BARRIER
1636
1637
1638	#ifdef CPP_IOIPSL
1639	if (ok_dyn_ins) then
1640	CALL writehist_loc(itau,vcov,ucov,teta,phi,q,
1641	& masse,ps,phis)
1642	endif ! of if (ok_dyn_ins)
1643	#endif
1644	ENDIF ! of IF(MOD(itau,iecri).EQ.0)
1645
1646
1647	IF(itau.EQ.itaufin) THEN
1648	! if (planet_type.eq."earth") then
1649	CALL dynredem1_loc("restart.nc",0.0,
1650	. vcov,ucov,teta,q,masse,ps)
1651	! endif ! of if (planet_type.eq."earth")
1652	ENDIF ! of IF(itau.EQ.itaufin)
1653
1654	forward = .TRUE.
1655	GO TO 1
1656
1657	ENDIF ! of IF (forward)
1658
1659	END IF ! of IF(.not.purmats)
1660	c$OMP MASTER
1661	call fin_getparam
1662	call finalize_parallel
1663	c$OMP END MASTER
1664	abort_message = 'Simulation finished'
1665	call abort_gcm(modname,abort_message,0)
1666	RETURN
1667	END

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

Download in other formats:

Original Format