Index: branches/2016/dev_v3_6_STABLE_r6506_AGRIF_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcana.F90
===================================================================
--- branches/2016/dev_v3_6_STABLE_r6506_AGRIF_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcana.F90	(revision 7069)
+++ branches/2016/dev_v3_6_STABLE_r6506_AGRIF_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcana.F90	(revision 7077)
@@ -15,4 +15,5 @@
    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 phycst          ! physical constants
    USE in_out_manager  ! I/O manager
@@ -20,19 +21,36 @@
    USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
    USE lib_fortran
-
+   USE wrk_nemo
+#if defined key_lim3
+   USE ice, ONLY       : pfrld, a_i_b
+   USE limthd_dh       ! for CALL lim_thd_snwblow
+#endif
+   
    IMPLICIT NONE
    PRIVATE
 
-   PUBLIC   sbc_ana    ! routine called in sbcmod module
-   PUBLIC   sbc_gyre   ! routine called in sbcmod module
+   PUBLIC   sbc_ana         ! routine called in sbcmod module
+   PUBLIC   sbc_gyre        ! routine called in sbcmod module
+#if defined key_lim3
+   PUBLIC   ana_ice_tau     ! routine called in sbc_ice_lim module
+   PUBLIC   ana_ice_flx     ! routine called in sbc_ice_lim module
+#endif
 
    !                       !!* Namelist namsbc_ana *
-   INTEGER  ::   nn_tau000  ! nb of time-step during which the surface stress
-   !                        ! increase from 0 to its nominal value 
-   REAL(wp) ::   rn_utau0   ! constant wind stress value in i-direction
-   REAL(wp) ::   rn_vtau0   ! constant wind stress value in j-direction
-   REAL(wp) ::   rn_qns0    ! non solar heat flux
-   REAL(wp) ::   rn_qsr0    !     solar heat flux
-   REAL(wp) ::   rn_emp0    ! net freshwater flux
+   ! --- oce variables --- !
+   INTEGER  ::   nn_tau000 ! nb of time-step during which the surface stress
+   !                       ! increase from 0 to its nominal value 
+   REAL(wp) ::   rn_utau0  ! constant wind stress value in i-direction
+   REAL(wp) ::   rn_vtau0  ! constant wind stress value in j-direction
+   REAL(wp) ::   rn_qns0   ! non solar heat flux
+   REAL(wp) ::   rn_qsr0   !     solar heat flux
+   REAL(wp) ::   rn_emp0   ! net freshwater flux
+   ! --- ice variables --- !
+   REAL(wp) ::   rn_iutau0 ! constant wind stress value in i-direction over ice
+   REAL(wp) ::   rn_ivtau0 ! constant wind stress value in j-direction over ice
+   REAL(wp) ::   rn_iqns0  ! non solar heat flux over ice
+   REAL(wp) ::   rn_iqsr0  !     solar heat flux over ice
+   REAL(wp) ::   rn_sprec0 ! snow precip
+   REAL(wp) ::   rn_ievap0 ! sublimation
    
    !! * Substitutions
@@ -69,5 +87,6 @@
       REAL(wp) ::   zcoef, zty, zmod      !   -      -
       !!
-      NAMELIST/namsbc_ana/ nn_tau000, rn_utau0, rn_vtau0, rn_qns0, rn_qsr0, rn_emp0
+      NAMELIST/namsbc_ana/ nn_tau000, rn_utau0, rn_vtau0, rn_qns0, rn_qsr0, rn_emp0,  &
+         &                 rn_iutau0, rn_ivtau0, rn_iqsr0, rn_iqns0, rn_sprec0, rn_ievap0
       !!---------------------------------------------------------------------
       !
@@ -86,10 +105,16 @@
          IF(lwp) WRITE(numout,*)' sbc_ana : Constant surface fluxes read in namsbc_ana namelist'
          IF(lwp) WRITE(numout,*)' ~~~~~~~ '
-         IF(lwp) WRITE(numout,*)'              spin up of the stress  nn_tau000 = ', nn_tau000, ' time-steps'
-         IF(lwp) WRITE(numout,*)'              constant i-stress      rn_utau0  = ', rn_utau0 , ' N/m2'
-         IF(lwp) WRITE(numout,*)'              constant j-stress      rn_vtau0  = ', rn_vtau0 , ' N/m2'
-         IF(lwp) WRITE(numout,*)'              non solar heat flux    rn_qns0   = ', rn_qns0  , ' W/m2'
-         IF(lwp) WRITE(numout,*)'              solar heat flux        rn_qsr0   = ', rn_qsr0  , ' W/m2'
-         IF(lwp) WRITE(numout,*)'              net heat flux          rn_emp0   = ', rn_emp0  , ' Kg/m2/s'
+         IF(lwp) WRITE(numout,*)'              spin up of the stress         nn_tau000 = ', nn_tau000 , ' time-steps'
+         IF(lwp) WRITE(numout,*)'              constant i-stress             rn_utau0  = ', rn_utau0  , ' N/m2'
+         IF(lwp) WRITE(numout,*)'              constant j-stress             rn_vtau0  = ', rn_vtau0  , ' N/m2'
+         IF(lwp) WRITE(numout,*)'              non solar heat flux           rn_qns0   = ', rn_qns0   , ' W/m2'
+         IF(lwp) WRITE(numout,*)'              solar heat flux               rn_qsr0   = ', rn_qsr0   , ' W/m2'
+         IF(lwp) WRITE(numout,*)'              net freshwater flux           rn_emp0   = ', rn_emp0   , ' Kg/m2/s'
+         IF(lwp) WRITE(numout,*)'              constant ice-atm stress       rn_iutau0 = ', rn_iutau0 , ' N/m2'
+         IF(lwp) WRITE(numout,*)'              constant ice-atm stress       rn_ivtau0 = ', rn_ivtau0 , ' N/m2'
+         IF(lwp) WRITE(numout,*)'              solar heat flux over ice      rn_iqsr0  = ', rn_iqsr0  , ' W/m2'
+         IF(lwp) WRITE(numout,*)'              non solar heat flux over ice  rn_iqns0  = ', rn_iqns0  , ' W/m2'
+         IF(lwp) WRITE(numout,*)'              snow precip                   rn_sprec0 = ', rn_sprec0 , ' Kg/m2/s'
+         IF(lwp) WRITE(numout,*)'              sublimation                   rn_ievap0 = ', rn_ievap0 , ' Kg/m2/s'
          !
          nn_tau000 = MAX( nn_tau000, 1 )     ! must be >= 1
@@ -133,5 +158,66 @@
    END SUBROUTINE sbc_ana
 
-
+#if defined key_lim3
+   SUBROUTINE ana_ice_tau
+      !!---------------------------------------------------------------------
+      !!                     ***  ROUTINE ana_ice_tau  ***
+      !!
+      !! ** Purpose :   provide the surface boundary (momentum) condition over sea-ice
+      !!---------------------------------------------------------------------
+      utau_ice(:,:) = rn_iutau0
+      vtau_ice(:,:) = rn_ivtau0
+     
+   END SUBROUTINE ana_ice_tau
+   
+   SUBROUTINE ana_ice_flx
+      !!---------------------------------------------------------------------
+      !!                     ***  ROUTINE ana_ice_flx  ***
+      !!
+      !! ** Purpose :   provide the surface boundary (flux) condition over sea-ice
+      !!---------------------------------------------------------------------
+      REAL(wp), DIMENSION(:,:), POINTER ::   zsnw       ! snw distribution after wind blowing
+      !!---------------------------------------------------------------------
+      CALL wrk_alloc( jpi,jpj, zsnw ) 
+
+      ! ocean variables (renaming)
+      emp_oce (:,:)   = rn_emp0
+      qsr_oce (:,:)   = rn_qsr0
+      qns_oce (:,:)   = rn_qns0
+      
+      ! ice variables
+      alb_ice (:,:,:) = 0.7_wp ! useless
+      qsr_ice (:,:,:) = rn_iqsr0
+      qns_ice (:,:,:) = rn_iqns0
+      sprecip (:,:)   = rn_sprec0
+      evap_ice(:,:,:) = rn_ievap0
+
+      ! ice variables deduced from above
+      zsnw(:,:) = 0._wp
+      CALL lim_thd_snwblow( pfrld, zsnw )  ! snow distribution over ice after wind blowing 
+      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(:,:)   =   rhosn * ( sst_m(:,:) * cpic - lfus ) * tmask(:,:,1) ! in J/m3
+      qemp_oce (:,:)   = - emp_oce(:,:) * sst_m(:,:) * rcp
+      qemp_ice (:,:)   =   sprecip(:,:) * zsnw * ( sst_m(:,:) * cpic - lfus ) * tmask(:,:,1) ! solid precip (only)
+
+      ! total fluxes
+      emp_tot (:,:) = emp_ice  + emp_oce 
+      qns_tot (:,:) = pfrld(:,:) * qns_oce(:,:) + SUM( a_i_b(:,:,:) * qns_ice(:,:,:), dim=3 ) + qemp_ice(:,:) + qemp_oce(:,:)
+      qsr_tot (:,:) = pfrld(:,:) * qsr_oce(:,:) + SUM( a_i_b(:,:,:) * qsr_ice(:,:,:), dim=3 )
+
+      !--------------------------------------------------------------------
+      ! FRACTIONs of net shortwave radiation which is not absorbed in the
+      ! thin surface layer and penetrates inside the ice cover
+      ! ( Maykut and Untersteiner, 1971 ; Ebert and Curry, 1993 )
+      fr1_i0(:,:) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice )
+      fr2_i0(:,:) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice )
+
+      CALL wrk_dealloc( jpi,jpj, zsnw ) 
+      
+   END SUBROUTINE ana_ice_flx
+#endif
+
+   
    SUBROUTINE sbc_gyre( kt )
       !!---------------------------------------------------------------------
Index: branches/2016/dev_v3_6_STABLE_r6506_AGRIF_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90
===================================================================
--- branches/2016/dev_v3_6_STABLE_r6506_AGRIF_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90	(revision 7069)
+++ branches/2016/dev_v3_6_STABLE_r6506_AGRIF_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90	(revision 7077)
@@ -30,4 +30,5 @@
    USE sbcblk_clio     ! Surface boundary condition: CLIO bulk
    USE sbccpl          ! Surface boundary condition: coupled interface
+   USE sbcana          ! Surface boundary condition: analytic formulation
    USE albedo          ! ocean & ice albedo
 
@@ -109,5 +110,5 @@
       !!---------------------------------------------------------------------
       INTEGER, INTENT(in) ::   kt      ! ocean time step
-      INTEGER, INTENT(in) ::   kblk    ! type of bulk (=3 CLIO, =4 CORE, =5 COUPLED)
+      INTEGER, INTENT(in) ::   kblk    ! type of bulk (=1 ANALYTIC, =3 CLIO, =4 CORE, =5 COUPLED)
       !!
       INTEGER  ::   jl                 ! dummy loop index
@@ -156,4 +157,5 @@
          !-----------------------------------------------------------------
          SELECT CASE( kblk )
+            CASE( jp_ana     )   ;    CALL ana_ice_tau                              ! analytic formulation            
             CASE( jp_clio    )   ;    CALL blk_ice_clio_tau                         ! CLIO bulk formulation            
             CASE( jp_core    )   ;    CALL blk_ice_core_tau                         ! CORE bulk formulation
@@ -228,4 +230,7 @@
          SELECT CASE( kblk )
 
+            CASE( jp_ana )                                        ! analytic formulation
+                                      CALL ana_ice_flx
+               
             CASE( jp_clio )                                       ! CLIO bulk formulation
                ! In CLIO the cloud fraction is read in the climatology and the all-sky albedo 
