Cantera  3.1.0a1

Class Reactor is a general-purpose class for stirred reactors. More...

#include <Reactor.h>

Inheritance diagram for Reactor:
[legend]

Detailed Description

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:

  • rate of change of the total volume (m^3/s)
  • surface heat loss rate (W)
  • species surface production rates (kmol/s)

Definition at line 43 of file Reactor.h.

Public Member Functions

string type () const override
 String indicating the reactor model implemented. More...
 
virtual bool isOde () const
 Indicate whether the governing equations for this reactor type are a system of ODEs or DAEs. More...
 
virtual bool timeIsIndependent () const
 Indicates whether the governing equations for this reactor are functions of time or a spatial variable. More...
 
template<class G >
void insert (G &contents)
 Insert something into the reactor. More...
 
void insert (shared_ptr< Solution > sol)
 
void setKineticsMgr (Kinetics &kin) override
 Specify chemical kinetics governing the reactor. More...
 
void setChemistry (bool cflag=true) override
 Enable or disable changes in reactor composition due to chemical reactions. More...
 
bool chemistryEnabled () const
 Returns true if changes in the reactor composition due to chemical reactions are enabled. More...
 
void setEnergy (int eflag=1) override
 Set the energy equation on or off. More...
 
bool energyEnabled () const
 Returns true if solution of the energy equation is enabled. More...
 
size_t neq ()
 Number of equations (state variables) for this reactor. More...
 
virtual void getState (double *y)
 Get the the current state of the reactor. More...
 
virtual void getStateDae (double *y, double *ydot)
 Get the current state and derivative vector of the reactor for a DAE solver. More...
 
void initialize (double t0=0.0) override
 Initialize the reactor. More...
 
virtual void eval (double t, double *LHS, double *RHS)
 Evaluate the reactor governing equations. More...
 
virtual void evalDae (double t, double *y, double *ydot, double *residual)
 Evaluate the reactor governing equations. More...
 
virtual void getConstraints (double *constraints)
 Given a vector of length neq(), mark which variables should be considered algebraic constraints. More...
 
void syncState () override
 Set the state of the reactor to correspond to the state of the associated ThermoPhase object. More...
 
virtual void updateState (double *y)
 Set the state of the reactor to correspond to the state vector y. More...
 
virtual size_t nSensParams () const
 Number of sensitivity parameters associated with this reactor (including walls) More...
 
virtual void addSensitivityReaction (size_t rxn)
 Add a sensitivity parameter associated with the reaction number rxn (in the homogeneous phase). More...
 
virtual void addSensitivitySpeciesEnthalpy (size_t k)
 Add a sensitivity parameter associated with the enthalpy formation of species k (in the homogeneous phase) More...
 
virtual size_t componentIndex (const string &nm) const
 Return the index in the solution vector for this reactor of the component named nm. More...
 
virtual string componentName (size_t k)
 Return the name of the solution component with index i. More...
 
void setAdvanceLimits (const double *limits)
 Set absolute step size limits during advance. More...
 
bool hasAdvanceLimits () const
 Check whether Reactor object uses advance limits. More...
 
bool getAdvanceLimits (double *limits) const
 Retrieve absolute step size limits during advance. More...
 
void setAdvanceLimit (const string &nm, const double limit)
 Set individual step size limit for component name nm More...
 
virtual Eigen::SparseMatrix< double > jacobian ()
 Calculate the Jacobian of a specific Reactor specialization. More...
 
Eigen::SparseMatrix< double > finiteDifferenceJacobian ()
 Calculate the reactor-specific Jacobian using a finite difference method. More...
 
virtual void setDerivativeSettings (AnyMap &settings)
 Use this to set the kinetics objects derivative settings. More...
 
virtual void applySensitivity (double *params)
 Set reaction rate multipliers based on the sensitivity variables in params. More...
 
virtual void resetSensitivity (double *params)
 Reset the reaction rate multipliers. More...
 
virtual bool preconditionerSupported () const
 Return a false if preconditioning is not supported or true otherwise. More...
 
- Public Member Functions inherited from ReactorBase
 ReactorBase (const string &name="(none)")
 
 ReactorBase (const ReactorBase &)=delete
 
ReactorBaseoperator= (const ReactorBase &)=delete
 
string name () const
 Return the name of this reactor. More...
 
void setName (const string &name)
 Set the name of this reactor. More...
 
void restoreState ()
 Set the state of the Phase object associated with this reactor to the reactor's current state. More...
 
ThermoPhasecontents ()
 return a reference to the contents. More...
 
const ThermoPhasecontents () const
 
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. More...
 
ReactorNetnetwork ()
 The ReactorNet that this reactor belongs to. More...
 
void setNetwork (ReactorNet *net)
 Set the ReactorNet that this reactor belongs to. More...
 
void setInitialVolume (double vol)
 Set the initial reactor volume. By default, the volume is 1.0 m^3. More...
 
virtual void setThermoMgr (ThermoPhase &thermo)
 Specify the mixture contained in the reactor. More...
 
void addInlet (FlowDevice &inlet)
 Connect an inlet FlowDevice to this reactor. More...
 
void addOutlet (FlowDevice &outlet)
 Connect an outlet FlowDevice to this reactor. More...
 
FlowDeviceinlet (size_t n=0)
 Return a reference to the n-th inlet FlowDevice connected to this reactor. More...
 
FlowDeviceoutlet (size_t n=0)
 Return a reference to the n-th outlet FlowDevice connected to this reactor. More...
 
size_t nInlets ()
 Return the number of inlet FlowDevice objects connected to this reactor. More...
 
size_t nOutlets ()
 Return the number of outlet FlowDevice objects connected to this reactor. More...
 
size_t nWalls ()
 Return the number of Wall objects connected to this reactor. More...
 
void addWall (WallBase &w, int lr)
 Insert a Wall between this reactor and another reactor. More...
 
WallBasewall (size_t n)
 Return a reference to the n-th Wall connected to this reactor. More...
 
virtual void addSurface (ReactorSurface *surf)
 
ReactorSurfacesurface (size_t n)
 Return a reference to the n-th ReactorSurface connected to this reactor. More...
 
virtual size_t nSurfs ()
 Return the number of surfaces in a reactor. More...
 
double volume () const
 Returns the current volume (m^3) of the reactor. More...
 
double density () const
 Returns the current density (kg/m^3) of the reactor's contents. More...
 
double temperature () const
 Returns the current temperature (K) of the reactor's contents. More...
 
double enthalpy_mass () const
 Returns the current enthalpy (J/kg) of the reactor's contents. More...
 
double intEnergy_mass () const
 Returns the current internal energy (J/kg) of the reactor's contents. More...
 
double pressure () const
 Returns the current pressure (Pa) of the reactor. More...
 
double mass () const
 Returns the mass (kg) of the reactor's contents. More...
 
const double * massFractions () const
 Return the vector of species mass fractions. More...
 
double massFraction (size_t k) const
 Return the mass fraction of the k-th species. More...
 

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. More...
 
virtual void evalWalls (double t)
 Evaluate terms related to Walls. More...
 
virtual void evalSurfaces (double *LHS, double *RHS, double *sdot)
 Evaluate terms related to surface reactions. More...
 
virtual void evalSurfaces (double *RHS, double *sdot)
 
virtual void updateSurfaceState (double *y)
 Update the state of SurfPhase objects attached to this reactor. More...
 
virtual void updateConnected (bool updatePressure)
 Update the state information needed by connected reactors, flow devices, and reactor walls. More...
 
virtual void getSurfaceInitialConditions (double *y)
 Get initial conditions for SurfPhase objects attached to this reactor. More...
 

Protected Attributes

Kineticsm_kin = nullptr
 Pointer to the homogeneous Kinetics object that handles the reactions. More...
 
double m_vdot = 0.0
 net rate of volume change from moving walls [m^3/s] More...
 
double m_Qdot = 0.0
 net heat transfer into the reactor, through walls [W] More...
 
double m_mass = 0.0
 total mass More...
 
vector< double > m_work
 
vector< double > m_sdot
 Production rates of gas phase species on surfaces [kmol/s]. More...
 
vector< double > m_wdot
 Species net molar production rates. More...
 
vector< double > m_uk
 Species molar internal energies. More...
 
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 More...
 
vector< SensitivityParameter > m_sensParams
 
vector< Eigen::Triplet< double > > m_jac_trips
 Vector of triplets representing the jacobian. More...
 
- Protected Attributes inherited from ReactorBase
size_t m_nsp = 0
 Number of homogeneous species in the mixture. More...
 
ThermoPhasem_thermo = nullptr
 
double m_vol = 1.0
 Current volume of the reactor [m^3]. More...
 
double m_enthalpy = 0.0
 Current specific enthalpy of the reactor [J/kg]. More...
 
double m_intEnergy = 0.0
 Current internal energy of the reactor [J/kg]. More...
 
double m_pressure = 0.0
 Current pressure in the reactor [Pa]. More...
 
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. More...
 
string m_name
 
ReactorNetm_net = nullptr
 The ReactorNet that this reactor is part of. More...
 

Member Function Documentation

◆ type()

string type ( ) const
inlineoverridevirtual

String indicating the reactor model implemented.

Usually corresponds to the name of the derived class.

Reimplemented from ReactorBase.

Definition at line 48 of file Reactor.h.

◆ isOde()

virtual bool isOde ( ) const
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.

Definition at line 55 of file Reactor.h.

◆ timeIsIndependent()

virtual bool timeIsIndependent ( ) const
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.

Definition at line 61 of file Reactor.h.

◆ insert()

void insert ( G &  contents)
inline

Insert something into the reactor.

The 'something' must belong to a class that is a subclass of both ThermoPhase and Kinetics.

Definition at line 70 of file Reactor.h.

◆ setKineticsMgr()

void setKineticsMgr ( Kinetics kin)
overridevirtual

Specify chemical kinetics governing the reactor.

Reimplemented from ReactorBase.

Definition at line 39 of file Reactor.cpp.

◆ setChemistry()

void setChemistry ( bool  cflag = true)
inlineoverridevirtual

Enable or disable changes in reactor composition due to chemical reactions.

Reimplemented from ReactorBase.

Definition at line 79 of file Reactor.h.

◆ chemistryEnabled()

bool chemistryEnabled ( ) const
inline

Returns true if changes in the reactor composition due to chemical reactions are enabled.

Definition at line 85 of file Reactor.h.

◆ setEnergy()

void setEnergy ( int  eflag = 1)
inlineoverridevirtual

Set the energy equation on or off.

Reimplemented from ReactorBase.

Definition at line 89 of file Reactor.h.

◆ energyEnabled()

bool energyEnabled ( ) const
inline

Returns true if solution of the energy equation is enabled.

Definition at line 98 of file Reactor.h.

◆ neq()

size_t neq ( )
inline

Number of equations (state variables) for this reactor.

Definition at line 103 of file Reactor.h.

◆ getState()

void getState ( double *  y)
virtual

Get the the current state of the reactor.

Parameters
[out]ystate vector representing the initial state of the reactor

Reimplemented in MoleReactor, IdealGasReactor, IdealGasMoleReactor, IdealGasConstPressureReactor, IdealGasConstPressureMoleReactor, FlowReactor, ConstPressureReactor, and ConstPressureMoleReactor.

Definition at line 49 of file Reactor.cpp.

◆ getStateDae()

virtual void getStateDae ( double *  y,
double *  ydot 
)
inlinevirtual

Get the current state and derivative vector of the reactor for a DAE solver.

Parameters
[out]ystate vector representing the initial state of the reactor
[out]ydotstate vector representing the initial derivatives of the reactor

Reimplemented in FlowReactor.

Definition at line 122 of file Reactor.h.

◆ initialize()

void initialize ( double  t0 = 0.0)
overridevirtual

Initialize the reactor.

Called automatically by ReactorNet::initialize.

Reimplemented from ReactorBase.

Definition at line 84 of file Reactor.cpp.

◆ eval()

void eval ( double  t,
double *  LHS,
double *  RHS 
)
virtual

Evaluate the reactor governing equations.

Called by ReactorNet::eval.

Parameters
[in]ttime.
[out]LHSpointer to start of vector of left-hand side coefficients for governing equations, length m_nv, default values 1
[out]RHSpointer to start of vector of right-hand side coefficients for governing equations, length m_nv, default values 0

Reimplemented in MoleReactor, IdealGasReactor, IdealGasMoleReactor, IdealGasConstPressureReactor, IdealGasConstPressureMoleReactor, FlowReactor, ConstPressureReactor, and ConstPressureMoleReactor.

Definition at line 211 of file Reactor.cpp.

◆ evalDae()

virtual void evalDae ( double  t,
double *  y,
double *  ydot,
double *  residual 
)
inlinevirtual

Evaluate the reactor governing equations.

Called by ReactorNet::eval.

Parameters
[in]ttime.
[in]ysolution vector, length neq()
[in]ydotrate of change of solution vector, length neq()
[out]residualresiduals vector, length neq()

Reimplemented in FlowReactor.

Definition at line 143 of file Reactor.h.

◆ getConstraints()

virtual void getConstraints ( double *  constraints)
inlinevirtual

Given a vector of length neq(), mark which variables should be considered algebraic constraints.

Reimplemented in FlowReactor.

Definition at line 149 of file Reactor.h.

◆ syncState()

void syncState ( )
overridevirtual

Set the state of the reactor to correspond to the state of the associated ThermoPhase object.

This is the inverse of restoreState(). Calling this will trigger integrator reinitialization.

Reimplemented from ReactorBase.

Definition at line 121 of file Reactor.cpp.

◆ updateState()

void updateState ( double *  y)
virtual

Set the state of the reactor to correspond to the state vector y.

Reimplemented in MoleReactor, IdealGasReactor, IdealGasMoleReactor, IdealGasConstPressureReactor, IdealGasConstPressureMoleReactor, FlowReactor, ConstPressureReactor, and ConstPressureMoleReactor.

Definition at line 127 of file Reactor.cpp.

◆ nSensParams()

size_t nSensParams ( ) const
virtual

Number of sensitivity parameters associated with this reactor (including walls)

Definition at line 112 of file Reactor.cpp.

◆ addSensitivityReaction()

void addSensitivityReaction ( size_t  rxn)
virtual

Add a sensitivity parameter associated with the reaction number rxn (in the homogeneous phase).

Definition at line 398 of file Reactor.cpp.

◆ addSensitivitySpeciesEnthalpy()

void addSensitivitySpeciesEnthalpy ( size_t  k)
virtual

Add a sensitivity parameter associated with the enthalpy formation of species k (in the homogeneous phase)

Definition at line 411 of file Reactor.cpp.

◆ componentIndex()

size_t componentIndex ( const string &  nm) const
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 MoleReactor, IdealGasReactor, IdealGasMoleReactor, IdealGasConstPressureReactor, IdealGasConstPressureMoleReactor, FlowReactor, ConstPressureReactor, and ConstPressureMoleReactor.

Definition at line 448 of file Reactor.cpp.

◆ componentName()

string componentName ( size_t  k)
virtual

Return the name of the solution component with index i.

See also
componentIndex()

Reimplemented in MoleReactor, IdealGasReactor, IdealGasMoleReactor, IdealGasConstPressureReactor, IdealGasConstPressureMoleReactor, FlowReactor, ConstPressureReactor, and ConstPressureMoleReactor.

Definition at line 464 of file Reactor.cpp.

◆ setAdvanceLimits()

void setAdvanceLimits ( const double *  limits)

Set absolute step size limits during advance.

Parameters
limitsarray of step size limits with length neq

Definition at line 533 of file Reactor.cpp.

◆ hasAdvanceLimits()

bool hasAdvanceLimits ( ) const
inline

Check whether Reactor object uses advance limits.

Returns
True if at least one limit is set, False otherwise

Definition at line 186 of file Reactor.h.

◆ getAdvanceLimits()

bool getAdvanceLimits ( double *  limits) const

Retrieve absolute step size limits during advance.

Parameters
[out]limitsarray of step size limits with length neq
Returns
True if at least one limit is set, False otherwise

Definition at line 548 of file Reactor.cpp.

◆ setAdvanceLimit()

void setAdvanceLimit ( const string &  nm,
const double  limit 
)

Set individual step size limit for component name nm

Parameters
nmcomponent name
limitvalue for step size limit

Definition at line 559 of file Reactor.cpp.

◆ finiteDifferenceJacobian()

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.

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Definition at line 316 of file Reactor.cpp.

◆ setDerivativeSettings()

void setDerivativeSettings ( AnyMap settings)
virtual

Use this to set the kinetics objects derivative settings.

Definition at line 30 of file Reactor.cpp.

◆ applySensitivity()

void applySensitivity ( double *  params)
virtual

Set reaction rate multipliers based on the sensitivity variables in params.

Definition at line 490 of file Reactor.cpp.

◆ resetSensitivity()

void resetSensitivity ( double *  params)
virtual

Reset the reaction rate multipliers.

Definition at line 512 of file Reactor.cpp.

◆ preconditionerSupported()

virtual bool preconditionerSupported ( ) const
inlinevirtual

Return a false if preconditioning is not supported or true otherwise.

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.
Since
New in Cantera 3.0

Reimplemented in IdealGasMoleReactor, and IdealGasConstPressureMoleReactor.

Definition at line 238 of file Reactor.h.

◆ speciesIndex()

size_t speciesIndex ( const string &  nm) const
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 426 of file Reactor.cpp.

◆ evalWalls()

void evalWalls ( double  t)
protectedvirtual

Evaluate terms related to Walls.

Calculates m_vdot and m_Qdot based on wall movement and heat transfer.

Parameters
tthe current time

Definition at line 275 of file Reactor.cpp.

◆ evalSurfaces()

void evalSurfaces ( double *  LHS,
double *  RHS,
double *  sdot 
)
protectedvirtual

Evaluate terms related to surface reactions.

Parameters
[out]LHSMultiplicative factor on the left hand side of ODE for surface species coverages
[out]RHSRight hand side of ODE for surface species coverages
[out]sdotarray of production rates of bulk phase species on surfaces [kmol/s]

Reimplemented in MoleReactor.

Definition at line 287 of file Reactor.cpp.

◆ updateSurfaceState()

void updateSurfaceState ( double *  y)
protectedvirtual

Update the state of SurfPhase objects attached to this reactor.

Reimplemented in MoleReactor, and FlowReactor.

Definition at line 175 of file Reactor.cpp.

◆ updateConnected()

void updateConnected ( bool  updatePressure)
protectedvirtual

Update the state information needed by connected reactors, flow devices, and reactor walls.

Called from updateState().

Parameters
updatePressureIndicates 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 184 of file Reactor.cpp.

◆ getSurfaceInitialConditions()

void getSurfaceInitialConditions ( double *  y)
protectedvirtual

Get initial conditions for SurfPhase objects attached to this reactor.

Reimplemented in MoleReactor.

Definition at line 75 of file Reactor.cpp.

Member Data Documentation

◆ m_kin

Kinetics* m_kin = nullptr
protected

Pointer to the homogeneous Kinetics object that handles the reactions.

Definition at line 277 of file Reactor.h.

◆ m_vdot

double m_vdot = 0.0
protected

net rate of volume change from moving walls [m^3/s]

Definition at line 279 of file Reactor.h.

◆ m_Qdot

double m_Qdot = 0.0
protected

net heat transfer into the reactor, through walls [W]

Definition at line 281 of file Reactor.h.

◆ m_mass

double m_mass = 0.0
protected

total mass

Definition at line 283 of file Reactor.h.

◆ m_sdot

vector<double> m_sdot
protected

Production rates of gas phase species on surfaces [kmol/s].

Definition at line 287 of file Reactor.h.

◆ m_wdot

vector<double> m_wdot
protected

Species net molar production rates.

Definition at line 289 of file Reactor.h.

◆ m_uk

vector<double> m_uk
protected

Species molar internal energies.

Definition at line 290 of file Reactor.h.

◆ m_advancelimits

vector<double> m_advancelimits
protected

!< Number of variables associated with reactor surfaces

Advance step limit

Definition at line 296 of file Reactor.h.

◆ m_jac_trips

vector<Eigen::Triplet<double> > m_jac_trips
protected

Vector of triplets representing the jacobian.

Definition at line 302 of file Reactor.h.


The documentation for this class was generated from the following files: