diff options
Diffstat (limited to 'Source/Particles')
| -rw-r--r-- | Source/Particles/MultiParticleContainer.H | 16 | ||||
| -rw-r--r-- | Source/Particles/MultiParticleContainer.cpp | 87 | ||||
| -rw-r--r-- | Source/Particles/PhysicalParticleContainer.H | 21 | ||||
| -rw-r--r-- | Source/Particles/PhysicalParticleContainer.cpp | 107 | ||||
| -rw-r--r-- | Source/Particles/RigidInjectedParticleContainer.cpp | 7 | ||||
| -rw-r--r-- | Source/Particles/WarpXParticleContainer.cpp | 1 |
6 files changed, 208 insertions, 31 deletions
diff --git a/Source/Particles/MultiParticleContainer.H b/Source/Particles/MultiParticleContainer.H index 9db129b05..ed1c2f371 100644 --- a/Source/Particles/MultiParticleContainer.H +++ b/Source/Particles/MultiParticleContainer.H @@ -8,6 +8,7 @@ #include <RigidInjectedParticleContainer.H> #include <PhotonParticleContainer.H> #include <LaserParticleContainer.H> +#include <WarpXParserWrapper.H> #include <AMReX_Particles.H> #ifdef WARPX_QED @@ -215,6 +216,21 @@ public: IonizationProcess ionization_process; + std::string m_B_ext_particle_s = "default"; + std::string m_E_ext_particle_s = "default"; + // External fields added to particle fields. + amrex::Vector<amrex::Real> m_B_external_particle; + amrex::Vector<amrex::Real> m_E_external_particle; + // ParserWrapper for B_external on the particle + std::unique_ptr<ParserWrapper> m_Bx_particle_parser; + std::unique_ptr<ParserWrapper> m_By_particle_parser; + std::unique_ptr<ParserWrapper> m_Bz_particle_parser; + // ParserWrapper for E_external on the particle + std::unique_ptr<ParserWrapper> m_Ex_particle_parser; + std::unique_ptr<ParserWrapper> m_Ey_particle_parser; + std::unique_ptr<ParserWrapper> m_Ez_particle_parser; + + protected: // Particle container types diff --git a/Source/Particles/MultiParticleContainer.cpp b/Source/Particles/MultiParticleContainer.cpp index d84bc1afa..ab836ce9d 100644 --- a/Source/Particles/MultiParticleContainer.cpp +++ b/Source/Particles/MultiParticleContainer.cpp @@ -72,6 +72,93 @@ MultiParticleContainer::ReadParameters () { ParmParse pp("particles"); + // allocating and initializing default values of external fields for particles + m_E_external_particle.resize(3); + m_B_external_particle.resize(3); + // initialize E and B fields to 0.0 + for (int idim = 0; idim < 3; ++idim) { + m_E_external_particle[idim] = 0.0; + m_B_external_particle[idim] = 0.0; + } + // default values of E_external_particle and B_external_particle + // are used to set the E and B field when "constant" or "parser" + // is not explicitly used in the input + pp.query("B_ext_particle_init_style", m_B_ext_particle_s); + std::transform(m_B_ext_particle_s.begin(), + m_B_ext_particle_s.end(), + m_B_ext_particle_s.begin(), + ::tolower); + pp.query("E_ext_particle_init_style", m_E_ext_particle_s); + std::transform(m_E_ext_particle_s.begin(), + m_E_ext_particle_s.end(), + m_E_ext_particle_s.begin(), + ::tolower); + // if the input string for B_external on particles is "constant" + // then the values for the external B on particles must + // be provided in the input file. + if (m_B_ext_particle_s == "constant") + pp.getarr("B_external_particle", m_B_external_particle); + + // if the input string for E_external on particles is "constant" + // then the values for the external E on particles must + // be provided in the input file. + if (m_E_ext_particle_s == "constant") + pp.getarr("E_external_particle", m_E_external_particle); + + // if the input string for B_ext_particle_s is + // "parse_b_ext_particle_function" then the mathematical expression + // for the Bx_, By_, Bz_external_particle_function(x,y,z) + // must be provided in the input file. + if (m_B_ext_particle_s == "parse_b_ext_particle_function") { + // store the mathematical expression as string + std::string str_Bx_ext_particle_function; + std::string str_By_ext_particle_function; + std::string str_Bz_ext_particle_function; + Store_parserString(pp, "Bx_external_particle_function(x,y,z,t)", + str_Bx_ext_particle_function); + Store_parserString(pp, "By_external_particle_function(x,y,z,t)", + str_By_ext_particle_function); + Store_parserString(pp, "Bz_external_particle_function(x,y,z,t)", + str_Bz_ext_particle_function); + + // Parser for B_external on the particle + m_Bx_particle_parser.reset(new ParserWrapper( + makeParser(str_Bx_ext_particle_function))); + m_By_particle_parser.reset(new ParserWrapper( + makeParser(str_By_ext_particle_function))); + m_Bz_particle_parser.reset(new ParserWrapper( + makeParser(str_Bz_ext_particle_function))); + + } + + // if the input string for E_ext_particle_s is + // "parse_e_ext_particle_function" then the mathematical expression + // for the Ex_, Ey_, Ez_external_particle_function(x,y,z) + // must be provided in the input file. + if (m_E_ext_particle_s == "parse_e_ext_particle_function") { + // store the mathematical expression as string + std::string str_Ex_ext_particle_function; + std::string str_Ey_ext_particle_function; + std::string str_Ez_ext_particle_function; + Store_parserString(pp, "Ex_external_particle_function(x,y,z,t)", + str_Ex_ext_particle_function); + Store_parserString(pp, "Ey_external_particle_function(x,y,z,t)", + str_Ey_ext_particle_function); + Store_parserString(pp, "Ez_external_particle_function(x,y,z,t)", + str_Ez_ext_particle_function); + // Parser for E_external on the particle + m_Ex_particle_parser.reset(new ParserWrapper( + makeParser(str_Ex_ext_particle_function))); + m_Ey_particle_parser.reset(new ParserWrapper( + makeParser(str_Ey_ext_particle_function))); + m_Ez_particle_parser.reset(new ParserWrapper( + makeParser(str_Ez_ext_particle_function))); + + } + + + + pp.query("nspecies", nspecies); BL_ASSERT(nspecies >= 0); diff --git a/Source/Particles/PhysicalParticleContainer.H b/Source/Particles/PhysicalParticleContainer.H index 0192ffb37..74d1a0f62 100644 --- a/Source/Particles/PhysicalParticleContainer.H +++ b/Source/Particles/PhysicalParticleContainer.H @@ -47,6 +47,26 @@ public: amrex::Real t, amrex::Real dt) override; #endif // WARPX_DO_ELECTROSTATIC + /** + * \brief Apply external E and B fields on the particles. The E and B + * fields could be defined as a constant or using a parser for reading + * in a mathematical expression. The default value for the E- and B-fields + * is (0.0,0.0,0.0). + * + * \param[in,out] Exp-Bzp pointer to fields on particles modified based + * on external E and B + * \param[in] xp,yp,zp arrays of particle positions required to compute + * mathematical expression for the external fields + * using parser. + */ + void AssignExternalFieldOnParticles ( WarpXParIter& pti, + RealVector& Exp, RealVector& Eyp, + RealVector& Ezp, RealVector& Bxp, + RealVector& Byp, RealVector& Bzp, + amrex::Gpu::ManagedDeviceVector<amrex::ParticleReal> xp, + amrex::Gpu::ManagedDeviceVector<amrex::ParticleReal> yp, + amrex::Gpu::ManagedDeviceVector<amrex::ParticleReal> zp, int lev); + virtual void FieldGather (int lev, const amrex::MultiFab& Ex, const amrex::MultiFab& Ey, @@ -318,6 +338,7 @@ protected: //radiation reaction bool do_classical_radiation_reaction = false; + #ifdef WARPX_QED // A flag to enable quantum_synchrotron process for leptons bool m_do_qed_quantum_sync = false; diff --git a/Source/Particles/PhysicalParticleContainer.cpp b/Source/Particles/PhysicalParticleContainer.cpp index 936a64c9c..21c8ddd31 100644 --- a/Source/Particles/PhysicalParticleContainer.cpp +++ b/Source/Particles/PhysicalParticleContainer.cpp @@ -967,6 +967,80 @@ PhysicalParticleContainer::EvolveES (const Vector<std::array<std::unique_ptr<Mul #endif // WARPX_DO_ELECTROSTATIC void +PhysicalParticleContainer::AssignExternalFieldOnParticles(WarpXParIter& pti, + RealVector& Exp, RealVector& Eyp, RealVector& Ezp, + RealVector& Bxp, RealVector& Byp, RealVector& Bzp, + Gpu::ManagedDeviceVector<ParticleReal> xp, + Gpu::ManagedDeviceVector<ParticleReal> yp, + Gpu::ManagedDeviceVector<ParticleReal> zp, int lev) +{ + const long np = pti.numParticles(); + /// get WarpX class object + auto & warpx = WarpX::GetInstance(); + /// get MultiParticleContainer class object + auto & mypc = warpx.GetPartContainer(); + if (mypc.m_E_ext_particle_s=="constant" || + mypc.m_E_ext_particle_s=="default") { + Exp.assign(np,mypc.m_E_external_particle[0]); + Eyp.assign(np,mypc.m_E_external_particle[1]); + Ezp.assign(np,mypc.m_E_external_particle[2]); + } + if (mypc.m_B_ext_particle_s=="constant" || + mypc.m_B_ext_particle_s=="default") { + Bxp.assign(np,mypc.m_B_external_particle[0]); + Byp.assign(np,mypc.m_B_external_particle[1]); + Bzp.assign(np,mypc.m_B_external_particle[2]); + } + if (mypc.m_E_ext_particle_s=="parse_e_ext_particle_function") { + Real* const AMREX_RESTRICT xp_data = xp.dataPtr(); + Real* const AMREX_RESTRICT yp_data = yp.dataPtr(); + Real* const AMREX_RESTRICT zp_data = zp.dataPtr(); + Real* const AMREX_RESTRICT Exp_data = Exp.dataPtr(); + Real* const AMREX_RESTRICT Eyp_data = Eyp.dataPtr(); + Real* const AMREX_RESTRICT Ezp_data = Ezp.dataPtr(); + ParserWrapper *xfield_partparser = mypc.m_Ex_particle_parser.get(); + ParserWrapper *yfield_partparser = mypc.m_Ey_particle_parser.get(); + ParserWrapper *zfield_partparser = mypc.m_Ez_particle_parser.get(); + Real time = warpx.gett_new(lev); + amrex::ParallelFor(pti.numParticles(), + [=] AMREX_GPU_DEVICE (long i) { + Exp_data[i] = xfield_partparser->getField(xp_data[i],yp_data[i],zp_data[i],time); + Eyp_data[i] = yfield_partparser->getField(xp_data[i],yp_data[i],zp_data[i],time); + Ezp_data[i] = zfield_partparser->getField(xp_data[i],yp_data[i],zp_data[i],time); + }, + /* To allocate shared memory for the GPU threads. */ + /* But, for now only 4 doubles (x,y,z,t) are allocated. */ + amrex::Gpu::numThreadsPerBlockParallelFor() * sizeof(double) * 4 + ); + } + if (mypc.m_B_ext_particle_s=="parse_b_ext_particle_function") { + Real* const AMREX_RESTRICT xp_data = xp.dataPtr(); + Real* const AMREX_RESTRICT yp_data = yp.dataPtr(); + Real* const AMREX_RESTRICT zp_data = zp.dataPtr(); + Real* const AMREX_RESTRICT Bxp_data = Bxp.dataPtr(); + Real* const AMREX_RESTRICT Byp_data = Byp.dataPtr(); + Real* const AMREX_RESTRICT Bzp_data = Bzp.dataPtr(); + ParserWrapper *xfield_partparser = mypc.m_Bx_particle_parser.get(); + ParserWrapper *yfield_partparser = mypc.m_By_particle_parser.get(); + ParserWrapper *zfield_partparser = mypc.m_Bz_particle_parser.get(); + Real time = warpx.gett_new(lev); + amrex::ParallelFor(pti.numParticles(), + [=] AMREX_GPU_DEVICE (long i) { + Bxp_data[i] = xfield_partparser->getField(xp_data[i],yp_data[i],zp_data[i],time); + Byp_data[i] = yfield_partparser->getField(xp_data[i],yp_data[i],zp_data[i],time); + Bzp_data[i] = zfield_partparser->getField(xp_data[i],yp_data[i],zp_data[i],time); + }, + /* To allocate shared memory for the GPU threads. */ + /* But, for now only 4 doubles (x,y,z,t) are allocated. */ + amrex::Gpu::numThreadsPerBlockParallelFor() * sizeof(double) * 4 + ); + } + +} + + + +void PhysicalParticleContainer::FieldGather (int lev, const amrex::MultiFab& Ex, const amrex::MultiFab& Ey, @@ -1015,13 +1089,6 @@ PhysicalParticleContainer::FieldGather (int lev, const FArrayBox& byfab = By[pti]; const FArrayBox& bzfab = Bz[pti]; - Exp.assign(np,WarpX::E_external_particle[0]); - Eyp.assign(np,WarpX::E_external_particle[1]); - Ezp.assign(np,WarpX::E_external_particle[2]); - Bxp.assign(np,WarpX::B_external_particle[0]); - Byp.assign(np,WarpX::B_external_particle[1]); - Bzp.assign(np,WarpX::B_external_particle[2]); - // // copy data from particle container to temp arrays // @@ -1152,14 +1219,6 @@ PhysicalParticleContainer::Evolve (int lev, exfab, eyfab, ezfab, bxfab, byfab, bzfab); } - Exp.assign(np,WarpX::E_external_particle[0]); - Eyp.assign(np,WarpX::E_external_particle[1]); - Ezp.assign(np,WarpX::E_external_particle[2]); - - Bxp.assign(np,WarpX::B_external_particle[0]); - Byp.assign(np,WarpX::B_external_particle[1]); - Bzp.assign(np,WarpX::B_external_particle[2]); - // Determine which particles deposit/gather in the buffer, and // which particles deposit/gather in the fine patch long nfine_current = np; @@ -1799,14 +1858,6 @@ PhysicalParticleContainer::PushP (int lev, Real dt, const FArrayBox& byfab = By[pti]; const FArrayBox& bzfab = Bz[pti]; - Exp.assign(np,WarpX::E_external_particle[0]); - Eyp.assign(np,WarpX::E_external_particle[1]); - Ezp.assign(np,WarpX::E_external_particle[2]); - - Bxp.assign(np,WarpX::B_external_particle[0]); - Byp.assign(np,WarpX::B_external_particle[1]); - Bzp.assign(np,WarpX::B_external_particle[2]); - // // copy data from particle container to temp arrays // @@ -2176,9 +2227,16 @@ PhysicalParticleContainer::FieldGather (WarpXParIter& pti, AMREX_ALWAYS_ASSERT_WITH_MESSAGE((gather_lev==(lev-1)) || (gather_lev==(lev )), "Gather buffers only work for lev-1"); - // If no particles, do not do anything if (np_to_gather == 0) return; + + // initializing the field value to the externally applied field before + // gathering fields from the grid to the particles. + AssignExternalFieldOnParticles(pti, Exp, Eyp, Ezp, Bxp, Byp, Bzp, + m_xp[thread_num], m_yp[thread_num], + m_zp[thread_num], lev); + + // Get cell size on gather_lev const std::array<Real,3>& dx = WarpX::CellSize(std::max(gather_lev,0)); @@ -2424,4 +2482,5 @@ set_quantum_sync_engine_ptr(std::shared_ptr<QuantumSynchrotronEngine> ptr) { m_shr_p_qs_engine = ptr; } + #endif diff --git a/Source/Particles/RigidInjectedParticleContainer.cpp b/Source/Particles/RigidInjectedParticleContainer.cpp index bee71fba1..f3b502a0a 100644 --- a/Source/Particles/RigidInjectedParticleContainer.cpp +++ b/Source/Particles/RigidInjectedParticleContainer.cpp @@ -392,13 +392,6 @@ RigidInjectedParticleContainer::PushP (int lev, Real dt, const FArrayBox& byfab = By[pti]; const FArrayBox& bzfab = Bz[pti]; - Exp.assign(np,WarpX::E_external_particle[0]); - Eyp.assign(np,WarpX::E_external_particle[1]); - Ezp.assign(np,WarpX::E_external_particle[2]); - Bxp.assign(np,WarpX::B_external_particle[0]); - Byp.assign(np,WarpX::B_external_particle[1]); - Bzp.assign(np,WarpX::B_external_particle[2]); - // // copy data from particle container to temp arrays // diff --git a/Source/Particles/WarpXParticleContainer.cpp b/Source/Particles/WarpXParticleContainer.cpp index 6abc02139..85c4dc0d5 100644 --- a/Source/Particles/WarpXParticleContainer.cpp +++ b/Source/Particles/WarpXParticleContainer.cpp @@ -113,6 +113,7 @@ WarpXParticleContainer::WarpXParticleContainer (AmrCore* amr_core, int ispecies) m_xp.resize(num_threads); m_yp.resize(num_threads); m_zp.resize(num_threads); + } void |