Index: /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/CONFIG/SHARED/namelist_ref
===================================================================
--- /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/CONFIG/SHARED/namelist_ref	(revision 5822)
+++ /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/CONFIG/SHARED/namelist_ref	(revision 5823)
@@ -474,5 +474,5 @@
 !-----------------------------------------------------------------------
    rn_fiscpl = 43800    ! (number of time step) conservation period (maybe should be fix to the coupling frequencey of restart frequency)
-   ln_hfb    = .false.  ! activate conservation module (conservation exact after a time of rn_fiscpl)
+   ln_hsb    = .false.  ! activate conservation module (conservation exact after a time of rn_fiscpl)
 /
 !-----------------------------------------------------------------------
@@ -810,12 +810,4 @@
    rn_smsh          =     1.    !  Smagorinsky diffusivity: = 0 - use only sheer
    rn_aht_m         =  2000.    !  upper limit or stability criteria for lateral eddy diffusivity (m2/s)
-/
-!-----------------------------------------------------------------------
-&namtra_dmpfile    !   tracer: T & S newtonian damping
-!-----------------------------------------------------------------------
-!          !  file name                            ! frequency (hours) ! variable  ! time interp. !  clim  ! 'yearly'/ ! weights  ! rotation ! land/sea mask !
-!          !                                       !  (if <0  months)  !   name    !   (logical)  !  (T/F) ! 'monthly' ! filename ! pairing  ! filename      !
-   sn_dmpt  = 'resto',         -1        ,'Tinit' ,    .true.    , .true. , 'yearly'   , ''       ,   ''    ,    ''
-   sn_dmps  = 'resto',         -1        ,'Sinit' ,    .true.    , .true. , 'yearly'   , ''       ,   ''    ,    ''
 /
 !-----------------------------------------------------------------------
Index: /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90
===================================================================
--- /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90	(revision 5822)
+++ /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90	(revision 5823)
@@ -140,4 +140,5 @@
       ! 2 -  Content variations !
       ! ------------------------ !
+      ! glob_sum_full is needed because you keep the full interior domain to compute the sum (iscpl)
       zdiff_v2 = 0._wp
       zdiff_hc = 0._wp
@@ -160,6 +161,6 @@
             z2d1(:,:) = surf(:,:) * ( tsn(:,:,1,jp_sal) * sshn(:,:) - ssh_sc_loc_ini(:,:) ) 
          END IF
-         z_ssh_hc = glob_sum( z2d0 ) 
-         z_ssh_sc = glob_sum( z2d1 ) 
+         z_ssh_hc = glob_sum_full( z2d0 ) 
+         z_ssh_sc = glob_sum_full( z2d1 ) 
       ENDIF
 
@@ -201,11 +202,4 @@
 !      ENDIF
 !!gm end
-      IF ( lk_vvl ) THEN
-         IF (lwp) PRINT *, 'cons heat : ', kt, zdiff_hc / zvol_tot, zdiff_sc / zvol_tot 
-         IF (lwp) PRINT *, 'cons volu : ', kt, zdiff_v2 * 1.e-9      
-      ELSE
-         IF (lwp) PRINT *, 'cons heat : ', kt, zdiff_hc1 * 1.e-20 * rau0 * rcp, zdiff_sc1 * 1.e-9
-         IF (lwp) PRINT *, 'cons vol  : ', kt, zdiff_v1 * 1.e-9      
-      END IF
       IF( lk_vvl ) THEN
         CALL iom_put( 'bgtemper' , zdiff_hc / zvol_tot )              ! Temperature variation (C) 
@@ -269,5 +263,5 @@
               CALL iom_get( numror, 'frc_wn_s', frc_wn_s )
            ENDIF
-           CALL iom_get( numror, jpdom_autoglo, 'surf_ini', surf_ini )
+           CALL iom_get( numror, jpdom_autoglo, 'surf_ini', surf_ini ) ! ice sheet coupling
            CALL iom_get( numror, jpdom_autoglo, 'ssh_ini', ssh_ini )
            CALL iom_get( numror, jpdom_autoglo, 'e3t_ini', e3t_ini )
@@ -323,5 +317,5 @@
            CALL iom_rstput( kt, nitrst, numrow, 'frc_wn_s', frc_wn_s )
         ENDIF
-        CALL iom_rstput( kt, nitrst, numrow, 'surf_ini', surf_ini )
+        CALL iom_rstput( kt, nitrst, numrow, 'surf_ini', surf_ini )      ! ice sheet coupling
         CALL iom_rstput( kt, nitrst, numrow, 'ssh_ini', ssh_ini )
         CALL iom_rstput( kt, nitrst, numrow, 'e3t_ini', e3t_ini )
Index: /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90
===================================================================
--- /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90	(revision 5822)
+++ /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90	(revision 5823)
@@ -529,6 +529,7 @@
                WHERE( bathy(:,:) <= 0._wp )  risfdep(:,:) = 0._wp
 
-               ! set grounded point to 0 (treshold at 1cm, have to be update after first coupling experience)
-               WHERE (bathy(:,:) .LE. risfdep(:,:)+1e-2 )
+               ! set grounded point to 0 
+               ! (a treshold could be set here if needed, or set it offline based on the grounded fraction)
+               WHERE ( bathy(:,:) .LE. risfdep(:,:) )
                   misfdep(:,:) = 0 ; risfdep(:,:) = 0._wp
                   mbathy (:,:) = 0 ; bathy  (:,:) = 0._wp
@@ -575,5 +576,5 @@
          ! patch to avoid case bathy = ice shelf draft and bathy between 0 and zhmin
          IF ( ln_isfcav ) THEN
-            WHERE (bathy == risfdep)
+            WHERE ( bathy == risfdep )
                bathy   = 0.0_wp ; risfdep = 0.0_wp
             END WHERE
Index: /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/iscplhsb.F90
===================================================================
--- /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/iscplhsb.F90	(revision 5822)
+++ /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/iscplhsb.F90	(revision 5823)
@@ -103,5 +103,5 @@
                   ! heat diff
                   zdtem(ji,jj) = tsn(ji,jj,jk,jp_tem) * fse3t_n(ji,jj,jk) *  tmask  (ji,jj,jk)   &
-                            - tsb(ji,jj,jk,jp_tem) * pe3t_b (ji,jj,jk) * ptmask_b(ji,jj,jk)
+                               - tsb(ji,jj,jk,jp_tem) * pe3t_b (ji,jj,jk) * ptmask_b(ji,jj,jk)
                   ! salt diff
                   zdsal(ji,jj) = tsn(ji,jj,jk,jp_sal) * fse3t_n(ji,jj,jk) *  tmask  (ji,jj,jk)   &
Index: /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/iscplini.F90
===================================================================
--- /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/iscplini.F90	(revision 5822)
+++ /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/iscplini.F90	(revision 5823)
@@ -13,8 +13,7 @@
    USE dom_oce         ! ocean space and time domain
    USE oce             ! global tra/dyn variable
-   USE in_out_manager  ! I/O manager
    USE lib_mpp         ! MPP library
    USE lib_fortran     ! MPP library
-   USE iom
+   USE in_out_manager  ! I/O manager
 
    IMPLICIT NONE
@@ -45,5 +44,5 @@
       ALLOCATE( htsc_iscpl(jpi,jpj,jpk,jpts) , hdiv_iscpl(jpi,jpj,jpk) , STAT=iscpl_alloc )
          !
-      IF( lk_mpp              )   CALL mpp_sum ( iscpl_alloc )
+      IF( lk_mpp          )   CALL mpp_sum ( iscpl_alloc )
       IF( iscpl_alloc > 0 )   CALL ctl_warn('iscpl_alloc: allocation of arrays failed')
    END FUNCTION iscpl_alloc
Index: /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/iscplrst.F90
===================================================================
--- /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/iscplrst.F90	(revision 5822)
+++ /branches/NERC/dev_r5589_is_oce_cpl/NEMOGCM/NEMO/OPA_SRC/DOM/iscplrst.F90	(revision 5823)
@@ -21,14 +21,12 @@
    USE lib_fortran     ! MPP library
    USE wrk_nemo        ! Memory allocation
-   USE lbclnk
-   USE domngb
-   USE sbc_ice, ONLY : lk_lim3
-   USE iscplini
-   USE iscplhsb
+   USE lbclnk          ! communication
+   USE iscplini        ! ice sheet coupling: initialisation
+   USE iscplhsb        ! ice sheet coupling: conservation
 
    IMPLICIT NONE
    PRIVATE
    
-   PUBLIC   iscpl_stp       
+   PUBLIC   iscpl_stp          ! step management 
    PUBLIC   iscpl_rst_interpol ! routine to wet and dry
    !!
@@ -55,5 +53,5 @@
       REAL(wp), DIMENSION(:,:,:), POINTER :: ze3t_b  , ze3u_b  , ze3v_b  
       REAL(wp), DIMENSION(:,:,:), POINTER :: zdepw_b
-!!----------------------------------------------------------------------
+      !!----------------------------------------------------------------------
       INTEGER  ::   inum0
       CHARACTER(20) :: cfile
@@ -141,5 +139,5 @@
       REAL(wp), DIMENSION(:,:,:  ), INTENT(in ) :: ptmask_b, pumask_b, pvmask_b    !! mask before
       REAL(wp), DIMENSION(:,:,:  ), INTENT(in ) :: pe3t_b  , pe3u_b  , pe3v_b      !! scale factor before
-      REAL(wp), DIMENSION(:,:,:  ), INTENT(in ) :: pdepw_b
+      REAL(wp), DIMENSION(:,:,:  ), INTENT(in ) :: pdepw_b                         !! depth w before
       REAL(wp), DIMENSION(:,:    ), INTENT(in ) :: psmask_b                        !! mask before
       !!
@@ -157,5 +155,6 @@
       REAL(wp), DIMENSION(:,:,:  ), POINTER :: zwmaskn, zwmaskb, ztmp3d
       REAL(wp), DIMENSION(:,:,:,:), POINTER :: zts0
-      !!
+
+      !! allocate variables
       CALL wrk_alloc(jpi,jpj,jpk,2, zts0                                   )
       CALL wrk_alloc(jpi,jpj,jpk,   ztmask0, ztmask1 , ztrp, ztmp3d         ) 
@@ -165,4 +164,5 @@
       CALL wrk_alloc(jpi,jpj,       zbub   , zbvb    , zbun , zbvn         ) 
       CALL wrk_alloc(jpi,jpj,       zssh0  , zssh1, hu1, hv1               ) 
+
       !! mask value to be sure
       tsn(:,:,:,jp_tem) = tsn(:,:,:,jp_tem) * ptmask_b(:,:,:)
@@ -182,5 +182,5 @@
       zsmask0(:,:) = psmask_b(:,:)
       zsmask1(:,:) = psmask_b(:,:)
-      DO iz = 1,10    ! need to be tuned (configuration dependent)
+      DO iz = 1,10    ! need to be tuned (configuration dependent) (OK for ISOMIP+)
          zdsmask(:,:) = ssmask(:,:)-zsmask0(:,:)
          DO ji = 2,jpi-1
@@ -189,15 +189,15 @@
                jjp1=jj+1; jjm1=jj-1;
                summsk=(zsmask0(jip1,jj)+zsmask0(jim1,jj)+zsmask0(ji,jjp1)+zsmask0(ji,jjm1))
-               IF (zdsmask(ji,jj)==1 .AND. summsk .NE. 0) THEN
+               IF (zdsmask(ji,jj)==1._wp .AND. summsk .NE. 0._wp) THEN
                   sshn(ji,jj)=( zssh0(jip1,jj)*zsmask0(jip1,jj)     &
                   &           + zssh0(jim1,jj)*zsmask0(jim1,jj)     &
                   &           + zssh0(ji,jjp1)*zsmask0(ji,jjp1)     &
                   &           + zssh0(ji,jjm1)*zsmask0(ji,jjm1))/summsk
-                  zsmask1(ji,jj)=1
+                  zsmask1(ji,jj)=1._wp
                END IF
             END DO
          END DO
-         CALL lbc_lnk(sshn,'T',1.)
-         CALL lbc_lnk(zsmask1,'T',1.)
+         CALL lbc_lnk(sshn,'T',1._wp)
+         CALL lbc_lnk(zsmask1,'T',1._wp)
          zssh0   = sshn
          zsmask0 = zsmask1
@@ -257,7 +257,7 @@
          END DO
 
-         hu(:,:) = 0._wp                          ! Ocean depth at U-points
-         hv(:,:) = 0._wp                          ! Ocean depth at V-points
-         ht(:,:) = 0._wp                          ! Ocean depth at T-points
+      ! t-, u- and v- water column thickness
+      ! ------------------------------------
+         ht(:,:) = 0._wp ; hu(:,:) = 0._wp ; hv(:,:) = 0._wp
          DO jk = 1, jpkm1
             hu(:,:) = hu(:,:) + fse3u_n(:,:,jk) * umask(:,:,jk)
@@ -272,5 +272,4 @@
 
 !=============================================================================
-
 ! compute velocity
 ! compute velocity in order to conserve barotropic velocity (modification by poderation of the scale factor).
@@ -332,5 +331,5 @@
       ztmask1(:,:,:) = ptmask_b(:,:,:)
       ztmask0(:,:,:) = ptmask_b(:,:,:)
-      DO iz = 1,10
+      DO iz = 1,10 ! resolution dependent (OK for ISOMIP+ case)
           DO jk = 1,jpk-1
               zdmask=tmask(:,:,jk)-ztmask0(:,:,jk);
@@ -341,5 +340,5 @@
                       summsk= (ztmask0(jip1,jj  ,jk)+ztmask0(jim1,jj  ,jk)+ztmask0(ji  ,jjp1,jk)+ztmask0(ji  ,jjm1,jk))
                       !! horizontal basic extrapolation
-                      IF (zdmask(ji,jj)==1 .AND. summsk .NE. 0) THEN
+                      IF (zdmask(ji,jj)==1._wp .AND. summsk .NE. 0._wp) THEN
                          tsn(ji,jj,jk,1)=( zts0(jip1,jj  ,jk,1)*ztmask0(jip1,jj  ,jk) &
                          &                +zts0(jim1,jj  ,jk,1)*ztmask0(jim1,jj  ,jk) &
@@ -353,13 +352,13 @@
                       END IF
                       !! vertical extrapolation if horizontal extra failed
-                      IF (zdmask(ji,jj)==1 .AND. summsk==0) THEN
+                      IF (zdmask(ji,jj)==1._wp .AND. summsk==0._wp) THEN
                          jkm1=max(1,jk-1) ; jkp1=min(jpk,jk+1)
                          summsk=(ztmask0(ji,jj,jkm1)+ztmask0(ji,jj,jkp1))
-                         IF (zdmask(ji,jj)==1 .AND. summsk .NE. 0 ) THEN
+                         IF (zdmask(ji,jj)==1._wp .AND. summsk .NE. 0._wp ) THEN
                             tsn(ji,jj,jk,1)=( zts0(ji,jj,jkp1,1)*ztmask0(ji,jj,jkp1)     &
                             &                +zts0(ji,jj,jkm1,1)*ztmask0(ji,jj,jkm1))/summsk
                             tsn(ji,jj,jk,2)=( zts0(ji,jj,jkp1,2)*ztmask0(ji,jj,jkp1)     &
                             &                +zts0(ji,jj,jkm1,2)*ztmask0(ji,jj,jkm1))/summsk
-                            ztmask1(ji,jj,jk)=1
+                            ztmask1(ji,jj,jk)=1._wp
                          END IF
                       END IF
@@ -367,10 +366,16 @@
               END DO
           END DO
-          CALL lbc_lnk(tsn(:,:,:,1),'T',1.)
-          CALL lbc_lnk(tsn(:,:,:,2),'T',1.)
-          CALL lbc_lnk(ztmask1,'T',1.)
+          
+          CALL lbc_lnk(tsn(:,:,:,1),'T',1._wp)
+          CALL lbc_lnk(tsn(:,:,:,2),'T',1._wp)
+          CALL lbc_lnk(ztmask1,'T',1._wp)
+
+          ! update
           zts0(:,:,:,:) = tsn(:,:,:,:)
           ztmask0 = ztmask1
-      END DO
+
+      END DO
+
+      ! mask new tsn field
       tsn(:,:,:,jp_tem) = tsn(:,:,:,jp_tem) * tmask(:,:,:)
       tsn(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) * tmask(:,:,:)
@@ -381,10 +386,10 @@
             DO jj = 1,jpj
                DO ji = 1,jpi
-                  IF (zwmaskn(ji,jj,jk) * zwmaskb(ji,jj,jk) == 1.0_wp .AND. (tmask(ji,jj,1)==0 .OR. ptmask_b(ji,jj,1)==0) ) THEN
+                  IF (zwmaskn(ji,jj,jk) * zwmaskb(ji,jj,jk) == 1._wp .AND. (tmask(ji,jj,1)==0._wp .OR. ptmask_b(ji,jj,1)==0._wp) ) THEN
                      !compute weight
                      zdzp1 = MAX(0._wp,fsdepw_n(ji,jj,jk+1) - pdepw_b(ji,jj,jk+1))
                      zdz   =           fsdepw_n(ji,jj,jk+1) - pdepw_b(ji,jj,jk  ) 
                      zdzm1 = MAX(0._wp,pdepw_b(ji,jj,jk  ) - fsdepw_n(ji,jj,jk  ))
-                     IF (zdz .LT. 0.0_wp) THEN 
+                     IF (zdz .LT. 0._wp) THEN 
                         WRITE(numout,*) 'ERROR dz n ', ji,jj,jk,zdz,fsdepw_n(ji,jj,jk+1),fsdepw_n(ji,jj,jk),fsdepw_n(ji,jj,jk-1)
                         WRITE(numout,*) 'ERROR dz n             = ',fse3t_n (ji,jj,jk+1),fse3t_n (ji,jj,jk),fse3t_n (ji,jj,jk-1), sshn(ji,jj)
@@ -411,17 +416,20 @@
       END IF
 
-      ! Special value for closed pool and set the mask to 0 to run
-      WHERE (tmask(:,:,:) == 1.0 .AND. tsn(:,:,:,2) == 0._wp) 
-         tsn(:,:,:,2)=  -99._wp
-         tmask(:,:,:) = 0.0
-         umask(:,:,:) = 0.0
-         vmask(:,:,:) = 0.0
+      ! closed pool
+      ! -----------------------------------------------------------------------------------------
+      ! case we open a cell but no neigbour cells available to get an estimate of T and S
+      WHERE (tmask(:,:,:) == 1._wp .AND. tsn(:,:,:,2) == 0._wp) 
+         tsn(:,:,:,2)=  -99._wp  ! Special value for closed pool (checking purpose in output.init)
+         tmask(:,:,:) = 0._wp    ! set mask to 0 to run
+         umask(:,:,:) = 0._wp
+         vmask(:,:,:) = 0._wp
       END WHERE
-
+      
+      ! set mbkt and mikt to 1 in thiese location
       WHERE (SUM(tmask,dim=3) == 0)
          mbkt(:,:)=1 ; mbku(:,:)=1 ; mbkv(:,:)=1
          mikt(:,:)=1 ; miku(:,:)=1 ; mikv(:,:)=1
       END WHERE
-
+      ! -------------------------------------------------------------------------------------------
       ! compute new tn and sn if we close cell 
       ! nothing to do
