/[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 68 by guez, Wed Nov 14 16:59:30 2012 UTC trunk/Sources/phylmd/physiq.f revision 186 by guez, Mon Mar 21 15:36:26 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 (SVN revision 678)      ! From phylmd/physiq.F, version 1.22 2006/02/20 09:38:28
11      ! Author: Z.X. Li (LMD/CNRS) 1993      ! (subversion revision 678)
12    
13        ! 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    
# Line 16  contains Line 18  contains
18      USE abort_gcm_m, ONLY: abort_gcm      USE abort_gcm_m, ONLY: abort_gcm
19      use aeropt_m, only: aeropt      use aeropt_m, only: aeropt
20      use ajsec_m, only: ajsec      use ajsec_m, only: ajsec
     USE calendar, ONLY: ymds2ju  
21      use calltherm_m, only: calltherm      use calltherm_m, only: calltherm
22      USE clesphys, ONLY: cdhmax, cdmmax, co2_ppm, ecrit_hf, ecrit_ins, &      USE clesphys, ONLY: cdhmax, cdmmax, ecrit_hf, ecrit_ins, ecrit_mth, &
23           ecrit_mth, ecrit_reg, ecrit_tra, ksta, ksta_ter, ok_kzmin           ecrit_reg, ecrit_tra, ksta, ksta_ter, ok_kzmin
24      USE clesphys2, ONLY: cycle_diurne, iflag_con, nbapp_rad, new_oliq, &      USE clesphys2, ONLY: cycle_diurne, conv_emanuel, nbapp_rad, new_oliq, &
25           ok_orodr, ok_orolf, soil_model           ok_orodr, ok_orolf
26      USE clmain_m, ONLY: clmain      USE clmain_m, ONLY: clmain
27      USE comgeomphy, ONLY: airephy, cuphy, cvphy      use clouds_gno_m, only: clouds_gno
28        use comconst, only: dtphys
29        USE comgeomphy, ONLY: airephy
30      USE concvl_m, ONLY: concvl      USE concvl_m, ONLY: concvl
31      USE conf_gcm_m, ONLY: offline, raz_date      USE conf_gcm_m, ONLY: offline, raz_date, day_step, iphysiq
32      USE conf_phys_m, ONLY: conf_phys      USE conf_phys_m, ONLY: conf_phys
33      use conflx_m, only: conflx      use conflx_m, only: conflx
34      USE ctherm, ONLY: iflag_thermals, nsplit_thermals      USE ctherm, ONLY: iflag_thermals, nsplit_thermals
35      use diagcld2_m, only: diagcld2      use diagcld2_m, only: diagcld2
36      use diagetpq_m, only: diagetpq      use diagetpq_m, only: diagetpq
37      use diagphy_m, only: diagphy      use diagphy_m, only: diagphy
38      USE dimens_m, ONLY: iim, jjm, llm, nqmx      USE dimens_m, ONLY: llm, nqmx
39      USE dimphy, ONLY: klon, nbtr      USE dimphy, ONLY: klon
40      USE dimsoil, ONLY: nsoilmx      USE dimsoil, ONLY: nsoilmx
41      use drag_noro_m, only: drag_noro      use drag_noro_m, only: drag_noro
42        use dynetat0_m, only: day_ref, annee_ref
43      USE fcttre, ONLY: foeew, qsatl, qsats, thermcep      USE fcttre, ONLY: foeew, qsatl, qsats, thermcep
44      use fisrtilp_m, only: fisrtilp      use fisrtilp_m, only: fisrtilp
45      USE hgardfou_m, ONLY: hgardfou      USE hgardfou_m, ONLY: hgardfou
     USE histsync_m, ONLY: histsync  
     USE histwrite_m, ONLY: histwrite  
46      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, &
47           nbsrf           nbsrf
     USE ini_histhf_m, ONLY: ini_histhf  
     USE ini_histday_m, ONLY: ini_histday  
48      USE ini_histins_m, ONLY: ini_histins      USE ini_histins_m, ONLY: ini_histins
49        use netcdf95, only: NF95_CLOSE
50      use newmicro_m, only: newmicro      use newmicro_m, only: newmicro
51      USE oasis_m, ONLY: ok_oasis      use nuage_m, only: nuage
52      USE orbite_m, ONLY: orbite, zenang      USE orbite_m, ONLY: orbite
53      USE ozonecm_m, ONLY: ozonecm      USE ozonecm_m, ONLY: ozonecm
54      USE phyetat0_m, ONLY: phyetat0, rlat, rlon      USE phyetat0_m, ONLY: phyetat0, rlat, rlon
55      USE phyredem_m, ONLY: phyredem      USE phyredem_m, ONLY: phyredem
56        USE phyredem0_m, ONLY: phyredem0
57      USE phystokenc_m, ONLY: phystokenc      USE phystokenc_m, ONLY: phystokenc
58      USE phytrac_m, ONLY: phytrac      USE phytrac_m, ONLY: phytrac
59      USE qcheck_m, ONLY: qcheck      USE qcheck_m, ONLY: qcheck
60      use radlwsw_m, only: radlwsw      use radlwsw_m, only: radlwsw
61      use sugwd_m, only: sugwd      use readsulfate_m, only: readsulfate
62      USE suphec_m, ONLY: ra, rcpd, retv, rg, rlvtt, romega, rsigma, rtt      use readsulfate_preind_m, only: readsulfate_preind
63      USE temps, ONLY: annee_ref, day_ref, itau_phy      use yoegwd, only: sugwd
64        USE suphec_m, ONLY: rcpd, retv, rg, rlvtt, romega, rsigma, rtt
65        use transp_m, only: transp
66        use transp_lay_m, only: transp_lay
67      use unit_nml_m, only: unit_nml      use unit_nml_m, only: unit_nml
68        USE ymds2ju_m, ONLY: ymds2ju
69      USE yoethf_m, ONLY: r2es, rvtmp2      USE yoethf_m, ONLY: r2es, rvtmp2
70        use zenang_m, only: zenang
71    
72      ! Arguments:      logical, intent(in):: lafin ! dernier passage
73    
74      REAL, intent(in):: rdayvrai      integer, intent(in):: dayvrai
75      ! (elapsed time since January 1st 0h of the starting year, in days)      ! current day number, based at value 1 on January 1st of annee_ref
76    
77      REAL, intent(in):: time ! heure de la journée en fraction de jour      REAL, intent(in):: time ! heure de la journ\'ee en fraction de jour
     REAL, intent(in):: dtphys ! pas d'integration pour la physique (seconde)  
     logical, intent(in):: lafin ! dernier passage  
78    
79      REAL, intent(in):: paprs(klon, llm + 1)      REAL, intent(in):: paprs(:, :) ! (klon, llm + 1)
80      ! (pression pour chaque inter-couche, en Pa)      ! pression pour chaque inter-couche, en Pa
81    
82      REAL, intent(in):: play(klon, llm)      REAL, intent(in):: play(:, :) ! (klon, llm)
83      ! (input pression pour le mileu de chaque couche (en Pa))      ! pression pour le mileu de chaque couche (en Pa)
84    
85      REAL, intent(in):: pphi(klon, llm)      REAL, intent(in):: pphi(:, :) ! (klon, llm)
86      ! (input geopotentiel de chaque couche (g z) (reference sol))      ! géopotentiel de chaque couche (référence sol)
87    
88      REAL, intent(in):: pphis(klon) ! input geopotentiel du sol      REAL, intent(in):: pphis(:) ! (klon) géopotentiel du sol
89    
90      REAL, intent(in):: u(klon, llm)      REAL, intent(in):: u(:, :) ! (klon, llm)
91      ! vitesse dans la direction X (de O a E) en m/s      ! vitesse dans la direction X (de O a E) en m/s
92    
93      REAL, intent(in):: v(klon, llm) ! vitesse Y (de S a N) en m/s      REAL, intent(in):: v(:, :) ! (klon, llm) vitesse Y (de S a N) en m/s
94      REAL, intent(in):: t(klon, llm) ! input temperature (K)      REAL, intent(in):: t(:, :) ! (klon, llm) temperature (K)
95    
96      REAL, intent(in):: qx(klon, llm, nqmx)      REAL, intent(in):: qx(:, :, :) ! (klon, llm, nqmx)
97      ! (humidité spécifique et fractions massiques des autres traceurs)      ! (humidit\'e sp\'ecifique et fractions massiques des autres traceurs)
98    
99      REAL omega(klon, llm) ! input vitesse verticale en Pa/s      REAL, intent(in):: omega(:, :) ! (klon, llm) vitesse verticale en Pa/s
100      REAL, intent(out):: d_u(klon, llm) ! tendance physique de "u" (m/s/s)      REAL, intent(out):: d_u(:, :) ! (klon, llm) tendance physique de "u" (m s-2)
101      REAL, intent(out):: d_v(klon, llm) ! tendance physique de "v" (m/s/s)      REAL, intent(out):: d_v(:, :) ! (klon, llm) tendance physique de "v" (m s-2)
102      REAL, intent(out):: d_t(klon, llm) ! tendance physique de "t" (K/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  
103    
104      LOGICAL:: firstcal = .true.      REAL, intent(out):: d_qx(:, :, :) ! (klon, llm, nqmx)
105        ! tendance physique de "qx" (s-1)
106    
107      INTEGER nbteta      ! Local:
     PARAMETER(nbteta = 3)  
108    
109      REAL PVteta(klon, nbteta)      LOGICAL:: firstcal = .true.
     ! (output vorticite potentielle a des thetas constantes)  
110    
111      LOGICAL ok_gust ! pour activer l'effet des gust sur flux surface      LOGICAL ok_gust ! pour activer l'effet des gust sur flux surface
112      PARAMETER (ok_gust = .FALSE.)      PARAMETER (ok_gust = .FALSE.)
113    
114      LOGICAL check ! Verifier la conservation du modele en eau      LOGICAL, PARAMETER:: check = .FALSE.
115      PARAMETER (check = .FALSE.)      ! Verifier la conservation du modele en eau
116    
117      LOGICAL, PARAMETER:: ok_stratus = .FALSE.      LOGICAL, PARAMETER:: ok_stratus = .FALSE.
118      ! Ajouter artificiellement les stratus      ! Ajouter artificiellement les stratus
119    
     ! Parametres lies au coupleur OASIS:  
     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")  
   
     logical ok_ocean  
     SAVE ok_ocean  
   
     ! "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  
   
120      logical:: ok_journe = .false., ok_mensuel = .true., ok_instan = .false.      logical:: ok_journe = .false., ok_mensuel = .true., ok_instan = .false.
121      ! sorties journalieres, mensuelles et instantanees dans les      ! sorties journalieres, mensuelles et instantanees dans les
122      ! fichiers histday, histmth et histins      ! fichiers histday, histmth et histins
# Line 147  contains Line 129  contains
129      REAL entr_therm(klon, llm)      REAL entr_therm(klon, llm)
130      real, save:: q2(klon, llm + 1, nbsrf)      real, save:: q2(klon, llm + 1, nbsrf)
131    
132      INTEGER ivap ! indice de traceurs pour vapeur d'eau      INTEGER, PARAMETER:: ivap = 1 ! indice de traceur pour vapeur d'eau
133      PARAMETER (ivap = 1)      INTEGER, PARAMETER:: iliq = 2 ! indice de traceur pour eau liquide
     INTEGER iliq ! indice de traceurs pour eau liquide  
     PARAMETER (iliq = 2)  
134    
135      REAL, save:: t_ancien(klon, llm), q_ancien(klon, llm)      REAL, save:: t_ancien(klon, llm), q_ancien(klon, llm)
136      LOGICAL, save:: ancien_ok      LOGICAL, save:: ancien_ok
# Line 160  contains Line 140  contains
140    
141      real da(klon, llm), phi(klon, llm, llm), mp(klon, llm)      real da(klon, llm), phi(klon, llm, llm), mp(klon, llm)
142    
143      !IM Amip2 PV a theta constante      REAL swdn0(klon, llm + 1), swdn(klon, llm + 1)
144        REAL swup0(klon, llm + 1), swup(klon, llm + 1)
     CHARACTER(LEN = 3) ctetaSTD(nbteta)  
     DATA ctetaSTD/'350', '380', '405'/  
     REAL rtetaSTD(nbteta)  
     DATA rtetaSTD/350., 380., 405./  
   
     !MI Amip2 PV a theta constante  
   
     INTEGER klevp1  
     PARAMETER(klevp1 = llm + 1)  
   
     REAL swdn0(klon, klevp1), swdn(klon, klevp1)  
     REAL swup0(klon, klevp1), swup(klon, klevp1)  
145      SAVE swdn0, swdn, swup0, swup      SAVE swdn0, swdn, swup0, swup
146    
147      REAL lwdn0(klon, klevp1), lwdn(klon, klevp1)      REAL lwdn0(klon, llm + 1), lwdn(klon, llm + 1)
148      REAL lwup0(klon, klevp1), lwup(klon, klevp1)      REAL lwup0(klon, llm + 1), lwup(klon, llm + 1)
149      SAVE lwdn0, lwdn, lwup0, lwup      SAVE lwdn0, lwdn, lwup0, lwup
150    
151      !IM Amip2      ! Amip2
152      ! variables a une pression donnee      ! variables a une pression donnee
153    
154      integer nlevSTD      integer nlevSTD
155      PARAMETER(nlevSTD = 17)      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 '/  
156    
157      ! prw: precipitable water      ! prw: precipitable water
158      real prw(klon)      real prw(klon)
# Line 205  contains Line 165  contains
165      INTEGER kmax, lmax      INTEGER kmax, lmax
166      PARAMETER(kmax = 8, lmax = 8)      PARAMETER(kmax = 8, lmax = 8)
167      INTEGER kmaxm1, lmaxm1      INTEGER kmaxm1, lmaxm1
168      PARAMETER(kmaxm1 = kmax-1, lmaxm1 = lmax-1)      PARAMETER(kmaxm1 = kmax - 1, lmaxm1 = lmax - 1)
   
     REAL zx_tau(kmaxm1), zx_pc(lmaxm1)  
     DATA zx_tau/0.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  
169    
170      ! Variables propres a la physique      ! Variables propres a la physique
171    
172      INTEGER, save:: radpas      INTEGER, save:: radpas
173      ! (Radiative transfer computations are made every "radpas" call to      ! Radiative transfer computations are made every "radpas" call to
174      ! "physiq".)      ! "physiq".
175    
176      REAL radsol(klon)      REAL radsol(klon)
177      SAVE radsol ! bilan radiatif au sol calcule par code radiatif      SAVE radsol ! bilan radiatif au sol calcule par code radiatif
178    
179      INTEGER, SAVE:: itap ! number of calls to "physiq"      INTEGER:: itap = 0 ! number of calls to "physiq"
180    
181      REAL, save:: ftsol(klon, nbsrf) ! skin temperature of surface fraction      REAL, save:: ftsol(klon, nbsrf) ! skin temperature of surface fraction
182    
183      REAL, save:: ftsoil(klon, nsoilmx, nbsrf)      REAL, save:: ftsoil(klon, nsoilmx, nbsrf)
184      ! soil temperature of surface fraction      ! soil temperature of surface fraction
185    
186      REAL fevap(klon, nbsrf)      REAL, save:: fevap(klon, nbsrf) ! evaporation
     SAVE fevap ! evaporation  
187      REAL fluxlat(klon, nbsrf)      REAL fluxlat(klon, nbsrf)
188      SAVE fluxlat      SAVE fluxlat
189    
190      REAL fqsurf(klon, nbsrf)      REAL, save:: fqsurf(klon, nbsrf)
191      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  
192    
193      REAL fsnow(klon, nbsrf)      REAL, save:: qsol(klon)
194      SAVE fsnow ! epaisseur neigeuse      ! column-density of water in soil, in kg m-2
195    
196      REAL falbe(klon, nbsrf)      REAL, save:: fsnow(klon, nbsrf) ! epaisseur neigeuse
197      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  
198    
199      ! Paramètres de l'orographie à l'échelle sous-maille (OESM) :      ! Param\`etres de l'orographie \`a l'\'echelle sous-maille (OESM) :
200      REAL, save:: zmea(klon) ! orographie moyenne      REAL, save:: zmea(klon) ! orographie moyenne
201      REAL, save:: zstd(klon) ! deviation standard de l'OESM      REAL, save:: zstd(klon) ! deviation standard de l'OESM
202      REAL, save:: zsig(klon) ! pente de l'OESM      REAL, save:: zsig(klon) ! pente de l'OESM
# Line 297  contains Line 205  contains
205      REAL, save:: zpic(klon) ! Maximum de l'OESM      REAL, save:: zpic(klon) ! Maximum de l'OESM
206      REAL, save:: zval(klon) ! Minimum de l'OESM      REAL, save:: zval(klon) ! Minimum de l'OESM
207      REAL, save:: rugoro(klon) ! longueur de rugosite de l'OESM      REAL, save:: rugoro(klon) ! longueur de rugosite de l'OESM
   
208      REAL zulow(klon), zvlow(klon)      REAL zulow(klon), zvlow(klon)
209        INTEGER igwd, itest(klon)
     INTEGER igwd, idx(klon), itest(klon)  
210    
211      REAL agesno(klon, nbsrf)      REAL agesno(klon, nbsrf)
212      SAVE agesno ! age de la neige      SAVE agesno ! age de la neige
# Line 310  contains Line 216  contains
216      !KE43      !KE43
217      ! Variables liees a la convection de K. Emanuel (sb):      ! Variables liees a la convection de K. Emanuel (sb):
218    
     REAL bas, top ! cloud base and top levels  
     SAVE bas  
     SAVE top  
   
219      REAL Ma(klon, llm) ! undilute upward mass flux      REAL Ma(klon, llm) ! undilute upward mass flux
220      SAVE Ma      SAVE Ma
221      REAL qcondc(klon, llm) ! in-cld water content from convect      REAL qcondc(klon, llm) ! in-cld water content from convect
222      SAVE qcondc      SAVE qcondc
223      REAL ema_work1(klon, llm), ema_work2(klon, llm)      REAL, save:: sig1(klon, llm), w01(klon, llm)
224      SAVE ema_work1, ema_work2      REAL, save:: wd(klon)
225    
226      REAL wd(klon) ! sb      ! Variables pour la couche limite (al1):
     SAVE wd ! sb  
   
     ! Variables locales pour la couche limite (al1):  
   
     ! Variables locales:  
227    
228      REAL cdragh(klon) ! drag coefficient pour T and Q      REAL cdragh(klon) ! drag coefficient pour T and Q
229      REAL cdragm(klon) ! drag coefficient pour vent      REAL cdragm(klon) ! drag coefficient pour vent
230    
231      !AA Pour phytrac      ! Pour phytrac :
232      REAL ycoefh(klon, llm) ! coef d'echange pour phytrac      REAL ycoefh(klon, llm) ! coef d'echange pour phytrac
233      REAL yu1(klon) ! vents dans la premiere couche U      REAL yu1(klon) ! vents dans la premiere couche U
234      REAL yv1(klon) ! vents dans la premiere couche V      REAL yv1(klon) ! vents dans la premiere couche V
# Line 350  contains Line 247  contains
247      REAL frac_impa(klon, llm) ! fractions d'aerosols lessivees (impaction)      REAL frac_impa(klon, llm) ! fractions d'aerosols lessivees (impaction)
248      REAL frac_nucl(klon, llm) ! idem (nucleation)      REAL frac_nucl(klon, llm) ! idem (nucleation)
249    
250      REAL, save:: rain_fall(klon) ! pluie      REAL, save:: rain_fall(klon)
251      REAL, save:: snow_fall(klon) ! neige      ! liquid water mass flux (kg/m2/s), positive down
252    
253        REAL, save:: snow_fall(klon)
254        ! solid water mass flux (kg/m2/s), positive down
255    
256      REAL rain_tiedtke(klon), snow_tiedtke(klon)      REAL rain_tiedtke(klon), snow_tiedtke(klon)
257    
258      REAL evap(klon), devap(klon) ! evaporation et sa derivee      REAL evap(klon), devap(klon) ! evaporation and its derivative
259      REAL sens(klon), dsens(klon) ! chaleur sensible et sa derivee      REAL sens(klon), dsens(klon) ! chaleur sensible et sa derivee
260      REAL dlw(klon) ! derivee infra rouge      REAL dlw(klon) ! derivee infra rouge
261      SAVE dlw      SAVE dlw
262      REAL bils(klon) ! bilan de chaleur au sol      REAL bils(klon) ! bilan de chaleur au sol
263      REAL fder(klon) ! Derive de flux (sensible et latente)      REAL, save:: fder(klon) ! Derive de flux (sensible et latente)
     save fder  
264      REAL ve(klon) ! integr. verticale du transport meri. de l'energie      REAL ve(klon) ! integr. verticale du transport meri. de l'energie
265      REAL vq(klon) ! integr. verticale du transport meri. de l'eau      REAL vq(klon) ! integr. verticale du transport meri. de l'eau
266      REAL ue(klon) ! integr. verticale du transport zonal de l'energie      REAL ue(klon) ! integr. verticale du transport zonal de l'energie
267      REAL uq(klon) ! integr. verticale du transport zonal de l'eau      REAL uq(klon) ! integr. verticale du transport zonal de l'eau
268    
269      REAL frugs(klon, nbsrf) ! longueur de rugosite      REAL, save:: frugs(klon, nbsrf) ! longueur de rugosite
     save frugs  
270      REAL zxrugs(klon) ! longueur de rugosite      REAL zxrugs(klon) ! longueur de rugosite
271    
272      ! Conditions aux limites      ! Conditions aux limites
273    
274      INTEGER julien      INTEGER julien
   
275      INTEGER, SAVE:: lmt_pas ! number of time steps of "physics" per day      INTEGER, SAVE:: lmt_pas ! number of time steps of "physics" per day
276      REAL pctsrf(klon, nbsrf)      REAL, save:: pctsrf(klon, nbsrf) ! percentage of surface
277      !IM      REAL pctsrf_new(klon, nbsrf) ! pourcentage surfaces issus d'ORCHIDEE
278      REAL pctsrf_new(klon, nbsrf) !pourcentage surfaces issus d'ORCHIDEE      REAL, save:: albsol(klon) ! albedo du sol total visible
   
     SAVE pctsrf ! sous-fraction du sol  
     REAL albsol(klon)  
     SAVE albsol ! albedo du sol total  
     REAL albsollw(klon)  
     SAVE albsollw ! albedo du sol total  
   
279      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
280    
281      ! Declaration des procedures appelees      real, save:: clwcon(klon, llm), rnebcon(klon, llm)
282        real, save:: clwcon0(klon, llm), rnebcon0(klon, llm)
     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  
   
     real clwcon(klon, llm), rnebcon(klon, llm)  
     real clwcon0(klon, llm), rnebcon0(klon, llm)  
   
     save rnebcon, clwcon  
283    
284      REAL rhcl(klon, llm) ! humiditi relative ciel clair      REAL rhcl(klon, llm) ! humiditi relative ciel clair
285      REAL dialiq(klon, llm) ! eau liquide nuageuse      REAL dialiq(klon, llm) ! eau liquide nuageuse
# Line 421  contains Line 299  contains
299      REAL zxfluxu(klon, llm)      REAL zxfluxu(klon, llm)
300      REAL zxfluxv(klon, llm)      REAL zxfluxv(klon, llm)
301    
302      ! 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
303      ! les variables soient rémanentes.      ! les variables soient r\'emanentes.
304      REAL, save:: heat(klon, llm) ! chauffage solaire      REAL, save:: heat(klon, llm) ! chauffage solaire
305      REAL heat0(klon, llm) ! chauffage solaire ciel clair      REAL, save:: heat0(klon, llm) ! chauffage solaire ciel clair
306      REAL, save:: cool(klon, llm) ! refroidissement infrarouge      REAL, save:: cool(klon, llm) ! refroidissement infrarouge
307      REAL cool0(klon, llm) ! refroidissement infrarouge ciel clair      REAL, save:: cool0(klon, llm) ! refroidissement infrarouge ciel clair
308      REAL, save:: topsw(klon), toplw(klon), solsw(klon), sollw(klon)      REAL, save:: topsw(klon), toplw(klon), solsw(klon)
309      real sollwdown(klon) ! downward LW flux at surface      REAL, save:: sollw(klon) ! rayonnement infrarouge montant \`a la surface
310        real, save:: sollwdown(klon) ! downward LW flux at surface
311      REAL, save:: topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)      REAL, save:: topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)
312      REAL albpla(klon)      REAL, save:: albpla(klon)
313      REAL fsollw(klon, nbsrf) ! bilan flux IR pour chaque sous surface      REAL fsollw(klon, nbsrf) ! bilan flux IR pour chaque sous surface
314      REAL fsolsw(klon, nbsrf) ! flux solaire absorb. pour chaque sous surface      REAL fsolsw(klon, nbsrf) ! flux solaire absorb. pour chaque sous surface
     SAVE albpla, sollwdown  
     SAVE heat0, cool0  
   
     INTEGER itaprad  
     SAVE itaprad  
315    
316      REAL conv_q(klon, llm) ! convergence de l'humidite (kg/kg/s)      REAL conv_q(klon, llm) ! convergence de l'humidite (kg/kg/s)
317      REAL conv_t(klon, llm) ! convergence of temperature (K/s)      REAL conv_t(klon, llm) ! convergence of temperature (K/s)
# Line 447  contains Line 321  contains
321    
322      REAL zxtsol(klon), zxqsurf(klon), zxsnow(klon), zxfluxlat(klon)      REAL zxtsol(klon), zxqsurf(klon), zxsnow(klon), zxfluxlat(klon)
323    
324      REAL dist, rmu0(klon), fract(klon)      REAL dist, mu0(klon), fract(klon)
325      REAL zdtime ! pas de temps du rayonnement (s)      real longi
     real zlongi  
   
326      REAL z_avant(klon), z_apres(klon), z_factor(klon)      REAL z_avant(klon), z_apres(klon), z_factor(klon)
     LOGICAL zx_ajustq  
   
327      REAL za, zb      REAL za, zb
328      REAL zx_t, zx_qs, zdelta, zcor      REAL zx_t, zx_qs, zcor
329      real zqsat(klon, llm)      real zqsat(klon, llm)
330      INTEGER i, k, iq, nsrf      INTEGER i, k, iq, nsrf
331      REAL t_coup      REAL, PARAMETER:: t_coup = 234.
     PARAMETER (t_coup = 234.0)  
   
332      REAL zphi(klon, llm)      REAL zphi(klon, llm)
333    
334      !IM cf. AM Variables locales pour la CLA (hbtm2)      ! cf. Anne Mathieu variables pour la couche limite atmosphérique (hbtm)
335    
336      REAL, SAVE:: pblh(klon, nbsrf) ! Hauteur de couche limite      REAL, SAVE:: pblh(klon, nbsrf) ! Hauteur de couche limite
337      REAL, SAVE:: plcl(klon, nbsrf) ! Niveau de condensation de la CLA      REAL, SAVE:: plcl(klon, nbsrf) ! Niveau de condensation de la CLA
# Line 475  contains Line 343  contains
343      REAL, SAVE:: trmb1(klon, nbsrf) ! deep_cape      REAL, SAVE:: trmb1(klon, nbsrf) ! deep_cape
344      REAL, SAVE:: trmb2(klon, nbsrf) ! inhibition      REAL, SAVE:: trmb2(klon, nbsrf) ! inhibition
345      REAL, SAVE:: trmb3(klon, nbsrf) ! Point Omega      REAL, SAVE:: trmb3(klon, nbsrf) ! Point Omega
346      ! Grdeurs de sorties      ! Grandeurs de sorties
347      REAL s_pblh(klon), s_lcl(klon), s_capCL(klon)      REAL s_pblh(klon), s_lcl(klon), s_capCL(klon)
348      REAL s_oliqCL(klon), s_cteiCL(klon), s_pblt(klon)      REAL s_oliqCL(klon), s_cteiCL(klon), s_pblt(klon)
349      REAL s_therm(klon), s_trmb1(klon), s_trmb2(klon)      REAL s_therm(klon), s_trmb1(klon), s_trmb2(klon)
350      REAL s_trmb3(klon)      REAL s_trmb3(klon)
351    
352      ! Variables locales pour la convection de K. Emanuel :      ! Variables pour la convection de K. Emanuel :
353    
354      REAL upwd(klon, llm) ! saturated updraft mass flux      REAL upwd(klon, llm) ! saturated updraft mass flux
355      REAL dnwd(klon, llm) ! saturated downdraft mass flux      REAL dnwd(klon, llm) ! saturated downdraft mass flux
356      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  
357      REAL cape(klon) ! CAPE      REAL cape(klon) ! CAPE
358      SAVE cape      SAVE cape
359    
     REAL pbase(klon) ! cloud base pressure  
     SAVE pbase  
     REAL bbase(klon) ! cloud base buoyancy  
     SAVE bbase  
     REAL rflag(klon) ! flag fonctionnement de convect  
360      INTEGER iflagctrl(klon) ! flag fonctionnement de convect      INTEGER iflagctrl(klon) ! flag fonctionnement de convect
     ! -- convect43:  
     INTEGER ntra ! nb traceurs pour convect4.3  
     REAL dtvpdt1(klon, llm), dtvpdq1(klon, llm)  
     REAL dplcldt(klon), dplcldr(klon)  
361    
362      ! Variables du changement      ! Variables du changement
363    
364      ! con: convection      ! con: convection
365      ! lsc: large scale condensation      ! lsc: large scale condensation
366      ! ajs: ajustement sec      ! ajs: ajustement sec
367      ! eva: évaporation de l'eau liquide nuageuse      ! eva: \'evaporation de l'eau liquide nuageuse
368      ! vdf: vertical diffusion in boundary layer      ! vdf: vertical diffusion in boundary layer
369      REAL d_t_con(klon, llm), d_q_con(klon, llm)      REAL d_t_con(klon, llm), d_q_con(klon, llm)
370      REAL d_u_con(klon, llm), d_v_con(klon, llm)      REAL d_u_con(klon, llm), d_v_con(klon, llm)
# Line 515  contains Line 373  contains
373      REAL d_u_ajs(klon, llm), d_v_ajs(klon, llm)      REAL d_u_ajs(klon, llm), d_v_ajs(klon, llm)
374      REAL rneb(klon, llm)      REAL rneb(klon, llm)
375    
376      REAL pmfu(klon, llm), pmfd(klon, llm)      REAL mfu(klon, llm), mfd(klon, llm)
377      REAL pen_u(klon, llm), pen_d(klon, llm)      REAL pen_u(klon, llm), pen_d(klon, llm)
378      REAL pde_u(klon, llm), pde_d(klon, llm)      REAL pde_u(klon, llm), pde_d(klon, llm)
379      INTEGER kcbot(klon), kctop(klon), kdtop(klon)      INTEGER kcbot(klon), kctop(klon), kdtop(klon)
# Line 523  contains Line 381  contains
381      REAL prfl(klon, llm + 1), psfl(klon, llm + 1)      REAL prfl(klon, llm + 1), psfl(klon, llm + 1)
382    
383      INTEGER, save:: ibas_con(klon), itop_con(klon)      INTEGER, save:: ibas_con(klon), itop_con(klon)
384        real ema_pct(klon) ! Emanuel pressure at cloud top, in Pa
385    
386      REAL rain_con(klon), rain_lsc(klon)      REAL rain_con(klon), rain_lsc(klon)
387      REAL snow_con(klon), snow_lsc(klon)      REAL, save:: snow_con(klon) ! neige (mm / s)
388        real snow_lsc(klon)
389      REAL d_ts(klon, nbsrf)      REAL d_ts(klon, nbsrf)
390    
391      REAL d_u_vdf(klon, llm), d_v_vdf(klon, llm)      REAL d_u_vdf(klon, llm), d_v_vdf(klon, llm)
# Line 549  contains Line 409  contains
409      integer:: iflag_cldcon = 1      integer:: iflag_cldcon = 1
410      logical ptconv(klon, llm)      logical ptconv(klon, llm)
411    
412      ! Variables locales pour effectuer les appels en série :      ! Variables pour effectuer les appels en s\'erie :
413    
414      REAL t_seri(klon, llm), q_seri(klon, llm)      REAL t_seri(klon, llm), q_seri(klon, llm)
415      REAL ql_seri(klon, llm), qs_seri(klon, llm)      REAL ql_seri(klon, llm)
416      REAL u_seri(klon, llm), v_seri(klon, llm)      REAL u_seri(klon, llm), v_seri(klon, llm)
417        REAL tr_seri(klon, llm, nqmx - 2)
     REAL tr_seri(klon, llm, nbtr)  
     REAL d_tr(klon, llm, nbtr)  
418    
419      REAL zx_rh(klon, llm)      REAL zx_rh(klon, llm)
420    
# Line 565  contains Line 423  contains
423      REAL zustrph(klon), zvstrph(klon)      REAL zustrph(klon), zvstrph(klon)
424      REAL aam, torsfc      REAL aam, torsfc
425    
     REAL dudyn(iim + 1, jjm + 1, llm)  
   
426      REAL zx_tmp_fi2d(klon) ! variable temporaire grille physique      REAL zx_tmp_fi2d(klon) ! variable temporaire grille physique
     REAL zx_tmp_2d(iim, jjm + 1), zx_tmp_3d(iim, jjm + 1, llm)  
427    
428      INTEGER, SAVE:: nid_day, nid_ins      INTEGER, SAVE:: nid_ins
429    
430      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.
431      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.
432      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.
433      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.
434    
     REAL zsto  
   
     character(len = 20) modname  
     character(len = 80) abort_message  
     logical ok_sync  
435      real date0      real date0
436    
437      ! Variables liées au bilan d'énergie et d'enthalpie :      ! Variables li\'ees au bilan d'\'energie et d'enthalpie :
438      REAL ztsol(klon)      REAL ztsol(klon)
439      REAL d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec      REAL d_h_vcol, d_qt, d_ec
440      REAL, SAVE:: d_h_vcol_phy      REAL, SAVE:: d_h_vcol_phy
     REAL fs_bound, fq_bound  
441      REAL zero_v(klon)      REAL zero_v(klon)
442      CHARACTER(LEN = 15) tit      CHARACTER(LEN = 20) tit
443      INTEGER:: ip_ebil = 0 ! print level for energy conservation diagnostics      INTEGER:: ip_ebil = 0 ! print level for energy conservation diagnostics
444      INTEGER:: if_ebil = 0 ! verbosity for diagnostics of energy conservation      INTEGER:: if_ebil = 0 ! verbosity for diagnostics of energy conservation
445    
446      REAL d_t_ec(klon, llm) ! tendance due à la conversion Ec -> E thermique      REAL d_t_ec(klon, llm) ! tendance due \`a la conversion Ec -> E thermique
447      REAL ZRCPD      REAL ZRCPD
448    
449      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
450      REAL u10m(klon, nbsrf), v10m(klon, nbsrf) !vents a 10m      REAL u10m(klon, nbsrf), v10m(klon, nbsrf) ! vents a 10 m
451      REAL zt2m(klon), zq2m(klon) !temp., hum. 2m moyenne s/ 1 maille      REAL zt2m(klon), zq2m(klon) ! temp., hum. 2 m moyenne s/ 1 maille
452      REAL zu10m(klon), zv10m(klon) !vents a 10m moyennes s/1 maille      REAL zu10m(klon), zv10m(klon) ! vents a 10 m moyennes s/1 maille
453      !jq Aerosol effects (Johannes Quaas, 27/11/2003)  
454      REAL sulfate(klon, llm) ! SO4 aerosol concentration [ug/m3]      ! Aerosol effects:
455    
456        REAL sulfate(klon, llm) ! SO4 aerosol concentration (micro g/m3)
457    
458      REAL, save:: sulfate_pi(klon, llm)      REAL, save:: sulfate_pi(klon, llm)
459      ! (SO4 aerosol concentration, in ug/m3, pre-industrial value)      ! SO4 aerosol concentration, in \mu g/m3, pre-industrial value
460    
461      REAL cldtaupi(klon, llm)      REAL cldtaupi(klon, llm)
462      ! (Cloud optical thickness for pre-industrial (pi) aerosols)      ! cloud optical thickness for pre-industrial (pi) aerosols
463    
464      REAL re(klon, llm) ! Cloud droplet effective radius      REAL re(klon, llm) ! Cloud droplet effective radius
465      REAL fl(klon, llm) ! denominator of re      REAL fl(klon, llm) ! denominator of re
# Line 618  contains Line 469  contains
469      REAL, save:: cg_ae(klon, llm, 2)      REAL, save:: cg_ae(klon, llm, 2)
470    
471      REAL topswad(klon), solswad(klon) ! aerosol direct effect      REAL topswad(klon), solswad(klon) ! aerosol direct effect
     ! ok_ade --> ADE = topswad - topsw  
   
472      REAL topswai(klon), solswai(klon) ! aerosol indirect effect      REAL topswai(klon), solswai(klon) ! aerosol indirect effect
     ! ok_aie .and. ok_ade --> AIE = topswai - topswad  
     ! ok_aie .and. .not. ok_ade --> AIE = topswai - topsw  
473    
474      REAL aerindex(klon) ! POLDER aerosol index      REAL aerindex(klon) ! POLDER aerosol index
475    
# Line 630  contains Line 477  contains
477      LOGICAL:: ok_aie = .false. ! apply aerosol indirect effect      LOGICAL:: ok_aie = .false. ! apply aerosol indirect effect
478    
479      REAL:: bl95_b0 = 2., bl95_b1 = 0.2      REAL:: bl95_b0 = 2., bl95_b1 = 0.2
480      ! Parameters in the formula to link CDNC to aerosol mass conc      ! Parameters in equation (D) of Boucher and Lohmann (1995, Tellus
481      ! (Boucher and Lohmann, 1995), used in nuage.F      ! B). They link cloud droplet number concentration to aerosol mass
482        ! concentration.
483    
484      SAVE u10m      SAVE u10m
485      SAVE v10m      SAVE v10m
# Line 640  contains Line 488  contains
488      SAVE ffonte      SAVE ffonte
489      SAVE fqcalving      SAVE fqcalving
490      SAVE rain_con      SAVE rain_con
     SAVE snow_con  
491      SAVE topswai      SAVE topswai
492      SAVE topswad      SAVE topswad
493      SAVE solswai      SAVE solswai
494      SAVE solswad      SAVE solswad
495      SAVE d_u_con      SAVE d_u_con
496      SAVE d_v_con      SAVE d_v_con
     SAVE rnebcon0  
     SAVE clwcon0  
497    
498      real zmasse(klon, llm)      real zmasse(klon, llm)
499      ! (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)
500    
501      real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2      real, parameter:: dobson_u = 2.1415e-05 ! Dobson unit, in kg m-2
502        integer, save:: ncid_startphy, itau_phy
503    
504      namelist /physiq_nml/ ocean, ok_veget, ok_journe, ok_mensuel, ok_instan, &      namelist /physiq_nml/ ok_journe, ok_mensuel, ok_instan, fact_cldcon, &
505           fact_cldcon, facttemps, ok_newmicro, iflag_cldcon, ratqsbas, &           facttemps, ok_newmicro, iflag_cldcon, ratqsbas, ratqshaut, if_ebil, &
506           ratqshaut, if_ebil, ok_ade, ok_aie, bl95_b0, bl95_b1, iflag_thermals, &           ok_ade, ok_aie, bl95_b0, bl95_b1, iflag_thermals, nsplit_thermals
          nsplit_thermals  
507    
508      !----------------------------------------------------------------      !----------------------------------------------------------------
509    
510      modname = 'physiq'      IF (if_ebil >= 1) zero_v = 0.
511      IF (if_ebil >= 1) THEN      IF (nqmx < 2) CALL abort_gcm('physiq', &
512         DO i = 1, klon           'eaux vapeur et liquide sont indispensables')
           zero_v(i) = 0.  
        END DO  
     END IF  
     ok_sync = .TRUE.  
     IF (nqmx < 2) THEN  
        abort_message = 'eaux vapeur et liquide sont indispensables'  
        CALL abort_gcm(modname, abort_message, 1)  
     ENDIF  
513    
514      test_firstcal: IF (firstcal) THEN      test_firstcal: IF (firstcal) THEN
515         ! initialiser         ! initialiser
# Line 685  contains Line 522  contains
522         piz_ae = 0.         piz_ae = 0.
523         tau_ae = 0.         tau_ae = 0.
524         cg_ae = 0.         cg_ae = 0.
525         rain_con(:) = 0.         rain_con = 0.
526         snow_con(:) = 0.         snow_con = 0.
527         bl95_b0 = 0.         topswai = 0.
528         bl95_b1 = 0.         topswad = 0.
529         topswai(:) = 0.         solswai = 0.
530         topswad(:) = 0.         solswad = 0.
531         solswai(:) = 0.  
532         solswad(:) = 0.         d_u_con = 0.
533           d_v_con = 0.
534         d_u_con = 0.0         rnebcon0 = 0.
535         d_v_con = 0.0         clwcon0 = 0.
536         rnebcon0 = 0.0         rnebcon = 0.
537         clwcon0 = 0.0         clwcon = 0.
        rnebcon = 0.0  
        clwcon = 0.0  
538    
539         pblh =0. ! Hauteur de couche limite         pblh =0. ! Hauteur de couche limite
540         plcl =0. ! Niveau de condensation de la CLA         plcl =0. ! Niveau de condensation de la CLA
# Line 720  contains Line 555  contains
555         read(unit=*, nml=physiq_nml)         read(unit=*, nml=physiq_nml)
556         write(unit_nml, nml=physiq_nml)         write(unit_nml, nml=physiq_nml)
557    
        ! Appel à la lecture du run.def physique  
558         call conf_phys         call conf_phys
559    
560         ! Initialiser les compteurs:         ! Initialiser les compteurs:
561    
562         frugs = 0.         frugs = 0.
563         itap = 0         CALL phyetat0(pctsrf, ftsol, ftsoil, fqsurf, qsol, &
564         itaprad = 0              fsnow, falbe, fevap, rain_fall, snow_fall, solsw, sollw, dlw, &
565         CALL phyetat0("startphy.nc", pctsrf, ftsol, ftsoil, ocean, tslab, &              radsol, frugs, agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, &
566              seaice, fqsurf, qsol, fsnow, falbe, falblw, fevap, rain_fall, &              t_ancien, q_ancien, ancien_ok, rnebcon, ratqs, clwcon, &
567              snow_fall, solsw, sollwdown, dlw, radsol, frugs, agesno, zmea, &              run_off_lic_0, sig1, w01, ncid_startphy, itau_phy)
             zstd, zsig, zgam, zthe, zpic, zval, t_ancien, q_ancien, &  
             ancien_ok, rnebcon, ratqs, clwcon, run_off_lic_0)  
568    
569         ! ATTENTION : il faudra a terme relire q2 dans l'etat initial         ! ATTENTION : il faudra a terme relire q2 dans l'etat initial
570         q2 = 1.e-8         q2 = 1e-8
571    
572         radpas = NINT(86400. / dtphys / nbapp_rad)         lmt_pas = day_step / iphysiq
573           print *, 'Number of time steps of "physics" per day: ', lmt_pas
        ! on remet le calendrier a zero  
        IF (raz_date) itau_phy = 0  
   
        PRINT *, 'cycle_diurne = ', cycle_diurne  
   
        IF(ocean.NE.'force ') THEN  
           ok_ocean = .TRUE.  
        ENDIF  
574    
575         CALL printflag(radpas, ok_ocean, ok_oasis, ok_journe, ok_instan, &         radpas = lmt_pas / nbapp_rad
             ok_region)  
576    
577         IF (dtphys*REAL(radpas) > 21600..AND.cycle_diurne) THEN         ! On remet le calendrier a zero
578            print *, 'Nbre d appels au rayonnement insuffisant'         IF (raz_date) itau_phy = 0
           print *, "Au minimum 4 appels par jour si cycle diurne"  
           abort_message = 'Nbre d appels au rayonnement insuffisant'  
           call abort_gcm(modname, abort_message, 1)  
        ENDIF  
        print *, "Clef pour la convection, iflag_con = ", iflag_con  
579    
580         ! Initialisation pour la convection de K.E. (sb):         CALL printflag(radpas, ok_journe, ok_instan, ok_region)
        IF (iflag_con >= 3) THEN  
           print *, "Convection de Kerry Emanuel 4.3"  
581    
582            DO i = 1, klon         ! Initialisation pour le sch\'ema de convection d'Emanuel :
583               ibas_con(i) = 1         IF (conv_emanuel) THEN
584               itop_con(i) = 1            ibas_con = 1
585            ENDDO            itop_con = 1
586         ENDIF         ENDIF
587    
588         IF (ok_orodr) THEN         IF (ok_orodr) THEN
# Line 776  contains Line 592  contains
592            rugoro = 0.            rugoro = 0.
593         ENDIF         ENDIF
594    
        lmt_pas = NINT(86400. / dtphys) ! tous les jours  
        print *, 'Number of time steps of "physics" per day: ', lmt_pas  
   
595         ecrit_ins = NINT(ecrit_ins/dtphys)         ecrit_ins = NINT(ecrit_ins/dtphys)
596         ecrit_hf = NINT(ecrit_hf/dtphys)         ecrit_hf = NINT(ecrit_hf/dtphys)
597         ecrit_mth = NINT(ecrit_mth/dtphys)         ecrit_mth = NINT(ecrit_mth/dtphys)
598         ecrit_tra = NINT(86400.*ecrit_tra/dtphys)         ecrit_tra = NINT(86400.*ecrit_tra/dtphys)
599         ecrit_reg = NINT(ecrit_reg/dtphys)         ecrit_reg = NINT(ecrit_reg/dtphys)
600    
        ! Initialiser le couplage si necessaire  
   
        npas = 0  
        nexca = 0  
   
601         ! Initialisation des sorties         ! Initialisation des sorties
602    
603         call ini_histhf(dtphys, nid_hf, nid_hf3d)         call ini_histins(dtphys, ok_instan, nid_ins, itau_phy)
604         call ini_histday(dtphys, ok_journe, nid_day, nqmx)         CALL ymds2ju(annee_ref, 1, day_ref, 0., date0)
605         call ini_histins(dtphys, ok_instan, nid_ins)         ! Positionner date0 pour initialisation de ORCHIDEE
606         CALL ymds2ju(annee_ref, 1, int(day_ref), 0., date0)         print *, 'physiq date0: ', date0
607         !XXXPB Positionner date0 pour initialisation de ORCHIDEE         CALL phyredem0(lmt_pas, itau_phy)
        WRITE(*, *) 'physiq date0: ', date0  
608      ENDIF test_firstcal      ENDIF test_firstcal
609    
610      ! Mettre a zero des variables de sortie (pour securite)      ! We will modify variables *_seri and we will not touch variables
611        ! u, v, t, qx:
612      DO i = 1, klon      t_seri = t
613         d_ps(i) = 0.0      u_seri = u
614      ENDDO      v_seri = v
615      DO iq = 1, nqmx      q_seri = qx(:, :, ivap)
616         DO k = 1, llm      ql_seri = qx(:, :, iliq)
617            DO i = 1, klon      tr_seri = qx(:, :, 3:nqmx)
              d_qx(i, k, iq) = 0.0  
           ENDDO  
        ENDDO  
     ENDDO  
     da = 0.  
     mp = 0.  
     phi = 0.  
   
     ! Ne pas affecter les valeurs entrées de u, v, h, et q :  
618    
619      DO k = 1, llm      ztsol = sum(ftsol * pctsrf, dim = 2)
        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  
   
     DO i = 1, klon  
        ztsol(i) = 0.  
     ENDDO  
     DO nsrf = 1, nbsrf  
        DO i = 1, klon  
           ztsol(i) = ztsol(i) + ftsol(i, nsrf)*pctsrf(i, nsrf)  
        ENDDO  
     ENDDO  
620    
621      IF (if_ebil >= 1) THEN      IF (if_ebil >= 1) THEN
622         tit = 'after dynamics'         tit = 'after dynamics'
623         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, &
624              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)
625              d_ql, d_qs, d_ec)         ! Comme les tendances de la physique sont ajout\'es dans la
        ! Comme les tendances de la physique sont ajoutés dans la  
626         !  dynamique, la variation d'enthalpie par la dynamique devrait         !  dynamique, la variation d'enthalpie par la dynamique devrait
627         !  être égale à la variation de la physique au pas de temps         !  \^etre \'egale \`a la variation de la physique au pas de temps
628         !  précédent.  Donc la somme de ces 2 variations devrait être         !  pr\'ec\'edent.  Donc la somme de ces 2 variations devrait \^etre
629         !  nulle.         !  nulle.
630         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, &
631              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, &
632              d_qt, 0., fs_bound, fq_bound)              d_qt, 0.)
633      END IF      END IF
634    
635      ! Diagnostic de la tendance dynamique :      ! Diagnostic de la tendance dynamique :
# Line 869  contains Line 643  contains
643      ELSE      ELSE
644         DO k = 1, llm         DO k = 1, llm
645            DO i = 1, klon            DO i = 1, klon
646               d_t_dyn(i, k) = 0.0               d_t_dyn(i, k) = 0.
647               d_q_dyn(i, k) = 0.0               d_q_dyn(i, k) = 0.
648            ENDDO            ENDDO
649         ENDDO         ENDDO
650         ancien_ok = .TRUE.         ancien_ok = .TRUE.
# Line 886  contains Line 660  contains
660      ! Check temperatures:      ! Check temperatures:
661      CALL hgardfou(t_seri, ftsol)      CALL hgardfou(t_seri, ftsol)
662    
663      ! Incrementer le compteur de la physique      ! Incrémenter le compteur de la physique
664      itap = itap + 1      itap = itap + 1
665      julien = MOD(NINT(rdayvrai), 360)      julien = MOD(dayvrai, 360)
666      if (julien == 0) julien = 360      if (julien == 0) julien = 360
667    
668      forall (k = 1: llm) zmasse(:, k) = (paprs(:, k)-paprs(:, k + 1)) / rg      forall (k = 1: llm) zmasse(:, k) = (paprs(:, k) - paprs(:, k + 1)) / rg
669    
670      ! Mettre en action les conditions aux limites (albedo, sst, etc.).      ! Prescrire l'ozone :
   
     ! Prescrire l'ozone et calculer l'albedo sur l'ocean.  
671      wo = ozonecm(REAL(julien), paprs)      wo = ozonecm(REAL(julien), paprs)
672    
673      ! Évaporation de l'eau liquide nuageuse :      ! \'Evaporation de l'eau liquide nuageuse :
674      DO k = 1, llm      DO k = 1, llm
675         DO i = 1, klon         DO i = 1, klon
676            zb = MAX(0., ql_seri(i, k))            zb = MAX(0., ql_seri(i, k))
# Line 912  contains Line 684  contains
684      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
685         tit = 'after reevap'         tit = 'after reevap'
686         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, &
687              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)  
688         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, &
689              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)  
   
690      END IF      END IF
691    
692      ! Appeler la diffusion verticale (programme de couche limite)      frugs = MAX(frugs, 0.000015)
693        zxrugs = sum(frugs * pctsrf, dim = 2)
694    
695      DO i = 1, klon      ! Calculs nécessaires au calcul de l'albedo dans l'interface avec
696         zxrugs(i) = 0.0      ! la surface.
     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  
697    
698      ! calculs necessaires au calcul de l'albedo dans l'interface      CALL orbite(REAL(julien), longi, dist)
   
     CALL orbite(REAL(julien), zlongi, dist)  
699      IF (cycle_diurne) THEN      IF (cycle_diurne) THEN
700         zdtime = dtphys * REAL(radpas)         CALL zenang(longi, time, dtphys * radpas, mu0, fract)
        CALL zenang(zlongi, time, zdtime, rmu0, fract)  
701      ELSE      ELSE
702         rmu0 = -999.999         mu0 = - 999.999
703      ENDIF      ENDIF
704    
705      ! Calcul de l'abedo moyen par maille      ! Calcul de l'abedo moyen par maille
706      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  
707    
708      ! Repartition sous maille des flux LW et SW      ! R\'epartition sous maille des flux longwave et shortwave
709      ! Repartition du longwave par sous-surface linearisee      ! R\'epartition du longwave par sous-surface lin\'earis\'ee
710    
711      DO nsrf = 1, nbsrf      forall (nsrf = 1: nbsrf)
712         DO i = 1, klon         fsollw(:, nsrf) = sollw + 4. * RSIGMA * ztsol**3 &
713            fsollw(i, nsrf) = sollw(i) &              * (ztsol - ftsol(:, nsrf))
714                 + 4.0*RSIGMA*ztsol(i)**3 * (ztsol(i)-ftsol(i, nsrf))         fsolsw(:, nsrf) = solsw * (1. - falbe(:, nsrf)) / (1. - albsol)
715            fsolsw(i, nsrf) = solsw(i)*(1.-falbe(i, nsrf))/(1.-albsol(i))      END forall
        ENDDO  
     ENDDO  
716    
717      fder = dlw      fder = dlw
718    
719      ! Couche limite:      ! Couche limite:
720    
721      CALL clmain(dtphys, itap, date0, pctsrf, pctsrf_new, t_seri, q_seri, &      CALL clmain(dtphys, itap, pctsrf, pctsrf_new, t_seri, q_seri, u_seri, &
722           u_seri, v_seri, julien, rmu0, co2_ppm, ok_veget, ocean, npas, nexca, &           v_seri, julien, mu0, ftsol, cdmmax, cdhmax, ksta, ksta_ter, &
723           ftsol, soil_model, cdmmax, cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, &           ok_kzmin, ftsoil, qsol, paprs, play, fsnow, fqsurf, fevap, falbe, &
724           qsol, paprs, play, fsnow, fqsurf, fevap, falbe, falblw, fluxlat, &           fluxlat, rain_fall, snow_fall, fsolsw, fsollw, fder, rlat, frugs, &
725           rain_fall, snow_fall, fsolsw, fsollw, sollwdown, fder, rlon, rlat, &           firstcal, agesno, rugoro, d_t_vdf, d_q_vdf, d_u_vdf, d_v_vdf, d_ts, &
726           cuphy, cvphy, frugs, firstcal, lafin, agesno, rugoro, d_t_vdf, &           fluxt, fluxq, fluxu, fluxv, cdragh, cdragm, q2, dsens, devap, &
727           d_q_vdf, d_u_vdf, d_v_vdf, d_ts, fluxt, fluxq, fluxu, fluxv, cdragh, &           ycoefh, yu1, yv1, t2m, q2m, u10m, v10m, pblh, capCL, oliqCL, cteiCL, &
728           cdragm, q2, dsens, devap, ycoefh, yu1, yv1, t2m, q2m, u10m, v10m, &           pblT, therm, trmb1, trmb2, trmb3, plcl, fqcalving, ffonte, &
729           pblh, capCL, oliqCL, cteiCL, pblT, therm, trmb1, trmb2, trmb3, plcl, &           run_off_lic_0)
          fqcalving, ffonte, run_off_lic_0, fluxo, fluxg, tslab, seaice)  
730    
731      ! Incrémentation des flux      ! Incr\'ementation des flux
732    
733      zxfluxt = 0.      zxfluxt = 0.
734      zxfluxq = 0.      zxfluxq = 0.
# Line 991  contains Line 737  contains
737      DO nsrf = 1, nbsrf      DO nsrf = 1, nbsrf
738         DO k = 1, llm         DO k = 1, llm
739            DO i = 1, klon            DO i = 1, klon
740               zxfluxt(i, k) = zxfluxt(i, k) + &               zxfluxt(i, k) = zxfluxt(i, k) + fluxt(i, k, nsrf) * pctsrf(i, nsrf)
741                    fluxt(i, k, nsrf) * pctsrf(i, nsrf)               zxfluxq(i, k) = zxfluxq(i, k) + fluxq(i, k, nsrf) * pctsrf(i, nsrf)
742               zxfluxq(i, k) = zxfluxq(i, k) + &               zxfluxu(i, k) = zxfluxu(i, k) + fluxu(i, k, nsrf) * pctsrf(i, nsrf)
743                    fluxq(i, k, nsrf) * pctsrf(i, nsrf)               zxfluxv(i, k) = zxfluxv(i, k) + fluxv(i, k, nsrf) * pctsrf(i, nsrf)
              zxfluxu(i, k) = zxfluxu(i, k) + &  
                   fluxu(i, k, nsrf) * pctsrf(i, nsrf)  
              zxfluxv(i, k) = zxfluxv(i, k) + &  
                   fluxv(i, k, nsrf) * pctsrf(i, nsrf)  
744            END DO            END DO
745         END DO         END DO
746      END DO      END DO
747      DO i = 1, klon      DO i = 1, klon
748         sens(i) = - zxfluxt(i, 1) ! flux de chaleur sensible au sol         sens(i) = - zxfluxt(i, 1) ! flux de chaleur sensible au sol
749         evap(i) = - zxfluxq(i, 1) ! flux d'evaporation au sol         evap(i) = - zxfluxq(i, 1) ! flux d'\'evaporation au sol
750         fder(i) = dlw(i) + dsens(i) + devap(i)         fder(i) = dlw(i) + dsens(i) + devap(i)
751      ENDDO      ENDDO
752    
# Line 1020  contains Line 762  contains
762      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
763         tit = 'after clmain'         tit = 'after clmain'
764         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, &
765              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)  
766         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, &
767              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)  
768      END IF      END IF
769    
770      ! Update surface temperature:      ! Update surface temperature:
771    
772      DO i = 1, klon      DO i = 1, klon
773         zxtsol(i) = 0.0         zxtsol(i) = 0.
774         zxfluxlat(i) = 0.0         zxfluxlat(i) = 0.
775    
776         zt2m(i) = 0.0         zt2m(i) = 0.
777         zq2m(i) = 0.0         zq2m(i) = 0.
778         zu10m(i) = 0.0         zu10m(i) = 0.
779         zv10m(i) = 0.0         zv10m(i) = 0.
780         zxffonte(i) = 0.0         zxffonte(i) = 0.
781         zxfqcalving(i) = 0.0         zxfqcalving(i) = 0.
782    
783         s_pblh(i) = 0.0         s_pblh(i) = 0.
784         s_lcl(i) = 0.0         s_lcl(i) = 0.
785         s_capCL(i) = 0.0         s_capCL(i) = 0.
786         s_oliqCL(i) = 0.0         s_oliqCL(i) = 0.
787         s_cteiCL(i) = 0.0         s_cteiCL(i) = 0.
788         s_pblT(i) = 0.0         s_pblT(i) = 0.
789         s_therm(i) = 0.0         s_therm(i) = 0.
790         s_trmb1(i) = 0.0         s_trmb1(i) = 0.
791         s_trmb2(i) = 0.0         s_trmb2(i) = 0.
792         s_trmb3(i) = 0.0         s_trmb3(i) = 0.
793    
794         IF (abs(pctsrf(i, is_ter) + pctsrf(i, is_lic) + &         IF (abs(pctsrf(i, is_ter) + pctsrf(i, is_lic) + pctsrf(i, is_oce) &
795              pctsrf(i, is_oce) + pctsrf(i, is_sic) - 1.)  >  EPSFRA) &              + pctsrf(i, is_sic) - 1.)  >  EPSFRA) print *, &
796              THEN              'physiq : probl\`eme sous surface au point ', i, &
797            WRITE(*, *) 'physiq : pb sous surface au point ', i, &              pctsrf(i, 1 : nbsrf)
                pctsrf(i, 1 : nbsrf)  
        ENDIF  
798      ENDDO      ENDDO
799      DO nsrf = 1, nbsrf      DO nsrf = 1, nbsrf
800         DO i = 1, klon         DO i = 1, klon
# Line 1084  contains Line 822  contains
822         ENDDO         ENDDO
823      ENDDO      ENDDO
824    
825      ! 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 :
   
826      DO nsrf = 1, nbsrf      DO nsrf = 1, nbsrf
827         DO i = 1, klon         DO i = 1, klon
828            IF (pctsrf(i, nsrf) < epsfra) ftsol(i, nsrf) = zxtsol(i)            IF (pctsrf(i, nsrf) < epsfra) ftsol(i, nsrf) = zxtsol(i)
# Line 1110  contains Line 847  contains
847         ENDDO         ENDDO
848      ENDDO      ENDDO
849    
850      ! Calculer la derive du flux infrarouge      ! Calculer la dérive du flux infrarouge
851    
852      DO i = 1, klon      DO i = 1, klon
853         dlw(i) = - 4.0*RSIGMA*zxtsol(i)**3         dlw(i) = - 4. * RSIGMA * zxtsol(i)**3
854      ENDDO      ENDDO
855    
856      ! Appeler la convection (au choix)      IF (check) print *, "avantcon = ", qcheck(paprs, q_seri, ql_seri)
857    
858      DO k = 1, llm      ! Appeler la convection (au choix)
        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  
     zx_ajustq = iflag_con == 2  
     IF (zx_ajustq) THEN  
        DO i = 1, klon  
           z_avant(i) = 0.0  
        ENDDO  
        DO k = 1, llm  
           DO i = 1, klon  
              z_avant(i) = z_avant(i) + (q_seri(i, k) + ql_seri(i, k)) &  
                   *zmasse(i, k)  
           ENDDO  
        ENDDO  
     ENDIF  
859    
860      select case (iflag_con)      if (conv_emanuel) then
861      case (2)         da = 0.
862         CALL conflx(dtphys, paprs, play, t_seri, q_seri, conv_t, conv_q, &         mp = 0.
863              zxfluxq(1, 1), omega, d_t_con, d_q_con, rain_con, snow_con, pmfu, &         phi = 0.
864              pmfd, pen_u, pde_u, pen_d, pde_d, kcbot, kctop, kdtop, pmflxr, &         CALL concvl(dtphys, paprs, play, t_seri, q_seri, u_seri, v_seri, sig1, &
865              pmflxs)              w01, d_t_con, d_q_con, d_u_con, d_v_con, rain_con, ibas_con, &
866         WHERE (rain_con < 0.) rain_con = 0.              itop_con, upwd, dnwd, dnwd0, Ma, cape, iflagctrl, qcondc, wd, &
867         WHERE (snow_con < 0.) snow_con = 0.              pmflxr, da, phi, mp)
868         DO i = 1, klon         snow_con = 0.
           ibas_con(i) = llm + 1 - kcbot(i)  
           itop_con(i) = llm + 1 - kctop(i)  
        ENDDO  
     case (3:)  
        ! 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.  
        ntra = 1  
        ! Schéma de convection modularisé et vectorisé :  
        ! (driver commun aux versions 3 et 4)  
   
        CALL concvl(iflag_con, dtphys, paprs, play, t_seri, q_seri, u_seri, &  
             v_seri, tr_seri, ntra, ema_work1, ema_work2, 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)  
869         clwcon0 = qcondc         clwcon0 = qcondc
870         pmfu = upwd + dnwd         mfu = upwd + dnwd
871           IF (.NOT. ok_gust) wd = 0.
872    
873         IF (.NOT. ok_gust) THEN         IF (thermcep) THEN
874            do i = 1, klon            zqsat = MIN(0.5, r2es * FOEEW(t_seri, rtt >= t_seri) / play)
875               wd(i) = 0.0            zqsat = zqsat / (1. - retv * zqsat)
876            enddo         ELSE
877              zqsat = merge(qsats(t_seri), qsatl(t_seri), t_seri < t_coup) / play
878         ENDIF         ENDIF
879    
880         ! Calcul des propriétés des nuages convectifs         ! Properties of convective clouds
881           clwcon0 = fact_cldcon * clwcon0
        DO k = 1, llm  
           DO i = 1, klon  
              zx_t = t_seri(i, k)  
              IF (thermcep) THEN  
                 zdelta = MAX(0., SIGN(1., rtt-zx_t))  
                 zx_qs = r2es * FOEEW(zx_t, zdelta)/play(i, k)  
                 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  
   
        ! calcul des proprietes des nuages convectifs  
        clwcon0 = fact_cldcon*clwcon0  
882         call clouds_gno(klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, &         call clouds_gno(klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, &
883              rnebcon0)              rnebcon0)
884      case default  
885         print *, "iflag_con non-prevu", iflag_con         forall (i = 1:klon) ema_pct(i) = paprs(i,itop_con(i) + 1)
886         stop 1         mfd = 0.
887      END select         pen_u = 0.
888           pen_d = 0.
889           pde_d = 0.
890           pde_u = 0.
891        else
892           conv_q = d_q_dyn + d_q_vdf / dtphys
893           conv_t = d_t_dyn + d_t_vdf / dtphys
894           z_avant = sum((q_seri + ql_seri) * zmasse, dim=2)
895           CALL conflx(dtphys, paprs, play, t_seri(:, llm:1:- 1), &
896                q_seri(:, llm:1:- 1), conv_t, conv_q, zxfluxq(:, 1), omega, &
897                d_t_con, d_q_con, rain_con, snow_con, mfu(:, llm:1:- 1), &
898                mfd(:, llm:1:- 1), pen_u, pde_u, pen_d, pde_d, kcbot, kctop, &
899                kdtop, pmflxr, pmflxs)
900           WHERE (rain_con < 0.) rain_con = 0.
901           WHERE (snow_con < 0.) snow_con = 0.
902           ibas_con = llm + 1 - kcbot
903           itop_con = llm + 1 - kctop
904        END if
905    
906      DO k = 1, llm      DO k = 1, llm
907         DO i = 1, klon         DO i = 1, klon
# Line 1222  contains Line 915  contains
915      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
916         tit = 'after convect'         tit = 'after convect'
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         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, &
920              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)  
921      END IF      END IF
922    
923      IF (check) THEN      IF (check) THEN
924         za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)         za = qcheck(paprs, q_seri, ql_seri)
925         print *, "aprescon = ", za         print *, "aprescon = ", za
926         zx_t = 0.0         zx_t = 0.
927         za = 0.0         za = 0.
928         DO i = 1, klon         DO i = 1, klon
929            za = za + airephy(i)/REAL(klon)            za = za + airephy(i)/REAL(klon)
930            zx_t = zx_t + (rain_con(i)+ &            zx_t = zx_t + (rain_con(i)+ &
# Line 1242  contains Line 933  contains
933         zx_t = zx_t/za*dtphys         zx_t = zx_t/za*dtphys
934         print *, "Precip = ", zx_t         print *, "Precip = ", zx_t
935      ENDIF      ENDIF
936      IF (zx_ajustq) THEN  
937         DO i = 1, klon      IF (.not. conv_emanuel) THEN
938            z_apres(i) = 0.0         z_apres = sum((q_seri + ql_seri) * zmasse, dim=2)
939         ENDDO         z_factor = (z_avant - (rain_con + snow_con) * dtphys) / z_apres
        DO k = 1, llm  
           DO i = 1, klon  
              z_apres(i) = z_apres(i) + (q_seri(i, k) + ql_seri(i, k)) &  
                   *zmasse(i, k)  
           ENDDO  
        ENDDO  
        DO i = 1, klon  
           z_factor(i) = (z_avant(i)-(rain_con(i) + snow_con(i))*dtphys) &  
                /z_apres(i)  
        ENDDO  
940         DO k = 1, llm         DO k = 1, llm
941            DO i = 1, klon            DO i = 1, klon
942               IF (z_factor(i) > 1. + 1E-8 .OR. z_factor(i) < 1. - 1E-8) THEN               IF (z_factor(i) > 1. + 1E-8 .OR. z_factor(i) < 1. - 1E-8) THEN
# Line 1264  contains Line 945  contains
945            ENDDO            ENDDO
946         ENDDO         ENDDO
947      ENDIF      ENDIF
     zx_ajustq = .FALSE.  
948    
949      ! Convection sèche (thermiques ou ajustement)      ! Convection s\`eche (thermiques ou ajustement)
950    
951      d_t_ajs = 0.      d_t_ajs = 0.
952      d_u_ajs = 0.      d_u_ajs = 0.
# Line 1289  contains Line 969  contains
969      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
970         tit = 'after dry_adjust'         tit = 'after dry_adjust'
971         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, &
972              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)  
973      END IF      END IF
974    
975      ! Caclul des ratqs      ! Caclul des ratqs
976    
977      ! ratqs convectifs a l'ancienne en fonction de q(z = 0)-q / q      ! ratqs convectifs \`a l'ancienne en fonction de (q(z = 0) - q) / q
978      ! on ecrase le tableau ratqsc calcule par clouds_gno      ! on \'ecrase le tableau ratqsc calcul\'e par clouds_gno
979      if (iflag_cldcon == 1) then      if (iflag_cldcon == 1) then
980         do k = 1, llm         do k = 1, llm
981            do i = 1, klon            do i = 1, klon
982               if(ptconv(i, k)) then               if(ptconv(i, k)) then
983                  ratqsc(i, k) = ratqsbas &                  ratqsc(i, k) = ratqsbas + fact_cldcon &
984                       +fact_cldcon*(q_seri(i, 1)-q_seri(i, k))/q_seri(i, k)                       * (q_seri(i, 1) - q_seri(i, k)) / q_seri(i, k)
985               else               else
986                  ratqsc(i, k) = 0.                  ratqsc(i, k) = 0.
987               endif               endif
# Line 1313  contains Line 992  contains
992      ! ratqs stables      ! ratqs stables
993      do k = 1, llm      do k = 1, llm
994         do i = 1, klon         do i = 1, klon
995            ratqss(i, k) = ratqsbas + (ratqshaut-ratqsbas)* &            ratqss(i, k) = ratqsbas + (ratqshaut - ratqsbas) &
996                 min((paprs(i, 1)-play(i, k))/(paprs(i, 1)-30000.), 1.)                 * min((paprs(i, 1) - play(i, k)) / (paprs(i, 1) - 3e4), 1.)
997         enddo         enddo
998      enddo      enddo
999    
1000      ! ratqs final      ! ratqs final
1001      if (iflag_cldcon == 1 .or.iflag_cldcon == 2) then      if (iflag_cldcon == 1 .or. iflag_cldcon == 2) then
1002         ! les ratqs sont une conbinaison de ratqss et ratqsc         ! les ratqs sont une conbinaison de ratqss et ratqsc
1003         ! ratqs final         ! ratqs final
1004         ! 1e4 (en gros 3 heures), en dur pour le moment, est le temps de         ! 1e4 (en gros 3 heures), en dur pour le moment, est le temps de
1005         ! relaxation des ratqs         ! relaxation des ratqs
1006         facteur = exp(-dtphys*facttemps)         ratqs = max(ratqs * exp(- dtphys * facttemps), ratqss)
        ratqs = max(ratqs*facteur, ratqss)  
1007         ratqs = max(ratqs, ratqsc)         ratqs = max(ratqs, ratqsc)
1008      else      else
1009         ! on ne prend que le ratqs stable pour fisrtilp         ! on ne prend que le ratqs stable pour fisrtilp
1010         ratqs = ratqss         ratqs = ratqss
1011      endif      endif
1012    
     ! Processus de condensation à grande echelle et processus de  
     ! précipitation :  
1013      CALL fisrtilp(dtphys, paprs, play, t_seri, q_seri, ptconv, ratqs, &      CALL fisrtilp(dtphys, paprs, play, t_seri, q_seri, ptconv, ratqs, &
1014           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, &
1015           pfrac_impa, pfrac_nucl, pfrac_1nucl, frac_impa, frac_nucl, prfl, &           pfrac_impa, pfrac_nucl, pfrac_1nucl, frac_impa, frac_nucl, prfl, &
# Line 1351  contains Line 1027  contains
1027         ENDDO         ENDDO
1028      ENDDO      ENDDO
1029      IF (check) THEN      IF (check) THEN
1030         za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)         za = qcheck(paprs, q_seri, ql_seri)
1031         print *, "apresilp = ", za         print *, "apresilp = ", za
1032         zx_t = 0.0         zx_t = 0.
1033         za = 0.0         za = 0.
1034         DO i = 1, klon         DO i = 1, klon
1035            za = za + airephy(i)/REAL(klon)            za = za + airephy(i)/REAL(klon)
1036            zx_t = zx_t + (rain_lsc(i) &            zx_t = zx_t + (rain_lsc(i) &
# Line 1367  contains Line 1043  contains
1043      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
1044         tit = 'after fisrt'         tit = 'after fisrt'
1045         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, &
1046              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)  
1047         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, &
1048              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)  
1049      END IF      END IF
1050    
1051      ! PRESCRIPTION DES NUAGES POUR LE RAYONNEMENT      ! PRESCRIPTION DES NUAGES POUR LE RAYONNEMENT
1052    
1053      ! 1. NUAGES CONVECTIFS      ! 1. NUAGES CONVECTIFS
1054    
1055      IF (iflag_cldcon <= -1) THEN      IF (iflag_cldcon <= - 1) THEN
1056         ! seulement pour Tiedtke         ! seulement pour Tiedtke
1057         snow_tiedtke = 0.         snow_tiedtke = 0.
1058         if (iflag_cldcon == -1) then         if (iflag_cldcon == - 1) then
1059            rain_tiedtke = rain_con            rain_tiedtke = rain_con
1060         else         else
1061            rain_tiedtke = 0.            rain_tiedtke = 0.
1062            do k = 1, llm            do k = 1, llm
1063               do i = 1, klon               do i = 1, klon
1064                  if (d_q_con(i, k) < 0.) then                  if (d_q_con(i, k) < 0.) then
1065                     rain_tiedtke(i) = rain_tiedtke(i)-d_q_con(i, k)/dtphys &                     rain_tiedtke(i) = rain_tiedtke(i) - d_q_con(i, k)/dtphys &
1066                          *zmasse(i, k)                          *zmasse(i, k)
1067                  endif                  endif
1068               enddo               enddo
# Line 1396  contains Line 1070  contains
1070         endif         endif
1071    
1072         ! Nuages diagnostiques pour Tiedtke         ! Nuages diagnostiques pour Tiedtke
1073         CALL diagcld1(paprs, play, &         CALL diagcld1(paprs, play, rain_tiedtke, snow_tiedtke, ibas_con, &
1074              rain_tiedtke, snow_tiedtke, ibas_con, itop_con, &              itop_con, diafra, dialiq)
             diafra, dialiq)  
1075         DO k = 1, llm         DO k = 1, llm
1076            DO i = 1, klon            DO i = 1, klon
1077               IF (diafra(i, k) > cldfra(i, k)) THEN               IF (diafra(i, k) > cldfra(i, k)) THEN
# Line 1408  contains Line 1081  contains
1081            ENDDO            ENDDO
1082         ENDDO         ENDDO
1083      ELSE IF (iflag_cldcon == 3) THEN      ELSE IF (iflag_cldcon == 3) THEN
1084         ! On prend pour les nuages convectifs le max du calcul de la         ! On prend pour les nuages convectifs le maximum du calcul de
1085         ! convection et du calcul du pas de temps précédent diminué d'un facteur         ! la convection et du calcul du pas de temps pr\'ec\'edent diminu\'e
1086         ! facttemps         ! d'un facteur facttemps.
1087         facteur = dtphys *facttemps         facteur = dtphys * facttemps
1088         do k = 1, llm         do k = 1, llm
1089            do i = 1, klon            do i = 1, klon
1090               rnebcon(i, k) = rnebcon(i, k)*facteur               rnebcon(i, k) = rnebcon(i, k) * facteur
1091               if (rnebcon0(i, k)*clwcon0(i, k) > rnebcon(i, k)*clwcon(i, k)) &               if (rnebcon0(i, k) * clwcon0(i, k) &
1092                    then                    > rnebcon(i, k) * clwcon(i, k)) then
1093                  rnebcon(i, k) = rnebcon0(i, k)                  rnebcon(i, k) = rnebcon0(i, k)
1094                  clwcon(i, k) = clwcon0(i, k)                  clwcon(i, k) = clwcon0(i, k)
1095               endif               endif
# Line 1449  contains Line 1122  contains
1122      ENDDO      ENDDO
1123    
1124      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, &
1125           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, &
1126           d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)           d_qt, d_ec)
1127    
1128      ! Humidité relative pour diagnostic :      ! Humidit\'e relative pour diagnostic :
1129      DO k = 1, llm      DO k = 1, llm
1130         DO i = 1, klon         DO i = 1, klon
1131            zx_t = t_seri(i, k)            zx_t = t_seri(i, k)
1132            IF (thermcep) THEN            IF (thermcep) THEN
1133               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)  
1134               zx_qs = MIN(0.5, zx_qs)               zx_qs = MIN(0.5, zx_qs)
1135               zcor = 1./(1.-retv*zx_qs)               zcor = 1./(1. - retv*zx_qs)
1136               zx_qs = zx_qs*zcor               zx_qs = zx_qs*zcor
1137            ELSE            ELSE
1138               IF (zx_t < t_coup) THEN               IF (zx_t < t_coup) THEN
# Line 1477  contains Line 1149  contains
1149      ! Introduce the aerosol direct and first indirect radiative forcings:      ! Introduce the aerosol direct and first indirect radiative forcings:
1150      IF (ok_ade .OR. ok_aie) THEN      IF (ok_ade .OR. ok_aie) THEN
1151         ! Get sulfate aerosol distribution :         ! Get sulfate aerosol distribution :
1152         CALL readsulfate(rdayvrai, firstcal, sulfate)         CALL readsulfate(dayvrai, time, firstcal, sulfate)
1153         CALL readsulfate_preind(rdayvrai, firstcal, sulfate_pi)         CALL readsulfate_preind(dayvrai, time, firstcal, sulfate_pi)
1154    
1155         CALL aeropt(play, paprs, t_seri, sulfate, rhcl, tau_ae, piz_ae, cg_ae, &         CALL aeropt(play, paprs, t_seri, sulfate, rhcl, tau_ae, piz_ae, cg_ae, &
1156              aerindex)              aerindex)
# Line 1488  contains Line 1160  contains
1160         cg_ae = 0.         cg_ae = 0.
1161      ENDIF      ENDIF
1162    
1163      ! Paramètres optiques des nuages et quelques paramètres pour diagnostics :      ! Param\`etres optiques des nuages et quelques param\`etres pour
1164        ! diagnostics :
1165      if (ok_newmicro) then      if (ok_newmicro) then
1166         CALL newmicro(paprs, play, ok_newmicro, t_seri, cldliq, cldfra, &         CALL newmicro(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, &
1167              cldtau, cldemi, cldh, cldl, cldm, cldt, cldq, flwp, fiwp, flwc, &              cldh, cldl, cldm, cldt, cldq, flwp, fiwp, flwc, fiwc, ok_aie, &
1168              fiwc, ok_aie, sulfate, sulfate_pi, bl95_b0, bl95_b1, cldtaupi, &              sulfate, sulfate_pi, bl95_b0, bl95_b1, cldtaupi, re, fl)
             re, fl)  
1169      else      else
1170         CALL nuage(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, cldh, &         CALL nuage(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, cldh, &
1171              cldl, cldm, cldt, cldq, ok_aie, sulfate, sulfate_pi, bl95_b0, &              cldl, cldm, cldt, cldq, ok_aie, sulfate, sulfate_pi, bl95_b0, &
1172              bl95_b1, cldtaupi, re, fl)              bl95_b1, cldtaupi, re, fl)
1173      endif      endif
1174    
1175      ! Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.      IF (MOD(itap - 1, radpas) == 0) THEN
1176      IF (MOD(itaprad, radpas) == 0) THEN         ! Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.
1177         DO i = 1, klon         ! Calcul de l'abedo moyen par maille
1178            albsol(i) = falbe(i, is_oce) * pctsrf(i, is_oce) &         albsol = sum(falbe * pctsrf, dim = 2)
1179                 + 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  
1180         ! Rayonnement (compatible Arpege-IFS) :         ! Rayonnement (compatible Arpege-IFS) :
1181         CALL radlwsw(dist, rmu0, fract, paprs, play, zxtsol, albsol, &         CALL radlwsw(dist, mu0, fract, paprs, play, zxtsol, albsol, t_seri, &
1182              albsollw, t_seri, q_seri, wo, cldfra, cldemi, cldtau, heat, &              q_seri, wo, cldfra, cldemi, cldtau, heat, heat0, cool, cool0, &
1183              heat0, cool, cool0, radsol, albpla, topsw, toplw, solsw, sollw, &              radsol, albpla, topsw, toplw, solsw, sollw, sollwdown, topsw0, &
1184              sollwdown, topsw0, toplw0, solsw0, sollw0, lwdn0, lwdn, lwup0, &              toplw0, solsw0, sollw0, lwdn0, lwdn, lwup0, lwup, swdn0, swdn, &
1185              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, &
1186              cg_ae, topswad, solswad, cldtaupi, topswai, solswai)              solswad, cldtaupi, topswai, solswai)
        itaprad = 0  
1187      ENDIF      ENDIF
     itaprad = itaprad + 1  
1188    
1189      ! Ajouter la tendance des rayonnements (tous les pas)      ! Ajouter la tendance des rayonnements (tous les pas)
1190    
1191      DO k = 1, llm      DO k = 1, llm
1192         DO i = 1, klon         DO i = 1, klon
1193            t_seri(i, k) = t_seri(i, k) + (heat(i, k)-cool(i, k)) * dtphys/86400.            t_seri(i, k) = t_seri(i, k) + (heat(i, k) - cool(i, k)) * dtphys/86400.
1194         ENDDO         ENDDO
1195      ENDDO      ENDDO
1196    
1197      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
1198         tit = 'after rad'         tit = 'after rad'
1199         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, &
1200              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)  
1201         call diagphy(airephy, tit, ip_ebil, topsw, toplw, solsw, sollw, &         call diagphy(airephy, tit, ip_ebil, topsw, toplw, solsw, sollw, &
1202              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)  
1203      END IF      END IF
1204    
1205      ! Calculer l'hydrologie de la surface      ! Calculer l'hydrologie de la surface
1206      DO i = 1, klon      DO i = 1, klon
1207         zxqsurf(i) = 0.0         zxqsurf(i) = 0.
1208         zxsnow(i) = 0.0         zxsnow(i) = 0.
1209      ENDDO      ENDDO
1210      DO nsrf = 1, nbsrf      DO nsrf = 1, nbsrf
1211         DO i = 1, klon         DO i = 1, klon
# Line 1553  contains Line 1214  contains
1214         ENDDO         ENDDO
1215      ENDDO      ENDDO
1216    
1217      ! 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)
1218    
1219      DO i = 1, klon      DO i = 1, klon
1220         bils(i) = radsol(i) - sens(i) + zxfluxlat(i)         bils(i) = radsol(i) - sens(i) + zxfluxlat(i)
1221      ENDDO      ENDDO
1222    
1223      ! Paramétrisation de l'orographie à l'échelle sous-maille :      ! Param\'etrisation de l'orographie \`a l'\'echelle sous-maille :
1224    
1225      IF (ok_orodr) THEN      IF (ok_orodr) THEN
1226         ! selection des points pour lesquels le shema est actif:         ! S\'election des points pour lesquels le sch\'ema est actif :
1227         igwd = 0         igwd = 0
1228         DO i = 1, klon         DO i = 1, klon
1229            itest(i) = 0            itest(i) = 0
1230            IF (((zpic(i)-zmea(i)) > 100.).AND.(zstd(i) > 10.0)) THEN            IF (zpic(i) - zmea(i) > 100. .AND. zstd(i) > 10.) THEN
1231               itest(i) = 1               itest(i) = 1
1232               igwd = igwd + 1               igwd = igwd + 1
              idx(igwd) = i  
1233            ENDIF            ENDIF
1234         ENDDO         ENDDO
1235    
1236         CALL drag_noro(klon, llm, dtphys, paprs, play, zmea, zstd, zsig, zgam, &         CALL drag_noro(klon, llm, dtphys, paprs, play, zmea, zstd, zsig, zgam, &
1237              zthe, zpic, zval, igwd, idx, itest, t_seri, u_seri, v_seri, &              zthe, zpic, zval, itest, t_seri, u_seri, v_seri, zulow, zvlow, &
1238              zulow, zvlow, zustrdr, zvstrdr, d_t_oro, d_u_oro, d_v_oro)              zustrdr, zvstrdr, d_t_oro, d_u_oro, d_v_oro)
1239    
1240         ! ajout des tendances         ! ajout des tendances
1241         DO k = 1, llm         DO k = 1, llm
# Line 1588  contains Line 1248  contains
1248      ENDIF      ENDIF
1249    
1250      IF (ok_orolf) THEN      IF (ok_orolf) THEN
1251         ! Sélection des points pour lesquels le schéma est actif :         ! S\'election des points pour lesquels le sch\'ema est actif :
1252         igwd = 0         igwd = 0
1253         DO i = 1, klon         DO i = 1, klon
1254            itest(i) = 0            itest(i) = 0
1255            IF ((zpic(i) - zmea(i)) > 100.) THEN            IF (zpic(i) - zmea(i) > 100.) THEN
1256               itest(i) = 1               itest(i) = 1
1257               igwd = igwd + 1               igwd = igwd + 1
              idx(igwd) = i  
1258            ENDIF            ENDIF
1259         ENDDO         ENDDO
1260    
# Line 1613  contains Line 1272  contains
1272         ENDDO         ENDDO
1273      ENDIF      ENDIF
1274    
1275      ! Stress nécessaires : toute la physique      ! Stress n\'ecessaires : toute la physique
1276    
1277      DO i = 1, klon      DO i = 1, klon
1278         zustrph(i) = 0.         zustrph(i) = 0.
# Line 1628  contains Line 1287  contains
1287         ENDDO         ENDDO
1288      ENDDO      ENDDO
1289    
1290      CALL aaam_bud(ra, rg, romega, rlat, rlon, pphis, zustrdr, zustrli, &      CALL aaam_bud(rg, romega, rlat, rlon, pphis, zustrdr, zustrli, zustrph, &
1291           zustrph, zvstrdr, zvstrli, zvstrph, paprs, u, v, aam, torsfc)           zvstrdr, zvstrli, zvstrph, paprs, u, v, aam, torsfc)
1292    
1293      IF (if_ebil >= 2) CALL diagetpq(airephy, 'after orography', ip_ebil, 2, &      IF (if_ebil >= 2) CALL diagetpq(airephy, 'after orography', ip_ebil, 2, &
1294           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, &
1295           d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec)           d_qt, d_ec)
1296    
1297      ! Calcul des tendances traceurs      ! Calcul des tendances traceurs
1298      call phytrac(rnpb, itap, lmt_pas, julien, time, firstcal, lafin, nqmx-2, &      call phytrac(itap, lmt_pas, julien, time, firstcal, lafin, dtphys, t, &
1299           dtphys, u, t, paprs, play, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &           paprs, play, mfu, mfd, pde_u, pen_d, ycoefh, fm_therm, entr_therm, &
1300           ycoefh, fm_therm, entr_therm, yu1, yv1, ftsol, pctsrf, frac_impa, &           yu1, yv1, ftsol, pctsrf, frac_impa, frac_nucl, da, phi, mp, upwd, &
1301           frac_nucl, pphis, albsol, rhcl, cldfra, rneb, diafra, cldliq, &           dnwd, tr_seri, zmasse, ncid_startphy, nid_ins, itau_phy)
1302           pmflxr, pmflxs, prfl, psfl, da, phi, mp, upwd, dnwd, tr_seri, zmasse)  
1303        IF (offline) call phystokenc(dtphys, rlon, rlat, t, mfu, mfd, pen_u, &
1304      IF (offline) THEN           pde_u, pen_d, pde_d, fm_therm, entr_therm, ycoefh, yu1, yv1, ftsol, &
1305         call phystokenc(dtphys, rlon, rlat, t, pmfu, pmfd, pen_u, pde_u, &           pctsrf, frac_impa, frac_nucl, pphis, airephy, dtphys, itap)
             pen_d, pde_d, fm_therm, entr_therm, ycoefh, yu1, yv1, ftsol, &  
             pctsrf, frac_impa, frac_nucl, pphis, airephy, dtphys, itap)  
     ENDIF  
1306    
1307      ! Calculer le transport de l'eau et de l'energie (diagnostique)      ! Calculer le transport de l'eau et de l'energie (diagnostique)
1308      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)  
1309    
1310      ! diag. bilKP      ! diag. bilKP
1311    
1312      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, &
1313           ve_lay, vq_lay, ue_lay, uq_lay)           ve_lay, vq_lay, ue_lay, uq_lay)
1314    
1315      ! Accumuler les variables a stocker dans les fichiers histoire:      ! Accumuler les variables a stocker dans les fichiers histoire:
# Line 1673  contains Line 1328  contains
1328      IF (if_ebil >= 1) THEN      IF (if_ebil >= 1) THEN
1329         tit = 'after physic'         tit = 'after physic'
1330         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, &
1331              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)  
1332         ! Comme les tendances de la physique sont ajoute dans la dynamique,         ! Comme les tendances de la physique sont ajoute dans la dynamique,
1333         ! on devrait avoir que la variation d'entalpie par la dynamique         ! on devrait avoir que la variation d'entalpie par la dynamique
1334         ! 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.
1335         ! Donc la somme de ces 2 variations devrait etre nulle.         ! Donc la somme de ces 2 variations devrait etre nulle.
1336         call diagphy(airephy, tit, ip_ebil, topsw, toplw, solsw, sollw, sens, &         call diagphy(airephy, tit, ip_ebil, topsw, toplw, solsw, sollw, sens, &
1337              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)  
   
1338         d_h_vcol_phy = d_h_vcol         d_h_vcol_phy = d_h_vcol
   
1339      END IF      END IF
1340    
1341      ! SORTIES      ! SORTIES
1342    
1343      !cc prw = eau precipitable      ! prw = eau precipitable
1344      DO i = 1, klon      DO i = 1, klon
1345         prw(i) = 0.         prw(i) = 0.
1346         DO k = 1, llm         DO k = 1, llm
# Line 1709  contains Line 1360  contains
1360         ENDDO         ENDDO
1361      ENDDO      ENDDO
1362    
1363      IF (nqmx >= 3) THEN      DO iq = 3, nqmx
1364         DO iq = 3, nqmx         DO k = 1, llm
1365            DO k = 1, llm            DO i = 1, klon
1366               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  
1367            ENDDO            ENDDO
1368         ENDDO         ENDDO
1369      ENDIF      ENDDO
1370    
1371      ! 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:
1372      DO k = 1, llm      DO k = 1, llm
# Line 1727  contains Line 1376  contains
1376         ENDDO         ENDDO
1377      ENDDO      ENDDO
1378    
     ! Ecriture des sorties  
     call write_histhf  
     call write_histday  
1379      call write_histins      call write_histins
1380    
1381      ! Si c'est la fin, il faut conserver l'etat de redemarrage      IF (lafin) then
1382      IF (lafin) THEN         call NF95_CLOSE(ncid_startphy)
1383         itau_phy = itau_phy + itap         CALL phyredem(pctsrf, ftsol, ftsoil, fqsurf, qsol, &
1384         CALL phyredem("restartphy.nc", rlat, rlon, pctsrf, ftsol, ftsoil, &              fsnow, falbe, fevap, rain_fall, snow_fall, solsw, sollw, dlw, &
1385              tslab, seaice, fqsurf, qsol, fsnow, falbe, falblw, fevap, &              radsol, frugs, agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, &
1386              rain_fall, snow_fall, solsw, sollwdown, dlw, radsol, frugs, &              t_ancien, q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, &
1387              agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, t_ancien, &              w01)
1388              q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0)      end IF
     ENDIF  
1389    
1390      firstcal = .FALSE.      firstcal = .FALSE.
1391    
1392    contains    contains
1393    
     subroutine write_histday  
   
       use gr_phy_write_3d_m, only: gr_phy_write_3d  
       integer itau_w ! pas de temps ecriture  
   
       !------------------------------------------------  
   
       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  
   
     !***************************************************************  
   
1394      subroutine write_histins      subroutine write_histins
1395    
1396        ! From phylmd/write_histins.h, version 1.2 2005/05/25 13:10:09        ! From phylmd/write_histins.h, version 1.2 2005/05/25 13:10:09
1397    
1398        real zout        ! Ecriture des sorties
1399        integer itau_w ! pas de temps ecriture  
1400          use dimens_m, only: iim, jjm
1401          use gr_fi_ecrit_m, only: gr_fi_ecrit
1402          USE histsync_m, ONLY: histsync
1403          USE histwrite_m, ONLY: histwrite
1404    
1405          integer i, itau_w ! pas de temps ecriture
1406          REAL zx_tmp_2d(iim, jjm + 1), zx_tmp_3d(iim, jjm + 1, llm)
1407    
1408        !--------------------------------------------------        !--------------------------------------------------
1409    
1410        IF (ok_instan) THEN        IF (ok_instan) THEN
1411           ! Champs 2D:           ! Champs 2D:
1412    
          zsto = dtphys * ecrit_ins  
          zout = dtphys * ecrit_ins  
1413           itau_w = itau_phy + itap           itau_w = itau_phy + itap
1414    
          i = NINT(zout/zsto)  
1415           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, pphis, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, pphis, zx_tmp_2d)
1416           CALL histwrite(nid_ins, "phis", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "phis", itau_w, zx_tmp_2d)
1417    
          i = NINT(zout/zsto)  
1418           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, airephy, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, airephy, zx_tmp_2d)
1419           CALL histwrite(nid_ins, "aire", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "aire", itau_w, zx_tmp_2d)
1420    
# Line 1875  contains Line 1484  contains
1484           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, bils, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, bils, zx_tmp_2d)
1485           CALL histwrite(nid_ins, "bils", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "bils", itau_w, zx_tmp_2d)
1486    
1487           zx_tmp_fi2d(1:klon) = -1*sens(1:klon)           zx_tmp_fi2d(1:klon) = - sens(1:klon)
1488           ! CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sens, zx_tmp_2d)           ! CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sens, zx_tmp_2d)
1489           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1490           CALL histwrite(nid_ins, "sens", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "sens", itau_w, zx_tmp_2d)
# Line 1937  contains Line 1546  contains
1546              CALL histwrite(nid_ins, "rugs_"//clnsurf(nsrf), itau_w, &              CALL histwrite(nid_ins, "rugs_"//clnsurf(nsrf), itau_w, &
1547                   zx_tmp_2d)                   zx_tmp_2d)
1548    
1549              zx_tmp_fi2d(1 : klon) = falbe(1 : klon, nsrf)              zx_tmp_fi2d(1 : klon) = falbe(:, nsrf)
1550              CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)              CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zx_tmp_fi2d, zx_tmp_2d)
1551              CALL histwrite(nid_ins, "albe_"//clnsurf(nsrf), itau_w, &              CALL histwrite(nid_ins, "albe_"//clnsurf(nsrf), itau_w, &
1552                   zx_tmp_2d)                   zx_tmp_2d)
# Line 1945  contains Line 1554  contains
1554           END DO           END DO
1555           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, albsol, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, albsol, zx_tmp_2d)
1556           CALL histwrite(nid_ins, "albs", itau_w, 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)  
1557    
1558           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zxrugs, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zxrugs, zx_tmp_2d)
1559           CALL histwrite(nid_ins, "rugs", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "rugs", itau_w, zx_tmp_2d)
# Line 1983  contains Line 1590  contains
1590           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb3, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb3, zx_tmp_2d)
1591           CALL histwrite(nid_ins, "s_trmb3", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "s_trmb3", itau_w, zx_tmp_2d)
1592    
1593             if (conv_emanuel) then
1594                CALL gr_fi_ecrit(1, klon, iim, jjm + 1, ema_pct, zx_tmp_2d)
1595                CALL histwrite(nid_ins, "ptop", itau_w, zx_tmp_2d)
1596             end if
1597    
1598           ! Champs 3D:           ! Champs 3D:
1599    
1600           CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, t_seri, zx_tmp_3d)           CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, t_seri, zx_tmp_3d)
# Line 2006  contains Line 1618  contains
1618           CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, d_q_vdf, zx_tmp_3d)           CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, d_q_vdf, zx_tmp_3d)
1619           CALL histwrite(nid_ins, "dqvdf", itau_w, zx_tmp_3d)           CALL histwrite(nid_ins, "dqvdf", itau_w, zx_tmp_3d)
1620    
1621           if (ok_sync) then           CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, zx_rh, zx_tmp_3d)
1622              call histsync(nid_ins)           CALL histwrite(nid_ins, "rhum", itau_w, zx_tmp_3d)
1623           endif  
1624             call histsync(nid_ins)
1625        ENDIF        ENDIF
1626    
1627      end subroutine write_histins      end subroutine write_histins
1628    
     !****************************************************  
   
     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  
   
       if (ok_sync) then  
          call histsync(nid_hf3d)  
       endif  
   
     end subroutine write_histhf3d  
   
1629    END SUBROUTINE physiq    END SUBROUTINE physiq
1630    
1631  end module physiq_m  end module physiq_m
Legend:
Removed from v.68  
changed lines
  Added in v.186
  ViewVC Help
Powered by ViewVC 1.1.21