Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/cfgs/SHARED/field_def_nemo-oce.xml
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/cfgs/SHARED/field_def_nemo-oce.xml	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/cfgs/SHARED/field_def_nemo-oce.xml	(revision 11521)
@@ -264,4 +264,6 @@
 
           <!-- * variable related to ice shelf forcing * -->
+          <field id="isfdebug3d"                                                                                grid_ref="grid_T_3D" />
+          <field id="isfdebug2d"                                                                                />
           <field id="fwfisf_cav"      long_name="Ice shelf melting"                             unit="kg/m2/s"  />
           <field id="fwfisf_par"      long_name="Ice shelf melting"                             unit="kg/m2/s"  />
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/cfgs/SHARED/namelist_ref
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/cfgs/SHARED/namelist_ref	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/cfgs/SHARED/namelist_ref	(revision 11521)
@@ -435,53 +435,55 @@
 !-----------------------------------------------------------------------
    !
-   ! ---------------- general -------------------------------
+   ! ---------------- ice shelf load -------------------------------
+   !
+   cn_isfload = 'isomip'      ! scheme to compute ice shelf load (ln_isfcav = .true. in domain_cfg.nc)
+   !
+   ! ---------------- ice shelf melt formulation -------------------------------
    !
    ln_isf = .false.           ! activate ice shelf module
-   cn_isfdir = './'           ! directory for all ice shelf input file
-   !
-   ! ---------------- cavities opened -------------------------------
-   !
-   cn_isfload = 'isomip'
-   !
-   ln_isfcav_mlt = .false.    ! ice shelf melting into the cavity
-      cn_isfcav_mlt = '3eq'   ! ice shelf melting formulation (spe/2eq/3eq/oasis)
-      !                       ! spe = fwfisf is read from a forcing field
-      !                       ! 2eq = ISOMIP  like: 2 equations formulation (Hunter et al., 2006)
-      !                       ! 3eq = ISOMIP+ like: 3 equations formulation (Asay-Davis et al., 2015)
-      !                       ! oasis = fwfisf is given by oasis and pattern by file sn_isfcav_fwf
-      !              !  cn_isfcav_mlt = 2eq or 3eq cases:
-      cn_gammablk = 'ad15'    ! scheme to compute gammat/s (spe,ad15,hj99)
-      !                       ! ad15 = velocity dependend Gamma (u* * gammat/s)  (Jenkins et al. 2010)
-      !                       ! hj99 = velocity and stability dependent Gamma    (Holland et al. 1999)
-      rn_gammat0  = 1.e-4     ! gammat coefficient used in blk formula
-      rn_gammas0  = 1.e-4     ! gammas coefficient used in blk formula
+      cn_isfdir = './'           ! directory for all ice shelf input file
       !
-      rn_htbl     =  30.      ! thickness of the top boundary layer    (Losh et al. 2008)
-      !                       ! 0 => thickness of the tbl = thickness of the first wet cell
+      ! ---------------- cavities opened -------------------------------
       !
-      !* 'spe' and 'oasis' case
-      !___________!_____________!___________________!___________!_____________!_________!___________!__________!__________!_______________!
-      !           !  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_isfcav_fwf = 'isfmlt_cav',      -12       , 'fwflisf' ,  .false.    , .true.  , 'yearly'  ,    ''    ,   ''     ,    ''
-   !
-   ! ---------------- cavities parametrised -------------------------------
-   !
-   ln_isfpar_mlt = .false.   ! ice shelf melting parametrised
-      cn_isfpar_mlt = 'spe'  ! ice shelf melting parametrisation (spe/bg03/oasis)
-      !                      ! spe   = fwfisf is read from a forcing field
-      !                      ! bg03  = melt computed using Beckmann and Goosse parametrisation
-      !                      ! oasis = fwfisf is given by oasis and pattern by file sn_isfpar_fwf
+      ln_isfcav_mlt = .false.    ! ice shelf melting into the cavity (need ln_isfcav = .true. in domain_cfg.nc)
+         cn_isfcav_mlt = '3eq'   ! ice shelf melting formulation (spe/2eq/3eq/oasis)
+         !                       ! spe = fwfisf is read from a forcing field
+         !                       ! 2eq = ISOMIP  like: 2 equations formulation (Hunter et al., 2006)
+         !                       ! 3eq = ISOMIP+ like: 3 equations formulation (Asay-Davis et al., 2015)
+         !                       ! oasis = fwfisf is given by oasis and pattern by file sn_isfcav_fwf
+         !              !  cn_isfcav_mlt = 2eq or 3eq cases:
+         cn_gammablk = 'ad15'    ! scheme to compute gammat/s (spe,ad15,hj99)
+         !                       ! ad15 = velocity dependend Gamma (u* * gammat/s)  (Jenkins et al. 2010)
+         !                       ! hj99 = velocity and stability dependent Gamma    (Holland et al. 1999)
+         rn_gammat0  = 1.e-4     ! gammat coefficient used in blk formula
+         rn_gammas0  = 1.e-4     ! gammas coefficient used in blk formula
+         !
+         rn_htbl     =  30.      ! thickness of the top boundary layer    (Losh et al. 2008)
+         !                       ! 0 => thickness of the tbl = thickness of the first wet cell
+         !
+         !* 'spe' and 'oasis' case
+         !___________!_____________!___________________!___________!_____________!_________!___________!__________!__________!_______________!
+         !           !  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_isfcav_fwf = 'isfmlt_cav',      -12       , 'fwflisf' ,  .false.    , .true.  , 'yearly'  ,    ''    ,   ''     ,    ''
       !
-      !* all cases
-      !___________!_____________!___________________!___________!_____________!_________!___________!__________!__________!_______________!
-      !           !  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_isfpar_zmax = 'isfmlt_par',       0        ,'sozisfmax',  .false.    , .true.  , 'yearly'  ,    ''    ,   ''     ,    ''
-      sn_isfpar_zmin = 'isfmlt_par',       0        ,'sozisfmin',  .false.    , .true.  , 'yearly'  ,    ''    ,   ''     ,    ''
-      !* 'spe' and 'oasis' case
-      sn_isfpar_fwf = 'isfmlt_par' ,      -12       ,'fwfisf'   ,  .false.    , .true.  , 'yearly'  ,    ''    ,   ''     ,    ''
-      !* 'bg03' case
-      sn_isfpar_Leff = 'isfmlt_par',       0        ,'Leff'     ,  .false.    , .true.  , 'yearly'  ,    ''    ,   ''     ,    ''
+      ! ---------------- cavities parametrised -------------------------------
+      !
+      ln_isfpar_mlt = .false.   ! ice shelf melting parametrised
+         cn_isfpar_mlt = 'spe'  ! ice shelf melting parametrisation (spe/bg03/oasis)
+         !                      ! spe   = fwfisf is read from a forcing field
+         !                      ! bg03  = melt computed using Beckmann and Goosse parametrisation
+         !                      ! oasis = fwfisf is given by oasis and pattern by file sn_isfpar_fwf
+         !
+         !* all cases
+         !___________!_____________!___________________!___________!_____________!_________!___________!__________!__________!_______________!
+         !           !  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_isfpar_zmax = 'isfmlt_par',       0        ,'sozisfmax',  .false.    , .true.  , 'yearly'  ,    ''    ,   ''     ,    ''
+         sn_isfpar_zmin = 'isfmlt_par',       0        ,'sozisfmin',  .false.    , .true.  , 'yearly'  ,    ''    ,   ''     ,    ''
+         !* 'spe' and 'oasis' case
+         sn_isfpar_fwf = 'isfmlt_par' ,      -12       ,'fwfisf'   ,  .false.    , .true.  , 'yearly'  ,    ''    ,   ''     ,    ''
+         !* 'bg03' case
+         sn_isfpar_Leff = 'isfmlt_par',       0        ,'Leff'     ,  .false.    , .true.  , 'yearly'  ,    ''    ,   ''     ,    ''
    !
    ! ---------------- ice sheet coupling -------------------------------
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DIA/diahsb.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DIA/diahsb.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DIA/diahsb.F90	(revision 11521)
@@ -48,4 +48,5 @@
    REAL(wp), DIMENSION(:,:)  , ALLOCATABLE ::   ssh_hc_loc_ini, ssh_sc_loc_ini   !
    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   hc_loc_ini, sc_loc_ini, e3t_ini  !
+   REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   tmask_ini
 
    !! * Substitutions
@@ -156,13 +157,13 @@
       !
       DO jk = 1, jpkm1           ! volume variation (calculated with scale factors)
-         zwrk(:,:,jk) = ( surf(:,:)*e3t_n(:,:,jk) - surf_ini(:,:)*e3t_ini(:,:,jk) ) * tmask(:,:,jk)
+         zwrk(:,:,jk) = surf(:,:)*e3t_n(:,:,jk)*tmask(:,:,jk) - surf_ini(:,:)*e3t_ini(:,:,jk)*tmask_ini(:,:,jk)
       END DO
-      zdiff_v2 = glob_sum_full( 'diahsb', zwrk(:,:,:) )
+      zdiff_v2 = glob_sum_full( 'diahsb', zwrk(:,:,:) )     ! glob_sum_full needed as tmask and tmask_ini could be different
       DO jk = 1, jpkm1           ! heat content variation
-         zwrk(:,:,jk) = ( surf(:,:)*e3t_n(:,:,jk)*tsn(:,:,jk,jp_tem) - surf_ini(:,:)*hc_loc_ini(:,:,jk) ) * tmask(:,:,jk)
+         zwrk(:,:,jk) = ( surf(:,:)*e3t_n(:,:,jk)*tsn(:,:,jk,jp_tem) - surf_ini(:,:)*hc_loc_ini(:,:,jk) )
       END DO
       zdiff_hc = glob_sum_full( 'diahsb', zwrk(:,:,:) )
       DO jk = 1, jpkm1           ! salt content variation
-         zwrk(:,:,jk) = ( surf(:,:)*e3t_n(:,:,jk)*tsn(:,:,jk,jp_sal) - surf_ini(:,:)*sc_loc_ini(:,:,jk) ) * tmask(:,:,jk)
+         zwrk(:,:,jk) = ( surf(:,:)*e3t_n(:,:,jk)*tsn(:,:,jk,jp_sal) - surf_ini(:,:)*sc_loc_ini(:,:,jk) )
       END DO
       zdiff_sc = glob_sum_full( 'diahsb', zwrk(:,:,:) )
@@ -188,5 +189,5 @@
          zwrk(:,:,jk) = surf(:,:) * e3t_n(:,:,jk) * tmask(:,:,jk)
       END DO
-      zvol_tot = glob_sum_full( 'diahsb', zwrk(:,:,:) )
+      zvol_tot = glob_sum( 'diahsb', zwrk(:,:,:) )
 
 !!gm to be added ?
@@ -271,4 +272,5 @@
             CALL iom_get( numror, jpdom_autoglo, 'ssh_ini'   , ssh_ini   , ldxios = lrxios )
             CALL iom_get( numror, jpdom_autoglo, 'e3t_ini'   , e3t_ini   , ldxios = lrxios )
+            CALL iom_get( numror, jpdom_autoglo, 'tmask_ini' , tmask_ini , ldxios = lrxios )
             CALL iom_get( numror, jpdom_autoglo, 'hc_loc_ini', hc_loc_ini, ldxios = lrxios )
             CALL iom_get( numror, jpdom_autoglo, 'sc_loc_ini', sc_loc_ini, ldxios = lrxios )
@@ -326,4 +328,5 @@
          CALL iom_rstput( kt, nitrst, numrow, 'ssh_ini'   , ssh_ini   , ldxios = lwxios )
          CALL iom_rstput( kt, nitrst, numrow, 'e3t_ini'   , e3t_ini   , ldxios = lwxios )
+         CALL iom_rstput( kt, nitrst, numrow, 'tmask_ini' , tmask_ini , ldxios = lwxios )
          CALL iom_rstput( kt, nitrst, numrow, 'hc_loc_ini', hc_loc_ini, ldxios = lwxios )
          CALL iom_rstput( kt, nitrst, numrow, 'sc_loc_ini', sc_loc_ini, ldxios = lwxios )
@@ -397,5 +400,5 @@
       ! ------------------- !
       ALLOCATE( hc_loc_ini(jpi,jpj,jpk), sc_loc_ini(jpi,jpj,jpk), surf_ini(jpi,jpj), &
-         &      e3t_ini(jpi,jpj,jpk), surf(jpi,jpj),  ssh_ini(jpi,jpj), STAT=ierror  )
+         &      e3t_ini(jpi,jpj,jpk), surf(jpi,jpj),  ssh_ini(jpi,jpj), tmask_ini(jpi,jpj,jpk),STAT=ierror  )
       IF( ierror > 0 ) THEN
          CALL ctl_stop( 'dia_hsb_init: unable to allocate hc_loc_ini' )   ;   RETURN
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DIA/diawri.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DIA/diawri.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DIA/diawri.F90	(revision 11521)
@@ -27,4 +27,5 @@
    USE oce            ! ocean dynamics and tracers 
    USE isf
+   USE isfcpl
    USE dom_oce        ! ocean space and time domain
    USE phycst         ! physical constants
@@ -884,6 +885,7 @@
       !!----------------------------------------------------------------------
       CHARACTER (len=* ), INTENT( in ) ::   cdfile_name      ! name of the file created
+      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zisfdebug
       !!
-      INTEGER :: inum
+      INTEGER :: inum, jk
       !!----------------------------------------------------------------------
       ! 
@@ -923,4 +925,8 @@
          END IF
       END IF
+      DO jk = 1,jpk
+         zisfdebug(:,:,jk) = risfcpl_cons_vol(:,:,jk) * rdt_iscpl / e1e2t(:,:)
+      END DO
+      CALL iom_rstput( 0, 0, inum, 'isfdebug', zisfdebug)
       IF( ALLOCATED(ahtu) ) THEN
          CALL iom_rstput( 0, 0, inum,  'ahtu', ahtu              )    ! aht at u-point
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DOM/dom_oce.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DOM/dom_oce.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DOM/dom_oce.F90	(revision 11521)
@@ -35,5 +35,4 @@
    REAL(wp), PUBLIC ::   rn_atfp        !: asselin time filter parameter
    INTEGER , PUBLIC ::   nn_euler       !: =0 start with forward time step or not (=1)
-   LOGICAL , PUBLIC ::   ln_iscpl       !: coupling with ice sheet
    LOGICAL , PUBLIC ::   ln_crs         !: Apply grid coarsening to dynamical model output or online passive tracers
 
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DOM/domain.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DOM/domain.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DOM/domain.F90	(revision 11521)
@@ -185,9 +185,8 @@
       IF( lk_c1d         )   CALL cor_c1d       ! 1D configuration: Coriolis set at T-point
       !
-      IF( ln_meshmask .AND. .NOT.ln_iscpl )                        CALL dom_wri     ! Create a domain file
-      IF( ln_meshmask .AND.      ln_iscpl .AND. .NOT.ln_rstart )   CALL dom_wri     ! Create a domain file
-      IF(                                       .NOT.ln_rstart )   CALL dom_ctl     ! Domain control
-      !
-      IF( ln_write_cfg )   CALL cfg_write         ! create the configuration file
+      IF( ln_meshmask    )   CALL dom_wri       ! Create a domain file
+      IF( .NOT.ln_rstart )   CALL dom_ctl       ! Domain control
+      !
+      IF( ln_write_cfg   )   CALL cfg_write     ! create the configuration file
       !
       IF(lwp) THEN
@@ -285,5 +284,5 @@
          &             nn_it000, nn_itend , nn_date0    , nn_time0     , nn_leapy  , nn_istate ,     &
          &             nn_stock, nn_write , ln_mskland  , ln_clobber   , nn_chunksz, nn_euler  ,     &
-         &             ln_cfmeta, ln_iscpl, ln_xios_read, nn_wxios
+         &             ln_cfmeta, ln_xios_read, nn_wxios
       NAMELIST/namdom/ ln_linssh, rn_rdt, rn_atfp, ln_crs, ln_meshmask
 #if defined key_netcdf4
@@ -334,5 +333,4 @@
          WRITE(numout,*) '      overwrite an existing file      ln_clobber      = ', ln_clobber
          WRITE(numout,*) '      NetCDF chunksize (bytes)        nn_chunksz      = ', nn_chunksz
-         WRITE(numout,*) '      IS coupling at the restart step ln_iscpl        = ', ln_iscpl
          IF( TRIM(Agrif_CFixed()) == '0' ) THEN
             WRITE(numout,*) '      READ restart for a single file using XIOS ln_xios_read =', ln_xios_read
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DOM/domvvl.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DOM/domvvl.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DOM/domvvl.F90	(revision 11521)
@@ -20,4 +20,5 @@
    USE oce             ! ocean dynamics and tracers
    USE isfcpl
+   USE isf
    USE phycst          ! physical constant
    USE dom_oce         ! ocean space and time domain
@@ -816,5 +817,5 @@
             !
             ! coupling with an ice sheet model (grounding/calving line migration)
-            IF ( ln_iscpl ) THEN
+            IF ( ln_isfcpl ) THEN
                !
                ! extrapolate ssh over new water column if needed
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/divhor.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/divhor.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/divhor.F90	(revision 11521)
@@ -19,6 +19,8 @@
    !!----------------------------------------------------------------------
    USE oce             ! ocean dynamics and tracers
+   USE isf
+   USE isfutils
    USE dom_oce         ! ocean space and time domain
-   USE sbc_oce, ONLY : ln_rnf, ln_isf ! surface boundary condition: ocean
+   USE sbc_oce, ONLY : ln_rnf ! surface boundary condition: ocean
    USE sbcrnf          ! river runoff 
    USE isfhdiv         ! ice shelf
@@ -98,5 +100,7 @@
       ! 
 #endif
-      IF( ln_isf )   CALL isf_hdiv( kt, hdivn )      !==  ice shelf  ==!   (update hdivn field)
+      CALL debug('isfdivdebug 0',hdivn)
+      IF( ln_isf .OR. ll_isfcpl_cons )   CALL isf_hdiv( kt, hdivn )      !==  ice shelf  ==!   (update hdivn field)
+      CALL debug('isfdivdebug 3',hdivn)
       !
       CALL lbc_lnk( 'divhor', hdivn, 'T', 1. )   !   (no sign change)
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/dynnxt.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/dynnxt.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/dynnxt.F90	(revision 11521)
@@ -26,4 +26,5 @@
    !!-------------------------------------------------------------------------
    USE oce            ! ocean dynamics and tracers
+   USE isf
    USE dom_oce        ! ocean space and time domain
    USE sbc_oce        ! Surface boundary condition: ocean fields
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/dynspg_ts.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/dynspg_ts.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/dynspg_ts.F90	(revision 11521)
@@ -34,4 +34,5 @@
    USE zdfdrg          ! vertical physics: top/bottom drag coef.
    USE isf             ! ice shelf variable (fwfisf)
+   USE isfutils
    USE sbcapr          ! surface boundary condition: atmospheric pressure
    USE dynadv    , ONLY: ln_dynadv_vec
@@ -643,4 +644,8 @@
          zssh_frc(:,:) = zssh_frc(:,:) + div_sd(:,:)
       ENDIF
+      !
+      IF( ll_isfcpl_cons ) THEN
+         zssh_frc(:,:) = zssh_frc(:,:) + risfcpl_cons_ssh(:,:)
+      END IF
       !
 #if defined key_asminc
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/sshwzv.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/sshwzv.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/DYN/sshwzv.F90	(revision 11521)
@@ -18,4 +18,5 @@
    USE oce            ! ocean dynamics and tracers variables
    USE isf            ! ice shelf
+   USE isfutils
    USE dom_oce        ! ocean space and time domain variables 
    USE sbc_oce        ! surface boundary condition: ocean
@@ -255,6 +256,6 @@
          IF( .NOT.ln_linssh ) THEN                          ! before <-- with forcing removed
             zcoef = atfp * rdt * r1_rau0
-            sshb(:,:) = sshb(:,:) - zcoef * (     emp_b(:,:) - emp   (:,:)   &
-               &                             -    rnf_b(:,:) + rnf   (:,:)   &
+            sshb(:,:) = sshb(:,:) - zcoef * (  emp_b(:,:)        - emp   (:,:)       &
+               &                             - rnf_b(:,:)        + rnf   (:,:)       &
                &                             + fwfisf_cav_b(:,:) - fwfisf_cav(:,:)   &
                &                             + fwfisf_par_b(:,:) - fwfisf_par(:,:)   ) * ssmask(:,:)
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/IOM/iom.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/IOM/iom.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/IOM/iom.F90	(revision 11521)
@@ -378,15 +378,4 @@
         CALL iom_set_rstw_var_active('sshn')
         CALL iom_set_rstw_var_active('rhop')
-     ! extra variable needed for the ice sheet coupling
-        IF ( ln_iscpl ) THEN
-             CALL iom_set_rstw_var_active('tmask')
-             CALL iom_set_rstw_var_active('umask')
-             CALL iom_set_rstw_var_active('vmask')
-             CALL iom_set_rstw_var_active('smask')
-             CALL iom_set_rstw_var_active('e3t_n')
-             CALL iom_set_rstw_var_active('e3u_n')
-             CALL iom_set_rstw_var_active('e3v_n')
-             CALL iom_set_rstw_var_active('gdepw_n')
-        END IF
       ENDIF
       IF(ln_diurnal) CALL iom_set_rstw_var_active('Dsst')
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isf.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isf.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isf.F90	(revision 11521)
@@ -30,4 +30,5 @@
    ! 0.1 -------- ice shelf cavity parameter --------------
    CHARACTER(LEN=256), PUBLIC :: cn_isfdir
+   LOGICAL           , PUBLIC :: ln_isf
    !
    ! 0.2 -------- ice shelf cavity opened namelist parameter -------------
@@ -100,6 +101,7 @@
    INTEGER , PUBLIC                                        ::   nstp_iscpl   !:
    REAL(wp), PUBLIC                                        ::   rdt_iscpl    !: 
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   risfcpl_vol, risfcpl_cons_vol  !:
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   risfcpl_tsc, risfcpl_cons_tsc  !:
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   risfcpl_cons_ssh, risfcpl_cons_ssh_b               !:
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   risfcpl_vol, risfcpl_cons_vol, risfcpl_cons_vol_b  !:
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   risfcpl_tsc, risfcpl_cons_tsc, risfcpl_cons_tsc_b  !:
    !
    !!----------------------------------------------------------------------
@@ -184,5 +186,7 @@
       !
       ALLOCATE( risfcpl_tsc(jpi,jpj,jpk,jpts)      , risfcpl_vol(jpi,jpj,jpk)      ,             &
-         &      risfcpl_cons_tsc(jpi,jpj,jpk,jpts) , risfcpl_cons_vol(jpi,jpj,jpk) , STAT=ialloc )
+         &      risfcpl_cons_tsc(jpi,jpj,jpk,jpts) , risfcpl_cons_vol(jpi,jpj,jpk) ,             &
+         &      risfcpl_cons_tsc_b(jpi,jpj,jpk,jpts) , risfcpl_cons_vol_b(jpi,jpj,jpk),          &
+                risfcpl_cons_ssh_b(jpi,jpj), risfcpl_cons_ssh(jpi,jpj), STAT=ialloc )
       ierr = ierr + ialloc
       !
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcav.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcav.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcav.F90	(revision 11521)
@@ -78,6 +78,6 @@
       ! 
       ! output T/S/U/V for the top boundary layer
-      CALL iom_put('ttbl_cav',zttbl(:,:))
-      CALL iom_put('stbl'    ,zstbl(:,:))
+      CALL iom_put('ttbl_cav',zttbl(:,:) * mskisf_cav(:,:))
+      CALL iom_put('stbl'    ,zstbl(:,:) * mskisf_cav(:,:))
       !
       ! initialisation
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcavgam.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcavgam.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfcavgam.F90	(revision 11521)
@@ -64,6 +64,6 @@
          !
          ! mask velocity in tbl with ice shelf mask
-         !zutbl(:,:) = zutbl(:,:) * mskisf_cav(:,:)
-         !zvtbl(:,:) = zvtbl(:,:) * mskisf_cav(:,:)
+         zutbl(:,:) = zutbl(:,:) * mskisf_cav(:,:)
+         zvtbl(:,:) = zvtbl(:,:) * mskisf_cav(:,:)
          !
          ! output
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfdiags.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfdiags.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfdiags.F90	(revision 11521)
@@ -17,4 +17,5 @@
    USE dom_oce
    USE isf            ! ice shelf variable
+   USE isfutils
    USE iom            ! 
 
@@ -93,4 +94,6 @@
       ! compute 3d output
       zvar2d(:,:) = pvar2d(:,:) / phtbl(:,:)
+      zvar3d(:,:,:) = 0._wp
+      !
       DO jj = 1,jpj
          DO ji = 1,jpi
@@ -98,12 +101,10 @@
             ikb = kbot(ji,jj)
             DO jk = ikt, ikb - 1
-               zvar3d(ji,jj,jk) = zvar3d(ji,jj,jk) + zvar2d(ji,jj) * e3t_n(ji,jj,jk)
+               zvar3d(ji,jj,jk) = zvar2d(ji,jj) * e3t_n(ji,jj,jk)
             END DO
-            ! phtbl(ji,jj) * pfrac(ji,jj) = thickness of bottom affected by tbl
-            zvar3d(ji,jj,ikb) = zvar3d(ji,jj,ikb) + zvar2d(ji,jj) * e3t_n(ji,jj,ikb) * pfrac(ji,jj)
+            zvar3d(ji,jj,ikb) = zvar2d(ji,jj) * e3t_n(ji,jj,ikb) * pfrac(ji,jj)
          END DO
       END DO
       !
-      ! output data
       CALL iom_put( TRIM(cdvar) , zvar3d(:,:,:))
       !
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfhdiv.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfhdiv.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfhdiv.F90	(revision 11521)
@@ -2,5 +2,7 @@
 
    USE dom_oce
+   USE iom
    USE isf
+   USE isfutils
    USE phycst
    USE in_out_manager
@@ -28,9 +30,16 @@
       INTEGER, INTENT(in) :: kt
       !
-      ! ice shelf cavity contribution
-      IF ( ln_isfcav_mlt ) CALL isf_hdiv_mlt(misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav, fwfisf_cav, fwfisf_cav_b, phdiv)
+      IF (lwp) WRITE(numout,*) ' isfdiv start'
+      FLUSH(numout)
       !
-      ! ice shelf parametrisation contribution
-      IF ( ln_isfpar_mlt ) CALL isf_hdiv_mlt(misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par, fwfisf_par, fwfisf_par_b, phdiv)
+      IF ( ln_isf ) THEN
+         !
+         ! ice shelf cavity contribution
+         IF ( ln_isfcav_mlt ) CALL isf_hdiv_mlt(misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav, fwfisf_cav, fwfisf_cav_b, phdiv)
+         !
+         ! ice shelf parametrisation contribution
+         IF ( ln_isfpar_mlt ) CALL isf_hdiv_mlt(misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par, fwfisf_par, fwfisf_par_b, phdiv)
+         !
+      END IF
       !
       ! ice sheet coupling contribution (if conservation needed)
@@ -38,10 +47,12 @@
          !
          ! correct divergence only for the first time step
-         IF ( kt == nit000 ) CALL isf_hdiv_cpl(risfcpl_vol, phdiv)
+         !IF ( kt == nit000 ) CALL isf_hdiv_cpl(risfcpl_vol, phdiv)
          !
          ! correct divergence every time step to remove any trend due to coupling
-         IF ( ll_isfcpl_cons ) CALL isf_hdiv_cpl(risfcpl_cons_vol, phdiv)
+         IF ( ll_isfcpl_cons ) CALL isf_hdiv_cpl(risfcpl_cons_vol, risfcpl_cons_vol_b, phdiv)
          !
       END IF
+      IF (lwp) WRITE(numout,*) ' isfdiv end'
+      FLUSH(numout)
       !
    END SUBROUTINE isf_hdiv
@@ -78,12 +89,12 @@
       DO jj = 1,jpj
          DO ji = 1,jpi
-               ikt = ktop(ji,jj)
-               ikb = kbot(ji,jj)
-               ! level fully include in the ice shelf boundary layer
-               DO jk = ikt, ikb - 1
-                  phdiv(ji,jj,jk) = phdiv(ji,jj,jk) + zqvol(ji,jj)
-               END DO
-               ! level partially include in ice shelf boundary layer 
-               phdiv(ji,jj,ikb) = phdiv(ji,jj,ikb) + zqvol(ji,jj) * pfrac(ji,jj)
+            ikt = ktop(ji,jj)
+            ikb = kbot(ji,jj)
+            ! level fully include in the ice shelf boundary layer
+            DO jk = ikt, ikb - 1
+               phdiv(ji,jj,jk) = phdiv(ji,jj,jk) + zqvol(ji,jj)
+            END DO
+            ! level partially include in ice shelf boundary layer 
+            phdiv(ji,jj,ikb) = phdiv(ji,jj,ikb) + zqvol(ji,jj) * pfrac(ji,jj)
          END DO
       END DO
@@ -91,15 +102,26 @@
    END SUBROUTINE isf_hdiv_mlt
 
-   SUBROUTINE isf_hdiv_cpl(pqvol, phdiv)
+   SUBROUTINE isf_hdiv_cpl(pqvol, pqvol_b, phdiv)
+      !!----------------------------------------------------------------------
+      !!                  ***  SUBROUTINE isf_hdiv_cpl  ***
+      !!       
+      !! ** Purpose :   
+      !!
+      !! ** Method  :   
+      !!
+      !! ** Action  :  
+      !!----------------------------------------------------------------------
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: phdiv
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) :: pqvol
-   
-   INTEGER :: jk
-      
+      !!----------------------------------------------------------------------
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in   ) :: pqvol, pqvol_b
+      !!----------------------------------------------------------------------
+      INTEGER :: jk
+      !!----------------------------------------------------------------------
+      !
       DO jk=1,jpk 
          phdiv(:,:,jk) =  phdiv(:,:,jk) + pqvol(:,:,jk) * r1_e1e2t(:,:) / e3t_n(:,:,jk)
       END DO
-
-   END SUBROUTINE
+      !
+   END SUBROUTINE isf_hdiv_cpl
 
 END MODULE isfhdiv
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfnxt.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfnxt.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfnxt.F90	(revision 11521)
@@ -42,7 +42,4 @@
       IF ( ln_isfpar_mlt ) CALL isf_dynnxt_mlt(misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par, fwfisf_par, fwfisf_par_b, pcoef)
       !
-      ! ice sheet coupling
-      !IF ( ln_iscpl_hsb  ) CALL isf_dynnxt_cpl(...zcoef)
-      !
    END SUBROUTINE isf_dynnxt
 
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfparmlt.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfparmlt.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfparmlt.F90	(revision 11521)
@@ -148,8 +148,8 @@
       ! 3. ------------BG03 output
       ! output ttbl
-      CALL iom_put('ttbl_par', ztavg(:,:) )
+      CALL iom_put('ttbl_par', ztavg(:,:) * mskisf_par(:,:) )
       !
       ! output thermal driving
-      CALL iom_put('isfthermald_par',( ztfrz(:,:) - ztavg(:,:) ))
+      CALL iom_put('isfthermald_par',( ztfrz(:,:) - ztavg(:,:) ) * mskisf_par(:,:))
       !
       !
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfstp.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfstp.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfstp.F90	(revision 11521)
@@ -33,4 +33,5 @@
    USE isfcpl         ! isf variables
    USE isf            ! isf variables
+   USE isfutils
 
    IMPLICIT NONE
@@ -103,6 +104,8 @@
       END IF
 
-      IF ( ln_isfcpl ) THEN
+      IF ( ll_isfcpl ) THEN
+
          IF (lrst_oce) CALL isfcpl_rst_write(kt)
+
       END IF
       !
@@ -181,11 +184,4 @@
             WRITE(numout,*) ''
             !
-            WRITE(numout,*) '      Coupling to an ice sheet model          ln_isfcpl         = ', ln_isfcpl
-            IF ( ln_isfcpl ) THEN
-               WRITE(numout,*) '         conservation activated ln_isfcpl_cons           = ', ln_isfcpl_cons
-               WRITE(numout,*) '            number of call of the extrapolation loop = ', nn_drown
-            ENDIF
-            !
-            WRITE(numout,*) '      Ice shelf load method                   cn_isfload        = ', TRIM(cn_isfload)
          ELSE
             IF ( ln_isfcav ) THEN
@@ -196,4 +192,13 @@
          END IF
 
+         WRITE(numout,*) '      Coupling to an ice sheet model          ln_isfcpl         = ', ln_isfcpl
+         IF ( ln_isfcpl ) THEN
+            WRITE(numout,*) '         conservation activated ln_isfcpl_cons           = ', ln_isfcpl_cons
+            WRITE(numout,*) '            number of call of the extrapolation loop = ', nn_drown
+         ENDIF
+         !
+         IF (ln_isfcav) WRITE(numout,*) '      Ice shelf load method                   cn_isfload        = ', TRIM(cn_isfload)
+         WRITE(numout,*) ''
+
       END IF
       !
@@ -244,37 +249,39 @@
       !
       ! terminate routine now if no ice shelf melt formulation specify
-      IF ( .NOT. ln_isf ) RETURN
-      !
-      ! initialisation useful variable
-      r1_Lfusisf =  1._wp / rLfusisf
+      IF ( ln_isf ) THEN
+         !
+         ! initialisation useful variable
+         r1_Lfusisf =  1._wp / rLfusisf
+         !
+         ! initialisation melt in the cavity
+         IF ( ln_isfcav_mlt ) THEN
+            !
+            ! initialisation  of cav variable
+            CALL isf_cav_init()
+            !
+            ! read cav variable from restart
+            IF ( ln_rstart ) CALL isfrst_read('cav', risf_cav_tsc, fwfisf_cav, risf_cav_tsc_b, fwfisf_cav_b)
+            !
+         END IF
+         !
+         !---------------------------------------------------------------------------------------------------------------------
+         ! initialisation parametrised melt
+         IF ( ln_isfpar_mlt ) THEN
+            !
+            ! initialisation  of par variable
+            CALL isf_par_init()
+            !
+            ! read par variable from restart
+            IF ( ln_rstart ) CALL isfrst_read('par', risf_par_tsc, fwfisf_par, risf_par_tsc_b, fwfisf_par_b)
+            !
+         END IF
+      END IF
+      !
+      !---------------------------------------------------------------------------------------------------------------------
+      ! initialisation ice sheet coupling
       !
       ll_isfcpl     = .FALSE.
       ll_isfcpl_cons= .FALSE.
       !
-      ! initialisation melt in the cavity
-      IF ( ln_isfcav_mlt ) THEN
-         !
-         ! initialisation  of cav variable
-         CALL isf_cav_init()
-         !
-         ! read cav variable from restart
-         IF ( ln_rstart ) CALL isfrst_read('cav', risf_cav_tsc, fwfisf_cav, risf_cav_tsc_b, fwfisf_cav_b)
-         !
-      END IF
-      !
-      !---------------------------------------------------------------------------------------------------------------------
-      ! initialisation parametrised melt
-      IF ( ln_isfpar_mlt ) THEN
-         !
-         ! initialisation  of par variable
-         CALL isf_par_init()
-         !
-         ! read par variable from restart
-         IF ( ln_rstart ) CALL isfrst_read('par', risf_par_tsc, fwfisf_par, risf_par_tsc_b, fwfisf_par_b)
-         !
-      END IF
-      !
-      !---------------------------------------------------------------------------------------------------------------------
-      ! initialisation ice sheet coupling
       IF( ln_isfcpl ) THEN
 
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isftbl.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isftbl.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isftbl.F90	(revision 11521)
@@ -98,7 +98,4 @@
          !
       END SELECT
-      !
-      ! mask mean tbl value
-      pvarout(:,:) = pvarout(:,:) * ssmask(:,:)
       !
    END SUBROUTINE isf_tbl
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfutils.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfutils.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/ISF/isfutils.F90	(revision 11521)
@@ -65,5 +65,5 @@
       !!--------------------------------------------------------------------
       !
-      zsum = glob_sum( 'debug', pvar(:,:) )
+      zsum = glob_sum_full( 'debug', pvar(:,:) )
       zmin = MINVAL( pvar(:,:) ) ; CALL mpp_min( 'debug', zmin ) ! min over the global domain
       zmax = MAXVAL( pvar(:,:) ) ; CALL mpp_max( 'debug', zmax ) ! max over the global domain
@@ -89,7 +89,7 @@
       !!--------------------------------------------------------------------
       !
-      zsum = glob_sum( 'debug', pvar(:,:) )
-      zmin = MINVAL( pvar(:,:) ) ; CALL mpp_min( 'debug', zmin ) ! min over the global domain
-      zmax = MAXVAL( pvar(:,:) ) ; CALL mpp_max( 'debug', zmax ) ! max over the global domain
+      zsum = glob_sum_full( 'debug', pvar(:,:,:) )
+      zmin = MINVAL( pvar(:,:,:) ) ; CALL mpp_min( 'debug', zmin ) ! min over the global domain
+      zmax = MAXVAL( pvar(:,:,:) ) ; CALL mpp_max( 'debug', zmax ) ! max over the global domain
       !
       WRITE(numout,*) TRIM(cdtxt),' (min, max, sum) : ',zmin, zmax, zsum
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/SBC/sbc_oce.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/SBC/sbc_oce.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/SBC/sbc_oce.F90	(revision 11521)
@@ -43,5 +43,4 @@
    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)
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/SBC/sbccpl.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/SBC/sbccpl.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/SBC/sbccpl.F90	(revision 11521)
@@ -36,5 +36,5 @@
    USE eosbn2         ! 
    USE sbcrnf  , ONLY : l_rnfcpl
-   USE isf     , ONLY : l_isfoasis, fwfisf_oasis
+   USE isf     , ONLY : ln_isf, l_isfoasis, fwfisf_oasis
 #if defined key_cice
    USE ice_domain_size, only: ncat
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/TRA/traisf.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/TRA/traisf.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/TRA/traisf.F90	(revision 11521)
@@ -48,15 +48,23 @@
       IF( ln_timing )   CALL timing_start('tra_sbc')
       !
-      ! cavity case
-      IF ( ln_isfcav_mlt ) CALL tra_isf_mlt(misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav, risf_cav_tsc, risf_cav_tsc_b, tsa)
-      !
-      ! parametrisation case
-      IF ( ln_isfpar_mlt ) CALL tra_isf_mlt(misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par, risf_par_tsc, risf_par_tsc_b, tsa)
+      IF (ln_isf) THEN
+         !
+         ! cavity case
+         IF ( ln_isfcav_mlt ) CALL tra_isf_mlt(misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav, risf_cav_tsc, risf_cav_tsc_b, tsa)
+         !
+         ! parametrisation case
+         IF ( ln_isfpar_mlt ) CALL tra_isf_mlt(misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par, risf_par_tsc, risf_par_tsc_b, tsa)
+         !
+      END IF
       !
       ! ice sheet coupling case
-      IF ( ll_isfcpl .AND. kt == nit000 ) CALL tra_isf_cpl(risfcpl_tsc, tsa)
-      !
-      ! ice sheet coupling case (
-      IF ( ll_isfcpl_cons ) CALL tra_isf_cpl(risfcpl_cons_tsc, tsa)
+      IF ( ll_isfcpl ) THEN
+         !
+         IF ( kt == nit000 ) CALL tra_isf_cpl(risfcpl_tsc, tsa)
+         !
+         ! ensure 0 trend due to unconservation of the ice shelf coupling
+         IF ( ll_isfcpl_cons ) CALL tra_isf_cpl(risfcpl_cons_tsc, tsa)
+         !
+      END IF
       !
       IF( ln_timing )   CALL timing_stop('tra_isf')
@@ -95,5 +103,4 @@
             !
             ! level fully include in the ice shelf boundary layer
-            ! sign - because fwf sign of evapo (rnf sign of precip)
             DO jk = ikt, ikb - 1
                pts(ji,jj,jk,jp_tem) = pts(ji,jj,jk,jp_tem) + ztc(ji,jj)
@@ -117,5 +124,5 @@
       REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(inout) :: ptsa
       !!----------------------------------------------------------------------
-      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in   ) :: ptsc  !( >0 out )
+      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT(in   ) :: ptsc
       !!----------------------------------------------------------------------
       INTEGER :: jk
Index: NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/step.F90
===================================================================
--- NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/step.F90	(revision 11495)
+++ NEMO/branches/2019/ENHANCE-02_ISF_nemo/src/OCE/step.F90	(revision 11521)
@@ -39,4 +39,5 @@
    !
    USE iom              ! xIOs server
+   USE isf
 
    IMPLICIT NONE
@@ -113,5 +114,5 @@
       IF( ln_apr_dyn )   CALL sbc_apr ( kstp )                   ! atmospheric pressure (NB: call before bdy_dta which needs ssh_ib) 
       IF( ln_bdy     )   CALL bdy_dta ( kstp, time_offset=+1 )   ! update dynamic & tracer data at open boundaries
-      IF( ln_isf     )   CALL isf_stp ( kstp )
+      IF( ln_isf     .OR. ll_isfcpl )   CALL isf_stp ( kstp )
                          CALL sbc     ( kstp )                   ! Sea Boundary Condition (including sea-ice)
 
@@ -241,5 +242,5 @@
                          CALL tra_sbc       ( kstp )  ! surface boundary condition
       IF( ln_traqsr  )   CALL tra_qsr       ( kstp )  ! penetrative solar radiation qsr
-      IF( ln_isf     )   CALL tra_isf       ( kstp )  ! ice shelf heat flux
+      IF( ln_isf .OR. ll_isfcpl ) CALL tra_isf       ( kstp )  ! ice shelf heat flux
       IF( ln_trabbc  )   CALL tra_bbc       ( kstp )  ! bottom heat flux
       IF( ln_trabbl  )   CALL tra_bbl       ( kstp )  ! advective (and/or diffusive) bottom boundary layer scheme
