Class Reactor is a general-purpose class for stirred reactors. More...
#include <Reactor.h>
Class Reactor is a general-purpose class for stirred reactors.
The reactor may have an arbitrary number of inlets and outlets, each of which may be connected to a "flow device" such as a mass flow controller, a pressure regulator, etc. Additional reactors may be connected to the other end of the flow device, allowing construction of arbitrary reactor networks.
The reactor class integrates the same governing equations no matter what type of reactor is simulated. The differences among reactor types are completely specified by the attached flow devices and the time-dependent user-specified boundary conditions.
If an instance of class Reactor is used directly, it will simulate an adiabatic, constant volume reactor with gas-phase chemistry but no surface chemistry. Other reactor types may be simulated by deriving a class from Reactor. This method allows specifying the following in terms of the instantaneous reactor state:
See the Science Reference for the governing equations of class Reactor.
Public Member Functions | |
| Reactor (shared_ptr< Solution > sol, const string &name="(none)") | |
| Reactor (shared_ptr< Solution > sol, bool clone, const string &name="(none)") | |
| string | type () const override |
| String indicating the reactor model implemented. | |
| virtual bool | isOde () const |
| Indicate whether the governing equations for this reactor type are a system of ODEs or DAEs. | |
| virtual bool | timeIsIndependent () const |
| Indicates whether the governing equations for this reactor are functions of time or a spatial variable. | |
| void | setInitialVolume (double vol) override |
| Set the initial reactor volume. | |
| void | setChemistryEnabled (bool cflag=true) override |
| Enable or disable changes in reactor composition due to chemical reactions. | |
| bool | chemistryEnabled () const override |
Returns true if changes in the reactor composition due to chemical reactions are enabled. | |
| void | setEnergyEnabled (bool eflag=true) override |
| Set the energy equation on or off. | |
| bool | energyEnabled () const override |
Returns true if solution of the energy equation is enabled. | |
| size_t | neq () |
| Number of equations (state variables) for this reactor. | |
| virtual void | getState (double *y) |
| Get the the current state of the reactor. | |
| virtual void | getStateDae (double *y, double *ydot) |
| Get the current state and derivative vector of the reactor for a DAE solver. | |
| void | initialize (double t0=0.0) override |
| Initialize the reactor. | |
| virtual void | eval (double t, double *LHS, double *RHS) |
| Evaluate the reactor governing equations. | |
| virtual void | evalDae (double t, double *y, double *ydot, double *residual) |
| Evaluate the reactor governing equations. | |
| virtual void | getConstraints (double *constraints) |
| Given a vector of length neq(), mark which variables should be considered algebraic constraints. | |
| virtual vector< size_t > | steadyConstraints () const |
| Get the indices of equations that are algebraic constraints when solving the steady-state problem. | |
| virtual void | updateState (double *y) |
| Set the state of the reactor to correspond to the state vector y. | |
| size_t | nSensParams () const override |
| Number of sensitivity parameters associated with this reactor. | |
| void | addSensitivityReaction (size_t rxn) override |
| Add a sensitivity parameter associated with the reaction number rxn | |
| virtual void | addSensitivitySpeciesEnthalpy (size_t k) |
| Add a sensitivity parameter associated with the enthalpy formation of species k (in the homogeneous phase) | |
| virtual size_t | componentIndex (const string &nm) const |
| Return the index in the solution vector for this reactor of the component named nm. | |
| virtual string | componentName (size_t k) |
| Return the name of the solution component with index i. | |
| virtual double | upperBound (size_t k) const |
| Get the upper bound on the k-th component of the local state vector. | |
| virtual double | lowerBound (size_t k) const |
| Get the lower bound on the k-th component of the local state vector. | |
| virtual void | resetBadValues (double *y) |
| Reset physically or mathematically problematic values, such as negative species concentrations. | |
| void | setAdvanceLimits (const double *limits) |
| Set absolute step size limits during advance. | |
| bool | hasAdvanceLimits () const |
| Check whether Reactor object uses advance limits. | |
| bool | getAdvanceLimits (double *limits) const |
| Retrieve absolute step size limits during advance. | |
| void | setAdvanceLimit (const string &nm, const double limit) |
| Set individual step size limit for component name nm | |
| virtual Eigen::SparseMatrix< double > | jacobian () |
| Calculate the Jacobian of a specific Reactor specialization. | |
| Eigen::SparseMatrix< double > | finiteDifferenceJacobian () |
| Calculate the reactor-specific Jacobian using a finite difference method. | |
| virtual void | setDerivativeSettings (AnyMap &settings) |
| Use this to set the kinetics objects derivative settings. | |
| virtual void | applySensitivity (double *params) |
| Set reaction rate multipliers based on the sensitivity variables in params. | |
| virtual void | resetSensitivity (double *params) |
| Reset the reaction rate multipliers. | |
| virtual bool | preconditionerSupported () const |
| Return a false if preconditioning is not supported or true otherwise. | |
Public Member Functions inherited from ReactorBase | |
| ReactorBase (shared_ptr< Solution > sol, const string &name="(none)") | |
| Instantiate a ReactorBase object with Solution contents. | |
| ReactorBase (shared_ptr< Solution > sol, bool clone, const string &name="(none)") | |
| Instantiate a ReactorBase object with Solution contents. | |
| ReactorBase (const ReactorBase &)=delete | |
| ReactorBase & | operator= (const ReactorBase &)=delete |
| virtual string | type () const |
| String indicating the reactor model implemented. | |
| string | name () const |
| Return the name of this reactor. | |
| void | setName (const string &name) |
| Set the name of this reactor. | |
| bool | setDefaultName (map< string, int > &counts) |
Set the default name of a reactor. Returns false if it was previously set. | |
| shared_ptr< Solution > | phase () |
| Access the Solution object used to represent the contents of this reactor. | |
| shared_ptr< const Solution > | phase () const |
| Access the Solution object used to represent the contents of this reactor. | |
| virtual void | restoreState () |
| Set the state of the Phase object associated with this reactor to the reactor's current state. | |
| virtual void | syncState () |
| Set the state of the reactor to the associated ThermoPhase object. | |
| double | residenceTime () |
| Return the residence time (s) of the contents of this reactor, based on the outlet mass flow rates and the mass of the reactor contents. | |
| ReactorNet & | network () |
| The ReactorNet that this reactor belongs to. | |
| void | setNetwork (ReactorNet *net) |
| Set the ReactorNet that this reactor belongs to. | |
| virtual void | addSensitivityReaction (size_t rxn) |
| Add a sensitivity parameter associated with the reaction number rxn | |
| virtual size_t | nSensParams () const |
| Number of sensitivity parameters associated with this reactor. | |
| virtual double | area () const |
| Returns an area associated with a reactor [m²]. | |
| virtual void | setArea (double a) |
| Set an area associated with a reactor [m²]. | |
| virtual void | addInlet (FlowDevice &inlet) |
| Connect an inlet FlowDevice to this reactor. | |
| virtual void | addOutlet (FlowDevice &outlet) |
| Connect an outlet FlowDevice to this reactor. | |
| FlowDevice & | inlet (size_t n=0) |
| Return a reference to the n-th inlet FlowDevice connected to this reactor. | |
| FlowDevice & | outlet (size_t n=0) |
| Return a reference to the n-th outlet FlowDevice connected to this reactor. | |
| size_t | nInlets () |
| Return the number of inlet FlowDevice objects connected to this reactor. | |
| size_t | nOutlets () |
| Return the number of outlet FlowDevice objects connected to this reactor. | |
| size_t | nWalls () |
| Return the number of Wall objects connected to this reactor. | |
| virtual void | addWall (WallBase &w, int lr) |
| Insert a Wall between this reactor and another reactor. | |
| WallBase & | wall (size_t n) |
| Return a reference to the n-th Wall connected to this reactor. | |
| virtual void | addSurface (ReactorSurface *surf) |
| Add a ReactorSurface object to a Reactor object. | |
| ReactorSurface * | surface (size_t n) |
| Return a reference to the n-th ReactorSurface connected to this reactor. | |
| virtual size_t | nSurfs () const |
| Return the number of surfaces in a reactor. | |
| double | volume () const |
| Returns the current volume (m^3) of the reactor. | |
| double | density () const |
| Returns the current density (kg/m^3) of the reactor's contents. | |
| double | temperature () const |
| Returns the current temperature (K) of the reactor's contents. | |
| double | enthalpy_mass () const |
| Returns the current enthalpy (J/kg) of the reactor's contents. | |
| double | pressure () const |
| Returns the current pressure (Pa) of the reactor. | |
| double | mass () const |
| Returns the mass (kg) of the reactor's contents. | |
| const double * | massFractions () const |
| Return the vector of species mass fractions. | |
| double | massFraction (size_t k) const |
| Return the mass fraction of the k-th species. | |
Protected Member Functions | |
| virtual size_t | speciesIndex (const string &nm) const |
| Return the index in the solution vector for this reactor of the species named nm, in either the homogeneous phase or a surface phase, relative to the start of the species terms. | |
| virtual void | evalWalls (double t) |
| Evaluate terms related to Walls. | |
| virtual void | evalSurfaces (double *LHS, double *RHS, double *sdot) |
| Evaluate terms related to surface reactions. | |
| virtual void | evalSurfaces (double *RHS, double *sdot) |
| virtual void | updateSurfaceState (double *y) |
| Update the state of SurfPhase objects attached to this reactor. | |
| virtual void | updateConnected (bool updatePressure) |
| Update the state information needed by connected reactors, flow devices, and reactor walls. | |
| virtual void | getSurfaceInitialConditions (double *y) |
| Get initial conditions for SurfPhase objects attached to this reactor. | |
Protected Member Functions inherited from ReactorBase | |
| ReactorBase (const string &name="(none)") | |
Protected Attributes | |
| Kinetics * | m_kin = nullptr |
| Pointer to the homogeneous Kinetics object that handles the reactions. | |
| double | m_vdot = 0.0 |
| net rate of volume change from moving walls [m^3/s] | |
| double | m_Qdot = 0.0 |
| net heat transfer into the reactor, through walls [W] | |
| vector< double > | m_work |
| vector< double > | m_sdot |
| Production rates of gas phase species on surfaces [kmol/s]. | |
| vector< double > | m_wdot |
| Species net molar production rates. | |
| vector< double > | m_uk |
| Species molar internal energies. | |
| bool | m_chem = false |
| bool | m_energy = true |
| size_t | m_nv = 0 |
| size_t | m_nv_surf |
| vector< double > | m_advancelimits |
| !< Number of variables associated with reactor surfaces | |
| vector< Eigen::Triplet< double > > | m_jac_trips |
| Vector of triplets representing the jacobian. | |
Protected Attributes inherited from ReactorBase | |
| size_t | m_nsp = 0 |
| Number of homogeneous species in the mixture. | |
| ThermoPhase * | m_thermo = nullptr |
| double | m_vol = 0.0 |
| Current volume of the reactor [m^3]. | |
| double | m_mass = 0.0 |
| Current mass of the reactor [kg]. | |
| double | m_enthalpy = 0.0 |
| Current specific enthalpy of the reactor [J/kg]. | |
| double | m_pressure = 0.0 |
| Current pressure in the reactor [Pa]. | |
| vector< double > | m_state |
| vector< FlowDevice * > | m_inlet |
| vector< FlowDevice * > | m_outlet |
| vector< WallBase * > | m_wall |
| vector< ReactorSurface * > | m_surfaces |
| vector< int > | m_lr |
| Vector of length nWalls(), indicating whether this reactor is on the left (0) or right (1) of each wall. | |
| string | m_name |
| Reactor name. | |
| bool | m_defaultNameSet = false |
true if default name has been previously set. | |
| ReactorNet * | m_net = nullptr |
| The ReactorNet that this reactor is part of. | |
| shared_ptr< Solution > | m_solution |
| Composite thermo/kinetics/transport handler. | |
| vector< SensitivityParameter > | m_sensParams |
Definition at line 25 of file Reactor.cpp.
Definition at line 30 of file Reactor.cpp.
|
inlineoverridevirtual |
String indicating the reactor model implemented.
Usually corresponds to the name of the derived class.
Reimplemented from ReactorBase.
|
inlinevirtual |
Indicate whether the governing equations for this reactor type are a system of ODEs or DAEs.
In the first case, this class implements the eval() method. In the second case, this class implements the evalDae() method.
Reimplemented in FlowReactor.
|
inlinevirtual |
Indicates whether the governing equations for this reactor are functions of time or a spatial variable.
All reactors in a network must have the same value.
Reimplemented in FlowReactor.
|
inlineoverridevirtual |
Set the initial reactor volume.
Reimplemented from ReactorBase.
|
inlineoverridevirtual |
Enable or disable changes in reactor composition due to chemical reactions.
Reimplemented from ReactorBase.
|
inlineoverridevirtual |
Returns true if changes in the reactor composition due to chemical reactions are enabled.
Reimplemented from ReactorBase.
|
inlineoverridevirtual |
|
inlineoverridevirtual |
Returns true if solution of the energy equation is enabled.
Reimplemented from ReactorBase.
|
inline |
|
virtual |
Get the the current state of the reactor.
| [out] | y | state vector representing the initial state of the reactor |
Reimplemented in ConstPressureMoleReactor, ConstPressureReactor, FlowReactor, IdealGasConstPressureMoleReactor, IdealGasConstPressureReactor, IdealGasMoleReactor, IdealGasReactor, and MoleReactor.
Definition at line 47 of file Reactor.cpp.
|
inlinevirtual |
Get the current state and derivative vector of the reactor for a DAE solver.
| [out] | y | state vector representing the initial state of the reactor |
| [out] | ydot | state vector representing the initial derivatives of the reactor |
Reimplemented in FlowReactor.
|
overridevirtual |
Initialize the reactor.
Called automatically by ReactorNet::initialize.
Reimplemented from ReactorBase.
Definition at line 82 of file Reactor.cpp.
|
virtual |
Evaluate the reactor governing equations.
Called by ReactorNet::eval.
| [in] | t | time. |
| [out] | LHS | pointer to start of vector of left-hand side coefficients for governing equations, length m_nv, default values 1 |
| [out] | RHS | pointer to start of vector of right-hand side coefficients for governing equations, length m_nv, default values 0 |
Reimplemented in ConstPressureMoleReactor, ConstPressureReactor, FlowReactor, IdealGasConstPressureMoleReactor, IdealGasConstPressureReactor, IdealGasMoleReactor, IdealGasReactor, and MoleReactor.
Definition at line 202 of file Reactor.cpp.
|
inlinevirtual |
Evaluate the reactor governing equations.
Called by ReactorNet::eval.
| [in] | t | time. |
| [in] | y | solution vector, length neq() |
| [in] | ydot | rate of change of solution vector, length neq() |
| [out] | residual | residuals vector, length neq() |
Reimplemented in FlowReactor.
|
inlinevirtual |
Given a vector of length neq(), mark which variables should be considered algebraic constraints.
Reimplemented in FlowReactor.
|
virtual |
Get the indices of equations that are algebraic constraints when solving the steady-state problem.
Reimplemented in ConstPressureMoleReactor, ConstPressureReactor, FlowReactor, IdealGasConstPressureReactor, IdealGasMoleReactor, and IdealGasReactor.
Definition at line 307 of file Reactor.cpp.
|
virtual |
Set the state of the reactor to correspond to the state vector y.
Reimplemented in ConstPressureMoleReactor, ConstPressureReactor, FlowReactor, IdealGasConstPressureMoleReactor, IdealGasConstPressureReactor, IdealGasMoleReactor, IdealGasReactor, and MoleReactor.
Definition at line 119 of file Reactor.cpp.
|
overridevirtual |
Number of sensitivity parameters associated with this reactor.
(including walls)
Reimplemented from ReactorBase.
Definition at line 110 of file Reactor.cpp.
|
overridevirtual |
Add a sensitivity parameter associated with the reaction number rxn
Reimplemented from ReactorBase.
Definition at line 405 of file Reactor.cpp.
|
virtual |
Add a sensitivity parameter associated with the enthalpy formation of species k (in the homogeneous phase)
Definition at line 418 of file Reactor.cpp.
|
virtual |
Return the index in the solution vector for this reactor of the component named nm.
Possible values for nm are "mass", "volume", "int_energy", the name of a homogeneous phase species, or the name of a surface species.
Reimplemented in ConstPressureMoleReactor, ConstPressureReactor, FlowReactor, IdealGasConstPressureMoleReactor, IdealGasConstPressureReactor, IdealGasMoleReactor, IdealGasReactor, and MoleReactor.
Definition at line 456 of file Reactor.cpp.
|
virtual |
Return the name of the solution component with index i.
Reimplemented in ConstPressureMoleReactor, ConstPressureReactor, FlowReactor, IdealGasConstPressureMoleReactor, IdealGasConstPressureReactor, IdealGasMoleReactor, IdealGasReactor, and MoleReactor.
Definition at line 475 of file Reactor.cpp.
|
virtual |
Get the upper bound on the k-th component of the local state vector.
Reimplemented in ConstPressureMoleReactor, ConstPressureReactor, IdealGasConstPressureMoleReactor, IdealGasConstPressureReactor, IdealGasMoleReactor, IdealGasReactor, and MoleReactor.
Definition at line 501 of file Reactor.cpp.
|
virtual |
Get the lower bound on the k-th component of the local state vector.
Reimplemented in ConstPressureMoleReactor, ConstPressureReactor, IdealGasConstPressureMoleReactor, IdealGasConstPressureReactor, IdealGasMoleReactor, IdealGasReactor, and MoleReactor.
Definition at line 515 of file Reactor.cpp.
|
virtual |
Reset physically or mathematically problematic values, such as negative species concentrations.
| [in,out] | y | current state vector, to be updated; length neq() |
Reimplemented in ConstPressureMoleReactor, ConstPressureReactor, and MoleReactor.
Definition at line 529 of file Reactor.cpp.
| void setAdvanceLimits | ( | const double * | limits | ) |
Set absolute step size limits during advance.
| limits | array of step size limits with length neq |
Definition at line 578 of file Reactor.cpp.
|
inline |
| bool getAdvanceLimits | ( | double * | limits | ) | const |
Retrieve absolute step size limits during advance.
| [out] | limits | array of step size limits with length neq |
Definition at line 593 of file Reactor.cpp.
| void setAdvanceLimit | ( | const string & | nm, |
| const double | limit | ||
| ) |
Set individual step size limit for component name nm
| nm | component name |
| limit | value for step size limit |
Definition at line 604 of file Reactor.cpp.
| Eigen::SparseMatrix< double > finiteDifferenceJacobian | ( | ) |
Calculate the reactor-specific Jacobian using a finite difference method.
This method is used only for informational purposes. Jacobian calculations for the full reactor system are handled internally by CVODES.
Definition at line 323 of file Reactor.cpp.
|
virtual |
Use this to set the kinetics objects derivative settings.
Definition at line 38 of file Reactor.cpp.
|
virtual |
Set reaction rate multipliers based on the sensitivity variables in params.
Definition at line 535 of file Reactor.cpp.
|
virtual |
Reset the reaction rate multipliers.
Definition at line 557 of file Reactor.cpp.
|
inlinevirtual |
Return a false if preconditioning is not supported or true otherwise.
Reimplemented in IdealGasConstPressureMoleReactor, and IdealGasMoleReactor.
|
protectedvirtual |
Return the index in the solution vector for this reactor of the species named nm, in either the homogeneous phase or a surface phase, relative to the start of the species terms.
Used to implement componentIndex for specific reactor implementations.
Definition at line 433 of file Reactor.cpp.
|
protectedvirtual |
Evaluate terms related to Walls.
Calculates m_vdot and m_Qdot based on wall movement and heat transfer.
| t | the current time |
Definition at line 266 of file Reactor.cpp.
|
protectedvirtual |
Evaluate terms related to surface reactions.
| [out] | LHS | Multiplicative factor on the left hand side of ODE for surface species coverages |
| [out] | RHS | Right hand side of ODE for surface species coverages |
| [out] | sdot | array of production rates of bulk phase species on surfaces [kmol/s] |
Reimplemented in MoleReactor.
Definition at line 278 of file Reactor.cpp.
|
protectedvirtual |
Definition at line 376 of file Reactor.cpp.
|
protectedvirtual |
Update the state of SurfPhase objects attached to this reactor.
Reimplemented in FlowReactor, and MoleReactor.
Definition at line 167 of file Reactor.cpp.
|
protectedvirtual |
Update the state information needed by connected reactors, flow devices, and reactor walls.
Called from updateState().
| updatePressure | Indicates whether to update m_pressure. Should true for reactors where the pressure is a dependent property, calculated from the state, and false when the pressure is constant or an independent variable. |
Definition at line 176 of file Reactor.cpp.
|
protectedvirtual |
Get initial conditions for SurfPhase objects attached to this reactor.
Reimplemented in MoleReactor.
Definition at line 73 of file Reactor.cpp.
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |