#include <CHybridMethod.h>

Inheritance diagram for CHybridMethod:

Collaboration diagram for CHybridMethod:

Public Member Functions
	CHybridMethod (const CHybridMethod &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)

	~CHybridMethod ()

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::SubType &	getSubType () const

const CCopasiTask::Type &	getType () 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)

CCopasiParameterGroup *	assertGroup (const std::string &name)

template<class CType >
CCopasiParameter *	assertParameter (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

CCopasiParameterGroup *	getGroup (const std::string &name)

const CCopasiParameterGroup *	getGroup (const std::string &name) const

CCopasiParameterGroup *	getGroup (const size_t &index)

const CCopasiParameterGroup *	getGroup (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 CObjectInterface *	getObject (const CCopasiObjectName &cn) const

CCopasiParameter *	getParameter (const std::string &name)

const CCopasiParameter *	getParameter (const std::string &name) const

CCopasiParameter *	getParameter (const size_t &index)

const CCopasiParameter *	getParameter (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::Value &	getValue (const std::string &name) const

const CCopasiParameter::Value &	getValue (const size_t &index) const

CCopasiParameter::Value &	getValue (const std::string &name)

CCopasiParameter::Value &	getValue (const size_t &index)

CCopasiParameterGroup &	operator= (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::Type &	getType () const

const Value &	getValue () const

Value &	getValue ()

virtual void *	getValuePointer () const

CCopasiObject *	getValueReference () 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

CCopasiParameter &	operator= (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 objectMap &	getObjects () const

virtual std::string	getUnits () const

virtual const CCopasiObject *	getValueObject () 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 DataObjectSet &	getDirectDependencies (const DataObjectSet &context=DataObjectSet()) const

CCopasiContainer *	getObjectAncestor (const std::string &type) const

CCopasiDataModel *	getObjectDataModel ()

const CCopasiDataModel *	getObjectDataModel () const

const std::string &	getObjectName () const

CCopasiContainer *	getObjectParent () const

const std::string &	getObjectType () const

virtual const CObjectInterface::ObjectSet &	getPrerequisites () const

virtual Refresh *	getRefresh () const

UpdateMethod *	getUpdateMethod () 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 ()

Static Public Member Functions
static CHybridMethod *	createHybridMethod ()

Static Public Member Functions inherited from CTrajectoryMethod
static CTrajectoryMethod *	createMethod (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)

Protected Types
enum	metabStatus { LOW = 0, HIGH }

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 }

Protected Member Functions
void	calculateAmu (size_t rIndex)

void	calculateDerivative (std::vector< C_FLOAT64 > &deriv)

	CHybridMethod (const CCopasiContainer *pParent=NULL)

void	cleanup ()

virtual C_FLOAT64	doSingleStep (C_FLOAT64 currentTime, C_FLOAT64 endTime)=0

void	fireReaction (size_t rIndex)

C_FLOAT64	generateReactionTime (size_t rIndex)

std::set< std::string > *	getAffects (size_t rIndex)

std::set< std::string > *	getDependsOn (size_t rIndex)

std::set< size_t > *	getParticipatesIn (size_t rIndex)

void	getState (std::vector< C_FLOAT64 > &target)

void	getStochTimeAndIndex (C_FLOAT64 &ds, size_t &rIndex)

void	initMethod (C_FLOAT64 time)

void	insertDeterministicReaction (size_t rIndex)

void	integrateDeterministicPart (C_FLOAT64 ds)

void	integrateDeterministicPartEuler (C_FLOAT64 ds)

void	outputData (std::ostream &os, C_INT32 mode)

void	outputDebug (std::ostream &os, size_t level)

void	partitionSystem ()

void	removeDeterministicReaction (size_t rIndex)

void	rungeKutta (C_FLOAT64 dt)

void	setState (std::vector< C_FLOAT64 > &source)

void	setupBalances ()

void	setupDependencyGraph ()

void	setupMetab2React ()

void	setupMetab2ReactComplete ()

void	setupMetab2ReactPlusModifier ()

void	setupPartition ()

void	setupPriorityQueue (C_FLOAT64 startTime=0.0)

void	updatePriorityQueue (size_t rIndex, C_FLOAT64 time)

void	updateTauMu (size_t rIndex, C_FLOAT64 time)

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 >
CCopasiObject *	addMatrixReference (const std::string &name, CType &reference, const unsigned C_INT32 &flag=0)

template<class CType >
CCopasiObject *	addObjectReference (const std::string &name, CType &reference, const unsigned C_INT32 &flag=0)

template<class CType >
CCopasiObject *	addVectorReference (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)

Static Protected Member Functions
static C_INT32	checkModel (CModel *model)

static bool	modelHasAssignments (const CModel *pModel)

Protected Attributes
std::vector< C_FLOAT64 >	currentState

std::vector< C_FLOAT64 >	k1

std::vector< C_FLOAT64 >	k2

std::vector< C_FLOAT64 >	k3

std::vector< C_FLOAT64 >	k4

std::vector< C_FLOAT64 >	mAmu

std::vector< C_FLOAT64 >	mAmuOld

CDependencyGraph	mDG

bool	mDoCorrection

size_t	mFirstMetabIndex

CHybridStochFlag *	mFirstReactionFlag

bool	mHasAssignments

std::vector< std::vector < CHybridBalance > >	mLocalBalances

std::vector< std::vector < CHybridBalance > >	mLocalSubstrates

C_FLOAT64	mLowerStochLimit

size_t	mMaxBalance

size_t	mMaxIntBeforeStep

unsigned C_INT32	mMaxSteps

bool	mMaxStepsReached

std::vector< std::set< size_t > >	mMetab2React

std::vector< metabStatus >	mMetabFlags

size_t	mNumVariableMetabs

size_t	mOutputCounter

std::ofstream	mOutputFile

std::string	mOutputFileName

unsigned C_INT32	mPartitioningInterval

CCopasiVector< CMetab > *	mpMetabolites

CModel *	mpModel

CIndexedPriorityQueue	mPQ

CRandom *	mpRandomGenerator

const CCopasiVectorNS < CReaction > *	mpReactions

unsigned C_INT32	mRandomSeed

std::vector< CHybridStochFlag >	mReactionFlags

size_t	mStepsAfterPartitionSystem

C_FLOAT64	mStepsize

CMatrix< C_FLOAT64 >	mStoi

std::set< size_t >	mUpdateSet

C_FLOAT64	mUpperStochLimit

bool	mUseRandomSeed

std::vector< C_FLOAT64 >	temp

Protected Attributes inherited from CTrajectoryMethod
CState *	mpCurrentState

CTrajectoryProblem *	mpProblem

CVector< C_INT >	mRoots

Protected Attributes inherited from CCopasiMethod
CProcessReport *	mpCallBack

Protected Attributes inherited from CCopasiParameter
std::string	mKey

CCopasiObject *	mpValueReference

size_t	mSize

Value	mValue

Protected Attributes inherited from CCopasiContainer
objectMap	mObjects

Private Member Functions
void	initializeParameter ()

Friends
CTrajectoryMethod *	CTrajectoryMethod::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 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

Static Protected Attributes inherited from CCopasiObject
static CRenameHandler *	smpRenameHandler = NULL

Detailed Description

Definition at line 113 of file CHybridMethod.h.

Member Enumeration Documentation

enum CHybridMethod::metabStatus

protected

Status of the metabolites. Particle number can be high or low.

Enumerator
LOW
HIGH

Definition at line 462 of file CHybridMethod.h.

462 {LOW = 0, HIGH};

CHybridMethod::LOW

Definition: CHybridMethod.h:462

CHybridMethod::HIGH

Definition: CHybridMethod.h:462

Constructor & Destructor Documentation

CHybridMethod::CHybridMethod	(	const CHybridMethod &	src,
		const CCopasiContainer *	pParent = `NULL`
	)

Copy constructor

Parameters

const	CHybridMethod & src
const	CCopasiContainer * pParent (default: NULL)

Definition at line 70 of file CHybridMethod.cpp.

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

                                                               :
   CTrajectoryMethod(src, pParent)
 {
   mpRandomGenerator = CRandom::createGenerator(CRandom::mt19937);
   initializeParameter();
 }

CHybridMethod::~CHybridMethod ( )

Destructor.

Definition at line 81 of file CHybridMethod.cpp.

References cleanup(), and DESTRUCTOR_TRACE.

 {
   cleanup();
   DESTRUCTOR_TRACE;
 }

CHybridMethod::CHybridMethod ( const CCopasiContainer * pParent = NULL )

protected

Default constructor.

Parameters

const CCopasiContainer * pParent (default: NULL)

CHybridMethod

This class implements an hybrid algorithm for the simulation of a biochemical system over time.

File name: CHybridMethod.cpp Author: Juergen Pahle Email: juerg.nosp@m.en.p.nosp@m.ahle@.nosp@m.eml-.nosp@m.r.vil.nosp@m.la-b.nosp@m.osch..nosp@m.de

Last change: 14, December 2004

(C) European Media Lab 2003. Default constructor.

Definition at line 63 of file CHybridMethod.cpp.

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

                                                             :
   CTrajectoryMethod(CCopasiMethod::hybrid, pParent)
 {
   mpRandomGenerator = CRandom::createGenerator(CRandom::mt19937);
   initializeParameter();
 }

Member Function Documentation

void CHybridMethod::calculateAmu ( size_t rIndex )

protected

Calculates an amu value for a given reaction.

Parameters

rIndex A size_t specifying the reaction to be updated

Definition at line 691 of file CHybridMethod.cpp.

References C_FLOAT64, C_INT32, C_INT64, CCopasiParameter::getValue(), mAmu, mDoCorrection, mLocalSubstrates, and mpMetabolites.

Referenced by partitionSystem(), setupPriorityQueue(), and updatePriorityQueue().

 {
   if (!mDoCorrection)
     {
       mAmu[rIndex] = (*mpReactions)[rIndex]->calculateParticleFlux();
       return;
     }
 
   // We need the product of the cmu and hmu for this step.
   // We calculate this in one go, as there are fewer steps to
   // perform and we eliminate some possible rounding errors.
   C_FLOAT64 amu = 1; // initially
   //size_t total_substrates = 0;
   C_INT32 num_ident = 0;
   C_INT64 number = 0;
   C_INT64 lower_bound;
   // substrate_factor - The substrates, raised to their multiplicities,
   // multiplied with one another. If there are, e.g. m substrates of type m,
   // and n of type N, then substrate_factor = M^m * N^n.
   C_FLOAT64 substrate_factor = 1;
   // First, find the reaction associated with this index.
   // Keep a pointer to this.
   // Iterate through each substrate in the reaction
   const std::vector<CHybridBalance> & substrates = mLocalSubstrates[rIndex];
 
   int flag = 0;
 
   for (size_t i = 0; i < substrates.size(); i++)
     {
       num_ident = substrates[i].mMultiplicity;
 
       if (num_ident > 1)
         {
           flag = 1;
           number = static_cast<C_INT64>(floor((*mpMetabolites)[substrates[i].mIndex]->getValue()));
           lower_bound = number - num_ident;
           substrate_factor = substrate_factor * pow((double) number, (int) num_ident - 1); //optimization
 
           number--; // optimization
 
           while (number > lower_bound)
             {
               amu *= number;
               number--;
             }
         }
     }
 
   if ((amu == 0) || (substrate_factor == 0))  // at least one substrate particle number is zero
     {
       mAmu[rIndex] = 0;
       return;
     }
 
   // rate_factor is the rate function divided by substrate_factor.
   // It would be more efficient if this was generated directly, since in effect we
   // are multiplying and then dividing by the same thing (substrate_factor)!
   C_FLOAT64 rate_factor = (*mpReactions)[rIndex]->calculateParticleFlux();
 
   if (flag)
     {
       amu *= rate_factor / substrate_factor;;
       mAmu[rIndex] = amu;
     }
   else
     {
       mAmu[rIndex] = rate_factor;
     }
 
   return;
 
   // a more efficient way to calculate mass action kinetics could be included
 }

void CHybridMethod::calculateDerivative ( std::vector< C_FLOAT64 > & deriv )

protected

Calculates the derivative of the system and writes it into the vector deriv. Length of deriv must be mNumVariableMetabs. CAUTION: Only deterministic reactions are taken into account. That is, this is only the derivative of the deterministic part of the system.

Parameters

deriv A vector reference of length mNumVariableMetabs, into which the derivative is written

Definition at line 489 of file CHybridMethod.cpp.

References C_INT32, mFirstReactionFlag, CHybridStochFlag::mIndex, CHybridStochFlag::mpNext, mStoi, and CMatrix< CType >::numRows().

Referenced by integrateDeterministicPartEuler(), and rungeKutta().

 {
   size_t i;
   C_INT32 bal = 0;
   CHybridStochFlag * j;
 
   // Calculate all the needed kinetic functions of the deterministic reactions
   for (j = mFirstReactionFlag; j != NULL; j = j->mpNext)
     {
       (*mpReactions)[j->mIndex]->calculateParticleFlux();
     }
 
   // For each metabolite add up the contributions of the det. reactions
   // the number of rows in mStoi equals the number of non-fixed metabolites!
   //  for (i=0; i<mNumVariableMetabs; i++)
   for (i = 0; i < mStoi.numRows(); i++)
     {
       deriv[i] = 0.0;
 
       for (j = mFirstReactionFlag; j != NULL; j = j->mpNext)
         {
           // juergen: +0.5 to get a rounding out of the static_cast
           bal = static_cast< C_INT32 >(floor(mStoi[i][j->mIndex] + 0.5));
           deriv[i] += bal * (*mpReactions)[j->mIndex]->getParticleFlux(); //  balance * flux;
         }
     }
 
   /*
     for (; i < mNumVariableMetabs; i++) deriv[i] = 0.0; // important to get a correct deriv vector, because mStoi doesn't cover fixed metabolites
   */
   return;
 }

C_INT32 CHybridMethod::checkModel ( CModel * model )

staticprotected

Test the model if it is proper to perform stochastic simulations on. Several properties are tested (e.g. integer stoichometry, all reactions take place in one compartment only, irreversibility...).

Parameters

model The model to be checked

Returns: 1, if hybrid simulation is possible; <0, if an error occured.

Test the model if it is proper to perform stochastic simulations on. Several properties are tested (e.g. integer stoichometry, all reactions take place in one compartment only, irreversibility...).

Returns: 0, if everything is ok; <0, if an error occured.

Definition at line 804 of file CHybridMethod.cpp.

References C_FLOAT64, C_INT32, CModel::getReactions(), CModel::getStoi(), INT_EPSILON, mpReactions, mStoi, CMatrix< CType >::numRows(), and CCopasiVector< T >::size().

 {
   CCopasiVectorNS <CReaction> * mpReactions = &model->getReactions();
   CMatrix <C_FLOAT64> mStoi = model->getStoi();
   size_t i, numReactions = mpReactions->size();
   size_t j;
   C_INT32 multInt;
 
   C_FLOAT64 multFloat;
   //  size_t metabSize = mpMetabolites->size();
 
   for (i = 0; i < numReactions; i++) // for every reaction
     {
       // TEST getCompartmentNumber() == 1
       if ((*mpReactions)[i]->getCompartmentNumber() != 1) return - 1;
 
       // TEST isReversible() == 0
       if ((*mpReactions)[i]->isReversible() != 0) return - 2;
 
       // TEST integer stoichometry
       // Iterate through each the metabolites
       // juergen: the number of rows of mStoi equals the number of non-fixed metabs!
       //  for (j=0; i<metabSize; j++)
       for (j = 0; j < mStoi.numRows(); j++)
         {
           multFloat = mStoi[j][i];
           multInt = static_cast<C_INT32>(floor(multFloat + 0.5)); // +0.5 to get a rounding out of the static_cast to int!
 
           if ((multFloat - multInt) > INT_EPSILON) return - 3; // INT_EPSILON in CHybridMethod.h
         }
     }
 
   return 1; // Model is appropriate for hybrid simulation
 }

void CHybridMethod::cleanup ( )

protected

Cleans up memory, etc.

Definition at line 319 of file CHybridMethod.cpp.

References mpModel, and mpRandomGenerator.

Referenced by ~CHybridMethod().

 {
   delete mpRandomGenerator;
   mpRandomGenerator = NULL;
   mpModel = NULL;
   return;
 }

CHybridMethod * CHybridMethod::createHybridMethod ( )

static

Creates a HybridMethod adequate for the problem. (only CHybridNextReactionRKMethod so far)

Definition at line 153 of file CHybridMethod.cpp.

References C_INT32.

Referenced by CTrajectoryMethod::createMethod().

 {
   C_INT32 result = 1; // hybrid NextReactionRungeKutta method as default
   /*  if (pProblem && pProblem->getModel())
       {
       result = checkModel(pProblem->getModel());
       }*/
   CHybridMethod * method = NULL;
 
   switch (result)
     {
         /*    case - 3:        // non-integer stoichometry
         CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 1);
         break;
         case - 2:        // reversible reaction exists
         CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 2);
         break;
 
         case - 1:        // more than one compartment involved
         CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 3);
         break;*/
       case 1:
       default:
         // Everything alright: Hybrid simulation possible
         method = new CHybridNextReactionRKMethod();
         break;
     }
 
   return method;
 }

virtual C_FLOAT64 CHybridMethod::doSingleStep	(	C_FLOAT64	currentTime,
		C_FLOAT64	endTime
	)

protectedpure virtual

Simulates the system over the next interval of time. The current time and the end time of the current step() are given as arguments.

Parameters

currentTime	A C_FLOAT64 specifying the current time
endTime	A C_FLOAT64 specifying the end time of the step()

Returns: A C_FLOAT giving the new time

Implemented in CHybridNextReactionRKMethod.

Referenced by step().

bool CHybridMethod::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 143 of file CHybridMethod.cpp.

References initializeParameter().

 {
   initializeParameter();
   return true;
 }

void CHybridMethod::fireReaction ( size_t rIndex )

protected

Executes the specified reaction in the system once.

Parameters

rIndex A size_t specifying the index of the reaction, which will be fired.

Executes the specified reaction in the system once.

Parameters

rIndex	A size_t specifying the index of the reaction, which will be fired.
time	The current time

Definition at line 588 of file CHybridMethod.cpp.

References C_FLOAT64, CDependencyGraph::getDependents(), CModelEntity::getValue(), mDG, mLocalBalances, mUpdateSet, CMetab::refreshConcentration(), and CModelEntity::setValue().

Referenced by CHybridNextReactionRKMethod::doSingleStep().

 {
   // Change the particle numbers according to which step took place.
   // First, get the vector of balances in the reaction we've got.
   // (This vector expresses the number change of each metabolite
   // in the reaction.) Then step through each balance, using its
   // multiplicity to calculate a new value for the associated
   // metabolite. Finally, update the metabolite.
 
   size_t i;
   C_FLOAT64 newNumber;
   CMetab * pMetab;
 
   for (i = 0; i < mLocalBalances[rIndex].size(); i++)
     {
       pMetab = mLocalBalances[rIndex][i].mpMetabolite;
       newNumber = pMetab->getValue() + mLocalBalances[rIndex][i].mMultiplicity;
 
       pMetab->setValue(newNumber);
       pMetab->refreshConcentration();
     }
 
   // insert all dependent reactions into the mUpdateSet
   const std::set <size_t> & dependents = mDG.getDependents(rIndex);
   std::copy(dependents.begin(), dependents.end(),
             std::inserter(mUpdateSet, mUpdateSet.begin()));
 
   return;
 }

C_FLOAT64 CHybridMethod::generateReactionTime ( size_t rIndex )

protected

Generates a putative reaction time for the given reaction.

Parameters

rIndex A size_t specifying the index of the reaction

Returns: A C_FLOAT64 holding the calculated reaction time

Definition at line 678 of file CHybridMethod.cpp.

References C_FLOAT64, CRandom::getRandomOO(), mAmu, and mpRandomGenerator.

Referenced by partitionSystem(), setupPriorityQueue(), updatePriorityQueue(), and updateTauMu().

 {
   if (mAmu[rIndex] == 0) return std::numeric_limits<C_FLOAT64>::infinity();
 
   C_FLOAT64 rand2 = mpRandomGenerator->getRandomOO();
   return - 1 * log(rand2) / mAmu[rIndex];
 }

std::set< std::string > * CHybridMethod::getAffects ( size_t rIndex )

protected

Gets the set of metabolites which change number when a given reaction is executed.

Parameters

rIndex The index of the reaction being executed.

Returns: The set of affected metabolites.

Definition at line 1331 of file CHybridMethod.cpp.

References mLocalBalances.

Referenced by setupDependencyGraph().

 {
   size_t i;
   std::set<std::string> *retset = new std::set<std::string>;
 
   // Get the balances  associated with the reaction at this index
   // XXX We first get the chemical equation, then the balances, since the getBalances method in CReaction is unimplemented!
 
   for (i = 0; i < mLocalBalances[rIndex].size(); i++)
     {
       if (mLocalBalances[rIndex][i].mMultiplicity != 0)
         {
           retset->insert(mLocalBalances[rIndex][i].mpMetabolite->getKey());
         }
     }
 
   return retset;
 }

std::set< std::string > * CHybridMethod::getDependsOn ( size_t rIndex )

protected

Gets the set of metabolites on which a given reaction depends.

Parameters

rIndex The index of the reaction being executed.

Returns: The set of metabolites depended on.

Definition at line 1298 of file CHybridMethod.cpp.

References mpReactions, and CFunctionParameter::PARAMETER.

Referenced by setupDependencyGraph().

 {
   std::set<std::string> *retset = new std::set<std::string>;
 
   size_t i, imax = (*mpReactions)[rIndex]->getFunctionParameters().size();
   size_t j, jmax;
 
   for (i = 0; i < imax; ++i)
     {
       if ((*mpReactions)[rIndex]->getFunctionParameters()[i]->getUsage() == CFunctionParameter::PARAMETER)
         continue;
 
       //metablist = (*mpReactions)[rIndex]->getParameterMappingMetab(i);
       const std::vector <std::string> & metabKeylist =
         (*mpReactions)[rIndex]->getParameterMappings()[i];
       jmax = metabKeylist.size();
 
       for (j = 0; j < jmax; ++j)
         {
           retset->insert(metabKeylist[j]);
         }
     }
 
   return retset;
 }

std::set< size_t > * CHybridMethod::getParticipatesIn ( size_t rIndex )

protected

Gets the set of metabolites, which participate in the given reaction either as substrate, product or modifier.

Parameters

rIndex The index of the reaction being executed.

Returns: The set of participating metabolites.

Definition at line 1357 of file CHybridMethod.cpp.

Referenced by setupMetab2ReactPlusModifier().

 {
   std::set<size_t> *retset = new std::set<size_t>;
   return retset;
 }

void CHybridMethod::getState ( std::vector< C_FLOAT64 > & target )

protected

Gathers the state of the system into the array target. Later on CState should be used for this. Length of the array target must be mNumVariableMetabs.

Parameters

target A vector reference of length mNumVariableMetabs, into which the state of the system is written

Gathers the state of the system into the array target. Later on CState should be used for this. Length of the array target must be mNumVariableMetabs.

Parameters

target An array of C_FLOAT64s with length mNumVariableMetabs, into which the state of the system is written

Definition at line 529 of file CHybridMethod.cpp.

References mNumVariableMetabs.

Referenced by integrateDeterministicPartEuler(), and rungeKutta().

 {
   size_t i;
 
   for (i = 0; i < mNumVariableMetabs; i++)
     {
       target[i] = (*mpMetabolites)[i]->getValue();
     }
 
   return;
 }

void CHybridMethod::getStochTimeAndIndex	(	C_FLOAT64 &	ds,
		size_t &	rIndex
	)

protected

Find the reaction index and the reaction time of the stochastic (!) reaction with the lowest reaction time.

Parameters

ds	A reference to a C_FLOAT64. The putative reaction time for the first stochastic reaction is written into this variable.
rIndex	A reference to a size_t. The index of the first stochastic reaction is written into this variable.

Definition at line 574 of file CHybridMethod.cpp.

References mPQ, CIndexedPriorityQueue::topIndex(), and CIndexedPriorityQueue::topKey().

Referenced by CHybridNextReactionRKMethod::doSingleStep().

 {
   ds = mPQ.topKey();
   rIndex = mPQ.topIndex();
   return;
 }

void CHybridMethod::initializeParameter ( )

private

Intialize the method parameter

Definition at line 87 of file CHybridMethod.cpp.

Referenced by CHybridMethod(), and elevateChildren().

 {
   CCopasiParameter *pParm;
 
   assertParameter("Max Internal Steps", CCopasiParameter::INT, (C_INT32) MAX_STEPS);
   assertParameter("Lower Limit", CCopasiParameter::DOUBLE, (C_FLOAT64) LOWER_STOCH_LIMIT);
   assertParameter("Upper Limit", CCopasiParameter::DOUBLE, (C_FLOAT64) UPPER_STOCH_LIMIT);
   assertParameter("Runge Kutta Stepsize", CCopasiParameter::DOUBLE, (C_FLOAT64) RUNGE_KUTTA_STEPSIZE);
   assertParameter("Partitioning Interval", CCopasiParameter::UINT, (unsigned C_INT32) PARTITIONING_INTERVAL);
   assertParameter("Use Random Seed", CCopasiParameter::BOOL, (bool) USE_RANDOM_SEED);
   assertParameter("Random Seed", CCopasiParameter::UINT, (unsigned C_INT32) RANDOM_SEED);
 
   // Check whether we have a method with the old parameter names
   if ((pParm = getParameter("HYBRID.MaxSteps")) != NULL)
     {
       setValue("Max Internal Steps", *pParm->getValue().pINT);
       removeParameter("HYBRID.MaxSteps");
 
       if ((pParm = getParameter("HYBRID.LowerStochLimit")) != NULL)
         {
           setValue("Lower Limit", *pParm->getValue().pDOUBLE);
           removeParameter("HYBRID.LowerStochLimit");
         }
 
       if ((pParm = getParameter("HYBRID.UpperStochLimit")) != NULL)
         {
           setValue("Upper Limit", *pParm->getValue().pDOUBLE);
           removeParameter("HYBRID.UpperStochLimit");
         }
 
       if ((pParm = getParameter("HYBRID.RungeKuttaStepsize")) != NULL)
         {
           setValue("Runge Kutta Stepsize", *pParm->getValue().pDOUBLE);
           removeParameter("HYBRID.RungeKuttaStepsize");
         }
 
       if ((pParm = getParameter("HYBRID.PartitioningInterval")) != NULL)
         {
           setValue("Partitioning Interval", *pParm->getValue().pUINT);
           removeParameter("HYBRID.PartitioningInterval");
         }
 
       if ((pParm = getParameter("UseRandomSeed")) != NULL)
         {
           setValue("Use Random Seed", *pParm->getValue().pBOOL);
           removeParameter("UseRandomSeed");
         }
 
       if ((pParm = getParameter("")) != NULL)
         {
           setValue("Random Seed", *pParm->getValue().pUINT);
           removeParameter("");
         }
     }
 }

void CHybridMethod::initMethod ( C_FLOAT64 start_time )

protected

Initializes the solver.

Parameters

time	the current time

Initializes the solver and sets the model to be used.

Parameters

model A reference to an instance of a CModel

Definition at line 264 of file CHybridMethod.cpp.

References currentState, CModel::getMetabolitesX(), CModel::getNumDependentReactionMetabs(), CModel::getNumIndependentReactionMetabs(), CModel::getReactions(), CModel::getStoi(), CCopasiParameter::getValue(), CRandom::initialize(), k1, k2, k3, k4, mAmu, mAmuOld, mLowerStochLimit, mMaxSteps, mMaxStepsReached, mNumVariableMetabs, mPartitioningInterval, mpMetabolites, mpModel, mpRandomGenerator, mpReactions, mRandomSeed, mStepsAfterPartitionSystem, mStepsize, mStoi, mUpdateSet, mUpperStochLimit, mUseRandomSeed, CCopasiParameter::Value::pBOOL, CCopasiParameter::Value::pDOUBLE, CCopasiParameter::Value::pINT, CCopasiParameter::Value::pUINT, setupBalances(), setupDependencyGraph(), setupMetab2React(), setupPartition(), setupPriorityQueue(), CCopasiVector< T >::size(), and temp.

Referenced by start().

 {
   mpReactions = &mpModel->getReactions();
   mAmu.clear();
   mAmu.resize(mpReactions->size());
   mAmuOld.clear();
   mAmuOld.resize(mpReactions->size());
   mpMetabolites = &(const_cast < CCopasiVector < CMetab > & >(mpModel->getMetabolitesX()));
   //mNumVariableMetabs = mpModel->getNumVariableMetabs(); // ind + dep metabs, without fixed metabs
   mNumVariableMetabs = mpModel->getNumIndependentReactionMetabs() + mpModel->getNumDependentReactionMetabs(); // ind + dep metabs, without fixed metabs
   //  mNumVariableMetabs = mpCurrentState->getNumVariable(); // mpBeginFixed - mpBeginIndependent
 
   temp.clear();
   temp.resize(mNumVariableMetabs);
   currentState.clear();
   currentState.resize(mNumVariableMetabs);
 
   k1.clear();
   k1.resize(mNumVariableMetabs);
   k2.clear();
   k2.resize(mNumVariableMetabs);
   k3.clear();
   k3.resize(mNumVariableMetabs);
   k4.clear();
   k4.resize(mNumVariableMetabs);
 
   /* get configuration data */
   mMaxSteps = * getValue("Max Internal Steps").pINT;
   mLowerStochLimit = * getValue("Lower Limit").pDOUBLE;
   mUpperStochLimit = * getValue("Upper Limit").pDOUBLE;
   mStepsize = * getValue("Runge Kutta Stepsize").pDOUBLE;
   mPartitioningInterval = * getValue("Partitioning Interval").pUINT;
   mUseRandomSeed = * getValue("Use Random Seed").pBOOL;
   mRandomSeed = * getValue("Random Seed").pUINT;
 
   if (mUseRandomSeed) mpRandomGenerator->initialize(mRandomSeed);
 
   mStoi = mpModel->getStoi();
   mStepsAfterPartitionSystem = 0;
   mUpdateSet.clear();
 
   setupBalances(); // initialize mLocalBalances and mLocalSubstrates (has to be called first!)
   setupDependencyGraph(); // initialize mDG
   setupMetab2React(); // initialize mMetab2React
   setupPartition(); // initialize mReactionFlags
   setupPriorityQueue(start_time); // initialize mPQ
 
   mMaxStepsReached = false;
 
   return;
 }

void CHybridMethod::insertDeterministicReaction ( size_t rIndex )

protected

Inserts a new deterministic reaction into the linked list in the vector mReactionFlags.

Parameters

rIndex A size_t giving the index of the reaction to be inserted into the list of deterministic reactions.

Inserts a new deterministic reaction into the linked list in the array mReactionFlags.

Parameters

rIndex A size_t giving the index of the reaction to be inserted into the list of deterministic reactions.

Definition at line 1224 of file CHybridMethod.cpp.

References mAmu, mAmuOld, mFirstReactionFlag, CHybridStochFlag::mpPrev, and mReactionFlags.

Referenced by partitionSystem().

 {
   if (mReactionFlags[rIndex].mpPrev == NULL)
     // reaction is stochastic (avoids double insertions)
     {
       if (mFirstReactionFlag != NULL)
         // there are deterministic reactions already
         {
           mFirstReactionFlag->mpPrev = &mReactionFlags[rIndex];
           mReactionFlags[rIndex].mpNext = mFirstReactionFlag;
           mFirstReactionFlag = &mReactionFlags[rIndex];
           mFirstReactionFlag->mpPrev = mFirstReactionFlag;
         }
       else
         {
           // there are no deterministic reactions
           // Important to distinguish between stochastic (prev == NULL) and deterministic (prev != NULL) reactions
           mReactionFlags[rIndex].mpPrev = &mReactionFlags[rIndex];
           mFirstReactionFlag = &mReactionFlags[rIndex];
         }
 
       mAmu[rIndex] = 0.0;
       mAmuOld[rIndex] = 0.0;
     }
 
   return;
 }

void CHybridMethod::integrateDeterministicPart ( C_FLOAT64 dt )

protected

Integrates the deterministic reactions of the system over the specified time interval.

Parameters

ds	A C_FLOAT64 specifying the stepsize.

Definition at line 335 of file CHybridMethod.cpp.

References C_FLOAT64, CDependencyGraph::getDependents(), mDG, mFirstReactionFlag, CHybridStochFlag::mIndex, CHybridStochFlag::mpNext, mStepsize, mUpdateSet, and rungeKutta().

Referenced by CHybridNextReactionRKMethod::doSingleStep().

 {
   C_FLOAT64 integrationTime = 0.0;
   CHybridStochFlag * react = NULL;
 
   // This method uses a 4th order RungeKutta-method to integrate the deterministic part of the system. Maybe a better numerical method (adaptive stepsize, lsoda, ...) should be introduced here later on
 
   while ((dt - integrationTime) > mStepsize)
     {
       rungeKutta(mStepsize); // for the deterministic part of the system
       integrationTime += mStepsize;
     }
 
   rungeKutta(dt - integrationTime);
 
   // find the set union of all reactions, which depend on one of the deterministic reactions. The propensities of the stochastic reactions in this set union will be updated later in the method updatePriorityQueue().
   for (react = mFirstReactionFlag; react != NULL; react = react->mpNext)
     {
       const std::set <size_t> & dependents = mDG.getDependents(react->mIndex);
       std::copy(dependents.begin(), dependents.end(),
                 std::inserter(mUpdateSet, mUpdateSet.begin()));
     }
 
   return;
 }

void CHybridMethod::integrateDeterministicPartEuler ( C_FLOAT64 dt )

protected

Integrates the deterministic reactions of the system over the specified time interval.

Parameters

ds	A C_FLOAT64 specifying the stepsize.

Definition at line 367 of file CHybridMethod.cpp.

References C_FLOAT64, calculateDerivative(), currentState, CDependencyGraph::getDependents(), getState(), mDG, mFirstReactionFlag, CHybridStochFlag::mIndex, mNumVariableMetabs, CHybridStochFlag::mpNext, mStepsize, mUpdateSet, setState(), and temp.

 {
   C_FLOAT64 integrationTime = 0.0;
   CHybridStochFlag * react = NULL;
   size_t i;
 
   while ((dt - integrationTime) > mStepsize)
     {
       getState(currentState);
       calculateDerivative(temp);
 
       for (i = 0; i < mNumVariableMetabs; i++)
         temp[i] = currentState[i] + (temp[i] * mStepsize);
 
       setState(temp);
       integrationTime += mStepsize;
     }
 
   getState(currentState);
   calculateDerivative(temp);
 
   for (i = 0; i < mNumVariableMetabs; i++)
     temp[i] = currentState[i] + (temp[i] * (dt - integrationTime));
 
   setState(temp);
 
   // find the set union of all reactions, which depend on one of the deterministic reactions. The propensities of the stochastic reactions in this set union will be updated later in the method updatePriorityQueue().
   for (react = mFirstReactionFlag; react != NULL; react = react->mpNext)
     {
       const std::set <size_t> & dependents = mDG.getDependents(react->mIndex);
       std::copy(dependents.begin(), dependents.end(),
                 std::inserter(mUpdateSet, mUpdateSet.begin()));
     }
 
   return;
 }

bool CHybridMethod::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 1714 of file CHybridMethod.cpp.

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

 {
   if (!CTrajectoryMethod::isValidProblem(pProblem)) return false;
 
   const CTrajectoryProblem * pTP = dynamic_cast<const CTrajectoryProblem *>(pProblem);
 
   if (pTP->getDuration() < 0.0)
     {
       //back integration not possible
       CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 9);
       return false;
     }
 
   if (pTP->getModel()->getTotSteps() < 1)
     {
       //at least one reaction necessary
       CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 17);
       return false;
     }
 
   //check for rules
   size_t i, imax = pTP->getModel()->getNumModelValues();
 
   for (i = 0; i < imax; ++i)
     {
       if (pTP->getModel()->getModelValues()[i]->getStatus() == CModelEntity::ODE)
         {
           //ode rule found
           CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 18);
           return false;
         }
 
       /*      if (pTP->getModel()->getModelValues()[i]->getStatus()==CModelEntity::ASSIGNMENT)
               if (pTP->getModel()->getModelValues()[i]->isUsed())
                 {
                   //used assignment found
                   CCopasiMessage(CCopasiMessage::EXCEPTION, MCTrajectoryMethod + 19);
                   return false;
                 }*/
     }
 
   imax = pTP->getModel()->getNumMetabs();
 
   for (i = 0; i < imax; ++i)
     {
       if (pTP->getModel()->getMetabolites()[i]->getStatus() == CModelEntity::ODE)
         {
           //ode rule found
           CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 20);
           return false;
         }
 
       if (pTP->getModel()->getMetabolites()[i]->getStatus() == CModelEntity::ASSIGNMENT)
         if (pTP->getModel()->getMetabolites()[i]->isUsed())
           {
             //used assignment for species found
             CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 24);
             return false;
           }
     }
 
   imax = pTP->getModel()->getCompartments().size();
 
   for (i = 0; i < imax; ++i)
     {
       if (pTP->getModel()->getCompartments()[i]->getStatus() == CModelEntity::ODE)
         {
           //ode rule found
           CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 21);
           return false;
         }
     }
 
   //TODO: rewrite CModel::suitableForStochasticSimulation() to use
   //      CCopasiMessage
   std::string message = pTP->getModel()->suitableForStochasticSimulation();
 
   if (message != "")
     {
       //model not suitable, message describes the problem
       CCopasiMessage(CCopasiMessage::ERROR, message.c_str());
       return false;
     }
 
   /* Max Internal Steps */
   if (* getValue("Max Internal Steps").pINT <= 0)
     {
       //max steps should be at least 1
       CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 15);
       return false;
     }
 
   /* Lower Limit, Upper Limit */
   mLowerStochLimit = * getValue("Lower Limit").pDOUBLE;
   mUpperStochLimit = * getValue("Upper Limit").pDOUBLE;
 
   if (mLowerStochLimit > mUpperStochLimit)
     {
       CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 4, mLowerStochLimit, mUpperStochLimit);
       return false;
     }
 
   /* Runge Kutta Stepsize */
   if (* getValue("Runge Kutta Stepsize").pDOUBLE <= 0.0)
     {
       // Runge Kutta Stepsize must be positive
       CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 13);
       return false;
     }
 
   /* Partitioning Interval */
   // nothing to be done here so far
 
   /* Use Random Seed */
   // should be checked in the widget later on
 
   /* Random Seed */
   // nothing to be done here
 
   //events are not supported at the moment
   if (pTP->getModel()->getEvents().size() > 0)
     {
       CCopasiMessage(CCopasiMessage::ERROR, MCTrajectoryMethod + 23);
       return false;
     }
 
   return true;
 }

bool CHybridMethod::modelHasAssignments ( const CModel * pModel )

staticprotected

tests if the model contains a global value with an assignment rule that is used in calculations

Definition at line 1844 of file CHybridMethod.cpp.

References CModelEntity::ASSIGNMENT, CModel::getCompartments(), CModel::getMetabolites(), CModel::getModelValues(), CModel::getNumMetabs(), CModel::getNumModelValues(), and CCopasiVector< T >::size().

Referenced by start().

 {
   size_t i, imax = pModel->getNumModelValues();
 
   for (i = 0; i < imax; ++i)
     {
       if (pModel->getModelValues()[i]->getStatus() == CModelEntity::ASSIGNMENT)
         if (pModel->getModelValues()[i]->isUsed())
           {
             //used assignment found
             return true;
           }
     }
 
   imax = pModel->getNumMetabs();
 
   for (i = 0; i < imax; ++i)
     {
       if (pModel->getMetabolites()[i]->getStatus() == CModelEntity::ASSIGNMENT)
         if (pModel->getMetabolites()[i]->isUsed())
           {
             //used assignment found
             return true;
           }
     }
 
   imax = pModel->getCompartments().size();
 
   for (i = 0; i < imax; ++i)
     {
       if (pModel->getCompartments()[i]->getStatus() == CModelEntity::ASSIGNMENT)
         if (pModel->getCompartments()[i]->isUsed())
           {
             //used assignment found
             return true;
           }
     }
 
   return false;
 }

void CHybridMethod::outputData	(	std::ostream &	os,
		C_INT32	mode
	)

protected

Prints out data on standard output.

Prints out data on standard output. Deprecated.

Definition at line 1366 of file CHybridMethod.cpp.

References C_INT32, CState::getTime(), mOutputCounter, CTrajectoryMethod::mpCurrentState, mpMetabolites, and CCopasiVector< T >::size().

 {
   static unsigned C_INT32 counter = 0;
   size_t i;
 
   switch (mode)
     {
       case 0:
 
         if (mOutputCounter == (counter++))
           {
             counter = 0;
             os << mpCurrentState->getTime() << " : ";
 
             for (i = 0; i < mpMetabolites->size(); i++)
               {
                 os << (*mpMetabolites)[i]->getValue() << " ";
               }
 
             os << std::endl;
           }
 
         break;
 
       case 1:
         os << mpCurrentState->getTime() << " : ";
 
         for (i = 0; i < mpMetabolites->size(); i++)
           {
             os << (*mpMetabolites)[i]->getValue() << " ";
           }
 
         os << std::endl;
         break;
 
       default:
         ;
     }
 
   return;
 }

void CHybridMethod::outputDebug	(	std::ostream &	os,
		size_t	level
	)

protected

Prints out various data on standard output for debugging purposes.

Definition at line 1411 of file CHybridMethod.cpp.

References currentState, CCopasiObject::getObjectName(), CState::getTime(), k1, k2, k3, k4, LOW, mAmu, mAmuOld, mDG, mFirstReactionFlag, mLocalBalances, mLocalSubstrates, mLowerStochLimit, mMaxBalance, mMaxIntBeforeStep, mMaxSteps, mMetab2React, mMetabFlags, mNumVariableMetabs, mPartitioningInterval, CTrajectoryMethod::mpCurrentState, mpMetabolites, mPQ, mpRandomGenerator, mpReactions, mReactionFlags, mStepsAfterPartitionSystem, mStepsize, mStoi, mUpdateSet, mUpperStochLimit, CMatrix< CType >::numCols(), CMatrix< CType >::numRows(), CCopasiVector< T >::size(), and temp.

 {
   size_t i, j;
   std::set <size_t>::iterator iter, iterEnd;
 
   os << "outputDebug(" << level << ") *********************************************** BEGIN" << std::endl;
 
   switch (level)
     {
       case 0:                              // Everything !!!
         os << " Name: " << CCopasiParameter::getObjectName() << std::endl;
         os << "current time: " << mpCurrentState->getTime() << std::endl;
         os << "mNumVariableMetabs: " << mNumVariableMetabs << std::endl;
         os << "mMaxSteps: " << mMaxSteps << std::endl;
         os << "mMaxBalance: " << mMaxBalance << std::endl;
         os << "mMaxIntBeforeStep: " << mMaxIntBeforeStep << std::endl;
         os << "mpReactions.size(): " << mpReactions->size() << std::endl;
 
         for (i = 0; i < mpReactions->size(); i++)
           os << *(*mpReactions)[i] << std::endl;
 
         os << "mpMetabolites.size(): " << mpMetabolites->size() << std::endl;
 
         for (i = 0; i < mpMetabolites->size(); i++)
           os << *(*mpMetabolites)[i] << std::endl;
 
         os << "mStoi: " << std::endl;
 
         for (i = 0; i < mStoi.numRows(); i++)
           {
             for (j = 0; j < mStoi.numCols(); j++)
               os << mStoi[i][j] << " ";
 
             os << std::endl;
           }
 
         os << "temp: ";
 
         for (i = 0; i < mNumVariableMetabs; i++)
           os << temp[i] << " ";
 
         os << std::endl;
         os << "curentState: ";
 
         for (i = 0; i < mNumVariableMetabs; i++)
           os << currentState[i] << " ";
 
         os << std::endl;
         os << "k1: ";
 
         for (i = 0; i < mNumVariableMetabs; i++)
           os << k1[i] << " ";
 
         os << std::endl;
         os << "k2: ";
 
         for (i = 0; i < mNumVariableMetabs; i++)
           os << k2[i] << " ";
 
         os << std::endl;
         os << "k3: ";
 
         for (i = 0; i < mNumVariableMetabs; i++)
           os << k3[i] << " ";
 
         os << std::endl;
         os << "k4: ";
 
         for (i = 0; i < mNumVariableMetabs; i++)
           os << k4[i] << " ";
 
         os << std::endl;
         os << "mReactionFlags: " << std::endl;
 
         for (i = 0; i < mLocalBalances.size(); i++)
           os << mReactionFlags[i];
 
         os << "mFirstReactionFlag: " << std::endl;
         if (mFirstReactionFlag == NULL) os << "NULL" << std::endl; else os << *mFirstReactionFlag;
 
         os << "mMetabFlags: " << std::endl;
 
         for (i = 0; i < mMetabFlags.size(); i++)
           {
             if (mMetabFlags[i] == LOW)
               os << "LOW ";
             else
               os << "HIGH ";
           }
 
         os << std::endl;
         os << "mLocalBalances: " << std::endl;
 
         for (i = 0; i < mLocalBalances.size(); i++)
           {
             for (j = 0; j < mLocalBalances[i].size(); j++)
               os << mLocalBalances[i][j];
 
             os << std::endl;
           }
 
         os << "mLocalSubstrates: " << std::endl;
 
         for (i = 0; i < mLocalSubstrates.size(); i++)
           {
             for (j = 0; j < mLocalSubstrates[i].size(); j++)
               os << mLocalSubstrates[i][j];
 
             os << std::endl;
           }
 
         os << "mLowerStochLimit: " << mLowerStochLimit << std::endl;
         os << "mUpperStochLimit: " << mUpperStochLimit << std::endl;
         //deprecated:      os << "mOutputCounter: " << mOutputCounter << endl;
         os << "mPartitioningInterval: " << mPartitioningInterval << std::endl;
         os << "mStepsAfterPartitionSystem: " << mStepsAfterPartitionSystem << std::endl;
         os << "mStepsize: " << mStepsize << std::endl;
         os << "mMetab2React: " << std::endl;
 
         for (i = 0; i < mMetab2React.size(); i++)
           {
             os << i << ": ";
 
             for (iter = mMetab2React[i].begin(), iterEnd = mMetab2React[i].end(); iter != iterEnd; ++iter)
               os << *iter << " ";
 
             os << std::endl;
           }
 
         os << "mAmu: " << std::endl;
 
         for (i = 0; i < mpReactions->size(); i++)
           os << mAmu[i] << " ";
 
         os << std::endl;
         os << "mAmuOld: " << std::endl;
 
         for (i = 0; i < mpReactions->size(); i++)
           os << mAmuOld[i] << " ";
 
         os << std::endl;
         os << "mUpdateSet: " << std::endl;
 
         for (iter = mUpdateSet.begin(), iterEnd = mUpdateSet.end(); iter != iterEnd; iter++)
           os << *iter;
 
         os << std::endl;
         os << "mpRandomGenerator: " << mpRandomGenerator << std::endl;
         os << "mDG: " << std::endl << mDG;
         os << "mPQ: " << std::endl << mPQ;
         os << "Particle numbers: " << std::endl;
 
         for (i = 0; i < mpMetabolites->size(); i++)
           {
             os << (*mpMetabolites)[i]->getValue() << " ";
           }
 
         os << std::endl;
         break;
 
       case 1:                               // Variable values only
         os << "current time: " << mpCurrentState->getTime() << std::endl;
         /*
         case 1:
         os << "mTime: " << mpCurrentState->getTime() << std::endl;
         os << "oldState: ";
         for (i = 0; i < mDim; i++)
           os << oldState[i] << " ";
         os << std::endl;
         os << "x: ";
         for (i = 0; i < mDim; i++)
           os << x[i] << " ";
         os << std::endl;
         os << "y: ";
         for (i = 0; i < mDim; i++)
           os << y[i] << " ";
         os << std::endl;
         os << "increment: ";
         for (i = 0; i < mDim; i++)
           os << increment[i] << " ";
         os << std::endl;*/
         os << "temp: ";
 
         for (i = 0; i < mNumVariableMetabs; i++)
           os << temp[i] << " ";
 
         os << std::endl;
         os << "currentState: ";
 
         for (i = 0; i < mNumVariableMetabs; i++)
           os << currentState[i] << " ";
 
         os << std::endl;
         os << "k1: ";
 
         for (i = 0; i < mNumVariableMetabs; i++)
           os << k1[i] << " ";
 
         os << std::endl;
         os << "k2: ";
 
         for (i = 0; i < mNumVariableMetabs; i++)
           os << k2[i] << " ";
 
         os << std::endl;
         os << "k3: ";
 
         for (i = 0; i < mNumVariableMetabs; i++)
           os << k3[i] << " ";
 
         os << std::endl;
         os << "k4: ";
 
         for (i = 0; i < mNumVariableMetabs; i++)
           os << k4[i] << " ";
 
         os << std::endl;
         os << "mReactionFlags: " << std::endl;
 
         for (i = 0; i < mLocalBalances.size(); i++)
           os << mReactionFlags[i];
 
         os << "mFirstReactionFlag: " << std::endl;
         if (mFirstReactionFlag == NULL) os << "NULL" << std::endl; else os << *mFirstReactionFlag;
 
         os << "mMetabFlags: " << std::endl;
 
         for (i = 0; i < mMetabFlags.size(); i++)
           {
             if (mMetabFlags[i] == LOW)
               os << "LOW ";
             else
               os << "HIGH ";
           }
 
         os << std::endl;
         os << "mAmu: " << std::endl;
 
         for (i = 0; i < mpReactions->size(); i++)
           os << mAmu[i] << " ";
 
         os << std::endl;
         os << "mAmuOld: " << std::endl;
 
         for (i = 0; i < mpReactions->size(); i++)
           os << mAmuOld[i] << " ";
 
         os << std::endl;
         os << "mUpdateSet: " << std::endl;
 
         for (iter = mUpdateSet.begin(), iterEnd = mUpdateSet.end(); iter != iterEnd; iter++)
           os << *iter;
 
         os << std::endl;
         os << "mPQ: " << std::endl << mPQ;
         os << "Particle numbers: " << std::endl;
 
         for (i = 0; i < mpMetabolites->size(); i++)
           {
             os << (*mpMetabolites)[i]->getValue() << " ";
           }
 
         os << std::endl;
         break;
 
       case 2:
         break;
 
       default:
         ;
     }
 
   os << "outputDebug(" << level << ") ************************************************* END" << std::endl;
   return;
 }

void CHybridMethod::partitionSystem ( )

protected

Updates the partitioning of the system depending on the particle numbers present.

Definition at line 1162 of file CHybridMethod.cpp.

References C_FLOAT64, calculateAmu(), generateReactionTime(), CState::getTime(), CCopasiParameter::getValue(), HIGH, insertDeterministicReaction(), CIndexedPriorityQueue::insertStochReaction(), LOW, mAmu, mAmuOld, mLowerStochLimit, mMetab2React, mMetabFlags, mNumVariableMetabs, CTrajectoryMethod::mpCurrentState, mpMetabolites, mPQ, mReactionFlags, mUpperStochLimit, CCopasiParameter::mValue, removeDeterministicReaction(), and CIndexedPriorityQueue::removeStochReaction().

Referenced by CHybridNextReactionRKMethod::doSingleStep().

 {
   size_t i;
   std::set <size_t>::iterator iter, iterEnd;
   C_FLOAT64 key;
 
   for (i = 0; i < mNumVariableMetabs; i++)
     {
       if ((mMetabFlags[i] == LOW) && ((*mpMetabolites)[i]->getValue() >= mUpperStochLimit))
         {
           mMetabFlags[i] = HIGH;
 
           // go through all corresponding reactions and update flags
           for (iter = mMetab2React[i].begin(), iterEnd = mMetab2React[i].end(); iter != iterEnd; iter++)
             {
               mReactionFlags[*iter].mValue--;
 
               // if reaction gets deterministic, insert it into the linked list of deterministic reactions
               if (mReactionFlags[*iter].mValue == 0)
                 {
                   insertDeterministicReaction(*iter);
                   mPQ.removeStochReaction(*iter);
                 }
             }
         }
 
       if ((mMetabFlags[i] == HIGH) && ((*mpMetabolites)[i]->getValue() < mLowerStochLimit))
         {
           mMetabFlags[i] = LOW;
           (*mpMetabolites)[i]->setValue(floor((*mpMetabolites)[i]->getValue()));
           (*mpMetabolites)[i]->refreshConcentration();
 
           // go through all corresponding reactions and update flags
           for (iter = mMetab2React[i].begin(), iterEnd = mMetab2React[i].end(); iter != iterEnd; iter++)
             {
               if (mReactionFlags[*iter].mValue == 0)
                 {
                   removeDeterministicReaction(*iter);
                   /*
                     mPQ.insertStochReaction(*iter, 1234567.8);  // juergen: have to beautify this, number has to be the biggest C_FLOAT64 !!!
                   */
                   calculateAmu(*iter);
                   mAmuOld[*iter] = mAmu[*iter];
                   key = mpCurrentState->getTime() + generateReactionTime(*iter);
                   mPQ.insertStochReaction(*iter, key);
                 }
 
               mReactionFlags[*iter].mValue++;
             }
         }
     }
 
   return;
 }

void CHybridMethod::removeDeterministicReaction ( size_t rIndex )

protected

Removes a deterministic reaction from the linked list in the vector mReactionFlags.

Parameters

rIndex A size_t giving the index of the reaction to be removed from the list of deterministic reactions.

Removes a deterministic reaction from the linked list in the array mReactionFlags.

Parameters

rIndex A size_t giving the index of the reaction to be removed from the list of deterministic reactions.

Definition at line 1259 of file CHybridMethod.cpp.

References mFirstReactionFlag, CHybridStochFlag::mpPrev, and mReactionFlags.

Referenced by partitionSystem().

 {
   if (mReactionFlags[rIndex].mpPrev != NULL)
     // reaction is deterministic
     {
       if (&mReactionFlags[rIndex] != mFirstReactionFlag)
         // reactionFlag is not the first in the linked list
         {
           mReactionFlags[rIndex].mpPrev->mpNext = mReactionFlags[rIndex].mpNext;
 
           if (mReactionFlags[rIndex].mpNext != NULL)
             mReactionFlags[rIndex].mpNext->mpPrev = mReactionFlags[rIndex].mpPrev;
         }
       else
         // reactionFlag is the first in the linked list
         {
           if (mReactionFlags[rIndex].mpNext != NULL) // reactionFlag is not the only one in the linked list
             {
               mFirstReactionFlag = mReactionFlags[rIndex].mpNext;
               mFirstReactionFlag->mpPrev = mFirstReactionFlag;
             }
           else // reactionFlag is the only one in the linked list
             {
               mFirstReactionFlag = NULL;
             }
         }
     }
 
   mReactionFlags[rIndex].mpPrev = NULL;
   mReactionFlags[rIndex].mpNext = NULL;
   return;
 }

void CHybridMethod::rungeKutta ( C_FLOAT64 dt )

protected

Does one 4th order RungeKutta step to integrate the system numerically.

Parameters

dt	A C_FLOAT64 specifying the stepsize
result	A reference to a vector, into which the result, that is the increment vector, will be written

Definition at line 412 of file CHybridMethod.cpp.

References calculateDerivative(), currentState, getState(), k1, k2, k3, k4, mNumVariableMetabs, setState(), and temp.

Referenced by integrateDeterministicPart().

 {
   size_t i;
 
   /* save current state */
   getState(currentState);
 
   /* k1 step: k1 = dt*f(x(n)) */
   calculateDerivative(temp); // systemState == x(n)
 
   for (i = 0; i < mNumVariableMetabs; i++)
     {
       k1[i] = temp[i] * dt;
     }
 
   /* k2 step: k2 = dt*f(x(n) + k1/2) */
   for (i = 0; i < mNumVariableMetabs; i++)
     {
       temp[i] = k1[i] / 2.0 + currentState[i];
     }
 
   setState(temp);
   calculateDerivative(temp); // systemState == x(n) + k1/2
 
   for (i = 0; i < mNumVariableMetabs; i++)
     {
       k2[i] = temp[i] * dt;
     }
 
   /* k3 step: k3 = dt*f(x(n) + k2/2) */
   for (i = 0; i < mNumVariableMetabs; i++)
     {
       temp[i] = k2[i] / 2.0 + currentState[i];
     }
 
   setState(temp);
   calculateDerivative(temp); // systemState == x(n) + k2/2
 
   for (i = 0; i < mNumVariableMetabs; i++)
     {
       k3[i] = temp[i] * dt;
     }
 
   /* k4 step: k4 = dt*f(x(n) + k3); */
   for (i = 0; i < mNumVariableMetabs; i++)
     {
       temp[i] = k3[i] + currentState[i];
     }
 
   setState(temp);
   calculateDerivative(temp); // systemState == x(n) + k3
 
   for (i = 0; i < mNumVariableMetabs; i++)
     {
       k4[i] = temp[i] * dt;
     }
 
   /* Find next position: x(n+1) = x(n) + 1/6*(k1 + 2*k2 + 2*k3 + k4)  */
   for (i = 0; i < mNumVariableMetabs; i++)
     {
       temp[i] = currentState[i] + (1.0 / 6.0) * (k1[i] + 2.0 * k2[i] + 2.0 * k3[i] + k4[i]);
     }
 
   setState(temp);
 
   return;
 }

void CHybridMethod::setState ( std::vector< C_FLOAT64 > & source )

protected

Writes the state specified in the array source into the model. Length of the vector source must be mNumVariableMetabs. (Number of non-fixed metabolites in the model).

Parameters

source A vector reference with length mNumVariableMetabs, holding the state of the system to be set in the model

Writes the state specified in the vector source into the model. Length of the vector source must be mNumVariableMetabs. (Number of non-fixed metabolites in the model).

Parameters

source A vector reference with length mNumVariableMetabs, holding the state of the system to be set in the model

Definition at line 549 of file CHybridMethod.cpp.

References mNumVariableMetabs, mpModel, and CModel::updateSimulatedValues().

Referenced by integrateDeterministicPartEuler(), and rungeKutta().

 {
   size_t i;
 
   for (i = 0; i < mNumVariableMetabs; i++)
     {
       (*mpMetabolites)[i]->setValue(source[i]);
       (*mpMetabolites)[i]->refreshConcentration();
     }
 
   mpModel->updateSimulatedValues(false); //for assignments
   return;
 }

void CHybridMethod::setupBalances ( )

protected

Sets up an internal representation of the balances for each reaction. This is done in order to be able to deal with fixed metabolites and to avoid a time consuming search for the indices of metabolites in the model.

Definition at line 845 of file CHybridMethod.cpp.

References C_INT32, CModelEntity::FIXED, CCopasiVector< T >::getIndex(), CModel::getMetabolitesX(), CModelEntity::getStatus(), CHybridBalance::mIndex, mLocalBalances, mLocalSubstrates, mMaxBalance, mMaxIntBeforeStep, mMaxSteps, CHybridBalance::mMultiplicity, CHybridBalance::mpMetabolite, mpModel, mpReactions, and CCopasiVector< T >::size().

Referenced by initMethod().

 {
   size_t i, j;
   CHybridBalance newElement;
   C_INT32 maxBalance = 0;
   size_t numReactions;
 
   numReactions = mpReactions->size();
   mLocalBalances.clear();
   mLocalBalances.resize(numReactions);
   mLocalSubstrates.clear();
   mLocalSubstrates.resize(numReactions);
 
   for (i = 0; i < numReactions; i++)
     {
       const CCopasiVector <CChemEqElement> * balances =
         &(*mpReactions)[i]->getChemEq().getBalances();
 
       for (j = 0; j < balances->size(); j++)
         {
           newElement.mpMetabolite = const_cast < CMetab* >((*balances)[j]->getMetabolite());
           newElement.mIndex = mpModel->getMetabolitesX().getIndex(newElement.mpMetabolite);
           // + 0.5 to get a rounding out of the static_cast to size_t!
           newElement.mMultiplicity = static_cast<C_INT32>(floor((*balances)[j]->getMultiplicity() + 0.5));
 
           if ((newElement.mpMetabolite->getStatus()) != CModelEntity::FIXED)
             {
               if (newElement.mMultiplicity > maxBalance) maxBalance = newElement.mMultiplicity;
 
               mLocalBalances[i].push_back(newElement); // element is copied for the push_back
             }
         }
 
       balances = &(*mpReactions)[i]->getChemEq().getSubstrates();
 
       for (j = 0; j < balances->size(); j++)
         {
           newElement.mpMetabolite = const_cast < CMetab* >((*balances)[j]->getMetabolite());
           newElement.mIndex = mpModel->getMetabolitesX().getIndex(newElement.mpMetabolite);
           // + 0.5 to get a rounding out of the static_cast to size_t!
           newElement.mMultiplicity = static_cast<C_INT32>(floor((*balances)[j]->getMultiplicity() + 0.5));
 
           mLocalSubstrates[i].push_back(newElement); // element is copied for the push_back
         }
     }
 
   mMaxBalance = maxBalance;
   mMaxIntBeforeStep = INT_MAX - 1 - mMaxSteps * mMaxBalance;
 
   return;
 }

void CHybridMethod::setupDependencyGraph ( )

protected

Sets up the dependency graph

Sets up the dependency graph.

Definition at line 900 of file CHybridMethod.cpp.

References CDependencyGraph::addDependent(), CDependencyGraph::clear(), getAffects(), getDependsOn(), mDG, mpReactions, and CCopasiVector< T >::size().

Referenced by initMethod().

 {
   mDG.clear();
   std::vector< std::set<std::string>* > DependsOn;
   std::vector< std::set<std::string>* > Affects;
   size_t numReactions = mpReactions->size();
   size_t i, j;
 
   // Do for each reaction:
   for (i = 0; i < numReactions; i++)
     {
       // Get the set of metabolites  which affect the value of amu for this
       // reaction i.e. the set on which amu depends. This may be  more than
       // the set of substrates, since the kinetics can involve other
       // reactants, e.g. catalysts. We thus need to step through the
       // rate function and pick out every reactant which can vary.
       DependsOn.push_back(getDependsOn(i));
       // Get the set of metabolites which are affected when this reaction takes place
       Affects.push_back(getAffects(i));
     }
 
   // For each possible pair of reactions i and j, if the intersection of
   // Affects(i) with DependsOn(j) is non-empty, add a dependency edge from i to j.
   for (i = 0; i < numReactions; i++)
     {
       for (j = 0; j < numReactions; j++)
         {
           // Determine whether the intersection of these two sets is non-empty
           // Could also do this with set_intersection generic algorithm, but that
           // would require operator<() to be defined on the set elements.
 
           std::set<std::string>::iterator iter = Affects[i]->begin();
 
           for (; iter != Affects[i]->end(); iter++)
             {
               if (DependsOn[j]->count(*iter))
                 {
                   // The set intersection is non-empty
                   mDG.addDependent(i, j);
                   break;
                 }
             }
         }
 
       // Ensure that self edges are included
       //mDG.addDependent(i, i);
     }
 
   // Delete the memory allocated in getDependsOn() and getAffects()
   // since this is allocated in other functions.
   for (i = 0; i < numReactions; i++)
     {
       delete DependsOn[i];
       delete Affects[i];
     }
 
   return;
 }

void CHybridMethod::setupMetab2React ( )

protected

Creates for each metabolite a set of reaction indices. If the metabolite participates in a reaction as substrate or product this reaction is added to the corresponding set.

Creates for each metabolite a set of reaction indices. If the metabolite participates in a reaction as substrate or product (that means: balance != 0) this reaction is added to the corresponding set.

Definition at line 964 of file CHybridMethod.cpp.

References mLocalBalances, mMetab2React, mpMetabolites, and CCopasiVector< T >::size().

Referenced by initMethod().

 {
   size_t i, j;
   size_t metaboliteIndex;
 
   // Resize mMetab2React and create an initial set for each metabolite
   mMetab2React.clear();
   mMetab2React.resize(mpMetabolites->size());
 
   // Iterate over all reactions
   for (i = 0; i < mLocalBalances.size(); i++)
     {
       // Get the set of metabolites which take part in this reaction
       for (j = 0; j < mLocalBalances[i].size(); j++)
         {
           // find metaboliteIndex and insert the reaction into the set
           metaboliteIndex = mLocalBalances[i][j].mIndex;
           mMetab2React[metaboliteIndex].insert(i);
         }
     }
 
   return;
 }

void CHybridMethod::setupMetab2ReactComplete ( )

protected

Creates for each metabolite a set of reaction indices. Each reaction is dependent on each metabolite resulting in a complete switch.

Definition at line 1030 of file CHybridMethod.cpp.

References mMetab2React, mpMetabolites, mpReactions, and CCopasiVector< T >::size().

 {
   size_t i, j;
 
   // Resize mMetab2React and create an initial set for each metabolite
   mMetab2React.resize(mpMetabolites->size());
 
   // Iterate over all metabolites
   for (i = 0; i < mpMetabolites->size(); i++)
     {
       // Iterate over all reactions
       for (j = 0; j < mpReactions->size(); j++)
         {
           mMetab2React[i].insert(j);
         }
     }
 
   return;
 }

void CHybridMethod::setupMetab2ReactPlusModifier ( )

protected

Creates for each metabolite a set of reaction indices. If the metabolite participates in a reaction as substrate, product or modifier this reaction is added to the corresponding set.

Definition at line 993 of file CHybridMethod.cpp.

References getParticipatesIn(), mMetab2React, mpMetabolites, mpReactions, and CCopasiVector< T >::size().

 {
   std::vector< std::set<size_t>* > participatesIn;
   size_t numReactions = mpReactions->size();
   size_t i;
 
   // Resize mMetab2React and create an initial set for each metabolite
   mMetab2React.resize(mpMetabolites->size());
 
   // Do for each reaction:
   for (i = 0; i < numReactions; i++)
     {
       participatesIn.push_back(getParticipatesIn(i));
     }
 
   // Iterate over all reactions
   for (i = 0; i < numReactions; i++)
     {
       // Get the set of metabolites which take part in this reaction
       std::set<size_t>::iterator iter = participatesIn[i]->begin();
 
       for (; iter != participatesIn[i]->end(); iter++)
         mMetab2React[*iter].insert(i);
     }
 
   for (i = 0; i < numReactions; i++)
     {
       delete participatesIn[i];
     }
 
   return;
 }

void CHybridMethod::setupPartition ( )

protected

Creates an initial partitioning of the system. Deterministic and stochastic reactions are determined. The array mStochReactions is initialized.

Creates an initial partitioning of the system. Deterministic and stochastic reactions are determined. The vector mReactionFlags and the vector mMetabFlags are initialized.

Definition at line 1055 of file CHybridMethod.cpp.

References C_FLOAT64, CCopasiParameter::getValue(), HIGH, LOW, mFirstReactionFlag, mLocalBalances, mLowerStochLimit, mMetabFlags, mNumVariableMetabs, mpMetabolites, CHybridStochFlag::mpNext, mReactionFlags, mUpperStochLimit, and CCopasiParameter::mValue.

Referenced by initMethod().

 {
   size_t i, j;
   CHybridStochFlag * prevFlag;
   C_FLOAT64 averageStochLimit = (mUpperStochLimit + mLowerStochLimit) / 2;
 
   // initialize vector mMetabFlags
   mMetabFlags.clear();
   mMetabFlags.resize(mNumVariableMetabs);
 
   for (i = 0; i < mMetabFlags.size(); i++)
     {
       if ((*mpMetabolites)[i]->getValue() < averageStochLimit)
         {
           mMetabFlags[i] = LOW;
           (*mpMetabolites)[i]->setValue(floor((*mpMetabolites)[i]->getValue()));
           (*mpMetabolites)[i]->refreshConcentration();
         }
       else
         mMetabFlags[i] = HIGH;
     }
 
   // initialize vector mReactionFlags
   mReactionFlags.clear();
   mReactionFlags.resize(mLocalBalances.size());
 
   for (i = 0; i < mLocalBalances.size(); i++)
     {
       mReactionFlags[i].mIndex = i;
       mReactionFlags[i].mValue = 0;
 
       for (j = 0; j < mLocalBalances[i].size(); j++)
         {
           if (mMetabFlags[mLocalBalances[i][j].mIndex] == LOW)
             {
               mReactionFlags[i].mValue++;
             }
         }
     }
 
   mFirstReactionFlag = NULL;
   prevFlag = NULL;
 
   for (i = 0; i < mLocalBalances.size(); i++)
     {
       if (mReactionFlags[i].mValue == 0)
         {
           if (mFirstReactionFlag != NULL)
             {
               prevFlag->mpNext = &mReactionFlags[i];
               mReactionFlags[i].mpPrev = prevFlag;
               prevFlag = &mReactionFlags[i];
             }
           else
             {
               mFirstReactionFlag = &mReactionFlags[i];
               mReactionFlags[i].mpPrev = &mReactionFlags[i]; // Important to distinguish between stochastic (prev == NULL) and deterministic (prev != NULL) reactions
               prevFlag = &mReactionFlags[i];
             }
         }
       else
         {
           mReactionFlags[i].mpPrev = NULL;
           mReactionFlags[i].mpNext = NULL;
         }
     }
 
   if (prevFlag != NULL)
     {
       prevFlag->mpNext = NULL;
     }
 
   return;
 }

void CHybridMethod::setupPriorityQueue ( C_FLOAT64 startTime = 0.0 )

protected

Sets up the priority queue.

Parameters

startTime The time at which the simulation starts.

Definition at line 1135 of file CHybridMethod.cpp.

References C_FLOAT64, calculateAmu(), CIndexedPriorityQueue::clear(), generateReactionTime(), CIndexedPriorityQueue::initializeIndexPointer(), CIndexedPriorityQueue::insertStochReaction(), mPQ, mpReactions, mReactionFlags, and CCopasiVector< T >::size().

Referenced by initMethod().

 {
   size_t i;
   C_FLOAT64 time;
 
   mPQ.clear();
   mPQ.initializeIndexPointer(mpReactions->size());
 
   for (i = 0; i < mpReactions->size(); i++)
     {
       if (mReactionFlags[i].mpPrev == NULL) // Reaction is stochastic!
         {
           calculateAmu(i);
           time = startTime + generateReactionTime(i);
           mPQ.insertStochReaction(i, time);
         }
     }
 
   return;
 }

void CHybridMethod::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 230 of file CHybridMethod.cpp.

References CModel::deterministic, CStateTemplate::getIndex(), CModel::getMetabolitesX(), CCopasiProblem::getModel(), CModel::getModelType(), CModel::getStateTemplate(), CState::getTime(), initMethod(), mDoCorrection, mFirstMetabIndex, mHasAssignments, modelHasAssignments(), CTrajectoryMethod::mpCurrentState, mpModel, CTrajectoryMethod::mpProblem, CModel::setState(), and CModel::updateSimulatedValues().

 {
   *mpCurrentState = *initialState;
 
   mpModel = mpProblem->getModel();
   assert(mpModel);
 
   if (mpModel->getModelType() == CModel::deterministic)
     mDoCorrection = true;
   else
     mDoCorrection = false;
 
   mHasAssignments = modelHasAssignments(mpModel);
 
   mFirstMetabIndex = mpModel->getStateTemplate().getIndex(mpModel->getMetabolitesX()[0]);
 
   mpProblem->getModel()->setState(*mpCurrentState);
 
   mpModel->updateSimulatedValues(false); //for assignments
   //mpModel->updateNonSimulatedValues(); //for assignments
 
   // call init of the simulation method, can be overloaded in derived classes
   initMethod(mpCurrentState->getTime());
 
   return;
 }

CTrajectoryMethod::Status CHybridMethod::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 184 of file CHybridMethod.cpp.

References CState::beginIndependent(), C_FLOAT64, doSingleStep(), CModel::getMetabolitesX(), CCopasiProblem::getModel(), CModel::getNumVariableMetabs(), CState::getTime(), mFirstMetabIndex, mMaxIntBeforeStep, mMaxSteps, mMaxStepsReached, CTrajectoryMethod::mpCurrentState, CTrajectoryMethod::mpProblem, CTrajectoryMethod::NORMAL, CState::setTime(), and CCopasiMessage::WARNING.

 {
   // write the current state to the model
   //  mpProblem->getModel()->setState(mpCurrentState); // is that correct?
 
   // check for possible overflows
   size_t i;
   size_t imax;
 
   // :TODO: Bug 774: This assumes that the number of variable metabs is the number
   // of metabs determined by reaction. In addition they are expected at the beginning of the
   // MetabolitesX which is not the case if we have metabolites of type ODE.
   for (i = 0, imax = mpProblem->getModel()->getNumVariableMetabs(); i < imax; i++)
     if (mpProblem->getModel()->getMetabolitesX()[i]->getValue() >= mMaxIntBeforeStep)
       {
         // throw exception or something like that
       }
 
   // do several steps
   C_FLOAT64 time = mpCurrentState->getTime();
   C_FLOAT64 endTime = time + deltaT;
 
   for (i = 0; ((i < mMaxSteps) && (time < endTime)); i++)
     {
       time = doSingleStep(time, endTime);
     }
 
   mpCurrentState->setTime(time);
 
   if ((i >= mMaxSteps) && (!mMaxStepsReached))
     {
       mMaxStepsReached = true; //only report this message once
       CCopasiMessage(CCopasiMessage::WARNING, "maximum number of reaction events was reached in at least one simulation step.\nThat means time intervals in the output may not be what you requested.");
     }
 
   // get back the particle numbers
 
   /* Set the variable metabolites */
   C_FLOAT64 * Dbl = mpCurrentState->beginIndependent() + mFirstMetabIndex - 1;
 
   for (i = 0, imax = mpProblem->getModel()->getNumVariableMetabs(); i < imax; i++, Dbl++)
     *Dbl = mpProblem->getModel()->getMetabolitesX()[i]->getValue();
 
   return NORMAL;
 }

void CHybridMethod::updatePriorityQueue	(	size_t	rIndex,
		C_FLOAT64	time
	)

protected

Updates the priority queue.

Parameters

rIndex	A size_t giving the index of the fired reaction
time	A C_FLOAT64 holding the time taken by this reaction

Updates the priority queue.

Parameters

rIndex	A size_t giving the index of the fired reaction (-1, if no stochastic reaction has fired)
time	A C_FLOAT64 holding the current time

Definition at line 625 of file CHybridMethod.cpp.

References C_FLOAT64, C_INVALID_INDEX, calculateAmu(), generateReactionTime(), mAmu, mAmuOld, mHasAssignments, mpModel, mPQ, mpReactions, mReactionFlags, mUpdateSet, CCopasiVector< T >::size(), CIndexedPriorityQueue::updateNode(), CModel::updateSimulatedValues(), and updateTauMu().

Referenced by CHybridNextReactionRKMethod::doSingleStep().

 {
   C_FLOAT64 newTime;
   size_t index;
   std::set <size_t>::iterator iter, iterEnd;
 
   //if the model contains assignments we use a less efficient loop over all (stochastic) reactions to capture all changes
   // we do not know the exact dependencies. TODO: this should be changed later in order to get a more efficient update scheme
   if (mHasAssignments)
     {
       mpModel->updateSimulatedValues(false);
 
       for (index = 0; index < mpReactions->size(); index++)
         {
           if (mReactionFlags[index].mpPrev == NULL) // Reaction is stochastic!
             {
               mAmuOld[index] = mAmu[index];
               calculateAmu(index);
 
               if (mAmuOld[index] != mAmu[index])
                 if (index != rIndex) updateTauMu(index, time);
             }
         }
     }
   else
     {
       // iterate through the set of affected reactions and update the stochastic ones in the priority queue
       for (iter = mUpdateSet.begin(), iterEnd = mUpdateSet.end(); iter != iterEnd; iter++)
         {
           if (mReactionFlags[*iter].mpPrev == NULL) // reaction is stochastic!
             {
               index = *iter;
               mAmuOld[index] = mAmu[index];
               calculateAmu(index);
 
               if (*iter != rIndex) updateTauMu(index, time);
             }
         }
     }
 
   // draw new random number and update the reaction just fired
   if ((rIndex != C_INVALID_INDEX) && (mReactionFlags[rIndex].mpPrev == NULL))
     {
       // reaction is stochastic
       newTime = time + generateReactionTime(rIndex);
       mPQ.updateNode(rIndex, newTime);
     }
 
   // empty the mUpdateSet
   mUpdateSet.clear();
   return;
 }

void CHybridMethod::updateTauMu	(	size_t	rIndex,
		C_FLOAT64	time
	)

protected

Updates the putative reaction time of a stochastic reaction in the priority queue. The corresponding amu and amu_old must be set prior to the call of this method.

Parameters

rIndex	A size_t specifying the index of the reaction
time	A C_FLOAT64 specifying the current time

Updates the putative reaction time of a stochastic reaction in the priority queue. The corresponding amu and amu_old must be set prior to the call of this method.

Parameters

rIndex A size_t specifying the index of the reaction

Definition at line 772 of file CHybridMethod.cpp.

References C_FLOAT64, generateReactionTime(), CIndexedPriorityQueue::getKey(), mAmu, mAmuOld, mPQ, and CIndexedPriorityQueue::updateNode().

Referenced by updatePriorityQueue().

 {
   C_FLOAT64 newTime;
 
   // One must make sure that the calculation yields reasonable results even in the cases where mAmu=0 or mAmuOld=0 or both =0. Therefore mAmuOld=0 is checked. If mAmuOld equals 0, then a new random number has to be drawn, because tau equals inf and the stoch. information is lost. If both values equal 0, then tau should remain inf and the update of the queue can be skipped!
 
   if (mAmuOld[rIndex] == 0.0)
     {
       if (mAmu[rIndex] != 0.0)
         {
           newTime = time + generateReactionTime(rIndex);
           mPQ.updateNode(rIndex, newTime);
         }
     }
   else
     {
       newTime = time + (mAmuOld[rIndex] / mAmu[rIndex]) * (mPQ.getKey(rIndex) - time);
       mPQ.updateNode(rIndex, newTime);
     }
 
   return;
 }

Friends And Related Function Documentation

CTrajectoryMethod* CTrajectoryMethod::createMethod ( CCopasiMethod::SubType subType )

friend

Member Data Documentation

std::vector<C_FLOAT64> CHybridMethod::currentState

protected

Definition at line 531 of file CHybridMethod.h.

Referenced by initMethod(), integrateDeterministicPartEuler(), outputDebug(), and rungeKutta().

std::vector<C_FLOAT64> CHybridMethod::k1

protected

Definition at line 532 of file CHybridMethod.h.

Referenced by initMethod(), outputDebug(), and rungeKutta().

std::vector<C_FLOAT64> CHybridMethod::k2

protected

Definition at line 533 of file CHybridMethod.h.

Referenced by initMethod(), outputDebug(), and rungeKutta().

std::vector<C_FLOAT64> CHybridMethod::k3

protected

Definition at line 534 of file CHybridMethod.h.

Referenced by initMethod(), outputDebug(), and rungeKutta().

std::vector<C_FLOAT64> CHybridMethod::k4

protected

Definition at line 535 of file CHybridMethod.h.

Referenced by initMethod(), outputDebug(), and rungeKutta().

std::vector<C_FLOAT64> CHybridMethod::mAmu

protected

The propensities of the stochastic reactions.

Definition at line 593 of file CHybridMethod.h.

Referenced by calculateAmu(), generateReactionTime(), initMethod(), insertDeterministicReaction(), outputDebug(), partitionSystem(), updatePriorityQueue(), and updateTauMu().

std::vector<C_FLOAT64> CHybridMethod::mAmuOld

protected

Definition at line 594 of file CHybridMethod.h.

Referenced by initMethod(), insertDeterministicReaction(), outputDebug(), partitionSystem(), updatePriorityQueue(), and updateTauMu().

CDependencyGraph CHybridMethod::mDG

protected

The graph of reactions and their dependent reactions. When a reaction is executed, the propensities for each of its dependents must be updated.

Definition at line 610 of file CHybridMethod.h.

Referenced by fireReaction(), integrateDeterministicPart(), integrateDeterministicPartEuler(), outputDebug(), and setupDependencyGraph().

bool CHybridMethod::mDoCorrection

protected

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

Definition at line 510 of file CHybridMethod.h.

Referenced by calculateAmu(), and start().

size_t CHybridMethod::mFirstMetabIndex

protected

index of the first metab in CState

Definition at line 477 of file CHybridMethod.h.

Referenced by start(), and step().

CHybridStochFlag* CHybridMethod::mFirstReactionFlag

protected

Definition at line 543 of file CHybridMethod.h.

Referenced by calculateDerivative(), CHybridNextReactionRKMethod::doSingleStep(), insertDeterministicReaction(), integrateDeterministicPart(), integrateDeterministicPartEuler(), outputDebug(), removeDeterministicReaction(), and setupPartition().

bool CHybridMethod::mHasAssignments

protected

Indicates whether the model has global quantities with assignment rules. If it has, we will use a less efficient way to update the model state to handle this.

Definition at line 642 of file CHybridMethod.h.

Referenced by start(), and updatePriorityQueue().

std::vector<std::vector <CHybridBalance> > CHybridMethod::mLocalBalances

protected

Internal representation of the balances of each reaction. The index of each metabolite in the reaction is provided.

Definition at line 555 of file CHybridMethod.h.

Referenced by fireReaction(), getAffects(), outputDebug(), setupBalances(), setupMetab2React(), and setupPartition().

std::vector<std::vector <CHybridBalance> > CHybridMethod::mLocalSubstrates

protected

Internal representation of the substrates of each reaction. The index of each substrate-metabolite is provided.

Definition at line 561 of file CHybridMethod.h.

Referenced by calculateAmu(), outputDebug(), and setupBalances().

C_FLOAT64 CHybridMethod::mLowerStochLimit

protected

Limit for particle numbers. If a particle number is < StochLimit the corresponding reactions must be simulated stochastically.

Definition at line 567 of file CHybridMethod.h.

Referenced by initMethod(), isValidProblem(), outputDebug(), partitionSystem(), and setupPartition().

size_t CHybridMethod::mMaxBalance

protected

maximal increase of a particle number in one step.

Definition at line 500 of file CHybridMethod.h.

Referenced by outputDebug(), and setupBalances().

size_t CHybridMethod::mMaxIntBeforeStep

protected

This is set to maxint - mMaxSteps*mMaxBalance

Definition at line 505 of file CHybridMethod.h.

Referenced by outputDebug(), setupBalances(), and step().

unsigned C_INT32 CHybridMethod::mMaxSteps

protected

Max number of doSingleStep() per step()

Definition at line 482 of file CHybridMethod.h.

Referenced by initMethod(), outputDebug(), setupBalances(), and step().

bool CHybridMethod::mMaxStepsReached

protected

Definition at line 484 of file CHybridMethod.h.

Referenced by initMethod(), and step().

std::vector<std::set <size_t> > CHybridMethod::mMetab2React

protected

Vector of relations between metabolites to reactions.

Definition at line 588 of file CHybridMethod.h.

Referenced by outputDebug(), partitionSystem(), setupMetab2React(), setupMetab2ReactComplete(), and setupMetab2ReactPlusModifier().

std::vector<metabStatus> CHybridMethod::mMetabFlags

protected

Vector holding information on the status of metabolites. They can have low or high particle numbers.

Definition at line 549 of file CHybridMethod.h.

Referenced by outputDebug(), partitionSystem(), and setupPartition().

size_t CHybridMethod::mNumVariableMetabs

protected

Dimension of the system. Total number of metabolites.

Definition at line 472 of file CHybridMethod.h.

Referenced by getState(), initMethod(), integrateDeterministicPartEuler(), outputDebug(), partitionSystem(), rungeKutta(), setState(), and setupPartition().

size_t CHybridMethod::mOutputCounter

protected

Output counter.

Definition at line 635 of file CHybridMethod.h.

Referenced by outputData().

std::ofstream CHybridMethod::mOutputFile

protected

File output stream to write data.

Definition at line 625 of file CHybridMethod.h.

std::string CHybridMethod::mOutputFileName

protected

Output filename.

Definition at line 630 of file CHybridMethod.h.

unsigned C_INT32 CHybridMethod::mPartitioningInterval

protected

The system gets repartitioned after this number of elementary steps.

Definition at line 573 of file CHybridMethod.h.

Referenced by CHybridNextReactionRKMethod::doSingleStep(), initMethod(), and outputDebug().

CCopasiVector<CMetab>* CHybridMethod::mpMetabolites

protected

A pointer to the metabolites of the model.

Definition at line 520 of file CHybridMethod.h.

Referenced by calculateAmu(), initMethod(), outputData(), outputDebug(), partitionSystem(), setupMetab2React(), setupMetab2ReactComplete(), setupMetab2ReactPlusModifier(), and setupPartition().

CModel* CHybridMethod::mpModel

protected

Pointer to the model

Definition at line 467 of file CHybridMethod.h.

Referenced by cleanup(), initMethod(), setState(), setupBalances(), start(), and updatePriorityQueue().

CIndexedPriorityQueue CHybridMethod::mPQ

protected

The set of putative stochastic (!) reactions and associated times at which each reaction occurs. This is represented as a priority queue, sorted by the reaction time. This heap changes dynamically as stochastic reactions become deterministic (delete this reaction from the queue) or vice versa (insert a new reaction into the queue).

Definition at line 619 of file CHybridMethod.h.

Referenced by CHybridNextReactionRKMethod::doSingleStep(), getStochTimeAndIndex(), outputDebug(), partitionSystem(), setupPriorityQueue(), updatePriorityQueue(), and updateTauMu().

CRandom* CHybridMethod::mpRandomGenerator

protected

The random number generator.

Definition at line 604 of file CHybridMethod.h.

Referenced by CHybridMethod(), cleanup(), generateReactionTime(), initMethod(), and outputDebug().

const CCopasiVectorNS<CReaction>* CHybridMethod::mpReactions

protected

A pointer to the reactions of the model.

Definition at line 515 of file CHybridMethod.h.

Referenced by checkModel(), getDependsOn(), initMethod(), outputDebug(), setupBalances(), setupDependencyGraph(), setupMetab2ReactComplete(), setupMetab2ReactPlusModifier(), setupPriorityQueue(), and updatePriorityQueue().

unsigned C_INT32 CHybridMethod::mRandomSeed

protected

The random seed to use.

Definition at line 495 of file CHybridMethod.h.

Referenced by initMethod().

std::vector<CHybridStochFlag> CHybridMethod::mReactionFlags

protected

Vector to hold information about how many metabolites of a reaction have low particle numbers. If no metabolite of one reaction has low particle numbers this reaction will be simulated det.

Definition at line 542 of file CHybridMethod.h.

Referenced by insertDeterministicReaction(), outputDebug(), partitionSystem(), removeDeterministicReaction(), setupPartition(), setupPriorityQueue(), and updatePriorityQueue().