14 #ifndef rrCvodeInterfaceH
15 #define rrCvodeInterfaceH
18 #include "rrRoadRunnerOptions.h"
22 #include <sundials/sundials_linearsolver.h>
23 #include <sundials/sundials_nonlinearsolver.h>
24 #include <cvodes/cvodes.h>
25 #include <sundials/sundials_nvector.h>
26 #include "ForwardSensitivitySolver.h"
32 class ExecutableModel;
36 class ForwardSensitivitySolver;
52 using Integrator::Integrator;
95 std::string
getName()
const override;
109 std::string
getHint()
const override;
171 double integrate(
double t0,
double hstep)
override;
179 void restart(
double timeStart)
override;
218 std::string
ToString(
int val)
const;
224 std::string
ToString(
size_t val)
const;
251 static const int mDefaultMaxNumSteps;
252 static const int mDefaultMaxAdamsOrder;
253 static const int mDefaultMaxBDFOrder;
256 void *mCVODE_Memory =
nullptr;
257 N_Vector mStateVector =
nullptr;
258 SUNMatrix jac =
nullptr;
259 SUNNonlinearSolver nonLinSolver =
nullptr;
260 SUNLinearSolver linSolver =
nullptr;
262 IntegratorListenerPtr listener;
263 double mLastEventTime;
264 bool variableStepPendingEvent;
265 bool variableStepTimeEndEvent;
266 std::vector<double> variableStepPostEventState;
267 std::vector<unsigned char> eventStatus;
269 void testRootsAtInitialTime();
271 bool haveVariables()
const;
273 void assignResultsToModel()
const;
281 void setCVODETolerances();
286 void reInit(
double t0);
300 void applyPendingEvents(
double timeEnd);
302 void applyEvents(
double timeEnd, std::vector<unsigned char> &previousEventStatus);
304 double applyVariableStepPendingEvents();
308 void freeSundialsMemory();
310 bool stateVectorVariables;
312 unsigned long typecode_;
314 friend int cvodeDyDtFcn(
double t, N_Vector cv_y, N_Vector cv_ydot,
void *f_data);
316 friend int cvodeEventAndPiecewiseRootFcn(
double t, N_Vector y,
double *gout,
void *g_data);
320 template <
class SundialsType = CVODEIntegrator>
321 std::string decodeSundialsError(SundialsType* solver,
int cvodeError,
bool exInfo) {
323 std::stringstream ss;
324 ss << (int) solver->getValue(
"maximum_num_steps");
325 std::string max_steps = ss.str();
327 switch (cvodeError) {
328 case CV_TOO_MUCH_WORK:
329 result =
"CV_TOO_MUCH_WORK";
331 result +=
": The solver took mxstep (" + max_steps +
") internal steps but " +
332 "could not reach tout.";
335 case CV_TOO_MUCH_ACC:
336 result =
"CV_TOO_MUCH_ACC";
338 result +=
": The solver could not satisfy the accuracy "
339 "demanded by the user for some internal step.";
343 result =
"CV_ERR_FAILURE";
345 result +=
": Error test failures occurred too many times "
346 "(= MXNEF = 7) during one internal time step or"
347 "occurred with |h| = hmin.";
350 case CV_CONV_FAILURE:
351 result =
"CV_CONV_FAILURE";
353 result +=
": Convergence test failures occurred too many "
354 "times (= MXNCF = 10) during one internal time"
355 "step or occurred with |h| = hmin.";
359 result =
"CV_LINIT_FAIL";
361 result +=
": The linear solver's initialization function "
366 result =
"CV_LSETUP_FAIL";
368 result +=
": The linear solver's setup routine failed in an "
369 "unrecoverable manner.";
373 result =
"CV_LSOLVE_FAIL";
375 result +=
": The linear solver's solve routine failed in an "
376 "unrecoverable manner.";
379 case CV_RHSFUNC_FAIL:
380 result =
"CV_RHSFUNC_FAIL";
382 case CV_FIRST_RHSFUNC_ERR:
383 result =
"CV_FIRST_RHSFUNC_ERR";
385 case CV_REPTD_RHSFUNC_ERR:
386 result =
"CV_REPTD_RHSFUNC_ERR";
388 case CV_UNREC_RHSFUNC_ERR:
389 result =
"CV_UNREC_RHSFUNC_ERR";
392 result =
"CV_RTFUNC_FAIL";
395 result =
"CV_MEM_FAIL";
398 result =
"CV_MEM_NULL";
400 result +=
": The cvode_mem argument was NULL.";
404 result =
"CV_ILL_INPUT";
406 result +=
": One of the inputs to CVode is illegal. This "
407 "includes the situation when a component of the "
408 "error weight vectors becomes < 0 during "
409 "internal time-stepping. It also includes the "
410 "situation where a root of one of the root "
411 "functions was found both at t0 and very near t0. "
412 "The ILL_INPUT flag will also be returned if the "
413 "linear solver routine CV--- (called by the user "
414 "after calling CVodeCreate) failed to set one of "
415 "the linear solver-related fields in cvode_mem or "
416 "if the linear solver's init routine failed. In "
417 "any case, the user should see the printed "
418 "error message for more details.";
422 result =
"CV_NO_MALLOC";
424 result +=
": indicating that cvode_mem has not been "
425 "allocated (i.e., CVodeInit has not been "
432 result +=
": k is not in the range 0, 1, ..., qu.";
438 result +=
": t is not in the interval [tn-hu,tn].";
442 result =
"CV_BAD_DKY";
444 result +=
": The dky argument was NULL.";
448 result =
"CV_TOO_CLOSE:";
451 result =
"UNKNOWN_CODE";
RoadRunner's Gillespie SSA integrator.
A RoadRunner integrator based on CVODE; serves as RoadRunner's main integrator for ODEs.
Definition: CVODEIntegrator.h:46
std::vector< double > getAbsoluteToleranceVector() override
Get the absolute tolerance vector for the solver.
Definition: CVODEIntegrator.cpp:797
void tweakTolerances() override
Fix tolerances for SBML tests.
Definition: CVODEIntegrator.cpp:578
std::string getDescription() const override
Get the description for this integrator.
Definition: CVODEIntegrator.cpp:294
double integrate(double t0, double hstep) override
Main integration routine.
Definition: CVODEIntegrator.cpp:436
N_Vector getStateVector() const
getter for the internal state std::vector
Definition: CVODEIntegrator.cpp:1092
std::string ToString(int val) const
Converts integer to std::string for error print.
Definition: CVODEIntegrator.cpp:141
IntegratorListenerPtr getListener() override
Gets the integrator listener.
Definition: CVODEIntegrator.cpp:118
void setIndividualTolerance(std::string sid, double value) override
Sets tolerance for individual species.
Definition: CVODEIntegrator.cpp:314
void setMaxOrder(int newValue)
sets the value of maximum order, which defaults to 12 for Adams (non-stiff) and 5 for BDF (Stiff).
Definition: CVODEIntegrator.cpp:346
void syncWithModel(ExecutableModel *m) override
Called whenever a new model is loaded to allow integrator to reset internal state.
Definition: CVODEIntegrator.cpp:198
SUNNonlinearSolver getSolver() const
getter for the internal Sundials linear solver object
Definition: CVODEIntegrator.cpp:1096
IntegrationMethod getIntegrationMethod() const override
Always deterministic for CVODE.
Definition: CVODEIntegrator.cpp:310
void checkType() const
Does a RT type check which throws if it fails, EVEN IF RTTI IS DISABLED.
Definition: CVODEIntegrator.cpp:122
void checkIndex(int index, int size) const
Does a index check which throws if it is out of bound.
Definition: CVODEIntegrator.cpp:134
Solver * construct(ExecutableModel *executableModel) const override
construct an instance of type CVODEIntegrator.
Definition: CVODEIntegrator.cpp:305
std::string getName() const override
Get the name for this integrator.
Definition: CVODEIntegrator.cpp:290
void resetSettings() override
Reset all integrator settings to their respective default values.
Definition: CVODEIntegrator.cpp:154
void loadConfigSettings() override
It looks like this method only get used inside resetSettings.
Definition: CVODEIntegrator.cpp:221
void setListener(IntegratorListenerPtr) override
Sets the integrator listener.
Definition: CVODEIntegrator.cpp:114
void loadSBMLSettings(const std::string &filename) override
Load an SBML settings file and apply the configuration options.
Definition: CVODEIntegrator.cpp:243
void restart(double timeStart) override
Reset time to zero and reinitialize model.
Definition: CVODEIntegrator.cpp:920
void * getCvodeMemory() const
getter for the internal CVode memory buffer
Definition: CVODEIntegrator.cpp:1088
~CVODEIntegrator() override
Destructor.
Definition: CVODEIntegrator.cpp:108
CVODEIntegrator(ExecutableModel *oModel)
Constructor: takes an executable model, does not own the pointer.
Definition: CVODEIntegrator.cpp:83
std::string getHint() const override
Get the hint for this integrator.
Definition: CVODEIntegrator.cpp:301
void checkVectorSize(int expected, size_t real) const
Does a size check which throws if it fails.
Definition: CVODEIntegrator.cpp:127
Base class for all code generation systems; allows compiling and evaluating the model.
Definition: rrExecutableModel.h:118
Time based sensivitity solver.
Definition: ForwardSensitivitySolver.h:31
Integrator is an abstract base class that provides an interface to specific integrator class implemen...
Definition: Integrator.h:60
virtual void setValue(const std::string &key, Setting value)
Pull down the setValue from superclass.
Definition: Solver.cpp:125
Store a roadrunner option (or setting) as a Variant type.
Definition: Setting.h:78
Base class for all integrators and steady state solvers.
Definition: Solver.h:39