diff options
Diffstat (limited to 'Source')
| -rw-r--r-- | Source/Laser/LaserParticleContainer.H | 12 | ||||
| -rw-r--r-- | Source/Laser/LaserParticleContainer.cpp | 98 | ||||
| -rw-r--r-- | Source/Particles/MultiParticleContainer.H | 10 | ||||
| -rw-r--r-- | Source/Particles/MultiParticleContainer.cpp | 45 | ||||
| -rw-r--r-- | Source/Particles/PhysicalParticleContainer.H | 5 | ||||
| -rw-r--r-- | Source/Particles/PhysicalParticleContainer.cpp | 18 | ||||
| -rw-r--r-- | Source/Particles/WarpXParticleContainer.H | 20 | ||||
| -rw-r--r-- | Source/Utils/WarpXMovingWindow.cpp | 25 | ||||
| -rw-r--r-- | Source/WarpX.H | 4 | ||||
| -rw-r--r-- | Source/WarpX.cpp | 30 |
10 files changed, 218 insertions, 49 deletions
diff --git a/Source/Laser/LaserParticleContainer.H b/Source/Laser/LaserParticleContainer.H index 5dd94d0c7..516e368ab 100644 --- a/Source/Laser/LaserParticleContainer.H +++ b/Source/Laser/LaserParticleContainer.H @@ -84,11 +84,19 @@ private: std::string field_function; // laser particle domain - amrex::RealBox prob_domain; - + amrex::RealBox laser_injection_box; + // Theoretical position of the antenna. Used if do_continuous_injection=1. + // Track the position of the antenna until it enters the simulation domain. + amrex::Vector<amrex::Real> updated_position; + void ComputeSpacing (int lev, amrex::Real& Sx, amrex::Real& Sy) const; void ComputeWeightMobility (amrex::Real Sx, amrex::Real Sy); void InitData (int lev); + // Inject the laser antenna during the simulation, if it started + // outside of the simulation domain and enters it. + void ContinuousInjection(const amrex::RealBox& injection_box) override; + // Update position of the antenna + void UpdateContinuousInjectionPosition(amrex::Real dt) override; }; #endif diff --git a/Source/Laser/LaserParticleContainer.cpp b/Source/Laser/LaserParticleContainer.cpp index 2b56c3cfd..8a2590e5e 100644 --- a/Source/Laser/LaserParticleContainer.cpp +++ b/Source/Laser/LaserParticleContainer.cpp @@ -49,6 +49,7 @@ LaserParticleContainer::LaserParticleContainer (AmrCore* amr_core, int ispecies, pp.query("pusher_algo", pusher_algo); pp.get("wavelength", wavelength); pp.get("e_max", e_max); + pp.query("do_continuous_injection", do_continuous_injection); if ( profile == laser_t::Gaussian ) { // Parse the properties of the Gaussian profile @@ -148,16 +149,94 @@ LaserParticleContainer::LaserParticleContainer (AmrCore* amr_core, int ispecies, u_Y = {0., 1., 0.}; #endif - prob_domain = Geometry::ProbDomain(); + laser_injection_box= Geometry::ProbDomain(); { Vector<Real> lo, hi; if (pp.queryarr("prob_lo", lo)) { - prob_domain.setLo(lo); + laser_injection_box.setLo(lo); } if (pp.queryarr("prob_hi", hi)) { - prob_domain.setHi(hi); + laser_injection_box.setHi(hi); } } + + if (do_continuous_injection){ + // If laser antenna initially outside of the box, store its theoretical + // position in z_antenna_th + updated_position = position; + + // Sanity checks + int dir = WarpX::moving_window_dir; + std::vector<Real> windir(3, 0.0); +#if (AMREX_SPACEDIM==2) + windir[2*dir] = 1.0; +#else + windir[dir] = 1.0; +#endif + AMREX_ALWAYS_ASSERT_WITH_MESSAGE( + (nvec[0]-windir[0]) + (nvec[1]-windir[1]) + (nvec[2]-windir[2]) + < 1.e-12, "do_continous_injection for laser particle only works" + + " if moving window direction and laser propagation direction are the same"); + if ( WarpX::gamma_boost>1 ){ + AMREX_ALWAYS_ASSERT_WITH_MESSAGE( + (WarpX::boost_direction[0]-0)*(WarpX::boost_direction[0]-0) + + (WarpX::boost_direction[1]-0)*(WarpX::boost_direction[1]-0) + + (WarpX::boost_direction[2]-1)*(WarpX::boost_direction[2]-1) < 1.e-12, + "do_continous_injection for laser particle only works if " + + "warpx.boost_direction = z. TODO: all directions."); + } + } +} + +/* \brief Check if laser particles enter the box, and inject if necessary. + * \param injection_box: a RealBox where particles should be injected. + */ +void +LaserParticleContainer::ContinuousInjection (const RealBox& injection_box) +{ + // Input parameter injection_box contains small box where injection + // should occur. + // So far, LaserParticleContainer::laser_injection_box contains the + // outdated full problem domain at t=0. + + // Convert updated_position to Real* to use RealBox::contains(). +#if (AMREX_SPACEDIM == 3) + const Real* p_pos = updated_position.dataPtr(); +#else + const Real p_pos[2] = {updated_position[0], updated_position[2]}; +#endif + if ( injection_box.contains(p_pos) ){ + // Update laser_injection_box with current value + laser_injection_box = injection_box; + // Inject laser particles. LaserParticleContainer::InitData + // is called only once, when the antenna enters the simulation + // domain. + InitData(); + } +} + +/* \brief update position of the antenna if running in boosted frame. + * \param dt time step (level 0). + * The up-to-date antenna position is stored in updated_position. + */ +void +LaserParticleContainer::UpdateContinuousInjectionPosition(Real dt) +{ + int dir = WarpX::moving_window_dir; + if (do_continuous_injection and (WarpX::gamma_boost > 1)){ + // In boosted-frame simulations, the antenna has moved since the last + // call to this function, and injection position needs to be updated +#if ( AMREX_SPACEDIM == 3 ) + updated_position[dir] -= WarpX::beta_boost * + WarpX::boost_direction[dir] * PhysConst::c * dt; +#elif ( AMREX_SPACEDIM == 2 ) + // In 2D, dir=0 corresponds to x and dir=1 corresponds to z + // This needs to be converted in order to index `boost_direction` + // which has 3 components, for both 2D and 3D simulations. + updated_position[2*dir] -= WarpX::beta_boost * + WarpX::boost_direction[2*dir] * PhysConst::c * dt; +#endif + } } void @@ -175,6 +254,13 @@ LaserParticleContainer::InitData (int lev) ComputeSpacing(lev, S_X, S_Y); ComputeWeightMobility(S_X, S_Y); + // LaserParticleContainer::position contains the initial position of the + // laser antenna. In the boosted frame, the antenna is moving. + // Update its position with updated_position. + if (do_continuous_injection){ + position = updated_position; + } + auto Transform = [&](int i, int j) -> Vector<Real>{ #if (AMREX_SPACEDIM == 3) return { position[0] + (S_X*(i+0.5))*u_X[0] + (S_Y*(j+0.5))*u_Y[0], @@ -210,8 +296,8 @@ LaserParticleContainer::InitData (int lev) plane_hi[1] = std::max(plane_hi[1], j); }; - const Real* prob_lo = prob_domain.lo(); - const Real* prob_hi = prob_domain.hi(); + const Real* prob_lo = laser_injection_box.lo(); + const Real* prob_hi = laser_injection_box.hi(); #if (AMREX_SPACEDIM == 3) compute_min_max(prob_lo[0], prob_lo[1], prob_lo[2]); compute_min_max(prob_hi[0], prob_lo[1], prob_lo[2]); @@ -272,7 +358,7 @@ LaserParticleContainer::InitData (int lev) #else const Real x[2] = {pos[0], pos[2]}; #endif - if (prob_domain.contains(x)) + if (laser_injection_box.contains(x)) { for (int k = 0; k<2; ++k) { particle_x.push_back(pos[0]); diff --git a/Source/Particles/MultiParticleContainer.H b/Source/Particles/MultiParticleContainer.H index f3fd522a9..0c5e49c04 100644 --- a/Source/Particles/MultiParticleContainer.H +++ b/Source/Particles/MultiParticleContainer.H @@ -169,7 +169,15 @@ public: const amrex::Real z_old, const amrex::Real z_new, const amrex::Real t_boost, const amrex::Real t_lab, const amrex::Real dt, amrex::Vector<WarpXParticleContainer::DiagnosticParticleData>& parts) const; - + + // Inject particles during the simulation (for particles entering the + // simulation domain after some iterations, due to flowing plasma and/or + // moving window). + void ContinuousInjection(const amrex::RealBox& injection_box) const; + // Update injection position for continuously-injected species. + void UpdateContinuousInjectionPosition(amrex::Real dt) const; + int doContinuousInjection() const; + // // Parameters for the Cherenkov corrector in the FDTD solver. // Both stencils are calculated ar runtime. diff --git a/Source/Particles/MultiParticleContainer.cpp b/Source/Particles/MultiParticleContainer.cpp index a4df1f83a..440906348 100644 --- a/Source/Particles/MultiParticleContainer.cpp +++ b/Source/Particles/MultiParticleContainer.cpp @@ -414,3 +414,48 @@ MultiParticleContainer } } } + +/* \brief Continuous injection for particles initially outside of the domain. + * \param injection_box: Domain where new particles should be injected. + * Loop over all WarpXParticleContainer in MultiParticleContainer and + * calls virtual function ContinuousInjection. + */ +void +MultiParticleContainer::ContinuousInjection(const RealBox& injection_box) const +{ + for (int i=0; i<nspecies+nlasers; i++){ + auto& pc = allcontainers[i]; + if (pc->do_continuous_injection){ + pc->ContinuousInjection(injection_box); + } + } +} + +/* \brief Update position of continuous injection parameters. + * \param dt: simulation time step (level 0) + * All classes inherited from WarpXParticleContainer do not have + * a position to update (PhysicalParticleContainer does not do anything). + */ +void +MultiParticleContainer::UpdateContinuousInjectionPosition(Real dt) const +{ + for (int i=0; i<nspecies+nlasers; i++){ + auto& pc = allcontainers[i]; + if (pc->do_continuous_injection){ + pc->UpdateContinuousInjectionPosition(dt); + } + } +} + +int +MultiParticleContainer::doContinuousInjection() const +{ + int warpx_do_continuous_injection = 0; + for (int i=0; i<nspecies+nlasers; i++){ + auto& pc = allcontainers[i]; + if (pc->do_continuous_injection){ + warpx_do_continuous_injection = 1; + } + } + return warpx_do_continuous_injection; +} diff --git a/Source/Particles/PhysicalParticleContainer.H b/Source/Particles/PhysicalParticleContainer.H index 362683879..4f365768b 100644 --- a/Source/Particles/PhysicalParticleContainer.H +++ b/Source/Particles/PhysicalParticleContainer.H @@ -101,8 +101,6 @@ public: const amrex::Real t_lab, const amrex::Real dt, DiagnosticParticles& diagnostic_particles) final; - bool injected = false; - protected: std::string species_name; @@ -122,6 +120,9 @@ protected: int GetRefineFac(const amrex::Real x, const amrex::Real y, const amrex::Real z); std::unique_ptr<amrex::IArrayBox> m_refined_injection_mask = nullptr; + // Inject particles during the whole simulation + void ContinuousInjection(const amrex::RealBox& injection_box) override; + }; #endif diff --git a/Source/Particles/PhysicalParticleContainer.cpp b/Source/Particles/PhysicalParticleContainer.cpp index 17e6d98d9..2fa39d87d 100644 --- a/Source/Particles/PhysicalParticleContainer.cpp +++ b/Source/Particles/PhysicalParticleContainer.cpp @@ -24,7 +24,7 @@ NumParticlesToAdd(const Box& overlap_box, const RealBox& overlap_realbox, for (IntVect iv = overlap_box.smallEnd(); iv <= overlap_box.bigEnd(); overlap_box.next(iv)) { int fac; - if (injected) { + if (do_continuous_injection) { #if ( AMREX_SPACEDIM == 3 ) Real x = overlap_corner[0] + (iv[0] + 0.5)*dx[0]; Real y = overlap_corner[1] + (iv[1] + 0.5)*dx[1]; @@ -81,6 +81,7 @@ PhysicalParticleContainer::PhysicalParticleContainer (AmrCore* amr_core, int isp pp.query("do_backward_propagation", do_backward_propagation); pp.query("do_splitting", do_splitting); pp.query("split_type", split_type); + pp.query("do_continuous_injection", do_continuous_injection); } PhysicalParticleContainer::PhysicalParticleContainer (AmrCore* amr_core) @@ -361,7 +362,7 @@ PhysicalParticleContainer::AddPlasmaCPU (int lev, RealBox part_realbox) for (IntVect iv = overlap_box.smallEnd(); iv <= overlap_box.bigEnd(); overlap_box.next(iv)) { int fac; - if (injected) { + if (do_continuous_injection) { #if ( AMREX_SPACEDIM == 3 ) Real x = overlap_corner[0] + (iv[0] + 0.5)*dx[0]; Real y = overlap_corner[1] + (iv[1] + 0.5)*dx[1]; @@ -602,7 +603,7 @@ PhysicalParticleContainer::AddPlasmaGPU (int lev, RealBox part_realbox) for (IntVect iv = overlap_box.smallEnd(); iv <= overlap_box.bigEnd(); overlap_box.next(iv)) { int fac; - if (injected) { + if (do_continuous_injection) { #if ( AMREX_SPACEDIM == 3 ) Real x = overlap_corner[0] + (iv[0] + 0.5)*dx[0]; Real y = overlap_corner[1] + (iv[1] + 0.5)*dx[1]; @@ -2004,3 +2005,14 @@ int PhysicalParticleContainer::GetRefineFac(const Real x, const Real y, const Re return ref_fac; } + +/* \brief Inject particles during the simulation + * \param injection_box: domain where particles should be injected. + */ +void +PhysicalParticleContainer::ContinuousInjection(const RealBox& injection_box) +{ + // Inject plasma on level 0. Paticles will be redistributed. + const int lev=0; + AddPlasma(lev, injection_box); +} diff --git a/Source/Particles/WarpXParticleContainer.H b/Source/Particles/WarpXParticleContainer.H index 275554cd8..6fa02b824 100644 --- a/Source/Particles/WarpXParticleContainer.H +++ b/Source/Particles/WarpXParticleContainer.H @@ -183,7 +183,18 @@ public: int thread_num, int lev, amrex::Real dt ); - + + // If particles start outside of the domain, ContinuousInjection + // makes sure that they are initialized when they enter the domain, and + // NOT before. Virtual function, overriden by derived classes. + // Current status: + // PhysicalParticleContainer: implemented. + // LaserParticleContainer: implemented. + // RigidInjectedParticleContainer: not implemented. + virtual void ContinuousInjection(const amrex::RealBox& injection_box) {} + // Update optional sub-class-specific injection location. + virtual void UpdateContinuousInjectionPosition(amrex::Real dt) {} + /// /// This returns the total charge for all the particles in this ParticleContainer. /// This is needed when solving Poisson's equation with periodic boundary conditions. @@ -248,6 +259,12 @@ protected: static int do_not_push; + // Whether to allow particles outside of the simulation domain to be + // initialized when they enter the domain. + // This is currently required because continuous injection does not + // support all features allowed by direct injection. + int do_continuous_injection = 0; + amrex::Vector<amrex::FArrayBox> local_rho; amrex::Vector<amrex::FArrayBox> local_jx; amrex::Vector<amrex::FArrayBox> local_jy; @@ -258,6 +275,7 @@ protected: private: virtual void particlePostLocate(ParticleType& p, const amrex::ParticleLocData& pld, const int lev) override; + }; #endif diff --git a/Source/Utils/WarpXMovingWindow.cpp b/Source/Utils/WarpXMovingWindow.cpp index 05e171f22..18d89951d 100644 --- a/Source/Utils/WarpXMovingWindow.cpp +++ b/Source/Utils/WarpXMovingWindow.cpp @@ -9,7 +9,7 @@ WarpX::UpdatePlasmaInjectionPosition (Real dt) { int dir = moving_window_dir; // Continuously inject plasma in new cells (by default only on level 0) - if (WarpX::do_plasma_injection and (WarpX::gamma_boost > 1)){ + if (WarpX::warpx_do_continuous_injection and (WarpX::gamma_boost > 1)){ // In boosted-frame simulations, the plasma has moved since the last // call to this function, and injection position needs to be updated current_injection_position -= WarpX::beta_boost * @@ -33,7 +33,16 @@ WarpX::MoveWindow (bool move_j) // and of the plasma injection moving_window_x += moving_window_v * dt[0]; int dir = moving_window_dir; + + // Update warpx.current_injection_position + // PhysicalParticleContainer uses this injection position UpdatePlasmaInjectionPosition( dt[0] ); + if (WarpX::warpx_do_continuous_injection){ + // Update injection position for WarpXParticleContainer in mypc. + // Nothing to do for PhysicalParticleContainers + // For LaserParticleContainer, need to update the antenna position. + mypc->UpdateContinuousInjectionPosition( dt[0] ); + } // compute the number of cells to shift on the base level Real new_lo[AMREX_SPACEDIM]; @@ -134,7 +143,7 @@ WarpX::MoveWindow (bool move_j) } // Continuously inject plasma in new cells (by default only on level 0) - if (WarpX::do_plasma_injection) { + if (WarpX::warpx_do_continuous_injection) { const int lev = 0; @@ -162,15 +171,11 @@ WarpX::MoveWindow (bool move_j) particleBox.setLo( dir, new_injection_position ); particleBox.setHi( dir, current_injection_position ); } - // Perform the injection of new particles in particleBox + if (particleBox.ok() and (current_injection_position != new_injection_position)){ - for (int i = 0; i < num_injected_species; ++i) { - int ispecies = injected_plasma_species[i]; - WarpXParticleContainer& pc = mypc->GetParticleContainer(ispecies); - auto& ppc = dynamic_cast<PhysicalParticleContainer&>(pc); - ppc.AddPlasma(lev, particleBox); - } - // Update the injection position + // Performs continuous injection of all WarpXParticleContainer + // in mypc. + mypc->ContinuousInjection(particleBox); current_injection_position = new_injection_position; } } diff --git a/Source/WarpX.H b/Source/WarpX.H index 6cec8e5be..44e84033f 100644 --- a/Source/WarpX.H +++ b/Source/WarpX.H @@ -476,10 +476,10 @@ private: amrex::Real current_injection_position = 0; // Plasma injection parameters - int do_plasma_injection = 0; + int warpx_do_continuous_injection = 0; int num_injected_species = -1; amrex::Vector<int> injected_plasma_species; - + int do_electrostatic = 0; int n_buffer = 4; amrex::Real const_dt = 0.5e-11; diff --git a/Source/WarpX.cpp b/Source/WarpX.cpp index 47ead98df..a3a24897a 100644 --- a/Source/WarpX.cpp +++ b/Source/WarpX.cpp @@ -145,13 +145,14 @@ WarpX::WarpX () // Particle Container mypc = std::unique_ptr<MultiParticleContainer> (new MultiParticleContainer(this)); - - if (do_plasma_injection) { - for (int i = 0; i < num_injected_species; ++i) { - int ispecies = injected_plasma_species[i]; - WarpXParticleContainer& pc = mypc->GetParticleContainer(ispecies); - auto& ppc = dynamic_cast<PhysicalParticleContainer&>(pc); - ppc.injected = true; + warpx_do_continuous_injection = mypc->doContinuousInjection(); + if (warpx_do_continuous_injection){ + if (moving_window_v >= 0){ + // Inject particles continuously from the right end of the box + current_injection_position = geom[0].ProbHi(moving_window_dir); + } else { + // Inject particles continuously from the left end of the box + current_injection_position = geom[0].ProbLo(moving_window_dir); } } @@ -300,21 +301,6 @@ WarpX::ReadParameters () moving_window_v *= PhysConst::c; } - pp.query("do_plasma_injection", do_plasma_injection); - if (do_plasma_injection) { - pp.get("num_injected_species", num_injected_species); - injected_plasma_species.resize(num_injected_species); - pp.getarr("injected_plasma_species", injected_plasma_species, - 0, num_injected_species); - if (moving_window_v >= 0){ - // Inject particles continuously from the right end of the box - current_injection_position = geom[0].ProbHi(moving_window_dir); - } else { - // Inject particles continuously from the left end of the box - current_injection_position = geom[0].ProbLo(moving_window_dir); - } - } - pp.query("do_boosted_frame_diagnostic", do_boosted_frame_diagnostic); if (do_boosted_frame_diagnostic) { |