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

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

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

trunk/phylmd/clmain.f revision 82 by guez, Wed Mar 5 14:57:53 2014 UTC trunk/Sources/phylmd/clmain.f revision 228 by guez, Fri Nov 3 12:38:47 2017 UTC
# Line 4  module clmain_m Line 4  module clmain_m
4    
5  contains  contains
6    
7    SUBROUTINE clmain(dtime, itap, pctsrf, pctsrf_new, t, q, u, v, &    SUBROUTINE clmain(dtime, pctsrf, t, q, u, v, julien, mu0, ftsol, cdmmax, &
8         jour, rmu0, co2_ppm, ok_veget, ocean, ts, &         cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, qsol, paprs, pplay, fsnow, &
9         soil_model, cdmmax, cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, &         qsurf, evap, falbe, fluxlat, rain_fall, snow_f, fsolsw, fsollw, frugs, &
10         qsol, paprs, pplay, snow, qsurf, evap, albe, alblw, fluxlat, &         agesno, rugoro, d_t, d_q, d_u, d_v, d_ts, flux_t, flux_q, flux_u, &
11         rain_fall, snow_f, solsw, sollw, fder, rlon, rlat, &         flux_v, cdragh, cdragm, q2, dflux_t, dflux_q, ycoefh, t2m, q2m, &
12         rugos, debut, agesno, rugoro, d_t, d_q, d_u, d_v, &         u10m_srf, v10m_srf, pblh, capcl, oliqcl, cteicl, pblt, therm, trmb1, &
13         d_ts, flux_t, flux_q, flux_u, flux_v, cdragh, cdragm, q2, &         trmb2, trmb3, plcl, fqcalving, ffonte, run_off_lic_0)
        dflux_t, dflux_q, ycoefh, zu1, zv1, t2m, q2m, u10m, v10m, pblh, &  
        capcl, oliqcl, cteicl, pblt, therm, trmb1, trmb2, trmb3, plcl, &  
        fqcalving, ffonte, run_off_lic_0, flux_o, flux_g, tslab, seaice)  
14    
15      ! From phylmd/clmain.F, version 1.6 2005/11/16 14:47:19      ! From phylmd/clmain.F, version 1.6, 2005/11/16 14:47:19
16      ! Author: Z. X. Li (LMD/CNRS), date: 1993/08/18      ! Author: Z. X. Li (LMD/CNRS), date: 1993/08/18
17      ! Objet : interface de couche limite (diffusion verticale)      ! Objet : interface de couche limite (diffusion verticale)
18    
19      ! Tout ce qui a trait aux traceurs est dans "phytrac". Le calcul      ! Tout ce qui a trait aux traceurs est dans "phytrac". Le calcul
20      ! de la couche limite pour les traceurs se fait avec "cltrac" et      ! de la couche limite pour les traceurs se fait avec "cltrac" et
21      ! ne tient pas compte de la différentiation des sous-fractions de      ! ne tient pas compte de la diff\'erentiation des sous-fractions
22      ! sol.      ! de sol.
23    
     ! Pour pouvoir extraire les coefficients d'échanges et le vent  
     ! dans la première couche, trois champs ont été créés : "ycoefh",  
     ! "zu1" et "zv1". Nous avons moyenné les valeurs de ces trois  
     ! champs sur les quatre sous-surfaces du modèle.  
   
     use calendar, ONLY: ymds2ju  
24      use clqh_m, only: clqh      use clqh_m, only: clqh
25      use clvent_m, only: clvent      use clvent_m, only: clvent
26      use coefkz_m, only: coefkz      use coefkz_m, only: coefkz
27      use coefkzmin_m, only: coefkzmin      use coefkzmin_m, only: coefkzmin
28      USE conf_gcm_m, ONLY: prt_level      USE conf_gcm_m, ONLY: lmt_pas
29      USE conf_phys_m, ONLY: iflag_pbl      USE conf_phys_m, ONLY: iflag_pbl
     USE dimens_m, ONLY: iim, jjm  
30      USE dimphy, ONLY: klev, klon, zmasq      USE dimphy, ONLY: klev, klon, zmasq
31      USE dimsoil, ONLY: nsoilmx      USE dimsoil, ONLY: nsoilmx
     USE dynetat0_m, ONLY: day_ini  
     USE gath_cpl, ONLY: gath2cpl  
32      use hbtm_m, only: hbtm      use hbtm_m, only: hbtm
     USE histbeg_totreg_m, ONLY: histbeg_totreg  
     USE histdef_m, ONLY: histdef  
     USE histend_m, ONLY: histend  
     USE histsync_m, ONLY: histsync  
     use histwrite_m, only: histwrite  
33      USE indicesol, ONLY: epsfra, is_lic, is_oce, is_sic, is_ter, nbsrf      USE indicesol, ONLY: epsfra, is_lic, is_oce, is_sic, is_ter, nbsrf
34        USE interfoce_lim_m, ONLY: interfoce_lim
35        use stdlevvar_m, only: stdlevvar
36      USE suphec_m, ONLY: rd, rg, rkappa      USE suphec_m, ONLY: rd, rg, rkappa
37      USE temps, ONLY: annee_ref, itau_phy      use time_phylmdz, only: itap
38      use ustarhb_m, only: ustarhb      use ustarhb_m, only: ustarhb
     use vdif_kcay_m, only: vdif_kcay  
39      use yamada4_m, only: yamada4      use yamada4_m, only: yamada4
40    
     ! Arguments:  
   
41      REAL, INTENT(IN):: dtime ! interval du temps (secondes)      REAL, INTENT(IN):: dtime ! interval du temps (secondes)
     INTEGER, INTENT(IN):: itap ! numero du pas de temps  
     REAL, INTENT(inout):: pctsrf(klon, nbsrf)  
42    
43      ! la nouvelle repartition des surfaces sortie de l'interface      REAL, INTENT(inout):: pctsrf(klon, nbsrf)
44      REAL, INTENT(out):: pctsrf_new(klon, nbsrf)      ! tableau des pourcentages de surface de chaque maille
45    
46      REAL, INTENT(IN):: t(klon, klev) ! temperature (K)      REAL, INTENT(IN):: t(klon, klev) ! temperature (K)
47      REAL, INTENT(IN):: q(klon, klev) ! vapeur d'eau (kg/kg)      REAL, INTENT(IN):: q(klon, klev) ! vapeur d'eau (kg / kg)
48      REAL, INTENT(IN):: u(klon, klev), v(klon, klev) ! vitesse      REAL, INTENT(IN):: u(klon, klev), v(klon, klev) ! vitesse
49      INTEGER, INTENT(IN):: jour ! jour de l'annee en cours      INTEGER, INTENT(IN):: julien ! jour de l'annee en cours
50      REAL, intent(in):: rmu0(klon) ! cosinus de l'angle solaire zenithal          REAL, intent(in):: mu0(klon) ! cosinus de l'angle solaire zenithal    
51      REAL co2_ppm ! taux CO2 atmosphere      REAL, INTENT(IN):: ftsol(:, :) ! (klon, nbsrf) temp\'erature du sol (en K)
     LOGICAL ok_veget  
     CHARACTER(len=*), INTENT(IN):: ocean  
     REAL ts(klon, nbsrf) ! input-R- temperature du sol (en Kelvin)  
     LOGICAL, INTENT(IN):: soil_model  
52      REAL, INTENT(IN):: cdmmax, cdhmax ! seuils cdrm, cdrh      REAL, INTENT(IN):: cdmmax, cdhmax ! seuils cdrm, cdrh
53      REAL ksta, ksta_ter      REAL, INTENT(IN):: ksta, ksta_ter
54      LOGICAL ok_kzmin      LOGICAL, INTENT(IN):: ok_kzmin
55      REAL ftsoil(klon, nsoilmx, nbsrf)  
56      REAL qsol(klon)      REAL, INTENT(inout):: ftsoil(klon, nsoilmx, nbsrf)
57      REAL, INTENT(IN):: paprs(klon, klev+1) ! pression a intercouche (Pa)      ! soil temperature of surface fraction
58    
59        REAL, INTENT(inout):: qsol(:) ! (klon)
60        ! column-density of water in soil, in kg m-2
61    
62        REAL, INTENT(IN):: paprs(klon, klev + 1) ! pression a intercouche (Pa)
63      REAL, INTENT(IN):: pplay(klon, klev) ! pression au milieu de couche (Pa)      REAL, INTENT(IN):: pplay(klon, klev) ! pression au milieu de couche (Pa)
64      REAL snow(klon, nbsrf)      REAL, INTENT(inout):: fsnow(:, :) ! (klon, nbsrf) \'epaisseur neigeuse
65      REAL qsurf(klon, nbsrf)      REAL qsurf(klon, nbsrf)
66      REAL evap(klon, nbsrf)      REAL evap(klon, nbsrf)
67      REAL albe(klon, nbsrf)      REAL, intent(inout):: falbe(klon, nbsrf)
68      REAL alblw(klon, nbsrf)      REAL, intent(out):: fluxlat(:, :) ! (klon, nbsrf)
   
     REAL fluxlat(klon, nbsrf)  
69    
70      REAL, intent(in):: rain_fall(klon), snow_f(klon)      REAL, intent(in):: rain_fall(klon)
71      REAL, INTENT(IN):: solsw(klon, nbsrf), sollw(klon, nbsrf)      ! liquid water mass flux (kg / m2 / s), positive down
     REAL fder(klon)  
     REAL, INTENT(IN):: rlon(klon)  
     REAL, INTENT(IN):: rlat(klon) ! latitude en degrés  
72    
73      REAL rugos(klon, nbsrf)      REAL, intent(in):: snow_f(klon)
74      ! rugos----input-R- longeur de rugosite (en m)      ! solid water mass flux (kg / m2 / s), positive down
75    
76      LOGICAL, INTENT(IN):: debut      REAL, INTENT(IN):: fsolsw(klon, nbsrf), fsollw(klon, nbsrf)
77        REAL, intent(inout):: frugs(klon, nbsrf) ! longueur de rugosit\'e (en m)
78      real agesno(klon, nbsrf)      real agesno(klon, nbsrf)
79      REAL, INTENT(IN):: rugoro(klon)      REAL, INTENT(IN):: rugoro(klon)
80    
# Line 108  contains Line 85  contains
85      REAL, intent(out):: d_u(klon, klev), d_v(klon, klev)      REAL, intent(out):: d_u(klon, klev), d_v(klon, klev)
86      ! changement pour "u" et "v"      ! changement pour "u" et "v"
87    
88      REAL d_ts(klon, nbsrf)      REAL, intent(out):: d_ts(:, :) ! (klon, nbsrf) variation of ftsol
89      ! d_ts-----output-R- le changement pour "ts"  
90        REAL, intent(out):: flux_t(klon, nbsrf)
91        ! flux de chaleur sensible (Cp T) (W / m2) (orientation positive vers
92        ! le bas) à la surface
93    
94        REAL, intent(out):: flux_q(klon, nbsrf)
95        ! flux de vapeur d'eau (kg / m2 / s) à la surface
96    
97      REAL flux_t(klon, klev, nbsrf), flux_q(klon, klev, nbsrf)      REAL, intent(out):: flux_u(klon, nbsrf), flux_v(klon, nbsrf)
98      ! flux_t---output-R- flux de chaleur sensible (CpT) J/m**2/s (W/m**2)      ! tension du vent à la surface, en Pa
     !                    (orientation positive vers le bas)  
     ! flux_q---output-R- flux de vapeur d'eau (kg/m**2/s)  
   
     REAL flux_u(klon, klev, nbsrf), flux_v(klon, klev, nbsrf)  
     ! flux_u---output-R- tension du vent X: (kg m/s)/(m**2 s) ou Pascal  
     ! flux_v---output-R- tension du vent Y: (kg m/s)/(m**2 s) ou Pascal  
99    
100      REAL, INTENT(out):: cdragh(klon), cdragm(klon)      REAL, INTENT(out):: cdragh(klon), cdragm(klon)
101      real q2(klon, klev+1, nbsrf)      real q2(klon, klev + 1, nbsrf)
102    
103      REAL dflux_t(klon), dflux_q(klon)      REAL, INTENT(out):: dflux_t(klon), dflux_q(klon)
104      ! dflux_t derive du flux sensible      ! dflux_t derive du flux sensible
105      ! dflux_q derive du flux latent      ! dflux_q derive du flux latent
106      !IM "slab" ocean      ! IM "slab" ocean
107    
108      REAL, intent(out):: ycoefh(klon, klev)      REAL, intent(out):: ycoefh(klon, klev)
109      REAL, intent(out):: zu1(klon)      ! Pour pouvoir extraire les coefficients d'\'echange, le champ
110      REAL zv1(klon)      ! "ycoefh" a \'et\'e cr\'e\'e. Nous avons moyenn\'e les valeurs de
111      REAL t2m(klon, nbsrf), q2m(klon, nbsrf)      ! ce champ sur les quatre sous-surfaces du mod\`ele.
112      REAL u10m(klon, nbsrf), v10m(klon, nbsrf)  
113        REAL, INTENT(inout):: t2m(klon, nbsrf), q2m(klon, nbsrf)
114      !IM cf. AM : pbl, hbtm (Comme les autres diagnostics on cumule ds  
115      ! physiq ce qui permet de sortir les grdeurs par sous surface)      REAL, INTENT(inout):: u10m_srf(:, :), v10m_srf(:, :) ! (klon, nbsrf)
116      REAL pblh(klon, nbsrf)      ! composantes du vent \`a 10m sans spirale d'Ekman
117      ! pblh------- HCL  
118        ! Ionela Musat. Cf. Anne Mathieu : planetary boundary layer, hbtm.
119        ! Comme les autres diagnostics on cumule dans physiq ce qui permet
120        ! de sortir les grandeurs par sous-surface.
121        REAL pblh(klon, nbsrf) ! height of planetary boundary layer
122      REAL capcl(klon, nbsrf)      REAL capcl(klon, nbsrf)
123      REAL oliqcl(klon, nbsrf)      REAL oliqcl(klon, nbsrf)
124      REAL cteicl(klon, nbsrf)      REAL cteicl(klon, nbsrf)
125      REAL pblt(klon, nbsrf)      REAL, INTENT(inout):: pblt(klon, nbsrf) ! T au nveau HCL
     ! pblT------- T au nveau HCL  
126      REAL therm(klon, nbsrf)      REAL therm(klon, nbsrf)
127      REAL trmb1(klon, nbsrf)      REAL trmb1(klon, nbsrf)
128      ! trmb1-------deep_cape      ! trmb1-------deep_cape
# Line 154  contains Line 134  contains
134      REAL fqcalving(klon, nbsrf), ffonte(klon, nbsrf)      REAL fqcalving(klon, nbsrf), ffonte(klon, nbsrf)
135      ! ffonte----Flux thermique utilise pour fondre la neige      ! ffonte----Flux thermique utilise pour fondre la neige
136      ! fqcalving-Flux d'eau "perdue" par la surface et necessaire pour limiter la      ! fqcalving-Flux d'eau "perdue" par la surface et necessaire pour limiter la
137      !           hauteur de neige, en kg/m2/s      !           hauteur de neige, en kg / m2 / s
138      REAL run_off_lic_0(klon)      REAL run_off_lic_0(klon)
139    
140      REAL flux_o(klon), flux_g(klon)      ! Local:
     !IM "slab" ocean  
     ! flux_g---output-R-  flux glace (pour OCEAN='slab  ')  
     ! flux_o---output-R-  flux ocean (pour OCEAN='slab  ')  
   
     REAL tslab(klon)  
     ! tslab-in/output-R temperature du slab ocean (en Kelvin)  
     ! uniqmnt pour slab  
141    
142      REAL seaice(klon)      LOGICAL:: firstcal = .true.
     ! seaice---output-R-  glace de mer (kg/m2) (pour OCEAN='slab  ')  
143    
144      ! Local:      ! la nouvelle repartition des surfaces sortie de l'interface
145        REAL, save:: pctsrf_new_oce(klon)
146        REAL, save:: pctsrf_new_sic(klon)
147    
     REAL y_flux_o(klon), y_flux_g(klon)  
     real ytslab(klon)  
     real y_seaice(klon)  
148      REAL y_fqcalving(klon), y_ffonte(klon)      REAL y_fqcalving(klon), y_ffonte(klon)
149      real y_run_off_lic_0(klon)      real y_run_off_lic_0(klon)
   
150      REAL rugmer(klon)      REAL rugmer(klon)
   
151      REAL ytsoil(klon, nsoilmx)      REAL ytsoil(klon, nsoilmx)
   
152      REAL yts(klon), yrugos(klon), ypct(klon), yz0_new(klon)      REAL yts(klon), yrugos(klon), ypct(klon), yz0_new(klon)
153      REAL yalb(klon)      REAL yalb(klon)
154      REAL yalblw(klon)      REAL snow(klon), yqsurf(klon), yagesno(klon)
155      REAL yu1(klon), yv1(klon)      real yqsol(klon) ! column-density of water in soil, in kg m-2
156      ! on rajoute en output yu1 et yv1 qui sont les vents dans      REAL yrain_f(klon) ! liquid water mass flux (kg / m2 / s), positive down
157      ! la premiere couche      REAL ysnow_f(klon) ! solid water mass flux (kg / m2 / s), positive down
     REAL ysnow(klon), yqsurf(klon), yagesno(klon), yqsol(klon)  
     REAL yrain_f(klon), ysnow_f(klon)  
     REAL ysollw(klon), ysolsw(klon)  
     REAL yfder(klon), ytaux(klon), ytauy(klon)  
158      REAL yrugm(klon), yrads(klon), yrugoro(klon)      REAL yrugm(klon), yrads(klon), yrugoro(klon)
   
159      REAL yfluxlat(klon)      REAL yfluxlat(klon)
   
160      REAL y_d_ts(klon)      REAL y_d_ts(klon)
161      REAL y_d_t(klon, klev), y_d_q(klon, klev)      REAL y_d_t(klon, klev), y_d_q(klon, klev)
162      REAL y_d_u(klon, klev), y_d_v(klon, klev)      REAL y_d_u(klon, klev), y_d_v(klon, klev)
163      REAL y_flux_t(klon, klev), y_flux_q(klon, klev)      REAL y_flux_t(klon), y_flux_q(klon)
164      REAL y_flux_u(klon, klev), y_flux_v(klon, klev)      REAL y_flux_u(klon), y_flux_v(klon)
165      REAL y_dflux_t(klon), y_dflux_q(klon)      REAL y_dflux_t(klon), y_dflux_q(klon)
166      REAL coefh(klon, klev), coefm(klon, klev)      REAL coefh(klon, klev), coefm(klon, klev)
167      REAL yu(klon, klev), yv(klon, klev)      REAL yu(klon, klev), yv(klon, klev)
168      REAL yt(klon, klev), yq(klon, klev)      REAL yt(klon, klev), yq(klon, klev)
169      REAL ypaprs(klon, klev+1), ypplay(klon, klev), ydelp(klon, klev)      REAL ypaprs(klon, klev + 1), ypplay(klon, klev), ydelp(klon, klev)
   
170      REAL ycoefm0(klon, klev), ycoefh0(klon, klev)      REAL ycoefm0(klon, klev), ycoefh0(klon, klev)
171        REAL yzlay(klon, klev), zlev(klon, klev + 1), yteta(klon, klev)
172      REAL yzlay(klon, klev), yzlev(klon, klev+1), yteta(klon, klev)      REAL ykmm(klon, klev + 1), ykmn(klon, klev + 1)
173      REAL ykmm(klon, klev+1), ykmn(klon, klev+1)      REAL ykmq(klon, klev + 1)
174      REAL ykmq(klon, klev+1)      REAL yq2(klon, klev + 1)
     REAL yq2(klon, klev+1)  
     REAL q2diag(klon, klev+1)  
   
     REAL u1lay(klon), v1lay(klon)  
175      REAL delp(klon, klev)      REAL delp(klon, klev)
176      INTEGER i, k, nsrf      INTEGER i, k, nsrf
   
177      INTEGER ni(klon), knon, j      INTEGER ni(klon), knon, j
178    
179      REAL pctsrf_pot(klon, nbsrf)      REAL pctsrf_pot(klon, nbsrf)
180      ! "pourcentage potentiel" pour tenir compte des éventuelles      ! "pourcentage potentiel" pour tenir compte des \'eventuelles
181      ! apparitions ou disparitions de la glace de mer      ! apparitions ou disparitions de la glace de mer
182    
183      REAL zx_alf1, zx_alf2 !valeur ambiante par extrapola.      REAL yt2m(klon), yq2m(klon), wind10m(klon)
184        REAL ustar(klon)
     ! maf pour sorties IOISPL en cas de debugagage  
   
     CHARACTER(80) cldebug  
     SAVE cldebug  
     CHARACTER(8) cl_surf(nbsrf)  
     SAVE cl_surf  
     INTEGER nhoridbg, nidbg  
     SAVE nhoridbg, nidbg  
     INTEGER ndexbg(iim*(jjm+1))  
     REAL zx_lon(iim, jjm+1), zx_lat(iim, jjm+1), zjulian  
     REAL tabindx(klon)  
     REAL debugtab(iim, jjm+1)  
     LOGICAL first_appel  
     SAVE first_appel  
     DATA first_appel/ .TRUE./  
     LOGICAL:: debugindex = .FALSE.  
     INTEGER idayref  
   
     REAL yt2m(klon), yq2m(klon), yu10m(klon)  
     REAL yustar(klon)  
     ! -- LOOP  
     REAL yu10mx(klon)  
     REAL yu10my(klon)  
     REAL ywindsp(klon)  
     ! -- LOOP  
185    
186      REAL yt10m(klon), yq10m(klon)      REAL yt10m(klon), yq10m(klon)
187      REAL ypblh(klon)      REAL ypblh(klon)
# Line 263  contains Line 194  contains
194      REAL ytrmb1(klon)      REAL ytrmb1(klon)
195      REAL ytrmb2(klon)      REAL ytrmb2(klon)
196      REAL ytrmb3(klon)      REAL ytrmb3(klon)
197      REAL uzon(klon), vmer(klon)      REAL u1(klon), v1(klon)
198      REAL tair1(klon), qair1(klon), tairsol(klon)      REAL tair1(klon), qair1(klon), tairsol(klon)
199      REAL psfce(klon), patm(klon)      REAL psfce(klon), patm(klon)
200    
201      REAL qairsol(klon), zgeo1(klon)      REAL qairsol(klon), zgeo1(klon)
202      REAL rugo1(klon)      REAL rugo1(klon)
203    
     ! utiliser un jeu de fonctions simples                
     LOGICAL zxli  
     PARAMETER (zxli=.FALSE.)  
   
204      !------------------------------------------------------------      !------------------------------------------------------------
205    
206      ytherm = 0.      ytherm = 0.
207    
     IF (debugindex .AND. first_appel) THEN  
        first_appel = .FALSE.  
   
        ! initialisation sorties netcdf  
   
        idayref = day_ini  
        CALL ymds2ju(annee_ref, 1, idayref, 0., zjulian)  
        CALL gr_fi_ecrit(1, klon, iim, jjm+1, rlon, zx_lon)  
        DO i = 1, iim  
           zx_lon(i, 1) = rlon(i+1)  
           zx_lon(i, jjm+1) = rlon(i+1)  
        END DO  
        CALL gr_fi_ecrit(1, klon, iim, jjm+1, rlat, zx_lat)  
        cldebug = 'sous_index'  
        CALL histbeg_totreg(cldebug, zx_lon(:, 1), zx_lat(1, :), 1, &  
             iim, 1, jjm+1, itau_phy, zjulian, dtime, nhoridbg, nidbg)  
        ! no vertical axis  
        cl_surf(1) = 'ter'  
        cl_surf(2) = 'lic'  
        cl_surf(3) = 'oce'  
        cl_surf(4) = 'sic'  
        DO nsrf = 1, nbsrf  
           CALL histdef(nidbg, cl_surf(nsrf), cl_surf(nsrf), '-', iim, jjm+1, &  
                nhoridbg, 1, 1, 1, -99, 'inst', dtime, dtime)  
        END DO  
        CALL histend(nidbg)  
        CALL histsync(nidbg)  
     END IF  
   
208      DO k = 1, klev ! epaisseur de couche      DO k = 1, klev ! epaisseur de couche
209         DO i = 1, klon         DO i = 1, klon
210            delp(i, k) = paprs(i, k) - paprs(i, k+1)            delp(i, k) = paprs(i, k) - paprs(i, k + 1)
211         END DO         END DO
212      END DO      END DO
     DO i = 1, klon ! vent de la premiere couche  
        zx_alf1 = 1.0  
        zx_alf2 = 1.0 - zx_alf1  
        u1lay(i) = u(i, 1)*zx_alf1 + u(i, 2)*zx_alf2  
        v1lay(i) = v(i, 1)*zx_alf1 + v(i, 2)*zx_alf2  
     END DO  
213    
214      ! Initialization:      ! Initialization:
215      rugmer = 0.      rugmer = 0.
# Line 325  contains Line 217  contains
217      cdragm = 0.      cdragm = 0.
218      dflux_t = 0.      dflux_t = 0.
219      dflux_q = 0.      dflux_q = 0.
     zu1 = 0.  
     zv1 = 0.  
220      ypct = 0.      ypct = 0.
     yts = 0.  
     ysnow = 0.  
221      yqsurf = 0.      yqsurf = 0.
     yalb = 0.  
     yalblw = 0.  
222      yrain_f = 0.      yrain_f = 0.
223      ysnow_f = 0.      ysnow_f = 0.
     yfder = 0.  
     ytaux = 0.  
     ytauy = 0.  
     ysolsw = 0.  
     ysollw = 0.  
224      yrugos = 0.      yrugos = 0.
     yu1 = 0.  
     yv1 = 0.  
     yrads = 0.  
225      ypaprs = 0.      ypaprs = 0.
226      ypplay = 0.      ypplay = 0.
227      ydelp = 0.      ydelp = 0.
# Line 351  contains Line 229  contains
229      yv = 0.      yv = 0.
230      yt = 0.      yt = 0.
231      yq = 0.      yq = 0.
     pctsrf_new = 0.  
     y_flux_u = 0.  
     y_flux_v = 0.  
     !$$ PB  
232      y_dflux_t = 0.      y_dflux_t = 0.
233      y_dflux_q = 0.      y_dflux_q = 0.
     ytsoil = 999999.  
234      yrugoro = 0.      yrugoro = 0.
     ! -- LOOP  
     yu10mx = 0.  
     yu10my = 0.  
     ywindsp = 0.  
     ! -- LOOP  
235      d_ts = 0.      d_ts = 0.
     !§§§ PB  
     yfluxlat = 0.  
236      flux_t = 0.      flux_t = 0.
237      flux_q = 0.      flux_q = 0.
238      flux_u = 0.      flux_u = 0.
239      flux_v = 0.      flux_v = 0.
240        fluxlat = 0.
241      d_t = 0.      d_t = 0.
242      d_q = 0.      d_q = 0.
243      d_u = 0.      d_u = 0.
244      d_v = 0.      d_v = 0.
245      ycoefh = 0.      ycoefh = 0.
246    
247      ! Boucler sur toutes les sous-fractions du sol:      ! Initialisation des "pourcentages potentiels". On consid\`ere ici qu'on
248        ! peut avoir potentiellement de la glace sur tout le domaine oc\'eanique
249      ! Initialisation des "pourcentages potentiels". On considère ici qu'on      ! (\`a affiner)
     ! peut avoir potentiellement de la glace sur tout le domaine océanique  
     ! (à affiner)  
250    
251      pctsrf_pot = pctsrf      pctsrf_pot(:, is_ter) = pctsrf(:, is_ter)
252        pctsrf_pot(:, is_lic) = pctsrf(:, is_lic)
253      pctsrf_pot(:, is_oce) = 1. - zmasq      pctsrf_pot(:, is_oce) = 1. - zmasq
254      pctsrf_pot(:, is_sic) = 1. - zmasq      pctsrf_pot(:, is_sic) = 1. - zmasq
255    
256        ! Tester si c'est le moment de lire le fichier:
257        if (mod(itap - 1, lmt_pas) == 0) then
258           CALL interfoce_lim(julien, pctsrf_new_oce, pctsrf_new_sic)
259        endif
260    
261        ! Boucler sur toutes les sous-fractions du sol:
262    
263      loop_surface: DO nsrf = 1, nbsrf      loop_surface: DO nsrf = 1, nbsrf
264         ! Chercher les indices :         ! Chercher les indices :
265         ni = 0         ni = 0
266         knon = 0         knon = 0
267         DO i = 1, klon         DO i = 1, klon
268            ! Pour déterminer le domaine à traiter, on utilise les surfaces            ! Pour d\'eterminer le domaine \`a traiter, on utilise les surfaces
269            ! "potentielles"            ! "potentielles"
270            IF (pctsrf_pot(i, nsrf) > epsfra) THEN            IF (pctsrf_pot(i, nsrf) > epsfra) THEN
271               knon = knon + 1               knon = knon + 1
# Line 400  contains Line 273  contains
273            END IF            END IF
274         END DO         END DO
275    
        ! variables pour avoir une sortie IOIPSL des INDEX  
        IF (debugindex) THEN  
           tabindx = 0.  
           DO i = 1, knon  
              tabindx(i) = real(i)  
           END DO  
           debugtab = 0.  
           ndexbg = 0  
           CALL gath2cpl(tabindx, debugtab, klon, knon, iim, jjm, ni)  
           CALL histwrite(nidbg, cl_surf(nsrf), itap, debugtab)  
        END IF  
   
276         if_knon: IF (knon /= 0) then         if_knon: IF (knon /= 0) then
277            DO j = 1, knon            DO j = 1, knon
278               i = ni(j)               i = ni(j)
279               ypct(j) = pctsrf(i, nsrf)               ypct(j) = pctsrf(i, nsrf)
280               yts(j) = ts(i, nsrf)               yts(j) = ftsol(i, nsrf)
281               ytslab(i) = tslab(i)               snow(j) = fsnow(i, nsrf)
              ysnow(j) = snow(i, nsrf)  
282               yqsurf(j) = qsurf(i, nsrf)               yqsurf(j) = qsurf(i, nsrf)
283               yalb(j) = albe(i, nsrf)               yalb(j) = falbe(i, nsrf)
              yalblw(j) = alblw(i, nsrf)  
284               yrain_f(j) = rain_fall(i)               yrain_f(j) = rain_fall(i)
285               ysnow_f(j) = snow_f(i)               ysnow_f(j) = snow_f(i)
286               yagesno(j) = agesno(i, nsrf)               yagesno(j) = agesno(i, nsrf)
287               yfder(j) = fder(i)               yrugos(j) = frugs(i, nsrf)
              ytaux(j) = flux_u(i, 1, nsrf)  
              ytauy(j) = flux_v(i, 1, nsrf)  
              ysolsw(j) = solsw(i, nsrf)  
              ysollw(j) = sollw(i, nsrf)  
              yrugos(j) = rugos(i, nsrf)  
288               yrugoro(j) = rugoro(i)               yrugoro(j) = rugoro(i)
289               yu1(j) = u1lay(i)               yrads(j) = fsolsw(i, nsrf) + fsollw(i, nsrf)
290               yv1(j) = v1lay(i)               ypaprs(j, klev + 1) = paprs(i, klev + 1)
              yrads(j) = ysolsw(j) + ysollw(j)  
              ypaprs(j, klev+1) = paprs(i, klev+1)  
291               y_run_off_lic_0(j) = run_off_lic_0(i)               y_run_off_lic_0(j) = run_off_lic_0(i)
              yu10mx(j) = u10m(i, nsrf)  
              yu10my(j) = v10m(i, nsrf)  
              ywindsp(j) = sqrt(yu10mx(j)*yu10mx(j)+yu10my(j)*yu10my(j))  
292            END DO            END DO
293    
294            ! IF bucket model for continent, copy soil water content            ! For continent, copy soil water content
295            IF (nsrf == is_ter .AND. .NOT. ok_veget) THEN            IF (nsrf == is_ter) yqsol(:knon) = qsol(ni(:knon))
              DO j = 1, knon  
                 i = ni(j)  
                 yqsol(j) = qsol(i)  
              END DO  
           ELSE  
              yqsol = 0.  
           END IF  
296    
297            DO k = 1, nsoilmx            ytsoil(:knon, :) = ftsoil(ni(:knon), :, nsrf)
              DO j = 1, knon  
                 i = ni(j)  
                 ytsoil(j, k) = ftsoil(i, k, nsrf)  
              END DO  
           END DO  
298    
299            DO k = 1, klev            DO k = 1, klev
300               DO j = 1, knon               DO j = 1, knon
# Line 473  contains Line 310  contains
310            END DO            END DO
311    
312            ! calculer Cdrag et les coefficients d'echange            ! calculer Cdrag et les coefficients d'echange
313            CALL coefkz(nsrf, knon, ypaprs, ypplay, ksta, ksta_ter, yts, yrugos, &            CALL coefkz(nsrf, ypaprs, ypplay, ksta, ksta_ter, yts(:knon), &
314                 yu, yv, yt, yq, yqsurf, coefm(:knon, :), coefh(:knon, :))                 yrugos, yu, yv, yt, yq, yqsurf(:knon), coefm(:knon, :), &
315                   coefh(:knon, :))
316    
317            IF (iflag_pbl == 1) THEN            IF (iflag_pbl == 1) THEN
318               CALL coefkz2(nsrf, knon, ypaprs, ypplay, yt, ycoefm0, ycoefh0)               CALL coefkz2(nsrf, knon, ypaprs, ypplay, yt, ycoefm0, ycoefh0)
319               coefm(:knon, :) = max(coefm(:knon, :), ycoefm0(:knon, :))               coefm(:knon, :) = max(coefm(:knon, :), ycoefm0(:knon, :))
# Line 493  contains Line 332  contains
332                    coefm(:knon, 1), ycoefm0, ycoefh0)                    coefm(:knon, 1), ycoefm0, ycoefh0)
333               coefm(:knon, :) = max(coefm(:knon, :), ycoefm0(:knon, :))               coefm(:knon, :) = max(coefm(:knon, :), ycoefm0(:knon, :))
334               coefh(:knon, :) = max(coefh(:knon, :), ycoefh0(:knon, :))               coefh(:knon, :) = max(coefh(:knon, :), ycoefh0(:knon, :))
335             END IF            END IF
336    
337            IF (iflag_pbl >= 3) THEN            IF (iflag_pbl >= 6) THEN
338               ! MELLOR ET YAMADA adapté à Mars, Richard Fournier et               ! Mellor et Yamada adapt\'e \`a Mars, Richard Fournier et
339               ! Frédéric Hourdin               ! Fr\'ed\'eric Hourdin
340               yzlay(:knon, 1) = rd * yt(:knon, 1) / (0.5 * (ypaprs(:knon, 1) &               yzlay(:knon, 1) = rd * yt(:knon, 1) / (0.5 * (ypaprs(:knon, 1) &
341                    + ypplay(:knon, 1))) &                    + ypplay(:knon, 1))) &
342                    * (ypaprs(:knon, 1) - ypplay(:knon, 1)) / rg                    * (ypaprs(:knon, 1) - ypplay(:knon, 1)) / rg
343    
344               DO k = 2, klev               DO k = 2, klev
345                  yzlay(1:knon, k) = yzlay(1:knon, k-1) &                  yzlay(:knon, k) = yzlay(:knon, k-1) &
346                       + rd * 0.5 * (yt(1:knon, k-1) + yt(1:knon, k)) &                       + rd * 0.5 * (yt(1:knon, k-1) + yt(1:knon, k)) &
347                       / ypaprs(1:knon, k) &                       / ypaprs(1:knon, k) &
348                       * (ypplay(1:knon, k-1) - ypplay(1:knon, k)) / rg                       * (ypplay(1:knon, k-1) - ypplay(1:knon, k)) / rg
349               END DO               END DO
350    
351               DO k = 1, klev               DO k = 1, klev
352                  yteta(1:knon, k) = yt(1:knon, k)*(ypaprs(1:knon, 1) &                  yteta(1:knon, k) = yt(1:knon, k) * (ypaprs(1:knon, 1) &
353                       / ypplay(1:knon, k))**rkappa * (1.+0.61*yq(1:knon, k))                       / ypplay(1:knon, k))**rkappa * (1. + 0.61 * yq(1:knon, k))
354               END DO               END DO
355               yzlev(1:knon, 1) = 0.  
356               yzlev(:knon, klev+1) = 2. * yzlay(:knon, klev) &               zlev(:knon, 1) = 0.
357                 zlev(:knon, klev + 1) = 2. * yzlay(:knon, klev) &
358                    - yzlay(:knon, klev - 1)                    - yzlay(:knon, klev - 1)
359    
360               DO k = 2, klev               DO k = 2, klev
361                  yzlev(1:knon, k) = 0.5*(yzlay(1:knon, k)+yzlay(1:knon, k-1))                  zlev(:knon, k) = 0.5 * (yzlay(:knon, k) + yzlay(:knon, k-1))
362               END DO               END DO
363    
364               DO k = 1, klev + 1               DO k = 1, klev + 1
365                  DO j = 1, knon                  DO j = 1, knon
366                     i = ni(j)                     i = ni(j)
# Line 524  contains Line 368  contains
368                  END DO                  END DO
369               END DO               END DO
370    
371               CALL ustarhb(knon, yu, yv, coefm(:knon, 1), yustar)               ustar(:knon) = ustarhb(yu(:knon, 1), yv(:knon, 1), coefm(:knon, 1))
372    
373               IF (prt_level > 9) THEN               ! iflag_pbl peut \^etre utilis\'e comme longueur de m\'elange
                 PRINT *, 'USTAR = ', yustar  
              END IF  
374    
375               ! iflag_pbl peut être utilisé comme longueur de mélange               CALL yamada4(dtime, rg, zlev(:knon, :), yzlay(:knon, :), &
376                      yu(:knon, :), yv(:knon, :), yteta(:knon, :), &
377               IF (iflag_pbl >= 11) THEN                    coefm(:knon, 1), yq2(:knon, :), ykmm(:knon, :), &
378                  CALL vdif_kcay(knon, dtime, rg, rd, ypaprs, yt, yzlev, yzlay, &                    ykmn(:knon, :), ykmq(:knon, :), ustar(:knon), iflag_pbl)
                      yu, yv, yteta, coefm(:knon, 1), yq2, q2diag, ykmm, ykmn, &  
                      yustar, iflag_pbl)  
              ELSE  
                 CALL yamada4(knon, dtime, rg, yzlev, yzlay, yu, yv, yteta, &  
                      coefm(:knon, 1), yq2, ykmm, ykmn, ykmq, yustar, iflag_pbl)  
              END IF  
379    
380               coefm(:knon, 2:) = ykmm(:knon, 2:klev)               coefm(:knon, 2:) = ykmm(:knon, 2:klev)
381               coefh(:knon, 2:) = ykmn(:knon, 2:klev)               coefh(:knon, 2:) = ykmn(:knon, 2:klev)
382            END IF            END IF
383    
384            ! calculer la diffusion des vitesses "u" et "v"            ! calculer la diffusion des vitesses "u" et "v"
385            CALL clvent(knon, dtime, yu1, yv1, coefm(:knon, :), yt, yu, ypaprs, &            CALL clvent(knon, dtime, yu(:knon, 1), yv(:knon, 1), &
386                 ypplay, ydelp, y_d_u, y_flux_u)                 coefm(:knon, :), yt, yu, ypaprs, ypplay, ydelp, y_d_u, &
387            CALL clvent(knon, dtime, yu1, yv1, coefm(:knon, :), yt, yv, ypaprs, &                 y_flux_u(:knon))
388                 ypplay, ydelp, y_d_v, y_flux_v)            CALL clvent(knon, dtime, yu(:knon, 1), yv(:knon, 1), &
389                   coefm(:knon, :), yt, yv, ypaprs, ypplay, ydelp, y_d_v, &
390            ! pour le couplage                 y_flux_v(:knon))
           ytaux = y_flux_u(:, 1)  
           ytauy = y_flux_v(:, 1)  
391    
392            ! calculer la diffusion de "q" et de "h"            ! calculer la diffusion de "q" et de "h"
393            CALL clqh(dtime, itap, jour, debut, rlat, knon, nsrf, ni, pctsrf, &            CALL clqh(dtime, julien, firstcal, nsrf, ni(:knon), &
394                 soil_model, ytsoil, yqsol, ok_veget, ocean, rmu0, co2_ppm, &                 ytsoil(:knon, :), yqsol(:knon), mu0, yrugos, yrugoro, &
395                 yrugos, yrugoro, yu1, yv1, coefh(:knon, :), yt, yq, yts, &                 yu(:knon, 1), yv(:knon, 1), coefh(:knon, :), yt, yq, &
396                 ypaprs, ypplay, ydelp, yrads, yalb, yalblw, ysnow, yqsurf, &                 yts(:knon), ypaprs, ypplay, ydelp, yrads(:knon), yalb(:knon), &
397                 yrain_f, ysnow_f, yfder, ysolsw, yfluxlat, pctsrf_new, &                 snow(:knon), yqsurf, yrain_f, ysnow_f, yfluxlat(:knon), &
398                 yagesno, y_d_t, y_d_q, y_d_ts, yz0_new, y_flux_t, y_flux_q, &                 pctsrf_new_sic, yagesno(:knon), y_d_t, y_d_q, y_d_ts(:knon), &
399                 y_dflux_t, y_dflux_q, y_fqcalving, y_ffonte, y_run_off_lic_0, &                 yz0_new, y_flux_t(:knon), y_flux_q(:knon), y_dflux_t(:knon), &
400                 y_flux_o, y_flux_g, ytslab, y_seaice)                 y_dflux_q(:knon), y_fqcalving, y_ffonte, y_run_off_lic_0)
401    
402            ! calculer la longueur de rugosite sur ocean            ! calculer la longueur de rugosite sur ocean
403            yrugm = 0.            yrugm = 0.
404            IF (nsrf == is_oce) THEN            IF (nsrf == is_oce) THEN
405               DO j = 1, knon               DO j = 1, knon
406                  yrugm(j) = 0.018*coefm(j, 1)*(yu1(j)**2+yv1(j)**2)/rg + &                  yrugm(j) = 0.018 * coefm(j, 1) * (yu(j, 1)**2 + yv(j, 1)**2) &
407                       0.11*14E-6/sqrt(coefm(j, 1)*(yu1(j)**2+yv1(j)**2))                       / rg + 0.11 * 14E-6 &
408                         / sqrt(coefm(j, 1) * (yu(j, 1)**2 + yv(j, 1)**2))
409                  yrugm(j) = max(1.5E-05, yrugm(j))                  yrugm(j) = max(1.5E-05, yrugm(j))
410               END DO               END DO
411            END IF            END IF
412            DO j = 1, knon            DO j = 1, knon
413               y_dflux_t(j) = y_dflux_t(j)*ypct(j)               y_dflux_t(j) = y_dflux_t(j) * ypct(j)
414               y_dflux_q(j) = y_dflux_q(j)*ypct(j)               y_dflux_q(j) = y_dflux_q(j) * ypct(j)
              yu1(j) = yu1(j)*ypct(j)  
              yv1(j) = yv1(j)*ypct(j)  
415            END DO            END DO
416    
417            DO k = 1, klev            DO k = 1, klev
418               DO j = 1, knon               DO j = 1, knon
419                  i = ni(j)                  i = ni(j)
420                  coefh(j, k) = coefh(j, k)*ypct(j)                  coefh(j, k) = coefh(j, k) * ypct(j)
421                  coefm(j, k) = coefm(j, k)*ypct(j)                  coefm(j, k) = coefm(j, k) * ypct(j)
422                  y_d_t(j, k) = y_d_t(j, k)*ypct(j)                  y_d_t(j, k) = y_d_t(j, k) * ypct(j)
423                  y_d_q(j, k) = y_d_q(j, k)*ypct(j)                  y_d_q(j, k) = y_d_q(j, k) * ypct(j)
424                  flux_t(i, k, nsrf) = y_flux_t(j, k)                  y_d_u(j, k) = y_d_u(j, k) * ypct(j)
425                  flux_q(i, k, nsrf) = y_flux_q(j, k)                  y_d_v(j, k) = y_d_v(j, k) * ypct(j)
                 flux_u(i, k, nsrf) = y_flux_u(j, k)  
                 flux_v(i, k, nsrf) = y_flux_v(j, k)  
                 y_d_u(j, k) = y_d_u(j, k)*ypct(j)  
                 y_d_v(j, k) = y_d_v(j, k)*ypct(j)  
426               END DO               END DO
427            END DO            END DO
428    
429            evap(:, nsrf) = -flux_q(:, 1, nsrf)            flux_t(ni(:knon), nsrf) = y_flux_t(:knon)
430              flux_q(ni(:knon), nsrf) = y_flux_q(:knon)
431              flux_u(ni(:knon), nsrf) = y_flux_u(:knon)
432              flux_v(ni(:knon), nsrf) = y_flux_v(:knon)
433    
434            albe(:, nsrf) = 0.            evap(:, nsrf) = -flux_q(:, nsrf)
435            alblw(:, nsrf) = 0.  
436            snow(:, nsrf) = 0.            falbe(:, nsrf) = 0.
437              fsnow(:, nsrf) = 0.
438            qsurf(:, nsrf) = 0.            qsurf(:, nsrf) = 0.
439            rugos(:, nsrf) = 0.            frugs(:, nsrf) = 0.
           fluxlat(:, nsrf) = 0.  
440            DO j = 1, knon            DO j = 1, knon
441               i = ni(j)               i = ni(j)
442               d_ts(i, nsrf) = y_d_ts(j)               d_ts(i, nsrf) = y_d_ts(j)
443               albe(i, nsrf) = yalb(j)               falbe(i, nsrf) = yalb(j)
444               alblw(i, nsrf) = yalblw(j)               fsnow(i, nsrf) = snow(j)
              snow(i, nsrf) = ysnow(j)  
445               qsurf(i, nsrf) = yqsurf(j)               qsurf(i, nsrf) = yqsurf(j)
446               rugos(i, nsrf) = yz0_new(j)               frugs(i, nsrf) = yz0_new(j)
447               fluxlat(i, nsrf) = yfluxlat(j)               fluxlat(i, nsrf) = yfluxlat(j)
448               IF (nsrf == is_oce) THEN               IF (nsrf == is_oce) THEN
449                  rugmer(i) = yrugm(j)                  rugmer(i) = yrugm(j)
450                  rugos(i, nsrf) = yrugm(j)                  frugs(i, nsrf) = yrugm(j)
451               END IF               END IF
452               agesno(i, nsrf) = yagesno(j)               agesno(i, nsrf) = yagesno(j)
453               fqcalving(i, nsrf) = y_fqcalving(j)               fqcalving(i, nsrf) = y_fqcalving(j)
# Line 625  contains Line 456  contains
456               cdragm(i) = cdragm(i) + coefm(j, 1)               cdragm(i) = cdragm(i) + coefm(j, 1)
457               dflux_t(i) = dflux_t(i) + y_dflux_t(j)               dflux_t(i) = dflux_t(i) + y_dflux_t(j)
458               dflux_q(i) = dflux_q(i) + y_dflux_q(j)               dflux_q(i) = dflux_q(i) + y_dflux_q(j)
              zu1(i) = zu1(i) + yu1(j)  
              zv1(i) = zv1(i) + yv1(j)  
459            END DO            END DO
460            IF (nsrf == is_ter) THEN            IF (nsrf == is_ter) THEN
461               DO j = 1, knon               qsol(ni(:knon)) = yqsol(:knon)
462                  i = ni(j)            else IF (nsrf == is_lic) THEN
                 qsol(i) = yqsol(j)  
              END DO  
           END IF  
           IF (nsrf == is_lic) THEN  
463               DO j = 1, knon               DO j = 1, knon
464                  i = ni(j)                  i = ni(j)
465                  run_off_lic_0(i) = y_run_off_lic_0(j)                  run_off_lic_0(i) = y_run_off_lic_0(j)
466               END DO               END DO
467            END IF            END IF
468            !$$$ PB ajout pour soil  
469            ftsoil(:, :, nsrf) = 0.            ftsoil(:, :, nsrf) = 0.
470            DO k = 1, nsoilmx            ftsoil(ni(:knon), :, nsrf) = ytsoil(:knon, :)
              DO j = 1, knon  
                 i = ni(j)  
                 ftsoil(i, k, nsrf) = ytsoil(j, k)  
              END DO  
           END DO  
471    
472            DO j = 1, knon            DO j = 1, knon
473               i = ni(j)               i = ni(j)
# Line 660  contains Line 480  contains
480               END DO               END DO
481            END DO            END DO
482    
483            !cc diagnostic t, q a 2m et u, v a 10m            ! diagnostic t, q a 2m et u, v a 10m
484    
485            DO j = 1, knon            DO j = 1, knon
486               i = ni(j)               i = ni(j)
487               uzon(j) = yu(j, 1) + y_d_u(j, 1)               u1(j) = yu(j, 1) + y_d_u(j, 1)
488               vmer(j) = yv(j, 1) + y_d_v(j, 1)               v1(j) = yv(j, 1) + y_d_v(j, 1)
489               tair1(j) = yt(j, 1) + y_d_t(j, 1)               tair1(j) = yt(j, 1) + y_d_t(j, 1)
490               qair1(j) = yq(j, 1) + y_d_q(j, 1)               qair1(j) = yq(j, 1) + y_d_q(j, 1)
491               zgeo1(j) = rd*tair1(j)/(0.5*(ypaprs(j, 1)+ypplay(j, &               zgeo1(j) = rd * tair1(j) / (0.5 * (ypaprs(j, 1) + ypplay(j, &
492                    1)))*(ypaprs(j, 1)-ypplay(j, 1))                    1))) * (ypaprs(j, 1)-ypplay(j, 1))
493               tairsol(j) = yts(j) + y_d_ts(j)               tairsol(j) = yts(j) + y_d_ts(j)
494               rugo1(j) = yrugos(j)               rugo1(j) = yrugos(j)
495               IF (nsrf == is_oce) THEN               IF (nsrf == is_oce) THEN
496                  rugo1(j) = rugos(i, nsrf)                  rugo1(j) = frugs(i, nsrf)
497               END IF               END IF
498               psfce(j) = ypaprs(j, 1)               psfce(j) = ypaprs(j, 1)
499               patm(j) = ypplay(j, 1)               patm(j) = ypplay(j, 1)
# Line 681  contains Line 501  contains
501               qairsol(j) = yqsurf(j)               qairsol(j) = yqsurf(j)
502            END DO            END DO
503    
504            CALL stdlevvar(klon, knon, nsrf, zxli, uzon, vmer, tair1, qair1, &            CALL stdlevvar(klon, knon, nsrf, u1(:knon), v1(:knon), tair1(:knon), &
505                 zgeo1, tairsol, qairsol, rugo1, psfce, patm, yt2m, yq2m, &                 qair1, zgeo1, tairsol, qairsol, rugo1, psfce, patm, yt2m, &
506                 yt10m, yq10m, yu10m, yustar)                 yq2m, yt10m, yq10m, wind10m(:knon), ustar)
507    
508            DO j = 1, knon            DO j = 1, knon
509               i = ni(j)               i = ni(j)
510               t2m(i, nsrf) = yt2m(j)               t2m(i, nsrf) = yt2m(j)
511               q2m(i, nsrf) = yq2m(j)               q2m(i, nsrf) = yq2m(j)
512    
513               ! u10m, v10m : composantes du vent a 10m sans spirale de Ekman               u10m_srf(i, nsrf) = (wind10m(j) * u1(j)) &
514               u10m(i, nsrf) = (yu10m(j)*uzon(j))/sqrt(uzon(j)**2+vmer(j)**2)                    / sqrt(u1(j)**2 + v1(j)**2)
515               v10m(i, nsrf) = (yu10m(j)*vmer(j))/sqrt(uzon(j)**2+vmer(j)**2)               v10m_srf(i, nsrf) = (wind10m(j) * v1(j)) &
516                      / sqrt(u1(j)**2 + v1(j)**2)
517            END DO            END DO
518    
519            CALL hbtm(knon, ypaprs, ypplay, yt2m, yt10m, yq2m, yq10m, yustar, &            CALL hbtm(ypaprs, ypplay, yt2m, yq2m, ustar(:knon), y_flux_t(:knon), &
520                 y_flux_t, y_flux_q, yu, yv, yt, yq, ypblh, ycapcl, yoliqcl, &                 y_flux_q(:knon), yu, yv, yt, yq, ypblh(:knon), ycapcl, &
521                 ycteicl, ypblt, ytherm, ytrmb1, ytrmb2, ytrmb3, ylcl)                 yoliqcl, ycteicl, ypblt, ytherm, ytrmb1, ytrmb2, ytrmb3, ylcl)
522    
523            DO j = 1, knon            DO j = 1, knon
524               i = ni(j)               i = ni(j)
# Line 720  contains Line 540  contains
540                  q2(i, k, nsrf) = yq2(j, k)                  q2(i, k, nsrf) = yq2(j, k)
541               END DO               END DO
542            END DO            END DO
543            !IM "slab" ocean         else
544            IF (nsrf == is_oce) THEN            fsnow(:, nsrf) = 0.
              DO j = 1, knon  
                 ! on projette sur la grille globale  
                 i = ni(j)  
                 IF (pctsrf_new(i, is_oce)>epsfra) THEN  
                    flux_o(i) = y_flux_o(j)  
                 ELSE  
                    flux_o(i) = 0.  
                 END IF  
              END DO  
           END IF  
   
           IF (nsrf == is_sic) THEN  
              DO j = 1, knon  
                 i = ni(j)  
                 ! On pondère lorsque l'on fait le bilan au sol :  
                 IF (pctsrf_new(i, is_sic)>epsfra) THEN  
                    flux_g(i) = y_flux_g(j)  
                 ELSE  
                    flux_g(i) = 0.  
                 END IF  
              END DO  
   
           END IF  
           IF (ocean == 'slab  ') THEN  
              IF (nsrf == is_oce) THEN  
                 tslab(1:klon) = ytslab(1:klon)  
                 seaice(1:klon) = y_seaice(1:klon)  
              END IF  
           END IF  
545         end IF if_knon         end IF if_knon
546      END DO loop_surface      END DO loop_surface
547    
548      ! On utilise les nouvelles surfaces      ! On utilise les nouvelles surfaces
549        frugs(:, is_oce) = rugmer
550        pctsrf(:, is_oce) = pctsrf_new_oce
551        pctsrf(:, is_sic) = pctsrf_new_sic
552    
553      rugos(:, is_oce) = rugmer      firstcal = .false.
     pctsrf = pctsrf_new  
554    
555    END SUBROUTINE clmain    END SUBROUTINE clmain
556    
Legend:
Removed from v.82  
changed lines
  Added in v.228
  ViewVC Help
Powered by ViewVC 1.1.21