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 13295 for NEMO/trunk/src/OCE/ZDF/zdfgls.F90 – NEMO

Ignore:
Timestamp:
2020-07-10T20:24:21+02:00 (4 years ago)
Author:
acc
Message:

Replace do-loop macros in the trunk with alternative forms with greater flexibility for extra halo applications. This alters a lot of routines but does not change any behaviour or results. do_loop_substitute.h90 is greatly simplified by this change. SETTE results are identical to those with the previous revision

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk/src/OCE/ZDF/zdfgls.F90

    r13286 r13295  
    168168 
    169169      ! Compute surface, top and bottom friction at T-points 
    170       DO_2D_00_00 
     170      DO_2D( 0, 0, 0, 0 ) 
    171171         ! 
    172172         ! surface friction 
     
    181181      END_2D 
    182182      IF( ln_isfcav ) THEN       !top friction 
    183          DO_2D_00_00 
     183         DO_2D( 0, 0, 0, 0 ) 
    184184            zmsku = ( 2. - umask(ji-1,jj,mikt(ji,jj)) * umask(ji,jj,mikt(ji,jj)) ) 
    185185            zmskv = ( 2. - vmask(ji,jj-1,mikt(ji,jj)) * vmask(ji,jj,mikt(ji,jj)) )     ! (CAUTION: CdU<0) 
     
    204204      END SELECT 
    205205      ! 
    206       DO_3D_10_10( 2, jpkm1 ) 
     206      DO_3D( 1, 0, 1, 0, 2, jpkm1 ) 
    207207         eps(ji,jj,jk)  = rc03 * en(ji,jj,jk) * SQRT( en(ji,jj,jk) ) / hmxl_n(ji,jj,jk) 
    208208      END_3D 
     
    213213 
    214214      IF( nn_clos == 0 ) THEN    ! Mellor-Yamada 
    215          DO_3D_00_00( 2, jpkm1 ) 
     215         DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    216216            zup   = hmxl_n(ji,jj,jk) * gdepw(ji,jj,mbkt(ji,jj)+1,Kmm) 
    217217            zdown = vkarmn * gdepw(ji,jj,jk,Kmm) * ( -gdepw(ji,jj,jk,Kmm) + gdepw(ji,jj,mbkt(ji,jj)+1,Kmm) ) 
     
    234234      ! Warning : after this step, en : right hand side of the matrix 
    235235 
    236       DO_3D_00_00( 2, jpkm1 ) 
     236      DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    237237         ! 
    238238         buoy = - p_avt(ji,jj,jk) * rn2(ji,jj,jk)     ! stratif. destruction 
     
    330330         !                      ! en(ibot) = u*^2 / Co2 and hmxl_n(ibot) = rn_lmin 
    331331         !                      ! Balance between the production and the dissipation terms 
    332          DO_2D_00_00 
     332         DO_2D( 0, 0, 0, 0 ) 
    333333!!gm This means that bottom and ocean w-level above have a specified "en" value.   Sure ???? 
    334334!!   With thick deep ocean level thickness, this may be quite large, no ??? 
     
    348348         ! 
    349349         IF( ln_isfcav) THEN     ! top boundary   (ocean cavity) 
    350             DO_2D_00_00 
     350            DO_2D( 0, 0, 0, 0 ) 
    351351               itop   = mikt(ji,jj)       ! k   top w-point 
    352352               itopp1 = mikt(ji,jj) + 1   ! k+1 1st w-point below the top one 
     
    366366      CASE ( 1 )             ! Neumman boundary condition 
    367367         !                       
    368          DO_2D_00_00 
     368         DO_2D( 0, 0, 0, 0 ) 
    369369            ibot   = mbkt(ji,jj) + 1      ! k   bottom level of w-point 
    370370            ibotm1 = mbkt(ji,jj)          ! k-1 bottom level of w-point but >=1 
     
    380380         END_2D 
    381381         IF( ln_isfcav) THEN     ! top boundary   (ocean cavity) 
    382             DO_2D_00_00 
     382            DO_2D( 0, 0, 0, 0 ) 
    383383               itop   = mikt(ji,jj)       ! k   top w-point 
    384384               itopp1 = mikt(ji,jj) + 1   ! k+1 1st w-point below the top one 
     
    400400      ! ---------------------------------------------------------- 
    401401      ! 
    402       DO_3D_00_00( 2, jpkm1 ) 
     402      DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    403403         zdiag(ji,jj,jk) = zdiag(ji,jj,jk) - zd_lw(ji,jj,jk) * zd_up(ji,jj,jk-1) / zdiag(ji,jj,jk-1) 
    404404      END_3D 
    405       DO_3D_00_00( 2, jpk ) 
     405      DO_3D( 0, 0, 0, 0, 2, jpk ) 
    406406         zd_lw(ji,jj,jk) = en(ji,jj,jk) - zd_lw(ji,jj,jk) / zdiag(ji,jj,jk-1) * zd_lw(ji,jj,jk-1) 
    407407      END_3D 
    408       DO_3DS_00_00( jpk-1, 2, -1 ) 
     408      DO_3DS( 0, 0, 0, 0, jpk-1, 2, -1 ) 
    409409         en(ji,jj,jk) = ( zd_lw(ji,jj,jk) - zd_up(ji,jj,jk) * en(ji,jj,jk+1) ) / zdiag(ji,jj,jk) 
    410410      END_3D 
     
    421421      ! 
    422422      CASE( 0 )               ! k-kl  (Mellor-Yamada) 
    423          DO_3D_00_00( 2, jpkm1 ) 
     423         DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    424424            psi(ji,jj,jk)  = eb(ji,jj,jk) * hmxl_b(ji,jj,jk) 
    425425         END_3D 
    426426         ! 
    427427      CASE( 1 )               ! k-eps 
    428          DO_3D_00_00( 2, jpkm1 ) 
     428         DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    429429            psi(ji,jj,jk)  = eps(ji,jj,jk) 
    430430         END_3D 
    431431         ! 
    432432      CASE( 2 )               ! k-w 
    433          DO_3D_00_00( 2, jpkm1 ) 
     433         DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    434434            psi(ji,jj,jk)  = SQRT( eb(ji,jj,jk) ) / ( rc0 * hmxl_b(ji,jj,jk) ) 
    435435         END_3D 
    436436         ! 
    437437      CASE( 3 )               ! generic 
    438          DO_3D_00_00( 2, jpkm1 ) 
     438         DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    439439            psi(ji,jj,jk)  = rc02 * eb(ji,jj,jk) * hmxl_b(ji,jj,jk)**rnn  
    440440         END_3D 
     
    449449      ! Warning : after this step, en : right hand side of the matrix 
    450450 
    451       DO_3D_00_00( 2, jpkm1 ) 
     451      DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    452452         ! 
    453453         ! psi / k 
     
    546546         !                      ! en(ibot) = u*^2 / Co2 and hmxl_n(ibot) = vkarmn * r_z0_bot 
    547547         !                      ! Balance between the production and the dissipation terms 
    548          DO_2D_00_00 
     548         DO_2D( 0, 0, 0, 0 ) 
    549549            ibot   = mbkt(ji,jj) + 1      ! k   bottom level of w-point 
    550550            ibotm1 = mbkt(ji,jj)          ! k-1 bottom level of w-point but >=1 
     
    565565      CASE ( 1 )             ! Neumman boundary condition 
    566566         !                       
    567          DO_2D_00_00 
     567         DO_2D( 0, 0, 0, 0 ) 
    568568            ibot   = mbkt(ji,jj) + 1      ! k   bottom level of w-point 
    569569            ibotm1 = mbkt(ji,jj)          ! k-1 bottom level of w-point but >=1 
     
    593593      ! ---------------- 
    594594      ! 
    595       DO_3D_00_00( 2, jpkm1 ) 
     595      DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    596596         zdiag(ji,jj,jk) = zdiag(ji,jj,jk) - zd_lw(ji,jj,jk) * zd_up(ji,jj,jk-1) / zdiag(ji,jj,jk-1) 
    597597      END_3D 
    598       DO_3D_00_00( 2, jpk ) 
     598      DO_3D( 0, 0, 0, 0, 2, jpk ) 
    599599         zd_lw(ji,jj,jk) = psi(ji,jj,jk) - zd_lw(ji,jj,jk) / zdiag(ji,jj,jk-1) * zd_lw(ji,jj,jk-1) 
    600600      END_3D 
    601       DO_3DS_00_00( jpk-1, 2, -1 ) 
     601      DO_3DS( 0, 0, 0, 0, jpk-1, 2, -1 ) 
    602602         psi(ji,jj,jk) = ( zd_lw(ji,jj,jk) - zd_up(ji,jj,jk) * psi(ji,jj,jk+1) ) / zdiag(ji,jj,jk) 
    603603      END_3D 
     
    609609      ! 
    610610      CASE( 0 )               ! k-kl  (Mellor-Yamada) 
    611          DO_3D_00_00( 1, jpkm1 ) 
     611         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
    612612            eps(ji,jj,jk) = rc03 * en(ji,jj,jk) * en(ji,jj,jk) * SQRT( en(ji,jj,jk) ) / MAX( psi(ji,jj,jk), rn_epsmin) 
    613613         END_3D 
    614614         ! 
    615615      CASE( 1 )               ! k-eps 
    616          DO_3D_00_00( 1, jpkm1 ) 
     616         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
    617617            eps(ji,jj,jk) = psi(ji,jj,jk) 
    618618         END_3D 
    619619         ! 
    620620      CASE( 2 )               ! k-w 
    621          DO_3D_00_00( 1, jpkm1 ) 
     621         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
    622622            eps(ji,jj,jk) = rc04 * en(ji,jj,jk) * psi(ji,jj,jk)  
    623623         END_3D 
     
    627627         zex1  =      ( 1.5_wp + rmm/rnn ) 
    628628         zex2  = -1._wp / rnn 
    629          DO_3D_00_00( 1, jpkm1 ) 
     629         DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
    630630            eps(ji,jj,jk) = zcoef * en(ji,jj,jk)**zex1 * psi(ji,jj,jk)**zex2 
    631631         END_3D 
     
    635635      ! Limit dissipation rate under stable stratification 
    636636      ! -------------------------------------------------- 
    637       DO_3D_00_00( 1, jpkm1 ) 
     637      DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 
    638638         ! limitation 
    639639         eps   (ji,jj,jk)  = MAX( eps(ji,jj,jk), rn_epsmin ) 
     
    651651      ! 
    652652      CASE ( 0 , 1 )             ! Galperin or Kantha-Clayson stability functions 
    653          DO_3D_00_00( 2, jpkm1 ) 
     653         DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    654654            ! zcof =  l²/q² 
    655655            zcof = hmxl_b(ji,jj,jk) * hmxl_b(ji,jj,jk) / ( 2._wp*eb(ji,jj,jk) ) 
     
    668668         ! 
    669669      CASE ( 2, 3 )               ! Canuto stability functions 
    670          DO_3D_00_00( 2, jpkm1 ) 
     670         DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 
    671671            ! zcof =  l²/q² 
    672672            zcof = hmxl_b(ji,jj,jk)*hmxl_b(ji,jj,jk) / ( 2._wp * eb(ji,jj,jk) ) 
     
    700700      ! default value, in case jpk > mbkt(ji,jj)+1. Not needed but avoid a bug when looking for undefined values (-fpe0) 
    701701      zstm(:,:,jpk) = 0.   
    702       DO_2D_00_00 
     702      DO_2D( 0, 0, 0, 0 ) 
    703703         zstm(ji,jj,mbkt(ji,jj)+1) = zstm(ji,jj,mbkt(ji,jj)) 
    704704      END_2D 
     
    715715      !     later overwritten by surface/bottom boundaries conditions, so we don't really care of p_avm(:,:1) and p_avm(:,:jpk) 
    716716      !     for zd_lw and zd_up but they have to be defined to avoid a bug when looking for undefined values (-fpe0) 
    717       DO_3D_00_00( 1, jpk ) 
     717      DO_3D( 0, 0, 0, 0, 1, jpk ) 
    718718         zsqen = SQRT( 2._wp * en(ji,jj,jk) ) * hmxl_n(ji,jj,jk) 
    719719         zavt  = zsqen * zstt(ji,jj,jk) 
Note: See TracChangeset for help on using the changeset viewer.