aboutsummaryrefslogtreecommitdiff
path: root/Source/Particles
diff options
context:
space:
mode:
Diffstat (limited to 'Source/Particles')
-rw-r--r--Source/Particles/MultiParticleContainer.cpp2
-rw-r--r--Source/Particles/PhysicalParticleContainer.cpp81
-rw-r--r--Source/Particles/RigidInjectedParticleContainer.cpp1
-rw-r--r--Source/Particles/WarpXParticleContainer.H3
-rw-r--r--Source/Particles/WarpXParticleContainer.cpp30
5 files changed, 75 insertions, 42 deletions
diff --git a/Source/Particles/MultiParticleContainer.cpp b/Source/Particles/MultiParticleContainer.cpp
index 9d39ec2f9..6475d1463 100644
--- a/Source/Particles/MultiParticleContainer.cpp
+++ b/Source/Particles/MultiParticleContainer.cpp
@@ -297,7 +297,7 @@ MultiParticleContainer::GetChargeDensity (int lev, bool local)
std::unique_ptr<MultiFab> rho = allcontainers[0]->GetChargeDensity(lev, true);
for (unsigned i = 1, n = allcontainers.size(); i < n; ++i) {
std::unique_ptr<MultiFab> rhoi = allcontainers[i]->GetChargeDensity(lev, true);
- MultiFab::Add(*rho, *rhoi, 0, 0, 1, rho->nGrow());
+ MultiFab::Add(*rho, *rhoi, 0, 0, rho->nComp(), rho->nGrow());
}
if (!local) {
const Geometry& gm = allcontainers[0]->Geom(lev);
diff --git a/Source/Particles/PhysicalParticleContainer.cpp b/Source/Particles/PhysicalParticleContainer.cpp
index 43b46ec49..93a0ad7ea 100644
--- a/Source/Particles/PhysicalParticleContainer.cpp
+++ b/Source/Particles/PhysicalParticleContainer.cpp
@@ -41,16 +41,19 @@ NumParticlesToAdd(const Box& overlap_box, const RealBox& overlap_realbox,
int ref_num_ppc = num_ppc * AMREX_D_TERM(fac, *fac, *fac);
for (int i_part=0; i_part<ref_num_ppc;i_part++) {
- std::array<Real, 3> r;
- plasma_injector->getPositionUnitBox(r, i_part, fac);
+ std::array<Real, 3> point;
+ plasma_injector->getPositionUnitBox(point, i_part, fac);
+ Real x = overlap_corner[0] + (iv[0] + point[0])*dx[0];
#if ( AMREX_SPACEDIM == 3 )
- Real x = overlap_corner[0] + (iv[0] + r[0])*dx[0];
- Real y = overlap_corner[1] + (iv[1] + r[1])*dx[1];
- Real z = overlap_corner[2] + (iv[2] + r[2])*dx[2];
+ Real y = overlap_corner[1] + (iv[1] + point[1])*dx[1];
+ Real z = overlap_corner[2] + (iv[2] + point[2])*dx[2];
#elif ( AMREX_SPACEDIM == 2 )
- Real x = overlap_corner[0] + (iv[0] + r[0])*dx[0];
Real y = 0;
- Real z = overlap_corner[1] + (iv[1] + r[1])*dx[1];
+#ifdef WARPX_RZ
+ Real z = overlap_corner[1] + (iv[1] + point[2])*dx[1];
+#else
+ Real z = overlap_corner[1] + (iv[1] + point[1])*dx[1];
+#endif
#endif
// If the new particle is not inside the tile box,
// go to the next generated particle.
@@ -448,16 +451,19 @@ PhysicalParticleContainer::AddPlasmaCPU (int lev, RealBox part_realbox)
int ref_num_ppc = num_ppc * AMREX_D_TERM(fac, *fac, *fac);
for (int i_part=0; i_part<ref_num_ppc;i_part++) {
- std::array<Real, 3> r;
- plasma_injector->getPositionUnitBox(r, i_part, fac);
+ std::array<Real, 3> point;
+ plasma_injector->getPositionUnitBox(point, i_part, fac);
+ Real x = overlap_corner[0] + (iv[0] + point[0])*dx[0];
#if ( AMREX_SPACEDIM == 3 )
- Real x = overlap_corner[0] + (iv[0] + r[0])*dx[0];
- Real y = overlap_corner[1] + (iv[1] + r[1])*dx[1];
- Real z = overlap_corner[2] + (iv[2] + r[2])*dx[2];
+ Real y = overlap_corner[1] + (iv[1] + point[1])*dx[1];
+ Real z = overlap_corner[2] + (iv[2] + point[2])*dx[2];
#elif ( AMREX_SPACEDIM == 2 )
- Real x = overlap_corner[0] + (iv[0] + r[0])*dx[0];
Real y = 0;
- Real z = overlap_corner[1] + (iv[1] + r[1])*dx[1];
+#ifdef WARPX_RZ
+ Real z = overlap_corner[1] + (iv[1] + point[2])*dx[1];
+#else
+ Real z = overlap_corner[1] + (iv[1] + point[1])*dx[1];
+#endif
#endif
// If the new particle is not inside the tile box,
// go to the next generated particle.
@@ -473,11 +479,18 @@ PhysicalParticleContainer::AddPlasmaCPU (int lev, RealBox part_realbox)
Real yb = y;
#ifdef WARPX_RZ
- // Replace the x and y, choosing the angle randomly.
+ // Replace the x and y, setting an angle theta.
// These x and y are used to get the momentum and density
- Real theta = 2.*MathConst::pi*amrex::Random();
- y = x*std::sin(theta);
- x = x*std::cos(theta);
+ Real theta;
+ if (WarpX::nmodes == 1) {
+ // With only 1 mode, the angle doesn't matter so
+ // choose it randomly.
+ theta = 2.*MathConst::pi*amrex::Random();
+ } else {
+ theta = 2.*MathConst::pi*point[1];
+ }
+ x = xb*std::cos(theta);
+ y = xb*std::sin(theta);
#endif
Real dens;
@@ -690,16 +703,19 @@ PhysicalParticleContainer::AddPlasmaGPU (int lev, RealBox part_realbox)
int ref_num_ppc = num_ppc * AMREX_D_TERM(fac, *fac, *fac);
for (int i_part=0; i_part<ref_num_ppc;i_part++) {
- std::array<Real, 3> r;
- plasma_injector->getPositionUnitBox(r, i_part, fac);
+ std::array<Real, 3> point;
+ plasma_injector->getPositionUnitBox(point, i_part, fac);
+ Real x = overlap_corner[0] + (iv[0] + point[0])*dx[0];
#if ( AMREX_SPACEDIM == 3 )
- Real x = overlap_corner[0] + (iv[0] + r[0])*dx[0];
- Real y = overlap_corner[1] + (iv[1] + r[1])*dx[1];
- Real z = overlap_corner[2] + (iv[2] + r[2])*dx[2];
+ Real y = overlap_corner[1] + (iv[1] + point[1])*dx[1];
+ Real z = overlap_corner[2] + (iv[2] + point[2])*dx[2];
#elif ( AMREX_SPACEDIM == 2 )
- Real x = overlap_corner[0] + (iv[0] + r[0])*dx[0];
Real y = 0;
- Real z = overlap_corner[1] + (iv[1] + r[1])*dx[1];
+#ifdef WARPX_RZ
+ Real z = overlap_corner[1] + (iv[1] + point[2])*dx[1];
+#else
+ Real z = overlap_corner[1] + (iv[1] + point[1])*dx[1];
+#endif
#endif
// If the new particle is not inside the tile box,
// go to the next generated particle.
@@ -715,9 +731,16 @@ PhysicalParticleContainer::AddPlasmaGPU (int lev, RealBox part_realbox)
Real yb = y;
#ifdef WARPX_RZ
- // Replace the x and y, choosing the angle randomly.
+ // Replace the x and y, setting an angle theta.
// These x and y are used to get the momentum and density
- Real theta = 2.*MathConst::pi*amrex::Random();
+ Real theta;
+ if (WarpX::nmodes == 1) {
+ // With only 1 mode, the angle doesn't matter so
+ // choose it randomly.
+ theta = 2.*MathConst::pi*amrex::Random();
+ } else {
+ theta = 2.*MathConst::pi*point[1];
+ }
x = xb*std::cos(theta);
y = xb*std::sin(theta);
#endif
@@ -1133,6 +1156,7 @@ PhysicalParticleContainer::FieldGather (int lev,
BL_TO_FORTRAN_ANYD(bxfab),
BL_TO_FORTRAN_ANYD(byfab),
BL_TO_FORTRAN_ANYD(bzfab),
+ &WarpX::nmodes,
&ll4symtry, &WarpX::l_lower_order_in_v, &WarpX::do_nodal,
&lvect_fieldgathe, &WarpX::field_gathering_algo);
@@ -1425,6 +1449,7 @@ PhysicalParticleContainer::Evolve (int lev,
BL_TO_FORTRAN_ANYD(*bxfab),
BL_TO_FORTRAN_ANYD(*byfab),
BL_TO_FORTRAN_ANYD(*bzfab),
+ &WarpX::nmodes,
&ll4symtry, &WarpX::l_lower_order_in_v, &WarpX::do_nodal,
&lvect_fieldgathe, &WarpX::field_gathering_algo);
@@ -1512,6 +1537,7 @@ PhysicalParticleContainer::Evolve (int lev,
BL_TO_FORTRAN_ANYD(*cbxfab),
BL_TO_FORTRAN_ANYD(*cbyfab),
BL_TO_FORTRAN_ANYD(*cbzfab),
+ &WarpX::nmodes,
&ll4symtry, &WarpX::l_lower_order_in_v, &WarpX::do_nodal,
&lvect_fieldgathe, &WarpX::field_gathering_algo);
}
@@ -1853,6 +1879,7 @@ PhysicalParticleContainer::PushP (int lev, Real dt,
BL_TO_FORTRAN_ANYD(bxfab),
BL_TO_FORTRAN_ANYD(byfab),
BL_TO_FORTRAN_ANYD(bzfab),
+ &WarpX::nmodes,
&ll4symtry, &WarpX::l_lower_order_in_v, &WarpX::do_nodal,
&lvect_fieldgathe, &WarpX::field_gathering_algo);
diff --git a/Source/Particles/RigidInjectedParticleContainer.cpp b/Source/Particles/RigidInjectedParticleContainer.cpp
index 2a3e8dd0d..d659b3854 100644
--- a/Source/Particles/RigidInjectedParticleContainer.cpp
+++ b/Source/Particles/RigidInjectedParticleContainer.cpp
@@ -426,6 +426,7 @@ RigidInjectedParticleContainer::PushP (int lev, Real dt,
BL_TO_FORTRAN_ANYD(bxfab),
BL_TO_FORTRAN_ANYD(byfab),
BL_TO_FORTRAN_ANYD(bzfab),
+ &WarpX::nmodes,
&ll4symtry, &l_lower_order_in_v, &WarpX::do_nodal,
&lvect_fieldgathe, &WarpX::field_gathering_algo);
diff --git a/Source/Particles/WarpXParticleContainer.H b/Source/Particles/WarpXParticleContainer.H
index 0e800bf1d..b2d86c587 100644
--- a/Source/Particles/WarpXParticleContainer.H
+++ b/Source/Particles/WarpXParticleContainer.H
@@ -42,6 +42,9 @@ namespace ParticleStringNames
{"Bx", PIdx::Bx },
{"By", PIdx::By },
{"Bz", PIdx::Bz }
+#ifdef WARPX_RZ
+ ,{"theta", PIdx::theta}
+#endif
};
}
diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp
index 317f46fd4..a1517cb44 100644
--- a/Source/Particles/WarpXParticleContainer.cpp
+++ b/Source/Particles/WarpXParticleContainer.cpp
@@ -354,9 +354,9 @@ WarpXParticleContainer::DepositCurrent(WarpXParIter& pti,
tby.grow(ngJ);
tbz.grow(ngJ);
- local_jx[thread_num].resize(tbx);
- local_jy[thread_num].resize(tby);
- local_jz[thread_num].resize(tbz);
+ local_jx[thread_num].resize(tbx, jx->nComp());
+ local_jy[thread_num].resize(tby, jy->nComp());
+ local_jz[thread_num].resize(tbz, jz->nComp());
Real* jx_ptr = local_jx[thread_num].dataPtr();
Real* jy_ptr = local_jy[thread_num].dataPtr();
@@ -379,6 +379,7 @@ WarpXParticleContainer::DepositCurrent(WarpXParIter& pti,
jx_ptr, &ngJ, jxntot.getVect(),
jy_ptr, &ngJ, jyntot.getVect(),
jz_ptr, &ngJ, jzntot.getVect(),
+ &WarpX::nmodes,
&np_to_depose,
m_xp[thread_num].dataPtr() + offset,
m_yp[thread_num].dataPtr() + offset,
@@ -399,6 +400,7 @@ WarpXParticleContainer::DepositCurrent(WarpXParIter& pti,
jx_ptr, &ngJ, jxntot.getVect(),
jy_ptr, &ngJ, jyntot.getVect(),
jz_ptr, &ngJ, jzntot.getVect(),
+ &WarpX::nmodes,
&xyzmin[0], &dx[0]);
#endif
BL_PROFILE_VAR_STOP(blp_pxr_cd);
@@ -407,9 +409,9 @@ WarpXParticleContainer::DepositCurrent(WarpXParIter& pti,
BL_PROFILE_VAR_START(blp_accumulate);
// CPU, tiling: atomicAdd local_jx into jx
// (same for jx and jz)
- (*jx)[pti].atomicAdd(local_jx[thread_num], tbx, tbx, 0, 0, 1);
- (*jy)[pti].atomicAdd(local_jy[thread_num], tby, tby, 0, 0, 1);
- (*jz)[pti].atomicAdd(local_jz[thread_num], tbz, tbz, 0, 0, 1);
+ (*jx)[pti].atomicAdd(local_jx[thread_num], tbx, tbx, 0, 0, local_jx[thread_num].nComp());
+ (*jy)[pti].atomicAdd(local_jy[thread_num], tby, tby, 0, 0, local_jy[thread_num].nComp());
+ (*jz)[pti].atomicAdd(local_jz[thread_num], tbz, tbz, 0, 0, local_jz[thread_num].nComp());
BL_PROFILE_VAR_STOP(blp_accumulate);
#endif
}
@@ -440,11 +442,11 @@ WarpXParticleContainer::DepositCharge ( WarpXParIter& pti, RealVector& wp,
const std::array<Real, 3>& xyzmin = xyzmin_tile;
#ifdef AMREX_USE_GPU
- data_ptr = (*rhomf)[pti].dataPtr(icomp);
+ data_ptr = (*rhomf)[pti].dataPtr(icomp*(rhomf->nComp()/2));
auto rholen = (*rhomf)[pti].length();
#else
tile_box.grow(ngRho);
- local_rho[thread_num].resize(tile_box);
+ local_rho[thread_num].resize(tile_box, rhomf->nComp());
data_ptr = local_rho[thread_num].dataPtr();
auto rholen = local_rho[thread_num].length();
@@ -483,7 +485,7 @@ WarpXParticleContainer::DepositCharge ( WarpXParIter& pti, RealVector& wp,
#ifndef AMREX_USE_GPU
BL_PROFILE_VAR_START(blp_accumulate);
- (*rhomf)[pti].atomicAdd(local_rho[thread_num], tile_box, tile_box, 0, icomp, 1);
+ (*rhomf)[pti].atomicAdd(local_rho[thread_num], tile_box, tile_box, 0, icomp*(rhomf->nComp()/2), (rhomf->nComp()/2));
BL_PROFILE_VAR_STOP(blp_accumulate);
#endif
@@ -502,7 +504,7 @@ WarpXParticleContainer::DepositCharge ( WarpXParIter& pti, RealVector& wp,
#else
tile_box = amrex::convert(ctilebox, IntVect::TheUnitVector());
tile_box.grow(ngRho);
- local_rho[thread_num].resize(tile_box);
+ local_rho[thread_num].resize(tile_box, crhomf->nComp());
data_ptr = local_rho[thread_num].dataPtr();
auto rholen = local_rho[thread_num].length();
@@ -543,7 +545,7 @@ WarpXParticleContainer::DepositCharge ( WarpXParIter& pti, RealVector& wp,
#ifndef AMREX_USE_GPU
BL_PROFILE_VAR_START(blp_accumulate);
- (*crhomf)[pti].atomicAdd(local_rho[thread_num], tile_box, tile_box, 0, icomp, 1);
+ (*crhomf)[pti].atomicAdd(local_rho[thread_num], tile_box, tile_box, 0, icomp*(crhomf->nComp()/2), (crhomf->nComp()/2));
BL_PROFILE_VAR_STOP(blp_accumulate);
#endif
@@ -598,7 +600,7 @@ WarpXParticleContainer::DepositCharge (Vector<std::unique_ptr<MultiFab> >& rho,
BoxArray coarsened_fine_BA = fine_BA;
coarsened_fine_BA.coarsen(m_gdb->refRatio(lev));
- MultiFab coarsened_fine_data(coarsened_fine_BA, fine_dm, 1, 0);
+ MultiFab coarsened_fine_data(coarsened_fine_BA, fine_dm, rho[lev+1]->nComp(), 0);
coarsened_fine_data.setVal(0.0);
IntVect ratio(AMREX_D_DECL(2, 2, 2)); // FIXME
@@ -629,7 +631,7 @@ WarpXParticleContainer::GetChargeDensity (int lev, bool local)
const int ng = WarpX::nox;
- auto rho = std::unique_ptr<MultiFab>(new MultiFab(nba,dm,1,ng));
+ auto rho = std::unique_ptr<MultiFab>(new MultiFab(nba,dm,WarpX::ncomps,ng));
rho->setVal(0.0);
#ifdef _OPENMP
@@ -657,7 +659,7 @@ WarpXParticleContainer::GetChargeDensity (int lev, bool local)
Box tile_box = convert(pti.tilebox(), IntVect::TheUnitVector());
const std::array<Real, 3>& xyzmin = xyzmin_tile;
tile_box.grow(ng);
- rho_loc.resize(tile_box);
+ rho_loc.resize(tile_box, rho->nComp());
rho_loc = 0.0;
data_ptr = rho_loc.dataPtr();
auto rholen = rho_loc.length();