Ignore:
Timestamp:
Dec 18, 2018 12:31:28 PM (3 years ago)
Author:
knoop
Message:

OpenACC port for SPEC

File:
1 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/advec_ws.f90

    r3589 r3634  
    2525! -----------------
    2626! $Id$
     27! OpenACC port for SPEC
     28!
     29! 3589 2018-11-30 15:09:51Z suehring
    2730! Move the control parameter "salsa" from salsa_mod to control_parameters
    2831! (M. Kurppa)
     
    11431146!--    beginning of prognostic_equations.
    11441147       IF ( ws_scheme_mom )  THEN
     1148          !$ACC KERNELS PRESENT(sums_wsus_ws_l, sums_wsvs_ws_l) &
     1149          !$ACC PRESENT(sums_us2_ws_l, sums_vs2_ws_l, sums_ws2_ws_l)
    11451150          sums_wsus_ws_l = 0.0_wp
    11461151          sums_wsvs_ws_l = 0.0_wp
     
    11481153          sums_vs2_ws_l  = 0.0_wp
    11491154          sums_ws2_ws_l  = 0.0_wp
     1155          !$ACC END KERNELS
    11501156       ENDIF
    11511157
    11521158       IF ( ws_scheme_sca )  THEN
     1159          !$ACC KERNELS PRESENT(sums_wspts_ws_l)
    11531160          sums_wspts_ws_l = 0.0_wp
     1161          !$ACC END KERNELS
    11541162          IF ( humidity       )  sums_wsqs_ws_l = 0.0_wp
    11551163          IF ( passive_scalar )  sums_wsss_ws_l = 0.0_wp
     
    33003308
    33013309       CHARACTER (LEN = *), INTENT(IN)    ::  sk_char !< string identifier, used for assign fluxes to the correct dimension in the analysis array
     3310       INTEGER(iwp) ::  sk_num !< integer identifier, used for assign fluxes to the correct dimension in the analysis array
    33023311       
    33033312       INTEGER(iwp) ::  i      !< grid index along x-direction
     
    33203329       REAL(wp) ::  ibit1  !< flag indicating 3rd-order scheme along x-direction
    33213330       REAL(wp) ::  ibit2  !< flag indicating 5th-order scheme along x-direction
     3331#ifdef _OPENACC
     3332       REAL(wp) ::  ibit0_l  !< flag indicating 1st-order scheme along x-direction
     3333       REAL(wp) ::  ibit1_l  !< flag indicating 3rd-order scheme along x-direction
     3334       REAL(wp) ::  ibit2_l  !< flag indicating 5th-order scheme along x-direction
     3335#endif
    33223336       REAL(wp) ::  ibit3  !< flag indicating 1st-order scheme along y-direction
    33233337       REAL(wp) ::  ibit4  !< flag indicating 3rd-order scheme along y-direction
    33243338       REAL(wp) ::  ibit5  !< flag indicating 5th-order scheme along y-direction
     3339#ifdef _OPENACC
     3340       REAL(wp) ::  ibit3_s  !< flag indicating 1st-order scheme along y-direction
     3341       REAL(wp) ::  ibit4_s  !< flag indicating 3rd-order scheme along y-direction
     3342       REAL(wp) ::  ibit5_s  !< flag indicating 5th-order scheme along y-direction
     3343#endif
    33253344       REAL(wp) ::  ibit6  !< flag indicating 1st-order scheme along z-direction
    33263345       REAL(wp) ::  ibit7  !< flag indicating 3rd-order scheme along z-direction
     
    33303349       REAL(wp) ::  flux_d !< 6th-order flux at grid box bottom
    33313350       REAL(wp) ::  u_comp !< advection velocity along x-direction
     3351#ifdef _OPENACC
     3352       REAL(wp) ::  u_comp_l !< advection velocity along x-direction
     3353#endif
    33323354       REAL(wp) ::  v_comp !< advection velocity along y-direction
     3355#ifdef _OPENACC
     3356       REAL(wp) ::  v_comp_s !< advection velocity along y-direction
     3357#endif
    33333358       
    3334        REAL(wp), DIMENSION(nzb:nzt)  ::  diss_n !< discretized artificial dissipation at northward-side of the grid box
    3335        REAL(wp), DIMENSION(nzb:nzt)  ::  diss_r !< discretized artificial dissipation at rightward-side of the grid box
    3336        REAL(wp), DIMENSION(nzb:nzt)  ::  diss_t !< discretized artificial dissipation at rightward-side of the grid box
    3337        REAL(wp), DIMENSION(nzb:nzt)  ::  flux_n !< discretized 6th-order flux at northward-side of the grid box
    3338        REAL(wp), DIMENSION(nzb:nzt)  ::  flux_r !< discretized 6th-order flux at rightward-side of the grid box
    3339        REAL(wp), DIMENSION(nzb:nzt)  ::  flux_t !< discretized 6th-order flux at rightward-side of the grid box
     3359       REAL(wp) ::  diss_n !< discretized artificial dissipation at northward-side of the grid box
     3360       REAL(wp) ::  diss_r !< discretized artificial dissipation at rightward-side of the grid box
     3361       REAL(wp) ::  diss_t !< discretized artificial dissipation at rightward-side of the grid box
     3362       REAL(wp) ::  flux_n !< discretized 6th-order flux at northward-side of the grid box
     3363       REAL(wp) ::  flux_r !< discretized 6th-order flux at rightward-side of the grid box
     3364       REAL(wp) ::  flux_t !< discretized 6th-order flux at rightward-side of the grid box
    33403365       
     3366       REAL(wp) ::  diss_s !< discretized artificial dissipation term at southward-side of the grid box
     3367       REAL(wp) ::  flux_s !< discretized 6th-order flux at northward-side of the grid box
     3368#ifndef _OPENACC
    33413369       REAL(wp), DIMENSION(nzb+1:nzt) ::  swap_diss_y_local !< discretized artificial dissipation term at southward-side of the grid box
    33423370       REAL(wp), DIMENSION(nzb+1:nzt) ::  swap_flux_y_local !< discretized 6th-order flux at northward-side of the grid box
     3371#endif
    33433372       
     3373       REAL(wp) ::  diss_l !< discretized artificial dissipation term at leftward-side of the grid box
     3374       REAL(wp) ::  flux_l !< discretized 6th-order flux at leftward-side of the grid box
     3375#ifndef _OPENACC
    33443376       REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  swap_diss_x_local !< discretized artificial dissipation term at leftward-side of the grid box
    33453377       REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  swap_flux_x_local !< discretized 6th-order flux at leftward-side of the grid box
     3378#endif
    33463379       
    3347 
     3380       SELECT CASE ( sk_char )
     3381
     3382           CASE ( 'pt' )
     3383              sk_num = 1
     3384           CASE ( 'sa' )
     3385              sk_num = 2
     3386           CASE ( 'q' )
     3387              sk_num = 3
     3388           CASE ( 'qc' )
     3389              sk_num = 4
     3390           CASE ( 'qr' )
     3391              sk_num = 5
     3392           CASE ( 'nc' )
     3393              sk_num = 6
     3394           CASE ( 'nr' )
     3395              sk_num = 7
     3396           CASE ( 's' )
     3397              sk_num = 8
     3398           CASE ( 'aerosol_mass', 'aerosol_number', 'salsa_gas' )
     3399              sk_num = 9
     3400
     3401       END SELECT
     3402
     3403#ifndef _OPENACC
    33483404!
    33493405!--    Compute the fluxes for the whole left boundary of the processor domain.
     
    34083464
    34093465       ENDDO
    3410 
     3466#endif
     3467
     3468       !$ACC PARALLEL LOOP COLLAPSE(2) FIRSTPRIVATE(tn, sk_num) &
     3469       !$ACC PRIVATE(i, j, k, k_mm, k_pp, k_ppp) &
     3470       !$ACC PRIVATE(ibit0, ibit1, ibit2, ibit3, ibit4, ibit5) &
     3471       !$ACC PRIVATE(ibit0_l, ibit1_l, ibit2_l) &
     3472       !$ACC PRIVATE(ibit3_s, ibit4_s, ibit5_s) &
     3473       !$ACC PRIVATE(ibit6, ibit7, ibit8) &
     3474       !$ACC PRIVATE(flux_r, diss_r, flux_l, diss_l) &
     3475       !$ACC PRIVATE(flux_n, diss_n, flux_s, diss_s) &
     3476       !$ACC PRIVATE(flux_t, diss_t, flux_d, diss_d) &
     3477       !$ACC PRIVATE(div, u_comp, u_comp_l, v_comp, v_comp_s) &
     3478       !$ACC PRESENT(advc_flags_1) &
     3479       !$ACC PRESENT(sk, u, v, w, u_stokes_zu, v_stokes_zu) &
     3480       !$ACC PRESENT(drho_air, rho_air_zw, ddzw) &
     3481       !$ACC PRESENT(tend) &
     3482       !$ACC PRESENT(hom(nzb+1:nzb_max,1,1:3,0)) &
     3483       !$ACC PRESENT(weight_substep(intermediate_timestep_count)) &
     3484       !$ACC PRESENT(sums_wspts_ws_l, sums_wssas_ws_l) &
     3485       !$ACC PRESENT(sums_wsqs_ws_l, sums_wsqcs_ws_l) &
     3486       !$ACC PRESENT(sums_wsqrs_ws_l, sums_wsncs_ws_l) &
     3487       !$ACC PRESENT(sums_wsnrs_ws_l, sums_wsss_ws_l) &
     3488       !$ACC PRESENT(sums_salsa_ws_l)
    34113489       DO  i = nxl, nxr
    34123490
     3491#ifndef _OPENACC
    34133492          j = nys
    34143493          DO  k = nzb+1, nzb_max
     
    34673546
    34683547          ENDDO
     3548#endif
    34693549
    34703550          DO  j = nys, nyn
    34713551
    3472              flux_t(0) = 0.0_wp
    3473              diss_t(0) = 0.0_wp
    34743552             flux_d    = 0.0_wp
    34753553             diss_d    = 0.0_wp
     
    34813559                ibit0 = REAL( IBITS(advc_flags_1(k,j,i),0,1), KIND = wp )
    34823560
    3483                 u_comp    = u(k,j,i+1) - u_gtrans + u_stokes_zu(k)
    3484                 flux_r(k) = u_comp * (                                        &
     3561                u_comp = u(k,j,i+1) - u_gtrans + u_stokes_zu(k)
     3562                flux_r = u_comp * (                                           &
    34853563                          ( 37.0_wp * ibit2 * adv_sca_5                       &
    34863564                      +      7.0_wp * ibit1 * adv_sca_3                       &
     
    34973575                                     )
    34983576
    3499                 diss_r(k) = -ABS( u_comp ) * (                                &
     3577                diss_r = -ABS( u_comp ) * (                                   &
    35003578                          ( 10.0_wp * ibit2 * adv_sca_5                       &
    35013579                       +     3.0_wp * ibit1 * adv_sca_3                       &
     
    35123590                                             )
    35133591
     3592#ifdef _OPENACC
     3593!
     3594!--             Recompute the left fluxes.
     3595                ibit2_l = REAL( IBITS(advc_flags_1(k,j,i-1),2,1), KIND = wp )
     3596                ibit1_l = REAL( IBITS(advc_flags_1(k,j,i-1),1,1), KIND = wp )
     3597                ibit0_l = REAL( IBITS(advc_flags_1(k,j,i-1),0,1), KIND = wp )
     3598
     3599                u_comp_l = u(k,j,i) - u_gtrans + u_stokes_zu(k)
     3600                flux_l = u_comp_l * (                                         &
     3601                                             ( 37.0_wp * ibit2_l * adv_sca_5  &
     3602                                          +     7.0_wp * ibit1_l * adv_sca_3  &
     3603                                          +              ibit0_l * adv_sca_1  &
     3604                                             ) *                              &
     3605                                          ( sk(k,j,i)   + sk(k,j,i-1)    )    &
     3606                                      -      (  8.0_wp * ibit2_l * adv_sca_5  &
     3607                                          +              ibit1_l * adv_sca_3  &
     3608                                             ) *                              &
     3609                                          ( sk(k,j,i+1) + sk(k,j,i-2)    )    &
     3610                                      +      (           ibit2_l * adv_sca_5  &
     3611                                             ) *                              &
     3612                                          ( sk(k,j,i+2) + sk(k,j,i-3)    )    &
     3613                                               )
     3614
     3615                 diss_l = -ABS( u_comp_l ) * (                                &
     3616                                             ( 10.0_wp * ibit2_l * adv_sca_5  &
     3617                                          +     3.0_wp * ibit1_l * adv_sca_3  &
     3618                                          +              ibit0_l * adv_sca_1  &
     3619                                             ) *                              &
     3620                                          ( sk(k,j,i)   - sk(k,j,i-1) )       &
     3621                                      -      (  5.0_wp * ibit2_l * adv_sca_5  &
     3622                                          +              ibit1_l * adv_sca_3  &
     3623                                             ) *                              &
     3624                                         ( sk(k,j,i+1) - sk(k,j,i-2)  )       &
     3625                                      +      (           ibit2_l * adv_sca_5  &
     3626                                             ) *                              &
     3627                                          ( sk(k,j,i+2) - sk(k,j,i-3) )       &
     3628                                                        )
     3629#else
     3630                flux_l = swap_flux_x_local(k,j)
     3631                diss_l = swap_diss_x_local(k,j)
     3632#endif
     3633
    35143634                ibit5 = REAL( IBITS(advc_flags_1(k,j,i),5,1), KIND = wp )
    35153635                ibit4 = REAL( IBITS(advc_flags_1(k,j,i),4,1), KIND = wp )
    35163636                ibit3 = REAL( IBITS(advc_flags_1(k,j,i),3,1), KIND = wp )
    35173637
    3518                 v_comp    = v(k,j+1,i) - v_gtrans + v_stokes_zu(k)
    3519                 flux_n(k) = v_comp * (                                        &
     3638                v_comp = v(k,j+1,i) - v_gtrans + v_stokes_zu(k)
     3639                flux_n = v_comp * (                                           &
    35203640                          ( 37.0_wp * ibit5 * adv_sca_5                       &
    35213641                       +     7.0_wp * ibit4 * adv_sca_3                       &
     
    35323652                                     )
    35333653
    3534                 diss_n(k) = -ABS( v_comp ) * (                                &
     3654                diss_n = -ABS( v_comp ) * (                                   &
    35353655                          ( 10.0_wp * ibit5 * adv_sca_5                       &
    35363656                       +     3.0_wp * ibit4 * adv_sca_3                       &
     
    35463666                             ( sk(k,j+3,i) - sk(k,j-2,i) )                    &
    35473667                                             )
     3668
     3669#ifdef _OPENACC
     3670!
     3671!--             Recompute the south fluxes.
     3672                ibit5_s = REAL( IBITS(advc_flags_1(k,j-1,i),5,1), KIND = wp )
     3673                ibit4_s = REAL( IBITS(advc_flags_1(k,j-1,i),4,1), KIND = wp )
     3674                ibit3_s = REAL( IBITS(advc_flags_1(k,j-1,i),3,1), KIND = wp )
     3675
     3676                v_comp_s = v(k,j,i) - v_gtrans + v_stokes_zu(k)
     3677                flux_s = v_comp_s * (                                         &
     3678                                             ( 37.0_wp * ibit5_s * adv_sca_5  &
     3679                                          +     7.0_wp * ibit4_s * adv_sca_3  &
     3680                                          +              ibit3_s * adv_sca_1  &
     3681                                             ) *                              &
     3682                                         ( sk(k,j,i)  + sk(k,j-1,i)     )     &
     3683                                       -     (  8.0_wp * ibit5_s * adv_sca_5  &
     3684                                          +              ibit4_s * adv_sca_3  &
     3685                                              ) *                             &
     3686                                         ( sk(k,j+1,i) + sk(k,j-2,i)    )     &
     3687                                      +      (           ibit5_s * adv_sca_5  &
     3688                                             ) *                              &
     3689                                        ( sk(k,j+2,i) + sk(k,j-3,i)     )     &
     3690                                             )
     3691
     3692                diss_s = -ABS( v_comp_s ) * (                                 &
     3693                                             ( 10.0_wp * ibit5_s * adv_sca_5  &
     3694                                          +     3.0_wp * ibit4_s * adv_sca_3  &
     3695                                          +              ibit3_s * adv_sca_1  &
     3696                                             ) *                              &
     3697                                          ( sk(k,j,i)   - sk(k,j-1,i)    )    &
     3698                                      -      (  5.0_wp * ibit5_s * adv_sca_5  &
     3699                                          +              ibit4_s * adv_sca_3  &
     3700                                             ) *                              &
     3701                                          ( sk(k,j+1,i) - sk(k,j-2,i)    )    &
     3702                                      +      (           ibit5_s * adv_sca_5  &
     3703                                             ) *                              &
     3704                                          ( sk(k,j+2,i) - sk(k,j-3,i)    )    &
     3705                                                     )
     3706#else
     3707                flux_s = swap_flux_y_local(k)
     3708                diss_s = swap_diss_y_local(k)
     3709#endif
     3710
    35483711!
    35493712!--             k index has to be modified near bottom and top, else array
     
    35583721
    35593722
    3560                 flux_t(k) = w(k,j,i) * rho_air_zw(k) * (                      &
     3723                flux_t = w(k,j,i) * rho_air_zw(k) * (                         &
    35613724                           ( 37.0_wp * ibit8 * adv_sca_5                      &
    35623725                        +     7.0_wp * ibit7 * adv_sca_3                      &
     
    35723735                                       )
    35733736
    3574                 diss_t(k) = -ABS( w(k,j,i) ) * rho_air_zw(k) * (              &
     3737                diss_t = -ABS( w(k,j,i) ) * rho_air_zw(k) * (                 &
    35753738                           ( 10.0_wp * ibit8 * adv_sca_5                      &
    35763739                        +     3.0_wp * ibit7 * adv_sca_3                      &
     
    36163779
    36173780                tend(k,j,i) = tend(k,j,i) - (                                 &
    3618                         ( flux_r(k) + diss_r(k) - swap_flux_x_local(k,j) -    &
    3619                           swap_diss_x_local(k,j)            ) * ddx           &
    3620                       + ( flux_n(k) + diss_n(k) - swap_flux_y_local(k)   -    &
    3621                           swap_diss_y_local(k)              ) * ddy           &
    3622                       + ( ( flux_t(k) + diss_t(k) ) -                         &
    3623                           ( flux_d    + diss_d    )                           &
    3624                                                     ) * drho_air(k) * ddzw(k) &
     3781                        ( ( flux_r + diss_r )                                 &
     3782                      -   ( flux_l + diss_l ) ) * ddx                         &
     3783                      + ( ( flux_n + diss_n )                                 &
     3784                      -   ( flux_s + diss_s ) ) * ddy                         &
     3785                      + ( ( flux_t + diss_t )                                 &
     3786                      -   ( flux_d + diss_d ) ) * drho_air(k) * ddzw(k)       &
    36253787                                            ) + sk(k,j,i) * div
    36263788
    3627                 swap_flux_y_local(k)   = flux_n(k)
    3628                 swap_diss_y_local(k)   = diss_n(k)
    3629                 swap_flux_x_local(k,j) = flux_r(k)
    3630                 swap_diss_x_local(k,j) = diss_r(k)
    3631                 flux_d                 = flux_t(k)
    3632                 diss_d                 = diss_t(k)
     3789#ifndef _OPENACC
     3790                swap_flux_y_local(k)   = flux_n
     3791                swap_diss_y_local(k)   = diss_n
     3792                swap_flux_x_local(k,j) = flux_r
     3793                swap_diss_x_local(k,j) = diss_r
     3794#endif
     3795                flux_d                 = flux_t
     3796                diss_d                 = diss_t
     3797
     3798!
     3799!--             Evaluation of statistics.
     3800                SELECT CASE ( sk_num )
     3801
     3802                    CASE ( 1 )
     3803                       !$ACC ATOMIC
     3804                       sums_wspts_ws_l(k,tn) = sums_wspts_ws_l(k,tn)           &
     3805                          + ( flux_t                                           &
     3806                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
     3807                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
     3808                            + diss_t                                           &
     3809                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
     3810                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
     3811                            ) * weight_substep(intermediate_timestep_count)
     3812                    CASE ( 2 )
     3813                       !$ACC ATOMIC
     3814                       sums_wssas_ws_l(k,tn) = sums_wssas_ws_l(k,tn)           &
     3815                          + ( flux_t                                           &
     3816                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
     3817                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
     3818                            + diss_t                                           &
     3819                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
     3820                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
     3821                            ) * weight_substep(intermediate_timestep_count)
     3822                    CASE ( 3 )
     3823                       !$ACC ATOMIC
     3824                       sums_wsqs_ws_l(k,tn)  = sums_wsqs_ws_l(k,tn)            &
     3825                          + ( flux_t                                           &
     3826                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
     3827                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
     3828                            + diss_t                                           &
     3829                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
     3830                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
     3831                            ) * weight_substep(intermediate_timestep_count)
     3832                    CASE ( 4 )
     3833                       !$ACC ATOMIC
     3834                       sums_wsqcs_ws_l(k,tn)  = sums_wsqcs_ws_l(k,tn)          &
     3835                          + ( flux_t                                           &
     3836                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
     3837                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
     3838                            + diss_t                                           &
     3839                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
     3840                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
     3841                            ) * weight_substep(intermediate_timestep_count)
     3842                    CASE ( 5 )
     3843                       !$ACC ATOMIC
     3844                       sums_wsqrs_ws_l(k,tn)  = sums_wsqrs_ws_l(k,tn)          &
     3845                          + ( flux_t                                           &
     3846                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
     3847                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
     3848                            + diss_t                                           &
     3849                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
     3850                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
     3851                            ) * weight_substep(intermediate_timestep_count)
     3852                    CASE ( 6 )
     3853                       !$ACC ATOMIC
     3854                       sums_wsncs_ws_l(k,tn)  = sums_wsncs_ws_l(k,tn)          &
     3855                          + ( flux_t                                           &
     3856                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
     3857                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
     3858                            + diss_t                                           &
     3859                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
     3860                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
     3861                            ) * weight_substep(intermediate_timestep_count)
     3862                    CASE ( 7 )
     3863                       !$ACC ATOMIC
     3864                       sums_wsnrs_ws_l(k,tn)  = sums_wsnrs_ws_l(k,tn)          &
     3865                          + ( flux_t                                           &
     3866                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
     3867                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
     3868                            + diss_t                                           &
     3869                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
     3870                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
     3871                            ) * weight_substep(intermediate_timestep_count)
     3872                    CASE ( 8 )
     3873                       !$ACC ATOMIC
     3874                       sums_wsss_ws_l(k,tn)  = sums_wsss_ws_l(k,tn)            &
     3875                          + ( flux_t                                           &
     3876                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
     3877                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
     3878                            + diss_t                                           &
     3879                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
     3880                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
     3881                            ) * weight_substep(intermediate_timestep_count)
     3882                    CASE ( 9 )
     3883                        !$ACC ATOMIC
     3884                        sums_salsa_ws_l(k,tn)  = sums_salsa_ws_l(k,tn)         &
     3885                          + ( flux_t                                           &
     3886                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
     3887                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
     3888                            + diss_t                                           &
     3889                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
     3890                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
     3891                            ) * weight_substep(intermediate_timestep_count)
     3892
     3893                END SELECT
    36333894
    36343895             ENDDO
     
    36363897             DO  k = nzb_max+1, nzt
    36373898
    3638                 u_comp    = u(k,j,i+1) - u_gtrans + u_stokes_zu(k)
    3639                 flux_r(k) = u_comp * (                                        &
     3899                u_comp = u(k,j,i+1) - u_gtrans + u_stokes_zu(k)
     3900                flux_r = u_comp * (                                           &
    36403901                      37.0_wp * ( sk(k,j,i+1) + sk(k,j,i)   )                 &
    36413902                    -  8.0_wp * ( sk(k,j,i+2) + sk(k,j,i-1) )                 &
    36423903                    +           ( sk(k,j,i+3) + sk(k,j,i-2) ) ) * adv_sca_5
    3643                 diss_r(k) = -ABS( u_comp ) * (                                &
     3904                diss_r = -ABS( u_comp ) * (                                   &
    36443905                      10.0_wp * ( sk(k,j,i+1) - sk(k,j,i)   )                 &
    36453906                    -  5.0_wp * ( sk(k,j,i+2) - sk(k,j,i-1) )                 &
    36463907                    +           ( sk(k,j,i+3) - sk(k,j,i-2) ) ) * adv_sca_5
    36473908
    3648                 v_comp    = v(k,j+1,i) - v_gtrans + v_stokes_zu(k)
    3649                 flux_n(k) = v_comp * (                                        &
     3909#ifdef _OPENACC
     3910!
     3911!--             Recompute the left fluxes.
     3912                u_comp_l = u(k,j,i) - u_gtrans + u_stokes_zu(k)
     3913                flux_l = u_comp_l * (                                         &
     3914                                      37.0_wp * ( sk(k,j,i)   + sk(k,j,i-1) ) &
     3915                                    -  8.0_wp * ( sk(k,j,i+1) + sk(k,j,i-2) ) &
     3916                                    +           ( sk(k,j,i+2) + sk(k,j,i-3) ) &
     3917                                               ) * adv_sca_5
     3918
     3919                diss_l = -ABS( u_comp_l ) * (                                 &
     3920                                      10.0_wp * ( sk(k,j,i)   - sk(k,j,i-1) ) &
     3921                                    -  5.0_wp * ( sk(k,j,i+1) - sk(k,j,i-2) ) &
     3922                                    +           ( sk(k,j,i+2) - sk(k,j,i-3) ) &
     3923                                                       ) * adv_sca_5
     3924#else
     3925                flux_l = swap_flux_x_local(k,j)
     3926                diss_l = swap_diss_x_local(k,j)
     3927#endif
     3928
     3929                v_comp = v(k,j+1,i) - v_gtrans + v_stokes_zu(k)
     3930                flux_n = v_comp * (                                           &
    36503931                      37.0_wp * ( sk(k,j+1,i) + sk(k,j,i)   )                 &
    36513932                    -  8.0_wp * ( sk(k,j+2,i) + sk(k,j-1,i) )                 &
    36523933                    +           ( sk(k,j+3,i) + sk(k,j-2,i) ) ) * adv_sca_5
    3653                 diss_n(k) = -ABS( v_comp ) * (                                &
     3934                diss_n = -ABS( v_comp ) * (                                   &
    36543935                      10.0_wp * ( sk(k,j+1,i) - sk(k,j,i)   )                 &
    36553936                    -  5.0_wp * ( sk(k,j+2,i) - sk(k,j-1,i) )                 &
    36563937                    +           ( sk(k,j+3,i) - sk(k,j-2,i) ) ) * adv_sca_5
     3938
     3939#ifdef _OPENACC
     3940!
     3941!--             Recompute the south fluxes.
     3942                v_comp_s = v(k,j,i) - v_gtrans + v_stokes_zu(k)
     3943                flux_s = v_comp_s * (                                        &
     3944                                    37.0_wp * ( sk(k,j,i)   + sk(k,j-1,i) )  &
     3945                                  -  8.0_wp * ( sk(k,j+1,i) + sk(k,j-2,i) )  &
     3946                                  +           ( sk(k,j+2,i) + sk(k,j-3,i) )  &
     3947                                             ) * adv_sca_5
     3948                diss_s = -ABS( v_comp_s ) * (                                &
     3949                                    10.0_wp * ( sk(k,j,i)   - sk(k,j-1,i) )  &
     3950                                  -  5.0_wp * ( sk(k,j+1,i) - sk(k,j-2,i) )  &
     3951                                  +             sk(k,j+2,i) - sk(k,j-3,i)    &
     3952                                                      ) * adv_sca_5
     3953#else
     3954                flux_s = swap_flux_y_local(k)
     3955                diss_s = swap_diss_y_local(k)
     3956#endif
     3957
    36573958!
    36583959!--             k index has to be modified near bottom and top, else array
     
    36673968
    36683969
    3669                 flux_t(k) = w(k,j,i) * rho_air_zw(k) * (                      &
     3970                flux_t = w(k,j,i) * rho_air_zw(k) * (                      &
    36703971                           ( 37.0_wp * ibit8 * adv_sca_5                      &
    36713972                        +     7.0_wp * ibit7 * adv_sca_3                      &
     
    36813982                                       )
    36823983
    3683                 diss_t(k) = -ABS( w(k,j,i) ) * rho_air_zw(k) * (              &
     3984                diss_t = -ABS( w(k,j,i) ) * rho_air_zw(k) * (              &
    36843985                           ( 10.0_wp * ibit8 * adv_sca_5                      &
    36853986                        +     3.0_wp * ibit7 * adv_sca_3                      &
     
    37064007
    37074008                tend(k,j,i) = tend(k,j,i) - (                                 &
    3708                         ( flux_r(k) + diss_r(k) - swap_flux_x_local(k,j) -    &
    3709                           swap_diss_x_local(k,j)            ) * ddx           &
    3710                       + ( flux_n(k) + diss_n(k) - swap_flux_y_local(k)   -    &
    3711                           swap_diss_y_local(k)              ) * ddy           &
    3712                       + ( ( flux_t(k) + diss_t(k) ) -                         &
    3713                           ( flux_d    + diss_d    )                           &
    3714                                                     ) * drho_air(k) * ddzw(k) &
     4009                        ( ( flux_r + diss_r )                                 &
     4010                      -   ( flux_l + diss_l ) ) * ddx                         &
     4011                      + ( ( flux_n + diss_n )                                 &
     4012                      -   ( flux_s + diss_s ) ) * ddy                         &
     4013                      + ( ( flux_t + diss_t )                                 &
     4014                      -   ( flux_d + diss_d ) ) * drho_air(k) * ddzw(k)       &
    37154015                                            ) + sk(k,j,i) * div
    37164016
    3717                 swap_flux_y_local(k)   = flux_n(k)
    3718                 swap_diss_y_local(k)   = diss_n(k)
    3719                 swap_flux_x_local(k,j) = flux_r(k)
    3720                 swap_diss_x_local(k,j) = diss_r(k)
    3721                 flux_d                 = flux_t(k)
    3722                 diss_d                 = diss_t(k)
    3723 
    3724              ENDDO
    3725 !
    3726 !--          Evaluation of statistics.
    3727              SELECT CASE ( sk_char )
    3728 
    3729                  CASE ( 'pt' )
    3730                     DO  k = nzb, nzt
     4017#ifndef _OPENACC
     4018                swap_flux_y_local(k)   = flux_n
     4019                swap_diss_y_local(k)   = diss_n
     4020                swap_flux_x_local(k,j) = flux_r
     4021                swap_diss_x_local(k,j) = diss_r
     4022#endif
     4023                flux_d                 = flux_t
     4024                diss_d                 = diss_t
     4025
     4026!
     4027!--             Evaluation of statistics.
     4028                SELECT CASE ( sk_num )
     4029
     4030                    CASE ( 1 )
     4031                       !$ACC ATOMIC
    37314032                       sums_wspts_ws_l(k,tn) = sums_wspts_ws_l(k,tn)           &
    3732                           + ( flux_t(k)                                        &
     4033                          + ( flux_t                                           &
    37334034                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
    37344035                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
    3735                             + diss_t(k)                                        &
     4036                            + diss_t                                           &
    37364037                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
    37374038                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
    37384039                            ) * weight_substep(intermediate_timestep_count)
    3739                     ENDDO
    3740                  CASE ( 'sa' )
    3741                     DO  k = nzb, nzt
     4040                    CASE ( 2 )
     4041                       !$ACC ATOMIC
    37424042                       sums_wssas_ws_l(k,tn) = sums_wssas_ws_l(k,tn)           &
    3743                           + ( flux_t(k)                                        &
     4043                          + ( flux_t                                           &
    37444044                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
    37454045                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
    3746                             + diss_t(k)                                        &
     4046                            + diss_t                                           &
    37474047                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
    37484048                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
    37494049                            ) * weight_substep(intermediate_timestep_count)
    3750                     ENDDO
    3751                  CASE ( 'q' )
    3752                     DO  k = nzb, nzt
     4050                    CASE ( 3 )
     4051                       !$ACC ATOMIC
    37534052                       sums_wsqs_ws_l(k,tn)  = sums_wsqs_ws_l(k,tn)            &
    3754                           + ( flux_t(k)                                        &
     4053                          + ( flux_t                                           &
    37554054                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
    37564055                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
    3757                             + diss_t(k)                                        &
     4056                            + diss_t                                           &
    37584057                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
    37594058                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
    37604059                            ) * weight_substep(intermediate_timestep_count)
    3761                     ENDDO
    3762                  CASE ( 'qc' )
    3763                     DO  k = nzb, nzt
     4060                    CASE ( 4 )
     4061                       !$ACC ATOMIC
    37644062                       sums_wsqcs_ws_l(k,tn)  = sums_wsqcs_ws_l(k,tn)          &
    3765                           + ( flux_t(k)                                        &
     4063                          + ( flux_t                                           &
    37664064                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
    37674065                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
    3768                             + diss_t(k)                                        &
     4066                            + diss_t                                           &
    37694067                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
    37704068                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
    37714069                            ) * weight_substep(intermediate_timestep_count)
    3772                     ENDDO
    3773                  CASE ( 'qr' )
    3774                     DO  k = nzb, nzt
     4070                    CASE ( 5 )
     4071                       !$ACC ATOMIC
    37754072                       sums_wsqrs_ws_l(k,tn)  = sums_wsqrs_ws_l(k,tn)          &
    3776                           + ( flux_t(k)                                        &
     4073                          + ( flux_t                                           &
    37774074                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
    37784075                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
    3779                             + diss_t(k)                                        &
     4076                            + diss_t                                           &
    37804077                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
    37814078                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
    37824079                            ) * weight_substep(intermediate_timestep_count)
    3783                     ENDDO
    3784                  CASE ( 'nc' )
    3785                     DO  k = nzb, nzt
     4080                    CASE ( 6 )
     4081                       !$ACC ATOMIC
    37864082                       sums_wsncs_ws_l(k,tn)  = sums_wsncs_ws_l(k,tn)          &
    3787                           + ( flux_t(k)                                        &
     4083                          + ( flux_t                                           &
    37884084                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
    37894085                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
    3790                             + diss_t(k)                                        &
     4086                            + diss_t                                           &
    37914087                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
    37924088                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
    37934089                            ) * weight_substep(intermediate_timestep_count)
    3794                     ENDDO
    3795                  CASE ( 'nr' )
    3796                     DO  k = nzb, nzt
     4090                    CASE ( 7 )
     4091                       !$ACC ATOMIC
    37974092                       sums_wsnrs_ws_l(k,tn)  = sums_wsnrs_ws_l(k,tn)          &
    3798                           + ( flux_t(k)                                        &
     4093                          + ( flux_t                                           &
    37994094                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
    38004095                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
    3801                             + diss_t(k)                                        &
     4096                            + diss_t                                           &
    38024097                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
    38034098                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
    38044099                            ) * weight_substep(intermediate_timestep_count)
    3805                     ENDDO
    3806                  CASE ( 's' )
    3807                     DO  k = nzb, nzt
     4100                    CASE ( 8 )
     4101                       !$ACC ATOMIC
    38084102                       sums_wsss_ws_l(k,tn)  = sums_wsss_ws_l(k,tn)            &
    3809                           + ( flux_t(k)                                        &
     4103                          + ( flux_t                                           &
    38104104                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
    38114105                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
    3812                             + diss_t(k)                                        &
     4106                            + diss_t                                           &
    38134107                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
    38144108                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
    38154109                            ) * weight_substep(intermediate_timestep_count)
    3816                     ENDDO 
    3817      
    3818                  CASE ( 'aerosol_mass', 'aerosol_number', 'salsa_gas' )
    3819                      DO  k = nzb, nzt
     4110                    CASE ( 9 )
     4111                        !$ACC ATOMIC
    38204112                        sums_salsa_ws_l(k,tn)  = sums_salsa_ws_l(k,tn)         &
    3821                           + ( flux_t(k)                                        &
     4113                          + ( flux_t                                           &
    38224114                                / ( w(k,j,i) + SIGN( 1.0E-20_wp, w(k,j,i) ) )  &
    38234115                                * ( w(k,j,i) - hom(k,1,3,0)                 )  &
    3824                             + diss_t(k)                                        &
     4116                            + diss_t                                           &
    38254117                                / ( ABS(w(k,j,i)) + 1.0E-20_wp              )  &
    38264118                                *   ABS(w(k,j,i) - hom(k,1,3,0)             )  &
    38274119                            ) * weight_substep(intermediate_timestep_count)
    3828                      ENDDO                   
    3829                                    
    3830 
    3831               END SELECT
     4120
     4121                END SELECT
     4122
     4123             ENDDO
    38324124
    38334125         ENDDO
     
    38744166       REAL(wp)    ::  ibit10 !< flag indicating 3rd-order scheme along x-direction
    38754167       REAL(wp)    ::  ibit11 !< flag indicating 5th-order scheme along x-direction
     4168#ifdef _OPENACC
     4169       REAL(wp)    ::  ibit9_l  !< flag indicating 1st-order scheme along x-direction
     4170       REAL(wp)    ::  ibit10_l !< flag indicating 3rd-order scheme along x-direction
     4171       REAL(wp)    ::  ibit11_l !< flag indicating 5th-order scheme along x-direction
     4172#endif
    38764173       REAL(wp)    ::  ibit12 !< flag indicating 1st-order scheme along y-direction
    38774174       REAL(wp)    ::  ibit13 !< flag indicating 3rd-order scheme along y-direction
    38784175       REAL(wp)    ::  ibit14 !< flag indicating 5th-order scheme along y-direction
     4176#ifdef _OPENACC
     4177       REAL(wp)    ::  ibit12_s !< flag indicating 1st-order scheme along y-direction
     4178       REAL(wp)    ::  ibit13_s !< flag indicating 3rd-order scheme along y-direction
     4179       REAL(wp)    ::  ibit14_s !< flag indicating 5th-order scheme along y-direction
     4180#endif
    38794181       REAL(wp)    ::  ibit15 !< flag indicating 1st-order scheme along z-direction
    38804182       REAL(wp)    ::  ibit16 !< flag indicating 3rd-order scheme along z-direction
     
    38864188       REAL(wp)    ::  gv     !< Galilei-transformation velocity along y
    38874189       REAL(wp)    ::  v_comp !< advection velocity along y
     4190#ifdef _OPENACC
     4191       REAL(wp)    ::  v_comp_s !< advection velocity along y
     4192#endif
    38884193       REAL(wp)    ::  w_comp !< advection velocity along z
    38894194       
     4195       REAL(wp)    :: diss_s  !< discretized artificial dissipation at southward-side of the grid box
     4196       REAL(wp)    :: flux_s  !< discretized 6th-order flux at southward-side of the grid box
     4197#ifndef _OPENACC
    38904198       REAL(wp), DIMENSION(nzb+1:nzt) ::  swap_diss_y_local_u !< discretized artificial dissipation at southward-side of the grid box
    38914199       REAL(wp), DIMENSION(nzb+1:nzt) ::  swap_flux_y_local_u !< discretized 6th-order flux at southward-side of the grid box
     4200#endif
    38924201       
     4202       REAL(wp)    :: diss_l  !< discretized artificial dissipation at leftward-side of the grid box
     4203       REAL(wp)    :: flux_l  !< discretized 6th-order flux at leftward-side of the grid box
     4204#ifndef _OPENACC
    38934205       REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  swap_diss_x_local_u !< discretized artificial dissipation at leftward-side of the grid box
    38944206       REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  swap_flux_x_local_u !< discretized 6th-order flux at leftward-side of the grid box
     4207#endif
    38954208       
    3896        REAL(wp), DIMENSION(nzb:nzt) ::  diss_n !< discretized artificial dissipation at northward-side of the grid box
    3897        REAL(wp), DIMENSION(nzb:nzt) ::  diss_r !< discretized artificial dissipation at leftward-side of the grid box
    3898        REAL(wp), DIMENSION(nzb:nzt) ::  diss_t !< discretized artificial dissipation at top of the grid box
    3899        REAL(wp), DIMENSION(nzb:nzt) ::  flux_n !< discretized 6th-order flux at northward-side of the grid box
    3900        REAL(wp), DIMENSION(nzb:nzt) ::  flux_r !< discretized 6th-order flux at rightward-side of the grid box
    3901        REAL(wp), DIMENSION(nzb:nzt) ::  flux_t !< discretized 6th-order flux at top of the grid box
    3902        REAL(wp), DIMENSION(nzb:nzt) ::  u_comp !< advection velocity along x
     4209       REAL(wp) ::  diss_n !< discretized artificial dissipation at northward-side of the grid box
     4210       REAL(wp) ::  diss_r !< discretized artificial dissipation at leftward-side of the grid box
     4211       REAL(wp) ::  diss_t !< discretized artificial dissipation at top of the grid box
     4212       REAL(wp) ::  flux_n !< discretized 6th-order flux at northward-side of the grid box
     4213       REAL(wp) ::  flux_r !< discretized 6th-order flux at rightward-side of the grid box
     4214       REAL(wp) ::  flux_t !< discretized 6th-order flux at top of the grid box
     4215       REAL(wp) ::  u_comp !< advection velocity along x
     4216#ifdef _OPENACC
     4217       REAL(wp)    ::  u_comp_l !<
     4218#endif
    39034219 
    39044220       gu = 2.0_wp * u_gtrans
    39054221       gv = 2.0_wp * v_gtrans
    39064222
     4223#ifndef _OPENACC
    39074224!
    39084225!--    Compute the fluxes for the whole left boundary of the processor domain.
     
    39154232             ibit9  = REAL( IBITS(advc_flags_1(k,j,i-1),9,1),  KIND = wp )
    39164233
    3917              u_comp(k)                = u(k,j,i) + u(k,j,i-1) - gu
    3918              swap_flux_x_local_u(k,j) = u_comp(k) * (                          &
     4234             u_comp                   = u(k,j,i) + u(k,j,i-1) - gu
     4235             swap_flux_x_local_u(k,j) = u_comp * (                             &
    39194236                                       ( 37.0_wp * ibit11 * adv_mom_5          &
    39204237                                    +     7.0_wp * ibit10 * adv_mom_3          &
     
    39314248                                                   )
    39324249
    3933               swap_diss_x_local_u(k,j) = - ABS( u_comp(k) ) * (                &
     4250              swap_diss_x_local_u(k,j) = - ABS( u_comp ) * (                   &
    39344251                                       ( 10.0_wp * ibit11 * adv_mom_5          &
    39354252                                    +     3.0_wp * ibit10 * adv_mom_3          &
     
    39504267          DO  k = nzb_max+1, nzt
    39514268
    3952              u_comp(k)         = u(k,j,i) + u(k,j,i-1) - gu
    3953              swap_flux_x_local_u(k,j) = u_comp(k) * (                          &
     4269             u_comp            = u(k,j,i) + u(k,j,i-1) - gu
     4270             swap_flux_x_local_u(k,j) = u_comp * (                             &
    39544271                             37.0_wp * ( u(k,j,i) + u(k,j,i-1)   )             &
    39554272                           -  8.0_wp * ( u(k,j,i+1) + u(k,j,i-2) )             &
    39564273                           +           ( u(k,j,i+2) + u(k,j,i-3) ) ) * adv_mom_5
    3957              swap_diss_x_local_u(k,j) = - ABS(u_comp(k)) * (                   &
     4274             swap_diss_x_local_u(k,j) = - ABS(u_comp) * (                      &
    39584275                             10.0_wp * ( u(k,j,i) - u(k,j,i-1)   )             &
    39594276                           -  5.0_wp * ( u(k,j,i+1) - u(k,j,i-2) )             &
     
    39624279          ENDDO
    39634280       ENDDO
    3964 
     4281#endif
     4282
     4283       !$ACC PARALLEL LOOP COLLAPSE(2) FIRSTPRIVATE(tn, gu, gv) &
     4284       !$ACC PRIVATE(i, j, k, k_mm, k_pp, k_ppp) &
     4285       !$ACC PRIVATE(ibit9, ibit10, ibit11, ibit12, ibit13, ibit14) &
     4286       !$ACC PRIVATE(ibit9_l, ibit10_l, ibit11_l) &
     4287       !$ACC PRIVATE(ibit12_s, ibit13_s, ibit14_s) &
     4288       !$ACC PRIVATE(ibit15, ibit16, ibit17) &
     4289       !$ACC PRIVATE(flux_r, diss_r, flux_l, diss_l) &
     4290       !$ACC PRIVATE(flux_n, diss_n, flux_s, diss_s) &
     4291       !$ACC PRIVATE(flux_t, diss_t, flux_d, diss_d) &
     4292       !$ACC PRIVATE(div, u_comp, u_comp_l, v_comp, v_comp_s, w_comp) &
     4293       !$ACC PRESENT(advc_flags_1) &
     4294       !$ACC PRESENT(u, v, w) &
     4295       !$ACC PRESENT(drho_air, rho_air_zw, ddzw) &
     4296       !$ACC PRESENT(tend) &
     4297       !$ACC PRESENT(hom(nzb+1:nzb_max,1,1:3,0)) &
     4298       !$ACC PRESENT(weight_substep(intermediate_timestep_count)) &
     4299       !$ACC PRESENT(sums_us2_ws_l, sums_wsus_ws_l)
    39654300       DO i = nxlu, nxr
     4301#ifndef _OPENACC
    39664302!       
    39674303!--       The following loop computes the fluxes for the south boundary points
     
    40194355
    40204356          ENDDO
     4357#endif
     4358
    40214359!
    40224360!--       Computation of interior fluxes and tendency terms
    40234361          DO  j = nys, nyn
    40244362
    4025              flux_t(0) = 0.0_wp
    4026              diss_t(0) = 0.0_wp
    40274363             flux_d    = 0.0_wp
    40284364             diss_d    = 0.0_wp
     
    40344370                ibit9  = REAL( IBITS(advc_flags_1(k,j,i),9,1),  KIND = wp )
    40354371
    4036                 u_comp(k) = u(k,j,i+1) + u(k,j,i)
    4037                 flux_r(k) = ( u_comp(k) - gu ) * (                           &
     4372                u_comp = u(k,j,i+1) + u(k,j,i)
     4373                flux_r = ( u_comp - gu ) * (                                 &
    40384374                          ( 37.0_wp * ibit11 * adv_mom_5                     &
    40394375                       +     7.0_wp * ibit10 * adv_mom_3                     &
     
    40504386                                                 )
    40514387
    4052                 diss_r(k) = - ABS( u_comp(k) - gu ) * (                      &
     4388                diss_r = - ABS( u_comp - gu ) * (                            &
    40534389                          ( 10.0_wp * ibit11 * adv_mom_5                     &
    40544390                       +     3.0_wp * ibit10 * adv_mom_3                     &
     
    40654401                                                     )
    40664402
     4403#ifdef _OPENACC
     4404!
     4405!--             Recompute the left fluxes.
     4406                ibit11_l = REAL( IBITS(advc_flags_1(k,j,i-1),11,1), KIND = wp )
     4407                ibit10_l = REAL( IBITS(advc_flags_1(k,j,i-1),10,1), KIND = wp )
     4408                ibit9_l  = REAL( IBITS(advc_flags_1(k,j,i-1),9,1),  KIND = wp )
     4409
     4410                u_comp_l = u(k,j,i) + u(k,j,i-1) - gu
     4411                flux_l   = u_comp_l * (                                        &
     4412                                       ( 37.0_wp * ibit11_l * adv_mom_5           &
     4413                                    +     7.0_wp * ibit10_l * adv_mom_3           &
     4414                                    +              ibit9_l  * adv_mom_1           &
     4415                                       ) *                                     &
     4416                                     ( u(k,j,i)   + u(k,j,i-1) )               &
     4417                                -      (  8.0_wp * ibit11_l * adv_mom_5           &
     4418                                    +              ibit10_l * adv_mom_3           &
     4419                                       ) *                                     &
     4420                                     ( u(k,j,i+1) + u(k,j,i-2) )               &
     4421                                +      (           ibit11_l * adv_mom_5           &
     4422                                       ) *                                     &
     4423                                     ( u(k,j,i+2) + u(k,j,i-3) )               &
     4424                                                   )
     4425
     4426                diss_l   = - ABS( u_comp_l ) * (                               &
     4427                                       ( 10.0_wp * ibit11_l * adv_mom_5           &
     4428                                    +     3.0_wp * ibit10_l * adv_mom_3           &
     4429                                    +              ibit9_l  * adv_mom_1           &
     4430                                       ) *                                     &
     4431                                     ( u(k,j,i)   - u(k,j,i-1) )               &
     4432                                -      (  5.0_wp * ibit11_l * adv_mom_5           &
     4433                                    +              ibit10_l * adv_mom_3           &
     4434                                       ) *                                     &
     4435                                     ( u(k,j,i+1) - u(k,j,i-2) )               &
     4436                                +      (           ibit11_l * adv_mom_5           &
     4437                                       ) *                                     &
     4438                                     ( u(k,j,i+2) - u(k,j,i-3) )               &
     4439                                                             )
     4440#else
     4441                flux_l = swap_flux_x_local_u(k,j)
     4442                diss_l = swap_diss_x_local_u(k,j)
     4443#endif
     4444
    40674445                ibit14 = REAL( IBITS(advc_flags_1(k,j,i),14,1), KIND = wp )
    40684446                ibit13 = REAL( IBITS(advc_flags_1(k,j,i),13,1), KIND = wp )
    40694447                ibit12 = REAL( IBITS(advc_flags_1(k,j,i),12,1), KIND = wp )
    40704448
    4071                 v_comp    = v(k,j+1,i) + v(k,j+1,i-1) - gv
    4072                 flux_n(k) = v_comp * (                                       &
     4449                v_comp = v(k,j+1,i) + v(k,j+1,i-1) - gv
     4450                flux_n = v_comp * (                                          &
    40734451                          ( 37.0_wp * ibit14 * adv_mom_5                     &
    40744452                       +     7.0_wp * ibit13 * adv_mom_3                     &
     
    40854463                                                 )
    40864464
    4087                 diss_n(k) = - ABS ( v_comp ) * (                             &
     4465                diss_n = - ABS ( v_comp ) * (                                &
    40884466                          ( 10.0_wp * ibit14 * adv_mom_5                     &
    40894467                       +     3.0_wp * ibit13 * adv_mom_3                     &
     
    40994477                                 ( u(k,j+3,i) - u(k,j-2,i) )                 &
    41004478                                                      )
     4479
     4480#ifdef _OPENACC
     4481!
     4482!--             Recompute the south fluxes.
     4483                ibit14_s = REAL( IBITS(advc_flags_1(k,j-1,i),14,1), KIND = wp )
     4484                ibit13_s = REAL( IBITS(advc_flags_1(k,j-1,i),13,1), KIND = wp )
     4485                ibit12_s = REAL( IBITS(advc_flags_1(k,j-1,i),12,1), KIND = wp )
     4486
     4487                v_comp_s = v(k,j,i) + v(k,j,i-1) - gv
     4488                flux_s   = v_comp_s * (                                      &
     4489                                   ( 37.0_wp * ibit14_s * adv_mom_5             &
     4490                                +     7.0_wp * ibit13_s * adv_mom_3             &
     4491                                +              ibit12_s * adv_mom_1             &
     4492                                   ) *                                       &
     4493                                     ( u(k,j,i)   + u(k,j-1,i) )             &
     4494                            -      (  8.0_wp * ibit14_s * adv_mom_5             &
     4495                            +                  ibit13_s * adv_mom_3             &
     4496                                   ) *                                       &
     4497                                     ( u(k,j+1,i) + u(k,j-2,i) )             &
     4498                        +      (               ibit14_s * adv_mom_5             &
     4499                               ) *                                           &
     4500                                     ( u(k,j+2,i) + u(k,j-3,i) )             &
     4501                                               )
     4502
     4503                diss_s   = - ABS ( v_comp_s ) * (                            &
     4504                                   ( 10.0_wp * ibit14_s * adv_mom_5              &
     4505                                +     3.0_wp * ibit13_s * adv_mom_3              &
     4506                                +              ibit12_s * adv_mom_1              &
     4507                                   ) *                                        &
     4508                                     ( u(k,j,i)   - u(k,j-1,i) )              &
     4509                            -      (  5.0_wp * ibit14_s * adv_mom_5              &
     4510                                +              ibit13_s * adv_mom_3              &
     4511                                   ) *                                        &
     4512                                     ( u(k,j+1,i) - u(k,j-2,i) )              &
     4513                            +      (           ibit14_s * adv_mom_5              &
     4514                                   ) *                                        &
     4515                                     ( u(k,j+2,i) - u(k,j-3,i) )              &
     4516                                                         )
     4517#else
     4518                flux_s = swap_flux_y_local_u(k)
     4519                diss_s = swap_diss_y_local_u(k)
     4520#endif
     4521
    41014522!
    41024523!--             k index has to be modified near bottom and top, else array
     
    41104531                k_mm  = k - 2 * ibit17
    41114532
    4112                 w_comp    = w(k,j,i) + w(k,j,i-1)
    4113                 flux_t(k) = w_comp * rho_air_zw(k) * (                       &
     4533                w_comp = w(k,j,i) + w(k,j,i-1)
     4534                flux_t = w_comp * rho_air_zw(k) * (                          &
    41144535                          ( 37.0_wp * ibit17 * adv_mom_5                     &
    41154536                       +     7.0_wp * ibit16 * adv_mom_3                     &
     
    41264547                                      )
    41274548
    4128                 diss_t(k) = - ABS( w_comp ) * rho_air_zw(k) * (              &
     4549                diss_t = - ABS( w_comp ) * rho_air_zw(k) * (                 &
    41294550                          ( 10.0_wp * ibit17 * adv_mom_5                     &
    41304551                       +     3.0_wp * ibit16 * adv_mom_3                     &
     
    41444565!--             correction is needed to overcome numerical instabilities caused
    41454566!--             by a not sufficient reduction of divergences near topography.
    4146                 div = ( ( u_comp(k) * ( ibit9 + ibit10 + ibit11 )             &
     4567                div = ( ( u_comp * ( ibit9 + ibit10 + ibit11 )                &
    41474568                - ( u(k,j,i)   + u(k,j,i-1)   )                               &
    41484569                                    * (                                       &
     
    41734594
    41744595                tend(k,j,i) = tend(k,j,i) - (                                  &
    4175                  ( flux_r(k) + diss_r(k)                                       &
    4176                -   swap_flux_x_local_u(k,j) - swap_diss_x_local_u(k,j) ) * ddx &
    4177                + ( flux_n(k) + diss_n(k)                                       &
    4178                -   swap_flux_y_local_u(k)   - swap_diss_y_local_u(k)   ) * ddy &
    4179                + ( ( flux_t(k) + diss_t(k) )                                   &
    4180                -   ( flux_d    + diss_d    )                                   &
    4181                                                     ) * drho_air(k) * ddzw(k)  &
     4596                 ( ( flux_r + diss_r )                                         &
     4597               -   ( flux_l + diss_l ) ) * ddx                                 &
     4598               + ( ( flux_n + diss_n )                                         &
     4599               -   ( flux_s + diss_s ) ) * ddy                                 &
     4600               + ( ( flux_t + diss_t )                                         &
     4601               -   ( flux_d + diss_d ) ) * drho_air(k) * ddzw(k)               &
    41824602                                           ) + div * u(k,j,i)
    41834603
    4184                 swap_flux_x_local_u(k,j) = flux_r(k)
    4185                 swap_diss_x_local_u(k,j) = diss_r(k)
    4186                 swap_flux_y_local_u(k)   = flux_n(k)
    4187                 swap_diss_y_local_u(k)   = diss_n(k)
    4188                 flux_d                   = flux_t(k)
    4189                 diss_d                   = diss_t(k)
     4604#ifndef _OPENACC
     4605                swap_flux_x_local_u(k,j) = flux_r
     4606                swap_diss_x_local_u(k,j) = diss_r
     4607                swap_flux_y_local_u(k)   = flux_n
     4608                swap_diss_y_local_u(k)   = diss_n
     4609#endif
     4610                flux_d                   = flux_t
     4611                diss_d                   = diss_t
    41904612!
    41914613!--             Statistical Evaluation of u'u'. The factor has to be applied
    41924614!--             for right evaluation when gallilei_trans = .T. .
     4615                !$ACC ATOMIC
    41934616                sums_us2_ws_l(k,tn) = sums_us2_ws_l(k,tn)                      &
    4194                 + ( flux_r(k)                                                  &
    4195                     * ( u_comp(k) - 2.0_wp * hom(k,1,1,0)                   )  &
    4196                     / ( u_comp(k) - gu + SIGN( 1.0E-20_wp, u_comp(k) - gu ) )  &
    4197                   + diss_r(k)                                                  &
    4198                     *   ABS( u_comp(k) - 2.0_wp * hom(k,1,1,0)              )  &
    4199                     / ( ABS( u_comp(k) - gu ) + 1.0E-20_wp                  )  &
     4617                + ( flux_r                                                     &
     4618                    * ( u_comp - 2.0_wp * hom(k,1,1,0)                   )     &
     4619                    / ( u_comp - gu + SIGN( 1.0E-20_wp, u_comp - gu )    )     &
     4620                  + diss_r                                                     &
     4621                    *   ABS( u_comp - 2.0_wp * hom(k,1,1,0)              )     &
     4622                    / ( ABS( u_comp - gu ) + 1.0E-20_wp                  )     &
    42004623                  ) *   weight_substep(intermediate_timestep_count)
    42014624!
    42024625!--             Statistical Evaluation of w'u'.
     4626                !$ACC ATOMIC
    42034627                sums_wsus_ws_l(k,tn) = sums_wsus_ws_l(k,tn)                    &
    4204                 + ( flux_t(k)                                                  &
     4628                + ( flux_t                                                     &
    42054629                    * ( w_comp - 2.0_wp * hom(k,1,3,0)                   )     &
    42064630                    / ( w_comp + SIGN( 1.0E-20_wp, w_comp )              )     &
    4207                   + diss_t(k)                                                  &
     4631                  + diss_t                                                     &
    42084632                    *   ABS( w_comp - 2.0_wp * hom(k,1,3,0)              )     &
    42094633                    / ( ABS( w_comp ) + 1.0E-20_wp                       )     &
     
    42144638             DO  k = nzb_max+1, nzt
    42154639
    4216                 u_comp(k) = u(k,j,i+1) + u(k,j,i)
    4217                 flux_r(k) = ( u_comp(k) - gu ) * (                            &
     4640                u_comp = u(k,j,i+1) + u(k,j,i)
     4641                flux_r = ( u_comp - gu ) * (                                  &
    42184642                         37.0_wp * ( u(k,j,i+1) + u(k,j,i)   )                   &
    42194643                       -  8.0_wp * ( u(k,j,i+2) + u(k,j,i-1) )                   &
    42204644                       +           ( u(k,j,i+3) + u(k,j,i-2) ) ) * adv_mom_5
    4221                 diss_r(k) = - ABS( u_comp(k) - gu ) * (                       &
     4645                diss_r = - ABS( u_comp - gu ) * (                             &
    42224646                         10.0_wp * ( u(k,j,i+1) - u(k,j,i)   )                   &
    42234647                       -  5.0_wp * ( u(k,j,i+2) - u(k,j,i-1) )                   &
    42244648                       +           ( u(k,j,i+3) - u(k,j,i-2) ) ) * adv_mom_5
    42254649
    4226                 v_comp    = v(k,j+1,i) + v(k,j+1,i-1) - gv
    4227                 flux_n(k) = v_comp * (                                        &
     4650#ifdef _OPENACC
     4651!
     4652!--             Recompute the left fluxes.
     4653                u_comp_l = u(k,j,i) + u(k,j,i-1) - gu
     4654                flux_l   = u_comp_l * (                                       &
     4655                             37.0_wp * ( u(k,j,i) + u(k,j,i-1)   )                &
     4656                           -  8.0_wp * ( u(k,j,i+1) + u(k,j,i-2) )                &
     4657                           +           ( u(k,j,i+2) + u(k,j,i-3) ) ) * adv_mom_5
     4658                diss_l   = - ABS(u_comp_l) * (                                &
     4659                             10.0_wp * ( u(k,j,i) - u(k,j,i-1)   )                &
     4660                           -  5.0_wp * ( u(k,j,i+1) - u(k,j,i-2) )                &
     4661                           +           ( u(k,j,i+2) - u(k,j,i-3) ) ) * adv_mom_5
     4662#else
     4663                flux_l = swap_flux_x_local_u(k,j)
     4664                diss_l = swap_diss_x_local_u(k,j)
     4665#endif
     4666
     4667                v_comp = v(k,j+1,i) + v(k,j+1,i-1) - gv
     4668                flux_n = v_comp * (                                           &
    42284669                         37.0_wp * ( u(k,j+1,i) + u(k,j,i)   )                   &
    42294670                       -  8.0_wp * ( u(k,j+2,i) + u(k,j-1,i) )                   &
    42304671                       +           ( u(k,j+3,i) + u(k,j-2,i) ) ) * adv_mom_5
    4231                 diss_n(k) = - ABS( v_comp ) * (                               &
     4672                diss_n = - ABS( v_comp ) * (                                  &
    42324673                         10.0_wp * ( u(k,j+1,i) - u(k,j,i)   )                   &
    42334674                       -  5.0_wp * ( u(k,j+2,i) - u(k,j-1,i) )                   &
    42344675                       +           ( u(k,j+3,i) - u(k,j-2,i) ) ) * adv_mom_5
     4676
     4677#ifdef _OPENACC
     4678!
     4679!--             Recompute the south fluxes.
     4680                v_comp_s = v(k,j,i) + v(k,j,i-1) - gv
     4681                flux_s   = v_comp_s * (                                       &
     4682                           37.0_wp * ( u(k,j,i) + u(k,j-1,i)   )                 &
     4683                         -  8.0_wp * ( u(k,j+1,i) + u(k,j-2,i) )                 &
     4684                         +           ( u(k,j+2,i) + u(k,j-3,i) ) ) * adv_mom_5
     4685                diss_s   = - ABS( v_comp_s ) * (                              &
     4686                           10.0_wp * ( u(k,j,i) - u(k,j-1,i)   )                 &
     4687                         -  5.0_wp * ( u(k,j+1,i) - u(k,j-2,i) )                 &
     4688                         +           ( u(k,j+2,i) - u(k,j-3,i) ) ) * adv_mom_5
     4689#else
     4690                flux_s = swap_flux_y_local_u(k)
     4691                diss_s = swap_diss_y_local_u(k)
     4692#endif
     4693
    42354694!
    42364695!--             k index has to be modified near bottom and top, else array
     
    42444703                k_mm  = k - 2 * ibit17
    42454704
    4246                 w_comp    = w(k,j,i) + w(k,j,i-1)
    4247                 flux_t(k) = w_comp * rho_air_zw(k) * (                       &
     4705                w_comp = w(k,j,i) + w(k,j,i-1)
     4706                flux_t = w_comp * rho_air_zw(k) * (                          &
    42484707                          ( 37.0_wp * ibit17 * adv_mom_5                        &
    42494708                       +     7.0_wp * ibit16 * adv_mom_3                        &
     
    42604719                                      )
    42614720
    4262                 diss_t(k) = - ABS( w_comp ) * rho_air_zw(k) * (              &
     4721                diss_t = - ABS( w_comp ) * rho_air_zw(k) * (                 &
    42634722                          ( 10.0_wp * ibit17 * adv_mom_5                        &
    42644723                       +     3.0_wp * ibit16 * adv_mom_3                        &
     
    42784737!--             correction is needed to overcome numerical instabilities caused
    42794738!--             by a not sufficient reduction of divergences near topography.
    4280                 div = ( ( u_comp(k)   - ( u(k,j,i)   + u(k,j,i-1)   ) ) * ddx &
     4739                div = ( ( u_comp      - ( u(k,j,i)   + u(k,j,i-1)   ) ) * ddx &
    42814740                     +  ( v_comp + gv - ( v(k,j,i)   + v(k,j,i-1 )  ) ) * ddy &
    42824741                     +  (   w_comp                      * rho_air_zw(k) -     &
     
    42864745
    42874746                tend(k,j,i) = tend(k,j,i) - (                                  &
    4288                  ( flux_r(k) + diss_r(k)                                       &
    4289                -   swap_flux_x_local_u(k,j) - swap_diss_x_local_u(k,j) ) * ddx &
    4290                + ( flux_n(k) + diss_n(k)                                       &
    4291                -   swap_flux_y_local_u(k)   - swap_diss_y_local_u(k)   ) * ddy &
    4292                + ( ( flux_t(k) + diss_t(k) )                                   &
    4293                -   ( flux_d    + diss_d    )                                   &
    4294                                                     ) * drho_air(k) * ddzw(k)  &
     4747                 ( ( flux_r + diss_r )                                         &
     4748               -   ( flux_l + diss_l ) ) * ddx                                 &
     4749               + ( ( flux_n + diss_n )                                         &
     4750               -   ( flux_s + diss_s ) ) * ddy                                 &
     4751               + ( ( flux_t + diss_t )                                         &
     4752               -   ( flux_d + diss_d ) ) * drho_air(k) * ddzw(k)               &
    42954753                                           ) + div * u(k,j,i)
    42964754
    4297                 swap_flux_x_local_u(k,j) = flux_r(k)
    4298                 swap_diss_x_local_u(k,j) = diss_r(k)
    4299                 swap_flux_y_local_u(k)   = flux_n(k)
    4300                 swap_diss_y_local_u(k)   = diss_n(k)
    4301                 flux_d                   = flux_t(k)
    4302                 diss_d                   = diss_t(k)
     4755#ifndef _OPENACC
     4756                swap_flux_x_local_u(k,j) = flux_r
     4757                swap_diss_x_local_u(k,j) = diss_r
     4758                swap_flux_y_local_u(k)   = flux_n
     4759                swap_diss_y_local_u(k)   = diss_n
     4760#endif
     4761                flux_d                   = flux_t
     4762                diss_d                   = diss_t
    43034763!
    43044764!--             Statistical Evaluation of u'u'. The factor has to be applied
    43054765!--             for right evaluation when gallilei_trans = .T. .
     4766                !$ACC ATOMIC
    43064767                sums_us2_ws_l(k,tn) = sums_us2_ws_l(k,tn)                      &
    4307                 + ( flux_r(k)                                                  &
    4308                     * ( u_comp(k) - 2.0_wp * hom(k,1,1,0)                   )  &
    4309                     / ( u_comp(k) - gu + SIGN( 1.0E-20_wp, u_comp(k) - gu ) )  &
    4310                   + diss_r(k)                                                  &
    4311                     *   ABS( u_comp(k) - 2.0_wp * hom(k,1,1,0)              )  &
    4312                     / ( ABS( u_comp(k) - gu ) + 1.0E-20_wp                  )  &
     4768                + ( flux_r                                                     &
     4769                    * ( u_comp - 2.0_wp * hom(k,1,1,0)                   )     &
     4770                    / ( u_comp - gu + SIGN( 1.0E-20_wp, u_comp - gu )    )     &
     4771                  + diss_r                                                     &
     4772                    *   ABS( u_comp - 2.0_wp * hom(k,1,1,0)              )     &
     4773                    / ( ABS( u_comp - gu ) + 1.0E-20_wp                  )     &
    43134774                  ) *   weight_substep(intermediate_timestep_count)
    43144775!
    43154776!--             Statistical Evaluation of w'u'.
     4777                !$ACC ATOMIC
    43164778                sums_wsus_ws_l(k,tn) = sums_wsus_ws_l(k,tn)                    &
    4317                 + ( flux_t(k)                                                  &
     4779                + ( flux_t                                                     &
    43184780                    * ( w_comp - 2.0_wp * hom(k,1,3,0)                   )     &
    43194781                    / ( w_comp + SIGN( 1.0E-20_wp, w_comp )              )     &
    4320                   + diss_t(k)                                                  &
     4782                  + diss_t                                                     &
    43214783                    *   ABS( w_comp - 2.0_wp * hom(k,1,3,0)              )     &
    43224784                    / ( ABS( w_comp ) + 1.0E-20_wp                       )     &
     
    43694831       REAL(wp)    ::  ibit19 !< flag indicating 3rd-order scheme along x-direction
    43704832       REAL(wp)    ::  ibit20 !< flag indicating 5th-order scheme along x-direction
     4833#ifdef _OPENACC
     4834       REAL(wp)    ::  ibit18_l !< flag indicating 1st-order scheme along x-direction
     4835       REAL(wp)    ::  ibit19_l !< flag indicating 3rd-order scheme along x-direction
     4836       REAL(wp)    ::  ibit20_l !< flag indicating 5th-order scheme along x-direction
     4837#endif
    43714838       REAL(wp)    ::  ibit21 !< flag indicating 1st-order scheme along y-direction
    43724839       REAL(wp)    ::  ibit22 !< flag indicating 3rd-order scheme along y-direction
    43734840       REAL(wp)    ::  ibit23 !< flag indicating 5th-order scheme along y-direction
     4841#ifdef _OPENACC
     4842       REAL(wp)    ::  ibit21_s !< flag indicating 1st-order scheme along y-direction
     4843       REAL(wp)    ::  ibit22_s !< flag indicating 3rd-order scheme along y-direction
     4844       REAL(wp)    ::  ibit23_s !< flag indicating 5th-order scheme along y-direction
     4845#endif
    43744846       REAL(wp)    ::  ibit24 !< flag indicating 1st-order scheme along z-direction
    43754847       REAL(wp)    ::  ibit25 !< flag indicating 3rd-order scheme along z-direction
     
    43814853       REAL(wp)    ::  gv     !< Galilei-transformation velocity along y
    43824854       REAL(wp)    ::  u_comp !< advection velocity along x
     4855#ifdef _OPENACC
     4856       REAL(wp)    ::  u_comp_l !< advection velocity along x
     4857#endif
    43834858       REAL(wp)    ::  w_comp !< advection velocity along z
    43844859       
     4860       REAL(wp)    ::  diss_s !< discretized artificial dissipation at southward-side of the grid box
     4861       REAL(wp)    ::  flux_s !< discretized 6th-order flux at southward-side of the grid box
     4862#ifndef _OPENACC
    43854863       REAL(wp), DIMENSION(nzb+1:nzt) ::  swap_diss_y_local_v !< discretized artificial dissipation at southward-side of the grid box
    43864864       REAL(wp), DIMENSION(nzb+1:nzt) ::  swap_flux_y_local_v !< discretized 6th-order flux at southward-side of the grid box
     4865#endif
    43874866       
     4867       REAL(wp)    ::  diss_l !< discretized artificial dissipation at leftward-side of the grid box
     4868       REAL(wp)    ::  flux_l !< discretized 6th-order flux at leftward-side of the grid box
     4869#ifndef _OPENACC
    43884870       REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  swap_diss_x_local_v !< discretized artificial dissipation at leftward-side of the grid box
    43894871       REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  swap_flux_x_local_v !< discretized 6th-order flux at leftward-side of the grid box
     4872#endif
    43904873       
    4391        REAL(wp), DIMENSION(nzb:nzt) ::  diss_n !< discretized artificial dissipation at northward-side of the grid box
    4392        REAL(wp), DIMENSION(nzb:nzt) ::  diss_r !< discretized artificial dissipation at rightward-side of the grid box
    4393        REAL(wp), DIMENSION(nzb:nzt) ::  diss_t !< discretized artificial dissipation at top of the grid box
    4394        REAL(wp), DIMENSION(nzb:nzt) ::  flux_n !< discretized 6th-order flux at northward-side of the grid box
    4395        REAL(wp), DIMENSION(nzb:nzt) ::  flux_r !< discretized 6th-order flux at rightward-side of the grid box
    4396        REAL(wp), DIMENSION(nzb:nzt) ::  flux_t !< discretized 6th-order flux at top of the grid box
    4397        REAL(wp), DIMENSION(nzb:nzt) ::  v_comp !< advection velocity along y
     4874       REAL(wp)    ::  diss_n !< discretized artificial dissipation at northward-side of the grid box
     4875       REAL(wp)    ::  diss_r !< discretized artificial dissipation at rightward-side of the grid box
     4876       REAL(wp)    ::  diss_t !< discretized artificial dissipation at top of the grid box
     4877       REAL(wp)    ::  flux_n !< discretized 6th-order flux at northward-side of the grid box
     4878       REAL(wp)    ::  flux_r !< discretized 6th-order flux at rightward-side of the grid box
     4879       REAL(wp)    ::  flux_t !< discretized 6th-order flux at top of the grid box
     4880       REAL(wp)    ::  v_comp !< advection velocity along y
     4881#ifdef _OPENACC
     4882       REAL(wp)    ::  v_comp_s !<
     4883#endif
    43984884
    43994885       gu = 2.0_wp * u_gtrans
    44004886       gv = 2.0_wp * v_gtrans
     4887
     4888#ifndef _OPENACC
    44014889!
    44024890!--    First compute the whole left boundary of the processor domain
     
    44574945
    44584946       ENDDO
    4459 
     4947#endif
     4948
     4949       !$ACC PARALLEL LOOP COLLAPSE(2) FIRSTPRIVATE(tn, gu, gv) &
     4950       !$ACC PRIVATE(i, j, k, k_mm, k_pp, k_ppp) &
     4951       !$ACC PRIVATE(ibit18, ibit19, ibit20, ibit21, ibit22, ibit23) &
     4952       !$ACC PRIVATE(ibit18_l, ibit19_l, ibit20_l) &
     4953       !$ACC PRIVATE(ibit21_s, ibit22_s, ibit23_s) &
     4954       !$ACC PRIVATE(ibit24, ibit25, ibit26) &
     4955       !$ACC PRIVATE(flux_r, diss_r, flux_l, diss_l) &
     4956       !$ACC PRIVATE(flux_n, diss_n, flux_s, diss_s) &
     4957       !$ACC PRIVATE(flux_t, diss_t, flux_d, diss_d) &
     4958       !$ACC PRIVATE(div, u_comp, u_comp_l, v_comp, v_comp_s, w_comp) &
     4959       !$ACC PRESENT(advc_flags_1) &
     4960       !$ACC PRESENT(u, v, w) &
     4961       !$ACC PRESENT(drho_air, rho_air_zw, ddzw) &
     4962       !$ACC PRESENT(tend) &
     4963       !$ACC PRESENT(hom(nzb+1:nzb_max,1,2:3,0)) &
     4964       !$ACC PRESENT(weight_substep(intermediate_timestep_count)) &
     4965       !$ACC PRESENT(sums_vs2_ws_l, sums_wsvs_ws_l)
    44604966       DO i = nxl, nxr
    44614967
     4968#ifndef _OPENACC
    44624969          j = nysv
    44634970          DO  k = nzb+1, nzb_max
     
    44674974             ibit21 = REAL( IBITS(advc_flags_1(k,j-1,i),21,1), KIND = wp )
    44684975
    4469              v_comp(k)              = v(k,j,i) + v(k,j-1,i) - gv
    4470              swap_flux_y_local_v(k) = v_comp(k) * (                           &
     4976             v_comp                 = v(k,j,i) + v(k,j-1,i) - gv
     4977             swap_flux_y_local_v(k) = v_comp * (                              &
    44714978                                   ( 37.0_wp * ibit23 * adv_mom_5                &
    44724979                                +     7.0_wp * ibit22 * adv_mom_3                &
     
    44834990                                                 )
    44844991
    4485              swap_diss_y_local_v(k) = - ABS( v_comp(k) ) * (                  &
     4992             swap_diss_y_local_v(k) = - ABS( v_comp ) * (                     &
    44864993                                   ( 10.0_wp * ibit23 * adv_mom_5                &
    44874994                                +     3.0_wp * ibit22 * adv_mom_3                &
     
    45025009          DO  k = nzb_max+1, nzt
    45035010
    4504              v_comp(k)              = v(k,j,i) + v(k,j-1,i) - gv
    4505              swap_flux_y_local_v(k) = v_comp(k) * (                           &
     5011             v_comp                 = v(k,j,i) + v(k,j-1,i) - gv
     5012             swap_flux_y_local_v(k) = v_comp * (                              &
    45065013                           37.0_wp * ( v(k,j,i) + v(k,j-1,i)   )                 &
    45075014                         -  8.0_wp * ( v(k,j+1,i) + v(k,j-2,i) )                 &
    45085015                         +           ( v(k,j+2,i) + v(k,j-3,i) ) ) * adv_mom_5
    4509              swap_diss_y_local_v(k) = - ABS( v_comp(k) ) * (                  &
     5016             swap_diss_y_local_v(k) = - ABS( v_comp ) * (                     &
    45105017                           10.0_wp * ( v(k,j,i) - v(k,j-1,i)   )                 &
    45115018                         -  5.0_wp * ( v(k,j+1,i) - v(k,j-2,i) )                 &
     
    45135020
    45145021          ENDDO
     5022#endif
    45155023
    45165024          DO  j = nysv, nyn
    45175025
    4518              flux_t(0) = 0.0_wp
    4519              diss_t(0) = 0.0_wp
    45205026             flux_d    = 0.0_wp
    45215027             diss_d    = 0.0_wp
     
    45275033                ibit18 = REAL( IBITS(advc_flags_1(k,j,i),18,1), KIND = wp )
    45285034
    4529                 u_comp    = u(k,j-1,i+1) + u(k,j,i+1) - gu
    4530                 flux_r(k) = u_comp * (                                       &
     5035                u_comp = u(k,j-1,i+1) + u(k,j,i+1) - gu
     5036                flux_r = u_comp * (                                          &
    45315037                          ( 37.0_wp * ibit20 * adv_mom_5                        &
    45325038                       +     7.0_wp * ibit19 * adv_mom_3                        &
     
    45435049                                     )
    45445050
    4545                 diss_r(k) = - ABS( u_comp ) * (                              &
     5051                diss_r = - ABS( u_comp ) * (                                 &
    45465052                          ( 10.0_wp * ibit20 * adv_mom_5                        &
    45475053                       +     3.0_wp * ibit19 * adv_mom_3                        &
     
    45585064                                              )
    45595065
     5066#ifdef _OPENACC
     5067!
     5068!--             Recompute the left fluxes.
     5069                ibit20_l = REAL( IBITS(advc_flags_1(k,j,i-1),20,1), KIND = wp )
     5070                ibit19_l = REAL( IBITS(advc_flags_1(k,j,i-1),19,1), KIND = wp )
     5071                ibit18_l = REAL( IBITS(advc_flags_1(k,j,i-1),18,1), KIND = wp )
     5072
     5073                u_comp_l  = u(k,j-1,i) + u(k,j,i) - gu
     5074                flux_l    = u_comp_l * (                                     &
     5075                                      ( 37.0_wp * ibit20_l * adv_mom_5          &
     5076                                   +     7.0_wp * ibit19_l * adv_mom_3          &
     5077                                   +              ibit18_l * adv_mom_1          &
     5078                                      ) *                                    &
     5079                                     ( v(k,j,i)   + v(k,j,i-1) )             &
     5080                               -      (  8.0_wp * ibit20_l * adv_mom_5          &
     5081                                   +              ibit19_l * adv_mom_3          &
     5082                                      ) *                                    &
     5083                                     ( v(k,j,i+1) + v(k,j,i-2) )             &
     5084                               +      (           ibit20_l * adv_mom_5          &
     5085                                      ) *                                    &
     5086                                     ( v(k,j,i+2) + v(k,j,i-3) )             &
     5087                                                 )
     5088
     5089                 diss_l   = - ABS( u_comp_l ) * (                            &
     5090                                      ( 10.0_wp * ibit20_l * adv_mom_5          &
     5091                                   +     3.0_wp * ibit19_l * adv_mom_3          &
     5092                                   +              ibit18_l * adv_mom_1          &
     5093                                      ) *                                    &
     5094                                     ( v(k,j,i)   - v(k,j,i-1) )             &
     5095                               -      (  5.0_wp * ibit20_l * adv_mom_5          &
     5096                                   +              ibit19_l * adv_mom_3          &
     5097                                      ) *                                    &
     5098                                     ( v(k,j,i+1) - v(k,j,i-2) )             &
     5099                               +      (           ibit20_l * adv_mom_5          &
     5100                                      ) *                                    &
     5101                                     ( v(k,j,i+2) - v(k,j,i-3) )             &
     5102                                                           )
     5103#else
     5104                flux_l = swap_flux_x_local_v(k,j)
     5105                diss_l = swap_diss_x_local_v(k,j)
     5106#endif
     5107
    45605108                ibit23 = REAL( IBITS(advc_flags_1(k,j,i),23,1), KIND = wp )
    45615109                ibit22 = REAL( IBITS(advc_flags_1(k,j,i),22,1), KIND = wp )
    45625110                ibit21 = REAL( IBITS(advc_flags_1(k,j,i),21,1), KIND = wp )
    45635111
    4564                 v_comp(k) = v(k,j+1,i) + v(k,j,i)
    4565                 flux_n(k) = ( v_comp(k) - gv ) * (                           &
     5112                v_comp = v(k,j+1,i) + v(k,j,i)
     5113                flux_n = ( v_comp - gv ) * (                                 &
    45665114                          ( 37.0_wp * ibit23 * adv_mom_5                        &
    45675115                       +     7.0_wp * ibit22 * adv_mom_3                        &
     
    45785126                                     )
    45795127
    4580                 diss_n(k) = - ABS( v_comp(k) - gv ) * (                      &
     5128                diss_n = - ABS( v_comp - gv ) * (                            &
    45815129                          ( 10.0_wp * ibit23 * adv_mom_5                        &
    45825130                       +     3.0_wp * ibit22 * adv_mom_3                        &
     
    45925140                                 ( v(k,j+3,i) - v(k,j-2,i) )                 &
    45935141                                                      )
     5142
     5143#ifdef _OPENACC
     5144!
     5145!--             Recompute the south fluxes.
     5146                ibit23_s = REAL( IBITS(advc_flags_1(k,j-1,i),23,1), KIND = wp )
     5147                ibit22_s = REAL( IBITS(advc_flags_1(k,j-1,i),22,1), KIND = wp )
     5148                ibit21_s = REAL( IBITS(advc_flags_1(k,j-1,i),21,1), KIND = wp )
     5149
     5150                v_comp_s = v(k,j,i) + v(k,j-1,i) - gv
     5151                flux_s   = v_comp_s * (                                      &
     5152                                   ( 37.0_wp * ibit23_s * adv_mom_5             &
     5153                                +     7.0_wp * ibit22_s * adv_mom_3             &
     5154                                +              ibit21_s * adv_mom_1             &
     5155                                   ) *                                       &
     5156                                     ( v(k,j,i)   + v(k,j-1,i) )             &
     5157                            -      (  8.0_wp * ibit23_s * adv_mom_5             &
     5158                                +              ibit22_s * adv_mom_3             &
     5159                                   ) *                                       &
     5160                                     ( v(k,j+1,i) + v(k,j-2,i) )             &
     5161                            +      (           ibit23_s * adv_mom_5             &
     5162                                   ) *                                       &
     5163                                     ( v(k,j+2,i) + v(k,j-3,i) )             &
     5164                                                 )
     5165
     5166                diss_s   = - ABS( v_comp_s ) * (                             &
     5167                                   ( 10.0_wp * ibit23_s * adv_mom_5             &
     5168                                +     3.0_wp * ibit22_s * adv_mom_3             &
     5169                                +              ibit21_s * adv_mom_1             &
     5170                                   ) *                                       &
     5171                                     ( v(k,j,i)   - v(k,j-1,i) )             &
     5172                            -      (  5.0_wp * ibit23_s * adv_mom_5             &
     5173                                +              ibit22_s * adv_mom_3             &
     5174                                   ) *                                       &
     5175                                     ( v(k,j+1,i) - v(k,j-2,i) )             &
     5176                            +      (           ibit23_s * adv_mom_5             &
     5177                                   ) *                                       &
     5178                                     ( v(k,j+2,i) - v(k,j-3,i) )             &
     5179                                                          )
     5180#else
     5181               flux_s = swap_flux_y_local_v(k)
     5182               diss_s = swap_diss_y_local_v(k)
     5183#endif
     5184
    45945185!
    45955186!--             k index has to be modified near bottom and top, else array
     
    46035194                k_mm  = k - 2 * ibit26
    46045195
    4605                 w_comp    = w(k,j-1,i) + w(k,j,i)
    4606                 flux_t(k) = w_comp * rho_air_zw(k) * (                       &
     5196                w_comp = w(k,j-1,i) + w(k,j,i)
     5197                flux_t = w_comp * rho_air_zw(k) * (                          &
    46075198                          ( 37.0_wp * ibit26 * adv_mom_5                        &
    46085199                       +     7.0_wp * ibit25 * adv_mom_3                        &
     
    46195210                                      )
    46205211
    4621                 diss_t(k) = - ABS( w_comp ) * rho_air_zw(k) * (              &
     5212                diss_t = - ABS( w_comp ) * rho_air_zw(k) * (                 &
    46225213                          ( 10.0_wp * ibit26 * adv_mom_5                        &
    46235214                       +     3.0_wp * ibit25 * adv_mom_3                        &
     
    46465237                                         )                                    &
    46475238                  ) * ddx                                                     &
    4648                +  ( v_comp(k)                                                 &
     5239               +  ( v_comp                                                    &
    46495240                                       * ( ibit21 + ibit22 + ibit23 )         &
    46505241                - ( v(k,j,i)     + v(k,j-1,i) )                               &
     
    46685259
    46695260                tend(k,j,i) = tend(k,j,i) - (                                 &
    4670                        ( flux_r(k) + diss_r(k)                                &
    4671                      -   swap_flux_x_local_v(k,j) - swap_diss_x_local_v(k,j)  &
    4672                        ) * ddx                                                &
    4673                      + ( flux_n(k) + diss_n(k)                                &
    4674                      -   swap_flux_y_local_v(k) - swap_diss_y_local_v(k)      &
    4675                        ) * ddy                                                &
    4676                      + ( ( flux_t(k) + diss_t(k) )                            &
    4677                      -   ( flux_d    + diss_d    )                            &
    4678                        ) * drho_air(k) * ddzw(k)                              &
     5261                       ( ( flux_r + diss_r )                                  &
     5262                     -   ( flux_l + diss_l ) ) * ddx                          &
     5263                     + ( ( flux_n + diss_n )                                  &
     5264                     -   ( flux_s + diss_s ) ) * ddy                          &
     5265                     + ( ( flux_t + diss_t )                                  &
     5266                     -   ( flux_d + diss_d ) ) * drho_air(k) * ddzw(k)        &
    46795267                                            )  + v(k,j,i) * div
    46805268
    4681                 swap_flux_x_local_v(k,j) = flux_r(k)
    4682                 swap_diss_x_local_v(k,j) = diss_r(k)
    4683                 swap_flux_y_local_v(k)   = flux_n(k)
    4684                 swap_diss_y_local_v(k)   = diss_n(k)
    4685                 flux_d                   = flux_t(k)
    4686                 diss_d                   = diss_t(k)
     5269#ifndef _OPENACC
     5270                swap_flux_x_local_v(k,j) = flux_r
     5271                swap_diss_x_local_v(k,j) = diss_r
     5272                swap_flux_y_local_v(k)   = flux_n
     5273                swap_diss_y_local_v(k)   = diss_n
     5274#endif
     5275                flux_d                   = flux_t
     5276                diss_d                   = diss_t
    46875277
    46885278!
    46895279!--             Statistical Evaluation of v'v'. The factor has to be applied
    46905280!--             for right evaluation when gallilei_trans = .T. .
     5281                !$ACC ATOMIC
    46915282                sums_vs2_ws_l(k,tn) = sums_vs2_ws_l(k,tn)                      &
    4692                 + ( flux_n(k)                                                  &
    4693                     * ( v_comp(k) - 2.0_wp * hom(k,1,2,0)                   )  &
    4694                     / ( v_comp(k) - gv + SIGN( 1.0E-20_wp, v_comp(k) - gv ) )  &
    4695                +   diss_n(k)                                                   &
    4696                     *   ABS( v_comp(k) - 2.0_wp * hom(k,1,2,0)              )  &
    4697                     / ( ABS( v_comp(k) - gv ) + 1.0E-20_wp                  )  &
     5283                + ( flux_n                                                     &
     5284                    * ( v_comp - 2.0_wp * hom(k,1,2,0)                   )     &
     5285                    / ( v_comp - gv + SIGN( 1.0E-20_wp, v_comp - gv )    )     &
     5286               +   diss_n                                                      &
     5287                    *   ABS( v_comp - 2.0_wp * hom(k,1,2,0)              )     &
     5288                    / ( ABS( v_comp - gv ) + 1.0E-20_wp                  )     &
    46985289                  ) *   weight_substep(intermediate_timestep_count)
    46995290!
    47005291!--             Statistical Evaluation of w'u'.
     5292                !$ACC ATOMIC
    47015293                sums_wsvs_ws_l(k,tn) = sums_wsvs_ws_l(k,tn)                    &
    4702                 + ( flux_t(k)                                                  &
     5294                + ( flux_t                                                     &
    47035295                    * ( w_comp - 2.0_wp * hom(k,1,3,0)                   )     &
    47045296                    / ( w_comp + SIGN( 1.0E-20_wp, w_comp )              )     &
    4705                +   diss_t(k)                                                   &
     5297               +   diss_t                                                      &
    47065298                    *   ABS( w_comp - 2.0_wp * hom(k,1,3,0)              )     &
    47075299                    / ( ABS( w_comp ) + 1.0E-20_wp                       )     &
     
    47125304             DO  k = nzb_max+1, nzt
    47135305
    4714                 u_comp    = u(k,j-1,i+1) + u(k,j,i+1) - gu
    4715                 flux_r(k) = u_comp * (                                        &
     5306                u_comp = u(k,j-1,i+1) + u(k,j,i+1) - gu
     5307                flux_r = u_comp * (                                           &
    47165308                      37.0_wp * ( v(k,j,i+1) + v(k,j,i)   )                      &
    47175309                    -  8.0_wp * ( v(k,j,i+2) + v(k,j,i-1) )                      &
    47185310                    +           ( v(k,j,i+3) + v(k,j,i-2) ) ) * adv_mom_5
    47195311
    4720                 diss_r(k) = - ABS( u_comp ) * (                               &
     5312                diss_r = - ABS( u_comp ) * (                                  &
    47215313                      10.0_wp * ( v(k,j,i+1) - v(k,j,i) )                        &
    47225314                    -  5.0_wp * ( v(k,j,i+2) - v(k,j,i-1) )                      &
    47235315                    +           ( v(k,j,i+3) - v(k,j,i-2) ) ) * adv_mom_5
    47245316
    4725 
    4726                 v_comp(k) = v(k,j+1,i) + v(k,j,i)
    4727                 flux_n(k) = ( v_comp(k) - gv ) * (                            &
     5317#ifdef _OPENACC
     5318!
     5319!--             Recompute the left fluxes.
     5320                u_comp_l = u(k,j-1,i) + u(k,j,i) - gu
     5321                flux_l   = u_comp_l * (                                       &
     5322                             37.0_wp * ( v(k,j,i) + v(k,j,i-1)   )               &
     5323                           -  8.0_wp * ( v(k,j,i+1) + v(k,j,i-2) )               &
     5324                           +           ( v(k,j,i+2) + v(k,j,i-3) ) ) * adv_mom_5
     5325                diss_l   = - ABS( u_comp_l ) * (                              &
     5326                             10.0_wp * ( v(k,j,i) - v(k,j,i-1)   )               &
     5327                           -  5.0_wp * ( v(k,j,i+1) - v(k,j,i-2) )               &
     5328                           +           ( v(k,j,i+2) - v(k,j,i-3) ) ) * adv_mom_5
     5329#else
     5330                flux_l = swap_flux_x_local_v(k,j)
     5331                diss_l = swap_diss_x_local_v(k,j)
     5332#endif
     5333
     5334                v_comp = v(k,j+1,i) + v(k,j,i)
     5335                flux_n = ( v_comp - gv ) * (                                  &
    47285336                      37.0_wp * ( v(k,j+1,i) + v(k,j,i)   )                      &
    47295337                    -  8.0_wp * ( v(k,j+2,i) + v(k,j-1,i) )                      &
    47305338                      +         ( v(k,j+3,i) + v(k,j-2,i) ) ) * adv_mom_5
    47315339
    4732                 diss_n(k) = - ABS( v_comp(k) - gv ) * (                       &
     5340                diss_n = - ABS( v_comp - gv ) * (                             &
    47335341                      10.0_wp * ( v(k,j+1,i) - v(k,j,i)   )                      &
    47345342                    -  5.0_wp * ( v(k,j+2,i) - v(k,j-1,i) )                      &
    47355343                    +           ( v(k,j+3,i) - v(k,j-2,i) ) ) * adv_mom_5
     5344
     5345#ifdef _OPENACC
     5346!
     5347!--             Recompute the south fluxes.
     5348                v_comp_s = v(k,j,i) + v(k,j-1,i) - gv
     5349                flux_s   = v_comp_s * (                                       &
     5350                           37.0_wp * ( v(k,j,i) + v(k,j-1,i)   )                 &
     5351                         -  8.0_wp * ( v(k,j+1,i) + v(k,j-2,i) )                 &
     5352                         +           ( v(k,j+2,i) + v(k,j-3,i) ) ) * adv_mom_5
     5353                diss_s   = - ABS( v_comp_s ) * (                              &
     5354                           10.0_wp * ( v(k,j,i) - v(k,j-1,i)   )                 &
     5355                         -  5.0_wp * ( v(k,j+1,i) - v(k,j-2,i) )                 &
     5356                         +           ( v(k,j+2,i) - v(k,j-3,i) ) ) * adv_mom_5
     5357#else
     5358                flux_s = swap_flux_y_local_v(k)
     5359                diss_s = swap_diss_y_local_v(k)
     5360#endif
     5361
    47365362!
    47375363!--             k index has to be modified near bottom and top, else array
     
    47455371                k_mm  = k - 2 * ibit26
    47465372
    4747                 w_comp    = w(k,j-1,i) + w(k,j,i)
    4748                 flux_t(k) = w_comp * rho_air_zw(k) * (                       &
     5373                w_comp = w(k,j-1,i) + w(k,j,i)
     5374                flux_t = w_comp * rho_air_zw(k) * (                          &
    47495375                          ( 37.0_wp * ibit26 * adv_mom_5                        &
    47505376                       +     7.0_wp * ibit25 * adv_mom_3                        &
     
    47615387                                      )
    47625388
    4763                 diss_t(k) = - ABS( w_comp ) * rho_air_zw(k) * (              &
     5389                diss_t = - ABS( w_comp ) * rho_air_zw(k) * (                 &
    47645390                          ( 10.0_wp * ibit26 * adv_mom_5                        &
    47655391                       +     3.0_wp * ibit25 * adv_mom_3                        &
     
    47805406!--             by a not sufficient reduction of divergences near topography.
    47815407                div = ( ( u_comp + gu - ( u(k,j-1,i)   + u(k,j,i)   ) ) * ddx &
    4782                      +  ( v_comp(k)   - ( v(k,j,i)     + v(k,j-1,i) ) ) * ddy &
     5408                     +  ( v_comp      - ( v(k,j,i)     + v(k,j-1,i) ) ) * ddy &
    47835409                     +  (   w_comp                      * rho_air_zw(k) -     &
    47845410                          ( w(k-1,j-1,i) + w(k-1,j,i) ) * rho_air_zw(k-1)     &
     
    47875413 
    47885414                tend(k,j,i) = tend(k,j,i) - (                                 &
    4789                        ( flux_r(k) + diss_r(k)                                &
    4790                      -   swap_flux_x_local_v(k,j) - swap_diss_x_local_v(k,j)  &
    4791                        ) * ddx                                                &
    4792                      + ( flux_n(k) + diss_n(k)                                &
    4793                      -   swap_flux_y_local_v(k) - swap_diss_y_local_v(k)      &
    4794                        ) * ddy                                                &
    4795                      + ( ( flux_t(k) + diss_t(k) )                            &
    4796                      -   ( flux_d    + diss_d    )                            &
    4797                        ) * drho_air(k) * ddzw(k)                              &
     5415                       ( ( flux_r + diss_r )                                  &
     5416                     -   ( flux_l + diss_l ) ) * ddx                          &
     5417                     + ( ( flux_n + diss_n )                                  &
     5418                     -   ( flux_s + diss_s ) ) * ddy                          &
     5419                     + ( ( flux_t + diss_t )                                  &
     5420                     -   ( flux_d + diss_d ) ) * drho_air(k) * ddzw(k)        &
    47985421                                            )  + v(k,j,i) * div
    47995422
    4800                 swap_flux_x_local_v(k,j) = flux_r(k)
    4801                 swap_diss_x_local_v(k,j) = diss_r(k)
    4802                 swap_flux_y_local_v(k)   = flux_n(k)
    4803                 swap_diss_y_local_v(k)   = diss_n(k)
    4804                 flux_d                   = flux_t(k)
    4805                 diss_d                   = diss_t(k)
     5423#ifndef _OPENACC
     5424                swap_flux_x_local_v(k,j) = flux_r
     5425                swap_diss_x_local_v(k,j) = diss_r
     5426                swap_flux_y_local_v(k)   = flux_n
     5427                swap_diss_y_local_v(k)   = diss_n
     5428#endif
     5429                flux_d                   = flux_t
     5430                diss_d                   = diss_t
    48065431
    48075432!
    48085433!--             Statistical Evaluation of v'v'. The factor has to be applied
    48095434!--             for right evaluation when gallilei_trans = .T. .
     5435                !$ACC ATOMIC
    48105436                sums_vs2_ws_l(k,tn) = sums_vs2_ws_l(k,tn)                      &
    4811                 + ( flux_n(k)                                                  &
    4812                     * ( v_comp(k) - 2.0_wp * hom(k,1,2,0)                   )  &
    4813                     / ( v_comp(k) - gv + SIGN( 1.0E-20_wp, v_comp(k) - gv ) )  &
    4814                +   diss_n(k)                                                   &
    4815                     *   ABS( v_comp(k) - 2.0_wp * hom(k,1,2,0)              )  &
    4816                     / ( ABS( v_comp(k) - gv ) + 1.0E-20_wp                  )  &
     5437                + ( flux_n                                                     &
     5438                    * ( v_comp - 2.0_wp * hom(k,1,2,0)                   )     &
     5439                    / ( v_comp - gv + SIGN( 1.0E-20_wp, v_comp - gv )    )     &
     5440               +   diss_n                                                      &
     5441                    *   ABS( v_comp - 2.0_wp * hom(k,1,2,0)              )     &
     5442                    / ( ABS( v_comp - gv ) + 1.0E-20_wp                  )     &
    48175443                  ) *   weight_substep(intermediate_timestep_count)
    48185444!
    48195445!--             Statistical Evaluation of w'u'.
     5446                !$ACC ATOMIC
    48205447                sums_wsvs_ws_l(k,tn) = sums_wsvs_ws_l(k,tn)                    &
    4821                 + ( flux_t(k)                                                  &
     5448                + ( flux_t                                                     &
    48225449                    * ( w_comp - 2.0_wp * hom(k,1,3,0)                   )     &
    48235450                    / ( w_comp + SIGN( 1.0E-20_wp, w_comp )              )     &
    4824                +   diss_t(k)                                                   &
     5451               +   diss_t                                                      &
    48255452                    *   ABS( w_comp - 2.0_wp * hom(k,1,3,0)              )     &
    48265453                    / ( ABS( w_comp ) + 1.0E-20_wp                       )     &
     
    48305457          ENDDO
    48315458       ENDDO
     5459!$ACC UPDATE HOST(sums_vs2_ws_l(nzb+1,tn))
    48325460       sums_vs2_ws_l(nzb,tn) = sums_vs2_ws_l(nzb+1,tn)
     5461!$ACC UPDATE DEVICE(sums_vs2_ws_l(nzb,tn))
    48335462
    48345463
     
    48745503       REAL(wp)    ::  ibit28 !< flag indicating 3rd-order scheme along x-direction
    48755504       REAL(wp)    ::  ibit29 !< flag indicating 5th-order scheme along x-direction
     5505#ifdef _OPENACC
     5506       REAL(wp)    ::  ibit27_l !< flag indicating 1st-order scheme along x-direction
     5507       REAL(wp)    ::  ibit28_l !< flag indicating 3rd-order scheme along x-direction
     5508       REAL(wp)    ::  ibit29_l !< flag indicating 5th-order scheme along x-direction
     5509#endif
    48765510       REAL(wp)    ::  ibit30 !< flag indicating 1st-order scheme along y-direction
    48775511       REAL(wp)    ::  ibit31 !< flag indicating 3rd-order scheme along y-direction
    48785512       REAL(wp)    ::  ibit32 !< flag indicating 5th-order scheme along y-direction
     5513#ifdef _OPENACC
     5514       REAL(wp)    ::  ibit30_s !< flag indicating 1st-order scheme along y-direction
     5515       REAL(wp)    ::  ibit31_s !< flag indicating 3rd-order scheme along y-direction
     5516       REAL(wp)    ::  ibit32_s !< flag indicating 5th-order scheme along y-direction
     5517#endif
    48795518       REAL(wp)    ::  ibit33 !< flag indicating 1st-order scheme along z-direction
    48805519       REAL(wp)    ::  ibit34 !< flag indicating 3rd-order scheme along z-direction
     
    48865525       REAL(wp)    ::  gv     !< Galilei-transformation velocity along y
    48875526       REAL(wp)    ::  u_comp !< advection velocity along x
     5527#ifdef _OPENACC
     5528       REAL(wp)    ::  u_comp_l !< advection velocity along x
     5529#endif
    48885530       REAL(wp)    ::  v_comp !< advection velocity along y
     5531#ifdef _OPENACC
     5532       REAL(wp)    ::  v_comp_s !< advection velocity along y
     5533#endif
    48895534       REAL(wp)    ::  w_comp !< advection velocity along z
    48905535       
    4891        REAL(wp), DIMENSION(nzb:nzt)    ::  diss_t !< discretized artificial dissipation at top of the grid box
    4892        REAL(wp), DIMENSION(nzb:nzt)    ::  flux_t !< discretized 6th-order flux at top of the grid box
     5536       REAL(wp)    ::  diss_t !< discretized artificial dissipation at top of the grid box
     5537       REAL(wp)    ::  flux_t !< discretized 6th-order flux at top of the grid box
    48935538       
    4894        REAL(wp), DIMENSION(nzb+1:nzt)  ::  diss_n !< discretized artificial dissipation at northward-side of the grid box
    4895        REAL(wp), DIMENSION(nzb+1:nzt)  ::  diss_r !< discretized artificial dissipation at rightward-side of the grid box
    4896        REAL(wp), DIMENSION(nzb+1:nzt)  ::  flux_n !< discretized 6th-order flux at northward-side of the grid box
    4897        REAL(wp), DIMENSION(nzb+1:nzt)  ::  flux_r !< discretized 6th-order flux at rightward-side of the grid box
     5539       REAL(wp)    ::  diss_n !< discretized artificial dissipation at northward-side of the grid box
     5540       REAL(wp)    ::  diss_r !< discretized artificial dissipation at rightward-side of the grid box
     5541       REAL(wp)    ::  flux_n !< discretized 6th-order flux at northward-side of the grid box
     5542       REAL(wp)    ::  flux_r !< discretized 6th-order flux at rightward-side of the grid box
     5543
     5544       REAL(wp)    ::  diss_s !< discretized artificial dissipation at southward-side of the grid box
     5545       REAL(wp)    ::  flux_s !< discretized 6th-order flux at southward-side of the grid box
     5546#ifndef _OPENACC
    48985547       REAL(wp), DIMENSION(nzb+1:nzt)  ::  swap_diss_y_local_w !< discretized artificial dissipation at southward-side of the grid box
    48995548       REAL(wp), DIMENSION(nzb+1:nzt)  ::  swap_flux_y_local_w !< discretized 6th-order flux at southward-side of the grid box
     5549#endif
    49005550       
     5551       REAL(wp)    ::  diss_l !< discretized artificial dissipation at leftward-side of the grid box
     5552       REAL(wp)    ::  flux_l !< discretized 6th-order flux at leftward-side of the grid box
     5553#ifndef _OPENACC
    49015554       REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  swap_diss_x_local_w !< discretized artificial dissipation at leftward-side of the grid box
    49025555       REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) ::  swap_flux_x_local_w !< discretized 6th-order flux at leftward-side of the grid box
     5556#endif
    49035557 
    49045558       gu = 2.0_wp * u_gtrans
    49055559       gv = 2.0_wp * v_gtrans
     5560
     5561#ifndef _OPENACC
    49065562!
    49075563!--   compute the whole left boundary of the processor domain
     
    49625618
    49635619       ENDDO
    4964 
     5620#endif
     5621
     5622       !$ACC PARALLEL LOOP COLLAPSE(2) FIRSTPRIVATE(tn, gu, gv) &
     5623       !$ACC PRIVATE(i, j, k, k_mm, k_pp, k_ppp) &
     5624       !$ACC PRIVATE(ibit27, ibit28, ibit29, ibit30, ibit31, ibit32) &
     5625       !$ACC PRIVATE(ibit27_l, ibit28_l, ibit29_l) &
     5626       !$ACC PRIVATE(ibit30_s, ibit31_s, ibit32_s) &
     5627       !$ACC PRIVATE(ibit33, ibit34, ibit35) &
     5628       !$ACC PRIVATE(flux_r, diss_r, flux_l, diss_l) &
     5629       !$ACC PRIVATE(flux_n, diss_n, flux_s, diss_s) &
     5630       !$ACC PRIVATE(flux_t, diss_t, flux_d, diss_d) &
     5631       !$ACC PRIVATE(div, u_comp, u_comp_l, v_comp, v_comp_s, w_comp) &
     5632       !$ACC PRESENT(advc_flags_1, advc_flags_2) &
     5633       !$ACC PRESENT(u, v, w) &
     5634       !$ACC PRESENT(rho_air, drho_air_zw, ddzu) &
     5635       !$ACC PRESENT(tend) &
     5636       !$ACC PRESENT(hom(nzb+1:nzb_max,1,3,0)) &
     5637       !$ACC PRESENT(weight_substep(intermediate_timestep_count)) &
     5638       !$ACC PRESENT(sums_ws2_ws_l(nzb+1:nzb_max,0))
    49655639       DO i = nxl, nxr
    49665640
     5641#ifndef _OPENACC
    49675642          j = nys
    49685643          DO  k = nzb+1, nzb_max
     
    50185693
    50195694          ENDDO
     5695#endif
    50205696
    50215697          DO  j = nys, nyn
     
    50255701!--          at the first w-level. For topography wall this is done implicitely
    50265702!--          by advc_flags_1.
    5027              k         = nzb + 1
    5028              w_comp    = w(k,j,i) + w(k-1,j,i)
    5029              flux_t(0) = w_comp       * ( w(k,j,i) + w(k-1,j,i) ) * adv_mom_1
    5030              diss_t(0) = -ABS(w_comp) * ( w(k,j,i) - w(k-1,j,i) ) * adv_mom_1
    5031              flux_d    = flux_t(0)
    5032              diss_d    = diss_t(0)
     5703             k      = nzb + 1
     5704             w_comp = w(k,j,i) + w(k-1,j,i)
     5705             flux_d = w_comp       * ( w(k,j,i) + w(k-1,j,i) ) * adv_mom_1
     5706             diss_d = -ABS(w_comp) * ( w(k,j,i) - w(k-1,j,i) ) * adv_mom_1
    50335707
    50345708             DO  k = nzb+1, nzb_max
     
    50385712                ibit27 = REAL( IBITS(advc_flags_1(k,j,i),27,1), KIND = wp )
    50395713
    5040                 u_comp    = u(k+1,j,i+1) + u(k,j,i+1) - gu
    5041                 flux_r(k) = u_comp * (                                       &
     5714                u_comp = u(k+1,j,i+1) + u(k,j,i+1) - gu
     5715                flux_r = u_comp * (                                          &
    50425716                          ( 37.0_wp * ibit29 * adv_mom_5                        &
    50435717                       +     7.0_wp * ibit28 * adv_mom_3                        &
     
    50545728                                     )
    50555729
    5056                 diss_r(k) = - ABS( u_comp ) * (                              &
     5730                diss_r = - ABS( u_comp ) * (                                 &
    50575731                          ( 10.0_wp * ibit29 * adv_mom_5                        &
    50585732                       +     3.0_wp * ibit28 * adv_mom_3                        &
     
    50695743                                              )
    50705744
     5745#ifdef _OPENACC
     5746!
     5747!--             Recompute the left fluxes.
     5748                ibit29_l = REAL( IBITS(advc_flags_1(k,j,i-1),29,1), KIND = wp )
     5749                ibit28_l = REAL( IBITS(advc_flags_1(k,j,i-1),28,1), KIND = wp )
     5750                ibit27_l = REAL( IBITS(advc_flags_1(k,j,i-1),27,1), KIND = wp )
     5751
     5752                u_comp_l = u(k+1,j,i) + u(k,j,i) - gu
     5753                flux_l   = u_comp_l * (                                      &
     5754                                      ( 37.0_wp * ibit29_l * adv_mom_5          &
     5755                                   +     7.0_wp * ibit28_l * adv_mom_3          &
     5756                                   +              ibit27_l * adv_mom_1          &
     5757                                      ) *                                    &
     5758                                     ( w(k,j,i)   + w(k,j,i-1) )             &
     5759                               -      (  8.0_wp * ibit29_l * adv_mom_5          &
     5760                                   +              ibit28_l * adv_mom_3          &
     5761                                      ) *                                    &
     5762                                     ( w(k,j,i+1) + w(k,j,i-2) )             &
     5763                               +      (           ibit29_l * adv_mom_5          &
     5764                                      ) *                                    &
     5765                                     ( w(k,j,i+2) + w(k,j,i-3) )             &
     5766                                                 )
     5767
     5768                diss_l   = - ABS( u_comp_l ) * (                             &
     5769                                        ( 10.0_wp * ibit29_l * adv_mom_5        &
     5770                                     +     3.0_wp * ibit28_l * adv_mom_3        &
     5771                                     +              ibit27_l * adv_mom_1        &
     5772                                        ) *                                  &
     5773                                     ( w(k,j,i)   - w(k,j,i-1) )             &
     5774                                 -      (  5.0_wp * ibit29_l * adv_mom_5        &
     5775                                     +              ibit28_l * adv_mom_3        &
     5776                                        ) *                                  &
     5777                                     ( w(k,j,i+1) - w(k,j,i-2) )             &
     5778                                 +      (           ibit29_l * adv_mom_5        &
     5779                                        ) *                                  &
     5780                                     ( w(k,j,i+2) - w(k,j,i-3) )             &
     5781                                                            )
     5782#else
     5783                flux_l = swap_flux_x_local_w(k,j)
     5784                diss_l = swap_diss_x_local_w(k,j)
     5785#endif
     5786
     5787
    50715788                ibit32 = REAL( IBITS(advc_flags_2(k,j,i),0,1),  KIND = wp )
    50725789                ibit31 = REAL( IBITS(advc_flags_1(k,j,i),31,1), KIND = wp )
    50735790                ibit30 = REAL( IBITS(advc_flags_1(k,j,i),30,1), KIND = wp )
    50745791
    5075                 v_comp    = v(k+1,j+1,i) + v(k,j+1,i) - gv
    5076                 flux_n(k) = v_comp * (                                       &
     5792                v_comp = v(k+1,j+1,i) + v(k,j+1,i) - gv
     5793                flux_n = v_comp * (                                          &
    50775794                          ( 37.0_wp * ibit32 * adv_mom_5                        &
    50785795                       +     7.0_wp * ibit31 * adv_mom_3                        &
     
    50895806                                     )
    50905807
    5091                 diss_n(k) = - ABS( v_comp ) * (                              &
     5808                diss_n = - ABS( v_comp ) * (                                 &
    50925809                          ( 10.0_wp * ibit32 * adv_mom_5                        &
    50935810                       +     3.0_wp * ibit31 * adv_mom_3                        &
     
    51035820                                 ( w(k,j+3,i) - w(k,j-2,i) )                 &
    51045821                                              )
     5822
     5823#ifdef _OPENACC
     5824!
     5825!--             Recompute the south fluxes.
     5826                ibit32_s = REAL( IBITS(advc_flags_2(k,j-1,i),0,1),  KIND = wp )
     5827                ibit31_s = REAL( IBITS(advc_flags_1(k,j-1,i),31,1), KIND = wp )
     5828                ibit30_s = REAL( IBITS(advc_flags_1(k,j-1,i),30,1), KIND = wp )
     5829
     5830                v_comp_s = v(k+1,j,i) + v(k,j,i) - gv
     5831                flux_s   = v_comp_s * (                                      &
     5832                                    ( 37.0_wp * ibit32_s * adv_mom_5            &
     5833                                 +     7.0_wp * ibit31_s * adv_mom_3            &
     5834                                 +              ibit30_s * adv_mom_1            &
     5835                                    ) *                                      &
     5836                                     ( w(k,j,i)   + w(k,j-1,i) )             &
     5837                             -      (  8.0_wp * ibit32_s * adv_mom_5            &
     5838                                 +              ibit31_s * adv_mom_3            &
     5839                                    ) *                                      &
     5840                                     ( w(k,j+1,i) + w(k,j-2,i) )             &
     5841                             +      (           ibit32_s * adv_mom_5            &
     5842                                    ) *                                      &
     5843                                     ( w(k,j+2,i) + w(k,j-3,i) )             &
     5844                                               )
     5845
     5846                diss_s   = - ABS( v_comp_s ) * (                             &
     5847                                    ( 10.0_wp * ibit32_s * adv_mom_5            &
     5848                                 +     3.0_wp * ibit31_s * adv_mom_3            &
     5849                                 +              ibit30_s * adv_mom_1            &
     5850                                    ) *                                      &
     5851                                     ( w(k,j,i)   - w(k,j-1,i) )             &
     5852                             -      (  5.0_wp * ibit32_s * adv_mom_5            &
     5853                                 +              ibit31_s * adv_mom_3            &
     5854                                    ) *                                      &
     5855                                     ( w(k,j+1,i) - w(k,j-2,i) )             &
     5856                             +      (           ibit32_s * adv_mom_5            &
     5857                                    ) *                                      &
     5858                                     ( w(k,j+2,i) - w(k,j-3,i) )             &
     5859                                                        )
     5860#else
     5861                flux_s = swap_flux_y_local_w(k)
     5862                diss_s = swap_diss_y_local_w(k)
     5863#endif
     5864
    51055865!
    51065866!--             k index has to be modified near bottom and top, else array
     
    51145874                k_mm  = k - 2 * ibit35
    51155875
    5116                 w_comp    = w(k+1,j,i) + w(k,j,i)
    5117                 flux_t(k) = w_comp * rho_air(k+1) * (                        &
     5876                w_comp = w(k+1,j,i) + w(k,j,i)
     5877                flux_t = w_comp * rho_air(k+1) * (                           &
    51185878                          ( 37.0_wp * ibit35 * adv_mom_5                        &
    51195879                       +     7.0_wp * ibit34 * adv_mom_3                        &
     
    51305890                                       )
    51315891
    5132                 diss_t(k) = - ABS( w_comp ) * rho_air(k+1) * (               &
     5892                diss_t = - ABS( w_comp ) * rho_air(k+1) * (                  &
    51335893                          ( 10.0_wp * ibit35 * adv_mom_5                        &
    51345894                       +     3.0_wp * ibit34 * adv_mom_3                        &
     
    51775937
    51785938                tend(k,j,i) = tend(k,j,i) - (                                 &
    5179                       ( flux_r(k) + diss_r(k)                                 &
    5180                     -   swap_flux_x_local_w(k,j) - swap_diss_x_local_w(k,j)   &
    5181                       ) * ddx                                                 &
    5182                     + ( flux_n(k) + diss_n(k)                                 &
    5183                     -   swap_flux_y_local_w(k)   - swap_diss_y_local_w(k)     &
    5184                       ) * ddy                                                 &
    5185                     + ( ( flux_t(k) + diss_t(k) )                             &
    5186                     -   ( flux_d    + diss_d    )                             &
    5187                       ) * drho_air_zw(k) * ddzu(k+1)                          &
     5939                      ( ( flux_r + diss_r )                                   &
     5940                    -   ( flux_l + diss_l ) ) * ddx                           &
     5941                    + ( ( flux_n + diss_n )                                   &
     5942                    -   ( flux_s + diss_s ) ) * ddy                           &
     5943                    + ( ( flux_t + diss_t )                                   &
     5944                    -   ( flux_d + diss_d ) ) * drho_air_zw(k) * ddzu(k+1)    &
    51885945                                            )  + div * w(k,j,i)
    51895946
    5190                 swap_flux_x_local_w(k,j) = flux_r(k)
    5191                 swap_diss_x_local_w(k,j) = diss_r(k)
    5192                 swap_flux_y_local_w(k)   = flux_n(k)
    5193                 swap_diss_y_local_w(k)   = diss_n(k)
    5194                 flux_d                   = flux_t(k)
    5195                 diss_d                   = diss_t(k)
    5196 
     5947#ifndef _OPENACC
     5948                swap_flux_x_local_w(k,j) = flux_r
     5949                swap_diss_x_local_w(k,j) = diss_r
     5950                swap_flux_y_local_w(k)   = flux_n
     5951                swap_diss_y_local_w(k)   = diss_n
     5952#endif
     5953                flux_d                   = flux_t
     5954                diss_d                   = diss_t
     5955
     5956                !$ACC ATOMIC
    51975957                sums_ws2_ws_l(k,tn)  = sums_ws2_ws_l(k,tn)                    &
    5198                       + ( flux_t(k)                                           &
     5958                      + ( flux_t                                              &
    51995959                       * ( w_comp - 2.0_wp * hom(k,1,3,0)                   ) &
    52005960                       / ( w_comp + SIGN( 1.0E-20_wp, w_comp )              ) &
    5201                         + diss_t(k)                                           &
     5961                        + diss_t                                              &
    52025962                       *   ABS( w_comp - 2.0_wp * hom(k,1,3,0)              ) &
    52035963                       / ( ABS( w_comp ) + 1.0E-20_wp                       ) &
     
    52085968             DO  k = nzb_max+1, nzt
    52095969
    5210                 u_comp    = u(k+1,j,i+1) + u(k,j,i+1) - gu
    5211                 flux_r(k) = u_comp * (                                      &
     5970                u_comp = u(k+1,j,i+1) + u(k,j,i+1) - gu
     5971                flux_r = u_comp * (                                         &
    52125972                      37.0_wp * ( w(k,j,i+1) + w(k,j,i)   )                    &
    52135973                    -  8.0_wp * ( w(k,j,i+2) + w(k,j,i-1) )                    &
    52145974                    +           ( w(k,j,i+3) + w(k,j,i-2) ) ) * adv_mom_5
    52155975
    5216                 diss_r(k) = - ABS( u_comp ) * (                             &
     5976                diss_r = - ABS( u_comp ) * (                                &
    52175977                      10.0_wp * ( w(k,j,i+1) - w(k,j,i)   )                    &
    52185978                    -  5.0_wp * ( w(k,j,i+2) - w(k,j,i-1) )                    &
    52195979                    +           ( w(k,j,i+3) - w(k,j,i-2) ) ) * adv_mom_5
    52205980
    5221                 v_comp    = v(k+1,j+1,i) + v(k,j+1,i) - gv
    5222                 flux_n(k) = v_comp * (                                      &
     5981#ifdef _OPENACC
     5982!
     5983!--             Recompute the left fluxes.
     5984                u_comp_l = u(k+1,j,i) + u(k,j,i) - gu
     5985                flux_l   = u_comp_l * (                                     &
     5986                            37.0_wp * ( w(k,j,i) + w(k,j,i-1)   )              &
     5987                          -  8.0_wp * ( w(k,j,i+1) + w(k,j,i-2) )              &
     5988                          +           ( w(k,j,i+2) + w(k,j,i-3) ) ) * adv_mom_5
     5989                diss_l   = - ABS( u_comp_l ) * (                            &
     5990                            10.0_wp * ( w(k,j,i) - w(k,j,i-1)   )              &
     5991                          -  5.0_wp * ( w(k,j,i+1) - w(k,j,i-2) )              &
     5992                          +           ( w(k,j,i+2) - w(k,j,i-3) ) ) * adv_mom_5
     5993#else
     5994                flux_l = swap_flux_x_local_w(k,j)
     5995                diss_l = swap_diss_x_local_w(k,j)
     5996#endif
     5997
     5998                v_comp = v(k+1,j+1,i) + v(k,j+1,i) - gv
     5999                flux_n = v_comp * (                                         &
    52236000                      37.0_wp * ( w(k,j+1,i) + w(k,j,i)   )                    &
    52246001                    -  8.0_wp * ( w(k,j+2,i) + w(k,j-1,i) )                    &
    52256002                    +           ( w(k,j+3,i) + w(k,j-2,i) ) ) * adv_mom_5
    52266003
    5227                 diss_n(k) = - ABS( v_comp ) * (                             &
     6004                diss_n = - ABS( v_comp ) * (                                &
    52286005                      10.0_wp * ( w(k,j+1,i) - w(k,j,i)   )                    &
    52296006                    -  5.0_wp * ( w(k,j+2,i) - w(k,j-1,i) )                    &
    52306007                    +           ( w(k,j+3,i) - w(k,j-2,i) ) ) * adv_mom_5
     6008
     6009#ifdef _OPENACC
     6010!
     6011!--             Recompute the south fluxes.
     6012                v_comp_s = v(k+1,j,i) + v(k,j,i) - gv
     6013                flux_s   = v_comp_s * (                                     &
     6014                           37.0_wp * ( w(k,j,i) + w(k,j-1,i)   )               &
     6015                         -  8.0_wp * ( w(k,j+1,i) +w(k,j-2,i)  )               &
     6016                         +           ( w(k,j+2,i) + w(k,j-3,i) ) ) * adv_mom_5
     6017                diss_s   = - ABS( v_comp_s ) * (                            &
     6018                           10.0_wp * ( w(k,j,i) - w(k,j-1,i)   )               &
     6019                         -  5.0_wp * ( w(k,j+1,i) - w(k,j-2,i) )               &
     6020                         +           ( w(k,j+2,i) - w(k,j-3,i) ) ) * adv_mom_5
     6021#else
     6022                flux_s = swap_flux_y_local_w(k)
     6023                diss_s = swap_diss_y_local_w(k)
     6024#endif
     6025
    52316026!
    52326027!--             k index has to be modified near bottom and top, else array
     
    52406035                k_mm  = k - 2 * ibit35
    52416036
    5242                 w_comp    = w(k+1,j,i) + w(k,j,i)
    5243                 flux_t(k) = w_comp * rho_air(k+1) * (                        &
     6037                w_comp = w(k+1,j,i) + w(k,j,i)
     6038                flux_t = w_comp * rho_air(k+1) * (                           &
    52446039                          ( 37.0_wp * ibit35 * adv_mom_5                        &
    52456040                       +     7.0_wp * ibit34 * adv_mom_3                        &
     
    52566051                                       )
    52576052
    5258                 diss_t(k) = - ABS( w_comp ) * rho_air(k+1) * (               &
     6053                diss_t = - ABS( w_comp ) * rho_air(k+1) * (                  &
    52596054                          ( 10.0_wp * ibit35 * adv_mom_5                        &
    52606055                       +     3.0_wp * ibit34 * adv_mom_3                        &
     
    52826077
    52836078                tend(k,j,i) = tend(k,j,i) - (                                 &
    5284                       ( flux_r(k) + diss_r(k)                                 &
    5285                     -   swap_flux_x_local_w(k,j) - swap_diss_x_local_w(k,j)   &
    5286                       ) * ddx                                                 &
    5287                     + ( flux_n(k) + diss_n(k)                                 &
    5288                     -   swap_flux_y_local_w(k)   - swap_diss_y_local_w(k)     &
    5289                       ) * ddy                                                 &
    5290                     + ( ( flux_t(k) + diss_t(k) )                             &
    5291                     -   ( flux_d    + diss_d    )                             &
    5292                       ) * drho_air_zw(k) * ddzu(k+1)                          &
     6079                      ( ( flux_r + diss_r )                                   &
     6080                    -   ( flux_l + diss_l ) ) * ddx                           &
     6081                    + ( ( flux_n + diss_n )                                   &
     6082                    -   ( flux_s + diss_s ) ) * ddy                           &
     6083                    + ( ( flux_t + diss_t )                                   &
     6084                    -   ( flux_d + diss_d ) ) * drho_air_zw(k) * ddzu(k+1)    &
    52936085                                            )  + div * w(k,j,i)
    52946086
    5295                 swap_flux_x_local_w(k,j) = flux_r(k)
    5296                 swap_diss_x_local_w(k,j) = diss_r(k)
    5297                 swap_flux_y_local_w(k)   = flux_n(k)
    5298                 swap_diss_y_local_w(k)   = diss_n(k)
    5299                 flux_d                   = flux_t(k)
    5300                 diss_d                   = diss_t(k)
    5301 
     6087#ifndef _OPENACC
     6088                swap_flux_x_local_w(k,j) = flux_r
     6089                swap_diss_x_local_w(k,j) = diss_r
     6090                swap_flux_y_local_w(k)   = flux_n
     6091                swap_diss_y_local_w(k)   = diss_n
     6092#endif
     6093                flux_d                   = flux_t
     6094                diss_d                   = diss_t
     6095
     6096                !$ACC ATOMIC
    53026097                sums_ws2_ws_l(k,tn)  = sums_ws2_ws_l(k,tn)                    &
    5303                       + ( flux_t(k)                                           &
     6098                      + ( flux_t                                              &
    53046099                       * ( w_comp - 2.0_wp * hom(k,1,3,0)                   ) &
    53056100                       / ( w_comp + SIGN( 1.0E-20_wp, w_comp )              ) &
    5306                         + diss_t(k)                                           &
     6101                        + diss_t                                              &
    53076102                       *   ABS( w_comp - 2.0_wp * hom(k,1,3,0)              ) &
    53086103                       / ( ABS( w_comp ) + 1.0E-20_wp                       ) &
Note: See TracChangeset for help on using the changeset viewer.