New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 14667 – NEMO

Changeset 14667


Ignore:
Timestamp:
2021-04-01T10:43:32+02:00 (3 years ago)
Author:
francesca
Message:

TRA additional changes + first set of DYN files

Location:
NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE
Files:
18 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/DOM/domqco.F90

    r14511 r14667  
    154154#if ! defined key_qcoTest_FluxForm 
    155155      !                                ! no 'key_qcoTest_FluxForm' : surface weighted ssh average 
    156          DO_2D( 0, 0, 0, 0 ) 
     156         ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) 
     157         DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
    157158            pr3u(ji,jj) = 0.5_wp * (  e1e2t(ji  ,jj) * pssh(ji  ,jj)  & 
    158159               &                    + e1e2t(ji+1,jj) * pssh(ji+1,jj)  ) * r1_hu_0(ji,jj) * r1_e1e2u(ji,jj) 
     
    162163!!st      ELSE                                         !- Flux Form   (simple averaging) 
    163164#else 
    164          DO_2D( 0, 0, 0, 0 ) 
     165         ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) 
     166         DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
    165167            pr3u(ji,jj) = 0.5_wp * (  pssh(ji,jj) + pssh(ji+1,jj  )  ) * r1_hu_0(ji,jj) 
    166168            pr3v(ji,jj) = 0.5_wp * (  pssh(ji,jj) + pssh(ji  ,jj+1)  ) * r1_hv_0(ji,jj) 
     
    170172      ! 
    171173      IF( .NOT.PRESENT( pr3f ) ) THEN              !- lbc on ratio at u-, v-points only 
    172          CALL lbc_lnk( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp ) 
     174         IF (nn_hls.eq.1) CALL lbc_lnk( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp ) 
    173175         ! 
    174176         ! 
     
    179181         !                                ! no 'key_qcoTest_FluxForm' : surface weighted ssh average 
    180182 
    181             DO_2D( 0, 0, 0, 0 )                               ! start from 1 since lbc_lnk('F') doesn't update the 1st row/line 
    182                pr3f(ji,jj) = 0.25_wp * (  e1e2t(ji  ,jj  ) * pssh(ji  ,jj  )  & 
    183                   &                     + e1e2t(ji+1,jj  ) * pssh(ji+1,jj  )  & 
    184                   &                     + e1e2t(ji  ,jj+1) * pssh(ji  ,jj+1)  & 
    185                   &                     + e1e2t(ji+1,jj+1) * pssh(ji+1,jj+1)  ) * r1_hf_0(ji,jj) * r1_e1e2f(ji,jj) 
     183            ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )                               ! start from 1 since lbc_lnk('F') doesn't update the 1st row/line 
     184            DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )                               ! start from 1 since lbc_lnk('F') doesn't update the 1st row/line 
     185               ! round brackets added to fix the order of floating point operations 
     186               ! needed to ensure halo 1 - halo 2 compatibility 
     187               pr3f(ji,jj) = 0.25_wp * ( ( e1e2t(ji  ,jj  ) * pssh(ji  ,jj  )   & 
     188                  &                      + e1e2t(ji+1,jj  ) * pssh(ji+1,jj  )   & 
     189                  &                      )                                      & ! bracket for halo 1 - halo 2 compatibility 
     190                  &                     + ( e1e2t(ji  ,jj+1) * pssh(ji  ,jj+1)  & 
     191                  &                       + e1e2t(ji+1,jj+1) * pssh(ji+1,jj+1)  & 
     192                  &                       )                                     & ! bracket for halo 1 - halo 2 compatibility 
     193                  &                    ) * r1_hf_0(ji,jj) * r1_e1e2f(ji,jj) 
    186194            END_2D 
    187195!!st         ELSE                                      !- Flux Form   (simple averaging) 
    188196#else 
    189             DO_2D( 0, 0, 0, 0 )                               ! start from 1 since lbc_lnk('F') doesn't update the 1st row/line 
    190                pr3f(ji,jj) = 0.25_wp * (  pssh(ji,jj  ) + pssh(ji+1,jj  )  & 
    191                   &                     + pssh(ji,jj+1) + pssh(ji+1,jj+1)  ) * r1_hf_0(ji,jj) 
     197            ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )                               ! start from 1 since lbc_lnk('F') doesn't update the 1st row/line 
     198            DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
     199               ! round brackets added to fix the order of floating point operations 
     200               ! needed to ensure halo 1 - halo 2 compatibility 
     201               pr3f(ji,jj) = 0.25_wp * ( ( pssh(ji,jj  ) + pssh(ji+1,jj  ) ) & 
     202                  &                     + ( pssh(ji,jj+1) + pssh(ji+1,jj+1)  &  
     203                  &                       )                                  & ! bracket for halo 1 - halo 2 compatibility 
     204                  &                    ) * r1_hf_0(ji,jj) 
    192205            END_2D 
    193206!!st         ENDIF 
    194207#endif 
    195208         !                                                 ! lbc on ratio at u-,v-,f-points 
    196          CALL lbc_lnk( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp, pr3f, 'F', 1._wp ) 
     209         IF (nn_hls.eq.1) CALL lbc_lnk( 'dom_qco_r3c', pr3u, 'U', 1._wp, pr3v, 'V', 1._wp, pr3f, 'F', 1._wp ) 
    197210         ! 
    198211      ENDIF 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/DYN/divhor.F90

    r13558 r14667  
    6464      ! 
    6565      INTEGER  ::   ji, jj, jk    ! dummy loop indices 
    66       REAL(wp) ::   zraur, zdep   ! local scalars 
    67       REAL(wp), DIMENSION(jpi,jpj) :: ztmp 
    6866      !!---------------------------------------------------------------------- 
    6967      ! 
     
    7775      ENDIF 
    7876      ! 
    79       DO_3D( 0, 0, 0, 0, 1, jpkm1 )                                    !==  Horizontal divergence  ==! 
    80          hdiv(ji,jj,jk) = (   e2u(ji  ,jj) * e3u(ji  ,jj,jk,Kmm) * uu(ji  ,jj,jk,Kmm)      & 
    81             &               - e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm) * uu(ji-1,jj,jk,Kmm)      & 
    82             &               + e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kmm) * vv(ji,jj  ,jk,Kmm)      & 
    83             &               - e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) * vv(ji,jj-1,jk,Kmm)  )   & 
    84             &            * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     77      ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 )                                    !==  Horizontal divergence  ==! 
     78      DO_3D( nn_hls-1, nn_hls, nn_hls-1, nn_hls, 1, jpkm1 )                                    !==  Horizontal divergence  ==! 
     79         ! round brackets added to fix the order of floating point operations 
     80         ! needed to ensure halo 1 - halo 2 compatibility 
     81         hdiv(ji,jj,jk) = (  ( e2u(ji  ,jj) * e3u(ji  ,jj,jk,Kmm) * uu(ji  ,jj,jk,Kmm)     & 
     82            &                - e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm) * uu(ji-1,jj,jk,Kmm)     & 
     83            &                )                                                             & ! bracket for halo 1 - halo 2 compatibility 
     84            &              + ( e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kmm) * vv(ji,jj  ,jk,Kmm)     & 
     85            &                - e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) * vv(ji,jj-1,jk,Kmm)     &  
     86            &                )                                                             & ! bracket for halo 1 - halo 2 compatibility 
     87            &             )  * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    8588      END_3D 
    8689      ! 
     
    9194      !  
    9295#endif 
    93       ! 
     96      ! WED025 + isomip true  
    9497      IF( ln_isf )                      CALL isf_hdiv( kt, Kmm, hdiv )           !==  ice shelf         ==!   (update hdiv field) 
    9598      ! 
    96       CALL lbc_lnk( 'divhor', hdiv, 'T', 1.0_wp )   !   (no sign change) 
     99      IF (nn_hls.eq.1) CALL lbc_lnk( 'divhor', hdiv, 'T', 1.0_wp )   !   (no sign change) 
    97100      ! 
    98101      IF( ln_timing )   CALL timing_stop('div_hor') 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/DYN/dynadv.F90

    r14053 r14667  
    7676         CALL dyn_zad     ( kt,                 Kmm, puu, pvv, Krhs )    ! vector form : vertical advection 
    7777      CASE( np_FLX_c2  )  
     78         ! [comm_cleanup: dyn_adv_cen2 NOT TESTED]  
    7879         CALL dyn_adv_cen2( kt,                 Kmm, puu, pvv, Krhs )    ! 2nd order centered scheme 
    7980      CASE( np_FLX_ubs )    
     81         ! [comm_cleanup: dyn_adv_ubs NOT TESTED]  
    8082         CALL dyn_adv_ubs ( kt,            Kbb, Kmm, puu, pvv, Krhs )    ! 3rd order UBS      scheme (UP3) 
    8183      END SELECT 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/DYN/dynadv_cen2.F90

    r13497 r14667  
    7272         zfu(:,:,jk) = 0.25_wp * e2u(:,:) * e3u(:,:,jk,Kmm) * puu(:,:,jk,Kmm) 
    7373         zfv(:,:,jk) = 0.25_wp * e1v(:,:) * e3v(:,:,jk,Kmm) * pvv(:,:,jk,Kmm) 
    74          DO_2D( 1, 0, 1, 0 )              ! horizontal momentum fluxes (at T- and F-point) 
     74         ! [comm_cleanup] ! DO_2D( 1, 0, 1, 0 )              ! horizontal momentum fluxes (at T- and F-point) 
     75         DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )              ! horizontal momentum fluxes (at T- and F-point) 
    7576            zfu_t(ji+1,jj  ,jk) = ( zfu(ji,jj,jk) + zfu(ji+1,jj,jk) ) * ( puu(ji,jj,jk,Kmm) + puu(ji+1,jj  ,jk,Kmm) ) 
    7677            zfv_f(ji  ,jj  ,jk) = ( zfv(ji,jj,jk) + zfv(ji+1,jj,jk) ) * ( puu(ji,jj,jk,Kmm) + puu(ji  ,jj+1,jk,Kmm) ) 
     
    7879            zfv_t(ji  ,jj+1,jk) = ( zfv(ji,jj,jk) + zfv(ji,jj+1,jk) ) * ( pvv(ji,jj,jk,Kmm) + pvv(ji  ,jj+1,jk,Kmm) ) 
    7980         END_2D 
    80          DO_2D( 0, 0, 0, 0 )              ! divergence of horizontal momentum fluxes 
     81         ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )              ! divergence of horizontal momentum fluxes 
     82         DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )              ! divergence of horizontal momentum fluxes 
    8183            puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - (  zfu_t(ji+1,jj,jk) - zfu_t(ji,jj  ,jk)    & 
    8284               &                           + zfv_f(ji  ,jj,jk) - zfv_f(ji,jj-1,jk)  ) * r1_e1e2u(ji,jj)   & 
     
    98100      !                             !==  Vertical advection  ==! 
    99101      ! 
    100       DO_2D( 0, 0, 0, 0 )                 ! surface/bottom advective fluxes set to zero 
     102      ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )                          ! surface/bottom advective fluxes set to zero 
     103      DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )                          ! surface/bottom advective fluxes set to zero 
    101104         zfu_uw(ji,jj,jpk) = 0._wp   ;   zfv_vw(ji,jj,jpk) = 0._wp 
    102105         zfu_uw(ji,jj, 1 ) = 0._wp   ;   zfv_vw(ji,jj, 1 ) = 0._wp 
    103106      END_2D 
    104107      IF( ln_linssh ) THEN                ! linear free surface: advection through the surface 
    105          DO_2D( 0, 0, 0, 0 ) 
     108         ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) 
     109         DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 
    106110            zfu_uw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * ww(ji,jj,1) + e1e2t(ji+1,jj) * ww(ji+1,jj,1) ) * puu(ji,jj,1,Kmm) 
    107111            zfv_vw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * ww(ji,jj,1) + e1e2t(ji,jj+1) * ww(ji,jj+1,1) ) * pvv(ji,jj,1,Kmm) 
     
    109113      ENDIF 
    110114      DO jk = 2, jpkm1                    ! interior advective fluxes 
    111          DO_2D( 0, 1, 0, 1 )                  ! 1/4 * Vertical transport 
     115         ! [comm_cleanup] ! DO_2D( 0, 1, 0, 1 )                  ! 1/4 * Vertical transport 
     116         DO_2D( nn_hls-1, nn_hls, nn_hls-1, nn_hls )                  ! 1/4 * Vertical transport  
    112117            zfw(ji,jj,jk) = 0.25_wp * e1e2t(ji,jj) * ww(ji,jj,jk) 
    113118         END_2D 
    114          DO_2D( 0, 0, 0, 0 ) 
     119         ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) 
     120         DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 
    115121            zfu_uw(ji,jj,jk) = ( zfw(ji,jj,jk) + zfw(ji+1,jj  ,jk) ) * ( puu(ji,jj,jk,Kmm) + puu(ji,jj,jk-1,Kmm) ) 
    116122            zfv_vw(ji,jj,jk) = ( zfw(ji,jj,jk) + zfw(ji  ,jj+1,jk) ) * ( pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk-1,Kmm) ) 
    117123         END_2D 
    118124      END DO 
    119       DO_3D( 0, 0, 0, 0, 1, jpkm1 )       ! divergence of vertical momentum flux divergence 
     125      ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 )       ! divergence of vertical momentum flux divergence 
     126      DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 )       ! divergence of vertical momentum flux divergence 
    120127         puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj)   & 
    121128            &                                      / e3u(ji,jj,jk,Kmm) 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/DYN/dynadv_ubs.F90

    r14511 r14667  
    108108         zfv(:,:,jk) = e1v(:,:) * e3v(:,:,jk,Kmm) * pvv(:,:,jk,Kmm) 
    109109         !             
    110          DO_2D( 0, 0, 0, 0 )                       ! laplacian 
     110         ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )                       ! laplacian 
     111         DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )                       ! laplacian 
    111112            zlu_uu(ji,jj,jk,1) = ( puu (ji+1,jj  ,jk,Kbb) - 2.*puu (ji,jj,jk,Kbb) + puu (ji-1,jj  ,jk,Kbb) ) * umask(ji,jj,jk) 
    112113            zlv_vv(ji,jj,jk,1) = ( pvv (ji  ,jj+1,jk,Kbb) - 2.*pvv (ji,jj,jk,Kbb) + pvv (ji  ,jj-1,jk,Kbb) ) * vmask(ji,jj,jk) 
     
    124125         END_2D 
    125126      END DO 
    126       CALL lbc_lnk( 'dynadv_ubs', zlu_uu(:,:,:,1), 'U', 1.0_wp , zlu_uv(:,:,:,1), 'U', 1.0_wp,  & 
    127          &                        zlu_uu(:,:,:,2), 'U', 1.0_wp , zlu_uv(:,:,:,2), 'U', 1.0_wp,  &  
    128          &                        zlv_vv(:,:,:,1), 'V', 1.0_wp , zlv_vu(:,:,:,1), 'V', 1.0_wp,  & 
    129          &                        zlv_vv(:,:,:,2), 'V', 1.0_wp , zlv_vu(:,:,:,2), 'V', 1.0_wp   ) 
     127      IF (nn_hls.eq.1) CALL lbc_lnk( 'dynadv_ubs', zlu_uu(:,:,:,1), 'U', 1.0_wp , zlu_uv(:,:,:,1), 'U', 1.0_wp,  & 
     128                          &                        zlu_uu(:,:,:,2), 'U', 1.0_wp , zlu_uv(:,:,:,2), 'U', 1.0_wp,  &  
     129                          &                        zlv_vv(:,:,:,1), 'V', 1.0_wp , zlv_vu(:,:,:,1), 'V', 1.0_wp,  & 
     130                          &                        zlv_vv(:,:,:,2), 'V', 1.0_wp , zlv_vu(:,:,:,2), 'V', 1.0_wp   ) 
    130131      ! 
    131132      !                                      ! ====================== ! 
     
    136137         zfv(:,:,jk) = 0.25_wp * e1v(:,:) * e3v(:,:,jk,Kmm) * pvv(:,:,jk,Kmm) 
    137138         ! 
    138          DO_2D( 1, 0, 1, 0 )                       ! horizontal momentum fluxes at T- and F-point 
     139         ! [comm_cleanup] ! DO_2D( 1, 0, 1, 0 )              ! horizontal momentum fluxes at T- and F-point 
     140         DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )              ! horizontal momentum fluxes at T- and F-point 
    139141            zui = ( puu(ji,jj,jk,Kmm) + puu(ji+1,jj  ,jk,Kmm) ) 
    140142            zvj = ( pvv(ji,jj,jk,Kmm) + pvv(ji  ,jj+1,jk,Kmm) ) 
     
    168170               &                * ( pvv(ji,jj,jk,Kmm) + pvv(ji+1,jj  ,jk,Kmm) - gamma1 * zl_v ) 
    169171         END_2D 
    170          DO_2D( 0, 0, 0, 0 )                       ! divergence of horizontal momentum fluxes 
     172         ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )              ! divergence of horizontal momentum fluxes 
     173         DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )              ! divergence of horizontal momentum fluxes 
    171174            puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - (  zfu_t(ji+1,jj,jk) - zfu_t(ji,jj  ,jk)    & 
    172175               &                           + zfv_f(ji  ,jj,jk) - zfv_f(ji,jj-1,jk)  ) * r1_e1e2u(ji,jj)   & 
     
    187190      !                                      !  Vertical advection  ! 
    188191      !                                      ! ==================== ! 
    189       DO_2D( 0, 0, 0, 0 )                          ! surface/bottom advective fluxes set to zero 
     192      ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )                          ! surface/bottom advective fluxes set to zero 
     193      DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )                          ! surface/bottom advective fluxes set to zero 
    190194         zfu_uw(ji,jj,jpk) = 0._wp 
    191195         zfv_vw(ji,jj,jpk) = 0._wp 
     
    194198      END_2D 
    195199      IF( ln_linssh ) THEN                         ! constant volume : advection through the surface 
    196          DO_2D( 0, 0, 0, 0 ) 
     200         ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) 
     201         DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 
    197202            zfu_uw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * ww(ji,jj,1) + e1e2t(ji+1,jj) * ww(ji+1,jj,1) ) * puu(ji,jj,1,Kmm) 
    198203            zfv_vw(ji,jj,1) = 0.5_wp * ( e1e2t(ji,jj) * ww(ji,jj,1) + e1e2t(ji,jj+1) * ww(ji,jj+1,1) ) * pvv(ji,jj,1,Kmm) 
     
    200205      ENDIF 
    201206      DO jk = 2, jpkm1                          ! interior fluxes 
    202          DO_2D( 0, 1, 0, 1 ) 
     207         ! [comm_cleanup] ! DO_2D( 0, 1, 0, 1 ) 
     208         DO_2D( nn_hls-1, nn_hls, nn_hls-1, nn_hls )  
    203209            zfw(ji,jj,jk) = 0.25_wp * e1e2t(ji,jj) * ww(ji,jj,jk) 
    204210         END_2D 
    205          DO_2D( 0, 0, 0, 0 ) 
     211         ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) 
     212         DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 
    206213            zfu_uw(ji,jj,jk) = ( zfw(ji,jj,jk)+ zfw(ji+1,jj,jk) ) * ( puu(ji,jj,jk,Kmm) + puu(ji,jj,jk-1,Kmm) ) 
    207214            zfv_vw(ji,jj,jk) = ( zfw(ji,jj,jk)+ zfw(ji,jj+1,jk) ) * ( pvv(ji,jj,jk,Kmm) + pvv(ji,jj,jk-1,Kmm) ) 
    208215         END_2D 
    209216      END DO 
    210       DO_3D( 0, 0, 0, 0, 1, jpkm1 )             ! divergence of vertical momentum flux divergence 
     217      ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 )       ! divergence of vertical momentum flux divergence 
     218      DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 )       ! divergence of vertical momentum flux divergence 
    211219         puu(ji,jj,jk,Krhs) =  puu(ji,jj,jk,Krhs) - ( zfu_uw(ji,jj,jk) - zfu_uw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj)   & 
    212220            &                                       / e3u(ji,jj,jk,Kmm) 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/DYN/dynkeg.F90

    r13497 r14667  
    101101      ! 
    102102      CASE ( nkeg_C2 )                          !--  Standard scheme  --! 
    103          DO_3D( 0, 1, 0, 1, 1, jpkm1 ) 
     103         ! [comm_cleanup] ! DO_3D( 0, 1, 0, 1, 1, jpkm1 ) 
     104         DO_3D( nn_hls-1, nn_hls, nn_hls-1, nn_hls, 1, jpkm1 ) 
    104105            zu =    puu(ji-1,jj  ,jk,Kmm) * puu(ji-1,jj  ,jk,Kmm)   & 
    105106               &  + puu(ji  ,jj  ,jk,Kmm) * puu(ji  ,jj  ,jk,Kmm) 
     
    109110         END_3D 
    110111      CASE ( nkeg_HW )                          !--  Hollingsworth scheme  --! 
    111          DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     112         ! [comm_cleanup: Hollingsworth scheme NOT TESTED] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     113         DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 ) 
    112114            zu = 8._wp * ( puu(ji-1,jj  ,jk,Kmm) * puu(ji-1,jj  ,jk,Kmm)    & 
    113115               &         + puu(ji  ,jj  ,jk,Kmm) * puu(ji  ,jj  ,jk,Kmm) )  & 
     
    121123            zhke(ji,jj,jk) = r1_48 * ( zv + zu ) 
    122124         END_3D 
    123          CALL lbc_lnk( 'dynkeg', zhke, 'T', 1.0_wp ) 
     125         IF (nn_hls.eq.1) CALL lbc_lnk( 'dynkeg', zhke, 'T', 1.0_wp ) 
    124126         ! 
    125127      END SELECT  
    126128      ! 
    127       DO_3D( 0, 0, 0, 0, 1, jpkm1 )       !==  grad( KE ) added to the general momentum trends  ==! 
     129      ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 )       !==  grad( KE ) added to the general momentum trends  ==! 
     130      DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 )       !==  grad( KE ) added to the general momentum trends  ==! 
    128131         puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zhke(ji+1,jj  ,jk) - zhke(ji,jj,jk) ) / e1u(ji,jj) 
    129132         pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) - ( zhke(ji  ,jj+1,jk) - zhke(ji,jj,jk) ) / e2v(ji,jj) 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/DYN/dynzad.F90

    r14072 r14667  
    7979 
    8080      DO jk = 2, jpkm1                ! Vertical momentum advection at level w and u- and v- vertical 
    81          DO_2D( 0, 1, 0, 1 )              ! vertical fluxes 
     81         ! [comm_cleanup] ! DO_2D( 0, 1, 0, 1 )              ! vertical fluxes  
     82         DO_2D( nn_hls-1, nn_hls, nn_hls-1, nn_hls )              ! vertical fluxes 
    8283          IF( ln_vortex_force ) THEN 
    8384            zww(ji,jj) = 0.25_wp * e1e2t(ji,jj) * ( ww(ji,jj,jk) + wsd(ji,jj,jk) ) 
     
    8687          ENDIF 
    8788         END_2D 
    88          DO_2D( 0, 0, 0, 0 )              ! vertical momentum advection at w-point 
     89         ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )              ! vertical momentum advection at w-point 
     90         DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )              ! vertical momentum advection at w-point 
    8991            zwuw(ji,jj,jk) = ( zww(ji+1,jj  ) + zww(ji,jj) ) * ( puu(ji,jj,jk-1,Kmm) - puu(ji,jj,jk,Kmm) ) 
    9092            zwvw(ji,jj,jk) = ( zww(ji  ,jj+1) + zww(ji,jj) ) * ( pvv(ji,jj,jk-1,Kmm) - pvv(ji,jj,jk,Kmm) ) 
     
    9395      ! 
    9496      ! Surface and bottom advective fluxes set to zero 
    95       DO_2D( 0, 0, 0, 0 ) 
     97      ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )  
     98      DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 
    9699         zwuw(ji,jj, 1 ) = 0._wp 
    97100         zwvw(ji,jj, 1 ) = 0._wp 
     
    100103      END_2D 
    101104      ! 
    102       DO_3D( 0, 0, 0, 0, 1, jpkm1 )   ! Vertical momentum advection at u- and v-points 
     105      ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 )   ! Vertical momentum advection at u- and v-points 
     106      DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 )   ! Vertical momentum advection at u- and v-points 
    103107         puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) - ( zwuw(ji,jj,jk) + zwuw(ji,jj,jk+1) ) * r1_e1e2u(ji,jj)   & 
    104108            &                                      / e3u(ji,jj,jk,Kmm) 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/DYN/sshwzv.F90

    r14205 r14667  
    7878      REAL(wp), DIMENSION(jpi,jpj,jpt), INTENT(inout) ::   pssh           ! sea-surface height 
    7979      !  
    80       INTEGER  ::   jk      ! dummy loop index 
     80      INTEGER  ::   ji, jj, jk      ! dummy loop index 
    8181      REAL(wp) ::   zcoef   ! local scalar 
    8282      REAL(wp), DIMENSION(jpi,jpj) ::   zhdiv   ! 2D workspace 
     
    103103      ! 
    104104      zhdiv(:,:) = 0._wp 
    105       DO jk = 1, jpkm1                                 ! Horizontal divergence of barotropic transports 
    106         zhdiv(:,:) = zhdiv(:,:) + e3t(:,:,jk,Kmm) * hdiv(:,:,jk) 
    107       END DO 
     105      ! [comm_cleanup] ! DO jk = 1, jpkm1                                 ! Horizontal divergence of barotropic transports  
     106      DO_3D( nn_hls-1, nn_hls, nn_hls-1, nn_hls, 1, jpkm1 )                                 ! Horizontal divergence of barotropic transports  
     107        zhdiv(ji,jj) = zhdiv(ji,jj) + e3t(ji,jj,jk,Kmm) * hdiv(ji,jj,jk) 
     108      END_3D 
    108109      !                                                ! Sea surface elevation time stepping 
    109110      ! In time-split case we need a first guess of the ssh after (using the baroclinic timestep) in order to 
    110111      ! compute the vertical velocity which can be used to compute the non-linear terms of the momentum equations. 
    111112      !  
    112       pssh(:,:,Kaa) = (  pssh(:,:,Kbb) - rDt * ( zcoef * ( emp_b(:,:) + emp(:,:) ) + zhdiv(:,:) )  ) * ssmask(:,:) 
     113      ! [comm_cleanup]  
     114      DO_2D( nn_hls-1, nn_hls, nn_hls-1, nn_hls ) 
     115         pssh(ji,jj,Kaa) = (  pssh(ji,jj,Kbb) - rDt * ( zcoef * ( emp_b(ji,jj) + emp(ji,jj) ) + zhdiv(ji,jj) )  ) * ssmask(ji,jj) 
     116      END_2D 
    113117      ! 
    114118#if defined key_agrif 
     
    119123      IF ( .NOT.ln_dynspg_ts ) THEN 
    120124         IF( ln_bdy ) THEN 
    121             CALL lbc_lnk( 'sshwzv', pssh(:,:,Kaa), 'T', 1.0_wp )    ! Not sure that's necessary 
     125            ! [comm_cleanup]  
     126            IF (nn_hls.eq.1) CALL lbc_lnk( 'sshwzv', pssh(:,:,Kaa), 'T', 1.0_wp )    ! Not sure that's necessary 
    122127            CALL bdy_ssh( pssh(:,:,Kaa) )             ! Duplicate sea level across open boundaries 
    123128         ENDIF 
     
    178183            ! horizontal divergence of thickness diffusion transport ( velocity multiplied by e3t) 
    179184            ! - ML - note: computation already done in dom_vvl_sf_nxt. Could be optimized (not critical and clearer this way) 
    180             DO_2D( 0, 0, 0, 0 ) 
     185            ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) 
     186            DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 
    181187               zhdiv(ji,jj,jk) = r1_e1e2t(ji,jj) * ( un_td(ji,jj,jk) - un_td(ji-1,jj,jk) + vn_td(ji,jj,jk) - vn_td(ji,jj-1,jk) ) 
    182188            END_2D 
    183189         END DO 
    184          CALL lbc_lnk('sshwzv', zhdiv, 'T', 1.0_wp)  ! - ML - Perhaps not necessary: not used for horizontal "connexions" 
     190         IF (nn_hls.eq.1) CALL lbc_lnk('sshwzv', zhdiv, 'T', 1.0_wp)  ! - ML - Perhaps not necessary: not used for horizontal "connexions" 
    185191         !                             ! Is it problematic to have a wrong vertical velocity in boundary cells? 
    186192         !                             ! Same question holds for hdiv. Perhaps just for security 
     
    357363      zdt = 2._wp * rn_Dt                            ! 2*rn_Dt and not rDt (for restartability) 
    358364      IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN 
    359          DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     365         ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     366         DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 ) 
    360367            z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm) 
    361368            Cu_adv(ji,jj,jk) =   zdt *                                                         & 
     
    374381         END_3D 
    375382      ELSE 
    376          DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     383         ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     384         DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 ) 
    377385            z1_e3t = 1._wp / e3t(ji,jj,jk,Kmm) 
    378386            Cu_adv(ji,jj,jk) =   zdt *                                                      & 
     
    387395         END_3D 
    388396      ENDIF 
    389       CALL lbc_lnk( 'sshwzv', Cu_adv, 'T', 1.0_wp ) 
     397      IF (nn_hls.eq.1) CALL lbc_lnk( 'sshwzv', Cu_adv, 'T', 1.0_wp ) 
    390398      ! 
    391399      CALL iom_put("Courant",Cu_adv) 
    392400      ! 
    393401      IF( MAXVAL( Cu_adv(:,:,:) ) > Cu_min ) THEN       ! Quick check if any breaches anywhere 
    394          DO_3DS( 1, 1, 1, 1, jpkm1, 2, -1 )             ! or scan Courant criterion and partition ! w where necessary 
     402         ! [comm_cleanup] ! DO_3DS( 1, 1, 1, 1, jpkm1, 2, -1 )             ! or scan Courant criterion and partition ! w where necessary 
     403         DO_3DS( nn_hls, nn_hls, nn_hls, nn_hls, jpkm1, 2, -1 )             ! or scan Courant criterion and partition ! w where necessary 
    395404            ! 
    396405            zCu = MAX( Cu_adv(ji,jj,jk) , Cu_adv(ji,jj,jk-1) ) 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/ISF/isfhdiv.F90

    r13295 r14667  
    100100      ! 
    101101      ! update divergence at each level affected by ice shelf top boundary layer 
    102       DO_2D( 1, 1, 1, 1 ) 
     102      ! [comm_cleanup] ! DO_2D( 1, 1, 1, 1 ) 
     103      DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )  
    103104         ikt = ktop(ji,jj) 
    104105         ikb = kbot(ji,jj) 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/LDF/ldfc1d_c2d.F90

    r14511 r14667  
    135135      ! 
    136136      CASE( 'DYN' )                       ! T- and F-points 
    137          DO_2D( 1, 1, 1, 1 ) 
     137         DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) 
    138138            pah1(ji,jj,1) = pUfac * MAX( e1t(ji,jj) , e2t(ji,jj) )**knn 
    139139            pah2(ji,jj,1) = pUfac * MAX( e1f(ji,jj) , e2f(ji,jj) )**knn 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/LDF/ldftra.F90

    r14538 r14667  
    633633      INTEGER                         , INTENT(in   ) ::   kt             ! ocean time-step index 
    634634      INTEGER                         , INTENT(in   ) ::   Kmm            ! ocean time level indices 
    635       REAL(wp)                        , INTENT(inout) ::   paei0          ! max value            [m2/s] 
     635      REAL(wp)                        , INTENT(in   ) ::   paei0          ! max value            [m2/s] 
    636636      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) ::   paeiu, paeiv   ! eiv coefficient      [m2/s] 
    637637      ! 
    638638      INTEGER  ::   ji, jj, jk    ! dummy loop indices 
    639       REAL(wp) ::   zfw, ze3w, zn2, z1_f20, zaht, zaht_min, zzaei    ! local scalars 
     639      REAL(wp) ::   zfw, ze3w, zn2, z1_f20, zzaei    ! local scalars 
    640640      REAL(wp), DIMENSION(jpi,jpj) ::   zn, zah, zhw, zRo, zaeiw   ! 2D workspace 
    641641      !!---------------------------------------------------------------------- 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/SBC/sbcrnf.F90

    r14072 r14667  
    206206      IF( ln_rnf_depth .OR. ln_rnf_depth_ini ) THEN      !==   runoff distributed over several levels   ==! 
    207207         IF( ln_linssh ) THEN    !* constant volume case : just apply the runoff input flow 
    208             DO_2D( 1, 1, 1, 1 ) 
     208            ! [comm_cleanup] ! DO_2D( 1, 1, 1, 1 ) 
     209            DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) 
    209210               DO jk = 1, nk_rnf(ji,jj) 
    210211                  phdivn(ji,jj,jk) = phdivn(ji,jj,jk) - ( rnf(ji,jj) + rnf_b(ji,jj) ) * zfact * r1_rho0 / h_rnf(ji,jj) 
     
    212213            END_2D 
    213214         ELSE                    !* variable volume case 
    214             DO_2D( 1, 1, 1, 1 )              ! update the depth over which runoffs are distributed 
     215            ! [comm_cleanup] ! DO_2D( 1, 1, 1, 1 )              ! update the depth over which runoffs are distributed 
     216            DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )              ! update the depth over which runoffs are distributed 
    215217               h_rnf(ji,jj) = 0._wp 
    216218               DO jk = 1, nk_rnf(ji,jj)                             ! recalculates h_rnf to be the depth in metres 
     
    358360         ! 
    359361         nk_rnf(:,:) = 0                               ! set the number of level over which river runoffs are applied 
    360          DO_2D( 1, 1, 1, 1 ) 
     362         ! [comm_cleanup] ! DO_2D( 1, 1, 1, 1 ) 
     363         DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) 
    361364            IF( h_rnf(ji,jj) > 0._wp ) THEN 
    362365               jk = 2 
     
    371374            ENDIF 
    372375         END_2D 
    373          DO_2D( 1, 1, 1, 1 )                           ! set the associated depth 
     376         ! [comm_cleanup] ! DO_2D( 1, 1, 1, 1 )                           ! set the associated depth 
     377         DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )                           ! set the associated depth 
    374378            h_rnf(ji,jj) = 0._wp 
    375379            DO jk = 1, nk_rnf(ji,jj) 
     
    401405         WHERE( zrnfcl(:,:,1) > 0._wp )  h_rnf(:,:) = zacoef * zrnfcl(:,:,1)   ! compute depth for all runoffs 
    402406         ! 
    403          DO_2D( 1, 1, 1, 1 )                ! take in account min depth of ocean rn_hmin 
     407         ! [comm_cleanup] ! DO_2D( 1, 1, 1, 1 )                ! take in account min depth of ocean rn_hmin 
     408         DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )                ! take in account min depth of ocean rn_hmin 
    404409            IF( zrnfcl(ji,jj,1) > 0._wp ) THEN 
    405410               jk = mbkt(ji,jj) 
     
    409414         ! 
    410415         nk_rnf(:,:) = 0                       ! number of levels on which runoffs are distributed 
    411          DO_2D( 1, 1, 1, 1 ) 
     416         ! [comm_cleanup] ! DO_2D( 1, 1, 1, 1 ) 
     417         DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) 
    412418            IF( zrnfcl(ji,jj,1) > 0._wp ) THEN 
    413419               jk = 2 
     
    420426         END_2D 
    421427         ! 
    422          DO_2D( 1, 1, 1, 1 )                          ! set the associated depth 
     428         ! [comm_cleanup] ! DO_2D( 1, 1, 1, 1 )                          ! set the associated depth 
     429         DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )                          ! set the associated depth 
    423430            h_rnf(ji,jj) = 0._wp 
    424431            DO jk = 1, nk_rnf(ji,jj) 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/TRA/traldf.F90

    r14538 r14667  
    9292            CALL tra_ldf_triad( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs), jpts,  1 ) 
    9393         CASE ( np_blp , np_blp_i , np_blp_it )             ! bilaplacian: iso-level & iso-neutral operators 
    94             ! [comm_cleanup] 
    95             ! IF (nn_hls.EQ.2) CALL lbc_lnk( 'tra_ldf', pts(:,:,:,:,Kbb), 'T',1.) 
    9694            CALL tra_ldf_blp  ( kt, Kmm, nit000,'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, pts(:,:,:,:,Kbb), pts(:,:,:,:,Krhs),             jpts, nldf_tra ) 
    9795         END SELECT 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/TRA/traldf_iso.F90

    r14538 r14667  
    180180            zmskv = wmask(ji,jj,jk) / MAX(   vmask(ji,jj  ,jk-1) + vmask(ji,jj-1,jk)          & 
    181181               &                           + vmask(ji,jj-1,jk-1) + vmask(ji,jj  ,jk) , 1._wp  ) 
    182                ! 
    183             zahu_w = (   pahu(ji  ,jj,jk-1) + pahu(ji-1,jj,jk)    & 
    184                &       + pahu(ji-1,jj,jk-1) + pahu(ji  ,jj,jk)  ) * zmsku 
    185             zahv_w = (   pahv(ji,jj  ,jk-1) + pahv(ji,jj-1,jk)    & 
    186                &       + pahv(ji,jj-1,jk-1) + pahv(ji,jj  ,jk)  ) * zmskv 
     182            ! round brackets added to fix the order of floating point operations 
     183            ! needed to ensure halo 1 - halo 2 compatibility 
     184            zahu_w = ( (  pahu(ji  ,jj,jk-1) + pahu(ji-1,jj,jk)                    & 
     185               &       )                                                           & ! bracket for halo 1 - halo 2 compatibility 
     186               &       + ( pahu(ji-1,jj,jk-1) + pahu(ji  ,jj,jk)                   & 
     187               &         ) ) * zmsku                                                 ! bracket for halo 1 - halo 2 compatibility 
     188            zahv_w = ( (  pahv(ji,jj  ,jk-1) + pahv(ji,jj-1,jk)                    & 
     189               &       )                                                           & ! bracket for halo 1 - halo 2 compatibility 
     190               &       + ( pahv(ji,jj-1,jk-1) + pahv(ji,jj  ,jk)                   &  
     191               &         ) ) * zmskv                                                 ! bracket for halo 1 - halo 2 compatibility 
    187192               ! 
    188193            ah_wslp2(ji,jj,jk) = zahu_w * wslpi(ji,jj,jk) * wslpi(ji,jj,jk)   & 
     
    193198            ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    194199            DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 )  
     200               ! round brackets added to fix the order of floating point operations 
     201               ! needed to ensure halo 1 - halo 2 compatibility 
    195202               akz(ji,jj,jk) = 0.25_wp * (                                                                     & 
    196                   &              ( pahu(ji  ,jj,jk) + pahu(ji  ,jj,jk-1) ) / ( e1u(ji  ,jj) * e1u(ji  ,jj) )   & 
     203                  &            ( ( pahu(ji  ,jj,jk) + pahu(ji  ,jj,jk-1) ) / ( e1u(ji  ,jj) * e1u(ji  ,jj) )   & 
    197204                  &            + ( pahu(ji-1,jj,jk) + pahu(ji-1,jj,jk-1) ) / ( e1u(ji-1,jj) * e1u(ji-1,jj) )   & 
    198                   &            + ( pahv(ji,jj  ,jk) + pahv(ji,jj  ,jk-1) ) / ( e2v(ji,jj  ) * e2v(ji,jj  ) )   & 
    199                   &            + ( pahv(ji,jj-1,jk) + pahv(ji,jj-1,jk-1) ) / ( e2v(ji,jj-1) * e2v(ji,jj-1) )   ) 
     205                  &            )                                                                               & ! bracket for halo 1 - halo 2 compatibility 
     206                  &            + ( ( pahv(ji,jj  ,jk) + pahv(ji,jj  ,jk-1) ) / ( e2v(ji,jj  ) * e2v(ji,jj  ) ) & 
     207                  &              + ( pahv(ji,jj-1,jk) + pahv(ji,jj-1,jk-1) ) / ( e2v(ji,jj-1) * e2v(ji,jj-1) ) &   
     208                  &              ) )                                                                             ! bracket for halo 1 - halo 2 compatibility 
    200209            END_3D 
    201210            ! 
     
    289298               zcof2 = - pahv(ji,jj,jk) * e1v(ji,jj) * vslp(ji,jj,jk) * zmskv 
    290299               ! 
    291                zftu(ji,jj,jk ) = (  zabe1 * zdit(ji,jj,jk)   & 
    292                   &               + zcof1 * (  zdkt (ji+1,jj) + zdk1t(ji,jj)      & 
    293                   &                          + zdk1t(ji+1,jj) + zdkt (ji,jj)  )  ) * umask(ji,jj,jk) 
    294                zftv(ji,jj,jk) = (  zabe2 * zdjt(ji,jj,jk)   & 
    295                   &               + zcof2 * (  zdkt (ji,jj+1) + zdk1t(ji,jj)      & 
    296                   &                          + zdk1t(ji,jj+1) + zdkt (ji,jj)  )  ) * vmask(ji,jj,jk) 
     300               ! round brackets added to fix the order of floating point operations 
     301               ! needed to ensure halo 1 - halo 2 compatibility 
     302               zftu(ji,jj,jk ) = (  zabe1 * zdit(ji,jj,jk)                       & 
     303                  &               + zcof1 * ( ( zdkt (ji+1,jj) + zdk1t(ji,jj)    & 
     304                  &                           )                                  & ! bracket for halo 1 - halo 2 compatibility 
     305                  &                         + ( zdk1t(ji+1,jj) + zdkt (ji,jj)    &  
     306                  &                           )                                  & ! bracket for halo 1 - halo 2 compatibility 
     307                  &                         ) ) * umask(ji,jj,jk) 
     308               zftv(ji,jj,jk) = (  zabe2 * zdjt(ji,jj,jk)                        & 
     309                  &              + zcof2 * ( ( zdkt (ji,jj+1) + zdk1t(ji,jj)     & 
     310                  &                           )                                  & ! bracket for halo 1 - halo 2 compatibility 
     311                  &                         + ( zdk1t(ji,jj+1) + zdkt (ji,jj)    & 
     312                  &                           )                                  & ! bracket for halo 1 - halo 2 compatibility 
     313                  &                         ) ) * vmask(ji,jj,jk)       
    297314            END_2D 
    298315            ! 
    299316            ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 )           !== horizontal divergence and add to pta 
    300317            DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )           !== horizontal divergence and add to pta 
    301                pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)    & 
    302                   &       + zsign * (  zftu(ji,jj,jk) - zftu(ji-1,jj,jk) + zftv(ji,jj,jk) - zftv(ji,jj-1,jk)  )   & 
    303                   &                                              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     318               ! round brackets added to fix the order of floating point operations 
     319               ! needed to ensure halo 1 - halo 2 compatibility 
     320               pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)                         & 
     321                  &       + zsign * ( ( zftu(ji,jj,jk) - zftu(ji-1,jj,jk)        & 
     322                  &                   )                                          & ! bracket for halo 1 - halo 2 compatibility 
     323                  &                 + ( zftv(ji,jj,jk) - zftv(ji,jj-1,jk)        & 
     324                  &                   )                                          & ! bracket for halo 1 - halo 2 compatibility 
     325                  &                 ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    304326            END_2D 
    305327         END DO                                        !   End of slab 
     
    330352            zcoef4 = - zahv_w * e1t(ji,jj) * zmskv * wslpj (ji,jj,jk) 
    331353            ! 
    332             ztfw(ji,jj,jk) = zcoef3 * (   zdit(ji  ,jj  ,jk-1) + zdit(ji-1,jj  ,jk)      & 
    333                &                        + zdit(ji-1,jj  ,jk-1) + zdit(ji  ,jj  ,jk)  )   & 
    334                &           + zcoef4 * (   zdjt(ji  ,jj  ,jk-1) + zdjt(ji  ,jj-1,jk)      & 
    335                &                        + zdjt(ji  ,jj-1,jk-1) + zdjt(ji  ,jj  ,jk)  ) 
     354            ! round brackets added to fix the order of floating point operations 
     355            ! needed to ensure halo 1 - halo 2 compatibility 
     356            ztfw(ji,jj,jk) = zcoef3 * ( ( zdit(ji  ,jj  ,jk-1) + zdit(ji-1,jj  ,jk)    & 
     357                  &                     )                                              & ! bracket for halo 1 - halo 2 compatibility 
     358                  &                   + ( zdit(ji-1,jj  ,jk-1) + zdit(ji  ,jj  ,jk)    & 
     359                  &                     )                                              & ! bracket for halo 1 - halo 2 compatibility  
     360                  &                   )                                                & 
     361                  &        + zcoef4 * ( ( zdjt(ji  ,jj  ,jk-1) + zdjt(ji  ,jj-1,jk)    & 
     362                  &                     )                                              & ! bracket for halo 1 - halo 2 compatibility 
     363                  &                   + ( zdjt(ji  ,jj-1,jk-1) + zdjt(ji  ,jj  ,jk)    & 
     364                  &                     )                                              & ! bracket for halo 1 - halo 2 compatibility 
     365                  &                   ) 
    336366         END_3D 
    337367         !                                !==  add the vertical 33 flux  ==! 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/TRA/traldf_lap_blp.F90

    r14609 r14667  
    159159         ! 
    160160         DO_3D( isi, iei, isj, iej, 1, jpkm1 )            !== Second derivative (divergence) added to the general tracer trends  ==! 
    161             pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + (  ztu(ji,jj,jk) - ztu(ji-1,jj,jk)     & 
    162                &                                      +    ztv(ji,jj,jk) - ztv(ji,jj-1,jk) )   & 
    163                &                                      / ( e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) ) 
     161            ! round brackets added to fix the order of floating point operations 
     162            ! needed to ensure halo 1 - halo 2 compatibility 
     163            pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + ( ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk)    & 
     164               &                                          )                                    & ! bracket for halo 1 - halo 2 compatibility 
     165               &                                      +   ( ztv(ji,jj,jk) - ztv(ji,jj-1,jk)    & 
     166               &                                          )                                    & ! bracket for halo 1 - halo 2 compatibility 
     167               &                                        ) / ( e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm) ) 
    164168         END_3D 
    165169         ! 
     
    239243      !                                               ! Partial top/bottom cell: GRADh( zlap ) 
    240244      IF( ln_zps ) THEN 
    241          ! [comm_cleanup] ! halo > 2 needed to delete this comm 
    242          IF (nn_hls.EQ.2) CALL lbc_lnk( 'traldf_lap_blp', zlap(:,:,:,:) , 'T', 1.0_wp )     ! Lateral boundary conditions (unchanged sign) 
    243245         IF( ln_isfcav ) THEN   ;   CALL zps_hde_isf( kt, Kmm, kjpt, zlap, zglu, zglv, zgui, zgvi )  ! both top & bottom 
    244246         ELSE   ;   CALL zps_hde    ( kt, Kmm, kjpt, zlap, zglu, zglv )              ! only bottom 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/TRA/traldf_triad.F90

    r14538 r14667  
    109109      REAL(wp), DIMENSION(A2D_T(ktt_rhs),JPK,KJPT), INTENT(inout) ::   pt_rhs     ! tracer trend 
    110110      ! 
    111       INTEGER  ::  ji, jj, jk, jn   ! dummy loop indices 
    112       INTEGER  ::  ip,jp,kp         ! dummy loop indices 
    113       INTEGER  ::  ierr            ! local integer 
    114       REAL(wp) ::  zmsku, zabe1, zcof1, zcoef3    ! local scalars 
    115       REAL(wp) ::  zmskv, zabe2, zcof2, zcoef4    !   -      - 
     111      INTEGER  ::  ji, jj, jk, jn, kp   ! dummy loop indices 
    116112      REAL(wp) ::  zcoef0, ze3w_2, zsign          !   -      - 
    117113      ! 
    118       REAL(wp) ::   zslope_skew, zslope_iso, zslope2, zbu, zbv 
    119       REAL(wp) ::   ze1ur, ze2vr, ze3wr, zdxt, zdyt, zdzt 
    120       REAL(wp) ::   zah, zah_slp, zaei_slp 
     114      REAL(wp) ::   zslope2, zbu, zbv, zbu1, zbv1, zslope21, zah, zah1, zah_ip1, zah_jp1, zbu_ip1, zbv_jp1 
     115      REAL(wp) ::   ze1ur, ze2vr, ze3wr, zdxt, zdyt, zdzt, zdyt_jp1, ze3wr_jp1, zdzt_jp1, zah_slp1, zah_slp_jp1, zaei_slp_jp1 
     116      REAL(wp) ::   zah_slp, zaei_slp, zdxt_ip1, ze3wr_ip1, zdzt_ip1, zah_slp_ip1, zaei_slp_ip1, zaei_slp1 
    121117      REAL(wp), DIMENSION(A2D(nn_hls),0:1)     ::   zdkt3d                         ! vertical tracer gradient at 2 levels 
    122118      REAL(wp), DIMENSION(A2D(nn_hls)        ) ::   z2d                            ! 2D workspace 
     
    153149         ! 
    154150         ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpk ) 
    155          DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpk )  
     151         DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpk ) 
    156152            akz     (ji,jj,jk) = 0._wp 
    157153            ah_wslp2(ji,jj,jk) = 0._wp 
    158154         END_3D 
    159155         ! 
    160          DO ip = 0, 1                            ! i-k triads 
    161             DO kp = 0, 1 
    162                ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
    163                DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 )  
    164                   ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
    165                   zbu   = e1e2u(ji-ip,jj) * e3u(ji-ip,jj,jk,Kmm) 
    166                   zah   = 0.25_wp * pahu(ji-ip,jj,jk) 
    167                   zslope_skew = triadi_g(ji,jj,jk,1-ip,kp) 
    168                   ! Subtract s-coordinate slope at t-points to give slope rel to s-surfaces (do this by *adding* gradient of depth) 
    169                   zslope2 = zslope_skew + ( gdept(ji-ip+1,jj,jk,Kmm) - gdept(ji-ip,jj,jk,Kmm) ) * r1_e1u(ji-ip,jj) * umask(ji-ip,jj,jk+kp) 
    170                   zslope2 = zslope2 *zslope2 
    171                   ah_wslp2(ji,jj,jk+kp) = ah_wslp2(ji,jj,jk+kp) + zah * zbu * ze3wr * r1_e1e2t(ji,jj) * zslope2 
    172                   akz     (ji,jj,jk+kp) = akz     (ji,jj,jk+kp) + zah * r1_e1u(ji-ip,jj)       & 
    173                      &                                                      * r1_e1u(ji-ip,jj) * umask(ji-ip,jj,jk+kp) 
    174                      ! 
    175                END_3D 
    176             END DO 
     156         DO kp = 0, 1                            ! i-k triads 
     157            ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     158            DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 ) 
     159               ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     160               zbu   = e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
     161               zbu1  = e1e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm) 
     162               zah   = 0.25_wp * pahu(ji,jj,jk) 
     163               zah1  = 0.25_wp * pahu(ji-1,jj,jk) 
     164               ! Subtract s-coordinate slope at t-points to give slope rel to s-surfaces (do this by *adding* gradient of depth) 
     165               zslope2 = triadi_g(ji,jj,jk,1,kp) + ( gdept(ji+1,jj,jk,Kmm) - gdept(ji,jj,jk,Kmm) ) * r1_e1u(ji,jj) * umask(ji,jj,jk+kp) 
     166               zslope2 = zslope2 *zslope2 
     167               zslope21 = triadi_g(ji,jj,jk,0,kp) + ( gdept(ji,jj,jk,Kmm) - gdept(ji-1,jj,jk,Kmm) ) * r1_e1u(ji-1,jj) * umask(ji-1,jj,jk+kp) 
     168               zslope21 = zslope21 *zslope21 
     169               ! round brackets added to fix the order of floating point operations 
     170               ! needed to ensure halo 1 - halo 2 compatibility 
     171               ah_wslp2(ji,jj,jk+kp) =  ah_wslp2(ji,jj,jk+kp) + ( zah * zbu * ze3wr * r1_e1e2t(ji,jj) * zslope2                    & 
     172                        &                                       + zah1 * zbu1 * ze3wr * r1_e1e2t(ji,jj) * zslope21                 & 
     173                        &                                       )                                                                  ! bracket for halo 1 - halo 2 compatibility 
     174               akz     (ji,jj,jk+kp) =  akz     (ji,jj,jk+kp) + ( zah * r1_e1u(ji,jj) * r1_e1u(ji,jj) * umask(ji,jj,jk+kp)         & 
     175                                                                + zah1 * r1_e1u(ji-1,jj) * r1_e1u(ji-1,jj) * umask(ji-1,jj,jk+kp)  & 
     176                        &                                       )                                                                  ! bracket for halo 1 - halo 2 compatibility 
     177            END_3D 
    177178         END DO 
    178179         ! 
    179          DO jp = 0, 1                            ! j-k triads 
    180             DO kp = 0, 1 
    181                ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
    182                DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 )  
    183                   ze3wr = 1.0_wp / e3w(ji,jj,jk+kp,Kmm) 
    184                   zbv   = e1e2v(ji,jj-jp) * e3v(ji,jj-jp,jk,Kmm) 
    185                   zah   = 0.25_wp * pahv(ji,jj-jp,jk) 
    186                   zslope_skew = triadj_g(ji,jj,jk,1-jp,kp) 
    187                   ! Subtract s-coordinate slope at t-points to give slope rel to s surfaces 
    188                   !    (do this by *adding* gradient of depth) 
    189                   zslope2 = zslope_skew + ( gdept(ji,jj-jp+1,jk,Kmm) - gdept(ji,jj-jp,jk,Kmm) ) * r1_e2v(ji,jj-jp) * vmask(ji,jj-jp,jk+kp) 
    190                   zslope2 = zslope2 * zslope2 
    191                   ah_wslp2(ji,jj,jk+kp) = ah_wslp2(ji,jj,jk+kp) + zah * zbv * ze3wr * r1_e1e2t(ji,jj) * zslope2 
    192                   akz     (ji,jj,jk+kp) = akz     (ji,jj,jk+kp) + zah * r1_e2v(ji,jj-jp)     & 
    193                      &                                                      * r1_e2v(ji,jj-jp) * vmask(ji,jj-jp,jk+kp) 
    194                   ! 
    195                END_3D 
    196             END DO 
     180         DO kp = 0, 1                            ! j-k triads 
     181            ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     182            DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 ) 
     183               ze3wr = 1.0_wp / e3w(ji,jj,jk+kp,Kmm) 
     184               zbv   = e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
     185               zbv1   = e1e2v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) 
     186               zah   = 0.25_wp * pahv(ji,jj,jk) 
     187               zah1   = 0.25_wp * pahv(ji,jj-1,jk) 
     188               ! Subtract s-coordinate slope at t-points to give slope rel to s surfaces 
     189               !    (do this by *adding* gradient of depth) 
     190               zslope2 = triadj_g(ji,jj,jk,1,kp) + ( gdept(ji,jj+1,jk,Kmm) - gdept(ji,jj,jk,Kmm) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk+kp) 
     191               zslope2 = zslope2 * zslope2 
     192               zslope21 = triadj_g(ji,jj,jk,0,kp) + ( gdept(ji,jj,jk,Kmm) - gdept(ji,jj-1,jk,Kmm) ) * r1_e2v(ji,jj-1) * vmask(ji,jj-1,jk+kp) 
     193               zslope21 = zslope21 * zslope21 
     194               ! round brackets added to fix the order of floating point operations 
     195               ! needed to ensure halo 1 - halo 2 compatibility 
     196               ah_wslp2(ji,jj,jk+kp) = ah_wslp2(ji,jj,jk+kp) + ( zah * zbv * ze3wr * r1_e1e2t(ji,jj) * zslope2                     & 
     197                        &                                      + zah1 * zbv1 * ze3wr * r1_e1e2t(ji,jj) * zslope21                  & 
     198                        &                                      )                                                                   ! bracket for halo 1 - halo 2 compatibility 
     199               akz     (ji,jj,jk+kp) = akz     (ji,jj,jk+kp) + ( zah * r1_e2v(ji,jj) * r1_e2v(ji,jj) * vmask(ji,jj,jk+kp)          & 
     200                        &                                      + zah1 * r1_e2v(ji,jj-1) * r1_e2v(ji,jj-1) * vmask(ji,jj-1,jk+kp)   & 
     201                        &                                      )                                                                   ! bracket for halo 1 - halo 2 compatibility 
     202               ! 
     203            END_3D 
    197204         END DO 
    198205         ! 
     
    201208            IF( ln_traldf_blp ) THEN                ! bilaplacian operator 
    202209               ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    203                DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 )  
     210               DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 
    204211                  akz(ji,jj,jk) = 16._wp           & 
    205212                     &   * ah_wslp2   (ji,jj,jk)   & 
     
    210217            ELSEIF( ln_traldf_lap ) THEN              ! laplacian operator 
    211218               ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    212                DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 )  
     219               DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 
    213220                  ze3w_2 = e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) 
    214221                  zcoef0 = rDt * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2  ) 
     
    219226         ELSE                                    ! 33 flux set to zero with akz=ah_wslp2 ==>> computed in full implicit 
    220227            ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpk ) 
    221             DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpk )  
     228            DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpk ) 
    222229               akz(ji,jj,jk) = ah_wslp2(ji,jj,jk) 
    223230            END_3D 
     
    229236               IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) 
    230237 
    231                zpsi_uw(:,:,:) = 0._wp 
    232                zpsi_vw(:,:,:) = 0._wp 
    233  
    234                DO jp = 0, 1 
     238                  zpsi_uw(:,:,:) = 0._wp 
     239                  zpsi_vw(:,:,:) = 0._wp 
     240                   
    235241                  DO kp = 0, 1 
    236                      ! [comm_cleanup] ! DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 
    237                      DO_3D( nn_hls, nn_hls-1, nn_hls, nn_hls-1, 1, jpkm1 )  
    238                         zpsi_uw(ji,jj,jk+kp) = zpsi_uw(ji,jj,jk+kp) & 
    239                            & + 0.25_wp * aeiu(ji,jj,jk) * e2u(ji,jj) * triadi_g(ji+jp,jj,jk,1-jp,kp) 
    240                         zpsi_vw(ji,jj,jk+kp) = zpsi_vw(ji,jj,jk+kp) & 
    241                            & + 0.25_wp * aeiv(ji,jj,jk) * e1v(ji,jj) * triadj_g(ji,jj+jp,jk,1-jp,kp) 
     242                     DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 
     243                        ! round brackets added to fix the order of floating point operations 
     244                        ! needed to ensure halo 1 - halo 2 compatibility 
     245                        zpsi_uw(ji,jj,jk+kp) = zpsi_uw(ji,jj,jk+kp)                                     & 
     246                           & + ( 0.25_wp * aeiu(ji,jj,jk) * e2u(ji,jj) * triadi_g(ji,jj,jk,1,kp)        &  
     247                           &   + 0.25_wp * aeiu(ji,jj,jk) * e2u(ji,jj) * triadi_g(ji+1,jj,jk,0,kp)      & 
     248                           &   )                                                                        ! bracket for halo 1 - halo 2 compatibility 
     249                        zpsi_vw(ji,jj,jk+kp) = zpsi_vw(ji,jj,jk+kp)                                     & 
     250                           & + ( 0.25_wp * aeiv(ji,jj,jk) * e1v(ji,jj) * triadj_g(ji,jj,jk,1,kp)        & 
     251                           &   + 0.25_wp * aeiv(ji,jj,jk) * e1v(ji,jj) * triadj_g(ji,jj+1,jk,0,kp)      & 
     252                           &   )                                                                        ! bracket for halo 1 - halo 2 compatibility 
    242253                     END_3D 
    243254                  END DO 
    244                END DO 
    245255               CALL ldf_eiv_dia( zpsi_uw, zpsi_vw, Kmm ) 
    246256 
     
    261271         ! 
    262272         ! [comm_cleanup] ! DO_3D( 1, 0, 1, 0, 1, jpkm1 )    !==  before lateral T & S gradients at T-level jk  ==! 
    263          DO_3D( nn_hls, nn_hls-1, nn_hls, nn_hls-1, 1, jpkm1 )  
     273         DO_3D( nn_hls, nn_hls-1, nn_hls, nn_hls-1, 1, jpkm1 )    !==  before lateral T & S gradients at T-level jk  ==! 
    264274            zdit(ji,jj,jk) = ( pt(ji+1,jj  ,jk,jn) - pt(ji,jj,jk,jn) ) * umask(ji,jj,jk) 
    265275            zdjt(ji,jj,jk) = ( pt(ji  ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 
     
    267277         IF( ln_zps .AND. l_grad_zps ) THEN    ! partial steps: correction at top/bottom ocean level 
    268278            ! [comm_cleanup] ! DO_2D( 1, 0, 1, 0 )                    ! bottom level 
    269             DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )  
     279            DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )                    ! bottom level 
    270280               zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn) 
    271281               zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn) 
     
    273283            IF( ln_isfcav ) THEN                   ! top level (ocean cavities only) 
    274284               ! [comm_cleanup] ! DO_2D( 1, 0, 1, 0 ) 
    275                DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )  
     285               DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
    276286                  IF( miku(ji,jj)  > 1 )   zdit(ji,jj,miku(ji,jj) ) = pgui(ji,jj,jn) 
    277287                  IF( mikv(ji,jj)  > 1 )   zdjt(ji,jj,mikv(ji,jj) ) = pgvi(ji,jj,jn) 
     
    287297            !                    !==  Vertical tracer gradient at level jk and jk+1 
    288298            ! [comm_cleanup] ! DO_2D( 1, 1, 1, 1 ) 
    289             DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )  
     299            DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) 
    290300               zdkt3d(ji,jj,1) = ( pt(ji,jj,jk,jn) - pt(ji,jj,jk+1,jn) ) * tmask(ji,jj,jk+1) 
    291301            END_2D 
     
    295305            ELSE 
    296306               ! [comm_cleanup] ! DO_2D( 1, 1, 1, 1 ) 
    297                DO_2D( nn_hls, nn_hls, nn_hls, nn_hls )  
     307               DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) 
    298308                  zdkt3d(ji,jj,0) = ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) * tmask(ji,jj,jk) 
    299309               END_2D 
     
    303313            ! 
    304314            IF( ln_botmix_triad ) THEN 
    305                DO ip = 0, 1              !==  Horizontal & vertical fluxes 
    306                   DO kp = 0, 1 
    307                      ! [comm_cleanup] ! DO_2D( 1, 0, 1, 0 ) 
    308                      DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )  
    309                         ze1ur = r1_e1u(ji,jj) 
    310                         zdxt  = zdit(ji,jj,jk) * ze1ur 
    311                         ze3wr = 1._wp / e3w(ji+ip,jj,jk+kp,Kmm) 
    312                         zdzt  = zdkt3d(ji+ip,jj,kp) * ze3wr 
    313                         zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 
    314                         zslope_iso  = triadi  (ji+ip,jj,jk,1-ip,kp) 
    315                         ! 
    316                         zbu = 0.25_wp * e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
    317                         ! ln_botmix_triad is .T. don't mask zah for bottom half cells    !!gm ?????   ahu is masked.... 
    318                         zah = pahu(ji,jj,jk) 
    319                         zah_slp  = zah * zslope_iso 
    320                         IF( ln_ldfeiv )   zaei_slp = aeiu(ji,jj,jk) * zslope_skew 
    321                         zftu(ji   ,jj,jk   ) = zftu(ji   ,jj,jk   ) - ( zah * zdxt + (zah_slp - zaei_slp) * zdzt ) * zbu * ze1ur 
    322                         ztfw(ji+ip,jj,jk+kp) = ztfw(ji+ip,jj,jk+kp) - ( zah_slp + zaei_slp) * zdxt                 * zbu * ze3wr 
    323                      END_2D 
    324                   END DO 
     315               DO kp = 0, 1              !==  Horizontal & vertical fluxes 
     316                  ! [comm_cleanup] ! DO_2D( 1, 0, 1, 0 ) 
     317                  DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
     318                     ze1ur = r1_e1u(ji,jj) 
     319                     zdxt  = zdit(ji,jj,jk) * ze1ur 
     320                     zdxt_ip1  = zdit(ji+1,jj,jk) * r1_e1u(ji+1,jj) 
     321                     ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     322                     ze3wr_ip1 = 1._wp / e3w(ji+1,jj,jk+kp,Kmm) 
     323                     zdzt  = zdkt3d(ji,jj,kp) * ze3wr 
     324                     zdzt_ip1  = zdkt3d(ji+1,jj,kp) * ze3wr_ip1 
     325                     ! 
     326                     zbu = 0.25_wp * e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
     327                     zbu_ip1 = 0.25_wp * e1e2u(ji+1,jj) * e3u(ji+1,jj,jk,Kmm) 
     328                     ! ln_botmix_triad is .T. don't mask zah for bottom half cells    !!gm ?????   ahu is masked.... 
     329                     zah = pahu(ji,jj,jk)  
     330                     zah_ip1 = pahu(ji+1,jj,jk)  
     331                     zah_slp  = zah * triadi(ji,jj,jk,1,kp) 
     332                     zah_slp_ip1  = zah_ip1 * triadi(ji+1,jj,jk,1,kp) 
     333                     zah_slp1  = zah * triadi(ji+1,jj,jk,0,kp) 
     334                     IF( ln_ldfeiv )   THEN 
     335                        zaei_slp = aeiu(ji,jj,jk) * triadi_g(ji,jj,jk,1,kp) 
     336                        zaei_slp_ip1 = aeiu(ji+1,jj,jk) * triadi_g(ji+1,jj,jk,1,kp) 
     337                        zaei_slp1 = aeiu(ji,jj,jk) * triadi_g(ji+1,jj,jk,0,kp) 
     338                     ENDIF 
     339                     ! round brackets added to fix the order of floating point operations 
     340                     ! needed to ensure halo 1 - halo 2 compatibility 
     341                     zftu(ji   ,jj,jk  ) =  zftu(ji   ,jj,jk )                                                               & 
     342                                         &    - ( ( zah * zdxt + ( zah_slp - zaei_slp ) * zdzt ) * zbu * ze1ur               & 
     343                                         &      + ( zah * zdxt + zah_slp1 * zdzt_ip1 - zaei_slp1 * zdzt_ip1 ) * zbu * ze1ur  & 
     344                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     345                     ztfw(ji+1,jj,jk+kp) =  ztfw(ji+1,jj,jk+kp)                                                              &  
     346                                         &    - ( (zah_slp_ip1 + zaei_slp_ip1) * zdxt_ip1 * zbu_ip1 * ze3wr_ip1              & 
     347                                         &      + ( zah_slp1 + zaei_slp1) * zdxt * zbu * ze3wr_ip1                           &  
     348                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     349                  END_2D 
    325350               END DO 
    326351               ! 
    327                DO jp = 0, 1 
    328                   DO kp = 0, 1 
    329                      ! [comm_cleanup] ! DO_2D( 1, 0, 1, 0 ) 
    330                      DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )  
    331                         ze2vr = r1_e2v(ji,jj) 
    332                         zdyt  = zdjt(ji,jj,jk) * ze2vr 
    333                         ze3wr = 1._wp / e3w(ji,jj+jp,jk+kp,Kmm) 
    334                         zdzt  = zdkt3d(ji,jj+jp,kp) * ze3wr 
    335                         zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 
    336                         zslope_iso  = triadj(ji,jj+jp,jk,1-jp,kp) 
    337                         zbv = 0.25_wp * e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
    338                         ! ln_botmix_triad is .T. don't mask zah for bottom half cells    !!gm ?????  ahv is masked... 
    339                         zah = pahv(ji,jj,jk) 
    340                         zah_slp = zah * zslope_iso 
    341                         IF( ln_ldfeiv )   zaei_slp = aeiv(ji,jj,jk) * zslope_skew 
    342                         zftv(ji,jj   ,jk   ) = zftv(ji,jj   ,jk   ) - ( zah * zdyt + (zah_slp - zaei_slp) * zdzt ) * zbv * ze2vr 
    343                         ztfw(ji,jj+jp,jk+kp) = ztfw(ji,jj+jp,jk+kp) - ( zah_slp + zaei_slp ) * zdyt                * zbv * ze3wr 
    344                      END_2D 
    345                   END DO 
     352               DO kp = 0, 1 
     353                  ! [comm_cleanup] ! DO_2D( 1, 0, 1, 0 ) 
     354                  DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
     355                     ze2vr = r1_e2v(ji,jj) 
     356                     zdyt  = zdjt(ji,jj,jk) * ze2vr 
     357                     zdyt_jp1  = zdjt(ji,jj+1,jk) * r1_e2v(ji,jj+1) 
     358                     ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     359                     ze3wr_jp1 = 1._wp / e3w(ji,jj+1,jk+kp,Kmm) 
     360                     zdzt  = zdkt3d(ji,jj,kp) * ze3wr 
     361                     zdzt_jp1  = zdkt3d(ji,jj+1,kp) * ze3wr_jp1 
     362                     zbv = 0.25_wp * e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
     363                     zbv_jp1 = 0.25_wp * e1e2v(ji,jj+1) * e3v(ji,jj+1,jk,Kmm) 
     364                     ! ln_botmix_triad is .T. don't mask zah for bottom half cells    !!gm ?????   ahu is masked.... 
     365                     zah = pahv(ji,jj,jk)          ! pahv(ji,jj+jp,jk)  ???? 
     366                     zah_jp1 = pahv(ji,jj+1,jk)  
     367                     zah_slp = zah * triadj(ji,jj,jk,1,kp) 
     368                     zah_slp1 = zah * triadj(ji,jj+1,jk,0,kp) 
     369                     zah_slp_jp1 = zah_jp1 * triadj(ji,jj+1,jk,1,kp) 
     370                     IF( ln_ldfeiv )   THEN 
     371                        zaei_slp = aeiv(ji,jj,jk) * triadj_g(ji,jj,jk,1,kp)   
     372                        zaei_slp_jp1 = aeiv(ji,jj+1,jk) * triadj_g(ji,jj+1,jk,1,kp)   
     373                        zaei_slp1 = aeiv(ji,jj,jk) * triadj_g(ji,jj+1,jk,0,kp) 
     374                     ENDIF 
     375                     ! round brackets added to fix the order of floating point operations 
     376                     ! needed to ensure halo 1 - halo 2 compatibility 
     377                     zftv(ji,jj  ,jk   ) =  zftv(ji,jj  ,jk   )                                                              & 
     378                                         &    - ( ( zah * zdyt + ( zah_slp - zaei_slp ) * zdzt ) * zbv * ze2vr               & 
     379                                         &      + ( zah * zdyt + zah_slp1 * zdzt_jp1 - zaei_slp1 * zdzt_jp1 ) * zbv * ze2vr  & 
     380                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     381                     ztfw(ji,jj+1,jk+kp) =  ztfw(ji,jj+1,jk+kp)                                                              & 
     382                                         &    - ( ( zah_slp_jp1 + zaei_slp_jp1) * zdyt_jp1 * zbv_jp1 * ze3wr_jp1             & 
     383                                         &      + ( zah_slp1 + zaei_slp1) * zdyt * zbv * ze3wr_jp1                           & 
     384                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     385                  END_2D 
    346386               END DO 
    347387               ! 
    348388            ELSE 
    349389               ! 
    350                DO ip = 0, 1               !==  Horizontal & vertical fluxes 
    351                   DO kp = 0, 1 
    352                      ! [comm_cleanup] ! DO_2D( 1, 0, 1, 0 ) 
    353                      DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )  
    354                         ze1ur = r1_e1u(ji,jj) 
    355                         zdxt  = zdit(ji,jj,jk) * ze1ur 
    356                         ze3wr = 1._wp / e3w(ji+ip,jj,jk+kp,Kmm) 
    357                         zdzt  = zdkt3d(ji+ip,jj,kp) * ze3wr 
    358                         zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 
    359                         zslope_iso  = triadi(ji+ip,jj,jk,1-ip,kp) 
    360                         ! 
    361                         zbu = 0.25_wp * e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
    362                         ! ln_botmix_triad is .F. mask zah for bottom half cells 
    363                         zah = pahu(ji,jj,jk) * umask(ji,jj,jk+kp)         ! pahu(ji+ip,jj,jk)   ===>>  ???? 
    364                         zah_slp  = zah * zslope_iso 
    365                         IF( ln_ldfeiv )   zaei_slp = aeiu(ji,jj,jk) * zslope_skew        ! aeit(ji+ip,jj,jk)*zslope_skew 
    366                         zftu(ji   ,jj,jk   ) = zftu(ji   ,jj,jk   ) - ( zah * zdxt + (zah_slp - zaei_slp) * zdzt ) * zbu * ze1ur 
    367                         ztfw(ji+ip,jj,jk+kp) = ztfw(ji+ip,jj,jk+kp) - (zah_slp + zaei_slp) * zdxt * zbu * ze3wr 
    368                      END_2D 
    369                   END DO 
     390               DO kp = 0, 1               !==  Horizontal & vertical fluxes 
     391                  ! [comm_cleanup] ! DO_2D( 1, 0, 1, 0 ) 
     392                  DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
     393                     ze1ur = r1_e1u(ji,jj) 
     394                     zdxt  = zdit(ji,jj,jk) * ze1ur 
     395                     zdxt_ip1  = zdit(ji+1,jj,jk) * r1_e1u(ji+1,jj) 
     396                     ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     397                     ze3wr_ip1 = 1._wp / e3w(ji+1,jj,jk+kp,Kmm) 
     398                     zdzt  = zdkt3d(ji,jj,kp) * ze3wr 
     399                     zdzt_ip1  = zdkt3d(ji+1,jj,kp) * ze3wr_ip1 
     400                     ! 
     401                     zbu = 0.25_wp * e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
     402                     zbu_ip1 = 0.25_wp * e1e2u(ji+1,jj) * e3u(ji+1,jj,jk,Kmm) 
     403                     ! ln_botmix_triad is .F. mask zah for bottom half cells 
     404                     zah = pahu(ji,jj,jk) * umask(ji,jj,jk+kp)         ! pahu(ji+ip,jj,jk)   ===>>  ???? 
     405                     zah_ip1 = pahu(ji+1,jj,jk) * umask(ji+1,jj,jk+kp)  
     406                     zah_slp  = zah * triadi(ji,jj,jk,1,kp) 
     407                     zah_slp_ip1  = zah_ip1 * triadi(ji+1,jj,jk,1,kp) 
     408                     zah_slp1  = zah * triadi(ji+1,jj,jk,0,kp) 
     409                     IF( ln_ldfeiv )   THEN 
     410                        zaei_slp = aeiu(ji,jj,jk) * triadi_g(ji,jj,jk,1,kp)  
     411                        zaei_slp_ip1 = aeiu(ji+1,jj,jk) * triadi_g(ji+1,jj,jk,1,kp) 
     412                        zaei_slp1 = aeiu(ji,jj,jk) * triadi_g(ji+1,jj,jk,0,kp) 
     413                     ENDIF 
     414!                     zftu(ji   ,jj,jk  ) = zftu(ji   ,jj,jk ) - ( zah * zdxt + (zah_slp - zaei_slp) * zdzt ) * zbu * ze1ur - ( zah * zdxt + (zah_slp1 - zaei_slp1) * zdzt_ip1 ) * zbu * ze1ur 
     415!                     ztfw(ji+1,jj,jk+kp) = ztfw(ji+1,jj,jk+kp) - (zah_slp_ip1 + zaei_slp_ip1) * zdxt_ip1 * zbu_ip1 * ze3wr_ip1 - (zah_slp1 + zaei_slp1) * zdxt * zbu * ze3wr_ip1 
     416                     ! round brackets added to fix the order of floating point operations 
     417                     ! needed to ensure halo 1 - halo 2 compatibility 
     418                     zftu(ji   ,jj,jk  ) =  zftu(ji   ,jj,jk )                                                               & 
     419                                         &    - ( ( zah * zdxt + ( zah_slp - zaei_slp ) * zdzt ) * zbu * ze1ur               & 
     420                                         &      + ( zah * zdxt + zah_slp1 * zdzt_ip1 - zaei_slp1 * zdzt_ip1 ) * zbu * ze1ur  & 
     421                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     422                     ztfw(ji+1,jj,jk+kp) =  ztfw(ji+1,jj,jk+kp)                                                              &  
     423                                         &    - ( (zah_slp_ip1 + zaei_slp_ip1) * zdxt_ip1 * zbu_ip1 * ze3wr_ip1              & 
     424                                         &      + ( zah_slp1 + zaei_slp1) * zdxt * zbu * ze3wr_ip1                           & 
     425                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     426                  END_2D 
    370427               END DO 
    371428               ! 
    372                DO jp = 0, 1 
    373                   DO kp = 0, 1 
    374                      ! [comm_cleanup] ! DO_2D( 1, 0, 1, 0 ) 
    375                      DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )  
    376                         ze2vr = r1_e2v(ji,jj) 
    377                         zdyt  = zdjt(ji,jj,jk) * ze2vr 
    378                         ze3wr = 1._wp / e3w(ji,jj+jp,jk+kp,Kmm) 
    379                         zdzt  = zdkt3d(ji,jj+jp,kp) * ze3wr 
    380                         zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 
    381                         zslope_iso  = triadj(ji,jj+jp,jk,1-jp,kp) 
    382                         zbv = 0.25_wp * e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
    383                         ! ln_botmix_triad is .F. mask zah for bottom half cells 
    384                         zah = pahv(ji,jj,jk) * vmask(ji,jj,jk+kp)         ! pahv(ji,jj+jp,jk)  ???? 
    385                         zah_slp = zah * zslope_iso 
    386                         IF( ln_ldfeiv )   zaei_slp = aeiv(ji,jj,jk) * zslope_skew        ! aeit(ji,jj+jp,jk)*zslope_skew 
    387                         zftv(ji,jj,jk) = zftv(ji,jj,jk) - ( zah * zdyt + (zah_slp - zaei_slp) * zdzt ) * zbv * ze2vr 
    388                         ztfw(ji,jj+jp,jk+kp) = ztfw(ji,jj+jp,jk+kp) - (zah_slp + zaei_slp) * zdyt * zbv * ze3wr 
    389                      END_2D 
    390                   END DO 
     429               DO kp = 0, 1 
     430                  ! [comm_cleanup] ! DO_2D( 1, 0, 1, 0 ) 
     431                  DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
     432                     ze2vr = r1_e2v(ji,jj) 
     433                     zdyt  = zdjt(ji,jj,jk) * ze2vr 
     434                     zdyt_jp1  = zdjt(ji,jj+1,jk) * r1_e2v(ji,jj+1) 
     435                     ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     436                     ze3wr_jp1 = 1._wp / e3w(ji,jj+1,jk+kp,Kmm) 
     437                     zdzt  = zdkt3d(ji,jj,kp) * ze3wr 
     438                     zdzt_jp1  = zdkt3d(ji,jj+1,kp) * ze3wr_jp1 
     439                     zbv = 0.25_wp * e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
     440                     zbv_jp1 = 0.25_wp * e1e2v(ji,jj+1) * e3v(ji,jj+1,jk,Kmm) 
     441                     ! ln_botmix_triad is .F. mask zah for bottom half cells 
     442                     zah = pahv(ji,jj,jk) * vmask(ji,jj,jk+kp)         ! pahv(ji,jj+jp,jk)  ???? 
     443                     zah_jp1 = pahv(ji,jj+1,jk) * vmask(ji,jj+1,jk+kp) 
     444                     zah_slp = zah * triadj(ji,jj,jk,1,kp) 
     445                     zah_slp1 = zah * triadj(ji,jj+1,jk,0,kp) 
     446                     zah_slp_jp1 = zah_jp1 * triadj(ji,jj+1,jk,1,kp) 
     447                     IF( ln_ldfeiv )   THEN 
     448                        zaei_slp = aeiv(ji,jj,jk) * triadj_g(ji,jj,jk,1,kp)  
     449                        zaei_slp_jp1 = aeiv(ji,jj+1,jk) * triadj_g(ji,jj+1,jk,1,kp) 
     450                        zaei_slp1 = aeiv(ji,jj,jk) * triadj_g(ji,jj+1,jk,0,kp) 
     451                     ENDIF 
     452!                     zftv(ji,jj  ,jk   ) = zftv(ji,jj  ,jk   ) - ( zah * zdyt + (zah_slp - zaei_slp) * zdzt ) * zbv * ze2vr - ( zah * zdyt + (zah_slp1 - zaei_slp1) * zdzt_jp1 ) * zbv * ze2vr 
     453!                     ztfw(ji,jj+1,jk+kp) = ztfw(ji,jj+1,jk+kp) - ( zah_slp_jp1 + zaei_slp_jp1) * zdyt_jp1 * zbv_jp1 * ze3wr_jp1 - (zah_slp1 + zaei_slp1) * zdyt * zbv * ze3wr_jp1 
     454                     ! round brackets added to fix the order of floating point operations 
     455                     ! needed to ensure halo 1 - halo 2 compatibility 
     456                     zftv(ji,jj  ,jk   ) =  zftv(ji,jj  ,jk   )                                                              & 
     457                                         &    - ( ( zah * zdyt + ( zah_slp - zaei_slp ) * zdzt ) * zbv * ze2vr               & 
     458                                         &      + ( zah * zdyt + zah_slp1 * zdzt_jp1 - zaei_slp1 * zdzt_jp1 ) * zbv * ze2vr  & 
     459                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     460                     ztfw(ji,jj+1,jk+kp) =  ztfw(ji,jj+1,jk+kp)                                                              & 
     461                                         &    - ( ( zah_slp_jp1 + zaei_slp_jp1) * zdyt_jp1 * zbv_jp1 * ze3wr_jp1             & 
     462                                         &      + ( zah_slp1 + zaei_slp1) * zdyt * zbv * ze3wr_jp1                           & 
     463                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     464                  END_2D 
    391465               END DO 
    392466            ENDIF 
    393467            !                             !==  horizontal divergence and add to the general trend  ==! 
    394468            ! [comm_cleanup] ! DO_2D( 0, 0, 0, 0 ) 
    395             DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )  
    396                pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)    & 
    397                   &                       + zsign * (  zftu(ji-1,jj  ,jk) - zftu(ji,jj,jk)       & 
    398                   &                                           + zftv(ji,jj-1,jk) - zftv(ji,jj,jk)   )   & 
    399                   &                                        / (  e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm)  ) 
     469            DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 
     470               ! round brackets added to fix the order of floating point operations 
     471               ! needed to ensure halo 1 - halo 2 compatibility 
     472               pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)                                                & 
     473                  &                       + zsign * ( ( zftu(ji-1,jj  ,jk) - zftu(ji,jj,jk)             & 
     474                  &                                   )                                                 & ! bracket for halo 1 - halo 2 compatibility 
     475                  &                                 + ( zftv(ji,jj-1,jk) - zftv(ji,jj,jk)               & 
     476                  &                                   )                                                 & ! bracket for halo 1 - halo 2 compatibility 
     477                  &                                 ) / (  e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm)  ) 
    400478            END_2D 
    401479            ! 
     
    405483         IF( ln_traldf_lap ) THEN               ! laplacian case: eddy coef = ah_wslp2 - akz 
    406484            ! [comm_cleanup] ! DO_3D( 0, 0, 1, 0, 2, jpkm1 ) 
    407             DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 )  
     485            DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 
    408486               ztfw(ji,jj,jk) = ztfw(ji,jj,jk) - e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * tmask(ji,jj,jk)   & 
    409487                  &                            * ( ah_wslp2(ji,jj,jk) - akz(ji,jj,jk) )             & 
     
    414492            CASE(  1  )                            ! 1st pass : eddy coef = ah_wslp2 
    415493               ! [comm_cleanup] ! DO_3D( 0, 0, 1, 0, 2, jpkm1 ) 
    416                DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 )  
     494               DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 
    417495                  ztfw(ji,jj,jk) = ztfw(ji,jj,jk) - e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * tmask(ji,jj,jk)             & 
    418496                     &                            * ah_wslp2(ji,jj,jk) * ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) 
     
    420498            CASE(  2  )                            ! 2nd pass : eddy flux = ah_wslp2 and akz applied on pt  and pt2 gradients, resp. 
    421499               ! [comm_cleanup] ! DO_3D( 0, 0, 1, 0, 2, jpkm1 ) 
    422                DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 )  
     500               DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 2, jpkm1 ) 
    423501                  ztfw(ji,jj,jk) = ztfw(ji,jj,jk) - e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * tmask(ji,jj,jk)                      & 
    424502                     &                            * (  ah_wslp2(ji,jj,jk) * ( pt (ji,jj,jk-1,jn) - pt (ji,jj,jk,jn) )   & 
     
    429507         ! 
    430508         ! [comm_cleanup] ! DO_3D( 0, 0, 0, 0, 1, jpkm1 )      !==  Divergence of vertical fluxes added to pta  ==! 
    431          DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 )  
     509         DO_3D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1, 1, jpkm1 )      !==  Divergence of vertical fluxes added to pta  ==! 
    432510            pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)    & 
    433511            &                                  + zsign * (  ztfw(ji,jj,jk+1) - ztfw(ji,jj,jk)  )   & 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/TRA/zpshde.F90

    r14609 r14667  
    130130      DO jn = 1, kjpt      !==   Interpolation of tracers at the last ocean level   ==! 
    131131         ! 
    132          DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )              ! Gradient of density at the last level 
     132         DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )              ! Gradient of density at the last level 
    133133            iku = mbku(ji,jj)   ;   ikum1 = MAX( iku - 1 , 1 )    ! last and before last ocean level at u- & v-points 
    134134            ikv = mbkv(ji,jj)   ;   ikvm1 = MAX( ikv - 1 , 1 )    ! if level first is a p-step, ik.m1=1 
     
    174174         pgru(:,:) = 0._wp 
    175175         pgrv(:,:) = 0._wp                ! depth of the partial step level 
    176          DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 ) 
     176         DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 ) 
    177177            iku = mbku(ji,jj) 
    178178            ikv = mbkv(ji,jj) 
     
    190190         CALL eos( ztj, zhj, zrj )        ! at the partial step depth output in  zri, zrj 
    191191         ! 
    192          DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )              ! Gradient of density at the last level 
     192         DO_2D( nn_hls, nn_hls-1, nn_hls, nn_hls-1 )              ! Gradient of density at the last level 
    193193            iku = mbku(ji,jj) 
    194194            ikv = mbkv(ji,jj) 
  • NEMO/branches/2021/dev_r14393_HPC-03_Mele_Comm_Cleanup/src/OCE/stpmlf.F90

    r14609 r14667  
    187187            CALL lbc_lnk( 'stp_MLF', avm  , 'W', 1.0_wp , avt  , 'W', 1.0_wp , avs , 'W', 1.0_wp ) 
    188188         ENDIF 
    189          ! CALL lbc_lnk( 'stp_MLF', uu(:,:,:,Nnn), 'U', -1.0_wp, vv(:,:,:,Nnn), 'V', -1.0_wp ) 
    190          ! IF(.NOT.lk_linssh) CALL lbc_lnk( 'stp_MLF', r3u(:,:,Nnn), 'U', 1.0_wp, r3v(:,:,Nnn), 'V', 1.0_wp, r3t(:,:,Nnn), 'T', 1.0_wp ) 
    191189      ENDIF 
    192190                         CALL zdf_phy( kstp, Nbb, Nnn, Nrhs )   ! vertical physics update (top/bot drag, avt, avs, avm + MLD) 
     
    286284                         CALL ssh_atf    ( kstp, Nbb, Nnn, Naa, ssh )            ! time filtering of "now" sea surface height 
    287285      IF(.NOT.lk_linssh) CALL dom_qco_r3c( ssh(:,:,Nnn), r3t_f, r3u_f, r3v_f )   ! "now" ssh/h_0 ratio from filtrered ssh 
     286         ! [comm_cleanup] this should not be needed 
     287         IF(nn_hls.eq.2.AND..NOT.lk_linssh) CALL lbc_lnk( 'stp_MLF', r3u_f, 'U', 1.0_wp, r3v_f, 'V', 1.0_wp ) 
    288288#if defined key_top 
    289289      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 
Note: See TracChangeset for help on using the changeset viewer.