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

Diff of /trunk/Sources/phylmd/physiq.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

trunk/phylmd/physiq.f revision 101 by guez, Mon Jul 7 17:45:21 2014 UTC trunk/Sources/phylmd/physiq.f revision 195 by guez, Wed May 18 17:56:44 2016 UTC
# Line 4  module physiq_m Line 4  module physiq_m
4    
5  contains  contains
6    
7    SUBROUTINE physiq(lafin, rdayvrai, time, dtphys, paprs, play, pphi, pphis, &    SUBROUTINE physiq(lafin, dayvrai, time, paprs, play, pphi, pphis, u, v, t, &
8         u, v, t, qx, omega, d_u, d_v, d_t, d_qx)         qx, omega, d_u, d_v, d_t, d_qx)
9    
10      ! From phylmd/physiq.F, version 1.22 2006/02/20 09:38:28      ! From phylmd/physiq.F, version 1.22 2006/02/20 09:38:28
11      ! (subversion revision 678)      ! (subversion revision 678)
12    
13      ! Author: Z.X. Li (LMD/CNRS) 1993      ! Author: Z. X. Li (LMD/CNRS) 1993
14    
15      ! This is the main procedure for the "physics" part of the program.      ! This is the main procedure for the "physics" part of the program.
16    
17      use aaam_bud_m, only: aaam_bud      use aaam_bud_m, only: aaam_bud
18      USE abort_gcm_m, ONLY: abort_gcm      USE abort_gcm_m, ONLY: abort_gcm
     use aeropt_m, only: aeropt  
19      use ajsec_m, only: ajsec      use ajsec_m, only: ajsec
20      use calltherm_m, only: calltherm      use calltherm_m, only: calltherm
21      USE clesphys, ONLY: cdhmax, cdmmax, co2_ppm, ecrit_hf, ecrit_ins, &      USE clesphys, ONLY: cdhmax, cdmmax, ecrit_hf, ecrit_ins, ecrit_mth, &
22           ecrit_mth, ecrit_reg, ecrit_tra, ksta, ksta_ter, ok_kzmin           ecrit_reg, ecrit_tra, ksta, ksta_ter, ok_kzmin, ok_instan
23      USE clesphys2, ONLY: cycle_diurne, iflag_con, nbapp_rad, new_oliq, &      USE clesphys2, ONLY: cycle_diurne, conv_emanuel, nbapp_rad, new_oliq, &
24           ok_orodr, ok_orolf           ok_orodr, ok_orolf
25      USE clmain_m, ONLY: clmain      USE clmain_m, ONLY: clmain
26      use clouds_gno_m, only: clouds_gno      use clouds_gno_m, only: clouds_gno
27      USE comgeomphy, ONLY: airephy, cuphy, cvphy      use comconst, only: dtphys
28        USE comgeomphy, ONLY: airephy
29      USE concvl_m, ONLY: concvl      USE concvl_m, ONLY: concvl
30      USE conf_gcm_m, ONLY: offline, raz_date      USE conf_gcm_m, ONLY: offline, day_step, iphysiq
31      USE conf_phys_m, ONLY: conf_phys      USE conf_phys_m, ONLY: conf_phys
32      use conflx_m, only: conflx      use conflx_m, only: conflx
33      USE ctherm, ONLY: iflag_thermals, nsplit_thermals      USE ctherm, ONLY: iflag_thermals, nsplit_thermals
# Line 38  contains Line 38  contains
38      USE dimphy, ONLY: klon      USE dimphy, ONLY: klon
39      USE dimsoil, ONLY: nsoilmx      USE dimsoil, ONLY: nsoilmx
40      use drag_noro_m, only: drag_noro      use drag_noro_m, only: drag_noro
41        use dynetat0_m, only: day_ref, annee_ref
42      USE fcttre, ONLY: foeew, qsatl, qsats, thermcep      USE fcttre, ONLY: foeew, qsatl, qsats, thermcep
43      use fisrtilp_m, only: fisrtilp      use fisrtilp_m, only: fisrtilp
44      USE hgardfou_m, ONLY: hgardfou      USE hgardfou_m, ONLY: hgardfou
45        USE histsync_m, ONLY: histsync
46        USE histwrite_phy_m, ONLY: histwrite_phy
47      USE indicesol, ONLY: clnsurf, epsfra, is_lic, is_oce, is_sic, is_ter, &      USE indicesol, ONLY: clnsurf, epsfra, is_lic, is_oce, is_sic, is_ter, &
48           nbsrf           nbsrf
49      USE ini_histins_m, ONLY: ini_histins      USE ini_histins_m, ONLY: ini_histins, nid_ins
50        use netcdf95, only: NF95_CLOSE
51      use newmicro_m, only: newmicro      use newmicro_m, only: newmicro
52        use nr_util, only: assert
53        use nuage_m, only: nuage
54      USE orbite_m, ONLY: orbite      USE orbite_m, ONLY: orbite
55      USE ozonecm_m, ONLY: ozonecm      USE ozonecm_m, ONLY: ozonecm
56      USE phyetat0_m, ONLY: phyetat0, rlat, rlon      USE phyetat0_m, ONLY: phyetat0, rlat, rlon
57      USE phyredem_m, ONLY: phyredem      USE phyredem_m, ONLY: phyredem
58        USE phyredem0_m, ONLY: phyredem0
59      USE phystokenc_m, ONLY: phystokenc      USE phystokenc_m, ONLY: phystokenc
60      USE phytrac_m, ONLY: phytrac      USE phytrac_m, ONLY: phytrac
61      USE qcheck_m, ONLY: qcheck      USE qcheck_m, ONLY: qcheck
62      use radlwsw_m, only: radlwsw      use radlwsw_m, only: radlwsw
63      use readsulfate_m, only: readsulfate      use yoegwd, only: sugwd
64      use sugwd_m, only: sugwd      USE suphec_m, ONLY: rcpd, retv, rg, rlvtt, romega, rsigma, rtt
65      USE suphec_m, ONLY: ra, rcpd, retv, rg, rlvtt, romega, rsigma, rtt      use time_phylmdz, only: itap, increment_itap
66      USE temps, ONLY: annee_ref, day_ref, itau_phy      use transp_m, only: transp
67        use transp_lay_m, only: transp_lay
68      use unit_nml_m, only: unit_nml      use unit_nml_m, only: unit_nml
69      USE ymds2ju_m, ONLY: ymds2ju      USE ymds2ju_m, ONLY: ymds2ju
70      USE yoethf_m, ONLY: r2es, rvtmp2      USE yoethf_m, ONLY: r2es, rvtmp2
# Line 64  contains Line 72  contains
72    
73      logical, intent(in):: lafin ! dernier passage      logical, intent(in):: lafin ! dernier passage
74    
75      REAL, intent(in):: rdayvrai      integer, intent(in):: dayvrai
76      ! (elapsed time since January 1st 0h of the starting year, in days)      ! current day number, based at value 1 on January 1st of annee_ref
77    
78      REAL, intent(in):: time ! heure de la journ\'ee en fraction de jour      REAL, intent(in):: time ! heure de la journ\'ee en fraction de jour
     REAL, intent(in):: dtphys ! pas d'integration pour la physique (seconde)  
79    
80      REAL, intent(in):: paprs(:, :) ! (klon, llm + 1)      REAL, intent(in):: paprs(:, :) ! (klon, llm + 1)
81      ! pression pour chaque inter-couche, en Pa      ! pression pour chaque inter-couche, en Pa
# Line 76  contains Line 83  contains
83      REAL, intent(in):: play(:, :) ! (klon, llm)      REAL, intent(in):: play(:, :) ! (klon, llm)
84      ! pression pour le mileu de chaque couche (en Pa)      ! pression pour le mileu de chaque couche (en Pa)
85    
86      REAL, intent(in):: pphi(:, :) ! (klon, llm)      REAL, intent(in):: pphi(:, :) ! (klon, llm)
87      ! géopotentiel de chaque couche (référence sol)      ! géopotentiel de chaque couche (référence sol)
88    
89      REAL, intent(in):: pphis(:) ! (klon) géopotentiel du sol      REAL, intent(in):: pphis(:) ! (klon) géopotentiel du sol
# Line 102  contains Line 109  contains
109    
110      LOGICAL:: firstcal = .true.      LOGICAL:: firstcal = .true.
111    
112      LOGICAL ok_gust ! pour activer l'effet des gust sur flux surface      LOGICAL, PARAMETER:: check = .FALSE.
     PARAMETER (ok_gust = .FALSE.)  
   
     LOGICAL, PARAMETER:: check = .FALSE.  
113      ! Verifier la conservation du modele en eau      ! Verifier la conservation du modele en eau
114    
115      LOGICAL, PARAMETER:: ok_stratus = .FALSE.      LOGICAL, PARAMETER:: ok_stratus = .FALSE.
116      ! Ajouter artificiellement les stratus      ! Ajouter artificiellement les stratus
117    
     ! "slab" ocean  
     REAL, save:: tslab(klon) ! temperature of ocean slab  
     REAL, save:: seaice(klon) ! glace de mer (kg/m2)  
     REAL fluxo(klon) ! flux turbulents ocean-glace de mer  
     REAL fluxg(klon) ! flux turbulents ocean-atmosphere  
   
     logical:: ok_journe = .false., ok_mensuel = .true., ok_instan = .false.  
     ! sorties journalieres, mensuelles et instantanees dans les  
     ! fichiers histday, histmth et histins  
   
     LOGICAL ok_region ! sortir le fichier regional  
     PARAMETER (ok_region = .FALSE.)  
   
118      ! pour phsystoke avec thermiques      ! pour phsystoke avec thermiques
119      REAL fm_therm(klon, llm + 1)      REAL fm_therm(klon, llm + 1)
120      REAL entr_therm(klon, llm)      REAL entr_therm(klon, llm)
# Line 146  contains Line 137  contains
137    
138      REAL lwdn0(klon, llm + 1), lwdn(klon, llm + 1)      REAL lwdn0(klon, llm + 1), lwdn(klon, llm + 1)
139      REAL lwup0(klon, llm + 1), lwup(klon, llm + 1)      REAL lwup0(klon, llm + 1), lwup(klon, llm + 1)
140      SAVE lwdn0, lwdn, lwup0, lwup      SAVE lwdn0, lwdn, lwup0, lwup
   
     ! Amip2  
     ! variables a une pression donnee  
   
     integer nlevSTD  
     PARAMETER(nlevSTD = 17)  
     real rlevSTD(nlevSTD)  
     DATA rlevSTD/100000., 92500., 85000., 70000., &  
          60000., 50000., 40000., 30000., 25000., 20000., &  
          15000., 10000., 7000., 5000., 3000., 2000., 1000./  
     CHARACTER(LEN = 4) clevSTD(nlevSTD)  
     DATA clevSTD/'1000', '925 ', '850 ', '700 ', '600 ', &  
          '500 ', '400 ', '300 ', '250 ', '200 ', '150 ', '100 ', &  
          '70 ', '50 ', '30 ', '20 ', '10 '/  
141    
142      ! prw: precipitable water      ! prw: precipitable water
143      real prw(klon)      real prw(klon)
# Line 170  contains Line 147  contains
147      REAL flwp(klon), fiwp(klon)      REAL flwp(klon), fiwp(klon)
148      REAL flwc(klon, llm), fiwc(klon, llm)      REAL flwc(klon, llm), fiwc(klon, llm)
149    
     INTEGER kmax, lmax  
     PARAMETER(kmax = 8, lmax = 8)  
     INTEGER kmaxm1, lmaxm1  
     PARAMETER(kmaxm1 = kmax-1, lmaxm1 = lmax-1)  
   
     REAL zx_tau(kmaxm1), zx_pc(lmaxm1)  
     DATA zx_tau/0., 0.3, 1.3, 3.6, 9.4, 23., 60./  
     DATA zx_pc/50., 180., 310., 440., 560., 680., 800./  
   
     ! cldtopres pression au sommet des nuages  
     REAL cldtopres(lmaxm1)  
     DATA cldtopres/50., 180., 310., 440., 560., 680., 800./  
   
     ! taulev: numero du niveau de tau dans les sorties ISCCP  
     CHARACTER(LEN = 4) taulev(kmaxm1)  
   
     DATA taulev/'tau0', 'tau1', 'tau2', 'tau3', 'tau4', 'tau5', 'tau6'/  
     CHARACTER(LEN = 3) pclev(lmaxm1)  
     DATA pclev/'pc1', 'pc2', 'pc3', 'pc4', 'pc5', 'pc6', 'pc7'/  
   
     CHARACTER(LEN = 28) cnameisccp(lmaxm1, kmaxm1)  
     DATA cnameisccp/'pc< 50hPa, tau< 0.3', 'pc= 50-180hPa, tau< 0.3', &  
          'pc= 180-310hPa, tau< 0.3', 'pc= 310-440hPa, tau< 0.3', &  
          'pc= 440-560hPa, tau< 0.3', 'pc= 560-680hPa, tau< 0.3', &  
          'pc= 680-800hPa, tau< 0.3', 'pc< 50hPa, tau= 0.3-1.3', &  
          'pc= 50-180hPa, tau= 0.3-1.3', 'pc= 180-310hPa, tau= 0.3-1.3', &  
          'pc= 310-440hPa, tau= 0.3-1.3', 'pc= 440-560hPa, tau= 0.3-1.3', &  
          'pc= 560-680hPa, tau= 0.3-1.3', 'pc= 680-800hPa, tau= 0.3-1.3', &  
          'pc< 50hPa, tau= 1.3-3.6', 'pc= 50-180hPa, tau= 1.3-3.6', &  
          'pc= 180-310hPa, tau= 1.3-3.6', 'pc= 310-440hPa, tau= 1.3-3.6', &  
          'pc= 440-560hPa, tau= 1.3-3.6', 'pc= 560-680hPa, tau= 1.3-3.6', &  
          'pc= 680-800hPa, tau= 1.3-3.6', 'pc< 50hPa, tau= 3.6-9.4', &  
          'pc= 50-180hPa, tau= 3.6-9.4', 'pc= 180-310hPa, tau= 3.6-9.4', &  
          'pc= 310-440hPa, tau= 3.6-9.4', 'pc= 440-560hPa, tau= 3.6-9.4', &  
          'pc= 560-680hPa, tau= 3.6-9.4', 'pc= 680-800hPa, tau= 3.6-9.4', &  
          'pc< 50hPa, tau= 9.4-23', 'pc= 50-180hPa, tau= 9.4-23', &  
          'pc= 180-310hPa, tau= 9.4-23', 'pc= 310-440hPa, tau= 9.4-23', &  
          'pc= 440-560hPa, tau= 9.4-23', 'pc= 560-680hPa, tau= 9.4-23', &  
          'pc= 680-800hPa, tau= 9.4-23', 'pc< 50hPa, tau= 23-60', &  
          'pc= 50-180hPa, tau= 23-60', 'pc= 180-310hPa, tau= 23-60', &  
          'pc= 310-440hPa, tau= 23-60', 'pc= 440-560hPa, tau= 23-60', &  
          'pc= 560-680hPa, tau= 23-60', 'pc= 680-800hPa, tau= 23-60', &  
          'pc< 50hPa, tau> 60.', 'pc= 50-180hPa, tau> 60.', &  
          'pc= 180-310hPa, tau> 60.', 'pc= 310-440hPa, tau> 60.', &  
          'pc= 440-560hPa, tau> 60.', 'pc= 560-680hPa, tau> 60.', &  
          'pc= 680-800hPa, tau> 60.'/  
   
     ! ISCCP simulator v3.4  
   
150      ! Variables propres a la physique      ! Variables propres a la physique
151    
152      INTEGER, save:: radpas      INTEGER, save:: radpas
153      ! (Radiative transfer computations are made every "radpas" call to      ! Radiative transfer computations are made every "radpas" call to
154      ! "physiq".)      ! "physiq".
155    
156      REAL radsol(klon)      REAL radsol(klon)
157      SAVE radsol ! bilan radiatif au sol calcule par code radiatif      SAVE radsol ! bilan radiatif au sol calcule par code radiatif
158    
     INTEGER, SAVE:: itap ! number of calls to "physiq"  
   
159      REAL, save:: ftsol(klon, nbsrf) ! skin temperature of surface fraction      REAL, save:: ftsol(klon, nbsrf) ! skin temperature of surface fraction
160    
161      REAL, save:: ftsoil(klon, nsoilmx, nbsrf)      REAL, save:: ftsoil(klon, nsoilmx, nbsrf)
# Line 246  contains Line 172  contains
172      ! column-density of water in soil, in kg m-2      ! column-density of water in soil, in kg m-2
173    
174      REAL, save:: fsnow(klon, nbsrf) ! epaisseur neigeuse      REAL, save:: fsnow(klon, nbsrf) ! epaisseur neigeuse
175      REAL, save:: falbe(klon, nbsrf) ! albedo par type de surface      REAL, save:: falbe(klon, nbsrf) ! albedo visible par type de surface
     REAL, save:: falblw(klon, nbsrf) ! albedo par type de surface  
176    
177      ! Param\`etres de l'orographie \`a l'\'echelle sous-maille (OESM) :      ! Param\`etres de l'orographie \`a l'\'echelle sous-maille (OESM) :
178      REAL, save:: zmea(klon) ! orographie moyenne      REAL, save:: zmea(klon) ! orographie moyenne
# Line 258  contains Line 183  contains
183      REAL, save:: zpic(klon) ! Maximum de l'OESM      REAL, save:: zpic(klon) ! Maximum de l'OESM
184      REAL, save:: zval(klon) ! Minimum de l'OESM      REAL, save:: zval(klon) ! Minimum de l'OESM
185      REAL, save:: rugoro(klon) ! longueur de rugosite de l'OESM      REAL, save:: rugoro(klon) ! longueur de rugosite de l'OESM
   
186      REAL zulow(klon), zvlow(klon)      REAL zulow(klon), zvlow(klon)
187        INTEGER igwd, itest(klon)
188    
189      INTEGER igwd, idx(klon), itest(klon)      REAL, save:: agesno(klon, nbsrf) ! age de la neige
190        REAL, save:: run_off_lic_0(klon)
191    
192      REAL agesno(klon, nbsrf)      ! Variables li\'ees \`a la convection d'Emanuel :
193      SAVE agesno ! age de la neige      REAL, save:: Ma(klon, llm) ! undilute upward mass flux
194        REAL, save:: qcondc(klon, llm) ! in-cld water content from convect
     REAL run_off_lic_0(klon)  
     SAVE run_off_lic_0  
     !KE43  
     ! Variables liees a la convection de K. Emanuel (sb):  
   
     REAL Ma(klon, llm) ! undilute upward mass flux  
     SAVE Ma  
     REAL qcondc(klon, llm) ! in-cld water content from convect  
     SAVE qcondc  
195      REAL, save:: sig1(klon, llm), w01(klon, llm)      REAL, save:: sig1(klon, llm), w01(klon, llm)
     REAL, save:: wd(klon)  
   
     ! Variables locales pour la couche limite (al1):  
   
     ! Variables locales:  
196    
197        ! Variables pour la couche limite (Alain Lahellec) :
198      REAL cdragh(klon) ! drag coefficient pour T and Q      REAL cdragh(klon) ! drag coefficient pour T and Q
199      REAL cdragm(klon) ! drag coefficient pour vent      REAL cdragm(klon) ! drag coefficient pour vent
200    
# Line 289  contains Line 202  contains
202      REAL ycoefh(klon, llm) ! coef d'echange pour phytrac      REAL ycoefh(klon, llm) ! coef d'echange pour phytrac
203      REAL yu1(klon) ! vents dans la premiere couche U      REAL yu1(klon) ! vents dans la premiere couche U
204      REAL yv1(klon) ! vents dans la premiere couche V      REAL yv1(klon) ! vents dans la premiere couche V
205      REAL ffonte(klon, nbsrf) !Flux thermique utilise pour fondre la neige      REAL ffonte(klon, nbsrf) ! flux thermique utilise pour fondre la neige
206      REAL fqcalving(klon, nbsrf) !Flux d'eau "perdue" par la surface  
207      ! !et necessaire pour limiter la      REAL fqcalving(klon, nbsrf)
208      ! !hauteur de neige, en kg/m2/s      ! flux d'eau "perdue" par la surface et necessaire pour limiter la
209        ! hauteur de neige, en kg/m2/s
210    
211      REAL zxffonte(klon), zxfqcalving(klon)      REAL zxffonte(klon), zxfqcalving(klon)
212    
213      REAL pfrac_impa(klon, llm)! Produits des coefs lessivage impaction      REAL pfrac_impa(klon, llm)! Produits des coefs lessivage impaction
# Line 317  contains Line 232  contains
232      REAL dlw(klon) ! derivee infra rouge      REAL dlw(klon) ! derivee infra rouge
233      SAVE dlw      SAVE dlw
234      REAL bils(klon) ! bilan de chaleur au sol      REAL bils(klon) ! bilan de chaleur au sol
235      REAL fder(klon) ! Derive de flux (sensible et latente)      REAL, save:: fder(klon) ! Derive de flux (sensible et latente)
     save fder  
236      REAL ve(klon) ! integr. verticale du transport meri. de l'energie      REAL ve(klon) ! integr. verticale du transport meri. de l'energie
237      REAL vq(klon) ! integr. verticale du transport meri. de l'eau      REAL vq(klon) ! integr. verticale du transport meri. de l'eau
238      REAL ue(klon) ! integr. verticale du transport zonal de l'energie      REAL ue(klon) ! integr. verticale du transport zonal de l'energie
# Line 333  contains Line 247  contains
247      INTEGER, SAVE:: lmt_pas ! number of time steps of "physics" per day      INTEGER, SAVE:: lmt_pas ! number of time steps of "physics" per day
248      REAL, save:: pctsrf(klon, nbsrf) ! percentage of surface      REAL, save:: pctsrf(klon, nbsrf) ! percentage of surface
249      REAL pctsrf_new(klon, nbsrf) ! pourcentage surfaces issus d'ORCHIDEE      REAL pctsrf_new(klon, nbsrf) ! pourcentage surfaces issus d'ORCHIDEE
250      REAL, save:: albsol(klon) ! albedo du sol total      REAL, save:: albsol(klon) ! albedo du sol total visible
     REAL, save:: albsollw(klon) ! albedo du sol total  
251      REAL, SAVE:: wo(klon, llm) ! column density of ozone in a cell, in kDU      REAL, SAVE:: wo(klon, llm) ! column density of ozone in a cell, in kDU
252    
     ! Declaration des procedures appelees  
   
     EXTERNAL nuage ! calculer les proprietes radiatives  
     EXTERNAL transp ! transport total de l'eau et de l'energie  
   
     ! Variables locales  
   
253      real, save:: clwcon(klon, llm), rnebcon(klon, llm)      real, save:: clwcon(klon, llm), rnebcon(klon, llm)
254      real, save:: clwcon0(klon, llm), rnebcon0(klon, llm)      real, save:: clwcon0(klon, llm), rnebcon0(klon, llm)
255    
# Line 368  contains Line 274  contains
274      ! Le rayonnement n'est pas calcul\'e tous les pas, il faut donc que      ! Le rayonnement n'est pas calcul\'e tous les pas, il faut donc que
275      ! les variables soient r\'emanentes.      ! les variables soient r\'emanentes.
276      REAL, save:: heat(klon, llm) ! chauffage solaire      REAL, save:: heat(klon, llm) ! chauffage solaire
277      REAL heat0(klon, llm) ! chauffage solaire ciel clair      REAL, save:: heat0(klon, llm) ! chauffage solaire ciel clair
278      REAL, save:: cool(klon, llm) ! refroidissement infrarouge      REAL, save:: cool(klon, llm) ! refroidissement infrarouge
279      REAL cool0(klon, llm) ! refroidissement infrarouge ciel clair      REAL, save:: cool0(klon, llm) ! refroidissement infrarouge ciel clair
280      REAL, save:: topsw(klon), toplw(klon), solsw(klon)      REAL, save:: topsw(klon), toplw(klon), solsw(klon)
281      REAL, save:: sollw(klon) ! rayonnement infrarouge montant \`a la surface      REAL, save:: sollw(klon) ! rayonnement infrarouge montant \`a la surface
282      real, save:: sollwdown(klon) ! downward LW flux at surface      real, save:: sollwdown(klon) ! downward LW flux at surface
283      REAL, save:: topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)      REAL, save:: topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)
284      REAL albpla(klon)      REAL, save:: albpla(klon)
285      REAL fsollw(klon, nbsrf) ! bilan flux IR pour chaque sous surface      REAL fsollw(klon, nbsrf) ! bilan flux IR pour chaque sous-surface
286      REAL fsolsw(klon, nbsrf) ! flux solaire absorb. pour chaque sous surface      REAL fsolsw(klon, nbsrf) ! flux solaire absorb\'e pour chaque sous-surface
     SAVE albpla  
     SAVE heat0, cool0  
   
     INTEGER itaprad  
     SAVE itaprad  
287    
288      REAL conv_q(klon, llm) ! convergence de l'humidite (kg/kg/s)      REAL conv_q(klon, llm) ! convergence de l'humidite (kg/kg/s)
289      REAL conv_t(klon, llm) ! convergence of temperature (K/s)      REAL conv_t(klon, llm) ! convergence of temperature (K/s)
290    
291      REAL cldl(klon), cldm(klon), cldh(klon) !nuages bas, moyen et haut      REAL cldl(klon), cldm(klon), cldh(klon) ! nuages bas, moyen et haut
292      REAL cldt(klon), cldq(klon) !nuage total, eau liquide integree      REAL cldt(klon), cldq(klon) ! nuage total, eau liquide integree
293    
294      REAL zxtsol(klon), zxqsurf(klon), zxsnow(klon), zxfluxlat(klon)      REAL zxtsol(klon), zxqsurf(klon), zxsnow(klon), zxfluxlat(klon)
295    
296      REAL dist, rmu0(klon), fract(klon)      REAL dist, mu0(klon), fract(klon)
297      real zlongi      real longi
298      REAL z_avant(klon), z_apres(klon), z_factor(klon)      REAL z_avant(klon), z_apres(klon), z_factor(klon)
299      REAL za, zb      REAL za, zb
300      REAL zx_t, zx_qs, zdelta, zcor      REAL zx_t, zx_qs, zcor
301      real zqsat(klon, llm)      real zqsat(klon, llm)
302      INTEGER i, k, iq, nsrf      INTEGER i, k, iq, nsrf
303      REAL, PARAMETER:: t_coup = 234.      REAL, PARAMETER:: t_coup = 234.
304      REAL zphi(klon, llm)      REAL zphi(klon, llm)
305    
306      ! cf. AM Variables locales pour la CLA (hbtm2)      ! cf. Anne Mathieu variables pour la couche limite atmosphérique (hbtm)
307    
308      REAL, SAVE:: pblh(klon, nbsrf) ! Hauteur de couche limite      REAL, SAVE:: pblh(klon, nbsrf) ! Hauteur de couche limite
309      REAL, SAVE:: plcl(klon, nbsrf) ! Niveau de condensation de la CLA      REAL, SAVE:: plcl(klon, nbsrf) ! Niveau de condensation de la CLA
# Line 412  contains Line 313  contains
313      REAL, SAVE:: pblt(klon, nbsrf) ! T a la Hauteur de couche limite      REAL, SAVE:: pblt(klon, nbsrf) ! T a la Hauteur de couche limite
314      REAL, SAVE:: therm(klon, nbsrf)      REAL, SAVE:: therm(klon, nbsrf)
315      REAL, SAVE:: trmb1(klon, nbsrf) ! deep_cape      REAL, SAVE:: trmb1(klon, nbsrf) ! deep_cape
316      REAL, SAVE:: trmb2(klon, nbsrf) ! inhibition      REAL, SAVE:: trmb2(klon, nbsrf) ! inhibition
317      REAL, SAVE:: trmb3(klon, nbsrf) ! Point Omega      REAL, SAVE:: trmb3(klon, nbsrf) ! Point Omega
318      ! Grdeurs de sorties      ! Grandeurs de sorties
319      REAL s_pblh(klon), s_lcl(klon), s_capCL(klon)      REAL s_pblh(klon), s_lcl(klon), s_capCL(klon)
320      REAL s_oliqCL(klon), s_cteiCL(klon), s_pblt(klon)      REAL s_oliqCL(klon), s_cteiCL(klon), s_pblt(klon)
321      REAL s_therm(klon), s_trmb1(klon), s_trmb2(klon)      REAL s_therm(klon), s_trmb1(klon), s_trmb2(klon)
322      REAL s_trmb3(klon)      REAL s_trmb3(klon)
323    
324      ! Variables locales pour la convection de K. Emanuel :      ! Variables pour la convection de K. Emanuel :
325    
326      REAL upwd(klon, llm) ! saturated updraft mass flux      REAL upwd(klon, llm) ! saturated updraft mass flux
327      REAL dnwd(klon, llm) ! saturated downdraft mass flux      REAL dnwd(klon, llm) ! saturated downdraft mass flux
# Line 452  contains Line 353  contains
353      REAL prfl(klon, llm + 1), psfl(klon, llm + 1)      REAL prfl(klon, llm + 1), psfl(klon, llm + 1)
354    
355      INTEGER, save:: ibas_con(klon), itop_con(klon)      INTEGER, save:: ibas_con(klon), itop_con(klon)
356        real ema_pct(klon) ! Emanuel pressure at cloud top, in Pa
357    
358      REAL rain_con(klon), rain_lsc(klon)      REAL rain_con(klon), rain_lsc(klon)
359      REAL snow_con(klon), snow_lsc(klon)      REAL, save:: snow_con(klon) ! neige (mm / s)
360        real snow_lsc(klon)
361      REAL d_ts(klon, nbsrf)      REAL d_ts(klon, nbsrf)
362    
363      REAL d_u_vdf(klon, llm), d_v_vdf(klon, llm)      REAL d_u_vdf(klon, llm), d_v_vdf(klon, llm)
# Line 478  contains Line 381  contains
381      integer:: iflag_cldcon = 1      integer:: iflag_cldcon = 1
382      logical ptconv(klon, llm)      logical ptconv(klon, llm)
383    
384      ! Variables locales pour effectuer les appels en s\'erie :      ! Variables pour effectuer les appels en s\'erie :
385    
386      REAL t_seri(klon, llm), q_seri(klon, llm)      REAL t_seri(klon, llm), q_seri(klon, llm)
387      REAL ql_seri(klon, llm)      REAL ql_seri(klon, llm)
# Line 492  contains Line 395  contains
395      REAL zustrph(klon), zvstrph(klon)      REAL zustrph(klon), zvstrph(klon)
396      REAL aam, torsfc      REAL aam, torsfc
397    
     REAL zx_tmp_fi2d(klon) ! variable temporaire grille physique  
   
     INTEGER, SAVE:: nid_ins  
   
398      REAL ve_lay(klon, llm) ! transport meri. de l'energie a chaque niveau vert.      REAL ve_lay(klon, llm) ! transport meri. de l'energie a chaque niveau vert.
399      REAL vq_lay(klon, llm) ! transport meri. de l'eau a chaque niveau vert.      REAL vq_lay(klon, llm) ! transport meri. de l'eau a chaque niveau vert.
400      REAL ue_lay(klon, llm) ! transport zonal de l'energie a chaque niveau vert.      REAL ue_lay(klon, llm) ! transport zonal de l'energie a chaque niveau vert.
401      REAL uq_lay(klon, llm) ! transport zonal de l'eau a chaque niveau vert.      REAL uq_lay(klon, llm) ! transport zonal de l'eau a chaque niveau vert.
402    
     REAL zsto  
403      real date0      real date0
404    
405      ! Variables li\'ees au bilan d'\'energie et d'enthalpie :      ! Variables li\'ees au bilan d'\'energie et d'enthalpie :
# Line 511  contains Line 409  contains
409      REAL zero_v(klon)      REAL zero_v(klon)
410      CHARACTER(LEN = 20) tit      CHARACTER(LEN = 20) tit
411      INTEGER:: ip_ebil = 0 ! print level for energy conservation diagnostics      INTEGER:: ip_ebil = 0 ! print level for energy conservation diagnostics
412      INTEGER:: if_ebil = 0 ! verbosity for diagnostics of energy conservation      INTEGER:: if_ebil = 0 ! verbosity for diagnostics of energy conservation
413    
414      REAL d_t_ec(klon, llm) ! tendance due \`a la conversion Ec -> E thermique      REAL d_t_ec(klon, llm) ! tendance due \`a la conversion Ec -> E thermique
415      REAL ZRCPD      REAL ZRCPD
# Line 526  contains Line 424  contains
424      REAL sulfate(klon, llm) ! SO4 aerosol concentration (micro g/m3)      REAL sulfate(klon, llm) ! SO4 aerosol concentration (micro g/m3)
425    
426      REAL, save:: sulfate_pi(klon, llm)      REAL, save:: sulfate_pi(klon, llm)
427      ! SO4 aerosol concentration, in micro g/m3, pre-industrial value      ! SO4 aerosol concentration, in \mu g/m3, pre-industrial value
428    
429      REAL cldtaupi(klon, llm)      REAL cldtaupi(klon, llm)
430      ! cloud optical thickness for pre-industrial (pi) aerosols      ! cloud optical thickness for pre-industrial aerosols
431    
432      REAL re(klon, llm) ! Cloud droplet effective radius      REAL re(klon, llm) ! Cloud droplet effective radius
433      REAL fl(klon, llm) ! denominator of re      REAL fl(klon, llm) ! denominator of re
# Line 541  contains Line 439  contains
439      REAL topswad(klon), solswad(klon) ! aerosol direct effect      REAL topswad(klon), solswad(klon) ! aerosol direct effect
440      REAL topswai(klon), solswai(klon) ! aerosol indirect effect      REAL topswai(klon), solswai(klon) ! aerosol indirect effect
441    
     REAL aerindex(klon) ! POLDER aerosol index  
   
442      LOGICAL:: ok_ade = .false. ! apply aerosol direct effect      LOGICAL:: ok_ade = .false. ! apply aerosol direct effect
443      LOGICAL:: ok_aie = .false. ! apply aerosol indirect effect      LOGICAL:: ok_aie = .false. ! apply aerosol indirect effect
444    
# Line 558  contains Line 454  contains
454      SAVE ffonte      SAVE ffonte
455      SAVE fqcalving      SAVE fqcalving
456      SAVE rain_con      SAVE rain_con
     SAVE snow_con  
457      SAVE topswai      SAVE topswai
458      SAVE topswad      SAVE topswad
459      SAVE solswai      SAVE solswai
# Line 566  contains Line 461  contains
461      SAVE d_u_con      SAVE d_u_con
462      SAVE d_v_con      SAVE d_v_con
463    
464      real zmasse(klon, llm)      real zmasse(klon, llm)
465      ! (column-density of mass of air in a cell, in kg m-2)      ! (column-density of mass of air in a cell, in kg m-2)
466    
467      real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2      integer, save:: ncid_startphy
468    
469      namelist /physiq_nml/ ok_journe, ok_mensuel, ok_instan, fact_cldcon, &      namelist /physiq_nml/ fact_cldcon, facttemps, ok_newmicro, &
470           facttemps, ok_newmicro, iflag_cldcon, ratqsbas, ratqshaut, if_ebil, &           iflag_cldcon, ratqsbas, ratqshaut, if_ebil, ok_ade, ok_aie, bl95_b0, &
471           ok_ade, ok_aie, bl95_b0, bl95_b1, iflag_thermals, nsplit_thermals           bl95_b1, iflag_thermals, nsplit_thermals
472    
473      !----------------------------------------------------------------      !----------------------------------------------------------------
474    
475      IF (if_ebil >= 1) zero_v = 0.      IF (if_ebil >= 1) zero_v = 0.
476      IF (nqmx < 2) CALL abort_gcm('physiq', &      IF (nqmx < 2) CALL abort_gcm('physiq', &
477           'eaux vapeur et liquide sont indispensables', 1)           'eaux vapeur et liquide sont indispensables')
478    
479      test_firstcal: IF (firstcal) THEN      test_firstcal: IF (firstcal) THEN
480         ! initialiser         ! initialiser
# Line 614  contains Line 509  contains
509         pblt =0. ! T a la Hauteur de couche limite         pblt =0. ! T a la Hauteur de couche limite
510         therm =0.         therm =0.
511         trmb1 =0. ! deep_cape         trmb1 =0. ! deep_cape
512         trmb2 =0. ! inhibition         trmb2 =0. ! inhibition
513         trmb3 =0. ! Point Omega         trmb3 =0. ! Point Omega
514    
515         IF (if_ebil >= 1) d_h_vcol_phy = 0.         IF (if_ebil >= 1) d_h_vcol_phy = 0.
# Line 630  contains Line 525  contains
525         ! Initialiser les compteurs:         ! Initialiser les compteurs:
526    
527         frugs = 0.         frugs = 0.
528         itap = 0         CALL phyetat0(pctsrf, ftsol, ftsoil, fqsurf, qsol, fsnow, falbe, &
529         itaprad = 0              fevap, rain_fall, snow_fall, solsw, sollw, dlw, radsol, frugs, &
530         CALL phyetat0(pctsrf, ftsol, ftsoil, tslab, seaice, fqsurf, qsol, &              agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, t_ancien, &
531              fsnow, falbe, falblw, fevap, rain_fall, snow_fall, solsw, sollw, &              q_ancien, ancien_ok, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, &
532              dlw, radsol, frugs, agesno, zmea, zstd, zsig, zgam, zthe, zpic, &              w01, ncid_startphy)
             zval, t_ancien, q_ancien, ancien_ok, rnebcon, ratqs, clwcon, &  
             run_off_lic_0, sig1, w01)  
533    
534         ! ATTENTION : il faudra a terme relire q2 dans l'etat initial         ! ATTENTION : il faudra a terme relire q2 dans l'etat initial
535         q2 = 1e-8         q2 = 1e-8
536    
537         radpas = NINT(86400. / dtphys / nbapp_rad)         lmt_pas = day_step / iphysiq
538           print *, 'Number of time steps of "physics" per day: ', lmt_pas
        ! on remet le calendrier a zero  
        IF (raz_date) itau_phy = 0  
   
        PRINT *, 'cycle_diurne = ', cycle_diurne  
        CALL printflag(radpas, ok_journe, ok_instan, ok_region)  
539    
540         IF (dtphys * REAL(radpas) > 21600. .AND. cycle_diurne) THEN         radpas = lmt_pas / nbapp_rad
541            print *, "Au minimum 4 appels par jour si cycle diurne"         print *, "radpas = ", radpas
           call abort_gcm('physiq', &  
                "Nombre d'appels au rayonnement insuffisant", 1)  
        ENDIF  
542    
543         ! Initialisation pour le sch\'ema de convection d'Emanuel :         ! Initialisation pour le sch\'ema de convection d'Emanuel :
544         IF (iflag_con >= 3) THEN         IF (conv_emanuel) THEN
545            ibas_con = 1            ibas_con = 1
546            itop_con = 1            itop_con = 1
547         ENDIF         ENDIF
# Line 668  contains Line 553  contains
553            rugoro = 0.            rugoro = 0.
554         ENDIF         ENDIF
555    
        lmt_pas = NINT(86400. / dtphys) ! tous les jours  
        print *, 'Number of time steps of "physics" per day: ', lmt_pas  
   
556         ecrit_ins = NINT(ecrit_ins/dtphys)         ecrit_ins = NINT(ecrit_ins/dtphys)
557         ecrit_hf = NINT(ecrit_hf/dtphys)         ecrit_hf = NINT(ecrit_hf/dtphys)
558         ecrit_mth = NINT(ecrit_mth/dtphys)         ecrit_mth = NINT(ecrit_mth/dtphys)
# Line 679  contains Line 561  contains
561    
562         ! Initialisation des sorties         ! Initialisation des sorties
563    
564         call ini_histins(dtphys, ok_instan, nid_ins)         call ini_histins(dtphys)
565         CALL ymds2ju(annee_ref, 1, int(day_ref), 0., date0)         CALL ymds2ju(annee_ref, 1, day_ref, 0., date0)
566         ! Positionner date0 pour initialisation de ORCHIDEE         ! Positionner date0 pour initialisation de ORCHIDEE
567         print *, 'physiq date0: ', date0         print *, 'physiq date0: ', date0
568           CALL phyredem0(lmt_pas)
569      ENDIF test_firstcal      ENDIF test_firstcal
570    
571      ! We will modify variables *_seri and we will not touch variables      ! We will modify variables *_seri and we will not touch variables
# Line 692  contains Line 575  contains
575      v_seri = v      v_seri = v
576      q_seri = qx(:, :, ivap)      q_seri = qx(:, :, ivap)
577      ql_seri = qx(:, :, iliq)      ql_seri = qx(:, :, iliq)
578      tr_seri = qx(:, :, 3: nqmx)      tr_seri = qx(:, :, 3:nqmx)
579    
580      ztsol = sum(ftsol * pctsrf, dim = 2)      ztsol = sum(ftsol * pctsrf, dim = 2)
581    
582      IF (if_ebil >= 1) THEN      IF (if_ebil >= 1) THEN
583         tit = 'after dynamics'         tit = 'after dynamics'
584         CALL diagetpq(airephy, tit, ip_ebil, 1, 1, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 1, 1, dtphys, t_seri, q_seri, &
585              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
586         ! Comme les tendances de la physique sont ajout\'es dans la         ! Comme les tendances de la physique sont ajout\'es dans la
587         !  dynamique, la variation d'enthalpie par la dynamique devrait         ! dynamique, la variation d'enthalpie par la dynamique devrait
588         !  \^etre \'egale \`a la variation de la physique au pas de temps         ! \^etre \'egale \`a la variation de la physique au pas de temps
589         !  pr\'ec\'edent.  Donc la somme de ces 2 variations devrait \^etre         ! pr\'ec\'edent. Donc la somme de ces 2 variations devrait \^etre
590         !  nulle.         ! nulle.
591         call diagphy(airephy, tit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &         call diagphy(airephy, tit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
592              zero_v, zero_v, zero_v, zero_v, ztsol, d_h_vcol + d_h_vcol_phy, &              zero_v, zero_v, zero_v, zero_v, ztsol, d_h_vcol + d_h_vcol_phy, &
593              d_qt, 0.)              d_qt, 0.)
# Line 738  contains Line 621  contains
621      ! Check temperatures:      ! Check temperatures:
622      CALL hgardfou(t_seri, ftsol)      CALL hgardfou(t_seri, ftsol)
623    
624      ! Incrémenter le compteur de la physique      call increment_itap
625      itap = itap + 1      julien = MOD(dayvrai, 360)
     julien = MOD(NINT(rdayvrai), 360)  
626      if (julien == 0) julien = 360      if (julien == 0) julien = 360
627    
628      forall (k = 1: llm) zmasse(:, k) = (paprs(:, k)-paprs(:, k + 1)) / rg      forall (k = 1: llm) zmasse(:, k) = (paprs(:, k) - paprs(:, k + 1)) / rg
629    
630      ! Prescrire l'ozone :      ! Prescrire l'ozone :
631      wo = ozonecm(REAL(julien), paprs)      wo = ozonecm(REAL(julien), paprs)
# Line 759  contains Line 641  contains
641      ENDDO      ENDDO
642      ql_seri = 0.      ql_seri = 0.
643    
644      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
645         tit = 'after reevap'         tit = 'after reevap'
646         CALL diagetpq(airephy, tit, ip_ebil, 2, 1, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 2, 1, dtphys, t_seri, q_seri, &
647              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
# Line 770  contains Line 652  contains
652      frugs = MAX(frugs, 0.000015)      frugs = MAX(frugs, 0.000015)
653      zxrugs = sum(frugs * pctsrf, dim = 2)      zxrugs = sum(frugs * pctsrf, dim = 2)
654    
655      ! Calculs nécessaires au calcul de l'albedo dans l'interface      ! Calculs n\'ecessaires au calcul de l'albedo dans l'interface avec
656        ! la surface.
657    
658      CALL orbite(REAL(julien), zlongi, dist)      CALL orbite(REAL(julien), longi, dist)
659      IF (cycle_diurne) THEN      IF (cycle_diurne) THEN
660         CALL zenang(zlongi, time, dtphys * REAL(radpas), rmu0, fract)         CALL zenang(longi, time, dtphys * radpas, mu0, fract)
661      ELSE      ELSE
662         rmu0 = -999.999         mu0 = - 999.999
663      ENDIF      ENDIF
664    
665      ! Calcul de l'abedo moyen par maille      ! Calcul de l'abedo moyen par maille
666      albsol = sum(falbe * pctsrf, dim = 2)      albsol = sum(falbe * pctsrf, dim = 2)
     albsollw = sum(falblw * pctsrf, dim = 2)  
667    
668      ! R\'epartition sous maille des flux longwave et shortwave      ! R\'epartition sous maille des flux longwave et shortwave
669      ! R\'epartition du longwave par sous-surface lin\'earis\'ee      ! R\'epartition du longwave par sous-surface lin\'earis\'ee
# Line 796  contains Line 678  contains
678    
679      ! Couche limite:      ! Couche limite:
680    
681      CALL clmain(dtphys, itap, pctsrf, pctsrf_new, t_seri, q_seri, u_seri, &      CALL clmain(dtphys, pctsrf, pctsrf_new, t_seri, q_seri, u_seri, v_seri, &
682           v_seri, julien, rmu0, co2_ppm, ftsol, cdmmax, cdhmax, &           julien, mu0, ftsol, cdmmax, cdhmax, ksta, ksta_ter, ok_kzmin, &
683           ksta, ksta_ter, ok_kzmin, ftsoil, qsol, paprs, play, fsnow, fqsurf, &           ftsoil, qsol, paprs, play, fsnow, fqsurf, fevap, falbe, fluxlat, &
684           fevap, falbe, falblw, fluxlat, rain_fall, snow_fall, fsolsw, fsollw, &           rain_fall, snow_fall, fsolsw, fsollw, fder, rlat, frugs, firstcal, &
685           fder, rlat, frugs, firstcal, agesno, rugoro, d_t_vdf, d_q_vdf, &           agesno, rugoro, d_t_vdf, d_q_vdf, d_u_vdf, d_v_vdf, d_ts, fluxt, &
686           d_u_vdf, d_v_vdf, d_ts, fluxt, fluxq, fluxu, fluxv, cdragh, cdragm, &           fluxq, fluxu, fluxv, cdragh, cdragm, q2, dsens, devap, ycoefh, yu1, &
687           q2, dsens, devap, ycoefh, yu1, yv1, t2m, q2m, u10m, v10m, pblh, &           yv1, t2m, q2m, u10m, v10m, pblh, capCL, oliqCL, cteiCL, pblT, therm, &
688           capCL, oliqCL, cteiCL, pblT, therm, trmb1, trmb2, trmb3, plcl, &           trmb1, trmb2, trmb3, plcl, fqcalving, ffonte, run_off_lic_0)
          fqcalving, ffonte, run_off_lic_0, fluxo, fluxg, tslab)  
689    
690      ! Incr\'ementation des flux      ! Incr\'ementation des flux
691    
# Line 837  contains Line 718  contains
718         ENDDO         ENDDO
719      ENDDO      ENDDO
720    
721      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
722         tit = 'after clmain'         tit = 'after clmain'
723         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
724              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
# Line 858  contains Line 739  contains
739         zxffonte(i) = 0.         zxffonte(i) = 0.
740         zxfqcalving(i) = 0.         zxfqcalving(i) = 0.
741    
742         s_pblh(i) = 0.         s_pblh(i) = 0.
743         s_lcl(i) = 0.         s_lcl(i) = 0.
744         s_capCL(i) = 0.         s_capCL(i) = 0.
745         s_oliqCL(i) = 0.         s_oliqCL(i) = 0.
746         s_cteiCL(i) = 0.         s_cteiCL(i) = 0.
# Line 868  contains Line 749  contains
749         s_trmb1(i) = 0.         s_trmb1(i) = 0.
750         s_trmb2(i) = 0.         s_trmb2(i) = 0.
751         s_trmb3(i) = 0.         s_trmb3(i) = 0.
   
        IF (abs(pctsrf(i, is_ter) + pctsrf(i, is_lic) + pctsrf(i, is_oce) &  
             + pctsrf(i, is_sic) - 1.)  >  EPSFRA) print *, &  
             'physiq : probl\`eme sous surface au point ', i, &  
             pctsrf(i, 1 : nbsrf)  
752      ENDDO      ENDDO
753    
754        call assert(abs(sum(pctsrf, dim = 2) - 1.) <= EPSFRA, 'physiq: pctsrf')
755    
756      DO nsrf = 1, nbsrf      DO nsrf = 1, nbsrf
757         DO i = 1, klon         DO i = 1, klon
758            ftsol(i, nsrf) = ftsol(i, nsrf) + d_ts(i, nsrf)            ftsol(i, nsrf) = ftsol(i, nsrf) + d_ts(i, nsrf)
# Line 928  contains Line 807  contains
807      ! Calculer la dérive du flux infrarouge      ! Calculer la dérive du flux infrarouge
808    
809      DO i = 1, klon      DO i = 1, klon
810         dlw(i) = - 4. * RSIGMA * zxtsol(i)**3         dlw(i) = - 4. * RSIGMA * zxtsol(i)**3
     ENDDO  
   
     ! Appeler la convection (au choix)  
   
     DO k = 1, llm  
        DO i = 1, klon  
           conv_q(i, k) = d_q_dyn(i, k) + d_q_vdf(i, k) / dtphys  
           conv_t(i, k) = d_t_dyn(i, k) + d_t_vdf(i, k) / dtphys  
        ENDDO  
811      ENDDO      ENDDO
812    
813      IF (check) print *, "avantcon = ", qcheck(paprs, q_seri, ql_seri)      IF (check) print *, "avantcon = ", qcheck(paprs, q_seri, ql_seri)
814    
815      if (iflag_con == 2) then      ! Appeler la convection
        z_avant = sum((q_seri + ql_seri) * zmasse, dim=2)  
        CALL conflx(dtphys, paprs, play, t_seri(:, llm:1:-1), &  
             q_seri(:, llm:1:-1), conv_t, conv_q, zxfluxq(:, 1), omega, &  
             d_t_con, d_q_con, rain_con, snow_con, mfu(:, llm:1:-1), &  
             mfd(:, llm:1:-1), pen_u, pde_u, pen_d, pde_d, kcbot, kctop, &  
             kdtop, pmflxr, pmflxs)  
        WHERE (rain_con < 0.) rain_con = 0.  
        WHERE (snow_con < 0.) snow_con = 0.  
        ibas_con = llm + 1 - kcbot  
        itop_con = llm + 1 - kctop  
     else  
        ! iflag_con >= 3  
816    
817        if (conv_emanuel) then
818         da = 0.         da = 0.
819         mp = 0.         mp = 0.
820         phi = 0.         phi = 0.
821         CALL concvl(dtphys, paprs, play, t_seri, q_seri, u_seri, v_seri, sig1, &         CALL concvl(paprs, play, t_seri, q_seri, u_seri, v_seri, sig1, w01, &
822              w01, d_t_con, d_q_con, d_u_con, d_v_con, rain_con, snow_con, &              d_t_con, d_q_con, d_u_con, d_v_con, rain_con, ibas_con, itop_con, &
823              ibas_con, itop_con, upwd, dnwd, dnwd0, Ma, cape, iflagctrl, &              upwd, dnwd, dnwd0, Ma, cape, iflagctrl, qcondc, pmflxr, da, phi, mp)
824              qcondc, wd, pmflxr, pmflxs, da, phi, mp)         snow_con = 0.
825         clwcon0 = qcondc         clwcon0 = qcondc
826         mfu = upwd + dnwd         mfu = upwd + dnwd
        IF (.NOT. ok_gust) wd = 0.  
   
        ! Calcul des propri\'et\'es des nuages convectifs  
827    
828         DO k = 1, llm         IF (thermcep) THEN
829            DO i = 1, klon            zqsat = MIN(0.5, r2es * FOEEW(t_seri, rtt >= t_seri) / play)
830               IF (thermcep) THEN            zqsat = zqsat / (1. - retv * zqsat)
831                  zdelta = MAX(0., SIGN(1., rtt - t_seri(i, k)))         ELSE
832                  zqsat(i, k) = r2es * FOEEW(t_seri(i, k), zdelta) / play(i, k)            zqsat = merge(qsats(t_seri), qsatl(t_seri), t_seri < t_coup) / play
833                  zqsat(i, k) = MIN(0.5, zqsat(i, k))         ENDIF
                 zqsat(i, k) = zqsat(i, k) / (1.-retv*zqsat(i, k))  
              ELSE  
                 IF (t_seri(i, k) < t_coup) THEN  
                    zqsat(i, k) = qsats(t_seri(i, k))/play(i, k)  
                 ELSE  
                    zqsat(i, k) = qsatl(t_seri(i, k))/play(i, k)  
                 ENDIF  
              ENDIF  
           ENDDO  
        ENDDO  
834    
835         ! calcul des proprietes des nuages convectifs         ! Properties of convective clouds
836         clwcon0 = fact_cldcon * clwcon0         clwcon0 = fact_cldcon * clwcon0
837         call clouds_gno(klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, &         call clouds_gno(klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, &
838              rnebcon0)              rnebcon0)
839    
840           forall (i = 1:klon) ema_pct(i) = paprs(i, itop_con(i) + 1)
841         mfd = 0.         mfd = 0.
842         pen_u = 0.         pen_u = 0.
843         pen_d = 0.         pen_d = 0.
844         pde_d = 0.         pde_d = 0.
845         pde_u = 0.         pde_u = 0.
846        else
847           conv_q = d_q_dyn + d_q_vdf / dtphys
848           conv_t = d_t_dyn + d_t_vdf / dtphys
849           z_avant = sum((q_seri + ql_seri) * zmasse, dim=2)
850           CALL conflx(dtphys, paprs, play, t_seri(:, llm:1:- 1), &
851                q_seri(:, llm:1:- 1), conv_t, conv_q, zxfluxq(:, 1), omega, &
852                d_t_con, d_q_con, rain_con, snow_con, mfu(:, llm:1:- 1), &
853                mfd(:, llm:1:- 1), pen_u, pde_u, pen_d, pde_d, kcbot, kctop, &
854                kdtop, pmflxr, pmflxs)
855           WHERE (rain_con < 0.) rain_con = 0.
856           WHERE (snow_con < 0.) snow_con = 0.
857           ibas_con = llm + 1 - kcbot
858           itop_con = llm + 1 - kctop
859      END if      END if
860    
861      DO k = 1, llm      DO k = 1, llm
# Line 1007  contains Line 867  contains
867         ENDDO         ENDDO
868      ENDDO      ENDDO
869    
870      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
871         tit = 'after convect'         tit = 'after convect'
872         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
873              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
# Line 1029  contains Line 889  contains
889         print *, "Precip = ", zx_t         print *, "Precip = ", zx_t
890      ENDIF      ENDIF
891    
892      IF (iflag_con == 2) THEN      IF (.not. conv_emanuel) THEN
893         z_apres = sum((q_seri + ql_seri) * zmasse, dim=2)         z_apres = sum((q_seri + ql_seri) * zmasse, dim=2)
894         z_factor = (z_avant - (rain_con + snow_con) * dtphys) / z_apres         z_factor = (z_avant - (rain_con + snow_con) * dtphys) / z_apres
895         DO k = 1, llm         DO k = 1, llm
# Line 1061  contains Line 921  contains
921              q_seri, d_u_ajs, d_v_ajs, d_t_ajs, d_q_ajs, fm_therm, entr_therm)              q_seri, d_u_ajs, d_v_ajs, d_t_ajs, d_q_ajs, fm_therm, entr_therm)
922      endif      endif
923    
924      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
925         tit = 'after dry_adjust'         tit = 'after dry_adjust'
926         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
927              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
# Line 1088  contains Line 948  contains
948      do k = 1, llm      do k = 1, llm
949         do i = 1, klon         do i = 1, klon
950            ratqss(i, k) = ratqsbas + (ratqshaut - ratqsbas) &            ratqss(i, k) = ratqsbas + (ratqshaut - ratqsbas) &
951                 * min((paprs(i, 1) - play(i, k)) / (paprs(i, 1) - 3e4), 1.)                 * min((paprs(i, 1) - play(i, k)) / (paprs(i, 1) - 3e4), 1.)
952         enddo         enddo
953      enddo      enddo
954    
# Line 1135  contains Line 995  contains
995         print *, "Precip = ", zx_t         print *, "Precip = ", zx_t
996      ENDIF      ENDIF
997    
998      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
999         tit = 'after fisrt'         tit = 'after fisrt'
1000         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1001              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
# Line 1147  contains Line 1007  contains
1007    
1008      ! 1. NUAGES CONVECTIFS      ! 1. NUAGES CONVECTIFS
1009    
1010      IF (iflag_cldcon <= -1) THEN      IF (iflag_cldcon <= - 1) THEN
1011         ! seulement pour Tiedtke         ! seulement pour Tiedtke
1012         snow_tiedtke = 0.         snow_tiedtke = 0.
1013         if (iflag_cldcon == -1) then         if (iflag_cldcon == - 1) then
1014            rain_tiedtke = rain_con            rain_tiedtke = rain_con
1015         else         else
1016            rain_tiedtke = 0.            rain_tiedtke = 0.
1017            do k = 1, llm            do k = 1, llm
1018               do i = 1, klon               do i = 1, klon
1019                  if (d_q_con(i, k) < 0.) then                  if (d_q_con(i, k) < 0.) then
1020                     rain_tiedtke(i) = rain_tiedtke(i)-d_q_con(i, k)/dtphys &                     rain_tiedtke(i) = rain_tiedtke(i) - d_q_con(i, k)/dtphys &
1021                          *zmasse(i, k)                          *zmasse(i, k)
1022                  endif                  endif
1023               enddo               enddo
# Line 1225  contains Line 1085  contains
1085         DO i = 1, klon         DO i = 1, klon
1086            zx_t = t_seri(i, k)            zx_t = t_seri(i, k)
1087            IF (thermcep) THEN            IF (thermcep) THEN
1088               zdelta = MAX(0., SIGN(1., rtt-zx_t))               zx_qs = r2es * FOEEW(zx_t, rtt >= zx_t)/play(i, k)
              zx_qs = r2es * FOEEW(zx_t, zdelta)/play(i, k)  
1089               zx_qs = MIN(0.5, zx_qs)               zx_qs = MIN(0.5, zx_qs)
1090               zcor = 1./(1.-retv*zx_qs)               zcor = 1./(1. - retv*zx_qs)
1091               zx_qs = zx_qs*zcor               zx_qs = zx_qs*zcor
1092            ELSE            ELSE
1093               IF (zx_t < t_coup) THEN               IF (zx_t < t_coup) THEN
# Line 1243  contains Line 1102  contains
1102      ENDDO      ENDDO
1103    
1104      ! Introduce the aerosol direct and first indirect radiative forcings:      ! Introduce the aerosol direct and first indirect radiative forcings:
1105      IF (ok_ade .OR. ok_aie) THEN      tau_ae = 0.
1106         ! Get sulfate aerosol distribution :      piz_ae = 0.
1107         CALL readsulfate(rdayvrai, firstcal, sulfate)      cg_ae = 0.
        CALL readsulfate_preind(rdayvrai, firstcal, sulfate_pi)  
   
        CALL aeropt(play, paprs, t_seri, sulfate, rhcl, tau_ae, piz_ae, cg_ae, &  
             aerindex)  
     ELSE  
        tau_ae = 0.  
        piz_ae = 0.  
        cg_ae = 0.  
     ENDIF  
1108    
1109      ! Param\`etres optiques des nuages et quelques param\`etres pour      ! Param\`etres optiques des nuages et quelques param\`etres pour
1110      ! diagnostics :      ! diagnostics :
# Line 1268  contains Line 1118  contains
1118              bl95_b1, cldtaupi, re, fl)              bl95_b1, cldtaupi, re, fl)
1119      endif      endif
1120    
1121      ! Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.      IF (MOD(itap - 1, radpas) == 0) THEN
1122      IF (MOD(itaprad, radpas) == 0) THEN         ! Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.
1123         DO i = 1, klon         ! Calcul de l'abedo moyen par maille
1124            albsol(i) = falbe(i, is_oce) * pctsrf(i, is_oce) &         albsol = sum(falbe * pctsrf, dim = 2)
1125                 + falbe(i, is_lic) * pctsrf(i, is_lic) &  
                + falbe(i, is_ter) * pctsrf(i, is_ter) &  
                + falbe(i, is_sic) * pctsrf(i, is_sic)  
           albsollw(i) = falblw(i, is_oce) * pctsrf(i, is_oce) &  
                + falblw(i, is_lic) * pctsrf(i, is_lic) &  
                + falblw(i, is_ter) * pctsrf(i, is_ter) &  
                + falblw(i, is_sic) * pctsrf(i, is_sic)  
        ENDDO  
1126         ! Rayonnement (compatible Arpege-IFS) :         ! Rayonnement (compatible Arpege-IFS) :
1127         CALL radlwsw(dist, rmu0, fract, paprs, play, zxtsol, albsol, &         CALL radlwsw(dist, mu0, fract, paprs, play, zxtsol, albsol, t_seri, &
1128              albsollw, t_seri, q_seri, wo, cldfra, cldemi, cldtau, heat, &              q_seri, wo, cldfra, cldemi, cldtau, heat, heat0, cool, cool0, &
1129              heat0, cool, cool0, radsol, albpla, topsw, toplw, solsw, sollw, &              radsol, albpla, topsw, toplw, solsw, sollw, sollwdown, topsw0, &
1130              sollwdown, topsw0, toplw0, solsw0, sollw0, lwdn0, lwdn, lwup0, &              toplw0, solsw0, sollw0, lwdn0, lwdn, lwup0, lwup, swdn0, swdn, &
1131              lwup, swdn0, swdn, swup0, swup, ok_ade, ok_aie, tau_ae, piz_ae, &              swup0, swup, ok_ade, ok_aie, tau_ae, piz_ae, cg_ae, topswad, &
1132              cg_ae, topswad, solswad, cldtaupi, topswai, solswai)              solswad, cldtaupi, topswai, solswai)
        itaprad = 0  
1133      ENDIF      ENDIF
     itaprad = itaprad + 1  
1134    
1135      ! Ajouter la tendance des rayonnements (tous les pas)      ! Ajouter la tendance des rayonnements (tous les pas)
1136    
1137      DO k = 1, llm      DO k = 1, llm
1138         DO i = 1, klon         DO i = 1, klon
1139            t_seri(i, k) = t_seri(i, k) + (heat(i, k)-cool(i, k)) * dtphys/86400.            t_seri(i, k) = t_seri(i, k) + (heat(i, k) - cool(i, k)) * dtphys/86400.
1140         ENDDO         ENDDO
1141      ENDDO      ENDDO
1142    
1143      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
1144         tit = 'after rad'         tit = 'after rad'
1145         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1146              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
# Line 1328  contains Line 1169  contains
1169      ! Param\'etrisation de l'orographie \`a l'\'echelle sous-maille :      ! Param\'etrisation de l'orographie \`a l'\'echelle sous-maille :
1170    
1171      IF (ok_orodr) THEN      IF (ok_orodr) THEN
1172         ! selection des points pour lesquels le shema est actif:         ! S\'election des points pour lesquels le sch\'ema est actif :
1173         igwd = 0         igwd = 0
1174         DO i = 1, klon         DO i = 1, klon
1175            itest(i) = 0            itest(i) = 0
1176            IF (((zpic(i)-zmea(i)) > 100.).AND.(zstd(i) > 10.)) THEN            IF (zpic(i) - zmea(i) > 100. .AND. zstd(i) > 10.) THEN
1177               itest(i) = 1               itest(i) = 1
1178               igwd = igwd + 1               igwd = igwd + 1
              idx(igwd) = i  
1179            ENDIF            ENDIF
1180         ENDDO         ENDDO
1181    
1182         CALL drag_noro(klon, llm, dtphys, paprs, play, zmea, zstd, zsig, zgam, &         CALL drag_noro(klon, llm, dtphys, paprs, play, zmea, zstd, zsig, zgam, &
1183              zthe, zpic, zval, igwd, idx, itest, t_seri, u_seri, v_seri, &              zthe, zpic, zval, itest, t_seri, u_seri, v_seri, zulow, zvlow, &
1184              zulow, zvlow, zustrdr, zvstrdr, d_t_oro, d_u_oro, d_v_oro)              zustrdr, zvstrdr, d_t_oro, d_u_oro, d_v_oro)
1185    
1186         ! ajout des tendances         ! ajout des tendances
1187         DO k = 1, llm         DO k = 1, llm
# Line 1358  contains Line 1198  contains
1198         igwd = 0         igwd = 0
1199         DO i = 1, klon         DO i = 1, klon
1200            itest(i) = 0            itest(i) = 0
1201            IF ((zpic(i) - zmea(i)) > 100.) THEN            IF (zpic(i) - zmea(i) > 100.) THEN
1202               itest(i) = 1               itest(i) = 1
1203               igwd = igwd + 1               igwd = igwd + 1
              idx(igwd) = i  
1204            ENDIF            ENDIF
1205         ENDDO         ENDDO
1206    
# Line 1394  contains Line 1233  contains
1233         ENDDO         ENDDO
1234      ENDDO      ENDDO
1235    
1236      CALL aaam_bud(ra, rg, romega, rlat, rlon, pphis, zustrdr, zustrli, &      CALL aaam_bud(rg, romega, rlat, rlon, pphis, zustrdr, zustrli, zustrph, &
1237           zustrph, zvstrdr, zvstrli, zvstrph, paprs, u, v, aam, torsfc)           zvstrdr, zvstrli, zvstrph, paprs, u, v, aam, torsfc)
1238    
1239      IF (if_ebil >= 2) CALL diagetpq(airephy, 'after orography', ip_ebil, 2, &      IF (if_ebil >= 2) CALL diagetpq(airephy, 'after orography', ip_ebil, 2, &
1240           2, dtphys, t_seri, q_seri, ql_seri, u_seri, v_seri, paprs, d_h_vcol, &           2, dtphys, t_seri, q_seri, ql_seri, u_seri, v_seri, paprs, d_h_vcol, &
1241           d_qt, d_ec)           d_qt, d_ec)
1242    
1243      ! Calcul des tendances traceurs      ! Calcul des tendances traceurs
1244      call phytrac(itap, lmt_pas, julien, time, firstcal, lafin, dtphys, u, t, &      call phytrac(lmt_pas, julien, time, firstcal, lafin, dtphys, t, paprs, &
1245           paprs, play, mfu, mfd, pde_u, pen_d, ycoefh, fm_therm, entr_therm, &           play, mfu, mfd, pde_u, pen_d, ycoefh, fm_therm, entr_therm, yu1, &
1246           yu1, yv1, ftsol, pctsrf, frac_impa, frac_nucl, pphis, albsol, rhcl, &           yv1, ftsol, pctsrf, frac_impa, frac_nucl, da, phi, mp, upwd, dnwd, &
1247           cldfra, rneb, diafra, cldliq, pmflxr, pmflxs, prfl, psfl, da, phi, &           tr_seri, zmasse, ncid_startphy)
1248           mp, upwd, dnwd, tr_seri, zmasse)  
1249        IF (offline) call phystokenc(dtphys, t, mfu, mfd, pen_u, pde_u, pen_d, &
1250      IF (offline) call phystokenc(dtphys, rlon, rlat, t, mfu, mfd, pen_u, &           pde_d, fm_therm, entr_therm, ycoefh, yu1, yv1, ftsol, pctsrf, &
1251           pde_u, pen_d, pde_d, fm_therm, entr_therm, ycoefh, yu1, yv1, ftsol, &           frac_impa, frac_nucl, pphis, airephy, dtphys)
          pctsrf, frac_impa, frac_nucl, pphis, airephy, dtphys, itap)  
1252    
1253      ! Calculer le transport de l'eau et de l'energie (diagnostique)      ! Calculer le transport de l'eau et de l'energie (diagnostique)
1254      CALL transp(paprs, zxtsol, t_seri, q_seri, u_seri, v_seri, zphi, ve, vq, &      CALL transp(paprs, t_seri, q_seri, u_seri, v_seri, zphi, ve, vq, ue, uq)
          ue, uq)  
1255    
1256      ! diag. bilKP      ! diag. bilKP
1257    
1258      CALL transp_lay(paprs, zxtsol, t_seri, q_seri, u_seri, v_seri, zphi, &      CALL transp_lay(paprs, t_seri, q_seri, u_seri, v_seri, zphi, &
1259           ve_lay, vq_lay, ue_lay, uq_lay)           ve_lay, vq_lay, ue_lay, uq_lay)
1260    
1261      ! Accumuler les variables a stocker dans les fichiers histoire:      ! Accumuler les variables a stocker dans les fichiers histoire:
# Line 1434  contains Line 1271  contains
1271         END DO         END DO
1272      END DO      END DO
1273    
1274      IF (if_ebil >= 1) THEN      IF (if_ebil >= 1) THEN
1275         tit = 'after physic'         tit = 'after physic'
1276         CALL diagetpq(airephy, tit, ip_ebil, 1, 1, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 1, 1, dtphys, t_seri, q_seri, &
1277              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
1278         ! Comme les tendances de la physique sont ajoute dans la dynamique,         ! Comme les tendances de la physique sont ajoute dans la dynamique,
1279         ! on devrait avoir que la variation d'entalpie par la dynamique         ! on devrait avoir que la variation d'entalpie par la dynamique
1280         ! est egale a la variation de la physique au pas de temps precedent.         ! est egale a la variation de la physique au pas de temps precedent.
1281         ! Donc la somme de ces 2 variations devrait etre nulle.         ! Donc la somme de ces 2 variations devrait etre nulle.
# Line 1472  contains Line 1309  contains
1309      DO iq = 3, nqmx      DO iq = 3, nqmx
1310         DO k = 1, llm         DO k = 1, llm
1311            DO i = 1, klon            DO i = 1, klon
1312               d_qx(i, k, iq) = (tr_seri(i, k, iq-2) - qx(i, k, iq)) / dtphys               d_qx(i, k, iq) = (tr_seri(i, k, iq - 2) - qx(i, k, iq)) / dtphys
1313            ENDDO            ENDDO
1314         ENDDO         ENDDO
1315      ENDDO      ENDDO
# Line 1485  contains Line 1322  contains
1322         ENDDO         ENDDO
1323      ENDDO      ENDDO
1324    
1325      ! Ecriture des sorties      CALL histwrite_phy("phis", pphis)
1326      call write_histins      CALL histwrite_phy("aire", airephy)
1327        CALL histwrite_phy("psol", paprs(:, 1))
1328        CALL histwrite_phy("precip", rain_fall + snow_fall)
1329        CALL histwrite_phy("plul", rain_lsc + snow_lsc)
1330        CALL histwrite_phy("pluc", rain_con + snow_con)
1331        CALL histwrite_phy("tsol", zxtsol)
1332        CALL histwrite_phy("t2m", zt2m)
1333        CALL histwrite_phy("q2m", zq2m)
1334        CALL histwrite_phy("u10m", zu10m)
1335        CALL histwrite_phy("v10m", zv10m)
1336        CALL histwrite_phy("snow", snow_fall)
1337        CALL histwrite_phy("cdrm", cdragm)
1338        CALL histwrite_phy("cdrh", cdragh)
1339        CALL histwrite_phy("topl", toplw)
1340        CALL histwrite_phy("evap", evap)
1341        CALL histwrite_phy("sols", solsw)
1342        CALL histwrite_phy("soll", sollw)
1343        CALL histwrite_phy("solldown", sollwdown)
1344        CALL histwrite_phy("bils", bils)
1345        CALL histwrite_phy("sens", - sens)
1346        CALL histwrite_phy("fder", fder)
1347        CALL histwrite_phy("dtsvdfo", d_ts(:, is_oce))
1348        CALL histwrite_phy("dtsvdft", d_ts(:, is_ter))
1349        CALL histwrite_phy("dtsvdfg", d_ts(:, is_lic))
1350        CALL histwrite_phy("dtsvdfi", d_ts(:, is_sic))
1351    
1352      ! Si c'est la fin, il faut conserver l'etat de redemarrage      DO nsrf = 1, nbsrf
1353      IF (lafin) THEN         CALL histwrite_phy("pourc_"//clnsurf(nsrf), pctsrf(:, nsrf)*100.)
1354         itau_phy = itau_phy + itap         CALL histwrite_phy("fract_"//clnsurf(nsrf), pctsrf(:, nsrf))
1355         CALL phyredem("restartphy.nc", rlat, rlon, pctsrf, ftsol, ftsoil, &         CALL histwrite_phy("sens_"//clnsurf(nsrf), fluxt(:, 1, nsrf))
1356              tslab, seaice, fqsurf, qsol, fsnow, falbe, falblw, fevap, &         CALL histwrite_phy("lat_"//clnsurf(nsrf), fluxlat(:, nsrf))
1357              rain_fall, snow_fall, solsw, sollw, dlw, radsol, frugs, &         CALL histwrite_phy("tsol_"//clnsurf(nsrf), ftsol(:, nsrf))
1358              agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, t_ancien, &         CALL histwrite_phy("taux_"//clnsurf(nsrf), fluxu(:, 1, nsrf))
1359              q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, w01)         CALL histwrite_phy("tauy_"//clnsurf(nsrf), fluxv(:, 1, nsrf))
1360      ENDIF         CALL histwrite_phy("rugs_"//clnsurf(nsrf), frugs(:, nsrf))
1361           CALL histwrite_phy("albe_"//clnsurf(nsrf), falbe(:, nsrf))
1362      firstcal = .FALSE.      END DO
   
   contains  
   
     subroutine write_histins  
   
       ! From phylmd/write_histins.h, version 1.2 2005/05/25 13:10:09  
   
       use dimens_m, only: iim, jjm  
       USE histsync_m, ONLY: histsync  
       USE histwrite_m, ONLY: histwrite  
   
       real zout  
       integer itau_w ! pas de temps ecriture  
       REAL zx_tmp_2d(iim, jjm + 1), zx_tmp_3d(iim, jjm + 1, llm)  
   
       !--------------------------------------------------  
   
       IF (ok_instan) THEN  
          ! Champs 2D:  
   
          zsto = dtphys * ecrit_ins  
          zout = dtphys * ecrit_ins  
          itau_w = itau_phy + itap  
   
          i = NINT(zout/zsto)  
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, pphis, zx_tmp_2d)  
          CALL histwrite(nid_ins, "phis", itau_w, zx_tmp_2d)  
   
          i = NINT(zout/zsto)  
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, airephy, zx_tmp_2d)  
          CALL histwrite(nid_ins, "aire", itau_w, zx_tmp_2d)  
   
          DO i = 1, klon  
             zx_tmp_fi2d(i) = paprs(i, 1)  
          ENDDO  
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
          CALL histwrite(nid_ins, "psol", itau_w, zx_tmp_2d)  
   
          DO i = 1, klon  
             zx_tmp_fi2d(i) = rain_fall(i) + snow_fall(i)  
          ENDDO  
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
          CALL histwrite(nid_ins, "precip", itau_w, zx_tmp_2d)  
   
          DO i = 1, klon  
             zx_tmp_fi2d(i) = rain_lsc(i) + snow_lsc(i)  
          ENDDO  
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
          CALL histwrite(nid_ins, "plul", itau_w, zx_tmp_2d)  
   
          DO i = 1, klon  
             zx_tmp_fi2d(i) = rain_con(i) + snow_con(i)  
          ENDDO  
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
          CALL histwrite(nid_ins, "pluc", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zxtsol, zx_tmp_2d)  
          CALL histwrite(nid_ins, "tsol", itau_w, zx_tmp_2d)  
          !ccIM  
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zt2m, zx_tmp_2d)  
          CALL histwrite(nid_ins, "t2m", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zq2m, zx_tmp_2d)  
          CALL histwrite(nid_ins, "q2m", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zu10m, zx_tmp_2d)  
          CALL histwrite(nid_ins, "u10m", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zv10m, zx_tmp_2d)  
          CALL histwrite(nid_ins, "v10m", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, snow_fall, zx_tmp_2d)  
          CALL histwrite(nid_ins, "snow", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, cdragm, zx_tmp_2d)  
          CALL histwrite(nid_ins, "cdrm", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, cdragh, zx_tmp_2d)  
          CALL histwrite(nid_ins, "cdrh", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, toplw, zx_tmp_2d)  
          CALL histwrite(nid_ins, "topl", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, evap, zx_tmp_2d)  
          CALL histwrite(nid_ins, "evap", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, solsw, zx_tmp_2d)  
          CALL histwrite(nid_ins, "sols", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sollw, zx_tmp_2d)  
          CALL histwrite(nid_ins, "soll", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sollwdown, zx_tmp_2d)  
          CALL histwrite(nid_ins, "solldown", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, bils, zx_tmp_2d)  
          CALL histwrite(nid_ins, "bils", itau_w, zx_tmp_2d)  
   
          zx_tmp_fi2d(1:klon) = -1*sens(1:klon)  
          ! CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sens, zx_tmp_2d)  
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
          CALL histwrite(nid_ins, "sens", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, fder, zx_tmp_2d)  
          CALL histwrite(nid_ins, "fder", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_oce), zx_tmp_2d)  
          CALL histwrite(nid_ins, "dtsvdfo", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_ter), zx_tmp_2d)  
          CALL histwrite(nid_ins, "dtsvdft", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_lic), zx_tmp_2d)  
          CALL histwrite(nid_ins, "dtsvdfg", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_sic), zx_tmp_2d)  
          CALL histwrite(nid_ins, "dtsvdfi", itau_w, zx_tmp_2d)  
   
          DO nsrf = 1, nbsrf  
             !XXX  
             zx_tmp_fi2d(1 : klon) = pctsrf(1 : klon, nsrf)*100.  
             CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
             CALL histwrite(nid_ins, "pourc_"//clnsurf(nsrf), itau_w, &  
                  zx_tmp_2d)  
   
             zx_tmp_fi2d(1 : klon) = pctsrf(1 : klon, nsrf)  
             CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
             CALL histwrite(nid_ins, "fract_"//clnsurf(nsrf), itau_w, &  
                  zx_tmp_2d)  
   
             zx_tmp_fi2d(1 : klon) = fluxt(1 : klon, 1, nsrf)  
             CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
             CALL histwrite(nid_ins, "sens_"//clnsurf(nsrf), itau_w, &  
                  zx_tmp_2d)  
   
             zx_tmp_fi2d(1 : klon) = fluxlat(1 : klon, nsrf)  
             CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
             CALL histwrite(nid_ins, "lat_"//clnsurf(nsrf), itau_w, &  
                  zx_tmp_2d)  
   
             zx_tmp_fi2d(1 : klon) = ftsol(1 : klon, nsrf)  
             CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
             CALL histwrite(nid_ins, "tsol_"//clnsurf(nsrf), itau_w, &  
                  zx_tmp_2d)  
   
             zx_tmp_fi2d(1 : klon) = fluxu(1 : klon, 1, nsrf)  
             CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
             CALL histwrite(nid_ins, "taux_"//clnsurf(nsrf), itau_w, &  
                  zx_tmp_2d)  
   
             zx_tmp_fi2d(1 : klon) = fluxv(1 : klon, 1, nsrf)  
             CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
             CALL histwrite(nid_ins, "tauy_"//clnsurf(nsrf), itau_w, &  
                  zx_tmp_2d)  
   
             zx_tmp_fi2d(1 : klon) = frugs(1 : klon, nsrf)  
             CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
             CALL histwrite(nid_ins, "rugs_"//clnsurf(nsrf), itau_w, &  
                  zx_tmp_2d)  
   
             zx_tmp_fi2d(1 : klon) = falbe(1 : klon, nsrf)  
             CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)  
             CALL histwrite(nid_ins, "albe_"//clnsurf(nsrf), itau_w, &  
                  zx_tmp_2d)  
   
          END DO  
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, albsol, zx_tmp_2d)  
          CALL histwrite(nid_ins, "albs", itau_w, zx_tmp_2d)  
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, albsollw, zx_tmp_2d)  
          CALL histwrite(nid_ins, "albslw", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zxrugs, zx_tmp_2d)  
          CALL histwrite(nid_ins, "rugs", itau_w, zx_tmp_2d)  
   
          !HBTM2  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_pblh, zx_tmp_2d)  
          CALL histwrite(nid_ins, "s_pblh", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_pblt, zx_tmp_2d)  
          CALL histwrite(nid_ins, "s_pblt", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_lcl, zx_tmp_2d)  
          CALL histwrite(nid_ins, "s_lcl", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_capCL, zx_tmp_2d)  
          CALL histwrite(nid_ins, "s_capCL", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_oliqCL, zx_tmp_2d)  
          CALL histwrite(nid_ins, "s_oliqCL", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_cteiCL, zx_tmp_2d)  
          CALL histwrite(nid_ins, "s_cteiCL", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_therm, zx_tmp_2d)  
          CALL histwrite(nid_ins, "s_therm", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb1, zx_tmp_2d)  
          CALL histwrite(nid_ins, "s_trmb1", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb2, zx_tmp_2d)  
          CALL histwrite(nid_ins, "s_trmb2", itau_w, zx_tmp_2d)  
   
          CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb3, zx_tmp_2d)  
          CALL histwrite(nid_ins, "s_trmb3", itau_w, zx_tmp_2d)  
   
          ! Champs 3D:  
   
          CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, t_seri, zx_tmp_3d)  
          CALL histwrite(nid_ins, "temp", itau_w, zx_tmp_3d)  
   
          CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, u_seri, zx_tmp_3d)  
          CALL histwrite(nid_ins, "vitu", itau_w, zx_tmp_3d)  
   
          CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, v_seri, zx_tmp_3d)  
          CALL histwrite(nid_ins, "vitv", itau_w, zx_tmp_3d)  
   
          CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, zphi, zx_tmp_3d)  
          CALL histwrite(nid_ins, "geop", itau_w, zx_tmp_3d)  
   
          CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, play, zx_tmp_3d)  
          CALL histwrite(nid_ins, "pres", itau_w, zx_tmp_3d)  
   
          CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, d_t_vdf, zx_tmp_3d)  
          CALL histwrite(nid_ins, "dtvdf", itau_w, zx_tmp_3d)  
   
          CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, d_q_vdf, zx_tmp_3d)  
          CALL histwrite(nid_ins, "dqvdf", itau_w, zx_tmp_3d)  
1363    
1364           call histsync(nid_ins)      CALL histwrite_phy("albs", albsol)
1365        ENDIF      CALL histwrite_phy("rugs", zxrugs)
1366        CALL histwrite_phy("s_pblh", s_pblh)
1367        CALL histwrite_phy("s_pblt", s_pblt)
1368        CALL histwrite_phy("s_lcl", s_lcl)
1369        CALL histwrite_phy("s_capCL", s_capCL)
1370        CALL histwrite_phy("s_oliqCL", s_oliqCL)
1371        CALL histwrite_phy("s_cteiCL", s_cteiCL)
1372        CALL histwrite_phy("s_therm", s_therm)
1373        CALL histwrite_phy("s_trmb1", s_trmb1)
1374        CALL histwrite_phy("s_trmb2", s_trmb2)
1375        CALL histwrite_phy("s_trmb3", s_trmb3)
1376        if (conv_emanuel) CALL histwrite_phy("ptop", ema_pct)
1377        CALL histwrite_phy("temp", t_seri)
1378        CALL histwrite_phy("vitu", u_seri)
1379        CALL histwrite_phy("vitv", v_seri)
1380        CALL histwrite_phy("geop", zphi)
1381        CALL histwrite_phy("pres", play)
1382        CALL histwrite_phy("dtvdf", d_t_vdf)
1383        CALL histwrite_phy("dqvdf", d_q_vdf)
1384        CALL histwrite_phy("rhum", zx_rh)
1385    
1386        if (ok_instan) call histsync(nid_ins)
1387    
1388        IF (lafin) then
1389           call NF95_CLOSE(ncid_startphy)
1390           CALL phyredem(pctsrf, ftsol, ftsoil, fqsurf, qsol, &
1391                fsnow, falbe, fevap, rain_fall, snow_fall, solsw, sollw, dlw, &
1392                radsol, frugs, agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, &
1393                t_ancien, q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, &
1394                w01)
1395        end IF
1396    
1397      end subroutine write_histins      firstcal = .FALSE.
1398    
1399    END SUBROUTINE physiq    END SUBROUTINE physiq
1400    

Legend:
Removed from v.101  
changed lines
  Added in v.195

  ViewVC Help
Powered by ViewVC 1.1.21