libRoadRunner C++ API  1.3
 All Classes Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
rrExecutableModel.h
Go to the documentation of this file.
1 
2 // == PREAMBLE ================================================
3 
4 // * Licensed under the Apache License, Version 2.0; see README
5 
6 // == FILEDOC =================================================
7 
8 /** @file rrExecutableModel.h
9 * @author MTK, ETS, JKM
10 * @copyright Apache License, Version 2.0
11 * @brief Base class for all code generators in RoadRunner
12 **/
13 
14 # ifndef rrExecutableModelH
15 # define rrExecutableModelH
16 
17 // == INCLUDES ================================================
18 
19 # include "rrOSSpecifics.h"
20 # include "rrException.h"
21 
22 # include <stdint.h>
23 # include <string>
24 # include <vector>
25 # include <list>
26 # include <ostream>
27 
28 
29 # include "tr1proxy/rr_memory.h"
30 
31 // == CODE ====================================================
32 
33 namespace rr
34 {
35 
36 class ExecutableModel;
37 
38 /**
39  * @brief Notifies the user of SBML events
40  * @details RoadRunner has the capatiblity to notify user objects of any sbml event.
41  *
42  * In order to listen to sbml events, one simply implements the EventHandler
43  * interface and resgisters it with the ExecutableModel::setEventHandler method.
44  *
45  * To remove it, just pass in a 0.
46  *
47  * EventHanders are free to change any model parameters.
48  *
49  * They may return a result value specified by the Result enum. Currently, we
50  * we only have the HALT_SIMULATION which will result in RoadRunner::oneStep
51  * or RoadRunner::simulate to stop at the current time and return. This may be usefull
52  * if someone wants to run a simulation up until some threshold or state is reached.
53  */
55 {
56 public:
57  enum Result
58  {
59  HALT_SIMULATION = (0x1 << 0), // => 0x00000001
60  };
61 
62  virtual uint onTrigger(ExecutableModel* model, int eventIndex, const std::string& eventId) = 0;
63  virtual uint onAssignment(ExecutableModel* model, int eventIndex, const std::string& eventId) = 0;
64 
65 protected:
66  ~EventListener() {};
67 };
68 
69 /**
70  * listeners are shared objects, so use std smart pointers
71  * to manage them.
72  */
73 typedef cxx11_ns::shared_ptr<EventListener> EventListenerPtr;
74 
75 class EventListenerException: public std::exception
76 {
77 public:
78  explicit EventListenerException(uint resultCode) :
79  resultCode(resultCode)
80  {
81  msg = "EventHandlerException, resultCode: ";
82 
83  switch (resultCode)
84  {
85  case EventListener::HALT_SIMULATION:
86  msg += "HALT_SIMULATION";
87  break;
88  }
89  }
90 
91  virtual ~EventListenerException() throw()
92  {
93  };
94 
95  virtual const char* what() const throw()
96  {
97  return msg.c_str();
98  }
99 
100  uint getResultCode() const
101  {
102  return resultCode;
103  }
104 
105 private:
106  uint resultCode;
107  std::string msg;
108 };
109 
110 
111 
112 /**
113  * @brief Base class for all code generation systems; allows compiling
114  * and evaluating the model
115  * @details The ExecutableModel interface provides a way to access an
116  * sbml model that was compiled, JIT'd or interpreted
117  * as executable (runnable) module.
118  *
119  * An ExecutableModel holds a ModelData structure, all the simulation
120  * values are stored in the ModelData struct, i.e. the dynamic state
121  * of the model is fully contained in the ModelData structure.
122  *
123  * An ExecutableModel shoud also contain all of the initial condisions,
124  * rules, functions and whatever other semantic information that was
125  * specified in the sbml model.
126  */
127 class RR_DECLSPEC ExecutableModel
128 {
129 public:
130  /**
131  * @author JKM
132  * @date 07/31/2015
133  * @brief Returns a human-readable description of the code generation backend,
134  * e.g. LLVM, legacy C, etc.
135  */
136  virtual std::string getExecutableModelDesc() const = 0;
137 
138  /**
139  * get the name of the model
140  */
141  virtual std::string getModelName() = 0;
142 
143  virtual void setTime(double _time) = 0;
144  virtual double getTime() = 0;
145 
146 
147  /**
148  * Loads the initial conditions into the current model state.
149  *
150  * Initial conditions may have been updated at any time externally.
151  */
152  virtual void reset() = 0;
153 
154 
155 
156  /************************ Floating Species Section ****************************/
157  #if (1) /**********************************************************************/
158  /******************************************************************************/
159 
160  /**
161  * dependent species are defined by rules and the only way to change them
162  * is by changing the values on which they depend.
163  */
164  virtual int getNumDepFloatingSpecies() = 0;
165 
166  /**
167  * total number of floating species.
168  */
169  virtual int getNumFloatingSpecies() = 0;
170 
171  virtual int getFloatingSpeciesIndex(const std::string& eid) = 0;
172  virtual std::string getFloatingSpeciesId(int index) = 0;
173 
174  /**
175  * independent species do are not defined by rules, they typically participate
176  * in reactions and can have thier values set at any time.
177  */
178  virtual int getNumIndFloatingSpecies() = 0;
179 
180  /**
181  * get the floating species amounts
182  *
183  * @param[in] len the length of the indx and values arrays.
184  * @param[in] indx an array of length len of boundary species to return.
185  * @param[out] values an array of at least length len which will store the
186  * returned boundary species amounts.
187  */
188  virtual int getFloatingSpeciesAmounts(int len, int const *indx,
189  double *values) = 0;
190 
191  virtual int setFloatingSpeciesAmounts(int len, int const *indx,
192  const double *values) = 0;
193 
194  virtual int getFloatingSpeciesAmountRates(int len, int const *indx,
195  double *values) = 0;
196 
197 
198  virtual int getFloatingSpeciesConcentrationRates(int len, int const *indx,
199  double *values) = 0;
200 
201  /**
202  * get the floating species concentrations
203  *
204  * @param[in] len the length of the indx and values arrays.
205  * @param[in] indx an array of length len of boundary species to return.
206  * @param[out] values an array of at least length len which will store the
207  * returned boundary species amounts.
208  */
209  virtual int getFloatingSpeciesConcentrations(int len, int const *indx,
210  double *values) = 0;
211 
212  /**
213  * set the floating species concentrations
214  *
215  * @param[in] len the length of the indx and values arrays.
216  * @param[in] indx an array of length len of boundary species to return.
217  * @param[in] values an array of at least length len which will store the
218  * returned boundary species amounts.
219  */
220  virtual int setFloatingSpeciesConcentrations(int len, int const *indx,
221  double const *values) = 0;
222 
223  /**
224  * Set the initial concentrations of the floating species.
225  *
226  * Takes the same indices as the other floating species methods.
227  *
228  * Note, if a floating species has an initial assignment rule,
229  * than the initial conditions value can only be set by
230  * updating the values on which it depends, it can not be set
231  * directly.
232  */
233  virtual int setFloatingSpeciesInitConcentrations(int len, int const *indx,
234  double const *values) = 0;
235 
236  /**
237  * Get the initial concentrations of the floating species,
238  * uses the same indexing as the other floating species methods.
239  */
240  virtual int getFloatingSpeciesInitConcentrations(int len, int const *indx,
241  double *values) = 0;
242 
243  /**
244  * Set the initial amounts of the floating species.
245  *
246  * Takes the same indices as the other floating species methods.
247  *
248  * Note, if a floating species has an initial assignment rule,
249  * than the initial conditions value can only be set by
250  * updating the values on which it depends, it can not be set
251  * directly.
252  */
253  virtual int setFloatingSpeciesInitAmounts(int len, int const *indx,
254  double const *values) = 0;
255 
256  /**
257  * Get the initial amounts of the floating species,
258  * uses the same indexing as the other floating species methods.
259  */
260  virtual int getFloatingSpeciesInitAmounts(int len, int const *indx,
261  double *values) = 0;
262 
263  /************************ End Floating Species Section ************************/
264  #endif /***********************************************************************/
265  /******************************************************************************/
266 
267 
268 
269  /************************ Boundary Species Section ****************************/
270  #if (1) /**********************************************************************/
271  /******************************************************************************/
272 
273 
274  /**
275  * get the number of boundary species.
276  */
277  virtual int getNumBoundarySpecies() = 0;
278  virtual int getBoundarySpeciesIndex(const std::string &eid) = 0;
279  virtual std::string getBoundarySpeciesId(int index) = 0;
280 
281  /**
282  * get the boundary species amounts
283  *
284  * @param[in] len the length of the indx and values arrays.
285  * @param[in] indx an array of length len of boundary species to return.
286  * @param[out] values an array of at least length len which will store the
287  * returned boundary species amounts.
288  */
289  virtual int getBoundarySpeciesAmounts(int len, int const *indx,
290  double *values) = 0;
291 
292 
293  /**
294  * get the boundary species concentrations
295  *
296  * @param[in] len the length of the indx and values arrays.
297  * @param[in] indx an array of length len of boundary species to return.
298  * @param[out] values an array of at least length len which will store the
299  * returned boundary species amounts.
300  */
301  virtual int getBoundarySpeciesConcentrations(int len, int const *indx,
302  double *values) = 0;
303 
304  /**
305  * set the boundary species concentrations
306  *
307  * @param[in] len the length of the indx and values arrays.
308  * @param[in] indx an array of length len of boundary species to return.
309  * @param[in] values an array of at least length len which will store the
310  * returned boundary species amounts.
311  */
312  virtual int setBoundarySpeciesConcentrations(int len, int const *indx,
313  double const *values) = 0;
314 
315 
316  /************************ End Boundary Species Section ************************/
317  #endif /***********************************************************************/
318  /******************************************************************************/
319 
320 
321  /************************ Global Parameters Section ***************************/
322  #if (1) /**********************************************************************/
323  /******************************************************************************/
324 
325  /**
326  * get the number of global parameters
327  */
328  virtual int getNumGlobalParameters() = 0;
329 
330  /**
331  * index of the global parameter id, -1 if it does not exist.
332  */
333  virtual int getGlobalParameterIndex(const std::string& eid) = 0;
334 
335  /**
336  * id of the indexed global parameter.
337  */
338  virtual std::string getGlobalParameterId(int index) = 0;
339 
340  /**
341  * get the global parameter values
342  *
343  * @param[in] len the length of the indx and values arrays.
344  * @param[in] indx an array of length len of boundary species to return.
345  * @param[out] values an array of at least length len which will store the
346  * returned boundary species amounts.
347  */
348  virtual int getGlobalParameterValues(int len, int const *indx,
349  double *values) = 0;
350 
351  virtual int setGlobalParameterValues(int len, int const *indx,
352  const double *values) = 0;
353 
354 
355  /************************ Global Parameters Species Section *******************/
356  #endif /***********************************************************************/
357  /******************************************************************************/
358 
359 
360  /************************ Compartments Section ********************************/
361  #if (1) /**********************************************************************/
362  /******************************************************************************/
363 
364  virtual int getNumCompartments() = 0;
365  virtual int getCompartmentIndex(const std::string& eid) = 0;
366  virtual std::string getCompartmentId(int index) = 0;
367 
368  /**
369  * get the compartment volumes
370  *
371  * @param[in] len the length of the indx and values arrays.
372  * @param[in] indx an array of length len of boundary species to return.
373  * @param[out] values an array of at least length len which will store the
374  * returned boundary species amounts.
375  */
376  virtual int getCompartmentVolumes(int len, int const *indx,
377  double *values) = 0;
378 
379  virtual int setCompartmentVolumes(int len, int const *indx,
380  const double *values) = 0;
381 
382  /**
383  * Set the initial volumes of the compartments.
384  *
385  * Takes the same indices as the other compartment methods.
386  *
387  * Note, if a compartment has an initial assignment rule,
388  * than the initial conditions value can only be set by
389  * updating the values on which it depends, it can not be set
390  * directly.
391  */
392  virtual int setCompartmentInitVolumes(int len, int const *indx,
393  double const *values) = 0;
394 
395  /**
396  * Get the initial volume of the compartments,
397  * uses the same indexing as the other compartment methods.
398  */
399  virtual int getCompartmentInitVolumes(int len, int const *indx,
400  double *values) = 0;
401 
402 
403  /************************ End Compartments Species Section ********************/
404  #endif /***********************************************************************/
405  /******************************************************************************/
406 
407 
408  /************************ Selection Ids Species Section ***********************/
409  #if (1) /**********************************************************************/
410  /******************************************************************************/
411 
412  /**
413  * populates a given list with all the ids that this class can accept.
414  *
415  * @param ids: a list of strings that will be filled by this class.
416  * @param types: the types of ids that are requested. Can be set to
417  * 0xffffffff to request all the ids that this class supports.
418  * This should by a bitwise OR of the filelds in SelectionRecord::SelectionType
419  */
420  virtual void getIds(int types, std::list<std::string> &ids) = 0;
421 
422  /**
423  * returns a bit field of the ids that this class supports.
424  */
425  virtual int getSupportedIdTypes() = 0;
426 
427  /**
428  * gets the value for the given id string. The string must be a SelectionRecord
429  * string that is accepted by this class.
430  */
431  virtual double getValue(const std::string& id) = 0;
432 
433  /**
434  * sets the value coresponding to the given selection stringl
435  */
436  virtual void setValue(const std::string& id, double value) = 0;
437 
438 
439  /************************ End Selection Ids Species Section *******************/
440  #endif /***********************************************************************/
441  /******************************************************************************/
442 
443  /**
444  * allocate a block of memory and copy the stochiometric values into it,
445  * and return it.
446  *
447  * The caller is responsible for freeing the memory that is referenced by data.
448  *
449  * @param[out] rows will hold the number of rows in the matrix.
450  * @param[out] cols will hold the number of columns in the matrix.
451  * @param[out] data a pointer which will hold a newly allocated memory block.
452  */
453  virtual int getStoichiometryMatrix(int* rows, int* cols, double** data) = 0;
454 
455  /**
456  * Get the current stiochiometry value for the given species / reaction.
457  *
458  * If either are not valid, NaN is returned.
459  */
460  virtual double getStoichiometry(int speciesIndex, int reactionIndex) = 0;
461 
462 
463  virtual int getNumConservedMoieties() = 0;
464  virtual int getConservedMoietyIndex(const std::string& eid) = 0;
465  virtual std::string getConservedMoietyId(int index) = 0;
466  virtual int getConservedMoietyValues(int len, int const *indx, double *values) = 0;
467  virtual int setConservedMoietyValues(int len, int const *indx,
468  const double *values) = 0;
469 
470  virtual int getNumRateRules() = 0;
471 
472  /**
473  * @author JKM
474  * @date 07/31/2015
475  * @brief Gets the symbols defined by rate rules, i.e.
476  * returns all x such that x' = f(x) is a rule which defines parameter x.
477  */
478  virtual std::vector<std::string> getRateRuleSymbols() const {
479  throw NotImplementedException("getRateRuleSymbols not implemented in " + getExecutableModelDesc());
480  }
481 
482  /**
483  * get the number of reactions the model has
484  */
485  virtual int getNumReactions() = 0;
486 
487  /**
488  * get the index of a named reaction
489  * @returns >= 0 on success, < 0 on failure.
490  */
491  virtual int getReactionIndex(const std::string& eid) = 0;
492 
493  /**
494  * get the name of the specified reaction
495  */
496  virtual std::string getReactionId(int index) = 0;
497 
498  /**
499  * get the vector of reaction rates.
500  *
501  * @param len: the length of the suplied buffer, must be >= reaction rates size.
502  * @param indx: pointer to index array. If NULL, then it is ignored and the
503  * reaction rates are copied directly into the suplied buffer.
504  * @param values: pointer to user suplied buffer where rates will be stored.
505  */
506  virtual int getReactionRates(int len, int const *indx,
507  double *values) = 0;
508 
509  /**
510  * get the 'values' i.e. the what the rate rule integrates to, and
511  * store it in the given array.
512  *
513  * The length of rateRuleValues obviously must be the number of
514  * rate rules we have.
515  */
516  virtual void getRateRuleValues(double *rateRuleValues) = 0;
517 
518  /**
519  * get the id of an element of the state vector.
520  */
521  virtual std::string getStateVectorId(int index) = 0;
522 
523  /**
524  * The state vector is a vector of elements that are defined by
525  * differential equations (rate rules) or independent floating species
526  * are defined by reactions.
527  *
528  * To get the ids of the state vector elements, use getStateVectorId.
529  *
530  * copies the internal model state vector into the provided
531  * buffer.
532  *
533  * @param[out] stateVector a buffer to copy the state vector into, if NULL,
534  * return the size required.
535  *
536  * @return the number of items coppied into the provided buffer, if
537  * stateVector is NULL, returns the length of the state vector.
538  */
539  virtual int getStateVector(double *stateVector) = 0;
540 
541  /**
542  * sets the internal model state to the provided packed state vector.
543  *
544  * @param[in] an array which holds the packed state vector, must be
545  * at least the size returned by getStateVector.
546  *
547  * @return the number of items copied from the state vector, negative
548  * on failure.
549  */
550  virtual int setStateVector(const double *stateVector) = 0;
551 
552  /**
553  * the state vector y is the rate rule values and floating species
554  * concentrations concatenated. y is of length numFloatingSpecies + numRateRules.
555  *
556  * The state vector is packed such that the first n raterule elements are the
557  * values of the rate rules, and the last n floatingspecies are the floating
558  * species values.
559  *
560  * @param[in] time current simulator time
561  * @param[in] y state vector, must be either null, or have a size of that
562  * speciefied by getStateVector. If y is null, then the model is
563  * evaluated using its current state. If y is not null, then the
564  * y is considered the state vector.
565  * @param[out] dydt calculated rate of change of the state vector, if null,
566  * it is ignored.
567  */
568  virtual void getStateVectorRate(double time, const double *y, double* dydt=0) = 0;
569 
570  virtual void testConstraints() = 0;
571 
572  virtual std::string getInfo() = 0;
573 
574  virtual void print(std::ostream &stream) = 0;
575 
576  /******************************* Events Section *******************************/
577  #if (1) /**********************************************************************/
578  /******************************************************************************/
579 
580  virtual int getNumEvents() = 0;
581 
582  /**
583  * get the event status, false if the even is not triggered, true if it is.
584  *
585  * The reason this returns an unsigned char instead of a bool array is this
586  * array is typically stuffed into an std::vector, and std::vector<bool> is
587  * well, weird as it's actually implemented as a bitfield, and can not be
588  * used as a C array.
589  *
590  * So, on every modern system I'm aware of, bool is an unsigned char, so
591  * use that data type here.
592  */
593  virtual int getEventTriggers(int len, const int *indx, unsigned char *values) = 0;
594 
595 
596  /**
597  * Itterate through all of the current and pending events and apply them. If any
598  * events trigger a state change which triggers any additional events, these
599  * are applied as well. After this method finishes, all events are processed.
600  *
601  * @param timeEnd: model time when the event occured.
602  * @param previousEventStatus: array of previous event triggered states.
603  * @param initialState (optional): initial state vector, may be NULL, in which
604  * the current state is used.
605  * @param finalState (optional): final state vector, where the final state is
606  * coppied to. May be NULL, in which case, ignored.
607  */
608  virtual int applyEvents(double timeEnd, const unsigned char* previousEventStatus,
609  const double *initialState, double* finalState) = 0;
610 
611 
612  /**
613  * evaluate the event 'roots' -- when events transition form triggered - non-triggered
614  * or triggered to non-triggered state.
615  *
616  * Simplest method is to return 1 for triggered, -1 for not-triggered, so long
617  * as there is a zero crossing.
618  *
619  * @param time[in] current time
620  * @param y[in] the state vector
621  * @param gdot[out] result event roots, this is of length numEvents.
622  */
623  virtual void getEventRoots(double time, const double* y, double* gdot) = 0;
624 
625  virtual double getNextPendingEventTime(bool pop) = 0;
626 
627  virtual int getPendingEventSize() = 0;
628 
629  virtual void resetEvents() = 0;
630 
631  /**
632  * need a virtual destructor as object implementing this interface
633  * can be deleted directly, i.e.
634  * ExecutableModel *p = createModel(...);
635  * delete p;
636  */
637  virtual ~ExecutableModel() {};
638 
639  /******************************* Events Section *******************************/
640  #endif /**********************************************************************/
641  /******************************************************************************/
642 
643 
644  /**
645  * Gets the index for an event id.
646  * If there is no event with this id, returns -1.
647  */
648  virtual int getEventIndex(const std::string& eid) = 0;
649  virtual std::string getEventId(int index) = 0;
650  virtual void setEventListener(int index, EventListenerPtr eventHandler) = 0;
651  virtual EventListenerPtr getEventListener(int index) = 0;
652 
653  /**
654  * Get the amount rate of change for the i'th floating species
655  * given a reaction rates vector.
656  *
657  * TODO: This should be merged with getFloatingSpeciesAmountRates, but that will
658  * break inteface, will do in next point release.
659  *
660  * TODO: If the conversion factor changes in between getting the
661  * reaction rates vector via getReactionRates
662  *
663  * @param index: index of the desired floating speceis rate.
664  * @param reactionRates: pointer to buffer of reaction rates.
665  */
666  virtual double getFloatingSpeciesAmountRate(int index,
667  const double *reactionRates) = 0;
668 
669  /**
670  * reset the model accordign to a bitfield specified by the
671  * SelectionRecord::SelectionType values.
672  */
673  virtual void reset(int options) = 0;
674 
675  /**
676  * set the seed used by the random number generator. This will by definition
677  * reset the RNG.
678  */
679  virtual void setRandomSeed(int64_t) = 0;
680 
681  /**
682  * get the seed used by the RNG.
683  */
684  virtual int64_t getRandomSeed() = 0;
685 
686  /**
687  * Get a uniform random number between 0 and 1 created by the RNG.
688  *
689  * The type of RNG can be specified in the config file with the
690  * RANDOM_GENERATOR_TYPE key.
691  */
692  virtual double getRandom() = 0;
693 
694  /**
695  * Get the current set of flags
696  */
697  virtual uint32_t getFlags() const = 0;
698 
699  /**
700  * Set certain options that determine the state of the ExecutableModel,
701  * these are listed in
702  */
703  virtual void setFlags(uint32_t) = 0;
704 
706  /**
707  * A simulation is currently running. This means that the model
708  * should not have to re-calculate the reaction rate vector
709  * as it was calculated in the previous integration step.
710  */
711  INTEGRATION = (0x1 << 0), // => 0x00000001
712 
713  /**
714  * optimize not-recalculating the reaction rates during selection.
715  */
716  OPTIMIZE_REACTION_RATE_SELECTION = (0x1 << 1), // => 0x00000010
717  };
718 
719  /**
720  * for source compatability
721  */
723 
724  friend class RoadRunner;
725 
726 protected:
727 
728  /**
729  * is integration is currently proceeding.
730  */
731  void setIntegration(bool value) {
732  uint32_t flags = getFlags();
733  if(value) {
734  flags |= INTEGRATION;
735  } else {
736  flags &= ~INTEGRATION;
737  }
738  setFlags(flags);
739  }
740 
741 
742 };
743 
744 
745 
746 
747 
748 /**
749  * dump the model to a stream convenience func
750  */
751 RR_DECLSPEC std::ostream& operator << (std::ostream &stream, ExecutableModel* model);
752 
753 
754 }
755 #endif
Definition: CVODEIntegrator.cpp:26
Notifies the user of SBML events.
Definition: rrExecutableModel.h:54
virtual ~ExecutableModel()
Definition: rrExecutableModel.h:637
Base class for all code generation systems; allows compiling and evaluating the model.
Definition: rrExecutableModel.h:127
ExecutableModelFlags
Definition: rrExecutableModel.h:705
void setIntegration(bool value)
Definition: rrExecutableModel.h:731
virtual std::vector< std::string > getRateRuleSymbols() const
Gets the symbols defined by rate rules, i.e. returns all x such that x' = f(x) is a rule which define...
Definition: rrExecutableModel.h:478
void computeAllRatesOfChange()
Definition: rrExecutableModel.h:722
Definition: rrRoadRunner.h:35