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

Contents of /trunk/phylmd/physiq.f

Parent Directory Parent Directory | Revision Log Revision Log

Revision 17 - (show annotations)
Tue Aug 5 13:31:32 2008 UTC (15 years, 9 months ago) by guez
Original Path: trunk/libf/phylmd/physiq.f90
File size: 76848 byte(s) 
Created rule for "compare_sampl_*" files in
"Documentation/Manuel_LMDZE.texfol/Graphiques/GNUmakefile".

Extracted "qcheck", "radiornpb", "minmaxqfi" into separate files.

Read pressure coordinate of ozone coefficients once per run instead of
every day.

Added some "intent" attributes.

Added argument "nq" to "ini_histday". Replaced calls to "gr_fi_ecrit"
by calls to "gr_phy_write_2d". "Sigma_O3_Royer" is written to
"histday.nc" only if "nq >= 4". Moved "ini_histrac" to module
"ini_hist".

Compute "zmasse" in "physiq", pass it to "phytrac".

Removed computations of "pftsol*" and "ppsrf*" in "phytrac".

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