trunk/dyn3d/gcm.f

PROGRAM gcm

  ! General circulation model of LMD

  ! From "gcm.F", version 1.4, 2006/04/04 15:05:16

  ! Avec coordonnées verticales hybrides, avec nouveaux opérateurs de
  ! dissipation * (gradiv2, divgrad2, nxgraro2)

  ! Authors: P. Le Van, L. Fairhead, F. Hourdin

  ! Possibilité de choisir le schéma pour l'advection de "q", en
  ! modifiant "iadv" dans "traceur.def".

  ! Pour Van-Leer plus vapeur d'eau saturée : iadv(1)=4
  ! Pour Van-Leer : iadv=10

  USE IOIPSL, only: ioconf_calendar, histclo
  use dimens_m, only: iim, jjm, llm, nqmx
  use dimphy, only: klon
  use paramet_m, only: ip1jm, ip1jmp1
  use comconst, only: daysec, cpp, dtvr, g, rad, r, initialize
  use comdissnew, only: lstardis, nitergdiv, nitergrot, niterh, tetagdiv, &
       tetagrot, tetatemp
  use conf_gcm_m, only: day_step, iperiod, anneeref, dayref, iecri, iphysiq, &
       nday, raz_date, periodav, conf_gcm
  use logic, only: iflag_phys
  use comgeom, only: rlatu, aire_2d, cu_2d, cv_2d, rlonv
  use temps, only: day_ref, annee_ref, day_end, itau_dyn
  use com_io_dyn, only: histid, histvid, histaveid
  use tracstoke, only: istdyn, istphy
  use abort_gcm_m, only: abort_gcm
  use inithist_m, only: inithist
  use initdynav_m, only: initdynav
  use dynetat0_m, only: dynetat0, day_ini
  use grid_change, only: dyn_phy, init_dyn_phy
  use iniadvtrac_m, only: iniadvtrac
  use leapfrog_m, only: leapfrog
  use dynredem0_m, only: dynredem0
  use clesphys2, only: read_clesphys2
  use inigeom_m, only: inigeom

  IMPLICIT NONE

  REAL zdtvr ! time step for dynamics, in s

  ! Variables dynamiques :
  REAL vcov(ip1jm, llm), ucov(ip1jmp1, llm) ! vents covariants
  REAL teta(ip1jmp1, llm) ! température potentielle 
  REAL q(ip1jmp1, llm, nqmx) ! champs advectés
  REAL ps(ip1jmp1) ! pression au sol (Pa)

  REAL masse(ip1jmp1, llm) ! masse d'air
  REAL phis(iim + 1, jjm + 1) ! géopotentiel au sol

  ! Variables pour le fichier histoire :
  REAL time_0 ! time in day, as a fraction of day, in [0, 1[

  ! Calendrier :
  LOGICAL:: true_calendar = .false. ! default value

  ! Variables pour l'initialisation de la physique :
  REAL zcufi(klon), zcvfi(klon) ! "cu" and "cv" values on the scalar grid
  REAL latfi(klon), lonfi(klon)
  REAL airefi(klon)

  logical mask_v(iim + 1, jjm) 
  ! (mask for points in the "v" grid, first index is for longitude,
  ! second index is for latitude)

  namelist /main_nml/true_calendar

  !------------------------------------------------------------

  print *, "Enter namelist 'main_nml'."
  read (unit=*, nml=main_nml)
  write(unit=*, nml=main_nml)

  ! Initialisations:
  call initialize

  ! Choix du calendrier :
  if (true_calendar) then
     call ioconf_calendar('gregorian')
  else
     call ioconf_calendar('360d')
  endif

  ! Lecture des fichiers "gcm.def" ou "run.def" :
  call read_clesphys2
  CALL conf_gcm

  ! Initialisation des traceurs
  ! Choix du schéma pour l'advection dans le fichier "traceur.def" ou via INCA
  call iniadvtrac

  ! Lecture du fichier "start.nc" :
  CALL dynetat0(vcov, ucov, teta, q, masse, ps, phis, time_0)

  ! Lecture des paramètres de contrôle pour la simulation :
  ! on recalcule éventuellement le pas de temps
  IF (MOD(day_step, iperiod) /= 0) THEN
     call abort_gcm(modname = "gcm", message = &
          'Il faut choisir un nombre de pas par jour multiple de "iperiod".', &
          ierr = 1)
  ENDIF

  IF (MOD(day_step,iphysiq)/=0) THEN
     call abort_gcm(modname = "gcm", message = &
          'Il faut choisir un nombre de pas par jour multiple de "iphysiq".', &
          ierr = 1)
  ENDIF

  ! On remet le calendrier à zero si demandé:
  if (annee_ref /= anneeref .or. day_ref /= dayref) then
     print *, 'Attention : les dates initiales lues dans le fichier ' // &
          '"start" ne correspondent pas à celles lues dans "gcm.def".'
     if (raz_date) then
        print *, 'On réinitialise à la date lue dans "gcm.def".'
        annee_ref = anneeref
        day_ref = dayref
        day_ini = dayref
        itau_dyn = 0
        time_0 = 0.
     else
        print *, 'On garde les dates du fichier "start".'
     endif
  ELSE
     raz_date = .false.
  endif

  ! Initialisation des constantes dynamiques :
  zdtvr = daysec / REAL(day_step)
  IF (dtvr /= zdtvr) THEN
     print *, 'Warning: the time steps in the ".def" file and in ' // &
          '"start.nc" are different'
     print *, 'dtvr (from "start.nc") = ', dtvr
     print *, 'zdtvr (from ".def") = ', zdtvr
     print *, 'Using the value from the ".def" file.'
     dtvr = zdtvr
  ENDIF
  CALL iniconst

  ! Initialisation de la géometrie :
  CALL inigeom

  ! Initialisation du filtre :
  CALL inifilr

  ! Initialisation de la dissipation :
  CALL inidissip(lstardis, nitergdiv, nitergrot, niterh, tetagdiv, tetagrot, &
       tetatemp)

  call init_dyn_phy

  ! Initialisation de la physique :
  IF (iflag_phys == 1) THEN
     latfi(1)=rlatu(1)
     latfi(2:klon-1) = pack(spread(rlatu(2:jjm), 1, iim), .true.)
     latfi(klon)= rlatu(jjm + 1)

     lonfi(1)=0.
     lonfi(2:klon-1) = pack(spread(rlonv(:iim), 2, jjm - 1), .true.)
     lonfi(klon)= 0.

     zcufi = pack(cu_2d, dyn_phy)

     ! Construct a mask for points in the "v" grid:
     mask_v = .true.
     mask_v(2:, 1) = .false.
     mask_v(iim + 1, 2:) = .false.

     zcvfi(:klon - 1) = pack(cv_2d, mask_v)
     zcvfi(klon) = cv_2d(1, jjm)
     ! (that value of "cv_2d" is used twice in "zcvfi")

     airefi = pack(aire_2d, dyn_phy)
     print *, 'Attention : vitesse verticale nulle dans la physique.'
     CALL iniphysiq(klon, llm, latfi, lonfi, airefi, &
          zcufi, zcvfi, rad, g, r, cpp)
  ENDIF

  ! Initialisation des entrées-sorties :
  day_end = day_ini + nday
  print *, "day_ini = ", day_ini
  print *, "day_end = ", day_end

  CALL dynredem0("restart.nc", day_end, phis)
  CALL inithist(day_ref, annee_ref, zdtvr, nqmx, histid, histvid, &
       infile="dyn_hist.nc", t_ops = iecri * daysec, t_wrt = iecri * daysec)
  CALL initdynav(day_ref, annee_ref, zdtvr, nqmx, histaveid, &
       infile='dyn_hist_ave.nc', t_ops = iperiod * zdtvr, &
       t_wrt = periodav * daysec)

  ! Choix des fréquences de stockage pour le hors-ligne :
  istdyn = day_step / 4     ! stockage toutes les 6 h = 1 jour / 4
  istphy = istdyn / iphysiq     

  ! Intégration temporelle du modèle :
  CALL leapfrog(ucov, vcov, teta, ps, masse, phis, q, time_0)

  call histclo
  print *, 'Simulation finished'
  print *, 'Everything is cool'

END PROGRAM gcm
1	guez	3	PROGRAM gcm
2
3			! General circulation model of LMD
4
5			! From "gcm.F", version 1.4, 2006/04/04 15:05:16
6
7			! Avec coordonnées verticales hybrides, avec nouveaux opérateurs de
8			! dissipation * (gradiv2, divgrad2, nxgraro2)
9
10			! Authors: P. Le Van, L. Fairhead, F. Hourdin
11
12			! Possibilité de choisir le schéma pour l'advection de "q", en
13			! modifiant "iadv" dans "traceur.def".
14
15			! Pour Van-Leer plus vapeur d'eau saturée : iadv(1)=4
16			! Pour Van-Leer : iadv=10
17
18	guez	10	USE IOIPSL, only: ioconf_calendar, histclo
19	guez	3	use dimens_m, only: iim, jjm, llm, nqmx
20			use dimphy, only: klon
21			use paramet_m, only: ip1jm, ip1jmp1
22	guez	23	use comconst, only: daysec, cpp, dtvr, g, rad, r, initialize
23	guez	3	use comdissnew, only: lstardis, nitergdiv, nitergrot, niterh, tetagdiv, &
24			tetagrot, tetatemp
25			use conf_gcm_m, only: day_step, iperiod, anneeref, dayref, iecri, iphysiq, &
26			nday, raz_date, periodav, conf_gcm
27			use logic, only: iflag_phys
28			use comgeom, only: rlatu, aire_2d, cu_2d, cv_2d, rlonv
29	guez	25	use temps, only: day_ref, annee_ref, day_end, itau_dyn
30	guez	3	use com_io_dyn, only: histid, histvid, histaveid
31			use tracstoke, only: istdyn, istphy
32			use abort_gcm_m, only: abort_gcm
33			use inithist_m, only: inithist
34			use initdynav_m, only: initdynav
35	guez	25	use dynetat0_m, only: dynetat0, day_ini
36	guez	3	use grid_change, only: dyn_phy, init_dyn_phy
37	guez	18	use iniadvtrac_m, only: iniadvtrac
38	guez	3	use leapfrog_m, only: leapfrog
39	guez	7	use dynredem0_m, only: dynredem0
40	guez	13	use clesphys2, only: read_clesphys2
41	guez	25	use inigeom_m, only: inigeom
42	guez	3
43			IMPLICIT NONE
44
45			REAL zdtvr ! time step for dynamics, in s
46
47			! Variables dynamiques :
48			REAL vcov(ip1jm, llm), ucov(ip1jmp1, llm) ! vents covariants
49			REAL teta(ip1jmp1, llm) ! température potentielle
50			REAL q(ip1jmp1, llm, nqmx) ! champs advectés
51			REAL ps(ip1jmp1) ! pression au sol (Pa)
52
53			REAL masse(ip1jmp1, llm) ! masse d'air
54	guez	7	REAL phis(iim + 1, jjm + 1) ! géopotentiel au sol
55	guez	3
56			! Variables pour le fichier histoire :
57	guez	25	REAL time_0 ! time in day, as a fraction of day, in [0, 1[
58	guez	3
59			! Calendrier :
60			LOGICAL:: true_calendar = .false. ! default value
61
62			! Variables pour l'initialisation de la physique :
63			REAL zcufi(klon), zcvfi(klon) ! "cu" and "cv" values on the scalar grid
64			REAL latfi(klon), lonfi(klon)
65			REAL airefi(klon)
66
67			logical mask_v(iim + 1, jjm)
68			! (mask for points in the "v" grid, first index is for longitude,
69			! second index is for latitude)
70
71			namelist /main_nml/true_calendar
72
73			!------------------------------------------------------------
74
75			print *, "Enter namelist 'main_nml'."
76			read (unit=*, nml=main_nml)
77			write(unit=*, nml=main_nml)
78
79			! Initialisations:
80			call initialize
81
82			! Choix du calendrier :
83			if (true_calendar) then
84			call ioconf_calendar('gregorian')
85			else
86			call ioconf_calendar('360d')
87			endif
88
89			! Lecture des fichiers "gcm.def" ou "run.def" :
90	guez	13	call read_clesphys2
91			CALL conf_gcm
92	guez	3
93			! Initialisation des traceurs
94			! Choix du schéma pour l'advection dans le fichier "traceur.def" ou via INCA
95	guez	23	call iniadvtrac
96	guez	3
97			! Lecture du fichier "start.nc" :
98			CALL dynetat0(vcov, ucov, teta, q, masse, ps, phis, time_0)
99
100			! Lecture des paramètres de contrôle pour la simulation :
101			! on recalcule éventuellement le pas de temps
102			IF (MOD(day_step, iperiod) /= 0) THEN
103			call abort_gcm(modname = "gcm", message = &
104			'Il faut choisir un nombre de pas par jour multiple de "iperiod".', &
105			ierr = 1)
106			ENDIF
107
108			IF (MOD(day_step,iphysiq)/=0) THEN
109			call abort_gcm(modname = "gcm", message = &
110			'Il faut choisir un nombre de pas par jour multiple de "iphysiq".', &
111			ierr = 1)
112			ENDIF
113
114			! On remet le calendrier à zero si demandé:
115			if (annee_ref /= anneeref .or. day_ref /= dayref) then
116			print *, 'Attention : les dates initiales lues dans le fichier ' // &
117			'"start" ne correspondent pas à celles lues dans "gcm.def".'
118	guez	15	if (raz_date) then
119	guez	3	print *, 'On réinitialise à la date lue dans "gcm.def".'
120			annee_ref = anneeref
121			day_ref = dayref
122			day_ini = dayref
123			itau_dyn = 0
124			time_0 = 0.
125	guez	15	else
126			print *, 'On garde les dates du fichier "start".'
127	guez	3	endif
128			ELSE
129	guez	15	raz_date = .false.
130	guez	3	endif
131
132			! Initialisation des constantes dynamiques :
133			zdtvr = daysec / REAL(day_step)
134			IF (dtvr /= zdtvr) THEN
135			print *, 'Warning: the time steps in the ".def" file and in ' // &
136			'"start.nc" are different'
137			print *, 'dtvr (from "start.nc") = ', dtvr
138			print *, 'zdtvr (from ".def") = ', zdtvr
139			print *, 'Using the value from the ".def" file.'
140			dtvr = zdtvr
141			ENDIF
142			CALL iniconst
143
144			! Initialisation de la géometrie :
145			CALL inigeom
146
147			! Initialisation du filtre :
148			CALL inifilr
149
150			! Initialisation de la dissipation :
151			CALL inidissip(lstardis, nitergdiv, nitergrot, niterh, tetagdiv, tetagrot, &
152			tetatemp)
153
154			call init_dyn_phy
155
156			! Initialisation de la physique :
157			IF (iflag_phys == 1) THEN
158			latfi(1)=rlatu(1)
159			latfi(2:klon-1) = pack(spread(rlatu(2:jjm), 1, iim), .true.)
160			latfi(klon)= rlatu(jjm + 1)
161
162			lonfi(1)=0.
163			lonfi(2:klon-1) = pack(spread(rlonv(:iim), 2, jjm - 1), .true.)
164			lonfi(klon)= 0.
165
166			zcufi = pack(cu_2d, dyn_phy)
167
168			! Construct a mask for points in the "v" grid:
169			mask_v = .true.
170			mask_v(2:, 1) = .false.
171			mask_v(iim + 1, 2:) = .false.
172
173			zcvfi(:klon - 1) = pack(cv_2d, mask_v)
174			zcvfi(klon) = cv_2d(1, jjm)
175			! (that value of "cv_2d" is used twice in "zcvfi")
176
177			airefi = pack(aire_2d, dyn_phy)
178			print *, 'Attention : vitesse verticale nulle dans la physique.'
179	guez	15	CALL iniphysiq(klon, llm, latfi, lonfi, airefi, &
180	guez	3	zcufi, zcvfi, rad, g, r, cpp)
181			ENDIF
182
183			! Initialisation des entrées-sorties :
184			day_end = day_ini + nday
185			print *, "day_ini = ", day_ini
186			print *, "day_end = ", day_end
187
188	guez	5	CALL dynredem0("restart.nc", day_end, phis)
189	guez	3	CALL inithist(day_ref, annee_ref, zdtvr, nqmx, histid, histvid, &
190			infile="dyn_hist.nc", t_ops = iecri * daysec, t_wrt = iecri * daysec)
191			CALL initdynav(day_ref, annee_ref, zdtvr, nqmx, histaveid, &
192			infile='dyn_hist_ave.nc', t_ops = iperiod * zdtvr, &
193			t_wrt = periodav * daysec)
194
195			! Choix des fréquences de stockage pour le hors-ligne :
196			istdyn = day_step / 4 ! stockage toutes les 6 h = 1 jour / 4
197			istphy = istdyn / iphysiq
198
199			! Intégration temporelle du modèle :
200	guez	23	CALL leapfrog(ucov, vcov, teta, ps, masse, phis, q, time_0)
201	guez	3
202	guez	10	call histclo
203			print *, 'Simulation finished'
204			print *, 'Everything is cool'
205
206	guez	3	END PROGRAM gcm