source: LMDZ5/trunk/libf/dyn3dmem/gcm.F90 @ 2402

Last change on this file since 2402 was 2372, checked in by acozic, 9 years ago

Change call to inca initialisation to fit with new sections dynamique/physic

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
File size: 14.1 KB
RevLine 
[2225]1! $Id: $
[1632]2
[2263]3PROGRAM gcm
4
[1632]5#ifdef CPP_IOIPSL
[2263]6  USE IOIPSL
[1632]7#endif
8
[2263]9  USE mod_const_mpi, ONLY: init_const_mpi
10  USE parallel_lmdz
11  USE infotrac
[2351]12!#ifdef CPP_PHYS
13!  USE mod_interface_dyn_phys, ONLY: init_interface_dyn_phys
14!#endif
[2263]15  USE mod_hallo
16  USE Bands
17  USE filtreg_mod
18  USE control_mod
[1632]19
[1673]20#ifdef CPP_PHYS
[2347]21  USE iniphysiq_mod, ONLY: iniphysiq
[1632]22#endif
[2263]23  IMPLICIT NONE
[1632]24
[2263]25  !      ......   Version  du 10/01/98    ..........
[1632]26
[2263]27  !             avec  coordonnees  verticales hybrides
28  !   avec nouveaux operat. dissipation * ( gradiv2,divgrad2,nxgraro2 )
[1632]29
[2263]30  !=======================================================================
31  !
32  !   Auteur:  P. Le Van /L. Fairhead/F.Hourdin
33  !   -------
34  !
35  !   Objet:
36  !   ------
37  !
38  !   GCM LMD nouvelle grille
39  !
40  !=======================================================================
41  !
42  !  ... Dans inigeom , nouveaux calculs pour les elongations  cu , cv
43  !      et possibilite d'appeler une fonction f(y)  a derivee tangente
44  !      hyperbolique a la  place de la fonction a derivee sinusoidale.
45  !  ... Possibilite de choisir le schema pour l'advection de
46  !        q  , en modifiant iadv dans traceur.def  (MAF,10/02) .
47  !
48  !      Pour Van-Leer + Vapeur d'eau saturee, iadv(1)=4. (F.Codron,10/99)
49  !      Pour Van-Leer iadv=10
50  !
51  !-----------------------------------------------------------------------
52  !   Declarations:
53  !   -------------
54  include "dimensions.h"
55  include "paramet.h"
56  include "comconst.h"
57  include "comdissnew.h"
58  include "comvert.h"
59  include "comgeom.h"
60  include "logic.h"
61  include "temps.h"
62  include "ener.h"
63  include "description.h"
64  include "serre.h"
65  !include "com_io_dyn.h"
66  include "iniprint.h"
67  include "tracstoke.h"
[1658]68
[1632]69
[2263]70  REAL zdtvr
[1632]71
[2263]72  !   variables dynamiques
73  REAL,ALLOCATABLE,SAVE  :: vcov(:,:),ucov(:,:) ! vents covariants
74  REAL,ALLOCATABLE,SAVE  :: teta(:,:)     ! temperature potentielle
75  REAL, ALLOCATABLE,SAVE :: q(:,:,:)      ! champs advectes
76  REAL,ALLOCATABLE,SAVE  :: ps(:)         ! pression  au sol
77  !      REAL p (ip1jmp1,llmp1  )               ! pression aux interfac.des couches
78  REAL,ALLOCATABLE,SAVE  :: masse(:,:)    ! masse d'air
79  REAL,ALLOCATABLE,SAVE  :: phis(:)       ! geopotentiel au sol
80  !      REAL phi(ip1jmp1,llm)                  ! geopotentiel
81  !      REAL w(ip1jmp1,llm)                    ! vitesse verticale
[1632]82
[2263]83  ! variables dynamiques intermediaire pour le transport
[1632]84
[2263]85  !   variables pour le fichier histoire
86  REAL dtav      ! intervalle de temps elementaire
[1632]87
[2263]88  REAL time_0
[1632]89
[2263]90  LOGICAL lafin
[1632]91
[2263]92  real time_step, t_wrt, t_ops
[1632]93
[2263]94  !+jld variables test conservation energie
95  !      REAL ecin(ip1jmp1,llm),ecin0(ip1jmp1,llm)
96  !     Tendance de la temp. potentiel d (theta)/ d t due a la
97  !     tansformation d'energie cinetique en energie thermique
98  !     cree par la dissipation
99  !      REAL dhecdt(ip1jmp1,llm)
100  !      REAL vcont(ip1jm,llm),ucont(ip1jmp1,llm)
101  !      REAL      d_h_vcol, d_qt, d_qw, d_ql, d_ec
102  !      CHARACTER (len=15) :: ztit
103  !-jld
[1632]104
105
[2263]106  character (len=80) :: dynhist_file, dynhistave_file
107  character (len=20) :: modname
108  character (len=80) :: abort_message
109  ! locales pour gestion du temps
110  INTEGER :: an, mois, jour
111  REAL :: heure
[1632]112
113
[2263]114  !-----------------------------------------------------------------------
115  !   Initialisations:
116  !   ----------------
[1632]117
[2263]118  abort_message = 'last timestep reached'
119  modname = 'gcm'
120  descript = 'Run GCM LMDZ'
121  lafin    = .FALSE.
122  dynhist_file = 'dyn_hist'
123  dynhistave_file = 'dyn_hist_ave'
[1632]124
125
126
[2263]127  !----------------------------------------------------------------------
128  !  lecture des fichiers gcm.def ou run.def
129  !  ---------------------------------------
130  !
131  CALL conf_gcm( 99, .TRUE. )
132  if (mod(iphysiq, iperiod) /= 0) call abort_gcm("conf_gcm", &
133       "iphysiq must be a multiple of iperiod", 1)
134  !
135  !
136  !------------------------------------
137  !   Initialisation partie parallele
138  !------------------------------------
139  CALL init_const_mpi
[1632]140
[2263]141  call init_parallel
142  call Read_Distrib
[1673]143
[2351]144!#ifdef CPP_PHYS
145!  CALL Init_Phys_lmdz(iim,jjp1,llm,mpi_size,distrib_phys)
[2263]146  !#endif
147  !      CALL set_bands
148  !#ifdef CPP_PHYS
[2351]149!  CALL Init_interface_dyn_phys
150!#endif
[2263]151  CALL barrier
[1632]152
[2263]153  CALL set_bands
154  if (mpi_rank==0) call WriteBands
155  call Set_Distrib(distrib_caldyn)
[1632]156
[2263]157  !$OMP PARALLEL
158  call Init_Mod_hallo
159  !$OMP END PARALLEL
[1632]160
[2263]161  !#ifdef CPP_PHYS
162  !c$OMP PARALLEL
163  !      call InitComgeomphy ! now done in iniphysiq
164  !c$OMP END PARALLEL
165  !#endif
[1632]166
[2263]167  !-----------------------------------------------------------------------
168  !   Choix du calendrier
169  !   -------------------
[1632]170
[2263]171  !      calend = 'earth_365d'
[1632]172
173#ifdef CPP_IOIPSL
[2263]174  if (calend == 'earth_360d') then
175     call ioconf_calendar('360d')
176     write(lunout,*)'CALENDRIER CHOISI: Terrestre a 360 jours/an'
177  else if (calend == 'earth_365d') then
178     call ioconf_calendar('noleap')
179     write(lunout,*)'CALENDRIER CHOISI: Terrestre a 365 jours/an'
180  else if (calend == 'gregorian') then
181     call ioconf_calendar('gregorian')
182     write(lunout,*)'CALENDRIER CHOISI: Terrestre bissextile'
183  else
184     abort_message = 'Mauvais choix de calendrier'
185     call abort_gcm(modname,abort_message,1)
186  endif
[1632]187#endif
188
189
[2263]190  !-----------------------------------------------------------------------
191  !   Initialisation des traceurs
192  !   ---------------------------
193  !  Choix du nombre de traceurs et du schema pour l'advection
194  !  dans fichier traceur.def, par default ou via INCA
195  call infotrac_init
[1632]196
[2263]197  ! Allocation de la tableau q : champs advectes   
198  ALLOCATE(ucov(ijb_u:ije_u,llm))
199  ALLOCATE(vcov(ijb_v:ije_v,llm))
200  ALLOCATE(teta(ijb_u:ije_u,llm))
201  ALLOCATE(masse(ijb_u:ije_u,llm))
202  ALLOCATE(ps(ijb_u:ije_u))
203  ALLOCATE(phis(ijb_u:ije_u))
204  ALLOCATE(q(ijb_u:ije_u,llm,nqtot))
[1632]205
[2263]206  !-----------------------------------------------------------------------
207  !   Lecture de l'etat initial :
208  !   ---------------------------
[1632]209
[2263]210  !  lecture du fichier start.nc
211  if (read_start) then
212     ! we still need to run iniacademic to initialize some
213     ! constants & fields, if we run the 'newtonian' or 'SW' cases:
214     if (iflag_phys.ne.1) then
215        CALL iniacademic_loc(vcov,ucov,teta,q,masse,ps,phis,time_0)
216     endif
[1657]217
[2263]218     !        if (planet_type.eq."earth") then
219     ! Load an Earth-format start file
220     CALL dynetat0_loc("start.nc",vcov,ucov, &
221          teta,q,masse,ps,phis, time_0)
222     !        endif ! of if (planet_type.eq."earth")
[1673]223
[2263]224     !       write(73,*) 'ucov',ucov
225     !       write(74,*) 'vcov',vcov
226     !       write(75,*) 'teta',teta
227     !       write(76,*) 'ps',ps
228     !       write(77,*) 'q',q
[1632]229
[2263]230  endif ! of if (read_start)
[1632]231
[2263]232  ! le cas echeant, creation d un etat initial
233  IF (prt_level > 9) WRITE(lunout,*) &
234       'GCM: AVANT iniacademic AVANT AVANT AVANT AVANT'
235  if (.not.read_start) then
236     CALL iniacademic_loc(vcov,ucov,teta,q,masse,ps,phis,time_0)
237  endif
[1632]238
[2263]239  !-----------------------------------------------------------------------
240  !   Lecture des parametres de controle pour la simulation :
241  !   -------------------------------------------------------
242  !  on recalcule eventuellement le pas de temps
[1632]243
[2263]244  IF(MOD(day_step,iperiod).NE.0) THEN
245     abort_message =  &
246          'Il faut choisir un nb de pas par jour multiple de iperiod'
247     call abort_gcm(modname,abort_message,1)
248  ENDIF
[1632]249
[2263]250  IF(MOD(day_step,iphysiq).NE.0) THEN
251     abort_message =  &
252          'Il faut choisir un nb de pas par jour multiple de iphysiq'
253     call abort_gcm(modname,abort_message,1)
254  ENDIF
[1632]255
[2263]256  zdtvr    = daysec/REAL(day_step)
257  IF(dtvr.NE.zdtvr) THEN
258     WRITE(lunout,*) &
259          'WARNING!!! changement de pas de temps',dtvr,'>',zdtvr
260  ENDIF
[1632]261
[2263]262  !
263  ! on remet le calendrier \`a zero si demande
264  !
265  IF (start_time /= starttime) then
266     WRITE(lunout,*)' GCM: Attention l''heure de depart lue dans le' &
267          ,' fichier restart ne correspond pas a celle lue dans le run.def'
268     IF (raz_date == 1) then
269        WRITE(lunout,*)'Je prends l''heure lue dans run.def'
270        start_time = starttime
271     ELSE
272        WRITE(lunout,*)'Je m''arrete'
273        CALL abort
274     ENDIF
275  ENDIF
276  IF (raz_date == 1) THEN
277     annee_ref = anneeref
278     day_ref = dayref
279     day_ini = dayref
280     itau_dyn = 0
281     itau_phy = 0
282     time_0 = 0.
283     write(lunout,*) &
284          'GCM: On reinitialise a la date lue dans gcm.def'
285  ELSE IF (annee_ref .ne. anneeref .or. day_ref .ne. dayref) THEN
286     write(lunout,*) &
287          'GCM: Attention les dates initiales lues dans le fichier'
288     write(lunout,*) &
289          ' restart ne correspondent pas a celles lues dans '
290     write(lunout,*)' gcm.def'
291     write(lunout,*)' annee_ref=',annee_ref," anneeref=",anneeref
292     write(lunout,*)' day_ref=',day_ref," dayref=",dayref
293     write(lunout,*)' Pas de remise a zero'
294  ENDIF
295  !      if (annee_ref .ne. anneeref .or. day_ref .ne. dayref) then
296  !        write(lunout,*)
297  !     .  'GCM: Attention les dates initiales lues dans le fichier'
298  !        write(lunout,*)
299  !     .  ' restart ne correspondent pas a celles lues dans '
300  !        write(lunout,*)' gcm.def'
301  !        write(lunout,*)' annee_ref=',annee_ref," anneeref=",anneeref
302  !        write(lunout,*)' day_ref=',day_ref," dayref=",dayref
303  !        if (raz_date .ne. 1) then
304  !          write(lunout,*)
305  !     .    'GCM: On garde les dates du fichier restart'
306  !        else
307  !          annee_ref = anneeref
308  !          day_ref = dayref
309  !          day_ini = dayref
310  !          itau_dyn = 0
311  !          itau_phy = 0
312  !          time_0 = 0.
313  !          write(lunout,*)
314  !     .   'GCM: On reinitialise a la date lue dans gcm.def'
315  !        endif
316  !      ELSE
317  !        raz_date = 0
318  !      endif
[1632]319
320#ifdef CPP_IOIPSL
[2263]321  mois = 1
322  heure = 0.
323  call ymds2ju(annee_ref, mois, day_ref, heure, jD_ref)
324  jH_ref = jD_ref - int(jD_ref)
325  jD_ref = int(jD_ref)
[1632]326
[2263]327  call ioconf_startdate(INT(jD_ref), jH_ref)
[1632]328
[2263]329  write(lunout,*)'DEBUG'
330  write(lunout,*)'annee_ref, mois, day_ref, heure, jD_ref'
331  write(lunout,*)annee_ref, mois, day_ref, heure, jD_ref
332  call ju2ymds(jD_ref+jH_ref,an, mois, jour, heure)
333  write(lunout,*)'jD_ref+jH_ref,an, mois, jour, heure'
334  write(lunout,*)jD_ref+jH_ref,an, mois, jour, heure
[1632]335#else
[2263]336  ! Ehouarn: we still need to define JD_ref and JH_ref
337  ! and since we don't know how many days there are in a year
338  ! we set JD_ref to 0 (this should be improved ...)
339  jD_ref=0
340  jH_ref=0
[1632]341#endif
342
[2263]343  if (iflag_phys.eq.1) then
344     ! these initialisations have already been done (via iniacademic)
345     ! if running in SW or Newtonian mode
346     !-----------------------------------------------------------------------
347     !   Initialisation des constantes dynamiques :
348     !   ------------------------------------------
349     dtvr = zdtvr
350     CALL iniconst
[1632]351
[2263]352     !-----------------------------------------------------------------------
353     !   Initialisation de la geometrie :
354     !   --------------------------------
355     CALL inigeom
[1632]356
[2263]357     !-----------------------------------------------------------------------
358     !   Initialisation du filtre :
359     !   --------------------------
360     CALL inifilr
361  endif ! of if (iflag_phys.eq.1)
362  !
363  !-----------------------------------------------------------------------
364  !   Initialisation de la dissipation :
365  !   ----------------------------------
[1632]366
[2263]367  CALL inidissip( lstardis, nitergdiv, nitergrot, niterh   , &
368       tetagdiv, tetagrot , tetatemp, vert_prof_dissip)
[1632]369
[2263]370  !-----------------------------------------------------------------------
371  !   Initialisation de la physique :
372  !   -------------------------------
373  IF ((iflag_phys==1).or.(iflag_phys>=100)) THEN
374     ! Physics:
[1673]375#ifdef CPP_PHYS
[2351]376     CALL iniphysiq(iim,jjm,llm, &
377          distrib_phys(mpi_rank),comm_lmdz, &
378          daysec,day_ini,dtphys/nsplit_phys, &
[2346]379          rlatu,rlatv,rlonu,rlonv,aire,cu,cv,rad,g,r,cpp, &
[2263]380          iflag_phys)
[1632]381#endif
[2263]382  ENDIF ! of IF ((iflag_phys==1).or.(iflag_phys>=100))
[1632]383
384
[2263]385  !-----------------------------------------------------------------------
386  !   Initialisation des I/O :
387  !   ------------------------
[1632]388
389
[2263]390  if (nday>=0) then
391     day_end = day_ini + nday
392  else
393     day_end = day_ini - nday/day_step
394  endif
[1632]395
[2263]396  WRITE(lunout,300)day_ini,day_end
397300 FORMAT('1'/,15x,'run du jour',i7,2x,'au jour',i7//)
398
[1632]399#ifdef CPP_IOIPSL
[2263]400  call ju2ymds(jD_ref + day_ini - day_ref, an, mois, jour, heure)
401  write (lunout,301)jour, mois, an
402  call ju2ymds(jD_ref + day_end - day_ref, an, mois, jour, heure)
403  write (lunout,302)jour, mois, an
404301 FORMAT('1'/,15x,'run du ', i2,'/',i2,'/',i4)
405302 FORMAT('1'/,15x,'    au ', i2,'/',i2,'/',i4)
[1632]406#endif
407
[2263]408  !      if (planet_type.eq."earth") then
409  ! Write an Earth-format restart file
410  CALL dynredem0_loc("restart.nc", day_end, phis)
411  !      endif
[1632]412
[2263]413  ecripar = .TRUE.
[1632]414
415#ifdef CPP_IOIPSL
[2263]416  time_step = zdtvr
417  IF (mpi_rank==0) then
418     if (ok_dyn_ins) then
419        ! initialize output file for instantaneous outputs
420        ! t_ops = iecri * daysec ! do operations every t_ops
421        t_ops =((1.0*iecri)/day_step) * daysec 
422        t_wrt = daysec ! iecri * daysec ! write output every t_wrt
423        t_wrt = daysec ! iecri * daysec ! write output every t_wrt
424        CALL inithist(day_ref,annee_ref,time_step, &
425             t_ops,t_wrt)
426     endif
[1632]427
[2263]428     IF (ok_dyn_ave) THEN
429        ! initialize output file for averaged outputs
430        t_ops = iperiod * time_step ! do operations every t_ops
431        t_wrt = periodav * daysec   ! write output every t_wrt
432        CALL initdynav_loc(day_ref,annee_ref,time_step,t_ops,t_wrt)
433     END IF
434  ENDIF
435  dtav = iperiod*dtvr/daysec
[1632]436#endif
[2263]437  ! #endif of #ifdef CPP_IOIPSL
[1632]438
[2263]439  !  Choix des frequences de stokage pour le offline
440  !      istdyn=day_step/4     ! stockage toutes les 6h=1jour/4
441  !      istdyn=day_step/12     ! stockage toutes les 2h=1jour/12
442  istdyn=day_step/4     ! stockage toutes les 6h=1jour/12
443  istphy=istdyn/iphysiq     
[1632]444
445
[2263]446  !
447  !-----------------------------------------------------------------------
448  !   Integration temporelle du modele :
449  !   ----------------------------------
[1632]450
[2263]451  !       write(78,*) 'ucov',ucov
452  !       write(78,*) 'vcov',vcov
453  !       write(78,*) 'teta',teta
454  !       write(78,*) 'ps',ps
455  !       write(78,*) 'q',q
[1632]456
[2263]457  !$OMP PARALLEL DEFAULT(SHARED) COPYIN(/temps/,/logici/,/logicl/)
458  CALL leapfrog_loc(ucov,vcov,teta,ps,masse,phis,q,time_0)
459  !$OMP END PARALLEL
[1632]460
[2263]461  !      OPEN(unit=5487,file='ok_lmdz',status='replace')
462  !      WRITE(5487,*) 'ok_lmdz'
463  !      CLOSE(5487)
464END PROGRAM gcm
Note: See TracBrowser for help on using the repository browser.