Index: NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/DIA/diaar5.F90
===================================================================
--- NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/DIA/diaar5.F90	(revision 13053)
+++ NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/DIA/diaar5.F90	(revision 13054)
@@ -34,4 +34,5 @@
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:  ) ::   thick0       ! ocean thickness (interior domain)
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   sn0          ! initial salinity
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   hstr_adv, hstr_ldf
 
    LOGICAL  :: l_ar5
@@ -53,5 +54,6 @@
       !!----------------------------------------------------------------------
       !
-      ALLOCATE( thick0(jpi,jpj) , sn0(jpi,jpj,jpk) , STAT=dia_ar5_alloc )
+      ALLOCATE( thick0(jpi,jpj) , sn0(jpi,jpj,jpk) , &
+         &      hstr_adv(jpi,jpj,jpts,2), hstr_ldf(jpi,jpj,jpts,2), STAT=dia_ar5_alloc )
       !
       CALL mpp_sum ( 'diaar5', dia_ar5_alloc )
@@ -301,6 +303,6 @@
    END SUBROUTINE dia_ar5
 
-
-   SUBROUTINE dia_ar5_hst( ktra, cptr, puflx, pvflx ) 
+   ! TODO: These changes and lbc_lnk not necessary if using XIOS (subdomain support, will not output haloes)
+   SUBROUTINE dia_ar5_hst( ktra, cptr, puflx, pvflx )
       !!----------------------------------------------------------------------
       !!                    ***  ROUTINE dia_ar5_htr ***
@@ -311,39 +313,50 @@
       INTEGER                         , INTENT(in )  :: ktra  ! tracer index
       CHARACTER(len=3)                , INTENT(in)   :: cptr  ! transport type  'adv'/'ldf'
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)   :: puflx  ! u-flux of advection/diffusion
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)   :: pvflx  ! v-flux of advection/diffusion
+      REAL(wp), DIMENSION(A2D,jpk)    , INTENT(in)   :: puflx  ! u-flux of advection/diffusion
+      REAL(wp), DIMENSION(A2D,jpk)    , INTENT(in)   :: pvflx  ! v-flux of advection/diffusion
       !
       INTEGER    ::  ji, jj, jk
-      REAL(wp), DIMENSION(jpi,jpj)  :: z2d
-
-    
-      z2d(:,:) = puflx(:,:,1) 
-      DO_3D_00_00( 1, jpkm1 )
-         z2d(ji,jj) = z2d(ji,jj) + puflx(ji,jj,jk) 
-      END_3D
-       CALL lbc_lnk( 'diaar5', z2d, 'U', -1. )
-       IF( cptr == 'adv' ) THEN
-          IF( ktra == jp_tem ) CALL iom_put( 'uadv_heattr' , rho0_rcp * z2d )  ! advective heat transport in i-direction
-          IF( ktra == jp_sal ) CALL iom_put( 'uadv_salttr' , rho0     * z2d )  ! advective salt transport in i-direction
-       ENDIF
-       IF( cptr == 'ldf' ) THEN
-          IF( ktra == jp_tem ) CALL iom_put( 'udiff_heattr' , rho0_rcp * z2d ) ! diffusive heat transport in i-direction
-          IF( ktra == jp_sal ) CALL iom_put( 'udiff_salttr' , rho0     * z2d ) ! diffusive salt transport in i-direction
-       ENDIF
-       !
-       z2d(:,:) = pvflx(:,:,1) 
-       DO_3D_00_00( 1, jpkm1 )
-          z2d(ji,jj) = z2d(ji,jj) + pvflx(ji,jj,jk) 
-       END_3D
-       CALL lbc_lnk( 'diaar5', z2d, 'V', -1. )
-       IF( cptr == 'adv' ) THEN
-          IF( ktra == jp_tem ) CALL iom_put( 'vadv_heattr' , rho0_rcp * z2d )  ! advective heat transport in j-direction
-          IF( ktra == jp_sal ) CALL iom_put( 'vadv_salttr' , rho0     * z2d )  ! advective salt transport in j-direction
-       ENDIF
-       IF( cptr == 'ldf' ) THEN
-          IF( ktra == jp_tem ) CALL iom_put( 'vdiff_heattr' , rho0_rcp * z2d ) ! diffusive heat transport in j-direction
-          IF( ktra == jp_sal ) CALL iom_put( 'vdiff_salttr' , rho0     * z2d ) ! diffusive salt transport in j-direction
-       ENDIF
-          
+
+      IF( cptr /= 'adv' .AND. cptr /= 'ldf' ) RETURN
+      IF( ktra /= jp_tem .AND. ktra /= jp_sal ) RETURN
+
+      IF( cptr == 'adv' ) THEN
+         DO_2D_00_00
+            hstr_adv(ji,jj,ktra,1) = puflx(ji,jj,1)
+            hstr_adv(ji,jj,ktra,2) = pvflx(ji,jj,1)
+         END_2D
+         DO_3D_00_00( 1, jpkm1 )
+            hstr_adv(ji,jj,ktra,1) = hstr_adv(ji,jj,ktra,1) + puflx(ji,jj,jk)
+            hstr_adv(ji,jj,ktra,2) = hstr_adv(ji,jj,ktra,2) + pvflx(ji,jj,jk)
+         END_3D
+      ELSE IF( cptr == 'ldf' ) THEN
+         DO_2D_00_00
+            hstr_ldf(ji,jj,ktra,1) = puflx(ji,jj,1)
+            hstr_ldf(ji,jj,ktra,2) = pvflx(ji,jj,1)
+         END_2D
+         DO_3D_00_00( 1, jpkm1 )
+            hstr_ldf(ji,jj,ktra,1) = hstr_ldf(ji,jj,ktra,1) + puflx(ji,jj,jk)
+            hstr_ldf(ji,jj,ktra,2) = hstr_ldf(ji,jj,ktra,2) + pvflx(ji,jj,jk)
+         END_3D
+      ENDIF
+
+      IF( ntile == 0 .OR. ntile == nijtile ) THEN
+         IF( cptr == 'adv' ) THEN
+            CALL lbc_lnk( 'diaar5', hstr_adv(:,:,ktra,1), 'U', -1. )
+            IF( ktra == jp_tem ) CALL iom_put( 'uadv_heattr' , rho0_rcp * hstr_adv(:,:,ktra,1) )  ! advective heat transport in i-direction
+            IF( ktra == jp_sal ) CALL iom_put( 'uadv_salttr' , rho0     * hstr_adv(:,:,ktra,1) )  ! advective salt transport in i-direction
+            CALL lbc_lnk( 'diaar5', hstr_adv(:,:,ktra,2), 'V', -1. )
+            IF( ktra == jp_tem ) CALL iom_put( 'vadv_heattr' , rho0_rcp * hstr_adv(:,:,ktra,2) )  ! advective heat transport in j-direction
+            IF( ktra == jp_sal ) CALL iom_put( 'vadv_salttr' , rho0     * hstr_adv(:,:,ktra,2) )  ! advective salt transport in j-direction
+         ENDIF
+         IF( cptr == 'ldf' ) THEN
+            CALL lbc_lnk( 'diaar5', hstr_ldf(:,:,ktra,1), 'U', -1. )
+            IF( ktra == jp_tem ) CALL iom_put( 'udiff_heattr' , rho0_rcp * hstr_ldf(:,:,ktra,1) ) ! diffusive heat transport in i-direction
+            IF( ktra == jp_sal ) CALL iom_put( 'udiff_salttr' , rho0     * hstr_ldf(:,:,ktra,1) ) ! diffusive salt transport in i-direction
+            CALL lbc_lnk( 'diaar5', hstr_ldf(:,:,ktra,2), 'V', -1. )
+            IF( ktra == jp_tem ) CALL iom_put( 'vdiff_heattr' , rho0_rcp * hstr_ldf(:,:,ktra,2) ) ! diffusive heat transport in j-direction
+            IF( ktra == jp_sal ) CALL iom_put( 'vdiff_salttr' , rho0     * hstr_ldf(:,:,ktra,2) ) ! diffusive salt transport in j-direction
+         ENDIF
+      ENDIF
    END SUBROUTINE dia_ar5_hst
 
@@ -367,5 +380,9 @@
       IF(   iom_use( 'voltot'  ) .OR. iom_use( 'sshtot'    )  .OR. iom_use( 'sshdyn' )  .OR.  & 
          &  iom_use( 'masstot' ) .OR. iom_use( 'temptot'   )  .OR. iom_use( 'saltot' ) .OR.  &    
-         &  iom_use( 'botpres' ) .OR. iom_use( 'sshthster' )  .OR. iom_use( 'sshsteric' )  ) L_ar5 = .TRUE.
+         &  iom_use( 'botpres' ) .OR. iom_use( 'sshthster' )  .OR. iom_use( 'sshsteric' ) .OR. &
+         &  iom_use( 'uadv_heattr' ) .OR. iom_use( 'udiff_heattr' ) .OR. &
+         &  iom_use( 'uadv_salttr' ) .OR. iom_use( 'udiff_salttr' ) .OR. &
+         &  iom_use( 'vadv_heattr' ) .OR. iom_use( 'vdiff_heattr' ) .OR. &
+         &  iom_use( 'vadv_salttr' ) .OR. iom_use( 'vdiff_salttr' ) ) L_ar5 = .TRUE.
   
       IF( l_ar5 ) THEN
Index: NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/DIA/diaptr.F90
===================================================================
--- NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/DIA/diaptr.F90	(revision 13053)
+++ NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/DIA/diaptr.F90	(revision 13054)
@@ -22,4 +22,5 @@
    USE oce              ! ocean dynamics and active tracers
    USE dom_oce          ! ocean space and time domain
+   USE domain, ONLY : dom_tile
    USE phycst           ! physical constants
    !
@@ -32,4 +33,8 @@
    PRIVATE
 
+   INTERFACE ptr_sum
+      MODULE PROCEDURE ptr_sum_3d, ptr_sum_2d
+   END INTERFACE
+
    INTERFACE ptr_sj
       MODULE PROCEDURE ptr_sj_3d, ptr_sj_2d
@@ -43,9 +48,12 @@
 
    !                                  !!** namelist  namptr  **
-   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hstr_adv, hstr_ldf, hstr_eiv   !: Heat/Salt TRansports(adv, diff, Bolus.)
-   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hstr_ove, hstr_btr, hstr_vtr   !: heat Salt TRansports(overturn, baro, merional)
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   hstr_adv, hstr_ldf, hstr_eiv   !: Heat/Salt TRansports(adv, diff, Bolus.)
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   hstr_ove, hstr_btr, hstr_vtr   !: heat Salt TRansports(overturn, baro, merional)
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) ::   pvtr_int, pzon_int             !: Other zonal integrals
 
    LOGICAL , PUBLIC ::   l_diaptr        !: tracers  trend flag (set from namelist in trdini)
    INTEGER, PARAMETER, PUBLIC ::   nptr = 5  ! (glo, atl, pac, ind, ipc)
+   INTEGER, PARAMETER         ::   jp_msk = 3
+   INTEGER, PARAMETER         ::   jp_vtr = 4
 
    REAL(wp) ::   rc_sv    = 1.e-6_wp   ! conversion from m3/s to Sverdrup
@@ -55,7 +63,4 @@
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: btmsk   ! T-point basin interior masks
    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: btmsk34 ! mask out Southern Ocean (=0 south of 34°S)
-
-   REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:)   :: p_fval1d
-   REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:,:) :: p_fval2d
 
    LOGICAL ::   ll_init = .TRUE.        !: tracers  trend flag (set from namelist in trdini)
@@ -69,4 +74,5 @@
 CONTAINS
 
+   ! TODO: Most changes and some code in this module not necessary if using XIOS (subdomain support, axis operations)
    SUBROUTINE dia_ptr( kt, Kmm, pvtr )
       !!----------------------------------------------------------------------
@@ -78,9 +84,6 @@
       !
       INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
-      REAL(wp) ::   zsfc,zvfc               ! local scalar
+      INTEGER  ::   itile
       REAL(wp), DIMENSION(jpi,jpj)     ::  z2d   ! 2D workspace
-      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zmask   ! 3D workspace
-      REAL(wp), DIMENSION(jpi,jpj,jpk) ::  z3d    ! 3D workspace
-      REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) ::  zts   ! 3D workspace
       REAL(wp), DIMENSION(jpj)      ::  zvsum, ztsum, zssum   ! 1D workspace
       !
@@ -99,9 +102,19 @@
       IF( .NOT. l_diaptr )   RETURN
 
+      ! Calculate zonal integrals
+      IF( PRESENT( pvtr ) ) THEN
+         CALL dia_ptr_zint( Kmm, pvtr )
+      ELSE
+         CALL dia_ptr_zint( Kmm )
+      ENDIF
+
+      ! Calculate diagnostics only when zonal integrals have finished
+      IF( ntile /= 0 .AND. ntile /= nijtile )  RETURN
+
       IF( PRESENT( pvtr ) ) THEN
          IF( iom_use( 'zomsf' ) ) THEN    ! effective MSF
             DO jn = 1, nptr                                    ! by sub-basins
-               z4d1(1,:,:,jn) =  ptr_sjk( pvtr(:,:,:), btmsk34(:,:,jn) )  ! zonal cumulative effective transport excluding closed seas
-               DO jk = jpkm1, 1, -1 
+               z4d1(1,:,:,jn) =  pvtr_int(:,:,jp_vtr,jn)                  ! zonal cumulative effective transport excluding closed seas
+               DO jk = jpkm1, 1, -1
                   z4d1(1,:,jk,jn) = z4d1(1,:,jk+1,jn) - z4d1(1,:,jk,jn)    ! effective j-Stream-Function (MSF)
                END DO
@@ -112,25 +125,13 @@
             CALL iom_put( 'zomsf', z4d1 * rc_sv )
          ENDIF
-         IF(  iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR.   &
-            & iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN
-            ! define fields multiplied by scalar
-            zmask(:,:,:) = 0._wp
-            zts(:,:,:,:) = 0._wp
-            DO_3D_10_11( 1, jpkm1 )
-               zvfc = e1v(ji,jj) * e3v(ji,jj,jk,Kmm)
-               zmask(ji,jj,jk)      = vmask(ji,jj,jk)      * zvfc
-               zts(ji,jj,jk,jp_tem) = (ts(ji,jj,jk,jp_tem,Kmm)+ts(ji,jj+1,jk,jp_tem,Kmm)) * 0.5 * zvfc  !Tracers averaged onto V grid
-               zts(ji,jj,jk,jp_sal) = (ts(ji,jj,jk,jp_sal,Kmm)+ts(ji,jj+1,jk,jp_sal,Kmm)) * 0.5 * zvfc
-            END_3D
-         ENDIF
          IF( iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) ) THEN
             DO jn = 1, nptr
-               sjk(:,:,jn) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) )
+               sjk(:,:,jn) = pvtr_int(:,:,jp_msk,jn)
                r1_sjk(:,:,jn) = 0._wp
                WHERE( sjk(:,:,jn) /= 0._wp )   r1_sjk(:,:,jn) = 1._wp / sjk(:,:,jn)
                ! i-mean T and S, j-Stream-Function, basin
-               zt_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) * r1_sjk(:,:,jn)
-               zs_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) * r1_sjk(:,:,jn)
-               v_msf(:,:,jn) = ptr_sjk( pvtr(:,:,:), btmsk34(:,:,jn) ) 
+               zt_jk(:,:,jn) = pvtr_int(:,:,jp_tem,jn) * r1_sjk(:,:,jn)
+               zs_jk(:,:,jn) = pvtr_int(:,:,jp_sal,jn) * r1_sjk(:,:,jn)
+               v_msf(:,:,jn) = pvtr_int(:,:,jp_vtr,jn)
                hstr_ove(:,jp_tem,jn) = SUM( v_msf(:,:,jn)*zt_jk(:,:,jn), 2 )
                hstr_ove(:,jp_sal,jn) = SUM( v_msf(:,:,jn)*zs_jk(:,:,jn), 2 )
@@ -156,11 +157,11 @@
             ! Calculate barotropic heat and salt transport here 
             DO jn = 1, nptr
-               sjk(:,1,jn) = ptr_sj( zmask(:,:,:), btmsk(:,:,jn) )
+               sjk(:,1,jn) = SUM( pvtr_int(:,:,jp_msk,jn), 2 )
                r1_sjk(:,1,jn) = 0._wp
                WHERE( sjk(:,1,jn) /= 0._wp )   r1_sjk(:,1,jn) = 1._wp / sjk(:,1,jn)
                !
-               zvsum(:) = ptr_sj( pvtr(:,:,:), btmsk34(:,:,jn) )
-               ztsum(:) = ptr_sj( zts(:,:,:,jp_tem), btmsk(:,:,jn) )
-               zssum(:) = ptr_sj( zts(:,:,:,jp_sal), btmsk(:,:,jn) )
+               zvsum(:) =    SUM( pvtr_int(:,:,jp_vtr,jn), 2 )
+               ztsum(:) =    SUM( pvtr_int(:,:,jp_tem,jn), 2 )
+               zssum(:) =    SUM( pvtr_int(:,:,jp_sal,jn), 2 )
                hstr_btr(:,jp_tem,jn) = zvsum(:) * ztsum(:) * r1_sjk(:,1,jn)
                hstr_btr(:,jp_sal,jn) = zvsum(:) * zssum(:) * r1_sjk(:,1,jn)
@@ -183,26 +184,21 @@
          ENDIF 
          !
+         hstr_ove(:,:,:) = 0._wp       ! Zero before next timestep
+         hstr_btr(:,:,:) = 0._wp
+         pvtr_int(:,:,:,:) = 0._wp
       ELSE
-         !
-         zmask(:,:,:) = 0._wp
-         zts(:,:,:,:) = 0._wp
-         IF( iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. iom_use( 'zosrf' )  ) THEN    ! i-mean i-k-surface 
-            DO_3D_11_11( 1, jpkm1 )
-               zsfc = e1t(ji,jj) * e3t(ji,jj,jk,Kmm)
-               zmask(ji,jj,jk)      = tmask(ji,jj,jk)      * zsfc
-               zts(ji,jj,jk,jp_tem) = ts(ji,jj,jk,jp_tem,Kmm) * zsfc
-               zts(ji,jj,jk,jp_sal) = ts(ji,jj,jk,jp_sal,Kmm) * zsfc
-            END_3D
+         IF( iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. iom_use( 'zosrf' )  ) THEN    ! i-mean i-k-surface
             !
             DO jn = 1, nptr
-               zmask(1,:,:) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) )
-               z4d1(:,:,:,jn) = zmask(:,:,:)
+               z4d1(1,:,:,jn) = pzon_int(:,:,jp_msk,jn)
+               DO ji = 2, jpi
+                  z4d1(ji,:,:,jn) = z4d1(1,:,:,jn)
+               ENDDO
             ENDDO
             CALL iom_put( 'zosrf', z4d1 )
             !
             DO jn = 1, nptr
-               z4d2(1,:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) &
-                  &            / MAX( z4d1(1,:,:,jn), 10.e-15 )
-               DO ji = 1, jpi
+               z4d2(1,:,:,jn) = pzon_int(:,:,jp_tem,jn) / MAX( z4d1(1,:,:,jn), 10.e-15 )
+               DO ji = 2, jpi
                   z4d2(ji,:,:,jn) = z4d2(1,:,:,jn)
                ENDDO
@@ -211,7 +207,6 @@
             !
             DO jn = 1, nptr
-               z4d2(1,:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) &
-                  &            / MAX( z4d1(1,:,:,jn), 10.e-15 )
-               DO ji = 1, jpi
+               z4d2(1,:,:,jn) = pzon_int(:,:,jp_sal,jn) / MAX( z4d1(1,:,:,jn), 10.e-15 )
+               DO ji = 2, jpi
                   z4d2(ji,:,:,jn) = z4d2(1,:,:,jn)
                ENDDO
@@ -277,12 +272,4 @@
          !
          IF( iom_use( 'sopstvtr' ) .OR. iom_use( 'sophtvtr' ) ) THEN
-            zts(:,:,:,:) = 0._wp
-            DO_3D_10_11( 1, jpkm1 )
-               zvfc = e1v(ji,jj) * e3v(ji,jj,jk,Kmm)
-               zts(ji,jj,jk,jp_tem) = (ts(ji,jj,jk,jp_tem,Kmm)+ts(ji,jj+1,jk,jp_tem,Kmm)) * 0.5 * zvfc  !Tracers averaged onto V grid
-               zts(ji,jj,jk,jp_sal) = (ts(ji,jj,jk,jp_sal,Kmm)+ts(ji,jj+1,jk,jp_sal,Kmm)) * 0.5 * zvfc
-            END_3D
-             CALL dia_ptr_hst( jp_tem, 'vtr', zts(:,:,:,jp_tem) )
-             CALL dia_ptr_hst( jp_sal, 'vtr', zts(:,:,:,jp_sal) )
              DO jn = 1, nptr
                 z3dtr(1,:,jn) = hstr_vtr(:,jp_tem,jn) * rc_pwatt  !  (conversion in PW)
@@ -301,10 +288,19 @@
          ENDIF
          !
+         ! TODO: Possibly not necessary if using XIOS (if cumulative axis operations are possible)
          IF( iom_use( 'uocetr_vsum_cumul' ) ) THEN
+            itile = ntile
+            CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 )                            ! Use full domain
             CALL iom_get_var(  'uocetr_vsum_op', z2d ) ! get uocetr_vsum_op from xml
             z2d(:,:) = ptr_ci_2d( z2d(:,:) )  
             CALL iom_put( 'uocetr_vsum_cumul', z2d )
-         ENDIF
-         !
+            IF( ntile /= itile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = itile )   ! Revert to tile domain
+         ENDIF
+         !
+         hstr_adv(:,:,:) = 0._wp       ! Zero before next timestep
+         hstr_ldf(:,:,:) = 0._wp
+         hstr_eiv(:,:,:) = 0._wp
+         hstr_vtr(:,:,:) = 0._wp
+         pzon_int(:,:,:,:) = 0._wp
       ENDIF
       !
@@ -312,4 +308,99 @@
       !
    END SUBROUTINE dia_ptr
+
+
+   SUBROUTINE dia_ptr_zint( Kmm, pvtr )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE dia_ptr_zint ***
+      !!----------------------------------------------------------------------
+      !! ** Purpose : i and i-k sum operations on arrays
+      !!
+      !! ** Method  : - Call ptr_sjk (i sum) or ptr_sj (i-k sum) to perform the sum operation
+      !!              - Call ptr_sum to add this result to the sum over tiles
+      !!
+      !! ** Action  : pvtr_int - terms for volume streamfunction, heat/salt transport barotropic/overturning terms
+      !!              pzon_int - terms for i mean temperature/salinity
+      !!----------------------------------------------------------------------
+      INTEGER                     , INTENT(in)           :: Kmm          ! time level index
+      REAL(wp), DIMENSION(A2D,jpk), INTENT(in), OPTIONAL :: pvtr         ! j-effective transport
+      REAL(wp), DIMENSION(A2D,jpk)                       :: zmask        ! 3D workspace
+      REAL(wp), DIMENSION(A2D,jpk,jpts)                  :: zts          ! 4D workspace
+      REAL(wp), DIMENSION(A1Dj,jpk,nptr)                 :: sjk, v_msf   ! Zonal sum: i-k surface area, j-effective transport
+      REAL(wp), DIMENSION(A1Dj,jpk,nptr)                 :: zt_jk, zs_jk ! Zonal sum: i-k surface area * (T, S)
+      REAL(wp)                                           :: zsfc, zvfc   ! i-k surface area
+      INTEGER  ::   ji, jj, jk, jn                                       ! dummy loop indices
+
+      IF( PRESENT( pvtr ) ) THEN
+         ! i sum of effective j transport excluding closed seas
+         IF( iom_use( 'zomsf' ) .OR. iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) ) THEN
+            DO jn = 1, nptr
+               v_msf(:,:,jn) = ptr_sjk( pvtr(:,:,:), btmsk34(:,:,jn) )
+            ENDDO
+
+            CALL ptr_sum( pvtr_int(:,:,jp_vtr,:), v_msf(:,:,:) )
+         ENDIF
+
+         ! i sum of j surface area, j surface area - temperature/salinity product on V grid
+         IF(  iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR.   &
+            & iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN
+            zmask(:,:,:) = 0._wp
+            zts(:,:,:,:) = 0._wp
+
+            DO_3D_10_11( 1, jpkm1 )
+               zvfc = e1v(ji,jj) * e3v(ji,jj,jk,Kmm)
+               zmask(ji,jj,jk)      = vmask(ji,jj,jk)      * zvfc
+               zts(ji,jj,jk,jp_tem) = (ts(ji,jj,jk,jp_tem,Kmm)+ts(ji,jj+1,jk,jp_tem,Kmm)) * 0.5 * zvfc !Tracers averaged onto V grid
+               zts(ji,jj,jk,jp_sal) = (ts(ji,jj,jk,jp_sal,Kmm)+ts(ji,jj+1,jk,jp_sal,Kmm)) * 0.5 * zvfc
+            END_3D
+
+            DO jn = 1, nptr
+               sjk(:,:,jn)   = ptr_sjk( zmask(:,:,:)     , btmsk(:,:,jn) )
+               zt_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) )
+               zs_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) )
+            ENDDO
+
+            CALL ptr_sum( pvtr_int(:,:,jp_msk,:), sjk(:,:,:)   )
+            CALL ptr_sum( pvtr_int(:,:,jp_tem,:), zt_jk(:,:,:) )
+            CALL ptr_sum( pvtr_int(:,:,jp_sal,:), zs_jk(:,:,:) )
+         ENDIF
+      ELSE
+         ! i sum of j surface area - temperature/salinity product on T grid
+         IF( iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. iom_use( 'zosrf' )  ) THEN
+            zmask(:,:,:) = 0._wp
+            zts(:,:,:,:) = 0._wp
+
+            DO_3D_11_11( 1, jpkm1 )
+               zsfc = e1t(ji,jj) * e3t(ji,jj,jk,Kmm)
+               zmask(ji,jj,jk)      = tmask(ji,jj,jk)      * zsfc
+               zts(ji,jj,jk,jp_tem) = ts(ji,jj,jk,jp_tem,Kmm) * zsfc
+               zts(ji,jj,jk,jp_sal) = ts(ji,jj,jk,jp_sal,Kmm) * zsfc
+            END_3D
+
+            DO jn = 1, nptr
+               sjk(:,:,jn)   = ptr_sjk( zmask(:,:,:)     , btmsk(:,:,jn) )
+               zt_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) )
+               zs_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) )
+            ENDDO
+
+            CALL ptr_sum( pzon_int(:,:,jp_msk,:), sjk(:,:,:)   )
+            CALL ptr_sum( pzon_int(:,:,jp_tem,:), zt_jk(:,:,:) )
+            CALL ptr_sum( pzon_int(:,:,jp_sal,:), zs_jk(:,:,:) )
+         ENDIF
+
+         ! i-k sum of j surface area - temperature/salinity product on V grid
+         IF( iom_use( 'sopstvtr' ) .OR. iom_use( 'sophtvtr' ) ) THEN
+            zts(:,:,:,:) = 0._wp
+
+            DO_3D_10_11( 1, jpkm1 )
+               zvfc = e1v(ji,jj) * e3v(ji,jj,jk,Kmm)
+               zts(ji,jj,jk,jp_tem) = (ts(ji,jj,jk,jp_tem,Kmm)+ts(ji,jj+1,jk,jp_tem,Kmm)) * 0.5 * zvfc  !Tracers averaged onto V grid
+               zts(ji,jj,jk,jp_sal) = (ts(ji,jj,jk,jp_sal,Kmm)+ts(ji,jj+1,jk,jp_sal,Kmm)) * 0.5 * zvfc
+            END_3D
+
+            CALL dia_ptr_hst( jp_tem, 'vtr', zts(:,:,:,jp_tem) )
+            CALL dia_ptr_hst( jp_sal, 'vtr', zts(:,:,:,jp_sal) )
+         ENDIF
+      ENDIF
+   END SUBROUTINE dia_ptr_zint
 
 
@@ -380,4 +471,6 @@
          hstr_btr(:,:,:) = 0._wp           !
          hstr_vtr(:,:,:) = 0._wp           !
+         pvtr_int(:,:,:,:) = 0._wp
+         pzon_int(:,:,:,:) = 0._wp
          !
          ll_init = .FALSE.
@@ -397,61 +490,101 @@
       INTEGER                         , INTENT(in )  :: ktra  ! tracer index
       CHARACTER(len=3)                , INTENT(in)   :: cptr  ! transport type  'adv'/'ldf'/'eiv'
-      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in)   :: pvflx   ! 3D input array of advection/diffusion
+      REAL(wp), DIMENSION(A2D,jpk)    , INTENT(in)   :: pvflx ! 3D input array of advection/diffusion
+      REAL(wp), DIMENSION(A1Dj,nptr)                 :: zsj   !
       INTEGER                                        :: jn    !
 
+      DO jn = 1, nptr
+         zsj(:,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
+      ENDDO
       !
       IF( cptr == 'adv' ) THEN
-         IF( ktra == jp_tem )  THEN
-             DO jn = 1, nptr
-                hstr_adv(:,jp_tem,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
-         IF( ktra == jp_sal )  THEN
-             DO jn = 1, nptr
-                hstr_adv(:,jp_sal,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
+         IF( ktra == jp_tem )  CALL ptr_sum( hstr_adv(:,jp_tem,:), zsj(:,:) )
+         IF( ktra == jp_sal )  CALL ptr_sum( hstr_adv(:,jp_sal,:), zsj(:,:) )
+      ELSE IF( cptr == 'ldf' ) THEN
+         IF( ktra == jp_tem )  CALL ptr_sum( hstr_ldf(:,jp_tem,:), zsj(:,:) )
+         IF( ktra == jp_sal )  CALL ptr_sum( hstr_ldf(:,jp_sal,:), zsj(:,:) )
+      ELSE IF( cptr == 'eiv' ) THEN
+         IF( ktra == jp_tem )  CALL ptr_sum( hstr_eiv(:,jp_tem,:), zsj(:,:) )
+         IF( ktra == jp_sal )  CALL ptr_sum( hstr_eiv(:,jp_sal,:), zsj(:,:) )
+      ELSE IF( cptr == 'vtr' ) THEN
+         IF( ktra == jp_tem )  CALL ptr_sum( hstr_vtr(:,jp_tem,:), zsj(:,:) )
+         IF( ktra == jp_sal )  CALL ptr_sum( hstr_vtr(:,jp_sal,:), zsj(:,:) )
       ENDIF
       !
-      IF( cptr == 'ldf' ) THEN
-         IF( ktra == jp_tem )  THEN
-             DO jn = 1, nptr
-                hstr_ldf(:,jp_tem,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
-         IF( ktra == jp_sal )  THEN
-             DO jn = 1, nptr
-                hstr_ldf(:,jp_sal,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
+   END SUBROUTINE dia_ptr_hst
+
+
+   SUBROUTINE ptr_sum_2d( phstr, pva )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE ptr_sum_2d ***
+      !!----------------------------------------------------------------------
+      !! ** Purpose : Add two 2D arrays with (j,nptr) dimensions
+      !!
+      !! ** Method  : - phstr = phstr + pva
+      !!              - Call mpp_sum if the final tile
+      !!
+      !! ** Action  : phstr
+      !!----------------------------------------------------------------------
+      REAL(wp), DIMENSION(jpj,nptr) , INTENT(inout)         ::  phstr  !
+      REAL(wp), DIMENSION(A1Dj,nptr), INTENT(in)            ::  pva    !
+      INTEGER                                               ::  jj
+#if defined key_mpp_mpi
+      INTEGER, DIMENSION(1)          ::  ish1d
+      INTEGER, DIMENSION(2)          ::  ish2d
+      REAL(wp), DIMENSION(jpj*nptr)  ::  zwork
+#endif
+
+      DO jj = Ntjs0, Ntje0
+         phstr(jj,:) = phstr(jj,:)  + pva(jj,:)
+      END DO
+
+#if defined key_mpp_mpi
+      IF( ntile == 0 .OR. ntile == nijtile ) THEN
+         ish1d(1) = jpj*nptr
+         ish2d(1) = jpj ; ish2d(2) = nptr
+         zwork(:) = RESHAPE( phstr(:,:), ish1d )
+         CALL mpp_sum( 'diaptr', zwork, ish1d(1), ncomm_znl )
+         phstr(:,:) = RESHAPE( zwork, ish2d )
       ENDIF
-      !
-      IF( cptr == 'eiv' ) THEN
-         IF( ktra == jp_tem )  THEN
-             DO jn = 1, nptr
-                hstr_eiv(:,jp_tem,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
-         IF( ktra == jp_sal )  THEN
-             DO jn = 1, nptr
-                hstr_eiv(:,jp_sal,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
+#endif
+   END SUBROUTINE ptr_sum_2d
+
+
+   SUBROUTINE ptr_sum_3d( phstr, pva )
+      !!----------------------------------------------------------------------
+      !!                    ***  ROUTINE ptr_sum_3d ***
+      !!----------------------------------------------------------------------
+      !! ** Purpose : Add two 3D arrays with (j,k,nptr) dimensions
+      !!
+      !! ** Method  : - phstr = phstr + pva
+      !!              - Call mpp_sum if the final tile
+      !!
+      !! ** Action  : phstr
+      !!----------------------------------------------------------------------
+      REAL(wp), DIMENSION(jpj,jpk,nptr) , INTENT(inout)     ::  phstr  !
+      REAL(wp), DIMENSION(A1Dj,jpk,nptr), INTENT(in)        ::  pva    !
+      INTEGER                                               ::  jj, jk
+#if defined key_mpp_mpi
+      INTEGER, DIMENSION(1)              ::  ish1d
+      INTEGER, DIMENSION(3)              ::  ish3d
+      REAL(wp), DIMENSION(jpj*jpk*nptr)  ::  zwork
+#endif
+
+      DO jk = 1, jpk
+         DO jj = Ntjs0, Ntje0
+            phstr(jj,jk,:) = phstr(jj,jk,:)  + pva(jj,jk,:)
+         END DO
+      END DO
+
+#if defined key_mpp_mpi
+      IF( ntile == 0 .OR. ntile == nijtile ) THEN
+         ish1d(1) = jpj*jpk*nptr
+         ish3d(1) = jpj ; ish3d(2) = jpk ; ish3d(3) = nptr
+         zwork(:) = RESHAPE( phstr(:,:,:), ish1d )
+         CALL mpp_sum( 'diaptr', zwork, ish1d(1), ncomm_znl )
+         phstr(:,:,:) = RESHAPE( zwork, ish3d )
       ENDIF
-      !
-      IF( cptr == 'vtr' ) THEN
-         IF( ktra == jp_tem )  THEN
-             DO jn = 1, nptr
-                hstr_vtr(:,jp_tem,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
-         IF( ktra == jp_sal )  THEN
-             DO jn = 1, nptr
-                hstr_vtr(:,jp_sal,jn) = ptr_sj( pvflx(:,:,:), btmsk(:,:,jn) )
-             ENDDO
-         ENDIF
-      ENDIF
-      !
-   END SUBROUTINE dia_ptr_hst
+#endif
+   END SUBROUTINE ptr_sum_3d
 
 
@@ -461,5 +594,5 @@
       !!----------------------------------------------------------------------
       INTEGER               ::   dia_ptr_alloc   ! return value
-      INTEGER, DIMENSION(3) ::   ierr
+      INTEGER, DIMENSION(2) ::   ierr
       !!----------------------------------------------------------------------
       ierr(:) = 0
@@ -471,5 +604,6 @@
             &      hstr_ldf(jpj,jpts,nptr), hstr_vtr(jpj,jpts,nptr), STAT=ierr(1)  )
             !
-         ALLOCATE( p_fval1d(jpj), p_fval2d(jpj,jpk), Stat=ierr(2))
+         ALLOCATE( pvtr_int(jpj,jpk,jpts+2,nptr), &
+            &      pzon_int(jpj,jpk,jpts+1,nptr), STAT=ierr(2) )
          !
          dia_ptr_alloc = MAXVAL( ierr )
@@ -491,23 +625,15 @@
       !! ** Action  : - p_fval: i-k-mean poleward flux of pvflx
       !!----------------------------------------------------------------------
-      REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk)  ::   pvflx  ! mask flux array at V-point
-      REAL(wp), INTENT(in), DIMENSION(jpi,jpj)      ::   pmsk   ! Optional 2D basin mask
+      REAL(wp), INTENT(in), DIMENSION(A2D,jpk)  ::   pvflx  ! mask flux array at V-point
+      REAL(wp), INTENT(in), DIMENSION(jpi,jpj)  ::   pmsk   ! Optional 2D basin mask
       !
       INTEGER                  ::   ji, jj, jk   ! dummy loop arguments
-      INTEGER                  ::   ijpj         ! ???
-      REAL(wp), POINTER, DIMENSION(:) :: p_fval  ! function value
+      REAL(wp), DIMENSION(A1Dj) :: p_fval  ! function value
       !!--------------------------------------------------------------------
       !
-      p_fval => p_fval1d
-
-      ijpj = jpj
       p_fval(:) = 0._wp
       DO_3D_00_00( 1, jpkm1 )
          p_fval(jj) = p_fval(jj) + pvflx(ji,jj,jk) * pmsk(ji,jj) * tmask_i(ji,jj)
       END_3D
-#if defined key_mpp_mpi
-      CALL mpp_sum( 'diaptr', p_fval, ijpj, ncomm_znl)
-#endif
-      !
    END FUNCTION ptr_sj_3d
 
@@ -524,23 +650,15 @@
       !! ** Action  : - p_fval: i-k-mean poleward flux of pvflx
       !!----------------------------------------------------------------------
-      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   pvflx  ! mask flux array at V-point
+      REAL(wp) , INTENT(in), DIMENSION(A2D)     ::   pvflx  ! mask flux array at V-point
       REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   pmsk   ! Optional 2D basin mask
       !
       INTEGER                  ::   ji,jj       ! dummy loop arguments
-      INTEGER                  ::   ijpj        ! ???
-      REAL(wp), POINTER, DIMENSION(:) :: p_fval ! function value
+      REAL(wp), DIMENSION(A1Dj) :: p_fval ! function value
       !!--------------------------------------------------------------------
       ! 
-      p_fval => p_fval1d
-
-      ijpj = jpj
       p_fval(:) = 0._wp
       DO_2D_00_00
          p_fval(jj) = p_fval(jj) + pvflx(ji,jj) * pmsk(ji,jj) * tmask_i(ji,jj)
       END_2D
-#if defined key_mpp_mpi
-      CALL mpp_sum( 'diaptr', p_fval, ijpj, ncomm_znl )
-#endif
-      ! 
    END FUNCTION ptr_sj_2d
 
@@ -587,33 +705,15 @@
       !!
       IMPLICIT none
-      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) ::   pta    ! mask flux array at V-point
-      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj)     ::   pmsk   ! Optional 2D basin mask
+      REAL(wp) , INTENT(in), DIMENSION(A2D,jpk) ::   pta    ! mask flux array at V-point
+      REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) ::   pmsk   ! Optional 2D basin mask
       !!
       INTEGER                           :: ji, jj, jk ! dummy loop arguments
-      REAL(wp), POINTER, DIMENSION(:,:) :: p_fval     ! return function value
-#if defined key_mpp_mpi
-      INTEGER, DIMENSION(1) ::   ish
-      INTEGER, DIMENSION(2) ::   ish2
-      INTEGER               ::   ijpjjpk
-      REAL(wp), DIMENSION(jpj*jpk) ::   zwork    ! mask flux array at V-point
-#endif
+      REAL(wp), DIMENSION(A1Dj,jpk) :: p_fval     ! return function value
       !!--------------------------------------------------------------------
       !
-      p_fval => p_fval2d
-
       p_fval(:,:) = 0._wp
-      !
       DO_3D_00_00( 1, jpkm1 )
          p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk) * pmsk(ji,jj) * tmask_i(ji,jj)
       END_3D
-      !
-#if defined key_mpp_mpi
-      ijpjjpk = jpj*jpk
-      ish(1) = ijpjjpk  ;   ish2(1) = jpj   ;   ish2(2) = jpk
-      zwork(1:ijpjjpk) = RESHAPE( p_fval, ish )
-      CALL mpp_sum( 'diaptr', zwork, ijpjjpk, ncomm_znl )
-      p_fval(:,:) = RESHAPE( zwork, ish2 )
-#endif
-      !
    END FUNCTION ptr_sjk
 
Index: NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/do_loop_substitute.h90
===================================================================
--- NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/do_loop_substitute.h90	(revision 13053)
+++ NEMO/branches/UKMO/dev_r12866_HPC-02_Daley_Tiling_trial_extra_halo/src/OCE/do_loop_substitute.h90	(revision 13054)
@@ -114,5 +114,7 @@
 #define Ntie2      Ntie0 + nn_hls
 #define Ntje2      Ntje0 + nn_hls
-#define A2D        Ntis2:Ntie2,Ntjs2:Ntje2
+#define A1Di       Ntis2:Ntie2
+#define A1Dj       Ntjs2:Ntje2
+#define A2D        A1Di,A1Dj
 
 #define DO_2D_00_00   DO jj = Ntjs0, Ntje0   ;   DO ji = Ntis0, Ntie0
