Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/CONFIG/SHARED/namelist_ref
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/CONFIG/SHARED/namelist_ref	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/CONFIG/SHARED/namelist_ref	(revision 5624)
@@ -255,9 +255,5 @@
    ln_dm2dc    = .false.   !  daily mean to diurnal cycle on short wave
    ln_rnf      = .true.    !  runoffs                                   (T   => fill namsbc_rnf)
-   nn_isf      = 0         !  ice shelf melting/freezing                (/=0 => fill namsbc_isf)
-                           !  0 =no isf                  1 = presence of ISF
-                           !  2 = bg03 parametrisation   3 = rnf file for isf
-                           !  4 = ISF fwf specified
-                           !  option 1 and 4 need ln_isfcav = .true. (domzgr)
+   ln_isf      = .false.   !  ice shelf                                 (T   => fill namsbc_isf)
    ln_ssr      = .true.    !  Sea Surface Restoring on T and/or S       (T => fill namsbc_ssr)
    nn_fwb      = 2         !  FreshWater Budget: =0 unchecked
@@ -454,5 +450,8 @@
    sn_Leff_isf = 'rnfisf' ,       0          ,'Leff'         ,   .false.  , .true.  , 'yearly'  ,  ''      ,   ''
 ! for all case
-   ln_divisf   = .true.  ! apply isf melting as a mass flux or in the salinity trend. (maybe I should remove this option as for runoff?)
+   nn_isf      = 1         !  ice shelf melting/freezing
+                           !  1 = presence of ISF    2 = bg03 parametrisation 
+                           !  3 = rnf file for isf   4 = ISF fwf specified
+                           !  option 1 and 4 need ln_isfcav = .true. (domzgr)
 ! only for nn_isf = 1 or 2
    rn_gammat0  = 1.0e-4   ! gammat coefficient used in blk formula
@@ -462,5 +461,4 @@
    rn_hisf_tbl =  30.      ! thickness of the top boundary layer           (Losh et al. 2008)
                           ! 0 => thickness of the tbl = thickness of the first wet cell
-   ln_conserve = .true.   ! conservative case (take into account meltwater advection)
    nn_gammablk = 1        ! 0 = cst Gammat (= gammat/s)
                           ! 1 = velocity dependend Gamma (u* * gammat/s)  (Jenkins et al. 2010)
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/BDY/bdyvol.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/BDY/bdyvol.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/BDY/bdyvol.F90	(revision 5624)
@@ -91,5 +91,5 @@
       ! Calculate the cumulate surface Flux z_cflxemp (m3/s) over all the domain
       ! -----------------------------------------------------------------------
-      z_cflxemp = SUM ( ( emp(:,:)-rnf(:,:)+rdivisf*fwfisf(:,:) ) * bdytmask(:,:) * e1t(:,:) * e2t(:,:) ) / rau0
+      z_cflxemp = SUM ( ( emp(:,:) - rnf(:,:) + fwfisf(:,:) ) * bdytmask(:,:) * e1t(:,:) * e2t(:,:) ) / rau0
       IF( lk_mpp )   CALL mpp_sum( z_cflxemp )     ! sum over the global domain
 
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90	(revision 5624)
@@ -93,5 +93,5 @@
       ! 1 - Trends due to forcing !
       ! ------------------------- !
-      z_frc_trd_v = r1_rau0 * glob_sum( - ( emp(:,:) - rnf(:,:) + rdivisf * fwfisf(:,:) ) * surf(:,:) ) ! volume fluxes
+      z_frc_trd_v = r1_rau0 * glob_sum( - ( emp(:,:) - rnf(:,:) + fwfisf(:,:) ) * surf(:,:) ) ! volume fluxes
       z_frc_trd_t =           glob_sum( sbc_tsc(:,:,jp_tem) * surf(:,:) )                               ! heat fluxes
       z_frc_trd_s =           glob_sum( sbc_tsc(:,:,jp_sal) * surf(:,:) )                               ! salt fluxes
@@ -100,9 +100,5 @@
       IF( ln_rnf_sal)   z_frc_trd_s = z_frc_trd_s + glob_sum( rnf_tsc(:,:,jp_sal) * surf(:,:) )
       ! Add ice shelf heat & salt input
-      IF( nn_isf .GE. 1 )  THEN
-          z_frc_trd_t = z_frc_trd_t &
-              &   + glob_sum( ( risf_tsc(:,:,jp_tem) - rdivisf * fwfisf(:,:) * (-1.9) * r1_rau0 ) * surf(:,:) )
-          z_frc_trd_s = z_frc_trd_s + (1.0_wp - rdivisf) * glob_sum( risf_tsc(:,:,jp_sal) * surf(:,:) )
-      ENDIF
+      IF( ln_isf    ) z_frc_trd_t = z_frc_trd_t + glob_sum( risf_tsc(:,:,jp_tem) * surf(:,:) )
 
       ! Add penetrative solar radiation
@@ -200,5 +196,4 @@
 !      ENDIF
 !!gm end
-
 
       IF( lk_vvl ) THEN
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90	(revision 5624)
@@ -1002,5 +1002,5 @@
 
       IF( ioptio /= 1 )   CALL ctl_stop( 'Choose ONE vertical coordinate in namelist nam_vvl' )
-      IF( .NOT. ln_vvl_zstar .AND. nn_isf .NE. 0) CALL ctl_stop( 'Only vvl_zstar has been tested with ice shelf cavity' )
+      IF( .NOT. ln_vvl_zstar .AND. ln_isf ) CALL ctl_stop( 'Only vvl_zstar has been tested with ice shelf cavity' )
 
       IF(lwp) THEN                   ! Print the choice
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DYN/divcur.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DYN/divcur.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DYN/divcur.F90	(revision 5624)
@@ -26,5 +26,5 @@
    USE oce             ! ocean dynamics and tracers
    USE dom_oce         ! ocean space and time domain
-   USE sbc_oce, ONLY : ln_rnf, nn_isf ! surface boundary condition: ocean
+   USE sbc_oce, ONLY : ln_rnf, ln_isf ! surface boundary condition: ocean
    USE sbcrnf          ! river runoff 
    USE sbcisf          ! ice shelf 
@@ -229,5 +229,5 @@
 
       IF( ln_rnf      )   CALL sbc_rnf_div( hdivn )          ! runoffs   (update hdivn field)
-      IF( ln_divisf .AND. (nn_isf /= 0) )   CALL sbc_isf_div( hdivn )          ! ice shelf (update hdivn field)
+      IF( ln_isf      )   CALL sbc_isf_div( hdivn )          ! ice shelf (update hdivn field)
       IF( nn_cla == 1 )   CALL cla_div    ( kt )             ! Cross Land Advection (Update Hor. divergence)
       
@@ -328,5 +328,5 @@
 
       IF( ln_rnf      )   CALL sbc_rnf_div( hdivn )                            ! runoffs (update hdivn field)
-      IF( ln_divisf .AND. (nn_isf .GT. 0) )   CALL sbc_isf_div( hdivn )          ! ice shelf (update hdivn field)
+      IF( ln_isf      )   CALL sbc_isf_div( hdivn )          ! ice shelf (update hdivn field)
       IF( nn_cla == 1 )   CALL cla_div    ( kt )             ! Cross Land Advection (update hdivn field)
       !
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90	(revision 5624)
@@ -266,15 +266,13 @@
                ! Add volume filter correction: compatibility with tracer advection scheme
                ! => time filter + conservation correction (only at the first level)
-               fse3t_b(:,:,1) = fse3t_b(:,:,1) - atfp * rdt * r1_rau0 * ( emp_b(:,:) - emp(:,:) &
-                              &                                          -rnf_b(:,:) + rnf(:,:) ) * tmask(:,:,1)
-               !DO jj = 1,jpj
-               !   DO ji = 1,jpi
-               !      jk = mikt(ji,jj)
-               !      fse3t_b(ji,jj,jk) = fse3t_b(ji,jj,jk) - atfp * rdt * r1_rau0 
-               !                        &                          * ( (emp_b(ji,jj) - emp(ji,jj) ) &
-               !                        &                            - (rnf_b(ji,jj) - rnf(ji,jj) ) &    
-               !                        &                            + (fwfisf_b(ji,jj) - fwfisf(ji,jj) ) ) * tmask(ji,jj,jk) ! ?????
-               !   END DO
-               !END DO
+               DO jj = 1,jpj
+                  DO ji = 1,jpi
+                     jk = mikt(ji,jj)
+                     fse3t_b(ji,jj,jk) = fse3t_b(ji,jj,jk) - atfp * rdt * r1_rau0                  & 
+                                       &                          * ( (emp_b(ji,jj) - emp(ji,jj) ) &
+                                       &                            - (rnf_b(ji,jj) - rnf(ji,jj) ) &    
+                                       &                            + (fwfisf_b(ji,jj) - fwfisf(ji,jj) ) ) * tmask(ji,jj,jk)
+                  END DO
+               END DO
             ENDIF
             !
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90	(revision 5624)
@@ -454,8 +454,8 @@
       !                                         ! Surface net water flux and rivers
       IF (ln_bt_fw) THEN
-         zssh_frc(:,:) = zraur * ( emp(:,:) - rnf(:,:) + rdivisf * fwfisf(:,:) )
+         zssh_frc(:,:) = zraur * ( emp(:,:) - rnf(:,:) + fwfisf(:,:) )
       ELSE
          zssh_frc(:,:) = zraur * z1_2 * (  emp(:,:) + emp_b(:,:) - rnf(:,:) - rnf_b(:,:)   &
-                &                        + rdivisf * ( fwfisf(:,:) + fwfisf_b(:,:) )       )
+                &                        + fwfisf(:,:) + fwfisf_b(:,:)                     )
       ENDIF
 #if defined key_asminc
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90	(revision 5624)
@@ -268,7 +268,7 @@
       ELSE                                         !** Leap-Frog time-stepping: Asselin filter + swap
          sshb(:,:) = sshn(:,:) + atfp * ( sshb(:,:) - 2 * sshn(:,:) + ssha(:,:) )     ! before <-- now filtered
-         IF( lk_vvl ) sshb(:,:) = sshb(:,:) - atfp * rdt / rau0 * ( emp_b(:,:) - emp(:,:) - rnf_b(:,:) + rnf(:,:) ) * ssmask(:,:)
-         !IF( lk_vvl ) sshb(:,:) = sshb(:,:) - atfp * rdt / rau0 * ( (emp_b(:,:) - emp(:,:)) - (rnf_b(:,:) - rnf(:,:))  &
-         !                 &                                       + (fwfisf_b(:,:) - fwfisf(:,:))) * ssmask(:,:)
+         IF( lk_vvl ) sshb(:,:) = sshb(:,:) - atfp * rdt / rau0 * ( (emp_b(:,:) - emp(:,:))      &
+                          &                                       - (rnf_b(:,:) - rnf(:,:))      &
+                          &                                       + (fwfisf_b(:,:) - fwfisf(:,:))) * ssmask(:,:)
          sshn(:,:) = ssha(:,:)                           ! now <-- after
       ENDIF
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90	(revision 5624)
@@ -44,8 +44,8 @@
    LOGICAL , PUBLIC ::   ln_dm2dc       !: Daily mean to Diurnal Cycle short wave (qsr)
    LOGICAL , PUBLIC ::   ln_rnf         !: runoffs / runoff mouths
+   LOGICAL , PUBLIC ::   ln_isf         !: ice shelf melting
    LOGICAL , PUBLIC ::   ln_ssr         !: Sea Surface restoring on SST and/or SSS      
    LOGICAL , PUBLIC ::   ln_apr_dyn     !: Atmospheric pressure forcing used on dynamics (ocean & ice)
    INTEGER , PUBLIC ::   nn_ice         !: flag for ice in the surface boundary condition (=0/1/2/3)
-   INTEGER , PUBLIC ::   nn_isf         !: flag for isf in the surface boundary condition (=0/1/2/3/4) 
    INTEGER , PUBLIC ::   nn_ice_embd    !: flag for levitating/embedding sea-ice in the ocean
    !                                             !: =0 levitating ice (no mass exchange, concentration/dilution effect)
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcfwb.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcfwb.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcfwb.F90	(revision 5624)
@@ -108,5 +108,5 @@
          !
          IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN
-            z_fwf = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) + rdivisf * fwfisf(:,:) -  snwice_fmass(:,:) ) ) / area   ! sum over the global domain
+            z_fwf = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) + fwfisf(:,:) -  snwice_fmass(:,:) ) ) / area   ! sum over the global domain
             zcoef = z_fwf * rcp
             emp(:,:) = emp(:,:) - z_fwf              * tmask(:,:,1)
@@ -162,5 +162,5 @@
             zsurf_pos = glob_sum( e1e2t(:,:)*ztmsk_pos(:,:) )
             !                                                  ! fwf global mean (excluding ocean to ice/snow exchanges) 
-            z_fwf     = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) + rdivisf * fwfisf(:,:) - snwice_fmass(:,:) ) ) / area
+            z_fwf     = glob_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) + fwfisf(:,:) - snwice_fmass(:,:) ) ) / area
             !            
             IF( z_fwf < 0._wp ) THEN         ! spread out over >0 erp area to increase evaporation
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcisf.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcisf.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcisf.F90	(revision 5624)
@@ -39,11 +39,9 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   qisf              !: net heat flux from ice shelf
    REAL(wp), PUBLIC ::   rn_hisf_tbl                 !: thickness of top boundary layer [m]
-   LOGICAL , PUBLIC ::   ln_divisf                   !: flag to correct divergence 
+   INTEGER , PUBLIC ::   nn_isf                      !: flag to choose between explicit/param/specified  
    INTEGER , PUBLIC ::   nn_isfblk                   !: 
    INTEGER , PUBLIC ::   nn_gammablk                 !:
-   LOGICAL , PUBLIC ::   ln_conserve                 !:
    REAL(wp), PUBLIC ::   rn_gammat0                  !: temperature exchange coeficient
    REAL(wp), PUBLIC ::   rn_gammas0                  !: salinity    exchange coeficient 
-   REAL(wp), PUBLIC ::   rdivisf                     !: flag to test if fwf apply on divergence
 
    REAL(wp)   , PUBLIC, ALLOCATABLE, SAVE, DIMENSION (:,:)     ::  rzisf_tbl              !:depth of calving front (shallowest point) nn_isf ==2/3
@@ -91,4 +89,5 @@
     INTEGER             :: ikt, ikb   ! top and bottom level of the isf boundary layer
     REAL(wp)            :: zhk
+    REAL(wp)            ::   zt_frz, zpress
     CHARACTER(len=256)  :: cvarzisf, cvarhisf   ! name for isf file
     CHARACTER(LEN=32 )  :: cvarLeff                    ! variable name for efficient Length scale
@@ -96,5 +95,5 @@
       !
       !!---------------------------------------------------------------------
-      NAMELIST/namsbc_isf/ nn_isfblk, rn_hisf_tbl, ln_divisf, ln_conserve, rn_gammat0, rn_gammas0, nn_gammablk, &
+      NAMELIST/namsbc_isf/ nn_isfblk, rn_hisf_tbl, rn_gammat0, rn_gammas0, nn_gammablk, nn_isf      , &
                          & sn_fwfisf, sn_qisf, sn_rnfisf, sn_depmax_isf, sn_depmin_isf, sn_Leff_isf
       !
@@ -121,14 +120,8 @@
          IF ( lwp ) WRITE(numout,*) '        nn_isfblk   = ', nn_isfblk
          IF ( lwp ) WRITE(numout,*) '        rn_hisf_tbl = ', rn_hisf_tbl
-         IF ( lwp ) WRITE(numout,*) '        ln_divisf   = ', ln_divisf 
          IF ( lwp ) WRITE(numout,*) '        nn_gammablk = ', nn_gammablk 
          IF ( lwp ) WRITE(numout,*) '        rn_gammat0  = ', rn_gammat0  
          IF ( lwp ) WRITE(numout,*) '        rn_gammas0  = ', rn_gammas0  
          IF ( lwp ) WRITE(numout,*) '        rn_tfri2    = ', rn_tfri2 
-         IF (ln_divisf) THEN       ! keep it in the namelist ??? used true anyway as for runoff ? (PM)
-            rdivisf = 1._wp
-         ELSE
-            rdivisf = 0._wp
-         END IF
          !
          ! Allocate public variable
@@ -270,9 +263,11 @@
          ! compute tsc due to isf
          ! WARNING water add at temp = 0C, correction term is added in trasbc, maybe better here but need a 3D variable).
-         risf_tsc(:,:,jp_tem) = qisf(:,:) * r1_rau0_rcp !
+!         zpress = grav*rau0*fsdept(ji,jj,jk)*1.e-04
+         zt_frz = -1.9_wp !CALL eos_fzp( tsn(ji,jj,jk,jp_sal), zt_frz, zpress )
+         risf_tsc(:,:,jp_tem) = qisf(:,:) * r1_rau0_rcp - fwfisf(:,:) * zt_frz * r1_rau0 !
          
          ! salt effect already take into account in vertical advection
-         risf_tsc(:,:,jp_sal) = (1.0_wp-rdivisf) * fwfisf(:,:) * stbl(:,:) * r1_rau0
-          
+         risf_tsc(:,:,jp_sal) = 0.0_wp
+
          ! lbclnk
          CALL lbc_lnk(risf_tsc(:,:,jp_tem),'T',1.)
@@ -296,5 +291,5 @@
          ! output
          IF( iom_use('qisf'  ) )   CALL iom_put('qisf'  , qisf)
-         IF( iom_use('fwfisf') )   CALL iom_put('fwfisf', fwfisf * stbl / soce )
+         IF( iom_use('fwfisf') )   CALL iom_put('fwfisf', fwfisf)
       END IF
   
@@ -426,7 +421,4 @@
 
       ! initialisation
-      IF (ln_conserve) THEN ;  zcfac=1.0_wp
-      ELSE                  ;  zcfac=0.0_wp
-      ENDIF
       zpress (:,:)=0.0_wp
       zgammat(:,:)=rn_gammat0 ; zgammas (:,:)=rn_gammas0;
@@ -500,10 +492,7 @@
                   ! zfwflx is upward water flux
                   ! zhtflx is upward heat flux (out of ocean)
-! If non conservative we have zcfac=0.0 so zhtflx is as ISOMIP but with different zfrz value
-! If non conservative we have zcfac=0.0 so what follows is then zfwflx*sss_m/zsfrz
                   ! compute the upward water and heat flux
                   zfwflx(ji,jj) = rau0 * zgammas(ji,jj) * (zsfrz-stbl(ji,jj)) / MAX(zsfrz,zeps)
-                  zhtflx(ji,jj) = ( zgammat(ji,jj) * rau0 - zcfac * zfwflx(ji,jj) ) * rcp * (zti(ji,jj) - zfrz(ji,jj) ) 
-                  zfwflx(ji,jj) = ( zgammas(ji,jj) * rau0 - zcfac * zfwflx(ji,jj) ) * (zsfrz - stbl(ji,jj)) / MAX(stbl(ji,jj),zeps)
+                  zhtflx(ji,jj) = zgammat(ji,jj) * rau0 * rcp * (zti(ji,jj) - zfrz(ji,jj) ) 
                END DO
             END DO
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90	(revision 5624)
@@ -88,5 +88,5 @@
       NAMELIST/namsbc/ nn_fsbc   , ln_ana    , ln_flx, ln_blk_clio, ln_blk_core, ln_mixcpl,   &
          &             ln_blk_mfs, ln_apr_dyn, nn_ice, nn_ice_embd, ln_dm2dc   , ln_rnf   ,   &
-         &             ln_ssr    , nn_isf    , nn_fwb, ln_cdgw    , ln_wave    , ln_sdw   ,   &
+         &             ln_ssr    , ln_isf    , nn_fwb, ln_cdgw    , ln_wave    , ln_sdw   ,   &
          &             nn_lsm    , nn_limflx , nn_components, ln_cpl
       INTEGER  ::   ios
@@ -141,5 +141,5 @@
          WRITE(numout,*) '              daily mean to diurnal cycle qsr            ln_dm2dc    = ', ln_dm2dc 
          WRITE(numout,*) '              runoff / runoff mouths                     ln_rnf      = ', ln_rnf
-         WRITE(numout,*) '              iceshelf formulation                       nn_isf      = ', nn_isf
+         WRITE(numout,*) '              iceshelf melting                           ln_isf      = ', ln_isf
          WRITE(numout,*) '              Sea Surface Restoring on SST and/or SSS    ln_ssr      = ', ln_ssr
          WRITE(numout,*) '              FreshWater Budget control  (=0/1/2)        nn_fwb      = ', nn_fwb
@@ -179,9 +179,8 @@
 
       !                          ! Checks:
-      IF( nn_isf .EQ. 0 ) THEN                      ! variable initialisation if no ice shelf 
+      IF( ln_isf ) THEN                      ! variable initialisation if no ice shelf 
          IF( sbc_isf_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_init : unable to allocate sbc_isf arrays' )
          fwfisf  (:,:) = 0.0_wp
          fwfisf_b(:,:) = 0.0_wp
-         rdivisf       = 0.0_wp
       END IF
       IF( nn_ice == 0 .AND. nn_components /= jp_iam_opa )   fr_i(:,:) = 0.e0 ! no ice in the domain, ice fraction is always zero
@@ -396,5 +395,5 @@
       IF( ln_icebergs    )   CALL icb_stp( kt )                   ! compute icebergs
 
-      IF( nn_isf   /= 0  )   CALL sbc_isf( kt )                    ! compute iceshelves
+      IF( ln_isf         )   CALL sbc_isf( kt )                   ! compute iceshelves
 
       IF( ln_rnf         )   CALL sbc_rnf( kt )                   ! add runoffs to fresh water fluxes
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/TRA/tranxt.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/TRA/tranxt.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/TRA/tranxt.F90	(revision 5624)
@@ -28,4 +28,5 @@
    USE sbc_oce         ! surface boundary condition: ocean
    USE sbcrnf          ! river runoffs
+   USE sbcisf          ! ice shelf melting
    USE zdf_oce         ! ocean vertical mixing
    USE domvvl          ! variable volume
@@ -279,5 +280,5 @@
 
       !!     
-      LOGICAL  ::   ll_tra_hpg, ll_traqsr, ll_rnf   ! local logical
+      LOGICAL  ::   ll_tra_hpg, ll_traqsr, ll_rnf, ll_isf   ! local logical
       INTEGER  ::   ji, jj, jk, jn              ! dummy loop indices
       REAL(wp) ::   zfact1, ztc_a , ztc_n , ztc_b , ztc_f , ztc_d    ! local scalar
@@ -295,8 +296,10 @@
          ll_traqsr  = ln_traqsr        ! active  tracers case  and  solar penetration
          ll_rnf     = ln_rnf           ! active  tracers case  and  river runoffs
+         ll_isf     = ln_isf           ! active  tracers case  and  ice shelf melting
       ELSE                          
          ll_tra_hpg = .FALSE.          ! passive tracers case or NO semi-implicit hpg
          ll_traqsr  = .FALSE.          ! active  tracers case and NO solar penetration
          ll_rnf     = .FALSE.          ! passive tracers or NO river runoffs
+         ll_isf     = .FALSE.          ! passive tracers or NO ice shelf melting
       ENDIF
       !
@@ -321,10 +324,8 @@
                   ztc_f  = ztc_n  + atfp * ztc_d
                   !
-                  IF( jk == 1 ) THEN           ! first level 
-                !  IF( jk == mikt(ji,jj) ) THEN           ! first level 
-                     ze3t_f = ze3t_f - zfact2 * ( emp_b(ji,jj) - emp(ji,jj) + rnf(ji,jj) - rnf_b(ji,jj))
-                !     ze3t_f = ze3t_f - zfact2 * ( (emp_b(ji,jj) - emp(ji,jj))          &
-                !            &                   - (rnf_b(ji,jj) - rnf(ji,jj))          &
-                !            &                   + (fwfisf_b(ji,jj) - fwfisf(ji,jj)) )
+                  IF( jk == mikt(ji,jj) ) THEN           ! first level 
+                     ze3t_f = ze3t_f - zfact2 * ( (emp_b(ji,jj) - emp(ji,jj))          &
+                            &                   - (rnf_b(ji,jj) - rnf(ji,jj))          &
+                            &                   + (fwfisf_b(ji,jj) - fwfisf(ji,jj)) )
                             
                      ztc_f  = ztc_f  - zfact1 * ( psbc_tc(ji,jj,jn) - psbc_tc_b(ji,jj,jn) )
@@ -334,11 +335,20 @@
                      &     ztc_f  = ztc_f  - zfact1 * ( qsr_hc(ji,jj,jk) - qsr_hc_b(ji,jj,jk) ) 
 
-                  IF( ll_rnf .AND. jk <= nk_rnf(ji,jj) )   & 		      ! river runoffs
+                  ! river runoffs
+                  IF( ll_rnf .AND. jk <= nk_rnf(ji,jj) )                                          & 
                      &     ztc_f  = ztc_f  - zfact1 * ( rnf_tsc(ji,jj,jn) - rnf_tsc_b(ji,jj,jn) ) & 
                      &                              * fse3t_n(ji,jj,jk) / h_rnf(ji,jj)
 
-                  !IF( nn_isf .NE. 0 .AND. jk == mikt(ji,jj) )   & 		      ! ice shelf
-                  !   &     ztc_f  = ztc_f  + zfact1 * ( risf_tsc(ji,jj,jn) - risf_tsc_b(ji,jj,jn) ) & 
-                  !   &                              * fse3t_n(ji,jj,jk) * r1_hisf_tbl (ji,jj)
+                  ! ice shelf
+                  IF( ll_isf ) THEN
+                     ! level fully include in the Losch_2008 ice shelf boundary layer
+                     IF ( jk >= misfkt(ji,jj) .AND. jk < misfkb(ji,jj) )                             &
+                        ztc_f  = ztc_f  - zfact1 * ( risf_tsc(ji,jj,jn) - risf_tsc_b(ji,jj,jn) ) & 
+                               &                 * fse3t_n(ji,jj,jk) * r1_hisf_tbl (ji,jj)
+                     ! level partially include in Losch_2008 ice shelf boundary layer 
+                     IF ( jk == misfkb(ji,jj) )                                                    &
+                        ztc_f  = ztc_f  - zfact1 * ( risf_tsc(ji,jj,jn) - risf_tsc_b(ji,jj,jn) ) & 
+                               &                 * fse3t_n(ji,jj,jk) * r1_hisf_tbl (ji,jj) * ralpha(ji,jj)
+                  END IF
 
                   ze3t_f = 1.e0 / ze3t_f
Index: /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/TRA/trasbc.F90
===================================================================
--- /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/TRA/trasbc.F90	(revision 5623)
+++ /branches/2015/dev_r5151_UKMO_ISF/NEMOGCM/NEMO/OPA_SRC/TRA/trasbc.F90	(revision 5624)
@@ -217,5 +217,5 @@
       !----------------------------------------
       !
-      IF( nn_isf > 0 ) THEN
+      IF( ln_isf ) THEN
          zfact = 0.5_wp
          DO jj = 2, jpj
@@ -226,32 +226,17 @@
    
                ! level fully include in the ice shelf boundary layer
-               ! if isfdiv, we have to remove heat flux due to inflow at 0oC (as in rnf when you add rnf at sst)
                ! sign - because fwf sign of evapo (rnf sign of precip)
                DO jk = ikt, ikb - 1
-               ! compute tfreez for the temperature correction (we add water at freezing temperature)
-!                  zpress = grav*rau0*fsdept(ji,jj,jk)*1.e-04
-                  zt_frz = -1.9_wp !CALL eos_fzp( tsn(ji,jj,jk,jp_sal), zt_frz, zpress )
                ! compute trend
                   tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem)                                                &
-                     &           + zfact * ( risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)               &
-                     &                     - rdivisf * (fwfisf(ji,jj) + fwfisf_b(ji,jj)) * zt_frz * r1_rau0) &
-                     &           * r1_hisf_tbl(ji,jj)
-                  tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal)                                                &
-                     &           + zfact * (risf_tsc_b(ji,jj,jp_sal) + risf_tsc(ji,jj,jp_sal))               &
+                     &           + zfact * ( risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem) )             &
                      &           * r1_hisf_tbl(ji,jj)
                END DO
    
                ! level partially include in ice shelf boundary layer 
-               ! compute tfreez for the temperature correction (we add water at freezing temperature)
-!               zpress = grav*rau0*fsdept(ji,jj,ikb)*1.e-04
-               zt_frz = -1.9_wp !CALL eos_fzp( tsn(ji,jj,ikb,jp_sal), zt_frz, zpress )
                ! compute trend
                tsa(ji,jj,ikb,jp_tem) = tsa(ji,jj,ikb,jp_tem)                                                 &
-                  &              + zfact * ( risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem)               &
-                  &                        - rdivisf * (fwfisf(ji,jj) + fwfisf_b(ji,jj)) * zt_frz * r1_rau0) & 
+                  &              + zfact * ( risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem) )             &
                   &              * r1_hisf_tbl(ji,jj) * ralpha(ji,jj)
-               tsa(ji,jj,ikb,jp_sal) = tsa(ji,jj,ikb,jp_sal)                                                 &
-                  &              + zfact * (risf_tsc_b(ji,jj,jp_sal) + risf_tsc(ji,jj,jp_sal))               &
-                  &              * r1_hisf_tbl(ji,jj) * ralpha(ji,jj) 
 
             END DO
