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 12340 for NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/TRA/traldf_iso.F90 – NEMO

Ignore:
Timestamp:
2020-01-27T15:31:53+01:00 (4 years ago)
Author:
acc
Message:

Branch 2019/dev_r11943_MERGE_2019. This commit introduces basic do loop macro
substitution to the 2019 option 1, merge branch. These changes have been SETTE
tested. The only addition is the do_loop_substitute.h90 file in the OCE directory but
the macros defined therein are used throughout the code to replace identifiable, 2D-
and 3D- nested loop opening and closing statements with single-line alternatives. Code
indents are also adjusted accordingly.

The following explanation is taken from comments in the new header file:

This header file contains preprocessor definitions and macros used in the do-loop
substitutions introduced between version 4.0 and 4.2. The primary aim of these macros
is to assist in future applications of tiling to improve performance. This is expected
to be achieved by alternative versions of these macros in selected locations. The
initial introduction of these macros simply replaces all identifiable nested 2D- and
3D-loops with single line statements (and adjusts indenting accordingly). Do loops
are identifiable if they comform to either:

DO jk = ....

DO jj = .... DO jj = ...

DO ji = .... DO ji = ...
. OR .
. .

END DO END DO

END DO END DO

END DO

and white-space variants thereof.

Additionally, only loops with recognised jj and ji loops limits are treated; these are:
Lower limits of 1, 2 or fs_2
Upper limits of jpi, jpim1 or fs_jpim1 (for ji) or jpj, jpjm1 or fs_jpjm1 (for jj)

The macro naming convention takes the form: DO_2D_BT_LR where:

B is the Bottom offset from the PE's inner domain;
T is the Top offset from the PE's inner domain;
L is the Left offset from the PE's inner domain;
R is the Right offset from the PE's inner domain

So, given an inner domain of 2,jpim1 and 2,jpjm1, a typical example would replace:

DO jj = 2, jpj

DO ji = 1, jpim1
.
.

END DO

END DO

with:

DO_2D_01_10
.
.
END_2D

similar conventions apply to the 3D loops macros. jk loop limits are retained
through macro arguments and are not restricted. This includes the possibility of
strides for which an extra set of DO_3DS macros are defined.

In the example definition below the inner PE domain is defined by start indices of
(kIs, kJs) and end indices of (kIe, KJe)

#define DO_2D_00_00 DO jj = kJs, kJe ; DO ji = kIs, kIe
#define END_2D END DO ; END DO

TO DO:


Only conventional nested loops have been identified and replaced by this step. There are constructs such as:

DO jk = 2, jpkm1

z2d(:,:) = z2d(:,:) + e3w(:,:,jk,Kmm) * z3d(:,:,jk) * wmask(:,:,jk)

END DO

which may need to be considered.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/TRA/traldf_iso.F90

    r12193 r12340  
    4141   !! * Substitutions 
    4242#  include "vectopt_loop_substitute.h90" 
     43#  include "do_loop_substitute.h90" 
    4344   !!---------------------------------------------------------------------- 
    4445   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    145146      IF( kpass == 1 ) THEN                  !==  first pass only  ==! 
    146147         ! 
    147          DO jk = 2, jpkm1 
    148             DO jj = 2, jpjm1 
    149                DO ji = fs_2, fs_jpim1   ! vector opt. 
    150                   ! 
    151                   zmsku = wmask(ji,jj,jk) / MAX(   umask(ji  ,jj,jk-1) + umask(ji-1,jj,jk)          & 
    152                      &                           + umask(ji-1,jj,jk-1) + umask(ji  ,jj,jk) , 1._wp  ) 
    153                   zmskv = wmask(ji,jj,jk) / MAX(   vmask(ji,jj  ,jk-1) + vmask(ji,jj-1,jk)          & 
    154                      &                           + vmask(ji,jj-1,jk-1) + vmask(ji,jj  ,jk) , 1._wp  ) 
    155                      ! 
    156                   zahu_w = (   pahu(ji  ,jj,jk-1) + pahu(ji-1,jj,jk)    & 
    157                      &       + pahu(ji-1,jj,jk-1) + pahu(ji  ,jj,jk)  ) * zmsku 
    158                   zahv_w = (   pahv(ji,jj  ,jk-1) + pahv(ji,jj-1,jk)    & 
    159                      &       + pahv(ji,jj-1,jk-1) + pahv(ji,jj  ,jk)  ) * zmskv 
    160                      ! 
    161                   ah_wslp2(ji,jj,jk) = zahu_w * wslpi(ji,jj,jk) * wslpi(ji,jj,jk)   & 
    162                      &               + zahv_w * wslpj(ji,jj,jk) * wslpj(ji,jj,jk) 
    163                END DO 
    164             END DO 
    165          END DO 
     148         DO_3D_00_00( 2, jpkm1 ) 
     149            ! 
     150            zmsku = wmask(ji,jj,jk) / MAX(   umask(ji  ,jj,jk-1) + umask(ji-1,jj,jk)          & 
     151               &                           + umask(ji-1,jj,jk-1) + umask(ji  ,jj,jk) , 1._wp  ) 
     152            zmskv = wmask(ji,jj,jk) / MAX(   vmask(ji,jj  ,jk-1) + vmask(ji,jj-1,jk)          & 
     153               &                           + vmask(ji,jj-1,jk-1) + vmask(ji,jj  ,jk) , 1._wp  ) 
     154               ! 
     155            zahu_w = (   pahu(ji  ,jj,jk-1) + pahu(ji-1,jj,jk)    & 
     156               &       + pahu(ji-1,jj,jk-1) + pahu(ji  ,jj,jk)  ) * zmsku 
     157            zahv_w = (   pahv(ji,jj  ,jk-1) + pahv(ji,jj-1,jk)    & 
     158               &       + pahv(ji,jj-1,jk-1) + pahv(ji,jj  ,jk)  ) * zmskv 
     159               ! 
     160            ah_wslp2(ji,jj,jk) = zahu_w * wslpi(ji,jj,jk) * wslpi(ji,jj,jk)   & 
     161               &               + zahv_w * wslpj(ji,jj,jk) * wslpj(ji,jj,jk) 
     162         END_3D 
    166163         ! 
    167164         IF( ln_traldf_msc ) THEN                ! stabilizing vertical diffusivity coefficient 
    168             DO jk = 2, jpkm1 
    169                DO jj = 2, jpjm1 
    170                   DO ji = fs_2, fs_jpim1 
    171                      akz(ji,jj,jk) = 0.25_wp * (                                                                     & 
    172                         &              ( pahu(ji  ,jj,jk) + pahu(ji  ,jj,jk-1) ) / ( e1u(ji  ,jj) * e1u(ji  ,jj) )   & 
    173                         &            + ( pahu(ji-1,jj,jk) + pahu(ji-1,jj,jk-1) ) / ( e1u(ji-1,jj) * e1u(ji-1,jj) )   & 
    174                         &            + ( pahv(ji,jj  ,jk) + pahv(ji,jj  ,jk-1) ) / ( e2v(ji,jj  ) * e2v(ji,jj  ) )   & 
    175                         &            + ( pahv(ji,jj-1,jk) + pahv(ji,jj-1,jk-1) ) / ( e2v(ji,jj-1) * e2v(ji,jj-1) )   ) 
    176                   END DO 
    177                END DO 
    178             END DO 
     165            DO_3D_00_00( 2, jpkm1 ) 
     166               akz(ji,jj,jk) = 0.25_wp * (                                                                     & 
     167                  &              ( pahu(ji  ,jj,jk) + pahu(ji  ,jj,jk-1) ) / ( e1u(ji  ,jj) * e1u(ji  ,jj) )   & 
     168                  &            + ( pahu(ji-1,jj,jk) + pahu(ji-1,jj,jk-1) ) / ( e1u(ji-1,jj) * e1u(ji-1,jj) )   & 
     169                  &            + ( pahv(ji,jj  ,jk) + pahv(ji,jj  ,jk-1) ) / ( e2v(ji,jj  ) * e2v(ji,jj  ) )   & 
     170                  &            + ( pahv(ji,jj-1,jk) + pahv(ji,jj-1,jk-1) ) / ( e2v(ji,jj-1) * e2v(ji,jj-1) )   ) 
     171            END_3D 
    179172            ! 
    180173            IF( ln_traldf_blp ) THEN                ! bilaplacian operator 
    181                DO jk = 2, jpkm1 
    182                   DO jj = 1, jpjm1 
    183                      DO ji = 1, fs_jpim1 
    184                         akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk)   & 
    185                            &          * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) )  ) 
    186                      END DO 
    187                   END DO 
    188                END DO 
     174               DO_3D_10_10( 2, jpkm1 ) 
     175                  akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk)   & 
     176                     &          * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) )  ) 
     177               END_3D 
    189178            ELSEIF( ln_traldf_lap ) THEN              ! laplacian operator 
    190                DO jk = 2, jpkm1 
    191                   DO jj = 1, jpjm1 
    192                      DO ji = 1, fs_jpim1 
    193                         ze3w_2 = e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) 
    194                         zcoef0 = z2dt * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2  ) 
    195                         akz(ji,jj,jk) = MAX( zcoef0 - 0.5_wp , 0._wp ) * ze3w_2 * z1_2dt 
    196                      END DO 
    197                   END DO 
    198                END DO 
     179               DO_3D_10_10( 2, jpkm1 ) 
     180                  ze3w_2 = e3w(ji,jj,jk,Kmm) * e3w(ji,jj,jk,Kmm) 
     181                  zcoef0 = z2dt * (  akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2  ) 
     182                  akz(ji,jj,jk) = MAX( zcoef0 - 0.5_wp , 0._wp ) * ze3w_2 * z1_2dt 
     183               END_3D 
    199184           ENDIF 
    200185           ! 
     
    217202 
    218203         ! Horizontal tracer gradient  
    219          DO jk = 1, jpkm1 
    220             DO jj = 1, jpjm1 
    221                DO ji = 1, fs_jpim1   ! vector opt. 
    222                   zdit(ji,jj,jk) = ( pt(ji+1,jj  ,jk,jn) - pt(ji,jj,jk,jn) ) * umask(ji,jj,jk) 
    223                   zdjt(ji,jj,jk) = ( pt(ji  ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 
    224                END DO 
    225             END DO 
    226          END DO 
     204         DO_3D_10_10( 1, jpkm1 ) 
     205            zdit(ji,jj,jk) = ( pt(ji+1,jj  ,jk,jn) - pt(ji,jj,jk,jn) ) * umask(ji,jj,jk) 
     206            zdjt(ji,jj,jk) = ( pt(ji  ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 
     207         END_3D 
    227208         IF( ln_zps ) THEN      ! botton and surface ocean correction of the horizontal gradient 
    228             DO jj = 1, jpjm1              ! bottom correction (partial bottom cell) 
    229                DO ji = 1, fs_jpim1   ! vector opt. 
    230                   zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn)           
    231                   zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn) 
    232                END DO 
    233             END DO 
     209            DO_2D_10_10 
     210               zdit(ji,jj,mbku(ji,jj)) = pgu(ji,jj,jn)           
     211               zdjt(ji,jj,mbkv(ji,jj)) = pgv(ji,jj,jn) 
     212            END_2D 
    234213            IF( ln_isfcav ) THEN      ! first wet level beneath a cavity 
    235                DO jj = 1, jpjm1 
    236                   DO ji = 1, fs_jpim1   ! vector opt. 
    237                      IF( miku(ji,jj) > 1 )   zdit(ji,jj,miku(ji,jj)) = pgui(ji,jj,jn)           
    238                      IF( mikv(ji,jj) > 1 )   zdjt(ji,jj,mikv(ji,jj)) = pgvi(ji,jj,jn)      
    239                   END DO 
    240                END DO 
     214               DO_2D_10_10 
     215                  IF( miku(ji,jj) > 1 )   zdit(ji,jj,miku(ji,jj)) = pgui(ji,jj,jn)           
     216                  IF( mikv(ji,jj) > 1 )   zdjt(ji,jj,mikv(ji,jj)) = pgvi(ji,jj,jn)      
     217               END_2D 
    241218            ENDIF 
    242219         ENDIF 
     
    254231            ELSE                 ;   zdkt(:,:) = ( pt(:,:,jk-1,jn) - pt(:,:,jk,jn) ) * wmask(:,:,jk) 
    255232            ENDIF 
    256             DO jj = 1 , jpjm1            !==  Horizontal fluxes 
    257                DO ji = 1, fs_jpim1   ! vector opt. 
    258                   zabe1 = pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm) 
    259                   zabe2 = pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
    260                   ! 
    261                   zmsku = 1. / MAX(  wmask(ji+1,jj,jk  ) + wmask(ji,jj,jk+1)   & 
    262                      &             + wmask(ji+1,jj,jk+1) + wmask(ji,jj,jk  ), 1. ) 
    263                   ! 
    264                   zmskv = 1. / MAX(  wmask(ji,jj+1,jk  ) + wmask(ji,jj,jk+1)   & 
    265                      &             + wmask(ji,jj+1,jk+1) + wmask(ji,jj,jk  ), 1. ) 
    266                   ! 
    267                   zcof1 = - pahu(ji,jj,jk) * e2u(ji,jj) * uslp(ji,jj,jk) * zmsku 
    268                   zcof2 = - pahv(ji,jj,jk) * e1v(ji,jj) * vslp(ji,jj,jk) * zmskv 
    269                   ! 
    270                   zftu(ji,jj,jk ) = (  zabe1 * zdit(ji,jj,jk)   & 
    271                      &               + zcof1 * (  zdkt (ji+1,jj) + zdk1t(ji,jj)      & 
    272                      &                          + zdk1t(ji+1,jj) + zdkt (ji,jj)  )  ) * umask(ji,jj,jk) 
    273                   zftv(ji,jj,jk) = (  zabe2 * zdjt(ji,jj,jk)   & 
    274                      &               + zcof2 * (  zdkt (ji,jj+1) + zdk1t(ji,jj)      & 
    275                      &                          + zdk1t(ji,jj+1) + zdkt (ji,jj)  )  ) * vmask(ji,jj,jk)                   
    276                END DO 
    277             END DO 
    278             ! 
    279             DO jj = 2 , jpjm1          !== horizontal divergence and add to pt_rhs 
    280                DO ji = fs_2, fs_jpim1   ! vector opt. 
    281                   pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  zftu(ji,jj,jk) - zftu(ji-1,jj,jk)      & 
    282                      &                                                 + zftv(ji,jj,jk) - zftv(ji,jj-1,jk)  )   & 
    283                      &                                              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    284                END DO 
    285             END DO 
     233            DO_2D_10_10 
     234               zabe1 = pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm) 
     235               zabe2 = pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm) 
     236               ! 
     237               zmsku = 1. / MAX(  wmask(ji+1,jj,jk  ) + wmask(ji,jj,jk+1)   & 
     238                  &             + wmask(ji+1,jj,jk+1) + wmask(ji,jj,jk  ), 1. ) 
     239               ! 
     240               zmskv = 1. / MAX(  wmask(ji,jj+1,jk  ) + wmask(ji,jj,jk+1)   & 
     241                  &             + wmask(ji,jj+1,jk+1) + wmask(ji,jj,jk  ), 1. ) 
     242               ! 
     243               zcof1 = - pahu(ji,jj,jk) * e2u(ji,jj) * uslp(ji,jj,jk) * zmsku 
     244               zcof2 = - pahv(ji,jj,jk) * e1v(ji,jj) * vslp(ji,jj,jk) * zmskv 
     245               ! 
     246               zftu(ji,jj,jk ) = (  zabe1 * zdit(ji,jj,jk)   & 
     247                  &               + zcof1 * (  zdkt (ji+1,jj) + zdk1t(ji,jj)      & 
     248                  &                          + zdk1t(ji+1,jj) + zdkt (ji,jj)  )  ) * umask(ji,jj,jk) 
     249               zftv(ji,jj,jk) = (  zabe2 * zdjt(ji,jj,jk)   & 
     250                  &               + zcof2 * (  zdkt (ji,jj+1) + zdk1t(ji,jj)      & 
     251                  &                          + zdk1t(ji,jj+1) + zdkt (ji,jj)  )  ) * vmask(ji,jj,jk)                   
     252            END_2D 
     253            ! 
     254            DO_2D_00_00 
     255               pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  zftu(ji,jj,jk) - zftu(ji-1,jj,jk)      & 
     256                  &                                                 + zftv(ji,jj,jk) - zftv(ji,jj-1,jk)  )   & 
     257                  &                                              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     258            END_2D 
    286259         END DO                                        !   End of slab   
    287260 
     
    297270         ztfw(:,:, 1 ) = 0._wp      ;      ztfw(:,:,jpk) = 0._wp 
    298271          
    299          DO jk = 2, jpkm1           ! interior (2=<jk=<jpk-1) 
    300             DO jj = 2, jpjm1 
    301                DO ji = fs_2, fs_jpim1   ! vector opt. 
    302                   ! 
    303                   zmsku = wmask(ji,jj,jk) / MAX(   umask(ji  ,jj,jk-1) + umask(ji-1,jj,jk)          & 
    304                      &                           + umask(ji-1,jj,jk-1) + umask(ji  ,jj,jk) , 1._wp  ) 
    305                   zmskv = wmask(ji,jj,jk) / MAX(   vmask(ji,jj  ,jk-1) + vmask(ji,jj-1,jk)          & 
    306                      &                           + vmask(ji,jj-1,jk-1) + vmask(ji,jj  ,jk) , 1._wp  ) 
    307                      ! 
    308                   zahu_w = (   pahu(ji  ,jj,jk-1) + pahu(ji-1,jj,jk)    & 
    309                      &       + pahu(ji-1,jj,jk-1) + pahu(ji  ,jj,jk)  ) * zmsku 
    310                   zahv_w = (   pahv(ji,jj  ,jk-1) + pahv(ji,jj-1,jk)    & 
    311                      &       + pahv(ji,jj-1,jk-1) + pahv(ji,jj  ,jk)  ) * zmskv 
    312                      ! 
    313                   zcoef3 = - zahu_w * e2t(ji,jj) * zmsku * wslpi (ji,jj,jk)   !wslpi & j are already w-masked 
    314                   zcoef4 = - zahv_w * e1t(ji,jj) * zmskv * wslpj (ji,jj,jk) 
    315                   ! 
    316                   ztfw(ji,jj,jk) = zcoef3 * (   zdit(ji  ,jj  ,jk-1) + zdit(ji-1,jj  ,jk)      & 
    317                      &                        + zdit(ji-1,jj  ,jk-1) + zdit(ji  ,jj  ,jk)  )   & 
    318                      &           + zcoef4 * (   zdjt(ji  ,jj  ,jk-1) + zdjt(ji  ,jj-1,jk)      & 
    319                      &                        + zdjt(ji  ,jj-1,jk-1) + zdjt(ji  ,jj  ,jk)  ) 
    320                END DO 
    321             END DO 
    322          END DO 
     272         DO_3D_00_00( 2, jpkm1 ) 
     273            ! 
     274            zmsku = wmask(ji,jj,jk) / MAX(   umask(ji  ,jj,jk-1) + umask(ji-1,jj,jk)          & 
     275               &                           + umask(ji-1,jj,jk-1) + umask(ji  ,jj,jk) , 1._wp  ) 
     276            zmskv = wmask(ji,jj,jk) / MAX(   vmask(ji,jj  ,jk-1) + vmask(ji,jj-1,jk)          & 
     277               &                           + vmask(ji,jj-1,jk-1) + vmask(ji,jj  ,jk) , 1._wp  ) 
     278               ! 
     279            zahu_w = (   pahu(ji  ,jj,jk-1) + pahu(ji-1,jj,jk)    & 
     280               &       + pahu(ji-1,jj,jk-1) + pahu(ji  ,jj,jk)  ) * zmsku 
     281            zahv_w = (   pahv(ji,jj  ,jk-1) + pahv(ji,jj-1,jk)    & 
     282               &       + pahv(ji,jj-1,jk-1) + pahv(ji,jj  ,jk)  ) * zmskv 
     283               ! 
     284            zcoef3 = - zahu_w * e2t(ji,jj) * zmsku * wslpi (ji,jj,jk)   !wslpi & j are already w-masked 
     285            zcoef4 = - zahv_w * e1t(ji,jj) * zmskv * wslpj (ji,jj,jk) 
     286            ! 
     287            ztfw(ji,jj,jk) = zcoef3 * (   zdit(ji  ,jj  ,jk-1) + zdit(ji-1,jj  ,jk)      & 
     288               &                        + zdit(ji-1,jj  ,jk-1) + zdit(ji  ,jj  ,jk)  )   & 
     289               &           + zcoef4 * (   zdjt(ji  ,jj  ,jk-1) + zdjt(ji  ,jj-1,jk)      & 
     290               &                        + zdjt(ji  ,jj-1,jk-1) + zdjt(ji  ,jj  ,jk)  ) 
     291         END_3D 
    323292         !                                !==  add the vertical 33 flux  ==! 
    324293         IF( ln_traldf_lap ) THEN               ! laplacian case: eddy coef = ah_wslp2 - akz 
    325             DO jk = 2, jpkm1        
    326                DO jj = 2, jpjm1 
    327                   DO ji = fs_2, fs_jpim1 
    328                      ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * wmask(ji,jj,jk)   & 
    329                         &                            * ( ah_wslp2(ji,jj,jk) - akz(ji,jj,jk) )               & 
    330                         &                            * (  pt(ji,jj,jk-1,jn) -  pt(ji,jj,jk,jn) ) 
    331                   END DO 
    332                END DO 
    333             END DO 
     294            DO_3D_00_00( 2, jpkm1 ) 
     295               ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * wmask(ji,jj,jk)   & 
     296                  &                            * ( ah_wslp2(ji,jj,jk) - akz(ji,jj,jk) )               & 
     297                  &                            * (  pt(ji,jj,jk-1,jn) -  pt(ji,jj,jk,jn) ) 
     298            END_3D 
    334299            ! 
    335300         ELSE                                   ! bilaplacian  
    336301            SELECT CASE( kpass ) 
    337302            CASE(  1  )                            ! 1st pass : eddy coef = ah_wslp2 
    338                DO jk = 2, jpkm1  
    339                   DO jj = 2, jpjm1 
    340                      DO ji = fs_2, fs_jpim1 
    341                         ztfw(ji,jj,jk) = ztfw(ji,jj,jk)                       & 
    342                            &           + ah_wslp2(ji,jj,jk)  * e1e2t(ji,jj)   & 
    343                            &           * ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) / e3w(ji,jj,jk,Kmm) * wmask(ji,jj,jk) 
    344                      END DO 
    345                   END DO 
    346                END DO  
     303               DO_3D_00_00( 2, jpkm1 ) 
     304                  ztfw(ji,jj,jk) = ztfw(ji,jj,jk)                       & 
     305                     &           + ah_wslp2(ji,jj,jk)  * e1e2t(ji,jj)   & 
     306                     &           * ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) / e3w(ji,jj,jk,Kmm) * wmask(ji,jj,jk) 
     307               END_3D 
    347308            CASE(  2  )                         ! 2nd pass : eddy flux = ah_wslp2 and akz applied on pt  and pt2 gradients, resp. 
    348                DO jk = 2, jpkm1  
    349                   DO jj = 2, jpjm1 
    350                      DO ji = fs_2, fs_jpim1 
    351                         ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * wmask(ji,jj,jk)                  & 
    352                            &                            * (  ah_wslp2(ji,jj,jk) * ( pt (ji,jj,jk-1,jn) - pt (ji,jj,jk,jn) )   & 
    353                            &                            +         akz(ji,jj,jk) * ( pt2(ji,jj,jk-1,jn) - pt2(ji,jj,jk,jn) )   ) 
    354                      END DO 
    355                   END DO 
    356                END DO 
     309               DO_3D_00_00( 2, jpkm1 ) 
     310                  ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w(ji,jj,jk,Kmm) * wmask(ji,jj,jk)                  & 
     311                     &                            * (  ah_wslp2(ji,jj,jk) * ( pt (ji,jj,jk-1,jn) - pt (ji,jj,jk,jn) )   & 
     312                     &                            +         akz(ji,jj,jk) * ( pt2(ji,jj,jk-1,jn) - pt2(ji,jj,jk,jn) )   ) 
     313               END_3D 
    357314            END SELECT 
    358315         ENDIF 
    359316         !          
    360          DO jk = 1, jpkm1                 !==  Divergence of vertical fluxes added to pt_rhs  ==! 
    361             DO jj = 2, jpjm1 
    362                DO ji = fs_2, fs_jpim1   ! vector opt. 
    363                   pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  ztfw (ji,jj,jk) - ztfw(ji,jj,jk+1)  )   & 
    364                      &                                              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
    365                END DO 
    366             END DO 
    367          END DO 
     317         DO_3D_00_00( 1, jpkm1 ) 
     318            pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * (  ztfw (ji,jj,jk) - ztfw(ji,jj,jk+1)  )   & 
     319               &                                              * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kmm) 
     320         END_3D 
    368321         ! 
    369322         IF( ( kpass == 1 .AND. ln_traldf_lap ) .OR.  &     !==  first pass only (  laplacian)  ==! 
Note: See TracChangeset for help on using the changeset viewer.