/[lmdze]/trunk/Sources/dyn3d/Guide/guide.f
ViewVC logotype

Contents of /trunk/Sources/dyn3d/Guide/guide.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 37 - (show annotations)
Tue Dec 21 15:45:48 2010 UTC (13 years, 5 months ago) by guez
Original Path: trunk/libf/dyn3d/guide.f90
File size: 10740 byte(s)
Inlined procedure "pression".

Split "guide.f90" into "guide.f90" and "tau2alpha.f90". Split
"read_reanalyse.f" into single-procedure files in directory
"Read_reanalyse".

Useless copy of variables in "iniphysiq". Directly define module
variables in "gcm" and remove procedure "iniphysiq".

Added "pressure-altitude" in "test_disvert".

1 MODULE guide_m
2
3 ! From dyn3d/guide.F, version 1.3 2005/05/25 13:10:09
4 ! and dyn3d/guide.h, version 1.1.1.1 2004/05/19 12:53:06
5
6 IMPLICIT NONE
7
8 REAL tau_min_u, tau_max_u
9 REAL tau_min_v, tau_max_v
10 REAL tau_min_t, tau_max_t
11 REAL tau_min_q, tau_max_q
12 REAL tau_min_p, tau_max_p
13 REAL aire_min, aire_max
14
15
16 LOGICAL guide_u, guide_v, guide_t, guide_q, guide_p
17 LOGICAL ncep, ini_anal
18 INTEGER online
19
20 CONTAINS
21
22 SUBROUTINE guide(itau, ucov, vcov, teta, q, masse, ps)
23
24 ! Author: F.Hourdin
25
26 USE dimens_m, ONLY : jjm, llm
27 USE paramet_m, ONLY : iip1, ip1jm, ip1jmp1, jjp1, llmp1
28 USE comconst, ONLY : cpp, daysec, dtvr, kappa, pi
29 USE comvert, ONLY : ap, bp, preff, presnivs
30 USE conf_gcm_m, ONLY : day_step, iperiod
31 USE comgeom, ONLY : aire, rlatu, rlonv
32 USE serre, ONLY : clat, clon
33 USE q_sat_m, ONLY : q_sat
34 USE exner_hyb_m, ONLY : exner_hyb
35 USE inigrads_m, ONLY : inigrads
36 use netcdf, only: nf90_nowrite, nf90_open, nf90_close
37 use tau2alpha_m, only: tau2alpha
38
39 INCLUDE 'netcdf.inc'
40
41 ! variables dynamiques
42 REAL vcov(ip1jm, llm), ucov(ip1jmp1, llm) ! vents covariants
43 REAL, intent(inout):: teta(ip1jmp1, llm) ! temperature potentielle
44 REAL q(ip1jmp1, llm) ! temperature potentielle
45 REAL ps(ip1jmp1) ! pression au sol
46 REAL masse(ip1jmp1, llm) ! masse d'air
47
48 ! common passe pour des sorties
49 REAL dxdys(iip1, jjp1), dxdyu(iip1, jjp1), dxdyv(iip1, jjm)
50 COMMON /comdxdy/dxdys, dxdyu, dxdyv
51
52 ! variables dynamiques pour les reanalyses.
53 REAL ucovrea1(ip1jmp1, llm), vcovrea1(ip1jm, llm) !vts cov reas
54 REAL tetarea1(ip1jmp1, llm) ! temp pot reales
55 REAL qrea1(ip1jmp1, llm) ! temp pot reales
56 REAL psrea1(ip1jmp1) ! ps
57 REAL ucovrea2(ip1jmp1, llm), vcovrea2(ip1jm, llm) !vts cov reas
58 REAL tetarea2(ip1jmp1, llm) ! temp pot reales
59 REAL qrea2(ip1jmp1, llm) ! temp pot reales
60 REAL masserea2(ip1jmp1, llm) ! masse
61 REAL psrea2(ip1jmp1) ! ps
62
63 REAL alpha_q(ip1jmp1)
64 REAL alpha_t(ip1jmp1), alpha_p(ip1jmp1)
65 REAL alpha_u(ip1jmp1), alpha_v(ip1jm)
66 REAL dday_step, toto, reste, itau_test
67 INTEGER step_rea, count_no_rea
68
69 INTEGER ilon, ilat
70 REAL factt, ztau(ip1jmp1)
71
72 INTEGER, INTENT (IN) :: itau
73 INTEGER ij, l
74 INTEGER ncidpl, varidpl, nlev, status
75 INTEGER rcod, rid
76 REAL ditau, tau, a
77 SAVE nlev
78
79 ! TEST SUR QSAT
80 REAL p(ip1jmp1, llmp1), pk(ip1jmp1, llm), pks(ip1jmp1)
81 REAL pkf(ip1jmp1, llm)
82 REAL pres(ip1jmp1, llm)
83
84 REAL qsat(ip1jmp1, llm)
85 REAL unskap
86 REAL tnat(ip1jmp1, llm)
87
88 LOGICAL:: first = .TRUE.
89
90 SAVE ucovrea1, vcovrea1, tetarea1, psrea1, qrea1
91 SAVE ucovrea2, vcovrea2, tetarea2, masserea2, psrea2, qrea2
92
93 SAVE alpha_t, alpha_q, alpha_u, alpha_v, alpha_p, itau_test
94 SAVE step_rea, count_no_rea
95
96 CHARACTER (10) file
97 INTEGER igrads
98 REAL dtgrads
99 SAVE igrads, dtgrads
100 DATA igrads, dtgrads/2, 100./
101
102 !-----------------------------------------------------------------------
103
104 PRINT *, 'Call sequence information: guide'
105
106 ! calcul de l'humidite saturante
107
108 forall (l = 1: llm + 1) p(:, l) = ap(l) + bp(l) * ps
109 CALL massdair(p, masse)
110 PRINT *, 'OK1'
111 CALL exner_hyb(ps, p, pks, pk, pkf)
112 PRINT *, 'OK2'
113 tnat(:, :) = pk(:, :)*teta(:, :)/cpp
114 PRINT *, 'OK3'
115 unskap = 1./kappa
116 pres(:, :) = preff*(pk(:, :)/cpp)**unskap
117 PRINT *, 'OK4'
118 qsat = q_sat(tnat, pres)
119
120 ! initialisations pour la lecture des reanalyses.
121 ! alpha determine la part des injections de donnees a chaque etape
122 ! alpha=1 signifie pas d'injection
123 ! alpha=0 signifie injection totale
124
125 PRINT *, 'ONLINE=', online
126 IF (online==-1) THEN
127 RETURN
128 END IF
129
130 IF (first) THEN
131
132 PRINT *, 'initialisation du guide '
133 CALL conf_guide
134 PRINT *, 'apres conf_guide'
135
136 file = 'guide'
137 CALL inigrads(igrads, rlonv, 180./pi, -180., 180., rlatu, -90., 90., &
138 180./pi, presnivs, 1., dtgrads, file, 'dyn_zon ')
139
140 PRINT *, &
141 '1: en-ligne, 0: hors-ligne (x=x_rea), -1: climat (x=x_gcm)'
142
143 IF (online==-1) RETURN
144 IF (online==1) THEN
145
146 ! Constantes de temps de rappel en jour
147 ! 0.1 c'est en gros 2h30.
148 ! 1e10 est une constante infinie donc en gros pas de guidage
149
150 ! coordonnees du centre du zoom
151 CALL coordij(clon, clat, ilon, ilat)
152 ! aire de la maille au centre du zoom
153 aire_min = aire(ilon+(ilat-1)*iip1)
154 ! aire maximale de la maille
155 aire_max = 0.
156 DO ij = 1, ip1jmp1
157 aire_max = max(aire_max, aire(ij))
158 END DO
159 ! factt = pas de temps en fraction de jour
160 factt = dtvr*iperiod/daysec
161
162 CALL tau2alpha(3, iip1, jjm, factt, tau_min_v, tau_max_v, alpha_v)
163 CALL tau2alpha(2, iip1, jjp1, factt, tau_min_u, tau_max_u, alpha_u)
164 CALL tau2alpha(1, iip1, jjp1, factt, tau_min_t, tau_max_t, alpha_t)
165 CALL tau2alpha(1, iip1, jjp1, factt, tau_min_p, tau_max_p, alpha_p)
166 CALL tau2alpha(1, iip1, jjp1, factt, tau_min_q, tau_max_q, alpha_q)
167
168 CALL dump2d(iip1, jjp1, aire, 'AIRE MAILLe ')
169 CALL dump2d(iip1, jjp1, alpha_u, 'COEFF U ')
170 CALL dump2d(iip1, jjp1, alpha_t, 'COEFF T ')
171
172 ! Cas ou on force exactement par les variables analysees
173 ELSE
174 alpha_t = 0.
175 alpha_u = 0.
176 alpha_v = 0.
177 alpha_p = 0.
178 ! physic=.false.
179 END IF
180
181 itau_test = 1001
182 step_rea = 1
183 count_no_rea = 0
184 ncidpl = -99
185
186 ! itau_test montre si l'importation a deja ete faite au rang itau
187 ! lecture d'un fichier netcdf pour determiner le nombre de niveaux
188 if (guide_u) then
189 if (ncidpl.eq.-99) rcod=nf90_open('u.nc',Nf90_NOWRITe,ncidpl)
190 endif
191
192 if (guide_v) then
193 if (ncidpl.eq.-99) rcod=nf90_open('v.nc',nf90_nowrite,ncidpl)
194 endif
195
196 if (guide_T) then
197 if (ncidpl.eq.-99) rcod=nf90_open('T.nc',nf90_nowrite,ncidpl)
198 endif
199
200 if (guide_Q) then
201 if (ncidpl.eq.-99) rcod=nf90_open('hur.nc',nf90_nowrite, ncidpl)
202 endif
203
204 IF (ncep) THEN
205 status = nf_inq_dimid(ncidpl, 'LEVEL', rid)
206 ELSE
207 status = nf_inq_dimid(ncidpl, 'PRESSURE', rid)
208 END IF
209 status = nf_inq_dimlen(ncidpl, rid, nlev)
210 PRINT *, 'nlev', nlev
211 rcod = nf90_close(ncidpl)
212 ! Lecture du premier etat des reanalyses.
213 CALL read_reanalyse(1, ps, ucovrea2, vcovrea2, tetarea2, qrea2, &
214 masserea2, psrea2, 1, nlev)
215 qrea2(:, :) = max(qrea2(:, :), 0.1)
216
217
218 ! Debut de l'integration temporelle:
219 END IF ! first
220
221 ! IMPORTATION DES VENTS, PRESSION ET TEMPERATURE REELS:
222
223 ditau = real(itau)
224 dday_step = real(day_step)
225 WRITE (*, *) 'ditau, dday_step'
226 WRITE (*, *) ditau, dday_step
227 toto = 4*ditau/dday_step
228 reste = toto - aint(toto)
229
230 IF (reste==0.) THEN
231 IF (itau_test==itau) THEN
232 WRITE (*, *) 'deuxieme passage de advreel a itau=', itau
233 STOP
234 ELSE
235 vcovrea1(:, :) = vcovrea2(:, :)
236 ucovrea1(:, :) = ucovrea2(:, :)
237 tetarea1(:, :) = tetarea2(:, :)
238 qrea1(:, :) = qrea2(:, :)
239
240 PRINT *, 'LECTURE REANALYSES, pas ', step_rea, 'apres ', &
241 count_no_rea, ' non lectures'
242 step_rea = step_rea + 1
243 itau_test = itau
244 CALL read_reanalyse(step_rea, ps, ucovrea2, vcovrea2, tetarea2, &
245 qrea2, masserea2, psrea2, 1, nlev)
246 qrea2(:, :) = max(qrea2(:, :), 0.1)
247 factt = dtvr*iperiod/daysec
248 ztau(:) = factt/max(alpha_t(:), 1.E-10)
249 CALL wrgrads(igrads, 1, aire, 'aire ', 'aire ')
250 CALL wrgrads(igrads, 1, dxdys, 'dxdy ', 'dxdy ')
251 CALL wrgrads(igrads, 1, alpha_u, 'au ', 'au ')
252 CALL wrgrads(igrads, 1, alpha_t, 'at ', 'at ')
253 CALL wrgrads(igrads, 1, ztau, 'taut ', 'taut ')
254 CALL wrgrads(igrads, llm, ucov, 'u ', 'u ')
255 CALL wrgrads(igrads, llm, ucovrea2, 'ua ', 'ua ')
256 CALL wrgrads(igrads, llm, teta, 'T ', 'T ')
257 CALL wrgrads(igrads, llm, tetarea2, 'Ta ', 'Ta ')
258 CALL wrgrads(igrads, llm, qrea2, 'Qa ', 'Qa ')
259 CALL wrgrads(igrads, llm, q, 'Q ', 'Q ')
260
261 CALL wrgrads(igrads, llm, qsat, 'QSAT ', 'QSAT ')
262
263 END IF
264 ELSE
265 count_no_rea = count_no_rea + 1
266 END IF
267
268 ! Guidage
269 ! x_gcm = a * x_gcm + (1-a) * x_reanalyses
270
271 IF (ini_anal) PRINT *, 'ATTENTION !!! ON PART DU GUIDAGE'
272
273 ditau = real(itau)
274 dday_step = real(day_step)
275
276
277 tau = 4*ditau/dday_step
278 tau = tau - aint(tau)
279
280 ! ucov
281 IF (guide_u) THEN
282 DO l = 1, llm
283 DO ij = 1, ip1jmp1
284 a = (1.-tau)*ucovrea1(ij, l) + tau*ucovrea2(ij, l)
285 ucov(ij, l) = (1.-alpha_u(ij))*ucov(ij, l) + alpha_u(ij)*a
286 IF (first .AND. ini_anal) ucov(ij, l) = a
287 END DO
288 END DO
289 END IF
290
291 IF (guide_t) THEN
292 DO l = 1, llm
293 DO ij = 1, ip1jmp1
294 a = (1.-tau)*tetarea1(ij, l) + tau*tetarea2(ij, l)
295 teta(ij, l) = (1.-alpha_t(ij))*teta(ij, l) + alpha_t(ij)*a
296 IF (first .AND. ini_anal) teta(ij, l) = a
297 END DO
298 END DO
299 END IF
300
301 ! P
302 IF (guide_p) THEN
303 DO ij = 1, ip1jmp1
304 a = (1.-tau)*psrea1(ij) + tau*psrea2(ij)
305 ps(ij) = (1.-alpha_p(ij))*ps(ij) + alpha_p(ij)*a
306 IF (first .AND. ini_anal) ps(ij) = a
307 END DO
308 forall (l = 1: llm + 1) p(:, l) = ap(l) + bp(l) * ps
309 CALL massdair(p, masse)
310 END IF
311
312
313 ! q
314 IF (guide_q) THEN
315 DO l = 1, llm
316 DO ij = 1, ip1jmp1
317 a = (1.-tau)*qrea1(ij, l) + tau*qrea2(ij, l)
318 ! hum relative en % -> hum specif
319 a = qsat(ij, l)*a*0.01
320 q(ij, l) = (1.-alpha_q(ij))*q(ij, l) + alpha_q(ij)*a
321 IF (first .AND. ini_anal) q(ij, l) = a
322 END DO
323 END DO
324 END IF
325
326 ! vcov
327 IF (guide_v) THEN
328 DO l = 1, llm
329 DO ij = 1, ip1jm
330 a = (1.-tau)*vcovrea1(ij, l) + tau*vcovrea2(ij, l)
331 vcov(ij, l) = (1.-alpha_v(ij))*vcov(ij, l) + alpha_v(ij)*a
332 IF (first .AND. ini_anal) vcov(ij, l) = a
333 END DO
334 IF (first .AND. ini_anal) vcov(ij, l) = a
335 END DO
336 END IF
337
338 first = .FALSE.
339
340 END SUBROUTINE guide
341
342 END MODULE guide_m

  ViewVC Help
Powered by ViewVC 1.1.21