Cantera  2.1.2
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
ReactorNet Class Reference

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

#include <ReactorNet.h>

Inheritance diagram for ReactorNet:
[legend]
Collaboration diagram for ReactorNet:
[legend]

Public Member Functions

void addReactor (ReactorBase *r, bool iown=false)
 Add the reactor r to this reactor network. More...
 
ReactorBasereactor (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 (doublereal *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 std::string &species, size_t p, int reactor=0)
 Return the sensitivity of the species named species with respect to the p-th sensitivity parameter. More...
 
void evalJacobian (doublereal t, doublereal *y, doublereal *ydot, doublereal *p, Array2D *j)
 Evaluate the Jacobian matrix for the reactor network. More...
 
virtual size_t neq ()
 Number of equations. More...
 
virtual void eval (doublereal t, doublereal *y, doublereal *ydot, doublereal *p)
 Evaluate the right-hand-side function. More...
 
virtual void getInitialConditions (doublereal t0, size_t leny, doublereal *y)
 Fill the solution vector with the initial conditions at initial time t0. More...
 
virtual size_t nparams ()
 Number of sensitivity parameters. More...
 
size_t globalComponentIndex (const std::string &species, size_t reactor=0)
 Return the index corresponding to the species named species in the reactor with index reactor in the global state vector for the reactor network. More...
 
void registerSensitivityReaction (void *reactor, size_t reactionIndex, const std::string &name, int leftright=0)
 Used by Reactor and Wall objects to register the addition of sensitivity reactions so that the ReactorNet can keep track of the order in which sensitivity parameters are added. More...
 
const std::string & sensitivityParameterName (size_t p)
 The name of the p-th sensitivity parameter added to this ReactorNet. More...
 
void connect (size_t i, size_t j)
 
bool connected (size_t i, size_t j)
 
Methods to set up a simulation.
void setInitialTime (doublereal time)
 Set initial time. More...
 
void setMaxTimeStep (double maxstep)
 Set the maximum time step. More...
 
void setMaxErrTestFails (int nmax)
 Set the maximum number of error test failures permitted by the CVODES integrator in a single time step. More...
 
void setTolerances (doublereal rtol, doublereal atol)
 Set the relative and absolute tolerances for the integrator. More...
 
void setSensitivityTolerances (doublereal rtol, doublereal atol)
 Set the relative and absolute tolerances for integrating the sensitivity equations. More...
 
doublereal time ()
 Current value of the simulation time. More...
 
doublereal rtol ()
 Relative tolerance. More...
 
doublereal atol ()
 Absolute integration tolerance. More...
 
doublereal rtolSensitivity () const
 Relative sensitivity tolerance. More...
 
doublereal atolSensitivity () const
 Absolute sensitivity tolerance. More...
 
void advance (doublereal time)
 Advance the state of all reactors in time. More...
 
double step (doublereal time)
 Advance the state of all reactors in time. More...
 

Protected Member Functions

void initialize ()
 Initialize the reactor network. More...
 

Protected Attributes

std::vector< ReactorBase * > m_r
 
std::vector< Reactor * > m_reactors
 
size_t m_nr
 
size_t m_nreactors
 
Integratorm_integ
 
doublereal m_time
 
bool m_init
 
size_t m_nv
 
std::vector< size_t > m_size
 
vector_fp m_atol
 
doublereal m_rtol
 
doublereal m_rtolsens
 
doublereal m_atols
 
doublereal m_atolsens
 
doublereal m_maxstep
 
int m_maxErrTestFails
 
bool m_verbose
 
size_t m_ntotpar
 
std::vector< size_t > m_nparams
 
std::vector< std::string > m_paramNames
 Names corresponding to each sensitivity parameter. More...
 
std::map< std::pair< void
*, int >, std::map< size_t,
size_t > > 
m_sensOrder
 Structure used to determine the order of sensitivity parameters m_sensOrder[Reactor or Wall, leftright][reaction number] = parameter index. More...
 
std::vector< size_t > m_sensIndex
 Mapping from the order in which sensitivity parameters were added to the ReactorNet to the order in which they occur in the integrator output. More...
 
vector_int m_connect
 
vector_fp m_ydot
 
std::vector< bool > m_iown
 

Detailed Description

A class representing a network of connected reactors.

This class is used to integrate the time-dependent governing equations for a network of reactors (Reactor, ConstPressureReactor) connected by various means, e.g. Wall, MassFlowController, Valve, PressureController.

Definition at line 23 of file ReactorNet.h.

Member Function Documentation

void setInitialTime ( doublereal  time)
inline

Set initial time.

Default = 0.0 s. Restarts integration from this time using the current mixture state as the initial condition.

Definition at line 37 of file ReactorNet.h.

References ReactorNet::time().

void setMaxTimeStep ( double  maxstep)
inline

Set the maximum time step.

Definition at line 43 of file ReactorNet.h.

void setMaxErrTestFails ( int  nmax)
inline

Set the maximum number of error test failures permitted by the CVODES integrator in a single time step.

Definition at line 50 of file ReactorNet.h.

void setTolerances ( doublereal  rtol,
doublereal  atol 
)
inline

Set the relative and absolute tolerances for the integrator.

Definition at line 56 of file ReactorNet.h.

References ReactorNet::atol(), and ReactorNet::rtol().

void setSensitivityTolerances ( doublereal  rtol,
doublereal  atol 
)
inline

Set the relative and absolute tolerances for integrating the sensitivity equations.

Definition at line 68 of file ReactorNet.h.

References ReactorNet::atol(), and ReactorNet::rtol().

doublereal time ( )
inline

Current value of the simulation time.

Definition at line 79 of file ReactorNet.h.

Referenced by ReactorNet::advance(), and ReactorNet::setInitialTime().

doublereal rtol ( )
inline

Relative tolerance.

Definition at line 84 of file ReactorNet.h.

Referenced by ReactorNet::setSensitivityTolerances(), and ReactorNet::setTolerances().

doublereal atol ( )
inline

Absolute integration tolerance.

Definition at line 89 of file ReactorNet.h.

Referenced by ReactorNet::setSensitivityTolerances(), and ReactorNet::setTolerances().

doublereal rtolSensitivity ( ) const
inline

Relative sensitivity tolerance.

Definition at line 94 of file ReactorNet.h.

doublereal atolSensitivity ( ) const
inline

Absolute sensitivity tolerance.

Definition at line 99 of file ReactorNet.h.

void advance ( doublereal  time)

Advance the state of all reactors in time.

Take as many internal timesteps as necessary to reach time.

Parameters
timeTime to advance to (s).

Definition at line 148 of file ReactorNet.cpp.

References ReactorNet::initialize(), Integrator::integrate(), Integrator::solution(), ReactorNet::time(), and ReactorNet::updateState().

double step ( doublereal  time)

Advance the state of all reactors in time.

Take a single timestep toward time.

Definition at line 161 of file ReactorNet.cpp.

References ReactorNet::initialize(), Integrator::solution(), Integrator::step(), and ReactorNet::updateState().

void addReactor ( ReactorBase r,
bool  iown = false 
)

Add the reactor r to this reactor network.

Definition at line 174 of file ReactorNet.cpp.

References Cantera::int2str(), ReactorBase::name(), ReactorBase::setNetwork(), ReactorBase::type(), and Cantera::writelog().

ReactorBase& 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 122 of file ReactorNet.h.

Referenced by ReactorNet::globalComponentIndex(), and ReactorNet::sensitivity().

bool verbose ( ) const
inline

Returns true if verbose logging output is enabled.

Definition at line 127 of file ReactorNet.h.

void setVerbose ( bool  v = true)
inline

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

Definition at line 133 of file ReactorNet.h.

Integrator& integrator ( )
inline

Return a reference to the integrator.

Definition at line 138 of file ReactorNet.h.

void updateState ( doublereal *  y)

Update the state of all the reactors in the network to correspond to the values in the solution vector y.

Definition at line 238 of file ReactorNet.cpp.

Referenced by ReactorNet::advance(), ReactorNet::eval(), and ReactorNet::step().

double sensitivity ( size_t  k,
size_t  p 
)
inline

Return the sensitivity of the k-th solution component with respect to the p-th sensitivity parameter.

Definition at line 148 of file ReactorNet.h.

References ReactorNet::initialize(), ReactorNet::m_sensIndex, and Integrator::solution().

Referenced by ReactorNet::sensitivity().

double sensitivity ( const std::string &  species,
size_t  p,
int  reactor = 0 
)
inline

Return the sensitivity of the species named species with respect to the p-th sensitivity parameter.

Definition at line 157 of file ReactorNet.h.

References ReactorNet::globalComponentIndex(), ReactorNet::reactor(), and ReactorNet::sensitivity().

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

Evaluate the Jacobian matrix for the reactor network.

Parameters
[in]tTime at which to evaluate the Jacobian
[in]yGlobal state vector at time t
[out]ydotTime derivative of the state vector evaluated at t.
[in]psensitivity parameter vector (unused?)
[out]jJacobian matrix, size neq() by neq().

Definition at line 204 of file ReactorNet.cpp.

References DATA_PTR, Cantera::error(), ReactorNet::eval(), and CanteraError::what().

virtual size_t neq ( )
inlinevirtual

Number of equations.

Implements FuncEval.

Definition at line 174 of file ReactorNet.h.

Referenced by ReactorNet::initialize().

void eval ( doublereal  t,
doublereal *  y,
doublereal *  ydot,
doublereal *  p 
)
virtual

Evaluate the right-hand-side 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()

Implements FuncEval.

Definition at line 188 of file ReactorNet.cpp.

References ReactorNet::updateState().

Referenced by ReactorNet::evalJacobian().

void getInitialConditions ( doublereal  t0,
size_t  leny,
doublereal *  y 
)
virtual

Fill the solution vector with the initial conditions at initial time t0.

Implements FuncEval.

Definition at line 247 of file ReactorNet.cpp.

virtual size_t nparams ( )
inlinevirtual

Number of sensitivity parameters.

Reimplemented from FuncEval.

Definition at line 181 of file ReactorNet.h.

size_t globalComponentIndex ( const std::string &  species,
size_t  reactor = 0 
)

Return the index corresponding to the species named species in the reactor with index reactor in the global state vector for the reactor network.

Definition at line 257 of file ReactorNet.cpp.

References ReactorNet::initialize(), and ReactorNet::reactor().

Referenced by ReactorNet::sensitivity().

void registerSensitivityReaction ( void *  reactor,
size_t  reactionIndex,
const std::string &  name,
int  leftright = 0 
)

Used by Reactor and Wall objects to register the addition of sensitivity reactions so that the ReactorNet can keep track of the order in which sensitivity parameters are added.

Definition at line 270 of file ReactorNet.cpp.

References ReactorNet::m_paramNames, and ReactorNet::m_sensOrder.

Referenced by Reactor::addSensitivityReaction().

const std::string& sensitivityParameterName ( size_t  p)
inline

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

Definition at line 197 of file ReactorNet.h.

References ReactorNet::m_paramNames.

void initialize ( )
protected

Member Data Documentation

std::vector<std::string> m_paramNames
protected

Names corresponding to each sensitivity parameter.

Definition at line 236 of file ReactorNet.h.

Referenced by ReactorNet::registerSensitivityReaction(), and ReactorNet::sensitivityParameterName().

std::map<std::pair<void*, int>, std::map<size_t, size_t> > m_sensOrder
protected

Structure used to determine the order of sensitivity parameters m_sensOrder[Reactor or Wall, leftright][reaction number] = parameter index.

Definition at line 240 of file ReactorNet.h.

Referenced by ReactorNet::initialize(), and ReactorNet::registerSensitivityReaction().

std::vector<size_t> m_sensIndex
protected

Mapping from the order in which sensitivity parameters were added to the ReactorNet to the order in which they occur in the integrator output.

Definition at line 245 of file ReactorNet.h.

Referenced by ReactorNet::initialize(), and ReactorNet::sensitivity().


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