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

Contents of /trunk/phylmd/physiq.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 34 - (show annotations)
Wed Jun 2 11:01:12 2010 UTC (13 years, 11 months ago) by guez
Original Path: trunk/libf/phylmd/physiq.f90
File size: 75217 byte(s)
Split "ini_hist.f90" into single-procedure files.

In "calfis" and "physiq", removed dummy argument "nq" since "nq" must
be equal to "nqmx".

In "calfis", renamed dummy argument "pq" to "q", same name as actual
argument in "leapfrog". Renamed local variable "zqfi" to "qx", same
name as dummy argument in "physiq".

Removed arguments "itop_con" and "ibas_con" of "phytrac", which were
not used.

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

  ViewVC Help
Powered by ViewVC 1.1.21