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

Contents of /trunk/phylmd/physiq.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 188 - (show annotations)
Tue Mar 22 16:31:39 2016 UTC (8 years, 2 months ago) by guez
Original Path: trunk/Sources/phylmd/physiq.f
File size: 57315 byte(s)
Removed argument ncum of cv30_unsat, arguments nloc, ncum, nd, na of cv30_yield.

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

  ViewVC Help
Powered by ViewVC 1.1.21