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

Annotation of /trunk/phylmd/physiq.f90

Parent Directory Parent Directory | Revision Log Revision Log


Revision 40 - (hide annotations)
Tue Feb 22 13:49:36 2011 UTC (13 years, 2 months ago) by guez
Original Path: trunk/libf/phylmd/physiq.f90
File size: 75122 byte(s)
"alpha" useless, always 0, in "exner_hyb".

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 38 use clmain_m, only: clmain
35 guez 3 use comgeomphy
36 guez 32 use conf_gcm_m, only: raz_date, offline
37     use conf_phys_m, only: conf_phys
38 guez 3 use ctherm
39 guez 34 use dimens_m, only: jjm, iim, llm, nqmx
40 guez 32 use dimphy, only: klon, nbtr
41     use dimsoil, only: nsoilmx
42     use hgardfou_m, only: hgardfou
43     USE histcom, only: histsync
44     USE histwrite_m, only: histwrite
45     use indicesol, only: nbsrf, is_ter, is_lic, is_sic, is_oce, &
46     clnsurf, epsfra
47 guez 34 use ini_histhf_m, only: ini_histhf
48     use ini_histday_m, only: ini_histday
49     use ini_histins_m, only: ini_histins
50 guez 32 use iniprint, only: prt_level
51 guez 3 use oasis_m
52     use orbite_m, only: orbite, zenang
53 guez 32 use ozonecm_m, only: ozonecm
54 guez 3 use phyetat0_m, only: phyetat0, rlat, rlon
55 guez 15 use phyredem_m, only: phyredem
56 guez 32 use phystokenc_m, only: phystokenc
57     use phytrac_m, only: phytrac
58 guez 17 use qcheck_m, only: qcheck
59 guez 32 use radepsi
60     use radopt
61     use temps, only: itau_phy, day_ref, annee_ref
62 guez 38 use yoethf_m
63     use SUPHEC_M, only: rcpd, rtt, rlvtt, rg, ra, rsigma, retv, romega
64 guez 3
65     ! Declaration des constantes et des fonctions thermodynamiques :
66     use fcttre, only: thermcep, foeew, qsats, qsatl
67    
68     ! Variables argument:
69    
70 guez 15 REAL, intent(in):: rdayvrai
71     ! (elapsed time since January 1st 0h of the starting year, in days)
72    
73 guez 3 REAL, intent(in):: gmtime ! heure de la journée en fraction de jour
74 guez 12 REAL, intent(in):: pdtphys ! pas d'integration pour la physique (seconde)
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 guez 35 REAL, intent(in):: v(klon, llm) ! vitesse Y (de S a N) en m/s
90 guez 3 REAL t(klon, llm) ! input temperature (K)
91    
92 guez 34 REAL, intent(in):: qx(klon, llm, nqmx)
93     ! (humidité spécifique 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 guez 34 REAL d_qx(klon, llm, nqmx) ! output tendance physique de "qx" (kg/kg/s)
100 guez 3 REAL d_ps(klon) ! output tendance physique de la pression au sol
101    
102 guez 35 LOGICAL:: firstcal = .true.
103    
104 guez 3 INTEGER nbteta
105     PARAMETER(nbteta=3)
106    
107     REAL PVteta(klon, nbteta)
108     ! (output vorticite potentielle a des thetas constantes)
109    
110     LOGICAL ok_cvl ! pour activer le nouveau driver pour convection KE
111     PARAMETER (ok_cvl=.TRUE.)
112     LOGICAL ok_gust ! pour activer l'effet des gust sur flux surface
113     PARAMETER (ok_gust=.FALSE.)
114    
115     LOGICAL check ! Verifier la conservation du modele en eau
116     PARAMETER (check=.FALSE.)
117     LOGICAL ok_stratus ! Ajouter artificiellement les stratus
118     PARAMETER (ok_stratus=.FALSE.)
119    
120     ! Parametres lies au coupleur OASIS:
121     INTEGER, SAVE :: npas, nexca
122     logical rnpb
123     parameter(rnpb=.true.)
124    
125 guez 12 character(len=6), save:: ocean
126     ! (type de modèle océan à utiliser: "force" ou "slab" mais pas "couple")
127    
128 guez 3 logical ok_ocean
129     SAVE ok_ocean
130    
131     !IM "slab" ocean
132     REAL tslab(klon) !Temperature du slab-ocean
133     SAVE tslab
134     REAL seaice(klon) !glace de mer (kg/m2)
135     SAVE seaice
136     REAL fluxo(klon) !flux turbulents ocean-glace de mer
137     REAL fluxg(klon) !flux turbulents ocean-atmosphere
138    
139     ! Modele thermique du sol, a activer pour le cycle diurne:
140 guez 17 logical, save:: ok_veget
141     LOGICAL, save:: ok_journe ! sortir le fichier journalier
142 guez 3
143     LOGICAL ok_mensuel ! sortir le fichier mensuel
144    
145     LOGICAL ok_instan ! sortir le fichier instantane
146     save ok_instan
147    
148     LOGICAL ok_region ! sortir le fichier regional
149     PARAMETER (ok_region=.FALSE.)
150    
151     ! pour phsystoke avec thermiques
152     REAL fm_therm(klon, llm+1)
153     REAL entr_therm(klon, llm)
154     real q2(klon, llm+1, nbsrf)
155     save q2
156    
157     INTEGER ivap ! indice de traceurs pour vapeur d'eau
158     PARAMETER (ivap=1)
159     INTEGER iliq ! indice de traceurs pour eau liquide
160     PARAMETER (iliq=2)
161    
162     REAL t_ancien(klon, llm), q_ancien(klon, llm)
163     SAVE t_ancien, q_ancien
164     LOGICAL ancien_ok
165     SAVE ancien_ok
166    
167     REAL d_t_dyn(klon, llm) ! tendance dynamique pour "t" (K/s)
168     REAL d_q_dyn(klon, llm) ! tendance dynamique pour "q" (kg/kg/s)
169    
170     real da(klon, llm), phi(klon, llm, llm), mp(klon, llm)
171    
172     !IM Amip2 PV a theta constante
173    
174     CHARACTER(LEN=3) ctetaSTD(nbteta)
175     DATA ctetaSTD/'350', '380', '405'/
176     REAL rtetaSTD(nbteta)
177     DATA rtetaSTD/350., 380., 405./
178    
179     !MI Amip2 PV a theta constante
180    
181     INTEGER klevp1
182     PARAMETER(klevp1=llm+1)
183    
184     REAL swdn0(klon, klevp1), swdn(klon, klevp1)
185     REAL swup0(klon, klevp1), swup(klon, klevp1)
186     SAVE swdn0, swdn, swup0, swup
187    
188     REAL lwdn0(klon, klevp1), lwdn(klon, klevp1)
189     REAL lwup0(klon, klevp1), lwup(klon, klevp1)
190     SAVE lwdn0, lwdn, lwup0, lwup
191    
192     !IM Amip2
193     ! variables a une pression donnee
194    
195     integer nlevSTD
196     PARAMETER(nlevSTD=17)
197     real rlevSTD(nlevSTD)
198     DATA rlevSTD/100000., 92500., 85000., 70000., &
199     60000., 50000., 40000., 30000., 25000., 20000., &
200     15000., 10000., 7000., 5000., 3000., 2000., 1000./
201     CHARACTER(LEN=4) clevSTD(nlevSTD)
202     DATA clevSTD/'1000', '925 ', '850 ', '700 ', '600 ', &
203     '500 ', '400 ', '300 ', '250 ', '200 ', '150 ', '100 ', &
204     '70 ', '50 ', '30 ', '20 ', '10 '/
205    
206     ! prw: precipitable water
207     real prw(klon)
208    
209     ! flwp, fiwp = Liquid Water Path & Ice Water Path (kg/m2)
210     ! flwc, fiwc = Liquid Water Content & Ice Water Content (kg/kg)
211     REAL flwp(klon), fiwp(klon)
212     REAL flwc(klon, llm), fiwc(klon, llm)
213    
214 guez 23 INTEGER kmax, lmax
215 guez 3 PARAMETER(kmax=8, lmax=8)
216     INTEGER kmaxm1, lmaxm1
217     PARAMETER(kmaxm1=kmax-1, lmaxm1=lmax-1)
218    
219     REAL zx_tau(kmaxm1), zx_pc(lmaxm1)
220     DATA zx_tau/0.0, 0.3, 1.3, 3.6, 9.4, 23., 60./
221     DATA zx_pc/50., 180., 310., 440., 560., 680., 800./
222    
223     ! cldtopres pression au sommet des nuages
224     REAL cldtopres(lmaxm1)
225     DATA cldtopres/50., 180., 310., 440., 560., 680., 800./
226    
227     ! taulev: numero du niveau de tau dans les sorties ISCCP
228     CHARACTER(LEN=4) taulev(kmaxm1)
229    
230     DATA taulev/'tau0', 'tau1', 'tau2', 'tau3', 'tau4', 'tau5', 'tau6'/
231     CHARACTER(LEN=3) pclev(lmaxm1)
232     DATA pclev/'pc1', 'pc2', 'pc3', 'pc4', 'pc5', 'pc6', 'pc7'/
233    
234     CHARACTER(LEN=28) cnameisccp(lmaxm1, kmaxm1)
235     DATA cnameisccp/'pc< 50hPa, tau< 0.3', 'pc= 50-180hPa, tau< 0.3', &
236     'pc= 180-310hPa, tau< 0.3', 'pc= 310-440hPa, tau< 0.3', &
237     'pc= 440-560hPa, tau< 0.3', 'pc= 560-680hPa, tau< 0.3', &
238     'pc= 680-800hPa, tau< 0.3', 'pc< 50hPa, tau= 0.3-1.3', &
239     'pc= 50-180hPa, tau= 0.3-1.3', 'pc= 180-310hPa, tau= 0.3-1.3', &
240     'pc= 310-440hPa, tau= 0.3-1.3', 'pc= 440-560hPa, tau= 0.3-1.3', &
241     'pc= 560-680hPa, tau= 0.3-1.3', 'pc= 680-800hPa, tau= 0.3-1.3', &
242     'pc< 50hPa, tau= 1.3-3.6', 'pc= 50-180hPa, tau= 1.3-3.6', &
243     'pc= 180-310hPa, tau= 1.3-3.6', 'pc= 310-440hPa, tau= 1.3-3.6', &
244     'pc= 440-560hPa, tau= 1.3-3.6', 'pc= 560-680hPa, tau= 1.3-3.6', &
245     'pc= 680-800hPa, tau= 1.3-3.6', 'pc< 50hPa, tau= 3.6-9.4', &
246     'pc= 50-180hPa, tau= 3.6-9.4', 'pc= 180-310hPa, tau= 3.6-9.4', &
247     'pc= 310-440hPa, tau= 3.6-9.4', 'pc= 440-560hPa, tau= 3.6-9.4', &
248     'pc= 560-680hPa, tau= 3.6-9.4', 'pc= 680-800hPa, tau= 3.6-9.4', &
249     'pc< 50hPa, tau= 9.4-23', 'pc= 50-180hPa, tau= 9.4-23', &
250     'pc= 180-310hPa, tau= 9.4-23', 'pc= 310-440hPa, tau= 9.4-23', &
251     'pc= 440-560hPa, tau= 9.4-23', 'pc= 560-680hPa, tau= 9.4-23', &
252     'pc= 680-800hPa, tau= 9.4-23', 'pc< 50hPa, tau= 23-60', &
253     'pc= 50-180hPa, tau= 23-60', 'pc= 180-310hPa, tau= 23-60', &
254     'pc= 310-440hPa, tau= 23-60', 'pc= 440-560hPa, tau= 23-60', &
255     'pc= 560-680hPa, tau= 23-60', 'pc= 680-800hPa, tau= 23-60', &
256     'pc< 50hPa, tau> 60.', 'pc= 50-180hPa, tau> 60.', &
257     'pc= 180-310hPa, tau> 60.', 'pc= 310-440hPa, tau> 60.', &
258     'pc= 440-560hPa, tau> 60.', 'pc= 560-680hPa, tau> 60.', &
259     'pc= 680-800hPa, tau> 60.'/
260    
261     !IM ISCCP simulator v3.4
262    
263     integer nid_hf, nid_hf3d
264     save nid_hf, nid_hf3d
265    
266     ! Variables propres a la physique
267    
268     INTEGER, save:: radpas
269     ! (Radiative transfer computations are made every "radpas" call to
270     ! "physiq".)
271    
272     REAL radsol(klon)
273     SAVE radsol ! bilan radiatif au sol calcule par code radiatif
274    
275 guez 7 INTEGER, SAVE:: itap ! number of calls to "physiq"
276 guez 3
277     REAL ftsol(klon, nbsrf)
278     SAVE ftsol ! temperature du sol
279    
280     REAL ftsoil(klon, nsoilmx, nbsrf)
281     SAVE ftsoil ! temperature dans le sol
282    
283     REAL fevap(klon, nbsrf)
284     SAVE fevap ! evaporation
285     REAL fluxlat(klon, nbsrf)
286     SAVE fluxlat
287    
288     REAL fqsurf(klon, nbsrf)
289     SAVE fqsurf ! humidite de l'air au contact de la surface
290    
291     REAL qsol(klon)
292     SAVE qsol ! hauteur d'eau dans le sol
293    
294     REAL fsnow(klon, nbsrf)
295     SAVE fsnow ! epaisseur neigeuse
296    
297     REAL falbe(klon, nbsrf)
298     SAVE falbe ! albedo par type de surface
299     REAL falblw(klon, nbsrf)
300     SAVE falblw ! albedo par type de surface
301    
302 guez 13 ! Paramètres de l'orographie à l'échelle sous-maille (OESM) :
303     REAL, save:: zmea(klon) ! orographie moyenne
304     REAL, save:: zstd(klon) ! deviation standard de l'OESM
305     REAL, save:: zsig(klon) ! pente de l'OESM
306     REAL, save:: zgam(klon) ! anisotropie de l'OESM
307     REAL, save:: zthe(klon) ! orientation de l'OESM
308     REAL, save:: zpic(klon) ! Maximum de l'OESM
309     REAL, save:: zval(klon) ! Minimum de l'OESM
310     REAL, save:: rugoro(klon) ! longueur de rugosite de l'OESM
311 guez 3
312     REAL zulow(klon), zvlow(klon)
313    
314     INTEGER igwd, idx(klon), itest(klon)
315    
316     REAL agesno(klon, nbsrf)
317     SAVE agesno ! age de la neige
318    
319     REAL run_off_lic_0(klon)
320     SAVE run_off_lic_0
321     !KE43
322     ! Variables liees a la convection de K. Emanuel (sb):
323    
324     REAL bas, top ! cloud base and top levels
325     SAVE bas
326     SAVE top
327    
328     REAL Ma(klon, llm) ! undilute upward mass flux
329     SAVE Ma
330     REAL qcondc(klon, llm) ! in-cld water content from convect
331     SAVE qcondc
332     REAL ema_work1(klon, llm), ema_work2(klon, llm)
333     SAVE ema_work1, ema_work2
334    
335     REAL wd(klon) ! sb
336     SAVE wd ! sb
337    
338     ! Variables locales pour la couche limite (al1):
339    
340     ! Variables locales:
341    
342     REAL cdragh(klon) ! drag coefficient pour T and Q
343     REAL cdragm(klon) ! drag coefficient pour vent
344    
345     !AA Pour phytrac
346     REAL ycoefh(klon, llm) ! coef d'echange pour phytrac
347     REAL yu1(klon) ! vents dans la premiere couche U
348     REAL yv1(klon) ! vents dans la premiere couche V
349     REAL ffonte(klon, nbsrf) !Flux thermique utilise pour fondre la neige
350     REAL fqcalving(klon, nbsrf) !Flux d'eau "perdue" par la surface
351     ! !et necessaire pour limiter la
352     ! !hauteur de neige, en kg/m2/s
353     REAL zxffonte(klon), zxfqcalving(klon)
354    
355     REAL pfrac_impa(klon, llm)! Produits des coefs lessivage impaction
356     save pfrac_impa
357     REAL pfrac_nucl(klon, llm)! Produits des coefs lessivage nucleation
358     save pfrac_nucl
359     REAL pfrac_1nucl(klon, llm)! Produits des coefs lessi nucl (alpha = 1)
360     save pfrac_1nucl
361     REAL frac_impa(klon, llm) ! fractions d'aerosols lessivees (impaction)
362     REAL frac_nucl(klon, llm) ! idem (nucleation)
363    
364     !AA
365     REAL rain_fall(klon) ! pluie
366     REAL snow_fall(klon) ! neige
367     save snow_fall, rain_fall
368     !IM cf FH pour Tiedtke 080604
369     REAL rain_tiedtke(klon), snow_tiedtke(klon)
370    
371     REAL evap(klon), devap(klon) ! evaporation et sa derivee
372     REAL sens(klon), dsens(klon) ! chaleur sensible et sa derivee
373     REAL dlw(klon) ! derivee infra rouge
374     SAVE dlw
375     REAL bils(klon) ! bilan de chaleur au sol
376     REAL fder(klon) ! Derive de flux (sensible et latente)
377     save fder
378     REAL ve(klon) ! integr. verticale du transport meri. de l'energie
379     REAL vq(klon) ! integr. verticale du transport meri. de l'eau
380     REAL ue(klon) ! integr. verticale du transport zonal de l'energie
381     REAL uq(klon) ! integr. verticale du transport zonal de l'eau
382    
383     REAL frugs(klon, nbsrf) ! longueur de rugosite
384     save frugs
385     REAL zxrugs(klon) ! longueur de rugosite
386    
387     ! Conditions aux limites
388    
389     INTEGER julien
390    
391 guez 7 INTEGER, SAVE:: lmt_pas ! number of time steps of "physics" per day
392 guez 3 REAL pctsrf(klon, nbsrf)
393     !IM
394     REAL pctsrf_new(klon, nbsrf) !pourcentage surfaces issus d'ORCHIDEE
395    
396     SAVE pctsrf ! sous-fraction du sol
397     REAL albsol(klon)
398     SAVE albsol ! albedo du sol total
399     REAL albsollw(klon)
400     SAVE albsollw ! albedo du sol total
401    
402 guez 17 REAL, SAVE:: wo(klon, llm) ! column density of ozone in a cell, in kDU
403 guez 3
404     ! Declaration des procedures appelees
405    
406     EXTERNAL alboc ! calculer l'albedo sur ocean
407     EXTERNAL ajsec ! ajustement sec
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 38 ! Couche limite:
1039 guez 3
1040 guez 40 CALL clmain(pdtphys, itap, date0, pctsrf, pctsrf_new, t_seri, q_seri, &
1041     u_seri, v_seri, julien, rmu0, co2_ppm, ok_veget, ocean, npas, nexca, &
1042     ftsol, soil_model, cdmmax, cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, &
1043     qsol, paprs, pplay, fsnow, fqsurf, fevap, falbe, falblw, fluxlat, &
1044     rain_fall, snow_fall, fsolsw, fsollw, sollwdown, fder, rlon, rlat, &
1045     cuphy, cvphy, frugs, firstcal, lafin, agesno, rugoro, d_t_vdf, &
1046     d_q_vdf, d_u_vdf, d_v_vdf, d_ts, fluxt, fluxq, fluxu, fluxv, cdragh, &
1047     cdragm, q2, dsens, devap, ycoefh, yu1, yv1, t2m, q2m, u10m, v10m, &
1048     pblh, capCL, oliqCL, cteiCL, pblT, therm, trmb1, trmb2, trmb3, plcl, &
1049     fqcalving, ffonte, run_off_lic_0, fluxo, fluxg, tslab, seaice)
1050 guez 3
1051 guez 40 ! Incrémentation des flux
1052    
1053 guez 3 zxfluxt=0.
1054     zxfluxq=0.
1055     zxfluxu=0.
1056     zxfluxv=0.
1057     DO nsrf = 1, nbsrf
1058     DO k = 1, llm
1059     DO i = 1, klon
1060     zxfluxt(i, k) = zxfluxt(i, k) + &
1061     fluxt(i, k, nsrf) * pctsrf( i, nsrf)
1062     zxfluxq(i, k) = zxfluxq(i, k) + &
1063     fluxq(i, k, nsrf) * pctsrf( i, nsrf)
1064     zxfluxu(i, k) = zxfluxu(i, k) + &
1065     fluxu(i, k, nsrf) * pctsrf( i, nsrf)
1066     zxfluxv(i, k) = zxfluxv(i, k) + &
1067     fluxv(i, k, nsrf) * pctsrf( i, nsrf)
1068     END DO
1069     END DO
1070     END DO
1071     DO i = 1, klon
1072     sens(i) = - zxfluxt(i, 1) ! flux de chaleur sensible au sol
1073     evap(i) = - zxfluxq(i, 1) ! flux d'evaporation au sol
1074     fder(i) = dlw(i) + dsens(i) + devap(i)
1075     ENDDO
1076    
1077     DO k = 1, llm
1078     DO i = 1, klon
1079     t_seri(i, k) = t_seri(i, k) + d_t_vdf(i, k)
1080     q_seri(i, k) = q_seri(i, k) + d_q_vdf(i, k)
1081     u_seri(i, k) = u_seri(i, k) + d_u_vdf(i, k)
1082     v_seri(i, k) = v_seri(i, k) + d_v_vdf(i, k)
1083     ENDDO
1084     ENDDO
1085    
1086     IF (if_ebil >= 2) THEN
1087     ztit='after clmain'
1088 guez 12 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, pdtphys &
1089 guez 10 , t_seri, q_seri, ql_seri, qs_seri, u_seri, v_seri, paprs &
1090 guez 3 , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
1091     call diagphy(airephy, ztit, ip_ebil &
1092     , zero_v, zero_v, zero_v, zero_v, sens &
1093     , evap, zero_v, zero_v, ztsol &
1094     , d_h_vcol, d_qt, d_ec &
1095     , fs_bound, fq_bound )
1096     END IF
1097    
1098     ! Incrementer la temperature du sol
1099    
1100     DO i = 1, klon
1101     zxtsol(i) = 0.0
1102     zxfluxlat(i) = 0.0
1103    
1104     zt2m(i) = 0.0
1105     zq2m(i) = 0.0
1106     zu10m(i) = 0.0
1107     zv10m(i) = 0.0
1108     zxffonte(i) = 0.0
1109     zxfqcalving(i) = 0.0
1110    
1111     s_pblh(i) = 0.0
1112     s_lcl(i) = 0.0
1113     s_capCL(i) = 0.0
1114     s_oliqCL(i) = 0.0
1115     s_cteiCL(i) = 0.0
1116     s_pblT(i) = 0.0
1117     s_therm(i) = 0.0
1118     s_trmb1(i) = 0.0
1119     s_trmb2(i) = 0.0
1120     s_trmb3(i) = 0.0
1121    
1122     IF ( abs( pctsrf(i, is_ter) + pctsrf(i, is_lic) + &
1123     pctsrf(i, is_oce) + pctsrf(i, is_sic) - 1.) .GT. EPSFRA) &
1124     THEN
1125     WRITE(*, *) 'physiq : pb sous surface au point ', i, &
1126     pctsrf(i, 1 : nbsrf)
1127     ENDIF
1128     ENDDO
1129     DO nsrf = 1, nbsrf
1130     DO i = 1, klon
1131     ftsol(i, nsrf) = ftsol(i, nsrf) + d_ts(i, nsrf)
1132     zxtsol(i) = zxtsol(i) + ftsol(i, nsrf)*pctsrf(i, nsrf)
1133     zxfluxlat(i) = zxfluxlat(i) + fluxlat(i, nsrf)*pctsrf(i, nsrf)
1134    
1135     zt2m(i) = zt2m(i) + t2m(i, nsrf)*pctsrf(i, nsrf)
1136     zq2m(i) = zq2m(i) + q2m(i, nsrf)*pctsrf(i, nsrf)
1137     zu10m(i) = zu10m(i) + u10m(i, nsrf)*pctsrf(i, nsrf)
1138     zv10m(i) = zv10m(i) + v10m(i, nsrf)*pctsrf(i, nsrf)
1139     zxffonte(i) = zxffonte(i) + ffonte(i, nsrf)*pctsrf(i, nsrf)
1140     zxfqcalving(i) = zxfqcalving(i) + &
1141     fqcalving(i, nsrf)*pctsrf(i, nsrf)
1142     s_pblh(i) = s_pblh(i) + pblh(i, nsrf)*pctsrf(i, nsrf)
1143     s_lcl(i) = s_lcl(i) + plcl(i, nsrf)*pctsrf(i, nsrf)
1144     s_capCL(i) = s_capCL(i) + capCL(i, nsrf) *pctsrf(i, nsrf)
1145     s_oliqCL(i) = s_oliqCL(i) + oliqCL(i, nsrf) *pctsrf(i, nsrf)
1146     s_cteiCL(i) = s_cteiCL(i) + cteiCL(i, nsrf) *pctsrf(i, nsrf)
1147     s_pblT(i) = s_pblT(i) + pblT(i, nsrf) *pctsrf(i, nsrf)
1148     s_therm(i) = s_therm(i) + therm(i, nsrf) *pctsrf(i, nsrf)
1149     s_trmb1(i) = s_trmb1(i) + trmb1(i, nsrf) *pctsrf(i, nsrf)
1150     s_trmb2(i) = s_trmb2(i) + trmb2(i, nsrf) *pctsrf(i, nsrf)
1151     s_trmb3(i) = s_trmb3(i) + trmb3(i, nsrf) *pctsrf(i, nsrf)
1152     ENDDO
1153     ENDDO
1154    
1155     ! Si une sous-fraction n'existe pas, elle prend la temp. moyenne
1156    
1157     DO nsrf = 1, nbsrf
1158     DO i = 1, klon
1159 guez 7 IF (pctsrf(i, nsrf) < epsfra) ftsol(i, nsrf) = zxtsol(i)
1160 guez 3
1161 guez 7 IF (pctsrf(i, nsrf) < epsfra) t2m(i, nsrf) = zt2m(i)
1162     IF (pctsrf(i, nsrf) < epsfra) q2m(i, nsrf) = zq2m(i)
1163     IF (pctsrf(i, nsrf) < epsfra) u10m(i, nsrf) = zu10m(i)
1164     IF (pctsrf(i, nsrf) < epsfra) v10m(i, nsrf) = zv10m(i)
1165     IF (pctsrf(i, nsrf) < epsfra) ffonte(i, nsrf) = zxffonte(i)
1166     IF (pctsrf(i, nsrf) < epsfra) &
1167 guez 3 fqcalving(i, nsrf) = zxfqcalving(i)
1168 guez 7 IF (pctsrf(i, nsrf) < epsfra) pblh(i, nsrf)=s_pblh(i)
1169     IF (pctsrf(i, nsrf) < epsfra) plcl(i, nsrf)=s_lcl(i)
1170     IF (pctsrf(i, nsrf) < epsfra) capCL(i, nsrf)=s_capCL(i)
1171     IF (pctsrf(i, nsrf) < epsfra) oliqCL(i, nsrf)=s_oliqCL(i)
1172     IF (pctsrf(i, nsrf) < epsfra) cteiCL(i, nsrf)=s_cteiCL(i)
1173     IF (pctsrf(i, nsrf) < epsfra) pblT(i, nsrf)=s_pblT(i)
1174     IF (pctsrf(i, nsrf) < epsfra) therm(i, nsrf)=s_therm(i)
1175     IF (pctsrf(i, nsrf) < epsfra) trmb1(i, nsrf)=s_trmb1(i)
1176     IF (pctsrf(i, nsrf) < epsfra) trmb2(i, nsrf)=s_trmb2(i)
1177     IF (pctsrf(i, nsrf) < epsfra) trmb3(i, nsrf)=s_trmb3(i)
1178 guez 3 ENDDO
1179     ENDDO
1180    
1181     ! Calculer la derive du flux infrarouge
1182    
1183     DO i = 1, klon
1184     dlw(i) = - 4.0*RSIGMA*zxtsol(i)**3
1185     ENDDO
1186    
1187     ! Appeler la convection (au choix)
1188    
1189     DO k = 1, llm
1190     DO i = 1, klon
1191     conv_q(i, k) = d_q_dyn(i, k) &
1192 guez 12 + d_q_vdf(i, k)/pdtphys
1193 guez 3 conv_t(i, k) = d_t_dyn(i, k) &
1194 guez 12 + d_t_vdf(i, k)/pdtphys
1195 guez 3 ENDDO
1196     ENDDO
1197     IF (check) THEN
1198     za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)
1199 guez 12 print *, "avantcon=", za
1200 guez 3 ENDIF
1201     zx_ajustq = .FALSE.
1202     IF (iflag_con == 2) zx_ajustq=.TRUE.
1203     IF (zx_ajustq) THEN
1204     DO i = 1, klon
1205     z_avant(i) = 0.0
1206     ENDDO
1207     DO k = 1, llm
1208     DO i = 1, klon
1209     z_avant(i) = z_avant(i) + (q_seri(i, k)+ql_seri(i, k)) &
1210 guez 17 *zmasse(i, k)
1211 guez 3 ENDDO
1212     ENDDO
1213     ENDIF
1214     IF (iflag_con == 1) THEN
1215     stop 'reactiver le call conlmd dans physiq.F'
1216     ELSE IF (iflag_con == 2) THEN
1217 guez 12 CALL conflx(pdtphys, paprs, pplay, t_seri, q_seri, &
1218 guez 3 conv_t, conv_q, zxfluxq(1, 1), omega, &
1219     d_t_con, d_q_con, rain_con, snow_con, &
1220     pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &
1221     kcbot, kctop, kdtop, pmflxr, pmflxs)
1222     WHERE (rain_con < 0.) rain_con = 0.
1223     WHERE (snow_con < 0.) snow_con = 0.
1224     DO i = 1, klon
1225     ibas_con(i) = llm+1 - kcbot(i)
1226     itop_con(i) = llm+1 - kctop(i)
1227     ENDDO
1228     ELSE IF (iflag_con >= 3) THEN
1229     ! nb of tracers for the KE convection:
1230     ! MAF la partie traceurs est faite dans phytrac
1231     ! on met ntra=1 pour limiter les appels mais on peut
1232     ! supprimer les calculs / ftra.
1233     ntra = 1
1234     ! Schema de convection modularise et vectorise:
1235     ! (driver commun aux versions 3 et 4)
1236    
1237     IF (ok_cvl) THEN ! new driver for convectL
1238 guez 18 CALL concvl(iflag_con, pdtphys, paprs, pplay, t_seri, q_seri, &
1239 guez 3 u_seri, v_seri, tr_seri, ntra, &
1240     ema_work1, ema_work2, &
1241     d_t_con, d_q_con, d_u_con, d_v_con, d_tr, &
1242     rain_con, snow_con, ibas_con, itop_con, &
1243     upwd, dnwd, dnwd0, &
1244     Ma, cape, tvp, iflagctrl, &
1245     pbase, bbase, dtvpdt1, dtvpdq1, dplcldt, dplcldr, qcondc, wd, &
1246     pmflxr, pmflxs, &
1247     da, phi, mp)
1248    
1249     clwcon0=qcondc
1250 guez 13 pmfu=upwd+dnwd
1251 guez 3 ELSE ! ok_cvl
1252     ! MAF conema3 ne contient pas les traceurs
1253 guez 17 CALL conema3 (pdtphys, paprs, pplay, t_seri, q_seri, &
1254 guez 3 u_seri, v_seri, tr_seri, ntra, &
1255     ema_work1, ema_work2, &
1256     d_t_con, d_q_con, d_u_con, d_v_con, d_tr, &
1257     rain_con, snow_con, ibas_con, itop_con, &
1258     upwd, dnwd, dnwd0, bas, top, &
1259     Ma, cape, tvp, rflag, &
1260     pbase &
1261     , bbase, dtvpdt1, dtvpdq1, dplcldt, dplcldr &
1262     , clwcon0)
1263     ENDIF ! ok_cvl
1264    
1265     IF (.NOT. ok_gust) THEN
1266     do i = 1, klon
1267     wd(i)=0.0
1268     enddo
1269     ENDIF
1270    
1271     ! Calcul des proprietes des nuages convectifs
1272    
1273     DO k = 1, llm
1274     DO i = 1, klon
1275     zx_t = t_seri(i, k)
1276     IF (thermcep) THEN
1277     zdelta = MAX(0., SIGN(1., rtt-zx_t))
1278     zx_qs = r2es * FOEEW(zx_t, zdelta)/pplay(i, k)
1279     zx_qs = MIN(0.5, zx_qs)
1280     zcor = 1./(1.-retv*zx_qs)
1281     zx_qs = zx_qs*zcor
1282     ELSE
1283 guez 7 IF (zx_t < t_coup) THEN
1284 guez 3 zx_qs = qsats(zx_t)/pplay(i, k)
1285     ELSE
1286     zx_qs = qsatl(zx_t)/pplay(i, k)
1287     ENDIF
1288     ENDIF
1289     zqsat(i, k)=zx_qs
1290     ENDDO
1291     ENDDO
1292    
1293     ! calcul des proprietes des nuages convectifs
1294 guez 13 clwcon0=fact_cldcon*clwcon0
1295 guez 3 call clouds_gno &
1296     (klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, rnebcon0)
1297     ELSE
1298 guez 12 print *, "iflag_con non-prevu", iflag_con
1299 guez 3 stop 1
1300     ENDIF
1301    
1302     DO k = 1, llm
1303     DO i = 1, klon
1304     t_seri(i, k) = t_seri(i, k) + d_t_con(i, k)
1305     q_seri(i, k) = q_seri(i, k) + d_q_con(i, k)
1306     u_seri(i, k) = u_seri(i, k) + d_u_con(i, k)
1307     v_seri(i, k) = v_seri(i, k) + d_v_con(i, k)
1308     ENDDO
1309     ENDDO
1310    
1311     IF (if_ebil >= 2) THEN
1312     ztit='after convect'
1313 guez 12 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, pdtphys &
1314 guez 10 , t_seri, q_seri, ql_seri, qs_seri, u_seri, v_seri, paprs &
1315 guez 3 , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
1316     call diagphy(airephy, ztit, ip_ebil &
1317     , zero_v, zero_v, zero_v, zero_v, zero_v &
1318     , zero_v, rain_con, snow_con, ztsol &
1319     , d_h_vcol, d_qt, d_ec &
1320     , fs_bound, fq_bound )
1321     END IF
1322    
1323     IF (check) THEN
1324     za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)
1325 guez 12 print *,"aprescon=", za
1326 guez 3 zx_t = 0.0
1327     za = 0.0
1328     DO i = 1, klon
1329     za = za + airephy(i)/REAL(klon)
1330     zx_t = zx_t + (rain_con(i)+ &
1331     snow_con(i))*airephy(i)/REAL(klon)
1332     ENDDO
1333 guez 12 zx_t = zx_t/za*pdtphys
1334     print *,"Precip=", zx_t
1335 guez 3 ENDIF
1336     IF (zx_ajustq) THEN
1337     DO i = 1, klon
1338     z_apres(i) = 0.0
1339     ENDDO
1340     DO k = 1, llm
1341     DO i = 1, klon
1342     z_apres(i) = z_apres(i) + (q_seri(i, k)+ql_seri(i, k)) &
1343 guez 17 *zmasse(i, k)
1344 guez 3 ENDDO
1345     ENDDO
1346     DO i = 1, klon
1347 guez 12 z_factor(i) = (z_avant(i)-(rain_con(i)+snow_con(i))*pdtphys) &
1348 guez 3 /z_apres(i)
1349     ENDDO
1350     DO k = 1, llm
1351     DO i = 1, klon
1352     IF (z_factor(i).GT.(1.0+1.0E-08) .OR. &
1353 guez 7 z_factor(i) < (1.0-1.0E-08)) THEN
1354 guez 3 q_seri(i, k) = q_seri(i, k) * z_factor(i)
1355     ENDIF
1356     ENDDO
1357     ENDDO
1358     ENDIF
1359     zx_ajustq=.FALSE.
1360    
1361     ! Convection seche (thermiques ou ajustement)
1362    
1363 guez 13 d_t_ajs=0.
1364     d_u_ajs=0.
1365     d_v_ajs=0.
1366     d_q_ajs=0.
1367     fm_therm=0.
1368     entr_therm=0.
1369 guez 3
1370 guez 12 IF(prt_level>9)print *, &
1371 guez 3 'AVANT LA CONVECTION SECHE, iflag_thermals=' &
1372     , iflag_thermals, ' nsplit_thermals=', nsplit_thermals
1373 guez 7 if(iflag_thermals < 0) then
1374 guez 3 ! Rien
1375 guez 12 IF(prt_level>9)print *,'pas de convection'
1376 guez 3 else if(iflag_thermals == 0) then
1377     ! Ajustement sec
1378 guez 12 IF(prt_level>9)print *,'ajsec'
1379 guez 3 CALL ajsec(paprs, pplay, t_seri, q_seri, d_t_ajs, d_q_ajs)
1380 guez 13 t_seri = t_seri + d_t_ajs
1381     q_seri = q_seri + d_q_ajs
1382 guez 3 else
1383     ! Thermiques
1384 guez 12 IF(prt_level>9)print *,'JUSTE AVANT, iflag_thermals=' &
1385 guez 3 , iflag_thermals, ' nsplit_thermals=', nsplit_thermals
1386     call calltherm(pdtphys &
1387     , pplay, paprs, pphi &
1388     , u_seri, v_seri, t_seri, q_seri &
1389     , d_u_ajs, d_v_ajs, d_t_ajs, d_q_ajs &
1390     , fm_therm, entr_therm)
1391     endif
1392    
1393     IF (if_ebil >= 2) THEN
1394     ztit='after dry_adjust'
1395 guez 12 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, pdtphys &
1396 guez 10 , t_seri, q_seri, ql_seri, qs_seri, u_seri, v_seri, paprs &
1397 guez 3 , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
1398     END IF
1399    
1400     ! Caclul des ratqs
1401    
1402     ! ratqs convectifs a l'ancienne en fonction de q(z=0)-q / q
1403     ! on ecrase le tableau ratqsc calcule par clouds_gno
1404     if (iflag_cldcon == 1) then
1405     do k=1, llm
1406     do i=1, klon
1407     if(ptconv(i, k)) then
1408     ratqsc(i, k)=ratqsbas &
1409     +fact_cldcon*(q_seri(i, 1)-q_seri(i, k))/q_seri(i, k)
1410     else
1411     ratqsc(i, k)=0.
1412     endif
1413     enddo
1414     enddo
1415     endif
1416    
1417     ! ratqs stables
1418     do k=1, llm
1419     do i=1, klon
1420     ratqss(i, k)=ratqsbas+(ratqshaut-ratqsbas)* &
1421     min((paprs(i, 1)-pplay(i, k))/(paprs(i, 1)-30000.), 1.)
1422     enddo
1423     enddo
1424    
1425     ! ratqs final
1426     if (iflag_cldcon == 1 .or.iflag_cldcon == 2) then
1427     ! les ratqs sont une conbinaison de ratqss et ratqsc
1428     ! ratqs final
1429     ! 1e4 (en gros 3 heures), en dur pour le moment, est le temps de
1430     ! relaxation des ratqs
1431     facteur=exp(-pdtphys*facttemps)
1432 guez 13 ratqs=max(ratqs*facteur, ratqss)
1433     ratqs=max(ratqs, ratqsc)
1434 guez 3 else
1435     ! on ne prend que le ratqs stable pour fisrtilp
1436 guez 13 ratqs=ratqss
1437 guez 3 endif
1438    
1439     ! Appeler le processus de condensation a grande echelle
1440     ! et le processus de precipitation
1441 guez 12 CALL fisrtilp(pdtphys, paprs, pplay, &
1442 guez 3 t_seri, q_seri, ptconv, ratqs, &
1443     d_t_lsc, d_q_lsc, d_ql_lsc, rneb, cldliq, &
1444     rain_lsc, snow_lsc, &
1445     pfrac_impa, pfrac_nucl, pfrac_1nucl, &
1446     frac_impa, frac_nucl, &
1447     prfl, psfl, rhcl)
1448    
1449     WHERE (rain_lsc < 0) rain_lsc = 0.
1450     WHERE (snow_lsc < 0) snow_lsc = 0.
1451     DO k = 1, llm
1452     DO i = 1, klon
1453     t_seri(i, k) = t_seri(i, k) + d_t_lsc(i, k)
1454     q_seri(i, k) = q_seri(i, k) + d_q_lsc(i, k)
1455     ql_seri(i, k) = ql_seri(i, k) + d_ql_lsc(i, k)
1456     cldfra(i, k) = rneb(i, k)
1457     IF (.NOT.new_oliq) cldliq(i, k) = ql_seri(i, k)
1458     ENDDO
1459     ENDDO
1460     IF (check) THEN
1461     za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)
1462 guez 12 print *,"apresilp=", za
1463 guez 3 zx_t = 0.0
1464     za = 0.0
1465     DO i = 1, klon
1466     za = za + airephy(i)/REAL(klon)
1467     zx_t = zx_t + (rain_lsc(i) &
1468     + snow_lsc(i))*airephy(i)/REAL(klon)
1469     ENDDO
1470 guez 12 zx_t = zx_t/za*pdtphys
1471     print *,"Precip=", zx_t
1472 guez 3 ENDIF
1473    
1474     IF (if_ebil >= 2) THEN
1475     ztit='after fisrt'
1476 guez 12 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, pdtphys &
1477 guez 10 , t_seri, q_seri, ql_seri, qs_seri, u_seri, v_seri, paprs &
1478 guez 3 , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
1479     call diagphy(airephy, ztit, ip_ebil &
1480     , zero_v, zero_v, zero_v, zero_v, zero_v &
1481     , zero_v, rain_lsc, snow_lsc, ztsol &
1482     , d_h_vcol, d_qt, d_ec &
1483     , fs_bound, fq_bound )
1484     END IF
1485    
1486     ! PRESCRIPTION DES NUAGES POUR LE RAYONNEMENT
1487    
1488     ! 1. NUAGES CONVECTIFS
1489    
1490     IF (iflag_cldcon.le.-1) THEN ! seulement pour Tiedtke
1491     snow_tiedtke=0.
1492     if (iflag_cldcon == -1) then
1493     rain_tiedtke=rain_con
1494     else
1495     rain_tiedtke=0.
1496     do k=1, llm
1497     do i=1, klon
1498 guez 7 if (d_q_con(i, k) < 0.) then
1499 guez 3 rain_tiedtke(i)=rain_tiedtke(i)-d_q_con(i, k)/pdtphys &
1500 guez 17 *zmasse(i, k)
1501 guez 3 endif
1502     enddo
1503     enddo
1504     endif
1505    
1506     ! Nuages diagnostiques pour Tiedtke
1507     CALL diagcld1(paprs, pplay, &
1508     rain_tiedtke, snow_tiedtke, ibas_con, itop_con, &
1509     diafra, dialiq)
1510     DO k = 1, llm
1511     DO i = 1, klon
1512     IF (diafra(i, k).GT.cldfra(i, k)) THEN
1513     cldliq(i, k) = dialiq(i, k)
1514     cldfra(i, k) = diafra(i, k)
1515     ENDIF
1516     ENDDO
1517     ENDDO
1518    
1519     ELSE IF (iflag_cldcon == 3) THEN
1520 guez 7 ! On prend pour les nuages convectifs le max du calcul de la
1521     ! convection et du calcul du pas de temps précédent diminué d'un facteur
1522     ! facttemps
1523 guez 3 facteur = pdtphys *facttemps
1524     do k=1, llm
1525     do i=1, klon
1526     rnebcon(i, k)=rnebcon(i, k)*facteur
1527     if (rnebcon0(i, k)*clwcon0(i, k).gt.rnebcon(i, k)*clwcon(i, k)) &
1528     then
1529     rnebcon(i, k)=rnebcon0(i, k)
1530     clwcon(i, k)=clwcon0(i, k)
1531     endif
1532     enddo
1533     enddo
1534    
1535     ! On prend la somme des fractions nuageuses et des contenus en eau
1536 guez 13 cldfra=min(max(cldfra, rnebcon), 1.)
1537     cldliq=cldliq+rnebcon*clwcon
1538 guez 3
1539     ENDIF
1540    
1541     ! 2. NUAGES STARTIFORMES
1542    
1543     IF (ok_stratus) THEN
1544     CALL diagcld2(paprs, pplay, t_seri, q_seri, diafra, dialiq)
1545     DO k = 1, llm
1546     DO i = 1, klon
1547     IF (diafra(i, k).GT.cldfra(i, k)) THEN
1548     cldliq(i, k) = dialiq(i, k)
1549     cldfra(i, k) = diafra(i, k)
1550     ENDIF
1551     ENDDO
1552     ENDDO
1553     ENDIF
1554    
1555     ! Precipitation totale
1556    
1557     DO i = 1, klon
1558     rain_fall(i) = rain_con(i) + rain_lsc(i)
1559     snow_fall(i) = snow_con(i) + snow_lsc(i)
1560     ENDDO
1561    
1562     IF (if_ebil >= 2) THEN
1563     ztit="after diagcld"
1564 guez 12 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, pdtphys &
1565 guez 10 , t_seri, q_seri, ql_seri, qs_seri, u_seri, v_seri, paprs &
1566 guez 3 , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
1567     END IF
1568    
1569     ! Calculer l'humidite relative pour diagnostique
1570    
1571     DO k = 1, llm
1572     DO i = 1, klon
1573     zx_t = t_seri(i, k)
1574     IF (thermcep) THEN
1575     zdelta = MAX(0., SIGN(1., rtt-zx_t))
1576     zx_qs = r2es * FOEEW(zx_t, zdelta)/pplay(i, k)
1577     zx_qs = MIN(0.5, zx_qs)
1578     zcor = 1./(1.-retv*zx_qs)
1579     zx_qs = zx_qs*zcor
1580     ELSE
1581 guez 7 IF (zx_t < t_coup) THEN
1582 guez 3 zx_qs = qsats(zx_t)/pplay(i, k)
1583     ELSE
1584     zx_qs = qsatl(zx_t)/pplay(i, k)
1585     ENDIF
1586     ENDIF
1587     zx_rh(i, k) = q_seri(i, k)/zx_qs
1588     zqsat(i, k)=zx_qs
1589     ENDDO
1590     ENDDO
1591     !jq - introduce the aerosol direct and first indirect radiative forcings
1592     !jq - Johannes Quaas, 27/11/2003 (quaas@lmd.jussieu.fr)
1593     IF (ok_ade.OR.ok_aie) THEN
1594     ! Get sulfate aerosol distribution
1595 guez 7 CALL readsulfate(rdayvrai, firstcal, sulfate)
1596     CALL readsulfate_preind(rdayvrai, firstcal, sulfate_pi)
1597 guez 3
1598     ! Calculate aerosol optical properties (Olivier Boucher)
1599     CALL aeropt(pplay, paprs, t_seri, sulfate, rhcl, &
1600     tau_ae, piz_ae, cg_ae, aerindex)
1601     ELSE
1602     tau_ae(:, :, :)=0.0
1603     piz_ae(:, :, :)=0.0
1604     cg_ae(:, :, :)=0.0
1605     ENDIF
1606    
1607     ! Calculer les parametres optiques des nuages et quelques
1608     ! parametres pour diagnostiques:
1609    
1610     if (ok_newmicro) then
1611     CALL newmicro (paprs, pplay, ok_newmicro, &
1612     t_seri, cldliq, cldfra, cldtau, cldemi, &
1613     cldh, cldl, cldm, cldt, cldq, &
1614     flwp, fiwp, flwc, fiwc, &
1615     ok_aie, &
1616     sulfate, sulfate_pi, &
1617     bl95_b0, bl95_b1, &
1618     cldtaupi, re, fl)
1619     else
1620     CALL nuage (paprs, pplay, &
1621     t_seri, cldliq, cldfra, cldtau, cldemi, &
1622     cldh, cldl, cldm, cldt, cldq, &
1623     ok_aie, &
1624     sulfate, sulfate_pi, &
1625     bl95_b0, bl95_b1, &
1626     cldtaupi, re, fl)
1627    
1628     endif
1629    
1630     ! Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.
1631    
1632     IF (MOD(itaprad, radpas) == 0) THEN
1633     DO i = 1, klon
1634     albsol(i) = falbe(i, is_oce) * pctsrf(i, is_oce) &
1635     + falbe(i, is_lic) * pctsrf(i, is_lic) &
1636     + falbe(i, is_ter) * pctsrf(i, is_ter) &
1637     + falbe(i, is_sic) * pctsrf(i, is_sic)
1638     albsollw(i) = falblw(i, is_oce) * pctsrf(i, is_oce) &
1639     + falblw(i, is_lic) * pctsrf(i, is_lic) &
1640     + falblw(i, is_ter) * pctsrf(i, is_ter) &
1641     + falblw(i, is_sic) * pctsrf(i, is_sic)
1642     ENDDO
1643     ! nouveau rayonnement (compatible Arpege-IFS):
1644     CALL radlwsw(dist, rmu0, fract, &
1645     paprs, pplay, zxtsol, albsol, albsollw, t_seri, q_seri, &
1646     wo, &
1647     cldfra, cldemi, cldtau, &
1648     heat, heat0, cool, cool0, radsol, albpla, &
1649     topsw, toplw, solsw, sollw, &
1650     sollwdown, &
1651     topsw0, toplw0, solsw0, sollw0, &
1652     lwdn0, lwdn, lwup0, lwup, &
1653     swdn0, swdn, swup0, swup, &
1654     ok_ade, ok_aie, & ! new for aerosol radiative effects
1655     tau_ae, piz_ae, cg_ae, &
1656     topswad, solswad, &
1657     cldtaupi, &
1658     topswai, solswai)
1659     itaprad = 0
1660     ENDIF
1661     itaprad = itaprad + 1
1662    
1663     ! Ajouter la tendance des rayonnements (tous les pas)
1664    
1665     DO k = 1, llm
1666     DO i = 1, klon
1667     t_seri(i, k) = t_seri(i, k) &
1668 guez 12 + (heat(i, k)-cool(i, k)) * pdtphys/86400.
1669 guez 3 ENDDO
1670     ENDDO
1671    
1672     IF (if_ebil >= 2) THEN
1673     ztit='after rad'
1674 guez 12 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, pdtphys &
1675 guez 10 , t_seri, q_seri, ql_seri, qs_seri, u_seri, v_seri, paprs &
1676 guez 3 , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
1677     call diagphy(airephy, ztit, ip_ebil &
1678     , topsw, toplw, solsw, sollw, zero_v &
1679     , zero_v, zero_v, zero_v, ztsol &
1680     , d_h_vcol, d_qt, d_ec &
1681     , fs_bound, fq_bound )
1682     END IF
1683    
1684     ! Calculer l'hydrologie de la surface
1685    
1686     DO i = 1, klon
1687     zxqsurf(i) = 0.0
1688     zxsnow(i) = 0.0
1689     ENDDO
1690     DO nsrf = 1, nbsrf
1691     DO i = 1, klon
1692     zxqsurf(i) = zxqsurf(i) + fqsurf(i, nsrf)*pctsrf(i, nsrf)
1693     zxsnow(i) = zxsnow(i) + fsnow(i, nsrf)*pctsrf(i, nsrf)
1694     ENDDO
1695     ENDDO
1696    
1697     ! Calculer le bilan du sol et la derive de temperature (couplage)
1698    
1699     DO i = 1, klon
1700     bils(i) = radsol(i) - sens(i) + zxfluxlat(i)
1701     ENDDO
1702    
1703 guez 13 !mod deb lott(jan95)
1704 guez 3 ! Appeler le programme de parametrisation de l'orographie
1705     ! a l'echelle sous-maille:
1706    
1707     IF (ok_orodr) THEN
1708     ! selection des points pour lesquels le shema est actif:
1709     igwd=0
1710     DO i=1, klon
1711     itest(i)=0
1712     IF (((zpic(i)-zmea(i)).GT.100.).AND.(zstd(i).GT.10.0)) THEN
1713     itest(i)=1
1714     igwd=igwd+1
1715     idx(igwd)=i
1716     ENDIF
1717     ENDDO
1718    
1719 guez 12 CALL drag_noro(klon, llm, pdtphys, paprs, pplay, &
1720 guez 3 zmea, zstd, zsig, zgam, zthe, zpic, zval, &
1721     igwd, idx, itest, &
1722     t_seri, u_seri, v_seri, &
1723     zulow, zvlow, zustrdr, zvstrdr, &
1724     d_t_oro, d_u_oro, d_v_oro)
1725    
1726     ! ajout des tendances
1727     DO k = 1, llm
1728     DO i = 1, klon
1729     t_seri(i, k) = t_seri(i, k) + d_t_oro(i, k)
1730     u_seri(i, k) = u_seri(i, k) + d_u_oro(i, k)
1731     v_seri(i, k) = v_seri(i, k) + d_v_oro(i, k)
1732     ENDDO
1733     ENDDO
1734 guez 13 ENDIF
1735 guez 3
1736     IF (ok_orolf) THEN
1737    
1738     ! selection des points pour lesquels le shema est actif:
1739     igwd=0
1740     DO i=1, klon
1741     itest(i)=0
1742     IF ((zpic(i)-zmea(i)).GT.100.) THEN
1743     itest(i)=1
1744     igwd=igwd+1
1745     idx(igwd)=i
1746     ENDIF
1747     ENDDO
1748    
1749 guez 12 CALL lift_noro(klon, llm, pdtphys, paprs, pplay, &
1750 guez 3 rlat, zmea, zstd, zpic, &
1751     itest, &
1752     t_seri, u_seri, v_seri, &
1753     zulow, zvlow, zustrli, zvstrli, &
1754     d_t_lif, d_u_lif, d_v_lif)
1755    
1756     ! ajout des tendances
1757     DO k = 1, llm
1758     DO i = 1, klon
1759     t_seri(i, k) = t_seri(i, k) + d_t_lif(i, k)
1760     u_seri(i, k) = u_seri(i, k) + d_u_lif(i, k)
1761     v_seri(i, k) = v_seri(i, k) + d_v_lif(i, k)
1762     ENDDO
1763     ENDDO
1764    
1765     ENDIF ! fin de test sur ok_orolf
1766    
1767     ! STRESS NECESSAIRES: TOUTE LA PHYSIQUE
1768    
1769     DO i = 1, klon
1770     zustrph(i)=0.
1771     zvstrph(i)=0.
1772     ENDDO
1773     DO k = 1, llm
1774     DO i = 1, klon
1775 guez 17 zustrph(i)=zustrph(i)+(u_seri(i, k)-u(i, k))/pdtphys* zmasse(i, k)
1776     zvstrph(i)=zvstrph(i)+(v_seri(i, k)-v(i, k))/pdtphys* zmasse(i, k)
1777 guez 3 ENDDO
1778     ENDDO
1779    
1780     !IM calcul composantes axiales du moment angulaire et couple des montagnes
1781    
1782 guez 17 CALL aaam_bud(27, klon, llm, gmtime, &
1783 guez 3 ra, rg, romega, &
1784     rlat, rlon, pphis, &
1785     zustrdr, zustrli, zustrph, &
1786     zvstrdr, zvstrli, zvstrph, &
1787     paprs, u, v, &
1788     aam, torsfc)
1789    
1790     IF (if_ebil >= 2) THEN
1791     ztit='after orography'
1792 guez 12 CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, pdtphys &
1793 guez 10 , t_seri, q_seri, ql_seri, qs_seri, u_seri, v_seri, paprs &
1794 guez 3 , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
1795     END IF
1796    
1797 guez 31 ! Calcul des tendances traceurs
1798 guez 35 call phytrac(rnpb, itap, lmt_pas, julien, gmtime, firstcal, lafin, &
1799     nqmx-2, pdtphys, u, t, paprs, pplay, pmfu, pmfd, pen_u, pde_u, &
1800     pen_d, pde_d, ycoefh, fm_therm, entr_therm, yu1, yv1, ftsol, pctsrf, &
1801     frac_impa, frac_nucl, pphis, pphi, albsol, rhcl, cldfra, rneb, &
1802 guez 34 diafra, cldliq, pmflxr, pmflxs, prfl, psfl, da, phi, mp, upwd, dnwd, &
1803     tr_seri, zmasse)
1804 guez 3
1805     IF (offline) THEN
1806 guez 31 call phystokenc(pdtphys, rlon, rlat, t, pmfu, pmfd, pen_u, pde_u, &
1807     pen_d, pde_d, fm_therm, entr_therm, ycoefh, yu1, yv1, ftsol, &
1808     pctsrf, frac_impa, frac_nucl, pphis, airephy, pdtphys, itap)
1809 guez 3 ENDIF
1810    
1811     ! Calculer le transport de l'eau et de l'energie (diagnostique)
1812 guez 31 CALL transp(paprs, zxtsol, t_seri, q_seri, u_seri, v_seri, zphi, ve, vq, &
1813     ue, uq)
1814 guez 3
1815 guez 31 ! diag. bilKP
1816 guez 3
1817     CALL transp_lay (paprs, zxtsol, &
1818     t_seri, q_seri, u_seri, v_seri, zphi, &
1819     ve_lay, vq_lay, ue_lay, uq_lay)
1820    
1821     ! Accumuler les variables a stocker dans les fichiers histoire:
1822    
1823     !+jld ec_conser
1824     DO k = 1, llm
1825     DO i = 1, klon
1826     ZRCPD = RCPD*(1.0+RVTMP2*q_seri(i, k))
1827     d_t_ec(i, k)=0.5/ZRCPD &
1828     *(u(i, k)**2+v(i, k)**2-u_seri(i, k)**2-v_seri(i, k)**2)
1829     t_seri(i, k)=t_seri(i, k)+d_t_ec(i, k)
1830 guez 12 d_t_ec(i, k) = d_t_ec(i, k)/pdtphys
1831 guez 3 END DO
1832     END DO
1833     !-jld ec_conser
1834     IF (if_ebil >= 1) THEN
1835     ztit='after physic'
1836 guez 12 CALL diagetpq(airephy, ztit, ip_ebil, 1, 1, pdtphys &
1837 guez 10 , t_seri, q_seri, ql_seri, qs_seri, u_seri, v_seri, paprs &
1838 guez 3 , d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)
1839     ! Comme les tendances de la physique sont ajoute dans la dynamique,
1840     ! on devrait avoir que la variation d'entalpie par la dynamique
1841     ! est egale a la variation de la physique au pas de temps precedent.
1842     ! Donc la somme de ces 2 variations devrait etre nulle.
1843     call diagphy(airephy, ztit, ip_ebil &
1844     , topsw, toplw, solsw, sollw, sens &
1845     , evap, rain_fall, snow_fall, ztsol &
1846     , d_h_vcol, d_qt, d_ec &
1847     , fs_bound, fq_bound )
1848    
1849     d_h_vcol_phy=d_h_vcol
1850    
1851     END IF
1852    
1853     ! SORTIES
1854    
1855     !cc prw = eau precipitable
1856     DO i = 1, klon
1857     prw(i) = 0.
1858     DO k = 1, llm
1859 guez 17 prw(i) = prw(i) + q_seri(i, k)*zmasse(i, k)
1860 guez 3 ENDDO
1861     ENDDO
1862    
1863     ! Convertir les incrementations en tendances
1864    
1865     DO k = 1, llm
1866     DO i = 1, klon
1867 guez 12 d_u(i, k) = ( u_seri(i, k) - u(i, k) ) / pdtphys
1868     d_v(i, k) = ( v_seri(i, k) - v(i, k) ) / pdtphys
1869     d_t(i, k) = ( t_seri(i, k)-t(i, k) ) / pdtphys
1870     d_qx(i, k, ivap) = ( q_seri(i, k) - qx(i, k, ivap) ) / pdtphys
1871     d_qx(i, k, iliq) = ( ql_seri(i, k) - qx(i, k, iliq) ) / pdtphys
1872 guez 3 ENDDO
1873     ENDDO
1874    
1875 guez 34 IF (nqmx >= 3) THEN
1876     DO iq = 3, nqmx
1877 guez 3 DO k = 1, llm
1878     DO i = 1, klon
1879 guez 13 d_qx(i, k, iq) = (tr_seri(i, k, iq-2) - qx(i, k, iq)) / pdtphys
1880 guez 3 ENDDO
1881     ENDDO
1882     ENDDO
1883     ENDIF
1884    
1885     ! Sauvegarder les valeurs de t et q a la fin de la physique:
1886     DO k = 1, llm
1887     DO i = 1, klon
1888     t_ancien(i, k) = t_seri(i, k)
1889     q_ancien(i, k) = q_seri(i, k)
1890     ENDDO
1891     ENDDO
1892    
1893     ! Ecriture des sorties
1894     call write_histhf
1895     call write_histday
1896     call write_histins
1897    
1898     ! Si c'est la fin, il faut conserver l'etat de redemarrage
1899     IF (lafin) THEN
1900     itau_phy = itau_phy + itap
1901 guez 13 CALL phyredem("restartphy.nc", rlat, rlon, pctsrf, ftsol, &
1902 guez 12 ftsoil, tslab, seaice, fqsurf, qsol, &
1903 guez 3 fsnow, falbe, falblw, fevap, rain_fall, snow_fall, &
1904     solsw, sollwdown, dlw, &
1905     radsol, frugs, agesno, &
1906 guez 13 zmea, zstd, zsig, zgam, zthe, zpic, zval, &
1907 guez 3 t_ancien, q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0)
1908     ENDIF
1909    
1910 guez 35 firstcal = .FALSE.
1911    
1912 guez 3 contains
1913    
1914 guez 15 subroutine write_histday
1915 guez 3
1916 guez 32 use gr_phy_write_3d_m, only: gr_phy_write_3d
1917 guez 17 integer itau_w ! pas de temps ecriture
1918 guez 3
1919 guez 15 !------------------------------------------------
1920 guez 3
1921     if (ok_journe) THEN
1922     itau_w = itau_phy + itap
1923 guez 34 if (nqmx <= 4) then
1924 guez 17 call histwrite(nid_day, "Sigma_O3_Royer", itau_w, &
1925     gr_phy_write_3d(wo) * 1e3)
1926     ! (convert "wo" from kDU to DU)
1927     end if
1928 guez 3 if (ok_sync) then
1929     call histsync(nid_day)
1930     endif
1931     ENDIF
1932    
1933     End subroutine write_histday
1934    
1935     !****************************
1936    
1937     subroutine write_histhf
1938    
1939     ! From phylmd/write_histhf.h, v 1.5 2005/05/25 13:10:09
1940    
1941 guez 17 !------------------------------------------------
1942 guez 3
1943     call write_histhf3d
1944    
1945     IF (ok_sync) THEN
1946     call histsync(nid_hf)
1947     ENDIF
1948    
1949     end subroutine write_histhf
1950    
1951     !***************************************************************
1952    
1953     subroutine write_histins
1954    
1955     ! From phylmd/write_histins.h, v 1.2 2005/05/25 13:10:09
1956    
1957     real zout
1958 guez 17 integer itau_w ! pas de temps ecriture
1959 guez 3
1960     !--------------------------------------------------
1961    
1962     IF (ok_instan) THEN
1963     ! Champs 2D:
1964    
1965 guez 12 zsto = pdtphys * ecrit_ins
1966     zout = pdtphys * ecrit_ins
1967 guez 3 itau_w = itau_phy + itap
1968    
1969     i = NINT(zout/zsto)
1970     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), pphis, zx_tmp_2d)
1971 guez 15 CALL histwrite(nid_ins, "phis", itau_w, zx_tmp_2d)
1972 guez 3
1973     i = NINT(zout/zsto)
1974     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), airephy, zx_tmp_2d)
1975 guez 15 CALL histwrite(nid_ins, "aire", itau_w, zx_tmp_2d)
1976 guez 3
1977     DO i = 1, klon
1978     zx_tmp_fi2d(i) = paprs(i, 1)
1979     ENDDO
1980     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
1981 guez 15 CALL histwrite(nid_ins, "psol", itau_w, zx_tmp_2d)
1982 guez 3
1983     DO i = 1, klon
1984     zx_tmp_fi2d(i) = rain_fall(i) + snow_fall(i)
1985     ENDDO
1986     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
1987 guez 15 CALL histwrite(nid_ins, "precip", itau_w, zx_tmp_2d)
1988 guez 3
1989     DO i = 1, klon
1990     zx_tmp_fi2d(i) = rain_lsc(i) + snow_lsc(i)
1991     ENDDO
1992     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
1993 guez 15 CALL histwrite(nid_ins, "plul", itau_w, zx_tmp_2d)
1994 guez 3
1995     DO i = 1, klon
1996     zx_tmp_fi2d(i) = rain_con(i) + snow_con(i)
1997     ENDDO
1998     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
1999 guez 15 CALL histwrite(nid_ins, "pluc", itau_w, zx_tmp_2d)
2000 guez 3
2001     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zxtsol, zx_tmp_2d)
2002 guez 15 CALL histwrite(nid_ins, "tsol", itau_w, zx_tmp_2d)
2003 guez 3 !ccIM
2004     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zt2m, zx_tmp_2d)
2005 guez 15 CALL histwrite(nid_ins, "t2m", itau_w, zx_tmp_2d)
2006 guez 3
2007     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zq2m, zx_tmp_2d)
2008 guez 15 CALL histwrite(nid_ins, "q2m", itau_w, zx_tmp_2d)
2009 guez 3
2010     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zu10m, zx_tmp_2d)
2011 guez 15 CALL histwrite(nid_ins, "u10m", itau_w, zx_tmp_2d)
2012 guez 3
2013     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zv10m, zx_tmp_2d)
2014 guez 15 CALL histwrite(nid_ins, "v10m", itau_w, zx_tmp_2d)
2015 guez 3
2016     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), snow_fall, zx_tmp_2d)
2017 guez 15 CALL histwrite(nid_ins, "snow", itau_w, zx_tmp_2d)
2018 guez 3
2019     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), cdragm, zx_tmp_2d)
2020 guez 15 CALL histwrite(nid_ins, "cdrm", itau_w, zx_tmp_2d)
2021 guez 3
2022     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), cdragh, zx_tmp_2d)
2023 guez 15 CALL histwrite(nid_ins, "cdrh", itau_w, zx_tmp_2d)
2024 guez 3
2025     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), toplw, zx_tmp_2d)
2026 guez 15 CALL histwrite(nid_ins, "topl", itau_w, zx_tmp_2d)
2027 guez 3
2028     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), evap, zx_tmp_2d)
2029 guez 15 CALL histwrite(nid_ins, "evap", itau_w, zx_tmp_2d)
2030 guez 3
2031     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), solsw, zx_tmp_2d)
2032 guez 15 CALL histwrite(nid_ins, "sols", itau_w, zx_tmp_2d)
2033 guez 3
2034     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), sollw, zx_tmp_2d)
2035 guez 15 CALL histwrite(nid_ins, "soll", itau_w, zx_tmp_2d)
2036 guez 3
2037     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), sollwdown, zx_tmp_2d)
2038 guez 15 CALL histwrite(nid_ins, "solldown", itau_w, zx_tmp_2d)
2039 guez 3
2040     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), bils, zx_tmp_2d)
2041 guez 15 CALL histwrite(nid_ins, "bils", itau_w, zx_tmp_2d)
2042 guez 3
2043     zx_tmp_fi2d(1:klon)=-1*sens(1:klon)
2044     ! CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), sens, zx_tmp_2d)
2045     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
2046 guez 15 CALL histwrite(nid_ins, "sens", itau_w, zx_tmp_2d)
2047 guez 3
2048     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), fder, zx_tmp_2d)
2049 guez 15 CALL histwrite(nid_ins, "fder", itau_w, zx_tmp_2d)
2050 guez 3
2051     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), d_ts(1, is_oce), zx_tmp_2d)
2052 guez 15 CALL histwrite(nid_ins, "dtsvdfo", itau_w, zx_tmp_2d)
2053 guez 3
2054     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), d_ts(1, is_ter), zx_tmp_2d)
2055 guez 15 CALL histwrite(nid_ins, "dtsvdft", itau_w, zx_tmp_2d)
2056 guez 3
2057     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), d_ts(1, is_lic), zx_tmp_2d)
2058 guez 15 CALL histwrite(nid_ins, "dtsvdfg", itau_w, zx_tmp_2d)
2059 guez 3
2060     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), d_ts(1, is_sic), zx_tmp_2d)
2061 guez 15 CALL histwrite(nid_ins, "dtsvdfi", itau_w, zx_tmp_2d)
2062 guez 3
2063     DO nsrf = 1, nbsrf
2064     !XXX
2065     zx_tmp_fi2d(1 : klon) = pctsrf( 1 : klon, nsrf)*100.
2066     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
2067     CALL histwrite(nid_ins, "pourc_"//clnsurf(nsrf), itau_w, &
2068 guez 15 zx_tmp_2d)
2069 guez 3
2070     zx_tmp_fi2d(1 : klon) = pctsrf( 1 : klon, nsrf)
2071     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
2072     CALL histwrite(nid_ins, "fract_"//clnsurf(nsrf), itau_w, &
2073 guez 15 zx_tmp_2d)
2074 guez 3
2075     zx_tmp_fi2d(1 : klon) = fluxt( 1 : klon, 1, nsrf)
2076     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
2077     CALL histwrite(nid_ins, "sens_"//clnsurf(nsrf), itau_w, &
2078 guez 15 zx_tmp_2d)
2079 guez 3
2080     zx_tmp_fi2d(1 : klon) = fluxlat( 1 : klon, nsrf)
2081     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
2082     CALL histwrite(nid_ins, "lat_"//clnsurf(nsrf), itau_w, &
2083 guez 15 zx_tmp_2d)
2084 guez 3
2085     zx_tmp_fi2d(1 : klon) = ftsol( 1 : klon, nsrf)
2086     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
2087     CALL histwrite(nid_ins, "tsol_"//clnsurf(nsrf), itau_w, &
2088 guez 15 zx_tmp_2d)
2089 guez 3
2090     zx_tmp_fi2d(1 : klon) = fluxu( 1 : klon, 1, nsrf)
2091     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
2092     CALL histwrite(nid_ins, "taux_"//clnsurf(nsrf), itau_w, &
2093 guez 15 zx_tmp_2d)
2094 guez 3
2095     zx_tmp_fi2d(1 : klon) = fluxv( 1 : klon, 1, nsrf)
2096     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
2097     CALL histwrite(nid_ins, "tauy_"//clnsurf(nsrf), itau_w, &
2098 guez 15 zx_tmp_2d)
2099 guez 3
2100     zx_tmp_fi2d(1 : klon) = frugs( 1 : klon, nsrf)
2101     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
2102     CALL histwrite(nid_ins, "rugs_"//clnsurf(nsrf), itau_w, &
2103 guez 15 zx_tmp_2d)
2104 guez 3
2105     zx_tmp_fi2d(1 : klon) = falbe( 1 : klon, nsrf)
2106     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zx_tmp_fi2d, zx_tmp_2d)
2107     CALL histwrite(nid_ins, "albe_"//clnsurf(nsrf), itau_w, &
2108 guez 15 zx_tmp_2d)
2109 guez 3
2110     END DO
2111     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), albsol, zx_tmp_2d)
2112 guez 15 CALL histwrite(nid_ins, "albs", itau_w, zx_tmp_2d)
2113 guez 3 CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), albsollw, zx_tmp_2d)
2114 guez 15 CALL histwrite(nid_ins, "albslw", itau_w, zx_tmp_2d)
2115 guez 3
2116     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zxrugs, zx_tmp_2d)
2117 guez 15 CALL histwrite(nid_ins, "rugs", itau_w, zx_tmp_2d)
2118 guez 3
2119     !IM cf. AM 081204 BEG
2120    
2121     !HBTM2
2122    
2123     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_pblh, zx_tmp_2d)
2124 guez 15 CALL histwrite(nid_ins, "s_pblh", itau_w, zx_tmp_2d)
2125 guez 3
2126     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_pblt, zx_tmp_2d)
2127 guez 15 CALL histwrite(nid_ins, "s_pblt", itau_w, zx_tmp_2d)
2128 guez 3
2129     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_lcl, zx_tmp_2d)
2130 guez 15 CALL histwrite(nid_ins, "s_lcl", itau_w, zx_tmp_2d)
2131 guez 3
2132     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_capCL, zx_tmp_2d)
2133 guez 15 CALL histwrite(nid_ins, "s_capCL", itau_w, zx_tmp_2d)
2134 guez 3
2135     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_oliqCL, zx_tmp_2d)
2136 guez 15 CALL histwrite(nid_ins, "s_oliqCL", itau_w, zx_tmp_2d)
2137 guez 3
2138     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_cteiCL, zx_tmp_2d)
2139 guez 15 CALL histwrite(nid_ins, "s_cteiCL", itau_w, zx_tmp_2d)
2140 guez 3
2141     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_therm, zx_tmp_2d)
2142 guez 15 CALL histwrite(nid_ins, "s_therm", itau_w, zx_tmp_2d)
2143 guez 3
2144     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_trmb1, zx_tmp_2d)
2145 guez 15 CALL histwrite(nid_ins, "s_trmb1", itau_w, zx_tmp_2d)
2146 guez 3
2147     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_trmb2, zx_tmp_2d)
2148 guez 15 CALL histwrite(nid_ins, "s_trmb2", itau_w, zx_tmp_2d)
2149 guez 3
2150     CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_trmb3, zx_tmp_2d)
2151 guez 15 CALL histwrite(nid_ins, "s_trmb3", itau_w, zx_tmp_2d)
2152 guez 3
2153     !IM cf. AM 081204 END
2154    
2155     ! Champs 3D:
2156    
2157     CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), t_seri, zx_tmp_3d)
2158 guez 15 CALL histwrite(nid_ins, "temp", itau_w, zx_tmp_3d)
2159 guez 3
2160     CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), u_seri, zx_tmp_3d)
2161 guez 15 CALL histwrite(nid_ins, "vitu", itau_w, zx_tmp_3d)
2162 guez 3
2163     CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), v_seri, zx_tmp_3d)
2164 guez 15 CALL histwrite(nid_ins, "vitv", itau_w, zx_tmp_3d)
2165 guez 3
2166     CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), zphi, zx_tmp_3d)
2167 guez 15 CALL histwrite(nid_ins, "geop", itau_w, zx_tmp_3d)
2168 guez 3
2169     CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), pplay, zx_tmp_3d)
2170 guez 15 CALL histwrite(nid_ins, "pres", itau_w, zx_tmp_3d)
2171 guez 3
2172     CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), d_t_vdf, zx_tmp_3d)
2173 guez 15 CALL histwrite(nid_ins, "dtvdf", itau_w, zx_tmp_3d)
2174 guez 3
2175     CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), d_q_vdf, zx_tmp_3d)
2176 guez 15 CALL histwrite(nid_ins, "dqvdf", itau_w, zx_tmp_3d)
2177 guez 3
2178     if (ok_sync) then
2179     call histsync(nid_ins)
2180     endif
2181     ENDIF
2182    
2183     end subroutine write_histins
2184    
2185     !****************************************************
2186    
2187     subroutine write_histhf3d
2188    
2189     ! From phylmd/write_histhf3d.h, v 1.2 2005/05/25 13:10:09
2190    
2191 guez 17 integer itau_w ! pas de temps ecriture
2192 guez 3
2193 guez 17 !-------------------------------------------------------
2194    
2195 guez 3 itau_w = itau_phy + itap
2196    
2197     ! Champs 3D:
2198    
2199     CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), t_seri, zx_tmp_3d)
2200 guez 15 CALL histwrite(nid_hf3d, "temp", itau_w, zx_tmp_3d)
2201 guez 3
2202     CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), qx(1, 1, ivap), zx_tmp_3d)
2203 guez 15 CALL histwrite(nid_hf3d, "ovap", itau_w, zx_tmp_3d)
2204 guez 3
2205     CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), u_seri, zx_tmp_3d)
2206 guez 15 CALL histwrite(nid_hf3d, "vitu", itau_w, zx_tmp_3d)
2207 guez 3
2208     CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), v_seri, zx_tmp_3d)
2209 guez 15 CALL histwrite(nid_hf3d, "vitv", itau_w, zx_tmp_3d)
2210 guez 3
2211 guez 6 if (nbtr >= 3) then
2212     CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), tr_seri(1, 1, 3), &
2213     zx_tmp_3d)
2214 guez 15 CALL histwrite(nid_hf3d, "O3", itau_w, zx_tmp_3d)
2215 guez 6 end if
2216 guez 3
2217     if (ok_sync) then
2218     call histsync(nid_hf3d)
2219     endif
2220    
2221     end subroutine write_histhf3d
2222    
2223     END SUBROUTINE physiq
2224    
2225     end module physiq_m

  ViewVC Help
Powered by ViewVC 1.1.21