roadrunner  2.6.0
Fast simulator for SBML models
rrCompiledExecutableModel.h
1 #ifndef rrCompiledExecutableModel_H_
2 #define rrCompiledExecutableModel_H_
3 #include <list>
4 #include "rrExecutableModel.h"
5 #include "rrModelSharedLibrary.h"
6 #include "rrCModelGenerator.h"
7 #include "rrModelData.h"
8 #include "rrNOMSupport.h"
9 #include "rrPendingAssignment.h"
10 #include "rrEvent.h"
11 #include "rr-libstruct/lsLibStructural.h"
12 #include <stack>
13 
14 namespace rr
15 {
16 using namespace ls;
17 using Poco::SharedLibrary;
18 class CGenerator;
19 class CompiledModelState;
20 
21 
22 //Function pointer typedefs..
23 typedef void (RR_CDECL *c_void_MDS)(ModelData*); //MDS stands for ModelDataStructure
24 typedef int (RR_CDECL *c_int_MDS)(ModelData*);
25 typedef int (RR_CDECL *c_int_MDS_int)(ModelData*, int);
26 typedef char* (RR_CDECL *c_charStar_MDS)(ModelData*);
27 typedef void (RR_CDECL *c_void_MDS_doubleStar)(ModelData*, const double*);
28 typedef double (RR_CDECL *c_double_MDS_int)(ModelData*, int);
29 typedef double* (RR_CDECL *c_doubleStar_MDS)(ModelData*);
30 typedef void (RR_CDECL *c_void_MDS_double_doubleStar)(ModelData*, double, const double*);
31 typedef void (RR_CDECL *c_void_MDS_int_double)(ModelData*, int, double);
32 
33 typedef ComputeEventAssignmentHandler* (RR_CDECL *c_ComputeEventAssignmentHandlerStar)();
34 typedef EventDelayHandler* (RR_CDECL *c_GetEventDelayHandlerStar)();
35 
42 class RR_DECLSPEC CompiledExecutableModel: public ExecutableModel
43 {
44 
45 public:
46 
53  virtual std::string getExecutableModelDesc() const {
54  return "Compiled (legacy C) Executable Model";
55  }
56 
62  virtual ~CompiledExecutableModel();
63 
64  virtual bool getConservedSumChanged();
65  virtual void setConservedSumChanged(bool);
66 
67  virtual std::string getModelName();
68 
69  virtual void setTime(double _time);
70  virtual double getTime();
71 
72  virtual void evalInitialConditions();
73 
77  virtual void reset();
78 
83  virtual void reset(int options);
84 
85 
87  {
92  PopDiscard = 0x00000001
93  };
94 
95 
104  virtual int pushState(unsigned options = 0);
105 
112  virtual int popState(unsigned options = 0);
113 
114  // functions --------------------------------------------------------
115  virtual int getNumIndFloatingSpecies();
116  virtual int getNumDepFloatingSpecies();
117 
118  virtual int getNumFloatingSpecies();
119  virtual int getFloatingSpeciesIndex(const std::string& name);
120  virtual std::string getFloatingSpeciesId(int index);
121 
122  virtual int getNumBoundarySpecies();
123  virtual int getBoundarySpeciesIndex(const std::string &name);
124  virtual std::string getBoundarySpeciesId(int index);
125  virtual int getBoundarySpeciesCompartmentIndex(int index);
126 
135  virtual int getFloatingSpeciesAmounts(int len, int const *indx,
136  double *values);
137 
138  virtual int setFloatingSpeciesAmounts(int len, int const *indx,
139  const double *values);
140 
149  virtual int getFloatingSpeciesConcentrations(int len, int const *indx,
150  double *values);
151 
160  virtual int setFloatingSpeciesConcentrations(int len, int const *indx,
161  double const *values);
162 
163  virtual int getFloatingSpeciesAmountRates(int len, int const *indx,
164  double *values);
165 
166  virtual int getFloatingSpeciesConcentrationRates(int len, int const *indx,
167  double *values);
168 
169 
178  virtual int getBoundarySpeciesAmounts(int len, int const *indx,
179  double *values);
180 
189  virtual int getBoundarySpeciesConcentrations(int len, int const *indx,
190  double *values);
191 
200  virtual int setBoundarySpeciesConcentrations(int len, int const *indx,
201  double const *values);
202 
203  virtual int getNumGlobalParameters();
204  virtual int getGlobalParameterIndex(const std::string& name);
205  virtual std::string getGlobalParameterId(int index);
206 
215  virtual int getGlobalParameterValues(int len, int const *indx,
216  double *values);
217 
218  virtual int setGlobalParameterValues(int len, int const *indx,
219  const double *values);
220 
221  virtual int getNumCompartments();
222  virtual int getCompartmentIndex(const std::string& name);
223  virtual std::string getCompartmentId(int index);
224 
233  virtual int getCompartmentVolumes(int len, int const *indx,
234  double *values);
235 
236  virtual int getNumRateRules();
237 
241  virtual int getNumReactions();
242 
243  virtual int getReactionRates(int len, const int* indx, double* values);
244 
249  virtual int getReactionIndex(const std::string& reactionName);
250 
254  virtual std::string getReactionId(int index);
255 
256  virtual int getNumEvents();
257  virtual void computeEventPriorites();
258  void setConcentration(int index, double value);
259  virtual void evalReactionRates();
260  virtual void setCompartmentVolumes();
261  virtual int getNumLocalParameters(int reactionId);
262 
263  virtual void computeRules();
264 
265  virtual void initializeInitialConditions();
266  virtual void setParameterValues();
267  virtual void setBoundaryConditions();
268  virtual void setInitialConditions();
269  virtual void evalInitialAssignments();
270 
271  virtual void convertToAmounts();
272  virtual void computeConservedTotals();
273 
274 
275  //Access dll data
276  virtual void getRateRuleValues(double *rateRuleValues);
277 
278  double getAmount(const int i);
279  virtual void initializeRates();
280 
281  virtual std::string getStateVectorId(int index);
282 
283  virtual void setRateRuleValues(const double *rateRuleValues);
284 
295  virtual int getStateVector(double *stateVector);
296 
306  virtual int setStateVector(const double *stateVector);
307 
308  virtual void convertToConcentrations();
309  virtual void updateDependentSpeciesValues();
310  virtual void computeAllRatesOfChange();
311 
325  virtual void getStateVectorRate(double time, const double *y, double *dydt = 0);
326 
327  virtual void evalEvents(const double time, const double *y);
328  virtual void resetEvents();
329  virtual void testConstraints();
330  virtual void initializeRateRuleSymbols();
331  virtual std::string getInfo();
332 
333  virtual void print(std::ostream &stream);
334 
335  virtual void getIds(int types, std::list<std::string> &ids);
336 
337  virtual int getSupportedIdTypes();
338 
339  virtual double getValue(const std::string& id);
340 
341  virtual void setValue(const std::string& id, double value);
342 
343  virtual int getEventDelays(int len, int const *indx, double *values);
344 
345  virtual int getEventPriorities(int len, int const *indx, double *values);
346 
347  virtual void eventAssignment(int eventId);
348 
349  virtual double* evalEventAssignment(int eventId);
350 
351  virtual void applyEventAssignment(int eventId, double *values);
352 
353  virtual int getEventTriggers(int len, const int *indx, unsigned char *values);
354 
355  virtual int getNumConservedMoieties();
356  virtual int getConservedMoietyIndex(const std::string& name);
357  virtual std::string getConservedMoietyId(int index);
358  virtual int getConservedMoietyValues(int len, int const *indx, double *values);
359  virtual int setConservedMoietyValues(int len, int const *indx,
360  const double *values);
361 
362  virtual int setCompartmentVolumes(int len, int const *indx,
363  const double *values);
364  virtual int setFloatingSpeciesInitConcentrations(int len, int const *indx,
365  double const *values);
366  virtual int getFloatingSpeciesInitConcentrations(int len, int const *indx,
367  double *values);
368 
369  virtual double getStoichiometry(int speciesIndex, int reactionIndex);
370 
381  virtual int getStoichiometryMatrix(int* rows, int* cols, double** data);
382 
383 
384  double getFloatingSpeciesConcentration(int index);
385 
386  bool mConservedSumChanged;
387 
391  bool setupModelData();
392 
397  virtual bool setupDLLFunctions();
398 
399  //This structure holds data generated/used in the shared model lib..
400  //some of it could be made global in the dll later on, like modelName..
401  int mDummyInt;
402  int mDummyDouble;
403  double* mDummyDoubleArray;
404 
410  ModelSymbols ms;
411 
416  ModelSharedLibrary* mDLL;
417 
418  std::stack<CompiledModelState*> modelStates;
419 
420  //Function pointers...
421  c_int_MDS cInitModel;
422  c_int_MDS cInitModelData;
423  c_charStar_MDS cgetModelName;
424  c_void_MDS cinitializeInitialConditions;
425  c_void_MDS csetParameterValues;
426  c_void_MDS csetCompartmentVolumes;
427  c_int_MDS_int cgetNumLocalParameters;
428  c_void_MDS csetBoundaryConditions;
429  c_void_MDS csetInitialConditions;
430  c_void_MDS cevalInitialAssignments;
431  c_void_MDS cupdateDependentSpeciesValues;
432  c_void_MDS ccomputeRules;
433  c_void_MDS cconvertToAmounts;
434  c_void_MDS ccomputeConservedTotals;
435  c_double_MDS_int cgetConcentration;
436  c_doubleStar_MDS cGetCurrentValues;
437  c_void_MDS_double_doubleStar cevalModel;
438  c_void_MDS cconvertToConcentrations;
439  c_void_MDS_double_doubleStar cevalEvents;
440  c_void_MDS ccomputeAllRatesOfChange;
441  c_void_MDS cAssignRates_a;
442  c_void_MDS_doubleStar cAssignRates_b;
443  c_void_MDS ctestConstraints;
444  c_void_MDS cresetEvents;
445  c_void_MDS cInitializeRates;
446  c_void_MDS cInitializeRateRuleSymbols;
447  c_void_MDS_int_double csetConcentration;
448  c_void_MDS cComputeReactionRates;
449  c_void_MDS ccomputeEventPriorities;
450 
451  std::vector<PendingAssignment> mAssignments;
452 
453  std::vector<double> mAssignmentTimes;
454 
455  DoubleMatrix stoichiometryMatrix;
456 
457  std::vector<int> retestEvents(const double& timeEnd,
458  const std::vector<int>& handledEvents, std::vector<int>& removeEvents);
459 
460  std::vector<int> retestEvents(const double& timeEnd, std::vector<int>& handledEvents,
461  const bool& assignOldState);
462 
463  std::vector<int> retestEvents(const double& timeEnd,
464  const std::vector<int>& handledEvents,
465  const bool& assignOldState, std::vector<int>& removeEvents);
466 
467  virtual int applyPendingEvents(double timeEnd);
468 
469  virtual int applyEvents(double timeEnd, const unsigned char* previousEventStatus,
470  const double *initialState, double* finalState);
471 
472  virtual void getEventRoots(double time, const double* y, double* gdot);
473 
474  virtual double getNextPendingEventTime(bool pop);
475 
476  virtual int getPendingEventSize();
477 
478  void removePendingAssignmentForIndex(int eventIndex);
479 
480  void sortEventsByPriority(std::vector<Event>& firedEvents);
481 
482  void sortEventsByPriority(std::vector<int>& firedEvents);
483 
484 
485  virtual int setFloatingSpeciesInitAmounts(int len, int const *indx,
486  double const *values);
487 
488  virtual int getFloatingSpeciesInitAmounts(int len, int const *indx,
489  double *values);
490 
491  virtual int setCompartmentInitVolumes(int len, int const *indx,
492  double const *values);
493 
494  virtual int getCompartmentInitVolumes(int len, int const *indx,
495  double *values);
496 
497  virtual int getEventIndex(const std::string&);
498  virtual std::string getEventId(int);
499  virtual void setEventListener(int, rr::EventListenerPtr);
500  virtual rr::EventListenerPtr getEventListener(int);
501 
502  friend class CompiledModelState;
503 
504 
505 
506  virtual double getFloatingSpeciesAmountRate(int index,
507  const double *reactionRates);
508 
509  virtual void setRandomSeed(int64_t);
510 
511  virtual int64_t getRandomSeed();
512 
513  virtual double getRandom();
514 
515  virtual uint32_t getFlags() const
516  {
517  return 0;
518  }
519 
520  virtual void setFlags(uint32_t)
521  {
522 
523  }
524 };
525 }
526 #endif
Both the CModelGenerator and the CSharpModelGenerator use the same paradigm of producing source code,...
Definition: rrCompiledExecutableModel.h:43
ModelData mData
the data that is exchanged with the loaded shared lib, and all sorts of other routines such as CVODE.
Definition: rrCompiledExecutableModel.h:409
virtual uint32_t getFlags() const
Get the current set of flags.
Definition: rrCompiledExecutableModel.h:515
virtual void setFlags(uint32_t)
Set certain options that determine the state of the ExecutableModel, these are listed in.
Definition: rrCompiledExecutableModel.h:520
bool mIsInitialized
If all functions are found properly in the dll, this one is true.
Definition: rrCompiledExecutableModel.h:415
virtual std::string getExecutableModelDesc() const
Returns a human-readable description of the code generation backend, e.g. LLVM, legacy C,...
Definition: rrCompiledExecutableModel.h:53
StateStackOptions
Definition: rrCompiledExecutableModel.h:87
Saves the 'state' of a Model.
Definition: rrCompiledModelState.h:26
Base class for all code generation systems; allows compiling and evaluating the model.
Definition: rrExecutableModel.h:118
Access an actual compiled shared library (.so, .dll or .dylib) that was compiled by a ModelGenerator ...
Definition: rrModelSharedLibrary.h:22
Definition: rrModelSymbols.h:20
C_DECL_SPEC char *rrcCallConv getModelName(RRHandle handle)
Returns the name of currently loaded SBML model.
C_DECL_SPEC RRDoubleMatrixPtr rrcCallConv getStoichiometryMatrix(RRHandle handle)
Retrieve the stoichiometry matrix for the current model.
C_DECL_SPEC RRVectorPtr rrcCallConv getBoundarySpeciesAmounts(RRHandle handle)
Retrieve the amounts for all the boundary species in a vector.
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 RRVectorPtr rrcCallConv getFloatingSpeciesConcentrations(RRHandle handle)
Retrieve in a vector the concentrations for all the floating species.
C_DECL_SPEC RRVectorPtr rrcCallConv getFloatingSpeciesAmounts(RRHandle handle)
Retrieve in a vector the amounts for all the floating species.
C_DECL_SPEC bool rrcCallConv setFloatingSpeciesConcentrations(RRHandle handle, const RRVectorPtr vec)
Set the floating species concentration to the vector vec.
C_DECL_SPEC RRVectorPtr rrcCallConv getGlobalParameterValues(RRHandle handle)
Retrieve the values for all the global parameter values in a vector.
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 getValue(RRHandle handle, const char *symbolId, double *value)
Get the value for a given symbol, use getAvailableTimeCourseSymbols(void) for a list of symbols.
C_DECL_SPEC bool rrcCallConv setValue(RRHandle handle, const char *symbolId, const double value)
Set the value for a given symbol, use getAvailableTimeCourseSymbols(void) for a list of symbols.
C_DECL_SPEC char *rrcCallConv getInfo(RRHandle handle)
Retrieve info about current state of roadrunner, e.g. loaded model, conservationAnalysis etc.
Base class for all code generators in RoadRunner.
cxx11_ns::shared_ptr< EventListener > EventListenerPtr
listeners are shared objects, so use std smart pointers to manage them.
Definition: rrExecutableModel.h:71
A data structure that is that allows data to be exchanged with running SBML models.
Definition: rrModelData.h:49