Index: branches/2016/dev_CNRS_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90
===================================================================
--- branches/2016/dev_CNRS_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90	(revision 7278)
+++ branches/2016/dev_CNRS_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90	(revision 7280)
@@ -13,5 +13,6 @@
    !!            3.4  ! 2011-11  (C. Harris) CICE added as an option
    !!            3.5  ! 2012-11  (A. Coward, G. Madec) Rethink of heat, mass and salt surface fluxes
-   !!            3.6  ! 2014-11  (P. Mathiot, C. Harris) add ice shelves melting                    
+   !!            3.6  ! 2014-11  (P. Mathiot, C. Harris) add ice shelves melting
+   !!            4.0  ! 2016-06  (L. Brodeau) new general bulk formulation
    !!----------------------------------------------------------------------
 
@@ -30,12 +31,10 @@
    USE sbcssm         ! surface boundary condition: sea-surface mean variables
    USE sbcflx         ! surface boundary condition: flux formulation
-   USE sbcblk_clio    ! surface boundary condition: bulk formulation : CLIO
-   USE sbcblk_core    ! surface boundary condition: bulk formulation : CORE
-   USE sbcblk_mfs     ! surface boundary condition: bulk formulation : MFS
+   USE sbcblk         ! surface boundary condition: bulk formulation
    USE sbcice_if      ! surface boundary condition: ice-if sea-ice model
    USE sbcice_lim     ! surface boundary condition: LIM 3.0 sea-ice model
    USE sbcice_lim_2   ! surface boundary condition: LIM 2.0 sea-ice model
    USE sbcice_cice    ! surface boundary condition: CICE    sea-ice model
-   USE sbccpl         ! surface boundary condition: coupled florulation
+   USE sbccpl         ! surface boundary condition: coupled formulation
    USE cpl_oasis3     ! OASIS routines for coupling
    USE sbcssr         ! surface boundary condition: sea surface restoring
@@ -63,7 +62,7 @@
    PUBLIC   sbc        ! routine called by step.F90
    PUBLIC   sbc_init   ! routine called by opa.F90
-   
+
    INTEGER ::   nsbc   ! type of surface boundary condition (deduced from namsbc informations)
-      
+
    !!----------------------------------------------------------------------
    !! NEMO/OPA 4.0 , NEMO-consortium (2016) 
@@ -85,16 +84,15 @@
       !!              - nsbc: type of sbc
       !!----------------------------------------------------------------------
-      INTEGER ::   icpt   ! local integer
-      !!
-      NAMELIST/namsbc/ nn_fsbc  , ln_usr   , ln_flx, ln_blk_clio, ln_blk_core, ln_blk_mfs,   &
-         &             ln_cpl   , ln_mixcpl, nn_components      , nn_limflx  ,               &
-         &             ln_traqsr, ln_dm2dc ,                                                 &  
-         &             nn_ice   , nn_ice_embd,                                               &
-         &             ln_rnf   , ln_ssr   , ln_isf   , nn_fwb    , ln_apr_dyn,              &
-         &             ln_wave  ,                                                            &
-         &             nn_lsm   
-      INTEGER  ::   ios
-      INTEGER  ::   ierr, ierr0, ierr1, ierr2, ierr3, jpm
-      LOGICAL  ::   ll_purecpl
+      INTEGER ::   ios, icpt                         ! local integer
+      LOGICAL ::   ll_purecpl, ll_opa, ll_not_nemo   ! local logical
+      !!
+      NAMELIST/namsbc/ nn_fsbc  ,                                                    &
+         &             ln_usr   , ln_flx   , ln_blk       ,                          &
+         &             ln_cpl   , ln_mixcpl, nn_components, nn_limflx,               &
+         &             nn_ice   , nn_ice_embd,                                       &
+         &             ln_traqsr, ln_dm2dc ,                                         &
+         &             ln_rnf   , nn_fwb   , ln_ssr       , ln_isf   , ln_apr_dyn,   &
+         &             ln_wave  ,                                                    &
+         &             nn_lsm
       !!----------------------------------------------------------------------
       !
@@ -105,41 +103,41 @@
       ENDIF
       !
-      REWIND( numnam_ref )       ! Namelist namsbc in reference namelist : Surface boundary
+      !                       !**  read Surface Module namelist
+      REWIND( numnam_ref )          !* Namelist namsbc in reference namelist : Surface boundary
       READ  ( numnam_ref, namsbc, IOSTAT = ios, ERR = 901)
 901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namsbc in reference namelist', lwp )
       !
-      REWIND( numnam_cfg )       ! Namelist namsbc in configuration namelist : Parameters of the run
+      REWIND( numnam_cfg )          !* Namelist namsbc in configuration namelist : Parameters of the run
       READ  ( numnam_cfg, namsbc, IOSTAT = ios, ERR = 902 )
 902   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namsbc in configuration namelist', lwp )
       IF(lwm) WRITE( numond, namsbc )
       !
-      !                          ! overwrite namelist parameter using CPP key information
+      !                             !* overwrite namelist parameter using CPP key information
       IF( Agrif_Root() ) THEN                ! AGRIF zoom
-        IF( lk_lim2 )   nn_ice      = 2
-        IF( lk_lim3 )   nn_ice      = 3
-        IF( lk_cice )   nn_ice      = 4
-      ENDIF
-      !
-      IF(lwp) THEN               ! Control print
+         IF( lk_lim2 )   nn_ice      = 2
+         IF( lk_lim3 )   nn_ice      = 3
+         IF( lk_cice )   nn_ice      = 4
+      ENDIF
+      !
+      IF(lwp) THEN                  !* Control print
          WRITE(numout,*) '   Namelist namsbc (partly overwritten with CPP key setting)'
-         WRITE(numout,*) '      Frequency update of sbc (and ice)             nn_fsbc       = ', nn_fsbc
+         WRITE(numout,*) '      frequency update of sbc (and ice)             nn_fsbc       = ', nn_fsbc
          WRITE(numout,*) '      Type of air-sea fluxes : '
          WRITE(numout,*) '         user defined formulation                   ln_usr        = ', ln_usr
          WRITE(numout,*) '         flux         formulation                   ln_flx        = ', ln_flx
-         WRITE(numout,*) '         CLIO bulk    formulation                   ln_blk_clio   = ', ln_blk_clio
-         WRITE(numout,*) '         CORE bulk    formulation                   ln_blk_core   = ', ln_blk_core
-         WRITE(numout,*) '         MFS  bulk    formulation                   ln_blk_mfs    = ', ln_blk_mfs
+         WRITE(numout,*) '         bulk         formulation                   ln_blk        = ', ln_blk
          WRITE(numout,*) '      Type of coupling (Ocean/Ice/Atmosphere) : '
          WRITE(numout,*) '         ocean-atmosphere coupled formulation       ln_cpl        = ', ln_cpl
-         WRITE(numout,*) '         forced-coupled mixed formulation           ln_mixcpl     = ', ln_mixcpl
+         WRITE(numout,*) '         mixed forced-coupled     formulation       ln_mixcpl     = ', ln_mixcpl
+!!gm  lk_oasis is controlled by key_oasis3  ===>>>  It shoud be removed from the namelist 
          WRITE(numout,*) '         OASIS coupling (with atm or sas)           lk_oasis      = ', lk_oasis
          WRITE(numout,*) '         components of your executable              nn_components = ', nn_components
          WRITE(numout,*) '         Multicategory heat flux formulation (LIM3) nn_limflx     = ', nn_limflx
          WRITE(numout,*) '      Sea-ice : '
-         WRITE(numout,*) '         ice management in the sbc (=0/1/2/3)       nn_ice        = ', nn_ice 
+         WRITE(numout,*) '         ice management in the sbc (=0/1/2/3)       nn_ice        = ', nn_ice
          WRITE(numout,*) '         ice-ocean embedded/levitating (=0/1/2)     nn_ice_embd   = ', nn_ice_embd
          WRITE(numout,*) '      Misc. options of sbc : '
          WRITE(numout,*) '         Light penetration in temperature Eq.       ln_traqsr     = ', ln_traqsr
-         WRITE(numout,*) '            daily mean to diurnal cycle qsr            ln_dm2dc   = ', ln_dm2dc 
+         WRITE(numout,*) '            daily mean to diurnal cycle qsr            ln_dm2dc   = ', ln_dm2dc
          WRITE(numout,*) '         Sea Surface Restoring on SST and/or SSS    ln_ssr        = ', ln_ssr
          WRITE(numout,*) '         FreshWater Budget control  (=0/1/2)        nn_fwb        = ', nn_fwb
@@ -149,5 +147,5 @@
          WRITE(numout,*) '         closed sea (=0/1) (set in namdom)          nn_closea     = ', nn_closea
          WRITE(numout,*) '         nb of iterations if land-sea-mask applied  nn_lsm        = ', nn_lsm
-         WRITE(numout,*) '         surface wave                               ln_wave       = ', ln_wave  
+         WRITE(numout,*) '         surface wave                               ln_wave       = ', ln_wave
       ENDIF
       !
@@ -157,104 +155,133 @@
          IF( MOD( rdt  , 2.  ) /= 0. )   CALL ctl_stop( 'the time step (in second) must be an even number'           )
       ENDIF
-      !
-      IF(lwp) THEN
-         WRITE(numout,*)
-         SELECT CASE ( nn_limflx )        ! LIM3 Multi-category heat flux formulation
-         CASE ( -1 )   ;   WRITE(numout,*) '   LIM3: use per-category fluxes (nn_limflx = -1) '
-         CASE (  0 )   ;   WRITE(numout,*) '   LIM3: use average per-category fluxes (nn_limflx = 0) ' 
-         CASE (  1 )   ;   WRITE(numout,*) '   LIM3: use average then redistribute per-category fluxes (nn_limflx = 1) '
-         CASE (  2 )   ;   WRITE(numout,*) '   LIM3: Redistribute a single flux over categories (nn_limflx = 2) '
+      !                       !**  check option consistency
+      !
+      IF(lwp) WRITE(numout,*)       !* Single / Multi - executable (NEMO / OPA+SAS) 
+      SELECT CASE( nn_components )
+      CASE( jp_iam_nemo )
+         IF(lwp) WRITE(numout,*) '   NEMO configured as a single executable (i.e. including both OPA and Surface module'
+      CASE( jp_iam_opa  )
+         IF(lwp) WRITE(numout,*) '   Multi executable configuration. Here, OPA component'
+         IF( .NOT.lk_oasis )   CALL ctl_stop( 'sbc_init : OPA-SAS coupled via OASIS, but key_oasis3 disabled' )
+         IF( ln_cpl        )   CALL ctl_stop( 'sbc_init : OPA-SAS coupled via OASIS, but ln_cpl = T in OPA'   )
+         IF( ln_mixcpl     )   CALL ctl_stop( 'sbc_init : OPA-SAS coupled via OASIS, but ln_mixcpl = T in OPA' )
+      CASE( jp_iam_sas  )
+         IF(lwp) WRITE(numout,*) '   Multi executable configuration. Here, SAS component'
+         IF( .NOT.lk_oasis )   CALL ctl_stop( 'sbc_init : OPA-SAS coupled via OASIS, but key_oasis3 disabled' )
+         IF( ln_mixcpl     )   CALL ctl_stop( 'sbc_init : OPA-SAS coupled via OASIS, but ln_mixcpl = T in OPA' )
+      CASE DEFAULT
+         CALL ctl_stop( 'sbc_init : unsupported value for nn_components' )
+      END SELECT
+      !                             !* coupled options
+      IF( ln_cpl ) THEN
+         IF( .NOT. lk_oasis )   CALL ctl_stop( 'sbc_init : coupled mode with an atmosphere model (ln_cpl=T)',   &
+            &                                  '           required to defined key_oasis3' )
+      ENDIF
+      IF( ln_mixcpl ) THEN
+         IF( .NOT. lk_oasis )   CALL ctl_stop( 'sbc_init : mixed forced-coupled mode (ln_mixcpl=T) ',   &
+            &                                  '           required to defined key_oasis3' )
+         IF( .NOT.ln_cpl    )   CALL ctl_stop( 'sbc_init : mixed forced-coupled mode (ln_mixcpl=T) requires ln_cpl = T' )
+         IF( nn_components /= jp_iam_nemo )    &
+            &                   CALL ctl_stop( 'sbc_init : the mixed forced-coupled mode (ln_mixcpl=T) ',   &
+            &                                   '          not yet working with sas-opa coupling via oasis' )
+      ENDIF
+      !                             !* sea-ice
+      SELECT CASE( nn_ice )
+      CASE( 0 )                        !- no ice in the domain
+      CASE( 1 )                        !- Ice-cover climatology ("Ice-if" model)  
+      CASE( 2 )                        !- LIM2 ice model
+         IF( .NOT.( ln_blk .OR. ln_cpl ) )   CALL ctl_stop( 'sbc_init : LIM2 sea-ice model requires ln_blk or ln_cpl = T' )
+      CASE( 3 )                        !- LIM3 ice model
+         IF( .NOT.( ln_blk .OR. ln_cpl ) )   CALL ctl_stop( 'sbc_init : LIM3 sea-ice model requires ln_blk or ln_cpl = T' )
+         IF( nn_ice_embd == 0            )   CALL ctl_stop( 'sbc_init : LIM3 sea-ice models require nn_ice_embd = 1 or 2' )
+      CASE( 4 )                        !- CICE ice model
+         IF( .NOT.( ln_blk .OR. ln_cpl ) )   CALL ctl_stop( 'sbc_init : CICE sea-ice model requires ln_blk or ln_cpl = T' )
+         IF( nn_ice_embd == 0            )   CALL ctl_stop( 'sbc_init : CICE sea-ice models require nn_ice_embd = 1 or 2' )
+         IF( lk_agrif                    )   CALL ctl_stop( 'sbc_init : CICE sea-ice model not currently available with AGRIF' ) 
+      CASE DEFAULT                     !- not supported
+      END SELECT
+      !
+      IF( nn_ice == 3 ) THEN           !- LIM3 case: multi-category flux option
+         IF(lwp) WRITE(numout,*)
+         SELECT CASE( nn_limflx )         ! LIM3 Multi-category heat flux formulation
+         CASE ( -1 )
+            IF(lwp) WRITE(numout,*) '   LIM3: use per-category fluxes (nn_limflx = -1) '
+            IF( ln_cpl )   CALL ctl_stop( 'sbc_init : the chosen nn_limflx for LIM3 in coupled mode must be 0 or 2' )
+         CASE ( 0  )
+            IF(lwp) WRITE(numout,*) '   LIM3: use average per-category fluxes (nn_limflx = 0) '
+         CASE ( 1  )
+            IF(lwp) WRITE(numout,*) '   LIM3: use average then redistribute per-category fluxes (nn_limflx = 1) '
+            IF( ln_cpl )   CALL ctl_stop( 'sbc_init : the chosen nn_limflx for LIM3 in coupled mode must be 0 or 2' )
+         CASE ( 2  )
+            IF(lwp) WRITE(numout,*) '   LIM3: Redistribute a single flux over categories (nn_limflx = 2) '
+            IF( .NOT.ln_cpl )   CALL ctl_stop( 'sbc_init : the chosen nn_limflx for LIM3 in forced mode cannot be 2' )
+         CASE DEFAULT
+            CALL ctl_stop( 'sbcmod: LIM3 option, nn_limflx, should be between -1 and 2' )
          END SELECT
-      ENDIF
-      !
-      IF( nn_components /= jp_iam_nemo .AND. .NOT. lk_oasis )   &
-         &      CALL ctl_stop( 'sbc_init : OPA-SAS coupled via OASIS, but key_oasis3 disabled' )
-      IF( nn_components == jp_iam_opa .AND. ln_cpl )   &
-         &      CALL ctl_stop( 'sbc_init : OPA-SAS coupled via OASIS, but ln_cpl = T in OPA' )
-      IF( nn_components == jp_iam_opa .AND. ln_mixcpl )   &
-         &      CALL ctl_stop( 'sbc_init : OPA-SAS coupled via OASIS, but ln_mixcpl = T in OPA' )
-      IF( ln_cpl .AND. .NOT. lk_oasis )    &
-         &      CALL ctl_stop( 'sbc_init : OASIS-coupled atmosphere model, but key_oasis3 disabled' )
-      IF( ln_mixcpl .AND. .NOT. lk_oasis )    &
-         &      CALL ctl_stop( 'the forced-coupled mixed mode (ln_mixcpl) requires the cpp key key_oasis3' )
-      IF( ln_mixcpl .AND. .NOT. ln_cpl )    &
-         &      CALL ctl_stop( 'the forced-coupled mixed mode (ln_mixcpl) requires ln_cpl = T' )
-      IF( ln_mixcpl .AND. nn_components /= jp_iam_nemo )    &
-         &      CALL ctl_stop( 'the forced-coupled mixed mode (ln_mixcpl) is not yet working with sas-opa coupling via oasis' )
-
-      !                              ! allocate sbc arrays
+      ELSE                             ! other sea-ice model
+         IF( nn_limflx >= 0  )   CALL ctl_warn( 'sbc_init : multi-category flux option (nn_limflx) only available in LIM3' )
+      ENDIF
+      !
+      !                       !**  allocate and set required variables
+      !
+      !                             !* allocate sbc arrays
       IF( sbc_oce_alloc() /= 0 )   CALL ctl_stop( 'sbc_init : unable to allocate sbc_oce arrays' )
-
-      !                          ! Checks:
-      IF( .NOT. ln_isf ) THEN                      ! variable initialisation if no ice shelf 
+      !
+      IF( .NOT.ln_isf ) THEN        !* No ice-shelf in the domain : allocate and set to zero
          IF( sbc_isf_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_init : unable to allocate sbc_isf arrays' )
-         fwfisf  (:,:)   = 0.0_wp ; fwfisf_b  (:,:)   = 0.0_wp
-         risf_tsc(:,:,:) = 0.0_wp ; risf_tsc_b(:,:,:) = 0.0_wp
+         fwfisf  (:,:)   = 0._wp   ;   risf_tsc  (:,:,:) = 0._wp
+         fwfisf_b(:,:)   = 0._wp   ;   risf_tsc_b(:,:,:) = 0._wp
       END IF
-      IF( nn_ice == 0 .AND. nn_components /= jp_iam_opa )   fr_i(:,:) = 0._wp    ! no ice in the domain, ice fraction is always zero
-
-      sfx(:,:) = 0._wp                             ! the salt flux due to freezing/melting will be computed (i.e. will be non-zero) 
-                                                   ! only if sea-ice is present
- 
-      fmmflx(:,:) = 0._wp                          ! freezing-melting array initialisation
-      
-      taum(:,:) = 0._wp                            ! Initialise taum for use in gls in case of reduced restart
-
-      !                                            ! restartability   
-      IF( ( nn_ice == 2 .OR. nn_ice ==3 ) .AND. .NOT.( ln_blk_clio .OR. ln_blk_core .OR. ln_cpl ) )   &
-         &   CALL ctl_stop( 'LIM sea-ice model requires a bulk formulation or coupled configuration' )
-      IF( nn_ice == 4 .AND. .NOT.( ln_blk_core .OR. ln_cpl ) )   &
-         &   CALL ctl_stop( 'CICE sea-ice model requires ln_blk_core or ln_cpl' )
-      IF( nn_ice == 4 .AND. lk_agrif )   &
-         &   CALL ctl_stop( 'CICE sea-ice model not currently available with AGRIF' )
-      IF( ( nn_ice == 3 .OR. nn_ice == 4 ) .AND. nn_ice_embd == 0 )   &
-         &   CALL ctl_stop( 'LIM3 and CICE sea-ice models require nn_ice_embd = 1 or 2' )
-      IF( ( nn_ice /= 3 ) .AND. ( nn_limflx >= 0 ) )   &
-         &   WRITE(numout,*) 'The nn_limflx>=0 option has no effect if sea ice model is not LIM3'
-      IF( ( nn_ice == 3 ) .AND. ( ln_cpl ) .AND. ( ( nn_limflx == -1 ) .OR. ( nn_limflx == 1 ) ) )   &
-         &   CALL ctl_stop( 'The chosen nn_limflx for LIM3 in coupled mode must be 0 or 2' )
-      IF( ( nn_ice == 3 ) .AND. ( .NOT. ln_cpl ) .AND. ( nn_limflx == 2 ) )   &
-         &   CALL ctl_stop( 'The chosen nn_limflx for LIM3 in forced mode cannot be 2' )
-
-      IF( ln_dm2dc )   nday_qsr = -1   ! initialisation flag
-
-      IF( ln_dm2dc .AND. .NOT.( ln_flx .OR. ln_blk_core ) .AND. nn_components /= jp_iam_opa )   &
-         &   CALL ctl_stop( 'diurnal cycle into qsr field from daily values requires a flux or core-bulk formulation' )
-      
-      !                          ! Choice of the Surface Boudary Condition (set nsbc)
-      ll_purecpl = ln_cpl .AND. .NOT. ln_mixcpl
-      !
+      IF( nn_ice == 0 ) THEN        !* No sea-ice in the domain : ice fraction is always zero
+         IF( nn_components /= jp_iam_opa )   fr_i(:,:) = 0._wp    ! except for OPA in SAS-OPA coupled case
+      ENDIF
+      !
+      sfx   (:,:) = 0._wp           !* salt flux due to freezing/melting
+      fmmflx(:,:) = 0._wp           !* freezing minus melting flux
+
+      taum(:,:) = 0._wp             !* wind stress module (needed in GLS in case of reduced restart)
+
+      IF( ln_dm2dc ) THEN           !* daily mean to diurnal cycle
+         nday_qsr = -1   ! allow initialization at the 1st call
+         IF( .NOT.( ln_flx .OR. ln_blk ) .AND. nn_components /= jp_iam_opa )   &
+            &   CALL ctl_stop( 'qsr diurnal cycle from daily values requires a flux or bulk formulation' )
+      ENDIF
+
+      !                             !* Choice of the Surface Boudary Condition
+      !                             (set nsbc)
+      !
+      ll_purecpl  = ln_cpl .AND. .NOT.ln_mixcpl
+      ll_opa      = nn_components == jp_iam_opa
+      ll_not_nemo = nn_components /= jp_iam_nemo
       icpt = 0
+      !
       IF( ln_usr          ) THEN   ;   nsbc = jp_usr     ; icpt = icpt + 1   ;   ENDIF       ! user defined         formulation
       IF( ln_flx          ) THEN   ;   nsbc = jp_flx     ; icpt = icpt + 1   ;   ENDIF       ! flux                 formulation
-      IF( ln_blk_clio     ) THEN   ;   nsbc = jp_clio    ; icpt = icpt + 1   ;   ENDIF       ! CLIO bulk            formulation
-      IF( ln_blk_core     ) THEN   ;   nsbc = jp_core    ; icpt = icpt + 1   ;   ENDIF       ! CORE bulk            formulation
-      IF( ln_blk_mfs      ) THEN   ;   nsbc = jp_mfs     ; icpt = icpt + 1   ;   ENDIF       ! MFS  bulk            formulation
+      IF( ln_blk          ) THEN   ;   nsbc = jp_blk     ; icpt = icpt + 1   ;   ENDIF       ! bulk                 formulation
       IF( ll_purecpl      ) THEN   ;   nsbc = jp_purecpl ; icpt = icpt + 1   ;   ENDIF       ! Pure Coupled         formulation
-      IF( nn_components == jp_iam_opa )   &
-         &                  THEN   ;   nsbc = jp_none    ; icpt = icpt + 1   ;   ENDIF       ! opa coupling via SAS module
-      !
-      IF( icpt /= 1 )    CALL ctl_stop( 'sbc_init: choose ONE and only ONE sbc option' )
-      !
-      IF(lwp) THEN
+      IF( ll_opa          ) THEN   ;   nsbc = jp_none    ; icpt = icpt + 1   ;   ENDIF       ! opa coupling via SAS module
+      !
+      IF( icpt /= 1 )    CALL ctl_stop( 'sbc_init : choose ONE and only ONE sbc option' )
+      !
+      IF(lwp) THEN                     !- print the choice of surface flux formulation
          WRITE(numout,*)
          SELECT CASE( nsbc )
          CASE( jp_usr     )   ;   WRITE(numout,*) '   user defined formulation'
-         CASE( jp_flx     )   ;   WRITE(numout,*) '   flux formulation'
-         CASE( jp_clio    )   ;   WRITE(numout,*) '   CLIO bulk formulation'
-         CASE( jp_core    )   ;   WRITE(numout,*) '   CORE bulk formulation'
-         CASE( jp_purecpl )   ;   WRITE(numout,*) '   pure coupled formulation'
-         CASE( jp_mfs     )   ;   WRITE(numout,*) '   MFS Bulk formulation'
-         CASE( jp_none    )   ;   WRITE(numout,*) '   OPA coupled to SAS via oasis'
-            IF( ln_mixcpl )       WRITE(numout,*) '       + forced-coupled mixed formulation'
+         CASE( jp_flx     )   ;   WRITE(numout,*) '      ===>>   flux formulation'
+         CASE( jp_blk     )   ;   WRITE(numout,*) '      ===>>   bulk formulation'
+         CASE( jp_purecpl )   ;   WRITE(numout,*) '      ===>>   pure coupled formulation'
+!!gm abusive use of jp_none ??   ===>>> need to be check and changed by adding a jp_sas parameter
+         CASE( jp_none    )   ;   WRITE(numout,*) '      ===>>   OPA coupled to SAS via oasis'
+            IF( ln_mixcpl )       WRITE(numout,*) '                  + forced-coupled mixed formulation'
          END SELECT
-         IF( nn_components/= jp_iam_nemo )  &
-            &                     WRITE(numout,*) '       + OASIS coupled SAS'
-      ENDIF
-      !
-      IF( lk_oasis )   CALL sbc_cpl_init (nn_ice)   ! OASIS initialisation. must be done before: (1) first time step
-      !                                             !                                            (2) the use of nn_fsbc
+         IF( ll_not_nemo  )       WRITE(numout,*) '                  + OASIS coupled SAS'
+      ENDIF
+      !
+      !                             !* OASIS initialization
+      !
+      IF( lk_oasis )   CALL sbc_cpl_init( nn_ice )   ! Must be done before: (1) first time step
+      !                                              !                      (2) the use of nn_fsbc
       !     nn_fsbc initialization if OPA-SAS coupling via OASIS
-      !     sas model time step has to be declared in OASIS (mandatory) -> nn_fsbc has to be modified accordingly
+      !     SAS time-step has to be declared in OASIS (mandatory) -> nn_fsbc has to be modified accordingly
       IF( nn_components /= jp_iam_nemo ) THEN
          IF( nn_components == jp_iam_opa )   nn_fsbc = cpl_freq('O_SFLX') / NINT(rdt)
@@ -268,7 +295,8 @@
       ENDIF
       !
+      !                             !* check consistency between model timeline and nn_fsbc
       IF( MOD( nitend - nit000 + 1, nn_fsbc) /= 0 .OR.   &
-          MOD( nstock             , nn_fsbc) /= 0 ) THEN 
-         WRITE(ctmp1,*) 'experiment length (', nitend - nit000 + 1, ') or nstock (', nstock,   &
+          MOD( nstock             , nn_fsbc) /= 0 ) THEN
+         WRITE(ctmp1,*) 'sbc_init : experiment length (', nitend - nit000 + 1, ') or nstock (', nstock,   &
             &           ' is NOT a multiple of nn_fsbc (', nn_fsbc, ')'
          CALL ctl_stop( ctmp1, 'Impossible to properly do model restart' )
@@ -276,18 +304,23 @@
       !
       IF( MOD( rday, REAL(nn_fsbc, wp) * rdt ) /= 0 )   &
-         &  CALL ctl_warn( 'nn_fsbc is NOT a multiple of the number of time steps in a day' )
-      !
-      IF( ln_dm2dc .AND. ( ( NINT(rday) / ( nn_fsbc * NINT(rdt) ) )  < 8 ) )   &
-         &   CALL ctl_warn( 'diurnal cycle for qsr: the sampling of the diurnal cycle is too small...' )
-      !
-                          CALL sbc_ssm_init               ! Sea-surface mean fields initialisation
-      !
-      IF( ln_ssr      )   CALL sbc_ssr_init               ! Sea-Surface Restoring initialisation
-      !
-                          CALL sbc_rnf_init               ! Runof initialisation
-      !
-      IF( nn_ice == 3 )   CALL sbc_lim_init               ! LIM3 initialisation
-      !
-      IF( nn_ice == 4 )   CALL cice_sbc_init( nsbc )      ! CICE initialisation
+         &  CALL ctl_warn( 'sbc_init : nn_fsbc is NOT a multiple of the number of time steps in a day' )
+      !
+      IF( ln_dm2dc .AND. NINT(rday) / ( nn_fsbc * NINT(rdt) ) < 8  )   &
+         &   CALL ctl_warn( 'sbc_init : diurnal cycle for qsr: the sampling of the diurnal cycle is too small...' )
+      !
+   
+      !                       !**  associated modules : initialization
+      !
+                          CALL sbc_ssm_init            ! Sea-surface mean fields initialization
+      !
+      IF( ln_blk      )   CALL sbc_blk_init            ! bulk formulae initialization
+
+      IF( ln_ssr      )   CALL sbc_ssr_init            ! Sea-Surface Restoring initialization
+      !
+                          CALL sbc_rnf_init            ! Runof initialization
+      !
+      IF( nn_ice == 3 )   CALL sbc_lim_init            ! LIM3 initialization
+      !
+      IF( nn_ice == 4 )   CALL cice_sbc_init( nsbc )   ! CICE initialization
       !
    END SUBROUTINE sbc_init
@@ -297,18 +330,20 @@
       !!---------------------------------------------------------------------
       !!                    ***  ROUTINE sbc  ***
-      !!              
+      !!
       !! ** Purpose :   provide at each time-step the ocean surface boundary
       !!                condition (momentum, heat and freshwater fluxes)
       !!
-      !! ** Method  :   blah blah  to be written ????????? 
+      !! ** Method  :   blah blah  to be written ?????????
       !!                CAUTION : never mask the surface stress field (tke sbc)
       !!
-      !! ** Action  : - set the ocean surface boundary condition at before and now 
-      !!                time step, i.e.  
+      !! ** Action  : - set the ocean surface boundary condition at before and now
+      !!                time step, i.e.
       !!                utau_b, vtau_b, qns_b, qsr_b, emp_n, sfx_b, qrp_b, erp_b
       !!                utau  , vtau  , qns  , qsr  , emp  , sfx  , qrp  , erp
       !!              - updte the ice fraction : fr_i
       !!----------------------------------------------------------------------
-      INTEGER, INTENT(in) ::   kt       ! ocean time step
+      INTEGER, INTENT(in) ::   kt   ! ocean time step
+      !
+      LOGICAL ::   ll_sas, ll_opa   ! local logical
       !!---------------------------------------------------------------------
       !
@@ -332,30 +367,27 @@
       !                                            ! ---------------------------------------- !
       !
-      IF( nn_components /= jp_iam_sas )   CALL sbc_ssm ( kt )  ! ocean sea surface variables (sst_m, sss_m, ssu_m, ssv_m)
-      !                                                        ! averaged over nf_sbc time-step
-      IF( ln_wave                     )   CALL sbc_wave( kt )  ! surface waves
-      
-      
-                                                   !==  sbc formulation  ==!
-                                                            
+      ll_sas = nn_components == jp_iam_sas               ! component flags
+      ll_opa = nn_components == jp_iam_opa
+      !
+      IF( .NOT.ll_sas )   CALL sbc_ssm ( kt )            ! mean ocean sea surface variables (sst_m, sss_m, ssu_m, ssv_m)
+      IF( ln_wave     )   CALL sbc_wave( kt )            ! surface waves
+
+      !
+      !                                            !==  sbc formulation  ==!
+      !                                                   
       SELECT CASE( nsbc )                                ! Compute ocean surface boundary condition
       !                                                  ! (i.e. utau,vtau, qns, qsr, emp, sfx)
-      CASE( jp_usr   )   ;   CALL usr_def_sbc ( kt )                    ! user defined formulation 
-      CASE( jp_flx   )   ;   CALL sbc_flx     ( kt )                    ! flux formulation
-      CASE( jp_clio  )   ;   CALL sbc_blk_clio( kt )                    ! bulk formulation : CLIO for the ocean
-      CASE( jp_core  )   
-         IF( nn_components == jp_iam_sas ) &
-            &                CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice )   ! OPA-SAS coupling: SAS receiving fields from OPA 
-                             CALL sbc_blk_core( kt )                    ! bulk formulation : CORE for the ocean
-                                                                        ! from oce: sea surface variables (sst_m, sss_m,  ssu_m,  ssv_m)
-      CASE( jp_purecpl )  ;  CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice )   ! pure coupled formulation
-                                                                        !
-      CASE( jp_mfs   )   ;   CALL sbc_blk_mfs ( kt )                    ! bulk formulation : MFS for the ocean
-      CASE( jp_none  ) 
-         IF( nn_components == jp_iam_opa )   &
-            &                CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice )   ! OPA-SAS coupling: OPA receiving fields from SAS
+      CASE( jp_usr   )     ;   CALL usr_def_sbc( kt )                    ! user defined formulation 
+      CASE( jp_flx     )   ;   CALL sbc_flx    ( kt )                    ! flux formulation
+      CASE( jp_blk     )
+         IF( ll_sas    )       CALL sbc_cpl_rcv( kt, nn_fsbc, nn_ice )   ! OPA-SAS coupling: SAS receiving fields from OPA
+                               CALL sbc_blk    ( kt )                    ! bulk formulation for the ocean
+                               !
+      CASE( jp_purecpl )   ;   CALL sbc_cpl_rcv( kt, nn_fsbc, nn_ice )   ! pure coupled formulation
+      CASE( jp_none    )
+         IF( ll_opa    )       CALL sbc_cpl_rcv( kt, nn_fsbc, nn_ice )   ! OPA-SAS coupling: OPA receiving fields from SAS
       END SELECT
 
-      IF( ln_mixcpl )        CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice )   ! forced-coupled mixed formulation after forcing
+      IF( ln_mixcpl )          CALL sbc_cpl_rcv( kt, nn_fsbc, nn_ice )   ! forced-coupled mixed formulation after forcing
 
       !
@@ -367,5 +399,5 @@
       CASE(  3 )   ;         CALL sbc_ice_lim  ( kt, nsbc )          ! LIM-3 ice model
       CASE(  4 )   ;         CALL sbc_ice_cice ( kt, nsbc )          ! CICE ice model
-      END SELECT                                              
+      END SELECT
 
       IF( ln_icebergs    )   CALL icb_stp( kt )                   ! compute icebergs
@@ -374,11 +406,10 @@
 
       IF( ln_rnf         )   CALL sbc_rnf( kt )                   ! add runoffs to fresh water fluxes
- 
+
       IF( ln_ssr         )   CALL sbc_ssr( kt )                   ! add SST/SSS damping term
 
       IF( nn_fwb    /= 0 )   CALL sbc_fwb( kt, nn_fwb, nn_fsbc )  ! control the freshwater budget
 
-      ! treatment of closed sea in the model domain 
-      ! (update freshwater fluxes)
+      ! treatment of closed sea in the model domain   (update freshwater fluxes)
       ! Should not be ran if ln_diurnal_only
       IF( .NOT.ln_diurnal_only .AND. nn_closea == 1 )   CALL sbc_clo( kt, cn_cfg, nn_cfg )   
@@ -391,5 +422,5 @@
          !                                             ! ---------------------------------------- !
          IF( ln_rstart .AND.    &                               !* Restart: read in restart file
-            & iom_varid( numror, 'utau_b', ldstop = .FALSE. ) > 0 ) THEN 
+            & iom_varid( numror, 'utau_b', ldstop = .FALSE. ) > 0 ) THEN
             IF(lwp) WRITE(numout,*) '          nit000-1 surface forcing fields red in the restart file'
             CALL iom_get( numror, jpdom_autoglo, 'utau_b', utau_b )   ! before i-stress  (U-point)
@@ -407,9 +438,9 @@
          ELSE                                                   !* no restart: set from nit000 values
             IF(lwp) WRITE(numout,*) '          nit000-1 surface forcing fields set to nit000'
-            utau_b(:,:) = utau(:,:) 
+            utau_b(:,:) = utau(:,:)
             vtau_b(:,:) = vtau(:,:)
             qns_b (:,:) = qns (:,:)
-            emp_b (:,:) = emp(:,:)
-            sfx_b (:,:) = sfx(:,:)
+            emp_b (:,:) = emp (:,:)
+            sfx_b (:,:) = sfx (:,:)
          ENDIF
       ENDIF
@@ -435,18 +466,16 @@
          CALL iom_put( "empmr"  , emp    - rnf )                ! upward water flux
          CALL iom_put( "empbmr" , emp_b  - rnf )                ! before upward water flux ( needed to recalculate the time evolution of ssh in offline )
-         CALL iom_put( "saltflx", sfx  )                        ! downward salt flux  
-                                                                ! (includes virtual salt flux beneath ice 
-                                                                ! in linear free surface case)
+         CALL iom_put( "saltflx", sfx  )                        ! downward salt flux (includes virtual salt flux beneath ice in linear free surface case)
          CALL iom_put( "fmmflx", fmmflx  )                      ! Freezing-melting water flux
-         CALL iom_put( "qt"    , qns  + qsr )                   ! total heat flux 
+         CALL iom_put( "qt"    , qns  + qsr )                   ! total heat flux
          CALL iom_put( "qns"   , qns        )                   ! solar heat flux
          CALL iom_put( "qsr"   ,       qsr  )                   ! solar heat flux
-         IF( nn_ice > 0 .OR. nn_components == jp_iam_opa )   CALL iom_put( "ice_cover", fr_i )   ! ice fraction 
-         CALL iom_put( "taum"  , taum       )                   ! wind stress module 
+         IF( nn_ice > 0 .OR. ll_opa )   CALL iom_put( "ice_cover", fr_i )   ! ice fraction
+         CALL iom_put( "taum"  , taum       )                   ! wind stress module
          CALL iom_put( "wspd"  , wndm       )                   ! wind speed  module over free ocean or leads in presence of sea-ice
       ENDIF
       !
-      CALL iom_put( "utau", utau )   ! i-wind stress   (stress can be updated at 
-      CALL iom_put( "vtau", vtau )   ! j-wind stress    each time step in sea-ice)
+      CALL iom_put( "utau", utau )   ! i-wind stress   (stress can be updated at each time step in sea-ice)
+      CALL iom_put( "vtau", vtau )   ! j-wind stress
       !
       IF(ln_ctl) THEN         ! print mean trends (used for debugging)
