diff options
-rw-r--r-- | Docs/source/running_cpp/parameters.rst | 25 | ||||
-rw-r--r-- | Source/Initialization/PlasmaInjector.H | 1 | ||||
-rw-r--r-- | Source/Initialization/PlasmaInjector.cpp | 1 | ||||
-rw-r--r-- | Source/Particles/PhysicalParticleContainer.H | 5 | ||||
-rw-r--r-- | Source/Particles/PhysicalParticleContainer.cpp | 97 |
5 files changed, 89 insertions, 40 deletions
diff --git a/Docs/source/running_cpp/parameters.rst b/Docs/source/running_cpp/parameters.rst index ba64c1fcd..c38ef472d 100644 --- a/Docs/source/running_cpp/parameters.rst +++ b/Docs/source/running_cpp/parameters.rst @@ -187,6 +187,15 @@ Particle initialization * ``NRandomPerCell``: injection with a fixed number of randomly-distributed particles per cell. This requires the additional parameter ``<species_name>.num_particles_per_cell``. + * ``gaussian_beam``: Inject particle beam with gaussian distribution in + space in all directions. This requires additional parameters: + ``<species_name>.q_tot`` (beam charge), + ``<species_name>.npart`` (number of particles in the beam), + ``<species_name>.x/y/z_m`` (average position in `x/y/z`), + ``<species_name>.x/y/z_rms`` (standard deviation in `x/y/z`), + and optional argument ``<species_name>.do_symmetrize`` (whether to + symmetrize the beam in the x and y directions). + * ``<species_name>.do_continuous_injection`` (`0` or `1`) Whether to inject particles during the simulation, and not only at initialization. This can be required whith a moving window and/or when @@ -466,14 +475,14 @@ Laser initialization See definition in Akturk et al., Opt Express, vol 12, no 19 (2014). * ``<laser_name>.do_continuous_injection`` (`0` or `1`) optional (default `0`). - Whether or not to use continuous injection (`0` or not `0`). - If the antenna starts outside of the simulation domain but enters it - at some point (due to moving window or moving antenna in the boosted - frame), use this so that the laser antenna is injected when it reaches - the box boundary. If running in a boosted frame, this requires the - boost direction, moving window direction and laser propagation direction - to be along `z`. If not running in a boosted frame, this requires the - moving window and laser propagation directions to be the same (`x`, `y` + Whether or not to use continuous injection. + If the antenna starts outside of the simulation domain but enters it + at some point (due to moving window or moving antenna in the boosted + frame), use this so that the laser antenna is injected when it reaches + the box boundary. If running in a boosted frame, this requires the + boost direction, moving window direction and laser propagation direction + to be along `z`. If not running in a boosted frame, this requires the + moving window and laser propagation directions to be the same (`x`, `y` or `z`) * ``warpx.num_mirrors`` (`int`) optional (default `0`) 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; |