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

Diff of /trunk/phylmd/physiq.f

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

trunk/libf/phylmd/physiq.f90 revision 72 by guez, Tue Jul 23 13:00:07 2013 UTC trunk/Sources/phylmd/physiq.f revision 202 by guez, Wed Jun 8 12:23:41 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, d_ps, dudyn, PVteta)         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
     USE calendar, ONLY: ymds2ju  
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_ins, ksta, ksta_ter, ok_kzmin, &
22           ecrit_mth, 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, soil_model           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, lmt_pas
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
34      use diagcld2_m, only: diagcld2      use diagcld2_m, only: diagcld2
35      use diagetpq_m, only: diagetpq      use diagetpq_m, only: diagetpq
36      use diagphy_m, only: diagphy      use diagphy_m, only: diagphy
37      USE dimens_m, ONLY: iim, jjm, llm, nqmx      USE dimens_m, ONLY: llm, nqmx
38      USE dimphy, ONLY: klon, nbtr      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      USE histsync_m, ONLY: histsync
46      USE histwrite_m, ONLY: histwrite      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_histhf_m, ONLY: ini_histhf      USE ini_histins_m, ONLY: ini_histins, nid_ins
50      USE ini_histday_m, ONLY: ini_histday      use netcdf95, only: NF95_CLOSE
     USE ini_histins_m, ONLY: ini_histins  
51      use newmicro_m, only: newmicro      use newmicro_m, only: newmicro
52      USE oasis_m, ONLY: ok_oasis      use nr_util, only: assert
53      USE orbite_m, ONLY: orbite, zenang      use nuage_m, only: nuage
54        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
70      USE yoethf_m, ONLY: r2es, rvtmp2      USE yoethf_m, ONLY: r2es, rvtmp2
71        use zenang_m, only: zenang
72    
     ! Arguments:  
   
     REAL, intent(in):: rdayvrai  
     ! (elapsed time since January 1st 0h of the starting year, in days)  
   
     REAL, intent(in):: time ! heure de la journée en fraction de jour  
     REAL, intent(in):: dtphys ! pas d'integration pour la physique (seconde)  
73      logical, intent(in):: lafin ! dernier passage      logical, intent(in):: lafin ! dernier passage
74    
75      REAL, intent(in):: paprs(klon, llm + 1)      integer, intent(in):: dayvrai
76      ! (pression pour chaque inter-couche, en Pa)      ! current day number, based at value 1 on January 1st of annee_ref
77    
78      REAL, intent(in):: play(klon, llm)      REAL, intent(in):: time ! heure de la journ\'ee en fraction de jour
     ! (input pression pour le mileu de chaque couche (en Pa))  
79    
80      REAL, intent(in):: pphi(klon, llm)      REAL, intent(in):: paprs(:, :) ! (klon, llm + 1)
81      ! (input geopotentiel de chaque couche (g z) (reference sol))      ! pression pour chaque inter-couche, en Pa
82    
83      REAL, intent(in):: pphis(klon) ! input geopotentiel du sol      REAL, intent(in):: play(:, :) ! (klon, llm)
84        ! pression pour le mileu de chaque couche (en Pa)
85    
86      REAL, intent(in):: u(klon, llm)      REAL, intent(in):: pphi(:, :) ! (klon, llm)
87      ! vitesse dans la direction X (de O a E) en m/s      ! géopotentiel de chaque couche (référence sol)
88    
89      REAL, intent(in):: v(klon, llm) ! vitesse Y (de S a N) en m/s      REAL, intent(in):: pphis(:) ! (klon) géopotentiel du sol
     REAL, intent(in):: t(klon, llm) ! input temperature (K)  
90    
91      REAL, intent(in):: qx(klon, llm, nqmx)      REAL, intent(in):: u(:, :) ! (klon, llm)
92      ! (humidité spécifique et fractions massiques des autres traceurs)      ! vitesse dans la direction X (de O a E) en m / s
93    
94      REAL omega(klon, llm) ! input vitesse verticale en Pa/s      REAL, intent(in):: v(:, :) ! (klon, llm) vitesse Y (de S a N) en m / s
95      REAL, intent(out):: d_u(klon, llm) ! tendance physique de "u" (m/s/s)      REAL, intent(in):: t(:, :) ! (klon, llm) temperature (K)
     REAL, intent(out):: d_v(klon, llm) ! tendance physique de "v" (m/s/s)  
     REAL, intent(out):: d_t(klon, llm) ! tendance physique de "t" (K/s)  
     REAL d_qx(klon, llm, nqmx) ! output tendance physique de "qx" (kg/kg/s)  
     REAL d_ps(klon) ! output tendance physique de la pression au sol  
96    
97      LOGICAL:: firstcal = .true.      REAL, intent(in):: qx(:, :, :) ! (klon, llm, nqmx)
98        ! (humidit\'e sp\'ecifique et fractions massiques des autres traceurs)
99    
100        REAL, intent(in):: omega(:, :) ! (klon, llm) vitesse verticale en Pa / s
101        REAL, intent(out):: d_u(:, :) ! (klon, llm) tendance physique de "u" (m s-2)
102        REAL, intent(out):: d_v(:, :) ! (klon, llm) tendance physique de "v" (m s-2)
103        REAL, intent(out):: d_t(:, :) ! (klon, llm) tendance physique de "t" (K / s)
104    
105      INTEGER nbteta      REAL, intent(out):: d_qx(:, :, :) ! (klon, llm, nqmx)
106      PARAMETER(nbteta = 3)      ! tendance physique de "qx" (s-1)
107    
108      REAL PVteta(klon, nbteta)      ! Local:
     ! (output vorticite potentielle a des thetas constantes)  
109    
110      LOGICAL ok_gust ! pour activer l'effet des gust sur flux surface      LOGICAL:: firstcal = .true.
     PARAMETER (ok_gust = .FALSE.)  
111    
112      LOGICAL check ! Verifier la conservation du modele en eau      LOGICAL, PARAMETER:: check = .FALSE.
113      PARAMETER (check = .FALSE.)      ! 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    
118      ! Parametres lies au coupleur OASIS:      ! pour phystoke avec thermiques
     INTEGER, SAVE:: npas, nexca  
     logical rnpb  
     parameter(rnpb = .true.)  
   
     character(len = 6):: ocean = 'force '  
     ! (type de modèle océan à utiliser: "force" ou "slab" mais pas "couple")  
   
     ! "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  
   
     ! Modele thermique du sol, a activer pour le cycle diurne:  
     logical:: ok_veget = .false. ! type de modele de vegetation utilise  
   
     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.)  
   
     ! 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)
121      real, save:: q2(klon, llm + 1, nbsrf)      real, save:: q2(klon, llm + 1, nbsrf)
122    
123      INTEGER ivap ! indice de traceurs pour vapeur d'eau      INTEGER, PARAMETER:: ivap = 1 ! indice de traceur pour vapeur d'eau
124      PARAMETER (ivap = 1)      INTEGER, PARAMETER:: iliq = 2 ! indice de traceur pour eau liquide
     INTEGER iliq ! indice de traceurs pour eau liquide  
     PARAMETER (iliq = 2)  
125    
126      REAL, save:: t_ancien(klon, llm), q_ancien(klon, llm)      REAL, save:: t_ancien(klon, llm), q_ancien(klon, llm)
127      LOGICAL, save:: ancien_ok      LOGICAL, save:: ancien_ok
128    
129      REAL d_t_dyn(klon, llm) ! tendance dynamique pour "t" (K/s)      REAL d_t_dyn(klon, llm) ! tendance dynamique pour "t" (K / s)
130      REAL d_q_dyn(klon, llm) ! tendance dynamique pour "q" (kg/kg/s)      REAL d_q_dyn(klon, llm) ! tendance dynamique pour "q" (kg / kg / s)
131    
132      real da(klon, llm), phi(klon, llm, llm), mp(klon, llm)      real da(klon, llm), phi(klon, llm, llm), mp(klon, llm)
133    
     !IM Amip2 PV a theta constante  
   
     CHARACTER(LEN = 3) ctetaSTD(nbteta)  
     DATA ctetaSTD/'350', '380', '405'/  
     REAL rtetaSTD(nbteta)  
     DATA rtetaSTD/350., 380., 405./  
   
     !MI Amip2 PV a theta constante  
   
134      REAL swdn0(klon, llm + 1), swdn(klon, llm + 1)      REAL swdn0(klon, llm + 1), swdn(klon, llm + 1)
135      REAL swup0(klon, llm + 1), swup(klon, llm + 1)      REAL swup0(klon, llm + 1), swup(klon, llm + 1)
136      SAVE swdn0, swdn, swup0, swup      SAVE swdn0, swdn, swup0, swup
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
   
     !IM 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)
144    
145      ! flwp, fiwp = Liquid Water Path & Ice Water Path (kg/m2)      ! flwp, fiwp = Liquid Water Path & Ice Water Path (kg / m2)
146      ! flwc, fiwc = Liquid Water Content & Ice Water Content (kg/kg)      ! flwc, fiwc = Liquid Water Content & Ice Water Content (kg / kg)
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.'/  
   
     !IM ISCCP simulator v3.4  
   
     integer nid_hf, nid_hf3d  
     save nid_hf, nid_hf3d  
   
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 272  contains Line 165  contains
165      REAL fluxlat(klon, nbsrf)      REAL fluxlat(klon, nbsrf)
166      SAVE fluxlat      SAVE fluxlat
167    
168      REAL fqsurf(klon, nbsrf)      REAL, save:: fqsurf(klon, nbsrf)
169      SAVE fqsurf ! humidite de l'air au contact de la surface      ! humidite de l'air au contact de la surface
   
     REAL, save:: qsol(klon) ! hauteur d'eau dans le sol  
170    
171      REAL fsnow(klon, nbsrf)      REAL, save:: qsol(klon)
172      SAVE fsnow ! epaisseur neigeuse      ! column-density of water in soil, in kg m-2
173    
174      REAL falbe(klon, nbsrf)      REAL, save:: fsnow(klon, nbsrf) ! epaisseur neigeuse
175      SAVE falbe ! albedo par type de surface      REAL, save:: falbe(klon, nbsrf) ! albedo visible par type de surface
     REAL falblw(klon, nbsrf)  
     SAVE falblw ! albedo par type de surface  
176    
177      ! Paramètres de l'orographie à l'échelle 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
179      REAL, save:: zstd(klon) ! deviation standard de l'OESM      REAL, save:: zstd(klon) ! deviation standard de l'OESM
180      REAL, save:: zsig(klon) ! pente de l'OESM      REAL, save:: zsig(klon) ! pente de l'OESM
# Line 294  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)
     REAL agesno(klon, nbsrf)  
     SAVE agesno ! age de la neige  
191    
192      REAL run_off_lic_0(klon)      ! Variables li\'ees \`a la convection d'Emanuel :
193      SAVE run_off_lic_0      REAL, save:: Ma(klon, llm) ! undilute upward mass flux
194      !KE43      REAL, save:: qcondc(klon, llm) ! in-cld water content from convect
     ! Variables liees a la convection de K. Emanuel (sb):  
   
     REAL bas, top ! cloud base and top levels  
     SAVE bas  
     SAVE top  
   
     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 329  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 344  contains Line 219  contains
219      REAL frac_impa(klon, llm) ! fractions d'aerosols lessivees (impaction)      REAL frac_impa(klon, llm) ! fractions d'aerosols lessivees (impaction)
220      REAL frac_nucl(klon, llm) ! idem (nucleation)      REAL frac_nucl(klon, llm) ! idem (nucleation)
221    
222      REAL, save:: rain_fall(klon) ! pluie      REAL, save:: rain_fall(klon)
223      REAL, save:: snow_fall(klon) ! neige      ! liquid water mass flux (kg / m2 / s), positive down
224    
225        REAL, save:: snow_fall(klon)
226        ! solid water mass flux (kg / m2 / s), positive down
227    
228      REAL rain_tiedtke(klon), snow_tiedtke(klon)      REAL rain_tiedtke(klon), snow_tiedtke(klon)
229    
# Line 354  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
239      REAL uq(klon) ! integr. verticale du transport zonal de l'eau      REAL uq(klon) ! integr. verticale du transport zonal de l'eau
240    
241      REAL frugs(klon, nbsrf) ! longueur de rugosite      REAL, save:: frugs(klon, nbsrf) ! longueur de rugosite
     save frugs  
242      REAL zxrugs(klon) ! longueur de rugosite      REAL zxrugs(klon) ! longueur de rugosite
243    
244      ! Conditions aux limites      ! Conditions aux limites
245    
246      INTEGER julien      INTEGER julien
   
     INTEGER, SAVE:: lmt_pas ! number of time steps of "physics" per day  
247      REAL, save:: pctsrf(klon, nbsrf) ! percentage of surface      REAL, save:: pctsrf(klon, nbsrf) ! percentage of surface
248      REAL pctsrf_new(klon, nbsrf) ! pourcentage surfaces issus d'ORCHIDEE      REAL, save:: albsol(klon) ! albedo du sol total visible
   
     REAL albsol(klon)  
     SAVE albsol ! albedo du sol total  
     REAL albsollw(klon)  
     SAVE albsollw ! albedo du sol total  
   
249      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
250    
     ! Declaration des procedures appelees  
   
     EXTERNAL alboc ! calculer l'albedo sur ocean  
     !KE43  
     EXTERNAL conema3 ! convect4.3  
     EXTERNAL nuage ! calculer les proprietes radiatives  
     EXTERNAL transp ! transport total de l'eau et de l'energie  
   
     ! Variables locales  
   
251      real, save:: clwcon(klon, llm), rnebcon(klon, llm)      real, save:: clwcon(klon, llm), rnebcon(klon, llm)
252      real, save:: clwcon0(klon, llm), rnebcon0(klon, llm)      real, save:: clwcon0(klon, llm), rnebcon0(klon, llm)
253    
# Line 411  contains Line 269  contains
269      REAL zxfluxu(klon, llm)      REAL zxfluxu(klon, llm)
270      REAL zxfluxv(klon, llm)      REAL zxfluxv(klon, llm)
271    
272      ! Le rayonnement n'est pas calculé tous les pas, il faut donc que      ! Le rayonnement n'est pas calcul\'e tous les pas, il faut donc que
273      ! les variables soient rémanentes.      ! les variables soient r\'emanentes.
274      REAL, save:: heat(klon, llm) ! chauffage solaire      REAL, save:: heat(klon, llm) ! chauffage solaire
275      REAL heat0(klon, llm) ! chauffage solaire ciel clair      REAL, save:: heat0(klon, llm) ! chauffage solaire ciel clair
276      REAL, save:: cool(klon, llm) ! refroidissement infrarouge      REAL, save:: cool(klon, llm) ! refroidissement infrarouge
277      REAL cool0(klon, llm) ! refroidissement infrarouge ciel clair      REAL, save:: cool0(klon, llm) ! refroidissement infrarouge ciel clair
278      REAL, save:: topsw(klon), toplw(klon), solsw(klon)      REAL, save:: topsw(klon), toplw(klon), solsw(klon)
279      REAL, save:: sollw(klon) ! rayonnement infrarouge montant à la surface      REAL, save:: sollw(klon) ! rayonnement infrarouge montant \`a la surface
280      real, save:: sollwdown(klon) ! downward LW flux at surface      real, save:: sollwdown(klon) ! downward LW flux at surface
281      REAL, save:: topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)      REAL, save:: topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)
282      REAL albpla(klon)      REAL, save:: albpla(klon)
283      REAL fsollw(klon, nbsrf) ! bilan flux IR pour chaque sous surface      REAL fsollw(klon, nbsrf) ! bilan flux IR pour chaque sous-surface
284      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  
285    
286      INTEGER itaprad      REAL conv_q(klon, llm) ! convergence de l'humidite (kg / kg / s)
287      SAVE itaprad      REAL conv_t(klon, llm) ! convergence of temperature (K / s)
288    
289      REAL conv_q(klon, llm) ! convergence de l'humidite (kg/kg/s)      REAL cldl(klon), cldm(klon), cldh(klon) ! nuages bas, moyen et haut
290      REAL conv_t(klon, llm) ! convergence of temperature (K/s)      REAL cldt(klon), cldq(klon) ! nuage total, eau liquide integree
   
     REAL cldl(klon), cldm(klon), cldh(klon) !nuages bas, moyen et haut  
     REAL cldt(klon), cldq(klon) !nuage total, eau liquide integree  
291    
292      REAL zxtsol(klon), zxqsurf(klon), zxsnow(klon), zxfluxlat(klon)      REAL zxtsol(klon), zxqsurf(klon), zxsnow(klon), zxfluxlat(klon)
293    
294      REAL dist, rmu0(klon), fract(klon)      REAL dist, mu0(klon), fract(klon)
295      REAL zdtime ! pas de temps du rayonnement (s)      real longi
     real zlongi  
296      REAL z_avant(klon), z_apres(klon), z_factor(klon)      REAL z_avant(klon), z_apres(klon), z_factor(klon)
297      REAL za, zb      REAL za, zb
298      REAL zx_t, zx_qs, zdelta, zcor      REAL zx_t, zx_qs, zcor
299      real zqsat(klon, llm)      real zqsat(klon, llm)
300      INTEGER i, k, iq, nsrf      INTEGER i, k, iq, nsrf
301      REAL, PARAMETER:: t_coup = 234.      REAL, PARAMETER:: t_coup = 234.
302      REAL zphi(klon, llm)      REAL zphi(klon, llm)
303    
304      !IM cf. AM Variables locales pour la CLA (hbtm2)      ! cf. Anne Mathieu, variables pour la couche limite atmosphérique (hbtm)
305    
306      REAL, SAVE:: pblh(klon, nbsrf) ! Hauteur de couche limite      REAL, SAVE:: pblh(klon, nbsrf) ! Hauteur de couche limite
307      REAL, SAVE:: plcl(klon, nbsrf) ! Niveau de condensation de la CLA      REAL, SAVE:: plcl(klon, nbsrf) ! Niveau de condensation de la CLA
# Line 459  contains Line 311  contains
311      REAL, SAVE:: pblt(klon, nbsrf) ! T a la Hauteur de couche limite      REAL, SAVE:: pblt(klon, nbsrf) ! T a la Hauteur de couche limite
312      REAL, SAVE:: therm(klon, nbsrf)      REAL, SAVE:: therm(klon, nbsrf)
313      REAL, SAVE:: trmb1(klon, nbsrf) ! deep_cape      REAL, SAVE:: trmb1(klon, nbsrf) ! deep_cape
314      REAL, SAVE:: trmb2(klon, nbsrf) ! inhibition      REAL, SAVE:: trmb2(klon, nbsrf) ! inhibition
315      REAL, SAVE:: trmb3(klon, nbsrf) ! Point Omega      REAL, SAVE:: trmb3(klon, nbsrf) ! Point Omega
316      ! Grdeurs de sorties      ! Grandeurs de sorties
317      REAL s_pblh(klon), s_lcl(klon), s_capCL(klon)      REAL s_pblh(klon), s_lcl(klon), s_capCL(klon)
318      REAL s_oliqCL(klon), s_cteiCL(klon), s_pblt(klon)      REAL s_oliqCL(klon), s_cteiCL(klon), s_pblt(klon)
319      REAL s_therm(klon), s_trmb1(klon), s_trmb2(klon)      REAL s_therm(klon), s_trmb1(klon), s_trmb2(klon)
320      REAL s_trmb3(klon)      REAL s_trmb3(klon)
321    
322      ! Variables locales pour la convection de K. Emanuel :      ! Variables pour la convection de K. Emanuel :
323    
324      REAL upwd(klon, llm) ! saturated updraft mass flux      REAL upwd(klon, llm) ! saturated updraft mass flux
325      REAL dnwd(klon, llm) ! saturated downdraft mass flux      REAL dnwd(klon, llm) ! saturated downdraft mass flux
326      REAL dnwd0(klon, llm) ! unsaturated downdraft mass flux      REAL dnwd0(klon, llm) ! unsaturated downdraft mass flux
     REAL tvp(klon, llm) ! virtual temp of lifted parcel  
327      REAL cape(klon) ! CAPE      REAL cape(klon) ! CAPE
328      SAVE cape      SAVE cape
329    
     REAL pbase(klon) ! cloud base pressure  
     SAVE pbase  
     REAL bbase(klon) ! cloud base buoyancy  
     SAVE bbase  
     REAL rflag(klon) ! flag fonctionnement de convect  
330      INTEGER iflagctrl(klon) ! flag fonctionnement de convect      INTEGER iflagctrl(klon) ! flag fonctionnement de convect
     ! -- convect43:  
     REAL dtvpdt1(klon, llm), dtvpdq1(klon, llm)  
     REAL dplcldt(klon), dplcldr(klon)  
331    
332      ! Variables du changement      ! Variables du changement
333    
334      ! con: convection      ! con: convection
335      ! lsc: large scale condensation      ! lsc: large scale condensation
336      ! ajs: ajustement sec      ! ajs: ajustement sec
337      ! eva: évaporation de l'eau liquide nuageuse      ! eva: \'evaporation de l'eau liquide nuageuse
338      ! vdf: vertical diffusion in boundary layer      ! vdf: vertical diffusion in boundary layer
339      REAL d_t_con(klon, llm), d_q_con(klon, llm)      REAL d_t_con(klon, llm), d_q_con(klon, llm)
340      REAL d_u_con(klon, llm), d_v_con(klon, llm)      REAL d_u_con(klon, llm), d_v_con(klon, llm)
# Line 508  contains Line 351  contains
351      REAL prfl(klon, llm + 1), psfl(klon, llm + 1)      REAL prfl(klon, llm + 1), psfl(klon, llm + 1)
352    
353      INTEGER, save:: ibas_con(klon), itop_con(klon)      INTEGER, save:: ibas_con(klon), itop_con(klon)
354        real ema_pct(klon) ! Emanuel pressure at cloud top, in Pa
355    
356      REAL rain_con(klon), rain_lsc(klon)      REAL rain_con(klon), rain_lsc(klon)
357      REAL snow_con(klon), snow_lsc(klon)      REAL, save:: snow_con(klon) ! neige (mm / s)
358        real snow_lsc(klon)
359      REAL d_ts(klon, nbsrf)      REAL d_ts(klon, nbsrf)
360    
361      REAL d_u_vdf(klon, llm), d_v_vdf(klon, llm)      REAL d_u_vdf(klon, llm), d_v_vdf(klon, llm)
# Line 534  contains Line 379  contains
379      integer:: iflag_cldcon = 1      integer:: iflag_cldcon = 1
380      logical ptconv(klon, llm)      logical ptconv(klon, llm)
381    
382      ! Variables locales pour effectuer les appels en série :      ! Variables pour effectuer les appels en s\'erie :
383    
384      REAL t_seri(klon, llm), q_seri(klon, llm)      REAL t_seri(klon, llm), q_seri(klon, llm)
385      REAL ql_seri(klon, llm), qs_seri(klon, llm)      REAL ql_seri(klon, llm)
386      REAL u_seri(klon, llm), v_seri(klon, llm)      REAL u_seri(klon, llm), v_seri(klon, llm)
387        REAL tr_seri(klon, llm, nqmx - 2)
     REAL tr_seri(klon, llm, nbtr)  
     REAL d_tr(klon, llm, nbtr)  
388    
389      REAL zx_rh(klon, llm)      REAL zx_rh(klon, llm)
390    
# Line 550  contains Line 393  contains
393      REAL zustrph(klon), zvstrph(klon)      REAL zustrph(klon), zvstrph(klon)
394      REAL aam, torsfc      REAL aam, torsfc
395    
     REAL dudyn(iim + 1, jjm + 1, llm)  
   
     REAL zx_tmp_fi2d(klon) ! variable temporaire grille physique  
     REAL zx_tmp_2d(iim, jjm + 1), zx_tmp_3d(iim, jjm + 1, llm)  
   
     INTEGER, SAVE:: nid_day, nid_ins  
   
396      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.
397      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.
398      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.
399      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.
400    
     REAL zsto  
   
     logical ok_sync  
401      real date0      real date0
402    
403      ! Variables liées au bilan d'énergie et d'enthalpie :      ! Variables li\'ees au bilan d'\'energie et d'enthalpie :
404      REAL ztsol(klon)      REAL ztsol(klon)
405      REAL d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec      REAL d_h_vcol, d_qt, d_ec
406      REAL, SAVE:: d_h_vcol_phy      REAL, SAVE:: d_h_vcol_phy
     REAL fs_bound, fq_bound  
407      REAL zero_v(klon)      REAL zero_v(klon)
408      CHARACTER(LEN = 15) tit      CHARACTER(LEN = 20) tit
409      INTEGER:: ip_ebil = 0 ! print level for energy conservation diagnostics      INTEGER:: ip_ebil = 0 ! print level for energy conservation diagnostics
410      INTEGER:: if_ebil = 0 ! verbosity for diagnostics of energy conservation      INTEGER:: if_ebil = 0 ! verbosity for diagnostics of energy conservation
411    
412        REAL d_t_ec(klon, llm)
413        ! tendance due \`a la conversion Ec en énergie thermique
414    
     REAL d_t_ec(klon, llm) ! tendance due à la conversion Ec -> E thermique  
415      REAL ZRCPD      REAL ZRCPD
416    
417      REAL t2m(klon, nbsrf), q2m(klon, nbsrf) ! temperature and humidity at 2 m      REAL t2m(klon, nbsrf), q2m(klon, nbsrf) ! temperature and humidity at 2 m
# Line 587  contains Line 421  contains
421    
422      ! Aerosol effects:      ! Aerosol effects:
423    
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 605  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 622  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 630  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/ ocean, ok_veget, ok_journe, ok_mensuel, ok_instan, &      namelist /physiq_nml/ fact_cldcon, facttemps, ok_newmicro, &
470           fact_cldcon, facttemps, ok_newmicro, iflag_cldcon, ratqsbas, &           iflag_cldcon, ratqsbas, ratqshaut, if_ebil, ok_ade, ok_aie, bl95_b0, &
471           ratqshaut, if_ebil, ok_ade, ok_aie, bl95_b0, bl95_b1, iflag_thermals, &           bl95_b1, iflag_thermals, nsplit_thermals
          nsplit_thermals  
472    
473      !----------------------------------------------------------------      !----------------------------------------------------------------
474    
475      IF (if_ebil >= 1) zero_v = 0.      IF (if_ebil >= 1) zero_v = 0.
     ok_sync = .TRUE.  
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 658  contains Line 487  contains
487         piz_ae = 0.         piz_ae = 0.
488         tau_ae = 0.         tau_ae = 0.
489         cg_ae = 0.         cg_ae = 0.
490         rain_con(:) = 0.         rain_con = 0.
491         snow_con(:) = 0.         snow_con = 0.
492         topswai(:) = 0.         topswai = 0.
493         topswad(:) = 0.         topswad = 0.
494         solswai(:) = 0.         solswai = 0.
495         solswad(:) = 0.         solswad = 0.
496    
497         d_u_con = 0.         d_u_con = 0.
498         d_v_con = 0.         d_v_con = 0.
# Line 680  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 696  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("startphy.nc", pctsrf, ftsol, ftsoil, ocean, tslab, &              agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, t_ancien, &
531              seaice, fqsurf, qsol, fsnow, falbe, falblw, fevap, rain_fall, &              q_ancien, ancien_ok, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, &
532              snow_fall, solsw, sollw, dlw, radsol, frugs, agesno, zmea, &              w01, ncid_startphy)
             zstd, zsig, zgam, zthe, zpic, 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)         radpas = lmt_pas / nbapp_rad
538           print *, "radpas = ", radpas
        ! on remet le calendrier a zero  
        IF (raz_date) itau_phy = 0  
   
        PRINT *, 'cycle_diurne = ', cycle_diurne  
        CALL printflag(radpas, ocean /= 'force', ok_oasis, ok_journe, &  
             ok_instan, ok_region)  
   
        IF (dtphys * REAL(radpas) > 21600. .AND. cycle_diurne) THEN  
           print *, "Au minimum 4 appels par jour si cycle diurne"  
           call abort_gcm('physiq', &  
                "Nombre d'appels au rayonnement insuffisant", 1)  
        ENDIF  
539    
540         ! Initialisation pour le schéma de convection d'Emanuel :         ! Initialisation pour le sch\'ema de convection d'Emanuel :
541         IF (iflag_con >= 3) THEN         IF (conv_emanuel) THEN
542            ibas_con = 1            ibas_con = 1
543            itop_con = 1            itop_con = 1
544         ENDIF         ENDIF
# Line 735  contains Line 550  contains
550            rugoro = 0.            rugoro = 0.
551         ENDIF         ENDIF
552    
553         lmt_pas = NINT(86400. / dtphys) ! tous les jours         ecrit_ins = NINT(ecrit_ins / dtphys)
        print *, 'Number of time steps of "physics" per day: ', lmt_pas  
   
        ecrit_ins = NINT(ecrit_ins/dtphys)  
        ecrit_hf = NINT(ecrit_hf/dtphys)  
        ecrit_mth = NINT(ecrit_mth/dtphys)  
        ecrit_tra = NINT(86400.*ecrit_tra/dtphys)  
        ecrit_reg = NINT(ecrit_reg/dtphys)  
   
        ! Initialiser le couplage si necessaire  
   
        npas = 0  
        nexca = 0  
554    
555         ! Initialisation des sorties         ! Initialisation des sorties
556    
557         call ini_histhf(dtphys, nid_hf, nid_hf3d)         call ini_histins(dtphys)
558         call ini_histday(dtphys, ok_journe, nid_day, nqmx)         CALL ymds2ju(annee_ref, 1, day_ref, 0., date0)
        call ini_histins(dtphys, ok_instan, nid_ins)  
        CALL ymds2ju(annee_ref, 1, int(day_ref), 0., date0)  
559         ! Positionner date0 pour initialisation de ORCHIDEE         ! Positionner date0 pour initialisation de ORCHIDEE
560         print *, 'physiq date0: ', date0         print *, 'physiq date0: ', date0
561           CALL phyredem0
562      ENDIF test_firstcal      ENDIF test_firstcal
563    
564      ! Mettre a zero des variables de sortie (pour securite)      ! We will modify variables *_seri and we will not touch variables
565        ! u, v, t, qx:
566      DO i = 1, klon      t_seri = t
567         d_ps(i) = 0.      u_seri = u
568      ENDDO      v_seri = v
569      DO iq = 1, nqmx      q_seri = qx(:, :, ivap)
570         DO k = 1, llm      ql_seri = qx(:, :, iliq)
571            DO i = 1, klon      tr_seri = qx(:, :, 3:nqmx)
              d_qx(i, k, iq) = 0.  
           ENDDO  
        ENDDO  
     ENDDO  
     da = 0.  
     mp = 0.  
     phi = 0.  
   
     ! Ne pas affecter les valeurs entrées de u, v, h, et q :  
   
     DO k = 1, llm  
        DO i = 1, klon  
           t_seri(i, k) = t(i, k)  
           u_seri(i, k) = u(i, k)  
           v_seri(i, k) = v(i, k)  
           q_seri(i, k) = qx(i, k, ivap)  
           ql_seri(i, k) = qx(i, k, iliq)  
           qs_seri(i, k) = 0.  
        ENDDO  
     ENDDO  
     IF (nqmx >= 3) THEN  
        tr_seri(:, :, :nqmx-2) = qx(:, :, 3:nqmx)  
     ELSE  
        tr_seri(:, :, 1) = 0.  
     ENDIF  
572    
573      DO i = 1, klon      ztsol = sum(ftsol * pctsrf, dim = 2)
        ztsol(i) = 0.  
     ENDDO  
     DO nsrf = 1, nbsrf  
        DO i = 1, klon  
           ztsol(i) = ztsol(i) + ftsol(i, nsrf)*pctsrf(i, nsrf)  
        ENDDO  
     ENDDO  
574    
575      IF (if_ebil >= 1) THEN      IF (if_ebil >= 1) THEN
576         tit = 'after dynamics'         tit = 'after dynamics'
577         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, &
578              ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
579              d_ql, d_qs, d_ec)         ! Comme les tendances de la physique sont ajout\'es dans la
580         ! Comme les tendances de la physique sont ajoutés dans la         ! dynamique, la variation d'enthalpie par la dynamique devrait
581         !  dynamique, la variation d'enthalpie par la dynamique devrait         ! \^etre \'egale \`a la variation de la physique au pas de temps
582         !  être égale à la variation de la physique au pas de temps         ! pr\'ec\'edent. Donc la somme de ces 2 variations devrait \^etre
583         !  précédent.  Donc la somme de ces 2 variations devrait être         ! nulle.
        !  nulle.  
584         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, &
585              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, &
586              d_qt, 0., fs_bound, fq_bound)              d_qt, 0.)
587      END IF      END IF
588    
589      ! Diagnostic de la tendance dynamique :      ! Diagnostic de la tendance dynamique :
# Line 845  contains Line 614  contains
614      ! Check temperatures:      ! Check temperatures:
615      CALL hgardfou(t_seri, ftsol)      CALL hgardfou(t_seri, ftsol)
616    
617      ! Incrementer le compteur de la physique      call increment_itap
618      itap = itap + 1      julien = MOD(dayvrai, 360)
     julien = MOD(NINT(rdayvrai), 360)  
619      if (julien == 0) julien = 360      if (julien == 0) julien = 360
620    
621      forall (k = 1: llm) zmasse(:, k) = (paprs(:, k)-paprs(:, k + 1)) / rg      forall (k = 1: llm) zmasse(:, k) = (paprs(:, k) - paprs(:, k + 1)) / rg
622    
623      ! Mettre en action les conditions aux limites (albedo, sst etc.).      ! Prescrire l'ozone :
   
     ! Prescrire l'ozone et calculer l'albedo sur l'ocean.  
624      wo = ozonecm(REAL(julien), paprs)      wo = ozonecm(REAL(julien), paprs)
625    
626      ! Évaporation de l'eau liquide nuageuse :      ! \'Evaporation de l'eau liquide nuageuse :
627      DO k = 1, llm      DO k = 1, llm
628         DO i = 1, klon         DO i = 1, klon
629            zb = MAX(0., ql_seri(i, k))            zb = MAX(0., ql_seri(i, k))
# Line 868  contains Line 634  contains
634      ENDDO      ENDDO
635      ql_seri = 0.      ql_seri = 0.
636    
637      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
638         tit = 'after reevap'         tit = 'after reevap'
639         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, &
640              ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
             d_ql, d_qs, d_ec)  
641         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, &
642              zero_v, zero_v, zero_v, zero_v, ztsol, d_h_vcol, d_qt, d_ec, &              zero_v, zero_v, zero_v, zero_v, ztsol, d_h_vcol, d_qt, d_ec)
             fs_bound, fq_bound)  
   
643      END IF      END IF
644    
645      ! Appeler la diffusion verticale (programme de couche limite)      frugs = MAX(frugs, 0.000015)
646        zxrugs = sum(frugs * pctsrf, dim = 2)
     DO i = 1, klon  
        zxrugs(i) = 0.  
     ENDDO  
     DO nsrf = 1, nbsrf  
        DO i = 1, klon  
           frugs(i, nsrf) = MAX(frugs(i, nsrf), 0.000015)  
        ENDDO  
     ENDDO  
     DO nsrf = 1, nbsrf  
        DO i = 1, klon  
           zxrugs(i) = zxrugs(i) + frugs(i, nsrf)*pctsrf(i, nsrf)  
        ENDDO  
     ENDDO  
647    
648      ! calculs necessaires au calcul de l'albedo dans l'interface      ! Calculs n\'ecessaires au calcul de l'albedo dans l'interface avec
649        ! la surface.
650    
651      CALL orbite(REAL(julien), zlongi, dist)      CALL orbite(REAL(julien), longi, dist)
652      IF (cycle_diurne) THEN      IF (cycle_diurne) THEN
653         zdtime = dtphys * REAL(radpas)         CALL zenang(longi, time, dtphys * radpas, mu0, fract)
        CALL zenang(zlongi, time, zdtime, rmu0, fract)  
654      ELSE      ELSE
655         rmu0 = -999.999         mu0 = - 999.999
656      ENDIF      ENDIF
657    
658      ! Calcul de l'abedo moyen par maille      ! Calcul de l'abedo moyen par maille
659      albsol(:) = 0.      albsol = sum(falbe * pctsrf, dim = 2)
     albsollw(:) = 0.  
     DO nsrf = 1, nbsrf  
        DO i = 1, klon  
           albsol(i) = albsol(i) + falbe(i, nsrf) * pctsrf(i, nsrf)  
           albsollw(i) = albsollw(i) + falblw(i, nsrf) * pctsrf(i, nsrf)  
        ENDDO  
     ENDDO  
660    
661      ! Répartition sous maille des flux longwave et shortwave      ! R\'epartition sous maille des flux longwave et shortwave
662      ! Répartition du longwave par sous-surface linéarisée      ! R\'epartition du longwave par sous-surface lin\'earis\'ee
663    
664      DO nsrf = 1, nbsrf      forall (nsrf = 1: nbsrf)
665         DO i = 1, klon         fsollw(:, nsrf) = sollw + 4. * RSIGMA * ztsol**3 &
666            fsollw(i, nsrf) = sollw(i) &              * (ztsol - ftsol(:, nsrf))
667                 + 4. * RSIGMA * ztsol(i)**3 * (ztsol(i) - ftsol(i, nsrf))         fsolsw(:, nsrf) = solsw * (1. - falbe(:, nsrf)) / (1. - albsol)
668            fsolsw(i, nsrf) = solsw(i) * (1. - falbe(i, nsrf)) / (1. - albsol(i))      END forall
        ENDDO  
     ENDDO  
669    
670      fder = dlw      fder = dlw
671    
672      ! Couche limite:      ! Couche limite:
673    
674      CALL clmain(dtphys, itap, pctsrf, pctsrf_new, t_seri, q_seri, &      CALL clmain(dtphys, pctsrf, t_seri, q_seri, u_seri, v_seri, julien, mu0, &
675           u_seri, v_seri, julien, rmu0, co2_ppm, ok_veget, ocean, &           ftsol, cdmmax, cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, qsol, &
676           ftsol, soil_model, cdmmax, cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, &           paprs, play, fsnow, fqsurf, fevap, falbe, fluxlat, rain_fall, &
677           qsol, paprs, play, fsnow, fqsurf, fevap, falbe, falblw, fluxlat, &           snow_fall, fsolsw, fsollw, fder, rlat, frugs, agesno, rugoro, &
678           rain_fall, snow_fall, fsolsw, fsollw, fder, rlon, rlat, &           d_t_vdf, d_q_vdf, d_u_vdf, d_v_vdf, d_ts, fluxt, fluxq, fluxu, &
679           frugs, firstcal, agesno, rugoro, d_t_vdf, &           fluxv, cdragh, cdragm, q2, dsens, devap, ycoefh, yu1, yv1, t2m, q2m, &
680           d_q_vdf, d_u_vdf, d_v_vdf, d_ts, fluxt, fluxq, fluxu, fluxv, cdragh, &           u10m, v10m, pblh, capCL, oliqCL, cteiCL, pblT, therm, trmb1, trmb2, &
681           cdragm, q2, dsens, devap, ycoefh, yu1, yv1, t2m, q2m, u10m, v10m, &           trmb3, plcl, fqcalving, ffonte, run_off_lic_0)
          pblh, capCL, oliqCL, cteiCL, pblT, therm, trmb1, trmb2, trmb3, plcl, &  
          fqcalving, ffonte, run_off_lic_0, fluxo, fluxg, tslab, seaice)  
682    
683      ! Incrémentation des flux      ! Incr\'ementation des flux
684    
685      zxfluxt = 0.      zxfluxt = 0.
686      zxfluxq = 0.      zxfluxq = 0.
# Line 959  contains Line 698  contains
698      END DO      END DO
699      DO i = 1, klon      DO i = 1, klon
700         sens(i) = - zxfluxt(i, 1) ! flux de chaleur sensible au sol         sens(i) = - zxfluxt(i, 1) ! flux de chaleur sensible au sol
701         evap(i) = - zxfluxq(i, 1) ! flux d'évaporation au sol         evap(i) = - zxfluxq(i, 1) ! flux d'\'evaporation au sol
702         fder(i) = dlw(i) + dsens(i) + devap(i)         fder(i) = dlw(i) + dsens(i) + devap(i)
703      ENDDO      ENDDO
704    
# Line 972  contains Line 711  contains
711         ENDDO         ENDDO
712      ENDDO      ENDDO
713    
714      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
715         tit = 'after clmain'         tit = 'after clmain'
716         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, &
717              ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
             d_ql, d_qs, d_ec)  
718         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, &
719              sens, evap, zero_v, zero_v, ztsol, d_h_vcol, d_qt, d_ec, &              sens, evap, zero_v, zero_v, ztsol, d_h_vcol, d_qt, d_ec)
             fs_bound, fq_bound)  
720      END IF      END IF
721    
722      ! Update surface temperature:      ! Update surface temperature:
723    
724      DO i = 1, klon      DO i = 1, klon
        zxtsol(i) = 0.  
725         zxfluxlat(i) = 0.         zxfluxlat(i) = 0.
726    
727         zt2m(i) = 0.         zt2m(i) = 0.
# Line 995  contains Line 731  contains
731         zxffonte(i) = 0.         zxffonte(i) = 0.
732         zxfqcalving(i) = 0.         zxfqcalving(i) = 0.
733    
734         s_pblh(i) = 0.         s_pblh(i) = 0.
735         s_lcl(i) = 0.         s_lcl(i) = 0.
736         s_capCL(i) = 0.         s_capCL(i) = 0.
737         s_oliqCL(i) = 0.         s_oliqCL(i) = 0.
738         s_cteiCL(i) = 0.         s_cteiCL(i) = 0.
# Line 1005  contains Line 741  contains
741         s_trmb1(i) = 0.         s_trmb1(i) = 0.
742         s_trmb2(i) = 0.         s_trmb2(i) = 0.
743         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ème sous surface au point ', i, pctsrf(i, 1 : nbsrf)  
744      ENDDO      ENDDO
745    
746        call assert(abs(sum(pctsrf, dim = 2) - 1.) <= EPSFRA, 'physiq: pctsrf')
747    
748        ftsol = ftsol + d_ts
749        zxtsol = sum(ftsol * pctsrf, dim = 2)
750      DO nsrf = 1, nbsrf      DO nsrf = 1, nbsrf
751         DO i = 1, klon         DO i = 1, klon
752            ftsol(i, nsrf) = ftsol(i, nsrf) + d_ts(i, nsrf)            zxfluxlat(i) = zxfluxlat(i) + fluxlat(i, nsrf) * pctsrf(i, nsrf)
753            zxtsol(i) = zxtsol(i) + ftsol(i, nsrf)*pctsrf(i, nsrf)  
754            zxfluxlat(i) = zxfluxlat(i) + fluxlat(i, nsrf)*pctsrf(i, nsrf)            zt2m(i) = zt2m(i) + t2m(i, nsrf) * pctsrf(i, nsrf)
755              zq2m(i) = zq2m(i) + q2m(i, nsrf) * pctsrf(i, nsrf)
756            zt2m(i) = zt2m(i) + t2m(i, nsrf)*pctsrf(i, nsrf)            zu10m(i) = zu10m(i) + u10m(i, nsrf) * pctsrf(i, nsrf)
757            zq2m(i) = zq2m(i) + q2m(i, nsrf)*pctsrf(i, nsrf)            zv10m(i) = zv10m(i) + v10m(i, nsrf) * pctsrf(i, nsrf)
758            zu10m(i) = zu10m(i) + u10m(i, nsrf)*pctsrf(i, nsrf)            zxffonte(i) = zxffonte(i) + ffonte(i, nsrf) * pctsrf(i, nsrf)
           zv10m(i) = zv10m(i) + v10m(i, nsrf)*pctsrf(i, nsrf)  
           zxffonte(i) = zxffonte(i) + ffonte(i, nsrf)*pctsrf(i, nsrf)  
759            zxfqcalving(i) = zxfqcalving(i) + &            zxfqcalving(i) = zxfqcalving(i) + &
760                 fqcalving(i, nsrf)*pctsrf(i, nsrf)                 fqcalving(i, nsrf) * pctsrf(i, nsrf)
761            s_pblh(i) = s_pblh(i) + pblh(i, nsrf)*pctsrf(i, nsrf)            s_pblh(i) = s_pblh(i) + pblh(i, nsrf) * pctsrf(i, nsrf)
762            s_lcl(i) = s_lcl(i) + plcl(i, nsrf)*pctsrf(i, nsrf)            s_lcl(i) = s_lcl(i) + plcl(i, nsrf) * pctsrf(i, nsrf)
763            s_capCL(i) = s_capCL(i) + capCL(i, nsrf) *pctsrf(i, nsrf)            s_capCL(i) = s_capCL(i) + capCL(i, nsrf) * pctsrf(i, nsrf)
764            s_oliqCL(i) = s_oliqCL(i) + oliqCL(i, nsrf) *pctsrf(i, nsrf)            s_oliqCL(i) = s_oliqCL(i) + oliqCL(i, nsrf) * pctsrf(i, nsrf)
765            s_cteiCL(i) = s_cteiCL(i) + cteiCL(i, nsrf) *pctsrf(i, nsrf)            s_cteiCL(i) = s_cteiCL(i) + cteiCL(i, nsrf) * pctsrf(i, nsrf)
766            s_pblT(i) = s_pblT(i) + pblT(i, nsrf) *pctsrf(i, nsrf)            s_pblT(i) = s_pblT(i) + pblT(i, nsrf) * pctsrf(i, nsrf)
767            s_therm(i) = s_therm(i) + therm(i, nsrf) *pctsrf(i, nsrf)            s_therm(i) = s_therm(i) + therm(i, nsrf) * pctsrf(i, nsrf)
768            s_trmb1(i) = s_trmb1(i) + trmb1(i, nsrf) *pctsrf(i, nsrf)            s_trmb1(i) = s_trmb1(i) + trmb1(i, nsrf) * pctsrf(i, nsrf)
769            s_trmb2(i) = s_trmb2(i) + trmb2(i, nsrf) *pctsrf(i, nsrf)            s_trmb2(i) = s_trmb2(i) + trmb2(i, nsrf) * pctsrf(i, nsrf)
770            s_trmb3(i) = s_trmb3(i) + trmb3(i, nsrf) *pctsrf(i, nsrf)            s_trmb3(i) = s_trmb3(i) + trmb3(i, nsrf) * pctsrf(i, nsrf)
771         ENDDO         ENDDO
772      ENDDO      ENDDO
773    
774      ! Si une sous-fraction n'existe pas, elle prend la temp. moyenne      ! Si une sous-fraction n'existe pas, elle prend la température moyenne :
   
775      DO nsrf = 1, nbsrf      DO nsrf = 1, nbsrf
776         DO i = 1, klon         DO i = 1, klon
777            IF (pctsrf(i, nsrf) < epsfra) ftsol(i, nsrf) = zxtsol(i)            IF (pctsrf(i, nsrf) < epsfra) ftsol(i, nsrf) = zxtsol(i)
# Line 1062  contains Line 796  contains
796         ENDDO         ENDDO
797      ENDDO      ENDDO
798    
799      ! Calculer la derive du flux infrarouge      ! Calculer la dérive du flux infrarouge
800    
801      DO i = 1, klon      DO i = 1, klon
802         dlw(i) = - 4. * RSIGMA * zxtsol(i)**3         dlw(i) = - 4. * RSIGMA * zxtsol(i)**3
803      ENDDO      ENDDO
804    
805      ! Appeler la convection (au choix)      IF (check) print *, "avantcon = ", qcheck(paprs, q_seri, ql_seri)
806    
807      DO k = 1, llm      ! Appeler la convection
        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  
     ENDDO  
   
     IF (check) THEN  
        za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)  
        print *, "avantcon = ", za  
     ENDIF  
   
     if (iflag_con == 2) then  
        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  
   
        CALL concvl(dtphys, paprs, play, t_seri, q_seri, u_seri, &  
             v_seri, tr_seri, sig1, w01, d_t_con, d_q_con, &  
             d_u_con, d_v_con, d_tr, rain_con, snow_con, ibas_con, &  
             itop_con, upwd, dnwd, dnwd0, Ma, cape, tvp, iflagctrl, &  
             pbase, bbase, dtvpdt1, dtvpdq1, dplcldt, dplcldr, qcondc, &  
             wd, pmflxr, pmflxs, da, phi, mp, ntra=1)  
        ! (number of tracers for the convection scheme of Kerry Emanuel:  
        ! la partie traceurs est faite dans phytrac  
        ! on met ntra = 1 pour limiter les appels mais on peut  
        ! supprimer les calculs / ftra.)  
808    
809        if (conv_emanuel) then
810           da = 0.
811           mp = 0.
812           phi = 0.
813           CALL concvl(paprs, play, t_seri, q_seri, u_seri, v_seri, sig1, w01, &
814                d_t_con, d_q_con, d_u_con, d_v_con, rain_con, ibas_con, itop_con, &
815                upwd, dnwd, dnwd0, Ma, cape, iflagctrl, qcondc, pmflxr, da, phi, mp)
816           snow_con = 0.
817         clwcon0 = qcondc         clwcon0 = qcondc
818         mfu = upwd + dnwd         mfu = upwd + dnwd
        IF (.NOT. ok_gust) wd = 0.  
   
        ! Calcul des propriétés des nuages convectifs  
819    
820         DO k = 1, llm         IF (thermcep) THEN
821            DO i = 1, klon            zqsat = MIN(0.5, r2es * FOEEW(t_seri, rtt >= t_seri) / play)
822               zx_t = t_seri(i, k)            zqsat = zqsat / (1. - retv * zqsat)
823               IF (thermcep) THEN         ELSE
824                  zdelta = MAX(0., SIGN(1., rtt-zx_t))            zqsat = merge(qsats(t_seri), qsatl(t_seri), t_seri < t_coup) / play
825                  zx_qs = r2es * FOEEW(zx_t, zdelta) / play(i, k)         ENDIF
                 zx_qs = MIN(0.5, zx_qs)  
                 zcor = 1./(1.-retv*zx_qs)  
                 zx_qs = zx_qs*zcor  
              ELSE  
                 IF (zx_t < t_coup) THEN  
                    zx_qs = qsats(zx_t)/play(i, k)  
                 ELSE  
                    zx_qs = qsatl(zx_t)/play(i, k)  
                 ENDIF  
              ENDIF  
              zqsat(i, k) = zx_qs  
           ENDDO  
        ENDDO  
826    
827         ! calcul des proprietes des nuages convectifs         ! Properties of convective clouds
828         clwcon0 = fact_cldcon * clwcon0         clwcon0 = fact_cldcon * clwcon0
829         call clouds_gno(klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, &         call clouds_gno(klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, &
830              rnebcon0)              rnebcon0)
831    
832           forall (i = 1:klon) ema_pct(i) = paprs(i, itop_con(i) + 1)
833         mfd = 0.         mfd = 0.
834         pen_u = 0.         pen_u = 0.
835         pen_d = 0.         pen_d = 0.
836         pde_d = 0.         pde_d = 0.
837         pde_u = 0.         pde_u = 0.
838        else
839           conv_q = d_q_dyn + d_q_vdf / dtphys
840           conv_t = d_t_dyn + d_t_vdf / dtphys
841           z_avant = sum((q_seri + ql_seri) * zmasse, dim=2)
842           CALL conflx(dtphys, paprs, play, t_seri(:, llm:1:- 1), &
843                q_seri(:, llm:1:- 1), conv_t, conv_q, zxfluxq(:, 1), omega, &
844                d_t_con, d_q_con, rain_con, snow_con, mfu(:, llm:1:- 1), &
845                mfd(:, llm:1:- 1), pen_u, pde_u, pen_d, pde_d, kcbot, kctop, &
846                kdtop, pmflxr, pmflxs)
847           WHERE (rain_con < 0.) rain_con = 0.
848           WHERE (snow_con < 0.) snow_con = 0.
849           ibas_con = llm + 1 - kcbot
850           itop_con = llm + 1 - kctop
851      END if      END if
852    
853      DO k = 1, llm      DO k = 1, llm
# Line 1154  contains Line 859  contains
859         ENDDO         ENDDO
860      ENDDO      ENDDO
861    
862      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
863         tit = 'after convect'         tit = 'after convect'
864         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, &
865              ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
             d_ql, d_qs, d_ec)  
866         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, &
867              zero_v, zero_v, rain_con, snow_con, ztsol, d_h_vcol, d_qt, d_ec, &              zero_v, zero_v, rain_con, snow_con, ztsol, d_h_vcol, d_qt, d_ec)
             fs_bound, fq_bound)  
868      END IF      END IF
869    
870      IF (check) THEN      IF (check) THEN
871         za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)         za = qcheck(paprs, q_seri, ql_seri)
872         print *, "aprescon = ", za         print *, "aprescon = ", za
873         zx_t = 0.         zx_t = 0.
874         za = 0.         za = 0.
875         DO i = 1, klon         DO i = 1, klon
876            za = za + airephy(i)/REAL(klon)            za = za + airephy(i) / REAL(klon)
877            zx_t = zx_t + (rain_con(i)+ &            zx_t = zx_t + (rain_con(i)+ &
878                 snow_con(i))*airephy(i)/REAL(klon)                 snow_con(i)) * airephy(i) / REAL(klon)
879         ENDDO         ENDDO
880         zx_t = zx_t/za*dtphys         zx_t = zx_t / za * dtphys
881         print *, "Precip = ", zx_t         print *, "Precip = ", zx_t
882      ENDIF      ENDIF
883    
884      IF (iflag_con == 2) THEN      IF (.not. conv_emanuel) THEN
885         z_apres = sum((q_seri + ql_seri) * zmasse, dim=2)         z_apres = sum((q_seri + ql_seri) * zmasse, dim=2)
886         z_factor = (z_avant - (rain_con + snow_con) * dtphys) / z_apres         z_factor = (z_avant - (rain_con + snow_con) * dtphys) / z_apres
887         DO k = 1, llm         DO k = 1, llm
# Line 1190  contains Line 893  contains
893         ENDDO         ENDDO
894      ENDIF      ENDIF
895    
896      ! Convection sèche (thermiques ou ajustement)      ! Convection s\`eche (thermiques ou ajustement)
897    
898      d_t_ajs = 0.      d_t_ajs = 0.
899      d_u_ajs = 0.      d_u_ajs = 0.
# Line 1205  contains Line 908  contains
908         t_seri = t_seri + d_t_ajs         t_seri = t_seri + d_t_ajs
909         q_seri = q_seri + d_q_ajs         q_seri = q_seri + d_q_ajs
910      else      else
        ! Thermiques  
911         call calltherm(dtphys, play, paprs, pphi, u_seri, v_seri, t_seri, &         call calltherm(dtphys, play, paprs, pphi, u_seri, v_seri, t_seri, &
912              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)
913      endif      endif
914    
915      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
916         tit = 'after dry_adjust'         tit = 'after dry_adjust'
917         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, &
918              ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
             d_ql, d_qs, d_ec)  
919      END IF      END IF
920    
921      ! Caclul des ratqs      ! Caclul des ratqs
922    
923      ! ratqs convectifs à l'ancienne en fonction de (q(z = 0) - q) / q      ! ratqs convectifs \`a l'ancienne en fonction de (q(z = 0) - q) / q
924      ! on écrase le tableau ratqsc calculé par clouds_gno      ! on \'ecrase le tableau ratqsc calcul\'e par clouds_gno
925      if (iflag_cldcon == 1) then      if (iflag_cldcon == 1) then
926         do k = 1, llm         do k = 1, llm
927            do i = 1, klon            do i = 1, klon
# Line 1238  contains Line 939  contains
939      do k = 1, llm      do k = 1, llm
940         do i = 1, klon         do i = 1, klon
941            ratqss(i, k) = ratqsbas + (ratqshaut - ratqsbas) &            ratqss(i, k) = ratqsbas + (ratqshaut - ratqsbas) &
942                 * min((paprs(i, 1) - play(i, k)) / (paprs(i, 1) - 3e4), 1.)                 * min((paprs(i, 1) - play(i, k)) / (paprs(i, 1) - 3e4), 1.)
943         enddo         enddo
944      enddo      enddo
945    
# Line 1255  contains Line 956  contains
956         ratqs = ratqss         ratqs = ratqss
957      endif      endif
958    
     ! Processus de condensation à grande echelle et processus de  
     ! précipitation :  
959      CALL fisrtilp(dtphys, paprs, play, t_seri, q_seri, ptconv, ratqs, &      CALL fisrtilp(dtphys, paprs, play, t_seri, q_seri, ptconv, ratqs, &
960           d_t_lsc, d_q_lsc, d_ql_lsc, rneb, cldliq, rain_lsc, snow_lsc, &           d_t_lsc, d_q_lsc, d_ql_lsc, rneb, cldliq, rain_lsc, snow_lsc, &
961           pfrac_impa, pfrac_nucl, pfrac_1nucl, frac_impa, frac_nucl, prfl, &           pfrac_impa, pfrac_nucl, pfrac_1nucl, frac_impa, frac_nucl, prfl, &
# Line 1274  contains Line 973  contains
973         ENDDO         ENDDO
974      ENDDO      ENDDO
975      IF (check) THEN      IF (check) THEN
976         za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)         za = qcheck(paprs, q_seri, ql_seri)
977         print *, "apresilp = ", za         print *, "apresilp = ", za
978         zx_t = 0.         zx_t = 0.
979         za = 0.         za = 0.
980         DO i = 1, klon         DO i = 1, klon
981            za = za + airephy(i)/REAL(klon)            za = za + airephy(i) / REAL(klon)
982            zx_t = zx_t + (rain_lsc(i) &            zx_t = zx_t + (rain_lsc(i) &
983                 + snow_lsc(i))*airephy(i)/REAL(klon)                 + snow_lsc(i)) * airephy(i) / REAL(klon)
984         ENDDO         ENDDO
985         zx_t = zx_t/za*dtphys         zx_t = zx_t / za * dtphys
986         print *, "Precip = ", zx_t         print *, "Precip = ", zx_t
987      ENDIF      ENDIF
988    
989      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
990         tit = 'after fisrt'         tit = 'after fisrt'
991         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, &
992              ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
             d_ql, d_qs, d_ec)  
993         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, &
994              zero_v, zero_v, rain_lsc, snow_lsc, ztsol, d_h_vcol, d_qt, d_ec, &              zero_v, zero_v, rain_lsc, snow_lsc, ztsol, d_h_vcol, d_qt, d_ec)
             fs_bound, fq_bound)  
995      END IF      END IF
996    
997      ! PRESCRIPTION DES NUAGES POUR LE RAYONNEMENT      ! PRESCRIPTION DES NUAGES POUR LE RAYONNEMENT
998    
999      ! 1. NUAGES CONVECTIFS      ! 1. NUAGES CONVECTIFS
1000    
1001      IF (iflag_cldcon <= -1) THEN      IF (iflag_cldcon <= - 1) THEN
1002         ! seulement pour Tiedtke         ! seulement pour Tiedtke
1003         snow_tiedtke = 0.         snow_tiedtke = 0.
1004         if (iflag_cldcon == -1) then         if (iflag_cldcon == - 1) then
1005            rain_tiedtke = rain_con            rain_tiedtke = rain_con
1006         else         else
1007            rain_tiedtke = 0.            rain_tiedtke = 0.
1008            do k = 1, llm            do k = 1, llm
1009               do i = 1, klon               do i = 1, klon
1010                  if (d_q_con(i, k) < 0.) then                  if (d_q_con(i, k) < 0.) then
1011                     rain_tiedtke(i) = rain_tiedtke(i)-d_q_con(i, k)/dtphys &                     rain_tiedtke(i) = rain_tiedtke(i) - d_q_con(i, k) / dtphys &
1012                          *zmasse(i, k)                          * zmasse(i, k)
1013                  endif                  endif
1014               enddo               enddo
1015            enddo            enddo
# Line 1331  contains Line 1028  contains
1028         ENDDO         ENDDO
1029      ELSE IF (iflag_cldcon == 3) THEN      ELSE IF (iflag_cldcon == 3) THEN
1030         ! On prend pour les nuages convectifs le maximum du calcul de         ! On prend pour les nuages convectifs le maximum du calcul de
1031         ! la convection et du calcul du pas de temps précédent diminué         ! la convection et du calcul du pas de temps pr\'ec\'edent diminu\'e
1032         ! d'un facteur facttemps.         ! d'un facteur facttemps.
1033         facteur = dtphys * facttemps         facteur = dtphys * facttemps
1034         do k = 1, llm         do k = 1, llm
# Line 1347  contains Line 1044  contains
1044    
1045         ! On prend la somme des fractions nuageuses et des contenus en eau         ! On prend la somme des fractions nuageuses et des contenus en eau
1046         cldfra = min(max(cldfra, rnebcon), 1.)         cldfra = min(max(cldfra, rnebcon), 1.)
1047         cldliq = cldliq + rnebcon*clwcon         cldliq = cldliq + rnebcon * clwcon
1048      ENDIF      ENDIF
1049    
1050      ! 2. Nuages stratiformes      ! 2. Nuages stratiformes
# Line 1371  contains Line 1068  contains
1068      ENDDO      ENDDO
1069    
1070      IF (if_ebil >= 2) CALL diagetpq(airephy, "after diagcld", ip_ebil, 2, 2, &      IF (if_ebil >= 2) CALL diagetpq(airephy, "after diagcld", ip_ebil, 2, 2, &
1071           dtphys, t_seri, q_seri, ql_seri, qs_seri, u_seri, v_seri, paprs, &           dtphys, t_seri, q_seri, ql_seri, u_seri, v_seri, paprs, d_h_vcol, &
1072           d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)           d_qt, d_ec)
1073    
1074      ! Humidité relative pour diagnostic :      ! Humidit\'e relative pour diagnostic :
1075      DO k = 1, llm      DO k = 1, llm
1076         DO i = 1, klon         DO i = 1, klon
1077            zx_t = t_seri(i, k)            zx_t = t_seri(i, k)
1078            IF (thermcep) THEN            IF (thermcep) THEN
1079               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)  
1080               zx_qs = MIN(0.5, zx_qs)               zx_qs = MIN(0.5, zx_qs)
1081               zcor = 1./(1.-retv*zx_qs)               zcor = 1. / (1. - retv * zx_qs)
1082               zx_qs = zx_qs*zcor               zx_qs = zx_qs * zcor
1083            ELSE            ELSE
1084               IF (zx_t < t_coup) THEN               IF (zx_t < t_coup) THEN
1085                  zx_qs = qsats(zx_t)/play(i, k)                  zx_qs = qsats(zx_t) / play(i, k)
1086               ELSE               ELSE
1087                  zx_qs = qsatl(zx_t)/play(i, k)                  zx_qs = qsatl(zx_t) / play(i, k)
1088               ENDIF               ENDIF
1089            ENDIF            ENDIF
1090            zx_rh(i, k) = q_seri(i, k)/zx_qs            zx_rh(i, k) = q_seri(i, k) / zx_qs
1091            zqsat(i, k) = zx_qs            zqsat(i, k) = zx_qs
1092         ENDDO         ENDDO
1093      ENDDO      ENDDO
1094    
1095      ! Introduce the aerosol direct and first indirect radiative forcings:      ! Introduce the aerosol direct and first indirect radiative forcings:
1096      IF (ok_ade .OR. ok_aie) THEN      tau_ae = 0.
1097         ! Get sulfate aerosol distribution :      piz_ae = 0.
1098         CALL readsulfate(rdayvrai, firstcal, sulfate)      cg_ae = 0.
        CALL readsulfate_preind(rdayvrai, firstcal, sulfate_pi)  
1099    
1100         CALL aeropt(play, paprs, t_seri, sulfate, rhcl, tau_ae, piz_ae, cg_ae, &      ! Param\`etres optiques des nuages et quelques param\`etres pour
1101              aerindex)      ! diagnostics :
     ELSE  
        tau_ae = 0.  
        piz_ae = 0.  
        cg_ae = 0.  
     ENDIF  
   
     ! Paramètres optiques des nuages et quelques paramètres pour diagnostics :  
1102      if (ok_newmicro) then      if (ok_newmicro) then
1103         CALL newmicro(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, &         CALL newmicro(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, &
1104              cldh, cldl, cldm, cldt, cldq, flwp, fiwp, flwc, fiwc, ok_aie, &              cldh, cldl, cldm, cldt, cldq, flwp, fiwp, flwc, fiwc, ok_aie, &
# Line 1421  contains Line 1109  contains
1109              bl95_b1, cldtaupi, re, fl)              bl95_b1, cldtaupi, re, fl)
1110      endif      endif
1111    
1112      ! Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.      IF (MOD(itap - 1, radpas) == 0) THEN
1113      IF (MOD(itaprad, radpas) == 0) THEN         ! Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.
1114         DO i = 1, klon         ! Calcul de l'abedo moyen par maille
1115            albsol(i) = falbe(i, is_oce) * pctsrf(i, is_oce) &         albsol = sum(falbe * pctsrf, dim = 2)
1116                 + 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  
1117         ! Rayonnement (compatible Arpege-IFS) :         ! Rayonnement (compatible Arpege-IFS) :
1118         CALL radlwsw(dist, rmu0, fract, paprs, play, zxtsol, albsol, &         CALL radlwsw(dist, mu0, fract, paprs, play, zxtsol, albsol, t_seri, &
1119              albsollw, t_seri, q_seri, wo, cldfra, cldemi, cldtau, heat, &              q_seri, wo, cldfra, cldemi, cldtau, heat, heat0, cool, cool0, &
1120              heat0, cool, cool0, radsol, albpla, topsw, toplw, solsw, sollw, &              radsol, albpla, topsw, toplw, solsw, sollw, sollwdown, topsw0, &
1121              sollwdown, topsw0, toplw0, solsw0, sollw0, lwdn0, lwdn, lwup0, &              toplw0, solsw0, sollw0, lwdn0, lwdn, lwup0, lwup, swdn0, swdn, &
1122              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, &
1123              cg_ae, topswad, solswad, cldtaupi, topswai, solswai)              solswad, cldtaupi, topswai, solswai)
        itaprad = 0  
1124      ENDIF      ENDIF
     itaprad = itaprad + 1  
1125    
1126      ! Ajouter la tendance des rayonnements (tous les pas)      ! Ajouter la tendance des rayonnements (tous les pas)
1127    
1128      DO k = 1, llm      DO k = 1, llm
1129         DO i = 1, klon         DO i = 1, klon
1130            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 &
1131                   / 86400.
1132         ENDDO         ENDDO
1133      ENDDO      ENDDO
1134    
1135      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
1136         tit = 'after rad'         tit = 'after rad'
1137         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, &
1138              ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
             d_ql, d_qs, d_ec)  
1139         call diagphy(airephy, tit, ip_ebil, topsw, toplw, solsw, sollw, &         call diagphy(airephy, tit, ip_ebil, topsw, toplw, solsw, sollw, &
1140              zero_v, zero_v, zero_v, zero_v, ztsol, d_h_vcol, d_qt, d_ec, &              zero_v, zero_v, zero_v, zero_v, ztsol, d_h_vcol, d_qt, d_ec)
             fs_bound, fq_bound)  
1141      END IF      END IF
1142    
1143      ! Calculer l'hydrologie de la surface      ! Calculer l'hydrologie de la surface
# Line 1469  contains Line 1147  contains
1147      ENDDO      ENDDO
1148      DO nsrf = 1, nbsrf      DO nsrf = 1, nbsrf
1149         DO i = 1, klon         DO i = 1, klon
1150            zxqsurf(i) = zxqsurf(i) + fqsurf(i, nsrf)*pctsrf(i, nsrf)            zxqsurf(i) = zxqsurf(i) + fqsurf(i, nsrf) * pctsrf(i, nsrf)
1151            zxsnow(i) = zxsnow(i) + fsnow(i, nsrf)*pctsrf(i, nsrf)            zxsnow(i) = zxsnow(i) + fsnow(i, nsrf) * pctsrf(i, nsrf)
1152         ENDDO         ENDDO
1153      ENDDO      ENDDO
1154    
1155      ! Calculer le bilan du sol et la dérive de température (couplage)      ! Calculer le bilan du sol et la d\'erive de temp\'erature (couplage)
1156    
1157      DO i = 1, klon      DO i = 1, klon
1158         bils(i) = radsol(i) - sens(i) + zxfluxlat(i)         bils(i) = radsol(i) - sens(i) + zxfluxlat(i)
1159      ENDDO      ENDDO
1160    
1161      ! Paramétrisation de l'orographie à l'échelle sous-maille :      ! Param\'etrisation de l'orographie \`a l'\'echelle sous-maille :
1162    
1163      IF (ok_orodr) THEN      IF (ok_orodr) THEN
1164         ! selection des points pour lesquels le shema est actif:         ! S\'election des points pour lesquels le sch\'ema est actif :
1165         igwd = 0         igwd = 0
1166         DO i = 1, klon         DO i = 1, klon
1167            itest(i) = 0            itest(i) = 0
1168            IF (((zpic(i)-zmea(i)) > 100.).AND.(zstd(i) > 10.)) THEN            IF (zpic(i) - zmea(i) > 100. .AND. zstd(i) > 10.) THEN
1169               itest(i) = 1               itest(i) = 1
1170               igwd = igwd + 1               igwd = igwd + 1
              idx(igwd) = i  
1171            ENDIF            ENDIF
1172         ENDDO         ENDDO
1173    
1174         CALL drag_noro(klon, llm, dtphys, paprs, play, zmea, zstd, zsig, zgam, &         CALL drag_noro(klon, llm, dtphys, paprs, play, zmea, zstd, zsig, zgam, &
1175              zthe, zpic, zval, igwd, idx, itest, t_seri, u_seri, v_seri, &              zthe, zpic, zval, itest, t_seri, u_seri, v_seri, zulow, zvlow, &
1176              zulow, zvlow, zustrdr, zvstrdr, d_t_oro, d_u_oro, d_v_oro)              zustrdr, zvstrdr, d_t_oro, d_u_oro, d_v_oro)
1177    
1178         ! ajout des tendances         ! ajout des tendances
1179         DO k = 1, llm         DO k = 1, llm
# Line 1509  contains Line 1186  contains
1186      ENDIF      ENDIF
1187    
1188      IF (ok_orolf) THEN      IF (ok_orolf) THEN
1189         ! Sélection des points pour lesquels le schéma est actif :         ! S\'election des points pour lesquels le sch\'ema est actif :
1190         igwd = 0         igwd = 0
1191         DO i = 1, klon         DO i = 1, klon
1192            itest(i) = 0            itest(i) = 0
1193            IF ((zpic(i) - zmea(i)) > 100.) THEN            IF (zpic(i) - zmea(i) > 100.) THEN
1194               itest(i) = 1               itest(i) = 1
1195               igwd = igwd + 1               igwd = igwd + 1
              idx(igwd) = i  
1196            ENDIF            ENDIF
1197         ENDDO         ENDDO
1198    
# Line 1534  contains Line 1210  contains
1210         ENDDO         ENDDO
1211      ENDIF      ENDIF
1212    
1213      ! Stress nécessaires : toute la physique      ! Stress n\'ecessaires : toute la physique
1214    
1215      DO i = 1, klon      DO i = 1, klon
1216         zustrph(i) = 0.         zustrph(i) = 0.
# Line 1549  contains Line 1225  contains
1225         ENDDO         ENDDO
1226      ENDDO      ENDDO
1227    
1228      CALL aaam_bud(ra, rg, romega, rlat, rlon, pphis, zustrdr, zustrli, &      CALL aaam_bud(rg, romega, rlat, rlon, pphis, zustrdr, zustrli, zustrph, &
1229           zustrph, zvstrdr, zvstrli, zvstrph, paprs, u, v, aam, torsfc)           zvstrdr, zvstrli, zvstrph, paprs, u, v, aam, torsfc)
1230    
1231      IF (if_ebil >= 2) CALL diagetpq(airephy, 'after orography', ip_ebil, 2, &      IF (if_ebil >= 2) CALL diagetpq(airephy, 'after orography', ip_ebil, 2, &
1232           2, dtphys, t_seri, q_seri, ql_seri, qs_seri, u_seri, v_seri, paprs, &           2, dtphys, t_seri, q_seri, ql_seri, u_seri, v_seri, paprs, d_h_vcol, &
1233           d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)           d_qt, d_ec)
1234    
1235      ! Calcul des tendances traceurs      ! Calcul des tendances traceurs
1236      call phytrac(rnpb, itap, lmt_pas, julien, time, firstcal, lafin, nqmx-2, &      call phytrac(julien, time, firstcal, lafin, dtphys, t, paprs, play, mfu, &
1237           dtphys, u, t, paprs, play, mfu, mfd, pen_u, pde_u, pen_d, pde_d, &           mfd, pde_u, pen_d, ycoefh, fm_therm, entr_therm, yu1, yv1, ftsol, &
1238           ycoefh, fm_therm, entr_therm, yu1, yv1, ftsol, pctsrf, frac_impa, &           pctsrf, frac_impa, frac_nucl, da, phi, mp, upwd, dnwd, tr_seri, &
1239           frac_nucl, pphis, albsol, rhcl, cldfra, rneb, diafra, cldliq, &           zmasse, ncid_startphy)
1240           pmflxr, pmflxs, prfl, psfl, da, phi, mp, upwd, dnwd, tr_seri, zmasse)  
1241        IF (offline) call phystokenc(dtphys, t, mfu, mfd, pen_u, pde_u, pen_d, &
1242      IF (offline) THEN           pde_d, fm_therm, entr_therm, ycoefh, yu1, yv1, ftsol, pctsrf, &
1243         call phystokenc(dtphys, rlon, rlat, t, mfu, mfd, pen_u, pde_u, &           frac_impa, frac_nucl, pphis, airephy, dtphys)
             pen_d, pde_d, fm_therm, entr_therm, ycoefh, yu1, yv1, ftsol, &  
             pctsrf, frac_impa, frac_nucl, pphis, airephy, dtphys, itap)  
     ENDIF  
1244    
1245      ! Calculer le transport de l'eau et de l'energie (diagnostique)      ! Calculer le transport de l'eau et de l'energie (diagnostique)
1246      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)  
1247    
1248      ! diag. bilKP      ! diag. bilKP
1249    
1250      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, &
1251           ve_lay, vq_lay, ue_lay, uq_lay)           ve_lay, vq_lay, ue_lay, uq_lay)
1252    
1253      ! Accumuler les variables a stocker dans les fichiers histoire:      ! Accumuler les variables a stocker dans les fichiers histoire:
1254    
1255      ! conversion Ec -> E thermique      ! conversion Ec en énergie thermique
1256      DO k = 1, llm      DO k = 1, llm
1257         DO i = 1, klon         DO i = 1, klon
1258            ZRCPD = RCPD * (1. + RVTMP2 * q_seri(i, k))            ZRCPD = RCPD * (1. + RVTMP2 * q_seri(i, k))
# Line 1591  contains Line 1263  contains
1263         END DO         END DO
1264      END DO      END DO
1265    
1266      IF (if_ebil >= 1) THEN      IF (if_ebil >= 1) THEN
1267         tit = 'after physic'         tit = 'after physic'
1268         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, &
1269              ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &              ql_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_ec)
1270              d_ql, d_qs, d_ec)         ! Comme les tendances de la physique sont ajoute dans la dynamique,
        ! Comme les tendances de la physique sont ajoute dans la dynamique,  
1271         ! on devrait avoir que la variation d'entalpie par la dynamique         ! on devrait avoir que la variation d'entalpie par la dynamique
1272         ! 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.
1273         ! Donc la somme de ces 2 variations devrait etre nulle.         ! Donc la somme de ces 2 variations devrait etre nulle.
1274         call diagphy(airephy, tit, ip_ebil, topsw, toplw, solsw, sollw, sens, &         call diagphy(airephy, tit, ip_ebil, topsw, toplw, solsw, sollw, sens, &
1275              evap, rain_fall, snow_fall, ztsol, d_h_vcol, d_qt, d_ec, &              evap, rain_fall, snow_fall, ztsol, d_h_vcol, d_qt, d_ec)
             fs_bound, fq_bound)  
   
1276         d_h_vcol_phy = d_h_vcol         d_h_vcol_phy = d_h_vcol
   
1277      END IF      END IF
1278    
1279      ! SORTIES      ! SORTIES
# Line 1614  contains Line 1282  contains
1282      DO i = 1, klon      DO i = 1, klon
1283         prw(i) = 0.         prw(i) = 0.
1284         DO k = 1, llm         DO k = 1, llm
1285            prw(i) = prw(i) + q_seri(i, k)*zmasse(i, k)            prw(i) = prw(i) + q_seri(i, k) * zmasse(i, k)
1286         ENDDO         ENDDO
1287      ENDDO      ENDDO
1288    
# Line 1630  contains Line 1298  contains
1298         ENDDO         ENDDO
1299      ENDDO      ENDDO
1300    
1301      IF (nqmx >= 3) THEN      DO iq = 3, nqmx
1302         DO iq = 3, nqmx         DO k = 1, llm
1303            DO k = 1, llm            DO i = 1, klon
1304               DO i = 1, klon               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  
              ENDDO  
1305            ENDDO            ENDDO
1306         ENDDO         ENDDO
1307      ENDIF      ENDDO
1308    
1309      ! Sauvegarder les valeurs de t et q a la fin de la physique:      ! Sauvegarder les valeurs de t et q a la fin de la physique:
1310      DO k = 1, llm      DO k = 1, llm
# Line 1648  contains Line 1314  contains
1314         ENDDO         ENDDO
1315      ENDDO      ENDDO
1316    
1317      ! Ecriture des sorties      CALL histwrite_phy("phis", pphis)
1318      call write_histhf      CALL histwrite_phy("aire", airephy)
1319      call write_histday      CALL histwrite_phy("psol", paprs(:, 1))
1320      call write_histins      CALL histwrite_phy("precip", rain_fall + snow_fall)
1321        CALL histwrite_phy("plul", rain_lsc + snow_lsc)
1322      ! Si c'est la fin, il faut conserver l'etat de redemarrage      CALL histwrite_phy("pluc", rain_con + snow_con)
1323      IF (lafin) THEN      CALL histwrite_phy("tsol", zxtsol)
1324         itau_phy = itau_phy + itap      CALL histwrite_phy("t2m", zt2m)
1325         CALL phyredem("restartphy.nc", rlat, rlon, pctsrf, ftsol, ftsoil, &      CALL histwrite_phy("q2m", zq2m)
1326              tslab, seaice, fqsurf, qsol, fsnow, falbe, falblw, fevap, &      CALL histwrite_phy("u10m", zu10m)
1327              rain_fall, snow_fall, solsw, sollw, dlw, radsol, frugs, &      CALL histwrite_phy("v10m", zv10m)
1328              agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, t_ancien, &      CALL histwrite_phy("snow", snow_fall)
1329              q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, w01)      CALL histwrite_phy("cdrm", cdragm)
1330      ENDIF      CALL histwrite_phy("cdrh", cdragh)
1331        CALL histwrite_phy("topl", toplw)
1332      firstcal = .FALSE.      CALL histwrite_phy("evap", evap)
1333        CALL histwrite_phy("sols", solsw)
1334    contains      CALL histwrite_phy("soll", sollw)
1335        CALL histwrite_phy("solldown", sollwdown)
1336      subroutine write_histday      CALL histwrite_phy("bils", bils)
1337        CALL histwrite_phy("sens", - sens)
1338        use gr_phy_write_3d_m, only: gr_phy_write_3d      CALL histwrite_phy("fder", fder)
1339        integer itau_w ! pas de temps ecriture      CALL histwrite_phy("dtsvdfo", d_ts(:, is_oce))
1340        CALL histwrite_phy("dtsvdft", d_ts(:, is_ter))
1341        !------------------------------------------------      CALL histwrite_phy("dtsvdfg", d_ts(:, is_lic))
1342        CALL histwrite_phy("dtsvdfi", d_ts(:, is_sic))
       if (ok_journe) THEN  
          itau_w = itau_phy + itap  
          if (nqmx <= 4) then  
             call histwrite(nid_day, "Sigma_O3_Royer", itau_w, &  
                  gr_phy_write_3d(wo) * 1e3)  
             ! (convert "wo" from kDU to DU)  
          end if  
          if (ok_sync) then  
             call histsync(nid_day)  
          endif  
       ENDIF  
   
     End subroutine write_histday  
   
     !****************************  
   
     subroutine write_histhf  
   
       ! From phylmd/write_histhf.h, version 1.5 2005/05/25 13:10:09  
   
       !------------------------------------------------  
   
       call write_histhf3d  
   
       IF (ok_sync) THEN  
          call histsync(nid_hf)  
       ENDIF  
   
     end subroutine write_histhf  
   
     !***************************************************************  
   
     subroutine write_histins  
   
       ! From phylmd/write_histins.h, version 1.2 2005/05/25 13:10:09  
   
       real zout  
       integer itau_w ! pas de temps ecriture  
   
       !--------------------------------------------------  
   
       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)  
1343    
1344           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, toplw, zx_tmp_2d)      DO nsrf = 1, nbsrf
1345           CALL histwrite(nid_ins, "topl", itau_w, zx_tmp_2d)         CALL histwrite_phy("pourc_"//clnsurf(nsrf), pctsrf(:, nsrf) * 100.)
1346           CALL histwrite_phy("fract_"//clnsurf(nsrf), pctsrf(:, nsrf))
1347           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, evap, zx_tmp_2d)         CALL histwrite_phy("sens_"//clnsurf(nsrf), fluxt(:, 1, nsrf))
1348           CALL histwrite(nid_ins, "evap", itau_w, zx_tmp_2d)         CALL histwrite_phy("lat_"//clnsurf(nsrf), fluxlat(:, nsrf))
1349           CALL histwrite_phy("tsol_"//clnsurf(nsrf), ftsol(:, nsrf))
1350           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, solsw, zx_tmp_2d)         CALL histwrite_phy("taux_"//clnsurf(nsrf), fluxu(:, 1, nsrf))
1351           CALL histwrite(nid_ins, "sols", itau_w, zx_tmp_2d)         CALL histwrite_phy("tauy_"//clnsurf(nsrf), fluxv(:, 1, nsrf))
1352           CALL histwrite_phy("rugs_"//clnsurf(nsrf), frugs(:, nsrf))
1353           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sollw, zx_tmp_2d)         CALL histwrite_phy("albe_"//clnsurf(nsrf), falbe(:, nsrf))
1354           CALL histwrite(nid_ins, "soll", itau_w, zx_tmp_2d)      END DO
   
          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)  
   
          if (ok_sync) then  
             call histsync(nid_ins)  
          endif  
       ENDIF  
   
     end subroutine write_histins  
   
     !****************************************************  
   
     subroutine write_histhf3d  
   
       ! From phylmd/write_histhf3d.h, version 1.2 2005/05/25 13:10:09  
   
       integer itau_w ! pas de temps ecriture  
   
       !-------------------------------------------------------  
   
       itau_w = itau_phy + itap  
   
       ! Champs 3D:  
   
       CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, t_seri, zx_tmp_3d)  
       CALL histwrite(nid_hf3d, "temp", itau_w, zx_tmp_3d)  
   
       CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, qx(1, 1, ivap), zx_tmp_3d)  
       CALL histwrite(nid_hf3d, "ovap", itau_w, zx_tmp_3d)  
   
       CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, u_seri, zx_tmp_3d)  
       CALL histwrite(nid_hf3d, "vitu", itau_w, zx_tmp_3d)  
   
       CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, v_seri, zx_tmp_3d)  
       CALL histwrite(nid_hf3d, "vitv", itau_w, zx_tmp_3d)  
   
       if (nbtr >= 3) then  
          CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, tr_seri(1, 1, 3), &  
               zx_tmp_3d)  
          CALL histwrite(nid_hf3d, "O3", itau_w, zx_tmp_3d)  
       end if  
1355    
1356        if (ok_sync) then      CALL histwrite_phy("albs", albsol)
1357           call histsync(nid_hf3d)      CALL histwrite_phy("rugs", zxrugs)
1358        endif      CALL histwrite_phy("s_pblh", s_pblh)
1359        CALL histwrite_phy("s_pblt", s_pblt)
1360        CALL histwrite_phy("s_lcl", s_lcl)
1361        CALL histwrite_phy("s_capCL", s_capCL)
1362        CALL histwrite_phy("s_oliqCL", s_oliqCL)
1363        CALL histwrite_phy("s_cteiCL", s_cteiCL)
1364        CALL histwrite_phy("s_therm", s_therm)
1365        CALL histwrite_phy("s_trmb1", s_trmb1)
1366        CALL histwrite_phy("s_trmb2", s_trmb2)
1367        CALL histwrite_phy("s_trmb3", s_trmb3)
1368        if (conv_emanuel) CALL histwrite_phy("ptop", ema_pct)
1369        CALL histwrite_phy("temp", t_seri)
1370        CALL histwrite_phy("vitu", u_seri)
1371        CALL histwrite_phy("vitv", v_seri)
1372        CALL histwrite_phy("geop", zphi)
1373        CALL histwrite_phy("pres", play)
1374        CALL histwrite_phy("dtvdf", d_t_vdf)
1375        CALL histwrite_phy("dqvdf", d_q_vdf)
1376        CALL histwrite_phy("rhum", zx_rh)
1377    
1378        if (ok_instan) call histsync(nid_ins)
1379    
1380        IF (lafin) then
1381           call NF95_CLOSE(ncid_startphy)
1382           CALL phyredem(pctsrf, ftsol, ftsoil, fqsurf, qsol, &
1383                fsnow, falbe, fevap, rain_fall, snow_fall, solsw, sollw, dlw, &
1384                radsol, frugs, agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, &
1385                t_ancien, q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, &
1386                w01)
1387        end IF
1388    
1389      end subroutine write_histhf3d      firstcal = .FALSE.
1390    
1391    END SUBROUTINE physiq    END SUBROUTINE physiq
1392    

Legend:
Removed from v.72  
changed lines
  Added in v.202

  ViewVC Help
Powered by ViewVC 1.1.21