Cantera  3.1.0a1
ReactorNet Class Reference

A class representing a network of connected reactors. More...

#include <ReactorNet.h>

Inheritance diagram for ReactorNet:
[legend]

Detailed Description

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
 
ReactorNetoperator= (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< Integratorm_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< PreconditionerBasem_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...
 
Reactorreactor (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...
 
Integratorintegrator ()
 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...
 

Member Function Documentation

◆ setLinearSolverType()

void setLinearSolverType ( const string &  linSolverType = "DENSE")

Set the type of linear solver used in the integration.

Parameters
linSolverTypetype of linear solver. Default type: "DENSE" Other options include: "DIAG", "DENSE", "GMRES", "BAND"

Definition at line 151 of file ReactorNet.cpp.

◆ setPreconditioner()

void setPreconditioner ( shared_ptr< PreconditionerBase preconditioner)

Set preconditioner used by the linear solver.

Parameters
preconditionerpreconditioner object used for the linear solver

Definition at line 157 of file ReactorNet.cpp.

◆ setInitialTime()

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.

◆ getInitialTime()

double getInitialTime ( ) const
inline

Get the initial value of the independent variable (typically time).

Since
New in Cantera 3.0.

Definition at line 58 of file ReactorNet.h.

◆ maxTimeStep()

double maxTimeStep ( ) const
inline

Get the maximum integrator step.

Definition at line 63 of file ReactorNet.h.

◆ setMaxTimeStep()

void setMaxTimeStep ( double  maxstep)

Set the maximum integrator step.

Definition at line 35 of file ReactorNet.cpp.

◆ setMaxErrTestFails()

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.

◆ setTolerances()

void setTolerances ( double  rtol,
double  atol 
)

Set the relative and absolute tolerances for the integrator.

Definition at line 46 of file ReactorNet.cpp.

◆ setSensitivityTolerances()

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.

◆ time()

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.

◆ distance()

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.

◆ rtol()

double rtol ( )
inline

Relative tolerance.

Definition at line 90 of file ReactorNet.h.

◆ atol()

double atol ( )
inline

Absolute integration tolerance.

Definition at line 95 of file ReactorNet.h.

◆ rtolSensitivity()

double rtolSensitivity ( ) const
inline

Relative sensitivity tolerance.

Definition at line 100 of file ReactorNet.h.

◆ atolSensitivity()

double atolSensitivity ( ) const
inline

Absolute sensitivity tolerance.

Definition at line 105 of file ReactorNet.h.

◆ linearSolverType()

string linearSolverType ( ) const

Problem type of integrator.

Definition at line 510 of file ReactorNet.cpp.

◆ maxSteps()

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.

◆ advance() [1/2]

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.

Parameters
tTime/distance to advance to (s or m).

Definition at line 173 of file ReactorNet.cpp.

◆ advance() [2/2]

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.

Parameters
tTime/distance to advance to (s or m).
applylimitLimit advance step (boolean).

Definition at line 185 of file ReactorNet.cpp.

◆ step()

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.

◆ addReactor()

void addReactor ( Reactor r)

Add the reactor r to this reactor network.

Definition at line 284 of file ReactorNet.cpp.

◆ reactor()

Reactor& reactor ( int  n)
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.

◆ verbose()

bool verbose ( ) const
inline

Returns true if verbose logging output is enabled.

Definition at line 149 of file ReactorNet.h.

◆ setVerbose()

void setVerbose ( bool  v = true)
inline

Enable or disable verbose logging while setting up and integrating the reactor network.

Definition at line 155 of file ReactorNet.h.

◆ integrator()

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.

◆ updateState()

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.

◆ sensitivity() [1/2]

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.

◆ sensitivity() [2/2]

double sensitivity ( const string &  component,
size_t  p,
int  reactor = 0 
)
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.

◆ evalJacobian()

void evalJacobian ( double  t,
double *  y,
double *  ydot,
double *  p,
Array2D j 
)

Evaluate the Jacobian matrix for the reactor network.

Parameters
[in]tTime/distance at which to evaluate the Jacobian
[in]yGlobal state vector at t
[out]ydotDerivative of the state vector evaluated at t, with respect to t.
[in]psensitivity parameter vector (unused?)
[out]jJacobian matrix, size neq() by neq().

Definition at line 376 of file ReactorNet.cpp.

◆ neq()

size_t neq ( ) const
inlineoverridevirtual

Number of equations.

Implements FuncEval.

Definition at line 208 of file ReactorNet.h.

◆ eval()

void eval ( double  t,
double *  y,
double *  ydot,
double *  p 
)
overridevirtual

Evaluate the right-hand-side ODE function.

Called by the integrator.

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

Reimplemented from FuncEval.

Definition at line 318 of file ReactorNet.cpp.

◆ evalDae()

void evalDae ( double  t,
double *  y,
double *  ydot,
double *  p,
double *  residual 
)
overridevirtual

eval coupling for IDA / DAEs

Reimplemented from FuncEval.

Definition at line 341 of file ReactorNet.cpp.

◆ getState()

void getState ( double *  y)
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.

◆ getStateDae()

void getStateDae ( double *  y,
double *  ydot 
)
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.

◆ getDerivative()

void getDerivative ( int  k,
double *  dky 
)
virtual

Return k-th derivative at the current state of the system.

Definition at line 406 of file ReactorNet.cpp.

◆ getConstraints()

void getConstraints ( double *  constraints)
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.

◆ nparams()

size_t nparams ( ) const
inlineoverridevirtual

Number of sensitivity parameters.

Reimplemented from FuncEval.

Definition at line 230 of file ReactorNet.h.

◆ globalComponentIndex()

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.

◆ componentName()

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.

◆ registerSensitivityParameter()

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.

Parameters
nameA name describing the parameter, for example the reaction string
valueThe nominal value of the parameter
scaleA scaling factor to be applied to the sensitivity coefficient
Returns
the index of this parameter in the vector of sensitivity parameters (global across all reactors)

Definition at line 479 of file ReactorNet.cpp.

◆ sensitivityParameterName()

const string& sensitivityParameterName ( size_t  p) const
inline

The name of the p-th sensitivity parameter added to this ReactorNet.

Definition at line 256 of file ReactorNet.h.

◆ initialize()

void initialize ( )

Initialize the reactor network.

Called automatically the first time advance or step is called.

Definition at line 86 of file ReactorNet.cpp.

◆ reinitialize()

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.

◆ setNeedsReinit()

void setNeedsReinit ( )
inline

Called to trigger integrator reinitialization before further integration.

Definition at line 272 of file ReactorNet.h.

◆ setMaxSteps()

void setMaxSteps ( int  nmax)
virtual

Set the maximum number of internal integration steps the integrator will take before reaching the next output point.

Parameters
nmaxThe maximum number of steps, setting this value to zero disables this option.

Definition at line 163 of file ReactorNet.cpp.

◆ setAdvanceLimits()

void setAdvanceLimits ( const double *  limits)

Set absolute step size limits during advance.

Definition at line 417 of file ReactorNet.cpp.

◆ hasAdvanceLimits()

bool hasAdvanceLimits ( ) const

Check whether ReactorNet object uses advance limits.

Definition at line 427 of file ReactorNet.cpp.

◆ getAdvanceLimits()

bool getAdvanceLimits ( double *  limits) const

Retrieve absolute step size limits during advance.

Definition at line 436 of file ReactorNet.cpp.

◆ preconditionerSetup()

void preconditionerSetup ( double  t,
double *  y,
double  gamma 
)
overridevirtual

Evaluate the setup processes for the Jacobian preconditioner.

Parameters
[in]ttime.
[in]ysolution vector, length neq()
gammathe gamma in M=I-gamma*J
Warning
This function is an experimental part of the Cantera API and may be changed or removed without notice.

Reimplemented from FuncEval.

Definition at line 528 of file ReactorNet.cpp.

◆ preconditionerSolve()

void preconditionerSolve ( double *  rhs,
double *  output 
)
overridevirtual

Evaluate the linear system Ax=b where A is the preconditioner.

Parameters
[in]rhsright hand side vector used in linear system
[out]outputoutput vector for solution
Warning
This function is an experimental part of the Cantera API and may be changed or removed without notice.

Reimplemented from FuncEval.

Definition at line 519 of file ReactorNet.cpp.

◆ solverStats()

AnyMap solverStats ( ) const

Get solver stats from integrator.

Definition at line 501 of file ReactorNet.cpp.

◆ setDerivativeSettings()

void setDerivativeSettings ( AnyMap settings)
virtual

Set derivative settings of all reactors.

Parameters
settingsthe settings map propagated to all reactors and kinetics objects

Definition at line 493 of file ReactorNet.cpp.

◆ checkPreconditionerSupported()

void checkPreconditionerSupported ( ) const
protectedvirtual

Check that preconditioning is supported by all reactors in the network.

Definition at line 571 of file ReactorNet.cpp.

◆ updatePreconditioner()

void updatePreconditioner ( double  gamma)
overrideprotectedvirtual

Update the preconditioner based on already computed jacobian values.

Reimplemented from FuncEval.

Definition at line 560 of file ReactorNet.cpp.

◆ getEstimate()

void getEstimate ( double  time,
int  k,
double *  yest 
)
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.

◆ lastOrder()

int lastOrder ( ) const
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.

Member Data Documentation

◆ m_time

double m_time = 0.0
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.

◆ m_initial_time

double m_initial_time = 0.0
protected

The initial value of the independent variable in the system.

Definition at line 326 of file ReactorNet.h.

◆ m_integrator_init

bool m_integrator_init = false
protected

True if integrator initialization is current.

Definition at line 329 of file ReactorNet.h.

◆ m_start

vector<size_t> m_start
protected

m_start[n] is the starting point in the state vector for reactor n

Definition at line 333 of file ReactorNet.h.

◆ m_maxstep

double m_maxstep = 0.0
protected

Maximum integrator internal timestep. Default of 0.0 means infinity.

Definition at line 344 of file ReactorNet.h.

◆ m_timeIsIndependent

bool m_timeIsIndependent = true
protected

Indicates whether time or space is the independent variable.

Definition at line 349 of file ReactorNet.h.

◆ m_paramNames

vector<string> m_paramNames
protected

Names corresponding to each sensitivity parameter.

Definition at line 352 of file ReactorNet.h.

◆ m_LHS

vector<double> m_LHS
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.


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