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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 49 - (show annotations)
Wed Aug 24 11:43:14 2011 UTC (12 years, 8 months ago) by guez
File size: 72229 byte(s)
LMDZE now uses library Jumble.

Removed all calls to "flinget". Replaced calls to "flinget",
"flininfo", "flinopen_nozoom" by calls to NetCDF95 and Jumble.

Split file "cv_driver.f" into "cv_driver.f90", "cv_flag.f90" and
"cv_thermo.f90".

Bug fix: "QANCIEN" was read twice in "phyeytat0".

In "physiq", initialization of "d_t", "d_u", "d_v" was useless.

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

  ViewVC Help
Powered by ViewVC 1.1.21