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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 40 - (show annotations)
Tue Feb 22 13:49:36 2011 UTC (13 years, 2 months ago) by guez
File size: 75122 byte(s)
"alpha" useless, always 0, in "exner_hyb".

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

  ViewVC Help
Powered by ViewVC 1.1.21