MODULE dynstcor
   !!======================================================================
   !!                       ***  MODULE  dynstcor  ***
   !! Ocean dynamics: Stokes-Coriolis effect.
   !!
   !!======================================================================
   !! History :  0.1  !  2012-10  (Oyvind Breivik)
   !!----------------------------------------------------------------------
   !! dyn_stcor        : Add the Stokes-Coriolis forcing to the momentum equation
   !!----------------------------------------------------------------------
   USE oce              ! ocean dynamics and tracers
   USE dom_oce          ! ocean space and time domain
   USE in_out_manager   ! I/O manager
   USE lib_mpp          ! distributed memory computing
   USE prtctl           ! Print control
   USE phycst
   USE lbclnk
   USE wrk_nemo         ! Memory Allocation
   USE sbcmod           ! Access to ln_stcor (sbc_oce) and wave parameters (sbc_wave)
   USE sbcwave          ! Wave module
   USE timing          ! Timing

   IMPLICIT NONE
   PRIVATE

   !! * Routine accessibility
   PUBLIC dyn_stcor      ! routine called by step.F90

   !! * Substitutions
#  include "vectopt_loop_substitute.h90"
   !!----------------------------------------------------------------------
   !! NEMO/OPA 3.6 , NEMO Consortium (2010)
   !! $Id:$
   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
   !!----------------------------------------------------------------------
 CONTAINS

   SUBROUTINE dyn_stcor( kt )
      !!----------------------------------------------------------------------
      !!                  ***  ROUTINE dyn_stcor ***
      !!
      !! ** Purpose:  Add Stokes-Coriolis forcing to horizontal momentum equation.
      !!
      !! ** History:  3.4  !   2012-10  (O. Breivik)   Initial version
      !!              3.6  !   2014-10  (R. Benshila) 
      !!----------------------------------------------------------------------
      INTEGER, INTENT( in ) ::   kt  ! ocean time-step index
      !!
      INTEGER  ::  ji, jj, jk        ! dummy loop indices
      REAL(wp) :: zx1, zx2, zy1, zy2
      !!----------------------------------------------------------------------
      !
      IF( nn_timing == 1 )  CALL timing_start('dyn_stcor')
      !
      IF( kt == nit000 ) THEN
         IF(lwp) WRITE(numout,*)
         IF(lwp) WRITE(numout,*) 'dyn_stcor : time stepping'
         IF(lwp) WRITE(numout,*) '~~~~~~~'
      ENDIF

      !
      ! Update velocity tendencies ua, va by adding the Stokes-Coriolis velocities 
      !
      DO jk = 1, jpkm1
         DO jj = 2, jpjm1
            DO ji = 2, jpim1
               zy1 = ff(ji  ,jj-1) * ( vsd3d(ji  ,jj-1,jk) + vsd3d(ji+1,jj-1,jk) )
               zy2 = ff(ji  ,jj  ) * ( vsd3d(ji  ,jj  ,jk) + vsd3d(ji+1,jj  ,jk) )
               zx1 = ff(ji-1,jj  ) * ( usd3d(ji-1,jj  ,jk) + usd3d(ji-1,jj+1,jk) )
               zx2 = ff(ji  ,jj  ) * ( usd3d(ji  ,jj  ,jk) + usd3d(ji  ,jj+1,jk) )

               ua(ji,jj,jk) = ua(ji,jj,jk) + 0.25 * ( zy1 + zy2 ) *umask(ji,jj,jk)
               va(ji,jj,jk) = va(ji,jj,jk) - 0.25 * ( zx1 + zx2 ) *vmask(ji,jj,jk)
            ENDDO
         ENDDO
      ENDDO
      !
      IF( nn_timing == 1 )  CALL timing_stop('dyn_stcor')
      !
   END SUBROUTINE dyn_stcor

END MODULE dynstcor