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

Contents of /trunk/Sources/phylmd/physiq.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 226 - (show annotations)
Mon Oct 16 13:04:05 2017 UTC (4 months ago) by guez
File size: 39032 byte(s)
In clmain, u1lay and v1lay do not depend on the subsurface. So output
values of zu1 and zv1 are simply u(:, 1) and v(:, 1). Just remove the
corresponding computations from clmain and define yu1 and yv1 in physiq.

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 ajsec_m, only: ajsec
20 use calltherm_m, only: calltherm
21 USE clesphys, ONLY: cdhmax, cdmmax, ecrit_ins, ksta, ksta_ter, ok_kzmin, &
22 ok_instan
23 USE clesphys2, ONLY: conv_emanuel, nbapp_rad, new_oliq, ok_orodr, ok_orolf
24 USE clmain_m, ONLY: clmain
25 use clouds_gno_m, only: clouds_gno
26 use comconst, only: dtphys
27 USE comgeomphy, ONLY: airephy
28 USE concvl_m, ONLY: concvl
29 USE conf_gcm_m, ONLY: lmt_pas
30 USE conf_phys_m, ONLY: conf_phys
31 use conflx_m, only: conflx
32 USE ctherm, ONLY: iflag_thermals, nsplit_thermals
33 use diagcld2_m, only: diagcld2
34 USE dimens_m, ONLY: llm, nqmx
35 USE dimphy, ONLY: klon
36 USE dimsoil, ONLY: nsoilmx
37 use drag_noro_m, only: drag_noro
38 use dynetat0_m, only: day_ref, annee_ref
39 USE fcttre, ONLY: foeew
40 use fisrtilp_m, only: fisrtilp
41 USE hgardfou_m, ONLY: hgardfou
42 USE histsync_m, ONLY: histsync
43 USE histwrite_phy_m, ONLY: histwrite_phy
44 USE indicesol, ONLY: clnsurf, epsfra, is_lic, is_oce, is_sic, is_ter, &
45 nbsrf
46 USE ini_histins_m, ONLY: ini_histins, nid_ins
47 use netcdf95, only: NF95_CLOSE
48 use newmicro_m, only: newmicro
49 use nr_util, only: assert
50 use nuage_m, only: nuage
51 USE orbite_m, ONLY: orbite
52 USE ozonecm_m, ONLY: ozonecm
53 USE phyetat0_m, ONLY: phyetat0, rlat, rlon
54 USE phyredem_m, ONLY: phyredem
55 USE phyredem0_m, ONLY: phyredem0
56 USE phytrac_m, ONLY: phytrac
57 use radlwsw_m, only: radlwsw
58 use yoegwd, only: sugwd
59 USE suphec_m, ONLY: rcpd, retv, rg, rlvtt, romega, rsigma, rtt, rmo3, md
60 use time_phylmdz, only: itap, increment_itap
61 use transp_m, only: transp
62 use transp_lay_m, only: transp_lay
63 use unit_nml_m, only: unit_nml
64 USE ymds2ju_m, ONLY: ymds2ju
65 USE yoethf_m, ONLY: r2es, rvtmp2
66 use zenang_m, only: zenang
67
68 logical, intent(in):: lafin ! dernier passage
69
70 integer, intent(in):: dayvrai
71 ! current day number, based at value 1 on January 1st of annee_ref
72
73 REAL, intent(in):: time ! heure de la journ\'ee en fraction de jour
74
75 REAL, intent(in):: paprs(:, :) ! (klon, llm + 1)
76 ! pression pour chaque inter-couche, en Pa
77
78 REAL, intent(in):: play(:, :) ! (klon, llm)
79 ! pression pour le mileu de chaque couche (en Pa)
80
81 REAL, intent(in):: pphi(:, :) ! (klon, llm)
82 ! géopotentiel de chaque couche (référence sol)
83
84 REAL, intent(in):: pphis(:) ! (klon) géopotentiel du sol
85
86 REAL, intent(in):: u(:, :) ! (klon, llm)
87 ! vitesse dans la direction X (de O a E) en m / s
88
89 REAL, intent(in):: v(:, :) ! (klon, llm) vitesse Y (de S a N) en m / s
90 REAL, intent(in):: t(:, :) ! (klon, llm) temperature (K)
91
92 REAL, intent(in):: qx(:, :, :) ! (klon, llm, nqmx)
93 ! (humidit\'e sp\'ecifique et fractions massiques des autres traceurs)
94
95 REAL, intent(in):: omega(:, :) ! (klon, llm) vitesse verticale en Pa / s
96 REAL, intent(out):: d_u(:, :) ! (klon, llm) tendance physique de "u" (m s-2)
97 REAL, intent(out):: d_v(:, :) ! (klon, llm) tendance physique de "v" (m s-2)
98 REAL, intent(out):: d_t(:, :) ! (klon, llm) tendance physique de "t" (K / s)
99
100 REAL, intent(out):: d_qx(:, :, :) ! (klon, llm, nqmx)
101 ! tendance physique de "qx" (s-1)
102
103 ! Local:
104
105 LOGICAL:: firstcal = .true.
106
107 LOGICAL, PARAMETER:: ok_stratus = .FALSE.
108 ! Ajouter artificiellement les stratus
109
110 ! pour phystoke avec thermiques
111 REAL fm_therm(klon, llm + 1)
112 REAL entr_therm(klon, llm)
113 real, save:: q2(klon, llm + 1, nbsrf)
114
115 INTEGER, PARAMETER:: ivap = 1 ! indice de traceur pour vapeur d'eau
116 INTEGER, PARAMETER:: iliq = 2 ! indice de traceur pour eau liquide
117
118 REAL, save:: t_ancien(klon, llm), q_ancien(klon, llm)
119 LOGICAL, save:: ancien_ok
120
121 REAL d_t_dyn(klon, llm) ! tendance dynamique pour "t" (K / s)
122 REAL d_q_dyn(klon, llm) ! tendance dynamique pour "q" (kg / kg / s)
123
124 real da(klon, llm), phi(klon, llm, llm), mp(klon, llm)
125
126 REAL, save:: swdn0(klon, llm + 1), swdn(klon, llm + 1)
127 REAL, save:: swup0(klon, llm + 1), swup(klon, llm + 1)
128
129 REAL, save:: lwdn0(klon, llm + 1), lwdn(klon, llm + 1)
130 REAL, save:: lwup0(klon, llm + 1), lwup(klon, llm + 1)
131
132 ! prw: precipitable water
133 real prw(klon)
134
135 ! flwp, fiwp = Liquid Water Path & Ice Water Path (kg / m2)
136 ! flwc, fiwc = Liquid Water Content & Ice Water Content (kg / kg)
137 REAL flwp(klon), fiwp(klon)
138 REAL flwc(klon, llm), fiwc(klon, llm)
139
140 ! Variables propres a la physique
141
142 INTEGER, save:: radpas
143 ! Radiative transfer computations are made every "radpas" call to
144 ! "physiq".
145
146 REAL, save:: radsol(klon) ! bilan radiatif au sol calcule par code radiatif
147 REAL, save:: ftsol(klon, nbsrf) ! skin temperature of surface fraction
148
149 REAL, save:: ftsoil(klon, nsoilmx, nbsrf)
150 ! soil temperature of surface fraction
151
152 REAL, save:: fevap(klon, nbsrf) ! evaporation
153 REAL fluxlat(klon, nbsrf)
154
155 REAL, save:: fqsurf(klon, nbsrf)
156 ! humidite de l'air au contact de la surface
157
158 REAL, save:: qsol(klon) ! column-density of water in soil, in kg m-2
159 REAL, save:: fsnow(klon, nbsrf) ! \'epaisseur neigeuse
160 REAL, save:: falbe(klon, nbsrf) ! albedo visible par type de surface
161
162 ! Param\`etres de l'orographie \`a l'\'echelle sous-maille (OESM) :
163 REAL, save:: zmea(klon) ! orographie moyenne
164 REAL, save:: zstd(klon) ! deviation standard de l'OESM
165 REAL, save:: zsig(klon) ! pente de l'OESM
166 REAL, save:: zgam(klon) ! anisotropie de l'OESM
167 REAL, save:: zthe(klon) ! orientation de l'OESM
168 REAL, save:: zpic(klon) ! Maximum de l'OESM
169 REAL, save:: zval(klon) ! Minimum de l'OESM
170 REAL, save:: rugoro(klon) ! longueur de rugosite de l'OESM
171 REAL zulow(klon), zvlow(klon)
172 INTEGER igwd, itest(klon)
173
174 REAL, save:: agesno(klon, nbsrf) ! age de la neige
175 REAL, save:: run_off_lic_0(klon)
176
177 ! Variables li\'ees \`a la convection d'Emanuel :
178 REAL, save:: Ma(klon, llm) ! undilute upward mass flux
179 REAL, save:: qcondc(klon, llm) ! in-cld water content from convect
180 REAL, save:: sig1(klon, llm), w01(klon, llm)
181
182 ! Variables pour la couche limite (Alain Lahellec) :
183 REAL cdragh(klon) ! drag coefficient pour T and Q
184 REAL cdragm(klon) ! drag coefficient pour vent
185
186 ! Pour phytrac :
187 REAL ycoefh(klon, llm) ! coef d'echange pour phytrac
188 REAL yu1(klon), yv1(klon) ! vent dans la premi\`ere couche
189
190 REAL, save:: ffonte(klon, nbsrf)
191 ! flux thermique utilise pour fondre la neige
192
193 REAL, save:: fqcalving(klon, nbsrf)
194 ! flux d'eau "perdue" par la surface et necessaire pour limiter la
195 ! hauteur de neige, en kg / m2 / s
196
197 REAL zxffonte(klon), zxfqcalving(klon)
198
199 REAL, save:: pfrac_impa(klon, llm)! Produits des coefs lessivage impaction
200 REAL, save:: pfrac_nucl(klon, llm)! Produits des coefs lessivage nucleation
201
202 REAL, save:: pfrac_1nucl(klon, llm)
203 ! Produits des coefs lessi nucl (alpha = 1)
204
205 REAL frac_impa(klon, llm) ! fraction d'a\'erosols lessiv\'es (impaction)
206 REAL frac_nucl(klon, llm) ! idem (nucleation)
207
208 REAL, save:: rain_fall(klon)
209 ! liquid water mass flux (kg / m2 / s), positive down
210
211 REAL, save:: snow_fall(klon)
212 ! solid water mass flux (kg / m2 / s), positive down
213
214 REAL rain_tiedtke(klon), snow_tiedtke(klon)
215
216 REAL evap(klon) ! flux d'\'evaporation au sol
217 real devap(klon) ! derivative of the evaporation flux at the surface
218 REAL sens(klon) ! flux de chaleur sensible au sol
219 real dsens(klon) ! derivee du flux de chaleur sensible au sol
220 REAL, save:: dlw(klon) ! derivative of infra-red flux
221 REAL bils(klon) ! bilan de chaleur au sol
222 REAL fder(klon) ! Derive de flux (sensible et latente)
223 REAL ve(klon) ! integr. verticale du transport meri. de l'energie
224 REAL vq(klon) ! integr. verticale du transport meri. de l'eau
225 REAL ue(klon) ! integr. verticale du transport zonal de l'energie
226 REAL uq(klon) ! integr. verticale du transport zonal de l'eau
227
228 REAL, save:: frugs(klon, nbsrf) ! longueur de rugosite
229 REAL zxrugs(klon) ! longueur de rugosite
230
231 ! Conditions aux limites
232
233 INTEGER julien
234 REAL, save:: pctsrf(klon, nbsrf) ! percentage of surface
235 REAL, save:: albsol(klon) ! albedo du sol total, visible, moyen par maille
236 REAL, SAVE:: wo(klon, llm) ! column density of ozone in a cell, in kDU
237 real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2
238
239 real, save:: clwcon(klon, llm), rnebcon(klon, llm)
240 real, save:: clwcon0(klon, llm), rnebcon0(klon, llm)
241
242 REAL rhcl(klon, llm) ! humiditi relative ciel clair
243 REAL dialiq(klon, llm) ! eau liquide nuageuse
244 REAL diafra(klon, llm) ! fraction nuageuse
245 REAL cldliq(klon, llm) ! eau liquide nuageuse
246 REAL cldfra(klon, llm) ! fraction nuageuse
247 REAL cldtau(klon, llm) ! epaisseur optique
248 REAL cldemi(klon, llm) ! emissivite infrarouge
249
250 REAL flux_q(klon, nbsrf) ! flux turbulent d'humidite à la surface
251 REAL flux_t(klon, nbsrf) ! flux turbulent de chaleur à la surface
252 REAL flux_u(klon, nbsrf) ! flux turbulent de vitesse u à la surface
253 REAL flux_v(klon, nbsrf) ! flux turbulent de vitesse v à la surface
254
255 ! Le rayonnement n'est pas calcul\'e tous les pas, il faut donc que
256 ! les variables soient r\'emanentes.
257 REAL, save:: heat(klon, llm) ! chauffage solaire
258 REAL, save:: heat0(klon, llm) ! chauffage solaire ciel clair
259 REAL, save:: cool(klon, llm) ! refroidissement infrarouge
260 REAL, save:: cool0(klon, llm) ! refroidissement infrarouge ciel clair
261 REAL, save:: topsw(klon), toplw(klon), solsw(klon)
262 REAL, save:: sollw(klon) ! rayonnement infrarouge montant \`a la surface
263 real, save:: sollwdown(klon) ! downward LW flux at surface
264 REAL, save:: topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)
265 REAL, save:: albpla(klon)
266 REAL fsollw(klon, nbsrf) ! bilan flux IR pour chaque sous-surface
267 REAL fsolsw(klon, nbsrf) ! flux solaire absorb\'e pour chaque sous-surface
268
269 REAL conv_q(klon, llm) ! convergence de l'humidite (kg / kg / s)
270 REAL conv_t(klon, llm) ! convergence of temperature (K / s)
271
272 REAL cldl(klon), cldm(klon), cldh(klon) ! nuages bas, moyen et haut
273 REAL cldt(klon), cldq(klon) ! nuage total, eau liquide integree
274
275 REAL zxfluxlat(klon)
276 REAL dist, mu0(klon), fract(klon)
277 real longi
278 REAL z_avant(klon), z_apres(klon), z_factor(klon)
279 REAL zb
280 REAL zx_t, zx_qs, zcor
281 real zqsat(klon, llm)
282 INTEGER i, k, iq, nsrf
283 REAL zphi(klon, llm)
284
285 ! cf. Anne Mathieu, variables pour la couche limite atmosphérique (hbtm)
286
287 REAL, SAVE:: pblh(klon, nbsrf) ! Hauteur de couche limite
288 REAL, SAVE:: plcl(klon, nbsrf) ! Niveau de condensation de la CLA
289 REAL, SAVE:: capCL(klon, nbsrf) ! CAPE de couche limite
290 REAL, SAVE:: oliqCL(klon, nbsrf) ! eau_liqu integree de couche limite
291 REAL, SAVE:: cteiCL(klon, nbsrf) ! cloud top instab. crit. couche limite
292 REAL, SAVE:: pblt(klon, nbsrf) ! T \`a la hauteur de couche limite
293 REAL, SAVE:: therm(klon, nbsrf)
294 REAL, SAVE:: trmb1(klon, nbsrf) ! deep_cape
295 REAL, SAVE:: trmb2(klon, nbsrf) ! inhibition
296 REAL, SAVE:: trmb3(klon, nbsrf) ! Point Omega
297 ! Grandeurs de sorties
298 REAL s_pblh(klon), s_lcl(klon), s_capCL(klon)
299 REAL s_oliqCL(klon), s_cteiCL(klon), s_pblt(klon)
300 REAL s_therm(klon), s_trmb1(klon), s_trmb2(klon)
301 REAL s_trmb3(klon)
302
303 ! Variables pour la convection de K. Emanuel :
304
305 REAL upwd(klon, llm) ! saturated updraft mass flux
306 REAL dnwd(klon, llm) ! saturated downdraft mass flux
307 REAL, save:: cape(klon)
308
309 INTEGER iflagctrl(klon) ! flag fonctionnement de convect
310
311 ! Variables du changement
312
313 ! con: convection
314 ! lsc: large scale condensation
315 ! ajs: ajustement sec
316 ! eva: \'evaporation de l'eau liquide nuageuse
317 ! vdf: vertical diffusion in boundary layer
318 REAL d_t_con(klon, llm), d_q_con(klon, llm)
319 REAL, save:: d_u_con(klon, llm), d_v_con(klon, llm)
320 REAL d_t_lsc(klon, llm), d_q_lsc(klon, llm), d_ql_lsc(klon, llm)
321 REAL d_t_ajs(klon, llm), d_q_ajs(klon, llm)
322 REAL d_u_ajs(klon, llm), d_v_ajs(klon, llm)
323 REAL rneb(klon, llm)
324
325 REAL mfu(klon, llm), mfd(klon, llm)
326 REAL pen_u(klon, llm), pen_d(klon, llm)
327 REAL pde_u(klon, llm), pde_d(klon, llm)
328 INTEGER kcbot(klon), kctop(klon), kdtop(klon)
329 REAL pmflxr(klon, llm + 1), pmflxs(klon, llm + 1)
330 REAL prfl(klon, llm + 1), psfl(klon, llm + 1)
331
332 INTEGER, save:: ibas_con(klon), itop_con(klon)
333 real ema_pct(klon) ! Emanuel pressure at cloud top, in Pa
334
335 REAL, save:: rain_con(klon)
336 real rain_lsc(klon)
337 REAL, save:: snow_con(klon) ! neige (mm / s)
338 real snow_lsc(klon)
339 REAL d_ts(klon, nbsrf) ! variation of ftsol
340
341 REAL d_u_vdf(klon, llm), d_v_vdf(klon, llm)
342 REAL d_t_vdf(klon, llm), d_q_vdf(klon, llm)
343
344 REAL d_u_oro(klon, llm), d_v_oro(klon, llm)
345 REAL d_t_oro(klon, llm)
346 REAL d_u_lif(klon, llm), d_v_lif(klon, llm)
347 REAL d_t_lif(klon, llm)
348
349 REAL, save:: ratqs(klon, llm)
350 real ratqss(klon, llm), ratqsc(klon, llm)
351 real:: ratqsbas = 0.01, ratqshaut = 0.3
352
353 ! Parametres lies au nouveau schema de nuages (SB, PDF)
354 real:: fact_cldcon = 0.375
355 real:: facttemps = 1.e-4
356 logical:: ok_newmicro = .true.
357 real facteur
358
359 integer:: iflag_cldcon = 1
360 logical ptconv(klon, llm)
361
362 ! Variables pour effectuer les appels en s\'erie :
363
364 REAL t_seri(klon, llm), q_seri(klon, llm)
365 REAL ql_seri(klon, llm)
366 REAL u_seri(klon, llm), v_seri(klon, llm)
367 REAL tr_seri(klon, llm, nqmx - 2)
368
369 REAL zx_rh(klon, llm)
370
371 REAL zustrdr(klon), zvstrdr(klon)
372 REAL zustrli(klon), zvstrli(klon)
373 REAL zustrph(klon), zvstrph(klon)
374 REAL aam, torsfc
375
376 REAL ve_lay(klon, llm) ! transport meri. de l'energie a chaque niveau vert.
377 REAL vq_lay(klon, llm) ! transport meri. de l'eau a chaque niveau vert.
378 REAL ue_lay(klon, llm) ! transport zonal de l'energie a chaque niveau vert.
379 REAL uq_lay(klon, llm) ! transport zonal de l'eau a chaque niveau vert.
380
381 real date0
382 REAL tsol(klon)
383
384 REAL d_t_ec(klon, llm)
385 ! tendance due \`a la conversion d'\'energie cin\'etique en
386 ! énergie thermique
387
388 REAL, save:: t2m(klon, nbsrf), q2m(klon, nbsrf)
389 ! temperature and humidity at 2 m
390
391 REAL, save:: u10m_srf(klon, nbsrf), v10m_srf(klon, nbsrf)
392 ! composantes du vent \`a 10 m
393
394 REAL zt2m(klon), zq2m(klon) ! température, humidité 2 m moyenne sur 1 maille
395 REAL u10m(klon), v10m(klon) ! vent \`a 10 m moyenn\' sur les sous-surfaces
396
397 ! Aerosol effects:
398
399 REAL, save:: topswad(klon), solswad(klon) ! aerosol direct effect
400 LOGICAL:: ok_ade = .false. ! apply aerosol direct effect
401
402 REAL:: bl95_b0 = 2., bl95_b1 = 0.2
403 ! Parameters in equation (D) of Boucher and Lohmann (1995, Tellus
404 ! B). They link cloud droplet number concentration to aerosol mass
405 ! concentration.
406
407 real zmasse(klon, llm)
408 ! (column-density of mass of air in a cell, in kg m-2)
409
410 integer, save:: ncid_startphy
411
412 namelist /physiq_nml/ fact_cldcon, facttemps, ok_newmicro, iflag_cldcon, &
413 ratqsbas, ratqshaut, ok_ade, bl95_b0, bl95_b1, iflag_thermals, &
414 nsplit_thermals
415
416 !----------------------------------------------------------------
417
418 IF (nqmx < 2) CALL abort_gcm('physiq', &
419 'eaux vapeur et liquide sont indispensables')
420
421 test_firstcal: IF (firstcal) THEN
422 ! initialiser
423 u10m_srf = 0.
424 v10m_srf = 0.
425 t2m = 0.
426 q2m = 0.
427 ffonte = 0.
428 fqcalving = 0.
429 rain_con = 0.
430 snow_con = 0.
431 d_u_con = 0.
432 d_v_con = 0.
433 rnebcon0 = 0.
434 clwcon0 = 0.
435 rnebcon = 0.
436 clwcon = 0.
437 pblh =0. ! Hauteur de couche limite
438 plcl =0. ! Niveau de condensation de la CLA
439 capCL =0. ! CAPE de couche limite
440 oliqCL =0. ! eau_liqu integree de couche limite
441 cteiCL =0. ! cloud top instab. crit. couche limite
442 pblt =0.
443 therm =0.
444 trmb1 =0. ! deep_cape
445 trmb2 =0. ! inhibition
446 trmb3 =0. ! Point Omega
447
448 iflag_thermals = 0
449 nsplit_thermals = 1
450 print *, "Enter namelist 'physiq_nml'."
451 read(unit=*, nml=physiq_nml)
452 write(unit_nml, nml=physiq_nml)
453
454 call conf_phys
455
456 ! Initialiser les compteurs:
457
458 frugs = 0.
459 CALL phyetat0(pctsrf, ftsol, ftsoil, fqsurf, qsol, fsnow, falbe, &
460 fevap, rain_fall, snow_fall, solsw, sollw, dlw, radsol, frugs, &
461 agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, t_ancien, &
462 q_ancien, ancien_ok, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, &
463 w01, ncid_startphy)
464
465 ! ATTENTION : il faudra a terme relire q2 dans l'etat initial
466 q2 = 1e-8
467
468 radpas = lmt_pas / nbapp_rad
469 print *, "radpas = ", radpas
470
471 ! Initialisation pour le sch\'ema de convection d'Emanuel :
472 IF (conv_emanuel) THEN
473 ibas_con = 1
474 itop_con = 1
475 ENDIF
476
477 IF (ok_orodr) THEN
478 rugoro = MAX(1e-5, zstd * zsig / 2)
479 CALL SUGWD(paprs, play)
480 else
481 rugoro = 0.
482 ENDIF
483
484 ecrit_ins = NINT(ecrit_ins / dtphys)
485
486 ! Initialisation des sorties
487
488 call ini_histins(dtphys, ok_newmicro)
489 CALL ymds2ju(annee_ref, 1, day_ref, 0., date0)
490 ! Positionner date0 pour initialisation de ORCHIDEE
491 print *, 'physiq date0: ', date0
492 CALL phyredem0
493 ENDIF test_firstcal
494
495 ! We will modify variables *_seri and we will not touch variables
496 ! u, v, t, qx:
497 t_seri = t
498 u_seri = u
499 v_seri = v
500 q_seri = qx(:, :, ivap)
501 ql_seri = qx(:, :, iliq)
502 tr_seri = qx(:, :, 3:nqmx)
503
504 tsol = sum(ftsol * pctsrf, dim = 2)
505
506 ! Diagnostic de la tendance dynamique :
507 IF (ancien_ok) THEN
508 DO k = 1, llm
509 DO i = 1, klon
510 d_t_dyn(i, k) = (t_seri(i, k) - t_ancien(i, k)) / dtphys
511 d_q_dyn(i, k) = (q_seri(i, k) - q_ancien(i, k)) / dtphys
512 ENDDO
513 ENDDO
514 ELSE
515 DO k = 1, llm
516 DO i = 1, klon
517 d_t_dyn(i, k) = 0.
518 d_q_dyn(i, k) = 0.
519 ENDDO
520 ENDDO
521 ancien_ok = .TRUE.
522 ENDIF
523
524 ! Ajouter le geopotentiel du sol:
525 DO k = 1, llm
526 DO i = 1, klon
527 zphi(i, k) = pphi(i, k) + pphis(i)
528 ENDDO
529 ENDDO
530
531 ! Check temperatures:
532 CALL hgardfou(t_seri, ftsol)
533
534 call increment_itap
535 julien = MOD(dayvrai, 360)
536 if (julien == 0) julien = 360
537
538 forall (k = 1: llm) zmasse(:, k) = (paprs(:, k) - paprs(:, k + 1)) / rg
539
540 ! \'Evaporation de l'eau liquide nuageuse :
541 DO k = 1, llm
542 DO i = 1, klon
543 zb = MAX(0., ql_seri(i, k))
544 t_seri(i, k) = t_seri(i, k) &
545 - zb * RLVTT / RCPD / (1. + RVTMP2 * q_seri(i, k))
546 q_seri(i, k) = q_seri(i, k) + zb
547 ENDDO
548 ENDDO
549 ql_seri = 0.
550
551 frugs = MAX(frugs, 0.000015)
552 zxrugs = sum(frugs * pctsrf, dim = 2)
553
554 ! Calculs n\'ecessaires au calcul de l'albedo dans l'interface avec
555 ! la surface.
556
557 CALL orbite(REAL(julien), longi, dist)
558 CALL zenang(longi, time, dtphys * radpas, mu0, fract)
559 albsol = sum(falbe * pctsrf, dim = 2)
560
561 ! R\'epartition sous maille des flux longwave et shortwave
562 ! R\'epartition du longwave par sous-surface lin\'earis\'ee
563
564 forall (nsrf = 1: nbsrf)
565 fsollw(:, nsrf) = sollw + 4. * RSIGMA * tsol**3 &
566 * (tsol - ftsol(:, nsrf))
567 fsolsw(:, nsrf) = solsw * (1. - falbe(:, nsrf)) / (1. - albsol)
568 END forall
569
570 CALL clmain(dtphys, pctsrf, t_seri, q_seri, u_seri, v_seri, julien, mu0, &
571 ftsol, cdmmax, cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, qsol, &
572 paprs, play, fsnow, fqsurf, fevap, falbe, fluxlat, rain_fall, &
573 snow_fall, fsolsw, fsollw, frugs, agesno, rugoro, d_t_vdf, d_q_vdf, &
574 d_u_vdf, d_v_vdf, d_ts, flux_t, flux_q, flux_u, flux_v, cdragh, &
575 cdragm, q2, dsens, devap, ycoefh, t2m, q2m, u10m_srf, v10m_srf, &
576 pblh, capCL, oliqCL, cteiCL, pblT, therm, trmb1, trmb2, trmb3, plcl, &
577 fqcalving, ffonte, run_off_lic_0)
578
579 ! Incr\'ementation des flux
580
581 sens = - sum(flux_t * pctsrf, dim = 2)
582 evap = - sum(flux_q * pctsrf, dim = 2)
583 fder = dlw + dsens + devap
584
585 yu1 = u_seri(:, 1)
586 yv1 = v_seri(:, 1)
587
588 DO k = 1, llm
589 DO i = 1, klon
590 t_seri(i, k) = t_seri(i, k) + d_t_vdf(i, k)
591 q_seri(i, k) = q_seri(i, k) + d_q_vdf(i, k)
592 u_seri(i, k) = u_seri(i, k) + d_u_vdf(i, k)
593 v_seri(i, k) = v_seri(i, k) + d_v_vdf(i, k)
594 ENDDO
595 ENDDO
596
597 ! Update surface temperature:
598
599 call assert(abs(sum(pctsrf, dim = 2) - 1.) <= EPSFRA, 'physiq: pctsrf')
600 ftsol = ftsol + d_ts
601 tsol = sum(ftsol * pctsrf, dim = 2)
602 zxfluxlat = sum(fluxlat * pctsrf, dim = 2)
603 zt2m = sum(t2m * pctsrf, dim = 2)
604 zq2m = sum(q2m * pctsrf, dim = 2)
605 u10m = sum(u10m_srf * pctsrf, dim = 2)
606 v10m = sum(v10m_srf * pctsrf, dim = 2)
607 zxffonte = sum(ffonte * pctsrf, dim = 2)
608 zxfqcalving = sum(fqcalving * pctsrf, dim = 2)
609 s_pblh = sum(pblh * pctsrf, dim = 2)
610 s_lcl = sum(plcl * pctsrf, dim = 2)
611 s_capCL = sum(capCL * pctsrf, dim = 2)
612 s_oliqCL = sum(oliqCL * pctsrf, dim = 2)
613 s_cteiCL = sum(cteiCL * pctsrf, dim = 2)
614 s_pblT = sum(pblT * pctsrf, dim = 2)
615 s_therm = sum(therm * pctsrf, dim = 2)
616 s_trmb1 = sum(trmb1 * pctsrf, dim = 2)
617 s_trmb2 = sum(trmb2 * pctsrf, dim = 2)
618 s_trmb3 = sum(trmb3 * pctsrf, dim = 2)
619
620 ! Si une sous-fraction n'existe pas, elle prend la valeur moyenne :
621 DO nsrf = 1, nbsrf
622 DO i = 1, klon
623 IF (pctsrf(i, nsrf) < epsfra) then
624 ftsol(i, nsrf) = tsol(i)
625 t2m(i, nsrf) = zt2m(i)
626 q2m(i, nsrf) = zq2m(i)
627 u10m_srf(i, nsrf) = u10m(i)
628 v10m_srf(i, nsrf) = v10m(i)
629 ffonte(i, nsrf) = zxffonte(i)
630 fqcalving(i, nsrf) = zxfqcalving(i)
631 pblh(i, nsrf) = s_pblh(i)
632 plcl(i, nsrf) = s_lcl(i)
633 capCL(i, nsrf) = s_capCL(i)
634 oliqCL(i, nsrf) = s_oliqCL(i)
635 cteiCL(i, nsrf) = s_cteiCL(i)
636 pblT(i, nsrf) = s_pblT(i)
637 therm(i, nsrf) = s_therm(i)
638 trmb1(i, nsrf) = s_trmb1(i)
639 trmb2(i, nsrf) = s_trmb2(i)
640 trmb3(i, nsrf) = s_trmb3(i)
641 end IF
642 ENDDO
643 ENDDO
644
645 dlw = - 4. * RSIGMA * tsol**3
646
647 ! Appeler la convection
648
649 if (conv_emanuel) then
650 CALL concvl(paprs, play, t_seri, q_seri, u_seri, v_seri, sig1, w01, &
651 d_t_con, d_q_con, d_u_con, d_v_con, rain_con, ibas_con, itop_con, &
652 upwd, dnwd, Ma, cape, iflagctrl, qcondc, pmflxr, da, phi, mp)
653 snow_con = 0.
654 clwcon0 = qcondc
655 mfu = upwd + dnwd
656
657 zqsat = MIN(0.5, r2es * FOEEW(t_seri, rtt >= t_seri) / play)
658 zqsat = zqsat / (1. - retv * zqsat)
659
660 ! Properties of convective clouds
661 clwcon0 = fact_cldcon * clwcon0
662 call clouds_gno(klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, &
663 rnebcon0)
664
665 forall (i = 1:klon) ema_pct(i) = paprs(i, itop_con(i) + 1)
666 mfd = 0.
667 pen_u = 0.
668 pen_d = 0.
669 pde_d = 0.
670 pde_u = 0.
671 else
672 conv_q = d_q_dyn + d_q_vdf / dtphys
673 conv_t = d_t_dyn + d_t_vdf / dtphys
674 z_avant = sum((q_seri + ql_seri) * zmasse, dim=2)
675 CALL conflx(dtphys, paprs, play, t_seri(:, llm:1:- 1), &
676 q_seri(:, llm:1:- 1), conv_t, conv_q, - evap, omega, &
677 d_t_con, d_q_con, rain_con, snow_con, mfu(:, llm:1:- 1), &
678 mfd(:, llm:1:- 1), pen_u, pde_u, pen_d, pde_d, kcbot, kctop, &
679 kdtop, pmflxr, pmflxs)
680 WHERE (rain_con < 0.) rain_con = 0.
681 WHERE (snow_con < 0.) snow_con = 0.
682 ibas_con = llm + 1 - kcbot
683 itop_con = llm + 1 - kctop
684 END if
685
686 DO k = 1, llm
687 DO i = 1, klon
688 t_seri(i, k) = t_seri(i, k) + d_t_con(i, k)
689 q_seri(i, k) = q_seri(i, k) + d_q_con(i, k)
690 u_seri(i, k) = u_seri(i, k) + d_u_con(i, k)
691 v_seri(i, k) = v_seri(i, k) + d_v_con(i, k)
692 ENDDO
693 ENDDO
694
695 IF (.not. conv_emanuel) THEN
696 z_apres = sum((q_seri + ql_seri) * zmasse, dim=2)
697 z_factor = (z_avant - (rain_con + snow_con) * dtphys) / z_apres
698 DO k = 1, llm
699 DO i = 1, klon
700 IF (z_factor(i) > 1. + 1E-8 .OR. z_factor(i) < 1. - 1E-8) THEN
701 q_seri(i, k) = q_seri(i, k) * z_factor(i)
702 ENDIF
703 ENDDO
704 ENDDO
705 ENDIF
706
707 ! Convection s\`eche (thermiques ou ajustement)
708
709 d_t_ajs = 0.
710 d_u_ajs = 0.
711 d_v_ajs = 0.
712 d_q_ajs = 0.
713 fm_therm = 0.
714 entr_therm = 0.
715
716 if (iflag_thermals == 0) then
717 ! Ajustement sec
718 CALL ajsec(paprs, play, t_seri, q_seri, d_t_ajs, d_q_ajs)
719 t_seri = t_seri + d_t_ajs
720 q_seri = q_seri + d_q_ajs
721 else
722 call calltherm(dtphys, play, paprs, pphi, u_seri, v_seri, t_seri, &
723 q_seri, d_u_ajs, d_v_ajs, d_t_ajs, d_q_ajs, fm_therm, entr_therm)
724 endif
725
726 ! Caclul des ratqs
727
728 ! ratqs convectifs \`a l'ancienne en fonction de (q(z = 0) - q) / q
729 ! on \'ecrase le tableau ratqsc calcul\'e par clouds_gno
730 if (iflag_cldcon == 1) then
731 do k = 1, llm
732 do i = 1, klon
733 if(ptconv(i, k)) then
734 ratqsc(i, k) = ratqsbas + fact_cldcon &
735 * (q_seri(i, 1) - q_seri(i, k)) / q_seri(i, k)
736 else
737 ratqsc(i, k) = 0.
738 endif
739 enddo
740 enddo
741 endif
742
743 ! ratqs stables
744 do k = 1, llm
745 do i = 1, klon
746 ratqss(i, k) = ratqsbas + (ratqshaut - ratqsbas) &
747 * min((paprs(i, 1) - play(i, k)) / (paprs(i, 1) - 3e4), 1.)
748 enddo
749 enddo
750
751 ! ratqs final
752 if (iflag_cldcon == 1 .or. iflag_cldcon == 2) then
753 ! les ratqs sont une conbinaison de ratqss et ratqsc
754 ! ratqs final
755 ! 1e4 (en gros 3 heures), en dur pour le moment, est le temps de
756 ! relaxation des ratqs
757 ratqs = max(ratqs * exp(- dtphys * facttemps), ratqss)
758 ratqs = max(ratqs, ratqsc)
759 else
760 ! on ne prend que le ratqs stable pour fisrtilp
761 ratqs = ratqss
762 endif
763
764 CALL fisrtilp(dtphys, paprs, play, t_seri, q_seri, ptconv, ratqs, &
765 d_t_lsc, d_q_lsc, d_ql_lsc, rneb, cldliq, rain_lsc, snow_lsc, &
766 pfrac_impa, pfrac_nucl, pfrac_1nucl, frac_impa, frac_nucl, prfl, &
767 psfl, rhcl)
768
769 WHERE (rain_lsc < 0) rain_lsc = 0.
770 WHERE (snow_lsc < 0) snow_lsc = 0.
771 DO k = 1, llm
772 DO i = 1, klon
773 t_seri(i, k) = t_seri(i, k) + d_t_lsc(i, k)
774 q_seri(i, k) = q_seri(i, k) + d_q_lsc(i, k)
775 ql_seri(i, k) = ql_seri(i, k) + d_ql_lsc(i, k)
776 cldfra(i, k) = rneb(i, k)
777 IF (.NOT.new_oliq) cldliq(i, k) = ql_seri(i, k)
778 ENDDO
779 ENDDO
780
781 ! PRESCRIPTION DES NUAGES POUR LE RAYONNEMENT
782
783 ! 1. NUAGES CONVECTIFS
784
785 IF (iflag_cldcon <= - 1) THEN
786 ! seulement pour Tiedtke
787 snow_tiedtke = 0.
788 if (iflag_cldcon == - 1) then
789 rain_tiedtke = rain_con
790 else
791 rain_tiedtke = 0.
792 do k = 1, llm
793 do i = 1, klon
794 if (d_q_con(i, k) < 0.) then
795 rain_tiedtke(i) = rain_tiedtke(i) - d_q_con(i, k) / dtphys &
796 * zmasse(i, k)
797 endif
798 enddo
799 enddo
800 endif
801
802 ! Nuages diagnostiques pour Tiedtke
803 CALL diagcld1(paprs, play, rain_tiedtke, snow_tiedtke, ibas_con, &
804 itop_con, diafra, dialiq)
805 DO k = 1, llm
806 DO i = 1, klon
807 IF (diafra(i, k) > cldfra(i, k)) THEN
808 cldliq(i, k) = dialiq(i, k)
809 cldfra(i, k) = diafra(i, k)
810 ENDIF
811 ENDDO
812 ENDDO
813 ELSE IF (iflag_cldcon == 3) THEN
814 ! On prend pour les nuages convectifs le maximum du calcul de
815 ! la convection et du calcul du pas de temps pr\'ec\'edent diminu\'e
816 ! d'un facteur facttemps.
817 facteur = dtphys * facttemps
818 do k = 1, llm
819 do i = 1, klon
820 rnebcon(i, k) = rnebcon(i, k) * facteur
821 if (rnebcon0(i, k) * clwcon0(i, k) &
822 > rnebcon(i, k) * clwcon(i, k)) then
823 rnebcon(i, k) = rnebcon0(i, k)
824 clwcon(i, k) = clwcon0(i, k)
825 endif
826 enddo
827 enddo
828
829 ! On prend la somme des fractions nuageuses et des contenus en eau
830 cldfra = min(max(cldfra, rnebcon), 1.)
831 cldliq = cldliq + rnebcon * clwcon
832 ENDIF
833
834 ! 2. Nuages stratiformes
835
836 IF (ok_stratus) THEN
837 CALL diagcld2(paprs, play, t_seri, q_seri, diafra, dialiq)
838 DO k = 1, llm
839 DO i = 1, klon
840 IF (diafra(i, k) > cldfra(i, k)) THEN
841 cldliq(i, k) = dialiq(i, k)
842 cldfra(i, k) = diafra(i, k)
843 ENDIF
844 ENDDO
845 ENDDO
846 ENDIF
847
848 ! Precipitation totale
849 DO i = 1, klon
850 rain_fall(i) = rain_con(i) + rain_lsc(i)
851 snow_fall(i) = snow_con(i) + snow_lsc(i)
852 ENDDO
853
854 ! Humidit\'e relative pour diagnostic :
855 DO k = 1, llm
856 DO i = 1, klon
857 zx_t = t_seri(i, k)
858 zx_qs = r2es * FOEEW(zx_t, rtt >= zx_t) / play(i, k)
859 zx_qs = MIN(0.5, zx_qs)
860 zcor = 1. / (1. - retv * zx_qs)
861 zx_qs = zx_qs * zcor
862 zx_rh(i, k) = q_seri(i, k) / zx_qs
863 zqsat(i, k) = zx_qs
864 ENDDO
865 ENDDO
866
867 ! Param\`etres optiques des nuages et quelques param\`etres pour
868 ! diagnostics :
869 if (ok_newmicro) then
870 CALL newmicro(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, &
871 cldh, cldl, cldm, cldt, cldq, flwp, fiwp, flwc, fiwc)
872 else
873 CALL nuage(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, cldh, &
874 cldl, cldm, cldt, cldq)
875 endif
876
877 IF (MOD(itap - 1, radpas) == 0) THEN
878 wo = ozonecm(REAL(julien), paprs)
879 albsol = sum(falbe * pctsrf, dim = 2)
880 CALL radlwsw(dist, mu0, fract, paprs, play, tsol, albsol, t_seri, &
881 q_seri, wo, cldfra, cldemi, cldtau, heat, heat0, cool, cool0, &
882 radsol, albpla, topsw, toplw, solsw, sollw, sollwdown, topsw0, &
883 toplw0, solsw0, sollw0, lwdn0, lwdn, lwup0, lwup, swdn0, swdn, &
884 swup0, swup, ok_ade, topswad, solswad)
885 ENDIF
886
887 ! Ajouter la tendance des rayonnements (tous les pas)
888 DO k = 1, llm
889 DO i = 1, klon
890 t_seri(i, k) = t_seri(i, k) + (heat(i, k) - cool(i, k)) * dtphys &
891 / 86400.
892 ENDDO
893 ENDDO
894
895 ! Calculer le bilan du sol et la d\'erive de temp\'erature (couplage)
896 DO i = 1, klon
897 bils(i) = radsol(i) - sens(i) + zxfluxlat(i)
898 ENDDO
899
900 ! Param\'etrisation de l'orographie \`a l'\'echelle sous-maille :
901
902 IF (ok_orodr) THEN
903 ! S\'election des points pour lesquels le sch\'ema est actif :
904 igwd = 0
905 DO i = 1, klon
906 itest(i) = 0
907 IF (zpic(i) - zmea(i) > 100. .AND. zstd(i) > 10.) THEN
908 itest(i) = 1
909 igwd = igwd + 1
910 ENDIF
911 ENDDO
912
913 CALL drag_noro(klon, llm, dtphys, paprs, play, zmea, zstd, zsig, zgam, &
914 zthe, zpic, zval, itest, t_seri, u_seri, v_seri, zulow, zvlow, &
915 zustrdr, zvstrdr, d_t_oro, d_u_oro, d_v_oro)
916
917 ! ajout des tendances
918 DO k = 1, llm
919 DO i = 1, klon
920 t_seri(i, k) = t_seri(i, k) + d_t_oro(i, k)
921 u_seri(i, k) = u_seri(i, k) + d_u_oro(i, k)
922 v_seri(i, k) = v_seri(i, k) + d_v_oro(i, k)
923 ENDDO
924 ENDDO
925 ENDIF
926
927 IF (ok_orolf) THEN
928 ! S\'election des points pour lesquels le sch\'ema est actif :
929 igwd = 0
930 DO i = 1, klon
931 itest(i) = 0
932 IF (zpic(i) - zmea(i) > 100.) THEN
933 itest(i) = 1
934 igwd = igwd + 1
935 ENDIF
936 ENDDO
937
938 CALL lift_noro(klon, llm, dtphys, paprs, play, rlat, zmea, zstd, zpic, &
939 itest, t_seri, u_seri, v_seri, zulow, zvlow, zustrli, zvstrli, &
940 d_t_lif, d_u_lif, d_v_lif)
941
942 ! Ajout des tendances :
943 DO k = 1, llm
944 DO i = 1, klon
945 t_seri(i, k) = t_seri(i, k) + d_t_lif(i, k)
946 u_seri(i, k) = u_seri(i, k) + d_u_lif(i, k)
947 v_seri(i, k) = v_seri(i, k) + d_v_lif(i, k)
948 ENDDO
949 ENDDO
950 ENDIF
951
952 ! Stress n\'ecessaires : toute la physique
953
954 DO i = 1, klon
955 zustrph(i) = 0.
956 zvstrph(i) = 0.
957 ENDDO
958 DO k = 1, llm
959 DO i = 1, klon
960 zustrph(i) = zustrph(i) + (u_seri(i, k) - u(i, k)) / dtphys &
961 * zmasse(i, k)
962 zvstrph(i) = zvstrph(i) + (v_seri(i, k) - v(i, k)) / dtphys &
963 * zmasse(i, k)
964 ENDDO
965 ENDDO
966
967 CALL aaam_bud(rg, romega, rlat, rlon, pphis, zustrdr, zustrli, zustrph, &
968 zvstrdr, zvstrli, zvstrph, paprs, u, v, aam, torsfc)
969
970 ! Calcul des tendances traceurs
971 call phytrac(julien, time, firstcal, lafin, dtphys, t, paprs, play, mfu, &
972 mfd, pde_u, pen_d, ycoefh, fm_therm, entr_therm, yu1, yv1, ftsol, &
973 pctsrf, frac_impa, frac_nucl, da, phi, mp, upwd, dnwd, tr_seri, &
974 zmasse, ncid_startphy)
975
976 ! Calculer le transport de l'eau et de l'energie (diagnostique)
977 CALL transp(paprs, t_seri, q_seri, u_seri, v_seri, zphi, ve, vq, ue, uq)
978
979 ! diag. bilKP
980
981 CALL transp_lay(paprs, t_seri, q_seri, u_seri, v_seri, zphi, &
982 ve_lay, vq_lay, ue_lay, uq_lay)
983
984 ! Accumuler les variables a stocker dans les fichiers histoire:
985
986 ! conversion Ec en énergie thermique
987 DO k = 1, llm
988 DO i = 1, klon
989 d_t_ec(i, k) = 0.5 / (RCPD * (1. + RVTMP2 * q_seri(i, k))) &
990 * (u(i, k)**2 + v(i, k)**2 - u_seri(i, k)**2 - v_seri(i, k)**2)
991 t_seri(i, k) = t_seri(i, k) + d_t_ec(i, k)
992 d_t_ec(i, k) = d_t_ec(i, k) / dtphys
993 END DO
994 END DO
995
996 ! SORTIES
997
998 ! prw = eau precipitable
999 DO i = 1, klon
1000 prw(i) = 0.
1001 DO k = 1, llm
1002 prw(i) = prw(i) + q_seri(i, k) * zmasse(i, k)
1003 ENDDO
1004 ENDDO
1005
1006 ! Convertir les incrementations en tendances
1007
1008 DO k = 1, llm
1009 DO i = 1, klon
1010 d_u(i, k) = (u_seri(i, k) - u(i, k)) / dtphys
1011 d_v(i, k) = (v_seri(i, k) - v(i, k)) / dtphys
1012 d_t(i, k) = (t_seri(i, k) - t(i, k)) / dtphys
1013 d_qx(i, k, ivap) = (q_seri(i, k) - qx(i, k, ivap)) / dtphys
1014 d_qx(i, k, iliq) = (ql_seri(i, k) - qx(i, k, iliq)) / dtphys
1015 ENDDO
1016 ENDDO
1017
1018 DO iq = 3, nqmx
1019 DO k = 1, llm
1020 DO i = 1, klon
1021 d_qx(i, k, iq) = (tr_seri(i, k, iq - 2) - qx(i, k, iq)) / dtphys
1022 ENDDO
1023 ENDDO
1024 ENDDO
1025
1026 ! Sauvegarder les valeurs de t et q a la fin de la physique:
1027 DO k = 1, llm
1028 DO i = 1, klon
1029 t_ancien(i, k) = t_seri(i, k)
1030 q_ancien(i, k) = q_seri(i, k)
1031 ENDDO
1032 ENDDO
1033
1034 CALL histwrite_phy("phis", pphis)
1035 CALL histwrite_phy("aire", airephy)
1036 CALL histwrite_phy("psol", paprs(:, 1))
1037 CALL histwrite_phy("precip", rain_fall + snow_fall)
1038 CALL histwrite_phy("plul", rain_lsc + snow_lsc)
1039 CALL histwrite_phy("pluc", rain_con + snow_con)
1040 CALL histwrite_phy("tsol", tsol)
1041 CALL histwrite_phy("t2m", zt2m)
1042 CALL histwrite_phy("q2m", zq2m)
1043 CALL histwrite_phy("u10m", u10m)
1044 CALL histwrite_phy("v10m", v10m)
1045 CALL histwrite_phy("snow", snow_fall)
1046 CALL histwrite_phy("cdrm", cdragm)
1047 CALL histwrite_phy("cdrh", cdragh)
1048 CALL histwrite_phy("topl", toplw)
1049 CALL histwrite_phy("evap", evap)
1050 CALL histwrite_phy("sols", solsw)
1051 CALL histwrite_phy("soll", sollw)
1052 CALL histwrite_phy("solldown", sollwdown)
1053 CALL histwrite_phy("bils", bils)
1054 CALL histwrite_phy("sens", - sens)
1055 CALL histwrite_phy("fder", fder)
1056 CALL histwrite_phy("dtsvdfo", d_ts(:, is_oce))
1057 CALL histwrite_phy("dtsvdft", d_ts(:, is_ter))
1058 CALL histwrite_phy("dtsvdfg", d_ts(:, is_lic))
1059 CALL histwrite_phy("dtsvdfi", d_ts(:, is_sic))
1060
1061 DO nsrf = 1, nbsrf
1062 CALL histwrite_phy("pourc_"//clnsurf(nsrf), pctsrf(:, nsrf) * 100.)
1063 CALL histwrite_phy("fract_"//clnsurf(nsrf), pctsrf(:, nsrf))
1064 CALL histwrite_phy("sens_"//clnsurf(nsrf), flux_t(:, nsrf))
1065 CALL histwrite_phy("lat_"//clnsurf(nsrf), fluxlat(:, nsrf))
1066 CALL histwrite_phy("tsol_"//clnsurf(nsrf), ftsol(:, nsrf))
1067 CALL histwrite_phy("taux_"//clnsurf(nsrf), flux_u(:, nsrf))
1068 CALL histwrite_phy("tauy_"//clnsurf(nsrf), flux_v(:, nsrf))
1069 CALL histwrite_phy("rugs_"//clnsurf(nsrf), frugs(:, nsrf))
1070 CALL histwrite_phy("albe_"//clnsurf(nsrf), falbe(:, nsrf))
1071 CALL histwrite_phy("u10m_"//clnsurf(nsrf), u10m_srf(:, nsrf))
1072 CALL histwrite_phy("v10m_"//clnsurf(nsrf), v10m_srf(:, nsrf))
1073 END DO
1074
1075 CALL histwrite_phy("albs", albsol)
1076 CALL histwrite_phy("tro3", wo * dobson_u * 1e3 / zmasse / rmo3 * md)
1077 CALL histwrite_phy("rugs", zxrugs)
1078 CALL histwrite_phy("s_pblh", s_pblh)
1079 CALL histwrite_phy("s_pblt", s_pblt)
1080 CALL histwrite_phy("s_lcl", s_lcl)
1081 CALL histwrite_phy("s_capCL", s_capCL)
1082 CALL histwrite_phy("s_oliqCL", s_oliqCL)
1083 CALL histwrite_phy("s_cteiCL", s_cteiCL)
1084 CALL histwrite_phy("s_therm", s_therm)
1085 CALL histwrite_phy("s_trmb1", s_trmb1)
1086 CALL histwrite_phy("s_trmb2", s_trmb2)
1087 CALL histwrite_phy("s_trmb3", s_trmb3)
1088
1089 if (conv_emanuel) then
1090 CALL histwrite_phy("ptop", ema_pct)
1091 CALL histwrite_phy("dnwd0", - mp)
1092 end if
1093
1094 CALL histwrite_phy("temp", t_seri)
1095 CALL histwrite_phy("vitu", u_seri)
1096 CALL histwrite_phy("vitv", v_seri)
1097 CALL histwrite_phy("geop", zphi)
1098 CALL histwrite_phy("pres", play)
1099 CALL histwrite_phy("dtvdf", d_t_vdf)
1100 CALL histwrite_phy("dqvdf", d_q_vdf)
1101 CALL histwrite_phy("rhum", zx_rh)
1102 CALL histwrite_phy("d_t_ec", d_t_ec)
1103 CALL histwrite_phy("dtsw0", heat0 / 86400.)
1104 CALL histwrite_phy("dtlw0", - cool0 / 86400.)
1105 CALL histwrite_phy("msnow", sum(fsnow * pctsrf, dim = 2))
1106 call histwrite_phy("qsurf", sum(fqsurf * pctsrf, dim = 2))
1107
1108 if (ok_instan) call histsync(nid_ins)
1109
1110 IF (lafin) then
1111 call NF95_CLOSE(ncid_startphy)
1112 CALL phyredem(pctsrf, ftsol, ftsoil, fqsurf, qsol, &
1113 fsnow, falbe, fevap, rain_fall, snow_fall, solsw, sollw, dlw, &
1114 radsol, frugs, agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, &
1115 t_ancien, q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, &
1116 w01)
1117 end IF
1118
1119 firstcal = .FALSE.
1120
1121 END SUBROUTINE physiq
1122
1123 end module physiq_m

  ViewVC Help
Powered by ViewVC 1.1.17