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

Contents of /trunk/phylmd/physiq.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 30 - (show annotations)
Thu Apr 1 09:07:28 2010 UTC (14 years, 1 month ago) by guez
Original Path: trunk/libf/phylmd/physiq.f90
File size: 75385 byte(s)
Imported Source files of the external library "IOIPSL_Lionel" into
"libf/IOIPSL".

Split "cray.f90" into "scopy.f90" and "ssum.f90".

Rewrote "leapfrog" in order to have a clearer algorithmic structure.

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

  ViewVC Help
Powered by ViewVC 1.1.21