/[lmdze]/trunk/Sources/phylmd/phyetat0.f
ViewVC logotype

Contents of /trunk/Sources/phylmd/phyetat0.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 191 - (show annotations)
Mon May 9 19:56:28 2016 UTC (8 years ago) by guez
File size: 11400 byte(s)
Extracted the call to read_comdissnew out of conf_gcm.

Made ok_instan a variable of module clesphys, itau_phy a variable of
module phyetat0_m, nid_ins a variable of module ini_histins_m, itap a
variable of new module time_phylmdz, so that histwrite_phy can be
called from any procedure without the need to cascade those variables
into that procedure. Made itau_w a variable of module time_phylmdz so
that it is computed only once per time step of physics.

Extracted variables of module clesphys which were in namelist
conf_phys_nml into their own namelist, clesphys_nml, and created
procedure read_clesphys reading clesphys_nml, to avoid side effect.

No need for double precision in procedure getso4fromfile. Assume there
is a single variable for the whole year in the NetCDF file instead of
one variable per month.

Created generic procedure histwrite_phy and removed procedure
write_histins, following LMDZ. histwrite_phy has only two arguments,
can be called from anywhere, and should manage the logic of writing or
not writing into various history files with various operations. So the
test on ok_instan goes inside histwrite_phy.

Test for raz_date in phyetat0 instead of physiq to avoid side effect.

Created procedure increment_itap to avoid side effect.

Removed unnecessary differences between procedures readsulfate and
readsulfate_pi.

1 module phyetat0_m
2
3 use dimphy, only: klon
4
5 IMPLICIT none
6
7 REAL, save:: rlat(klon), rlon(klon)
8 ! latitude and longitude of a point of the scalar grid identified
9 ! by a simple index, in degrees
10
11 integer, save:: itau_phy
12
13 private klon
14
15 contains
16
17 SUBROUTINE phyetat0(pctsrf, tsol, tsoil, qsurf, qsol, snow, albe, evap, &
18 rain_fall, snow_fall, solsw, sollw, fder, radsol, frugs, agesno, zmea, &
19 zstd, zsig, zgam, zthe, zpic, zval, t_ancien, q_ancien, ancien_ok, &
20 rnebcon, ratqs, clwcon, run_off_lic_0, sig1, w01, ncid_startphy)
21
22 ! From phylmd/phyetat0.F, version 1.4 2005/06/03 10:03:07
23 ! Author: Z.X. Li (LMD/CNRS)
24 ! Date: 1993/08/18
25 ! Objet : lecture de l'état initial pour la physique
26
27 USE conf_gcm_m, ONLY: raz_date
28 use dimphy, only: zmasq, klev
29 USE dimsoil, ONLY : nsoilmx
30 USE indicesol, ONLY : epsfra, is_lic, is_oce, is_sic, is_ter, nbsrf
31 use netcdf, only: nf90_global, nf90_inq_varid, NF90_NOERR, NF90_NOWRITE
32 use netcdf95, only: nf95_get_att, nf95_get_var, nf95_inq_varid, &
33 nf95_inquire_variable, NF95_OPEN
34
35 REAL, intent(out):: pctsrf(klon, nbsrf)
36 REAL, intent(out):: tsol(klon, nbsrf)
37 REAL, intent(out):: tsoil(klon, nsoilmx, nbsrf)
38 REAL, intent(out):: qsurf(klon, nbsrf)
39 REAL, intent(out):: qsol(:) ! (klon)
40 REAL, intent(out):: snow(klon, nbsrf)
41 REAL, intent(out):: albe(klon, nbsrf)
42 REAL, intent(out):: evap(klon, nbsrf)
43 REAL, intent(out):: rain_fall(klon)
44 REAL, intent(out):: snow_fall(klon)
45 real, intent(out):: solsw(klon)
46 REAL, intent(out):: sollw(klon)
47 real, intent(out):: fder(klon)
48 REAL, intent(out):: radsol(klon)
49 REAL, intent(out):: frugs(klon, nbsrf)
50 REAL, intent(out):: agesno(klon, nbsrf)
51 REAL, intent(out):: zmea(klon)
52 REAL, intent(out):: zstd(klon)
53 REAL, intent(out):: zsig(klon)
54 REAL, intent(out):: zgam(klon)
55 REAL, intent(out):: zthe(klon)
56 REAL, intent(out):: zpic(klon)
57 REAL, intent(out):: zval(klon)
58 REAL, intent(out):: t_ancien(klon, klev), q_ancien(klon, klev)
59 LOGICAL, intent(out):: ancien_ok
60 real, intent(out):: rnebcon(klon, klev), ratqs(klon, klev)
61 REAL, intent(out):: clwcon(klon, klev), run_off_lic_0(klon)
62 real, intent(out):: sig1(klon, klev) ! section adiabatic updraft
63
64 real, intent(out):: w01(klon, klev)
65 ! vertical velocity within adiabatic updraft
66
67 integer, intent(out):: ncid_startphy
68
69 ! Local:
70 REAL fractint(klon)
71 INTEGER varid, ndims
72 INTEGER ierr, i
73
74 !---------------------------------------------------------------
75
76 print *, "Call sequence information: phyetat0"
77
78 ! Fichier contenant l'état initial :
79 call NF95_OPEN("startphy.nc", NF90_NOWRITE, ncid_startphy)
80
81 IF (raz_date) then
82 itau_phy = 0
83 else
84 call nf95_get_att(ncid_startphy, nf90_global, "itau_phy", itau_phy)
85 end IF
86
87 ! Lecture des latitudes (coordonnees):
88
89 call NF95_INQ_VARID(ncid_startphy, "latitude", varid)
90 call NF95_GET_VAR(ncid_startphy, varid, rlat)
91
92 ! Lecture des longitudes (coordonnees):
93
94 call NF95_INQ_VARID(ncid_startphy, "longitude", varid)
95 call NF95_GET_VAR(ncid_startphy, varid, rlon)
96
97 ! Lecture du masque terre mer
98
99 call NF95_INQ_VARID(ncid_startphy, "masque", varid)
100 call nf95_get_var(ncid_startphy, varid, zmasq)
101
102 ! Lecture des fractions pour chaque sous-surface
103
104 ! initialisation des sous-surfaces
105
106 pctsrf = 0.
107
108 ! fraction de terre
109
110 ierr = NF90_INQ_VARID(ncid_startphy, "FTER", varid)
111 IF (ierr == NF90_NOERR) THEN
112 call nf95_get_var(ncid_startphy, varid, pctsrf(:, is_ter))
113 else
114 PRINT *, 'phyetat0: Le champ <FTER> est absent'
115 ENDIF
116
117 ! fraction de glace de terre
118
119 ierr = NF90_INQ_VARID(ncid_startphy, "FLIC", varid)
120 IF (ierr == NF90_NOERR) THEN
121 call nf95_get_var(ncid_startphy, varid, pctsrf(:, is_lic))
122 else
123 PRINT *, 'phyetat0: Le champ <FLIC> est absent'
124 ENDIF
125
126 ! fraction d'ocean
127
128 ierr = NF90_INQ_VARID(ncid_startphy, "FOCE", varid)
129 IF (ierr == NF90_NOERR) THEN
130 call nf95_get_var(ncid_startphy, varid, pctsrf(:, is_oce))
131 else
132 PRINT *, 'phyetat0: Le champ <FOCE> est absent'
133 ENDIF
134
135 ! fraction glace de mer
136
137 ierr = NF90_INQ_VARID(ncid_startphy, "FSIC", varid)
138 IF (ierr == NF90_NOERR) THEN
139 call nf95_get_var(ncid_startphy, varid, pctsrf(:, is_sic))
140 else
141 PRINT *, 'phyetat0: Le champ <FSIC> est absent'
142 ENDIF
143
144 ! Verification de l'adequation entre le masque et les sous-surfaces
145
146 fractint = pctsrf(:, is_ter) + pctsrf(:, is_lic)
147 DO i = 1 , klon
148 IF ( abs(fractint(i) - zmasq(i) ) > EPSFRA ) THEN
149 WRITE(*, *) 'phyetat0: attention fraction terre pas ', &
150 'coherente ', i, zmasq(i), pctsrf(i, is_ter) &
151 , pctsrf(i, is_lic)
152 ENDIF
153 END DO
154 fractint = pctsrf(:, is_oce) + pctsrf(:, is_sic)
155 DO i = 1 , klon
156 IF ( abs( fractint(i) - (1. - zmasq(i))) > EPSFRA ) THEN
157 WRITE(*, *) 'phyetat0 attention fraction ocean pas ', &
158 'coherente ', i, zmasq(i) , pctsrf(i, is_oce) &
159 , pctsrf(i, is_sic)
160 ENDIF
161 END DO
162
163 ! Lecture des temperatures du sol:
164 call NF95_INQ_VARID(ncid_startphy, "TS", varid)
165 call nf95_inquire_variable(ncid_startphy, varid, ndims = ndims)
166 if (ndims == 2) then
167 call NF95_GET_VAR(ncid_startphy, varid, tsol)
168 else
169 print *, "Found only one surface type for soil temperature."
170 call nf95_get_var(ncid_startphy, varid, tsol(:, 1))
171 tsol(:, 2:nbsrf) = spread(tsol(:, 1), dim = 2, ncopies = nbsrf - 1)
172 end if
173
174 ! Lecture des temperatures du sol profond:
175
176 call NF95_INQ_VARID(ncid_startphy, 'Tsoil', varid)
177 call NF95_GET_VAR(ncid_startphy, varid, tsoil)
178
179 ! Lecture de l'humidite de l'air juste au dessus du sol:
180
181 call NF95_INQ_VARID(ncid_startphy, "QS", varid)
182 call nf95_get_var(ncid_startphy, varid, qsurf)
183
184 ! Eau dans le sol (pour le modele de sol "bucket")
185
186 ierr = NF90_INQ_VARID(ncid_startphy, "QSOL", varid)
187 IF (ierr == NF90_NOERR) THEN
188 call nf95_get_var(ncid_startphy, varid, qsol)
189 else
190 PRINT *, 'phyetat0: Le champ <QSOL> est absent'
191 PRINT *, ' Valeur par defaut nulle'
192 qsol = 0.
193 ENDIF
194
195 ! Lecture de neige au sol:
196
197 call NF95_INQ_VARID(ncid_startphy, "SNOW", varid)
198 call nf95_get_var(ncid_startphy, varid, snow)
199
200 ! Lecture de albedo au sol:
201
202 call NF95_INQ_VARID(ncid_startphy, "ALBE", varid)
203 call nf95_get_var(ncid_startphy, varid, albe)
204
205 ! Lecture de evaporation:
206
207 call NF95_INQ_VARID(ncid_startphy, "EVAP", varid)
208 call nf95_get_var(ncid_startphy, varid, evap)
209
210 ! Lecture precipitation liquide:
211
212 call NF95_INQ_VARID(ncid_startphy, "rain_f", varid)
213 call NF95_GET_VAR(ncid_startphy, varid, rain_fall)
214
215 ! Lecture precipitation solide:
216
217 call NF95_INQ_VARID(ncid_startphy, "snow_f", varid)
218 call NF95_GET_VAR(ncid_startphy, varid, snow_fall)
219
220 ! Lecture rayonnement solaire au sol:
221
222 ierr = NF90_INQ_VARID(ncid_startphy, "solsw", varid)
223 IF (ierr /= NF90_NOERR) THEN
224 PRINT *, 'phyetat0: Le champ <solsw> est absent'
225 PRINT *, 'mis a zero'
226 solsw = 0.
227 ELSE
228 call nf95_get_var(ncid_startphy, varid, solsw)
229 ENDIF
230
231 ! Lecture rayonnement IF au sol:
232
233 ierr = NF90_INQ_VARID(ncid_startphy, "sollw", varid)
234 IF (ierr /= NF90_NOERR) THEN
235 PRINT *, 'phyetat0: Le champ <sollw> est absent'
236 PRINT *, 'mis a zero'
237 sollw = 0.
238 ELSE
239 call nf95_get_var(ncid_startphy, varid, sollw)
240 ENDIF
241
242 ! Lecture derive des flux:
243
244 ierr = NF90_INQ_VARID(ncid_startphy, "fder", varid)
245 IF (ierr /= NF90_NOERR) THEN
246 PRINT *, 'phyetat0: Le champ <fder> est absent'
247 PRINT *, 'mis a zero'
248 fder = 0.
249 ELSE
250 call nf95_get_var(ncid_startphy, varid, fder)
251 ENDIF
252
253 ! Lecture du rayonnement net au sol:
254
255 call NF95_INQ_VARID(ncid_startphy, "RADS", varid)
256 call NF95_GET_VAR(ncid_startphy, varid, radsol)
257
258 ! Lecture de la longueur de rugosite
259
260 call NF95_INQ_VARID(ncid_startphy, "RUG", varid)
261 call nf95_get_var(ncid_startphy, varid, frugs)
262
263 ! Lecture de l'age de la neige:
264
265 call NF95_INQ_VARID(ncid_startphy, "AGESNO", varid)
266 call nf95_get_var(ncid_startphy, varid, agesno)
267
268 call NF95_INQ_VARID(ncid_startphy, "ZMEA", varid)
269 call NF95_GET_VAR(ncid_startphy, varid, zmea)
270
271 call NF95_INQ_VARID(ncid_startphy, "ZSTD", varid)
272 call NF95_GET_VAR(ncid_startphy, varid, zstd)
273
274 call NF95_INQ_VARID(ncid_startphy, "ZSIG", varid)
275 call NF95_GET_VAR(ncid_startphy, varid, zsig)
276
277 call NF95_INQ_VARID(ncid_startphy, "ZGAM", varid)
278 call NF95_GET_VAR(ncid_startphy, varid, zgam)
279
280 call NF95_INQ_VARID(ncid_startphy, "ZTHE", varid)
281 call NF95_GET_VAR(ncid_startphy, varid, zthe)
282
283 call NF95_INQ_VARID(ncid_startphy, "ZPIC", varid)
284 call NF95_GET_VAR(ncid_startphy, varid, zpic)
285
286 call NF95_INQ_VARID(ncid_startphy, "ZVAL", varid)
287 call NF95_GET_VAR(ncid_startphy, varid, zval)
288
289 ancien_ok = .TRUE.
290
291 ierr = NF90_INQ_VARID(ncid_startphy, "TANCIEN", varid)
292 IF (ierr /= NF90_NOERR) THEN
293 PRINT *, "phyetat0: Le champ <TANCIEN> est absent"
294 PRINT *, "Depart legerement fausse. Mais je continue"
295 ancien_ok = .FALSE.
296 ELSE
297 call nf95_get_var(ncid_startphy, varid, t_ancien)
298 ENDIF
299
300 ierr = NF90_INQ_VARID(ncid_startphy, "QANCIEN", varid)
301 IF (ierr /= NF90_NOERR) THEN
302 PRINT *, "phyetat0: Le champ <QANCIEN> est absent"
303 PRINT *, "Depart legerement fausse. Mais je continue"
304 ancien_ok = .FALSE.
305 ELSE
306 call nf95_get_var(ncid_startphy, varid, q_ancien)
307 ENDIF
308
309 ierr = NF90_INQ_VARID(ncid_startphy, "CLWCON", varid)
310 IF (ierr /= NF90_NOERR) THEN
311 PRINT *, "phyetat0: Le champ CLWCON est absent"
312 PRINT *, "Depart legerement fausse. Mais je continue"
313 clwcon = 0.
314 ELSE
315 call nf95_get_var(ncid_startphy, varid, clwcon(:, 1))
316 clwcon(:, 2:) = 0.
317 ENDIF
318
319 ierr = NF90_INQ_VARID(ncid_startphy, "RNEBCON", varid)
320 IF (ierr /= NF90_NOERR) THEN
321 PRINT *, "phyetat0: Le champ RNEBCON est absent"
322 PRINT *, "Depart legerement fausse. Mais je continue"
323 rnebcon = 0.
324 ELSE
325 call nf95_get_var(ncid_startphy, varid, rnebcon(:, 1))
326 rnebcon(:, 2:) = 0.
327 ENDIF
328
329 ! Lecture ratqs
330
331 ierr = NF90_INQ_VARID(ncid_startphy, "RATQS", varid)
332 IF (ierr /= NF90_NOERR) THEN
333 PRINT *, "phyetat0: Le champ <RATQS> est absent"
334 PRINT *, "Depart legerement fausse. Mais je continue"
335 ratqs = 0.
336 ELSE
337 call nf95_get_var(ncid_startphy, varid, ratqs(:, 1))
338 ratqs(:, 2:) = 0.
339 ENDIF
340
341 ! Lecture run_off_lic_0
342
343 ierr = NF90_INQ_VARID(ncid_startphy, "RUNOFFLIC0", varid)
344 IF (ierr /= NF90_NOERR) THEN
345 PRINT *, "phyetat0: Le champ <RUNOFFLIC0> est absent"
346 PRINT *, "Depart legerement fausse. Mais je continue"
347 run_off_lic_0 = 0.
348 ELSE
349 call nf95_get_var(ncid_startphy, varid, run_off_lic_0)
350 ENDIF
351
352 call nf95_inq_varid(ncid_startphy, "sig1", varid)
353 call nf95_get_var(ncid_startphy, varid, sig1)
354
355 call nf95_inq_varid(ncid_startphy, "w01", varid)
356 call nf95_get_var(ncid_startphy, varid, w01)
357
358 END SUBROUTINE phyetat0
359
360 end module phyetat0_m

  ViewVC Help
Powered by ViewVC 1.1.21