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

Changeset 14820


Ignore:
Timestamp:
2021-05-10T10:26:13+02:00 (3 years ago)
Author:
francesca
Message:

merge ticket2607_r14608_halo1_halo2_compatibility into trunk

Location:
NEMO/trunk/src/OCE
Files:
10 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk/src/OCE/DOM/domqco.F90

    r14433 r14820  
    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/trunk/src/OCE/DYN/divhor.F90

    r13558 r14820  
    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/trunk/src/OCE/DYN/dynadv_ubs.F90

    r14433 r14820  
    109109         !             
    110110         DO_2D( 0, 0, 0, 0 )                       ! laplacian 
    111             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) 
    112             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) 
     111            ! round brackets added to fix the order of floating point operations 
     112            ! needed to ensure halo 1 - halo 2 compatibility 
     113            zlu_uu(ji,jj,jk,1) = ( (puu (ji+1,jj  ,jk,Kbb) - puu (ji  ,jj  ,jk,Kbb)) +                         & 
     114               &                   (puu (ji-1,jj  ,jk,Kbb) - puu (ji  ,jj  ,jk,Kbb)) ) * umask(ji  ,jj  ,jk) 
     115            zlv_vv(ji,jj,jk,1) = ( (pvv (ji  ,jj+1,jk,Kbb) - pvv (ji  ,jj  ,jk,Kbb)) +                         & 
     116               &                   (pvv (ji  ,jj-1,jk,Kbb) - pvv (ji  ,jj  ,jk,Kbb)) ) * vmask(ji  ,jj  ,jk) 
    113117            zlu_uv(ji,jj,jk,1) = ( puu (ji  ,jj+1,jk,Kbb) - puu (ji  ,jj  ,jk,Kbb) ) * fmask(ji  ,jj  ,jk)   & 
    114118               &               - ( puu (ji  ,jj  ,jk,Kbb) - puu (ji  ,jj-1,jk,Kbb) ) * fmask(ji  ,jj-1,jk) 
     
    116120               &               - ( pvv (ji  ,jj  ,jk,Kbb) - pvv (ji-1,jj  ,jk,Kbb) ) * fmask(ji-1,jj  ,jk) 
    117121            ! 
    118             zlu_uu(ji,jj,jk,2) = ( zfu(ji+1,jj  ,jk) - 2.*zfu(ji,jj,jk) + zfu(ji-1,jj  ,jk) ) * umask(ji,jj,jk) 
    119             zlv_vv(ji,jj,jk,2) = ( zfv(ji  ,jj+1,jk) - 2.*zfv(ji,jj,jk) + zfv(ji  ,jj-1,jk) ) * vmask(ji,jj,jk) 
     122            zlu_uu(ji,jj,jk,2) = ( (zfu(ji+1,jj  ,jk) - zfu(ji  ,jj  ,jk)) +                                   & 
     123               &                   (zfu(ji-1,jj  ,jk) - zfu(ji  ,jj  ,jk)) ) * umask(ji  ,jj  ,jk) 
     124            zlv_vv(ji,jj,jk,2) = ( (zfv(ji  ,jj+1,jk) - zfv(ji  ,jj  ,jk)) +                                   & 
     125               &                   (zfv(ji  ,jj-1,jk) - zfv(ji  ,jj  ,jk)) ) * vmask(ji  ,jj  ,jk) 
    120126            zlu_uv(ji,jj,jk,2) = ( zfu(ji  ,jj+1,jk) - zfu(ji  ,jj  ,jk) ) * fmask(ji  ,jj  ,jk)   & 
    121127               &               - ( zfu(ji  ,jj  ,jk) - zfu(ji  ,jj-1,jk) ) * fmask(ji  ,jj-1,jk) 
     
    124130         END_2D 
    125131      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   ) 
     132      CALL lbc_lnk( 'dynadv_ubs', zlu_uu(:,:,:,1), 'U', -1.0_wp , zlu_uv(:,:,:,1), 'U', -1.0_wp,  & 
     133         &                        zlu_uu(:,:,:,2), 'U', -1.0_wp , zlu_uv(:,:,:,2), 'U', -1.0_wp,  & 
     134         &                        zlv_vv(:,:,:,1), 'V', -1.0_wp , zlv_vu(:,:,:,1), 'V', -1.0_wp,  & 
     135         &                        zlv_vv(:,:,:,2), 'V', -1.0_wp , zlv_vu(:,:,:,2), 'V', -1.0_wp   ) 
    130136      ! 
    131137      !                                      ! ====================== ! 
  • NEMO/trunk/src/OCE/DYN/dynkeg.F90

    r13497 r14820  
    110110      CASE ( nkeg_HW )                          !--  Hollingsworth scheme  --! 
    111111         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
     112            ! round brackets added to fix the order of floating point operations 
     113            ! needed to ensure halo 1 - halo 2 compatibility 
    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) )  & 
    114                &   +     ( puu(ji-1,jj-1,jk,Kmm) + puu(ji-1,jj+1,jk,Kmm) ) * ( puu(ji-1,jj-1,jk,Kmm) + puu(ji-1,jj+1,jk,Kmm) )   & 
    115                &   +     ( puu(ji  ,jj-1,jk,Kmm) + puu(ji  ,jj+1,jk,Kmm) ) * ( puu(ji  ,jj-1,jk,Kmm) + puu(ji  ,jj+1,jk,Kmm) ) 
     116               &   +     ( ( puu(ji-1,jj-1,jk,Kmm) + puu(ji-1,jj+1,jk,Kmm) ) * ( puu(ji-1,jj-1,jk,Kmm) + puu(ji-1,jj+1,jk,Kmm) )   & 
     117               &   +       ( puu(ji  ,jj-1,jk,Kmm) + puu(ji  ,jj+1,jk,Kmm) ) * ( puu(ji  ,jj-1,jk,Kmm) + puu(ji  ,jj+1,jk,Kmm) )   & 
     118               &         )                                                               ! bracket for halo 1 - halo 2 compatibility 
    116119               ! 
    117120            zv = 8._wp * ( pvv(ji  ,jj-1,jk,Kmm) * pvv(ji  ,jj-1,jk,Kmm)    & 
    118121               &         + pvv(ji  ,jj  ,jk,Kmm) * pvv(ji  ,jj  ,jk,Kmm) )  & 
    119                &  +      ( pvv(ji-1,jj-1,jk,Kmm) + pvv(ji+1,jj-1,jk,Kmm) ) * ( pvv(ji-1,jj-1,jk,Kmm) + pvv(ji+1,jj-1,jk,Kmm) )   & 
    120                &  +      ( pvv(ji-1,jj  ,jk,Kmm) + pvv(ji+1,jj  ,jk,Kmm) ) * ( pvv(ji-1,jj  ,jk,Kmm) + pvv(ji+1,jj  ,jk,Kmm) ) 
     122               &  +      ( ( pvv(ji-1,jj-1,jk,Kmm) + pvv(ji+1,jj-1,jk,Kmm) ) * ( pvv(ji-1,jj-1,jk,Kmm) + pvv(ji+1,jj-1,jk,Kmm) )  & 
     123               &  +        ( pvv(ji-1,jj  ,jk,Kmm) + pvv(ji+1,jj  ,jk,Kmm) ) * ( pvv(ji-1,jj  ,jk,Kmm) + pvv(ji+1,jj  ,jk,Kmm) )  &  
     124               &         )                                                               ! bracket for halo 1 - halo 2 compatibility 
    121125            zhke(ji,jj,jk) = r1_48 * ( zv + zu ) 
    122126         END_3D 
  • NEMO/trunk/src/OCE/DYN/dynvor.F90

    r14433 r14820  
    632632         CASE ( 0 )                                   ! original formulation  (masked averaging of e3t divided by 4) 
    633633            DO_2D( 1, 0, 1, 0 ) 
    634                ze3f = (  e3t(ji  ,jj+1,jk,Kmm)*tmask(ji  ,jj+1,jk)   & 
    635                   &    + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk)   & 
    636                   &    + e3t(ji  ,jj  ,jk,Kmm)*tmask(ji  ,jj  ,jk)   & 
    637                   &    + e3t(ji+1,jj  ,jk,Kmm)*tmask(ji+1,jj  ,jk)  ) 
     634               ! round brackets added to fix the order of floating point operations 
     635               ! needed to ensure halo 1 - halo 2 compatibility 
     636               ze3f = (  (e3t(ji  ,jj+1,jk,Kmm)*tmask(ji  ,jj+1,jk)    & 
     637                  &    +  e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk))   & 
     638                  &    + (e3t(ji  ,jj  ,jk,Kmm)*tmask(ji  ,jj  ,jk)    & 
     639                  &    +  e3t(ji+1,jj  ,jk,Kmm)*tmask(ji+1,jj  ,jk))  ) 
    638640               IF( ze3f /= 0._wp ) THEN   ;   z1_e3f(ji,jj) = 4._wp / ze3f 
    639641               ELSE                       ;   z1_e3f(ji,jj) = 0._wp 
     
    642644         CASE ( 1 )                                   ! new formulation  (masked averaging of e3t divided by the sum of mask) 
    643645            DO_2D( 1, 0, 1, 0 ) 
    644                ze3f = (  e3t(ji  ,jj+1,jk,Kmm)*tmask(ji  ,jj+1,jk)   & 
    645                   &    + e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk)   & 
    646                   &    + e3t(ji  ,jj  ,jk,Kmm)*tmask(ji  ,jj  ,jk)   & 
    647                   &    + e3t(ji+1,jj  ,jk,Kmm)*tmask(ji+1,jj  ,jk)  ) 
     646               ! round brackets added to fix the order of floating point operations 
     647               ! needed to ensure halo 1 - halo 2 compatibility 
     648               ze3f = (  (e3t(ji  ,jj+1,jk,Kmm)*tmask(ji  ,jj+1,jk)    & 
     649                  &    +  e3t(ji+1,jj+1,jk,Kmm)*tmask(ji+1,jj+1,jk))   & 
     650                  &    + (e3t(ji  ,jj  ,jk,Kmm)*tmask(ji  ,jj  ,jk)    & 
     651                  &    +  e3t(ji+1,jj  ,jk,Kmm)*tmask(ji+1,jj  ,jk))  ) 
    648652               zmsk = (                    tmask(ji,jj+1,jk) +                     tmask(ji+1,jj+1,jk)   & 
    649653                  &                      + tmask(ji,jj  ,jk) +                     tmask(ji+1,jj  ,jk)  ) 
     
    678682         CASE ( np_CRV )                           !* Coriolis + relative vorticity 
    679683            DO_2D( 1, 0, 1, 0 ) 
    680                zwz(ji,jj,jk) = (  ff_f(ji,jj) + (  e2v(ji+1,jj  ) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk)      & 
    681                   &                              - e1u(ji  ,jj+1) * pu(ji,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk)  )   & 
    682                   &                           * r1_e1e2f(ji,jj)   ) * z1_e3f(ji,jj) 
     684               ! round brackets added to fix the order of floating point operations 
     685               ! needed to ensure halo 1 - halo 2 compatibility 
     686               zwz(ji,jj,jk) = (  ff_f(ji,jj) + ( ( e2v(ji+1,jj  ) * pv(ji+1,jj,jk) - e2v(ji,jj) * pv(ji,jj,jk)      & 
     687                  &                               )                                                                  & ! bracket for halo 1 - halo 2 compatibility 
     688                  &                             - ( e1u(ji  ,jj+1) * pu(ji,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk)      &  
     689                  &                               )                                                                  & ! bracket for halo 1 - halo 2 compatibility 
     690                  &                             ) * r1_e1e2f(ji,jj)   ) * z1_e3f(ji,jj) 
    683691            END_2D 
    684692            IF( ln_dynvor_msk ) THEN                     ! mask the relative vorticity 
     
    781789         CASE ( np_RVO )                           !* relative vorticity 
    782790            DO_2D( 1, 0, 1, 0 ) 
    783                zwz(ji,jj,jk) = (  e2v(ji+1,jj  ) * pv(ji+1,jj  ,jk) - e2v(ji,jj) * pv(ji,jj,jk)    & 
    784                   &             - e1u(ji  ,jj+1) * pu(ji  ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk)  ) & 
     791               ! round brackets added to fix the order of floating point operations 
     792               ! needed to ensure halo 1 - halo 2 compatibility 
     793               zwz(ji,jj,jk) = (  (e2v(ji+1,jj  ) * pv(ji+1,jj  ,jk) - e2v(ji,jj) * pv(ji,jj,jk))    & 
     794                  &             - (e1u(ji  ,jj+1) * pu(ji  ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk))  ) & 
    785795                  &          * r1_e1e2f(ji,jj) 
    786796            END_2D 
     
    797807         CASE ( np_CRV )                           !* Coriolis + relative vorticity 
    798808            DO_2D( 1, 0, 1, 0 ) 
    799                zwz(ji,jj,jk) = (  ff_f(ji,jj) + (  e2v(ji+1,jj  ) * pv(ji+1,jj  ,jk) - e2v(ji,jj) * pv(ji,jj,jk)    & 
    800                   &                              - e1u(ji  ,jj+1) * pu(ji  ,jj+1,jk) + e1u(ji,jj) * pu(ji,jj,jk)  ) & 
     809               ! round brackets added to fix the order of floating point operations 
     810               ! needed to ensure halo 1 - halo 2 compatibility 
     811               zwz(ji,jj,jk) = (  ff_f(ji,jj) + (  (e2v(ji+1,jj  ) * pv(ji+1,jj  ,jk) - e2v(ji,jj) * pv(ji,jj,jk))    & 
     812                  &                              - (e1u(ji  ,jj+1) * pu(ji  ,jj+1,jk) - e1u(ji,jj) * pu(ji,jj,jk))  ) & 
    801813                  &                         * r1_e1e2f(ji,jj)    ) 
    802814            END_2D 
  • NEMO/trunk/src/OCE/TRA/traadv_fct.F90

    r14433 r14820  
    238238               END_2D 
    239239            END DO 
    240             CALL lbc_lnk( 'traadv_fct', zltu, 'T', 1.0_wp , zltv, 'T', 1.0_wp )   ! Lateral boundary cond. (unchanged sgn) 
     240            ! NOTE [ comm_cleanup ] : need to change sign to ensure halo 1 - halo 2 compatibility 
     241            CALL lbc_lnk( 'traadv_fct', zltu, 'T', -1.0_wp , zltv, 'T', -1.0_wp )   ! Lateral boundary cond. (unchanged sgn) 
    241242            ! 
    242243            DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 
     
    244245               zC2t_v = pt(ji,jj,jk,jn,Kmm) + pt(ji  ,jj+1,jk,jn,Kmm) 
    245246               !                                                        ! C4 minus upstream advective fluxes 
    246                zwx(ji,jj,jk) =  0.5_wp * pU(ji,jj,jk) * ( zC2t_u + zltu(ji,jj,jk) - zltu(ji+1,jj,jk) ) - zwx(ji,jj,jk) 
    247                zwy(ji,jj,jk) =  0.5_wp * pV(ji,jj,jk) * ( zC2t_v + zltv(ji,jj,jk) - zltv(ji,jj+1,jk) ) - zwy(ji,jj,jk) 
     247               ! round brackets added to fix the order of floating point operations 
     248               ! needed to ensure halo 1 - halo 2 compatibility 
     249               zwx(ji,jj,jk) =  0.5_wp * pU(ji,jj,jk) * ( zC2t_u + ( zltu(ji,jj,jk) - zltu(ji+1,jj,jk)   & 
     250                             &                                     )                                     & ! bracket for halo 1 - halo 2 compatibility 
     251                             &                          ) - zwx(ji,jj,jk) 
     252               zwy(ji,jj,jk) =  0.5_wp * pV(ji,jj,jk) * ( zC2t_v + ( zltv(ji,jj,jk) - zltv(ji,jj+1,jk)   & 
     253                             &                                     )                                     & ! bracket for halo 1 - halo 2 compatibility 
     254                             &                          ) - zwy(ji,jj,jk) 
    248255            END_3D 
    249256            IF (nn_hls.EQ.2) CALL lbc_lnk( 'traadv_fct', zwx, 'U', -1.0_wp, zwy, 'V', -1.0_wp )   ! Lateral boundary cond. (unchanged sgn) 
  • NEMO/trunk/src/OCE/TRA/trabbl.F90

    r14433 r14820  
    141141         IF( ntile == 0 .OR. ntile == nijtile ) THEN                       ! Do only on the last tile 
    142142            ! lateral boundary conditions ; just need for outputs 
    143             CALL lbc_lnk( 'trabbl', utr_bbl, 'U', 1.0_wp , vtr_bbl, 'V', 1.0_wp ) 
     143            CALL lbc_lnk( 'trabbl', utr_bbl, 'U', -1.0_wp , vtr_bbl, 'V', -1.0_wp ) 
    144144            CALL iom_put( "uoce_bbl", utr_bbl )  ! bbl i-transport 
    145145            CALL iom_put( "voce_bbl", vtr_bbl )  ! bbl j-transport 
  • NEMO/trunk/src/OCE/TRA/traldf_iso.F90

    r14072 r14820  
    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/trunk/src/OCE/TRA/traldf_lap_blp.F90

    r14215 r14820  
    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/trunk/src/OCE/TRA/traldf_triad.F90

    r14215 r14820  
    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 
     
    289299            ! 
    290300            zaei_slp = 0._wp 
     301            zaei_slp_ip1 = 0._wp 
     302            zaei_slp_jp1 = 0._wp 
     303            zaei_slp1 = 0._wp 
    291304            ! 
    292305            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 
     306               DO kp = 0, 1              !==  Horizontal & vertical fluxes 
     307                  DO_2D( 1, 0, 1, 0 ) 
     308                     ze1ur = r1_e1u(ji,jj) 
     309                     zdxt  = zdit(ji,jj,jk) * ze1ur 
     310                     zdxt_ip1  = zdit(ji+1,jj,jk) * r1_e1u(ji+1,jj) 
     311                     ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     312                     ze3wr_ip1 = 1._wp / e3w(ji+1,jj,jk+kp,Kmm) 
     313                     zdzt  = zdkt3d(ji,jj,kp) * ze3wr 
     314                     zdzt_ip1  = zdkt3d(ji+1,jj,kp) * ze3wr_ip1 
     315                     ! 
     316                     zbu = 0.25_wp * e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
     317                     zbu_ip1 = 0.25_wp * e1e2u(ji+1,jj) * e3u(ji+1,jj,jk,Kmm) 
     318                     ! ln_botmix_triad is .T. don't mask zah for bottom half cells    !!gm ?????   ahu is masked.... 
     319                     zah = pahu(ji,jj,jk)  
     320                     zah_ip1 = pahu(ji+1,jj,jk)  
     321                     zah_slp  = zah * triadi(ji,jj,jk,1,kp) 
     322                     zah_slp_ip1  = zah_ip1 * triadi(ji+1,jj,jk,1,kp) 
     323                     zah_slp1  = zah * triadi(ji+1,jj,jk,0,kp) 
     324                     IF( ln_ldfeiv )   THEN 
     325                        zaei_slp = aeiu(ji,jj,jk) * triadi_g(ji,jj,jk,1,kp) 
     326                        zaei_slp_ip1 = aeiu(ji+1,jj,jk) * triadi_g(ji+1,jj,jk,1,kp) 
     327                        zaei_slp1 = aeiu(ji,jj,jk) * triadi_g(ji+1,jj,jk,0,kp) 
     328                     ENDIF 
     329                     ! round brackets added to fix the order of floating point operations 
     330                     ! needed to ensure halo 1 - halo 2 compatibility 
     331                     zftu(ji   ,jj,jk  ) =  zftu(ji   ,jj,jk )                                                               & 
     332                                         &    - ( ( zah * zdxt + ( zah_slp - zaei_slp ) * zdzt ) * zbu * ze1ur               & 
     333                                         &      + ( zah * zdxt + zah_slp1 * zdzt_ip1 - zaei_slp1 * zdzt_ip1 ) * zbu * ze1ur  & 
     334                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     335                     ztfw(ji+1,jj,jk+kp) =  ztfw(ji+1,jj,jk+kp)                                                              &  
     336                                         &    - ( (zah_slp_ip1 + zaei_slp_ip1) * zdxt_ip1 * zbu_ip1 * ze3wr_ip1              & 
     337                                         &      + ( zah_slp1 + zaei_slp1) * zdxt * zbu * ze3wr_ip1                           &  
     338                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
    310339                     END_2D 
    311340                  END DO 
    312                END DO 
    313341               ! 
    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 
     342               DO kp = 0, 1 
     343                  DO_2D( 1, 0, 1, 0 ) 
     344                     ze2vr = r1_e2v(ji,jj) 
     345                     zdyt  = zdjt(ji,jj,jk) * ze2vr 
     346                     zdyt_jp1  = zdjt(ji,jj+1,jk) * r1_e2v(ji,jj+1) 
     347                     ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     348                     ze3wr_jp1 = 1._wp / e3w(ji,jj+1,jk+kp,Kmm) 
     349                     zdzt  = zdkt3d(ji,jj,kp) * ze3wr 
     350                     zdzt_jp1  = zdkt3d(ji,jj+1,kp) * ze3wr_jp1 
     351                     zbv = 0.25_wp * e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
     352                     zbv_jp1 = 0.25_wp * e1e2v(ji,jj+1) * e3v(ji,jj+1,jk,Kmm) 
     353                     ! ln_botmix_triad is .T. don't mask zah for bottom half cells    !!gm ?????   ahu is masked.... 
     354                     zah = pahv(ji,jj,jk)          ! pahv(ji,jj+jp,jk)  ???? 
     355                     zah_jp1 = pahv(ji,jj+1,jk)  
     356                     zah_slp = zah * triadj(ji,jj,jk,1,kp) 
     357                     zah_slp1 = zah * triadj(ji,jj+1,jk,0,kp) 
     358                     zah_slp_jp1 = zah_jp1 * triadj(ji,jj+1,jk,1,kp) 
     359                     IF( ln_ldfeiv )   THEN 
     360                        zaei_slp = aeiv(ji,jj,jk) * triadj_g(ji,jj,jk,1,kp)   
     361                        zaei_slp_jp1 = aeiv(ji,jj+1,jk) * triadj_g(ji,jj+1,jk,1,kp)   
     362                        zaei_slp1 = aeiv(ji,jj,jk) * triadj_g(ji,jj+1,jk,0,kp) 
     363                     ENDIF 
     364                     ! round brackets added to fix the order of floating point operations 
     365                     ! needed to ensure halo 1 - halo 2 compatibility 
     366                     zftv(ji,jj  ,jk   ) =  zftv(ji,jj  ,jk   )                                                              & 
     367                                         &    - ( ( zah * zdyt + ( zah_slp - zaei_slp ) * zdzt ) * zbv * ze2vr               & 
     368                                         &      + ( zah * zdyt + zah_slp1 * zdzt_jp1 - zaei_slp1 * zdzt_jp1 ) * zbv * ze2vr  & 
     369                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     370                     ztfw(ji,jj+1,jk+kp) =  ztfw(ji,jj+1,jk+kp)                                                              & 
     371                                         &    - ( ( zah_slp_jp1 + zaei_slp_jp1) * zdyt_jp1 * zbv_jp1 * ze3wr_jp1             & 
     372                                         &      + ( zah_slp1 + zaei_slp1) * zdyt * zbv * ze3wr_jp1                           & 
     373                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     374                  END_2D 
    332375               END DO 
    333376               ! 
    334377            ELSE 
    335378               ! 
    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 
     379               DO kp = 0, 1 
     380                  DO_2D( 1, 0, 1, 0 ) 
     381                     ze1ur = r1_e1u(ji,jj) 
     382                     zdxt  = zdit(ji,jj,jk) * ze1ur 
     383                     zdxt_ip1  = zdit(ji+1,jj,jk) * r1_e1u(ji+1,jj) 
     384                     ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     385                     ze3wr_ip1 = 1._wp / e3w(ji+1,jj,jk+kp,Kmm) 
     386                     zdzt  = zdkt3d(ji,jj,kp) * ze3wr 
     387                     zdzt_ip1  = zdkt3d(ji+1,jj,kp) * ze3wr_ip1 
     388                     ! 
     389                     zbu = 0.25_wp * e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm) 
     390                     zbu_ip1 = 0.25_wp * e1e2u(ji+1,jj) * e3u(ji+1,jj,jk,Kmm) 
     391                     ! ln_botmix_triad is .F. mask zah for bottom half cells 
     392                     zah = pahu(ji,jj,jk) * umask(ji,jj,jk+kp)         ! pahu(ji+ip,jj,jk)   ===>>  ???? 
     393                     zah_ip1 = pahu(ji+1,jj,jk) * umask(ji+1,jj,jk+kp)  
     394                     zah_slp  = zah * triadi(ji,jj,jk,1,kp) 
     395                     zah_slp_ip1  = zah_ip1 * triadi(ji+1,jj,jk,1,kp) 
     396                     zah_slp1  = zah * triadi(ji+1,jj,jk,0,kp) 
     397                     IF( ln_ldfeiv )   THEN 
     398                        zaei_slp = aeiu(ji,jj,jk) * triadi_g(ji,jj,jk,1,kp)  
     399                        zaei_slp_ip1 = aeiu(ji+1,jj,jk) * triadi_g(ji+1,jj,jk,1,kp) 
     400                        zaei_slp1 = aeiu(ji,jj,jk) * triadi_g(ji+1,jj,jk,0,kp) 
     401                     ENDIF 
     402                     ! round brackets added to fix the order of floating point operations 
     403                     ! needed to ensure halo 1 - halo 2 compatibility 
     404                     zftu(ji   ,jj,jk  ) =  zftu(ji   ,jj,jk )                                                               & 
     405                                         &    - ( ( zah * zdxt + ( zah_slp - zaei_slp ) * zdzt ) * zbu * ze1ur               & 
     406                                         &      + ( zah * zdxt + zah_slp1 * zdzt_ip1 - zaei_slp1 * zdzt_ip1 ) * zbu * ze1ur  & 
     407                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     408                     ztfw(ji+1,jj,jk+kp) =  ztfw(ji+1,jj,jk+kp)                                                              &  
     409                                         &    - ( (zah_slp_ip1 + zaei_slp_ip1) * zdxt_ip1 * zbu_ip1 * ze3wr_ip1              & 
     410                                         &      + ( zah_slp1 + zaei_slp1) * zdxt * zbu * ze3wr_ip1                           & 
     411                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     412                  END_2D 
    355413               END DO 
    356414               ! 
    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 
     415               DO kp = 0, 1 
     416                  DO_2D( 1, 0, 1, 0 ) 
     417                     ze2vr = r1_e2v(ji,jj) 
     418                     zdyt  = zdjt(ji,jj,jk) * ze2vr 
     419                     zdyt_jp1  = zdjt(ji,jj+1,jk) * r1_e2v(ji,jj+1) 
     420                     ze3wr = 1._wp / e3w(ji,jj,jk+kp,Kmm) 
     421                     ze3wr_jp1 = 1._wp / e3w(ji,jj+1,jk+kp,Kmm) 
     422                     zdzt  = zdkt3d(ji,jj,kp) * ze3wr 
     423                     zdzt_jp1  = zdkt3d(ji,jj+1,kp) * ze3wr_jp1 
     424                     zbv = 0.25_wp * e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
     425                     zbv_jp1 = 0.25_wp * e1e2v(ji,jj+1) * e3v(ji,jj+1,jk,Kmm) 
     426                     ! ln_botmix_triad is .F. mask zah for bottom half cells 
     427                     zah = pahv(ji,jj,jk) * vmask(ji,jj,jk+kp)         ! pahv(ji,jj+jp,jk)  ???? 
     428                     zah_jp1 = pahv(ji,jj+1,jk) * vmask(ji,jj+1,jk+kp) 
     429                     zah_slp = zah * triadj(ji,jj,jk,1,kp) 
     430                     zah_slp1 = zah * triadj(ji,jj+1,jk,0,kp) 
     431                     zah_slp_jp1 = zah_jp1 * triadj(ji,jj+1,jk,1,kp) 
     432                     IF( ln_ldfeiv )   THEN 
     433                        zaei_slp = aeiv(ji,jj,jk) * triadj_g(ji,jj,jk,1,kp)  
     434                        zaei_slp_jp1 = aeiv(ji,jj+1,jk) * triadj_g(ji,jj+1,jk,1,kp) 
     435                        zaei_slp1 = aeiv(ji,jj,jk) * triadj_g(ji,jj+1,jk,0,kp) 
     436                     ENDIF 
     437                     ! round brackets added to fix the order of floating point operations 
     438                     ! needed to ensure halo 1 - halo 2 compatibility 
     439                     zftv(ji,jj  ,jk   ) =  zftv(ji,jj  ,jk   )                                                              & 
     440                                         &    - ( ( zah * zdyt + ( zah_slp - zaei_slp ) * zdzt ) * zbv * ze2vr               & 
     441                                         &      + ( zah * zdyt + zah_slp1 * zdzt_jp1 - zaei_slp1 * zdzt_jp1 ) * zbv * ze2vr  & 
     442                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     443                     ztfw(ji,jj+1,jk+kp) =  ztfw(ji,jj+1,jk+kp)                                                              & 
     444                                         &    - ( ( zah_slp_jp1 + zaei_slp_jp1) * zdyt_jp1 * zbv_jp1 * ze3wr_jp1             & 
     445                                         &      + ( zah_slp1 + zaei_slp1) * zdyt * zbv * ze3wr_jp1                           & 
     446                                         &      )                                                                            ! bracket for halo 1 - halo 2 compatibility 
     447                  END_2D 
    375448               END DO 
    376449            ENDIF 
    377450            !                             !==  horizontal divergence and add to the general trend  ==! 
    378451            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)  ) 
     452               ! round brackets added to fix the order of floating point operations 
     453               ! needed to ensure halo 1 - halo 2 compatibility 
     454               pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn)                                                & 
     455                  &                       + zsign * ( ( zftu(ji-1,jj  ,jk) - zftu(ji,jj,jk)             & 
     456                  &                                   )                                                 & ! bracket for halo 1 - halo 2 compatibility 
     457                  &                                 + ( zftv(ji,jj-1,jk) - zftv(ji,jj,jk)               & 
     458                  &                                   )                                                 & ! bracket for halo 1 - halo 2 compatibility 
     459                  &                                 ) / (  e1e2t(ji,jj) * e3t(ji,jj,jk,Kmm)  ) 
    383460            END_2D 
    384461            ! 
Note: See TracChangeset for help on using the changeset viewer.