diff --git a/tests/BENCH/MY_SRC/usrdef_sbc.F90 b/tests/BENCH/MY_SRC/usrdef_sbc.F90
index fb690fe8e436f550a69bf82751ca49db6a1ef512..899c20bb3d6cfbdc8b273b08a7efb771a9fd4fc2 100644
--- a/tests/BENCH/MY_SRC/usrdef_sbc.F90
+++ b/tests/BENCH/MY_SRC/usrdef_sbc.F90
@@ -17,14 +17,14 @@ MODULE usrdef_sbc
    USE oce             ! ocean dynamics and tracers
    USE sbc_oce         ! Surface boundary condition: ocean fields
    USE sbc_ice         ! Surface boundary condition: ocean fields
-   USE in_out_manager  ! I/O manager
+   USE sbc_phy, ONLY   : pp_cldf
    USE phycst          ! physical constants
-   USE lib_mpp         ! MPP library
-   USE lbclnk          ! lateral boundary conditions - mpp exchanges
-
 #if defined key_si3
-   USE ice, ONLY       : at_i_b, a_i_b
+   USE ice, ONLY       : jpl, at_i_b, a_i_b
 #endif
+   USE in_out_manager  ! I/O manager
+   USE lib_mpp         ! MPP library
+   USE lbclnk          ! lateral boundary conditions - mpp exchanges
 
    IMPLICIT NONE
    PRIVATE
@@ -132,7 +132,10 @@ CONTAINS
       REAL(wp), DIMENSION(:,:,:), INTENT(in)  ::   phs    ! snow thickness
       REAL(wp), DIMENSION(:,:,:), INTENT(in)  ::   phi    ! ice thickness
       !!
+      INTEGER  ::   jl
+      REAL(wp) ::   zfr1, zfr2                ! local variables
       REAL(wp), DIMENSION(A2D(0)) ::   zsnw   ! snw distribution after wind blowing
+      REAL(wp), DIMENSION(A2D(0)) ::   ztri
       !!---------------------------------------------------------------------
 #if defined key_si3
       !
@@ -166,8 +169,19 @@ CONTAINS
       qns_tot (:,:) = at_i_b(:,:) * qns_oce(:,:) + SUM( a_i_b(:,:,:) * qns_ice(:,:,:), dim=3 ) + qemp_ice(:,:) + qemp_oce(:,:)
       qsr_tot (:,:) = at_i_b(:,:) * qsr_oce(:,:) + SUM( a_i_b(:,:,:) * qsr_ice(:,:,:), dim=3 )
 
-      ! --- shortwave radiation transmitted thru the surface scattering layer (W/m2) --- !
-      qtr_ice_top(:,:,:) = 0._wp
+      ! --- shortwave radiation transmitted below the surface (W/m2, see Grenfell Maykut 77) --- !
+      cloud_fra(:,:) = pp_cldf
+      ztri(:,:) = 0.18 * ( 1.0 - cloud_fra(:,:) ) + 0.35 * cloud_fra(:,:)  ! surface transmission when hi>10cm
+      !
+      DO jl = 1, jpl
+         WHERE    ( phs(A2D(0),jl) <= 0._wp .AND. phi(A2D(0),jl) <  0.1_wp )     ! linear decrease from hi=0 to 10cm  
+            qtr_ice_top(:,:,jl) = qsr_ice(:,:,jl) * ( ztri(:,:) + ( 1._wp - ztri(:,:) ) * ( 1._wp - phi(A2D(0),jl) * 10._wp ) )
+         ELSEWHERE( phs(A2D(0),jl) <= 0._wp .AND. phi(A2D(0),jl) >= 0.1_wp )     ! constant (ztri) when hi>10cm
+            qtr_ice_top(:,:,jl) = qsr_ice(:,:,jl) * ztri(:,:)
+         ELSEWHERE                                                         ! zero when hs>0
+            qtr_ice_top(:,:,jl) = 0._wp
+         END WHERE
+      ENDDO
 #endif
 
    END SUBROUTINE usrdef_sbc_ice_flx
diff --git a/tests/DONUT/MY_SRC/usrdef_sbc.F90 b/tests/DONUT/MY_SRC/usrdef_sbc.F90
index aee05bbde386bdb19b4ab20addfeca1f187010c1..07fa6c770e5f420b32fa1ffc0e61b11e351e02e3 100644
--- a/tests/DONUT/MY_SRC/usrdef_sbc.F90
+++ b/tests/DONUT/MY_SRC/usrdef_sbc.F90
@@ -16,11 +16,12 @@ MODULE usrdef_sbc
    USE dom_oce         ! ocean space and time domain
    USE sbc_oce         ! Surface boundary condition: ocean fields
    USE sbc_ice         ! Surface boundary condition: ice fields
-   USE in_out_manager  ! I/O manager
+   USE sbc_phy, ONLY   : pp_cldf
    USE phycst          ! physical constants
 #if defined key_si3
-   USE ice, ONLY       : at_i_b, a_i_b
+   USE ice, ONLY       : jpl, at_i_b, a_i_b
 #endif
+   USE in_out_manager  ! I/O manager
 
    IMPLICIT NONE
    PRIVATE
@@ -29,6 +30,8 @@ MODULE usrdef_sbc
    PUBLIC   usrdef_sbc_ice_tau  ! routine called by icestp.F90 for ice dynamics
    PUBLIC   usrdef_sbc_ice_flx  ! routine called by icestp.F90 for ice thermo
 
+   !! * Substitutions
+#  include "do_loop_substitute.h90"
    !!----------------------------------------------------------------------
    !! NEMO/OCE 4.0 , NEMO Consortium (2018)
    !! $Id: usrdef_sbc.F90 10074 2018-08-28 16:15:49Z nicolasmartin $
@@ -107,7 +110,10 @@ CONTAINS
       REAL(wp), DIMENSION(:,:,:), INTENT(in) ::   phs    ! snow thickness
       REAL(wp), DIMENSION(:,:,:), INTENT(in) ::   phi    ! ice thickness
       !!
-      REAL(wp), DIMENSION(jpi,jpj) ::   zsnw   ! snw distribution after wind blowing
+      INTEGER  ::   jl
+      REAL(wp) ::   zfr1, zfr2                ! local variables
+      REAL(wp), DIMENSION(A2D(0)) ::   zsnw   ! snw distribution after wind blowing
+      REAL(wp), DIMENSION(A2D(0)) ::   ztri
       !!---------------------------------------------------------------------
 #if defined key_si3
       !
@@ -132,17 +138,28 @@ CONTAINS
       emp_ice  (:,:)   = SUM( a_i_b(:,:,:) * evap_ice(:,:,:), dim=3 ) - sprecip(:,:) * zsnw(:,:)
       emp_oce  (:,:)   = emp_oce(:,:) - sprecip(:,:) * (1._wp - zsnw(:,:) )
       qevap_ice(:,:,:) =   0._wp
-      qprec_ice(:,:)   =   rhos * ( sst_m(:,:) * rcpi - rLfus ) * tmask(:,:,1) !  in J/m3
-      qemp_oce (:,:)   = - emp_oce(:,:) * sst_m(:,:) * rcp
-      qemp_ice (:,:)   =   sprecip(:,:) * zsnw * ( sst_m(:,:) * rcpi - rLfus ) * tmask(:,:,1) ! solid precip (only)
+      qprec_ice(:,:)   =   rhos * ( sst_m(A2D(0)) * rcpi - rLfus ) * smask0(:,:) !  in J/m3
+      qemp_oce (:,:)   = - emp_oce(:,:) * sst_m(A2D(0)) * rcp
+      qemp_ice (:,:)   =   sprecip(:,:) * zsnw * ( sst_m(A2D(0)) * rcpi - rLfus ) * smask0(:,:) ! solid precip (only)
 
       ! total fluxes
       emp_tot (:,:) = emp_ice  + emp_oce
       qns_tot (:,:) = at_i_b(:,:) * qns_oce(:,:) + SUM( a_i_b(:,:,:) * qns_ice(:,:,:), dim=3 ) + qemp_ice(:,:) + qemp_oce(:,:)
       qsr_tot (:,:) = at_i_b(:,:) * qsr_oce(:,:) + SUM( a_i_b(:,:,:) * qsr_ice(:,:,:), dim=3 )
 
-      ! --- shortwave radiation transmitted thru the surface scattering layer (W/m2) --- !
-      qtr_ice_top(:,:,:) = 0._wp
+      ! --- shortwave radiation transmitted below the surface (W/m2, see Grenfell Maykut 77) --- !
+      cloud_fra(:,:) = pp_cldf
+      ztri(:,:) = 0.18 * ( 1.0 - cloud_fra(:,:) ) + 0.35 * cloud_fra(:,:)  ! surface transmission when hi>10cm
+      !
+      DO jl = 1, jpl
+         WHERE    ( phs(A2D(0),jl) <= 0._wp .AND. phi(A2D(0),jl) <  0.1_wp )     ! linear decrease from hi=0 to 10cm  
+            qtr_ice_top(:,:,jl) = qsr_ice(:,:,jl) * ( ztri(:,:) + ( 1._wp - ztri(:,:) ) * ( 1._wp - phi(A2D(0),jl) * 10._wp ) )
+         ELSEWHERE( phs(A2D(0),jl) <= 0._wp .AND. phi(A2D(0),jl) >= 0.1_wp )     ! constant (ztri) when hi>10cm
+            qtr_ice_top(:,:,jl) = qsr_ice(:,:,jl) * ztri(:,:)
+         ELSEWHERE                                                         ! zero when hs>0
+            qtr_ice_top(:,:,jl) = 0._wp
+         END WHERE
+      ENDDO
 #endif
 
    END SUBROUTINE usrdef_sbc_ice_flx
diff --git a/tests/DONUT/cpp_DONUT.fcm b/tests/DONUT/cpp_DONUT.fcm
index 647c4d5f591339110d6c9bf1de8efb89a605f8cd..676164c5a7a20f8628911ec7bbb0694ac3aa0321 100644
--- a/tests/DONUT/cpp_DONUT.fcm
+++ b/tests/DONUT/cpp_DONUT.fcm
@@ -1 +1 @@
- bld::tool::fppkeys key_si3 key_qco key_xios key_zco
+ bld::tool::fppkeys key_si3 key_qco key_xios key_vco_1d
diff --git a/tests/ICE_ADV1D/MY_SRC/usrdef_zgr.F90 b/tests/ICE_ADV1D/MY_SRC/usrdef_zgr.F90
index e9e25e4a25eee8e1fd25706ca646e6ae40d41ec0..559fdecfae5d5e5b4b74e940fab0cfbd533ff106 100644
--- a/tests/ICE_ADV1D/MY_SRC/usrdef_zgr.F90
+++ b/tests/ICE_ADV1D/MY_SRC/usrdef_zgr.F90
@@ -32,11 +32,11 @@ MODULE usrdef_zgr
 CONTAINS             
 
    SUBROUTINE usr_def_zgr( ld_zco  , ld_zps  , ld_sco  , ld_isfcav,    &   ! type of vertical coordinate
+      &                    k_top   , k_bot                        ,    &   ! top & bottom ocean level
       &                    pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d  ,    &   ! 1D reference vertical coordinate
+      &                    pe3t  , pe3u  , pe3v   , pe3f ,             &   ! vertical scale factors
       &                    pdept , pdepw ,                             &   ! 3D t & w-points depth
-      &                    pe3t  , pe3u  , pe3v , pe3f ,               &   ! vertical scale factors
-      &                    pe3w  , pe3uw , pe3vw,                      &   !     -      -      -
-      &                    k_top  , k_bot    )                             ! top & bottom ocean level
+      &                    pe3w  , pe3uw , pe3vw                       )   ! vertical scale factors
       !!---------------------------------------------------------------------
       !!              ***  ROUTINE usr_def_zgr  ***
       !!
@@ -45,12 +45,12 @@ CONTAINS
       !!----------------------------------------------------------------------
       LOGICAL                   , INTENT(out) ::   ld_zco, ld_zps, ld_sco      ! vertical coordinate flags
       LOGICAL                   , INTENT(out) ::   ld_isfcav                   ! under iceshelf cavity flag
+      INTEGER , DIMENSION(:,:)  , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
       REAL(wp), DIMENSION(:)    , INTENT(out) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth     [m]
       REAL(wp), DIMENSION(:)    , INTENT(out) ::   pe3t_1d , pe3w_1d           ! 1D grid-point depth     [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pdept, pdepw                ! grid-point depth        [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors  [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pe3w , pe3uw, pe3vw         ! i-scale factors 
-      INTEGER , DIMENSION(:,:)  , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
+      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pdept, pdepw                ! grid-point depth        [m]
+      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors  [m]
+      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pe3w , pe3uw, pe3vw         ! i-scale factors 
       !
       INTEGER  ::   jk, k_dz  ! dummy indices
       !!----------------------------------------------------------------------
diff --git a/tests/ICE_ADV1D/cpp_ICE_ADV1D.fcm b/tests/ICE_ADV1D/cpp_ICE_ADV1D.fcm
index 21a0f66fd4bc8b4d87e9d98cbb3c188a68301f66..90a47cd951089c9bb1a48ea9d0600b1ccaec872d 100644
--- a/tests/ICE_ADV1D/cpp_ICE_ADV1D.fcm
+++ b/tests/ICE_ADV1D/cpp_ICE_ADV1D.fcm
@@ -1 +1 @@
-bld::tool::fppkeys key_si3 key_xios key_linssh
+bld::tool::fppkeys key_si3 key_xios key_linssh key_vco_1d
diff --git a/tests/ICE_ADV2D/MY_SRC/usrdef_sbc.F90 b/tests/ICE_ADV2D/MY_SRC/usrdef_sbc.F90
index 508e0eb4b852730f042845d3e9eaf3819fcc3563..f0804e33c32a3fced05a771ae7e0ae19f16059b1 100644
--- a/tests/ICE_ADV2D/MY_SRC/usrdef_sbc.F90
+++ b/tests/ICE_ADV2D/MY_SRC/usrdef_sbc.F90
@@ -21,9 +21,6 @@ MODULE usrdef_sbc
    USE ice, ONLY       : jpl, at_i_b, a_i_b
    !
    USE in_out_manager  ! I/O manager
-   USE lib_mpp         ! distribued memory computing library
-   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
-   USE lib_fortran     ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 
 
    IMPLICIT NONE
    PRIVATE
diff --git a/tests/ICE_ADV2D/MY_SRC/usrdef_zgr.F90 b/tests/ICE_ADV2D/MY_SRC/usrdef_zgr.F90
index a4ce102ab18fff2faf845ece0bc214d18b40effc..5f64ec94bb3d468901932a353753e9a0b7aeeb78 100644
--- a/tests/ICE_ADV2D/MY_SRC/usrdef_zgr.F90
+++ b/tests/ICE_ADV2D/MY_SRC/usrdef_zgr.F90
@@ -32,11 +32,11 @@ MODULE usrdef_zgr
 CONTAINS             
 
    SUBROUTINE usr_def_zgr( ld_zco  , ld_zps  , ld_sco  , ld_isfcav,    &   ! type of vertical coordinate
+      &                    k_top   , k_bot                        ,    &   ! top & bottom ocean level
       &                    pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d  ,    &   ! 1D reference vertical coordinate
+      &                    pe3t  , pe3u  , pe3v   , pe3f ,             &   ! vertical scale factors
       &                    pdept , pdepw ,                             &   ! 3D t & w-points depth
-      &                    pe3t  , pe3u  , pe3v , pe3f ,               &   ! vertical scale factors
-      &                    pe3w  , pe3uw , pe3vw,                      &   !     -      -      -
-      &                    k_top  , k_bot    )                             ! top & bottom ocean level
+      &                    pe3w  , pe3uw , pe3vw                       )   ! vertical scale factors
       !!---------------------------------------------------------------------
       !!              ***  ROUTINE usr_def_zgr  ***
       !!
@@ -45,12 +45,12 @@ CONTAINS
       !!----------------------------------------------------------------------
       LOGICAL                   , INTENT(out) ::   ld_zco, ld_zps, ld_sco      ! vertical coordinate flags
       LOGICAL                   , INTENT(out) ::   ld_isfcav                   ! under iceshelf cavity flag
+      INTEGER , DIMENSION(:,:)  , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
       REAL(wp), DIMENSION(:)    , INTENT(out) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth     [m]
       REAL(wp), DIMENSION(:)    , INTENT(out) ::   pe3t_1d , pe3w_1d           ! 1D grid-point depth     [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pdept, pdepw                ! grid-point depth        [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors  [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pe3w , pe3uw, pe3vw         ! i-scale factors 
-      INTEGER , DIMENSION(:,:)  , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
+      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pdept, pdepw                ! grid-point depth        [m]
+      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors  [m]
+      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pe3w , pe3uw, pe3vw         ! i-scale factors 
       !
       INTEGER  ::   jk, k_dz  ! dummy indices
       !!----------------------------------------------------------------------
diff --git a/tests/ICE_ADV2D/cpp_ICE_ADV2D.fcm b/tests/ICE_ADV2D/cpp_ICE_ADV2D.fcm
index 47e40c89b444e4a2267e0fd795bf074100de7148..36ee7c29d7228c625d9669e5a05eca4e9ef95612 100644
--- a/tests/ICE_ADV2D/cpp_ICE_ADV2D.fcm
+++ b/tests/ICE_ADV2D/cpp_ICE_ADV2D.fcm
@@ -1 +1 @@
-bld::tool::fppkeys key_si3 key_linssh key_xios
+bld::tool::fppkeys key_si3 key_linssh key_xios key_vco_1d
diff --git a/tests/ICE_RHEO/MY_SRC/usrdef_zgr.F90 b/tests/ICE_RHEO/MY_SRC/usrdef_zgr.F90
index 096e524a2b51bc2f9bff953f82c462c8722ae995..e62e295a42e703f271d04900f2f3ca4495abd6b5 100644
--- a/tests/ICE_RHEO/MY_SRC/usrdef_zgr.F90
+++ b/tests/ICE_RHEO/MY_SRC/usrdef_zgr.F90
@@ -32,11 +32,11 @@ MODULE usrdef_zgr
 CONTAINS             
 
    SUBROUTINE usr_def_zgr( ld_zco  , ld_zps  , ld_sco  , ld_isfcav,    &   ! type of vertical coordinate
+      &                    k_top   , k_bot                        ,    &   ! top & bottom ocean level
       &                    pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d  ,    &   ! 1D reference vertical coordinate
+      &                    pe3t  , pe3u  , pe3v   , pe3f ,             &   ! vertical scale factors
       &                    pdept , pdepw ,                             &   ! 3D t & w-points depth
-      &                    pe3t  , pe3u  , pe3v , pe3f ,               &   ! vertical scale factors
-      &                    pe3w  , pe3uw , pe3vw,                      &   !     -      -      -
-      &                    k_top  , k_bot    )                             ! top & bottom ocean level
+      &                    pe3w  , pe3uw , pe3vw                       )   ! vertical scale factors
       !!---------------------------------------------------------------------
       !!              ***  ROUTINE usr_def_zgr  ***
       !!
@@ -45,12 +45,12 @@ CONTAINS
       !!----------------------------------------------------------------------
       LOGICAL                   , INTENT(out) ::   ld_zco, ld_zps, ld_sco      ! vertical coordinate flags
       LOGICAL                   , INTENT(out) ::   ld_isfcav                   ! under iceshelf cavity flag
+      INTEGER , DIMENSION(:,:)  , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
       REAL(wp), DIMENSION(:)    , INTENT(out) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth     [m]
       REAL(wp), DIMENSION(:)    , INTENT(out) ::   pe3t_1d , pe3w_1d           ! 1D grid-point depth     [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pdept, pdepw                ! grid-point depth        [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors  [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pe3w , pe3uw, pe3vw         ! i-scale factors 
-      INTEGER , DIMENSION(:,:)  , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
+      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pdept, pdepw                ! grid-point depth        [m]
+      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors  [m]
+      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pe3w , pe3uw, pe3vw         ! i-scale factors 
       !
       INTEGER  ::   jk, k_dz  ! dummy indices
       !!----------------------------------------------------------------------
@@ -88,7 +88,7 @@ CONTAINS
       !
       !                       !==  z-coordinate  ==!   (step-like topography)
       !                                !* bottom ocean compute from the depth of grid-points
-      jpkm1 = jpk
+      !!jpkm1 = jpk
       k_bot(:,:) = 1    ! here use k_top as a land mask
       !                                !* horizontally uniform coordinate (reference z-co everywhere)
       DO jk = 1, jpk
diff --git a/tests/ICE_RHEO/cpp_ICE_RHEO.fcm b/tests/ICE_RHEO/cpp_ICE_RHEO.fcm
index 47e40c89b444e4a2267e0fd795bf074100de7148..36ee7c29d7228c625d9669e5a05eca4e9ef95612 100644
--- a/tests/ICE_RHEO/cpp_ICE_RHEO.fcm
+++ b/tests/ICE_RHEO/cpp_ICE_RHEO.fcm
@@ -1 +1 @@
-bld::tool::fppkeys key_si3 key_linssh key_xios
+bld::tool::fppkeys key_si3 key_linssh key_xios key_vco_1d