From 2107beb1926c66fee40d7bceca98eafafcf9e680 Mon Sep 17 00:00:00 2001 From: Edoardo Zoni <59625522+EZoni@users.noreply.github.com> Date: Mon, 12 Jul 2021 12:37:52 -0700 Subject: Add Cost Calculations for Cartesian/RZ Filtering (#2074) * Add Cost Calculations to ApplyFilter * Add Cost Calculations to ApplyStencil * Update Doxygen --- .../SpectralSolver/SpectralFieldDataRZ.cpp | 35 ++++++++++++++++++++-- 1 file changed, 33 insertions(+), 2 deletions(-) (limited to 'Source/FieldSolver/SpectralSolver/SpectralFieldDataRZ.cpp') 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* 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* 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); + } } } -- cgit v1.2.3