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

Legend:
Removed from v.34  
changed lines
  Added in v.169

  ViewVC Help
Powered by ViewVC 1.1.21