C API Documentation
RK45Integrator.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 
14 #ifndef RK45Integrator_H_
15 #define RK45Integrator_H_
16 
17 #include <Integrator.h>
18 #include <rrRoadRunnerOptions.h>
19 
20 namespace rr
21 {
22 
23 
30  class RK45Integrator: public Integrator
31  {
32  public:
38 
43  virtual ~RK45Integrator();
44 
50  virtual void syncWithModel(ExecutableModel* m);
51 
52 
56  public:
57 
69  virtual double integrate(double t, double h);
70 
74  virtual void restart(double t0);
75 
76  // ** Meta Info ********************************************************
77 
83  std::string getName() const;
84 
89  static std::string getRK45Name();
90 
96  std::string getDescription() const;
97 
102  static std::string getRK45Description();
103 
109  std::string getHint() const;
110 
115  static std::string getRK45Hint();
116 
117  // ** Getters / Setters ************************************************
118 
119  virtual Variant getValue(std::string key);
120 
125  IntegrationMethod getIntegrationMethod() const;
126 
132 
133  // ** Listeners ********************************************************
134 
139  virtual void setListener(IntegratorListenerPtr);
140 
144  virtual IntegratorListenerPtr getListener();
145 
146  public:
147 
148 
149  private:
150  ExecutableModel *model;
151 
152  unsigned stateVectorSize;
153 
157  double *k1, *k2, *k3, *k4, *y, *ytmp;
158 
159  // new scalars for RK45
160  double hCurrent, hmin, hmax;
161 
162  double *k5, *k6, *err;
163 
164  void testRootsAtInitialTime();
165  void applyEvents(double timeEnd, std::vector<unsigned char> &previousEventStatus);
166 
167  };
168 
169 
170  // ** Registration *********************************************************
171 
172 
174  public:
179  virtual std::string getName() const {
181  }
182 
187  virtual std::string getDescription() const {
189  }
190 
195  virtual std::string getHint() const {
197  }
198 
203  virtual Integrator* construct(ExecutableModel *model) const {
204  return new RK45Integrator(model);
205  }
206  };
207 
208 } /* namespace rr */
209 
210 #endif /* RK45Integrator_H_ */
RoadRunner's Gillespie SSA integrator.
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
A Runge-Kutta Fehlberg method for roadrunner.
Definition: RK45Integrator.h:31
virtual double integrate(double t, double h)
Integrates the model from t to t + h.
virtual void syncWithModel(ExecutableModel *m)
Called whenever a new model is loaded to allow integrator to reset internal state.
std::string getDescription() const
Get the description for this integrator.
RK45Integrator(ExecutableModel *m)
Constructor: takes an executable model, does not own the pointer.
virtual ~RK45Integrator()
Destructor.
virtual void setListener(IntegratorListenerPtr)
void resetSettings()
Reset all integrator settings to their respective default values.
static std::string getRK45Name()
Get the name for this integrator.
static std::string getRK45Hint()
Get the hint for this integrator.
static std::string getRK45Description()
Get the description for this integrator.
IntegrationMethod getIntegrationMethod() const
Always deterministic for RK45.
virtual void restart(double t0)
Restarts the integrator.
std::string getHint() const
Get the hint for this integrator.
std::string getName() const
Get the name for this integrator.
virtual IntegratorListenerPtr getListener()
Definition: RK45Integrator.h:173
virtual std::string getDescription() const
Gets the description associated with this integrator type.
Definition: RK45Integrator.h:187
virtual std::string getName() const
Gets the name associated with this integrator type.
Definition: RK45Integrator.h:179
virtual std::string getHint() const
Gets the hint associated with this integrator type.
Definition: RK45Integrator.h:195
virtual Integrator * construct(ExecutableModel *model) const
Constructs a new integrator of a given type.
Definition: RK45Integrator.h:203
Definition: Variant.h:75