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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 158 - (show annotations)
Tue Jul 21 14:44:45 2015 UTC (8 years, 9 months ago) by guez
File size: 57632 byte(s)
Subroutine sugwd sets variables of module yoegwd. Better to put it
into module yoegwd.

Variables of module yoegwd other than NKTOPG, NSTRA can be symbolic
constants. sugwd now only sets NKTOPG, NSTRA. Simplified the
computation of NKTOPG, NSTRA by making the local variable zpm1r an
array instead of a scalar and calling ifirstloc.

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

  ViewVC Help
Powered by ViewVC 1.1.21