Index: NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traldf_lap_blp.F90
===================================================================
--- NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traldf_lap_blp.F90	(revision 13899)
+++ NEMO/branches/2020/tickets_icb_1900/src/OCE/TRA/traldf_lap_blp.F90	(revision 14012)
@@ -13,4 +13,5 @@
    USE oce            ! ocean dynamics and active tracers
    USE dom_oce        ! ocean space and time domain
+   USE domutl, ONLY : is_tile
    USE ldftra         ! lateral physics: eddy diffusivity
    USE traldf_iso     ! iso-neutral lateral diffusion (standard operator)     (tra_ldf_iso   routine)
@@ -46,7 +47,29 @@
 CONTAINS
 
-   SUBROUTINE tra_ldf_lap( kt, Kmm, kit000, cdtype, pahu, pahv  ,               &
-      &                                             pgu , pgv   , pgui, pgvi,   &
-      &                                             pt  , pt_rhs, kjpt, kpass ) 
+   SUBROUTINE tra_ldf_lap( kt, Kmm, kit000, cdtype, pahu, pahv,             &
+      &                                             pgu , pgv , pgui, pgvi, &
+      &                                             pt, pt_rhs, kjpt, kpass )
+      !!
+      INTEGER                     , INTENT(in   ) ::   kt         ! ocean time-step index
+      INTEGER                     , INTENT(in   ) ::   kit000     ! first time step index
+      CHARACTER(len=3)            , INTENT(in   ) ::   cdtype     ! =TRA or TRC (tracer indicator)
+      INTEGER                     , INTENT(in   ) ::   kjpt       ! number of tracers
+      INTEGER                     , INTENT(in   ) ::   kpass      ! =1/2 first or second passage
+      INTEGER                     , INTENT(in   ) ::   Kmm        ! ocean time level index
+      REAL(wp), DIMENSION(:,:,:)  , INTENT(in   ) ::   pahu, pahv ! eddy diffusivity at u- and v-points  [m2/s]
+      REAL(wp), DIMENSION(:,:,:)  , INTENT(in   ) ::   pgu, pgv   ! tracer gradient at pstep levels
+      REAL(wp), DIMENSION(:,:,:)  , INTENT(in   ) ::   pgui, pgvi ! tracer gradient at top   levels
+      REAL(wp), DIMENSION(:,:,:,:), INTENT(in   ) ::   pt         ! before tracer fields
+      REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) ::   pt_rhs     ! tracer trend
+      !!
+      CALL tra_ldf_lap_t( kt, Kmm, kit000, cdtype, pahu, pahv, is_tile(pahu),                            &
+      &                                            pgu , pgv , is_tile(pgu) , pgui, pgvi, is_tile(pgui), &
+      &                                            pt, is_tile(pt), pt_rhs, is_tile(pt_rhs), kjpt, kpass )
+   END SUBROUTINE tra_ldf_lap
+
+
+   SUBROUTINE tra_ldf_lap_t( kt, Kmm, kit000, cdtype, pahu, pahv, ktah,                   &
+      &                                               pgu , pgv , ktg , pgui, pgvi, ktgi, &
+      &                                               pt, ktt, pt_rhs, ktt_rhs, kjpt, kpass )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE tra_ldf_lap  ***
@@ -72,26 +95,30 @@
       INTEGER                              , INTENT(in   ) ::   kpass      ! =1/2 first or second passage
       INTEGER                              , INTENT(in   ) ::   Kmm        ! ocean time level index
-      REAL(wp), DIMENSION(jpi,jpj,jpk)     , INTENT(in   ) ::   pahu, pahv ! eddy diffusivity at u- and v-points  [m2/s]
-      REAL(wp), DIMENSION(jpi,jpj    ,kjpt), INTENT(in   ) ::   pgu, pgv   ! tracer gradient at pstep levels
-      REAL(wp), DIMENSION(jpi,jpj,    kjpt), INTENT(in   ) ::   pgui, pgvi ! tracer gradient at top   levels
-      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   pt         ! before tracer fields
-      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pt_rhs     ! tracer trend 
-      !
-      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
-      REAL(wp) ::   zsign            ! local scalars
-      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   ztu, ztv, zaheeu, zaheev
-      !!----------------------------------------------------------------------
-      !
-      IF( kt == nit000 .AND. lwp )  THEN
-         WRITE(numout,*)
-         WRITE(numout,*) 'tra_ldf_lap : iso-level laplacian diffusion on ', cdtype, ', pass=', kpass
-         WRITE(numout,*) '~~~~~~~~~~~ '
-      ENDIF
-      !
-      l_hst = .FALSE.
-      l_ptr = .FALSE.
-      IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) )     l_ptr = .TRUE.
-      IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. &
-         &                        iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )  l_hst = .TRUE.
+      INTEGER                              , INTENT(in   ) ::   ktah, ktg, ktgi, ktt, ktt_rhs
+      REAL(wp), DIMENSION(A2D_T(ktah),   JPK)     , INTENT(in   ) ::   pahu, pahv ! eddy diffusivity at u- and v-points  [m2/s]
+      REAL(wp), DIMENSION(A2D_T(ktg),        KJPT), INTENT(in   ) ::   pgu, pgv   ! tracer gradient at pstep levels
+      REAL(wp), DIMENSION(A2D_T(ktgi),       KJPT), INTENT(in   ) ::   pgui, pgvi ! tracer gradient at top   levels
+      REAL(wp), DIMENSION(A2D_T(ktt),    JPK,KJPT), INTENT(in   ) ::   pt         ! before tracer fields
+      REAL(wp), DIMENSION(A2D_T(ktt_rhs),JPK,KJPT), INTENT(inout) ::   pt_rhs     ! tracer trend
+      !
+      INTEGER  ::   ji, jj, jk, jn      ! dummy loop indices
+      INTEGER  ::   isi, iei, isj, iej  ! local integers
+      REAL(wp) ::   zsign               ! local scalars
+      REAL(wp), DIMENSION(A2D(nn_hls),jpk) ::   ztu, ztv, zaheeu, zaheev
+      !!----------------------------------------------------------------------
+      !
+      IF( ntile == 0 .OR. ntile == 1 )  THEN                       ! Do only on the first tile
+         IF( kt == nit000 .AND. lwp )  THEN
+            WRITE(numout,*)
+            WRITE(numout,*) 'tra_ldf_lap : iso-level laplacian diffusion on ', cdtype, ', pass=', kpass
+            WRITE(numout,*) '~~~~~~~~~~~ '
+         ENDIF
+         !
+         l_hst = .FALSE.
+         l_ptr = .FALSE.
+         IF( cdtype == 'TRA' .AND. ( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) )     l_ptr = .TRUE.
+         IF( cdtype == 'TRA' .AND. ( iom_use("uadv_heattr") .OR. iom_use("vadv_heattr") .OR. &
+            &                        iom_use("uadv_salttr") .OR. iom_use("vadv_salttr")  ) )  l_hst = .TRUE.
+      ENDIF
       !
       !                                !==  Initialization of metric arrays used for all tracers  ==!
@@ -99,5 +126,11 @@
       ELSE                    ;   zsign = -1._wp
       ENDIF
-      DO_3D( 1, 0, 1, 0, 1, jpkm1 )
+
+      IF( ntsi == Nis0 ) THEN ; isi = nn_hls - 1 ; ELSE ; isi = 0 ; ENDIF    ! Avoid double-counting when using tiling
+      IF( ntsj == Njs0 ) THEN ; isj = nn_hls - 1 ; ELSE ; isj = 0 ; ENDIF
+      IF( ntei == Nie0 ) THEN ; iei = nn_hls - 1 ; ELSE ; iei = 0 ; ENDIF
+      IF( ntej == Nje0 ) THEN ; iej = nn_hls - 1 ; ELSE ; iej = 0 ; ENDIF
+
+      DO_3D( nn_hls, nn_hls-1, nn_hls, nn_hls-1, 1, jpkm1 )            !== First derivative (gradient)  ==!
          zaheeu(ji,jj,jk) = zsign * pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm)   !!gm   * umask(ji,jj,jk) pah masked!
          zaheev(ji,jj,jk) = zsign * pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm)   !!gm   * vmask(ji,jj,jk)
@@ -108,15 +141,15 @@
          !                          ! =========== !    
          !                               
-         DO_3D( 1, 0, 1, 0, 1, jpkm1 )            !== First derivative (gradient)  ==!
+         DO_3D( nn_hls, nn_hls-1, nn_hls, nn_hls-1, 1, jpkm1 )            !== First derivative (gradient)  ==!
             ztu(ji,jj,jk) = zaheeu(ji,jj,jk) * ( pt(ji+1,jj  ,jk,jn) - pt(ji,jj,jk,jn) )
             ztv(ji,jj,jk) = zaheev(ji,jj,jk) * ( pt(ji  ,jj+1,jk,jn) - pt(ji,jj,jk,jn) )
          END_3D
          IF( ln_zps ) THEN                             ! set gradient at bottom/top ocean level
-            DO_2D( 1, 0, 1, 0 )                              ! bottom
+            DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )                              ! bottom
                ztu(ji,jj,mbku(ji,jj)) = zaheeu(ji,jj,mbku(ji,jj)) * pgu(ji,jj,jn)
                ztv(ji,jj,mbkv(ji,jj)) = zaheev(ji,jj,mbkv(ji,jj)) * pgv(ji,jj,jn)
             END_2D
             IF( ln_isfcav ) THEN                             ! top in ocean cavities only
-               DO_2D( 1, 0, 1, 0 )
+               DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )
                   IF( miku(ji,jj) > 1 )   ztu(ji,jj,miku(ji,jj)) = zaheeu(ji,jj,miku(ji,jj)) * pgui(ji,jj,jn) 
                   IF( mikv(ji,jj) > 1 )   ztv(ji,jj,mikv(ji,jj)) = zaheev(ji,jj,mikv(ji,jj)) * pgvi(ji,jj,jn) 
@@ -125,5 +158,6 @@
          ENDIF
          !
-         DO_3D( 0, 0, 0, 0, 1, jpkm1 )            !== Second derivative (divergence) added to the general tracer trends  ==!
+         ! NOTE: [tiling-comms-merge] Bounds changed to avoid repeating this calculation for overlapping rows when using tiling
+         DO_3D( isj, iej, isi, iei, 1, jpkm1 )            !== Second derivative (divergence) added to the general tracer trends  ==!
             pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + (  ztu(ji,jj,jk) - ztu(ji-1,jj,jk)     &
                &                                      +    ztv(ji,jj,jk) - ztv(ji,jj-1,jk) )   &
@@ -142,5 +176,5 @@
       !                             ! ==================
       !
-   END SUBROUTINE tra_ldf_lap
+   END SUBROUTINE tra_ldf_lap_t
    
 
@@ -173,17 +207,19 @@
       !
       INTEGER ::   ji, jj, jk, jn   ! dummy loop indices
-      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt) :: zlap         ! laplacian at t-point
-      REAL(wp), DIMENSION(jpi,jpj,    kjpt) :: zglu, zglv   ! bottom GRADh of the laplacian (u- and v-points)
-      REAL(wp), DIMENSION(jpi,jpj,    kjpt) :: zgui, zgvi   ! top    GRADh of the laplacian (u- and v-points)
+      REAL(wp), DIMENSION(A2D(nn_hls),jpk,kjpt) :: zlap         ! laplacian at t-point
+      REAL(wp), DIMENSION(A2D(nn_hls),    kjpt) :: zglu, zglv   ! bottom GRADh of the laplacian (u- and v-points)
+      REAL(wp), DIMENSION(A2D(nn_hls),    kjpt) :: zgui, zgvi   ! top    GRADh of the laplacian (u- and v-points)
       !!---------------------------------------------------------------------
       !
-      IF( kt == kit000 .AND. lwp )  THEN
-         WRITE(numout,*)
-         SELECT CASE ( kldf )
-         CASE ( np_blp    )   ;   WRITE(numout,*) 'tra_ldf_blp : iso-level   bilaplacian operator on ', cdtype
-         CASE ( np_blp_i  )   ;   WRITE(numout,*) 'tra_ldf_blp : iso-neutral bilaplacian operator on ', cdtype, ' (Standard)'
-         CASE ( np_blp_it )   ;   WRITE(numout,*) 'tra_ldf_blp : iso-neutral bilaplacian operator on ', cdtype, ' (triad)'
-         END SELECT
-         WRITE(numout,*) '~~~~~~~~~~~'
+      IF( ntile == 0 .OR. ntile == 1 )  THEN                       ! Do only on the first tile
+         IF( kt == kit000 .AND. lwp )  THEN
+            WRITE(numout,*)
+            SELECT CASE ( kldf )
+            CASE ( np_blp    )   ;   WRITE(numout,*) 'tra_ldf_blp : iso-level   bilaplacian operator on ', cdtype
+            CASE ( np_blp_i  )   ;   WRITE(numout,*) 'tra_ldf_blp : iso-neutral bilaplacian operator on ', cdtype, ' (Standard)'
+            CASE ( np_blp_it )   ;   WRITE(numout,*) 'tra_ldf_blp : iso-neutral bilaplacian operator on ', cdtype, ' (triad)'
+            END SELECT
+            WRITE(numout,*) '~~~~~~~~~~~'
+         ENDIF
       ENDIF
 
@@ -200,4 +236,5 @@
       END SELECT
       !
+      ! NOTE: [tiling-comms-merge] Needed for both nn_hls as tra_ldf_iso and tra_ldf_triad have not yet been adjusted to work with nn_hls = 2. In the zps case the lbc_lnk in zps_hde handles this, but in the zco case zlap always needs this lbc_lnk. I did try adjusting the bounds in tra_ldf_iso and tra_ldf_triad so this lbc_lnk was only needed for nn_hls = 1, but this was not correct and I did not have time to figure out why
       CALL lbc_lnk( 'traldf_lap_blp', zlap(:,:,:,:) , 'T', 1.0_wp )     ! Lateral boundary conditions (unchanged sign)
       !                                               ! Partial top/bottom cell: GRADh( zlap )  
