  MODULE traobc2dapply

   !!======================================================================
   !!                   ***  MODULE  traobc2dapply  ***
   !! Ocean dynamics : apply2d open boundary conditions of variables
   !!                  temperature and salinity
   !!======================================================================

#include "toolcpp.h"

   !! * Modules used

   IMPLICIT NONE
   PRIVATE

   !! * Accessibility
   PUBLIC ts_obc2d_apply                ! routine called by stepu.F90 stepv.F90

   !! * Shared module variables

   !! * Private variables

  CONTAINS
   !!======================================================================

  SUBROUTINE ts_obc2d_apply(u,v,sal,temp)

   !&E---------------------------------------------------------------------
   !&E                 ***  ROUTINE ts_obc2d_apply  ***
   !&E
   !&E ** Purpose : apply2d open boundary conditions for temperature
   !&E              and salinity field
   !&E
   !&E ** Description : apply2d external value when water mass enters the domain
   !&E
   !&E ** Called by : stepu, stepv
   !&E
   !&E ** External calls :
   !&E
   !&E ** Used ij-arrays : sal,temp,hx,u,obc_sw,obc_tw,dx_u,obc_se,obc_te,v,
   !&E                     obc_ss,obc_ts,dy_v,obc_sn,obc_tn
   !&E
   !&E ** Modified variables : sal temp
   !&E
   !&E ** Reference :
   !&E
   !&E ** History :
   !&E       !                        Original code
   !&E       !  2007-05  (V. Garnier) Mise en forme - introduction de obc_* variables
   !&E       !  2007-09  (V. Garnier) homogeneneity of dx dy...
   !&E       !  2009-04  (A. Mano) MPI2D
   !&E       !  2010-04  (T. Odaka) MPI2D cyclic
   !&E
   !&E---------------------------------------------------------------------
   !! * Modules used
   USE parameters, ONLY : imin,imax,jmin,jmax,kmax,                      &
                          ljminp2,ljmaxm1,ljmin,ljmax,ljminm1,ljmaxp1,   &
                          limin,limax,liminm1,limaxp2,liminp2,limaxm1,   &
                          ljmaxp2,liminm2,ljminm2,rsh,rlg,               &
                          ishift_obc_west,ishift_obc_east,jshift_obc_south,jshift_obc_north
   USE comvars2d,  ONLY : dt,ig,hx,h0fond,dx_u,id,jb,hy,dy_v,jh,l_obc_cycl_x,l_obc_cycl_y
   USE comvars3d,  ONLY : l_testcase,valmanq
   USE obcogcm,    ONLY : obc_tw,obc_se,obc_te,obc_ss,obc_ts,obc_sn,obc_tn,obc_sw,obc_coefrel
   USE obccombine, ONLY : l_obc_south,l_obc_north,l_obc_east,l_obc_west
   USE_MPI toolmpi,ONLY : ex_i_rsh,ex_j_rsh,ex_peri_x_rsh,ex_peri_y_rsh

   !! * Arguments
   REAL(KIND=rsh),DIMENSION(kmax,liminm1:limaxp2,ljminm1:ljmaxp2),INTENT(inout) :: sal,temp
   REAL(KIND=rsh),DIMENSION(liminm2:limaxp2,ljminm2:ljmaxp2),INTENT(in) :: u,v

   !! * Local declarations
   INTEGER                   :: i,j,k
   REAL(KIND=rsh)            :: uu,vv,cou,dtrelax
   REAL(KIND=rlg)            :: dtm

   !!----------------------------------------------------------------------
   !! * Executable part

   PRINT_DBG*, 'enter TS_OBC2D_APPLY routine'

   IF (.NOT.l_testcase) THEN

     dtm=0.5_rlg*dt
     dtrelax=dtm*obc_coefrel  ! dt / time scale of relaxation
     IF (obc_coefrel==0.0_rsh) dtrelax=1.0_rsh

!$OMP SINGLE
   ! periodic conditions along x-axis
     IF (l_obc_cycl_x) THEN

#ifdef key_MPI_2D
!!!!!!#if defined key_tssub_adv_fourthorder  || defined key_tssub_adv_fifthorder
#ifdef key_shiftobc
       CALL_MPI ex_peri_x_rsh(3,1,liminm1,limaxp2,ljminm1,ljmaxp2,sal(1,liminm1:limaxp2,ljminm1:ljmaxp2))
       CALL_MPI ex_peri_x_rsh(3,1,liminm1,limaxp2,ljminm1,ljmaxp2,temp(1,liminm1:limaxp2,ljminm1:ljmaxp2))
#else
       CALL_MPI ex_peri_x_rsh(1,1,liminm1,limaxp2,ljminm1,ljmaxp2,sal(1,liminm1:limaxp2,ljminm1:ljmaxp2))
       CALL_MPI ex_peri_x_rsh(1,1,liminm1,limaxp2,ljminm1,ljmaxp2,temp(1,liminm1:limaxp2,ljminm1:ljmaxp2))
#endif
#else
       DO j=ljmin,ljmax
   ! western boundary
         sal(1,ig(j),j) = sal(1,id(j)-1,j)
         sal(1,ig(j)-1,j) = sal(1,id(j)-2,j)
         temp(1,ig(j),j) = temp(1,id(j)-1,j)
         temp(1,ig(j)-1,j) = temp(1,id(j)-2,j)
!!!!!!#if defined key_tssub_adv_fourthorder  || defined key_tssub_adv_fifthorder
#ifdef key_shiftobc
         sal(1,ig(j)-2,j) = sal(1,id(j)-3,j)
         temp(1,ig(j)-2,j) = temp(1,id(j)-3,j)
#endif
   ! eastern boundary
         sal(1,id(j),j) = sal(1,ig(j)+1,j)
         temp(1,id(j),j) = temp(1,ig(j)+1,j)
!!!!#if defined key_tssub_adv_fourthorder  || defined key_tssub_adv_fifthorder
#ifdef key_shiftobc
         sal(1,id(j)+1,j) = sal(1,ig(j)+2,j)
         sal(1,id(j)+2,j) = sal(1,ig(j)+3,j)
         temp(1,id(j)+1,j) = temp(1,ig(j)+2,j)
         temp(1,id(j)+2,j) = temp(1,ig(j)+3,j)
#endif
       END DO
#endif
     END IF  ! ends periodic conditions along x-axis

   ! periodic conditions along y-axis
     IF (l_obc_cycl_y) THEN

#ifdef key_MPI_2D
!!!!!#if defined key_tssub_adv_fourthorder  || defined key_tssub_adv_fifthorder
#ifdef key_shiftobc
       CALL_MPI ex_peri_y_rsh(3,1,liminm1,limaxp2,ljminm1,ljmaxp2,sal(1,liminm1:limaxp2,ljminm1:ljmaxp2))
       CALL_MPI ex_peri_y_rsh(3,1,liminm1,limaxp2,ljminm1,ljmaxp2,temp(1,liminm1:limaxp2,ljminm1:ljmaxp2))
#else
       CALL_MPI ex_peri_y_rsh(1,1,liminm1,limaxp2,ljminm1,ljmaxp2,sal(1,liminm1:limaxp2,ljminm1:ljmaxp2))
       CALL_MPI ex_peri_y_rsh(1,1,liminm1,limaxp2,ljminm1,ljmaxp2,temp(1,liminm1:limaxp2,ljminm1:ljmaxp2))
#endif
#else
       DO i=limin,limax
   ! southern boundary
         sal(1,i,jb(i)) = sal(1,i,jh(i)-1)
         sal(1,i,jb(i)-1) = sal(1,i,jh(i)-2)
         temp(1,i,jb(i)) = temp(1,i,jh(i)-1)
         temp(1,i,jb(i)-1) = temp(1,i,jh(i)-2)
!!!#if defined key_tssub_adv_fourthorder  || defined key_tssub_adv_fifthorder
#ifdef key_shiftobc
         sal(1,i,jb(i)-2) = sal(1,i,jh(i)-3)
         temp(1,i,jb(i)-2) = temp(1,i,jh(i)-3)
#endif
   ! northern boundary
         sal(1,i,jh(i)) = sal(1,i,jb(i)+1)
         temp(1,i,jh(i)) = temp(1,i,jb(i)+1)
!!!!#if defined key_tssub_adv_fourthorder  || defined key_tssub_adv_fifthorder
#ifdef key_shiftobc
         sal(1,i,jh(i)+1) = sal(1,i,jb(i)+2)
         sal(1,i,jh(i)+2) = sal(1,i,jb(i)+3)
         temp(1,i,jh(i)+1) = temp(1,i,jb(i)+2)
         temp(1,i,jh(i)+2) = temp(1,i,jb(i)+3)
#endif
       END DO
#endif

     END IF  ! ends periodic conditions along y-axis
!$OMP END SINGLE


     IF (.NOT.l_obc_cycl_x .AND. (l_obc_west .OR. l_obc_east)) THEN

       DO j=ljminp2,ljmaxm1

   ! Western boundary
         i=ig(j)
         IF_MPI (i==imin+1+ishift_obc_west .AND. limin<=i .AND. i<=limax) THEN
         IF (hx(i,j) /= -valmanq-h0fond) THEN
           uu=u(i,j)
           IF (uu > 0.0_rsh) THEN
             sal(1,i,j)=sal(1,i,j)*(1.0_rsh-dtrelax)+obc_sw(1,1,j)*dtrelax
             temp(1,i,j)=temp(1,i,j)*(1.0_rsh-dtrelax)+obc_tw(1,1,j)*dtrelax
           ELSE
             cou=uu*dtm/dx_u(i,j)
             sal(1,i,j)=sal(1,i,j)*(1.0_rsh+cou)-cou*sal(1,i+1,j)
             temp(1,i,j)=temp(1,i,j)*(1.0_rsh+cou)-cou*temp(1,i+1,j)
           END IF
         END IF
         ENDIF_MPI

   ! Eastern boundary
         i=id(j)
         IF_MPI (i==imax-ishift_obc_east .AND. limin<=i .AND. i<=limax) THEN
         IF (hx(i-1,j) /= -valmanq-h0fond) THEN
           uu=u(i-1,j)
           IF (uu < 0.0_rsh) THEN
             sal(1,i,j)=sal(1,i,j)*(1.0_rsh-dtrelax)+obc_se(1,imax,j)*dtrelax
             temp(1,i,j)=temp(1,i,j)*(1.0_rsh-dtrelax)+obc_te(1,imax,j)*dtrelax
           ELSE
             cou=uu*dtm/dx_u(i-1,j)
             sal(1,i,j)=sal(1,i,j)*(1.0_rsh-cou)+cou*sal(1,i-1,j)
             temp(1,i,j)=temp(1,i,j)*(1.0_rsh-cou)+cou*temp(1,i-1,j)
           END IF
         END IF
         ENDIF_MPI

       END DO

     END IF  ! ends test on .not.l_obc_cycl_x

     IF (.NOT.l_obc_cycl_y .AND. (l_obc_north .OR. l_obc_south)) THEN

       DO j=ljmin,ljmax
         DO i=liminp2,limaxm1

   ! Southern boundary
           IF (j == jb(i) .AND. j==jmin+1+jshift_obc_south) THEN
             IF (hy(i,j) /= -valmanq-h0fond) THEN
               vv=v(i,j)
               IF (vv > 0.0_rsh) THEN
                 sal(1,i,j)=sal(1,i,j)*(1.0_rsh-dtrelax)+obc_ss(1,i,1)*dtrelax
                 temp(1,i,j)=temp(1,i,j)*(1.0_rsh-dtrelax)+obc_ts(1,i,1)*dtrelax
               ELSE
                 cou=vv*dtm/dy_v(i,j)
                 sal(1,i,j)=sal(1,i,j)*(1.0_rsh+cou)-cou*sal(1,i,j+1)
                 temp(1,i,j)=temp(1,i,j)*(1.0_rsh+cou)-cou*temp(1,i,j+1)
               END IF
             END IF

   ! Northern boundary
           ELSE IF (j == jh(i) .AND. j==jmax-jshift_obc_north) THEN
             IF (hy(i,j-1) /= -valmanq-h0fond) THEN
               vv=v(i,j-1)
               IF (vv < 0.0_rsh) THEN
                 sal(1,i,j)=sal(1,i,j)*(1.0_rsh-dtrelax)+obc_sn(1,i,jmax)*dtrelax
                 temp(1,i,j)=temp(1,i,j)*(1.0_rsh-dtrelax)+obc_tn(1,i,jmax)*dtrelax
               ELSE
                 cou=vv*dtm/dy_v(i,j-1)
                 sal(1,i,j)=sal(1,i,j)*(1.0_rsh-cou)+cou*sal(1,i,j-1)
                 temp(1,i,j)=temp(1,i,j)*(1.0_rsh-cou)+cou*temp(1,i,j-1)
               END IF
             END IF
           END IF

         END DO
       END DO

     END IF  ! ends test on .not.l_obc_cycl_y

   END IF   ! ends test on l_testcase

OMPMPI barrier
OMPMPI master
   CALL_MPI ex_i_rsh(-1,2,kmax,liminm1,limaxp2,ljminm1,ljmaxp2,sal(:,liminm1:limaxp2,ljminm1:ljmaxp2))
   CALL_MPI ex_j_rsh(-1,2,kmax,liminm1,limaxp2,ljminm1,ljmaxp2,sal(:,liminm1:limaxp2,ljminm1:ljmaxp2))
   CALL_MPI ex_i_rsh(-1,2,kmax,liminm1,limaxp2,ljminm1,ljmaxp2,temp(:,liminm1:limaxp2,ljminm1:ljmaxp2))
   CALL_MPI ex_j_rsh(-1,2,kmax,liminm1,limaxp2,ljminm1,ljmaxp2,temp(:,liminm1:limaxp2,ljminm1:ljmaxp2))
OMPMPI end master
OMPMPI barrier
OMPMPI flush(sal,temp)

   PRINT_DBG*, 'exit TS_OBC2D_APPLY routine'

  END SUBROUTINE ts_obc2d_apply

   !!======================================================================

  END MODULE traobc2dapply
