Changeset 11191


Ignore:
Timestamp:
2019-06-27T10:14:39+02:00 (14 months ago)
Author:
girrmann
Message:

dev_r10984_HPC-13 : bdy treatment can now handel a rim 0 and a rim 1, results are unchanged when only rim 1 is provided, see #2288 and #2285

Location:
NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdy_oce.F90

    r11071 r11191  
    2222      INTEGER ,          DIMENSION(jpbgrd) ::  nblen 
    2323      INTEGER ,          DIMENSION(jpbgrd) ::  nblenrim 
     24      INTEGER ,          DIMENSION(jpbgrd) ::  nblenrim0 
    2425      INTEGER , POINTER, DIMENSION(:,:)    ::  nbi 
    2526      INTEGER , POINTER, DIMENSION(:,:)    ::  nbj 
     
    139140   TYPE(OBC_DATA) , DIMENSION(jp_bdy), TARGET      ::   dta_bdy           !: bdy external data (local process) 
    140141!$AGRIF_END_DO_NOT_TREAT 
    141    LOGICAL, ALLOCATABLE, DIMENSION(:,:,:) ::   lsend_bdy      !: mark needed communication for given boundary, grid and neighbour 
    142    LOGICAL, ALLOCATABLE, DIMENSION(:,:,:) ::   lrecv_bdy      !:  when searching in any direction 
    143    LOGICAL, ALLOCATABLE, DIMENSION(:,:,:) ::   lsend_bdyint   !: mark needed communication for given boundary, grid and neighbour 
    144    LOGICAL, ALLOCATABLE, DIMENSION(:,:,:) ::   lrecv_bdyint   !:  when searching towards the interior of the computational domain 
    145    LOGICAL, ALLOCATABLE, DIMENSION(:,:,:) ::   lsend_bdyext   !: mark needed communication for given boundary, grid and neighbour 
    146    LOGICAL, ALLOCATABLE, DIMENSION(:,:,:) ::   lrecv_bdyext   !:  when searching towards the exterior of the computational domain 
     142   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) ::   lsend_bdy      !: mark needed communication for given boundary, grid and neighbour 
     143   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) ::   lrecv_bdy      !:  when searching in any direction 
     144   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) ::   lsend_bdyint   !: mark needed communication for given boundary, grid and neighbour 
     145   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) ::   lrecv_bdyint   !:  when searching towards the interior of the computational domain 
     146   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) ::   lsend_bdyext   !: mark needed communication for given boundary, grid and neighbour 
     147   LOGICAL, ALLOCATABLE, DIMENSION(:,:,:,:) ::   lrecv_bdyext   !:  when searching towards the exterior of the computational domain 
    147148   !!---------------------------------------------------------------------- 
    148149   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
  • NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdydyn2d.F90

    r11149 r11191  
    5050      REAL(wp), DIMENSION(jpi,jpj), INTENT(in   ) :: pssh 
    5151      !! 
    52       INTEGER                                  ::   ib_bdy     ! Loop counter 
     52      INTEGER  ::   ib_bdy, ir     ! BDY set index, rim index 
     53      LOGICAL  ::   llrim0         ! indicate if rim 0 is treated 
    5354      LOGICAL, DIMENSION(4) :: llsend2, llrecv2, llsend3, llrecv3  ! indicate how communications are to be carried out 
    5455       
    55       DO ib_bdy=1, nb_bdy 
    56          SELECT CASE( cn_dyn2d(ib_bdy) ) 
    57          CASE('none') 
    58             CYCLE 
    59          CASE('frs') 
    60             CALL bdy_dyn2d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, pua2d, pva2d ) 
    61          CASE('flather') 
    62             CALL bdy_dyn2d_fla( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, pua2d, pva2d, pssh, phur, phvr ) 
    63          CASE('orlanski') 
    64             CALL bdy_dyn2d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, & 
    65                                      & pua2d, pva2d, pub2d, pvb2d, ll_npo=.false.) 
    66          CASE('orlanski_npo') 
    67             CALL bdy_dyn2d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, & 
    68                                      & pua2d, pva2d, pub2d, pvb2d, ll_npo=.true. ) 
    69          CASE DEFAULT 
    70             CALL ctl_stop( 'bdy_dyn2d : unrecognised option for open boundaries for barotropic variables' ) 
    71          END SELECT 
    72       ENDDO 
    73       ! 
    74       llsend2(:) = .false. 
    75       llrecv2(:) = .false. 
    76       llsend3(:) = .false. 
    77       llrecv3(:) = .false. 
    78       DO ib_bdy=1, nb_bdy 
    79          SELECT CASE( cn_dyn2d(ib_bdy) ) 
    80          CASE('flather') 
    81             llsend2(1:2) = llsend2(1:2) .OR. lsend_bdyint(ib_bdy,2,1:2)   ! west/east, U points 
    82             llsend2(1)   = llsend2(1)   .OR. lsend_bdyext(ib_bdy,2,1)     ! neighbour might search point towards bdy on its east 
    83             llrecv2(1:2) = llrecv2(1:2) .OR. lrecv_bdyint(ib_bdy,2,1:2)   ! west/east, U points 
    84             llrecv2(2)   = llrecv2(2)   .OR. lrecv_bdyext(ib_bdy,2,2)     ! might search point towards bdy on the east 
    85             llsend3(3:4) = llsend3(3:4) .OR. lsend_bdyint(ib_bdy,3,3:4)   ! north/south, V points 
    86             llsend3(3)   = llsend3(3)   .OR. lsend_bdyext(ib_bdy,3,3)     ! neighbour might search point towards bdy on its north 
    87             llrecv3(3:4) = llrecv3(3:4) .OR. lrecv_bdyint(ib_bdy,3,3:4)   ! north/south, V points 
    88             llrecv3(4)   = llrecv3(4)   .OR. lrecv_bdyext(ib_bdy,3,4)     ! might search point towards bdy on the north 
    89          CASE('orlanski', 'orlanski_npo') 
    90             llsend2(:) = llsend2(:) .OR. lsend_bdy(ib_bdy,2,:)   ! possibly every direction, U points 
    91             llrecv2(:) = llrecv2(:) .OR. lrecv_bdy(ib_bdy,2,:)   ! possibly every direction, U points 
    92             llsend3(:) = llsend3(:) .OR. lsend_bdy(ib_bdy,3,:)   ! possibly every direction, V points 
    93             llrecv3(:) = llrecv3(:) .OR. lrecv_bdy(ib_bdy,3,:)   ! possibly every direction, V points 
    94          END SELECT 
    95       END DO 
    96       IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN   ! if need to send/recv in at least one direction 
    97          CALL lbc_bdy_lnk( 'bdydyn2d', llsend2, llrecv2, pua2d, 'U', -1. ) 
    98       END IF 
    99       IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction 
    100          CALL lbc_bdy_lnk( 'bdydyn2d', llsend3, llrecv3, pva2d, 'V', -1. ) 
    101       END IF 
     56      DO ir = 1, 0, -1   ! treat rim 1 before rim 0 
     57         IF( ir == 0 ) THEN   ;   llrim0 = .TRUE. 
     58         ELSE                 ;   llrim0 = .FALSE. 
     59         END IF 
     60         DO ib_bdy=1, nb_bdy 
     61            SELECT CASE( cn_dyn2d(ib_bdy) ) 
     62            CASE('none') 
     63               CYCLE 
     64            CASE('frs')   ! treat the whole boundary at once 
     65               IF( llrim0 )   CALL bdy_dyn2d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, pua2d, pva2d ) 
     66            CASE('flather') 
     67               CALL bdy_dyn2d_fla( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, pua2d, pva2d, pssh, phur, phvr, llrim0 ) 
     68            CASE('orlanski') 
     69               CALL bdy_dyn2d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, & 
     70                    & pua2d, pva2d, pub2d, pvb2d, llrim0, ll_npo=.false. ) 
     71            CASE('orlanski_npo') 
     72               CALL bdy_dyn2d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, & 
     73                    & pua2d, pva2d, pub2d, pvb2d, llrim0, ll_npo=.true.  ) 
     74            CASE DEFAULT 
     75               CALL ctl_stop( 'bdy_dyn2d : unrecognised option for open boundaries for barotropic variables' ) 
     76            END SELECT 
     77         ENDDO 
     78         ! 
     79         llsend2(:) = .false. 
     80         llrecv2(:) = .false. 
     81         llsend3(:) = .false. 
     82         llrecv3(:) = .false. 
     83         DO ib_bdy=1, nb_bdy 
     84            SELECT CASE( cn_dyn2d(ib_bdy) ) 
     85            CASE('flather') 
     86               llsend2(1:2) = llsend2(1:2) .OR. lsend_bdyint(ib_bdy,2,1:2,ir)   ! west/east, U points 
     87               llsend2(1)   = llsend2(1)   .OR. lsend_bdyext(ib_bdy,2,1,ir)     ! neighbour might search point towards its east bdy 
     88               llrecv2(1:2) = llrecv2(1:2) .OR. lrecv_bdyint(ib_bdy,2,1:2,ir)   ! west/east, U points 
     89               llrecv2(2)   = llrecv2(2)   .OR. lrecv_bdyext(ib_bdy,2,2,ir)     ! might search point towards bdy on the east 
     90               llsend3(3:4) = llsend3(3:4) .OR. lsend_bdyint(ib_bdy,3,3:4,ir)   ! north/south, V points 
     91               llsend3(3)   = llsend3(3)   .OR. lsend_bdyext(ib_bdy,3,3,ir)     ! neighbour might search point towards its north bdy  
     92               llrecv3(3:4) = llrecv3(3:4) .OR. lrecv_bdyint(ib_bdy,3,3:4,ir)   ! north/south, V points 
     93               llrecv3(4)   = llrecv3(4)   .OR. lrecv_bdyext(ib_bdy,3,4,ir)     ! might search point towards bdy on the north 
     94            CASE('orlanski', 'orlanski_npo') 
     95               llsend2(:) = llsend2(:) .OR. lsend_bdy(ib_bdy,2,:,ir)   ! possibly every direction, U points 
     96               llrecv2(:) = llrecv2(:) .OR. lrecv_bdy(ib_bdy,2,:,ir)   ! possibly every direction, U points 
     97               llsend3(:) = llsend3(:) .OR. lsend_bdy(ib_bdy,3,:,ir)   ! possibly every direction, V points 
     98               llrecv3(:) = llrecv3(:) .OR. lrecv_bdy(ib_bdy,3,:,ir)   ! possibly every direction, V points 
     99            END SELECT 
     100         END DO 
     101         IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN   ! if need to send/recv in at least one direction 
     102            CALL lbc_bdy_lnk( 'bdydyn2d', llsend2, llrecv2, pua2d, 'U', -1. ) 
     103         END IF 
     104         IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction 
     105            CALL lbc_bdy_lnk( 'bdydyn2d', llsend3, llrecv3, pva2d, 'V', -1. ) 
     106         END IF 
     107         ! 
     108      END DO 
    102109      ! 
    103110   END SUBROUTINE bdy_dyn2d 
     
    119126      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pua2d, pva2d  
    120127      !! 
    121       INTEGER  ::   jb, jk         ! dummy loop indices 
     128      INTEGER  ::   jb             ! dummy loop indices 
    122129      INTEGER  ::   ii, ij, igrd   ! local integers 
    123130      REAL(wp) ::   zwgt           ! boundary weight 
     
    143150 
    144151 
    145    SUBROUTINE bdy_dyn2d_fla( idx, dta, ib_bdy, pua2d, pva2d, pssh, phur, phvr ) 
     152   SUBROUTINE bdy_dyn2d_fla( idx, dta, ib_bdy, pua2d, pva2d, pssh, phur, phvr, llrim0 ) 
    146153      !!---------------------------------------------------------------------- 
    147154      !!                 ***  SUBROUTINE bdy_dyn2d_fla  *** 
     
    166173      INTEGER,                      INTENT(in) ::   ib_bdy  ! BDY set index 
    167174      REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: pua2d, pva2d 
    168       REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   pssh, phur, phvr  
    169  
     175      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   pssh, phur, phvr 
     176      LOGICAL                     , INTENT(in) ::   llrim0   ! indicate if rim 0 is treated 
     177      INTEGER  ::   ibeg, iend                       ! length of rim to be treated (rim 0 or rim 1) 
    170178      INTEGER  ::   jb, igrd                         ! dummy loop indices 
    171179      INTEGER  ::   ii, ij                           ! 2D addresses 
    172180      INTEGER  ::   iiTrim, ijTrim                   ! T pts i/j-indice on the rim 
    173181      INTEGER  ::   iiToce, ijToce, iiUoce, ijVoce   ! T, U and V pts i/j-indice of the ocean next to the rim 
    174       REAL(wp), POINTER :: flagu, flagv              ! short cuts 
     182      REAL(wp) ::   flagu, flagv                     ! short cuts 
    175183      REAL(wp) ::   zfla                             ! Flather correction 
    176184      REAL(wp) ::   z1_2                             !  
     
    185193      ! Fill temporary array with ssh data (here we use spgu with the alias sshdta): 
    186194      igrd = 1 
     195      IF( llrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd) 
     196      ELSE                ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd) 
     197      END IF 
    187198      sshdta(:,:) = 0.0 
    188       DO jb = 1, idx%nblenrim(igrd) 
     199      DO jb = ibeg, iend 
    189200         ii = idx%nbi(jb,igrd) 
    190201         ij = idx%nbj(jb,igrd) 
     
    196207      END DO 
    197208      ! 
    198       igrd = 2      ! Flather bc on u-velocity;  
     209      igrd = 2      ! Flather bc on u-velocity 
    199210      !             ! remember that flagu=-1 if normal velocity direction is outward 
    200211      !             ! I think we should rather use after ssh ? 
    201       DO jb = 1, idx%nblenrim(igrd) 
    202          ii = idx%nbi(jb,igrd) 
    203          ij = idx%nbj(jb,igrd) 
    204          flagu => idx%flagu(jb,igrd) 
     212      IF( llrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd) 
     213      ELSE                ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd) 
     214      END IF 
     215      DO jb = ibeg, iend 
     216         ii    = idx%nbi(jb,igrd) 
     217         ij    = idx%nbj(jb,igrd) 
     218         flagu = idx%flagu(jb,igrd) 
    205219         IF( flagu == 0. ) THEN 
    206220            pua2d(ii,ij) = dta%u2d(jb) 
     
    222236      igrd = 3      ! Flather bc on v-velocity 
    223237      !             ! remember that flagv=-1 if normal velocity direction is outward 
    224       DO jb = 1, idx%nblenrim(igrd) 
    225          ii = idx%nbi(jb,igrd) 
    226          ij = idx%nbj(jb,igrd) 
    227          flagv => idx%flagv(jb,igrd) 
     238      IF( llrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd) 
     239      ELSE                ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd) 
     240      END IF 
     241      DO jb = ibeg, iend 
     242         ii    = idx%nbi(jb,igrd) 
     243         ij    = idx%nbj(jb,igrd) 
     244         flagv = idx%flagv(jb,igrd) 
    228245         IF( flagv == 0. ) THEN 
    229246            pva2d(ii,ij) = dta%v2d(jb) 
     
    246263 
    247264 
    248    SUBROUTINE bdy_dyn2d_orlanski( idx, dta, ib_bdy, pua2d, pva2d, pub2d, pvb2d, ll_npo ) 
     265   SUBROUTINE bdy_dyn2d_orlanski( idx, dta, ib_bdy, pua2d, pva2d, pub2d, pvb2d, llrim0, ll_npo ) 
    249266      !!---------------------------------------------------------------------- 
    250267      !!                 ***  SUBROUTINE bdy_dyn2d_orlanski  *** 
     
    263280      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pub2d, pvb2d  
    264281      LOGICAL,                      INTENT(in) ::   ll_npo  ! flag for NPO version 
    265  
     282      LOGICAL,                      INTENT(in) ::   llrim0   ! indicate if rim 0 is treated 
    266283      INTEGER  ::   ib, igrd                               ! dummy loop indices 
    267284      INTEGER  ::   ii, ij, iibm1, ijbm1                   ! indices 
     
    270287      igrd = 2      ! Orlanski bc on u-velocity;  
    271288      !             
    272       CALL bdy_orlanski_2d( idx, igrd, pub2d, pua2d, dta%u2d, ll_npo ) 
     289      CALL bdy_orlanski_2d( idx, igrd, pub2d, pua2d, dta%u2d, llrim0, ll_npo ) 
    273290 
    274291      igrd = 3      ! Orlanski bc on v-velocity 
    275292      !   
    276       CALL bdy_orlanski_2d( idx, igrd, pvb2d, pva2d, dta%v2d, ll_npo ) 
     293      CALL bdy_orlanski_2d( idx, igrd, pvb2d, pva2d, dta%v2d, llrim0, ll_npo ) 
    277294      ! 
    278295   END SUBROUTINE bdy_dyn2d_orlanski 
     
    288305      REAL(wp), DIMENSION(jpi,jpj,1), INTENT(inout) ::   zssh ! Sea level, need 3 dimensions to be used by bdy_nmn 
    289306      !! 
    290       INTEGER  ::   ib_bdy          ! bdy index 
     307      INTEGER ::   ib_bdy, ir      ! bdy index, rim index 
     308      INTEGER ::   ibeg, iend      ! length of rim to be treated (rim 0 or rim 1) 
     309      LOGICAL ::   llrim0          ! indicate if rim 0 is treated 
    291310      LOGICAL, DIMENSION(4) :: llsend1, llrecv1  ! indicate how communications are to be carried out 
    292311      !!---------------------------------------------------------------------- 
    293       llsend1(:) = .false. 
    294       llrecv1(:) = .false. 
    295       DO ib_bdy = 1, nb_bdy 
    296          CALL bdy_nmn( idx_bdy(ib_bdy), 1, zssh )   ! zssh is masked 
    297          llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:)   ! possibly every direction, T points 
    298          llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:)   ! possibly every direction, T points 
    299       END DO 
    300       IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction 
    301          CALL lbc_bdy_lnk( 'bdydyn2d', llsend1, llrecv1, zssh(:,:,1), 'T',  1. ) 
    302       END IF 
     312      DO ir = 1, 0, -1   ! treat rim 1 before rim 0 
     313         llsend1(:) = .false. 
     314         llrecv1(:) = .false. 
     315         IF( ir == 0 ) THEN   ;   llrim0 = .TRUE. 
     316         ELSE                 ;   llrim0 = .FALSE. 
     317         END IF 
     318         DO ib_bdy = 1, nb_bdy 
     319            CALL bdy_nmn( idx_bdy(ib_bdy), 1, zssh, llrim0 )   ! zssh is masked 
     320            llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:,ir)   ! possibly every direction, T points 
     321            llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:,ir)   ! possibly every direction, T points 
     322         END DO 
     323         IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction 
     324            CALL lbc_bdy_lnk( 'bdydyn2d', llsend1, llrecv1, zssh(:,:,1), 'T',  1. ) 
     325         END IF 
     326      END DO 
    303327      ! 
    304328   END SUBROUTINE bdy_ssh 
  • NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdydyn3d.F90

    r11071 r11191  
    4242      INTEGER, INTENT(in) ::   kt   ! Main time step counter 
    4343      ! 
    44       INTEGER ::   ib_bdy   ! loop index 
     44      INTEGER  ::   ib_bdy, ir     ! BDY set index, rim index 
     45      LOGICAL  ::   llrim0         ! indicate if rim 0 is treated 
    4546      LOGICAL, DIMENSION(4) :: llsend2, llrecv2, llsend3, llrecv3  ! indicate how communications are to be carried out 
    4647 
    4748      !!---------------------------------------------------------------------- 
    48       DO ib_bdy=1, nb_bdy 
     49      DO ir = 1, 0, -1   ! treat rim 1 before rim 0 
     50         IF( ir == 0 ) THEN   ;   llrim0 = .TRUE. 
     51         ELSE                 ;   llrim0 = .FALSE. 
     52         END IF 
     53         DO ib_bdy=1, nb_bdy 
     54            ! 
     55            SELECT CASE( cn_dyn3d(ib_bdy) ) 
     56            CASE('none')        ;   CYCLE 
     57            CASE('frs' )        ! treat the whole boundary at once 
     58               IF( ir == 0) CALL bdy_dyn3d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 
     59            CASE('specified')   ! treat the whole rim      at once 
     60               IF( ir == 0) CALL bdy_dyn3d_spe( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 
     61            CASE('zero')        ! treat the whole rim      at once 
     62               IF( ir == 0) CALL bdy_dyn3d_zro( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 
     63            CASE('orlanski' )   ;   CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, llrim0, ll_npo=.false. ) 
     64            CASE('orlanski_npo');   CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, llrim0, ll_npo=.true.  ) 
     65            CASE('zerograd')    ;   CALL bdy_dyn3d_zgrad( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy, llrim0 ) 
     66            CASE('neumann')     ;   CALL bdy_dyn3d_nmn( idx_bdy(ib_bdy), ib_bdy, llrim0 ) 
     67            CASE DEFAULT        ;   CALL ctl_stop( 'bdy_dyn3d : unrecognised option for open boundaries for baroclinic velocities' ) 
     68            END SELECT 
     69         END DO 
    4970         ! 
    50          SELECT CASE( cn_dyn3d(ib_bdy) ) 
    51          CASE('none')        ;   CYCLE 
    52          CASE('frs' )        ;   CALL bdy_dyn3d_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 
    53          CASE('specified')   ;   CALL bdy_dyn3d_spe( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 
    54          CASE('zero')        ;   CALL bdy_dyn3d_zro( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 
    55          CASE('orlanski' )   ;   CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, ll_npo=.false. ) 
    56          CASE('orlanski_npo');   CALL bdy_dyn3d_orlanski( idx_bdy(ib_bdy), dta_bdy(ib_bdy), ib_bdy, ll_npo=.true. ) 
    57          CASE('zerograd')    ;   CALL bdy_dyn3d_zgrad( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 
    58          CASE('neumann')     ;   CALL bdy_dyn3d_nmn( idx_bdy(ib_bdy), ib_bdy ) 
    59          CASE DEFAULT        ;   CALL ctl_stop( 'bdy_dyn3d : unrecognised option for open boundaries for baroclinic velocities' ) 
    60          END SELECT 
    61       END DO 
    62       ! 
    63       llsend2(:) = .false. 
    64       llrecv2(:) = .false. 
    65       llsend3(:) = .false. 
    66       llrecv3(:) = .false. 
    67       DO ib_bdy=1, nb_bdy 
    68          SELECT CASE( cn_dyn3d(ib_bdy) ) 
    69          CASE('orlanski', 'orlanski_npo') 
    70             llsend2(:) = llsend2(:) .OR. lsend_bdy(ib_bdy,2,:)   ! possibly every direction, U points 
    71             llrecv2(:) = llrecv2(:) .OR. lrecv_bdy(ib_bdy,2,:)   ! possibly every direction, U points 
    72             llsend3(:) = llsend3(:) .OR. lsend_bdy(ib_bdy,3,:)   ! possibly every direction, V points 
    73             llrecv3(:) = llrecv3(:) .OR. lrecv_bdy(ib_bdy,3,:)   ! possibly every direction, V points 
    74          CASE('zerograd') 
    75             llsend2(3:4) = llsend2(3:4) .OR. lsend_bdyint(ib_bdy,2,3:4)   ! north/south, U points 
    76             llrecv2(3:4) = llrecv2(3:4) .OR. lrecv_bdyint(ib_bdy,2,3:4)   ! north/south, U points 
    77             llsend3(1:2) = llsend3(1:2) .OR. lsend_bdyint(ib_bdy,3,1:2)   ! west/east, V points 
    78             llrecv3(1:2) = llrecv3(1:2) .OR. lrecv_bdyint(ib_bdy,3,1:2)   ! west/east, V points 
    79          CASE('neumann') 
    80             llsend2(:) = llsend2(:) .OR. lsend_bdyint(ib_bdy,2,:)   ! possibly every direction, U points 
    81             llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(ib_bdy,2,:)   ! possibly every direction, U points 
    82             llsend3(:) = llsend3(:) .OR. lsend_bdyint(ib_bdy,3,:)   ! possibly every direction, V points 
    83             llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(ib_bdy,3,:)   ! possibly every direction, V points 
    84          END SELECT 
    85       END DO 
    86       ! 
    87       IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN   ! if need to send/recv in at least one direction 
    88          CALL lbc_bdy_lnk( 'bdydyn2d', llsend2, llrecv2, ua, 'U', -1. ) 
    89       END IF 
    90       IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction 
    91          CALL lbc_bdy_lnk( 'bdydyn2d', llsend3, llrecv3, va, 'V', -1. ) 
    92       END IF 
     71         llsend2(:) = .false. 
     72         llrecv2(:) = .false. 
     73         llsend3(:) = .false. 
     74         llrecv3(:) = .false. 
     75         DO ib_bdy=1, nb_bdy 
     76            SELECT CASE( cn_dyn3d(ib_bdy) ) 
     77            CASE('orlanski', 'orlanski_npo') 
     78               llsend2(:) = llsend2(:) .OR. lsend_bdy(ib_bdy,2,:,ir)   ! possibly every direction, U points 
     79               llrecv2(:) = llrecv2(:) .OR. lrecv_bdy(ib_bdy,2,:,ir)   ! possibly every direction, U points 
     80               llsend3(:) = llsend3(:) .OR. lsend_bdy(ib_bdy,3,:,ir)   ! possibly every direction, V points 
     81               llrecv3(:) = llrecv3(:) .OR. lrecv_bdy(ib_bdy,3,:,ir)   ! possibly every direction, V points 
     82            CASE('zerograd') 
     83               llsend2(3:4) = llsend2(3:4) .OR. lsend_bdyint(ib_bdy,2,3:4,ir)   ! north/south, U points 
     84               llrecv2(3:4) = llrecv2(3:4) .OR. lrecv_bdyint(ib_bdy,2,3:4,ir)   ! north/south, U points 
     85               llsend3(1:2) = llsend3(1:2) .OR. lsend_bdyint(ib_bdy,3,1:2,ir)   ! west/east, V points 
     86               llrecv3(1:2) = llrecv3(1:2) .OR. lrecv_bdyint(ib_bdy,3,1:2,ir)   ! west/east, V points 
     87            CASE('neumann') 
     88               llsend2(:) = llsend2(:) .OR. lsend_bdyint(ib_bdy,2,:,ir)   ! possibly every direction, U points 
     89               llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(ib_bdy,2,:,ir)   ! possibly every direction, U points 
     90               llsend3(:) = llsend3(:) .OR. lsend_bdyint(ib_bdy,3,:,ir)   ! possibly every direction, V points 
     91               llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(ib_bdy,3,:,ir)   ! possibly every direction, V points 
     92            END SELECT 
     93         END DO 
     94         ! 
     95         IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN   ! if need to send/recv in at least one direction 
     96            CALL lbc_bdy_lnk( 'bdydyn2d', llsend2, llrecv2, ua, 'U', -1. ) 
     97         END IF 
     98         IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction 
     99            CALL lbc_bdy_lnk( 'bdydyn2d', llsend3, llrecv3, va, 'V', -1. ) 
     100         END IF 
     101      END DO 
    93102      ! 
    94103   END SUBROUTINE bdy_dyn3d 
     
    133142 
    134143 
    135    SUBROUTINE bdy_dyn3d_zgrad( idx, dta, kt , ib_bdy ) 
     144   SUBROUTINE bdy_dyn3d_zgrad( idx, dta, kt, ib_bdy, llrim0 ) 
    136145      !!---------------------------------------------------------------------- 
    137146      !!                  ***  SUBROUTINE bdy_dyn3d_zgrad  *** 
     
    141150      !!---------------------------------------------------------------------- 
    142151      INTEGER                     ::   kt 
    143       TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices 
    144       TYPE(OBC_DATA),  INTENT(in) ::   dta  ! OBC external data 
    145       INTEGER,         INTENT(in) ::   ib_bdy  ! BDY set index 
     152      TYPE(OBC_INDEX), INTENT(in) ::   idx      ! OBC indices 
     153      TYPE(OBC_DATA),  INTENT(in) ::   dta      ! OBC external data 
     154      INTEGER,         INTENT(in) ::   ib_bdy   ! BDY set index 
     155      LOGICAL,         INTENT(in) ::   llrim0   ! indicate if rim 0 is treated 
    146156      !! 
    147157      INTEGER  ::   jb, jk         ! dummy loop indices 
    148158      INTEGER  ::   ii, ij, igrd   ! local integers 
    149159      INTEGER  ::   flagu, flagv           ! short cuts 
     160      INTEGER  ::   ibeg, iend     ! length of rim to be treated (rim 0 or rim 1 or both) 
    150161      !!---------------------------------------------------------------------- 
    151162      ! 
    152163      igrd = 2                      ! Copying tangential velocity into bdy points 
    153       DO jb = 1, idx%nblenrim(igrd) 
     164      IF( llrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd) 
     165      ELSE                ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd) 
     166      ENDIF 
     167      DO jb = ibeg, iend 
    154168         ii    = idx%nbi(jb,igrd) 
    155169         ij    = idx%nbj(jb,igrd) 
     
    168182      ! 
    169183      igrd = 3                      ! Copying tangential velocity into bdy points 
    170       DO jb = 1, idx%nblenrim(igrd) 
     184      IF( llrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd) 
     185      ELSE               ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd) 
     186      ENDIF 
     187      DO jb = ibeg, iend 
    171188         ii    = idx%nbi(jb,igrd) 
    172189         ij    = idx%nbj(jb,igrd) 
     
    268285 
    269286 
    270    SUBROUTINE bdy_dyn3d_orlanski( idx, dta, ib_bdy, ll_npo ) 
     287   SUBROUTINE bdy_dyn3d_orlanski( idx, dta, ib_bdy, llrim0, ll_npo ) 
    271288      !!---------------------------------------------------------------------- 
    272289      !!                 ***  SUBROUTINE bdy_dyn3d_orlanski  *** 
     
    280297      TYPE(OBC_INDEX),              INTENT(in) ::   idx  ! OBC indices 
    281298      TYPE(OBC_DATA),               INTENT(in) ::   dta  ! OBC external data 
    282       INTEGER,                      INTENT(in) ::   ib_bdy  ! BDY set index 
    283       LOGICAL,                      INTENT(in) ::   ll_npo  ! switch for NPO version 
     299      INTEGER,                      INTENT(in) ::   ib_bdy   ! BDY set index 
     300      LOGICAL,                      INTENT(in) ::   llrim0   ! indicate if rim 0 is treated 
     301      LOGICAL,                      INTENT(in) ::   ll_npo   ! switch for NPO version 
    284302 
    285303      INTEGER  ::   jb, igrd                               ! dummy loop indices 
     
    290308      igrd = 2      ! Orlanski bc on u-velocity;  
    291309      !             
    292       CALL bdy_orlanski_3d( idx, igrd, ub, ua, dta%u3d, ll_npo ) 
     310      CALL bdy_orlanski_3d( idx, igrd, ub, ua, dta%u3d, ll_npo, llrim0 ) 
    293311 
    294312      igrd = 3      ! Orlanski bc on v-velocity 
    295313      !   
    296       CALL bdy_orlanski_3d( idx, igrd, vb, va, dta%v3d, ll_npo ) 
     314      CALL bdy_orlanski_3d( idx, igrd, vb, va, dta%v3d, ll_npo, llrim0 ) 
    297315      ! 
    298316   END SUBROUTINE bdy_dyn3d_orlanski 
     
    347365 
    348366 
    349    SUBROUTINE bdy_dyn3d_nmn( idx, ib_bdy ) 
     367   SUBROUTINE bdy_dyn3d_nmn( idx, ib_bdy, llrim0 ) 
    350368      !!---------------------------------------------------------------------- 
    351369      !!                 ***  SUBROUTINE bdy_dyn3d_nmn  *** 
     
    356374      !! 
    357375      !!---------------------------------------------------------------------- 
    358       TYPE(OBC_INDEX),              INTENT(in) ::   idx  ! OBC indices 
    359       INTEGER,                      INTENT(in) ::   ib_bdy  ! BDY set index 
    360       INTEGER  ::   igrd                                    ! dummy indice 
     376      TYPE(OBC_INDEX), INTENT(in) ::   idx      ! OBC indices 
     377      INTEGER,         INTENT(in) ::   ib_bdy   ! BDY set index 
     378      LOGICAL,         INTENT(in) ::   llrim0   ! indicate if rim 0 is treated 
     379      INTEGER  ::   igrd                        ! dummy indice 
    361380      !!---------------------------------------------------------------------- 
    362381      ! 
     
    365384      igrd = 2      ! Neumann bc on u-velocity;  
    366385      !             
    367       CALL bdy_nmn( idx, igrd, ua )   ! ua is masked 
     386      CALL bdy_nmn( idx, igrd, ua, llrim0 )   ! ua is masked 
    368387 
    369388      igrd = 3      ! Neumann bc on v-velocity 
    370389      !   
    371       CALL bdy_nmn( idx, igrd, va )   ! va is masked 
     390      CALL bdy_nmn( idx, igrd, va, llrim0 )   ! va is masked 
    372391      ! 
    373392   END SUBROUTINE bdy_dyn3d_nmn 
  • NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdyice.F90

    r11071 r11191  
    5555      INTEGER, INTENT(in) ::   kt   ! Main time step counter 
    5656      ! 
    57       INTEGER ::   jbdy                               ! BDY set index 
    58       LOGICAL, DIMENSION(4)       :: llsend1, llrecv1   ! indicate how communications are to be carried out 
     57      INTEGER ::   jbdy, ir                             ! BDY set index, rim index 
     58      INTEGER ::   ibeg, iend                           ! length of rim to be treated (rim 0 or rim 1) 
     59      LOGICAL ::   llrim0                               ! indicate if rim 0 is treated 
     60      LOGICAL, DIMENSION(4)  :: llsend1, llrecv1        ! indicate how communications are to be carried out 
    5961      !!---------------------------------------------------------------------- 
    6062      ! controls 
     
    6466      CALL ice_var_glo2eqv 
    6567      ! 
    66       DO jbdy = 1, nb_bdy 
     68      DO ir = 1, 0, -1   ! treat rim 1 before rim 0 
     69         IF( ir == 0 ) THEN   ;   llrim0 = .TRUE. 
     70         ELSE                 ;   llrim0 = .FALSE. 
     71         END IF 
     72         DO jbdy = 1, nb_bdy 
     73            ! 
     74            SELECT CASE( cn_ice(jbdy) ) 
     75            CASE('none')   ;   CYCLE 
     76            CASE('frs' )   ;   CALL bdy_ice_frs( idx_bdy(jbdy), dta_bdy(jbdy), kt, jbdy, llrim0 ) 
     77            CASE DEFAULT 
     78               CALL ctl_stop( 'bdy_ice : unrecognised option for open boundaries for ice fields' ) 
     79            END SELECT 
     80            ! 
     81         END DO 
    6782         ! 
    68          SELECT CASE( cn_ice(jbdy) ) 
    69          CASE('none')   ;   CYCLE 
    70          CASE('frs' )   ;   CALL bdy_ice_frs( idx_bdy(jbdy), dta_bdy(jbdy), kt, jbdy ) 
    71          CASE DEFAULT 
    72             CALL ctl_stop( 'bdy_ice : unrecognised option for open boundaries for ice fields' ) 
    73          END SELECT 
    74          ! 
    75       END DO 
    76       ! 
    77       ! Update bdy points 
    78       llsend1(:) = .false. 
    79       llrecv1(:) = .false. 
    80       DO jbdy = 1, nb_bdy 
    81          IF( cn_ice(jbdy) == 'frs' ) THEN 
    82             llsend1(:) = llsend1(:) .OR. lsend_bdyint(jbdy,1,:)   ! possibly every direction, T points 
    83             llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(jbdy,1,:)   ! possibly every direction, T points 
     83         ! Update bdy points 
     84         llsend1(:) = .false. 
     85         llrecv1(:) = .false. 
     86         DO jbdy = 1, nb_bdy 
     87            IF( cn_ice(jbdy) == 'frs' ) THEN 
     88               llsend1(:) = llsend1(:) .OR. lsend_bdyint(jbdy,1,:,ir)   ! possibly every direction, T points 
     89               llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(jbdy,1,:,ir)   ! possibly every direction, T points 
     90            END IF 
     91         END DO   ! jbdy 
     92         IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction 
     93            ! exchange 3d arrays 
     94            CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, a_i , 'T', 1., h_i , 'T', 1., h_s , 'T', 1. & 
     95                 &                                            , oa_i, 'T', 1., a_ip, 'T', 1., v_ip, 'T', 1. & 
     96                 &                                            , s_i , 'T', 1., t_su, 'T', 1., v_i , 'T', 1. & 
     97                 &                                            , v_s , 'T', 1., sv_i, 'T', 1.                ) 
     98            ! exchange 4d arrays 
     99            CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, t_s , 'T', 1., e_s , 'T', 1. )   ! third dimension = 1 
     100            CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, t_i , 'T', 1., e_i , 'T', 1. )   ! third dimension = jpk 
    84101         END IF 
    85       END DO 
    86       IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction 
    87          ! exchange 3d arrays 
    88          CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, a_i , 'T', 1., h_i , 'T', 1., h_s , 'T', 1. & 
    89               &                                            , oa_i, 'T', 1., a_ip, 'T', 1., v_ip, 'T', 1. & 
    90               &                                            , s_i , 'T', 1., t_su, 'T', 1., v_i , 'T', 1. & 
    91               &                                            , v_s , 'T', 1., sv_i, 'T', 1.                ) 
    92          ! exchange 4d arrays 
    93          CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, t_s , 'T', 1., e_s , 'T', 1. )   ! third dimension = 1 
    94          CALL lbc_bdy_lnk_multi( 'bdyice', llsend1, llrecv1, t_i , 'T', 1., e_i , 'T', 1. )   ! third dimension = jpk 
    95       END IF 
     102      END DO   ! ir 
    96103      ! 
    97104      CALL ice_cor( kt , 0 )      ! -- In case categories are out of bounds, do a remapping 
     
    108115 
    109116 
    110    SUBROUTINE bdy_ice_frs( idx, dta, kt, jbdy ) 
     117   SUBROUTINE bdy_ice_frs( idx, dta, kt, jbdy, llrim0 ) 
    111118      !!------------------------------------------------------------------------------ 
    112119      !!                 ***  SUBROUTINE bdy_ice_frs  *** 
     
    117124      !!             dimensional baroclinic ocean model with realistic topography. Tellus, 365-382. 
    118125      !!------------------------------------------------------------------------------ 
    119       TYPE(OBC_INDEX), INTENT(in) ::   idx     ! OBC indices 
    120       TYPE(OBC_DATA),  INTENT(in) ::   dta     ! OBC external data 
    121       INTEGER,         INTENT(in) ::   kt      ! main time-step counter 
    122       INTEGER,         INTENT(in) ::   jbdy    ! BDY set index 
     126      TYPE(OBC_INDEX), INTENT(in) ::   idx      ! OBC indices 
     127      TYPE(OBC_DATA),  INTENT(in) ::   dta      ! OBC external data 
     128      INTEGER,         INTENT(in) ::   kt       ! main time-step counter 
     129      INTEGER,         INTENT(in) ::   jbdy     ! BDY set index 
     130      LOGICAL,         INTENT(in) ::   llrim0   ! indicate if rim 0 is treated 
    123131      ! 
    124132      INTEGER  ::   jpbound            ! 0 = incoming ice 
    125133      !                                ! 1 = outgoing ice 
     134      INTEGER  ::   ibeg, iend         ! length of rim to be treated (rim 0 or rim 1) 
    126135      INTEGER  ::   i_bdy, jgrd        ! dummy loop indices 
    127136      INTEGER  ::   ji, jj, jk, jl, ib, jb 
     
    132141      ! 
    133142      jgrd = 1      ! Everything is at T-points here 
     143      IF( llrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(jgrd) 
     144      ELSE                ;   ibeg = idx%nblenrim0(jgrd)+1   ;   iend = idx%nblenrim(jgrd) 
     145      END IF 
    134146      ! 
    135147      DO jl = 1, jpl 
    136          DO i_bdy = 1, idx%nblenrim(jgrd) 
     148         DO i_bdy = ibeg, iend 
    137149            ji    = idx%nbi(i_bdy,jgrd) 
    138150            jj    = idx%nbj(i_bdy,jgrd) 
     
    162174 
    163175      DO jl = 1, jpl 
    164          DO i_bdy = 1, idx%nblenrim(jgrd) 
     176         DO i_bdy = ibeg, iend 
    165177            ji = idx%nbi(i_bdy,jgrd) 
    166178            jj = idx%nbj(i_bdy,jgrd) 
     
    295307      CHARACTER(len=1), INTENT(in)  ::   cd_type   ! nature of velocity grid-points 
    296308      ! 
    297       INTEGER  ::   i_bdy, jgrd      ! dummy loop indices 
    298       INTEGER  ::   ji, jj           ! local scalar 
    299       INTEGER  ::   jbdy             ! BDY set index 
     309      INTEGER  ::   i_bdy, jgrd       ! dummy loop indices 
     310      INTEGER  ::   ji, jj            ! local scalar 
     311      INTEGER  ::   jbdy, ir     ! BDY set index, rim index 
     312      INTEGER  ::   ibeg, iend   ! length of rim to be treated (rim 0 or rim 1) 
    300313      REAL(wp) ::   zmsk1, zmsk2, zflag 
    301314      LOGICAL, DIMENSION(4) :: llsend2, llrecv2, llsend3, llrecv3  ! indicate how communications are to be carried out 
     
    303316      IF( ln_timing )   CALL timing_start('bdy_ice_dyn') 
    304317      ! 
    305       DO jbdy=1, nb_bdy 
     318      DO ir = 1, 0, -1 
     319         DO jbdy = 1, nb_bdy 
     320            ! 
     321            SELECT CASE( cn_ice(jbdy) ) 
     322               ! 
     323            CASE('none') 
     324               CYCLE 
     325               ! 
     326            CASE('frs') 
     327               ! 
     328               IF( nn_ice_dta(jbdy) == 0 ) CYCLE            ! case ice boundaries = initial conditions  
     329               !                                            !      do not change ice velocity (it is only computed by rheology) 
     330               SELECT CASE ( cd_type ) 
     331                  !      
     332               CASE ( 'U' )   
     333                  jgrd = 2      ! u velocity 
     334                  IF( ir == 0 ) THEN   ;   ibeg = 1                                 ;   iend = idx_bdy(jbdy)%nblenrim0(jgrd) 
     335                  ELSE                 ;   ibeg = idx_bdy(jbdy)%nblenrim0(jgrd)+1   ;   iend = idx_bdy(jbdy)%nblenrim(jgrd) 
     336                  END IF 
     337                  DO i_bdy = ibeg, iend 
     338                     ji    = idx_bdy(jbdy)%nbi(i_bdy,jgrd) 
     339                     jj    = idx_bdy(jbdy)%nbj(i_bdy,jgrd) 
     340                     zflag = idx_bdy(jbdy)%flagu(i_bdy,jgrd) 
     341                     !     i-1  i   i    |  !        i  i i+1 |  !          i  i i+1 | 
     342                     !      >  ice  >    |  !        o  > ice |  !          o  >  o  |       
     343                     ! => set at u_ice(i-1) !  => set to O       !  => unchanged 
     344                     IF( zflag == -1. .AND. ji > 1 .AND. ji < jpi )   THEN   
     345                        IF    ( vt_i(ji  ,jj) > 0. )   THEN   ;   u_ice(ji,jj) = u_ice(ji-1,jj)  
     346                        ELSEIF( vt_i(ji+1,jj) > 0. )   THEN   ;   u_ice(ji,jj) = 0._wp 
     347                        END IF 
     348                     END IF 
     349                     ! |    i  i+1 i+1        !  |  i   i i+1        !  | i  i i+1 
     350                     ! |    >  ice  >         !  | ice  >  o         !  | o  >  o    
     351                     ! => set at u_ice(i+1)   !     => set to O      !     =>  unchanged 
     352                     IF( zflag ==  1. .AND. ji+1 < jpi+1 )   THEN 
     353                        IF    ( vt_i(ji+1,jj) > 0. )   THEN   ;   u_ice(ji,jj) = u_ice(ji+1,jj) 
     354                        ELSEIF( vt_i(ji  ,jj) > 0. )   THEN   ;   u_ice(ji,jj) = 0._wp 
     355                        END IF 
     356                     END IF 
     357                     ! 
     358                     IF( zflag ==  0. )   u_ice(ji,jj) = 0._wp   ! u_ice = 0  if north/south bdy   
     359                     ! 
     360                  END DO 
     361                  ! 
     362               CASE ( 'V' ) 
     363                  jgrd = 3      ! v velocity 
     364                  IF( ir == 0 ) THEN   ;   ibeg = 1                                 ;   iend = idx_bdy(jbdy)%nblenrim0(jgrd) 
     365                  ELSE                 ;   ibeg = idx_bdy(jbdy)%nblenrim0(jgrd)+1   ;   iend = idx_bdy(jbdy)%nblenrim(jgrd) 
     366                  END IF 
     367                  DO i_bdy = ibeg, iend 
     368                     ji    = idx_bdy(jbdy)%nbi(i_bdy,jgrd) 
     369                     jj    = idx_bdy(jbdy)%nbj(i_bdy,jgrd) 
     370                     zflag = idx_bdy(jbdy)%flagv(i_bdy,jgrd) 
     371                     !    ¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨     !  ¨¨¨¨ïce¨¨¨(jj+1)¨¨     ! ¨¨¨¨¨¨ö¨¨¨¨(jj+1)        
     372                     !       ^    (jj  )       !       ^    (jj  )       !       ^    (jj  )        
     373                     !      ice   (jj  )       !       o    (jj  )       !       o    (jj  )        
     374                     !       ^    (jj-1)       !                         ! 
     375                     ! => set to u_ice(jj-1)   !  =>   set to 0          !   => unchanged         
     376                     IF( zflag == -1. .AND. jj > 1 .AND. jj < jpj )   THEN                  
     377                        IF    ( vt_i(ji,jj  ) > 0. )   THEN   ;   v_ice(ji,jj) = v_ice(ji,jj-1) 
     378                        ELSEIF( vt_i(ji,jj+1) > 0. )   THEN   ;   v_ice(ji,jj) = 0._wp 
     379                        END IF 
     380                     END IF 
     381                     !       ^    (jj+1)       !                         !               
     382                     !      ice   (jj+1)       !       o    (jj+1)       !       o    (jj+1)        
     383                     !       ^    (jj  )       !       ^    (jj  )       !       ^    (jj  ) 
     384                     !   ________________      !  ____ice___(jj  )_      !  _____o____(jj  )  
     385                     ! => set to u_ice(jj+1)   !    => set to 0          !    => unchanged   
     386                     IF( zflag ==  1. .AND. jj < jpj )   THEN               
     387                        IF    ( vt_i(ji,jj+1) > 0. )   THEN   ;   v_ice(ji,jj) = v_ice(ji,jj+1) 
     388                        ELSEIF( vt_i(ji,jj  ) > 0. )   THEN   ;   v_ice(ji,jj) = 0._wp 
     389                        END IF 
     390                     END IF 
     391                     ! 
     392                     IF( zflag ==  0. )   v_ice(ji,jj) = 0._wp   ! v_ice = 0  if west/east bdy   
     393                     ! 
     394                  END DO 
     395                  ! 
     396               END SELECT 
     397               ! 
     398            CASE DEFAULT 
     399               CALL ctl_stop( 'bdy_ice_dyn : unrecognised option for open boundaries for ice fields' ) 
     400            END SELECT 
     401            ! 
     402         END DO    ! jbdy 
    306403         ! 
    307          SELECT CASE( cn_ice(jbdy) ) 
    308          ! 
    309          CASE('none') 
    310             CYCLE 
    311             ! 
    312          CASE('frs') 
    313             ! 
    314             IF( nn_ice_dta(jbdy) == 0 ) CYCLE            ! case ice boundaries = initial conditions  
    315             !                                            !      do not change ice velocity (it is only computed by rheology) 
    316             SELECT CASE ( cd_type ) 
    317             !      
    318             CASE ( 'U' )   
    319                jgrd = 2      ! u velocity 
    320                DO i_bdy = 1, idx_bdy(jbdy)%nblenrim(jgrd) 
    321                   ji    = idx_bdy(jbdy)%nbi(i_bdy,jgrd) 
    322                   jj    = idx_bdy(jbdy)%nbj(i_bdy,jgrd) 
    323                   zflag = idx_bdy(jbdy)%flagu(i_bdy,jgrd) 
    324                   !     i-1  i   i    |  !        i  i i+1 |  !          i  i i+1 | 
    325                   !      >  ice  >    |  !        o  > ice |  !          o  >  o  |       
    326                   ! => set at u_ice(i-1) !  => set to O       !  => unchanged 
    327                   IF( zflag == -1. .AND. ji > 1 .AND. ji < jpi )   THEN   
    328                      IF    ( vt_i(ji  ,jj) > 0. )   THEN   ;   u_ice(ji,jj) = u_ice(ji-1,jj)  
    329                      ELSEIF( vt_i(ji+1,jj) > 0. )   THEN   ;   u_ice(ji,jj) = 0._wp 
    330                      END IF 
    331                   END IF 
    332                   ! |    i  i+1 i+1        !  |  i   i i+1        !  | i  i i+1 
    333                   ! |    >  ice  >         !  | ice  >  o         !  | o  >  o    
    334                   ! => set at u_ice(i+1)   !     => set to O      !     =>  unchanged 
    335                   IF( zflag ==  1. .AND. ji+1 < jpi+1 )   THEN 
    336                      IF    ( vt_i(ji+1,jj) > 0. )   THEN   ;   u_ice(ji,jj) = u_ice(ji+1,jj) 
    337                      ELSEIF( vt_i(ji  ,jj) > 0. )   THEN   ;   u_ice(ji,jj) = 0._wp 
    338                      END IF 
    339                   END IF 
    340                   ! 
    341                   IF( zflag ==  0. )   u_ice(ji,jj) = 0._wp   ! u_ice = 0  if north/south bdy   
    342                   ! 
    343                END DO 
    344                ! 
    345             CASE ( 'V' ) 
    346                jgrd = 3      ! v velocity 
    347                DO i_bdy = 1, idx_bdy(jbdy)%nblenrim(jgrd) 
    348                   ji    = idx_bdy(jbdy)%nbi(i_bdy,jgrd) 
    349                   jj    = idx_bdy(jbdy)%nbj(i_bdy,jgrd) 
    350                   zflag = idx_bdy(jbdy)%flagv(i_bdy,jgrd) 
    351                   !    ¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨     !  ¨¨¨¨ïce¨¨¨(jj+1)¨¨     ! ¨¨¨¨¨¨ö¨¨¨¨(jj+1)        
    352                   !       ^    (jj  )       !       ^    (jj  )       !       ^    (jj  )        
    353                   !      ice   (jj  )       !       o    (jj  )       !       o    (jj  )        
    354                   !       ^    (jj-1)       !                         ! 
    355                   ! => set to u_ice(jj-1)   !  =>   set to 0          !   => unchanged         
    356                   IF( zflag == -1. .AND. jj > 1 .AND. jj < jpj )   THEN                  
    357                     IF    ( vt_i(ji,jj  ) > 0. )   THEN   ;   v_ice(ji,jj) = v_ice(ji,jj-1) 
    358                     ELSEIF( vt_i(ji,jj+1) > 0. )   THEN   ;   v_ice(ji,jj) = 0._wp 
    359                     END IF 
    360                   END IF  
    361                   !       ^    (jj+1)       !                         !               
    362                   !      ice   (jj+1)       !       o    (jj+1)       !       o    (jj+1)        
    363                   !       ^    (jj  )       !       ^    (jj  )       !       ^    (jj  ) 
    364                   !   ________________      !  ____ice___(jj  )_      !  _____o____(jj  )  
    365                   ! => set to u_ice(jj+1)   !    => set to 0          !    => unchanged   
    366                   IF( zflag ==  1. .AND. jj < jpj )   THEN               
    367                      IF    ( vt_i(ji,jj+1) > 0. )   THEN   ;   v_ice(ji,jj) = v_ice(ji,jj+1) 
    368                      ELSEIF( vt_i(ji,jj  ) > 0. )   THEN   ;   v_ice(ji,jj) = 0._wp 
    369                      END IF 
    370                   END IF                                           
    371                   ! 
    372                   IF( zflag ==  0. )   v_ice(ji,jj) = 0._wp   ! v_ice = 0  if west/east bdy   
    373                   ! 
    374                END DO 
    375                ! 
    376             END SELECT 
    377             ! 
    378          CASE DEFAULT 
    379             CALL ctl_stop( 'bdy_ice_dyn : unrecognised option for open boundaries for ice fields' ) 
     404         SELECT CASE ( cd_type )         
     405         CASE ( 'U' )  
     406            llsend2(:) = .false.   ;   llrecv2(:) = .false. 
     407            DO jbdy = 1, nb_bdy 
     408               IF( cn_ice(jbdy) == 'frs' .AND. nn_ice_dta(jbdy) /= 0 ) THEN 
     409                  llsend2(:) = llsend2(:) .OR. lsend_bdyint(jbdy,2,:,ir)   ! possibly every direction, U points 
     410                  llsend2(1) = llsend2(1) .OR. lsend_bdyext(jbdy,2,1,ir)   ! neighbour might search point towards its west bdy 
     411                  llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(jbdy,2,:,ir)   ! possibly every direction, U points 
     412                  llrecv2(2) = llrecv2(2) .OR. lrecv_bdyext(jbdy,2,2,ir)   ! might search point towards east bdy 
     413               END IF 
     414            END DO 
     415            IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN   ! if need to send/recv in at least one direction 
     416               CALL lbc_bdy_lnk( 'bdyice', llsend2, llrecv2, u_ice, 'U', -1. ) 
     417            END IF 
     418         CASE ( 'V' ) 
     419            llsend3(:) = .false.   ;   llrecv3(:) = .false. 
     420            DO jbdy = 1, nb_bdy 
     421               IF( cn_ice(jbdy) == 'frs' .AND. nn_ice_dta(jbdy) /= 0 ) THEN 
     422                  llsend3(:) = llsend3(:) .OR. lsend_bdyint(jbdy,3,:,ir)   ! possibly every direction, V points 
     423                  llsend3(3) = llsend3(3) .OR. lsend_bdyext(jbdy,3,3,ir)   ! neighbour might search point towards its south bdy 
     424                  llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(jbdy,3,:,ir)   ! possibly every direction, V points 
     425                  llrecv3(4) = llrecv3(4) .OR. lrecv_bdyext(jbdy,3,4,ir)   ! might search point towards north bdy 
     426               END IF 
     427            END DO 
     428            IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction 
     429               CALL lbc_bdy_lnk( 'bdyice', llsend3, llrecv3, v_ice, 'V', -1. ) 
     430            END IF 
    380431         END SELECT 
    381          ! 
    382       END DO   ! jbdy 
    383       ! 
    384       SELECT CASE ( cd_type )         
    385       CASE ( 'U' )  
    386          llsend2(:) = .false.   ;   llrecv2(:) = .false. 
    387          DO jbdy = 1, nb_bdy 
    388             IF( cn_ice(jbdy) == 'frs' .AND. nn_ice_dta(jbdy) /= 0 ) THEN 
    389                llsend2(:) = llsend2(:) .OR. lsend_bdyint(jbdy,2,:)   ! possibly every direction, U points 
    390                llsend2(1) = llsend2(1) .OR. lsend_bdyext(jbdy,2,1)   ! neighbour might search point towards bdy on its east 
    391                llrecv2(:) = llrecv2(:) .OR. lrecv_bdyint(jbdy,2,:)   ! possibly every direction, U points 
    392                llrecv2(2) = llrecv2(2) .OR. lrecv_bdyext(jbdy,2,2)   ! might search point towards bdy on the east 
    393             END IF 
    394          END DO 
    395          IF( ANY(llsend2) .OR. ANY(llrecv2) ) THEN   ! if need to send/recv in at least one direction 
    396             CALL lbc_bdy_lnk( 'bdyice', llsend2, llrecv2, u_ice, 'U', -1. ) 
    397          END IF 
    398       CASE ( 'V' ) 
    399          llsend3(:) = .false.   ;   llrecv3(:) = .false. 
    400          DO jbdy = 1, nb_bdy 
    401             IF( cn_ice(jbdy) == 'frs' .AND. nn_ice_dta(jbdy) /= 0 ) THEN 
    402                llsend3(:) = llsend3(:) .OR. lsend_bdyint(jbdy,3,:)   ! possibly every direction, V points 
    403                llsend3(3) = llsend3(3) .OR. lsend_bdyext(jbdy,3,3)   ! neighbour might search point towards bdy on its north 
    404                llrecv3(:) = llrecv3(:) .OR. lrecv_bdyint(jbdy,3,:)   ! possibly every direction, V points 
    405                llrecv3(4) = llrecv3(4) .OR. lrecv_bdyext(jbdy,3,4)   ! might search point towards bdy on the north 
    406             END IF 
    407          END DO 
    408          IF( ANY(llsend3) .OR. ANY(llrecv3) ) THEN   ! if need to send/recv in at least one direction 
    409             CALL lbc_bdy_lnk( 'bdyice', llsend3, llrecv3, v_ice, 'V', -1. ) 
    410          END IF 
    411       END SELECT 
     432      END DO   ! ir 
    412433      ! 
    413434      IF( ln_timing )   CALL timing_stop('bdy_ice_dyn') 
  • NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdyini.F90

    r11183 r11191  
    3333   PRIVATE 
    3434 
    35    PUBLIC   bdy_init   ! routine called in nemo_init 
     35   PUBLIC   bdy_init    ! routine called in nemo_init 
    3636   PUBLIC   find_neib   ! routine called in bdy_nmn 
    3737 
     
    126126      !! ** Input   :  bdy_init.nc, input file for unstructured open boundaries 
    127127      !!----------------------------------------------------------------------       
    128       INTEGER  ::   ib_bdy, ii, ij, ik, igrd, ib, ir, iseg ! dummy loop indices 
    129       INTEGER  ::   icount, icountr, ibr_max, ilen1, ibm1  ! local integers 
     128      INTEGER  ::   ib_bdy, ii, ij, igrd, ib, ir, iseg     ! dummy loop indices 
     129      INTEGER  ::   icount, icountr, icountr0, ibr_max     ! local integers 
     130      INTEGER  ::   ilen1, ibm1                            !   -       - 
    130131      INTEGER  ::   iwe, ies, iso, ino, inum, id_dummy     !   -       - 
    131       INTEGER  ::   igrd_start, igrd_end, jpbdta           !   -       - 
    132       INTEGER  ::   jpbdtau, jpbdtas                       !   -       - 
     132      INTEGER  ::   jpbdta, jpbdtau, jpbdtas               !   -       - 
    133133      INTEGER  ::   ib_bdy1, ib_bdy2, ib1, ib2             !   -       - 
    134       INTEGER  ::   i_offset, j_offset, inn                !   -       - 
    135134      INTEGER  ::   ii1, ii2, ii3, ij1, ij2, ij3           !   -       - 
    136135      INTEGER  ::   iibe, ijbe, iibi, ijbi                 !   -       - 
    137136      INTEGER  ::   flagu, flagv                           ! short cuts 
    138       INTEGER , POINTER  ::  nbi, nbj, nbr                 ! short cuts 
    139137      REAL(wp), POINTER, DIMENSION(:,:)       ::   pmask   ! pointer to 2D mask fields 
    140138      REAL(wp) ::   zefl, zwfl, znfl, zsfl                 ! local scalars 
     
    145143      CHARACTER(LEN=1),DIMENSION(jpbgrd)      ::   cgrid 
    146144      INTEGER :: jpk_max 
    147       LOGICAL :: llsend_ea, llsend_we, llsend_so, llsend_no  ! Flags for boundaries sending 
    148       LOGICAL :: llrecv_ea, llrecv_we, llrecv_so, llrecv_no  ! Flags for boundaries receiving 
    149       REAL(wp), TARGET, DIMENSION(jpi,jpj) ::   zfmask  ! temporary fmask array excluding coastal boundary condition (shlat) 
    150       REAL(wp)        , DIMENSION(jpi,jpj) ::   ztmp 
    151       LOGICAL  ::   llnon, llson, llean, llwen          ! local logicals 
     145      REAL(wp), DIMENSION(jpi,jpj) ::   zfmask                  ! temporary fmask array excluding coastal boundary condition (shlat) 
     146      REAL(wp), DIMENSION(jpi,jpj) ::   ztmask, zumask, zvmask  ! temporary u/v mask array 
    152147      !! 
    153148      CHARACTER(LEN=1)                     ::   ctypebdy   !     -        -  
     
    798793      ! 
    799794      iwe = mig(1) 
    800       ies = mig(nlci) 
     795      ies = mig(jpi) 
    801796      iso = mjg(1)  
    802       ino = mjg(nlcj)  
     797      ino = mjg(jpj)  
    803798      ! 
    804799      DO ib_bdy = 1, nb_bdy 
    805800         DO igrd = 1, jpbgrd 
    806             icount  = 0   ! initialization of local bdy length 
    807             icountr = 0   ! initialization of local rim bdy length 
    808             idx_bdy(ib_bdy)%nblen(igrd)    = 0 
    809             idx_bdy(ib_bdy)%nblenrim(igrd) = 0 
     801            icount   = 0   ! initialization of local bdy length 
     802            icountr  = 0   ! initialization of local rim 0 and rim 1 bdy length 
     803            icountr0 = 0   ! initialization of local rim 0 bdy length 
     804            idx_bdy(ib_bdy)%nblen(igrd)     = 0 
     805            idx_bdy(ib_bdy)%nblenrim(igrd)  = 0 
     806            idx_bdy(ib_bdy)%nblenrim0(igrd) = 0 
    810807            DO ib = 1, nblendta(igrd,ib_bdy) 
    811808               ! check that data is in correct order in file 
     
    823820                  ! 
    824821                  icount = icount + 1 
    825                   IF( nbrdta(ib,igrd,ib_bdy) == 1 )   icountr = icountr + 1 
     822                  IF( nbrdta(ib,igrd,ib_bdy) == 1 .OR. nbrdta(ib,igrd,ib_bdy) == 0 )   icountr = icountr + 1 
     823                  IF( nbrdta(ib,igrd,ib_bdy) == 0 )   icountr0 = icountr0 + 1 
    826824               ENDIF 
    827825            END DO 
    828             idx_bdy(ib_bdy)%nblenrim(igrd) = icountr !: length of rim boundary data on each proc 
    829             idx_bdy(ib_bdy)%nblen   (igrd) = icount  !: length of boundary data on each proc         
     826            idx_bdy(ib_bdy)%nblen    (igrd) = icount   !: length of boundary data on each proc 
     827            idx_bdy(ib_bdy)%nblenrim (igrd) = icountr  !: length of rim 0 and rim 1 boundary data on each proc    
     828            idx_bdy(ib_bdy)%nblenrim0(igrd) = icountr0 !: length of rim 0 boundary data on each proc      
    830829         END DO   ! igrd 
    831830 
     
    849848            icount  = 0 
    850849            ! Outer loop on rimwidth to ensure outermost points come first in the local arrays. 
    851             DO ir=1, nn_rimwidth(ib_bdy) 
     850            DO ir = 0, nn_rimwidth(ib_bdy) 
    852851               DO ib = 1, nblendta(igrd,ib_bdy) 
    853852                  ! check if point is in local domain and equals ir 
     
    870869      ! Initialize array indicating communications in bdy 
    871870      ! ------------------------------------------------- 
    872       ALLOCATE( lsend_bdy(nb_bdy,jpbgrd,4), lrecv_bdy(nb_bdy,jpbgrd,4) ) 
    873       lsend_bdy(:,:,:) = .false. 
    874       lrecv_bdy(:,:,:) = .false.  
     871      ALLOCATE( lsend_bdy(nb_bdy,jpbgrd,4,0:1), lrecv_bdy(nb_bdy,jpbgrd,4,0:1) ) 
     872      lsend_bdy(:,:,:,:) = .false. 
     873      lrecv_bdy(:,:,:,:) = .false.  
    875874 
    876875      DO ib_bdy = 1, nb_bdy 
     
    879878               ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
    880879               ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 
     880               IF( ib .LE. idx_bdy(ib_bdy)%nblenrim0(igrd) ) THEN   ;   ir = 0 
     881               ELSE                                                 ;   ir = 1 
     882               END IF 
    881883               ! 
    882884               ! check if point has to be sent     to   a neighbour 
    883885               ! W neighbour and on the inner left  side 
    884                IF( ii == 2     .and. (nbondi == 0 .or. nbondi ==  1) )   lsend_bdy(ib_bdy,igrd,1) = .true. 
     886               IF( ii == 2     .and. (nbondi == 0 .or. nbondi ==  1) )   lsend_bdy(ib_bdy,igrd,1,ir) = .true. 
    885887               ! E neighbour and on the inner right side 
    886                IF( ii == jpi-1 .and. (nbondi == 0 .or. nbondi == -1) )   lsend_bdy(ib_bdy,igrd,2) = .true. 
     888               IF( ii == jpi-1 .and. (nbondi == 0 .or. nbondi == -1) )   lsend_bdy(ib_bdy,igrd,2,ir) = .true. 
    887889               ! S neighbour and on the inner down side 
    888                IF( ij == 2     .and. (nbondj == 0 .or. nbondj ==  1) )   lsend_bdy(ib_bdy,igrd,3) = .true. 
     890               IF( ij == 2     .and. (nbondj == 0 .or. nbondj ==  1) )   lsend_bdy(ib_bdy,igrd,3,ir) = .true. 
    889891               ! N neighbour and on the inner up   side 
    890                IF( ij == jpj-1 .and. (nbondj == 0 .or. nbondj == -1) )   lsend_bdy(ib_bdy,igrd,4) = .true. 
     892               IF( ij == jpj-1 .and. (nbondj == 0 .or. nbondj == -1) )   lsend_bdy(ib_bdy,igrd,4,ir) = .true. 
    891893               ! 
    892894               ! check if point has to be received from a neighbour 
    893895               ! W neighbour and on the outter left  side 
    894                IF( ii == 1   .and. (nbondi == 0 .or. nbondi ==  1) )    lrecv_bdy(ib_bdy,igrd,1) = .true. 
     896               IF( ii == 1     .and. (nbondi == 0 .or. nbondi ==  1) )   lrecv_bdy(ib_bdy,igrd,1,ir) = .true. 
    895897               ! E neighbour and on the outter right side 
    896                IF( ii == jpi .and. (nbondi == 0 .or. nbondi == -1) )    lrecv_bdy(ib_bdy,igrd,2) = .true. 
     898               IF( ii == jpi   .and. (nbondi == 0 .or. nbondi == -1) )   lrecv_bdy(ib_bdy,igrd,2,ir) = .true. 
    897899               ! S neighbour and on the outter down side 
    898                IF( ij == 1   .and. (nbondj == 0 .or. nbondj ==  1) )    lrecv_bdy(ib_bdy,igrd,3) = .true. 
     900               IF( ij == 1     .and. (nbondj == 0 .or. nbondj ==  1) )   lrecv_bdy(ib_bdy,igrd,3,ir) = .true. 
    899901               ! N neighbour and on the outter up   side 
    900                IF( ij == jpj .and. (nbondj == 0 .or. nbondj == -1) )    lrecv_bdy(ib_bdy,igrd,4) = .true. 
     902               IF( ij == jpj   .and. (nbondj == 0 .or. nbondj == -1) )   lrecv_bdy(ib_bdy,igrd,4,ir) = .true. 
    901903               ! 
    902904            END DO 
     
    907909         DO igrd = 1, jpbgrd 
    908910            DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd) 
    909                nbr => idx_bdy(ib_bdy)%nbr(ib,igrd) 
    910                idx_bdy(ib_bdy)%nbw(ib,igrd) = 1.- TANH( REAL( nbr - 1 ) *0.5 )      ! tanh formulation 
    911 !               idx_bdy(ib_bdy)%nbw(ib,igrd) = (REAL(nn_rimwidth(ib_bdy)+1-nbr)/REAL(nn_rimwidth(ib_bdy)))**2.  ! quadratic 
    912 !               idx_bdy(ib_bdy)%nbw(ib,igrd) =  REAL(nn_rimwidth(ib_bdy)+1-nbr)/REAL(nn_rimwidth(ib_bdy))       ! linear 
     911               ir = MAX( 1, idx_bdy(ib_bdy)%nbr(ib,igrd) )   ! both rim 0 and rim 1 have the same weights 
     912               idx_bdy(ib_bdy)%nbw(ib,igrd) = 1.- TANH( REAL( ir - 1 ) *0.5 )      ! tanh formulation 
     913!               idx_bdy(ib_bdy)%nbw(ib,igrd) = (REAL(nn_rimwidth(ib_bdy)+1-ir)/REAL(nn_rimwidth(ib_bdy)))**2.  ! quadratic 
     914!               idx_bdy(ib_bdy)%nbw(ib,igrd) =  REAL(nn_rimwidth(ib_bdy)+1-ir)/REAL(nn_rimwidth(ib_bdy))       ! linear 
    913915            END DO 
    914916         END DO  
     
    918920         DO igrd = 1, jpbgrd 
    919921            DO ib = 1, idx_bdy(ib_bdy)%nblen(igrd) 
    920                nbr => idx_bdy(ib_bdy)%nbr(ib,igrd) 
     922               ir = MAX( 1, idx_bdy(ib_bdy)%nbr(ib,igrd) )   ! both rim 0 and rim 1 have the same damping coefficients 
    921923               idx_bdy(ib_bdy)%nbd(ib,igrd) = 1. / ( rn_time_dmp(ib_bdy) * rday ) &  
    922                & *(REAL(nn_rimwidth(ib_bdy)+1-nbr)/REAL(nn_rimwidth(ib_bdy)))**2.   ! quadratic 
     924               & *(REAL(nn_rimwidth(ib_bdy)+1-ir)/REAL(nn_rimwidth(ib_bdy)))**2.   ! quadratic 
    923925               idx_bdy(ib_bdy)%nbdout(ib,igrd) = 1. / ( rn_time_dmp_out(ib_bdy) * rday ) &  
    924                & *(REAL(nn_rimwidth(ib_bdy)+1-nbr)/REAL(nn_rimwidth(ib_bdy)))**2.   ! quadratic 
     926               & *(REAL(nn_rimwidth(ib_bdy)+1-ir)/REAL(nn_rimwidth(ib_bdy)))**2.   ! quadratic 
    925927            END DO 
    926928         END DO  
     
    931933      ! Initialise masks and find normal/tangential directions 
    932934      ! ------------------------------------------------------ 
     935 
     936      ! ------------------------------------------ 
     937      ! handel rim0, do as if rim 1 was free ocean 
     938      ! ------------------------------------------ 
     939 
     940      ztmask(:,:) = tmask(:,:,1)   ;   zumask(:,:) = umask(:,:,1)   ;   zvmask(:,:) = vmask(:,:,1) 
     941      ! For the flagu/flagv calculation below we require a version of fmask without 
     942      ! the land boundary condition (shlat) included: 
     943      DO ij = 2, jpjm1 
     944         DO ii = 2, jpim1 
     945            zfmask(ii,ij) = ztmask(ii,ij  ) * ztmask(ii+1,ij  )   & 
     946           &              * ztmask(ii,ij+1) * ztmask(ii+1,ij+1) 
     947         END DO       
     948      END DO 
     949      CALL lbc_lnk( 'bdyini', zfmask, 'F', 1. ) 
    933950 
    934951      ! Read global 2D mask at T-points: bdytmask 
     
    940957 
    941958      ! Derive mask on U and V grid from mask on T grid 
    942  
    943       bdyumask(:,:) = 0._wp 
    944       bdyvmask(:,:) = 0._wp 
    945959      DO ij = 1, jpjm1 
    946960         DO ii = 1, jpim1 
    947             bdyumask(ii,ij) = bdytmask(ii,ij) * bdytmask(ii+1, ij ) 
     961            bdyumask(ii,ij) = bdytmask(ii,ij) * bdytmask(ii+1,ij ) 
    948962            bdyvmask(ii,ij) = bdytmask(ii,ij) * bdytmask(ii  ,ij+1)   
    949963         END DO 
    950964      END DO 
    951       CALL lbc_lnk_multi( 'bdyini', bdyumask, 'U', 1. , bdyvmask, 'V', 1. )   ! Lateral boundary cond.  
    952  
    953       ! bdy masks are now set to zero on boundary points: 
    954       ! 
    955       igrd = 1 
     965      CALL lbc_lnk_multi( 'bdyini', bdyumask, 'U', 1., bdyvmask, 'V', 1. )   ! Lateral boundary cond. 
     966 
     967      ! bdy masks are now set to zero on rim 0 points: 
    956968      DO ib_bdy = 1, nb_bdy 
    957         DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)       
    958           bdytmask(idx_bdy(ib_bdy)%nbi(ib,igrd), idx_bdy(ib_bdy)%nbj(ib,igrd)) = 0._wp 
    959         END DO 
     969         DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(1)   ! extent of rim 0 
     970            bdytmask(idx_bdy(ib_bdy)%nbi(ib,1), idx_bdy(ib_bdy)%nbj(ib,1)) = 0._wp 
     971         END DO 
     972         DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(2)   ! extent of rim 0 
     973            bdyumask(idx_bdy(ib_bdy)%nbi(ib,2), idx_bdy(ib_bdy)%nbj(ib,2)) = 0._wp 
     974         END DO 
     975         DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(3)   ! extent of rim 0 
     976            bdyvmask(idx_bdy(ib_bdy)%nbi(ib,3), idx_bdy(ib_bdy)%nbj(ib,3)) = 0._wp 
     977         END DO 
    960978      END DO 
    961       ! 
    962       igrd = 2 
     979 
     980      CALL bdy_rim_treat( zumask, zvmask, zfmask, .true. )   ! compute flagu, flagv, ntreat on rim 0 
     981 
     982      ! ------------------------------------ 
     983      ! handel rim1, do as if rim 0 was land 
     984      ! ------------------------------------ 
     985 
    963986      DO ib_bdy = 1, nb_bdy 
    964         DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
    965           bdyumask(idx_bdy(ib_bdy)%nbi(ib,igrd), idx_bdy(ib_bdy)%nbj(ib,igrd)) = 0._wp 
    966         END DO 
     987         DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(1)   ! extent of rim 0 
     988            ztmask(idx_bdy(ib_bdy)%nbi(ib,1), idx_bdy(ib_bdy)%nbj(ib,1)) = 0._wp 
     989         END DO 
     990         DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(2)   ! extent of rim 0 
     991            zumask(idx_bdy(ib_bdy)%nbi(ib,2), idx_bdy(ib_bdy)%nbj(ib,2)) = 0._wp 
     992         END DO 
     993         DO ib = 1, idx_bdy(ib_bdy)%nblenrim0(3)   ! extent of rim 0 
     994            zvmask(idx_bdy(ib_bdy)%nbi(ib,3), idx_bdy(ib_bdy)%nbj(ib,3)) = 0._wp 
     995         END DO 
    967996      END DO 
    968       ! 
    969       igrd = 3 
    970       DO ib_bdy = 1, nb_bdy 
    971         DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
    972           bdyvmask(idx_bdy(ib_bdy)%nbi(ib,igrd), idx_bdy(ib_bdy)%nbj(ib,igrd)) = 0._wp 
    973         END DO 
    974       END DO 
    975  
    976       ! For the flagu/flagv calculation below we require a version of fmask without 
    977       ! the land boundary condition (shlat) included: 
    978       zfmask(:,:) = 0 
     997 
     998      ! Recompute zfmask 
    979999      DO ij = 2, jpjm1 
    9801000         DO ii = 2, jpim1 
    981             zfmask(ii,ij) = tmask(ii,ij  ,1) * tmask(ii+1,ij  ,1)   & 
    982            &              * tmask(ii,ij+1,1) * tmask(ii+1,ij+1,1) 
     1001            zfmask(ii,ij) = ztmask(ii,ij  ) * ztmask(ii+1,ij  )   & 
     1002           &              * ztmask(ii,ij+1) * ztmask(ii+1,ij+1) 
    9831003         END DO       
    9841004      END DO 
    985  
    986       ! Lateral boundary conditions 
    987       CALL lbc_lnk( 'bdyini', zfmask, 'F', 1. )  
    988       CALL lbc_lnk_multi( 'bdyini', bdyumask, 'U', 1. , bdyvmask, 'V', 1., bdytmask, 'T', 1. ) 
     1005      CALL lbc_lnk( 'bdyini', zfmask, 'F', 1. ) 
     1006 
     1007      ! bdy masks are now set to zero on rim1 points: 
     1008      DO ib_bdy = 1, nb_bdy 
     1009         DO ib = idx_bdy(ib_bdy)%nblenrim0(1) + 1,  idx_bdy(ib_bdy)%nblenrim(1)   ! extent of rim 1 
     1010            bdytmask(idx_bdy(ib_bdy)%nbi(ib,1), idx_bdy(ib_bdy)%nbj(ib,1)) = 0._wp 
     1011         END DO 
     1012         DO ib = idx_bdy(ib_bdy)%nblenrim0(2) + 1,  idx_bdy(ib_bdy)%nblenrim(2)   ! extent of rim 1 
     1013            bdyumask(idx_bdy(ib_bdy)%nbi(ib,2), idx_bdy(ib_bdy)%nbj(ib,2)) = 0._wp 
     1014         END DO 
     1015         DO ib = idx_bdy(ib_bdy)%nblenrim0(3) + 1,  idx_bdy(ib_bdy)%nblenrim(3)   ! extent of rim 1 
     1016            bdyvmask(idx_bdy(ib_bdy)%nbi(ib,3), idx_bdy(ib_bdy)%nbj(ib,3)) = 0._wp 
     1017         END DO 
     1018      END DO 
     1019 
     1020      CALL bdy_rim_treat( zumask, zvmask, zfmask, .false. )   ! compute flagu, flagv, ntreat on rim 1 
     1021 
     1022 
     1023      ! 
     1024      ! Check which boundaries might need communication 
     1025      ALLOCATE( lsend_bdyint(nb_bdy,jpbgrd,4,0:1), lrecv_bdyint(nb_bdy,jpbgrd,4,0:1) ) 
     1026      lsend_bdyint(:,:,:,:) = .false. 
     1027      lrecv_bdyint(:,:,:,:) = .false.  
     1028      ALLOCATE( lsend_bdyext(nb_bdy,jpbgrd,4,0:1), lrecv_bdyext(nb_bdy,jpbgrd,4,0:1) ) 
     1029      lsend_bdyext(:,:,:,:) = .false. 
     1030      lrecv_bdyext(:,:,:,:) = .false. 
     1031      ! 
     1032      DO igrd = 1, jpbgrd 
     1033         DO ib_bdy = 1, nb_bdy 
     1034            DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
     1035               ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
     1036               ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 
     1037               IF( ib .LE. idx_bdy(ib_bdy)%nblenrim0(igrd) ) THEN   ;   ir = 0 
     1038               ELSE                                                 ;   ir = 1 
     1039               END IF 
     1040               flagu = NINT(idx_bdy(ib_bdy)%flagu(ib,igrd)) 
     1041               flagv = NINT(idx_bdy(ib_bdy)%flagv(ib,igrd)) 
     1042               iibe = ii - flagu   ! neighbouring point towards the exterior of the computational domain 
     1043               ijbe = ij - flagv 
     1044               iibi = ii + flagu   ! neighbouring point towards the interior of the computational domain 
     1045               ijbi = ij + flagv 
     1046               CALL find_neib( ii, ij, idx_bdy(ib_bdy)%ntreat(ib,igrd), ii1, ij1, ii2, ij2, ii3, ij3 )   ! free ocean neighbours 
     1047               ! 
     1048               ! search neighbour in the  west/east  direction 
     1049               ! Rim is on the halo and computed ocean is towards exterior of mpi domain   
     1050               !      <--    (o exterior)     -->   
     1051               ! (1)  o|x         OR    (2)   x|o 
     1052               !       |___                 ___|  
     1053               IF( iibi == 0     .OR. ii1 == 0     .OR. ii2 == 0     .OR. ii3 == 0     )   lrecv_bdyint(ib_bdy,igrd,1,ir) = .true. 
     1054               IF( iibi == jpi+1 .OR. ii1 == jpi+1 .OR. ii2 == jpi+1 .OR. ii3 == jpi+1 )   lrecv_bdyint(ib_bdy,igrd,2,ir) = .true.   
     1055               IF( iibe == 0     )   lrecv_bdyext(ib_bdy,igrd,1,ir) = .true. 
     1056               IF( iibe == jpi+1 )   lrecv_bdyext(ib_bdy,igrd,2,ir) = .true.   
     1057               ! Check if neighbour has its rim parallel to its mpi subdomain border and located next to its halo 
     1058               ! :¨¨¨¨¨|¨¨-->    |                                             |    <--¨¨|¨¨¨¨¨:  
     1059               ! :     |  x:o    |    neighbour limited by ... would need o    |    o:x  |     : 
     1060               ! :.....|_._:_____|   (1) W neighbour         E neighbour (2)   |_____:_._|.....: 
     1061               IF( ii == 2     .AND. ( nbondi ==  1 .OR. nbondi == 0 ) .AND. & 
     1062                 & ( iibi == 3     .OR. ii1 == 3     .OR. ii2 == 3     .OR. ii3 == 3    ) )   lsend_bdyint(ib_bdy,igrd,1,ir) = .true. 
     1063               IF( ii == jpi-1 .AND. ( nbondi == -1 .OR. nbondi == 0 ) .AND. & 
     1064                 & ( iibi == jpi-2 .OR. ii1 == jpi-2 .OR. ii2 == jpi-2 .OR. ii3 == jpi-2) )   lsend_bdyint(ib_bdy,igrd,2,ir) = .true. 
     1065               IF( ii == 2     .AND. ( nbondi ==  1 .OR. nbondi == 0 ) .AND. iibe == 3    )   lsend_bdyext(ib_bdy,igrd,1,ir) = .true. 
     1066               IF( ii == jpi-1 .AND. ( nbondi == -1 .OR. nbondi == 0 ) .AND. iibe == jpi-2)   lsend_bdyext(ib_bdy,igrd,2,ir) = .true. 
     1067               ! 
     1068               ! search neighbour in the north/south direction    
     1069               ! Rim is on the halo and computed ocean is towards exterior of mpi domain 
     1070               !(3)   |       |         ^   ___o___      
     1071               !  |   |___x___|   OR    |  |   x   | 
     1072               !  v       o           (4)  |       | 
     1073               IF( ijbi == 0     .OR. ij1 == 0     .OR. ij2 == 0     .OR. ij3 == 0     )   lrecv_bdyint(ib_bdy,igrd,3,ir) = .true. 
     1074               IF( ijbi == jpj+1 .OR. ij1 == jpj+1 .OR. ij2 == jpj+1 .OR. ij3 == jpj+1 )   lrecv_bdyint(ib_bdy,igrd,4,ir) = .true. 
     1075               IF( ijbe == 0     )   lrecv_bdyext(ib_bdy,igrd,3,ir) = .true. 
     1076               IF( ijbe == jpj+1 )   lrecv_bdyext(ib_bdy,igrd,4,ir) = .true. 
     1077               ! Check if neighbour has its rim parallel to its mpi subdomain     _________  border and next to its halo 
     1078               !   ^  |    o    |                                                :         :  
     1079               !   |  |¨¨¨¨x¨¨¨¨|   neighbour limited by ... would need o     |  |....x....| 
     1080               !      :_________:  (3) S neighbour          N neighbour (4)   v  |    o    |    
     1081               IF( ij == 2     .AND. ( nbondj ==  1 .OR. nbondj == 0 ) .AND. & 
     1082                 & ( ijbi == 3     .OR. ij1 == 3     .OR. ij2 == 3     .OR. ij3 == 3    ) )   lsend_bdyint(ib_bdy,igrd,3,ir) = .true. 
     1083               IF( ij == jpj-1 .AND. ( nbondj == -1 .OR. nbondj == 0 ) .AND. & 
     1084                 & ( ijbi == jpj-2 .OR. ij1 == jpj-2 .OR. ij2 == jpj-2 .OR. ij3 == jpj-2) )   lsend_bdyint(ib_bdy,igrd,4,ir) = .true. 
     1085               IF( ij == 2     .AND. ( nbondj ==  1 .OR. nbondj == 0 ) .AND. ijbe == 3    )   lsend_bdyext(ib_bdy,igrd,3,ir) = .true. 
     1086               IF( ij == jpj-1 .AND. ( nbondj == -1 .OR. nbondj == 0 ) .AND. ijbe == jpj-2)   lsend_bdyext(ib_bdy,igrd,4,ir) = .true. 
     1087            END DO 
     1088         END DO 
     1089      END DO 
     1090 
     1091      DO ib_bdy = 1,nb_bdy 
     1092         IF( cn_dyn2d(ib_bdy) == 'orlanski' .OR. cn_dyn2d(ib_bdy) == 'orlanski_npo' .OR. & 
     1093           & cn_dyn3d(ib_bdy) == 'orlanski' .OR. cn_dyn3d(ib_bdy) == 'orlanski_npo' .OR. & 
     1094           & cn_tra(ib_bdy)   == 'orlanski' .OR. cn_tra(ib_bdy)   == 'orlanski_npo'      ) THEN 
     1095            DO igrd = 1, jpbgrd 
     1096               DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
     1097                  ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
     1098                  ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 
     1099                  IF(  mig(ii) > 2 .AND. mig(ii) < jpiglo-2 .AND. mjg(ij) > 2 .AND. mjg(ij) < jpjglo-2  ) THEN 
     1100                     WRITE(ctmp1,*) ' Orlanski is not safe when the open boundaries are on the interior of the computational domain' 
     1101                     WRITE(ctmp2,*) ' ========== ' 
     1102                     CALL ctl_stop( ' ', ctmp1, ctmp2, ' ' ) 
     1103                  END IF 
     1104               END DO 
     1105            END DO 
     1106         END IF 
     1107      END DO 
     1108 
     1109      ! 
     1110      ! Tidy up 
     1111      !-------- 
     1112      IF( nb_bdy>0 )   DEALLOCATE( nbidta, nbjdta, nbrdta ) 
     1113      ! 
     1114   END SUBROUTINE bdy_segs 
     1115 
     1116 
     1117   SUBROUTINE bdy_rim_treat( zumask, zvmask, zfmask, lrim0 ) 
     1118      !!---------------------------------------------------------------------- 
     1119      !!                 ***  ROUTINE bdy_rim_treat  *** 
     1120      !! 
     1121      !! ** Purpose :   Initialize structures ( flagu, flagv, ntreat ) indicating how rim points 
     1122      !!                  are to be handeled in the boundary condition treatment 
     1123      !! 
     1124      !! ** Method  :   - to handel rim 0 zmasks must indicate ocean points      (set at one on rim 0 and rim 1 and interior) 
     1125      !!                            and bdymasks must be set at 0 on rim 0       (set at one on rim 1 and interior) 
     1126      !!                    (as if rim 1 was free ocean) 
     1127      !!                - to handel rim 1 zmasks must be set at 0 on rim 0       (set at one on rim 1 and interior) 
     1128      !!                            and bdymasks must indicate free ocean points (set at one on interior) 
     1129      !!                    (as if rim 0 was land) 
     1130      !!                - we can then check in which direction the interior of the computational domain is with the difference 
     1131      !!                         mask array values on both sides to compute flagu and flagv 
     1132      !!                - and look at the ocean neighbours to compute ntreat 
     1133      !!---------------------------------------------------------------------- 
     1134      REAL(wp), TARGET, DIMENSION(jpi,jpj), INTENT (in   ) :: zfmask   ! temporary fmask excluding coastal boundary condition (shlat) 
     1135      REAL(wp), TARGET, DIMENSION(jpi,jpj), INTENT (in   ) :: zumask, zvmask   ! temporary t/u/v mask array 
     1136      LOGICAL                             , INTENT (in   ) :: lrim0    ! .true. -> rim 0   .false. -> rim 1 
     1137      INTEGER  ::   ib_bdy, ii, ij, igrd, ib, icount       ! dummy loop indices 
     1138      INTEGER  ::   i_offset, j_offset, inn                ! local integer 
     1139      INTEGER  ::   ibeg, iend                             ! local integer 
     1140      LOGICAL  ::   llnon, llson, llean, llwen             ! local logicals indicating the presence of a ocean neighbour 
     1141      REAL(wp), POINTER, DIMENSION(:,:)       ::   pmask   ! pointer to 2D mask fields 
     1142      REAL(wp) ::   zefl, zwfl, znfl, zsfl                 ! local scalars 
     1143      CHARACTER(LEN=1), DIMENSION(jpbgrd)     ::   cgrid 
     1144      REAL(wp)        , DIMENSION(jpi,jpj)    ::   ztmp 
     1145      !!---------------------------------------------------------------------- 
     1146 
     1147      cgrid = (/'t','u','v'/) 
     1148 
    9891149      DO ib_bdy = 1, nb_bdy       ! Indices and directions of rim velocity components 
    990          icount = 0  
    9911150 
    9921151         ! Calculate relationship of U direction to the local orientation of the boundary 
     
    9941153         ! flagu =  0 : u is tangential 
    9951154         ! flagu =  1 : u is normal to the boundary and is direction is inward 
    996    
    9971155         DO igrd = 1, jpbgrd  
    9981156            SELECT CASE( igrd ) 
    999                CASE( 1 )   ;   pmask => umask   (:,:,1)   ;   i_offset = 0 
    1000                CASE( 2 )   ;   pmask => bdytmask(:,:)     ;   i_offset = 1 
    1001                CASE( 3 )   ;   pmask => zfmask  (:,:)     ;   i_offset = 0 
     1157               CASE( 1 )   ;   pmask => zumask     ;   i_offset = 0 
     1158               CASE( 2 )   ;   pmask => bdytmask   ;   i_offset = 1 
     1159               CASE( 3 )   ;   pmask => zfmask     ;   i_offset = 0 
    10021160            END SELECT  
    10031161            icount = 0 
    10041162            ztmp(:,:) = 0._wp 
    1005             DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)   
     1163            IF( lrim0 ) THEN   ! extent of rim 0 
     1164               ibeg = 1                                     ;   iend = idx_bdy(ib_bdy)%nblenrim0(igrd) 
     1165            ELSE               ! extent of rim 1 
     1166               ibeg = idx_bdy(ib_bdy)%nblenrim0(igrd) + 1   ;   iend = idx_bdy(ib_bdy)%nblenrim(igrd) 
     1167            END IF 
     1168            DO ib = ibeg, iend  
    10061169               ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
    10071170               ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 
    10081171               IF( ii == 1 .OR. ii == jpi .OR. ij == 1 .OR. ij == jpj )  CYCLE 
    1009                zefl = pmask(ii+i_offset-1,ij) 
    1010                zwfl = pmask(ii+i_offset  ,ij) 
     1172               zwfl = pmask(ii+i_offset-1,ij) 
     1173               zefl = pmask(ii+i_offset  ,ij) 
    10111174               ! This error check only works if you are using the bdyXmask arrays 
    10121175               IF( i_offset == 1 .and. zefl + zwfl == 2. ) THEN 
     
    10141177                  IF(lwp) WRITE(numout,*) 'Problem with igrd = ',igrd,' at (global) nbi, nbj : ',mig(ii),mjg(ij) 
    10151178               ELSE 
    1016                   ztmp(ii,ij) = -zefl + zwfl 
     1179                  ztmp(ii,ij) = -zwfl + zefl 
    10171180               ENDIF 
    10181181            END DO 
     
    10241187            ENDIF  
    10251188            CALL lbc_lnk( 'bdyini', ztmp, 'T', 1. )  
    1026             DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)   
     1189            DO ib = ibeg, iend 
    10271190               ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
    10281191               ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 
     
    10351198         ! flagv =  0 : v is tangential 
    10361199         ! flagv =  1 : v is normal to the boundary and is direction is inward 
    1037  
    10381200         DO igrd = 1, jpbgrd  
    10391201            SELECT CASE( igrd ) 
    1040                CASE( 1 )   ;   pmask => vmask (:,:,1)   ;   j_offset = 0 
    1041                CASE( 2 )   ;   pmask => zfmask(:,:)     ;   j_offset = 0 
    1042                CASE( 3 )   ;   pmask => bdytmask        ;   j_offset = 1 
     1202               CASE( 1 )   ;   pmask => zvmask     ;   j_offset = 0 
     1203               CASE( 2 )   ;   pmask => zfmask     ;   j_offset = 0 
     1204               CASE( 3 )   ;   pmask => bdytmask   ;   j_offset = 1 
    10431205            END SELECT  
    10441206            icount = 0 
    10451207            ztmp(:,:) = 0._wp 
    1046             DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)   
     1208            IF( lrim0 ) THEN   ! extent of rim 0 
     1209               ibeg = 1                                     ;   iend = idx_bdy(ib_bdy)%nblenrim0(igrd) 
     1210            ELSE               ! extent of rim 1 
     1211               ibeg = idx_bdy(ib_bdy)%nblenrim0(igrd) + 1   ;   iend = idx_bdy(ib_bdy)%nblenrim(igrd) 
     1212            END IF 
     1213            DO ib = ibeg, iend 
    10471214               ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
    10481215               ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 
    10491216               IF( ii == 1 .OR. ii == jpi .OR. ij == 1 .OR. ij == jpj )  CYCLE 
    1050                znfl = pmask(ii,ij+j_offset-1) 
    1051                zsfl = pmask(ii,ij+j_offset  ) 
     1217               zsfl = pmask(ii,ij+j_offset-1) 
     1218               znfl = pmask(ii,ij+j_offset  ) 
    10521219               ! This error check only works if you are using the bdyXmask arrays 
    10531220               IF( j_offset == 1 .and. znfl + zsfl == 2. ) THEN 
     
    10551222                  icount = icount + 1 
    10561223               ELSE 
    1057                   ztmp(ii,ij) = -znfl + zsfl 
     1224                  ztmp(ii,ij) = -zsfl + znfl 
    10581225               END IF 
    10591226            END DO 
     
    10651232            ENDIF 
    10661233            CALL lbc_lnk( 'bdyini', ztmp, 'T', 1. )  
    1067             DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd)   
     1234            DO ib = ibeg, iend 
    10681235               ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
    10691236               ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 
     
    10721239         END DO 
    10731240         ! 
    1074       END DO 
     1241      END DO ! ib_bdy 
    10751242       
    10761243      DO ib_bdy = 1, nb_bdy 
     
    10821249            END SELECT 
    10831250            ztmp(:,:) = 0._wp 
    1084             DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
     1251            IF( lrim0 ) THEN   ! extent of rim 0 
     1252               ibeg = 1                                     ;   iend = idx_bdy(ib_bdy)%nblenrim0(igrd) 
     1253            ELSE               ! extent of rim 1 
     1254               ibeg = idx_bdy(ib_bdy)%nblenrim0(igrd) + 1   ;   iend = idx_bdy(ib_bdy)%nblenrim(igrd) 
     1255            END IF 
     1256            DO ib = ibeg, iend 
    10851257               ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
    10861258               ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 
    1087                IF( ii == 1 .OR. ii == jpi .OR. ij == 1 .OR. ij == jpj )   CYCLE  
    1088                llean = pmask(ii+1,ij  ) == 1. 
    1089                llwen = pmask(ii-1,ij  ) == 1. 
    1090                llnon = pmask(ii  ,ij+1) == 1. 
    1091                llson = pmask(ii  ,ij-1) == 1. 
     1259               IF( ii == 1 .OR. ii == jpi .OR. ij == 1 .OR. ij == jpj )   CYCLE 
     1260               llnon = pmask(ii  ,ij+1) == 1.   
     1261               llson = pmask(ii  ,ij-1) == 1.  
     1262               llean = pmask(ii+1,ij  ) == 1.  
     1263               llwen = pmask(ii-1,ij  ) == 1.  
    10921264               inn  = COUNT( (/ llnon, llson, llean, llwen /) ) 
    10931265               IF( inn == 0 ) THEN   ! no neighbours -> interior of a corner  or  cluster of rim points 
     
    11011273                  ELSE 
    11021274                     WRITE(ctmp1,*) ' E R R O R : Problem with  ',cgrid(igrd) ,' grid point', ii, ij,   & 
    1103                           ' on boundary set ', ib_bdy, ' has no free ocean neighbour' 
    1104                      WRITE(ctmp2,*) ' There seems to be a cluster of rim points.' 
     1275                       ' on boundary set ', ib_bdy, ' has no free ocean neighbour' 
     1276                     IF( lrim0 ) THEN 
     1277                        WRITE(ctmp2,*) ' There seems to be a cluster of rim 0 points.' 
     1278                     ELSE 
     1279                        WRITE(ctmp2,*) ' There seems to be a cluster of rim 1 points.' 
     1280                     END IF 
    11051281                     WRITE(ctmp3,*) ' ========== ' 
    11061282                     CALL ctl_stop( ' ', ctmp1, ctmp2, ctmp3, ' ' ) 
     
    11421318            END DO 
    11431319            CALL lbc_lnk( 'bdyini', ztmp, 'T', 1. ) 
    1144             DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
     1320            DO ib = ibeg, iend 
    11451321               ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
    11461322               ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 
     
    11501326      END DO 
    11511327 
    1152  
    1153       ALLOCATE( lsend_bdyint(nb_bdy,jpbgrd,4), lrecv_bdyint(nb_bdy,jpbgrd,4) ) 
    1154       lsend_bdyint(:,:,:) = .false. 
    1155       lrecv_bdyint(:,:,:) = .false.  
    1156       ALLOCATE( lsend_bdyext(nb_bdy,jpbgrd,4), lrecv_bdyext(nb_bdy,jpbgrd,4) ) 
    1157       lsend_bdyext(:,:,:) = .false. 
    1158       lrecv_bdyext(:,:,:) = .false. 
    1159       ! 
    1160       ! Check which boundaries might need communication 
    1161       DO igrd = 1, jpbgrd 
    1162          DO ib_bdy = 1, nb_bdy 
    1163             DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
    1164                ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
    1165                ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 
    1166                flagu = NINT(idx_bdy(ib_bdy)%flagu(ib,igrd)) 
    1167                flagv = NINT(idx_bdy(ib_bdy)%flagv(ib,igrd)) 
    1168                iibe = ii - flagu   ! neighbouring point towards the exterior of the computational domain 
    1169                ijbe = ij - flagv 
    1170                iibi = ii + flagu   ! neighbouring point towards the interior of the computational domain 
    1171                ijbi = ij + flagv 
    1172                CALL find_neib( ii, ij, idx_bdy(ib_bdy)%ntreat(ib,igrd), ii1, ij1, ii2, ij2, ii3, ij3 ) 
    1173                ! 
    1174                ! search neighbour in the  west/east  direction 
    1175                ! Rim is on the halo and computed ocean is towards exterior of mpi domain   
    1176                !      <--    (o exterior)     -->   
    1177                ! (1)  o|x         OR    (2)   x|o 
    1178                !       |___                 ___|  
    1179                IF( iibi == 0     .OR. ii1 == 0     .OR. ii2 == 0     .OR. ii3 == 0     )   lrecv_bdyint(ib_bdy,igrd,1) = .true. 
    1180                IF( iibi == jpi+1 .OR. ii1 == jpi+1 .OR. ii2 == jpi+1 .OR. ii3 == jpi+1 )   lrecv_bdyint(ib_bdy,igrd,2) = .true.   
    1181                IF( iibe == 0     )   lrecv_bdyext(ib_bdy,igrd,1) = .true. 
    1182                IF( iibe == jpi+1 )   lrecv_bdyext(ib_bdy,igrd,2) = .true.   
    1183                ! Check if neighbour has its rim parallel to its mpi subdomain border and located next to its halo 
    1184                ! :¨¨¨¨¨|¨¨-->    |                                             |    <--¨¨|¨¨¨¨¨:  
    1185                ! :     |  x:o    |    neighbour limited by ... would need o    |    o:x  |     : 
    1186                ! :.....|_._:_____|   (1) W neighbour         E neighbour (2)   |_____:_._|.....: 
    1187                IF( ii == 2     .AND. ( nbondi ==  1 .OR. nbondi == 0 ) .AND. & 
    1188                     & ( iibi == 3     .OR. ii1 == 3     .OR. ii2 == 3     .OR. ii3 == 3    ) )   lsend_bdyint(ib_bdy,igrd,1) = .true. 
    1189                IF( ii == jpi-1 .AND. ( nbondi == -1 .OR. nbondi == 0 ) .AND. & 
    1190                     & ( iibi == jpi-2 .OR. ii1 == jpi-2 .OR. ii2 == jpi-2 .OR. ii3 == jpi-2) )   lsend_bdyint(ib_bdy,igrd,2) = .true. 
    1191                IF( ii == 2     .AND. ( nbondi ==  1 .OR. nbondi == 0 ) .AND. iibe == 3     )   lsend_bdyext(ib_bdy,igrd,1) = .true. 
    1192                IF( ii == jpi-1 .AND. ( nbondi == -1 .OR. nbondi == 0 ) .AND. iibe == jpi-2 )   lsend_bdyext(ib_bdy,igrd,2) = .true. 
    1193                ! 
    1194                ! search neighbour in the north/south direction    
    1195                ! Rim is on the halo and computed ocean is towards exterior of mpi domain 
    1196                !(3)   |       |         ^   ___o___      
    1197                !  |   |___x___|   OR    |  |   x   | 
    1198                !  v       o           (4)  |       | 
    1199                IF( ijbi == 0     .OR. ij1 == 0     .OR. ij2 == 0     .OR. ij3 == 0     )   lrecv_bdyint(ib_bdy,igrd,3) = .true. 
    1200                IF( ijbi == jpj+1 .OR. ij1 == jpj+1 .OR. ij2 == jpj+1 .OR. ij3 == jpj+1 )   lrecv_bdyint(ib_bdy,igrd,4) = .true. 
    1201                IF( ijbe == 0     )   lrecv_bdyext(ib_bdy,igrd,3) = .true. 
    1202                IF( ijbe == jpj+1 )   lrecv_bdyext(ib_bdy,igrd,4) = .true. 
    1203                ! Check if neighbour has its rim parallel to its mpi subdomain     _________  border and next to its halo 
    1204                !   ^  |    o    |                                                :         :  
    1205                !   |  |¨¨¨¨x¨¨¨¨|   neighbour limited by ... would need o     |  |....x....| 
    1206                !      :_________:  (3) S neighbour          N neighbour (4)   v  |    o    |    
    1207                IF( ij == 2     .AND. ( nbondj ==  1 .OR. nbondj == 0 ) .AND. & 
    1208                     & ( ijbi == 3     .OR. ij1 == 3     .OR. ij2 == 3     .OR. ij3 == 3    ) )   lsend_bdyint(ib_bdy,igrd,3) = .true. 
    1209                IF( ij == jpj-1 .AND. ( nbondj == -1 .OR. nbondj == 0 ) .AND. & 
    1210                     & ( ijbi == jpj-2 .OR. ij1 == jpj-2 .OR. ij2 == jpj-2 .OR. ij3 == jpj-2) )   lsend_bdyint(ib_bdy,igrd,4) = .true. 
    1211                IF( ij == 2     .AND. ( nbondj ==  1 .OR. nbondj == 0 ) .AND. ijbe == 3     )   lsend_bdyext(ib_bdy,igrd,3) = .true. 
    1212                IF( ij == jpj-1 .AND. ( nbondj == -1 .OR. nbondj == 0 ) .AND. ijbe == jpj-2 )   lsend_bdyext(ib_bdy,igrd,4) = .true. 
    1213             END DO 
    1214          END DO 
    1215       END DO 
    1216  
    1217       DO ib_bdy = 1,nb_bdy 
    1218          IF( cn_dyn2d(ib_bdy) == 'orlanski' .OR. cn_dyn2d(ib_bdy) == 'orlanski_npo' .OR. & 
    1219            & cn_dyn3d(ib_bdy) == 'orlanski' .OR. cn_dyn3d(ib_bdy) == 'orlanski_npo' .OR. & 
    1220            & cn_tra(ib_bdy)   == 'orlanski' .OR. cn_tra(ib_bdy)   == 'orlanski_npo'      ) THEN 
    1221             DO igrd = 1, jpbgrd 
    1222                DO ib = 1, idx_bdy(ib_bdy)%nblenrim(igrd) 
    1223                   ii = idx_bdy(ib_bdy)%nbi(ib,igrd) 
    1224                   ij = idx_bdy(ib_bdy)%nbj(ib,igrd) 
    1225                   IF(  mig(ii) > 2 .AND. mig(ii) < jpiglo-2 .AND. mjg(ij) > 2 .AND. mjg(ij) < jpjglo-2  ) THEN 
    1226                      WRITE(ctmp1,*) ' Orlanski can not be used when the open boundaries are on the interior of the computational domain' 
    1227                      WRITE(ctmp2,*) ' ========== ' 
    1228                      CALL ctl_stop( ' ', ctmp1, ctmp2, ' ' ) 
    1229                   END IF 
    1230                END DO 
    1231             END DO 
    1232          END IF 
    1233       END DO 
    1234  
    1235       ! 
    1236       ! Tidy up 
    1237       !-------- 
    1238       IF( nb_bdy>0 )   DEALLOCATE( nbidta, nbjdta, nbrdta ) 
    1239       ! 
    1240    END SUBROUTINE bdy_segs 
    1241  
    1242  
    1243    SUBROUTINE find_neib( ii, ij, itreat, ii1, ij1, ii2, ij2, ii3, ij3 ) 
    1244      !!---------------------------------------------------------------------- 
    1245      !!                 ***  ROUTINE find_neib  *** 
    1246      !! 
    1247      !! ** Purpose :   get ii1, ij1, ii2, ij2, ii3, ij3, the indices of 
    1248      !!               the free ocean neighbours of (ii,ij) for bdy treatment 
    1249      !! 
    1250      !!---------------------------------------------------------------------- 
    1251      INTEGER, INTENT(in   )      ::   ii, ij, itreat 
    1252      INTEGER, INTENT(  out)      ::   ii1, ij1, ii2, ij2, ii3, ij3 
    1253      !!---------------------------------------------------------------------- 
    1254      SELECT CASE( itreat )   ! points that will be used by bdy routines, -1 will be discarded 
    1255         !               !               !     _____     !     _____      
    1256         !  1 |   o      !  2  o   |     !  3 | x        !  4     x |     
    1257         !    |_x_ _     !    _ _x_|     !    |   o      !      o   | 
    1258      CASE( 1 )    ;   ii1 = ii+1   ;   ij1 = ij+1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
    1259      CASE( 2 )    ;   ii1 = ii-1   ;   ij1 = ij+1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
    1260      CASE( 3 )    ;   ii1 = ii+1   ;   ij1 = ij-1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
    1261      CASE( 4 )    ;   ii1 = ii-1   ;   ij1 = ij-1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
    1262         !    |          !         |     !      o        !    ______                   ! or incomplete corner 
    1263         ! 5  | x o      ! 6   o x |     ! 7  __x__      ! 8    x                      !  7  ____ o 
    1264         !    |          !         |     !               !      o                      !         |x___ 
    1265      CASE( 5 )    ;   ii1 = ii+1   ;   ij1 = ij     ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
    1266      CASE( 6 )    ;   ii1 = ii-1   ;   ij1 = ij     ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
    1267      CASE( 7 )    ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
    1268      CASE( 8 )    ;   ii1 = ii     ;   ij1 = ij-1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
    1269         !        o      !        o      !    _____|     !       |_____   
    1270         !  9 ____x o    ! 10   o x___   ! 11     x o    ! 12   o x       
    1271         !         |     !       |       !        o      !        o       
    1272      CASE( 9  )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii+1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1  
    1273      CASE( 10 )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii-1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1 
    1274      CASE( 11 )   ;   ii1 = ii     ;   ij1 = ij-1   ;   ii2 = ii+1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1 
    1275      CASE( 12 )   ;   ii1 = ii     ;   ij1 = ij-1   ;   ii2 = ii-1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1 
    1276         !    |_  o      !        o  _|  !     ¨¨|_|¨¨   !       o          
    1277         ! 13  _| x o    !  14  o x |_   !  15  o x o    ! 16  o x o        
    1278         !    |   o      !        o   |  !        o      !    __|¨|__  
    1279      CASE( 13 )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii+1   ;   ij2 = ij     ;   ii3 = ii     ;   ij3 = ij-1    
    1280      CASE( 14 )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii-1   ;   ij2 = ij     ;   ii3 = ii     ;   ij3 = ij-1  
    1281      CASE( 15 )   ;   ii1 = ii-1   ;   ij1 = ij     ;   ii2 = ii     ;   ij2 = ij-1   ;   ii3 = ii+1   ;   ij3 = ij    
    1282      CASE( 16 )   ;   ii1 = ii-1   ;   ij1 = ij     ;   ii2 = ii     ;   ij2 = ij+1   ;   ii3 = ii+1   ;   ij3 = ij 
    1283      END SELECT 
    1284    END SUBROUTINE find_neib 
    1285  
     1328    END SUBROUTINE bdy_rim_treat 
     1329 
     1330    
     1331    SUBROUTINE find_neib( ii, ij, itreat, ii1, ij1, ii2, ij2, ii3, ij3 ) 
     1332      !!---------------------------------------------------------------------- 
     1333      !!                 ***  ROUTINE find_neib  *** 
     1334      !! 
     1335      !! ** Purpose :   get ii1, ij1, ii2, ij2, ii3, ij3, the indices of 
     1336      !!               the free ocean neighbours of (ii,ij) for bdy treatment 
     1337      !! 
     1338      !! ** Method  :  use itreat input to select a case 
     1339      !!               N.B. ntreat is defined for all bdy points in routine bdy_segs 
     1340      !! 
     1341      !!---------------------------------------------------------------------- 
     1342      INTEGER, INTENT(in   )      ::   ii, ij, itreat 
     1343      INTEGER, INTENT(  out)      ::   ii1, ij1, ii2, ij2, ii3, ij3 
     1344      !!---------------------------------------------------------------------- 
     1345      SELECT CASE( itreat )   ! points that will be used by bdy routines, -1 will be discarded 
     1346         !               !               !     _____     !     _____      
     1347         !  1 |   o      !  2  o   |     !  3 | x        !  4     x |     
     1348         !    |_x_ _     !    _ _x_|     !    |   o      !      o   | 
     1349      CASE( 1 )    ;   ii1 = ii+1   ;   ij1 = ij+1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
     1350      CASE( 2 )    ;   ii1 = ii-1   ;   ij1 = ij+1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
     1351      CASE( 3 )    ;   ii1 = ii+1   ;   ij1 = ij-1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
     1352      CASE( 4 )    ;   ii1 = ii-1   ;   ij1 = ij-1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
     1353         !    |          !         |     !      o        !    ______                   ! or incomplete corner 
     1354         ! 5  | x o      ! 6   o x |     ! 7  __x__      ! 8    x                      !  7  ____ o 
     1355         !    |          !         |     !               !      o                      !         |x___ 
     1356      CASE( 5 )    ;   ii1 = ii+1   ;   ij1 = ij     ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
     1357      CASE( 6 )    ;   ii1 = ii-1   ;   ij1 = ij     ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
     1358      CASE( 7 )    ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
     1359      CASE( 8 )    ;   ii1 = ii     ;   ij1 = ij-1   ;   ii2 = -1     ;   ij2 = -1     ;   ii3 = -1     ;   ij3 = -1 
     1360         !        o      !        o      !    _____|     !       |_____   
     1361         !  9 ____x o    ! 10   o x___   ! 11     x o    ! 12   o x       
     1362         !         |     !       |       !        o      !        o       
     1363      CASE( 9  )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii+1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1  
     1364      CASE( 10 )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii-1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1 
     1365      CASE( 11 )   ;   ii1 = ii     ;   ij1 = ij-1   ;   ii2 = ii+1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1 
     1366      CASE( 12 )   ;   ii1 = ii     ;   ij1 = ij-1   ;   ii2 = ii-1   ;   ij2 = ij     ;   ii3 = -1     ;   ij3 = -1 
     1367         !    |_  o      !        o  _|  !     ¨¨|_|¨¨   !       o          
     1368         ! 13  _| x o    !  14  o x |_   !  15  o x o    ! 16  o x o        
     1369         !    |   o      !        o   |  !        o      !    __|¨|__  
     1370      CASE( 13 )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii+1   ;   ij2 = ij     ;   ii3 = ii     ;   ij3 = ij-1    
     1371      CASE( 14 )   ;   ii1 = ii     ;   ij1 = ij+1   ;   ii2 = ii-1   ;   ij2 = ij     ;   ii3 = ii     ;   ij3 = ij-1  
     1372      CASE( 15 )   ;   ii1 = ii-1   ;   ij1 = ij     ;   ii2 = ii     ;   ij2 = ij-1   ;   ii3 = ii+1   ;   ij3 = ij    
     1373      CASE( 16 )   ;   ii1 = ii-1   ;   ij1 = ij     ;   ii2 = ii     ;   ij2 = ij+1   ;   ii3 = ii+1   ;   ij3 = ij 
     1374      END SELECT 
     1375    END SUBROUTINE find_neib 
     1376     
    12861377 
    12871378   SUBROUTINE bdy_ctl_seg 
  • NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdylib.F90

    r11183 r11191  
    9292 
    9393 
    94    SUBROUTINE bdy_orl( idx, ptb, pta, dta, ll_npo ) 
     94   SUBROUTINE bdy_orl( idx, ptb, pta, dta, lrim0, ll_npo ) 
    9595      !!---------------------------------------------------------------------- 
    9696      !!                 ***  SUBROUTINE bdy_orl  *** 
     
    104104      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   ptb  ! before tracer field 
    105105      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pta  ! tracer trend 
     106      LOGICAL                 , OPTIONAL,  INTENT(in) ::   lrim0   ! indicate if rim 0 is treated 
    106107      LOGICAL,                             INTENT(in) ::   ll_npo  ! switch for NPO version 
    107108      !! 
     
    111112      igrd = 1                       ! Everything is at T-points here 
    112113      ! 
    113       CALL bdy_orlanski_3d( idx, igrd, ptb(:,:,:), pta(:,:,:), dta, ll_npo ) 
     114      CALL bdy_orlanski_3d( idx, igrd, ptb(:,:,:), pta(:,:,:), dta, lrim0, ll_npo ) 
    114115      ! 
    115116   END SUBROUTINE bdy_orl 
    116117 
    117118 
    118    SUBROUTINE bdy_orlanski_2d( idx, igrd, phib, phia, phi_ext, ll_npo ) 
     119   SUBROUTINE bdy_orlanski_2d( idx, igrd, phib, phia, phi_ext, lrim0, ll_npo ) 
    119120      !!---------------------------------------------------------------------- 
    120121      !!                 ***  SUBROUTINE bdy_orlanski_2d  *** 
     
    132133      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   phia     ! model after 2D field (to be updated) 
    133134      REAL(wp), DIMENSION(:)  , INTENT(in   ) ::   phi_ext  ! external forcing data 
     135      LOGICAL, OPTIONAL,        INTENT(in   ) ::   lrim0    ! indicate if rim 0 is treated 
    134136      LOGICAL ,                 INTENT(in   ) ::   ll_npo   ! switch for NPO version 
    135137      ! 
     
    140142      INTEGER  ::   ii_offset, ij_offset                   ! offsets for mask indices 
    141143      INTEGER  ::   flagu, flagv                           ! short cuts 
     144      INTEGER  ::   ibeg, iend                             ! length of rim to be treated (rim 0 or rim 1 or both) 
    142145      REAL(wp) ::   zmask_x, zmask_y1, zmask_y2 
    143146      REAL(wp) ::   zex1, zex2, zey, zey1, zey2 
     
    185188      END SELECT 
    186189      ! 
    187       DO jb = 1, idx%nblenrim(igrd) 
     190      IF( PRESENT(lrim0) ) THEN 
     191         IF( lrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd)   ! rim 0 
     192         ELSE               ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd)    ! rim 1 
     193         END IF 
     194      ELSE                  ;   ibeg = 1                       ;   iend = idx%nblenrim(igrd)    ! both 
     195      END IF 
     196      ! 
     197      DO jb = ibeg, iend 
    188198         ii  = idx%nbi(jb,igrd) 
    189199         ij  = idx%nbj(jb,igrd)  
     
    272282 
    273283 
    274    SUBROUTINE bdy_orlanski_3d( idx, igrd, phib, phia, phi_ext, ll_npo ) 
     284   SUBROUTINE bdy_orlanski_3d( idx, igrd, phib, phia, phi_ext, lrim0, ll_npo ) 
    275285      !!---------------------------------------------------------------------- 
    276286      !!                 ***  SUBROUTINE bdy_orlanski_3d  *** 
     
    288298      REAL(wp), DIMENSION(:,:,:), INTENT(inout) ::   phia     ! model after 3D field (to be updated) 
    289299      REAL(wp), DIMENSION(:,:)  , INTENT(in   ) ::   phi_ext  ! external forcing data 
     300      LOGICAL, OPTIONAL,          INTENT(in   ) ::   lrim0    ! indicate if rim 0 is treated 
    290301      LOGICAL ,                   INTENT(in   ) ::   ll_npo   ! switch for NPO version 
    291302      ! 
     
    296307      INTEGER  ::   ii_offset, ij_offset                   ! offsets for mask indices 
    297308      INTEGER  ::   flagu, flagv                           ! short cuts 
     309      INTEGER  ::   ibeg, iend                             ! length of rim to be treated (rim 0 or rim 1 or both) 
    298310      REAL(wp) ::   zmask_x, zmask_y1, zmask_y2 
    299311      REAL(wp) ::   zex1, zex2, zey, zey1, zey2 
     
    340352         CASE DEFAULT ;   CALL ctl_stop( 'unrecognised value for igrd in bdy_orlanksi_2d' ) 
    341353      END SELECT 
    342  
     354      ! 
     355      IF( PRESENT(lrim0) ) THEN 
     356         IF( lrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd)   ! rim 0 
     357         ELSE               ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd)    ! rim 1 
     358         END IF 
     359      ELSE                  ;   ibeg = 1                       ;   iend = idx%nblenrim(igrd)    ! both 
     360      END IF 
     361      ! 
    343362      DO jk = 1, jpk 
    344363         !             
    345          DO jb = 1, idx%nblenrim(igrd) 
     364         DO jb = ibeg, iend 
    346365            ii  = idx%nbi(jb,igrd) 
    347366            ij  = idx%nbj(jb,igrd)  
     
    430449   END SUBROUTINE bdy_orlanski_3d 
    431450 
    432    SUBROUTINE bdy_nmn( idx, igrd, phia ) 
     451   SUBROUTINE bdy_nmn( idx, igrd, phia, lrim0 ) 
    433452      !!---------------------------------------------------------------------- 
    434453      !!                 ***  SUBROUTINE bdy_nmn  *** 
     
    444463      !!                                                   !      o       
    445464      !!---------------------------------------------------------------------- 
    446       INTEGER,                    INTENT(in)     ::   igrd     ! grid index 
     465      INTEGER,                    INTENT(in   )  ::   igrd     ! grid index 
    447466      REAL(wp), DIMENSION(:,:,:), INTENT(inout)  ::   phia     ! model after 3D field (to be updated), must be masked 
    448       TYPE(OBC_INDEX), INTENT(in) ::   idx  ! OBC indices 
     467      TYPE(OBC_INDEX),            INTENT(in   )  ::   idx      ! OBC indices 
     468      LOGICAL, OPTIONAL,          INTENT(in   )  ::   lrim0    ! indicate if rim 0 is treated 
    449469      !!  
    450470      REAL(wp) ::   zweight 
     
    453473      INTEGER  ::   ii, ij   ! 2D addresses 
    454474      INTEGER  ::   ipkm1    ! size of phia third dimension minus 1 
     475      INTEGER  ::   ibeg, iend                          ! length of rim to be treated (rim 0 or rim 1 or both) 
    455476      INTEGER  ::   ii1, ii2, ii3, ij1, ij2, ij3, itreat 
    456477      !!---------------------------------------------------------------------- 
     
    465486      END SELECT 
    466487      ! 
    467       DO ib = 1, idx%nblenrim(igrd) 
     488      IF( PRESENT(lrim0) ) THEN 
     489         IF( lrim0 ) THEN   ;   ibeg = 1                       ;   iend = idx%nblenrim0(igrd)   ! rim 0 
     490         ELSE               ;   ibeg = idx%nblenrim0(igrd)+1   ;   iend = idx%nblenrim(igrd)    ! rim 1 
     491         END IF 
     492      ELSE                  ;   ibeg = 1                       ;   iend = idx%nblenrim(igrd)    ! both 
     493      END IF 
     494      ! 
     495      DO ib = ibeg, iend 
    468496         ii = idx%nbi(ib,igrd) 
    469497         ij = idx%nbj(ib,igrd) 
    470498         itreat = idx%ntreat(ib,igrd) 
    471          CALL find_neib( ii, ij, itreat, ii1, ij1, ii2, ij2, ii3, ij3 ) 
     499         CALL find_neib( ii, ij, itreat, ii1, ij1, ii2, ij2, ii3, ij3 )   ! find free ocean neighbours 
    472500         SELECT CASE( itreat ) 
    473501         CASE( 1:8 ) 
  • NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/BDY/bdytra.F90

    r11071 r11191  
    4949      INTEGER, INTENT(in) ::   kt   ! Main time step counter 
    5050      ! 
    51       INTEGER                        :: ib_bdy, jn, igrd   ! Loop indeces 
    52       TYPE(ztrabdy), DIMENSION(jpts) :: zdta               ! Temporary data structure 
    53       LOGICAL, DIMENSION(4)          :: llsend1, llrecv1     ! indicate how communications are to be carried out 
     51      INTEGER                        :: ib_bdy, jn, igrd, ir   ! Loop indeces 
     52      TYPE(ztrabdy), DIMENSION(jpts) :: zdta                   ! Temporary data structure 
     53      LOGICAL                        :: llrim0                 ! indicate if rim 0 is treated 
     54      LOGICAL, DIMENSION(4)          :: llsend1, llrecv1       ! indicate how communications are to be carried out 
    5455      !!---------------------------------------------------------------------- 
    5556      igrd = 1  
    5657 
    57       DO ib_bdy=1, nb_bdy 
     58      DO ir = 1, 0, -1   ! treat rim 1 before rim 0 
     59         IF( ir == 0 ) THEN   ;   llrim0 = .TRUE. 
     60         ELSE                 ;   llrim0 = .FALSE. 
     61         END IF 
     62         DO ib_bdy=1, nb_bdy 
     63            ! 
     64            zdta(1)%tra => dta_bdy(ib_bdy)%tem 
     65            zdta(2)%tra => dta_bdy(ib_bdy)%sal 
     66            ! 
     67            DO jn = 1, jpts 
     68               ! 
     69               SELECT CASE( TRIM(cn_tra(ib_bdy)) ) 
     70               CASE('none'        )   ;   CYCLE 
     71               CASE('frs'         )   ! treat the whole boundary at once 
     72                  IF( ir == 0 ) CALL bdy_frs ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), zdta(jn)%tra ) 
     73               CASE('specified'   )   ! treat the whole rim      at once 
     74                  IF( ir == 0 ) CALL bdy_spe ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), zdta(jn)%tra ) 
     75               CASE('neumann'     )   ;   CALL bdy_nmn ( idx_bdy(ib_bdy), igrd         , tsa(:,:,:,jn), llrim0 )   ! tsa masked 
     76               CASE('orlanski'    )   ;   CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), & 
     77                    & zdta(jn)%tra, llrim0, ll_npo=.false. ) 
     78               CASE('orlanski_npo')   ;   CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), & 
     79                    & zdta(jn)%tra, llrim0, ll_npo=.true.  ) 
     80               CASE('runoff'      )   ;   CALL bdy_rnf ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), jn, llrim0 ) 
     81               CASE DEFAULT           ;   CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' ) 
     82               END SELECT 
     83               !  
     84            END DO 
     85         END DO 
    5886         ! 
    59          zdta(1)%tra => dta_bdy(ib_bdy)%tem 
    60          zdta(2)%tra => dta_bdy(ib_bdy)%sal 
    61          ! 
    62          DO jn = 1, jpts 
    63             ! 
     87         llsend1(:) = .false. 
     88         llrecv1(:) = .false. 
     89         DO ib_bdy=1, nb_bdy 
    6490            SELECT CASE( TRIM(cn_tra(ib_bdy)) ) 
    65             CASE('none'        )   ;   CYCLE 
    66             CASE('frs'         )   ;   CALL bdy_frs ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), zdta(jn)%tra ) 
    67             CASE('specified'   )   ;   CALL bdy_spe ( idx_bdy(ib_bdy),                tsa(:,:,:,jn), zdta(jn)%tra ) 
    68             CASE('neumann'     )   ;   CALL bdy_nmn ( idx_bdy(ib_bdy), igrd         , tsa(:,:,:,jn)               )   ! tsa masked 
    69             CASE('orlanski'    )   ;   CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), zdta(jn)%tra, ll_npo=.false. ) 
    70             CASE('orlanski_npo')   ;   CALL bdy_orl ( idx_bdy(ib_bdy), tsb(:,:,:,jn), tsa(:,:,:,jn), zdta(jn)%tra, ll_npo=.true. ) 
    71             CASE('runoff'      )   ;   CALL bdy_rnf ( idx_bdy(ib_bdy),                tsa(:,:,:,jn),               jn ) 
    72             CASE DEFAULT           ;   CALL ctl_stop( 'bdy_tra : unrecognised option for open boundaries for T and S' ) 
     91            CASE('neumann','runoff') 
     92               llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:,ir)   ! possibly every direction, T points 
     93               llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:,ir)   ! possibly every direction, T points 
     94            CASE('orlanski', 'orlanski_npo') 
     95               llsend1(:) = llsend1(:) .OR. lsend_bdy(ib_bdy,1,:,ir)   ! possibly every direction, T points 
     96               llrecv1(:) = llrecv1(:) .OR. lrecv_bdy(ib_bdy,1,:,ir)   ! possibly every direction, T points 
    7397            END SELECT 
    74             !  
    7598         END DO 
     99         IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction 
     100            CALL lbc_bdy_lnk( 'bdytra', llsend1, llrecv1, tsa, 'T',  1. ) 
     101         END IF 
    76102      END DO 
    77       ! 
    78       llsend1(:) = .false. 
    79       llrecv1(:) = .false. 
    80       DO ib_bdy=1, nb_bdy 
    81          SELECT CASE( TRIM(cn_tra(ib_bdy)) ) 
    82          CASE('neumann','runoff') 
    83             llsend1(:) = llsend1(:) .OR. lsend_bdyint(ib_bdy,1,:)   ! possibly every direction, T points 
    84             llrecv1(:) = llrecv1(:) .OR. lrecv_bdyint(ib_bdy,1,:)   ! possibly every direction, T points 
    85          CASE('orlanski', 'orlanski_npo') 
    86             llsend1(:) = llsend1(:) .OR. lsend_bdy(ib_bdy,1,:)   ! possibly every direction, T points 
    87             llrecv1(:) = llrecv1(:) .OR. lrecv_bdy(ib_bdy,1,:)   ! possibly every direction, T points 
    88          END SELECT 
    89       END DO 
    90       IF( ANY(llsend1) .OR. ANY(llrecv1) ) THEN   ! if need to send/recv in at least one direction 
    91          CALL lbc_bdy_lnk( 'bdytra', llsend1, llrecv1, tsa, 'T',  1. ) 
    92       END IF 
    93103      ! 
    94104   END SUBROUTINE bdy_tra 
    95105 
    96106 
    97    SUBROUTINE bdy_rnf( idx, pta, jpa ) 
     107   SUBROUTINE bdy_rnf( idx, pta, jpa, llrim0 ) 
    98108      !!---------------------------------------------------------------------- 
    99109      !!                 ***  SUBROUTINE bdy_rnf  *** 
     
    104114      !!  
    105115      !!---------------------------------------------------------------------- 
    106       TYPE(OBC_INDEX),                     INTENT(in) ::   idx  ! OBC indices 
    107       REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pta  ! tracer trend 
    108       INTEGER,                             INTENT(in) ::   jpa  ! TRA index 
     116      TYPE(OBC_INDEX),                     INTENT(in) ::   idx      ! OBC indices 
     117      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   pta      ! tracer trend 
     118      INTEGER,                             INTENT(in) ::   jpa      ! TRA index 
     119      LOGICAL,                             INTENT(in) ::   llrim0   ! indicate if rim 0 is treated 
    109120      ! 
    110121      INTEGER  ::   ib, ii, ij, igrd   ! dummy loop indices 
     
    113124      ! 
    114125      igrd = 1                       ! Everything is at T-points here 
    115       IF(      jpa == jp_tem )   THEN 
    116          CALL bdy_nmn( idx, igrd, pta ) 
    117       ELSE IF( jpa == jp_sal )   THEN 
    118          DO ib = 1, idx%nblenrim(igrd) 
     126      IF(      jpa == jp_tem ) THEN 
     127         CALL bdy_nmn( idx, igrd, pta, llrim0 ) 
     128      ELSE IF( jpa == jp_sal ) THEN 
     129         IF( .NOT. llrim0 )   RETURN 
     130         DO ib = 1, idx%nblenrim(igrd)   ! if llrim0 then treat the whole rim 
    119131            ii = idx%nbi(ib,igrd) 
    120132            ij = idx%nbj(ib,igrd) 
  • NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/OCE/DYN/dynkeg.F90

    r11071 r11191  
    140140            llrecv1(:) = .false. 
    141141            DO ib_bdy = 1, nb_bdy 
    142                llsend1(2) = llsend1(2) .OR. lsend_bdy(ib_bdy,igrd,2)   ! send east 
    143                llsend1(4) = llsend1(4) .OR. lsend_bdy(ib_bdy,igrd,4)   ! send north 
    144                llrecv1(1) = llrecv1(1) .OR. lrecv_bdy(ib_bdy,igrd,1)   ! receive west  
    145                llrecv1(3) = llrecv1(3) .OR. lrecv_bdy(ib_bdy,igrd,3)   ! receive south 
     142               llsend1(2) = llsend1(2) .OR. ANY(lsend_bdy(ib_bdy,igrd,2,:))   ! send east 
     143               llsend1(4) = llsend1(4) .OR. ANY(lsend_bdy(ib_bdy,igrd,4,:))   ! send north 
     144               llrecv1(1) = llrecv1(1) .OR. ANY(lrecv_bdy(ib_bdy,igrd,1,:))   ! receive west  
     145               llrecv1(3) = llrecv1(3) .OR. ANY(lrecv_bdy(ib_bdy,igrd,3,:))   ! receive south 
    146146            END DO 
    147147    
Note: See TracChangeset for help on using the changeset viewer.