Ignore:
Timestamp:
Sep 27, 2012 9:23:24 AM (9 years ago)
Author:
raasch
Message:

Starting with changes required for GPU optimization. OpenACC statements for using NVIDIA GPUs added.
Adjustment of mixing length to the Prandtl mixing length at first grid point above ground removed.
mask array is set zero for ghost boundaries

File:
1 edited

Legend:

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

    r668 r1015  
    44! Current revisions:
    55! -----------------
     6! OpenACC statements added + code changes required for GPU optimization,
     7! adjustment of mixing length to the Prandtl mixing length at first grid point
     8! above ground removed
    69!
    710! Former revisions:
     
    5255    INTEGER ::  i, j, k, omp_get_thread_num, sr, tn
    5356
    54     REAL    ::  dvar_dz, l_stable, var_reference
    55 
    56     REAL, SAVE ::  phi_m = 1.0
     57    REAL    ::  dvar_dz, l, ll, l_stable, sqrt_e, var_reference
    5758
    5859    REAL    ::  var(nzb:nzt+1,nysg:nyng,nxlg:nxrg)
    59 
    60     REAL, DIMENSION(1:nzt) ::  l, ll, sqrt_e
    6160
    6261
     
    7170!
    7271!-- Compute the turbulent diffusion coefficient for momentum
    73     !$OMP PARALLEL PRIVATE (dvar_dz,i,j,k,l,ll,l_stable,phi_m,sqrt_e,sr,tn)
     72    !$OMP PARALLEL PRIVATE (dvar_dz,i,j,k,l,ll,l_stable,sqrt_e,sr,tn)
    7473!$  tn = omp_get_thread_num()
    7574
     75!
     76!-- Data declerations for accelerators
     77    !$acc data present( dd2zu, e, km, kh, l_grid, l_wall, nzb_s_inner, rif, var )
     78    !$acc kernels
     79
     80!
     81!-- Introduce an optional minimum tke
     82    IF ( e_min > 0.0 )  THEN
     83       !$OMP DO
     84       !$acc loop
     85       DO  i = nxlg, nxrg
     86          DO  j = nysg, nyng
     87             !$acc loop vector( 32 )
     88             DO  k = 1, nzt
     89                IF ( k > nzb_s_inner(j,i) )  THEN
     90                   e(k,j,i) = MAX( e(k,j,i), e_min )
     91                ENDIF
     92             ENDDO
     93          ENDDO
     94       ENDDO
     95    ENDIF
     96
    7697    !$OMP DO
     98    !$acc loop
    7799    DO  i = nxlg, nxrg
    78100       DO  j = nysg, nyng
     101          !$acc loop vector( 32 )
     102          DO  k = 1, nzt
    79103
     104             IF ( k > nzb_s_inner(j,i) )  THEN
     105
     106                sqrt_e = SQRT( e(k,j,i) )
    80107!
    81 !--       Compute the Phi-function for a possible adaption of the mixing length
    82 !--       to the Prandtl mixing length
    83           IF ( adjust_mixing_length  .AND.  prandtl_layer )  THEN
    84              IF ( rif(j,i) >= 0.0 )  THEN
    85                 phi_m = 1.0 + 5.0 * rif(j,i)
    86              ELSE
    87                 phi_m = 1.0 / SQRT( SQRT( 1.0 - 16.0 * rif(j,i) ) )
    88              ENDIF
    89           ENDIF
    90          
     108!--             Determine the mixing length
     109                dvar_dz = atmos_ocean_sign * &  ! inverse effect of pt/rho gradient
     110                          ( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k)
     111                IF ( dvar_dz > 0.0 ) THEN
     112                   IF ( use_reference )  THEN
     113                      l_stable = 0.76 * sqrt_e / &
     114                                 SQRT( g / var_reference * dvar_dz ) + 1E-5
     115                   ELSE
     116                      l_stable = 0.76 * sqrt_e / &
     117                                 SQRT( g / var(k,j,i) * dvar_dz ) + 1E-5
     118                   ENDIF
     119                ELSE
     120                   l_stable = l_grid(k)
     121                ENDIF
    91122!
    92 !--       Introduce an optional minimum tke
    93           IF ( e_min > 0.0 )  THEN
    94              DO  k = nzb_s_inner(j,i)+1, nzt
    95                 e(k,j,i) = MAX( e(k,j,i), e_min )
    96              ENDDO
    97           ENDIF
     123!--             Adjustment of the mixing length
     124                IF ( wall_adjustment )  THEN
     125                   l  = MIN( l_wall(k,j,i), l_grid(k), l_stable )
     126                   ll = MIN( l_wall(k,j,i), l_grid(k) )
     127                ELSE
     128                   l  = MIN( l_grid(k), l_stable )
     129                   ll = l_grid(k)
     130                ENDIF
    98131
    99 !
    100 !--       Calculate square root of e in a seperate loop, because it is used
    101 !--       twice in the next loop (better vectorization)
    102           DO  k = nzb_s_inner(j,i)+1, nzt
    103              sqrt_e(k) = SQRT( e(k,j,i) )
    104           ENDDO
     132      !
     133      !--       Compute diffusion coefficients for momentum and heat
     134                km(k,j,i) = 0.1 * l * sqrt_e
     135                kh(k,j,i) = ( 1.0 + 2.0 * l / ll ) * km(k,j,i)
    105136
    106 !
    107 !--       Determine the mixing length
    108           DO  k = nzb_s_inner(j,i)+1, nzt
    109              dvar_dz = atmos_ocean_sign * &  ! inverse effect of pt/rho gradient
    110                        ( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k)
    111              IF ( dvar_dz > 0.0 ) THEN
    112                 IF ( use_reference )  THEN
    113                    l_stable = 0.76 * sqrt_e(k) / &
    114                                      SQRT( g / var_reference * dvar_dz ) + 1E-5
    115                 ELSE
    116                    l_stable = 0.76 * sqrt_e(k) / &
    117                                      SQRT( g / var(k,j,i) * dvar_dz ) + 1E-5
    118                 ENDIF
    119              ELSE
    120                 l_stable = l_grid(k)
    121              ENDIF
    122 !
    123 !--          Adjustment of the mixing length
    124              IF ( wall_adjustment )  THEN
    125                 l(k)  = MIN( l_wall(k,j,i), l_grid(k), l_stable )
    126                 ll(k) = MIN( l_wall(k,j,i), l_grid(k) )
    127              ELSE
    128                 l(k)  = MIN( l_grid(k), l_stable )
    129                 ll(k) = l_grid(k)
    130              ENDIF
    131              IF ( adjust_mixing_length  .AND.  prandtl_layer )  THEN
    132                 l(k)  = MIN( l(k),  kappa * &
    133                                     ( zu(k) - zw(nzb_s_inner(j,i)) ) / phi_m )
    134                 ll(k) = MIN( ll(k), kappa * &
    135                                     ( zu(k) - zw(nzb_s_inner(j,i)) ) / phi_m )
    136137             ENDIF
    137138
     139#if ! defined( __openacc )
    138140!
    139 !--          Compute diffusion coefficients for momentum and heat
    140              km(k,j,i) = 0.1 * l(k) * sqrt_e(k)
    141              kh(k,j,i) = ( 1.0 + 2.0 * l(k) / ll(k) ) * km(k,j,i)
     141!++          Statistics still have to be realized for accelerators
     142!--          Summation for averaged profile (cf. flow_statistics)
     143             DO  sr = 0, statistic_regions
     144                sums_l_l(k,sr,tn) = sums_l_l(k,sr,tn) + l * rmask(j,i,sr)
     145             ENDDO
     146#endif
    142147
    143148          ENDDO
    144 
    145 !
    146 !--       Summation for averaged profile (cf. flow_statistics)
    147 !--       (the IF statement still requires a performance check on NEC machines)
    148           DO  sr = 0, statistic_regions
    149              IF ( rmask(j,i,sr) /= 0.0 .AND.  &
    150                   i >= nxl .AND. i <= nxr .AND. j >= nys .AND. j <= nyn )  THEN
    151                 DO  k = nzb_s_inner(j,i)+1, nzt
    152                    sums_l_l(k,sr,tn) = sums_l_l(k,sr,tn) + l(k)
    153                 ENDDO
    154              ENDIF
    155           ENDDO
    156 
    157149       ENDDO
    158150    ENDDO
    159151
     152#if ! defined( __openacc )
     153!
     154!++ Statistics still have to be realized for accelerators
    160155    sums_l_l(nzt+1,:,tn) = sums_l_l(nzt,:,tn)   ! quasi boundary-condition for
    161156                                                  ! data output
    162 
     157#endif
    163158    !$OMP END PARALLEL
    164159
     
    169164!-- values of the diffusivities are not needed
    170165    !$OMP PARALLEL DO
     166    !$acc loop
    171167    DO  i = nxlg, nxrg
    172168       DO  j = nysg, nyng
     
    198194    ENDIF
    199195
     196    !$acc end kernels
     197    !$acc end data
    200198
    201199 END SUBROUTINE diffusivities
Note: See TracChangeset for help on using the changeset viewer.