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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 12 - (show annotations)
Mon Jul 21 16:05:07 2008 UTC (15 years, 9 months ago) by guez
File size: 92356 byte(s)
-- Minor modification of input/output:

Created procedure "read_logic". Variables of module "logic" are read
by "read_logic" instead of "conf_gcm". Variable "offline" of module
"conf_gcm" is read from namelist instead of "*.def".

Deleted arguments "dtime", "co2_ppm_etat0", "solaire_etat0",
"tabcntr0" and local variables "radpas", "tab_cntrl" of
"phyetat0". "phyetat0" does not read "controle" in "startphy.nc" any
longer. "phyetat0" now reads global attribute "itau_phy" from
"startphy.nc". "phyredem" does not create variable "controle" in
"startphy.nc" any longer. "phyredem" now writes global attribute
"itau_phy" of "startphy.nc". Deleted argument "tabcntr0" of
"printflag". Removed diagnostic messages written by "printflag" for
comparison of the variable "controle" of "startphy.nc" and the
variables read from "*.def" or namelist input.

-- Removing unwanted functionality:

Removed variable "lunout" from module "iniprint", replaced everywhere
by standard output.

Removed case "ocean == 'couple'" in "clmain", "interfsurf_hq" and
"physiq". Removed procedure "interfoce_cpl".

-- Should not change anything at run time:

Automated creation of graphs in documentation. More documentation on
input files.

Converted Fortran files to free format: "phyredem.f90", "printflag.f90".

Split module "clesphy" into "clesphys" and "clesphys2".

Removed variables "conser", "leapf", "forward", "apphys", "apdiss" and
"statcl" from module "logic". Added arguments "conser" to "advect",
"leapf" to "integrd". Added local variables "forward", "leapf",
"apphys", "conser", "apdiss" in "leapfrog".

Added intent attributes.

Deleted arguments "dtime" of "phyredem", "pdtime" of "flxdtdq", "sh"
of "phytrac", "dt" of "yamada".

Deleted local variables "dtime", "co2_ppm_etat0", "solaire_etat0",
"length", "tabcntr0" in "physiq". Replaced all references to "dtime"
by references to "pdtphys".

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

  ViewVC Help
Powered by ViewVC 1.1.21