7#include <boost/math/special_functions/gamma.hpp>
35PlasmaPhase::~PlasmaPhase()
37 if (shared_ptr<Solution> soln =
m_soln.lock()) {
38 soln->removeChangedCallback(
this);
39 soln->kinetics()->removeReactionAddedCallback(
this);
50 throw CanteraError(
"PlasmaPhase::updateElectronEnergyDistribution",
51 "Invalid for discretized electron energy distribution.");
64 throw CanteraError(
"PlasmaPhase::normalizeElectronEnergyDistribution",
65 "The norm is negative. This might be caused by bad "
66 "electron energy distribution");
73 if (
type ==
"discretized" ||
74 type ==
"isotropic") {
77 throw CanteraError(
"PlasmaPhase::setElectronEnergyDistributionType",
78 "Unknown type for electron energy distribution.");
85 double x = m_isotropicShapeFactor;
86 double gamma1 = boost::math::tgamma(3.0 / 2.0 * x);
87 double gamma2 = boost::math::tgamma(5.0 / 2.0 * x);
88 double c1 = x * std::pow(gamma2, 1.5) / std::pow(gamma1, 2.5);
89 double c2 = x * std::pow(gamma2 / gamma1, x);
108 bool updateEnergyDist)
140 throw CanteraError(
"PlasmaPhase::checkElectronEnergyLevels",
141 "Values of electron energy levels need to be positive and "
142 "monotonically increasing.");
151 throw CanteraError(
"PlasmaPhase::checkElectronEnergyDistribution",
152 "Values of electron energy distribution cannot be negative.");
155 warn_user(
"PlasmaPhase::checkElectronEnergyDistribution",
156 "The value of the last element of electron energy distribution exceed 0.01. "
157 "This indicates that the value of electron energy level is not high enough "
158 "to contain the isotropic distribution at mean electron energy of "
170 Eigen::Map<const Eigen::ArrayXd>(dist, length);
194 if (epsilon_m < 0.0) {
195 throw CanteraError(
"PlasmaPhase::updateElectronTemperatureFromEnergyDist",
196 "The electron energy distribution produces negative electron temperature.");
203 m_isotropicShapeFactor = x;
214 eedf[
"energy-levels"] = levels;
216 eedf[
"shape-factor"] = m_isotropicShapeFactor;
221 eedf[
"distribution"] = dist;
224 phaseNode[
"electron-energy-distribution"] = std::move(eedf);
231 if (phaseNode.
hasKey(
"electron-energy-distribution")) {
232 const AnyMap eedf = phaseNode[
"electron-energy-distribution"].as<
AnyMap>();
235 if (eedf.
hasKey(
"shape-factor")) {
236 m_isotropicShapeFactor = eedf[
"shape-factor"].asDouble();
239 "isotropic type requires shape-factor key.");
241 if (eedf.
hasKey(
"energy-levels")) {
243 eedf[
"energy-levels"].asVector<double>().size(),
246 if (eedf.
hasKey(
"mean-electron-energy")) {
247 double energy = eedf.
convert(
"mean-electron-energy",
"eV");
252 "isotropic type requires electron-temperature key.");
255 if (!eedf.
hasKey(
"energy-levels")) {
257 "Cannot find key energy-levels.");
259 if (!eedf.
hasKey(
"distribution")) {
261 "Cannot find key distribution.");
263 if (eedf.
hasKey(
"normalize")) {
267 eedf[
"distribution"].asVector<double>().data(),
268 eedf[
"energy-levels"].asVector<double>().size());
278 if (spec->composition.find(
"E") != spec->composition.end() &&
279 spec->composition.size() == 1 &&
280 spec->composition[
"E"] == 1) {
285 "Cannot add species, {}. "
286 "Only one electron species is allowed.", spec->name);
298 "No electron species found.");
306 if (shared_ptr<Solution> soln =
m_soln.lock()) {
307 soln->registerChangedCallback(
this, [&]() {
319 if (shared_ptr<Solution> soln =
m_soln.lock()) {
320 shared_ptr<Kinetics> kin = soln->kinetics();
326 for (
size_t i = 0; i < kin->nReactions(); i++) {
327 std::shared_ptr<Reaction> R = kin->reaction(i);
328 if (R->rate()->type() !=
"electron-collision-plasma") {
336 kin->registerReactionAddedCallback(
this, [
this, kin]() {
337 size_t i = kin->nReactions() - 1;
338 if (kin->reaction(i)->type() ==
"electron-collision-plasma") {
339 addCollision(kin->reaction(i));
351 collision->registerSetRateCallback(
this, [
this, i, collision]() {
354 std::dynamic_pointer_cast<ElectronCollisionPlasmaRate>(collision->rate());
358 for (
const auto& [
name, _] : collision->reactants) {
368 std::dynamic_pointer_cast<ElectronCollisionPlasmaRate>(collision->rate()));
396 for (
size_t i = 1; i <
m_nPoints - 1; i++) {
402 (h1 * h0) / (h1 + h0);
431 if (interpChanged[i]) {
450 auto cs_array = Eigen::Map<const Eigen::ArrayXd>(
492 if (cached.
state1 != tempNow || cached.
state2 != electronTempNow) {
496 cached.
state2 = electronTempNow;
539 for (
size_t k = 0; k <
m_kk; k++) {
540 ubar[k] =
RT() * (_h[k] - 1.0);
543 ubar[k] =
RTe() * (_h[k] - 1.0);
551 mu[k] += (
RTe() -
RT()) * log(xx);
565 copy(_s.begin(), _s.end(), sr);
567 for (
size_t k = 0; k <
m_kk; k++) {
579 copy(gibbsrt.begin(), gibbsrt.end(), grt);
581 for (
size_t k = 0; k <
m_kk; k++) {
Header for plasma reaction rates parameterized by electron collision cross section and electron energ...
Base class for kinetics managers and also contains the kineticsmgr module documentation (see Kinetics...
Header file for class PlasmaPhase.
Declaration for class Cantera::Species.
A map of string keys to values whose type can vary at runtime.
bool hasKey(const string &key) const
Returns true if the map contains an item named key.
double convert(const string &key, const string &units) const
Convert the item stored by the given key to the units specified in units.
Base class for exceptions thrown by Cantera classes.
const vector< double > & entropy_R_ref() const
Returns a reference to the dimensionless reference state Entropy vector.
double enthalpy_mole() const override
Return the Molar enthalpy. Units: J/kmol.
void getPartialMolarEnthalpies(double *hbar) const override
Returns an array of partial molar enthalpies for the species in the mixture.
void getChemPotentials(double *mu) const override
Get the species chemical potentials. Units: J/kmol.
double pressure() const override
Pressure.
vector< double > m_g0_RT
Temporary storage for dimensionless reference state Gibbs energies.
vector< double > m_h0_RT
Temporary storage for dimensionless reference state enthalpies.
void getGibbs_ref(double *g) const override
Returns the vector of the Gibbs function of the reference state at the current temperature of the sol...
void getStandardChemPotentials(double *mu) const override
Get the array of chemical potentials at unit activity for the species at their standard states at the...
const vector< double > & gibbs_RT_ref() const
Returns a reference to the dimensionless reference state Gibbs free energy vector.
void getStandardVolumes_ref(double *vol) const override
Get the molar volumes of the species reference states at the current T and P_ref of the solution.
virtual void updateThermo() const
Update the species reference state thermodynamic functions.
vector< double > m_s0_R
Temporary storage for dimensionless reference state entropies.
vector< double > m_cp0_R
Temporary storage for dimensionless reference state heat capacities.
bool addSpecies(shared_ptr< Species > spec) override
Add a Species to this Phase.
void getPartialMolarEntropies(double *sbar) const override
Returns an array of partial molar entropies of the species in the solution.
const vector< double > & enthalpy_RT_ref() const
Returns a reference to the dimensionless reference state enthalpy vector.
virtual void update_single(size_t k, double T, double *cp_R, double *h_RT, double *s_R) const
Get reference-state properties for a single species.
ValueCache m_cache
Cached for saved calculations within each ThermoPhase.
size_t m_kk
Number of species in the phase.
double temperature() const
Temperature (K).
virtual double concentration(const size_t k) const
Concentration of species k.
size_t speciesIndex(const string &name) const
Returns the index of a species named 'name' within the Phase object.
double moleFraction(size_t k) const
Return the mole fraction of a single species.
double molecularWeight(size_t k) const
Molecular weight of species k.
string name() const
Return the name of the phase.
void addCollision(std::shared_ptr< Reaction > collision)
Add a collision and record the target species.
void checkElectronEnergyDistribution() const
Check the electron energy distribution.
void setCollisions()
Set collisions.
double meanElectronEnergy() const
Mean electron energy [eV].
double enthalpy_mole() const override
Return the Molar enthalpy. Units: J/kmol.
size_t m_nPoints
Number of points of electron energy levels.
virtual void setSolution(std::weak_ptr< Solution > soln) override
Set the link to the Solution object that owns this ThermoPhase.
void getPartialMolarEnthalpies(double *hbar) const override
Returns an array of partial molar enthalpies for the species in the mixture.
void getChemPotentials(double *mu) const override
Get the species chemical potentials. Units: J/kmol.
void normalizeElectronEnergyDistribution()
Electron energy distribution norm.
void getStandardChemPotentials(double *muStar) const override
Get the array of chemical potentials at unit activity for the species at their standard states at the...
void updateThermo() const override
Update the species reference state thermodynamic functions.
vector< size_t > m_targetSpeciesIndices
The collision-target species indices of m_collisions.
void setElectronTemperature(double Te) override
Set the internally stored electron temperature of the phase (K).
void electronEnergyLevelChanged()
When electron energy level changed, plasma properties such as electron-collision reaction rates need ...
void setElectronEnergyLevels(const double *levels, size_t length, bool updateEnergyDist=true)
Set electron energy levels.
double elasticPowerLoss()
The elastic power loss (J/s/m³)
int m_levelNum
Electron energy level change variable.
bool updateInterpolatedCrossSection(size_t k)
Update interpolated cross section of a collision.
void electronEnergyDistributionChanged()
When electron energy distribution changed, plasma properties such as electron-collision reaction rate...
void getEntropy_R(double *sr) const override
Get the array of nondimensional Entropy functions for the standard state species at the current T and...
size_t nElectronEnergyLevels() const
Number of electron levels.
void getGibbs_ref(double *g) const override
Returns the vector of the Gibbs function of the reference state at the current temperature of the sol...
size_t nCollisions() const
Number of collisions.
Eigen::ArrayXd m_electronEnergyDist
Normalized electron energy distribution vector [-] Length: m_nPoints.
Eigen::ArrayXd m_electronEnergyLevels
electron energy levels [ev]. Length: m_nPoints
void setDiscretizedElectronEnergyDist(const double *levels, const double *distrb, size_t length)
Set discretized electron energy distribution.
void getParameters(AnyMap &phaseNode) const override
Store the parameters of a ThermoPhase object such that an identical one could be reconstructed using ...
string type() const override
String indicating the thermodynamic model implemented.
void checkElectronEnergyLevels() const
Check the electron energy levels.
void initThermo() override
Initialize the ThermoPhase object after all species have been set up.
void updateElasticElectronEnergyLossCoefficients()
Update elastic electron energy loss coefficients.
void updateElectronTemperatureFromEnergyDist()
Update electron temperature (K) From energy distribution.
vector< double > m_interp_cs
Interpolated cross sections.
string m_distributionType
Electron energy distribution type.
void getStandardVolumes_ref(double *vol) const override
Get the molar volumes of the species reference states at the current T and P_ref of the solution.
void updateElectronEnergyDistribution()
Update electron energy distribution.
vector< double > m_elasticElectronEnergyLossCoefficients
Elastic electron energy loss coefficients (eV m3/s)
string m_quadratureMethod
Numerical quadrature method for electron energy distribution.
PlasmaPhase(const string &inputFile="", const string &id="")
Construct and initialize a PlasmaPhase object directly from an input file.
double m_electronTemp
Electron temperature [K].
double RTe() const
Return the Gas Constant multiplied by the current electron temperature.
void getGibbs_RT(double *grt) const override
Get the nondimensional Gibbs functions for the species in their standard states at the current T and ...
bool m_do_normalizeElectronEnergyDist
Flag of normalizing electron energy distribution.
AnyMap m_root
Data for initiate reaction.
void updateElectronEnergyDistDifference()
Update electron energy distribution difference.
void updateElasticElectronEnergyLossCoefficient(size_t i)
Updates the elastic electron energy loss coefficient for collision index i.
void getPartialMolarIntEnergies(double *ubar) const override
Return an array of partial molar internal energies for the species in the mixture.
string electronSpeciesName() const
Electron species name.
void setIsotropicElectronEnergyDistribution()
Set isotropic electron energy distribution.
Eigen::ArrayXd m_electronEnergyDistDiff
Electron energy distribution Difference dF/dε (V^-5/2)
bool addSpecies(shared_ptr< Species > spec) override
Add a Species to this Phase.
vector< bool > m_interp_cs_ready
The list of whether the interpolated cross sections is ready.
vector< shared_ptr< ElectronCollisionPlasmaRate > > m_collisionRates
The list of shared pointers of collision rates.
vector< shared_ptr< Reaction > > m_collisions
The list of shared pointers of plasma collision reactions.
void setParameters(const AnyMap &phaseNode, const AnyMap &rootNode=AnyMap()) override
Set equation of state parameters from an AnyMap phase description.
void setMeanElectronEnergy(double energy)
Set mean electron energy [eV].
size_t m_electronSpeciesIndex
Index of electron species.
void setElectronEnergyDistributionType(const string &type)
Set electron energy distribution type.
void getPartialMolarEntropies(double *sbar) const override
Returns an array of partial molar entropies of the species in the solution.
virtual double electronPressure() const
Electron pressure.
double electronTemperature() const override
Electron Temperature (K)
void setIsotropicShapeFactor(double x)
Set the shape factor of isotropic electron energy distribution.
void enableNormalizeElectronEnergyDist(bool enable)
Set flag of automatically normalize electron energy distribution Flag: m_do_normalizeElectronEnergyDi...
int m_distNum
Electron energy distribution change variable.
virtual void setParameters(const AnyMap &phaseNode, const AnyMap &rootNode=AnyMap())
Set equation of state parameters from an AnyMap phase description.
virtual void getParameters(AnyMap &phaseNode) const
Store the parameters of a ThermoPhase object such that an identical one could be reconstructed using ...
double RT() const
Return the Gas Constant multiplied by the current temperature.
virtual void setSolution(std::weak_ptr< Solution > soln)
Set the link to the Solution object that owns this ThermoPhase.
virtual void initThermo()
Initialize the ThermoPhase object after all species have been set up.
void initThermoFile(const string &inputFile, const string &id)
Initialize a ThermoPhase object using an input file.
std::weak_ptr< Solution > m_soln
reference to Solution
MultiSpeciesThermo m_spthermo
Pointer to the calculation manager for species reference-state thermodynamic properties.
virtual double refPressure() const
Returns the reference pressure in Pa.
CachedScalar getScalar(int id)
Get a reference to a CachedValue object representing a scalar (double) with the given id.
int getId()
Get a unique id for a cached value.
Header for a file containing miscellaneous numerical functions.
This file contains definitions for utility functions and text for modules, inputfiles and logging,...
double numericalQuadrature(const string &method, const Eigen::ArrayXd &f, const Eigen::ArrayXd &x)
Numerical integration of a function.
const double Boltzmann
Boltzmann constant [J/K].
const double Avogadro
Avogadro's Number [number/kmol].
const double GasConstant
Universal Gas Constant [J/kmol/K].
const double ElectronCharge
Elementary charge [C].
const double ElectronMass
Electron Mass [kg].
void warn_user(const string &method, const string &msg, const Args &... args)
Print a user warning raised from method as CanteraWarning.
Namespace for the Cantera kernel.
const size_t npos
index returned by functions to indicate "no position"
const double SmallNumber
smallest number to compare to zero.
A cached property value and the state at which it was evaluated.
double state2
Value of the second state variable for the state at which value was evaluated, for example density or...
bool validate(double state1New)
Check whether the currently cached value is valid based on a single state variable.
double state1
Value of the first state variable for the state at which value was evaluated, for example temperature...