Index: /branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icedyn_rhg_evp.F90
===================================================================
--- /branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icedyn_rhg_evp.F90	(revision 9481)
+++ /branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icedyn_rhg_evp.F90	(revision 9482)
@@ -176,9 +176,4 @@
 
       IF( kt == nit000 .AND. lwp )   WRITE(numout,*) '-- ice_dyn_rhg_evp: EVP sea-ice rheology'
-
-#if defined key_agrif 
-      CALL agrif_interp_lim3( 'U', 0, nn_nevp )   ! First interpolation of coarse values
-      CALL agrif_interp_lim3( 'V', 0, nn_nevp )
-#endif
       !
 !!gm for Clem:  OPTIMIZATION:  I think zfmask can be computed one for all at the initialization....
Index: /branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icestp.F90
===================================================================
--- /branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icestp.F90	(revision 9481)
+++ /branches/2017/dev_merge_2017/NEMOGCM/NEMO/LIM_SRC_3/icestp.F90	(revision 9482)
@@ -110,7 +110,4 @@
          kt_ice = kt                              ! -- Ice model time step
          !
-# if defined key_agrif
-         IF( .NOT. Agrif_Root() )  lim_nbstep = MOD( lim_nbstep, Agrif_irhot() * Agrif_Parent(nn_fsbc) / nn_fsbc ) + 1
-# endif
          u_oce(:,:) = ssu_m(:,:)                  ! -- mean surface ocean current
          v_oce(:,:) = ssv_m(:,:)
@@ -119,5 +116,14 @@
          t_bo(:,:) = ( t_bo(:,:) + rt0 ) * tmask(:,:,1) + rt0 * ( 1._wp - tmask(:,:,1) )
          !
-         CALL store_fields                        ! -- Store now ice values
+         !                          !==  AGRIF Parent to Child  ==!
+#if defined key_agrif
+         !                              ! lim_nbstep ranges from 1 to the nb of child ocean steps inside one parent ice step
+         IF( .NOT. Agrif_Root() )       lim_nbstep = MOD( lim_nbstep, Agrif_irhot() * Agrif_Parent(nn_fsbc) / nn_fsbc ) + 1
+         !                              ! these calls must remain here for restartability purposes
+                                        CALL agrif_interp_lim3( 'T' ) 
+                                        CALL agrif_interp_lim3( 'U' )
+                                        CALL agrif_interp_lim3( 'V' )
+#endif
+                                        CALL store_fields             ! Store now ice values
          !
          !------------------------------------------------!
@@ -131,6 +137,6 @@
          ! --- ice dynamics and advection  --- !
          !-------------------------------------!
-         CALL diag_set0             ! set diag of mass, heat and salt fluxes to 0
-         CALL ice_rst_opn( kt )     ! Open Ice restart file (if necessary) 
+                                        CALL diag_set0                ! set diag of mass, heat and salt fluxes to 0
+                                        CALL ice_rst_opn( kt )        ! Open Ice restart file (if necessary) 
          !
          IF( ln_icedyn .AND. .NOT.lk_c1d )   &
@@ -138,15 +144,11 @@
          !
          !                          !==  lateral boundary conditions  ==!
-#if defined key_agrif
-         IF( .NOT. Agrif_Root()     )   CALL agrif_interp_lim3('T')   ! -- AGRIF 
-#endif
          IF( ln_icethd .AND. ln_bdy )   CALL bdy_ice( kt )            ! -- bdy ice thermo
          !
-         !
          !                          !==  previous lead fraction and ice volume for flux calculations
-         CALL ice_var_glo2eqv            ! h_i and h_s for ice albedo calculation
-         CALL ice_var_agg(1)             ! at_i for coupling 
-         CALL store_fields               ! Store now ice values
-
+                                        CALL ice_var_glo2eqv          ! h_i and h_s for ice albedo calculation
+                                        CALL ice_var_agg(1)           ! at_i for coupling 
+                                        CALL store_fields             ! Store now ice values
+         !
          !------------------------------------------------------!
          ! --- Thermodynamical coupling with the atmosphere --- !
@@ -170,19 +172,21 @@
          !
          IF( ln_icethd )                CALL ice_cor( kt , 2 )        ! -- Corrections
+         !
+         !                          !==  AGRIF Child to Parent  ==!
 # if defined key_agrif
-         IF( .NOT. Agrif_Root() )       CALL agrif_update_lim3( kt )
+                                        CALL agrif_update_lim3( kt )
 # endif
-         CALL ice_var_glo2eqv        ! necessary calls (at least for coupling)
-         CALL ice_var_agg( 2 )       ! necessary calls (at least for coupling)
-                                     !
+                                        CALL ice_var_glo2eqv          ! necessary calls (at least for coupling)
+                                        CALL ice_var_agg( 2 )         ! necessary calls (at least for coupling)
+         !
 !! clem: one should switch the calculation of the fluxes onto the parent grid but the following calls do not work
 !!       moreover it should only be called at the update frequency (as in agrif_lim3_update.F90)
-!!# if defined key_agrif
-!!         IF( .NOT. Agrif_Root() )     CALL Agrif_ChildGrid_To_ParentGrid()
-!!# endif
+!# if defined key_agrif
+!         IF( .NOT. Agrif_Root() )     CALL Agrif_ChildGrid_To_ParentGrid()
+!# endif
                                         CALL ice_update_flx( kt )     ! -- Update ocean surface mass, heat and salt fluxes
-!!# if defined key_agrif
-!!         IF( .NOT. Agrif_Root() )     CALL Agrif_ParentGrid_To_ChildGrid()
-!!# endif
+!# if defined key_agrif
+!         IF( .NOT. Agrif_Root() )     CALL Agrif_ParentGrid_To_ChildGrid()
+!# endif
          IF( ln_icediahsb )             CALL ice_dia( kt )            ! -- Diagnostics outputs 
          !
Index: /branches/2017/dev_merge_2017/NEMOGCM/NEMO/NST_SRC/agrif_lim3_interp.F90
===================================================================
--- /branches/2017/dev_merge_2017/NEMOGCM/NEMO/NST_SRC/agrif_lim3_interp.F90	(revision 9481)
+++ /branches/2017/dev_merge_2017/NEMOGCM/NEMO/NST_SRC/agrif_lim3_interp.F90	(revision 9482)
@@ -53,12 +53,12 @@
       !!-----------------------------------------------------------------------
       !
-      IF( Agrif_Root() .OR. nn_ice==0 )  RETURN   ! do not interpolate if inside Parent domain or if child domain does not have ice
+      IF( Agrif_Root() .OR. nn_ice==0 )  RETURN   ! do not interpolate if inside Parent Grid or if child domain does not have ice
       !
       SELECT CASE( cd_type )
       CASE('U','V')
-         IF( PRESENT( kiter ) ) THEN  ! interpolation at the child sub-time step (only for ice rheology)
+         IF( PRESENT( kiter ) ) THEN  ! interpolation at the child ice sub-time step (only for ice rheology)
             zbeta = ( REAL(lim_nbstep) - REAL(kitermax - kiter) / REAL(kitermax) ) /  &
                &    ( Agrif_Rhot() * REAL(Agrif_Parent(nn_fsbc)) / REAL(nn_fsbc) )
-         ELSE                         ! interpolation at the child time step
+         ELSE                         ! interpolation at the child ice time step
             zbeta = REAL(lim_nbstep) / ( Agrif_Rhot() * REAL(Agrif_Parent(nn_fsbc)) / REAL(nn_fsbc) )
          ENDIF
@@ -204,10 +204,10 @@
                !
                DO jk = 1, nlay_s
-                  e_s(i1:i2,j1:j2,jk,jl) = ptab(:,:,jm) * tmask(i1:i2,j1:j2,1)
+                  e_s(i1:i2,j1:j2,jk,jl) = ptab(i1:i2,j1:j2,jm) * tmask(i1:i2,j1:j2,1)
                   jm = jm + 1
                END DO
                !
                DO jk = 1, nlay_i
-                  e_i(i1:i2,j1:j2,jk,jl) = ptab(:,:,jm) * tmask(i1:i2,j1:j2,1)
+                  e_i(i1:i2,j1:j2,jk,jl) = ptab(i1:i2,j1:j2,jm) * tmask(i1:i2,j1:j2,1)
                   jm = jm + 1
                END DO
@@ -216,11 +216,9 @@
             !
 !!==> clem: this interpolation does not work because it creates negative values, due
-!           to negative coefficients when mixing points (for ex. z7)
+!!          to negative coefficients when mixing points (for ex. z7)
 !!
 !         ELSE                          ! ==> complex interpolation (only one ghost cell available)
 !            !! Use a more complex interpolation since we mix solutions over a couple of grid points
 !            !! it is advised to use it for fields modified by high order schemes (e.g. advection UM5...)
-!            !! clem: for some reason (I don't know why), the following lines do not work 
-!            !        I think there is an issue with Agrif_SpecialValue here (not taken into account properly)
 !            ! record ztab
 !            jm = 1
@@ -358,16 +356,4 @@
 !         ENDIF  ! nbghostcells=1
          
-         ! integrated values
-         vt_i (i1:i2,j1:j2) = SUM(      v_i(i1:i2,j1:j2,:),            dim=3 )
-         vt_s (i1:i2,j1:j2) = SUM(      v_s(i1:i2,j1:j2,:),            dim=3 )
-         at_i (i1:i2,j1:j2) = SUM(      a_i(i1:i2,j1:j2,:),            dim=3 )
-         et_s(i1:i2,j1:j2)  = SUM( SUM( e_s(i1:i2,j1:j2,:,:), dim=4 ), dim=3 )
-         et_i(i1:i2,j1:j2)  = SUM( SUM( e_i(i1:i2,j1:j2,:,:), dim=4 ), dim=3 )
-
-         at_ip(i1:i2,j1:j2) = SUM( a_ip(i1:i2,j1:j2,:), dim=3 ) ! melt ponds
-         vt_ip(i1:i2,j1:j2) = SUM( v_ip(i1:i2,j1:j2,:), dim=3 )
-         !
-         ato_i(i1:i2,j1:j2) = 1._wp - at_i(i1:i2,j1:j2)         ! open water fraction  
-
          DO jl = 1, jpl
             WHERE( tmask(i1:i2,j1:j2,1) == 0._wp )   t_su(i1:i2,j1:j2,jl) = rt0   ! to avoid a division by 0 in sbcblk.F90
Index: /branches/2017/dev_merge_2017/NEMOGCM/NEMO/NST_SRC/agrif_lim3_update.F90
===================================================================
--- /branches/2017/dev_merge_2017/NEMOGCM/NEMO/NST_SRC/agrif_lim3_update.F90	(revision 9481)
+++ /branches/2017/dev_merge_2017/NEMOGCM/NEMO/NST_SRC/agrif_lim3_update.F90	(revision 9482)
@@ -1,3 +1,4 @@
 #define TWO_WAY
+!!#undef TWO_WAY
 
 MODULE agrif_lim3_update
@@ -32,5 +33,5 @@
    PRIVATE
 
-   PUBLIC   agrif_update_lim3   ! called by agrif_user.F90
+   PUBLIC   agrif_update_lim3   ! called by agrif_user.F90 and icestp.F90
 
    !!----------------------------------------------------------------------
@@ -51,9 +52,7 @@
       !!----------------------------------------------------------------------
       !
-      !! clem: I think the update should take place each time the ocean sees the surface forcings
-      !!       (but maybe I am wrong and we should update every rhot time steps) 
-      IF( ( MOD( (kt-nit000)/nn_fsbc + 1, Agrif_irhot() * Agrif_Parent(nn_fsbc) / nn_fsbc ) /=0 ) .AND. (kt /= 0) ) RETURN ! do not update if nb of child time steps differ from time refinement
-                                                                                                                           ! i.e. update only at the parent time step
-      IF( nn_ice == 0 ) RETURN   ! do not update if child domain does not have ice
+      IF( Agrif_Root() .OR. nn_ice == 0 ) RETURN   ! do not update if inside Parent Grid or if child domain does not have ice
+      !
+      IF( ( MOD( (kt-nit000)/nn_fsbc + 1, Agrif_irhot() * Agrif_Parent(nn_fsbc)/nn_fsbc ) /=0 ) .AND. (kt /= 0) ) RETURN   ! update only at the parent ice time step
       !
       Agrif_SpecialValueFineGrid    = -9999.
@@ -148,16 +147,4 @@
          END DO
          !
-         ! integrated values
-         vt_i(i1:i2,j1:j2) = SUM(      v_i(i1:i2,j1:j2,:)           , dim=3 )
-         vt_s(i1:i2,j1:j2) = SUM(      v_s(i1:i2,j1:j2,:)           , dim=3 )
-         at_i(i1:i2,j1:j2) = SUM(      a_i(i1:i2,j1:j2,:)           , dim=3 )
-         et_s(i1:i2,j1:j2) = SUM( SUM( e_s(i1:i2,j1:j2,:,:), dim=4 ), dim=3 )
-         et_i(i1:i2,j1:j2) = SUM( SUM( e_i(i1:i2,j1:j2,:,:), dim=4 ), dim=3 )
-
-         at_ip(i1:i2,j1:j2) = SUM( a_ip(i1:i2,j1:j2,:), dim=3 ) ! melt ponds
-         vt_ip(i1:i2,j1:j2) = SUM( v_ip(i1:i2,j1:j2,:), dim=3 )
-         !
-         ato_i(i1:i2,j1:j2) = 1._wp - at_i(i1:i2,j1:j2)         ! open water fraction  
-
          DO jl = 1, jpl
             WHERE( tmask(i1:i2,j1:j2,1) == 0._wp )   t_su(i1:i2,j1:j2,jl) = rt0   ! to avoid a division by 0 in sbcblk.F90
Index: /branches/2017/dev_merge_2017/NEMOGCM/NEMO/NST_SRC/agrif_user.F90
===================================================================
--- /branches/2017/dev_merge_2017/NEMOGCM/NEMO/NST_SRC/agrif_user.F90	(revision 9481)
+++ /branches/2017/dev_merge_2017/NEMOGCM/NEMO/NST_SRC/agrif_user.F90	(revision 9482)
@@ -307,4 +307,7 @@
    USE agrif_top_update
 #endif
+#if defined key_lim3
+   USE agrif_lim3_update
+#endif
    !
    IMPLICIT NONE
@@ -323,4 +326,8 @@
 !!   CALL Agrif_Update_tke(0)
 
+#if defined key_lim3
+   CALL agrif_update_lim3(0)
+#endif
+   
    CALL Agrif_ChildGrid_To_ParentGrid()
    CALL Agrif_Update_ini()
@@ -497,5 +504,4 @@
    USE agrif_ice
    USE in_out_manager
-   USE agrif_lim3_update
    USE agrif_lim3_interp
    USE lib_mpp
@@ -510,9 +516,9 @@
    ! Controls
 
-   ! clem: For some reason, nn_fsbc(child)/=1 does not work properly (signal is largely degraded by the agrif zoom)
+   ! clem: For some reason, nn_fsbc(child)/=1 does not work properly (signal can be largely degraded by the agrif zoom)
    !          the run must satisfy CFL=Uice/(dx/dt) < 0.6/nn_fsbc(child)
    !          therefore, if nn_fsbc(child)>1 one must reduce the time-step in proportion to nn_fsbc(child), which is not acceptable
-   !       If a solution is found, the following stop could be removed
-   IF( nn_fsbc > 1 )  CALL ctl_stop('nn_fsbc(child) must be set to 1 otherwise agrif and lim3 do not work properly')
+   !       If a solution is found, the following stop could be removed because the rest of the code take nn_fsbc(child) into account
+   IF( nn_fsbc > 1 )  CALL ctl_stop('nn_fsbc(child) must be set to 1 otherwise agrif and sea-ice may not work properly')
 
    ! stop if rhot * nn_fsbc(parent) /= N * nn_fsbc(child) with N being integer
@@ -528,8 +534,4 @@
    lim_nbstep = 0
    
-   ! Update in case 2 ways
-   !----------------------
-   CALL agrif_update_lim3(0)
-
    !
 END SUBROUTINE Agrif_InitValues_cont_lim3
