aboutsummaryrefslogtreecommitdiff
path: root/Source
diff options
context:
space:
mode:
Diffstat (limited to 'Source')
-rw-r--r--Source/Initialization/PlasmaInjector.H1
-rw-r--r--Source/Initialization/PlasmaInjector.cpp1
-rw-r--r--Source/Particles/PhysicalParticleContainer.H5
-rw-r--r--Source/Particles/PhysicalParticleContainer.cpp97
4 files changed, 72 insertions, 32 deletions
diff --git a/Source/Initialization/PlasmaInjector.H b/Source/Initialization/PlasmaInjector.H
index cd5802a55..f998e217e 100644
--- a/Source/Initialization/PlasmaInjector.H
+++ b/Source/Initialization/PlasmaInjector.H
@@ -293,6 +293,7 @@ public:
amrex::Real z_rms;
amrex::Real q_tot;
long npart;
+ int do_symmetrize = 0;
bool radially_weighted = true;
diff --git a/Source/Initialization/PlasmaInjector.cpp b/Source/Initialization/PlasmaInjector.cpp
index 52b5d8b61..f9642d1b6 100644
--- a/Source/Initialization/PlasmaInjector.cpp
+++ b/Source/Initialization/PlasmaInjector.cpp
@@ -286,6 +286,7 @@ PlasmaInjector::PlasmaInjector(int ispecies, const std::string& name)
pp.get("z_rms", z_rms);
pp.get("q_tot", q_tot);
pp.get("npart", npart);
+ pp.query("do_symmetrize", do_symmetrize);
gaussian_beam = true;
parseMomentum(pp);
}
diff --git a/Source/Particles/PhysicalParticleContainer.H b/Source/Particles/PhysicalParticleContainer.H
index df1e0296a..bd8cfb47e 100644
--- a/Source/Particles/PhysicalParticleContainer.H
+++ b/Source/Particles/PhysicalParticleContainer.H
@@ -94,7 +94,10 @@ public:
void AddGaussianBeam(amrex::Real x_m, amrex::Real y_m, amrex::Real z_m,
amrex::Real x_rms, amrex::Real y_rms, amrex::Real z_rms,
- amrex::Real q_tot, long npart);
+ amrex::Real q_tot, long npart, int do_symmetrize);
+
+ void CheckAndAddParticle(amrex::Real x, amrex::Real y, amrex::Real z,
+ std::array<amrex::Real, 3> u, amrex::Real weight);
virtual void GetParticleSlice(const int direction, const amrex::Real z_old,
const amrex::Real z_new, const amrex::Real t_boost,
diff --git a/Source/Particles/PhysicalParticleContainer.cpp b/Source/Particles/PhysicalParticleContainer.cpp
index ab37bf6ff..1f517fccb 100644
--- a/Source/Particles/PhysicalParticleContainer.cpp
+++ b/Source/Particles/PhysicalParticleContainer.cpp
@@ -181,7 +181,8 @@ void PhysicalParticleContainer::MapParticletoBoostedFrame(Real& x, Real& y, Real
void
PhysicalParticleContainer::AddGaussianBeam(Real x_m, Real y_m, Real z_m,
Real x_rms, Real y_rms, Real z_rms,
- Real q_tot, long npart) {
+ Real q_tot, long npart,
+ int do_symmetrize) {
const Geometry& geom = m_gdb->Geom(0);
RealBox containing_bx = geom.ProbDomain();
@@ -191,13 +192,15 @@ PhysicalParticleContainer::AddGaussianBeam(Real x_m, Real y_m, Real z_m,
std::normal_distribution<double> disty(y_m, y_rms);
std::normal_distribution<double> distz(z_m, z_rms);
- std::array<Real,PIdx::nattribs> attribs;
- attribs.fill(0.0);
-
if (ParallelDescriptor::IOProcessor()) {
- std::array<Real, 3> u;
- Real weight;
- for (long i = 0; i < npart; ++i) {
+ std::array<Real, 3> u;
+ Real weight;
+ // If do_symmetrize, create 4x fewer particles, and
+ // Replicate each particle 4 times (x,y) (-x,y) (x,-y) (-x,-y)
+ if (do_symmetrize){
+ npart /= 4;
+ }
+ for (long i = 0; i < npart; ++i) {
#if ( AMREX_SPACEDIM == 3 | WARPX_RZ)
weight = q_tot/npart/charge;
Real x = distx(mt);
@@ -209,29 +212,27 @@ PhysicalParticleContainer::AddGaussianBeam(Real x_m, Real y_m, Real z_m,
Real y = 0.;
Real z = distz(mt);
#endif
- if (plasma_injector->insideBounds(x, y, z)) {
- plasma_injector->getMomentum(u, x, y, z);
- if (WarpX::gamma_boost > 1.) {
- MapParticletoBoostedFrame(x, y, z, u);
- }
- attribs[PIdx::ux] = u[0];
- attribs[PIdx::uy] = u[1];
- attribs[PIdx::uz] = u[2];
- attribs[PIdx::w ] = weight;
-
- if (WarpX::do_boosted_frame_diagnostic && do_boosted_frame_diags)
- {
- auto& particle_tile = DefineAndReturnParticleTile(0, 0, 0);
- particle_tile.push_back_real(particle_comps["xold"], x);
- particle_tile.push_back_real(particle_comps["yold"], y);
- particle_tile.push_back_real(particle_comps["zold"], z);
-
- particle_tile.push_back_real(particle_comps["uxold"], u[0]);
- particle_tile.push_back_real(particle_comps["uyold"], u[1]);
- particle_tile.push_back_real(particle_comps["uzold"], u[2]);
- }
-
- AddOneParticle(0, 0, 0, x, y, z, attribs);
+ if (plasma_injector->insideBounds(x, y, z)) {
+ plasma_injector->getMomentum(u, x, y, z);
+ if (do_symmetrize){
+ std::array<Real, 3> u_tmp;
+ Real x_tmp, y_tmp;
+ // Add four particles to the beam:
+ // (x,ux,y,uy) (-x,-ux,y,uy) (x,ux,-y,-uy) (-x,-ux,-y,-uy)
+ for (int ix=0; ix<2; ix++){
+ for (int iy=0; iy<2; iy++){
+ u_tmp = u;
+ x_tmp = x*std::pow(-1,ix);
+ u_tmp[0] *= std::pow(-1,ix);
+ y_tmp = y*std::pow(-1,iy);
+ u_tmp[1] *= std::pow(-1,iy);
+ CheckAndAddParticle(x_tmp, y_tmp, z,
+ u_tmp, weight/4);
+ }
+ }
+ } else {
+ CheckAndAddParticle(x, y, z, u, weight);
+ }
}
}
}
@@ -239,6 +240,39 @@ PhysicalParticleContainer::AddGaussianBeam(Real x_m, Real y_m, Real z_m,
}
void
+PhysicalParticleContainer::CheckAndAddParticle(Real x, Real y, Real z,
+ std::array<Real, 3> u,
+ Real weight)
+{
+ std::array<Real,PIdx::nattribs> attribs;
+ attribs.fill(0.0);
+
+ // update attribs with input arguments
+ if (WarpX::gamma_boost > 1.) {
+ MapParticletoBoostedFrame(x, y, z, u);
+ }
+ attribs[PIdx::ux] = u[0];
+ attribs[PIdx::uy] = u[1];
+ attribs[PIdx::uz] = u[2];
+ attribs[PIdx::w ] = weight;
+
+ if (WarpX::do_boosted_frame_diagnostic && do_boosted_frame_diags)
+ {
+ // need to create old values
+ auto& particle_tile = DefineAndReturnParticleTile(0, 0, 0);
+ particle_tile.push_back_real(particle_comps["xold"], x);
+ particle_tile.push_back_real(particle_comps["yold"], y);
+ particle_tile.push_back_real(particle_comps["zold"], z);
+
+ particle_tile.push_back_real(particle_comps["uxold"], u[0]);
+ particle_tile.push_back_real(particle_comps["uyold"], u[1]);
+ particle_tile.push_back_real(particle_comps["uzold"], u[2]);
+ }
+ // add particle
+ AddOneParticle(0, 0, 0, x, y, z, attribs);
+}
+
+void
PhysicalParticleContainer::AddParticles (int lev)
{
BL_PROFILE("PhysicalParticleContainer::AddParticles()");
@@ -263,7 +297,8 @@ PhysicalParticleContainer::AddParticles (int lev)
plasma_injector->y_rms,
plasma_injector->z_rms,
plasma_injector->q_tot,
- plasma_injector->npart);
+ plasma_injector->npart,
+ plasma_injector->do_symmetrize);
return;