ReactorNet.h Source File#
ReactorNet.h
Go to the documentation of this file.
A map of string keys to values whose type can vary at runtime.
Definition AnyMap.h:427
A class for 2D arrays stored in column-major (Fortran-compatible) form.
Definition Array.h:32
Virtual base class for ODE/DAE right-hand-side function evaluators.
Definition FuncEval.h:32
vector< double > m_sens_params
Values for the problem parameters for which sensitivities are computed This is the array which is per...
Definition FuncEval.h:181
void setPreconditioner(shared_ptr< PreconditionerBase > preconditioner)
Set preconditioner used by the linear solver.
Definition ReactorNet.cpp:157
void setLinearSolverType(const string &linSolverType="DENSE")
Set the type of linear solver used in the integration.
Definition ReactorNet.cpp:151
void preconditionerSetup(double t, double *y, double gamma) override
Evaluate the setup processes for the Jacobian preconditioner.
Definition ReactorNet.cpp:528
double step()
Advance the state of all reactors with respect to the independent variable (time or space).
Definition ReactorNet.cpp:240
virtual int lastOrder() const
Returns the order used for last solution step of the ODE integrator The function is intended for inte...
Definition ReactorNet.cpp:275
const string & sensitivityParameterName(size_t p) const
The name of the p-th sensitivity parameter added to this ReactorNet.
Definition ReactorNet.h:256
void eval(double t, double *y, double *ydot, double *p) override
Evaluate the right-hand-side ODE function.
Definition ReactorNet.cpp:318
void advance(double t)
Advance the state of all reactors in the independent variable (time or space).
Definition ReactorNet.cpp:173
vector< size_t > m_start
m_start[n] is the starting point in the state vector for reactor n
Definition ReactorNet.h:333
vector< double > m_LHS
m_LHS is a vector representing the coefficients on the "left hand side" of each governing equation
Definition ReactorNet.h:359
double m_initial_time
The initial value of the independent variable in the system.
Definition ReactorNet.h:326
void evalJacobian(double t, double *y, double *ydot, double *p, Array2D *j)
Evaluate the Jacobian matrix for the reactor network.
Definition ReactorNet.cpp:376
double time()
Current value of the simulation time [s], for reactor networks that are solved in the time domain.
Definition ReactorNet.cpp:68
double getInitialTime() const
Get the initial value of the independent variable (typically time).
Definition ReactorNet.h:58
void setNeedsReinit()
Called to trigger integrator reinitialization before further integration.
Definition ReactorNet.h:272
void getConstraints(double *constraints) override
Given a vector of length neq(), mark which variables should be considered algebraic constraints.
Definition ReactorNet.cpp:353
Reactor & reactor(int n)
Return a reference to the n-th reactor in this network.
Definition ReactorNet.h:144
virtual void setMaxSteps(int nmax)
Set the maximum number of internal integration steps the integrator will take before reaching the nex...
Definition ReactorNet.cpp:163
string componentName(size_t i) const
Return the name of the i-th component of the global state vector.
Definition ReactorNet.cpp:467
void getStateDae(double *y, double *ydot) override
Fill in the vectors y and ydot with the current state of the system.
Definition ReactorNet.cpp:452
void setInitialTime(double time)
Set the initial value of the independent variable (typically time).
Definition ReactorNet.cpp:28
virtual void getDerivative(int k, double *dky)
Return k-th derivative at the current state of the system.
Definition ReactorNet.cpp:406
void setMaxErrTestFails(int nmax)
Set the maximum number of error test failures permitted by the CVODES integrator in a single step.
Definition ReactorNet.cpp:41
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 React...
Definition ReactorNet.cpp:479
double m_maxstep
Maximum integrator internal timestep. Default of 0.0 means infinity.
Definition ReactorNet.h:344
double distance()
Current position [m] along the length of the reactor network, for reactors that are solved as a funct...
Definition ReactorNet.cpp:77
void setSensitivityTolerances(double rtol, double atol)
Set the relative and absolute tolerances for integrating the sensitivity equations.
Definition ReactorNet.cpp:57
int maxSteps()
Returns the maximum number of internal integration steps the integrator will take before reaching the...
Definition ReactorNet.cpp:168
virtual void setDerivativeSettings(AnyMap &settings)
Set derivative settings of all reactors.
Definition ReactorNet.cpp:493
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.
Definition ReactorNet.cpp:360
void updateState(double *y)
Update the state of all the reactors in the network to correspond to the values in the solution vecto...
Definition ReactorNet.cpp:398
void getState(double *y) override
Fill in the vector y with the current state of the system.
Definition ReactorNet.cpp:445
void setAdvanceLimits(const double *limits)
Set absolute step size limits during advance.
Definition ReactorNet.cpp:417
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 ...
Definition ReactorNet.cpp:459
bool m_timeIsIndependent
Indicates whether time or space is the independent variable.
Definition ReactorNet.h:349
bool hasAdvanceLimits() const
Check whether ReactorNet object uses advance limits.
Definition ReactorNet.cpp:427
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 paramete...
Definition ReactorNet.h:190
void evalDae(double t, double *y, double *ydot, double *p, double *residual) override
eval coupling for IDA / DAEs
Definition ReactorNet.cpp:341
virtual void checkPreconditionerSupported() const
Check that preconditioning is supported by all reactors in the network.
Definition ReactorNet.cpp:571
bool getAdvanceLimits(double *limits) const
Retrieve absolute step size limits during advance.
Definition ReactorNet.cpp:436
void setVerbose(bool v=true)
Enable or disable verbose logging while setting up and integrating the reactor network.
Definition ReactorNet.h:155
void updatePreconditioner(double gamma) override
Update the preconditioner based on already computed jacobian values.
Definition ReactorNet.cpp:560
void preconditionerSolve(double *rhs, double *output) override
Evaluate the linear system Ax=b where A is the preconditioner.
Definition ReactorNet.cpp:519
vector< string > m_paramNames
Names corresponding to each sensitivity parameter.
Definition ReactorNet.h:352
void setTolerances(double rtol, double atol)
Set the relative and absolute tolerances for the integrator.
Definition ReactorNet.cpp:46
virtual void getEstimate(double time, int k, double *yest)
Estimate a future state based on current derivatives.
Definition ReactorNet.cpp:252
Class Reactor is a general-purpose class for stirred reactors.
Definition Reactor.h:44
void scale(InputIter begin, InputIter end, OutputIter out, S scale_factor)
Multiply elements of an array by a scale factor.
Definition utilities.h:104
Generated by