WarpX
MultiParticleContainer.H
Go to the documentation of this file.
1 /* Copyright 2019-2020 Andrew Myers, Ann Almgren, Axel Huebl
2  * David Grote, Jean-Luc Vay, Junmin Gu
3  * Luca Fedeli, Mathieu Lobet, Maxence Thevenet
4  * Remi Lehe, Revathi Jambunathan, Weiqun Zhang
5  * Yinjian Zhao
6  *
7  * This file is part of WarpX.
8  *
9  * License: BSD-3-Clause-LBNL
10  */
11 #ifndef WARPX_ParticleContainer_H_
12 #define WARPX_ParticleContainer_H_
13 
15 
16 #include "Evolve/WarpXDtType.H"
17 #include "Evolve/WarpXPushType.H"
19 #ifdef WARPX_QED
22 #endif
24 #include "Utils/TextMsg.H"
25 #include "Utils/WarpXConst.H"
26 #include "WarpXParticleContainer.H"
27 #include "ParticleBoundaries.H"
28 
29 #include <AMReX_BLassert.H>
30 #include <AMReX_Box.H>
31 #include <AMReX_Config.H>
32 #include <AMReX_GpuControl.H>
33 #include <AMReX_INT.H>
34 #include <AMReX_MFIter.H>
35 #include <AMReX_REAL.H>
36 #include <AMReX_RealBox.H>
37 #include <AMReX_Vector.H>
38 
39 #include <AMReX_BaseFwd.H>
40 #include <AMReX_AmrCoreFwd.H>
41 
42 #include <algorithm>
43 #include <array>
44 #include <iosfwd>
45 #include <iterator>
46 #include <limits>
47 #include <memory>
48 #include <string>
49 #include <vector>
50 
66 {
67 
68 public:
69 
71 
73 
78 
79  [[nodiscard]] WarpXParticleContainer&
81 
82  [[nodiscard]] WarpXParticleContainer*
83  GetParticleContainerPtr (int index) const {return allcontainers[index].get();}
84 
85  [[nodiscard]] WarpXParticleContainer&
86  GetParticleContainerFromName (const std::string& name) const;
87 
88  std::array<amrex::ParticleReal, 3> meanParticleVelocity(int index) {
89  return allcontainers[index]->meanParticleVelocity();
90  }
91 
92  void AllocData ();
93 
94  void InitData ();
95 
97 
103  void Evolve (int lev,
104  const amrex::MultiFab& Ex, const amrex::MultiFab& Ey, const amrex::MultiFab& Ez,
105  const amrex::MultiFab& Bx, const amrex::MultiFab& By, const amrex::MultiFab& Bz,
108  amrex::MultiFab* rho, amrex::MultiFab* crho,
109  const amrex::MultiFab* cEx, const amrex::MultiFab* cEy, const amrex::MultiFab* cEz,
110  const amrex::MultiFab* cBx, const amrex::MultiFab* cBy, const amrex::MultiFab* cBz,
111  amrex::Real t, amrex::Real dt, DtType a_dt_type=DtType::Full, bool skip_deposition=false,
112  PushType push_type=PushType::Explicit);
113 
118  void PushX (amrex::Real dt);
119 
126  void PushP (int lev, amrex::Real dt,
127  const amrex::MultiFab& Ex, const amrex::MultiFab& Ey, const amrex::MultiFab& Ez,
128  const amrex::MultiFab& Bx, const amrex::MultiFab& By, const amrex::MultiFab& Bz);
129 
136  std::unique_ptr<amrex::MultiFab> GetZeroChargeDensity(int lev);
137 
147  void
148  DepositCharge (amrex::Vector<std::unique_ptr<amrex::MultiFab> >& rho,
149  amrex::Real relative_time);
150 
162  void
163  DepositCurrent (amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3 > >& J,
164  amrex::Real dt, amrex::Real relative_time);
165 
169 
170  std::unique_ptr<amrex::MultiFab> GetChargeDensity(int lev, bool local = false);
171 
172  void doFieldIonization (int lev,
173  const amrex::MultiFab& Ex, const amrex::MultiFab& Ey, const amrex::MultiFab& Ez,
174  const amrex::MultiFab& Bx, const amrex::MultiFab& By, const amrex::MultiFab& Bz);
175 
176  void doCollisions (amrex::Real cur_time, amrex::Real dt);
177 
183  void doResampling (int timestep, bool verbose);
184 
185 #ifdef WARPX_QED
193  void doQEDSchwinger ();
194 
199  [[nodiscard]] amrex::Box ComputeSchwingerGlobalBox () const;
200 #endif
201 
202  void Restart (const std::string& dir);
203 
204  void PostRestart ();
205 
206  void ReadHeader (std::istream& is);
207 
208  void WriteHeader (std::ostream& os) const;
209 
210  void SortParticlesByBin (amrex::IntVect bin_size);
211 
212  void Redistribute ();
213 
214  void defineAllParticleTiles ();
215 
216  void deleteInvalidParticles ();
217 
218  void RedistributeLocal (int num_ghost);
219 
222  void ApplyBoundaryConditions ();
223 
231  [[nodiscard]] amrex::Vector<amrex::Long> GetZeroParticlesInGrid(int lev) const;
232 
233  [[nodiscard]] amrex::Vector<amrex::Long> NumberOfParticlesInGrid(int lev) const;
234 
235  void Increment (amrex::MultiFab& mf, int lev);
236 
237  void SetParticleBoxArray (int lev, amrex::BoxArray& new_ba);
239 
240  [[nodiscard]] int nSpecies () const {return static_cast<int>(species_names.size());}
241  [[nodiscard]] int nLasers () const {return static_cast<int>(lasers_names.size());}
242  [[nodiscard]] int nContainers () const {return static_cast<int>(allcontainers.size());}
243 
248  void SetDoBackTransformedParticles (bool do_back_transformed_particles);
254  void SetDoBackTransformedParticles (const std::string& species_name, bool do_back_transformed_particles);
255 
256  [[nodiscard]] int nSpeciesDepositOnMainGrid () const
257  {
258  bool const onMainGrid = true;
259  auto const & v = m_deposit_on_main_grid;
260  return static_cast<int>(std::count( v.begin(), v.end(), onMainGrid ));
261  }
262 
263  [[nodiscard]] int nSpeciesGatherFromMainGrid() const
264  {
265  bool const fromMainGrid = true;
266  auto const & v = m_gather_from_main_grid;
267  return static_cast<int>(std::count( v.begin(), v.end(), fromMainGrid ));
268  }
269 
270  // Inject particles during the simulation (for particles entering the
271  // simulation domain after some iterations, due to flowing plasma and/or
272  // moving window).
273  void ContinuousInjection(const amrex::RealBox& injection_box) const;
274 
279  void UpdateAntennaPosition(amrex::Real dt) const;
280 
281  [[nodiscard]] int doContinuousInjection() const;
282 
283  // Inject particles from a surface during the simulation
284  void ContinuousFluxInjection(amrex::Real t, amrex::Real dt) const;
285 
286  [[nodiscard]] std::vector<std::string> GetSpeciesNames() const { return species_names; }
287 
288  [[nodiscard]] std::vector<std::string> GetLasersNames() const { return lasers_names; }
289 
290  [[nodiscard]] std::vector<std::string> GetSpeciesAndLasersNames() const
291  {
292  std::vector<std::string> tmp = species_names;
293  tmp.insert(tmp.end(), lasers_names.begin(), lasers_names.end());
294  return tmp;
295  }
296 
298 
299  void ScrapeParticlesAtEB (const amrex::Vector<const amrex::MultiFab*>& distance_to_eb);
300 
301  std::string m_B_ext_particle_s = "none";
302  std::string m_E_ext_particle_s = "none";
303  // Parser for B_external on the particle
304  std::unique_ptr<amrex::Parser> m_Bx_particle_parser;
305  std::unique_ptr<amrex::Parser> m_By_particle_parser;
306  std::unique_ptr<amrex::Parser> m_Bz_particle_parser;
307  // Parser for E_external on the particle
308  std::unique_ptr<amrex::Parser> m_Ex_particle_parser;
309  std::unique_ptr<amrex::Parser> m_Ey_particle_parser;
310  std::unique_ptr<amrex::Parser> m_Ez_particle_parser;
311 
312  amrex::ParticleReal m_repeated_plasma_lens_period;
321 
322 #ifdef WARPX_QED
326  void doQedEvents (int lev,
327  const amrex::MultiFab& Ex,
328  const amrex::MultiFab& Ey,
329  const amrex::MultiFab& Ez,
330  const amrex::MultiFab& Bx,
331  const amrex::MultiFab& By,
332  const amrex::MultiFab& Bz);
333 #endif
334 
335  [[nodiscard]] int getSpeciesID (const std::string& product_str) const;
336 
339 
340 protected:
341 
342 #ifdef WARPX_QED
346  void doQedBreitWheeler (int lev,
347  const amrex::MultiFab& Ex,
348  const amrex::MultiFab& Ey,
349  const amrex::MultiFab& Ez,
350  const amrex::MultiFab& Bx,
351  const amrex::MultiFab& By,
352  const amrex::MultiFab& Bz);
353 
357  void doQedQuantumSync (int lev,
358  const amrex::MultiFab& Ex,
359  const amrex::MultiFab& Ey,
360  const amrex::MultiFab& Ez,
361  const amrex::MultiFab& Bx,
362  const amrex::MultiFab& By,
363  const amrex::MultiFab& Bz);
364 #endif
365 
366  // Particle container types
368 
369  std::vector<std::string> species_names;
370 
371  std::vector<std::string> lasers_names;
372 
373  std::unique_ptr<CollisionHandler> collisionhandler;
374 
376  std::vector<bool> m_deposit_on_main_grid;
377  std::vector<bool> m_laser_deposit_on_main_grid;
378 
380  std::vector<bool> m_gather_from_main_grid;
381 
382  std::vector<PCTypes> species_types;
383 
384  template<typename ...Args>
385  [[nodiscard]]
387  Args const&... pc_dsts) const noexcept
388  {
389  amrex::MFItInfo info;
390 
391  MFItInfoCheckTiling(pc_src, pc_dsts...);
392 
395  }
396 
397 #ifdef AMREX_USE_OMP
398  info.SetDynamic(true);
399 #endif
400 
401  return info;
402  }
403 
404 
405 #ifdef WARPX_QED
406  // The QED engines
407  std::shared_ptr<BreitWheelerEngine> m_shr_p_bw_engine;
408  std::shared_ptr<QuantumSynchrotronEngine> m_shr_p_qs_engine;
409  //_______________________________
410 
415  void InitQED ();
416 
417  //Variables to store how many species need a QED process
420  //________
421 
423  static_cast<amrex::ParticleReal>(
433  [[nodiscard]] int NSpeciesQuantumSync() const { return m_nspecies_quantum_sync;}
434 
438  [[nodiscard]] int NSpeciesBreitWheeler() const { return m_nspecies_breit_wheeler;}
439 
443  void InitQuantumSync ();
444 
448  void InitBreitWheeler ();
449 
455 
461 
463  bool m_do_qed_schwinger = false;
482  amrex::Real m_qed_schwinger_xmin = std::numeric_limits<amrex::Real>::lowest();
483  amrex::Real m_qed_schwinger_xmax = std::numeric_limits<amrex::Real>::max();
484  amrex::Real m_qed_schwinger_ymin = std::numeric_limits<amrex::Real>::lowest();
485  amrex::Real m_qed_schwinger_ymax = std::numeric_limits<amrex::Real>::max();
486  amrex::Real m_qed_schwinger_zmin = std::numeric_limits<amrex::Real>::lowest();
487  amrex::Real m_qed_schwinger_zmax = std::numeric_limits<amrex::Real>::max();
488 
489 #endif
490 
491 private:
492 
493  // physical particles (+ laser)
495  // Temporary particle container, used e.g. for particle splitting.
496  std::unique_ptr<PhysicalParticleContainer> pc_tmp;
497 
498  void ReadParameters ();
499 
500  void mapSpeciesProduct ();
501 
503 
504  void MFItInfoCheckTiling(const WarpXParticleContainer& /*pc_src*/) const noexcept
505  {}
506 
507  template<typename First, typename ...Args>
509  First const& pc_dst, Args const&... others) const noexcept
510  {
512  WARPX_ALWAYS_ASSERT_WITH_MESSAGE(pc_dst.do_tiling,
513  "For particle creation processes, either all or none of the "
514  "particle species must use tiling.");
515  }
516 
517  MFItInfoCheckTiling(pc_src, others...);
518  }
519 
526 
527 #ifdef WARPX_QED
533  void CheckQEDProductSpecies();
534 #endif
535 
536 
537 };
538 #endif /*WARPX_ParticleContainer_H_*/
#define WARPX_ALWAYS_ASSERT_WITH_MESSAGE(EX, MSG)
Definition: TextMsg.H:13
DtType
Definition: WarpXDtType.H:11
PushType
Definition: WarpXPushType.H:12
Definition: MultiParticleContainer.H:66
void RedistributeLocal(int num_ghost)
Definition: MultiParticleContainer.cpp:647
void deleteInvalidParticles()
Definition: MultiParticleContainer.cpp:639
std::unique_ptr< amrex::Parser > m_Bz_particle_parser
Definition: MultiParticleContainer.H:306
void InitQED()
Definition: MultiParticleContainer.cpp:971
std::unique_ptr< amrex::Parser > m_Ey_particle_parser
Definition: MultiParticleContainer.H:309
std::unique_ptr< amrex::MultiFab > GetChargeDensity(int lev, bool local=false)
Definition: MultiParticleContainer.cpp:589
int nLasers() const
Definition: MultiParticleContainer.H:241
amrex::Real m_qed_schwinger_zmax
Definition: MultiParticleContainer.H:487
amrex::Real m_qed_schwinger_zmin
Definition: MultiParticleContainer.H:486
void InitMultiPhysicsModules()
Definition: MultiParticleContainer.cpp:435
std::string m_B_ext_particle_s
Definition: MultiParticleContainer.H:301
amrex::Vector< std::unique_ptr< WarpXParticleContainer > >::iterator begin()
Definition: MultiParticleContainer.H:337
void QuantumSyncGenerateTable()
Definition: MultiParticleContainer.cpp:1134
amrex::Vector< amrex::ParticleReal > h_repeated_plasma_lens_strengths_E
Definition: MultiParticleContainer.H:315
void AllocData()
Definition: MultiParticleContainer.cpp:403
void SetParticleDistributionMap(int lev, amrex::DistributionMapping &new_dm)
Definition: MultiParticleContainer.cpp:711
void InitData()
Definition: MultiParticleContainer.cpp:412
amrex::Real m_qed_schwinger_ymin
Definition: MultiParticleContainer.H:484
std::vector< bool > m_laser_deposit_on_main_grid
Definition: MultiParticleContainer.H:377
void Increment(amrex::MultiFab &mf, int lev)
Definition: MultiParticleContainer.cpp:695
void mapSpeciesProduct()
Definition: MultiParticleContainer.cpp:772
void ContinuousInjection(const amrex::RealBox &injection_box) const
Definition: MultiParticleContainer.cpp:724
amrex::Gpu::DeviceVector< amrex::ParticleReal > d_repeated_plasma_lens_lengths
Definition: MultiParticleContainer.H:318
MultiParticleContainer(MultiParticleContainer const &)=delete
int NSpeciesBreitWheeler() const
Definition: MultiParticleContainer.H:438
void ReadHeader(std::istream &is)
Definition: ParticleIO.cpp:220
static constexpr auto m_default_quantum_sync_photon_creation_energy_threshold
Definition: MultiParticleContainer.H:422
void ApplyBoundaryConditions()
Definition: MultiParticleContainer.cpp:655
std::vector< std::string > lasers_names
Definition: MultiParticleContainer.H:371
amrex::Gpu::DeviceVector< amrex::ParticleReal > d_repeated_plasma_lens_starts
Definition: MultiParticleContainer.H:317
bool m_do_qed_schwinger
Definition: MultiParticleContainer.H:463
std::vector< std::string > GetLasersNames() const
Definition: MultiParticleContainer.H:288
std::unique_ptr< amrex::Parser > m_By_particle_parser
Definition: MultiParticleContainer.H:305
void PushX(amrex::Real dt)
This pushes the particle positions by one time step for all the species in the MultiParticleContainer...
Definition: MultiParticleContainer.cpp:481
amrex::Vector< amrex::ParticleReal > h_repeated_plasma_lens_lengths
Definition: MultiParticleContainer.H:314
std::string m_qed_schwinger_ele_product_name
Definition: MultiParticleContainer.H:465
int m_nspecies_quantum_sync
Definition: MultiParticleContainer.H:418
int m_qed_schwinger_ele_product
Definition: MultiParticleContainer.H:469
int nSpecies() const
Definition: MultiParticleContainer.H:240
int m_nspecies_breit_wheeler
Definition: MultiParticleContainer.H:419
void MFItInfoCheckTiling(const WarpXParticleContainer &) const noexcept
Definition: MultiParticleContainer.H:504
void InitQuantumSync()
Definition: MultiParticleContainer.cpp:1002
void DepositCurrent(amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 > > &J, amrex::Real dt, amrex::Real relative_time)
Deposit current density.
Definition: MultiParticleContainer.cpp:524
std::shared_ptr< BreitWheelerEngine > m_shr_p_bw_engine
Definition: MultiParticleContainer.H:407
void Redistribute()
Definition: MultiParticleContainer.cpp:623
std::vector< PCTypes > species_types
Definition: MultiParticleContainer.H:382
void PostRestart()
Definition: MultiParticleContainer.cpp:424
int NSpeciesQuantumSync() const
Definition: MultiParticleContainer.H:433
amrex::ParticleReal m_quantum_sync_photon_creation_energy_threshold
Definition: MultiParticleContainer.H:427
WarpXParticleContainer & GetParticleContainerFromName(const std::string &name) const
Definition: MultiParticleContainer.cpp:392
MultiParticleContainer & operator=(MultiParticleContainer &&)=default
std::vector< std::string > species_names
Definition: MultiParticleContainer.H:369
amrex::Box ComputeSchwingerGlobalBox() const
Definition: MultiParticleContainer.cpp:1416
amrex::Gpu::DeviceVector< amrex::ParticleReal > d_repeated_plasma_lens_strengths_B
Definition: MultiParticleContainer.H:320
amrex::Real m_qed_schwinger_ymax
Definition: MultiParticleContainer.H:485
void doQEDSchwinger()
Definition: MultiParticleContainer.cpp:1309
std::string m_E_ext_particle_s
Definition: MultiParticleContainer.H:302
amrex::Real m_qed_schwinger_xmin
Definition: MultiParticleContainer.H:482
int nSpeciesGatherFromMainGrid() const
Definition: MultiParticleContainer.H:263
void defineAllParticleTiles()
Definition: MultiParticleContainer.cpp:631
void BreitWheelerGenerateTable()
Definition: MultiParticleContainer.cpp:1224
std::unique_ptr< PhysicalParticleContainer > pc_tmp
Definition: MultiParticleContainer.H:496
std::unique_ptr< amrex::MultiFab > GetZeroChargeDensity(int lev)
This returns a MultiFAB filled with zeros. It is used to return the charge density when there is no p...
Definition: MultiParticleContainer.cpp:499
void ReadParameters()
Definition: MultiParticleContainer.cpp:132
int doContinuousInjection() const
Definition: MultiParticleContainer.cpp:744
amrex::Vector< amrex::ParticleReal > h_repeated_plasma_lens_strengths_B
Definition: MultiParticleContainer.H:316
int nSpeciesDepositOnMainGrid() const
Definition: MultiParticleContainer.H:256
void DepositCharge(amrex::Vector< std::unique_ptr< amrex::MultiFab > > &rho, amrex::Real relative_time)
Deposit charge density.
Definition: MultiParticleContainer.cpp:552
std::unique_ptr< CollisionHandler > collisionhandler
Definition: MultiParticleContainer.H:373
amrex::ParticleReal m_repeated_plasma_lens_period
Definition: MultiParticleContainer.H:312
amrex::Real m_qed_schwinger_xmax
Definition: MultiParticleContainer.H:483
WarpXParticleContainer * GetParticleContainerPtr(int index) const
Definition: MultiParticleContainer.H:83
void PushP(int lev, amrex::Real dt, const amrex::MultiFab &Ex, const amrex::MultiFab &Ey, const amrex::MultiFab &Ez, const amrex::MultiFab &Bx, const amrex::MultiFab &By, const amrex::MultiFab &Bz)
Definition: MultiParticleContainer.cpp:489
void MFItInfoCheckTiling(const WarpXParticleContainer &pc_src, First const &pc_dst, Args const &... others) const noexcept
Definition: MultiParticleContainer.H:508
int getSpeciesID(const std::string &product_str) const
Definition: MultiParticleContainer.cpp:817
std::vector< std::string > GetSpeciesAndLasersNames() const
Definition: MultiParticleContainer.H:290
amrex::Vector< amrex::Long > NumberOfParticlesInGrid(int lev) const
Definition: MultiParticleContainer.cpp:672
std::vector< std::string > GetSpeciesNames() const
Definition: MultiParticleContainer.H:286
void InitBreitWheeler()
Definition: MultiParticleContainer.cpp:1074
void doFieldIonization(int lev, const amrex::MultiFab &Ex, const amrex::MultiFab &Ey, const amrex::MultiFab &Ez, const amrex::MultiFab &Bx, const amrex::MultiFab &By, const amrex::MultiFab &Bz)
Definition: MultiParticleContainer.cpp:865
amrex::Vector< amrex::Long > GetZeroParticlesInGrid(int lev) const
This returns a vector filled with zeros whose size is the number of boxes in the simulation boxarray....
Definition: MultiParticleContainer.cpp:663
bool m_do_back_transformed_particles
Definition: MultiParticleContainer.H:502
void SortParticlesByBin(amrex::IntVect bin_size)
Definition: MultiParticleContainer.cpp:611
amrex::Vector< std::unique_ptr< WarpXParticleContainer > > allcontainers
Definition: MultiParticleContainer.H:494
void SetParticleBoxArray(int lev, amrex::BoxArray &new_ba)
Definition: MultiParticleContainer.cpp:703
void doQedQuantumSync(int lev, const amrex::MultiFab &Ex, const amrex::MultiFab &Ey, const amrex::MultiFab &Ez, const amrex::MultiFab &Bx, const amrex::MultiFab &By, const amrex::MultiFab &Bz)
Performs QED photon emission for the species for which it is enabled.
Definition: MultiParticleContainer.cpp:1575
void ScrapeParticlesAtEB(const amrex::Vector< const amrex::MultiFab * > &distance_to_eb)
Definition: MultiParticleContainer.cpp:959
void UpdateAntennaPosition(amrex::Real dt) const
Update antenna position for continuous injection of lasers in a boosted frame. Empty function for con...
Definition: MultiParticleContainer.cpp:734
MultiParticleContainer(MultiParticleContainer &&)=default
void doResampling(int timestep, bool verbose)
This function loops over all species and performs resampling if appropriate.
Definition: MultiParticleContainer.cpp:937
void CheckIonizationProductSpecies()
Definition: MultiParticleContainer.cpp:948
std::unique_ptr< amrex::Parser > m_Ex_particle_parser
Definition: MultiParticleContainer.H:308
std::array< amrex::ParticleReal, 3 > meanParticleVelocity(int index)
Definition: MultiParticleContainer.H:88
amrex::Vector< amrex::ParticleReal > h_repeated_plasma_lens_starts
Definition: MultiParticleContainer.H:313
int m_qed_schwinger_threshold_poisson_gaussian
Definition: MultiParticleContainer.H:478
std::vector< bool > m_gather_from_main_grid
instead of gathering fields from the finest patch level, gather from the coarsest
Definition: MultiParticleContainer.H:380
void doQedEvents(int lev, const amrex::MultiFab &Ex, const amrex::MultiFab &Ey, const amrex::MultiFab &Ez, const amrex::MultiFab &Bx, const amrex::MultiFab &By, const amrex::MultiFab &Bz)
Performs QED events (Breit-Wheeler process and photon emission)
Definition: MultiParticleContainer.cpp:1478
void CheckQEDProductSpecies()
Definition: MultiParticleContainer.cpp:1656
MultiParticleContainer & operator=(MultiParticleContainer const &)=delete
amrex::Real m_qed_schwinger_y_size
Definition: MultiParticleContainer.H:473
amrex::Gpu::DeviceVector< amrex::ParticleReal > d_repeated_plasma_lens_strengths_E
Definition: MultiParticleContainer.H:319
PhysicalParticleContainer & GetPCtmp()
Definition: MultiParticleContainer.H:297
int nContainers() const
Definition: MultiParticleContainer.H:242
void WriteHeader(std::ostream &os) const
Definition: ParticleIO.cpp:231
std::string m_qed_schwinger_pos_product_name
Definition: MultiParticleContainer.H:467
void Restart(const std::string &dir)
Definition: ParticleIO.cpp:124
void ContinuousFluxInjection(amrex::Real t, amrex::Real dt) const
Definition: MultiParticleContainer.cpp:760
int m_qed_schwinger_pos_product
Definition: MultiParticleContainer.H:471
amrex::MFItInfo getMFItInfo(const WarpXParticleContainer &pc_src, Args const &... pc_dsts) const noexcept
Definition: MultiParticleContainer.H:386
MultiParticleContainer(amrex::AmrCore *amr_core)
Definition: MultiParticleContainer.cpp:96
WarpXParticleContainer & GetParticleContainer(int index) const
Definition: MultiParticleContainer.H:80
void doCollisions(amrex::Real cur_time, amrex::Real dt)
Definition: MultiParticleContainer.cpp:931
void Evolve(int lev, const amrex::MultiFab &Ex, const amrex::MultiFab &Ey, const amrex::MultiFab &Ez, const amrex::MultiFab &Bx, const amrex::MultiFab &By, const amrex::MultiFab &Bz, amrex::MultiFab &jx, amrex::MultiFab &jy, amrex::MultiFab &jz, amrex::MultiFab *cjx, amrex::MultiFab *cjy, amrex::MultiFab *cjz, amrex::MultiFab *rho, amrex::MultiFab *crho, const amrex::MultiFab *cEx, const amrex::MultiFab *cEy, const amrex::MultiFab *cEz, const amrex::MultiFab *cBx, const amrex::MultiFab *cBy, const amrex::MultiFab *cBz, amrex::Real t, amrex::Real dt, DtType a_dt_type=DtType::Full, bool skip_deposition=false, PushType push_type=PushType::Explicit)
This evolves all the particles by one PIC time step, including current deposition,...
Definition: MultiParticleContainer.cpp:453
std::unique_ptr< amrex::Parser > m_Ez_particle_parser
Definition: MultiParticleContainer.H:310
std::shared_ptr< QuantumSynchrotronEngine > m_shr_p_qs_engine
Definition: MultiParticleContainer.H:408
~MultiParticleContainer()=default
std::vector< bool > m_deposit_on_main_grid
instead of depositing (current, charge) on the finest patch level, deposit to the coarsest grid
Definition: MultiParticleContainer.H:376
void SetDoBackTransformedParticles(bool do_back_transformed_particles)
Definition: MultiParticleContainer.cpp:841
PCTypes
Definition: MultiParticleContainer.H:367
amrex::Vector< std::unique_ptr< WarpXParticleContainer > >::iterator end()
Definition: MultiParticleContainer.H:338
std::unique_ptr< amrex::Parser > m_Bx_particle_parser
Definition: MultiParticleContainer.H:304
void doQedBreitWheeler(int lev, const amrex::MultiFab &Ex, const amrex::MultiFab &Ey, const amrex::MultiFab &Ez, const amrex::MultiFab &Bx, const amrex::MultiFab &By, const amrex::MultiFab &Bz)
Performs Breit-Wheeler process for the species for which it is enabled.
Definition: MultiParticleContainer.cpp:1492
Definition: PhysicalParticleContainer.H:47
Definition: WarpXParticleContainer.H:111
static AMREX_EXPORT bool do_tiling
static AMREX_EXPORT IntVect tile_size
Long size() const noexcept
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
bool notInLaunchRegion() noexcept
int count
Definition: run_alltests.py:322
name
Definition: run_automated.py:229
index
Definition: run_automated.py:328
float dt
Definition: stencil.py:442
MFItInfo & EnableTiling(const IntVect &ts=FabArrayBase::mfiter_tile_size) noexcept
MFItInfo & SetDynamic(bool f) noexcept