Cantera
3.1.0a1
|
A class representing a network of connected reactors. More...
#include <ReactorNet.h>
A class representing a network of connected reactors.
This class is used to integrate the governing equations for a network of reactors that are time dependent (Reactor, ConstPressureReactor) connected by various means, for example Wall, MassFlowController, Valve, or PressureController; or reactors dependent on a single spatial variable (FlowReactor).
Definition at line 29 of file ReactorNet.h.
Public Member Functions | |
ReactorNet (const ReactorNet &)=delete | |
ReactorNet & | operator= (const ReactorNet &)=delete |
Public Member Functions inherited from FuncEval | |
int | evalNoThrow (double t, double *y, double *ydot) |
Evaluate the right-hand side using return code to indicate status. More... | |
int | evalDaeNoThrow (double t, double *y, double *ydot, double *residual) |
Evaluate the right-hand side using return code to indicate status. More... | |
int | preconditioner_setup_nothrow (double t, double *y, double gamma) |
Preconditioner setup that doesn't throw an error but returns a CVODES flag. More... | |
int | preconditioner_solve_nothrow (double *rhs, double *output) |
Preconditioner solve that doesn't throw an error but returns a CVODES flag. More... | |
void | suppressErrors (bool suppress) |
Enable or disable suppression of errors when calling eval() More... | |
bool | suppressErrors () const |
Get current state of error suppression. More... | |
string | getErrors () const |
Return a string containing the text of any suppressed errors. More... | |
void | clearErrors () |
Clear any previously-stored suppressed errors. More... | |
Methods to set up a simulation | |
vector< Reactor * > | m_reactors |
unique_ptr< Integrator > | m_integ |
double | m_time = 0.0 |
The independent variable in the system. More... | |
double | m_initial_time = 0.0 |
The initial value of the independent variable in the system. More... | |
bool | m_init = false |
bool | m_integrator_init = false |
True if integrator initialization is current. More... | |
size_t | m_nv = 0 |
vector< size_t > | m_start |
m_start[n] is the starting point in the state vector for reactor n More... | |
vector< double > | m_atol |
double | m_rtol = 1.0e-9 |
double | m_rtolsens = 1.0e-4 |
double | m_atols = 1.0e-15 |
double | m_atolsens = 1.0e-6 |
shared_ptr< PreconditionerBase > | m_precon |
string | m_linearSolverType |
double | m_maxstep = 0.0 |
Maximum integrator internal timestep. Default of 0.0 means infinity. More... | |
bool | m_verbose = false |
bool | m_timeIsIndependent = true |
Indicates whether time or space is the independent variable. More... | |
vector< string > | m_paramNames |
Names corresponding to each sensitivity parameter. More... | |
vector< double > | m_ydot |
vector< double > | m_yest |
vector< double > | m_advancelimits |
vector< double > | m_LHS |
m_LHS is a vector representing the coefficients on the "left hand side" of each governing equation More... | |
vector< double > | m_RHS |
void | setLinearSolverType (const string &linSolverType="DENSE") |
Set the type of linear solver used in the integration. More... | |
void | setPreconditioner (shared_ptr< PreconditionerBase > preconditioner) |
Set preconditioner used by the linear solver. More... | |
void | setInitialTime (double time) |
Set the initial value of the independent variable (typically time). More... | |
double | getInitialTime () const |
Get the initial value of the independent variable (typically time). More... | |
double | maxTimeStep () const |
Get the maximum integrator step. More... | |
void | setMaxTimeStep (double maxstep) |
Set the maximum integrator step. More... | |
void | setMaxErrTestFails (int nmax) |
Set the maximum number of error test failures permitted by the CVODES integrator in a single step. More... | |
void | setTolerances (double rtol, double atol) |
Set the relative and absolute tolerances for the integrator. More... | |
void | setSensitivityTolerances (double rtol, double atol) |
Set the relative and absolute tolerances for integrating the sensitivity equations. More... | |
double | time () |
Current value of the simulation time [s], for reactor networks that are solved in the time domain. More... | |
double | distance () |
Current position [m] along the length of the reactor network, for reactors that are solved as a function of space. More... | |
double | rtol () |
Relative tolerance. More... | |
double | atol () |
Absolute integration tolerance. More... | |
double | rtolSensitivity () const |
Relative sensitivity tolerance. More... | |
double | atolSensitivity () const |
Absolute sensitivity tolerance. More... | |
string | linearSolverType () const |
Problem type of integrator. More... | |
int | maxSteps () |
Returns the maximum number of internal integration steps the integrator will take before reaching the next output point. More... | |
void | advance (double t) |
Advance the state of all reactors in the independent variable (time or space). More... | |
double | advance (double t, bool applylimit) |
Advance the state of all reactors in the independent variable (time or space). More... | |
double | step () |
Advance the state of all reactors with respect to the independent variable (time or space). More... | |
void | addReactor (Reactor &r) |
Add the reactor r to this reactor network. More... | |
Reactor & | reactor (int n) |
Return a reference to the n-th reactor in this network. More... | |
bool | verbose () const |
Returns true if verbose logging output is enabled. More... | |
void | setVerbose (bool v=true) |
Enable or disable verbose logging while setting up and integrating the reactor network. More... | |
Integrator & | integrator () |
Return a reference to the integrator. More... | |
void | updateState (double *y) |
Update the state of all the reactors in the network to correspond to the values in the solution vector y. More... | |
double | sensitivity (size_t k, size_t p) |
Return the sensitivity of the k-th solution component with respect to the p-th sensitivity parameter. More... | |
double | sensitivity (const string &component, size_t p, int reactor=0) |
Return the sensitivity of the component named component with respect to the p-th sensitivity parameter. More... | |
void | evalJacobian (double t, double *y, double *ydot, double *p, Array2D *j) |
Evaluate the Jacobian matrix for the reactor network. More... | |
size_t | neq () const override |
Number of equations. More... | |
size_t | nReactors () const |
void | eval (double t, double *y, double *ydot, double *p) override |
Evaluate the right-hand-side ODE function. More... | |
void | evalDae (double t, double *y, double *ydot, double *p, double *residual) override |
eval coupling for IDA / DAEs More... | |
void | getState (double *y) override |
Fill in the vector y with the current state of the system. More... | |
void | getStateDae (double *y, double *ydot) override |
Fill in the vectors y and ydot with the current state of the system. More... | |
virtual void | getDerivative (int k, double *dky) |
Return k-th derivative at the current state of the system. More... | |
void | getConstraints (double *constraints) override |
Given a vector of length neq(), mark which variables should be considered algebraic constraints. More... | |
size_t | nparams () const override |
Number of sensitivity parameters. More... | |
size_t | globalComponentIndex (const string &component, size_t reactor=0) |
Return the index corresponding to the component named component in the reactor with index reactor in the global state vector for the reactor network. More... | |
string | componentName (size_t i) const |
Return the name of the i-th component of the global state vector. More... | |
size_t | registerSensitivityParameter (const string &name, double value, double scale) |
Used by Reactor and Wall objects to register the addition of sensitivity parameters so that the ReactorNet can keep track of the order in which sensitivity parameters are added. More... | |
const string & | sensitivityParameterName (size_t p) const |
The name of the p-th sensitivity parameter added to this ReactorNet. More... | |
void | initialize () |
Initialize the reactor network. More... | |
void | reinitialize () |
Reinitialize the integrator. More... | |
void | setNeedsReinit () |
Called to trigger integrator reinitialization before further integration. More... | |
virtual void | setMaxSteps (int nmax) |
Set the maximum number of internal integration steps the integrator will take before reaching the next output point. More... | |
void | setAdvanceLimits (const double *limits) |
Set absolute step size limits during advance. More... | |
bool | hasAdvanceLimits () const |
Check whether ReactorNet object uses advance limits. More... | |
bool | getAdvanceLimits (double *limits) const |
Retrieve absolute step size limits during advance. More... | |
void | preconditionerSetup (double t, double *y, double gamma) override |
Evaluate the setup processes for the Jacobian preconditioner. More... | |
void | preconditionerSolve (double *rhs, double *output) override |
Evaluate the linear system Ax=b where A is the preconditioner. More... | |
AnyMap | solverStats () const |
Get solver stats from integrator. More... | |
virtual void | setDerivativeSettings (AnyMap &settings) |
Set derivative settings of all reactors. More... | |
virtual void | checkPreconditionerSupported () const |
Check that preconditioning is supported by all reactors in the network. More... | |
void | updatePreconditioner (double gamma) override |
Update the preconditioner based on already computed jacobian values. More... | |
virtual void | getEstimate (double time, int k, double *yest) |
Estimate a future state based on current derivatives. More... | |
virtual int | lastOrder () const |
Returns the order used for last solution step of the ODE integrator The function is intended for internal use by ReactorNet::advance and deliberately not exposed in external interfaces. More... | |
Additional Inherited Members | |
Public Attributes inherited from FuncEval | |
vector< double > | m_sens_params |
Values for the problem parameters for which sensitivities are computed This is the array which is perturbed and passed back as the fourth argument to eval(). More... | |
vector< double > | m_paramScales |
Scaling factors for each sensitivity parameter. More... | |
Protected Attributes inherited from FuncEval | |
bool | m_suppress_errors = false |
vector< string > | m_errors |
Errors occurring during function evaluations. More... | |
void setLinearSolverType | ( | const string & | linSolverType = "DENSE" | ) |
Set the type of linear solver used in the integration.
linSolverType | type of linear solver. Default type: "DENSE" Other options include: "DIAG", "DENSE", "GMRES", "BAND" |
Definition at line 151 of file ReactorNet.cpp.
void setPreconditioner | ( | shared_ptr< PreconditionerBase > | preconditioner | ) |
Set preconditioner used by the linear solver.
preconditioner | preconditioner object used for the linear solver |
Definition at line 157 of file ReactorNet.cpp.
void setInitialTime | ( | double | time | ) |
Set the initial value of the independent variable (typically time).
Default = 0.0 s. Restarts integration from this value using the current mixture state as the initial condition.
Definition at line 28 of file ReactorNet.cpp.
|
inline |
Get the initial value of the independent variable (typically time).
Definition at line 58 of file ReactorNet.h.
|
inline |
Get the maximum integrator step.
Definition at line 63 of file ReactorNet.h.
void setMaxTimeStep | ( | double | maxstep | ) |
Set the maximum integrator step.
Definition at line 35 of file ReactorNet.cpp.
void setMaxErrTestFails | ( | int | nmax | ) |
Set the maximum number of error test failures permitted by the CVODES integrator in a single step.
Definition at line 41 of file ReactorNet.cpp.
void setTolerances | ( | double | rtol, |
double | atol | ||
) |
Set the relative and absolute tolerances for the integrator.
Definition at line 46 of file ReactorNet.cpp.
void setSensitivityTolerances | ( | double | rtol, |
double | atol | ||
) |
Set the relative and absolute tolerances for integrating the sensitivity equations.
Definition at line 57 of file ReactorNet.cpp.
double time | ( | ) |
Current value of the simulation time [s], for reactor networks that are solved in the time domain.
Definition at line 68 of file ReactorNet.cpp.
double distance | ( | ) |
Current position [m] along the length of the reactor network, for reactors that are solved as a function of space.
Definition at line 77 of file ReactorNet.cpp.
|
inline |
Relative tolerance.
Definition at line 90 of file ReactorNet.h.
|
inline |
Absolute integration tolerance.
Definition at line 95 of file ReactorNet.h.
|
inline |
Relative sensitivity tolerance.
Definition at line 100 of file ReactorNet.h.
|
inline |
Absolute sensitivity tolerance.
Definition at line 105 of file ReactorNet.h.
string linearSolverType | ( | ) | const |
Problem type of integrator.
Definition at line 510 of file ReactorNet.cpp.
int maxSteps | ( | ) |
Returns the maximum number of internal integration steps the integrator will take before reaching the next output point.
Definition at line 168 of file ReactorNet.cpp.
void advance | ( | double | t | ) |
Advance the state of all reactors in the independent variable (time or space).
Take as many internal steps as necessary to reach t.
t | Time/distance to advance to (s or m). |
Definition at line 173 of file ReactorNet.cpp.
double advance | ( | double | t, |
bool | applylimit | ||
) |
Advance the state of all reactors in the independent variable (time or space).
Take as many internal steps as necessary towards t. If applylimit is true, the advance step will be automatically reduced if needed to stay within limits (set by setAdvanceLimit). Returns the time/distance at the end of integration.
t | Time/distance to advance to (s or m). |
applylimit | Limit advance step (boolean). |
Definition at line 185 of file ReactorNet.cpp.
double step | ( | ) |
Advance the state of all reactors with respect to the independent variable (time or space).
Returns the new value of the independent variable [s or m].
Definition at line 240 of file ReactorNet.cpp.
void addReactor | ( | Reactor & | r | ) |
Add the reactor r to this reactor network.
Definition at line 284 of file ReactorNet.cpp.
|
inline |
Return a reference to the n-th reactor in this network.
The reactor indices are determined by the order in which the reactors were added to the reactor network.
Definition at line 144 of file ReactorNet.h.
|
inline |
Returns true
if verbose logging output is enabled.
Definition at line 149 of file ReactorNet.h.
|
inline |
Enable or disable verbose logging while setting up and integrating the reactor network.
Definition at line 155 of file ReactorNet.h.
Integrator & integrator | ( | ) |
Return a reference to the integrator.
Only valid after adding at least one reactor to the network.
Definition at line 310 of file ReactorNet.cpp.
void updateState | ( | double * | y | ) |
Update the state of all the reactors in the network to correspond to the values in the solution vector y.
Definition at line 398 of file ReactorNet.cpp.
double sensitivity | ( | size_t | k, |
size_t | p | ||
) |
Return the sensitivity of the k-th solution component with respect to the p-th sensitivity parameter.
The sensitivity coefficient \( S_{ki} \) of solution variable \( y_k \) with respect to sensitivity parameter \( p_i \) is defined as:
\[ S_{ki} = \frac{1}{y_k} \frac{\partial y_k}{\partial p_i} \]
For reaction sensitivities, the parameter is a multiplier on the forward rate constant (and implicitly on the reverse rate constant for reversible reactions) which has a nominal value of 1.0, and the sensitivity is nondimensional.
For species enthalpy sensitivities, the parameter is a perturbation to the molar enthalpy of formation, such that the dimensions of the sensitivity are kmol/J.
Definition at line 360 of file ReactorNet.cpp.
|
inline |
Return the sensitivity of the component named component with respect to the p-th sensitivity parameter.
The sensitivity coefficient \( S_{ki} \) of solution variable \( y_k \) with respect to sensitivity parameter \( p_i \) is defined as:
\[ S_{ki} = \frac{1}{y_k} \frac{\partial y_k}{\partial p_i} \]
For reaction sensitivities, the parameter is a multiplier on the forward rate constant (and implicitly on the reverse rate constant for reversible reactions) which has a nominal value of 1.0, and the sensitivity is nondimensional.
For species enthalpy sensitivities, the parameter is a perturbation to the molar enthalpy of formation, such that the dimensions of the sensitivity are kmol/J.
Definition at line 190 of file ReactorNet.h.
void evalJacobian | ( | double | t, |
double * | y, | ||
double * | ydot, | ||
double * | p, | ||
Array2D * | j | ||
) |
Evaluate the Jacobian matrix for the reactor network.
[in] | t | Time/distance at which to evaluate the Jacobian |
[in] | y | Global state vector at t |
[out] | ydot | Derivative of the state vector evaluated at t, with respect to t. |
[in] | p | sensitivity parameter vector (unused?) |
[out] | j | Jacobian matrix, size neq() by neq(). |
Definition at line 376 of file ReactorNet.cpp.
|
inlineoverridevirtual |
|
overridevirtual |
Evaluate the right-hand-side ODE function.
Called by the integrator.
[in] | t | time. |
[in] | y | solution vector, length neq() |
[out] | ydot | rate of change of solution vector, length neq() |
[in] | p | sensitivity parameter vector, length nparams() |
Reimplemented from FuncEval.
Definition at line 318 of file ReactorNet.cpp.
|
overridevirtual |
eval coupling for IDA / DAEs
Reimplemented from FuncEval.
Definition at line 341 of file ReactorNet.cpp.
|
overridevirtual |
Fill in the vector y with the current state of the system.
Used for getting the initial state for ODE systems.
Reimplemented from FuncEval.
Definition at line 445 of file ReactorNet.cpp.
|
overridevirtual |
Fill in the vectors y and ydot with the current state of the system.
Used for getting the initial state for DAE systems.
Reimplemented from FuncEval.
Definition at line 452 of file ReactorNet.cpp.
|
virtual |
Return k-th derivative at the current state of the system.
Definition at line 406 of file ReactorNet.cpp.
|
overridevirtual |
Given a vector of length neq(), mark which variables should be considered algebraic constraints.
Reimplemented from FuncEval.
Definition at line 353 of file ReactorNet.cpp.
|
inlineoverridevirtual |
Number of sensitivity parameters.
Reimplemented from FuncEval.
Definition at line 230 of file ReactorNet.h.
size_t globalComponentIndex | ( | const string & | component, |
size_t | reactor = 0 |
||
) |
Return the index corresponding to the component named component in the reactor with index reactor in the global state vector for the reactor network.
Definition at line 459 of file ReactorNet.cpp.
string componentName | ( | size_t | i | ) | const |
Return the name of the i-th component of the global state vector.
The name returned includes both the name of the reactor and the specific component, for example ‘'reactor1: CH4’`.
Definition at line 467 of file ReactorNet.cpp.
size_t registerSensitivityParameter | ( | const string & | name, |
double | value, | ||
double | scale | ||
) |
Used by Reactor and Wall objects to register the addition of sensitivity parameters so that the ReactorNet can keep track of the order in which sensitivity parameters are added.
name | A name describing the parameter, for example the reaction string |
value | The nominal value of the parameter |
scale | A scaling factor to be applied to the sensitivity coefficient |
Definition at line 479 of file ReactorNet.cpp.
|
inline |
The name of the p-th sensitivity parameter added to this ReactorNet.
Definition at line 256 of file ReactorNet.h.
void initialize | ( | ) |
Initialize the reactor network.
Called automatically the first time advance or step is called.
Definition at line 86 of file ReactorNet.cpp.
void reinitialize | ( | ) |
Reinitialize the integrator.
Used to solve a new problem (different initial conditions) but with the same configuration of the reactor network. Can be called manually, or automatically after calling setInitialTime or modifying a reactor's contents.
Definition at line 137 of file ReactorNet.cpp.
|
inline |
Called to trigger integrator reinitialization before further integration.
Definition at line 272 of file ReactorNet.h.
|
virtual |
Set the maximum number of internal integration steps the integrator will take before reaching the next output point.
nmax | The maximum number of steps, setting this value to zero disables this option. |
Definition at line 163 of file ReactorNet.cpp.
void setAdvanceLimits | ( | const double * | limits | ) |
Set absolute step size limits during advance.
Definition at line 417 of file ReactorNet.cpp.
bool hasAdvanceLimits | ( | ) | const |
Check whether ReactorNet object uses advance limits.
Definition at line 427 of file ReactorNet.cpp.
bool getAdvanceLimits | ( | double * | limits | ) | const |
Retrieve absolute step size limits during advance.
Definition at line 436 of file ReactorNet.cpp.
|
overridevirtual |
Evaluate the setup processes for the Jacobian preconditioner.
[in] | t | time. |
[in] | y | solution vector, length neq() |
gamma | the gamma in M=I-gamma*J |
Reimplemented from FuncEval.
Definition at line 528 of file ReactorNet.cpp.
|
overridevirtual |
Evaluate the linear system Ax=b where A is the preconditioner.
[in] | rhs | right hand side vector used in linear system |
[out] | output | output vector for solution |
Reimplemented from FuncEval.
Definition at line 519 of file ReactorNet.cpp.
AnyMap solverStats | ( | ) | const |
Get solver stats from integrator.
Definition at line 501 of file ReactorNet.cpp.
|
virtual |
Set derivative settings of all reactors.
settings | the settings map propagated to all reactors and kinetics objects |
Definition at line 493 of file ReactorNet.cpp.
|
protectedvirtual |
Check that preconditioning is supported by all reactors in the network.
Definition at line 571 of file ReactorNet.cpp.
|
overrideprotectedvirtual |
Update the preconditioner based on already computed jacobian values.
Reimplemented from FuncEval.
Definition at line 560 of file ReactorNet.cpp.
|
protectedvirtual |
Estimate a future state based on current derivatives.
The function is intended for internal use by ReactorNet::advance and deliberately not exposed in external interfaces.
Definition at line 252 of file ReactorNet.cpp.
|
protectedvirtual |
Returns the order used for last solution step of the ODE integrator The function is intended for internal use by ReactorNet::advance and deliberately not exposed in external interfaces.
Definition at line 275 of file ReactorNet.cpp.
|
protected |
The independent variable in the system.
May be either time or space depending on the type of reactors in the network.
Definition at line 323 of file ReactorNet.h.
|
protected |
The initial value of the independent variable in the system.
Definition at line 326 of file ReactorNet.h.
|
protected |
True if integrator initialization is current.
Definition at line 329 of file ReactorNet.h.
|
protected |
m_start[n] is the starting point in the state vector for reactor n
Definition at line 333 of file ReactorNet.h.
|
protected |
Maximum integrator internal timestep. Default of 0.0 means infinity.
Definition at line 344 of file ReactorNet.h.
|
protected |
Indicates whether time or space is the independent variable.
Definition at line 349 of file ReactorNet.h.
|
protected |
Names corresponding to each sensitivity parameter.
Definition at line 352 of file ReactorNet.h.
|
protected |
m_LHS is a vector representing the coefficients on the "left hand side" of each governing equation
Definition at line 359 of file ReactorNet.h.