source: LMDZ5/trunk/libf/dyn3d/gcm.F90 @ 2343

Last change on this file since 2343 was 2247, checked in by lguez, 10 years ago

We want to keep the same "*.def" files with programs ce0l and gcm. But
there is only a sequential version of program ce0l, which has not the
FFT filter. So we have to ignore the setting of use_filtre_fft in
program ce0l. Moved the test on use_filtre_fft from procedure conf_gcm
to main units ce0l and gcm.

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