Index: trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90
===================================================================
--- trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90	(revision 5341)
+++ trunk/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90	(revision 5407)
@@ -110,4 +110,5 @@
       REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_os, zalb_cs  ! ice albedo under overcast/clear sky
       REAL(wp), POINTER, DIMENSION(:,:,:)   ::   zalb_ice          ! mean ice albedo (for coupled)
+      REAL(wp), POINTER, DIMENSION(:,:  )   ::   zutau_ice, zvtau_ice 
       !!----------------------------------------------------------------------
 
@@ -115,4 +116,5 @@
 
       IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !  Ice time-step only
+
          !-----------------------!                                           
          ! --- Bulk Formulae --- !                                           
@@ -124,18 +126,16 @@
          t_bo(:,:) = ( eos_fzp( sss_m ) + rt0 ) * tmask(:,:,1) + rt0 * ( 1._wp - tmask(:,:,1) )  
          !                                                                                      
-         ! Ice albedo
-         CALL wrk_alloc( jpi,jpj,jpl, zalb_os, zalb_cs, zalb_ice )
-         CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os )  ! cloud-sky and overcast-sky ice albedos
-
-         ! CORE and COUPLED bulk formulations
-         SELECT CASE( kblk )
-         CASE( jp_core , jp_cpl )
-
-            ! albedo depends on cloud fraction because of non-linear spectral effects
-            zalb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:)
-            ! In CLIO the cloud fraction is read in the climatology and the all-sky albedo 
-            ! (zalb_ice) is computed within the bulk routine
-            
-         END SELECT
+!!clem         ! Ice albedo
+!!clem         CALL wrk_@lloc( jpi,jpj,jpl, zalb_os, zalb_cs, zalb_ice )
+!!clem         CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os )  ! cloud-sky and overcast-sky ice albedos
+!!
+!!         ! CORE and COUPLED bulk formulations
+!!         SELECT CASE( kblk )
+!!         CASE( jp_core , jp_purecpl )
+!!            ! albedo depends on cloud fraction because of non-linear spectral effects
+!!            zalb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:)
+!!            ! In CLIO the cloud fraction is read in the climatology and the all-sky albedo 
+!!            ! (zalb_ice) is computed within the bulk routine
+!!clem         END SELECT
          
          ! Mask sea ice surface temperature (set to rt0 over land)
@@ -154,24 +154,25 @@
          SELECT CASE( kblk )
          CASE( jp_clio )                                       ! CLIO bulk formulation
-            CALL blk_ice_clio( t_su , zalb_cs    , zalb_os    , zalb_ice  ,               &
-               &                      utau_ice   , vtau_ice   , qns_ice   , qsr_ice   ,   &
-               &                      qla_ice    , dqns_ice   , dqla_ice  ,               &
-               &                      tprecip    , sprecip    ,                           &
-               &                      fr1_i0     , fr2_i0     , cp_ice_msh, jpl  )
-            !         
-            IF( nn_limflx /= 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   &
-               &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx )
+!!clem            CALL blk_ice_clio( t_su , zalb_cs    , zalb_os    , zalb_ice  ,               &
+!!               &                      utau_ice   , vtau_ice   , qns_ice   , qsr_ice   ,   &
+!!               &                      qla_ice    , dqns_ice   , dqla_ice  ,               &
+!!               &                      tprecip    , sprecip    ,                           &
+!!               &                      fr1_i0     , fr2_i0     , cp_ice_msh, jpl  )
+!!            !         
+!!            IF( nn_limflx /= 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   &
+!!               &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx )
+            CALL blk_ice_clio_tau
 
          CASE( jp_core )                                       ! CORE bulk formulation
-            CALL blk_ice_core( t_su , u_ice     , v_ice     , zalb_ice   ,               &
-               &                      utau_ice  , vtau_ice  , qns_ice    , qsr_ice   ,   &
-               &                      qla_ice   , dqns_ice  , dqla_ice   ,               &
-               &                      tprecip   , sprecip   ,                            &
-               &                      fr1_i0    , fr2_i0    , cp_ice_msh, jpl  )
-               !
-            IF( nn_limflx /= 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   &
-               &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx )
+!!clem            CALL blk_ice_core( t_su , u_ice     , v_ice     , zalb_ice   ,               &
+!!clem               &                      utau_ice  , vtau_ice  , qns_ice    , qsr_ice   ,   &
+!!clem               &                      qla_ice   , dqns_ice  , dqla_ice   ,               &
+!!clem               &                      tprecip   , sprecip   ,                            &
+!!clem               &                      fr1_i0    , fr2_i0    , cp_ice_msh, jpl  )
+!!clem            IF( nn_limflx /= 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   &
+!!clem               &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx )
+            CALL blk_ice_core_tau
             !
-         CASE ( jp_cpl )
+         CASE ( jp_purecpl )
             
             CALL sbc_cpl_ice_tau( utau_ice , vtau_ice )
@@ -179,7 +180,16 @@
          END SELECT
          
-         !------------------------------!
-         ! --- LIM-3 main time-step --- !
-         !------------------------------!
+         IF( ln_mixcpl) THEN
+            CALL wrk_alloc( jpi,jpj    , zutau_ice, zvtau_ice)
+            CALL sbc_cpl_ice_tau( zutau_ice , zvtau_ice )
+            utau_ice(:,:) = utau_ice(:,:) * xcplmask(:,:,0) + zutau_ice(:,:) * ( 1. - xcplmask(:,:,0) )
+            vtau_ice(:,:) = vtau_ice(:,:) * xcplmask(:,:,0) + zvtau_ice(:,:) * ( 1. - xcplmask(:,:,0) )
+            CALL wrk_dealloc( jpi,jpj  , zutau_ice, zvtau_ice)
+         ENDIF
+
+         !                                           !----------------------!
+         !                                           ! LIM-3  time-stepping !
+         !                                           !----------------------!
+         ! 
          numit = numit + nn_fsbc                     ! Ice model time step
          !                                                   
@@ -220,13 +230,42 @@
          phicif(:,:)  = vt_i(:,:)
          
+         ! Ice albedo
+         CALL wrk_alloc( jpi,jpj,jpl, zalb_os, zalb_cs, zalb_ice )
+         CALL albedo_ice( t_su, ht_i, ht_s, zalb_cs, zalb_os )  ! cloud-sky and overcast-sky ice albedos
+ 
          SELECT CASE( kblk )
-         CASE ( jp_cpl )
-            CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=zalb_ice, psst=sst_m, pist=t_su    )
+         CASE( jp_clio )                                       ! CLIO bulk formulation
+            ! In CLIO the cloud fraction is read in the climatology and the all-sky albedo 
+            ! (zalb_ice) is computed within the bulk routine
+!           CALL blk_ice_clio_flx( t_su , zalb_cs, zalb_os  , zalb_ice, qns_ice   , qsr_ice   ,    &
+!              &                      qla_ice, dqns_ice   , dqla_ice  , tprecip, sprecip    ,  &
+!              &                      fr1_i0     , fr2_i0     , jpl  )
+!           !         
+            CALL blk_ice_clio_flx( t_su, zalb_cs, zalb_os, zalb_ice )
+            IF( ln_mixcpl ) CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=zalb_ice, psst=sst_m, pist=t_su )
+            IF( nn_limflx /= 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   &
+               &                                           dqns_ice, evap_ice, devap_ice, nn_limflx )
+
+         CASE( jp_core )                                       ! CORE bulk formulation
+            ! albedo depends on cloud fraction because of non-linear spectral effects
+            zalb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:)
+            CALL blk_ice_core_flx( t_su, zalb_ice )
+            IF( ln_mixcpl ) CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=zalb_ice, psst=sst_m, pist=t_su )
+            IF( nn_limflx /= 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   &
+               &                                           dqns_ice, evap_ice, devap_ice, nn_limflx )
+
+         CASE ( jp_purecpl )
+            ! albedo depends on cloud fraction because of non-linear spectral effects
+            zalb_ice(:,:,:) = ( 1. - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:)
+            CALL sbc_cpl_ice_flx( p_frld=pfrld, palbi=zalb_ice, psst=sst_m, pist=t_su )
             IF( nn_limflx == 2 )   CALL ice_lim_flx( t_su, zalb_ice, qns_ice, qsr_ice ,   &
-               &                                           dqns_ice, qla_ice, dqla_ice, nn_limflx )
+               &                                           dqns_ice, evap_ice, devap_ice, nn_limflx )
             ! Latent heat flux is forced to 0 in coupled: it is included in qns (non-solar heat flux)
-            qla_ice  (:,:,:) = 0._wp
-            dqla_ice (:,:,:) = 0._wp
+            evap_ice  (:,:,:) = 0._wp
+            devap_ice (:,:,:) = 0._wp
+
          END SELECT
+         CALL wrk_dealloc( jpi,jpj,jpl, zalb_os, zalb_cs, zalb_ice )
+
          !
          CALL lim_thd( kt )                         ! Ice thermodynamics 
@@ -247,5 +286,4 @@
          IF( ln_icectl )  CALL lim_ctl( kt )        ! alerts in case of model crash
          !
-         CALL wrk_dealloc( jpi,jpj,jpl, zalb_os, zalb_cs, zalb_ice )
          !
       ENDIF   ! End sea-ice time step only
@@ -476,5 +514,5 @@
    
    SUBROUTINE ice_lim_flx( ptn_ice, palb_ice, pqns_ice, pqsr_ice,   &
-         &                          pdqn_ice, pqla_ice, pdql_ice, k_limflx )
+         &                          pdqn_ice, pevap_ice, pdevap_ice, k_limflx )
       !!---------------------------------------------------------------------
       !!                  ***  ROUTINE ice_lim_flx  ***
@@ -494,6 +532,6 @@
       REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pqsr_ice   ! net solar flux
       REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pdqn_ice   ! non solar flux sensitivity
-      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pqla_ice   ! latent heat flux
-      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pdql_ice   ! latent heat flux sensitivity
+      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pevap_ice  ! sublimation
+      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   pdevap_ice ! sublimation sensitivity
       !
       INTEGER  ::   jl      ! dummy loop index
@@ -504,7 +542,7 @@
       REAL(wp), POINTER, DIMENSION(:,:) :: z_qsr_m   ! Mean solar heat flux over all categories
       REAL(wp), POINTER, DIMENSION(:,:) :: z_qns_m   ! Mean non solar heat flux over all categories
-      REAL(wp), POINTER, DIMENSION(:,:) :: z_qla_m   ! Mean latent heat flux over all categories
+      REAL(wp), POINTER, DIMENSION(:,:) :: z_evap_m  ! Mean sublimation over all categories
       REAL(wp), POINTER, DIMENSION(:,:) :: z_dqn_m   ! Mean d(qns)/dT over all categories
-      REAL(wp), POINTER, DIMENSION(:,:) :: z_dql_m   ! Mean d(qla)/dT over all categories
+      REAL(wp), POINTER, DIMENSION(:,:) :: z_devap_m ! Mean d(evap)/dT over all categories
       !!----------------------------------------------------------------------
 
@@ -514,14 +552,14 @@
       SELECT CASE( k_limflx )                              !==  averaged on all ice categories  ==!
       CASE( 0 , 1 )
-         CALL wrk_alloc( jpi,jpj, z_qsr_m, z_qns_m, z_qla_m, z_dqn_m, z_dql_m)
+         CALL wrk_alloc( jpi,jpj, z_qsr_m, z_qns_m, z_evap_m, z_dqn_m, z_devap_m)
          !
          z_qns_m(:,:) = fice_ice_ave ( pqns_ice (:,:,:) )
          z_qsr_m(:,:) = fice_ice_ave ( pqsr_ice (:,:,:) )
          z_dqn_m(:,:) = fice_ice_ave ( pdqn_ice (:,:,:) )
-         z_qla_m(:,:) = fice_ice_ave ( pqla_ice (:,:,:) )
-         z_dql_m(:,:) = fice_ice_ave ( pdql_ice (:,:,:) )
+         z_evap_m(:,:) = fice_ice_ave ( pevap_ice (:,:,:) )
+         z_devap_m(:,:) = fice_ice_ave ( pdevap_ice (:,:,:) )
          DO jl = 1, jpl
             pdqn_ice(:,:,jl) = z_dqn_m(:,:)
-            pdql_ice(:,:,jl) = z_dql_m(:,:)
+            pdevap_ice(:,:,jl) = z_devap_m(:,:)
          END DO
          !
@@ -529,8 +567,8 @@
             pqns_ice(:,:,jl) = z_qns_m(:,:)
             pqsr_ice(:,:,jl) = z_qsr_m(:,:)
-            pqla_ice(:,:,jl) = z_qla_m(:,:)
+            pevap_ice(:,:,jl) = z_evap_m(:,:)
          END DO
          !
-         CALL wrk_dealloc( jpi,jpj, z_qsr_m, z_qns_m, z_qla_m, z_dqn_m, z_dql_m)
+         CALL wrk_dealloc( jpi,jpj, z_qsr_m, z_qns_m, z_evap_m, z_dqn_m, z_devap_m)
       END SELECT
 
@@ -543,5 +581,5 @@
          DO jl = 1, jpl
             pqns_ice(:,:,jl) = pqns_ice(:,:,jl) + pdqn_ice(:,:,jl) * (ptn_ice(:,:,jl) - ztem_m(:,:))
-            pqla_ice(:,:,jl) = pqla_ice(:,:,jl) + pdql_ice(:,:,jl) * (ptn_ice(:,:,jl) - ztem_m(:,:))
+            pevap_ice(:,:,jl) = pevap_ice(:,:,jl) + pdevap_ice(:,:,jl) * (ptn_ice(:,:,jl) - ztem_m(:,:))
             pqsr_ice(:,:,jl) = pqsr_ice(:,:,jl) * ( 1._wp - palb_ice(:,:,jl) ) / ( 1._wp - zalb_m(:,:) ) 
          END DO
@@ -593,5 +631,4 @@
       wfx_spr(:,:) = 0._wp   ;   
       
-      hfx_in (:,:) = 0._wp   ;   hfx_out(:,:) = 0._wp
       hfx_thd(:,:) = 0._wp   ;   
       hfx_snw(:,:) = 0._wp   ;   hfx_opw(:,:) = 0._wp
@@ -610,5 +647,6 @@
       
    END SUBROUTINE sbc_lim_diag0
-      
+
+     
    FUNCTION fice_cell_ave ( ptab )
       !!--------------------------------------------------------------------------
