Context Navigation

leapfrog_loc.F @ 5452

Last change on this file since 5452 was 2381, checked in by acozic, 9 years ago

Make some commit to fit with INCA coupling

merge with rev 2180 on trunk
merge with rev 2185 on trunk
merge with rev 2200 on trunk
change gregorian calendar in wxios

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

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

Download in other formats:

Original Format