C API Documentation
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 
9 #ifdef _MSC_VER
10 #pragma warning(disable: 26812)
11 #pragma warning(disable: 26451)
12 #endif
13 #include "rr-libstruct/lsMatrix.h"
14 #ifdef _MSC_VER
15 #pragma warning(disable: 26812)
16 #pragma warning(disable: 26451)
17 #endif
18 
19 #include <string>
20 #include <vector>
21 #include <list>
22 #include <set>
23 
24 namespace ls
25 {
26  class LibStructural;
27 }
28 
29 namespace rr
30 {
31 
32  class ModelGenerator;
33  class SBMLModelSimulation;
34  class ExecutableModel;
35  class Integrator;
36  class SteadyStateSolver;
37 
46  class RR_DECLSPEC RoadRunner
47  {
48 
49  public:
50 
54  RoadRunner(unsigned int level = 3, unsigned int version = 2);
55 
67  RoadRunner(const std::string& uriOrSBML, const Dictionary* options = 0);
68 
79  RoadRunner(const std::string& compiler, const std::string& tempDir,
80  const std::string& supportCodeDir);
81 
86  RoadRunner(const RoadRunner& rr);
87 
91  virtual ~RoadRunner();
92 
97 
102 
107  static std::string getParamPromotedSBML(const std::string& sArg);
108 
112  std::string getInfo();
113 
117  class Compiler* getCompiler();
118 
125  void setCompiler(const std::string& compiler);
126 
132 
136  Integrator* getIntegratorByName(const std::string& name);
137 
141  Integrator* makeIntegrator(std::string name);
142 
147 
148  /* Return a list of the names of all existing integrators. */
149  std::vector<std::string> getExistingIntegratorNames();
150 
154  static std::vector<std::string> getRegisteredIntegratorNames();
155 
159  static std::vector<std::string> getRegisteredSteadyStateSolverNames();
160 
164  static void ensureSolversRegistered();
165 
166  // DEPRECATED
167  //Integrator* getIntegrator(std::string name);
168 
169  void setIntegrator(std::string name);
170 
171  bool integratorExists(std::string name);
172 
173  void setSteadyStateSolver(std::string name);
174 
175  bool steadyStateSolverExists(std::string name);
176 
177  bool isModelLoaded();
178 
182  std::string getModelName();
183 
190  bool clearModel();
191 
200  double oneStep(double currentTime, double stepSize, bool reset = true);
201 
209  double internalOneStep(double currentTime, double stepSize, bool reset = true);
210 
266  const ls::DoubleMatrix *simulate(const Dictionary* options = 0);
267 
268  /*
269  * Saves this roadrunner instance to a file so it can be reloaded later
270  * If opt == 'b' (the default value), this function will output a platform-specific
271  * binary file which can be reloaded later
272  * If opt == 'r', this function will output a human readable file which cannot be reloaded later
273  */
274  void saveState(std::string filename, char opt = 'b');
275 
276  /*
277  * Loads a roadrunner instance saved by saveState with the 'b' option
278  */
279  void loadState(std::string filename);
280 
285  const ls::DoubleMatrix* getSimulationData() const;
286 
287 #ifndef SWIG // deprecated methods not SWIG'ed
288 
289 #endif
290 
291  void setSimulateOptions(const SimulateOptions& settings);
292 
298 
306 
307 
308  void setOptions(const RoadRunnerOptions&);
309 
318  std::string getSBML(int level = 0, int version = 0);
319 
320 
331  std::string getCurrentSBML(int level = 0, int version = 0);
332 
339  void reset();
340 
355  void reset(int options);
356 
361 
370  void changeInitialConditions(const std::vector<double>& ic);
371 
372 
377 
389  void load(const std::string& uriOrSBML,
390  const Dictionary* options = 0);
391 
392 
393 /************************ Selection Ids Species Section ***********************/
394 #if (1) /**********************************************************************/
395 /******************************************************************************/
396 
400  rr::SelectionRecord createSelection(const std::string& str);
401 
406  std::vector<rr::SelectionRecord>& getSelections();
407 
412  double getValue(const std::string& sel);
413 
414  double getValue(const SelectionRecord& record);
415 
416 
417  void setSelections(const std::vector<std::string>& selections);
418 
419  void setSelections(const std::vector<rr::SelectionRecord>& selections);
420 
424  std::vector<double> getSelectedValues();
425 
429  void getIds(int types, std::list<std::string> &ids);
430 
439  std::vector<std::string> getIndependentFloatingSpeciesIds();
440 
447  std::vector<std::string> getDependentFloatingSpeciesIds();
448 
453  std::vector<std::string> getFloatingSpeciesConcentrationIds();
454 
460  std::vector<std::string> getFloatingSpeciesInitialConcentrationIds();
461 
466 
467 
473  void setValue(const std::string& id, double value);
474 
475 /************************ End Selection Ids Species Section *******************/
476 #endif /***********************************************************************/
477 /******************************************************************************/
478 
484 
490 
496 
502 
507  std::vector<double> getRatesOfChange();
508 
513  ls::DoubleMatrix getRatesOfChangeNamedArray();
514 
519  std::vector<double> getIndependentRatesOfChange();
520 
526 
531  std::vector<double> getDependentRatesOfChange();
532 
538 
542  ls::DoubleMatrix getFullJacobian();
543 
544  ls::DoubleMatrix getFullReorderedJacobian();
545 
551  ls::DoubleMatrix getReducedJacobian(double h = -1.0);
552 
560  std::vector<ls::Complex> getFullEigenValues();
561 
569  std::vector<ls::Complex> getReducedEigenValues();
570 
571 
572  ls::DoubleMatrix getLinkMatrix();
573 
580  ls::DoubleMatrix getNrMatrix();
581 
582 
587  ls::DoubleMatrix getKMatrix();
588 
595  ls::DoubleMatrix getReducedStoichiometryMatrix();
596 
601  ls::DoubleMatrix getFullStoichiometryMatrix();
602 
603  ls::DoubleMatrix getExtendedStoichiometryMatrix();
604 
605 
606  ls::DoubleMatrix getL0Matrix();
607 
608 
609  ls::DoubleMatrix getConservationMatrix();
612  ls::DoubleMatrix getUnscaledFluxControlCoefficientMatrix();
613  ls::DoubleMatrix getScaledFluxControlCoefficientMatrix();
614 
615 
620  std::vector<std::string> getEigenValueIds();
621 
626  double getUnscaledParameterElasticity(const string& reactionName,
627  const string& parameterName);
628 
629 
630  ls::DoubleMatrix getFrequencyResponse(double startFrequency,
631  int numberOfDecades, int numberOfPoints,
632  const std::string& parameterName, const std::string& variableName,
633  bool useDB, bool useHz);
634 
638  void setConservedMoietyAnalysis(bool value);
639 
644 
645 
649  static std::string getExtendedVersionInfo();
650 
651 
656  double getDiffStepSize() const;
657 
662  void setDiffStepSize(double val);
663 
670  double getSteadyStateThreshold() const;
671 
678  void setSteadyStateThreshold(double val);
679 
688  double getuCC(const std::string& variableName, const std::string& parameterName);
689 
698  double getCC(const std::string& variableName, const std::string& parameterName);
699 
703  double getuEE(const std::string& reactionName, const std::string& parameterName);
704 
709  double getuEE(const std::string& reactionName, const std::string& parameterName,
710  bool computeSteadystate);
711 
715  double getEE(const std::string& reactionName, const std::string& parameterName);
716 
721  double getEE(const std::string& reactionName, const std::string& parameterName,
722  bool computeSteadyState);
723 
727  ls::DoubleMatrix getUnscaledElasticityMatrix();
728 
732  ls::DoubleMatrix getScaledElasticityMatrix();
733 
737  double getScaledFloatingSpeciesElasticity(const std::string& reactionName,
738  const std::string& speciesName);
739 
745  double getUnscaledSpeciesElasticity(int reactionId, int speciesIndex);
746 
760  void addSpecies(const std::string& sid, const std::string& compartment, double initAmount = 0, bool hasOnlySubstanceUnits=false, bool boundaryCondition=false, const std::string& substanceUnits = "", bool forceRegenerate = true);
761 
762 
763  /*
764  * Remove a species from the current model. Note that all reactions related to this species(as reactants,
765  * products or modifiers will be removed as well.
766  * @param sid: the ID of the species to be removed
767  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
768  * after this function call
769  * default value is true to regenerate model after each call
770  * of editing function
771  * to save time for editing for multiple times, one could
772  * set this flag to true only in the last call of editing
773  */
774  void removeSpecies(const std::string& sid, bool forceRegenerate = true);
775 
788  void setBoundary(const std::string& sid, bool boundaryCondition, bool forceRegenerate = true);
789 
801  void setHasOnlySubstanceUnits(const std::string& sid, bool hasOnlySubstanceUnits, bool forceRegenerate = true);
802 
803 
815  void setInitAmount(const std::string& sid, double initAmount, bool forceRegenerate = true);
816 
817 
829  void setInitConcentration(const std::string& sid, double initConcentration, bool forceRegenerate = true);
830 
831 
844  void setConstant(const std::string& sid, bool constant, bool forceRegenerate = true);
845 
846 
847 
848  /*
849  * Add a reaction to the current model
850  * @param sbmlRep: the SBML representation (i.e. a reaction tag) describing the reaction to be added
851  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
852  * after this function call
853  * default value is true to regenerate model after each call
854  * of editing function
855  * to save time for editing for multiple times, one could
856  * set this flag to true only in the last call of editing
857  */
858  void addReaction(const std::string& sbmlRep, bool forceRegenerate = true);
859 
860 
861 
862  /*
863  * Add a reaction to the current model
864  * By our default, the reaction is not reversible and not fast.
865  * @param rid: the ID of reaction to be added
866  * @param reactants: the list of reactant ID, double value could be inserted before ID as stoichiometry
867  e.g, [2S1] or [1.5S1]
868  * @param products: the list of product stoichiometry and ID, double value could be inserted before ID as stoichiometry
869  e.g, [2S1] or [1.5S1]
870  * @param kineticLaw: the kinetic formula of reaction to be added
871  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
872  * after this function call
873  * default value is true to regenerate model after each call
874  * of editing function
875  * to save time for editing for multiple times, one could
876  * set this flag to true only in the last call of editing
877  */
878  void addReaction(const std::string& rid, std::vector<string> reactants, std::vector<string> products, const std::string& kineticLaw, bool forceRegenerate = true);
879 
890  void removeReaction(const std::string& rid, bool deleteUnusedParameters = false, bool forceRegenerate = true);
891 
892  /*
893  * Set the reversible attribut for an existing reaction in the current model
894  * @param rid: the ID of reaction to be modified
895  * @param reversible: the reversible attribute to be set
896  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
897  * after this function call
898  * default value is true to regenerate model after each call
899  * of editing function
900  * to save time for editing for multiple times, one could
901  * set this flag to true only in the last call of editing
902  */
903  void setReversible(const std::string& rid, bool reversible, bool forceRegenerate = true);
904 
905 
906  /*
907  * Set the kinetic law for an existing reaction in the current model
908  * @param rid: the ID of reaction to be modified
909  * @param kineticLaw: the kinetic formular of reaction
910  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
911  * after this function call
912  * default value is true to regenerate model after each call
913  * of editing function
914  * to save time for editing for multiple times, one could
915  * set this flag to true only in the last call of editing
916  */
917  void setKineticLaw(const std::string& rid, const std::string& kineticLaw, bool forceRegenerate = true);
918 
923  std::string getKineticLaw(const std::string& rid);
924 
925 
937  void addParameter(const std::string& pid, double value, bool forceRegenerate = true);
938 
949  void removeParameter(const std::string& pid, bool forceRegenerate = true);
950 
951 
952 
964  void addCompartment(const std::string& cid, double initVolume, bool forceRegenerate = true);
965 
976  void removeCompartment(const std::string& cid, bool forceRegenerate = true);
977 
978 
979  /*
980  * Add an assignment rule to the current model
981  * @param vid: ID of variable that the new rule assigns formula to
982  * @param formula: the math formula of assignment rule to be added
983  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
984  * after this function call
985  * default value is true to regenerate model after each call
986  * of editing function
987  * to save time for editing for multiple times, one could
988  * set this flag to true only in the last call of editing
989  */
990  void addAssignmentRule(const std::string& vid, const std::string& formula, bool forceRegenerate = true);
991 
992 
993  /*
994  * Add a rate rule to the current model
995  * @param vid: ID of variable that rules assigns formula to
996  * @param formula: the math formula of rate rule to be added
997  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
998  * after this function call
999  * default value is true to regenerate model after each call
1000  * of editing function
1001  * to save time for editing for multiple times, one could
1002  * set this flag to true only in the last call of editing
1003  */
1004  void addRateRule(const std::string& vid, const std::string& formula, bool forceRegenerate = true);
1005 
1006 
1020  void removeRules(const std::string& vid, bool useInitialValue = false, bool forceRegenerate = true);
1021 
1022 
1023  /*
1024  * Add an initial assignment to an exsiting symbol of the current model
1025  * @param vid: ID of symbol
1026  * @param formula: the math formula of the initial assignment
1027  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1028  * after this function call
1029  * default value is true to regenerate model after each call
1030  * of editing function
1031  * to save time for editing for multiple times, one could
1032  * set this flag to true only in the last call of editing
1033  */
1034  void addInitialAssignment(const std::string& vid, const std::string& formula, bool forceRegenerate = true);
1035 
1036 
1047  void removeInitialAssignment(const std::string& vid, bool forceRegenerate = true);
1048 
1049 
1050  /*
1051  * Add an event to the current model
1052  * @param eid: the ID of the event to be added
1053  * @param useValuesFromTriggerTime: indicate the moment at which the eventÂ’s assignments are to be evaluated
1054  * @param trigger: the math formula of event trigger
1055  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1056  * after this function call
1057  * default value is true to regenerate model after each call
1058  * of editing function
1059  * to save time for editing for multiple times, one could
1060  * set this flag to true only in the last call of editing
1061  */
1062  void addEvent(const std::string& eid, bool useValuesFromTriggerTime, const std::string& trigger, bool forceRegenerate = true);
1063 
1064  /*
1065  * Add trigger to an existing event in the model
1066  * If the given event already has a trigger object, the given trigger will replace the old trigger of the event
1067  * By default, the peresistent attribute is false and the initial value attribute is false
1068  * @param eid: the ID of the event to add trigger
1069  * @param trigger: the math formula of event trigger
1070  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1071  * after this function call
1072  * default value is true to regenerate model after each call
1073  * of editing function
1074  * to save time for editing for multiple times, one could
1075  * set this flag to true only in the last call of editing
1076  */
1077  void addTrigger(const std::string& eid, const std::string& trigger, bool forceRegenerate = true);
1078 
1079  /*
1080  * Set the persistent attribute of the trigger of given event
1081  * @param eid: the ID of the event of the trigger
1082  * @param persistent: the persistent attribute to be set
1083  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1084  * after this function call
1085  * default value is true to regenerate model after each call
1086  * of editing function
1087  * to save time for editing for multiple times, one could
1088  * set this flag to true only in the last call of editing
1089  */
1090  void setPersistent(const std::string& eid, bool persistent, bool forceRegenerate = true);
1091 
1092  /*
1093  * Set the initial value attribute of the trigger of given event
1094  * @param eid: the ID of the event of the trigger
1095  * @param initValue: the initial value attribute to be set
1096  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1097  * after this function call
1098  * default value is true to regenerate model after each call
1099  * of editing function
1100  * to save time for editing for multiple times, one could
1101  * set this flag to true only in the last call of editing
1102  */
1103  void setTriggerInitialValue(const std::string& eid, bool initValue, bool forceRegenerate = true);
1104 
1105  /*
1106  * Add priority to an existing event in the model
1107  * If the given event already has a priority object, the given priority will replace the old priority of the event
1108  * @param eid: the ID of the event to add priority
1109  * @param priority: the math formula of event priority
1110  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1111  * after this function call
1112  * default value is true to regenerate model after each call
1113  * of editing function
1114  * to save time for editing for multiple times, one could
1115  * set this flag to true only in the last call of editing
1116  */
1117  void addPriority(const std::string& eid, const std::string& priority, bool forceRegenerate = true);
1118 
1119  /*
1120  * Add delay to an existing event in the model
1121  * If the given event already has a delay object, the given delay will replace the old delay of the event
1122  * @param eid: the ID of the event to add priority
1123  * @param delay: the math formula of event delay
1124  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1125  * after this function call
1126  * default value is true to regenerate model after each call
1127  * of editing function
1128  * to save time for editing for multiple times, one could
1129  * set this flag to true only in the last call of editing
1130  */
1131  void addDelay(const std::string& eid, const std::string& delay, bool forceRegenerate = true);
1132 
1133 
1134 
1135  /*
1136  * Add an event assignment to an existing event in the current model
1137  * @param eid: the ID of the event to add the assignment
1138  * @param vid: the ID of the variable to assign formula
1139  * @param fomula: the math formula to assign
1140  * @param forceRegenerate: a boolean value to indicate if the model is regenerated
1141  * after this function call
1142  * default value is true to regenerate model after each call
1143  * of editing function
1144  * to save time for editing for multiple times, one could
1145  * set this flag to true only in the last call of editing
1146  */
1147  void addEventAssignment(const std::string& eid, const std::string& vid, const std::string& formula, bool forceRegenerate = true);
1148 
1149 
1161  void removeEventAssignments(const std::string& eid, const std::string& vid, bool forceRegenerate = true);
1162 
1173  void removeEvent(const std::string& eid, bool forceRegenerate = true);
1174 
1179 
1180 
1181  /******************************* Steady State Section *************************/
1182 #if (1) /**********************************************************************/
1183  /******************************************************************************/
1184 
1185  double mcaSteadyState();
1186 
1198  double steadyState(const Dictionary* dict = 0);
1199 
1209  double steadyStateApproximate(const Dictionary* dict = 0);
1210 
1214  ls::DoubleMatrix steadyStateNamedArray(const Dictionary* dict = 0);
1215 
1219  std::vector<rr::SelectionRecord>& getSteadyStateSelections();
1220 
1225  void setSteadyStateSelections(const std::vector<std::string>&
1226  steadyStateSelections);
1227 
1232  void setSteadyStateSelections(const std::vector<rr::SelectionRecord>&
1233  steadyStateSelections);
1234 
1240  std::vector<double> getSteadyStateValues();
1241 
1245  std::vector<std::string> getSteadyStateSelectionStrings() const;
1246 
1251 
1255  void regenerate(bool forceRegenerate = true, bool reset = false);
1256 
1257  /******************************* End Steady State Section *********************/
1258 #endif /***********************************************************************/
1259  /******************************************************************************/
1260 
1261  /********* Used by rrplugins *************************/
1262 
1267  void setBoundarySpeciesByIndex(const int& index, const double& value);
1268 
1274 
1279  std::vector<std::string> getGlobalParameterIds();
1280 
1285  std::vector<std::string> getBoundarySpeciesIds();
1286 
1287 
1292  std::vector<std::string> getBoundarySpeciesConcentrationIds();
1293 
1298  double getBoundarySpeciesByIndex(const int& index);
1299 
1304  double getGlobalParameterByIndex(const int& index);
1305 
1310  void setGlobalParameterByName(const std::string& param, double value);
1311 
1317  double getGlobalParameterByName(const std::string& param);
1318 
1319 
1320 
1321  /******** !!! DEPRECATED INTERNAL METHODS * THESE WILL BE REMOVED!!! **********/
1322 #if (1) /**********************************************************************/
1323  /******************************************************************************/
1324 
1329  std::vector<double> getConservedMoietyValues();
1330 
1331  std::vector<std::string> getConservedMoietyIds();
1332 
1333 #ifndef SWIG // deprecated methods not SWIG'ed
1334 
1339  RR_DEPRECATED(int getNumberOfReactions());
1340 
1345  RR_DEPRECATED(double getReactionRate(const int& index));
1346 
1351  RR_DEPRECATED(double getRateOfChange(const int& index));
1352 
1357  RR_DEPRECATED(std::vector<std::string> getRateOfChangeIds());
1358 
1368  RR_DEPRECATED(int getNumberOfCompartments());
1369 
1374  RR_DEPRECATED(void setCompartmentByIndex(const int& index, const double& value));
1375 
1380  RR_DEPRECATED(double getCompartmentByIndex(const int& index));
1381 
1386  RR_DEPRECATED(std::vector<std::string> getCompartmentIds());
1387 
1392  RR_DEPRECATED(int getNumberOfBoundarySpecies());
1393 
1398  RR_DEPRECATED(std::vector<double> getBoundarySpeciesConcentrations());
1399 
1404  RR_DEPRECATED(void setBoundarySpeciesConcentrations(const std::vector<double>& values));
1405 
1410  RR_DEPRECATED(int getNumberOfFloatingSpecies());
1411 
1416  RR_DEPRECATED(double getFloatingSpeciesByIndex(int index));
1417 
1422  RR_DEPRECATED(void setFloatingSpeciesByIndex(int index, double value));
1423 
1428  RR_DEPRECATED(std::vector<double> getFloatingSpeciesConcentrationsV());
1429 
1434  RR_DEPRECATED(std::vector<double> getFloatingSpeciesAmountsV());
1435 
1440  RR_DEPRECATED(std::vector<double> getBoundarySpeciesConcentrationsV());
1441 
1446  RR_DEPRECATED(std::vector<double> getBoundarySpeciesAmountsV());
1447 
1452  RR_DEPRECATED(std::vector<double> getFloatingSpeciesInitialConcentrations());
1453 
1458  RR_DEPRECATED(void setFloatingSpeciesConcentrations(const std::vector<double>& values));
1459 
1464  RR_DEPRECATED(void setFloatingSpeciesInitialConcentrationByIndex(const int& index,
1465  const double& value));
1466 
1471  RR_DEPRECATED(void setFloatingSpeciesInitialConcentrations(const std::vector<double>& values));
1472 
1477  RR_DEPRECATED(std::vector<std::string> getFloatingSpeciesIds());
1478 
1483  RR_DEPRECATED(std::vector<std::string> getFloatingSpeciesInitialConditionIds());
1484 
1489  RR_DEPRECATED(size_t getNumberOfGlobalParameters());
1490 
1495  RR_DEPRECATED(void setGlobalParameterByIndex(const int index, const double value));
1496 
1497 
1502  RR_DEPRECATED(std::vector<double> getGlobalParameterValues());
1503 
1508  void evalModel();
1509 
1514  RR_DEPRECATED(int getNumberOfDependentSpecies());
1515 
1516 
1521  RR_DEPRECATED(std::vector<double> getReactionRates());
1522 
1529  RR_DEPRECATED(std::vector<std::string> getReactionIds());
1530 
1539  void setTempDir(const std::string& folder);
1540 
1548  std::string getTempDir();
1549 
1550 #endif // #ifndef SWIG
1551 
1552 
1553  /******** !!! DEPRECATED INTERNAL METHODS * THESE WILL BE REMOVED!!! **********/
1554 #endif /**********************************************************************/
1555  /******************************************************************************/
1556 
1557  private:
1558 
1559  void fixDependentSpeciesValues(int except, double* ref);
1560 
1561 
1562  size_t createDefaultSteadyStateSelectionList();
1563  size_t createDefaultTimeCourseSelectionList();
1564 
1569  void getSelectedValues(ls::DoubleMatrix& results, int nRow,
1570  double currentTime);
1571 
1575  void getSelectedValues(std::vector<double> &results, double currentTime);
1576 
1577  bool populateResult();
1578 
1579 
1580  double getNthSelectedOutput(size_t index, double currentTime);
1581 
1582  bool isParameterUsed(const std::string& sid);
1583 
1584  void getAllVariables(const libsbml::ASTNode* node, std::set<std::string>& ids);
1585 
1587  int getTimeRowIndex();
1588 
1589  enum VariableType
1590  {
1591  vtSpecies = 0, vtFlux
1592  };
1593 
1594  double getVariableValue(const VariableType variableType,
1595  const int variableIndex);
1596 
1600  ls::LibStructural* getLibStruct();
1601 
1606  //void updateIntegrator();
1607 
1608  bool createDefaultSelectionLists();
1609 
1614  size_t createTimeCourseSelectionList();
1615 
1616  std::vector<SelectionRecord> getSelectionList();
1617 
1623  void applySimulateOptions();
1624 
1625 
1626  enum JacobianMode {
1627  JACOBIAN_FULL, JACOBIAN_REDUCED
1628  };
1629 
1630  std::vector< std::complex<double> > getEigenValues(JacobianMode mode);
1631 
1636  class RoadRunnerImpl* impl;
1637 
1638  /*
1639  * Check if the id already existed in the model
1640  */
1641  void checkID(const std::string& functionName, const std::string& sid);
1642 
1643  /*
1644  * Parse a string with format stoichiometry + sID and return its stoichiometry value and sID
1645  */
1646  void parseSpecies(const string& species, double* stoichiometry, char** sid);
1647 
1648  /*
1649  * Remove a variable from the current model.
1650  */
1651  void removeVariable(const string& sid);
1652 
1653  /*
1654  * check recursively if a ASTnode or any of its child has the given variable
1655  */
1656  bool hasVariable(const libsbml::ASTNode* node, const string& sid);
1657 
1658  /*
1659  * Get the names of all the species involved in a given AST
1660  */
1661  void getSpeciesIdsFromAST(const libsbml::ASTNode* node, std::vector<string>& species);
1662  void getSpeciesIdsFromAST(const libsbml::ASTNode* node, std::vector<string>& species, std::vector<string>& speciesNames);
1663 
1664  /*
1665  * check and remove all parameter without any assignments
1666  */
1667  void checkGlobalParameters();
1668  void saveSelectionVector(std::ostream&, std::vector<SelectionRecord>&);
1669  void loadSelectionVector(std::istream&, std::vector<SelectionRecord>&);
1670  const int fileMagicNumber = 0xAD6F52;
1671  const int dataVersionNumber = 1;
1672  };
1673 
1674 }
1675 
1676 #endif
This class is frozen, no new features A dictionary interface that objects can implement....
Definition: Dictionary.h:31
Base class for all code generation systems; allows compiling and evaluating the model.
Definition: rrExecutableModel.h:118
Definition: Integrator.h:62
Definition: rrRoadRunner.h:47
ls::DoubleMatrix getBoundarySpeciesAmountsNamedArray()
Returns the boundary species amounts as a named array.
std::vector< double > getIndependentRatesOfChange()
Returns the rate of change of the independent floating species as an array.
RoadRunner(unsigned int level=3, unsigned int version=2)
virtual ~RoadRunner()
static void ensureSolversRegistered()
SimulateOptions & getSimulateOptions()
void removeReaction(const std::string &rid, bool deleteUnusedParameters=false, bool forceRegenerate=true)
void setSteadyStateSelections(const std::vector< std::string > &steadyStateSelections)
const ls::DoubleMatrix * simulate(const Dictionary *options=0)
void setGlobalParameterByName(const std::string &param, double value)
set the parameter with id
const ls::DoubleMatrix * getSimulationData() const
ls::DoubleMatrix getRatesOfChangeNamedArray()
Returns the rate of change of the floating species as a named array.
ls::DoubleMatrix getKMatrix()
std::vector< double > getConservedMoietyValues()
Returns the sum of each conserved cycle.
double oneStep(double currentTime, double stepSize, bool reset=true)
void setValue(const std::string &id, double value)
double internalOneStep(double currentTime, double stepSize, bool reset=true)
RoadRunnerOptions & getOptions()
void setBoundary(const std::string &sid, bool boundaryCondition, bool forceRegenerate=true)
std::string getInfo()
SteadyStateSolver * getSteadyStateSolver()
static std::string getParamPromotedSBML(const std::string &sArg)
std::vector< std::string > getIndependentFloatingSpeciesIds()
Gets the ids for all independent floating species.
double getCC(const std::string &variableName, const std::string &parameterName)
std::vector< rr::SelectionRecord > & getSteadyStateSelections()
std::vector< double > getRatesOfChange()
Returns the rate of change of the floating species as an array.
std::string getKineticLaw(const std::string &rid)
RoadRunner(const std::string &compiler, const std::string &tempDir, const std::string &supportCodeDir)
class Compiler * getCompiler()
ls::DoubleMatrix getFullStoichiometryMatrix()
rr::SelectionRecord createSelection(const std::string &str)
ls::DoubleMatrix getFloatingSpeciesConcentrationsNamedArray()
Returns the floating species concentrations as a named array.
int getInstanceCount()
void removeInitialAssignment(const std::string &vid, bool forceRegenerate=true)
void validateCurrentSBML()
double getEE(const std::string &reactionName, const std::string &parameterName)
void setSteadyStateSelections(const std::vector< rr::SelectionRecord > &steadyStateSelections)
ls::DoubleMatrix getScaledElasticityMatrix()
std::vector< std::string > getFloatingSpeciesConcentrationIds()
Gets the ids for all floating species concentrations.
ls::DoubleMatrix getFloatingSpeciesAmountsNamedArray()
Returns the floating species amounts as a named array.
ls::DoubleMatrix getReducedStoichiometryMatrix()
std::vector< double > getSteadyStateValues()
double getEE(const std::string &reactionName, const std::string &parameterName, bool computeSteadyState)
double getUnscaledSpeciesElasticity(int reactionId, int speciesIndex)
Integrator * getIntegratorByName(const std::string &name)
void resetSelectionLists()
std::vector< std::string > getFloatingSpeciesInitialConcentrationIds()
Gets the ids for all initial floating species concentrations.
ls::DoubleMatrix getReducedJacobian(double h=-1.0)
bool clearModel()
Clears the currently loaded model and all associated memory.
void setSteadyStateThreshold(double val)
Set the steady state threshold used in getCC.
int getSupportedIdTypes()
std::vector< std::string > getDependentFloatingSpeciesIds()
Gets the ids for all dependent floating species.
void addParameter(const std::string &pid, double value, bool forceRegenerate=true)
double getuEE(const std::string &reactionName, const std::string &parameterName, bool computeSteadystate)
double getScaledFloatingSpeciesElasticity(const std::string &reactionName, const std::string &speciesName)
std::vector< rr::SelectionRecord > & getSelections()
double steadyState(const Dictionary *dict=0)
void removeEventAssignments(const std::string &eid, const std::string &vid, bool forceRegenerate=true)
std::vector< double > getDependentRatesOfChange()
Returns the rate of change of the dependent floating species as an array.
static std::vector< std::string > getRegisteredIntegratorNames()
std::vector< std::string > getBoundarySpeciesConcentrationIds()
Gets the ids for all boundary species concentrations.
double getuCC(const std::string &variableName, const std::string &parameterName)
void addCompartment(const std::string &cid, double initVolume, bool forceRegenerate=true)
void reset(int options)
void setCompiler(const std::string &compiler)
ls::DoubleMatrix getFullJacobian()
ls::DoubleMatrix getIndependentRatesOfChangeNamedArray()
Returns the rate of change of the independent floating species as a named array.
ls::DoubleMatrix getNrMatrix()
void getIds(int types, std::list< std::string > &ids)
std::vector< std::string > getSteadyStateSelectionStrings() const
void setInitAmount(const std::string &sid, double initAmount, bool forceRegenerate=true)
std::vector< ls::Complex > getReducedEigenValues()
void setInitConcentration(const std::string &sid, double initConcentration, bool forceRegenerate=true)
void regenerateModel(bool forceRegenerate=true, bool reset=false)
void load(const std::string &uriOrSBML, const Dictionary *options=0)
double getValue(const std::string &sel)
void removeParameter(const std::string &pid, bool forceRegenerate=true)
std::vector< std::string > getEigenValueIds()
void setConstant(const std::string &sid, bool constant, bool forceRegenerate=true)
RoadRunner(const RoadRunner &rr)
Integrator * getIntegrator()
std::vector< ls::Complex > getFullEigenValues()
ls::DoubleMatrix getSteadyStateValuesNamedArray()
ExecutableModel * getModel()
ls::DoubleMatrix getBoundarySpeciesConcentrationsNamedArray()
Returns the boundary species concentrations as a named array.
double steadyStateApproximate(const Dictionary *dict=0)
std::string getCurrentSBML(int level=0, int version=0)
std::string getModelName()
static std::vector< std::string > getRegisteredSteadyStateSolverNames()
void addSpecies(const std::string &sid, const std::string &compartment, double initAmount=0, bool hasOnlySubstanceUnits=false, bool boundaryCondition=false, const std::string &substanceUnits="", bool forceRegenerate=true)
double getUnscaledParameterElasticity(const string &reactionName, const string &parameterName)
double getSteadyStateThreshold() const
Get the steady state threshold used in getCC.
static std::string getExtendedVersionInfo()
double getuEE(const std::string &reactionName, const std::string &parameterName)
Integrator * makeIntegrator(std::string name)
void removeCompartment(const std::string &cid, bool forceRegenerate=true)
void removeEvent(const std::string &eid, bool forceRegenerate=true)
RoadRunner(const std::string &uriOrSBML, const Dictionary *options=0)
void setConservedMoietyAnalysis(bool value)
void removeRules(const std::string &vid, bool useInitialValue=false, bool forceRegenerate=true)
double getDiffStepSize() const
Set the differential step size used in MCA routines like getCC.
std::string getSBML(int level=0, int version=0)
void setDiffStepSize(double val)
Set the differential step size used in MCA routines like getCC.
ls::DoubleMatrix getUnscaledElasticityMatrix()
ls::DoubleMatrix getDependentRatesOfChangeNamedArray()
Returns the rate of change of the dependent floating species as a named array.
void setHasOnlySubstanceUnits(const std::string &sid, bool hasOnlySubstanceUnits, bool forceRegenerate=true)
std::vector< double > getSelectedValues()
double getGlobalParameterByName(const std::string &param)
get the
ls::DoubleMatrix steadyStateNamedArray(const Dictionary *dict=0)
bool getConservedMoietyAnalysis()
Definition: rrSelectionRecord.h:16
This class is frozen, no new features RoadRunner simulation options.
Definition: rrRoadRunnerOptions.h:242
Definition: SteadyStateSolver.h:41
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 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 RRVectorPtr rrcCallConv getBoundarySpeciesConcentrations(RRHandle handle)
Retrieve the concentrations for all the boundary species in a vector.
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 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 setPersistent(RRHandle handle, const char *eid, bool persistent)
Set the persistent attribute of the trigger of given event.
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 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 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 setTriggerInitialValue(RRHandle handle, const char *eid, bool initValue)
Set the initial value attribute of the trigger of given event.
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 bool rrcCallConv removeSpecies(RRHandle handle, const char *sid)
Remove a species 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 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 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 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 bool rrcCallConv loadState(RRHandle handle, const char *filename)
Reload a road runner instance's state saved by saveState.
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 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 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 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 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.
Definition: rrRoadRunnerOptions.h:352