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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 175 - (show annotations)
Fri Feb 5 16:02:34 2016 UTC (8 years, 3 months ago) by guez
File size: 57159 byte(s)
Added argument itau_phy to ini_histins, phyetat0, phytrac and
phyredem0. Removed variable itau_phy of module temps. Avoiding side
effect in etat0 and phyetat0. The procedures ini_histins, phyetat0,
phytrac and phyredem0 are all called by physiq so there is no
cascading variable penalty.

In procedure inifilr, made the condition on colat0 weaker to allow for
rounding error.

Removed arguments flux_o, flux_g and t_slab of clmain, flux_o and
flux_g of clqh and interfsurf_hq, tslab and seaice of phyetat0 and
phyredem. NetCDF variables TSLAB and SEAICE no longer in
restartphy.nc. All these variables were related to the not-implemented
slab ocean. seaice and tslab were just set to 0 in phyetat0 and never
used nor changed. flux_o and flux_g were computed in clmain but never
used in physiq.

Removed argument swnet of clqh. Was used only to compute a local
variable, swdown, which was not used.

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

  ViewVC Help
Powered by ViewVC 1.1.21