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 14666 – NEMO

Changeset 14666


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

TRA & DYN changes for halo 1 - halo 2 compatibility

Location:
NEMO/branches/2021/ticket2607_r14608_halo1_halo2_compatibility/src/OCE
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2021/ticket2607_r14608_halo1_halo2_compatibility/src/OCE/DOM/domqco.F90

    r14433 r14666  
    180180 
    181181            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) 
     182               ! round brackets added to fix the order of floating point operations 
     183               ! needed to ensure halo 1 - halo 2 compatibility 
     184               pr3f(ji,jj) = 0.25_wp * ( ( e1e2t(ji  ,jj  ) * pssh(ji  ,jj  )   & 
     185                  &                      + e1e2t(ji+1,jj  ) * pssh(ji+1,jj  )   & 
     186                  &                      )                                      & ! bracket for halo 1 - halo 2 compatibility 
     187                  &                     + ( e1e2t(ji  ,jj+1) * pssh(ji  ,jj+1)  & 
     188                  &                       + e1e2t(ji+1,jj+1) * pssh(ji+1,jj+1)  & 
     189                  &                       )                                     & ! bracket for halo 1 - halo 2 compatibility 
     190                  &                    ) * r1_hf_0(ji,jj) * r1_e1e2f(ji,jj) 
    186191            END_2D 
    187192!!st         ELSE                                      !- Flux Form   (simple averaging) 
    188193#else 
    189194            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) 
     195               ! round brackets added to fix the order of floating point operations 
     196               ! needed to ensure halo 1 - halo 2 compatibility 
     197               pr3f(ji,jj) = 0.25_wp * ( ( pssh(ji,jj  ) + pssh(ji+1,jj  ) ) & 
     198                  &                     + ( pssh(ji,jj+1) + pssh(ji+1,jj+1)  &  
     199                  &                       )                                  & ! bracket for halo 1 - halo 2 compatibility 
     200                  &                    ) * r1_hf_0(ji,jj) 
    192201            END_2D 
    193202!!st         ENDIF 
  • NEMO/branches/2021/ticket2607_r14608_halo1_halo2_compatibility/src/OCE/DYN/divhor.F90

    r13558 r14666  
    7878      ! 
    7979      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) 
     80         ! round brackets added to fix the order of floating point operations 
     81         ! needed to ensure halo 1 - halo 2 compatibility 
     82         hdiv(ji,jj,jk) = (  ( e2u(ji  ,jj) * e3u(ji  ,jj,jk,Kmm) * uu(ji  ,jj,jk,Kmm)     & 
     83            &                - e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm) * uu(ji-1,jj,jk,Kmm)     & 
     84            &                )                                                             & ! bracket for halo 1 - halo 2 compatibility 
     85            &              + ( e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kmm) * vv(ji,jj  ,jk,Kmm)     & 
     86            &                - e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) * vv(ji,jj-1,jk,Kmm)     &  
     87            &                )                                                             & ! bracket for halo 1 - halo 2 compatibility 
     88            &             )  * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    8589      END_3D 
    8690      ! 
  • NEMO/branches/2021/ticket2607_r14608_halo1_halo2_compatibility/src/OCE/TRA/traldf_iso.F90

    r14072 r14666  
    179179               &                           + vmask(ji,jj-1,jk-1) + vmask(ji,jj  ,jk) , 1._wp  ) 
    180180               ! 
    181             zahu_w = (   pahu(ji  ,jj,jk-1) + pahu(ji-1,jj,jk)    & 
    182                &       + pahu(ji-1,jj,jk-1) + pahu(ji  ,jj,jk)  ) * zmsku 
    183             zahv_w = (   pahv(ji,jj  ,jk-1) + pahv(ji,jj-1,jk)    & 
    184                &       + pahv(ji,jj-1,jk-1) + pahv(ji,jj  ,jk)  ) * zmskv 
     181            ! round brackets added to fix the order of floating point operations 
     182            ! needed to ensure halo 1 - halo 2 compatibility 
     183            zahu_w = ( (  pahu(ji  ,jj,jk-1) + pahu(ji-1,jj,jk)                    & 
     184               &       )                                                           & ! bracket for halo 1 - halo 2 compatibility 
     185               &       + ( pahu(ji-1,jj,jk-1) + pahu(ji  ,jj,jk)                   & 
     186               &         )                                                         & ! bracket for halo 1 - halo 2 compatibility 
     187               &     ) * zmsku                                        
     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               &         )                                                         & ! bracket for halo 1 - halo 2 compatibility 
     192               &     ) * zmskv                                                
    185193               ! 
    186194            ah_wslp2(ji,jj,jk) = zahu_w * wslpi(ji,jj,jk) * wslpi(ji,jj,jk)   & 
     
    190198         IF( ln_traldf_msc ) THEN                ! stabilizing vertical diffusivity coefficient 
    191199            DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
     200               ! round brackets added to fix the order of floating point operations 
     201               ! needed to ensure halo 1 - halo 2 compatibility 
    192202               akz(ji,jj,jk) = 0.25_wp * (                                                                     & 
    193                   &              ( 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) )   & 
    194204                  &            + ( pahu(ji-1,jj,jk) + pahu(ji-1,jj,jk-1) ) / ( e1u(ji-1,jj) * e1u(ji-1,jj) )   & 
    195                   &            + ( pahv(ji,jj  ,jk) + pahv(ji,jj  ,jk-1) ) / ( e2v(ji,jj  ) * e2v(ji,jj  ) )   & 
    196                   &            + ( 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 
     209                  &                      )                                                                          
    197210            END_3D 
    198211            ! 
     
    278291               zcof2 = - pahv(ji,jj,jk) * e1v(ji,jj) * vslp(ji,jj,jk) * zmskv 
    279292               ! 
    280                zftu(ji,jj,jk ) = (  zabe1 * zdit(ji,jj,jk)   & 
    281                   &               + zcof1 * (  zdkt (ji+1,jj) + zdk1t(ji,jj)      & 
    282                   &                          + zdk1t(ji+1,jj) + zdkt (ji,jj)  )  ) * umask(ji,jj,jk) 
    283                zftv(ji,jj,jk) = (  zabe2 * zdjt(ji,jj,jk)   & 
    284                   &               + zcof2 * (  zdkt (ji,jj+1) + zdk1t(ji,jj)      & 
    285                   &                          + zdk1t(ji,jj+1) + zdkt (ji,jj)  )  ) * vmask(ji,jj,jk) 
     293               ! round brackets added to fix the order of floating point operations 
     294               ! needed to ensure halo 1 - halo 2 compatibility 
     295               zftu(ji,jj,jk ) = (  zabe1 * zdit(ji,jj,jk)                       & 
     296                  &               + zcof1 * ( ( zdkt (ji+1,jj) + zdk1t(ji,jj)    & 
     297                  &                           )                                  & ! bracket for halo 1 - halo 2 compatibility 
     298                  &                         + ( zdk1t(ji+1,jj) + zdkt (ji,jj)    &  
     299                  &                           )                                  & ! bracket for halo 1 - halo 2 compatibility 
     300                  &                         ) ) * umask(ji,jj,jk) 
     301               zftv(ji,jj,jk) = (  zabe2 * zdjt(ji,jj,jk)                        & 
     302                  &              + zcof2 * ( ( zdkt (ji,jj+1) + zdk1t(ji,jj)     & 
     303                  &                           )                                  & ! bracket for halo 1 - halo 2 compatibility 
     304                  &                         + ( zdk1t(ji,jj+1) + zdkt (ji,jj)    & 
     305                  &                           )                                  & ! bracket for halo 1 - halo 2 compatibility 
     306                  &                         ) ) * vmask(ji,jj,jk)       
    286307            END_2D 
    287308            ! 
    288309            DO_2D( 0, 0, 0, 0 )           !== horizontal divergence and add to pta 
    289                pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)    & 
    290                   &       + zsign * (  zftu(ji,jj,jk) - zftu(ji-1,jj,jk) + zftv(ji,jj,jk) - zftv(ji,jj-1,jk)  )   & 
    291                   &                                              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     310               ! round brackets added to fix the order of floating point operations 
     311               ! needed to ensure halo 1 - halo 2 compatibility 
     312               pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)                         & 
     313                  &       + zsign * ( ( zftu(ji,jj,jk) - zftu(ji-1,jj,jk)        & 
     314                  &                   )                                          & ! bracket for halo 1 - halo 2 compatibility 
     315                  &                 + ( zftv(ji,jj,jk) - zftv(ji,jj-1,jk)        & 
     316                  &                   )                                          & ! bracket for halo 1 - halo 2 compatibility 
     317                  &                 ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    292318            END_2D 
    293319         END DO                                        !   End of slab 
     
    317343            zcoef4 = - zahv_w * e1t(ji,jj) * zmskv * wslpj (ji,jj,jk) 
    318344            ! 
    319             ztfw(ji,jj,jk) = zcoef3 * (   zdit(ji  ,jj  ,jk-1) + zdit(ji-1,jj  ,jk)      & 
    320                &                        + zdit(ji-1,jj  ,jk-1) + zdit(ji  ,jj  ,jk)  )   & 
    321                &           + zcoef4 * (   zdjt(ji  ,jj  ,jk-1) + zdjt(ji  ,jj-1,jk)      & 
    322                &                        + zdjt(ji  ,jj-1,jk-1) + zdjt(ji  ,jj  ,jk)  ) 
     345            ! round brackets added to fix the order of floating point operations 
     346            ! needed to ensure halo 1 - halo 2 compatibility 
     347            ztfw(ji,jj,jk) = zcoef3 * ( ( zdit(ji  ,jj  ,jk-1) + zdit(ji-1,jj  ,jk)    & 
     348                  &                     )                                              & ! bracket for halo 1 - halo 2 compatibility 
     349                  &                   + ( zdit(ji-1,jj  ,jk-1) + zdit(ji  ,jj  ,jk)    & 
     350                  &                     )                                              & ! bracket for halo 1 - halo 2 compatibility  
     351                  &                   )                                                & 
     352                  &        + zcoef4 * ( ( zdjt(ji  ,jj  ,jk-1) + zdjt(ji  ,jj-1,jk)    & 
     353                  &                     )                                              & ! bracket for halo 1 - halo 2 compatibility 
     354                  &                   + ( zdjt(ji  ,jj-1,jk-1) + zdjt(ji  ,jj  ,jk)    & 
     355                  &                     )                                              & ! bracket for halo 1 - halo 2 compatibility 
     356                  &                   ) 
    323357         END_3D 
    324358         !                                !==  add the vertical 33 flux  ==! 
  • NEMO/branches/2021/ticket2607_r14608_halo1_halo2_compatibility/src/OCE/TRA/traldf_lap_blp.F90

    r14215 r14666  
    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         ! 
  • NEMO/branches/2021/ticket2607_r14608_halo1_halo2_compatibility/src/OCE/TRA/traldf_triad.F90

    r14215 r14666  
    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 
     
    157153         END_3D 
    158154         ! 
    159          DO ip = 0, 1                            ! i-k triads 
    160             DO kp = 0, 1 
    161                DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
    162                   ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
    163                   zbu   = e1e2u(ji-ip,jj) * e3u(ji-ip,jj,jk,Kmm) 
    164                   zah   = 0.25_wp * pahu(ji-ip,jj,jk) 
    165                   zslope_skew = triadi_g(ji,jj,jk,1-ip,kp) 
    166                   ! Subtract s-coordinate slope at t-points to give slope rel to s-surfaces (do this by *adding* gradient of depth) 
    167                   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) 
    168                   zslope2 = zslope2 *zslope2 
    169                   ah_wslp2(ji,jj,jk+kp) = ah_wslp2(ji,jj,jk+kp) + zah * zbu * ze3wr * r1_e1e2t(ji,jj) * zslope2 
    170                   akz     (ji,jj,jk+kp) = akz     (ji,jj,jk+kp) + zah * r1_e1u(ji-ip,jj)       & 
    171                      &                                                      * r1_e1u(ji-ip,jj) * umask(ji-ip,jj,jk+kp) 
    172                      ! 
    173                END_3D 
    174             END DO 
     155         DO kp = 0, 1                            ! i-k triads 
     156            DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     157               ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     158               zbu   = e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
     159               zbu1  = e1e2u(ji-1,jj) * e3u(ji-1,jj,jk,Kmm) 
     160               zah   = 0.25_wp * pahu(ji,jj,jk) 
     161               zah1  = 0.25_wp * pahu(ji-1,jj,jk) 
     162               ! Subtract s-coordinate slope at t-points to give slope rel to s-surfaces (do this by *adding* gradient of depth) 
     163               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) 
     164               zslope2 = zslope2 *zslope2 
     165               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) 
     166               zslope21 = zslope21 *zslope21 
     167               ! round brackets added to fix the order of floating point operations 
     168               ! needed to ensure halo 1 - halo 2 compatibility 
     169               ah_wslp2(ji,jj,jk+kp) =  ah_wslp2(ji,jj,jk+kp) + ( zah * zbu * ze3wr * r1_e1e2t(ji,jj) * zslope2                    & 
     170                        &                                       + zah1 * zbu1 * ze3wr * r1_e1e2t(ji,jj) * zslope21                 & 
     171                        &                                       )                                                                  ! bracket for halo 1 - halo 2 compatibility 
     172               akz     (ji,jj,jk+kp) =  akz     (ji,jj,jk+kp) + ( zah * r1_e1u(ji,jj) * r1_e1u(ji,jj) * umask(ji,jj,jk+kp)         & 
     173                                                                + zah1 * r1_e1u(ji-1,jj) * r1_e1u(ji-1,jj) * umask(ji-1,jj,jk+kp)  & 
     174                        &                                       )                                                                  ! bracket for halo 1 - halo 2 compatibility 
     175            END_3D 
    175176         END DO 
    176177         ! 
    177          DO jp = 0, 1                            ! j-k triads 
    178             DO kp = 0, 1 
    179                DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
    180                   ze3wr = 1.0_wp / e3w(ji,jj,jk+kp,Kmm) 
    181                   zbv   = e1e2v(ji,jj-jp) * e3v(ji,jj-jp,jk,Kmm) 
    182                   zah   = 0.25_wp * pahv(ji,jj-jp,jk) 
    183                   zslope_skew = triadj_g(ji,jj,jk,1-jp,kp) 
    184                   ! Subtract s-coordinate slope at t-points to give slope rel to s surfaces 
    185                   !    (do this by *adding* gradient of depth) 
    186                   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) 
    187                   zslope2 = zslope2 * zslope2 
    188                   ah_wslp2(ji,jj,jk+kp) = ah_wslp2(ji,jj,jk+kp) + zah * zbv * ze3wr * r1_e1e2t(ji,jj) * zslope2 
    189                   akz     (ji,jj,jk+kp) = akz     (ji,jj,jk+kp) + zah * r1_e2v(ji,jj-jp)     & 
    190                      &                                                      * r1_e2v(ji,jj-jp) * vmask(ji,jj-jp,jk+kp) 
    191                   ! 
    192                END_3D 
    193             END DO 
     178         DO kp = 0, 1                            ! j-k triads 
     179            DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     180               ze3wr = 1.0_wp / e3w(ji,jj,jk+kp,Kmm) 
     181               zbv   = e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
     182               zbv1   = e1e2v(ji,jj-1) * e3v(ji,jj-1,jk,Kmm) 
     183               zah   = 0.25_wp * pahv(ji,jj,jk) 
     184               zah1   = 0.25_wp * pahv(ji,jj-1,jk) 
     185               ! Subtract s-coordinate slope at t-points to give slope rel to s surfaces 
     186               !    (do this by *adding* gradient of depth) 
     187               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) 
     188               zslope2 = zslope2 * zslope2 
     189               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) 
     190               zslope21 = zslope21 * zslope21 
     191               ! round brackets added to fix the order of floating point operations 
     192               ! needed to ensure halo 1 - halo 2 compatibility 
     193               ah_wslp2(ji,jj,jk+kp) = ah_wslp2(ji,jj,jk+kp) + ( zah * zbv * ze3wr * r1_e1e2t(ji,jj) * zslope2                     & 
     194                        &                                      + zah1 * zbv1 * ze3wr * r1_e1e2t(ji,jj) * zslope21                  & 
     195                        &                                      )                                                                   ! bracket for halo 1 - halo 2 compatibility 
     196               akz     (ji,jj,jk+kp) = akz     (ji,jj,jk+kp) + ( zah * r1_e2v(ji,jj) * r1_e2v(ji,jj) * vmask(ji,jj,jk+kp)          & 
     197                        &                                      + zah1 * r1_e2v(ji,jj-1) * r1_e2v(ji,jj-1) * vmask(ji,jj-1,jk+kp)   & 
     198                        &                                      )                                                                   ! bracket for halo 1 - halo 2 compatibility 
     199            END_3D 
    194200         END DO 
    195201         ! 
     
    223229               IF( ln_tile ) CALL dom_tile( ntsi, ntsj, ntei, ntej, ktile = 0 ) 
    224230 
    225                zpsi_uw(:,:,:) = 0._wp 
    226                zpsi_vw(:,:,:) = 0._wp 
    227  
    228                DO jp = 0, 1 
     231                  zpsi_uw(:,:,:) = 0._wp 
     232                  zpsi_vw(:,:,:) = 0._wp 
     233 
    229234                  DO kp = 0, 1 
    230235                     DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 
    231                         zpsi_uw(ji,jj,jk+kp) = zpsi_uw(ji,jj,jk+kp) & 
    232                            & + 0.25_wp * aeiu(ji,jj,jk) * e2u(ji,jj) * triadi_g(ji+jp,jj,jk,1-jp,kp) 
    233                         zpsi_vw(ji,jj,jk+kp) = zpsi_vw(ji,jj,jk+kp) & 
    234                            & + 0.25_wp * aeiv(ji,jj,jk) * e1v(ji,jj) * triadj_g(ji,jj+jp,jk,1-jp,kp) 
     236                        ! round brackets added to fix the order of floating point operations 
     237                        ! needed to ensure halo 1 - halo 2 compatibility [ NOT TESTED ] 
     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,jj,jk,1,kp)        &  
     240                           &   + 0.25_wp * aeiu(ji,jj,jk) * e2u(ji,jj) * triadi_g(ji+1,jj,jk,0,kp)      & 
     241                           &   )                                                                        ! bracket for halo 1 - halo 2 compatibility 
     242                        zpsi_vw(ji,jj,jk+kp) = zpsi_vw(ji,jj,jk+kp)                                     & 
     243                           & + ( 0.25_wp * aeiv(ji,jj,jk) * e1v(ji,jj) * triadj_g(ji,jj,jk,1,kp)        & 
     244                           &   + 0.25_wp * aeiv(ji,jj,jk) * e1v(ji,jj) * triadj_g(ji,jj+1,jk,0,kp)      & 
     245                           &   )                                                                        ! bracket for halo 1 - halo 2 compatibility 
    235246                     END_3D 
    236247                  END DO 
    237                END DO 
    238248               CALL ldf_eiv_dia( zpsi_uw, zpsi_vw, Kmm ) 
    239249 
     
    291301            ! 
    292302            IF( ln_botmix_triad ) THEN 
    293                DO ip = 0, 1              !==  Horizontal & vertical fluxes 
    294                   DO kp = 0, 1 
    295                      DO_2D( 1, 0, 1, 0 ) 
    296                         ze1ur = r1_e1u(ji,jj) 
    297                         zdxt  = zdit(ji,jj,jk) * ze1ur 
    298                         ze3wr = 1._wp / e3w(ji+ip,jj,jk+kp,Kmm) 
    299                         zdzt  = zdkt3d(ji+ip,jj,kp) * ze3wr 
    300                         zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 
    301                         zslope_iso  = triadi  (ji+ip,jj,jk,1-ip,kp) 
    302                         ! 
    303                         zbu = 0.25_wp * e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
    304                         ! ln_botmix_triad is .T. don't mask zah for bottom half cells    !!gm ?????   ahu is masked.... 
    305                         zah = pahu(ji,jj,jk) 
    306                         zah_slp  = zah * zslope_iso 
    307                         IF( ln_ldfeiv )   zaei_slp = aeiu(ji,jj,jk) * zslope_skew 
    308                         zftu(ji   ,jj,jk   ) = zftu(ji   ,jj,jk   ) - ( zah * zdxt + (zah_slp - zaei_slp) * zdzt ) * zbu * ze1ur 
    309                         ztfw(ji+ip,jj,jk+kp) = ztfw(ji+ip,jj,jk+kp) - ( zah_slp + zaei_slp) * zdxt                 * zbu * ze3wr 
     303               DO kp = 0, 1              !==  Horizontal & vertical fluxes 
     304                  DO_2D( 1, 0, 1, 0 ) 
     305                     ze1ur = r1_e1u(ji,jj) 
     306                     zdxt  = zdit(ji,jj,jk) * ze1ur 
     307                     zdxt_ip1  = zdit(ji+1,jj,jk) * r1_e1u(ji+1,jj) 
     308                     ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     309                     ze3wr_ip1 = 1._wp / e3w(ji+1,jj,jk+kp,Kmm) 
     310                     zdzt  = zdkt3d(ji,jj,kp) * ze3wr 
     311                     zdzt_ip1  = zdkt3d(ji+1,jj,kp) * ze3wr_ip1 
     312                     ! 
     313                     zbu = 0.25_wp * e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
     314                     zbu_ip1 = 0.25_wp * e1e2u(ji+1,jj) * e3u(ji+1,jj,jk,Kmm) 
     315                     ! ln_botmix_triad is .T. don't mask zah for bottom half cells    !!gm ?????   ahu is masked.... 
     316                     zah = pahu(ji,jj,jk)  
     317                     zah_ip1 = pahu(ji+1,jj,jk)  
     318                     zah_slp  = zah * triadi(ji,jj,jk,1,kp) 
     319                     zah_slp_ip1  = zah_ip1 * triadi(ji+1,jj,jk,1,kp) 
     320                     zah_slp1  = zah * triadi(ji+1,jj,jk,0,kp) 
     321                     IF( ln_ldfeiv )   THEN 
     322                        zaei_slp = aeiu(ji,jj,jk) * triadi_g(ji,jj,jk,1,kp) 
     323                        zaei_slp_ip1 = aeiu(ji+1,jj,jk) * triadi_g(ji+1,jj,jk,1,kp) 
     324                        zaei_slp1 = aeiu(ji,jj,jk) * triadi_g(ji+1,jj,jk,0,kp) 
     325                     ENDIF 
     326                     ! round brackets added to fix the order of floating point operations 
     327                     ! needed to ensure halo 1 - halo 2 compatibility 
     328                     zftu(ji   ,jj,jk  ) =  zftu(ji   ,jj,jk )                                                               & 
     329                                         &    - ( ( zah * zdxt + ( zah_slp - zaei_slp ) * zdzt ) * zbu * ze1ur               & 
     330                                         &      + ( zah * zdxt + zah_slp1 * zdzt_ip1 - zaei_slp1 * zdzt_ip1 ) * zbu * ze1ur  & 
     331                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     332                     ztfw(ji+1,jj,jk+kp) =  ztfw(ji+1,jj,jk+kp)                                                              &  
     333                                         &    - ( (zah_slp_ip1 + zaei_slp_ip1) * zdxt_ip1 * zbu_ip1 * ze3wr_ip1              & 
     334                                         &      + ( zah_slp1 + zaei_slp1) * zdxt * zbu * ze3wr_ip1                           &  
     335                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
    310336                     END_2D 
    311337                  END DO 
    312                END DO 
    313338               ! 
    314                DO jp = 0, 1 
    315                   DO kp = 0, 1 
    316                      DO_2D( 1, 0, 1, 0 ) 
    317                         ze2vr = r1_e2v(ji,jj) 
    318                         zdyt  = zdjt(ji,jj,jk) * ze2vr 
    319                         ze3wr = 1._wp / e3w(ji,jj+jp,jk+kp,Kmm) 
    320                         zdzt  = zdkt3d(ji,jj+jp,kp) * ze3wr 
    321                         zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 
    322                         zslope_iso  = triadj(ji,jj+jp,jk,1-jp,kp) 
    323                         zbv = 0.25_wp * e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
    324                         ! ln_botmix_triad is .T. don't mask zah for bottom half cells    !!gm ?????  ahv is masked... 
    325                         zah = pahv(ji,jj,jk) 
    326                         zah_slp = zah * zslope_iso 
    327                         IF( ln_ldfeiv )   zaei_slp = aeiv(ji,jj,jk) * zslope_skew 
    328                         zftv(ji,jj   ,jk   ) = zftv(ji,jj   ,jk   ) - ( zah * zdyt + (zah_slp - zaei_slp) * zdzt ) * zbv * ze2vr 
    329                         ztfw(ji,jj+jp,jk+kp) = ztfw(ji,jj+jp,jk+kp) - ( zah_slp + zaei_slp ) * zdyt                * zbv * ze3wr 
    330                      END_2D 
    331                   END DO 
     339               DO kp = 0, 1 
     340                  DO_2D( 1, 0, 1, 0 ) 
     341                     ze2vr = r1_e2v(ji,jj) 
     342                     zdyt  = zdjt(ji,jj,jk) * ze2vr 
     343                     zdyt_jp1  = zdjt(ji,jj+1,jk) * r1_e2v(ji,jj+1) 
     344                     ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     345                     ze3wr_jp1 = 1._wp / e3w(ji,jj+1,jk+kp,Kmm) 
     346                     zdzt  = zdkt3d(ji,jj,kp) * ze3wr 
     347                     zdzt_jp1  = zdkt3d(ji,jj+1,kp) * ze3wr_jp1 
     348                     zbv = 0.25_wp * e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
     349                     zbv_jp1 = 0.25_wp * e1e2v(ji,jj+1) * e3v(ji,jj+1,jk,Kmm) 
     350                     ! ln_botmix_triad is .T. don't mask zah for bottom half cells    !!gm ?????   ahu is masked.... 
     351                     zah = pahv(ji,jj,jk)          ! pahv(ji,jj+jp,jk)  ???? 
     352                     zah_jp1 = pahv(ji,jj+1,jk)  
     353                     zah_slp = zah * triadj(ji,jj,jk,1,kp) 
     354                     zah_slp1 = zah * triadj(ji,jj+1,jk,0,kp) 
     355                     zah_slp_jp1 = zah_jp1 * triadj(ji,jj+1,jk,1,kp) 
     356                     IF( ln_ldfeiv )   THEN 
     357                        zaei_slp = aeiv(ji,jj,jk) * triadj_g(ji,jj,jk,1,kp)   
     358                        zaei_slp_jp1 = aeiv(ji,jj+1,jk) * triadj_g(ji,jj+1,jk,1,kp)   
     359                        zaei_slp1 = aeiv(ji,jj,jk) * triadj_g(ji,jj+1,jk,0,kp) 
     360                     ENDIF 
     361                     ! round brackets added to fix the order of floating point operations 
     362                     ! needed to ensure halo 1 - halo 2 compatibility 
     363                     zftv(ji,jj  ,jk   ) =  zftv(ji,jj  ,jk   )                                                              & 
     364                                         &    - ( ( zah * zdyt + ( zah_slp - zaei_slp ) * zdzt ) * zbv * ze2vr               & 
     365                                         &      + ( zah * zdyt + zah_slp1 * zdzt_jp1 - zaei_slp1 * zdzt_jp1 ) * zbv * ze2vr  & 
     366                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     367                     ztfw(ji,jj+1,jk+kp) =  ztfw(ji,jj+1,jk+kp)                                                              & 
     368                                         &    - ( ( zah_slp_jp1 + zaei_slp_jp1) * zdyt_jp1 * zbv_jp1 * ze3wr_jp1             & 
     369                                         &      + ( zah_slp1 + zaei_slp1) * zdyt * zbv * ze3wr_jp1                           & 
     370                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     371                  END_2D 
    332372               END DO 
    333373               ! 
    334374            ELSE 
    335375               ! 
    336                DO ip = 0, 1               !==  Horizontal & vertical fluxes 
    337                   DO kp = 0, 1 
    338                      DO_2D( 1, 0, 1, 0 ) 
    339                         ze1ur = r1_e1u(ji,jj) 
    340                         zdxt  = zdit(ji,jj,jk) * ze1ur 
    341                         ze3wr = 1._wp / e3w(ji+ip,jj,jk+kp,Kmm) 
    342                         zdzt  = zdkt3d(ji+ip,jj,kp) * ze3wr 
    343                         zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 
    344                         zslope_iso  = triadi(ji+ip,jj,jk,1-ip,kp) 
    345                         ! 
    346                         zbu = 0.25_wp * e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
    347                         ! ln_botmix_triad is .F. mask zah for bottom half cells 
    348                         zah = pahu(ji,jj,jk) * umask(ji,jj,jk+kp)         ! pahu(ji+ip,jj,jk)   ===>>  ???? 
    349                         zah_slp  = zah * zslope_iso 
    350                         IF( ln_ldfeiv )   zaei_slp = aeiu(ji,jj,jk) * zslope_skew        ! aeit(ji+ip,jj,jk)*zslope_skew 
    351                         zftu(ji   ,jj,jk   ) = zftu(ji   ,jj,jk   ) - ( zah * zdxt + (zah_slp - zaei_slp) * zdzt ) * zbu * ze1ur 
    352                         ztfw(ji+ip,jj,jk+kp) = ztfw(ji+ip,jj,jk+kp) - (zah_slp + zaei_slp) * zdxt * zbu * ze3wr 
    353                      END_2D 
    354                   END DO 
     376               DO kp = 0, 1 
     377                  DO_2D( 1, 0, 1, 0 ) 
     378                     ze1ur = r1_e1u(ji,jj) 
     379                     zdxt  = zdit(ji,jj,jk) * ze1ur 
     380                     zdxt_ip1  = zdit(ji+1,jj,jk) * r1_e1u(ji+1,jj) 
     381                     ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     382                     ze3wr_ip1 = 1._wp / e3w(ji+1,jj,jk+kp,Kmm) 
     383                     zdzt  = zdkt3d(ji,jj,kp) * ze3wr 
     384                     zdzt_ip1  = zdkt3d(ji+1,jj,kp) * ze3wr_ip1 
     385                     ! 
     386                     zbu = 0.25_wp * e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
     387                     zbu_ip1 = 0.25_wp * e1e2u(ji+1,jj) * e3u(ji+1,jj,jk,Kmm) 
     388                     ! ln_botmix_triad is .F. mask zah for bottom half cells 
     389                     zah = pahu(ji,jj,jk) * umask(ji,jj,jk+kp)         ! pahu(ji+ip,jj,jk)   ===>>  ???? 
     390                     zah_ip1 = pahu(ji+1,jj,jk) * umask(ji+1,jj,jk+kp)  
     391                     zah_slp  = zah * triadi(ji,jj,jk,1,kp) 
     392                     zah_slp_ip1  = zah_ip1 * triadi(ji+1,jj,jk,1,kp) 
     393                     zah_slp1  = zah * triadi(ji+1,jj,jk,0,kp) 
     394                     IF( ln_ldfeiv )   THEN 
     395                        zaei_slp = aeiu(ji,jj,jk) * triadi_g(ji,jj,jk,1,kp)  
     396                        zaei_slp_ip1 = aeiu(ji+1,jj,jk) * triadi_g(ji+1,jj,jk,1,kp) 
     397                        zaei_slp1 = aeiu(ji,jj,jk) * triadi_g(ji+1,jj,jk,0,kp) 
     398                     ENDIF 
     399                     ! round brackets added to fix the order of floating point operations 
     400                     ! needed to ensure halo 1 - halo 2 compatibility 
     401                     zftu(ji   ,jj,jk  ) =  zftu(ji   ,jj,jk )                                                               & 
     402                                         &    - ( ( zah * zdxt + ( zah_slp - zaei_slp ) * zdzt ) * zbu * ze1ur               & 
     403                                         &      + ( zah * zdxt + zah_slp1 * zdzt_ip1 - zaei_slp1 * zdzt_ip1 ) * zbu * ze1ur  & 
     404                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     405                     ztfw(ji+1,jj,jk+kp) =  ztfw(ji+1,jj,jk+kp)                                                              &  
     406                                         &    - ( (zah_slp_ip1 + zaei_slp_ip1) * zdxt_ip1 * zbu_ip1 * ze3wr_ip1              & 
     407                                         &      + ( zah_slp1 + zaei_slp1) * zdxt * zbu * ze3wr_ip1                           & 
     408                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     409                  END_2D 
    355410               END DO 
    356411               ! 
    357                DO jp = 0, 1 
    358                   DO kp = 0, 1 
    359                      DO_2D( 1, 0, 1, 0 ) 
    360                         ze2vr = r1_e2v(ji,jj) 
    361                         zdyt  = zdjt(ji,jj,jk) * ze2vr 
    362                         ze3wr = 1._wp / e3w(ji,jj+jp,jk+kp,Kmm) 
    363                         zdzt  = zdkt3d(ji,jj+jp,kp) * ze3wr 
    364                         zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 
    365                         zslope_iso  = triadj(ji,jj+jp,jk,1-jp,kp) 
    366                         zbv = 0.25_wp * e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
    367                         ! ln_botmix_triad is .F. mask zah for bottom half cells 
    368                         zah = pahv(ji,jj,jk) * vmask(ji,jj,jk+kp)         ! pahv(ji,jj+jp,jk)  ???? 
    369                         zah_slp = zah * zslope_iso 
    370                         IF( ln_ldfeiv )   zaei_slp = aeiv(ji,jj,jk) * zslope_skew        ! aeit(ji,jj+jp,jk)*zslope_skew 
    371                         zftv(ji,jj,jk) = zftv(ji,jj,jk) - ( zah * zdyt + (zah_slp - zaei_slp) * zdzt ) * zbv * ze2vr 
    372                         ztfw(ji,jj+jp,jk+kp) = ztfw(ji,jj+jp,jk+kp) - (zah_slp + zaei_slp) * zdyt * zbv * ze3wr 
    373                      END_2D 
    374                   END DO 
     412               DO kp = 0, 1 
     413                  DO_2D( 1, 0, 1, 0 ) 
     414                     ze2vr = r1_e2v(ji,jj) 
     415                     zdyt  = zdjt(ji,jj,jk) * ze2vr 
     416                     zdyt_jp1  = zdjt(ji,jj+1,jk) * r1_e2v(ji,jj+1) 
     417                     ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     418                     ze3wr_jp1 = 1._wp / e3w(ji,jj+1,jk+kp,Kmm) 
     419                     zdzt  = zdkt3d(ji,jj,kp) * ze3wr 
     420                     zdzt_jp1  = zdkt3d(ji,jj+1,kp) * ze3wr_jp1 
     421                     zbv = 0.25_wp * e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
     422                     zbv_jp1 = 0.25_wp * e1e2v(ji,jj+1) * e3v(ji,jj+1,jk,Kmm) 
     423                     ! ln_botmix_triad is .F. mask zah for bottom half cells 
     424                     zah = pahv(ji,jj,jk) * vmask(ji,jj,jk+kp)         ! pahv(ji,jj+jp,jk)  ???? 
     425                     zah_jp1 = pahv(ji,jj+1,jk) * vmask(ji,jj+1,jk+kp) 
     426                     zah_slp = zah * triadj(ji,jj,jk,1,kp) 
     427                     zah_slp1 = zah * triadj(ji,jj+1,jk,0,kp) 
     428                     zah_slp_jp1 = zah_jp1 * triadj(ji,jj+1,jk,1,kp) 
     429                     IF( ln_ldfeiv )   THEN 
     430                        zaei_slp = aeiv(ji,jj,jk) * triadj_g(ji,jj,jk,1,kp)  
     431                        zaei_slp_jp1 = aeiv(ji,jj+1,jk) * triadj_g(ji,jj+1,jk,1,kp) 
     432                        zaei_slp1 = aeiv(ji,jj,jk) * triadj_g(ji,jj+1,jk,0,kp) 
     433                     ENDIF 
     434                     ! round brackets added to fix the order of floating point operations 
     435                     ! needed to ensure halo 1 - halo 2 compatibility 
     436                     zftv(ji,jj  ,jk   ) =  zftv(ji,jj  ,jk   )                                                              & 
     437                                         &    - ( ( zah * zdyt + ( zah_slp - zaei_slp ) * zdzt ) * zbv * ze2vr               & 
     438                                         &      + ( zah * zdyt + zah_slp1 * zdzt_jp1 - zaei_slp1 * zdzt_jp1 ) * zbv * ze2vr  & 
     439                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     440                     ztfw(ji,jj+1,jk+kp) =  ztfw(ji,jj+1,jk+kp)                                                              & 
     441                                         &    - ( ( zah_slp_jp1 + zaei_slp_jp1) * zdyt_jp1 * zbv_jp1 * ze3wr_jp1             & 
     442                                         &      + ( zah_slp1 + zaei_slp1) * zdyt * zbv * ze3wr_jp1                           & 
     443                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     444                  END_2D 
    375445               END DO 
    376446            ENDIF 
    377447            !                             !==  horizontal divergence and add to the general trend  ==! 
    378448            DO_2D( 0, 0, 0, 0 ) 
    379                pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)    & 
    380                   &                       + zsign * (  zftu(ji-1,jj  ,jk) - zftu(ji,jj,jk)       & 
    381                   &                                           + zftv(ji,jj-1,jk) - zftv(ji,jj,jk)   )   & 
    382                   &                                        / (  e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm)  ) 
     449               ! round brackets added to fix the order of floating point operations 
     450               ! needed to ensure halo 1 - halo 2 compatibility 
     451               pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)                                                & 
     452                  &                       + zsign * ( ( zftu(ji-1,jj  ,jk) - zftu(ji,jj,jk)             & 
     453                  &                                   )                                                 & ! bracket for halo 1 - halo 2 compatibility 
     454                  &                                 + ( zftv(ji,jj-1,jk) - zftv(ji,jj,jk)               & 
     455                  &                                   )                                                 & ! bracket for halo 1 - halo 2 compatibility 
     456                  &                                 ) / (  e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm)  ) 
    383457            END_2D 
    384458            ! 
Note: See TracChangeset for help on using the changeset viewer.