aboutsummaryrefslogtreecommitdiff
path: root/Source
diff options
context:
space:
mode:
Diffstat (limited to 'Source')
-rw-r--r--Source/Particles/MultiParticleContainer.H1
-rw-r--r--Source/Particles/PhotonParticleContainer.H14
-rw-r--r--Source/Particles/PhotonParticleContainer.cpp51
-rw-r--r--Source/Particles/PhysicalParticleContainer.H38
-rw-r--r--Source/Particles/PhysicalParticleContainer.cpp119
-rw-r--r--Source/Particles/WarpXParticleContainer.H3
-rw-r--r--Source/QED/BreitWheelerEngineWrapper.H8
-rw-r--r--Source/QED/BreitWheelerEngineWrapper.cpp6
-rw-r--r--Source/QED/QuantumSyncEngineWrapper.H10
-rw-r--r--Source/QED/QuantumSyncEngineWrapper.cpp6
10 files changed, 235 insertions, 21 deletions
diff --git a/Source/Particles/MultiParticleContainer.H b/Source/Particles/MultiParticleContainer.H
index f9a0e51d7..06898a5c1 100644
--- a/Source/Particles/MultiParticleContainer.H
+++ b/Source/Particles/MultiParticleContainer.H
@@ -11,6 +11,7 @@
#include <AMReX_Particles.H>
#ifdef WARPX_QED
+ #include <QedChiFunctions.H>
#include <BreitWheelerEngineWrapper.H>
#include <QuantumSyncEngineWrapper.H>
#endif
diff --git a/Source/Particles/PhotonParticleContainer.H b/Source/Particles/PhotonParticleContainer.H
index 851726d57..580f13f86 100644
--- a/Source/Particles/PhotonParticleContainer.H
+++ b/Source/Particles/PhotonParticleContainer.H
@@ -46,7 +46,7 @@ public:
amrex::Gpu::ManagedDeviceVector<amrex::ParticleReal>& zp,
amrex::Real dt, DtType a_dt_type=DtType::Full) override;
- // Don't push momenta for photons
+ // Do nothing
virtual void PushP (int lev,
amrex::Real dt,
const amrex::MultiFab& Ex,
@@ -79,6 +79,18 @@ public:
return false;
};
+#ifdef WARPX_QED
+ /**
+ * This function evolves the optical depth of the photons if QED effects
+ * are enabled.
+ * @param[in,out] pti particle iterator (optical depth will be modified)
+ * @param[in] dt temporal step
+ */
+ virtual void EvolveOpticalDepth(WarpXParIter& pti,
+ amrex::Real dt) override;
+
+#endif
+
};
#endif // #ifndef WARPX_PhotonParticleContainer_H_
diff --git a/Source/Particles/PhotonParticleContainer.cpp b/Source/Particles/PhotonParticleContainer.cpp
index fd47ac8a0..424569c9a 100644
--- a/Source/Particles/PhotonParticleContainer.cpp
+++ b/Source/Particles/PhotonParticleContainer.cpp
@@ -78,8 +78,6 @@ PhotonParticleContainer::PushPX(WarpXParIter& pti,
copy_attribs(pti, x, y, z);
}
- //No need to update momentum for photons (for now)
-
amrex::ParallelFor(
pti.numParticles(),
[=] AMREX_GPU_DEVICE (long i) {
@@ -88,7 +86,6 @@ PhotonParticleContainer::PushPX(WarpXParIter& pti,
ux[i], uy[i], uz[i], dt );
}
);
-
}
void
@@ -116,3 +113,51 @@ PhotonParticleContainer::Evolve (int lev,
t, dt);
}
+
+#ifdef WARPX_QED
+
+void
+PhotonParticleContainer::EvolveOpticalDepth(
+ WarpXParIter& pti,amrex::Real dt)
+{
+ //m_shr_p_bw_engine->are_lookup_tables_initialized() is necessary here if we want
+ //to perform just initialization tests of the optical depth without actually
+ //enabling QED effects (this requires lookup tables).
+ if(!has_breit_wheeler() || !m_shr_p_bw_engine->are_lookup_tables_initialized())
+ return;
+
+ auto& attribs = pti.GetAttribs();
+ ParticleReal* const AMREX_RESTRICT ux = attribs[PIdx::ux].dataPtr();
+ ParticleReal* const AMREX_RESTRICT uy = attribs[PIdx::uy].dataPtr();
+ ParticleReal* const AMREX_RESTRICT uz = attribs[PIdx::uz].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT Ex = attribs[PIdx::Ex].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT Ey = attribs[PIdx::Ey].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT Ez = attribs[PIdx::Ez].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT Bx = attribs[PIdx::Bx].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT By = attribs[PIdx::By].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT Bz = attribs[PIdx::Bz].dataPtr();
+
+ BreitWheelerEvolveOpticalDepth evolve_opt =
+ m_shr_p_bw_engine->build_evolve_functor();
+
+ amrex::Real* AMREX_RESTRICT p_tau =
+ pti.GetAttribs(particle_comps["tau"]).dataPtr();
+
+ const auto me = PhysConst::m_e;
+
+ amrex::ParallelFor(
+ pti.numParticles(),
+ [=] AMREX_GPU_DEVICE (long i) {
+ const ParticleReal px = me * ux[i];
+ const ParticleReal py = me * uy[i];
+ const ParticleReal pz = me * uz[i];
+
+ bool has_event_happened = evolve_opt(
+ px, py, pz,
+ Ex[i], Ey[i], Ez[i],
+ Bx[i], By[i], Bz[i],
+ dt, p_tau[i]);
+ }
+ );
+}
+#endif
diff --git a/Source/Particles/PhysicalParticleContainer.H b/Source/Particles/PhysicalParticleContainer.H
index 174488eb9..e70b470b8 100644
--- a/Source/Particles/PhysicalParticleContainer.H
+++ b/Source/Particles/PhysicalParticleContainer.H
@@ -187,13 +187,47 @@ public:
amrex::FArrayBox const * & byfab, amrex::FArrayBox const * & bzfab);
#ifdef WARPX_QED
+ //Functions decleared in WarpXParticleContainer.H
+ //containers for which QED processes could be relevant
+ //are expected to override these functions
+
+ /**
+ * Tells if this PhysicalParticleContainer has Quantum
+ * Synchrotron process enabled
+ * @return true if process is enabled
+ */
bool has_quantum_sync() override;
+
+ /**
+ * Tells if this PhysicalParticleContainer has Breit
+ * Wheeler process enabled
+ * @return true if process is enabled
+ */
bool has_breit_wheeler() override;
+ /**
+ * Acquires a shared smart pointer to a BreitWheelerEngine
+ * @param[in] ptr the pointer
+ */
void set_breit_wheeler_engine_ptr
- (std::shared_ptr<BreitWheelerEngine>) override;
+ (std::shared_ptr<BreitWheelerEngine> ptr) override;
+
+ /**
+ * Acquires a shared smart pointer to a QuantumSynchrotronEngine
+ * @param[in] ptr the pointer
+ */
void set_quantum_sync_engine_ptr
- (std::shared_ptr<QuantumSynchrotronEngine>) override;
+ (std::shared_ptr<QuantumSynchrotronEngine> ptr) override;
+ //__________
+
+ /**
+ * This function evolves the optical depth of the particles if QED effects
+ * are enabled.
+ * @param[in,out] pti particle iterator (optical depth will be modified)
+ * @param[in] dt temporal step
+ */
+ virtual void EvolveOpticalDepth(WarpXParIter& pti,
+ amrex::Real dt);
#endif
protected:
diff --git a/Source/Particles/PhysicalParticleContainer.cpp b/Source/Particles/PhysicalParticleContainer.cpp
index 938c80de5..9c13d34e1 100644
--- a/Source/Particles/PhysicalParticleContainer.cpp
+++ b/Source/Particles/PhysicalParticleContainer.cpp
@@ -1052,6 +1052,7 @@ PhysicalParticleContainer::Evolve (int lev,
BL_PROFILE("PPC::Evolve()");
BL_PROFILE_VAR_NS("PPC::Evolve::Copy", blp_copy);
BL_PROFILE_VAR_NS("PPC::FieldGather", blp_fg);
+ BL_PROFILE_VAR_NS("PPC::EvolveOpticalDepth", blp_ppc_qed_ev);
BL_PROFILE_VAR_NS("PPC::ParticlePush", blp_ppc_pp);
const std::array<Real,3>& dx = WarpX::CellSize(lev);
@@ -1243,6 +1244,15 @@ PhysicalParticleContainer::Evolve (int lev,
BL_PROFILE_VAR_STOP(blp_fg);
+#ifdef WARPX_QED
+ //
+ //Evolve Optical Depth
+ //
+ BL_PROFILE_VAR_START(blp_ppc_qed_ev);
+ EvolveOpticalDepth(pti, dt);
+ BL_PROFILE_VAR_STOP(blp_ppc_qed_ev);
+#endif
+
//
// Particle Push
//
@@ -1590,8 +1600,45 @@ PhysicalParticleContainer::PushPX(WarpXParIter& pti,
const Real q = this->charge;
const Real m = this-> mass;
+#ifdef WARPX_QED
- //Assumes that all consistency checks have been done at initialization
+ auto t_chi_max = m_shr_p_qs_engine->get_ref_ctrl().chi_part_min;
+
+ if(do_classical_radiation_reaction){
+ if(m_do_qed_quantum_sync){
+ amrex::ParallelFor(
+ pti.numParticles(),
+ [=] AMREX_GPU_DEVICE (long i) {
+ auto chi = QedUtils::chi_lepton(m*ux[i], m*uy[i], m*uz[i],
+ Ex[i], Ey[i], Ez[i],
+ Bx[i], By[i], Bz[i]);
+ if(chi < t_chi_max){
+ UpdateMomentumBorisWithRadiationReaction( ux[i], uy[i], uz[i],
+ Ex[i], Ey[i], Ez[i], Bx[i],
+ By[i], Bz[i], q, m, dt);
+ }
+ else{
+ UpdateMomentumBoris( ux[i], uy[i], uz[i],
+ Ex[i], Ey[i], Ez[i], Bx[i],
+ By[i], Bz[i], q, m, dt);
+ }
+ UpdatePosition( x[i], y[i], z[i],
+ ux[i], uy[i], uz[i], dt );
+ }
+ );
+ }else{
+ amrex::ParallelFor(
+ pti.numParticles(),
+ [=] AMREX_GPU_DEVICE (long i) {
+ UpdateMomentumBorisWithRadiationReaction( ux[i], uy[i], uz[i],
+ Ex[i], Ey[i], Ez[i], Bx[i],
+ By[i], Bz[i], q, m, dt);
+ UpdatePosition( x[i], y[i], z[i],
+ ux[i], uy[i], uz[i], dt );
+ }
+ );
+ }
+#else
if(do_classical_radiation_reaction){
amrex::ParallelFor(
pti.numParticles(),
@@ -1605,6 +1652,7 @@ PhysicalParticleContainer::PushPX(WarpXParIter& pti,
ux[i], uy[i], uz[i], dt );
}
);
+#endif
} else if (WarpX::particle_pusher_algo == ParticlePusherAlgo::Boris){
amrex::ParallelFor(
pti.numParticles(),
@@ -1649,6 +1697,52 @@ PhysicalParticleContainer::PushPX(WarpXParIter& pti,
};
}
+#ifdef WARPX_QED
+void PhysicalParticleContainer::EvolveOpticalDepth(
+ WarpXParIter& pti, amrex::Real dt)
+{
+ //m_shr_p_qs_engine->are_lookup_tables_initialized() is necessary here if we want
+ //to perform just initialization tests of the optical depth without actually
+ //enabling QED effects (this requires lookup tables).
+ if(!has_quantum_sync() || !m_shr_p_qs_engine->are_lookup_tables_initialized())
+ return;
+
+ QuantumSynchrotronEvolveOpticalDepth evolve_opt =
+ m_shr_p_qs_engine->build_evolve_functor();
+
+ auto& attribs = pti.GetAttribs();
+ const ParticleReal* const AMREX_RESTRICT ux = attribs[PIdx::ux].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT uy = attribs[PIdx::uy].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT uz = attribs[PIdx::uz].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT Ex = attribs[PIdx::Ex].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT Ey = attribs[PIdx::Ey].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT Ez = attribs[PIdx::Ez].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT Bx = attribs[PIdx::Bx].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT By = attribs[PIdx::By].dataPtr();
+ const ParticleReal* const AMREX_RESTRICT Bz = attribs[PIdx::Bz].dataPtr();
+
+ ParticleReal* const AMREX_RESTRICT p_tau =
+ pti.GetAttribs(particle_comps["tau"]).dataPtr();
+
+ const ParticleReal m = this->mass;
+
+ amrex::ParallelFor(pti.numParticles(),
+ [=] AMREX_GPU_DEVICE (long i) {
+ const ParticleReal px = m * ux[i];
+ const ParticleReal py = m * uy[i];
+ const ParticleReal pz = m * uz[i];
+
+ bool has_event_happened = evolve_opt(
+ px, py, pz,
+ Ex[i], Ey[i], Ez[i],
+ Bx[i], By[i], Bz[i],
+ dt, p_tau[i]);
+ }
+ );
+
+}
+#endif
+
void
PhysicalParticleContainer::PushP (int lev, Real dt,
const MultiFab& Ex, const MultiFab& Ey, const MultiFab& Ez,
@@ -1732,12 +1826,10 @@ PhysicalParticleContainer::PushP (int lev, Real dt,
ion_lev = pti.GetiAttribs(particle_icomps["ionization_level"]).dataPtr();
}
-
//Assumes that all consistency checks have been done at initialization
if(do_classical_radiation_reaction){
- amrex::ParallelFor(
- pti.numParticles(),
- [=] AMREX_GPU_DEVICE (long i) {
+ amrex::ParallelFor(pti.numParticles(),
+ [=] AMREX_GPU_DEVICE (long i){
Real qp = q;
if (ion_lev){ qp *= ion_lev[i]; }
UpdateMomentumBorisWithRadiationReaction(
@@ -1748,7 +1840,7 @@ PhysicalParticleContainer::PushP (int lev, Real dt,
}
);
} else if (WarpX::particle_pusher_algo == ParticlePusherAlgo::Boris){
- amrex::ParallelFor( pti.numParticles(),
+ amrex::ParallelFor(pti.numParticles(),
[=] AMREX_GPU_DEVICE (long i) {
Real qp = q;
if (ion_lev){ qp *= ion_lev[i]; }
@@ -1759,8 +1851,8 @@ PhysicalParticleContainer::PushP (int lev, Real dt,
qp, m, dt);
}
);
- } else if (WarpX::particle_pusher_algo == ParticlePusherAlgo::Vay) {
- amrex::ParallelFor( pti.numParticles(),
+ } else if (WarpX::particle_pusher_algo == ParticlePusherAlgo::Vay){
+ amrex::ParallelFor(pti.numParticles(),
[=] AMREX_GPU_DEVICE (long i) {
Real qp = q;
if (ion_lev){ qp *= ion_lev[i]; }
@@ -1771,16 +1863,19 @@ PhysicalParticleContainer::PushP (int lev, Real dt,
qp, m, dt);
}
);
- } else if (WarpX::particle_pusher_algo == ParticlePusherAlgo::HigueraCary) {
- amrex::ParallelFor( pti.numParticles(),
+ } else if (WarpX::particle_pusher_algo == ParticlePusherAlgo::HigueraCary){
+ amrex::ParallelFor(pti.numParticles(),
[=] AMREX_GPU_DEVICE (long i) {
UpdateMomentumHigueraCary( ux[i], uy[i], uz[i],
- Expp[i], Eypp[i], Ezpp[i], Bxpp[i], Bypp[i], Bzpp[i], q, m, dt);
+ Expp[i], Eypp[i], Ezpp[i],
+ Bxpp[i], Bypp[i], Bzpp[i],
+ q, m, dt);
}
);
} else {
amrex::Abort("Unknown particle pusher");
- };
+ }
+
}
}
}
diff --git a/Source/Particles/WarpXParticleContainer.H b/Source/Particles/WarpXParticleContainer.H
index 7f86930b2..1a89ded02 100644
--- a/Source/Particles/WarpXParticleContainer.H
+++ b/Source/Particles/WarpXParticleContainer.H
@@ -322,9 +322,12 @@ protected:
#ifdef WARPX_QED
bool m_do_qed = false;
+ //Species for which QED effects are relevant should override these methods
virtual bool has_quantum_sync(){return false;};
virtual bool has_breit_wheeler(){return false;};
+ //Species can receive a shared pointer to a QED engine (species for
+ //which this is relevant should override these functions)
virtual void
set_breit_wheeler_engine_ptr(std::shared_ptr<BreitWheelerEngine>){};
virtual void
diff --git a/Source/QED/BreitWheelerEngineWrapper.H b/Source/QED/BreitWheelerEngineWrapper.H
index 1033ff7c9..3cdfdeae6 100644
--- a/Source/QED/BreitWheelerEngineWrapper.H
+++ b/Source/QED/BreitWheelerEngineWrapper.H
@@ -151,7 +151,7 @@ class BreitWheelerGeneratePairs
{
public:
/**
- * Constructor acquires a reference to control parameters and
+ * Constructor acquires pointers to control parameters and
* lookup tables data.
* lookup_table uses non-owning vectors under the hood. So no new data
* allocations should be triggered on GPU
@@ -288,6 +288,12 @@ public:
*/
PicsarBreitWheelerCtrl get_default_ctrl() const;
+ /**
+ * returns a constant reference to the control parameters
+ * @return const reference to control parameters
+ */
+ const PicsarBreitWheelerCtrl& get_ref_ctrl() const;
+
private:
bool m_lookup_tables_initialized = false;
diff --git a/Source/QED/BreitWheelerEngineWrapper.cpp b/Source/QED/BreitWheelerEngineWrapper.cpp
index 42953c97f..78ff13fc5 100644
--- a/Source/QED/BreitWheelerEngineWrapper.cpp
+++ b/Source/QED/BreitWheelerEngineWrapper.cpp
@@ -176,6 +176,12 @@ BreitWheelerEngine::get_default_ctrl() const
return PicsarBreitWheelerCtrl();
}
+const PicsarBreitWheelerCtrl&
+BreitWheelerEngine::get_ref_ctrl() const
+{
+ return m_innards.ctrl;
+}
+
void BreitWheelerEngine::compute_lookup_tables (
PicsarBreitWheelerCtrl ctrl)
{
diff --git a/Source/QED/QuantumSyncEngineWrapper.H b/Source/QED/QuantumSyncEngineWrapper.H
index 1a6ffe4f3..fd1571720 100644
--- a/Source/QED/QuantumSyncEngineWrapper.H
+++ b/Source/QED/QuantumSyncEngineWrapper.H
@@ -79,7 +79,7 @@ class QuantumSynchrotronEvolveOpticalDepth
{
public:
/**
- * Constructor acquires a reference to control parameters and
+ * Constructor acquires pointers to control parameters and
* lookup tables data.
* lookup_table uses non-owning vectors under the hood. So no new data
* allocations should be triggered on GPU
@@ -152,7 +152,7 @@ class QuantumSynchrotronGeneratePhotonAndUpdateMomentum
{
public:
/**
- * Constructor acquires a reference to control parameters and
+ * Constructor acquires pointers to control parameters and
* lookup tables data.
* lookup_table uses non-owning vectors under the hood. So no new data
* allocations should be triggered on GPU
@@ -285,6 +285,12 @@ public:
*/
PicsarQuantumSynchrotronCtrl get_default_ctrl() const;
+ /**
+ * returns a constant reference to the control parameters
+ * @return const reference to control parameters
+ */
+ const PicsarQuantumSynchrotronCtrl& get_ref_ctrl() const;
+
private:
bool m_lookup_tables_initialized = false;
diff --git a/Source/QED/QuantumSyncEngineWrapper.cpp b/Source/QED/QuantumSyncEngineWrapper.cpp
index b2630dc4d..c2d9f0abf 100644
--- a/Source/QED/QuantumSyncEngineWrapper.cpp
+++ b/Source/QED/QuantumSyncEngineWrapper.cpp
@@ -175,6 +175,12 @@ QuantumSynchrotronEngine::get_default_ctrl() const
return PicsarQuantumSynchrotronCtrl();
}
+const PicsarQuantumSynchrotronCtrl&
+QuantumSynchrotronEngine::get_ref_ctrl() const
+{
+ return m_innards.ctrl;
+}
+
void QuantumSynchrotronEngine::compute_lookup_tables (
PicsarQuantumSynchrotronCtrl ctrl)
{
generated by cgit v1.2.3 (git 2.46.0) at 2025-08-13 12:15:21 +0000