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

#include <CTauLeapMethod.h>

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

Classes

class  CReactionDependencies
 

Public Member Functions

 CTauLeapMethod (const CTauLeapMethod &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)
 
 ~CTauLeapMethod ()
 
- 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

const C_FLOAT64calculateAmu (const size_t &index)
 
void cleanup ()
 
 CTauLeapMethod (const CCopasiContainer *pParent=NULL)
 
C_FLOAT64 doSingleStep (C_FLOAT64 ds)
 
void updatePropensities ()
 
bool updateSystem ()
 
- 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
 
CVector< C_FLOAT64mAvgDX
 
bool mDoCorrection
 
C_FLOAT64 mEpsilon
 
size_t mFirstReactionSpeciesIndex
 
CVector< C_FLOAT64mK
 
unsigned C_INT32 mMaxSteps
 
CState mMethodState
 
size_t mNumReactions
 
size_t mNumReactionSpecies
 
CModelmpModel
 
CRandommpRandomGenerator
 
unsigned C_INT32 mRandomSeed
 
std::vector
< CReactionDependencies
mReactionDependencies
 
CVector< C_FLOAT64mSigDX
 
bool mUseRandomSeed
 
- 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 47 of file CTauLeapMethod.h.

Constructor & Destructor Documentation

CTauLeapMethod::CTauLeapMethod ( const CTauLeapMethod src,
const CCopasiContainer pParent = NULL 
)

Copy constructor

Parameters
constCTauLeapMethod & src
constCCopasiContainer * pParent (default: NULL)

Definition at line 112 of file CTauLeapMethod.cpp.

References CRandom::createGenerator(), initializeParameter(), mpRandomGenerator, and CRandom::mt19937.

113  :
114  CTrajectoryMethod(src, pParent),
116 {
119 }
void initializeParameter()
std::vector< CReactionDependencies > mReactionDependencies
static CRandom * createGenerator(CRandom::Type type=CRandom::mt19937, unsigned C_INT32 seed=0)
Definition: CRandom.cpp:49
CRandom * mpRandomGenerator
CTauLeapMethod::~CTauLeapMethod ( )

Destructor.

Definition at line 124 of file CTauLeapMethod.cpp.

References cleanup(), and DESTRUCTOR_TRACE.

125 {
126  cleanup();
128 }
#define DESTRUCTOR_TRACE
Definition: copasi.h:206
CTauLeapMethod::CTauLeapMethod ( const CCopasiContainer pParent = NULL)
protected

Default constructor.

Parameters
constCCopasiContainer * pParent (default: NULL)

Default constructor.

Definition at line 104 of file CTauLeapMethod.cpp.

References CRandom::createGenerator(), initializeParameter(), mpRandomGenerator, and CRandom::mt19937.

104  :
107 {
110 }
void initializeParameter()
std::vector< CReactionDependencies > mReactionDependencies
static CRandom * createGenerator(CRandom::Type type=CRandom::mt19937, unsigned C_INT32 seed=0)
Definition: CRandom.cpp:49
CRandom * mpRandomGenerator

Member Function Documentation

const C_FLOAT64 & CTauLeapMethod::calculateAmu ( const size_t &  index)
protected

Calculate one of the propensities

Parameters
constsize_t & index
Returns
const C_FLOAT64 & amu

Definition at line 473 of file CTauLeapMethod.cpp.

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

Referenced by updatePropensities().

474 {
475  const CReactionDependencies & Dependencies = mReactionDependencies[index];
476  C_FLOAT64 & Amu = mAmu[index];
477 
478  Amu = *Dependencies.mpParticleFlux;
479 
480  if (Amu < 0.0)
481  {
482  // TODO CRITICAL Create a warning message
483  Amu = 0.0;
484  }
485 
486  if (!mDoCorrection)
487  {
488  return Amu;
489  }
490 
491  C_FLOAT64 SubstrateMultiplier = 1.0;
492  C_FLOAT64 SubstrateDevisor = 1.0;
493  C_FLOAT64 Multiplicity;
494  C_FLOAT64 LowerBound;
495  C_FLOAT64 Number;
496 
497  bool ApplyCorrection = false;
498 
499  const C_FLOAT64 * pMultiplicity = Dependencies.mSubstrateMultiplier.array();
500  const C_FLOAT64 * pEndMultiplicity = pMultiplicity + Dependencies.mSubstrateMultiplier.size();
501  C_FLOAT64 *const* ppLocalSubstrate = Dependencies.mMethodSubstrates.array();
502  C_FLOAT64 *const* ppModelSubstrate = Dependencies.mModelSubstrates.array();
503 
504  for (; pMultiplicity != pEndMultiplicity; ++pMultiplicity, ++ppLocalSubstrate, ++ppModelSubstrate)
505  {
506  Multiplicity = *pMultiplicity;
507 
508  // TODO We should check the error introduced through rounding.
509  **ppLocalSubstrate = floor(**ppModelSubstrate + 0.5);
510 
511  if (Multiplicity > 1.01)
512  {
513  ApplyCorrection = true;
514 
515  Number = **ppLocalSubstrate;
516 
517  LowerBound = Number - Multiplicity;
518  SubstrateDevisor *= pow(Number, Multiplicity - 1.0); //optimization
519  Number -= 1.0;
520 
521  while (Number > LowerBound)
522  {
523  SubstrateMultiplier *= Number;
524  Number -= 1.0;
525  }
526  }
527  }
528 
529  // at least one substrate particle number is zero
530  if (SubstrateMultiplier < 0.5 || SubstrateDevisor < 0.5)
531  {
532  Amu = 0.0;
533  }
534  else if (ApplyCorrection)
535  {
536  Amu *= SubstrateMultiplier / SubstrateDevisor;
537  }
538 
539  return Amu;
540 }
std::vector< CReactionDependencies > mReactionDependencies
CVector< C_FLOAT64 > mAmu
#define C_FLOAT64
Definition: copasi.h:92
void CTauLeapMethod::cleanup ( )
protected

Cleans up memory, etc.

Initializes the solver and sets the model to be used.

Parameters
modelA reference to an instance of a CModel Cleans up memory, etc.

Definition at line 357 of file CTauLeapMethod.cpp.

References mpModel, and mpRandomGenerator.

Referenced by ~CTauLeapMethod().

358 {
359  delete mpRandomGenerator;
360  mpRandomGenerator = NULL;
361  mpModel = NULL;
362  return;
363 }
CRandom * mpRandomGenerator
C_FLOAT64 CTauLeapMethod::doSingleStep ( C_FLOAT64  ds)
protected

Simulates the system over the next interval of time. The timestep is given as argument.

Parameters
dsA C_FLOAT64 specifying the timestep

Definition at line 371 of file CTauLeapMethod.cpp.

References CVectorCore< CType >::array(), CState::beginIndependent(), C_FLOAT64, CCopasiMessage::EXCEPTION, CRandom::getRandomCC(), CRandom::getRandomPoisson(), mAmu, mAvgDX, MCTrajectoryMethod, mEpsilon, mFirstReactionSpeciesIndex, min, mK, mMethodState, mNumReactions, mNumReactionSpecies, mpRandomGenerator, mReactionDependencies, mSigDX, updatePropensities(), and updateSystem().

Referenced by step().

372 {
373  C_FLOAT64 Lambda, Tmp, Tau, Tau1, Tau2;
374 
376 
377  mAvgDX = 0.0;
378  mSigDX = 0.0;
379 
380  std::vector< CReactionDependencies >::const_iterator itReaction = mReactionDependencies.begin();
381  const C_FLOAT64 * pAmu = mAmu.array();
382  const C_FLOAT64 * pAmuEnd = pAmu + mNumReactions;
383 
384  for (; pAmu != pAmuEnd; ++pAmu, ++itReaction)
385  {
386  const C_FLOAT64 * pMultiplicity = itReaction->mSpeciesMultiplier.array();
387  const C_FLOAT64 * pMultiplicityEnd = pMultiplicity + itReaction->mSpeciesMultiplier.size();
388  const size_t * pIndex = itReaction->mMethodSpeciesIndex.array();
389 
390  for (; pMultiplicity != pMultiplicityEnd; ++pMultiplicity, ++pIndex)
391  {
392  mAvgDX[*pIndex] += *pMultiplicity * *pAmu;
393  mSigDX[*pIndex] += *pMultiplicity * *pMultiplicity * *pAmu;
394  }
395  }
396 
397  Tau1 = Tau2 = std::numeric_limits< C_FLOAT64 >::infinity();
398 
400  const C_FLOAT64 * pNumberEnd = pNumber + mNumReactionSpecies;
401  C_FLOAT64 * pAvgDX = mAvgDX.array();
402  C_FLOAT64 * pSigDX = mSigDX.array();
403 
404  for (; pNumber != pNumberEnd; ++pNumber, ++pAvgDX, ++pSigDX)
405  {
406  if ((Tmp = mEpsilon * fabs(*pNumber)) < 1.0)
407  Tmp = 1.0;
408 
409  *pAvgDX = Tmp / fabs(*pAvgDX);
410  *pSigDX = (Tmp * Tmp) / fabs(*pSigDX);
411 
412  if (Tau1 > *pAvgDX)
413  Tau1 = *pAvgDX;
414 
415  if (Tau2 > *pSigDX)
416  Tau2 = *pSigDX;
417  }
418 
419  Tau = std::min(Tau1, Tau2);
420 
421  if (ds < Tau)
422  Tau = ds;
423 
424  pAmu = mAmu.array();
425  C_FLOAT64 * pK = mK.array();
426  C_FLOAT64 * pKEnd = pK + mNumReactions;
427 
428  for (; pAmu != pAmuEnd; ++pAmu, ++pK)
429  {
430  Lambda = *pAmu * Tau;
431 
432  if (Lambda < 0.0)
434  else if (Lambda > 2.0e9)
436 
437  *pK = mpRandomGenerator->getRandomPoisson(Lambda);
438  }
439 
440  while (!updateSystem())
441  {
442  Tau *= 0.5;
443  pK = mK.array();
444 
445  for (; pK != pKEnd; ++pK)
446  {
447  *pK *= 0.5;
448 
449  if (*pK < floor(*pK + 0.75))
450  {
451  *pK += mpRandomGenerator->getRandomCC() < 0.5 ? - 0.5 : 0.5;
452  }
453  }
454  }
455 
456 
457  return Tau;
458 }
C_FLOAT64 mEpsilon
CVector< C_FLOAT64 > mK
size_t mFirstReactionSpeciesIndex
#define MCTrajectoryMethod
std::vector< CReactionDependencies > mReactionDependencies
virtual C_FLOAT64 getRandomPoisson(const C_FLOAT64 &mean)
Definition: CRandom.cpp:314
virtual C_FLOAT64 getRandomCC()
Definition: CRandom.cpp:235
CVector< C_FLOAT64 > mAmu
CRandom * mpRandomGenerator
CVector< C_FLOAT64 > mSigDX
#define C_FLOAT64
Definition: copasi.h:92
CType * array()
Definition: CVector.h:139
size_t mNumReactionSpecies
CVector< C_FLOAT64 > mAvgDX
#define min(a, b)
Definition: f2c.h:175
C_FLOAT64 * beginIndependent()
Definition: CState.cpp:328
bool CTauLeapMethod::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 158 of file CTauLeapMethod.cpp.

References initializeParameter().

159 {
161  return true;
162 }
void initializeParameter()
void CTauLeapMethod::initializeParameter ( )
private

Initialize the method parameter

Definition at line 130 of file CTauLeapMethod.cpp.

References CCopasiParameterGroup::assertParameter(), CCopasiParameter::BOOL, C_FLOAT64, C_INT32, CCopasiParameter::DOUBLE, CCopasiParameterGroup::getParameter(), CCopasiParameter::getValue(), CCopasiParameter::Value::pBOOL, CCopasiParameter::Value::pUINT, CCopasiParameterGroup::removeParameter(), CCopasiParameterGroup::setValue(), and CCopasiParameter::UINT.

Referenced by CTauLeapMethod(), and elevateChildren().

131 {
132  CCopasiParameter *pParm;
133 
135  assertParameter("Max Internal Steps", CCopasiParameter::UINT, (unsigned C_INT32) 10000);
136  assertParameter("Use Random Seed", CCopasiParameter::BOOL, false);
137  assertParameter("Random Seed", CCopasiParameter::UINT, (unsigned C_INT32) 1);
138 
139  // Check whether we have a method with the old parameter names
140  if ((pParm = getParameter("TAULEAP.Tau")) != NULL)
141  {
142  removeParameter("TAULEAP.Tau");
143 
144  if ((pParm = getParameter("TAULEAP.UseRandomSeed")) != NULL)
145  {
146  setValue("Use Random Seed", *pParm->getValue().pBOOL);
147  removeParameter("TAULEAP.UseRandomSeed");
148  }
149 
150  if ((pParm = getParameter("TAULEAP.RandomSeed")) != NULL)
151  {
152  setValue("Random Seed", *pParm->getValue().pUINT);
153  removeParameter("TAULEAP.RandomSeed");
154  }
155  }
156 }
#define C_INT32
Definition: copasi.h:90
bool removeParameter(const std::string &name)
const Value & getValue() const
bool setValue(const std::string &name, const CType &value)
unsigned C_INT32 * pUINT
CCopasiParameter * getParameter(const std::string &name)
#define C_FLOAT64
Definition: copasi.h:92
CCopasiParameter * assertParameter(const std::string &name, const CCopasiParameter::Type type, const CType &defaultValue)
bool CTauLeapMethod::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 585 of file CTauLeapMethod.cpp.

References CCopasiMessage::ERROR, CModel::getCompartments(), CTrajectoryProblem::getDuration(), CModel::getEvents(), CModel::getMetabolites(), CCopasiProblem::getModel(), CModel::getModelValues(), CModel::getNumMetabs(), CModel::getNumModelValues(), CModel::getTotSteps(), CTrajectoryMethod::isValidProblem(), MCTrajectoryMethod, CModelEntity::ODE, CCopasiVector< T >::size(), and CModel::suitableForStochasticSimulation().

586 {
587  if (!CTrajectoryMethod::isValidProblem(pProblem)) return false;
588 
589  const CTrajectoryProblem * pTP = dynamic_cast<const CTrajectoryProblem *>(pProblem);
590 
591  if (pTP->getDuration() < 0.0)
592  {
593  //back integration not possible
595  return false;
596  }
597 
598  if (pTP->getModel()->getTotSteps() < 1)
599  {
600  //at least one reaction necessary
602  return false;
603  }
604 
605  //check for ODE rules
606  size_t i, imax = pTP->getModel()->getNumModelValues();
607 
608  for (i = 0; i < imax; ++i)
609  {
610  if (pTP->getModel()->getModelValues()[i]->getStatus() == CModelEntity::ODE)
611  {
612  //ode rule found
614  return false;
615  }
616  }
617 
618  imax = pTP->getModel()->getNumMetabs();
619 
620  for (i = 0; i < imax; ++i)
621  {
622  if (pTP->getModel()->getMetabolites()[i]->getStatus() == CModelEntity::ODE)
623  {
624  //ode rule found
626  return false;
627  }
628  }
629 
630  imax = pTP->getModel()->getCompartments().size();
631 
632  for (i = 0; i < imax; ++i)
633  {
634  if (pTP->getModel()->getCompartments()[i]->getStatus() == CModelEntity::ODE)
635  {
636  //ode rule found
638  return false;
639  }
640  }
641 
642  //TODO: rewrite CModel::suitableForStochasticSimulation() to use
643  // CCopasiMessage
644  std::string message = pTP->getModel()->suitableForStochasticSimulation();
645 
646  if (message != "")
647  {
648  //model not suitable, message describes the problem
649  CCopasiMessage(CCopasiMessage::ERROR, message.c_str());
650  return false;
651  }
652 
653  //events are not supported at the moment
654  if (pTP->getModel()->getEvents().size() > 0)
655  {
657  return false;
658  }
659 
660  return true;
661 }
CCopasiVectorN< CEvent > & getEvents()
Definition: CModel.cpp:1110
const CCopasiVector< CMetab > & getMetabolites() const
Definition: CModel.cpp:1051
const CCopasiVectorN< CModelValue > & getModelValues() const
Definition: CModel.cpp:1060
virtual size_t size() const
size_t getNumMetabs() const
Definition: CModel.cpp:1118
virtual bool isValidProblem(const CCopasiProblem *pProblem)
size_t getTotSteps() const
Definition: CModel.cpp:1136
#define MCTrajectoryMethod
const C_FLOAT64 & getDuration() const
CCopasiVectorNS< CCompartment > & getCompartments()
Definition: CModel.cpp:1145
size_t getNumModelValues() const
Definition: CModel.cpp:1139
std::string suitableForStochasticSimulation() const
Definition: CModel.cpp:3901
CModel * getModel() const
void CTauLeapMethod::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 195 of file CTauLeapMethod.cpp.

References CCopasiVector< T >::begin(), CStateTemplate::beginIndependent(), CState::beginIndependent(), C_FLOAT64, C_INT32, CModel::deterministic, CCopasiVector< T >::end(), CStateTemplate::endFixed(), CStateTemplate::getIndex(), CCopasiProblem::getModel(), CModel::getModelType(), CModel::getNumDependentReactionMetabs(), CModel::getNumIndependentReactionMetabs(), CModel::getReactions(), CModel::getState(), CModel::getStateTemplate(), CModelEntity::getStatus(), CCopasiParameter::getValue(), CCopasiObject::getValuePointer(), CModelEntity::getValueReference(), CRandom::initialize(), mAmu, mAvgDX, mDoCorrection, mEpsilon, mFirstReactionSpeciesIndex, mK, mMaxSteps, mMethodState, mNumReactions, mNumReactionSpecies, CTrajectoryMethod::mpCurrentState, mpModel, CTrajectoryMethod::mpProblem, mpRandomGenerator, mRandomSeed, mReactionDependencies, mSigDX, mUseRandomSeed, CCopasiParameter::Value::pBOOL, CCopasiParameter::Value::pDOUBLE, CCopasiParameter::Value::pUINT, CModelEntity::REACTIONS, CVector< CType >::resize(), CModel::setState(), CCopasiVector< T >::size(), and CModel::updateSimulatedValues().

196 {
197  /* get configuration data */
198 
199  bool useRandomSeed = * getValue("Use Random Seed").pBOOL;
200  unsigned C_INT32 randomSeed = * getValue("Random Seed").pUINT;
201 
202  if (useRandomSeed) mpRandomGenerator->initialize(randomSeed);
203 
204  mEpsilon = * getValue("Epsilon").pDOUBLE;
205  mUseRandomSeed = * getValue("Use Random Seed").pBOOL;
206  mRandomSeed = * getValue("Random Seed").pUINT;
207  mMaxSteps = * getValue("Max Internal Steps").pUINT;
208 
209  *mpCurrentState = *initialState;
210 
212  assert(mpModel);
213 
215  mDoCorrection = true;
216  else
217  mDoCorrection = false;
218 
219  // Size the arrays
221 
225  mAmu = 0.0;
226 
228 
231 
232  // Create a local copy of the state where the particle number are rounded to integers
234 
235  const CStateTemplate & StateTemplate = mpModel->getStateTemplate();
236 
237  CModelEntity *const* ppEntity = StateTemplate.beginIndependent();
238  CModelEntity *const* endEntity = StateTemplate.endFixed();
240 
242  size_t Index = 1;
243 
244  for (; ppEntity != endEntity; ++ppEntity, ++pValue, ++Index)
245  {
246  if (dynamic_cast< const CMetab * >(*ppEntity) != NULL)
247  {
248  *pValue = floor(*pValue + 0.5);
249 
250  if (mFirstReactionSpeciesIndex == 0 &&
251  (*ppEntity)->getStatus() == CModelEntity::REACTIONS)
252  {
254  }
255  }
256  }
257 
258  // Update the model state so that the species are all represented by integers.
260  mpModel->updateSimulatedValues(false); //for assignments
261 
262  C_FLOAT64 * pMethodStateValue = mMethodState.beginIndependent() - 1;
263 
264  // Build the reaction dependencies
265  size_t NumReactions = 0;
266 
269  std::vector< CReactionDependencies >::iterator itDependencies = mReactionDependencies.begin();
270 
271  for (; it != end; ++it)
272  {
273  const CCopasiVector<CChemEqElement> & Balances = (*it)->getChemEq().getBalances();
274  const CCopasiVector<CChemEqElement> & Substrates = (*it)->getChemEq().getSubstrates();
275 
276  // This reactions does not change anything we ignore it
277  if (Balances.size() == 0 && Substrates.size() == 0)
278  {
279  continue;
280  }
281 
282  itDependencies->mpParticleFlux = (C_FLOAT64 *)(*it)->getParticleFluxReference()->getValuePointer();
283 
284  itDependencies->mMethodSpeciesIndex.resize(Balances.size());
285  itDependencies->mSpeciesMultiplier.resize(Balances.size());
286  itDependencies->mMethodSpecies.resize(Balances.size());
287  itDependencies->mModelSpecies.resize(Balances.size());
288 
291 
292  Index = 0;
293 
294  for (; itBalance != endBalance; ++itBalance)
295  {
296  const CMetab * pMetab = (*itBalance)->getMetabolite();
297 
298  if (pMetab->getStatus() == CModelEntity::REACTIONS)
299  {
300  itDependencies->mMethodSpeciesIndex[Index] = StateTemplate.getIndex(pMetab) - mFirstReactionSpeciesIndex;
301  itDependencies->mSpeciesMultiplier[Index] = floor((*itBalance)->getMultiplicity() + 0.5);
302  itDependencies->mMethodSpecies[Index] = pMethodStateValue + StateTemplate.getIndex(pMetab);
303  itDependencies->mModelSpecies[Index] = (C_FLOAT64 *) pMetab->getValueReference()->getValuePointer();
304 
305  Index++;
306  }
307  }
308 
309  // Correct allocation for metabolites which are not determined by reactions
310  itDependencies->mMethodSpeciesIndex.resize(Index, true);
311  itDependencies->mSpeciesMultiplier.resize(Index, true);
312  itDependencies->mMethodSpecies.resize(Index, true);
313  itDependencies->mModelSpecies.resize(Index, true);
314 
315  itDependencies->mSubstrateMultiplier.resize(Substrates.size());
316  itDependencies->mMethodSubstrates.resize(Substrates.size());
317  itDependencies->mModelSubstrates.resize(Substrates.size());
318 
319  CCopasiVector< CChemEqElement >::const_iterator itSubstrate = Substrates.begin();
320  CCopasiVector< CChemEqElement >::const_iterator endSubstrate = Substrates.end();
321 
322  Index = 0;
323 
324  for (; itSubstrate != endSubstrate; ++itSubstrate, ++Index)
325  {
326  const CMetab * pMetab = (*itSubstrate)->getMetabolite();
327 
328  itDependencies->mSubstrateMultiplier[Index] = floor((*itSubstrate)->getMultiplicity() + 0.5);
329  itDependencies->mMethodSubstrates[Index] = pMethodStateValue + StateTemplate.getIndex(pMetab);
330  itDependencies->mModelSubstrates[Index] = (C_FLOAT64 *) pMetab->getValueReference()->getValuePointer();
331  }
332 
333  ++itDependencies;
334  ++NumReactions;
335  }
336 
337  mNumReactions = NumReactions;
338 
340  mAmu.resize(mNumReactions, true);
341  mK.resize(mNumReactions, true);
342 
343  return;
344 }
C_FLOAT64 mEpsilon
virtual size_t size() const
iterator begin()
void updateSimulatedValues(const bool &updateMoieties)
Definition: CModel.cpp:1851
CVector< C_FLOAT64 > mK
CTrajectoryProblem * mpProblem
void resize(size_t size, const bool &copy=false)
Definition: CVector.h:301
size_t mFirstReactionSpeciesIndex
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
std::vector< CReactionDependencies > mReactionDependencies
size_t getNumDependentReactionMetabs() const
Definition: CModel.cpp:1133
CModelEntity ** endFixed()
Definition: CState.cpp:213
void setState(const CState &state)
Definition: CModel.cpp:1785
unsigned C_INT32 mMaxSteps
const Value & getValue() const
CVector< C_FLOAT64 > mAmu
unsigned C_INT32 * pUINT
CModelEntity ** beginIndependent()
Definition: CState.cpp:208
size_t getNumIndependentReactionMetabs() const
Definition: CModel.cpp:1130
CRandom * mpRandomGenerator
CVector< C_FLOAT64 > mSigDX
#define C_FLOAT64
Definition: copasi.h:92
const CStateTemplate & getStateTemplate() const
Definition: CModel.cpp:1172
const ModelType & getModelType() const
Definition: CModel.cpp:2339
size_t mNumReactionSpecies
CCopasiVectorNS< CReaction > & getReactions()
Definition: CModel.cpp:1039
virtual void * getValuePointer() const
CVector< C_FLOAT64 > mAvgDX
unsigned C_INT32 mRandomSeed
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 CTauLeapMethod::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 164 of file CTauLeapMethod.cpp.

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

165 {
166  // do several steps
167  C_FLOAT64 Time = mpCurrentState->getTime();
168  C_FLOAT64 EndTime = Time + deltaT;
169 
170  size_t Steps = 0;
171 
172  while (Time < EndTime)
173  {
174  mMethodState.setTime(Time);
177 
178  // We do not need to update the the method state since the only independent state
179  // values are species of type reaction which are all controlled by the method.
180 
181  Time += doSingleStep(EndTime - Time);
182 
183  if (++Steps > mMaxSteps)
184  {
186  }
187  }
188 
190  mpCurrentState->setTime(Time);
191 
192  return NORMAL;
193 }
void updateSimulatedValues(const bool &updateMoieties)
Definition: CModel.cpp:1851
CTrajectoryProblem * mpProblem
C_FLOAT64 doSingleStep(C_FLOAT64 ds)
void setTime(const C_FLOAT64 &time)
Definition: CState.cpp:326
#define MCTrajectoryMethod
void setState(const CState &state)
Definition: CModel.cpp:1785
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
void CTauLeapMethod::updatePropensities ( )
protected

Calculate the propensities for all reactions

Definition at line 460 of file CTauLeapMethod.cpp.

References calculateAmu(), mA0, and mNumReactions.

Referenced by doSingleStep().

461 {
462  //mA0Old = mA0;
463  mA0 = 0;
464 
465  for (size_t i = 0; i < mNumReactions; i++)
466  {
467  mA0 += calculateAmu(i);
468  }
469 
470  return;
471 }
const C_FLOAT64 & calculateAmu(const size_t &index)
bool CTauLeapMethod::updateSystem ( )
protected

Updates the system according to the probabilistic number of firings mK[i] of each reaction i

Definition at line 546 of file CTauLeapMethod.cpp.

References CVectorCore< CType >::array(), CState::beginIndependent(), C_FLOAT64, mFirstReactionSpeciesIndex, mK, mMethodState, mNumReactions, mNumReactionSpecies, and mReactionDependencies.

Referenced by doSingleStep().

547 {
548  std::vector< CReactionDependencies >::const_iterator itReaction = mReactionDependencies.begin();
549 
550  CState OldState(mMethodState);
551 
552  const C_FLOAT64 * pK = mK.array();
553  const C_FLOAT64 * pKEnd = pK + mNumReactions;
554 
555  for (; pK != pKEnd; ++pK, ++itReaction)
556  {
557  const C_FLOAT64 * pMultiplicity = itReaction->mSpeciesMultiplier.array();
558  const C_FLOAT64 * pMultiplicityEnd = pMultiplicity + itReaction->mSpeciesMultiplier.size();
559  C_FLOAT64 * const * pSpecies = itReaction->mMethodSpecies.array();
560 
561  for (; pMultiplicity != pMultiplicityEnd; ++pMultiplicity, ++pSpecies)
562  {
563  **pSpecies += *pK * *pMultiplicity;
564  }
565  }
566 
568 
569  const C_FLOAT64 * pSpeciesEnd = pSpecies + mNumReactionSpecies;
570 
571  for (; pSpecies != pSpeciesEnd; ++pSpecies)
572  {
573  if (*pSpecies < -0.5)
574  {
575  // We need to undo the changes
576  mMethodState = OldState;
577  return false;
578  }
579  }
580 
581  return true;
582 }
CVector< C_FLOAT64 > mK
Definition: CState.h:305
size_t mFirstReactionSpeciesIndex
std::vector< CReactionDependencies > mReactionDependencies
#define C_FLOAT64
Definition: copasi.h:92
CType * array()
Definition: CVector.h:139
size_t mNumReactionSpecies
C_FLOAT64 * beginIndependent()
Definition: CState.cpp:328

Friends And Related Function Documentation

Member Data Documentation

C_FLOAT64 CTauLeapMethod::mA0
protected

Total propensity (sum over mAmu[i])

Definition at line 252 of file CTauLeapMethod.h.

Referenced by updatePropensities().

CVector< C_FLOAT64 > CTauLeapMethod::mAmu
protected

A vector of reaction propensities

Definition at line 247 of file CTauLeapMethod.h.

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

CVector< C_FLOAT64> CTauLeapMethod::mAvgDX
protected

For tau-selection method

Definition at line 263 of file CTauLeapMethod.h.

Referenced by doSingleStep(), and start().

bool CTauLeapMethod::mDoCorrection
protected

indicates if the correction N^2 -> N*(N-1) should be performed

Definition at line 299 of file CTauLeapMethod.h.

Referenced by calculateAmu(), and start().

C_FLOAT64 CTauLeapMethod::mEpsilon
protected

The tolerance ratio x(t+t')< eps*x(t)

Definition at line 273 of file CTauLeapMethod.h.

Referenced by doSingleStep(), and start().

size_t CTauLeapMethod::mFirstReactionSpeciesIndex
protected

index of first species in a CState

Definition at line 304 of file CTauLeapMethod.h.

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

CVector< C_FLOAT64 > CTauLeapMethod::mK
protected

The k-values of the reactions, that is the probabilistic number of firings within one leap.

Definition at line 258 of file CTauLeapMethod.h.

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

unsigned C_INT32 CTauLeapMethod::mMaxSteps
protected

The maximum number of tau leap steps allowed for an integrations step

Definition at line 278 of file CTauLeapMethod.h.

Referenced by start(), and step().

CState CTauLeapMethod::mMethodState
protected

The method internal state which contains particle rounded particle numbers.

Definition at line 227 of file CTauLeapMethod.h.

Referenced by doSingleStep(), start(), step(), and updateSystem().

size_t CTauLeapMethod::mNumReactions
protected

Number of reactions.

Definition at line 232 of file CTauLeapMethod.h.

Referenced by doSingleStep(), start(), updatePropensities(), and updateSystem().

size_t CTauLeapMethod::mNumReactionSpecies
protected

Number of variable metabolites.

Definition at line 242 of file CTauLeapMethod.h.

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

CModel* CTauLeapMethod::mpModel
protected

Pointer to the model.

Definition at line 222 of file CTauLeapMethod.h.

Referenced by cleanup(), start(), and step().

CRandom* CTauLeapMethod::mpRandomGenerator
protected

The random number generator.

Definition at line 294 of file CTauLeapMethod.h.

Referenced by cleanup(), CTauLeapMethod(), doSingleStep(), and start().

unsigned C_INT32 CTauLeapMethod::mRandomSeed
protected

The random seed to use.

Definition at line 289 of file CTauLeapMethod.h.

Referenced by start().

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

A vector containing dependency information to minimize the required updates.

Definition at line 237 of file CTauLeapMethod.h.

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

CVector< C_FLOAT64> CTauLeapMethod::mSigDX
protected

For tau-selection method

Definition at line 268 of file CTauLeapMethod.h.

Referenced by doSingleStep(), and start().

bool CTauLeapMethod::mUseRandomSeed
protected

Specifies if the mRandomSeed should be used. otherwise a randomly chosen seed is used.

Definition at line 284 of file CTauLeapMethod.h.

Referenced by start().


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