/[lmdze]/trunk/phylmd/physiq.f
ViewVC logotype

Annotation of /trunk/phylmd/physiq.f

Parent Directory Parent Directory | Revision Log Revision Log

Revision 32 - (hide annotations)
Tue Apr 6 17:52:58 2010 UTC (14 years, 1 month ago) by guez
Original Path: trunk/libf/phylmd/physiq.f90
File size: 75235 byte(s) 
Split "stringop.f90" into single-procedure files. Gathered files in directory
"IOIPSL/Stringop".

Split "flincom.f90" into "flincom.f90" and "flinget.f90". Removed
unused procedures from module "flincom". Removed unused argument
"filename" of procedure "flinopen_nozoom".

Removed unused files.

Split "grid_change.f90" into "grid_change.f90" and
"gr_phy_write_3d.f90".

Removed unused procedures from modules "calendar", "ioipslmpp",
"grid_atob", "gath_cpl" and "getincom". Removed unused procedures in
files "ppm3d.f" and "thermcell.f".

Split "mathelp.f90" into "mathelp.f90" and "mathop.f90".

Removed unused variable "dpres" of module "comvert".

Use argument "itau" instead of local variables "iadvtr" and "first" to
control algorithm in procedure "fluxstokenc".

Removed unused arguments of procedure "integrd".

Removed useless computations at the end of procedure "leapfrog".

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