diff options
Diffstat (limited to 'Source/Particles/PhysicalParticleContainer.cpp')
| -rw-r--r-- | Source/Particles/PhysicalParticleContainer.cpp | 187 |
1 files changed, 153 insertions, 34 deletions
diff --git a/Source/Particles/PhysicalParticleContainer.cpp b/Source/Particles/PhysicalParticleContainer.cpp index 02dee1967..d10803320 100644 --- a/Source/Particles/PhysicalParticleContainer.cpp +++ b/Source/Particles/PhysicalParticleContainer.cpp @@ -1,6 +1,8 @@ #include <limits> #include <sstream> +#include <PhysicalParticleContainer.H> + #include <MultiParticleContainer.H> #include <WarpX_f.H> #include <WarpX.H> @@ -49,6 +51,33 @@ PhysicalParticleContainer::PhysicalParticleContainer (AmrCore* amr_core, int isp "version of Redistribute in AMReX does not work with runtime parameters"); #endif + +#ifdef WARPX_QED + //Add real component if QED is enabled + pp.query("do_qed", do_qed); + if(do_qed) + AddRealComp("tau"); + + //IF do_qed is enabled, find out if Quantum Synchrotron process is enabled + if(do_qed) + pp.query("do_qed_quantum_sync", do_qed_quantum_sync); + + //TODO: SHOULD CHECK IF SPECIES IS EITHER ELECTRONS OR POSITRONS!! +#endif + + //variable to set plot_flags size + int plot_flag_size = PIdx::nattribs; + if(WarpX::do_boosted_frame_diagnostic && do_boosted_frame_diags) + plot_flag_size += 6; + +#ifdef WARPX_QED + if(do_qed){ + // plot_flag will have an entry for the optical depth + plot_flag_size++; + } +#endif + //_______________________________ + pp.query("plot_species", plot_species); int do_user_plot_vars; do_user_plot_vars = pp.queryarr("plot_vars", plot_vars); @@ -56,18 +85,11 @@ PhysicalParticleContainer::PhysicalParticleContainer (AmrCore* amr_core, int isp // By default, all particle variables are dumped to plotfiles, // including {x,y,z,ux,uy,uz}old variables when running in a // boosted frame - if (WarpX::do_boosted_frame_diagnostic && do_boosted_frame_diags){ - plot_flags.resize(PIdx::nattribs + 6, 1); - } else { - plot_flags.resize(PIdx::nattribs, 1); - } + plot_flags.resize(plot_flag_size, 1); } else { // Set plot_flag to 0 for all attribs - if (WarpX::do_boosted_frame_diagnostic && do_boosted_frame_diags){ - plot_flags.resize(PIdx::nattribs + 6, 0); - } else { - plot_flags.resize(PIdx::nattribs, 0); - } + plot_flags.resize(plot_flag_size, 0); + // If not none, set plot_flags values to 1 for elements in plot_vars. if (plot_vars[0] != "none"){ for (const auto& var : plot_vars){ @@ -79,6 +101,13 @@ PhysicalParticleContainer::PhysicalParticleContainer (AmrCore* amr_core, int isp } } } + + #ifdef WARPX_QED + if(do_qed){ + //Optical depths is always plotted if QED is on + plot_flags[plot_flag_size-1] = 1; + } + #endif } PhysicalParticleContainer::PhysicalParticleContainer (AmrCore* amr_core) @@ -92,7 +121,9 @@ void PhysicalParticleContainer::InitData() // Init ionization module here instead of in the PhysicalParticleContainer // constructor because dt is required if (do_field_ionization) {InitIonizationModule();} + AddParticles(0); // Note - add on level 0 + Redistribute(); // We then redistribute } @@ -156,6 +187,16 @@ 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); + // Allocate temporary vectors on the CPU + Gpu::HostVector<ParticleReal> particle_x; + Gpu::HostVector<ParticleReal> particle_y; + Gpu::HostVector<ParticleReal> particle_z; + Gpu::HostVector<ParticleReal> particle_ux; + Gpu::HostVector<ParticleReal> particle_uy; + Gpu::HostVector<ParticleReal> particle_uz; + Gpu::HostVector<ParticleReal> particle_w; + int np = 0; + if (ParallelDescriptor::IOProcessor()) { // If do_symmetrize, create 4x fewer particles, and // Replicate each particle 4 times (x,y) (-x,y) (x,-y) (-x,-y) @@ -181,44 +222,61 @@ PhysicalParticleContainer::AddGaussianBeam(Real x_m, Real y_m, Real z_m, u.z *= PhysConst::c; if (do_symmetrize){ // Add four particles to the beam: - CheckAndAddParticle( x, y, z, { u.x, u.y, u.z}, weight/4. ); - CheckAndAddParticle( x,-y, z, { u.x,-u.y, u.z}, weight/4. ); - CheckAndAddParticle(-x, y, z, {-u.x, u.y, u.z}, weight/4. ); - CheckAndAddParticle(-x,-y, z, {-u.x,-u.y, u.z}, weight/4. ); + CheckAndAddParticle(x, y, z, { u.x, u.y, u.z}, weight/4., + particle_x, particle_y, particle_z, + particle_ux, particle_uy, particle_uz, + particle_w); + CheckAndAddParticle(x, -y, z, { u.x, -u.y, u.z}, weight/4., + particle_x, particle_y, particle_z, + particle_ux, particle_uy, particle_uz, + particle_w); + CheckAndAddParticle(-x, y, z, { -u.x, u.y, u.z}, weight/4., + particle_x, particle_y, particle_z, + particle_ux, particle_uy, particle_uz, + particle_w); + CheckAndAddParticle(-x, -y, z, { -u.x, -u.y, u.z}, weight/4., + particle_x, particle_y, particle_z, + particle_ux, particle_uy, particle_uz, + particle_w); } else { - CheckAndAddParticle(x, y, z, {u.x,u.y,u.z}, weight); + CheckAndAddParticle(x, y, z, { u.x, u.y, u.z}, weight, + particle_x, particle_y, particle_z, + particle_ux, particle_uy, particle_uz, + particle_w); } } } } - Redistribute(); + // Add the temporary CPU vectors to the particle structure + np = particle_z.size(); + AddNParticles(0,np, + particle_x.dataPtr(), particle_y.dataPtr(), particle_z.dataPtr(), + particle_ux.dataPtr(), particle_uy.dataPtr(), particle_uz.dataPtr(), + 1, particle_w.dataPtr(),1); } void PhysicalParticleContainer::CheckAndAddParticle(Real x, Real y, Real z, std::array<Real, 3> u, - Real weight) + Real weight, + Gpu::HostVector<ParticleReal>& particle_x, + Gpu::HostVector<ParticleReal>& particle_y, + Gpu::HostVector<ParticleReal>& particle_z, + Gpu::HostVector<ParticleReal>& particle_ux, + Gpu::HostVector<ParticleReal>& particle_uy, + Gpu::HostVector<ParticleReal>& particle_uz, + Gpu::HostVector<ParticleReal>& particle_w) { - std::array<ParticleReal,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 ( (NumRuntimeRealComps()>0) || (NumRuntimeIntComps()>0) ) - { - // need to create old values - auto& particle_tile = DefineAndReturnParticleTile(0, 0, 0); - } - - // add particle - AddOneParticle(0, 0, 0, x, y, z, attribs); + particle_x.push_back(x); + particle_y.push_back(y); + particle_z.push_back(z); + particle_ux.push_back(u[0]); + particle_uy.push_back(u[1]); + particle_uz.push_back(u[2]); + particle_w.push_back(weight); } void @@ -450,6 +508,30 @@ PhysicalParticleContainer::AddPlasma (int lev, RealBox part_realbox) pi = soa.GetIntData(particle_icomps["ionization_level"]).data() + old_size; } +#ifdef WARPX_QED + //Pointer to the optical depth component + amrex::Real* p_tau; + + //If a QED effect is enabled, tau has to be initialized + bool loc_has_quantum_sync = has_quantum_sync(); + bool loc_has_breit_wheeler = has_breit_wheeler(); + if(loc_has_quantum_sync || loc_has_breit_wheeler){ + p_tau = soa.GetRealData(particle_comps["tau"]).data() + old_size; + } + + //If needed, get the appropriate functors from the engines + QuantumSynchrotronGetOpticalDepth quantum_sync_get_opt; + BreitWheelerGetOpticalDepth breit_wheeler_get_opt; + if(loc_has_quantum_sync){ + quantum_sync_get_opt = + shr_ptr_qs_engine->build_optical_depth_functor(); + } + if(loc_has_breit_wheeler){ + breit_wheeler_get_opt = + shr_ptr_bw_engine->build_optical_depth_functor(); + } +#endif + const GpuArray<Real,AMREX_SPACEDIM> overlap_corner {AMREX_D_DECL(overlap_realbox.lo(0), overlap_realbox.lo(1), @@ -597,6 +679,16 @@ PhysicalParticleContainer::AddPlasma (int lev, RealBox part_realbox) pi[ip] = loc_ionization_initial_level; } +#ifdef WARPX_QED + if(loc_has_quantum_sync){ + p_tau[ip] = quantum_sync_get_opt(); + } + + if(loc_has_breit_wheeler){ + p_tau[ip] = breit_wheeler_get_opt(); + } +#endif + u.x *= PhysConst::c; u.y *= PhysConst::c; u.z *= PhysConst::c; @@ -2047,3 +2139,30 @@ PhysicalParticleContainer::buildIonizationMask (const amrex::MFIter& mfi, const } ); } + +#ifdef WARPX_QED + +bool PhysicalParticleContainer::has_quantum_sync() +{ + return do_qed_quantum_sync; +} + +bool PhysicalParticleContainer::has_breit_wheeler() +{ + return do_qed_breit_wheeler; +} + +void +PhysicalParticleContainer:: +set_breit_wheeler_engine_ptr(std::shared_ptr<BreitWheelerEngine> ptr) +{ + shr_ptr_bw_engine = ptr; +} + +void +PhysicalParticleContainer:: +set_quantum_sync_engine_ptr(std::shared_ptr<QuantumSynchrotronEngine> ptr) +{ + shr_ptr_qs_engine = ptr; +} +#endif |