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

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

Parent Directory Parent Directory | Revision Log Revision Log


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

  ViewVC Help
Powered by ViewVC 1.1.21