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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 71 - (show annotations)
Mon Jul 8 18:12:18 2013 UTC (10 years, 9 months ago) by guez
File size: 68722 byte(s)
No reason to call inidissip in ce0l.

In inidissip, set random seed to 1 beacuse PGI compiler does not
accept all zeros.

dq was computed needlessly in caladvtrac. Arguments masse and dq of
calfis not used.

Replaced real*8 by double precision.

Pass arrays with inverted order of vertical levels to conflx instead
of creating local variables for this inside conflx.

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

  ViewVC Help
Powered by ViewVC 1.1.21