Index: LMDZ.3.3/branches/rel-LF/libf/dyn3d/etat0_netcdf.F
===================================================================
--- LMDZ.3.3/branches/rel-LF/libf/dyn3d/etat0_netcdf.F	(revision 110)
+++ LMDZ.3.3/branches/rel-LF/libf/dyn3d/etat0_netcdf.F	(revision 112)
@@ -10,10 +10,13 @@
       !
       !
-      INTEGER, PARAMETER :: KIDIA=1, KFDIA=iim*(jjm-1)+2, 
-     .KLON=KFDIA-KIDIA+1,KLEV=llm
+c      INTEGER, PARAMETER :: KIDIA=1, KFDIA=iim*(jjm-1)+2, 
+c     .KLON=KFDIA-KIDIA+1,KLEV=llm
       !
 #include "comgeom2.h"
 #include "comvert.h"
 #include "comconst.h"
+#include "indicesol.h"
+#include "dimphy.h"
+#include "dimsoil.h"
       !
       REAL :: latfi(klon), lonfi(klon)
@@ -25,14 +28,31 @@
       REAL :: t3d(iip1, jjp1, llm), tpot(iip1, jjp1, llm)
       REAL :: q3d(iip1, jjp1, llm,nqmx), qsat(iip1, jjp1, llm)
-      REAL :: tsol(klon), qsol(klon), sn(klon), radsol(klon)
-      REAL :: deltat(klon), rugmer(klon), agesno(klon)
+      REAL :: tsol(klon), qsol(klon), sn(klon)
+      REAL :: tsolsrf(klon,nbsrf), qsolsrf(klon,nbsrf),snsrf(klon,nbsrf) 
+      REAL :: albe(klon,nbsrf), evap(klon,nbsrf)
+      REAL :: tsoil(klon,nsoilmx,nbsrf) 
+      REAL :: radsol(klon),rain_fall(klon), snow_fall(klon)
+      REAL :: solsw(klon), sollw(klon)
+      REAL :: deltat(klon), frugs(klon,nbsrf), agesno(klon),rugmer(klon)
       REAL :: zmea(iip1*jjp1), zstd(iip1*jjp1)
       REAL :: zsig(iip1*jjp1), zgam(iip1*jjp1), zthe(iip1*jjp1)
       REAL :: zpic(iip1*jjp1), zval(iip1*jjp1), rugsrel(iip1*jjp1)
       REAL :: qd(iip1, jjp1, llm)
-      !
+      REAL :: pctsrf(klon, nbsrf)
+      REAL :: t_ancien(klon,klev), q_ancien(klon,klev)      !
+      ! declarations pour lecture glace de mer
+      INTEGER :: iml_lic, jml_lic, llm_tmp, ttm_tmp, iret
+      INTEGER :: itaul(1), fid
+      REAL :: lev(1), date, dt
+      REAL, ALLOCATABLE, DIMENSION(:,:) :: lon_lic, lat_lic
+      REAL, ALLOCATABLE, DIMENSION(:)  :: dlon_lic, dlat_lic
+      REAL, ALLOCATABLE, DIMENSION (:,:) :: fraclic
+      REAL :: flic_tmp(iip1, jjp1)
+      REAL :: champint(iim, jjp1)
+      !
+
       CHARACTER*80 :: varname
       !
-      INTEGER :: i,j, ig, l
+      INTEGER :: i,j, ig, l, ji
       REAL :: xpi
       !
@@ -144,4 +164,22 @@
       !
       !
+
+
+C
+C on initialise les sous surfaces
+C
+      pctsrf=0.
+      !cree le masque a partir du fichier relief
+      varname = 'zmasq'
+      zmasq(:) = 0.
+      CALL startget(varname,iip1,jjp1,rlonv,rlatu,klon,zmasq,0.0)
+      WHERE (zmasq(1 : klon) .LE. EPSFRA)
+          zmasq(1 : klon) = 0.
+      END WHERE 
+      WRITE(*,*)zmasq
+
+
+
+
       varname = 'psol'
       psol(:,:) = 0.0
@@ -227,5 +265,5 @@
       ! This line needs to be replaced by a call to restget to get the values in the restart file
       tsol(:) = 0.0
-      CALL startget(varname, iip1, jjp1, rlonv, rlatu, klon, tsol, 0.0)
+      CALL startget(varname, iip1, jjp1, rlonv, rlatu, klon, tsol,0.0)
       !
       WRITE(*,*) 'TSOL construit :'
@@ -234,9 +272,9 @@
       varname = 'qsol'
       qsol(:) = 0.0
-      CALL startget(varname, iip1, jjp1, rlonv, rlatu, klon, qsol, 0.0)
+      CALL startget(varname, iip1, jjp1, rlonv, rlatu, klon, qsol,0.0)
       !
       varname = 'snow'
       sn(:) = 0.0
-      CALL startget(varname, iip1, jjp1, rlonv, rlatu, klon, sn, 0.0)
+      CALL startget(varname, iip1, jjp1, rlonv, rlatu, klon, sn,0.0)
       !
       varname = 'rads'
@@ -278,4 +316,88 @@
       CALL startget(varname,iip1,jjp1,rlonv,rlatu,klon,zval,0.0)
       rugsrel(:) = 0.0
+
+
+C
+C lecture du fichier glace de terre pour fixer la fraction de terre 
+C et de glace de terre
+C
+      CALL flininfo("landiceref.nc", iml_lic, jml_lic,llm_tmp, ttm_tmp
+     $    , fid)
+      ALLOCATE(lat_lic(iml_lic, jml_lic), stat=iret)
+      ALLOCATE(lon_lic(iml_lic, jml_lic), stat=iret)
+      ALLOCATE(dlon_lic(iml_lic), stat=iret)
+      ALLOCATE(dlat_lic(jml_lic), stat=iret)
+      ALLOCATE(fraclic(iml_lic, jml_lic), stat=iret)
+      CALL flinopen("landiceref.nc", .FALSE., iml_lic, jml_lic, llm_tmp
+     $    , lon_lic, lat_lic, lev, ttm_tmp, itaul, date, dt, fid)
+      CALL flinget(fid, 'landice', iml_lic, jml_lic, llm_tmp, ttm_tmp
+     $    , 1, 1, fraclic)
+      CALL flinclo(fid)
+C
+C interpolation sur la grille T du modele
+C
+      WRITE(*,*) 'dimensions de landice iml_lic, jml_lic : ', 
+     $    iml_lic, jml_lic
+c
+C sil les coordonnees sont en degres, on les transforme
+C
+      IF( MAXVAL( lon_lic(:,:) ) .GT. 2.0 * asin(1.0) )  THEN
+          lon_lic(:,:) = lon_lic(:,:) * 2.0* ASIN(1.0) / 180.
+      ENDIF 
+      IF( maxval( lat_lic(:,:) ) .GT. 2.0 * asin(1.0)) THEN 
+          lat_lic(:,:) = lat_lic(:,:) * 2.0 * asin(1.0) / 180.
+      ENDIF 
+
+      dlon_lic(1 : iml_lic) = lon_lic(1 : iml_lic, 1)
+      dlat_lic(1 : jml_lic) = lat_lic(1 , 1 : jml_lic) 
+C
+      CALL grille_m(iml_lic, jml_lic, dlon_lic, dlat_lic, fraclic
+     $    ,iim, jjp1,
+     $    rlonv, rlatu, flic_tmp(1 : iim, 1 : jjp1))
+c$$$      flic_tmp(1 : iim, 1 : jjp1) = champint(1: iim, 1 : jjp1)
+      flic_tmp(iip1, 1 : jjp1) = flic_tmp(1 , 1 : jjp1)
+C
+C passage sur la grille physique
+C
+      CALL gr_dyn_fi(1, iip1, jjp1, klon, flic_tmp,
+     $    pctsrf(1:klon, is_lic))
+C adequation avec le maque terre/mer
+      WHERE (pctsrf(1 : klon, is_lic) .LE. EPSFRA ) 
+          pctsrf(1 : klon, is_lic) = 0. 
+      END WHERE
+      WHERE (zmasq( 1 : klon) .LE. EPSFRA) 
+          pctsrf(1 : klon, is_lic) = 0.
+      END WHERE 
+      pctsrf(1 : klon, is_ter) = zmasq(1 : klon)
+      DO ji = 1, klon
+        IF (zmasq(ji) .GT. EPSFRA) THEN 
+            IF ( pctsrf(ji, is_lic) .GE. zmasq(ji)) THEN
+                pctsrf(ji, is_lic) = zmasq(ji)
+                pctsrf(ji, is_ter) = 0.
+            ELSE 
+                pctsrf(ji,is_ter) = zmasq(ji) - pctsrf(ji, is_lic)
+            ENDIF 
+        ENDIF 
+      END DO 
+C
+C sous surface ocean et glace de mer (pour demarrer on met glace de mer a 0)
+C
+      pctsrf(1 : klon, is_oce) = (1. - zmasq(1 : klon))
+      WHERE (pctsrf(1 : klon, is_oce) .LT. EPSFRA)
+          pctsrf(1 : klon, is_oce) = 0.
+      END WHERE 
+C
+C verif que somme des sous surface = 1
+C
+      ji=count( (abs( sum(pctsrf(1 : klon, 1 : nbsrf), dim = 2)) - 1.0 ) 
+     $    .GT. EPSFRA)
+      IF (ji .NE. 0) THEN
+          WRITE(*,*) 'pb repartition sous maille pour ',ji,' points'
+      ENDIF 
+
+
+
+
+
 C Calcul intermediaire
 c 
@@ -333,10 +455,55 @@
       solaire   = 1370.0
 
-      call physdem(lonfi, latfi, phystep,radpas,co2_ppm,
-     .                   solaire,tsol, qsol,
-     .                   sn, radsol, deltat, rugmer,
-     .                   agesno, zmea, zstd, zsig,
-     .                   zgam, zthe, zpic, zval,
-     .                   rugsrel)
+c      call physdem(lonfi, latfi, phystep,radpas,co2_ppm,
+c     .                   solaire,tsol, qsol,
+c     .                   sn, radsol, deltat, rugmer,
+c     .                   agesno, zmea, zstd, zsig,
+c     .                   zgam, zthe, zpic, zval,
+c     .                   rugsrel)
+
+c
+c Initialisation 
+c tsol, qsol, sn,albe, evap,tsoil,rain_fall, snow_fall,solsw, sollw,frugs
+c
+      tsolsrf(:,is_ter) = tsol
+      tsolsrf(:,is_lic) = tsol
+      tsolsrf(:,is_oce) = tsol
+      tsolsrf(:,is_sic) = tsol
+      snsrf(:,is_ter) = sn
+      snsrf(:,is_lic) = sn
+      snsrf(:,is_oce) = sn
+      snsrf(:,is_sic) = sn
+      albe(:,is_ter) = 0.08
+      albe(:,is_lic) = 0.6
+      albe(:,is_oce) = 0.5
+      albe(:,is_sic) = 0.6
+      evap(:,:) = 0.
+      qsolsrf(:,is_ter) = qsol
+      qsolsrf(:,is_lic) = qsol
+      qsolsrf(:,is_oce) = 150.
+      qsolsrf(:,is_sic) = 150.
+      do i = 1, nbsrf
+        do j = 1, nsoilmx
+          tsoil(:,j,i) = tsol
+        enddo
+      enddo
+      rain_fall = 0.; snow_fall = 0.
+      solsw = 165.
+      sollw = -53.
+      t_ancien = 273.15
+      q_ancien = 0.
+      agesno = 0.
+      deltat = 0.
+      frugs(:,is_oce) = rugmer
+      frugs(:,is_ter) = rugmer
+      frugs(:,is_lic) = rugmer
+      frugs(:,is_sic) = rugmer
+
+      call phyredem("startphy.nc",phystep,radpas, co2_ppm, solaire,
+     $    latfi, lonfi, pctsrf, tsolsrf, tsoil, deltat, qsolsrf, snsrf, 
+     $    albe, evap, rain_fall, snow_fall, solsw, sollw, 
+     $    radsol, frugs,  agesno, 
+     $    zmea, zstd, zsig, zgam, zthe, zpic, zval, rugsrel, 
+     $    t_ancien, q_ancien)
 
 C     Sortie Visu pour les champs dynamiques
Index: LMDZ.3.3/branches/rel-LF/libf/dyn3d/gcm.F
===================================================================
--- LMDZ.3.3/branches/rel-LF/libf/dyn3d/gcm.F	(revision 110)
+++ LMDZ.3.3/branches/rel-LF/libf/dyn3d/gcm.F	(revision 112)
@@ -437,5 +437,5 @@
 
 c rajout debug
-c        lafin = .true.
+        lafin = .true.
 
         CALL calfis( nqmx, lafin ,rdayvrai,rday_ecri,time  ,
@@ -503,8 +503,8 @@
 
 c ajout debug
-c               IF( lafin ) then  
-c                 abort_message = 'Simulation finished'
-c                 call abort_gcm(modname,abort_message,0)
-c               ENDIF
+               IF( lafin ) then  
+                 abort_message = 'Simulation finished'
+                 call abort_gcm(modname,abort_message,0)
+               ENDIF
         
 c   ********************************************************************
Index: LMDZ.3.3/branches/rel-LF/libf/dyn3d/modif_etat0.F
===================================================================
--- LMDZ.3.3/branches/rel-LF/libf/dyn3d/modif_etat0.F	(revision 110)
+++ LMDZ.3.3/branches/rel-LF/libf/dyn3d/modif_etat0.F	(revision 112)
@@ -245,5 +245,5 @@
 C Ecriture etat initial physique
 C
-      WRITE(*,*) 'zmasq avant phyredem'
+      WRITE(*,*)'zmasq avant phyredem'
       WRITE(*,*) zmasq
 
Index: LMDZ.3.3/branches/rel-LF/libf/phylmd/YOMCST.inc
===================================================================
--- LMDZ.3.3/branches/rel-LF/libf/phylmd/YOMCST.inc	(revision 110)
+++ LMDZ.3.3/branches/rel-LF/libf/phylmd/YOMCST.inc	(revision 112)
@@ -1,22 +1,22 @@
 ! A1.0 Fundamental constants
-      REAL RPI,RCLUM,RHPLA,RKBOL,RNAVO
+      REAL :: RPI,RCLUM,RHPLA,RKBOL,RNAVO
 ! A1.1 Astronomical constants
-      REAL RDAY,REA,REPSM,RSIYEA,RSIDAY,ROMEGA
+      REAL :: RDAY,REA,REPSM,RSIYEA,RSIDAY,ROMEGA
 ! A1.1.bis Constantes concernant l'orbite de la Terre:
-      REAL R_ecc, R_peri, R_incl
+      REAL :: R_ecc, R_peri, R_incl
 ! A1.2 Geoide
-      REAL RA,RG,R1SA
+      REAL :: RA,RG,R1SA
 ! A1.3 Radiation
-      REAL RSIGMA,RI0
+      REAL :: RSIGMA,RI0
 ! A1.4 Thermodynamic gas phase
-      REAL R,RMD,RMV,RD,RV,RCPD,RCPV,RCVD,RCVV
-      REAL RKAPPA,RETV
+      REAL :: R,RMD,RMV,RD,RV,RCPD,RCPV,RCVD,RCVV
+      REAL :: RKAPPA,RETV
 ! A1.5,6 Thermodynamic liquid,solid phases
-      REAL RCW,RCS
+      REAL :: RCW,RCS
 ! A1.7 Thermodynamic transition of phase
-      REAL RLVTT,RLSTT,RLMLT,RTT,RATM
+      REAL :: RLVTT,RLSTT,RLMLT,RTT,RATM
 ! A1.8 Curve of saturation
-      REAL RESTT,RALPW,RBETW,RGAMW,RALPS,RBETS,RGAMS
-      REAL RALPD,RBETD,RGAMD
+      REAL :: RESTT,RALPW,RBETW,RGAMW,RALPS,RBETS,RGAMS
+      REAL :: RALPD,RBETD,RGAMD
 !
       COMMON/YOMCST/RPI   ,RCLUM ,RHPLA ,RKBOL ,RNAVO &
Index: LMDZ.3.3/branches/rel-LF/libf/phylmd/clmain.F
===================================================================
--- LMDZ.3.3/branches/rel-LF/libf/phylmd/clmain.F	(revision 110)
+++ LMDZ.3.3/branches/rel-LF/libf/phylmd/clmain.F	(revision 112)
@@ -1,9 +1,10 @@
-      SUBROUTINE clmain(dtime,pctsrf,t,q,u,v,
+      SUBROUTINE clmain(dtime,itap,pctsrf,
+     .                  t,q,u,v,
      .                  jour, rmu0,
-     .                  ok_veget,ts,
+     .                  ok_veget, ocean, npas, nexca, ts,
      .                  paprs,pplay,radsol,snow,qsol,evap,albe,
      .                  rain_f, snow_f, solsw, sollw, fder,
      .                  rlon, rlat, rugos,
-     .                  debut, lafin,
+     .                  debut, lafin, agesno,
      .                  d_t,d_q,d_u,d_v,d_ts,
      .                  flux_t,flux_q,flux_u,flux_v,cdragh,cdragm,
@@ -33,4 +34,5 @@
 c Arguments:
 c dtime----input-R- interval du temps (secondes)
+c itap-----input-I- numero du pas de temps
 c t--------input-R- temperature (K)
 c q--------input-R- vapeur d'eau (kg/kg)
@@ -75,4 +77,5 @@
 c
       REAL dtime
+      integer itap
       REAL t(klon,klev), q(klon,klev)
       REAL u(klon,klev), v(klon,klev)
@@ -84,9 +87,11 @@
       REAL dflux_t(klon), dflux_q(klon)
       REAL flux_u(klon,klev, nbsrf), flux_v(klon,klev, nbsrf)
-      REAL rugmer(klon)
+      REAL rugmer(klon), agesno(klon)
       REAL cdragh(klon), cdragm(klon)
       integer jour            ! jour de l'annee en cours
       real rmu0(klon)         ! cosinus de l'angle solaire zenithal
       LOGICAL debut, lafin, ok_veget
+      character*6 ocean
+      integer npas, nexca
 cAA      INTEGER itr
 cAA      REAL tr(klon,klev,nbtr)
@@ -288,7 +293,4 @@
         yv1(j) = v1lay(i)
         yrads(j) = totalflu(i)
-c        ycal(j) = cal(i)
-c        ybeta(j) = beta(i)
-c        ydif(j) = dif_grnd(i)
         ypaprs(j,klev+1) = paprs(i,klev+1)
       ENDDO
@@ -331,6 +333,9 @@
       ytauy = y_flux_v(:,1)
 c calculer la diffusion de "q" et de "h"
-      CALL clqh(knon, dtime, nsrf, ni, pctsrf, rlon, rlat,
-     e          jour, rmu0,
+      CALL clqh(dtime, itap, jour, debut,lafin,
+     e          rlon, rlat,
+     e          knon, nsrf, ni, pctsrf,
+     e          ok_veget, ocean, npas, nexca,
+     e          rmu0,
      e          yu1, yv1, ycoefh,
      e          yt,yq,yts,ypaprs,ypplay,
@@ -338,5 +343,5 @@
      e          yrain_f, ysnow_f, yfder, ytaux, ytauy, 
      e          ysollw, ysolsw,
-     s          pctsrf_new,
+     s          pctsrf_new, agesno,
      s          y_d_t, y_d_q, y_d_ts, yz0_new,
      s          y_flux_t, y_flux_q, y_dflux_t, y_dflux_q)
@@ -349,15 +354,4 @@
       ENDDO
       ENDIF
-c
-cAA MAINTENANT DANS PHYTRAC
-cAAc calculer la diffusion des traceurs
-cAA      IF (itr.GE.1) THEN
-cAA      DO it = 1, itr
-cAA      CALL cltrac(knon,dtime,ycoefh, yt, ytr(1,1,it), yflxsrf(1,it),
-cAA     e            ypaprs, ypplay, ydelp,
-cAA     s            y_d_tr(1,1,it))
-cAA      ENDDO
-cAA      ENDIF
-c
       DO j = 1, knon
          y_dflux_t(j) = y_dflux_t(j) * ypct(j)
@@ -427,29 +421,19 @@
       ENDDO
 c
-cAA      IF (itr.GE.1) THEN
-cAA      DO it = 1, itr
-cAA      DO k = 1, klev
-cAA      DO j = 1, knon
-cAA         y_d_tr(j,k,it) = y_d_tr(j,k,it) * ypct(j)
-cAA      ENDDO
-cAA      ENDDO
-cAA      ENDDO
-cAA      DO j = 1, knon
-cAA      i = ni(j)
-cAA      DO it = 1, itr
-cAA      DO k = 1, klev
-cAA         d_tr(i,k,it) = d_tr(i,k,it) + y_d_tr(j,k,it)
-cAA      ENDDO
-cAA      ENDDO
-cAA      ENDDO
-cAA      ENDIF
-c
 99999 CONTINUE
 c
+C
+C On utilise les nouvelles surfaces
+C A rajouter: conservation de l'albedo
+C
+      pctsrf = pctsrf_new
 
       RETURN
       END
-      SUBROUTINE clqh(knon,dtime,nisurf,knindex,pctsrf, rlon, rlat,
-     e                jour, rmu0
+      SUBROUTINE clqh(dtime,itime, jour,debut,lafin,
+     e                rlon, rlat,
+     e                knon, nisurf, knindex, pctsrf,
+     e                ok_veget, ocean, npas, nexca,
+     e                rmu0,
      e                u1lay,v1lay,coef,
      e                t,q,ts,paprs,pplay,
@@ -457,5 +441,5 @@
      e                precip_rain, precip_snow, fder, taux, tauy,
      e                lwdown, swdown, 
-     s                pctsrf_new,
+     s                pctsrf_new, agesno,
      s                d_t, d_q, d_ts, z0_new, 
      s                flux_t, flux_q,dflux_s,dflux_l)
@@ -499,4 +483,5 @@
       REAL qsol(klon)         ! humidite de la surface
       real precip_rain(klon), precip_snow(klon)
+      REAL agesno(klon)
       integer jour            ! jour de l'annee en cours
       real rmu0(klon)         ! cosinus de l'angle solaire zenithal
@@ -504,4 +489,8 @@
       real pctsrf(klon,nbsrf)
       real rlon(klon), rlat(klon)
+      logical ok_veget 
+      character*6 ocean
+      integer npas, nexca
+
 c
       REAL d_t(klon,klev)     ! incrementation de "t"
@@ -565,5 +554,5 @@
       integer itime
       integer nisurf
-      logical debut, lafin, ok_veget
+      logical debut, lafin
       real zlev1(klon)
       real fder(klon), taux(klon), tauy(klon)
@@ -575,5 +564,4 @@
       real p1lay(klon)
       real coef1lay(klon)
-      character*6 ocean
 
 ! Parametres de sortie
@@ -691,7 +679,4 @@
 C Appel a interfsurf (appel generique) routine d'interface avec la surface
 
-      ok_veget = .false.
-      ocean = 'force '
-
       petAcoef=zx_ch(:,1)
       peqAcoef=zx_cq(:,1)
@@ -717,7 +702,8 @@
      e albedo, snow, qsol,
      e ts, p1lay, psref, radsol,
-     e ocean,zmasq,
+     e ocean, npas, nexca, zmasq,
      s evap, fluxsens, fluxlat, dflux_l, dflux_s,              
-     s tsol_rad, tsurf_new, alb_new, emis_new, z0_new, pctsrf_new)
+     s tsol_rad, tsurf_new, alb_new, emis_new, z0_new, 
+     s pctsrf_NEW, agesno)
 
       flux_t(:,1) = fluxsens
Index: LMDZ.3.3/branches/rel-LF/libf/phylmd/interface_surf.F90
===================================================================
--- LMDZ.3.3/branches/rel-LF/libf/phylmd/interface_surf.F90	(revision 110)
+++ LMDZ.3.3/branches/rel-LF/libf/phylmd/interface_surf.F90	(revision 112)
@@ -36,7 +36,8 @@
 
 #include "YOMCST.inc"
-
-
-  ! run_off      ruissellement total
+#include "indicesol.inc"
+
+
+! run_off      ruissellement total
   real, allocatable, dimension(:),save    :: run_off
 
@@ -56,7 +57,7 @@
       & albedo, snow, qsol, &
       & tsurf, p1lay, ps, radsol, &
-      & ocean, zmasq, &
+      & ocean, npas, nexca, zmasq, &
       & evap, fluxsens, fluxlat, dflux_l, dflux_s, &              
-      & tsol_rad, tsurf_new, alb_new, emis_new, z0_new, pctsrf_new)
+      & tsol_rad, tsurf_new, alb_new, emis_new, z0_new, pctsrf_new, agesno)
 
 
@@ -125,5 +126,4 @@
 !   pctsrf_new   nouvelle repartition des surfaces
 
-  include 'indicesol.h'
 
 ! Parametres d'entree
@@ -153,4 +153,5 @@
   real, dimension(klon), intent(IN) :: fder, taux, tauy
   character (len = 6)  :: ocean
+  integer              :: npas, nexca ! nombre et pas de temps couplage
   real, dimension(knon), intent(INOUT) :: evap, snow, qsol
 
@@ -161,4 +162,5 @@
   real, dimension(knon), intent(OUT):: dflux_l, dflux_s
   real, dimension(klon,nbsrf), intent(OUT) :: pctsrf_new
+  real, dimension(klon), intent(INOUT):: agesno
 
 ! Local
@@ -166,14 +168,15 @@
   character (len = 80) :: abort_message 
   logical, save        :: first_call = .true.
-  integer              :: error
+  INTEGER              :: error, ii
   logical              :: check = .true.
   real, dimension(knon):: cal, beta, dif_grnd, capsol
   real, parameter      :: calice=1.0/(5.1444e+06*0.15), tau_gl=1./86400.*5.
   real, parameter      :: calsno=1./(2.3867e+06*.15)
-  integer              :: nexca !pas de temps couplage
   real, dimension(knon):: alb_ice
   real, dimension(knon):: tsurf_temp
-  real, dimension(klon):: agesno, alb_neig_grid, alb_eau
+  real, dimension(klon):: alb_neig_grid, alb_eau
   real, dimension(knon):: alb_neig
+  REAL, DIMENSION(knon):: lmt_rug, lmt_alb
+  real, DIMENSION(knon):: zfra
 
   if (check) write(*,*) 'Entree ', modname
@@ -203,4 +206,5 @@
       abort_message='voir ci-dessus'
       call abort_gcm(modname,abort_message,1)
+    endif
   endif
   first_call = .false.
@@ -233,8 +237,8 @@
 !
 
-  CALL albsno(agesno,alb_neig_grid)  
-!
-!
-!
+  CALL albsno(klon,agesno,alb_neig_grid)  
+  
+  
+ 
     if (.not. ok_veget) then
 !
@@ -256,9 +260,19 @@
      & klon, nisurf, knon, knindex, debut,  &
      & lmt_alb, lmt_rug)
-       alb_neig = alb_neig_grid(knindex)
+!
+! Pb compilo sun
+!       alb_neig = alb_neig_grid(knindex)
+!      alb_new = alb_neig*zfra + lmt_alb(knindex)*(1.0-zfra)
+!      z0_new = lmt_rug(knindex)
+!
+       DO ii = 1, knon
+         alb_neig(ii) = alb_neig_grid(knindex(ii))
+         alb_new(ii) = lmt_alb(knindex(ii))
+       enddo
        zfra = MAX(0.0,MIN(1.0,snow/(snow+10.0)))
-       alb_new = alb_neig*zfra + lmt_alb(knindex)*(1.0-zfra)
-       z0_new = lmt_rug(knindex)
-   
+       alb_new = alb_neig*zfra + alb_new*(1.0-zfra)
+       DO ii = 1, knon
+         z0_new(ii) = lmt_rug(knindex(ii))
+       enddo   
     else
 !
@@ -293,10 +307,8 @@
       call interfoce(itime, dtime, &
       & klon, iim, jjm, nisurf, pctsrf, knon, knindex, rlon, rlat, &
-      & ocean, nexca, debut, lafin, &
+      & ocean, npas, nexca, debut, lafin, &
       & swdown, lwdown, precip_rain, precip_snow, evap, tsurf, &
       & fder, albedo, taux, tauy, zmasq, &
       & tsurf_new, alb_new, alb_ice, pctsrf_new)
-
-      tsurf_temp = tsurf_new
 
 !    else if (ocean == 'slab  ') then
@@ -310,4 +322,5 @@
     endif
 
+    tsurf_temp = tsurf_new
     cal = 0.
     beta = 1.
@@ -324,11 +337,12 @@
 !
 
-     if ( minval(rmu0) == maxval(rmu0) && minval(rmu0) = -999.999 ) then
-       CALL alboc(FLOAT(jour),rlat,alb_eau)
-     else  ! cycle diurne
-       CALL alboc_cd(rmu0,alb_eau)
-     endif
-     alb_new = alb_eau(knindex)
-
+    if ( minval(rmu0) == maxval(rmu0) .and. minval(rmu0) == -999.999 ) then
+      CALL alboc(FLOAT(jour),rlat,alb_eau)
+    else  ! cycle diurne
+      CALL alboc_cd(rmu0,alb_eau)
+    endif
+    DO ii =1, knon
+      alb_new(ii) = alb_eau(knindex(ii))
+    enddo
 !
   else if (nisurf == is_sic) then
@@ -341,9 +355,8 @@
 !
     if (ocean == 'couple') then
-      nexca = 0
 
       call interfoce(itime, dtime, &
       & klon, iim, jjm, nisurf, pctsrf, knon, knindex, rlon, rlat, &
-      & ocean, nexca, debut, lafin, &
+      & ocean, npas, nexca, debut, lafin, &
       & swdown, lwdown, precip_rain, precip_snow, evap, tsurf, &
       & fder, albedo, taux, tauy, zmasq, &
@@ -357,9 +370,8 @@
 !      call interfoce(nisurf)
     else                              ! lecture conditions limites
-!      call interfoce(itime, dtime, jour, & 
-!     &  klon, nisurf, knon, knindex, &
-!     &  debut, &
-!     &  tsurf_new, pctsrf_new)
-!   endif
+      call interfoce(itime, dtime, jour, & 
+     &  klon, nisurf, knon, knindex, &
+     &  debut, &
+     &  tsurf_new, pctsrf_new)
 
       cal = calice
@@ -381,5 +393,7 @@
 !
        zfra = MAX(0.0,MIN(1.0,snow/(snow+10.0)))
-       alb_neig = alb_neig_grid(knindex)
+       DO ii = 1, knon
+         alb_neig = alb_neig_grid(knindex(ii))
+       enddo
        alb_new = alb_neig*zfra + 0.6 * (1.0-zfra)
 
@@ -409,5 +423,7 @@
 !
        zfra = MAX(0.0,MIN(1.0,snow/(snow+10.0)))
-       alb_neig = alb_neig_grid(knindex)
+       DO ii =1, knon
+         alb_neig = alb_neig_grid(knindex(ii))
+       enddo
        alb_new = alb_neig*zfra + 0.6 * (1.0-zfra)
 
@@ -634,5 +650,5 @@
   SUBROUTINE interfoce_cpl(itime, dtime, &
       & klon, iim, jjm, nisurf, pctsrf, knon, knindex, rlon, rlat, &
-      & ocean, nexca, debut, lafin, &
+      & ocean, npas, nexca, debut, lafin, &
       & swdown, lwdown, precip_rain, precip_snow, evap, tsurf, &
       & fder, albsol, taux, tauy, zmasq, &
@@ -691,5 +707,4 @@
 !
 
-#include 'indicesol.h'
 
 ! Parametres d'entree
@@ -708,5 +723,5 @@
   real, dimension(knon), intent(IN) :: precip_rain, precip_snow
   real, dimension(knon), intent(IN) :: tsurf, fder, albsol, taux, tauy
-  integer              :: nexca
+  INTEGER              :: nexca, npas
   real, dimension(klon), intent(IN) :: zmasq
 
@@ -719,5 +734,4 @@
 ! Variables locales
   integer                    :: j, error, sum_error, ig
-  integer                    :: npas
   character (len = 20) :: modname = 'interfoce_cpl'
   character (len = 80) :: abort_message
@@ -792,5 +806,5 @@
 ! initialisation couplage
 !
-    call inicma(npas, nexca, dtime)
+    call inicma(npas , nexca, dtime)
 !
 ! 1ere lecture champs ocean
@@ -824,7 +838,7 @@
   endif ! fin if (debut)
 
-!! fichier restart et fichiers histoires
-
-!! calcul des fluxs a passer
+! fichier restart et fichiers histoires
+
+! calcul des fluxs a passer
 
   cpl_sols(:,nisurf) = cpl_sols(:,nisurf) + swdown      / FLOAT(nexca)
@@ -1037,5 +1051,4 @@
 !
 
-#include "indicesol.h"
 
 ! Parametres d'entree
@@ -1055,5 +1068,5 @@
 ! Variables locales
   integer     :: ii
-  integer     :: lmt_pas     ! frequence de lecture des conditions limites 
+  INTEGER,save :: lmt_pas     ! frequence de lecture des conditions limites 
                              ! (en pas de physique)
   logical,save :: deja_lu    ! pour indiquer que le jour a lire a deja
@@ -1063,6 +1076,6 @@
   character (len = 20) :: modname = 'interfoce_lim'
   character (len = 80) :: abort_message
-  character (len = 20) :: fich ='limit'
-  logical     :: newlmt = .false.
+  character (len = 20) :: fich ='limit.nc'
+  LOGICAL     :: newlmt = .TRUE.
   logical     :: check = .true.
 ! Champs lus dans le fichier de CL
@@ -1079,5 +1092,5 @@
 !
     
-  if (debut) then
+  if (debut .and. .not. allocated(sst_lu)) then
     lmt_pas = nint(86400./dtime * 1.0) ! pour une lecture une fois par jour
     jour_lu = jour - 1
@@ -1096,4 +1109,5 @@
 ! Ouverture du fichier
 !
+    fich = trim(fich)
     ierr = NF_OPEN (fich, NF_NOWRITE,nid)
     if (ierr.NE.NF_NOERR) then
@@ -1119,7 +1133,7 @@
       endif
 #ifdef NC_DOUBLE
-      ierr = NF_GET_VARA_DOUBLE(nid,nvarid,start,epais,pctsrf_new(1,is_oce))
+      ierr = NF_GET_VARA_DOUBLE(nid,nvarid,start,epais,pct_tmp(1,is_oce))
 #else
-      ierr = NF_GET_VARA_REAL(nid,nvarid,start,epais,pctsrf_new(1,is_oce))
+      ierr = NF_GET_VARA_REAL(nid,nvarid,start,epais,pct_tmp(1,is_oce))
 #endif
       if (ierr /= NF_NOERR) then
@@ -1136,7 +1150,7 @@
       endif
 #ifdef NC_DOUBLE
-      ierr = NF_GET_VARA_DOUBLE(nid,nvarid,start,epais,pctsrf_new(1,is_sic))
+      ierr = NF_GET_VARA_DOUBLE(nid,nvarid,start,epais,pct_tmp(1,is_sic))
 #else
-      ierr = NF_GET_VARA_REAL(nid,nvarid,start,epais,pctsrf_new(1,is_sic))
+      ierr = NF_GET_VARA_REAL(nid,nvarid,start,epais,pct_tmp(1,is_sic))
 #endif
       if (ierr /= NF_NOERR) then
@@ -1153,7 +1167,7 @@
       endif
 #ifdef NC_DOUBLE
-      ierr = NF_GET_VARA_DOUBLE(nid,nvarid,start,epais,pctsrf_new(1,is_ter))
+      ierr = NF_GET_VARA_DOUBLE(nid,nvarid,start,epais,pct_tmp(1,is_ter))
 #else
-      ierr = NF_GET_VARA_REAL(nid,nvarid,start,epais,pctsrf_new(1,is_ter))
+      ierr = NF_GET_VARA_REAL(nid,nvarid,start,epais,pct_tmp(1,is_ter))
 #endif
       if (ierr /= NF_NOERR) then
@@ -1170,7 +1184,7 @@
       endif
 #ifdef NC_DOUBLE
-      ierr = NF_GET_VARA_DOUBLE(nid,nvarid,start,epais,pctsrf_new(1,is_lic))
+      ierr = NF_GET_VARA_DOUBLE(nid,nvarid,start,epais,pct_tmp(1,is_lic))
 #else
-      ierr = NF_GET_VARA_REAL(nid,nvarid,start,epais,pctsrf_new(1,is_lic))
+      ierr = NF_GET_VARA_REAL(nid,nvarid,start,epais,pct_tmp(1,is_lic))
 #endif
       if (ierr /= NF_NOERR) then
@@ -1240,5 +1254,9 @@
 ! Recopie des variables dans les champs de sortie
 !
-  lmt_sst = sst_lu(knindex)
+  do ii = 1, knon
+    lmt_sst(ii) = sst_lu(knindex(ii))
+  enddo
+! je peux pas utiliser la ligne suivante a cause du compilo Sun
+!  lmt_sst = sst_lu(knindex)
   pctsrf_new = pct_tmp
 
@@ -1275,5 +1293,4 @@
 !
 
-#include "indicesol.h"
 
 ! Parametres d'entree
@@ -1301,5 +1318,5 @@
   character (len = 20) :: modname = 'interfoce_lim'
   character (len = 80) :: abort_message
-  character (len = 20) :: fich ='limit'
+  character (len = 20) :: fich ='limit.nc'
   logical     :: newlmt = .false.
   logical     :: check = .true.
@@ -1323,7 +1340,7 @@
   endif
 
-  if ((jour - jour_lu_sur) /= 0) deja_lu = .false.
+  if ((jour - jour_lu_sur) /= 0) deja_lu_sur = .false.
   
-  if (check) write(*,*)modname,':: jour_lu, deja_lu_sur', jour_lu, deja_lu_sur 
+  if (check) write(*,*)modname,':: jour_lu_sur, deja_lu_sur', jour_lu_sur, deja_lu_sur 
 
 ! Tester d'abord si c'est le moment de lire le fichier
@@ -1332,4 +1349,5 @@
 ! Ouverture du fichier
 !
+    fich = trim(fich)
     ierr = NF_OPEN (fich, NF_NOWRITE,nid)
     if (ierr.NE.NF_NOERR) then
@@ -1339,5 +1357,5 @@
 ! 
 ! La tranche de donnees a lire:
-!
+ 
     start(1) = 1
     start(2) = jour + 1
@@ -1389,6 +1407,8 @@
 ! Recopie des variables dans les champs de sortie
 !
-  lmt_alb = alb_lu(knindex)
-  lmt_rug = rug_lu(knindex)
+  DO ii = 1, knon
+    lmt_alb(ii) = alb_lu(knindex(ii))
+    lmt_rug(ii) = rug_lu(knindex(ii))
+  enddo
 
   END SUBROUTINE interfsur_lim
@@ -1441,5 +1461,4 @@
 #include "YOETHF.inc"
 #include "FCTTRE.inc"
-#include 'indicesol.h'
 
 ! Parametres d'entree
@@ -1606,118 +1625,4 @@
 !#########################################################################
 !
-
-  SUBROUTINE sol_dem_write(itime, klon, rlon, rlat, &
-  &                     pctsrf_new,tsurf_new,alb_new)
-
-! Routine d'ecriture de l'etat de redemarrage pour le sol
-!
-! L.Fairhead
-!
-! input:
-!   itime        numero du pas de temps
-!   klon         nombre total de points de grille
-!   rlon         longitudes
-!   rlat         latitudes
-!   tsurf_new    temperature au sol
-!   alb_new      albedo
-!   pctsrf_new   repartition des surfaces
-
-  include 'indicesol.h'
-#include 'temps.inc'
-  include 'netcdf.inc'
-
-! Parametres d'entree
-  integer, intent(IN) :: itime
-  integer, intent(IN) :: klon
-  real, dimension(klon), intent(IN) :: rlon, rlat
-  real, dimension(klon,nbsrf), intent(IN)  :: tsurf_new, alb_new
-  real, dimension(klon,nbsrf), intent(IN) :: pctsrf_new
-
-! Variables locales
-  integer             :: ierr, nid
-  integer             :: idim1, idim2, idim3
-  integer,parameter   :: length = 100
-  character (len = 20) :: modname = 'sol_dem_write'
-  character (len = 80) :: abort_message 
-  real, dimension(length) :: tab_cntrl = 0.
-  integer                 :: nvarid
-
-  ierr = NF_CREATE('restartsol', NF_CLOBBER, nid)
-  IF (ierr.NE.NF_NOERR) THEN
-    abort_message=' Pb d''ouverture du fichier restartsol'
-    CALL abort_gcm(modname,abort_message,ierr)
-  ENDIF
-
-  ierr = NF_PUT_ATT_TEXT (nid, NF_GLOBAL, "title", 23,    &
-     &                    "Fichier redemmarage sol")
-  ierr = NF_DEF_DIM (nid, "index", length, idim1)
-  ierr = NF_DEF_DIM (nid, "points_physiques", klon, idim2)
-  ierr = NF_DEF_DIM (nid, "nombre_surfaces", nbsrf, idim3)
-  ierr = NF_ENDDEF(nid)
-
-  tab_cntrl(13) = day_end
-  tab_cntrl(14) = anne_ini
-
-  ierr = NF_REDEF (nid)
-#ifdef NC_DOUBLE
-  ierr = NF_DEF_VAR (nid, "controle", NF_DOUBLE, 1, idim1,nvarid)
-#else
-  ierr = NF_DEF_VAR (nid, "controle", NF_FLOAT, 1, idim1,nvarid)
-#endif
-  ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 22,  &
-     &                        "Parametres de controle")
-  ierr = NF_ENDDEF(nid)
-#ifdef NC_DOUBLE
-  ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,tab_cntrl)
-#else
-  ierr = NF_PUT_VAR_REAL (nid,nvarid,tab_cntrl)
-#endif
-
-  ierr = NF_REDEF (nid)
-#ifdef NC_DOUBLE
-  ierr = NF_DEF_VAR (nid, "longitude", NF_DOUBLE, 1, idim2,nvarid)
-#else
-  ierr = NF_DEF_VAR (nid, "longitude", NF_FLOAT, 1, idim2,nvarid)
-#endif
-  ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 32,  &
-     &               "Longitudes de la grille physique")
-  ierr = NF_ENDDEF(nid)
-#ifdef NC_DOUBLE
-  ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,rlon)
-#else
-  ierr = NF_PUT_VAR_REAL (nid,nvarid,rlon)
-#endif
-!
-  ierr = NF_REDEF (nid)
-#ifdef NC_DOUBLE
-  ierr = NF_DEF_VAR (nid, "latitude", NF_DOUBLE, 1, idim2,nvarid)
-#else
-  ierr = NF_DEF_VAR (nid, "latitude", NF_FLOAT, 1, idim2,nvarid)
-#endif
-  ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 31,    &
-     &                        "Latitudes de la grille physique")
-  ierr = NF_ENDDEF(nid)
-#ifdef NC_DOUBLE
-  ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,rlat)
-#else
-  ierr = NF_PUT_VAR_REAL (nid,nvarid,rlat)
-#endif
-        ierr = NF_REDEF (nid)
-#ifdef NC_DOUBLE
-        ierr = NF_DEF_VAR (nid, "TS", NF_DOUBLE, 1, idim2,nvarid)
-#else
-        ierr = NF_DEF_VAR (nid, "TS", NF_FLOAT, 1, idim2,nvarid)
-#endif
-        ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 22, &
-     &                        "Temperature de surface")
-        ierr = NF_ENDDEF(nid)
-
-
-
-
-  END SUBROUTINE sol_dem_write
-!
-!#########################################################################
-!
   SUBROUTINE gath2cpl(champ_in, champ_out, klon, knon, iim, jjm, knindex)
 
@@ -1806,37 +1711,31 @@
 !#########################################################################
 !
-!
-!#########################################################################
-!
-  SUBROUTINE albsno(agesno,alb_neig_grid)
+  SUBROUTINE albsno(klon, agesno,alb_neig_grid)
   IMPLICIT none
-c
-#include "dimensions.h"
-#include "dimphy.h"
-  INTEGER nvm
-  PARAMETER (nvm=8)
-  REAL veget(klon,nvm)
-  REAL alb_neig(klon)
-  REAL agesno(klon)
-c
-  INTEGER i, nv
-c
-  REAL init(nvm), decay(nvm), as
-  SAVE init, decay
+ 
+  integer :: klon
+  INTEGER, PARAMETER :: nvm = 8
+  REAL, dimension(klon,nvm) :: veget
+  REAL, DIMENSION(klon) :: alb_neig_grid, agesno
+ 
+  INTEGER :: i, nv
+ 
+  REAL, DIMENSION(nvm),SAVE :: init, decay
+  REAL :: as
   DATA init /0.55, 0.14, 0.18, 0.29, 0.15, 0.15, 0.14, 0./
   DATA decay/0.30, 0.67, 0.63, 0.45, 0.40, 0.14, 0.06, 1./
-c
+ 
   veget = 0.
   veget(:,1) = 1.     ! desert partout
   DO i = 1, klon
-    alb_neig(i) = 0.0
+    alb_neig_grid(i) = 0.0
   ENDDO
   DO nv = 1, nvm
     DO i = 1, klon
       as = init(nv)+decay(nv)*EXP(-agesno(i)/5.)
-      alb_neig(i) = alb_neig(i) + veget(i,nv)*as
+      alb_neig_grid(i) = alb_neig_grid(i) + veget(i,nv)*as
     ENDDO
   ENDDO
-c
+ 
   END SUBROUTINE albsno
 !
Index: LMDZ.3.3/branches/rel-LF/libf/phylmd/phyetat0.F
===================================================================
--- LMDZ.3.3/branches/rel-LF/libf/phylmd/phyetat0.F	(revision 110)
+++ LMDZ.3.3/branches/rel-LF/libf/phylmd/phyetat0.F	(revision 112)
@@ -2,5 +2,5 @@
      .            rlat,rlon, pctsrf, tsol,tsoil,deltat,qsol,snow,
      .           albe, evap, rain_fall, snow_fall, solsw, sollw,
-     .           radsol,rugmer,agesno,clesphy0,
+     .           radsol,frugs,agesno,clesphy0,
      .           zmea,zstd,zsig,zgam,zthe,zpic,zval,rugsrel,tabcntr0,
      .           t_ancien,q_ancien,ancien_ok)
@@ -35,5 +35,5 @@
       REAL sollw(klon)
       real solsw(klon)
-      REAL rugmer(klon)
+      REAL frugs(klon,nbsrf)
       REAL agesno(klon)
       REAL zmea(klon)
@@ -640,5 +640,5 @@
               xmax = MAX(evap(i,nsrf),xmax)
            ENDDO
-           PRINT*,'Neige du sol EVAP**:', nsrf, xmin, xmax
+           PRINT*,'evap du sol EVAP**:', nsrf, xmin, xmax
          ENDDO
       ELSE
@@ -660,5 +660,5 @@
             xmax = MAX(evap(i,1),xmax)
          ENDDO
-         PRINT*,'Neige du sol <EVAP>', xmin, xmax
+         PRINT*,'Evap du sol <EVAP>', xmin, xmax
          DO nsrf = 2, nbsrf
          DO i = 1, klon
@@ -793,27 +793,65 @@
       PRINT*,'Rayonnement net au sol radsol:', xmin, xmax
 c
-c Lecture de la longueur de rugosite en mer:
-c
-      ierr = NF_INQ_VARID (nid, "RUGMER", nvarid)
-      IF (ierr.NE.NF_NOERR) THEN
-         PRINT*, 'phyetat0: Le champ <RUGMER> est absent'
-         CALL abort
-      ENDIF
-#ifdef NC_DOUBLE
-      ierr = NF_GET_VAR_DOUBLE(nid, nvarid, rugmer)
-#else
-      ierr = NF_GET_VAR_REAL(nid, nvarid, rugmer)
-#endif
-      IF (ierr.NE.NF_NOERR) THEN
-         PRINT*, 'phyetat0: Lecture echouee pour <RUGMER>'
-         CALL abort
-      ENDIF
-      xmin = 1.0E+20
-      xmax = -1.0E+20
-      DO i = 1, klon
-         xmin = MIN(rugmer(i),xmin)
-         xmax = MAX(rugmer(i),xmax)
-      ENDDO
-      PRINT*,'Rugosite sur la mer rugmer:', xmin, xmax
+c Lecture de la longueur de rugosite 
+c
+c
+      ierr = NF_INQ_VARID (nid, "RUG", nvarid)
+      IF (ierr.NE.NF_NOERR) THEN
+         PRINT*, 'phyetat0: Le champ <RUG> est absent'
+         PRINT*, '          Mais je vais essayer de lire RUG**'
+         DO nsrf = 1, nbsrf
+           IF (nsrf.GT.99) THEN
+             PRINT*, "Trop de sous-mailles"
+             CALL abort
+           ENDIF
+           WRITE(str2,'(i2.2)') nsrf
+           ierr = NF_INQ_VARID (nid, "RUG"//str2, nvarid)
+           IF (ierr.NE.NF_NOERR) THEN
+              PRINT*, "phyetat0: Le champ <RUG"//str2//"> est absent"
+              CALL abort
+           ENDIF
+#ifdef NC_DOUBLE
+           ierr = NF_GET_VAR_DOUBLE(nid, nvarid, frugs(1,nsrf))
+#else
+           ierr = NF_GET_VAR_REAL(nid, nvarid, frugs(1,nsrf))
+#endif
+           IF (ierr.NE.NF_NOERR) THEN
+             PRINT*, "phyetat0: Lecture echouee pour <RUG"//str2//">"
+             CALL abort
+           ENDIF
+           xmin = 1.0E+20
+           xmax = -1.0E+20
+           DO i = 1, klon
+              xmin = MIN(frugs(i,nsrf),xmin)
+              xmax = MAX(frugs(i,nsrf),xmax)
+           ENDDO
+           PRINT*,'evap du sol RUG**:', nsrf, xmin, xmax
+         ENDDO
+      ELSE
+         PRINT*, 'phyetat0: Le champ <RUG> est present'
+         PRINT*, '          J ignore donc les autres RUG**'
+#ifdef NC_DOUBLE
+         ierr = NF_GET_VAR_DOUBLE(nid, nvarid, frugs(1,1))
+#else
+         ierr = NF_GET_VAR_REAL(nid, nvarid, frugs(1,1))
+#endif
+         IF (ierr.NE.NF_NOERR) THEN
+            PRINT*, "phyetat0: Lecture echouee pour <RUG>"
+            CALL abort
+         ENDIF
+         xmin = 1.0E+20
+         xmax = -1.0E+20
+         DO i = 1, klon
+            xmin = MIN(frugs(i,1),xmin)
+            xmax = MAX(frugs(i,1),xmax)
+         ENDDO
+         PRINT*,'Neige du sol <RUG>', xmin, xmax
+         DO nsrf = 2, nbsrf
+         DO i = 1, klon
+            frugs(i,nsrf) = frugs(i,1)
+         ENDDO
+         ENDDO
+      ENDIF
+
 c
 c Lecture de l'age de la neige:
Index: LMDZ.3.3/branches/rel-LF/libf/phylmd/phyredem.F
===================================================================
--- LMDZ.3.3/branches/rel-LF/libf/phylmd/phyredem.F	(revision 110)
+++ LMDZ.3.3/branches/rel-LF/libf/phylmd/phyredem.F	(revision 112)
@@ -3,5 +3,5 @@
      .           albedo, evap, rain_fall, snow_fall,
      .           solsw, sollw,
-     .           radsol,rugmer,agesno,
+     .           radsol,frugs,agesno,
      .           zmea,zstd,zsig,zgam,zthe,zpic,zval,rugsrel,
      .           t_ancien, q_ancien)
@@ -38,5 +38,5 @@
       real sollw(klon)
       REAL radsol(klon)
-      REAL rugmer(klon)
+      REAL frugs(klon,nbsrf)
       REAL agesno(klon)
       REAL zmea(klon)
@@ -464,18 +464,26 @@
 #endif
 c
-      ierr = NF_REDEF (nid)
-#ifdef NC_DOUBLE
-      ierr = NF_DEF_VAR (nid, "RUGMER", NF_DOUBLE, 1, idim2,nvarid)
-#else
-      ierr = NF_DEF_VAR (nid, "RUGMER", NF_FLOAT, 1, idim2,nvarid)
-#endif
-      ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 28,
-     .                        "Longueur de rugosite sur mer")
-      ierr = NF_ENDDEF(nid)
-#ifdef NC_DOUBLE
-      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,rugmer)
-#else
-      ierr = NF_PUT_VAR_REAL (nid,nvarid,rugmer)
-#endif
+      DO nsrf = 1, nbsrf
+        IF (nsrf.LE.99) THEN
+        WRITE(str2,'(i2.2)') nsrf
+        ierr = NF_REDEF (nid)
+#ifdef NC_DOUBLE
+        ierr = NF_DEF_VAR (nid,"RUG"//str2,NF_DOUBLE,1,idim2,nvarid)
+#else
+        ierr = NF_DEF_VAR (nid,"RUG"//str2,NF_FLOAT,1,idim2,nvarid)
+#endif
+        ierr = NF_PUT_ATT_TEXT (nid,nvarid,"title", 23,
+     .                        "rugosite de surface No."//str2)
+        ierr = NF_ENDDEF(nid)
+        ELSE
+        PRINT*, "Trop de sous-mailles"
+        CALL abort
+        ENDIF
+#ifdef NC_DOUBLE
+      ierr = NF_PUT_VAR_DOUBLE (nid,nvarid,frugs(1,nsrf))
+#else
+      ierr = NF_PUT_VAR_REAL (nid,nvarid,frugs(1,nsrf))
+#endif
+      ENDDO
 c
       ierr = NF_REDEF (nid)
Index: LMDZ.3.3/branches/rel-LF/libf/phylmd/physiq.F
===================================================================
--- LMDZ.3.3/branches/rel-LF/libf/phylmd/physiq.F	(revision 110)
+++ LMDZ.3.3/branches/rel-LF/libf/phylmd/physiq.F	(revision 112)
@@ -70,10 +70,12 @@
       logical rnpb
       parameter(rnpb=.true.)
-      PARAMETER (npas=1440)
-      PARAMETER (nexca=48)
-      PARAMETER (itimestep=1800)
+c      PARAMETER (npas=1440)
+c      PARAMETER (nexca=48)
+c      PARAMETER (itimestep=1800)
       EXTERNAL fromcpl, intocpl, inicma
       REAL cpl_sst(iim,jjmp1), cpl_sic(iim,jjmp1)
       REAL cpl_alb_sst(iim,jjmp1), cpl_alb_sic(iim,jjmp1)
+      character *6 ocean
+      parameter (ocean = 'force ')
 c======================================================================
 c ok_ocean indique l'utilisation du modele oceanique "slab ocean",
@@ -669,4 +671,18 @@
          PRINT*, 'La frequence de sortie region est de ', ecrit_reg
          ENDIF
+
+c
+c Initialiser le couplage si necessaire
+c
+      npas = 0
+      nexca = 0
+      if (ocean == 'couple') then
+        npas = itaufin/ iphysiq
+        nexca = 86400 / dtime
+        write(*,*)' ##### Ocean couple #####'
+        write(*,*)' Valeurs des pas de temps'
+        write(*,*)' npas = ', npas
+        write(*,*)' nexca = ', nexca
+      endif        
 c
 c
@@ -1536,5 +1552,6 @@
       DO i = 1, klon
         if (.not. ok_veget) then
-         frugs(i,is_ter) = SQRT(frugs(i,is_ter)**2+rugoro(i)**2)
+          frugs(i,is_ter) = SQRT(frugs(i,is_ter)**2+rugoro(i)**2)
+        endif 
          frugs(i,is_lic) = rugoro(i)
          frugs(i,is_oce) = rugmer(i)
@@ -1563,12 +1580,12 @@
       ENDIF
 
-      CALL clmain(dtime,pctsrf,
+      CALL clmain(dtime,itap,pctsrf,
      e            t_seri,q_seri,u_seri,v_seri,
      e            julien, rmu0, 
-     e            ok_veget, ftsol,
+     e            ok_veget, ocean, npas, nexca, ftsol,
      e            paprs,pplay,radsol, fsnow,fqsol,fevap,falbe,
      e            rain_fall, snow_fall, solsw, sollw, fder,
      e            rlon, rlat, frugs,
-     e            debut, lafin,
+     e            debut, lafin, agesno,
      s            d_t_vdf,d_q_vdf,d_u_vdf,d_v_vdf,d_ts,
      s            fluxt,fluxq,fluxu,fluxv,cdragh,cdragm,rugmer,
@@ -2793,5 +2810,5 @@
      .      falbe, fevap, rain_fall, snow_fall,
      .      solsw, sollw,
-     .      radsol,rugmer,agesno,
+     .      radsol,frugs,agesno,
      .      zmea,zstd,zsig,zgam,zthe,zpic,zval,rugoro,
      .      t_ancien, q_ancien)
Index: LMDZ.3.3/branches/rel-LF/liste_des_sources_f90
===================================================================
--- LMDZ.3.3/branches/rel-LF/liste_des_sources_f90	(revision 110)
+++ LMDZ.3.3/branches/rel-LF/liste_des_sources_f90	(revision 112)
@@ -1,22 +1,41 @@
        USE IOIPSL
 libf/bibio/initdynav.F
-      USE ioipsl
-libf/bibio/writedynav.F
+       USE IOIPSL
+libf/bibio/initfluxsto.F
        USE IOIPSL
 libf/bibio/inithist.F
       USE ioipsl
+libf/bibio/writedynav.F
+      USE ioipsl
 libf/bibio/writehist.F
       USE IOIPSL
+libf/dyn3d/abort_gcm.F
+      USE IOIPSL
+libf/dyn3d/dynredem.F
+      USE ioipsl
+libf/dyn3d/etat0_netcdf.F
+       USE IOIPSL
+libf/dyn3d/fluxstokenc.F
+      USE IOIPSL
+libf/dyn3d/gcm.F
+      USE IOIPSL
 libf/dyn3d/gcmtest.F
-      USE IOIPSL
-libf/dyn3d/abort_gcm.F
       USE ioipsl
 libf/dyn3d/offline.F
       USE ioipsl
 libf/dyn3d/startvar.F
-      USE IOIPSL
-libf/dyn3d/etat0_netcdf.F
+      USE ioipsl
+libf/dyn3d/modif_etat0.F
+       USE IOIPSL
+libf/phylmd/initphysto.F
       USE ioipsl
 libf/phylmd/physiq.F
       USE ioipsl
+libf/phylmd/phystokenc.F
+      USE ioipsl
 libf/phylmd/phytrac.F
+      USE ioipsl
+libf/phylmd/donneesgrads.F
+      USE IOIPSL
+libf/phylmd/1DUTILS.h
+libf/phylmd/interface_surf.F90
Index: LMDZ.3.3/branches/rel-LF/tmp
===================================================================
--- LMDZ.3.3/branches/rel-LF/tmp	(revision 110)
+++ LMDZ.3.3/branches/rel-LF/tmp	(revision 112)
@@ -1,1 +1,2 @@
-fast
+g
+C
Index: LMDZ.3.3/branches/rel-LF/tmp90
===================================================================
--- LMDZ.3.3/branches/rel-LF/tmp90	(revision 110)
+++ LMDZ.3.3/branches/rel-LF/tmp90	(revision 112)
@@ -1,22 +1,41 @@
        USE IOIPSL
 libf/bibio/initdynav.F
-      USE ioipsl
-libf/bibio/writedynav.F
+       USE IOIPSL
+libf/bibio/initfluxsto.F
        USE IOIPSL
 libf/bibio/inithist.F
       USE ioipsl
+libf/bibio/writedynav.F
+      USE ioipsl
 libf/bibio/writehist.F
       USE IOIPSL
+libf/dyn3d/abort_gcm.F
+      USE IOIPSL
+libf/dyn3d/dynredem.F
+      USE ioipsl
+libf/dyn3d/etat0_netcdf.F
+       USE IOIPSL
+libf/dyn3d/fluxstokenc.F
+      USE IOIPSL
+libf/dyn3d/gcm.F
+      USE IOIPSL
 libf/dyn3d/gcmtest.F
-      USE IOIPSL
-libf/dyn3d/abort_gcm.F
       USE ioipsl
 libf/dyn3d/offline.F
       USE ioipsl
 libf/dyn3d/startvar.F
-      USE IOIPSL
-libf/dyn3d/etat0_netcdf.F
+      USE ioipsl
+libf/dyn3d/modif_etat0.F
+       USE IOIPSL
+libf/phylmd/initphysto.F
       USE ioipsl
 libf/phylmd/physiq.F
       USE ioipsl
+libf/phylmd/phystokenc.F
+      USE ioipsl
 libf/phylmd/phytrac.F
+      USE ioipsl
+libf/phylmd/donneesgrads.F
+      USE IOIPSL
+libf/phylmd/1DUTILS.h
+libf/phylmd/interface_surf.F90
