6#ifndef CT_STEADYSTATESYSTEM_H
7#define CT_STEADYSTATESYSTEM_H
38 virtual void eval(
double* x,
double* r,
double rdt=-1.0,
int count=1) = 0;
73 int solve(
double* x0,
double* x1,
int loglevel);
81 void solve(
int loglevel=0);
88 double ssnorm(
double* x,
double* r);
107 virtual string componentName(
size_t i)
const {
return fmt::format(
"{}", i); }
114 return {
"",
"component name"};
145 return (
m_rdt != 0.0);
150 return (
m_rdt == 0.0);
176 double timeStep(
int nsteps,
double dt,
double* x,
double* r,
int loglevel);
191 void setTimeStep(
double stepsize,
size_t n,
const int* tsteps);
227 void setJacAge(
int ss_age,
int ts_age=-1);
261 int loglevel,
int attempt_counter) {}
Declarations for class SystemJacobian.
Base class for 'functor' classes that evaluate a function of one variable.
Newton iterator for multi-domain, one-dimensional problems.
Base class for representing a system of differential-algebraic equations and solving for its steady-s...
int solve(double *x0, double *x1, int loglevel)
Solve , where is the residual function.
int m_nsteps
Number of time steps taken in the current call to solve()
size_t m_size
Solution vector size
int m_nsteps_max
Maximum number of timesteps allowed per call to solve()
virtual void resize()
Call to set the size of internal data structures after first defining the system or if the problem si...
shared_ptr< SystemJacobian > linearSolver() const
Get the the linear solver being used to hold the Jacobian matrix and solve linear systems as part of ...
vector< double > m_xnew
Work array used to hold the residual or the new solution.
unique_ptr< MultiNewton > m_newt
Newton iterator.
virtual void resetBadValues(double *x)
Reset values such as negative species concentrations.
double m_jacobianAbsPerturb
Absolute perturbation of each component in finite difference Jacobian.
size_t size() const
Total solution vector length;.
virtual double weightedNorm(const double *step) const =0
Compute the weighted norm of step.
double ssnorm(double *x, double *r)
Steady-state max norm (infinity norm) of the residual evaluated using solution x.
vector< int > & transientMask()
Access the vector indicating which equations contain a transient term.
virtual double upperBound(size_t i) const =0
Get the upper bound for global component i in the state vector.
double rdt() const
Reciprocal of the time step.
virtual void initTimeInteg(double dt, double *x)
Prepare for time stepping beginning with solution x and timestep dt.
virtual string componentName(size_t i) const
Get the name of the i-th component of the state vector.
virtual string componentTableLabel(size_t i) const
Get elements of the component name, aligned with the column headings given by componentTableHeader().
size_t bandwidth() const
Jacobian bandwidth.
double m_rdt
Reciprocal of time step.
double m_jacobianThreshold
Threshold for ignoring small elements in Jacobian.
void setMaxTimeStep(double tmax)
Set the maximum time step allowed during time stepping.
shared_ptr< SystemJacobian > m_jac
Jacobian evaluator.
shared_ptr< vector< double > > m_state
Solution vector.
void setMinTimeStep(double tmin)
Set the minimum time step allowed during time stepping.
vector< int > m_mask
Transient mask.
void setTimeStepCallback(Func1 *callback)
Set a function that will be called after each successful timestep.
Func1 * m_interrupt
Function called at the start of every call to eval.
virtual double lowerBound(size_t i) const =0
Get the lower bound for global component i in the state vector.
bool transient() const
True if transient mode.
void setMaxTimeStepCount(int nmax)
Set the maximum number of timeteps allowed before successful steady-state solve.
void setTimeStepFactor(double tfactor)
Sets a factor by which the time step is reduced if the time stepping fails.
vector< int > m_steps
Array of number of steps to take after each unsuccessful steady-state solve before re-attempting the ...
size_t m_bw
Jacobian bandwidth.
double m_tfactor
Factor time step is multiplied by if time stepping fails ( < 1 )
bool m_jac_ok
If true, Jacobian is current.
double m_tstep
Initial timestep.
int maxTimeStepCount() const
Get the maximum number of timeteps allowed before successful steady-state solve.
int m_ts_jac_age
Maximum age of the Jacobian in time-stepping mode.
virtual void eval(double *x, double *r, double rdt=-1.0, int count=1)=0
Evaluate the residual function.
double timeStep(int nsteps, double dt, double *x, double *r, int loglevel)
Take time steps using Backward Euler.
void setJacAge(int ss_age, int ts_age=-1)
Set the maximum number of steps that can be taken using the same Jacobian before it must be re-evalua...
void evalSSJacobian(double *x, double *rsd)
Evaluate the steady-state Jacobian, accessible via linearSolver()
virtual void writeDebugInfo(const string &header_suffix, const string &message, int loglevel, int attempt_counter)
Write solver debugging based on the specified log level.
double m_tmin
Minimum timestep size.
void setJacobianPerturbation(double relative, double absolute, double threshold)
Configure perturbations used to evaluate finite difference Jacobian.
double m_tmax
Maximum timestep size.
int m_ss_jac_age
Maximum age of the Jacobian in steady-state mode.
virtual void setSteadyMode()
Prepare to solve the steady-state problem.
virtual void clearDebugFile()
Delete debug output file that may be created by writeDebugInfo() when solving with high loglevel.
int m_attempt_counter
Counter used to manage the number of states stored in the debug log file generated by writeDebugInfo(...
void setLinearSolver(shared_ptr< SystemJacobian > solver)
Set the linear solver used to hold the Jacobian matrix and solve linear systems as part of each Newto...
void setTimeStep(double stepsize, size_t n, const int *tsteps)
Set the number of time steps to try when the steady Newton solver is unsuccessful.
int m_max_history
Constant that determines the maximum number of states stored in the debug log file generated by write...
void setInitialGuess(const double *x)
Set the initial guess. Should be called before solve().
void getState(double *x) const
Get the converged steady-state solution after calling solve().
Func1 * m_time_step_callback
User-supplied function called after each successful timestep.
virtual pair< string, string > componentTableHeader() const
Get header lines describing the column names included in a component label.
virtual void evalJacobian(double *x0)=0
Evaluates the Jacobian at x0 using finite differences.
double m_jacobianRelPerturb
Relative perturbation of each component in finite difference Jacobian.
vector< double > m_xlast_ts
State vector after the last successful set of time steps.
bool steady() const
True if steady mode.
vector< double > m_work1
Work arrays used during Jacobian evaluation.
MultiNewton & newton()
Return a reference to the Newton iterator.
void setInterrupt(Func1 *interrupt)
Set a function that will be called every time eval is called.
This file contains definitions of constants, types and terms that are used in internal routines and a...
Namespace for the Cantera kernel.