/[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

revision 61 by guez, Fri Apr 20 14:58:43 2012 UTC revision 62 by guez, Thu Jul 26 14:37:37 2012 UTC
# Line 4  module clmain_m Line 4  module clmain_m
4    
5  contains  contains
6    
7    SUBROUTINE clmain(dtime, itap, date0, pctsrf, pctsrf_new, t, q, u, v,&    SUBROUTINE clmain(dtime, itap, date0, pctsrf, pctsrf_new, t, q, u, v, &
8         jour, rmu0, co2_ppm, ok_veget, ocean, npas, nexca, ts,&         jour, rmu0, co2_ppm, ok_veget, ocean, npas, nexca, ts, &
9         soil_model, cdmmax, cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil,&         soil_model, cdmmax, cdhmax, ksta, ksta_ter, ok_kzmin, ftsoil, &
10         qsol, paprs, pplay, snow, qsurf, evap, albe, alblw, fluxlat,&         qsol, paprs, pplay, snow, qsurf, evap, albe, alblw, fluxlat, &
11         rain_f, snow_f, solsw, sollw, sollwdown, fder, rlon, rlat, cufi,&         rain_fall, snow_f, solsw, sollw, sollwdown, fder, rlon, rlat, cufi, &
12         cvfi, rugos, debut, lafin, agesno, rugoro, d_t, d_q, d_u, d_v,&         cvfi, rugos, debut, lafin, agesno, rugoro, d_t, d_q, d_u, d_v, &
13         d_ts, flux_t, flux_q, flux_u, flux_v, cdragh, cdragm, q2,&         d_ts, flux_t, flux_q, flux_u, flux_v, cdragh, cdragm, q2, &
14         dflux_t, dflux_q, zcoefh, zu1, zv1, t2m, q2m, u10m, v10m, pblh,&         dflux_t, dflux_q, zcoefh, zu1, zv1, t2m, q2m, u10m, v10m, pblh, &
15         capcl, oliqcl, cteicl, pblt, therm, trmb1, trmb2, trmb3, plcl,&         capcl, oliqcl, cteicl, pblt, therm, trmb1, trmb2, trmb3, plcl, &
16         fqcalving, ffonte, run_off_lic_0, flux_o, flux_g, tslab, seaice)         fqcalving, ffonte, run_off_lic_0, flux_o, flux_g, tslab, seaice)
17    
18      ! 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
19      ! Author: Z.X. Li (LMD/CNRS), date: 1993/08/18      ! Author: Z. X. Li (LMD/CNRS), date: 1993/08/18
20      ! Objet : interface de "couche limite" (diffusion verticale)      ! Objet : interface de couche limite (diffusion verticale)
21    
22      ! Tout ce qui a trait aux traceurs est dans "phytrac" maintenant.      ! Tout ce qui a trait aux traceurs est dans "phytrac". Le calcul
23      ! Pour l'instant le calcul de la couche limite pour les traceurs      ! de la couche limite pour les traceurs se fait avec "cltrac" et
24      ! se fait avec "cltrac" et ne tient pas compte de la différentiation      ! ne tient pas compte de la différentiation des sous-fractions de
25      ! des sous-fractions de sol.      ! sol.
26    
27      ! Pour pouvoir extraire les coefficients d'échanges et le vent      ! Pour pouvoir extraire les coefficients d'échanges et le vent
28      ! dans la première couche, trois champs supplémentaires ont été      ! dans la première couche, trois champs ont été créés : "zcoefh",
29      ! créés : "zcoefh", "zu1" et "zv1". Pour l'instant nous avons      ! "zu1" et "zv1". Nous avons moyenné les valeurs de ces trois
30      ! moyenné les valeurs de ces trois champs sur les 4 sous-surfaces      ! champs sur les quatre sous-surfaces du modèle.
     ! du modèle. Dans l'avenir, si les informations des sous-surfaces  
     ! doivent être prises en compte, il faudra sortir ces mêmes champs  
     ! en leur ajoutant une dimension, c'est-à-dire "nbsrf" (nombre de  
     ! sous-surfaces).  
31    
32      use calendar, ONLY : ymds2ju      use calendar, ONLY: ymds2ju
33      use clqh_m, only: clqh      use clqh_m, only: clqh
34        use clvent_m, only: clvent
35      use coefkz_m, only: coefkz      use coefkz_m, only: coefkz
36      use coefkzmin_m, only: coefkzmin      use coefkzmin_m, only: coefkzmin
37      USE conf_phys_m, ONLY : iflag_pbl      USE conf_gcm_m, ONLY: prt_level
38      USE dimens_m, ONLY : iim, jjm      USE conf_phys_m, ONLY: iflag_pbl
39      USE dimphy, ONLY : klev, klon, zmasq      USE dimens_m, ONLY: iim, jjm
40      USE dimsoil, ONLY : nsoilmx      USE dimphy, ONLY: klev, klon, zmasq
41      USE dynetat0_m, ONLY : day_ini      USE dimsoil, ONLY: nsoilmx
42      USE gath_cpl, ONLY : gath2cpl      USE dynetat0_m, ONLY: day_ini
43        USE gath_cpl, ONLY: gath2cpl
44      use hbtm_m, only: hbtm      use hbtm_m, only: hbtm
45      USE histsync_m, ONLY : histsync      USE histbeg_totreg_m, ONLY: histbeg_totreg
46      USE histbeg_totreg_m, ONLY : histbeg_totreg      USE histdef_m, ONLY: histdef
47      USE histend_m, ONLY : histend      USE histend_m, ONLY: histend
48      USE histdef_m, ONLY : histdef      USE histsync_m, ONLY: histsync
49      use histwrite_m, only: histwrite      use histwrite_m, only: histwrite
50      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
51      USE conf_gcm_m, ONLY : prt_level      USE suphec_m, ONLY: rd, rg, rkappa
52      USE suphec_m, ONLY : rd, rg, rkappa      USE temps, ONLY: annee_ref, itau_phy
53      USE temps, ONLY : annee_ref, itau_phy      use ustarhb_m, only: ustarhb
54        use vdif_kcay_m, only: vdif_kcay
55      use yamada4_m, only: yamada4      use yamada4_m, only: yamada4
56    
57      ! Arguments:      ! Arguments:
58    
59      REAL, INTENT (IN) :: dtime ! interval du temps (secondes)      REAL, INTENT(IN):: dtime ! interval du temps (secondes)
60      REAL date0      INTEGER, INTENT(IN):: itap ! numero du pas de temps
61      ! date0----input-R- jour initial      REAL, INTENT(IN):: date0 ! jour initial
62      INTEGER, INTENT (IN) :: itap      REAL, INTENT(inout):: pctsrf(klon, nbsrf)
63      ! itap-----input-I- numero du pas de temps  
64      REAL, INTENT(IN):: t(klon, klev), q(klon, klev)      ! la nouvelle repartition des surfaces sortie de l'interface
65      ! t--------input-R- temperature (K)      REAL, INTENT(out):: pctsrf_new(klon, nbsrf)
66      ! q--------input-R- vapeur d'eau (kg/kg)  
67      REAL, INTENT (IN):: u(klon, klev), v(klon, klev)      REAL, INTENT(IN):: t(klon, klev) ! temperature (K)
68      ! u--------input-R- vitesse u      REAL, INTENT(IN):: q(klon, klev) ! vapeur d'eau (kg/kg)
69      ! v--------input-R- vitesse v      REAL, INTENT(IN):: u(klon, klev), v(klon, klev) ! vitesse
70      REAL, INTENT (IN):: paprs(klon, klev+1)      INTEGER, INTENT(IN):: jour ! jour de l'annee en cours
71      ! paprs----input-R- pression a intercouche (Pa)      REAL, intent(in):: rmu0(klon) ! cosinus de l'angle solaire zenithal    
72      REAL, INTENT (IN):: pplay(klon, klev)      REAL, INTENT(IN):: paprs(klon, klev+1) ! pression a intercouche (Pa)
73      ! pplay----input-R- pression au milieu de couche (Pa)      REAL, INTENT(IN):: pplay(klon, klev) ! pression au milieu de couche (Pa)
74      REAL, INTENT (IN):: rlon(klon), rlat(klon)      REAL, INTENT(IN):: rlon(klon)
75      ! rlat-----input-R- latitude en degree      REAL, INTENT(IN):: rlat(klon) ! latitude en degrés
76      REAL cufi(klon), cvfi(klon)      REAL cufi(klon), cvfi(klon)
77      ! cufi-----input-R- resolution des mailles en x (m)      ! cufi-----input-R- resolution des mailles en x (m)
78      ! cvfi-----input-R- resolution des mailles en y (m)      ! cvfi-----input-R- resolution des mailles en y (m)
79      REAL d_t(klon, klev), d_q(klon, klev)      REAL d_t(klon, klev), d_q(klon, klev)
80      ! d_t------output-R- le changement pour "t"      ! d_t------output-R- le changement pour "t"
81      ! d_q------output-R- le changement pour "q"      ! d_q------output-R- le changement pour "q"
82      REAL d_u(klon, klev), d_v(klon, klev)  
83      ! d_u------output-R- le changement pour "u"      REAL, intent(out):: d_u(klon, klev), d_v(klon, klev)
84      ! d_v------output-R- le changement pour "v"      ! changement pour "u" et "v"
85    
86      REAL flux_t(klon, klev, nbsrf), flux_q(klon, klev, nbsrf)      REAL flux_t(klon, klev, nbsrf), flux_q(klon, klev, nbsrf)
87      ! flux_t---output-R- flux de chaleur sensible (CpT) J/m**2/s (W/m**2)      ! flux_t---output-R- flux de chaleur sensible (CpT) J/m**2/s (W/m**2)
88      !                    (orientation positive vers le bas)      !                    (orientation positive vers le bas)
# Line 112  contains Line 112  contains
112      ! flux_u---output-R- tension du vent X: (kg m/s)/(m**2 s) ou Pascal      ! flux_u---output-R- tension du vent X: (kg m/s)/(m**2 s) ou Pascal
113      ! flux_v---output-R- tension du vent Y: (kg m/s)/(m**2 s) ou Pascal      ! flux_v---output-R- tension du vent Y: (kg m/s)/(m**2 s) ou Pascal
114      REAL rugmer(klon), agesno(klon, nbsrf)      REAL rugmer(klon), agesno(klon, nbsrf)
115      REAL, INTENT (IN) :: rugoro(klon)      REAL, INTENT(IN):: rugoro(klon)
116      REAL cdragh(klon), cdragm(klon)      REAL, INTENT(out):: cdragh(klon), cdragm(klon)
     ! jour de l'annee en cours                  
     INTEGER jour  
     REAL rmu0(klon) ! cosinus de l'angle solaire zenithal      
117      ! taux CO2 atmosphere                          ! taux CO2 atmosphere                    
118      REAL co2_ppm      REAL co2_ppm
119      LOGICAL, INTENT (IN) :: debut      LOGICAL, INTENT(IN):: debut
120      LOGICAL, INTENT (IN) :: lafin      LOGICAL, INTENT(IN):: lafin
121      LOGICAL ok_veget      LOGICAL ok_veget
122      CHARACTER (len=*), INTENT (IN) :: ocean      CHARACTER(len=*), INTENT(IN):: ocean
123      INTEGER npas, nexca      INTEGER npas, nexca
124    
     REAL pctsrf(klon, nbsrf)  
125      REAL ts(klon, nbsrf)      REAL ts(klon, nbsrf)
126      ! ts-------input-R- temperature du sol (en Kelvin)      ! ts-------input-R- temperature du sol (en Kelvin)
127      REAL d_ts(klon, nbsrf)      REAL d_ts(klon, nbsrf)
# Line 138  contains Line 134  contains
134    
135      REAL fluxlat(klon, nbsrf)      REAL fluxlat(klon, nbsrf)
136    
137      REAL rain_f(klon), snow_f(klon)      REAL, intent(in):: rain_fall(klon), snow_f(klon)
138      REAL fder(klon)      REAL fder(klon)
139    
140      REAL sollw(klon, nbsrf), solsw(klon, nbsrf), sollwdown(klon)      REAL sollw(klon, nbsrf), solsw(klon, nbsrf), sollwdown(klon)
141      REAL rugos(klon, nbsrf)      REAL rugos(klon, nbsrf)
142      ! rugos----input-R- longeur de rugosite (en m)      ! rugos----input-R- longeur de rugosite (en m)
     ! la nouvelle repartition des surfaces sortie de l'interface  
     REAL pctsrf_new(klon, nbsrf)  
143    
144      REAL zcoefh(klon, klev)      REAL zcoefh(klon, klev)
145      REAL zu1(klon)      REAL zu1(klon)
146      REAL zv1(klon)      REAL zv1(klon)
147    
148      !$$$ PB ajout pour soil      !$$$ PB ajout pour soil
149      LOGICAL, INTENT (IN) :: soil_model      LOGICAL, INTENT(IN):: soil_model
150      !IM ajout seuils cdrm, cdrh      !IM ajout seuils cdrm, cdrh
151      REAL cdmmax, cdhmax      REAL cdmmax, cdhmax
152    
# Line 163  contains Line 157  contains
157      REAL ytsoil(klon, nsoilmx)      REAL ytsoil(klon, nsoilmx)
158      REAL qsol(klon)      REAL qsol(klon)
159    
     EXTERNAL clvent, calbeta, cltrac  
   
160      REAL yts(klon), yrugos(klon), ypct(klon), yz0_new(klon)      REAL yts(klon), yrugos(klon), ypct(klon), yz0_new(klon)
161      REAL yalb(klon)      REAL yalb(klon)
162      REAL yalblw(klon)      REAL yalblw(klon)
# Line 185  contains Line 177  contains
177      REAL y_flux_t(klon, klev), y_flux_q(klon, klev)      REAL y_flux_t(klon, klev), y_flux_q(klon, klev)
178      REAL y_flux_u(klon, klev), y_flux_v(klon, klev)      REAL y_flux_u(klon, klev), y_flux_v(klon, klev)
179      REAL y_dflux_t(klon), y_dflux_q(klon)      REAL y_dflux_t(klon), y_dflux_q(klon)
180      REAL ycoefh(klon, klev), ycoefm(klon, klev)      REAL coefh(klon, klev), coefm(klon, klev)
181      REAL yu(klon, klev), yv(klon, klev)      REAL yu(klon, klev), yv(klon, klev)
182      REAL yt(klon, klev), yq(klon, klev)      REAL yt(klon, klev), yq(klon, klev)
183      REAL ypaprs(klon, klev+1), ypplay(klon, klev), ydelp(klon, klev)      REAL ypaprs(klon, klev+1), ypplay(klon, klev), ydelp(klon, klev)
# Line 214  contains Line 206  contains
206    
207      ! maf pour sorties IOISPL en cas de debugagage      ! maf pour sorties IOISPL en cas de debugagage
208    
209      CHARACTER (80) cldebug      CHARACTER(80) cldebug
210      SAVE cldebug      SAVE cldebug
211      CHARACTER (8) cl_surf(nbsrf)      CHARACTER(8) cl_surf(nbsrf)
212      SAVE cl_surf      SAVE cl_surf
213      INTEGER nhoridbg, nidbg      INTEGER nhoridbg, nidbg
214      SAVE nhoridbg, nidbg      SAVE nhoridbg, nidbg
# Line 227  contains Line 219  contains
219      LOGICAL first_appel      LOGICAL first_appel
220      SAVE first_appel      SAVE first_appel
221      DATA first_appel/ .TRUE./      DATA first_appel/ .TRUE./
222      LOGICAL :: debugindex = .FALSE.      LOGICAL:: debugindex = .FALSE.
223      INTEGER idayref      INTEGER idayref
224      REAL t2m(klon, nbsrf), q2m(klon, nbsrf)      REAL t2m(klon, nbsrf), q2m(klon, nbsrf)
225      REAL u10m(klon, nbsrf), v10m(klon, nbsrf)      REAL u10m(klon, nbsrf), v10m(klon, nbsrf)
# Line 268  contains Line 260  contains
260      REAL ytrmb1(klon)      REAL ytrmb1(klon)
261      REAL ytrmb2(klon)      REAL ytrmb2(klon)
262      REAL ytrmb3(klon)      REAL ytrmb3(klon)
     REAL y_cd_h(klon), y_cd_m(klon)  
263      REAL uzon(klon), vmer(klon)      REAL uzon(klon), vmer(klon)
264      REAL tair1(klon), qair1(klon), tairsol(klon)      REAL tair1(klon), qair1(klon), tairsol(klon)
265      REAL psfce(klon), patm(klon)      REAL psfce(klon), patm(klon)
# Line 284  contains Line 275  contains
275      REAL t_coup      REAL t_coup
276      PARAMETER (t_coup=273.15)      PARAMETER (t_coup=273.15)
277    
278      CHARACTER (len=20) :: modname = 'clmain'      CHARACTER(len=20):: modname = 'clmain'
279    
280      !------------------------------------------------------------      !------------------------------------------------------------
281    
# Line 425  contains Line 416  contains
416            CALL histwrite(nidbg, cl_surf(nsrf), itap, debugtab)            CALL histwrite(nidbg, cl_surf(nsrf), itap, debugtab)
417         END IF         END IF
418    
419         IF (knon == 0) CYCLE         if_knon: IF (knon /= 0) then
   
        DO j = 1, knon  
           i = ni(j)  
           ypct(j) = pctsrf(i, nsrf)  
           yts(j) = ts(i, nsrf)  
           ytslab(i) = tslab(i)  
           ysnow(j) = snow(i, nsrf)  
           yqsurf(j) = qsurf(i, nsrf)  
           yalb(j) = albe(i, nsrf)  
           yalblw(j) = alblw(i, nsrf)  
           yrain_f(j) = rain_f(i)  
           ysnow_f(j) = snow_f(i)  
           yagesno(j) = agesno(i, nsrf)  
           yfder(j) = fder(i)  
           ytaux(j) = flux_u(i, 1, nsrf)  
           ytauy(j) = flux_v(i, 1, nsrf)  
           ysolsw(j) = solsw(i, nsrf)  
           ysollw(j) = sollw(i, nsrf)  
           ysollwdown(j) = sollwdown(i)  
           yrugos(j) = rugos(i, nsrf)  
           yrugoro(j) = rugoro(i)  
           yu1(j) = u1lay(i)  
           yv1(j) = v1lay(i)  
           yrads(j) = ysolsw(j) + ysollw(j)  
           ypaprs(j, klev+1) = paprs(i, klev+1)  
           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))  
        END DO  
   
        ! IF bucket model for continent, copy soil water content  
        IF (nsrf == is_ter .AND. .NOT. ok_veget) THEN  
           DO j = 1, knon  
              i = ni(j)  
              yqsol(j) = qsol(i)  
           END DO  
        ELSE  
           yqsol = 0.  
        END IF  
        !$$$ PB ajour pour soil  
        DO k = 1, nsoilmx  
           DO j = 1, knon  
              i = ni(j)  
              ytsoil(j, k) = ftsoil(i, k, nsrf)  
           END DO  
        END DO  
        DO k = 1, klev  
420            DO j = 1, knon            DO j = 1, knon
421               i = ni(j)               i = ni(j)
422               ypaprs(j, k) = paprs(i, k)               ypct(j) = pctsrf(i, nsrf)
423               ypplay(j, k) = pplay(i, k)               yts(j) = ts(i, nsrf)
424               ydelp(j, k) = delp(i, k)               ytslab(i) = tslab(i)
425               yu(j, k) = u(i, k)               ysnow(j) = snow(i, nsrf)
426               yv(j, k) = v(i, k)               yqsurf(j) = qsurf(i, nsrf)
427               yt(j, k) = t(i, k)               yalb(j) = albe(i, nsrf)
428               yq(j, k) = q(i, k)               yalblw(j) = alblw(i, nsrf)
429                 yrain_f(j) = rain_fall(i)
430                 ysnow_f(j) = snow_f(i)
431                 yagesno(j) = agesno(i, nsrf)
432                 yfder(j) = fder(i)
433                 ytaux(j) = flux_u(i, 1, nsrf)
434                 ytauy(j) = flux_v(i, 1, nsrf)
435                 ysolsw(j) = solsw(i, nsrf)
436                 ysollw(j) = sollw(i, nsrf)
437                 ysollwdown(j) = sollwdown(i)
438                 yrugos(j) = rugos(i, nsrf)
439                 yrugoro(j) = rugoro(i)
440                 yu1(j) = u1lay(i)
441                 yv1(j) = v1lay(i)
442                 yrads(j) = ysolsw(j) + ysollw(j)
443                 ypaprs(j, klev+1) = paprs(i, klev+1)
444                 y_run_off_lic_0(j) = run_off_lic_0(i)
445                 yu10mx(j) = u10m(i, nsrf)
446                 yu10my(j) = v10m(i, nsrf)
447                 ywindsp(j) = sqrt(yu10mx(j)*yu10mx(j)+yu10my(j)*yu10my(j))
448            END DO            END DO
        END DO  
449    
450         ! calculer Cdrag et les coefficients d'echange            ! IF bucket model for continent, copy soil water content
451         CALL coefkz(nsrf, knon, ypaprs, ypplay, ksta, ksta_ter, yts,&            IF (nsrf == is_ter .AND. .NOT. ok_veget) THEN
452              yrugos, yu, yv, yt, yq, yqsurf, ycoefm, ycoefh)               DO j = 1, knon
453         IF (iflag_pbl == 1) THEN                  i = ni(j)
454            CALL coefkz2(nsrf, knon, ypaprs, ypplay, yt, ycoefm0, ycoefh0)                  yqsol(j) = qsol(i)
           DO k = 1, klev  
              DO i = 1, knon  
                 ycoefm(i, k) = max(ycoefm(i, k), ycoefm0(i, k))  
                 ycoefh(i, k) = max(ycoefh(i, k), ycoefh0(i, k))  
455               END DO               END DO
456            END DO            ELSE
457         END IF               yqsol = 0.
458              END IF
        ! on seuille ycoefm et ycoefh  
        IF (nsrf == is_oce) THEN  
           DO j = 1, knon  
              ycoefm(j, 1) = min(ycoefm(j, 1), cdmmax)  
              ycoefh(j, 1) = min(ycoefh(j, 1), cdhmax)  
           END DO  
        END IF  
   
        IF (ok_kzmin) THEN  
           ! Calcul d'une diffusion minimale pour les conditions tres stables  
           CALL coefkzmin(knon, ypaprs, ypplay, yu, yv, yt, yq, ycoefm(:, 1), &  
                ycoefm0, ycoefh0)  
459    
460            DO k = 1, klev            DO k = 1, nsoilmx
461               DO i = 1, knon               DO j = 1, knon
462                  ycoefm(i, k) = max(ycoefm(i, k), ycoefm0(i, k))                  i = ni(j)
463                  ycoefh(i, k) = max(ycoefh(i, k), ycoefh0(i, k))                  ytsoil(j, k) = ftsoil(i, k, nsrf)
464               END DO               END DO
465            END DO            END DO
        END IF  
466    
        IF (iflag_pbl >= 3) THEN  
           ! MELLOR ET YAMADA adapté à Mars, Richard Fournier et Frédéric Hourdin  
           yzlay(1:knon, 1) = rd*yt(1:knon, 1)/(0.5*(ypaprs(1:knon, &  
                1)+ypplay(1:knon, 1)))*(ypaprs(1:knon, 1)-ypplay(1:knon, 1))/rg  
           DO k = 2, klev  
              yzlay(1:knon, k) = yzlay(1:knon, k-1) &  
                   + rd * 0.5 * (yt(1:knon, k-1) + yt(1:knon, k)) &  
                   / ypaprs(1:knon, k) &  
                   * (ypplay(1:knon, k-1) - ypplay(1:knon, k)) / rg  
           END DO  
467            DO k = 1, klev            DO k = 1, klev
              yteta(1:knon, k) = yt(1:knon, k)*(ypaprs(1:knon, 1) &  
                   / ypplay(1:knon, k))**rkappa * (1.+0.61*yq(1:knon, k))  
           END DO  
           yzlev(1:knon, 1) = 0.  
           yzlev(1:knon, klev+1) = 2.*yzlay(1:knon, klev) - yzlay(1:knon, klev-1)  
           DO k = 2, klev  
              yzlev(1:knon, k) = 0.5*(yzlay(1:knon, k)+yzlay(1:knon, k-1))  
           END DO  
           DO k = 1, klev + 1  
468               DO j = 1, knon               DO j = 1, knon
469                  i = ni(j)                  i = ni(j)
470                  yq2(j, k) = q2(i, k, nsrf)                  ypaprs(j, k) = paprs(i, k)
471                    ypplay(j, k) = pplay(i, k)
472                    ydelp(j, k) = delp(i, k)
473                    yu(j, k) = u(i, k)
474                    yv(j, k) = v(i, k)
475                    yt(j, k) = t(i, k)
476                    yq(j, k) = q(i, k)
477               END DO               END DO
478            END DO            END DO
479    
480            y_cd_m(1:knon) = ycoefm(1:knon, 1)            ! calculer Cdrag et les coefficients d'echange
481            y_cd_h(1:knon) = ycoefh(1:knon, 1)            CALL coefkz(nsrf, knon, ypaprs, ypplay, ksta, ksta_ter, yts, yrugos, &
482            CALL ustarhb(knon, yu, yv, y_cd_m, yustar)                 yu, yv, yt, yq, yqsurf, coefm(:knon, :), coefh(:knon, :))
483              IF (iflag_pbl == 1) THEN
484                 CALL coefkz2(nsrf, knon, ypaprs, ypplay, yt, ycoefm0, ycoefh0)
485                 coefm(:knon, :) = max(coefm(:knon, :), ycoefm0(:knon, :))
486                 coefh(:knon, :) = max(coefh(:knon, :), ycoefh0(:knon, :))
487              END IF
488    
489            IF (prt_level>9) THEN            ! on seuille coefm et coefh
490               PRINT *, 'USTAR = ', yustar            IF (nsrf == is_oce) THEN
491                 coefm(:knon, 1) = min(coefm(:knon, 1), cdmmax)
492                 coefh(:knon, 1) = min(coefh(:knon, 1), cdhmax)
493            END IF            END IF
494    
495            ! iflag_pbl peut être utilisé comme longueur de mélange            IF (ok_kzmin) THEN
496                 ! Calcul d'une diffusion minimale pour les conditions tres stables
497                 CALL coefkzmin(knon, ypaprs, ypplay, yu, yv, yt, yq, &
498                      coefm(:, 1), ycoefm0, ycoefh0)
499                 coefm(:knon, :) = max(coefm(:knon, :), ycoefm0(:knon, :))
500                 coefh(:knon, :) = max(coefh(:knon, :), ycoefh0(:knon, :))
501               END IF
502    
503              IF (iflag_pbl >= 3) THEN
504                 ! MELLOR ET YAMADA adapté à Mars, Richard Fournier et
505                 ! Frédéric Hourdin
506                 yzlay(:knon, 1) = rd * yt(:knon, 1) / (0.5 * (ypaprs(:knon, 1) &
507                      + ypplay(:knon, 1))) &
508                      * (ypaprs(:knon, 1) - ypplay(:knon, 1)) / rg
509                 DO k = 2, klev
510                    yzlay(1:knon, k) = yzlay(1:knon, k-1) &
511                         + rd * 0.5 * (yt(1:knon, k-1) + yt(1:knon, k)) &
512                         / ypaprs(1:knon, k) &
513                         * (ypplay(1:knon, k-1) - ypplay(1:knon, k)) / rg
514                 END DO
515                 DO k = 1, klev
516                    yteta(1:knon, k) = yt(1:knon, k)*(ypaprs(1:knon, 1) &
517                         / ypplay(1:knon, k))**rkappa * (1.+0.61*yq(1:knon, k))
518                 END DO
519                 yzlev(1:knon, 1) = 0.
520                 yzlev(:knon, klev+1) = 2. * yzlay(:knon, klev) &
521                      - yzlay(:knon, klev - 1)
522                 DO k = 2, klev
523                    yzlev(1:knon, k) = 0.5*(yzlay(1:knon, k)+yzlay(1:knon, k-1))
524                 END DO
525                 DO k = 1, klev + 1
526                    DO j = 1, knon
527                       i = ni(j)
528                       yq2(j, k) = q2(i, k, nsrf)
529                    END DO
530                 END DO
531    
532            IF (iflag_pbl >= 11) THEN               CALL ustarhb(knon, yu, yv, coefm(:knon, 1), yustar)
533               CALL vdif_kcay(knon, dtime, rg, rd, ypaprs, yt, yzlev, yzlay, &  
534                    yu, yv, yteta, y_cd_m, yq2, q2diag, ykmm, ykmn, yustar, &               IF (prt_level > 9) THEN
535                    iflag_pbl)                  PRINT *, 'USTAR = ', yustar
536            ELSE               END IF
537               CALL yamada4(knon, dtime, rg, yzlev, yzlay, yu, yv, yteta, &  
538                    y_cd_m, yq2, ykmm, ykmn, ykmq, yustar, iflag_pbl)               ! iflag_pbl peut être utilisé comme longueur de mélange
539    
540                 IF (iflag_pbl >= 11) THEN
541                    CALL vdif_kcay(knon, dtime, rg, rd, ypaprs, yt, yzlev, yzlay, &
542                         yu, yv, yteta, coefm(:knon, 1), yq2, q2diag, ykmm, ykmn, &
543                         yustar, iflag_pbl)
544                 ELSE
545                    CALL yamada4(knon, dtime, rg, yzlev, yzlay, yu, yv, yteta, &
546                         coefm(:knon, 1), yq2, ykmm, ykmn, ykmq, yustar, iflag_pbl)
547                 END IF
548    
549                 coefm(:knon, 2:) = ykmm(:knon, 2:klev)
550                 coefh(:knon, 2:) = ykmn(:knon, 2:klev)
551            END IF            END IF
552    
553            ycoefm(1:knon, 1) = y_cd_m(1:knon)            ! calculer la diffusion des vitesses "u" et "v"
554            ycoefh(1:knon, 1) = y_cd_h(1:knon)            CALL clvent(knon, dtime, yu1, yv1, coefm, yt, yu, ypaprs, ypplay, &
555            ycoefm(1:knon, 2:klev) = ykmm(1:knon, 2:klev)                 ydelp, y_d_u, y_flux_u)
556            ycoefh(1:knon, 2:klev) = ykmn(1:knon, 2:klev)            CALL clvent(knon, dtime, yu1, yv1, coefm, yt, yv, ypaprs, ypplay, &
557         END IF                 ydelp, y_d_v, y_flux_v)
558    
559              ! pour le couplage
560              ytaux = y_flux_u(:, 1)
561              ytauy = y_flux_v(:, 1)
562    
563              ! calculer la diffusion de "q" et de "h"
564              CALL clqh(dtime, itap, date0, jour, debut, lafin, rlon, rlat, &
565                   cufi, cvfi, knon, nsrf, ni, pctsrf, soil_model, ytsoil, &
566                   yqsol, ok_veget, ocean, npas, nexca, rmu0, co2_ppm, yrugos, &
567                   yrugoro, yu1, yv1, coefh, yt, yq, yts, ypaprs, ypplay, &
568                   ydelp, yrads, yalb, yalblw, ysnow, yqsurf, yrain_f, ysnow_f, &
569                   yfder, ytaux, ytauy, ywindsp, ysollw, ysollwdown, ysolsw, &
570                   yfluxlat, pctsrf_new, yagesno, y_d_t, y_d_q, y_d_ts, &
571                   yz0_new, y_flux_t, y_flux_q, y_dflux_t, y_dflux_q, &
572                   y_fqcalving, y_ffonte, y_run_off_lic_0, y_flux_o, y_flux_g, &
573                   ytslab, y_seaice)
574    
575         ! calculer la diffusion des vitesses "u" et "v"            ! calculer la longueur de rugosite sur ocean
576         CALL clvent(knon, dtime, yu1, yv1, ycoefm, yt, yu, ypaprs, ypplay, &            yrugm = 0.
577              ydelp, y_d_u, y_flux_u)            IF (nsrf == is_oce) THEN
578         CALL clvent(knon, dtime, yu1, yv1, ycoefm, yt, yv, ypaprs, ypplay, &               DO j = 1, knon
579              ydelp, y_d_v, y_flux_v)                  yrugm(j) = 0.018*coefm(j, 1)*(yu1(j)**2+yv1(j)**2)/rg + &
580                         0.11*14E-6/sqrt(coefm(j, 1)*(yu1(j)**2+yv1(j)**2))
581         ! pour le couplage                  yrugm(j) = max(1.5E-05, yrugm(j))
582         ytaux = y_flux_u(:, 1)               END DO
583         ytauy = y_flux_v(:, 1)            END IF
   
        ! calculer la diffusion de "q" et de "h"  
        CALL clqh(dtime, itap, date0, jour, debut, lafin, rlon, rlat,&  
             cufi, cvfi, knon, nsrf, ni, pctsrf, soil_model, ytsoil,&  
             yqsol, ok_veget, ocean, npas, nexca, rmu0, co2_ppm, yrugos,&  
             yrugoro, yu1, yv1, ycoefh, yt, yq, yts, ypaprs, ypplay,&  
             ydelp, yrads, yalb, yalblw, ysnow, yqsurf, yrain_f, ysnow_f, &  
             yfder, ytaux, ytauy, ywindsp, ysollw, ysollwdown, ysolsw,&  
             yfluxlat, pctsrf_new, yagesno, y_d_t, y_d_q, y_d_ts,&  
             yz0_new, y_flux_t, y_flux_q, y_dflux_t, y_dflux_q,&  
             y_fqcalving, y_ffonte, y_run_off_lic_0, y_flux_o, y_flux_g,&  
             ytslab, y_seaice)  
   
        ! calculer la longueur de rugosite sur ocean  
        yrugm = 0.  
        IF (nsrf == is_oce) THEN  
584            DO j = 1, knon            DO j = 1, knon
585               yrugm(j) = 0.018*ycoefm(j, 1)*(yu1(j)**2+yv1(j)**2)/rg + &               y_dflux_t(j) = y_dflux_t(j)*ypct(j)
586                    0.11*14E-6/sqrt(ycoefm(j, 1)*(yu1(j)**2+yv1(j)**2))               y_dflux_q(j) = y_dflux_q(j)*ypct(j)
587               yrugm(j) = max(1.5E-05, yrugm(j))               yu1(j) = yu1(j)*ypct(j)
588                 yv1(j) = yv1(j)*ypct(j)
589            END DO            END DO
        END IF  
        DO j = 1, knon  
           y_dflux_t(j) = y_dflux_t(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)  
        END DO  
590    
591         DO k = 1, klev            DO k = 1, klev
592            DO j = 1, knon               DO j = 1, knon
593               i = ni(j)                  i = ni(j)
594               ycoefh(j, k) = ycoefh(j, k)*ypct(j)                  coefh(j, k) = coefh(j, k)*ypct(j)
595               ycoefm(j, k) = ycoefm(j, k)*ypct(j)                  coefm(j, k) = coefm(j, k)*ypct(j)
596               y_d_t(j, k) = y_d_t(j, k)*ypct(j)                  y_d_t(j, k) = y_d_t(j, k)*ypct(j)
597               y_d_q(j, k) = y_d_q(j, k)*ypct(j)                  y_d_q(j, k) = y_d_q(j, k)*ypct(j)
598               flux_t(i, k, nsrf) = y_flux_t(j, k)                  flux_t(i, k, nsrf) = y_flux_t(j, k)
599               flux_q(i, k, nsrf) = y_flux_q(j, k)                  flux_q(i, k, nsrf) = y_flux_q(j, k)
600               flux_u(i, k, nsrf) = y_flux_u(j, k)                  flux_u(i, k, nsrf) = y_flux_u(j, k)
601               flux_v(i, k, nsrf) = y_flux_v(j, k)                  flux_v(i, k, nsrf) = y_flux_v(j, k)
602               y_d_u(j, k) = y_d_u(j, k)*ypct(j)                  y_d_u(j, k) = y_d_u(j, k)*ypct(j)
603               y_d_v(j, k) = y_d_v(j, k)*ypct(j)                  y_d_v(j, k) = y_d_v(j, k)*ypct(j)
604                 END DO
605            END DO            END DO
        END DO  
606    
607         evap(:, nsrf) = -flux_q(:, 1, nsrf)            evap(:, nsrf) = -flux_q(:, 1, nsrf)
608    
609         albe(:, nsrf) = 0.            albe(:, nsrf) = 0.
610         alblw(:, nsrf) = 0.            alblw(:, nsrf) = 0.
611         snow(:, nsrf) = 0.            snow(:, nsrf) = 0.
612         qsurf(:, nsrf) = 0.            qsurf(:, nsrf) = 0.
613         rugos(:, nsrf) = 0.            rugos(:, nsrf) = 0.
614         fluxlat(:, nsrf) = 0.            fluxlat(:, nsrf) = 0.
        DO j = 1, knon  
           i = ni(j)  
           d_ts(i, nsrf) = y_d_ts(j)  
           albe(i, nsrf) = yalb(j)  
           alblw(i, nsrf) = yalblw(j)  
           snow(i, nsrf) = ysnow(j)  
           qsurf(i, nsrf) = yqsurf(j)  
           rugos(i, nsrf) = yz0_new(j)  
           fluxlat(i, nsrf) = yfluxlat(j)  
           IF (nsrf == is_oce) THEN  
              rugmer(i) = yrugm(j)  
              rugos(i, nsrf) = yrugm(j)  
           END IF  
           agesno(i, nsrf) = yagesno(j)  
           fqcalving(i, nsrf) = y_fqcalving(j)  
           ffonte(i, nsrf) = y_ffonte(j)  
           cdragh(i) = cdragh(i) + ycoefh(j, 1)  
           cdragm(i) = cdragm(i) + ycoefm(j, 1)  
           dflux_t(i) = dflux_t(i) + y_dflux_t(j)  
           dflux_q(i) = dflux_q(i) + y_dflux_q(j)  
           zu1(i) = zu1(i) + yu1(j)  
           zv1(i) = zv1(i) + yv1(j)  
        END DO  
        IF (nsrf == is_ter) THEN  
615            DO j = 1, knon            DO j = 1, knon
616               i = ni(j)               i = ni(j)
617               qsol(i) = yqsol(j)               d_ts(i, nsrf) = y_d_ts(j)
618                 albe(i, nsrf) = yalb(j)
619                 alblw(i, nsrf) = yalblw(j)
620                 snow(i, nsrf) = ysnow(j)
621                 qsurf(i, nsrf) = yqsurf(j)
622                 rugos(i, nsrf) = yz0_new(j)
623                 fluxlat(i, nsrf) = yfluxlat(j)
624                 IF (nsrf == is_oce) THEN
625                    rugmer(i) = yrugm(j)
626                    rugos(i, nsrf) = yrugm(j)
627                 END IF
628                 agesno(i, nsrf) = yagesno(j)
629                 fqcalving(i, nsrf) = y_fqcalving(j)
630                 ffonte(i, nsrf) = y_ffonte(j)
631                 cdragh(i) = cdragh(i) + coefh(j, 1)
632                 cdragm(i) = cdragm(i) + coefm(j, 1)
633                 dflux_t(i) = dflux_t(i) + y_dflux_t(j)
634                 dflux_q(i) = dflux_q(i) + y_dflux_q(j)
635                 zu1(i) = zu1(i) + yu1(j)
636                 zv1(i) = zv1(i) + yv1(j)
637            END DO            END DO
638         END IF            IF (nsrf == is_ter) THEN
639         IF (nsrf == is_lic) THEN               DO j = 1, knon
640                    i = ni(j)
641                    qsol(i) = yqsol(j)
642                 END DO
643              END IF
644              IF (nsrf == is_lic) THEN
645                 DO j = 1, knon
646                    i = ni(j)
647                    run_off_lic_0(i) = y_run_off_lic_0(j)
648                 END DO
649              END IF
650              !$$$ PB ajout pour soil
651              ftsoil(:, :, nsrf) = 0.
652              DO k = 1, nsoilmx
653                 DO j = 1, knon
654                    i = ni(j)
655                    ftsoil(i, k, nsrf) = ytsoil(j, k)
656                 END DO
657              END DO
658    
659            DO j = 1, knon            DO j = 1, knon
660               i = ni(j)               i = ni(j)
661               run_off_lic_0(i) = y_run_off_lic_0(j)               DO k = 1, klev
662                    d_t(i, k) = d_t(i, k) + y_d_t(j, k)
663                    d_q(i, k) = d_q(i, k) + y_d_q(j, k)
664                    d_u(i, k) = d_u(i, k) + y_d_u(j, k)
665                    d_v(i, k) = d_v(i, k) + y_d_v(j, k)
666                    zcoefh(i, k) = zcoefh(i, k) + coefh(j, k)
667                 END DO
668            END DO            END DO
669         END IF  
670         !$$$ PB ajout pour soil            !cc diagnostic t, q a 2m et u, v a 10m
671         ftsoil(:, :, nsrf) = 0.  
        DO k = 1, nsoilmx  
672            DO j = 1, knon            DO j = 1, knon
673               i = ni(j)               i = ni(j)
674               ftsoil(i, k, nsrf) = ytsoil(j, k)               uzon(j) = yu(j, 1) + y_d_u(j, 1)
675            END DO               vmer(j) = yv(j, 1) + y_d_v(j, 1)
676         END DO               tair1(j) = yt(j, 1) + y_d_t(j, 1)
677                 qair1(j) = yq(j, 1) + y_d_q(j, 1)
678                 zgeo1(j) = rd*tair1(j)/(0.5*(ypaprs(j, 1)+ypplay(j, &
679                      1)))*(ypaprs(j, 1)-ypplay(j, 1))
680                 tairsol(j) = yts(j) + y_d_ts(j)
681                 rugo1(j) = yrugos(j)
682                 IF (nsrf == is_oce) THEN
683                    rugo1(j) = rugos(i, nsrf)
684                 END IF
685                 psfce(j) = ypaprs(j, 1)
686                 patm(j) = ypplay(j, 1)
687    
688         DO j = 1, knon               qairsol(j) = yqsurf(j)
           i = ni(j)  
           DO k = 1, klev  
              d_t(i, k) = d_t(i, k) + y_d_t(j, k)  
              d_q(i, k) = d_q(i, k) + y_d_q(j, k)  
              d_u(i, k) = d_u(i, k) + y_d_u(j, k)  
              d_v(i, k) = d_v(i, k) + y_d_v(j, k)  
              zcoefh(i, k) = zcoefh(i, k) + ycoefh(j, k)  
689            END DO            END DO
        END DO  
   
        !cc diagnostic t, q a 2m et u, v a 10m  
690    
691         DO j = 1, knon            CALL stdlevvar(klon, knon, nsrf, zxli, uzon, vmer, tair1, qair1, &
692            i = ni(j)                 zgeo1, tairsol, qairsol, rugo1, psfce, patm, yt2m, yq2m, &
693            uzon(j) = yu(j, 1) + y_d_u(j, 1)                 yt10m, yq10m, yu10m, yustar)
           vmer(j) = yv(j, 1) + y_d_v(j, 1)  
           tair1(j) = yt(j, 1) + y_d_t(j, 1)  
           qair1(j) = yq(j, 1) + y_d_q(j, 1)  
           zgeo1(j) = rd*tair1(j)/(0.5*(ypaprs(j, 1)+ypplay(j, &  
                1)))*(ypaprs(j, 1)-ypplay(j, 1))  
           tairsol(j) = yts(j) + y_d_ts(j)  
           rugo1(j) = yrugos(j)  
           IF (nsrf == is_oce) THEN  
              rugo1(j) = rugos(i, nsrf)  
           END IF  
           psfce(j) = ypaprs(j, 1)  
           patm(j) = ypplay(j, 1)  
694    
695            qairsol(j) = yqsurf(j)            DO j = 1, knon
696         END DO               i = ni(j)
697                 t2m(i, nsrf) = yt2m(j)
698                 q2m(i, nsrf) = yq2m(j)
699    
700         CALL stdlevvar(klon, knon, nsrf, zxli, uzon, vmer, tair1, qair1, zgeo1, &               ! u10m, v10m : composantes du vent a 10m sans spirale de Ekman
701              tairsol, qairsol, rugo1, psfce, patm, yt2m, yq2m, yt10m, yq10m, &               u10m(i, nsrf) = (yu10m(j)*uzon(j))/sqrt(uzon(j)**2+vmer(j)**2)
702              yu10m, yustar)               v10m(i, nsrf) = (yu10m(j)*vmer(j))/sqrt(uzon(j)**2+vmer(j)**2)
   
        DO j = 1, knon  
           i = ni(j)  
           t2m(i, nsrf) = yt2m(j)  
           q2m(i, nsrf) = yq2m(j)  
   
           ! u10m, v10m : composantes du vent a 10m sans spirale de Ekman  
           u10m(i, nsrf) = (yu10m(j)*uzon(j))/sqrt(uzon(j)**2+vmer(j)**2)  
           v10m(i, nsrf) = (yu10m(j)*vmer(j))/sqrt(uzon(j)**2+vmer(j)**2)  
703    
704         END DO            END DO
705    
706         DO i = 1, knon            CALL hbtm(knon, ypaprs, ypplay, yt2m, yt10m, yq2m, yq10m, yustar, &
707            y_cd_h(i) = ycoefh(i, 1)                 y_flux_t, y_flux_q, yu, yv, yt, yq, ypblh, ycapcl, yoliqcl, &
708            y_cd_m(i) = ycoefm(i, 1)                 ycteicl, ypblt, ytherm, ytrmb1, ytrmb2, ytrmb3, ylcl)
        END DO  
        CALL hbtm(knon, ypaprs, ypplay, yt2m, yt10m, yq2m, yq10m, yustar, &  
             y_flux_t, y_flux_q, yu, yv, yt, yq, ypblh, ycapcl, yoliqcl, &  
             ycteicl, ypblt, ytherm, ytrmb1, ytrmb2, ytrmb3, ylcl)  
   
        DO j = 1, knon  
           i = ni(j)  
           pblh(i, nsrf) = ypblh(j)  
           plcl(i, nsrf) = ylcl(j)  
           capcl(i, nsrf) = ycapcl(j)  
           oliqcl(i, nsrf) = yoliqcl(j)  
           cteicl(i, nsrf) = ycteicl(j)  
           pblt(i, nsrf) = ypblt(j)  
           therm(i, nsrf) = ytherm(j)  
           trmb1(i, nsrf) = ytrmb1(j)  
           trmb2(i, nsrf) = ytrmb2(j)  
           trmb3(i, nsrf) = ytrmb3(j)  
        END DO  
709    
        DO j = 1, knon  
           DO k = 1, klev + 1  
              i = ni(j)  
              q2(i, k, nsrf) = yq2(j, k)  
           END DO  
        END DO  
        !IM "slab" ocean  
        IF (nsrf == is_oce) THEN  
710            DO j = 1, knon            DO j = 1, knon
              ! on projette sur la grille globale  
711               i = ni(j)               i = ni(j)
712               IF (pctsrf_new(i, is_oce)>epsfra) THEN               pblh(i, nsrf) = ypblh(j)
713                  flux_o(i) = y_flux_o(j)               plcl(i, nsrf) = ylcl(j)
714               ELSE               capcl(i, nsrf) = ycapcl(j)
715                  flux_o(i) = 0.               oliqcl(i, nsrf) = yoliqcl(j)
716               END IF               cteicl(i, nsrf) = ycteicl(j)
717                 pblt(i, nsrf) = ypblt(j)
718                 therm(i, nsrf) = ytherm(j)
719                 trmb1(i, nsrf) = ytrmb1(j)
720                 trmb2(i, nsrf) = ytrmb2(j)
721                 trmb3(i, nsrf) = ytrmb3(j)
722            END DO            END DO
        END IF  
723    
        IF (nsrf == is_sic) THEN  
724            DO j = 1, knon            DO j = 1, knon
725               i = ni(j)               DO k = 1, klev + 1
726               ! On pondère lorsque l'on fait le bilan au sol :                  i = ni(j)
727               IF (pctsrf_new(i, is_sic)>epsfra) THEN                  q2(i, k, nsrf) = yq2(j, k)
728                  flux_g(i) = y_flux_g(j)               END DO
              ELSE  
                 flux_g(i) = 0.  
              END IF  
729            END DO            END DO
730              !IM "slab" ocean
        END IF  
        IF (ocean == 'slab  ') THEN  
731            IF (nsrf == is_oce) THEN            IF (nsrf == is_oce) THEN
732               tslab(1:klon) = ytslab(1:klon)               DO j = 1, knon
733               seaice(1:klon) = y_seaice(1:klon)                  ! on projette sur la grille globale
734                    i = ni(j)
735                    IF (pctsrf_new(i, is_oce)>epsfra) THEN
736                       flux_o(i) = y_flux_o(j)
737                    ELSE
738                       flux_o(i) = 0.
739                    END IF
740                 END DO
741            END IF            END IF
742         END IF  
743              IF (nsrf == is_sic) THEN
744                 DO j = 1, knon
745                    i = ni(j)
746                    ! On pondère lorsque l'on fait le bilan au sol :
747                    IF (pctsrf_new(i, is_sic)>epsfra) THEN
748                       flux_g(i) = y_flux_g(j)
749                    ELSE
750                       flux_g(i) = 0.
751                    END IF
752                 END DO
753    
754              END IF
755              IF (ocean == 'slab  ') THEN
756                 IF (nsrf == is_oce) THEN
757                    tslab(1:klon) = ytslab(1:klon)
758                    seaice(1:klon) = y_seaice(1:klon)
759                 END IF
760              END IF
761           end IF if_knon
762      END DO loop_surface      END DO loop_surface
763    
764      ! On utilise les nouvelles surfaces      ! On utilise les nouvelles surfaces

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

  ViewVC Help
Powered by ViewVC 1.1.21