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

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

Parent Directory Parent Directory | Revision Log Revision Log

Revision 35 - (hide annotations)
Tue Jun 8 15:37:21 2010 UTC (13 years, 11 months ago) by guez
File size: 75209 byte(s) 
Created intermediary variable for meridional wind in "calfis". Removed
unused variables.

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