roadrunner  2.6.0
Fast simulator for SBML models
SensitivitySolver.h
1 //
2 // Created by Ciaran on 08/07/2021.
3 //
4 
5 #ifndef ROADRUNNER_SENSITIVITYSOLVER_H
6 #define ROADRUNNER_SENSITIVITYSOLVER_H
7 
8 #include "Solver.h"
9 #include "Matrix3D.h"
10 
11 
12 namespace rr {
13 
14  template <typename T>
15  class Matrix;
16 
17  class ForwardSensitivitySolver;
18 
22  class SensitivitySolver : public Solver {
23  public:
24  using Solver::Solver;
25 
26  virtual ~SensitivitySolver() = default;
27 
33  virtual rr::Matrix<double> getSensitivityMatrix(int k = 0) = 0;
34 
35  };
36 
43  public:
44  using SensitivitySolver::SensitivitySolver;
45 
46  ~TimeSeriesSensitivitySolver() override = default;
47 
55  virtual double integrate(double t0, double hstep) = 0;
56 
76  double start, double stop, int num,
77  std::vector<std::string> params,
78  std::vector<std::string> species,
79  int k) = 0;
80 
81  };
82 
88  public:
89  using SensitivitySolver::SensitivitySolver;
90 
91  ~SteadyStateSensitivitySolver() override = default;
92 
98  virtual double solveSteadyState() = 0;
99 
112  std::vector<std::string> params = std::vector<std::string>(),
113  int k = 0) = 0;
114 
115  };
116 
117 }
118 
119 #endif //ROADRUNNER_SENSITIVITYSOLVER_H
Contains the base class for RoadRunner solvers.
A basic local 3D version of the Matrix class, based on initializer_list.
Definition: Matrix3D.h:17
A basic local matrix class, based on the libstruct version.
Definition: Matrix.h:18
generic interface for all SensitivitySolvers
Definition: SensitivitySolver.h:22
virtual rr::Matrix< double > getSensitivityMatrix(int k=0)=0
get current values of sensitivities of model variables to parameters.
Base class for all integrators and steady state solvers.
Definition: Solver.h:39
genetic interface for sensitivity solvers that solve for steady state before computing model sensitiv...
Definition: SensitivitySolver.h:87
virtual double solveSteadyState()=0
solves the model for steady state.
virtual rr::Matrix< double > solveSensitivities(std::vector< std::string > params=std::vector< std::string >(), int k=0)=0
compute sensitivities at steady state
generic interface for sensitivity solvers that integrate the model and compute sensitivities at each ...
Definition: SensitivitySolver.h:42
virtual Matrix3D< double, double > solveSensitivities(double start, double stop, int num, std::vector< std::string > params, std::vector< std::string > species, int k)=0
simulate a timeseries with sensitivities from start to step with num data points.
virtual double integrate(double t0, double hstep)=0
integrate the model from t0 to t0 + hstep.