diff options
Diffstat (limited to 'Source/FieldSolver/SpectralSolver')
3 files changed, 41 insertions, 8 deletions
diff --git a/Source/FieldSolver/SpectralSolver/SpectralFieldDataRZ.H b/Source/FieldSolver/SpectralSolver/SpectralFieldDataRZ.H index 85aa41d3b..87158314d 100644 --- a/Source/FieldSolver/SpectralSolver/SpectralFieldDataRZ.H +++ b/Source/FieldSolver/SpectralSolver/SpectralFieldDataRZ.H @@ -67,8 +67,9 @@ class SpectralFieldDataRZ void InitFilter (amrex::IntVect const & filter_npass_each_dir, bool const compensation, SpectralKSpaceRZ const & k_space); - void ApplyFilter (int const field_index); - void ApplyFilter (int const field_index1, int const field_index2, int const field_index3); + void ApplyFilter (const int lev, int const field_index); + void ApplyFilter (const int lev, int const field_index1, + int const field_index2, int const field_index3); // Returns an array that holds the kr for all of the modes HankelTransform::RealVector const & getKrArray (amrex::MFIter const & mfi) const { diff --git a/Source/FieldSolver/SpectralSolver/SpectralFieldDataRZ.cpp b/Source/FieldSolver/SpectralSolver/SpectralFieldDataRZ.cpp index a87bfdb54..80760afb3 100644 --- a/Source/FieldSolver/SpectralSolver/SpectralFieldDataRZ.cpp +++ b/Source/FieldSolver/SpectralSolver/SpectralFieldDataRZ.cpp @@ -715,10 +715,18 @@ SpectralFieldDataRZ::InitFilter (amrex::IntVect const & filter_npass_each_dir, b /* \brief Apply K-space filtering on a scalar */ void -SpectralFieldDataRZ::ApplyFilter (int const field_index) +SpectralFieldDataRZ::ApplyFilter (const int lev, int const field_index) { + amrex::LayoutData<amrex::Real>* cost = WarpX::getCosts(lev); for (amrex::MFIter mfi(binomialfilter); mfi.isValid(); ++mfi){ + + if (cost && WarpX::load_balance_costs_update_algo == LoadBalanceCostsUpdateAlgo::Timers) + { + amrex::Gpu::synchronize(); + } + amrex::Real wt = amrex::second(); + auto const & filter_r = binomialfilter[mfi].getFilterArrayR(); auto const & filter_z = binomialfilter[mfi].getFilterArrayZ(); auto const & filter_r_arr = filter_r.dataPtr(); @@ -738,15 +746,31 @@ SpectralFieldDataRZ::ApplyFilter (int const field_index) int const ir = i + nr*mode; fields_arr(i,j,k,ic) *= filter_r_arr[ir]*filter_z_arr[j]; }); + + if (cost && WarpX::load_balance_costs_update_algo == LoadBalanceCostsUpdateAlgo::Timers) + { + amrex::Gpu::synchronize(); + wt = amrex::second() - wt; + amrex::HostDevice::Atomic::Add( &(*cost)[mfi.index()], wt); + } } } /* \brief Apply K-space filtering on a vector */ void -SpectralFieldDataRZ::ApplyFilter (int const field_index1, int const field_index2, int const field_index3) +SpectralFieldDataRZ::ApplyFilter (const int lev, int const field_index1, + int const field_index2, int const field_index3) { + amrex::LayoutData<amrex::Real>* cost = WarpX::getCosts(lev); for (amrex::MFIter mfi(binomialfilter); mfi.isValid(); ++mfi){ + + if (cost && WarpX::load_balance_costs_update_algo == LoadBalanceCostsUpdateAlgo::Timers) + { + amrex::Gpu::synchronize(); + } + amrex::Real wt = amrex::second(); + auto const & filter_r = binomialfilter[mfi].getFilterArrayR(); auto const & filter_z = binomialfilter[mfi].getFilterArrayZ(); auto const & filter_r_arr = filter_r.dataPtr(); @@ -770,5 +794,12 @@ SpectralFieldDataRZ::ApplyFilter (int const field_index1, int const field_index2 fields_arr(i,j,k,ic2) *= filter_r_arr[ir]*filter_z_arr[j]; fields_arr(i,j,k,ic3) *= filter_r_arr[ir]*filter_z_arr[j]; }); + + if (cost && WarpX::load_balance_costs_update_algo == LoadBalanceCostsUpdateAlgo::Timers) + { + amrex::Gpu::synchronize(); + wt = amrex::second() - wt; + amrex::HostDevice::Atomic::Add( &(*cost)[mfi.index()], wt); + } } } diff --git a/Source/FieldSolver/SpectralSolver/SpectralSolverRZ.H b/Source/FieldSolver/SpectralSolver/SpectralSolverRZ.H index 114d199c5..1c6b3d613 100644 --- a/Source/FieldSolver/SpectralSolver/SpectralSolverRZ.H +++ b/Source/FieldSolver/SpectralSolver/SpectralSolverRZ.H @@ -69,15 +69,16 @@ class SpectralSolverRZ } /* \brief Apply K space filtering for a scalar */ - void ApplyFilter (int const field_index) + void ApplyFilter (const int lev, int const field_index) { - field_data.ApplyFilter(field_index); + field_data.ApplyFilter(lev, field_index); } /* \brief Apply K space filtering for a vector */ - void ApplyFilter (int const field_index1, int const field_index2, int const field_index3) + void ApplyFilter (const int lev, int const field_index1, + int const field_index2, int const field_index3) { - field_data.ApplyFilter(field_index1, field_index2, field_index3); + field_data.ApplyFilter(lev, field_index1, field_index2, field_index3); } /** |