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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 47 - (hide annotations)
Fri Jul 1 15:00:48 2011 UTC (12 years, 10 months ago) by guez
File size: 72649 byte(s)
Split "thermcell.f" and "cv3_routines.f".
Removed copies of files that are now in "L_util".
Moved "mva9" and "diagetpq" to their own files.
Unified variable names across procedures.

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

  ViewVC Help
Powered by ViewVC 1.1.21