14 #ifndef rrCvodeInterfaceH
15 #define rrCvodeInterfaceH
20 #include "rrRoadRunnerOptions.h"
24 #include <sundials/sundials_linearsolver.h>
25 #include <sundials/sundials_nonlinearsolver.h>
38 class ExecutableModel;
75 void loadConfigSettings()
override;
277 static const int mDefaultMaxNumSteps;
278 static const int mDefaultMaxAdamsOrder;
279 static const int mDefaultMaxBDFOrder;
286 SUNMatrix jac =
nullptr;
287 SUNNonlinearSolver nonLinSolver =
nullptr;
288 SUNLinearSolver linSolver =
nullptr;
290 IntegratorListenerPtr listener;
291 double lastEventTime;
292 bool variableStepPendingEvent;
293 bool variableStepTimeEndEvent;
294 std::vector<double> variableStepPostEventState;
295 std::vector<unsigned char> eventStatus;
297 void testRootsAtInitialTime();
298 bool haveVariables()
const;
299 void assignResultsToModel()
const;
305 void setCVODETolerances();
310 void reInit(
double t0);
324 void applyPendingEvents(
double timeEnd);
326 void applyEvents(
double timeEnd, std::vector<unsigned char> &previousEventStatus);
328 double applyVariableStepPendingEvents();
334 bool stateVectorVariables;
336 unsigned long typecode_;
338 friend int cvodeDyDtFcn(
double t,
N_Vector cv_y,
N_Vector cv_ydot,
void *f_data);
340 friend int cvodeRootFcn(
double t,
N_Vector y,
double *gout,
void *g_data);
struct _generic_N_Vector * N_Vector
Definition: CVODEIntegrator.h:32
RoadRunner's Gillespie SSA integrator.
A RoadRunner integrator based on CVODE; serves as RoadRunner's main integrator for ODEs.
Definition: CVODEIntegrator.h:50
void tweakTolerances() override
Fix tolerances for SBML tests.
static std::string getCVODEIntegratorName()
Get the name for this integrator.
std::string getDescription() const override
Get the description for this integrator.
N_Vector getStateVector() const
getter for the internal state vector
void setValue(std::string setting, const Variant &value) override
Sets the value of an integrator setting (e.g. absolute_tolerance)
IntegrationMethod getIntegrationMethod() const override
Always deterministic for CVODE.
std::string ToString(int val) const
Converts integer to string for error print.
static std::string getCVODEIntegratorHint()
Get the hint for this integrator.
IntegratorListenerPtr getListener() override
Gets the integrator listener.
void * getCvodeMemory() const
getter for the internal CVode memory buffer
std::string ToString(size_t val) const
Converts size_t to string for error print.
void syncWithModel(ExecutableModel *m) override
Called whenever a new model is loaded to allow integrator to reset internal state.
SUNNonlinearSolver getSolver() const
getter for the internal Sundials linear solver object
void checkType() const
Does a RT type check which throws if it fails, EVEN IF RTTI IS DISABLED.
void checkIndex(int index, int size) const
Does a index check which throws if it is out of bound.
std::string getName() const override
Get the name for this integrator.
std::string cvodeDecodeError(int cvodeError, bool exInfo=true)
decode the cvode error code to a string
double integrate(double t0, double tf) override
Main integration routine.
void resetSettings() override
Reset all integrator settings to their respective default values.
std::vector< double > getConcentrationTolerance() override
Gets tolerance based on concentration of species.
void setListener(IntegratorListenerPtr) override
Sets the integrator listener.
void loadSBMLSettings(const std::string &filename) override
Load an SBML settings file and apply the configuration options.
void restart(double timeStart) override
Reset time to zero and reinitialize model.
void setIndividualTolerance(string sid, double value) override
Sets tolerance for individual species.
~CVODEIntegrator() override
Destructor.
CVODEIntegrator(ExecutableModel *oModel)
Constructor: takes an executable model, does not own the pointer.
std::string getHint() const override
Get the hint for this integrator.
static std::string getCVODEIntegratorDescription()
Get the description for this integrator.
void setConcentrationTolerance(const Variant &value) override
Sets tolerance based on concentration of species.
void checkVectorSize(int expected, size_t real) const
Does a size check which throws if it fails.
Definition: CVODEIntegrator.h:350
std::string getName() const override
Gets the name associated with this integrator type.
Definition: CVODEIntegrator.h:356
std::string getHint() const override
Gets the hint associated with this integrator type.
Definition: CVODEIntegrator.h:372
Integrator * construct(ExecutableModel *model) const override
Constructs a new integrator of a given type.
Definition: CVODEIntegrator.h:380
std::string getDescription() const override
Gets the description associated with this integrator type.
Definition: CVODEIntegrator.h:364
Base class for all code generation systems; allows compiling and evaluating the model.
Definition: rrExecutableModel.h:118
Definition: Integrator.h:62
Handles constructing an integrator and contains meta information about it.
Definition: Integrator.h:161