Index: /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/CONFIG/SHARED/namelist_ref
===================================================================
--- /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/CONFIG/SHARED/namelist_ref	(revision 8755)
+++ /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/CONFIG/SHARED/namelist_ref	(revision 8756)
@@ -393,4 +393,5 @@
    sn_rcv_wnum   =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
    sn_rcv_tauoc  =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
+   sn_rcv_tauw   =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
    sn_rcv_wdrag  =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
 !
@@ -1295,4 +1296,5 @@
    ln_stcor    = .false.   !  Activate Stokes Coriolis term
    ln_tauoc    = .false.   !  Activate ocean stress modified by external wave induced stress
+   ln_tauw     = .false.   !  Activate ocean stress components from wave model
    ln_phioc    = .false.   !  Activate wave to ocean energy
    ln_rough    = .false.   !  Wave roughness equals the significant wave height
@@ -1312,4 +1314,6 @@
    sn_wnum     =  'sdw_wave' ,        1          , 'wave_num'   ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
    sn_tauoc    =  'sdw_wave' ,        1          , 'wave_stress',     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
+   sn_tauwx    =  'sdw_wave' ,        1          , 'wave_stress',     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
+   sn_tauwy    =  'sdw_wave' ,        1          , 'wave_stress',     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
    sn_phioc    =  'sdw_wave' ,        1          , 'wave_energy',     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
    cn_dir      = './'  !  root directory for the location of wave forcing files
Index: /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90
===================================================================
--- /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90	(revision 8755)
+++ /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90	(revision 8756)
@@ -69,4 +69,5 @@
    LOGICAL , PUBLIC ::   ln_stcor       !: true if Stokes-Coriolis term is used 
    LOGICAL , PUBLIC ::   ln_tauoc       !: true if normalized stress from wave is used 
+   LOGICAL , PUBLIC ::   ln_tauw        !: true if ocean stress components from wave is used 
    LOGICAL , PUBLIC ::   ln_phioc       !: true if wave energy to ocean is used 
    LOGICAL , PUBLIC ::   ln_rough       !: true if wave roughness equals significant wave height
Index: /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90
===================================================================
--- /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90	(revision 8755)
+++ /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90	(revision 8756)
@@ -116,5 +116,7 @@
    INTEGER, PARAMETER ::   jpr_wdrag  = 51            ! Neutral surface drag coefficient 
    INTEGER, PARAMETER ::   jpr_wfreq  = 52            ! Wave peak frequency 
-   INTEGER, PARAMETER ::   jprcv      = 52            ! total number of fields received
+   INTEGER, PARAMETER ::   jpr_tauwx  = 53            ! x component of the ocean stress from waves
+   INTEGER, PARAMETER ::   jpr_tauwy  = 54            ! y component of the ocean stress from waves
+   INTEGER, PARAMETER ::   jprcv      = 54            ! total number of fields received
 
    INTEGER, PARAMETER ::   jps_fice   =  1            ! ice fraction sent to the atmosphere
@@ -169,5 +171,6 @@
    ! Received from waves  
    TYPE(FLD_C) ::   sn_rcv_hsig,sn_rcv_phioc,sn_rcv_sdrft,sn_rcv_wper, &
-                    sn_rcv_wfreq,sn_rcv_wnum,sn_rcv_tauoc,sn_rcv_wdrag
+                    sn_rcv_wfreq,sn_rcv_wnum,sn_rcv_tauoc,sn_rcv_tauw, &
+                    sn_rcv_wdrag
    ! Other namelist parameters                        !
    INTEGER     ::   nn_cplmodel            ! Maximum number of models to/from which NEMO is potentialy sending/receiving data
@@ -246,5 +249,6 @@
          &                  sn_rcv_sdrft, sn_rcv_wper  , sn_rcv_wnum  , sn_rcv_wfreq, sn_rcv_tauoc,    &
          &                  sn_rcv_wdrag, sn_rcv_qns   , sn_rcv_emp   , sn_rcv_rnf  , sn_rcv_cal ,     &
-         &                  sn_rcv_iceflx, sn_rcv_co2   , sn_rcv_mslp  , nn_cplmodel, ln_usecplmask
+         &                  sn_rcv_iceflx, sn_rcv_co2   , sn_rcv_mslp , sn_rcv_tauw ,                  &
+         &                  nn_cplmodel, ln_usecplmask
       !!---------------------------------------------------------------------
       !
@@ -295,4 +299,5 @@
          WRITE(numout,*)'      Wave peak frequency             = ', TRIM(sn_rcv_wfreq%cldes ), ' (', TRIM(sn_rcv_wfreq%clcat ), ')'  
          WRITE(numout,*)'      Stress frac adsorbed by waves   = ', TRIM(sn_rcv_tauoc%cldes ), ' (', TRIM(sn_rcv_tauoc%clcat ), ')'  
+         WRITE(numout,*)'      Stress components by waves      = ', TRIM(sn_rcv_tauw%cldes  ), ' (', TRIM(sn_rcv_tauw%clcat  ), ')'  
          WRITE(numout,*)'      Neutral surf drag coefficient   = ', TRIM(sn_rcv_wdrag%cldes ), ' (', TRIM(sn_rcv_wdrag%clcat ), ')'
          WRITE(numout,*)'  sent fields (multiple ice categories)'
@@ -573,4 +578,11 @@
          cpl_tauoc = .TRUE. 
       ENDIF 
+      srcv(jpr_tauwx)%clname = 'O_Tauwx'      ! ocean stress from wave in the x direction
+      srcv(jpr_tauwy)%clname = 'O_Tauwy'      ! ocean stress from wave in the y direction
+      IF( TRIM(sn_rcv_tauw%cldes ) == 'coupled' )  THEN 
+         srcv(jpr_tauwx)%laction = .TRUE. 
+         srcv(jpr_tauwy)%laction = .TRUE. 
+         cpl_tauw = .TRUE. 
+      ENDIF 
       srcv(jpr_wdrag)%clname = 'O_WDrag'     ! neutral surface drag coefficient 
       IF( TRIM(sn_rcv_wdrag%cldes ) == 'coupled' )  THEN 
@@ -578,4 +590,9 @@
          cpl_wdrag = .TRUE. 
       ENDIF 
+      !
+      IF( srcv(jpr_tauoc)%laction .AND. srcv(jpr_tauwx)%laction .AND. srcv(jpr_tauwy)%laction ) &
+            CALL ctl_stop( 'More than one method for modifying the ocean stress has been selected ', &
+                                     '(sn_rcv_tauoc=coupled and sn_rcv_tauw=coupled)' )
+      !
       !  
       !                                                      ! ------------------------------- !
@@ -1205,4 +1222,16 @@
          WHERE(tauoc_wave <   0.0 ) tauoc_wave = 1.0
          WHERE(tauoc_wave > 100.0 ) tauoc_wave = 1.0
+      ENDIF
+      !                                                      ! ========================= !  
+      !                                                      ! Stress component by waves ! 
+      !                                                      ! ========================= !  
+      IF( srcv(jpr_tauwx)%laction .AND. srcv(jpr_tauwy)%laction .AND. ln_tauw ) THEN
+         tauw_x(:,:) = frcv(jpr_tauwx)%z3(:,:,1)
+         tauw_y(:,:) = frcv(jpr_tauwy)%z3(:,:,1)
+         ! cap the value of tauoc
+         WHERE(tauw_x < -100.0 ) tauw_x = 0.0
+         WHERE(tauw_x >  100.0 ) tauw_x = 0.0
+         WHERE(tauw_y < -100.0 ) tauw_y = 0.0
+         WHERE(tauw_y >  100.0 ) tauw_y = 0.0
       ENDIF
       
Index: /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcflx.F90
===================================================================
--- /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcflx.F90	(revision 8755)
+++ /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcflx.F90	(revision 8756)
@@ -30,10 +30,11 @@
    PUBLIC sbc_flx       ! routine called by step.F90
 
-   INTEGER , PARAMETER ::   jpfld   = 5   ! maximum number of files to read 
-   INTEGER , PARAMETER ::   jp_utau = 1   ! index of wind stress (i-component) file
-   INTEGER , PARAMETER ::   jp_vtau = 2   ! index of wind stress (j-component) file
-   INTEGER , PARAMETER ::   jp_qtot = 3   ! index of total (non solar+solar) heat file
-   INTEGER , PARAMETER ::   jp_qsr  = 4   ! index of solar heat file
-   INTEGER , PARAMETER ::   jp_emp  = 5   ! index of evaporation-precipation file
+   INTEGER             ::   jpfld   = 6   ! maximum number of files to read 
+   INTEGER             ::   jp_utau       ! index of wind stress (i-component) file
+   INTEGER             ::   jp_vtau       ! index of wind stress (j-component) file
+   INTEGER             ::   jp_qtot       ! index of total (non solar+solar) heat file
+   INTEGER             ::   jp_qsr        ! index of solar heat file
+   INTEGER             ::   jp_emp        ! index of evaporation-precipation file
+   INTEGER             ::   jp_press      ! index of pressure for UKMO shelf fluxes
    TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf    ! structure of input fields (file informations, fields read)
    LOGICAL , PUBLIC    ::   ln_shelf_flx = .FALSE. ! UKMO SHELF specific flux flag
@@ -92,4 +93,5 @@
       CHARACTER(len=100) ::  cn_dir                               ! Root directory for location of flx files
       NAMELIST/namsbc_flx/ ln_shelf_flx, ln_rel_wind, rn_wfac     ! Put here to allow merging with another UKMO branch
+      LOGICAL  :: ln_readtau                                      ! Is it necessary to read tau from file?
       TYPE(FLD_N), DIMENSION(jpfld) ::   slf_i                    ! array of namelist information structures
       TYPE(FLD_N) ::   sn_utau, sn_vtau, sn_qtot, sn_qsr, sn_emp  ! informations about the fields to be read
@@ -97,4 +99,21 @@
       !!---------------------------------------------------------------------
       !
+      ln_readtau = .NOT. (ln_wave .AND. ln_tauw )
+
+      ! prepare the index of the fields that have to be read
+      jpfld = 0
+      IF( ln_readtau ) THEN
+         jp_utau = jpfld+1
+         jp_vtau = jpfld+2
+         jpfld = jpfld+2
+      ELSE
+         jp_utau = 0   ;  jp_vtau = 0
+      ENDIF
+      jp_qtot = jpfld+1
+      jp_qsr = jpfld+2
+      jp_emp = jpfld+3
+      jp_press = jpfld+4
+      jpfld = jpfld+4
+
       IF( kt == nit000 ) THEN                ! First call kt=nit000  
          ! set file information
@@ -122,5 +141,7 @@
          !
          !                                         ! store namelist information in an array
+         IF( ln_readtau ) THEN
          slf_i(jp_utau) = sn_utau   ;   slf_i(jp_vtau) = sn_vtau
+         ENDIF
          slf_i(jp_qtot) = sn_qtot   ;   slf_i(jp_qsr ) = sn_qsr 
          slf_i(jp_emp ) = sn_emp
@@ -146,5 +167,5 @@
 
          !!UKMO SHELF wind speed relative to surface currents - put here to allow merging with coupling branch
-         IF( ln_shelf_flx ) THEN
+         IF( ln_shelf_flx .AND. ln_readtau ) THEN
             CALL wrk_alloc( jpi,jpj, zwnd_i, zwnd_j )
 
@@ -178,5 +199,5 @@
          !                                                        ! add modification due to drag coefficient read from wave forcing
          !                                                        ! this code is inefficient but put here to allow merging with another UKMO branch
-         IF( ln_shelf_flx ) THEN
+         IF( ln_shelf_flx .AND. ln_readtau ) THEN
             IF( ln_cdgw .AND. nn_drag == jp_std ) THEN
                IF( cpl_wdrag ) THEN 
@@ -212,4 +233,5 @@
          !
          !                                                        ! module of wind stress and wind speed at T-point
+         IF( ln_readtau ) THEN
          zcoef = 1. / ( zrhoa * zcdrag )
 !CDIR NOVERRCHK
@@ -233,4 +255,5 @@
          taum(:,:) = taum(:,:) * tmask(:,:,1) ; wndm(:,:) = wndm(:,:) * tmask(:,:,1)
          CALL lbc_lnk( taum(:,:), 'T', 1. )   ;   CALL lbc_lnk( wndm(:,:), 'T', 1. )
+         ENDIF
 
          IF( nitend-nit000 <= 100 .AND. lwp ) THEN                ! control print (if less than 100 time-step asked)
@@ -248,5 +271,5 @@
          ENDIF
          !
-         IF( ln_shelf_flx ) THEN
+         IF( ln_shelf_flx .AND. ln_readtau ) THEN
             CALL wrk_dealloc( jpi,jpj, zwnd_i, zwnd_j )
          ENDIF
Index: /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90
===================================================================
--- /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90	(revision 8755)
+++ /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90	(revision 8756)
@@ -376,20 +376,5 @@
       IF( ln_mixcpl .OR. ln_wavcpl )  CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice )   ! forced-coupled mixed formulation after forcing
 
-      IF ( ln_wave .AND. ln_tauoc) THEN                 ! Wave stress subctracted 
-            utau(:,:) = utau(:,:)*tauoc_wave(:,:) 
-            vtau(:,:) = vtau(:,:)*tauoc_wave(:,:) 
-            taum(:,:) = taum(:,:)*tauoc_wave(:,:) 
-      ! 
-            SELECT CASE( nsbc ) 
-            CASE(  0,1,2,3,5,-1 )  ; 
-                IF(lwp .AND. kt == nit000 ) THEN
-                   write(numout,*)
-                   WRITE(numout,*) 'WARNING: You are subtracting the wave stress to the ocean. & 
-                                                      & If not requested select ln_tauoc=.false' 
-                   write(numout,*)
-                END IF
-            END SELECT 
-      ! 
-      END IF
+      IF ( ln_wave .AND. (ln_tauoc .OR. ln_tauw) ) CALL sbc_stress( )   ! Wave stress update 
 
       !                                            !==  Misc. Options  ==!
Index: /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90
===================================================================
--- /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90	(revision 8755)
+++ /branches/UKMO/r6232_HZG_WAVE-coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90	(revision 8756)
@@ -32,4 +32,5 @@
 
    PUBLIC   sbc_stokes      ! routine called in sbccpl 
+   PUBLIC   sbc_stress      ! routine called in sbcmod
    PUBLIC   sbc_wave        ! routine called in sbcmod 
    PUBLIC   sbc_wave_init   ! routine called in sbcmod
@@ -43,4 +44,5 @@
    LOGICAL, PUBLIC ::   cpl_wnum   = .FALSE.
    LOGICAL, PUBLIC ::   cpl_tauoc  = .FALSE.
+   LOGICAL, PUBLIC ::   cpl_tauw   = .FALSE.
    LOGICAL, PUBLIC ::   cpl_wdrag  = .FALSE.
 
@@ -61,4 +63,5 @@
    TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::  sf_wn    ! structure of input fields (file informations, fields read) wave number for Qiao
    TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::  sf_tauoc ! structure of input fields (file informations, fields read) normalized wave stress into the ocean
+   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::  sf_tauw ! structure of input fields (file informations, fields read) ocean stress components from wave model
    TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::  sf_phioc ! structure of input fields (file informations, fields read) wave to ocean energy 
    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   cdn_wave            !:
@@ -67,4 +70,6 @@
    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   rn_crban            !: Craig and Banner constant for surface breaking waves mixing
    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   tauoc_wave          !:
+   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   tauw_x              !:
+   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   tauw_y              !:
    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   tsd2d               !:
    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)   ::   div_sd              !: barotropic stokes drift divergence
@@ -232,4 +237,41 @@
 
 
+   SUBROUTINE sbc_stress( )
+      !!---------------------------------------------------------------------
+      !!                     ***  ROUTINE sbc_stress  ***
+      !!
+      !! ** Purpose :   Updates the ocean momentum modified by waves
+      !!
+      !! ** Method  : - Calculate u,v components of stress depending on stress
+      !!                model 
+      !!              - Calculate the stress module
+      !!              - The wind module is not modified by waves 
+      !! ** action  
+      !!---------------------------------------------------------------------
+      INTEGER  ::   jj, ji   ! dummy loop argument
+      !
+      IF( ln_tauoc ) THEN
+         utau(:,:) = utau(:,:)*tauoc_wave(:,:)
+         vtau(:,:) = vtau(:,:)*tauoc_wave(:,:)
+         taum(:,:) = taum(:,:)*tauoc_wave(:,:)
+      ENDIF
+      !
+      IF( ln_tauw ) THEN
+         DO jj = 1, jpjm1
+            DO ji = 1, jpim1
+               ! Stress components at u- & v-points
+               utau(ji,jj) = 0.5_wp * ( tauw_x(ji,jj) + tauw_x(ji+1,jj) )
+               vtau(ji,jj) = 0.5_wp * ( tauw_y(ji,jj) + tauw_y(ji,jj+1) )
+               !
+               ! Stress module at t points
+               taum(ji,jj) = SQRT( tauw_x(ji,jj)*tauw_x(ji,jj) + tauw_y(ji,jj)*tauw_y(ji,jj) )
+            END DO
+         END DO
+
+      ENDIF
+      !
+   END SUBROUTINE sbc_stress
+
+
    SUBROUTINE sbc_wave( kt )
       !!---------------------------------------------------------------------
@@ -263,4 +305,15 @@
          WHERE( tauoc_wave <   0.0 ) tauoc_wave = 1.0
          WHERE( tauoc_wave > 100.0 ) tauoc_wave = 1.0
+      ENDIF
+
+      IF( ln_tauw .AND. .NOT. cpl_tauw ) THEN      !==  Wave induced stress  ==!
+         CALL fld_read( kt, nn_fsbc, sf_tauw )           ! read ocean stress components from external forcing (T grid)
+         tauw_x(:,:) = sf_tauw(1)%fnow(:,:,1)
+         WHERE( tauw_x < -100.0 ) tauw_x = 0.0
+         WHERE( tauw_x >  100.0 ) tauw_x = 0.0
+
+         tauw_y(:,:) = sf_tauw(2)%fnow(:,:,1)
+         WHERE( tauw_y < -100.0 ) tauw_y = 0.0
+         WHERE( tauw_y >  100.0 ) tauw_y = 0.0
       ENDIF
 
@@ -343,11 +396,12 @@
       !!
       CHARACTER(len=100)     ::  cn_dir                          ! Root directory for location of drag coefficient files
-      TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) ::   slf_i     ! array of namelist informations on the fields to read
+      TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) ::   slf_i, slf_j     ! array of namelist informations on the fields to read
       TYPE(FLD_N)            ::  sn_cdg, sn_usd, sn_vsd, sn_phioc, &
                              &   sn_hsw, sn_wmp, sn_wfr, sn_wnum , &
-                             &   sn_tauoc      ! informations about the fields to be read
-      !
-      NAMELIST/namsbc_wave/  sn_cdg, cn_dir, sn_usd, sn_vsd, sn_hsw, sn_wmp, sn_wfr, sn_wnum, sn_tauoc, sn_phioc, &
-                             ln_cdgw, ln_sdw, ln_stcor, ln_phioc, ln_tauoc, ln_zdfqiao, ln_rough,                 &
+                             &   sn_tauoc, sn_tauwx, sn_tauwy      ! information about the fields to be read
+      !
+      NAMELIST/namsbc_wave/  sn_cdg, cn_dir, sn_usd, sn_vsd, sn_hsw, sn_wmp, sn_wfr, sn_wnum, sn_phioc,          &
+                             sn_tauoc, sn_tauwx, sn_tauwy,                                                       &
+                             ln_cdgw, ln_sdw, ln_stcor, ln_phioc, ln_tauoc, ln_tauw, ln_zdfqiao, ln_rough,       &
                              nn_sdrift, nn_wmix
       !!---------------------------------------------------------------------
@@ -371,4 +425,5 @@
          WRITE(numout,*) '      Stokes coriolis term                    ln_stcor    = ', ln_stcor 
          WRITE(numout,*) '      wave modified ocean stress              ln_tauoc    = ', ln_tauoc 
+         WRITE(numout,*) '      wave modified ocean stress components   ln_tauw     = ', ln_tauw 
          WRITE(numout,*) '      wave to ocean energy                    ln_phioc    = ', ln_phioc
          WRITE(numout,*) '        vertical mixing parametrization       nn_wmix     = ', nn_wmix 
@@ -380,10 +435,11 @@
       IF ( ln_wave ) THEN
          ! Activated wave physics but no wave physics components activated 
-         IF ( .NOT.(ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_stcor .OR. ln_phioc .OR. ln_rough .OR. ln_zdfqiao) )   THEN  
-             CALL ctl_warn( 'Ask for wave coupling but ln_cdgw=F, ln_sdw=F, ln_tauoc=F, ln_stcor=F, ln_phioc=F ', &
-                                                      'ln_rough=F, ln_zdfqiao=F' ) 
-         ELSEIF (ln_stcor .AND. .NOT. ln_sdw) THEN  
+         IF ( .NOT.(ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_tauw .OR. ln_stcor .OR. ln_phioc &
+                                                                    .OR. ln_rough .OR. ln_zdfqiao) )   THEN  
+             CALL ctl_warn( 'Ask for wave coupling but ln_cdgw=F, ln_sdw=F, ln_tauoc=F, ln_tauw=F, ln_stcor=F ', &
+                                                      'ln_phioc=F, ln_rough=F, ln_zdfqiao=F' ) 
+         ELSE
+         IF (ln_stcor .AND. .NOT. ln_sdw) & 
              CALL ctl_stop( 'Stokes-Coriolis term calculated only if activated Stokes Drift ln_sdw=T')
-         ENDIF
          IF ( ln_cdgw .AND. .NOT.(nn_drag==jp_ukmo .OR. nn_drag==jp_std .OR. nn_drag==jp_const .OR. nn_drag==jp_mcore) ) & 
              CALL ctl_stop( 'The chosen nn_drag for momentum calculation must be 0, 1, 2, or 3')
@@ -398,6 +454,15 @@
          IF( ln_zdfqiao .AND. .NOT.ln_sdw ) & 
             CALL ctl_stop( 'Qiao vertical mixing can not be used without Stokes drift (ln_sdw)' )
+         IF( ln_tauoc .AND. ln_tauw ) & 
+            CALL ctl_stop( 'More than one method for modifying the ocean stress has been selected ', &
+                                     '(ln_tauoc=.true. and ln_tauw=.true.)' )
+         IF( ln_tauoc ) &
+             CALL ctl_warn( 'You are subtracting the wave stress to the ocean (ln_tauoc=.true.)' )
+         IF( ln_tauw ) &
+             CALL ctl_warn( 'The wave modified ocean stress components are used (ln_tauw=.true.) ', &
+                                  'This will override any other specification of the ocean stress' ) 
+         ENDIF
       ELSE
-         IF ( ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_stcor .OR. ln_phioc .OR. ln_rough .OR. ln_zdfqiao ) &  
+         IF ( ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_tauw .OR. ln_stcor .OR. ln_phioc .OR. ln_rough .OR. ln_zdfqiao ) &  
             &   CALL ctl_stop( 'Not Activated Wave Module (ln_wave=F) but asked coupling ',    & 
             &                  'with drag coefficient (ln_cdgw =T) '  ,                        & 
@@ -405,4 +470,5 @@
             &                  'or Stokes-Coriolis term (ln_stcor=T)',                         &
             &                  'or ocean stress modification due to waves (ln_tauoc=T) ',      & 
+            &                  'or ocean stress components from waves (ln_tauw=T) ',          & 
             &                  'or wave to ocean energy modification (ln_phioc=T) ',           & 
             &                  'or wave surface roughness (ln_rough=T) ',                      & 
@@ -432,4 +498,22 @@
          ENDIF
          ALLOCATE( tauoc_wave(jpi,jpj) )
+      ENDIF
+
+      IF( ln_tauw ) THEN
+         IF( .NOT. cpl_tauw ) THEN
+            ALLOCATE( sf_tauw(2), STAT=ierror )           !* allocate and fill sf_wave with sn_tauwx/y
+            IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave_init: unable to allocate sf_tauw structure' )
+            !
+            ALLOCATE( slf_j(2) )
+            slf_j(1) = sn_tauwx
+            slf_j(2) = sn_tauwy
+                                    ALLOCATE( sf_tauw(1)%fnow(jpi,jpj,1)   )
+                                    ALLOCATE( sf_tauw(2)%fnow(jpi,jpj,1)   )
+            IF( slf_j(1)%ln_tint )  ALLOCATE( sf_tauw(1)%fdta(jpi,jpj,1,2) )
+            IF( slf_j(2)%ln_tint )  ALLOCATE( sf_tauw(2)%fdta(jpi,jpj,1,2) )
+            CALL fld_fill( sf_tauw, (/ slf_j /), cn_dir, 'sbc_wave_init', 'read wave input', 'namsbc_wave' )
+         ENDIF
+         ALLOCATE( tauw_x(jpi,jpj) )
+         ALLOCATE( tauw_y(jpi,jpj) )
       ENDIF
 
