/[lmdze]/trunk/phylmd/physiq.f
ViewVC logotype

Contents of /trunk/phylmd/physiq.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 182 - (show annotations)
Wed Mar 16 11:11:27 2016 UTC (8 years, 2 months ago) by guez
Original Path: trunk/Sources/phylmd/physiq.f
File size: 57203 byte(s)
Replaced integer variable iflag_con of module clesphys2 by logical
variable conv_emanuel.

1 module physiq_m
2
3 IMPLICIT none
4
5 contains
6
7 SUBROUTINE physiq(lafin, dayvrai, time, paprs, play, pphi, pphis, u, v, t, &
8 qx, omega, d_u, d_v, d_t, d_qx)
9
10 ! From phylmd/physiq.F, version 1.22 2006/02/20 09:38:28
11 ! (subversion revision 678)
12
13 ! Author: Z. X. Li (LMD/CNRS) 1993
14
15 ! This is the main procedure for the "physics" part of the program.
16
17 use aaam_bud_m, only: aaam_bud
18 USE abort_gcm_m, ONLY: abort_gcm
19 use aeropt_m, only: aeropt
20 use ajsec_m, only: ajsec
21 use calltherm_m, only: calltherm
22 USE clesphys, ONLY: cdhmax, cdmmax, ecrit_hf, ecrit_ins, ecrit_mth, &
23 ecrit_reg, ecrit_tra, ksta, ksta_ter, ok_kzmin
24 USE clesphys2, ONLY: cycle_diurne, conv_emanuel, nbapp_rad, new_oliq, &
25 ok_orodr, ok_orolf
26 USE clmain_m, ONLY: clmain
27 use clouds_gno_m, only: clouds_gno
28 use comconst, only: dtphys
29 USE comgeomphy, ONLY: airephy
30 USE concvl_m, ONLY: concvl
31 USE conf_gcm_m, ONLY: offline, raz_date, day_step, iphysiq
32 USE conf_phys_m, ONLY: conf_phys
33 use conflx_m, only: conflx
34 USE ctherm, ONLY: iflag_thermals, nsplit_thermals
35 use diagcld2_m, only: diagcld2
36 use diagetpq_m, only: diagetpq
37 use diagphy_m, only: diagphy
38 USE dimens_m, ONLY: llm, nqmx
39 USE dimphy, ONLY: klon
40 USE dimsoil, ONLY: nsoilmx
41 use drag_noro_m, only: drag_noro
42 use dynetat0_m, only: day_ref, annee_ref
43 USE fcttre, ONLY: foeew, qsatl, qsats, thermcep
44 use fisrtilp_m, only: fisrtilp
45 USE hgardfou_m, ONLY: hgardfou
46 USE indicesol, ONLY: clnsurf, epsfra, is_lic, is_oce, is_sic, is_ter, &
47 nbsrf
48 USE ini_histins_m, ONLY: ini_histins
49 use netcdf95, only: NF95_CLOSE
50 use newmicro_m, only: newmicro
51 use nuage_m, only: nuage
52 USE orbite_m, ONLY: orbite
53 USE ozonecm_m, ONLY: ozonecm
54 USE phyetat0_m, ONLY: phyetat0, rlat, rlon
55 USE phyredem_m, ONLY: phyredem
56 USE phyredem0_m, ONLY: phyredem0
57 USE phystokenc_m, ONLY: phystokenc
58 USE phytrac_m, ONLY: phytrac
59 USE qcheck_m, ONLY: qcheck
60 use radlwsw_m, only: radlwsw
61 use readsulfate_m, only: readsulfate
62 use readsulfate_preind_m, only: readsulfate_preind
63 use yoegwd, only: sugwd
64 USE suphec_m, ONLY: rcpd, retv, rg, rlvtt, romega, rsigma, rtt
65 use transp_m, only: transp
66 use transp_lay_m, only: transp_lay
67 use unit_nml_m, only: unit_nml
68 USE ymds2ju_m, ONLY: ymds2ju
69 USE yoethf_m, ONLY: r2es, rvtmp2
70 use zenang_m, only: zenang
71
72 logical, intent(in):: lafin ! dernier passage
73
74 integer, intent(in):: dayvrai
75 ! current day number, based at value 1 on January 1st of annee_ref
76
77 REAL, intent(in):: time ! heure de la journ\'ee en fraction de jour
78
79 REAL, intent(in):: paprs(:, :) ! (klon, llm + 1)
80 ! pression pour chaque inter-couche, en Pa
81
82 REAL, intent(in):: play(:, :) ! (klon, llm)
83 ! pression pour le mileu de chaque couche (en Pa)
84
85 REAL, intent(in):: pphi(:, :) ! (klon, llm)
86 ! géopotentiel de chaque couche (référence sol)
87
88 REAL, intent(in):: pphis(:) ! (klon) géopotentiel du sol
89
90 REAL, intent(in):: u(:, :) ! (klon, llm)
91 ! vitesse dans la direction X (de O a E) en m/s
92
93 REAL, intent(in):: v(:, :) ! (klon, llm) vitesse Y (de S a N) en m/s
94 REAL, intent(in):: t(:, :) ! (klon, llm) temperature (K)
95
96 REAL, intent(in):: qx(:, :, :) ! (klon, llm, nqmx)
97 ! (humidit\'e sp\'ecifique et fractions massiques des autres traceurs)
98
99 REAL, intent(in):: omega(:, :) ! (klon, llm) vitesse verticale en Pa/s
100 REAL, intent(out):: d_u(:, :) ! (klon, llm) tendance physique de "u" (m s-2)
101 REAL, intent(out):: d_v(:, :) ! (klon, llm) tendance physique de "v" (m s-2)
102 REAL, intent(out):: d_t(:, :) ! (klon, llm) tendance physique de "t" (K/s)
103
104 REAL, intent(out):: d_qx(:, :, :) ! (klon, llm, nqmx)
105 ! tendance physique de "qx" (s-1)
106
107 ! Local:
108
109 LOGICAL:: firstcal = .true.
110
111 LOGICAL ok_gust ! pour activer l'effet des gust sur flux surface
112 PARAMETER (ok_gust = .FALSE.)
113
114 LOGICAL, PARAMETER:: check = .FALSE.
115 ! Verifier la conservation du modele en eau
116
117 LOGICAL, PARAMETER:: ok_stratus = .FALSE.
118 ! Ajouter artificiellement les stratus
119
120 logical:: ok_journe = .false., ok_mensuel = .true., ok_instan = .false.
121 ! sorties journalieres, mensuelles et instantanees dans les
122 ! fichiers histday, histmth et histins
123
124 LOGICAL ok_region ! sortir le fichier regional
125 PARAMETER (ok_region = .FALSE.)
126
127 ! pour phsystoke avec thermiques
128 REAL fm_therm(klon, llm + 1)
129 REAL entr_therm(klon, llm)
130 real, save:: q2(klon, llm + 1, nbsrf)
131
132 INTEGER, PARAMETER:: ivap = 1 ! indice de traceur pour vapeur d'eau
133 INTEGER, PARAMETER:: iliq = 2 ! indice de traceur pour eau liquide
134
135 REAL, save:: t_ancien(klon, llm), q_ancien(klon, llm)
136 LOGICAL, save:: ancien_ok
137
138 REAL d_t_dyn(klon, llm) ! tendance dynamique pour "t" (K/s)
139 REAL d_q_dyn(klon, llm) ! tendance dynamique pour "q" (kg/kg/s)
140
141 real da(klon, llm), phi(klon, llm, llm), mp(klon, llm)
142
143 REAL swdn0(klon, llm + 1), swdn(klon, llm + 1)
144 REAL swup0(klon, llm + 1), swup(klon, llm + 1)
145 SAVE swdn0, swdn, swup0, swup
146
147 REAL lwdn0(klon, llm + 1), lwdn(klon, llm + 1)
148 REAL lwup0(klon, llm + 1), lwup(klon, llm + 1)
149 SAVE lwdn0, lwdn, lwup0, lwup
150
151 ! Amip2
152 ! variables a une pression donnee
153
154 integer nlevSTD
155 PARAMETER(nlevSTD = 17)
156
157 ! prw: precipitable water
158 real prw(klon)
159
160 ! flwp, fiwp = Liquid Water Path & Ice Water Path (kg/m2)
161 ! flwc, fiwc = Liquid Water Content & Ice Water Content (kg/kg)
162 REAL flwp(klon), fiwp(klon)
163 REAL flwc(klon, llm), fiwc(klon, llm)
164
165 INTEGER kmax, lmax
166 PARAMETER(kmax = 8, lmax = 8)
167 INTEGER kmaxm1, lmaxm1
168 PARAMETER(kmaxm1 = kmax - 1, lmaxm1 = lmax - 1)
169
170 ! Variables propres a la physique
171
172 INTEGER, save:: radpas
173 ! Radiative transfer computations are made every "radpas" call to
174 ! "physiq".
175
176 REAL radsol(klon)
177 SAVE radsol ! bilan radiatif au sol calcule par code radiatif
178
179 INTEGER:: itap = 0 ! number of calls to "physiq"
180
181 REAL, save:: ftsol(klon, nbsrf) ! skin temperature of surface fraction
182
183 REAL, save:: ftsoil(klon, nsoilmx, nbsrf)
184 ! soil temperature of surface fraction
185
186 REAL, save:: fevap(klon, nbsrf) ! evaporation
187 REAL fluxlat(klon, nbsrf)
188 SAVE fluxlat
189
190 REAL, save:: fqsurf(klon, nbsrf)
191 ! humidite de l'air au contact de la surface
192
193 REAL, save:: qsol(klon)
194 ! column-density of water in soil, in kg m-2
195
196 REAL, save:: fsnow(klon, nbsrf) ! epaisseur neigeuse
197 REAL, save:: falbe(klon, nbsrf) ! albedo visible par type de surface
198
199 ! Param\`etres de l'orographie \`a l'\'echelle sous-maille (OESM) :
200 REAL, save:: zmea(klon) ! orographie moyenne
201 REAL, save:: zstd(klon) ! deviation standard de l'OESM
202 REAL, save:: zsig(klon) ! pente de l'OESM
203 REAL, save:: zgam(klon) ! anisotropie de l'OESM
204 REAL, save:: zthe(klon) ! orientation de l'OESM
205 REAL, save:: zpic(klon) ! Maximum de l'OESM
206 REAL, save:: zval(klon) ! Minimum de l'OESM
207 REAL, save:: rugoro(klon) ! longueur de rugosite de l'OESM
208 REAL zulow(klon), zvlow(klon)
209 INTEGER igwd, itest(klon)
210
211 REAL agesno(klon, nbsrf)
212 SAVE agesno ! age de la neige
213
214 REAL run_off_lic_0(klon)
215 SAVE run_off_lic_0
216 !KE43
217 ! Variables liees a la convection de K. Emanuel (sb):
218
219 REAL Ma(klon, llm) ! undilute upward mass flux
220 SAVE Ma
221 REAL qcondc(klon, llm) ! in-cld water content from convect
222 SAVE qcondc
223 REAL, save:: sig1(klon, llm), w01(klon, llm)
224 REAL, save:: wd(klon)
225
226 ! Variables pour la couche limite (al1):
227
228 REAL cdragh(klon) ! drag coefficient pour T and Q
229 REAL cdragm(klon) ! drag coefficient pour vent
230
231 ! Pour phytrac :
232 REAL ycoefh(klon, llm) ! coef d'echange pour phytrac
233 REAL yu1(klon) ! vents dans la premiere couche U
234 REAL yv1(klon) ! vents dans la premiere couche V
235 REAL ffonte(klon, nbsrf) !Flux thermique utilise pour fondre la neige
236 REAL fqcalving(klon, nbsrf) !Flux d'eau "perdue" par la surface
237 ! !et necessaire pour limiter la
238 ! !hauteur de neige, en kg/m2/s
239 REAL zxffonte(klon), zxfqcalving(klon)
240
241 REAL pfrac_impa(klon, llm)! Produits des coefs lessivage impaction
242 save pfrac_impa
243 REAL pfrac_nucl(klon, llm)! Produits des coefs lessivage nucleation
244 save pfrac_nucl
245 REAL pfrac_1nucl(klon, llm)! Produits des coefs lessi nucl (alpha = 1)
246 save pfrac_1nucl
247 REAL frac_impa(klon, llm) ! fractions d'aerosols lessivees (impaction)
248 REAL frac_nucl(klon, llm) ! idem (nucleation)
249
250 REAL, save:: rain_fall(klon)
251 ! liquid water mass flux (kg/m2/s), positive down
252
253 REAL, save:: snow_fall(klon)
254 ! solid water mass flux (kg/m2/s), positive down
255
256 REAL rain_tiedtke(klon), snow_tiedtke(klon)
257
258 REAL evap(klon), devap(klon) ! evaporation and its derivative
259 REAL sens(klon), dsens(klon) ! chaleur sensible et sa derivee
260 REAL dlw(klon) ! derivee infra rouge
261 SAVE dlw
262 REAL bils(klon) ! bilan de chaleur au sol
263 REAL, save:: fder(klon) ! Derive de flux (sensible et latente)
264 REAL ve(klon) ! integr. verticale du transport meri. de l'energie
265 REAL vq(klon) ! integr. verticale du transport meri. de l'eau
266 REAL ue(klon) ! integr. verticale du transport zonal de l'energie
267 REAL uq(klon) ! integr. verticale du transport zonal de l'eau
268
269 REAL, save:: frugs(klon, nbsrf) ! longueur de rugosite
270 REAL zxrugs(klon) ! longueur de rugosite
271
272 ! Conditions aux limites
273
274 INTEGER julien
275 INTEGER, SAVE:: lmt_pas ! number of time steps of "physics" per day
276 REAL, save:: pctsrf(klon, nbsrf) ! percentage of surface
277 REAL pctsrf_new(klon, nbsrf) ! pourcentage surfaces issus d'ORCHIDEE
278 REAL, save:: albsol(klon) ! albedo du sol total visible
279 REAL, SAVE:: wo(klon, llm) ! column density of ozone in a cell, in kDU
280
281 real, save:: clwcon(klon, llm), rnebcon(klon, llm)
282 real, save:: clwcon0(klon, llm), rnebcon0(klon, llm)
283
284 REAL rhcl(klon, llm) ! humiditi relative ciel clair
285 REAL dialiq(klon, llm) ! eau liquide nuageuse
286 REAL diafra(klon, llm) ! fraction nuageuse
287 REAL cldliq(klon, llm) ! eau liquide nuageuse
288 REAL cldfra(klon, llm) ! fraction nuageuse
289 REAL cldtau(klon, llm) ! epaisseur optique
290 REAL cldemi(klon, llm) ! emissivite infrarouge
291
292 REAL fluxq(klon, llm, nbsrf) ! flux turbulent d'humidite
293 REAL fluxt(klon, llm, nbsrf) ! flux turbulent de chaleur
294 REAL fluxu(klon, llm, nbsrf) ! flux turbulent de vitesse u
295 REAL fluxv(klon, llm, nbsrf) ! flux turbulent de vitesse v
296
297 REAL zxfluxt(klon, llm)
298 REAL zxfluxq(klon, llm)
299 REAL zxfluxu(klon, llm)
300 REAL zxfluxv(klon, llm)
301
302 ! Le rayonnement n'est pas calcul\'e tous les pas, il faut donc que
303 ! les variables soient r\'emanentes.
304 REAL, save:: heat(klon, llm) ! chauffage solaire
305 REAL, save:: heat0(klon, llm) ! chauffage solaire ciel clair
306 REAL, save:: cool(klon, llm) ! refroidissement infrarouge
307 REAL, save:: cool0(klon, llm) ! refroidissement infrarouge ciel clair
308 REAL, save:: topsw(klon), toplw(klon), solsw(klon)
309 REAL, save:: sollw(klon) ! rayonnement infrarouge montant \`a la surface
310 real, save:: sollwdown(klon) ! downward LW flux at surface
311 REAL, save:: topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)
312 REAL, save:: albpla(klon)
313 REAL fsollw(klon, nbsrf) ! bilan flux IR pour chaque sous surface
314 REAL fsolsw(klon, nbsrf) ! flux solaire absorb. pour chaque sous surface
315
316 REAL conv_q(klon, llm) ! convergence de l'humidite (kg/kg/s)
317 REAL conv_t(klon, llm) ! convergence of temperature (K/s)
318
319 REAL cldl(klon), cldm(klon), cldh(klon) !nuages bas, moyen et haut
320 REAL cldt(klon), cldq(klon) !nuage total, eau liquide integree
321
322 REAL zxtsol(klon), zxqsurf(klon), zxsnow(klon), zxfluxlat(klon)
323
324 REAL dist, mu0(klon), fract(klon)
325 real longi
326 REAL z_avant(klon), z_apres(klon), z_factor(klon)
327 REAL za, zb
328 REAL zx_t, zx_qs, zcor
329 real zqsat(klon, llm)
330 INTEGER i, k, iq, nsrf
331 REAL, PARAMETER:: t_coup = 234.
332 REAL zphi(klon, llm)
333
334 ! cf. AM Variables pour la CLA (hbtm2)
335
336 REAL, SAVE:: pblh(klon, nbsrf) ! Hauteur de couche limite
337 REAL, SAVE:: plcl(klon, nbsrf) ! Niveau de condensation de la CLA
338 REAL, SAVE:: capCL(klon, nbsrf) ! CAPE de couche limite
339 REAL, SAVE:: oliqCL(klon, nbsrf) ! eau_liqu integree de couche limite
340 REAL, SAVE:: cteiCL(klon, nbsrf) ! cloud top instab. crit. couche limite
341 REAL, SAVE:: pblt(klon, nbsrf) ! T a la Hauteur de couche limite
342 REAL, SAVE:: therm(klon, nbsrf)
343 REAL, SAVE:: trmb1(klon, nbsrf) ! deep_cape
344 REAL, SAVE:: trmb2(klon, nbsrf) ! inhibition
345 REAL, SAVE:: trmb3(klon, nbsrf) ! Point Omega
346 ! Grdeurs de sorties
347 REAL s_pblh(klon), s_lcl(klon), s_capCL(klon)
348 REAL s_oliqCL(klon), s_cteiCL(klon), s_pblt(klon)
349 REAL s_therm(klon), s_trmb1(klon), s_trmb2(klon)
350 REAL s_trmb3(klon)
351
352 ! Variables pour la convection de K. Emanuel :
353
354 REAL upwd(klon, llm) ! saturated updraft mass flux
355 REAL dnwd(klon, llm) ! saturated downdraft mass flux
356 REAL dnwd0(klon, llm) ! unsaturated downdraft mass flux
357 REAL cape(klon) ! CAPE
358 SAVE cape
359
360 INTEGER iflagctrl(klon) ! flag fonctionnement de convect
361
362 ! Variables du changement
363
364 ! con: convection
365 ! lsc: large scale condensation
366 ! ajs: ajustement sec
367 ! eva: \'evaporation de l'eau liquide nuageuse
368 ! vdf: vertical diffusion in boundary layer
369 REAL d_t_con(klon, llm), d_q_con(klon, llm)
370 REAL d_u_con(klon, llm), d_v_con(klon, llm)
371 REAL d_t_lsc(klon, llm), d_q_lsc(klon, llm), d_ql_lsc(klon, llm)
372 REAL d_t_ajs(klon, llm), d_q_ajs(klon, llm)
373 REAL d_u_ajs(klon, llm), d_v_ajs(klon, llm)
374 REAL rneb(klon, llm)
375
376 REAL mfu(klon, llm), mfd(klon, llm)
377 REAL pen_u(klon, llm), pen_d(klon, llm)
378 REAL pde_u(klon, llm), pde_d(klon, llm)
379 INTEGER kcbot(klon), kctop(klon), kdtop(klon)
380 REAL pmflxr(klon, llm + 1), pmflxs(klon, llm + 1)
381 REAL prfl(klon, llm + 1), psfl(klon, llm + 1)
382
383 INTEGER, save:: ibas_con(klon), itop_con(klon)
384
385 REAL rain_con(klon), rain_lsc(klon)
386 REAL, save:: snow_con(klon)
387 real snow_lsc(klon)
388 REAL d_ts(klon, nbsrf)
389
390 REAL d_u_vdf(klon, llm), d_v_vdf(klon, llm)
391 REAL d_t_vdf(klon, llm), d_q_vdf(klon, llm)
392
393 REAL d_u_oro(klon, llm), d_v_oro(klon, llm)
394 REAL d_t_oro(klon, llm)
395 REAL d_u_lif(klon, llm), d_v_lif(klon, llm)
396 REAL d_t_lif(klon, llm)
397
398 REAL, save:: ratqs(klon, llm)
399 real ratqss(klon, llm), ratqsc(klon, llm)
400 real:: ratqsbas = 0.01, ratqshaut = 0.3
401
402 ! Parametres lies au nouveau schema de nuages (SB, PDF)
403 real:: fact_cldcon = 0.375
404 real:: facttemps = 1.e-4
405 logical:: ok_newmicro = .true.
406 real facteur
407
408 integer:: iflag_cldcon = 1
409 logical ptconv(klon, llm)
410
411 ! Variables pour effectuer les appels en s\'erie :
412
413 REAL t_seri(klon, llm), q_seri(klon, llm)
414 REAL ql_seri(klon, llm)
415 REAL u_seri(klon, llm), v_seri(klon, llm)
416 REAL tr_seri(klon, llm, nqmx - 2)
417
418 REAL zx_rh(klon, llm)
419
420 REAL zustrdr(klon), zvstrdr(klon)
421 REAL zustrli(klon), zvstrli(klon)
422 REAL zustrph(klon), zvstrph(klon)
423 REAL aam, torsfc
424
425 REAL zx_tmp_fi2d(klon) ! variable temporaire grille physique
426
427 INTEGER, SAVE:: nid_ins
428
429 REAL ve_lay(klon, llm) ! transport meri. de l'energie a chaque niveau vert.
430 REAL vq_lay(klon, llm) ! transport meri. de l'eau a chaque niveau vert.
431 REAL ue_lay(klon, llm) ! transport zonal de l'energie a chaque niveau vert.
432 REAL uq_lay(klon, llm) ! transport zonal de l'eau a chaque niveau vert.
433
434 real date0
435
436 ! Variables li\'ees au bilan d'\'energie et d'enthalpie :
437 REAL ztsol(klon)
438 REAL d_h_vcol, d_qt, d_ec
439 REAL, SAVE:: d_h_vcol_phy
440 REAL zero_v(klon)
441 CHARACTER(LEN = 20) tit
442 INTEGER:: ip_ebil = 0 ! print level for energy conservation diagnostics
443 INTEGER:: if_ebil = 0 ! verbosity for diagnostics of energy conservation
444
445 REAL d_t_ec(klon, llm) ! tendance due \`a la conversion Ec -> E thermique
446 REAL ZRCPD
447
448 REAL t2m(klon, nbsrf), q2m(klon, nbsrf) ! temperature and humidity at 2 m
449 REAL u10m(klon, nbsrf), v10m(klon, nbsrf) ! vents a 10 m
450 REAL zt2m(klon), zq2m(klon) ! temp., hum. 2 m moyenne s/ 1 maille
451 REAL zu10m(klon), zv10m(klon) ! vents a 10 m moyennes s/1 maille
452
453 ! Aerosol effects:
454
455 REAL sulfate(klon, llm) ! SO4 aerosol concentration (micro g/m3)
456
457 REAL, save:: sulfate_pi(klon, llm)
458 ! SO4 aerosol concentration, in \mu g/m3, pre-industrial value
459
460 REAL cldtaupi(klon, llm)
461 ! cloud optical thickness for pre-industrial (pi) aerosols
462
463 REAL re(klon, llm) ! Cloud droplet effective radius
464 REAL fl(klon, llm) ! denominator of re
465
466 ! Aerosol optical properties
467 REAL, save:: tau_ae(klon, llm, 2), piz_ae(klon, llm, 2)
468 REAL, save:: cg_ae(klon, llm, 2)
469
470 REAL topswad(klon), solswad(klon) ! aerosol direct effect
471 REAL topswai(klon), solswai(klon) ! aerosol indirect effect
472
473 REAL aerindex(klon) ! POLDER aerosol index
474
475 LOGICAL:: ok_ade = .false. ! apply aerosol direct effect
476 LOGICAL:: ok_aie = .false. ! apply aerosol indirect effect
477
478 REAL:: bl95_b0 = 2., bl95_b1 = 0.2
479 ! Parameters in equation (D) of Boucher and Lohmann (1995, Tellus
480 ! B). They link cloud droplet number concentration to aerosol mass
481 ! concentration.
482
483 SAVE u10m
484 SAVE v10m
485 SAVE t2m
486 SAVE q2m
487 SAVE ffonte
488 SAVE fqcalving
489 SAVE rain_con
490 SAVE topswai
491 SAVE topswad
492 SAVE solswai
493 SAVE solswad
494 SAVE d_u_con
495 SAVE d_v_con
496
497 real zmasse(klon, llm)
498 ! (column-density of mass of air in a cell, in kg m-2)
499
500 real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2
501 integer, save:: ncid_startphy, itau_phy
502
503 namelist /physiq_nml/ ok_journe, ok_mensuel, ok_instan, fact_cldcon, &
504 facttemps, ok_newmicro, iflag_cldcon, ratqsbas, ratqshaut, if_ebil, &
505 ok_ade, ok_aie, bl95_b0, bl95_b1, iflag_thermals, nsplit_thermals
506
507 !----------------------------------------------------------------
508
509 IF (if_ebil >= 1) zero_v = 0.
510 IF (nqmx < 2) CALL abort_gcm('physiq', &
511 'eaux vapeur et liquide sont indispensables')
512
513 test_firstcal: IF (firstcal) THEN
514 ! initialiser
515 u10m = 0.
516 v10m = 0.
517 t2m = 0.
518 q2m = 0.
519 ffonte = 0.
520 fqcalving = 0.
521 piz_ae = 0.
522 tau_ae = 0.
523 cg_ae = 0.
524 rain_con = 0.
525 snow_con = 0.
526 topswai = 0.
527 topswad = 0.
528 solswai = 0.
529 solswad = 0.
530
531 d_u_con = 0.
532 d_v_con = 0.
533 rnebcon0 = 0.
534 clwcon0 = 0.
535 rnebcon = 0.
536 clwcon = 0.
537
538 pblh =0. ! Hauteur de couche limite
539 plcl =0. ! Niveau de condensation de la CLA
540 capCL =0. ! CAPE de couche limite
541 oliqCL =0. ! eau_liqu integree de couche limite
542 cteiCL =0. ! cloud top instab. crit. couche limite
543 pblt =0. ! T a la Hauteur de couche limite
544 therm =0.
545 trmb1 =0. ! deep_cape
546 trmb2 =0. ! inhibition
547 trmb3 =0. ! Point Omega
548
549 IF (if_ebil >= 1) d_h_vcol_phy = 0.
550
551 iflag_thermals = 0
552 nsplit_thermals = 1
553 print *, "Enter namelist 'physiq_nml'."
554 read(unit=*, nml=physiq_nml)
555 write(unit_nml, nml=physiq_nml)
556
557 call conf_phys
558
559 ! Initialiser les compteurs:
560
561 frugs = 0.
562 CALL phyetat0(pctsrf, ftsol, ftsoil, fqsurf, qsol, &
563 fsnow, falbe, fevap, rain_fall, snow_fall, solsw, sollw, dlw, &
564 radsol, frugs, agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, &
565 t_ancien, q_ancien, ancien_ok, rnebcon, ratqs, clwcon, &
566 run_off_lic_0, sig1, w01, ncid_startphy, itau_phy)
567
568 ! ATTENTION : il faudra a terme relire q2 dans l'etat initial
569 q2 = 1e-8
570
571 lmt_pas = day_step / iphysiq
572 print *, 'Number of time steps of "physics" per day: ', lmt_pas
573
574 radpas = lmt_pas / nbapp_rad
575
576 ! On remet le calendrier a zero
577 IF (raz_date) itau_phy = 0
578
579 CALL printflag(radpas, ok_journe, ok_instan, ok_region)
580
581 ! Initialisation pour le sch\'ema de convection d'Emanuel :
582 IF (conv_emanuel) THEN
583 ibas_con = 1
584 itop_con = 1
585 ENDIF
586
587 IF (ok_orodr) THEN
588 rugoro = MAX(1e-5, zstd * zsig / 2)
589 CALL SUGWD(paprs, play)
590 else
591 rugoro = 0.
592 ENDIF
593
594 ecrit_ins = NINT(ecrit_ins/dtphys)
595 ecrit_hf = NINT(ecrit_hf/dtphys)
596 ecrit_mth = NINT(ecrit_mth/dtphys)
597 ecrit_tra = NINT(86400.*ecrit_tra/dtphys)
598 ecrit_reg = NINT(ecrit_reg/dtphys)
599
600 ! Initialisation des sorties
601
602 call ini_histins(dtphys, ok_instan, nid_ins, itau_phy)
603 CALL ymds2ju(annee_ref, 1, day_ref, 0., date0)
604 ! Positionner date0 pour initialisation de ORCHIDEE
605 print *, 'physiq date0: ', date0
606 CALL phyredem0(lmt_pas, itau_phy)
607 ENDIF test_firstcal
608
609 ! We will modify variables *_seri and we will not touch variables
610 ! u, v, t, qx:
611 t_seri = t
612 u_seri = u
613 v_seri = v
614 q_seri = qx(:, :, ivap)
615 ql_seri = qx(:, :, iliq)
616 tr_seri = qx(:, :, 3:nqmx)
617
618 ztsol = sum(ftsol * pctsrf, dim = 2)
619
620 IF (if_ebil >= 1) THEN
621 tit = 'after dynamics'
622 CALL diagetpq(airephy, tit, ip_ebil, 1, 1, dtphys, t_seri, q_seri, &
623 ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
624 ! Comme les tendances de la physique sont ajout\'es dans la
625 ! dynamique, la variation d'enthalpie par la dynamique devrait
626 ! \^etre \'egale \`a la variation de la physique au pas de temps
627 ! pr\'ec\'edent. Donc la somme de ces 2 variations devrait \^etre
628 ! nulle.
629 call diagphy(airephy, tit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
630 zero_v, zero_v, zero_v, zero_v, ztsol, d_h_vcol + d_h_vcol_phy, &
631 d_qt, 0.)
632 END IF
633
634 ! Diagnostic de la tendance dynamique :
635 IF (ancien_ok) THEN
636 DO k = 1, llm
637 DO i = 1, klon
638 d_t_dyn(i, k) = (t_seri(i, k) - t_ancien(i, k)) / dtphys
639 d_q_dyn(i, k) = (q_seri(i, k) - q_ancien(i, k)) / dtphys
640 ENDDO
641 ENDDO
642 ELSE
643 DO k = 1, llm
644 DO i = 1, klon
645 d_t_dyn(i, k) = 0.
646 d_q_dyn(i, k) = 0.
647 ENDDO
648 ENDDO
649 ancien_ok = .TRUE.
650 ENDIF
651
652 ! Ajouter le geopotentiel du sol:
653 DO k = 1, llm
654 DO i = 1, klon
655 zphi(i, k) = pphi(i, k) + pphis(i)
656 ENDDO
657 ENDDO
658
659 ! Check temperatures:
660 CALL hgardfou(t_seri, ftsol)
661
662 ! Incrémenter le compteur de la physique
663 itap = itap + 1
664 julien = MOD(dayvrai, 360)
665 if (julien == 0) julien = 360
666
667 forall (k = 1: llm) zmasse(:, k) = (paprs(:, k) - paprs(:, k + 1)) / rg
668
669 ! Prescrire l'ozone :
670 wo = ozonecm(REAL(julien), paprs)
671
672 ! \'Evaporation de l'eau liquide nuageuse :
673 DO k = 1, llm
674 DO i = 1, klon
675 zb = MAX(0., ql_seri(i, k))
676 t_seri(i, k) = t_seri(i, k) &
677 - zb * RLVTT / RCPD / (1. + RVTMP2 * q_seri(i, k))
678 q_seri(i, k) = q_seri(i, k) + zb
679 ENDDO
680 ENDDO
681 ql_seri = 0.
682
683 IF (if_ebil >= 2) THEN
684 tit = 'after reevap'
685 CALL diagetpq(airephy, tit, ip_ebil, 2, 1, dtphys, t_seri, q_seri, &
686 ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
687 call diagphy(airephy, tit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
688 zero_v, zero_v, zero_v, zero_v, ztsol, d_h_vcol, d_qt, d_ec)
689 END IF
690
691 frugs = MAX(frugs, 0.000015)
692 zxrugs = sum(frugs * pctsrf, dim = 2)
693
694 ! Calculs nécessaires au calcul de l'albedo dans l'interface avec
695 ! la surface.
696
697 CALL orbite(REAL(julien), longi, dist)
698 IF (cycle_diurne) THEN
699 CALL zenang(longi, time, dtphys * radpas, mu0, fract)
700 ELSE
701 mu0 = - 999.999
702 ENDIF
703
704 ! Calcul de l'abedo moyen par maille
705 albsol = sum(falbe * pctsrf, dim = 2)
706
707 ! R\'epartition sous maille des flux longwave et shortwave
708 ! R\'epartition du longwave par sous-surface lin\'earis\'ee
709
710 forall (nsrf = 1: nbsrf)
711 fsollw(:, nsrf) = sollw + 4. * RSIGMA * ztsol**3 &
712 * (ztsol - ftsol(:, nsrf))
713 fsolsw(:, nsrf) = solsw * (1. - falbe(:, nsrf)) / (1. - albsol)
714 END forall
715
716 fder = dlw
717
718 ! Couche limite:
719
720 CALL clmain(dtphys, itap, pctsrf, pctsrf_new, t_seri, q_seri, u_seri, &
721 v_seri, julien, mu0, ftsol, cdmmax, cdhmax, ksta, ksta_ter, &
722 ok_kzmin, ftsoil, qsol, paprs, play, fsnow, fqsurf, fevap, falbe, &
723 fluxlat, rain_fall, snow_fall, fsolsw, fsollw, fder, rlat, frugs, &
724 firstcal, agesno, rugoro, d_t_vdf, d_q_vdf, d_u_vdf, d_v_vdf, d_ts, &
725 fluxt, fluxq, fluxu, fluxv, cdragh, cdragm, q2, dsens, devap, &
726 ycoefh, yu1, yv1, t2m, q2m, u10m, v10m, pblh, capCL, oliqCL, cteiCL, &
727 pblT, therm, trmb1, trmb2, trmb3, plcl, fqcalving, ffonte, &
728 run_off_lic_0)
729
730 ! Incr\'ementation des flux
731
732 zxfluxt = 0.
733 zxfluxq = 0.
734 zxfluxu = 0.
735 zxfluxv = 0.
736 DO nsrf = 1, nbsrf
737 DO k = 1, llm
738 DO i = 1, klon
739 zxfluxt(i, k) = zxfluxt(i, k) + fluxt(i, k, nsrf) * pctsrf(i, nsrf)
740 zxfluxq(i, k) = zxfluxq(i, k) + fluxq(i, k, nsrf) * pctsrf(i, nsrf)
741 zxfluxu(i, k) = zxfluxu(i, k) + fluxu(i, k, nsrf) * pctsrf(i, nsrf)
742 zxfluxv(i, k) = zxfluxv(i, k) + fluxv(i, k, nsrf) * pctsrf(i, nsrf)
743 END DO
744 END DO
745 END DO
746 DO i = 1, klon
747 sens(i) = - zxfluxt(i, 1) ! flux de chaleur sensible au sol
748 evap(i) = - zxfluxq(i, 1) ! flux d'\'evaporation au sol
749 fder(i) = dlw(i) + dsens(i) + devap(i)
750 ENDDO
751
752 DO k = 1, llm
753 DO i = 1, klon
754 t_seri(i, k) = t_seri(i, k) + d_t_vdf(i, k)
755 q_seri(i, k) = q_seri(i, k) + d_q_vdf(i, k)
756 u_seri(i, k) = u_seri(i, k) + d_u_vdf(i, k)
757 v_seri(i, k) = v_seri(i, k) + d_v_vdf(i, k)
758 ENDDO
759 ENDDO
760
761 IF (if_ebil >= 2) THEN
762 tit = 'after clmain'
763 CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
764 ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
765 call diagphy(airephy, tit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
766 sens, evap, zero_v, zero_v, ztsol, d_h_vcol, d_qt, d_ec)
767 END IF
768
769 ! Update surface temperature:
770
771 DO i = 1, klon
772 zxtsol(i) = 0.
773 zxfluxlat(i) = 0.
774
775 zt2m(i) = 0.
776 zq2m(i) = 0.
777 zu10m(i) = 0.
778 zv10m(i) = 0.
779 zxffonte(i) = 0.
780 zxfqcalving(i) = 0.
781
782 s_pblh(i) = 0.
783 s_lcl(i) = 0.
784 s_capCL(i) = 0.
785 s_oliqCL(i) = 0.
786 s_cteiCL(i) = 0.
787 s_pblT(i) = 0.
788 s_therm(i) = 0.
789 s_trmb1(i) = 0.
790 s_trmb2(i) = 0.
791 s_trmb3(i) = 0.
792
793 IF (abs(pctsrf(i, is_ter) + pctsrf(i, is_lic) + pctsrf(i, is_oce) &
794 + pctsrf(i, is_sic) - 1.) > EPSFRA) print *, &
795 'physiq : probl\`eme sous surface au point ', i, &
796 pctsrf(i, 1 : nbsrf)
797 ENDDO
798 DO nsrf = 1, nbsrf
799 DO i = 1, klon
800 ftsol(i, nsrf) = ftsol(i, nsrf) + d_ts(i, nsrf)
801 zxtsol(i) = zxtsol(i) + ftsol(i, nsrf)*pctsrf(i, nsrf)
802 zxfluxlat(i) = zxfluxlat(i) + fluxlat(i, nsrf)*pctsrf(i, nsrf)
803
804 zt2m(i) = zt2m(i) + t2m(i, nsrf)*pctsrf(i, nsrf)
805 zq2m(i) = zq2m(i) + q2m(i, nsrf)*pctsrf(i, nsrf)
806 zu10m(i) = zu10m(i) + u10m(i, nsrf)*pctsrf(i, nsrf)
807 zv10m(i) = zv10m(i) + v10m(i, nsrf)*pctsrf(i, nsrf)
808 zxffonte(i) = zxffonte(i) + ffonte(i, nsrf)*pctsrf(i, nsrf)
809 zxfqcalving(i) = zxfqcalving(i) + &
810 fqcalving(i, nsrf)*pctsrf(i, nsrf)
811 s_pblh(i) = s_pblh(i) + pblh(i, nsrf)*pctsrf(i, nsrf)
812 s_lcl(i) = s_lcl(i) + plcl(i, nsrf)*pctsrf(i, nsrf)
813 s_capCL(i) = s_capCL(i) + capCL(i, nsrf) *pctsrf(i, nsrf)
814 s_oliqCL(i) = s_oliqCL(i) + oliqCL(i, nsrf) *pctsrf(i, nsrf)
815 s_cteiCL(i) = s_cteiCL(i) + cteiCL(i, nsrf) *pctsrf(i, nsrf)
816 s_pblT(i) = s_pblT(i) + pblT(i, nsrf) *pctsrf(i, nsrf)
817 s_therm(i) = s_therm(i) + therm(i, nsrf) *pctsrf(i, nsrf)
818 s_trmb1(i) = s_trmb1(i) + trmb1(i, nsrf) *pctsrf(i, nsrf)
819 s_trmb2(i) = s_trmb2(i) + trmb2(i, nsrf) *pctsrf(i, nsrf)
820 s_trmb3(i) = s_trmb3(i) + trmb3(i, nsrf) *pctsrf(i, nsrf)
821 ENDDO
822 ENDDO
823
824 ! Si une sous-fraction n'existe pas, elle prend la température moyenne :
825 DO nsrf = 1, nbsrf
826 DO i = 1, klon
827 IF (pctsrf(i, nsrf) < epsfra) ftsol(i, nsrf) = zxtsol(i)
828
829 IF (pctsrf(i, nsrf) < epsfra) t2m(i, nsrf) = zt2m(i)
830 IF (pctsrf(i, nsrf) < epsfra) q2m(i, nsrf) = zq2m(i)
831 IF (pctsrf(i, nsrf) < epsfra) u10m(i, nsrf) = zu10m(i)
832 IF (pctsrf(i, nsrf) < epsfra) v10m(i, nsrf) = zv10m(i)
833 IF (pctsrf(i, nsrf) < epsfra) ffonte(i, nsrf) = zxffonte(i)
834 IF (pctsrf(i, nsrf) < epsfra) &
835 fqcalving(i, nsrf) = zxfqcalving(i)
836 IF (pctsrf(i, nsrf) < epsfra) pblh(i, nsrf) = s_pblh(i)
837 IF (pctsrf(i, nsrf) < epsfra) plcl(i, nsrf) = s_lcl(i)
838 IF (pctsrf(i, nsrf) < epsfra) capCL(i, nsrf) = s_capCL(i)
839 IF (pctsrf(i, nsrf) < epsfra) oliqCL(i, nsrf) = s_oliqCL(i)
840 IF (pctsrf(i, nsrf) < epsfra) cteiCL(i, nsrf) = s_cteiCL(i)
841 IF (pctsrf(i, nsrf) < epsfra) pblT(i, nsrf) = s_pblT(i)
842 IF (pctsrf(i, nsrf) < epsfra) therm(i, nsrf) = s_therm(i)
843 IF (pctsrf(i, nsrf) < epsfra) trmb1(i, nsrf) = s_trmb1(i)
844 IF (pctsrf(i, nsrf) < epsfra) trmb2(i, nsrf) = s_trmb2(i)
845 IF (pctsrf(i, nsrf) < epsfra) trmb3(i, nsrf) = s_trmb3(i)
846 ENDDO
847 ENDDO
848
849 ! Calculer la dérive du flux infrarouge
850
851 DO i = 1, klon
852 dlw(i) = - 4. * RSIGMA * zxtsol(i)**3
853 ENDDO
854
855 IF (check) print *, "avantcon = ", qcheck(paprs, q_seri, ql_seri)
856
857 ! Appeler la convection (au choix)
858
859 if (conv_emanuel) then
860 da = 0.
861 mp = 0.
862 phi = 0.
863 CALL concvl(dtphys, paprs, play, t_seri, q_seri, u_seri, v_seri, sig1, &
864 w01, d_t_con, d_q_con, d_u_con, d_v_con, rain_con, snow_con, &
865 ibas_con, itop_con, upwd, dnwd, dnwd0, Ma, cape, iflagctrl, &
866 qcondc, wd, pmflxr, da, phi, mp)
867 clwcon0 = qcondc
868 mfu = upwd + dnwd
869 IF (.NOT. ok_gust) wd = 0.
870
871 IF (thermcep) THEN
872 zqsat = MIN(0.5, r2es * FOEEW(t_seri, rtt >= t_seri) / play)
873 zqsat = zqsat / (1. - retv * zqsat)
874 ELSE
875 zqsat = merge(qsats(t_seri), qsatl(t_seri), t_seri < t_coup) / play
876 ENDIF
877
878 ! Properties of convective clouds
879 clwcon0 = fact_cldcon * clwcon0
880 call clouds_gno(klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, &
881 rnebcon0)
882
883 mfd = 0.
884 pen_u = 0.
885 pen_d = 0.
886 pde_d = 0.
887 pde_u = 0.
888 else
889 conv_q = d_q_dyn + d_q_vdf / dtphys
890 conv_t = d_t_dyn + d_t_vdf / dtphys
891 z_avant = sum((q_seri + ql_seri) * zmasse, dim=2)
892 CALL conflx(dtphys, paprs, play, t_seri(:, llm:1:- 1), &
893 q_seri(:, llm:1:- 1), conv_t, conv_q, zxfluxq(:, 1), omega, &
894 d_t_con, d_q_con, rain_con, snow_con, mfu(:, llm:1:- 1), &
895 mfd(:, llm:1:- 1), pen_u, pde_u, pen_d, pde_d, kcbot, kctop, &
896 kdtop, pmflxr, pmflxs)
897 WHERE (rain_con < 0.) rain_con = 0.
898 WHERE (snow_con < 0.) snow_con = 0.
899 ibas_con = llm + 1 - kcbot
900 itop_con = llm + 1 - kctop
901 END if
902
903 DO k = 1, llm
904 DO i = 1, klon
905 t_seri(i, k) = t_seri(i, k) + d_t_con(i, k)
906 q_seri(i, k) = q_seri(i, k) + d_q_con(i, k)
907 u_seri(i, k) = u_seri(i, k) + d_u_con(i, k)
908 v_seri(i, k) = v_seri(i, k) + d_v_con(i, k)
909 ENDDO
910 ENDDO
911
912 IF (if_ebil >= 2) THEN
913 tit = 'after convect'
914 CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
915 ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
916 call diagphy(airephy, tit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
917 zero_v, zero_v, rain_con, snow_con, ztsol, d_h_vcol, d_qt, d_ec)
918 END IF
919
920 IF (check) THEN
921 za = qcheck(paprs, q_seri, ql_seri)
922 print *, "aprescon = ", za
923 zx_t = 0.
924 za = 0.
925 DO i = 1, klon
926 za = za + airephy(i)/REAL(klon)
927 zx_t = zx_t + (rain_con(i)+ &
928 snow_con(i))*airephy(i)/REAL(klon)
929 ENDDO
930 zx_t = zx_t/za*dtphys
931 print *, "Precip = ", zx_t
932 ENDIF
933
934 IF (.not. conv_emanuel) THEN
935 z_apres = sum((q_seri + ql_seri) * zmasse, dim=2)
936 z_factor = (z_avant - (rain_con + snow_con) * dtphys) / z_apres
937 DO k = 1, llm
938 DO i = 1, klon
939 IF (z_factor(i) > 1. + 1E-8 .OR. z_factor(i) < 1. - 1E-8) THEN
940 q_seri(i, k) = q_seri(i, k) * z_factor(i)
941 ENDIF
942 ENDDO
943 ENDDO
944 ENDIF
945
946 ! Convection s\`eche (thermiques ou ajustement)
947
948 d_t_ajs = 0.
949 d_u_ajs = 0.
950 d_v_ajs = 0.
951 d_q_ajs = 0.
952 fm_therm = 0.
953 entr_therm = 0.
954
955 if (iflag_thermals == 0) then
956 ! Ajustement sec
957 CALL ajsec(paprs, play, t_seri, q_seri, d_t_ajs, d_q_ajs)
958 t_seri = t_seri + d_t_ajs
959 q_seri = q_seri + d_q_ajs
960 else
961 ! Thermiques
962 call calltherm(dtphys, play, paprs, pphi, u_seri, v_seri, t_seri, &
963 q_seri, d_u_ajs, d_v_ajs, d_t_ajs, d_q_ajs, fm_therm, entr_therm)
964 endif
965
966 IF (if_ebil >= 2) THEN
967 tit = 'after dry_adjust'
968 CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
969 ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
970 END IF
971
972 ! Caclul des ratqs
973
974 ! ratqs convectifs \`a l'ancienne en fonction de (q(z = 0) - q) / q
975 ! on \'ecrase le tableau ratqsc calcul\'e par clouds_gno
976 if (iflag_cldcon == 1) then
977 do k = 1, llm
978 do i = 1, klon
979 if(ptconv(i, k)) then
980 ratqsc(i, k) = ratqsbas + fact_cldcon &
981 * (q_seri(i, 1) - q_seri(i, k)) / q_seri(i, k)
982 else
983 ratqsc(i, k) = 0.
984 endif
985 enddo
986 enddo
987 endif
988
989 ! ratqs stables
990 do k = 1, llm
991 do i = 1, klon
992 ratqss(i, k) = ratqsbas + (ratqshaut - ratqsbas) &
993 * min((paprs(i, 1) - play(i, k)) / (paprs(i, 1) - 3e4), 1.)
994 enddo
995 enddo
996
997 ! ratqs final
998 if (iflag_cldcon == 1 .or. iflag_cldcon == 2) then
999 ! les ratqs sont une conbinaison de ratqss et ratqsc
1000 ! ratqs final
1001 ! 1e4 (en gros 3 heures), en dur pour le moment, est le temps de
1002 ! relaxation des ratqs
1003 ratqs = max(ratqs * exp(- dtphys * facttemps), ratqss)
1004 ratqs = max(ratqs, ratqsc)
1005 else
1006 ! on ne prend que le ratqs stable pour fisrtilp
1007 ratqs = ratqss
1008 endif
1009
1010 CALL fisrtilp(dtphys, paprs, play, t_seri, q_seri, ptconv, ratqs, &
1011 d_t_lsc, d_q_lsc, d_ql_lsc, rneb, cldliq, rain_lsc, snow_lsc, &
1012 pfrac_impa, pfrac_nucl, pfrac_1nucl, frac_impa, frac_nucl, prfl, &
1013 psfl, rhcl)
1014
1015 WHERE (rain_lsc < 0) rain_lsc = 0.
1016 WHERE (snow_lsc < 0) snow_lsc = 0.
1017 DO k = 1, llm
1018 DO i = 1, klon
1019 t_seri(i, k) = t_seri(i, k) + d_t_lsc(i, k)
1020 q_seri(i, k) = q_seri(i, k) + d_q_lsc(i, k)
1021 ql_seri(i, k) = ql_seri(i, k) + d_ql_lsc(i, k)
1022 cldfra(i, k) = rneb(i, k)
1023 IF (.NOT.new_oliq) cldliq(i, k) = ql_seri(i, k)
1024 ENDDO
1025 ENDDO
1026 IF (check) THEN
1027 za = qcheck(paprs, q_seri, ql_seri)
1028 print *, "apresilp = ", za
1029 zx_t = 0.
1030 za = 0.
1031 DO i = 1, klon
1032 za = za + airephy(i)/REAL(klon)
1033 zx_t = zx_t + (rain_lsc(i) &
1034 + snow_lsc(i))*airephy(i)/REAL(klon)
1035 ENDDO
1036 zx_t = zx_t/za*dtphys
1037 print *, "Precip = ", zx_t
1038 ENDIF
1039
1040 IF (if_ebil >= 2) THEN
1041 tit = 'after fisrt'
1042 CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1043 ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
1044 call diagphy(airephy, tit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
1045 zero_v, zero_v, rain_lsc, snow_lsc, ztsol, d_h_vcol, d_qt, d_ec)
1046 END IF
1047
1048 ! PRESCRIPTION DES NUAGES POUR LE RAYONNEMENT
1049
1050 ! 1. NUAGES CONVECTIFS
1051
1052 IF (iflag_cldcon <= - 1) THEN
1053 ! seulement pour Tiedtke
1054 snow_tiedtke = 0.
1055 if (iflag_cldcon == - 1) then
1056 rain_tiedtke = rain_con
1057 else
1058 rain_tiedtke = 0.
1059 do k = 1, llm
1060 do i = 1, klon
1061 if (d_q_con(i, k) < 0.) then
1062 rain_tiedtke(i) = rain_tiedtke(i) - d_q_con(i, k)/dtphys &
1063 *zmasse(i, k)
1064 endif
1065 enddo
1066 enddo
1067 endif
1068
1069 ! Nuages diagnostiques pour Tiedtke
1070 CALL diagcld1(paprs, play, rain_tiedtke, snow_tiedtke, ibas_con, &
1071 itop_con, diafra, dialiq)
1072 DO k = 1, llm
1073 DO i = 1, klon
1074 IF (diafra(i, k) > cldfra(i, k)) THEN
1075 cldliq(i, k) = dialiq(i, k)
1076 cldfra(i, k) = diafra(i, k)
1077 ENDIF
1078 ENDDO
1079 ENDDO
1080 ELSE IF (iflag_cldcon == 3) THEN
1081 ! On prend pour les nuages convectifs le maximum du calcul de
1082 ! la convection et du calcul du pas de temps pr\'ec\'edent diminu\'e
1083 ! d'un facteur facttemps.
1084 facteur = dtphys * facttemps
1085 do k = 1, llm
1086 do i = 1, klon
1087 rnebcon(i, k) = rnebcon(i, k) * facteur
1088 if (rnebcon0(i, k) * clwcon0(i, k) &
1089 > rnebcon(i, k) * clwcon(i, k)) then
1090 rnebcon(i, k) = rnebcon0(i, k)
1091 clwcon(i, k) = clwcon0(i, k)
1092 endif
1093 enddo
1094 enddo
1095
1096 ! On prend la somme des fractions nuageuses et des contenus en eau
1097 cldfra = min(max(cldfra, rnebcon), 1.)
1098 cldliq = cldliq + rnebcon*clwcon
1099 ENDIF
1100
1101 ! 2. Nuages stratiformes
1102
1103 IF (ok_stratus) THEN
1104 CALL diagcld2(paprs, play, t_seri, q_seri, diafra, dialiq)
1105 DO k = 1, llm
1106 DO i = 1, klon
1107 IF (diafra(i, k) > cldfra(i, k)) THEN
1108 cldliq(i, k) = dialiq(i, k)
1109 cldfra(i, k) = diafra(i, k)
1110 ENDIF
1111 ENDDO
1112 ENDDO
1113 ENDIF
1114
1115 ! Precipitation totale
1116 DO i = 1, klon
1117 rain_fall(i) = rain_con(i) + rain_lsc(i)
1118 snow_fall(i) = snow_con(i) + snow_lsc(i)
1119 ENDDO
1120
1121 IF (if_ebil >= 2) CALL diagetpq(airephy, "after diagcld", ip_ebil, 2, 2, &
1122 dtphys, t_seri, q_seri, ql_seri, u_seri, v_seri, paprs, d_h_vcol, &
1123 d_qt, d_ec)
1124
1125 ! Humidit\'e relative pour diagnostic :
1126 DO k = 1, llm
1127 DO i = 1, klon
1128 zx_t = t_seri(i, k)
1129 IF (thermcep) THEN
1130 zx_qs = r2es * FOEEW(zx_t, rtt >= zx_t)/play(i, k)
1131 zx_qs = MIN(0.5, zx_qs)
1132 zcor = 1./(1. - retv*zx_qs)
1133 zx_qs = zx_qs*zcor
1134 ELSE
1135 IF (zx_t < t_coup) THEN
1136 zx_qs = qsats(zx_t)/play(i, k)
1137 ELSE
1138 zx_qs = qsatl(zx_t)/play(i, k)
1139 ENDIF
1140 ENDIF
1141 zx_rh(i, k) = q_seri(i, k)/zx_qs
1142 zqsat(i, k) = zx_qs
1143 ENDDO
1144 ENDDO
1145
1146 ! Introduce the aerosol direct and first indirect radiative forcings:
1147 IF (ok_ade .OR. ok_aie) THEN
1148 ! Get sulfate aerosol distribution :
1149 CALL readsulfate(dayvrai, time, firstcal, sulfate)
1150 CALL readsulfate_preind(dayvrai, time, firstcal, sulfate_pi)
1151
1152 CALL aeropt(play, paprs, t_seri, sulfate, rhcl, tau_ae, piz_ae, cg_ae, &
1153 aerindex)
1154 ELSE
1155 tau_ae = 0.
1156 piz_ae = 0.
1157 cg_ae = 0.
1158 ENDIF
1159
1160 ! Param\`etres optiques des nuages et quelques param\`etres pour
1161 ! diagnostics :
1162 if (ok_newmicro) then
1163 CALL newmicro(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, &
1164 cldh, cldl, cldm, cldt, cldq, flwp, fiwp, flwc, fiwc, ok_aie, &
1165 sulfate, sulfate_pi, bl95_b0, bl95_b1, cldtaupi, re, fl)
1166 else
1167 CALL nuage(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, cldh, &
1168 cldl, cldm, cldt, cldq, ok_aie, sulfate, sulfate_pi, bl95_b0, &
1169 bl95_b1, cldtaupi, re, fl)
1170 endif
1171
1172 IF (MOD(itap - 1, radpas) == 0) THEN
1173 ! Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.
1174 ! Calcul de l'abedo moyen par maille
1175 albsol = sum(falbe * pctsrf, dim = 2)
1176
1177 ! Rayonnement (compatible Arpege-IFS) :
1178 CALL radlwsw(dist, mu0, fract, paprs, play, zxtsol, albsol, t_seri, &
1179 q_seri, wo, cldfra, cldemi, cldtau, heat, heat0, cool, cool0, &
1180 radsol, albpla, topsw, toplw, solsw, sollw, sollwdown, topsw0, &
1181 toplw0, solsw0, sollw0, lwdn0, lwdn, lwup0, lwup, swdn0, swdn, &
1182 swup0, swup, ok_ade, ok_aie, tau_ae, piz_ae, cg_ae, topswad, &
1183 solswad, cldtaupi, topswai, solswai)
1184 ENDIF
1185
1186 ! Ajouter la tendance des rayonnements (tous les pas)
1187
1188 DO k = 1, llm
1189 DO i = 1, klon
1190 t_seri(i, k) = t_seri(i, k) + (heat(i, k) - cool(i, k)) * dtphys/86400.
1191 ENDDO
1192 ENDDO
1193
1194 IF (if_ebil >= 2) THEN
1195 tit = 'after rad'
1196 CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1197 ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
1198 call diagphy(airephy, tit, ip_ebil, topsw, toplw, solsw, sollw, &
1199 zero_v, zero_v, zero_v, zero_v, ztsol, d_h_vcol, d_qt, d_ec)
1200 END IF
1201
1202 ! Calculer l'hydrologie de la surface
1203 DO i = 1, klon
1204 zxqsurf(i) = 0.
1205 zxsnow(i) = 0.
1206 ENDDO
1207 DO nsrf = 1, nbsrf
1208 DO i = 1, klon
1209 zxqsurf(i) = zxqsurf(i) + fqsurf(i, nsrf)*pctsrf(i, nsrf)
1210 zxsnow(i) = zxsnow(i) + fsnow(i, nsrf)*pctsrf(i, nsrf)
1211 ENDDO
1212 ENDDO
1213
1214 ! Calculer le bilan du sol et la d\'erive de temp\'erature (couplage)
1215
1216 DO i = 1, klon
1217 bils(i) = radsol(i) - sens(i) + zxfluxlat(i)
1218 ENDDO
1219
1220 ! Param\'etrisation de l'orographie \`a l'\'echelle sous-maille :
1221
1222 IF (ok_orodr) THEN
1223 ! S\'election des points pour lesquels le sch\'ema est actif :
1224 igwd = 0
1225 DO i = 1, klon
1226 itest(i) = 0
1227 IF (zpic(i) - zmea(i) > 100. .AND. zstd(i) > 10.) THEN
1228 itest(i) = 1
1229 igwd = igwd + 1
1230 ENDIF
1231 ENDDO
1232
1233 CALL drag_noro(klon, llm, dtphys, paprs, play, zmea, zstd, zsig, zgam, &
1234 zthe, zpic, zval, itest, t_seri, u_seri, v_seri, zulow, zvlow, &
1235 zustrdr, zvstrdr, d_t_oro, d_u_oro, d_v_oro)
1236
1237 ! ajout des tendances
1238 DO k = 1, llm
1239 DO i = 1, klon
1240 t_seri(i, k) = t_seri(i, k) + d_t_oro(i, k)
1241 u_seri(i, k) = u_seri(i, k) + d_u_oro(i, k)
1242 v_seri(i, k) = v_seri(i, k) + d_v_oro(i, k)
1243 ENDDO
1244 ENDDO
1245 ENDIF
1246
1247 IF (ok_orolf) THEN
1248 ! S\'election des points pour lesquels le sch\'ema est actif :
1249 igwd = 0
1250 DO i = 1, klon
1251 itest(i) = 0
1252 IF (zpic(i) - zmea(i) > 100.) THEN
1253 itest(i) = 1
1254 igwd = igwd + 1
1255 ENDIF
1256 ENDDO
1257
1258 CALL lift_noro(klon, llm, dtphys, paprs, play, rlat, zmea, zstd, zpic, &
1259 itest, t_seri, u_seri, v_seri, zulow, zvlow, zustrli, zvstrli, &
1260 d_t_lif, d_u_lif, d_v_lif)
1261
1262 ! Ajout des tendances :
1263 DO k = 1, llm
1264 DO i = 1, klon
1265 t_seri(i, k) = t_seri(i, k) + d_t_lif(i, k)
1266 u_seri(i, k) = u_seri(i, k) + d_u_lif(i, k)
1267 v_seri(i, k) = v_seri(i, k) + d_v_lif(i, k)
1268 ENDDO
1269 ENDDO
1270 ENDIF
1271
1272 ! Stress n\'ecessaires : toute la physique
1273
1274 DO i = 1, klon
1275 zustrph(i) = 0.
1276 zvstrph(i) = 0.
1277 ENDDO
1278 DO k = 1, llm
1279 DO i = 1, klon
1280 zustrph(i) = zustrph(i) + (u_seri(i, k) - u(i, k)) / dtphys &
1281 * zmasse(i, k)
1282 zvstrph(i) = zvstrph(i) + (v_seri(i, k) - v(i, k)) / dtphys &
1283 * zmasse(i, k)
1284 ENDDO
1285 ENDDO
1286
1287 CALL aaam_bud(rg, romega, rlat, rlon, pphis, zustrdr, zustrli, zustrph, &
1288 zvstrdr, zvstrli, zvstrph, paprs, u, v, aam, torsfc)
1289
1290 IF (if_ebil >= 2) CALL diagetpq(airephy, 'after orography', ip_ebil, 2, &
1291 2, dtphys, t_seri, q_seri, ql_seri, u_seri, v_seri, paprs, d_h_vcol, &
1292 d_qt, d_ec)
1293
1294 ! Calcul des tendances traceurs
1295 call phytrac(itap, lmt_pas, julien, time, firstcal, lafin, dtphys, t, &
1296 paprs, play, mfu, mfd, pde_u, pen_d, ycoefh, fm_therm, entr_therm, &
1297 yu1, yv1, ftsol, pctsrf, frac_impa, frac_nucl, da, phi, mp, upwd, &
1298 dnwd, tr_seri, zmasse, ncid_startphy, nid_ins, itau_phy)
1299
1300 IF (offline) call phystokenc(dtphys, rlon, rlat, t, mfu, mfd, pen_u, &
1301 pde_u, pen_d, pde_d, fm_therm, entr_therm, ycoefh, yu1, yv1, ftsol, &
1302 pctsrf, frac_impa, frac_nucl, pphis, airephy, dtphys, itap)
1303
1304 ! Calculer le transport de l'eau et de l'energie (diagnostique)
1305 CALL transp(paprs, t_seri, q_seri, u_seri, v_seri, zphi, ve, vq, ue, uq)
1306
1307 ! diag. bilKP
1308
1309 CALL transp_lay(paprs, t_seri, q_seri, u_seri, v_seri, zphi, &
1310 ve_lay, vq_lay, ue_lay, uq_lay)
1311
1312 ! Accumuler les variables a stocker dans les fichiers histoire:
1313
1314 ! conversion Ec -> E thermique
1315 DO k = 1, llm
1316 DO i = 1, klon
1317 ZRCPD = RCPD * (1. + RVTMP2 * q_seri(i, k))
1318 d_t_ec(i, k) = 0.5 / ZRCPD &
1319 * (u(i, k)**2 + v(i, k)**2 - u_seri(i, k)**2 - v_seri(i, k)**2)
1320 t_seri(i, k) = t_seri(i, k) + d_t_ec(i, k)
1321 d_t_ec(i, k) = d_t_ec(i, k) / dtphys
1322 END DO
1323 END DO
1324
1325 IF (if_ebil >= 1) THEN
1326 tit = 'after physic'
1327 CALL diagetpq(airephy, tit, ip_ebil, 1, 1, dtphys, t_seri, q_seri, &
1328 ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
1329 ! Comme les tendances de la physique sont ajoute dans la dynamique,
1330 ! on devrait avoir que la variation d'entalpie par la dynamique
1331 ! est egale a la variation de la physique au pas de temps precedent.
1332 ! Donc la somme de ces 2 variations devrait etre nulle.
1333 call diagphy(airephy, tit, ip_ebil, topsw, toplw, solsw, sollw, sens, &
1334 evap, rain_fall, snow_fall, ztsol, d_h_vcol, d_qt, d_ec)
1335 d_h_vcol_phy = d_h_vcol
1336 END IF
1337
1338 ! SORTIES
1339
1340 ! prw = eau precipitable
1341 DO i = 1, klon
1342 prw(i) = 0.
1343 DO k = 1, llm
1344 prw(i) = prw(i) + q_seri(i, k)*zmasse(i, k)
1345 ENDDO
1346 ENDDO
1347
1348 ! Convertir les incrementations en tendances
1349
1350 DO k = 1, llm
1351 DO i = 1, klon
1352 d_u(i, k) = (u_seri(i, k) - u(i, k)) / dtphys
1353 d_v(i, k) = (v_seri(i, k) - v(i, k)) / dtphys
1354 d_t(i, k) = (t_seri(i, k) - t(i, k)) / dtphys
1355 d_qx(i, k, ivap) = (q_seri(i, k) - qx(i, k, ivap)) / dtphys
1356 d_qx(i, k, iliq) = (ql_seri(i, k) - qx(i, k, iliq)) / dtphys
1357 ENDDO
1358 ENDDO
1359
1360 DO iq = 3, nqmx
1361 DO k = 1, llm
1362 DO i = 1, klon
1363 d_qx(i, k, iq) = (tr_seri(i, k, iq - 2) - qx(i, k, iq)) / dtphys
1364 ENDDO
1365 ENDDO
1366 ENDDO
1367
1368 ! Sauvegarder les valeurs de t et q a la fin de la physique:
1369 DO k = 1, llm
1370 DO i = 1, klon
1371 t_ancien(i, k) = t_seri(i, k)
1372 q_ancien(i, k) = q_seri(i, k)
1373 ENDDO
1374 ENDDO
1375
1376 call write_histins
1377
1378 IF (lafin) then
1379 call NF95_CLOSE(ncid_startphy)
1380 CALL phyredem(pctsrf, ftsol, ftsoil, fqsurf, qsol, &
1381 fsnow, falbe, fevap, rain_fall, snow_fall, solsw, sollw, dlw, &
1382 radsol, frugs, agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, &
1383 t_ancien, q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, &
1384 w01)
1385 end IF
1386
1387 firstcal = .FALSE.
1388
1389 contains
1390
1391 subroutine write_histins
1392
1393 ! From phylmd/write_histins.h, version 1.2 2005/05/25 13:10:09
1394
1395 ! Ecriture des sorties
1396
1397 use dimens_m, only: iim, jjm
1398 USE histsync_m, ONLY: histsync
1399 USE histwrite_m, ONLY: histwrite
1400
1401 integer i, itau_w ! pas de temps ecriture
1402 REAL zx_tmp_2d(iim, jjm + 1), zx_tmp_3d(iim, jjm + 1, llm)
1403
1404 !--------------------------------------------------
1405
1406 IF (ok_instan) THEN
1407 ! Champs 2D:
1408
1409 itau_w = itau_phy + itap
1410
1411 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, pphis, zx_tmp_2d)
1412 CALL histwrite(nid_ins, "phis", itau_w, zx_tmp_2d)
1413
1414 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, airephy, zx_tmp_2d)
1415 CALL histwrite(nid_ins, "aire", itau_w, zx_tmp_2d)
1416
1417 DO i = 1, klon
1418 zx_tmp_fi2d(i) = paprs(i, 1)
1419 ENDDO
1420 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1421 CALL histwrite(nid_ins, "psol", itau_w, zx_tmp_2d)
1422
1423 DO i = 1, klon
1424 zx_tmp_fi2d(i) = rain_fall(i) + snow_fall(i)
1425 ENDDO
1426 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1427 CALL histwrite(nid_ins, "precip", itau_w, zx_tmp_2d)
1428
1429 DO i = 1, klon
1430 zx_tmp_fi2d(i) = rain_lsc(i) + snow_lsc(i)
1431 ENDDO
1432 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1433 CALL histwrite(nid_ins, "plul", itau_w, zx_tmp_2d)
1434
1435 DO i = 1, klon
1436 zx_tmp_fi2d(i) = rain_con(i) + snow_con(i)
1437 ENDDO
1438 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1439 CALL histwrite(nid_ins, "pluc", itau_w, zx_tmp_2d)
1440
1441 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zxtsol, zx_tmp_2d)
1442 CALL histwrite(nid_ins, "tsol", itau_w, zx_tmp_2d)
1443 !ccIM
1444 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zt2m, zx_tmp_2d)
1445 CALL histwrite(nid_ins, "t2m", itau_w, zx_tmp_2d)
1446
1447 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zq2m, zx_tmp_2d)
1448 CALL histwrite(nid_ins, "q2m", itau_w, zx_tmp_2d)
1449
1450 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zu10m, zx_tmp_2d)
1451 CALL histwrite(nid_ins, "u10m", itau_w, zx_tmp_2d)
1452
1453 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zv10m, zx_tmp_2d)
1454 CALL histwrite(nid_ins, "v10m", itau_w, zx_tmp_2d)
1455
1456 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, snow_fall, zx_tmp_2d)
1457 CALL histwrite(nid_ins, "snow", itau_w, zx_tmp_2d)
1458
1459 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, cdragm, zx_tmp_2d)
1460 CALL histwrite(nid_ins, "cdrm", itau_w, zx_tmp_2d)
1461
1462 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, cdragh, zx_tmp_2d)
1463 CALL histwrite(nid_ins, "cdrh", itau_w, zx_tmp_2d)
1464
1465 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, toplw, zx_tmp_2d)
1466 CALL histwrite(nid_ins, "topl", itau_w, zx_tmp_2d)
1467
1468 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, evap, zx_tmp_2d)
1469 CALL histwrite(nid_ins, "evap", itau_w, zx_tmp_2d)
1470
1471 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, solsw, zx_tmp_2d)
1472 CALL histwrite(nid_ins, "sols", itau_w, zx_tmp_2d)
1473
1474 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sollw, zx_tmp_2d)
1475 CALL histwrite(nid_ins, "soll", itau_w, zx_tmp_2d)
1476
1477 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sollwdown, zx_tmp_2d)
1478 CALL histwrite(nid_ins, "solldown", itau_w, zx_tmp_2d)
1479
1480 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, bils, zx_tmp_2d)
1481 CALL histwrite(nid_ins, "bils", itau_w, zx_tmp_2d)
1482
1483 zx_tmp_fi2d(1:klon) = - sens(1:klon)
1484 ! CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sens, zx_tmp_2d)
1485 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1486 CALL histwrite(nid_ins, "sens", itau_w, zx_tmp_2d)
1487
1488 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, fder, zx_tmp_2d)
1489 CALL histwrite(nid_ins, "fder", itau_w, zx_tmp_2d)
1490
1491 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_oce), zx_tmp_2d)
1492 CALL histwrite(nid_ins, "dtsvdfo", itau_w, zx_tmp_2d)
1493
1494 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_ter), zx_tmp_2d)
1495 CALL histwrite(nid_ins, "dtsvdft", itau_w, zx_tmp_2d)
1496
1497 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_lic), zx_tmp_2d)
1498 CALL histwrite(nid_ins, "dtsvdfg", itau_w, zx_tmp_2d)
1499
1500 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_sic), zx_tmp_2d)
1501 CALL histwrite(nid_ins, "dtsvdfi", itau_w, zx_tmp_2d)
1502
1503 DO nsrf = 1, nbsrf
1504 !XXX
1505 zx_tmp_fi2d(1 : klon) = pctsrf(1 : klon, nsrf)*100.
1506 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1507 CALL histwrite(nid_ins, "pourc_"//clnsurf(nsrf), itau_w, &
1508 zx_tmp_2d)
1509
1510 zx_tmp_fi2d(1 : klon) = pctsrf(1 : klon, nsrf)
1511 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1512 CALL histwrite(nid_ins, "fract_"//clnsurf(nsrf), itau_w, &
1513 zx_tmp_2d)
1514
1515 zx_tmp_fi2d(1 : klon) = fluxt(1 : klon, 1, nsrf)
1516 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1517 CALL histwrite(nid_ins, "sens_"//clnsurf(nsrf), itau_w, &
1518 zx_tmp_2d)
1519
1520 zx_tmp_fi2d(1 : klon) = fluxlat(1 : klon, nsrf)
1521 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1522 CALL histwrite(nid_ins, "lat_"//clnsurf(nsrf), itau_w, &
1523 zx_tmp_2d)
1524
1525 zx_tmp_fi2d(1 : klon) = ftsol(1 : klon, nsrf)
1526 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1527 CALL histwrite(nid_ins, "tsol_"//clnsurf(nsrf), itau_w, &
1528 zx_tmp_2d)
1529
1530 zx_tmp_fi2d(1 : klon) = fluxu(1 : klon, 1, nsrf)
1531 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1532 CALL histwrite(nid_ins, "taux_"//clnsurf(nsrf), itau_w, &
1533 zx_tmp_2d)
1534
1535 zx_tmp_fi2d(1 : klon) = fluxv(1 : klon, 1, nsrf)
1536 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1537 CALL histwrite(nid_ins, "tauy_"//clnsurf(nsrf), itau_w, &
1538 zx_tmp_2d)
1539
1540 zx_tmp_fi2d(1 : klon) = frugs(1 : klon, nsrf)
1541 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1542 CALL histwrite(nid_ins, "rugs_"//clnsurf(nsrf), itau_w, &
1543 zx_tmp_2d)
1544
1545 zx_tmp_fi2d(1 : klon) = falbe(:, nsrf)
1546 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1547 CALL histwrite(nid_ins, "albe_"//clnsurf(nsrf), itau_w, &
1548 zx_tmp_2d)
1549
1550 END DO
1551 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, albsol, zx_tmp_2d)
1552 CALL histwrite(nid_ins, "albs", itau_w, zx_tmp_2d)
1553
1554 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zxrugs, zx_tmp_2d)
1555 CALL histwrite(nid_ins, "rugs", itau_w, zx_tmp_2d)
1556
1557 !HBTM2
1558
1559 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_pblh, zx_tmp_2d)
1560 CALL histwrite(nid_ins, "s_pblh", itau_w, zx_tmp_2d)
1561
1562 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_pblt, zx_tmp_2d)
1563 CALL histwrite(nid_ins, "s_pblt", itau_w, zx_tmp_2d)
1564
1565 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_lcl, zx_tmp_2d)
1566 CALL histwrite(nid_ins, "s_lcl", itau_w, zx_tmp_2d)
1567
1568 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_capCL, zx_tmp_2d)
1569 CALL histwrite(nid_ins, "s_capCL", itau_w, zx_tmp_2d)
1570
1571 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_oliqCL, zx_tmp_2d)
1572 CALL histwrite(nid_ins, "s_oliqCL", itau_w, zx_tmp_2d)
1573
1574 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_cteiCL, zx_tmp_2d)
1575 CALL histwrite(nid_ins, "s_cteiCL", itau_w, zx_tmp_2d)
1576
1577 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_therm, zx_tmp_2d)
1578 CALL histwrite(nid_ins, "s_therm", itau_w, zx_tmp_2d)
1579
1580 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb1, zx_tmp_2d)
1581 CALL histwrite(nid_ins, "s_trmb1", itau_w, zx_tmp_2d)
1582
1583 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb2, zx_tmp_2d)
1584 CALL histwrite(nid_ins, "s_trmb2", itau_w, zx_tmp_2d)
1585
1586 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb3, zx_tmp_2d)
1587 CALL histwrite(nid_ins, "s_trmb3", itau_w, zx_tmp_2d)
1588
1589 ! Champs 3D:
1590
1591 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, t_seri, zx_tmp_3d)
1592 CALL histwrite(nid_ins, "temp", itau_w, zx_tmp_3d)
1593
1594 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, u_seri, zx_tmp_3d)
1595 CALL histwrite(nid_ins, "vitu", itau_w, zx_tmp_3d)
1596
1597 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, v_seri, zx_tmp_3d)
1598 CALL histwrite(nid_ins, "vitv", itau_w, zx_tmp_3d)
1599
1600 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, zphi, zx_tmp_3d)
1601 CALL histwrite(nid_ins, "geop", itau_w, zx_tmp_3d)
1602
1603 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, play, zx_tmp_3d)
1604 CALL histwrite(nid_ins, "pres", itau_w, zx_tmp_3d)
1605
1606 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, d_t_vdf, zx_tmp_3d)
1607 CALL histwrite(nid_ins, "dtvdf", itau_w, zx_tmp_3d)
1608
1609 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, d_q_vdf, zx_tmp_3d)
1610 CALL histwrite(nid_ins, "dqvdf", itau_w, zx_tmp_3d)
1611
1612 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, zx_rh, zx_tmp_3d)
1613 CALL histwrite(nid_ins, "rhum", itau_w, zx_tmp_3d)
1614
1615 call histsync(nid_ins)
1616 ENDIF
1617
1618 end subroutine write_histins
1619
1620 END SUBROUTINE physiq
1621
1622 end module physiq_m

  ViewVC Help
Powered by ViewVC 1.1.21