diff options
author | 2023-04-28 12:52:14 -0700 | |
---|---|---|
committer | 2023-04-28 12:52:14 -0700 | |
commit | 05f8542ab93c63ef09eedf57c72c86a86d304216 (patch) | |
tree | 2a065eed40af8ea754a3baa6982290432f07eb14 /Source/WarpX.cpp | |
parent | d7da74aff0d61fa4d753fb6459eec8cfcab28ac9 (diff) | |
download | WarpX-05f8542ab93c63ef09eedf57c72c86a86d304216.tar.gz WarpX-05f8542ab93c63ef09eedf57c72c86a86d304216.tar.zst WarpX-05f8542ab93c63ef09eedf57c72c86a86d304216.zip |
Build buffer masks on GPU (#3880)
Diffstat (limited to '')
-rw-r--r-- | Source/WarpX.cpp | 60 |
1 files changed, 13 insertions, 47 deletions
diff --git a/Source/WarpX.cpp b/Source/WarpX.cpp index 4afbb8087..b0d4adce6 100644 --- a/Source/WarpX.cpp +++ b/Source/WarpX.cpp @@ -2893,57 +2893,23 @@ void WarpX::BuildBufferMasksInBox ( const amrex::Box tbx, amrex::IArrayBox &buffer_mask, const amrex::IArrayBox &guard_mask, const int ng ) { - bool setnull; - const amrex::Dim3 lo = amrex::lbound( tbx ); - const amrex::Dim3 hi = amrex::ubound( tbx ); - Array4<int> msk = buffer_mask.array(); - Array4<int const> gmsk = guard_mask.array(); -#if defined(WARPX_DIM_1D_Z) - int k = lo.z; - int j = lo.y; - for (int i = lo.x; i <= hi.x; ++i) { - setnull = false; - // If gmsk=0 for any neighbor within ng cells, current cell is in the buffer region - for (int ii = i-ng; ii <= i+ng; ++ii) { - if ( gmsk(ii,j,k) == 0 ) setnull = true; - } - if ( setnull ) msk(i,j,k) = 0; - else msk(i,j,k) = 1; - } -#elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ) - int k = lo.z; - for (int j = lo.y; j <= hi.y; ++j) { - for (int i = lo.x; i <= hi.x; ++i) { - setnull = false; - // If gmsk=0 for any neighbor within ng cells, current cell is in the buffer region - for (int jj = j-ng; jj <= j+ng; ++jj) { - for (int ii = i-ng; ii <= i+ng; ++ii) { - if ( gmsk(ii,jj,k) == 0 ) setnull = true; - } - } - if ( setnull ) msk(i,j,k) = 0; - else msk(i,j,k) = 1; - } - } -#elif defined(WARPX_DIM_3D) - for (int k = lo.z; k <= hi.z; ++k) { - for (int j = lo.y; j <= hi.y; ++j) { - for (int i = lo.x; i <= hi.x; ++i) { - setnull = false; - // If gmsk=0 for any neighbor within ng cells, current cell is in the buffer region - for (int kk = k-ng; kk <= k+ng; ++kk) { - for (int jj = j-ng; jj <= j+ng; ++jj) { - for (int ii = i-ng; ii <= i+ng; ++ii) { - if ( gmsk(ii,jj,kk) == 0 ) setnull = true; - } + auto const& msk = buffer_mask.array(); + auto const& gmsk = guard_mask.const_array(); + amrex::Dim3 ng3 = amrex::IntVect(ng).dim3(); + amrex::ParallelFor(tbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) + { + for (int kk = k-ng3.z; kk <= k+ng3.z; ++kk) { + for (int jj = j-ng3.y; jj <= j+ng3.y; ++jj) { + for (int ii = i-ng3.x; ii <= i+ng3.x; ++ii) { + if ( gmsk(ii,jj,kk) == 0 ) { + msk(i,j,k) = 0; + return; } } - if ( setnull ) msk(i,j,k) = 0; - else msk(i,j,k) = 1; } } - } -#endif + msk(i,j,k) = 1; + }); } amrex::Vector<amrex::Real> WarpX::getFornbergStencilCoefficients(const int n_order, const short a_grid_type) |