Index: NEMO/trunk/src/OCE/BDY/bdydyn3d.F90
===================================================================
--- NEMO/trunk/src/OCE/BDY/bdydyn3d.F90	(revision 10529)
+++ NEMO/trunk/src/OCE/BDY/bdydyn3d.F90	(revision 11536)
@@ -42,21 +42,65 @@
       INTEGER, INTENT(in) ::   kt   ! Main time step counter
       !
-      INTEGER ::   ib_bdy   ! loop index
-      !!----------------------------------------------------------------------
-      !
-      DO ib_bdy=1, nb_bdy
+      INTEGER  ::   ib_bdy, ir     ! BDY set index, rim index
+      LOGICAL  ::   llrim0         ! indicate if rim 0 is treated
+      LOGICAL, DIMENSION(4) :: llsend2, llrecv2, llsend3, llrecv3  ! indicate how communications are to be carried out
+
+      !!----------------------------------------------------------------------
+      llsend2(:) = .false.   ;   llrecv2(:) = .false.
+      llsend3(:) = .false.   ;   llrecv3(:) = .false.
+      DO ir = 1, 0, -1   ! treat rim 1 before rim 0
+         IF( ir == 0 ) THEN   ;   llrim0 = .TRUE.
+         ELSE                 ;   llrim0 = .FALSE.
+         END IF
+         DO ib_bdy=1, nb_bdy
+            !
+            SELECT CASE( cn_dyn3d(ib_bdy) )
+            CASE('none')        ;   CYCLE
+            CASE('frs' )        ! treat the whole boundary at once
+               IF( ir == 0) CALL bdy_dyn3d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
+            CASE('specified')   ! treat the whole rim      at once
+               IF( ir == 0) CALL bdy_dyn3d_spe( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
+            CASE('zero')        ! treat the whole rim      at once
+               IF( ir == 0) CALL bdy_dyn3d_zro( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
+            CASE('orlanski' )   ;   CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, llrim0, ll_npo=.false. )
+            CASE('orlanski_npo');   CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, llrim0, ll_npo=.true.  )
+            CASE('zerograd')    ;   CALL bdy_dyn3d_zgrad( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy, llrim0 )
+            CASE('neumann')     ;   CALL bdy_dyn3d_nmn( idx_bdy(ib_bdy), ib_bdy, llrim0 )
+            CASE DEFAULT        ;   CALL ctl_stop( 'bdy_dyn3d : unrecognised option for open boundaries for baroclinic velocities' )
+            END SELECT
+         END DO
          !
-         SELECT CASE( cn_dyn3d(ib_bdy) )
-         CASE('none')        ;   CYCLE
-         CASE('frs' )        ;   CALL bdy_dyn3d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
-         CASE('specified')   ;   CALL bdy_dyn3d_spe( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
-         CASE('zero')        ;   CALL bdy_dyn3d_zro( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
-         CASE('orlanski' )   ;   CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, ll_npo=.false. )
-         CASE('orlanski_npo');   CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, ll_npo=.true. )
-         CASE('zerograd')    ;   CALL bdy_dyn3d_zgrad( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy )
-         CASE('neumann')     ;   CALL bdy_dyn3d_nmn( idx_bdy(ib_bdy), ib_bdy )
-         CASE DEFAULT        ;   CALL ctl_stop( 'bdy_dyn3d : unrecognised option for open boundaries for baroclinic velocities' )
-         END SELECT
-      END DO
+         IF( nn_hls > 1 .AND. ir == 1 ) CYCLE   ! at least 2 halos will be corrected -> no need to correct rim 1 before rim 0
+         IF( nn_hls == 1 ) THEN
+            llsend2(:) = .false.   ;   llrecv2(:) = .false.
+            llsend3(:) = .false.   ;   llrecv3(:) = .false.
+         END IF
+         DO ib_bdy=1, nb_bdy
+            SELECT CASE( cn_dyn3d(ib_bdy) )
+            CASE('orlanski', 'orlanski_npo')
+               llsend2(:) = llsend2(:) .OR. lsend_bdy(ib_bdy,2,:,ir)   ! possibly every direction, U points
+               llrecv2(:) = llrecv2(:) .OR. lrecv_bdy(ib_bdy,2,:,ir)   ! possibly every direction, U points
+               llsend3(:) = llsend3(:) .OR. lsend_bdy(ib_bdy,3,:,ir)   ! possibly every direction, V points
+               llrecv3(:) = llrecv3(:) .OR. lrecv_bdy(ib_bdy,3,:,ir)   ! possibly every direction, V points
+            CASE('zerograd')
+               llsend2(3:4) = llsend2(3:4) .OR. lsend_bdyint(ib_bdy,2,3:4,ir)   ! north/south, U points
+               llrecv2(3:4) = llrecv2(3:4) .OR. lrecv_bdyint(ib_bdy,2,3:4,ir)   ! north/south, U points
+               llsend3(1:2) = llsend3(1:2) .OR. lsend_bdyint(ib_bdy,3,1:2,ir)   ! west/east, V points
+               llrecv3(1:2) = llrecv3(1:2) .OR. lrecv_bdyint(ib_bdy,3,1:2,ir)   ! west/east, V points
+            CASE('neumann')
+               llsend2(:) = llsend2(:) .OR. lsend_bdyint(ib_bdy,2,:,ir)   ! possibly every direction, U points
+               llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(ib_bdy,2,:,ir)   ! possibly every direction, U points
+               llsend3(:) = llsend3(:) .OR. lsend_bdyint(ib_bdy,3,:,ir)   ! possibly every direction, V points
+               llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(ib_bdy,3,:,ir)   ! possibly every direction, V points
+            END SELECT
+         END DO
+         !
+         IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN   ! if need to send/recv in at least one direction
+            CALL lbc_lnk( 'bdydyn2d', ua, 'U', -1., kfillmode=jpfillnothing ,lsend=llsend2, lrecv=llrecv2 )
+         END IF
+         IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction
+            CALL lbc_lnk( 'bdydyn2d', va, 'V', -1., kfillmode=jpfillnothing ,lsend=llsend3, lrecv=llrecv3 )
+         END IF
+      END DO   ! ir
       !
    END SUBROUTINE bdy_dyn3d
@@ -78,5 +122,4 @@
       INTEGER  ::   jb, jk         ! dummy loop indices
       INTEGER  ::   ii, ij, igrd   ! local integers
-      REAL(wp) ::   zwgt           ! boundary weight
       !!----------------------------------------------------------------------
       !
@@ -98,13 +141,9 @@
          END DO
       END DO
-      CALL lbc_bdy_lnk( 'bdydyn3d', ua, 'U', -1., ib_bdy )   ! Boundary points should be updated  
-      CALL lbc_bdy_lnk( 'bdydyn3d', va, 'V', -1., ib_bdy )   
-      !
-      IF( kt == nit000 )   CLOSE( unit = 102 )
       !
    END SUBROUTINE bdy_dyn3d_spe
 
 
-   SUBROUTINE bdy_dyn3d_zgrad( idx, dta, kt , ib_bdy )
+   SUBROUTINE bdy_dyn3d_zgrad( idx, dta, kt, ib_bdy, llrim0 )
       !!----------------------------------------------------------------------
       !!                  ***  SUBROUTINE bdy_dyn3d_zgrad  ***
@@ -114,39 +153,54 @@
       !!----------------------------------------------------------------------
       INTEGER                     ::   kt
-      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
-      TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data
-      INTEGER,         INTENT(in) ::   ib_bdy  ! BDY set index
+      TYPE(OBC_INDEX), INTENT(in) ::   idx      ! OBC indices
+      TYPE(OBC_DATA),  INTENT(in) ::   dta      ! OBC external data
+      INTEGER,         INTENT(in) ::   ib_bdy   ! BDY set index
+      LOGICAL,         INTENT(in) ::   llrim0   ! indicate if rim 0 is treated
       !!
       INTEGER  ::   jb, jk         ! dummy loop indices
       INTEGER  ::   ii, ij, igrd   ! local integers
-      REAL(wp) ::   zwgt           ! boundary weight
-      INTEGER  ::   fu, fv
+      INTEGER  ::   flagu, flagv           ! short cuts
+      INTEGER  ::   ibeg, iend     ! length of rim to be treated (rim 0 or rim 1 or both)
       !!----------------------------------------------------------------------
       !
       igrd = 2                      ! Copying tangential velocity into bdy points
-      DO jb = 1, idx%nblenrim(igrd)
-         DO jk = 1, jpkm1
-            ii   = idx%nbi(jb,igrd)
-            ij   = idx%nbj(jb,igrd)
-            fu   = ABS( ABS (NINT( idx%flagu(jb,igrd) ) ) - 1 )
-            ua(ii,ij,jk) = ua(ii,ij,jk) * REAL( 1 - fu) + ( ua(ii,ij+fu,jk) * umask(ii,ij+fu,jk) &
-                        &+ ua(ii,ij-fu,jk) * umask(ii,ij-fu,jk) ) * umask(ii,ij,jk) * REAL( fu )
-         END DO
+      IF( llrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd)
+      ELSE                ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd)
+      ENDIF
+      DO jb = ibeg, iend
+         ii    = idx%nbi(jb,igrd)
+         ij    = idx%nbj(jb,igrd)
+         flagu = NINT(idx%flagu(jb,igrd))
+         flagv = NINT(idx%flagv(jb,igrd))
+         !
+         IF( flagu == 0 )   THEN              ! north/south bdy
+            IF( ij+flagv > jpj .OR. ij+flagv < 1 )   CYCLE      
+            !
+            DO jk = 1, jpkm1
+               ua(ii,ij,jk) = ua(ii,ij+flagv,jk) * umask(ii,ij+flagv,jk)
+            END DO
+            !
+         END IF
       END DO
       !
       igrd = 3                      ! Copying tangential velocity into bdy points
-      DO jb = 1, idx%nblenrim(igrd)
-         DO jk = 1, jpkm1
-            ii   = idx%nbi(jb,igrd)
-            ij   = idx%nbj(jb,igrd)
-            fv   = ABS( ABS (NINT( idx%flagv(jb,igrd) ) ) - 1 )
-            va(ii,ij,jk) = va(ii,ij,jk) * REAL( 1 - fv ) + ( va(ii+fv,ij,jk) * vmask(ii+fv,ij,jk) &
-                        &+ va(ii-fv,ij,jk) * vmask(ii-fv,ij,jk) ) * vmask(ii,ij,jk) * REAL( fv )
-         END DO
-      END DO
-      CALL lbc_bdy_lnk( 'bdydyn3d', ua, 'U', -1., ib_bdy )   ! Boundary points should be updated  
-      CALL lbc_bdy_lnk( 'bdydyn3d', va, 'V', -1., ib_bdy )   
-      !
-      IF( kt == nit000 )   CLOSE( unit = 102 )
+      IF( llrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd)
+      ELSE                ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd)
+      ENDIF
+      DO jb = ibeg, iend
+         ii    = idx%nbi(jb,igrd)
+         ij    = idx%nbj(jb,igrd)
+         flagu = NINT(idx%flagu(jb,igrd))
+         flagv = NINT(idx%flagv(jb,igrd))
+         !
+         IF( flagv == 0 )   THEN              !  west/east  bdy
+            IF( ii+flagu > jpi .OR. ii+flagu < 1 )   CYCLE      
+            !
+            DO jk = 1, jpkm1
+               va(ii,ij,jk) = va(ii+flagu,ij,jk) * vmask(ii+flagu,ij,jk)
+            END DO
+            !
+         END IF
+      END DO
       !
    END SUBROUTINE bdy_dyn3d_zgrad
@@ -167,5 +221,4 @@
       INTEGER  ::   ib, ik         ! dummy loop indices
       INTEGER  ::   ii, ij, igrd   ! local integers
-      REAL(wp) ::   zwgt           ! boundary weight
       !!----------------------------------------------------------------------
       !
@@ -178,5 +231,5 @@
          END DO
       END DO
-
+      !
       igrd = 3                       ! Everything is at T-points here
       DO ib = 1, idx%nblenrim(igrd)
@@ -187,8 +240,4 @@
          END DO
       END DO
-      !
-      CALL lbc_bdy_lnk( 'bdydyn3d', ua, 'U', -1., ib_bdy )   ;   CALL lbc_bdy_lnk( 'bdydyn3d', va, 'V', -1.,ib_bdy )   ! Boundary points should be updated
-      !
-      IF( kt == nit000 )   CLOSE( unit = 102 )
       !
    END SUBROUTINE bdy_dyn3d_zro
@@ -234,14 +283,10 @@
             va(ii,ij,jk) = ( va(ii,ij,jk) + zwgt * ( dta%v3d(jb,jk) - va(ii,ij,jk) ) ) * vmask(ii,ij,jk)
          END DO
-      END DO 
-      CALL lbc_bdy_lnk( 'bdydyn3d', ua, 'U', -1., ib_bdy )    ! Boundary points should be updated
-      CALL lbc_bdy_lnk( 'bdydyn3d', va, 'V', -1., ib_bdy )   
-      !
-      IF( kt == nit000 )   CLOSE( unit = 102 )
+      END DO   
       !
    END SUBROUTINE bdy_dyn3d_frs
 
 
-   SUBROUTINE bdy_dyn3d_orlanski( idx, dta, ib_bdy, ll_npo )
+   SUBROUTINE bdy_dyn3d_orlanski( idx, dta, ib_bdy, llrim0, ll_npo )
       !!----------------------------------------------------------------------
       !!                 ***  SUBROUTINE bdy_dyn3d_orlanski  ***
@@ -255,6 +300,7 @@
       TYPE(OBC_INDEX),              INTENT(in) ::   idx  ! OBC indices
       TYPE(OBC_DATA),               INTENT(in) ::   dta  ! OBC external data
-      INTEGER,                      INTENT(in) ::   ib_bdy  ! BDY set index
-      LOGICAL,                      INTENT(in) ::   ll_npo  ! switch for NPO version
+      INTEGER,                      INTENT(in) ::   ib_bdy   ! BDY set index
+      LOGICAL,                      INTENT(in) ::   llrim0   ! indicate if rim 0 is treated
+      LOGICAL,                      INTENT(in) ::   ll_npo   ! switch for NPO version
 
       INTEGER  ::   jb, igrd                               ! dummy loop indices
@@ -265,12 +311,9 @@
       igrd = 2      ! Orlanski bc on u-velocity; 
       !            
-      CALL bdy_orlanski_3d( idx, igrd, ub, ua, dta%u3d, ll_npo )
+      CALL bdy_orlanski_3d( idx, igrd, ub, ua, dta%u3d, ll_npo, llrim0 )
 
       igrd = 3      ! Orlanski bc on v-velocity
       !  
-      CALL bdy_orlanski_3d( idx, igrd, vb, va, dta%v3d, ll_npo )
-      !
-      CALL lbc_bdy_lnk( 'bdydyn3d', ua, 'U', -1., ib_bdy )    ! Boundary points should be updated
-      CALL lbc_bdy_lnk( 'bdydyn3d', va, 'V', -1., ib_bdy )   
+      CALL bdy_orlanski_3d( idx, igrd, vb, va, dta%v3d, ll_npo, llrim0 )
       !
    END SUBROUTINE bdy_dyn3d_orlanski
@@ -320,6 +363,4 @@
       END DO
       !
-      CALL lbc_lnk_multi( 'bdydyn3d', ua, 'U', -1.,  va, 'V', -1. )   ! Boundary points should be updated
-      !
       IF( ln_timing )   CALL timing_stop('bdy_dyn3d_dmp')
       !
@@ -327,5 +368,5 @@
 
 
-   SUBROUTINE bdy_dyn3d_nmn( idx, ib_bdy )
+   SUBROUTINE bdy_dyn3d_nmn( idx, ib_bdy, llrim0 )
       !!----------------------------------------------------------------------
       !!                 ***  SUBROUTINE bdy_dyn3d_nmn  ***
@@ -336,8 +377,8 @@
       !!
       !!----------------------------------------------------------------------
-      TYPE(OBC_INDEX),              INTENT(in) ::   idx  ! OBC indices
-      INTEGER,                      INTENT(in) ::   ib_bdy  ! BDY set index
-
-      INTEGER  ::   jb, igrd                               ! dummy loop indices
+      TYPE(OBC_INDEX), INTENT(in) ::   idx      ! OBC indices
+      INTEGER,         INTENT(in) ::   ib_bdy   ! BDY set index
+      LOGICAL,         INTENT(in) ::   llrim0   ! indicate if rim 0 is treated
+      INTEGER  ::   igrd                        ! dummy indice
       !!----------------------------------------------------------------------
       !
@@ -346,12 +387,9 @@
       igrd = 2      ! Neumann bc on u-velocity; 
       !            
-      CALL bdy_nmn( idx, igrd, ua )
+      CALL bdy_nmn( idx, igrd, ua, llrim0 )   ! ua is masked
 
       igrd = 3      ! Neumann bc on v-velocity
       !  
-      CALL bdy_nmn( idx, igrd, va )
-      !
-      CALL lbc_bdy_lnk( 'bdydyn3d', ua, 'U', -1., ib_bdy )    ! Boundary points should be updated
-      CALL lbc_bdy_lnk( 'bdydyn3d', va, 'V', -1., ib_bdy )
+      CALL bdy_nmn( idx, igrd, va, llrim0 )   ! va is masked
       !
    END SUBROUTINE bdy_dyn3d_nmn
