chombo-discharge
|
Class which implements CdrPlasmaPhysics and parses plasma chemistry from a JSON input file. More...
#include <CD_CdrPlasmaJSON.H>
Public Types | |
using | InitialDataFunction = std::function< Real(const RealVect a_position, const Real a_time)> |
Function alias for initial data function. More... | |
using | FunctionEN = std::function< Real(const Real a_E, const Real a_N)> |
Function for encapsulating operations f = f(E,N). field in Townsend units. More... | |
using | FunctionX = std::function< Real(const RealVect a_position)> |
Function for encapsulating a function f = f(x) where x is the physical coordinates. More... | |
using | FunctionEX = std::function< Real(const Real E, const RealVect x)> |
Function for encapsulating a function f = f(E, x) where E is the electric field at physical coordinates x. More... | |
using | FunctionT = std::function< Real(const Real a_T)> |
Function for encapsulating a function f = f(T) where T is the temperature of some species. More... | |
using | FunctionTT = std::function< Real(const Real a_T1, const Real a_T2)> |
Function for encapsulating a function f = f(T1, T2) where T1/T2 are temperatures of two species. More... | |
Public Member Functions | |
CdrPlasmaJSON () | |
Default constructor. Puts object in usable state. | |
CdrPlasmaJSON (const int a_dummy) | |
Dummy constructor. Leaves object in undefined state. More... | |
virtual | ~CdrPlasmaJSON () |
Destructor. | |
virtual void | parseRuntimeOptions () override |
Parse run-time class options. | |
virtual int | getNumberOfPlotVariables () const override |
Get number of plot variables for this physics class. More... | |
virtual Vector< std::string > | getPlotVariableNames () const override |
Get plot variable names. The names and positions between this routine and getPlotVariables must be consistent! | |
virtual Vector< Real > | getPlotVariables (const Vector< Real > a_cdrDensities, const Vector< RealVect > a_cdrGradients, const Vector< Real > a_rteDensities, const RealVect a_E, const RealVect a_position, const Real a_dx, const Real a_dt, const Real a_time, const Real a_kappa) const override |
Provide plot variables. This is used by CdrPlasmaStepper when writing plot files. More... | |
virtual Real | computeAlpha (const Real E, const RealVect a_position) const override |
Compute alpha. Should return Townsend ionization coefficient. More... | |
virtual Real | computeEta (const Real a_E, const RealVect a_position) const override |
Compute eta. Should return Townsend attachment coefficient. More... | |
virtual void | advanceReactionNetwork (Vector< Real > &a_cdrSources, Vector< Real > &a_rteSources, const Vector< Real > a_cdrDensities, const Vector< RealVect > a_cdrGradients, const Vector< Real > a_rteDensities, const RealVect a_E, const RealVect a_pos, const Real a_dx, const Real a_dt, const Real a_time, const Real a_kappa) const override |
Routine intended for advancing a reaction network over a time a_dt. More... | |
virtual Vector< RealVect > | computeCdrDriftVelocities (const Real a_time, const RealVect a_pos, const RealVect a_E, const Vector< Real > a_cdrDensities) const override |
Compute velocities for the CDR equations. More... | |
virtual Vector< Real > | computeCdrDiffusionCoefficients (const Real a_time, const RealVect a_pos, const RealVect a_E, const Vector< Real > a_cdrDensities) const override |
Compute diffusion coefficients for the CDR equations. More... | |
virtual Vector< Real > | computeCdrElectrodeFluxes (const Real a_time, const RealVect a_pos, const RealVect a_normal, const RealVect a_E, const Vector< Real > a_cdrDensities, const Vector< Real > a_cdrVelocities, const Vector< Real > a_cdrGradients, const Vector< Real > a_rteFluxes, const Vector< Real > a_extrapCdrFluxes) const override |
Compute CDR fluxes on electrode-gas interfaces. This is used as a boundary condition in the CDR equations. More... | |
virtual Vector< Real > | computeCdrDielectricFluxes (const Real a_time, const RealVect a_pos, const RealVect a_normal, const RealVect a_E, const Vector< Real > a_cdrDensities, const Vector< Real > a_cdrVelocities, const Vector< Real > a_cdrGradients, const Vector< Real > a_rteFluxes, const Vector< Real > a_extrapCdrFluxes) const override |
Compute CDR fluxes on dielectric-gas interfaces. This is used as a boundary condition in the CDR equations. More... | |
virtual Vector< Real > | computeCdrDomainFluxes (const Real a_time, const RealVect a_pos, const int a_dir, const Side::LoHiSide a_side, const RealVect a_E, const Vector< Real > a_cdrDensities, const Vector< Real > a_cdrVelocities, const Vector< Real > a_cdrGradients, const Vector< Real > a_rteFluxes, const Vector< Real > a_extrapCdrFluxes) const override |
Compute CDR fluxes through domain sides. This is used as a boundary condition in the CDR equations. More... | |
virtual Real | initialSigma (const Real a_time, const RealVect a_pos) const override |
Set the initial surface charge. More... | |
Public Member Functions inherited from Physics::CdrPlasma::CdrPlasmaPhysics | |
CdrPlasmaPhysics () | |
Default constructor. Does nothing. | |
virtual | ~CdrPlasmaPhysics () |
Base destructor. Does nothing. | |
const Vector< RefCountedPtr< CdrSpecies > > & | getCdrSpecies () const |
Get all CDR species. | |
const Vector< RefCountedPtr< RtSpecies > > & | getRtSpecies () const |
Get all RTE species. | |
int | getNumCdrSpecies () const |
Return number of CDR species that we solve for. | |
int | getNumRtSpecies () const |
Return number of RTE equations that we solve for. | |
Protected Types | |
enum class | LookupMethod { Constant , FunctionX , FunctionT , FunctionTT , FunctionEN , FunctionEX , TableEN , TableEnergy , AlphaV , EtaV } |
Enum class for distinguishing types of computation methods when computing transport data stuff. | |
enum class | ReactionIntegrator { None , ExplicitEuler , ExplicitTrapezoidal , ExplicitMidpoint , ExplicitRK4 } |
Enum for distinguishing integration methods. Note that 'None' just fills directly with source terms. | |
enum class | ReactiveEnergyLoss { AddMean , SubtractMean , AddDirect , SubtractDirect , External } |
Enum class for distinguishing how we add/lose energy when running LEA-based models. The user will specify how to do this in the chemistry file, where 'AddMean' means that we add the mean energy to the equation. Likewise, 'SubtractMean' implies removing the mean energy, and 'External' is a user-specified energy to be added/removed. | |
Protected Member Functions | |
virtual void | parseOptions () |
Parse class options. | |
virtual void | parseIntegrator () |
Parse the reactive integrator. | |
virtual void | parseJSON () |
Parse the JSON file. | |
virtual void | initializeSigma () |
Initialize surface charge. More... | |
virtual void | initializeNeutralSpecies () |
Initialize neutral species. | |
virtual void | initializePlasmaSpecies () |
Initialize species. More... | |
virtual void | initializePhotonSpecies () |
Initialize photon species. More... | |
virtual void | parseMobilities () |
Initialize species mobilities. | |
virtual void | parseDiffusion () |
Initialize species diffusion coefficients. | |
virtual void | parseTemperatures () |
Initialize species temperatures. | |
virtual void | parseAlpha () |
Parse the Townsend ionization coefficient. | |
virtual void | parseEta () |
Parse the Townsend attachment coefficient. | |
virtual void | parsePlasmaReactions () |
Parse plasma reactions. | |
virtual InitialDataFunction | parsePlasmaSpeciesInitialData (const json &a_json) const |
Generate an initial data function for a given plasma species. More... | |
virtual List< PointParticle > | parsePlasmaSpeciesInitialParticles (const json &a_json) const |
Generate initial particles for a given plasma species. More... | |
virtual std::list< std::tuple< std::string, std::vector< std::string >, std::vector< std::string > > > | parseReactionWildcards (const std::vector< std::string > &a_reactants, const std::vector< std::string > &a_products, const json &a_reaction) |
Make a reaction set into a superset. This parses wildcards '@' in reaction string. More... | |
virtual void | parseReactionString (std::vector< std::string > &a_reactants, std::vector< std::string > &a_products, const std::string &a_reaction) const |
Parses a reaction string into reactangs and products. More... | |
virtual void | sanctifyPlasmaReaction (const std::vector< std::string > &a_reactants, const std::vector< std::string > &a_products, const std::string a_reaction) const |
Check if a plasma-reaction makes sense in terms of the species that have been defined. More... | |
virtual void | getReactionSpecies (std::list< int > &a_plasmaReactants, std::list< int > &a_neutralReactants, std::list< int > &a_photonReactants, std::list< int > &a_plasmaProducts, std::list< int > &a_neutralProducts, std::list< int > &a_photonProducts, const std::vector< std::string > &a_reactants, const std::vector< std::string > &a_products) const |
Get the int-encoding corresponding to species involved in some reaction. More... | |
virtual void | parsePlasmaReactionRate (const int a_reactionIndex, const json &a_reactionJSON) |
Parse reaction rate for plasma reaction. More... | |
virtual void | parsePlasmaReactionScaling (const int a_reactionIndex, const json &a_reactionJSON) |
Parse scaling factors for reactions. More... | |
virtual void | parsePlasmaReactionPlot (const int a_reactionIndex, const json &a_reactionJSON) |
Parse reaction plotting. More... | |
virtual void | parsePlasmaReactionDescription (const int a_reactionIndex, const json &a_reactionJSON, const std::string a_wildcard) |
Parse plasma reaction descriptions. More... | |
virtual void | parsePlasmaReactionSoloviev (const int a_reactionIndex, const json &a_reactionJSON) |
Parse plasma reaction energy correction. More... | |
virtual void | parsePlasmaReactionEnergyLosses (const int a_reactionIndex, const json &a_reactionJSON) |
Parse plasma reaction energy losses. More... | |
virtual void | parsePhotoReactions () |
Parse photo-reactions. | |
virtual void | parsePhotoReactionScaling (const int a_reactionIndex, const json &a_reactionJSON) |
Parse scaling for photo-reactions. Includes Helmholtz corrections if doing Helmholtz reconstruction of photoionization profiles. More... | |
virtual void | parsePhotoReactionEnergyLosses (const int a_reactionIndex, const json &a_reactionJSON) |
Parse photo-reaction energy losses. More... | |
virtual void | sanityCheckSpecies () const |
Do a species sanity check. More... | |
virtual void | sanctifyPhotoReaction (const std::vector< std::string > &a_reactants, const std::vector< std::string > &a_products, const std::string a_reaction) const |
Check if a photo-reaction makes sense in terms of the species that have been defined. More... | |
virtual void | parseElectrodeReactions () |
Parse secondary emission on electrodes. | |
virtual void | parseElectrodeReactionRate (const int a_reactionIndex, const json &a_reactionJSON) |
Parse reaction rate for electrode surface reactions. More... | |
virtual void | parseElectrodeReactionScaling (const int a_reactionIndex, const json &a_reactionJSON) |
Parse electrode reaction scaling for a specific reaction. More... | |
virtual void | parseElectrodeReactionEnergyLosses (const int a_reactionIndex, const json &a_reactionJSON) |
Parse electrode-reaction energy losses. More... | |
virtual void | parseDielectricReactions () |
Parse secondary emission on dielectrics. | |
virtual void | parseDielectricReactionRate (const int a_reactionIndex, const json &a_reactionJSON) |
Parse reaction rate for dielectric surface reactions. More... | |
virtual void | parseDielectricReactionScaling (const int a_reactionIndex, const json &a_reactionJSON) |
Parse dielectric electrode reaction scaling for a specific reaction. More... | |
virtual void | parseDielectricReactionEnergyLosses (const int a_reactionIndex, const json &a_reactionJSON) |
Parse dielectric-reaction energy losses. More... | |
virtual void | parseDomainReactions () |
Parse secondary emission on domain. | |
virtual void | parseDomainReactionRate (const int a_reactionIndex, const json &a_reactionJSON, const std::vector< std::string > &a_sides) |
Parse reaction rate for domain reactions. More... | |
virtual void | parseDomainReactionScaling (const int a_reactionIndex, const json &a_reactionJSON, const std::vector< std::string > &a_sides) |
Parse domain reaction scaling for a specific reaction. More... | |
virtual void | sanctifySurfaceReaction (const std::vector< std::string > &a_reactants, const std::vector< std::string > &a_products, const std::string a_reaction) const |
Check if a surface-reaction makes sense in terms of the species that have been defined. More... | |
virtual std::vector< Real > | computePlasmaSpeciesTemperatures (const RealVect &a_position, const RealVect &a_E, const std::vector< Real > &a_cdrDensities) const |
Compute the various plasma species temperatures. More... | |
virtual std::vector< Real > | computePlasmaSpeciesEnergies (const RealVect &a_position, const RealVect &a_E, const std::vector< Real > &a_cdrDensities) const |
Compute the various plasma species energies. Returns a list of energies in electron-volts. We assume that temperature-energy relations are e = 3/2 * kB * T. More... | |
virtual std::vector< Real > | computePlasmaSpeciesMobilities (const RealVect &a_position, const RealVect &a_E, const std::vector< Real > &a_cdrDensities) const |
Compute the various plasma species mobilities. More... | |
virtual std::vector< Real > | computePlasmaSpeciesDiffusion (const RealVect a_position, const RealVect a_E, const std::vector< Real > a_cdrDensities) const |
Compute the various plasma species diffusion coefficients. More... | |
virtual Real | computePlasmaReactionRate (const int &a_reactionIndex, const std::vector< Real > &a_cdrDensities, const std::vector< Real > &a_cdrMobilities, const std::vector< Real > &a_cdrDiffusionCoefficients, const std::vector< Real > &a_cdrTemperatures, const std::vector< Real > &a_cdrEnergies, const std::vector< RealVect > &a_cdrGradients, const RealVect &a_pos, const RealVect &a_vectorE, const Real &a_E, const Real &a_Etd, const Real &a_N, const Real &a_alpha, const Real &a_eta, const Real &a_time) const |
Compute the reaction rate for a plasma reaction. More... | |
void | throwParserError (const std::string a_error) const |
Throw a parser error. More... | |
void | throwParserWarning (const std::string a_warning) const |
Throw a parser wearning. More... | |
bool | containsWildcard (const std::string a_str) const |
Protect the @ character in a string. More... | |
bool | containsBracket (const std::string a_str) const |
Protect all kinds of brackets in a string. More... | |
bool | isBracketed (const std::string a_str) const |
Return true if string starts and ends with a paranthesis. More... | |
std::string | trim (const std::string &a_string) const |
Remove whitespace from string. | |
bool | isNeutralSpecies (const std::string &a_name) const |
Return true if species exists in map and false otherwise. More... | |
bool | isPlasmaSpecies (const std::string &a_name) const |
Return true if species exists in map and false otherwise. More... | |
bool | isPhotonSpecies (const std::string &a_name) const |
Return true if species exists in map and false otherwise. More... | |
bool | doesFileExist (const std::string a_filename) const |
Check if file exists. More... | |
virtual void | addPhotoIonization (std::vector< Real > &a_cdrSources, const std::vector< Real > &a_rteDensities, const RealVect a_position, const Real a_E, const Real a_dt, const Real a_dx) const |
Add photoionization products to transport equations source terms. More... | |
virtual void | integrateReactions (std::vector< Real > &a_cdrDensities, std::vector< Real > &a_photonProduction, const std::vector< RealVect > a_cdrGradients, const RealVect a_E, const RealVect a_pos, const Real a_dx, const Real a_dt, const Real a_time, const Real a_kappa) const |
Routine for integrating the reactive-only problem using various algorithms. More... | |
void | fillSourceTerms (std::vector< Real > &a_cdrSources, std::vector< Real > &a_rteSources, const std::vector< Real > a_cdrDensities, const std::vector< RealVect > a_cdrGradients, const RealVect a_E, const RealVect a_pos, const Real a_dx, const Real a_time, const Real a_kappa) const |
Routine for filling the source terms in the reactive problem. More... | |
void | integrateReactionsExplicitEuler (std::vector< Real > &a_cdrDensities, std::vector< Real > &a_photonProduction, const std::vector< RealVect > a_cdrGradients, const RealVect a_E, const RealVect a_pos, const Real a_dx, const Real a_dt, const Real a_time, const Real a_kappa) const |
Routine for integrating the reactive-only problem using the explicit Euler rule. More... | |
void | integrateReactionsImplicitEuler (std::vector< Real > &a_cdrDensities, std::vector< Real > &a_photonProduction, const std::vector< RealVect > a_cdrGradients, const RealVect a_E, const RealVect a_pos, const Real a_dx, const Real a_dt, const Real a_time, const Real a_kappa) const |
Routine for integrating the reactive-only problem using the implicit Euler rule. More... | |
void | integrateReactionsExplicitRK2 (std::vector< Real > &a_cdrDensities, std::vector< Real > &a_photonProduction, const std::vector< RealVect > a_cdrGradients, const RealVect a_E, const RealVect a_pos, const Real a_dx, const Real a_dt, const Real a_time, const Real a_kappa, const Real a_tableuAlpha) const |
Routine for integrating the reactive-only problem using a second order Runge-Kutta method. More... | |
void | integrateReactionsExplicitRK4 (std::vector< Real > &a_cdrDensities, std::vector< Real > &a_photonProduction, const std::vector< RealVect > a_cdrGradients, const RealVect a_E, const RealVect a_pos, const Real a_dx, const Real a_dt, const Real a_time, const Real a_kappa) const |
Routine for integrating the reactive-only problem using the foruth order Runge-Kutta method. More... | |
Protected Attributes | |
bool | m_verbose |
Verbose or not. | |
bool | m_plotGas |
Plot gas pressure, density, and temperature. | |
bool | m_plotAlpha |
Plot Townsend ionization coefficient. | |
bool | m_plotEta |
Plot Townsend attachment coefficient. | |
bool | m_discretePhotons |
Using discrete photons or not. | |
bool | m_skipReactions |
A flag for skipping reactions completely. | |
json | m_json |
JSON definition. This is populated when calling parseJSON. | |
ReactionIntegrator | m_reactionIntegrator |
Reaction integrator. | |
std::string | m_jsonFile |
Input JSON file name. | |
std::vector< json > | m_cdrSpeciesJSON |
JSON entries for species in the defined field 'plasma species'. | |
std::vector< json > | m_rteSpeciesJSON |
JSON entries for species in photon_species. | |
std::function< Real(const RealVect a_position, const Real a_time)> | m_initialSigma |
Initial surface charge. | |
FunctionX | m_gasPressure |
Gas pressure (in Pascal). | |
FunctionX | m_gasTemperature |
Gas temperature (in Kelvin) | |
FunctionX | m_gasDensity |
Gas number density (in m^(-3)) | |
Real | m_chemistryDt |
Chemistry time step. | |
std::vector< FunctionX > | m_neutralSpeciesDensities |
Neutral species densities. | |
std::vector< std::shared_ptr< NeutralSpeciesJSON > > | m_neutralSpecies |
These are the neutral species. | |
std::map< std::string, int > | m_neutralSpeciesMap |
Map for figuring out which where in m_neutralSpecies a neutral species is found. | |
std::map< int, std::string > | m_neutralSpeciesInverseMap |
Inverse of m_neutralSpeciesMap. | |
std::map< std::string, int > | m_cdrSpeciesMap |
string-int encoding of the CDr species. More... | |
std::map< int, std::string > | m_cdrSpeciesInverseMap |
int-string encoding of the CDR species. More... | |
std::map< int, bool > | m_cdrIsEnergySolver |
int-bool encoding for determining if a solver is an energy solver. | |
std::map< int, bool > | m_cdrHasEnergySolver |
int-bool encoding for determining if a CDR solver HAS an associated energy solver. | |
std::map< int, std::tuple< Real, Real, Real > > | m_cdrEnergyComputation |
Parameters for computing the mean energy from energy density and density. More... | |
std::map< int, int > | m_cdrTransportEnergyMap |
int-int encoding for associating a transport solver with an energy solver. More... | |
std::map< int, Real > | m_cdrMasses |
Map of the species masses. This is needed for imposing BCs on the energy equations. More... | |
std::map< std::string, int > | m_rteSpeciesMap |
string-int encoding of the RTE species. More... | |
std::map< int, std::string > | m_rteSpeciesInverseMap |
int-string encoding of the RTE species. More... | |
LookupMethod | m_alphaLookup |
Lookup method for Townsend ionization coefficient. | |
LookupMethod | m_etaLookup |
Lookup method for Townsend attachment coefficient. | |
Real | m_alphaConstant |
For when we can use alpha = constant. | |
Real | m_etaConstant |
For when we can use eta = constant. | |
FunctionEN | m_alphaFunctionEN |
For when we can put alpha = alpha(E,N) as an analytic function. | |
FunctionEN | m_etaFunctionEN |
For when we can put eta = eta(E,N) as an analytic function. | |
LookupTable1D< Real, 1 > | m_alphaTableEN |
For when we can put alpha = table(E,N) | |
LookupTable1D< Real, 1 > | m_etaTableEN |
For when we can put eta = table(E,N) | |
std::map< int, LookupMethod > | m_mobilityLookup |
Mobility lookup method for each species. | |
std::map< int, Real > | m_mobilityConstants |
Map for constant mobilities. | |
std::map< int, FunctionEN > | m_mobilityFunctionsEN |
Map for function-based mobilities mu = mu(E,N) | |
std::map< int, FunctionEX > | m_mobilityFunctionsEX |
Map for function-based mobilities mu = mu(E,x) | |
std::map< int, LookupTable1D< Real, 1 > > | m_mobilityTablesEN |
Map for table-based mobilities. Stored as tables (E/N, mu*N) | |
std::map< int, LookupTable1D< Real, 1 > > | m_mobilityTablesEnergy |
Map for table-based mobilities as function of energy. More... | |
std::map< int, LookupMethod > | m_diffusionLookup |
Diffusion lookup method. | |
std::map< int, Real > | m_diffusionConstants |
Map for constant diffusion coefficients. | |
std::map< int, FunctionEN > | m_diffusionFunctionsEN |
Map for function-based diffusion coefficients. . | |
std::map< int, LookupTable1D< Real, 1 > > | m_diffusionTablesEN |
Map for table-based diffusion coefficients D = D(E,N). More... | |
std::map< int, LookupTable1D< Real, 1 > > | m_diffusionTablesEnergy |
Map for table-based diffusion coefficients as function of energy. More... | |
std::map< int, LookupMethod > | m_temperatureLookup |
Temperature lookup method. | |
std::map< int, FunctionX > | m_temperatureConstants |
Constant temperatures. | |
std::map< int, LookupTable1D< Real, 1 > > | m_temperatureTablesEN |
Temperatures as functions of E/N. | |
std::map< int, std::string > | m_plasmaReactionDescriptions |
Description of plasma reactions. Only used for I/O. | |
std::map< int, LookupMethod > | m_plasmaReactionLookup |
Map for figuring out how to look up the rate for a certain plasma reaction. | |
std::map< int, Real > | m_plasmaReactionConstants |
Constant plasma reaction rates. | |
std::map< int, int > | m_plasmaReactionAlphaV |
Plasma reaction rates that are alpha*|v|. | |
std::map< int, int > | m_plasmaReactionEtaV |
Plasma reaction rates that are eta*|v|. | |
std::map< int, std::pair< int, FunctionT > > | m_plasmaReactionFunctionsT |
Maps for functions of the type k = f(T) where T is the temperature of some species. More... | |
std::map< int, std::tuple< int, int, FunctionTT > > | m_plasmaReactionFunctionsTT |
Maps for functions of the type k = f(T1,T2) where T1 and T2 are the temperatures of some species. More... | |
std::map< int, FunctionEN > | m_plasmaReactionFunctionsEN |
Function-based plasma reaction rates. | |
std::map< int, LookupTable1D< Real, 1 > > | m_plasmaReactionTablesEN |
Map for table-based reaction coefficients, where k = k(E,N). | |
std::map< int, std::pair< int, LookupTable1D< Real, 1 > > > | m_plasmaReactionTablesEnergy |
Map for table-based reaction coefficients where k = k(energy). More... | |
std::map< int, FunctionEX > | m_plasmaReactionEfficiencies |
Scaled plasma reactions. These account for e.g. reaction efficiencies, collisional quenching, etc. | |
std::vector< CdrPlasmaReactionJSON > | m_plasmaReactions |
Plasma reactions. | |
std::map< int, bool > | m_plasmaReactionPlot |
Plot plasma reaction or not. | |
std::map< int, std::pair< bool, int > > | m_plasmaReactionSolovievCorrection |
Flag for whether or not reaction includes Soloviev energy correction. More... | |
std::map< int, std::map< int, std::pair< ReactiveEnergyLoss, Real > > > | m_plasmaReactionEnergyLosses |
For mapping reactive energy losses for all reactions. More... | |
std::map< int, bool > | m_plasmaReactionHasEnergyLoss |
Associative container for determining if a reaction is associated with an energy loss/gain. | |
std::map< int, FunctionEX > | m_photoReactionEfficiencies |
Flag for photo-reaction efficiencies. Includes Helmholtz corrections, if present. | |
std::map< int, bool > | m_photoReactionUseHelmholtz |
Map over the Helmholtz reconstructions. | |
std::vector< CdrPlasmaPhotoReactionJSON > | m_photoReactions |
Photo-reactions. | |
std::map< int, std::list< std::pair< int, Real > > > | m_photoReactionEnergyLosses |
Associated energy losses for a photo-reaction. More... | |
std::map< int, bool > | m_photoReactionHasEnergyLoss |
Associative container for determining if a reaction is associated with an energy loss/gain. | |
std::map< int, LookupMethod > | m_electrodeReactionLookup |
Lookup method for the electrode surface reaction rates. | |
std::map< int, Real > | m_electrodeReactionConstants |
Constant electrode reaction rate. | |
std::map< int, FunctionEX > | m_electrodeReactionEfficiencies |
Electrode reaction effiencies. Used for scaling reactions on electrodes in a "generic" way. | |
std::vector< CdrPlasmaSurfaceReactionJSON > | m_electrodeReactions |
List of electrode reactions. | |
std::map< int, std::list< std::pair< int, Real > > > | m_electrodeReactionEnergyLosses |
Associated energy losses for a surface reaction on electrodes. More... | |
std::map< int, bool > | m_electrodeReactionHasEnergyLoss |
Associative container for determining if a reaction is associated with an energy loss/gain. | |
std::map< int, bool > | m_electrodeExtrapBC |
A container which determines if we should add the extrapolated flux as an inflow condition. | |
std::map< int, LookupMethod > | m_dielectricReactionLookup |
Lookup method for the dielectric surface reaction rates. | |
std::map< int, Real > | m_dielectricReactionConstants |
Constant dielectric reaction rate. | |
std::map< int, FunctionEX > | m_dielectricReactionEfficiencies |
Dielectric reaction effiencies. Used for scaling reactions on dielectrics in a "generic" way. | |
std::vector< CdrPlasmaSurfaceReactionJSON > | m_dielectricReactions |
List of dielectric reactions. | |
std::map< int, std::list< std::pair< int, Real > > > | m_dielectricReactionEnergyLosses |
Associated energy losses for a surface reaction on dielectrics. More... | |
std::map< int, bool > | m_dielectricReactionHasEnergyLoss |
Associative container for determining if a reaction is associated with an energy loss/gain. | |
std::map< int, bool > | m_dielectricExtrapBC |
A container which determines if we should add the extrapolated flux as an inflow condition. | |
std::map< std::pair< int, Side::LoHiSide >, std::map< int, LookupMethod > > | m_domainReactionLookup |
Lookup method for the domain reaction rates. The pair is made up of an int representing direction (0=x, 1=y, 2=z) and a Side::LoHiSide representing side (Side::Lo, Side::Hi) | |
std::map< std::pair< int, Side::LoHiSide >, std::map< int, Real > > | m_domainReactionConstants |
Constant domain reaction rate. The pair is made up of an int representing direction (0=x, 1=y, 2=z) and a Side::LoHiSide representing side (Side::Lo, Side::Hi) | |
std::map< std::pair< int, Side::LoHiSide >, std::map< int, FunctionEX > > | m_domainReactionEfficiencies |
Domain reaction effiencies. Used for scaling reactions on domains in a "generic" way. The pair is made up of an int representing direction (0=x, 1=y, 2=z) and a Side::LoHiSide representing side (Side::Lo, Side::Hi) | |
std::map< std::pair< int, Side::LoHiSide >, std::vector< CdrPlasmaSurfaceReactionJSON > > | m_domainReactions |
List of domain reactions. The pair is made up of an int representing direction (0=x, 1=y, 2=z) and a Side::LoHiSide representing side (Side::Lo, Side::Hi) | |
const std::map< char, int > | m_dirCharToInt {{'x', 0}, {'y', 1}, {'z', 2}} |
map to translate dir from char to int | |
const std::map< std::string, Side::LoHiSide > | m_sideStringToSide {{"lo", Side::Lo}, {"hi", Side::Hi}} |
map to translate side from std::string to Side::LoHiSide | |
std::map< std::tuple< int, Side::LoHiSide, int >, bool > | m_domainExtrapBC |
A container which determines if we should add the extrapolated flux as an inflow condition. More... | |
Protected Attributes inherited from Physics::CdrPlasma::CdrPlasmaPhysics | |
Vector< RefCountedPtr< CdrSpecies > > | m_cdrSpecies |
List of species. | |
Vector< RefCountedPtr< RtSpecies > > | m_rtSpecies |
List of optical transitions between species. | |
int | m_numCdrSpecies |
Number of species. | |
int | m_numRtSpecies |
Number of RTE species. | |
Class which implements CdrPlasmaPhysics and parses plasma chemistry from a JSON input file.
using Physics::CdrPlasma::CdrPlasmaJSON::FunctionEN = std::function<Real(const Real a_E, const Real a_N)> |
Function for encapsulating operations f = f(E,N). field in Townsend units.
[in] | a_E | Electric field in SI units. |
[in] | a_N | Neutral density. |
using Physics::CdrPlasma::CdrPlasmaJSON::FunctionEX = std::function<Real(const Real E, const RealVect x)> |
Function for encapsulating a function f = f(E, x) where E is the electric field at physical coordinates x.
[in] | E | Electric field magnitude in SI units. |
[in] | x | Physical coordinates |
using Physics::CdrPlasma::CdrPlasmaJSON::FunctionT = std::function<Real(const Real a_T)> |
Function for encapsulating a function f = f(T) where T is the temperature of some species.
[in] | a_T | Some temperature. |
using Physics::CdrPlasma::CdrPlasmaJSON::FunctionTT = std::function<Real(const Real a_T1, const Real a_T2)> |
Function for encapsulating a function f = f(T1, T2) where T1/T2 are temperatures of two species.
[in] | a_T1 | Some species temperature. |
[in] | a_T2 | Some other species temperature. |
using Physics::CdrPlasma::CdrPlasmaJSON::FunctionX = std::function<Real(const RealVect a_position)> |
Function for encapsulating a function f = f(x) where x is the physical coordinates.
[in] | a_position | Physical coordinates |
using Physics::CdrPlasma::CdrPlasmaJSON::InitialDataFunction = std::function<Real(const RealVect a_position, const Real a_time)> |
Function alias for initial data function.
[in] | a_position | Physical coordinates |
[in] | a_time | Time |
CdrPlasmaJSON::CdrPlasmaJSON | ( | const int | a_dummy | ) |
Dummy constructor. Leaves object in undefined state.
This constructor is present in order to allow new derived class of CdrPlasmaJSON that instantiate without calling the all-defining base constructor.
[in] | a_dummy | Dummy argument. |
|
protectedvirtual |
Add photoionization products to transport equations source terms.
[in,out] | a_cdrSources | Source terms for CDR densities. |
[in] | a_rteDensities | RTE mesh densities. |
[in] | a_position | Physical coordinates |
[in] | a_E | Electric field (SI units) |
[in] | a_dt | Time step |
[in] | a_dx | Grid resolution. |
|
overridevirtual |
Routine intended for advancing a reaction network over a time a_dt.
This routine assumes that the subsequent advance is in the form phi^(k+1) = phi^k + S*a_dt. Thus, this routine exists such that users can EITHER fill a_cdrSources and a_rteSources directly with an explicit rule, OR they can perform a fully implicit advance within this routine and set S from that.
[out] | a_cdrSources | Source terms for CDR equations. |
[out] | a_rteSources | Source terms for RTE equations. |
[in] | a_cdrDensities | Grid-based density for particle species. |
[in] | a_cdrGradients | Grid-based gradients for particle species. |
[in] | a_rteDensities | Grid-based densities for photons. |
[in] | a_E | Electric field. |
[in] | a_pos | Position in space. |
[in] | a_dx | Grid resolution. |
[in] | a_dt | Advanced time. |
[in] | a_time | Current time. |
[in] | a_kappa | Grid cell unit volume. |
Implements Physics::CdrPlasma::CdrPlasmaPhysics.
|
overridevirtual |
Compute alpha. Should return Townsend ionization coefficient.
This function is mostly used for the cell tagging classes, but can be used for reactions as well.
[in] | a_E | Electric field. |
Implements Physics::CdrPlasma::CdrPlasmaPhysics.
|
overridevirtual |
Compute CDR fluxes on dielectric-gas interfaces. This is used as a boundary condition in the CDR equations.
[in] | a_time | Time |
[in] | a_pos | Position |
[in] | a_normal | Normal vector. This points into the gas phase. |
[in] | a_E | Electric field |
[in] | a_cdrDensities | CDR densities (on the EB) |
[in] | a_cdrVelocities | Normal component of CDR velocities (on the EB). |
[in] | a_cdrGradients | Normal gradients of cdr densities |
[in] | a_rteFluxes | RTE fluxes (normal component only) |
[in] | a_extrapCdrFluxes | Extrapolated fluxes from the gas side. |
Implements Physics::CdrPlasma::CdrPlasmaPhysics.
|
overridevirtual |
Compute diffusion coefficients for the CDR equations.
[in] | a_time | Time |
[in] | a_pos | Position |
[in] | a_E | Electric field |
[in] | a_cdrDensities | CDR densities |
Implements Physics::CdrPlasma::CdrPlasmaPhysics.
|
overridevirtual |
Compute CDR fluxes through domain sides. This is used as a boundary condition in the CDR equations.
[in] | a_time | Time |
[in] | a_pos | Position |
[in] | a_dir | Direction (0 = x, 1=y etc) |
[in] | a_side | Side (low or high side) |
[in] | a_E | Electric field |
[in] | a_cdrDensities | CDR densities. |
[in] | a_cdrVelocities | CDR velocities (normal component only). |
[in] | a_cdrGradients | CDR gradients (normal component only) |
[in] | a_rteFluxes | RTE fluxes (normal component only) |
[in] | a_extrapCdrFluxes | Extrapolated fluxes from the gas side. |
Implements Physics::CdrPlasma::CdrPlasmaPhysics.
|
overridevirtual |
Compute velocities for the CDR equations.
[in] | a_time | Time |
[in] | a_pos | Position |
[in] | a_E | Electric field |
[in] | a_cdrDensities | CDR densities |
Implements Physics::CdrPlasma::CdrPlasmaPhysics.
|
overridevirtual |
Compute CDR fluxes on electrode-gas interfaces. This is used as a boundary condition in the CDR equations.
[in] | a_time | Time |
[in] | a_pos | Position |
[in] | a_normal | Boundary normal vector. This points into the gas phase. |
[in] | a_E | Electric field |
[in] | a_cdrVelocities | CDR velocities. Normal component only. |
[in] | a_cdrDensities | CDR densities. |
[in] | a_cdrGradients | Normal gradients of cdr densities |
[in] | a_rteFluxes | RTE fluxes (normal component only) |
[in] | a_extrapCdrFluxes | Extrapolated fluxes from the gas side. |
Implements Physics::CdrPlasma::CdrPlasmaPhysics.
|
overridevirtual |
Compute eta. Should return Townsend attachment coefficient.
This function is mostly used for the cell tagging classes, but can be used for reactions as well.
[in] | a_E | Electric field. |
[in] | a_position | Position |
Implements Physics::CdrPlasma::CdrPlasmaPhysics.
|
protectedvirtual |
Compute the reaction rate for a plasma reaction.
This routine exists because we need to compute the rates both in advanceReactionNetwork and getPlotVariables. This function reduces code duplication.
[in] | a_reactionIndex | Reaction index |
[in] | a_cdrDensities | Plasma species densities. |
[in] | a_cdrMobilities | Plasma species mobilities. |
[in] | a_cdrDiffusionCoefficients | Plasma species diffusion coefficients. |
[in] | a_cdrTemperatures | Plasma species temperatures. |
[in] | a_cdrEnergies | Plasma species energies. |
[in] | a_cdrGradients | Plasma species gradients. |
[in] | a_pos | Position (physical coordinates) |
[in] | a_vectorE | Electric field (vector) |
[in] | a_E | Electric field magnitude (SI units) |
[in] | a_Etd | Electric field magnitude (Townsend units) |
[in] | a_N | Neutral density |
[in] | a_alpha | Townsend ionization coefficient |
[in] | a_eta | Townsend attachment coefficient |
[in] | a_time | Time |
|
protectedvirtual |
Compute the various plasma species diffusion coefficients.
[in] | a_position | Physical coordinates |
[in] | a_E | Electric field (SI units) |
[in] | a_cdrDensities | List of plasma species densities |
|
protectedvirtual |
Compute the various plasma species energies. Returns a list of energies in electron-volts. We assume that temperature-energy relations are e = 3/2 * kB * T.
If a transported species is associated with an energy solver the energy is computed from e = (n_eps)/ne. Otherwise, we compute the energy and returns 3/2 * kB * T.
[in] | a_position | Physical coordinates |
[in] | a_E | Electric field (SI units) |
[in] | a_cdrDensities | List of plasma species densities. |
|
protectedvirtual |
Compute the various plasma species mobilities.
[in] | a_position | Physical coordinates |
[in] | a_E | Electric field (SI units) |
[in] | a_cdrDensities | List of plasma species densities. |
|
protectedvirtual |
Compute the various plasma species temperatures.
[in] | a_position | Physical coordinates |
[in] | a_E | Electric field (SI units) |
[in] | a_cdrDensities | List of plasma species densities. |
|
protected |
Protect all kinds of brackets in a string.
[in] | a_str | Input string. If it contains any bracket we throw an error. |
|
protected |
Protect the @ character in a string.
[in] | a_str | Input string. If it contains the at character we throw an error. |
|
protected |
Check if file exists.
[in] | a_filename | File name |
|
protected |
Routine for filling the source terms in the reactive problem.
[out] | a_cdrSources | Contains source term for CDR equations. |
[out] | a_rteSources | Contains source terms for RTE equations. |
[in] | a_cdrDensities | CDR densities |
[in] | a_cdrGradients | CDR gradients at time a_time |
[in] | a_E | Electric field |
[in] | a_pos | Physical coordinates |
[in] | a_dx | Grid resolution |
[in] | a_time | Time |
[in] | a_kappa | Volume fraction |
|
overridevirtual |
Get number of plot variables for this physics class.
This is used by CdrPlasmaStepper for pre-allocating data that will be put in a plot file. The current implementation plots the gas pressure, temperature, and density
Reimplemented from Physics::CdrPlasma::CdrPlasmaPhysics.
|
overridevirtual |
Provide plot variables. This is used by CdrPlasmaStepper when writing plot files.
The plot variables should have length this->getNumberOfPlotVariables (as should getPlotVariableNames)
[in] | a_cdrDensities | Grid-based density for particle species. |
[in] | a_cdrGradients | Grid-based gradients for particle species. |
[in] | a_rteDensities | Grid-based densities for photons. |
[in] | a_E | Electric field. |
[in] | a_pos | Position in space. |
[in] | a_dx | Grid resolution. |
[in] | a_dt | Advanced time. |
[in] | a_time | Current time. |
[in] | a_kappa | Grid cell unit volume. |
Reimplemented from Physics::CdrPlasma::CdrPlasmaPhysics.
|
protectedvirtual |
Get the int-encoding corresponding to species involved in some reaction.
[out] | a_plasmaReactants | Plasma reactants |
[out] | a_neutralReactants | Neutral reactants |
[out] | a_photonReactants | Photon reactants |
[out] | a_plasmaProducts | Plasma products |
[out] | a_neutralProducts | Neutral products |
[out] | a_photonProducts | Photon products |
[in] | a_reactants | Reactant names |
[in] | a_products | Reaction names |
|
protectedvirtual |
Initialize photon species.
This will initialize the radiative transfer species absed on the "photon_species" field in the JSON file.
|
protectedvirtual |
Initialize species.
This will initialize the species based on the "plasma_species" field in the JSON file.
|
protectedvirtual |
Initialize surface charge.
This will initialize the surface charge based on the "sigma" field in the JSON file. If you want more complex initial surface charges, you will have to either extend this routine or overwrite it.
|
overridevirtual |
Set the initial surface charge.
[in] | a_time | Time |
[in] | a_pos | Position |
Implements Physics::CdrPlasma::CdrPlasmaPhysics.
|
protectedvirtual |
Routine for integrating the reactive-only problem using various algorithms.
[in,out] | a_cdrDensities | On input, contains n(t). On output it contains n(t+dt). |
[out] | a_photonProduction | On input, should be equal to zero. On output it will contain the number of photons produced during the time step. |
[in] | a_cdrGradients | CDR gradients at time a_time |
[in] | a_E | Electric field |
[in] | a_pos | Physical coordinates |
[in] | a_dx | Grid resolution |
[in] | a_dt | Time step |
[in] | a_kappa | Volume fraction |
|
protected |
Routine for integrating the reactive-only problem using the explicit Euler rule.
[in,out] | a_cdrDensities | On input, contains n(t). On output it contains n(t+dt). |
[out] | a_photonProduction | On input, should be equal to zero. On output it will contain the number of photons produced during the time step. |
[in] | a_cdrGradients | CDR gradients at time a_time |
[in] | a_E | Electric field |
[in] | a_pos | Physical coordinates |
[in] | a_dx | Grid resolution |
[in] | a_dt | Time step |
[in] | a_time | Time |
[in] | a_kappa | Volume fraction |
|
protected |
Routine for integrating the reactive-only problem using a second order Runge-Kutta method.
This includes the tableu through the a_tableuAlpha parameter.
[in,out] | a_cdrDensities | On input, contains n(t). On output it contains n(t+dt). |
[out] | a_photonProduction | On input, should be equal to zero. On output it will contain the number of photons produced during the time step. |
[in] | a_cdrGradients | CDR gradients at time a_time |
[in] | a_E | Electric field |
[in] | a_pos | Physical coordinates |
[in] | a_dx | Grid resolution |
[in] | a_dt | Time step |
[in] | a_time | Time |
[in] | a_kappa | Volume fraction |
[in] | a_tableuAlpha | RK2 tableu alpha. Use 0.5 for midpoint and 1.0 for trapezoidal (Heun's method) |
|
protected |
Routine for integrating the reactive-only problem using the foruth order Runge-Kutta method.
[in,out] | a_cdrDensities | On input, contains n(t). On output it contains n(t+dt). |
[out] | a_photonProduction | On input, should be equal to zero. On output it will contain the number of photons produced during the time step. |
[in] | a_cdrGradients | CDR gradients at time a_time |
[in] | a_E | Electric field |
[in] | a_pos | Physical coordinates |
[in] | a_dx | Grid resolution |
[in] | a_dt | Time step |
[in] | a_time | Time |
[in] | a_kappa | Volume fraction |
|
protected |
Routine for integrating the reactive-only problem using the implicit Euler rule.
[in,out] | a_cdrDensities | On input, contains n(t). On output it contains n(t+dt). |
[out] | a_photonProduction | On input, should be equal to zero. On output it will contain the number of photons produced during the time step. |
[in] | a_cdrGradients | CDR gradients at time a_time |
[in] | a_E | Electric field |
[in] | a_pos | Physical coordinates |
[in] | a_dx | Grid resolution |
[in] | a_dt | Time step |
[in] | a_time | Time |
[in] | a_kappa | Volume fraction |
|
protected |
Return true if string starts and ends with a paranthesis.
[in] | a_str | Input string. |
|
protected |
Return true if species exists in map and false otherwise.
[in] | a_name | Neutral species name |
|
protected |
Return true if species exists in map and false otherwise.
[in] | a_name | Rte species name |
|
protected |
Return true if species exists in map and false otherwise.
[in] | a_name | Cdr species name |
|
protectedvirtual |
Parse dielectric-reaction energy losses.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'dielectric reactions' array. |
|
protectedvirtual |
Parse reaction rate for dielectric surface reactions.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'electrode reactions' array. |
|
protectedvirtual |
Parse dielectric electrode reaction scaling for a specific reaction.
This can be used to adjust for e.g. positional or field dependence.
[in] | a_reactionIndex | Reaction index (this is the array index in the dielectric reaction array supplied in the JSON file) |
[in] | a_reactionJSON | JSON entry for the reaction in the JSON input file. |
|
protectedvirtual |
Parse reaction rate for domain reactions.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'domain reactions' array. |
[in] | a_sides | List of the domain sides that this occurs at |
|
protectedvirtual |
Parse domain reaction scaling for a specific reaction.
This can be used to adjust for e.g. positional or field dependence.
[in] | a_reactionIndex | Reaction index (this is the array index in the domain reaction array supplied in the JSON file) |
[in] | a_reactionJSON | JSON entry for the reaction in the JSON input file. |
[in] | a_sides | List of the domain sides that this occurs at |
|
protectedvirtual |
Parse electrode-reaction energy losses.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'electrode reactions' array. |
|
protectedvirtual |
Parse reaction rate for electrode surface reactions.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'electrode reactions' array. |
|
protectedvirtual |
Parse electrode reaction scaling for a specific reaction.
This can be used to adjust for e.g. positional or field dependence.
[in] | a_reactionIndex | Reaction index (this is the array index in the reaction array supplied in the JSON file) |
[in] | a_reactionJSON | JSON entry for the reaction in the JSON input file. |
|
protectedvirtual |
Parse photo-reaction energy losses.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'photo reactions' array. |
|
protectedvirtual |
Parse scaling for photo-reactions. Includes Helmholtz corrections if doing Helmholtz reconstruction of photoionization profiles.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'photo reactions' array. |
|
protectedvirtual |
Parse plasma reaction descriptions.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'plasma reactions' array. |
|
protectedvirtual |
Parse plasma reaction energy losses.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'plasma reactions' array. |
|
protectedvirtual |
Parse reaction plotting.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'plasma reactions' array. |
|
protectedvirtual |
Parse reaction rate for plasma reaction.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'plasma reactions' array. |
Add the rate and lookup method.
|
protectedvirtual |
Parse scaling factors for reactions.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'plasma reactions' array. |
|
protectedvirtual |
Parse plasma reaction energy correction.
[in] | a_reactionIndex | Reaction index. |
[in] | a_reactionJSON | Input reaction. Must be one of the entries in the 'plasma reactions' array. |
|
protectedvirtual |
Generate an initial data function for a given plasma species.
[in] | a_json | JSON field, usually (always?) describing one of the objects in the 'plasma species' field. |
|
protectedvirtual |
Generate initial particles for a given plasma species.
[in] | a_json | JSON field, usually (always?) describing one of the objects in the 'plasma species' field. |
|
protectedvirtual |
Parses a reaction string into reactangs and products.
[out] | a_reactants | Left-hand side of reaction |
[out] | a_products | Right-hand side of reaction |
[in] | a_reaction | Reaction string. Must be in format "a + b + c -> e + f + g". |
|
protectedvirtual |
Make a reaction set into a superset. This parses wildcards '@' in reaction string.
[in] | a_reactants | List of reactants. Can contain wildcard. |
[in] | a_products | List of products. Can contain wildcard. |
[in] | a_reaction | JSON reaction entry. |
|
protectedvirtual |
Check if a photo-reaction makes sense in terms of the species that have been defined.
This will throw errors if the species do not exist or charge is not conserved.
[in] | a_reactants | Reactants |
[in] | a_products | Reaction products |
[in] | a_reaction | Reaction string |
|
protectedvirtual |
Check if a plasma-reaction makes sense in terms of the species that have been defined.
This will throw errors if the species do not exist or charge is not conserved.
[in] | a_reactants | Reactants |
[in] | a_products | Reaction products |
[in] | a_reaction | Reaction string |
|
protectedvirtual |
Check if a surface-reaction makes sense in terms of the species that have been defined.
[in] | a_reactants | Reactants |
[in] | a_products | Reaction products |
[in] | a_reaction | Reaction string |
|
protectedvirtual |
Do a species sanity check.
This will make sure that species names are not duplicates.
|
protected |
Throw a parser error.
[in] | a_error | Error code. |
|
protected |
Throw a parser wearning.
[in] | a_error | Warning |
|
protected |
Parameters for computing the mean energy from energy density and density.
The energy is computed as e = max(Emin, min(Emax, E)) where E = n_energy/(min(n_density, safety));
|
protected |
Map of the species masses. This is needed for imposing BCs on the energy equations.
This is populated in initializePlasmaSpecies.
|
protected |
int-string encoding of the CDR species.
This is needed because we sometimes need to use the species name for indexing in the vector.
|
protected |
string-int encoding of the CDr species.
This is needed because we sometimes need to use the species name for indexing in the vector.
|
protected |
int-int encoding for associating a transport solver with an energy solver.
This is only defined for the species that has/is an energy solver. The first index is the transport solver index and the second index is the energy solver index.
|
protected |
Associated energy losses for a surface reaction on dielectrics.
The list is a list of species and the corresponding energy losses for a reaction. Note that although we write 'loss', this can also be an energy 'gain'.
|
protected |
Map for table-based diffusion coefficients D = D(E,N).
Tables stored as (E/N, D*N)
|
protected |
Map for table-based diffusion coefficients as function of energy.
Stored as tables (eV, D*N)
|
protected |
A container which determines if we should add the extrapolated flux as an inflow condition.
The first entry in the map is the species integer index. The tuple indicates (coordinate_direction, side, doExtrap).
|
protected |
Associated energy losses for a surface reaction on electrodes.
The list is a list of species and the corresponding energy losses for a reaction. Note that although we write 'loss', this can also be an energy 'gain'.
|
protected |
Map for table-based mobilities as function of energy.
Stored as tables (eV, mu*N)
|
protected |
Associated energy losses for a photo-reaction.
The list is a list of species and the corresponding energy losses for a reaction. Note that although we write 'loss', this can also be an energy 'gain'.
|
protected |
For mapping reactive energy losses for all reactions.
The index in the first map is the reaction (i.e., index in m_plasmaReactions). The second map determines how reactive energy losses/gains are computed for each species. The first index in the second map indicates the species which will add/lose energy, and the std::pair indicates how this loss is computed.
|
protected |
Maps for functions of the type k = f(T) where T is the temperature of some species.
In the tuple the first index in the pair indicates the species we are talking about here. A special rule is enforced when this index is < 0 in which case the temperature is replaced by the background gas temperature.
|
protected |
Maps for functions of the type k = f(T1,T2) where T1 and T2 are the temperatures of some species.
This signature is absolute horrific – what it means is that we have a reaction which should be evaluated as f(T1, T2), but we need to know which species are involved. By design, this should be general so that T1 and T2 can be the temperatures of any species, including neutral species. So, we make a tuple for indicating which species we are talking about. The first index in the tuple is the first species, the second is the second species and the third entry in the tuple is the actual function. A special rule occurs if one of the first two indices is < 0 in which case the temperature is replaced by the background gas temperature.
|
protected |
Flag for whether or not reaction includes Soloviev energy correction.
If this is true, the rate for a reaction (in the local field approximation) will be modified as k * (1 + E.(D * grad(phi))/n *
|
protected |
Map for table-based reaction coefficients where k = k(energy).
The first index is the reaction index, while the pair describes which species energy and the tabulated data.
|
protected |
int-string encoding of the RTE species.
This is needed because we sometimes need to use the species name for indexing in the vector.
|
protected |
string-int encoding of the RTE species.
This is needed because we sometimes need to use the species name for indexing in the vector.