aboutsummaryrefslogtreecommitdiff
path: root/Source/Particles/WarpXParticleContainer.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'Source/Particles/WarpXParticleContainer.cpp')
-rw-r--r--Source/Particles/WarpXParticleContainer.cpp195
1 files changed, 56 insertions, 139 deletions
diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp
index 5905ccdb7..b52e31e42 100644
--- a/Source/Particles/WarpXParticleContainer.cpp
+++ b/Source/Particles/WarpXParticleContainer.cpp
@@ -9,6 +9,7 @@
*/
#include "WarpXParticleContainer.H"
+#include "ABLASTR/DepositCharge.H"
#include "Deposition/ChargeDeposition.H"
#include "Deposition/CurrentDeposition.H"
#include "Pusher/GetAndSetPosition.H"
@@ -598,146 +599,62 @@ WarpXParticleContainer::DepositCharge (WarpXParIter& pti, RealVector& wp,
const long offset, const long np_to_depose,
int thread_num, int lev, int depos_lev)
{
- AMREX_ALWAYS_ASSERT_WITH_MESSAGE((depos_lev==(lev-1)) ||
- (depos_lev==(lev )),
- "Deposition buffers only work for lev-1");
-
- // If no particles, do not do anything
- if (np_to_depose == 0) return;
-
- // If user decides not to deposit
- if (do_not_deposit) return;
-
- // Number of guard cells for local deposition of rho
- WarpX& warpx = WarpX::GetInstance();
- const amrex::IntVect& ng_rho = warpx.get_ng_depos_rho();
-
- // Extract deposition order and check that particles shape fits within the guard cells.
- // NOTE: In specific situations where the staggering of rho and the charge deposition algorithm
- // are not trivial, this check might be too strict and we might need to relax it, as currently
- // done for the current deposition.
-
-#if defined(WARPX_DIM_1D_Z)
- const amrex::IntVect shape_extent = amrex::IntVect(static_cast<int>(WarpX::noz/2+1));
-#elif defined(WARPX_DIM_XZ) || defined(WARPX_DIM_RZ)
- const amrex::IntVect shape_extent = amrex::IntVect(static_cast<int>(WarpX::nox/2+1),
- static_cast<int>(WarpX::noz/2+1));
-#elif defined(WARPX_DIM_3D)
- const amrex::IntVect shape_extent = amrex::IntVect(static_cast<int>(WarpX::nox/2+1),
- static_cast<int>(WarpX::noy/2+1),
- static_cast<int>(WarpX::noz/2+1));
-#endif
-
- // On CPU: particles deposit on tile arrays, which have a small number of guard cells ng_rho
- // On GPU: particles deposit directly on the rho array, which usually have a larger number of guard cells
-#ifndef AMREX_USE_GPU
- const amrex::IntVect range = ng_rho - shape_extent;
-#else
- const amrex::IntVect range = rho->nGrowVect() - shape_extent;
-#endif
-
- AMREX_ALWAYS_ASSERT_WITH_MESSAGE(
- amrex::numParticlesOutOfRange(pti, range) == 0,
- "Particles shape does not fit within tile (CPU) or guard cells (GPU) used for charge deposition");
-
- const std::array<Real,3>& dx = WarpX::CellSize(std::max(depos_lev,0));
- const Real q = this->charge;
-
- WARPX_PROFILE_VAR_NS("WarpXParticleContainer::DepositCharge::ChargeDeposition", blp_ppc_chd);
- WARPX_PROFILE_VAR_NS("WarpXParticleContainer::DepositCharge::Accumulate", blp_accumulate);
-
- // Get tile box where charge is deposited.
- // The tile box is different when depositing in the buffers (depos_lev<lev)
- // or when depositing inside the level (depos_lev=lev)
- Box tilebox;
- if (lev == depos_lev) {
- tilebox = pti.tilebox();
- } else {
- const IntVect& ref_ratio = WarpX::RefRatio(depos_lev);
- tilebox = amrex::coarsen(pti.tilebox(),ref_ratio);
- }
-
-#ifndef AMREX_USE_GPU
- // Staggered tile box
- Box tb = amrex::convert( tilebox, rho->ixType().toIntVect() );
-#endif
-
- tilebox.grow(ng_rho);
-
- const int nc = WarpX::ncomps;
-
-#ifdef AMREX_USE_GPU
- amrex::ignore_unused(thread_num);
- // GPU, no tiling: rho_fab points to the full rho array
- MultiFab rhoi(*rho, amrex::make_alias, icomp*nc, nc);
- auto & rho_fab = rhoi.get(pti);
-#else
- tb.grow(ng_rho);
-
- // CPU, tiling: rho_fab points to local_rho[thread_num]
- local_rho[thread_num].resize(tb, nc);
-
- // local_rho[thread_num] is set to zero
- local_rho[thread_num].setVal(0.0);
-
- auto & rho_fab = local_rho[thread_num];
-#endif
-
- const auto GetPosition = GetParticlePosition(pti, offset);
-
- // Lower corner of tile box physical domain
- // Note that this includes guard cells since it is after tilebox.ngrow
- Real cur_time = warpx.gett_new(lev);
- Real dt = warpx.getdt(lev);
- const auto& time_of_last_gal_shift = warpx.time_of_last_gal_shift;
- // Take into account Galilean shift
- Real time_shift_rho_old = (cur_time - time_of_last_gal_shift);
- Real time_shift_rho_new = (cur_time + dt - time_of_last_gal_shift);
- amrex::Array<amrex::Real,3> galilean_shift;
- if (icomp==0){
- galilean_shift = {
- m_v_galilean[0]*time_shift_rho_old,
- m_v_galilean[1]*time_shift_rho_old,
- m_v_galilean[2]*time_shift_rho_old };
- } else{
- galilean_shift = {
- m_v_galilean[0]*time_shift_rho_new,
- m_v_galilean[1]*time_shift_rho_new,
- m_v_galilean[2]*time_shift_rho_new };
- }
- const std::array<Real, 3>& xyzmin = WarpX::LowerCorner(tilebox, galilean_shift, depos_lev);
-
- // Indices of the lower bound
- const Dim3 lo = lbound(tilebox);
-
- WARPX_PROFILE_VAR_START(blp_ppc_chd);
- amrex::LayoutData<amrex::Real>* costs = WarpX::getCosts(lev);
- amrex::Real* cost = costs ? &((*costs)[pti.index()]) : nullptr;
-
- if (WarpX::nox == 1){
- doChargeDepositionShapeN<1>(GetPosition, wp.dataPtr()+offset, ion_lev,
- rho_fab, np_to_depose, dx, xyzmin, lo, q,
- WarpX::n_rz_azimuthal_modes, cost,
- WarpX::load_balance_costs_update_algo);
- } else if (WarpX::nox == 2){
- doChargeDepositionShapeN<2>(GetPosition, wp.dataPtr()+offset, ion_lev,
- rho_fab, np_to_depose, dx, xyzmin, lo, q,
- WarpX::n_rz_azimuthal_modes, cost,
- WarpX::load_balance_costs_update_algo);
- } else if (WarpX::nox == 3){
- doChargeDepositionShapeN<3>(GetPosition, wp.dataPtr()+offset, ion_lev,
- rho_fab, np_to_depose, dx, xyzmin, lo, q,
- WarpX::n_rz_azimuthal_modes, cost,
- WarpX::load_balance_costs_update_algo);
+ if (!do_not_deposit) {
+ WarpX& warpx = WarpX::GetInstance();
+ const amrex::IntVect& ng_rho = warpx.get_ng_depos_rho();
+ const std::array<amrex::Real,3>& dx = WarpX::CellSize(std::max(depos_lev,0));
+ amrex::IntVect ref_ratio;
+ if (lev == depos_lev) {
+ ref_ratio = IntVect(AMREX_D_DECL(1, 1, 1 ));
+ } else {
+ ref_ratio = WarpX::RefRatio(depos_lev);
+ }
+ const int nc = WarpX::ncomps;
+
+ // Get tile box where charge is deposited.
+ // The tile box is different when depositing in the buffers (depos_lev<lev)
+ // or when depositing inside the level (depos_lev=lev)
+ amrex::Box tilebox;
+ if (lev == depos_lev) {
+ tilebox = pti.tilebox();
+ } else {
+ tilebox = amrex::coarsen(pti.tilebox(), ref_ratio);
+ }
+ tilebox.grow(ng_rho);
+
+ // Lower corner of tile box physical domain
+ // Note that this includes guard cells since it is after tilebox.ngrow
+ // Take into account Galilean shift
+ const amrex::Real cur_time = warpx.gett_new(lev);
+ const amrex::Real dt = warpx.getdt(lev);
+ const amrex::Real time_of_last_gal_shift = warpx.time_of_last_gal_shift;
+ const amrex::Real time_shift_rho_old = (cur_time - time_of_last_gal_shift);
+ const amrex::Real time_shift_rho_new = (cur_time + dt - time_of_last_gal_shift);
+ amrex::Array<amrex::Real,3> galilean_shift;
+ if (icomp==0){
+ galilean_shift = {
+ m_v_galilean[0]*time_shift_rho_old,
+ m_v_galilean[1]*time_shift_rho_old,
+ m_v_galilean[2]*time_shift_rho_old };
+ } else{
+ galilean_shift = {
+ m_v_galilean[0]*time_shift_rho_new,
+ m_v_galilean[1]*time_shift_rho_new,
+ m_v_galilean[2]*time_shift_rho_new };
+ }
+ const auto& xyzmin = WarpX::LowerCorner(tilebox, galilean_shift, depos_lev);
+
+ // pointer to costs data
+ amrex::LayoutData<amrex::Real>* costs = WarpX::getCosts(lev);
+ amrex::Real* cost = costs ? &((*costs)[pti.index()]) : nullptr;
+
+ ablastr::DepositCharge<WarpXParticleContainer>
+ (pti, wp, ion_lev, rho, icomp, nc, offset, np_to_depose,
+ local_rho[thread_num], lev, depos_lev, this->charge,
+ WarpX::nox, WarpX::noy, WarpX::noz, ng_rho, dx, xyzmin, ref_ratio,
+ cost, WarpX::n_rz_azimuthal_modes, WarpX::load_balance_costs_update_algo,
+ WarpX::do_device_synchronize);
}
- WARPX_PROFILE_VAR_STOP(blp_ppc_chd);
-
-#ifndef AMREX_USE_GPU
- // CPU, tiling: atomicAdd local_rho into rho
- WARPX_PROFILE_VAR_START(blp_accumulate);
- (*rho)[pti].atomicAdd(local_rho[thread_num], tb, tb, 0, icomp*nc, nc);
- WARPX_PROFILE_VAR_STOP(blp_accumulate);
-#endif
}
void
generated by cgit v1.2.3 (git 2.46.0) at 2025-08-18 08:06:16 +0000