WarpX
QEDPhotonEmission.H
Go to the documentation of this file.
1 /* Copyright 2019 Luca Fedeli
2  *
3  * This file is part of WarpX.
4  *
5  * License: BSD-3-Clause-LBNL
6  */
7 
8 #ifndef WARPX_QED_PHOTON_EMISSION_H_
9 #define WARPX_QED_PHOTON_EMISSION_H_
10 
16 #include "Utils/WarpXConst.H"
17 
18 #include <AMReX_Array.H>
19 #include <AMReX_Array4.H>
20 #include <AMReX_Dim3.H>
21 #include <AMReX_Extension.H>
22 #include <AMReX_GpuLaunch.H>
23 #include <AMReX_GpuQualifiers.H>
24 #include <AMReX_IndexType.H>
25 #include <AMReX_ParticleTile.H>
26 #include <AMReX_REAL.H>
27 
28 #include <AMReX_BaseFwd.H>
29 
30 #include <algorithm>
31 #include <limits>
32 
44 {
45 public:
46 
52  PhotonEmissionFilterFunc(int const opt_depth_runtime_comp)
53  : m_opt_depth_runtime_comp(opt_depth_runtime_comp)
54  {}
55 
64  template <typename PData>
66  bool operator() (const PData& ptd, int const i, amrex::RandomEngine const&) const noexcept
67  {
68  using namespace amrex;
69 
70  const amrex::ParticleReal opt_depth =
71  ptd.m_runtime_rdata[m_opt_depth_runtime_comp][i];
72  return (opt_depth < 0.0_rt);
73  }
74 
75 private:
77 };
78 
83 {
84 
85 public:
86 
112  QuantumSynchrotronGetOpticalDepth opt_depth_functor,
113  int opt_depth_runtime_comp,
114  QuantumSynchrotronPhotonEmission emission_functor,
115  const WarpXParIter& a_pti, int lev, amrex::IntVect ngEB,
116  amrex::FArrayBox const& exfab,
117  amrex::FArrayBox const& eyfab,
118  amrex::FArrayBox const& ezfab,
119  amrex::FArrayBox const& bxfab,
120  amrex::FArrayBox const& byfab,
121  amrex::FArrayBox const& bzfab,
122  amrex::Vector<amrex::ParticleReal>& E_external_particle,
123  amrex::Vector<amrex::ParticleReal>& B_external_particle,
124  int a_offset = 0);
125 
136  template <typename DstData, typename SrcData>
138  void operator() (DstData& dst, SrcData& src, int i_src, int i_dst,
139  amrex::RandomEngine const& engine) const noexcept
140  {
141  using namespace amrex;
142 
143  // gather E and B
144  amrex::ParticleReal xp, yp, zp;
145  m_get_position(i_src, xp, yp, zp);
146 
147  amrex::ParticleReal ex = m_Ex_external_particle;
148  amrex::ParticleReal ey = m_Ey_external_particle;
149  amrex::ParticleReal ez = m_Ez_external_particle;
150  amrex::ParticleReal bx = m_Bx_external_particle;
151  amrex::ParticleReal by = m_By_external_particle;
152  amrex::ParticleReal bz = m_Bz_external_particle;
153 
154  m_get_externalEB(i_src, ex, ey, ez, bx, by, bz);
155 
156  doGatherShapeN(xp, yp, zp, ex, ey, ez, bx, by, bz,
157  m_ex_arr, m_ey_arr, m_ez_arr, m_bx_arr, m_by_arr, m_bz_arr,
158  m_ex_type, m_ey_type, m_ez_type, m_bx_type, m_by_type, m_bz_type,
159  m_dx_arr, m_xyzmin_arr, m_lo, m_n_rz_azimuthal_modes,
160  m_nox, m_galerkin_interpolation);
161 
162  auto& ux = src.m_rdata[PIdx::ux][i_src];
163  auto& uy = src.m_rdata[PIdx::uy][i_src];
164  auto& uz = src.m_rdata[PIdx::uz][i_src];
165  auto& g_ux = dst.m_rdata[PIdx::ux][i_dst];
166  auto& g_uy = dst.m_rdata[PIdx::uy][i_dst];
167  auto& g_uz = dst.m_rdata[PIdx::uz][i_dst];
168  m_emission_functor(
169  ux, uy, uz,
170  ex, ey, ez,
171  bx, by, bz,
172  g_ux, g_uy, g_uz,
173  engine);
174 
175  //Initialize the optical depth component of the source species.
176  src.m_runtime_rdata[m_opt_depth_runtime_comp][i_src] =
177  m_opt_depth_functor(engine);
178  }
179 
180 private:
184  const int m_opt_depth_runtime_comp = 0;
191  amrex::ParticleReal m_Ex_external_particle;
192  amrex::ParticleReal m_Ey_external_particle;
193  amrex::ParticleReal m_Ez_external_particle;
194  amrex::ParticleReal m_Bx_external_particle;
195  amrex::ParticleReal m_By_external_particle;
196  amrex::ParticleReal m_Bz_external_particle;
197 
204 
211 
214 
216  int m_nox;
218 
220 };
221 
222 
235 template <typename PTile>
237  PTile& ptile,
238  const int old_size, const int num_added,
239  const amrex::ParticleReal energy_threshold)
240 {
241  auto& soa = ptile.GetStructOfArrays();
242  auto p_idcpu = soa.GetIdCPUData().data() + old_size;
243  const auto p_ux = soa.GetRealData(PIdx::ux).data() + old_size;
244  const auto p_uy = soa.GetRealData(PIdx::uy).data() + old_size;
245  const auto p_uz = soa.GetRealData(PIdx::uz).data() + old_size;
246 
247  //The square of the energy threshold
248  const auto energy_threshold2 = std::max(
249  energy_threshold*energy_threshold,
250  std::numeric_limits<amrex::ParticleReal>::min());
251 
252  amrex::ParallelFor(num_added, [=] AMREX_GPU_DEVICE (int ip) noexcept
253  {
254  const auto ux = p_ux[ip];
255  const auto uy = p_uy[ip];
256  const auto uz = p_uz[ip];
257 
258  //The square of the photon energy (in SI units)
259  // ( Particle momentum is stored as gamma * velocity.)
260  constexpr amrex::ParticleReal me_c = PhysConst::m_e*PhysConst::c;
261  const auto phot_energy2 = (ux*ux + uy*uy + uz*uz)*me_c*me_c;
262 
263  if (phot_energy2 < energy_threshold2) {
264  p_idcpu[ip] = amrex::ParticleIdCpus::Invalid;
265  }
266  });
267 }
268 
269 
270 #endif //WARPX_QED_PHOTON_EMISSION_H_
#define AMREX_FORCE_INLINE
#define AMREX_GPU_DEVICE
#define AMREX_GPU_HOST_DEVICE
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void doGatherShapeN(const amrex::ParticleReal xp, const amrex::ParticleReal yp, const amrex::ParticleReal zp, amrex::ParticleReal &Exp, amrex::ParticleReal &Eyp, amrex::ParticleReal &Ezp, amrex::ParticleReal &Bxp, amrex::ParticleReal &Byp, amrex::ParticleReal &Bzp, amrex::Array4< amrex::Real const > const &ex_arr, amrex::Array4< amrex::Real const > const &ey_arr, amrex::Array4< amrex::Real const > const &ez_arr, amrex::Array4< amrex::Real const > const &bx_arr, amrex::Array4< amrex::Real const > const &by_arr, amrex::Array4< amrex::Real const > const &bz_arr, const amrex::IndexType ex_type, const amrex::IndexType ey_type, const amrex::IndexType ez_type, const amrex::IndexType bx_type, const amrex::IndexType by_type, const amrex::IndexType bz_type, const amrex::GpuArray< amrex::Real, 3 > &dx, const amrex::GpuArray< amrex::Real, 3 > &xyzmin, const amrex::Dim3 &lo, const int n_rz_azimuthal_modes)
Field gather for a single particle.
Definition: FieldGather.H:38
void cleanLowEnergyPhotons(PTile &ptile, const int old_size, const int num_added, const amrex::ParticleReal energy_threshold)
Free function to call to remove immediately low energy photons by setting their ID to -1....
Definition: QEDPhotonEmission.H:236
Filter functor for the QED photon emission process.
Definition: QEDPhotonEmission.H:44
int m_opt_depth_runtime_comp
Definition: QEDPhotonEmission.H:76
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool operator()(const PData &ptd, int const i, amrex::RandomEngine const &) const noexcept
Functor call. This method determines if a given (electron or positron) particle should undergo QED ph...
Definition: QEDPhotonEmission.H:66
PhotonEmissionFilterFunc(int const opt_depth_runtime_comp)
Constructor of the PhotonEmissionFilterFunc functor.
Definition: QEDPhotonEmission.H:52
Transform functor for the QED photon emission process.
Definition: QEDPhotonEmission.H:83
GetParticlePosition< PIdx > m_get_position
Definition: QEDPhotonEmission.H:189
amrex::ParticleReal m_Ey_external_particle
Definition: QEDPhotonEmission.H:192
amrex::Array4< const amrex::Real > m_ey_arr
Definition: QEDPhotonEmission.H:199
amrex::GpuArray< amrex::Real, 3 > m_xyzmin_arr
Definition: QEDPhotonEmission.H:213
amrex::ParticleReal m_By_external_particle
Definition: QEDPhotonEmission.H:195
bool m_galerkin_interpolation
Definition: QEDPhotonEmission.H:215
amrex::IndexType m_bz_type
Definition: QEDPhotonEmission.H:210
int m_nox
Definition: QEDPhotonEmission.H:216
amrex::ParticleReal m_Ez_external_particle
Definition: QEDPhotonEmission.H:193
const QuantumSynchrotronGetOpticalDepth m_opt_depth_functor
Definition: QEDPhotonEmission.H:182
amrex::IndexType m_by_type
Definition: QEDPhotonEmission.H:209
amrex::Array4< const amrex::Real > m_ex_arr
Definition: QEDPhotonEmission.H:198
amrex::IndexType m_ex_type
Definition: QEDPhotonEmission.H:205
amrex::IndexType m_bx_type
Definition: QEDPhotonEmission.H:208
amrex::IndexType m_ey_type
Definition: QEDPhotonEmission.H:206
amrex::IndexType m_ez_type
Definition: QEDPhotonEmission.H:207
amrex::Array4< const amrex::Real > m_bx_arr
Definition: QEDPhotonEmission.H:201
amrex::GpuArray< amrex::Real, 3 > m_dx_arr
Definition: QEDPhotonEmission.H:212
int m_n_rz_azimuthal_modes
Definition: QEDPhotonEmission.H:217
amrex::Array4< const amrex::Real > m_ez_arr
Definition: QEDPhotonEmission.H:200
amrex::ParticleReal m_Bz_external_particle
Definition: QEDPhotonEmission.H:196
amrex::ParticleReal m_Ex_external_particle
Definition: QEDPhotonEmission.H:191
amrex::Array4< const amrex::Real > m_by_arr
Definition: QEDPhotonEmission.H:202
const QuantumSynchrotronPhotonEmission m_emission_functor
Definition: QEDPhotonEmission.H:187
GetExternalEBField m_get_externalEB
Definition: QEDPhotonEmission.H:190
amrex::ParticleReal m_Bx_external_particle
Definition: QEDPhotonEmission.H:194
amrex::Dim3 m_lo
Definition: QEDPhotonEmission.H:219
amrex::Array4< const amrex::Real > m_bz_arr
Definition: QEDPhotonEmission.H:203
Definition: QuantumSyncEngineWrapper.H:76
Definition: QuantumSyncEngineWrapper.H:179
Definition: WarpXParticleContainer.H:53
static constexpr auto c
vacuum speed of light [m/s]
Definition: constant.H:44
static constexpr auto m_e
electron mass [kg]
Definition: constant.H:52
constexpr std::uint64_t Invalid
std::enable_if_t< std::is_integral_v< T > > ParallelFor(TypeList< CTOs... >, std::array< int, sizeof...(CTOs)> const &runtime_options, T N, F &&f)
i
Definition: check_interp_points_and_weights.py:174
Functor class that assigns external field values (E and B) to particles.
Definition: GetExternalFields.H:25
@ uz
Definition: NamedComponentParticleContainer.H:34
@ uy
Definition: NamedComponentParticleContainer.H:34
@ ux
Definition: NamedComponentParticleContainer.H:34