Hipace
MultiPlasma.H
Go to the documentation of this file.
1 /* Copyright 2021-2022
2  *
3  * This file is part of HiPACE++.
4  *
5  * Authors: AlexanderSinn, MaxThevenet, Severin Diederichs
6  * License: BSD-3-Clause-LBNL
7  */
8 #ifndef MULTIPLASMA_H_
9 #define MULTIPLASMA_H_
10 
13 #include "fields/Fields.H"
14 #include "laser/MultiLaser.H"
16 
18 {
19 
20 public:
21 
24  MultiPlasma ();
25 
28 
39 
40 
54  void DepositCurrent (
55  Fields & fields, const MultiLaser & multi_laser, int which_slice, bool deposit_jx_jy,
56  bool deposit_jz, bool deposit_rho, bool deposit_chi, bool deposit_rhomjz,
57  amrex::Vector<amrex::Geometry> const& gm, int const lev);
58 
66  void ExplicitDeposition (Fields& fields, const MultiLaser& multi_laser,
67  amrex::Vector<amrex::Geometry> const& gm, int const lev);
68 
74  amrex::Real maxDensity (amrex::Real z);
75 
84  void AdvanceParticles (
85  const Fields & fields, const MultiLaser & multi_laser, amrex::Vector<amrex::Geometry> const& gm,
86  bool temp_slice, int lev);
87 
97  Fields & fields, const MultiLaser & multi_laser, int which_slice,
98  amrex::Vector<amrex::Geometry> const& gm, int const lev);
99 
106  void DoFieldIonization (const int lev, const amrex::Geometry& geom, const Fields& fields);
107 
108  bool IonizationOn () const;
110  bool AnySpeciesNeutralizeBackground () const;
111 
117  void TileSort (amrex::Box bx, amrex::Geometry geom);
118 
120  const amrex::Vector<std::string>& GetNames() const {return m_names;}
121 
123  int GetNPlasmas() const {return m_nplasmas;}
124 
128  void ReorderParticles (const int islice);
129 
135  void TagByLevel (const int current_N_level, amrex::Vector<amrex::Geometry> const& geom3D,
136  const bool to_prev=false);
137 
145  void InSituComputeDiags (int step, int islice, int max_step,
146  amrex::Real physical_time, amrex::Real max_time);
154  void InSituWriteToFile (int step, amrex::Real time, const amrex::Geometry& geom,
155  int max_step, amrex::Real max_time);
156 
157  int m_sort_bin_size {32};
160  int m_nplasmas = 0;
164 private:
166  amrex::Real m_adaptive_density = 0.;
167 
168 };
169 
170 #endif // MULTIPLASMA_H_
Main class handling all field data structures and operations.
Definition: Fields.H:85
Definition: MultiLaser.H:99
Definition: MultiPlasma.H:18
const amrex::Vector< std::string > & GetNames() const
Definition: MultiPlasma.H:120
int GetNPlasmas() const
Definition: MultiPlasma.H:123
amrex::Vector< std::string > m_names
Definition: MultiPlasma.H:163
void TagByLevel(const int current_N_level, amrex::Vector< amrex::Geometry > const &geom3D, const bool to_prev=false)
Store the finest level of every plasma particle in the cpu() attribute.
Definition: MultiPlasma.cpp:172
bool AnySpeciesNeutralizeBackground() const
whether any plasma species uses a neutralizing background, e.g. no ion motion
Definition: MultiPlasma.cpp:144
amrex::Vector< PlasmaParticleContainer > m_all_plasmas
Definition: MultiPlasma.H:159
void InSituWriteToFile(int step, amrex::Real time, const amrex::Geometry &geom, int max_step, amrex::Real max_time)
Definition: MultiPlasma.cpp:193
void DepositNeutralizingBackground(Fields &fields, const MultiLaser &multi_laser, int which_slice, amrex::Vector< amrex::Geometry > const &gm, int const lev)
Loop over plasma species and deposit their neutralizing background, if needed.
Definition: MultiPlasma.cpp:111
void InSituComputeDiags(int step, int islice, int max_step, amrex::Real physical_time, amrex::Real max_time)
Definition: MultiPlasma.cpp:181
amrex::Vector< PlasmaBins > m_all_bins
Definition: MultiPlasma.H:162
void TileSort(amrex::Box bx, amrex::Geometry geom)
sort particles of all containers by tile logically, and store results in m_all_bins
Definition: MultiPlasma.cpp:154
void ReorderParticles(const int islice)
Definition: MultiPlasma.cpp:164
void AdvanceParticles(const Fields &fields, const MultiLaser &multi_laser, amrex::Vector< amrex::Geometry > const &gm, bool temp_slice, int lev)
Gather field values and push particles.
Definition: MultiPlasma.cpp:100
amrex::Real maxDensity(amrex::Real z)
Return max density, to compute the adaptive time step.
Definition: MultiPlasma.cpp:67
bool IonizationOn() const
Definition: MultiPlasma.cpp:134
void DoFieldIonization(const int lev, const amrex::Geometry &geom, const Fields &fields)
Definition: MultiPlasma.cpp:125
void InitData(amrex::Vector< amrex::BoxArray > slice_ba, amrex::Vector< amrex::DistributionMapping > slice_dm, amrex::Vector< amrex::Geometry > slice_gm, amrex::Vector< amrex::Geometry > gm)
Loop over plasma species and initialize them.
Definition: MultiPlasma.cpp:40
int m_nplasmas
Definition: MultiPlasma.H:160
amrex::Real m_adaptive_density
Definition: MultiPlasma.H:166
MultiPlasma()
Definition: MultiPlasma.cpp:18
int m_sort_bin_size
Definition: MultiPlasma.H:157
void ExplicitDeposition(Fields &fields, const MultiLaser &multi_laser, amrex::Vector< amrex::Geometry > const &gm, int const lev)
Definition: MultiPlasma.cpp:91
~MultiPlasma()
Definition: MultiPlasma.H:27
void DepositCurrent(Fields &fields, const MultiLaser &multi_laser, int which_slice, bool deposit_jx_jy, bool deposit_jz, bool deposit_rho, bool deposit_chi, bool deposit_rhomjz, amrex::Vector< amrex::Geometry > const &gm, int const lev)
Definition: MultiPlasma.cpp:78