source: LMDZ6/trunk/libf/dyn3dmem/leapfrog_loc.F @ 4379

Last change on this file since 4379 was 4187, checked in by acozic, 2 years ago

Add a syntax correction in leapfrog_loc for compilation with inca / add a condition on inca in phys_output_write to be use without inca

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