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

Annotation of /trunk/phylmd/physiq.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 327 - (hide annotations)
Thu Jun 13 13:59:19 2019 UTC (4 years, 11 months ago) by guez
File size: 36331 byte(s)
Move the update of ftsol by `d_ts` inside `pbl_surface` (following
LMDZ). This makes the procedure physiq lighter. And it is clearer now
that the update of ftsol does come from `pbl_surface`.

1 guez 3 module physiq_m
2    
3     IMPLICIT none
4    
5     contains
6    
7 guez 154 SUBROUTINE physiq(lafin, dayvrai, time, paprs, play, pphi, pphis, u, v, t, &
8     qx, omega, d_u, d_v, d_t, d_qx)
9 guez 3
10 guez 72 ! From phylmd/physiq.F, version 1.22 2006/02/20 09:38:28
11     ! (subversion revision 678)
12    
13 guez 154 ! Author: Z. X. Li (LMD/CNRS) 1993
14 guez 3
15 guez 49 ! This is the main procedure for the "physics" part of the program.
16 guez 3
17 guez 56 use aaam_bud_m, only: aaam_bud
18 guez 51 USE abort_gcm_m, ONLY: abort_gcm
19 guez 53 use ajsec_m, only: ajsec
20 guez 52 use calltherm_m, only: calltherm
21 guez 250 USE clesphys, ONLY: cdhmax, cdmmax, ecrit_ins, ok_instan
22 guez 209 USE clesphys2, ONLY: conv_emanuel, nbapp_rad, new_oliq, ok_orodr, ok_orolf
23 guez 301 USE conf_interface_m, ONLY: conf_interface
24 guez 267 USE pbl_surface_m, ONLY: pbl_surface
25 guez 72 use clouds_gno_m, only: clouds_gno
26 guez 154 use comconst, only: dtphys
27 guez 137 USE comgeomphy, ONLY: airephy
28 guez 51 USE concvl_m, ONLY: concvl
29 guez 224 USE conf_gcm_m, ONLY: lmt_pas
30 guez 51 USE conf_phys_m, ONLY: conf_phys
31 guez 62 use conflx_m, only: conflx
32 guez 323 USE ctherm_m, ONLY: iflag_thermals, ctherm
33 guez 52 use diagcld2_m, only: diagcld2
34 guez 265 USE dimensions, ONLY: llm, nqmx
35 guez 98 USE dimphy, ONLY: klon
36 guez 51 USE dimsoil, ONLY: nsoilmx
37 guez 52 use drag_noro_m, only: drag_noro
38 guez 313 use dynetat0_chosen_m, only: day_ref, annee_ref
39 guez 221 USE fcttre, ONLY: foeew
40 guez 68 use fisrtilp_m, only: fisrtilp
41 guez 51 USE hgardfou_m, ONLY: hgardfou
42 guez 191 USE histsync_m, ONLY: histsync
43     USE histwrite_phy_m, ONLY: histwrite_phy
44 guez 327 USE indicesol, ONLY: clnsurf, epsfra, nbsrf
45 guez 191 USE ini_histins_m, ONLY: ini_histins, nid_ins
46 guez 227 use lift_noro_m, only: lift_noro
47 guez 157 use netcdf95, only: NF95_CLOSE
48 guez 68 use newmicro_m, only: newmicro
49 guez 191 use nr_util, only: assert
50 guez 175 use nuage_m, only: nuage
51 guez 98 USE orbite_m, ONLY: orbite
52 guez 51 USE ozonecm_m, ONLY: ozonecm
53 guez 227 USE phyetat0_m, ONLY: phyetat0
54 guez 51 USE phyredem_m, ONLY: phyredem
55 guez 157 USE phyredem0_m, ONLY: phyredem0
56 guez 51 USE phytrac_m, ONLY: phytrac
57 guez 53 use radlwsw_m, only: radlwsw
58 guez 158 use yoegwd, only: sugwd
59 guez 212 USE suphec_m, ONLY: rcpd, retv, rg, rlvtt, romega, rsigma, rtt, rmo3, md
60 guez 191 use time_phylmdz, only: itap, increment_itap
61 guez 169 use transp_m, only: transp
62 guez 178 use transp_lay_m, only: transp_lay
63 guez 68 use unit_nml_m, only: unit_nml
64 guez 92 USE ymds2ju_m, ONLY: ymds2ju
65 guez 51 USE yoethf_m, ONLY: r2es, rvtmp2
66 guez 98 use zenang_m, only: zenang
67 guez 3
68 guez 91 logical, intent(in):: lafin ! dernier passage
69 guez 3
70 guez 130 integer, intent(in):: dayvrai
71     ! current day number, based at value 1 on January 1st of annee_ref
72 guez 15
73 guez 90 REAL, intent(in):: time ! heure de la journ\'ee en fraction de jour
74 guez 3
75 guez 98 REAL, intent(in):: paprs(:, :) ! (klon, llm + 1)
76     ! pression pour chaque inter-couche, en Pa
77 guez 15
78 guez 98 REAL, intent(in):: play(:, :) ! (klon, llm)
79     ! pression pour le mileu de chaque couche (en Pa)
80 guez 3
81 guez 190 REAL, intent(in):: pphi(:, :) ! (klon, llm)
82 guez 97 ! géopotentiel de chaque couche (référence sol)
83 guez 3
84 guez 98 REAL, intent(in):: pphis(:) ! (klon) géopotentiel du sol
85 guez 3
86 guez 98 REAL, intent(in):: u(:, :) ! (klon, llm)
87 guez 202 ! vitesse dans la direction X (de O a E) en m / s
88 guez 51
89 guez 202 REAL, intent(in):: v(:, :) ! (klon, llm) vitesse Y (de S a N) en m / s
90 guez 98 REAL, intent(in):: t(:, :) ! (klon, llm) temperature (K)
91 guez 3
92 guez 98 REAL, intent(in):: qx(:, :, :) ! (klon, llm, nqmx)
93 guez 90 ! (humidit\'e sp\'ecifique et fractions massiques des autres traceurs)
94 guez 3
95 guez 202 REAL, intent(in):: omega(:, :) ! (klon, llm) vitesse verticale en Pa / s
96 guez 98 REAL, intent(out):: d_u(:, :) ! (klon, llm) tendance physique de "u" (m s-2)
97     REAL, intent(out):: d_v(:, :) ! (klon, llm) tendance physique de "v" (m s-2)
98 guez 202 REAL, intent(out):: d_t(:, :) ! (klon, llm) tendance physique de "t" (K / s)
99 guez 3
100 guez 98 REAL, intent(out):: d_qx(:, :, :) ! (klon, llm, nqmx)
101     ! tendance physique de "qx" (s-1)
102    
103 guez 91 ! Local:
104    
105 guez 35 LOGICAL:: firstcal = .true.
106    
107 guez 51 LOGICAL, PARAMETER:: ok_stratus = .FALSE.
108 guez 49 ! Ajouter artificiellement les stratus
109    
110 guez 200 ! pour phystoke avec thermiques
111 guez 51 REAL fm_therm(klon, llm + 1)
112 guez 3 REAL entr_therm(klon, llm)
113 guez 51 real, save:: q2(klon, llm + 1, nbsrf)
114 guez 3
115 guez 98 INTEGER, PARAMETER:: ivap = 1 ! indice de traceur pour vapeur d'eau
116     INTEGER, PARAMETER:: iliq = 2 ! indice de traceur pour eau liquide
117 guez 3
118 guez 49 REAL, save:: t_ancien(klon, llm), q_ancien(klon, llm)
119     LOGICAL, save:: ancien_ok
120 guez 3
121 guez 202 REAL d_t_dyn(klon, llm) ! tendance dynamique pour "t" (K / s)
122     REAL d_q_dyn(klon, llm) ! tendance dynamique pour "q" (kg / kg / s)
123 guez 3
124     real da(klon, llm), phi(klon, llm, llm), mp(klon, llm)
125    
126 guez 205 REAL, save:: swdn0(klon, llm + 1), swdn(klon, llm + 1)
127     REAL, save:: swup0(klon, llm + 1), swup(klon, llm + 1)
128 guez 3
129 guez 205 REAL, save:: lwdn0(klon, llm + 1), lwdn(klon, llm + 1)
130     REAL, save:: lwup0(klon, llm + 1), lwup(klon, llm + 1)
131 guez 3
132     ! prw: precipitable water
133     real prw(klon)
134    
135 guez 202 ! flwp, fiwp = Liquid Water Path & Ice Water Path (kg / m2)
136     ! flwc, fiwc = Liquid Water Content & Ice Water Content (kg / kg)
137 guez 3 REAL flwp(klon), fiwp(klon)
138     REAL flwc(klon, llm), fiwc(klon, llm)
139    
140     ! Variables propres a la physique
141    
142     INTEGER, save:: radpas
143 guez 125 ! Radiative transfer computations are made every "radpas" call to
144     ! "physiq".
145 guez 3
146 guez 308 REAL, save:: radsol(klon)
147 guez 324 ! Bilan radiatif net au sol (W/m2), positif vers le bas. Must be
148     ! saved because radlwsw is not called at every time step.
149 guez 308
150 guez 310 REAL, save:: ftsol(klon, nbsrf) ! skin temperature of surface fraction, in K
151 guez 3
152 guez 49 REAL, save:: ftsoil(klon, nsoilmx, nbsrf)
153     ! soil temperature of surface fraction
154 guez 3
155 guez 311 REAL fluxlat(klon, nbsrf) ! flux de chaleur latente, en W m-2
156 guez 3
157 guez 98 REAL, save:: fqsurf(klon, nbsrf)
158     ! humidite de l'air au contact de la surface
159 guez 3
160 guez 215 REAL, save:: qsol(klon) ! column-density of water in soil, in kg m-2
161     REAL, save:: fsnow(klon, nbsrf) ! \'epaisseur neigeuse
162 guez 155 REAL, save:: falbe(klon, nbsrf) ! albedo visible par type de surface
163 guez 3
164 guez 90 ! Param\`etres de l'orographie \`a l'\'echelle sous-maille (OESM) :
165 guez 13 REAL, save:: zmea(klon) ! orographie moyenne
166     REAL, save:: zstd(klon) ! deviation standard de l'OESM
167     REAL, save:: zsig(klon) ! pente de l'OESM
168     REAL, save:: zgam(klon) ! anisotropie de l'OESM
169     REAL, save:: zthe(klon) ! orientation de l'OESM
170     REAL, save:: zpic(klon) ! Maximum de l'OESM
171     REAL, save:: zval(klon) ! Minimum de l'OESM
172     REAL, save:: rugoro(klon) ! longueur de rugosite de l'OESM
173 guez 3 REAL zulow(klon), zvlow(klon)
174 guez 247 INTEGER ktest(klon)
175 guez 3
176 guez 189 REAL, save:: agesno(klon, nbsrf) ! age de la neige
177     REAL, save:: run_off_lic_0(klon)
178 guez 3
179 guez 189 ! Variables li\'ees \`a la convection d'Emanuel :
180     REAL, save:: Ma(klon, llm) ! undilute upward mass flux
181 guez 72 REAL, save:: sig1(klon, llm), w01(klon, llm)
182 guez 3
183 guez 189 ! Variables pour la couche limite (Alain Lahellec) :
184 guez 3 REAL cdragh(klon) ! drag coefficient pour T and Q
185     REAL cdragm(klon) ! drag coefficient pour vent
186    
187 guez 244 REAL coefh(klon, 2:llm) ! coef d'echange pour phytrac
188 guez 191
189 guez 205 REAL, save:: ffonte(klon, nbsrf)
190     ! flux thermique utilise pour fondre la neige
191    
192 guez 279 REAL fqcalving(klon, nbsrf)
193     ! flux d'eau "perdue" par la surface et n\'ecessaire pour limiter
194     ! la hauteur de neige, en kg / m2 / s
195 guez 191
196 guez 279 REAL zxffonte(klon)
197 guez 3
198 guez 205 REAL, save:: pfrac_impa(klon, llm)! Produits des coefs lessivage impaction
199     REAL, save:: pfrac_nucl(klon, llm)! Produits des coefs lessivage nucleation
200    
201     REAL, save:: pfrac_1nucl(klon, llm)
202     ! Produits des coefs lessi nucl (alpha = 1)
203    
204 guez 213 REAL frac_impa(klon, llm) ! fraction d'a\'erosols lessiv\'es (impaction)
205 guez 3 REAL frac_nucl(klon, llm) ! idem (nucleation)
206    
207 guez 101 REAL, save:: rain_fall(klon)
208 guez 202 ! liquid water mass flux (kg / m2 / s), positive down
209 guez 62
210 guez 101 REAL, save:: snow_fall(klon)
211 guez 202 ! solid water mass flux (kg / m2 / s), positive down
212 guez 101
213 guez 3 REAL rain_tiedtke(klon), snow_tiedtke(klon)
214    
215 guez 206 REAL evap(klon) ! flux d'\'evaporation au sol
216 guez 302 real dflux_q(klon) ! derivative of the evaporation flux at the surface
217 guez 206 REAL sens(klon) ! flux de chaleur sensible au sol
218 guez 302 real dflux_t(klon) ! derivee du flux de chaleur sensible au sol
219 guez 223 REAL, save:: dlw(klon) ! derivative of infra-red flux
220 guez 3 REAL bils(klon) ! bilan de chaleur au sol
221 guez 325 REAL fder(klon) ! d\'erive de flux (sensible et latente)
222 guez 3 REAL ve(klon) ! integr. verticale du transport meri. de l'energie
223     REAL vq(klon) ! integr. verticale du transport meri. de l'eau
224     REAL ue(klon) ! integr. verticale du transport zonal de l'energie
225     REAL uq(klon) ! integr. verticale du transport zonal de l'eau
226    
227 guez 98 REAL, save:: frugs(klon, nbsrf) ! longueur de rugosite
228 guez 3 REAL zxrugs(klon) ! longueur de rugosite
229    
230     ! Conditions aux limites
231    
232     INTEGER julien
233 guez 70 REAL, save:: pctsrf(klon, nbsrf) ! percentage of surface
234 guez 217 REAL, save:: albsol(klon) ! albedo du sol total, visible, moyen par maille
235 guez 17 REAL, SAVE:: wo(klon, llm) ! column density of ozone in a cell, in kDU
236 guez 212 real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2
237 guez 3
238 guez 72 real, save:: clwcon(klon, llm), rnebcon(klon, llm)
239     real, save:: clwcon0(klon, llm), rnebcon0(klon, llm)
240 guez 3
241 guez 266 REAL rhcl(klon, llm) ! humidit\'e relative ciel clair
242 guez 47 REAL dialiq(klon, llm) ! eau liquide nuageuse
243     REAL diafra(klon, llm) ! fraction nuageuse
244     REAL cldliq(klon, llm) ! eau liquide nuageuse
245     REAL cldfra(klon, llm) ! fraction nuageuse
246     REAL cldtau(klon, llm) ! epaisseur optique
247     REAL cldemi(klon, llm) ! emissivite infrarouge
248 guez 3
249 guez 206 REAL flux_q(klon, nbsrf) ! flux turbulent d'humidite à la surface
250 guez 3
251 guez 309 REAL flux_t(klon, nbsrf)
252     ! flux de chaleur sensible (c_p T) (W / m2) (orientation positive
253     ! vers le bas) à la surface
254    
255 guez 229 REAL flux_u(klon, nbsrf), flux_v(klon, nbsrf)
256     ! tension du vent (flux turbulent de vent) à la surface, en Pa
257    
258 guez 90 ! Le rayonnement n'est pas calcul\'e tous les pas, il faut donc que
259     ! les variables soient r\'emanentes.
260 guez 53 REAL, save:: heat(klon, llm) ! chauffage solaire
261 guez 154 REAL, save:: heat0(klon, llm) ! chauffage solaire ciel clair
262 guez 62 REAL, save:: cool(klon, llm) ! refroidissement infrarouge
263 guez 154 REAL, save:: cool0(klon, llm) ! refroidissement infrarouge ciel clair
264 guez 72 REAL, save:: topsw(klon), toplw(klon), solsw(klon)
265 guez 308
266     REAL, save:: sollw(klon) ! surface net downward longwave flux, in W m-2
267 guez 325 real, save:: sollwdown(klon) ! downwelling longwave flux at surface
268 guez 62 REAL, save:: topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)
269 guez 154 REAL, save:: albpla(klon)
270 guez 3
271 guez 202 REAL conv_q(klon, llm) ! convergence de l'humidite (kg / kg / s)
272     REAL conv_t(klon, llm) ! convergence of temperature (K / s)
273 guez 3
274 guez 191 REAL cldl(klon), cldm(klon), cldh(klon) ! nuages bas, moyen et haut
275     REAL cldt(klon), cldq(klon) ! nuage total, eau liquide integree
276 guez 3
277 guez 221 REAL zxfluxlat(klon)
278 guez 118 REAL dist, mu0(klon), fract(klon)
279     real longi
280 guez 3 REAL z_avant(klon), z_apres(klon), z_factor(klon)
281 guez 205 REAL zb
282 guez 323 REAL zx_qs, zcor
283 guez 3 real zqsat(klon, llm)
284     INTEGER i, k, iq, nsrf
285     REAL zphi(klon, llm)
286    
287 guez 200 ! cf. Anne Mathieu, variables pour la couche limite atmosphérique (hbtm)
288 guez 3
289 guez 49 REAL, SAVE:: pblh(klon, nbsrf) ! Hauteur de couche limite
290     REAL, SAVE:: plcl(klon, nbsrf) ! Niveau de condensation de la CLA
291     REAL, SAVE:: capCL(klon, nbsrf) ! CAPE de couche limite
292     REAL, SAVE:: oliqCL(klon, nbsrf) ! eau_liqu integree de couche limite
293     REAL, SAVE:: cteiCL(klon, nbsrf) ! cloud top instab. crit. couche limite
294 guez 207 REAL, SAVE:: pblt(klon, nbsrf) ! T \`a la hauteur de couche limite
295 guez 49 REAL, SAVE:: therm(klon, nbsrf)
296 guez 186 ! Grandeurs de sorties
297 guez 3 REAL s_pblh(klon), s_lcl(klon), s_capCL(klon)
298     REAL s_oliqCL(klon), s_cteiCL(klon), s_pblt(klon)
299 guez 252 REAL s_therm(klon)
300 guez 3
301 guez 175 ! Variables pour la convection de K. Emanuel :
302 guez 3
303 guez 47 REAL upwd(klon, llm) ! saturated updraft mass flux
304     REAL dnwd(klon, llm) ! saturated downdraft mass flux
305 guez 205 REAL, save:: cape(klon)
306 guez 3
307 guez 47 INTEGER iflagctrl(klon) ! flag fonctionnement de convect
308 guez 3
309     ! Variables du changement
310    
311     ! con: convection
312 guez 51 ! lsc: large scale condensation
313 guez 3 ! ajs: ajustement sec
314 guez 90 ! eva: \'evaporation de l'eau liquide nuageuse
315 guez 51 ! vdf: vertical diffusion in boundary layer
316 guez 3 REAL d_t_con(klon, llm), d_q_con(klon, llm)
317 guez 205 REAL, save:: d_u_con(klon, llm), d_v_con(klon, llm)
318 guez 3 REAL d_t_lsc(klon, llm), d_q_lsc(klon, llm), d_ql_lsc(klon, llm)
319     REAL d_t_ajs(klon, llm), d_q_ajs(klon, llm)
320     REAL d_u_ajs(klon, llm), d_v_ajs(klon, llm)
321     REAL rneb(klon, llm)
322    
323 guez 71 REAL mfu(klon, llm), mfd(klon, llm)
324 guez 3 REAL pen_u(klon, llm), pen_d(klon, llm)
325     REAL pde_u(klon, llm), pde_d(klon, llm)
326     INTEGER kcbot(klon), kctop(klon), kdtop(klon)
327 guez 51 REAL pmflxr(klon, llm + 1), pmflxs(klon, llm + 1)
328     REAL prfl(klon, llm + 1), psfl(klon, llm + 1)
329 guez 3
330 guez 62 INTEGER, save:: ibas_con(klon), itop_con(klon)
331 guez 183 real ema_pct(klon) ! Emanuel pressure at cloud top, in Pa
332 guez 3
333 guez 305 REAL rain_con(klon)
334 guez 205 real rain_lsc(klon)
335 guez 305 REAL snow_con(klon) ! neige (mm / s)
336 guez 180 real snow_lsc(klon)
337 guez 3
338     REAL d_u_vdf(klon, llm), d_v_vdf(klon, llm)
339     REAL d_t_vdf(klon, llm), d_q_vdf(klon, llm)
340    
341     REAL d_u_oro(klon, llm), d_v_oro(klon, llm)
342     REAL d_t_oro(klon, llm)
343     REAL d_u_lif(klon, llm), d_v_lif(klon, llm)
344     REAL d_t_lif(klon, llm)
345    
346 guez 68 REAL, save:: ratqs(klon, llm)
347     real ratqss(klon, llm), ratqsc(klon, llm)
348     real:: ratqsbas = 0.01, ratqshaut = 0.3
349 guez 3
350     ! Parametres lies au nouveau schema de nuages (SB, PDF)
351 guez 68 real:: fact_cldcon = 0.375
352     real:: facttemps = 1.e-4
353     logical:: ok_newmicro = .true.
354 guez 3 real facteur
355    
356 guez 68 integer:: iflag_cldcon = 1
357 guez 3 logical ptconv(klon, llm)
358    
359 guez 175 ! Variables pour effectuer les appels en s\'erie :
360 guez 3
361     REAL t_seri(klon, llm), q_seri(klon, llm)
362 guez 98 REAL ql_seri(klon, llm)
363 guez 3 REAL u_seri(klon, llm), v_seri(klon, llm)
364 guez 98 REAL tr_seri(klon, llm, nqmx - 2)
365 guez 3
366     REAL zx_rh(klon, llm)
367    
368     REAL zustrdr(klon), zvstrdr(klon)
369     REAL zustrli(klon), zvstrli(klon)
370     REAL aam, torsfc
371    
372     REAL ve_lay(klon, llm) ! transport meri. de l'energie a chaque niveau vert.
373     REAL vq_lay(klon, llm) ! transport meri. de l'eau a chaque niveau vert.
374     REAL ue_lay(klon, llm) ! transport zonal de l'energie a chaque niveau vert.
375     REAL uq_lay(klon, llm) ! transport zonal de l'eau a chaque niveau vert.
376    
377 guez 221 REAL tsol(klon)
378 guez 51
379 guez 202 REAL d_t_ec(klon, llm)
380 guez 213 ! tendance due \`a la conversion d'\'energie cin\'etique en
381     ! énergie thermique
382 guez 200
383 guez 205 REAL, save:: t2m(klon, nbsrf), q2m(klon, nbsrf)
384     ! temperature and humidity at 2 m
385    
386 guez 225 REAL, save:: u10m_srf(klon, nbsrf), v10m_srf(klon, nbsrf)
387     ! composantes du vent \`a 10 m
388    
389 guez 207 REAL zt2m(klon), zq2m(klon) ! température, humidité 2 m moyenne sur 1 maille
390 guez 225 REAL u10m(klon), v10m(klon) ! vent \`a 10 m moyenn\' sur les sous-surfaces
391 guez 3
392 guez 69 ! Aerosol effects:
393    
394 guez 205 REAL, save:: topswad(klon), solswad(klon) ! aerosol direct effect
395 guez 68 LOGICAL:: ok_ade = .false. ! apply aerosol direct effect
396 guez 3
397 guez 68 REAL:: bl95_b0 = 2., bl95_b1 = 0.2
398 guez 69 ! Parameters in equation (D) of Boucher and Lohmann (1995, Tellus
399     ! B). They link cloud droplet number concentration to aerosol mass
400     ! concentration.
401 guez 68
402 guez 190 real zmasse(klon, llm)
403 guez 17 ! (column-density of mass of air in a cell, in kg m-2)
404    
405 guez 191 integer, save:: ncid_startphy
406 guez 17
407 guez 204 namelist /physiq_nml/ fact_cldcon, facttemps, ok_newmicro, iflag_cldcon, &
408 guez 323 ratqsbas, ratqshaut, ok_ade, bl95_b0, bl95_b1
409 guez 68
410 guez 3 !----------------------------------------------------------------
411    
412 guez 69 IF (nqmx < 2) CALL abort_gcm('physiq', &
413 guez 171 'eaux vapeur et liquide sont indispensables')
414 guez 3
415 guez 7 test_firstcal: IF (firstcal) THEN
416 guez 47 ! initialiser
417 guez 225 u10m_srf = 0.
418     v10m_srf = 0.
419 guez 51 t2m = 0.
420     q2m = 0.
421     ffonte = 0.
422 guez 72 d_u_con = 0.
423     d_v_con = 0.
424     rnebcon0 = 0.
425     clwcon0 = 0.
426     rnebcon = 0.
427     clwcon = 0.
428 guez 47 pblh =0. ! Hauteur de couche limite
429     plcl =0. ! Niveau de condensation de la CLA
430     capCL =0. ! CAPE de couche limite
431     oliqCL =0. ! eau_liqu integree de couche limite
432     cteiCL =0. ! cloud top instab. crit. couche limite
433 guez 207 pblt =0.
434 guez 47 therm =0.
435 guez 3
436 guez 68 print *, "Enter namelist 'physiq_nml'."
437     read(unit=*, nml=physiq_nml)
438     write(unit_nml, nml=physiq_nml)
439    
440 guez 323 call ctherm
441 guez 68 call conf_phys
442 guez 3
443     ! Initialiser les compteurs:
444    
445     frugs = 0.
446 guez 191 CALL phyetat0(pctsrf, ftsol, ftsoil, fqsurf, qsol, fsnow, falbe, &
447 guez 304 rain_fall, snow_fall, solsw, sollw, dlw, radsol, frugs, agesno, &
448     zmea, zstd, zsig, zgam, zthe, zpic, zval, t_ancien, q_ancien, &
449     ancien_ok, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, w01, &
450     ncid_startphy)
451 guez 3
452 guez 47 ! ATTENTION : il faudra a terme relire q2 dans l'etat initial
453 guez 69 q2 = 1e-8
454 guez 3
455 guez 154 radpas = lmt_pas / nbapp_rad
456 guez 191 print *, "radpas = ", radpas
457 guez 154
458 guez 90 ! Initialisation pour le sch\'ema de convection d'Emanuel :
459 guez 182 IF (conv_emanuel) THEN
460 guez 69 ibas_con = 1
461     itop_con = 1
462 guez 3 ENDIF
463    
464     IF (ok_orodr) THEN
465 guez 13 rugoro = MAX(1e-5, zstd * zsig / 2)
466 guez 54 CALL SUGWD(paprs, play)
467 guez 13 else
468     rugoro = 0.
469 guez 3 ENDIF
470    
471 guez 47 ! Initialisation des sorties
472 guez 298 call ini_histins(ok_newmicro)
473 guez 202 CALL phyredem0
474 guez 301 call conf_interface
475 guez 7 ENDIF test_firstcal
476 guez 3
477 guez 91 ! We will modify variables *_seri and we will not touch variables
478 guez 98 ! u, v, t, qx:
479     t_seri = t
480     u_seri = u
481     v_seri = v
482     q_seri = qx(:, :, ivap)
483     ql_seri = qx(:, :, iliq)
484 guez 157 tr_seri = qx(:, :, 3:nqmx)
485 guez 3
486 guez 221 tsol = sum(ftsol * pctsrf, dim = 2)
487 guez 3
488 guez 51 ! Diagnostic de la tendance dynamique :
489 guez 3 IF (ancien_ok) THEN
490     DO k = 1, llm
491     DO i = 1, klon
492 guez 49 d_t_dyn(i, k) = (t_seri(i, k) - t_ancien(i, k)) / dtphys
493     d_q_dyn(i, k) = (q_seri(i, k) - q_ancien(i, k)) / dtphys
494 guez 3 ENDDO
495     ENDDO
496     ELSE
497     DO k = 1, llm
498     DO i = 1, klon
499 guez 72 d_t_dyn(i, k) = 0.
500     d_q_dyn(i, k) = 0.
501 guez 3 ENDDO
502     ENDDO
503     ancien_ok = .TRUE.
504     ENDIF
505    
506     ! Ajouter le geopotentiel du sol:
507     DO k = 1, llm
508     DO i = 1, klon
509     zphi(i, k) = pphi(i, k) + pphis(i)
510     ENDDO
511     ENDDO
512    
513 guez 49 ! Check temperatures:
514 guez 3 CALL hgardfou(t_seri, ftsol)
515    
516 guez 191 call increment_itap
517 guez 130 julien = MOD(dayvrai, 360)
518 guez 3 if (julien == 0) julien = 360
519    
520 guez 103 forall (k = 1: llm) zmasse(:, k) = (paprs(:, k) - paprs(:, k + 1)) / rg
521 guez 17
522 guez 90 ! \'Evaporation de l'eau liquide nuageuse :
523 guez 51 DO k = 1, llm
524 guez 3 DO i = 1, klon
525 guez 51 zb = MAX(0., ql_seri(i, k))
526     t_seri(i, k) = t_seri(i, k) &
527     - zb * RLVTT / RCPD / (1. + RVTMP2 * q_seri(i, k))
528 guez 3 q_seri(i, k) = q_seri(i, k) + zb
529     ENDDO
530     ENDDO
531 guez 51 ql_seri = 0.
532 guez 3
533 guez 98 frugs = MAX(frugs, 0.000015)
534     zxrugs = sum(frugs * pctsrf, dim = 2)
535 guez 3
536 guez 191 ! Calculs n\'ecessaires au calcul de l'albedo dans l'interface avec
537 guez 118 ! la surface.
538 guez 3
539 guez 118 CALL orbite(REAL(julien), longi, dist)
540 guez 209 CALL zenang(longi, time, dtphys * radpas, mu0, fract)
541 guez 3
542 guez 298 CALL pbl_surface(pctsrf, t_seri, q_seri, u_seri, v_seri, julien, mu0, &
543     ftsol, cdmmax, cdhmax, ftsoil, qsol, paprs, play, fsnow, fqsurf, &
544 guez 307 falbe, fluxlat, rain_fall, snow_fall, frugs, agesno, rugoro, d_t_vdf, &
545 guez 327 d_q_vdf, d_u_vdf, d_v_vdf, flux_t, flux_q, flux_u, flux_v, cdragh, &
546     cdragm, q2, dflux_t, dflux_q, coefh, t2m, q2m, u10m_srf, v10m_srf, &
547     pblh, capCL, oliqCL, cteiCL, pblT, therm, plcl, fqcalving, ffonte, &
548     run_off_lic_0, albsol, sollw, solsw, tsol)
549 guez 3
550 guez 90 ! Incr\'ementation des flux
551 guez 40
552 guez 309 sens = sum(flux_t * pctsrf, dim = 2)
553 guez 206 evap = - sum(flux_q * pctsrf, dim = 2)
554 guez 302 fder = dlw + dflux_t + dflux_q
555 guez 3
556     DO k = 1, llm
557     DO i = 1, klon
558     t_seri(i, k) = t_seri(i, k) + d_t_vdf(i, k)
559     q_seri(i, k) = q_seri(i, k) + d_q_vdf(i, k)
560     u_seri(i, k) = u_seri(i, k) + d_u_vdf(i, k)
561     v_seri(i, k) = v_seri(i, k) + d_v_vdf(i, k)
562     ENDDO
563     ENDDO
564    
565 guez 191 call assert(abs(sum(pctsrf, dim = 2) - 1.) <= EPSFRA, 'physiq: pctsrf')
566 guez 221 tsol = sum(ftsol * pctsrf, dim = 2)
567 guez 208 zxfluxlat = sum(fluxlat * pctsrf, dim = 2)
568     zt2m = sum(t2m * pctsrf, dim = 2)
569     zq2m = sum(q2m * pctsrf, dim = 2)
570 guez 225 u10m = sum(u10m_srf * pctsrf, dim = 2)
571     v10m = sum(v10m_srf * pctsrf, dim = 2)
572 guez 208 zxffonte = sum(ffonte * pctsrf, dim = 2)
573     s_pblh = sum(pblh * pctsrf, dim = 2)
574     s_lcl = sum(plcl * pctsrf, dim = 2)
575     s_capCL = sum(capCL * pctsrf, dim = 2)
576     s_oliqCL = sum(oliqCL * pctsrf, dim = 2)
577     s_cteiCL = sum(cteiCL * pctsrf, dim = 2)
578     s_pblT = sum(pblT * pctsrf, dim = 2)
579     s_therm = sum(therm * pctsrf, dim = 2)
580 guez 3
581 guez 205 ! Si une sous-fraction n'existe pas, elle prend la valeur moyenne :
582 guez 3 DO nsrf = 1, nbsrf
583     DO i = 1, klon
584 guez 205 IF (pctsrf(i, nsrf) < epsfra) then
585 guez 221 ftsol(i, nsrf) = tsol(i)
586 guez 205 t2m(i, nsrf) = zt2m(i)
587     q2m(i, nsrf) = zq2m(i)
588 guez 225 u10m_srf(i, nsrf) = u10m(i)
589     v10m_srf(i, nsrf) = v10m(i)
590 guez 205 ffonte(i, nsrf) = zxffonte(i)
591     pblh(i, nsrf) = s_pblh(i)
592     plcl(i, nsrf) = s_lcl(i)
593     capCL(i, nsrf) = s_capCL(i)
594     oliqCL(i, nsrf) = s_oliqCL(i)
595     cteiCL(i, nsrf) = s_cteiCL(i)
596     pblT(i, nsrf) = s_pblT(i)
597     therm(i, nsrf) = s_therm(i)
598     end IF
599 guez 3 ENDDO
600     ENDDO
601    
602 guez 223 dlw = - 4. * RSIGMA * tsol**3
603 guez 3
604 guez 190 ! Appeler la convection
605 guez 3
606 guez 182 if (conv_emanuel) then
607 guez 195 CALL concvl(paprs, play, t_seri, q_seri, u_seri, v_seri, sig1, w01, &
608     d_t_con, d_q_con, d_u_con, d_v_con, rain_con, ibas_con, itop_con, &
609 guez 266 upwd, dnwd, Ma, cape, iflagctrl, clwcon0, pmflxr, da, phi, mp)
610 guez 183 snow_con = 0.
611 guez 71 mfu = upwd + dnwd
612 guez 3
613 guez 207 zqsat = MIN(0.5, r2es * FOEEW(t_seri, rtt >= t_seri) / play)
614     zqsat = zqsat / (1. - retv * zqsat)
615 guez 3
616 guez 103 ! Properties of convective clouds
617 guez 71 clwcon0 = fact_cldcon * clwcon0
618 guez 62 call clouds_gno(klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, &
619     rnebcon0)
620 guez 72
621 guez 190 forall (i = 1:klon) ema_pct(i) = paprs(i, itop_con(i) + 1)
622 guez 72 mfd = 0.
623     pen_u = 0.
624     pen_d = 0.
625     pde_d = 0.
626     pde_u = 0.
627 guez 182 else
628     conv_q = d_q_dyn + d_q_vdf / dtphys
629     conv_t = d_t_dyn + d_t_vdf / dtphys
630     z_avant = sum((q_seri + ql_seri) * zmasse, dim=2)
631 guez 298 CALL conflx(paprs, play, t_seri(:, llm:1:- 1), q_seri(:, llm:1:- 1), &
632     conv_t, conv_q, - evap, omega, d_t_con, d_q_con, rain_con, &
633     snow_con, mfu(:, llm:1:- 1), mfd(:, llm:1:- 1), pen_u, pde_u, &
634     pen_d, pde_d, kcbot, kctop, kdtop, pmflxr, pmflxs)
635 guez 182 WHERE (rain_con < 0.) rain_con = 0.
636     WHERE (snow_con < 0.) snow_con = 0.
637     ibas_con = llm + 1 - kcbot
638     itop_con = llm + 1 - kctop
639 guez 69 END if
640 guez 3
641     DO k = 1, llm
642     DO i = 1, klon
643     t_seri(i, k) = t_seri(i, k) + d_t_con(i, k)
644     q_seri(i, k) = q_seri(i, k) + d_q_con(i, k)
645     u_seri(i, k) = u_seri(i, k) + d_u_con(i, k)
646     v_seri(i, k) = v_seri(i, k) + d_v_con(i, k)
647     ENDDO
648     ENDDO
649    
650 guez 182 IF (.not. conv_emanuel) THEN
651 guez 69 z_apres = sum((q_seri + ql_seri) * zmasse, dim=2)
652     z_factor = (z_avant - (rain_con + snow_con) * dtphys) / z_apres
653 guez 3 DO k = 1, llm
654     DO i = 1, klon
655 guez 325 IF (z_factor(i) /= 1.) THEN
656 guez 3 q_seri(i, k) = q_seri(i, k) * z_factor(i)
657     ENDIF
658     ENDDO
659     ENDDO
660     ENDIF
661    
662 guez 90 ! Convection s\`eche (thermiques ou ajustement)
663 guez 3
664 guez 51 d_t_ajs = 0.
665     d_u_ajs = 0.
666     d_v_ajs = 0.
667     d_q_ajs = 0.
668     fm_therm = 0.
669     entr_therm = 0.
670 guez 3
671 guez 323 if (iflag_thermals) then
672     call calltherm(play, paprs, pphi, u_seri, v_seri, t_seri, q_seri, &
673     d_u_ajs, d_v_ajs, d_t_ajs, d_q_ajs, fm_therm, entr_therm)
674     else
675 guez 47 CALL ajsec(paprs, play, t_seri, q_seri, d_t_ajs, d_q_ajs)
676 guez 13 t_seri = t_seri + d_t_ajs
677     q_seri = q_seri + d_q_ajs
678 guez 3 endif
679    
680 guez 47 ! Caclul des ratqs
681 guez 3
682     if (iflag_cldcon == 1) then
683 guez 266 ! ratqs convectifs \`a l'ancienne en fonction de (q(z = 0) - q) / q
684     ! on \'ecrase le tableau ratqsc calcul\'e par clouds_gno
685 guez 51 do k = 1, llm
686     do i = 1, klon
687 guez 3 if(ptconv(i, k)) then
688 guez 70 ratqsc(i, k) = ratqsbas + fact_cldcon &
689     * (q_seri(i, 1) - q_seri(i, k)) / q_seri(i, k)
690 guez 3 else
691 guez 51 ratqsc(i, k) = 0.
692 guez 3 endif
693     enddo
694     enddo
695     endif
696    
697 guez 47 ! ratqs stables
698 guez 51 do k = 1, llm
699     do i = 1, klon
700 guez 70 ratqss(i, k) = ratqsbas + (ratqshaut - ratqsbas) &
701 guez 190 * min((paprs(i, 1) - play(i, k)) / (paprs(i, 1) - 3e4), 1.)
702 guez 3 enddo
703     enddo
704    
705 guez 47 ! ratqs final
706 guez 69 if (iflag_cldcon == 1 .or. iflag_cldcon == 2) then
707 guez 47 ! les ratqs sont une conbinaison de ratqss et ratqsc
708     ! ratqs final
709     ! 1e4 (en gros 3 heures), en dur pour le moment, est le temps de
710     ! relaxation des ratqs
711 guez 70 ratqs = max(ratqs * exp(- dtphys * facttemps), ratqss)
712 guez 51 ratqs = max(ratqs, ratqsc)
713 guez 3 else
714 guez 47 ! on ne prend que le ratqs stable pour fisrtilp
715 guez 51 ratqs = ratqss
716 guez 3 endif
717    
718 guez 298 CALL fisrtilp(paprs, play, t_seri, q_seri, ptconv, ratqs, d_t_lsc, &
719 guez 266 d_q_lsc, d_ql_lsc, rneb, cldliq, rain_lsc, snow_lsc, pfrac_impa, &
720     pfrac_nucl, pfrac_1nucl, frac_impa, frac_nucl, prfl, psfl, rhcl)
721 guez 3
722     WHERE (rain_lsc < 0) rain_lsc = 0.
723     WHERE (snow_lsc < 0) snow_lsc = 0.
724     DO k = 1, llm
725     DO i = 1, klon
726     t_seri(i, k) = t_seri(i, k) + d_t_lsc(i, k)
727     q_seri(i, k) = q_seri(i, k) + d_q_lsc(i, k)
728     ql_seri(i, k) = ql_seri(i, k) + d_ql_lsc(i, k)
729     cldfra(i, k) = rneb(i, k)
730     IF (.NOT.new_oliq) cldliq(i, k) = ql_seri(i, k)
731     ENDDO
732     ENDDO
733    
734 guez 47 ! PRESCRIPTION DES NUAGES POUR LE RAYONNEMENT
735 guez 3
736     ! 1. NUAGES CONVECTIFS
737    
738 guez 174 IF (iflag_cldcon <= - 1) THEN
739 guez 62 ! seulement pour Tiedtke
740 guez 51 snow_tiedtke = 0.
741 guez 174 if (iflag_cldcon == - 1) then
742 guez 51 rain_tiedtke = rain_con
743 guez 3 else
744 guez 51 rain_tiedtke = 0.
745     do k = 1, llm
746     do i = 1, klon
747 guez 7 if (d_q_con(i, k) < 0.) then
748 guez 202 rain_tiedtke(i) = rain_tiedtke(i) - d_q_con(i, k) / dtphys &
749     * zmasse(i, k)
750 guez 3 endif
751     enddo
752     enddo
753     endif
754    
755     ! Nuages diagnostiques pour Tiedtke
756 guez 69 CALL diagcld1(paprs, play, rain_tiedtke, snow_tiedtke, ibas_con, &
757     itop_con, diafra, dialiq)
758 guez 3 DO k = 1, llm
759     DO i = 1, klon
760 guez 51 IF (diafra(i, k) > cldfra(i, k)) THEN
761 guez 3 cldliq(i, k) = dialiq(i, k)
762     cldfra(i, k) = diafra(i, k)
763     ENDIF
764     ENDDO
765     ENDDO
766     ELSE IF (iflag_cldcon == 3) THEN
767 guez 72 ! On prend pour les nuages convectifs le maximum du calcul de
768 guez 90 ! la convection et du calcul du pas de temps pr\'ec\'edent diminu\'e
769 guez 72 ! d'un facteur facttemps.
770     facteur = dtphys * facttemps
771 guez 51 do k = 1, llm
772     do i = 1, klon
773 guez 70 rnebcon(i, k) = rnebcon(i, k) * facteur
774 guez 72 if (rnebcon0(i, k) * clwcon0(i, k) &
775     > rnebcon(i, k) * clwcon(i, k)) then
776 guez 51 rnebcon(i, k) = rnebcon0(i, k)
777     clwcon(i, k) = clwcon0(i, k)
778 guez 3 endif
779     enddo
780     enddo
781    
782 guez 47 ! On prend la somme des fractions nuageuses et des contenus en eau
783 guez 51 cldfra = min(max(cldfra, rnebcon), 1.)
784 guez 202 cldliq = cldliq + rnebcon * clwcon
785 guez 3 ENDIF
786    
787 guez 51 ! 2. Nuages stratiformes
788 guez 3
789     IF (ok_stratus) THEN
790 guez 47 CALL diagcld2(paprs, play, t_seri, q_seri, diafra, dialiq)
791 guez 3 DO k = 1, llm
792     DO i = 1, klon
793 guez 51 IF (diafra(i, k) > cldfra(i, k)) THEN
794 guez 3 cldliq(i, k) = dialiq(i, k)
795     cldfra(i, k) = diafra(i, k)
796     ENDIF
797     ENDDO
798     ENDDO
799     ENDIF
800    
801     ! Precipitation totale
802     DO i = 1, klon
803     rain_fall(i) = rain_con(i) + rain_lsc(i)
804     snow_fall(i) = snow_con(i) + snow_lsc(i)
805     ENDDO
806    
807 guez 90 ! Humidit\'e relative pour diagnostic :
808 guez 3 DO k = 1, llm
809     DO i = 1, klon
810 guez 323 zx_qs = r2es * FOEEW(t_seri(i, k), rtt >= t_seri(i, k)) / play(i, k)
811 guez 207 zx_qs = MIN(0.5, zx_qs)
812     zcor = 1. / (1. - retv * zx_qs)
813     zx_qs = zx_qs * zcor
814 guez 202 zx_rh(i, k) = q_seri(i, k) / zx_qs
815 guez 51 zqsat(i, k) = zx_qs
816 guez 3 ENDDO
817     ENDDO
818 guez 52
819 guez 97 ! Param\`etres optiques des nuages et quelques param\`etres pour
820     ! diagnostics :
821 guez 3 if (ok_newmicro) then
822 guez 69 CALL newmicro(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, &
823 guez 217 cldh, cldl, cldm, cldt, cldq, flwp, fiwp, flwc, fiwc)
824 guez 3 else
825 guez 52 CALL nuage(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, cldh, &
826 guez 217 cldl, cldm, cldt, cldq)
827 guez 3 endif
828    
829 guez 154 IF (MOD(itap - 1, radpas) == 0) THEN
830 guez 212 wo = ozonecm(REAL(julien), paprs)
831 guez 155 albsol = sum(falbe * pctsrf, dim = 2)
832 guez 221 CALL radlwsw(dist, mu0, fract, paprs, play, tsol, albsol, t_seri, &
833 guez 155 q_seri, wo, cldfra, cldemi, cldtau, heat, heat0, cool, cool0, &
834     radsol, albpla, topsw, toplw, solsw, sollw, sollwdown, topsw0, &
835     toplw0, solsw0, sollw0, lwdn0, lwdn, lwup0, lwup, swdn0, swdn, &
836 guez 217 swup0, swup, ok_ade, topswad, solswad)
837 guez 3 ENDIF
838 guez 118
839 guez 3 ! Ajouter la tendance des rayonnements (tous les pas)
840     DO k = 1, llm
841     DO i = 1, klon
842 guez 202 t_seri(i, k) = t_seri(i, k) + (heat(i, k) - cool(i, k)) * dtphys &
843     / 86400.
844 guez 3 ENDDO
845     ENDDO
846    
847 guez 90 ! Calculer le bilan du sol et la d\'erive de temp\'erature (couplage)
848 guez 3 DO i = 1, klon
849 guez 309 bils(i) = radsol(i) + sens(i) + zxfluxlat(i)
850 guez 3 ENDDO
851    
852 guez 90 ! Param\'etrisation de l'orographie \`a l'\'echelle sous-maille :
853 guez 3
854     IF (ok_orodr) THEN
855 guez 174 ! S\'election des points pour lesquels le sch\'ema est actif :
856 guez 51 DO i = 1, klon
857 guez 247 ktest(i) = 0
858 guez 174 IF (zpic(i) - zmea(i) > 100. .AND. zstd(i) > 10.) THEN
859 guez 247 ktest(i) = 1
860 guez 3 ENDIF
861     ENDDO
862    
863 guez 298 CALL drag_noro(paprs, play, zmea, zstd, zsig, zgam, zthe, zpic, zval, &
864     ktest, t_seri, u_seri, v_seri, zulow, zvlow, zustrdr, zvstrdr, &
865     d_t_oro, d_u_oro, d_v_oro)
866 guez 3
867 guez 47 ! ajout des tendances
868 guez 3 DO k = 1, llm
869     DO i = 1, klon
870     t_seri(i, k) = t_seri(i, k) + d_t_oro(i, k)
871     u_seri(i, k) = u_seri(i, k) + d_u_oro(i, k)
872     v_seri(i, k) = v_seri(i, k) + d_v_oro(i, k)
873     ENDDO
874     ENDDO
875 guez 13 ENDIF
876 guez 3
877     IF (ok_orolf) THEN
878 guez 90 ! S\'election des points pour lesquels le sch\'ema est actif :
879 guez 51 DO i = 1, klon
880 guez 247 ktest(i) = 0
881 guez 174 IF (zpic(i) - zmea(i) > 100.) THEN
882 guez 247 ktest(i) = 1
883 guez 3 ENDIF
884     ENDDO
885    
886 guez 298 CALL lift_noro(paprs, play, zmea, zstd, zpic, ktest, t_seri, u_seri, &
887     v_seri, zulow, zvlow, zustrli, zvstrli, d_t_lif, d_u_lif, d_v_lif)
888 guez 3
889 guez 51 ! Ajout des tendances :
890 guez 3 DO k = 1, llm
891     DO i = 1, klon
892     t_seri(i, k) = t_seri(i, k) + d_t_lif(i, k)
893     u_seri(i, k) = u_seri(i, k) + d_u_lif(i, k)
894     v_seri(i, k) = v_seri(i, k) + d_v_lif(i, k)
895     ENDDO
896     ENDDO
897 guez 49 ENDIF
898 guez 3
899 guez 229 CALL aaam_bud(rg, romega, pphis, zustrdr, zustrli, &
900     sum((u_seri - u) / dtphys * zmasse, dim = 2), zvstrdr, &
901     zvstrli, sum((v_seri - v) / dtphys * zmasse, dim = 2), paprs, u, v, &
902     aam, torsfc)
903 guez 3
904 guez 47 ! Calcul des tendances traceurs
905 guez 298 call phytrac(julien, time, firstcal, lafin, t, paprs, play, mfu, mfd, &
906     pde_u, pen_d, coefh, cdragh, fm_therm, entr_therm, u(:, 1), v(:, 1), &
907     ftsol, pctsrf, frac_impa, frac_nucl, da, phi, mp, upwd, dnwd, &
908     tr_seri, zmasse, ncid_startphy)
909 guez 3
910     ! Calculer le transport de l'eau et de l'energie (diagnostique)
911 guez 171 CALL transp(paprs, t_seri, q_seri, u_seri, v_seri, zphi, ve, vq, ue, uq)
912 guez 3
913 guez 31 ! diag. bilKP
914 guez 3
915 guez 324 CALL transp_lay(paprs, t_seri, q_seri, u_seri, v_seri, zphi, ve_lay, &
916     vq_lay, ue_lay, uq_lay)
917 guez 3
918     ! Accumuler les variables a stocker dans les fichiers histoire:
919    
920 guez 200 ! conversion Ec en énergie thermique
921 guez 3 DO k = 1, llm
922     DO i = 1, klon
923 guez 213 d_t_ec(i, k) = 0.5 / (RCPD * (1. + RVTMP2 * q_seri(i, k))) &
924 guez 51 * (u(i, k)**2 + v(i, k)**2 - u_seri(i, k)**2 - v_seri(i, k)**2)
925     t_seri(i, k) = t_seri(i, k) + d_t_ec(i, k)
926     d_t_ec(i, k) = d_t_ec(i, k) / dtphys
927 guez 3 END DO
928     END DO
929 guez 51
930 guez 47 ! SORTIES
931 guez 3
932 guez 69 ! prw = eau precipitable
933 guez 3 DO i = 1, klon
934     prw(i) = 0.
935     DO k = 1, llm
936 guez 202 prw(i) = prw(i) + q_seri(i, k) * zmasse(i, k)
937 guez 3 ENDDO
938     ENDDO
939    
940     ! Convertir les incrementations en tendances
941    
942     DO k = 1, llm
943     DO i = 1, klon
944 guez 49 d_u(i, k) = (u_seri(i, k) - u(i, k)) / dtphys
945     d_v(i, k) = (v_seri(i, k) - v(i, k)) / dtphys
946     d_t(i, k) = (t_seri(i, k) - t(i, k)) / dtphys
947     d_qx(i, k, ivap) = (q_seri(i, k) - qx(i, k, ivap)) / dtphys
948     d_qx(i, k, iliq) = (ql_seri(i, k) - qx(i, k, iliq)) / dtphys
949 guez 3 ENDDO
950     ENDDO
951    
952 guez 98 DO iq = 3, nqmx
953     DO k = 1, llm
954     DO i = 1, klon
955 guez 174 d_qx(i, k, iq) = (tr_seri(i, k, iq - 2) - qx(i, k, iq)) / dtphys
956 guez 3 ENDDO
957     ENDDO
958 guez 98 ENDDO
959 guez 3
960     ! Sauvegarder les valeurs de t et q a la fin de la physique:
961     DO k = 1, llm
962     DO i = 1, klon
963     t_ancien(i, k) = t_seri(i, k)
964     q_ancien(i, k) = q_seri(i, k)
965     ENDDO
966     ENDDO
967    
968 guez 191 CALL histwrite_phy("phis", pphis)
969     CALL histwrite_phy("aire", airephy)
970     CALL histwrite_phy("psol", paprs(:, 1))
971     CALL histwrite_phy("precip", rain_fall + snow_fall)
972     CALL histwrite_phy("plul", rain_lsc + snow_lsc)
973     CALL histwrite_phy("pluc", rain_con + snow_con)
974 guez 221 CALL histwrite_phy("tsol", tsol)
975 guez 191 CALL histwrite_phy("t2m", zt2m)
976     CALL histwrite_phy("q2m", zq2m)
977 guez 225 CALL histwrite_phy("u10m", u10m)
978     CALL histwrite_phy("v10m", v10m)
979 guez 191 CALL histwrite_phy("snow", snow_fall)
980     CALL histwrite_phy("cdrm", cdragm)
981     CALL histwrite_phy("cdrh", cdragh)
982     CALL histwrite_phy("topl", toplw)
983     CALL histwrite_phy("evap", evap)
984     CALL histwrite_phy("sols", solsw)
985 guez 308 CALL histwrite_phy("rls", sollw)
986 guez 191 CALL histwrite_phy("solldown", sollwdown)
987     CALL histwrite_phy("bils", bils)
988 guez 309 CALL histwrite_phy("sens", sens)
989 guez 191 CALL histwrite_phy("fder", fder)
990 guez 279 CALL histwrite_phy("zxfqcalving", sum(fqcalving * pctsrf, dim = 2))
991 guez 191 CALL histwrite_phy("albs", albsol)
992 guez 212 CALL histwrite_phy("tro3", wo * dobson_u * 1e3 / zmasse / rmo3 * md)
993 guez 191 CALL histwrite_phy("rugs", zxrugs)
994     CALL histwrite_phy("s_pblh", s_pblh)
995     CALL histwrite_phy("s_pblt", s_pblt)
996     CALL histwrite_phy("s_lcl", s_lcl)
997     CALL histwrite_phy("s_capCL", s_capCL)
998     CALL histwrite_phy("s_oliqCL", s_oliqCL)
999     CALL histwrite_phy("s_cteiCL", s_cteiCL)
1000     CALL histwrite_phy("s_therm", s_therm)
1001     CALL histwrite_phy("temp", t_seri)
1002     CALL histwrite_phy("vitu", u_seri)
1003     CALL histwrite_phy("vitv", v_seri)
1004     CALL histwrite_phy("geop", zphi)
1005     CALL histwrite_phy("pres", play)
1006     CALL histwrite_phy("dtvdf", d_t_vdf)
1007     CALL histwrite_phy("dqvdf", d_q_vdf)
1008     CALL histwrite_phy("rhum", zx_rh)
1009 guez 213 CALL histwrite_phy("d_t_ec", d_t_ec)
1010     CALL histwrite_phy("dtsw0", heat0 / 86400.)
1011     CALL histwrite_phy("dtlw0", - cool0 / 86400.)
1012 guez 326 call histwrite_phy("pmflxr", pmflxr(:, :llm))
1013 guez 215 CALL histwrite_phy("msnow", sum(fsnow * pctsrf, dim = 2))
1014 guez 221 call histwrite_phy("qsurf", sum(fqsurf * pctsrf, dim = 2))
1015 guez 306 call histwrite_phy("flat", zxfluxlat)
1016 guez 191
1017 guez 306 DO nsrf = 1, nbsrf
1018     CALL histwrite_phy("fract_"//clnsurf(nsrf), pctsrf(:, nsrf))
1019     CALL histwrite_phy("sens_"//clnsurf(nsrf), flux_t(:, nsrf))
1020     CALL histwrite_phy("lat_"//clnsurf(nsrf), fluxlat(:, nsrf))
1021     CALL histwrite_phy("tsol_"//clnsurf(nsrf), ftsol(:, nsrf))
1022     CALL histwrite_phy("taux_"//clnsurf(nsrf), flux_u(:, nsrf))
1023     CALL histwrite_phy("tauy_"//clnsurf(nsrf), flux_v(:, nsrf))
1024     CALL histwrite_phy("rugs_"//clnsurf(nsrf), frugs(:, nsrf))
1025     CALL histwrite_phy("albe_"//clnsurf(nsrf), falbe(:, nsrf))
1026     CALL histwrite_phy("u10m_"//clnsurf(nsrf), u10m_srf(:, nsrf))
1027     CALL histwrite_phy("v10m_"//clnsurf(nsrf), v10m_srf(:, nsrf))
1028     END DO
1029    
1030     if (conv_emanuel) then
1031     CALL histwrite_phy("ptop", ema_pct)
1032     CALL histwrite_phy("dnwd0", - mp)
1033     end if
1034    
1035 guez 191 if (ok_instan) call histsync(nid_ins)
1036    
1037 guez 157 IF (lafin) then
1038     call NF95_CLOSE(ncid_startphy)
1039 guez 304 CALL phyredem(pctsrf, ftsol, ftsoil, fqsurf, qsol, fsnow, falbe, &
1040     rain_fall, snow_fall, solsw, sollw, dlw, radsol, frugs, agesno, &
1041     zmea, zstd, zsig, zgam, zthe, zpic, zval, t_ancien, q_ancien, &
1042     rnebcon, ratqs, clwcon, run_off_lic_0, sig1, w01)
1043 guez 157 end IF
1044 guez 3
1045 guez 35 firstcal = .FALSE.
1046    
1047 guez 3 END SUBROUTINE physiq
1048    
1049     end module physiq_m

  ViewVC Help
Powered by ViewVC 1.1.21