COPASI API  4.16.103
Classes | Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | Friends | List of all members
CStochDirectMethod Class Reference

#include <CStochDirectMethod.h>

Inheritance diagram for CStochDirectMethod:
Inheritance graph
[legend]
Collaboration diagram for CStochDirectMethod:
Collaboration graph
[legend]

Classes

class  CReactionDependencies
 

Public Member Functions

 CStochDirectMethod (const CStochDirectMethod &src, const CCopasiContainer *pParent=NULL)
 
virtual bool elevateChildren ()
 
virtual bool isValidProblem (const CCopasiProblem *pProblem)
 
virtual void start (const CState *initialState)
 
virtual Status step (const double &deltaT)
 
 ~CStochDirectMethod ()
 
- Public Member Functions inherited from CTrajectoryMethod
 CTrajectoryMethod (const CTrajectoryMethod &src, const CCopasiContainer *pParent=NULL)
 
const CVector< C_INT > & getRoots () const
 
void setCurrentState (CState *currentState)
 
void setProblem (CTrajectoryProblem *problem)
 
virtual void stateChanged ()
 
 ~CTrajectoryMethod ()
 
- Public Member Functions inherited from CCopasiMethod
 CCopasiMethod (const CCopasiMethod &src, const CCopasiContainer *pParent=NULL)
 
const CCopasiMethod::SubTypegetSubType () const
 
const CCopasiTask::TypegetType () const
 
virtual void load (CReadConfig &configBuffer, CReadConfig::Mode mode=CReadConfig::SEARCH)
 
virtual void print (std::ostream *ostream) const
 
virtual void printResult (std::ostream *ostream) const
 
virtual bool setCallBack (CProcessReport *pCallBack)
 
virtual ~CCopasiMethod ()
 
- Public Member Functions inherited from CCopasiParameterGroup
bool addGroup (const std::string &name)
 
bool addParameter (const CCopasiParameter &parameter)
 
bool addParameter (const std::string &name, const CCopasiParameter::Type type)
 
template<class CType >
bool addParameter (const std::string &name, const CCopasiParameter::Type type, const CType &value)
 
void addParameter (CCopasiParameter *pParameter)
 
CCopasiParameterGroupassertGroup (const std::string &name)
 
template<class CType >
CCopasiParameterassertParameter (const std::string &name, const CCopasiParameter::Type type, const CType &defaultValue)
 
index_iterator beginIndex () const
 
name_iterator beginName () const
 
 CCopasiParameterGroup (const CCopasiParameterGroup &src, const CCopasiContainer *pParent=NULL)
 
 CCopasiParameterGroup (const std::string &name, const CCopasiContainer *pParent=NULL, const std::string &objectType="ParameterGroup")
 
void clear ()
 
index_iterator endIndex () const
 
name_iterator endName () const
 
CCopasiParameterGroupgetGroup (const std::string &name)
 
const CCopasiParameterGroupgetGroup (const std::string &name) const
 
CCopasiParameterGroupgetGroup (const size_t &index)
 
const CCopasiParameterGroupgetGroup (const size_t &index) const
 
size_t getIndex (const std::string &name) const
 
std::string getKey (const std::string &name) const
 
std::string getKey (const size_t &index) const
 
virtual const std::string & getName (const size_t &index) const
 
virtual const CObjectInterfacegetObject (const CCopasiObjectName &cn) const
 
CCopasiParametergetParameter (const std::string &name)
 
const CCopasiParametergetParameter (const std::string &name) const
 
CCopasiParametergetParameter (const size_t &index)
 
const CCopasiParametergetParameter (const size_t &index) const
 
CCopasiParameter::Type getType (const std::string &name) const
 
CCopasiParameter::Type getType (const size_t &index) const
 
std::string getUniqueParameterName (const CCopasiParameter *pParameter) const
 
const CCopasiParameter::ValuegetValue (const std::string &name) const
 
const CCopasiParameter::ValuegetValue (const size_t &index) const
 
CCopasiParameter::ValuegetValue (const std::string &name)
 
CCopasiParameter::ValuegetValue (const size_t &index)
 
CCopasiParameterGroupoperator= (const CCopasiParameterGroup &rhs)
 
bool removeParameter (const std::string &name)
 
bool removeParameter (const size_t &index)
 
template<class CType >
bool setValue (const std::string &name, const CType &value)
 
template<class CType >
bool setValue (const size_t &index, const CType &value)
 
size_t size () const
 
bool swap (const size_t &iFrom, const size_t &iTo)
 
bool swap (index_iterator &from, index_iterator &to)
 
virtual ~CCopasiParameterGroup ()
 
- Public Member Functions inherited from CCopasiParameter
 CCopasiParameter (const CCopasiParameter &src, const CCopasiContainer *pParent=NULL)
 
 CCopasiParameter (const std::string &name, const Type &type, const void *pValue=NULL, const CCopasiContainer *pParent=NULL, const std::string &objectType="Parameter")
 
virtual CCopasiObjectName getCN () const
 
virtual const std::string & getKey () const
 
virtual std::string getObjectDisplayName (bool regular=true, bool richtext=false) const
 
const CCopasiParameter::TypegetType () const
 
const ValuegetValue () const
 
ValuegetValue ()
 
virtual voidgetValuePointer () const
 
CCopasiObjectgetValueReference () const
 
bool isValidValue (const C_FLOAT64 &value) const
 
bool isValidValue (const C_INT32 &value) const
 
bool isValidValue (const unsigned C_INT32 &value) const
 
bool isValidValue (const bool &value) const
 
bool isValidValue (const std::string &value) const
 
bool isValidValue (const CCopasiObjectName &value) const
 
bool isValidValue (const std::vector< CCopasiParameter * > &value) const
 
CCopasiParameteroperator= (const CCopasiParameter &rhs)
 
template<class CType >
bool setValue (const CType &value)
 
bool setValue (const std::vector< CCopasiParameter * > &value)
 
virtual ~CCopasiParameter ()
 
- Public Member Functions inherited from CCopasiContainer
virtual bool add (CCopasiObject *pObject, const bool &adopt=true)
 
 CCopasiContainer (const std::string &name, const CCopasiContainer *pParent=NULL, const std::string &type="CN", const unsigned C_INT32 &flag=CCopasiObject::Container)
 
 CCopasiContainer (const CCopasiContainer &src, const CCopasiContainer *pParent=NULL)
 
virtual std::string getChildObjectUnits (const CCopasiObject *pObject) const
 
virtual const objectMapgetObjects () const
 
virtual std::string getUnits () const
 
virtual const CCopasiObjectgetValueObject () const
 
virtual bool remove (CCopasiObject *pObject)
 
virtual ~CCopasiContainer ()
 
- Public Member Functions inherited from CCopasiObject
void addDirectDependency (const CCopasiObject *pObject)
 
 CCopasiObject (const CCopasiObject &src, const CCopasiContainer *pParent=NULL)
 
void clearDirectDependencies ()
 
void clearRefresh ()
 
bool dependsOn (DataObjectSet candidates, const DataObjectSet &context=DataObjectSet()) const
 
void getAllDependencies (DataObjectSet &dependencies, const DataObjectSet &context) const
 
virtual const DataObjectSetgetDirectDependencies (const DataObjectSet &context=DataObjectSet()) const
 
CCopasiContainergetObjectAncestor (const std::string &type) const
 
CCopasiDataModelgetObjectDataModel ()
 
const CCopasiDataModelgetObjectDataModel () const
 
const std::string & getObjectName () const
 
CCopasiContainergetObjectParent () const
 
const std::string & getObjectType () const
 
virtual const
CObjectInterface::ObjectSet
getPrerequisites () const
 
virtual RefreshgetRefresh () const
 
UpdateMethodgetUpdateMethod () const
 
bool hasCircularDependencies (DataObjectSet &candidates, DataObjectSet &verified, const DataObjectSet &context) const
 
bool hasUpdateMethod () const
 
bool isArray () const
 
bool isContainer () const
 
bool isDataModel () const
 
bool isMatrix () const
 
bool isNameVector () const
 
bool isNonUniqueName () const
 
virtual bool isPrerequisiteForContext (const CObjectInterface *pObject, const CMath::SimulationContextFlag &context, const CObjectInterface::ObjectSet &changedObjects) const
 
bool isReference () const
 
bool isRoot () const
 
bool isSeparator () const
 
bool isStaticString () const
 
bool isValueBool () const
 
bool isValueDbl () const
 
bool isValueInt () const
 
bool isValueInt64 () const
 
bool isValueString () const
 
bool isVector () const
 
virtual bool mustBeDeleted (const DataObjectSet &deletedObjects) const
 
void removeDirectDependency (const CCopasiObject *pObject)
 
void setDirectDependencies (const DataObjectSet &directDependencies)
 
bool setObjectName (const std::string &name)
 
virtual bool setObjectParent (const CCopasiContainer *pParent)
 
void setObjectValue (const C_FLOAT64 &value)
 
void setObjectValue (const C_INT32 &value)
 
void setObjectValue (const bool &value)
 
template<class CType >
void setRefresh (CType *pType, void(CType::*method)(void))
 
template<class CType >
void setUpdateMethod (CType *pType, void(CType::*method)(const C_FLOAT64 &))
 
template<class CType >
void setUpdateMethod (CType *pType, void(CType::*method)(const C_INT32 &))
 
template<class CType >
void setUpdateMethod (CType *pType, void(CType::*method)(const bool &))
 
virtual ~CCopasiObject ()
 
- Public Member Functions inherited from CObjectInterface
 CObjectInterface ()
 
virtual ~CObjectInterface ()
 

Protected Member Functions

void calculateAmu (const size_t &index)
 
 CStochDirectMethod (const CCopasiContainer *pParent=NULL)
 
C_FLOAT64 doSingleStep (const C_FLOAT64 &curTime, const C_FLOAT64 &endTime)
 
- Protected Member Functions inherited from CTrajectoryMethod
 CTrajectoryMethod (const CCopasiMethod::SubType &subType, const CCopasiContainer *pParent=NULL)
 
- Protected Member Functions inherited from CCopasiMethod
 CCopasiMethod (const CCopasiTask::Type &taskType, const SubType &subType, const CCopasiContainer *pParent=NULL)
 
- Protected Member Functions inherited from CCopasiParameterGroup
 CCopasiParameterGroup ()
 
- Protected Member Functions inherited from CCopasiContainer
template<class CType >
CCopasiObjectaddMatrixReference (const std::string &name, CType &reference, const unsigned C_INT32 &flag=0)
 
template<class CType >
CCopasiObjectaddObjectReference (const std::string &name, CType &reference, const unsigned C_INT32 &flag=0)
 
template<class CType >
CCopasiObjectaddVectorReference (const std::string &name, CType &reference, const unsigned C_INT32 &flag=0)
 
void initObjects ()
 
- Protected Member Functions inherited from CCopasiObject
 CCopasiObject ()
 
 CCopasiObject (const std::string &name, const CCopasiContainer *pParent=NULL, const std::string &type="CN", const unsigned C_INT32 &flag=0)
 

Protected Attributes

C_FLOAT64 mA0
 
CVector< C_FLOAT64mAmu
 
bool mDoCorrection
 
unsigned C_INT32 mMaxSteps
 
bool mMaxStepsReached
 
CState mMethodState
 
size_t mNextReactionIndex
 
C_FLOAT64 mNextReactionTime
 
size_t mNumReactions
 
CModelmpModel
 
CRandommpRandomGenerator
 
std::vector
< CReactionDependencies
mReactionDependencies
 
- Protected Attributes inherited from CTrajectoryMethod
CStatempCurrentState
 
CTrajectoryProblemmpProblem
 
CVector< C_INTmRoots
 
- Protected Attributes inherited from CCopasiMethod
CProcessReportmpCallBack
 
- Protected Attributes inherited from CCopasiParameter
std::string mKey
 
CCopasiObjectmpValueReference
 
size_t mSize
 
Value mValue
 
- Protected Attributes inherited from CCopasiContainer
objectMap mObjects
 

Private Member Functions

void initializeParameter ()
 

Friends

CTrajectoryMethodCTrajectoryMethod::createMethod (CCopasiMethod::SubType subType)
 

Additional Inherited Members

- Public Types inherited from CTrajectoryMethod
enum  Status { FAILURE = -1, NORMAL = 0, ROOT = 1 }
 
- Public Types inherited from CCopasiMethod
enum  SubType {
  unset = 0, RandomSearch, RandomSearchMaster, SimulatedAnnealing,
  CoranaWalk, DifferentialEvolution, ScatterSearch, GeneticAlgorithm,
  EvolutionaryProgram, SteepestDescent, HybridGASA, GeneticAlgorithmSR,
  HookeJeeves, LevenbergMarquardt, NelderMead, SRES,
  Statistics, ParticleSwarm, Praxis, TruncatedNewton,
  Newton, deterministic, LSODAR, directMethod,
  stochastic, tauLeap, adaptiveSA, hybrid,
  hybridLSODA, hybridODE45, DsaLsodar, tssILDM,
  tssILDMModified, tssCSP, mcaMethodReder, scanMethod,
  lyapWolf, sensMethod, EFMAlgorithm, EFMBitPatternTreeAlgorithm,
  EFMBitPatternAlgorithm, Householder, crossSectionMethod, linearNoiseApproximation
}
 
- Public Types inherited from CCopasiParameterGroup
typedef parameterGroup::iterator index_iterator
 
typedef
CCopasiContainer::objectMap::iterator 
name_iterator
 
typedef std::vector
< CCopasiParameter * > 
parameterGroup
 
- Public Types inherited from CCopasiParameter
enum  Type {
  DOUBLE = 0, UDOUBLE, INT, UINT,
  BOOL, GROUP, STRING, CN,
  KEY, FILE, EXPRESSION, INVALID
}
 
- Public Types inherited from CCopasiContainer
typedef std::multimap
< std::string, CCopasiObject * > 
objectMap
 
- Public Types inherited from CCopasiObject
typedef std::set< const
CCopasiObject * > 
DataObjectSet
 
typedef std::vector< Refresh * > DataUpdateSequence
 
- Public Types inherited from CObjectInterface
typedef std::set< const
CObjectInterface * > 
ObjectSet
 
typedef std::vector
< CObjectInterface * > 
UpdateSequence
 
- Static Public Member Functions inherited from CTrajectoryMethod
static CTrajectoryMethodcreateMethod (CCopasiMethod::SubType subType=CCopasiMethod::deterministic)
 
- Static Public Member Functions inherited from CCopasiObject
static std::vector< Refresh * > buildUpdateSequence (const DataObjectSet &objects, const DataObjectSet &uptoDateObjects, const DataObjectSet &context=DataObjectSet())
 
static void setRenameHandler (CRenameHandler *rh)
 
- Static Public Attributes inherited from CCopasiMethod
static const std::string SubTypeName []
 
static const char * XMLSubType []
 
- Static Public Attributes inherited from CCopasiParameter
static const std::string TypeName []
 
static const char * XMLType []
 
- Static Public Attributes inherited from CCopasiContainer
static const std::vector
< CCopasiContainer * > 
EmptyList
 
- Protected Types inherited from CCopasiObject
enum  Flag {
  Container = 0x1, Vector = 0x2, Matrix = 0x4, NameVector = 0x8,
  Reference = 0x10, ValueBool = 0x20, ValueInt = 0x40, ValueInt64 = 0x80,
  ValueDbl = 0x100, NonUniqueName = 0x200, StaticString = 0x400, ValueString = 0x800,
  Separator = 0x1000, ModelEntity = 0x2000, Array = 0x4000, DataModel = 0x8000,
  Root = 0x10000, Gui = 0x20000
}
 
- Static Protected Attributes inherited from CCopasiObject
static CRenameHandlersmpRenameHandler = NULL
 

Detailed Description

Definition at line 39 of file CStochDirectMethod.h.

Constructor & Destructor Documentation

CStochDirectMethod::CStochDirectMethod ( const CCopasiContainer pParent = NULL)
protected

Default constructor.

Parameters
constCCopasiContainer * pParent (default: NULL)

Definition at line 86 of file CStochDirectMethod.cpp.

References initializeParameter().

86  :
89  mpModel(NULL),
90  mNumReactions(0),
91  mMaxSteps(1000000),
92  mNextReactionTime(0.0),
94  mDoCorrection(true),
95  mAmu(0),
96  mA0(0.0),
97  mMethodState(),
99  mMaxStepsReached(false)
100 {
102 }
CVector< C_FLOAT64 > mAmu
std::vector< CReactionDependencies > mReactionDependencies
#define C_INVALID_INDEX
Definition: copasi.h:222
static CRandom * createGenerator(CRandom::Type type=CRandom::mt19937, unsigned C_INT32 seed=0)
Definition: CRandom.cpp:49
unsigned C_INT32 mMaxSteps
CStochDirectMethod::CStochDirectMethod ( const CStochDirectMethod src,
const CCopasiContainer pParent = NULL 
)

Copy constructor.

Parameters
constCStochDirectMethod & src,
constCCopasiContainer * pParent (Default: NULL)

Definition at line 104 of file CStochDirectMethod.cpp.

References initializeParameter().

105  :
106  CTrajectoryMethod(src, pParent),
108  mpModel(src.mpModel),
110  mMaxSteps(src.mMaxSteps),
114  mAmu(src.mAmu),
115  mA0(src.mA0),
119 {
121 }
CVector< C_FLOAT64 > mAmu
std::vector< CReactionDependencies > mReactionDependencies
static CRandom * createGenerator(CRandom::Type type=CRandom::mt19937, unsigned C_INT32 seed=0)
Definition: CRandom.cpp:49
unsigned C_INT32 mMaxSteps
CStochDirectMethod::~CStochDirectMethod ( )

Destructor.

Definition at line 123 of file CStochDirectMethod.cpp.

References mpRandomGenerator, and pdelete.

124 {
126 }
#define pdelete(p)
Definition: copasi.h:215

Member Function Documentation

void CStochDirectMethod::calculateAmu ( const size_t &  index)
protected

Calculate the propensity of the indexed reaction

Parameters
constsize_t & index

Definition at line 522 of file CStochDirectMethod.cpp.

References CVectorCore< CType >::array(), C_FLOAT64, mAmu, mDoCorrection, CStochDirectMethod::CReactionDependencies::mMethodSubstrates, CStochDirectMethod::CReactionDependencies::mModelSubstrates, CStochDirectMethod::CReactionDependencies::mpParticleFlux, mReactionDependencies, CStochDirectMethod::CReactionDependencies::mSubstrateMultiplier, and CVectorCore< CType >::size().

Referenced by doSingleStep(), and start().

523 {
524  CReactionDependencies & Dependencies = mReactionDependencies[index];
525  C_FLOAT64 & Amu = mAmu[index];
526 
527  Amu = *Dependencies.mpParticleFlux;
528 
529  if (Amu < 0.0)
530  {
531  // TODO CRITICAL Create a warning message
532  Amu = 0.0;
533  }
534 
535  if (!mDoCorrection)
536  {
537  return;
538  }
539 
540  C_FLOAT64 SubstrateMultiplier = 1.0;
541  C_FLOAT64 SubstrateDevisor = 1.0;
542  C_FLOAT64 Multiplicity;
543  C_FLOAT64 LowerBound;
544  C_FLOAT64 Number;
545 
546  bool ApplyCorrection = false;
547 
548  C_FLOAT64 * pMultiplier = Dependencies.mSubstrateMultiplier.array();
549  C_FLOAT64 * endMultiplier = pMultiplier + Dependencies.mSubstrateMultiplier.size();
550  C_FLOAT64 ** ppLocalSubstrate = Dependencies.mMethodSubstrates.array();
551  C_FLOAT64 ** ppModelSubstrate = Dependencies.mModelSubstrates.array();
552 
553  for (; pMultiplier != endMultiplier; ++pMultiplier, ++ppLocalSubstrate, ++ppModelSubstrate)
554  {
555  Multiplicity = *pMultiplier;
556 
557  // TODO We should check the error introduced through rounding.
558  **ppLocalSubstrate = floor(**ppModelSubstrate + 0.5);
559 
560  if (Multiplicity > 1.01)
561  {
562  ApplyCorrection = true;
563 
564  Number = **ppLocalSubstrate;
565 
566  LowerBound = Number - Multiplicity;
567  SubstrateDevisor *= pow(Number, Multiplicity - 1.0); //optimization
568  Number -= 1.0;
569 
570  while (Number > LowerBound)
571  {
572  SubstrateMultiplier *= Number;
573  Number -= 1.0;
574  }
575  }
576  }
577 
578  // at least one substrate particle number is zero
579  if (SubstrateMultiplier < 0.5 || SubstrateDevisor < 0.5)
580  {
581  Amu = 0.0;
582  return;
583  }
584 
585  // rate_factor is the rate function divided by substrate_factor.
586  // It would be more efficient if this was generated directly, since in effect we
587  // are multiplying and then dividing by the same thing (substrate_factor)!
588  //C_FLOAT64 rate_factor = mpModel->getReactions()[index]->calculateParticleFlux();
589  if (ApplyCorrection)
590  {
591  Amu *= SubstrateMultiplier / SubstrateDevisor;
592  }
593 
594  return;
595 }
CVector< C_FLOAT64 > mAmu
std::vector< CReactionDependencies > mReactionDependencies
#define C_FLOAT64
Definition: copasi.h:92
C_FLOAT64 CStochDirectMethod::doSingleStep ( const C_FLOAT64 curTime,
const C_FLOAT64 endTime 
)
protected

Fire the next reaction if it fire before the endTime

Parameters
constC_FLOAT64 & curTime
constC_FLOAT64 & endTime
Returns
C_FLOAT64 timeAfterStep

Definition at line 429 of file CStochDirectMethod.cpp.

References CVectorCore< CType >::array(), C_FLOAT64, C_INVALID_INDEX, calculateAmu(), CCopasiMessage::EXCEPTION, CRandom::getRandomOO(), mA0, mAmu, CStochDirectMethod::CReactionDependencies::mCalculations, MCTrajectoryMethod, CStochDirectMethod::CReactionDependencies::mDependentReactions, CStochDirectMethod::CReactionDependencies::mMethodSpecies, CStochDirectMethod::CReactionDependencies::mModelSpecies, mNextReactionIndex, mNextReactionTime, mNumReactions, mpModel, mpRandomGenerator, mReactionDependencies, CStochDirectMethod::CReactionDependencies::mSpeciesMultiplier, CModel::setTime(), and CVectorCore< CType >::size().

Referenced by step().

430 {
432  {
433  if (mA0 == 0)
434  {
435  return endTime - curTime;
436  }
437 
438  // We need to throw an exception if mA0 is NaN
439  if (isnan(mA0))
440  {
442  }
443 
444  mNextReactionTime = curTime - log(mpRandomGenerator->getRandomOO()) / mA0;
445 
446  // We are sure that we have at least 1 reaction
447  mNextReactionIndex = 0;
448 
449  C_FLOAT64 sum = 0.0;
451 
452  C_FLOAT64 * pAmu = mAmu.array();
453  C_FLOAT64 * endAmu = pAmu + mNumReactions;
454 
455  for (; (sum < rand) && (pAmu != endAmu); ++pAmu, ++mNextReactionIndex)
456  {
457  sum += *pAmu;
458  }
459 
460  mNextReactionIndex--;
461  }
462 
463  if (mNextReactionTime >= endTime)
464  {
465  return endTime - curTime;
466  }
467 
468  // Update the model time (for explicitly time dependent models)
470 
471  CReactionDependencies & Dependencies = mReactionDependencies[mNextReactionIndex];
472 
473  // Update the method internal and model species numbers
474  C_FLOAT64 ** ppModelSpecies = Dependencies.mModelSpecies.array();
475  C_FLOAT64 ** ppLocalSpecies = Dependencies.mMethodSpecies.array();
476  C_FLOAT64 * pMultiplier = Dependencies.mSpeciesMultiplier.array();
477  C_FLOAT64 * endMultiplier = pMultiplier + Dependencies.mSpeciesMultiplier.size();
478 
479  for (; pMultiplier != endMultiplier; ++pMultiplier, ++ppLocalSpecies, ++ppModelSpecies)
480  {
481  **ppLocalSpecies += *pMultiplier;
482  **ppModelSpecies = **ppLocalSpecies;
483  }
484 
485  // Calculate all values which depend on the firing reaction
486  std::vector< Refresh * >::const_iterator itCalcualtion = Dependencies.mCalculations.begin();
487  std::vector< Refresh * >::const_iterator endCalcualtion = Dependencies.mCalculations.end();
488 
489  while (itCalcualtion != endCalcualtion)
490  {
491  (**itCalcualtion++)();
492  }
493 
494  // We do not need to update the the method state since the only independent state
495  // values are species of type reaction which are all controlled by the method.
496 
497  // calculate the propensities which depend on the firing reaction
498  size_t * pDependentReaction = Dependencies.mDependentReactions.array();
499  size_t * endDependentReactions = pDependentReaction + Dependencies.mDependentReactions.size();
500 
501  for (; pDependentReaction != endDependentReactions; ++pDependentReaction)
502  {
503  calculateAmu(*pDependentReaction);
504  }
505 
506  // calculate the total propensity
507  C_FLOAT64 * pAmu = mAmu.array();
508  C_FLOAT64 * endAmu = pAmu + mNumReactions;
509 
510  mA0 = 0.0;
511 
512  for (; pAmu != endAmu; ++pAmu)
513  {
514  mA0 += *pAmu;
515  }
516 
517  mNextReactionIndex = C_INVALID_INDEX;
518 
519  return mNextReactionTime - curTime;
520 }
virtual C_FLOAT64 getRandomOO()
Definition: CRandom.cpp:254
void calculateAmu(const size_t &index)
CVector< C_FLOAT64 > mAmu
std::vector< CReactionDependencies > mReactionDependencies
#define C_INVALID_INDEX
Definition: copasi.h:222
#define MCTrajectoryMethod
void setTime(const C_FLOAT64 &time)
Definition: CModel.cpp:1187
#define C_FLOAT64
Definition: copasi.h:92
CType * array()
Definition: CVector.h:139
bool CStochDirectMethod::elevateChildren ( )
virtual

This methods must be called to elevate subgroups to derived objects. The default implementation does nothing.

Returns
bool success

Reimplemented from CCopasiParameterGroup.

Definition at line 135 of file CStochDirectMethod.cpp.

References initializeParameter().

136 {
138  return true;
139 }
void CStochDirectMethod::initializeParameter ( )
private

Initialize the method parameter

Definition at line 128 of file CStochDirectMethod.cpp.

References CCopasiParameterGroup::assertParameter(), CCopasiParameter::BOOL, C_INT32, CCopasiParameter::INT, and CCopasiParameter::UINT.

Referenced by CStochDirectMethod(), and elevateChildren().

129 {
130  assertParameter("Max Internal Steps", CCopasiParameter::INT, (C_INT32) 1000000);
131  assertParameter("Use Random Seed", CCopasiParameter::BOOL, false);
132  assertParameter("Random Seed", CCopasiParameter::UINT, (unsigned C_INT32) 1);
133 }
#define C_INT32
Definition: copasi.h:90
CCopasiParameter * assertParameter(const std::string &name, const CCopasiParameter::Type type, const CType &defaultValue)
bool CStochDirectMethod::isValidProblem ( const CCopasiProblem pProblem)
virtual

Check if the method is suitable for this problem

Returns
bool suitability of the method

Reimplemented from CTrajectoryMethod.

Definition at line 351 of file CStochDirectMethod.cpp.

References CStateTemplate::beginIndependent(), CStateTemplate::endIndependent(), CCopasiMessage::ERROR, CTrajectoryProblem::getDuration(), CModel::getEvents(), CCopasiProblem::getModel(), CModel::getStateTemplate(), CModel::getTotSteps(), CCopasiParameter::getValue(), CTrajectoryMethod::isValidProblem(), MCTrajectoryMethod, CModelEntity::ODE, CCopasiParameter::Value::pINT, CCopasiVector< T >::size(), and CModel::suitableForStochasticSimulation().

352 {
353  if (!CTrajectoryMethod::isValidProblem(pProblem)) return false;
354 
355  const CTrajectoryProblem * pTP = dynamic_cast<const CTrajectoryProblem *>(pProblem);
356 
357  if (pTP->getDuration() < 0.0)
358  {
359  //back integration not possible
361  return false;
362  }
363 
364  if (pTP->getModel()->getTotSteps() < 1)
365  {
366  //at least one reaction necessary
368  return false;
369  }
370 
371  // check for ODEs
372  const CStateTemplate & StateTemplate = pTP->getModel()->getStateTemplate();
373  CModelEntity *const* ppEntity = StateTemplate.beginIndependent();
374  CModelEntity *const* ppEntityEnd = StateTemplate.endIndependent();
375 
376  for (; ppEntity != ppEntityEnd; ++ppEntity)
377  {
378  if ((*ppEntity)->getStatus() == CModelEntity::ODE)
379  {
380  if (dynamic_cast<const CModelValue *>(*ppEntity) != NULL)
381  {
382  // global quantity ode rule found
384  return false;
385  }
386  else if (dynamic_cast<const CCompartment *>(*ppEntity) != NULL)
387  {
388  // compartment ode rule found
390  return false;
391  }
392  else
393  {
394  // species ode rule found
396  return false;
397  }
398  }
399  }
400 
401  //TODO: rewrite CModel::suitableForStochasticSimulation() to use
402  // CCopasiMessage
403  std::string message = pTP->getModel()->suitableForStochasticSimulation();
404 
405  if (message != "")
406  {
407  //model not suitable, message describes the problem
408  CCopasiMessage(CCopasiMessage::ERROR, message.c_str());
409  return false;
410  }
411 
412  if (* getValue("Max Internal Steps").pINT <= 0)
413  {
414  //max steps should be at least 1
416  return false;
417  }
418 
419  //events are not supported at the moment
420  if (pTP->getModel()->getEvents().size() > 0)
421  {
423  return false;
424  }
425 
426  return true;
427 }
CCopasiVectorN< CEvent > & getEvents()
Definition: CModel.cpp:1110
virtual size_t size() const
virtual bool isValidProblem(const CCopasiProblem *pProblem)
size_t getTotSteps() const
Definition: CModel.cpp:1136
#define MCTrajectoryMethod
const C_FLOAT64 & getDuration() const
const Value & getValue() const
CModelEntity ** beginIndependent()
Definition: CState.cpp:208
const CStateTemplate & getStateTemplate() const
Definition: CModel.cpp:1172
std::string suitableForStochasticSimulation() const
Definition: CModel.cpp:3901
CModel * getModel() const
CModelEntity ** endIndependent()
Definition: CState.cpp:209
void CStochDirectMethod::start ( const CState initialState)
virtual

This instructs the method to prepare for integration starting with the initialState given.

Parameters
const CState *initialState

Reimplemented from CTrajectoryMethod.

Definition at line 166 of file CStochDirectMethod.cpp.

References CCopasiVector< T >::begin(), CStateTemplate::beginIndependent(), CState::beginIndependent(), CCopasiObject::buildUpdateSequence(), C_FLOAT64, C_INT32, C_INVALID_INDEX, calculateAmu(), CModel::deterministic, CCopasiVector< T >::end(), CStateTemplate::endFixed(), CStateTemplate::getIndex(), CCopasiProblem::getModel(), CModel::getModelType(), CModel::getReactions(), CModel::getState(), CModel::getStateTemplate(), CModelEntity::getStatus(), CState::getTime(), CCopasiParameter::getValue(), CCopasiObject::getValuePointer(), CModelEntity::getValueReference(), CRandom::initialize(), mA0, mAmu, mDoCorrection, mMaxSteps, mMaxStepsReached, mMethodState, mNextReactionIndex, mNextReactionTime, mNumReactions, CTrajectoryMethod::mpCurrentState, mpModel, CTrajectoryMethod::mpProblem, mpRandomGenerator, mReactionDependencies, CCopasiParameter::Value::pBOOL, CCopasiParameter::Value::pINT, CCopasiParameter::Value::pUINT, CModelEntity::REACTIONS, CVector< CType >::resize(), CModel::setState(), CCopasiVector< T >::size(), and CModel::updateSimulatedValues().

167 {
168  /* get configuration data */
169  mMaxSteps = * getValue("Max Internal Steps").pINT;
170 
171  bool useRandomSeed = * getValue("Use Random Seed").pBOOL;
172  unsigned C_INT32 randomSeed = * getValue("Random Seed").pUINT;
173 
174  if (useRandomSeed) mpRandomGenerator->initialize(randomSeed);
175 
176  //mpCurrentState is initialized. This state is not used internally in the
177  //stochastic solver, but it is used for returning the result after each step.
178  *mpCurrentState = *initialState;
179 
181  assert(mpModel);
182 
184  mDoCorrection = true;
185  else
186  mDoCorrection = false;
187 
188  //initialize the vector of ints that contains the particle numbers
189  //for the discrete simulation. This also floors all particle numbers in the model.
190 
192 
195  mAmu = 0.0;
196 
197  // Create a local copy of the state where the particle number are rounded to integers
199 
200  const CStateTemplate & StateTemplate = mpModel->getStateTemplate();
201 
202  CModelEntity *const* ppEntity = StateTemplate.beginIndependent();
203  CModelEntity *const* endEntity = StateTemplate.endFixed();
205 
206  for (; ppEntity != endEntity; ++ppEntity, ++pValue)
207  {
208  if (dynamic_cast< const CMetab * >(*ppEntity) != NULL)
209  {
210  *pValue = floor(*pValue + 0.5);
211  }
212  }
213 
214  // Update the model state so that the species are all represented by integers.
216  mpModel->updateSimulatedValues(false); //for assignments
217 
218  // TODO CRITICAL Handle species of type ASSIGNMENT.
219  // These need to be checked whether they are sufficiently close to an integer
220 
221  C_FLOAT64 * pMethodStateValue = mMethodState.beginIndependent() - 1;
222 
223  // Build the reaction dependencies
224  size_t NumReactions = 0;
225 
228  std::vector< CReactionDependencies >::iterator itDependencies = mReactionDependencies.begin();
229 
230  for (; it != end; ++it)
231  {
232  const CCopasiVector<CChemEqElement> & Balances = (*it)->getChemEq().getBalances();
233  const CCopasiVector<CChemEqElement> & Substrates = (*it)->getChemEq().getSubstrates();
234 
235  // This reactions does not change anything we ignore it
236  if (Balances.size() == 0 && Substrates.size() == 0)
237  {
238  continue;
239  }
240 
241  itDependencies->mpParticleFlux = (C_FLOAT64 *)(*it)->getParticleFluxReference()->getValuePointer();
242 
243  itDependencies->mSpeciesMultiplier.resize(Balances.size());
244  itDependencies->mMethodSpecies.resize(Balances.size());
245  itDependencies->mModelSpecies.resize(Balances.size());
246 
247  std::set< const CCopasiObject * > changed;
248 
249  // The time is always updated
250  changed.insert(mpModel->getValueReference());
251 
254 
255  size_t Index = 0;
256 
257  for (; itBalance != endBalance; ++itBalance)
258  {
259  const CMetab * pMetab = (*itBalance)->getMetabolite();
260 
261  if (pMetab->getStatus() == CModelEntity::REACTIONS)
262  {
263  itDependencies->mSpeciesMultiplier[Index] = floor((*itBalance)->getMultiplicity() + 0.5);
264  itDependencies->mMethodSpecies[Index] = pMethodStateValue + StateTemplate.getIndex(pMetab);
265  itDependencies->mModelSpecies[Index] = (C_FLOAT64 *) pMetab->getValueReference()->getValuePointer();
266 
267  changed.insert(pMetab->getValueReference());
268 
269  Index++;
270  }
271  }
272 
273  // Correct allocation for metabolites which are not determined by reactions
274  itDependencies->mSpeciesMultiplier.resize(Index, true);
275  itDependencies->mMethodSpecies.resize(Index, true);
276  itDependencies->mModelSpecies.resize(Index, true);
277 
278  itDependencies->mSubstrateMultiplier.resize(Substrates.size());
279  itDependencies->mMethodSubstrates.resize(Substrates.size());
280  itDependencies->mModelSubstrates.resize(Substrates.size());
281 
282  CCopasiVector< CChemEqElement >::const_iterator itSubstrate = Substrates.begin();
283  CCopasiVector< CChemEqElement >::const_iterator endSubstrate = Substrates.end();
284 
285  Index = 0;
286 
287  for (; itSubstrate != endSubstrate; ++itSubstrate, ++Index)
288  {
289  const CMetab * pMetab = (*itSubstrate)->getMetabolite();
290 
291  itDependencies->mSubstrateMultiplier[Index] = floor((*itSubstrate)->getMultiplicity() + 0.5);
292  itDependencies->mMethodSubstrates[Index] = pMethodStateValue + StateTemplate.getIndex(pMetab);
293  itDependencies->mModelSubstrates[Index] = (C_FLOAT64 *) pMetab->getValueReference()->getValuePointer();
294  }
295 
296  std::set< const CCopasiObject * > dependend;
297 
299  itDependencies->mDependentReactions.resize(mNumReactions);
300 
301  Index = 0;
302  size_t Count = 0;
303 
304  for (; itReaction != end; ++itReaction, ++Index)
305  {
306  if ((*itReaction)->getParticleFluxReference()->dependsOn(changed))
307  {
308  dependend.insert((*itReaction)->getParticleFluxReference());
309  itDependencies->mDependentReactions[Count] = Index;
310 
311  Count++;
312  }
313  }
314 
315  itDependencies->mDependentReactions.resize(Count, true);
316 
317  itDependencies->mCalculations = CCopasiObject::buildUpdateSequence(dependend, changed);
318 
319  ++itDependencies;
320  ++NumReactions;
321  }
322 
323  mNumReactions = NumReactions;
324 
326  mAmu.resize(mNumReactions, true);
327 
328  size_t i;
329 
330  for (i = 0; i < mNumReactions; i++)
331  {
332  calculateAmu(i);
333  }
334 
335  mA0 = 0;
336 
337  for (i = 0; i < mNumReactions; i++)
338  {
339  mA0 += mAmu[i];
340  }
341 
342  mMaxStepsReached = false;
343 
346 
347  return;
348 }
virtual size_t size() const
void calculateAmu(const size_t &index)
CVector< C_FLOAT64 > mAmu
iterator begin()
void updateSimulatedValues(const bool &updateMoieties)
Definition: CModel.cpp:1851
std::vector< CReactionDependencies > mReactionDependencies
CTrajectoryProblem * mpProblem
void resize(size_t size, const bool &copy=false)
Definition: CVector.h:301
#define C_INVALID_INDEX
Definition: copasi.h:222
size_t getIndex(const CModelEntity *entity) const
Definition: CState.cpp:231
#define C_INT32
Definition: copasi.h:90
Definition: CMetab.h:178
iterator end()
std::vector< CType * >::const_iterator const_iterator
Definition: CCopasiVector.h:57
CModelEntity ** endFixed()
Definition: CState.cpp:213
void setState(const CState &state)
Definition: CModel.cpp:1785
const Value & getValue() const
unsigned C_INT32 * pUINT
CModelEntity ** beginIndependent()
Definition: CState.cpp:208
unsigned C_INT32 mMaxSteps
const C_FLOAT64 & getTime() const
Definition: CState.cpp:325
#define C_FLOAT64
Definition: copasi.h:92
const CStateTemplate & getStateTemplate() const
Definition: CModel.cpp:1172
const ModelType & getModelType() const
Definition: CModel.cpp:2339
CCopasiVectorNS< CReaction > & getReactions()
Definition: CModel.cpp:1039
static std::vector< Refresh * > buildUpdateSequence(const DataObjectSet &objects, const DataObjectSet &uptoDateObjects, const DataObjectSet &context=DataObjectSet())
virtual void * getValuePointer() const
const CState & getState() const
Definition: CModel.cpp:1771
const CModelEntity::Status & getStatus() const
CModel * getModel() const
C_FLOAT64 * beginIndependent()
Definition: CState.cpp:328
virtual void initialize(unsigned C_INT32 seed=CRandom::getSystemSeed())
Definition: CRandom.cpp:162
CCopasiObject * getValueReference() const
CTrajectoryMethod::Status CStochDirectMethod::step ( const double &  deltaT)
virtual

This instructs the method to calculate a time step of deltaT starting with the current state, i.e., the result of the previous step. The new state (after deltaT) is expected in the current state. The return value is the actual timestep taken.

Parameters
const double &deltaT
Returns
Status status

Reimplemented from CTrajectoryMethod.

Definition at line 141 of file CStochDirectMethod.cpp.

References C_FLOAT64, doSingleStep(), CCopasiMessage::EXCEPTION, CCopasiProblem::getModel(), CModel::getState(), CState::getTime(), MCTrajectoryMethod, mMaxSteps, mMethodState, CTrajectoryMethod::mpCurrentState, CTrajectoryMethod::mpProblem, CTrajectoryMethod::NORMAL, and CState::setTime().

142 {
143  // do several steps
144  C_FLOAT64 Time = mpCurrentState->getTime();
145  C_FLOAT64 EndTime = Time + deltaT;
146 
147  size_t Steps = 0;
148 
149  while (Time < EndTime)
150  {
151  mMethodState.setTime(Time);
152  Time += doSingleStep(Time, EndTime);
153 
154  if (++Steps > mMaxSteps)
155  {
157  }
158  }
159 
161  mpCurrentState->setTime(Time);
162 
163  return NORMAL;
164 }
CTrajectoryProblem * mpProblem
void setTime(const C_FLOAT64 &time)
Definition: CState.cpp:326
#define MCTrajectoryMethod
C_FLOAT64 doSingleStep(const C_FLOAT64 &curTime, const C_FLOAT64 &endTime)
unsigned C_INT32 mMaxSteps
const C_FLOAT64 & getTime() const
Definition: CState.cpp:325
#define C_FLOAT64
Definition: copasi.h:92
const CState & getState() const
Definition: CModel.cpp:1771
CModel * getModel() const

Friends And Related Function Documentation

Member Data Documentation

C_FLOAT64 CStochDirectMethod::mA0
protected

Total propensity (sum over mAmu[i])

Definition at line 237 of file CStochDirectMethod.h.

Referenced by doSingleStep(), and start().

CVector< C_FLOAT64 > CStochDirectMethod::mAmu
protected

A vector of reaction propensities

Definition at line 232 of file CStochDirectMethod.h.

Referenced by calculateAmu(), doSingleStep(), and start().

bool CStochDirectMethod::mDoCorrection
protected

A boolean flag indicating whether correction for higher order reactions need to be applied

Definition at line 227 of file CStochDirectMethod.h.

Referenced by calculateAmu(), and start().

unsigned C_INT32 CStochDirectMethod::mMaxSteps
protected

max number of single stochastic steps to do in one step()

Definition at line 212 of file CStochDirectMethod.h.

Referenced by start(), and step().

bool CStochDirectMethod::mMaxStepsReached
protected

A boolean flag indicating whether the maximum steps have been reached. This is used to avoid multiple messages.

Definition at line 253 of file CStochDirectMethod.h.

Referenced by start().

CState CStochDirectMethod::mMethodState
protected

The method internal state which contains particle rounded particle numbers.

Definition at line 242 of file CStochDirectMethod.h.

Referenced by start(), and step().

size_t CStochDirectMethod::mNextReactionIndex
protected

The index of the next reaction which fires

Definition at line 222 of file CStochDirectMethod.h.

Referenced by doSingleStep(), and start().

C_FLOAT64 CStochDirectMethod::mNextReactionTime
protected

The time the next reaction fires

Definition at line 217 of file CStochDirectMethod.h.

Referenced by doSingleStep(), and start().

size_t CStochDirectMethod::mNumReactions
protected

The particle and reaction numbers

Definition at line 207 of file CStochDirectMethod.h.

Referenced by doSingleStep(), and start().

CModel* CStochDirectMethod::mpModel
protected

A pointer to the instance of CModel being used.

Definition at line 202 of file CStochDirectMethod.h.

Referenced by doSingleStep(), and start().

CRandom* CStochDirectMethod::mpRandomGenerator
protected

The random number generator

Definition at line 197 of file CStochDirectMethod.h.

Referenced by doSingleStep(), start(), and ~CStochDirectMethod().

std::vector< CReactionDependencies > CStochDirectMethod::mReactionDependencies
protected

A vector containing dependency information to minimize the required updates.

Definition at line 247 of file CStochDirectMethod.h.

Referenced by calculateAmu(), doSingleStep(), and start().


The documentation for this class was generated from the following files: