libRoadRunner C++ API  1.3
 All Classes Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
rrRoadRunner.h
1 #ifndef rrRoadRunnerH
2 #define rrRoadRunnerH
3 
4 #include "rrOSSpecifics.h"
5 #include "rr-libstruct/lsMatrix.h"
6 #include "rrSelectionRecord.h"
7 #include "rrRoadRunnerOptions.h"
8 
9 #include <string>
10 #include <vector>
11 #include <list>
12 
13 namespace ls
14 {
15 class LibStructural;
16 }
17 
18 namespace rr
19 {
20 
21 class ModelGenerator;
22 class SBMLModelSimulation;
23 class ExecutableModel;
24 class Integrator;
25 class SteadyStateSolver;
26 
27 /**
28  * The main RoadRunner class.
29  *
30  * The RoadRunner class is responsible for loading and simulating SBML models.
31  *
32  * MemoryManagment: Any pointer returned by a get... method is owned by the
33  * RoadRunner object and does NOT have to be deleted.
34  */
35 class RR_DECLSPEC RoadRunner
36 {
37 
38 public:
39 
40  /**
41  * create an empty RoadRunner object.
42  */
43  RoadRunner();
44 
45  /**
46  * load an sbml document from anywhere.
47  *
48  * same arguments as load.
49  *
50  * If options is not null, then the RoadRunner::computeAndAssignConservationLaws
51  * flag is set to whatever value is specified in the options struct.
52  *
53  * @param uriOrSBML: a URI, local path or sbml document contents.
54  * @param options: an options struct, if null, default values are used.
55  */
56  RoadRunner(const std::string& uriOrSBML, const Dictionary* options = 0);
57 
58  /**
59  * All three of the RoadRunner options default to the empty string, in this
60  * case, the default values are used.
61  *
62  * @param compiler: If LLVM build is enabled, the compiler defaults to LLVM.
63  * @param tempDir: If the old external C compiler is used, this is the
64  * where the C files are written to.
65  * @param supportCodeDir: If the old external C compiler is used, this is
66  * the location where roadrunner C include files are.
67  */
68  RoadRunner(const std::string& compiler, const std::string& tempDir,
69  const std::string& supportCodeDir);
70 
71  /**
72  * free any memory this class allocated
73  */
74  virtual ~RoadRunner();
75 
76  /**
77  * When there are multiple instances of RoadRunner, this is the instance id.
78  */
79  int getInstanceID();
80 
81  /**
82  * Number of currently running RoadRunner instances.
83  */
84  int getInstanceCount();
85 
86  /**
87  * given an sbml document, this method moves all the local parameters
88  * to global parameters.
89  */
90  static std::string getParamPromotedSBML(const std::string& sArg);
91 
92  /**
93  * information about the current state of this object.
94  */
95  std::string getInfo();
96 
97  /**
98  * The Compiler that the ModelGenerator is using to compile / interpret sbml code.
99  */
100  class Compiler* getCompiler();
101 
102  /**
103  * Set the name of the externa compiler to use. Some ModelGenerators may have no use
104  * for this value.
105  *
106  * An exception is raised if the string is invalid.
107  */
108  void setCompiler(const std::string& compiler);
109 
110  /**
111  * get a pointer to the integrator which is currently being used to
112  * time evolve the system.
113  */
114  Integrator* getIntegrator();
115 
116  /**
117  * Get an integrator by name without switching the current integrator
118  */
119  Integrator* getIntegratorByName(const std::string& name);
120 
121  /**
122  * Make an integrator for the given name
123  */
124  Integrator* makeIntegrator(std::string name);
125 
126  /**
127  * get a pointer to the current steady state solver
128  */
129  SteadyStateSolver* getSteadyStateSolver();
130 
131  /* Return a list of the names of all existing integrators. */
132  std::vector<std::string> getExistingIntegratorNames();
133 
134  /**
135  * Return a vector of the names of all registered integrators
136  */
137  static std::vector<std::string> getRegisteredIntegratorNames();
138 
139  /**
140  * Return a vector of the names of all registered integrators
141  */
142  static std::vector<std::string> getRegisteredSteadyStateSolverNames();
143 
144  /**
145  * Ensures all integrators and steady state solvers are registered
146  */
147  static void ensureSolversRegistered();
148 
149  // DEPRECATED
150  //Integrator* getIntegrator(std::string name);
151 
152  void setIntegrator(std::string name);
153 
154  bool integratorExists(std::string name);
155 
156  void setSteadyStateSolver(std::string name);
157 
158  bool steadyStateSolverExists(std::string name);
159 
160  bool isModelLoaded();
161 
162  /**
163  * returns the model name if a model is loaded, empty string otherwise.
164  */
165  std::string getModelName();
166 
167  /**
168  * @brief Clears the currently loaded model and all associated memory
169  * @details Deletes jitted code and libStruct data
170  * @returns True if memory was freed, false if no model was loaded
171  * in the first place
172  */
173  bool clearModel();
174 
175 
176  /**
177  * Carry out a single integration step using a stepsize as indicated
178  * in the method call.
179  *
180  * @param t0 starting time
181  * @param tf final time
182  * @param options override current options.
183  */
184  double integrate(double t0, double tf, const SimulateOptions* options = 0);
185 
186  /**
187  * @deprecated, use integrate instead.
188  */
189  double oneStep(double currentTime, double stepSize, bool reset = true);
190 
191  /**
192  * Simulate the current SBML model.
193  *
194  * If options is null, then the current simulation settings (start time,
195  * end time, n steps) are used. If options is not null, then the
196  * current simulation settings are set to the values specified by
197  * options and they are used.
198  *
199  * The options Dictionary may contain a large number of options, for a complete
200  * list of all available options for each integrator type, @see IntegratorFactory,
201  * @see SimulateOptions.
202  *
203  * For example, to perform a simulation from time 0 to 10 with 1000 steps, using a
204  * stiff integtator, you would:
205  * @code
206  * RoadRunner r = RoadRunner("someFile.xml");
207  * BasicDictionary opt;
208  * opt.setItem("start", 0);
209  * opt.setItem("duration", 10);
210  * opt.setItem("steps", 1000);
211  * opt.setItem("stiff", true);
212  * const DoubleMatrix *result = r.simulate(&opt);
213  * @endcode
214  *
215  * Similarly, if one wants to use a stochastic integrator, such as the Gillespie
216  * integrator, this is set via the "integrator" key, i.e.
217  * @code
218  * RoadRunner r = RoadRunner("someFile.xml");
219  * BasicDictionary opt;
220  * opt.setItem("integrator", "gillespie");
221  * opt.setItem("start", 0);
222  * opt.setItem("duration", 10);
223  * opt.setItem("steps", 1000);
224  * opt.setItem("stiff", true);
225  * opt.setItem("seed", 12345);
226  * const DoubleMatrix *result = r.simulate(&opt);
227  * @endcode
228  * Here, the "integrator" specifies the integrator to use. The "stiff" key
229  * is only used by the deterministic solvers, and it is safely ignored by the
230  * stochastic solvers. Also, the "seed" sets the random seed that the integrator
231  * uses. For more information about all of the avaialble options for each integrator,
232  * @see IntegratorFactory::getIntegratorOptions".
233  *
234  * @throws an std::exception if any options are invalid.
235  * @returns a borrowed reference to a DoubleMatrix object if successfull.
236  */
237  const ls::DoubleMatrix *simulate(const Dictionary* options = 0);
238 
239  /**
240  * RoadRunner keeps a copy of the simulation data around until the
241  * next call to simulate. This matrix can be obtained here.
242  */
243  const ls::DoubleMatrix* getSimulationData() const;
244 
245  #ifndef SWIG // deprecated methods not SWIG'ed
246 
247  #endif
248 
249  void setSimulateOptions(const SimulateOptions& settings);
250 
251  /**
252  * get a reference to the SimulateOptions that were set either
253  * by setSimulateOptions or simulate.
254  */
255  SimulateOptions& getSimulateOptions();
256 
257  /**
258  * Get a reference to the options that determine how this class should behave.
259  *
260  * These are general options. For options specific for loading or simulationg,
261  * @see getSimulateOptions.
262  */
263  RoadRunnerOptions& getOptions();
264 
265 
266  void setOptions(const RoadRunnerOptions&);
267 
268  /**
269  * get the originally loaded sbml document as a string.
270  *
271  * This may optionally up or downconvert the document to a different version, if
272  * the level and version arguments are non-zero.
273  *
274  * If both arguments are zero, then the document is left alone and the
275  */
276  std::string getSBML(int level = 0, int version = 0);
277 
278 
279  /**
280  * Returns the SBML with the current model parameters. This is different than
281  * getSBML which returns the orginal SBML.
282  *
283  * This may optionally up or downconvert the document to a different version, if
284  * the level and version arguments are non-zero.
285  *
286  * If both arguments are zero, then the document is left alone and the original
287  * version is returned.
288  */
289  std::string getCurrentSBML(int level = 0, int version = 0);
290 
291  /**
292  * Picks up default options from config.
293  *
294  * The current default is
295  * SelectionRecord::TIME | SelectionRecord::RATE | SelectionRecord::FLOATING.
296  */
297  void reset();
298 
299  /**
300  * reset the model accordign to a bitfield specified by the
301  * SelectionRecord::SelectionType values.
302  *
303  * Note, this would make more sense as an unsigned, however SWIG has issues mapping
304  * an enum to an unsigned, but seems to map enums just fine to an int.
305  *
306  * For example, to reset the floating species, time and rate rule values:
307  * @code
308  * r.reset(SelectionRecord::TIME | SelectionRecord::RATE | SelectionRecord::FLOATING);
309  * @endcode
310  *
311  * @param options a bitmask made from the SelectionRecord::SelectionTypes values.
312  */
313  void reset(int options);
314 
315  /**
316  * @internal
317  * set the floating species initial concentrations.
318  *
319  * equivalent to ExecutableModel::reset, then ExecutableModel::setFloatingSpeciesConcentrations
320  *
321  * @deprecated
322  */
323  void changeInitialConditions(const std::vector<double>& ic);
324 
325 
326  /**
327  * get a pointer to the ExecutableModel owned by the RoadRunner object.
328  */
329  ExecutableModel* getModel();
330 
331  /**
332  * load an sbml document from anywhere.
333  *
334  * If options is not null, then the RoadRunner::computeAndAssignConservationLaws
335  * flag is set to whatever value is specified in the options struct.
336  *
337  * WARNING, will no longer return a value in next version.
338  *
339  * @param uriOrSBML: a URI, local path or sbml document contents.
340  * @param options: an options struct, if null, default values are used.
341  */
342  void load(const std::string& uriOrSBML,
343  const Dictionary* options = 0);
344 
345 
346 /************************ Selection Ids Species Section ***********************/
347 #if (1) /**********************************************************************/
348 /******************************************************************************/
349 
350  /**
351  * create a selection record. This record can be used to select values.
352  */
353  rr::SelectionRecord createSelection(const std::string& str);
354 
355  /**
356  * Returns the currently selected columns that will be returned by
357  * calls to simulate() or simulateEx(,,).
358  */
359  std::vector<rr::SelectionRecord>& getSelections();
360 
361  /**
362  * Creates a new selection based on the selection string,
363  * and returns the value it queries.
364  */
365  double getValue(const std::string& sel);
366 
367  double getValue(const SelectionRecord& record);
368 
369 
370  void setSelections(const std::vector<std::string>& selections);
371 
372  void setSelections(const std::vector<rr::SelectionRecord>& selections);
373 
374 
375  /**
376  * returns the values selected with SimulateOptions for the current model time / timestep")
377  */
378  std::vector<double> getSelectedValues();
379 
380  /**
381  * populates a given list with all the ids that this class can accept.
382  */
383  void getIds(int types, std::list<std::string> &ids);
384 
385  /**
386  * @author JKM
387  * @brief Gets the ids for all independent floating species
388  * @details Independent means (in this case) that the species is unique
389  * up to a conserved cycle, i.e. each conserved cycle counts for only
390  * one species. See Sauro, Systems Biology: Introduction to Pathway
391  * Modeling, 1st ed. pp. 60.
392  */
393  std::vector<std::string> getIndependentFloatingSpeciesIds();
394 
395  /**
396  * @author JKM
397  * @brief Gets the ids for all dependent floating species
398  * @details See @ref getIndependentFloatingSpeciesIds for an explanation
399  * of independent vs. dependent.
400  */
401  std::vector<std::string> getDependentFloatingSpeciesIds();
402 
403  /**
404  * @brief Gets the ids for all initial floating species concentrations
405  * @details See @ref getIndependentFloatingSpeciesIds for an explanation
406  * of independent vs. dependent.
407  */
408  std::vector<std::string> getFloatingSpeciesInitialConcentrationIds();
409 
410  /**
411  * returns a bit field of the ids that this class supports.
412  */
413  int getSupportedIdTypes();
414 
415 
416  /**
417  * sets the value coresponding to the given selection string
418  *
419  * raises an exception in the selection string is invalid.
420  */
421  void setValue(const std::string& id, double value);
422 
423 /************************ End Selection Ids Species Section *******************/
424 #endif /***********************************************************************/
425 /******************************************************************************/
426 
427  /**
428  * @author JKM
429  * @brief Returns the floating species amounts as a named array
430  */
431  ls::DoubleMatrix getFloatingSpeciesAmounts();
432 
433  /**
434  * @author JKM
435  * @brief Returns the floating species concentrations as a named array
436  */
437  ls::DoubleMatrix getFloatingSpeciesConcentrations();
438 
439  /**
440  * compute the full Jacobian at the current operating point
441  */
442  ls::DoubleMatrix getFullJacobian();
443 
444  ls::DoubleMatrix getFullReorderedJacobian();
445 
446  /**
447  * Compute the reduced Jacobian at the current operating point.
448  * @param h The step sized used for central difference method.
449  * If negative, the default value from the config file is used.
450  */
451  ls::DoubleMatrix getReducedJacobian(double h = -1.0);
452 
453  /**
454  * Returns the eigenvalues of the full jacobian.
455  *
456  * If the eigenvalues are all real, this returns a N x 1 matrix,
457  * if complex, returns an N x 2 matrix where the first column is the
458  * real values and the second is the imaginary part.
459  */
460  std::vector<ls::Complex> getFullEigenValues();
461 
462  /**
463  * Returns the eigenvalues of the reduced jacobian.
464  *
465  * If the eigenvalues are all real, this returns a N x 1 matrix,
466  * if complex, returns an N x 2 matrix where the first column is the
467  * real values and the second is the imaginary part.
468  */
469  std::vector<ls::Complex> getReducedEigenValues();
470 
471 
472  ls::DoubleMatrix getLinkMatrix();
473 
474  /**
475  * get the reduced stochiometry matrix. If conservation conversion is enabled,
476  * this is the matrix that coresponds to the independent species.
477  *
478  * A synonym for getReducedStoichiometryMatrix().
479  */
480  ls::DoubleMatrix getNrMatrix();
481 
482 
483  /**
484  * Returns the K matrix (right nullspace of Nr)
485  * The K matrix has the structure, [I K0]'
486  */
487  ls::DoubleMatrix getKMatrix();
488 
489  /**
490  * get the reduced stochiometry matrix. If conservation conversion is enabled,
491  * this is the matrix that coresponds to the independent species.
492  *
493  * A synonym for getNrMatrix().
494  */
495  ls::DoubleMatrix getReducedStoichiometryMatrix();
496 
497  /**
498  * Get the stoichiometry matrix that coresponds to the full model, even it
499  * it was converted via conservation conversion.
500  */
501  ls::DoubleMatrix getFullStoichiometryMatrix();
502 
503 
504  ls::DoubleMatrix getL0Matrix();
505 
506 
507  ls::DoubleMatrix getConservationMatrix();
508  ls::DoubleMatrix getUnscaledConcentrationControlCoefficientMatrix();
509  ls::DoubleMatrix getScaledConcentrationControlCoefficientMatrix();
510  ls::DoubleMatrix getUnscaledFluxControlCoefficientMatrix();
511  ls::DoubleMatrix getScaledFluxControlCoefficientMatrix();
512 
513 
514  /**
515  * returns the list of floating species, but with a "eigen(...)" string
516  * wrapped around them.
517  */
518  std::vector<std::string> getEigenValueIds();
519 
520  /**
521  * Returns the unscaled elasticity for a named reaction with respect to a
522  * named parameter
523  */
524  double getUnscaledParameterElasticity(const string& reactionName,
525  const string& parameterName);
526 
527 
528  ls::DoubleMatrix getFrequencyResponse(double startFrequency,
529  int numberOfDecades, int numberOfPoints,
530  const string& parameterName, const string& variableName,
531  bool useDB, bool useHz);
532 
533  /**
534  * This method turns on / off the computation and adherence to conservation laws.
535  */
536  void setConservedMoietyAnalysis(bool value);
537 
538  /**
539  * is conservation analysis enabled. This is set
540  */
541  bool getConservedMoietyAnalysis();
542 
543 
544  /**
545  * getVersion plus info about dependent libs versions..
546  */
547  static std::string getExtendedVersionInfo();
548 
549 
550  /**
551  * @author JKM
552  * @brief Set the differential step size used in MCA routines like @ref getCC
553  */
554  double getDiffStepSize() const;
555 
556  /**
557  * @author JKM
558  * @brief Set the differential step size used in MCA routines like @ref getCC
559  */
560  void setDiffStepSize(double val);
561 
562  /**
563  * @author JKM
564  * @brief Get the steady state threshold used in getCC
565  * @details In the MCA routines, RoadRunner will keep trying to
566  * converge to a steady state until this threshold is met
567  */
568  double getSteadyStateThreshold() const;
569 
570  /**
571  * @author JKM
572  * @brief Set the steady state threshold used in getCC
573  * @details In the MCA routines, RoadRunner will keep trying to
574  * converge to a steady state until this threshold is met
575  */
576  void setSteadyStateThreshold(double val);
577 
578  /**
579  * Get unscaled control coefficient with respect to a global parameter
580  *
581  * variableName must be either a reaction or floating species.
582  *
583  * parameterName must be eithe a global parameter, boundary species, or
584  * conserved sum.
585  */
586  double getuCC(const std::string& variableName, const std::string& parameterName);
587 
588  /**
589  * Get scaled control coefficient with respect to a global parameter
590  *
591  * The variableName must be either a reaction id, or a floating species id.
592  *
593  * The parameterName must be either a global parameter, boundary species,
594  * or conserved sum.
595  */
596  double getCC(const std::string& variableName, const std::string& parameterName);
597 
598  /**
599  * Get unscaled elasticity coefficient with respect to a global parameter or species
600  */
601  double getuEE(const std::string& reactionName, const std::string& parameterName);
602 
603  /**
604  * Get unscaled elasticity coefficient with respect to a global parameter or species.
605  * Optionally the model is brought to steady state after the computation.
606  */
607  double getuEE(const std::string& reactionName, const std::string& parameterName,
608  bool computeSteadystate);
609 
610  /**
611  * Get scaled elasticity coefficient with respect to a global parameter or species
612  */
613  double getEE(const std::string& reactionName, const std::string& parameterName);
614 
615  /**
616  * Get scaled elasticity coefficient with respect to a global parameter or species.
617  * Optionally the model is brought to steady state after the computation.
618  */
619  double getEE(const std::string& reactionName, const std::string& parameterName,
620  bool computeSteadyState);
621 
622  /**
623  * Compute the unscaled species elasticity matrix at the current operating point
624  */
625  ls::DoubleMatrix getUnscaledElasticityMatrix();
626 
627  /**
628  * Compute the unscaled elasticity matrix at the current operating point
629  */
630  ls::DoubleMatrix getScaledElasticityMatrix();
631 
632  /**
633  * Compute the scaled elasticity for a given reaction and given species
634  */
635  double getScaledFloatingSpeciesElasticity(const std::string& reactionName,
636  const std::string& speciesName);
637 
638  /**
639  * Get a single species elasticity value
640  * IMPORTANT:
641  * Assumes that the reaction rates have been precomputed at the operating point !!
642  */
643  double getUnscaledSpeciesElasticity(int reactionId, int speciesIndex);
644 
645 
646  /******************************* Steady State Section *************************/
647  #if (1) /**********************************************************************/
648  /******************************************************************************/
649 
650  /**
651  * Compute the steady state of the model, returns the sum of squares of the
652  * solution
653  *
654  * The steady state solver and whatever options it needs may be specified
655  * via the given dictionary. For a list of all available steady state solvers,
656  * @see SteadyStateSolverFactory.
657  *
658  * @param dict a pointer to a dictionary which has the steady state options.
659  * May be NULL, in this case the existing options are used.
660  */
661  double steadyState(const Dictionary* dict = 0);
662 
663  /**
664  * Like @ref steadyState but returns a named array of the steady state values
665  */
666  ls::DoubleMatrix steadyStateNamedArray(const Dictionary* dict = 0);
667 
668  /**
669  * returns the current set of steady state selections.
670  */
671  std::vector<rr::SelectionRecord>& getSteadyStateSelections();
672 
673  /**
674  * parses the given list of strings and generates selections records
675  * which will be used for the steady state selections.
676  */
677  void setSteadyStateSelections(const std::vector<std::string>&
678  steadyStateSelections);
679 
680  /**
681  * makes a copy of an existing list of selection records. These will be
682  * saved and used for selection values in getSteadyStateValues().
683  */
684  void setSteadyStateSelections(const std::vector<rr::SelectionRecord>&
685  steadyStateSelections);
686 
687  /**
688  * Performs a steady state calculation (evolves the system to a steady
689  * state), then calculates and returns the set of values specifed by
690  * the steady state selections.
691  */
692  std::vector<double> getSteadyStateValues();
693 
694  /**
695  * Returns a vector of the steady state selection strings
696  */
697  std::vector<std::string> getSteadyStateSelectionStrings() const;
698 
699  /**
700  * Like @ref getSteadyStateValues but returns a named array
701  */
702  ls::DoubleMatrix getSteadyStateValuesNamedArray();
703 
704  /******************************* End Steady State Section *********************/
705  #endif /***********************************************************************/
706  /******************************************************************************/
707 
708  /********* Used by rrplugins *************************/
709 
710  /**
711  * @internal
712  * @deprecated
713  */
714  void setBoundarySpeciesByIndex(const int& index, const double& value);
715 
716  /**
717  * @internal
718  * @deprecated
719  */
720  int getNumberOfIndependentSpecies();
721 
722  /**
723  * @internal
724  * @deprecated use ExecutableModel::getGlobalParameterIds
725  */
726  std::vector<std::string> getGlobalParameterIds();
727 
728  /**
729  * @internal
730  * @deprecated
731  */
732  std::vector<std::string> getBoundarySpeciesIds();
733 
734  /**
735  * @internal
736  * @deprecated
737  */
738  double getBoundarySpeciesByIndex(const int& index);
739 
740  /**
741  * @internal
742  * @deprecated use ExecutableModel::getGlobalParameterValues
743  */
744  double getGlobalParameterByIndex(const int& index);
745 
746 
747  /******** !!! DEPRECATED INTERNAL METHODS * THESE WILL BE REMOVED!!! **********/
748  #if (1) /**********************************************************************/
749  /******************************************************************************/
750 
751  /**
752  * @author MTK, JKM
753  * @brief Returns the sum of each conserved cycle
754  */
755  std::vector<double> getConservedMoietyValues();
756 
757  #ifndef SWIG // deprecated methods not SWIG'ed
758 
759  /**
760  * @internal
761  * @deprecated
762  */
763  RR_DEPRECATED(int getNumberOfReactions());
764 
765  /**
766  * @internal
767  * @deprecated
768  */
769  RR_DEPRECATED(double getReactionRate(const int& index));
770 
771  /**
772  * @internal
773  * @deprecated
774  */
775  RR_DEPRECATED(double getRateOfChange(const int& index));
776 
777  /**
778  * @internal
779  * @deprecated
780  */
781  RR_DEPRECATED(std::vector<std::string> getRateOfChangeIds());
782 
783  /**
784  * @internal
785  * @deprecated
786  */
787  RR_DEPRECATED(std::vector<std::string> getConservedMoietyIds());
788 
789  /**
790  * @internal
791  * @deprecated
792  */
793  RR_DEPRECATED(int getNumberOfCompartments());
794 
795  /**
796  * @internal
797  * @deprecated
798  */
799  RR_DEPRECATED(void setCompartmentByIndex(const int& index, const double& value));
800 
801  /**
802  * @internal
803  * @deprecated
804  */
805  RR_DEPRECATED(double getCompartmentByIndex(const int& index));
806 
807  /**
808  * @internal
809  * @deprecated
810  */
811  RR_DEPRECATED(std::vector<std::string> getCompartmentIds());
812 
813  /**
814  * @internal
815  * @deprecated
816  */
817  RR_DEPRECATED(int getNumberOfBoundarySpecies());
818 
819  /**
820  * @internal
821  * @deprecated
822  */
823  RR_DEPRECATED(std::vector<double> getBoundarySpeciesConcentrations());
824 
825  /**
826  * @internal
827  * @deprecated
828  */
829  RR_DEPRECATED(void setBoundarySpeciesConcentrations(const std::vector<double>& values));
830 
831  /**
832  * @internal
833  * @deprecated
834  */
835  RR_DEPRECATED(int getNumberOfFloatingSpecies());
836 
837  /**
838  * @internal
839  * @deprecated
840  */
841  RR_DEPRECATED(double getFloatingSpeciesByIndex(int index));
842 
843  /**
844  * @internal
845  * @deprecated
846  */
847  RR_DEPRECATED(void setFloatingSpeciesByIndex(int index, double value));
848 
849  /**
850  * @internal
851  * @deprecated
852  */
853  RR_DEPRECATED(std::vector<double> getFloatingSpeciesConcentrationsV());
854 
855  /**
856  * @internal
857  * @deprecated
858  */
859  RR_DEPRECATED(std::vector<double> getFloatingSpeciesInitialConcentrations());
860 
861  /**
862  * @internal
863  * @deprecated
864  */
865  RR_DEPRECATED(void setFloatingSpeciesConcentrations(const std::vector<double>& values));
866 
867  /**
868  * @internal
869  * @deprecated
870  */
871  RR_DEPRECATED(void setFloatingSpeciesInitialConcentrationByIndex(const int& index,
872  const double& value));
873 
874  /**
875  * @internal
876  * @deprecated
877  */
878  RR_DEPRECATED(void setFloatingSpeciesInitialConcentrations(const std::vector<double>& values));
879 
880  /**
881  * @internal
882  * @deprecated
883  */
884  RR_DEPRECATED(std::vector<std::string> getFloatingSpeciesIds());
885 
886  /**
887  * @internal
888  * @deprecated
889  */
890  RR_DEPRECATED(std::vector<std::string> getFloatingSpeciesInitialConditionIds());
891 
892  /**
893  * @internal
894  * @deprecated use ExecutableModel::getNumGlobalParameters
895  */
896  RR_DEPRECATED(int getNumberOfGlobalParameters());
897 
898  /**
899  * @internal
900  * @deprecated use ExecutableModel::setGlobalParameterValues
901  */
902  RR_DEPRECATED(void setGlobalParameterByIndex(const int index, const double value));
903 
904  /**
905  * @internal
906  * @deprecated use ExecutableModel::getGlobalParameterValues
907  */
908  RR_DEPRECATED(std::vector<double> getGlobalParameterValues());
909 
910  /**
911  * @internal
912  * @deprecated
913  */
914  void evalModel();
915 
916  /**
917  * @internal
918  * @deprecated
919  */
920  RR_DEPRECATED(int getNumberOfDependentSpecies());
921 
922 
923  /**
924  * @internal
925  * @deprecated, use ExecutableModel::getReactionRates
926  */
927  RR_DEPRECATED(std::vector<double> getReactionRates());
928 
929  /**
930  * @internal
931  * @deprecated
932  * returns a list of reaction ids obtained from
933  * ExecutableModel::getReactionId
934  */
935  RR_DEPRECATED(std::vector<std::string> getReactionIds());
936 
937  /**
938  * @internal
939  * @deprecated
940  *
941  * C backend only
942  * set the location where the ModelGenerator creates temporary files, such
943  * as shared libraries.
944  */
945  void setTempDir(const std::string& folder);
946 
947  /**
948  * @internal
949  * @deprecated
950  *
951  * legacy C backend only
952  * get the ModelGenerator's temporary file directory.
953  */
954  std::string getTempDir();
955 
956  #endif // #ifndef SWIG
957 
958 
959  /******** !!! DEPRECATED INTERNAL METHODS * THESE WILL BE REMOVED!!! **********/
960  #endif /**********************************************************************/
961  /******************************************************************************/
962 
963 private:
964 
965 
966  int createDefaultSteadyStateSelectionList();
967  int createDefaultTimeCourseSelectionList();
968 
969  /**
970  * copies the current selection values into the n'th row of the
971  * given matrix
972  */
973  void getSelectedValues(ls::DoubleMatrix& results, int nRow,
974  double currentTime);
975 
976  /**
977  * copies the current selection values into the given vector.
978  */
979  void getSelectedValues(std::vector<double> &results, double currentTime);
980 
981  bool populateResult();
982 
983 
984  double getNthSelectedOutput(unsigned index, double currentTime);
985 
986 
987  /// Get the row index of the time variable in the output array (returns -1 if time is not selected)
988  int getTimeRowIndex();
989 
990  enum VariableType
991  {
992  vtSpecies = 0, vtFlux
993  };
994 
995  double getVariableValue(const VariableType variableType,
996  const int variableIndex);
997 
998  /**
999  * the LibStruct is normally null, only created on demand here.
1000  */
1001  ls::LibStructural* getLibStruct();
1002 
1003  /**
1004  * If the specified integrator does not exist, create it, and point the
1005  * integrator pointer to it.
1006  */
1007  //void updateIntegrator();
1008 
1009  bool createDefaultSelectionLists();
1010 
1011  /**
1012  * creates a selection list from the amounts / conc / variables ivars of the
1013  * SimulationOptions struct.
1014  */
1015  int createTimeCourseSelectionList();
1016 
1017  std::vector<SelectionRecord> getSelectionList();
1018 
1019  /**
1020  * The simulateOptions may be changed at any time. This updates
1021  * the integrators just before they are used with the
1022  * potentially changed options.
1023  */
1024  void applySimulateOptions();
1025 
1026 
1027  enum JacobianMode {
1028  JACOBIAN_FULL, JACOBIAN_REDUCED
1029  };
1030 
1031  std::vector< std::complex<double> > getEigenValues(JacobianMode mode);
1032 
1033  /**
1034  * private implementation class, can only access if inside
1035  * the implementation file.
1036  */
1037  class RoadRunnerImpl* impl;
1038 };
1039 
1040 }
1041 
1042 #endif
Definition: rrRoadRunner.h:13
This class is frozen, no new features RoadRunner simulation options.
Definition: rrRoadRunnerOptions.h:231
Definition: rrCompiler.h:25
Definition: CVODEIntegrator.cpp:26
Definition: rrRoadRunnerOptions.h:332
Base class for all code generation systems; allows compiling and evaluating the model.
Definition: rrExecutableModel.h:127
Definition: rrRoadRunner.h:35
Definition: rrSelectionRecord.h:15
This class is frozen, no new features A dictionary interface that objects can implement. These are automatically mapped to python dictionaries.
Definition: Dictionary.h:30