Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/CONFIG/SHARED/namelist_ref
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/CONFIG/SHARED/namelist_ref	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/CONFIG/SHARED/namelist_ref	(revision 8563)
@@ -564,6 +564,6 @@
 ! for lim3
 !   bn_a_i     = 'amm12_bdyT_ice',         24        , 'ileadfra',    .true.   , .false. ,  'daily'  ,    ''    ,   ''     ,     ''
-!   bn_ht_i    = 'amm12_bdyT_ice',         24        , 'iicethic',    .true.   , .false. ,  'daily'  ,    ''    ,   ''     ,     ''
-!   bn_ht_s    = 'amm12_bdyT_ice',         24        , 'isnowthi',    .true.   , .false. ,  'daily'  ,    ''    ,   ''     ,     ''
+!   bn_h_i    = 'amm12_bdyT_ice',         24        , 'iicethic',    .true.   , .false. ,  'daily'  ,    ''    ,   ''     ,     ''
+!   bn_h_s    = 'amm12_bdyT_ice',         24        , 'isnowthi',    .true.   , .false. ,  'daily'  ,    ''    ,   ''     ,     ''
 
    cn_dir      = 'bdydta/' !  root directory for the location of the bulk files
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/ice.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/ice.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/ice.F90	(revision 8563)
@@ -118,6 +118,6 @@
    !!-------------|-------------|---------------------------------|-------|
    !!                                                                     |
-   !! ht_i        | ht_i_1d     |    Ice thickness                | m     |
-   !! ht_s        ! ht_s_1d     |    Snow depth                   | m     |
+   !! h_i        | h_i_1d     |    Ice thickness                | m     |
+   !! h_s        ! h_s_1d     |    Snow depth                   | m     |
    !! sm_i        ! sm_i_1d     |    Sea ice bulk salinity        ! ppt   |
    !! s_i         ! s_i_1d      |    Sea ice salinity profile     ! ppt   |
@@ -312,9 +312,9 @@
    !!----------------------------------------------------------------------
    !! Variables defined for each ice category
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   ht_i      !: Ice thickness (m)
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   h_i      !: Ice thickness (m)
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   a_i       !: Ice fractional areas (concentration)
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   v_i       !: Ice volume per unit area (m)
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   v_s       !: Snow volume per unit area(m)
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   ht_s      !: Snow thickness (m)
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   h_s      !: Snow thickness (m)
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   t_su      !: Sea-Ice Surface Temperature (K)
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sm_i      !: Sea-Ice Bulk salinity (ppt)
@@ -375,5 +375,5 @@
    !! * Old values of global variables
    !!----------------------------------------------------------------------
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   v_s_b, v_i_b, ht_s_b, ht_i_b  !: snow and ice volumes/thickness
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   v_s_b, v_i_b, h_s_b, h_i_b  !: snow and ice volumes/thickness
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   a_i_b, smv_i_b, oa_i_b        !:
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   e_s_b                         !: snow heat content
@@ -463,6 +463,6 @@
       ii = ii + 1
       ALLOCATE( ftr_ice(jpi,jpj,jpl) ,                                                 &
-         &      ht_i   (jpi,jpj,jpl) , a_i   (jpi,jpj,jpl) , v_i   (jpi,jpj,jpl) ,     &
-         &      v_s    (jpi,jpj,jpl) , ht_s  (jpi,jpj,jpl) , t_su  (jpi,jpj,jpl) ,     &
+         &      h_i   (jpi,jpj,jpl) , a_i   (jpi,jpj,jpl) , v_i   (jpi,jpj,jpl) ,     &
+         &      v_s    (jpi,jpj,jpl) , h_s  (jpi,jpj,jpl) , t_su  (jpi,jpj,jpl) ,     &
          &      sm_i   (jpi,jpj,jpl) , smv_i (jpi,jpj,jpl) , o_i   (jpi,jpj,jpl) ,     &
          &      oa_i   (jpi,jpj,jpl) , bv_i  (jpi,jpj,jpl) ,  STAT=ierr(ii) )
@@ -512,5 +512,5 @@
       ! * Old values of global variables
       ii = ii + 1
-      ALLOCATE( v_s_b  (jpi,jpj,jpl) , v_i_b  (jpi,jpj,jpl) , ht_s_b(jpi,jpj,jpl)        , ht_i_b(jpi,jpj,jpl)        ,   &
+      ALLOCATE( v_s_b  (jpi,jpj,jpl) , v_i_b  (jpi,jpj,jpl) , h_s_b(jpi,jpj,jpl)        , h_i_b(jpi,jpj,jpl)        ,   &
          &      a_i_b  (jpi,jpj,jpl) , smv_i_b(jpi,jpj,jpl) , e_i_b (jpi,jpj,nlay_i,jpl) , e_s_b (jpi,jpj,nlay_s,jpl) ,   &
          &      oa_i_b (jpi,jpj,jpl)                                                     , STAT=ierr(ii) )
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/ice1D.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/ice1D.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/ice1D.F90	(revision 8563)
@@ -105,7 +105,7 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   a_i_1d        !: <==> the 2D  a_i
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   a_ib_1d       !: <==> the 2D  a_i_b
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   ht_i_1d       !:
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   ht_ib_1d      !:
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   ht_s_1d       !:
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   h_i_1d       !:
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   h_ib_1d      !:
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   h_s_1d       !:
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   fc_su         !: Surface Conduction flux 
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) ::   fc_bo_i       !: Bottom  Conduction flux 
@@ -153,8 +153,8 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   v_ip_2d 
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   t_su_2d 
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ht_i_2d
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   h_i_2d
 
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   a_ib_2d
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ht_ib_2d
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   h_ib_2d
    
    !!----------------------------------------------------------------------
@@ -200,5 +200,5 @@
       ii = ii + 1
       ALLOCATE( t_su_1d  (jpij) , t_si_1d   (jpij) , a_i_1d  (jpij) , a_ib_1d(jpij) ,                  &
-         &      ht_i_1d  (jpij) , ht_ib_1d  (jpij) , ht_s_1d (jpij) , fc_su  (jpij) , fc_bo_i(jpij) ,  &    
+         &      h_i_1d  (jpij) , h_ib_1d  (jpij) , h_s_1d (jpij) , fc_su  (jpij) , fc_bo_i(jpij) ,  &    
          &      dh_s_tot (jpij) , dh_i_surf (jpij) , dh_i_sub(jpij) ,                                  &    
          &      dh_i_bott(jpij) , dh_snowice(jpij) , sm_i_1d (jpij) , s_i_new(jpij) , &
@@ -218,5 +218,5 @@
       !
       ii = ii + 1
-      ALLOCATE( a_i_2d(jpij,jpl) , a_ib_2d(jpij,jpl) , ht_i_2d(jpij,jpl) , ht_ib_2d(jpij,jpl) , &
+      ALLOCATE( a_i_2d(jpij,jpl) , a_ib_2d(jpij,jpl) , h_i_2d(jpij,jpl) , h_ib_2d(jpij,jpl) , &
          &      v_i_2d(jpij,jpl) ,v_s_2d(jpij,jpl) ,oa_i_2d(jpij,jpl) ,smv_i_2d(jpij,jpl) ,  &
          &      a_ip_2d(jpij,jpl) ,v_ip_2d(jpij,jpl) ,t_su_2d(jpij,jpl) ,  &
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icecor.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icecor.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icecor.F90	(revision 8563)
@@ -70,8 +70,8 @@
          !                          !  thickness of the smallest category above himin    !
          !                          !-----------------------------------------------------
-         WHERE( a_i(:,:,1) >= epsi20 )   ;   ht_i(:,:,1) = v_i (:,:,1) / a_i(:,:,1)
-         ELSEWHERE                       ;   ht_i(:,:,1) = 0._wp
+         WHERE( a_i(:,:,1) >= epsi20 )   ;   h_i(:,:,1) = v_i (:,:,1) / a_i(:,:,1)
+         ELSEWHERE                       ;   h_i(:,:,1) = 0._wp
          END WHERE
-         WHERE( ht_i(:,:,1) < rn_himin )     a_i (:,:,1) = a_i (:,:,1) * ht_i(:,:,1) / rn_himin
+         WHERE( h_i(:,:,1) < rn_himin )      a_i(:,:,1) = a_i (:,:,1) * h_i(:,:,1) / rn_himin
          !
       ENDIF
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icectl.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icectl.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icectl.F90	(revision 8563)
@@ -248,5 +248,5 @@
       DO jj = 1, jpj
          DO ji = 1, jpi
-            IF(   ht_i(ji,jj,jl)  >  50._wp   ) THEN
+            IF(   h_i(ji,jj,jl)  >  50._wp   ) THEN
                !CALL ice_prt( kt, ji, jj, 2, ' ALERTE 3 :   Very thick ice ' )
                inb_alp(ialert_id) = inb_alp(ialert_id) + 1
@@ -469,6 +469,6 @@
                   WRITE(numout,*) ' - Category (', jl,')'
                   WRITE(numout,*) ' a_i           : ', a_i(ji,jj,jl)
-                  WRITE(numout,*) ' ht_i          : ', ht_i(ji,jj,jl)
-                  WRITE(numout,*) ' ht_s          : ', ht_s(ji,jj,jl)
+                  WRITE(numout,*) ' h_i           : ', h_i(ji,jj,jl)
+                  WRITE(numout,*) ' h_s           : ', h_s(ji,jj,jl)
                   WRITE(numout,*) ' v_i           : ', v_i(ji,jj,jl)
                   WRITE(numout,*) ' v_s           : ', v_s(ji,jj,jl)
@@ -520,5 +520,5 @@
                   WRITE(numout,*) ' - Category (',jl,')'
                   WRITE(numout,*) '   ~~~~~~~~         ' 
-                  WRITE(numout,*) ' ht_i       : ', ht_i(ji,jj,jl)             , ' ht_s       : ', ht_s(ji,jj,jl)
+                  WRITE(numout,*) ' h_i        : ', h_i(ji,jj,jl)              , ' h_s        : ', h_s(ji,jj,jl)
                   WRITE(numout,*) ' t_i        : ', t_i(ji,jj,1:nlay_i,jl)
                   WRITE(numout,*) ' t_su       : ', t_su(ji,jj,jl)             , ' t_s        : ', t_s(ji,jj,1,jl)
@@ -630,6 +630,6 @@
          CALL prt_ctl_info(' - Category : ', ivar1=jl)
          CALL prt_ctl_info('   ~~~~~~~~~~')
-         CALL prt_ctl(tab2d_1=ht_i       (:,:,jl)        , clinfo1= ' ht_i        : ')
-         CALL prt_ctl(tab2d_1=ht_s       (:,:,jl)        , clinfo1= ' ht_s        : ')
+         CALL prt_ctl(tab2d_1=h_i        (:,:,jl)        , clinfo1= ' h_i        : ')
+         CALL prt_ctl(tab2d_1=h_s        (:,:,jl)        , clinfo1= ' h_s        : ')
          CALL prt_ctl(tab2d_1=t_su       (:,:,jl)        , clinfo1= ' t_su        : ')
          CALL prt_ctl(tab2d_1=t_s        (:,:,1,jl)      , clinfo1= ' t_snow      : ')
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icedyn.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icedyn.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icedyn.F90	(revision 8563)
@@ -85,15 +85,15 @@
          tau_icebfr(:,:) = 0._wp
          DO jl = 1, jpl
-            WHERE( ht_i(:,:,jl) > ht_n(:,:) * rn_gamma )   tau_icebfr(:,:) = tau_icebfr(:,:) + a_i(:,:,jl) * rn_icebfr
+            WHERE( h_i(:,:,jl) > ht_n(:,:) * rn_gamma )   tau_icebfr(:,:) = tau_icebfr(:,:) + a_i(:,:,jl) * rn_icebfr
          END DO
          IF( iom_use('tau_icebfr') )   CALL iom_put( 'tau_icebfr', tau_icebfr )  
       ENDIF
 
-      zhmax(:,:,:) = ht_i_b(:,:,:)      !-- Record max of the surrounding 9-pts ice thick. (for CALL Hbig)
+      zhmax(:,:,:) = h_i_b(:,:,:)      !-- Record max of the surrounding 9-pts ice thick. (for CALL Hbig)
       DO jl = 1, jpl
          DO jj = 2, jpjm1
             DO ji = 2, jpim1
 !!gm use of MAXVAL here is very probably less efficient than expending the 9 values
-               zhmax(ji,jj,jl) = MAX( epsi20, MAXVAL( ht_i_b(ji-1:ji+1,jj-1:jj+1,jl) ) )
+               zhmax(ji,jj,jl) = MAX( epsi20, MAXVAL( h_i_b(ji-1:ji+1,jj-1:jj+1,jl) ) )
             END DO
          END DO
@@ -166,5 +166,5 @@
                   IF ( ( zdv >  0.0 .AND. zh > phmax(ji,jj,jl) .AND. at_i_b(ji,jj) < 0.80 ) .OR. &
                      & ( zdv <= 0.0 .AND. zh > phmax(ji,jj,jl) ) ) THEN
-                     a_i (ji,jj,jl) = v_i(ji,jj,jl) / MIN( phmax(ji,jj,jl), hi_max(jpl) )   !-- bound ht_i to hi_max (99 m)
+                     a_i (ji,jj,jl) = v_i(ji,jj,jl) / MIN( phmax(ji,jj,jl), hi_max(jpl) )   !-- bound h_i to hi_max (99 m)
                   ENDIF
                   !
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceforcing.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceforcing.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceforcing.F90	(revision 8563)
@@ -130,5 +130,5 @@
 
       ! --- cloud-sky and overcast-sky ice albedos --- !
-      CALL ice_alb( t_su, ht_i, ht_s, a_ip_frac, h_ip, ln_pnd_rad, zalb_cs, zalb_os )
+      CALL ice_alb( t_su, h_i, h_s, a_ip_frac, h_ip, ln_pnd_rad, zalb_cs, zalb_os )
 
       ! albedo depends on cloud fraction because of non-linear spectral effects
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceistate.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceistate.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceistate.F90	(revision 8563)
@@ -285,5 +285,5 @@
                DO ji = 1, jpi
                   a_i(ji,jj,jl)   = zswitch(ji,jj) * za_i_ini(ji,jj,jl)                       ! concentration
-                  ht_i(ji,jj,jl)  = zswitch(ji,jj) * zh_i_ini(ji,jj,jl)                       ! ice thickness
+                  h_i(ji,jj,jl)  = zswitch(ji,jj) * zh_i_ini(ji,jj,jl)                       ! ice thickness
                   sm_i(ji,jj,jl)  = zswitch(ji,jj) * zsm_i_ini(ji,jj)                         ! salinity
                   o_i(ji,jj,jl)   = 0._wp                                                     ! age (0 day)
@@ -291,20 +291,20 @@
 
                   IF( zht_i_ini(ji,jj) > 0._wp )THEN
-                    ht_s(ji,jj,jl)= ht_i(ji,jj,jl) * ( zht_s_ini(ji,jj) / zht_i_ini(ji,jj) )  ! snow depth
+                    h_s(ji,jj,jl)= h_i(ji,jj,jl) * ( zht_s_ini(ji,jj) / zht_i_ini(ji,jj) )  ! snow depth
                   ELSE
-                    ht_s(ji,jj,jl)= 0._wp
+                    h_s(ji,jj,jl)= 0._wp
                   ENDIF
 
-                  ! This case below should not be used if (ht_s/ht_i) is ok in namelist
+                  ! This case below should not be used if (h_s/h_i) is ok in namelist
                   ! In case snow load is in excess that would lead to transformation from snow to ice
                   ! Then, transfer the snow excess into the ice (different from icethd_dh)
-                  zdh = MAX( 0._wp, ( rhosn * ht_s(ji,jj,jl) + ( rhoic - rau0 ) * ht_i(ji,jj,jl) ) * r1_rau0 ) 
-                  ! recompute ht_i, ht_s avoiding out of bounds values
-                  ht_i(ji,jj,jl) = MIN( hi_max(jl), ht_i(ji,jj,jl) + zdh )
-                  ht_s(ji,jj,jl) = MAX( 0._wp, ht_s(ji,jj,jl) - zdh * rhoic * r1_rhosn )
+                  zdh = MAX( 0._wp, ( rhosn * h_s(ji,jj,jl) + ( rhoic - rau0 ) * h_i(ji,jj,jl) ) * r1_rau0 ) 
+                  ! recompute h_i, h_s avoiding out of bounds values
+                  h_i(ji,jj,jl) = MIN( hi_max(jl), h_i(ji,jj,jl) + zdh )
+                  h_s(ji,jj,jl) = MAX( 0._wp, h_s(ji,jj,jl) - zdh * rhoic * r1_rhosn )
 
                   ! ice volume, salt content, age content
-                  v_i(ji,jj,jl)   = ht_i(ji,jj,jl) * a_i(ji,jj,jl)              ! ice volume
-                  v_s(ji,jj,jl)   = ht_s(ji,jj,jl) * a_i(ji,jj,jl)              ! snow volume
+                  v_i(ji,jj,jl)   = h_i(ji,jj,jl) * a_i(ji,jj,jl)              ! ice volume
+                  v_s(ji,jj,jl)   = h_s(ji,jj,jl) * a_i(ji,jj,jl)              ! snow volume
                   smv_i(ji,jj,jl) = MIN( sm_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl) ! salt content
                   oa_i(ji,jj,jl)  = o_i(ji,jj,jl) * a_i(ji,jj,jl)               ! age content
@@ -381,6 +381,6 @@
          smv_i(:,:,:) = 0._wp
          oa_i (:,:,:) = 0._wp
-         ht_i (:,:,:) = 0._wp
-         ht_s (:,:,:) = 0._wp
+         h_i  (:,:,:) = 0._wp
+         h_s  (:,:,:) = 0._wp
          sm_i (:,:,:) = 0._wp
          o_i  (:,:,:) = 0._wp
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceitd.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceitd.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceitd.F90	(revision 8563)
@@ -101,6 +101,6 @@
          zhbnew(:,:) = 0._wp
 
-         CALL tab_3d_2d( nidx, idxice(1:nidx), ht_i_2d (1:nidx,1:jpl), ht_i   )
-         CALL tab_3d_2d( nidx, idxice(1:nidx), ht_ib_2d(1:nidx,1:jpl), ht_i_b )
+         CALL tab_3d_2d( nidx, idxice(1:nidx), h_i_2d (1:nidx,1:jpl), h_i   )
+         CALL tab_3d_2d( nidx, idxice(1:nidx), h_ib_2d(1:nidx,1:jpl), h_i_b )
          CALL tab_3d_2d( nidx, idxice(1:nidx), a_i_2d  (1:nidx,1:jpl), a_i    )
          CALL tab_3d_2d( nidx, idxice(1:nidx), a_ib_2d (1:nidx,1:jpl), a_i_b  )
@@ -109,5 +109,5 @@
             ! Compute thickness change in each ice category
             DO ji = 1, nidx
-               zdhice(ji,jl) = ht_i_2d(ji,jl) - ht_ib_2d(ji,jl)
+               zdhice(ji,jl) = h_i_2d(ji,jl) - h_ib_2d(ji,jl)
             END DO
          END DO
@@ -119,9 +119,9 @@
                !
                ! --- New boundary: Hn* = Hn + Fn*dt --- !
-               !     Fn*dt = ( fn + (fn+1 - fn)/(hn+1 - hn) * (Hn - hn) ) * dt = zdhice + zslope * (Hmax - ht_i_b)
+               !     Fn*dt = ( fn + (fn+1 - fn)/(hn+1 - hn) * (Hn - hn) ) * dt = zdhice + zslope * (Hmax - h_i_b)
                !
                IF    ( a_ib_2d(ji,jl) >  epsi10 .AND. a_ib_2d(ji,jl+1) >  epsi10 ) THEN   ! a(jl+1) & a(jl) /= 0
-                  zslope           = ( zdhice(ji,jl+1) - zdhice(ji,jl) ) / ( ht_ib_2d(ji,jl+1) - ht_ib_2d(ji,jl) )
-                  zhbnew(ji,jl) = hi_max(jl) + zdhice(ji,jl) + zslope * ( hi_max(jl) - ht_ib_2d(ji,jl) )
+                  zslope           = ( zdhice(ji,jl+1) - zdhice(ji,jl) ) / ( h_ib_2d(ji,jl+1) - h_ib_2d(ji,jl) )
+                  zhbnew(ji,jl) = hi_max(jl) + zdhice(ji,jl) + zslope * ( hi_max(jl) - h_ib_2d(ji,jl) )
                ELSEIF( a_ib_2d(ji,jl) >  epsi10 .AND. a_ib_2d(ji,jl+1) <= epsi10 ) THEN   ! a(jl+1)=0 => Hn* = Hn + fn*dt
                   zhbnew(ji,jl) = hi_max(jl) + zdhice(ji,jl)
@@ -136,6 +136,6 @@
                !    Note: hn(t+1) must not be too close to either HR or HL otherwise a division by nearly 0 is possible 
                !          in ice_itd_glinear in the case (HR-HL) = 3(Hice - HL) or = 3(HR - Hice)
-               IF( a_i_2d(ji,jl  ) > epsi10 .AND. ht_i_2d(ji,jl  ) > ( zhbnew(ji,jl) - epsi10 ) )   idxice(ji) = 0
-               IF( a_i_2d(ji,jl+1) > epsi10 .AND. ht_i_2d(ji,jl+1) < ( zhbnew(ji,jl) + epsi10 ) )   idxice(ji) = 0
+               IF( a_i_2d(ji,jl  ) > epsi10 .AND. h_i_2d(ji,jl  ) > ( zhbnew(ji,jl) - epsi10 ) )   idxice(ji) = 0
+               IF( a_i_2d(ji,jl+1) > epsi10 .AND. h_i_2d(ji,jl+1) < ( zhbnew(ji,jl) + epsi10 ) )   idxice(ji) = 0
                
                ! 2) Hn-1 < Hn* < Hn+1  
@@ -149,5 +149,5 @@
          DO ji = 1, nidx
             IF( a_i_2d(ji,jpl) > epsi10 ) THEN
-               zhbnew(ji,jpl) = MAX( hi_max(jpl-1), 3._wp * ht_i_2d(ji,jpl) - 2._wp * zhbnew(ji,jpl-1) )
+               zhbnew(ji,jpl) = MAX( hi_max(jpl-1), 3._wp * h_i_2d(ji,jpl) - 2._wp * zhbnew(ji,jpl-1) )
             ELSE
                zhbnew(ji,jpl) = hi_max(jpl)  
@@ -158,6 +158,6 @@
             !    h1(t) must not be too close to either HR or HL otherwise a division by nearly 0 is possible 
             !    in ice_itd_glinear in the case (HR-HL) = 3(Hice - HL) or = 3(HR - Hice)
-            IF( ht_ib_2d(ji,1) < ( hi_max(0) + epsi10 ) )   idxice(ji) = 0
-            IF( ht_ib_2d(ji,1) > ( hi_max(1) - epsi10 ) )   idxice(ji) = 0
+            IF( h_ib_2d(ji,1) < ( hi_max(0) + epsi10 ) )   idxice(ji) = 0
+            IF( h_ib_2d(ji,1) > ( hi_max(1) - epsi10 ) )   idxice(ji) = 0
          END DO
          !
@@ -189,6 +189,6 @@
          DO jl = 1, jpl
             !
-            CALL tab_2d_1d( nidx, idxice(1:nidx), ht_ib_1d(1:nidx), ht_i_b(:,:,jl) )
-            CALL tab_2d_1d( nidx, idxice(1:nidx), ht_i_1d (1:nidx), ht_i(:,:,jl)   )
+            CALL tab_2d_1d( nidx, idxice(1:nidx), h_ib_1d(1:nidx), h_i_b(:,:,jl) )
+            CALL tab_2d_1d( nidx, idxice(1:nidx), h_i_1d (1:nidx), h_i(:,:,jl)   )
             CALL tab_2d_1d( nidx, idxice(1:nidx), a_i_1d  (1:nidx), a_i(:,:,jl)    )
             CALL tab_2d_1d( nidx, idxice(1:nidx), v_i_1d  (1:nidx), v_i(:,:,jl)    )
@@ -197,5 +197,5 @@
                !  
                ! --- g(h) for category 1 --- !
-               CALL ice_itd_glinear( zhb0(1:nidx)  , zhb1(1:nidx)  , ht_ib_1d(1:nidx)  , a_i_1d(1:nidx)  ,  &   ! in
+               CALL ice_itd_glinear( zhb0(1:nidx)  , zhb1(1:nidx)  , h_ib_1d(1:nidx)  , a_i_1d(1:nidx)  ,  &   ! in
                   &                  g0  (1:nidx,1), g1  (1:nidx,1), hL      (1:nidx,1), hR    (1:nidx,1)   )   ! out
                   !
@@ -205,5 +205,5 @@
                   IF( a_i_1d(ji) > epsi10 ) THEN
                      !
-                     zdh0 =  ht_i_1d(ji) - ht_ib_1d(ji)                
+                     zdh0 =  h_i_1d(ji) - h_ib_1d(ji)                
                      IF( zdh0 < 0.0 ) THEN      !remove area from category 1
                         zdh0 = MIN( -zdh0, hi_max(1) )
@@ -215,11 +215,11 @@
                            zx2    = 0.5 * zetamax * zetamax 
                            zda0   = g1(ji,1) * zx2 + g0(ji,1) * zx1                        ! ice area removed
-                           zdamax = a_i_1d(ji) * (1.0 - ht_i_1d(ji) / ht_ib_1d(ji) ) ! Constrain new thickness <= ht_i                
+                           zdamax = a_i_1d(ji) * (1.0 - h_i_1d(ji) / h_ib_1d(ji) ) ! Constrain new thickness <= h_i                
                            zda0   = MIN( zda0, zdamax )                                                  ! ice area lost due to melting 
                            !     of thin ice (zdamax > 0)
                            ! Remove area, conserving volume
-                           ht_i_1d(ji) = ht_i_1d(ji) * a_i_1d(ji) / ( a_i_1d(ji) - zda0 )
+                           h_i_1d(ji) = h_i_1d(ji) * a_i_1d(ji) / ( a_i_1d(ji) - zda0 )
                            a_i_1d(ji)  = a_i_1d(ji) - zda0
-                           v_i_1d(ji)  = a_i_1d(ji) * ht_i_1d(ji) ! clem-useless ?
+                           v_i_1d(ji)  = a_i_1d(ji) * h_i_1d(ji) ! clem-useless ?
                         ENDIF
                         !
@@ -233,5 +233,5 @@
                END DO
                !
-               CALL tab_1d_2d( nidx, idxice(1:nidx), ht_i_1d (1:nidx), ht_i(:,:,jl)   )
+               CALL tab_1d_2d( nidx, idxice(1:nidx), h_i_1d (1:nidx), h_i(:,:,jl)   )
                CALL tab_1d_2d( nidx, idxice(1:nidx), a_i_1d  (1:nidx), a_i(:,:,jl)    )
                CALL tab_1d_2d( nidx, idxice(1:nidx), v_i_1d  (1:nidx), v_i(:,:,jl)    )
@@ -240,5 +240,5 @@
             !
             ! --- g(h) for each thickness category --- !  
-            CALL ice_itd_glinear( zhbnew(1:nidx,jl-1), zhbnew(1:nidx,jl), ht_i_1d(1:nidx)   , a_i_1d(1:nidx)   ,  &   ! in
+            CALL ice_itd_glinear( zhbnew(1:nidx,jl-1), zhbnew(1:nidx,jl), h_i_1d(1:nidx)   , a_i_1d(1:nidx)   ,  &   ! in
                &                  g0    (1:nidx,jl  ), g1    (1:nidx,jl), hL     (1:nidx,jl), hR    (1:nidx,jl)   )   ! out
             !
@@ -284,23 +284,23 @@
          
          !----------------------------------------------------------------------------------------------
-         ! 7) Make sure ht_i >= minimum ice thickness hi_min
+         ! 7) Make sure h_i >= minimum ice thickness hi_min
          !----------------------------------------------------------------------------------------------
-         CALL tab_2d_1d( nidx, idxice(1:nidx), ht_i_1d (1:nidx), ht_i(:,:,1)   )
+         CALL tab_2d_1d( nidx, idxice(1:nidx), h_i_1d (1:nidx), h_i(:,:,1)   )
          CALL tab_2d_1d( nidx, idxice(1:nidx), a_i_1d  (1:nidx), a_i(:,:,1)    )
          CALL tab_2d_1d( nidx, idxice(1:nidx), a_ip_1d (1:nidx), a_ip(:,:,1)   )
          
          DO ji = 1, nidx
-            IF ( a_i_1d(ji) > epsi10 .AND. ht_i_1d(ji) < rn_himin ) THEN
-               a_i_1d (ji) = a_i_1d(ji) * ht_i_1d(ji) / rn_himin 
+            IF ( a_i_1d(ji) > epsi10 .AND. h_i_1d(ji) < rn_himin ) THEN
+               a_i_1d (ji) = a_i_1d(ji) * h_i_1d(ji) / rn_himin 
                ! MV MP 2016
                IF ( nn_pnd_scheme > 0 ) THEN
-                  a_ip_1d(ji) = a_ip_1d(ji) * ht_i_1d(ji) / rn_himin
+                  a_ip_1d(ji) = a_ip_1d(ji) * h_i_1d(ji) / rn_himin
                ENDIF
                ! END MV MP 2016
-               ht_i_1d(ji) = rn_himin
+               h_i_1d(ji) = rn_himin
             ENDIF
          END DO
          !
-         CALL tab_1d_2d( nidx, idxice(1:nidx), ht_i_1d (1:nidx), ht_i(:,:,1)   )
+         CALL tab_1d_2d( nidx, idxice(1:nidx), h_i_1d (1:nidx), h_i(:,:,1)   )
          CALL tab_1d_2d( nidx, idxice(1:nidx), a_i_1d  (1:nidx), a_i(:,:,1)    )
          CALL tab_1d_2d( nidx, idxice(1:nidx), a_ip_1d (1:nidx), a_ip(:,:,1)   )
@@ -392,5 +392,5 @@
       !!------------------------------------------------------------------
          
-      CALL tab_3d_2d( nidx, idxice(1:nidx), ht_i_2d (1:nidx,1:jpl), ht_i   )
+      CALL tab_3d_2d( nidx, idxice(1:nidx), h_i_2d (1:nidx,1:jpl), h_i   )
       CALL tab_3d_2d( nidx, idxice(1:nidx), a_i_2d  (1:nidx,1:jpl), a_i    )
       CALL tab_3d_2d( nidx, idxice(1:nidx), v_i_2d  (1:nidx,1:jpl), v_i    )
@@ -518,12 +518,12 @@
       !-------------------------------------------------------------------------------
       WHERE( a_i_2d(1:nidx,:) >= epsi20 )
-         ht_i_2d(1:nidx,:)  =  v_i_2d(1:nidx,:) / a_i_2d(1:nidx,:) 
+         h_i_2d(1:nidx,:)  =  v_i_2d(1:nidx,:) / a_i_2d(1:nidx,:) 
          t_su_2d(1:nidx,:)  =  zaTsfn(1:nidx,:) / a_i_2d(1:nidx,:) 
       ELSEWHERE
-         ht_i_2d(1:nidx,:)  = 0._wp
+         h_i_2d(1:nidx,:)  = 0._wp
          t_su_2d(1:nidx,:)  = rt0
       END WHERE
       !
-      CALL tab_2d_3d( nidx, idxice(1:nidx), ht_i_2d (1:nidx,1:jpl), ht_i  )
+      CALL tab_2d_3d( nidx, idxice(1:nidx), h_i_2d (1:nidx,1:jpl), h_i  )
       CALL tab_2d_3d( nidx, idxice(1:nidx), a_i_2d  (1:nidx,1:jpl), a_i   )
       CALL tab_2d_3d( nidx, idxice(1:nidx), v_i_2d  (1:nidx,1:jpl), v_i   )
@@ -574,5 +574,5 @@
          END DO
          !
-!!clem   CALL tab_2d_1d( nidx, idxice(1:nidx), ht_i_1d (1:nidx), ht_i(:,:,jl)   )
+!!clem   CALL tab_2d_1d( nidx, idxice(1:nidx), h_i_1d (1:nidx), h_i(:,:,jl)   )
          CALL tab_2d_1d( nidx, idxice(1:nidx), a_i_1d  (1:nidx), a_i(:,:,jl)    )
          CALL tab_2d_1d( nidx, idxice(1:nidx), v_i_1d  (1:nidx), v_i(:,:,jl)    )
@@ -582,5 +582,5 @@
             ! how much of a_i you send in cat sup is somewhat arbitrary
 !!clem: these do not work properly after a restart (I do not know why)
-!!          zdaice(ji,jl)  = a_i_1d(ji) * ( ht_i_1d(ji) - hi_max(jl) + epsi10 ) / ht_i_1d(ji)  
+!!          zdaice(ji,jl)  = a_i_1d(ji) * ( h_i_1d(ji) - hi_max(jl) + epsi10 ) / h_i_1d(ji)  
 !!          zdvice(ji,jl)  = v_i_1d(ji) - ( a_i_1d(ji) - zdaice(ji,jl) ) * ( hi_max(jl) - epsi10 )
 !!clem: these do not work properly after a restart (I do not know why)
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icestp.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icestp.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icestp.F90	(revision 8563)
@@ -146,5 +146,5 @@
          !
          !                          !==  previous lead fraction and ice volume for flux calculations
-         CALL ice_var_glo2eqv            ! ht_i and ht_s for ice albedo calculation
+         CALL ice_var_glo2eqv            ! h_i and h_s for ice albedo calculation
          CALL ice_var_agg(1)             ! at_i for coupling 
          CALL store_fields               ! Store now ice values
@@ -358,9 +358,9 @@
       e_i_b  (:,:,:,:) = e_i  (:,:,:,:)   ! ice thermal energy
       WHERE( a_i_b(:,:,:) >= epsi20 )
-         ht_i_b(:,:,:) = v_i_b (:,:,:) / a_i_b(:,:,:)   ! ice thickness
-         ht_s_b(:,:,:) = v_s_b (:,:,:) / a_i_b(:,:,:)   ! snw thickness
+         h_i_b(:,:,:) = v_i_b (:,:,:) / a_i_b(:,:,:)   ! ice thickness
+         h_s_b(:,:,:) = v_s_b (:,:,:) / a_i_b(:,:,:)   ! snw thickness
       ELSEWHERE
-         ht_i_b(:,:,:) = 0._wp
-         ht_s_b(:,:,:) = 0._wp
+         h_i_b(:,:,:) = 0._wp
+         h_s_b(:,:,:) = 0._wp
       END WHERE
       
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd.F90	(revision 8563)
@@ -302,5 +302,5 @@
             
             ! mask temperature
-            rswitch       = 1._wp - MAX( 0._wp , SIGN( 1._wp , - ht_i_1d(ji) ) ) 
+            rswitch       = 1._wp - MAX( 0._wp , SIGN( 1._wp , - h_i_1d(ji) ) ) 
             t_i_1d(ji,jk) = rswitch * t_i_1d(ji,jk) + ( 1._wp - rswitch ) * rt0
          END DO 
@@ -326,14 +326,14 @@
          zdh_mel = MIN( 0._wp, dh_i_surf(ji) + dh_i_bott(ji) + dh_snowice(ji) + dh_i_sub(ji) )
          IF( zdh_mel < 0._wp .AND. a_i_1d(ji) > 0._wp )  THEN
-            zvi          = a_i_1d(ji) * ht_i_1d(ji)
-            zvs          = a_i_1d(ji) * ht_s_1d(ji)
+            zvi          = a_i_1d(ji) * h_i_1d(ji)
+            zvs          = a_i_1d(ji) * h_s_1d(ji)
             ! lateral melting = concentration change
-            zhi_bef     = ht_i_1d(ji) - zdh_mel
+            zhi_bef     = h_i_1d(ji) - zdh_mel
             rswitch     = MAX( 0._wp , SIGN( 1._wp , zhi_bef - epsi20 ) )
             zda_mel     = rswitch * a_i_1d(ji) * zdh_mel / ( 2._wp * MAX( zhi_bef, epsi20 ) )
             a_i_1d(ji)  = MAX( epsi20, a_i_1d(ji) + zda_mel ) 
             ! adjust thickness
-            ht_i_1d(ji) = zvi / a_i_1d(ji)            
-            ht_s_1d(ji) = zvs / a_i_1d(ji)            
+            h_i_1d(ji) = zvi / a_i_1d(ji)            
+            h_s_1d(ji) = zvs / a_i_1d(ji)            
             ! retrieve total concentration
             at_i_1d(ji) = a_i_1d(ji)
@@ -362,6 +362,6 @@
          CALL tab_2d_1d( nidx, idxice(1:nidx), at_i_1d(1:nidx), at_i             )
          CALL tab_2d_1d( nidx, idxice(1:nidx), a_i_1d (1:nidx), a_i (:,:,kl)     )
-         CALL tab_2d_1d( nidx, idxice(1:nidx), ht_i_1d(1:nidx), ht_i(:,:,kl)     )
-         CALL tab_2d_1d( nidx, idxice(1:nidx), ht_s_1d(1:nidx), ht_s(:,:,kl)     )
+         CALL tab_2d_1d( nidx, idxice(1:nidx), h_i_1d(1:nidx), h_i(:,:,kl)     )
+         CALL tab_2d_1d( nidx, idxice(1:nidx), h_s_1d(1:nidx), h_s(:,:,kl)     )
          CALL tab_2d_1d( nidx, idxice(1:nidx), t_su_1d(1:nidx), t_su(:,:,kl)     )
          CALL tab_2d_1d( nidx, idxice(1:nidx), sm_i_1d(1:nidx), sm_i(:,:,kl)     )
@@ -436,8 +436,8 @@
          ! --- Change units of e_i, e_s from J/m2 to J/m3 --- !
          DO jk = 1, nlay_i
-            WHERE( ht_i_1d(1:nidx)>0._wp ) e_i_1d(1:nidx,jk) = e_i_1d(1:nidx,jk) / (ht_i_1d(1:nidx) * a_i_1d(1:nidx)) * nlay_i
+            WHERE( h_i_1d(1:nidx)>0._wp ) e_i_1d(1:nidx,jk) = e_i_1d(1:nidx,jk) / (h_i_1d(1:nidx) * a_i_1d(1:nidx)) * nlay_i
          END DO
          DO jk = 1, nlay_s
-            WHERE( ht_s_1d(1:nidx)>0._wp ) e_s_1d(1:nidx,jk) = e_s_1d(1:nidx,jk) / (ht_s_1d(1:nidx) * a_i_1d(1:nidx)) * nlay_s
+            WHERE( h_s_1d(1:nidx)>0._wp ) e_s_1d(1:nidx,jk) = e_s_1d(1:nidx,jk) / (h_s_1d(1:nidx) * a_i_1d(1:nidx)) * nlay_s
          END DO
          !
@@ -447,19 +447,19 @@
          ! --- Change units of e_i, e_s from J/m3 to J/m2 --- !
          DO jk = 1, nlay_i
-            e_i_1d(1:nidx,jk) = e_i_1d(1:nidx,jk) * ht_i_1d(1:nidx) * a_i_1d(1:nidx) * r1_nlay_i
+            e_i_1d(1:nidx,jk) = e_i_1d(1:nidx,jk) * h_i_1d(1:nidx) * a_i_1d(1:nidx) * r1_nlay_i
          END DO
          DO jk = 1, nlay_s
-            e_s_1d(1:nidx,jk) = e_s_1d(1:nidx,jk) * ht_s_1d(1:nidx) * a_i_1d(1:nidx) * r1_nlay_s
+            e_s_1d(1:nidx,jk) = e_s_1d(1:nidx,jk) * h_s_1d(1:nidx) * a_i_1d(1:nidx) * r1_nlay_s
          END DO
          !
          ! Change thickness to volume (replaces routine ice_var_eqv2glo)
-         v_i_1d(1:nidx)   = ht_i_1d(1:nidx) * a_i_1d(1:nidx)
-         v_s_1d(1:nidx)   = ht_s_1d(1:nidx) * a_i_1d(1:nidx)
+         v_i_1d(1:nidx)   = h_i_1d(1:nidx) * a_i_1d(1:nidx)
+         v_s_1d(1:nidx)   = h_s_1d(1:nidx) * a_i_1d(1:nidx)
          smv_i_1d(1:nidx) = sm_i_1d(1:nidx) * v_i_1d(1:nidx)
          
          CALL tab_1d_2d( nidx, idxice(1:nidx), at_i_1d(1:nidx), at_i             )
          CALL tab_1d_2d( nidx, idxice(1:nidx), a_i_1d (1:nidx), a_i (:,:,kl)     )
-         CALL tab_1d_2d( nidx, idxice(1:nidx), ht_i_1d(1:nidx), ht_i(:,:,kl)     )
-         CALL tab_1d_2d( nidx, idxice(1:nidx), ht_s_1d(1:nidx), ht_s(:,:,kl)     )
+         CALL tab_1d_2d( nidx, idxice(1:nidx), h_i_1d(1:nidx), h_i(:,:,kl)     )
+         CALL tab_1d_2d( nidx, idxice(1:nidx), h_s_1d(1:nidx), h_s(:,:,kl)     )
          CALL tab_1d_2d( nidx, idxice(1:nidx), t_su_1d(1:nidx), t_su(:,:,kl)     )
          CALL tab_1d_2d( nidx, idxice(1:nidx), sm_i_1d(1:nidx), sm_i(:,:,kl)     )
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_da.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_da.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_da.F90	(revision 8563)
@@ -135,20 +135,20 @@
             
             ! Contribution to salt flux
-            sfx_lam_1d(ji) = sfx_lam_1d(ji) + rhoic *  ht_i_1d(ji) * zda * sm_i_1d(ji) * r1_rdtice
+            sfx_lam_1d(ji) = sfx_lam_1d(ji) + rhoic *  h_i_1d(ji) * zda * sm_i_1d(ji) * r1_rdtice
             
             ! Contribution to heat flux into the ocean [W.m-2], (<0)  
-            hfx_thd_1d(ji) = hfx_thd_1d(ji) - zda * r1_rdtice * ( ht_i_1d(ji) * r1_nlay_i * SUM( e_i_1d(ji,1:nlay_i) )  &
-                                                                + ht_s_1d(ji) * r1_nlay_s * SUM( e_s_1d(ji,1:nlay_s) ) ) 
+            hfx_thd_1d(ji) = hfx_thd_1d(ji) - zda * r1_rdtice * ( h_i_1d(ji) * r1_nlay_i * SUM( e_i_1d(ji,1:nlay_i) )  &
+                                                                + h_s_1d(ji) * r1_nlay_s * SUM( e_s_1d(ji,1:nlay_s) ) ) 
             
             ! Contribution to mass flux
-            wfx_lam_1d(ji) =  wfx_lam_1d(ji) + zda * r1_rdtice * ( rhoic * ht_i_1d(ji) + rhosn * ht_s_1d(ji) )
+            wfx_lam_1d(ji) =  wfx_lam_1d(ji) + zda * r1_rdtice * ( rhoic * h_i_1d(ji) + rhosn * h_s_1d(ji) )
             
             ! new concentration
             a_i_1d(ji) = a_i_1d(ji) - zda
 
-            ! ensure that ht_i = 0 where a_i = 0
+            ! ensure that h_i = 0 where a_i = 0
             IF( a_i_1d(ji) == 0._wp ) THEN
-               ht_i_1d(ji) = 0._wp
-               ht_s_1d(ji) = 0._wp
+               h_i_1d(ji) = 0._wp
+               h_s_1d(ji) = 0._wp
             ENDIF
          ENDIF
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_dh.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_dh.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_dh.F90	(revision 8563)
@@ -121,5 +121,5 @@
       DO jk = 1, nlay_i
          DO ji = 1, nidx
-            h_i_old (ji,jk) = ht_i_1d(ji) * r1_nlay_i
+            h_i_old (ji,jk) = h_i_1d(ji) * r1_nlay_i
             eh_i_old(ji,jk) = e_i_1d(ji,jk) * h_i_old(ji,jk)
          END DO
@@ -145,9 +145,9 @@
          IF( t_s_1d(ji,1) > rt0 ) THEN !!! Internal melting
             ! Contribution to heat flux to the ocean [W.m-2], < 0  
-            hfx_res_1d(ji) = hfx_res_1d(ji) + e_s_1d(ji,1) * ht_s_1d(ji) * a_i_1d(ji) * r1_rdtice
+            hfx_res_1d(ji) = hfx_res_1d(ji) + e_s_1d(ji,1) * h_s_1d(ji) * a_i_1d(ji) * r1_rdtice
             ! Contribution to mass flux
-            wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) + rhosn * ht_s_1d(ji) * a_i_1d(ji) * r1_rdtice
+            wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) + rhosn * h_s_1d(ji) * a_i_1d(ji) * r1_rdtice
             ! updates
-            ht_s_1d(ji)   = 0._wp
+            h_s_1d(ji)   = 0._wp
             e_s_1d (ji,1) = 0._wp
             t_s_1d (ji,1) = rt0
@@ -160,5 +160,5 @@
       DO jk = 1, nlay_i
          DO ji = 1, nidx
-            zh_i(ji,jk) = ht_i_1d(ji) * r1_nlay_i
+            zh_i(ji,jk) = h_i_1d(ji) * r1_nlay_i
             zeh_i(ji)   = zeh_i(ji) + e_i_1d(ji,jk) * zh_i(ji,jk)
          END DO
@@ -216,5 +216,5 @@
 
          ! update thickness
-         ht_s_1d(ji) = MAX( 0._wp , ht_s_1d(ji) + zdh_s_pre(ji) )
+         h_s_1d(ji) = MAX( 0._wp , h_s_1d(ji) + zdh_s_pre(ji) )
       END DO
 
@@ -225,8 +225,8 @@
          DO ji = 1, nidx
             ! thickness change
-            rswitch          = 1._wp - MAX( 0._wp, SIGN( 1._wp, - ht_s_1d(ji) ) ) 
+            rswitch          = 1._wp - MAX( 0._wp, SIGN( 1._wp, - h_s_1d(ji) ) ) 
             rswitch          = rswitch * ( MAX( 0._wp, SIGN( 1._wp, e_s_1d(ji,jk) - epsi20 ) ) ) 
             zdeltah  (ji,jk) = - rswitch * zq_su(ji) / MAX( e_s_1d(ji,jk), epsi20 )
-            zdeltah  (ji,jk) = MAX( zdeltah(ji,jk) , - ht_s_1d(ji) ) ! bound melting
+            zdeltah  (ji,jk) = MAX( zdeltah(ji,jk) , - h_s_1d(ji) ) ! bound melting
             zdh_s_mel(ji)    = zdh_s_mel(ji) + zdeltah(ji,jk)    
             ! heat used to melt snow(W.m-2, >0)
@@ -236,5 +236,5 @@
             ! updates available heat + thickness
             zq_su (ji)  = MAX( 0._wp , zq_su (ji) + zdeltah(ji,jk) * e_s_1d(ji,jk) )
-            ht_s_1d(ji) = MAX( 0._wp , ht_s_1d(ji) + zdeltah(ji,jk) )
+            h_s_1d(ji) = MAX( 0._wp , h_s_1d(ji) + zdeltah(ji,jk) )
          END DO
       END DO
@@ -247,5 +247,5 @@
       zdeltah(1:nidx,:) = 0._wp
       DO ji = 1, nidx
-         zdh_s_sub(ji)  = MAX( - ht_s_1d(ji) , - evap_ice_1d(ji) * r1_rhosn * rdt_ice )
+         zdh_s_sub(ji)  = MAX( - h_s_1d(ji) , - evap_ice_1d(ji) * r1_rhosn * rdt_ice )
          ! remaining evap in kg.m-2 (used for ice melting later on)
          zevap_rema(ji)  = evap_ice_1d(ji) * rdt_ice + zdh_s_sub(ji) * rhosn
@@ -258,5 +258,5 @@
 
          ! new snow thickness
-         ht_s_1d(ji)    =  MAX( 0._wp , ht_s_1d(ji) + zdh_s_sub(ji) )
+         h_s_1d(ji)    =  MAX( 0._wp , h_s_1d(ji) + zdh_s_sub(ji) )
          ! update precipitations after sublimation and correct sublimation
          zdh_s_pre(ji) = zdh_s_pre(ji) + zdeltah(ji,1)
@@ -275,8 +275,8 @@
       DO jk = 1, nlay_s
          DO ji = 1,nidx
-            rswitch       = MAX( 0._wp , SIGN( 1._wp, ht_s_1d(ji) - epsi20 ) )
-            e_s_1d(ji,jk) = rswitch / MAX( ht_s_1d(ji), epsi20 ) *           &
+            rswitch       = MAX( 0._wp , SIGN( 1._wp, h_s_1d(ji) - epsi20 ) )
+            e_s_1d(ji,jk) = rswitch / MAX( h_s_1d(ji), epsi20 ) *           &
               &            ( ( zdh_s_pre(ji)               ) * zqprec(ji) +  &
-              &              ( ht_s_1d(ji) - zdh_s_pre(ji) ) * rhosn * ( cpic * ( rt0 - t_s_1d(ji,jk) ) + lfus ) )
+              &              ( h_s_1d(ji) - zdh_s_pre(ji) ) * rhosn * ( cpic * ( rt0 - t_s_1d(ji,jk) ) + lfus ) )
          END DO
       END DO
@@ -376,5 +376,5 @@
       ! update ice thickness
       DO ji = 1, nidx
-         ht_i_1d(ji) =  MAX( 0._wp , ht_i_1d(ji) + dh_i_surf(ji) + dh_i_sub(ji) )
+         h_i_1d(ji) =  MAX( 0._wp , h_i_1d(ji) + dh_i_surf(ji) + dh_i_sub(ji) )
       END DO
 
@@ -555,5 +555,5 @@
       !-------------------------------------------
       DO ji = 1, nidx
-         ht_i_1d(ji) =  MAX( 0._wp , ht_i_1d(ji) + dh_i_bott(ji) )
+         h_i_1d(ji) =  MAX( 0._wp , h_i_1d(ji) + dh_i_bott(ji) )
       END DO  
 
@@ -566,10 +566,10 @@
       DO ji = 1, nidx
          zq_rema(ji)     = zq_su(ji) + zq_bo(ji) 
-         rswitch         = 1._wp - MAX( 0._wp, SIGN( 1._wp, - ht_s_1d(ji) ) )   ! =1 if snow
+         rswitch         = 1._wp - MAX( 0._wp, SIGN( 1._wp, - h_s_1d(ji) ) )   ! =1 if snow
          rswitch         = rswitch * MAX( 0._wp, SIGN( 1._wp, e_s_1d(ji,1) - epsi20 ) )
          zdeltah  (ji,1) = - rswitch * zq_rema(ji) / MAX( e_s_1d(ji,1), epsi20 )
-         zdeltah  (ji,1) = MIN( 0._wp , MAX( zdeltah(ji,1) , - ht_s_1d(ji) ) ) ! bound melting
+         zdeltah  (ji,1) = MIN( 0._wp , MAX( zdeltah(ji,1) , - h_s_1d(ji) ) ) ! bound melting
          dh_s_tot (ji)   = dh_s_tot(ji)  + zdeltah(ji,1)
-         ht_s_1d   (ji)  = ht_s_1d(ji)   + zdeltah(ji,1)
+         h_s_1d   (ji)  = h_s_1d(ji)   + zdeltah(ji,1)
         
          zq_rema(ji)     = zq_rema(ji) + zdeltah(ji,1) * e_s_1d(ji,1)                ! update available heat (J.m-2)
@@ -594,8 +594,8 @@
       DO ji = 1, nidx
          !
-         dh_snowice(ji) = MAX(  0._wp , ( rhosn * ht_s_1d(ji) + (rhoic-rau0) * ht_i_1d(ji) ) * z1_rho )
-
-         ht_i_1d(ji)    = ht_i_1d(ji) + dh_snowice(ji)
-         ht_s_1d(ji)    = ht_s_1d(ji) - dh_snowice(ji)
+         dh_snowice(ji) = MAX(  0._wp , ( rhosn * h_s_1d(ji) + (rhoic-rau0) * h_i_1d(ji) ) * z1_rho )
+
+         h_i_1d(ji)    = h_i_1d(ji) + dh_snowice(ji)
+         h_s_1d(ji)    = h_s_1d(ji) - dh_snowice(ji)
 
          ! Contribution to energy flux to the ocean [J/m2], >0 (if sst<0)
@@ -632,5 +632,5 @@
       !-------------------------------------------
       DO ji = 1, nidx
-         rswitch     =  1.0 - MAX( 0._wp , SIGN( 1._wp , - ht_i_1d(ji) ) ) 
+         rswitch     =  1.0 - MAX( 0._wp , SIGN( 1._wp , - h_i_1d(ji) ) ) 
          t_su_1d(ji) =  rswitch * t_su_1d(ji) + ( 1.0 - rswitch ) * rt0
       END DO
@@ -639,5 +639,5 @@
          DO ji = 1,nidx
             ! mask enthalpy
-            rswitch       = 1._wp - MAX(  0._wp , SIGN( 1._wp, - ht_s_1d(ji) )  )
+            rswitch       = 1._wp - MAX(  0._wp , SIGN( 1._wp, - h_s_1d(ji) )  )
             e_s_1d(ji,jk) = rswitch * e_s_1d(ji,jk)
             ! recalculate t_s_1d from e_s_1d
@@ -646,6 +646,6 @@
       END DO
 
-      ! --- ensure that a_i = 0 where ht_i = 0 ---
-      WHERE( ht_i_1d(1:nidx) == 0._wp )   a_i_1d(1:nidx) = 0._wp
+      ! --- ensure that a_i = 0 where h_i = 0 ---
+      WHERE( h_i_1d(1:nidx) == 0._wp )   a_i_1d(1:nidx) = 0._wp
       !
    END SUBROUTINE ice_thd_dh
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_do.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_do.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_do.F90	(revision 8563)
@@ -76,5 +76,5 @@
       !!             - Computation of variation of ice volume and mass
       !!             - Computation of a_i after lateral accretion and 
-      !!               update ht_s_1d, ht_i_1d and tbif_1d(:,:)      
+      !!               update h_s_1d, h_i_1d      
       !!------------------------------------------------------------------------
       INTEGER  ::   ji,jj,jk,jl      ! dummy loop indices
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_sal.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_sal.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_sal.F90	(revision 8563)
@@ -75,8 +75,8 @@
             !  Update ice salinity from snow-ice and bottom growth
             !---------------------------------------------------------
-            IF( ht_i_1d(ji) > 0._wp ) THEN
+            IF( h_i_1d(ji) > 0._wp ) THEN
                zs_sni   = sss_1d(ji) * ( rhoic - rhosn ) * r1_rhoic                                 ! Salinity of snow ice
-               zsm_i_si = ( zs_sni      - sm_i_1d(ji) ) *             dh_snowice(ji)  / ht_i_1d(ji) ! snow-ice    
-               zsm_i_bg = ( s_i_new(ji) - sm_i_1d(ji) ) * MAX( 0._wp, dh_i_bott(ji) ) / ht_i_1d(ji) ! bottom growth
+               zsm_i_si = ( zs_sni      - sm_i_1d(ji) ) *             dh_snowice(ji)  / h_i_1d(ji) ! snow-ice    
+               zsm_i_bg = ( s_i_new(ji) - sm_i_1d(ji) ) * MAX( 0._wp, dh_i_bott(ji) ) / h_i_1d(ji) ! bottom growth
                ! Update salinity (nb: salt flux already included in icethd_dh)
                sm_i_1d(ji) = sm_i_1d(ji) + zsm_i_bg + zsm_i_si
@@ -97,5 +97,5 @@
                
                ! Salt flux
-               sfx_bri_1d(ji) = sfx_bri_1d(ji) - rhoic * a_i_1d(ji) * ht_i_1d(ji) * ( zsm_i_fl + zsm_i_gd ) * r1_rdtice
+               sfx_bri_1d(ji) = sfx_bri_1d(ji) - rhoic * a_i_1d(ji) * h_i_1d(ji) * ( zsm_i_fl + zsm_i_gd ) * r1_rdtice
             ENDIF
          END DO
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_zdf.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_zdf.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_zdf.F90	(revision 8563)
@@ -81,5 +81,5 @@
       !!           ice salinities          : s_i_1d
       !!           number of layers in the ice/snow: nlay_i, nlay_s
-      !!           total ice/snow thickness : ht_i_1d, ht_s_1d
+      !!           total ice/snow thickness : h_i_1d, h_s_1d
       !!-------------------------------------------------------------------
       INTEGER ::   ji, jk         ! spatial loop index
@@ -145,6 +145,6 @@
       ! --- diag error on heat diffusion - PART 1 --- !
       DO ji = 1, nidx
-         zq_ini(ji) = ( SUM( e_i_1d(ji,1:nlay_i) ) * ht_i_1d(ji) * r1_nlay_i +  &
-            &           SUM( e_s_1d(ji,1:nlay_s) ) * ht_s_1d(ji) * r1_nlay_s ) 
+         zq_ini(ji) = ( SUM( e_i_1d(ji,1:nlay_i) ) * h_i_1d(ji) * r1_nlay_i +  &
+            &           SUM( e_s_1d(ji,1:nlay_s) ) * h_s_1d(ji) * r1_nlay_s ) 
       END DO
 
@@ -153,8 +153,8 @@
       !------------------
       DO ji = 1, nidx
-         isnow(ji)= 1._wp - MAX( 0._wp , SIGN(1._wp, - ht_s_1d(ji) ) )  ! is there snow or not
+         isnow(ji)= 1._wp - MAX( 0._wp , SIGN(1._wp, - h_s_1d(ji) ) )  ! is there snow or not
          ! layer thickness
-         zh_i(ji) = ht_i_1d(ji) * r1_nlay_i
-         zh_s(ji) = ht_s_1d(ji) * r1_nlay_s
+         zh_i(ji) = h_i_1d(ji) * r1_nlay_i
+         zh_s(ji) = h_s_1d(ji) * r1_nlay_s
       END DO
       !
@@ -187,5 +187,5 @@
          ! ftr_ice = io.qsr_ice.exp(-k(h_i)) transmitted below the ice 
          ! fr1_i0_1d = i0 for a thin ice cover, fr1_i0_2d = i0 for a thick ice cover
-         zfac = MAX( 0._wp , 1._wp - ( ht_i_1d(ji) * z1_hsu ) )     
+         zfac = MAX( 0._wp , 1._wp - ( h_i_1d(ji) * z1_hsu ) )     
          i0(ji) = ( 1._wp - isnow(ji) ) * ( fr1_i0_1d(ji) + zfac * fr2_i0_1d(ji) )
 
@@ -279,5 +279,5 @@
             DO ji = 1, nidx
                zcnd_i = SUM( ztcond_i(ji,:) ) / REAL( nlay_i+1, wp )                             ! Mean sea ice thermal conductivity
-               zhe = ( rn_cnd_s * ht_i_1d(ji) + zcnd_i * ht_s_1d(ji) ) / ( rn_cnd_s + zcnd_i )   ! Effective thickness he (zhe)
+               zhe = ( rn_cnd_s * h_i_1d(ji) + zcnd_i * h_s_1d(ji) ) / ( rn_cnd_s + zcnd_i )   ! Effective thickness he (zhe)
                IF( zhe >=  zepsilon * 0.5_wp * EXP(1._wp) ) THEN
                   zghe(ji) = MIN( 2._wp, 0.5_wp * ( 1._wp + LOG( 2._wp * zhe / zepsilon ) ) )    ! G(he)
@@ -382,5 +382,5 @@
          DO ji = 1, nidx
             !                               !---------------------!
-            IF ( ht_s_1d(ji) > 0.0 ) THEN   !  snow-covered cells !
+            IF ( h_s_1d(ji) > 0.0 ) THEN   !  snow-covered cells !
                !                            !---------------------!
                ! snow interior terms (bottom equation has the same form as the others)
@@ -524,5 +524,5 @@
          DO ji = 1, nidx
             ! snow temperatures      
-            IF( ht_s_1d(ji) > 0._wp ) THEN
+            IF( h_s_1d(ji) > 0._wp ) THEN
                t_s_1d(ji,nlay_s) = ( zindtbis(ji,nlay_s+1) - ztrid(ji,nlay_s+1,3) * t_i_1d(ji,1) )  &
                   &                / zdiagbis(ji,nlay_s+1)
@@ -598,6 +598,6 @@
       !     zhfx_err = correction on the diagnosed heat flux due to non-convergence of the algorithm used to solve heat equation
       DO ji = 1, nidx
-         zdq = - zq_ini(ji) + ( SUM( e_i_1d(ji,1:nlay_i) ) * ht_i_1d(ji) * r1_nlay_i +  &
-            &                   SUM( e_s_1d(ji,1:nlay_s) ) * ht_s_1d(ji) * r1_nlay_s )
+         zdq = - zq_ini(ji) + ( SUM( e_i_1d(ji,1:nlay_i) ) * h_i_1d(ji) * r1_nlay_i +  &
+            &                   SUM( e_s_1d(ji,1:nlay_s) ) * h_s_1d(ji) * r1_nlay_s )
 
          IF( t_su_1d(ji) < rt0 ) THEN  ! case T_su < 0degC
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceupdate.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceupdate.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceupdate.F90	(revision 8563)
@@ -199,5 +199,5 @@
       ! Snow/ice albedo (only if sent to coupler, useless in forced mode)
       !------------------------------------------------------------------
-      CALL ice_alb( t_su, ht_i, ht_s, a_ip_frac, h_ip, ln_pnd_rad, zalb_cs, zalb_os ) ! cloud-sky and overcast-sky ice albedos
+      CALL ice_alb( t_su, h_i, h_s, a_ip_frac, h_ip, ln_pnd_rad, zalb_cs, zalb_os ) ! cloud-sky and overcast-sky ice albedos
       !
       alb_ice(:,:,:) = ( 1._wp - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:)
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icevar.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icevar.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icevar.F90	(revision 8563)
@@ -15,7 +15,7 @@
    !!                        - oa_i (jpi,jpj,jpl)
    !!                 VEQV : equivalent variables sometimes used in the model
-   !!                        - ht_i(jpi,jpj,jpl)
-   !!                        - ht_s(jpi,jpj,jpl)
-   !!                        - t_i (jpi,jpj,nlay_i,jpl)
+   !!                        - h_i(jpi,jpj,jpl)
+   !!                        - h_s(jpi,jpj,jpl)
+   !!                        - t_i(jpi,jpj,nlay_i,jpl)
    !!                        ...
    !!                 VAGG : aggregate variables, averaged/summed over all
@@ -174,15 +174,15 @@
       END WHERE
       !
-      ht_i(:,:,:) = v_i (:,:,:) * z1_a_i(:,:,:)    !--- ice thickness
+      h_i(:,:,:) = v_i (:,:,:) * z1_a_i(:,:,:)    !--- ice thickness
 
       zhmax    =          hi_max(jpl)
       z1_zhmax =  1._wp / hi_max(jpl)               
-      WHERE( ht_i(:,:,jpl) > zhmax )               !--- bound ht_i by hi_max (i.e. 99 m) with associated update of ice area
-         ht_i  (:,:,jpl) = zhmax
+      WHERE( h_i(:,:,jpl) > zhmax )               !--- bound h_i by hi_max (i.e. 99 m) with associated update of ice area
+         h_i  (:,:,jpl) = zhmax
          a_i   (:,:,jpl) = v_i(:,:,jpl) * z1_zhmax 
-         z1_a_i(:,:,jpl) = zhmax * z1_v_i(:,:,jpl)          ! NB: v_i always /=0 as ht_i > hi_max
+         z1_a_i(:,:,jpl) = zhmax * z1_v_i(:,:,jpl)          ! NB: v_i always /=0 as h_i > hi_max
       END WHERE
 
-      ht_s(:,:,:) = v_s (:,:,:) * z1_a_i(:,:,:)    !--- snow thickness
+      h_s(:,:,:) = v_s (:,:,:) * z1_a_i(:,:,:)    !--- snow thickness
       
       o_i(:,:,:)  = oa_i(:,:,:) * z1_a_i(:,:,:)    !--- ice age
@@ -253,6 +253,6 @@
       !!-------------------------------------------------------------------
       !
-      v_i  (:,:,:) = ht_i(:,:,:) * a_i(:,:,:)
-      v_s  (:,:,:) = ht_s(:,:,:) * a_i(:,:,:)
+      v_i  (:,:,:) = h_i(:,:,:) * a_i(:,:,:)
+      v_s  (:,:,:) = h_s(:,:,:) * a_i(:,:,:)
       smv_i(:,:,:) = sm_i(:,:,:) * v_i(:,:,:)
       !
@@ -305,5 +305,5 @@
          END DO
          !                                      ! Slope of the linear profile 
-         WHERE( ht_i(:,:,:) > epsi20 )   ;   z_slope_s(:,:,:) = 2._wp * sm_i(:,:,:) / ht_i(:,:,:)
+         WHERE( h_i(:,:,:) > epsi20 )   ;   z_slope_s(:,:,:) = 2._wp * sm_i(:,:,:) / h_i(:,:,:)
          ELSEWHERE                       ;   z_slope_s(:,:,:) = 0._wp
          END WHERE
@@ -326,5 +326,5 @@
                   DO ji = 1, jpi
                      !                          ! linear profile with 0 surface value
-                     zs0 = z_slope_s(ji,jj,jl) * ( REAL(jk,wp) - 0.5_wp ) * ht_i(ji,jj,jl) * r1_nlay_i
+                     zs0 = z_slope_s(ji,jj,jl) * ( REAL(jk,wp) - 0.5_wp ) * h_i(ji,jj,jl) * r1_nlay_i
                      zs  = zalpha(ji,jj,jl) * zs0 + ( 1._wp - zalpha(ji,jj,jl) ) * sm_i(ji,jj,jl)     ! weighting the profile
                      s_i(ji,jj,jk,jl) = MIN( rn_simax, MAX( zs, rn_simin ) )
@@ -389,5 +389,5 @@
          !
          !                                      ! Slope of the linear profile 
-         WHERE( ht_i_1d(1:nidx) > epsi20 )   ;   z_slope_s(1:nidx) = 2._wp * sm_i_1d(1:nidx) / ht_i_1d(1:nidx)
+         WHERE( h_i_1d(1:nidx) > epsi20 )   ;   z_slope_s(1:nidx) = 2._wp * sm_i_1d(1:nidx) / h_i_1d(1:nidx)
          ELSEWHERE                           ;   z_slope_s(1:nidx) = 0._wp
          END WHERE
@@ -404,5 +404,5 @@
             DO ji = 1, nidx
                !                          ! linear profile with 0 surface value
-               zs0 = z_slope_s(ji) * ( REAL(jk,wp) - 0.5_wp ) * ht_i_1d(ji) * r1_nlay_i
+               zs0 = z_slope_s(ji) * ( REAL(jk,wp) - 0.5_wp ) * h_i_1d(ji) * r1_nlay_i
                zs  = zalpha(ji) * zs0 + ( 1._wp - zalpha(ji) ) * sm_i_1d(ji)
                s_i_1d(ji,jk) = MIN( rn_simax , MAX( zs , rn_simin ) )
@@ -447,9 +447,9 @@
          ! Zap ice energy and use ocean heat to melt ice
          !-----------------------------------------------------------------
-         WHERE( a_i(:,:,jl) > epsi20 )   ;   ht_i(:,:,jl) = v_i(:,:,jl) / a_i(:,:,jl)
-         ELSEWHERE                       ;   ht_i(:,:,jl) = 0._wp
+         WHERE( a_i(:,:,jl) > epsi20 )   ;   h_i(:,:,jl) = v_i(:,:,jl) / a_i(:,:,jl)
+         ELSEWHERE                       ;   h_i(:,:,jl) = 0._wp
          END WHERE
          !
-         WHERE( a_i(:,:,jl) < epsi10 .OR. v_i(:,:,jl) < epsi10 .OR. ht_i(:,:,jl) < epsi10 )   ;   zswitch(:,:) = 0._wp
+         WHERE( a_i(:,:,jl) < epsi10 .OR. v_i(:,:,jl) < epsi10 .OR. h_i(:,:,jl) < epsi10 )   ;   zswitch(:,:) = 0._wp
          ELSEWHERE                                                                            ;   zswitch(:,:) = 1._wp
          END WHERE
@@ -490,6 +490,6 @@
                smv_i(ji,jj,jl) = smv_i(ji,jj,jl) * zswitch(ji,jj)
 
-               ht_i (ji,jj,jl) = ht_i (ji,jj,jl) * zswitch(ji,jj)
-               ht_s (ji,jj,jl) = ht_s (ji,jj,jl) * zswitch(ji,jj)
+               h_i (ji,jj,jl) = h_i (ji,jj,jl) * zswitch(ji,jj)
+               h_s (ji,jj,jl) = h_s (ji,jj,jl) * zswitch(ji,jj)
 
             END DO
@@ -519,5 +519,5 @@
 
 
-   SUBROUTINE ice_var_itd( zhti, zhts, zai, zht_i, zht_s, za_i )
+   SUBROUTINE ice_var_itd( zhti, zhts, zati, zh_i, zh_s, za_i )
       !!-------------------------------------------------------------------
       !!                ***  ROUTINE ice_var_itd   ***
@@ -543,5 +543,5 @@
       !! ** Arguments : zhti: 1-cat ice thickness
       !!                zhts: 1-cat snow depth
-      !!                zai : 1-cat ice concentration
+      !!                zati: 1-cat ice concentration
       !!
       !! ** Output    : jpl-cat 
@@ -552,6 +552,6 @@
       INTEGER  :: ijpij, i_fill, jl0  
       REAL(wp) :: zarg, zV, zconv, zdh, zdv
-      REAL(wp), DIMENSION(:),   INTENT(in)    ::   zhti, zhts, zai    ! input ice/snow variables
-      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   zht_i, zht_s, za_i ! output ice/snow variables
+      REAL(wp), DIMENSION(:),   INTENT(in)    ::   zhti, zhts, zati    ! input ice/snow variables
+      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   zh_i, zh_s, za_i ! output ice/snow variables
       INTEGER , DIMENSION(4)                  ::   itest
       !!-------------------------------------------------------------------
@@ -564,6 +564,6 @@
       ! volume and area conservation, positivity and ice categories bounds
       ijpij = SIZE( zhti , 1 )
-      zht_i(1:ijpij,1:jpl) = 0._wp
-      zht_s(1:ijpij,1:jpl) = 0._wp
+      zh_i(1:ijpij,1:jpl) = 0._wp
+      zh_s(1:ijpij,1:jpl) = 0._wp
       za_i (1:ijpij,1:jpl) = 0._wp
 
@@ -587,22 +587,22 @@
                i_fill = i_fill - 1
                !
-               zht_i(ji,1:jpl) = 0._wp
+               zh_i(ji,1:jpl) = 0._wp
                za_i (ji,1:jpl) = 0._wp
                itest(:)        = 0      
                
                IF ( i_fill == 1 ) THEN      !-- case very thin ice: fill only category 1
-                  zht_i(ji,1) = zhti(ji)
-                  za_i (ji,1) = zai (ji)
+                  zh_i(ji,1) = zhti(ji)
+                  za_i (ji,1) = zati (ji)
                ELSE                         !-- case ice is thicker: fill categories >1
                   ! thickness
                   DO jl = 1, i_fill - 1
-                     zht_i(ji,jl) = hi_mean(jl)
+                     zh_i(ji,jl) = hi_mean(jl)
                   END DO
                   
                   ! concentration
-                  za_i(ji,jl0) = zai(ji) / SQRT(REAL(jpl))
+                  za_i(ji,jl0) = zati(ji) / SQRT(REAL(jpl))
                   DO jl = 1, i_fill - 1
                      IF ( jl /= jl0 ) THEN
-                        zarg        = ( zht_i(ji,jl) - zhti(ji) ) / ( zhti(ji) * 0.5_wp )
+                        zarg        = ( zh_i(ji,jl) - zhti(ji) ) / ( zhti(ji) * 0.5_wp )
                         za_i(ji,jl) =   za_i (ji,jl0) * EXP(-zarg**2)
                      ENDIF
@@ -610,7 +610,7 @@
                   
                   ! last category
-                  za_i(ji,i_fill) = zai(ji) - SUM( za_i(ji,1:i_fill-1) )
-                  zV = SUM( za_i(ji,1:i_fill-1) * zht_i(ji,1:i_fill-1) )
-                  zht_i(ji,i_fill) = ( zhti(ji) * zai(ji) - zV ) / MAX( za_i(ji,i_fill), epsi10 ) 
+                  za_i(ji,i_fill) = zati(ji) - SUM( za_i(ji,1:i_fill-1) )
+                  zV = SUM( za_i(ji,1:i_fill-1) * zh_i(ji,1:i_fill-1) )
+                  zh_i(ji,i_fill) = ( zhti(ji) * zati(ji) - zV ) / MAX( za_i(ji,i_fill), epsi10 ) 
                   
                   ! clem: correction if concentration of upper cat is greater than lower cat
@@ -619,6 +619,6 @@
                      DO jl = jpl, jl0+1, -1
                         IF ( za_i(ji,jl) > za_i(ji,jl-1) ) THEN
-                           zdv = zht_i(ji,jl) * za_i(ji,jl)
-                           zht_i(ji,jl    ) = 0._wp
+                           zdv = zh_i(ji,jl) * za_i(ji,jl)
+                           zh_i(ji,jl    ) = 0._wp
                            za_i (ji,jl    ) = 0._wp
                            za_i (ji,1:jl-1) = za_i(ji,1:jl-1) + zdv / MAX( REAL(jl-1) * zhti(ji), epsi10 )
@@ -630,11 +630,11 @@
             
                ! Compatibility tests
-               zconv = ABS( zai(ji) - SUM( za_i(ji,1:jpl) ) ) 
+               zconv = ABS( zati(ji) - SUM( za_i(ji,1:jpl) ) ) 
                IF ( zconv < epsi06 ) itest(1) = 1                                        ! Test 1: area conservation
             
-               zconv = ABS( zhti(ji)*zai(ji) - SUM( za_i(ji,1:jpl)*zht_i(ji,1:jpl) ) )
+               zconv = ABS( zhti(ji)*zati(ji) - SUM( za_i(ji,1:jpl)*zh_i(ji,1:jpl) ) )
                IF ( zconv < epsi06 ) itest(2) = 1                                        ! Test 2: volume conservation
                
-               IF ( zht_i(ji,i_fill) >= hi_max(i_fill-1) ) itest(3) = 1                  ! Test 3: thickness of the last category is in-bounds ?
+               IF ( zh_i(ji,i_fill) >= hi_max(i_fill-1) ) itest(3) = 1                  ! Test 3: thickness of the last category is in-bounds ?
                
                itest(4) = 1
@@ -652,11 +652,11 @@
          DO ji = 1, ijpij
             IF( za_i(ji,jl) > 0._wp ) THEN
-               zht_s(ji,jl) = zht_i(ji,jl) * ( zhts(ji) / zhti(ji) )
+               zh_s(ji,jl) = zh_i(ji,jl) * ( zhts(ji) / zhti(ji) )
                ! In case snow load is in excess that would lead to transformation from snow to ice
                ! Then, transfer the snow excess into the ice (different from icethd_dh)
-               zdh = MAX( 0._wp, ( rhosn * zht_s(ji,jl) + ( rhoic - rau0 ) * zht_i(ji,jl) ) * r1_rau0 ) 
-               ! recompute ht_i, ht_s avoiding out of bounds values
-               zht_i(ji,jl) = MIN( hi_max(jl), zht_i(ji,jl) + zdh )
-               zht_s(ji,jl) = MAX( 0._wp, zht_s(ji,jl) - zdh * rhoic * r1_rhosn )
+               zdh = MAX( 0._wp, ( rhosn * zh_s(ji,jl) + ( rhoic - rau0 ) * zh_i(ji,jl) ) * r1_rau0 ) 
+               ! recompute h_i, h_s avoiding out of bounds values
+               zh_i(ji,jl) = MIN( hi_max(jl), zh_i(ji,jl) + zdh )
+               zh_s(ji,jl) = MAX( 0._wp, zh_s(ji,jl) - zdh * rhoic * r1_rhosn )
             ENDIF
          END DO
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icewri.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icewri.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icewri.F90	(revision 8563)
@@ -137,6 +137,6 @@
       !----------------------------------
       IF ( iom_use('iceconc_cat' ) )  CALL iom_put( "iceconc_cat"      , a_i   * zswi2   )        ! area for categories
-      IF ( iom_use('icethic_cat' ) )  CALL iom_put( "icethic_cat"      , ht_i  * zswi2   )        ! thickness for categories
-      IF ( iom_use('snowthic_cat') )  CALL iom_put( "snowthic_cat"     , ht_s  * zswi2   )        ! snow depth for categories
+      IF ( iom_use('icethic_cat' ) )  CALL iom_put( "icethic_cat"      , h_i  * zswi2   )        ! thickness for categories
+      IF ( iom_use('snowthic_cat') )  CALL iom_put( "snowthic_cat"     , h_s  * zswi2   )        ! snow depth for categories
       IF ( iom_use('salinity_cat') )  CALL iom_put( "salinity_cat"     , sm_i  * zswi2   )        ! salinity for categories
       IF ( iom_use('icetemp_cat' ) )  CALL iom_put( "icetemp_cat", ( SUM( t_i(:,:,:,:), dim=3 ) * r1_nlay_i - rt0 ) * zswi2 )    ! ice temperature
@@ -164,6 +164,6 @@
 
       IF ( iom_use( "iceconc_cat_mv"  ) )  CALL iom_put( "iceconc_cat_mv" , a_i(:,:,:)  * zswi2(:,:,:) + zmiss2(:,:,:) )        ! Area for categories
-      IF ( iom_use( "icethic_cat_mv"  ) )  CALL iom_put( "icethic_cat_mv" , ht_i(:,:,:) * zswi2(:,:,:) + zmiss2(:,:,:) )        ! Thickness for categories
-      IF ( iom_use( "snowthic_cat_mv" ) )  CALL iom_put( "snowthic_cat_mv", ht_s(:,:,:) * zswi2(:,:,:) + zmiss2(:,:,:) )        ! Snow depth for categories
+      IF ( iom_use( "icethic_cat_mv"  ) )  CALL iom_put( "icethic_cat_mv" , h_i(:,:,:) * zswi2(:,:,:) + zmiss2(:,:,:) )        ! Thickness for categories
+      IF ( iom_use( "snowthic_cat_mv" ) )  CALL iom_put( "snowthic_cat_mv", h_s(:,:,:) * zswi2(:,:,:) + zmiss2(:,:,:) )        ! Snow depth for categories
 
       IF ( iom_use( "icestK"   ) ) CALL iom_put( "icestK"      , tm_su(:,:)                * zswi(:,:) + zmiss(:,:) )           ! Ice surface temperature
@@ -387,9 +387,9 @@
          CALL histwrite( kid, "vfxpnd", kt, wfx_pnd     , jpi*jpj, (/1/) )
 
-      CALL histwrite( kid, "sithicat", kt, ht_i        , jpi*jpj*jpl, (/1/) )    
+      CALL histwrite( kid, "sithicat", kt, h_i         , jpi*jpj*jpl, (/1/) )    
       CALL histwrite( kid, "siconcat", kt, a_i         , jpi*jpj*jpl, (/1/) )    
       CALL histwrite( kid, "sisalcat", kt, sm_i        , jpi*jpj*jpl, (/1/) )    
       CALL histwrite( kid, "sitemcat", kt, tm_i - rt0  , jpi*jpj*jpl, (/1/) )    
-      CALL histwrite( kid, "snthicat", kt, ht_s        , jpi*jpj*jpl, (/1/) )    
+      CALL histwrite( kid, "snthicat", kt, h_s         , jpi*jpj*jpl, (/1/) )    
       CALL histwrite( kid, "sntemcat", kt, tm_su - rt0 , jpi*jpj*jpl, (/1/) )    
 
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdy_oce.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdy_oce.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdy_oce.F90	(revision 8563)
@@ -57,9 +57,9 @@
 #if defined key_lim3
       LOGICAL                           ::   ll_a_i
-      LOGICAL                           ::   ll_ht_i
-      LOGICAL                           ::   ll_ht_s
+      LOGICAL                           ::   ll_h_i
+      LOGICAL                           ::   ll_h_s
       REAL(wp), POINTER, DIMENSION(:,:) ::   a_i    !: now ice leads fraction climatology
-      REAL(wp), POINTER, DIMENSION(:,:) ::   ht_i   !: Now ice  thickness climatology
-      REAL(wp), POINTER, DIMENSION(:,:) ::   ht_s   !: now snow thickness
+      REAL(wp), POINTER, DIMENSION(:,:) ::   h_i    !: Now ice  thickness climatology
+      REAL(wp), POINTER, DIMENSION(:,:) ::   h_s    !: now snow thickness
 #endif
 #if defined key_top
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90	(revision 8563)
@@ -187,5 +187,5 @@
                   END DO
                ENDIF
-               IF( dta%ll_ht_i ) THEN
+               IF( dta%ll_h_i ) THEN
                   igrd = 1   
                   DO jl = 1, jpl
@@ -193,9 +193,9 @@
                         ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
                         ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
-                        dta_bdy(ib_bdy)%ht_i (ib,jl) =  ht_i(ii,ij,jl) * tmask(ii,ij,1) 
+                        dta_bdy(ib_bdy)%h_i (ib,jl) =  h_i(ii,ij,jl) * tmask(ii,ij,1) 
                      END DO
                   END DO
                ENDIF
-               IF( dta%ll_ht_s ) THEN
+               IF( dta%ll_h_s ) THEN
                   igrd = 1   
                   DO jl = 1, jpl
@@ -203,5 +203,5 @@
                         ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
                         ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
-                        dta_bdy(ib_bdy)%ht_s (ib,jl) =  ht_s(ii,ij,jl) * tmask(ii,ij,1) 
+                        dta_bdy(ib_bdy)%h_s (ib,jl) =  h_s(ii,ij,jl) * tmask(ii,ij,1) 
                      END DO
                   END DO
@@ -345,5 +345,5 @@
                IF( .NOT. ll_bdylim3 .AND. cn_ice_lim(ib_bdy) /= 'none' .AND. nn_ice_lim_dta(ib_bdy) == 1 ) THEN ! bdy ice input (case input is 1cat)
      	          CALL ice_var_itd ( bf(jfld_hti)%fnow(:,1,1), bf(jfld_hts)%fnow(:,1,1), bf(jfld_ai)%fnow(:,1,1), &
-                                  & dta_bdy(ib_bdy)%ht_i,     dta_bdy(ib_bdy)%ht_s,     dta_bdy(ib_bdy)%a_i     )
+                                  & dta_bdy(ib_bdy)%h_i,     dta_bdy(ib_bdy)%h_s,     dta_bdy(ib_bdy)%a_i     )
                ENDIF
 #endif
@@ -420,9 +420,9 @@
       TYPE(FLD_N) ::   bn_ssh, bn_u2d, bn_v2d           ! informations about the fields to be read
 #if defined key_lim3
-      TYPE(FLD_N) ::   bn_a_i, bn_ht_i, bn_ht_s      
+      TYPE(FLD_N) ::   bn_a_i, bn_h_i, bn_h_s      
 #endif
       NAMELIST/nambdy_dta/ cn_dir, bn_tem, bn_sal, bn_u3d, bn_v3d, bn_ssh, bn_u2d, bn_v2d 
 #if defined key_lim3
-      NAMELIST/nambdy_dta/ bn_a_i, bn_ht_i, bn_ht_s
+      NAMELIST/nambdy_dta/ bn_a_i, bn_h_i, bn_h_s
 #endif
       NAMELIST/nambdy_dta/ ln_full_vel, nb_jpk_bdy
@@ -638,7 +638,7 @@
                ENDIF
 
-               IF( dta%ll_ht_i ) THEN
-                  jfld = jfld + 1
-                  blf_i(jfld) = bn_ht_i
+               IF( dta%ll_h_i ) THEN
+                  jfld = jfld + 1
+                  blf_i(jfld) = bn_h_i
                   ibdy(jfld) = ib_bdy
                   igrid(jfld) = 1
@@ -647,7 +647,7 @@
                ENDIF
 
-               IF( dta%ll_ht_s ) THEN
-                  jfld = jfld + 1
-                   blf_i(jfld) = bn_ht_s
+               IF( dta%ll_h_s ) THEN
+                  jfld = jfld + 1
+                   blf_i(jfld) = bn_h_s
                   ibdy(jfld) = ib_bdy
                   igrid(jfld) = 1
@@ -785,15 +785,15 @@
          IF (cn_ice_lim(ib_bdy) /= 'none') THEN
             IF( nn_ice_lim_dta(ib_bdy) == 0 ) THEN
-               ALLOCATE( dta_bdy(ib_bdy)%a_i (nblen(1),jpl) )
-               ALLOCATE( dta_bdy(ib_bdy)%ht_i(nblen(1),jpl) )
-               ALLOCATE( dta_bdy(ib_bdy)%ht_s(nblen(1),jpl) )
+               ALLOCATE( dta_bdy(ib_bdy)%a_i(nblen(1),jpl) )
+               ALLOCATE( dta_bdy(ib_bdy)%h_i(nblen(1),jpl) )
+               ALLOCATE( dta_bdy(ib_bdy)%h_s(nblen(1),jpl) )
             ELSE
                IF ( ll_bdylim3 ) THEN ! case input is Xcat
                   jfld = jfld + 1
-                  dta_bdy(ib_bdy)%a_i  => bf(jfld)%fnow(:,1,:)
-                  jfld = jfld + 1
-                  dta_bdy(ib_bdy)%ht_i => bf(jfld)%fnow(:,1,:)
-                  jfld = jfld + 1
-                  dta_bdy(ib_bdy)%ht_s => bf(jfld)%fnow(:,1,:)
+                  dta_bdy(ib_bdy)%a_i => bf(jfld)%fnow(:,1,:)
+                  jfld = jfld + 1
+                  dta_bdy(ib_bdy)%h_i => bf(jfld)%fnow(:,1,:)
+                  jfld = jfld + 1
+                  dta_bdy(ib_bdy)%h_s => bf(jfld)%fnow(:,1,:)
                ELSE ! case input is 1cat
                   jfld_ai  = jfld + 1
@@ -801,10 +801,10 @@
                   jfld_hts = jfld + 3
                   jfld     = jfld + 3
-                  ALLOCATE( dta_bdy(ib_bdy)%a_i (nblen(1),jpl) )
-                  ALLOCATE( dta_bdy(ib_bdy)%ht_i(nblen(1),jpl) )
-                  ALLOCATE( dta_bdy(ib_bdy)%ht_s(nblen(1),jpl) )
-                  dta_bdy(ib_bdy)%a_i (:,:) = 0._wp
-                  dta_bdy(ib_bdy)%ht_i(:,:) = 0._wp
-                  dta_bdy(ib_bdy)%ht_s(:,:) = 0._wp
+                  ALLOCATE( dta_bdy(ib_bdy)%a_i(nblen(1),jpl) )
+                  ALLOCATE( dta_bdy(ib_bdy)%h_i(nblen(1),jpl) )
+                  ALLOCATE( dta_bdy(ib_bdy)%h_s(nblen(1),jpl) )
+                  dta_bdy(ib_bdy)%a_i(:,:) = 0._wp
+                  dta_bdy(ib_bdy)%h_i(:,:) = 0._wp
+                  dta_bdy(ib_bdy)%h_s(:,:) = 0._wp
                ENDIF
 
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdyice.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdyice.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdyice.F90	(revision 8563)
@@ -111,7 +111,7 @@
             zwgt  = idx%nbw(jb,jgrd)
             zwgt1 = 1.e0 - idx%nbw(jb,jgrd)
-            a_i (ji,jj,jl) = ( a_i (ji,jj,jl) * zwgt1 + dta%a_i (jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Leads fraction 
-            ht_i(ji,jj,jl) = ( ht_i(ji,jj,jl) * zwgt1 + dta%ht_i(jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Ice depth 
-            ht_s(ji,jj,jl) = ( ht_s(ji,jj,jl) * zwgt1 + dta%ht_s(jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Snow depth
+            a_i(ji,jj,jl) = ( a_i(ji,jj,jl) * zwgt1 + dta%a_i(jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Leads fraction 
+            h_i(ji,jj,jl) = ( h_i(ji,jj,jl) * zwgt1 + dta%h_i(jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Ice depth 
+            h_s(ji,jj,jl) = ( h_s(ji,jj,jl) * zwgt1 + dta%h_s(jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Snow depth
 
             ! -----------------
@@ -122,16 +122,16 @@
 
             ! Then, a) transfer the snow excess into the ice (different from icethd_dh)
-            zdh = MAX( 0._wp, ( rhosn * ht_s(ji,jj,jl) + ( rhoic - rau0 ) * ht_i(ji,jj,jl) ) * r1_rau0 )
+            zdh = MAX( 0._wp, ( rhosn * h_s(ji,jj,jl) + ( rhoic - rau0 ) * h_i(ji,jj,jl) ) * r1_rau0 )
             ! Or, b) transfer all the snow into ice (if incoming ice is likely to melt as it comes into a warmer environment)
-            !zdh = MAX( 0._wp, ht_s(ji,jj,jl) * rhosn / rhoic )
-
-            ! recompute ht_i, ht_s
-            ht_i(ji,jj,jl) = MIN( hi_max(jl), ht_i(ji,jj,jl) + zdh )
-            ht_s(ji,jj,jl) = MAX( 0._wp, ht_s(ji,jj,jl) - zdh * rhoic / rhosn ) 
+            !zdh = MAX( 0._wp, h_s(ji,jj,jl) * rhosn / rhoic )
+
+            ! recompute h_i, h_s
+            h_i(ji,jj,jl) = MIN( hi_max(jl), h_i(ji,jj,jl) + zdh )
+            h_s(ji,jj,jl) = MAX( 0._wp, h_s(ji,jj,jl) - zdh * rhoic / rhosn ) 
 
          ENDDO
-         CALL lbc_bdy_lnk(  a_i(:,:,jl), 'T', 1., ib_bdy )
-         CALL lbc_bdy_lnk( ht_i(:,:,jl), 'T', 1., ib_bdy )
-         CALL lbc_bdy_lnk( ht_s(:,:,jl), 'T', 1., ib_bdy )
+         CALL lbc_bdy_lnk( a_i(:,:,jl), 'T', 1., ib_bdy )
+         CALL lbc_bdy_lnk( h_i(:,:,jl), 'T', 1., ib_bdy )
+         CALL lbc_bdy_lnk( h_s(:,:,jl), 'T', 1., ib_bdy )
       ENDDO
       ! retrieve at_i
@@ -161,11 +161,11 @@
             !
             ! concentration and thickness
-            a_i (ji,jj,jl) = a_i (ii,ij,jl) * rswitch
-            ht_i(ji,jj,jl) = ht_i(ii,ij,jl) * rswitch
-            ht_s(ji,jj,jl) = ht_s(ii,ij,jl) * rswitch
+            a_i(ji,jj,jl) = a_i(ii,ij,jl) * rswitch
+            h_i(ji,jj,jl) = h_i(ii,ij,jl) * rswitch
+            h_s(ji,jj,jl) = h_s(ii,ij,jl) * rswitch
             !
             ! Ice and snow volumes
-            v_i(ji,jj,jl) = ht_i(ji,jj,jl) * a_i(ji,jj,jl)
-            v_s(ji,jj,jl) = ht_s(ji,jj,jl) * a_i(ji,jj,jl)
+            v_i(ji,jj,jl) = h_i(ji,jj,jl) * a_i(ji,jj,jl)
+            v_s(ji,jj,jl) = h_s(ji,jj,jl) * a_i(ji,jj,jl)
             !
             SELECT CASE( jpbound )
@@ -222,5 +222,5 @@
                   - rcp      * ( ztmelts - rt0 ) )
                ! Mutliply by ice volume, and divide by number of layers to get heat content in J/m2
-               e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * a_i(ji,jj,jl) * ht_i(ji,jj,jl) * r1_nlay_i
+               e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * a_i(ji,jj,jl) * h_i(ji,jj,jl) * r1_nlay_i
             END DO
             !
@@ -228,6 +228,6 @@
          !
          CALL lbc_bdy_lnk(  a_i(:,:,jl), 'T', 1., ib_bdy )
-         CALL lbc_bdy_lnk( ht_i(:,:,jl), 'T', 1., ib_bdy )
-         CALL lbc_bdy_lnk( ht_s(:,:,jl), 'T', 1., ib_bdy )
+         CALL lbc_bdy_lnk(  h_i(:,:,jl), 'T', 1., ib_bdy )
+         CALL lbc_bdy_lnk(  h_s(:,:,jl), 'T', 1., ib_bdy )
          CALL lbc_bdy_lnk(  v_i(:,:,jl), 'T', 1., ib_bdy )
          CALL lbc_bdy_lnk(  v_s(:,:,jl), 'T', 1., ib_bdy )
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdyini.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdyini.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/BDY/bdyini.F90	(revision 8563)
@@ -356,12 +356,12 @@
           CASE('none')
              IF(lwp) WRITE(numout,*) '      no open boundary condition'        
-             dta_bdy(ib_bdy)%ll_a_i  = .false.
-             dta_bdy(ib_bdy)%ll_ht_i = .false.
-             dta_bdy(ib_bdy)%ll_ht_s = .false.
+             dta_bdy(ib_bdy)%ll_a_i = .false.
+             dta_bdy(ib_bdy)%ll_h_i = .false.
+             dta_bdy(ib_bdy)%ll_h_s = .false.
           CASE('frs')
              IF(lwp) WRITE(numout,*) '      Flow Relaxation Scheme'
-             dta_bdy(ib_bdy)%ll_a_i  = .true.
-             dta_bdy(ib_bdy)%ll_ht_i = .true.
-             dta_bdy(ib_bdy)%ll_ht_s = .true.
+             dta_bdy(ib_bdy)%ll_a_i = .true.
+             dta_bdy(ib_bdy)%ll_h_i = .true.
+             dta_bdy(ib_bdy)%ll_h_s = .true.
           CASE DEFAULT   ;   CALL ctl_stop( 'unrecognised value for cn_ice_lim' )
         END SELECT
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/DIA/diadct.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/DIA/diadct.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/DIA/diadct.F90	(revision 8563)
@@ -768,11 +768,11 @@
 #if defined key_lim3
               DO jl=1,jpl
-                 transports_2d(1,jsec,jseg) = transports_2d(1,jsec,jseg) + (zTnorm)*     &
-                                   a_i(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J,jl) * &
-                                  ( ht_i(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J,jl) +  &
-                                    ht_s(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J,jl) )
+                 transports_2d(1,jsec,jseg) = transports_2d(1,jsec,jseg) + (zTnorm)*       &
+                                    a_i(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J,jl) *  &
+                                  ( h_i(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J,jl) +  &
+                                    h_s(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J,jl) )
                                    
                  transports_2d(2,jsec,jseg) = transports_2d(2,jsec,jseg) + (zTnorm)*   &
-                                   a_i(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J,jl)
+                                    a_i(sec%listPoint(jseg)%I,sec%listPoint(jseg)%J,jl)
               END DO
 #endif
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_ice.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_ice.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_ice.F90	(revision 8563)
@@ -91,5 +91,5 @@
    ! already defined in ice.F90 for LIM3
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  a_i
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  ht_i, ht_s
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::  h_i, h_s
 
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   tatm_ice       !: air temperature [K]
@@ -143,5 +143,5 @@
          &                     emp_ice(jpi,jpj)      , qns_ice(jpi,jpj,1)    , dqns_ice(jpi,jpj,1)   , &
          &                     STAT= ierr(3) )      
-      IF( ln_cpl )   ALLOCATE( ht_i(jpi,jpj,jpl) , ht_s(jpi,jpj,jpl) , STAT=ierr(4) )
+      IF( ln_cpl )   ALLOCATE( h_i(jpi,jpj,jpl) , h_s(jpi,jpj,jpl) , STAT=ierr(4) )
 #endif
 
@@ -173,5 +173,5 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   emp_ice
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   qsr_ice
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   ht_i, ht_s
+   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   h_i, h_s
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   topmelt, botmelt
    !
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90	(revision 8563)
@@ -1998,5 +1998,5 @@
             ! we must send the surface potential temperature 
             IF( l_useCT )  THEN    ;   ztmp1(:,:) = eos_pt_from_ct( tsn(:,:,1,jp_tem), tsn(:,:,1,jp_sal) )
-            ELSE                    ;   ztmp1(:,:) = tsn(:,:,1,jp_tem)
+            ELSE                   ;   ztmp1(:,:) = tsn(:,:,1,jp_tem)
             ENDIF
             !
@@ -2111,11 +2111,11 @@
             SELECT CASE( sn_snd_thick%clcat )
             CASE( 'yes' )   
-               ztmp3(:,:,1:jpl) =  ht_i(:,:,1:jpl) * a_i(:,:,1:jpl)
-               ztmp4(:,:,1:jpl) =  ht_s(:,:,1:jpl) * a_i(:,:,1:jpl)
+               ztmp3(:,:,1:jpl) =  h_i(:,:,1:jpl) * a_i(:,:,1:jpl)
+               ztmp4(:,:,1:jpl) =  h_s(:,:,1:jpl) * a_i(:,:,1:jpl)
             CASE( 'no' )
                ztmp3(:,:,:) = 0.0   ;  ztmp4(:,:,:) = 0.0
                DO jl=1,jpl
-                  ztmp3(:,:,1) = ztmp3(:,:,1) + ht_i(:,:,jl) * a_i(:,:,jl)
-                  ztmp4(:,:,1) = ztmp4(:,:,1) + ht_s(:,:,jl) * a_i(:,:,jl)
+                  ztmp3(:,:,1) = ztmp3(:,:,1) + h_i(:,:,jl) * a_i(:,:,jl)
+                  ztmp4(:,:,1) = ztmp4(:,:,1) + h_s(:,:,jl) * a_i(:,:,jl)
                ENDDO
             CASE default                  ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_thick%clcat' )
@@ -2124,10 +2124,10 @@
             SELECT CASE( sn_snd_thick%clcat )
             CASE( 'yes' )
-               ztmp3(:,:,1:jpl) = ht_i(:,:,1:jpl)
-               ztmp4(:,:,1:jpl) = ht_s(:,:,1:jpl)
+               ztmp3(:,:,1:jpl) = h_i(:,:,1:jpl)
+               ztmp4(:,:,1:jpl) = h_s(:,:,1:jpl)
             CASE( 'no' )
                WHERE( SUM( a_i, dim=3 ) /= 0. )
-                  ztmp3(:,:,1) = SUM( ht_i * a_i, dim=3 ) / SUM( a_i, dim=3 )
-                  ztmp4(:,:,1) = SUM( ht_s * a_i, dim=3 ) / SUM( a_i, dim=3 )
+                  ztmp3(:,:,1) = SUM( h_i * a_i, dim=3 ) / SUM( a_i, dim=3 )
+                  ztmp4(:,:,1) = SUM( h_s * a_i, dim=3 ) / SUM( a_i, dim=3 )
                ELSEWHERE
                  ztmp3(:,:,1) = 0.
Index: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90
===================================================================
--- branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90	(revision 8562)
+++ branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_cice.F90	(revision 8563)
@@ -721,6 +721,6 @@
 
       DO jl = 1,ncat
-         CALL cice2nemo(vsnon(:,:,jl,:),ht_s(:,:,jl),'T', 1. )
-         CALL cice2nemo(vicen(:,:,jl,:),ht_i(:,:,jl),'T', 1. )
+         CALL cice2nemo(vsnon(:,:,jl,:),h_s(:,:,jl),'T', 1. )
+         CALL cice2nemo(vicen(:,:,jl,:),h_i(:,:,jl),'T', 1. )
       ENDDO
       !
