source: trunk/LMDZ.COMMON/libf/dyn3d/nogcm_pluto.F90 @ 3553

Last change on this file since 3553 was 3229, checked in by tbertrand, 10 months ago

Pluto PCM:
Adding Nogcm.F : 2D surface model with 1 atm layer

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