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

Legend:
Removed from v.47  
changed lines
  Added in v.182

  ViewVC Help
Powered by ViewVC 1.1.21