Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/albedo.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/albedo.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/albedo.F90	(revision 4792)
@@ -212,5 +212,5 @@
       READ  ( numnam_cfg, namsbc_alb, IOSTAT = ios, ERR = 902 )
 902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_alb in configuration namelist', lwp )
-      WRITE ( numond, namsbc_alb )
+      IF(lwm) WRITE ( numond, namsbc_alb )
       !
       IF(lwp) THEN                      ! Control print
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90	(revision 4792)
@@ -40,5 +40,5 @@
       LOGICAL              ::   ln_clim     ! climatology or not (T/F)
       CHARACTER(len = 8)   ::   cltype      ! type of data file 'daily', 'monthly' or yearly'
-      CHARACTER(len = 34)  ::   wname       ! generic name of a NetCDF weights file to be used, blank if not
+      CHARACTER(len = 256) ::   wname       ! generic name of a NetCDF weights file to be used, blank if not
       CHARACTER(len = 34)  ::   vcomp       ! symbolic component name if a vector that needs rotation
       !                                     ! a string starting with "U" or "V" for each component   
@@ -489,9 +489,12 @@
             !       forcing record :    1 
             !                            
-            ztmp = REAL( nday, wp ) / REAL( nyear_len(1), wp ) + 0.5 + REAL( it_offset, wp )
+            ztmp = REAL( nsec_year, wp ) / ( REAL( nyear_len(1), wp ) * rday ) + 0.5 &
+           &       + REAL( it_offset, wp ) / ( REAL( nyear_len(1), wp ) * rday )
             sdjf%nrec_a(1) = 1 + INT( ztmp ) - COUNT((/llbefore/))
             ! swap at the middle of the year
-            IF( llbefore ) THEN   ;   sdjf%nrec_a(2) = nsec1jan000 - NINT(0.5 * rday) * nyear_len(0)
-            ELSE                  ;   sdjf%nrec_a(2) = nsec1jan000 + NINT(0.5 * rday) * nyear_len(1)   
+            IF( llbefore ) THEN   ;   sdjf%nrec_a(2) = nsec1jan000 - (1 - INT(ztmp)) * NINT(0.5 * rday) * nyear_len(0) + &
+                                    & INT(ztmp) * NINT( 0.5 * rday) * nyear_len(1) 
+            ELSE                  ;   sdjf%nrec_a(2) = nsec1jan000 + (1 - INT(ztmp)) * NINT(0.5 * rday) * nyear_len(1) + &
+                                    & INT(ztmp) * INT(rday) * nyear_len(1) + INT(ztmp) * NINT( 0.5 * rday) * nyear_len(2) 
             ENDIF
          ELSE                                    ! no time interpolation
@@ -517,5 +520,6 @@
             !       forcing record :  nmonth 
             !                            
-            ztmp = REAL( nday, wp ) / REAL( nmonth_len(nmonth), wp ) + 0.5 + REAL( it_offset, wp )
+            ztmp = REAL( nsec_month, wp ) / ( REAL( nmonth_len(nmonth), wp ) * rday ) + 0.5 &
+           &       + REAL( it_offset, wp ) / ( REAL( nmonth_len(nmonth), wp ) * rday )
             imth = nmonth + INT( ztmp ) - COUNT((/llbefore/))
             IF( sdjf%cltype == 'monthly' ) THEN   ;   sdjf%nrec_a(1) = 1 + INT( ztmp ) - COUNT((/llbefore/))
@@ -689,5 +693,5 @@
       !!----------------------------------------------------------------------
 #if defined key_bdy
-      USE bdy_oce, ONLY:  dta_global, dta_global2         ! workspace to read in global data arrays
+      USE bdy_oce, ONLY:  dta_global, dta_global_z, dta_global2, dta_global2_z         ! workspace to read in global data arrays
 #endif 
       INTEGER                   , INTENT(in ) ::   num     ! stream number
@@ -706,4 +710,5 @@
       INTEGER                                 ::   ib, ik, ji, jj   ! loop counters
       INTEGER                                 ::   ierr
+      REAL(wp)                                ::   fv ! fillvalue and alternative -ABS(fv)
       REAL(wp), POINTER, DIMENSION(:,:,:)     ::   dta_read    ! work space for global data
       REAL(wp), POINTER, DIMENSION(:,:,:)     ::   dta_read_z  ! work space for global data
@@ -753,5 +758,7 @@
          END SELECT
          CALL iom_getatt(num, '_FillValue', fv, cdvar=clvar )
+#if defined key_bdy
          CALL fld_bdy_interp(dta_read, dta_read_z, map, jpk_bdy, igrd, ibdy, fv, dta)
+#endif
       ELSE ! boundary data assumed to be on model grid
          CALL iom_get ( num, jpdom_unknown, clvar, dta_read(1:ilendta,1:ipj,1:ipk), nrec )                    
@@ -776,4 +783,5 @@
    END SUBROUTINE fld_map
    
+#if defined key_bdy
    SUBROUTINE fld_bdy_interp(dta_read, dta_read_z, map, jpk_bdy, igrd, ibdy, fv, dta)
 
@@ -784,7 +792,5 @@
       !!                boundary data from non-native vertical grid
       !!----------------------------------------------------------------------
-#if defined key_bdy
       USE bdy_oce, ONLY:  idx_bdy         ! indexing for map <-> ij transformation
-#endif 
 
       REAL(wp), POINTER, DIMENSION(:,:,:), INTENT(in )     ::   dta_read    ! work space for global data
@@ -792,6 +798,7 @@
       REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   dta        ! output field on model grid (2 dimensional)
       INTEGER,  DIMENSION(:)    , INTENT(in ) ::   map        ! global-to-local mapping indices
-      INTEGER  , INTENT(in)                   ::   igrd, ib_bdy, jpk_bdy      ! number of levels in bdy data
+      INTEGER  , INTENT(in)                   ::   igrd, ibdy, jpk_bdy      ! number of levels in bdy data
       INTEGER                                 ::   jpkm1_bdy    ! number of levels in bdy data minus 1
+      REAL(wp) , INTENT(in)                                ::   fv ! fillvalue and alternative -ABS(fv)
       !!
       INTEGER                                 ::   ipi        ! length of boundary data on local process
@@ -800,6 +807,7 @@
       INTEGER                                 ::   ilendta    ! length of data in file
       INTEGER                                 ::   ib, ik, ikk! loop counters
+      INTEGER                                 ::   ji, jj ! loop counters
       REAL(wp)                                ::   zl, zi     ! tmp variable for current depth and interpolation factor
-      REAL(wp)                                ::   fv, fv_alt ! fillvalue and alternative -ABS(fv)
+      REAL(wp)                                ::   fv_alt ! fillvalue and alternative -ABS(fv)
       !!---------------------------------------------------------------------
 
@@ -824,5 +832,5 @@
          DO ib = 1, ipi
             DO ik = 1, ipk                      
-               zl =  gdept_1(idx_bdy(ib_bdy)%nbi(ib,igrd),idx_bdy(ib_bdy)%nbj(ib,igrd),ik)   ! if using in step could use fsdept instead of gdept_1?
+               zl =  gdept_0(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd),ik)   ! if using in step could use fsdept instead of gdept_0?
                IF( zl < dta_read_z(map(ib),1,1) ) THEN                                         ! above the first level of external data
                   dta(ib,1,ik) =  dta_read(map(ib),1,1)
@@ -830,5 +838,5 @@
                   dta(ib,1,ik) =  dta_read(map(ib),1,MAXLOC(dta_read_z(map(ib),1,:),1))
                ELSE                                                                          ! inbetween : vertical interpolation between ikk & ikk+1
-                  DO ikk = 1, jpkm1_bdy                                                          ! when  gdept_1(ikk) < zl < gdept_1(ikk+1)
+                  DO ikk = 1, jpkm1_bdy                                                          ! when  gdept_0(ikk) < zl < gdept_0(ikk+1)
                      IF( ( (zl-dta_read_z(map(ib),1,ikk)) * (zl-dta_read_z(map(ib),1,ikk+1)) <= 0._wp)   &
                     &    .AND. (dta_read_z(map(ib),1,ikk+1) /= fv_alt)) THEN
@@ -857,16 +865,16 @@
             ji=map(ib)-(jj-1)*ilendta
             DO ik = 1, ipk                      
-               zl =  gdept_1(idx_bdy(ib_bdy)%nbi(ib,igrd),idx_bdy(ib_bdy)%nbj(ib,igrd),ik)   ! if using in step could use fsdept instead of gdept_1?
+               zl =  gdept_0(idx_bdy(ibdy)%nbi(ib,igrd),idx_bdy(ibdy)%nbj(ib,igrd),ik)   ! if using in step could use fsdept instead of gdept_0?
                IF( zl < dta_read_z(ji,jj,1) ) THEN                                         ! above the first level of external data
-                  dta(ib,1,ik) =  dta_read(ji,jj,1,1)
+                  dta(ib,1,ik) =  dta_read(ji,jj,1)
                ELSEIF( zl > MAXVAL(dta_read_z(ji,ji,:),1) ) THEN                           ! below the last level of external data 
                   dta(ib,1,ik) =  dta_read(ji,jj,MAXLOC(dta_read_z(ji,jj,:),1))
                ELSE                                                                          ! inbetween : vertical interpolation between ikk & ikk+1
-                  DO ikk = 1, jpkm1_bdy                                                          ! when  gdept_1(ikk) < zl < gdept_1(ikk+1)
+                  DO ikk = 1, jpkm1_bdy                                                          ! when  gdept_0(ikk) < zl < gdept_0(ikk+1)
                      IF( ( (zl-dta_read_z(ji,jj,ikk)) * (zl-dta_read_z(ji,jj,ikk+1)) <= 0._wp)   &
                     &    .AND. (dta_read_z(ji,jj,ikk+1) /= fv_alt)) THEN
                         zi = ( zl - dta_read_z(ji,jj,ikk) ) / (dta_read_z(ji,jj,ikk+1)-dta_read_z(ji,jj,ikk))
                         dta(ib,1,ik) = dta_read(ji,jj,ikk) + &
-                      &                ( dta_read(ji,jj,1,ikk+1) -  dta_read(ji,jj,ikk) ) * zi
+                      &                ( dta_read(ji,jj,ikk+1) -  dta_read(ji,jj,ikk) ) * zi
                      ENDIF
                   END DO
@@ -877,4 +885,5 @@
 
    END SUBROUTINE fld_bdy_interp
+#endif
 
    SUBROUTINE fld_rot( kt, sd )
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcana.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcana.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcana.F90	(revision 4792)
@@ -81,5 +81,5 @@
          READ  ( numnam_cfg, namsbc_ana, IOSTAT = ios, ERR = 902 )
 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_ana in configuration namelist', lwp )
-         WRITE ( numond, namsbc_ana )
+         IF(lwm) WRITE ( numond, namsbc_ana )
          !
          IF(lwp) WRITE(numout,*)' '
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcapr.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcapr.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcapr.F90	(revision 4792)
@@ -83,5 +83,5 @@
          READ  ( numnam_cfg, namsbc_apr, IOSTAT = ios, ERR = 902 )
 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_apr in configuration namelist', lwp )
-         WRITE ( numond, namsbc_apr )
+         IF(lwm) WRITE ( numond, namsbc_apr )
          !
          ALLOCATE( sf_apr(1), STAT=ierror )           !* allocate and fill sf_sst (forcing structure) with sn_sst
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_clio.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_clio.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_clio.F90	(revision 4792)
@@ -148,5 +148,5 @@
          READ  ( numnam_cfg, namsbc_clio, IOSTAT = ios, ERR = 902 )
 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_clio in configuration namelist', lwp )
-         WRITE ( numond, namsbc_clio )
+         IF(lwm) WRITE ( numond, namsbc_clio )
 
          ! store namelist information in an array
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90	(revision 4792)
@@ -154,5 +154,5 @@
 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_core in configuration namelist', lwp )
 
-         WRITE ( numond, namsbc_core )
+         IF(lwm) WRITE ( numond, namsbc_core )
          !                                         ! check: do we plan to use ln_dm2dc with non-daily forcing?
          IF( ln_dm2dc .AND. sn_qsr%nfreqh /= 24 )   & 
@@ -563,4 +563,6 @@
       zcoef_dqsb   = rhoa * cpa * Cice
       zcoef_frca   = 1.0  - 0.3
+      ! MV 2014 the proper cloud fraction (mean summer months from the CLIO climato, NH+SH) is 0.19
+      zcoef_frca   = 1.0  - 0.19
 
 !!gm brutal....
@@ -648,10 +650,5 @@
                p_qsr(ji,jj,jl) = zztmp * ( 1. - palb(ji,jj,jl) ) * qsr(ji,jj)
                ! Long  Wave (lw)
-               ! iovino
-               IF( ff(ji,jj) .GT. 0._wp ) THEN
-                  z_qlw(ji,jj,jl) = ( 0.95 * sf(jp_qlw)%fnow(ji,jj,1) - Stef * pst(ji,jj,jl) * zst3 ) * tmask(ji,jj,1)
-               ELSE
-                  z_qlw(ji,jj,jl) = 0.95 * ( sf(jp_qlw)%fnow(ji,jj,1) - Stef * pst(ji,jj,jl) * zst3 ) * tmask(ji,jj,1)
-               ENDIF
+               z_qlw(ji,jj,jl) = 0.95 * ( sf(jp_qlw)%fnow(ji,jj,1) - Stef * pst(ji,jj,jl) * zst3 ) * tmask(ji,jj,1)
                ! lw sensitivity
                z_dqlw(ji,jj,jl) = zcoef_dqlw * zst3                                               
@@ -668,5 +665,11 @@
                   &                         * (  11637800. * EXP( -5897.8 / pst(ji,jj,jl) ) / rhoa - sf(jp_humi)%fnow(ji,jj,1)  ) )
                ! Latent heat sensitivity for ice (Dqla/Dt)
-               p_dqla(ji,jj,jl) = rn_efac * zcoef_dqla * z_wnds_t(ji,jj) / ( zst2 ) * EXP( -5897.8 / pst(ji,jj,jl) )
+               ! MV we also have to cap the sensitivity if the flux is zero
+               IF ( p_qla(ji,jj,jl) .GT. 0.0 ) THEN
+                  p_dqla(ji,jj,jl) = rn_efac * zcoef_dqla * z_wnds_t(ji,jj) / ( zst2 ) * EXP( -5897.8 / pst(ji,jj,jl) )
+               ELSE
+                  p_dqla(ji,jj,jl) = 0.0
+               ENDIF
+                             
                ! Sensible heat sensitivity (Dqsb_ice/Dtn_ice)
                z_dqsb(ji,jj,jl) = zcoef_dqsb * z_wnds_t(ji,jj)
@@ -820,6 +823,9 @@
            sqrt_Cd=kappa/((kappa/sqrt_Cd_n10) - zpsi_m) ; Cd=sqrt_Cd*sqrt_Cd;
          ELSE
-           !! Shifting the wind speed to 10m and neutral stability :
-           U_n10 = dU10*1./(1. + sqrt_Cd_n10/kappa*(log(zu/10.) - zpsi_m)) !  L & Y eq. (9a)
+           !! Shifting the wind speed to 10m and neutral stability :  L & Y eq. (9a)
+           !   In very rare low-wind conditions, the old way of estimating the
+           !   neutral wind speed at 10m leads to a negative value that causes the code
+           !   to crash. To prevent this a threshold of 0.25m/s is now imposed.
+           U_n10 = MAX( 0.25 , dU10/(1. + sqrt_Cd_n10/kappa*(log(zu/10.) - zpsi_m)) )
 
            !! Updating the neutral 10m transfer coefficients :
@@ -956,7 +962,9 @@
          zpsi_m  = psi_m(zeta_u)
          !!
-         !! Shifting the wind speed to 10m and neutral stability : (L & Y eq.(9a))
-!        U_n10 = dU10/(1. + sqrt_Cd_n10/kappa*(log(zu/10.) - psi_m(zeta_u)))
-         U_n10 = dU10/(1. + sqrt_Cd_n10/kappa*(log(zu/10.) - zpsi_m))
+         !! Shifting the wind speed to 10m and neutral stability : L & Y eq.(9a)
+         !   In very rare low-wind conditions, the old way of estimating the
+         !   neutral wind speed at 10m leads to a negative value that causes the code
+         !   to crash. To prevent this a threshold of 0.25m/s is now imposed.
+         U_n10 = MAX( 0.25 , dU10/(1. + sqrt_Cd_n10/kappa*(log(zu/10.) - zpsi_m)) )
          !!
          !! Shifting temperature and humidity at zu :          (L & Y eq. (9b-9c))
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_mfs.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_mfs.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_mfs.F90	(revision 4792)
@@ -141,5 +141,5 @@
          READ  ( numnam_cfg, namsbc_mfs, IOSTAT = ios, ERR = 902 )
 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_mfs in configuration namelist', lwp )
-         WRITE ( numond, namsbc_mfs )
+         IF(lwm) WRITE ( numond, namsbc_mfs )
          !
          ! store namelist information in an array
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90	(revision 4792)
@@ -244,5 +244,5 @@
       READ  ( numnam_cfg, namsbc_cpl, IOSTAT = ios, ERR = 902 )
 902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_cpl in configuration namelist', lwp )
-      WRITE ( numond, namsbc_cpl )
+      IF(lwm) WRITE ( numond, namsbc_cpl )
 
       IF(lwp) THEN                        ! control print
@@ -1374,5 +1374,5 @@
             END SELECT
          CASE( 'mixed oce-ice'        )   
-            ztmp1(:,:) = ( tsn(:,:,1,1) + rt0 ) * zfr_l(:,:) 
+            ztmp1(:,:) = ( tsn(:,:,1,jp_tem) + rt0 ) * zfr_l(:,:) 
             DO jl=1,jpl
                ztmp1(:,:) = ztmp1(:,:) + tn_ice(:,:,jl) * a_i(:,:,jl)
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcflx.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcflx.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcflx.F90	(revision 4792)
@@ -98,5 +98,5 @@
          READ  ( numnam_cfg, namsbc_flx, IOSTAT = ios, ERR = 902 )
 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_flx in configuration namelist', lwp )
-         WRITE ( numond, namsbc_flx ) 
+         IF(lwm) WRITE ( numond, namsbc_flx ) 
          !
          !                                         ! check: do we plan to use ln_dm2dc with non-daily forcing?
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcfwb.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcfwb.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcfwb.F90	(revision 4792)
@@ -57,4 +57,6 @@
       !!                =1 global mean of emp set to zero at each nn_fsbc time step
       !!                =2 annual global mean corrected from previous year
+      !!                =3 global mean of emp set to zero at each nn_fsbc time step
+      !!                   & spread out over erp area depending its sign
       !! Note: if sea ice is embedded it is taken into account when computing the budget 
       !!----------------------------------------------------------------------
@@ -81,5 +83,8 @@
             IF( kn_fwb == 1 )   WRITE(numout,*) '          instantaneously set to zero'
             IF( kn_fwb == 2 )   WRITE(numout,*) '          adjusted from previous year budget'
-         ENDIF
+            IF( kn_fwb == 3 )   WRITE(numout,*) '          fwf set to zero and spread out over erp area'
+         ENDIF
+         !
+         IF( kn_fwb == 3 .AND. nn_sssr /= 2 )   CALL ctl_stop( 'sbc_fwb: nn_fwb = 3 requires nn_sssr = 2, we stop ' )
          !
          area = glob_sum( e1e2t(:,:) )           ! interior global domain surface
@@ -142,6 +147,58 @@
          ENDIF
          !
+      CASE ( 3 )                             !==  global fwf set to zero and spread out over erp area  ==!
+         !
+         IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN
+            ztmsk_pos(:,:) = tmask_i(:,:)                      ! Select <0 and >0 area of erp
+            WHERE( erp < 0._wp )   ztmsk_pos = 0._wp
+            ztmsk_neg(:,:) = tmask_i(:,:) - ztmsk_pos(:,:)
+            !
+            zsurf_neg = glob_sum( e1e2t(:,:)*ztmsk_neg(:,:) )  ! Area filled by <0 and >0 erp 
+            zsurf_pos = glob_sum( e1e2t(:,:)*ztmsk_pos(:,:) )
+            !                                                  ! fwf global mean (excluding ocean to ice/snow exchanges) 
+            z_fwf     = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) - snwice_fmass(:,:) ) ) / area
+            !            
+            IF( z_fwf < 0._wp ) THEN         ! spread out over >0 erp area to increase evaporation
+                zsurf_tospread      = zsurf_pos
+                ztmsk_tospread(:,:) = ztmsk_pos(:,:)
+            ELSE                             ! spread out over <0 erp area to increase precipitation
+                zsurf_tospread      = zsurf_neg
+                ztmsk_tospread(:,:) = ztmsk_neg(:,:)
+            ENDIF
+            !
+            zsum_fwf   = glob_sum( e1e2t(:,:) * z_fwf )         ! fwf global mean over <0 or >0 erp area
+!!gm :  zsum_fwf   = z_fwf * area   ???  it is right?  I think so....
+            z_fwf_nsrf =  zsum_fwf / ( zsurf_tospread + rsmall )
+            !                                                  ! weight to respect erp field 2D structure 
+            zsum_erp   = glob_sum( ztmsk_tospread(:,:) * erp(:,:) * e1e2t(:,:) )
+            z_wgt(:,:) = ztmsk_tospread(:,:) * erp(:,:) / ( zsum_erp + rsmall )
+            !                                                  ! final correction term to apply
+            zerp_cor(:,:) = -1. * z_fwf_nsrf * zsurf_tospread * z_wgt(:,:)
+            !
+!!gm   ===>>>>  lbc_lnk should be useless as all the computation is done over the whole domain !
+            CALL lbc_lnk( zerp_cor, 'T', 1. )
+            !
+            emp(:,:) = emp(:,:) + zerp_cor(:,:)
+            qns(:,:) = qns(:,:) - zerp_cor(:,:) * rcp * sst_m(:,:)  ! account for change to the heat budget due to fw correction
+            erp(:,:) = erp(:,:) + zerp_cor(:,:)
+            !
+            IF( nprint == 1 .AND. lwp ) THEN                   ! control print
+               IF( z_fwf < 0._wp ) THEN
+                  WRITE(numout,*)'   z_fwf < 0'
+                  WRITE(numout,*)'   SUM(erp+)     = ', SUM( ztmsk_tospread(:,:)*erp(:,:)*e1e2t(:,:) )*1.e-9,' Sv'
+               ELSE
+                  WRITE(numout,*)'   z_fwf >= 0'
+                  WRITE(numout,*)'   SUM(erp-)     = ', SUM( ztmsk_tospread(:,:)*erp(:,:)*e1e2t(:,:) )*1.e-9,' Sv'
+               ENDIF
+               WRITE(numout,*)'   SUM(empG)     = ', SUM( z_fwf*e1e2t(:,:) )*1.e-9,' Sv'
+               WRITE(numout,*)'   z_fwf         = ', z_fwf      ,' Kg/m2/s'
+               WRITE(numout,*)'   z_fwf_nsrf    = ', z_fwf_nsrf ,' Kg/m2/s'
+               WRITE(numout,*)'   MIN(zerp_cor) = ', MINVAL(zerp_cor) 
+               WRITE(numout,*)'   MAX(zerp_cor) = ', MAXVAL(zerp_cor) 
+            ENDIF
+         ENDIF
+         !
       CASE DEFAULT                           !==  you should never be there  ==!
-         CALL ctl_stop( 'sbc_fwb : wrong nn_fwb value for the FreshWater Budget correction, choose either 1 or 2' )
+         CALL ctl_stop( 'sbc_fwb : wrong nn_fwb value for the FreshWater Budget correction, choose either 1, 2 or 3' )
          !
       END SELECT
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90	(revision 4792)
@@ -60,6 +60,6 @@
    PUBLIC sbc_ice_cice    ! routine called by sbc
 
-   INTEGER , PARAMETER ::   ji_off = INT ( (jpiglo - nx_global) / 2 )
-   INTEGER , PARAMETER ::   jj_off = INT ( (jpjglo - ny_global) / 2 )
+   INTEGER             ::   ji_off
+   INTEGER             ::   jj_off
 
    INTEGER , PARAMETER ::   jpfld   = 13   ! maximum number of files to read 
@@ -158,4 +158,7 @@
       IF(lwp) WRITE(numout,*)'cice_sbc_init'
 
+      ji_off = INT ( (jpiglo - nx_global) / 2 )
+      jj_off = INT ( (jpjglo - ny_global) / 2 )
+
 ! Initialize CICE
       CALL CICE_Initialize
@@ -220,34 +223,4 @@
          sshb(:,:) = sshb(:,:) - snwice_mass(:,:) * r1_rau0
          !
-         ! Note: Changed the initial values of sshb and sshn=>  need to recompute ssh[u,v,f]_[b,n] 
-         !       which were previously set in domvvl
-         IF ( lk_vvl ) THEN            ! Is this necessary? embd 2 should be restricted to vvl only???
-            DO jj = 1, jpjm1
-               DO ji = 1, jpim1                    ! caution: use of Vector Opt. not possible
-                  zcoefu = 0.5  * umask(ji,jj,1) / ( e1u(ji,jj) * e2u(ji,jj) )
-                  zcoefv = 0.5  * vmask(ji,jj,1) / ( e1v(ji,jj) * e2v(ji,jj) )
-                  zcoeff = 0.25 * umask(ji,jj,1) * umask(ji,jj+1,1)
-                  sshu_b(ji,jj) = zcoefu * ( e1t(ji  ,jj) * e2t(ji  ,jj) * sshb(ji  ,jj)     &
-                     &                     + e1t(ji+1,jj) * e2t(ji+1,jj) * sshb(ji+1,jj) )
-                  sshv_b(ji,jj) = zcoefv * ( e1t(ji,jj  ) * e2t(ji,jj  ) * sshb(ji,jj  )     &
-                     &                     + e1t(ji,jj+1) * e2t(ji,jj+1) * sshb(ji,jj+1) )
-                  sshu_n(ji,jj) = zcoefu * ( e1t(ji  ,jj) * e2t(ji  ,jj) * sshn(ji  ,jj)     &
-                     &                     + e1t(ji+1,jj) * e2t(ji+1,jj) * sshn(ji+1,jj) )
-                  sshv_n(ji,jj) = zcoefv * ( e1t(ji,jj  ) * e2t(ji,jj  ) * sshn(ji,jj  )     &
-                     &                     + e1t(ji,jj+1) * e2t(ji,jj+1) * sshn(ji,jj+1) )
-               END DO
-            END DO
-            CALL lbc_lnk( sshu_b, 'U', 1. )   ;   CALL lbc_lnk( sshu_n, 'U', 1. )
-            CALL lbc_lnk( sshv_b, 'V', 1. )   ;   CALL lbc_lnk( sshv_n, 'V', 1. )
-            DO jj = 1, jpjm1
-               DO ji = 1, jpim1      ! NO Vector Opt.
-                  sshf_n(ji,jj) = 0.5  * umask(ji,jj,1) * umask(ji,jj+1,1)                   &
-                       &               / ( e1f(ji,jj  ) * e2f(ji,jj  ) )                     &
-                       &               * ( e1u(ji,jj  ) * e2u(ji,jj  ) * sshu_n(ji,jj  )     &
-                       &                 + e1u(ji,jj+1) * e2u(ji,jj+1) * sshu_n(ji,jj+1) )
-               END DO
-            END DO
-            CALL lbc_lnk( sshf_n, 'F', 1. )
-          ENDIF
       ENDIF
  
@@ -747,5 +720,5 @@
          READ  ( numnam_cfg, namsbc_cice, IOSTAT = ios, ERR = 902 )
 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_cice in configuration namelist', lwp )
-         WRITE ( numond, namsbc_cice )
+         IF(lwm) WRITE ( numond, namsbc_cice )
 
          ! store namelist information in an array
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_if.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_if.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_if.F90	(revision 4792)
@@ -78,5 +78,5 @@
          READ  ( numnam_cfg, namsbc_iif, IOSTAT = ios, ERR = 902 )
 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_iif in configuration namelist', lwp )
-         WRITE ( numond, namsbc_iif )
+         IF(lwm) WRITE ( numond, namsbc_iif )
 
          ALLOCATE( sf_ice(1), STAT=ierror )
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90	(revision 4792)
@@ -59,4 +59,5 @@
    USE prtctl          ! Print control
    USE lib_fortran     ! 
+   USE cpl_oasis3, ONLY : lk_cpl
 
 #if defined key_bdy 
@@ -68,4 +69,5 @@
 
    PUBLIC sbc_ice_lim  ! routine called by sbcmod.F90
+   PUBLIC lim_prt_state
    
    !! * Substitutions
@@ -133,5 +135,5 @@
       INTEGER, INTENT(in) ::   kblk    ! type of bulk (=3 CLIO, =4 CORE)
       !!
-      INTEGER  ::   jl      ! dummy loop index
+      INTEGER  ::   ji, jj, jl, jk      ! dummy loop index
       REAL(wp) ::   zcoef   ! local scalar
       REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_ice_os, zalb_ice_cs  ! albedo of the ice under overcast/clear sky
@@ -146,4 +148,6 @@
       REAL(wp), POINTER, DIMENSION(:,:) :: z_dqns_ice_all  ! Mean d(qns)/dT over all categories
       REAL(wp), POINTER, DIMENSION(:,:) :: z_dqla_ice_all  ! Mean d(qla)/dT over all categories
+      REAL(wp) ::   ztmelts           ! clem 2014: for HC diags
+      REAL(wp) ::   epsi20 = 1.e-20   !
       !!----------------------------------------------------------------------
 
@@ -152,11 +156,11 @@
       IF( nn_timing == 1 )  CALL timing_start('sbc_ice_lim')
 
-      CALL wrk_alloc( jpi,jpj,jpl, zalb_ice_os, zalb_ice_cs )
-
-#if defined key_coupled
-      IF ( ln_cpl .OR. ln_iceflx_ave .OR. ln_iceflx_linear ) CALL wrk_alloc( jpi,jpj,jpl, zalb_ice)
-      IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) &
-         &   CALL wrk_alloc( jpi,jpj, ztem_ice_all, zalb_ice_all, z_qsr_ice_all, z_qns_ice_all, z_qla_ice_all, z_dqns_ice_all, z_dqla_ice_all)
-#endif
+      CALL wrk_alloc( jpi,jpj,jpl, zalb_ice_os, zalb_ice_cs, zalb_ice )
+
+      IF( lk_cpl ) THEN
+         IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) &
+            &   CALL wrk_alloc( jpi, jpj, ztem_ice_all , zalb_ice_all  , z_qsr_ice_all, z_qns_ice_all,   &
+            &                             z_qla_ice_all, z_dqns_ice_all, z_dqla_ice_all)
+      ENDIF
 
       IF( kt == nit000 ) THEN
@@ -168,5 +172,5 @@
          !
          IF( ln_nicep ) THEN      ! control print at a given point
-            jiindx = 177   ;   jjindx = 112
+            jiindx = 15    ;   jjindx =  44
             IF(lwp) WRITE(numout,*) ' The debugging point is : jiindx : ',jiindx, ' jjindx : ',jjindx
          ENDIF
@@ -176,12 +180,13 @@
       IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !  Ice time-step only  !
          !                                     !----------------------!
-         !                                           !  Bulk Formulea !
+         !                                           !  Bulk Formulae !
          !                                           !----------------!
          !
          u_oce(:,:) = ssu_m(:,:)                     ! mean surface ocean current at ice velocity point
          v_oce(:,:) = ssv_m(:,:)                     ! (C-grid dynamics :  U- & V-points as the ocean)
-         !
-         t_bo(:,:) = tfreez( sss_m ) +  rt0          ! masked sea surface freezing temperature [Kelvin]
-         !                                           ! (set to rt0 over land)
+
+         ! masked sea surface freezing temperature [Kelvin]
+         t_bo(:,:) = ( tfreez( sss_m ) +  rt0 ) * tmask(:,:,1) + rt0 * ( 1. - tmask(:,:,1) )
+
          CALL albedo_ice( t_su, ht_i, ht_s, zalb_ice_cs, zalb_ice_os )  ! ... ice albedo
 
@@ -192,13 +197,13 @@
          IF ( ln_cpl ) zalb_ice (:,:,:) = 0.5 * ( zalb_ice_cs (:,:,:) +  zalb_ice_os (:,:,:) )
          
-#if defined key_coupled
-         IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) THEN
-            !
-            ! Compute mean albedo and temperature
-            zalb_ice_all (:,:) = fice_ice_ave ( zalb_ice (:,:,:) ) 
-            ztem_ice_all (:,:) = fice_ice_ave ( tn_ice   (:,:,:) ) 
-            !
+         IF( lk_cpl ) THEN
+            IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) THEN
+               !
+               ! Compute mean albedo and temperature
+               zalb_ice_all (:,:) = fice_ice_ave ( zalb_ice (:,:,:) ) 
+               ztem_ice_all (:,:) = fice_ice_ave ( tn_ice   (:,:,:) ) 
+               !
+            ENDIF
          ENDIF
-#endif
                                                ! Bulk formulea - provides the following fields:
          ! utau_ice, vtau_ice : surface ice stress                     (U- & V-points)   [N/m2]
@@ -218,5 +223,13 @@
             !         
          CASE( 4 )                                       ! CORE bulk formulation
-            CALL blk_ice_core( t_su , u_ice     , v_ice     , zalb_ice_cs,               &
+            ! MV 2014
+            ! We must account for cloud fraction in the computation of the albedo
+            ! The present ref just uses the clear sky value
+            ! The overcast sky value is 0.06 higher, and polar skies are mostly overcast
+            ! CORE has no cloud fraction, hence we must prescribe it
+            ! Mean summer cloud fraction computed from CLIO = 0.81
+            zalb_ice(:,:,:) = 0.19 * zalb_ice_cs(:,:,:) + 0.81 * zalb_ice_os(:,:,:)
+            ! Following line, we replace zalb_ice_cs by simply zalb_ice
+            CALL blk_ice_core( t_su , u_ice     , v_ice     , zalb_ice   ,               &
                &                      utau_ice  , vtau_ice  , qns_ice    , qsr_ice   ,   &
                &                      qla_ice   , dqns_ice  , dqla_ice   ,               &
@@ -239,5 +252,5 @@
 
          ! Average over all categories
-#if defined key_coupled
+         IF( lk_cpl ) THEN
          IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) THEN
 
@@ -269,5 +282,5 @@
             END IF
          END IF
-#endif
+         ENDIF
          !                                           !----------------------!
          !                                           ! LIM-3  time-stepping !
@@ -285,8 +298,8 @@
          old_smv_i(:,:,:)   = smv_i(:,:,:)     ! salt content
          old_oa_i (:,:,:)   = oa_i (:,:,:)     ! areal age content
-         !
-         old_u_ice(:,:) = u_ice(:,:)
-         old_v_ice(:,:) = v_ice(:,:)
-         !                                           ! intialisation to zero    !!gm is it truly necessary ???
+         old_u_ice(:,:)     = u_ice(:,:)
+         old_v_ice(:,:)     = v_ice(:,:)
+
+         ! trends    !!gm is it truly necessary ???
          d_a_i_thd  (:,:,:)   = 0._wp   ;   d_a_i_trp  (:,:,:)   = 0._wp
          d_v_i_thd  (:,:,:)   = 0._wp   ;   d_v_i_trp  (:,:,:)   = 0._wp
@@ -296,33 +309,40 @@
          d_smv_i_thd(:,:,:)   = 0._wp   ;   d_smv_i_trp(:,:,:)   = 0._wp
          d_oa_i_thd (:,:,:)   = 0._wp   ;   d_oa_i_trp (:,:,:)   = 0._wp
-         !
-         d_u_ice_dyn(:,:) = 0._wp
-         d_v_ice_dyn(:,:) = 0._wp
-         !
-         sfx    (:,:) = 0._wp   ;   sfx_thd  (:,:) = 0._wp
-         sfx_bri(:,:) = 0._wp   ;   sfx_mec  (:,:) = 0._wp   ;   sfx_res  (:,:) = 0._wp
-         fhbri  (:,:) = 0._wp   ;   fheat_mec(:,:) = 0._wp   ;   fheat_res(:,:) = 0._wp
-         fhmec  (:,:) = 0._wp   ;   
-         fmmec  (:,:) = 0._wp
-         fmmflx (:,:) = 0._wp     
-         focea2D(:,:) = 0._wp
-         fsup2D (:,:) = 0._wp
-
-         ! used in limthd.F90
-         rdvosif(:,:) = 0._wp   ! variation of ice volume at surface
-         rdvobif(:,:) = 0._wp   ! variation of ice volume at bottom
-         fdvolif(:,:) = 0._wp   ! total variation of ice volume
-         rdvonif(:,:) = 0._wp   ! lateral variation of ice volume
-         fstric (:,:) = 0._wp   ! part of solar radiation transmitted through the ice
-         ffltbif(:,:) = 0._wp   ! linked with fstric
-         qfvbq  (:,:) = 0._wp   ! linked with fstric
-         rdm_snw(:,:) = 0._wp   ! variation of snow mass per unit area
-         rdm_ice(:,:) = 0._wp   ! variation of ice mass per unit area
-         hicifp (:,:) = 0._wp   ! daily thermodynamic ice production. 
-         !
-         diag_sni_gr(:,:) = 0._wp   ;   diag_lat_gr(:,:) = 0._wp
-         diag_bot_gr(:,:) = 0._wp   ;   diag_dyn_gr(:,:) = 0._wp
-         diag_bot_me(:,:) = 0._wp   ;   diag_sur_me(:,:) = 0._wp
-         diag_res_pr(:,:) = 0._wp   ;   diag_trp_vi(:,:) = 0._wp
+         d_u_ice_dyn(:,:)     = 0._wp   ;   d_v_ice_dyn(:,:)     = 0._wp
+
+         ! salt, heat and mass fluxes
+         sfx    (:,:) = 0._wp   ;
+         sfx_bri(:,:) = 0._wp   ;   sfx_dyn(:,:) = 0._wp 
+         sfx_sni(:,:) = 0._wp   ;   sfx_opw(:,:) = 0._wp
+         sfx_bog(:,:) = 0._wp   ;   sfx_dyn(:,:) = 0._wp
+         sfx_bom(:,:) = 0._wp   ;   sfx_sum(:,:) = 0._wp
+         sfx_res(:,:) = 0._wp
+
+         wfx_snw(:,:) = 0._wp   ;   wfx_ice(:,:) = 0._wp
+         wfx_sni(:,:) = 0._wp   ;   wfx_opw(:,:) = 0._wp
+         wfx_bog(:,:) = 0._wp   ;   wfx_dyn(:,:) = 0._wp
+         wfx_bom(:,:) = 0._wp   ;   wfx_sum(:,:) = 0._wp
+         wfx_res(:,:) = 0._wp   ;   wfx_sub(:,:) = 0._wp
+         wfx_spr(:,:) = 0._wp   ;   
+
+         hfx_in (:,:) = 0._wp   ;   hfx_out(:,:) = 0._wp
+         hfx_thd(:,:) = 0._wp   ;   
+         hfx_snw(:,:) = 0._wp   ;   hfx_opw(:,:) = 0._wp
+         hfx_bog(:,:) = 0._wp   ;   hfx_dyn(:,:) = 0._wp
+         hfx_bom(:,:) = 0._wp   ;   hfx_sum(:,:) = 0._wp
+         hfx_res(:,:) = 0._wp   ;   hfx_sub(:,:) = 0._wp
+         hfx_spr(:,:) = 0._wp   ;   hfx_dif(:,:) = 0._wp 
+         hfx_err(:,:) = 0._wp   ;   hfx_err_rem(:,:) = 0._wp
+
+         !
+         fhld  (:,:)    = 0._wp 
+         fmmflx(:,:)    = 0._wp     
+         ! part of solar radiation transmitted through the ice
+         ftr_ice(:,:,:) = 0._wp
+
+         ! diags
+         diag_trp_vi  (:,:) = 0._wp  ; diag_trp_vs(:,:) = 0._wp  ;  diag_trp_ei(:,:) = 0._wp  ;  diag_trp_es(:,:) = 0._wp
+         diag_heat_dhc(:,:) = 0._wp  
+
          ! dynamical invariants
          delta_i(:,:) = 0._wp       ;   divu_i(:,:) = 0._wp       ;   shear_i(:,:) = 0._wp
@@ -375,5 +395,4 @@
                           zcoef = rdt_ice /rday           !  Ice natural aging
                           oa_i(:,:,:) = oa_i(:,:,:) + a_i(:,:,:) * zcoef
-                          CALL lim_var_glo2eqv            ! this CALL is maybe not necessary (Martin)
          IF( ln_nicep )   CALL lim_prt_state( kt, jiindx, jjindx, 1, ' - ice thermodyn. - ' )   ! control print
                           CALL lim_itd_th( kt )           !  Remap ice categories, lateral accretion  !
@@ -391,7 +410,7 @@
          !                                           ! Diagnostics and outputs 
          IF (ln_limdiaout) CALL lim_diahsb
-!clem # if ! defined key_iomput
+
                           CALL lim_wri( 1  )              ! Ice outputs 
-!clem # endif
+
          IF( kt == nit000 .AND. ln_rstart )   &
             &             CALL iom_close( numrir )        ! clem: close input ice restart file
@@ -413,12 +432,11 @@
       
 !!gm   remark, the ocean-ice stress is not saved in ice diag call above .....  find a solution!!!
-      !
-      CALL wrk_dealloc( jpi,jpj,jpl, zalb_ice_os, zalb_ice_cs )
-
-#if defined key_coupled
-      IF ( ln_cpl .OR. ln_iceflx_ave .OR. ln_iceflx_linear ) CALL wrk_dealloc( jpi,jpj,jpl, zalb_ice)
-      IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) &
-         &    CALL wrk_dealloc( jpi,jpj, ztem_ice_all, zalb_ice_all, z_qsr_ice_all, z_qns_ice_all, z_qla_ice_all, z_dqns_ice_all, z_dqla_ice_all)
-#endif
+      CALL wrk_dealloc( jpi,jpj,jpl, zalb_ice_os, zalb_ice_cs, zalb_ice )
+
+      IF( lk_cpl ) THEN
+         IF ( ln_iceflx_ave .OR. ln_iceflx_linear ) &
+            &    CALL wrk_dealloc( jpi, jpj, ztem_ice_all , zalb_ice_all , z_qsr_ice_all, z_qns_ice_all,   &
+            &                                z_qla_ice_all, z_dqns_ice_all, z_dqla_ice_all)
+      ENDIF
       !
       IF( nn_timing == 1 )  CALL timing_stop('sbc_ice_lim')
@@ -534,5 +552,4 @@
 !                 WRITE(numout,*) ' sst                  : ', sst_m(ji,jj)
 !                 WRITE(numout,*) ' sss                  : ', sss_m(ji,jj)
-!                 WRITE(numout,*) ' s_i_newice           : ', s_i_newice(ji,jj,1:jpl)
 !                 WRITE(numout,*) 
                   inb_alp(ialert_id) = inb_alp(ialert_id) + 1
@@ -591,17 +608,4 @@
                !WRITE(numout,*) ' sst       : ', sst_m(ji,jj)
                !WRITE(numout,*) ' sss       : ', sss_m(ji,jj)
-               !WRITE(numout,*) ' qcmif     : ', qcmif(ji,jj)
-               !WRITE(numout,*) ' qldif     : ', qldif(ji,jj)
-               !WRITE(numout,*) ' qcmif     : ', qcmif(ji,jj) / rdt_ice
-               !WRITE(numout,*) ' qldif     : ', qldif(ji,jj) / rdt_ice
-               !WRITE(numout,*) ' qfvbq     : ', qfvbq(ji,jj)
-               !WRITE(numout,*) ' qdtcn     : ', qdtcn(ji,jj)
-               !WRITE(numout,*) ' qfvbq / dt: ', qfvbq(ji,jj) / rdt_ice
-               !WRITE(numout,*) ' qdtcn / dt: ', qdtcn(ji,jj) / rdt_ice
-               !WRITE(numout,*) ' fdtcn     : ', fdtcn(ji,jj) 
-               !WRITE(numout,*) ' fhmec     : ', fhmec(ji,jj) 
-               !WRITE(numout,*) ' fheat_mec : ', fheat_mec(ji,jj) 
-               !WRITE(numout,*) ' fheat_res : ', fheat_res(ji,jj) 
-               !WRITE(numout,*) ' fhbri     : ', fhbri(ji,jj) 
                !
                !CALL lim_prt_state( kt, ji, jj, 2, '   ')
@@ -790,14 +794,8 @@
                WRITE(numout,*) ' - Heat / FW fluxes '
                WRITE(numout,*) '   ~~~~~~~~~~~~~~~~ '
-               WRITE(numout,*) ' emp        : ', emp      (ji,jj)
-               WRITE(numout,*) ' sfx        : ', sfx      (ji,jj)
-               WRITE(numout,*) ' sfx_thd    : ', sfx_thd(ji,jj)
-               WRITE(numout,*) ' sfx_bri    : ', sfx_bri  (ji,jj)
-               WRITE(numout,*) ' sfx_mec    : ', sfx_mec  (ji,jj)
-               WRITE(numout,*) ' sfx_res    : ', sfx_res(ji,jj)
-               WRITE(numout,*) ' fmmec      : ', fmmec    (ji,jj)
-               WRITE(numout,*) ' fhmec      : ', fhmec    (ji,jj)
-               WRITE(numout,*) ' fhbri      : ', fhbri    (ji,jj)
-               WRITE(numout,*) ' fheat_mec  : ', fheat_mec(ji,jj)
+               WRITE(numout,*) ' - Heat fluxes in and out the ice ***'
+               WRITE(numout,*) ' qsr_ini       : ', pfrld(ji,jj) * qsr(ji,jj) + SUM( old_a_i(ji,jj,:) * qsr_ice(ji,jj,:) )
+               WRITE(numout,*) ' qns_ini       : ', pfrld(ji,jj) * qns(ji,jj) + SUM( old_a_i(ji,jj,:) * qns_ice(ji,jj,:) )
+               WRITE(numout,*)
                WRITE(numout,*) 
                WRITE(numout,*) ' sst        : ', sst_m(ji,jj)  
@@ -829,21 +827,27 @@
                WRITE(numout,*) ' qsr       : ', qsr(ji,jj)
                WRITE(numout,*) ' qns       : ', qns(ji,jj)
-               WRITE(numout,*) ' fdtcn     : ', fdtcn(ji,jj)
-               WRITE(numout,*) ' qcmif     : ', qcmif(ji,jj) * r1_rdtice
-               WRITE(numout,*) ' qldif     : ', qldif(ji,jj) * r1_rdtice
+               WRITE(numout,*)
+               WRITE(numout,*) ' hfx_mass     : ', hfx_thd(ji,jj) + hfx_dyn(ji,jj) + hfx_snw(ji,jj) + hfx_res(ji,jj)
+               WRITE(numout,*) ' hfx_in       : ', hfx_in(ji,jj)
+               WRITE(numout,*) ' hfx_out      : ', hfx_out(ji,jj)
+               WRITE(numout,*) ' dhc          : ', diag_heat_dhc(ji,jj)              
+               WRITE(numout,*)
+               WRITE(numout,*) ' hfx_dyn      : ', hfx_dyn(ji,jj)
+               WRITE(numout,*) ' hfx_thd      : ', hfx_thd(ji,jj)
+               WRITE(numout,*) ' hfx_res      : ', hfx_res(ji,jj)
+               WRITE(numout,*) ' fhtur        : ', fhtur(ji,jj) 
+               WRITE(numout,*) ' qlead        : ', qlead(ji,jj) * r1_rdtice
                WRITE(numout,*)
                WRITE(numout,*) ' - Salt fluxes at bottom interface ***'
                WRITE(numout,*) ' emp       : ', emp    (ji,jj)
-               WRITE(numout,*) ' sfx_bri   : ', sfx_bri(ji,jj)
                WRITE(numout,*) ' sfx       : ', sfx    (ji,jj)
                WRITE(numout,*) ' sfx_res   : ', sfx_res(ji,jj)
-               WRITE(numout,*) ' sfx_mec   : ', sfx_mec(ji,jj)
-               WRITE(numout,*) ' - Heat fluxes at bottom interface ***'
-               WRITE(numout,*) ' fheat_res : ', fheat_res(ji,jj)
+               WRITE(numout,*) ' sfx_bri   : ', sfx_bri(ji,jj)
+               WRITE(numout,*) ' sfx_dyn   : ', sfx_dyn(ji,jj)
                WRITE(numout,*)
                WRITE(numout,*) ' - Momentum fluxes '
                WRITE(numout,*) ' utau      : ', utau(ji,jj) 
                WRITE(numout,*) ' vtau      : ', vtau(ji,jj)
-            ENDIF
+            ENDIF 
             WRITE(numout,*) ' '
             !
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim_2.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim_2.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim_2.F90	(revision 4792)
@@ -53,4 +53,8 @@
    USE agrif_lim2_update
 # endif
+
+#if defined key_bdy 
+   USE bdyice_lim       ! unstructured open boundary data  (bdy_ice_lim routine)
+#endif
 
    IMPLICIT NONE
@@ -205,4 +209,7 @@
                            CALL lim_trp_2      ( kt )      ! Ice transport   ( Advection/diffusion )
            IF( ln_limdmp ) CALL lim_dmp_2      ( kt )      ! Ice damping 
+#if defined key_bdy
+                           CALL bdy_ice_lim( kt ) ! bdy ice thermo
+#endif
          END IF
 #if defined key_coupled
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90	(revision 4792)
@@ -101,5 +101,5 @@
       READ  ( numnam_cfg, namsbc, IOSTAT = ios, ERR = 902 )
 902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc in configuration namelist', lwp )
-      WRITE ( numond, namsbc )
+      IF(lwm) WRITE ( numond, namsbc )
 
       !                          ! overwrite namelist parameter using CPP key information
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcrnf.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcrnf.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcrnf.F90	(revision 4792)
@@ -263,5 +263,5 @@
       READ  ( numnam_cfg, namsbc_rnf, IOSTAT = ios, ERR = 902 )
 902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_rnf in configuration namelist', lwp )
-      WRITE ( numond, namsbc_rnf )
+      IF(lwm) WRITE ( numond, namsbc_rnf )
       !
       !                                         ! Control print
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcssr.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcssr.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcssr.F90	(revision 4792)
@@ -174,5 +174,5 @@
       READ  ( numnam_cfg, namsbc_ssr, IOSTAT = ios, ERR = 902 )
 902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_ssr in configuration namelist', lwp )
-      WRITE ( numond, namsbc_ssr )
+      IF(lwm) WRITE ( numond, namsbc_ssr )
 
       IF(lwp) THEN                 !* control print
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90	(revision 4792)
@@ -90,5 +90,5 @@
          READ  ( numnam_cfg, namsbc_wave, IOSTAT = ios, ERR = 902 )
 902      IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_wave in configuration namelist', lwp )
-         WRITE ( numond, namsbc_wave )
+         IF(lwm) WRITE ( numond, namsbc_wave )
          !
 
Index: branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/tideini.F90
===================================================================
--- branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/tideini.F90	(revision 4694)
+++ branches/2014/dev_r4621_NOC4_BDY_VERT_INTERP/NEMOGCM/NEMO/OPA_SRC/SBC/tideini.F90	(revision 4792)
@@ -72,5 +72,5 @@
        READ  ( numnam_cfg, nam_tide, IOSTAT = ios, ERR = 902 )
 902    IF( ios /= 0 ) CALL ctl_nam ( ios , 'nam_tide in configuration namelist', lwp )
-       WRITE ( numond, nam_tide )
+       IF(lwm) WRITE ( numond, nam_tide )
        !
        nb_harmo=0
