Index: NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/OCE/TRA/traqsr.F90
===================================================================
--- NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/OCE/TRA/traqsr.F90	(revision 14384)
+++ NEMO/branches/2021/dev_r14383_PISCES_NEWDEV_PISCO/src/OCE/TRA/traqsr.F90	(revision 15127)
@@ -108,5 +108,5 @@
       !
       INTEGER  ::   ji, jj, jk               ! dummy loop indices
-      INTEGER  ::   irgb, isi, iei, isj, iej ! local integers
+      INTEGER  ::   irgb                     ! local integers
       REAL(wp) ::   zchl, zcoef, z1_2        ! local scalars
       REAL(wp) ::   zc0 , zc1 , zc2 , zc3    !    -         -
@@ -121,5 +121,5 @@
       IF( ln_timing )   CALL timing_start('tra_qsr')
       !
-      IF( ntile == 0 .OR. ntile == 1 )  THEN                       ! Do only on the first tile
+      IF( .NOT. l_istiled .OR. ntile == 1 )  THEN                       ! Do only on the first tile
          IF( kt == nit000 ) THEN
             IF(lwp) WRITE(numout,*)
@@ -137,13 +137,8 @@
       !                         !  before qsr induced heat content  !
       !                         !-----------------------------------!
-      IF( ntsi == Nis0 ) THEN ; isi = nn_hls ; ELSE ; isi = 0 ; ENDIF    ! Avoid double-counting when using tiling
-      IF( ntsj == Njs0 ) THEN ; isj = nn_hls ; ELSE ; isj = 0 ; ENDIF
-      IF( ntei == Nie0 ) THEN ; iei = nn_hls ; ELSE ; iei = 0 ; ENDIF
-      IF( ntej == Nje0 ) THEN ; iej = nn_hls ; ELSE ; iej = 0 ; ENDIF
-
       IF( kt == nit000 ) THEN          !==  1st time step  ==!
          IF( ln_rstart .AND. .NOT.l_1st_euler ) THEN    ! read in restart
             z1_2 = 0.5_wp
-            IF( ntile == 0 .OR. ntile == 1 )  THEN                        ! Do only on the first tile
+            IF( .NOT. l_istiled .OR. ntile == 1 )  THEN                        ! Do only on the first tile
                IF(lwp) WRITE(numout,*) '          nit000-1 qsr tracer content forcing field read in the restart file'
                CALL iom_get( numror, jpdom_auto, 'qsr_hc_b', qsr_hc_b )   ! before heat content trend due to Qsr flux
@@ -151,5 +146,5 @@
          ELSE                                           ! No restart or Euler forward at 1st time step
             z1_2 = 1._wp
-            DO_3D( isi, iei, isj, iej, 1, jpk )
+            DO_3D_OVR( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpk )
                qsr_hc_b(ji,jj,jk) = 0._wp
             END_3D
@@ -157,5 +152,5 @@
       ELSE                             !==  Swap of qsr heat content  ==!
          z1_2 = 0.5_wp
-         DO_3D( isi, iei, isj, iej, 1, jpk )
+         DO_3D_OVR( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpk )
             qsr_hc_b(ji,jj,jk) = qsr_hc(ji,jj,jk)
          END_3D
@@ -168,5 +163,5 @@
       CASE( np_BIO )                   !==  bio-model fluxes  ==!
          !
-         DO_3D( isi, iei, isj, iej, 1, nksr )
+         DO_3D_OVR( nn_hls, nn_hls, nn_hls, nn_hls, 1, nksr )
             qsr_hc(ji,jj,jk) = r1_rho0_rcp * ( etot3(ji,jj,jk) - etot3(ji,jj,jk+1) )
          END_3D
@@ -179,8 +174,8 @@
          !
          IF( nqsr == np_RGBc ) THEN          !*  Variable Chlorophyll
-            IF( ntile == 0 .OR. ntile == 1 )  THEN                                         ! Do only for the full domain
-               IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 )            ! Use full domain
+            IF( .NOT. l_istiled .OR. ntile == 1 )  THEN                                         ! Do only for the full domain
+               IF( ln_tile ) CALL dom_tile_stop( ldhold=.TRUE. )             ! Use full domain
                CALL fld_read( kt, 1, sf_chl )         ! Read Chl data and provides it at the current time step
-               IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 1 )            ! Revert to tile domain
+               IF( ln_tile ) CALL dom_tile_start( ldhold=.TRUE. )            ! Revert to tile domain
             ENDIF
             !
@@ -190,5 +185,5 @@
             ! most expensive calculations)
             !
-            DO_2D( isi, iei, isj, iej )
+            DO_2D_OVR( nn_hls, nn_hls, nn_hls, nn_hls )
                        ! zlogc = log(zchl)
                zlogc = LOG ( MIN( 10. , MAX( 0.03, sf_chl(1)%fnow(ji,jj,1) ) ) )
@@ -209,5 +204,5 @@
 
 !
-            DO_3D( isi, iei, isj, iej, 1, nksr + 1 )
+            DO_3D_OVR( nn_hls, nn_hls, nn_hls, nn_hls, 1, nksr + 1 )
                ! zchl    = ALOG( ze0(ji,jj) )
                zlogc = ze0(ji,jj)
@@ -239,5 +234,5 @@
          !
          zcoef  = ( 1. - rn_abs ) / 3._wp    !* surface equi-partition in R-G-B
-         DO_2D( isi, iei, isj, iej )
+         DO_2D_OVR( nn_hls, nn_hls, nn_hls, nn_hls )
             ze0(ji,jj) = rn_abs * qsr(ji,jj)
             ze1(ji,jj) = zcoef  * qsr(ji,jj)
@@ -250,5 +245,5 @@
          !
          !                                    !* interior equi-partition in R-G-B depending on vertical profile of Chl
-         DO_3D( isi, iei, isj, iej, 2, nksr + 1 )
+         DO_3D_OVR( nn_hls, nn_hls, nn_hls, nn_hls, 2, nksr + 1 )
             ze3t = e3t(ji,jj,jk-1,Kmm)
             irgb = NINT( ztmp3d(ji,jj,jk) )
@@ -264,5 +259,5 @@
          END_3D
          !
-         DO_3D( isi, iei, isj, iej, 1, nksr )          !* now qsr induced heat content
+         DO_3D_OVR( nn_hls, nn_hls, nn_hls, nn_hls, 1, nksr )          !* now qsr induced heat content
             qsr_hc(ji,jj,jk) = r1_rho0_rcp * ( ztmp3d(ji,jj,jk) - ztmp3d(ji,jj,jk+1) )
          END_3D
@@ -274,5 +269,5 @@
          zz0 =        rn_abs   * r1_rho0_rcp      ! surface equi-partition in 2-bands
          zz1 = ( 1. - rn_abs ) * r1_rho0_rcp
-         DO_3D( isi, iei, isj, iej, 1, nksr )          !* now qsr induced heat content
+         DO_3D_OVR( nn_hls, nn_hls, nn_hls, nn_hls, 1, nksr )          !* now qsr induced heat content
             zc0 = zz0 * EXP( -gdepw(ji,jj,jk  ,Kmm)*xsi0r ) + zz1 * EXP( -gdepw(ji,jj,jk  ,Kmm)*xsi1r )
             zc1 = zz0 * EXP( -gdepw(ji,jj,jk+1,Kmm)*xsi0r ) + zz1 * EXP( -gdepw(ji,jj,jk+1,Kmm)*xsi1r )
@@ -292,5 +287,5 @@
       !
       ! sea-ice: store the 1st ocean level attenuation coefficient
-      DO_2D( isi, iei, isj, iej )
+      DO_2D_OVR( nn_hls, nn_hls, nn_hls, nn_hls )
          IF( qsr(ji,jj) /= 0._wp ) THEN   ;   fraqsr_1lev(ji,jj) = qsr_hc(ji,jj,1) / ( r1_rho0_rcp * qsr(ji,jj) )
          ELSE                             ;   fraqsr_1lev(ji,jj) = 1._wp
@@ -298,18 +293,15 @@
       END_2D
       !
-      ! TEMP: [tiling] This change not necessary and working array can use A2D(nn_hls) if using XIOS (subdomain support)
-      IF( ntile == 0 .OR. ntile == nijtile )  THEN                ! Do only for the full domain
-         IF( iom_use('qsr3d') ) THEN      ! output the shortwave Radiation distribution
-            ALLOCATE( zetot(jpi,jpj,jpk) )
-            zetot(:,:,nksr+1:jpk) = 0._wp     ! below ~400m set to zero
-            DO jk = nksr, 1, -1
-               zetot(:,:,jk) = zetot(:,:,jk+1) + qsr_hc(:,:,jk) * rho0_rcp
-            END DO
-            CALL iom_put( 'qsr3d', zetot )   ! 3D distribution of shortwave Radiation
-            DEALLOCATE( zetot )
-         ENDIF
-      ENDIF
-      !
-      IF( ntile == 0 .OR. ntile == nijtile )  THEN                ! Do only on the last tile
+      IF( iom_use('qsr3d') ) THEN      ! output the shortwave Radiation distribution
+         ALLOCATE( zetot(A2D(nn_hls),jpk) )
+         zetot(:,:,nksr+1:jpk) = 0._wp     ! below ~400m set to zero
+         DO_3DS(0, 0, 0, 0, nksr, 1, -1)
+            zetot(ji,jj,jk) = zetot(ji,jj,jk+1) + qsr_hc(ji,jj,jk) * rho0_rcp
+         END_3D
+         CALL iom_put( 'qsr3d', zetot )   ! 3D distribution of shortwave Radiation
+         DEALLOCATE( zetot )
+      ENDIF
+      !
+      IF( .NOT. l_istiled .OR. ntile == nijtile )  THEN                ! Do only on the last tile
          IF( lrst_oce ) THEN     ! write in the ocean restart file
             CALL iom_rstput( kt, nitrst, numrow, 'qsr_hc_b'   , qsr_hc      )
