Cantera  2.5.1
DAE_Solver.h
Go to the documentation of this file.
1 /**
2  * @file DAE_Solver.h
3  * Header file for class DAE_Solver
4  */
5 
6 // This file is part of Cantera. See License.txt in the top-level directory or
7 // at https://cantera.org/license.txt for license and copyright information.
8 
9 #ifndef CT_DAE_Solver_H
10 #define CT_DAE_Solver_H
11 
12 #include "ResidJacEval.h"
13 #include "cantera/base/global.h"
14 
15 namespace Cantera
16 {
17 
18 class Jacobian
19 {
20 public:
21  Jacobian() {}
22  virtual ~Jacobian() {}
23  virtual bool supplied() {
24  return false;
25  }
26  virtual bool isBanded() {
27  return false;
28  }
29  virtual int lowerBandWidth() {
30  return 0;
31  }
32  virtual int upperBandWidth() {
33  return 0;
34  }
35 };
36 
37 class BandedJacobian : public Jacobian
38 {
39 public:
40  BandedJacobian(int ml, int mu) {
41  m_ml = ml;
42  m_mu = mu;
43  }
44  virtual bool supplied() {
45  return false;
46  }
47  virtual bool isBanded() {
48  return true;
49  }
50  virtual int lowerBandWidth() {
51  return m_ml;
52  }
53  virtual int upperBandWidth() {
54  return m_mu;
55  }
56 protected:
57  int m_ml, m_mu;
58 };
59 
60 const int cDirect = 0;
61 const int cKrylov = 1;
62 
63 
64 /**
65  * Wrapper for DAE solvers
66  *
67  * @attention This class currently does not have any test cases or examples. Its
68  * implementation may be incomplete, and future changes to Cantera may
69  * unexpectedly cause this class to stop working. If you use this class,
70  * please consider contributing examples or test cases. In the absence of
71  * new tests or examples, this class may be deprecated and removed in a
72  * future version of Cantera. See
73  * https://github.com/Cantera/cantera/issues/267 for additional information.
74  */
76 {
77 public:
79  m_resid(f),
80  m_neq(f.nEquations()),
81  m_time(0.0) {
82  }
83 
84  virtual ~DAE_Solver() {}
85 
86  /**
87  * Set error tolerances. This version specifies a scalar
88  * relative tolerance, and a vector absolute tolerance.
89  */
90  virtual void setTolerances(doublereal reltol,
91  doublereal* abstol) {
92  warn("setTolerances");
93  }
94 
95  /**
96  * Set error tolerances. This version specifies a scalar
97  * relative tolerance, and a scalar absolute tolerance.
98  */
99  virtual void setTolerances(doublereal reltol, doublereal abstol) {
100  warn("setTolerances");
101  }
102 
103  /**
104  * Specify a Jacobian evaluator. If this method is not called,
105  * the Jacobian will be computed by finite difference.
106  */
107  void setJacobian(Jacobian& jac) {
108  warn("setJacobian");
109  }
110 
111  virtual void setLinearSolverType(int solverType) {
112  warn("setLinearSolverType");
113  }
114 
115  virtual void setDenseLinearSolver() {
116  warn("setDenseLinearSolver");
117  }
118 
119  virtual void setBandedLinearSolver(int m_upper, int m_lower) {
120  warn("setBandedLinearSolver");
121  }
122  virtual void setMaxStepSize(doublereal dtmax) {
123  warn("setMaxStepSize");
124  }
125  virtual void setMaxOrder(int n) {
126  warn("setMaxOrder");
127  }
128  virtual void setMaxNumSteps(int n) {
129  warn("setMaxNumSteps");
130  }
131  virtual void setInitialStepSize(doublereal h0) {
132  warn("setInitialStepSize");
133  }
134  virtual void setStopTime(doublereal tstop) {
135  warn("setStopTime");
136  }
137  virtual void setMaxErrTestFailures(int n) {
138  warn("setMaxErrTestFailures");
139  }
140  virtual void setMaxNonlinIterations(int n) {
141  warn("setMaxNonlinIterations");
142  }
143  virtual void setMaxNonlinConvFailures(int n) {
144  warn("setMaxNonlinConvFailures");
145  }
146  virtual void inclAlgebraicInErrorTest(bool yesno) {
147  warn("inclAlgebraicInErrorTest");
148  }
149 
150  //! Calculate consistent value of the starting solution given the starting
151  //! solution derivatives
152  /**
153  * This method may be called if the initial conditions do not satisfy the
154  * residual equation F = 0. Given the derivatives of all variables, this
155  * method computes the initial y values.
156  */
157  virtual void correctInitial_Y_given_Yp(doublereal* y, doublereal* yp,
158  doublereal tout) {
159  warn("correctInitial_Y_given_Yp");
160  }
161 
162  //! Calculate consistent value of the algebraic constraints and derivatives
163  //! at the start of the problem
164  /**
165  * This method may be called if the initial conditions do not satisfy the
166  * residual equation F = 0. Given the initial values of all differential
167  * variables, it computes the initial values of all algebraic variables and
168  * the initial derivatives of all differential variables.
169  * @param y Calculated value of the solution vector after the procedure ends
170  * @param yp Calculated value of the solution derivative after the procedure
171  * @param tout The first value of t at which a soluton will be
172  * requested (from IDASolve). (This is needed here to
173  * determine the direction of integration and rough scale
174  * in the independent variable t.
175  */
176  virtual void correctInitial_YaYp_given_Yd(doublereal* y, doublereal* yp,
177  doublereal tout) {
178  warn("correctInitial_YaYp_given_Yd");
179  }
180 
181  /**
182  * Solve the system of equations up to time tout.
183  */
184  virtual int solve(doublereal tout) {
185  warn("solve");
186  return 0;
187  }
188 
189  /**
190  * Take one internal step.
191  */
192  virtual doublereal step(doublereal tout) {
193  warn("step");
194  return 0;
195  }
196 
197  /// Number of equations.
198  int nEquations() const {
199  return m_resid.nEquations();
200  }
201 
202  /**
203  * initialize. Base class method does nothing.
204  */
205  virtual void init(doublereal t0) {}
206 
207  /**
208  * Set a solver-specific input parameter.
209  */
210  virtual void setInputParameter(int flag, doublereal value) {
211  warn("setInputParameter");
212  }
213 
214  /**
215  * Get the value of a solver-specific output parameter.
216  */
217  virtual doublereal getOutputParameter(int flag) const {
218  warn("getOutputParameter");
219  return 0.0;
220  }
221 
222  /// the current value of solution component k.
223  virtual doublereal solution(int k) const {
224  warn("solution");
225  return 0.0;
226  }
227 
228  virtual const doublereal* solutionVector() const {
229  warn("solutionVector");
230  return &m_dummy;
231  }
232 
233  /// the current value of the derivative of solution component k.
234  virtual doublereal derivative(int k) const {
235  warn("derivative");
236  return 0.0;
237  }
238 
239  virtual const doublereal* derivativeVector() const {
240  warn("derivativeVector");
241  return &m_dummy;
242  }
243 
244 protected:
245  doublereal m_dummy;
246  ResidJacEval& m_resid;
247 
248  //! Number of total equations in the system
249  integer m_neq;
250  doublereal m_time;
251 
252 private:
253  void warn(const std::string& msg) const {
254  writelog(">>>> Warning: method "+msg+" of base class "
255  +"DAE_Solver called. Nothing done.\n");
256  }
257 };
258 
259 
260 //! Factor method for choosing a DAE solver
261 /*!
262  * @param itype String identifying the type (IDA is the only option)
263  * @param f Residual function to be solved by the DAE algorithm
264  * @returns a point to the instantiated DAE_Solver object
265  */
266 DAE_Solver* newDAE_Solver(const std::string& itype, ResidJacEval& f);
267 
268 }
269 
270 #endif
Dense, Square (not sparse) matrices.
Wrapper for DAE solvers.
Definition: DAE_Solver.h:76
virtual doublereal step(doublereal tout)
Take one internal step.
Definition: DAE_Solver.h:192
virtual doublereal derivative(int k) const
the current value of the derivative of solution component k.
Definition: DAE_Solver.h:234
virtual void setTolerances(doublereal reltol, doublereal abstol)
Set error tolerances.
Definition: DAE_Solver.h:99
virtual void init(doublereal t0)
initialize.
Definition: DAE_Solver.h:205
virtual void correctInitial_Y_given_Yp(doublereal *y, doublereal *yp, doublereal tout)
Calculate consistent value of the starting solution given the starting solution derivatives.
Definition: DAE_Solver.h:157
virtual void setTolerances(doublereal reltol, doublereal *abstol)
Set error tolerances.
Definition: DAE_Solver.h:90
int nEquations() const
Number of equations.
Definition: DAE_Solver.h:198
integer m_neq
Number of total equations in the system.
Definition: DAE_Solver.h:249
virtual int solve(doublereal tout)
Solve the system of equations up to time tout.
Definition: DAE_Solver.h:184
virtual doublereal solution(int k) const
the current value of solution component k.
Definition: DAE_Solver.h:223
virtual void setInputParameter(int flag, doublereal value)
Set a solver-specific input parameter.
Definition: DAE_Solver.h:210
void setJacobian(Jacobian &jac)
Specify a Jacobian evaluator.
Definition: DAE_Solver.h:107
virtual void correctInitial_YaYp_given_Yd(doublereal *y, doublereal *yp, doublereal tout)
Calculate consistent value of the algebraic constraints and derivatives at the start of the problem.
Definition: DAE_Solver.h:176
virtual doublereal getOutputParameter(int flag) const
Get the value of a solver-specific output parameter.
Definition: DAE_Solver.h:217
Wrappers for the function evaluators for Nonlinear solvers and Time steppers.
Definition: ResidJacEval.h:56
virtual int nEquations() const
Return the number of equations in the equation system.
This file contains definitions for utility functions and text for modules, inputfiles,...
void writelog(const std::string &fmt, const Args &... args)
Write a formatted message to the screen.
Definition: global.h:158
Namespace for the Cantera kernel.
Definition: AnyMap.cpp:264
DAE_Solver * newDAE_Solver(const std::string &itype, ResidJacEval &f)
Factor method for choosing a DAE solver.
Definition: DAE_solvers.cpp:16