Context Navigation

leapfrog_p.F @ 4359

Last change on this file since 4359 was 935, checked in by Laurent Fairhead, 17 years ago
Rajout ok_dynzon pour faire ou non les sorties bilan LF
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 44.9 KB

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

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

Download in other formats:

Original Format