1 | #include "conf_gcm.F90" |
---|
2 | |
---|
3 | ! |
---|
4 | ! $Id: 1DUTILS.h 3686 2020-05-27 12:59:10Z idelkadi $ |
---|
5 | ! |
---|
6 | ! |
---|
7 | ! |
---|
8 | SUBROUTINE conf_unicol |
---|
9 | ! |
---|
10 | #ifdef CPP_IOIPSL |
---|
11 | use IOIPSL |
---|
12 | #else |
---|
13 | ! if not using IOIPSL, we still need to use (a local version of) getin |
---|
14 | use ioipsl_getincom |
---|
15 | #endif |
---|
16 | USE print_control_mod, ONLY: lunout |
---|
17 | IMPLICIT NONE |
---|
18 | !----------------------------------------------------------------------- |
---|
19 | ! Auteurs : A. Lahellec . |
---|
20 | ! |
---|
21 | ! Declarations : |
---|
22 | ! -------------- |
---|
23 | |
---|
24 | #include "compar1d.h" |
---|
25 | #include "flux_arp.h" |
---|
26 | #include "tsoilnudge.h" |
---|
27 | #include "fcg_gcssold.h" |
---|
28 | #include "fcg_racmo.h" |
---|
29 | ! |
---|
30 | ! |
---|
31 | ! local: |
---|
32 | ! ------ |
---|
33 | |
---|
34 | ! CHARACTER ch1*72,ch2*72,ch3*72,ch4*12 |
---|
35 | |
---|
36 | ! |
---|
37 | ! ------------------------------------------------------------------- |
---|
38 | ! |
---|
39 | ! ......... Initilisation parametres du lmdz1D .......... |
---|
40 | ! |
---|
41 | !--------------------------------------------------------------------- |
---|
42 | ! initialisations: |
---|
43 | ! ---------------- |
---|
44 | |
---|
45 | !Config Key = lunout |
---|
46 | !Config Desc = unite de fichier pour les impressions |
---|
47 | !Config Def = 6 |
---|
48 | !Config Help = unite de fichier pour les impressions |
---|
49 | !Config (defaut sortie standard = 6) |
---|
50 | lunout=6 |
---|
51 | ! CALL getin('lunout', lunout) |
---|
52 | IF (lunout /= 5 .and. lunout /= 6) THEN |
---|
53 | OPEN(lunout,FILE='lmdz.out') |
---|
54 | ENDIF |
---|
55 | |
---|
56 | !Config Key = prt_level |
---|
57 | !Config Desc = niveau d'impressions de debogage |
---|
58 | !Config Def = 0 |
---|
59 | !Config Help = Niveau d'impression pour le debogage |
---|
60 | !Config (0 = minimum d'impression) |
---|
61 | ! prt_level = 0 |
---|
62 | ! CALL getin('prt_level',prt_level) |
---|
63 | |
---|
64 | !----------------------------------------------------------------------- |
---|
65 | ! Parametres de controle du run: |
---|
66 | !----------------------------------------------------------------------- |
---|
67 | |
---|
68 | !Config Key = restart |
---|
69 | !Config Desc = on repart des startphy et start1dyn |
---|
70 | !Config Def = false |
---|
71 | !Config Help = les fichiers restart doivent etre renomme en start |
---|
72 | restart =.false. |
---|
73 | CALL getin('restart',restart) |
---|
74 | |
---|
75 | !Config Key = forcing_type |
---|
76 | !Config Desc = defines the way the SCM is forced: |
---|
77 | !Config Def = 0 |
---|
78 | !!Config Help = 0 ==> forcing_les = .true. |
---|
79 | ! initial profiles from file prof.inp.001 |
---|
80 | ! no forcing by LS convergence ; |
---|
81 | ! surface temperature imposed ; |
---|
82 | ! radiative cooling may be imposed (iflag_radia=0 in physiq.def) |
---|
83 | ! = 1 ==> forcing_radconv = .true. |
---|
84 | ! idem forcing_type = 0, but the imposed radiative cooling |
---|
85 | ! is set to 0 (hence, if iflag_radia=0 in physiq.def, |
---|
86 | ! then there is no radiative cooling at all) |
---|
87 | ! = 2 ==> forcing_toga = .true. |
---|
88 | ! initial profiles from TOGA-COARE IFA files |
---|
89 | ! LS convergence and SST imposed from TOGA-COARE IFA files |
---|
90 | ! = 3 ==> forcing_GCM2SCM = .true. |
---|
91 | ! initial profiles from the GCM output |
---|
92 | ! LS convergence imposed from the GCM output |
---|
93 | ! = 4 ==> forcing_twpi = .true. |
---|
94 | ! initial profiles from TWPICE nc files |
---|
95 | ! LS convergence and SST imposed from TWPICE nc files |
---|
96 | ! = 5 ==> forcing_rico = .true. |
---|
97 | ! initial profiles from RICO idealized |
---|
98 | ! LS convergence imposed from RICO (cst) |
---|
99 | ! = 6 ==> forcing_amma = .true. |
---|
100 | ! = 10 ==> forcing_case = .true. |
---|
101 | ! initial profiles from case.nc file |
---|
102 | ! = 40 ==> forcing_GCSSold = .true. |
---|
103 | ! initial profile from GCSS file |
---|
104 | ! LS convergence imposed from GCSS file |
---|
105 | ! = 50 ==> forcing_fire = .true. |
---|
106 | ! = 59 ==> forcing_sandu = .true. |
---|
107 | ! initial profiles from sanduref file: see prof.inp.001 |
---|
108 | ! SST varying with time and divergence constante: see ifa_sanduref.txt file |
---|
109 | ! Radiation has to be computed interactively |
---|
110 | ! = 60 ==> forcing_astex = .true. |
---|
111 | ! initial profiles from file: see prof.inp.001 |
---|
112 | ! SST,divergence,ug,vg,ufa,vfa varying with time : see ifa_astex.txt file |
---|
113 | ! Radiation has to be computed interactively |
---|
114 | ! = 61 ==> forcing_armcu = .true. |
---|
115 | ! initial profiles from file: see prof.inp.001 |
---|
116 | ! sensible and latent heat flux imposed: see ifa_arm_cu_1.txt |
---|
117 | ! large scale advective forcing & radiative tendencies applied below 1000m: see ifa_arm_cu_2.txt |
---|
118 | ! use geostrophic wind ug=10m/s vg=0m/s. Duration of the case 53100s |
---|
119 | ! Radiation to be switched off |
---|
120 | ! > 100 ==> forcing_case = .true. or forcing_case2 = .true. |
---|
121 | ! initial profiles from case.nc file |
---|
122 | ! |
---|
123 | forcing_type = 0 |
---|
124 | CALL getin('forcing_type',forcing_type) |
---|
125 | imp_fcg_gcssold = .false. |
---|
126 | ts_fcg_gcssold = .false. |
---|
127 | Tp_fcg_gcssold = .false. |
---|
128 | Tp_ini_gcssold = .false. |
---|
129 | xTurb_fcg_gcssold = .false. |
---|
130 | IF (forcing_type .eq.40) THEN |
---|
131 | CALL getin('imp_fcg',imp_fcg_gcssold) |
---|
132 | CALL getin('ts_fcg',ts_fcg_gcssold) |
---|
133 | CALL getin('tp_fcg',Tp_fcg_gcssold) |
---|
134 | CALL getin('tp_ini',Tp_ini_gcssold) |
---|
135 | CALL getin('turb_fcg',xTurb_fcg_gcssold) |
---|
136 | ENDIF |
---|
137 | |
---|
138 | !Parametres de forcage |
---|
139 | !Config Key = tend_t |
---|
140 | !Config Desc = forcage ou non par advection de T |
---|
141 | !Config Def = false |
---|
142 | !Config Help = forcage ou non par advection de T |
---|
143 | tend_t =0 |
---|
144 | CALL getin('tend_t',tend_t) |
---|
145 | |
---|
146 | !Config Key = tend_q |
---|
147 | !Config Desc = forcage ou non par advection de q |
---|
148 | !Config Def = false |
---|
149 | !Config Help = forcage ou non par advection de q |
---|
150 | tend_q =0 |
---|
151 | CALL getin('tend_q',tend_q) |
---|
152 | |
---|
153 | !Config Key = tend_u |
---|
154 | !Config Desc = forcage ou non par advection de u |
---|
155 | !Config Def = false |
---|
156 | !Config Help = forcage ou non par advection de u |
---|
157 | tend_u =0 |
---|
158 | CALL getin('tend_u',tend_u) |
---|
159 | |
---|
160 | !Config Key = tend_v |
---|
161 | !Config Desc = forcage ou non par advection de v |
---|
162 | !Config Def = false |
---|
163 | !Config Help = forcage ou non par advection de v |
---|
164 | tend_v =0 |
---|
165 | CALL getin('tend_v',tend_v) |
---|
166 | |
---|
167 | !Config Key = tend_w |
---|
168 | !Config Desc = forcage ou non par vitesse verticale |
---|
169 | !Config Def = false |
---|
170 | !Config Help = forcage ou non par vitesse verticale |
---|
171 | tend_w =0 |
---|
172 | CALL getin('tend_w',tend_w) |
---|
173 | |
---|
174 | !Config Key = tend_rayo |
---|
175 | !Config Desc = forcage ou non par dtrad |
---|
176 | !Config Def = false |
---|
177 | !Config Help = forcage ou non par dtrad |
---|
178 | tend_rayo =0 |
---|
179 | CALL getin('tend_rayo',tend_rayo) |
---|
180 | |
---|
181 | |
---|
182 | !Config Key = nudge_t |
---|
183 | !Config Desc = constante de nudging de T |
---|
184 | !Config Def = false |
---|
185 | !Config Help = constante de nudging de T |
---|
186 | nudge_t =0. |
---|
187 | CALL getin('nudge_t',nudge_t) |
---|
188 | |
---|
189 | !Config Key = nudge_q |
---|
190 | !Config Desc = constante de nudging de q |
---|
191 | !Config Def = false |
---|
192 | !Config Help = constante de nudging de q |
---|
193 | nudge_q =0. |
---|
194 | CALL getin('nudge_q',nudge_q) |
---|
195 | |
---|
196 | !Config Key = nudge_u |
---|
197 | !Config Desc = constante de nudging de u |
---|
198 | !Config Def = false |
---|
199 | !Config Help = constante de nudging de u |
---|
200 | nudge_u =0. |
---|
201 | CALL getin('nudge_u',nudge_u) |
---|
202 | |
---|
203 | !Config Key = nudge_v |
---|
204 | !Config Desc = constante de nudging de v |
---|
205 | !Config Def = false |
---|
206 | !Config Help = constante de nudging de v |
---|
207 | nudge_v =0. |
---|
208 | CALL getin('nudge_v',nudge_v) |
---|
209 | |
---|
210 | !Config Key = nudge_w |
---|
211 | !Config Desc = constante de nudging de w |
---|
212 | !Config Def = false |
---|
213 | !Config Help = constante de nudging de w |
---|
214 | nudge_w =0. |
---|
215 | CALL getin('nudge_w',nudge_w) |
---|
216 | |
---|
217 | |
---|
218 | !Config Key = iflag_nudge |
---|
219 | !Config Desc = atmospheric nudging ttype (decimal code) |
---|
220 | !Config Def = 0 |
---|
221 | !Config Help = 0 ==> no nudging |
---|
222 | ! If digit number n of iflag_nudge is set, then nudging of type n is on |
---|
223 | ! If digit number n of iflag_nudge is not set, then nudging of type n is off |
---|
224 | ! (digits are numbered from the right) |
---|
225 | iflag_nudge = 0 |
---|
226 | CALL getin('iflag_nudge',iflag_nudge) |
---|
227 | |
---|
228 | !Config Key = ok_flux_surf |
---|
229 | !Config Desc = forcage ou non par les flux de surface |
---|
230 | !Config Def = false |
---|
231 | !Config Help = forcage ou non par les flux de surface |
---|
232 | ok_flux_surf =.false. |
---|
233 | CALL getin('ok_flux_surf',ok_flux_surf) |
---|
234 | |
---|
235 | !Config Key = ok_prescr_ust |
---|
236 | !Config Desc = ustar impose ou non |
---|
237 | !Config Def = false |
---|
238 | !Config Help = ustar impose ou non |
---|
239 | ok_prescr_ust = .false. |
---|
240 | CALL getin('ok_prescr_ust',ok_prescr_ust) |
---|
241 | |
---|
242 | !Config Key = ok_old_disvert |
---|
243 | !Config Desc = utilisation de l ancien programme disvert0 (dans 1DUTILS.h) |
---|
244 | !Config Def = false |
---|
245 | !Config Help = utilisation de l ancien programme disvert0 (dans 1DUTILS.h) |
---|
246 | ok_old_disvert = .FALSE. |
---|
247 | CALL getin('ok_old_disvert',ok_old_disvert) |
---|
248 | |
---|
249 | !Config Key = time_ini |
---|
250 | !Config Desc = meaningless in this case |
---|
251 | !Config Def = 0. |
---|
252 | !Config Help = |
---|
253 | tsurf = 0. |
---|
254 | CALL getin('time_ini',time_ini) |
---|
255 | |
---|
256 | !Config Key = rlat et rlon |
---|
257 | !Config Desc = latitude et longitude |
---|
258 | !Config Def = 0.0 0.0 |
---|
259 | !Config Help = fixe la position de la colonne |
---|
260 | xlat = 0. |
---|
261 | xlon = 0. |
---|
262 | CALL getin('rlat',xlat) |
---|
263 | CALL getin('rlon',xlon) |
---|
264 | |
---|
265 | !Config Key = airephy |
---|
266 | !Config Desc = Grid cell area |
---|
267 | !Config Def = 1.e11 |
---|
268 | !Config Help = |
---|
269 | airefi = 1.e11 |
---|
270 | CALL getin('airephy',airefi) |
---|
271 | |
---|
272 | !Config Key = nat_surf |
---|
273 | !Config Desc = surface type |
---|
274 | !Config Def = 0 (ocean) |
---|
275 | !Config Help = 0=ocean,1=land,2=glacier,3=banquise |
---|
276 | nat_surf = 0. |
---|
277 | CALL getin('nat_surf',nat_surf) |
---|
278 | |
---|
279 | !Config Key = tsurf |
---|
280 | !Config Desc = surface temperature |
---|
281 | !Config Def = 290. |
---|
282 | !Config Help = not used if type_ts_forcing=1 in lmdz1d.F |
---|
283 | tsurf = 290. |
---|
284 | CALL getin('tsurf',tsurf) |
---|
285 | |
---|
286 | !Config Key = psurf |
---|
287 | !Config Desc = surface pressure |
---|
288 | !Config Def = 102400. |
---|
289 | !Config Help = |
---|
290 | psurf = 102400. |
---|
291 | CALL getin('psurf',psurf) |
---|
292 | |
---|
293 | !Config Key = zsurf |
---|
294 | !Config Desc = surface altitude |
---|
295 | !Config Def = 0. |
---|
296 | !Config Help = |
---|
297 | zsurf = 0. |
---|
298 | CALL getin('zsurf',zsurf) |
---|
299 | |
---|
300 | !Config Key = rugos |
---|
301 | !Config Desc = coefficient de frottement |
---|
302 | !Config Def = 0.0001 |
---|
303 | !Config Help = calcul du Cdrag |
---|
304 | rugos = 0.0001 |
---|
305 | CALL getin('rugos',rugos) |
---|
306 | ! FH/2020/04/08/confinement: Pour le nouveau format standard, la rugosite s'appelle z0 |
---|
307 | CALL getin('z0',rugos) |
---|
308 | |
---|
309 | !Config Key = rugosh |
---|
310 | !Config Desc = coefficient de frottement |
---|
311 | !Config Def = rugos |
---|
312 | !Config Help = calcul du Cdrag |
---|
313 | rugosh = rugos |
---|
314 | CALL getin('rugosh',rugosh) |
---|
315 | |
---|
316 | |
---|
317 | |
---|
318 | !Config Key = snowmass |
---|
319 | !Config Desc = mass de neige de la surface en kg/m2 |
---|
320 | !Config Def = 0.0000 |
---|
321 | !Config Help = snowmass |
---|
322 | snowmass = 0.0000 |
---|
323 | CALL getin('snowmass',snowmass) |
---|
324 | |
---|
325 | !Config Key = wtsurf et wqsurf |
---|
326 | !Config Desc = ??? |
---|
327 | !Config Def = 0.0 0.0 |
---|
328 | !Config Help = |
---|
329 | wtsurf = 0.0 |
---|
330 | wqsurf = 0.0 |
---|
331 | CALL getin('wtsurf',wtsurf) |
---|
332 | CALL getin('wqsurf',wqsurf) |
---|
333 | |
---|
334 | !Config Key = albedo |
---|
335 | !Config Desc = albedo |
---|
336 | !Config Def = 0.09 |
---|
337 | !Config Help = |
---|
338 | albedo = 0.09 |
---|
339 | CALL getin('albedo',albedo) |
---|
340 | |
---|
341 | !Config Key = agesno |
---|
342 | !Config Desc = age de la neige |
---|
343 | !Config Def = 30.0 |
---|
344 | !Config Help = |
---|
345 | xagesno = 30.0 |
---|
346 | CALL getin('agesno',xagesno) |
---|
347 | |
---|
348 | !Config Key = restart_runoff |
---|
349 | !Config Desc = age de la neige |
---|
350 | !Config Def = 30.0 |
---|
351 | !Config Help = |
---|
352 | restart_runoff = 0.0 |
---|
353 | CALL getin('restart_runoff',restart_runoff) |
---|
354 | |
---|
355 | !Config Key = qsolinp |
---|
356 | !Config Desc = initial bucket water content (kg/m2) when land (5std) |
---|
357 | !Config Def = 30.0 |
---|
358 | !Config Help = |
---|
359 | qsolinp = 1. |
---|
360 | CALL getin('qsolinp',qsolinp) |
---|
361 | |
---|
362 | !Config Key = zpicinp |
---|
363 | !Config Desc = denivellation orographie |
---|
364 | !Config Def = 0. |
---|
365 | !Config Help = input brise |
---|
366 | zpicinp = 0. |
---|
367 | CALL getin('zpicinp',zpicinp) |
---|
368 | !Config key = nudge_tsoil |
---|
369 | !Config Desc = activation of soil temperature nudging |
---|
370 | !Config Def = .FALSE. |
---|
371 | !Config Help = ... |
---|
372 | |
---|
373 | nudge_tsoil=.FALSE. |
---|
374 | CALL getin('nudge_tsoil',nudge_tsoil) |
---|
375 | |
---|
376 | !Config key = isoil_nudge |
---|
377 | !Config Desc = level number where soil temperature is nudged |
---|
378 | !Config Def = 3 |
---|
379 | !Config Help = ... |
---|
380 | |
---|
381 | isoil_nudge=3 |
---|
382 | CALL getin('isoil_nudge',isoil_nudge) |
---|
383 | |
---|
384 | !Config key = Tsoil_nudge |
---|
385 | !Config Desc = target temperature for tsoil(isoil_nudge) |
---|
386 | !Config Def = 300. |
---|
387 | !Config Help = ... |
---|
388 | |
---|
389 | Tsoil_nudge=300. |
---|
390 | CALL getin('Tsoil_nudge',Tsoil_nudge) |
---|
391 | |
---|
392 | !Config key = tau_soil_nudge |
---|
393 | !Config Desc = nudging relaxation time for tsoil |
---|
394 | !Config Def = 3600. |
---|
395 | !Config Help = ... |
---|
396 | |
---|
397 | tau_soil_nudge=3600. |
---|
398 | CALL getin('tau_soil_nudge',tau_soil_nudge) |
---|
399 | |
---|
400 | !---------------------------------------------------------- |
---|
401 | ! Param??tres de for??age pour les forcages communs: |
---|
402 | ! Pour les forcages communs: ces entiers valent 0 ou 1 |
---|
403 | ! tadv= advection tempe, tadvv= adv tempe verticale, tadvh= adv tempe horizontale |
---|
404 | ! qadv= advection q, qadvv= adv q verticale, qadvh= adv q horizontale |
---|
405 | ! trad= 0 (rayonnement actif) ou 1 (prescrit par tend_rad) ou adv (prescir et contenu dans les tadv) |
---|
406 | ! forcages en omega, w, vent geostrophique ou ustar |
---|
407 | ! Parametres de nudging en u,v,t,q valent 0 ou 1 ou le temps de nudging |
---|
408 | !---------------------------------------------------------- |
---|
409 | |
---|
410 | !Config Key = tadv |
---|
411 | !Config Desc = forcage ou non par advection totale de T |
---|
412 | !Config Def = false |
---|
413 | !Config Help = forcage ou non par advection totale de T |
---|
414 | tadv =0 |
---|
415 | CALL getin('tadv',tadv) |
---|
416 | |
---|
417 | !Config Key = tadvv |
---|
418 | !Config Desc = forcage ou non par advection verticale de T |
---|
419 | !Config Def = false |
---|
420 | !Config Help = forcage ou non par advection verticale de T |
---|
421 | tadvv =0 |
---|
422 | CALL getin('tadvv',tadvv) |
---|
423 | |
---|
424 | !Config Key = tadvh |
---|
425 | !Config Desc = forcage ou non par advection horizontale de T |
---|
426 | !Config Def = false |
---|
427 | !Config Help = forcage ou non par advection horizontale de T |
---|
428 | tadvh =0 |
---|
429 | CALL getin('tadvh',tadvh) |
---|
430 | |
---|
431 | !Config Key = thadv |
---|
432 | !Config Desc = forcage ou non par advection totale de Theta |
---|
433 | !Config Def = false |
---|
434 | !Config Help = forcage ou non par advection totale de Theta |
---|
435 | thadv =0 |
---|
436 | CALL getin('thadv',thadv) |
---|
437 | |
---|
438 | !Config Key = thadvv |
---|
439 | !Config Desc = forcage ou non par advection verticale de Theta |
---|
440 | !Config Def = false |
---|
441 | !Config Help = forcage ou non par advection verticale de Theta |
---|
442 | thadvv =0 |
---|
443 | CALL getin('thadvv',thadvv) |
---|
444 | |
---|
445 | !Config Key = thadvh |
---|
446 | !Config Desc = forcage ou non par advection horizontale de Theta |
---|
447 | !Config Def = false |
---|
448 | !Config Help = forcage ou non par advection horizontale de Theta |
---|
449 | thadvh =0 |
---|
450 | CALL getin('thadvh',thadvh) |
---|
451 | |
---|
452 | !Config Key = qadv |
---|
453 | !Config Desc = forcage ou non par advection totale de Q |
---|
454 | !Config Def = false |
---|
455 | !Config Help = forcage ou non par advection totale de Q |
---|
456 | qadv =0 |
---|
457 | CALL getin('qadv',qadv) |
---|
458 | |
---|
459 | !Config Key = qadvv |
---|
460 | !Config Desc = forcage ou non par advection verticale de Q |
---|
461 | !Config Def = false |
---|
462 | !Config Help = forcage ou non par advection verticale de Q |
---|
463 | qadvv =0 |
---|
464 | CALL getin('qadvv',qadvv) |
---|
465 | |
---|
466 | !Config Key = qadvh |
---|
467 | !Config Desc = forcage ou non par advection horizontale de Q |
---|
468 | !Config Def = false |
---|
469 | !Config Help = forcage ou non par advection horizontale de Q |
---|
470 | qadvh =0 |
---|
471 | CALL getin('qadvh',qadvh) |
---|
472 | |
---|
473 | !Config Key = trad |
---|
474 | !Config Desc = forcage ou non par tendance radiative |
---|
475 | !Config Def = false |
---|
476 | !Config Help = forcage ou non par tendance radiative |
---|
477 | trad =0 |
---|
478 | CALL getin('trad',trad) |
---|
479 | |
---|
480 | !Config Key = forc_omega |
---|
481 | !Config Desc = forcage ou non par omega |
---|
482 | !Config Def = false |
---|
483 | !Config Help = forcage ou non par omega |
---|
484 | forc_omega =0 |
---|
485 | CALL getin('forc_omega',forc_omega) |
---|
486 | |
---|
487 | !Config Key = forc_u |
---|
488 | !Config Desc = forcage ou non par u |
---|
489 | !Config Def = false |
---|
490 | !Config Help = forcage ou non par u |
---|
491 | forc_u =0 |
---|
492 | CALL getin('forc_u',forc_u) |
---|
493 | |
---|
494 | !Config Key = forc_v |
---|
495 | !Config Desc = forcage ou non par v |
---|
496 | !Config Def = false |
---|
497 | !Config Help = forcage ou non par v |
---|
498 | forc_v =0 |
---|
499 | CALL getin('forc_v',forc_v) |
---|
500 | !Config Key = forc_w |
---|
501 | !Config Desc = forcage ou non par w |
---|
502 | !Config Def = false |
---|
503 | !Config Help = forcage ou non par w |
---|
504 | forc_w =0 |
---|
505 | CALL getin('forc_w',forc_w) |
---|
506 | |
---|
507 | !Config Key = forc_geo |
---|
508 | !Config Desc = forcage ou non par geo |
---|
509 | !Config Def = false |
---|
510 | !Config Help = forcage ou non par geo |
---|
511 | forc_geo =0 |
---|
512 | CALL getin('forc_geo',forc_geo) |
---|
513 | |
---|
514 | ! Meme chose que ok_precr_ust |
---|
515 | !Config Key = forc_ustar |
---|
516 | !Config Desc = forcage ou non par ustar |
---|
517 | !Config Def = false |
---|
518 | !Config Help = forcage ou non par ustar |
---|
519 | forc_ustar =0 |
---|
520 | CALL getin('forc_ustar',forc_ustar) |
---|
521 | IF (forc_ustar .EQ. 1) ok_prescr_ust=.true. |
---|
522 | |
---|
523 | !Config Key = nudging_u |
---|
524 | !Config Desc = forcage ou non par nudging sur u |
---|
525 | !Config Def = false |
---|
526 | !Config Help = forcage ou non par nudging sur u |
---|
527 | nudging_u =0 |
---|
528 | CALL getin('nudging_u',nudging_u) |
---|
529 | |
---|
530 | !Config Key = nudging_v |
---|
531 | !Config Desc = forcage ou non par nudging sur v |
---|
532 | !Config Def = false |
---|
533 | !Config Help = forcage ou non par nudging sur v |
---|
534 | nudging_v =0 |
---|
535 | CALL getin('nudging_v',nudging_v) |
---|
536 | |
---|
537 | !Config Key = nudging_w |
---|
538 | !Config Desc = forcage ou non par nudging sur w |
---|
539 | !Config Def = false |
---|
540 | !Config Help = forcage ou non par nudging sur w |
---|
541 | nudging_w =0 |
---|
542 | CALL getin('nudging_w',nudging_w) |
---|
543 | |
---|
544 | ! RELIQUE ANCIENS FORMAT. ECRASE PAR LE SUIVANT |
---|
545 | !Config Key = nudging_q |
---|
546 | !Config Desc = forcage ou non par nudging sur q |
---|
547 | !Config Def = false |
---|
548 | !Config Help = forcage ou non par nudging sur q |
---|
549 | nudging_qv =0 |
---|
550 | CALL getin('nudging_q',nudging_qv) |
---|
551 | CALL getin('nudging_qv',nudging_qv) |
---|
552 | |
---|
553 | p_nudging_u=11000. |
---|
554 | p_nudging_v=11000. |
---|
555 | p_nudging_t=11000. |
---|
556 | p_nudging_qv=11000. |
---|
557 | CALL getin('p_nudging_u',p_nudging_u) |
---|
558 | CALL getin('p_nudging_v',p_nudging_v) |
---|
559 | CALL getin('p_nudging_t',p_nudging_t) |
---|
560 | CALL getin('p_nudging_qv',p_nudging_qv) |
---|
561 | |
---|
562 | !Config Key = nudging_t |
---|
563 | !Config Desc = forcage ou non par nudging sur t |
---|
564 | !Config Def = false |
---|
565 | !Config Help = forcage ou non par nudging sur t |
---|
566 | nudging_t =0 |
---|
567 | CALL getin('nudging_t',nudging_t) |
---|
568 | |
---|
569 | |
---|
570 | |
---|
571 | write(lunout,*)' +++++++++++++++++++++++++++++++++++++++' |
---|
572 | write(lunout,*)' Configuration des parametres du gcm1D: ' |
---|
573 | write(lunout,*)' +++++++++++++++++++++++++++++++++++++++' |
---|
574 | write(lunout,*)' restart = ', restart |
---|
575 | write(lunout,*)' forcing_type = ', forcing_type |
---|
576 | write(lunout,*)' time_ini = ', time_ini |
---|
577 | write(lunout,*)' rlat = ', xlat |
---|
578 | write(lunout,*)' rlon = ', xlon |
---|
579 | write(lunout,*)' airephy = ', airefi |
---|
580 | write(lunout,*)' nat_surf = ', nat_surf |
---|
581 | write(lunout,*)' tsurf = ', tsurf |
---|
582 | write(lunout,*)' psurf = ', psurf |
---|
583 | write(lunout,*)' zsurf = ', zsurf |
---|
584 | write(lunout,*)' rugos = ', rugos |
---|
585 | write(lunout,*)' snowmass=', snowmass |
---|
586 | write(lunout,*)' wtsurf = ', wtsurf |
---|
587 | write(lunout,*)' wqsurf = ', wqsurf |
---|
588 | write(lunout,*)' albedo = ', albedo |
---|
589 | write(lunout,*)' xagesno = ', xagesno |
---|
590 | write(lunout,*)' restart_runoff = ', restart_runoff |
---|
591 | write(lunout,*)' qsolinp = ', qsolinp |
---|
592 | write(lunout,*)' zpicinp = ', zpicinp |
---|
593 | write(lunout,*)' nudge_tsoil = ', nudge_tsoil |
---|
594 | write(lunout,*)' isoil_nudge = ', isoil_nudge |
---|
595 | write(lunout,*)' Tsoil_nudge = ', Tsoil_nudge |
---|
596 | write(lunout,*)' tau_soil_nudge = ', tau_soil_nudge |
---|
597 | write(lunout,*)' tadv = ', tadv |
---|
598 | write(lunout,*)' tadvv = ', tadvv |
---|
599 | write(lunout,*)' tadvh = ', tadvh |
---|
600 | write(lunout,*)' thadv = ', thadv |
---|
601 | write(lunout,*)' thadvv = ', thadvv |
---|
602 | write(lunout,*)' thadvh = ', thadvh |
---|
603 | write(lunout,*)' qadv = ', qadv |
---|
604 | write(lunout,*)' qadvv = ', qadvv |
---|
605 | write(lunout,*)' qadvh = ', qadvh |
---|
606 | write(lunout,*)' trad = ', trad |
---|
607 | write(lunout,*)' forc_omega = ', forc_omega |
---|
608 | write(lunout,*)' forc_w = ', forc_w |
---|
609 | write(lunout,*)' forc_geo = ', forc_geo |
---|
610 | write(lunout,*)' forc_ustar = ', forc_ustar |
---|
611 | write(lunout,*)' nudging_u = ', nudging_u |
---|
612 | write(lunout,*)' nudging_v = ', nudging_v |
---|
613 | write(lunout,*)' nudging_t = ', nudging_t |
---|
614 | write(lunout,*)' nudging_qv = ', nudging_qv |
---|
615 | IF (forcing_type .eq.40) THEN |
---|
616 | write(lunout,*) '--- Forcing type GCSS Old --- with:' |
---|
617 | write(lunout,*)'imp_fcg',imp_fcg_gcssold |
---|
618 | write(lunout,*)'ts_fcg',ts_fcg_gcssold |
---|
619 | write(lunout,*)'tp_fcg',Tp_fcg_gcssold |
---|
620 | write(lunout,*)'tp_ini',Tp_ini_gcssold |
---|
621 | write(lunout,*)'xturb_fcg',xTurb_fcg_gcssold |
---|
622 | ENDIF |
---|
623 | |
---|
624 | write(lunout,*)' +++++++++++++++++++++++++++++++++++++++' |
---|
625 | write(lunout,*) |
---|
626 | ! |
---|
627 | RETURN |
---|
628 | END |
---|
629 | ! |
---|
630 | ! $Id: dyn1deta0.F 1279 2010/07/30 A Lahellec$ |
---|
631 | ! |
---|
632 | ! |
---|
633 | SUBROUTINE dyn1deta0(fichnom,plev,play,phi,phis,presnivs, & |
---|
634 | & ucov,vcov,temp,q,omega2) |
---|
635 | USE dimphy |
---|
636 | USE mod_grid_phy_lmdz |
---|
637 | USE mod_phys_lmdz_para |
---|
638 | USE iophy |
---|
639 | USE phys_state_var_mod |
---|
640 | USE iostart |
---|
641 | USE write_field_phy |
---|
642 | USE infotrac |
---|
643 | use control_mod |
---|
644 | USE comconst_mod, ONLY: im, jm, lllm |
---|
645 | USE logic_mod, ONLY: fxyhypb, ysinus |
---|
646 | USE temps_mod, ONLY: annee_ref, day_ini, day_ref, itau_dyn |
---|
647 | |
---|
648 | IMPLICIT NONE |
---|
649 | !======================================================= |
---|
650 | ! Ecriture du fichier de redemarrage sous format NetCDF |
---|
651 | !======================================================= |
---|
652 | ! Declarations: |
---|
653 | ! ------------- |
---|
654 | include "dimensions.h" |
---|
655 | !!#include "control.h" |
---|
656 | include "netcdf.inc" |
---|
657 | |
---|
658 | ! Arguments: |
---|
659 | ! ---------- |
---|
660 | CHARACTER*(*) fichnom |
---|
661 | !Al1 plev tronque pour .nc mais plev(klev+1):=0 |
---|
662 | real :: plev(klon,klev+1),play (klon,klev),phi(klon,klev) |
---|
663 | real :: presnivs(klon,klev) |
---|
664 | real :: ucov(klon,klev),vcov(klon,klev),temp(klon,klev) |
---|
665 | real :: q(klon,klev,nqtot),omega2(klon,klev) |
---|
666 | ! real :: ug(klev),vg(klev),fcoriolis |
---|
667 | real :: phis(klon) |
---|
668 | |
---|
669 | ! Variables locales pour NetCDF: |
---|
670 | ! ------------------------------ |
---|
671 | INTEGER iq |
---|
672 | INTEGER length |
---|
673 | PARAMETER (length = 100) |
---|
674 | REAL tab_cntrl(length) ! tableau des parametres du run |
---|
675 | character*4 nmq(nqtot) |
---|
676 | character*12 modname |
---|
677 | character*80 abort_message |
---|
678 | LOGICAL found |
---|
679 | |
---|
680 | modname = 'dyn1deta0 : ' |
---|
681 | !! nmq(1)="vap" |
---|
682 | !! nmq(2)="cond" |
---|
683 | !! do iq=3,nqtot |
---|
684 | !! write(nmq(iq),'("tra",i1)') iq-2 |
---|
685 | !! enddo |
---|
686 | DO iq = 1,nqtot |
---|
687 | nmq(iq) = trim(tname(iq)) |
---|
688 | ENDDO |
---|
689 | print*,'in dyn1deta0 ',fichnom,klon,klev,nqtot |
---|
690 | CALL open_startphy(fichnom) |
---|
691 | print*,'after open startphy ',fichnom,nmq |
---|
692 | |
---|
693 | ! |
---|
694 | ! Lecture des parametres de controle: |
---|
695 | ! |
---|
696 | CALL get_var("controle",tab_cntrl) |
---|
697 | |
---|
698 | |
---|
699 | im = tab_cntrl(1) |
---|
700 | jm = tab_cntrl(2) |
---|
701 | lllm = tab_cntrl(3) |
---|
702 | day_ref = tab_cntrl(4) |
---|
703 | annee_ref = tab_cntrl(5) |
---|
704 | ! rad = tab_cntrl(6) |
---|
705 | ! omeg = tab_cntrl(7) |
---|
706 | ! g = tab_cntrl(8) |
---|
707 | ! cpp = tab_cntrl(9) |
---|
708 | ! kappa = tab_cntrl(10) |
---|
709 | ! daysec = tab_cntrl(11) |
---|
710 | ! dtvr = tab_cntrl(12) |
---|
711 | ! etot0 = tab_cntrl(13) |
---|
712 | ! ptot0 = tab_cntrl(14) |
---|
713 | ! ztot0 = tab_cntrl(15) |
---|
714 | ! stot0 = tab_cntrl(16) |
---|
715 | ! ang0 = tab_cntrl(17) |
---|
716 | ! pa = tab_cntrl(18) |
---|
717 | ! preff = tab_cntrl(19) |
---|
718 | ! |
---|
719 | ! clon = tab_cntrl(20) |
---|
720 | ! clat = tab_cntrl(21) |
---|
721 | ! grossismx = tab_cntrl(22) |
---|
722 | ! grossismy = tab_cntrl(23) |
---|
723 | ! |
---|
724 | IF ( tab_cntrl(24).EQ.1. ) THEN |
---|
725 | fxyhypb =.true. |
---|
726 | ! dzoomx = tab_cntrl(25) |
---|
727 | ! dzoomy = tab_cntrl(26) |
---|
728 | ! taux = tab_cntrl(28) |
---|
729 | ! tauy = tab_cntrl(29) |
---|
730 | ELSE |
---|
731 | fxyhypb = .false. |
---|
732 | ysinus = .false. |
---|
733 | IF( tab_cntrl(27).EQ.1. ) ysinus =.true. |
---|
734 | ENDIF |
---|
735 | |
---|
736 | day_ini = tab_cntrl(30) |
---|
737 | itau_dyn = tab_cntrl(31) |
---|
738 | ! ................................................................. |
---|
739 | ! |
---|
740 | ! |
---|
741 | ! PRINT*,'rad,omeg,g,cpp,kappa',rad,omeg,g,cpp,kappa |
---|
742 | !Al1 |
---|
743 | Print*,'day_ref,annee_ref,day_ini,itau_dyn', & |
---|
744 | & day_ref,annee_ref,day_ini,itau_dyn |
---|
745 | |
---|
746 | ! Lecture des champs |
---|
747 | ! |
---|
748 | CALL get_field("play",play,found) |
---|
749 | IF (.NOT. found) PRINT*, modname//'Le champ <Play> est absent' |
---|
750 | CALL get_field("phi",phi,found) |
---|
751 | IF (.NOT. found) PRINT*, modname//'Le champ <Phi> est absent' |
---|
752 | CALL get_field("phis",phis,found) |
---|
753 | IF (.NOT. found) PRINT*, modname//'Le champ <Phis> est absent' |
---|
754 | CALL get_field("presnivs",presnivs,found) |
---|
755 | IF (.NOT. found) PRINT*, modname//'Le champ <Presnivs> est absent' |
---|
756 | CALL get_field("ucov",ucov,found) |
---|
757 | IF (.NOT. found) PRINT*, modname//'Le champ <ucov> est absent' |
---|
758 | CALL get_field("vcov",vcov,found) |
---|
759 | IF (.NOT. found) PRINT*, modname//'Le champ <vcov> est absent' |
---|
760 | CALL get_field("temp",temp,found) |
---|
761 | IF (.NOT. found) PRINT*, modname//'Le champ <temp> est absent' |
---|
762 | CALL get_field("omega2",omega2,found) |
---|
763 | IF (.NOT. found) PRINT*, modname//'Le champ <omega2> est absent' |
---|
764 | plev(1,klev+1)=0. |
---|
765 | CALL get_field("plev",plev(:,1:klev),found) |
---|
766 | IF (.NOT. found) PRINT*, modname//'Le champ <Plev> est absent' |
---|
767 | |
---|
768 | Do iq=1,nqtot |
---|
769 | CALL get_field("q"//nmq(iq),q(:,:,iq),found) |
---|
770 | IF (.NOT.found)PRINT*, modname//'Le champ <q'//nmq//'> est absent' |
---|
771 | EndDo |
---|
772 | |
---|
773 | CALL close_startphy |
---|
774 | print*,' close startphy',fichnom,play(1,1),play(1,klev),temp(1,klev) |
---|
775 | ! |
---|
776 | RETURN |
---|
777 | END |
---|
778 | ! |
---|
779 | ! $Id: dyn1dredem.F 1279 2010/07/29 A Lahellec$ |
---|
780 | ! |
---|
781 | ! |
---|
782 | SUBROUTINE dyn1dredem(fichnom,plev,play,phi,phis,presnivs, & |
---|
783 | & ucov,vcov,temp,q,omega2) |
---|
784 | USE dimphy |
---|
785 | USE mod_grid_phy_lmdz |
---|
786 | USE mod_phys_lmdz_para |
---|
787 | USE phys_state_var_mod |
---|
788 | USE iostart |
---|
789 | USE infotrac |
---|
790 | use control_mod |
---|
791 | USE comconst_mod, ONLY: cpp, daysec, dtvr, g, kappa, omeg, rad |
---|
792 | USE logic_mod, ONLY: fxyhypb, ysinus |
---|
793 | USE temps_mod, ONLY: annee_ref,day_end,day_ref,itau_dyn,itaufin |
---|
794 | |
---|
795 | IMPLICIT NONE |
---|
796 | !======================================================= |
---|
797 | ! Ecriture du fichier de redemarrage sous format NetCDF |
---|
798 | !======================================================= |
---|
799 | ! Declarations: |
---|
800 | ! ------------- |
---|
801 | include "dimensions.h" |
---|
802 | !!#include "control.h" |
---|
803 | include "netcdf.inc" |
---|
804 | |
---|
805 | ! Arguments: |
---|
806 | ! ---------- |
---|
807 | CHARACTER*(*) fichnom |
---|
808 | !Al1 plev tronque pour .nc mais plev(klev+1):=0 |
---|
809 | real :: plev(klon,klev),play (klon,klev),phi(klon,klev) |
---|
810 | real :: presnivs(klon,klev) |
---|
811 | real :: ucov(klon,klev),vcov(klon,klev),temp(klon,klev) |
---|
812 | real :: q(klon,klev,nqtot) |
---|
813 | real :: omega2(klon,klev),rho(klon,klev+1) |
---|
814 | ! real :: ug(klev),vg(klev),fcoriolis |
---|
815 | real :: phis(klon) |
---|
816 | |
---|
817 | ! Variables locales pour NetCDF: |
---|
818 | ! ------------------------------ |
---|
819 | INTEGER nid |
---|
820 | INTEGER ierr |
---|
821 | INTEGER iq,l |
---|
822 | INTEGER length |
---|
823 | PARAMETER (length = 100) |
---|
824 | REAL tab_cntrl(length) ! tableau des parametres du run |
---|
825 | character*4 nmq(nqtot) |
---|
826 | character*20 modname |
---|
827 | character*80 abort_message |
---|
828 | ! |
---|
829 | INTEGER pass |
---|
830 | |
---|
831 | CALL open_restartphy(fichnom) |
---|
832 | print*,'redm1 ',fichnom,klon,klev,nqtot |
---|
833 | !! nmq(1)="vap" |
---|
834 | !! nmq(2)="cond" |
---|
835 | !! nmq(3)="tra1" |
---|
836 | !! nmq(4)="tra2" |
---|
837 | DO iq = 1,nqtot |
---|
838 | nmq(iq) = trim(tname(iq)) |
---|
839 | ENDDO |
---|
840 | |
---|
841 | ! modname = 'dyn1dredem' |
---|
842 | ! ierr = NF_OPEN(fichnom, NF_WRITE, nid) |
---|
843 | ! IF (ierr .NE. NF_NOERR) THEN |
---|
844 | ! abort_message="Pb. d ouverture "//fichnom |
---|
845 | ! CALL abort_gcm('Modele 1D',abort_message,1) |
---|
846 | ! ENDIF |
---|
847 | |
---|
848 | DO l=1,length |
---|
849 | tab_cntrl(l) = 0. |
---|
850 | ENDDO |
---|
851 | tab_cntrl(1) = FLOAT(iim) |
---|
852 | tab_cntrl(2) = FLOAT(jjm) |
---|
853 | tab_cntrl(3) = FLOAT(llm) |
---|
854 | tab_cntrl(4) = FLOAT(day_ref) |
---|
855 | tab_cntrl(5) = FLOAT(annee_ref) |
---|
856 | tab_cntrl(6) = rad |
---|
857 | tab_cntrl(7) = omeg |
---|
858 | tab_cntrl(8) = g |
---|
859 | tab_cntrl(9) = cpp |
---|
860 | tab_cntrl(10) = kappa |
---|
861 | tab_cntrl(11) = daysec |
---|
862 | tab_cntrl(12) = dtvr |
---|
863 | ! tab_cntrl(13) = etot0 |
---|
864 | ! tab_cntrl(14) = ptot0 |
---|
865 | ! tab_cntrl(15) = ztot0 |
---|
866 | ! tab_cntrl(16) = stot0 |
---|
867 | ! tab_cntrl(17) = ang0 |
---|
868 | ! tab_cntrl(18) = pa |
---|
869 | ! tab_cntrl(19) = preff |
---|
870 | ! |
---|
871 | ! ..... parametres pour le zoom ...... |
---|
872 | |
---|
873 | ! tab_cntrl(20) = clon |
---|
874 | ! tab_cntrl(21) = clat |
---|
875 | ! tab_cntrl(22) = grossismx |
---|
876 | ! tab_cntrl(23) = grossismy |
---|
877 | ! |
---|
878 | IF ( fxyhypb ) THEN |
---|
879 | tab_cntrl(24) = 1. |
---|
880 | ! tab_cntrl(25) = dzoomx |
---|
881 | ! tab_cntrl(26) = dzoomy |
---|
882 | tab_cntrl(27) = 0. |
---|
883 | ! tab_cntrl(28) = taux |
---|
884 | ! tab_cntrl(29) = tauy |
---|
885 | ELSE |
---|
886 | tab_cntrl(24) = 0. |
---|
887 | ! tab_cntrl(25) = dzoomx |
---|
888 | ! tab_cntrl(26) = dzoomy |
---|
889 | tab_cntrl(27) = 0. |
---|
890 | tab_cntrl(28) = 0. |
---|
891 | tab_cntrl(29) = 0. |
---|
892 | IF( ysinus ) tab_cntrl(27) = 1. |
---|
893 | ENDIF |
---|
894 | !Al1 iday_end -> day_end |
---|
895 | tab_cntrl(30) = FLOAT(day_end) |
---|
896 | tab_cntrl(31) = FLOAT(itau_dyn + itaufin) |
---|
897 | ! |
---|
898 | DO pass=1,2 |
---|
899 | CALL put_var(pass,"controle","Param. de controle Dyn1D",tab_cntrl) |
---|
900 | ! |
---|
901 | |
---|
902 | ! Ecriture/extension de la coordonnee temps |
---|
903 | |
---|
904 | |
---|
905 | ! Ecriture des champs |
---|
906 | ! |
---|
907 | CALL put_field(pass,"plev","p interfaces sauf la nulle",plev) |
---|
908 | CALL put_field(pass,"play","",play) |
---|
909 | CALL put_field(pass,"phi","geopotentielle",phi) |
---|
910 | CALL put_field(pass,"phis","geopotentiell de surface",phis) |
---|
911 | CALL put_field(pass,"presnivs","",presnivs) |
---|
912 | CALL put_field(pass,"ucov","",ucov) |
---|
913 | CALL put_field(pass,"vcov","",vcov) |
---|
914 | CALL put_field(pass,"temp","",temp) |
---|
915 | CALL put_field(pass,"omega2","",omega2) |
---|
916 | |
---|
917 | Do iq=1,nqtot |
---|
918 | CALL put_field(pass,"q"//nmq(iq),"eau vap ou condens et traceurs", & |
---|
919 | & q(:,:,iq)) |
---|
920 | EndDo |
---|
921 | IF (pass==1) CALL enddef_restartphy |
---|
922 | IF (pass==2) CALL close_restartphy |
---|
923 | |
---|
924 | |
---|
925 | ENDDO |
---|
926 | |
---|
927 | ! |
---|
928 | RETURN |
---|
929 | END |
---|
930 | SUBROUTINE gr_fi_dyn(nfield,ngrid,im,jm,pfi,pdyn) |
---|
931 | IMPLICIT NONE |
---|
932 | !======================================================================= |
---|
933 | ! passage d'un champ de la grille scalaire a la grille physique |
---|
934 | !======================================================================= |
---|
935 | |
---|
936 | !----------------------------------------------------------------------- |
---|
937 | ! declarations: |
---|
938 | ! ------------- |
---|
939 | |
---|
940 | INTEGER im,jm,ngrid,nfield |
---|
941 | REAL pdyn(im,jm,nfield) |
---|
942 | REAL pfi(ngrid,nfield) |
---|
943 | |
---|
944 | INTEGER i,j,ifield,ig |
---|
945 | |
---|
946 | !----------------------------------------------------------------------- |
---|
947 | ! calcul: |
---|
948 | ! ------- |
---|
949 | |
---|
950 | DO ifield=1,nfield |
---|
951 | ! traitement des poles |
---|
952 | DO i=1,im |
---|
953 | pdyn(i,1,ifield)=pfi(1,ifield) |
---|
954 | pdyn(i,jm,ifield)=pfi(ngrid,ifield) |
---|
955 | ENDDO |
---|
956 | |
---|
957 | ! traitement des point normaux |
---|
958 | DO j=2,jm-1 |
---|
959 | ig=2+(j-2)*(im-1) |
---|
960 | CALL SCOPY(im-1,pfi(ig,ifield),1,pdyn(1,j,ifield),1) |
---|
961 | pdyn(im,j,ifield)=pdyn(1,j,ifield) |
---|
962 | ENDDO |
---|
963 | ENDDO |
---|
964 | |
---|
965 | RETURN |
---|
966 | END |
---|
967 | |
---|
968 | |
---|
969 | |
---|
970 | SUBROUTINE abort_gcm(modname, message, ierr) |
---|
971 | |
---|
972 | USE IOIPSL |
---|
973 | ! |
---|
974 | ! Stops the simulation cleanly, closing files and printing various |
---|
975 | ! comments |
---|
976 | ! |
---|
977 | ! Input: modname = name of calling program |
---|
978 | ! message = stuff to print |
---|
979 | ! ierr = severity of situation ( = 0 normal ) |
---|
980 | |
---|
981 | character(len=*) modname |
---|
982 | integer ierr |
---|
983 | character(len=*) message |
---|
984 | |
---|
985 | write(*,*) 'in abort_gcm' |
---|
986 | call histclo |
---|
987 | ! call histclo(2) |
---|
988 | ! call histclo(3) |
---|
989 | ! call histclo(4) |
---|
990 | ! call histclo(5) |
---|
991 | write(*,*) 'out of histclo' |
---|
992 | write(*,*) 'Stopping in ', modname |
---|
993 | write(*,*) 'Reason = ',message |
---|
994 | call getin_dump |
---|
995 | ! |
---|
996 | if (ierr .eq. 0) then |
---|
997 | write(*,*) 'Everything is cool' |
---|
998 | else |
---|
999 | write(*,*) 'Houston, we have a problem ', ierr |
---|
1000 | endif |
---|
1001 | STOP |
---|
1002 | END |
---|
1003 | REAL FUNCTION fq_sat(kelvin, millibar) |
---|
1004 | ! |
---|
1005 | IMPLICIT none |
---|
1006 | !====================================================================== |
---|
1007 | ! Autheur(s): Z.X. Li (LMD/CNRS) |
---|
1008 | ! Objet: calculer la vapeur d'eau saturante (formule Centre Euro.) |
---|
1009 | !====================================================================== |
---|
1010 | ! Arguments: |
---|
1011 | ! kelvin---input-R: temperature en Kelvin |
---|
1012 | ! millibar--input-R: pression en mb |
---|
1013 | ! |
---|
1014 | ! fq_sat----output-R: vapeur d'eau saturante en kg/kg |
---|
1015 | !====================================================================== |
---|
1016 | ! |
---|
1017 | REAL kelvin, millibar |
---|
1018 | ! |
---|
1019 | REAL r2es |
---|
1020 | PARAMETER (r2es=611.14 *18.0153/28.9644) |
---|
1021 | ! |
---|
1022 | REAL r3les, r3ies, r3es |
---|
1023 | PARAMETER (R3LES=17.269) |
---|
1024 | PARAMETER (R3IES=21.875) |
---|
1025 | ! |
---|
1026 | REAL r4les, r4ies, r4es |
---|
1027 | PARAMETER (R4LES=35.86) |
---|
1028 | PARAMETER (R4IES=7.66) |
---|
1029 | ! |
---|
1030 | REAL rtt |
---|
1031 | PARAMETER (rtt=273.16) |
---|
1032 | ! |
---|
1033 | REAL retv |
---|
1034 | PARAMETER (retv=28.9644/18.0153 - 1.0) |
---|
1035 | ! |
---|
1036 | REAL zqsat |
---|
1037 | REAL temp, pres |
---|
1038 | ! ------------------------------------------------------------------ |
---|
1039 | ! |
---|
1040 | ! |
---|
1041 | temp = kelvin |
---|
1042 | pres = millibar * 100.0 |
---|
1043 | ! write(*,*)'kelvin,millibar=',kelvin,millibar |
---|
1044 | ! write(*,*)'temp,pres=',temp,pres |
---|
1045 | ! |
---|
1046 | IF (temp .LE. rtt) THEN |
---|
1047 | r3es = r3ies |
---|
1048 | r4es = r4ies |
---|
1049 | ELSE |
---|
1050 | r3es = r3les |
---|
1051 | r4es = r4les |
---|
1052 | ENDIF |
---|
1053 | ! |
---|
1054 | zqsat=r2es/pres * EXP ( r3es*(temp-rtt) / (temp-r4es) ) |
---|
1055 | zqsat=MIN(0.5,ZQSAT) |
---|
1056 | zqsat=zqsat/(1.-retv *zqsat) |
---|
1057 | ! |
---|
1058 | fq_sat = zqsat |
---|
1059 | ! |
---|
1060 | RETURN |
---|
1061 | END |
---|
1062 | |
---|
1063 | SUBROUTINE gr_dyn_fi(nfield,im,jm,ngrid,pdyn,pfi) |
---|
1064 | IMPLICIT NONE |
---|
1065 | !======================================================================= |
---|
1066 | ! passage d'un champ de la grille scalaire a la grille physique |
---|
1067 | !======================================================================= |
---|
1068 | |
---|
1069 | !----------------------------------------------------------------------- |
---|
1070 | ! declarations: |
---|
1071 | ! ------------- |
---|
1072 | |
---|
1073 | INTEGER im,jm,ngrid,nfield |
---|
1074 | REAL pdyn(im,jm,nfield) |
---|
1075 | REAL pfi(ngrid,nfield) |
---|
1076 | |
---|
1077 | INTEGER j,ifield,ig |
---|
1078 | |
---|
1079 | !----------------------------------------------------------------------- |
---|
1080 | ! calcul: |
---|
1081 | ! ------- |
---|
1082 | |
---|
1083 | IF(ngrid.NE.2+(jm-2)*(im-1).AND.ngrid.NE.1) & |
---|
1084 | & STOP 'probleme de dim' |
---|
1085 | ! traitement des poles |
---|
1086 | CALL SCOPY(nfield,pdyn,im*jm,pfi,ngrid) |
---|
1087 | CALL SCOPY(nfield,pdyn(1,jm,1),im*jm,pfi(ngrid,1),ngrid) |
---|
1088 | |
---|
1089 | ! traitement des point normaux |
---|
1090 | DO ifield=1,nfield |
---|
1091 | DO j=2,jm-1 |
---|
1092 | ig=2+(j-2)*(im-1) |
---|
1093 | CALL SCOPY(im-1,pdyn(1,j,ifield),1,pfi(ig,ifield),1) |
---|
1094 | ENDDO |
---|
1095 | ENDDO |
---|
1096 | |
---|
1097 | RETURN |
---|
1098 | END |
---|
1099 | |
---|
1100 | SUBROUTINE disvert0(pa,preff,ap,bp,dpres,presnivs,nivsigs,nivsig) |
---|
1101 | |
---|
1102 | ! Ancienne version disvert dont on a modifie nom pour utiliser |
---|
1103 | ! le disvert de dyn3d (qui permet d'utiliser grille avec ab,bp imposes) |
---|
1104 | ! (MPL 18092012) |
---|
1105 | ! |
---|
1106 | ! Auteur : P. Le Van . |
---|
1107 | ! |
---|
1108 | IMPLICIT NONE |
---|
1109 | |
---|
1110 | include "dimensions.h" |
---|
1111 | include "paramet.h" |
---|
1112 | ! |
---|
1113 | !======================================================================= |
---|
1114 | ! |
---|
1115 | ! |
---|
1116 | ! s = sigma ** kappa : coordonnee verticale |
---|
1117 | ! dsig(l) : epaisseur de la couche l ds la coord. s |
---|
1118 | ! sig(l) : sigma a l'interface des couches l et l-1 |
---|
1119 | ! ds(l) : distance entre les couches l et l-1 en coord.s |
---|
1120 | ! |
---|
1121 | !======================================================================= |
---|
1122 | ! |
---|
1123 | REAL pa,preff |
---|
1124 | REAL ap(llmp1),bp(llmp1),dpres(llm),nivsigs(llm),nivsig(llmp1) |
---|
1125 | REAL presnivs(llm) |
---|
1126 | ! |
---|
1127 | ! declarations: |
---|
1128 | ! ------------- |
---|
1129 | ! |
---|
1130 | REAL sig(llm+1),dsig(llm) |
---|
1131 | ! |
---|
1132 | INTEGER l |
---|
1133 | REAL snorm |
---|
1134 | REAL alpha,beta,gama,delta,deltaz,h |
---|
1135 | INTEGER np,ierr |
---|
1136 | REAL pi,x |
---|
1137 | |
---|
1138 | !----------------------------------------------------------------------- |
---|
1139 | ! |
---|
1140 | pi=2.*ASIN(1.) |
---|
1141 | |
---|
1142 | OPEN(99,file='sigma.def',status='old',form='formatted', & |
---|
1143 | & iostat=ierr) |
---|
1144 | |
---|
1145 | !----------------------------------------------------------------------- |
---|
1146 | ! cas 1 on lit les options dans sigma.def: |
---|
1147 | ! ---------------------------------------- |
---|
1148 | |
---|
1149 | IF (ierr.eq.0) THEN |
---|
1150 | |
---|
1151 | print*,'WARNING!!! on lit les options dans sigma.def' |
---|
1152 | READ(99,*) deltaz |
---|
1153 | READ(99,*) h |
---|
1154 | READ(99,*) beta |
---|
1155 | READ(99,*) gama |
---|
1156 | READ(99,*) delta |
---|
1157 | READ(99,*) np |
---|
1158 | CLOSE(99) |
---|
1159 | alpha=deltaz/(llm*h) |
---|
1160 | ! |
---|
1161 | |
---|
1162 | DO 1 l = 1, llm |
---|
1163 | dsig(l) = (alpha+(1.-alpha)*exp(-beta*(llm-l)))* & |
---|
1164 | & ( (tanh(gama*l)/tanh(gama*llm))**np + & |
---|
1165 | & (1.-l/FLOAT(llm))*delta ) |
---|
1166 | 1 CONTINUE |
---|
1167 | |
---|
1168 | sig(1)=1. |
---|
1169 | DO 101 l=1,llm-1 |
---|
1170 | sig(l+1)=sig(l)*(1.-dsig(l))/(1.+dsig(l)) |
---|
1171 | 101 CONTINUE |
---|
1172 | sig(llm+1)=0. |
---|
1173 | |
---|
1174 | DO 2 l = 1, llm |
---|
1175 | dsig(l) = sig(l)-sig(l+1) |
---|
1176 | 2 CONTINUE |
---|
1177 | ! |
---|
1178 | |
---|
1179 | ELSE |
---|
1180 | !----------------------------------------------------------------------- |
---|
1181 | ! cas 2 ancienne discretisation (LMD5...): |
---|
1182 | ! ---------------------------------------- |
---|
1183 | |
---|
1184 | PRINT*,'WARNING!!! Ancienne discretisation verticale' |
---|
1185 | |
---|
1186 | h=7. |
---|
1187 | snorm = 0. |
---|
1188 | DO l = 1, llm |
---|
1189 | x = 2.*asin(1.) * (FLOAT(l)-0.5) / float(llm+1) |
---|
1190 | dsig(l) = 1.0 + 7.0 * SIN(x)**2 |
---|
1191 | snorm = snorm + dsig(l) |
---|
1192 | ENDDO |
---|
1193 | snorm = 1./snorm |
---|
1194 | DO l = 1, llm |
---|
1195 | dsig(l) = dsig(l)*snorm |
---|
1196 | ENDDO |
---|
1197 | sig(llm+1) = 0. |
---|
1198 | DO l = llm, 1, -1 |
---|
1199 | sig(l) = sig(l+1) + dsig(l) |
---|
1200 | ENDDO |
---|
1201 | |
---|
1202 | ENDIF |
---|
1203 | |
---|
1204 | |
---|
1205 | DO l=1,llm |
---|
1206 | nivsigs(l) = FLOAT(l) |
---|
1207 | ENDDO |
---|
1208 | |
---|
1209 | DO l=1,llmp1 |
---|
1210 | nivsig(l)= FLOAT(l) |
---|
1211 | ENDDO |
---|
1212 | |
---|
1213 | ! |
---|
1214 | ! .... Calculs de ap(l) et de bp(l) .... |
---|
1215 | ! ......................................... |
---|
1216 | ! |
---|
1217 | ! |
---|
1218 | ! ..... pa et preff sont lus sur les fichiers start par lectba ..... |
---|
1219 | ! |
---|
1220 | |
---|
1221 | bp(llmp1) = 0. |
---|
1222 | |
---|
1223 | DO l = 1, llm |
---|
1224 | !c |
---|
1225 | !cc ap(l) = 0. |
---|
1226 | !cc bp(l) = sig(l) |
---|
1227 | |
---|
1228 | bp(l) = EXP( 1. -1./( sig(l)*sig(l)) ) |
---|
1229 | ap(l) = pa * ( sig(l) - bp(l) ) |
---|
1230 | ! |
---|
1231 | ENDDO |
---|
1232 | ap(llmp1) = pa * ( sig(llmp1) - bp(llmp1) ) |
---|
1233 | |
---|
1234 | PRINT *,' BP ' |
---|
1235 | PRINT *, bp |
---|
1236 | PRINT *,' AP ' |
---|
1237 | PRINT *, ap |
---|
1238 | |
---|
1239 | DO l = 1, llm |
---|
1240 | dpres(l) = bp(l) - bp(l+1) |
---|
1241 | presnivs(l) = 0.5 *( ap(l)+bp(l)*preff + ap(l+1)+bp(l+1)*preff ) |
---|
1242 | ENDDO |
---|
1243 | |
---|
1244 | PRINT *,' PRESNIVS ' |
---|
1245 | PRINT *,presnivs |
---|
1246 | |
---|
1247 | RETURN |
---|
1248 | END |
---|
1249 | |
---|
1250 | !====================================================================== |
---|
1251 | SUBROUTINE read_tsurf1d(knon,sst_out) |
---|
1252 | |
---|
1253 | ! This subroutine specifies the surface temperature to be used in 1D simulations |
---|
1254 | |
---|
1255 | USE dimphy, ONLY : klon |
---|
1256 | |
---|
1257 | INTEGER, INTENT(IN) :: knon ! nomber of points on compressed grid |
---|
1258 | REAL, DIMENSION(klon), INTENT(OUT) :: sst_out ! tsurf used to force the single-column model |
---|
1259 | |
---|
1260 | INTEGER :: i |
---|
1261 | ! COMMON defined in lmdz1d.F: |
---|
1262 | real ts_cur |
---|
1263 | common /sst_forcing/ts_cur |
---|
1264 | |
---|
1265 | DO i = 1, knon |
---|
1266 | sst_out(i) = ts_cur |
---|
1267 | ENDDO |
---|
1268 | |
---|
1269 | END SUBROUTINE read_tsurf1d |
---|
1270 | |
---|
1271 | !=============================================================== |
---|
1272 | subroutine advect_vert(llm,w,dt,q,plev) |
---|
1273 | !=============================================================== |
---|
1274 | ! Schema amont pour l'advection verticale en 1D |
---|
1275 | ! w est la vitesse verticale dp/dt en Pa/s |
---|
1276 | ! Traitement en volumes finis |
---|
1277 | ! d / dt ( zm q ) = delta_z ( omega q ) |
---|
1278 | ! d / dt ( zm ) = delta_z ( omega ) |
---|
1279 | ! avec zm = delta_z ( p ) |
---|
1280 | ! si * designe la valeur au pas de temps t+dt |
---|
1281 | ! zm*(l) q*(l) - zm(l) q(l) = w(l+1) q(l+1) - w(l) q(l) |
---|
1282 | ! zm*(l) -zm(l) = w(l+1) - w(l) |
---|
1283 | ! avec w=omega * dt |
---|
1284 | !--------------------------------------------------------------- |
---|
1285 | implicit none |
---|
1286 | ! arguments |
---|
1287 | integer llm |
---|
1288 | real w(llm+1),q(llm),plev(llm+1),dt |
---|
1289 | |
---|
1290 | ! local |
---|
1291 | integer l |
---|
1292 | real zwq(llm+1),zm(llm+1),zw(llm+1) |
---|
1293 | real qold |
---|
1294 | |
---|
1295 | !--------------------------------------------------------------- |
---|
1296 | |
---|
1297 | do l=1,llm |
---|
1298 | zw(l)=dt*w(l) |
---|
1299 | zm(l)=plev(l)-plev(l+1) |
---|
1300 | zwq(l)=q(l)*zw(l) |
---|
1301 | enddo |
---|
1302 | zwq(llm+1)=0. |
---|
1303 | zw(llm+1)=0. |
---|
1304 | |
---|
1305 | do l=1,llm |
---|
1306 | qold=q(l) |
---|
1307 | q(l)=(q(l)*zm(l)+zwq(l+1)-zwq(l))/(zm(l)+zw(l+1)-zw(l)) |
---|
1308 | print*,'ADV Q ',zm(l),zw(l),zwq(l),qold,q(l) |
---|
1309 | enddo |
---|
1310 | |
---|
1311 | |
---|
1312 | return |
---|
1313 | end |
---|
1314 | |
---|
1315 | !=============================================================== |
---|
1316 | |
---|
1317 | |
---|
1318 | SUBROUTINE advect_va(llm,omega,d_t_va,d_q_va,d_u_va,d_v_va, & |
---|
1319 | & q,temp,u,v,play) |
---|
1320 | !itlmd |
---|
1321 | !---------------------------------------------------------------------- |
---|
1322 | ! Calcul de l'advection verticale (ascendance et subsidence) de |
---|
1323 | ! temperature et d'humidite. Hypothese : ce qui rentre de l'exterieur |
---|
1324 | ! a les memes caracteristiques que l'air de la colonne 1D (WTG) ou |
---|
1325 | ! sans WTG rajouter une advection horizontale |
---|
1326 | !---------------------------------------------------------------------- |
---|
1327 | implicit none |
---|
1328 | #include "YOMCST.h" |
---|
1329 | ! argument |
---|
1330 | integer llm |
---|
1331 | real omega(llm+1),d_t_va(llm), d_q_va(llm,3) |
---|
1332 | real d_u_va(llm), d_v_va(llm) |
---|
1333 | real q(llm,3),temp(llm) |
---|
1334 | real u(llm),v(llm) |
---|
1335 | real play(llm) |
---|
1336 | ! interne |
---|
1337 | integer l |
---|
1338 | real alpha,omgdown,omgup |
---|
1339 | |
---|
1340 | do l= 1,llm |
---|
1341 | if(l.eq.1) then |
---|
1342 | !si omgup pour la couche 1, alors tendance nulle |
---|
1343 | omgdown=max(omega(2),0.0) |
---|
1344 | alpha = rkappa*temp(l)*(1.+q(l,1)*rv/rd)/(play(l)*(1.+q(l,1))) |
---|
1345 | d_t_va(l)= alpha*(omgdown)-omgdown*(temp(l)-temp(l+1)) & |
---|
1346 | & /(play(l)-play(l+1)) |
---|
1347 | |
---|
1348 | d_q_va(l,:)= -omgdown*(q(l,:)-q(l+1,:))/(play(l)-play(l+1)) |
---|
1349 | |
---|
1350 | d_u_va(l)= -omgdown*(u(l)-u(l+1))/(play(l)-play(l+1)) |
---|
1351 | d_v_va(l)= -omgdown*(v(l)-v(l+1))/(play(l)-play(l+1)) |
---|
1352 | |
---|
1353 | |
---|
1354 | elseif(l.eq.llm) then |
---|
1355 | omgup=min(omega(l),0.0) |
---|
1356 | alpha = rkappa*temp(l)*(1.+q(l,1)*rv/rd)/(play(l)*(1.+q(l,1))) |
---|
1357 | d_t_va(l)= alpha*(omgup)- & |
---|
1358 | |
---|
1359 | !bug? & omgup*(temp(l-1)-temp(l))/(play(l-1)-plev(l)) |
---|
1360 | & omgup*(temp(l-1)-temp(l))/(play(l-1)-play(l)) |
---|
1361 | d_q_va(l,:)= -omgup*(q(l-1,:)-q(l,:))/(play(l-1)-play(l)) |
---|
1362 | d_u_va(l)= -omgup*(u(l-1)-u(l))/(play(l-1)-play(l)) |
---|
1363 | d_v_va(l)= -omgup*(v(l-1)-v(l))/(play(l-1)-play(l)) |
---|
1364 | |
---|
1365 | else |
---|
1366 | omgup=min(omega(l),0.0) |
---|
1367 | omgdown=max(omega(l+1),0.0) |
---|
1368 | alpha = rkappa*temp(l)*(1.+q(l,1)*rv/rd)/(play(l)*(1.+q(l,1))) |
---|
1369 | d_t_va(l)= alpha*(omgup+omgdown)-omgdown*(temp(l)-temp(l+1)) & |
---|
1370 | & /(play(l)-play(l+1))- & |
---|
1371 | !bug? & omgup*(temp(l-1)-temp(l))/(play(l-1)-plev(l)) |
---|
1372 | & omgup*(temp(l-1)-temp(l))/(play(l-1)-play(l)) |
---|
1373 | ! print*, ' ??? ' |
---|
1374 | |
---|
1375 | d_q_va(l,:)= -omgdown*(q(l,:)-q(l+1,:)) & |
---|
1376 | & /(play(l)-play(l+1))- & |
---|
1377 | & omgup*(q(l-1,:)-q(l,:))/(play(l-1)-play(l)) |
---|
1378 | d_u_va(l)= -omgdown*(u(l)-u(l+1)) & |
---|
1379 | & /(play(l)-play(l+1))- & |
---|
1380 | & omgup*(u(l-1)-u(l))/(play(l-1)-play(l)) |
---|
1381 | d_v_va(l)= -omgdown*(v(l)-v(l+1)) & |
---|
1382 | & /(play(l)-play(l+1))- & |
---|
1383 | & omgup*(v(l-1)-v(l))/(play(l-1)-play(l)) |
---|
1384 | |
---|
1385 | endif |
---|
1386 | |
---|
1387 | enddo |
---|
1388 | !fin itlmd |
---|
1389 | return |
---|
1390 | end |
---|
1391 | ! SUBROUTINE lstendH(llm,omega,d_t_va,d_q_va,d_u_va,d_v_va, |
---|
1392 | SUBROUTINE lstendH(llm,nqtot,omega,d_t_va,d_q_va, & |
---|
1393 | & q,temp,u,v,play) |
---|
1394 | !itlmd |
---|
1395 | !---------------------------------------------------------------------- |
---|
1396 | ! Calcul de l'advection verticale (ascendance et subsidence) de |
---|
1397 | ! temperature et d'humidite. Hypothese : ce qui rentre de l'exterieur |
---|
1398 | ! a les memes caracteristiques que l'air de la colonne 1D (WTG) ou |
---|
1399 | ! sans WTG rajouter une advection horizontale |
---|
1400 | !---------------------------------------------------------------------- |
---|
1401 | implicit none |
---|
1402 | #include "YOMCST.h" |
---|
1403 | ! argument |
---|
1404 | integer llm,nqtot |
---|
1405 | real omega(llm+1),d_t_va(llm), d_q_va(llm,nqtot) |
---|
1406 | ! real d_u_va(llm), d_v_va(llm) |
---|
1407 | real q(llm,nqtot),temp(llm) |
---|
1408 | real u(llm),v(llm) |
---|
1409 | real play(llm) |
---|
1410 | real cor(llm) |
---|
1411 | ! real dph(llm),dudp(llm),dvdp(llm),dqdp(llm),dtdp(llm) |
---|
1412 | real dph(llm),dqdp(llm),dtdp(llm) |
---|
1413 | ! interne |
---|
1414 | integer k |
---|
1415 | real omdn,omup |
---|
1416 | |
---|
1417 | ! dudp=0. |
---|
1418 | ! dvdp=0. |
---|
1419 | dqdp=0. |
---|
1420 | dtdp=0. |
---|
1421 | ! d_u_va=0. |
---|
1422 | ! d_v_va=0. |
---|
1423 | |
---|
1424 | cor(:) = rkappa*temp*(1.+q(:,1)*rv/rd)/(play*(1.+q(:,1))) |
---|
1425 | |
---|
1426 | |
---|
1427 | do k=2,llm-1 |
---|
1428 | |
---|
1429 | dph (k-1) = (play(k )- play(k-1 )) |
---|
1430 | ! dudp (k-1) = (u (k )- u (k-1 ))/dph(k-1) |
---|
1431 | ! dvdp (k-1) = (v (k )- v (k-1 ))/dph(k-1) |
---|
1432 | dqdp (k-1) = (q (k,1)- q (k-1,1))/dph(k-1) |
---|
1433 | dtdp (k-1) = (temp(k )- temp(k-1 ))/dph(k-1) |
---|
1434 | |
---|
1435 | enddo |
---|
1436 | |
---|
1437 | ! dudp ( llm ) = dudp ( llm-1 ) |
---|
1438 | ! dvdp ( llm ) = dvdp ( llm-1 ) |
---|
1439 | dqdp ( llm ) = dqdp ( llm-1 ) |
---|
1440 | dtdp ( llm ) = dtdp ( llm-1 ) |
---|
1441 | |
---|
1442 | do k=2,llm-1 |
---|
1443 | omdn=max(0.0,omega(k+1)) |
---|
1444 | omup=min(0.0,omega( k )) |
---|
1445 | |
---|
1446 | ! d_u_va(k) = -omdn*dudp(k)-omup*dudp(k-1) |
---|
1447 | ! d_v_va(k) = -omdn*dvdp(k)-omup*dvdp(k-1) |
---|
1448 | d_q_va(k,1)= -omdn*dqdp(k)-omup*dqdp(k-1) |
---|
1449 | d_t_va(k) = -omdn*dtdp(k)-omup*dtdp(k-1)+(omup+omdn)*cor(k) |
---|
1450 | enddo |
---|
1451 | |
---|
1452 | omdn=max(0.0,omega( 2 )) |
---|
1453 | omup=min(0.0,omega(llm)) |
---|
1454 | ! d_u_va( 1 ) = -omdn*dudp( 1 ) |
---|
1455 | ! d_u_va(llm) = -omup*dudp(llm) |
---|
1456 | ! d_v_va( 1 ) = -omdn*dvdp( 1 ) |
---|
1457 | ! d_v_va(llm) = -omup*dvdp(llm) |
---|
1458 | d_q_va( 1 ,1) = -omdn*dqdp( 1 ) |
---|
1459 | d_q_va(llm,1) = -omup*dqdp(llm) |
---|
1460 | d_t_va( 1 ) = -omdn*dtdp( 1 )+omdn*cor( 1 ) |
---|
1461 | d_t_va(llm) = -omup*dtdp(llm)!+omup*cor(llm) |
---|
1462 | |
---|
1463 | ! if(abs(rlat(1))>10.) then |
---|
1464 | ! Calculate the tendency due agestrophic motions |
---|
1465 | ! du_age = fcoriolis*(v-vg) |
---|
1466 | ! dv_age = fcoriolis*(ug-u) |
---|
1467 | ! endif |
---|
1468 | |
---|
1469 | ! call writefield_phy('d_t_va',d_t_va,llm) |
---|
1470 | |
---|
1471 | return |
---|
1472 | end |
---|
1473 | |
---|
1474 | !====================================================================== |
---|
1475 | |
---|
1476 | ! Subroutines for nudging |
---|
1477 | |
---|
1478 | Subroutine Nudge_RHT_init (paprs,pplay,t,q,t_targ,rh_targ) |
---|
1479 | ! ======================================================== |
---|
1480 | USE dimphy |
---|
1481 | |
---|
1482 | implicit none |
---|
1483 | |
---|
1484 | ! ======================================================== |
---|
1485 | REAL paprs(klon,klevp1) |
---|
1486 | REAL pplay(klon,klev) |
---|
1487 | ! |
---|
1488 | ! Variables d'etat |
---|
1489 | REAL t(klon,klev) |
---|
1490 | REAL q(klon,klev) |
---|
1491 | ! |
---|
1492 | ! Profiles cible |
---|
1493 | REAL t_targ(klon,klev) |
---|
1494 | REAL rh_targ(klon,klev) |
---|
1495 | ! |
---|
1496 | INTEGER k,i |
---|
1497 | REAL zx_qs |
---|
1498 | |
---|
1499 | ! Declaration des constantes et des fonctions thermodynamiques |
---|
1500 | ! |
---|
1501 | include "YOMCST.h" |
---|
1502 | include "YOETHF.h" |
---|
1503 | ! |
---|
1504 | ! ---------------------------------------- |
---|
1505 | ! Statement functions |
---|
1506 | include "FCTTRE.h" |
---|
1507 | ! ---------------------------------------- |
---|
1508 | ! |
---|
1509 | DO k = 1,klev |
---|
1510 | DO i = 1,klon |
---|
1511 | t_targ(i,k) = t(i,k) |
---|
1512 | IF (t(i,k).LT.RTT) THEN |
---|
1513 | zx_qs = qsats(t(i,k))/(pplay(i,k)) |
---|
1514 | ELSE |
---|
1515 | zx_qs = qsatl(t(i,k))/(pplay(i,k)) |
---|
1516 | ENDIF |
---|
1517 | rh_targ(i,k) = q(i,k)/zx_qs |
---|
1518 | ENDDO |
---|
1519 | ENDDO |
---|
1520 | print *, 't_targ',t_targ |
---|
1521 | print *, 'rh_targ',rh_targ |
---|
1522 | ! |
---|
1523 | ! |
---|
1524 | RETURN |
---|
1525 | END |
---|
1526 | |
---|
1527 | Subroutine Nudge_UV_init (paprs,pplay,u,v,u_targ,v_targ) |
---|
1528 | ! ======================================================== |
---|
1529 | USE dimphy |
---|
1530 | |
---|
1531 | implicit none |
---|
1532 | |
---|
1533 | ! ======================================================== |
---|
1534 | REAL paprs(klon,klevp1) |
---|
1535 | REAL pplay(klon,klev) |
---|
1536 | ! |
---|
1537 | ! Variables d'etat |
---|
1538 | REAL u(klon,klev) |
---|
1539 | REAL v(klon,klev) |
---|
1540 | ! |
---|
1541 | ! Profiles cible |
---|
1542 | REAL u_targ(klon,klev) |
---|
1543 | REAL v_targ(klon,klev) |
---|
1544 | ! |
---|
1545 | INTEGER k,i |
---|
1546 | ! |
---|
1547 | DO k = 1,klev |
---|
1548 | DO i = 1,klon |
---|
1549 | u_targ(i,k) = u(i,k) |
---|
1550 | v_targ(i,k) = v(i,k) |
---|
1551 | ENDDO |
---|
1552 | ENDDO |
---|
1553 | print *, 'u_targ',u_targ |
---|
1554 | print *, 'v_targ',v_targ |
---|
1555 | ! |
---|
1556 | ! |
---|
1557 | RETURN |
---|
1558 | END |
---|
1559 | |
---|
1560 | Subroutine Nudge_RHT (dtime,paprs,pplay,t_targ,rh_targ,t,q, & |
---|
1561 | & d_t,d_q) |
---|
1562 | ! ======================================================== |
---|
1563 | USE dimphy |
---|
1564 | |
---|
1565 | implicit none |
---|
1566 | |
---|
1567 | ! ======================================================== |
---|
1568 | REAL dtime |
---|
1569 | REAL paprs(klon,klevp1) |
---|
1570 | REAL pplay(klon,klev) |
---|
1571 | ! |
---|
1572 | ! Variables d'etat |
---|
1573 | REAL t(klon,klev) |
---|
1574 | REAL q(klon,klev) |
---|
1575 | ! |
---|
1576 | ! Tendances |
---|
1577 | REAL d_t(klon,klev) |
---|
1578 | REAL d_q(klon,klev) |
---|
1579 | ! |
---|
1580 | ! Profiles cible |
---|
1581 | REAL t_targ(klon,klev) |
---|
1582 | REAL rh_targ(klon,klev) |
---|
1583 | ! |
---|
1584 | ! Temps de relaxation |
---|
1585 | REAL tau |
---|
1586 | !c DATA tau /3600./ |
---|
1587 | !! DATA tau /5400./ |
---|
1588 | DATA tau /1800./ |
---|
1589 | ! |
---|
1590 | INTEGER k,i |
---|
1591 | REAL zx_qs, rh, tnew, d_rh, rhnew |
---|
1592 | |
---|
1593 | ! Declaration des constantes et des fonctions thermodynamiques |
---|
1594 | ! |
---|
1595 | include "YOMCST.h" |
---|
1596 | include "YOETHF.h" |
---|
1597 | ! |
---|
1598 | ! ---------------------------------------- |
---|
1599 | ! Statement functions |
---|
1600 | include "FCTTRE.h" |
---|
1601 | ! ---------------------------------------- |
---|
1602 | ! |
---|
1603 | print *,'dtime, tau ',dtime,tau |
---|
1604 | print *, 't_targ',t_targ |
---|
1605 | print *, 'rh_targ',rh_targ |
---|
1606 | print *,'temp ',t |
---|
1607 | print *,'hum ',q |
---|
1608 | ! |
---|
1609 | DO k = 1,klev |
---|
1610 | DO i = 1,klon |
---|
1611 | IF (paprs(i,1)-pplay(i,k) .GT. 10000.) THEN |
---|
1612 | IF (t(i,k).LT.RTT) THEN |
---|
1613 | zx_qs = qsats(t(i,k))/(pplay(i,k)) |
---|
1614 | ELSE |
---|
1615 | zx_qs = qsatl(t(i,k))/(pplay(i,k)) |
---|
1616 | ENDIF |
---|
1617 | rh = q(i,k)/zx_qs |
---|
1618 | ! |
---|
1619 | d_t(i,k) = d_t(i,k) + 1./tau*(t_targ(i,k)-t(i,k)) |
---|
1620 | d_rh = 1./tau*(rh_targ(i,k)-rh) |
---|
1621 | ! |
---|
1622 | tnew = t(i,k)+d_t(i,k)*dtime |
---|
1623 | !jyg< |
---|
1624 | ! Formule pour q : |
---|
1625 | ! d_q = (1/tau) [rh_targ*qsat(T_new) - q] |
---|
1626 | ! |
---|
1627 | ! Cette formule remplace d_q = (1/tau) [rh_targ - rh] qsat(T_new) |
---|
1628 | ! qui n'etait pas correcte. |
---|
1629 | ! |
---|
1630 | IF (tnew.LT.RTT) THEN |
---|
1631 | zx_qs = qsats(tnew)/(pplay(i,k)) |
---|
1632 | ELSE |
---|
1633 | zx_qs = qsatl(tnew)/(pplay(i,k)) |
---|
1634 | ENDIF |
---|
1635 | !! d_q(i,k) = d_q(i,k) + d_rh*zx_qs |
---|
1636 | d_q(i,k) = d_q(i,k) + (1./tau)*(rh_targ(i,k)*zx_qs - q(i,k)) |
---|
1637 | rhnew = (q(i,k)+d_q(i,k)*dtime)/zx_qs |
---|
1638 | ! |
---|
1639 | print *,' k,d_t,rh,d_rh,rhnew,d_q ', & |
---|
1640 | k,d_t(i,k),rh,d_rh,rhnew,d_q(i,k) |
---|
1641 | ENDIF |
---|
1642 | ! |
---|
1643 | ENDDO |
---|
1644 | ENDDO |
---|
1645 | ! |
---|
1646 | RETURN |
---|
1647 | END |
---|
1648 | |
---|
1649 | Subroutine Nudge_UV (dtime,paprs,pplay,u_targ,v_targ,u,v, & |
---|
1650 | & d_u,d_v) |
---|
1651 | ! ======================================================== |
---|
1652 | USE dimphy |
---|
1653 | |
---|
1654 | implicit none |
---|
1655 | |
---|
1656 | ! ======================================================== |
---|
1657 | REAL dtime |
---|
1658 | REAL paprs(klon,klevp1) |
---|
1659 | REAL pplay(klon,klev) |
---|
1660 | ! |
---|
1661 | ! Variables d'etat |
---|
1662 | REAL u(klon,klev) |
---|
1663 | REAL v(klon,klev) |
---|
1664 | ! |
---|
1665 | ! Tendances |
---|
1666 | REAL d_u(klon,klev) |
---|
1667 | REAL d_v(klon,klev) |
---|
1668 | ! |
---|
1669 | ! Profiles cible |
---|
1670 | REAL u_targ(klon,klev) |
---|
1671 | REAL v_targ(klon,klev) |
---|
1672 | ! |
---|
1673 | ! Temps de relaxation |
---|
1674 | REAL tau |
---|
1675 | !c DATA tau /3600./ |
---|
1676 | ! DATA tau /5400./ |
---|
1677 | DATA tau /43200./ |
---|
1678 | ! |
---|
1679 | INTEGER k,i |
---|
1680 | |
---|
1681 | ! |
---|
1682 | print *,'dtime, tau ',dtime,tau |
---|
1683 | print *, 'u_targ',u_targ |
---|
1684 | print *, 'v_targ',v_targ |
---|
1685 | print *,'zonal velocity ',u |
---|
1686 | print *,'meridional velocity ',v |
---|
1687 | DO k = 1,klev |
---|
1688 | DO i = 1,klon |
---|
1689 | !CR: nudging everywhere |
---|
1690 | ! IF (paprs(i,1)-pplay(i,k) .GT. 10000.) THEN |
---|
1691 | ! |
---|
1692 | d_u(i,k) = d_u(i,k) + 1./tau*(u_targ(i,k)-u(i,k)) |
---|
1693 | d_v(i,k) = d_v(i,k) + 1./tau*(v_targ(i,k)-v(i,k)) |
---|
1694 | ! |
---|
1695 | print *,' k,u,d_u,v,d_v ', & |
---|
1696 | k,u(i,k),d_u(i,k),v(i,k),d_v(i,k) |
---|
1697 | ! ENDIF |
---|
1698 | ! |
---|
1699 | ENDDO |
---|
1700 | ENDDO |
---|
1701 | ! |
---|
1702 | RETURN |
---|
1703 | END |
---|
1704 | |
---|
1705 | !===================================================================== |
---|
1706 | SUBROUTINE interp2_case_vertical(play,nlev_cas,plev_prof_cas & |
---|
1707 | & ,t_prof_cas,th_prof_cas,thv_prof_cas,thl_prof_cas & |
---|
1708 | & ,qv_prof_cas,ql_prof_cas,qi_prof_cas,u_prof_cas,v_prof_cas & |
---|
1709 | & ,ug_prof_cas,vg_prof_cas,vitw_prof_cas,omega_prof_cas & |
---|
1710 | & ,du_prof_cas,hu_prof_cas,vu_prof_cas,dv_prof_cas,hv_prof_cas,vv_prof_cas & |
---|
1711 | & ,dt_prof_cas,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas & |
---|
1712 | & ,dth_prof_cas,hth_prof_cas,vth_prof_cas & |
---|
1713 | ! |
---|
1714 | & ,t_mod_cas,theta_mod_cas,thv_mod_cas,thl_mod_cas & |
---|
1715 | & ,qv_mod_cas,ql_mod_cas,qi_mod_cas,u_mod_cas,v_mod_cas & |
---|
1716 | & ,ug_mod_cas,vg_mod_cas,w_mod_cas,omega_mod_cas & |
---|
1717 | & ,du_mod_cas,hu_mod_cas,vu_mod_cas,dv_mod_cas,hv_mod_cas,vv_mod_cas & |
---|
1718 | & ,dt_mod_cas,ht_mod_cas,vt_mod_cas,dtrad_mod_cas,dq_mod_cas,hq_mod_cas,vq_mod_cas & |
---|
1719 | & ,dth_mod_cas,hth_mod_cas,vth_mod_cas,mxcalc) |
---|
1720 | |
---|
1721 | implicit none |
---|
1722 | |
---|
1723 | #include "YOMCST.h" |
---|
1724 | #include "dimensions.h" |
---|
1725 | |
---|
1726 | !------------------------------------------------------------------------- |
---|
1727 | ! Vertical interpolation of generic case forcing data onto mod_casel levels |
---|
1728 | !------------------------------------------------------------------------- |
---|
1729 | |
---|
1730 | integer nlevmax |
---|
1731 | parameter (nlevmax=41) |
---|
1732 | integer nlev_cas,mxcalc |
---|
1733 | ! real play(llm), plev_prof(nlevmax) |
---|
1734 | ! real t_prof(nlevmax),q_prof(nlevmax) |
---|
1735 | ! real u_prof(nlevmax),v_prof(nlevmax), w_prof(nlevmax) |
---|
1736 | ! real ht_prof(nlevmax),vt_prof(nlevmax) |
---|
1737 | ! real hq_prof(nlevmax),vq_prof(nlevmax) |
---|
1738 | |
---|
1739 | real play(llm), plev_prof_cas(nlev_cas) |
---|
1740 | real t_prof_cas(nlev_cas),th_prof_cas(nlev_cas),thv_prof_cas(nlev_cas),thl_prof_cas(nlev_cas) |
---|
1741 | real qv_prof_cas(nlev_cas),ql_prof_cas(nlev_cas),qi_prof_cas(nlev_cas) |
---|
1742 | real u_prof_cas(nlev_cas),v_prof_cas(nlev_cas) |
---|
1743 | real ug_prof_cas(nlev_cas),vg_prof_cas(nlev_cas), vitw_prof_cas(nlev_cas),omega_prof_cas(nlev_cas) |
---|
1744 | real du_prof_cas(nlev_cas),hu_prof_cas(nlev_cas),vu_prof_cas(nlev_cas) |
---|
1745 | real dv_prof_cas(nlev_cas),hv_prof_cas(nlev_cas),vv_prof_cas(nlev_cas) |
---|
1746 | real dt_prof_cas(nlev_cas),ht_prof_cas(nlev_cas),vt_prof_cas(nlev_cas),dtrad_prof_cas(nlev_cas) |
---|
1747 | real dth_prof_cas(nlev_cas),hth_prof_cas(nlev_cas),vth_prof_cas(nlev_cas) |
---|
1748 | real dq_prof_cas(nlev_cas),hq_prof_cas(nlev_cas),vq_prof_cas(nlev_cas) |
---|
1749 | |
---|
1750 | real t_mod_cas(llm),theta_mod_cas(llm),thv_mod_cas(llm),thl_mod_cas(llm) |
---|
1751 | real qv_mod_cas(llm),ql_mod_cas(llm),qi_mod_cas(llm) |
---|
1752 | real u_mod_cas(llm),v_mod_cas(llm) |
---|
1753 | real ug_mod_cas(llm),vg_mod_cas(llm), w_mod_cas(llm),omega_mod_cas(llm) |
---|
1754 | real du_mod_cas(llm),hu_mod_cas(llm),vu_mod_cas(llm) |
---|
1755 | real dv_mod_cas(llm),hv_mod_cas(llm),vv_mod_cas(llm) |
---|
1756 | real dt_mod_cas(llm),ht_mod_cas(llm),vt_mod_cas(llm),dtrad_mod_cas(llm) |
---|
1757 | real dth_mod_cas(llm),hth_mod_cas(llm),vth_mod_cas(llm) |
---|
1758 | real dq_mod_cas(llm),hq_mod_cas(llm),vq_mod_cas(llm) |
---|
1759 | |
---|
1760 | integer l,k,k1,k2 |
---|
1761 | real frac,frac1,frac2,fact |
---|
1762 | |
---|
1763 | ! do l = 1, llm |
---|
1764 | ! print *,'debut interp2, play=',l,play(l) |
---|
1765 | ! enddo |
---|
1766 | ! do l = 1, nlev_cas |
---|
1767 | ! print *,'debut interp2, plev_prof_cas=',l,play(l),plev_prof_cas(l) |
---|
1768 | ! enddo |
---|
1769 | |
---|
1770 | do l = 1, llm |
---|
1771 | |
---|
1772 | if (play(l).ge.plev_prof_cas(nlev_cas)) then |
---|
1773 | |
---|
1774 | mxcalc=l |
---|
1775 | ! print *,'debut interp2, mxcalc=',mxcalc |
---|
1776 | k1=0 |
---|
1777 | k2=0 |
---|
1778 | |
---|
1779 | if (play(l).le.plev_prof_cas(1)) then |
---|
1780 | |
---|
1781 | do k = 1, nlev_cas-1 |
---|
1782 | if (play(l).le.plev_prof_cas(k).and. play(l).gt.plev_prof_cas(k+1)) then |
---|
1783 | k1=k |
---|
1784 | k2=k+1 |
---|
1785 | endif |
---|
1786 | enddo |
---|
1787 | |
---|
1788 | if (k1.eq.0 .or. k2.eq.0) then |
---|
1789 | write(*,*) 'PB! k1, k2 = ',k1,k2 |
---|
1790 | write(*,*) 'l,play(l) = ',l,play(l)/100 |
---|
1791 | do k = 1, nlev_cas-1 |
---|
1792 | write(*,*) 'k,plev_prof_cas(k) = ',k,plev_prof_cas(k)/100 |
---|
1793 | enddo |
---|
1794 | endif |
---|
1795 | |
---|
1796 | frac = (plev_prof_cas(k2)-play(l))/(plev_prof_cas(k2)-plev_prof_cas(k1)) |
---|
1797 | t_mod_cas(l)= t_prof_cas(k2) - frac*(t_prof_cas(k2)-t_prof_cas(k1)) |
---|
1798 | theta_mod_cas(l)= th_prof_cas(k2) - frac*(th_prof_cas(k2)-th_prof_cas(k1)) |
---|
1799 | if(theta_mod_cas(l).NE.0) t_mod_cas(l)= theta_mod_cas(l)*(play(l)/100000.)**(RD/RCPD) |
---|
1800 | thv_mod_cas(l)= thv_prof_cas(k2) - frac*(thv_prof_cas(k2)-thv_prof_cas(k1)) |
---|
1801 | thl_mod_cas(l)= thl_prof_cas(k2) - frac*(thl_prof_cas(k2)-thl_prof_cas(k1)) |
---|
1802 | qv_mod_cas(l)= qv_prof_cas(k2) - frac*(qv_prof_cas(k2)-qv_prof_cas(k1)) |
---|
1803 | ql_mod_cas(l)= ql_prof_cas(k2) - frac*(ql_prof_cas(k2)-ql_prof_cas(k1)) |
---|
1804 | qi_mod_cas(l)= qi_prof_cas(k2) - frac*(qi_prof_cas(k2)-qi_prof_cas(k1)) |
---|
1805 | u_mod_cas(l)= u_prof_cas(k2) - frac*(u_prof_cas(k2)-u_prof_cas(k1)) |
---|
1806 | v_mod_cas(l)= v_prof_cas(k2) - frac*(v_prof_cas(k2)-v_prof_cas(k1)) |
---|
1807 | ug_mod_cas(l)= ug_prof_cas(k2) - frac*(ug_prof_cas(k2)-ug_prof_cas(k1)) |
---|
1808 | vg_mod_cas(l)= vg_prof_cas(k2) - frac*(vg_prof_cas(k2)-vg_prof_cas(k1)) |
---|
1809 | w_mod_cas(l)= vitw_prof_cas(k2) - frac*(vitw_prof_cas(k2)-vitw_prof_cas(k1)) |
---|
1810 | omega_mod_cas(l)= omega_prof_cas(k2) - frac*(omega_prof_cas(k2)-omega_prof_cas(k1)) |
---|
1811 | du_mod_cas(l)= du_prof_cas(k2) - frac*(du_prof_cas(k2)-du_prof_cas(k1)) |
---|
1812 | hu_mod_cas(l)= hu_prof_cas(k2) - frac*(hu_prof_cas(k2)-hu_prof_cas(k1)) |
---|
1813 | vu_mod_cas(l)= vu_prof_cas(k2) - frac*(vu_prof_cas(k2)-vu_prof_cas(k1)) |
---|
1814 | dv_mod_cas(l)= dv_prof_cas(k2) - frac*(dv_prof_cas(k2)-dv_prof_cas(k1)) |
---|
1815 | hv_mod_cas(l)= hv_prof_cas(k2) - frac*(hv_prof_cas(k2)-hv_prof_cas(k1)) |
---|
1816 | vv_mod_cas(l)= vv_prof_cas(k2) - frac*(vv_prof_cas(k2)-vv_prof_cas(k1)) |
---|
1817 | dt_mod_cas(l)= dt_prof_cas(k2) - frac*(dt_prof_cas(k2)-dt_prof_cas(k1)) |
---|
1818 | ht_mod_cas(l)= ht_prof_cas(k2) - frac*(ht_prof_cas(k2)-ht_prof_cas(k1)) |
---|
1819 | vt_mod_cas(l)= vt_prof_cas(k2) - frac*(vt_prof_cas(k2)-vt_prof_cas(k1)) |
---|
1820 | dth_mod_cas(l)= dth_prof_cas(k2) - frac*(dth_prof_cas(k2)-dth_prof_cas(k1)) |
---|
1821 | hth_mod_cas(l)= hth_prof_cas(k2) - frac*(hth_prof_cas(k2)-hth_prof_cas(k1)) |
---|
1822 | vth_mod_cas(l)= vth_prof_cas(k2) - frac*(vth_prof_cas(k2)-vth_prof_cas(k1)) |
---|
1823 | dq_mod_cas(l)= dq_prof_cas(k2) - frac*(dq_prof_cas(k2)-dq_prof_cas(k1)) |
---|
1824 | hq_mod_cas(l)= hq_prof_cas(k2) - frac*(hq_prof_cas(k2)-hq_prof_cas(k1)) |
---|
1825 | vq_mod_cas(l)= vq_prof_cas(k2) - frac*(vq_prof_cas(k2)-vq_prof_cas(k1)) |
---|
1826 | dtrad_mod_cas(l)= dtrad_prof_cas(k2) - frac*(dtrad_prof_cas(k2)-dtrad_prof_cas(k1)) |
---|
1827 | |
---|
1828 | else !play>plev_prof_cas(1) |
---|
1829 | |
---|
1830 | k1=1 |
---|
1831 | k2=2 |
---|
1832 | print *,'interp2_vert, k1,k2=',plev_prof_cas(k1),plev_prof_cas(k2) |
---|
1833 | frac1 = (play(l)-plev_prof_cas(k2))/(plev_prof_cas(k1)-plev_prof_cas(k2)) |
---|
1834 | frac2 = (play(l)-plev_prof_cas(k1))/(plev_prof_cas(k1)-plev_prof_cas(k2)) |
---|
1835 | t_mod_cas(l)= frac1*t_prof_cas(k1) - frac2*t_prof_cas(k2) |
---|
1836 | theta_mod_cas(l)= frac1*th_prof_cas(k1) - frac2*th_prof_cas(k2) |
---|
1837 | if(theta_mod_cas(l).NE.0) t_mod_cas(l)= theta_mod_cas(l)*(play(l)/100000.)**(RD/RCPD) |
---|
1838 | thv_mod_cas(l)= frac1*thv_prof_cas(k1) - frac2*thv_prof_cas(k2) |
---|
1839 | thl_mod_cas(l)= frac1*thl_prof_cas(k1) - frac2*thl_prof_cas(k2) |
---|
1840 | qv_mod_cas(l)= frac1*qv_prof_cas(k1) - frac2*qv_prof_cas(k2) |
---|
1841 | ql_mod_cas(l)= frac1*ql_prof_cas(k1) - frac2*ql_prof_cas(k2) |
---|
1842 | qi_mod_cas(l)= frac1*qi_prof_cas(k1) - frac2*qi_prof_cas(k2) |
---|
1843 | u_mod_cas(l)= frac1*u_prof_cas(k1) - frac2*u_prof_cas(k2) |
---|
1844 | v_mod_cas(l)= frac1*v_prof_cas(k1) - frac2*v_prof_cas(k2) |
---|
1845 | ug_mod_cas(l)= frac1*ug_prof_cas(k1) - frac2*ug_prof_cas(k2) |
---|
1846 | vg_mod_cas(l)= frac1*vg_prof_cas(k1) - frac2*vg_prof_cas(k2) |
---|
1847 | w_mod_cas(l)= frac1*vitw_prof_cas(k1) - frac2*vitw_prof_cas(k2) |
---|
1848 | omega_mod_cas(l)= frac1*omega_prof_cas(k1) - frac2*omega_prof_cas(k2) |
---|
1849 | du_mod_cas(l)= frac1*du_prof_cas(k1) - frac2*du_prof_cas(k2) |
---|
1850 | hu_mod_cas(l)= frac1*hu_prof_cas(k1) - frac2*hu_prof_cas(k2) |
---|
1851 | vu_mod_cas(l)= frac1*vu_prof_cas(k1) - frac2*vu_prof_cas(k2) |
---|
1852 | dv_mod_cas(l)= frac1*dv_prof_cas(k1) - frac2*dv_prof_cas(k2) |
---|
1853 | hv_mod_cas(l)= frac1*hv_prof_cas(k1) - frac2*hv_prof_cas(k2) |
---|
1854 | vv_mod_cas(l)= frac1*vv_prof_cas(k1) - frac2*vv_prof_cas(k2) |
---|
1855 | dt_mod_cas(l)= frac1*dt_prof_cas(k1) - frac2*dt_prof_cas(k2) |
---|
1856 | ht_mod_cas(l)= frac1*ht_prof_cas(k1) - frac2*ht_prof_cas(k2) |
---|
1857 | vt_mod_cas(l)= frac1*vt_prof_cas(k1) - frac2*vt_prof_cas(k2) |
---|
1858 | dth_mod_cas(l)= frac1*dth_prof_cas(k1) - frac2*dth_prof_cas(k2) |
---|
1859 | hth_mod_cas(l)= frac1*hth_prof_cas(k1) - frac2*hth_prof_cas(k2) |
---|
1860 | vth_mod_cas(l)= frac1*vth_prof_cas(k1) - frac2*vth_prof_cas(k2) |
---|
1861 | dq_mod_cas(l)= frac1*dq_prof_cas(k1) - frac2*dq_prof_cas(k2) |
---|
1862 | hq_mod_cas(l)= frac1*hq_prof_cas(k1) - frac2*hq_prof_cas(k2) |
---|
1863 | vq_mod_cas(l)= frac1*vq_prof_cas(k1) - frac2*vq_prof_cas(k2) |
---|
1864 | dtrad_mod_cas(l)= frac1*dtrad_prof_cas(k1) - frac2*dtrad_prof_cas(k2) |
---|
1865 | |
---|
1866 | endif ! play.le.plev_prof_cas(1) |
---|
1867 | |
---|
1868 | else ! above max altitude of forcing file |
---|
1869 | |
---|
1870 | !jyg |
---|
1871 | fact=20.*(plev_prof_cas(nlev_cas)-play(l))/plev_prof_cas(nlev_cas) !jyg |
---|
1872 | fact = max(fact,0.) !jyg |
---|
1873 | fact = exp(-fact) !jyg |
---|
1874 | t_mod_cas(l)= t_prof_cas(nlev_cas) !jyg |
---|
1875 | theta_mod_cas(l)= th_prof_cas(nlev_cas) !jyg |
---|
1876 | thv_mod_cas(l)= thv_prof_cas(nlev_cas) !jyg |
---|
1877 | thl_mod_cas(l)= thl_prof_cas(nlev_cas) !jyg |
---|
1878 | qv_mod_cas(l)= qv_prof_cas(nlev_cas)*fact !jyg |
---|
1879 | ql_mod_cas(l)= ql_prof_cas(nlev_cas)*fact !jyg |
---|
1880 | qi_mod_cas(l)= qi_prof_cas(nlev_cas)*fact !jyg |
---|
1881 | u_mod_cas(l)= u_prof_cas(nlev_cas)*fact !jyg |
---|
1882 | v_mod_cas(l)= v_prof_cas(nlev_cas)*fact !jyg |
---|
1883 | ug_mod_cas(l)= ug_prof_cas(nlev_cas)*fact !jyg |
---|
1884 | vg_mod_cas(l)= vg_prof_cas(nlev_cas)*fact !jyg |
---|
1885 | w_mod_cas(l)= 0.0 !jyg |
---|
1886 | omega_mod_cas(l)= 0.0 !jyg |
---|
1887 | du_mod_cas(l)= du_prof_cas(nlev_cas)*fact |
---|
1888 | hu_mod_cas(l)= hu_prof_cas(nlev_cas)*fact !jyg |
---|
1889 | vu_mod_cas(l)= vu_prof_cas(nlev_cas)*fact !jyg |
---|
1890 | dv_mod_cas(l)= dv_prof_cas(nlev_cas)*fact |
---|
1891 | hv_mod_cas(l)= hv_prof_cas(nlev_cas)*fact !jyg |
---|
1892 | vv_mod_cas(l)= vv_prof_cas(nlev_cas)*fact !jyg |
---|
1893 | dt_mod_cas(l)= dt_prof_cas(nlev_cas) |
---|
1894 | ht_mod_cas(l)= ht_prof_cas(nlev_cas) !jyg |
---|
1895 | vt_mod_cas(l)= vt_prof_cas(nlev_cas) !jyg |
---|
1896 | dth_mod_cas(l)= dth_prof_cas(nlev_cas) |
---|
1897 | hth_mod_cas(l)= hth_prof_cas(nlev_cas) !jyg |
---|
1898 | vth_mod_cas(l)= vth_prof_cas(nlev_cas) !jyg |
---|
1899 | dq_mod_cas(l)= dq_prof_cas(nlev_cas)*fact |
---|
1900 | hq_mod_cas(l)= hq_prof_cas(nlev_cas)*fact !jyg |
---|
1901 | vq_mod_cas(l)= vq_prof_cas(nlev_cas)*fact !jyg |
---|
1902 | dtrad_mod_cas(l)= dtrad_prof_cas(nlev_cas)*fact !jyg |
---|
1903 | |
---|
1904 | endif ! play |
---|
1905 | |
---|
1906 | enddo ! l |
---|
1907 | |
---|
1908 | return |
---|
1909 | end |
---|
1910 | !***************************************************************************** |
---|
1911 | |
---|
1912 | |
---|
1913 | |
---|
1914 | |
---|