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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 49 - (hide annotations)
Wed Aug 24 11:43:14 2011 UTC (12 years, 8 months ago) by guez
File size: 72229 byte(s)
LMDZE now uses library Jumble.

Removed all calls to "flinget". Replaced calls to "flinget",
"flininfo", "flinopen_nozoom" by calls to NetCDF95 and Jumble.

Split file "cv_driver.f" into "cv_driver.f90", "cv_flag.f90" and
"cv_thermo.f90".

Bug fix: "QANCIEN" was read twice in "phyeytat0".

In "physiq", initialization of "d_t", "d_u", "d_v" was useless.

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

  ViewVC Help
Powered by ViewVC 1.1.21