Index: /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdy_oce.F90
===================================================================
--- /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdy_oce.F90	(revision 11190)
+++ /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdy_oce.F90	(revision 11191)
@@ -22,4 +22,5 @@
       INTEGER ,          DIMENSION(jpbgrd) ::  nblen
       INTEGER ,          DIMENSION(jpbgrd) ::  nblenrim
+      INTEGER ,          DIMENSION(jpbgrd) ::  nblenrim0
       INTEGER , POINTER, DIMENSION(:,:)    ::  nbi
       INTEGER , POINTER, DIMENSION(:,:)    ::  nbj
@@ -139,10 +140,10 @@
    TYPE(OBC_DATA) , DIMENSION(jp_bdy), TARGET      ::   dta_bdy           !: bdy external data (local process)
 !$AGRIF_END_DO_NOT_TREAT
-   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:) ::   lsend_bdy      !: mark needed communication for given boundary, grid and neighbour
-   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:) ::   lrecv_bdy      !:  when searching in any direction
-   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:) ::   lsend_bdyint   !: mark needed communication for given boundary, grid and neighbour
-   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:) ::   lrecv_bdyint   !:  when searching towards the interior of the computational domain
-   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:) ::   lsend_bdyext   !: mark needed communication for given boundary, grid and neighbour
-   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:) ::   lrecv_bdyext   !:  when searching towards the exterior of the computational domain
+   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) ::   lsend_bdy      !: mark needed communication for given boundary, grid and neighbour
+   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) ::   lrecv_bdy      !:  when searching in any direction
+   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) ::   lsend_bdyint   !: mark needed communication for given boundary, grid and neighbour
+   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) ::   lrecv_bdyint   !:  when searching towards the interior of the computational domain
+   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) ::   lsend_bdyext   !: mark needed communication for given boundary, grid and neighbour
+   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) ::   lrecv_bdyext   !:  when searching towards the exterior of the computational domain
    !!----------------------------------------------------------------------
    !! NEMO/OCE 4.0 , NEMO Consortium (2018)
Index: /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdydyn2d.F90
===================================================================
--- /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdydyn2d.F90	(revision 11190)
+++ /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdydyn2d.F90	(revision 11191)
@@ -50,54 +50,61 @@
       REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: pssh
       !!
-      INTEGER                                  ::   ib_bdy     ! Loop counter
+      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
       
-      DO ib_bdy=1, nb_bdy
-         SELECT CASE( cn_dyn2d(ib_bdy) )
-         CASE('none')
-            CYCLE
-         CASE('frs')
-            CALL bdy_dyn2d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, pua2d, pva2d )
-         CASE('flather')
-            CALL bdy_dyn2d_fla( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, pua2d, pva2d, pssh, phur, phvr )
-         CASE('orlanski')
-            CALL bdy_dyn2d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, &
-                                     & pua2d, pva2d, pub2d, pvb2d, ll_npo=.false.)
-         CASE('orlanski_npo')
-            CALL bdy_dyn2d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, &
-                                     & pua2d, pva2d, pub2d, pvb2d, ll_npo=.true. )
-         CASE DEFAULT
-            CALL ctl_stop( 'bdy_dyn2d : unrecognised option for open boundaries for barotropic variables' )
-         END SELECT
-      ENDDO
-      !
-      llsend2(:) = .false.
-      llrecv2(:) = .false.
-      llsend3(:) = .false.
-      llrecv3(:) = .false.
-      DO ib_bdy=1, nb_bdy
-         SELECT CASE( cn_dyn2d(ib_bdy) )
-         CASE('flather')
-            llsend2(1:2) = llsend2(1:2) .OR. lsend_bdyint(ib_bdy,2,1:2)   ! west/east, U points
-            llsend2(1)   = llsend2(1)   .OR. lsend_bdyext(ib_bdy,2,1)     ! neighbour might search point towards bdy on its east
-            llrecv2(1:2) = llrecv2(1:2) .OR. lrecv_bdyint(ib_bdy,2,1:2)   ! west/east, U points
-            llrecv2(2)   = llrecv2(2)   .OR. lrecv_bdyext(ib_bdy,2,2)     ! might search point towards bdy on the east
-            llsend3(3:4) = llsend3(3:4) .OR. lsend_bdyint(ib_bdy,3,3:4)   ! north/south, V points
-            llsend3(3)   = llsend3(3)   .OR. lsend_bdyext(ib_bdy,3,3)     ! neighbour might search point towards bdy on its north
-            llrecv3(3:4) = llrecv3(3:4) .OR. lrecv_bdyint(ib_bdy,3,3:4)   ! north/south, V points
-            llrecv3(4)   = llrecv3(4)   .OR. lrecv_bdyext(ib_bdy,3,4)     ! might search point towards bdy on the north
-         CASE('orlanski', 'orlanski_npo')
-            llsend2(:) = llsend2(:) .OR. lsend_bdy(ib_bdy,2,:)   ! possibly every direction, U points
-            llrecv2(:) = llrecv2(:) .OR. lrecv_bdy(ib_bdy,2,:)   ! possibly every direction, U points
-            llsend3(:) = llsend3(:) .OR. lsend_bdy(ib_bdy,3,:)   ! possibly every direction, V points
-            llrecv3(:) = llrecv3(:) .OR. lrecv_bdy(ib_bdy,3,:)   ! 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_bdy_lnk( 'bdydyn2d', llsend2, llrecv2, pua2d, 'U', -1. )
-      END IF
-      IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction
-         CALL lbc_bdy_lnk( 'bdydyn2d', llsend3, llrecv3, pva2d, 'V', -1. )
-      END IF
+      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_dyn2d(ib_bdy) )
+            CASE('none')
+               CYCLE
+            CASE('frs')   ! treat the whole boundary at once
+               IF( llrim0 )   CALL bdy_dyn2d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, pua2d, pva2d )
+            CASE('flather')
+               CALL bdy_dyn2d_fla( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, pua2d, pva2d, pssh, phur, phvr, llrim0 )
+            CASE('orlanski')
+               CALL bdy_dyn2d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, &
+                    & pua2d, pva2d, pub2d, pvb2d, llrim0, ll_npo=.false. )
+            CASE('orlanski_npo')
+               CALL bdy_dyn2d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, &
+                    & pua2d, pva2d, pub2d, pvb2d, llrim0, ll_npo=.true.  )
+            CASE DEFAULT
+               CALL ctl_stop( 'bdy_dyn2d : unrecognised option for open boundaries for barotropic variables' )
+            END SELECT
+         ENDDO
+         !
+         llsend2(:) = .false.
+         llrecv2(:) = .false.
+         llsend3(:) = .false.
+         llrecv3(:) = .false.
+         DO ib_bdy=1, nb_bdy
+            SELECT CASE( cn_dyn2d(ib_bdy) )
+            CASE('flather')
+               llsend2(1:2) = llsend2(1:2) .OR. lsend_bdyint(ib_bdy,2,1:2,ir)   ! west/east, U points
+               llsend2(1)   = llsend2(1)   .OR. lsend_bdyext(ib_bdy,2,1,ir)     ! neighbour might search point towards its east bdy
+               llrecv2(1:2) = llrecv2(1:2) .OR. lrecv_bdyint(ib_bdy,2,1:2,ir)   ! west/east, U points
+               llrecv2(2)   = llrecv2(2)   .OR. lrecv_bdyext(ib_bdy,2,2,ir)     ! might search point towards bdy on the east
+               llsend3(3:4) = llsend3(3:4) .OR. lsend_bdyint(ib_bdy,3,3:4,ir)   ! north/south, V points
+               llsend3(3)   = llsend3(3)   .OR. lsend_bdyext(ib_bdy,3,3,ir)     ! neighbour might search point towards its north bdy 
+               llrecv3(3:4) = llrecv3(3:4) .OR. lrecv_bdyint(ib_bdy,3,3:4,ir)   ! north/south, V points
+               llrecv3(4)   = llrecv3(4)   .OR. lrecv_bdyext(ib_bdy,3,4,ir)     ! might search point towards bdy on the north
+            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
+            END SELECT
+         END DO
+         IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN   ! if need to send/recv in at least one direction
+            CALL lbc_bdy_lnk( 'bdydyn2d', llsend2, llrecv2, pua2d, 'U', -1. )
+         END IF
+         IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction
+            CALL lbc_bdy_lnk( 'bdydyn2d', llsend3, llrecv3, pva2d, 'V', -1. )
+         END IF
+         !
+      END DO
       !
    END SUBROUTINE bdy_dyn2d
@@ -119,5 +126,5 @@
       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pua2d, pva2d 
       !!
-      INTEGER  ::   jb, jk         ! dummy loop indices
+      INTEGER  ::   jb             ! dummy loop indices
       INTEGER  ::   ii, ij, igrd   ! local integers
       REAL(wp) ::   zwgt           ! boundary weight
@@ -143,5 +150,5 @@
 
 
-   SUBROUTINE bdy_dyn2d_fla( idx, dta, ib_bdy, pua2d, pva2d, pssh, phur, phvr )
+   SUBROUTINE bdy_dyn2d_fla( idx, dta, ib_bdy, pua2d, pva2d, pssh, phur, phvr, llrim0 )
       !!----------------------------------------------------------------------
       !!                 ***  SUBROUTINE bdy_dyn2d_fla  ***
@@ -166,11 +173,12 @@
       INTEGER,                      INTENT(in) ::   ib_bdy  ! BDY set index
       REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pua2d, pva2d
-      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   pssh, phur, phvr 
-
+      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   pssh, phur, phvr
+      LOGICAL                     , INTENT(in) ::   llrim0   ! indicate if rim 0 is treated
+      INTEGER  ::   ibeg, iend                       ! length of rim to be treated (rim 0 or rim 1)
       INTEGER  ::   jb, igrd                         ! dummy loop indices
       INTEGER  ::   ii, ij                           ! 2D addresses
       INTEGER  ::   iiTrim, ijTrim                   ! T pts i/j-indice on the rim
       INTEGER  ::   iiToce, ijToce, iiUoce, ijVoce   ! T, U and V pts i/j-indice of the ocean next to the rim
-      REAL(wp), POINTER :: flagu, flagv              ! short cuts
+      REAL(wp) ::   flagu, flagv                     ! short cuts
       REAL(wp) ::   zfla                             ! Flather correction
       REAL(wp) ::   z1_2                             ! 
@@ -185,6 +193,9 @@
       ! Fill temporary array with ssh data (here we use spgu with the alias sshdta):
       igrd = 1
+      IF( llrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd)
+      ELSE                ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd)
+      END IF
       sshdta(:,:) = 0.0
-      DO jb = 1, idx%nblenrim(igrd)
+      DO jb = ibeg, iend
          ii = idx%nbi(jb,igrd)
          ij = idx%nbj(jb,igrd)
@@ -196,11 +207,14 @@
       END DO
       !
-      igrd = 2      ! Flather bc on u-velocity; 
+      igrd = 2      ! Flather bc on u-velocity
       !             ! remember that flagu=-1 if normal velocity direction is outward
       !             ! I think we should rather use after ssh ?
-      DO jb = 1, idx%nblenrim(igrd)
-         ii = idx%nbi(jb,igrd)
-         ij = idx%nbj(jb,igrd)
-         flagu => idx%flagu(jb,igrd)
+      IF( llrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd)
+      ELSE                ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd)
+      END IF
+      DO jb = ibeg, iend
+         ii    = idx%nbi(jb,igrd)
+         ij    = idx%nbj(jb,igrd)
+         flagu = idx%flagu(jb,igrd)
          IF( flagu == 0. ) THEN
             pua2d(ii,ij) = dta%u2d(jb)
@@ -222,8 +236,11 @@
       igrd = 3      ! Flather bc on v-velocity
       !             ! remember that flagv=-1 if normal velocity direction is outward
-      DO jb = 1, idx%nblenrim(igrd)
-         ii = idx%nbi(jb,igrd)
-         ij = idx%nbj(jb,igrd)
-         flagv => idx%flagv(jb,igrd)
+      IF( llrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd)
+      ELSE                ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd)
+      END IF
+      DO jb = ibeg, iend
+         ii    = idx%nbi(jb,igrd)
+         ij    = idx%nbj(jb,igrd)
+         flagv = idx%flagv(jb,igrd)
          IF( flagv == 0. ) THEN
             pva2d(ii,ij) = dta%v2d(jb)
@@ -246,5 +263,5 @@
 
 
-   SUBROUTINE bdy_dyn2d_orlanski( idx, dta, ib_bdy, pua2d, pva2d, pub2d, pvb2d, ll_npo )
+   SUBROUTINE bdy_dyn2d_orlanski( idx, dta, ib_bdy, pua2d, pva2d, pub2d, pvb2d, llrim0, ll_npo )
       !!----------------------------------------------------------------------
       !!                 ***  SUBROUTINE bdy_dyn2d_orlanski  ***
@@ -263,5 +280,5 @@
       REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pub2d, pvb2d 
       LOGICAL,                      INTENT(in) ::   ll_npo  ! flag for NPO version
-
+      LOGICAL,                      INTENT(in) ::   llrim0   ! indicate if rim 0 is treated
       INTEGER  ::   ib, igrd                               ! dummy loop indices
       INTEGER  ::   ii, ij, iibm1, ijbm1                   ! indices
@@ -270,9 +287,9 @@
       igrd = 2      ! Orlanski bc on u-velocity; 
       !            
-      CALL bdy_orlanski_2d( idx, igrd, pub2d, pua2d, dta%u2d, ll_npo )
+      CALL bdy_orlanski_2d( idx, igrd, pub2d, pua2d, dta%u2d, llrim0, ll_npo )
 
       igrd = 3      ! Orlanski bc on v-velocity
       !  
-      CALL bdy_orlanski_2d( idx, igrd, pvb2d, pva2d, dta%v2d, ll_npo )
+      CALL bdy_orlanski_2d( idx, igrd, pvb2d, pva2d, dta%v2d, llrim0, ll_npo )
       !
    END SUBROUTINE bdy_dyn2d_orlanski
@@ -288,17 +305,24 @@
       REAL(wp), DIMENSION(jpi,jpj,1), INTENT(inout) ::   zssh ! Sea level, need 3 dimensions to be used by bdy_nmn
       !!
-      INTEGER  ::   ib_bdy          ! bdy index
+      INTEGER ::   ib_bdy, ir      ! bdy index, rim index
+      INTEGER ::   ibeg, iend      ! length of rim to be treated (rim 0 or rim 1)
+      LOGICAL ::   llrim0          ! indicate if rim 0 is treated
       LOGICAL, DIMENSION(4) :: llsend1, llrecv1  ! indicate how communications are to be carried out
       !!----------------------------------------------------------------------
-      llsend1(:) = .false.
-      llrecv1(:) = .false.
-      DO ib_bdy = 1, nb_bdy
-         CALL bdy_nmn( idx_bdy(ib_bdy), 1, zssh )   ! zssh is masked
-         llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:)   ! possibly every direction, T points
-         llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:)   ! possibly every direction, T points
-      END DO
-      IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction
-         CALL lbc_bdy_lnk( 'bdydyn2d', llsend1, llrecv1, zssh(:,:,1), 'T',  1. )
-      END IF
+      DO ir = 1, 0, -1   ! treat rim 1 before rim 0
+         llsend1(:) = .false.
+         llrecv1(:) = .false.
+         IF( ir == 0 ) THEN   ;   llrim0 = .TRUE.
+         ELSE                 ;   llrim0 = .FALSE.
+         END IF
+         DO ib_bdy = 1, nb_bdy
+            CALL bdy_nmn( idx_bdy(ib_bdy), 1, zssh, llrim0 )   ! zssh is masked
+            llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:,ir)   ! possibly every direction, T points
+            llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:,ir)   ! possibly every direction, T points
+         END DO
+         IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction
+            CALL lbc_bdy_lnk( 'bdydyn2d', llsend1, llrecv1, zssh(:,:,1), 'T',  1. )
+         END IF
+      END DO
       !
    END SUBROUTINE bdy_ssh
Index: /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdydyn3d.F90
===================================================================
--- /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdydyn3d.F90	(revision 11190)
+++ /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdydyn3d.F90	(revision 11191)
@@ -42,53 +42,62 @@
       INTEGER, INTENT(in) ::   kt   ! Main time step counter
       !
-      INTEGER ::   ib_bdy   ! loop index
+      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
 
       !!----------------------------------------------------------------------
-      DO ib_bdy=1, nb_bdy
+      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
-      !
-      llsend2(:) = .false.
-      llrecv2(:) = .false.
-      llsend3(:) = .false.
-      llrecv3(:) = .false.
-      DO ib_bdy=1, nb_bdy
-         SELECT CASE( cn_dyn3d(ib_bdy) )
-         CASE('orlanski', 'orlanski_npo')
-            llsend2(:) = llsend2(:) .OR. lsend_bdy(ib_bdy,2,:)   ! possibly every direction, U points
-            llrecv2(:) = llrecv2(:) .OR. lrecv_bdy(ib_bdy,2,:)   ! possibly every direction, U points
-            llsend3(:) = llsend3(:) .OR. lsend_bdy(ib_bdy,3,:)   ! possibly every direction, V points
-            llrecv3(:) = llrecv3(:) .OR. lrecv_bdy(ib_bdy,3,:)   ! possibly every direction, V points
-         CASE('zerograd')
-            llsend2(3:4) = llsend2(3:4) .OR. lsend_bdyint(ib_bdy,2,3:4)   ! north/south, U points
-            llrecv2(3:4) = llrecv2(3:4) .OR. lrecv_bdyint(ib_bdy,2,3:4)   ! north/south, U points
-            llsend3(1:2) = llsend3(1:2) .OR. lsend_bdyint(ib_bdy,3,1:2)   ! west/east, V points
-            llrecv3(1:2) = llrecv3(1:2) .OR. lrecv_bdyint(ib_bdy,3,1:2)   ! west/east, V points
-         CASE('neumann')
-            llsend2(:) = llsend2(:) .OR. lsend_bdyint(ib_bdy,2,:)   ! possibly every direction, U points
-            llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(ib_bdy,2,:)   ! possibly every direction, U points
-            llsend3(:) = llsend3(:) .OR. lsend_bdyint(ib_bdy,3,:)   ! possibly every direction, V points
-            llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(ib_bdy,3,:)   ! 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_bdy_lnk( 'bdydyn2d', llsend2, llrecv2, ua, 'U', -1. )
-      END IF
-      IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction
-         CALL lbc_bdy_lnk( 'bdydyn2d', llsend3, llrecv3, va, 'V', -1. )
-      END IF
+         llsend2(:) = .false.
+         llrecv2(:) = .false.
+         llsend3(:) = .false.
+         llrecv3(:) = .false.
+         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_bdy_lnk( 'bdydyn2d', llsend2, llrecv2, ua, 'U', -1. )
+         END IF
+         IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction
+            CALL lbc_bdy_lnk( 'bdydyn2d', llsend3, llrecv3, va, 'V', -1. )
+         END IF
+      END DO
       !
    END SUBROUTINE bdy_dyn3d
@@ -133,5 +142,5 @@
 
 
-   SUBROUTINE bdy_dyn3d_zgrad( idx, dta, kt , ib_bdy )
+   SUBROUTINE bdy_dyn3d_zgrad( idx, dta, kt, ib_bdy, llrim0 )
       !!----------------------------------------------------------------------
       !!                  ***  SUBROUTINE bdy_dyn3d_zgrad  ***
@@ -141,15 +150,20 @@
       !!----------------------------------------------------------------------
       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
       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)
+      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)
@@ -168,5 +182,8 @@
       !
       igrd = 3                      ! Copying tangential velocity into bdy points
-      DO jb = 1, idx%nblenrim(igrd)
+      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)
@@ -268,5 +285,5 @@
 
 
-   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  ***
@@ -280,6 +297,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
@@ -290,9 +308,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 bdy_orlanski_3d( idx, igrd, vb, va, dta%v3d, ll_npo, llrim0 )
       !
    END SUBROUTINE bdy_dyn3d_orlanski
@@ -347,5 +365,5 @@
 
 
-   SUBROUTINE bdy_dyn3d_nmn( idx, ib_bdy )
+   SUBROUTINE bdy_dyn3d_nmn( idx, ib_bdy, llrim0 )
       !!----------------------------------------------------------------------
       !!                 ***  SUBROUTINE bdy_dyn3d_nmn  ***
@@ -356,7 +374,8 @@
       !!
       !!----------------------------------------------------------------------
-      TYPE(OBC_INDEX),              INTENT(in) ::   idx  ! OBC indices
-      INTEGER,                      INTENT(in) ::   ib_bdy  ! BDY set index
-      INTEGER  ::   igrd                                    ! dummy indice
+      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
       !!----------------------------------------------------------------------
       !
@@ -365,9 +384,9 @@
       igrd = 2      ! Neumann bc on u-velocity; 
       !            
-      CALL bdy_nmn( idx, igrd, ua )   ! ua is masked
+      CALL bdy_nmn( idx, igrd, ua, llrim0 )   ! ua is masked
 
       igrd = 3      ! Neumann bc on v-velocity
       !  
-      CALL bdy_nmn( idx, igrd, va )   ! va is masked
+      CALL bdy_nmn( idx, igrd, va, llrim0 )   ! va is masked
       !
    END SUBROUTINE bdy_dyn3d_nmn
Index: /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdyice.F90
===================================================================
--- /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdyice.F90	(revision 11190)
+++ /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdyice.F90	(revision 11191)
@@ -55,6 +55,8 @@
       INTEGER, INTENT(in) ::   kt   ! Main time step counter
       !
-      INTEGER ::   jbdy                               ! BDY set index
-      LOGICAL, DIMENSION(4)       :: llsend1, llrecv1   ! indicate how communications are to be carried out
+      INTEGER ::   jbdy, ir                             ! BDY set index, rim index
+      INTEGER ::   ibeg, iend                           ! length of rim to be treated (rim 0 or rim 1)
+      LOGICAL ::   llrim0                               ! indicate if rim 0 is treated
+      LOGICAL, DIMENSION(4)  :: llsend1, llrecv1        ! indicate how communications are to be carried out
       !!----------------------------------------------------------------------
       ! controls
@@ -64,34 +66,39 @@
       CALL ice_var_glo2eqv
       !
-      DO jbdy = 1, nb_bdy
+      DO ir = 1, 0, -1   ! treat rim 1 before rim 0
+         IF( ir == 0 ) THEN   ;   llrim0 = .TRUE.
+         ELSE                 ;   llrim0 = .FALSE.
+         END IF
+         DO jbdy = 1, nb_bdy
+            !
+            SELECT CASE( cn_ice(jbdy) )
+            CASE('none')   ;   CYCLE
+            CASE('frs' )   ;   CALL bdy_ice_frs( idx_bdy(jbdy), dta_bdy(jbdy), kt, jbdy, llrim0 )
+            CASE DEFAULT
+               CALL ctl_stop( 'bdy_ice : unrecognised option for open boundaries for ice fields' )
+            END SELECT
+            !
+         END DO
          !
-         SELECT CASE( cn_ice(jbdy) )
-         CASE('none')   ;   CYCLE
-         CASE('frs' )   ;   CALL bdy_ice_frs( idx_bdy(jbdy), dta_bdy(jbdy), kt, jbdy )
-         CASE DEFAULT
-            CALL ctl_stop( 'bdy_ice : unrecognised option for open boundaries for ice fields' )
-         END SELECT
-         !
-      END DO
-      !
-      ! Update bdy points
-      llsend1(:) = .false.
-      llrecv1(:) = .false.
-      DO jbdy = 1, nb_bdy
-         IF( cn_ice(jbdy) == 'frs' ) THEN
-            llsend1(:) = llsend1(:) .OR. lsend_bdyint(jbdy,1,:)   ! possibly every direction, T points
-            llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(jbdy,1,:)   ! possibly every direction, T points
+         ! Update bdy points
+         llsend1(:) = .false.
+         llrecv1(:) = .false.
+         DO jbdy = 1, nb_bdy
+            IF( cn_ice(jbdy) == 'frs' ) THEN
+               llsend1(:) = llsend1(:) .OR. lsend_bdyint(jbdy,1,:,ir)   ! possibly every direction, T points
+               llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(jbdy,1,:,ir)   ! possibly every direction, T points
+            END IF
+         END DO   ! jbdy
+         IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction
+            ! exchange 3d arrays
+            CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, a_i , 'T', 1., h_i , 'T', 1., h_s , 'T', 1. &
+                 &                                            , oa_i, 'T', 1., a_ip, 'T', 1., v_ip, 'T', 1. &
+                 &                                            , s_i , 'T', 1., t_su, 'T', 1., v_i , 'T', 1. &
+                 &                                            , v_s , 'T', 1., sv_i, 'T', 1.                )
+            ! exchange 4d arrays
+            CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, t_s , 'T', 1., e_s , 'T', 1. )   ! third dimension = 1
+            CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, t_i , 'T', 1., e_i , 'T', 1. )   ! third dimension = jpk
          END IF
-      END DO
-      IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction
-         ! exchange 3d arrays
-         CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, a_i , 'T', 1., h_i , 'T', 1., h_s , 'T', 1. &
-              &                                            , oa_i, 'T', 1., a_ip, 'T', 1., v_ip, 'T', 1. &
-              &                                            , s_i , 'T', 1., t_su, 'T', 1., v_i , 'T', 1. &
-              &                                            , v_s , 'T', 1., sv_i, 'T', 1.                )
-         ! exchange 4d arrays
-         CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, t_s , 'T', 1., e_s , 'T', 1. )   ! third dimension = 1
-         CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, t_i , 'T', 1., e_i , 'T', 1. )   ! third dimension = jpk
-      END IF
+      END DO   ! ir
       !
       CALL ice_cor( kt , 0 )      ! -- In case categories are out of bounds, do a remapping
@@ -108,5 +115,5 @@
 
 
-   SUBROUTINE bdy_ice_frs( idx, dta, kt, jbdy )
+   SUBROUTINE bdy_ice_frs( idx, dta, kt, jbdy, llrim0 )
       !!------------------------------------------------------------------------------
       !!                 ***  SUBROUTINE bdy_ice_frs  ***
@@ -117,11 +124,13 @@
       !!             dimensional baroclinic ocean model with realistic topography. Tellus, 365-382.
       !!------------------------------------------------------------------------------
-      TYPE(OBC_INDEX), INTENT(in) ::   idx     ! OBC indices
-      TYPE(OBC_DATA),  INTENT(in) ::   dta     ! OBC external data
-      INTEGER,         INTENT(in) ::   kt      ! main time-step counter
-      INTEGER,         INTENT(in) ::   jbdy    ! BDY set index
+      TYPE(OBC_INDEX), INTENT(in) ::   idx      ! OBC indices
+      TYPE(OBC_DATA),  INTENT(in) ::   dta      ! OBC external data
+      INTEGER,         INTENT(in) ::   kt       ! main time-step counter
+      INTEGER,         INTENT(in) ::   jbdy     ! BDY set index
+      LOGICAL,         INTENT(in) ::   llrim0   ! indicate if rim 0 is treated
       !
       INTEGER  ::   jpbound            ! 0 = incoming ice
       !                                ! 1 = outgoing ice
+      INTEGER  ::   ibeg, iend         ! length of rim to be treated (rim 0 or rim 1)
       INTEGER  ::   i_bdy, jgrd        ! dummy loop indices
       INTEGER  ::   ji, jj, jk, jl, ib, jb
@@ -132,7 +141,10 @@
       !
       jgrd = 1      ! Everything is at T-points here
+      IF( llrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(jgrd)
+      ELSE                ;   ibeg = idx%nblenrim0(jgrd)+1   ;   iend = idx%nblenrim(jgrd)
+      END IF
       !
       DO jl = 1, jpl
-         DO i_bdy = 1, idx%nblenrim(jgrd)
+         DO i_bdy = ibeg, iend
             ji    = idx%nbi(i_bdy,jgrd)
             jj    = idx%nbj(i_bdy,jgrd)
@@ -162,5 +174,5 @@
 
       DO jl = 1, jpl
-         DO i_bdy = 1, idx%nblenrim(jgrd)
+         DO i_bdy = ibeg, iend
             ji = idx%nbi(i_bdy,jgrd)
             jj = idx%nbj(i_bdy,jgrd)
@@ -295,7 +307,8 @@
       CHARACTER(len=1), INTENT(in)  ::   cd_type   ! nature of velocity grid-points
       !
-      INTEGER  ::   i_bdy, jgrd      ! dummy loop indices
-      INTEGER  ::   ji, jj           ! local scalar
-      INTEGER  ::   jbdy             ! BDY set index
+      INTEGER  ::   i_bdy, jgrd       ! dummy loop indices
+      INTEGER  ::   ji, jj            ! local scalar
+      INTEGER  ::   jbdy, ir     ! BDY set index, rim index
+      INTEGER  ::   ibeg, iend   ! length of rim to be treated (rim 0 or rim 1)
       REAL(wp) ::   zmsk1, zmsk2, zflag
       LOGICAL, DIMENSION(4) :: llsend2, llrecv2, llsend3, llrecv3  ! indicate how communications are to be carried out
@@ -303,111 +316,119 @@
       IF( ln_timing )   CALL timing_start('bdy_ice_dyn')
       !
-      DO jbdy=1, nb_bdy
+      DO ir = 1, 0, -1
+         DO jbdy = 1, nb_bdy
+            !
+            SELECT CASE( cn_ice(jbdy) )
+               !
+            CASE('none')
+               CYCLE
+               !
+            CASE('frs')
+               !
+               IF( nn_ice_dta(jbdy) == 0 ) CYCLE            ! case ice boundaries = initial conditions 
+               !                                            !      do not change ice velocity (it is only computed by rheology)
+               SELECT CASE ( cd_type )
+                  !     
+               CASE ( 'U' )  
+                  jgrd = 2      ! u velocity
+                  IF( ir == 0 ) THEN   ;   ibeg = 1                                 ;   iend = idx_bdy(jbdy)%nblenrim0(jgrd)
+                  ELSE                 ;   ibeg = idx_bdy(jbdy)%nblenrim0(jgrd)+1   ;   iend = idx_bdy(jbdy)%nblenrim(jgrd)
+                  END IF
+                  DO i_bdy = ibeg, iend
+                     ji    = idx_bdy(jbdy)%nbi(i_bdy,jgrd)
+                     jj    = idx_bdy(jbdy)%nbj(i_bdy,jgrd)
+                     zflag = idx_bdy(jbdy)%flagu(i_bdy,jgrd)
+                     !     i-1  i   i    |  !        i  i i+1 |  !          i  i i+1 |
+                     !      >  ice  >    |  !        o  > ice |  !          o  >  o  |      
+                     ! => set at u_ice(i-1) !  => set to O       !  => unchanged
+                     IF( zflag == -1. .AND. ji > 1 .AND. ji < jpi )   THEN  
+                        IF    ( vt_i(ji  ,jj) > 0. )   THEN   ;   u_ice(ji,jj) = u_ice(ji-1,jj) 
+                        ELSEIF( vt_i(ji+1,jj) > 0. )   THEN   ;   u_ice(ji,jj) = 0._wp
+                        END IF
+                     END IF
+                     ! |    i  i+1 i+1        !  |  i   i i+1        !  | i  i i+1
+                     ! |    >  ice  >         !  | ice  >  o         !  | o  >  o   
+                     ! => set at u_ice(i+1)   !     => set to O      !     =>  unchanged
+                     IF( zflag ==  1. .AND. ji+1 < jpi+1 )   THEN
+                        IF    ( vt_i(ji+1,jj) > 0. )   THEN   ;   u_ice(ji,jj) = u_ice(ji+1,jj)
+                        ELSEIF( vt_i(ji  ,jj) > 0. )   THEN   ;   u_ice(ji,jj) = 0._wp
+                        END IF
+                     END IF
+                     !
+                     IF( zflag ==  0. )   u_ice(ji,jj) = 0._wp   ! u_ice = 0  if north/south bdy  
+                     !
+                  END DO
+                  !
+               CASE ( 'V' )
+                  jgrd = 3      ! v velocity
+                  IF( ir == 0 ) THEN   ;   ibeg = 1                                 ;   iend = idx_bdy(jbdy)%nblenrim0(jgrd)
+                  ELSE                 ;   ibeg = idx_bdy(jbdy)%nblenrim0(jgrd)+1   ;   iend = idx_bdy(jbdy)%nblenrim(jgrd)
+                  END IF
+                  DO i_bdy = ibeg, iend
+                     ji    = idx_bdy(jbdy)%nbi(i_bdy,jgrd)
+                     jj    = idx_bdy(jbdy)%nbj(i_bdy,jgrd)
+                     zflag = idx_bdy(jbdy)%flagv(i_bdy,jgrd)
+                     !    ¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨     !  ¨¨¨¨ïce¨¨¨(jj+1)¨¨     ! ¨¨¨¨¨¨ö¨¨¨¨(jj+1)       
+                     !       ^    (jj  )       !       ^    (jj  )       !       ^    (jj  )       
+                     !      ice   (jj  )       !       o    (jj  )       !       o    (jj  )       
+                     !       ^    (jj-1)       !                         !
+                     ! => set to u_ice(jj-1)   !  =>   set to 0          !   => unchanged        
+                     IF( zflag == -1. .AND. jj > 1 .AND. jj < jpj )   THEN                 
+                        IF    ( vt_i(ji,jj  ) > 0. )   THEN   ;   v_ice(ji,jj) = v_ice(ji,jj-1)
+                        ELSEIF( vt_i(ji,jj+1) > 0. )   THEN   ;   v_ice(ji,jj) = 0._wp
+                        END IF
+                     END IF
+                     !       ^    (jj+1)       !                         !              
+                     !      ice   (jj+1)       !       o    (jj+1)       !       o    (jj+1)       
+                     !       ^    (jj  )       !       ^    (jj  )       !       ^    (jj  )
+                     !   ________________      !  ____ice___(jj  )_      !  _____o____(jj  ) 
+                     ! => set to u_ice(jj+1)   !    => set to 0          !    => unchanged  
+                     IF( zflag ==  1. .AND. jj < jpj )   THEN              
+                        IF    ( vt_i(ji,jj+1) > 0. )   THEN   ;   v_ice(ji,jj) = v_ice(ji,jj+1)
+                        ELSEIF( vt_i(ji,jj  ) > 0. )   THEN   ;   v_ice(ji,jj) = 0._wp
+                        END IF
+                     END IF
+                     !
+                     IF( zflag ==  0. )   v_ice(ji,jj) = 0._wp   ! v_ice = 0  if west/east bdy  
+                     !
+                  END DO
+                  !
+               END SELECT
+               !
+            CASE DEFAULT
+               CALL ctl_stop( 'bdy_ice_dyn : unrecognised option for open boundaries for ice fields' )
+            END SELECT
+            !
+         END DO    ! jbdy
          !
-         SELECT CASE( cn_ice(jbdy) )
-         !
-         CASE('none')
-            CYCLE
-            !
-         CASE('frs')
-            !
-            IF( nn_ice_dta(jbdy) == 0 ) CYCLE            ! case ice boundaries = initial conditions 
-            !                                            !      do not change ice velocity (it is only computed by rheology)
-            SELECT CASE ( cd_type )
-            !     
-            CASE ( 'U' )  
-               jgrd = 2      ! u velocity
-               DO i_bdy = 1, idx_bdy(jbdy)%nblenrim(jgrd)
-                  ji    = idx_bdy(jbdy)%nbi(i_bdy,jgrd)
-                  jj    = idx_bdy(jbdy)%nbj(i_bdy,jgrd)
-                  zflag = idx_bdy(jbdy)%flagu(i_bdy,jgrd)
-                  !     i-1  i   i    |  !        i  i i+1 |  !          i  i i+1 |
-                  !      >  ice  >    |  !        o  > ice |  !          o  >  o  |      
-                  ! => set at u_ice(i-1) !  => set to O       !  => unchanged
-                  IF( zflag == -1. .AND. ji > 1 .AND. ji < jpi )   THEN  
-                     IF    ( vt_i(ji  ,jj) > 0. )   THEN   ;   u_ice(ji,jj) = u_ice(ji-1,jj) 
-                     ELSEIF( vt_i(ji+1,jj) > 0. )   THEN   ;   u_ice(ji,jj) = 0._wp
-                     END IF
-                  END IF
-                  ! |    i  i+1 i+1        !  |  i   i i+1        !  | i  i i+1
-                  ! |    >  ice  >         !  | ice  >  o         !  | o  >  o   
-                  ! => set at u_ice(i+1)   !     => set to O      !     =>  unchanged
-                  IF( zflag ==  1. .AND. ji+1 < jpi+1 )   THEN
-                     IF    ( vt_i(ji+1,jj) > 0. )   THEN   ;   u_ice(ji,jj) = u_ice(ji+1,jj)
-                     ELSEIF( vt_i(ji  ,jj) > 0. )   THEN   ;   u_ice(ji,jj) = 0._wp
-                     END IF
-                  END IF
-                  !
-                  IF( zflag ==  0. )   u_ice(ji,jj) = 0._wp   ! u_ice = 0  if north/south bdy  
-                  !
-               END DO
-               !
-            CASE ( 'V' )
-               jgrd = 3      ! v velocity
-               DO i_bdy = 1, idx_bdy(jbdy)%nblenrim(jgrd)
-                  ji    = idx_bdy(jbdy)%nbi(i_bdy,jgrd)
-                  jj    = idx_bdy(jbdy)%nbj(i_bdy,jgrd)
-                  zflag = idx_bdy(jbdy)%flagv(i_bdy,jgrd)
-                  !    ¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨     !  ¨¨¨¨ïce¨¨¨(jj+1)¨¨     ! ¨¨¨¨¨¨ö¨¨¨¨(jj+1)       
-                  !       ^    (jj  )       !       ^    (jj  )       !       ^    (jj  )       
-                  !      ice   (jj  )       !       o    (jj  )       !       o    (jj  )       
-                  !       ^    (jj-1)       !                         !
-                  ! => set to u_ice(jj-1)   !  =>   set to 0          !   => unchanged        
-                  IF( zflag == -1. .AND. jj > 1 .AND. jj < jpj )   THEN                 
-                    IF    ( vt_i(ji,jj  ) > 0. )   THEN   ;   v_ice(ji,jj) = v_ice(ji,jj-1)
-                    ELSEIF( vt_i(ji,jj+1) > 0. )   THEN   ;   v_ice(ji,jj) = 0._wp
-                    END IF
-                  END IF 
-                  !       ^    (jj+1)       !                         !              
-                  !      ice   (jj+1)       !       o    (jj+1)       !       o    (jj+1)       
-                  !       ^    (jj  )       !       ^    (jj  )       !       ^    (jj  )
-                  !   ________________      !  ____ice___(jj  )_      !  _____o____(jj  ) 
-                  ! => set to u_ice(jj+1)   !    => set to 0          !    => unchanged  
-                  IF( zflag ==  1. .AND. jj < jpj )   THEN              
-                     IF    ( vt_i(ji,jj+1) > 0. )   THEN   ;   v_ice(ji,jj) = v_ice(ji,jj+1)
-                     ELSEIF( vt_i(ji,jj  ) > 0. )   THEN   ;   v_ice(ji,jj) = 0._wp
-                     END IF
-                  END IF                                          
-                  !
-                  IF( zflag ==  0. )   v_ice(ji,jj) = 0._wp   ! v_ice = 0  if west/east bdy  
-                  !
-               END DO
-               !
-            END SELECT
-            !
-         CASE DEFAULT
-            CALL ctl_stop( 'bdy_ice_dyn : unrecognised option for open boundaries for ice fields' )
+         SELECT CASE ( cd_type )        
+         CASE ( 'U' ) 
+            llsend2(:) = .false.   ;   llrecv2(:) = .false.
+            DO jbdy = 1, nb_bdy
+               IF( cn_ice(jbdy) == 'frs' .AND. nn_ice_dta(jbdy) /= 0 ) THEN
+                  llsend2(:) = llsend2(:) .OR. lsend_bdyint(jbdy,2,:,ir)   ! possibly every direction, U points
+                  llsend2(1) = llsend2(1) .OR. lsend_bdyext(jbdy,2,1,ir)   ! neighbour might search point towards its west bdy
+                  llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(jbdy,2,:,ir)   ! possibly every direction, U points
+                  llrecv2(2) = llrecv2(2) .OR. lrecv_bdyext(jbdy,2,2,ir)   ! might search point towards east bdy
+               END IF
+            END DO
+            IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN   ! if need to send/recv in at least one direction
+               CALL lbc_bdy_lnk( 'bdyice', llsend2, llrecv2, u_ice, 'U', -1. )
+            END IF
+         CASE ( 'V' )
+            llsend3(:) = .false.   ;   llrecv3(:) = .false.
+            DO jbdy = 1, nb_bdy
+               IF( cn_ice(jbdy) == 'frs' .AND. nn_ice_dta(jbdy) /= 0 ) THEN
+                  llsend3(:) = llsend3(:) .OR. lsend_bdyint(jbdy,3,:,ir)   ! possibly every direction, V points
+                  llsend3(3) = llsend3(3) .OR. lsend_bdyext(jbdy,3,3,ir)   ! neighbour might search point towards its south bdy
+                  llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(jbdy,3,:,ir)   ! possibly every direction, V points
+                  llrecv3(4) = llrecv3(4) .OR. lrecv_bdyext(jbdy,3,4,ir)   ! might search point towards north bdy
+               END IF
+            END DO
+            IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction
+               CALL lbc_bdy_lnk( 'bdyice', llsend3, llrecv3, v_ice, 'V', -1. )
+            END IF
          END SELECT
-         !
-      END DO   ! jbdy
-      !
-      SELECT CASE ( cd_type )        
-      CASE ( 'U' ) 
-         llsend2(:) = .false.   ;   llrecv2(:) = .false.
-         DO jbdy = 1, nb_bdy
-            IF( cn_ice(jbdy) == 'frs' .AND. nn_ice_dta(jbdy) /= 0 ) THEN
-               llsend2(:) = llsend2(:) .OR. lsend_bdyint(jbdy,2,:)   ! possibly every direction, U points
-               llsend2(1) = llsend2(1) .OR. lsend_bdyext(jbdy,2,1)   ! neighbour might search point towards bdy on its east
-               llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(jbdy,2,:)   ! possibly every direction, U points
-               llrecv2(2) = llrecv2(2) .OR. lrecv_bdyext(jbdy,2,2)   ! might search point towards bdy on the east
-            END IF
-         END DO
-         IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN   ! if need to send/recv in at least one direction
-            CALL lbc_bdy_lnk( 'bdyice', llsend2, llrecv2, u_ice, 'U', -1. )
-         END IF
-      CASE ( 'V' )
-         llsend3(:) = .false.   ;   llrecv3(:) = .false.
-         DO jbdy = 1, nb_bdy
-            IF( cn_ice(jbdy) == 'frs' .AND. nn_ice_dta(jbdy) /= 0 ) THEN
-               llsend3(:) = llsend3(:) .OR. lsend_bdyint(jbdy,3,:)   ! possibly every direction, V points
-               llsend3(3) = llsend3(3) .OR. lsend_bdyext(jbdy,3,3)   ! neighbour might search point towards bdy on its north
-               llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(jbdy,3,:)   ! possibly every direction, V points
-               llrecv3(4) = llrecv3(4) .OR. lrecv_bdyext(jbdy,3,4)   ! might search point towards bdy on the north
-            END IF
-         END DO
-         IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction
-            CALL lbc_bdy_lnk( 'bdyice', llsend3, llrecv3, v_ice, 'V', -1. )
-         END IF
-      END SELECT
+      END DO   ! ir
       !
       IF( ln_timing )   CALL timing_stop('bdy_ice_dyn')
Index: /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdyini.F90
===================================================================
--- /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdyini.F90	(revision 11190)
+++ /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdyini.F90	(revision 11191)
@@ -33,5 +33,5 @@
    PRIVATE
 
-   PUBLIC   bdy_init   ! routine called in nemo_init
+   PUBLIC   bdy_init    ! routine called in nemo_init
    PUBLIC   find_neib   ! routine called in bdy_nmn
 
@@ -126,15 +126,13 @@
       !! ** Input   :  bdy_init.nc, input file for unstructured open boundaries
       !!----------------------------------------------------------------------      
-      INTEGER  ::   ib_bdy, ii, ij, ik, igrd, ib, ir, iseg ! dummy loop indices
-      INTEGER  ::   icount, icountr, ibr_max, ilen1, ibm1  ! local integers
+      INTEGER  ::   ib_bdy, ii, ij, igrd, ib, ir, iseg     ! dummy loop indices
+      INTEGER  ::   icount, icountr, icountr0, ibr_max     ! local integers
+      INTEGER  ::   ilen1, ibm1                            !   -       -
       INTEGER  ::   iwe, ies, iso, ino, inum, id_dummy     !   -       -
-      INTEGER  ::   igrd_start, igrd_end, jpbdta           !   -       -
-      INTEGER  ::   jpbdtau, jpbdtas                       !   -       -
+      INTEGER  ::   jpbdta, jpbdtau, jpbdtas               !   -       -
       INTEGER  ::   ib_bdy1, ib_bdy2, ib1, ib2             !   -       -
-      INTEGER  ::   i_offset, j_offset, inn                !   -       -
       INTEGER  ::   ii1, ii2, ii3, ij1, ij2, ij3           !   -       -
       INTEGER  ::   iibe, ijbe, iibi, ijbi                 !   -       -
       INTEGER  ::   flagu, flagv                           ! short cuts
-      INTEGER , POINTER  ::  nbi, nbj, nbr                 ! short cuts
       REAL(wp), POINTER, DIMENSION(:,:)       ::   pmask   ! pointer to 2D mask fields
       REAL(wp) ::   zefl, zwfl, znfl, zsfl                 ! local scalars
@@ -145,9 +143,6 @@
       CHARACTER(LEN=1),DIMENSION(jpbgrd)      ::   cgrid
       INTEGER :: jpk_max
-      LOGICAL :: llsend_ea, llsend_we, llsend_so, llsend_no  ! Flags for boundaries sending
-      LOGICAL :: llrecv_ea, llrecv_we, llrecv_so, llrecv_no  ! Flags for boundaries receiving
-      REAL(wp), TARGET, DIMENSION(jpi,jpj) ::   zfmask  ! temporary fmask array excluding coastal boundary condition (shlat)
-      REAL(wp)        , DIMENSION(jpi,jpj) ::   ztmp
-      LOGICAL  ::   llnon, llson, llean, llwen          ! local logicals
+      REAL(wp), DIMENSION(jpi,jpj) ::   zfmask                  ! temporary fmask array excluding coastal boundary condition (shlat)
+      REAL(wp), DIMENSION(jpi,jpj) ::   ztmask, zumask, zvmask  ! temporary u/v mask array
       !!
       CHARACTER(LEN=1)                     ::   ctypebdy   !     -        - 
@@ -798,14 +793,16 @@
       !
       iwe = mig(1)
-      ies = mig(nlci)
+      ies = mig(jpi)
       iso = mjg(1) 
-      ino = mjg(nlcj) 
+      ino = mjg(jpj) 
       !
       DO ib_bdy = 1, nb_bdy
          DO igrd = 1, jpbgrd
-            icount  = 0   ! initialization of local bdy length
-            icountr = 0   ! initialization of local rim bdy length
-            idx_bdy(ib_bdy)%nblen(igrd)    = 0
-            idx_bdy(ib_bdy)%nblenrim(igrd) = 0
+            icount   = 0   ! initialization of local bdy length
+            icountr  = 0   ! initialization of local rim 0 and rim 1 bdy length
+            icountr0 = 0   ! initialization of local rim 0 bdy length
+            idx_bdy(ib_bdy)%nblen(igrd)     = 0
+            idx_bdy(ib_bdy)%nblenrim(igrd)  = 0
+            idx_bdy(ib_bdy)%nblenrim0(igrd) = 0
             DO ib = 1, nblendta(igrd,ib_bdy)
                ! check that data is in correct order in file
@@ -823,9 +820,11 @@
                   !
                   icount = icount + 1
-                  IF( nbrdta(ib,igrd,ib_bdy) == 1 )   icountr = icountr + 1
+                  IF( nbrdta(ib,igrd,ib_bdy) == 1 .OR. nbrdta(ib,igrd,ib_bdy) == 0 )   icountr = icountr + 1
+                  IF( nbrdta(ib,igrd,ib_bdy) == 0 )   icountr0 = icountr0 + 1
                ENDIF
             END DO
-            idx_bdy(ib_bdy)%nblenrim(igrd) = icountr !: length of rim boundary data on each proc
-            idx_bdy(ib_bdy)%nblen   (igrd) = icount  !: length of boundary data on each proc        
+            idx_bdy(ib_bdy)%nblen    (igrd) = icount   !: length of boundary data on each proc
+            idx_bdy(ib_bdy)%nblenrim (igrd) = icountr  !: length of rim 0 and rim 1 boundary data on each proc   
+            idx_bdy(ib_bdy)%nblenrim0(igrd) = icountr0 !: length of rim 0 boundary data on each proc     
          END DO   ! igrd
 
@@ -849,5 +848,5 @@
             icount  = 0
             ! Outer loop on rimwidth to ensure outermost points come first in the local arrays.
-            DO ir=1, nn_rimwidth(ib_bdy)
+            DO ir = 0, nn_rimwidth(ib_bdy)
                DO ib = 1, nblendta(igrd,ib_bdy)
                   ! check if point is in local domain and equals ir
@@ -870,7 +869,7 @@
       ! Initialize array indicating communications in bdy
       ! -------------------------------------------------
-      ALLOCATE( lsend_bdy(nb_bdy,jpbgrd,4), lrecv_bdy(nb_bdy,jpbgrd,4) )
-      lsend_bdy(:,:,:) = .false.
-      lrecv_bdy(:,:,:) = .false. 
+      ALLOCATE( lsend_bdy(nb_bdy,jpbgrd,4,0:1), lrecv_bdy(nb_bdy,jpbgrd,4,0:1) )
+      lsend_bdy(:,:,:,:) = .false.
+      lrecv_bdy(:,:,:,:) = .false. 
 
       DO ib_bdy = 1, nb_bdy
@@ -879,24 +878,27 @@
                ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
                ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
+               IF( ib .LE. idx_bdy(ib_bdy)%nblenrim0(igrd) ) THEN   ;   ir = 0
+               ELSE                                                 ;   ir = 1
+               END IF
                !
                ! check if point has to be sent     to   a neighbour
                ! W neighbour and on the inner left  side
-               IF( ii == 2     .and. (nbondi == 0 .or. nbondi ==  1) )   lsend_bdy(ib_bdy,igrd,1) = .true.
+               IF( ii == 2     .and. (nbondi == 0 .or. nbondi ==  1) )   lsend_bdy(ib_bdy,igrd,1,ir) = .true.
                ! E neighbour and on the inner right side
-               IF( ii == jpi-1 .and. (nbondi == 0 .or. nbondi == -1) )   lsend_bdy(ib_bdy,igrd,2) = .true.
+               IF( ii == jpi-1 .and. (nbondi == 0 .or. nbondi == -1) )   lsend_bdy(ib_bdy,igrd,2,ir) = .true.
                ! S neighbour and on the inner down side
-               IF( ij == 2     .and. (nbondj == 0 .or. nbondj ==  1) )   lsend_bdy(ib_bdy,igrd,3) = .true.
+               IF( ij == 2     .and. (nbondj == 0 .or. nbondj ==  1) )   lsend_bdy(ib_bdy,igrd,3,ir) = .true.
                ! N neighbour and on the inner up   side
-               IF( ij == jpj-1 .and. (nbondj == 0 .or. nbondj == -1) )   lsend_bdy(ib_bdy,igrd,4) = .true.
+               IF( ij == jpj-1 .and. (nbondj == 0 .or. nbondj == -1) )   lsend_bdy(ib_bdy,igrd,4,ir) = .true.
                !
                ! check if point has to be received from a neighbour
                ! W neighbour and on the outter left  side
-               IF( ii == 1   .and. (nbondi == 0 .or. nbondi ==  1) )    lrecv_bdy(ib_bdy,igrd,1) = .true.
+               IF( ii == 1     .and. (nbondi == 0 .or. nbondi ==  1) )   lrecv_bdy(ib_bdy,igrd,1,ir) = .true.
                ! E neighbour and on the outter right side
-               IF( ii == jpi .and. (nbondi == 0 .or. nbondi == -1) )    lrecv_bdy(ib_bdy,igrd,2) = .true.
+               IF( ii == jpi   .and. (nbondi == 0 .or. nbondi == -1) )   lrecv_bdy(ib_bdy,igrd,2,ir) = .true.
                ! S neighbour and on the outter down side
-               IF( ij == 1   .and. (nbondj == 0 .or. nbondj ==  1) )    lrecv_bdy(ib_bdy,igrd,3) = .true.
+               IF( ij == 1     .and. (nbondj == 0 .or. nbondj ==  1) )   lrecv_bdy(ib_bdy,igrd,3,ir) = .true.
                ! N neighbour and on the outter up   side
-               IF( ij == jpj .and. (nbondj == 0 .or. nbondj == -1) )    lrecv_bdy(ib_bdy,igrd,4) = .true.
+               IF( ij == jpj   .and. (nbondj == 0 .or. nbondj == -1) )   lrecv_bdy(ib_bdy,igrd,4,ir) = .true.
                !
             END DO
@@ -907,8 +909,8 @@
          DO igrd = 1, jpbgrd
             DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
-               nbr => idx_bdy(ib_bdy)%nbr(ib,igrd)
-               idx_bdy(ib_bdy)%nbw(ib,igrd) = 1.- TANH( REAL( nbr - 1 ) *0.5 )      ! tanh formulation
-!               idx_bdy(ib_bdy)%nbw(ib,igrd) = (REAL(nn_rimwidth(ib_bdy)+1-nbr)/REAL(nn_rimwidth(ib_bdy)))**2.  ! quadratic
-!               idx_bdy(ib_bdy)%nbw(ib,igrd) =  REAL(nn_rimwidth(ib_bdy)+1-nbr)/REAL(nn_rimwidth(ib_bdy))       ! linear
+               ir = MAX( 1, idx_bdy(ib_bdy)%nbr(ib,igrd) )   ! both rim 0 and rim 1 have the same weights
+               idx_bdy(ib_bdy)%nbw(ib,igrd) = 1.- TANH( REAL( ir - 1 ) *0.5 )      ! tanh formulation
+!               idx_bdy(ib_bdy)%nbw(ib,igrd) = (REAL(nn_rimwidth(ib_bdy)+1-ir)/REAL(nn_rimwidth(ib_bdy)))**2.  ! quadratic
+!               idx_bdy(ib_bdy)%nbw(ib,igrd) =  REAL(nn_rimwidth(ib_bdy)+1-ir)/REAL(nn_rimwidth(ib_bdy))       ! linear
             END DO
          END DO 
@@ -918,9 +920,9 @@
          DO igrd = 1, jpbgrd
             DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd)
-               nbr => idx_bdy(ib_bdy)%nbr(ib,igrd)
+               ir = MAX( 1, idx_bdy(ib_bdy)%nbr(ib,igrd) )   ! both rim 0 and rim 1 have the same damping coefficients
                idx_bdy(ib_bdy)%nbd(ib,igrd) = 1. / ( rn_time_dmp(ib_bdy) * rday ) & 
-               & *(REAL(nn_rimwidth(ib_bdy)+1-nbr)/REAL(nn_rimwidth(ib_bdy)))**2.   ! quadratic
+               & *(REAL(nn_rimwidth(ib_bdy)+1-ir)/REAL(nn_rimwidth(ib_bdy)))**2.   ! quadratic
                idx_bdy(ib_bdy)%nbdout(ib,igrd) = 1. / ( rn_time_dmp_out(ib_bdy) * rday ) & 
-               & *(REAL(nn_rimwidth(ib_bdy)+1-nbr)/REAL(nn_rimwidth(ib_bdy)))**2.   ! quadratic
+               & *(REAL(nn_rimwidth(ib_bdy)+1-ir)/REAL(nn_rimwidth(ib_bdy)))**2.   ! quadratic
             END DO
          END DO 
@@ -931,4 +933,19 @@
       ! Initialise masks and find normal/tangential directions
       ! ------------------------------------------------------
+
+      ! ------------------------------------------
+      ! handel rim0, do as if rim 1 was free ocean
+      ! ------------------------------------------
+
+      ztmask(:,:) = tmask(:,:,1)   ;   zumask(:,:) = umask(:,:,1)   ;   zvmask(:,:) = vmask(:,:,1)
+      ! For the flagu/flagv calculation below we require a version of fmask without
+      ! the land boundary condition (shlat) included:
+      DO ij = 2, jpjm1
+         DO ii = 2, jpim1
+            zfmask(ii,ij) = ztmask(ii,ij  ) * ztmask(ii+1,ij  )   &
+           &              * ztmask(ii,ij+1) * ztmask(ii+1,ij+1)
+         END DO      
+      END DO
+      CALL lbc_lnk( 'bdyini', zfmask, 'F', 1. )
 
       ! Read global 2D mask at T-points: bdytmask
@@ -940,53 +957,195 @@
 
       ! Derive mask on U and V grid from mask on T grid
-
-      bdyumask(:,:) = 0._wp
-      bdyvmask(:,:) = 0._wp
       DO ij = 1, jpjm1
          DO ii = 1, jpim1
-            bdyumask(ii,ij) = bdytmask(ii,ij) * bdytmask(ii+1, ij )
+            bdyumask(ii,ij) = bdytmask(ii,ij) * bdytmask(ii+1,ij  )
             bdyvmask(ii,ij) = bdytmask(ii,ij) * bdytmask(ii  ,ij+1)  
          END DO
       END DO
-      CALL lbc_lnk_multi( 'bdyini', bdyumask, 'U', 1. , bdyvmask, 'V', 1. )   ! Lateral boundary cond. 
-
-      ! bdy masks are now set to zero on boundary points:
-      !
-      igrd = 1
+      CALL lbc_lnk_multi( 'bdyini', bdyumask, 'U', 1., bdyvmask, 'V', 1. )   ! Lateral boundary cond.
+
+      ! bdy masks are now set to zero on rim 0 points:
       DO ib_bdy = 1, nb_bdy
-        DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)      
-          bdytmask(idx_bdy(ib_bdy)%nbi(ib,igrd), idx_bdy(ib_bdy)%nbj(ib,igrd)) = 0._wp
-        END DO
+         DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(1)   ! extent of rim 0
+            bdytmask(idx_bdy(ib_bdy)%nbi(ib,1), idx_bdy(ib_bdy)%nbj(ib,1)) = 0._wp
+         END DO
+         DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(2)   ! extent of rim 0
+            bdyumask(idx_bdy(ib_bdy)%nbi(ib,2), idx_bdy(ib_bdy)%nbj(ib,2)) = 0._wp
+         END DO
+         DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(3)   ! extent of rim 0
+            bdyvmask(idx_bdy(ib_bdy)%nbi(ib,3), idx_bdy(ib_bdy)%nbj(ib,3)) = 0._wp
+         END DO
       END DO
-      !
-      igrd = 2
+
+      CALL bdy_rim_treat( zumask, zvmask, zfmask, .true. )   ! compute flagu, flagv, ntreat on rim 0
+
+      ! ------------------------------------
+      ! handel rim1, do as if rim 0 was land
+      ! ------------------------------------
+
       DO ib_bdy = 1, nb_bdy
-        DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)
-          bdyumask(idx_bdy(ib_bdy)%nbi(ib,igrd), idx_bdy(ib_bdy)%nbj(ib,igrd)) = 0._wp
-        END DO
+         DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(1)   ! extent of rim 0
+            ztmask(idx_bdy(ib_bdy)%nbi(ib,1), idx_bdy(ib_bdy)%nbj(ib,1)) = 0._wp
+         END DO
+         DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(2)   ! extent of rim 0
+            zumask(idx_bdy(ib_bdy)%nbi(ib,2), idx_bdy(ib_bdy)%nbj(ib,2)) = 0._wp
+         END DO
+         DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(3)   ! extent of rim 0
+            zvmask(idx_bdy(ib_bdy)%nbi(ib,3), idx_bdy(ib_bdy)%nbj(ib,3)) = 0._wp
+         END DO
       END DO
-      !
-      igrd = 3
-      DO ib_bdy = 1, nb_bdy
-        DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)
-          bdyvmask(idx_bdy(ib_bdy)%nbi(ib,igrd), idx_bdy(ib_bdy)%nbj(ib,igrd)) = 0._wp
-        END DO
-      END DO
-
-      ! For the flagu/flagv calculation below we require a version of fmask without
-      ! the land boundary condition (shlat) included:
-      zfmask(:,:) = 0
+
+      ! Recompute zfmask
       DO ij = 2, jpjm1
          DO ii = 2, jpim1
-            zfmask(ii,ij) = tmask(ii,ij  ,1) * tmask(ii+1,ij  ,1)   &
-           &              * tmask(ii,ij+1,1) * tmask(ii+1,ij+1,1)
+            zfmask(ii,ij) = ztmask(ii,ij  ) * ztmask(ii+1,ij  )   &
+           &              * ztmask(ii,ij+1) * ztmask(ii+1,ij+1)
          END DO      
       END DO
-
-      ! Lateral boundary conditions
-      CALL lbc_lnk( 'bdyini', zfmask, 'F', 1. ) 
-      CALL lbc_lnk_multi( 'bdyini', bdyumask, 'U', 1. , bdyvmask, 'V', 1., bdytmask, 'T', 1. )
+      CALL lbc_lnk( 'bdyini', zfmask, 'F', 1. )
+
+      ! bdy masks are now set to zero on rim1 points:
+      DO ib_bdy = 1, nb_bdy
+         DO ib = idx_bdy(ib_bdy)%nblenrim0(1) + 1,  idx_bdy(ib_bdy)%nblenrim(1)   ! extent of rim 1
+            bdytmask(idx_bdy(ib_bdy)%nbi(ib,1), idx_bdy(ib_bdy)%nbj(ib,1)) = 0._wp
+         END DO
+         DO ib = idx_bdy(ib_bdy)%nblenrim0(2) + 1,  idx_bdy(ib_bdy)%nblenrim(2)   ! extent of rim 1
+            bdyumask(idx_bdy(ib_bdy)%nbi(ib,2), idx_bdy(ib_bdy)%nbj(ib,2)) = 0._wp
+         END DO
+         DO ib = idx_bdy(ib_bdy)%nblenrim0(3) + 1,  idx_bdy(ib_bdy)%nblenrim(3)   ! extent of rim 1
+            bdyvmask(idx_bdy(ib_bdy)%nbi(ib,3), idx_bdy(ib_bdy)%nbj(ib,3)) = 0._wp
+         END DO
+      END DO
+
+      CALL bdy_rim_treat( zumask, zvmask, zfmask, .false. )   ! compute flagu, flagv, ntreat on rim 1
+
+
+      !
+      ! Check which boundaries might need communication
+      ALLOCATE( lsend_bdyint(nb_bdy,jpbgrd,4,0:1), lrecv_bdyint(nb_bdy,jpbgrd,4,0:1) )
+      lsend_bdyint(:,:,:,:) = .false.
+      lrecv_bdyint(:,:,:,:) = .false. 
+      ALLOCATE( lsend_bdyext(nb_bdy,jpbgrd,4,0:1), lrecv_bdyext(nb_bdy,jpbgrd,4,0:1) )
+      lsend_bdyext(:,:,:,:) = .false.
+      lrecv_bdyext(:,:,:,:) = .false.
+      !
+      DO igrd = 1, jpbgrd
+         DO ib_bdy = 1, nb_bdy
+            DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)
+               ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
+               ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
+               IF( ib .LE. idx_bdy(ib_bdy)%nblenrim0(igrd) ) THEN   ;   ir = 0
+               ELSE                                                 ;   ir = 1
+               END IF
+               flagu = NINT(idx_bdy(ib_bdy)%flagu(ib,igrd))
+               flagv = NINT(idx_bdy(ib_bdy)%flagv(ib,igrd))
+               iibe = ii - flagu   ! neighbouring point towards the exterior of the computational domain
+               ijbe = ij - flagv
+               iibi = ii + flagu   ! neighbouring point towards the interior of the computational domain
+               ijbi = ij + flagv
+               CALL find_neib( ii, ij, idx_bdy(ib_bdy)%ntreat(ib,igrd), ii1, ij1, ii2, ij2, ii3, ij3 )   ! free ocean neighbours
+               !
+               ! search neighbour in the  west/east  direction
+               ! Rim is on the halo and computed ocean is towards exterior of mpi domain  
+               !      <--    (o exterior)     -->  
+               ! (1)  o|x         OR    (2)   x|o
+               !       |___                 ___| 
+               IF( iibi == 0     .OR. ii1 == 0     .OR. ii2 == 0     .OR. ii3 == 0     )   lrecv_bdyint(ib_bdy,igrd,1,ir) = .true.
+               IF( iibi == jpi+1 .OR. ii1 == jpi+1 .OR. ii2 == jpi+1 .OR. ii3 == jpi+1 )   lrecv_bdyint(ib_bdy,igrd,2,ir) = .true.  
+               IF( iibe == 0     )   lrecv_bdyext(ib_bdy,igrd,1,ir) = .true.
+               IF( iibe == jpi+1 )   lrecv_bdyext(ib_bdy,igrd,2,ir) = .true.  
+               ! Check if neighbour has its rim parallel to its mpi subdomain border and located next to its halo
+               ! :¨¨¨¨¨|¨¨-->    |                                             |    <--¨¨|¨¨¨¨¨: 
+               ! :     |  x:o    |    neighbour limited by ... would need o    |    o:x  |     :
+               ! :.....|_._:_____|   (1) W neighbour         E neighbour (2)   |_____:_._|.....:
+               IF( ii == 2     .AND. ( nbondi ==  1 .OR. nbondi == 0 ) .AND. &
+                 & ( iibi == 3     .OR. ii1 == 3     .OR. ii2 == 3     .OR. ii3 == 3    ) )   lsend_bdyint(ib_bdy,igrd,1,ir) = .true.
+               IF( ii == jpi-1 .AND. ( nbondi == -1 .OR. nbondi == 0 ) .AND. &
+                 & ( iibi == jpi-2 .OR. ii1 == jpi-2 .OR. ii2 == jpi-2 .OR. ii3 == jpi-2) )   lsend_bdyint(ib_bdy,igrd,2,ir) = .true.
+               IF( ii == 2     .AND. ( nbondi ==  1 .OR. nbondi == 0 ) .AND. iibe == 3    )   lsend_bdyext(ib_bdy,igrd,1,ir) = .true.
+               IF( ii == jpi-1 .AND. ( nbondi == -1 .OR. nbondi == 0 ) .AND. iibe == jpi-2)   lsend_bdyext(ib_bdy,igrd,2,ir) = .true.
+               !
+               ! search neighbour in the north/south direction   
+               ! Rim is on the halo and computed ocean is towards exterior of mpi domain
+               !(3)   |       |         ^   ___o___     
+               !  |   |___x___|   OR    |  |   x   |
+               !  v       o           (4)  |       |
+               IF( ijbi == 0     .OR. ij1 == 0     .OR. ij2 == 0     .OR. ij3 == 0     )   lrecv_bdyint(ib_bdy,igrd,3,ir) = .true.
+               IF( ijbi == jpj+1 .OR. ij1 == jpj+1 .OR. ij2 == jpj+1 .OR. ij3 == jpj+1 )   lrecv_bdyint(ib_bdy,igrd,4,ir) = .true.
+               IF( ijbe == 0     )   lrecv_bdyext(ib_bdy,igrd,3,ir) = .true.
+               IF( ijbe == jpj+1 )   lrecv_bdyext(ib_bdy,igrd,4,ir) = .true.
+               ! Check if neighbour has its rim parallel to its mpi subdomain     _________  border and next to its halo
+               !   ^  |    o    |                                                :         : 
+               !   |  |¨¨¨¨x¨¨¨¨|   neighbour limited by ... would need o     |  |....x....|
+               !      :_________:  (3) S neighbour          N neighbour (4)   v  |    o    |   
+               IF( ij == 2     .AND. ( nbondj ==  1 .OR. nbondj == 0 ) .AND. &
+                 & ( ijbi == 3     .OR. ij1 == 3     .OR. ij2 == 3     .OR. ij3 == 3    ) )   lsend_bdyint(ib_bdy,igrd,3,ir) = .true.
+               IF( ij == jpj-1 .AND. ( nbondj == -1 .OR. nbondj == 0 ) .AND. &
+                 & ( ijbi == jpj-2 .OR. ij1 == jpj-2 .OR. ij2 == jpj-2 .OR. ij3 == jpj-2) )   lsend_bdyint(ib_bdy,igrd,4,ir) = .true.
+               IF( ij == 2     .AND. ( nbondj ==  1 .OR. nbondj == 0 ) .AND. ijbe == 3    )   lsend_bdyext(ib_bdy,igrd,3,ir) = .true.
+               IF( ij == jpj-1 .AND. ( nbondj == -1 .OR. nbondj == 0 ) .AND. ijbe == jpj-2)   lsend_bdyext(ib_bdy,igrd,4,ir) = .true.
+            END DO
+         END DO
+      END DO
+
+      DO ib_bdy = 1,nb_bdy
+         IF( cn_dyn2d(ib_bdy) == 'orlanski' .OR. cn_dyn2d(ib_bdy) == 'orlanski_npo' .OR. &
+           & cn_dyn3d(ib_bdy) == 'orlanski' .OR. cn_dyn3d(ib_bdy) == 'orlanski_npo' .OR. &
+           & cn_tra(ib_bdy)   == 'orlanski' .OR. cn_tra(ib_bdy)   == 'orlanski_npo'      ) THEN
+            DO igrd = 1, jpbgrd
+               DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)
+                  ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
+                  ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
+                  IF(  mig(ii) > 2 .AND. mig(ii) < jpiglo-2 .AND. mjg(ij) > 2 .AND. mjg(ij) < jpjglo-2  ) THEN
+                     WRITE(ctmp1,*) ' Orlanski is not safe when the open boundaries are on the interior of the computational domain'
+                     WRITE(ctmp2,*) ' ========== '
+                     CALL ctl_stop( ' ', ctmp1, ctmp2, ' ' )
+                  END IF
+               END DO
+            END DO
+         END IF
+      END DO
+
+      !
+      ! Tidy up
+      !--------
+      IF( nb_bdy>0 )   DEALLOCATE( nbidta, nbjdta, nbrdta )
+      !
+   END SUBROUTINE bdy_segs
+
+
+   SUBROUTINE bdy_rim_treat( zumask, zvmask, zfmask, lrim0 )
+      !!----------------------------------------------------------------------
+      !!                 ***  ROUTINE bdy_rim_treat  ***
+      !!
+      !! ** Purpose :   Initialize structures ( flagu, flagv, ntreat ) indicating how rim points
+      !!                  are to be handeled in the boundary condition treatment
+      !!
+      !! ** Method  :   - to handel rim 0 zmasks must indicate ocean points      (set at one on rim 0 and rim 1 and interior)
+      !!                            and bdymasks must be set at 0 on rim 0       (set at one on rim 1 and interior)
+      !!                    (as if rim 1 was free ocean)
+      !!                - to handel rim 1 zmasks must be set at 0 on rim 0       (set at one on rim 1 and interior)
+      !!                            and bdymasks must indicate free ocean points (set at one on interior)
+      !!                    (as if rim 0 was land)
+      !!                - we can then check in which direction the interior of the computational domain is with the difference
+      !!                         mask array values on both sides to compute flagu and flagv
+      !!                - and look at the ocean neighbours to compute ntreat
+      !!----------------------------------------------------------------------
+      REAL(wp), TARGET, DIMENSION(jpi,jpj), INTENT (in   ) :: zfmask   ! temporary fmask excluding coastal boundary condition (shlat)
+      REAL(wp), TARGET, DIMENSION(jpi,jpj), INTENT (in   ) :: zumask, zvmask   ! temporary t/u/v mask array
+      LOGICAL                             , INTENT (in   ) :: lrim0    ! .true. -> rim 0   .false. -> rim 1
+      INTEGER  ::   ib_bdy, ii, ij, igrd, ib, icount       ! dummy loop indices
+      INTEGER  ::   i_offset, j_offset, inn                ! local integer
+      INTEGER  ::   ibeg, iend                             ! local integer
+      LOGICAL  ::   llnon, llson, llean, llwen             ! local logicals indicating the presence of a ocean neighbour
+      REAL(wp), POINTER, DIMENSION(:,:)       ::   pmask   ! pointer to 2D mask fields
+      REAL(wp) ::   zefl, zwfl, znfl, zsfl                 ! local scalars
+      CHARACTER(LEN=1), DIMENSION(jpbgrd)     ::   cgrid
+      REAL(wp)        , DIMENSION(jpi,jpj)    ::   ztmp
+      !!----------------------------------------------------------------------
+
+      cgrid = (/'t','u','v'/)
+
       DO ib_bdy = 1, nb_bdy       ! Indices and directions of rim velocity components
-         icount = 0 
 
          ! Calculate relationship of U direction to the local orientation of the boundary
@@ -994,19 +1153,23 @@
          ! flagu =  0 : u is tangential
          ! flagu =  1 : u is normal to the boundary and is direction is inward
-  
          DO igrd = 1, jpbgrd 
             SELECT CASE( igrd )
-               CASE( 1 )   ;   pmask => umask   (:,:,1)   ;   i_offset = 0
-               CASE( 2 )   ;   pmask => bdytmask(:,:)     ;   i_offset = 1
-               CASE( 3 )   ;   pmask => zfmask  (:,:)     ;   i_offset = 0
+               CASE( 1 )   ;   pmask => zumask     ;   i_offset = 0
+               CASE( 2 )   ;   pmask => bdytmask   ;   i_offset = 1
+               CASE( 3 )   ;   pmask => zfmask     ;   i_offset = 0
             END SELECT 
             icount = 0
             ztmp(:,:) = 0._wp
-            DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)  
+            IF( lrim0 ) THEN   ! extent of rim 0
+               ibeg = 1                                     ;   iend = idx_bdy(ib_bdy)%nblenrim0(igrd)
+            ELSE               ! extent of rim 1
+               ibeg = idx_bdy(ib_bdy)%nblenrim0(igrd) + 1   ;   iend = idx_bdy(ib_bdy)%nblenrim(igrd)
+            END IF
+            DO ib = ibeg, iend 
                ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
                ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
                IF( ii == 1 .OR. ii == jpi .OR. ij == 1 .OR. ij == jpj )  CYCLE
-               zefl = pmask(ii+i_offset-1,ij)
-               zwfl = pmask(ii+i_offset  ,ij)
+               zwfl = pmask(ii+i_offset-1,ij)
+               zefl = pmask(ii+i_offset  ,ij)
                ! This error check only works if you are using the bdyXmask arrays
                IF( i_offset == 1 .and. zefl + zwfl == 2. ) THEN
@@ -1014,5 +1177,5 @@
                   IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(ii),mjg(ij)
                ELSE
-                  ztmp(ii,ij) = -zefl + zwfl
+                  ztmp(ii,ij) = -zwfl + zefl
                ENDIF
             END DO
@@ -1024,5 +1187,5 @@
             ENDIF 
             CALL lbc_lnk( 'bdyini', ztmp, 'T', 1. ) 
-            DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)  
+            DO ib = ibeg, iend
                ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
                ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
@@ -1035,19 +1198,23 @@
          ! flagv =  0 : v is tangential
          ! flagv =  1 : v is normal to the boundary and is direction is inward
-
          DO igrd = 1, jpbgrd 
             SELECT CASE( igrd )
-               CASE( 1 )   ;   pmask => vmask (:,:,1)   ;   j_offset = 0
-               CASE( 2 )   ;   pmask => zfmask(:,:)     ;   j_offset = 0
-               CASE( 3 )   ;   pmask => bdytmask        ;   j_offset = 1
+               CASE( 1 )   ;   pmask => zvmask     ;   j_offset = 0
+               CASE( 2 )   ;   pmask => zfmask     ;   j_offset = 0
+               CASE( 3 )   ;   pmask => bdytmask   ;   j_offset = 1
             END SELECT 
             icount = 0
             ztmp(:,:) = 0._wp
-            DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)  
+            IF( lrim0 ) THEN   ! extent of rim 0
+               ibeg = 1                                     ;   iend = idx_bdy(ib_bdy)%nblenrim0(igrd)
+            ELSE               ! extent of rim 1
+               ibeg = idx_bdy(ib_bdy)%nblenrim0(igrd) + 1   ;   iend = idx_bdy(ib_bdy)%nblenrim(igrd)
+            END IF
+            DO ib = ibeg, iend
                ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
                ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
                IF( ii == 1 .OR. ii == jpi .OR. ij == 1 .OR. ij == jpj )  CYCLE
-               znfl = pmask(ii,ij+j_offset-1)
-               zsfl = pmask(ii,ij+j_offset  )
+               zsfl = pmask(ii,ij+j_offset-1)
+               znfl = pmask(ii,ij+j_offset  )
                ! This error check only works if you are using the bdyXmask arrays
                IF( j_offset == 1 .and. znfl + zsfl == 2. ) THEN
@@ -1055,5 +1222,5 @@
                   icount = icount + 1
                ELSE
-                  ztmp(ii,ij) = -znfl + zsfl
+                  ztmp(ii,ij) = -zsfl + znfl
                END IF
             END DO
@@ -1065,5 +1232,5 @@
             ENDIF
             CALL lbc_lnk( 'bdyini', ztmp, 'T', 1. ) 
-            DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)  
+            DO ib = ibeg, iend
                ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
                ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
@@ -1072,5 +1239,5 @@
          END DO
          !
-      END DO
+      END DO ! ib_bdy
       
       DO ib_bdy = 1, nb_bdy
@@ -1082,12 +1249,17 @@
             END SELECT
             ztmp(:,:) = 0._wp
-            DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)
+            IF( lrim0 ) THEN   ! extent of rim 0
+               ibeg = 1                                     ;   iend = idx_bdy(ib_bdy)%nblenrim0(igrd)
+            ELSE               ! extent of rim 1
+               ibeg = idx_bdy(ib_bdy)%nblenrim0(igrd) + 1   ;   iend = idx_bdy(ib_bdy)%nblenrim(igrd)
+            END IF
+            DO ib = ibeg, iend
                ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
                ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
-               IF( ii == 1 .OR. ii == jpi .OR. ij == 1 .OR. ij == jpj )   CYCLE 
-               llean = pmask(ii+1,ij  ) == 1.
-               llwen = pmask(ii-1,ij  ) == 1.
-               llnon = pmask(ii  ,ij+1) == 1.
-               llson = pmask(ii  ,ij-1) == 1.
+               IF( ii == 1 .OR. ii == jpi .OR. ij == 1 .OR. ij == jpj )   CYCLE
+               llnon = pmask(ii  ,ij+1) == 1.  
+               llson = pmask(ii  ,ij-1) == 1. 
+               llean = pmask(ii+1,ij  ) == 1. 
+               llwen = pmask(ii-1,ij  ) == 1. 
                inn  = COUNT( (/ llnon, llson, llean, llwen /) )
                IF( inn == 0 ) THEN   ! no neighbours -> interior of a corner  or  cluster of rim points
@@ -1101,6 +1273,10 @@
                   ELSE
                      WRITE(ctmp1,*) ' E R R O R : Problem with  ',cgrid(igrd) ,' grid point', ii, ij,   &
-                          ' on boundary set ', ib_bdy, ' has no free ocean neighbour'
-                     WRITE(ctmp2,*) ' There seems to be a cluster of rim points.'
+                       ' on boundary set ', ib_bdy, ' has no free ocean neighbour'
+                     IF( lrim0 ) THEN
+                        WRITE(ctmp2,*) ' There seems to be a cluster of rim 0 points.'
+                     ELSE
+                        WRITE(ctmp2,*) ' There seems to be a cluster of rim 1 points.'
+                     END IF
                      WRITE(ctmp3,*) ' ========== '
                      CALL ctl_stop( ' ', ctmp1, ctmp2, ctmp3, ' ' )
@@ -1142,5 +1318,5 @@
             END DO
             CALL lbc_lnk( 'bdyini', ztmp, 'T', 1. )
-            DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)
+            DO ib = ibeg, iend
                ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
                ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
@@ -1150,138 +1326,53 @@
       END DO
 
-
-      ALLOCATE( lsend_bdyint(nb_bdy,jpbgrd,4), lrecv_bdyint(nb_bdy,jpbgrd,4) )
-      lsend_bdyint(:,:,:) = .false.
-      lrecv_bdyint(:,:,:) = .false. 
-      ALLOCATE( lsend_bdyext(nb_bdy,jpbgrd,4), lrecv_bdyext(nb_bdy,jpbgrd,4) )
-      lsend_bdyext(:,:,:) = .false.
-      lrecv_bdyext(:,:,:) = .false.
-      !
-      ! Check which boundaries might need communication
-      DO igrd = 1, jpbgrd
-         DO ib_bdy = 1, nb_bdy
-            DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)
-               ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
-               ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
-               flagu = NINT(idx_bdy(ib_bdy)%flagu(ib,igrd))
-               flagv = NINT(idx_bdy(ib_bdy)%flagv(ib,igrd))
-               iibe = ii - flagu   ! neighbouring point towards the exterior of the computational domain
-               ijbe = ij - flagv
-               iibi = ii + flagu   ! neighbouring point towards the interior of the computational domain
-               ijbi = ij + flagv
-               CALL find_neib( ii, ij, idx_bdy(ib_bdy)%ntreat(ib,igrd), ii1, ij1, ii2, ij2, ii3, ij3 )
-               !
-               ! search neighbour in the  west/east  direction
-               ! Rim is on the halo and computed ocean is towards exterior of mpi domain  
-               !      <--    (o exterior)     -->  
-               ! (1)  o|x         OR    (2)   x|o
-               !       |___                 ___| 
-               IF( iibi == 0     .OR. ii1 == 0     .OR. ii2 == 0     .OR. ii3 == 0     )   lrecv_bdyint(ib_bdy,igrd,1) = .true.
-               IF( iibi == jpi+1 .OR. ii1 == jpi+1 .OR. ii2 == jpi+1 .OR. ii3 == jpi+1 )   lrecv_bdyint(ib_bdy,igrd,2) = .true.  
-               IF( iibe == 0     )   lrecv_bdyext(ib_bdy,igrd,1) = .true.
-               IF( iibe == jpi+1 )   lrecv_bdyext(ib_bdy,igrd,2) = .true.  
-               ! Check if neighbour has its rim parallel to its mpi subdomain border and located next to its halo
-               ! :¨¨¨¨¨|¨¨-->    |                                             |    <--¨¨|¨¨¨¨¨: 
-               ! :     |  x:o    |    neighbour limited by ... would need o    |    o:x  |     :
-               ! :.....|_._:_____|   (1) W neighbour         E neighbour (2)   |_____:_._|.....:
-               IF( ii == 2     .AND. ( nbondi ==  1 .OR. nbondi == 0 ) .AND. &
-                    & ( iibi == 3     .OR. ii1 == 3     .OR. ii2 == 3     .OR. ii3 == 3    ) )   lsend_bdyint(ib_bdy,igrd,1) = .true.
-               IF( ii == jpi-1 .AND. ( nbondi == -1 .OR. nbondi == 0 ) .AND. &
-                    & ( iibi == jpi-2 .OR. ii1 == jpi-2 .OR. ii2 == jpi-2 .OR. ii3 == jpi-2) )   lsend_bdyint(ib_bdy,igrd,2) = .true.
-               IF( ii == 2     .AND. ( nbondi ==  1 .OR. nbondi == 0 ) .AND. iibe == 3     )   lsend_bdyext(ib_bdy,igrd,1) = .true.
-               IF( ii == jpi-1 .AND. ( nbondi == -1 .OR. nbondi == 0 ) .AND. iibe == jpi-2 )   lsend_bdyext(ib_bdy,igrd,2) = .true.
-               !
-               ! search neighbour in the north/south direction   
-               ! Rim is on the halo and computed ocean is towards exterior of mpi domain
-               !(3)   |       |         ^   ___o___     
-               !  |   |___x___|   OR    |  |   x   |
-               !  v       o           (4)  |       |
-               IF( ijbi == 0     .OR. ij1 == 0     .OR. ij2 == 0     .OR. ij3 == 0     )   lrecv_bdyint(ib_bdy,igrd,3) = .true.
-               IF( ijbi == jpj+1 .OR. ij1 == jpj+1 .OR. ij2 == jpj+1 .OR. ij3 == jpj+1 )   lrecv_bdyint(ib_bdy,igrd,4) = .true.
-               IF( ijbe == 0     )   lrecv_bdyext(ib_bdy,igrd,3) = .true.
-               IF( ijbe == jpj+1 )   lrecv_bdyext(ib_bdy,igrd,4) = .true.
-               ! Check if neighbour has its rim parallel to its mpi subdomain     _________  border and next to its halo
-               !   ^  |    o    |                                                :         : 
-               !   |  |¨¨¨¨x¨¨¨¨|   neighbour limited by ... would need o     |  |....x....|
-               !      :_________:  (3) S neighbour          N neighbour (4)   v  |    o    |   
-               IF( ij == 2     .AND. ( nbondj ==  1 .OR. nbondj == 0 ) .AND. &
-                    & ( ijbi == 3     .OR. ij1 == 3     .OR. ij2 == 3     .OR. ij3 == 3    ) )   lsend_bdyint(ib_bdy,igrd,3) = .true.
-               IF( ij == jpj-1 .AND. ( nbondj == -1 .OR. nbondj == 0 ) .AND. &
-                    & ( ijbi == jpj-2 .OR. ij1 == jpj-2 .OR. ij2 == jpj-2 .OR. ij3 == jpj-2) )   lsend_bdyint(ib_bdy,igrd,4) = .true.
-               IF( ij == 2     .AND. ( nbondj ==  1 .OR. nbondj == 0 ) .AND. ijbe == 3     )   lsend_bdyext(ib_bdy,igrd,3) = .true.
-               IF( ij == jpj-1 .AND. ( nbondj == -1 .OR. nbondj == 0 ) .AND. ijbe == jpj-2 )   lsend_bdyext(ib_bdy,igrd,4) = .true.
-            END DO
-         END DO
-      END DO
-
-      DO ib_bdy = 1,nb_bdy
-         IF( cn_dyn2d(ib_bdy) == 'orlanski' .OR. cn_dyn2d(ib_bdy) == 'orlanski_npo' .OR. &
-           & cn_dyn3d(ib_bdy) == 'orlanski' .OR. cn_dyn3d(ib_bdy) == 'orlanski_npo' .OR. &
-           & cn_tra(ib_bdy)   == 'orlanski' .OR. cn_tra(ib_bdy)   == 'orlanski_npo'      ) THEN
-            DO igrd = 1, jpbgrd
-               DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)
-                  ii = idx_bdy(ib_bdy)%nbi(ib,igrd)
-                  ij = idx_bdy(ib_bdy)%nbj(ib,igrd)
-                  IF(  mig(ii) > 2 .AND. mig(ii) < jpiglo-2 .AND. mjg(ij) > 2 .AND. mjg(ij) < jpjglo-2  ) THEN
-                     WRITE(ctmp1,*) ' Orlanski can not be used when the open boundaries are on the interior of the computational domain'
-                     WRITE(ctmp2,*) ' ========== '
-                     CALL ctl_stop( ' ', ctmp1, ctmp2, ' ' )
-                  END IF
-               END DO
-            END DO
-         END IF
-      END DO
-
-      !
-      ! Tidy up
-      !--------
-      IF( nb_bdy>0 )   DEALLOCATE( nbidta, nbjdta, nbrdta )
-      !
-   END SUBROUTINE bdy_segs
-
-
-   SUBROUTINE find_neib( ii, ij, itreat, ii1, ij1, ii2, ij2, ii3, ij3 )
-     !!----------------------------------------------------------------------
-     !!                 ***  ROUTINE find_neib  ***
-     !!
-     !! ** Purpose :   get ii1, ij1, ii2, ij2, ii3, ij3, the indices of
-     !!               the free ocean neighbours of (ii,ij) for bdy treatment
-     !!
-     !!----------------------------------------------------------------------
-     INTEGER, INTENT(in   )      ::   ii, ij, itreat
-     INTEGER, INTENT(  out)      ::   ii1, ij1, ii2, ij2, ii3, ij3
-     !!----------------------------------------------------------------------
-     SELECT CASE( itreat )   ! points that will be used by bdy routines, -1 will be discarded
-        !               !               !     _____     !     _____     
-        !  1 |   o      !  2  o   |     !  3 | x        !  4     x |    
-        !    |_x_ _     !    _ _x_|     !    |   o      !      o   |
-     CASE( 1 )    ;   ii1 = ii+1   ;   ij1 = ij+1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
-     CASE( 2 )    ;   ii1 = ii-1   ;   ij1 = ij+1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
-     CASE( 3 )    ;   ii1 = ii+1   ;   ij1 = ij-1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
-     CASE( 4 )    ;   ii1 = ii-1   ;   ij1 = ij-1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
-        !    |          !         |     !      o        !    ______                   ! or incomplete corner
-        ! 5  | x o      ! 6   o x |     ! 7  __x__      ! 8    x                      !  7  ____ o
-        !    |          !         |     !               !      o                      !         |x___
-     CASE( 5 )    ;   ii1 = ii+1   ;   ij1 = ij     ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
-     CASE( 6 )    ;   ii1 = ii-1   ;   ij1 = ij     ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
-     CASE( 7 )    ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
-     CASE( 8 )    ;   ii1 = ii     ;   ij1 = ij-1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
-        !        o      !        o      !    _____|     !       |_____  
-        !  9 ____x o    ! 10   o x___   ! 11     x o    ! 12   o x      
-        !         |     !       |       !        o      !        o      
-     CASE( 9  )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii+1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1 
-     CASE( 10 )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii-1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1
-     CASE( 11 )   ;   ii1 = ii     ;   ij1 = ij-1   ;   ii2 = ii+1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1
-     CASE( 12 )   ;   ii1 = ii     ;   ij1 = ij-1   ;   ii2 = ii-1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1
-        !    |_  o      !        o  _|  !     ¨¨|_|¨¨   !       o         
-        ! 13  _| x o    !  14  o x |_   !  15  o x o    ! 16  o x o       
-        !    |   o      !        o   |  !        o      !    __|¨|__ 
-     CASE( 13 )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii+1   ;   ij2 = ij     ;   ii3 = ii     ;   ij3 = ij-1   
-     CASE( 14 )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii-1   ;   ij2 = ij     ;   ii3 = ii     ;   ij3 = ij-1 
-     CASE( 15 )   ;   ii1 = ii-1   ;   ij1 = ij     ;   ii2 = ii     ;   ij2 = ij-1   ;   ii3 = ii+1   ;   ij3 = ij   
-     CASE( 16 )   ;   ii1 = ii-1   ;   ij1 = ij     ;   ii2 = ii     ;   ij2 = ij+1   ;   ii3 = ii+1   ;   ij3 = ij
-     END SELECT
-   END SUBROUTINE find_neib
-
+    END SUBROUTINE bdy_rim_treat
+
+   
+    SUBROUTINE find_neib( ii, ij, itreat, ii1, ij1, ii2, ij2, ii3, ij3 )
+      !!----------------------------------------------------------------------
+      !!                 ***  ROUTINE find_neib  ***
+      !!
+      !! ** Purpose :   get ii1, ij1, ii2, ij2, ii3, ij3, the indices of
+      !!               the free ocean neighbours of (ii,ij) for bdy treatment
+      !!
+      !! ** Method  :  use itreat input to select a case
+      !!               N.B. ntreat is defined for all bdy points in routine bdy_segs
+      !!
+      !!----------------------------------------------------------------------
+      INTEGER, INTENT(in   )      ::   ii, ij, itreat
+      INTEGER, INTENT(  out)      ::   ii1, ij1, ii2, ij2, ii3, ij3
+      !!----------------------------------------------------------------------
+      SELECT CASE( itreat )   ! points that will be used by bdy routines, -1 will be discarded
+         !               !               !     _____     !     _____     
+         !  1 |   o      !  2  o   |     !  3 | x        !  4     x |    
+         !    |_x_ _     !    _ _x_|     !    |   o      !      o   |
+      CASE( 1 )    ;   ii1 = ii+1   ;   ij1 = ij+1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
+      CASE( 2 )    ;   ii1 = ii-1   ;   ij1 = ij+1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
+      CASE( 3 )    ;   ii1 = ii+1   ;   ij1 = ij-1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
+      CASE( 4 )    ;   ii1 = ii-1   ;   ij1 = ij-1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
+         !    |          !         |     !      o        !    ______                   ! or incomplete corner
+         ! 5  | x o      ! 6   o x |     ! 7  __x__      ! 8    x                      !  7  ____ o
+         !    |          !         |     !               !      o                      !         |x___
+      CASE( 5 )    ;   ii1 = ii+1   ;   ij1 = ij     ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
+      CASE( 6 )    ;   ii1 = ii-1   ;   ij1 = ij     ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
+      CASE( 7 )    ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
+      CASE( 8 )    ;   ii1 = ii     ;   ij1 = ij-1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1
+         !        o      !        o      !    _____|     !       |_____  
+         !  9 ____x o    ! 10   o x___   ! 11     x o    ! 12   o x      
+         !         |     !       |       !        o      !        o      
+      CASE( 9  )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii+1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1 
+      CASE( 10 )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii-1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1
+      CASE( 11 )   ;   ii1 = ii     ;   ij1 = ij-1   ;   ii2 = ii+1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1
+      CASE( 12 )   ;   ii1 = ii     ;   ij1 = ij-1   ;   ii2 = ii-1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1
+         !    |_  o      !        o  _|  !     ¨¨|_|¨¨   !       o         
+         ! 13  _| x o    !  14  o x |_   !  15  o x o    ! 16  o x o       
+         !    |   o      !        o   |  !        o      !    __|¨|__ 
+      CASE( 13 )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii+1   ;   ij2 = ij     ;   ii3 = ii     ;   ij3 = ij-1   
+      CASE( 14 )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii-1   ;   ij2 = ij     ;   ii3 = ii     ;   ij3 = ij-1 
+      CASE( 15 )   ;   ii1 = ii-1   ;   ij1 = ij     ;   ii2 = ii     ;   ij2 = ij-1   ;   ii3 = ii+1   ;   ij3 = ij   
+      CASE( 16 )   ;   ii1 = ii-1   ;   ij1 = ij     ;   ii2 = ii     ;   ij2 = ij+1   ;   ii3 = ii+1   ;   ij3 = ij
+      END SELECT
+    END SUBROUTINE find_neib
+    
 
    SUBROUTINE bdy_ctl_seg
Index: /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdylib.F90
===================================================================
--- /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdylib.F90	(revision 11190)
+++ /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdylib.F90	(revision 11191)
@@ -92,5 +92,5 @@
 
 
-   SUBROUTINE bdy_orl( idx, ptb, pta, dta, ll_npo )
+   SUBROUTINE bdy_orl( idx, ptb, pta, dta, lrim0, ll_npo )
       !!----------------------------------------------------------------------
       !!                 ***  SUBROUTINE bdy_orl  ***
@@ -104,4 +104,5 @@
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptb  ! before tracer field
       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pta  ! tracer trend
+      LOGICAL                 , OPTIONAL,  INTENT(in) ::   lrim0   ! indicate if rim 0 is treated
       LOGICAL,                             INTENT(in) ::   ll_npo  ! switch for NPO version
       !!
@@ -111,10 +112,10 @@
       igrd = 1                       ! Everything is at T-points here
       !
-      CALL bdy_orlanski_3d( idx, igrd, ptb(:,:,:), pta(:,:,:), dta, ll_npo )
+      CALL bdy_orlanski_3d( idx, igrd, ptb(:,:,:), pta(:,:,:), dta, lrim0, ll_npo )
       !
    END SUBROUTINE bdy_orl
 
 
-   SUBROUTINE bdy_orlanski_2d( idx, igrd, phib, phia, phi_ext, ll_npo )
+   SUBROUTINE bdy_orlanski_2d( idx, igrd, phib, phia, phi_ext, lrim0, ll_npo )
       !!----------------------------------------------------------------------
       !!                 ***  SUBROUTINE bdy_orlanski_2d  ***
@@ -132,4 +133,5 @@
       REAL(wp), DIMENSION(:,:), INTENT(inout) ::   phia     ! model after 2D field (to be updated)
       REAL(wp), DIMENSION(:)  , INTENT(in   ) ::   phi_ext  ! external forcing data
+      LOGICAL, OPTIONAL,        INTENT(in   ) ::   lrim0    ! indicate if rim 0 is treated
       LOGICAL ,                 INTENT(in   ) ::   ll_npo   ! switch for NPO version
       !
@@ -140,4 +142,5 @@
       INTEGER  ::   ii_offset, ij_offset                   ! offsets for mask indices
       INTEGER  ::   flagu, flagv                           ! short cuts
+      INTEGER  ::   ibeg, iend                             ! length of rim to be treated (rim 0 or rim 1 or both)
       REAL(wp) ::   zmask_x, zmask_y1, zmask_y2
       REAL(wp) ::   zex1, zex2, zey, zey1, zey2
@@ -185,5 +188,12 @@
       END SELECT
       !
-      DO jb = 1, idx%nblenrim(igrd)
+      IF( PRESENT(lrim0) ) THEN
+         IF( lrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd)   ! rim 0
+         ELSE               ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd)    ! rim 1
+         END IF
+      ELSE                  ;   ibeg = 1                       ;   iend = idx%nblenrim(igrd)    ! both
+      END IF
+      !
+      DO jb = ibeg, iend
          ii  = idx%nbi(jb,igrd)
          ij  = idx%nbj(jb,igrd) 
@@ -272,5 +282,5 @@
 
 
-   SUBROUTINE bdy_orlanski_3d( idx, igrd, phib, phia, phi_ext, ll_npo )
+   SUBROUTINE bdy_orlanski_3d( idx, igrd, phib, phia, phi_ext, lrim0, ll_npo )
       !!----------------------------------------------------------------------
       !!                 ***  SUBROUTINE bdy_orlanski_3d  ***
@@ -288,4 +298,5 @@
       REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   phia     ! model after 3D field (to be updated)
       REAL(wp), DIMENSION(:,:)  , INTENT(in   ) ::   phi_ext  ! external forcing data
+      LOGICAL, OPTIONAL,          INTENT(in   ) ::   lrim0    ! indicate if rim 0 is treated
       LOGICAL ,                   INTENT(in   ) ::   ll_npo   ! switch for NPO version
       !
@@ -296,4 +307,5 @@
       INTEGER  ::   ii_offset, ij_offset                   ! offsets for mask indices
       INTEGER  ::   flagu, flagv                           ! short cuts
+      INTEGER  ::   ibeg, iend                             ! length of rim to be treated (rim 0 or rim 1 or both)
       REAL(wp) ::   zmask_x, zmask_y1, zmask_y2
       REAL(wp) ::   zex1, zex2, zey, zey1, zey2
@@ -340,8 +352,15 @@
          CASE DEFAULT ;   CALL ctl_stop( 'unrecognised value for igrd in bdy_orlanksi_2d' )
       END SELECT
-
+      !
+      IF( PRESENT(lrim0) ) THEN
+         IF( lrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd)   ! rim 0
+         ELSE               ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd)    ! rim 1
+         END IF
+      ELSE                  ;   ibeg = 1                       ;   iend = idx%nblenrim(igrd)    ! both
+      END IF
+      !
       DO jk = 1, jpk
          !            
-         DO jb = 1, idx%nblenrim(igrd)
+         DO jb = ibeg, iend
             ii  = idx%nbi(jb,igrd)
             ij  = idx%nbj(jb,igrd) 
@@ -430,5 +449,5 @@
    END SUBROUTINE bdy_orlanski_3d
 
-   SUBROUTINE bdy_nmn( idx, igrd, phia )
+   SUBROUTINE bdy_nmn( idx, igrd, phia, lrim0 )
       !!----------------------------------------------------------------------
       !!                 ***  SUBROUTINE bdy_nmn  ***
@@ -444,7 +463,8 @@
       !!                                                   !      o      
       !!----------------------------------------------------------------------
-      INTEGER,                    INTENT(in)     ::   igrd     ! grid index
+      INTEGER,                    INTENT(in   )  ::   igrd     ! grid index
       REAL(wp), DIMENSION(:,:,:), INTENT(inout)  ::   phia     ! model after 3D field (to be updated), must be masked
-      TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices
+      TYPE(OBC_INDEX),            INTENT(in   )  ::   idx      ! OBC indices
+      LOGICAL, OPTIONAL,          INTENT(in   )  ::   lrim0    ! indicate if rim 0 is treated
       !! 
       REAL(wp) ::   zweight
@@ -453,4 +473,5 @@
       INTEGER  ::   ii, ij   ! 2D addresses
       INTEGER  ::   ipkm1    ! size of phia third dimension minus 1
+      INTEGER  ::   ibeg, iend                          ! length of rim to be treated (rim 0 or rim 1 or both)
       INTEGER  ::   ii1, ii2, ii3, ij1, ij2, ij3, itreat
       !!----------------------------------------------------------------------
@@ -465,9 +486,16 @@
       END SELECT
       !
-      DO ib = 1, idx%nblenrim(igrd)
+      IF( PRESENT(lrim0) ) THEN
+         IF( lrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd)   ! rim 0
+         ELSE               ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd)    ! rim 1
+         END IF
+      ELSE                  ;   ibeg = 1                       ;   iend = idx%nblenrim(igrd)    ! both
+      END IF
+      !
+      DO ib = ibeg, iend
          ii = idx%nbi(ib,igrd)
          ij = idx%nbj(ib,igrd)
          itreat = idx%ntreat(ib,igrd)
-         CALL find_neib( ii, ij, itreat, ii1, ij1, ii2, ij2, ii3, ij3 )
+         CALL find_neib( ii, ij, itreat, ii1, ij1, ii2, ij2, ii3, ij3 )   ! find free ocean neighbours
          SELECT CASE( itreat )
          CASE( 1:8 )
Index: /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdytra.F90
===================================================================
--- /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdytra.F90	(revision 11190)
+++ /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdytra.F90	(revision 11191)
@@ -49,51 +49,61 @@
       INTEGER, INTENT(in) ::   kt   ! Main time step counter
       !
-      INTEGER                        :: ib_bdy, jn, igrd   ! Loop indeces
-      TYPE(ztrabdy), DIMENSION(jpts) :: zdta               ! Temporary data structure
-      LOGICAL, DIMENSION(4)          :: llsend1, llrecv1     ! indicate how communications are to be carried out
+      INTEGER                        :: ib_bdy, jn, igrd, ir   ! Loop indeces
+      TYPE(ztrabdy), DIMENSION(jpts) :: zdta                   ! Temporary data structure
+      LOGICAL                        :: llrim0                 ! indicate if rim 0 is treated
+      LOGICAL, DIMENSION(4)          :: llsend1, llrecv1       ! indicate how communications are to be carried out
       !!----------------------------------------------------------------------
       igrd = 1 
 
-      DO ib_bdy=1, nb_bdy
+      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
+            !
+            zdta(1)%tra => dta_bdy(ib_bdy)%tem
+            zdta(2)%tra => dta_bdy(ib_bdy)%sal
+            !
+            DO jn = 1, jpts
+               !
+               SELECT CASE( TRIM(cn_tra(ib_bdy)) )
+               CASE('none'        )   ;   CYCLE
+               CASE('frs'         )   ! treat the whole boundary at once
+                  IF( ir == 0 ) CALL bdy_frs ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), zdta(jn)%tra )
+               CASE('specified'   )   ! treat the whole rim      at once
+                  IF( ir == 0 ) CALL bdy_spe ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), zdta(jn)%tra )
+               CASE('neumann'     )   ;   CALL bdy_nmn ( idx_bdy(ib_bdy), igrd         , tsa(:,:,:,jn), llrim0 )   ! tsa masked
+               CASE('orlanski'    )   ;   CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), &
+                    & zdta(jn)%tra, llrim0, ll_npo=.false. )
+               CASE('orlanski_npo')   ;   CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), &
+                    & zdta(jn)%tra, llrim0, ll_npo=.true.  )
+               CASE('runoff'      )   ;   CALL bdy_rnf ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), jn, llrim0 )
+               CASE DEFAULT           ;   CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' )
+               END SELECT
+               ! 
+            END DO
+         END DO
          !
-         zdta(1)%tra => dta_bdy(ib_bdy)%tem
-         zdta(2)%tra => dta_bdy(ib_bdy)%sal
-         !
-         DO jn = 1, jpts
-            !
+         llsend1(:) = .false.
+         llrecv1(:) = .false.
+         DO ib_bdy=1, nb_bdy
             SELECT CASE( TRIM(cn_tra(ib_bdy)) )
-            CASE('none'        )   ;   CYCLE
-            CASE('frs'         )   ;   CALL bdy_frs ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), zdta(jn)%tra )
-            CASE('specified'   )   ;   CALL bdy_spe ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), zdta(jn)%tra )
-            CASE('neumann'     )   ;   CALL bdy_nmn ( idx_bdy(ib_bdy), igrd         , tsa(:,:,:,jn)               )   ! tsa masked
-            CASE('orlanski'    )   ;   CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), zdta(jn)%tra, ll_npo=.false. )
-            CASE('orlanski_npo')   ;   CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), zdta(jn)%tra, ll_npo=.true. )
-            CASE('runoff'      )   ;   CALL bdy_rnf ( idx_bdy(ib_bdy),                tsa(:,:,:,jn),               jn )
-            CASE DEFAULT           ;   CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' )
+            CASE('neumann','runoff')
+               llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:,ir)   ! possibly every direction, T points
+               llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:,ir)   ! possibly every direction, T points
+            CASE('orlanski', 'orlanski_npo')
+               llsend1(:) = llsend1(:) .OR. lsend_bdy(ib_bdy,1,:,ir)   ! possibly every direction, T points
+               llrecv1(:) = llrecv1(:) .OR. lrecv_bdy(ib_bdy,1,:,ir)   ! possibly every direction, T points
             END SELECT
-            ! 
          END DO
+         IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction
+            CALL lbc_bdy_lnk( 'bdytra', llsend1, llrecv1, tsa, 'T',  1. )
+         END IF
       END DO
-      !
-      llsend1(:) = .false.
-      llrecv1(:) = .false.
-      DO ib_bdy=1, nb_bdy
-         SELECT CASE( TRIM(cn_tra(ib_bdy)) )
-         CASE('neumann','runoff')
-            llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:)   ! possibly every direction, T points
-            llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:)   ! possibly every direction, T points
-         CASE('orlanski', 'orlanski_npo')
-            llsend1(:) = llsend1(:) .OR. lsend_bdy(ib_bdy,1,:)   ! possibly every direction, T points
-            llrecv1(:) = llrecv1(:) .OR. lrecv_bdy(ib_bdy,1,:)   ! possibly every direction, T points
-         END SELECT
-      END DO
-      IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction
-         CALL lbc_bdy_lnk( 'bdytra', llsend1, llrecv1, tsa, 'T',  1. )
-      END IF
       !
    END SUBROUTINE bdy_tra
 
 
-   SUBROUTINE bdy_rnf( idx, pta, jpa )
+   SUBROUTINE bdy_rnf( idx, pta, jpa, llrim0 )
       !!----------------------------------------------------------------------
       !!                 ***  SUBROUTINE bdy_rnf  ***
@@ -104,7 +114,8 @@
       !! 
       !!----------------------------------------------------------------------
-      TYPE(OBC_INDEX),                     INTENT(in) ::   idx  ! OBC indices
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pta  ! tracer trend
-      INTEGER,                             INTENT(in) ::   jpa  ! TRA index
+      TYPE(OBC_INDEX),                     INTENT(in) ::   idx      ! OBC indices
+      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pta      ! tracer trend
+      INTEGER,                             INTENT(in) ::   jpa      ! TRA index
+      LOGICAL,                             INTENT(in) ::   llrim0   ! indicate if rim 0 is treated
       !
       INTEGER  ::   ib, ii, ij, igrd   ! dummy loop indices
@@ -113,8 +124,9 @@
       !
       igrd = 1                       ! Everything is at T-points here
-      IF(      jpa == jp_tem )   THEN
-         CALL bdy_nmn( idx, igrd, pta )
-      ELSE IF( jpa == jp_sal )   THEN
-         DO ib = 1, idx%nblenrim(igrd)
+      IF(      jpa == jp_tem ) THEN
+         CALL bdy_nmn( idx, igrd, pta, llrim0 )
+      ELSE IF( jpa == jp_sal ) THEN
+         IF( .NOT. llrim0 )   RETURN
+         DO ib = 1, idx%nblenrim(igrd)   ! if llrim0 then treat the whole rim
             ii = idx%nbi(ib,igrd)
             ij = idx%nbj(ib,igrd)
Index: /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/DYN/dynkeg.F90
===================================================================
--- /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/DYN/dynkeg.F90	(revision 11190)
+++ /NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/DYN/dynkeg.F90	(revision 11191)
@@ -140,8 +140,8 @@
             llrecv1(:) = .false.
             DO ib_bdy = 1, nb_bdy
-               llsend1(2) = llsend1(2) .OR. lsend_bdy(ib_bdy,igrd,2)   ! send east
-               llsend1(4) = llsend1(4) .OR. lsend_bdy(ib_bdy,igrd,4)   ! send north
-               llrecv1(1) = llrecv1(1) .OR. lrecv_bdy(ib_bdy,igrd,1)   ! receive west 
-               llrecv1(3) = llrecv1(3) .OR. lrecv_bdy(ib_bdy,igrd,3)   ! receive south
+               llsend1(2) = llsend1(2) .OR. ANY(lsend_bdy(ib_bdy,igrd,2,:))   ! send east
+               llsend1(4) = llsend1(4) .OR. ANY(lsend_bdy(ib_bdy,igrd,4,:))   ! send north
+               llrecv1(1) = llrecv1(1) .OR. ANY(lrecv_bdy(ib_bdy,igrd,1,:))   ! receive west 
+               llrecv1(3) = llrecv1(3) .OR. ANY(lrecv_bdy(ib_bdy,igrd,3,:))   ! receive south
             END DO
    
