Context Navigation

leapfrog_p.F @ 91

Last change on this file since 91 was 52, checked in by aslmd, 14 years ago

chantier principal du commit
--- version LMDZ5 qui fonctionne pour tests geantes
--- prochaine etape, tests sur GNOME

M libf/dyn3dpar/comconst.h
M libf/dyn3dpar/conf_planete.F90
ajout du flux de chaleur intrinseque: ihf
[par defaut il est nul]

M libf/dyn3dpar/gcm.F
changements cosmetiques
[pour diff plus efficace avec version non par]

M libf/dyn3dpar/iniacademic.F
possibilites de variations latitudinales
de temperature plus originales
[seulement pour planet_type.eq."giant"]

M libf/dyn3dpar/leapfrog_p.F

ajout d'une tendance causee par le flux de chaleur intrinseque

(seulement prise en compte si planet_type.eq."giant")

correction bugs problematiques a la compilation et au run

--> probleme dans les boucles (l'indice etait llm et non l)
--> ajout de SAVE pour les variables paralleles
--> correction des declarations de variables manquantes

M libf/dyn3dpar/calfis_p.F
correction d'une deuxieme parenthese manquante sur ALLOCATE(zteta(klon,llm))

M libf/phylmd/regr_lat_time_climoz_m.F90
erreur a la compilation avec FCM... il s'agit d'une routine terrestre
il y a visiblement un probleme avec o3_in
en attendant, les lignes sont commentees avec !AS

A deftanks/giant 8 fichiers
ajout de fichiers de configuration typiques pour les geantes gazeuses
[experimental pour le moment... on est loin de jupiter]

--> comparaisons entre un run ancien [avec LMDZ5-dev sur SVN ipsl sans cp var]
et run avec version sur ce SVN planeto donne des resultats similaires

pratique

A ioipsl
A ioipsl/compile_ioipsl.bash
A ioipsl/util 16 fichiers
script et utilitaire pour compiler IOIPSL de facon independante
il suffit d'executer ./compile_ioipsl.bash

M arch/arch-AMD64_CICLAD.path
si IOIPSL a ete compile avec la methode precedente, les bons
PATH sont definis dans ce fichier [le NETCDF est aussi OK]

M 000-README-svn
mise a jour options "svn status"

M mars/libf/phymars/meso_callkeys.h
mise a jour mineure du fichier
[ecri_phys etait defini mais pas dans la liste]

File size: 58.2 KB

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

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

Download in other formats:

Original Format