From a1b2a770f24c6469e0f773b7cd1d5231c7ed82a4 Mon Sep 17 00:00:00 2001 From: MaxThevenet Date: Tue, 30 Jul 2019 09:28:39 -0700 Subject: some cleaning in current deposition --- Source/Particles/WarpXParticleContainer.cpp | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index a20f0035e..89f233b2c 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -532,17 +532,17 @@ WarpXParticleContainer::DepositCurrent(WarpXParIter& pti, if (WarpX::nox == 1){ doDepositionShapeN<1>(xp, yp, zp, wp.dataPtr(), uxp.dataPtr(), uyp.dataPtr(), uzp.dataPtr(), jx_arr, jy_arr, - jz_arr, offset, np_to_depose, dt, dx, + jz_arr, np_to_depose, dt, dx, xyzmin, lo, stagger_shift, q); } else if (WarpX::nox == 2){ doDepositionShapeN<2>(xp, yp, zp, wp.dataPtr(), uxp.dataPtr(), uyp.dataPtr(), uzp.dataPtr(), jx_arr, jy_arr, - jz_arr, offset, np_to_depose, dt, dx, + jz_arr, np_to_depose, dt, dx, xyzmin, lo, stagger_shift, q); } else if (WarpX::nox == 3){ doDepositionShapeN<3>(xp, yp, zp, wp.dataPtr(), uxp.dataPtr(), uyp.dataPtr(), uzp.dataPtr(), jx_arr, jy_arr, - jz_arr, offset, np_to_depose, dt, dx, + jz_arr, np_to_depose, dt, dx, xyzmin, lo, stagger_shift, q); } } -- cgit v1.2.3 From 720f05cd0dcf4bc6abd775181846c22b070d4dd9 Mon Sep 17 00:00:00 2001 From: Dave Grote Date: Wed, 31 Jul 2019 11:50:41 -0700 Subject: Implemented RZ current density inverse volume scaling in c++ --- Source/Particles/WarpXParticleContainer.cpp | 7 +++++++ 1 file changed, 7 insertions(+) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index 89f233b2c..d78837cf5 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -503,6 +503,7 @@ WarpXParticleContainer::DepositCurrent(WarpXParIter& pti, Real* AMREX_RESTRICT yp = m_yp[thread_num].dataPtr() + offset; // Lower corner of tile box physical domain + // Note that this includes guard cells since it is after tilebox.ngrow const std::array& xyzmin = WarpX::LowerCorner(tilebox, depos_lev);; // xyzmin is built on pti.tilebox(), so it does // not include staggering, so the stagger_shift has to be done by hand. @@ -547,6 +548,12 @@ WarpXParticleContainer::DepositCurrent(WarpXParIter& pti, } } +#ifdef WARPX_DIM_RZ + ApplyInverseVolumeScalingToCurrentDensity (jx_arr, jy_arr, jz_arr, + tbx, tby, tbz, + xyzmin[0], lo.x, dx[0], ngJ); +#endif + #ifndef AMREX_USE_GPU BL_PROFILE_VAR_START(blp_accumulate); // CPU, tiling: atomicAdd local_jx into jx -- cgit v1.2.3 From 2ddcc226e87e2d91b887bce09a1dc07b99685692 Mon Sep 17 00:00:00 2001 From: Dave Grote Date: Wed, 31 Jul 2019 13:05:01 -0700 Subject: Included ApplyInverseVolumeScalingToCurrentDensity --- Source/Particles/Deposition/CurrentDeposition.H | 83 ++++++++++++++++++++++++- Source/Particles/WarpXParticleContainer.cpp | 6 +- 2 files changed, 85 insertions(+), 4 deletions(-) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Particles/Deposition/CurrentDeposition.H b/Source/Particles/Deposition/CurrentDeposition.H index 9a6978dbe..16092d513 100644 --- a/Source/Particles/Deposition/CurrentDeposition.H +++ b/Source/Particles/Deposition/CurrentDeposition.H @@ -40,7 +40,7 @@ void doDepositionShapeN(const amrex::Real * const xp, const amrex::Real dzi = 1.0/dx[2]; const amrex::Real dts2dx = 0.5*dt*dxi; const amrex::Real dts2dz = 0.5*dt*dzi; -#if (defined WARPX_DIM_2D) +#if (AMREX_SPACEDIM == 2) const amrex::Real invvol = dxi*dzi; #elif (defined WARPX_DIM_3D) const amrex::Real dyi = 1.0/dx[1]; @@ -333,4 +333,85 @@ void doEsirkepovDepositionShapeN (const amrex::Real * const xp, ); } +#ifdef WARPX_RZ +/* \brief Apply the inverse volume scaling on the current density + * \param J[rtz]_arr : The current density arrays + * /param tb[rtz] : Tileboxes of the J grids + * \param rmin : Minimum radius of the grid + * \param irmin : Grid cell index at rmin + * \param dr : Radial grid cell size + * \param ngJ : Number of guard cells + */ +void ApplyInverseVolumeScalingToCurrentDensity (amrex::Array4 const& Jr_arr, + amrex::Array4 const& Jt_arr, + amrex::Array4 const& Jz_arr, + const amrex::Box& tbr, + const amrex::Box& tbt, + const amrex::Box& tbz, + const amrex::Real rmin, + const int irmin, + const amrex::Real dr, + const long ngJ) +{ + // Rescale current in r-z mode since the inverse volume factor was not + // included in the current deposition. + amrex::ParallelFor(tbr, + [=] AMREX_GPU_DEVICE (int i, int j, int k) + { + // Wrap the current density deposited in the guard cells around + // to the cells above the axis. + // Note that Jr(i==0) is at 1/2 dr. + if (rmin == 0. && 0 <= i && i < ngJ) { + Jr_arr(i,j,0) -= Jr_arr(-1-i,j,0); + } + // Apply the inverse volume scaling + // Since Jr is not cell centered in r, no need for distinction + // between on axis and off-axis factors + const amrex::Real r = std::abs(rmin + (i - irmin + 0.5)*dr); + Jr_arr(i,j,0) /= (2.*MathConst::pi*r); + }); + amrex::ParallelFor(tbt, + [=] AMREX_GPU_DEVICE (int i, int j, int k) + { + // Wrap the current density deposited in the guard cells around + // to the cells above the axis. + // Jt is located on the boundary + if (rmin == 0. && 0 < i && i <= ngJ) { + Jt_arr(i,j,0) += Jt_arr(-i,j,0); + } + + // Apply the inverse volume scaling + // Jt is forced to zero on axis. + const amrex::Real r = std::abs(rmin + (i - irmin)*dr); + if (r == 0.) { + Jt_arr(i,j,0) = 0.; + } else { + Jt_arr(i,j,0) /= (2.*MathConst::pi*r); + } + }); + amrex::ParallelFor(tbz, + [=] AMREX_GPU_DEVICE (int i, int j, int k) + { + // Wrap the current density deposited in the guard cells around + // to the cells above the axis. + // Jz is located on the boundary + if (rmin == 0. && 0 < i && i <= ngJ) { + Jz_arr(i,j,0) += Jz_arr(-i,j,0); + } + + // Apply the inverse volume scaling + const amrex::Real r = std::abs(rmin + (i - irmin)*dr); + if (r == 0.) { + // Verboncoeur JCP 164, 421-427 (2001) : corrected volume on axis + Jz_arr(i,j,0) /= (MathConst::pi*dr/3.); + // Alternative (when order in r is limited one, which is not implemented) + // Jz_arr(i,j,0) /= (MathConst::pi*dr/4.); + } else { + Jz_arr(i,j,0) /= (2.*MathConst::pi*r); + } + }); +} +#endif + + #endif // CURRENTDEPOSITION_H_ diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index d78837cf5..787c0e397 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -549,9 +549,9 @@ WarpXParticleContainer::DepositCurrent(WarpXParIter& pti, } #ifdef WARPX_DIM_RZ - ApplyInverseVolumeScalingToCurrentDensity (jx_arr, jy_arr, jz_arr, - tbx, tby, tbz, - xyzmin[0], lo.x, dx[0], ngJ); + ApplyInverseVolumeScalingToCurrentDensity(jx_arr, jy_arr, jz_arr, + tbx, tby, tbz, + xyzmin[0], lo.x, dx[0], ngJ); #endif #ifndef AMREX_USE_GPU -- cgit v1.2.3 From f5207b1bb5bc9dd9b0f2b7d4170bcf67fc14c8a2 Mon Sep 17 00:00:00 2001 From: Dave Grote Date: Thu, 1 Aug 2019 17:12:29 -0700 Subject: Moved ApplyInverseVolumeScalingToCurrentDensity to WarpX and is called after all particles are deposited. --- Source/Evolve/WarpXEvolveEM.cpp | 7 +++ Source/FieldSolver/WarpXPushFieldsEM.cpp | 96 +++++++++++++++++++++++++++++ Source/Particles/WarpXParticleContainer.cpp | 8 +-- Source/WarpX.H | 7 +++ 4 files changed, 111 insertions(+), 7 deletions(-) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Evolve/WarpXEvolveEM.cpp b/Source/Evolve/WarpXEvolveEM.cpp index 32a4747db..57a0c44c0 100644 --- a/Source/Evolve/WarpXEvolveEM.cpp +++ b/Source/Evolve/WarpXEvolveEM.cpp @@ -481,6 +481,13 @@ WarpX::PushParticlesandDepose (int lev, Real cur_time) Efield_cax[lev][0].get(), Efield_cax[lev][1].get(), Efield_cax[lev][2].get(), Bfield_cax[lev][0].get(), Bfield_cax[lev][1].get(), Bfield_cax[lev][2].get(), cur_time, dt[lev]); +#ifdef WARPX_DIM_RZ + // This is called after all particles have deposited their current. + ApplyInverseVolumeScalingToCurrentDensity(current_fp[lev][0].get(), current_fp[lev][1].get(), current_fp[lev][2].get(), lev); + if (current_buf[lev][0].get()) { + ApplyInverseVolumeScalingToCurrentDensity(current_buf[lev][0].get(), current_buf[lev][1].get(), current_buf[lev][2].get(), lev-1); + } +#endif } void diff --git a/Source/FieldSolver/WarpXPushFieldsEM.cpp b/Source/FieldSolver/WarpXPushFieldsEM.cpp index 4fce4717b..892c08aa2 100644 --- a/Source/FieldSolver/WarpXPushFieldsEM.cpp +++ b/Source/FieldSolver/WarpXPushFieldsEM.cpp @@ -560,3 +560,99 @@ WarpX::EvolveF (int lev, PatchType patch_type, Real a_dt, DtType a_dt_type) } } +#ifdef WARPX_DIM_RZ +// This scales the current by the inverse volume and wraps around the depostion at negative radius. +// It is faster to apply this on the grid than to do it particle by particle. +// It is put here since there isn't another nice place for it. +void +WarpX::ApplyInverseVolumeScalingToCurrentDensity (MultiFab* Jx, MultiFab* Jy, MultiFab* Jz, int lev) +{ + const long ngJ = Jx->nGrow(); + const std::array& dx = WarpX::CellSize(lev); + const Real dr = dx[0]; + + Box tilebox; + + for ( MFIter mfi(*Jx, TilingIfNotGPU()); mfi.isValid(); ++mfi ) + { + + Array4 const& Jr_arr = Jx->array(mfi); + Array4 const& Jt_arr = Jy->array(mfi); + Array4 const& Jz_arr = Jz->array(mfi); + + tilebox = mfi.tilebox(); + Box tbr = convert(tilebox, WarpX::jx_nodal_flag); + Box tbt = convert(tilebox, WarpX::jy_nodal_flag); + Box tbz = convert(tilebox, WarpX::jz_nodal_flag); + + // Lower corner of tile box physical domain + // Note that this is done before the tilebox.grow so that + // these do not include the guard cells. + const std::array& xyzmin = WarpX::LowerCorner(tilebox, lev); + const Dim3 lo = lbound(tilebox); + const Real rmin = xyzmin[0]; + const int irmin = lo.x; + + tilebox.grow(ngJ); + + // Rescale current in r-z mode since the inverse volume factor was not + // included in the current deposition. + amrex::ParallelFor(tbr, + [=] AMREX_GPU_DEVICE (int i, int j, int k) + { + // Wrap the current density deposited in the guard cells around + // to the cells above the axis. + // Note that Jr(i==0) is at 1/2 dr. + if (rmin == 0. && 0 <= i && i < ngJ) { + Jr_arr(i,j,0) -= Jr_arr(-1-i,j,0); + } + // Apply the inverse volume scaling + // Since Jr is not cell centered in r, no need for distinction + // between on axis and off-axis factors + const amrex::Real r = std::abs(rmin + (i - irmin + 0.5)*dr); + Jr_arr(i,j,0) /= (2.*MathConst::pi*r); + }); + amrex::ParallelFor(tbt, + [=] AMREX_GPU_DEVICE (int i, int j, int k) + { + // Wrap the current density deposited in the guard cells around + // to the cells above the axis. + // Jt is located on the boundary + if (rmin == 0. && 0 < i && i <= ngJ) { + Jt_arr(i,j,0) += Jt_arr(-i,j,0); + } + + // Apply the inverse volume scaling + // Jt is forced to zero on axis. + const amrex::Real r = std::abs(rmin + (i - irmin)*dr); + if (r == 0.) { + Jt_arr(i,j,0) = 0.; + } else { + Jt_arr(i,j,0) /= (2.*MathConst::pi*r); + } + }); + amrex::ParallelFor(tbz, + [=] AMREX_GPU_DEVICE (int i, int j, int k) + { + /* if (j == 5) std::cout << "Jz " << i << " " << Jz_arr(i,j,0) << "\n"; */ + // Wrap the current density deposited in the guard cells around + // to the cells above the axis. + // Jz is located on the boundary + if (rmin == 0. && 0 < i && i <= ngJ) { + Jz_arr(i,j,0) += Jz_arr(-i,j,0); + } + + // Apply the inverse volume scaling + const amrex::Real r = std::abs(rmin + (i - irmin)*dr); + if (r == 0.) { + // Verboncoeur JCP 164, 421-427 (2001) : corrected volume on axis + Jz_arr(i,j,0) /= (MathConst::pi*dr/3.); + // Alternative (when order in r is limited one, which is not implemented) + // Jz_arr(i,j,0) /= (MathConst::pi*dr/4.); + } else { + Jz_arr(i,j,0) /= (2.*MathConst::pi*r); + } + }); + } +} +#endif diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index 787c0e397..0e76353a4 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -504,7 +504,7 @@ WarpXParticleContainer::DepositCurrent(WarpXParIter& pti, // Lower corner of tile box physical domain // Note that this includes guard cells since it is after tilebox.ngrow - const std::array& xyzmin = WarpX::LowerCorner(tilebox, depos_lev);; + const std::array& xyzmin = WarpX::LowerCorner(tilebox, depos_lev); // xyzmin is built on pti.tilebox(), so it does // not include staggering, so the stagger_shift has to be done by hand. // Alternatively, we could define xyzminx from tbx (and the same for 3 @@ -548,12 +548,6 @@ WarpXParticleContainer::DepositCurrent(WarpXParIter& pti, } } -#ifdef WARPX_DIM_RZ - ApplyInverseVolumeScalingToCurrentDensity(jx_arr, jy_arr, jz_arr, - tbx, tby, tbz, - xyzmin[0], lo.x, dx[0], ngJ); -#endif - #ifndef AMREX_USE_GPU BL_PROFILE_VAR_START(blp_accumulate); // CPU, tiling: atomicAdd local_jx into jx diff --git a/Source/WarpX.H b/Source/WarpX.H index a25eef9e4..dde2278d7 100644 --- a/Source/WarpX.H +++ b/Source/WarpX.H @@ -178,6 +178,13 @@ public: void EvolveE (int lev, PatchType patch_type, amrex::Real dt); void EvolveF (int lev, PatchType patch_type, amrex::Real dt, DtType dt_type); +#ifdef WARPX_DIM_RZ + void ApplyInverseVolumeScalingToCurrentDensity(amrex::MultiFab* Jx, + amrex::MultiFab* Jy, + amrex::MultiFab* Jz, + int lev); +#endif + void DampPML (); void DampPML (int lev); void DampPML (int lev, PatchType patch_type); -- cgit v1.2.3 From 9b813235d4f27e4b95e0448a9e43cb82b35c2c89 Mon Sep 17 00:00:00 2001 From: Dave Grote Date: Thu, 1 Aug 2019 17:15:20 -0700 Subject: Fixed SetPosition for RZ with CUDA --- Source/Particles/WarpXParticleContainer.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index 0e76353a4..1cbcee75e 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -47,7 +47,7 @@ WarpXParIter::SetPosition (const Cuda::ManagedDeviceVector& x, const Cuda: #ifdef WARPX_RZ auto& attribs = GetAttribs(); auto& theta = attribs[PIdx::theta]; - Cuda::DeviceVector r(x.size()); + Cuda::ManagedDeviceVector r(x.size()); for (unsigned int i=0 ; i < x.size() ; i++) { theta[i] = std::atan2(y[i], x[i]); r[i] = std::sqrt(x[i]*x[i] + y[i]*y[i]); -- cgit v1.2.3 From eb6f9caf956935425da5664eba07a86355f32761 Mon Sep 17 00:00:00 2001 From: Dave Grote Date: Fri, 2 Aug 2019 18:36:58 -0700 Subject: Removed obsolete warpx_current_deposition_rz_volume_scaling --- Source/FortranInterface/WarpX_picsar.F90 | 50 ----------------------------- Source/Particles/WarpXParticleContainer.cpp | 8 ----- 2 files changed, 58 deletions(-) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/FortranInterface/WarpX_picsar.F90 b/Source/FortranInterface/WarpX_picsar.F90 index 3b2fc97a7..14bd79ad4 100644 --- a/Source/FortranInterface/WarpX_picsar.F90 +++ b/Source/FortranInterface/WarpX_picsar.F90 @@ -10,7 +10,6 @@ #define WRPX_PXR_GETEB_ENERGY_CONSERVING geteb2drz_energy_conserving_generic #define WRPX_PXR_CURRENT_DEPOSITION depose_jrjtjz_generic_rz #define WRPX_PXR_RZ_VOLUME_SCALING_RHO apply_rz_volume_scaling_rho -#define WRPX_PXR_RZ_VOLUME_SCALING_J apply_rz_volume_scaling_j #else @@ -355,53 +354,4 @@ subroutine warpx_charge_deposition(rho,np,xp,yp,zp,w,q,xmin,ymin,zmin,dx,dy,dz,n end subroutine warpx_current_deposition - ! _________________________________________________________________ - !> - !> @brief - !> Applies the inverse volume scaling for RZ current deposition - !> - !> @details - !> The scaling is done for single mode only - ! - !> @param[inout] jx,jy,jz current arrays - !> @param[in] jx_ntot,jy_ntot,jz_ntot vectors with total number of - !> cells (including guard cells) along each axis for each current - !> @param[in] jx_ng,jy_ng,jz_ng vectors with number of guard cells along each - !> axis for each current - !> @param[in] rmin tile grid minimum radius - !> @param[in] dr radial space discretization steps - !> - subroutine warpx_current_deposition_rz_volume_scaling( & - jx,jx_ng,jx_ntot,jy,jy_ng,jy_ntot,jz,jz_ng,jz_ntot, & - rmin,dr) & - bind(C, name="warpx_current_deposition_rz_volume_scaling") - - integer, intent(in) :: jx_ntot(AMREX_SPACEDIM), jy_ntot(AMREX_SPACEDIM), jz_ntot(AMREX_SPACEDIM) - integer(c_long), intent(in) :: jx_ng, jy_ng, jz_ng - real(amrex_real), intent(IN OUT):: jx(*), jy(*), jz(*) - real(amrex_real), intent(IN) :: rmin, dr - -#ifdef WARPX_RZ - integer(c_long) :: type_rz_depose = 1 -#endif - ! Compute the number of valid cells and guard cells - integer(c_long) :: jx_nvalid(AMREX_SPACEDIM), jy_nvalid(AMREX_SPACEDIM), jz_nvalid(AMREX_SPACEDIM), & - jx_nguards(AMREX_SPACEDIM), jy_nguards(AMREX_SPACEDIM), jz_nguards(AMREX_SPACEDIM) - jx_nvalid = jx_ntot - 2*jx_ng - jy_nvalid = jy_ntot - 2*jy_ng - jz_nvalid = jz_ntot - 2*jz_ng - jx_nguards = jx_ng - jy_nguards = jy_ng - jz_nguards = jz_ng - -#ifdef WARPX_RZ - CALL WRPX_PXR_RZ_VOLUME_SCALING_J( & - jx,jx_nguards,jx_nvalid, & - jy,jy_nguards,jy_nvalid, & - jz,jz_nguards,jz_nvalid, & - rmin,dr,type_rz_depose) -#endif - - end subroutine warpx_current_deposition_rz_volume_scaling - end module warpx_to_pxr_module diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index 1cbcee75e..f6c7afed5 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -394,14 +394,6 @@ WarpXParticleContainer::DepositCurrentFortran(WarpXParIter& pti, &WarpX::nox,&WarpX::noy,&WarpX::noz, &j_is_nodal, &lvect,&WarpX::current_deposition_algo); -#ifdef WARPX_RZ - // Rescale current in r-z mode - warpx_current_deposition_rz_volume_scaling( - jx_ptr, &ngJ, jxntot.getVect(), - jy_ptr, &ngJ, jyntot.getVect(), - jz_ptr, &ngJ, jzntot.getVect(), - &xyzmin[0], &dx[0]); -#endif BL_PROFILE_VAR_STOP(blp_pxr_cd); #ifndef AMREX_USE_GPU -- cgit v1.2.3 From ce63dce8b7d71c12cafe6309f6afc4c5cfc3d5c5 Mon Sep 17 00:00:00 2001 From: Dave Grote Date: Mon, 5 Aug 2019 09:43:57 -0700 Subject: Replaced WARPX_RZ with WARPX_DIM_RZ --- Source/Diagnostics/ParticleIO.cpp | 2 +- Source/Evolve/WarpXEvolveEM.cpp | 2 +- Source/FortranInterface/WarpX_f.H | 4 ++-- Source/FortranInterface/WarpX_picsar.F90 | 8 ++++---- Source/Make.WarpX | 1 - Source/Particles/PhysicalParticleContainer.cpp | 14 +++++++------- Source/Particles/Pusher/GetAndSetPosition.H | 6 +++--- Source/Particles/Pusher/UpdatePosition.H | 2 +- Source/Particles/WarpXParticleContainer.H | 2 +- Source/Particles/WarpXParticleContainer.cpp | 26 +++++++++++++------------- Source/Utils/WarpXAlgorithmSelection.cpp | 2 +- Source/WarpX.cpp | 16 ++++++++-------- 12 files changed, 42 insertions(+), 43 deletions(-) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Diagnostics/ParticleIO.cpp b/Source/Diagnostics/ParticleIO.cpp index f2a543ed5..f159e5302 100644 --- a/Source/Diagnostics/ParticleIO.cpp +++ b/Source/Diagnostics/ParticleIO.cpp @@ -98,7 +98,7 @@ MultiParticleContainer::WritePlotFile (const std::string& dir) const real_names.push_back("By"); real_names.push_back("Bz"); -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ real_names.push_back("theta"); #endif diff --git a/Source/Evolve/WarpXEvolveEM.cpp b/Source/Evolve/WarpXEvolveEM.cpp index 57a0c44c0..e9bf98f81 100644 --- a/Source/Evolve/WarpXEvolveEM.cpp +++ b/Source/Evolve/WarpXEvolveEM.cpp @@ -498,7 +498,7 @@ WarpX::ComputeDt () if (maxwell_fdtd_solver_id == 0) { // CFL time step Yee solver -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ // Derived semi-analytically by R. Lehe deltat = cfl * 1./( std::sqrt((1+0.2105)/(dx[0]*dx[0]) + 1./(dx[1]*dx[1])) * PhysConst::c ); #else diff --git a/Source/FortranInterface/WarpX_f.H b/Source/FortranInterface/WarpX_f.H index 689029059..07cfcb42e 100644 --- a/Source/FortranInterface/WarpX_f.H +++ b/Source/FortranInterface/WarpX_f.H @@ -62,7 +62,7 @@ #define WRPX_PUSH_LEAPFROG warpx_push_leapfrog_2d #define WRPX_PUSH_LEAPFROG_POSITIONS warpx_push_leapfrog_positions_2d -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ #define WRPX_COMPUTE_DIVE warpx_compute_dive_rz #else #define WRPX_COMPUTE_DIVE warpx_compute_dive_2d @@ -314,7 +314,7 @@ extern "C" const BL_FORT_FAB_ARG_ANYD(ey), const BL_FORT_FAB_ARG_ANYD(ez), const amrex::Real* dx -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ ,const amrex::Real* rmin #endif ); diff --git a/Source/FortranInterface/WarpX_picsar.F90 b/Source/FortranInterface/WarpX_picsar.F90 index e65c30e42..0625d40f9 100644 --- a/Source/FortranInterface/WarpX_picsar.F90 +++ b/Source/FortranInterface/WarpX_picsar.F90 @@ -4,7 +4,7 @@ #elif (AMREX_SPACEDIM == 2) -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ #define WRPX_PXR_CURRENT_DEPOSITION depose_jrjtjz_generic_rz #define WRPX_PXR_RZ_VOLUME_SCALING_RHO apply_rz_volume_scaling_rho @@ -138,7 +138,7 @@ subroutine warpx_charge_deposition(rho,np,xp,yp,zp,w,q,xmin,ymin,zmin,dx,dy,dz,n ! Dimension 2 #elif (AMREX_SPACEDIM==2) -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ logical(pxr_logical) :: l_2drz = .TRUE._c_long #else logical(pxr_logical) :: l_2drz = .FALSE._c_long @@ -175,7 +175,7 @@ subroutine warpx_charge_deposition(rho,np,xp,yp,zp,w,q,xmin,ymin,zmin,dx,dy,dz,n real(amrex_real), intent(IN OUT):: rho(*) real(amrex_real), intent(IN) :: rmin, dr -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ integer(c_long) :: type_rz_depose = 1 #endif @@ -184,7 +184,7 @@ subroutine warpx_charge_deposition(rho,np,xp,yp,zp,w,q,xmin,ymin,zmin,dx,dy,dz,n rho_nvalid = rho_ntot - 2*rho_ng rho_nguards = rho_ng -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ CALL WRPX_PXR_RZ_VOLUME_SCALING_RHO( & rho,rho_nguards,rho_nvalid, & rmin,dr,type_rz_depose) diff --git a/Source/Make.WarpX b/Source/Make.WarpX index d2551e059..069623b50 100644 --- a/Source/Make.WarpX +++ b/Source/Make.WarpX @@ -141,7 +141,6 @@ endif ifeq ($(USE_RZ),TRUE) USERSuffix := $(USERSuffix).RZ - DEFINES += -DWARPX_RZ endif ifeq ($(DO_ELECTROSTATIC),TRUE) diff --git a/Source/Particles/PhysicalParticleContainer.cpp b/Source/Particles/PhysicalParticleContainer.cpp index 08f8a77b4..7d90857b4 100644 --- a/Source/Particles/PhysicalParticleContainer.cpp +++ b/Source/Particles/PhysicalParticleContainer.cpp @@ -148,7 +148,7 @@ PhysicalParticleContainer::AddGaussianBeam(Real x_m, Real y_m, Real z_m, npart /= 4; } for (long i = 0; i < npart; ++i) { -#if ( AMREX_SPACEDIM == 3 | WARPX_RZ) +#if ( AMREX_SPACEDIM == 3 | WARPX_DIM_RZ) Real weight = q_tot/npart/charge; Real x = distx(mt); Real y = disty(mt); @@ -269,7 +269,7 @@ PhysicalParticleContainer::AddPlasma (int lev, RealBox part_realbox) if (!part_realbox.ok()) part_realbox = geom.ProbDomain(); int num_ppc = plasma_injector->num_particles_per_cell; -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ Real rmax = std::min(plasma_injector->xmax, part_realbox.hi(0)); #endif @@ -323,7 +323,7 @@ PhysicalParticleContainer::AddPlasma (int lev, RealBox part_realbox) Real density_min = plasma_injector->density_min; Real density_max = plasma_injector->density_max; -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ bool radially_weighted = plasma_injector->radially_weighted; #endif @@ -496,11 +496,11 @@ PhysicalParticleContainer::AddPlasma (int lev, RealBox part_realbox) #endif // Save the x and y values to use in the insideBounds checks. - // This is needed with WARPX_RZ since x and y are modified. + // This is needed with WARPX_DIM_RZ since x and y are modified. Real xb = x; Real yb = y; -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ // Replace the x and y, choosing the angle randomly. // These x and y are used to get the momentum and density Real theta = 2.*MathConst::pi*amrex::Random(); @@ -574,7 +574,7 @@ PhysicalParticleContainer::AddPlasma (int lev, RealBox part_realbox) // Real weight = dens * scale_fac / (AMREX_D_TERM(fac, *fac, *fac)); Real weight = dens * scale_fac; -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ if (radially_weighted) { weight *= 2.*MathConst::pi*xb; } else { @@ -603,7 +603,7 @@ PhysicalParticleContainer::AddPlasma (int lev, RealBox part_realbox) p.pos(1) = y; p.pos(2) = z; #elif (AMREX_SPACEDIM == 2) -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ pa[PIdx::theta][ip] = theta; #endif p.pos(0) = xb; diff --git a/Source/Particles/Pusher/GetAndSetPosition.H b/Source/Particles/Pusher/GetAndSetPosition.H index 42c61343e..3c74baeb2 100644 --- a/Source/Particles/Pusher/GetAndSetPosition.H +++ b/Source/Particles/Pusher/GetAndSetPosition.H @@ -5,7 +5,7 @@ #include #include -#ifndef WARPX_RZ +#ifndef WARPX_DIM_RZ /* \brief Extract the particle's coordinates from the ParticleType struct `p`, * and stores them in the variables `x`, `y`, `z`. */ @@ -42,7 +42,7 @@ void SetPosition( #endif } -# else // if WARPX_RZ is True +# elif defined WARPX_DIM_RZ /* \brief Extract the particle's coordinates from `theta` and the attributes * of the ParticleType struct `p` (which contains the radius), @@ -71,6 +71,6 @@ void SetCylindricalPositionFromCartesian( p.pos(1) = z; } -#endif // WARPX_RZ +#endif // WARPX_DIM_RZ #endif // WARPX_PARTICLES_PUSHER_GETANDSETPOSITION_H_ diff --git a/Source/Particles/Pusher/UpdatePosition.H b/Source/Particles/Pusher/UpdatePosition.H index 0a4f579f4..a9df63a30 100644 --- a/Source/Particles/Pusher/UpdatePosition.H +++ b/Source/Particles/Pusher/UpdatePosition.H @@ -20,7 +20,7 @@ void UpdatePosition( const amrex::Real inv_gamma = 1./std::sqrt(1. + (ux*ux + uy*uy + uz*uz)*inv_c2); // Update positions over one time step x += ux * inv_gamma * dt; -#if (AMREX_SPACEDIM == 3) || (defined WARPX_RZ) // RZ pushes particles in 3D +#if (AMREX_SPACEDIM == 3) || (defined WARPX_DIM_RZ) // RZ pushes particles in 3D y += uy * inv_gamma * dt; #endif z += uz * inv_gamma * dt; diff --git a/Source/Particles/WarpXParticleContainer.H b/Source/Particles/WarpXParticleContainer.H index 7cf53260a..a609b4cb3 100644 --- a/Source/Particles/WarpXParticleContainer.H +++ b/Source/Particles/WarpXParticleContainer.H @@ -13,7 +13,7 @@ struct PIdx enum { // Particle Attributes stored in amrex::ParticleContainer's struct of array w = 0, // weight ux, uy, uz, Ex, Ey, Ez, Bx, By, Bz, -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ theta, // RZ needs all three position components #endif nattribs diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index f6c7afed5..c5c6afc19 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -27,7 +27,7 @@ void WarpXParIter::GetPosition (Cuda::ManagedDeviceVector& x, Cuda::ManagedDeviceVector& y, Cuda::ManagedDeviceVector& z) const { amrex::ParIter<0,0,PIdx::nattribs>::GetPosition(x, z); -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ const auto& attribs = GetAttribs(); const auto& theta = attribs[PIdx::theta]; y.resize(x.size()); @@ -44,7 +44,7 @@ WarpXParIter::GetPosition (Cuda::ManagedDeviceVector& x, Cuda::ManagedDevi void WarpXParIter::SetPosition (const Cuda::ManagedDeviceVector& x, const Cuda::ManagedDeviceVector& y, const Cuda::ManagedDeviceVector& z) { -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ auto& attribs = GetAttribs(); auto& theta = attribs[PIdx::theta]; Cuda::ManagedDeviceVector r(x.size()); @@ -80,7 +80,7 @@ WarpXParticleContainer::WarpXParticleContainer (AmrCore* amr_core, int ispecies) particle_comps["Bx"] = PIdx::Bx; particle_comps["By"] = PIdx::By; particle_comps["Bz"] = PIdx::Bz; -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ particle_comps["theta"] = PIdx::theta; #endif @@ -163,7 +163,7 @@ WarpXParticleContainer::AddOneParticle (ParticleTileType& particle_tile, p.pos(1) = y; p.pos(2) = z; #elif (AMREX_SPACEDIM == 2) -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ attribs[PIdx::theta] = std::atan2(y, x); x = std::sqrt(x*x + y*y); #endif @@ -209,7 +209,7 @@ WarpXParticleContainer::AddNParticles (int lev, std::size_t np = iend-ibegin; -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ Vector theta(np); #endif @@ -228,7 +228,7 @@ WarpXParticleContainer::AddNParticles (int lev, p.pos(1) = y[i]; p.pos(2) = z[i]; #elif (AMREX_SPACEDIM == 2) -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ theta[i-ibegin] = std::atan2(y[i], x[i]); p.pos(0) = std::sqrt(x[i]*x[i] + y[i]*y[i]); #else @@ -265,7 +265,7 @@ WarpXParticleContainer::AddNParticles (int lev, for (int comp = PIdx::uz+1; comp < PIdx::nattribs; ++comp) { -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ if (comp == PIdx::theta) { particle_tile.push_back_real(comp, theta.front(), theta.back()); } @@ -610,7 +610,7 @@ WarpXParticleContainer::DepositCharge ( WarpXParIter& pti, RealVector& wp, &ngRho, &ngRho, &ngRho, &WarpX::nox,&WarpX::noy,&WarpX::noz, &lvect, &WarpX::charge_deposition_algo); -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ warpx_charge_deposition_rz_volume_scaling( data_ptr, &ngRho, rholen.getVect(), &xyzmin[0], &dx[0]); @@ -670,7 +670,7 @@ WarpXParticleContainer::DepositCharge ( WarpXParIter& pti, RealVector& wp, &ngRho, &ngRho, &ngRho, &WarpX::nox,&WarpX::noy,&WarpX::noz, &lvect, &WarpX::charge_deposition_algo); -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ warpx_charge_deposition_rz_volume_scaling( data_ptr, &ngRho, rholen.getVect(), &cxyzmin_tile[0], &cdx[0]); @@ -830,7 +830,7 @@ WarpXParticleContainer::GetChargeDensity (int lev, bool local) &dx[0], &dx[1], &dx[2], &nx, &ny, &nz, &nxg, &nyg, &nzg, &WarpX::nox,&WarpX::noy,&WarpX::noz, &lvect, &WarpX::charge_deposition_algo); -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ long ngRho = WarpX::nox; warpx_charge_deposition_rz_volume_scaling( data_ptr, &ngRho, rholen.getVect(), @@ -1015,7 +1015,7 @@ WarpXParticleContainer::PushX (int lev, Real dt) Real* AMREX_RESTRICT ux = attribs[PIdx::ux].dataPtr(); Real* AMREX_RESTRICT uy = attribs[PIdx::uy].dataPtr(); Real* AMREX_RESTRICT uz = attribs[PIdx::uz].dataPtr(); -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ Real* AMREX_RESTRICT theta = attribs[PIdx::theta].dataPtr(); #endif // Loop over the particles and update their position @@ -1023,12 +1023,12 @@ WarpXParticleContainer::PushX (int lev, Real dt) [=] AMREX_GPU_DEVICE (long i) { ParticleType& p = pstructs[i]; // Particle object that gets updated Real x, y, z; // Temporary variables -#ifndef WARPX_RZ +#ifndef WARPX_DIM_RZ GetPosition( x, y, z, p ); // Initialize x, y, z UpdatePosition( x, y, z, ux[i], uy[i], uz[i], dt); SetPosition( p, x, y, z ); // Update the object p #else - // For WARPX_RZ, the particles are still pushed in 3D Cartesian + // For WARPX_DIM_RZ, the particles are still pushed in 3D Cartesian GetCartesianPositionFromCylindrical( x, y, z, p, theta[i] ); UpdatePosition( x, y, z, ux[i], uy[i], uz[i], dt); SetCylindricalPositionFromCartesian( p, theta[i], x, y, z ); diff --git a/Source/Utils/WarpXAlgorithmSelection.cpp b/Source/Utils/WarpXAlgorithmSelection.cpp index 2c8038ccd..3aa4eb7b7 100644 --- a/Source/Utils/WarpXAlgorithmSelection.cpp +++ b/Source/Utils/WarpXAlgorithmSelection.cpp @@ -8,7 +8,7 @@ const std::map maxwell_solver_algo_to_int = { {"yee", MaxwellSolverAlgo::Yee }, -#ifndef WARPX_RZ // Not available in RZ +#ifndef WARPX_DIM_RZ // Not available in RZ {"ckc", MaxwellSolverAlgo::CKC }, #endif {"default", MaxwellSolverAlgo::Yee } diff --git a/Source/WarpX.cpp b/Source/WarpX.cpp index 18252cb64..d45dd3a71 100644 --- a/Source/WarpX.cpp +++ b/Source/WarpX.cpp @@ -993,7 +993,7 @@ WarpX::ComputeDivB (MultiFab& divB, int dcomp, { Real dxinv = 1./dx[0], dyinv = 1./dx[1], dzinv = 1./dx[2]; -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ const Real rmin = GetInstance().Geom(0).ProbLo(0); #endif @@ -1012,7 +1012,7 @@ WarpX::ComputeDivB (MultiFab& divB, int dcomp, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept { warpx_computedivb(i, j, k, dcomp, divBfab, Bxfab, Byfab, Bzfab, dxinv, dyinv, dzinv -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ ,rmin #endif ); @@ -1027,7 +1027,7 @@ WarpX::ComputeDivB (MultiFab& divB, int dcomp, { Real dxinv = 1./dx[0], dyinv = 1./dx[1], dzinv = 1./dx[2]; -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ const Real rmin = GetInstance().Geom(0).ProbLo(0); #endif @@ -1046,7 +1046,7 @@ WarpX::ComputeDivB (MultiFab& divB, int dcomp, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept { warpx_computedivb(i, j, k, dcomp, divBfab, Bxfab, Byfab, Bzfab, dxinv, dyinv, dzinv -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ ,rmin #endif ); @@ -1061,7 +1061,7 @@ WarpX::ComputeDivE (MultiFab& divE, int dcomp, { Real dxinv = 1./dx[0], dyinv = 1./dx[1], dzinv = 1./dx[2]; -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ const Real rmin = GetInstance().Geom(0).ProbLo(0); #endif @@ -1080,7 +1080,7 @@ WarpX::ComputeDivE (MultiFab& divE, int dcomp, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept { warpx_computedive(i, j, k, dcomp, divEfab, Exfab, Eyfab, Ezfab, dxinv, dyinv, dzinv -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ ,rmin #endif ); @@ -1095,7 +1095,7 @@ WarpX::ComputeDivE (MultiFab& divE, int dcomp, { Real dxinv = 1./dx[0], dyinv = 1./dx[1], dzinv = 1./dx[2]; -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ const Real rmin = GetInstance().Geom(0).ProbLo(0); #endif @@ -1114,7 +1114,7 @@ WarpX::ComputeDivE (MultiFab& divE, int dcomp, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept { warpx_computedive(i, j, k, dcomp, divEfab, Exfab, Eyfab, Ezfab, dxinv, dyinv, dzinv -#ifdef WARPX_RZ +#ifdef WARPX_DIM_RZ ,rmin #endif ); -- cgit v1.2.3 From 5cc13f04eaf1a41c22245129bf083f38e7c596ef Mon Sep 17 00:00:00 2001 From: Dave Grote Date: Thu, 8 Aug 2019 15:48:08 -0700 Subject: Converted DepositCharge to c++ (and removed duplication of code) --- Source/Laser/LaserParticleContainer.cpp | 14 +- Source/Particles/Deposition/Make.package | 1 + Source/Particles/PhysicalParticleContainer.cpp | 14 +- Source/Particles/WarpXParticleContainer.H | 10 +- Source/Particles/WarpXParticleContainer.cpp | 243 +++++++------------------ 5 files changed, 100 insertions(+), 182 deletions(-) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Laser/LaserParticleContainer.cpp b/Source/Laser/LaserParticleContainer.cpp index 515aa1f5d..786ebc622 100644 --- a/Source/Laser/LaserParticleContainer.cpp +++ b/Source/Laser/LaserParticleContainer.cpp @@ -453,7 +453,12 @@ LaserParticleContainer::Evolve (int lev, pti.GetPosition(m_xp[thread_num], m_yp[thread_num], m_zp[thread_num]); BL_PROFILE_VAR_STOP(blp_copy); - if (rho) DepositCharge(pti, wp, rho, crho, 0, np_current, np, thread_num, lev); + if (rho) { + DepositCharge(pti, wp, rho, 0, 0, np_current, thread_num, lev, lev); + if (crho) { + DepositCharge(pti, wp, crho, 0, np_current, np-np_current, thread_num, lev, lev-1); + } + } // // Particle Push @@ -522,7 +527,12 @@ LaserParticleContainer::Evolve (int lev, pti.SetPosition(m_xp[thread_num], m_yp[thread_num], m_zp[thread_num]); BL_PROFILE_VAR_STOP(blp_copy); - if (rho) DepositCharge(pti, wp, rho, crho, 1, np_current, np, thread_num, lev); + if (rho) { + DepositCharge(pti, wp, rho, 1, 0, np_current, thread_num, lev, lev); + if (crho) { + DepositCharge(pti, wp, crho, 1, np_current, np-np_current, thread_num, lev, lev-1); + } + } if (cost) { const Box& tbx = pti.tilebox(); diff --git a/Source/Particles/Deposition/Make.package b/Source/Particles/Deposition/Make.package index 0d5ebe2a7..e1aace998 100644 --- a/Source/Particles/Deposition/Make.package +++ b/Source/Particles/Deposition/Make.package @@ -1,3 +1,4 @@ CEXE_headers += CurrentDeposition.H +CEXE_headers += ChargeDeposition.H INCLUDE_LOCATIONS += $(WARPX_HOME)/Source/Particles/Deposition VPATH_LOCATIONS += $(WARPX_HOME)/Source/Particles/Deposition diff --git a/Source/Particles/PhysicalParticleContainer.cpp b/Source/Particles/PhysicalParticleContainer.cpp index 59d3ce76f..d10390204 100644 --- a/Source/Particles/PhysicalParticleContainer.cpp +++ b/Source/Particles/PhysicalParticleContainer.cpp @@ -1156,7 +1156,12 @@ PhysicalParticleContainer::Evolve (int lev, pti.GetPosition(m_xp[thread_num], m_yp[thread_num], m_zp[thread_num]); BL_PROFILE_VAR_STOP(blp_copy); - if (rho) DepositCharge(pti, wp, rho, crho, 0, np_current, np, thread_num, lev); + if (rho) { + DepositCharge(pti, wp, rho, 0, 0, np_current, thread_num, lev, lev); + if (has_buffer){ + DepositCharge(pti, wp, crho, 0, np_current, np-np_current, thread_num, lev, lev-1); + } + } if (! do_not_push) { @@ -1292,7 +1297,12 @@ PhysicalParticleContainer::Evolve (int lev, BL_PROFILE_VAR_STOP(blp_copy); } - if (rho) DepositCharge(pti, wp, rho, crho, 1, np_current, np, thread_num, lev); + if (rho) { + DepositCharge(pti, wp, rho, 1, 0, np_current, thread_num, lev, lev); + if (has_buffer){ + DepositCharge(pti, wp, crho, 1, np_current, np-np_current, thread_num, lev, lev-1); + } + } if (cost) { const Box& tbx = pti.tilebox(); diff --git a/Source/Particles/WarpXParticleContainer.H b/Source/Particles/WarpXParticleContainer.H index a609b4cb3..ac5b47ada 100644 --- a/Source/Particles/WarpXParticleContainer.H +++ b/Source/Particles/WarpXParticleContainer.H @@ -167,13 +167,13 @@ public: virtual void DepositCharge(WarpXParIter& pti, RealVector& wp, - amrex::MultiFab* rhomf, - amrex::MultiFab* crhomf, + amrex::MultiFab* rho, int icomp, - const long np_current, - const long np, + const long offset, + const long np_to_depose, int thread_num, - int lev ); + int lev, + int depos_lev); virtual void DepositCurrent(WarpXParIter& pti, RealVector& wp, diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index c5c6afc19..ec0f78564 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -12,6 +12,7 @@ #include #include #include +#include using namespace amrex; @@ -552,140 +553,91 @@ WarpXParticleContainer::DepositCurrent(WarpXParIter& pti, } void -WarpXParticleContainer::DepositCharge ( WarpXParIter& pti, RealVector& wp, - MultiFab* rhomf, MultiFab* crhomf, int icomp, - const long np_current, - const long np, int thread_num, int lev ) +WarpXParticleContainer::DepositCharge (WarpXParIter& pti, RealVector& wp, + MultiFab* rho, int icomp, + 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"); - BL_PROFILE_VAR_NS("PICSAR::ChargeDeposition", blp_pxr_chd); - BL_PROFILE_VAR_NS("PPC::Evolve::Accumulate", blp_accumulate); - - const std::array& xyzmin_tile = WarpX::LowerCorner(pti.tilebox(), lev); - const long lvect = 8; + // If no particles, do not do anything + if (np_to_depose == 0) return; - long ngRho = rhomf->nGrow(); - Real* data_ptr; - Box tile_box = convert(pti.tilebox(), IntVect::TheUnitVector()); + const long ngRho = rho->nGrow(); + const std::array& dx = WarpX::CellSize(std::max(depos_lev,0)); + Real q = this->charge; - const std::array& dx = WarpX::CellSize(lev); - const std::array& cdx = WarpX::CellSize(std::max(lev-1,0)); + BL_PROFILE_VAR_NS("PPC::ChargeDeposition", blp_ppc_chd); + BL_PROFILE_VAR_NS("PPC::Evolve::Accumulate", blp_accumulate); - // Deposit charge for particles that are not in the current buffers - if (np_current > 0) - { - const std::array& xyzmin = xyzmin_tile; + // Get tile box where charge is deposited. + // The tile box is different when depositing in the buffers (depos_lev const& rho_arr = rho->array(pti); #else - tile_box.grow(ngRho); - local_rho[thread_num].resize(tile_box); + // Tiling is on: rho_ptr points to local_rho[thread_num] + Box tb = amrex::convert(tilebox, IntVect::TheUnitVector()); - data_ptr = local_rho[thread_num].dataPtr(); - auto rholen = local_rho[thread_num].length(); + local_rho[thread_num].resize(tb); - local_rho[thread_num].setVal(0.0); -#endif + // local_rho[thread_num] is set to zero + local_rho[thread_num].setVal(0.0); -#if (AMREX_SPACEDIM == 3) - const long nx = rholen[0]-1-2*ngRho; - const long ny = rholen[1]-1-2*ngRho; - const long nz = rholen[2]-1-2*ngRho; -#else - const long nx = rholen[0]-1-2*ngRho; - const long ny = 0; - const long nz = rholen[1]-1-2*ngRho; -#endif - BL_PROFILE_VAR_START(blp_pxr_chd); - warpx_charge_deposition(data_ptr, &np_current, - m_xp[thread_num].dataPtr(), - m_yp[thread_num].dataPtr(), - m_zp[thread_num].dataPtr(), - wp.dataPtr(), - &this->charge, - &xyzmin[0], &xyzmin[1], &xyzmin[2], - &dx[0], &dx[1], &dx[2], &nx, &ny, &nz, - &ngRho, &ngRho, &ngRho, - &WarpX::nox,&WarpX::noy,&WarpX::noz, - &lvect, &WarpX::charge_deposition_algo); -#ifdef WARPX_DIM_RZ - warpx_charge_deposition_rz_volume_scaling( - data_ptr, &ngRho, rholen.getVect(), - &xyzmin[0], &dx[0]); + Array4 const& rho_arr = local_rho[thread_num].array(); #endif - BL_PROFILE_VAR_STOP(blp_pxr_chd); - -#ifndef AMREX_USE_GPU - BL_PROFILE_VAR_START(blp_accumulate); - - (*rhomf)[pti].atomicAdd(local_rho[thread_num], tile_box, tile_box, 0, icomp, 1); - - BL_PROFILE_VAR_STOP(blp_accumulate); -#endif - } - - // Deposit charge for particles that are in the current buffers - if (np_current < np) - { - const IntVect& ref_ratio = WarpX::RefRatio(lev-1); - const Box& ctilebox = amrex::coarsen(pti.tilebox(), ref_ratio); - const std::array& cxyzmin_tile = WarpX::LowerCorner(ctilebox, lev-1); - -#ifdef AMREX_USE_GPU - data_ptr = (*crhomf)[pti].dataPtr(icomp); - auto rholen = (*crhomf)[pti].length(); -#else - tile_box = amrex::convert(ctilebox, IntVect::TheUnitVector()); - tile_box.grow(ngRho); - local_rho[thread_num].resize(tile_box); + // GPU, no tiling: deposit directly in rho + // CPU, tiling: deposit into local_rho - data_ptr = local_rho[thread_num].dataPtr(); - auto rholen = local_rho[thread_num].length(); - - local_rho[thread_num].setVal(0.0); -#endif + Real* AMREX_RESTRICT xp = m_xp[thread_num].dataPtr() + offset; + Real* AMREX_RESTRICT zp = m_zp[thread_num].dataPtr() + offset; + Real* AMREX_RESTRICT yp = m_yp[thread_num].dataPtr() + offset; -#if (AMREX_SPACEDIM == 3) - const long nx = rholen[0]-1-2*ngRho; - const long ny = rholen[1]-1-2*ngRho; - const long nz = rholen[2]-1-2*ngRho; -#else - const long nx = rholen[0]-1-2*ngRho; - const long ny = 0; - const long nz = rholen[1]-1-2*ngRho; -#endif + // Lower corner of tile box physical domain + // Note that this includes guard cells since it is after tilebox.ngrow + const std::array& xyzmin = WarpX::LowerCorner(tilebox, depos_lev); + // Indices of the lower bound + const Dim3 lo = lbound(tilebox); - long ncrse = np - np_current; - BL_PROFILE_VAR_START(blp_pxr_chd); - warpx_charge_deposition(data_ptr, &ncrse, - m_xp[thread_num].dataPtr() + np_current, - m_yp[thread_num].dataPtr() + np_current, - m_zp[thread_num].dataPtr() + np_current, - wp.dataPtr() + np_current, - &this->charge, - &cxyzmin_tile[0], &cxyzmin_tile[1], &cxyzmin_tile[2], - &cdx[0], &cdx[1], &cdx[2], &nx, &ny, &nz, - &ngRho, &ngRho, &ngRho, - &WarpX::nox,&WarpX::noy,&WarpX::noz, - &lvect, &WarpX::charge_deposition_algo); + BL_PROFILE_VAR_START(blp_ppc_chd); + if (WarpX::nox == 1){ + doChargeDepositionShapeN<1>(xp, yp, zp, wp.dataPtr()+offset, rho_arr, + np_to_depose, dx, xyzmin, lo, q); + } else if (WarpX::nox == 2){ + doChargeDepositionShapeN<2>(xp, yp, zp, wp.dataPtr()+offset, rho_arr, + np_to_depose, dx, xyzmin, lo, q); + } else if (WarpX::nox == 3){ + doChargeDepositionShapeN<3>(xp, yp, zp, wp.dataPtr()+offset, rho_arr, + np_to_depose, dx, xyzmin, lo, q); + } #ifdef WARPX_DIM_RZ - warpx_charge_deposition_rz_volume_scaling( - data_ptr, &ngRho, rholen.getVect(), - &cxyzmin_tile[0], &cdx[0]); + warpx_charge_deposition_rz_volume_scaling( + data_ptr, &ngRho, rholen.getVect(), + &xyzmin[0], &dx[0]); #endif - BL_PROFILE_VAR_STOP(blp_pxr_chd); + BL_PROFILE_VAR_STOP(blp_ppc_chd); #ifndef AMREX_USE_GPU - BL_PROFILE_VAR_START(blp_accumulate); + BL_PROFILE_VAR_START(blp_accumulate); - (*crhomf)[pti].atomicAdd(local_rho[thread_num], tile_box, tile_box, 0, icomp, 1); + (*rho)[pti].atomicAdd(local_rho[thread_num], tb, tb, 0, icomp, 1); - BL_PROFILE_VAR_STOP(blp_accumulate); + BL_PROFILE_VAR_STOP(blp_accumulate); #endif - } -}; +} void WarpXParticleContainer::DepositCharge (Vector >& rho, bool local) @@ -762,8 +714,6 @@ WarpXParticleContainer::GetChargeDensity (int lev, bool local) BoxArray nba = ba; nba.surroundingNodes(); - const std::array& dx = WarpX::CellSize(lev); - const int ng = WarpX::nox; auto rho = std::unique_ptr(new MultiFab(nba,dm,1,ng)); @@ -773,74 +723,21 @@ WarpXParticleContainer::GetChargeDensity (int lev, bool local) #pragma omp parallel { #endif - Cuda::ManagedDeviceVector xp, yp, zp; #ifdef _OPENMP - FArrayBox rho_loc; + int thread_num = omp_get_thread_num(); +#else + int thread_num = 0; #endif for (WarpXParIter pti(*this, lev); pti.isValid(); ++pti) { + const long np = pti.numParticles(); auto& wp = pti.GetAttribs(PIdx::w); - const long np = pti.numParticles(); + pti.GetPosition(m_xp[thread_num], m_yp[thread_num], m_zp[thread_num]); - pti.GetPosition(xp, yp, zp); - - // Data on the grid - Real* data_ptr; - FArrayBox& rhofab = (*rho)[pti]; -#ifdef _OPENMP - const std::array& xyzmin_tile = WarpX::LowerCorner(pti.tilebox(), lev); - Box tile_box = convert(pti.tilebox(), IntVect::TheUnitVector()); - const std::array& xyzmin = xyzmin_tile; - tile_box.grow(ng); - rho_loc.resize(tile_box); - rho_loc = 0.0; - data_ptr = rho_loc.dataPtr(); - auto rholen = rho_loc.length(); -#else - const Box& box = pti.validbox(); - const std::array& xyzmin_grid = WarpX::LowerCorner(box, lev); - const std::array& xyzmin = xyzmin_grid; - data_ptr = rhofab.dataPtr(); - auto rholen = rhofab.length(); -#endif - -#if (AMREX_SPACEDIM == 3) - const long nx = rholen[0]-1-2*ng; - const long ny = rholen[1]-1-2*ng; - const long nz = rholen[2]-1-2*ng; -#else - const long nx = rholen[0]-1-2*ng; - const long ny = 0; - const long nz = rholen[1]-1-2*ng; -#endif - - long nxg = ng; - long nyg = ng; - long nzg = ng; - long lvect = 8; - - warpx_charge_deposition(data_ptr, - &np, - xp.dataPtr(), - yp.dataPtr(), - zp.dataPtr(), wp.dataPtr(), - &this->charge, &xyzmin[0], &xyzmin[1], &xyzmin[2], - &dx[0], &dx[1], &dx[2], &nx, &ny, &nz, - &nxg, &nyg, &nzg, &WarpX::nox,&WarpX::noy,&WarpX::noz, - &lvect, &WarpX::charge_deposition_algo); -#ifdef WARPX_DIM_RZ - long ngRho = WarpX::nox; - warpx_charge_deposition_rz_volume_scaling( - data_ptr, &ngRho, rholen.getVect(), - &xyzmin[0], &dx[0]); -#endif - -#ifdef _OPENMP - rhofab.atomicAdd(rho_loc); + DepositCharge(pti, wp, rho.get(), 0, 0, np, thread_num, lev, lev); } -#endif } if (!local) rho->SumBoundary(gm.periodicity()); -- cgit v1.2.3 From 9fa3a16002c585720e1741b5fdbe5e6d5b974af9 Mon Sep 17 00:00:00 2001 From: Dave Grote Date: Thu, 8 Aug 2019 16:27:54 -0700 Subject: Implemented ApplyInverseVolumeScalingToChargeDensity --- Source/Evolve/WarpXEvolveEM.cpp | 10 +++++- Source/FieldSolver/WarpXPushFieldsEM.cpp | 49 +++++++++++++++++++++++++++++ Source/Particles/WarpXParticleContainer.cpp | 5 --- Source/WarpX.H | 4 +++ 4 files changed, 62 insertions(+), 6 deletions(-) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Evolve/WarpXEvolveEM.cpp b/Source/Evolve/WarpXEvolveEM.cpp index a097d8814..e05fde03e 100644 --- a/Source/Evolve/WarpXEvolveEM.cpp +++ b/Source/Evolve/WarpXEvolveEM.cpp @@ -483,11 +483,19 @@ WarpX::PushParticlesandDepose (int lev, Real cur_time) Bfield_cax[lev][0].get(), Bfield_cax[lev][1].get(), Bfield_cax[lev][2].get(), cur_time, dt[lev]); #ifdef WARPX_DIM_RZ - // This is called after all particles have deposited their current. + // This is called after all particles have deposited their current and charge. ApplyInverseVolumeScalingToCurrentDensity(current_fp[lev][0].get(), current_fp[lev][1].get(), current_fp[lev][2].get(), lev); if (current_buf[lev][0].get()) { ApplyInverseVolumeScalingToCurrentDensity(current_buf[lev][0].get(), current_buf[lev][1].get(), current_buf[lev][2].get(), lev-1); } + if (rho_fp[lev].get()) { + for (int icomp = 0 ; icomp < rho_fp[lev].get()->nComp() ; icomp++) { + ApplyInverseVolumeScalingToChargeDensity(rho_fp[lev].get(), icomp, lev); + if (charge_buf[lev].get()) { + ApplyInverseVolumeScalingToChargeDensity(charge_buf[lev].get(), icomp, lev-1); + } + } + } #endif } diff --git a/Source/FieldSolver/WarpXPushFieldsEM.cpp b/Source/FieldSolver/WarpXPushFieldsEM.cpp index 21000424d..f199e0660 100644 --- a/Source/FieldSolver/WarpXPushFieldsEM.cpp +++ b/Source/FieldSolver/WarpXPushFieldsEM.cpp @@ -671,4 +671,53 @@ WarpX::ApplyInverseVolumeScalingToCurrentDensity (MultiFab* Jx, MultiFab* Jy, Mu }); } } + +void +WarpX::ApplyInverseVolumeScalingToChargeDensity (MultiFab* Rho, int icomp, int lev) +{ + const long ngRho = Rho->nGrow(); + const std::array& dx = WarpX::CellSize(lev); + const Real dr = dx[0]; + + Box tilebox; + + for ( MFIter mfi(*Rho, TilingIfNotGPU()); mfi.isValid(); ++mfi ) + { + + Array4 const& Rho_arr = Rho->array(mfi); + + tilebox = mfi.tilebox(); + Box tb = convert(tilebox, IntVect::TheUnitVector()); + + // Lower corner of tile box physical domain + // Note that this is done before the tilebox.grow so that + // these do not include the guard cells. + const std::array& xyzmin = WarpX::LowerCorner(tilebox, lev); + const Dim3 lo = lbound(tilebox); + const Real rmin = xyzmin[0]; + const int irmin = lo.x; + + // Rescale charge in r-z mode since the inverse volume factor was not + // included in the charge deposition. + amrex::ParallelFor(tb, + [=] AMREX_GPU_DEVICE (int i, int j, int k) + { + // Wrap the charge density deposited in the guard cells around + // to the cells above the axis. + // Rho is located on the boundary + if (rmin == 0. && 0 < i && i <= ngRho) { + Rho_arr(i,j,0,icomp) += Rho_arr(-i,j,0,icomp); + } + + // Apply the inverse volume scaling + const amrex::Real r = std::abs(rmin + (i - irmin)*dr); + if (r == 0.) { + // Verboncoeur JCP 164, 421-427 (2001) : corrected volume on axis + Rho_arr(i,j,0,icomp) /= (MathConst::pi*dr/3.); + } else { + Rho_arr(i,j,0,icomp) /= (2.*MathConst::pi*r); + } + }); + } +} #endif diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index ec0f78564..32abae3ad 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -623,11 +623,6 @@ WarpXParticleContainer::DepositCharge (WarpXParIter& pti, RealVector& wp, doChargeDepositionShapeN<3>(xp, yp, zp, wp.dataPtr()+offset, rho_arr, np_to_depose, dx, xyzmin, lo, q); } -#ifdef WARPX_DIM_RZ - warpx_charge_deposition_rz_volume_scaling( - data_ptr, &ngRho, rholen.getVect(), - &xyzmin[0], &dx[0]); -#endif BL_PROFILE_VAR_STOP(blp_ppc_chd); #ifndef AMREX_USE_GPU diff --git a/Source/WarpX.H b/Source/WarpX.H index d8b1b9ce0..67d2366d4 100644 --- a/Source/WarpX.H +++ b/Source/WarpX.H @@ -183,6 +183,10 @@ public: amrex::MultiFab* Jy, amrex::MultiFab* Jz, int lev); + + void ApplyInverseVolumeScalingToChargeDensity(amrex::MultiFab* Rho, + int icomp, + int lev); #endif void DampPML (); -- cgit v1.2.3 From 7b72953b793af38045e2508302594c7fe83a42e9 Mon Sep 17 00:00:00 2001 From: MaxThevenet Date: Thu, 8 Aug 2019 16:51:22 -0700 Subject: fix current deposition in buffers --- Source/Particles/WarpXParticleContainer.cpp | 24 ++++++++++++------------ 1 file changed, 12 insertions(+), 12 deletions(-) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index c5c6afc19..8e4adc528 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -506,35 +506,35 @@ WarpXParticleContainer::DepositCurrent(WarpXParIter& pti, if (WarpX::current_deposition_algo == CurrentDepositionAlgo::Esirkepov) { if (WarpX::nox == 1){ - doEsirkepovDepositionShapeN<1>(xp, yp, zp, wp.dataPtr(), uxp.dataPtr(), - uyp.dataPtr(), uzp.dataPtr(), jx_arr, jy_arr, + doEsirkepovDepositionShapeN<1>(xp, yp, zp, wp.dataPtr() + offset, uxp.dataPtr() + offset, + uyp.dataPtr() + offset, uzp.dataPtr() + offset, jx_arr, jy_arr, jz_arr, np_to_depose, dt, dx, xyzmin, lo, q); } else if (WarpX::nox == 2){ - doEsirkepovDepositionShapeN<2>(xp, yp, zp, wp.dataPtr(), uxp.dataPtr(), - uyp.dataPtr(), uzp.dataPtr(), jx_arr, jy_arr, + doEsirkepovDepositionShapeN<2>(xp, yp, zp, wp.dataPtr() + offset, uxp.dataPtr() + offset, + uyp.dataPtr() + offset, uzp.dataPtr() + offset, jx_arr, jy_arr, jz_arr, np_to_depose, dt, dx, xyzmin, lo, q); } else if (WarpX::nox == 3){ - doEsirkepovDepositionShapeN<3>(xp, yp, zp, wp.dataPtr(), uxp.dataPtr(), - uyp.dataPtr(), uzp.dataPtr(), jx_arr, jy_arr, + doEsirkepovDepositionShapeN<3>(xp, yp, zp, wp.dataPtr() + offset, uxp.dataPtr() + offset, + uyp.dataPtr() + offset, uzp.dataPtr() + offset, jx_arr, jy_arr, jz_arr, np_to_depose, dt, dx, xyzmin, lo, q); } } else { if (WarpX::nox == 1){ - doDepositionShapeN<1>(xp, yp, zp, wp.dataPtr(), uxp.dataPtr(), - uyp.dataPtr(), uzp.dataPtr(), jx_arr, jy_arr, + doDepositionShapeN<1>(xp, yp, zp, wp.dataPtr() + offset, uxp.dataPtr() + offset, + uyp.dataPtr() + offset, uzp.dataPtr() + offset, jx_arr, jy_arr, jz_arr, np_to_depose, dt, dx, xyzmin, lo, stagger_shift, q); } else if (WarpX::nox == 2){ - doDepositionShapeN<2>(xp, yp, zp, wp.dataPtr(), uxp.dataPtr(), - uyp.dataPtr(), uzp.dataPtr(), jx_arr, jy_arr, + doDepositionShapeN<2>(xp, yp, zp, wp.dataPtr() + offset, uxp.dataPtr() + offset, + uyp.dataPtr() + offset, uzp.dataPtr() + offset, jx_arr, jy_arr, jz_arr, np_to_depose, dt, dx, xyzmin, lo, stagger_shift, q); } else if (WarpX::nox == 3){ - doDepositionShapeN<3>(xp, yp, zp, wp.dataPtr(), uxp.dataPtr(), - uyp.dataPtr(), uzp.dataPtr(), jx_arr, jy_arr, + doDepositionShapeN<3>(xp, yp, zp, wp.dataPtr() + offset, uxp.dataPtr() + offset, + uyp.dataPtr() + offset, uzp.dataPtr() + offset, jx_arr, jy_arr, jz_arr, np_to_depose, dt, dx, xyzmin, lo, stagger_shift, q); } -- cgit v1.2.3 From c45bde5edf520cd808a620e243e426ace445441d Mon Sep 17 00:00:00 2001 From: Dave Grote Date: Thu, 8 Aug 2019 17:07:04 -0700 Subject: Minor clean up in charge deposition conversion --- Source/Particles/WarpXParticleContainer.cpp | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index 32abae3ad..4ccddaedf 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -567,7 +567,7 @@ WarpXParticleContainer::DepositCharge (WarpXParIter& pti, RealVector& wp, const long ngRho = rho->nGrow(); const std::array& dx = WarpX::CellSize(std::max(depos_lev,0)); - Real q = this->charge; + const Real q = this->charge; BL_PROFILE_VAR_NS("PPC::ChargeDeposition", blp_ppc_chd); BL_PROFILE_VAR_NS("PPC::Evolve::Accumulate", blp_accumulate); @@ -586,11 +586,11 @@ WarpXParticleContainer::DepositCharge (WarpXParIter& pti, RealVector& wp, tilebox.grow(ngRho); #ifdef AMREX_USE_GPU - // No tiling on GPU: rho_ptr points to the full rho array. + // No tiling on GPU: rho_arr points to the full rho array. Array4 const& rho_arr = rho->array(pti); #else - // Tiling is on: rho_ptr points to local_rho[thread_num] - Box tb = amrex::convert(tilebox, IntVect::TheUnitVector()); + // Tiling is on: rho_arr points to local_rho[thread_num] + const Box tb = amrex::convert(tilebox, IntVect::TheUnitVector()); local_rho[thread_num].resize(tb); -- cgit v1.2.3 From 00cd6d4cee3c5f01b5e68b080019fcd85dc73714 Mon Sep 17 00:00:00 2001 From: Dave Grote Date: Thu, 8 Aug 2019 17:57:58 -0700 Subject: In GetChargeDensity, added call to ApplyInverseVolumeScalingToChargeDensity --- Source/Particles/WarpXParticleContainer.cpp | 4 ++++ 1 file changed, 4 insertions(+) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index 4ccddaedf..c4e670bf8 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -735,6 +735,10 @@ WarpXParticleContainer::GetChargeDensity (int lev, bool local) } } +#ifdef WARPX_DIM_RZ + WarpX::GetInstance().ApplyInverseVolumeScalingToChargeDensity(rho.get(), lev); +#endif + if (!local) rho->SumBoundary(gm.periodicity()); return rho; -- cgit v1.2.3 From 8b455ff9c46aa8d8a378d8c262eea8edcdf10e3f Mon Sep 17 00:00:00 2001 From: Dave Grote Date: Fri, 9 Aug 2019 16:37:29 -0700 Subject: Fix charge deposition for GPU --- Source/Particles/WarpXParticleContainer.cpp | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index c4e670bf8..aa74f9c4a 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -587,7 +587,8 @@ WarpXParticleContainer::DepositCharge (WarpXParIter& pti, RealVector& wp, #ifdef AMREX_USE_GPU // No tiling on GPU: rho_arr points to the full rho array. - Array4 const& rho_arr = rho->array(pti); + MultiFab rhoi(*rho, amrex::make_alias, icomp, 1); + Array4 const& rho_arr = rhoi.array(pti); #else // Tiling is on: rho_arr points to local_rho[thread_num] const Box tb = amrex::convert(tilebox, IntVect::TheUnitVector()); -- cgit v1.2.3 From a38e581301c7be825a6b7ca0d742fd714aef882d Mon Sep 17 00:00:00 2001 From: Remi Lehe Date: Fri, 9 Aug 2019 20:21:33 -0700 Subject: Correct syntax error in the absence of OPENMP --- Source/Particles/WarpXParticleContainer.cpp | 2 ++ 1 file changed, 2 insertions(+) (limited to 'Source/Particles/WarpXParticleContainer.cpp') diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index aa74f9c4a..3bbb65d1b 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -734,7 +734,9 @@ WarpXParticleContainer::GetChargeDensity (int lev, bool local) DepositCharge(pti, wp, rho.get(), 0, 0, np, thread_num, lev, lev); } +#ifdef _OPENMP } +#endif #ifdef WARPX_DIM_RZ WarpX::GetInstance().ApplyInverseVolumeScalingToChargeDensity(rho.get(), lev); -- cgit v1.2.3