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

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

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

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

Legend:
Removed from v.61  
changed lines
  Added in v.189

  ViewVC Help
Powered by ViewVC 1.1.21