roadrunner  2.6.0
Fast simulator for SBML models
Public Member Functions | Public Attributes | Friends | List of all members
rr::ForwardSensitivitySolver Class Reference

Time based sensivitity solver. More...

#include <ForwardSensitivitySolver.h>

Inheritance diagram for rr::ForwardSensitivitySolver:
rr::TimeSeriesSensitivitySolver rr::SensitivitySolver rr::Solver rr::Registrable

Public Member Functions

 ForwardSensitivitySolver (ExecutableModel *executableModel)
 
 ForwardSensitivitySolver (ExecutableModel *executableModel, std::vector< std::string > whichParameters)
 
double integrate (double tStart, double hstep) override
 integrate the model from t0 to t0 + hstep. More...
 
void create ()
 instantiate the code necessary to use cvodes More...
 
void freeSundialsMemory ()
 free sundials memory associated with sensitivities
 
std::string getName () const override
 Get the name of this solver.
 
std::string getDescription () const override
 Get the description of this solver.
 
std::string getHint () const override
 Get a (user-readable) hint for this solver.
 
Solverconstruct (ExecutableModel *executableModel) const override
 Constructs a new Solver of a given type. More...
 
void resetSettings () override
 resets all settings back to default values
 
void syncWithModel (ExecutableModel *executableModel) override
 Called whenever a new model is loaded to allow integrator to reset internal state. More...
 
std::string toRepr () const override
 Return std::string representation a la Python repr method. More...
 
void loadConfigSettings ()
 
std::string toString () const override
 Return a std::string representation of the solver. More...
 
std::vector< std::string > getGlobalParameterNames ()
 returns a vector of global parameter names extracted from the model
 
std::vector< std::string > getVariableNames ()
 get a vector of variable names in the order that they appear in the model. More...
 
ParameterMap getModelParametersAsMap ()
 get global parameters as an unordered map, strings as keys and parameter values as values
 
std::vector< double > getModelParametersAsVector ()
 return a std::vector<double> of model parameters in the order they appear in the model.
 
void deducePlist ()
 returns the indexes of parameters that user wants sensitivities for, based off of whichParameters.
 
N_Vector getStateVector ()
 retuns pointer to the state vector used by sundials for solving ODE.
 
N_Vector * getSensitivityNVectorPtr ()
 retuns pointer to the state vector used by sundials for storing sensitivity matrix.
 
rr::Matrix< double > getSensitivityMatrix (int k=0) override
 get current values of sensitivities of model variables to parameters. More...
 
Matrix3D< double, double > solveSensitivities (double start, double stop, int num, std::vector< std::string > params=std::vector< std::string >(), std::vector< std::string > species=std::vector< std::string >(), int k=0) override
 simulate a timeseries with sensitivities from start to step with num data points. More...
 
void setValue (const std::string &key, Setting val) override
 
- Public Member Functions inherited from rr::SensitivitySolver
 Solver ()=default
 
 Solver (ExecutableModel *model)
 
- Public Member Functions inherited from rr::Solver
 Solver (ExecutableModel *model)
 
void updateSettings (Dictionary *inputSettings)
 Update settings values. More...
 
std::vector< std::string > getSettings () const
 Get a list of all settings for this solver. More...
 
std::unordered_map< std::string, Setting > & getSettingsMap ()
 get settings for this solver More...
 
virtual Setting getValue (const std::string &key) const
 Get the value of an integrator setting. More...
 
virtual Setting hasValue (const std::string &key) const
 Return true if this setting is supported by the integrator. More...
 
virtual size_t getNumParams () const
 Get the number of parameters. More...
 
virtual std::string getParamName (size_t n) const
 Get the name of the parameter at index n. More...
 
virtual std::string getParamDisplayName (int n) const
 Get the display name of the parameter at index n. More...
 
virtual std::string getParamHint (int n) const
 Get the hint of the parameter at index n. More...
 
virtual std::string getParamDesc (int n) const
 Get the description of the parameter at index n. More...
 
virtual std::string getValueAsString (const std::string &key)
 Wrapper for getValue which converts output to a specific type. More...
 
virtual std::string getSettingsRepr () const
 Get the solver settings as a std::string. More...
 
virtual std::string settingsPyDictRepr () const
 Python dictionary-style std::string representation of settings. More...
 
const std::string & getDisplayName (const std::string &key) const
 Gets the hint associated with a given key. More...
 
const std::string & getHint (const std::string &key) const
 Gets the hint associated with a given key. More...
 
const std::string & getDescription (const std::string &key) const
 Gets the description associated with a given key. More...
 
Setting::TypeId getType (const std::string &key) const
 Gets the type associated with a given key. More...
 
virtual ExecutableModelgetModel () const
 returns the pointer to the ExecutableModel
 
virtual std::string getName () const=0
 Gets the name associated with this Solver type. More...
 
virtual std::string getHint () const=0
 Gets the hint associated with this Solver type. More...
 
virtual std::string getDescription () const=0
 Gets the description associated with this Solver type. More...
 

Public Attributes

std::vector< double > p
 parameters in the model as a member variable which enables passing the underlying data pointer to sundials for finite difference approx More...
 
std::vector< double > pbar
 scaling factors. More...
 
std::vector< int > plist
 which parameters to get sensitivity for as int? More...
 
std::vector< std::string > whichParameters
 which parameters to get sensitivities for, as strings
 
ParameterMap globalParameterMap
 a map containing model parameter names to values
 
int Np = 0
 number of global parameters in the model
 
int Ns = 0
 number of parameters we want to find sensitivities for
 
int numModelVariables = 0
 the number of state variables in the model More...
 
- Public Attributes inherited from rr::Solver
SettingsList sorted_settings
 
SettingsMap settings
 
DisplayNameMap display_names_
 
HintMap hints
 
DescriptionMap descriptions
 

Friends

int FFSDyDtFcn (realtype time, N_Vector cv_y, N_Vector cv_ydot, void *userData)
 
int FFSRootFcn (realtype time, N_Vector y_vector, realtype *gout, void *user_data)
 

Additional Inherited Members

- Public Types inherited from rr::Solver
using SettingsList = std::vector< std::string >
 
using SettingsMap = std::unordered_map< std::string, Setting >
 
using DisplayNameMap = std::unordered_map< std::string, std::string >
 
using HintMap = std::unordered_map< std::string, std::string >
 
using DescriptionMap = std::unordered_map< std::string, std::string >
 
- Protected Member Functions inherited from rr::Solver
void addSetting (const std::string &name, const Setting &val, const std::string &display_name, const std::string &hint, const std::string &description)
 
- Protected Attributes inherited from rr::Solver
ExecutableModelmModel = nullptr
 non-owning pointer to model
 

Detailed Description

Time based sensivitity solver.

Uses CVODEIntegrator to integrate the ExecutableModel and cvodes to compute sensitivity information at each time point. Implements the TimeSeriesSensitivitySolver interface

Member Function Documentation

◆ construct()

Solver * rr::ForwardSensitivitySolver::construct ( ExecutableModel model) const
overridevirtual

Constructs a new Solver of a given type.

Author
JKM, WBC

the caller is responsible for deleting memory associated with the returned Solver*.

Implements rr::Registrable.

◆ create()

void rr::ForwardSensitivitySolver::create ( )

instantiate the code necessary to use cvodes

note the cvDirectDemo.c example from sundials is useful to see the different cvode options in action

Sensitivity setup Only do this if there are model variables and parameters since its possible to have an empty RoadRunner object from which we create a model using the direct interface.

◆ getSensitivityMatrix()

rr::Matrix< double > rr::ForwardSensitivitySolver::getSensitivityMatrix ( int  k = 0)
overridevirtual

get current values of sensitivities of model variables to parameters.

Parameters
kthe kth derivative of the sensitivities.

Implements rr::SensitivitySolver.

◆ getVariableNames()

std::vector< std::string > rr::ForwardSensitivitySolver::getVariableNames ( )

get a vector of variable names in the order that they appear in the model.

variables refers to each equation in the model and are obtained with ExecutableModel::getStateVectorId.

◆ integrate()

double rr::ForwardSensitivitySolver::integrate ( double  t0,
double  hstep 
)
overridevirtual

integrate the model from t0 to t0 + hstep.

Note
this signature is the same as that found in Integrator().

integrating the model should update the sensitivities, which are available from getSensitivities()

Implements rr::TimeSeriesSensitivitySolver.

◆ solveSensitivities()

Matrix3D< double, double > rr::ForwardSensitivitySolver::solveSensitivities ( double  start,
double  stop,
int  num,
std::vector< std::string >  params = std::vector<std::string>(),
std::vector< std::string >  species = std::vector<std::string>(),
int  k = 0 
)
overridevirtual

simulate a timeseries with sensitivities from start to step with num data points.

Matrix3D indexed by time. Each element of the 3D matrix is a Matrix<double> with rows and columns parameters and model variables respectively. The parameter k determines the kth order derivative of the sensitivity information that will be returned

Parameters
startstarting time for time series simulation
stoplast time point for time series simulation
numnumber of data points to simulate. Determines Z of Matrix3D.
paramsvector of parameters that you want sensitivity for. When empty (default), compute sensitivities for all parameters vs all variables.
speciesvector of species that you want sensitivity for. When empty (default), compute sensitivities for parameters vs all variables.
k(default 0) return the kth other derivative of the sensitivity data.
Note
for developers. Pass by value for the params is easier to wrap with swig.

Note this is slightly contrived logic that is necessary for enabling users to supply whichParameters in the constructor and here via the params argument. This was an "accident of design" and it would be better to just have the whichParameters variable set in this method only.

Implements rr::TimeSeriesSensitivitySolver.

◆ syncWithModel()

void rr::ForwardSensitivitySolver::syncWithModel ( ExecutableModel m)
overridevirtual

Called whenever a new model is loaded to allow integrator to reset internal state.

Author
JKM

Implements rr::Solver.

◆ toRepr()

std::string rr::ForwardSensitivitySolver::toRepr ( ) const
overridevirtual

Return std::string representation a la Python repr method.

Author
JKM

Reimplemented from rr::Solver.

◆ toString()

std::string rr::ForwardSensitivitySolver::toString ( ) const
overridevirtual

Return a std::string representation of the solver.

Author
JKM

Reimplemented from rr::Solver.

Member Data Documentation

◆ numModelVariables

int rr::ForwardSensitivitySolver::numModelVariables = 0

the number of state variables in the model

aka the size of the mStateVector

◆ p

std::vector<double> rr::ForwardSensitivitySolver::p

parameters in the model as a member variable which enables passing the underlying data pointer to sundials for finite difference approx

This is the full vector of model global parameters in the order indexed in the model.

See also
plist

◆ pbar

std::vector<double> rr::ForwardSensitivitySolver::pbar

scaling factors.

Set to the values of the parameters for sensitivity. This should be identical to p as acceptible default.

◆ plist

std::vector<int> rr::ForwardSensitivitySolver::plist

which parameters to get sensitivity for as int?

indexes parameters in p and pbar.


The documentation for this class was generated from the following files: