19 if (thermo.
type() !=
"IdealGas") {
21 "Incompatible phase type provided");
23 Reactor::setThermoMgr(thermo);
26 void IdealGasReactor::getState(
double* y)
30 "Error: reactor is empty.");
32 m_thermo->restoreState(m_state);
35 m_mass = m_thermo->density() * m_vol;
42 y[2] = m_thermo->temperature();
45 m_thermo->getMassFractions(y+3);
49 getSurfaceInitialConditions(y + m_nsp + 3);
52 void IdealGasReactor::initialize(doublereal t0)
54 Reactor::initialize(t0);
55 m_uk.resize(m_nsp, 0.0);
58 void IdealGasReactor::updateState(doublereal* y)
65 m_thermo->setMassFractions_NoNorm(y+3);
66 m_thermo->setState_TR(y[2], m_mass / m_vol);
67 updateSurfaceState(y + m_nsp + 3);
70 m_enthalpy = m_thermo->enthalpy_mass();
71 m_pressure = m_thermo->pressure();
72 m_intEnergy = m_thermo->intEnergy_mass();
73 m_thermo->saveState(m_state);
76 void IdealGasReactor::evalEqs(doublereal time, doublereal* y,
77 doublereal* ydot, doublereal* params)
81 double* dYdt = ydot + 3;
83 m_thermo->restoreState(m_state);
84 applySensitivity(params);
85 m_thermo->getPartialMolarIntEnergies(&m_uk[0]);
86 const vector_fp& mw = m_thermo->molecularWeights();
87 const doublereal* Y = m_thermo->massFractions();
90 m_kin->getNetProductionRates(&m_wdot[0]);
94 double mdot_surf = evalSurfaces(time, ydot + m_nsp + 3);
98 mcvdTdt += - m_pressure * m_vdot - m_Q;
100 for (
size_t n = 0; n < m_nsp; n++) {
102 mcvdTdt -= m_wdot[n] * m_uk[n] * m_vol;
103 mcvdTdt -= m_sdot[n] * m_uk[n];
105 dYdt[n] = (m_wdot[n] * m_vol + m_sdot[n]) * mw[n] / m_mass;
107 dYdt[n] -= Y[n] * mdot_surf / m_mass;
111 for (
size_t i = 0; i < m_outlet.size(); i++) {
112 double mdot_out = m_outlet[i]->massFlowRate(time);
114 mcvdTdt -= mdot_out * m_pressure * m_vol / m_mass;
118 for (
size_t i = 0; i < m_inlet.size(); i++) {
119 double mdot_in = m_inlet[i]->massFlowRate(time);
121 mcvdTdt += m_inlet[i]->enthalpy_mass() * mdot_in;
122 for (
size_t n = 0; n < m_nsp; n++) {
123 double mdot_spec = m_inlet[i]->outletSpeciesMassFlowRate(n);
125 dYdt[n] += (mdot_spec - mdot_in * Y[n]) / m_mass;
129 mcvdTdt -= m_uk[n] / mw[n] * mdot_spec;
136 ydot[2] = mcvdTdt / (m_mass * m_thermo->cv_mass());
141 resetSensitivity(params);
144 size_t IdealGasReactor::componentIndex(
const string& nm)
const 146 size_t k = speciesIndex(nm);
149 }
else if (nm ==
"mass") {
151 }
else if (nm ==
"volume") {
153 }
else if (nm ==
"temperature") {
160 std::string IdealGasReactor::componentName(
size_t k) {
162 return "temperature";
164 return Reactor::componentName(k);
Header file for class Wall.
const size_t npos
index returned by functions to indicate "no position"
Base class for a phase with thermodynamic properties.
Base class for exceptions thrown by Cantera classes.
virtual std::string type() const
String indicating the thermodynamic model implemented.
std::vector< double > vector_fp
Turn on the use of stl vectors for the basic array type within cantera Vector of doubles.
Namespace for the Cantera kernel.