Index: NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/cfgs/SHARED/namelist_ref
===================================================================
--- NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/cfgs/SHARED/namelist_ref	(revision 12765)
+++ NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/cfgs/SHARED/namelist_ref	(revision 12766)
@@ -94,4 +94,11 @@
    ln_use_jattr = .false.    !  use (T) the file attribute: open_ocean_jstart, if present
    !                         !  in netcdf input files, as the start j-row for reading
+/
+!-----------------------------------------------------------------------
+&namtile        !   parameters of the tiling
+!-----------------------------------------------------------------------
+   ln_tile = .false.     !  Use tiling (T) or not (F)
+   nn_tile_i = 10        !  Length of tiles in i
+   nn_tile_j = 10        !  Length of tiles in j
 /
 !-----------------------------------------------------------------------
Index: NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/DOM/dom_oce.F90
===================================================================
--- NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/DOM/dom_oce.F90	(revision 12765)
+++ NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/DOM/dom_oce.F90	(revision 12766)
@@ -72,4 +72,13 @@
    !                                !  = 7 bi-cyclic East-West AND North-South
    LOGICAL, PUBLIC ::   l_Iperio, l_Jperio   !   should we explicitely take care I/J periodicity 
+
+   ! Tiling namelist
+   LOGICAL, PUBLIC ::   ln_tile
+   INTEGER         ::   nn_tile_i, nn_tile_j
+
+   ! Tile indices
+   TYPE, PUBLIC    ::   TILE
+      INTEGER      ::   ntsi, ntsj, ntei, ntej, ntsim1, ntsjm1, nteip1, ntejp1, ntile
+   END TYPE TILE
 
    !                                 !  domain MPP decomposition parameters
Index: NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/DOM/domain.F90
===================================================================
--- NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/DOM/domain.F90	(revision 12765)
+++ NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/DOM/domain.F90	(revision 12766)
@@ -49,4 +49,5 @@
 
    PUBLIC   dom_init     ! called by nemogcm.F90
+   PUBLIC   dom_tile     ! called by step.F90
    PUBLIC   domain_cfg   ! called by nemogcm.F90
 
@@ -270,4 +271,48 @@
 
 
+   SUBROUTINE dom_tile(kntile, ktile)
+      !!----------------------------------------------------------------------
+      !!                     ***  ROUTINE dom_tile  ***
+      !!
+      !! ** Purpose :   Set domain indices for specified tile
+      !!
+      !! ** Action  : - ntile          : current tile number
+      !!              - ntsi, ntsj     : start of internal part of domain
+      !!              - ntei, ntej     : end of internal part of domain
+      !!              - ntsim1, ntsjm1 : start of domain
+      !!              - nteip1, ntejp1 : end of domain
+      !!----------------------------------------------------------------------
+      INTEGER   , INTENT(in ) :: kntile               ! Tile number
+      TYPE(TILE), INTENT(out) :: ktile                ! Tile indices
+      INTEGER                 :: iitile, ijtile       ! Tile number in i and j
+      !!----------------------------------------------------------------------
+
+      IF( ln_tile .AND. kntile > 0 ) THEN          ! Tile domain
+         iitile = 1 + MOD( kntile - 1, jpnitile )
+         ijtile = 1 + (kntile - 1) / jpnitile
+
+         ktile % ntile = kntile
+         ktile % ntsi = 2 + (iitile - 1) * nn_tile_i
+         ktile % ntsj = 2 + (ijtile - 1) * nn_tile_j
+         ktile % ntei = MIN(ktile % ntsi + nn_tile_i - 1, jpim1)   ! Size of last tile limited by full domain
+         ktile % ntej = MIN(ktile % ntsj + nn_tile_j - 1, jpjm1)   !
+         ktile % ntsim1 = ktile % ntsi - 1
+         ktile % ntsjm1 = ktile % ntsj - 1
+         ktile % nteip1 = ktile % ntei + 1
+         ktile % ntejp1 = ktile % ntej + 1
+      ELSE                                         ! Full domain
+         ktile % ntile = 1
+         ktile % ntsi = 2
+         ktile % ntsj = 2
+         ktile % ntei = jpim1
+         ktile % ntej = jpjm1
+         ktile % ntsim1 = 1
+         ktile % ntsjm1 = 1
+         ktile % nteip1 = jpi
+         ktile % ntejp1 = jpj
+      ENDIF
+   END SUBROUTINE dom_tile
+
+
    SUBROUTINE dom_nam
       !!----------------------------------------------------------------------
@@ -278,4 +323,5 @@
       !! ** input   : - namrun namelist
       !!              - namdom namelist
+      !!              - namtile namelist
       !!              - namnc4 namelist   ! "key_netcdf4" only
       !!----------------------------------------------------------------------
@@ -290,4 +336,5 @@
          &             ln_cfmeta, ln_xios_read, nn_wxios
       NAMELIST/namdom/ ln_linssh, rn_Dt, rn_atfp, ln_crs, ln_meshmask
+      NAMELIST/namtile/ ln_tile, nn_tile_i, nn_tile_j
 #if defined key_netcdf4
       NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip
@@ -417,4 +464,36 @@
       r1_Dt = 1._wp / rDt
 
+      READ  ( numnam_ref, namtile, IOSTAT = ios, ERR = 905 )
+905   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namtile in reference namelist' )
+      READ  ( numnam_cfg, namtile, IOSTAT = ios, ERR = 906 )
+906   IF( ios >  0 )   CALL ctl_nam ( ios , 'namtile in configuration namelist' )
+      IF(lwm) WRITE( numond, namtile )
+
+      ! Set tile decomposition
+      IF( ln_tile ) THEN
+         jpnitile = (jpi - 2) / nn_tile_i
+         jpnjtile = (jpj - 2) / nn_tile_j
+         IF( MOD( jpi - 2, nn_tile_i ) /= 0 ) jpnitile = jpnitile + 1
+         IF( MOD( jpj - 2, nn_tile_j ) /= 0 ) jpnjtile = jpnjtile + 1
+      ELSE
+         jpnitile = 1
+         jpnjtile = 1
+      ENDIF
+      jpnijtile = jpnitile * jpnjtile
+
+      IF(lwp) THEN
+         WRITE(numout,*)
+         WRITE(numout,*) '   Namelist : namtile   ---   tiling decomposition'
+         WRITE(numout,*) '      Tiling (T) or not (F)                ln_tile   = ', ln_tile
+         WRITE(numout,*) '      Length of tile in i                  nn_tile_i = ', nn_tile_i
+         WRITE(numout,*) '      Length of tile in j                  nn_tile_j = ', nn_tile_j
+         WRITE(numout,*)
+         IF( ln_tile ) THEN
+            WRITE(numout,*) '      The domain will be decomposed into', jpnijtile, 'tiles of size', nn_tile_i, 'x', nn_tile_j
+         ELSE
+            WRITE(numout,*) '      Domain tiling will NOT be used'
+         ENDIF
+      ENDIF
+
       IF( TRIM(Agrif_CFixed()) == '0' ) THEN
          lrxios = ln_xios_read.AND.ln_rstart
Index: NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/TRA/traldf.F90
===================================================================
--- NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/TRA/traldf.F90	(revision 12765)
+++ NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/TRA/traldf.F90	(revision 12766)
@@ -45,5 +45,5 @@
 CONTAINS
 
-   SUBROUTINE tra_ldf( kt, Kbb, Kmm, pts, Krhs )
+   SUBROUTINE tra_ldf(ktile, kt, Kbb, Kmm, pts, Krhs )
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE tra_ldf  ***
@@ -51,4 +51,5 @@
       !! ** Purpose :   compute the lateral ocean tracer physics.
       !!----------------------------------------------------------------------
+      TYPE(TILE),                                INTENT(in   ) :: ktile           ! Tile indices
       INTEGER,                                   INTENT(in   ) :: kt              ! ocean time-step index
       INTEGER,                                   INTENT(in   ) :: Kbb, Kmm, Krhs  ! ocean time level indices
@@ -58,10 +59,16 @@
       !!----------------------------------------------------------------------
       !
-      IF( ln_timing )   CALL timing_start('tra_ldf')
+      IF( ktile % ntile == 1 )  THEN              ! Do only on the first tile
+         ! TODO: TO BE TILED
+         IF( ln_timing )   CALL timing_start('tra_ldf')
+      ENDIF
       !
-      IF( l_trdtra )   THEN                    !* Save ta and sa trends
-         ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) ) 
-         ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) 
-         ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs)
+      IF( ktile % ntile == jpnijtile )  THEN      ! Do only after all tiles finish
+         IF( l_trdtra )   THEN                    !* Save ta and sa trends
+            ! TODO: TO BE TILED
+            ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) )
+            ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs)
+            ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs)
+         ENDIF
       ENDIF
       !
@@ -70,5 +77,5 @@
          CALL tra_ldf_lap  ( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs),                   jpts,  1 )
       CASE ( np_lap_i )                                  ! laplacian: standard iso-neutral operator (Madec)
-         CALL tra_ldf_iso  ( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts,  1 )
+         CALL tra_ldf_iso  ( ktile, kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts,  1 )
       CASE ( np_lap_it )                                 ! laplacian: triad iso-neutral operator (griffies)
          CALL tra_ldf_triad( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts,  1 )
@@ -77,16 +84,23 @@
       END SELECT
       !
-      IF( l_trdtra )   THEN                    !* save the horizontal diffusive trends for further diagnostics
-         ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) - ztrdt(:,:,:)
-         ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs) - ztrds(:,:,:)
-         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_ldf, ztrdt )
-         CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_ldf, ztrds )
-         DEALLOCATE( ztrdt, ztrds ) 
+      IF( ktile % ntile == jpnijtile )  THEN      ! Do only after all tiles finish
+         IF( l_trdtra )   THEN                    !* save the horizontal diffusive trends for further diagnostics
+            ! TODO: TO BE TILED
+            ztrdt(:,:,:) = pts(:,:,:,jp_tem,Krhs) - ztrdt(:,:,:)
+            ztrds(:,:,:) = pts(:,:,:,jp_sal,Krhs) - ztrds(:,:,:)
+            ! TODO: TO BE TILED
+            CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_tem, jptra_ldf, ztrdt )
+            CALL trd_tra( kt, Kmm, Krhs, 'TRA', jp_sal, jptra_ldf, ztrds )
+            DEALLOCATE( ztrdt, ztrds )
+         ENDIF
+
+         !                                        !* print mean trends (used for debugging)
+         ! TODO: TO BE TILED
+         IF(sn_cfctl%l_prtctl)   CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' ldf  - Ta: ', mask1=tmask,               &
+            &                                  tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
+         !
+         ! TODO: TO BE TILED
+         IF( ln_timing )   CALL timing_stop('tra_ldf')
       ENDIF
-      !                                        !* print mean trends (used for debugging)
-      IF(sn_cfctl%l_prtctl)   CALL prt_ctl( tab3d_1=pts(:,:,:,jp_tem,Krhs), clinfo1=' ldf  - Ta: ', mask1=tmask,               &
-         &                                  tab3d_2=pts(:,:,:,jp_sal,Krhs), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
-      !
-      IF( ln_timing )   CALL timing_stop('tra_ldf')
       !
    END SUBROUTINE tra_ldf
Index: NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/TRA/traldf_iso.F90
===================================================================
--- NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/TRA/traldf_iso.F90	(revision 12765)
+++ NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/TRA/traldf_iso.F90	(revision 12766)
@@ -36,7 +36,4 @@
    PUBLIC   tra_ldf_iso   ! routine called by step.F90
 
-   LOGICAL  ::   l_ptr   ! flag to compute poleward transport
-   LOGICAL  ::   l_hst   ! flag to compute heat transport
-
    !! * Substitutions
 #  include "do_loop_substitute.h90"
@@ -48,5 +45,5 @@
 CONTAINS
 
-  SUBROUTINE tra_ldf_iso( kt, Kmm, kit000, cdtype, pahu, pahv,                    &
+  SUBROUTINE tra_ldf_iso( ktile, kt, Kmm, kit000, cdtype, pahu, pahv,                    &
       &                                            pgu , pgv    ,   pgui, pgvi,   &
       &                                       pt , pt2 , pt_rhs , kjpt  , kpass )
@@ -91,4 +88,5 @@
       !! ** Action :   Update pt_rhs arrays with the before rotated diffusion
       !!----------------------------------------------------------------------
+      TYPE(TILE)                           , INTENT(in   ) ::   ktile      ! Tile indices
       INTEGER                              , INTENT(in   ) ::   kt         ! ocean time-step index
       INTEGER                              , INTENT(in   ) ::   kit000     ! first time step index
@@ -104,4 +102,6 @@
       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pt_rhs     ! tracer trend
       !
+      LOGICAL  ::  l_ptr                                 ! flag to compute poleward transport
+      LOGICAL  ::  l_hst                                 ! flag to compute heat transport
       INTEGER  ::  ji, jj, jk, jn   ! dummy loop indices
       INTEGER  ::  ikt
@@ -110,15 +110,20 @@
       REAL(wp) ::  zmskv, zahv_w, zabe2, zcof2, zcoef4   !   -      -
       REAL(wp) ::  zcoef0, ze3w_2, zsign                 !   -      -
-      REAL(wp), DIMENSION(jpi,jpj)     ::   zdkt, zdk1t, z2d
-      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zdit, zdjt, zftu, zftv, ztfw 
+      REAL(wp), DIMENSION(IND_2D)     ::   zdkt, zdk1t, z2d
+      REAL(wp), DIMENSION(IND_2D,jpk) ::   zdit, zdjt, zftu, zftv, ztfw
       !!----------------------------------------------------------------------
       !
       IF( kpass == 1 .AND. kt == kit000 )  THEN
-         IF(lwp) WRITE(numout,*)
-         IF(lwp) WRITE(numout,*) 'tra_ldf_iso : rotated laplacian diffusion operator on ', cdtype
-         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
-         !
-         akz     (:,:,:) = 0._wp      
-         ah_wslp2(:,:,:) = 0._wp
+         IF( ktile % ntile == 1 )  THEN               ! Do only on the first tile
+            ! TODO: TO BE TILED
+            IF(lwp) WRITE(numout,*)
+            IF(lwp) WRITE(numout,*) 'tra_ldf_iso : rotated laplacian diffusion operator on ', cdtype
+            IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
+         ENDIF
+         !
+         DO_3D_11_11_T( 1, jpk )
+            akz     (ji,jj,jk) = 0._wp
+            ah_wslp2(ji,jj,jk) = 0._wp
+         END_3D
       ENDIF
       !   
@@ -140,5 +145,5 @@
       IF( kpass == 1 ) THEN                  !==  first pass only  ==!
          !
-         DO_3D_00_00( 2, jpkm1 )
+         DO_3D_00_00_T( 2, jpkm1 )
             !
             zmsku = wmask(ji,jj,jk) / MAX(   umask(ji  ,jj,jk-1) + umask(ji-1,jj,jk)          &
@@ -157,5 +162,5 @@
          !
          IF( ln_traldf_msc ) THEN                ! stabilizing vertical diffusivity coefficient
-            DO_3D_00_00( 2, jpkm1 )
+            DO_3D_00_00_T( 2, jpkm1 )
                akz(ji,jj,jk) = 0.25_wp * (                                                                     &
                   &              ( pahu(ji  ,jj,jk) + pahu(ji  ,jj,jk-1) ) / ( e1u(ji  ,jj) * e1u(ji  ,jj) )   &
@@ -166,10 +171,10 @@
             !
             IF( ln_traldf_blp ) THEN                ! bilaplacian operator
-               DO_3D_10_10( 2, jpkm1 )
+               DO_3D_10_10_T( 2, jpkm1 )
                   akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk)   &
                      &          * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) )  )
                END_3D
             ELSEIF( ln_traldf_lap ) THEN              ! laplacian operator
-               DO_3D_10_10( 2, jpkm1 )
+               DO_3D_10_10_T( 2, jpkm1 )
                   ze3w_2 = e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm)
                   zcoef0 = rDt * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2  )
@@ -179,5 +184,7 @@
            !
          ELSE                                    ! 33 flux set to zero with akz=ah_wslp2 ==>> computed in full implicit
-            akz(:,:,:) = ah_wslp2(:,:,:)      
+            DO_3D_11_11_T( 1, jpk )
+               akz(ji,jj,jk) = ah_wslp2(ji,jj,jk)
+            END_3D
          ENDIF
       ENDIF
@@ -196,15 +203,15 @@
 
          ! Horizontal tracer gradient 
-         DO_3D_10_10( 1, jpkm1 )
+         DO_3D_10_10_T( 1, jpkm1 )
             zdit(ji,jj,jk) = ( pt(ji+1,jj  ,jk,jn) - pt(ji,jj,jk,jn) ) * umask(ji,jj,jk)
             zdjt(ji,jj,jk) = ( pt(ji  ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk)
          END_3D
          IF( ln_zps ) THEN      ! botton and surface ocean correction of the horizontal gradient
-            DO_2D_10_10
+            DO_2D_10_10_T
                zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn)          
                zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn)
             END_2D
             IF( ln_isfcav ) THEN      ! first wet level beneath a cavity
-               DO_2D_10_10
+               DO_2D_10_10_T
                   IF( miku(ji,jj) > 1 )   zdit(ji,jj,miku(ji,jj)) = pgui(ji,jj,jn)          
                   IF( mikv(ji,jj) > 1 )   zdjt(ji,jj,mikv(ji,jj)) = pgvi(ji,jj,jn)     
@@ -219,11 +226,14 @@
          DO jk = 1, jpkm1                                 ! Horizontal slab
             !
-            !                             !== Vertical tracer gradient
-            zdk1t(:,:) = ( pt(:,:,jk,jn) - pt(:,:,jk+1,jn) ) * wmask(:,:,jk+1)     ! level jk+1
-            !
-            IF( jk == 1 ) THEN   ;   zdkt(:,:) = zdk1t(:,:)                          ! surface: zdkt(jk=1)=zdkt(jk=2)
-            ELSE                 ;   zdkt(:,:) = ( pt(:,:,jk-1,jn) - pt(:,:,jk,jn) ) * wmask(:,:,jk)
-            ENDIF
-            DO_2D_10_10
+            DO_2D_11_11_T
+               !                             !== Vertical tracer gradient
+               zdk1t(ji,jj) = ( pt(ji,jj,jk,jn) - pt(ji,jj,jk+1,jn) ) * wmask(ji,jj,jk+1)     ! level jk+1
+               !
+               IF( jk == 1 ) THEN   ;   zdkt(ji,jj) = zdk1t(ji,jj)                          ! surface: zdkt(jk=1)=zdkt(jk=2)
+               ELSE                 ;   zdkt(ji,jj) = ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) * wmask(ji,jj,jk)
+               ENDIF
+            END_2D
+            !
+            DO_2D_10_10_T
                zabe1 = pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm)
                zabe2 = pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm)
@@ -246,5 +256,5 @@
             END_2D
             !
-            DO_2D_00_00
+            DO_2D_00_00_T
                pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  zftu(ji,jj,jk) - zftu(ji-1,jj,jk)      &
                   &                                                 + zftv(ji,jj,jk) - zftv(ji,jj-1,jk)  )   &
@@ -262,5 +272,5 @@
          ztfw(:,:, 1 ) = 0._wp      ;      ztfw(:,:,jpk) = 0._wp
          
-         DO_3D_00_00( 2, jpkm1 )
+         DO_3D_00_00_T( 2, jpkm1 )
             !
             zmsku = wmask(ji,jj,jk) / MAX(   umask(ji  ,jj,jk-1) + umask(ji-1,jj,jk)          &
@@ -284,5 +294,5 @@
          !                                !==  add the vertical 33 flux  ==!
          IF( ln_traldf_lap ) THEN               ! laplacian case: eddy coef = ah_wslp2 - akz
-            DO_3D_00_00( 2, jpkm1 )
+            DO_3D_00_00_T( 2, jpkm1 )
                ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * wmask(ji,jj,jk)   &
                   &                            * ( ah_wslp2(ji,jj,jk) - akz(ji,jj,jk) )               &
@@ -293,5 +303,5 @@
             SELECT CASE( kpass )
             CASE(  1  )                            ! 1st pass : eddy coef = ah_wslp2
-               DO_3D_00_00( 2, jpkm1 )
+               DO_3D_00_00_T( 2, jpkm1 )
                   ztfw(ji,jj,jk) = ztfw(ji,jj,jk)                       &
                      &           + ah_wslp2(ji,jj,jk)  * e1e2t(ji,jj)   &
@@ -299,5 +309,5 @@
                END_3D
             CASE(  2  )                         ! 2nd pass : eddy flux = ah_wslp2 and akz applied on pt  and pt2 gradients, resp.
-               DO_3D_00_00( 2, jpkm1 )
+               DO_3D_00_00_T( 2, jpkm1 )
                   ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * wmask(ji,jj,jk)                  &
                      &                            * (  ah_wslp2(ji,jj,jk) * ( pt (ji,jj,jk-1,jn) - pt (ji,jj,jk,jn) )   &
@@ -307,19 +317,23 @@
          ENDIF
          !         
-         DO_3D_00_00( 1, jpkm1 )
+         DO_3D_00_00_T( 1, jpkm1 )
             pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  ztfw (ji,jj,jk) - ztfw(ji,jj,jk+1)  )   &
                &                                              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm)
          END_3D
          !
-         IF( ( kpass == 1 .AND. ln_traldf_lap ) .OR.  &     !==  first pass only (  laplacian)  ==!
-             ( kpass == 2 .AND. ln_traldf_blp ) ) THEN      !==  2nd   pass      (bilaplacian)  ==!
-            !
-            !                             ! "Poleward" diffusive heat or salt transports (T-S case only)
-               ! note sign is reversed to give down-gradient diffusive transports )
-            IF( l_ptr )  CALL dia_ptr_hst( jn, 'ldf', -zftv(:,:,:)  )
-            !                          ! Diffusive heat transports
-            IF( l_hst )  CALL dia_ar5_hst( jn, 'ldf', -zftu(:,:,:), -zftv(:,:,:) )
-            !
-         ENDIF                                                    !== end pass selection  ==!
+         IF( ktile % ntile == jpnijtile )  THEN                ! Do only after all tiles finish
+            IF( ( kpass == 1 .AND. ln_traldf_lap ) .OR.  &     !==  first pass only (  laplacian)  ==!
+                ( kpass == 2 .AND. ln_traldf_blp ) ) THEN      !==  2nd   pass      (bilaplacian)  ==!
+               !
+               !                             ! "Poleward" diffusive heat or salt transports (T-S case only)
+                  ! note sign is reversed to give down-gradient diffusive transports )
+               ! TODO: TO BE TILED
+               IF( l_ptr )  CALL dia_ptr_hst( jn, 'ldf', -zftv(:,:,:)  )
+               !                          ! Diffusive heat transports
+               ! TODO: TO BE TILED
+               IF( l_hst )  CALL dia_ar5_hst( jn, 'ldf', -zftu(:,:,:), -zftv(:,:,:) )
+               !
+            ENDIF                                                    !== end pass selection  ==!
+         ENDIF
          !
          !                                                        ! ===============
Index: NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/TRA/traldf_lap_blp.F90
===================================================================
--- NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/TRA/traldf_lap_blp.F90	(revision 12765)
+++ NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/TRA/traldf_lap_blp.F90	(revision 12766)
@@ -171,4 +171,5 @@
       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pt_rhs     ! tracer trend
       !
+      TYPE(TILE)                            :: ktile        ! Tile indices
       INTEGER ::   ji, jj, jk, jn   ! dummy loop indices
       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt) :: zlap         ! laplacian at t-point
@@ -194,5 +195,5 @@
          CALL tra_ldf_lap  ( kt, Kmm, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, pt,     zlap, kjpt, 1 )
       CASE ( np_blp_i  )               ! rotated   bilaplacian : standard operator (Madec)
-         CALL tra_ldf_iso  ( kt, Kmm, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, pt, pt, zlap, kjpt, 1 )
+         CALL tra_ldf_iso  ( ktile, kt, Kmm, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, pt, pt, zlap, kjpt, 1 )
       CASE ( np_blp_it )               ! rotated  bilaplacian : triad operator (griffies)
          CALL tra_ldf_triad( kt, Kmm, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, pt, pt, zlap, kjpt, 1 )
@@ -210,5 +211,5 @@
          CALL tra_ldf_lap  ( kt, Kmm, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, pt_rhs,         kjpt, 2 )
       CASE ( np_blp_i  )               ! rotated   bilaplacian : standard operator (Madec)
-         CALL tra_ldf_iso  ( kt, Kmm, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, pt    , pt_rhs, kjpt, 2 )
+         CALL tra_ldf_iso  ( ktile, kt, Kmm, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, pt    , pt_rhs, kjpt, 2 )
       CASE ( np_blp_it )               ! rotated  bilaplacian : triad operator (griffies)
          CALL tra_ldf_triad( kt, Kmm, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, pt    , pt_rhs, kjpt, 2 )
Index: NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/do_loop_substitute.h90
===================================================================
--- NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/do_loop_substitute.h90	(revision 12765)
+++ NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/do_loop_substitute.h90	(revision 12766)
@@ -102,2 +102,49 @@
 #define END_2D   END DO   ;   END DO
 #define END_3D   END DO   ;   END DO   ;   END DO
+
+#define __kIs_T     ktile % ntsi
+#define __kJs_T     ktile % ntsj
+#define __kIsm1_T   ktile % ntsim1
+#define __kJsm1_T   ktile % ntsjm1
+
+#define __kIe_T     ktile % ntei
+#define __kJe_T     ktile % ntej
+#define __kIep1_T   ktile % nteip1
+#define __kJep1_T   ktile % ntejp1
+
+#define IND_2D          __kIsm1_T:__kIep1_T,__kJsm1_T:__kJep1_T
+
+#define DO_2D_00_00_T   DO jj = __kJs_T, __kJe_T   ;   DO ji = __kIs_T, __kIe_T
+#define DO_2D_00_01_T   DO jj = __kJs_T, __kJe_T   ;   DO ji = __kIs_T, __kIep1_T
+#define DO_2D_00_10_T   DO jj = __kJs_T, __kJe_T   ;   DO ji = __kIsm1_T, __kIe_T
+#define DO_2D_00_11_T   DO jj = __kJs_T, __kJe_T   ;   DO ji = __kIsm1_T, __kIep1_T
+
+#define DO_2D_01_00_T   DO jj = __kJs_T, __kJep1_T   ;   DO ji = __kIs_T, __kIe_T
+#define DO_2D_01_01_T   DO jj = __kJs_T, __kJep1_T   ;   DO ji = __kIs_T, __kIep1_T
+#define DO_2D_01_10_T   DO jj = __kJs_T, __kJep1_T   ;   DO ji = __kIsm1_T, __kIe_T
+#define DO_2D_01_11_T   DO jj = __kJs_T, __kJep1_T   ;   DO ji = __kIsm1_T, __kIep1_T
+
+#define DO_2D_10_00_T   DO jj = __kJsm1_T, __kJe_T   ;   DO ji = __kIs_T, __kIe_T
+#define DO_2D_10_10_T   DO jj = __kJsm1_T, __kJe_T   ;   DO ji = __kIsm1_T, __kIe_T
+#define DO_2D_10_11_T   DO jj = __kJsm1_T, __kJe_T   ;   DO ji = __kIsm1_T, __kIep1_T
+
+#define DO_2D_11_00_T   DO jj = __kJsm1_T, __kJep1_T   ;   DO ji = __kIs_T, __kIe_T
+#define DO_2D_11_01_T   DO jj = __kJsm1_T, __kJep1_T   ;   DO ji = __kIs_T, __kIep1_T
+#define DO_2D_11_10_T   DO jj = __kJsm1_T, __kJep1_T   ;   DO ji = __kIsm1_T, __kIe_T
+#define DO_2D_11_11_T   DO jj = __kJsm1_T, __kJep1_T   ;   DO ji = __kIsm1_T, __kIep1_T
+
+#define DO_3D_00_00_T(ks,ke)   DO jk = ks, ke   ;   DO_2D_00_00_T
+#define DO_3D_00_10_T(ks,ke)   DO jk = ks, ke   ;   DO_2D_00_10_T
+
+#define DO_3D_01_01_T(ks,ke)   DO jk = ks, ke   ;   DO_2D_01_01_T
+
+#define DO_3D_10_00_T(ks,ke)   DO jk = ks, ke   ;   DO_2D_10_00_T
+#define DO_3D_10_10_T(ks,ke)   DO jk = ks, ke   ;   DO_2D_10_10_T
+#define DO_3D_10_11_T(ks,ke)   DO jk = ks, ke   ;   DO_2D_10_11_T
+
+#define DO_3D_11_11_T(ks,ke)   DO jk = ks, ke   ;   DO_2D_11_11_T
+
+#define DO_3DS_00_00_T(ks,ke,ki)   DO jk = ks, ke, ki   ;   DO_2D_00_00_T
+#define DO_3DS_01_01_T(ks,ke,ki)   DO jk = ks, ke, ki   ;   DO_2D_01_01_T
+#define DO_3DS_10_10_T(ks,ke,ki)   DO jk = ks, ke, ki   ;   DO_2D_10_10_T
+#define DO_3DS_11_11_T(ks,ke,ki)   DO jk = ks, ke, ki   ;   DO_2D_11_11_T
Index: NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/par_oce.F90
===================================================================
--- NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/par_oce.F90	(revision 12765)
+++ NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/par_oce.F90	(revision 12766)
@@ -62,4 +62,9 @@
    INTEGER, PUBLIC ::   jpjmax! = ( jpjglo-2*nn_hls + (jpnj-1) ) / jpnj + 2*nn_hls !: maximum jpj
 
+   ! Tiling decomposition
+   INTEGER, PUBLIC ::   jpnitile                !: number of tiles following i
+   INTEGER, PUBLIC ::   jpnjtile                !: number of tiles following j
+   INTEGER, PUBLIC ::   jpnijtile               !: number of tiles in total (jpnitile x jpnjtile)
+
    !!---------------------------------------------------------------------
    !! Active tracer parameters
Index: NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/step.F90
===================================================================
--- NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/step.F90	(revision 12765)
+++ NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/step.F90	(revision 12766)
@@ -81,5 +81,6 @@
       !!              -8- Outputs and diagnostics
       !!----------------------------------------------------------------------
-      INTEGER ::   ji, jj, jk   ! dummy loop indice
+      TYPE(TILE)   :: stile     ! Tile indices
+      INTEGER ::   ji, jj, jk, jtile   ! dummy loop indice
       INTEGER ::   indic        ! error indicator if < 0
 !!gm kcall can be removed, I guess
@@ -264,6 +265,10 @@
       IF( lrst_oce .AND. ln_zdfosm ) &
            &             CALL osm_rst    ( kstp,      Nnn, 'WRITE'  )  ! write OSMOSIS outputs + ww (so must do here) to restarts
-                         CALL tra_ldf    ( kstp, Nbb, Nnn, ts, Nrhs )  ! lateral mixing
-
+
+      ! Loop over tile domains
+      DO jtile = 1, jpnijtile
+                         CALL dom_tile( jtile, stile )
+                         CALL tra_ldf    ( stile, kstp, Nbb, Nnn, ts, Nrhs ) ! lateral mixing
+      END DO
                          CALL tra_zdf    ( kstp, Nbb, Nnn, Nrhs, ts, Naa  )  ! vertical mixing and after tracer fields
       IF( ln_zdfnpc  )   CALL tra_npc    ( kstp,      Nnn, Nrhs, ts, Naa  )  ! update after fields by non-penetrative convection
Index: NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/step_oce.F90
===================================================================
--- NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/step_oce.F90	(revision 12765)
+++ NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/step_oce.F90	(revision 12766)
@@ -62,4 +62,5 @@
    USE domvvl          ! variable vertical scale factors  (dom_vvl_sf_nxt routine)
    !                                                       (dom_vvl_sf_swp routine)
+   USE domain   , ONLY : dom_tile
 
    USE ldfslp          ! iso-neutral slopes               (ldf_slp routine)
