Hipace
Public Member Functions | Public Attributes | Private Attributes | List of all members
MultiPlasma Class Reference

#include <MultiPlasma.H>

Public Member Functions

 MultiPlasma ()
 
 ~MultiPlasma ()
 
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. More...
 
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)
 
void ExplicitDeposition (Fields &fields, const MultiLaser &multi_laser, amrex::Vector< amrex::Geometry > const &gm, int const lev)
 
amrex::Real maxDensity (amrex::Real z)
 Return max density, to compute the adaptive time step. More...
 
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. More...
 
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. More...
 
void DoFieldIonization (const int lev, const amrex::Geometry &geom, const Fields &fields)
 
bool IonizationOn () const
 
bool AnySpeciesNeutralizeBackground () const
 whether any plasma species uses a neutralizing background, e.g. no ion motion More...
 
void TileSort (amrex::Box bx, amrex::Geometry geom)
 sort particles of all containers by tile logically, and store results in m_all_bins More...
 
const amrex::Vector< std::string > & GetNames () const
 
int GetNPlasmas () const
 
void ReorderParticles (const int islice)
 
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. More...
 
void InSituComputeDiags (int step, int islice, int max_step, amrex::Real physical_time, amrex::Real max_time)
 
void InSituWriteToFile (int step, amrex::Real time, const amrex::Geometry &geom, int max_step, amrex::Real max_time)
 

Public Attributes

int m_sort_bin_size {32}
 
amrex::Vector< PlasmaParticleContainerm_all_plasmas
 
int m_nplasmas = 0
 
amrex::Vector< PlasmaBinsm_all_bins
 
amrex::Vector< std::string > m_names {"no_plasma"}
 

Private Attributes

amrex::Real m_adaptive_density = 0.
 

Constructor & Destructor Documentation

◆ MultiPlasma()

MultiPlasma::MultiPlasma ( )

Constructor

◆ ~MultiPlasma()

MultiPlasma::~MultiPlasma ( )
inline

Destructor

Member Function Documentation

◆ AdvanceParticles()

void MultiPlasma::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.

Parameters
[in,out]fieldsthe general field class, modified by this function
[in]multi_laserLasers that affects the plasma during the deposition
[in]gmGeometry of the simulation, to get the cell size etc.
[in]temp_sliceif true, the temporary data (x_temp, ...) will be used
[in]levMR level

◆ AnySpeciesNeutralizeBackground()

bool MultiPlasma::AnySpeciesNeutralizeBackground ( ) const

whether any plasma species uses a neutralizing background, e.g. no ion motion

◆ DepositCurrent()

void MultiPlasma::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 
)

Loop over plasma species and depose their currents into the current 2D slice in fields

Parameters
[in,out]fieldsthe general field class, modified by this function
[in]multi_laserLasers that affects the plasma during the deposition
[in]which_slicedefines if this or the next slice is handled
[in]deposit_jx_jyif true, deposit to jx and jy
[in]deposit_jzif true, deposit to jz
[in]deposit_rhoif true, deposit to rho
[in]deposit_chiif true, deposit chi
[in]deposit_rhomjzif true, deposit rhomjz
[in]gmGeometry of the simulation, to get the cell size etc.
[in]levMR level

◆ DepositNeutralizingBackground()

void MultiPlasma::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.

Parameters
[in,out]fieldsthe general field class, modified by this function
[in]multi_laserthat affects the plasma during the deposition
[in]which_sliceslice in which the densities are deposited
[in]gmGeometry of the simulation, to get the cell size etc.
[in]levMR level

◆ DoFieldIonization()

void MultiPlasma::DoFieldIonization ( const int  lev,
const amrex::Geometry geom,
const Fields fields 
)

Calculates Ionization Probability and makes new Plasma Particles

Parameters
[in]levMR level
[in]geomGeometry of the simulation, to get the cell size
[in]fieldsthe general field class

◆ ExplicitDeposition()

void MultiPlasma::ExplicitDeposition ( Fields fields,
const MultiLaser multi_laser,
amrex::Vector< amrex::Geometry > const &  gm,
int const  lev 
)

Loop over plasma species and depose Sx and Sy into the current 2D slice in fields

Parameters
[in,out]fieldsthe general field class, modified by this function
[in]multi_laserLasers that affects the plasma during the deposition
[in]gmGeometry of the simulation, to get the cell size etc.
[in]levMR level

◆ GetNames()

const amrex::Vector<std::string>& MultiPlasma::GetNames ( ) const
inline

returns a Vector of names of the plasmas

◆ GetNPlasmas()

int MultiPlasma::GetNPlasmas ( ) const
inline

returns number of plasma species

◆ InitData()

void MultiPlasma::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.

Parameters
[in]slice_baslice boxarray, on which plasma particles are defined
[in]slice_dmDistributionMapping of the transverse slice domain
[in]slice_gmslice geometry
[in]gmGeometry of the simulation, to get the cell size

◆ InSituComputeDiags()

void MultiPlasma::InSituComputeDiags ( int  step,
int  islice,
int  max_step,
amrex::Real  physical_time,
amrex::Real  max_time 
)

Compute reduced plasma diagnostics of current slice, store in member variable.

Parameters
[in]steptime step of simulation
[in]islicecurrent slice, on which diags are computed.
[in]max_stepmaximum time step of simulation
[in]physical_timephysical time at the given step
[in]max_timemaximum time of simulation

◆ InSituWriteToFile()

void MultiPlasma::InSituWriteToFile ( int  step,
amrex::Real  time,
const amrex::Geometry geom,
int  max_step,
amrex::Real  max_time 
)

Write reduced beam diagnostics to file

Parameters
[in]steptime step of simulation
[in]timephysical time at the given step
[in]geomSimulation geometry
[in]max_stepmaximum time step of simulation
[in]max_timemaximum time of simulation

◆ IonizationOn()

bool MultiPlasma::IonizationOn ( ) const

◆ maxDensity()

amrex::Real MultiPlasma::maxDensity ( amrex::Real  z)

Return max density, to compute the adaptive time step.

the max is taken across species AND include m_adaptive_density, giving a way to specify a density to the adaptive time step calculator even with no plasma species.

◆ ReorderParticles()

void MultiPlasma::ReorderParticles ( const int  islice)

Reorder particles to speed-up current deposition

Parameters
[in]islicezeta slice index

◆ TagByLevel()

void MultiPlasma::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.

Parameters
[in]current_N_levelnumber of MR levels active on the current slice
[in]geom3DGeometry object for the whole domain
[in]to_previf particles should be tagged to x_prev and y_prev

◆ TileSort()

void MultiPlasma::TileSort ( amrex::Box  bx,
amrex::Geometry  geom 
)

sort particles of all containers by tile logically, and store results in m_all_bins

Parameters
[in]bxtransverse box on which the particles are sorted
[in]geomGeometry object

Member Data Documentation

◆ m_adaptive_density

amrex::Real MultiPlasma::m_adaptive_density = 0.
private

Background (hypothetical) density, used to compute the adaptive time step

◆ m_all_bins

amrex::Vector<PlasmaBins> MultiPlasma::m_all_bins

Logical tile bins for all plasma containers

◆ m_all_plasmas

amrex::Vector<PlasmaParticleContainer> MultiPlasma::m_all_plasmas

contains all plasma containers

◆ m_names

amrex::Vector<std::string> MultiPlasma::m_names {"no_plasma"}

names of all plasma containers

◆ m_nplasmas

int MultiPlasma::m_nplasmas = 0

number of plasma containers

◆ m_sort_bin_size

int MultiPlasma::m_sort_bin_size {32}

Tile size to sort plasma particles


The documentation for this class was generated from the following files: