C API Documentation
CVODEIntegrator.h
Go to the documentation of this file.
1 // == PREAMBLE ================================================
2 
3 // * Licensed under the Apache License, Version 2.0; see README
4 
5 // == FILEDOC =================================================
6 
13 //todo replace the include guards to match the filename
14 #ifndef rrCvodeInterfaceH
15 #define rrCvodeInterfaceH
16 
17 // == INCLUDES ================================================
18 
19 #include "Integrator.h"
20 #include "rrRoadRunnerOptions.h"
21 
22 #include <string>
23 #include <vector>
24 #include <sundials/sundials_linearsolver.h>
25 #include <sundials/sundials_nonlinearsolver.h>
26 
27 // == CODE ====================================================
28 
32 typedef struct _generic_N_Vector *N_Vector;
33 
34 namespace rr
35 {
36  using std::string;
37 
38  class ExecutableModel;
39  class RoadRunner;
40 
49  class CVODEIntegrator : public Integrator
50  {
51 
52  public:
57  explicit CVODEIntegrator(ExecutableModel* oModel);
58 
63  ~CVODEIntegrator() override;
64 
65 
71  void syncWithModel(ExecutableModel* m) override;
72 
73  // ** Loading Settings *************************************************
74 
75  void loadConfigSettings() override;
76 
82  void loadSBMLSettings(const std::string& filename) override;
83 
84  // ** Meta Info ********************************************************
85 
91  std::string getName() const override;
92 
97  static std::string getCVODEIntegratorName();
98 
104  std::string getDescription() const override;
105 
110  static std::string getCVODEIntegratorDescription();
111 
117  std::string getHint() const override;
118 
123  static std::string getCVODEIntegratorHint();
124 
125  // ** Getters / Setters ************************************************
126 
131  IntegrationMethod getIntegrationMethod() const override;
132 
137  void setValue(std::string setting, const Variant& value) override;
138 
139 
144  void setIndividualTolerance(string sid, double value) override;
145 
153  void setConcentrationTolerance(const Variant& value) override;
154 
159  std::vector<double> getConcentrationTolerance() override;
160 
165  void resetSettings() override;
166 
167 
176  void tweakTolerances() override;
177 
178 
179  // ** Integration Routines *********************************************
180 
185  double integrate(double t0, double tf) override;
186 
193  void restart(double timeStart) override;
194 
195  // ** Listeners ********************************************************
196 
201  IntegratorListenerPtr getListener() override;
202 
207  void setListener(IntegratorListenerPtr) override;
208 
213  void checkType() const;
214 
215 
220  void checkVectorSize(int expected, size_t real) const;
221 
226  void checkIndex(int index, int size) const;
227 
232  std::string ToString(int val) const;
233 
238  std::string ToString(size_t val) const;
239 
243  std::string cvodeDecodeError(int cvodeError, bool exInfo = true);
244 
245 
255 
264  SUNNonlinearSolver getSolver() const;
265 
274  void *getCvodeMemory() const;
275 
276  private:
277  static const int mDefaultMaxNumSteps;
278  static const int mDefaultMaxAdamsOrder;
279  static const int mDefaultMaxBDFOrder;
280 
281  ExecutableModel* mModel;
282 
283  // cvode components
284  void* mCVODE_Memory;
285  N_Vector mStateVector;
286  SUNMatrix jac = nullptr;
287  SUNNonlinearSolver nonLinSolver = nullptr;
288  SUNLinearSolver linSolver = nullptr;
289 
290  IntegratorListenerPtr listener;
291  double lastEventTime;
292  bool variableStepPendingEvent;
293  bool variableStepTimeEndEvent;
294  std::vector<double> variableStepPostEventState;
295  std::vector<unsigned char> eventStatus;
296 
297  void testRootsAtInitialTime();
298  bool haveVariables() const;
299  void assignResultsToModel() const;
305  void setCVODETolerances();
306 
310  void reInit(double t0);
311 
322  void updateCVODE();
323 
324  void applyPendingEvents(double timeEnd);
325 
326  void applyEvents(double timeEnd, std::vector<unsigned char> &previousEventStatus);
327 
328  double applyVariableStepPendingEvents();
329 
330  void createCVode();
331 
332  void freeCVode();
333 
334  bool stateVectorVariables;
335 
336  unsigned long typecode_;
337 
338  friend int cvodeDyDtFcn(double t, N_Vector cv_y, N_Vector cv_ydot, void *f_data);
339 
340  friend int cvodeRootFcn(double t, N_Vector y, double *gout, void *g_data);
341 
342 
343 
344  };
345 
346 
347  // ** Registration *********************************************************
348 
349 
351  public:
356  std::string getName() const override {
358  }
359 
364  std::string getDescription() const override {
366  }
367 
372  std::string getHint() const override {
374  }
375 
380  Integrator* construct(ExecutableModel *model) const override {
381  return new CVODEIntegrator(model);
382  }
383  };
384 }
385 
386 #endif
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
Definition: Variant.h:75