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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 178 - (show annotations)
Fri Mar 11 18:47:26 2016 UTC (8 years, 2 months ago) by guez
File size: 57235 byte(s)
Moved variables date0, deltat, datasz_max, ncvar_ids, point, buff_pos,
buffer, regular from module histcom_var to modules where they are
defined.

Removed procedure ioipslmpp, useless for a sequential program.

Added argument datasz_max to histwrite_real (to avoid circular
dependency with histwrite).

Removed useless variables and computations everywhere.

Changed real litteral constants from default kind to double precision
in lwb, lwu, lwvn, sw1s, swtt, swtt1, swu.

Removed unused arguments: paer of sw, sw1s, sw2s, swclr; pcldsw of
sw1s, sw2s; pdsig, prayl of swr; co2_ppm of clmain, clqh; tsol of
transp_lay; nsrf of screenp; kcrit and kknu of gwstress; pstd of
orosetup.

Added output of relative humidity.

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