Cantera 2.6.0
IdealGasReactor.cpp
Go to the documentation of this file.
1//! @file IdealGasReactor.cpp A zero-dimensional reactor
2
3// This file is part of Cantera. See License.txt in the top-level directory or
4// at https://cantera.org/license.txt for license and copyright information.
5
12
13using namespace std;
14
15namespace Cantera
16{
17
19{
20 //! @todo: Add a method to ThermoPhase that indicates whether a given
21 //! subclass is compatible with this reactor model
22 if (thermo.type() != "IdealGas") {
23 throw CanteraError("IdealGasReactor::setThermoMgr",
24 "Incompatible phase type provided");
25 }
27}
28
30{
31 if (m_thermo == 0) {
32 throw CanteraError("IdealGasReactor::getState",
33 "Error: reactor is empty.");
34 }
35 m_thermo->restoreState(m_state);
36
37 // set the first component to the total mass
38 m_mass = m_thermo->density() * m_vol;
39 y[0] = m_mass;
40
41 // set the second component to the total volume
42 y[1] = m_vol;
43
44 // Set the third component to the temperature
45 y[2] = m_thermo->temperature();
46
47 // set components y+3 ... y+K+2 to the mass fractions of each species
48 m_thermo->getMassFractions(y+3);
49
50 // set the remaining components to the surface species
51 // coverages on the walls
53}
54
55void IdealGasReactor::initialize(doublereal t0)
56{
58 m_uk.resize(m_nsp, 0.0);
59}
60
62{
63 // The components of y are [0] the total mass, [1] the total volume,
64 // [2] the temperature, [3...K+3] are the mass fractions of each species,
65 // and [K+3...] are the coverages of surface species on each wall.
66 m_mass = y[0];
67 m_vol = y[1];
68 m_thermo->setMassFractions_NoNorm(y+3);
69 m_thermo->setState_TR(y[2], m_mass / m_vol);
70 updateConnected(true);
72}
73
74void IdealGasReactor::eval(double time, double* LHS, double* RHS)
75{
76 double& dmdt = RHS[0]; // dm/dt (gas phase)
77 double& mcvdTdt = RHS[2]; // m * c_v * dT/dt
78 double* mdYdt = RHS + 3; // mass * dY/dt
79
80 evalWalls(time);
81 m_thermo->restoreState(m_state);
82 m_thermo->getPartialMolarIntEnergies(&m_uk[0]);
83 const vector_fp& mw = m_thermo->molecularWeights();
84 const doublereal* Y = m_thermo->massFractions();
85
86 if (m_chem) {
87 m_kin->getNetProductionRates(&m_wdot[0]); // "omega dot"
88 }
89
90 evalSurfaces(LHS + m_nsp + 3, RHS + m_nsp + 3, m_sdot.data());
91 double mdot_surf = dot(m_sdot.begin(), m_sdot.end(), mw.begin());
92 dmdt += mdot_surf;
93
94 // compression work and external heat transfer
95 mcvdTdt += - m_pressure * m_vdot + m_Qdot;
96
97 for (size_t n = 0; n < m_nsp; n++) {
98 // heat release from gas phase and surface reactions
99 mcvdTdt -= m_wdot[n] * m_uk[n] * m_vol;
100 mcvdTdt -= m_sdot[n] * m_uk[n];
101 // production in gas phase and from surfaces
102 mdYdt[n] = (m_wdot[n] * m_vol + m_sdot[n]) * mw[n];
103 // dilution by net surface mass flux
104 mdYdt[n] -= Y[n] * mdot_surf;
105 //Assign left-hand side of dYdt ODE as total mass
106 LHS[n+3] = m_mass;
107 }
108
109 // add terms for outlets
110 for (auto outlet : m_outlet) {
111 double mdot = outlet->massFlowRate();
112 dmdt -= mdot; // mass flow out of system
113 mcvdTdt -= mdot * m_pressure * m_vol / m_mass; // flow work
114 }
115
116 // add terms for inlets
117 for (auto inlet : m_inlet) {
118 double mdot = inlet->massFlowRate();
119 dmdt += mdot; // mass flow into system
120 mcvdTdt += inlet->enthalpy_mass() * mdot;
121 for (size_t n = 0; n < m_nsp; n++) {
122 double mdot_spec = inlet->outletSpeciesMassFlowRate(n);
123 // flow of species into system and dilution by other species
124 mdYdt[n] += mdot_spec - mdot * Y[n];
125
126 // In combination with h_in*mdot_in, flow work plus thermal
127 // energy carried with the species
128 mcvdTdt -= m_uk[n] / mw[n] * mdot_spec;
129 }
130 }
131
132 RHS[1] = m_vdot;
133 if (m_energy) {
134 LHS[2] = m_mass * m_thermo->cv_mass();
135 } else {
136 RHS[2] = 0;
137 }
138}
139
140size_t IdealGasReactor::componentIndex(const string& nm) const
141{
142 size_t k = speciesIndex(nm);
143 if (k != npos) {
144 return k + 3;
145 } else if (nm == "mass") {
146 return 0;
147 } else if (nm == "volume") {
148 return 1;
149 } else if (nm == "temperature") {
150 return 2;
151 } else {
152 return npos;
153 }
154}
155
156std::string IdealGasReactor::componentName(size_t k) {
157 if (k == 2) {
158 return "temperature";
159 } else {
160 return Reactor::componentName(k);
161 }
162}
163
164
165}
Base class for kinetics managers and also contains the kineticsmgr module documentation (see Kinetics...
Header file for class ThermoPhase, the base class for phases with thermodynamic properties,...
Header file for base class WallBase.
Base class for exceptions thrown by Cantera classes.
Definition: ctexceptions.h:61
double outletSpeciesMassFlowRate(size_t k)
Mass flow rate (kg/s) of outlet species k.
Definition: FlowDevice.cpp:59
double enthalpy_mass()
specific enthalpy
Definition: FlowDevice.cpp:71
double massFlowRate()
Mass flow rate (kg/s).
Definition: FlowDevice.h:48
virtual size_t componentIndex(const std::string &nm) const
Return the index in the solution vector for this reactor of the component named nm.
virtual void updateState(doublereal *y)
Set the state of the reactor to correspond to the state vector y.
virtual void eval(double t, double *LHS, double *RHS)
Evaluate the reactor governing equations.
virtual void getState(doublereal *y)
Get the the current state of the reactor.
std::string componentName(size_t k)
Return the name of the solution component with index i.
virtual void initialize(doublereal t0=0.0)
Initialize the reactor.
virtual void setThermoMgr(ThermoPhase &thermo)
Specify the mixture contained in the reactor.
vector_fp m_uk
Species molar internal energies.
virtual void getNetProductionRates(doublereal *wdot)
Species net production rates [kmol/m^3/s or kmol/m^2/s].
Definition: Kinetics.cpp:484
virtual void setMassFractions_NoNorm(const double *const y)
Set the mass fractions to the specified values without normalizing.
Definition: Phase.cpp:405
const vector_fp & molecularWeights() const
Return a const reference to the internal vector of molecular weights.
Definition: Phase.cpp:509
const double * massFractions() const
Return a const pointer to the mass fraction array.
Definition: Phase.h:532
void setState_TR(doublereal t, doublereal rho)
Set the internally stored temperature (K) and density (kg/m^3)
Definition: Phase.cpp:462
virtual double density() const
Density (kg/m^3).
Definition: Phase.h:679
doublereal temperature() const
Temperature (K).
Definition: Phase.h:654
void restoreState(const vector_fp &state)
Restore a state saved on a previous call to saveState.
Definition: Phase.cpp:310
void getMassFractions(double *const y) const
Get the species mass fractions.
Definition: Phase.cpp:585
FlowDevice & outlet(size_t n=0)
Return a reference to the n-th outlet FlowDevice connected to this reactor.
double m_pressure
Current pressure in the reactor [Pa].
Definition: ReactorBase.h:263
FlowDevice & inlet(size_t n=0)
Return a reference to the n-th inlet FlowDevice connected to this reactor.
double m_vol
Current volume of the reactor [m^3].
Definition: ReactorBase.h:260
size_t m_nsp
Number of homogeneous species in the mixture.
Definition: ReactorBase.h:257
virtual void setThermoMgr(ThermoPhase &thermo)
Specify the mixture contained in the reactor.
Definition: ReactorBase.cpp:27
virtual void evalSurfaces(double *LHS, double *RHS, double *sdot)
Evaluate terms related to surface reactions.
Definition: Reactor.cpp:279
virtual void updateSurfaceState(double *y)
Update the state of SurfPhase objects attached to this reactor.
Definition: Reactor.cpp:175
Kinetics * m_kin
Pointer to the homogeneous Kinetics object that handles the reactions.
Definition: Reactor.h:216
virtual void evalWalls(double t)
Evaluate terms related to Walls.
Definition: Reactor.cpp:267
double m_Qdot
net heat transfer into the reactor, through walls [W]
Definition: Reactor.h:225
doublereal m_mass
total mass
Definition: Reactor.h:227
virtual std::string componentName(size_t k)
Return the name of the solution component with index i.
Definition: Reactor.cpp:378
vector_fp m_sdot
Production rates of gas phase species on surfaces [kmol/s].
Definition: Reactor.h:231
virtual void initialize(doublereal t0=0.0)
Initialize the reactor.
Definition: Reactor.cpp:84
virtual void getSurfaceInitialConditions(double *y)
Get initial conditions for SurfPhase objects attached to this reactor.
Definition: Reactor.cpp:75
virtual size_t speciesIndex(const std::string &nm) const
Return the index in the solution vector for this reactor of the species named nm, in either the homog...
Definition: Reactor.cpp:340
vector_fp m_wdot
Species net molar production rates.
Definition: Reactor.h:233
virtual void updateConnected(bool updatePressure)
Update the state information needed by connected reactors and flow devices.
Definition: Reactor.cpp:184
doublereal m_vdot
net rate of volume change from moving walls [m^3/s]
Definition: Reactor.h:218
Base class for a phase with thermodynamic properties.
Definition: ThermoPhase.h:102
virtual void getPartialMolarIntEnergies(doublereal *ubar) const
Return an array of partial molar internal energies for the species in the mixture.
Definition: ThermoPhase.h:541
virtual std::string type() const
String indicating the thermodynamic model implemented.
Definition: ThermoPhase.h:113
doublereal cv_mass() const
Specific heat at constant volume. Units: J/kg/K.
Definition: ThermoPhase.h:773
Namespace for the Cantera kernel.
Definition: AnyMap.h:29
const size_t npos
index returned by functions to indicate "no position"
Definition: ct_defs.h:192
doublereal dot(InputIter x_begin, InputIter x_end, InputIter2 y_begin)
Function that calculates a templated inner product.
Definition: utilities.h:77
std::vector< double > vector_fp
Turn on the use of stl vectors for the basic array type within cantera Vector of doubles.
Definition: ct_defs.h:184
Various templated functions that carry out common vector operations (see Templated Utility Functions)...