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

#include <CLsodaMethod.h>

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

Classes

struct  Data
 

Public Member Functions

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

Static Public Member Functions

static void EvalF (const C_INT *n, const C_FLOAT64 *t, const C_FLOAT64 *y, C_FLOAT64 *ydot)
 
static void EvalJ (const C_INT *n, const C_FLOAT64 *t, const C_FLOAT64 *y, const C_INT *ml, const C_INT *mu, C_FLOAT64 *pd, const C_INT *nRowPD)
 
static void EvalR (const C_INT *n, const C_FLOAT64 *t, const C_FLOAT64 *y, const C_INT *nr, C_FLOAT64 *r)
 
- Static Public Member Functions inherited from CTrajectoryMethod
static CTrajectoryMethodcreateMethod (CCopasiMethod::SubType subType=CCopasiMethod::deterministic)
 
- Static Public Member Functions inherited from CCopasiObject
static std::vector< Refresh * > buildUpdateSequence (const DataObjectSet &objects, const DataObjectSet &uptoDateObjects, const DataObjectSet &context=DataObjectSet())
 
static void setRenameHandler (CRenameHandler *rh)
 

Protected Types

enum  RootMasking { NONE = 0, ALL, DISCRETE }
 
- 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

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

Protected Attributes

C_INT mLsodaStatus
 
CState mMethodState
 
bool mNoODE
 
CModelmpModel
 
bool * mpReducedModel
 
RootMasking mRootMasking
 
C_FLOAT64 mTime
 
- Protected Attributes inherited from CTrajectoryMethod
CStatempCurrentState
 
CTrajectoryProblemmpProblem
 
CVector< C_INTmRoots
 
- Protected Attributes inherited from CCopasiMethod
CProcessReportmpCallBack
 
- Protected Attributes inherited from CCopasiParameter
std::string mKey
 
CCopasiObjectmpValueReference
 
size_t mSize
 
Value mValue
 
- Protected Attributes inherited from CCopasiContainer
objectMap mObjects
 

Private Member Functions

void createRootMask ()
 
void initializeParameter ()
 
void maskRoots (CVectorCore< C_FLOAT64 > &rootValues)
 
Status peekAhead ()
 

Private Attributes

CVector< C_FLOAT64mAtol
 
Data mData
 
CVector< bool > mDiscreteRoots
 
C_FLOAT64 mDummy
 
CVector< C_FLOAT64mDWork
 
std::ostringstream mErrorMsg
 
CVector< C_INTmIWork
 
C_INT mJType
 
CLSODA mLSODA
 
CLSODAR mLSODAR
 
C_INT mNumRoots
 
C_FLOAT64mpAbsoluteTolerance
 
bool mPeekAheadMode
 
unsigned C_INT32mpMaxInternalSteps
 
C_FLOAT64mpRelativeTolerance
 
unsigned C_INT32 mRootCounter
 
CVector< bool > mRootMask
 
C_FLOAT64 mRtol
 
C_INT mState
 
C_FLOAT64 mTargetTime
 
C_FLOAT64mY
 
CVector< C_FLOAT64mYdot
 

Friends

CTrajectoryMethodCTrajectoryMethod::createMethod (CCopasiMethod::SubType subType)
 

Additional Inherited Members

- Public Types inherited from CTrajectoryMethod
enum  Status { FAILURE = -1, NORMAL = 0, ROOT = 1 }
 
- Public Types inherited from CCopasiMethod
enum  SubType {
  unset = 0, RandomSearch, RandomSearchMaster, SimulatedAnnealing,
  CoranaWalk, DifferentialEvolution, ScatterSearch, GeneticAlgorithm,
  EvolutionaryProgram, SteepestDescent, HybridGASA, GeneticAlgorithmSR,
  HookeJeeves, LevenbergMarquardt, NelderMead, SRES,
  Statistics, ParticleSwarm, Praxis, TruncatedNewton,
  Newton, deterministic, LSODAR, directMethod,
  stochastic, tauLeap, adaptiveSA, hybrid,
  hybridLSODA, hybridODE45, DsaLsodar, tssILDM,
  tssILDMModified, tssCSP, mcaMethodReder, scanMethod,
  lyapWolf, sensMethod, EFMAlgorithm, EFMBitPatternTreeAlgorithm,
  EFMBitPatternAlgorithm, Householder, crossSectionMethod, linearNoiseApproximation
}
 
- Public Types inherited from CCopasiParameterGroup
typedef parameterGroup::iterator index_iterator
 
typedef
CCopasiContainer::objectMap::iterator 
name_iterator
 
typedef std::vector
< CCopasiParameter * > 
parameterGroup
 
- Public Types inherited from CCopasiParameter
enum  Type {
  DOUBLE = 0, UDOUBLE, INT, UINT,
  BOOL, GROUP, STRING, CN,
  KEY, FILE, EXPRESSION, INVALID
}
 
- Public Types inherited from CCopasiContainer
typedef std::multimap
< std::string, CCopasiObject * > 
objectMap
 
- Public Types inherited from CCopasiObject
typedef std::set< const
CCopasiObject * > 
DataObjectSet
 
typedef std::vector< Refresh * > DataUpdateSequence
 
- Public Types inherited from CObjectInterface
typedef std::set< const
CObjectInterface * > 
ObjectSet
 
typedef std::vector
< CObjectInterface * > 
UpdateSequence
 
- Static Public 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 CRenameHandlersmpRenameHandler = NULL
 

Detailed Description

Definition at line 28 of file CLsodaMethod.h.

Member Enumeration Documentation

enum CLsodaMethod::RootMasking
protected
Enumerator
NONE 
ALL 
DISCRETE 

Definition at line 42 of file CLsodaMethod.h.

Constructor & Destructor Documentation

CLsodaMethod::CLsodaMethod ( const CCopasiMethod::SubType subType = deterministic,
const CCopasiContainer pParent = NULL 
)
protected

Default constructor.

Parameters
constCCopasiMethod::SubType & subType (default: deterministic)
constCCopasiContainer * pParent (default: NULL)

Definition at line 25 of file CLsodaMethod.cpp.

References CLsodaMethod::Data::dim, initializeParameter(), mData, and CLsodaMethod::Data::pMethod.

26  :
27  CTrajectoryMethod(subType, pParent),
28  mMethodState(),
29  mY(NULL),
30  mRootMask(),
31  mTargetTime(0.0),
32  mRootCounter(0),
33  mPeekAheadMode(false)
34 {
35  assert((void *) &mData == (void *) &mData.dim);
36 
37  mData.pMethod = this;
39 }
CVector< bool > mRootMask
Definition: CLsodaMethod.h:175
C_FLOAT64 * mY
Definition: CLsodaMethod.h:86
unsigned C_INT32 mRootCounter
Definition: CLsodaMethod.h:199
CState mMethodState
Definition: CLsodaMethod.h:74
C_FLOAT64 mTargetTime
Definition: CLsodaMethod.h:193
CLsodaMethod * pMethod
Definition: CLsodaMethod.h:37
bool mPeekAheadMode
Definition: CLsodaMethod.h:204
void initializeParameter()
CLsodaMethod::CLsodaMethod ( const CLsodaMethod src,
const CCopasiContainer pParent = NULL 
)

Copy constructor.

Parameters
const CLsodaMethod &src
constCCopasiContainer * pParent (default: NULL)

Definition at line 41 of file CLsodaMethod.cpp.

References CLsodaMethod::Data::dim, initializeParameter(), mData, and CLsodaMethod::Data::pMethod.

42  :
43  CTrajectoryMethod(src, pParent),
44  mMethodState(),
45  mY(NULL),
46  mRootMask(src.mRootMask),
50 {
51  assert((void *) &mData == (void *) &mData.dim);
52 
53  mData.pMethod = this;
55 }
CVector< bool > mRootMask
Definition: CLsodaMethod.h:175
C_FLOAT64 * mY
Definition: CLsodaMethod.h:86
unsigned C_INT32 mRootCounter
Definition: CLsodaMethod.h:199
CState mMethodState
Definition: CLsodaMethod.h:74
C_FLOAT64 mTargetTime
Definition: CLsodaMethod.h:193
CLsodaMethod * pMethod
Definition: CLsodaMethod.h:37
bool mPeekAheadMode
Definition: CLsodaMethod.h:204
void initializeParameter()
CLsodaMethod::~CLsodaMethod ( )

Destructor.

Definition at line 57 of file CLsodaMethod.cpp.

58 {}

Member Function Documentation

void CLsodaMethod::createRootMask ( )
private

Create a mask which hides all roots being constant and zero.

Definition at line 608 of file CLsodaMethod.cpp.

References ALL, CVectorCore< CType >::array(), C_FLOAT64, CModel::calculateRootDerivatives(), CModel::evaluateRoots(), min, mMethodState, mpAbsoluteTolerance, mpModel, mpReducedModel, mRootMask, mRootMasking, CTrajectoryMethod::mRoots, CVector< CType >::resize(), CModel::setState(), CVectorCore< CType >::size(), and CModel::updateSimulatedValues().

Referenced by step().

609 {
610  size_t NumRoots = mRoots.size();
611  mRootMask.resize(NumRoots);
612  CVector< C_FLOAT64 > RootValues;
613  RootValues.resize(NumRoots);
614  CVector< C_FLOAT64 > RootDerivatives;
615  RootDerivatives.resize(NumRoots);
616 
618 
619  if (*mpReducedModel)
620  {
622  }
623 
624  mpModel->evaluateRoots(RootValues, true);
625  mpModel->calculateRootDerivatives(RootDerivatives);
626 
627  bool *pMask = mRootMask.array();
628  bool *pMaskEnd = pMask + mRootMask.size();
629  C_FLOAT64 * pRootValue = RootValues.array();
630  C_FLOAT64 * pRootDerivative = RootDerivatives.array();
631 
632  for (; pMask != pMaskEnd; ++pMask, ++pRootValue, ++pRootDerivative)
633  {
634  *pMask = (fabs(*pRootDerivative) < *mpAbsoluteTolerance ||
635  fabs(*pRootValue) < 1e3 * std::numeric_limits< C_FLOAT64 >::min()) ? true : false;
636  }
637 
638  mRootMasking = ALL;
639 }
C_FLOAT64 * mpAbsoluteTolerance
Definition: CLsodaMethod.h:63
CVector< bool > mRootMask
Definition: CLsodaMethod.h:175
void evaluateRoots(CVectorCore< C_FLOAT64 > &rootValues, const bool &ignoreDiscrete)
Definition: CModel.cpp:4679
void updateSimulatedValues(const bool &updateMoieties)
Definition: CModel.cpp:1851
void resize(size_t size, const bool &copy=false)
Definition: CVector.h:301
CModel * mpModel
Definition: CLsodaMethod.h:159
void calculateRootDerivatives(CVector< C_FLOAT64 > &rootDerivatives)
Definition: CModel.cpp:4712
CState mMethodState
Definition: CLsodaMethod.h:74
CVector< C_INT > mRoots
void setState(const CState &state)
Definition: CModel.cpp:1785
bool * mpReducedModel
Definition: CLsodaMethod.h:52
size_t size() const
Definition: CVector.h:100
#define C_FLOAT64
Definition: copasi.h:92
RootMasking mRootMasking
Definition: CLsodaMethod.h:186
CType * array()
Definition: CVector.h:139
#define min(a, b)
Definition: f2c.h:175
void CLsodaMethod::destroyRootMask ( )
protected

Destroy the mask which hides all roots being constant and zero.

Definition at line 641 of file CLsodaMethod.cpp.

References mRootMask, mRootMasking, NONE, and CVector< CType >::resize().

Referenced by CTrajectoryMethodDsaLsodar::fireReaction(), start(), stateChanged(), and step().

642 {
643  mRootMask.resize(0);
644  mRootMasking = NONE;
645 }
CVector< bool > mRootMask
Definition: CLsodaMethod.h:175
void resize(size_t size, const bool &copy=false)
Definition: CVector.h:301
RootMasking mRootMasking
Definition: CLsodaMethod.h:186
bool CLsodaMethod::elevateChildren ( )
virtual

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

Returns
bool success

Reimplemented from CCopasiParameterGroup.

Reimplemented in CTrajectoryMethodDsaLsodar.

Definition at line 150 of file CLsodaMethod.cpp.

References initializeParameter().

Referenced by CTrajectoryMethodDsaLsodar::elevateChildren().

151 {
153  return true;
154 }
void initializeParameter()
void CLsodaMethod::EvalF ( const C_INT n,
const C_FLOAT64 t,
const C_FLOAT64 y,
C_FLOAT64 ydot 
)
static

This evaluates the derivatives

Definition at line 528 of file CLsodaMethod.cpp.

References evalF(), and CLsodaMethod::Data::pMethod.

Referenced by step().

529 {static_cast<Data *>((void *) n)->pMethod->evalF(t, y, ydot);}
void CLsodaMethod::evalF ( const C_FLOAT64 t,
const C_FLOAT64 y,
C_FLOAT64 ydot 
)
virtual

Reimplemented in CTrajectoryMethodDsaLsodar.

Definition at line 531 of file CLsodaMethod.cpp.

References CModel::calculateDerivatives(), CModel::calculateDerivativesX(), mMethodState, mNoODE, mpModel, mpReducedModel, CModel::setState(), CState::setTime(), and CModel::updateSimulatedValues().

Referenced by EvalF(), and stateChanged().

532 {
533  // If we have no ODEs add a constant one.
534  if (mNoODE)
535  {
536  *ydot = 1.0;
537  return;
538  }
539 
540  mMethodState.setTime(*t);
541 
544 
545  if (*mpReducedModel)
547  else
549 
550  return;
551 }
void updateSimulatedValues(const bool &updateMoieties)
Definition: CModel.cpp:1851
CModel * mpModel
Definition: CLsodaMethod.h:159
void setTime(const C_FLOAT64 &time)
Definition: CState.cpp:326
CState mMethodState
Definition: CLsodaMethod.h:74
void setState(const CState &state)
Definition: CModel.cpp:1785
bool * mpReducedModel
Definition: CLsodaMethod.h:52
void calculateDerivativesX(C_FLOAT64 *derivativesX)
Definition: CModel.cpp:1941
void calculateDerivatives(C_FLOAT64 *derivatives)
Definition: CModel.cpp:1903
void CLsodaMethod::EvalJ ( const C_INT n,
const C_FLOAT64 t,
const C_FLOAT64 y,
const C_INT ml,
const C_INT mu,
C_FLOAT64 pd,
const C_INT nRowPD 
)
static

This evaluates the Jacobian

Definition at line 582 of file CLsodaMethod.cpp.

References evalJ(), and CLsodaMethod::Data::pMethod.

Referenced by step().

584 {static_cast<Data *>((void *) n)->pMethod->evalJ(t, y, ml, mu, pd, nRowPD);}
void CLsodaMethod::evalJ ( const C_FLOAT64 t,
const C_FLOAT64 y,
const C_INT ml,
const C_INT mu,
C_FLOAT64 pd,
const C_INT nRowPD 
)
virtual

Definition at line 587 of file CLsodaMethod.cpp.

Referenced by EvalJ().

589 {
590  // TODO Implement me.
591 }
void CLsodaMethod::EvalR ( const C_INT n,
const C_FLOAT64 t,
const C_FLOAT64 y,
const C_INT nr,
C_FLOAT64 r 
)
static

This evaluates the roots

Definition at line 553 of file CLsodaMethod.cpp.

References evalR(), and CLsodaMethod::Data::pMethod.

Referenced by step().

555 {static_cast<Data *>((void *) n)->pMethod->evalR(t, y, nr, r);}
void CLsodaMethod::evalR ( const C_FLOAT64 t,
const C_FLOAT64 y,
const C_INT nr,
C_FLOAT64 r 
)
virtual

Reimplemented in CTrajectoryMethodDsaLsodar.

Definition at line 557 of file CLsodaMethod.cpp.

References C_INT, CModel::evaluateRoots(), maskRoots(), mMethodState, mpModel, mpReducedModel, mRootMasking, CTrajectoryMethod::mRoots, NONE, CModel::setState(), CState::setTime(), CVectorCore< CType >::size(), and CModel::updateSimulatedValues().

Referenced by EvalR().

559 {
560  assert(*nr == (C_INT) mRoots.size());
561 
562  mMethodState.setTime(*t);
563 
565 
566  if (*mpReducedModel)
567  {
569  }
570 
571  CVectorCore< C_FLOAT64 > RootValues(*nr, r);
572 
573  mpModel->evaluateRoots(RootValues, true);
574 
575  if (mRootMasking != NONE)
576  {
577  maskRoots(RootValues);
578  }
579 };
#define C_INT
Definition: copasi.h:115
void evaluateRoots(CVectorCore< C_FLOAT64 > &rootValues, const bool &ignoreDiscrete)
Definition: CModel.cpp:4679
void updateSimulatedValues(const bool &updateMoieties)
Definition: CModel.cpp:1851
CModel * mpModel
Definition: CLsodaMethod.h:159
void maskRoots(CVectorCore< C_FLOAT64 > &rootValues)
void setTime(const C_FLOAT64 &time)
Definition: CState.cpp:326
CState mMethodState
Definition: CLsodaMethod.h:74
CVector< C_INT > mRoots
void setState(const CState &state)
Definition: CModel.cpp:1785
bool * mpReducedModel
Definition: CLsodaMethod.h:52
size_t size() const
Definition: CVector.h:100
RootMasking mRootMasking
Definition: CLsodaMethod.h:186
void CLsodaMethod::initializeParameter ( )
private

Initialize the method parameter

Definition at line 60 of file CLsodaMethod.cpp.

References CCopasiParameterGroup::assertParameter(), CCopasiParameter::BOOL, C_FLOAT64, C_INT32, CModel::getCompartments(), CCopasiDataModel::getModel(), CCopasiObject::getObjectDataModel(), CCopasiParameterGroup::getParameter(), CModel::getQuantity2NumberFactor(), CCopasiParameter::getValue(), max, mpAbsoluteTolerance, mpMaxInternalSteps, mpReducedModel, mpRelativeTolerance, CCopasiParameter::Value::pBOOL, CCopasiParameter::Value::pUDOUBLE, CCopasiParameter::Value::pUINT, CCopasiParameterGroup::removeParameter(), CCopasiVector< T >::size(), CCopasiParameter::UDOUBLE, and CCopasiParameter::UINT.

Referenced by CLsodaMethod(), and elevateChildren().

61 {
62  CCopasiParameter *pParm;
63 
65  assertParameter("Integrate Reduced Model", CCopasiParameter::BOOL, (bool) false)->getValue().pBOOL;
67  assertParameter("Relative Tolerance", CCopasiParameter::UDOUBLE, (C_FLOAT64) 1.0e-6)->getValue().pUDOUBLE;
69  assertParameter("Absolute Tolerance", CCopasiParameter::UDOUBLE, (C_FLOAT64) 1.0e-12)->getValue().pUDOUBLE;
71  assertParameter("Max Internal Steps", CCopasiParameter::UINT, (unsigned C_INT32) 10000)->getValue().pUINT;
72 
73  // Check whether we have a method with the old parameter names
74  if ((pParm = getParameter("LSODA.RelativeTolerance")) != NULL)
75  {
77  removeParameter("LSODA.RelativeTolerance");
78 
79  if ((pParm = getParameter("LSODA.AbsoluteTolerance")) != NULL)
80  {
82  removeParameter("LSODA.AbsoluteTolerance");
83  }
84 
85  if ((pParm = getParameter("LSODA.AdamsMaxOrder")) != NULL)
86  {
87  removeParameter("LSODA.AdamsMaxOrder");
88  }
89 
90  if ((pParm = getParameter("LSODA.BDFMaxOrder")) != NULL)
91  {
92  removeParameter("LSODA.BDFMaxOrder");
93  }
94 
95  if ((pParm = getParameter("LSODA.MaxStepsInternal")) != NULL)
96  {
97  *mpMaxInternalSteps = *pParm->getValue().pUINT;
98  removeParameter("LSODA.MaxStepsInternal");
99  }
100  }
101 
102  // Check whether we have a method with "Use Default Absolute Tolerance"
103  if ((pParm = getParameter("Use Default Absolute Tolerance")) != NULL)
104  {
105  C_FLOAT64 NewValue;
106 
107  if (*pParm->getValue().pBOOL)
108  {
109  // The default
110  NewValue = 1.e-12;
111  }
112  else
113  {
114  C_FLOAT64 OldValue = *mpAbsoluteTolerance;
115  CCopasiDataModel* pDataModel = getObjectDataModel();
116  assert(pDataModel != NULL);
117  CModel * pModel = pDataModel->getModel();
118 
119  if (pModel == NULL)
120  // The default
121  NewValue = 1.e-12;
122  else
123  {
124  const CCopasiVectorNS< CCompartment > & Compartment = pModel->getCompartments();
125  size_t i, imax;
127 
128  for (i = 0, imax = Compartment.size(); i < imax; i++)
129  if (Compartment[i]->getValue() < Volume)
130  Volume = Compartment[i]->getValue();
131 
133  // The default
134  NewValue = 1.e-12;
135  else
136  // Invert the scaling as best as we can
137  NewValue = OldValue / (Volume * pModel->getQuantity2NumberFactor());
138  }
139  }
140 
141  *mpAbsoluteTolerance = NewValue;
142  removeParameter("Use Default Absolute Tolerance");
143  }
144 
145  // These parameters are no longer supported.
146  removeParameter("Adams Max Order");
147  removeParameter("BDF Max Order");
148 }
CCopasiDataModel * getObjectDataModel()
virtual size_t size() const
C_FLOAT64 * mpAbsoluteTolerance
Definition: CLsodaMethod.h:63
#define C_INT32
Definition: copasi.h:90
bool removeParameter(const std::string &name)
const C_FLOAT64 & getQuantity2NumberFactor() const
Definition: CModel.cpp:2354
bool * mpReducedModel
Definition: CLsodaMethod.h:52
const Value & getValue() const
unsigned C_INT32 * pUINT
CCopasiParameter * getParameter(const std::string &name)
unsigned C_INT32 * mpMaxInternalSteps
Definition: CLsodaMethod.h:68
CCopasiVectorNS< CCompartment > & getCompartments()
Definition: CModel.cpp:1145
#define C_FLOAT64
Definition: copasi.h:92
CCopasiParameter * assertParameter(const std::string &name, const CCopasiParameter::Type type, const CType &defaultValue)
Definition: CModel.h:50
C_FLOAT64 * mpRelativeTolerance
Definition: CLsodaMethod.h:58
#define max(a, b)
Definition: f2c.h:176
void CLsodaMethod::maskRoots ( CVectorCore< C_FLOAT64 > &  rootValues)
private

Mask roots which are constant and zero.

Parameters
CVectorCore<C_FLOAT64 > & rootValues

Definition at line 593 of file CLsodaMethod.cpp.

References CVectorCore< CType >::array(), C_FLOAT64, mRootMask, and CVectorCore< CType >::size().

Referenced by evalR().

594 {
595  const bool *pMask = mRootMask.array();
596  const bool *pMaskEnd = pMask + mRootMask.size();
597  C_FLOAT64 * pRoot = rootValues.array();
598 
599  for (; pMask != pMaskEnd; ++pMask, ++pRoot)
600  {
601  if (*pMask)
602  {
603  *pRoot = 1.0;
604  }
605  }
606 }
CVector< bool > mRootMask
Definition: CLsodaMethod.h:175
size_t size() const
Definition: CVector.h:100
#define C_FLOAT64
Definition: copasi.h:92
CType * array()
Definition: CVector.h:139
CTrajectoryMethod::Status CLsodaMethod::peekAhead ( )
private

Peek ahead to detect simultaneous roots.

Definition at line 647 of file CLsodaMethod.cpp.

References ALL, CVectorCore< CType >::array(), CState::beginIndependent(), C_FLOAT64, C_INT, CLsodaMethod::Data::dim, CTrajectoryMethod::FAILURE, CState::getTime(), mAtol, max, mData, mDWork, mIWork, mLSODAR, mMethodState, mpAbsoluteTolerance, mPeekAheadMode, mpRelativeTolerance, mRootMasking, CTrajectoryMethod::mRoots, mTargetTime, mTime, CTrajectoryMethod::NORMAL, CInternalSolver::resetState(), CTrajectoryMethod::ROOT, CInternalSolver::saveState(), CState::setTime(), CVectorCore< CType >::size(), and step().

Referenced by step().

648 {
649  // Save the current state
651 
652  CState StartState = mMethodState;
653 
654  CState ResetState = mMethodState;
655  mLSODAR.saveState();
656  CVector< C_FLOAT64 > ResetDWork = mDWork;
657  CVector< C_INT > ResetIWork = mIWork;
658 
659  mPeekAheadMode = true;
660  Status PeekAheadStatus = ROOT;
661 
662  CVector< C_INT > CombinedRoots = mRoots;
663 
664  C_FLOAT64 MaxPeekAheadTime = std::max(mTargetTime, mTime * (1.0 + 2.0 * *mpRelativeTolerance));
665 
666  while (mPeekAheadMode)
667  {
668  switch (step(MaxPeekAheadTime - mTime))
669  {
670  case ROOT:
671  {
672  // Check whether the new state is within the tolerances
673  C_FLOAT64 * pOld = StartState.beginIndependent();
674  C_FLOAT64 * pOldEnd = pOld + mData.dim;
676  C_FLOAT64 * pAtol = mAtol.array();
677 
678  for (; pOld != pOldEnd; ++pOld, ++pNew, ++pAtol)
679  {
680  if ((2.0 * fabs(*pNew - *pOld) > (fabs(*pNew) + fabs(*pOld)) * *mpRelativeTolerance) &&
681  fabs(*pNew) > *pAtol &&
682  fabs(*pOld) > *pAtol)
683  {
684  break;
685  }
686  }
687 
688  if (pOld != pOldEnd ||
689  (mTime > StartState.getTime() * (1.0 + *mpRelativeTolerance) &&
690  fabs(mTime) > *mpAbsoluteTolerance &&
691  fabs(StartState.getTime()) > *mpAbsoluteTolerance))
692  {
693  mPeekAheadMode = false;
694  }
695  else
696  {
697  ResetState = mMethodState;
698  mLSODAR.saveState();
699  ResetDWork = mDWork;
700  ResetIWork = mIWork;
701 
702  // Combine all the roots
703  C_INT * pRoot = mRoots.array();
704  C_INT * pRootEnd = pRoot + mRoots.size();
705  C_INT * pCombinedRoot = CombinedRoots.array();
706 
707  for (; pRoot != pRootEnd; ++pRoot, ++pCombinedRoot)
708  {
709  if (*pRoot > 0)
710  {
711  *pCombinedRoot = 1;
712  }
713  }
714  }
715  }
716 
717  break;
718 
719  case FAILURE:
720  mPeekAheadMode = false;
721  PeekAheadStatus = FAILURE;
722  break;
723 
724  case NORMAL:
725 
726  if (mRootMasking != ALL)
727  {
728  mPeekAheadMode = false;
729  }
730 
731  break;
732  }
733  }
734 
735  // Reset the integrator to the saved state
736  mMethodState = ResetState;
739  mDWork = ResetDWork;
740  mIWork = ResetIWork;
741 
742  mPeekAheadMode = false;
743 
744  mRoots = CombinedRoots;
745 
746  return PeekAheadStatus;
747 }
#define C_INT
Definition: copasi.h:115
CVector< C_FLOAT64 > mAtol
Definition: CLsodaMethod.h:118
C_FLOAT64 * mpAbsoluteTolerance
Definition: CLsodaMethod.h:63
Definition: CState.h:305
CVector< C_FLOAT64 > mDWork
Definition: CLsodaMethod.h:143
void setTime(const C_FLOAT64 &time)
Definition: CState.cpp:326
virtual Status step(const double &deltaT)
CState mMethodState
Definition: CLsodaMethod.h:74
C_FLOAT64 mTargetTime
Definition: CLsodaMethod.h:193
CVector< C_INT > mRoots
C_FLOAT64 mTime
Definition: CLsodaMethod.h:102
size_t size() const
Definition: CVector.h:100
const C_FLOAT64 & getTime() const
Definition: CState.cpp:325
#define C_FLOAT64
Definition: copasi.h:92
RootMasking mRootMasking
Definition: CLsodaMethod.h:186
bool mPeekAheadMode
Definition: CLsodaMethod.h:204
CType * array()
Definition: CVector.h:139
CLSODAR mLSODAR
Definition: CLsodaMethod.h:133
CVector< C_INT > mIWork
Definition: CLsodaMethod.h:148
C_FLOAT64 * mpRelativeTolerance
Definition: CLsodaMethod.h:58
C_FLOAT64 * beginIndependent()
Definition: CState.cpp:328
#define max(a, b)
Definition: f2c.h:176
void CLsodaMethod::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.

Reimplemented in CTrajectoryMethodDsaLsodar.

Definition at line 448 of file CLsodaMethod.cpp.

References CVectorCore< CType >::array(), CState::beginIndependent(), C_INT, destroyRootMask(), CLsodaMethod::Data::dim, CCopasiProblem::getModel(), CModel::getNumDependentReactionMetabs(), CState::getNumIndependent(), CModel::getNumRoots(), CModel::getRootFinders(), CState::getTime(), CModel::initializeAtolVector(), mAtol, mData, mDiscreteRoots, mDummy, mDWork, mErrorMsg, mIWork, mJType, mLSODA, mLSODAR, mLsodaStatus, mMethodState, mNoODE, mNumRoots, mpAbsoluteTolerance, mPeekAheadMode, mpMaxInternalSteps, mpModel, CTrajectoryMethod::mpProblem, mpReducedModel, mpRelativeTolerance, mRootCounter, CTrajectoryMethod::mRoots, mRtol, mState, mTargetTime, mTime, mY, mYdot, CVector< CType >::resize(), and CInternalSolver::setOstream().

Referenced by CTrajectoryMethodDsaLsodar::start().

449 {
450  /* Retrieve the model to calculate */
452 
453  /* Reset lsoda */
454  mLsodaStatus = 1;
455  mState = 1;
456  mJType = 2;
457  mErrorMsg.str("");
458 
459  /* Release previous state and make the initialState the current */
460  mMethodState = *initialState;
462  mTargetTime = mTime;
463  mRootCounter = 0;
464  mPeekAheadMode = false;
465 
468  destroyRootMask();
469 
470  if (*mpReducedModel)
472  else
474 
475  // When we have roots we need to add an artificial ODE dDummy/dt = 1
476  if (mData.dim == 0 && mNumRoots != 0)
477  {
478  mData.dim = 1;
479  mNoODE = true;
480  mAtol.resize(1);
482  mDummy = 0;
483  mY = &mDummy;
484  }
485  else
486  {
487  mNoODE = false;
490  }
491 
493 
494  /* Configure lsoda(r) */
496 
497  mDWork.resize(22 + mData.dim * std::max<C_INT>(16, mData.dim + 9) + 3 * mNumRoots);
498  mDWork[4] = mDWork[5] = mDWork[6] = mDWork[7] = mDWork[8] = mDWork[9] = 0.0;
499  mIWork.resize(20 + mData.dim);
500  mIWork[4] = mIWork[6] = mIWork[9] = 0;
501 
503  mIWork[7] = 12;
504  mIWork[8] = 5;
505 
506  if (mNumRoots > 0)
507  {
510 
511  CMathTrigger::CRootFinder * const* ppRootFinder = mpModel->getRootFinders().array();
512  CMathTrigger::CRootFinder * const* ppRootFinderEnd = ppRootFinder + mNumRoots;
513  bool * pDiscrete = mDiscreteRoots.array();
514 
515  for (; ppRootFinder != ppRootFinderEnd; ++ppRootFinder, ++pDiscrete)
516  {
517  *pDiscrete = (*ppRootFinder)->isDiscrete();
518  }
519  }
520  else
521  {
523  }
524 
525  return;
526 }
#define C_INT
Definition: copasi.h:115
C_INT mLsodaStatus
Definition: CLsodaMethod.h:107
CVector< C_FLOAT64 > mAtol
Definition: CLsodaMethod.h:118
C_FLOAT64 mDummy
Definition: CLsodaMethod.h:170
const CVector< CMathTrigger::CRootFinder * > & getRootFinders() const
Definition: CModel.cpp:4717
C_FLOAT64 * mpAbsoluteTolerance
Definition: CLsodaMethod.h:63
void setOstream(std::ostream &os)
C_FLOAT64 mRtol
Definition: CLsodaMethod.h:113
CTrajectoryProblem * mpProblem
CVector< C_FLOAT64 > mDWork
Definition: CLsodaMethod.h:143
void resize(size_t size, const bool &copy=false)
Definition: CVector.h:301
CVector< bool > mDiscreteRoots
Definition: CLsodaMethod.h:180
C_FLOAT64 * mY
Definition: CLsodaMethod.h:86
CModel * mpModel
Definition: CLsodaMethod.h:159
unsigned C_INT32 mRootCounter
Definition: CLsodaMethod.h:199
CState mMethodState
Definition: CLsodaMethod.h:74
C_FLOAT64 mTargetTime
Definition: CLsodaMethod.h:193
std::ostringstream mErrorMsg
Definition: CLsodaMethod.h:123
CVector< C_INT > mRoots
C_INT mNumRoots
Definition: CLsodaMethod.h:96
size_t getNumDependentReactionMetabs() const
Definition: CModel.cpp:1133
bool * mpReducedModel
Definition: CLsodaMethod.h:52
size_t getNumIndependent() const
Definition: CState.cpp:342
void destroyRootMask()
C_FLOAT64 mTime
Definition: CLsodaMethod.h:102
unsigned C_INT32 * mpMaxInternalSteps
Definition: CLsodaMethod.h:68
const C_FLOAT64 & getTime() const
Definition: CState.cpp:325
bool mPeekAheadMode
Definition: CLsodaMethod.h:204
CType * array()
Definition: CVector.h:139
CLSODAR mLSODAR
Definition: CLsodaMethod.h:133
CVector< C_FLOAT64 > initializeAtolVector(const C_FLOAT64 &baseTolerance, const bool &reducedModel) const
Definition: CModel.cpp:4368
CVector< C_FLOAT64 > mYdot
Definition: CLsodaMethod.h:91
CVector< C_INT > mIWork
Definition: CLsodaMethod.h:148
size_t getNumRoots() const
Definition: CModel.cpp:4707
CModel * getModel() const
C_FLOAT64 * mpRelativeTolerance
Definition: CLsodaMethod.h:58
C_FLOAT64 * beginIndependent()
Definition: CState.cpp:328
void CLsodaMethod::stateChanged ( )
virtual

Inform the trajectory method that the state has changed outside its control

Reimplemented from CTrajectoryMethod.

Reimplemented in CTrajectoryMethodDsaLsodar.

Definition at line 157 of file CLsodaMethod.cpp.

References CVectorCore< CType >::array(), CState::beginFixed(), CState::beginIndependent(), C_FLOAT64, destroyRootMask(), CLsodaMethod::Data::dim, CState::endFixed(), evalF(), CState::getTime(), mData, mLsodaStatus, mMethodState, mNoODE, CTrajectoryMethod::mpCurrentState, mPeekAheadMode, mTime, and mY.

Referenced by CTrajectoryMethodDsaLsodar::stateChanged().

158 {
159  if (!mNoODE && mLsodaStatus != 1)
160  {
161  // Compare the independent state variables
162  // This an be done directly by comparing mMethodState and *mpCurrentState
163  C_FLOAT64 *pMethod = mY;
164  C_FLOAT64 *pMethodEnd = pMethod + mData.dim;
166 
167  for (; pMethod != pMethodEnd; ++pMethod, ++pCurrent)
168  {
169  // We may need to use the absolute and relative tolerance
170  if (*pMethod != *pCurrent)
171  {
172  mLsodaStatus = 1;
173 
174  break;
175  }
176  }
177 
178  // Bug 2106 It is not sufficient to check for the state variables we also need
179  // to check for fixed values which are event targets.
180  pMethod = mMethodState.beginFixed();
181  pMethodEnd = mMethodState.endFixed();
182  pCurrent = mpCurrentState->beginFixed();
183 
184  for (; pMethod != pMethodEnd; ++pMethod, ++pCurrent)
185  {
186  // We may need to use the absolute and relative tolerance
187  if (*pMethod != *pCurrent)
188  {
189  mLsodaStatus = 1;
190 
191  break;
192  }
193  }
194 
195  if (mLsodaStatus != 1)
196  {
197  // Compare the rates of the independent state variables
198  // we need to call evalF for mMethodState and *mpCurrentState and compare the returned rates.
199  CVector< C_FLOAT64 > MethodRate(mData.dim);
200  CVector< C_FLOAT64 > CurrentRate(mData.dim);
201 
202  evalF(&mTime, mY, MethodRate.array());
203 
205  mTime = mMethodState.getTime();
206 
207  evalF(&mTime, mY, CurrentRate.array());
208 
209  pMethod = MethodRate.array();
210  pMethodEnd = pMethod + mData.dim;
211  pCurrent = CurrentRate.array();
212 
213  for (; pMethod != pMethodEnd; ++pMethod, ++pCurrent)
214  {
215  // We may need to use the absolute and relative tolerance
216  if (*pMethod != *pCurrent)
217  {
218  mLsodaStatus = 1;
219  break;
220  }
221  }
222  }
223  }
224 
225  mMethodState = *mpCurrentState;
226  mTime = mMethodState.getTime();
227 
228  mPeekAheadMode = false;
229  destroyRootMask();
230 }
C_INT mLsodaStatus
Definition: CLsodaMethod.h:107
C_FLOAT64 * mY
Definition: CLsodaMethod.h:86
CState mMethodState
Definition: CLsodaMethod.h:74
C_FLOAT64 * endFixed()
Definition: CState.cpp:333
void destroyRootMask()
C_FLOAT64 mTime
Definition: CLsodaMethod.h:102
#define C_FLOAT64
Definition: copasi.h:92
bool mPeekAheadMode
Definition: CLsodaMethod.h:204
C_FLOAT64 * beginFixed()
Definition: CState.cpp:332
virtual void evalF(const C_FLOAT64 *t, const C_FLOAT64 *y, C_FLOAT64 *ydot)
C_FLOAT64 * beginIndependent()
Definition: CState.cpp:328
CTrajectoryMethod::Status CLsodaMethod::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 time step taken.

Parameters
const double &deltaT
Returns
Status status

Reimplemented from CTrajectoryMethod.

Reimplemented in CTrajectoryMethodDsaLsodar.

Definition at line 232 of file CLsodaMethod.cpp.

References ALL, CVectorCore< CType >::array(), C_FLOAT64, C_INT, createRootMask(), destroyRootMask(), CLsodaMethod::Data::dim, DISCRETE, EvalF(), EvalJ(), EvalR(), CCopasiMessage::EXCEPTION, CTrajectoryMethod::FAILURE, CState::getTime(), CState::isValid(), mAtol, MCTrajectoryMethod, mData, mDiscreteRoots, mDWork, mErrorMsg, mIWork, mJType, mLSODA, mLSODAR, mLsodaStatus, mMethodState, mNumRoots, CTrajectoryMethod::mpCurrentState, mPeekAheadMode, mpMaxInternalSteps, mRootCounter, mRootMask, mRootMasking, CTrajectoryMethod::mRoots, mRtol, mState, mTargetTime, mTime, mY, NONE, CTrajectoryMethod::NORMAL, peekAhead(), CTrajectoryMethod::ROOT, CState::setTime(), and CVectorCore< CType >::size().

Referenced by CTrajectoryMethodDsaLsodar::integrateDeterministicPart(), and peekAhead().

233 {
234  if (mData.dim == 0 && mNumRoots == 0) //just do nothing if there are no variables
235  {
236  mTime = mTime + deltaT;
239 
240  return NORMAL;
241  }
242 
243  C_FLOAT64 EndTime = mTime + deltaT;
244 
245  if (mTargetTime != EndTime)
246  {
247  // We have a new end time and reset the root counter.
248  mTargetTime = EndTime;
249  mRootCounter = 0;
250  }
251  else
252  {
253  // We are called with the same end time which means a root has previously been
254  // found. We increase the root counter and check whether the limit is reached.
255  mRootCounter++;
256 
258  {
259  return FAILURE;
260  }
261  }
262 
263  C_INT ITOL = 2; // mRtol scalar, mAtol vector
264  C_INT one = 1;
265  C_INT DSize = (C_INT) mDWork.size();
266  C_INT ISize = (C_INT) mIWork.size();
267 
268  // The return status of the integrator.
269  Status Status = NORMAL;
270 
271  if (mRoots.size() > 0)
272  {
273  mLSODAR(&EvalF, // 1. evaluate F
274  &mData.dim, // 2. number of variables
275  mY, // 3. the array of current concentrations
276  &mTime, // 4. the current time
277  &EndTime, // 5. the final time
278  &ITOL, // 6. error control
279  &mRtol, // 7. relative tolerance array
280  mAtol.array(), // 8. absolute tolerance array
281  &mState, // 9. output by overshoot & interpolation
282  &mLsodaStatus, // 10. the state control variable
283  &one, // 11. further options (one)
284  mDWork.array(), // 12. the double work array
285  &DSize, // 13. the double work array size
286  mIWork.array(), // 14. the int work array
287  &ISize, // 15. the int work array size
288  &EvalJ, // 16. evaluate J (not given)
289  &mJType, // 17. type of j evaluation 2 internal full matrix
290  &EvalR, // 18. evaluate constraint functions
291  &mNumRoots, // 19. number of constraint functions g(i)
292  mRoots.array()); // 20. integer array of length NG for output of root information
293 
294  // There exist situations where LSODAR reports status = 3, which are actually status = -33
295  // Obviously the trivial case is where LSODAR did not advance at all, i.e, the start time
296  // equals the current time. It may also happen that a very small steps has been taken in
297  // we reset short before we reach the internal step limit.
298  if (mLsodaStatus == 3 &&
299  (mRootCounter > 0.99 * *mpMaxInternalSteps ||
301  {
302  mLsodaStatus = -33;
303  mRootCounter = 0;
304  }
305 
306  switch (mLsodaStatus)
307  {
308  case -33:
309 
310  switch (mRootMasking)
311  {
312  case NONE:
313  case DISCRETE:
314  // Reset the integrator to the state before the failed integration.
315 
318  mLsodaStatus = 1;
319 
320  // Create a mask which hides all roots being constant and zero.
321  createRootMask();
322  break;
323 
324  case ALL:
325  break;
326  }
327 
328  break;
329 
330  case 3:
331 
332  // If mLsodaStatus == 3 we have found a root. This needs to be indicated to
333  // the caller as it is not sufficient to rely on the fact that T < TOUT
334  if (mLsodaStatus == 3)
335  {
336  // It is sufficient to switch to 2. Eventual state changes due to events
337  // are indicated via the method stateChanged()
338  mLsodaStatus = 2;
339  Status = ROOT;
340  }
341 
342  // To detect simultaneous roots we have to peek ahead, i.e., continue
343  // integration until the state changes are larger that the relative
344  if (!mPeekAheadMode)
345  {
346  Status = peekAhead();
347  }
348 
349  // The break statement is intentionally missing since we
350  // have to continue to check the root masking state.
351  default:
352 
353  switch (mRootMasking)
354  {
355  case NONE:
356  case DISCRETE:
357  break;
358 
359  case ALL:
360  {
361  const bool * pDiscrete = mDiscreteRoots.array();
362  bool * pMask = mRootMask.array();
363  bool * pMaskEnd = pMask + mNumRoots;
364  bool Destroy = true;
365 
366  for (; pMask != pMaskEnd; ++pMask, ++pDiscrete)
367  {
368  if (*pMask)
369  {
370  if (*pDiscrete)
371  {
372  Destroy = false;
373  }
374  else
375  {
376  *pMask = false;
377  }
378  }
379  }
380 
381  if (Destroy)
382  {
383  destroyRootMask();
384  }
385  else
386  {
388  }
389 
390  // We have to restart the integrator
391  mLsodaStatus = 1;
392  }
393 
394  break;
395  }
396 
397  break;
398  }
399  }
400  else
401  {
402  mLSODA(&EvalF, // 1. evaluate F
403  &mData.dim, // 2. number of variables
404  mY, // 3. the array of current concentrations
405  &mTime, // 4. the current time
406  &EndTime, // 5. the final time
407  &ITOL, // 6. error control
408  &mRtol, // 7. relative tolerance array
409  mAtol.array(), // 8. absolute tolerance array
410  &mState, // 9. output by overshoot & interpolation
411  &mLsodaStatus, // 10. the state control variable
412  &one, // 11. further options (one)
413  mDWork.array(), // 12. the double work array
414  &DSize, // 13. the double work array size
415  mIWork.array(), // 14. the int work array
416  &ISize, // 15. the int work array size
417  EvalJ, // 16. evaluate J (not given)
418  &mJType); // 17. the type of jacobian calculate (2)
419  }
420 
421  // Why did we ignore this error?
422  // if (mLsodaStatus == -1) mLsodaStatus = 2;
423 
425 
426  if (!mPeekAheadMode)
427  {
429  }
430 
431  if ((mLsodaStatus <= 0))
432  {
433  Status = FAILURE;
434  mPeekAheadMode = false;
436  }
437 
438  if (!mMethodState.isValid())
439  {
440  Status = FAILURE;
441  mPeekAheadMode = false;
443  }
444 
445  return Status;
446 }
#define C_INT
Definition: copasi.h:115
C_INT mLsodaStatus
Definition: CLsodaMethod.h:107
static void EvalR(const C_INT *n, const C_FLOAT64 *t, const C_FLOAT64 *y, const C_INT *nr, C_FLOAT64 *r)
CVector< C_FLOAT64 > mAtol
Definition: CLsodaMethod.h:118
static void EvalF(const C_INT *n, const C_FLOAT64 *t, const C_FLOAT64 *y, C_FLOAT64 *ydot)
CVector< bool > mRootMask
Definition: CLsodaMethod.h:175
C_FLOAT64 mRtol
Definition: CLsodaMethod.h:113
CVector< C_FLOAT64 > mDWork
Definition: CLsodaMethod.h:143
CVector< bool > mDiscreteRoots
Definition: CLsodaMethod.h:180
C_FLOAT64 * mY
Definition: CLsodaMethod.h:86
unsigned C_INT32 mRootCounter
Definition: CLsodaMethod.h:199
static void EvalJ(const C_INT *n, const C_FLOAT64 *t, const C_FLOAT64 *y, const C_INT *ml, const C_INT *mu, C_FLOAT64 *pd, const C_INT *nRowPD)
void setTime(const C_FLOAT64 &time)
Definition: CState.cpp:326
CState mMethodState
Definition: CLsodaMethod.h:74
C_FLOAT64 mTargetTime
Definition: CLsodaMethod.h:193
std::ostringstream mErrorMsg
Definition: CLsodaMethod.h:123
CVector< C_INT > mRoots
#define MCTrajectoryMethod
C_INT mNumRoots
Definition: CLsodaMethod.h:96
void destroyRootMask()
bool isValid() const
Definition: CState.cpp:351
C_FLOAT64 mTime
Definition: CLsodaMethod.h:102
unsigned C_INT32 * mpMaxInternalSteps
Definition: CLsodaMethod.h:68
size_t size() const
Definition: CVector.h:100
const C_FLOAT64 & getTime() const
Definition: CState.cpp:325
#define C_FLOAT64
Definition: copasi.h:92
RootMasking mRootMasking
Definition: CLsodaMethod.h:186
bool mPeekAheadMode
Definition: CLsodaMethod.h:204
CType * array()
Definition: CVector.h:139
Status peekAhead()
void createRootMask()
CLSODAR mLSODAR
Definition: CLsodaMethod.h:133
CVector< C_INT > mIWork
Definition: CLsodaMethod.h:148

Friends And Related Function Documentation

Member Data Documentation

CVector< C_FLOAT64 > CLsodaMethod::mAtol
private

A vector of absolute tolerances.

Definition at line 118 of file CLsodaMethod.h.

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

Data CLsodaMethod::mData
private

mData.dim is the dimension of the ODE system. mData.pMethod contains CLsodaMethod * this to be used in the static method EvalF

Definition at line 81 of file CLsodaMethod.h.

Referenced by CLsodaMethod(), peekAhead(), start(), stateChanged(), and step().

CVector< bool > CLsodaMethod::mDiscreteRoots
private

A which indicates whether roots change only discretely.

Definition at line 180 of file CLsodaMethod.h.

Referenced by start(), and step().

C_FLOAT64 CLsodaMethod::mDummy
private

A dummy variable if we do not have any ODEs

Definition at line 170 of file CLsodaMethod.h.

Referenced by start().

CVector< C_FLOAT64 > CLsodaMethod::mDWork
private

LSODA C_FLOAT64 work area

Definition at line 143 of file CLsodaMethod.h.

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

std::ostringstream CLsodaMethod::mErrorMsg
private

Stream to capture LSODA error messages

Definition at line 123 of file CLsodaMethod.h.

Referenced by start(), and step().

CVector< C_INT > CLsodaMethod::mIWork
private

LSODA C_INT work area

Definition at line 148 of file CLsodaMethod.h.

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

C_INT CLsodaMethod::mJType
private

The way LSODA calculates the jacobian

Definition at line 153 of file CLsodaMethod.h.

Referenced by start(), and step().

CLSODA CLsodaMethod::mLSODA
private

The LSODA integrator

Definition at line 128 of file CLsodaMethod.h.

Referenced by start(), and step().

CLSODAR CLsodaMethod::mLSODAR
private

The LSODA integrator

Definition at line 133 of file CLsodaMethod.h.

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

C_INT CLsodaMethod::mLsodaStatus
protected

LSODA state.

Definition at line 107 of file CLsodaMethod.h.

Referenced by CTrajectoryMethodDsaLsodar::fireReaction(), start(), stateChanged(), and step().

CState CLsodaMethod::mMethodState
protected
bool CLsodaMethod::mNoODE
protected

A Boolean value indication whether we have no ODEs

Definition at line 164 of file CLsodaMethod.h.

Referenced by evalF(), CTrajectoryMethodDsaLsodar::evalF(), start(), and stateChanged().

C_INT CLsodaMethod::mNumRoots
private

Number of roots

Definition at line 96 of file CLsodaMethod.h.

Referenced by start(), and step().

C_FLOAT64* CLsodaMethod::mpAbsoluteTolerance
private

A pointer to the value of "Absolute Tolerance"

Definition at line 63 of file CLsodaMethod.h.

Referenced by createRootMask(), initializeParameter(), peekAhead(), and start().

bool CLsodaMethod::mPeekAheadMode
private

A Boolean indicating whether we are in peekAhead mode

Definition at line 204 of file CLsodaMethod.h.

Referenced by peekAhead(), start(), stateChanged(), and step().

unsigned C_INT32* CLsodaMethod::mpMaxInternalSteps
private

A pointer to the value of "Max Internal Steps"

Definition at line 68 of file CLsodaMethod.h.

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

CModel* CLsodaMethod::mpModel
protected
bool* CLsodaMethod::mpReducedModel
protected

A pointer to the value of "Integrate Reduced Model"

Definition at line 52 of file CLsodaMethod.h.

Referenced by createRootMask(), evalF(), CTrajectoryMethodDsaLsodar::evalF(), evalR(), initializeParameter(), CTrajectoryMethodDsaLsodar::integrateDeterministicPart(), and start().

C_FLOAT64* CLsodaMethod::mpRelativeTolerance
private

A pointer to the value of "Relative Tolerance"

Definition at line 58 of file CLsodaMethod.h.

Referenced by initializeParameter(), peekAhead(), and start().

unsigned C_INT32 CLsodaMethod::mRootCounter
private

Root counter to determine whether the internal step limit is exceeded.

Definition at line 199 of file CLsodaMethod.h.

Referenced by start(), and step().

CVector< bool > CLsodaMethod::mRootMask
private

A mask which hides all roots being constant and zero.

Definition at line 175 of file CLsodaMethod.h.

Referenced by createRootMask(), destroyRootMask(), maskRoots(), and step().

RootMasking CLsodaMethod::mRootMasking
protected

A Boolean flag indicating whether we should try masking roots

Definition at line 186 of file CLsodaMethod.h.

Referenced by createRootMask(), destroyRootMask(), evalR(), peekAhead(), and step().

C_FLOAT64 CLsodaMethod::mRtol
private

Relative tolerance.

Definition at line 113 of file CLsodaMethod.h.

Referenced by start(), and step().

C_INT CLsodaMethod::mState
private

The state of the integrator

Definition at line 138 of file CLsodaMethod.h.

Referenced by start(), and step().

C_FLOAT64 CLsodaMethod::mTargetTime
private

Store the targeted end time to determine whether the internal step limit is exceeded.

Definition at line 193 of file CLsodaMethod.h.

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

C_FLOAT64 CLsodaMethod::mTime
protected

Current time.

Definition at line 102 of file CLsodaMethod.h.

Referenced by CTrajectoryMethodDsaLsodar::fireReaction(), peekAhead(), start(), stateChanged(), and step().

C_FLOAT64* CLsodaMethod::mY
private

Pointer to the array with left hand side values.

Definition at line 86 of file CLsodaMethod.h.

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

CVector< C_FLOAT64 > CLsodaMethod::mYdot
private

Vector containing the derivatives after calling eval

Definition at line 91 of file CLsodaMethod.h.

Referenced by start().


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