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

Contents of /trunk/phylmd/physiq.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 189 - (show annotations)
Tue Mar 29 15:20:23 2016 UTC (8 years, 1 month ago) by guez
Original Path: trunk/Sources/phylmd/physiq.f
File size: 55611 byte(s)
There was a function gr_phy_write_3d in dyn3d and a function
gr_phy_write_2d in module grid_change. Moved them into a new module
gr_phy_write_m under a generic interface gr_phy_write. Replaced calls
to gr_fi_ecrit by calls to gr_phy_write.

Removed arguments len, nloc and nd of cv30_compress.

Removed arguments wd and wd1 of cv30_uncompress, wd of cv30_yield, wd
of concvl, wd1 of cv_driver. Was just filled with 0. Removed option
ok_gust in physiq, never used.

In cv30_unsat, cv30_yield and cv_driver, we only need to define b to
level nl - 1.

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

  ViewVC Help
Powered by ViewVC 1.1.21