roadrunner  2.6.0
Fast simulator for SBML models
rrRoadRunner.h
1 #ifndef rrRoadRunnerH
2 #define rrRoadRunnerH
3 
4 #include "rrOSSpecifics.h"
5 #include "rrSelectionRecord.h"
6 #include "rrRoadRunnerOptions.h"
7 #include "sbml/Species.h"
8 #include <thread>
9 
10 #ifdef _MSC_VER
11 #pragma warning(disable: 26812)
12 #pragma warning(disable: 26451)
13 #endif
14 
15 #include "rr-libstruct/lsMatrix.h"
16 
17 #ifdef _MSC_VER
18 #pragma warning(disable: 26812)
19 #pragma warning(disable: 26451)
20 #endif
21 
22 #include <string>
23 #include <vector>
24 #include <list>
25 #include <set>
26 
27 
28 namespace ls {
29  class LibStructural;
30 }
31 
32 namespace rr {
33 
34  class ModelGenerator;
35 
36  class SBMLModelSimulation;
37 
38  class ExecutableModel;
39 
40  class Integrator;
41 
42  class SteadyStateSolver;
43 
44  class SensitivitySolver;
45 
46  template<class IndexType, class DataType>
47  class Matrix3D;
48 
49  static std::mutex rrMtx;
50 
59  class RR_DECLSPEC RoadRunner {
60 
61  public:
62 
66  RoadRunner(unsigned int level = 3, unsigned int version = 2);
67 
79  RoadRunner(const std::string &uriOrSBML, const Dictionary *options = 0);
80 
91  RoadRunner(const std::string &compiler, const std::string &tempDir,
92  const std::string &supportCodeDir);
93 
97  RoadRunner(const RoadRunner &rr);
98 
102  void operator=(const RoadRunner& rr);
103 
107  virtual ~RoadRunner();
108 
112  int getInstanceID();
113 
117  int getInstanceCount();
118 
123  static std::string getParamPromotedSBML(const std::string &sArg);
124 
128  std::string getInfo();
129 
133  double getCurrentTime();
134 
138  class Compiler *getCompiler();
139 
146  void setCompiler(const std::string &compiler);
147 
152  Integrator *getIntegrator();
153 
157  SteadyStateSolver *getSteadyStateSolver();
158 
162  SensitivitySolver *getSensitivitySolver();
163 
167  Integrator *getIntegratorByName(const std::string &name);
168 
172  SteadyStateSolver *getSteadyStateSolverByName(const std::string &name);
173 
177  SensitivitySolver *getSensitivitySolverByName(const std::string &name);
178 
182  Integrator *makeIntegrator(const std::string &name);
183 
187  SteadyStateSolver *makeSteadyStateSolver(const std::string &name);
188 
192  SensitivitySolver *makeSensitivitySolver(const std::string &name);
193 
197  std::vector<std::string> getExistingIntegratorNames();
198 
202  std::vector<std::string> getExistingSteadyStateSolverNames();
203 
207  std::vector<std::string> getExistingSensitivitySolverNames();
208 
212  static std::vector<std::string> getRegisteredIntegratorNames();
213 
217  static std::vector<std::string> getRegisteredSteadyStateSolverNames();
218 
222  static std::vector<std::string> getRegisteredSensitivitySolverNames();
223 
224 
228  void setIntegrator(const std::string &name);
229 
235  void setSteadyStateSolver(const std::string &name);
236 
242  void setSensitivitySolver(const std::string &name);
243 
249  bool integratorExists(const std::string &name);
250 
254  bool steadyStateSolverExists(const std::string &name);
255 
259  bool sensitivitySolverExists(const std::string &name);
260 
264  static void registerSolvers();
265 
266  bool isModelLoaded();
267 
271  std::string getModelName();
272 
276  void setModelName(const std::string& name);
277 
281  std::string getModelId();
282 
286  void setModelId(const std::string& id);
287 
294  bool clearModel();
295 
304  double oneStep(double currentTime, double stepSize, bool reset = true);
305 
313  double internalOneStep(double currentTime, double stepSize, bool reset = true);
314 
369  const ls::DoubleMatrix *simulate(const SimulateOptions *options = 0);
370 
377  const ls::DoubleMatrix *simulate(double start, double stop, int points);
378 
383  const ls::DoubleMatrix *simulate(const std::vector<double> &times);
384 
403  Matrix3D<double, double> timeSeriesSensitivities(
404  double start, double stop, int num,
405  std::vector<std::string> params = std::vector<std::string>(),
406  std::vector<std::string> species = std::vector<std::string>(),
407  int k = 0);
408 
416  void saveState(std::string filename, char opt = 'b');
417 
428  std::stringstream* saveStateS(char opt = 'b');
429 
436  void loadState(const std::string& filename);
437 
447  void loadStateS(std::stringstream* state) ;
448 
453  const ls::DoubleMatrix *getSimulationData() const;
454 
455  void setSimulateOptions(const SimulateOptions &settings);
456 
461  SimulateOptions &getSimulateOptions();
462 
469  RoadRunnerOptions &getOptions();
470 
471 
472  void setOptions(const RoadRunnerOptions &);
473 
482  std::string getSBML(int level = 0, int version = 0);
483 
494  std::string getCurrentSBML(int level = 0, int version = 0);
495 
502  void reset();
503 
518  void reset(int options);
519 
523  void resetSelectionLists();
524 
530  void changeInitialConditions(const std::vector<double> &ic);
531 
535  ExecutableModel *getModel();
536 
548  void load(const std::string &uriOrSBML,
549  const Dictionary *options = 0);
550 
551 
552 
557  rr::SelectionRecord createSelection(const std::string &str);
558 
563  std::vector<rr::SelectionRecord> &getSelections();
564 
569  double getValue(const std::string &sel);
570 
571  double getValue(const SelectionRecord &record);
572 
573 
574  void setSelections(const std::vector<std::string> &selections);
575 
576  void setSelections(const std::vector<rr::SelectionRecord> &selections);
577 
581  std::vector<double> getSelectedValues();
582 
586  void getIds(int types, std::list<std::string> &ids);
587 
596  std::vector<std::string> getIndependentFloatingSpeciesIds();
597 
604  std::vector<std::string> getDependentFloatingSpeciesIds();
605 
610  std::vector<std::string> getFloatingSpeciesConcentrationIds();
611 
617  std::vector<std::string> getFloatingSpeciesInitialConcentrationIds();
618 
622  int getSupportedIdTypes();
623 
624 
630  void setValue(const std::string &id, double value);
631 
632 
638  ls::DoubleMatrix getFloatingSpeciesAmountsNamedArray();
639 
644  ls::DoubleMatrix getFloatingSpeciesConcentrationsNamedArray();
645 
650  ls::DoubleMatrix getBoundarySpeciesAmountsNamedArray();
651 
656  ls::DoubleMatrix getBoundarySpeciesConcentrationsNamedArray();
657 
662  std::vector<double> getIndependentFloatingSpeciesAmountsV();
663 
668  std::vector<double> getDependentFloatingSpeciesAmountsV();
669 
674  std::vector<double> getIndependentFloatingSpeciesConcentrationsV();
675 
680  std::vector<double> getDependentFloatingSpeciesConcentrationsV();
681 
686  ls::DoubleMatrix getIndependentFloatingSpeciesAmountsNamedArray();
687 
692  ls::DoubleMatrix getDependentFloatingSpeciesAmountsNamedArray();
693 
698  ls::DoubleMatrix getIndependentFloatingSpeciesConcentrationsNamedArray();
699 
704  ls::DoubleMatrix getDependentFloatingSpeciesConcentrationsNamedArray();
705 
710  std::vector<double> getRatesOfChange();
711 
716  ls::DoubleMatrix getRatesOfChangeNamedArray();
717 
722  std::vector<double> getIndependentRatesOfChange();
723 
728  ls::DoubleMatrix getIndependentRatesOfChangeNamedArray();
729 
734  std::vector<double> getDependentRatesOfChange();
735 
740  ls::DoubleMatrix getDependentRatesOfChangeNamedArray();
741 
745  ls::DoubleMatrix getFullJacobian();
746 
747  ls::DoubleMatrix getFullReorderedJacobian();
748 
754  ls::DoubleMatrix getReducedJacobian(double h = -1.0);
755 
763  std::vector<std::complex<double>> getFullEigenValues();
764 
772  std::vector<std::complex<double>> getReducedEigenValues();
773 
774 
783  ls::DoubleMatrix getFullEigenValuesNamedArray();
784 
793  ls::DoubleMatrix getReducedEigenValuesNamedArray();
794 
795 
796  ls::DoubleMatrix getLinkMatrix();
797 
804  ls::DoubleMatrix getNrMatrix();
805 
806 
811  ls::DoubleMatrix getKMatrix();
812 
819  ls::DoubleMatrix getReducedStoichiometryMatrix();
820 
825  ls::DoubleMatrix getFullStoichiometryMatrix();
826 
827  ls::DoubleMatrix getExtendedStoichiometryMatrix();
828 
829 
830  ls::DoubleMatrix getL0Matrix();
831 
832 
833  ls::DoubleMatrix getConservationMatrix();
834 
836 
838 
839  ls::DoubleMatrix getUnscaledFluxControlCoefficientMatrix();
840 
841  ls::DoubleMatrix getScaledFluxControlCoefficientMatrix();
842 
843 
848  std::vector<std::string> getEigenValueIds();
849 
854  double getUnscaledParameterElasticity(const std::string &reactionName,
855  const std::string &parameterName);
856 
857 
858  ls::DoubleMatrix getFrequencyResponse(double startFrequency,
859  int numberOfDecades, int numberOfPoints,
860  const std::string &parameterName, const std::string &variableName,
861  bool useDB, bool useHz);
862 
866  void setConservedMoietyAnalysis(bool value);
867 
871  bool getConservedMoietyAnalysis();
872 
873 
877  static std::string getExtendedVersionInfo();
878 
879 
884  double getDiffStepSize() const;
885 
890  void setDiffStepSize(double val);
891 
898  double getSteadyStateThreshold() const;
899 
906  void setSteadyStateThreshold(double val);
907 
913  double getFluxThreshold() const;
914 
919  void setFluxThreshold(double val);
920 
929  double getuCC(const std::string &variableName, const std::string &parameterName);
930 
939  double getCC(const std::string &variableName, const std::string &parameterName);
940 
944  double getuEE(const std::string &reactionName, const std::string &parameterName);
945 
950  double getuEE(const std::string &reactionName, const std::string &parameterName,
951  bool computeSteadystate);
952 
956  double getEE(const std::string &reactionName, const std::string &parameterName);
957 
962  double getEE(const std::string &reactionName, const std::string &parameterName,
963  bool computeSteadyState);
964 
968  ls::DoubleMatrix getUnscaledElasticityMatrix();
969 
973  ls::DoubleMatrix getScaledElasticityMatrix();
974 
978  double getScaledFloatingSpeciesElasticity(const std::string &reactionName,
979  const std::string &speciesName);
980 
986  double getUnscaledSpeciesElasticity(int reactionId, int speciesIndex);
987 
1003  void addSpeciesConcentration(const std::string &sid, const std::string &compartment, double initConcentration,
1004  bool hasOnlySubstanceUnits = false, bool boundaryCondition = false,
1005  const std::string &substanceUnits = "", bool forceRegenerate = true);
1006 
1022  void addSpeciesAmount(const std::string &sid, const std::string &compartment, double initAmount = 0,
1023  bool hasOnlySubstanceUnits = false, bool boundaryCondition = false,
1024  const std::string &substanceUnits = "", bool forceRegenerate = true);
1025 
1026 
1027  /*
1028  * Remove a species from the current model. Note that all reactions related to this species(as reactants,
1029  * products or modifiers will be removed as well.
1030  * @param sid: the ID of the species to be removed
1031  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1032  * after this function call
1033  * default value is true to regenerate model after each call
1034  * of editing function
1035  * to save time for editing for multiple times, one could
1036  * set this flag to true only in the last call of editing
1037  */
1038  void removeSpecies(const std::string &sid, bool forceRegenerate = true);
1039 
1052  void setBoundary(const std::string &sid, bool boundaryCondition, bool forceRegenerate = true);
1053 
1065  void setHasOnlySubstanceUnits(const std::string &sid, bool hasOnlySubstanceUnits, bool forceRegenerate = true);
1066 
1071  bool getHasOnlySubstanceUnits(const std::string &sid);
1072 
1084  void setInitAmount(const std::string &sid, double initAmount, bool forceRegenerate = true);
1085 
1086 
1098  void setInitConcentration(const std::string &sid, double initConcentration, bool forceRegenerate = true);
1099 
1100 
1113  void setConstant(const std::string &sid, bool constant, bool forceRegenerate = true);
1114 
1115 
1116  /*
1117  * Add a reaction to the current model
1118  * @param sbmlRep: the SBML representation (i.e. a reaction tag) describing the reaction to be added
1119  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1120  * after this function call
1121  * default value is true to regenerate model after each call
1122  * of editing function
1123  * to save time for editing for multiple times, one could
1124  * set this flag to true only in the last call of editing
1125  */
1126  void addReaction(const std::string &sbmlRep, bool forceRegenerate = true);
1127 
1128 
1129  /*
1130  * Add a reaction to the current model
1131  * By our default, the reaction is not reversible and not fast.
1132  * @param rid: the ID of reaction to be added
1133  * @param reactants: the list of reactant ID, double value could be inserted before ID as stoichiometry
1134  e.g, [2S1] or [1.5S1]
1135  * @param products: the list of product stoichiometry and ID, double value could be inserted before ID as stoichiometry
1136  e.g, [2S1] or [1.5S1]
1137  * @param kineticLaw: the kinetic formula of reaction to be added
1138  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1139  * after this function call
1140  * default value is true to regenerate model after each call
1141  * of editing function
1142  * to save time for editing for multiple times, one could
1143  * set this flag to true only in the last call of editing
1144  */
1145  void addReaction(const std::string &rid, std::vector<std::string> reactants, std::vector<std::string> products,
1146  const std::string &kineticLaw, bool forceRegenerate = true);
1147 
1159  void removeReaction(const std::string &rid, bool deleteUnusedParameters = false, bool forceRegenerate = true);
1160 
1161  /*
1162  * Set the reversible attribut for an existing reaction in the current model
1163  * @param rid: the ID of reaction to be modified
1164  * @param reversible: the reversible attribute to be set
1165  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1166  * after this function call
1167  * default value is true to regenerate model after each call
1168  * of editing function
1169  * to save time for editing for multiple times, one could
1170  * set this flag to true only in the last call of editing
1171  */
1172  void setReversible(const std::string &rid, bool reversible, bool forceRegenerate = true);
1173 
1174 
1175  /*
1176  * Set the kinetic law for an existing reaction in the current model
1177  * @param rid: the ID of reaction to be modified
1178  * @param kineticLaw: the kinetic formular of reaction
1179  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1180  * after this function call
1181  * default value is true to regenerate model after each call
1182  * of editing function
1183  * to save time for editing for multiple times, one could
1184  * set this flag to true only in the last call of editing
1185  */
1186  void setKineticLaw(const std::string &rid, const std::string &kineticLaw, bool forceRegenerate = true);
1187 
1192  std::string getKineticLaw(const std::string &rid);
1193 
1194 
1206  void addParameter(const std::string &pid, double value, bool forceRegenerate = true);
1207 
1218  void removeParameter(const std::string &pid, bool forceRegenerate = true);
1219 
1220 
1232  void addCompartment(const std::string &cid, double initVolume, bool forceRegenerate = true);
1233 
1244  void removeCompartment(const std::string &cid, bool forceRegenerate = true);
1245 
1246 
1247  /*
1248  * Add an assignment rule to the current model
1249  * @param vid: ID of variable that the new rule assigns formula to
1250  * @param formula: the math formula of assignment rule to be added
1251  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1252  * after this function call
1253  * default value is true to regenerate model after each call
1254  * of editing function
1255  * to save time for editing for multiple times, one could
1256  * set this flag to true only in the last call of editing
1257  */
1258  void addAssignmentRule(const std::string &vid, const std::string &formula, bool forceRegenerate = true);
1259 
1260 
1261  /*
1262  * Add a rate rule to the current model
1263  * @param vid: ID of variable that rules assigns formula to
1264  * @param formula: the math formula of rate rule to be added
1265  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1266  * after this function call
1267  * default value is true to regenerate model after each call
1268  * of editing function
1269  * to save time for editing for multiple times, one could
1270  * set this flag to true only in the last call of editing
1271  */
1272  void addRateRule(const std::string &vid, const std::string &formula, bool forceRegenerate = true);
1273 
1274 
1288  void removeRules(const std::string &vid, bool useInitialValue = false, bool forceRegenerate = true);
1289 
1290 
1291  /*
1292  * Add an initial assignment to an exsiting symbol of the current model
1293  * @param vid: ID of symbol
1294  * @param formula: the math formula of the initial assignment
1295  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1296  * after this function call
1297  * default value is true to regenerate model after each call
1298  * of editing function
1299  * to save time for editing for multiple times, one could
1300  * set this flag to true only in the last call of editing
1301  */
1302  void addInitialAssignment(const std::string &vid, const std::string &formula, bool forceRegenerate = true);
1303 
1304 
1315  void removeInitialAssignment(const std::string &vid, bool forceRegenerate = true, bool errIfNotExist = true);
1316 
1317 
1318  /*
1319  * Add an event to the current model
1320  * @param eid: the ID of the event to be added
1321  * @param useValuesFromTriggerTime: indicate the moment at which the eventÂ’s assignments are to be evaluated
1322  * @param trigger: the math formula of event trigger
1323  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1324  * after this function call
1325  * default value is true to regenerate model after each call
1326  * of editing function
1327  * to save time for editing for multiple times, one could
1328  * set this flag to true only in the last call of editing
1329  */
1330  void addEvent(const std::string &eid, bool useValuesFromTriggerTime, const std::string &trigger,
1331  bool forceRegenerate = true);
1332 
1333  /*
1334  * Add trigger to an existing event in the model
1335  * If the given event already has a trigger object, the given trigger will replace the old trigger of the event
1336  * By default, the peresistent attribute is false and the initial value attribute is false
1337  * @param eid: the ID of the event to add trigger
1338  * @param trigger: the math formula of event trigger
1339  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1340  * after this function call
1341  * default value is true to regenerate model after each call
1342  * of editing function
1343  * to save time for editing for multiple times, one could
1344  * set this flag to true only in the last call of editing
1345  */
1346  void addTrigger(const std::string &eid, const std::string &trigger, bool forceRegenerate = true);
1347 
1348  /*
1349  * Set the persistent attribute of the trigger of given event
1350  * @param eid: the ID of the event of the trigger
1351  * @param persistent: the persistent attribute to be set
1352  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1353  * after this function call
1354  * default value is true to regenerate model after each call
1355  * of editing function
1356  * to save time for editing for multiple times, one could
1357  * set this flag to true only in the last call of editing
1358  */
1359  void setPersistent(const std::string &eid, bool persistent, bool forceRegenerate = true);
1360 
1361  /*
1362  * Set the initial value attribute of the trigger of given event
1363  * @param eid: the ID of the event of the trigger
1364  * @param initValue: the initial value attribute to be set
1365  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1366  * after this function call
1367  * default value is true to regenerate model after each call
1368  * of editing function
1369  * to save time for editing for multiple times, one could
1370  * set this flag to true only in the last call of editing
1371  */
1372  void setTriggerInitialValue(const std::string &eid, bool initValue, bool forceRegenerate = true);
1373 
1374  /*
1375  * Add priority to an existing event in the model
1376  * If the given event already has a priority object, the given priority will replace the old priority of the event
1377  * @param eid: the ID of the event to add priority
1378  * @param priority: the math formula of event priority
1379  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1380  * after this function call
1381  * default value is true to regenerate model after each call
1382  * of editing function
1383  * to save time for editing for multiple times, one could
1384  * set this flag to true only in the last call of editing
1385  */
1386  void addPriority(const std::string &eid, const std::string &priority, bool forceRegenerate = true);
1387 
1388  /*
1389  * Add delay to an existing event in the model
1390  * If the given event already has a delay object, the given delay will replace the old delay of the event
1391  * @param eid: the ID of the event to add priority
1392  * @param delay: the math formula of event delay
1393  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1394  * after this function call
1395  * default value is true to regenerate model after each call
1396  * of editing function
1397  * to save time for editing for multiple times, one could
1398  * set this flag to true only in the last call of editing
1399  */
1400  void addDelay(const std::string &eid, const std::string &delay, bool forceRegenerate = true);
1401 
1402 
1403  /*
1404  * Add an event assignment to an existing event in the current model
1405  * @param eid: the ID of the event to add the assignment
1406  * @param vid: the ID of the variable to assign formula
1407  * @param fomula: the math formula to assign
1408  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1409  * after this function call
1410  * default value is true to regenerate model after each call
1411  * of editing function
1412  * to save time for editing for multiple times, one could
1413  * set this flag to true only in the last call of editing
1414  */
1415  void addEventAssignment(const std::string &eid, const std::string &vid, const std::string &formula,
1416  bool forceRegenerate = true);
1417 
1418 
1430  void removeEventAssignments(const std::string &eid, const std::string &vid, bool forceRegenerate = true);
1431 
1442  void removeEvent(const std::string &eid, bool forceRegenerate = true);
1443 
1447  void validateCurrentSBML();
1448 
1449 
1450 
1452  double mcaSteadyState();
1453 
1465  double steadyState(Dictionary *dict = 0);
1466 
1467 
1471  ls::DoubleMatrix steadyStateNamedArray(const Dictionary *dict = 0);
1472 
1476  std::vector<rr::SelectionRecord> &getSteadyStateSelections();
1477 
1482  void setSteadyStateSelections(const std::vector<std::string> &
1483  steadyStateSelections);
1484 
1489  void setSteadyStateSelections(const std::vector<rr::SelectionRecord> &
1490  steadyStateSelections);
1491 
1497  std::vector<double> getSteadyStateValues();
1498 
1502  std::vector<std::string> getSteadyStateSelectionStrings() const;
1503 
1507  ls::DoubleMatrix getSteadyStateValuesNamedArray();
1508 
1514  void regenerateModel(bool forceRegenerate = true, bool reset = false);
1515 
1516 
1517 
1519 
1524  void setBoundarySpeciesByIndex(const int &index, const double &value);
1525 
1530 
1535 
1540  std::vector<std::string> getGlobalParameterIds();
1541 
1545  std::vector<std::string> getBoundarySpeciesIds();
1546 
1550  std::vector<std::string> getAssignmentRuleIds();
1551 
1555  std::vector<std::string> getRateRuleIds();
1556 
1560  std::vector<std::string> getInitialAssignmentIds();
1561 
1565  std::vector<std::string> getBoundarySpeciesConcentrationIds();
1566 
1570  double getBoundarySpeciesByIndex(const int &index);
1571 
1576  double getGlobalParameterByIndex(const int &index);
1577 
1582  void setGlobalParameterByName(const std::string &param, double value);
1583 
1589  double getGlobalParameterByName(const std::string &param);
1590 
1591 
1592 
1593 
1599  std::vector<double> getConservedMoietyValues();
1600 
1601  std::vector<std::string> getConservedMoietyIds();
1602 
1606  void setSeed(long int seed, bool resetModel = true);
1607 
1611  int64_t getSeed(const std::string &integratorName="");
1612 
1613 
1619  void resetSeed();
1620 
1621 #ifndef SWIG // deprecated methods not SWIG'ed
1622 
1626  int getNumberOfReactions();
1627 
1631  double getReactionRate(const int &index);
1632 
1636  double getRateOfChange(const int &index);
1637 
1641  std::vector<std::string> getRateOfChangeIds();
1642 
1647 
1651  void setCompartmentByIndex(const int &index, const double &value);
1652 
1656  double getCompartmentByIndex(const int &index);
1657 
1661  std::vector<std::string> getCompartmentIds();
1662 
1667 
1671  void setBoundarySpeciesConcentrations(const std::vector<double> &values);
1672 
1676  void setBoundarySpeciesAmounts(const std::vector<double> &values);
1677 
1682 
1686  double getFloatingSpeciesByIndex(int index);
1687 
1691  void setFloatingSpeciesByIndex(int index, double value);
1692 
1696  std::vector<double> getFloatingSpeciesConcentrationsV();
1697 
1701  std::vector<double> getFloatingSpeciesAmountsV();
1702 
1706  std::vector<double> getBoundarySpeciesConcentrationsV();
1707 
1711  std::vector<double> getBoundarySpeciesAmountsV();
1712 
1716  std::vector<double> getFloatingSpeciesInitialConcentrations();
1717 
1721  void setFloatingSpeciesConcentrations(const std::vector<double> &values);
1722 
1726  void setFloatingSpeciesInitialConcentrationByIndex(const int &index,
1727  const double &value);
1728 
1732  void setFloatingSpeciesInitialConcentrations(const std::vector<double> &values);
1733 
1737  std::vector<std::string> getFloatingSpeciesIds();
1738 
1742  std::vector<std::string> getFloatingSpeciesInitialConditionIds();
1743 
1747  size_t getNumberOfGlobalParameters();
1748 
1752  void setGlobalParameterByIndex(const int index, const double value);
1753 
1754 
1759  std::vector<double> getGlobalParameterValues();
1760 
1764  void evalModel();
1765 
1766 
1771  std::vector<double> getReactionRates();
1772 
1778  std::vector<std::string> getReactionIds();
1779 
1788  void setTempDir(const std::string &folder);
1789 
1796  std::string getTempDir();
1797 
1801  static const libsbml::SBase* getElementWithMathematicalMeaning(const libsbml::Model* model, const std::string& id);
1804 #endif // #ifndef SWIG
1805 
1806  private:
1807 
1812  static bool llvmInitialized;
1813  static bool solversRegistered;
1814 
1819  class RoadRunnerImpl* impl;
1820 
1821  const int fileMagicNumber = 0xAD6F52;
1822  const int dataVersionNumber;
1823 
1824  void fixDependentSpeciesValues(int except, double *ref);
1825 
1829  void initLLVM();
1830 
1831  size_t createDefaultSteadyStateSelectionList();
1832 
1833  size_t createDefaultTimeCourseSelectionList();
1834 
1839  void getSelectedValues(ls::DoubleMatrix &results, int nRow,
1840  double currentTime);
1841 
1845  void getSelectedValues(std::vector<double> &results, double currentTime);
1846 
1847  bool populateResult();
1848 
1849 
1850  double getNthSelectedOutput(size_t index, double currentTime);
1851 
1852  bool isParameterUsed(const std::string &sid);
1853 
1854  void getAllVariables(const libsbml::ASTNode *node, std::set<std::string> &ids);
1855 
1857  int getTimeRowIndex();
1858 
1859  enum VariableType {
1860  vtSpecies = 0, vtFlux
1861  };
1862 
1863  double getVariableValue(const VariableType variableType,
1864  const int variableIndex);
1865 
1869  ls::LibStructural *getLibStruct();
1870 
1875  //void updateIntegrator();
1876 
1877  bool createDefaultSelectionLists();
1878 
1883  size_t createTimeCourseSelectionList();
1884 
1885  std::vector<SelectionRecord> getSelectionList();
1886 
1892  void applySimulateOptions();
1893 
1894 
1895  enum JacobianMode {
1896  JACOBIAN_FULL, JACOBIAN_REDUCED
1897  };
1898 
1899  std::vector<std::complex<double> > getEigenValues(JacobianMode mode);
1900 
1901  ls::DoubleMatrix getEigenValuesNamedArray(JacobianMode mode);
1902 
1903  /*
1904  * Check if the id already existed in the model
1905  */
1906  void checkID(const std::string &functionName, const std::string &sid);
1907 
1908  /*
1909  * Parse a std::string with format stoichiometry + sID and return its stoichiometry value and sID
1910  */
1911  void parseSpecies(const std::string &species, double *stoichiometry, char **sid);
1912 
1913  /*
1914  * Remove a variable from the current model.
1915  */
1916  void removeVariable(const std::string &sid);
1917 
1918  /*
1919  * check recursively if a ASTnode or any of its child has the given variable
1920  */
1921  bool hasVariable(const libsbml::ASTNode *node, const std::string &sid);
1922 
1923  /*
1924  * Get the names of all the species involved in a given AST
1925  */
1926  void getSpeciesIdsFromAST(const libsbml::ASTNode *node, std::vector<std::string> &species);
1927 
1928  void getSpeciesIdsFromAST(const libsbml::ASTNode *node, std::vector<std::string> &species,
1929  std::vector<std::string> &speciesNames);
1930 
1931  void saveSelectionVector(std::ostream &, std::vector<SelectionRecord> &);
1932 
1933  void loadSelectionVector(std::istream &, std::vector<SelectionRecord> &);
1934  };
1935 
1936 }
1937 
1938 #endif // rrRoadRunnerH
interface to manipulate 'compiler' settings.
Definition: rrCompiler.h:26
This class is frozen, no new features A dictionary interface that objects can implement.
Definition: Dictionary.h:30
Base class for all code generation systems; allows compiling and evaluating the model.
Definition: rrExecutableModel.h:118
Integrator is an abstract base class that provides an interface to specific integrator class implemen...
Definition: Integrator.h:60
A basic local 3D version of the Matrix class, based on initializer_list.
Definition: Matrix3D.h:17
implemention class, hide all details here.
Definition: rrRoadRunner.cpp:184
The main RoadRunner class.
Definition: rrRoadRunner.h:59
a way to find sbml model elements using the RoadRunner syntax.
Definition: rrSelectionRecord.h:18
generic interface for all SensitivitySolvers
Definition: SensitivitySolver.h:22
This class is frozen, no new features RoadRunner simulation options.
Definition: rrRoadRunnerOptions.h:294
SteadyStateSolver is an abstract base class that provides an interface to specific steady-state solve...
Definition: SteadyStateSolver.h:38
C_DECL_SPEC char *rrcCallConv getModelId(RRHandle handle)
Returns the id of currently loaded SBML model.
C_DECL_SPEC void rrcCallConv setModelName(RRHandle handle, char *name)
Sets the name of currently loaded SBML model.
C_DECL_SPEC void rrcCallConv setModelId(RRHandle handle, char *id)
Sets the id of currently loaded SBML model.
C_DECL_SPEC char *rrcCallConv getModelName(RRHandle handle)
Returns the name of currently loaded SBML model.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getReducedJacobian(RRHandle handle)
Retrieve the reduced Jacobian for the current model.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getFullJacobian(RRHandle handle)
Retrieve the full Jacobian for the current model.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getNrMatrix(RRHandle handle)
Retrieve the reduced stoichiometry matrix for the current model.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getConservationMatrix(RRHandle handle)
Retrieve the conservation matrix for the current model.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getL0Matrix(RRHandle handle)
Returns the L0 Matrix.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getLinkMatrix(RRHandle handle)
Retrieve the Link matrix for the current model.
C_DECL_SPEC bool rrcCallConv setBoundarySpeciesByIndex(RRHandle handle, const unsigned int index, const double value)
Set the concentration for a particular boundary species.
C_DECL_SPEC RRStringArrayPtr rrcCallConv getBoundarySpeciesConcentrationIds(RRHandle handle)
Obtain the list of boundary species concentration Ids.
C_DECL_SPEC RRStringArrayPtr rrcCallConv getBoundarySpeciesIds(RRHandle handle)
Obtain the list of boundary species Ids.
C_DECL_SPEC bool rrcCallConv getBoundarySpeciesByIndex(RRHandle handle, const int index, double *value)
Retrieve the concentration for a particular floating species.
C_DECL_SPEC int rrcCallConv getNumberOfBoundarySpecies(RRHandle handle)
Returns the number of boundary species in the model.
C_DECL_SPEC bool rrcCallConv setBoundarySpeciesConcentrations(RRHandle handle, const RRVectorPtr vec)
Set the boundary species concentration to the vector vec.
C_DECL_SPEC RRStringArrayPtr rrcCallConv getCompartmentIds(RRHandle handle)
Obtain the list of compartment Ids.
C_DECL_SPEC bool rrcCallConv getCompartmentByIndex(RRHandle handle, const int index, double *value)
Retrieve the compartment volume for a particular compartment.
C_DECL_SPEC bool rrcCallConv setCompartmentByIndex(RRHandle handle, const int index, const double value)
Set the volume for a particular compartment.
C_DECL_SPEC int rrcCallConv getNumberOfCompartments(RRHandle handle)
Returns the number of compartments in the model.
C_DECL_SPEC bool rrcCallConv removeEventAssignments(RRHandle handle, const char *eid, const char *vid)
Remove event assignments for given variable from an existing event.
C_DECL_SPEC bool rrcCallConv removeReaction(RRHandle handle, const char *rid)
Remove a reaction from the current model.
C_DECL_SPEC bool rrcCallConv addSpeciesConcentration(RRHandle handle, const char *sid, const char *compartment, double initialConcentration, bool hasOnlySubstanceUnits, bool boundaryCondition)
Add a species to the current model.
C_DECL_SPEC bool rrcCallConv setHasOnlySubstanceUnits(RRHandle handle, const char *sid, bool hasOnlySubstanceUnits)
Set the hasOnlySubstanceUnits attribute for an existing species.
C_DECL_SPEC bool rrcCallConv setReversible(RRHandle handle, const char *rid, bool reversible)
Set the reversible attribut for an existing reaction in the current model.
C_DECL_SPEC bool rrcCallConv addEvent(RRHandle handle, const char *eid, bool useValuesFromTriggerTime, const char *trigger)
Add an event to the current model.
C_DECL_SPEC bool rrcCallConv setInitConcentration(RRHandle handle, const char *sid, double initConcentration)
Set initial concentration for an existing species. Initial amount/concentration set before will be un...
C_DECL_SPEC bool rrcCallConv setPersistent(RRHandle handle, const char *eid, bool persistent)
Set the persistent attribute of the trigger of given event.
C_DECL_SPEC bool rrcCallConv validateCurrentSBML(RRHandle handle)
Validate the current SBML file.
C_DECL_SPEC bool rrcCallConv setKineticLaw(RRHandle handle, const char *rid, const char *kineticLaw)
Set the kinetic law for a existing reaction in the current model.
C_DECL_SPEC bool rrcCallConv setConstant(RRHandle handle, const char *sid, bool constant)
Set the constant attribute for an existing species/ parameter/ compartment.
C_DECL_SPEC bool rrcCallConv addRateRule(RRHandle handle, const char *vid, const char *formula)
Add a rate rule to the current model.
C_DECL_SPEC bool rrcCallConv addPriority(RRHandle handle, const char *eid, const char *priority)
Add priority to an existing event in the model If the given event already has a priority object,...
C_DECL_SPEC bool rrcCallConv addCompartment(RRHandle handle, const char *cid, double initVolume)
Add a compartment to the current model.
C_DECL_SPEC bool rrcCallConv setBoundary(RRHandle handle, const char *sid, bool boundaryCondition)
Set the boundary condition of an existing species.
C_DECL_SPEC bool rrcCallConv setInitAmount(RRHandle handle, const char *sid, double initAmount)
Set initial amount for an existing species. Initial amount/concentration set before will be unset.
C_DECL_SPEC bool rrcCallConv addEventAssignment(RRHandle handle, const char *eid, const char *vid, const char *formula)
Add an event assignment to an existing event in the current model.
C_DECL_SPEC bool rrcCallConv addReaction(RRHandle handle, const char *rid, const char **reactants, int numReactants, const char **products, int numProducts, const char *kineticLaw)
Add a reaction to the current model.
C_DECL_SPEC bool rrcCallConv addTrigger(RRHandle handle, const char *eid, const char *trigger)
Add trigger to an existing event in the model If the given event already has a trigger object,...
C_DECL_SPEC bool rrcCallConv addAssignmentRule(RRHandle handle, const char *vid, const char *formula)
Add an assignment rule to the current model.
C_DECL_SPEC bool rrcCallConv removeParameter(RRHandle handle, const char *pid)
Remove a parameter from the current model.
C_DECL_SPEC bool rrcCallConv removeRules(RRHandle handle, const char *vid)
Remove rules related to given variable from the current model.
C_DECL_SPEC bool rrcCallConv setTriggerInitialValue(RRHandle handle, const char *eid, bool initValue)
Set the initial value attribute of the trigger of given event.
C_DECL_SPEC bool rrcCallConv addParameter(RRHandle handle, const char *pid, double value)
Add a parameter to the current model.
C_DECL_SPEC bool rrcCallConv addDelay(RRHandle handle, const char *eid, const char *delay)
Add delay to an existing event in the model If the given event already has a delay object,...
C_DECL_SPEC char *rrcCallConv getKineticLaw(RRHandle handle, const char *rid)
Get the kinetic law for a existing reaction in the current model.
C_DECL_SPEC bool rrcCallConv addSpeciesAmount(RRHandle handle, const char *sid, const char *compartment, double initialAmount, bool hasOnlySubstanceUnits, bool boundaryCondition)
Add a species to the current model.
C_DECL_SPEC bool rrcCallConv removeSpecies(RRHandle handle, const char *sid)
Remove a species from the current model.
C_DECL_SPEC bool rrcCallConv removeCompartment(RRHandle handle, const char *cid)
Remove a compartment from the current model.
C_DECL_SPEC bool rrcCallConv removeEvent(RRHandle handle, const char *eid)
Remove an event from the current model.
C_DECL_SPEC bool rrcCallConv setFloatingSpeciesInitialConcentrationByIndex(RRHandle handle, int index, double value)
Set the initial concentration for a particular floating species.
C_DECL_SPEC int rrcCallConv getNumberOfDependentSpecies(RRHandle handle)
Returns the number of dependent species in the model.
C_DECL_SPEC RRStringArrayPtr rrcCallConv getDependentFloatingSpeciesIds(RRHandle handle)
Obtain the list of dependent floating species Id.
C_DECL_SPEC bool rrcCallConv setFloatingSpeciesByIndex(RRHandle handle, const int index, const double value)
Set the concentration for a particular floating species.
C_DECL_SPEC bool rrcCallConv getFloatingSpeciesByIndex(RRHandle handle, const int index, double *value)
Retrieve the concentration for a particular floating species.
C_DECL_SPEC int rrcCallConv getNumberOfFloatingSpecies(RRHandle handle)
Returns the number of floating species in the model.
C_DECL_SPEC bool rrcCallConv setFloatingSpeciesConcentrations(RRHandle handle, const RRVectorPtr vec)
Set the floating species concentration to the vector vec.
C_DECL_SPEC RRStringArrayPtr rrcCallConv getFloatingSpeciesIds(RRHandle handle)
Obtain the list of floating species Id.
C_DECL_SPEC RRStringArrayPtr rrcCallConv getIndependentFloatingSpeciesIds(RRHandle handle)
Obtain the list of independent floating species Id.
C_DECL_SPEC RRStringArrayPtr rrcCallConv getFloatingSpeciesConcentrationIds(RRHandle handle)
Obtain the list of floating species concentrations Id.
C_DECL_SPEC int rrcCallConv getNumberOfIndependentSpecies(RRHandle handle)
Returns the number of independent species in the model.
C_DECL_SPEC RRStringArrayPtr rrcCallConv getFloatingSpeciesInitialConditionIds(RRHandle handle)
Get the initial floating species Ids.
C_DECL_SPEC RRStringArrayPtr rrcCallConv getFloatingSpeciesInitialConcentrationIds(RRHandle handle)
Get the initial floating species Ids.
C_DECL_SPEC bool rrcCallConv setFloatingSpeciesInitialConcentrations(RRHandle handle, const RRVectorPtr vec)
Set the initial floating species concentrations.
C_DECL_SPEC RRVectorPtr rrcCallConv getFloatingSpeciesInitialConcentrations(RRHandle handle)
Get the initial floating species concentrations.
C_DECL_SPEC bool rrcCallConv isModelLoaded(RRHandle handle)
check if a model is loaded
C_DECL_SPEC char *rrcCallConv getSBML(RRHandle handle)
Retrieve the SBML model that was last loaded into roadRunner.
C_DECL_SPEC bool rrcCallConv loadState(RRHandle handle, const char *filename)
Reload a road runner instance's state saved by saveState.
C_DECL_SPEC bool rrcCallConv clearModel(RRHandle handle)
Unload current model.
C_DECL_SPEC char *rrcCallConv getCurrentSBML(RRHandle handle)
Retrieve the current state of the model in the form of an SBML string.
C_DECL_SPEC bool rrcCallConv saveState(RRHandle handle, const char *filename)
Save a road runner instance's state to a platform-specific binary file.
C_DECL_SPEC bool rrcCallConv getuEE(RRHandle handle, const char *name, const char *species, double *value)
Retrieve a single unscaled elasticity coefficient.
C_DECL_SPEC bool rrcCallConv setDiffStepSize(RRHandle handle, const double value)
Sets the differential step size used in routines such as getCC.
C_DECL_SPEC bool rrcCallConv getCC(RRHandle handle, const char *variable, const char *parameter, double *value)
Retrieve a single control coefficient.
C_DECL_SPEC bool rrcCallConv getEE(RRHandle handle, const char *name, const char *species, double *value)
Retrieve a single elasticity coefficient.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getUnscaledElasticityMatrix(RRHandle handle)
Retrieve the unscaled elasticity matrix for the current model.
C_DECL_SPEC bool rrcCallConv getDiffStepSize(RRHandle handle, double *value)
Retrieve the differential step size used in routines such as getCC.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getScaledElasticityMatrix(RRHandle handle)
Retrieve the scaled elasticity matrix for the current model.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getUnscaledFluxControlCoefficientMatrix(RRHandle handle)
Retrieve the matrix of unscaled flux control coefficients for the current model.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getScaledFluxControlCoefficientMatrix(RRHandle handle)
Retrieve the matrix of scaled flux control coefficients for the current model.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getUnscaledConcentrationControlCoefficientMatrix(RRHandle handle)
Retrieve the matrix of unscaled concentration control coefficients for the current model.
C_DECL_SPEC bool rrcCallConv getScaledFloatingSpeciesElasticity(RRHandle handle, const char *reactionId, const char *speciesId, double *value)
Retrieve the scaled elasticity matrix for the current model.
C_DECL_SPEC bool rrcCallConv getuCC(RRHandle handle, const char *variable, const char *parameter, double *value)
Retrieve a single unscaled control coefficient.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getScaledConcentrationControlCoefficientMatrix(RRHandle handle)
Retrieve the matrix of scaled concentration control coefficients for the current model.
C_DECL_SPEC RRStringArrayPtr rrcCallConv getGlobalParameterIds(RRHandle handle)
Obtain the list of global parameter Ids.
C_DECL_SPEC bool rrcCallConv setGlobalParameterByIndex(RRHandle handle, const int index, const double value)
Set the value for a particular global parameter.
C_DECL_SPEC int rrcCallConv getNumberOfGlobalParameters(RRHandle handle)
Returns the number of global parameters in the model.
C_DECL_SPEC RRVectorPtr rrcCallConv getGlobalParameterValues(RRHandle handle)
Retrieve the values for all the global parameter values in a vector.
C_DECL_SPEC bool rrcCallConv getGlobalParameterByIndex(RRHandle handle, const int index, double *value)
Retrieve the global parameter value.
C_DECL_SPEC bool rrcCallConv getRateOfChange(RRHandle handle, const int, double *value)
Retrieve the rate of change for a given floating species.
C_DECL_SPEC RRVectorPtr rrcCallConv getRatesOfChange(RRHandle handle)
Retrieve the vector of rates of change as determined by the current state of the model.
C_DECL_SPEC bool rrcCallConv getReactionRate(RRHandle handle, const int index, double *rate)
Retrieve a give reaction rate as indicated by the index parameter.
C_DECL_SPEC RRStringArrayPtr rrcCallConv getReactionIds(RRHandle handle)
Obtain the list of reaction Ids.
C_DECL_SPEC int rrcCallConv getNumberOfReactions(RRHandle handle)
Obtain the number of reactions in the loaded model.
C_DECL_SPEC RRVectorPtr rrcCallConv getReactionRates(RRHandle handle)
Retrieve a vector of reaction rates as determined by the current state of the model.
C_DECL_SPEC bool rrcCallConv reset(RRHandle handle)
Resets all variables of the model to their current initial values. Does not change the parameters.
C_DECL_SPEC bool rrcCallConv oneStep(RRHandle handle, const double currentTime, const double stepSize, double *value)
Carry out a one step integration of the model.
C_DECL_SPEC RRCDataPtr rrcCallConv simulate(RRHandle handle)
Carry out a time-course simulation. setTimeStart, setTimeEnd, setNumPoints, etc are used to set the s...
C_DECL_SPEC bool rrcCallConv evalModel(RRHandle handle)
Evaluate the current model, that it update all assignments and rates of change. Do not carry out an i...
C_DECL_SPEC bool rrcCallConv getValue(RRHandle handle, const char *symbolId, double *value)
Get the value for a given symbol, use getAvailableTimeCourseSymbols(void) for a list of symbols.
C_DECL_SPEC bool rrcCallConv setValue(RRHandle handle, const char *symbolId, const double value)
Set the value for a given symbol, use getAvailableTimeCourseSymbols(void) for a list of symbols.
C_DECL_SPEC bool rrcCallConv steadyState(RRHandle handle, double *value)
Compute the steady state of the current model.
C_DECL_SPEC bool rrcCallConv setCompiler(RRHandle handle, const char *fNameWithPath)
Set the path and filename to the compiler to be used by roadrunner.
C_DECL_SPEC char *rrcCallConv getInfo(RRHandle handle)
Retrieve info about current state of roadrunner, e.g. loaded model, conservationAnalysis etc.
C_DECL_SPEC char *rrcCallConv getCompiler(RRHandle handle)
Get the name of the compiler currently being used by roadrunner.
C_DECL_SPEC char *rrcCallConv getParamPromotedSBML(RRHandle handle, const char *sArg)
Promote any local parameters to global status.
C_DECL_SPEC bool rrcCallConv addInitialAssignment(RRHandle handle, const char *vid, const char *formula, bool forceRegenerate)
Add an initial assignment to an exsiting symbol of the current model.
A set of options that determine how the top level RoadRunner class should behave.
Definition: rrRoadRunnerOptions.h:424