/[lmdze]/trunk/libf/phylmd/physiq.f90
ViewVC logotype

Contents of /trunk/libf/phylmd/physiq.f90

Parent Directory Parent Directory | Revision Log Revision Log


Revision 61 - (show annotations)
Fri Apr 20 14:58:43 2012 UTC (12 years ago) by guez
File size: 71581 byte(s)
No more included file in LMDZE, not even "netcdf.inc".

Created a variable containing the list of common source files in
GNUmakefile. So we now also see clearly files that are specific to
each program.

Split module "histcom". Assembled resulting files in directory
"Histcom".

Removed aliasing in calls to "laplacien".

1 module physiq_m
2
3 IMPLICIT none
4
5 contains
6
7 SUBROUTINE physiq(lafin, rdayvrai, time, dtphys, paprs, play, pphi, pphis, &
8 u, v, t, qx, omega, d_u, d_v, d_t, d_qx, d_ps, dudyn, PVteta)
9
10 ! From phylmd/physiq.F, version 1.22 2006/02/20 09:38:28 (SVN revision 678)
11 ! Author: Z.X. Li (LMD/CNRS) 1993
12
13 ! This is the main procedure for the "physics" part of the program.
14
15 use aaam_bud_m, only: aaam_bud
16 USE abort_gcm_m, ONLY: abort_gcm
17 use ajsec_m, only: ajsec
18 USE calendar, ONLY: ymds2ju
19 use calltherm_m, only: calltherm
20 USE clesphys, ONLY: cdhmax, cdmmax, co2_ppm, ecrit_hf, ecrit_ins, &
21 ecrit_mth, ecrit_reg, ecrit_tra, ksta, ksta_ter, ok_kzmin
22 USE clesphys2, ONLY: cycle_diurne, iflag_con, nbapp_rad, new_oliq, &
23 ok_orodr, ok_orolf, soil_model
24 USE clmain_m, ONLY: clmain
25 USE comgeomphy, ONLY: airephy, cuphy, cvphy
26 USE concvl_m, ONLY: concvl
27 USE conf_gcm_m, ONLY: offline, raz_date
28 USE conf_phys_m, ONLY: conf_phys
29 USE ctherm, ONLY: iflag_thermals, nsplit_thermals
30 use diagcld2_m, only: diagcld2
31 use diagetpq_m, only: diagetpq
32 USE dimens_m, ONLY: iim, jjm, llm, nqmx
33 USE dimphy, ONLY: klon, nbtr
34 USE dimsoil, ONLY: nsoilmx
35 use drag_noro_m, only: drag_noro
36 USE fcttre, ONLY: foeew, qsatl, qsats, thermcep
37 USE hgardfou_m, ONLY: hgardfou
38 USE histsync_m, ONLY: histsync
39 USE histwrite_m, ONLY: histwrite
40 USE indicesol, ONLY: clnsurf, epsfra, is_lic, is_oce, is_sic, is_ter, &
41 nbsrf
42 USE ini_histhf_m, ONLY: ini_histhf
43 USE ini_histday_m, ONLY: ini_histday
44 USE ini_histins_m, ONLY: ini_histins
45 USE oasis_m, ONLY: ok_oasis
46 USE orbite_m, ONLY: orbite, zenang
47 USE ozonecm_m, ONLY: ozonecm
48 USE phyetat0_m, ONLY: phyetat0, rlat, rlon
49 USE phyredem_m, ONLY: phyredem
50 USE phystokenc_m, ONLY: phystokenc
51 USE phytrac_m, ONLY: phytrac
52 USE qcheck_m, ONLY: qcheck
53 use radlwsw_m, only: radlwsw
54 use sugwd_m, only: sugwd
55 USE suphec_m, ONLY: ra, rcpd, retv, rg, rlvtt, romega, rsigma, rtt
56 USE temps, ONLY: annee_ref, day_ref, itau_phy
57 USE yoethf_m, ONLY: r2es, rvtmp2
58
59 ! Arguments:
60
61 REAL, intent(in):: rdayvrai
62 ! (elapsed time since January 1st 0h of the starting year, in days)
63
64 REAL, intent(in):: time ! heure de la journée en fraction de jour
65 REAL, intent(in):: dtphys ! pas d'integration pour la physique (seconde)
66 logical, intent(in):: lafin ! dernier passage
67
68 REAL, intent(in):: paprs(klon, llm + 1)
69 ! (pression pour chaque inter-couche, en Pa)
70
71 REAL, intent(in):: play(klon, llm)
72 ! (input pression pour le mileu de chaque couche (en Pa))
73
74 REAL, intent(in):: pphi(klon, llm)
75 ! (input geopotentiel de chaque couche (g z) (reference sol))
76
77 REAL, intent(in):: pphis(klon) ! input geopotentiel du sol
78
79 REAL, intent(in):: u(klon, llm)
80 ! vitesse dans la direction X (de O a E) en m/s
81
82 REAL, intent(in):: v(klon, llm) ! vitesse Y (de S a N) en m/s
83 REAL, intent(in):: t(klon, llm) ! input temperature (K)
84
85 REAL, intent(in):: qx(klon, llm, nqmx)
86 ! (humidité spécifique et fractions massiques des autres traceurs)
87
88 REAL omega(klon, llm) ! input vitesse verticale en Pa/s
89 REAL, intent(out):: d_u(klon, llm) ! tendance physique de "u" (m/s/s)
90 REAL, intent(out):: d_v(klon, llm) ! tendance physique de "v" (m/s/s)
91 REAL, intent(out):: d_t(klon, llm) ! tendance physique de "t" (K/s)
92 REAL d_qx(klon, llm, nqmx) ! output tendance physique de "qx" (kg/kg/s)
93 REAL d_ps(klon) ! output tendance physique de la pression au sol
94
95 LOGICAL:: firstcal = .true.
96
97 INTEGER nbteta
98 PARAMETER(nbteta = 3)
99
100 REAL PVteta(klon, nbteta)
101 ! (output vorticite potentielle a des thetas constantes)
102
103 LOGICAL ok_cvl ! pour activer le nouveau driver pour convection KE
104 PARAMETER (ok_cvl = .TRUE.)
105 LOGICAL ok_gust ! pour activer l'effet des gust sur flux surface
106 PARAMETER (ok_gust = .FALSE.)
107
108 LOGICAL check ! Verifier la conservation du modele en eau
109 PARAMETER (check = .FALSE.)
110
111 LOGICAL, PARAMETER:: ok_stratus = .FALSE.
112 ! Ajouter artificiellement les stratus
113
114 ! Parametres lies au coupleur OASIS:
115 INTEGER, SAVE:: npas, nexca
116 logical rnpb
117 parameter(rnpb = .true.)
118
119 character(len = 6), save:: ocean
120 ! (type de modèle océan à utiliser: "force" ou "slab" mais pas "couple")
121
122 logical ok_ocean
123 SAVE ok_ocean
124
125 ! "slab" ocean
126 REAL, save:: tslab(klon) ! temperature of ocean slab
127 REAL, save:: seaice(klon) ! glace de mer (kg/m2)
128 REAL fluxo(klon) ! flux turbulents ocean-glace de mer
129 REAL fluxg(klon) ! flux turbulents ocean-atmosphere
130
131 ! Modele thermique du sol, a activer pour le cycle diurne:
132 logical, save:: ok_veget
133 LOGICAL, save:: ok_journe ! sortir le fichier journalier
134
135 LOGICAL ok_mensuel ! sortir le fichier mensuel
136
137 LOGICAL ok_instan ! sortir le fichier instantane
138 save ok_instan
139
140 LOGICAL ok_region ! sortir le fichier regional
141 PARAMETER (ok_region = .FALSE.)
142
143 ! pour phsystoke avec thermiques
144 REAL fm_therm(klon, llm + 1)
145 REAL entr_therm(klon, llm)
146 real, save:: q2(klon, llm + 1, nbsrf)
147
148 INTEGER ivap ! indice de traceurs pour vapeur d'eau
149 PARAMETER (ivap = 1)
150 INTEGER iliq ! indice de traceurs pour eau liquide
151 PARAMETER (iliq = 2)
152
153 REAL, save:: t_ancien(klon, llm), q_ancien(klon, llm)
154 LOGICAL, save:: ancien_ok
155
156 REAL d_t_dyn(klon, llm) ! tendance dynamique pour "t" (K/s)
157 REAL d_q_dyn(klon, llm) ! tendance dynamique pour "q" (kg/kg/s)
158
159 real da(klon, llm), phi(klon, llm, llm), mp(klon, llm)
160
161 !IM Amip2 PV a theta constante
162
163 CHARACTER(LEN = 3) ctetaSTD(nbteta)
164 DATA ctetaSTD/'350', '380', '405'/
165 REAL rtetaSTD(nbteta)
166 DATA rtetaSTD/350., 380., 405./
167
168 !MI Amip2 PV a theta constante
169
170 INTEGER klevp1
171 PARAMETER(klevp1 = llm + 1)
172
173 REAL swdn0(klon, klevp1), swdn(klon, klevp1)
174 REAL swup0(klon, klevp1), swup(klon, klevp1)
175 SAVE swdn0, swdn, swup0, swup
176
177 REAL lwdn0(klon, klevp1), lwdn(klon, klevp1)
178 REAL lwup0(klon, klevp1), lwup(klon, klevp1)
179 SAVE lwdn0, lwdn, lwup0, lwup
180
181 !IM Amip2
182 ! variables a une pression donnee
183
184 integer nlevSTD
185 PARAMETER(nlevSTD = 17)
186 real rlevSTD(nlevSTD)
187 DATA rlevSTD/100000., 92500., 85000., 70000., &
188 60000., 50000., 40000., 30000., 25000., 20000., &
189 15000., 10000., 7000., 5000., 3000., 2000., 1000./
190 CHARACTER(LEN = 4) clevSTD(nlevSTD)
191 DATA clevSTD/'1000', '925 ', '850 ', '700 ', '600 ', &
192 '500 ', '400 ', '300 ', '250 ', '200 ', '150 ', '100 ', &
193 '70 ', '50 ', '30 ', '20 ', '10 '/
194
195 ! prw: precipitable water
196 real prw(klon)
197
198 ! flwp, fiwp = Liquid Water Path & Ice Water Path (kg/m2)
199 ! flwc, fiwc = Liquid Water Content & Ice Water Content (kg/kg)
200 REAL flwp(klon), fiwp(klon)
201 REAL flwc(klon, llm), fiwc(klon, llm)
202
203 INTEGER kmax, lmax
204 PARAMETER(kmax = 8, lmax = 8)
205 INTEGER kmaxm1, lmaxm1
206 PARAMETER(kmaxm1 = kmax-1, lmaxm1 = lmax-1)
207
208 REAL zx_tau(kmaxm1), zx_pc(lmaxm1)
209 DATA zx_tau/0.0, 0.3, 1.3, 3.6, 9.4, 23., 60./
210 DATA zx_pc/50., 180., 310., 440., 560., 680., 800./
211
212 ! cldtopres pression au sommet des nuages
213 REAL cldtopres(lmaxm1)
214 DATA cldtopres/50., 180., 310., 440., 560., 680., 800./
215
216 ! taulev: numero du niveau de tau dans les sorties ISCCP
217 CHARACTER(LEN = 4) taulev(kmaxm1)
218
219 DATA taulev/'tau0', 'tau1', 'tau2', 'tau3', 'tau4', 'tau5', 'tau6'/
220 CHARACTER(LEN = 3) pclev(lmaxm1)
221 DATA pclev/'pc1', 'pc2', 'pc3', 'pc4', 'pc5', 'pc6', 'pc7'/
222
223 CHARACTER(LEN = 28) cnameisccp(lmaxm1, kmaxm1)
224 DATA cnameisccp/'pc< 50hPa, tau< 0.3', 'pc= 50-180hPa, tau< 0.3', &
225 'pc= 180-310hPa, tau< 0.3', 'pc= 310-440hPa, tau< 0.3', &
226 'pc= 440-560hPa, tau< 0.3', 'pc= 560-680hPa, tau< 0.3', &
227 'pc= 680-800hPa, tau< 0.3', 'pc< 50hPa, tau= 0.3-1.3', &
228 'pc= 50-180hPa, tau= 0.3-1.3', 'pc= 180-310hPa, tau= 0.3-1.3', &
229 'pc= 310-440hPa, tau= 0.3-1.3', 'pc= 440-560hPa, tau= 0.3-1.3', &
230 'pc= 560-680hPa, tau= 0.3-1.3', 'pc= 680-800hPa, tau= 0.3-1.3', &
231 'pc< 50hPa, tau= 1.3-3.6', 'pc= 50-180hPa, tau= 1.3-3.6', &
232 'pc= 180-310hPa, tau= 1.3-3.6', 'pc= 310-440hPa, tau= 1.3-3.6', &
233 'pc= 440-560hPa, tau= 1.3-3.6', 'pc= 560-680hPa, tau= 1.3-3.6', &
234 'pc= 680-800hPa, tau= 1.3-3.6', 'pc< 50hPa, tau= 3.6-9.4', &
235 'pc= 50-180hPa, tau= 3.6-9.4', 'pc= 180-310hPa, tau= 3.6-9.4', &
236 'pc= 310-440hPa, tau= 3.6-9.4', 'pc= 440-560hPa, tau= 3.6-9.4', &
237 'pc= 560-680hPa, tau= 3.6-9.4', 'pc= 680-800hPa, tau= 3.6-9.4', &
238 'pc< 50hPa, tau= 9.4-23', 'pc= 50-180hPa, tau= 9.4-23', &
239 'pc= 180-310hPa, tau= 9.4-23', 'pc= 310-440hPa, tau= 9.4-23', &
240 'pc= 440-560hPa, tau= 9.4-23', 'pc= 560-680hPa, tau= 9.4-23', &
241 'pc= 680-800hPa, tau= 9.4-23', 'pc< 50hPa, tau= 23-60', &
242 'pc= 50-180hPa, tau= 23-60', 'pc= 180-310hPa, tau= 23-60', &
243 'pc= 310-440hPa, tau= 23-60', 'pc= 440-560hPa, tau= 23-60', &
244 'pc= 560-680hPa, tau= 23-60', 'pc= 680-800hPa, tau= 23-60', &
245 'pc< 50hPa, tau> 60.', 'pc= 50-180hPa, tau> 60.', &
246 'pc= 180-310hPa, tau> 60.', 'pc= 310-440hPa, tau> 60.', &
247 'pc= 440-560hPa, tau> 60.', 'pc= 560-680hPa, tau> 60.', &
248 'pc= 680-800hPa, tau> 60.'/
249
250 !IM ISCCP simulator v3.4
251
252 integer nid_hf, nid_hf3d
253 save nid_hf, nid_hf3d
254
255 ! Variables propres a la physique
256
257 INTEGER, save:: radpas
258 ! (Radiative transfer computations are made every "radpas" call to
259 ! "physiq".)
260
261 REAL radsol(klon)
262 SAVE radsol ! bilan radiatif au sol calcule par code radiatif
263
264 INTEGER, SAVE:: itap ! number of calls to "physiq"
265
266 REAL, save:: ftsol(klon, nbsrf) ! skin temperature of surface fraction
267
268 REAL, save:: ftsoil(klon, nsoilmx, nbsrf)
269 ! soil temperature of surface fraction
270
271 REAL fevap(klon, nbsrf)
272 SAVE fevap ! evaporation
273 REAL fluxlat(klon, nbsrf)
274 SAVE fluxlat
275
276 REAL fqsurf(klon, nbsrf)
277 SAVE fqsurf ! humidite de l'air au contact de la surface
278
279 REAL, save:: qsol(klon) ! hauteur d'eau dans le sol
280
281 REAL fsnow(klon, nbsrf)
282 SAVE fsnow ! epaisseur neigeuse
283
284 REAL falbe(klon, nbsrf)
285 SAVE falbe ! albedo par type de surface
286 REAL falblw(klon, nbsrf)
287 SAVE falblw ! albedo par type de surface
288
289 ! Paramètres de l'orographie à l'échelle sous-maille (OESM) :
290 REAL, save:: zmea(klon) ! orographie moyenne
291 REAL, save:: zstd(klon) ! deviation standard de l'OESM
292 REAL, save:: zsig(klon) ! pente de l'OESM
293 REAL, save:: zgam(klon) ! anisotropie de l'OESM
294 REAL, save:: zthe(klon) ! orientation de l'OESM
295 REAL, save:: zpic(klon) ! Maximum de l'OESM
296 REAL, save:: zval(klon) ! Minimum de l'OESM
297 REAL, save:: rugoro(klon) ! longueur de rugosite de l'OESM
298
299 REAL zulow(klon), zvlow(klon)
300
301 INTEGER igwd, idx(klon), itest(klon)
302
303 REAL agesno(klon, nbsrf)
304 SAVE agesno ! age de la neige
305
306 REAL run_off_lic_0(klon)
307 SAVE run_off_lic_0
308 !KE43
309 ! Variables liees a la convection de K. Emanuel (sb):
310
311 REAL bas, top ! cloud base and top levels
312 SAVE bas
313 SAVE top
314
315 REAL Ma(klon, llm) ! undilute upward mass flux
316 SAVE Ma
317 REAL qcondc(klon, llm) ! in-cld water content from convect
318 SAVE qcondc
319 REAL ema_work1(klon, llm), ema_work2(klon, llm)
320 SAVE ema_work1, ema_work2
321
322 REAL wd(klon) ! sb
323 SAVE wd ! sb
324
325 ! Variables locales pour la couche limite (al1):
326
327 ! Variables locales:
328
329 REAL cdragh(klon) ! drag coefficient pour T and Q
330 REAL cdragm(klon) ! drag coefficient pour vent
331
332 !AA Pour phytrac
333 REAL ycoefh(klon, llm) ! coef d'echange pour phytrac
334 REAL yu1(klon) ! vents dans la premiere couche U
335 REAL yv1(klon) ! vents dans la premiere couche V
336 REAL ffonte(klon, nbsrf) !Flux thermique utilise pour fondre la neige
337 REAL fqcalving(klon, nbsrf) !Flux d'eau "perdue" par la surface
338 ! !et necessaire pour limiter la
339 ! !hauteur de neige, en kg/m2/s
340 REAL zxffonte(klon), zxfqcalving(klon)
341
342 REAL pfrac_impa(klon, llm)! Produits des coefs lessivage impaction
343 save pfrac_impa
344 REAL pfrac_nucl(klon, llm)! Produits des coefs lessivage nucleation
345 save pfrac_nucl
346 REAL pfrac_1nucl(klon, llm)! Produits des coefs lessi nucl (alpha = 1)
347 save pfrac_1nucl
348 REAL frac_impa(klon, llm) ! fractions d'aerosols lessivees (impaction)
349 REAL frac_nucl(klon, llm) ! idem (nucleation)
350
351 !AA
352 REAL rain_fall(klon) ! pluie
353 REAL snow_fall(klon) ! neige
354 save snow_fall, rain_fall
355 !IM cf FH pour Tiedtke 080604
356 REAL rain_tiedtke(klon), snow_tiedtke(klon)
357
358 REAL evap(klon), devap(klon) ! evaporation et sa derivee
359 REAL sens(klon), dsens(klon) ! chaleur sensible et sa derivee
360 REAL dlw(klon) ! derivee infra rouge
361 SAVE dlw
362 REAL bils(klon) ! bilan de chaleur au sol
363 REAL fder(klon) ! Derive de flux (sensible et latente)
364 save fder
365 REAL ve(klon) ! integr. verticale du transport meri. de l'energie
366 REAL vq(klon) ! integr. verticale du transport meri. de l'eau
367 REAL ue(klon) ! integr. verticale du transport zonal de l'energie
368 REAL uq(klon) ! integr. verticale du transport zonal de l'eau
369
370 REAL frugs(klon, nbsrf) ! longueur de rugosite
371 save frugs
372 REAL zxrugs(klon) ! longueur de rugosite
373
374 ! Conditions aux limites
375
376 INTEGER julien
377
378 INTEGER, SAVE:: lmt_pas ! number of time steps of "physics" per day
379 REAL pctsrf(klon, nbsrf)
380 !IM
381 REAL pctsrf_new(klon, nbsrf) !pourcentage surfaces issus d'ORCHIDEE
382
383 SAVE pctsrf ! sous-fraction du sol
384 REAL albsol(klon)
385 SAVE albsol ! albedo du sol total
386 REAL albsollw(klon)
387 SAVE albsollw ! albedo du sol total
388
389 REAL, SAVE:: wo(klon, llm) ! column density of ozone in a cell, in kDU
390
391 ! Declaration des procedures appelees
392
393 EXTERNAL alboc ! calculer l'albedo sur ocean
394 !KE43
395 EXTERNAL conema3 ! convect4.3
396 EXTERNAL fisrtilp ! schema de condensation a grande echelle (pluie)
397 EXTERNAL nuage ! calculer les proprietes radiatives
398 EXTERNAL transp ! transport total de l'eau et de l'energie
399
400 ! Variables locales
401
402 real clwcon(klon, llm), rnebcon(klon, llm)
403 real clwcon0(klon, llm), rnebcon0(klon, llm)
404
405 save rnebcon, clwcon
406
407 REAL rhcl(klon, llm) ! humiditi relative ciel clair
408 REAL dialiq(klon, llm) ! eau liquide nuageuse
409 REAL diafra(klon, llm) ! fraction nuageuse
410 REAL cldliq(klon, llm) ! eau liquide nuageuse
411 REAL cldfra(klon, llm) ! fraction nuageuse
412 REAL cldtau(klon, llm) ! epaisseur optique
413 REAL cldemi(klon, llm) ! emissivite infrarouge
414
415 REAL fluxq(klon, llm, nbsrf) ! flux turbulent d'humidite
416 REAL fluxt(klon, llm, nbsrf) ! flux turbulent de chaleur
417 REAL fluxu(klon, llm, nbsrf) ! flux turbulent de vitesse u
418 REAL fluxv(klon, llm, nbsrf) ! flux turbulent de vitesse v
419
420 REAL zxfluxt(klon, llm)
421 REAL zxfluxq(klon, llm)
422 REAL zxfluxu(klon, llm)
423 REAL zxfluxv(klon, llm)
424
425 ! Le rayonnement n'est pas calcule tous les pas, il faut donc
426 ! que les variables soient rémanentes
427 REAL, save:: heat(klon, llm) ! chauffage solaire
428 REAL heat0(klon, llm) ! chauffage solaire ciel clair
429 REAL cool(klon, llm) ! refroidissement infrarouge
430 REAL cool0(klon, llm) ! refroidissement infrarouge ciel clair
431 REAL topsw(klon), toplw(klon), solsw(klon), sollw(klon)
432 real sollwdown(klon) ! downward LW flux at surface
433 REAL topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)
434 REAL albpla(klon)
435 REAL fsollw(klon, nbsrf) ! bilan flux IR pour chaque sous surface
436 REAL fsolsw(klon, nbsrf) ! flux solaire absorb. pour chaque sous surface
437 SAVE cool, albpla, topsw, toplw, solsw, sollw, sollwdown
438 SAVE topsw0, toplw0, solsw0, sollw0, heat0, cool0
439
440 INTEGER itaprad
441 SAVE itaprad
442
443 REAL conv_q(klon, llm) ! convergence de l'humidite (kg/kg/s)
444 REAL conv_t(klon, llm) ! convergence of temperature (K/s)
445
446 REAL cldl(klon), cldm(klon), cldh(klon) !nuages bas, moyen et haut
447 REAL cldt(klon), cldq(klon) !nuage total, eau liquide integree
448
449 REAL zxtsol(klon), zxqsurf(klon), zxsnow(klon), zxfluxlat(klon)
450
451 REAL dist, rmu0(klon), fract(klon)
452 REAL zdtime ! pas de temps du rayonnement (s)
453 real zlongi
454
455 REAL z_avant(klon), z_apres(klon), z_factor(klon)
456 LOGICAL zx_ajustq
457
458 REAL za, zb
459 REAL zx_t, zx_qs, zdelta, zcor
460 real zqsat(klon, llm)
461 INTEGER i, k, iq, nsrf
462 REAL t_coup
463 PARAMETER (t_coup = 234.0)
464
465 REAL zphi(klon, llm)
466
467 !IM cf. AM Variables locales pour la CLA (hbtm2)
468
469 REAL, SAVE:: pblh(klon, nbsrf) ! Hauteur de couche limite
470 REAL, SAVE:: plcl(klon, nbsrf) ! Niveau de condensation de la CLA
471 REAL, SAVE:: capCL(klon, nbsrf) ! CAPE de couche limite
472 REAL, SAVE:: oliqCL(klon, nbsrf) ! eau_liqu integree de couche limite
473 REAL, SAVE:: cteiCL(klon, nbsrf) ! cloud top instab. crit. couche limite
474 REAL, SAVE:: pblt(klon, nbsrf) ! T a la Hauteur de couche limite
475 REAL, SAVE:: therm(klon, nbsrf)
476 REAL, SAVE:: trmb1(klon, nbsrf) ! deep_cape
477 REAL, SAVE:: trmb2(klon, nbsrf) ! inhibition
478 REAL, SAVE:: trmb3(klon, nbsrf) ! Point Omega
479 ! Grdeurs de sorties
480 REAL s_pblh(klon), s_lcl(klon), s_capCL(klon)
481 REAL s_oliqCL(klon), s_cteiCL(klon), s_pblt(klon)
482 REAL s_therm(klon), s_trmb1(klon), s_trmb2(klon)
483 REAL s_trmb3(klon)
484
485 ! Variables locales pour la convection de K. Emanuel (sb):
486
487 REAL upwd(klon, llm) ! saturated updraft mass flux
488 REAL dnwd(klon, llm) ! saturated downdraft mass flux
489 REAL dnwd0(klon, llm) ! unsaturated downdraft mass flux
490 REAL tvp(klon, llm) ! virtual temp of lifted parcel
491 REAL cape(klon) ! CAPE
492 SAVE cape
493
494 REAL pbase(klon) ! cloud base pressure
495 SAVE pbase
496 REAL bbase(klon) ! cloud base buoyancy
497 SAVE bbase
498 REAL rflag(klon) ! flag fonctionnement de convect
499 INTEGER iflagctrl(klon) ! flag fonctionnement de convect
500 ! -- convect43:
501 INTEGER ntra ! nb traceurs pour convect4.3
502 REAL dtvpdt1(klon, llm), dtvpdq1(klon, llm)
503 REAL dplcldt(klon), dplcldr(klon)
504
505 ! Variables du changement
506
507 ! con: convection
508 ! lsc: large scale condensation
509 ! ajs: ajustement sec
510 ! eva: évaporation de l'eau liquide nuageuse
511 ! vdf: vertical diffusion in boundary layer
512 REAL d_t_con(klon, llm), d_q_con(klon, llm)
513 REAL d_u_con(klon, llm), d_v_con(klon, llm)
514 REAL d_t_lsc(klon, llm), d_q_lsc(klon, llm), d_ql_lsc(klon, llm)
515 REAL d_t_ajs(klon, llm), d_q_ajs(klon, llm)
516 REAL d_u_ajs(klon, llm), d_v_ajs(klon, llm)
517 REAL rneb(klon, llm)
518
519 REAL pmfu(klon, llm), pmfd(klon, llm)
520 REAL pen_u(klon, llm), pen_d(klon, llm)
521 REAL pde_u(klon, llm), pde_d(klon, llm)
522 INTEGER kcbot(klon), kctop(klon), kdtop(klon)
523 REAL pmflxr(klon, llm + 1), pmflxs(klon, llm + 1)
524 REAL prfl(klon, llm + 1), psfl(klon, llm + 1)
525
526 INTEGER,save:: ibas_con(klon), itop_con(klon)
527
528 REAL rain_con(klon), rain_lsc(klon)
529 REAL snow_con(klon), snow_lsc(klon)
530 REAL d_ts(klon, nbsrf)
531
532 REAL d_u_vdf(klon, llm), d_v_vdf(klon, llm)
533 REAL d_t_vdf(klon, llm), d_q_vdf(klon, llm)
534
535 REAL d_u_oro(klon, llm), d_v_oro(klon, llm)
536 REAL d_t_oro(klon, llm)
537 REAL d_u_lif(klon, llm), d_v_lif(klon, llm)
538 REAL d_t_lif(klon, llm)
539
540 REAL ratqs(klon, llm), ratqss(klon, llm), ratqsc(klon, llm)
541 real ratqsbas, ratqshaut
542 save ratqsbas, ratqshaut, ratqs
543
544 ! Parametres lies au nouveau schema de nuages (SB, PDF)
545 real, save:: fact_cldcon
546 real, save:: facttemps
547 logical ok_newmicro
548 save ok_newmicro
549 real facteur
550
551 integer iflag_cldcon
552 save iflag_cldcon
553
554 logical ptconv(klon, llm)
555
556 ! Variables locales pour effectuer les appels en série :
557
558 REAL t_seri(klon, llm), q_seri(klon, llm)
559 REAL ql_seri(klon, llm), qs_seri(klon, llm)
560 REAL u_seri(klon, llm), v_seri(klon, llm)
561
562 REAL tr_seri(klon, llm, nbtr)
563 REAL d_tr(klon, llm, nbtr)
564
565 REAL zx_rh(klon, llm)
566
567 REAL zustrdr(klon), zvstrdr(klon)
568 REAL zustrli(klon), zvstrli(klon)
569 REAL zustrph(klon), zvstrph(klon)
570 REAL aam, torsfc
571
572 REAL dudyn(iim + 1, jjm + 1, llm)
573
574 REAL zx_tmp_fi2d(klon) ! variable temporaire grille physique
575 REAL zx_tmp_2d(iim, jjm + 1), zx_tmp_3d(iim, jjm + 1, llm)
576
577 INTEGER, SAVE:: nid_day, nid_ins
578
579 REAL ve_lay(klon, llm) ! transport meri. de l'energie a chaque niveau vert.
580 REAL vq_lay(klon, llm) ! transport meri. de l'eau a chaque niveau vert.
581 REAL ue_lay(klon, llm) ! transport zonal de l'energie a chaque niveau vert.
582 REAL uq_lay(klon, llm) ! transport zonal de l'eau a chaque niveau vert.
583
584 REAL zsto
585
586 character(len = 20) modname
587 character(len = 80) abort_message
588 logical ok_sync
589 real date0
590
591 ! Variables liées au bilan d'énergie et d'enthalpie :
592 REAL ztsol(klon)
593 REAL d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec
594 REAL, SAVE:: d_h_vcol_phy
595 REAL fs_bound, fq_bound
596 REAL zero_v(klon)
597 CHARACTER(LEN = 15) ztit
598 INTEGER:: ip_ebil = 0 ! print level for energy conservation diagnostics
599 INTEGER, SAVE:: if_ebil ! level for energy conservation diagnostics
600
601 REAL d_t_ec(klon, llm) ! tendance due à la conversion Ec -> E thermique
602 REAL ZRCPD
603
604 REAL t2m(klon, nbsrf), q2m(klon, nbsrf) ! temperature and humidity at 2 m
605 REAL u10m(klon, nbsrf), v10m(klon, nbsrf) !vents a 10m
606 REAL zt2m(klon), zq2m(klon) !temp., hum. 2m moyenne s/ 1 maille
607 REAL zu10m(klon), zv10m(klon) !vents a 10m moyennes s/1 maille
608 !jq Aerosol effects (Johannes Quaas, 27/11/2003)
609 REAL sulfate(klon, llm) ! SO4 aerosol concentration [ug/m3]
610
611 REAL, save:: sulfate_pi(klon, llm)
612 ! (SO4 aerosol concentration, in ug/m3, pre-industrial value)
613
614 REAL cldtaupi(klon, llm)
615 ! (Cloud optical thickness for pre-industrial (pi) aerosols)
616
617 REAL re(klon, llm) ! Cloud droplet effective radius
618 REAL fl(klon, llm) ! denominator of re
619
620 ! Aerosol optical properties
621 REAL tau_ae(klon, llm, 2), piz_ae(klon, llm, 2)
622 REAL cg_ae(klon, llm, 2)
623
624 REAL topswad(klon), solswad(klon) ! Aerosol direct effect.
625 ! ok_ade = True -ADE = topswad-topsw
626
627 REAL topswai(klon), solswai(klon) ! Aerosol indirect effect.
628 ! ok_aie = True ->
629 ! ok_ade = True -AIE = topswai-topswad
630 ! ok_ade = F -AIE = topswai-topsw
631
632 REAL aerindex(klon) ! POLDER aerosol index
633
634 ! Parameters
635 LOGICAL ok_ade, ok_aie ! Apply aerosol (in)direct effects or not
636 REAL bl95_b0, bl95_b1 ! Parameter in Boucher and Lohmann (1995)
637
638 SAVE ok_ade, ok_aie, bl95_b0, bl95_b1
639 SAVE u10m
640 SAVE v10m
641 SAVE t2m
642 SAVE q2m
643 SAVE ffonte
644 SAVE fqcalving
645 SAVE piz_ae
646 SAVE tau_ae
647 SAVE cg_ae
648 SAVE rain_con
649 SAVE snow_con
650 SAVE topswai
651 SAVE topswad
652 SAVE solswai
653 SAVE solswad
654 SAVE d_u_con
655 SAVE d_v_con
656 SAVE rnebcon0
657 SAVE clwcon0
658
659 real zmasse(klon, llm)
660 ! (column-density of mass of air in a cell, in kg m-2)
661
662 real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2
663
664 !----------------------------------------------------------------
665
666 modname = 'physiq'
667 IF (if_ebil >= 1) THEN
668 DO i = 1, klon
669 zero_v(i) = 0.
670 END DO
671 END IF
672 ok_sync = .TRUE.
673 IF (nqmx < 2) THEN
674 abort_message = 'eaux vapeur et liquide sont indispensables'
675 CALL abort_gcm(modname, abort_message, 1)
676 ENDIF
677
678 test_firstcal: IF (firstcal) THEN
679 ! initialiser
680 u10m = 0.
681 v10m = 0.
682 t2m = 0.
683 q2m = 0.
684 ffonte = 0.
685 fqcalving = 0.
686 piz_ae = 0.
687 tau_ae = 0.
688 cg_ae = 0.
689 rain_con(:) = 0.
690 snow_con(:) = 0.
691 bl95_b0 = 0.
692 bl95_b1 = 0.
693 topswai(:) = 0.
694 topswad(:) = 0.
695 solswai(:) = 0.
696 solswad(:) = 0.
697
698 d_u_con = 0.0
699 d_v_con = 0.0
700 rnebcon0 = 0.0
701 clwcon0 = 0.0
702 rnebcon = 0.0
703 clwcon = 0.0
704
705 pblh =0. ! Hauteur de couche limite
706 plcl =0. ! Niveau de condensation de la CLA
707 capCL =0. ! CAPE de couche limite
708 oliqCL =0. ! eau_liqu integree de couche limite
709 cteiCL =0. ! cloud top instab. crit. couche limite
710 pblt =0. ! T a la Hauteur de couche limite
711 therm =0.
712 trmb1 =0. ! deep_cape
713 trmb2 =0. ! inhibition
714 trmb3 =0. ! Point Omega
715
716 IF (if_ebil >= 1) d_h_vcol_phy = 0.
717
718 ! appel a la lecture du run.def physique
719
720 call conf_phys(ocean, ok_veget, ok_journe, ok_mensuel, &
721 ok_instan, fact_cldcon, facttemps, ok_newmicro, &
722 iflag_cldcon, ratqsbas, ratqshaut, if_ebil, &
723 ok_ade, ok_aie, &
724 bl95_b0, bl95_b1, &
725 iflag_thermals, nsplit_thermals)
726
727 ! Initialiser les compteurs:
728
729 frugs = 0.
730 itap = 0
731 itaprad = 0
732 CALL phyetat0("startphy.nc", pctsrf, ftsol, ftsoil, ocean, tslab, &
733 seaice, fqsurf, qsol, fsnow, falbe, falblw, fevap, rain_fall, &
734 snow_fall, solsw, sollwdown, dlw, radsol, frugs, agesno, zmea, &
735 zstd, zsig, zgam, zthe, zpic, zval, t_ancien, q_ancien, &
736 ancien_ok, rnebcon, ratqs, clwcon, run_off_lic_0)
737
738 ! ATTENTION : il faudra a terme relire q2 dans l'etat initial
739 q2 = 1.e-8
740
741 radpas = NINT(86400. / dtphys / nbapp_rad)
742
743 ! on remet le calendrier a zero
744 IF (raz_date) itau_phy = 0
745
746 PRINT *, 'cycle_diurne = ', cycle_diurne
747
748 IF(ocean.NE.'force ') THEN
749 ok_ocean = .TRUE.
750 ENDIF
751
752 CALL printflag(radpas, ok_ocean, ok_oasis, ok_journe, ok_instan, &
753 ok_region)
754
755 IF (dtphys*REAL(radpas) > 21600..AND.cycle_diurne) THEN
756 print *,'Nbre d appels au rayonnement insuffisant'
757 print *,"Au minimum 4 appels par jour si cycle diurne"
758 abort_message = 'Nbre d appels au rayonnement insuffisant'
759 call abort_gcm(modname, abort_message, 1)
760 ENDIF
761 print *,"Clef pour la convection, iflag_con = ", iflag_con
762 print *,"Clef pour le driver de la convection, ok_cvl = ", &
763 ok_cvl
764
765 ! Initialisation pour la convection de K.E. (sb):
766 IF (iflag_con >= 3) THEN
767
768 print *,"*** Convection de Kerry Emanuel 4.3 "
769
770 !IM15/11/02 rajout initialisation ibas_con, itop_con cf. SB =>BEG
771 DO i = 1, klon
772 ibas_con(i) = 1
773 itop_con(i) = 1
774 ENDDO
775 !IM15/11/02 rajout initialisation ibas_con, itop_con cf. SB =>END
776
777 ENDIF
778
779 IF (ok_orodr) THEN
780 rugoro = MAX(1e-5, zstd * zsig / 2)
781 CALL SUGWD(paprs, play)
782 else
783 rugoro = 0.
784 ENDIF
785
786 lmt_pas = NINT(86400. / dtphys) ! tous les jours
787 print *, 'Number of time steps of "physics" per day: ', lmt_pas
788
789 ecrit_ins = NINT(ecrit_ins/dtphys)
790 ecrit_hf = NINT(ecrit_hf/dtphys)
791 ecrit_mth = NINT(ecrit_mth/dtphys)
792 ecrit_tra = NINT(86400.*ecrit_tra/dtphys)
793 ecrit_reg = NINT(ecrit_reg/dtphys)
794
795 ! Initialiser le couplage si necessaire
796
797 npas = 0
798 nexca = 0
799
800 print *,'AVANT HIST IFLAG_CON = ', iflag_con
801
802 ! Initialisation des sorties
803
804 call ini_histhf(dtphys, nid_hf, nid_hf3d)
805 call ini_histday(dtphys, ok_journe, nid_day, nqmx)
806 call ini_histins(dtphys, ok_instan, nid_ins)
807 CALL ymds2ju(annee_ref, 1, int(day_ref), 0., date0)
808 !XXXPB Positionner date0 pour initialisation de ORCHIDEE
809 WRITE(*, *) 'physiq date0: ', date0
810 ENDIF test_firstcal
811
812 ! Mettre a zero des variables de sortie (pour securite)
813
814 DO i = 1, klon
815 d_ps(i) = 0.0
816 ENDDO
817 DO iq = 1, nqmx
818 DO k = 1, llm
819 DO i = 1, klon
820 d_qx(i, k, iq) = 0.0
821 ENDDO
822 ENDDO
823 ENDDO
824 da = 0.
825 mp = 0.
826 phi = 0.
827
828 ! Ne pas affecter les valeurs entrées de u, v, h, et q :
829
830 DO k = 1, llm
831 DO i = 1, klon
832 t_seri(i, k) = t(i, k)
833 u_seri(i, k) = u(i, k)
834 v_seri(i, k) = v(i, k)
835 q_seri(i, k) = qx(i, k, ivap)
836 ql_seri(i, k) = qx(i, k, iliq)
837 qs_seri(i, k) = 0.
838 ENDDO
839 ENDDO
840 IF (nqmx >= 3) THEN
841 tr_seri(:, :, :nqmx-2) = qx(:, :, 3:nqmx)
842 ELSE
843 tr_seri(:, :, 1) = 0.
844 ENDIF
845
846 DO i = 1, klon
847 ztsol(i) = 0.
848 ENDDO
849 DO nsrf = 1, nbsrf
850 DO i = 1, klon
851 ztsol(i) = ztsol(i) + ftsol(i, nsrf)*pctsrf(i, nsrf)
852 ENDDO
853 ENDDO
854
855 IF (if_ebil >= 1) THEN
856 ztit = 'after dynamics'
857 CALL diagetpq(airephy, ztit, ip_ebil, 1, 1, dtphys, t_seri, q_seri, &
858 ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &
859 d_ql, d_qs, d_ec)
860 ! Comme les tendances de la physique sont ajoutés dans la
861 ! dynamique, la variation d'enthalpie par la dynamique devrait
862 ! être égale à la variation de la physique au pas de temps
863 ! précédent. Donc la somme de ces 2 variations devrait être
864 ! nulle.
865 call diagphy(airephy, ztit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
866 zero_v, zero_v, zero_v, zero_v, ztsol, d_h_vcol + d_h_vcol_phy, &
867 d_qt, 0., fs_bound, fq_bound)
868 END IF
869
870 ! Diagnostic de la tendance dynamique :
871 IF (ancien_ok) THEN
872 DO k = 1, llm
873 DO i = 1, klon
874 d_t_dyn(i, k) = (t_seri(i, k) - t_ancien(i, k)) / dtphys
875 d_q_dyn(i, k) = (q_seri(i, k) - q_ancien(i, k)) / dtphys
876 ENDDO
877 ENDDO
878 ELSE
879 DO k = 1, llm
880 DO i = 1, klon
881 d_t_dyn(i, k) = 0.0
882 d_q_dyn(i, k) = 0.0
883 ENDDO
884 ENDDO
885 ancien_ok = .TRUE.
886 ENDIF
887
888 ! Ajouter le geopotentiel du sol:
889 DO k = 1, llm
890 DO i = 1, klon
891 zphi(i, k) = pphi(i, k) + pphis(i)
892 ENDDO
893 ENDDO
894
895 ! Check temperatures:
896 CALL hgardfou(t_seri, ftsol)
897
898 ! Incrementer le compteur de la physique
899 itap = itap + 1
900 julien = MOD(NINT(rdayvrai), 360)
901 if (julien == 0) julien = 360
902
903 forall (k = 1: llm) zmasse(:, k) = (paprs(:, k)-paprs(:, k + 1)) / rg
904
905 ! Mettre en action les conditions aux limites (albedo, sst, etc.).
906
907 ! Prescrire l'ozone et calculer l'albedo sur l'ocean.
908 wo = ozonecm(REAL(julien), paprs)
909
910 ! Évaporation de l'eau liquide nuageuse :
911 DO k = 1, llm
912 DO i = 1, klon
913 zb = MAX(0., ql_seri(i, k))
914 t_seri(i, k) = t_seri(i, k) &
915 - zb * RLVTT / RCPD / (1. + RVTMP2 * q_seri(i, k))
916 q_seri(i, k) = q_seri(i, k) + zb
917 ENDDO
918 ENDDO
919 ql_seri = 0.
920
921 IF (if_ebil >= 2) THEN
922 ztit = 'after reevap'
923 CALL diagetpq(airephy, ztit, ip_ebil, 2, 1, dtphys, t_seri, q_seri, &
924 ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &
925 d_ql, d_qs, d_ec)
926 call diagphy(airephy, ztit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
927 zero_v, zero_v, zero_v, zero_v, ztsol, d_h_vcol, d_qt, d_ec, &
928 fs_bound, fq_bound)
929
930 END IF
931
932 ! Appeler la diffusion verticale (programme de couche limite)
933
934 DO i = 1, klon
935 zxrugs(i) = 0.0
936 ENDDO
937 DO nsrf = 1, nbsrf
938 DO i = 1, klon
939 frugs(i, nsrf) = MAX(frugs(i, nsrf), 0.000015)
940 ENDDO
941 ENDDO
942 DO nsrf = 1, nbsrf
943 DO i = 1, klon
944 zxrugs(i) = zxrugs(i) + frugs(i, nsrf)*pctsrf(i, nsrf)
945 ENDDO
946 ENDDO
947
948 ! calculs necessaires au calcul de l'albedo dans l'interface
949
950 CALL orbite(REAL(julien), zlongi, dist)
951 IF (cycle_diurne) THEN
952 zdtime = dtphys * REAL(radpas)
953 CALL zenang(zlongi, time, zdtime, rmu0, fract)
954 ELSE
955 rmu0 = -999.999
956 ENDIF
957
958 ! Calcul de l'abedo moyen par maille
959 albsol(:) = 0.
960 albsollw(:) = 0.
961 DO nsrf = 1, nbsrf
962 DO i = 1, klon
963 albsol(i) = albsol(i) + falbe(i, nsrf) * pctsrf(i, nsrf)
964 albsollw(i) = albsollw(i) + falblw(i, nsrf) * pctsrf(i, nsrf)
965 ENDDO
966 ENDDO
967
968 ! Repartition sous maille des flux LW et SW
969 ! Repartition du longwave par sous-surface linearisee
970
971 DO nsrf = 1, nbsrf
972 DO i = 1, klon
973 fsollw(i, nsrf) = sollw(i) &
974 + 4.0*RSIGMA*ztsol(i)**3 * (ztsol(i)-ftsol(i, nsrf))
975 fsolsw(i, nsrf) = solsw(i)*(1.-falbe(i, nsrf))/(1.-albsol(i))
976 ENDDO
977 ENDDO
978
979 fder = dlw
980
981 ! Couche limite:
982
983 CALL clmain(dtphys, itap, date0, pctsrf, pctsrf_new, t_seri, q_seri, &
984 u_seri, v_seri, julien, rmu0, co2_ppm, ok_veget, ocean, npas, nexca, &
985 ftsol, soil_model, cdmmax, cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, &
986 qsol, paprs, play, fsnow, fqsurf, fevap, falbe, falblw, fluxlat, &
987 rain_fall, snow_fall, fsolsw, fsollw, sollwdown, fder, rlon, rlat, &
988 cuphy, cvphy, frugs, firstcal, lafin, agesno, rugoro, d_t_vdf, &
989 d_q_vdf, d_u_vdf, d_v_vdf, d_ts, fluxt, fluxq, fluxu, fluxv, cdragh, &
990 cdragm, q2, dsens, devap, ycoefh, yu1, yv1, t2m, q2m, u10m, v10m, &
991 pblh, capCL, oliqCL, cteiCL, pblT, therm, trmb1, trmb2, trmb3, plcl, &
992 fqcalving, ffonte, run_off_lic_0, fluxo, fluxg, tslab, seaice)
993
994 ! Incrémentation des flux
995
996 zxfluxt = 0.
997 zxfluxq = 0.
998 zxfluxu = 0.
999 zxfluxv = 0.
1000 DO nsrf = 1, nbsrf
1001 DO k = 1, llm
1002 DO i = 1, klon
1003 zxfluxt(i, k) = zxfluxt(i, k) + &
1004 fluxt(i, k, nsrf) * pctsrf(i, nsrf)
1005 zxfluxq(i, k) = zxfluxq(i, k) + &
1006 fluxq(i, k, nsrf) * pctsrf(i, nsrf)
1007 zxfluxu(i, k) = zxfluxu(i, k) + &
1008 fluxu(i, k, nsrf) * pctsrf(i, nsrf)
1009 zxfluxv(i, k) = zxfluxv(i, k) + &
1010 fluxv(i, k, nsrf) * pctsrf(i, nsrf)
1011 END DO
1012 END DO
1013 END DO
1014 DO i = 1, klon
1015 sens(i) = - zxfluxt(i, 1) ! flux de chaleur sensible au sol
1016 evap(i) = - zxfluxq(i, 1) ! flux d'evaporation au sol
1017 fder(i) = dlw(i) + dsens(i) + devap(i)
1018 ENDDO
1019
1020 DO k = 1, llm
1021 DO i = 1, klon
1022 t_seri(i, k) = t_seri(i, k) + d_t_vdf(i, k)
1023 q_seri(i, k) = q_seri(i, k) + d_q_vdf(i, k)
1024 u_seri(i, k) = u_seri(i, k) + d_u_vdf(i, k)
1025 v_seri(i, k) = v_seri(i, k) + d_v_vdf(i, k)
1026 ENDDO
1027 ENDDO
1028
1029 IF (if_ebil >= 2) THEN
1030 ztit = 'after clmain'
1031 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1032 ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &
1033 d_ql, d_qs, d_ec)
1034 call diagphy(airephy, ztit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
1035 sens, evap, zero_v, zero_v, ztsol, d_h_vcol, d_qt, d_ec, &
1036 fs_bound, fq_bound)
1037 END IF
1038
1039 ! Update surface temperature:
1040
1041 DO i = 1, klon
1042 zxtsol(i) = 0.0
1043 zxfluxlat(i) = 0.0
1044
1045 zt2m(i) = 0.0
1046 zq2m(i) = 0.0
1047 zu10m(i) = 0.0
1048 zv10m(i) = 0.0
1049 zxffonte(i) = 0.0
1050 zxfqcalving(i) = 0.0
1051
1052 s_pblh(i) = 0.0
1053 s_lcl(i) = 0.0
1054 s_capCL(i) = 0.0
1055 s_oliqCL(i) = 0.0
1056 s_cteiCL(i) = 0.0
1057 s_pblT(i) = 0.0
1058 s_therm(i) = 0.0
1059 s_trmb1(i) = 0.0
1060 s_trmb2(i) = 0.0
1061 s_trmb3(i) = 0.0
1062
1063 IF (abs(pctsrf(i, is_ter) + pctsrf(i, is_lic) + &
1064 pctsrf(i, is_oce) + pctsrf(i, is_sic) - 1.) > EPSFRA) &
1065 THEN
1066 WRITE(*, *) 'physiq : pb sous surface au point ', i, &
1067 pctsrf(i, 1 : nbsrf)
1068 ENDIF
1069 ENDDO
1070 DO nsrf = 1, nbsrf
1071 DO i = 1, klon
1072 ftsol(i, nsrf) = ftsol(i, nsrf) + d_ts(i, nsrf)
1073 zxtsol(i) = zxtsol(i) + ftsol(i, nsrf)*pctsrf(i, nsrf)
1074 zxfluxlat(i) = zxfluxlat(i) + fluxlat(i, nsrf)*pctsrf(i, nsrf)
1075
1076 zt2m(i) = zt2m(i) + t2m(i, nsrf)*pctsrf(i, nsrf)
1077 zq2m(i) = zq2m(i) + q2m(i, nsrf)*pctsrf(i, nsrf)
1078 zu10m(i) = zu10m(i) + u10m(i, nsrf)*pctsrf(i, nsrf)
1079 zv10m(i) = zv10m(i) + v10m(i, nsrf)*pctsrf(i, nsrf)
1080 zxffonte(i) = zxffonte(i) + ffonte(i, nsrf)*pctsrf(i, nsrf)
1081 zxfqcalving(i) = zxfqcalving(i) + &
1082 fqcalving(i, nsrf)*pctsrf(i, nsrf)
1083 s_pblh(i) = s_pblh(i) + pblh(i, nsrf)*pctsrf(i, nsrf)
1084 s_lcl(i) = s_lcl(i) + plcl(i, nsrf)*pctsrf(i, nsrf)
1085 s_capCL(i) = s_capCL(i) + capCL(i, nsrf) *pctsrf(i, nsrf)
1086 s_oliqCL(i) = s_oliqCL(i) + oliqCL(i, nsrf) *pctsrf(i, nsrf)
1087 s_cteiCL(i) = s_cteiCL(i) + cteiCL(i, nsrf) *pctsrf(i, nsrf)
1088 s_pblT(i) = s_pblT(i) + pblT(i, nsrf) *pctsrf(i, nsrf)
1089 s_therm(i) = s_therm(i) + therm(i, nsrf) *pctsrf(i, nsrf)
1090 s_trmb1(i) = s_trmb1(i) + trmb1(i, nsrf) *pctsrf(i, nsrf)
1091 s_trmb2(i) = s_trmb2(i) + trmb2(i, nsrf) *pctsrf(i, nsrf)
1092 s_trmb3(i) = s_trmb3(i) + trmb3(i, nsrf) *pctsrf(i, nsrf)
1093 ENDDO
1094 ENDDO
1095
1096 ! Si une sous-fraction n'existe pas, elle prend la temp. moyenne
1097
1098 DO nsrf = 1, nbsrf
1099 DO i = 1, klon
1100 IF (pctsrf(i, nsrf) < epsfra) ftsol(i, nsrf) = zxtsol(i)
1101
1102 IF (pctsrf(i, nsrf) < epsfra) t2m(i, nsrf) = zt2m(i)
1103 IF (pctsrf(i, nsrf) < epsfra) q2m(i, nsrf) = zq2m(i)
1104 IF (pctsrf(i, nsrf) < epsfra) u10m(i, nsrf) = zu10m(i)
1105 IF (pctsrf(i, nsrf) < epsfra) v10m(i, nsrf) = zv10m(i)
1106 IF (pctsrf(i, nsrf) < epsfra) ffonte(i, nsrf) = zxffonte(i)
1107 IF (pctsrf(i, nsrf) < epsfra) &
1108 fqcalving(i, nsrf) = zxfqcalving(i)
1109 IF (pctsrf(i, nsrf) < epsfra) pblh(i, nsrf) = s_pblh(i)
1110 IF (pctsrf(i, nsrf) < epsfra) plcl(i, nsrf) = s_lcl(i)
1111 IF (pctsrf(i, nsrf) < epsfra) capCL(i, nsrf) = s_capCL(i)
1112 IF (pctsrf(i, nsrf) < epsfra) oliqCL(i, nsrf) = s_oliqCL(i)
1113 IF (pctsrf(i, nsrf) < epsfra) cteiCL(i, nsrf) = s_cteiCL(i)
1114 IF (pctsrf(i, nsrf) < epsfra) pblT(i, nsrf) = s_pblT(i)
1115 IF (pctsrf(i, nsrf) < epsfra) therm(i, nsrf) = s_therm(i)
1116 IF (pctsrf(i, nsrf) < epsfra) trmb1(i, nsrf) = s_trmb1(i)
1117 IF (pctsrf(i, nsrf) < epsfra) trmb2(i, nsrf) = s_trmb2(i)
1118 IF (pctsrf(i, nsrf) < epsfra) trmb3(i, nsrf) = s_trmb3(i)
1119 ENDDO
1120 ENDDO
1121
1122 ! Calculer la derive du flux infrarouge
1123
1124 DO i = 1, klon
1125 dlw(i) = - 4.0*RSIGMA*zxtsol(i)**3
1126 ENDDO
1127
1128 ! Appeler la convection (au choix)
1129
1130 DO k = 1, llm
1131 DO i = 1, klon
1132 conv_q(i, k) = d_q_dyn(i, k) &
1133 + d_q_vdf(i, k)/dtphys
1134 conv_t(i, k) = d_t_dyn(i, k) &
1135 + d_t_vdf(i, k)/dtphys
1136 ENDDO
1137 ENDDO
1138 IF (check) THEN
1139 za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)
1140 print *, "avantcon = ", za
1141 ENDIF
1142 zx_ajustq = .FALSE.
1143 IF (iflag_con == 2) zx_ajustq = .TRUE.
1144 IF (zx_ajustq) THEN
1145 DO i = 1, klon
1146 z_avant(i) = 0.0
1147 ENDDO
1148 DO k = 1, llm
1149 DO i = 1, klon
1150 z_avant(i) = z_avant(i) + (q_seri(i, k) + ql_seri(i, k)) &
1151 *zmasse(i, k)
1152 ENDDO
1153 ENDDO
1154 ENDIF
1155
1156 select case (iflag_con)
1157 case (1)
1158 print *, 'Réactiver l''appel à "conlmd" dans "physiq.F".'
1159 stop 1
1160 case (2)
1161 CALL conflx(dtphys, paprs, play, t_seri, q_seri, conv_t, conv_q, &
1162 zxfluxq(1, 1), omega, d_t_con, d_q_con, rain_con, snow_con, pmfu, &
1163 pmfd, pen_u, pde_u, pen_d, pde_d, kcbot, kctop, kdtop, pmflxr, &
1164 pmflxs)
1165 WHERE (rain_con < 0.) rain_con = 0.
1166 WHERE (snow_con < 0.) snow_con = 0.
1167 DO i = 1, klon
1168 ibas_con(i) = llm + 1 - kcbot(i)
1169 itop_con(i) = llm + 1 - kctop(i)
1170 ENDDO
1171 case (3:)
1172 ! number of tracers for the convection scheme of Kerry Emanuel:
1173 ! la partie traceurs est faite dans phytrac
1174 ! on met ntra = 1 pour limiter les appels mais on peut
1175 ! supprimer les calculs / ftra.
1176 ntra = 1
1177 ! Schéma de convection modularisé et vectorisé :
1178 ! (driver commun aux versions 3 et 4)
1179
1180 IF (ok_cvl) THEN
1181 ! new driver for convectL
1182 CALL concvl(iflag_con, dtphys, paprs, play, t_seri, q_seri, &
1183 u_seri, v_seri, tr_seri, ntra, ema_work1, ema_work2, d_t_con, &
1184 d_q_con, d_u_con, d_v_con, d_tr, rain_con, snow_con, ibas_con, &
1185 itop_con, upwd, dnwd, dnwd0, Ma, cape, tvp, iflagctrl, pbase, &
1186 bbase, dtvpdt1, dtvpdq1, dplcldt, dplcldr, qcondc, wd, pmflxr, &
1187 pmflxs, da, phi, mp)
1188 clwcon0 = qcondc
1189 pmfu = upwd + dnwd
1190 ELSE
1191 ! conema3 ne contient pas les traceurs
1192 CALL conema3(dtphys, paprs, play, t_seri, q_seri, u_seri, v_seri, &
1193 tr_seri, ntra, ema_work1, ema_work2, d_t_con, d_q_con, &
1194 d_u_con, d_v_con, d_tr, rain_con, snow_con, ibas_con, &
1195 itop_con, upwd, dnwd, dnwd0, bas, top, Ma, cape, tvp, rflag, &
1196 pbase, bbase, dtvpdt1, dtvpdq1, dplcldt, dplcldr, clwcon0)
1197 ENDIF
1198
1199 IF (.NOT. ok_gust) THEN
1200 do i = 1, klon
1201 wd(i) = 0.0
1202 enddo
1203 ENDIF
1204
1205 ! Calcul des propriétés des nuages convectifs
1206
1207 DO k = 1, llm
1208 DO i = 1, klon
1209 zx_t = t_seri(i, k)
1210 IF (thermcep) THEN
1211 zdelta = MAX(0., SIGN(1., rtt-zx_t))
1212 zx_qs = r2es * FOEEW(zx_t, zdelta)/play(i, k)
1213 zx_qs = MIN(0.5, zx_qs)
1214 zcor = 1./(1.-retv*zx_qs)
1215 zx_qs = zx_qs*zcor
1216 ELSE
1217 IF (zx_t < t_coup) THEN
1218 zx_qs = qsats(zx_t)/play(i, k)
1219 ELSE
1220 zx_qs = qsatl(zx_t)/play(i, k)
1221 ENDIF
1222 ENDIF
1223 zqsat(i, k) = zx_qs
1224 ENDDO
1225 ENDDO
1226
1227 ! calcul des proprietes des nuages convectifs
1228 clwcon0 = fact_cldcon*clwcon0
1229 call clouds_gno &
1230 (klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, rnebcon0)
1231 case default
1232 print *, "iflag_con non-prevu", iflag_con
1233 stop 1
1234 END select
1235
1236 DO k = 1, llm
1237 DO i = 1, klon
1238 t_seri(i, k) = t_seri(i, k) + d_t_con(i, k)
1239 q_seri(i, k) = q_seri(i, k) + d_q_con(i, k)
1240 u_seri(i, k) = u_seri(i, k) + d_u_con(i, k)
1241 v_seri(i, k) = v_seri(i, k) + d_v_con(i, k)
1242 ENDDO
1243 ENDDO
1244
1245 IF (if_ebil >= 2) THEN
1246 ztit = 'after convect'
1247 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1248 ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &
1249 d_ql, d_qs, d_ec)
1250 call diagphy(airephy, ztit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
1251 zero_v, zero_v, rain_con, snow_con, ztsol, d_h_vcol, d_qt, d_ec, &
1252 fs_bound, fq_bound)
1253 END IF
1254
1255 IF (check) THEN
1256 za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)
1257 print *,"aprescon = ", za
1258 zx_t = 0.0
1259 za = 0.0
1260 DO i = 1, klon
1261 za = za + airephy(i)/REAL(klon)
1262 zx_t = zx_t + (rain_con(i)+ &
1263 snow_con(i))*airephy(i)/REAL(klon)
1264 ENDDO
1265 zx_t = zx_t/za*dtphys
1266 print *,"Precip = ", zx_t
1267 ENDIF
1268 IF (zx_ajustq) THEN
1269 DO i = 1, klon
1270 z_apres(i) = 0.0
1271 ENDDO
1272 DO k = 1, llm
1273 DO i = 1, klon
1274 z_apres(i) = z_apres(i) + (q_seri(i, k) + ql_seri(i, k)) &
1275 *zmasse(i, k)
1276 ENDDO
1277 ENDDO
1278 DO i = 1, klon
1279 z_factor(i) = (z_avant(i)-(rain_con(i) + snow_con(i))*dtphys) &
1280 /z_apres(i)
1281 ENDDO
1282 DO k = 1, llm
1283 DO i = 1, klon
1284 IF (z_factor(i) > 1. + 1E-8 .OR. z_factor(i) < 1. - 1E-8) THEN
1285 q_seri(i, k) = q_seri(i, k) * z_factor(i)
1286 ENDIF
1287 ENDDO
1288 ENDDO
1289 ENDIF
1290 zx_ajustq = .FALSE.
1291
1292 ! Convection sèche (thermiques ou ajustement)
1293
1294 d_t_ajs = 0.
1295 d_u_ajs = 0.
1296 d_v_ajs = 0.
1297 d_q_ajs = 0.
1298 fm_therm = 0.
1299 entr_therm = 0.
1300
1301 if (iflag_thermals == 0) then
1302 ! Ajustement sec
1303 CALL ajsec(paprs, play, t_seri, q_seri, d_t_ajs, d_q_ajs)
1304 t_seri = t_seri + d_t_ajs
1305 q_seri = q_seri + d_q_ajs
1306 else
1307 ! Thermiques
1308 call calltherm(dtphys, play, paprs, pphi, u_seri, v_seri, t_seri, &
1309 q_seri, d_u_ajs, d_v_ajs, d_t_ajs, d_q_ajs, fm_therm, entr_therm)
1310 endif
1311
1312 IF (if_ebil >= 2) THEN
1313 ztit = 'after dry_adjust'
1314 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1315 ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &
1316 d_ql, d_qs, d_ec)
1317 END IF
1318
1319 ! Caclul des ratqs
1320
1321 ! ratqs convectifs a l'ancienne en fonction de q(z = 0)-q / q
1322 ! on ecrase le tableau ratqsc calcule par clouds_gno
1323 if (iflag_cldcon == 1) then
1324 do k = 1, llm
1325 do i = 1, klon
1326 if(ptconv(i, k)) then
1327 ratqsc(i, k) = ratqsbas &
1328 +fact_cldcon*(q_seri(i, 1)-q_seri(i, k))/q_seri(i, k)
1329 else
1330 ratqsc(i, k) = 0.
1331 endif
1332 enddo
1333 enddo
1334 endif
1335
1336 ! ratqs stables
1337 do k = 1, llm
1338 do i = 1, klon
1339 ratqss(i, k) = ratqsbas + (ratqshaut-ratqsbas)* &
1340 min((paprs(i, 1)-play(i, k))/(paprs(i, 1)-30000.), 1.)
1341 enddo
1342 enddo
1343
1344 ! ratqs final
1345 if (iflag_cldcon == 1 .or.iflag_cldcon == 2) then
1346 ! les ratqs sont une conbinaison de ratqss et ratqsc
1347 ! ratqs final
1348 ! 1e4 (en gros 3 heures), en dur pour le moment, est le temps de
1349 ! relaxation des ratqs
1350 facteur = exp(-dtphys*facttemps)
1351 ratqs = max(ratqs*facteur, ratqss)
1352 ratqs = max(ratqs, ratqsc)
1353 else
1354 ! on ne prend que le ratqs stable pour fisrtilp
1355 ratqs = ratqss
1356 endif
1357
1358 ! Processus de condensation à grande echelle et processus de
1359 ! précipitation :
1360 CALL fisrtilp(dtphys, paprs, play, t_seri, q_seri, ptconv, ratqs, &
1361 d_t_lsc, d_q_lsc, d_ql_lsc, rneb, cldliq, rain_lsc, snow_lsc, &
1362 pfrac_impa, pfrac_nucl, pfrac_1nucl, frac_impa, frac_nucl, prfl, &
1363 psfl, rhcl)
1364
1365 WHERE (rain_lsc < 0) rain_lsc = 0.
1366 WHERE (snow_lsc < 0) snow_lsc = 0.
1367 DO k = 1, llm
1368 DO i = 1, klon
1369 t_seri(i, k) = t_seri(i, k) + d_t_lsc(i, k)
1370 q_seri(i, k) = q_seri(i, k) + d_q_lsc(i, k)
1371 ql_seri(i, k) = ql_seri(i, k) + d_ql_lsc(i, k)
1372 cldfra(i, k) = rneb(i, k)
1373 IF (.NOT.new_oliq) cldliq(i, k) = ql_seri(i, k)
1374 ENDDO
1375 ENDDO
1376 IF (check) THEN
1377 za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)
1378 print *,"apresilp = ", za
1379 zx_t = 0.0
1380 za = 0.0
1381 DO i = 1, klon
1382 za = za + airephy(i)/REAL(klon)
1383 zx_t = zx_t + (rain_lsc(i) &
1384 + snow_lsc(i))*airephy(i)/REAL(klon)
1385 ENDDO
1386 zx_t = zx_t/za*dtphys
1387 print *,"Precip = ", zx_t
1388 ENDIF
1389
1390 IF (if_ebil >= 2) THEN
1391 ztit = 'after fisrt'
1392 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1393 ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &
1394 d_ql, d_qs, d_ec)
1395 call diagphy(airephy, ztit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
1396 zero_v, zero_v, rain_lsc, snow_lsc, ztsol, d_h_vcol, d_qt, d_ec, &
1397 fs_bound, fq_bound)
1398 END IF
1399
1400 ! PRESCRIPTION DES NUAGES POUR LE RAYONNEMENT
1401
1402 ! 1. NUAGES CONVECTIFS
1403
1404 IF (iflag_cldcon.le.-1) THEN ! seulement pour Tiedtke
1405 snow_tiedtke = 0.
1406 if (iflag_cldcon == -1) then
1407 rain_tiedtke = rain_con
1408 else
1409 rain_tiedtke = 0.
1410 do k = 1, llm
1411 do i = 1, klon
1412 if (d_q_con(i, k) < 0.) then
1413 rain_tiedtke(i) = rain_tiedtke(i)-d_q_con(i, k)/dtphys &
1414 *zmasse(i, k)
1415 endif
1416 enddo
1417 enddo
1418 endif
1419
1420 ! Nuages diagnostiques pour Tiedtke
1421 CALL diagcld1(paprs, play, &
1422 rain_tiedtke, snow_tiedtke, ibas_con, itop_con, &
1423 diafra, dialiq)
1424 DO k = 1, llm
1425 DO i = 1, klon
1426 IF (diafra(i, k) > cldfra(i, k)) THEN
1427 cldliq(i, k) = dialiq(i, k)
1428 cldfra(i, k) = diafra(i, k)
1429 ENDIF
1430 ENDDO
1431 ENDDO
1432 ELSE IF (iflag_cldcon == 3) THEN
1433 ! On prend pour les nuages convectifs le max du calcul de la
1434 ! convection et du calcul du pas de temps précédent diminué d'un facteur
1435 ! facttemps
1436 facteur = dtphys *facttemps
1437 do k = 1, llm
1438 do i = 1, klon
1439 rnebcon(i, k) = rnebcon(i, k)*facteur
1440 if (rnebcon0(i, k)*clwcon0(i, k) > rnebcon(i, k)*clwcon(i, k)) &
1441 then
1442 rnebcon(i, k) = rnebcon0(i, k)
1443 clwcon(i, k) = clwcon0(i, k)
1444 endif
1445 enddo
1446 enddo
1447
1448 ! On prend la somme des fractions nuageuses et des contenus en eau
1449 cldfra = min(max(cldfra, rnebcon), 1.)
1450 cldliq = cldliq + rnebcon*clwcon
1451 ENDIF
1452
1453 ! 2. Nuages stratiformes
1454
1455 IF (ok_stratus) THEN
1456 CALL diagcld2(paprs, play, t_seri, q_seri, diafra, dialiq)
1457 DO k = 1, llm
1458 DO i = 1, klon
1459 IF (diafra(i, k) > cldfra(i, k)) THEN
1460 cldliq(i, k) = dialiq(i, k)
1461 cldfra(i, k) = diafra(i, k)
1462 ENDIF
1463 ENDDO
1464 ENDDO
1465 ENDIF
1466
1467 ! Precipitation totale
1468
1469 DO i = 1, klon
1470 rain_fall(i) = rain_con(i) + rain_lsc(i)
1471 snow_fall(i) = snow_con(i) + snow_lsc(i)
1472 ENDDO
1473
1474 IF (if_ebil >= 2) THEN
1475 ztit = "after diagcld"
1476 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1477 ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &
1478 d_ql, d_qs, d_ec)
1479 END IF
1480
1481 ! Humidité relative pour diagnostic:
1482 DO k = 1, llm
1483 DO i = 1, klon
1484 zx_t = t_seri(i, k)
1485 IF (thermcep) THEN
1486 zdelta = MAX(0., SIGN(1., rtt-zx_t))
1487 zx_qs = r2es * FOEEW(zx_t, zdelta)/play(i, k)
1488 zx_qs = MIN(0.5, zx_qs)
1489 zcor = 1./(1.-retv*zx_qs)
1490 zx_qs = zx_qs*zcor
1491 ELSE
1492 IF (zx_t < t_coup) THEN
1493 zx_qs = qsats(zx_t)/play(i, k)
1494 ELSE
1495 zx_qs = qsatl(zx_t)/play(i, k)
1496 ENDIF
1497 ENDIF
1498 zx_rh(i, k) = q_seri(i, k)/zx_qs
1499 zqsat(i, k) = zx_qs
1500 ENDDO
1501 ENDDO
1502
1503 ! Introduce the aerosol direct and first indirect radiative forcings:
1504 ! Johannes Quaas, 27/11/2003 (quaas@lmd.jussieu.fr)
1505 IF (ok_ade .OR. ok_aie) THEN
1506 ! Get sulfate aerosol distribution
1507 CALL readsulfate(rdayvrai, firstcal, sulfate)
1508 CALL readsulfate_preind(rdayvrai, firstcal, sulfate_pi)
1509
1510 ! Calculate aerosol optical properties (Olivier Boucher)
1511 CALL aeropt(play, paprs, t_seri, sulfate, rhcl, tau_ae, piz_ae, cg_ae, &
1512 aerindex)
1513 ELSE
1514 tau_ae = 0.
1515 piz_ae = 0.
1516 cg_ae = 0.
1517 ENDIF
1518
1519 ! Paramètres optiques des nuages et quelques paramètres pour
1520 ! diagnostics :
1521 if (ok_newmicro) then
1522 CALL newmicro(paprs, play, ok_newmicro, t_seri, cldliq, cldfra, &
1523 cldtau, cldemi, cldh, cldl, cldm, cldt, cldq, flwp, fiwp, flwc, &
1524 fiwc, ok_aie, sulfate, sulfate_pi, bl95_b0, bl95_b1, cldtaupi, &
1525 re, fl)
1526 else
1527 CALL nuage(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, cldh, &
1528 cldl, cldm, cldt, cldq, ok_aie, sulfate, sulfate_pi, bl95_b0, &
1529 bl95_b1, cldtaupi, re, fl)
1530 endif
1531
1532 ! Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.
1533 IF (MOD(itaprad, radpas) == 0) THEN
1534 DO i = 1, klon
1535 albsol(i) = falbe(i, is_oce) * pctsrf(i, is_oce) &
1536 + falbe(i, is_lic) * pctsrf(i, is_lic) &
1537 + falbe(i, is_ter) * pctsrf(i, is_ter) &
1538 + falbe(i, is_sic) * pctsrf(i, is_sic)
1539 albsollw(i) = falblw(i, is_oce) * pctsrf(i, is_oce) &
1540 + falblw(i, is_lic) * pctsrf(i, is_lic) &
1541 + falblw(i, is_ter) * pctsrf(i, is_ter) &
1542 + falblw(i, is_sic) * pctsrf(i, is_sic)
1543 ENDDO
1544 ! nouveau rayonnement (compatible Arpege-IFS):
1545 CALL radlwsw(dist, rmu0, fract, paprs, play, zxtsol, albsol, &
1546 albsollw, t_seri, q_seri, wo, cldfra, cldemi, cldtau, heat, &
1547 heat0, cool, cool0, radsol, albpla, topsw, toplw, solsw, sollw, &
1548 sollwdown, topsw0, toplw0, solsw0, sollw0, lwdn0, lwdn, lwup0, &
1549 lwup, swdn0, swdn, swup0, swup, ok_ade, ok_aie, tau_ae, piz_ae, &
1550 cg_ae, topswad, solswad, cldtaupi, topswai, solswai)
1551 itaprad = 0
1552 ENDIF
1553 itaprad = itaprad + 1
1554
1555 ! Ajouter la tendance des rayonnements (tous les pas)
1556
1557 DO k = 1, llm
1558 DO i = 1, klon
1559 t_seri(i, k) = t_seri(i, k) + (heat(i, k)-cool(i, k)) * dtphys/86400.
1560 ENDDO
1561 ENDDO
1562
1563 IF (if_ebil >= 2) THEN
1564 ztit = 'after rad'
1565 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1566 ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &
1567 d_ql, d_qs, d_ec)
1568 call diagphy(airephy, ztit, ip_ebil, topsw, toplw, solsw, sollw, &
1569 zero_v, zero_v, zero_v, zero_v, ztsol, d_h_vcol, d_qt, d_ec, &
1570 fs_bound, fq_bound)
1571 END IF
1572
1573 ! Calculer l'hydrologie de la surface
1574 DO i = 1, klon
1575 zxqsurf(i) = 0.0
1576 zxsnow(i) = 0.0
1577 ENDDO
1578 DO nsrf = 1, nbsrf
1579 DO i = 1, klon
1580 zxqsurf(i) = zxqsurf(i) + fqsurf(i, nsrf)*pctsrf(i, nsrf)
1581 zxsnow(i) = zxsnow(i) + fsnow(i, nsrf)*pctsrf(i, nsrf)
1582 ENDDO
1583 ENDDO
1584
1585 ! Calculer le bilan du sol et la dérive de température (couplage)
1586
1587 DO i = 1, klon
1588 bils(i) = radsol(i) - sens(i) + zxfluxlat(i)
1589 ENDDO
1590
1591 ! Paramétrisation de l'orographie à l'échelle sous-maille :
1592
1593 IF (ok_orodr) THEN
1594 ! selection des points pour lesquels le shema est actif:
1595 igwd = 0
1596 DO i = 1, klon
1597 itest(i) = 0
1598 IF (((zpic(i)-zmea(i)) > 100.).AND.(zstd(i) > 10.0)) THEN
1599 itest(i) = 1
1600 igwd = igwd + 1
1601 idx(igwd) = i
1602 ENDIF
1603 ENDDO
1604
1605 CALL drag_noro(klon, llm, dtphys, paprs, play, zmea, zstd, zsig, zgam, &
1606 zthe, zpic, zval, igwd, idx, itest, t_seri, u_seri, v_seri, &
1607 zulow, zvlow, zustrdr, zvstrdr, d_t_oro, d_u_oro, d_v_oro)
1608
1609 ! ajout des tendances
1610 DO k = 1, llm
1611 DO i = 1, klon
1612 t_seri(i, k) = t_seri(i, k) + d_t_oro(i, k)
1613 u_seri(i, k) = u_seri(i, k) + d_u_oro(i, k)
1614 v_seri(i, k) = v_seri(i, k) + d_v_oro(i, k)
1615 ENDDO
1616 ENDDO
1617 ENDIF
1618
1619 IF (ok_orolf) THEN
1620 ! Sélection des points pour lesquels le schéma est actif :
1621 igwd = 0
1622 DO i = 1, klon
1623 itest(i) = 0
1624 IF ((zpic(i) - zmea(i)) > 100.) THEN
1625 itest(i) = 1
1626 igwd = igwd + 1
1627 idx(igwd) = i
1628 ENDIF
1629 ENDDO
1630
1631 CALL lift_noro(klon, llm, dtphys, paprs, play, rlat, zmea, zstd, zpic, &
1632 itest, t_seri, u_seri, v_seri, zulow, zvlow, zustrli, zvstrli, &
1633 d_t_lif, d_u_lif, d_v_lif)
1634
1635 ! Ajout des tendances :
1636 DO k = 1, llm
1637 DO i = 1, klon
1638 t_seri(i, k) = t_seri(i, k) + d_t_lif(i, k)
1639 u_seri(i, k) = u_seri(i, k) + d_u_lif(i, k)
1640 v_seri(i, k) = v_seri(i, k) + d_v_lif(i, k)
1641 ENDDO
1642 ENDDO
1643 ENDIF
1644
1645 ! STRESS NECESSAIRES: TOUTE LA PHYSIQUE
1646
1647 DO i = 1, klon
1648 zustrph(i) = 0.
1649 zvstrph(i) = 0.
1650 ENDDO
1651 DO k = 1, llm
1652 DO i = 1, klon
1653 zustrph(i) = zustrph(i) + (u_seri(i, k)-u(i, k))/dtphys* zmasse(i, k)
1654 zvstrph(i) = zvstrph(i) + (v_seri(i, k)-v(i, k))/dtphys* zmasse(i, k)
1655 ENDDO
1656 ENDDO
1657
1658 CALL aaam_bud(ra, rg, romega, rlat, rlon, pphis, zustrdr, zustrli, &
1659 zustrph, zvstrdr, zvstrli, zvstrph, paprs, u, v, aam, torsfc)
1660
1661 IF (if_ebil >= 2) THEN
1662 ztit = 'after orography'
1663 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1664 ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &
1665 d_ql, d_qs, d_ec)
1666 END IF
1667
1668 ! Calcul des tendances traceurs
1669 call phytrac(rnpb, itap, lmt_pas, julien, time, firstcal, lafin, &
1670 nqmx-2, dtphys, u, t, paprs, play, pmfu, pmfd, pen_u, pde_u, &
1671 pen_d, pde_d, ycoefh, fm_therm, entr_therm, yu1, yv1, ftsol, pctsrf, &
1672 frac_impa, frac_nucl, pphis, albsol, rhcl, cldfra, rneb, &
1673 diafra, cldliq, pmflxr, pmflxs, prfl, psfl, da, phi, mp, upwd, dnwd, &
1674 tr_seri, zmasse)
1675
1676 IF (offline) THEN
1677 call phystokenc(dtphys, rlon, rlat, t, pmfu, pmfd, pen_u, pde_u, &
1678 pen_d, pde_d, fm_therm, entr_therm, ycoefh, yu1, yv1, ftsol, &
1679 pctsrf, frac_impa, frac_nucl, pphis, airephy, dtphys, itap)
1680 ENDIF
1681
1682 ! Calculer le transport de l'eau et de l'energie (diagnostique)
1683 CALL transp(paprs, zxtsol, t_seri, q_seri, u_seri, v_seri, zphi, ve, vq, &
1684 ue, uq)
1685
1686 ! diag. bilKP
1687
1688 CALL transp_lay(paprs, zxtsol, t_seri, q_seri, u_seri, v_seri, zphi, &
1689 ve_lay, vq_lay, ue_lay, uq_lay)
1690
1691 ! Accumuler les variables a stocker dans les fichiers histoire:
1692
1693 ! conversion Ec -> E thermique
1694 DO k = 1, llm
1695 DO i = 1, klon
1696 ZRCPD = RCPD * (1. + RVTMP2 * q_seri(i, k))
1697 d_t_ec(i, k) = 0.5 / ZRCPD &
1698 * (u(i, k)**2 + v(i, k)**2 - u_seri(i, k)**2 - v_seri(i, k)**2)
1699 t_seri(i, k) = t_seri(i, k) + d_t_ec(i, k)
1700 d_t_ec(i, k) = d_t_ec(i, k) / dtphys
1701 END DO
1702 END DO
1703
1704 IF (if_ebil >= 1) THEN
1705 ztit = 'after physic'
1706 CALL diagetpq(airephy, ztit, ip_ebil, 1, 1, dtphys, t_seri, q_seri, &
1707 ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &
1708 d_ql, d_qs, d_ec)
1709 ! Comme les tendances de la physique sont ajoute dans la dynamique,
1710 ! on devrait avoir que la variation d'entalpie par la dynamique
1711 ! est egale a la variation de la physique au pas de temps precedent.
1712 ! Donc la somme de ces 2 variations devrait etre nulle.
1713 call diagphy(airephy, ztit, ip_ebil, topsw, toplw, solsw, sollw, sens, &
1714 evap, rain_fall, snow_fall, ztsol, d_h_vcol, d_qt, d_ec, &
1715 fs_bound, fq_bound)
1716
1717 d_h_vcol_phy = d_h_vcol
1718
1719 END IF
1720
1721 ! SORTIES
1722
1723 !cc prw = eau precipitable
1724 DO i = 1, klon
1725 prw(i) = 0.
1726 DO k = 1, llm
1727 prw(i) = prw(i) + q_seri(i, k)*zmasse(i, k)
1728 ENDDO
1729 ENDDO
1730
1731 ! Convertir les incrementations en tendances
1732
1733 DO k = 1, llm
1734 DO i = 1, klon
1735 d_u(i, k) = (u_seri(i, k) - u(i, k)) / dtphys
1736 d_v(i, k) = (v_seri(i, k) - v(i, k)) / dtphys
1737 d_t(i, k) = (t_seri(i, k) - t(i, k)) / dtphys
1738 d_qx(i, k, ivap) = (q_seri(i, k) - qx(i, k, ivap)) / dtphys
1739 d_qx(i, k, iliq) = (ql_seri(i, k) - qx(i, k, iliq)) / dtphys
1740 ENDDO
1741 ENDDO
1742
1743 IF (nqmx >= 3) THEN
1744 DO iq = 3, nqmx
1745 DO k = 1, llm
1746 DO i = 1, klon
1747 d_qx(i, k, iq) = (tr_seri(i, k, iq-2) - qx(i, k, iq)) / dtphys
1748 ENDDO
1749 ENDDO
1750 ENDDO
1751 ENDIF
1752
1753 ! Sauvegarder les valeurs de t et q a la fin de la physique:
1754 DO k = 1, llm
1755 DO i = 1, klon
1756 t_ancien(i, k) = t_seri(i, k)
1757 q_ancien(i, k) = q_seri(i, k)
1758 ENDDO
1759 ENDDO
1760
1761 ! Ecriture des sorties
1762 call write_histhf
1763 call write_histday
1764 call write_histins
1765
1766 ! Si c'est la fin, il faut conserver l'etat de redemarrage
1767 IF (lafin) THEN
1768 itau_phy = itau_phy + itap
1769 CALL phyredem("restartphy.nc", rlat, rlon, pctsrf, ftsol, ftsoil, &
1770 tslab, seaice, fqsurf, qsol, fsnow, falbe, falblw, fevap, &
1771 rain_fall, snow_fall, solsw, sollwdown, dlw, radsol, frugs, &
1772 agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, t_ancien, &
1773 q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0)
1774 ENDIF
1775
1776 firstcal = .FALSE.
1777
1778 contains
1779
1780 subroutine write_histday
1781
1782 use gr_phy_write_3d_m, only: gr_phy_write_3d
1783 integer itau_w ! pas de temps ecriture
1784
1785 !------------------------------------------------
1786
1787 if (ok_journe) THEN
1788 itau_w = itau_phy + itap
1789 if (nqmx <= 4) then
1790 call histwrite(nid_day, "Sigma_O3_Royer", itau_w, &
1791 gr_phy_write_3d(wo) * 1e3)
1792 ! (convert "wo" from kDU to DU)
1793 end if
1794 if (ok_sync) then
1795 call histsync(nid_day)
1796 endif
1797 ENDIF
1798
1799 End subroutine write_histday
1800
1801 !****************************
1802
1803 subroutine write_histhf
1804
1805 ! From phylmd/write_histhf.h, version 1.5 2005/05/25 13:10:09
1806
1807 !------------------------------------------------
1808
1809 call write_histhf3d
1810
1811 IF (ok_sync) THEN
1812 call histsync(nid_hf)
1813 ENDIF
1814
1815 end subroutine write_histhf
1816
1817 !***************************************************************
1818
1819 subroutine write_histins
1820
1821 ! From phylmd/write_histins.h, version 1.2 2005/05/25 13:10:09
1822
1823 real zout
1824 integer itau_w ! pas de temps ecriture
1825
1826 !--------------------------------------------------
1827
1828 IF (ok_instan) THEN
1829 ! Champs 2D:
1830
1831 zsto = dtphys * ecrit_ins
1832 zout = dtphys * ecrit_ins
1833 itau_w = itau_phy + itap
1834
1835 i = NINT(zout/zsto)
1836 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, pphis, zx_tmp_2d)
1837 CALL histwrite(nid_ins, "phis", itau_w, zx_tmp_2d)
1838
1839 i = NINT(zout/zsto)
1840 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, airephy, zx_tmp_2d)
1841 CALL histwrite(nid_ins, "aire", itau_w, zx_tmp_2d)
1842
1843 DO i = 1, klon
1844 zx_tmp_fi2d(i) = paprs(i, 1)
1845 ENDDO
1846 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1847 CALL histwrite(nid_ins, "psol", itau_w, zx_tmp_2d)
1848
1849 DO i = 1, klon
1850 zx_tmp_fi2d(i) = rain_fall(i) + snow_fall(i)
1851 ENDDO
1852 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1853 CALL histwrite(nid_ins, "precip", itau_w, zx_tmp_2d)
1854
1855 DO i = 1, klon
1856 zx_tmp_fi2d(i) = rain_lsc(i) + snow_lsc(i)
1857 ENDDO
1858 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1859 CALL histwrite(nid_ins, "plul", itau_w, zx_tmp_2d)
1860
1861 DO i = 1, klon
1862 zx_tmp_fi2d(i) = rain_con(i) + snow_con(i)
1863 ENDDO
1864 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1865 CALL histwrite(nid_ins, "pluc", itau_w, zx_tmp_2d)
1866
1867 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zxtsol, zx_tmp_2d)
1868 CALL histwrite(nid_ins, "tsol", itau_w, zx_tmp_2d)
1869 !ccIM
1870 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zt2m, zx_tmp_2d)
1871 CALL histwrite(nid_ins, "t2m", itau_w, zx_tmp_2d)
1872
1873 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zq2m, zx_tmp_2d)
1874 CALL histwrite(nid_ins, "q2m", itau_w, zx_tmp_2d)
1875
1876 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zu10m, zx_tmp_2d)
1877 CALL histwrite(nid_ins, "u10m", itau_w, zx_tmp_2d)
1878
1879 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zv10m, zx_tmp_2d)
1880 CALL histwrite(nid_ins, "v10m", itau_w, zx_tmp_2d)
1881
1882 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, snow_fall, zx_tmp_2d)
1883 CALL histwrite(nid_ins, "snow", itau_w, zx_tmp_2d)
1884
1885 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, cdragm, zx_tmp_2d)
1886 CALL histwrite(nid_ins, "cdrm", itau_w, zx_tmp_2d)
1887
1888 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, cdragh, zx_tmp_2d)
1889 CALL histwrite(nid_ins, "cdrh", itau_w, zx_tmp_2d)
1890
1891 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, toplw, zx_tmp_2d)
1892 CALL histwrite(nid_ins, "topl", itau_w, zx_tmp_2d)
1893
1894 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, evap, zx_tmp_2d)
1895 CALL histwrite(nid_ins, "evap", itau_w, zx_tmp_2d)
1896
1897 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, solsw, zx_tmp_2d)
1898 CALL histwrite(nid_ins, "sols", itau_w, zx_tmp_2d)
1899
1900 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sollw, zx_tmp_2d)
1901 CALL histwrite(nid_ins, "soll", itau_w, zx_tmp_2d)
1902
1903 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sollwdown, zx_tmp_2d)
1904 CALL histwrite(nid_ins, "solldown", itau_w, zx_tmp_2d)
1905
1906 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, bils, zx_tmp_2d)
1907 CALL histwrite(nid_ins, "bils", itau_w, zx_tmp_2d)
1908
1909 zx_tmp_fi2d(1:klon) = -1*sens(1:klon)
1910 ! CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sens, zx_tmp_2d)
1911 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1912 CALL histwrite(nid_ins, "sens", itau_w, zx_tmp_2d)
1913
1914 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, fder, zx_tmp_2d)
1915 CALL histwrite(nid_ins, "fder", itau_w, zx_tmp_2d)
1916
1917 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_oce), zx_tmp_2d)
1918 CALL histwrite(nid_ins, "dtsvdfo", itau_w, zx_tmp_2d)
1919
1920 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_ter), zx_tmp_2d)
1921 CALL histwrite(nid_ins, "dtsvdft", itau_w, zx_tmp_2d)
1922
1923 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_lic), zx_tmp_2d)
1924 CALL histwrite(nid_ins, "dtsvdfg", itau_w, zx_tmp_2d)
1925
1926 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_sic), zx_tmp_2d)
1927 CALL histwrite(nid_ins, "dtsvdfi", itau_w, zx_tmp_2d)
1928
1929 DO nsrf = 1, nbsrf
1930 !XXX
1931 zx_tmp_fi2d(1 : klon) = pctsrf(1 : klon, nsrf)*100.
1932 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1933 CALL histwrite(nid_ins, "pourc_"//clnsurf(nsrf), itau_w, &
1934 zx_tmp_2d)
1935
1936 zx_tmp_fi2d(1 : klon) = pctsrf(1 : klon, nsrf)
1937 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1938 CALL histwrite(nid_ins, "fract_"//clnsurf(nsrf), itau_w, &
1939 zx_tmp_2d)
1940
1941 zx_tmp_fi2d(1 : klon) = fluxt(1 : klon, 1, nsrf)
1942 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1943 CALL histwrite(nid_ins, "sens_"//clnsurf(nsrf), itau_w, &
1944 zx_tmp_2d)
1945
1946 zx_tmp_fi2d(1 : klon) = fluxlat(1 : klon, nsrf)
1947 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1948 CALL histwrite(nid_ins, "lat_"//clnsurf(nsrf), itau_w, &
1949 zx_tmp_2d)
1950
1951 zx_tmp_fi2d(1 : klon) = ftsol(1 : klon, nsrf)
1952 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1953 CALL histwrite(nid_ins, "tsol_"//clnsurf(nsrf), itau_w, &
1954 zx_tmp_2d)
1955
1956 zx_tmp_fi2d(1 : klon) = fluxu(1 : klon, 1, nsrf)
1957 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1958 CALL histwrite(nid_ins, "taux_"//clnsurf(nsrf), itau_w, &
1959 zx_tmp_2d)
1960
1961 zx_tmp_fi2d(1 : klon) = fluxv(1 : klon, 1, nsrf)
1962 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1963 CALL histwrite(nid_ins, "tauy_"//clnsurf(nsrf), itau_w, &
1964 zx_tmp_2d)
1965
1966 zx_tmp_fi2d(1 : klon) = frugs(1 : klon, nsrf)
1967 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1968 CALL histwrite(nid_ins, "rugs_"//clnsurf(nsrf), itau_w, &
1969 zx_tmp_2d)
1970
1971 zx_tmp_fi2d(1 : klon) = falbe(1 : klon, nsrf)
1972 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1973 CALL histwrite(nid_ins, "albe_"//clnsurf(nsrf), itau_w, &
1974 zx_tmp_2d)
1975
1976 END DO
1977 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, albsol, zx_tmp_2d)
1978 CALL histwrite(nid_ins, "albs", itau_w, zx_tmp_2d)
1979 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, albsollw, zx_tmp_2d)
1980 CALL histwrite(nid_ins, "albslw", itau_w, zx_tmp_2d)
1981
1982 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zxrugs, zx_tmp_2d)
1983 CALL histwrite(nid_ins, "rugs", itau_w, zx_tmp_2d)
1984
1985 !HBTM2
1986
1987 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_pblh, zx_tmp_2d)
1988 CALL histwrite(nid_ins, "s_pblh", itau_w, zx_tmp_2d)
1989
1990 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_pblt, zx_tmp_2d)
1991 CALL histwrite(nid_ins, "s_pblt", itau_w, zx_tmp_2d)
1992
1993 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_lcl, zx_tmp_2d)
1994 CALL histwrite(nid_ins, "s_lcl", itau_w, zx_tmp_2d)
1995
1996 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_capCL, zx_tmp_2d)
1997 CALL histwrite(nid_ins, "s_capCL", itau_w, zx_tmp_2d)
1998
1999 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_oliqCL, zx_tmp_2d)
2000 CALL histwrite(nid_ins, "s_oliqCL", itau_w, zx_tmp_2d)
2001
2002 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_cteiCL, zx_tmp_2d)
2003 CALL histwrite(nid_ins, "s_cteiCL", itau_w, zx_tmp_2d)
2004
2005 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_therm, zx_tmp_2d)
2006 CALL histwrite(nid_ins, "s_therm", itau_w, zx_tmp_2d)
2007
2008 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb1, zx_tmp_2d)
2009 CALL histwrite(nid_ins, "s_trmb1", itau_w, zx_tmp_2d)
2010
2011 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb2, zx_tmp_2d)
2012 CALL histwrite(nid_ins, "s_trmb2", itau_w, zx_tmp_2d)
2013
2014 CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb3, zx_tmp_2d)
2015 CALL histwrite(nid_ins, "s_trmb3", itau_w, zx_tmp_2d)
2016
2017 ! Champs 3D:
2018
2019 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, t_seri, zx_tmp_3d)
2020 CALL histwrite(nid_ins, "temp", itau_w, zx_tmp_3d)
2021
2022 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, u_seri, zx_tmp_3d)
2023 CALL histwrite(nid_ins, "vitu", itau_w, zx_tmp_3d)
2024
2025 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, v_seri, zx_tmp_3d)
2026 CALL histwrite(nid_ins, "vitv", itau_w, zx_tmp_3d)
2027
2028 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, zphi, zx_tmp_3d)
2029 CALL histwrite(nid_ins, "geop", itau_w, zx_tmp_3d)
2030
2031 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, play, zx_tmp_3d)
2032 CALL histwrite(nid_ins, "pres", itau_w, zx_tmp_3d)
2033
2034 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, d_t_vdf, zx_tmp_3d)
2035 CALL histwrite(nid_ins, "dtvdf", itau_w, zx_tmp_3d)
2036
2037 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, d_q_vdf, zx_tmp_3d)
2038 CALL histwrite(nid_ins, "dqvdf", itau_w, zx_tmp_3d)
2039
2040 if (ok_sync) then
2041 call histsync(nid_ins)
2042 endif
2043 ENDIF
2044
2045 end subroutine write_histins
2046
2047 !****************************************************
2048
2049 subroutine write_histhf3d
2050
2051 ! From phylmd/write_histhf3d.h, version 1.2 2005/05/25 13:10:09
2052
2053 integer itau_w ! pas de temps ecriture
2054
2055 !-------------------------------------------------------
2056
2057 itau_w = itau_phy + itap
2058
2059 ! Champs 3D:
2060
2061 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, t_seri, zx_tmp_3d)
2062 CALL histwrite(nid_hf3d, "temp", itau_w, zx_tmp_3d)
2063
2064 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, qx(1, 1, ivap), zx_tmp_3d)
2065 CALL histwrite(nid_hf3d, "ovap", itau_w, zx_tmp_3d)
2066
2067 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, u_seri, zx_tmp_3d)
2068 CALL histwrite(nid_hf3d, "vitu", itau_w, zx_tmp_3d)
2069
2070 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, v_seri, zx_tmp_3d)
2071 CALL histwrite(nid_hf3d, "vitv", itau_w, zx_tmp_3d)
2072
2073 if (nbtr >= 3) then
2074 CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, tr_seri(1, 1, 3), &
2075 zx_tmp_3d)
2076 CALL histwrite(nid_hf3d, "O3", itau_w, zx_tmp_3d)
2077 end if
2078
2079 if (ok_sync) then
2080 call histsync(nid_hf3d)
2081 endif
2082
2083 end subroutine write_histhf3d
2084
2085 END SUBROUTINE physiq
2086
2087 end module physiq_m

  ViewVC Help
Powered by ViewVC 1.1.21