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

Diff of /trunk/phylmd/physiq.f

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

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

Legend:
Removed from v.18  
changed lines
  Added in v.175

  ViewVC Help
Powered by ViewVC 1.1.21