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/DOM/domvvl.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/DOM/domvvl.F90

    r12150 r12340  
    6565   !! * Substitutions 
    6666#  include "vectopt_loop_substitute.h90" 
     67#  include "do_loop_substitute.h90" 
    6768   !!---------------------------------------------------------------------- 
    6869   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
    … …  
    190191      gdept(:,:,1,Kbb) = 0.5_wp * e3w(:,:,1,Kbb) 
    191192      gdepw(:,:,1,Kbb) = 0.0_wp 
    192       DO jk = 2, jpk                               ! vertical sum 
    193          DO jj = 1,jpj 
    194             DO ji = 1,jpi 
    195                !    zcoef = tmask - wmask    ! 0 everywhere tmask = wmask, ie everywhere expect at jk = mikt 
    196                !                             ! 1 everywhere from mbkt to mikt + 1 or 1 (if no isf) 
    197                !                             ! 0.5 where jk = mikt      
     193      DO_3D_11_11( 2, jpk ) 
     194         !    zcoef = tmask - wmask    ! 0 everywhere tmask = wmask, ie everywhere expect at jk = mikt 
     195         !                             ! 1 everywhere from mbkt to mikt + 1 or 1 (if no isf) 
     196         !                             ! 0.5 where jk = mikt      
    198197!!gm ???????   BUG ?  gdept(:,:,:,Kmm) as well as gde3w  does not include the thickness of ISF ?? 
    199                zcoef = ( tmask(ji,jj,jk) - wmask(ji,jj,jk) ) 
    200                gdepw(ji,jj,jk,Kmm) = gdepw(ji,jj,jk-1,Kmm) + e3t(ji,jj,jk-1,Kmm) 
    201                gdept(ji,jj,jk,Kmm) =      zcoef  * ( gdepw(ji,jj,jk  ,Kmm) + 0.5 * e3w(ji,jj,jk,Kmm))  & 
    202                   &                + (1-zcoef) * ( gdept(ji,jj,jk-1,Kmm) +       e3w(ji,jj,jk,Kmm))  
    203                gde3w(ji,jj,jk) = gdept(ji,jj,jk,Kmm) - ssh(ji,jj,Kmm) 
    204                gdepw(ji,jj,jk,Kbb) = gdepw(ji,jj,jk-1,Kbb) + e3t(ji,jj,jk-1,Kbb) 
    205                gdept(ji,jj,jk,Kbb) =      zcoef  * ( gdepw(ji,jj,jk  ,Kbb) + 0.5 * e3w(ji,jj,jk,Kbb))  & 
    206                   &                + (1-zcoef) * ( gdept(ji,jj,jk-1,Kbb) +       e3w(ji,jj,jk,Kbb))  
    207             END DO 
    208          END DO 
    209       END DO 
     198         zcoef = ( tmask(ji,jj,jk) - wmask(ji,jj,jk) ) 
     199         gdepw(ji,jj,jk,Kmm) = gdepw(ji,jj,jk-1,Kmm) + e3t(ji,jj,jk-1,Kmm) 
     200         gdept(ji,jj,jk,Kmm) =      zcoef  * ( gdepw(ji,jj,jk  ,Kmm) + 0.5 * e3w(ji,jj,jk,Kmm))  & 
     201            &                + (1-zcoef) * ( gdept(ji,jj,jk-1,Kmm) +       e3w(ji,jj,jk,Kmm))  
     202         gde3w(ji,jj,jk) = gdept(ji,jj,jk,Kmm) - ssh(ji,jj,Kmm) 
     203         gdepw(ji,jj,jk,Kbb) = gdepw(ji,jj,jk-1,Kbb) + e3t(ji,jj,jk-1,Kbb) 
     204         gdept(ji,jj,jk,Kbb) =      zcoef  * ( gdepw(ji,jj,jk  ,Kbb) + 0.5 * e3w(ji,jj,jk,Kbb))  & 
     205            &                + (1-zcoef) * ( gdept(ji,jj,jk-1,Kbb) +       e3w(ji,jj,jk,Kbb))  
     206      END_3D 
    210207      ! 
    211208      !                    !==  thickness of the water column  !!   (ocean portion only) 
    … …  
    242239         ENDIF 
    243240         IF ( ln_vvl_zstar_at_eqtor ) THEN   ! use z-star in vicinity of the Equator 
    244             DO jj = 1, jpj 
    245                DO ji = 1, jpi 
     241            DO_2D_11_11 
    246242!!gm  case |gphi| >= 6 degrees is useless   initialized just above by default 
    247                   IF( ABS(gphit(ji,jj)) >= 6.) THEN 
    248                      ! values outside the equatorial band and transition zone (ztilde) 
    249                      frq_rst_e3t(ji,jj) =  2.0_wp * rpi / ( MAX( rn_rst_e3t  , rsmall ) * 86400.e0_wp ) 
    250                      frq_rst_hdv(ji,jj) =  2.0_wp * rpi / ( MAX( rn_lf_cutoff, rsmall ) * 86400.e0_wp ) 
    251                   ELSEIF( ABS(gphit(ji,jj)) <= 2.5) THEN    ! Equator strip ==> z-star 
    252                      ! values inside the equatorial band (ztilde as zstar) 
    253                      frq_rst_e3t(ji,jj) =  0.0_wp 
    254                      frq_rst_hdv(ji,jj) =  1.0_wp / rdt 
    255                   ELSE                                      ! transition band (2.5 to 6 degrees N/S) 
    256                      !                                      ! (linearly transition from z-tilde to z-star) 
    257                      frq_rst_e3t(ji,jj) = 0.0_wp + (frq_rst_e3t(ji,jj)-0.0_wp)*0.5_wp   & 
    258                         &            * (  1.0_wp - COS( rad*(ABS(gphit(ji,jj))-2.5_wp)  & 
    259                         &                                          * 180._wp / 3.5_wp ) ) 
    260                      frq_rst_hdv(ji,jj) = (1.0_wp / rdt)                                & 
    261                         &            + (  frq_rst_hdv(ji,jj)-(1.e0_wp / rdt) )*0.5_wp   & 
    262                         &            * (  1._wp  - COS( rad*(ABS(gphit(ji,jj))-2.5_wp)  & 
    263                         &                                          * 180._wp / 3.5_wp ) ) 
    264                   ENDIF 
    265                END DO 
    266             END DO 
     243               IF( ABS(gphit(ji,jj)) >= 6.) THEN 
     244                  ! values outside the equatorial band and transition zone (ztilde) 
     245                  frq_rst_e3t(ji,jj) =  2.0_wp * rpi / ( MAX( rn_rst_e3t  , rsmall ) * 86400.e0_wp ) 
     246                  frq_rst_hdv(ji,jj) =  2.0_wp * rpi / ( MAX( rn_lf_cutoff, rsmall ) * 86400.e0_wp ) 
     247               ELSEIF( ABS(gphit(ji,jj)) <= 2.5) THEN    ! Equator strip ==> z-star 
     248                  ! values inside the equatorial band (ztilde as zstar) 
     249                  frq_rst_e3t(ji,jj) =  0.0_wp 
     250                  frq_rst_hdv(ji,jj) =  1.0_wp / rdt 
     251               ELSE                                      ! transition band (2.5 to 6 degrees N/S) 
     252                  !                                      ! (linearly transition from z-tilde to z-star) 
     253                  frq_rst_e3t(ji,jj) = 0.0_wp + (frq_rst_e3t(ji,jj)-0.0_wp)*0.5_wp   & 
     254                     &            * (  1.0_wp - COS( rad*(ABS(gphit(ji,jj))-2.5_wp)  & 
     255                     &                                          * 180._wp / 3.5_wp ) ) 
     256                  frq_rst_hdv(ji,jj) = (1.0_wp / rdt)                                & 
     257                     &            + (  frq_rst_hdv(ji,jj)-(1.e0_wp / rdt) )*0.5_wp   & 
     258                     &            * (  1._wp  - COS( rad*(ABS(gphit(ji,jj))-2.5_wp)  & 
     259                     &                                          * 180._wp / 3.5_wp ) ) 
     260               ENDIF 
     261            END_2D 
    267262            IF( cn_cfg == "orca" .OR. cn_cfg == "ORCA" ) THEN 
    268263               IF( nn_cfg == 3 ) THEN   ! ORCA2: Suppress ztilde in the Foxe Basin for ORCA2 
    … …  
    413408         zwu(:,:) = 0._wp 
    414409         zwv(:,:) = 0._wp 
    415          DO jk = 1, jpkm1        ! a - first derivative: diffusive fluxes 
    416             DO jj = 1, jpjm1 
    417                DO ji = 1, fs_jpim1   ! vector opt. 
    418                   un_td(ji,jj,jk) = rn_ahe3 * umask(ji,jj,jk) * e2_e1u(ji,jj)           & 
    419                      &            * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji+1,jj  ,jk) ) 
    420                   vn_td(ji,jj,jk) = rn_ahe3 * vmask(ji,jj,jk) * e1_e2v(ji,jj)           &  
    421                      &            * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji  ,jj+1,jk) ) 
    422                   zwu(ji,jj) = zwu(ji,jj) + un_td(ji,jj,jk) 
    423                   zwv(ji,jj) = zwv(ji,jj) + vn_td(ji,jj,jk) 
    424                END DO 
    425             END DO 
    426          END DO 
    427          DO jj = 1, jpj          ! b - correction for last oceanic u-v points 
    428             DO ji = 1, jpi 
    429                un_td(ji,jj,mbku(ji,jj)) = un_td(ji,jj,mbku(ji,jj)) - zwu(ji,jj) 
    430                vn_td(ji,jj,mbkv(ji,jj)) = vn_td(ji,jj,mbkv(ji,jj)) - zwv(ji,jj) 
    431             END DO 
    432          END DO 
    433          DO jk = 1, jpkm1        ! c - second derivative: divergence of diffusive fluxes 
    434             DO jj = 2, jpjm1 
    435                DO ji = fs_2, fs_jpim1   ! vector opt. 
    436                   tilde_e3t_a(ji,jj,jk) = tilde_e3t_a(ji,jj,jk) + (   un_td(ji-1,jj  ,jk) - un_td(ji,jj,jk)    & 
    437                      &                                          +     vn_td(ji  ,jj-1,jk) - vn_td(ji,jj,jk)    & 
    438                      &                                            ) * r1_e1e2t(ji,jj) 
    439                END DO 
    440             END DO 
    441          END DO 
     410         DO_3D_10_10( 1, jpkm1 ) 
     411            un_td(ji,jj,jk) = rn_ahe3 * umask(ji,jj,jk) * e2_e1u(ji,jj)           & 
     412               &            * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji+1,jj  ,jk) ) 
     413            vn_td(ji,jj,jk) = rn_ahe3 * vmask(ji,jj,jk) * e1_e2v(ji,jj)           &  
     414               &            * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji  ,jj+1,jk) ) 
     415            zwu(ji,jj) = zwu(ji,jj) + un_td(ji,jj,jk) 
     416            zwv(ji,jj) = zwv(ji,jj) + vn_td(ji,jj,jk) 
     417         END_3D 
     418         DO_2D_11_11 
     419            un_td(ji,jj,mbku(ji,jj)) = un_td(ji,jj,mbku(ji,jj)) - zwu(ji,jj) 
     420            vn_td(ji,jj,mbkv(ji,jj)) = vn_td(ji,jj,mbkv(ji,jj)) - zwv(ji,jj) 
     421         END_2D 
     422         DO_3D_00_00( 1, jpkm1 ) 
     423            tilde_e3t_a(ji,jj,jk) = tilde_e3t_a(ji,jj,jk) + (   un_td(ji-1,jj  ,jk) - un_td(ji,jj,jk)    & 
     424               &                                          +     vn_td(ji  ,jj-1,jk) - vn_td(ji,jj,jk)    & 
     425               &                                            ) * r1_e1e2t(ji,jj) 
     426         END_3D 
    442427         !                       ! d - thickness diffusion transport: boundary conditions 
    443428         !                             (stored for tracer advction and continuity equation) 
    … …  
    670655      gdepw(:,:,1,Kmm) = 0.0_wp 
    671656      gde3w(:,:,1) = gdept(:,:,1,Kmm) - ssh(:,:,Kmm) 
    672       DO jk = 2, jpk 
    673          DO jj = 1,jpj 
    674             DO ji = 1,jpi 
    675               !    zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))   ! 0 everywhere tmask = wmask, ie everywhere expect at jk = mikt 
    676                                                                  ! 1 for jk = mikt 
    677                zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) 
    678                gdepw(ji,jj,jk,Kmm) = gdepw(ji,jj,jk-1,Kmm) + e3t(ji,jj,jk-1,Kmm) 
    679                gdept(ji,jj,jk,Kmm) =    zcoef  * ( gdepw(ji,jj,jk  ,Kmm) + 0.5 * e3w(ji,jj,jk,Kmm) )  & 
    680                    &             + (1-zcoef) * ( gdept(ji,jj,jk-1,Kmm) +       e3w(ji,jj,jk,Kmm) )  
    681                gde3w(ji,jj,jk) = gdept(ji,jj,jk,Kmm) - ssh(ji,jj,Kmm) 
    682             END DO 
    683          END DO 
    684       END DO 
     657      DO_3D_11_11( 2, jpk ) 
     658        !    zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk))   ! 0 everywhere tmask = wmask, ie everywhere expect at jk = mikt 
     659                                                           ! 1 for jk = mikt 
     660         zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) 
     661         gdepw(ji,jj,jk,Kmm) = gdepw(ji,jj,jk-1,Kmm) + e3t(ji,jj,jk-1,Kmm) 
     662         gdept(ji,jj,jk,Kmm) =    zcoef  * ( gdepw(ji,jj,jk  ,Kmm) + 0.5 * e3w(ji,jj,jk,Kmm) )  & 
     663             &             + (1-zcoef) * ( gdept(ji,jj,jk-1,Kmm) +       e3w(ji,jj,jk,Kmm) )  
     664         gde3w(ji,jj,jk) = gdept(ji,jj,jk,Kmm) - ssh(ji,jj,Kmm) 
     665      END_3D 
    685666 
    686667      ! Local depth and Inverse of the local depth of the water 
    … …  
    729710         ! 
    730711      CASE( 'U' )                   !* from T- to U-point : hor. surface weighted mean 
    731          DO jk = 1, jpk 
    732             DO jj = 1, jpjm1 
    733                DO ji = 1, fs_jpim1   ! vector opt. 
    734                   pe3_out(ji,jj,jk) = 0.5_wp * (  umask(ji,jj,jk) * (1.0_wp - zlnwd) + zlnwd ) * r1_e1e2u(ji,jj)   & 
    735                      &                       * (   e1e2t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) )     & 
    736                      &                           + e1e2t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) ) 
    737                END DO 
    738             END DO 
    739          END DO 
     712         DO_3D_10_10( 1, jpk ) 
     713            pe3_out(ji,jj,jk) = 0.5_wp * (  umask(ji,jj,jk) * (1.0_wp - zlnwd) + zlnwd ) * r1_e1e2u(ji,jj)   & 
     714               &                       * (   e1e2t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) )     & 
     715               &                           + e1e2t(ji+1,jj) * ( pe3_in(ji+1,jj,jk) - e3t_0(ji+1,jj,jk) ) ) 
     716         END_3D 
    740717         CALL lbc_lnk( 'domvvl', pe3_out(:,:,:), 'U', 1._wp ) 
    741718         pe3_out(:,:,:) = pe3_out(:,:,:) + e3u_0(:,:,:) 
    742719         ! 
    743720      CASE( 'V' )                   !* from T- to V-point : hor. surface weighted mean 
    744          DO jk = 1, jpk 
    745             DO jj = 1, jpjm1 
    746                DO ji = 1, fs_jpim1   ! vector opt. 
    747                   pe3_out(ji,jj,jk) = 0.5_wp * ( vmask(ji,jj,jk)  * (1.0_wp - zlnwd) + zlnwd ) * r1_e1e2v(ji,jj)   & 
    748                      &                       * (   e1e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) )     & 
    749                      &                           + e1e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) ) 
    750                END DO 
    751             END DO 
    752          END DO 
     721         DO_3D_10_10( 1, jpk ) 
     722            pe3_out(ji,jj,jk) = 0.5_wp * ( vmask(ji,jj,jk)  * (1.0_wp - zlnwd) + zlnwd ) * r1_e1e2v(ji,jj)   & 
     723               &                       * (   e1e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) )     & 
     724               &                           + e1e2t(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3t_0(ji,jj+1,jk) ) ) 
     725         END_3D 
    753726         CALL lbc_lnk( 'domvvl', pe3_out(:,:,:), 'V', 1._wp ) 
    754727         pe3_out(:,:,:) = pe3_out(:,:,:) + e3v_0(:,:,:) 
    755728         ! 
    756729      CASE( 'F' )                   !* from U-point to F-point : hor. surface weighted mean 
    757          DO jk = 1, jpk 
    758             DO jj = 1, jpjm1 
    759                DO ji = 1, fs_jpim1   ! vector opt. 
    760                   pe3_out(ji,jj,jk) = 0.5_wp * (  umask(ji,jj,jk) * umask(ji,jj+1,jk) * (1.0_wp - zlnwd) + zlnwd ) & 
    761                      &                       *    r1_e1e2f(ji,jj)                                                  & 
    762                      &                       * (   e1e2u(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3u_0(ji,jj  ,jk) )     & 
    763                      &                           + e1e2u(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3u_0(ji,jj+1,jk) ) ) 
    764                END DO 
    765             END DO 
    766          END DO 
     730         DO_3D_10_10( 1, jpk ) 
     731            pe3_out(ji,jj,jk) = 0.5_wp * (  umask(ji,jj,jk) * umask(ji,jj+1,jk) * (1.0_wp - zlnwd) + zlnwd ) & 
     732               &                       *    r1_e1e2f(ji,jj)                                                  & 
     733               &                       * (   e1e2u(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3u_0(ji,jj  ,jk) )     & 
     734               &                           + e1e2u(ji,jj+1) * ( pe3_in(ji,jj+1,jk) - e3u_0(ji,jj+1,jk) ) ) 
     735         END_3D 
    767736         CALL lbc_lnk( 'domvvl', pe3_out(:,:,:), 'F', 1._wp ) 
    768737         pe3_out(:,:,:) = pe3_out(:,:,:) + e3f_0(:,:,:) 
    … …  
    926895                  ssh(:,:,Kbb) = -ssh_ref 
    927896 
    928                   DO jj = 1, jpj 
    929                      DO ji = 1, jpi 
    930                         IF( ht_0(ji,jj)-ssh_ref <  rn_wdmin1 ) THEN ! if total depth is less than min depth 
    931                            ssh(ji,jj,Kbb) = rn_wdmin1 - (ht_0(ji,jj) ) 
    932                            ssh(ji,jj,Kmm) = rn_wdmin1 - (ht_0(ji,jj) ) 
    933                         ENDIF 
    934                      ENDDO 
    935                   ENDDO 
     897                  DO_2D_11_11 
     898                     IF( ht_0(ji,jj)-ssh_ref <  rn_wdmin1 ) THEN ! if total depth is less than min depth 
     899                        ssh(ji,jj,Kbb) = rn_wdmin1 - (ht_0(ji,jj) ) 
     900                        ssh(ji,jj,Kmm) = rn_wdmin1 - (ht_0(ji,jj) ) 
     901                     ENDIF 
     902                  END_2D 
    936903               ENDIF !If test case else 
    937904 
Note: See TracChangeset for help on using the changeset viewer.