/[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 49 by guez, Wed Aug 24 11:43:14 2011 UTC trunk/Sources/phylmd/physiq.f revision 174 by guez, Wed Nov 25 20:14:19 2015 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    
17      use abort_gcm_m, only: abort_gcm      use aaam_bud_m, only: aaam_bud
18      USE calendar, only: ymds2ju      USE abort_gcm_m, ONLY: abort_gcm
19      use clesphys, only: ecrit_hf, ecrit_ins, ecrit_mth, cdmmax, cdhmax, &      use aeropt_m, only: aeropt
20           co2_ppm, ecrit_reg, ecrit_tra, ksta, ksta_ter, ok_kzmin      use ajsec_m, only: ajsec
21      use clesphys2, only: iflag_con, ok_orolf, ok_orodr, nbapp_rad, &      use calltherm_m, only: calltherm
22           cycle_diurne, new_oliq, soil_model      USE clesphys, ONLY: cdhmax, cdmmax, co2_ppm, ecrit_hf, ecrit_ins, &
23      use clmain_m, only: clmain           ecrit_mth, ecrit_reg, ecrit_tra, ksta, ksta_ter, ok_kzmin
24      use comgeomphy      USE clesphys2, ONLY: cycle_diurne, iflag_con, nbapp_rad, new_oliq, &
25      use concvl_m, only: concvl           ok_orodr, ok_orolf
26      use conf_gcm_m, only: raz_date, offline      USE clmain_m, ONLY: clmain
27      use conf_phys_m, only: conf_phys      use clouds_gno_m, only: clouds_gno
28      use ctherm      use comconst, only: dtphys
29      use dimens_m, only: jjm, iim, llm, nqmx      USE comgeomphy, ONLY: airephy
30      use dimphy, only: klon, nbtr      USE concvl_m, ONLY: concvl
31      use dimsoil, only: nsoilmx      USE conf_gcm_m, ONLY: offline, raz_date, day_step, iphysiq
32      use fcttre, only: thermcep, foeew, qsats, qsatl      USE conf_phys_m, ONLY: conf_phys
33      use hgardfou_m, only: hgardfou      use conflx_m, only: conflx
34      USE histcom, only: histsync      USE ctherm, ONLY: iflag_thermals, nsplit_thermals
35      USE histwrite_m, only: histwrite      use diagcld2_m, only: diagcld2
36      use indicesol, only: nbsrf, is_ter, is_lic, is_sic, is_oce, clnsurf, epsfra      use diagetpq_m, only: diagetpq
37      use ini_histhf_m, only: ini_histhf      use diagphy_m, only: diagphy
38      use ini_histday_m, only: ini_histday      USE dimens_m, ONLY: llm, nqmx
39      use ini_histins_m, only: ini_histins      USE dimphy, ONLY: klon
40      use iniprint, only: prt_level      USE dimsoil, ONLY: nsoilmx
41      use oasis_m      use drag_noro_m, only: drag_noro
42      use orbite_m, only: orbite, zenang      use dynetat0_m, only: day_ref, annee_ref
43      use ozonecm_m, only: ozonecm      USE fcttre, ONLY: foeew, qsatl, qsats, thermcep
44      use phyetat0_m, only: phyetat0, rlat, rlon      use fisrtilp_m, only: fisrtilp
45      use phyredem_m, only: phyredem      USE hgardfou_m, ONLY: hgardfou
46      use phystokenc_m, only: phystokenc      USE indicesol, ONLY: clnsurf, epsfra, is_lic, is_oce, is_sic, is_ter, &
47      use phytrac_m, only: phytrac           nbsrf
48      use qcheck_m, only: qcheck      USE ini_histins_m, ONLY: ini_histins
49      use radepsi      use netcdf95, only: NF95_CLOSE
50      use radopt      use newmicro_m, only: newmicro
51      use SUPHEC_M, only: rcpd, rtt, rlvtt, rg, ra, rsigma, retv, romega      USE orbite_m, ONLY: orbite
52      use temps, only: itau_phy, day_ref, annee_ref      USE ozonecm_m, ONLY: ozonecm
53      use yoethf_m      USE phyetat0_m, ONLY: phyetat0, rlat, rlon
54        USE phyredem_m, ONLY: phyredem
55        USE phyredem0_m, ONLY: phyredem0
56        USE phystokenc_m, ONLY: phystokenc
57        USE phytrac_m, ONLY: phytrac
58        USE qcheck_m, ONLY: qcheck
59        use radlwsw_m, only: radlwsw
60        use readsulfate_m, only: readsulfate
61        use readsulfate_preind_m, only: readsulfate_preind
62        use yoegwd, only: sugwd
63        USE suphec_m, ONLY: rcpd, retv, rg, rlvtt, romega, rsigma, rtt
64        USE temps, ONLY: itau_phy
65        use transp_m, only: transp
66        use unit_nml_m, only: unit_nml
67        USE ymds2ju_m, ONLY: ymds2ju
68        USE yoethf_m, ONLY: r2es, rvtmp2
69        use zenang_m, only: zenang
70    
71      ! Variables argument:      logical, intent(in):: lafin ! dernier passage
72    
73      REAL, intent(in):: rdayvrai      integer, intent(in):: dayvrai
74      ! (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
75    
76      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  
77    
78      REAL, intent(in):: paprs(klon, llm+1)      REAL, intent(in):: paprs(:, :) ! (klon, llm + 1)
79      ! (pression pour chaque inter-couche, en Pa)      ! pression pour chaque inter-couche, en Pa
80    
81      REAL, intent(in):: play(klon, llm)      REAL, intent(in):: play(:, :) ! (klon, llm)
82      ! (input pression pour le mileu de chaque couche (en Pa))      ! pression pour le mileu de chaque couche (en Pa)
83    
84      REAL, intent(in):: pphi(klon, llm)      REAL, intent(in):: pphi(:, :) ! (klon, llm)
85      ! (input geopotentiel de chaque couche (g z) (reference sol))      ! géopotentiel de chaque couche (référence sol)
86    
87      REAL pphis(klon) ! input geopotentiel du sol      REAL, intent(in):: pphis(:) ! (klon) géopotentiel du sol
88    
89      REAL, intent(in):: u(klon, llm)      REAL, intent(in):: u(:, :) ! (klon, llm)
90      ! vitesse dans la direction X (de O a E) en m/s      ! vitesse dans la direction X (de O a E) en m/s
       
     REAL, intent(in):: v(klon, llm) ! vitesse Y (de S a N) en m/s  
     REAL t(klon, llm) ! input temperature (K)  
   
     REAL, intent(in):: qx(klon, llm, nqmx)  
     ! (humidité spécifique et fractions massiques des autres traceurs)  
   
     REAL omega(klon, llm) ! input vitesse verticale en Pa/s  
     REAL, intent(out):: d_u(klon, llm) ! tendance physique de "u" (m/s/s)  
     REAL, intent(out):: d_v(klon, llm) ! tendance physique de "v" (m/s/s)  
     REAL, intent(out):: d_t(klon, llm) ! tendance physique de "t" (K/s)  
     REAL d_qx(klon, llm, nqmx) ! output tendance physique de "qx" (kg/kg/s)  
     REAL d_ps(klon) ! output tendance physique de la pression au sol  
91    
92      LOGICAL:: firstcal = .true.      REAL, intent(in):: v(:, :) ! (klon, llm) vitesse Y (de S a N) en m/s
93        REAL, intent(in):: t(:, :) ! (klon, llm) temperature (K)
94    
95      INTEGER nbteta      REAL, intent(in):: qx(:, :, :) ! (klon, llm, nqmx)
96      PARAMETER(nbteta=3)      ! (humidit\'e sp\'ecifique et fractions massiques des autres traceurs)
97    
98      REAL PVteta(klon, nbteta)      REAL, intent(in):: omega(:, :) ! (klon, llm) vitesse verticale en Pa/s
99      ! (output vorticite potentielle a des thetas constantes)      REAL, intent(out):: d_u(:, :) ! (klon, llm) tendance physique de "u" (m s-2)
100        REAL, intent(out):: d_v(:, :) ! (klon, llm) tendance physique de "v" (m s-2)
101        REAL, intent(out):: d_t(:, :) ! (klon, llm) tendance physique de "t" (K/s)
102    
103      LOGICAL ok_cvl ! pour activer le nouveau driver pour convection KE      REAL, intent(out):: d_qx(:, :, :) ! (klon, llm, nqmx)
104      PARAMETER (ok_cvl=.TRUE.)      ! tendance physique de "qx" (s-1)
     LOGICAL ok_gust ! pour activer l'effet des gust sur flux surface  
     PARAMETER (ok_gust=.FALSE.)  
105    
106      LOGICAL check ! Verifier la conservation du modele en eau      ! Local:
     PARAMETER (check=.FALSE.)  
107    
108      LOGICAL, PARAMETER:: ok_stratus=.FALSE.      LOGICAL:: firstcal = .true.
     ! Ajouter artificiellement les stratus  
109    
110      ! Parametres lies au coupleur OASIS:      LOGICAL ok_gust ! pour activer l'effet des gust sur flux surface
111      INTEGER, SAVE :: npas, nexca      PARAMETER (ok_gust = .FALSE.)
     logical rnpb  
     parameter(rnpb=.true.)  
112    
113      character(len=6), save:: ocean      LOGICAL, PARAMETER:: check = .FALSE.
114      ! (type de modèle océan à utiliser: "force" ou "slab" mais pas "couple")      ! Verifier la conservation du modele en eau
115    
116      logical ok_ocean      LOGICAL, PARAMETER:: ok_stratus = .FALSE.
117      SAVE ok_ocean      ! Ajouter artificiellement les stratus
118    
119      ! "slab" ocean      ! "slab" ocean
120      REAL, save:: tslab(klon) ! temperature of ocean slab      REAL, save:: tslab(klon) ! temperature of ocean slab
# Line 122  contains Line 122  contains
122      REAL fluxo(klon) ! flux turbulents ocean-glace de mer      REAL fluxo(klon) ! flux turbulents ocean-glace de mer
123      REAL fluxg(klon) ! flux turbulents ocean-atmosphere      REAL fluxg(klon) ! flux turbulents ocean-atmosphere
124    
125      ! Modele thermique du sol, a activer pour le cycle diurne:      logical:: ok_journe = .false., ok_mensuel = .true., ok_instan = .false.
126      logical, save:: ok_veget      ! sorties journalieres, mensuelles et instantanees dans les
127      LOGICAL, save:: ok_journe ! sortir le fichier journalier      ! fichiers histday, histmth et histins
   
     LOGICAL ok_mensuel ! sortir le fichier mensuel  
   
     LOGICAL ok_instan ! sortir le fichier instantane  
     save ok_instan  
128    
129      LOGICAL ok_region ! sortir le fichier regional      LOGICAL ok_region ! sortir le fichier regional
130      PARAMETER (ok_region=.FALSE.)      PARAMETER (ok_region = .FALSE.)
131    
132      ! pour phsystoke avec thermiques      ! pour phsystoke avec thermiques
133      REAL fm_therm(klon, llm+1)      REAL fm_therm(klon, llm + 1)
134      REAL entr_therm(klon, llm)      REAL entr_therm(klon, llm)
135      real, save:: q2(klon, llm+1, nbsrf)      real, save:: q2(klon, llm + 1, nbsrf)
136    
137      INTEGER ivap ! indice de traceurs pour vapeur d'eau      INTEGER, PARAMETER:: ivap = 1 ! indice de traceur pour vapeur d'eau
138      PARAMETER (ivap=1)      INTEGER, PARAMETER:: iliq = 2 ! indice de traceur pour eau liquide
     INTEGER iliq ! indice de traceurs pour eau liquide  
     PARAMETER (iliq=2)  
139    
140      REAL, save:: t_ancien(klon, llm), q_ancien(klon, llm)      REAL, save:: t_ancien(klon, llm), q_ancien(klon, llm)
141      LOGICAL, save:: ancien_ok      LOGICAL, save:: ancien_ok
# Line 152  contains Line 145  contains
145    
146      real da(klon, llm), phi(klon, llm, llm), mp(klon, llm)      real da(klon, llm), phi(klon, llm, llm), mp(klon, llm)
147    
148      !IM Amip2 PV a theta constante      REAL swdn0(klon, llm + 1), swdn(klon, llm + 1)
149        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)  
150      SAVE swdn0, swdn, swup0, swup      SAVE swdn0, swdn, swup0, swup
151    
152      REAL lwdn0(klon, klevp1), lwdn(klon, klevp1)      REAL lwdn0(klon, llm + 1), lwdn(klon, llm + 1)
153      REAL lwup0(klon, klevp1), lwup(klon, klevp1)      REAL lwup0(klon, llm + 1), lwup(klon, llm + 1)
154      SAVE lwdn0, lwdn, lwup0, lwup      SAVE lwdn0, lwdn, lwup0, lwup
155    
156      !IM Amip2      ! Amip2
157      ! variables a une pression donnee      ! variables a une pression donnee
158    
159      integer nlevSTD      integer nlevSTD
160      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 '/  
161    
162      ! prw: precipitable water      ! prw: precipitable water
163      real prw(klon)      real prw(klon)
# Line 195  contains Line 168  contains
168      REAL flwc(klon, llm), fiwc(klon, llm)      REAL flwc(klon, llm), fiwc(klon, llm)
169    
170      INTEGER kmax, lmax      INTEGER kmax, lmax
171      PARAMETER(kmax=8, lmax=8)      PARAMETER(kmax = 8, lmax = 8)
172      INTEGER kmaxm1, lmaxm1      INTEGER kmaxm1, lmaxm1
173      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  
174    
175      ! Variables propres a la physique      ! Variables propres a la physique
176    
177      INTEGER, save:: radpas      INTEGER, save:: radpas
178      ! (Radiative transfer computations are made every "radpas" call to      ! Radiative transfer computations are made every "radpas" call to
179      ! "physiq".)      ! "physiq".
180    
181      REAL radsol(klon)      REAL radsol(klon)
182      SAVE radsol ! bilan radiatif au sol calcule par code radiatif      SAVE radsol ! bilan radiatif au sol calcule par code radiatif
183    
184      INTEGER, SAVE:: itap ! number of calls to "physiq"      INTEGER:: itap = 0 ! number of calls to "physiq"
185    
186      REAL, save:: ftsol(klon, nbsrf) ! skin temperature of surface fraction      REAL, save:: ftsol(klon, nbsrf) ! skin temperature of surface fraction
187    
188      REAL, save:: ftsoil(klon, nsoilmx, nbsrf)      REAL, save:: ftsoil(klon, nsoilmx, nbsrf)
189      ! soil temperature of surface fraction      ! soil temperature of surface fraction
190    
191      REAL fevap(klon, nbsrf)      REAL, save:: fevap(klon, nbsrf) ! evaporation
     SAVE fevap ! evaporation  
192      REAL fluxlat(klon, nbsrf)      REAL fluxlat(klon, nbsrf)
193      SAVE fluxlat      SAVE fluxlat
194    
195      REAL fqsurf(klon, nbsrf)      REAL, save:: fqsurf(klon, nbsrf)
196      SAVE fqsurf ! humidite de l'air au contact de la surface      ! humidite de l'air au contact de la surface
197    
198      REAL, save:: qsol(klon) ! hauteur d'eau dans le sol      REAL, save:: qsol(klon)
199        ! column-density of water in soil, in kg m-2
200    
201      REAL fsnow(klon, nbsrf)      REAL, save:: fsnow(klon, nbsrf) ! epaisseur neigeuse
202      SAVE fsnow ! epaisseur neigeuse      REAL, save:: falbe(klon, nbsrf) ! albedo visible par type de surface
203    
204      REAL falbe(klon, nbsrf)      ! Param\`etres de l'orographie \`a l'\'echelle sous-maille (OESM) :
     SAVE falbe ! albedo par type de surface  
     REAL falblw(klon, nbsrf)  
     SAVE falblw ! albedo par type de surface  
   
     ! Paramètres de l'orographie à l'échelle sous-maille (OESM) :  
205      REAL, save:: zmea(klon) ! orographie moyenne      REAL, save:: zmea(klon) ! orographie moyenne
206      REAL, save:: zstd(klon) ! deviation standard de l'OESM      REAL, save:: zstd(klon) ! deviation standard de l'OESM
207      REAL, save:: zsig(klon) ! pente de l'OESM      REAL, save:: zsig(klon) ! pente de l'OESM
# Line 302  contains Line 223  contains
223      !KE43      !KE43
224      ! Variables liees a la convection de K. Emanuel (sb):      ! Variables liees a la convection de K. Emanuel (sb):
225    
     REAL bas, top ! cloud base and top levels  
     SAVE bas  
     SAVE top  
   
226      REAL Ma(klon, llm) ! undilute upward mass flux      REAL Ma(klon, llm) ! undilute upward mass flux
227      SAVE Ma      SAVE Ma
228      REAL qcondc(klon, llm) ! in-cld water content from convect      REAL qcondc(klon, llm) ! in-cld water content from convect
229      SAVE qcondc      SAVE qcondc
230      REAL ema_work1(klon, llm), ema_work2(klon, llm)      REAL, save:: sig1(klon, llm), w01(klon, llm)
231      SAVE ema_work1, ema_work2      REAL, save:: wd(klon)
   
     REAL wd(klon) ! sb  
     SAVE wd ! sb  
232    
233      ! Variables locales pour la couche limite (al1):      ! Variables locales pour la couche limite (al1):
234    
# Line 323  contains Line 237  contains
237      REAL cdragh(klon) ! drag coefficient pour T and Q      REAL cdragh(klon) ! drag coefficient pour T and Q
238      REAL cdragm(klon) ! drag coefficient pour vent      REAL cdragm(klon) ! drag coefficient pour vent
239    
240      !AA Pour phytrac      ! Pour phytrac :
241      REAL ycoefh(klon, llm) ! coef d'echange pour phytrac      REAL ycoefh(klon, llm) ! coef d'echange pour phytrac
242      REAL yu1(klon) ! vents dans la premiere couche U      REAL yu1(klon) ! vents dans la premiere couche U
243      REAL yv1(klon) ! vents dans la premiere couche V      REAL yv1(klon) ! vents dans la premiere couche V
# Line 342  contains Line 256  contains
256      REAL frac_impa(klon, llm) ! fractions d'aerosols lessivees (impaction)      REAL frac_impa(klon, llm) ! fractions d'aerosols lessivees (impaction)
257      REAL frac_nucl(klon, llm) ! idem (nucleation)      REAL frac_nucl(klon, llm) ! idem (nucleation)
258    
259      !AA      REAL, save:: rain_fall(klon)
260      REAL rain_fall(klon) ! pluie      ! liquid water mass flux (kg/m2/s), positive down
261      REAL snow_fall(klon) ! neige  
262      save snow_fall, rain_fall      REAL, save:: snow_fall(klon)
263      !IM cf FH pour Tiedtke 080604      ! solid water mass flux (kg/m2/s), positive down
264    
265      REAL rain_tiedtke(klon), snow_tiedtke(klon)      REAL rain_tiedtke(klon), snow_tiedtke(klon)
266    
267      REAL evap(klon), devap(klon) ! evaporation et sa derivee      REAL evap(klon), devap(klon) ! evaporation and its derivative
268      REAL sens(klon), dsens(klon) ! chaleur sensible et sa derivee      REAL sens(klon), dsens(klon) ! chaleur sensible et sa derivee
269      REAL dlw(klon) ! derivee infra rouge      REAL dlw(klon) ! derivee infra rouge
270      SAVE dlw      SAVE dlw
271      REAL bils(klon) ! bilan de chaleur au sol      REAL bils(klon) ! bilan de chaleur au sol
272      REAL fder(klon) ! Derive de flux (sensible et latente)      REAL, save:: fder(klon) ! Derive de flux (sensible et latente)
     save fder  
273      REAL ve(klon) ! integr. verticale du transport meri. de l'energie      REAL ve(klon) ! integr. verticale du transport meri. de l'energie
274      REAL vq(klon) ! integr. verticale du transport meri. de l'eau      REAL vq(klon) ! integr. verticale du transport meri. de l'eau
275      REAL ue(klon) ! integr. verticale du transport zonal de l'energie      REAL ue(klon) ! integr. verticale du transport zonal de l'energie
276      REAL uq(klon) ! integr. verticale du transport zonal de l'eau      REAL uq(klon) ! integr. verticale du transport zonal de l'eau
277    
278      REAL frugs(klon, nbsrf) ! longueur de rugosite      REAL, save:: frugs(klon, nbsrf) ! longueur de rugosite
     save frugs  
279      REAL zxrugs(klon) ! longueur de rugosite      REAL zxrugs(klon) ! longueur de rugosite
280    
281      ! Conditions aux limites      ! Conditions aux limites
282    
283      INTEGER julien      INTEGER julien
   
284      INTEGER, SAVE:: lmt_pas ! number of time steps of "physics" per day      INTEGER, SAVE:: lmt_pas ! number of time steps of "physics" per day
285      REAL pctsrf(klon, nbsrf)      REAL, save:: pctsrf(klon, nbsrf) ! percentage of surface
286      !IM      REAL pctsrf_new(klon, nbsrf) ! pourcentage surfaces issus d'ORCHIDEE
287      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  
   
288      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
289    
290      ! Declaration des procedures appelees      ! Declaration des procedures appelees
291    
     EXTERNAL alboc ! calculer l'albedo sur ocean  
     EXTERNAL ajsec ! ajustement sec  
     !KE43  
     EXTERNAL conema3 ! convect4.3  
     EXTERNAL fisrtilp ! schema de condensation a grande echelle (pluie)  
292      EXTERNAL nuage ! calculer les proprietes radiatives      EXTERNAL nuage ! calculer les proprietes radiatives
     EXTERNAL radlwsw ! rayonnements solaire et infrarouge  
     EXTERNAL transp ! transport total de l'eau et de l'energie  
293    
294      ! Variables locales      ! Variables locales
295    
296      real clwcon(klon, llm), rnebcon(klon, llm)      real, save:: clwcon(klon, llm), rnebcon(klon, llm)
297      real clwcon0(klon, llm), rnebcon0(klon, llm)      real, save:: clwcon0(klon, llm), rnebcon0(klon, llm)
   
     save rnebcon, clwcon  
298    
299      REAL rhcl(klon, llm) ! humiditi relative ciel clair      REAL rhcl(klon, llm) ! humiditi relative ciel clair
300      REAL dialiq(klon, llm) ! eau liquide nuageuse      REAL dialiq(klon, llm) ! eau liquide nuageuse
# Line 418  contains Line 314  contains
314      REAL zxfluxu(klon, llm)      REAL zxfluxu(klon, llm)
315      REAL zxfluxv(klon, llm)      REAL zxfluxv(klon, llm)
316    
317      REAL heat(klon, llm) ! chauffage solaire      ! Le rayonnement n'est pas calcul\'e tous les pas, il faut donc que
318      REAL heat0(klon, llm) ! chauffage solaire ciel clair      ! les variables soient r\'emanentes.
319      REAL cool(klon, llm) ! refroidissement infrarouge      REAL, save:: heat(klon, llm) ! chauffage solaire
320      REAL cool0(klon, llm) ! refroidissement infrarouge ciel clair      REAL, save:: heat0(klon, llm) ! chauffage solaire ciel clair
321      REAL topsw(klon), toplw(klon), solsw(klon), sollw(klon)      REAL, save:: cool(klon, llm) ! refroidissement infrarouge
322      real sollwdown(klon) ! downward LW flux at surface      REAL, save:: cool0(klon, llm) ! refroidissement infrarouge ciel clair
323      REAL topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)      REAL, save:: topsw(klon), toplw(klon), solsw(klon)
324      REAL albpla(klon)      REAL, save:: sollw(klon) ! rayonnement infrarouge montant \`a la surface
325        real, save:: sollwdown(klon) ! downward LW flux at surface
326        REAL, save:: topsw0(klon), toplw0(klon), solsw0(klon), sollw0(klon)
327        REAL, save:: albpla(klon)
328      REAL fsollw(klon, nbsrf) ! bilan flux IR pour chaque sous surface      REAL fsollw(klon, nbsrf) ! bilan flux IR pour chaque sous surface
329      REAL fsolsw(klon, nbsrf) ! flux solaire absorb. pour chaque sous surface      REAL fsolsw(klon, nbsrf) ! flux solaire absorb. pour chaque sous surface
     ! Le rayonnement n'est pas calcule tous les pas, il faut donc  
     ! sauvegarder les sorties du rayonnement  
     SAVE heat, cool, albpla, topsw, toplw, solsw, sollw, sollwdown  
     SAVE topsw0, toplw0, solsw0, sollw0, heat0, cool0  
   
     INTEGER itaprad  
     SAVE itaprad  
330    
331      REAL conv_q(klon, llm) ! convergence de l'humidite (kg/kg/s)      REAL conv_q(klon, llm) ! convergence de l'humidite (kg/kg/s)
332      REAL conv_t(klon, llm) ! convergence of temperature (K/s)      REAL conv_t(klon, llm) ! convergence of temperature (K/s)
# Line 444  contains Line 336  contains
336    
337      REAL zxtsol(klon), zxqsurf(klon), zxsnow(klon), zxfluxlat(klon)      REAL zxtsol(klon), zxqsurf(klon), zxsnow(klon), zxfluxlat(klon)
338    
339      REAL dist, rmu0(klon), fract(klon)      REAL dist, mu0(klon), fract(klon)
340      REAL zdtime ! pas de temps du rayonnement (s)      real longi
     real zlongi  
   
341      REAL z_avant(klon), z_apres(klon), z_factor(klon)      REAL z_avant(klon), z_apres(klon), z_factor(klon)
     LOGICAL zx_ajustq  
   
342      REAL za, zb      REAL za, zb
343      REAL zx_t, zx_qs, zdelta, zcor, zlvdcp, zlsdcp      REAL zx_t, zx_qs, zcor
344      real zqsat(klon, llm)      real zqsat(klon, llm)
345      INTEGER i, k, iq, nsrf      INTEGER i, k, iq, nsrf
346      REAL t_coup      REAL, PARAMETER:: t_coup = 234.
     PARAMETER (t_coup=234.0)  
   
347      REAL zphi(klon, llm)      REAL zphi(klon, llm)
348    
349      !IM cf. AM Variables locales pour la CLA (hbtm2)      ! cf. AM Variables locales pour la CLA (hbtm2)
350    
351      REAL, SAVE:: pblh(klon, nbsrf) ! Hauteur de couche limite      REAL, SAVE:: pblh(klon, nbsrf) ! Hauteur de couche limite
352      REAL, SAVE:: plcl(klon, nbsrf) ! Niveau de condensation de la CLA      REAL, SAVE:: plcl(klon, nbsrf) ! Niveau de condensation de la CLA
# Line 478  contains Line 364  contains
364      REAL s_therm(klon), s_trmb1(klon), s_trmb2(klon)      REAL s_therm(klon), s_trmb1(klon), s_trmb2(klon)
365      REAL s_trmb3(klon)      REAL s_trmb3(klon)
366    
367      ! Variables locales pour la convection de K. Emanuel (sb):      ! Variables locales pour la convection de K. Emanuel :
368    
369      REAL upwd(klon, llm) ! saturated updraft mass flux      REAL upwd(klon, llm) ! saturated updraft mass flux
370      REAL dnwd(klon, llm) ! saturated downdraft mass flux      REAL dnwd(klon, llm) ! saturated downdraft mass flux
371      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  
372      REAL cape(klon) ! CAPE      REAL cape(klon) ! CAPE
373      SAVE cape      SAVE cape
374    
     REAL pbase(klon) ! cloud base pressure  
     SAVE pbase  
     REAL bbase(klon) ! cloud base buoyancy  
     SAVE bbase  
     REAL rflag(klon) ! flag fonctionnement de convect  
375      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)  
376    
377      ! Variables du changement      ! Variables du changement
378    
379      ! con: convection      ! con: convection
380      ! lsc: condensation a grande echelle (Large-Scale-Condensation)      ! lsc: large scale condensation
381      ! ajs: ajustement sec      ! ajs: ajustement sec
382      ! eva: evaporation de l'eau liquide nuageuse      ! eva: \'evaporation de l'eau liquide nuageuse
383      ! vdf: couche limite (Vertical DiFfusion)      ! vdf: vertical diffusion in boundary layer
384      REAL d_t_con(klon, llm), d_q_con(klon, llm)      REAL d_t_con(klon, llm), d_q_con(klon, llm)
385      REAL d_u_con(klon, llm), d_v_con(klon, llm)      REAL d_u_con(klon, llm), d_v_con(klon, llm)
386      REAL d_t_lsc(klon, llm), d_q_lsc(klon, llm), d_ql_lsc(klon, llm)      REAL d_t_lsc(klon, llm), d_q_lsc(klon, llm), d_ql_lsc(klon, llm)
# Line 512  contains Line 388  contains
388      REAL d_u_ajs(klon, llm), d_v_ajs(klon, llm)      REAL d_u_ajs(klon, llm), d_v_ajs(klon, llm)
389      REAL rneb(klon, llm)      REAL rneb(klon, llm)
390    
391      REAL pmfu(klon, llm), pmfd(klon, llm)      REAL mfu(klon, llm), mfd(klon, llm)
392      REAL pen_u(klon, llm), pen_d(klon, llm)      REAL pen_u(klon, llm), pen_d(klon, llm)
393      REAL pde_u(klon, llm), pde_d(klon, llm)      REAL pde_u(klon, llm), pde_d(klon, llm)
394      INTEGER kcbot(klon), kctop(klon), kdtop(klon)      INTEGER kcbot(klon), kctop(klon), kdtop(klon)
395      REAL pmflxr(klon, llm+1), pmflxs(klon, llm+1)      REAL pmflxr(klon, llm + 1), pmflxs(klon, llm + 1)
396      REAL prfl(klon, llm+1), psfl(klon, llm+1)      REAL prfl(klon, llm + 1), psfl(klon, llm + 1)
   
     INTEGER ibas_con(klon), itop_con(klon)  
397    
398      SAVE ibas_con, itop_con      INTEGER, save:: ibas_con(klon), itop_con(klon)
399    
400      REAL rain_con(klon), rain_lsc(klon)      REAL rain_con(klon), rain_lsc(klon)
401      REAL snow_con(klon), snow_lsc(klon)      REAL snow_con(klon), snow_lsc(klon)
# Line 535  contains Line 409  contains
409      REAL d_u_lif(klon, llm), d_v_lif(klon, llm)      REAL d_u_lif(klon, llm), d_v_lif(klon, llm)
410      REAL d_t_lif(klon, llm)      REAL d_t_lif(klon, llm)
411    
412      REAL ratqs(klon, llm), ratqss(klon, llm), ratqsc(klon, llm)      REAL, save:: ratqs(klon, llm)
413      real ratqsbas, ratqshaut      real ratqss(klon, llm), ratqsc(klon, llm)
414      save ratqsbas, ratqshaut, ratqs      real:: ratqsbas = 0.01, ratqshaut = 0.3
415    
416      ! Parametres lies au nouveau schema de nuages (SB, PDF)      ! Parametres lies au nouveau schema de nuages (SB, PDF)
417      real, save:: fact_cldcon      real:: fact_cldcon = 0.375
418      real, save:: facttemps      real:: facttemps = 1.e-4
419      logical ok_newmicro      logical:: ok_newmicro = .true.
     save ok_newmicro  
420      real facteur      real facteur
421    
422      integer iflag_cldcon      integer:: iflag_cldcon = 1
     save iflag_cldcon  
   
423      logical ptconv(klon, llm)      logical ptconv(klon, llm)
424    
425      ! Variables locales pour effectuer les appels en série      ! Variables locales pour effectuer les appels en s\'erie :
426    
427      REAL t_seri(klon, llm), q_seri(klon, llm)      REAL t_seri(klon, llm), q_seri(klon, llm)
428      REAL ql_seri(klon, llm), qs_seri(klon, llm)      REAL ql_seri(klon, llm)
429      REAL u_seri(klon, llm), v_seri(klon, llm)      REAL u_seri(klon, llm), v_seri(klon, llm)
430        REAL tr_seri(klon, llm, nqmx - 2)
     REAL tr_seri(klon, llm, nbtr)  
     REAL d_tr(klon, llm, nbtr)  
431    
432      REAL zx_rh(klon, llm)      REAL zx_rh(klon, llm)
433    
# Line 567  contains Line 436  contains
436      REAL zustrph(klon), zvstrph(klon)      REAL zustrph(klon), zvstrph(klon)
437      REAL aam, torsfc      REAL aam, torsfc
438    
     REAL dudyn(iim+1, jjm + 1, llm)  
   
439      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)  
440    
441      INTEGER, SAVE:: nid_day, nid_ins      INTEGER, SAVE:: nid_ins
442    
443      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.
444      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.
445      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.
446      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.
447    
     REAL zsto  
   
     character(len=20) modname  
     character(len=80) abort_message  
     logical ok_sync  
448      real date0      real date0
449    
450      ! Variables liees au bilan d'energie et d'enthalpi      ! Variables li\'ees au bilan d'\'energie et d'enthalpie :
451      REAL ztsol(klon)      REAL ztsol(klon)
452      REAL d_h_vcol, d_qt, d_qw, d_ql, d_qs, d_ec      REAL d_h_vcol, d_qt, d_ec
453      REAL d_h_vcol_phy      REAL, SAVE:: d_h_vcol_phy
     REAL fs_bound, fq_bound  
     SAVE d_h_vcol_phy  
454      REAL zero_v(klon)      REAL zero_v(klon)
455      CHARACTER(LEN=15) ztit      CHARACTER(LEN = 20) tit
456      INTEGER ip_ebil ! PRINT level for energy conserv. diag.      INTEGER:: ip_ebil = 0 ! print level for energy conservation diagnostics
457      SAVE ip_ebil      INTEGER:: if_ebil = 0 ! verbosity for diagnostics of energy conservation
458      DATA ip_ebil/0/  
459      INTEGER, SAVE:: if_ebil ! level for energy conservation diagnostics      REAL d_t_ec(klon, llm) ! tendance due \`a la conversion Ec -> E thermique
     !+jld ec_conser  
     REAL d_t_ec(klon, llm) ! tendance du a la conersion Ec -> E thermique  
460      REAL ZRCPD      REAL ZRCPD
461      !-jld ec_conser  
     !IM: t2m, q2m, u10m, v10m  
462      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
463      REAL u10m(klon, nbsrf), v10m(klon, nbsrf) !vents a 10m      REAL u10m(klon, nbsrf), v10m(klon, nbsrf) ! vents a 10 m
464      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
465      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
466      !jq Aerosol effects (Johannes Quaas, 27/11/2003)  
467      REAL sulfate(klon, llm) ! SO4 aerosol concentration [ug/m3]      ! Aerosol effects:
468    
469        REAL sulfate(klon, llm) ! SO4 aerosol concentration (micro g/m3)
470    
471      REAL, save:: sulfate_pi(klon, llm)      REAL, save:: sulfate_pi(klon, llm)
472      ! (SO4 aerosol concentration, in ug/m3, pre-industrial value)      ! SO4 aerosol concentration, in micro g/m3, pre-industrial value
473    
474      REAL cldtaupi(klon, llm)      REAL cldtaupi(klon, llm)
475      ! (Cloud optical thickness for pre-industrial (pi) aerosols)      ! cloud optical thickness for pre-industrial (pi) aerosols
476    
477      REAL re(klon, llm) ! Cloud droplet effective radius      REAL re(klon, llm) ! Cloud droplet effective radius
478      REAL fl(klon, llm) ! denominator of re      REAL fl(klon, llm) ! denominator of re
479    
480      ! Aerosol optical properties      ! Aerosol optical properties
481      REAL tau_ae(klon, llm, 2), piz_ae(klon, llm, 2)      REAL, save:: tau_ae(klon, llm, 2), piz_ae(klon, llm, 2)
482      REAL cg_ae(klon, llm, 2)      REAL, save:: cg_ae(klon, llm, 2)
   
     REAL topswad(klon), solswad(klon) ! Aerosol direct effect.  
     ! ok_ade=True -ADE=topswad-topsw  
483    
484      REAL topswai(klon), solswai(klon) ! Aerosol indirect effect.      REAL topswad(klon), solswad(klon) ! aerosol direct effect
485      ! ok_aie=True ->      REAL topswai(klon), solswai(klon) ! aerosol indirect effect
     ! ok_ade=True -AIE=topswai-topswad  
     ! ok_ade=F -AIE=topswai-topsw  
486    
487      REAL aerindex(klon) ! POLDER aerosol index      REAL aerindex(klon) ! POLDER aerosol index
488    
489      ! Parameters      LOGICAL:: ok_ade = .false. ! apply aerosol direct effect
490      LOGICAL ok_ade, ok_aie ! Apply aerosol (in)direct effects or not      LOGICAL:: ok_aie = .false. ! apply aerosol indirect effect
491      REAL bl95_b0, bl95_b1 ! Parameter in Boucher and Lohmann (1995)  
492        REAL:: bl95_b0 = 2., bl95_b1 = 0.2
493        ! Parameters in equation (D) of Boucher and Lohmann (1995, Tellus
494        ! B). They link cloud droplet number concentration to aerosol mass
495        ! concentration.
496    
     SAVE ok_ade, ok_aie, bl95_b0, bl95_b1  
497      SAVE u10m      SAVE u10m
498      SAVE v10m      SAVE v10m
499      SAVE t2m      SAVE t2m
500      SAVE q2m      SAVE q2m
501      SAVE ffonte      SAVE ffonte
502      SAVE fqcalving      SAVE fqcalving
     SAVE piz_ae  
     SAVE tau_ae  
     SAVE cg_ae  
503      SAVE rain_con      SAVE rain_con
504      SAVE snow_con      SAVE snow_con
505      SAVE topswai      SAVE topswai
# Line 655  contains Line 508  contains
508      SAVE solswad      SAVE solswad
509      SAVE d_u_con      SAVE d_u_con
510      SAVE d_v_con      SAVE d_v_con
     SAVE rnebcon0  
     SAVE clwcon0  
511    
512      real zmasse(klon, llm)      real zmasse(klon, llm)
513      ! (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)
514    
515      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
516        integer, save:: ncid_startphy
517    
518        namelist /physiq_nml/ ok_journe, ok_mensuel, ok_instan, fact_cldcon, &
519             facttemps, ok_newmicro, iflag_cldcon, ratqsbas, ratqshaut, if_ebil, &
520             ok_ade, ok_aie, bl95_b0, bl95_b1, iflag_thermals, nsplit_thermals
521    
522      !----------------------------------------------------------------      !----------------------------------------------------------------
523    
524      modname = 'physiq'      IF (if_ebil >= 1) zero_v = 0.
525      IF (if_ebil >= 1) THEN      IF (nqmx < 2) CALL abort_gcm('physiq', &
526         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  
527    
528      test_firstcal: IF (firstcal) THEN      test_firstcal: IF (firstcal) THEN
529         ! initialiser         ! initialiser
530         u10m=0.         u10m = 0.
531         v10m=0.         v10m = 0.
532         t2m=0.         t2m = 0.
533         q2m=0.         q2m = 0.
534         ffonte=0.         ffonte = 0.
535         fqcalving=0.         fqcalving = 0.
536         piz_ae=0.         piz_ae = 0.
537         tau_ae=0.         tau_ae = 0.
538         cg_ae=0.         cg_ae = 0.
539         rain_con(:)=0.         rain_con = 0.
540         snow_con(:)=0.         snow_con = 0.
541         bl95_b0=0.         topswai = 0.
542         bl95_b1=0.         topswad = 0.
543         topswai(:)=0.         solswai = 0.
544         topswad(:)=0.         solswad = 0.
545         solswai(:)=0.  
546         solswad(:)=0.         d_u_con = 0.
547           d_v_con = 0.
548         d_u_con = 0.0         rnebcon0 = 0.
549         d_v_con = 0.0         clwcon0 = 0.
550         rnebcon0 = 0.0         rnebcon = 0.
551         clwcon0 = 0.0         clwcon = 0.
        rnebcon = 0.0  
        clwcon = 0.0  
552    
553         pblh =0. ! Hauteur de couche limite         pblh =0. ! Hauteur de couche limite
554         plcl =0. ! Niveau de condensation de la CLA         plcl =0. ! Niveau de condensation de la CLA
# Line 715  contains Line 561  contains
561         trmb2 =0. ! inhibition         trmb2 =0. ! inhibition
562         trmb3 =0. ! Point Omega         trmb3 =0. ! Point Omega
563    
564         IF (if_ebil >= 1) d_h_vcol_phy=0.         IF (if_ebil >= 1) d_h_vcol_phy = 0.
565    
566         ! appel a la lecture du run.def physique         iflag_thermals = 0
567           nsplit_thermals = 1
568           print *, "Enter namelist 'physiq_nml'."
569           read(unit=*, nml=physiq_nml)
570           write(unit_nml, nml=physiq_nml)
571    
572         call conf_phys(ocean, ok_veget, ok_journe, ok_mensuel, &         call conf_phys
             ok_instan, fact_cldcon, facttemps, ok_newmicro, &  
             iflag_cldcon, ratqsbas, ratqshaut, if_ebil, &  
             ok_ade, ok_aie, &  
             bl95_b0, bl95_b1, &  
             iflag_thermals, nsplit_thermals)  
573    
574         ! Initialiser les compteurs:         ! Initialiser les compteurs:
575    
576         frugs = 0.         frugs = 0.
577         itap = 0         CALL phyetat0(pctsrf, ftsol, ftsoil, tslab, seaice, fqsurf, qsol, &
578         itaprad = 0              fsnow, falbe, fevap, rain_fall, snow_fall, solsw, sollw, dlw, &
579         CALL phyetat0("startphy.nc", pctsrf, ftsol, ftsoil, ocean, tslab, &              radsol, frugs, agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, &
580              seaice, fqsurf, qsol, fsnow, falbe, falblw, fevap, rain_fall, &              t_ancien, q_ancien, ancien_ok, rnebcon, ratqs, clwcon, &
581              snow_fall, solsw, sollwdown, dlw, radsol, frugs, agesno, zmea, &              run_off_lic_0, sig1, w01, ncid_startphy)
             zstd, zsig, zgam, zthe, zpic, zval, t_ancien, q_ancien, &  
             ancien_ok, rnebcon, ratqs, clwcon, run_off_lic_0)  
582    
583         ! ATTENTION : il faudra a terme relire q2 dans l'etat initial         ! ATTENTION : il faudra a terme relire q2 dans l'etat initial
584         q2=1.e-8         q2 = 1e-8
   
        radpas = NINT(86400. / dtphys / nbapp_rad)  
   
        ! on remet le calendrier a zero  
        IF (raz_date) itau_phy = 0  
585    
586         PRINT *, 'cycle_diurne = ', cycle_diurne         lmt_pas = day_step / iphysiq
587           print *, 'Number of time steps of "physics" per day: ', lmt_pas
588    
589         IF(ocean.NE.'force ') THEN         radpas = lmt_pas / nbapp_rad
           ok_ocean=.TRUE.  
        ENDIF  
590    
591         CALL printflag(radpas, ok_ocean, ok_oasis, ok_journe, ok_instan, &         ! On remet le calendrier a zero
592              ok_region)         IF (raz_date) itau_phy = 0
593    
594         IF (dtphys*REAL(radpas).GT.21600..AND.cycle_diurne) THEN         CALL printflag(radpas, ok_journe, ok_instan, ok_region)
           print *,'Nbre d appels au rayonnement insuffisant'  
           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  
        print *,"Clef pour le driver de la convection, ok_cvl=", &  
             ok_cvl  
595    
596         ! Initialisation pour la convection de K.E. (sb):         ! Initialisation pour le sch\'ema de convection d'Emanuel :
597         IF (iflag_con >= 3) THEN         IF (iflag_con >= 3) THEN
598              ibas_con = 1
599            print *,"*** Convection de Kerry Emanuel 4.3 "            itop_con = 1
   
           !IM15/11/02 rajout initialisation ibas_con, itop_con cf. SB =>BEG  
           DO i = 1, klon  
              ibas_con(i) = 1  
              itop_con(i) = 1  
           ENDDO  
           !IM15/11/02 rajout initialisation ibas_con, itop_con cf. SB =>END  
   
600         ENDIF         ENDIF
601    
602         IF (ok_orodr) THEN         IF (ok_orodr) THEN
603            rugoro = MAX(1e-5, zstd * zsig / 2)            rugoro = MAX(1e-5, zstd * zsig / 2)
604            CALL SUGWD(klon, llm, paprs, play)            CALL SUGWD(paprs, play)
605         else         else
606            rugoro = 0.            rugoro = 0.
607         ENDIF         ENDIF
608    
        lmt_pas = NINT(86400. / dtphys) ! tous les jours  
        print *, 'Number of time steps of "physics" per day: ', lmt_pas  
   
609         ecrit_ins = NINT(ecrit_ins/dtphys)         ecrit_ins = NINT(ecrit_ins/dtphys)
610         ecrit_hf = NINT(ecrit_hf/dtphys)         ecrit_hf = NINT(ecrit_hf/dtphys)
611         ecrit_mth = NINT(ecrit_mth/dtphys)         ecrit_mth = NINT(ecrit_mth/dtphys)
612         ecrit_tra = NINT(86400.*ecrit_tra/dtphys)         ecrit_tra = NINT(86400.*ecrit_tra/dtphys)
613         ecrit_reg = NINT(ecrit_reg/dtphys)         ecrit_reg = NINT(ecrit_reg/dtphys)
614    
        ! Initialiser le couplage si necessaire  
   
        npas = 0  
        nexca = 0  
   
        print *,'AVANT HIST IFLAG_CON=', iflag_con  
   
615         ! Initialisation des sorties         ! Initialisation des sorties
616    
        call ini_histhf(dtphys, nid_hf, nid_hf3d)  
        call ini_histday(dtphys, ok_journe, nid_day, nqmx)  
617         call ini_histins(dtphys, ok_instan, nid_ins)         call ini_histins(dtphys, ok_instan, nid_ins)
618         CALL ymds2ju(annee_ref, 1, int(day_ref), 0., date0)         CALL ymds2ju(annee_ref, 1, day_ref, 0., date0)
619         !XXXPB Positionner date0 pour initialisation de ORCHIDEE         ! Positionner date0 pour initialisation de ORCHIDEE
620         WRITE(*, *) 'physiq date0 : ', date0         print *, 'physiq date0: ', date0
621           CALL phyredem0(lmt_pas)
622      ENDIF test_firstcal      ENDIF test_firstcal
623    
624      ! Mettre a zero des variables de sortie (pour securite)      ! We will modify variables *_seri and we will not touch variables
625        ! u, v, t, qx:
626        t_seri = t
627        u_seri = u
628        v_seri = v
629        q_seri = qx(:, :, ivap)
630        ql_seri = qx(:, :, iliq)
631        tr_seri = qx(:, :, 3:nqmx)
632    
633      DO i = 1, klon      ztsol = sum(ftsol * pctsrf, dim = 2)
        d_ps(i) = 0.0  
     ENDDO  
     DO iq = 1, nqmx  
        DO k = 1, llm  
           DO i = 1, klon  
              d_qx(i, k, iq) = 0.0  
           ENDDO  
        ENDDO  
     ENDDO  
     da=0.  
     mp=0.  
     phi=0.  
   
     ! Ne pas affecter les valeurs entrees de u, v, h, et q  
   
     DO k = 1, llm  
        DO i = 1, klon  
           t_seri(i, k) = t(i, k)  
           u_seri(i, k) = u(i, k)  
           v_seri(i, k) = v(i, k)  
           q_seri(i, k) = qx(i, k, ivap)  
           ql_seri(i, k) = qx(i, k, iliq)  
           qs_seri(i, k) = 0.  
        ENDDO  
     ENDDO  
     IF (nqmx >= 3) THEN  
        tr_seri(:, :, :nqmx-2) = qx(:, :, 3:nqmx)  
     ELSE  
        tr_seri(:, :, 1) = 0.  
     ENDIF  
   
     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  
634    
635      IF (if_ebil >= 1) THEN      IF (if_ebil >= 1) THEN
636         ztit='after dynamic'         tit = 'after dynamics'
637         CALL diagetpq(airephy, ztit, ip_ebil, 1, 1, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 1, 1, dtphys, t_seri, q_seri, &
638              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)
639              d_ql, d_qs, d_ec)         ! Comme les tendances de la physique sont ajout\'es dans la
640         ! Comme les tendances de la physique sont ajoute dans la dynamique,         !  dynamique, la variation d'enthalpie par la dynamique devrait
641         ! on devrait avoir que la variation d'entalpie par la dynamique         !  \^etre \'egale \`a la variation de la physique au pas de temps
642         ! est egale a la variation de la physique au pas de temps precedent.         !  pr\'ec\'edent.  Donc la somme de ces 2 variations devrait \^etre
643         ! Donc la somme de ces 2 variations devrait etre nulle.         !  nulle.
644         call diagphy(airephy, ztit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &         call diagphy(airephy, tit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
645              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, &
646              d_qt, 0., fs_bound, fq_bound)              d_qt, 0.)
647      END IF      END IF
648    
649      ! Diagnostiquer la tendance dynamique      ! Diagnostic de la tendance dynamique :
650      IF (ancien_ok) THEN      IF (ancien_ok) THEN
651         DO k = 1, llm         DO k = 1, llm
652            DO i = 1, klon            DO i = 1, klon
# Line 879  contains Line 657  contains
657      ELSE      ELSE
658         DO k = 1, llm         DO k = 1, llm
659            DO i = 1, klon            DO i = 1, klon
660               d_t_dyn(i, k) = 0.0               d_t_dyn(i, k) = 0.
661               d_q_dyn(i, k) = 0.0               d_q_dyn(i, k) = 0.
662            ENDDO            ENDDO
663         ENDDO         ENDDO
664         ancien_ok = .TRUE.         ancien_ok = .TRUE.
# Line 896  contains Line 674  contains
674      ! Check temperatures:      ! Check temperatures:
675      CALL hgardfou(t_seri, ftsol)      CALL hgardfou(t_seri, ftsol)
676    
677      ! Incrementer le compteur de la physique      ! Incrémenter le compteur de la physique
678      itap = itap + 1      itap = itap + 1
679      julien = MOD(NINT(rdayvrai), 360)      julien = MOD(dayvrai, 360)
680      if (julien == 0) julien = 360      if (julien == 0) julien = 360
681    
682      forall (k = 1: llm) zmasse(:, k) = (paprs(:, k)-paprs(:, k+1)) / rg      forall (k = 1: llm) zmasse(:, k) = (paprs(:, k) - paprs(:, k + 1)) / rg
683    
684      ! Mettre en action les conditions aux limites (albedo, sst, etc.).      ! Prescrire l'ozone :
685        wo = ozonecm(REAL(julien), paprs)
686    
687      ! Prescrire l'ozone et calculer l'albedo sur l'ocean.      ! \'Evaporation de l'eau liquide nuageuse :
688      if (nqmx >= 5) then      DO k = 1, llm
        wo = qx(:, :, 5) * zmasse / dobson_u / 1e3  
     else IF (MOD(itap - 1, lmt_pas) == 0) THEN  
        wo = ozonecm(REAL(julien), paprs)  
     ENDIF  
   
     ! Re-evaporer l'eau liquide nuageuse  
   
     DO k = 1, llm ! re-evaporation de l'eau liquide nuageuse  
689         DO i = 1, klon         DO i = 1, klon
690            zlvdcp=RLVTT/RCPD/(1.0+RVTMP2*q_seri(i, k))            zb = MAX(0., ql_seri(i, k))
691            zlsdcp=RLVTT/RCPD/(1.0+RVTMP2*q_seri(i, k))            t_seri(i, k) = t_seri(i, k) &
692            zdelta = MAX(0., SIGN(1., RTT-t_seri(i, k)))                 - zb * RLVTT / RCPD / (1. + RVTMP2 * q_seri(i, k))
           zb = MAX(0.0, ql_seri(i, k))  
           za = - MAX(0.0, ql_seri(i, k)) &  
                * (zlvdcp*(1.-zdelta)+zlsdcp*zdelta)  
           t_seri(i, k) = t_seri(i, k) + za  
693            q_seri(i, k) = q_seri(i, k) + zb            q_seri(i, k) = q_seri(i, k) + zb
           ql_seri(i, k) = 0.0  
694         ENDDO         ENDDO
695      ENDDO      ENDDO
696        ql_seri = 0.
697    
698      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
699         ztit='after reevap'         tit = 'after reevap'
700         CALL diagetpq(airephy, ztit, ip_ebil, 2, 1, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 2, 1, dtphys, t_seri, q_seri, &
701              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)
702              d_ql, d_qs, d_ec)         call diagphy(airephy, tit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
703         call diagphy(airephy, ztit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &              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)  
   
704      END IF      END IF
705    
706      ! Appeler la diffusion verticale (programme de couche limite)      frugs = MAX(frugs, 0.000015)
707        zxrugs = sum(frugs * pctsrf, dim = 2)
708    
709      DO i = 1, klon      ! Calculs nécessaires au calcul de l'albedo dans l'interface avec
710         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  
   
     ! calculs necessaires au calcul de l'albedo dans l'interface  
711    
712      CALL orbite(REAL(julien), zlongi, dist)      CALL orbite(REAL(julien), longi, dist)
713      IF (cycle_diurne) THEN      IF (cycle_diurne) THEN
714         zdtime = dtphys * REAL(radpas)         CALL zenang(longi, time, dtphys * radpas, mu0, fract)
        CALL zenang(zlongi, time, zdtime, rmu0, fract)  
715      ELSE      ELSE
716         rmu0 = -999.999         mu0 = - 999.999
717      ENDIF      ENDIF
718    
719      ! Calcul de l'abedo moyen par maille      ! Calcul de l'abedo moyen par maille
720      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  
721    
722      ! Repartition sous maille des flux LW et SW      ! R\'epartition sous maille des flux longwave et shortwave
723      ! Repartition du longwave par sous-surface linearisee      ! R\'epartition du longwave par sous-surface lin\'earis\'ee
724    
725      DO nsrf = 1, nbsrf      forall (nsrf = 1: nbsrf)
726         DO i = 1, klon         fsollw(:, nsrf) = sollw + 4. * RSIGMA * ztsol**3 &
727            fsollw(i, nsrf) = sollw(i) &              * (ztsol - ftsol(:, nsrf))
728                 + 4.0*RSIGMA*ztsol(i)**3 * (ztsol(i)-ftsol(i, nsrf))         fsolsw(:, nsrf) = solsw * (1. - falbe(:, nsrf)) / (1. - albsol)
729            fsolsw(i, nsrf) = solsw(i)*(1.-falbe(i, nsrf))/(1.-albsol(i))      END forall
        ENDDO  
     ENDDO  
730    
731      fder = dlw      fder = dlw
732    
733      ! Couche limite:      ! Couche limite:
734    
735      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, &
736           u_seri, v_seri, julien, rmu0, co2_ppm, ok_veget, ocean, npas, nexca, &           v_seri, julien, mu0, co2_ppm, ftsol, cdmmax, cdhmax, ksta, ksta_ter, &
737           ftsol, soil_model, cdmmax, cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, &           ok_kzmin, ftsoil, qsol, paprs, play, fsnow, fqsurf, fevap, falbe, &
738           qsol, paprs, play, fsnow, fqsurf, fevap, falbe, falblw, fluxlat, &           fluxlat, rain_fall, snow_fall, fsolsw, fsollw, fder, rlat, frugs, &
739           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, &
740           cuphy, cvphy, frugs, firstcal, lafin, agesno, rugoro, d_t_vdf, &           fluxt, fluxq, fluxu, fluxv, cdragh, cdragm, q2, dsens, devap, &
741           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, &
742           cdragm, q2, dsens, devap, ycoefh, yu1, yv1, t2m, q2m, u10m, v10m, &           pblT, therm, trmb1, trmb2, trmb3, plcl, fqcalving, ffonte, &
743           pblh, capCL, oliqCL, cteiCL, pblT, therm, trmb1, trmb2, trmb3, plcl, &           run_off_lic_0, fluxo, fluxg, tslab)
744           fqcalving, ffonte, run_off_lic_0, fluxo, fluxg, tslab, seaice)  
745        ! Incr\'ementation des flux
746      ! Incrémentation des flux  
747        zxfluxt = 0.
748      zxfluxt=0.      zxfluxq = 0.
749      zxfluxq=0.      zxfluxu = 0.
750      zxfluxu=0.      zxfluxv = 0.
     zxfluxv=0.  
751      DO nsrf = 1, nbsrf      DO nsrf = 1, nbsrf
752         DO k = 1, llm         DO k = 1, llm
753            DO i = 1, klon            DO i = 1, klon
754               zxfluxt(i, k) = zxfluxt(i, k) + &               zxfluxt(i, k) = zxfluxt(i, k) + fluxt(i, k, nsrf) * pctsrf(i, nsrf)
755                    fluxt(i, k, nsrf) * pctsrf(i, nsrf)               zxfluxq(i, k) = zxfluxq(i, k) + fluxq(i, k, nsrf) * pctsrf(i, nsrf)
756               zxfluxq(i, k) = zxfluxq(i, k) + &               zxfluxu(i, k) = zxfluxu(i, k) + fluxu(i, k, nsrf) * pctsrf(i, nsrf)
757                    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)  
758            END DO            END DO
759         END DO         END DO
760      END DO      END DO
761      DO i = 1, klon      DO i = 1, klon
762         sens(i) = - zxfluxt(i, 1) ! flux de chaleur sensible au sol         sens(i) = - zxfluxt(i, 1) ! flux de chaleur sensible au sol
763         evap(i) = - zxfluxq(i, 1) ! flux d'evaporation au sol         evap(i) = - zxfluxq(i, 1) ! flux d'\'evaporation au sol
764         fder(i) = dlw(i) + dsens(i) + devap(i)         fder(i) = dlw(i) + dsens(i) + devap(i)
765      ENDDO      ENDDO
766    
# Line 1037  contains Line 774  contains
774      ENDDO      ENDDO
775    
776      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
777         ztit='after clmain'         tit = 'after clmain'
778         CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
779              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)
780              d_ql, d_qs, d_ec)         call diagphy(airephy, tit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
781         call diagphy(airephy, ztit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &              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)  
782      END IF      END IF
783    
784      ! Update surface temperature:      ! Update surface temperature:
785    
786      DO i = 1, klon      DO i = 1, klon
787         zxtsol(i) = 0.0         zxtsol(i) = 0.
788         zxfluxlat(i) = 0.0         zxfluxlat(i) = 0.
789    
790         zt2m(i) = 0.0         zt2m(i) = 0.
791         zq2m(i) = 0.0         zq2m(i) = 0.
792         zu10m(i) = 0.0         zu10m(i) = 0.
793         zv10m(i) = 0.0         zv10m(i) = 0.
794         zxffonte(i) = 0.0         zxffonte(i) = 0.
795         zxfqcalving(i) = 0.0         zxfqcalving(i) = 0.
796    
797         s_pblh(i) = 0.0         s_pblh(i) = 0.
798         s_lcl(i) = 0.0         s_lcl(i) = 0.
799         s_capCL(i) = 0.0         s_capCL(i) = 0.
800         s_oliqCL(i) = 0.0         s_oliqCL(i) = 0.
801         s_cteiCL(i) = 0.0         s_cteiCL(i) = 0.
802         s_pblT(i) = 0.0         s_pblT(i) = 0.
803         s_therm(i) = 0.0         s_therm(i) = 0.
804         s_trmb1(i) = 0.0         s_trmb1(i) = 0.
805         s_trmb2(i) = 0.0         s_trmb2(i) = 0.
806         s_trmb3(i) = 0.0         s_trmb3(i) = 0.
807    
808         IF (abs(pctsrf(i, is_ter) + pctsrf(i, is_lic) + &         IF (abs(pctsrf(i, is_ter) + pctsrf(i, is_lic) + pctsrf(i, is_oce) &
809              pctsrf(i, is_oce) + pctsrf(i, is_sic) - 1.) .GT. EPSFRA) &              + pctsrf(i, is_sic) - 1.)  >  EPSFRA) print *, &
810              THEN              'physiq : probl\`eme sous surface au point ', i, &
811            WRITE(*, *) 'physiq : pb sous surface au point ', i, &              pctsrf(i, 1 : nbsrf)
                pctsrf(i, 1 : nbsrf)  
        ENDIF  
812      ENDDO      ENDDO
813      DO nsrf = 1, nbsrf      DO nsrf = 1, nbsrf
814         DO i = 1, klon         DO i = 1, klon
# Line 1103  contains Line 836  contains
836         ENDDO         ENDDO
837      ENDDO      ENDDO
838    
839      ! 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 :
   
840      DO nsrf = 1, nbsrf      DO nsrf = 1, nbsrf
841         DO i = 1, klon         DO i = 1, klon
842            IF (pctsrf(i, nsrf) < epsfra) ftsol(i, nsrf) = zxtsol(i)            IF (pctsrf(i, nsrf) < epsfra) ftsol(i, nsrf) = zxtsol(i)
# Line 1116  contains Line 848  contains
848            IF (pctsrf(i, nsrf) < epsfra) ffonte(i, nsrf) = zxffonte(i)            IF (pctsrf(i, nsrf) < epsfra) ffonte(i, nsrf) = zxffonte(i)
849            IF (pctsrf(i, nsrf) < epsfra) &            IF (pctsrf(i, nsrf) < epsfra) &
850                 fqcalving(i, nsrf) = zxfqcalving(i)                 fqcalving(i, nsrf) = zxfqcalving(i)
851            IF (pctsrf(i, nsrf) < epsfra) pblh(i, nsrf)=s_pblh(i)            IF (pctsrf(i, nsrf) < epsfra) pblh(i, nsrf) = s_pblh(i)
852            IF (pctsrf(i, nsrf) < epsfra) plcl(i, nsrf)=s_lcl(i)            IF (pctsrf(i, nsrf) < epsfra) plcl(i, nsrf) = s_lcl(i)
853            IF (pctsrf(i, nsrf) < epsfra) capCL(i, nsrf)=s_capCL(i)            IF (pctsrf(i, nsrf) < epsfra) capCL(i, nsrf) = s_capCL(i)
854            IF (pctsrf(i, nsrf) < epsfra) oliqCL(i, nsrf)=s_oliqCL(i)            IF (pctsrf(i, nsrf) < epsfra) oliqCL(i, nsrf) = s_oliqCL(i)
855            IF (pctsrf(i, nsrf) < epsfra) cteiCL(i, nsrf)=s_cteiCL(i)            IF (pctsrf(i, nsrf) < epsfra) cteiCL(i, nsrf) = s_cteiCL(i)
856            IF (pctsrf(i, nsrf) < epsfra) pblT(i, nsrf)=s_pblT(i)            IF (pctsrf(i, nsrf) < epsfra) pblT(i, nsrf) = s_pblT(i)
857            IF (pctsrf(i, nsrf) < epsfra) therm(i, nsrf)=s_therm(i)            IF (pctsrf(i, nsrf) < epsfra) therm(i, nsrf) = s_therm(i)
858            IF (pctsrf(i, nsrf) < epsfra) trmb1(i, nsrf)=s_trmb1(i)            IF (pctsrf(i, nsrf) < epsfra) trmb1(i, nsrf) = s_trmb1(i)
859            IF (pctsrf(i, nsrf) < epsfra) trmb2(i, nsrf)=s_trmb2(i)            IF (pctsrf(i, nsrf) < epsfra) trmb2(i, nsrf) = s_trmb2(i)
860            IF (pctsrf(i, nsrf) < epsfra) trmb3(i, nsrf)=s_trmb3(i)            IF (pctsrf(i, nsrf) < epsfra) trmb3(i, nsrf) = s_trmb3(i)
861         ENDDO         ENDDO
862      ENDDO      ENDDO
863    
864      ! Calculer la derive du flux infrarouge      ! Calculer la dérive du flux infrarouge
865    
866      DO i = 1, klon      DO i = 1, klon
867         dlw(i) = - 4.0*RSIGMA*zxtsol(i)**3         dlw(i) = - 4. * RSIGMA * zxtsol(i)**3
868      ENDDO      ENDDO
869    
870        IF (check) print *, "avantcon = ", qcheck(paprs, q_seri, ql_seri)
871    
872      ! Appeler la convection (au choix)      ! Appeler la convection (au choix)
873    
874      DO k = 1, llm      if (iflag_con == 2) then
875         DO i = 1, klon         conv_q = d_q_dyn + d_q_vdf / dtphys
876            conv_q(i, k) = d_q_dyn(i, k) &         conv_t = d_t_dyn + d_t_vdf / dtphys
877                 + d_q_vdf(i, k)/dtphys         z_avant = sum((q_seri + ql_seri) * zmasse, dim=2)
878            conv_t(i, k) = d_t_dyn(i, k) &         CALL conflx(dtphys, paprs, play, t_seri(:, llm:1:- 1), &
879                 + d_t_vdf(i, k)/dtphys              q_seri(:, llm:1:- 1), conv_t, conv_q, zxfluxq(:, 1), omega, &
880         ENDDO              d_t_con, d_q_con, rain_con, snow_con, mfu(:, llm:1:- 1), &
881      ENDDO              mfd(:, llm:1:- 1), pen_u, pde_u, pen_d, pde_d, kcbot, kctop, &
882      IF (check) THEN              kdtop, pmflxr, pmflxs)
        za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)  
        print *, "avantcon=", za  
     ENDIF  
     zx_ajustq = .FALSE.  
     IF (iflag_con == 2) zx_ajustq=.TRUE.  
     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  
     IF (iflag_con == 1) THEN  
        stop 'reactiver le call conlmd dans physiq.F'  
     ELSE IF (iflag_con == 2) THEN  
        CALL conflx(dtphys, paprs, play, t_seri, q_seri, &  
             conv_t, conv_q, zxfluxq(1, 1), omega, &  
             d_t_con, d_q_con, rain_con, snow_con, &  
             pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, &  
             kcbot, kctop, kdtop, pmflxr, pmflxs)  
883         WHERE (rain_con < 0.) rain_con = 0.         WHERE (rain_con < 0.) rain_con = 0.
884         WHERE (snow_con < 0.) snow_con = 0.         WHERE (snow_con < 0.) snow_con = 0.
885         DO i = 1, klon         ibas_con = llm + 1 - kcbot
886            ibas_con(i) = llm+1 - kcbot(i)         itop_con = llm + 1 - kctop
887            itop_con(i) = llm+1 - kctop(i)      else
888         ENDDO         ! iflag_con >= 3
     ELSE IF (iflag_con >= 3) THEN  
        ! nb of tracers for the KE convection:  
        ! MAF la partie traceurs est faite dans phytrac  
        ! on met ntra=1 pour limiter les appels mais on peut  
        ! supprimer les calculs / ftra.  
        ntra = 1  
        ! Schema de convection modularise et vectorise:  
        ! (driver commun aux versions 3 et 4)  
   
        IF (ok_cvl) THEN ! new driver for convectL  
           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)  
889    
890            clwcon0=qcondc         da = 0.
891            pmfu=upwd+dnwd         mp = 0.
892           phi = 0.
893           CALL concvl(dtphys, paprs, play, t_seri, q_seri, u_seri, v_seri, sig1, &
894                w01, d_t_con, d_q_con, d_u_con, d_v_con, rain_con, snow_con, &
895                ibas_con, itop_con, upwd, dnwd, dnwd0, Ma, cape, iflagctrl, &
896                qcondc, wd, pmflxr, pmflxs, da, phi, mp)
897           clwcon0 = qcondc
898           mfu = upwd + dnwd
899           IF (.NOT. ok_gust) wd = 0.
900    
901           IF (thermcep) THEN
902              zqsat = MIN(0.5, r2es * FOEEW(t_seri, rtt >= t_seri) / play)
903              zqsat = zqsat / (1. - retv * zqsat)
904         ELSE         ELSE
905            ! MAF conema3 ne contient pas les traceurs            zqsat = merge(qsats(t_seri), qsatl(t_seri), t_seri < t_coup) / play
           CALL conema3 (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, bas, top, &  
                Ma, cape, tvp, rflag, &  
                pbase &  
                , bbase, dtvpdt1, dtvpdq1, dplcldt, dplcldr &  
                , clwcon0)  
        ENDIF ! ok_cvl  
   
        IF (.NOT. ok_gust) THEN  
           do i = 1, klon  
              wd(i)=0.0  
           enddo  
906         ENDIF         ENDIF
907    
908         ! Calcul des proprietes des nuages convectifs         ! Properties of convective clouds
909           clwcon0 = fact_cldcon * clwcon0
910         DO k = 1, llm         call clouds_gno(klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, &
911            DO i = 1, klon              rnebcon0)
912               zx_t = t_seri(i, k)  
913               IF (thermcep) THEN         mfd = 0.
914                  zdelta = MAX(0., SIGN(1., rtt-zx_t))         pen_u = 0.
915                  zx_qs = r2es * FOEEW(zx_t, zdelta)/play(i, k)         pen_d = 0.
916                  zx_qs = MIN(0.5, zx_qs)         pde_d = 0.
917                  zcor = 1./(1.-retv*zx_qs)         pde_u = 0.
918                  zx_qs = zx_qs*zcor      END if
              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  
        call clouds_gno &  
             (klon, llm, q_seri, zqsat, clwcon0, ptconv, ratqsc, rnebcon0)  
     ELSE  
        print *, "iflag_con non-prevu", iflag_con  
        stop 1  
     ENDIF  
919    
920      DO k = 1, llm      DO k = 1, llm
921         DO i = 1, klon         DO i = 1, klon
# Line 1256  contains Line 927  contains
927      ENDDO      ENDDO
928    
929      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
930         ztit='after convect'         tit = 'after convect'
931         CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
932              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)
933              d_ql, d_qs, d_ec)         call diagphy(airephy, tit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
934         call diagphy(airephy, ztit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &              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)  
935      END IF      END IF
936    
937      IF (check) THEN      IF (check) THEN
938         za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)         za = qcheck(paprs, q_seri, ql_seri)
939         print *,"aprescon=", za         print *, "aprescon = ", za
940         zx_t = 0.0         zx_t = 0.
941         za = 0.0         za = 0.
942         DO i = 1, klon         DO i = 1, klon
943            za = za + airephy(i)/REAL(klon)            za = za + airephy(i)/REAL(klon)
944            zx_t = zx_t + (rain_con(i)+ &            zx_t = zx_t + (rain_con(i)+ &
945                 snow_con(i))*airephy(i)/REAL(klon)                 snow_con(i))*airephy(i)/REAL(klon)
946         ENDDO         ENDDO
947         zx_t = zx_t/za*dtphys         zx_t = zx_t/za*dtphys
948         print *,"Precip=", zx_t         print *, "Precip = ", zx_t
949      ENDIF      ENDIF
950      IF (zx_ajustq) THEN  
951         DO i = 1, klon      IF (iflag_con == 2) THEN
952            z_apres(i) = 0.0         z_apres = sum((q_seri + ql_seri) * zmasse, dim=2)
953         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  
954         DO k = 1, llm         DO k = 1, llm
955            DO i = 1, klon            DO i = 1, klon
956               IF (z_factor(i).GT.(1.0+1.0E-08) .OR. &               IF (z_factor(i) > 1. + 1E-8 .OR. z_factor(i) < 1. - 1E-8) THEN
                   z_factor(i) < (1.0-1.0E-08)) THEN  
957                  q_seri(i, k) = q_seri(i, k) * z_factor(i)                  q_seri(i, k) = q_seri(i, k) * z_factor(i)
958               ENDIF               ENDIF
959            ENDDO            ENDDO
960         ENDDO         ENDDO
961      ENDIF      ENDIF
     zx_ajustq=.FALSE.  
962    
963      ! Convection seche (thermiques ou ajustement)      ! Convection s\`eche (thermiques ou ajustement)
964    
965      d_t_ajs=0.      d_t_ajs = 0.
966      d_u_ajs=0.      d_u_ajs = 0.
967      d_v_ajs=0.      d_v_ajs = 0.
968      d_q_ajs=0.      d_q_ajs = 0.
969      fm_therm=0.      fm_therm = 0.
970      entr_therm=0.      entr_therm = 0.
971    
972      if (iflag_thermals == 0) then      if (iflag_thermals == 0) then
973         ! Ajustement sec         ! Ajustement sec
# Line 1324  contains Line 981  contains
981      endif      endif
982    
983      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
984         ztit='after dry_adjust'         tit = 'after dry_adjust'
985         CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
986              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)  
987      END IF      END IF
988    
989      ! Caclul des ratqs      ! Caclul des ratqs
990    
991      ! 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
992      ! on ecrase le tableau ratqsc calcule par clouds_gno      ! on \'ecrase le tableau ratqsc calcul\'e par clouds_gno
993      if (iflag_cldcon == 1) then      if (iflag_cldcon == 1) then
994         do k=1, llm         do k = 1, llm
995            do i=1, klon            do i = 1, klon
996               if(ptconv(i, k)) then               if(ptconv(i, k)) then
997                  ratqsc(i, k)=ratqsbas &                  ratqsc(i, k) = ratqsbas + fact_cldcon &
998                       +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)
999               else               else
1000                  ratqsc(i, k)=0.                  ratqsc(i, k) = 0.
1001               endif               endif
1002            enddo            enddo
1003         enddo         enddo
1004      endif      endif
1005    
1006      ! ratqs stables      ! ratqs stables
1007      do k=1, llm      do k = 1, llm
1008         do i=1, klon         do i = 1, klon
1009            ratqss(i, k)=ratqsbas+(ratqshaut-ratqsbas)* &            ratqss(i, k) = ratqsbas + (ratqshaut - ratqsbas) &
1010                 min((paprs(i, 1)-play(i, k))/(paprs(i, 1)-30000.), 1.)                 * min((paprs(i, 1) - play(i, k)) / (paprs(i, 1) - 3e4), 1.)
1011         enddo         enddo
1012      enddo      enddo
1013    
1014      ! ratqs final      ! ratqs final
1015      if (iflag_cldcon == 1 .or.iflag_cldcon == 2) then      if (iflag_cldcon == 1 .or. iflag_cldcon == 2) then
1016         ! les ratqs sont une conbinaison de ratqss et ratqsc         ! les ratqs sont une conbinaison de ratqss et ratqsc
1017         ! ratqs final         ! ratqs final
1018         ! 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
1019         ! relaxation des ratqs         ! relaxation des ratqs
1020         facteur=exp(-dtphys*facttemps)         ratqs = max(ratqs * exp(- dtphys * facttemps), ratqss)
1021         ratqs=max(ratqs*facteur, ratqss)         ratqs = max(ratqs, ratqsc)
        ratqs=max(ratqs, ratqsc)  
1022      else      else
1023         ! on ne prend que le ratqs stable pour fisrtilp         ! on ne prend que le ratqs stable pour fisrtilp
1024         ratqs=ratqss         ratqs = ratqss
1025      endif      endif
1026    
1027      ! Appeler le processus de condensation a grande echelle      CALL fisrtilp(dtphys, paprs, play, t_seri, q_seri, ptconv, ratqs, &
1028      ! et le processus de precipitation           d_t_lsc, d_q_lsc, d_ql_lsc, rneb, cldliq, rain_lsc, snow_lsc, &
1029      CALL fisrtilp(dtphys, paprs, play, &           pfrac_impa, pfrac_nucl, pfrac_1nucl, frac_impa, frac_nucl, prfl, &
1030           t_seri, q_seri, ptconv, ratqs, &           psfl, rhcl)
          d_t_lsc, d_q_lsc, d_ql_lsc, rneb, cldliq, &  
          rain_lsc, snow_lsc, &  
          pfrac_impa, pfrac_nucl, pfrac_1nucl, &  
          frac_impa, frac_nucl, &  
          prfl, psfl, rhcl)  
1031    
1032      WHERE (rain_lsc < 0) rain_lsc = 0.      WHERE (rain_lsc < 0) rain_lsc = 0.
1033      WHERE (snow_lsc < 0) snow_lsc = 0.      WHERE (snow_lsc < 0) snow_lsc = 0.
# Line 1391  contains Line 1041  contains
1041         ENDDO         ENDDO
1042      ENDDO      ENDDO
1043      IF (check) THEN      IF (check) THEN
1044         za = qcheck(klon, llm, paprs, q_seri, ql_seri, airephy)         za = qcheck(paprs, q_seri, ql_seri)
1045         print *,"apresilp=", za         print *, "apresilp = ", za
1046         zx_t = 0.0         zx_t = 0.
1047         za = 0.0         za = 0.
1048         DO i = 1, klon         DO i = 1, klon
1049            za = za + airephy(i)/REAL(klon)            za = za + airephy(i)/REAL(klon)
1050            zx_t = zx_t + (rain_lsc(i) &            zx_t = zx_t + (rain_lsc(i) &
1051                 + snow_lsc(i))*airephy(i)/REAL(klon)                 + snow_lsc(i))*airephy(i)/REAL(klon)
1052         ENDDO         ENDDO
1053         zx_t = zx_t/za*dtphys         zx_t = zx_t/za*dtphys
1054         print *,"Precip=", zx_t         print *, "Precip = ", zx_t
1055      ENDIF      ENDIF
1056    
1057      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
1058         ztit='after fisrt'         tit = 'after fisrt'
1059         CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1060              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)
1061              d_ql, d_qs, d_ec)         call diagphy(airephy, tit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &
1062         call diagphy(airephy, ztit, ip_ebil, zero_v, zero_v, zero_v, zero_v, &              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)  
1063      END IF      END IF
1064    
1065      ! PRESCRIPTION DES NUAGES POUR LE RAYONNEMENT      ! PRESCRIPTION DES NUAGES POUR LE RAYONNEMENT
1066    
1067      ! 1. NUAGES CONVECTIFS      ! 1. NUAGES CONVECTIFS
1068    
1069      IF (iflag_cldcon.le.-1) THEN ! seulement pour Tiedtke      IF (iflag_cldcon <= - 1) THEN
1070         snow_tiedtke=0.         ! seulement pour Tiedtke
1071         if (iflag_cldcon == -1) then         snow_tiedtke = 0.
1072            rain_tiedtke=rain_con         if (iflag_cldcon == - 1) then
1073              rain_tiedtke = rain_con
1074         else         else
1075            rain_tiedtke=0.            rain_tiedtke = 0.
1076            do k=1, llm            do k = 1, llm
1077               do i=1, klon               do i = 1, klon
1078                  if (d_q_con(i, k) < 0.) then                  if (d_q_con(i, k) < 0.) then
1079                     rain_tiedtke(i)=rain_tiedtke(i)-d_q_con(i, k)/dtphys &                     rain_tiedtke(i) = rain_tiedtke(i) - d_q_con(i, k)/dtphys &
1080                          *zmasse(i, k)                          *zmasse(i, k)
1081                  endif                  endif
1082               enddo               enddo
# Line 1435  contains Line 1084  contains
1084         endif         endif
1085    
1086         ! Nuages diagnostiques pour Tiedtke         ! Nuages diagnostiques pour Tiedtke
1087         CALL diagcld1(paprs, play, &         CALL diagcld1(paprs, play, rain_tiedtke, snow_tiedtke, ibas_con, &
1088              rain_tiedtke, snow_tiedtke, ibas_con, itop_con, &              itop_con, diafra, dialiq)
             diafra, dialiq)  
1089         DO k = 1, llm         DO k = 1, llm
1090            DO i = 1, klon            DO i = 1, klon
1091               IF (diafra(i, k).GT.cldfra(i, k)) THEN               IF (diafra(i, k) > cldfra(i, k)) THEN
1092                  cldliq(i, k) = dialiq(i, k)                  cldliq(i, k) = dialiq(i, k)
1093                  cldfra(i, k) = diafra(i, k)                  cldfra(i, k) = diafra(i, k)
1094               ENDIF               ENDIF
1095            ENDDO            ENDDO
1096         ENDDO         ENDDO
1097      ELSE IF (iflag_cldcon == 3) THEN      ELSE IF (iflag_cldcon == 3) THEN
1098         ! On prend pour les nuages convectifs le max du calcul de la         ! On prend pour les nuages convectifs le maximum du calcul de
1099         ! 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
1100         ! facttemps         ! d'un facteur facttemps.
1101         facteur = dtphys *facttemps         facteur = dtphys * facttemps
1102         do k=1, llm         do k = 1, llm
1103            do i=1, klon            do i = 1, klon
1104               rnebcon(i, k)=rnebcon(i, k)*facteur               rnebcon(i, k) = rnebcon(i, k) * facteur
1105               if (rnebcon0(i, k)*clwcon0(i, k).gt.rnebcon(i, k)*clwcon(i, k)) &               if (rnebcon0(i, k) * clwcon0(i, k) &
1106                    then                    > rnebcon(i, k) * clwcon(i, k)) then
1107                  rnebcon(i, k)=rnebcon0(i, k)                  rnebcon(i, k) = rnebcon0(i, k)
1108                  clwcon(i, k)=clwcon0(i, k)                  clwcon(i, k) = clwcon0(i, k)
1109               endif               endif
1110            enddo            enddo
1111         enddo         enddo
1112    
1113         ! On prend la somme des fractions nuageuses et des contenus en eau         ! On prend la somme des fractions nuageuses et des contenus en eau
1114         cldfra=min(max(cldfra, rnebcon), 1.)         cldfra = min(max(cldfra, rnebcon), 1.)
1115         cldliq=cldliq+rnebcon*clwcon         cldliq = cldliq + rnebcon*clwcon
1116      ENDIF      ENDIF
1117    
1118      ! 2. NUAGES STARTIFORMES      ! 2. Nuages stratiformes
1119    
1120      IF (ok_stratus) THEN      IF (ok_stratus) THEN
1121         CALL diagcld2(paprs, play, t_seri, q_seri, diafra, dialiq)         CALL diagcld2(paprs, play, t_seri, q_seri, diafra, dialiq)
1122         DO k = 1, llm         DO k = 1, llm
1123            DO i = 1, klon            DO i = 1, klon
1124               IF (diafra(i, k).GT.cldfra(i, k)) THEN               IF (diafra(i, k) > cldfra(i, k)) THEN
1125                  cldliq(i, k) = dialiq(i, k)                  cldliq(i, k) = dialiq(i, k)
1126                  cldfra(i, k) = diafra(i, k)                  cldfra(i, k) = diafra(i, k)
1127               ENDIF               ENDIF
# Line 1482  contains Line 1130  contains
1130      ENDIF      ENDIF
1131    
1132      ! Precipitation totale      ! Precipitation totale
   
1133      DO i = 1, klon      DO i = 1, klon
1134         rain_fall(i) = rain_con(i) + rain_lsc(i)         rain_fall(i) = rain_con(i) + rain_lsc(i)
1135         snow_fall(i) = snow_con(i) + snow_lsc(i)         snow_fall(i) = snow_con(i) + snow_lsc(i)
1136      ENDDO      ENDDO
1137    
1138      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) CALL diagetpq(airephy, "after diagcld", ip_ebil, 2, 2, &
1139         ztit="after diagcld"           dtphys, t_seri, q_seri, ql_seri, u_seri, v_seri, paprs, d_h_vcol, &
1140         CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &           d_qt, d_ec)
             ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &  
             d_ql, d_qs, d_ec)  
     END IF  
   
     ! Calculer l'humidite relative pour diagnostique  
1141    
1142        ! Humidit\'e relative pour diagnostic :
1143      DO k = 1, llm      DO k = 1, llm
1144         DO i = 1, klon         DO i = 1, klon
1145            zx_t = t_seri(i, k)            zx_t = t_seri(i, k)
1146            IF (thermcep) THEN            IF (thermcep) THEN
1147               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)  
1148               zx_qs = MIN(0.5, zx_qs)               zx_qs = MIN(0.5, zx_qs)
1149               zcor = 1./(1.-retv*zx_qs)               zcor = 1./(1. - retv*zx_qs)
1150               zx_qs = zx_qs*zcor               zx_qs = zx_qs*zcor
1151            ELSE            ELSE
1152               IF (zx_t < t_coup) THEN               IF (zx_t < t_coup) THEN
# Line 1514  contains Line 1156  contains
1156               ENDIF               ENDIF
1157            ENDIF            ENDIF
1158            zx_rh(i, k) = q_seri(i, k)/zx_qs            zx_rh(i, k) = q_seri(i, k)/zx_qs
1159            zqsat(i, k)=zx_qs            zqsat(i, k) = zx_qs
1160         ENDDO         ENDDO
1161      ENDDO      ENDDO
1162      !jq - introduce the aerosol direct and first indirect radiative forcings  
1163      !jq - Johannes Quaas, 27/11/2003 (quaas@lmd.jussieu.fr)      ! Introduce the aerosol direct and first indirect radiative forcings:
1164      IF (ok_ade.OR.ok_aie) THEN      IF (ok_ade .OR. ok_aie) THEN
1165         ! Get sulfate aerosol distribution         ! Get sulfate aerosol distribution :
1166         CALL readsulfate(rdayvrai, firstcal, sulfate)         CALL readsulfate(dayvrai, time, firstcal, sulfate)
1167         CALL readsulfate_preind(rdayvrai, firstcal, sulfate_pi)         CALL readsulfate_preind(dayvrai, time, firstcal, sulfate_pi)
1168    
1169         ! Calculate aerosol optical properties (Olivier Boucher)         CALL aeropt(play, paprs, t_seri, sulfate, rhcl, tau_ae, piz_ae, cg_ae, &
1170         CALL aeropt(play, paprs, t_seri, sulfate, rhcl, &              aerindex)
             tau_ae, piz_ae, cg_ae, aerindex)  
1171      ELSE      ELSE
1172         tau_ae=0.0         tau_ae = 0.
1173         piz_ae=0.0         piz_ae = 0.
1174         cg_ae=0.0         cg_ae = 0.
1175      ENDIF      ENDIF
1176    
1177      ! Calculer les parametres optiques des nuages et quelques      ! Param\`etres optiques des nuages et quelques param\`etres pour
1178      ! parametres pour diagnostiques:      ! diagnostics :
   
1179      if (ok_newmicro) then      if (ok_newmicro) then
1180         CALL newmicro (paprs, play, ok_newmicro, &         CALL newmicro(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, &
1181              t_seri, cldliq, cldfra, cldtau, cldemi, &              cldh, cldl, cldm, cldt, cldq, flwp, fiwp, flwc, fiwc, ok_aie, &
1182              cldh, cldl, cldm, cldt, cldq, &              sulfate, sulfate_pi, bl95_b0, bl95_b1, cldtaupi, re, fl)
             flwp, fiwp, flwc, fiwc, &  
             ok_aie, &  
             sulfate, sulfate_pi, &  
             bl95_b0, bl95_b1, &  
             cldtaupi, re, fl)  
1183      else      else
1184         CALL nuage (paprs, play, &         CALL nuage(paprs, play, t_seri, cldliq, cldfra, cldtau, cldemi, cldh, &
1185              t_seri, cldliq, cldfra, cldtau, cldemi, &              cldl, cldm, cldt, cldq, ok_aie, sulfate, sulfate_pi, bl95_b0, &
1186              cldh, cldl, cldm, cldt, cldq, &              bl95_b1, cldtaupi, re, fl)
             ok_aie, &  
             sulfate, sulfate_pi, &  
             bl95_b0, bl95_b1, &  
             cldtaupi, re, fl)  
   
1187      endif      endif
1188    
1189      ! Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.      IF (MOD(itap - 1, radpas) == 0) THEN
1190           ! Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.
1191      IF (MOD(itaprad, radpas) == 0) THEN         ! Calcul de l'abedo moyen par maille
1192         DO i = 1, klon         albsol = sum(falbe * pctsrf, dim = 2)
1193            albsol(i) = falbe(i, is_oce) * pctsrf(i, is_oce) &  
1194                 + falbe(i, is_lic) * pctsrf(i, is_lic) &         ! Rayonnement (compatible Arpege-IFS) :
1195                 + falbe(i, is_ter) * pctsrf(i, is_ter) &         CALL radlwsw(dist, mu0, fract, paprs, play, zxtsol, albsol, t_seri, &
1196                 + falbe(i, is_sic) * pctsrf(i, is_sic)              q_seri, wo, cldfra, cldemi, cldtau, heat, heat0, cool, cool0, &
1197            albsollw(i) = falblw(i, is_oce) * pctsrf(i, is_oce) &              radsol, albpla, topsw, toplw, solsw, sollw, sollwdown, topsw0, &
1198                 + falblw(i, is_lic) * pctsrf(i, is_lic) &              toplw0, solsw0, sollw0, lwdn0, lwdn, lwup0, lwup, swdn0, swdn, &
1199                 + falblw(i, is_ter) * pctsrf(i, is_ter) &              swup0, swup, ok_ade, ok_aie, tau_ae, piz_ae, cg_ae, topswad, &
1200                 + falblw(i, is_sic) * pctsrf(i, is_sic)              solswad, cldtaupi, topswai, solswai)
        ENDDO  
        ! nouveau rayonnement (compatible Arpege-IFS):  
        CALL radlwsw(dist, rmu0, fract, paprs, play, zxtsol, albsol, &  
             albsollw, t_seri, q_seri, wo, cldfra, cldemi, cldtau, heat, &  
             heat0, cool, cool0, radsol, albpla, topsw, toplw, solsw, sollw, &  
             sollwdown, topsw0, toplw0, solsw0, sollw0, lwdn0, lwdn, lwup0, &  
             lwup, swdn0, swdn, swup0, swup, ok_ade, ok_aie, tau_ae, piz_ae, &  
             cg_ae, topswad, solswad, cldtaupi, topswai, solswai)  
        itaprad = 0  
1201      ENDIF      ENDIF
     itaprad = itaprad + 1  
1202    
1203      ! Ajouter la tendance des rayonnements (tous les pas)      ! Ajouter la tendance des rayonnements (tous les pas)
1204    
1205      DO k = 1, llm      DO k = 1, llm
1206         DO i = 1, klon         DO i = 1, klon
1207            t_seri(i, k) = t_seri(i, k) &            t_seri(i, k) = t_seri(i, k) + (heat(i, k) - cool(i, k)) * dtphys/86400.
                + (heat(i, k)-cool(i, k)) * dtphys/86400.  
1208         ENDDO         ENDDO
1209      ENDDO      ENDDO
1210    
1211      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) THEN
1212         ztit='after rad'         tit = 'after rad'
1213         CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &
1214              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)
1215              d_ql, d_qs, d_ec)         call diagphy(airephy, tit, ip_ebil, topsw, toplw, solsw, sollw, &
1216         call diagphy(airephy, ztit, ip_ebil, topsw, toplw, solsw, sollw, &              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)  
1217      END IF      END IF
1218    
1219      ! Calculer l'hydrologie de la surface      ! Calculer l'hydrologie de la surface
1220      DO i = 1, klon      DO i = 1, klon
1221         zxqsurf(i) = 0.0         zxqsurf(i) = 0.
1222         zxsnow(i) = 0.0         zxsnow(i) = 0.
1223      ENDDO      ENDDO
1224      DO nsrf = 1, nbsrf      DO nsrf = 1, nbsrf
1225         DO i = 1, klon         DO i = 1, klon
# Line 1611  contains Line 1228  contains
1228         ENDDO         ENDDO
1229      ENDDO      ENDDO
1230    
1231      ! Calculer le bilan du sol et la derive de temperature (couplage)      ! Calculer le bilan du sol et la d\'erive de temp\'erature (couplage)
1232    
1233      DO i = 1, klon      DO i = 1, klon
1234         bils(i) = radsol(i) - sens(i) + zxfluxlat(i)         bils(i) = radsol(i) - sens(i) + zxfluxlat(i)
1235      ENDDO      ENDDO
1236    
1237      !mod deb lott(jan95)      ! Param\'etrisation de l'orographie \`a l'\'echelle sous-maille :
     ! Appeler le programme de parametrisation de l'orographie  
     ! a l'echelle sous-maille:  
1238    
1239      IF (ok_orodr) THEN      IF (ok_orodr) THEN
1240         ! selection des points pour lesquels le shema est actif:         ! S\'election des points pour lesquels le sch\'ema est actif :
1241         igwd=0         igwd = 0
1242         DO i=1, klon         DO i = 1, klon
1243            itest(i)=0            itest(i) = 0
1244            IF (((zpic(i)-zmea(i)).GT.100.).AND.(zstd(i).GT.10.0)) THEN            IF (zpic(i) - zmea(i) > 100. .AND. zstd(i) > 10.) THEN
1245               itest(i)=1               itest(i) = 1
1246               igwd=igwd+1               igwd = igwd + 1
1247               idx(igwd)=i               idx(igwd) = i
1248            ENDIF            ENDIF
1249         ENDDO         ENDDO
1250    
1251         CALL drag_noro(klon, llm, dtphys, paprs, play, &         CALL drag_noro(klon, llm, dtphys, paprs, play, zmea, zstd, zsig, zgam, &
1252              zmea, zstd, zsig, zgam, zthe, zpic, zval, &              zthe, zpic, zval, itest, t_seri, u_seri, v_seri, zulow, zvlow, &
1253              igwd, idx, itest, &              zustrdr, zvstrdr, d_t_oro, d_u_oro, d_v_oro)
             t_seri, u_seri, v_seri, &  
             zulow, zvlow, zustrdr, zvstrdr, &  
             d_t_oro, d_u_oro, d_v_oro)  
1254    
1255         ! ajout des tendances         ! ajout des tendances
1256         DO k = 1, llm         DO k = 1, llm
# Line 1651  contains Line 1263  contains
1263      ENDIF      ENDIF
1264    
1265      IF (ok_orolf) THEN      IF (ok_orolf) THEN
1266         ! selection des points pour lesquels le shema est actif:         ! S\'election des points pour lesquels le sch\'ema est actif :
1267         igwd=0         igwd = 0
1268         DO i=1, klon         DO i = 1, klon
1269            itest(i)=0            itest(i) = 0
1270            IF ((zpic(i)-zmea(i)).GT.100.) THEN            IF (zpic(i) - zmea(i) > 100.) THEN
1271               itest(i)=1               itest(i) = 1
1272               igwd=igwd+1               igwd = igwd + 1
1273               idx(igwd)=i               idx(igwd) = i
1274            ENDIF            ENDIF
1275         ENDDO         ENDDO
1276    
# Line 1666  contains Line 1278  contains
1278              itest, t_seri, u_seri, v_seri, zulow, zvlow, zustrli, zvstrli, &              itest, t_seri, u_seri, v_seri, zulow, zvlow, zustrli, zvstrli, &
1279              d_t_lif, d_u_lif, d_v_lif)              d_t_lif, d_u_lif, d_v_lif)
1280    
1281         ! ajout des tendances         ! Ajout des tendances :
1282         DO k = 1, llm         DO k = 1, llm
1283            DO i = 1, klon            DO i = 1, klon
1284               t_seri(i, k) = t_seri(i, k) + d_t_lif(i, k)               t_seri(i, k) = t_seri(i, k) + d_t_lif(i, k)
# Line 1676  contains Line 1288  contains
1288         ENDDO         ENDDO
1289      ENDIF      ENDIF
1290    
1291      ! STRESS NECESSAIRES: TOUTE LA PHYSIQUE      ! Stress n\'ecessaires : toute la physique
1292    
1293      DO i = 1, klon      DO i = 1, klon
1294         zustrph(i)=0.         zustrph(i) = 0.
1295         zvstrph(i)=0.         zvstrph(i) = 0.
1296      ENDDO      ENDDO
1297      DO k = 1, llm      DO k = 1, llm
1298         DO i = 1, klon         DO i = 1, klon
1299            zustrph(i)=zustrph(i)+(u_seri(i, k)-u(i, k))/dtphys* zmasse(i, k)            zustrph(i) = zustrph(i) + (u_seri(i, k) - u(i, k)) / dtphys &
1300            zvstrph(i)=zvstrph(i)+(v_seri(i, k)-v(i, k))/dtphys* zmasse(i, k)                 * zmasse(i, k)
1301              zvstrph(i) = zvstrph(i) + (v_seri(i, k) - v(i, k)) / dtphys &
1302                   * zmasse(i, k)
1303         ENDDO         ENDDO
1304      ENDDO      ENDDO
1305    
1306      !IM calcul composantes axiales du moment angulaire et couple des montagnes      CALL aaam_bud(rg, romega, rlat, rlon, pphis, zustrdr, zustrli, zustrph, &
1307             zvstrdr, zvstrli, zvstrph, paprs, u, v, aam, torsfc)
     CALL aaam_bud(27, klon, llm, time, ra, rg, romega, rlat, rlon, pphis, &  
          zustrdr, zustrli, zustrph, zvstrdr, zvstrli, zvstrph, paprs, u, v, &  
          aam, torsfc)  
1308    
1309      IF (if_ebil >= 2) THEN      IF (if_ebil >= 2) CALL diagetpq(airephy, 'after orography', ip_ebil, 2, &
1310         ztit='after orography'           2, dtphys, t_seri, q_seri, ql_seri, u_seri, v_seri, paprs, d_h_vcol, &
1311         CALL diagetpq(airephy, ztit, ip_ebil, 2, 2, dtphys, t_seri, q_seri, &           d_qt, d_ec)
             ql_seri, qs_seri, u_seri, v_seri, paprs, d_h_vcol, d_qt, d_qw, &  
             d_ql, d_qs, d_ec)  
     END IF  
1312    
1313      ! Calcul des tendances traceurs      ! Calcul des tendances traceurs
1314      call phytrac(rnpb, itap, lmt_pas, julien, time, firstcal, lafin, &      call phytrac(itap, lmt_pas, julien, time, firstcal, lafin, dtphys, t, &
1315           nqmx-2, dtphys, u, t, paprs, play, pmfu, pmfd, pen_u, pde_u, &           paprs, play, mfu, mfd, pde_u, pen_d, ycoefh, fm_therm, entr_therm, &
1316           pen_d, pde_d, ycoefh, fm_therm, entr_therm, yu1, yv1, ftsol, pctsrf, &           yu1, yv1, ftsol, pctsrf, frac_impa, frac_nucl, da, phi, mp, upwd, &
1317           frac_impa, frac_nucl, pphis, albsol, rhcl, cldfra, rneb, &           dnwd, tr_seri, zmasse, ncid_startphy, nid_ins)
1318           diafra, cldliq, pmflxr, pmflxs, prfl, psfl, da, phi, mp, upwd, dnwd, &  
1319           tr_seri, zmasse)      IF (offline) call phystokenc(dtphys, rlon, rlat, t, mfu, mfd, pen_u, &
1320             pde_u, pen_d, pde_d, fm_therm, entr_therm, ycoefh, yu1, yv1, ftsol, &
1321      IF (offline) THEN           pctsrf, frac_impa, frac_nucl, pphis, airephy, dtphys, itap)
        call phystokenc(dtphys, rlon, rlat, t, pmfu, pmfd, pen_u, pde_u, &  
             pen_d, pde_d, fm_therm, entr_therm, ycoefh, yu1, yv1, ftsol, &  
             pctsrf, frac_impa, frac_nucl, pphis, airephy, dtphys, itap)  
     ENDIF  
1322    
1323      ! Calculer le transport de l'eau et de l'energie (diagnostique)      ! Calculer le transport de l'eau et de l'energie (diagnostique)
1324      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)  
1325    
1326      ! diag. bilKP      ! diag. bilKP
1327    
1328      CALL transp_lay (paprs, zxtsol, t_seri, q_seri, u_seri, v_seri, zphi, &      CALL transp_lay(paprs, zxtsol, t_seri, q_seri, u_seri, v_seri, zphi, &
1329           ve_lay, vq_lay, ue_lay, uq_lay)           ve_lay, vq_lay, ue_lay, uq_lay)
1330    
1331      ! Accumuler les variables a stocker dans les fichiers histoire:      ! Accumuler les variables a stocker dans les fichiers histoire:
1332    
1333      !+jld ec_conser      ! conversion Ec -> E thermique
1334      DO k = 1, llm      DO k = 1, llm
1335         DO i = 1, klon         DO i = 1, klon
1336            ZRCPD = RCPD*(1.0+RVTMP2*q_seri(i, k))            ZRCPD = RCPD * (1. + RVTMP2 * q_seri(i, k))
1337            d_t_ec(i, k)=0.5/ZRCPD &            d_t_ec(i, k) = 0.5 / ZRCPD &
1338                 *(u(i, k)**2+v(i, k)**2-u_seri(i, k)**2-v_seri(i, k)**2)                 * (u(i, k)**2 + v(i, k)**2 - u_seri(i, k)**2 - v_seri(i, k)**2)
1339            t_seri(i, k)=t_seri(i, k)+d_t_ec(i, k)            t_seri(i, k) = t_seri(i, k) + d_t_ec(i, k)
1340            d_t_ec(i, k) = d_t_ec(i, k)/dtphys            d_t_ec(i, k) = d_t_ec(i, k) / dtphys
1341         END DO         END DO
1342      END DO      END DO
1343      !-jld ec_conser  
1344      IF (if_ebil >= 1) THEN      IF (if_ebil >= 1) THEN
1345         ztit='after physic'         tit = 'after physic'
1346         CALL diagetpq(airephy, ztit, ip_ebil, 1, 1, dtphys, t_seri, q_seri, &         CALL diagetpq(airephy, tit, ip_ebil, 1, 1, dtphys, t_seri, q_seri, &
1347              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)  
1348         ! Comme les tendances de la physique sont ajoute dans la dynamique,         ! Comme les tendances de la physique sont ajoute dans la dynamique,
1349         ! on devrait avoir que la variation d'entalpie par la dynamique         ! on devrait avoir que la variation d'entalpie par la dynamique
1350         ! 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.
1351         ! Donc la somme de ces 2 variations devrait etre nulle.         ! Donc la somme de ces 2 variations devrait etre nulle.
1352         call diagphy(airephy, ztit, ip_ebil, topsw, toplw, solsw, sollw, sens, &         call diagphy(airephy, tit, ip_ebil, topsw, toplw, solsw, sollw, sens, &
1353              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)
1354              fs_bound, fq_bound)         d_h_vcol_phy = d_h_vcol
   
        d_h_vcol_phy=d_h_vcol  
   
1355      END IF      END IF
1356    
1357      ! SORTIES      ! SORTIES
1358    
1359      !cc prw = eau precipitable      ! prw = eau precipitable
1360      DO i = 1, klon      DO i = 1, klon
1361         prw(i) = 0.         prw(i) = 0.
1362         DO k = 1, llm         DO k = 1, llm
# Line 1777  contains Line 1376  contains
1376         ENDDO         ENDDO
1377      ENDDO      ENDDO
1378    
1379      IF (nqmx >= 3) THEN      DO iq = 3, nqmx
1380         DO iq = 3, nqmx         DO k = 1, llm
1381            DO k = 1, llm            DO i = 1, klon
1382               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  
1383            ENDDO            ENDDO
1384         ENDDO         ENDDO
1385      ENDIF      ENDDO
1386    
1387      ! 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:
1388      DO k = 1, llm      DO k = 1, llm
# Line 1795  contains Line 1392  contains
1392         ENDDO         ENDDO
1393      ENDDO      ENDDO
1394    
     ! Ecriture des sorties  
     call write_histhf  
     call write_histday  
1395      call write_histins      call write_histins
1396    
1397      ! Si c'est la fin, il faut conserver l'etat de redemarrage      IF (lafin) then
1398      IF (lafin) THEN         call NF95_CLOSE(ncid_startphy)
1399         itau_phy = itau_phy + itap         CALL phyredem(pctsrf, ftsol, ftsoil, tslab, seaice, fqsurf, qsol, &
1400         CALL phyredem("restartphy.nc", rlat, rlon, pctsrf, ftsol, ftsoil, &              fsnow, falbe, fevap, rain_fall, snow_fall, solsw, sollw, dlw, &
1401              tslab, seaice, fqsurf, qsol, fsnow, falbe, falblw, fevap, &              radsol, frugs, agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, &
1402              rain_fall, snow_fall, solsw, sollwdown, dlw, radsol, frugs, &              t_ancien, q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0, sig1, &
1403              agesno, zmea, zstd, zsig, zgam, zthe, zpic, zval, t_ancien, &              w01)
1404              q_ancien, rnebcon, ratqs, clwcon, run_off_lic_0)      end IF
     ENDIF  
1405    
1406      firstcal = .FALSE.      firstcal = .FALSE.
1407    
1408    contains    contains
1409    
     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  
   
     !***************************************************************  
   
1410      subroutine write_histins      subroutine write_histins
1411    
1412        ! 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
1413    
1414        real zout        ! Ecriture des sorties
1415        integer itau_w ! pas de temps ecriture  
1416          use dimens_m, only: iim, jjm
1417          USE histsync_m, ONLY: histsync
1418          USE histwrite_m, ONLY: histwrite
1419    
1420          integer i, itau_w ! pas de temps ecriture
1421          REAL zx_tmp_2d(iim, jjm + 1), zx_tmp_3d(iim, jjm + 1, llm)
1422    
1423        !--------------------------------------------------        !--------------------------------------------------
1424    
1425        IF (ok_instan) THEN        IF (ok_instan) THEN
1426           ! Champs 2D:           ! Champs 2D:
1427    
          zsto = dtphys * ecrit_ins  
          zout = dtphys * ecrit_ins  
1428           itau_w = itau_phy + itap           itau_w = itau_phy + itap
1429    
1430           i = NINT(zout/zsto)           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)  
1431           CALL histwrite(nid_ins, "phis", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "phis", itau_w, zx_tmp_2d)
1432    
1433           i = NINT(zout/zsto)           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)  
1434           CALL histwrite(nid_ins, "aire", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "aire", itau_w, zx_tmp_2d)
1435    
1436           DO i = 1, klon           DO i = 1, klon
1437              zx_tmp_fi2d(i) = paprs(i, 1)              zx_tmp_fi2d(i) = paprs(i, 1)
1438           ENDDO           ENDDO
1439           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)
1440           CALL histwrite(nid_ins, "psol", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "psol", itau_w, zx_tmp_2d)
1441    
1442           DO i = 1, klon           DO i = 1, klon
1443              zx_tmp_fi2d(i) = rain_fall(i) + snow_fall(i)              zx_tmp_fi2d(i) = rain_fall(i) + snow_fall(i)
1444           ENDDO           ENDDO
1445           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)
1446           CALL histwrite(nid_ins, "precip", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "precip", itau_w, zx_tmp_2d)
1447    
1448           DO i = 1, klon           DO i = 1, klon
1449              zx_tmp_fi2d(i) = rain_lsc(i) + snow_lsc(i)              zx_tmp_fi2d(i) = rain_lsc(i) + snow_lsc(i)
1450           ENDDO           ENDDO
1451           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)
1452           CALL histwrite(nid_ins, "plul", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "plul", itau_w, zx_tmp_2d)
1453    
1454           DO i = 1, klon           DO i = 1, klon
1455              zx_tmp_fi2d(i) = rain_con(i) + snow_con(i)              zx_tmp_fi2d(i) = rain_con(i) + snow_con(i)
1456           ENDDO           ENDDO
1457           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)
1458           CALL histwrite(nid_ins, "pluc", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "pluc", itau_w, zx_tmp_2d)
1459    
1460           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zxtsol, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zxtsol, zx_tmp_2d)
1461           CALL histwrite(nid_ins, "tsol", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "tsol", itau_w, zx_tmp_2d)
1462           !ccIM           !ccIM
1463           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zt2m, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zt2m, zx_tmp_2d)
1464           CALL histwrite(nid_ins, "t2m", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "t2m", itau_w, zx_tmp_2d)
1465    
1466           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zq2m, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zq2m, zx_tmp_2d)
1467           CALL histwrite(nid_ins, "q2m", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "q2m", itau_w, zx_tmp_2d)
1468    
1469           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zu10m, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zu10m, zx_tmp_2d)
1470           CALL histwrite(nid_ins, "u10m", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "u10m", itau_w, zx_tmp_2d)
1471    
1472           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), zv10m, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, zv10m, zx_tmp_2d)
1473           CALL histwrite(nid_ins, "v10m", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "v10m", itau_w, zx_tmp_2d)
1474    
1475           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), snow_fall, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, snow_fall, zx_tmp_2d)
1476           CALL histwrite(nid_ins, "snow", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "snow", itau_w, zx_tmp_2d)
1477    
1478           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), cdragm, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, cdragm, zx_tmp_2d)
1479           CALL histwrite(nid_ins, "cdrm", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "cdrm", itau_w, zx_tmp_2d)
1480    
1481           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), cdragh, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, cdragh, zx_tmp_2d)
1482           CALL histwrite(nid_ins, "cdrh", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "cdrh", itau_w, zx_tmp_2d)
1483    
1484           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), toplw, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, toplw, zx_tmp_2d)
1485           CALL histwrite(nid_ins, "topl", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "topl", itau_w, zx_tmp_2d)
1486    
1487           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), evap, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, evap, zx_tmp_2d)
1488           CALL histwrite(nid_ins, "evap", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "evap", itau_w, zx_tmp_2d)
1489    
1490           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), solsw, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, solsw, zx_tmp_2d)
1491           CALL histwrite(nid_ins, "sols", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "sols", itau_w, zx_tmp_2d)
1492    
1493           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), sollw, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sollw, zx_tmp_2d)
1494           CALL histwrite(nid_ins, "soll", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "soll", itau_w, zx_tmp_2d)
1495    
1496           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), sollwdown, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, sollwdown, zx_tmp_2d)
1497           CALL histwrite(nid_ins, "solldown", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "solldown", itau_w, zx_tmp_2d)
1498    
1499           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)
1500           CALL histwrite(nid_ins, "bils", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "bils", itau_w, zx_tmp_2d)
1501    
1502           zx_tmp_fi2d(1:klon)=-1*sens(1:klon)           zx_tmp_fi2d(1:klon) = - sens(1:klon)
1503           ! 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)
1504           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)
1505           CALL histwrite(nid_ins, "sens", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "sens", itau_w, zx_tmp_2d)
1506    
1507           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), fder, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, fder, zx_tmp_2d)
1508           CALL histwrite(nid_ins, "fder", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "fder", itau_w, zx_tmp_2d)
1509    
1510           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), d_ts(1, is_oce), zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_oce), zx_tmp_2d)
1511           CALL histwrite(nid_ins, "dtsvdfo", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "dtsvdfo", itau_w, zx_tmp_2d)
1512    
1513           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), d_ts(1, is_ter), zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_ter), zx_tmp_2d)
1514           CALL histwrite(nid_ins, "dtsvdft", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "dtsvdft", itau_w, zx_tmp_2d)
1515    
1516           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), d_ts(1, is_lic), zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_lic), zx_tmp_2d)
1517           CALL histwrite(nid_ins, "dtsvdfg", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "dtsvdfg", itau_w, zx_tmp_2d)
1518    
1519           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), d_ts(1, is_sic), zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, d_ts(1, is_sic), zx_tmp_2d)
1520           CALL histwrite(nid_ins, "dtsvdfi", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "dtsvdfi", itau_w, zx_tmp_2d)
1521    
1522           DO nsrf = 1, nbsrf           DO nsrf = 1, nbsrf
1523              !XXX              !XXX
1524              zx_tmp_fi2d(1 : klon) = pctsrf(1 : klon, nsrf)*100.              zx_tmp_fi2d(1 : klon) = pctsrf(1 : klon, nsrf)*100.
1525              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)
1526              CALL histwrite(nid_ins, "pourc_"//clnsurf(nsrf), itau_w, &              CALL histwrite(nid_ins, "pourc_"//clnsurf(nsrf), itau_w, &
1527                   zx_tmp_2d)                   zx_tmp_2d)
1528    
1529              zx_tmp_fi2d(1 : klon) = pctsrf(1 : klon, nsrf)              zx_tmp_fi2d(1 : klon) = pctsrf(1 : klon, nsrf)
1530              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)
1531              CALL histwrite(nid_ins, "fract_"//clnsurf(nsrf), itau_w, &              CALL histwrite(nid_ins, "fract_"//clnsurf(nsrf), itau_w, &
1532                   zx_tmp_2d)                   zx_tmp_2d)
1533    
1534              zx_tmp_fi2d(1 : klon) = fluxt(1 : klon, 1, nsrf)              zx_tmp_fi2d(1 : klon) = fluxt(1 : klon, 1, nsrf)
1535              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)
1536              CALL histwrite(nid_ins, "sens_"//clnsurf(nsrf), itau_w, &              CALL histwrite(nid_ins, "sens_"//clnsurf(nsrf), itau_w, &
1537                   zx_tmp_2d)                   zx_tmp_2d)
1538    
1539              zx_tmp_fi2d(1 : klon) = fluxlat(1 : klon, nsrf)              zx_tmp_fi2d(1 : klon) = fluxlat(1 : klon, nsrf)
1540              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)
1541              CALL histwrite(nid_ins, "lat_"//clnsurf(nsrf), itau_w, &              CALL histwrite(nid_ins, "lat_"//clnsurf(nsrf), itau_w, &
1542                   zx_tmp_2d)                   zx_tmp_2d)
1543    
1544              zx_tmp_fi2d(1 : klon) = ftsol(1 : klon, nsrf)              zx_tmp_fi2d(1 : klon) = ftsol(1 : klon, nsrf)
1545              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)
1546              CALL histwrite(nid_ins, "tsol_"//clnsurf(nsrf), itau_w, &              CALL histwrite(nid_ins, "tsol_"//clnsurf(nsrf), itau_w, &
1547                   zx_tmp_2d)                   zx_tmp_2d)
1548    
1549              zx_tmp_fi2d(1 : klon) = fluxu(1 : klon, 1, nsrf)              zx_tmp_fi2d(1 : klon) = fluxu(1 : klon, 1, 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, "taux_"//clnsurf(nsrf), itau_w, &              CALL histwrite(nid_ins, "taux_"//clnsurf(nsrf), itau_w, &
1552                   zx_tmp_2d)                   zx_tmp_2d)
1553    
1554              zx_tmp_fi2d(1 : klon) = fluxv(1 : klon, 1, nsrf)              zx_tmp_fi2d(1 : klon) = fluxv(1 : klon, 1, nsrf)
1555              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)
1556              CALL histwrite(nid_ins, "tauy_"//clnsurf(nsrf), itau_w, &              CALL histwrite(nid_ins, "tauy_"//clnsurf(nsrf), itau_w, &
1557                   zx_tmp_2d)                   zx_tmp_2d)
1558    
1559              zx_tmp_fi2d(1 : klon) = frugs(1 : klon, nsrf)              zx_tmp_fi2d(1 : klon) = frugs(1 : klon, nsrf)
1560              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)
1561              CALL histwrite(nid_ins, "rugs_"//clnsurf(nsrf), itau_w, &              CALL histwrite(nid_ins, "rugs_"//clnsurf(nsrf), itau_w, &
1562                   zx_tmp_2d)                   zx_tmp_2d)
1563    
1564              zx_tmp_fi2d(1 : klon) = falbe(1 : klon, nsrf)              zx_tmp_fi2d(1 : klon) = falbe(:, nsrf)
1565              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)
1566              CALL histwrite(nid_ins, "albe_"//clnsurf(nsrf), itau_w, &              CALL histwrite(nid_ins, "albe_"//clnsurf(nsrf), itau_w, &
1567                   zx_tmp_2d)                   zx_tmp_2d)
1568    
1569           END DO           END DO
1570           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)
1571           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)  
1572    
1573           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)
1574           CALL histwrite(nid_ins, "rugs", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "rugs", itau_w, zx_tmp_2d)
1575    
          !IM cf. AM 081204 BEG  
   
1576           !HBTM2           !HBTM2
1577    
1578           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_pblh, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_pblh, zx_tmp_2d)
1579           CALL histwrite(nid_ins, "s_pblh", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "s_pblh", itau_w, zx_tmp_2d)
1580    
1581           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_pblt, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_pblt, zx_tmp_2d)
1582           CALL histwrite(nid_ins, "s_pblt", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "s_pblt", itau_w, zx_tmp_2d)
1583    
1584           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_lcl, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_lcl, zx_tmp_2d)
1585           CALL histwrite(nid_ins, "s_lcl", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "s_lcl", itau_w, zx_tmp_2d)
1586    
1587           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_capCL, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_capCL, zx_tmp_2d)
1588           CALL histwrite(nid_ins, "s_capCL", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "s_capCL", itau_w, zx_tmp_2d)
1589    
1590           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_oliqCL, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_oliqCL, zx_tmp_2d)
1591           CALL histwrite(nid_ins, "s_oliqCL", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "s_oliqCL", itau_w, zx_tmp_2d)
1592    
1593           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_cteiCL, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_cteiCL, zx_tmp_2d)
1594           CALL histwrite(nid_ins, "s_cteiCL", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "s_cteiCL", itau_w, zx_tmp_2d)
1595    
1596           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_therm, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_therm, zx_tmp_2d)
1597           CALL histwrite(nid_ins, "s_therm", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "s_therm", itau_w, zx_tmp_2d)
1598    
1599           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_trmb1, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb1, zx_tmp_2d)
1600           CALL histwrite(nid_ins, "s_trmb1", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "s_trmb1", itau_w, zx_tmp_2d)
1601    
1602           CALL gr_fi_ecrit(1, klon, iim, (jjm + 1), s_trmb2, zx_tmp_2d)           CALL gr_fi_ecrit(1, klon, iim, jjm + 1, s_trmb2, zx_tmp_2d)
1603           CALL histwrite(nid_ins, "s_trmb2", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "s_trmb2", itau_w, zx_tmp_2d)
1604    
1605           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)
1606           CALL histwrite(nid_ins, "s_trmb3", itau_w, zx_tmp_2d)           CALL histwrite(nid_ins, "s_trmb3", itau_w, zx_tmp_2d)
1607    
          !IM cf. AM 081204 END  
   
1608           ! Champs 3D:           ! Champs 3D:
1609    
1610           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)
1611           CALL histwrite(nid_ins, "temp", itau_w, zx_tmp_3d)           CALL histwrite(nid_ins, "temp", itau_w, zx_tmp_3d)
1612    
1613           CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), u_seri, zx_tmp_3d)           CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, u_seri, zx_tmp_3d)
1614           CALL histwrite(nid_ins, "vitu", itau_w, zx_tmp_3d)           CALL histwrite(nid_ins, "vitu", itau_w, zx_tmp_3d)
1615    
1616           CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), v_seri, zx_tmp_3d)           CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, v_seri, zx_tmp_3d)
1617           CALL histwrite(nid_ins, "vitv", itau_w, zx_tmp_3d)           CALL histwrite(nid_ins, "vitv", itau_w, zx_tmp_3d)
1618    
1619           CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), zphi, zx_tmp_3d)           CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, zphi, zx_tmp_3d)
1620           CALL histwrite(nid_ins, "geop", itau_w, zx_tmp_3d)           CALL histwrite(nid_ins, "geop", itau_w, zx_tmp_3d)
1621    
1622           CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), play, zx_tmp_3d)           CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, play, zx_tmp_3d)
1623           CALL histwrite(nid_ins, "pres", itau_w, zx_tmp_3d)           CALL histwrite(nid_ins, "pres", itau_w, zx_tmp_3d)
1624    
1625           CALL gr_fi_ecrit(llm, klon, iim, (jjm + 1), d_t_vdf, zx_tmp_3d)           CALL gr_fi_ecrit(llm, klon, iim, jjm + 1, d_t_vdf, zx_tmp_3d)
1626           CALL histwrite(nid_ins, "dtvdf", itau_w, zx_tmp_3d)           CALL histwrite(nid_ins, "dtvdf", itau_w, zx_tmp_3d)
1627    
1628           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)
1629           CALL histwrite(nid_ins, "dqvdf", itau_w, zx_tmp_3d)           CALL histwrite(nid_ins, "dqvdf", itau_w, zx_tmp_3d)
1630    
1631           if (ok_sync) then           call histsync(nid_ins)
             call histsync(nid_ins)  
          endif  
1632        ENDIF        ENDIF
1633    
1634      end subroutine write_histins      end subroutine write_histins
1635    
     !****************************************************  
   
     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  
   
1636    END SUBROUTINE physiq    END SUBROUTINE physiq
1637    
1638  end module physiq_m  end module physiq_m

Legend:
Removed from v.49  
changed lines
  Added in v.174

  ViewVC Help
Powered by ViewVC 1.1.21