COPASI API  4.16.103
CModel.cpp
Go to the documentation of this file.
1 // Copyright (C) 2010 - 2015 by Pedro Mendes, Virginia Tech Intellectual
2 // Properties, Inc., University of Heidelberg, and The University
3 // of Manchester.
4 // All rights reserved.
5 
6 // Copyright (C) 2008 - 2009 by Pedro Mendes, Virginia Tech Intellectual
7 // Properties, Inc., EML Research, gGmbH, University of Heidelberg,
8 // and The University of Manchester.
9 // All rights reserved.
10 
11 // Copyright (C) 2001 - 2007 by Pedro Mendes, Virginia Tech Intellectual
12 // Properties, Inc. and EML Research, gGmbH.
13 // All rights reserved.
14 
15 //
16 
17 #include <sstream>
18 #include <string>
19 #include <vector>
20 #include <limits>
21 #include <cmath>
22 #include <algorithm>
23 
24 #include "copasi.h"
25 
26 // #define DEBUG_MATRIX
27 
28 #include "CCompartment.h"
29 #include "CMetab.h"
30 #include "CModel.h"
31 #include "CState.h"
32 #include "CModelValue.h"
33 #include "CObjectLists.h"
34 #include "function/CFunctionDB.h"
36 #include "report/CKeyFactory.h"
42 #include "utilities/CVector.h"
43 #include "utilities/CluX.h"
44 #include "utilities/utility.h"
46 #include "CReactionInterface.h"
48 #include "CMetabNameInterface.h"
49 #include "CMathModel.h"
50 
51 #ifdef USE_MATH_CONTAINER
52 # include "math/CMathContainer.h"
53 #endif //USE_MATH_CONTAINER
54 
55 #include "lapack/blaswrap.h"
56 #include "lapack/lapackwrap.h"
57 
58 #ifdef COPASI_DEBUG
59 #define CCHECK {check();}
60 #else
61 #define CCHECK
62 #endif
63 
64 #define MNumMetabolitesReactionDependent (mNumMetabolitesReaction - mNumMetabolitesReactionIndependent)
65 
66 const char * CModel::VolumeUnitNames[] =
67 {"dimensionless", "m\xc2\xb3", "l", "ml", "\xc2\xb5l", "nl", "pl", "fl", NULL};
68 
69 const char * CModel::AreaUnitNames[] =
70 {"dimensionless", "m\xc2\xb2", "dm\xc2\xb2", "cm\xc2\xb2", "mm\xc2\xb2", "\xc2\xb5m\xc2\xb2", "nm\xc2\xb2", "pm\xc2\xb2", "fm\xc2\xb2", NULL};
71 
72 const char * CModel::LengthUnitNames[] =
73 {"dimensionless", "m", "dm", "cm", "mm", "\xc2\xb5m", "nm", "pm", "fm", NULL};
74 
75 const char * CModel::TimeUnitNames[] =
76 {"dimensionless", "d", "h", "min", "s", "ms", "\xc2\xb5s", "ns", "ps", "fs", NULL};
77 
78 // "mol" is the correct name, however in the COPASI XML files "Mol" is used
79 // up to build 18
80 
81 const char * CModel::QuantityUnitOldXMLNames[] =
82 {"dimensionless", "Mol", "mMol", "\xc2\xb5Mol", "nMol", "pMol", "fMol", "#", NULL};
83 
84 const char * CModel::QuantityUnitNames[] =
85 {"dimensionless", "mol", "mmol", "\xc2\xb5mol", "nmol", "pmol", "fmol", "#", NULL};
86 
87 const char * CModel::ModelTypeNames[] =
88 {"deterministic", "stochastic", NULL};
89 
91  CModelEntity("New Model", pParent, "Model"),
92  mInitialState(),
93  mCurrentState(),
94  mStateTemplate(*this, this->mInitialState, this->mCurrentState),
95  mSimulatedUpToDateObjects(),
96  mInitialDependencies(),
97  mTransientDependencies(),
98  mPhysicalDependencies(),
99 #ifdef USE_MATH_CONTAINER
100  mpMathContainer(NULL),
101 #endif // USE_MATH_CONTAINER
102  mVolumeUnit(ml),
103  mAreaUnit(m2),
104  mLengthUnit(m),
105  mTimeUnit(s),
106  mQuantityUnit(mMol),
107  mType(deterministic),
108  mCompartments("Compartments", this),
109  mMetabolites("Metabolites", this),
110  mMetabolitesX("Reduced Model Metabolites", this),
111  mSteps("Reactions", this),
112  mEvents("Events", this),
113  mParticleFluxes(),
114  mValues("Values", this),
115  mParameterSet("Initial State", this),
116  mParameterSets("ParameterSets", this),
117  mActiveParameterSetKey(""),
118  mMoieties("Moieties", this),
119  mStoiInternal(),
120  mpStoiAnnotation(NULL),
121  mStoi(),
122  mRedStoi(),
123  mpRedStoiAnnotation(NULL),
124  mNumMetabolitesUnused(0),
125  mNumMetabolitesODE(0),
126  mNumMetabolitesReaction(0),
127  mNumMetabolitesAssignment(0),
128  mNumMetabolitesReactionIndependent(0),
129  mL(),
130  mpLinkMatrixAnnotation(NULL),
131  mLView(mL),
132  mAvogadro(6.02214179e23),
133  mQuantity2NumberFactor(1.0),
134  mNumber2QuantityFactor(1.0),
135  mpCompileHandler(NULL),
136  mInitialRefreshes(),
137  mSimulatedRefreshes(),
138  mApplyInitialValuesRefreshes(),
139  mNonSimulatedRefreshes(),
140  mReorderNeeded(false),
141  mIsAutonomous(true),
142  mBuildInitialSequence(true),
143  mpMathModel(NULL)
144 {
145  initObjects();
146 
147  setStatus(TIME);
148  setUsed(true);
149 
150  *mpIValue = 0.0;
151  *mpValue = std::numeric_limits<C_FLOAT64>::quiet_NaN();
152 
153  size_t i, imax = mSteps.size();
154 
155  for (i = 0; i < imax; i++)
156  mSteps[i]->compile(/*mCompartments*/);
157 
159 
160  forceCompile(NULL);
161 
162  /* This following 2 lines added by Liang Xu
163  Because of the failure to initialize the parameter when creating a new models
164  */
165  setQuantityUnit(mQuantityUnit); // set the factors
166  setVolumeUnit(mVolumeUnit); // set the factors
167 
169 }
170 
171 // CModel::CModel(const CModel & src):
172 // CModelEntity(src),
173 // mInitialState(),
174 // mCurrentState(),
175 // mStateTemplate(*this, this->mInitialState, this->mCurrentState),
176 // mComments(src.mComments),
177 // mVolumeUnit(src.mVolumeUnit),
178 // mTimeUnit(src.mTimeUnit),
179 // mQuantityUnit(src.mQuantityUnit),
180 // mType(src.mType),
181 // mCompartments(src.mCompartments, this),
182 // mMetabolites(src.mMetabolites, this),
183 // mMetabolitesX(src.mMetabolitesX, this),
184 // mSteps(src.mSteps, this),
185 // mEvents(src.mEvents,this),
186 // mParticleFluxes(src.mParticleFluxes),
187 // mValues(src.mValues, this),
188 // mParameterSets(src.mParameterSets, this),
189 // mActiveParameterSetKey(src.mActiveParameterSetKey),
190 // mMoieties(src.mMoieties, this),
191 // mStoiInternal(src.mStoiInternal),
192 // mpStoiAnnotation(NULL),
193 // mStoi(src.mStoi),
194 // mRedStoi(src.mRedStoi),
195 // mpRedStoiAnnotation(NULL),
196 // mNumMetabolitesUnused(src.mNumMetabolitesUnused),
197 // mNumMetabolitesODE(src.mNumMetabolitesODE),
198 // mNumMetabolitesReaction(src.mNumMetabolitesReaction),
199 // mNumMetabolitesAssignment(src.mNumMetabolitesAssignment),
200 // mNumMetabolitesIndependent(src.mNumMetabolitesIndependent),
201 // mL(src.mL),
202 // mpLinkMatrixAnnotation(NULL),
203 // mLView(mL, mNumMetabolitesIndependent),
204 // mQuantity2NumberFactor(src.mQuantity2NumberFactor),
205 // mNumber2QuantityFactor(src.mNumber2QuantityFactor),
206 // mpCompileHandler(NULL),
207 // mInitialRefreshes(),
208 // mSimulatedRefreshes(),
209 // mConstantRefreshes(),
210 // mNonSimulatedRefreshes(),
211 // mReorderNeeded(false),
212 // mIsAutonomous(false)
213 // {
214 // CONSTRUCTOR_TRACE;
215 // initObjects();
216 //
217 // size_t i, imax = mSteps.size();
218 //
219 // for (i = 0; i < imax; i++)
220 // mSteps[i]->compile(/*mCompartments*/);
221 //
222 // initializeMetabolites();
223 //
224 // forceCompile(NULL);
225 //}
226 
228 {
229  mpIValue = NULL;
230  mpValue = NULL;
231 
235 
236 #ifdef USE_MATH_CONTAINER
237  pdelete(mpMathContainer);
238 #endif // USE_MATH_CONTAINER
239 
241 
243 }
244 
245 // virtual
246 std::string CModel::getChildObjectUnits(const CCopasiObject * pObject) const
247 {
248  if (pObject->getObjectName() == "Initial Time" ||
249  pObject->getObjectName() == "Time")
250  return getTimeUnitName();
251 
252  return "";
253 }
254 
256 {
257  C_INT32 Size = 0;
258  C_INT32 Fail = 0;
259  size_t i;
260  std::string tmp;
261 
262  // For old Versions we must read the list of Metabolites beforehand
263  if ((Fail = configBuffer.getVariable("TotalMetabolites", "C_INT32",
264  &Size, CReadConfig::LOOP)))
265  return Fail;
266 
267  // :TODO: Remove OldMetabolites as part of the data model.
268  CCopasiDataModel* pDataModel = getObjectDataModel();
269  assert(pDataModel != NULL);
270  pDataModel->pOldMetabolites->load(configBuffer, Size);
271 
272  if ((Fail = configBuffer.getVariable("Title", "string", &tmp,
274  return Fail;
275 
276  setObjectName(tmp);
277 
278  std::string Notes;
279 
280  try
281  {
282  Fail = configBuffer.getVariable("Comments", "multiline", &Notes,
284  }
285 
286  catch (CCopasiException & Exception)
287  {
288  if ((MCReadConfig + 1) == Exception.getMessage().getNumber())
289  Notes = "";
290  else
291  throw Exception;
292  }
293 
294  setNotes(Notes);
295 
296  try
297  {
298  Fail = configBuffer.getVariable("TimeUnit", "string", &tmp,
300  }
301  catch (CCopasiException & Exception)
302  {
303  if ((MCReadConfig + 1) == Exception.getMessage().getNumber())
304  tmp = ""; //unknown?
305  else
306  throw Exception;
307  }
308 
309  setTimeUnit(tmp); // set the factors
310 
311  try
312  {
313  Fail = configBuffer.getVariable("ConcentrationUnit", "string", &tmp,
315  }
316  catch (CCopasiException & Exception)
317  {
318  if ((MCReadConfig + 1) == Exception.getMessage().getNumber())
319  tmp = ""; //unknown?
320  else
321  throw Exception;
322  }
323 
324  setQuantityUnit(tmp); // set the factors
325 
326  try
327  {
328  Fail = configBuffer.getVariable("VolumeUnit", "string", &tmp,
330  }
331  catch (CCopasiException & Exception)
332  {
333  if ((MCReadConfig + 1) == Exception.getMessage().getNumber())
334  tmp = ""; //unknown?
335  else
336  throw Exception;
337  }
338 
339  setVolumeUnit(tmp); // set the factors
340 
341  *mpIValue = 0;
342 
343  if ((Fail = configBuffer.getVariable("TotalCompartments", "C_INT32", &Size,
345  return Fail;
346 
347  mCompartments.load(configBuffer, Size);
348 
349  // Create the correct compartment / metabolite relationships
350  CMetab *pMetabolite;
351 
352  for (i = 0; i < pDataModel->pOldMetabolites->size(); i++)
353  {
354  pMetabolite = new CMetab;
355  mCompartments[(*pDataModel->pOldMetabolites)[i]->getIndex()]->
356  addMetabolite(pMetabolite);
357 
358  (*pMetabolite) = *(*pDataModel->pOldMetabolites)[i];
359  }
360 
362 
363  if ((Fail = CCopasiRootContainer::getFunctionList()->load(configBuffer))) // slow
364  return Fail;
365 
366  if ((Fail = configBuffer.getVariable("TotalSteps", "C_INT32", &Size,
368  return Fail;
369 
370  mSteps.load(configBuffer, Size); // slow
371 
372  for (i = 0; i < mSteps.size(); i++)
373  mSteps[i]->compile(/*mCompartments*/);
374 
375  pDataModel->pOldMetabolites->cleanup();
376 
377  setCompileFlag();
378  return Fail;
379 }
380 
382 {
383  bool success = true;
384 
386  {
388  }
389 
391 
393 
394  unsigned C_INT32 CompileStep = 0;
395  size_t hCompileStep;
396 
397  if (mpCompileHandler)
398  {
399  mpCompileHandler->setName("Compiling model...");
400  unsigned C_INT32 totalSteps = 7;
401  hCompileStep = mpCompileHandler->addItem("Compile Process",
402  CompileStep,
403  &totalSteps);
404  }
405 
406  // To assure that we do not produce access violations we clear the refresh sequences
407  // first
408  mInitialRefreshes.clear();
409  mSimulatedRefreshes.clear();
411  mNonSimulatedRefreshes.clear();
412 
413  CompileStep = 0;
414 
415  if (mpCompileHandler && !mpCompileHandler->progressItem(hCompileStep))
416  {
417  success = false;
418  goto finish;
419  }
420 
421  buildStoi();
422  CompileStep = 1;
423 
424  if (mpCompileHandler && !mpCompileHandler->progressItem(hCompileStep))
425  {
426  success = false;
427  goto finish;
428  }
429 
430  buildLinkZero();
431  CompileStep = 2;
432 
433  if (mpCompileHandler && !mpCompileHandler->progressItem(hCompileStep))
434  {
435  success = false;
436  goto finish;
437  }
438 
439  buildRedStoi();
440  CompileStep = 3;
441 
442  if (mpCompileHandler && !mpCompileHandler->progressItem(hCompileStep))
443  {
444  success = false;
445  goto finish;
446  }
447 
448  buildMoieties();
449  CompileStep = 4;
450 
451  if (mpCompileHandler && !mpCompileHandler->progressItem(hCompileStep))
452  {
453  success = false;
454  goto finish;
455  }
456 
458  CompileStep = 5;
459 
460  if (mpCompileHandler && !mpCompileHandler->progressItem(hCompileStep))
461  {
462  success = false;
463  goto finish;
464  }
465 
466  try
467  {
468  success &= buildInitialSequence();
469  success &= buildApplyInitialValuesSequence();
470  success &= buildSimulatedSequence();
471  success &= buildNonSimulatedSequence();
472  }
473  catch (...)
474  {
475  success = false;
476  }
477 
478  CompileStep = 6;
479 
480  if (mpCompileHandler && !mpCompileHandler->progressItem(hCompileStep))
481  {
482  success = false;
483  goto finish;
484  }
485 
486  buildUserOrder();
487 
488  if (mpCompileHandler) mpCompileHandler->finishItem(hCompileStep);
489 
490  {
493 
494  for (; itSpecies != endSpecies; ++itSpecies)
495  {
496  (*itSpecies)->compileIsInitialConcentrationChangeAllowed();
497  }
498  }
499 
500  //update annotations
502 
503  success &= compileEvents();
504  success &= mpMathModel->compile(this);
505 
506  if (!success)
507  {
508  mIsAutonomous = false;
509  }
510  else
511  {
512  mCompileIsNecessary = false;
514  }
515 
516  //writeDependenciesToDotFile();
517 
519 
520 #ifdef USE_MATH_CONTAINER
521  pdelete(mpMathContainer);
522  mpMathContainer = new CMathContainer(*this);
523 
524  // CMathContainer CopyModel(MathModel);
525 #endif // USE_MATH_CONTAINER
526 
527  // Update the parameter set
529 
530 finish:
531 
532  // Since we have applied the pivot to the stoichiometry matrix and the species
533  // we do not need them any longer. In fact it is detrimental if other functions rely
534  // on consistency between the stoichiometry matrix, reduced stoichiometry matrix and the Link matrix..
535  mL.clearPivoting();
536 
538  {
540  }
541 
542  mCompileIsNecessary = !success;
543 
544  return success;
545 }
546 
548 {
552 
553  // The initial values of the model entities
554  // We need to add the time for non-autonomous models.
555  CModelEntity **ppEntity = mStateTemplate.beginIndependent() - 1;
556  CModelEntity **ppEntityEnd = mStateTemplate.endFixed();
557 
558  for (; ppEntity != ppEntityEnd; ++ppEntity)
559  {
560  CMetab * pSpecies = dynamic_cast< CMetab * >(*ppEntity);
561 
562  mInitialDependencies.addObject((*ppEntity)->getInitialValueReference());
563  mTransientDependencies.addObject((*ppEntity)->getValueReference());
564 
565  if (pSpecies != NULL)
566  {
569  }
570 
571  switch ((*ppEntity)->getStatus())
572  {
573  case CModelEntity::ODE:
575  mTransientDependencies.addObject((*ppEntity)->getRateReference());
576  break;
577 
578  default:
579  break;
580  }
581  }
582 
585 
586  for (; itMoiety != endMoiety; ++itMoiety)
587  {
588  mInitialDependencies.addObject((*itMoiety)->getInitialValueReference());
589  mTransientDependencies.addObject((*itMoiety)->getTotalNumberReference());
590  mTransientDependencies.addObject((*itMoiety)->getDependentNumberReference());
591  }
592 
593 #ifdef COPASI_DEBUG_TRACE
594  std::ofstream File;
595 
596  File.open("InitialDependencies.dot");
597  mInitialDependencies.exportDOTFormat(File, "InitialDependencies");
598  File.close();
599 
600  File.open("SimulationDependencies.dot");
601  mTransientDependencies.exportDOTFormat(File, "SimulationDependencies");
602  File.close();
603 #endif //COPASI_DEBUG_TRACE
604  return true;
605 }
606 
608 {
610 }
611 
613 {
614  bool success = true;
615 
617  {
618  mpCompileHandler = pProcessReport;
619 
620  try
621  {
622  success &= compile();
623  }
624 
625  catch (...)
626  {
627  success = false;
628  }
629 
630  mpCompileHandler = NULL;
631  }
632 
633  return success;
634 }
635 
637 {
638  setCompileFlag();
639  return compileIfNecessary(pProcessReport);
640 }
641 
643 {
644  unsigned C_INT32 i, numCols;
645 
647 
648  size_t numRows;
649  numRows = mNumMetabolitesReaction;
650  numCols = (unsigned C_INT32) mSteps.size();
651 
652  mParticleFluxes.resize(numCols);
653  mStoiInternal.resize(numRows, numCols);
654  mStoiInternal = 0.0;
655 
656  size_t hProcess;
657 
658  if (mpCompileHandler)
659  {
660  i = 0;
661  hProcess = mpCompileHandler->addItem("Building Stoichiometry",
662  i,
663  &numCols);
664  }
665 
666  C_FLOAT64 * pCol, *pColEnd;
667  pCol = mStoiInternal.array();
668  pColEnd = mStoiInternal.array() + numCols;
669 
670  C_FLOAT64 * pRow, *pRowEnd;
671  pRowEnd = mStoiInternal.array() + numRows * numCols;
672 
675 
676  for (; pCol < pColEnd; ++pCol, ++itStep)
677  {
678  if (mpCompileHandler && !mpCompileHandler->progressItem(hProcess)) return;
679 
680  // Since we are stepping through the reactions we can check whether
681  // the kinetic functions are usable.
682  if (!(*itStep)->getFunction()->isUsable())
684  (*itStep)->getObjectName().c_str(),
685  (*itStep)->getFunction()->getObjectName().c_str());
686 
687  const CCopasiVector< CChemEqElement > & Balance =
688  (*itStep)->getChemEq().getBalances();
691 
692  for (; itBalance != endBalance; ++itBalance)
693  {
694  const std::string & key = (*itBalance)->getMetaboliteKey();
695 
696  for (pRow = pCol, itMetab = mMetabolitesX.begin() + mNumMetabolitesODE;
697  pRow < pRowEnd;
698  pRow += numCols, ++itMetab)
699  if ((*itMetab)->getKey() == key)
700  {
701  *pRow = (*itBalance)->getMultiplicity();
702  break;
703  }
704  }
705  }
706 
707  // We need to have all unused and fixed metabolites at the end of mMetabolites.
708  // However we can only detect unused metabolites after building the
709  // stoichiometry matrix.
711 
712 #ifdef DEBUG_MATRIX
713  DebugFile << "Stoichiometry Matrix" << std::endl;
714  DebugFile << mStoiInternal << std::endl;
715 #endif
716 
717  if (mpCompileHandler)
718  mpCompileHandler->finishItem(hProcess);
719 
720  return;
721 }
722 
724 {
725  size_t numRows, numCols;
726  numRows = mStoiInternal.numRows();
727  numCols = mStoiInternal.numCols();
728 
729  C_FLOAT64 * pStoi, *pStoiEnd, *pRowEnd;
730  pStoi = mStoiInternal.array();
731  pStoiEnd = mStoiInternal.array() + numRows * numCols;
732 
733  size_t i, NumUnused;
734  C_FLOAT64 tmp;
735  std::vector< size_t > Unused;
736 
737  for (i = 0; i < numRows; i++)
738  {
739  tmp = 0;
740 
741  for (pRowEnd = pStoi + numCols; pStoi < pRowEnd; ++pStoi)
742  tmp += fabs(*pStoi);
743 
744  if (tmp < std::numeric_limits< C_FLOAT64 >::min()) Unused.push_back(i);
745  }
746 
747  NumUnused = Unused.size();
748 
749  if (NumUnused == 0) return false;
750 
751  // We treat unused variables in the same way as fixed, i.e.
752  // they will be sorted to the end of the metabolite list.
753 
754  numRows -= NumUnused;
755 
756  CMatrix< C_FLOAT64 > NewStoi(numRows, numCols);
757  C_FLOAT64 * pNewStoi = NewStoi.array();
758  std::vector< CMetab * > UsedMetabolites(numRows);
759  std::vector< CMetab * >::iterator itUsedMetabolites = UsedMetabolites.begin();
760  std::vector< CMetab * > UnusedMetabolites(NumUnused);
761  std::vector< CMetab * >::iterator itUnusedMetabolites = UnusedMetabolites.begin();
762  std::vector< size_t >::const_iterator itUnused = Unused.begin();
763  std::vector< size_t >::const_iterator endUnused = Unused.end();
764 
767 
768  // Build new stoichiometry Matrix
769  pStoi = mStoiInternal.array();
770 
771  for (i = 0; itMetab != endMetab; ++itMetab, i++, pStoi += numCols)
772  {
773  if (itUnused != endUnused && i == *itUnused)
774  {
775  (*itMetab)->setUsed(false);
776  *itUnusedMetabolites = (*itMetab);
777 
778  ++itUnusedMetabolites;
779  ++itUnused;
780  }
781  else
782  {
783  (*itMetab)->setUsed(true);
784  *itUsedMetabolites = (*itMetab);
785  ++itUsedMetabolites;
786 
787  // The row needs to be copied to the new stoichiometry matrix
788  memcpy(pNewStoi, pStoi, sizeof(C_FLOAT64) * numCols);
789  pNewStoi += numCols;
790  }
791  }
792 
793  // Reorder metabolites
794  // Skip the metabolites determined by ODE
795  itMetab = mMetabolitesX.begin() + mNumMetabolitesODE;
796 
797  // Handle the metabolites determined by actually determined by reactions
798  itUsedMetabolites = UsedMetabolites.begin();
799  std::vector< CMetab * >::iterator itMetabolitesEnd = UsedMetabolites.end();
800 
801  for (; itUsedMetabolites != itMetabolitesEnd; ++itUsedMetabolites, ++itMetab)
802  *itMetab = *itUsedMetabolites;
803 
804  // Handle metabolites determined by assignment and marked as fixed
805  // This is just a shift of NumUnused.
806  endMetab = itMetab + mNumMetabolitesAssignment + mNumMetabolitesUnused;
807 
808  for (; itMetab != endMetab; ++itMetab)
809  *itMetab = *(itMetab + NumUnused);
810 
811  // Handle newly marked unused metabolites
812  itUnusedMetabolites = UnusedMetabolites.begin();
813  itMetabolitesEnd = UnusedMetabolites.end();
814 
815  for (; itUnusedMetabolites != itMetabolitesEnd; ++itUnusedMetabolites, ++itMetab)
816  *itMetab = *itUnusedMetabolites;
817 
818  // Now its time to update the number of metabolites determined by reactions
819  // and the one being unused.
820  mNumMetabolitesReaction -= NumUnused;
821  mNumMetabolitesUnused += NumUnused;
822 
823  // Update stoichiometry matrix
824  mStoiInternal = NewStoi;
825 
826  return true;
827 }
828 
830 {
831  mRedStoi = mStoi;
833 
834  // The first metabolites are determined by ODEs we therefore cannot simply
835  // apply the pivot.
836 
837  // Create a temporary copy of metabolites determined by reactions.
838  CVector< CMetab * > ReactionMetabolites(mNumMetabolitesReaction);
839  CMetab ** ppMetab = ReactionMetabolites.array();
840  CMetab ** ppMetabEnd = ppMetab + mNumMetabolitesReaction;
842 
843  for (; ppMetab != ppMetabEnd; ++ppMetab, ++itMetabX)
844  {
845  *ppMetab = *itMetabX;
846  }
847 
848  // Apply the pivot on the temporary copy
849  mL.applyRowPivot(ReactionMetabolites);
850 
851  // Map the result on the actual metabolites
852  ppMetab = ReactionMetabolites.array();
853  itMetabX = mMetabolitesX.begin() + mNumMetabolitesODE;
854 
855  for (; ppMetab != ppMetabEnd; ++ppMetab, ++itMetabX)
856  {
857  *itMetabX = *ppMetab;
858  }
859 
860  return;
861 }
862 
864 {
868 
871 
873  size_t j;
874 
875  for (j = 0; it != end; ++it, j++)
876  {
877  CN = (*it)->getCN();
878 
883  }
884 
886 
887  for (; it != end; ++it, j++)
888  {
889  CN = (*it)->getCN();
890 
893  }
894 }
895 
897 {
900 
901  for (; it != end; ++it)
902  (*it)->refreshInitialValue();
903 }
904 
906 {
907  // Independent metabs
910 
911  // Dependent metabs
912  CCopasiVector< CMetab >::iterator itDependent = endIndependent;
914 
915  C_FLOAT64 * pFactor = mL.array();
916 
917  CMoiety *pMoiety;
918  mMoieties.cleanup();
919 
920  for (; itDependent != endDependent; ++itDependent)
921  {
922  pMoiety = new CMoiety((*itDependent)->getObjectName());
923  pMoiety->add(1.0, *itDependent);
924 
925  if (pFactor != NULL)
926  {
927  for (itIndependent = mMetabolitesX.begin() + mNumMetabolitesODE; itIndependent != endIndependent; ++itIndependent, ++pFactor)
928  if (fabs(*pFactor) > std::numeric_limits< C_FLOAT64 >::epsilon())
929  pMoiety->add(- *pFactor, *itIndependent);
930  }
931 
932  mMoieties.add(pMoiety, true);
933  }
934 
936  return;
937 }
938 
939 //this is supposed to be so fast it can be called often to be kept up to date
940 //all the time. At the moment it creates the mMetabolites and sorts the fixed
941 //metabs to the end
943 {
944  // Create a vector of pointers to all metabolites.
945  // Note, the metabolites physically exist in the compartments.
947 
952 
953  for (; itCompartment != endCompartment; ++itCompartment)
954  {
955  itMetab = (*itCompartment)->getMetabolites().begin();
956  endMetab = (*itCompartment)->getMetabolites().end();
957 
958  for (; itMetab != endMetab; ++itMetab)
959  {
960  // Reset all moieties
961  (*itMetab)->setDependentOn(NULL);
962  mMetabolites.add(*itMetab);
963  }
964  }
965 
966  // We sort the metabolites by type
967  itMetab = mMetabolites.begin();
968  endMetab = mMetabolites.end();
969 
970  std::vector< CMetab *> ODEMetabs;
971  std::vector< CMetab *> ReactionMetabs;
972  std::vector< CMetab *> AssignmentMetabs;
973  std::vector< CMetab *> FixedMetabs;
974 
975  for (; itMetab != endMetab; ++itMetab)
976  switch ((*itMetab)->getStatus())
977  {
978  case FIXED:
979  FixedMetabs.push_back(*itMetab);
980  (*itMetab)->setUsed(false);
981  break;
982 
983  case ASSIGNMENT:
984  AssignmentMetabs.push_back(*itMetab);
985  (*itMetab)->setUsed(true);
986  break;
987 
988  case ODE:
989  ODEMetabs.push_back(*itMetab);
990  (*itMetab)->setUsed(true);
991  break;
992 
993  case REACTIONS:
994  ReactionMetabs.push_back(*itMetab);
995  (*itMetab)->setUsed(true);
996  break;
997 
998  default:
999  fatalError();
1000  break;
1001  }
1002 
1003  // Update mMetabolitesX to reflect the reordering.
1004  // We need to to this to allow the use of the full model for simulation.
1006 
1007  mNumMetabolitesODE = ODEMetabs.size();
1008  itMetab = mMetabolitesX.begin();
1009  std::vector< CMetab *>::const_iterator itSorted = ODEMetabs.begin();
1010  std::vector< CMetab *>::const_iterator endSorted = ODEMetabs.end();
1011 
1012  for (; itSorted != endSorted; ++itSorted, ++itMetab)
1013  *itMetab = *itSorted;
1014 
1015  mNumMetabolitesReaction = ReactionMetabs.size();
1016  itSorted = ReactionMetabs.begin();
1017  endSorted = ReactionMetabs.end();
1018 
1019  for (; itSorted != endSorted; ++itSorted, ++itMetab)
1020  *itMetab = *itSorted;
1021 
1022  mNumMetabolitesAssignment = AssignmentMetabs.size();
1023  itSorted = AssignmentMetabs.begin();
1024  endSorted = AssignmentMetabs.end();
1025 
1026  for (; itSorted != endSorted; ++itSorted, ++itMetab)
1027  *itMetab = *itSorted;
1028 
1029  mNumMetabolitesUnused = FixedMetabs.size();
1030  itSorted = FixedMetabs.begin();
1031  endSorted = FixedMetabs.end();
1032 
1033  for (; itSorted != endSorted; ++itSorted, ++itMetab)
1034  *itMetab = *itSorted;
1035 }
1036 
1037 //**********************************************************************
1038 
1040 {return mSteps;}
1041 
1043 {return mSteps;}
1044 
1046 {return mParticleFluxes;}
1047 
1049 {return mMetabolites;}
1050 
1052 {return mMetabolites;}
1053 
1055 {CCHECK return mMetabolitesX;}
1056 
1058 {CCHECK return mMetabolitesX;}
1059 
1061 {return mValues;}
1062 
1064 {return mValues;}
1065 
1067 {return mParameterSets;}
1068 
1070 {return mParameterSets;}
1071 
1073 {
1074  return mParameterSet;
1075 }
1076 
1078 {
1079  return mParameterSet;
1080 }
1081 
1083 {
1084  CModelParameterSet * pParameterSet =
1086 
1087  if (pParameterSet != NULL)
1088  {
1089  pParameterSet->updateModel();
1090  }
1091  else
1092  {
1093  /*
1094  CModelParameterSet * pParameterSet = new CModelParameterSet(UTCTimeStamp());
1095  mParameterSets.add(pParameterSet, true);
1096  mActiveParameterSetKey = pParameterSet->getKey();
1097  pParameterSet->createFromModel();
1098  */
1099  }
1100 
1103 }
1104 
1106 {
1108 }
1109 
1111 {return mEvents;}
1112 
1114 {return mEvents;}
1115 
1116 //********
1117 
1118 size_t CModel::getNumMetabs() const
1119 {return mMetabolites.size();}
1120 
1123 
1125 {CCHECK return mNumMetabolitesODE;}
1126 
1129 
1132 
1135 
1136 size_t CModel::getTotSteps() const
1137 {return mSteps.size();}
1138 
1140 {return mValues.size();}
1141 
1142 const std::string & CModel::getKey() const
1143 {return mKey;}
1144 
1146 {return mCompartments;}
1147 
1149 {return mCompartments;}
1150 
1151 /**
1152  * Get the Reduced Stoichiometry Matrix of this Model
1153  */
1155 {CCHECK return mRedStoi;}
1156 
1157 /**
1158  * Get the reordered stoichiometry matrix of this model
1159  */
1161 {CCHECK return mStoi;}
1162 
1164 {return mMoieties;}
1165 
1167 {CCHECK return mLView;}
1168 
1170 {return mL;}
1171 
1173 {CCHECK return mStateTemplate;}
1174 
1176 {CCHECK return mStateTemplate;}
1177 
1178 const std::set< const CCopasiObject * > & CModel::getUptoDateObjects() const
1179 {return mSimulatedUpToDateObjects;}
1180 
1182 {*mpIValue = time;}
1183 
1185 {return *mpIValue;}
1186 
1187 void CModel::setTime(const C_FLOAT64 & time)
1188 {*mpValue = time;}
1189 
1191 {return *mpValue;}
1192 
1193 //**********************************************************************
1194 
1195 /**
1196  * Returns the index of the metab
1197  */
1198 size_t CModel::findMetabByName(const std::string & name) const
1199 {
1200  size_t i, imax = mMetabolites.size();
1202 
1203  std::string Name = unQuote(name);
1204 
1205  for (i = 0; i < imax; i++, Target++)
1206  if (*Target &&
1207  ((*Target)->getObjectName() == name ||
1208  (*Target)->getObjectName() == Name)) return i;
1209 
1210  return C_INVALID_INDEX;
1211 }
1212 
1213 /**
1214  * Returns the index of the Moiety
1215  */
1216 size_t CModel::findMoiety(const std::string &Target) const
1217 {
1218  size_t i, s;
1219  std::string name;
1220 
1221  s = mMoieties.size();
1222 
1223  for (i = 0; i < s; i++)
1224  {
1225  name = mMoieties[i]->getObjectName();
1226 
1227  if (name == Target)
1228  return i;
1229  }
1230 
1231  return C_INVALID_INDEX;
1232 }
1233 
1234 //**********************************************************************
1235 
1237 {
1238  // Copy the initial state to the current state,
1240 
1241  // Since the initial state is in itself consistent we should not need to
1242  // do anything further. However, for species of type ODE and ASSIGNMENT
1243  // the effective state variable is the concentration, i.e., we need to update
1244  // their concentration here.
1245  std::vector< Refresh * >::const_iterator itRefresh = mApplyInitialValuesRefreshes.begin();
1246  std::vector< Refresh * >::const_iterator endRefresh = mApplyInitialValuesRefreshes.end();
1247 
1248  while (itRefresh != endRefresh)
1249  (**itRefresh++)();
1250 
1251  // Update all dependent objects needed for simulation.
1252  updateSimulatedValues(false);
1253 
1255 }
1256 
1258 {
1259  mMoieties.clear();
1260 }
1261 
1263 {
1265  CModelEntity ** ppEntity = Entities.array();
1266 
1269 
1270  for (; itValue != endValue; ++itValue)
1271  if ((*itValue)->getStatus() == ODE)
1272  {
1273  *ppEntity = *itValue;
1274  (*ppEntity++)->setUsed(true);
1275  }
1276 
1279 
1280  for (; itCompartment != endCompartment; ++itCompartment)
1281  if ((*itCompartment)->getStatus() == ODE)
1282  {
1283  *ppEntity = *itCompartment;
1284  (*ppEntity++)->setUsed(true);
1285  }
1286 
1289 
1290  for (; itMetab != endMetab; ++itMetab)
1291  {
1292  if (!(*itMetab)->isUsed()) break;
1293 
1294  *ppEntity++ = *itMetab;
1295  }
1296 
1297  itCompartment = mCompartments.begin();
1298 
1299  for (; itCompartment != endCompartment; ++itCompartment)
1300  if ((*itCompartment)->getStatus() == ASSIGNMENT)
1301  {
1302  *ppEntity = *itCompartment;
1303  (*ppEntity++)->setUsed(true);
1304  }
1305 
1306  itValue = mValues.begin();
1307 
1308  for (; itValue != endValue; ++itValue)
1309  if ((*itValue)->getStatus() == ASSIGNMENT)
1310  {
1311  *ppEntity = *itValue;
1312  (*ppEntity++)->setUsed(true);
1313  }
1314 
1315  for (; itMetab != endMetab; ++itMetab)
1316  *ppEntity++ = *itMetab;
1317 
1318  itCompartment = mCompartments.begin();
1319 
1320  for (; itCompartment != endCompartment; ++itCompartment)
1321  if ((*itCompartment)->isFixed())
1322  *ppEntity++ = *itCompartment;
1323 
1324  itValue = mValues.begin();
1325 
1326  for (; itValue != endValue; ++itValue)
1327  if ((*itValue)->isFixed())
1328  *ppEntity++ = *itValue;
1329 
1330  mStateTemplate.reorder(Entities);
1331  mReorderNeeded = false;
1332 
1333  // Now all entities and reactions can be compiled
1334  ppEntity = mStateTemplate.beginIndependent();
1335  CModelEntity ** ppEntityEnd = mStateTemplate.endFixed();
1336 
1337  for (; ppEntity != ppEntityEnd; ++ppEntity)
1338  (*ppEntity)->compile();
1339 
1342 
1343  for (; itReaction != endReaction; ++itReaction)
1344  (*itReaction)->compile();
1345 
1346  return true;
1347 }
1348 
1350 {
1352  CModelEntity ** ppEntity = Entities.array();
1353 
1356 
1357  for (; itMetab != endMetab; ++itMetab)
1358  *ppEntity++ = *itMetab;;
1359 
1362 
1363  for (; itCompartment != endCompartment; ++itCompartment)
1364  *ppEntity++ = *itCompartment;
1365 
1368 
1369  for (; itValue != endValue; ++itValue)
1370  *ppEntity++ = *itValue;
1371 
1372  mStateTemplate.setUserOrder(Entities);
1373 
1375  //now sized to the number of entities with ODEs + all metabolites dependent on reactions
1376 
1377  const size_t * pUserOrder = mStateTemplate.getUserOrder().array();
1378  const size_t * pUserOrderEnd = pUserOrder + mStateTemplate.getUserOrder().size();
1379  ppEntity = mStateTemplate.getEntities();
1380 
1381  size_t i;
1382 
1383  for (i = 0; pUserOrder != pUserOrderEnd; ++pUserOrder)
1384  {
1385  const Status & Status = ppEntity[*pUserOrder]->getStatus();
1386 
1387  if (Status == ODE ||
1388  (Status == REACTIONS && ppEntity[*pUserOrder]->isUsed()))
1389  mJacobianPivot[i++] = *pUserOrder - 1;
1390  }
1391 
1392  return true;
1393 }
1394 
1396 {
1397  bool success = true;
1398 
1399  // The objects which are changed are all initial values of of all model entities including
1400  // fixed and unused once. Additionally, all kinetic parameters are possibly changed.
1401  // This is basically all the parameters in the parameter overview whose value is editable.
1402 
1403  // Issue 1170: We need to add elements of the stoichiometry, reduced stoichiometry,
1404  // and link matrix.
1405 
1407 
1408  Objects.insert(static_cast< const CCopasiObject * >(getObject(CCopasiObjectName("Reference=Avogadro Constant"))));
1409 
1410  // The initial values of the model entities
1411  CModelEntity **ppEntity = mStateTemplate.beginIndependent() - 1; // Offset for time
1412  CModelEntity **ppEntityEnd = mStateTemplate.endFixed();
1413 
1414  for (; ppEntity != ppEntityEnd; ++ppEntity)
1415  {
1416  // Assignments have no initial values
1417  if ((*ppEntity)->getStatus() != ASSIGNMENT ||
1418  (*ppEntity)->getInitialValueReference()->getDirectDependencies().size() == 0)
1419  Objects.insert((*ppEntity)->getInitialValueReference());
1420  }
1421 
1422  // The reaction parameters
1425  size_t i, imax;
1426 
1427  for (; itReaction != endReaction; ++itReaction)
1428  {
1429  const CCopasiParameterGroup & Group = (*itReaction)->getParameters();
1430 
1431  for (i = 0, imax = Group.size(); i < imax; i++)
1432  Objects.insert(Group.getParameter(i)->getValueReference());
1433  }
1434 
1435  // Fix for Issue 1170: We need to add elements of the stoichiometry, reduced stoichiometry,
1436  // and link matrices.
1437  if (mpStoiAnnotation != NULL)
1439 
1440  if (mpRedStoiAnnotation != NULL)
1442 
1443  if (mpLinkMatrixAnnotation != NULL)
1445 
1446  try
1447  {
1449  }
1450  catch (...)
1451  {
1452  mInitialRefreshes.clear();
1453  success = false;
1454  }
1455 
1456  mBuildInitialSequence = false;
1457 
1458  return success;
1459 }
1460 
1462 {
1463  if (mCompileIsNecessary)
1464  {
1465  compileIfNecessary(NULL);
1466  }
1467 
1469  {
1471  }
1472 
1473  // Update all initial values.
1474  std::vector< Refresh * >::const_iterator itRefresh = mInitialRefreshes.begin();
1475  std::vector< Refresh * >::const_iterator endRefresh = mInitialRefreshes.end();
1476 
1477  while (itRefresh != endRefresh)
1478  (**itRefresh++)();
1479 
1480  return true;
1481 }
1482 
1484 {
1485  bool success = true;
1486 
1487  // We need to add each used model entity to the objects which need to be updated.
1488  mSimulatedUpToDateObjects.clear();
1489 
1490  // For CModelValues and CCompartment ODEs we need to add the Rate
1491  // For CMetab ODEs we need to add the Particle Rate
1494 
1495  for (; ppEntity != ppEntityEnd; ++ppEntity) //loop over all entities with ODEs
1496  mSimulatedUpToDateObjects.insert((*ppEntity)->getRateReference());
1497 
1498  // We do not add the rates for metabolites of type REACTION. These are automatically calculated
1499  // with dgemm in calculate derivatives based on the reaction fluxes added below.
1500  // In the case that other simulated values depend on such a rate this is taken care by
1501  // calculating all dependencies.
1502  // This mechanism may lead occasionally to multiple calculations of rates of metabolites when used
1503  // in assignments or ODEs. However this is acceptable and more than compensated by the performance
1504  // gains of dgemm.
1505 
1506  // Furthermore all reaction fluxes have to be calculated too (see CMetab REACTION above)
1509 
1510  for (; itReaction != endReaction; ++itReaction)
1511  mSimulatedUpToDateObjects.insert((*itReaction)->getParticleFluxReference());
1512 
1513  // We now detect unused assignments, i.e., the result of an assignment is not
1514  // used during updateSimulatedValues except for itself or another unused assignment.
1515  bool UnusedFound = true;
1516 
1517  std::set<const CCopasiObject * > Candidate;
1518  std::set< const CCopasiObject * >::iterator it;
1519  std::set< const CCopasiObject * >::iterator end = mSimulatedUpToDateObjects.end();
1520  CCopasiObject * pObject;
1521  CMetab * pMetab;
1522  ppEntityEnd = mStateTemplate.endDependent();
1523 
1524  while (UnusedFound)
1525  {
1526  UnusedFound = false;
1528  //now points to the first entity with assignment
1529 
1530  for (; ppEntity != ppEntityEnd; ++ppEntity) //over all entities with assignments
1531  if ((*ppEntity)->isUsed())
1532  {
1533  if ((*ppEntity)->getStatus() != ASSIGNMENT)
1534  pObject = *ppEntity;
1535  else if ((pMetab = dynamic_cast< CMetab *>(*ppEntity)) != NULL)
1536  pObject = pMetab->getConcentrationReference();
1537  else
1538  pObject = (*ppEntity)->getValueReference();
1539 
1540  Candidate.insert(pObject);
1541 
1542  for (it = mSimulatedUpToDateObjects.begin(), end = mSimulatedUpToDateObjects.end(); it != end; ++it)
1543  if (*it != pObject &&
1544  (*it)->dependsOn(Candidate))
1545  break;
1546 
1547  if (it == end)
1548  {
1549  UnusedFound = true;
1550  mReorderNeeded = true;
1551  (*ppEntity)->setUsed(false);
1552  }
1553 
1554  Candidate.erase(pObject);
1555  }
1556  }
1557 
1558  if (mReorderNeeded)
1559  {
1561  CModelEntity ** ppReorder = Reorder.array();
1562 
1563  ppEntity = mStateTemplate.beginIndependent();
1565 
1566  //range is for all entities with ODEs + all metabs dependent on Reactions
1567  for (; ppEntity != ppEntityEnd; ++ppEntity, ++ppReorder)
1568  *ppReorder = *ppEntity;
1569 
1570  // :TODO: This must be enhanced as the mMetaboliteX and the state template may get out of sync
1571  // when we use assignments for metabolites.
1572 
1573  // the entities with assignments are reordered according to the isUsed() flag
1574  ppEntityEnd = mStateTemplate.endDependent();
1575 
1576  for (; ppEntity != ppEntityEnd; ++ppEntity) //over all entities with assignments
1577  if ((*ppEntity)->isUsed())
1578  *ppReorder++ = *ppEntity;
1579 
1581 
1582  for (; ppEntity != ppEntityEnd; ++ppEntity) //again over all entities with assignments
1583  if (!(*ppEntity)->isUsed())
1584  *ppReorder++ = *ppEntity;
1585 
1586  ppEntityEnd = mStateTemplate.endFixed();
1587 
1588  for (; ppEntity != ppEntityEnd; ++ppEntity, ++ppReorder) //over all fixed entities
1589  *ppReorder = *ppEntity;
1590 
1591  mStateTemplate.reorder(Reorder);
1592  mReorderNeeded = false;
1593 
1594  // We need to recompile as pointers to values may have changed
1595  ppEntity = mStateTemplate.beginIndependent();
1596  ppEntityEnd = mStateTemplate.endFixed();
1597 
1598  for (; ppEntity != ppEntityEnd; ++ppEntity)
1599  (*ppEntity)->compile();
1600 
1601  itReaction = mSteps.begin();
1602  endReaction = mSteps.end();
1603 
1604  for (; itReaction != endReaction; ++itReaction)
1605  (*itReaction)->compile();
1606 
1607  // The compile might have broken some refresh pointers we need to rebuild the constant sequence
1610  }
1611 
1612  std::set< const CCopasiObject * > UpToDate;
1613 
1614  try
1615  {
1617  }
1618  catch (...)
1619  {
1620  mSimulatedRefreshes.clear();
1621  success = false;
1622  }
1623 
1624  return success;
1625 }
1626 
1628 {
1629  bool success = true;
1630 
1632 
1633  std::set< const CCopasiObject * > Objects;
1634 
1635  const CMetab * pMetab;
1636 
1638  CModelEntity ** ppEntityEnd = mStateTemplate.endFixed();
1639 
1640  for (; ppEntity != ppEntityEnd; ++ppEntity)
1641  {
1642  if ((*ppEntity)->getStatus() == ODE)
1643  {
1644  // TODO We need only to calculate rates which are constant since the other will
1645  // be updated by the simulation request.
1646  Objects.insert((*ppEntity)->getRateReference());
1647  }
1648 
1649  // Species of type assignment have a second pseudo state value the concentration,
1650  // which always can be directly calculated.
1651  if ((*ppEntity)->getStatus() == ASSIGNMENT &&
1652  (pMetab = dynamic_cast< const CMetab * >(*ppEntity)) != NULL)
1653  {
1654  mApplyInitialValuesRefreshes.push_back(pMetab->getConcentrationReference()->getApplyInitialValueRefresh());
1655  }
1656  }
1657 
1658  std::set< const CCopasiObject * > UpToDate;
1659 
1660  try
1661  {
1662  std::vector< Refresh * > RateRefreshes =
1663  CCopasiObject::buildUpdateSequence(Objects, UpToDate);
1664 
1666  RateRefreshes.begin(),
1667  RateRefreshes.end());
1668  }
1669  catch (...)
1670  {
1672  success = false;
1673  }
1674 
1675  return success;
1676 }
1677 
1679 {
1680  bool success = true;
1681 
1683 
1684  // Compartments
1687 
1688  for (; itComp != endComp; ++itComp)
1689  {
1690  Objects.insert((*itComp)->getValueReference());
1691 
1692  switch ((*itComp)->getStatus())
1693  {
1694  case ODE:
1695  Objects.insert((*itComp)->getRateReference());
1696  break;
1697 
1698  default:
1699  break;
1700  }
1701  }
1702 
1703  // Metabolites
1706 
1707  for (; itMetab != endMetab; ++itMetab)
1708  {
1709  Objects.insert((*itMetab)->getConcentrationReference());
1710  Objects.insert((*itMetab)->getValueReference());
1711 
1712  switch ((*itMetab)->getStatus())
1713  {
1714  case REACTIONS:
1715  case ODE:
1716  Objects.insert(static_cast< const CCopasiObject *>((*itMetab)->getObject(CCopasiObjectName("Reference=TransitionTime"))));
1717  Objects.insert((*itMetab)->getConcentrationRateReference());
1718  Objects.insert((*itMetab)->getRateReference());
1719  break;
1720 
1721  default:
1722  break;
1723  }
1724  }
1725 
1726  // Reactions
1729 
1730  for (; itStep != endStep; ++itStep)
1731  {
1732  Objects.insert(static_cast< const CCopasiObject * >((*itStep)->getObject(CCopasiObjectName("Reference=Flux"))));
1733  Objects.insert((*itStep)->getParticleFluxReference());
1734  }
1735 
1736  // Model Values
1739 
1740  for (; itValue != endValue; ++itValue)
1741  {
1742  Objects.insert((*itValue)->getValueReference());
1743 
1744  switch ((*itValue)->getStatus())
1745  {
1746  case ODE:
1747  Objects.insert((*itValue)->getRateReference());
1748  break;
1749 
1750  default:
1751  break;
1752  }
1753  }
1754 
1755  try
1756  {
1758  }
1759  catch (...)
1760  {
1761  mNonSimulatedRefreshes.clear();
1762  success = false;
1763  }
1764 
1765  return success;
1766 }
1767 
1769 {return mInitialState;}
1770 
1771 const CState & CModel::getState() const
1772 {return mCurrentState;}
1773 
1774 void CModel::setInitialState(const CState & state)
1775 {
1776  mInitialState = state;
1777 
1778  if (mIsAutonomous &&
1780  mInitialState.setTime(0.0);
1781 
1782  return;
1783 }
1784 
1785 void CModel::setState(const CState & state)
1786 {
1787  mCurrentState = state;
1788 
1789  return;
1790 }
1791 
1793  const CCopasiObject::DataObjectSet & changedObjects,
1794  const CCopasiObject::DataObjectSet & requestedObjects,
1795  CCopasiObject::DataUpdateSequence & updateSequence) const
1796 {
1798  context,
1799  changedObjects,
1800  requestedObjects,
1801  updateSequence);
1802 }
1803 
1805  const CCopasiObject::DataObjectSet & changedObjects,
1806  const CCopasiObject::DataObjectSet & requestedObjects,
1807  CCopasiObject::DataUpdateSequence & updateSequence) const
1808 {
1810  context,
1811  changedObjects,
1812  requestedObjects,
1813  updateSequence);
1814 }
1815 
1817  const CMath::SimulationContextFlag & context,
1818  const CCopasiObject::DataObjectSet & changedObjects,
1819  const CCopasiObject::DataObjectSet & requestedObjects,
1820  CCopasiObject::DataUpdateSequence & updateSequence) const
1821 {
1822  updateSequence.clear();
1823 
1825 
1826  if (!dependencyGraph.getUpdateSequence(context,
1827  *reinterpret_cast< const CObjectInterface::ObjectSet *>(&changedObjects),
1828  *reinterpret_cast< const CObjectInterface::ObjectSet *>(&requestedObjects),
1829  UpdateSequence))
1830  {
1831  return false;
1832  }
1833 
1834  CObjectInterface::UpdateSequence::iterator it = UpdateSequence.begin();
1835  CObjectInterface::UpdateSequence::iterator end = UpdateSequence.end();
1836  Refresh * pRefresh = NULL;
1837 
1838  for (; it != end; ++it)
1839  {
1840  pRefresh = static_cast< CCopasiObject * >(*it)->getRefresh();
1841 
1842  if (pRefresh != NULL)
1843  {
1844  updateSequence.push_back(pRefresh);
1845  }
1846  }
1847 
1848  return true;
1849 }
1850 
1851 void CModel::updateSimulatedValues(const bool & updateMoieties)
1852 {
1853  // Depending on which model we are using we need to update
1854  // the particle numbers for the dependent metabolites.
1855  if (updateMoieties)
1856  {
1857  C_FLOAT64 * pDependent = mCurrentState.beginDependent();
1860 
1861  for (; itMoiety != endMoiety; ++itMoiety, ++pDependent)
1862  *pDependent = (*itMoiety)->dependentNumber();
1863  }
1864 
1865  std::vector< Refresh * >::const_iterator itRefresh = mSimulatedRefreshes.begin();
1866  std::vector< Refresh * >::const_iterator endRefresh = mSimulatedRefreshes.end();
1867 
1868  while (itRefresh != endRefresh)
1869  (**itRefresh++)();
1870 
1871  // Store the particle fluxes for further calculations
1874 
1875  C_FLOAT64 * pFlux = mParticleFluxes.array();
1876 
1877  for (; it != end; ++it, ++pFlux)
1878  {
1879  *pFlux = (*it)->getParticleFlux();
1880 
1881 #ifdef COPASI_DEBUG
1882 
1883  if (isnan(*pFlux))
1884  {
1885  std::cout << "Reaction Flux[" << it - mSteps.begin() << "] is NaN" << std::endl;
1886  }
1887 
1888 #endif // COPASI_DEBUG
1889  }
1890 }
1891 
1893 {
1894  std::vector< Refresh * >::const_iterator itRefresh = mNonSimulatedRefreshes.begin();
1895  std::vector< Refresh * >::const_iterator endRefresh = mNonSimulatedRefreshes.end();
1896 
1897  while (itRefresh != endRefresh)
1898  (**itRefresh++)();
1899 
1900  return;
1901 }
1902 
1904 {
1905  C_FLOAT64 * pTmp = derivatives;
1906 
1907  // First retrieve derivatives of quantities determined by ODE
1908  // The offset 1 is for the model time which is always the first
1909  // state variable.
1910  CModelEntity ** ppIt = mStateTemplate.getEntities() + 1;
1911  CModelEntity ** ppEnd =
1913 
1914  for (; ppIt != ppEnd; ++ppIt, ++pTmp)
1915  {
1916  *pTmp = (*ppIt)->getRate();
1917 
1918 #ifdef COPASI_DEBUG
1919 
1920  if (isnan(*pTmp))
1921  {
1922  std::cout << "Entity[" << ppIt - (mStateTemplate.getEntities() + 1) << "] is NaN" << std::endl;
1923  }
1924 
1925 #endif // COPASI_DEBUG
1926  }
1927 
1928  // Now calculate derivatives of all metabolites determined by reactions
1929  char T = 'N';
1930  C_INT M = 1;
1932  C_INT K = (C_INT) mSteps.size();
1933  C_FLOAT64 Alpha = 1.0;
1934  C_FLOAT64 Beta = 0.0;
1935 
1936  if (K != 0)
1937  dgemm_(&T, &T, &M, &N, &K, &Alpha, mParticleFluxes.array(), &M,
1938  mStoi.array(), &K, &Beta, pTmp, &M);
1939 }
1940 
1942 {
1943  C_FLOAT64 * pTmp = derivativesX;
1944 
1945  // First retrieve derivatives of quantities determined by ODE
1946  // The offset 1 is for the model time which is always the first
1947  // state variable.
1948  CModelEntity ** ppIt = mStateTemplate.getEntities() + 1;
1949  CModelEntity ** ppEnd =
1951 
1952  for (; ppIt != ppEnd; ++ppIt, ++pTmp)
1953  *pTmp = (*ppIt)->getRate();
1954 
1955  // Now calculate derivatives of the independent metabolites determined by reactions
1956  char T = 'N';
1957  C_INT M = 1;
1959  C_INT K = (C_INT) mSteps.size();
1960  C_FLOAT64 Alpha = 1.0;
1961  C_FLOAT64 Beta = 0.0;
1962 
1963  if (K != 0)
1964  dgemm_(&T, &T, &M, &N, &K, &Alpha, mParticleFluxes.array(), &M,
1965  mRedStoi.array(), &K, &Beta, pTmp, &M);
1966 }
1967 
1969  const C_FLOAT64 & resolution)
1970 {
1971  size_t Col;
1972  size_t nCol = mElasticities.numCols();
1973 
1974  C_FLOAT64 * itE;
1975  C_FLOAT64 * beginE = mElasticities.array();
1976 
1980 
1981  CModelEntity ** itEntity;
1982  CModelEntity ** beginEntity = mStateTemplate.beginIndependent();
1983  CModelEntity ** endEntity = mStateTemplate.endDependent();
1984 
1985  for (itEntity = beginEntity, Col = 0; itEntity != endEntity; ++itEntity, ++Col)
1986  {
1987  // :TODO: This only works for entities of type metabolites.
1988  // The scaling factor for other entities should be 1.
1989  const C_FLOAT64 invVolume =
1990  1.0 / static_cast<CMetab *>(*itEntity)->getCompartment()->getValue();
1991  C_FLOAT64 * pX =
1992  const_cast<C_FLOAT64 *>(&static_cast<CMetab *>(*itEntity)->getConcentration());
1993 
1994  for (itReaction = beginReaction, itE = beginE + Col;
1995  itReaction != endReaction;
1996  ++itReaction, itE += nCol)
1997  * itE = invVolume * (*itReaction)->calculatePartialDerivative(pX, factor, resolution);
1998  }
1999 
2000  // DebugFile << "CModel::calculateElasticityMatrix()" << std::endl;
2001  // DebugFile << mElasticities << std::endl;
2002 
2003  return;
2004 }
2005 
2007  const C_FLOAT64 & derivationFactor,
2008  const C_FLOAT64 & /* resolution */)
2009 {
2010  C_FLOAT64 DerivationFactor = std::max(derivationFactor, 100.0 * std::numeric_limits< C_FLOAT64 >::epsilon());
2011 
2012  size_t Dim =
2014  //Dim now contains the number of entities with ODEs + number of metabs depending on reactions.
2015 
2016  size_t Col;
2017 
2018  jacobian.resize(Dim, Dim);
2019 
2020  C_FLOAT64 Store;
2021  C_FLOAT64 X1;
2022  C_FLOAT64 X2;
2023  C_FLOAT64 InvDelta;
2024 
2025  CVector< C_FLOAT64 > Y1(Dim);
2026  CVector< C_FLOAT64 > Y2(Dim);
2027 
2028  C_FLOAT64 * pY1;
2029  C_FLOAT64 * pY2;
2030 
2032  C_FLOAT64 * pXEnd = pX + Dim;
2033 
2034  C_FLOAT64 * pJacobian;
2035  C_FLOAT64 * pJacobianEnd = jacobian.array() + Dim * Dim;
2036 
2037  for (Col = 0; pX != pXEnd; ++pX, ++Col)
2038  {
2039  Store = *pX;
2040 
2041  // We only need to make sure that we do not have an underflow problem
2042  if (fabs(Store) < DerivationFactor)
2043  {
2044  X1 = 0.0;
2045 
2046  if (Store < 0.0)
2047  X2 = -2.0 * DerivationFactor;
2048  else
2049  X2 = 2.0 * DerivationFactor;
2050  }
2051  else
2052  {
2053  X1 = Store * (1.0 + DerivationFactor);
2054  X2 = Store * (1.0 - DerivationFactor);
2055  }
2056 
2057  InvDelta = 1.0 / (X2 - X1);
2058 
2059  *pX = X1;
2060  updateSimulatedValues(false);
2062 
2063  *pX = X2;
2064  updateSimulatedValues(false);
2066 
2067  *pX = Store;
2068 
2069  pJacobian = jacobian.array() + Col;
2070  pY1 = Y1.array();
2071  pY2 = Y2.array();
2072 
2073  for (; pJacobian < pJacobianEnd; pJacobian += Dim, ++pY1, ++pY2)
2074  * pJacobian = (*pY2 - *pY1) * InvDelta;
2075  }
2076 
2077  updateSimulatedValues(false);
2078 
2079  return;
2080 }
2081 
2083  const C_FLOAT64 & derivationFactor,
2084  const C_FLOAT64 & /* resolution */)
2085 {
2086  C_FLOAT64 DerivationFactor = std::max(derivationFactor, 100.0 * std::numeric_limits< C_FLOAT64 >::epsilon());
2087 
2088  size_t Dim = mCurrentState.getNumIndependent();
2089  size_t Col;
2090 
2091  jacobianX.resize(Dim, Dim);
2092 
2093  C_FLOAT64 Store;
2094  C_FLOAT64 X1;
2095  C_FLOAT64 X2;
2096  C_FLOAT64 InvDelta;
2097 
2098  CVector< C_FLOAT64 > Y1(Dim);
2099  CVector< C_FLOAT64 > Y2(Dim);
2100 
2101  C_FLOAT64 * pY1;
2102  C_FLOAT64 * pY2;
2103 
2105  C_FLOAT64 * pXEnd = pX + Dim;
2106 
2107  C_FLOAT64 * pJacobian;
2108  C_FLOAT64 * pJacobianEnd = jacobianX.array() + Dim * Dim;
2109 
2110  for (Col = 0; pX != pXEnd; ++pX, ++Col)
2111  {
2112  Store = *pX;
2113 
2114  // We only need to make sure that we do not have an underflow problem
2115  if (fabs(Store) < DerivationFactor)
2116  {
2117  X1 = 0.0;
2118 
2119  if (Store < 0.0)
2120  X2 = -2.0 * DerivationFactor;
2121  else
2122  X2 = 2.0 * DerivationFactor;;
2123  }
2124  else
2125  {
2126  X1 = Store * (1.0 + DerivationFactor);
2127  X2 = Store * (1.0 - DerivationFactor);
2128  }
2129 
2130  InvDelta = 1.0 / (X2 - X1);
2131 
2132  *pX = X1;
2133  updateSimulatedValues(true);
2135 
2136  *pX = X2;
2137  updateSimulatedValues(true);
2139 
2140  *pX = Store;
2141 
2142  pJacobian = jacobianX.array() + Col;
2143  pY1 = Y1.array();
2144  pY2 = Y2.array();
2145 
2146  for (; pJacobian < pJacobianEnd; pJacobian += Dim, ++pY1, ++pY2)
2147  * pJacobian = (*pY2 - *pY1) * InvDelta;
2148  }
2149 
2150  updateSimulatedValues(true);
2151 }
2152 
2154 {
2155  fatalError(); //not yet implemented
2156  return 0.0;
2157 }
2158 
2159 bool CModel::setVolumeUnit(const std::string & name)
2160 {
2161  return setVolumeUnit(toEnum(name.c_str(), VolumeUnitNames, ml));
2162 }
2163 
2165 {
2166  mVolumeUnit = unit;
2167  return true;
2168 }
2169 
2170 std::string CModel::getVolumeUnitName() const
2171 {
2172  return VolumeUnitNames[mVolumeUnit];
2173 }
2174 
2176 {
2177  return mVolumeUnit;
2178 }
2179 
2180 //****
2181 
2182 bool CModel::setAreaUnit(const std::string & name)
2183 {
2184  return setAreaUnit(toEnum(name.c_str(), AreaUnitNames, m2));
2185 }
2186 
2188 {
2189  mAreaUnit = unit;
2190  return true;
2191 }
2192 
2193 std::string CModel::getAreaUnitName() const
2194 {
2195  return AreaUnitNames[mAreaUnit];
2196 }
2197 
2199 {
2200  return mAreaUnit;
2201 }
2202 
2203 //****
2204 bool CModel::setLengthUnit(const std::string & name)
2205 {
2206  return setLengthUnit(toEnum(name.c_str(), LengthUnitNames, m));
2207 }
2208 
2210 {
2211  mLengthUnit = unit;
2212  return true;
2213 }
2214 
2215 std::string CModel::getLengthUnitName() const
2216 {
2217  return LengthUnitNames[mLengthUnit];
2218 }
2219 
2221 {
2222  return mLengthUnit;
2223 }
2224 
2225 //****
2226 
2227 bool CModel::setTimeUnit(const std::string & name)
2228 {
2229  return setTimeUnit(toEnum(name.c_str(), TimeUnitNames, s));
2230 }
2231 
2233 {
2234  mTimeUnit = unit;
2235  return true;
2236 }
2237 
2238 std::string CModel::getTimeUnitName() const
2239 {
2240  return TimeUnitNames[mTimeUnit];
2241 }
2242 
2244 {
2245  return mTimeUnit;
2246 }
2247 
2248 //****
2249 
2250 bool CModel::setQuantityUnit(const std::string & name)
2251 {
2252  QuantityUnit unit = toEnum(name.c_str(), QuantityUnitNames, OldXML);
2253 
2254  if (unit == OldXML)
2255  unit = toEnum(name.c_str(), QuantityUnitOldXMLNames, mMol);
2256 
2257  return setQuantityUnit(unit);
2258 }
2259 
2261 {
2262  bool success = true;
2263  mQuantityUnit = unit;
2264 
2265  switch (unit)
2266  {
2267  case Mol:
2269  break;
2270 
2271  case mMol:
2273  break;
2274 
2275  case microMol:
2277  break;
2278 
2279  case nMol:
2281  break;
2282 
2283  case pMol:
2285  break;
2286 
2287  case fMol:
2289  break;
2290 
2291  case number:
2292  mQuantity2NumberFactor = 1.0;
2293  break;
2294 
2295  case dimensionlessQuantity:
2296  mQuantity2NumberFactor = 1.0;
2297  break;
2298 
2299  default:
2301  mQuantity2NumberFactor = 1.0;
2302  success = false;
2303  break;
2304  }
2305 
2307 
2308  //adapt particle numbers
2309  size_t i, imax = mMetabolites.size();
2310 
2311  for (i = 0; i < imax; ++i)
2312  {
2313  //update particle numbers
2314  mMetabolites[i]->setInitialConcentration(mMetabolites[i]->getInitialConcentration());
2315  mMetabolites[i]->setConcentration(mMetabolites[i]->getConcentration());
2316  }
2317 
2318  return success;
2319 }
2320 
2321 std::string CModel::getQuantityUnitName() const
2322 {
2324 }
2325 
2327 {
2329 }
2330 
2332 {
2333  return mQuantityUnit;
2334 }
2335 
2337 {mType = modelType;}
2338 
2340 {return mType;}
2341 
2342 void CModel::setAvogadro(const C_FLOAT64 & avogadro)
2343 {
2344  mAvogadro = avogadro;
2345 
2347 }
2348 
2350 {
2351  return mAvogadro;
2352 }
2353 
2355 {return mQuantity2NumberFactor;}
2356 
2358 {return mNumber2QuantityFactor;}
2359 
2360 //*****
2361 
2362 //**********************************************************************
2363 
2364 bool CModel::appendDependentModelObjects(const std::set< const CCopasiObject * > & deletedObjects,
2365  std::set< const CCopasiObject * > & dependentReactions,
2366  std::set< const CCopasiObject * > & dependentMetabolites,
2367  std::set< const CCopasiObject * > & dependentCompartments,
2368  std::set< const CCopasiObject * > & dependentModelValues,
2369  std::set< const CCopasiObject * > & dependentEvents) const
2370 {
2371  // We need a local copy since we recursively add deleted objects.
2372  std::set< const CCopasiObject * > DeletedObjects = deletedObjects;
2373 
2374  bool ObjectsAppended = false;
2375  bool DeleteObjects = DeletedObjects.size() > 0;
2376 
2377  // This is this implemented recursively. Since deleting a container may result
2378  // in the deletion of objects not dependent on the original set of deleted objects.
2379 
2380  while (DeleteObjects)
2381  {
2382  DeleteObjects = false;
2383 
2384  DeleteObjects |= appendDependentReactions(DeletedObjects, dependentReactions);
2385 
2386  if (dependentReactions.size() > 0)
2387  {
2388  std::set< const CCopasiObject * >::const_iterator it, itEnd = dependentReactions.end();
2389 
2390  for (it = dependentReactions.begin(); it != itEnd; ++it)
2391  if (DeletedObjects.find(*it) == DeletedObjects.end())
2392  {
2393  CCopasiObject::DataObjectSet AdditionalObjects =
2394  static_cast< const CReaction * >(*it)->getDeletedObjects();
2395 
2396  std::set< const CCopasiObject * >::const_iterator itDeleted = AdditionalObjects.begin();
2397  std::set< const CCopasiObject * >::const_iterator endDeleted = AdditionalObjects.end();
2398 
2399  for (; itDeleted != endDeleted; ++itDeleted)
2400  DeletedObjects.insert(*itDeleted);
2401  }
2402  }
2403 
2404  DeleteObjects |= appendDependentMetabolites(DeletedObjects, dependentMetabolites);
2405 
2406  if (dependentMetabolites.size() > 0)
2407  {
2408  std::set< const CCopasiObject * >::const_iterator it, itEnd = dependentMetabolites.end();
2409 
2410  for (it = dependentMetabolites.begin(); it != itEnd; ++it)
2411  if (DeletedObjects.find(*it) == DeletedObjects.end())
2412  {
2413  CCopasiObject::DataObjectSet AdditionalObjects =
2414  static_cast< const CMetab * >(*it)->getDeletedObjects();
2415 
2416  std::set< const CCopasiObject * >::const_iterator itDeleted = AdditionalObjects.begin();
2417  std::set< const CCopasiObject * >::const_iterator endDeleted = AdditionalObjects.end();
2418 
2419  for (; itDeleted != endDeleted; ++itDeleted)
2420  DeletedObjects.insert(*itDeleted);
2421  }
2422  }
2423 
2424  DeleteObjects |= appendDependentModelValues(DeletedObjects, dependentModelValues);
2425 
2426  if (dependentModelValues.size() > 0)
2427  {
2428  std::set< const CCopasiObject * >::const_iterator it, itEnd = dependentModelValues.end();
2429 
2430  for (it = dependentModelValues.begin(); it != itEnd; ++it)
2431  if (DeletedObjects.find(*it) == DeletedObjects.end())
2432  {
2433  CCopasiObject::DataObjectSet AdditionalObjects =
2434  static_cast< const CModelValue * >(*it)->getDeletedObjects();
2435 
2436  std::set< const CCopasiObject * >::const_iterator itDeleted = AdditionalObjects.begin();
2437  std::set< const CCopasiObject * >::const_iterator endDeleted = AdditionalObjects.end();
2438 
2439  for (; itDeleted != endDeleted; ++itDeleted)
2440  DeletedObjects.insert(*itDeleted);
2441  }
2442  }
2443 
2444  DeleteObjects |= appendDependentCompartments(DeletedObjects, dependentCompartments);
2445 
2446  if (dependentCompartments.size() > 0)
2447  {
2448  std::set< const CCopasiObject * >::const_iterator it, itEnd = dependentCompartments.end();
2449 
2450  for (it = dependentCompartments.begin(); it != itEnd; ++it)
2451  if (DeletedObjects.find(*it) == DeletedObjects.end())
2452  {
2453  CCopasiObject::DataObjectSet AdditionalObjects =
2454  static_cast< const CCompartment * >(*it)->getDeletedObjects();
2455 
2456  std::set< const CCopasiObject * >::const_iterator itDeleted = AdditionalObjects.begin();
2457  std::set< const CCopasiObject * >::const_iterator endDeleted = AdditionalObjects.end();
2458 
2459  for (; itDeleted != endDeleted; ++itDeleted)
2460  DeletedObjects.insert(*itDeleted);
2461  }
2462  }
2463 
2464  DeleteObjects |= appendDependentEvents(DeletedObjects, dependentEvents);
2465 
2466  ObjectsAppended |= DeleteObjects;
2467  }
2468 
2469  return ObjectsAppended;
2470 }
2471 
2472 bool CModel::appendDependentReactions(std::set< const CCopasiObject * > candidates,
2473  std::set< const CCopasiObject * > & dependents) const
2474 {
2475  const_cast< CModel * >(this)->compileIfNecessary(NULL);
2476 
2477  size_t Size = dependents.size();
2478 
2481 
2482  CCopasiObject::DataObjectSet::const_iterator itSet;
2483  CCopasiObject::DataObjectSet::const_iterator endSet;
2484 
2485  for (; it != end; ++it)
2486  if (candidates.find(*it) == candidates.end())
2487  {
2488  // Check whether the reaction is already in the list of deleted objects
2489  if (candidates.find(*it) == candidates.end())
2490  {
2491  if ((*it)->mustBeDeleted(candidates))
2492  {
2493  dependents.insert(*it);
2494  }
2495  }
2496  }
2497 
2498  return Size < dependents.size();
2499 }
2500 
2501 bool CModel::appendDependentMetabolites(std::set< const CCopasiObject * > candidates,
2502  std::set< const CCopasiObject * > & dependents) const
2503 {
2504  const_cast< CModel * >(this)->compileIfNecessary(NULL);
2505 
2506  size_t Size = dependents.size();
2507 
2510 
2513 
2514  CCopasiObject::DataObjectSet::const_iterator itSet;
2515  CCopasiObject::DataObjectSet::const_iterator endSet;
2516 
2517  for (; itComp != endComp; ++itComp)
2518  {
2519  it = (*itComp)->getMetabolites().begin();
2520  end = (*itComp)->getMetabolites().end();
2521 
2522  for (; it != end; ++it)
2523  {
2524  // Check whether the species is already in the list of deleted objects
2525  if (candidates.find(*it) == candidates.end())
2526  {
2527  if ((*it)->mustBeDeleted(candidates))
2528  {
2529  dependents.insert(*it);
2530  }
2531  }
2532  }
2533  }
2534 
2535  return Size < dependents.size();
2536 }
2537 
2538 bool CModel::appendDependentCompartments(std::set< const CCopasiObject * > candidates,
2539  std::set< const CCopasiObject * > & dependents) const
2540 {
2541  const_cast< CModel * >(this)->compileIfNecessary(NULL);
2542 
2543  size_t Size = dependents.size();
2544 
2547 
2548  std::set< const CCopasiObject * >::const_iterator itSet;
2549  std::set< const CCopasiObject * >::const_iterator endSet;
2550 
2551  for (; it != end; ++it)
2552  {
2553  // Check whether the compartment is already in the list of deleted objects
2554  if (candidates.find(*it) == candidates.end())
2555  {
2556  if ((*it)->mustBeDeleted(candidates))
2557  {
2558  dependents.insert(*it);
2559  }
2560  }
2561  }
2562 
2563  return Size < dependents.size();
2564 }
2565 
2566 bool CModel::appendDependentModelValues(std::set< const CCopasiObject * > candidates,
2567  std::set< const CCopasiObject * > & dependents) const
2568 {
2569  const_cast< CModel * >(this)->compileIfNecessary(NULL);
2570 
2571  size_t Size = dependents.size();
2572 
2575 
2576  std::set< const CCopasiObject * >::const_iterator itSet;
2577  std::set< const CCopasiObject * >::const_iterator endSet;
2578 
2579  for (; it != end; ++it)
2580 
2581  // Check whether the model value is already in the list of deleted objects
2582  if (candidates.find(*it) == candidates.end())
2583  {
2584  if ((*it)->mustBeDeleted(candidates))
2585  {
2586  dependents.insert(*it);
2587  }
2588  }
2589 
2590  return Size < dependents.size();
2591 }
2592 
2593 bool CModel::appendDependentEvents(std::set< const CCopasiObject * > candidates,
2594  std::set< const CCopasiObject * > & dependents) const
2595 {
2596  const_cast< CModel * >(this)->compileIfNecessary(NULL);
2597 
2598  size_t Size = dependents.size();
2599 
2602 
2603  std::set< const CCopasiObject * >::const_iterator itSet;
2604  std::set< const CCopasiObject * >::const_iterator endSet;
2605 
2606  for (; it != end; ++it)
2607 
2608  // Check whether the model value is already in the list of deleted objects
2609  if (candidates.find(*it) == candidates.end())
2610  {
2611  if ((*it)->mustBeDeleted(candidates))
2612  {
2613  dependents.insert(*it);
2614  }
2615  }
2616 
2617  return Size < dependents.size();
2618 }
2619 
2620 //**********************************************************************
2621 
2622 CMetab* CModel::createMetabolite(const std::string & name,
2623  const std::string & compartment,
2624  const C_FLOAT64 & iconc,
2625  const CMetab::Status & status)
2626 {
2627  size_t Index;
2628 
2629  if (mCompartments.size() == 0)
2630  return NULL;
2631 
2632  if (compartment == "")
2633  Index = 0;
2634  else if ((Index = mCompartments.getIndex(compartment)) == C_INVALID_INDEX)
2635  return NULL;
2636 
2637  if (mCompartments[Index]->getMetabolites().getIndex(name) != C_INVALID_INDEX)
2638  return NULL;
2639 
2640  CMetab * pMetab = new CMetab(name);
2641 
2642  if (!mCompartments[Index]->addMetabolite(pMetab))
2643  {
2644  delete pMetab;
2645  return NULL;
2646  }
2647 
2648  pMetab->setStatus(status);
2649  pMetab->setInitialConcentration(iconc);
2650  pMetab->refreshInitialValue();
2651 
2652  if (!mMetabolites.add(pMetab))
2653  return NULL;
2654 
2655  mCompileIsNecessary = true;
2656 
2657  return pMetab;
2658 }
2659 
2660 bool CModel::removeMetabolite(const size_t index,
2661  const bool & recursive)
2662 {
2663  const CMetab* pMetabolite = getMetabolites()[index];
2664  return removeMetabolite(pMetabolite, recursive);
2665 }
2666 
2667 bool CModel::removeMetabolite(const std::string & key,
2668  const bool & recursive)
2669 {
2670  CMetab* pMetabolite =
2671  dynamic_cast< CMetab * >(CCopasiRootContainer::getKeyFactory()->get(key));
2672  return removeMetabolite(pMetabolite, recursive);
2673 }
2674 
2675 bool CModel::removeMetabolite(const CMetab* pMetabolite,
2676  const bool & recursive)
2677 {
2678  if (!pMetabolite)
2679  return false;
2680 
2681  if (recursive)
2682  {
2684  }
2685 
2686  /* Assure that all references are removed */
2687  mMetabolites.remove((CMetab *)pMetabolite);
2688  mMetabolitesX.remove((CMetab *)pMetabolite);
2689 
2690  pdelete(pMetabolite);
2691 
2692  clearMoieties();
2693  mCompileIsNecessary = true;
2694 
2695  return true;
2696 }
2697 
2698 CCompartment* CModel::createCompartment(const std::string & name,
2699  const C_FLOAT64 & volume)
2700 {
2701  // check if there is already a volume with this name
2702  if (mCompartments.getIndex(name) != C_INVALID_INDEX)
2703  return NULL;
2704 
2705  CCompartment * cpt = new CCompartment(name);
2706 
2707  cpt->setInitialValue(volume);
2708  //cpt->setVolume(volume);
2709 
2710  if (!mCompartments.add(cpt, true))
2711  {
2712  delete cpt;
2713  return NULL;
2714  }
2715 
2716  mCompileIsNecessary = true;
2717  return cpt;
2718 }
2719 
2720 bool CModel::removeCompartment(const size_t index,
2721  const bool & recursive)
2722 {
2723  const CCompartment * pCompartment = getCompartments()[index];
2724  return removeCompartment(pCompartment, recursive);
2725 }
2726 
2727 bool CModel::removeCompartment(const std::string & key,
2728  const bool & recursive)
2729 {
2730  CCompartment *pCompartment =
2731  dynamic_cast< CCompartment * >(CCopasiRootContainer::getKeyFactory()->get(key));
2732  return removeCompartment(pCompartment, recursive);
2733 }
2734 
2735 bool CModel::removeCompartment(const CCompartment * pCompartment,
2736  const bool & recursive)
2737 {
2738  if (!pCompartment)
2739  return false;
2740 
2741  if (recursive)
2742  {
2744  }
2745 
2746  //Check if Compartment with that name exists
2747  size_t index =
2748  mCompartments.CCopasiVector< CCompartment >::getIndex(pCompartment);
2749 
2750  if (index == C_INVALID_INDEX)
2751  return false;
2752 
2753  mCompartments.CCopasiVector< CCompartment >::remove(index);
2754 
2755  mCompileIsNecessary = true;
2756 
2757  return true;
2758 }
2759 
2760 CReaction* CModel::createReaction(const std::string & name)
2761 {
2762  if (mSteps.getIndex(name) != C_INVALID_INDEX)
2763  return NULL;
2764 
2765  CReaction * pReaction = new CReaction(name);
2766 
2767  if (!mSteps.add(pReaction, true))
2768  {
2769  delete pReaction;
2770  return NULL;
2771  }
2772 
2773  mCompileIsNecessary = true;
2774  return pReaction;
2775 }
2776 
2777 bool CModel::removeReaction(const std::string & key,
2778  const bool & recursive)
2779 {
2780  CReaction * pReaction =
2781  dynamic_cast< CReaction * >(CCopasiRootContainer::getKeyFactory()->get(key));
2782  return removeReaction(pReaction, recursive);
2783 }
2784 
2785 bool CModel::removeReaction(const size_t index,
2786  const bool & recursive)
2787 {
2788  const CReaction * pReaction = getReactions()[index];
2789  return removeReaction(pReaction, recursive);
2790 }
2791 
2792 bool CModel::removeReaction(const CReaction * pReaction,
2793  const bool & recursive)
2794 {
2795  if (!pReaction)
2796  return false;
2797 
2798  if (recursive)
2799  {
2801  }
2802 
2803  //Check if Reaction exists
2804  size_t index =
2805  mSteps.CCopasiVector< CReaction >::getIndex(pReaction);
2806 
2807  if (index == C_INVALID_INDEX)
2808  return false;
2809 
2810  mSteps.CCopasiVector< CReaction >::remove(index);
2811 
2812  clearMoieties();
2813  mCompileIsNecessary = true;
2814 
2815  return true;
2816 }
2817 
2818 bool CModel::removeLocalReactionParameter(const std::string & key,
2819  const bool & recursive)
2820 {
2821  CCopasiParameter * pParameter =
2822  dynamic_cast< CCopasiParameter * >(CCopasiRootContainer::getKeyFactory()->get(key));
2823 
2824  if (pParameter == NULL)
2825  return false;
2826 
2827  if (recursive)
2828  {
2829  std::set< const CCopasiObject * > DeletedObjects;
2830  DeletedObjects.insert(pParameter->getValueReference());
2831 
2832  removeDependentModelObjects(DeletedObjects);
2833  }
2834 
2835  return true;
2836 }
2837 
2838 CModelValue* CModel::createModelValue(const std::string & name,
2839  const C_FLOAT64 & value)
2840 {
2841  // check if there is already a value with this name
2842  if (mValues.getIndex(name) != C_INVALID_INDEX)
2843  return NULL;
2844 
2845  CModelValue * cmv = new CModelValue(name);
2846 
2847  cmv->setInitialValue(value);
2848 
2849  if (!mValues.add(cmv, true))
2850  {
2851  delete cmv;
2852  return NULL;
2853  }
2854 
2855  mCompileIsNecessary = true;
2856  return cmv;
2857 }
2858 
2859 void CModel::removeDependentModelObjects(const std::set<const CCopasiObject*> & deletedObjects)
2860 {
2861  std::set<const CCopasiObject*> Reactions;
2862  std::set<const CCopasiObject*> Metabolites;
2863  std::set<const CCopasiObject*> Values;
2864  std::set<const CCopasiObject*> Compartments;
2865  std::set<const CCopasiObject*> Events;
2866 
2867  appendDependentModelObjects(deletedObjects, Reactions, Metabolites, Compartments, Values, Events);
2868 
2869  std::set<const CCopasiObject*>::const_iterator it, end;
2870 
2871  for (it = Reactions.begin(), end = Reactions.end(); it != end; ++it)
2872  removeReaction((*it)->getKey(), false);
2873 
2874  for (it = Metabolites.begin(), end = Metabolites.end(); it != end; ++it)
2875  removeMetabolite((*it)->getKey(), false);
2876 
2877  for (it = Compartments.begin(), end = Compartments.end(); it != end; ++it)
2878  removeCompartment((*it)->getKey(), false);
2879 
2880  for (it = Values.begin(), end = Values.end(); it != end; ++it)
2881  removeModelValue((*it)->getKey(), false);
2882 
2883  for (it = Events.begin(), end = Events.end(); it != end; ++it)
2884  removeEvent((*it)->getKey(), false);
2885 
2886  return;
2887 }
2888 
2889 bool CModel::removeModelValue(const size_t index,
2890  const bool & recursive)
2891 {
2892  const CModelValue * pMV = getModelValues()[index];
2893  return removeModelValue(pMV, recursive);
2894 }
2895 bool CModel::removeModelValue(const std::string & key,
2896  const bool & recursive)
2897 {
2898  CModelValue *pModelValue =
2899  dynamic_cast< CModelValue * >(CCopasiRootContainer::getKeyFactory()->get(key));
2900  return removeModelValue(pModelValue, recursive);
2901 }
2902 
2903 bool CModel::removeModelValue(const CModelValue * pModelValue,
2904  const bool & recursive)
2905 {
2906  if (!pModelValue)
2907  return false;
2908 
2909  if (recursive)
2910  {
2912  }
2913 
2914  //Check if Value with that name exists
2915  size_t index =
2916  mValues.CCopasiVector< CModelValue >::getIndex(pModelValue);
2917 
2918  if (index == C_INVALID_INDEX)
2919  return false;
2920 
2921  mValues.CCopasiVector< CModelValue >::remove(index);
2922 
2923  mCompileIsNecessary = true;
2924 
2925  return true;
2926 }
2927 
2928 CEvent* CModel::createEvent(const std::string & name)
2929 {
2930  if (mEvents.getIndex(name) != C_INVALID_INDEX)
2931  return NULL;
2932 
2933  CEvent * pEvent = new CEvent(name, this);
2934 
2935  if (!mEvents.add(pEvent, true))
2936  {
2937  delete pEvent;
2938  return NULL;
2939  }
2940 
2941  mCompileIsNecessary = true;
2942  return pEvent;
2943 }
2944 
2945 bool CModel::removeEvent(const size_t index,
2946  const bool & recursive)
2947 {
2948  const CEvent * pEvent = mEvents[index];
2949 
2950  return removeEvent(pEvent, recursive);
2951 }
2952 
2953 bool CModel::removeEvent(const std::string & key,
2954  const bool & recursive)
2955 {
2956  CEvent * pEvent = dynamic_cast< CEvent * >(CCopasiRootContainer::getKeyFactory()->get(key));
2957 
2958  return removeEvent(pEvent, recursive);
2959 }
2960 
2961 bool CModel::removeEvent(const CEvent * pEvent,
2962  const bool & /* recursive */)
2963 {
2964  if (!pEvent)
2965  return false;
2966 
2967  //Check if Event exists
2968  size_t index =
2969  mEvents.CCopasiVector< CEvent >::getIndex(pEvent);
2970 
2971  if (index == C_INVALID_INDEX)
2972  return false;
2973 
2974  mEvents.CCopasiVector< CEvent >::remove(index);
2975 
2976  clearMoieties();
2977 
2978  mCompileIsNecessary = true;
2979 
2980  return true;
2981 }
2982 
2983 #ifdef WITH_PE_EVENT_CREATION
2984 
2988 
2991 
2994 
2996 
2997 std::string getNextId(const std::string& base, int count)
2998 {
2999  std::stringstream str;
3000  str << base << count;
3001  return str.str();
3002 }
3003 
3004 const CCopasiObject*
3005 getDependentOrNull(const std::map< CCopasiObject *, size_t > & dependentMap, int index)
3006 {
3007 
3008  std::map< CCopasiObject *, size_t >::const_iterator it = dependentMap.begin();
3009 
3010  while (it != dependentMap.end())
3011  {
3012  if (it->second == index)
3013  return it->first;
3014 
3015  ++it;
3016  }
3017 
3018  return NULL;
3019 }
3020 
3021 bool
3022 CModel::createEventsForTimeseries(CExperiment* experiment/* = NULL*/)
3023 {
3024 
3025 #pragma region //find_experiment
3026 
3027  if (experiment == NULL)
3028  {
3029  // find experiment and invoke with it
3030  const CExperiment* theExperiment = NULL;
3031 
3032  const CCopasiDataModel* dataModel = getObjectDataModel();
3033  const CFitTask* task = dynamic_cast<const CFitTask*>((*dataModel->getTaskList())["Parameter Estimation"]);
3034 
3035  if (task == NULL)
3036  {
3038  "The parameter estimation was not correctly setup.");
3039  return false;
3040  }
3041 
3042  const CFitProblem *problem = static_cast<const CFitProblem*>(task->getProblem());
3043 
3044  const CExperimentSet& experiments = problem->getExperiementSet();
3045 
3046  // find first time course experiment
3047  for (size_t i = 0; i < experiments.size(); ++i)
3048  {
3049  const CExperiment* exp = experiments.getExperiment(i);
3050 
3052  {
3053  theExperiment = exp;
3054  break;
3055  }
3056  }
3057 
3058  // if still not found, bail
3059  if (theExperiment == NULL)
3060  {
3062  "No suitable experiment could be found, please first define a time course experiment and map the data.");
3063  return false;
3064  }
3065 
3066  if (experiments.size() > 1)
3067  {
3069  "You have defined multiple experiments, this function will only create events for the first time course experiment, unless another is specified.");
3070  }
3071 
3072  return createEventsForTimeseries(const_cast<CExperiment*>(theExperiment));
3073  }
3074 
3075 #pragma endregion //find_experiment
3076 
3077  if (experiment->getExperimentType() != CCopasiTask::timeCourse)
3078  {
3080  "The selected experiment, is not a time series experiment.");
3081  return false;
3082  }
3083 
3084  // need to get at the time course data
3085 
3086  std::ifstream File;
3087  File.open(CLocaleString::fromUtf8(experiment->getFileName()).c_str());
3088 
3089  size_t CurrentLine = 1;
3090 
3091  if (!experiment->read(File, CurrentLine))
3092  {
3094  "The data file could not be read.");
3095  return false;
3096  }
3097 
3098  if (!experiment->compile())
3099  {
3101  "The experiment could not be compiled.");
3102  return false;
3103  }
3104 
3105  // grab time column
3106  const CVector<double>& time = experiment->getTimeData();
3107  size_t numRows = experiment->getNumDataRows();
3108 
3109  if (numRows <= 1)
3110  {
3112  "Need at least 2 data rows in the experiment.");
3113  return false;
3114  }
3115 
3116  // remember whether we had events before, if so it could be that
3117  // we have two events triggering at the same time
3118  bool hadEvents = getEvents().size() > 0;
3119 
3120  size_t numCols = experiment->getNumColumns();
3121  const std::map< CCopasiObject *, size_t > & dependentMap = experiment->getDependentObjects();
3122 
3123  // then go through each time point
3124  for (size_t i = 0; i < numRows; ++i)
3125  {
3126  double current = time[i];
3127 
3128  // skip initial time
3129  if (current == 0) continue;
3130 
3131  CEvent* pEvent = createEvent(getNextId("pe_event_", i));
3132 
3133  if (pEvent == NULL)
3134  {
3136  "Could not create event, please verify that the events have not been created before.");
3137  return false;
3138  }
3139 
3140  std::stringstream trigger; trigger
3141  << "<" << getObject(CRegisteredObjectName("Reference=Time"))->getCN()
3142  << ">" << " > " << current;
3143  pEvent->setTriggerExpression(trigger.str());
3144  pEvent->getTriggerExpressionPtr()->compile();
3145 
3146  const CMatrix<double>& data = experiment->getDependentData();
3147 
3148  // create event and assignment for each mapping with its value
3149  for (size_t j = 0; j < data.numCols(); ++j)
3150  {
3151 
3152  const CCopasiObject* currentObject = getDependentOrNull(dependentMap, j); //objects[j + 1];
3153 
3154  if (currentObject == NULL || currentObject->getObjectParent() == NULL) continue;
3155 
3156  double value = data(i, j);
3157 
3158  // don't include missing data points
3159  if (value != value)
3160  {
3161  std::string displayName = currentObject->getObjectParent()->getObjectDisplayName();
3163  "At time %.2f: a missing data point was encountered for '%s', the value has been ignored."
3164  , current, displayName.c_str());
3165  continue;
3166  }
3167 
3168  CEventAssignment * pNewAssignment =
3169  new CEventAssignment(currentObject->getObjectParent()->getKey());
3170  std::stringstream assignmentStr; assignmentStr << value;
3171  pNewAssignment->setExpression(assignmentStr.str());
3172  pNewAssignment->getExpressionPtr()->compile();
3173  pEvent->getAssignments().add(pNewAssignment, true);
3174  }
3175  }
3176 
3177  // ensure that the 'continue on simultaneous events' option is enabled
3178  // for time course simulations, this needs to happen whenever we already
3179  // have events in the model, as then there is a chance that two events trigger at the same time
3180 
3181  if (!hadEvents) return true;
3182 
3183  CTrajectoryTask* task = dynamic_cast<CTrajectoryTask*>((*getObjectDataModel()->getTaskList())["Time-Course"]);
3184 
3185  if (task != NULL)
3186  {
3187  CTrajectoryProblem* problem = dynamic_cast<CTrajectoryProblem*>(task->getProblem());
3188 
3189  if (!problem->getContinueSimultaneousEvents())
3190  {
3191  problem->setContinueSimultaneousEvents(true);
3193  "Since the model contained events, the option 'Continue on Simultaneous Events' "
3194  "has been enabled in the 'Time 'Course' task to ensure that simulation continues "
3195  "if multiple events trigger at the same time.");
3196  }
3197  }
3198 
3199  return true;
3200 }
3201 #endif
3202 
3203 //*****************************************************************
3204 
3206 {
3207  // TODO check if there are any reversible reactions
3208  // TODO warn the user
3209  // TODO tell the GUI about changes -> not from here
3210  // TODO generate report ?
3211  // TODO check if newly generated reaction names are valid
3212  // TODO map, so that the same function is split only once
3213 
3214  bool success = true;
3215 
3216  std::vector< CReaction * > reactionsToDelete;
3217 
3218  CReaction *reac0, *reac1, *reac2;
3219  std::string fn, rn1, rn2;
3220 
3221  //CModel* model = dynamic_cast< CModel * >(CCopasiRootContainer::getKeyFactory()->get(objKey));
3222  //if (!model) return false;
3223 
3224  CCopasiVectorN< CReaction > & steps = this->getReactions();
3225 
3226  size_t i, imax = steps.size();
3227 
3228  for (i = 0; i < imax; ++i)
3229  if (steps[i]->isReversible())
3230  {
3231  std::vector< std::pair< const CCopasiObject * , const CCopasiObject * > > ParameterMap;
3232 
3233  reac0 = steps[i];
3234  rn1 = reac0->getObjectName() + " (forward)";
3235  rn2 = reac0->getObjectName() + " (backward)";
3236 
3237  fn = reac0->getFunction()->getObjectName();
3238 
3239  const CFunction* pFunc = reac0->getFunction();
3240 
3241  bool massaction = (fn == "Mass action (reversible)");
3242 
3243  std::pair<CFunction *, CFunction *> tmp;
3244 
3245  if (massaction)
3246  {
3247  //set functions to mass action (irrev)
3248  tmp.first = dynamic_cast<CFunction*>
3249  (CCopasiRootContainer::getFunctionList()-> findFunction("Mass action (irreversible)"));
3250  assert(tmp.first);
3251  tmp.second = tmp.first;
3252  }
3253  else //not mass action
3254  {
3255  //try splitting
3256  tmp = pFunc->splitFunction(NULL, pFunc->getObjectName() + " (forward part)",
3257  pFunc->getObjectName() + " (backward part)");
3258 
3259  if ((tmp.first == NULL) || (tmp.second == NULL))
3260  {
3261  // Create a message that the conversion for this reaction failed.
3263  reac0->getObjectName().c_str(), fn.c_str());
3264  success = false;
3265 
3266  pdelete(tmp.first);
3267  pdelete(tmp.second);
3268  continue;
3269  }
3270 
3271  tmp.first = CCopasiRootContainer::getFunctionList()->addAndAdaptName(tmp.first);
3272 
3273  tmp.second = CCopasiRootContainer::getFunctionList()->addAndAdaptName(tmp.second);
3274  }
3275 
3276  size_t i, imax;
3277 
3278  //**** create 1st reaction.
3279  reac1 = createReaction(rn1);
3280  reac1->setReversible(false);
3281  //substrates
3282  imax = reac0->getChemEq().getSubstrates().size();
3283 
3284  for (i = 0; i < imax; ++i)
3285  reac1->addSubstrate(reac0->getChemEq().getSubstrates()[i]->getMetaboliteKey(),
3286  reac0->getChemEq().getSubstrates()[i]->getMultiplicity());
3287 
3288  //products
3289  imax = reac0->getChemEq().getProducts().size();
3290 
3291  for (i = 0; i < imax; ++i)
3292  reac1->addProduct(reac0->getChemEq().getProducts()[i]->getMetaboliteKey(),
3293  reac0->getChemEq().getProducts()[i]->getMultiplicity());
3294 
3295  //function
3296  reac1->setFunction(tmp.first);
3297 
3298  //**** create 2nd reaction.
3299  reac2 = createReaction(rn2);
3300  reac2->setReversible(false);
3301  //substrates -> products
3302  imax = reac0->getChemEq().getSubstrates().size();
3303 
3304  for (i = 0; i < imax; ++i)
3305  reac2->addProduct(reac0->getChemEq().getSubstrates()[i]->getMetaboliteKey(),
3306  reac0->getChemEq().getSubstrates()[i]->getMultiplicity());
3307 
3308  //products -> substrates
3309  imax = reac0->getChemEq().getProducts().size();
3310 
3311  for (i = 0; i < imax; ++i)
3312  reac2->addSubstrate(reac0->getChemEq().getProducts()[i]->getMetaboliteKey(),
3313  reac0->getChemEq().getProducts()[i]->getMultiplicity());
3314 
3315  //function
3316  reac2->setFunction(tmp.second);
3317 
3318  //mapping for both reactions
3319  if (massaction)
3320  {
3321  // the parameter names of the massaction kinetics are hardcoded here.
3322  if (reac0->isLocalParameter("k1"))
3323  reac1->setParameterValue("k1", reac0->getParameterValue("k1"));
3324  else
3325  reac1->setParameterMapping("k1", reac0->getParameterMapping("k1")[0]);
3326 
3327  reac1->setParameterMappingVector("substrate", reac0->getParameterMapping("substrate"));
3328 
3329  if (reac0->isLocalParameter("k2"))
3330  reac2->setParameterValue("k1", reac0->getParameterValue("k2"));
3331  else
3332  reac2->setParameterMapping("k1", reac0->getParameterMapping("k2")[0]);
3333 
3334  reac2->setParameterMappingVector("substrate", reac0->getParameterMapping("product"));
3335 
3336  // Bug 2143: We need to replace any occurrence of the reaction parameter of the
3337  // reversible reaction with its replacement in the forward or backwards reaction
3338  ParameterMap.push_back(std::make_pair(reac0->getParameters().getParameter("k1")->getValueReference(),
3339  reac1->getParameters().getParameter("k1")->getValueReference()));
3340  ParameterMap.push_back(std::make_pair(reac0->getParameters().getParameter("k2")->getValueReference(),
3341  reac2->getParameters().getParameter("k1")->getValueReference()));
3342  }
3343  else //not mass action
3344  {
3345  const CFunctionParameters & fps = reac0->getFunctionParameters();
3346  imax = fps.size();
3347 
3348  for (i = 0; i < imax; ++i)
3349  {
3350  const CFunctionParameter * fp = fps[i];
3351  assert(fp);
3352  assert(fp->getType() == CFunctionParameter::FLOAT64);
3353 
3354  switch (fp->getUsage())
3355  {
3357  reac1->setParameterMapping(fp->getObjectName(),
3358  reac0->getParameterMapping(fp->getObjectName())[0]);
3359 
3360  // It is possible (see Bug 1830) that the split function may have additional modifier.
3361  // This will happen e.g. if the product is referenced in the forward part of the reaction.
3362  if (reac2->setParameterMapping(fp->getObjectName(),
3363  reac0->getParameterMapping(fp->getObjectName())[0]))
3364  {
3365  reac2->addModifier(reac0->getParameterMapping(fp->getObjectName())[0]);
3366  }
3367 
3368  break;
3369 
3371 
3372  // It is possible (see Bug 1830) that the split function may have additional modifier.
3373  // This will happen e.g. if the product is referenced in the forward part of the reaction.
3374  if (reac1->setParameterMapping(fp->getObjectName(),
3375  reac0->getParameterMapping(fp->getObjectName())[0]))
3376  {
3377  reac1->addModifier(reac0->getParameterMapping(fp->getObjectName())[0]);
3378  }
3379 
3380  reac2->setParameterMapping(fp->getObjectName(),
3381  reac0->getParameterMapping(fp->getObjectName())[0]);
3382  break;
3383 
3385 
3386  if (reac1->setParameterMapping(fp->getObjectName(),
3387  reac0->getParameterMapping(fp->getObjectName())[0]))
3388  {
3389  // Add the modifier
3390  reac1->addModifier(reac0->getParameterMapping(fp->getObjectName())[0]);
3391  }
3392 
3393  if (reac2->setParameterMapping(fp->getObjectName(),
3394  reac0->getParameterMapping(fp->getObjectName())[0]))
3395  {
3396  // Add the modifier
3397  reac2->addModifier(reac0->getParameterMapping(fp->getObjectName())[0]);
3398  }
3399 
3400  break;
3401 
3403 
3404  if (reac0->isLocalParameter(fp->getObjectName()))
3405  {
3406 
3407  reac1->setParameterValue(fp->getObjectName(),
3408  reac0->getParameterValue(fp->getObjectName()));
3409  reac2->setParameterValue(fp->getObjectName(),
3410  reac0->getParameterValue(fp->getObjectName()));
3411  }
3412  else
3413  {
3414  reac1->setParameterMapping(fp->getObjectName(),
3415  reac0->getParameterMapping(fp->getObjectName())[0]);
3416  reac2->setParameterMapping(fp->getObjectName(),
3417  reac0->getParameterMapping(fp->getObjectName())[0]);
3418  }
3419 
3420  // Bug 2143: We need to replace any occurrence of the reaction parameter of the
3421  // reversible reaction with its replacement in the forward or backwards reaction
3422  {
3423  CCopasiParameter * pOldParameter = reac0->getParameters().getParameter(fp->getObjectName());
3424  CCopasiParameter * pNewParameter = reac1->getParameters().getParameter(fp->getObjectName());
3425 
3426  if (pNewParameter == NULL)
3427  {
3428  pNewParameter = reac2->getParameters().getParameter(fp->getObjectName());
3429  }
3430 
3431  ParameterMap.push_back(std::make_pair(pOldParameter->getValueReference(), pNewParameter->getValueReference()));
3432  }
3433 
3434  break;
3435 
3436  default:
3437  reac1->setParameterMapping(fp->getObjectName(),
3438  reac0->getParameterMapping(fp->getObjectName())[0]);
3439  reac2->setParameterMapping(fp->getObjectName(),
3440  reac0->getParameterMapping(fp->getObjectName())[0]);
3441  break;
3442  }
3443  }
3444  }
3445 
3446  reac1->compile();
3447  reac2->compile();
3448 
3449  std::string Old;
3450  std::string New;
3451  std::string Infix;
3452 
3453  // Bug 2143: We need to replace any occurrence of the reaction parameter of the
3454  // reversible reaction with its replacement in the forward or backwards reaction
3455  std::vector< std::pair< const CCopasiObject *, const CCopasiObject * > >::const_iterator itParameter = ParameterMap.begin();
3456  std::vector< std::pair< const CCopasiObject *, const CCopasiObject * > >::const_iterator endParameter = ParameterMap.end();
3457 
3458  for (; itParameter != endParameter; ++itParameter)
3459  {
3460  Old = "<" + itParameter->first->getCN() + ">";
3461  New = "<" + itParameter->second->getCN() + ">";
3462 
3463  // We need to find all expressions which reference the old parameter (first value of the itParameter)
3464  // Compartments (expression, initial expression)
3467 
3468  for (; itComp != endComp; ++itComp)
3469  {
3470  const CExpression * pExpression = (*itComp)->getExpressionPtr();
3471 
3472  if (pExpression != NULL)
3473  {
3474  const CCopasiObject::DataObjectSet & Dependencies = pExpression->getDirectDependencies();
3475 
3476  if (Dependencies.find(itParameter->first) != Dependencies.end())
3477  {
3478  Infix = (*itComp)->getExpression();
3479  stringReplace(Infix, Old, New);
3480  (*itComp)->setExpression(Infix);
3481  }
3482  }
3483 
3484  pExpression = (*itComp)->getInitialExpressionPtr();
3485 
3486  if (pExpression != NULL)
3487  {
3488  const CCopasiObject::DataObjectSet & Dependencies = pExpression->getDirectDependencies();
3489 
3490  if (Dependencies.find(itParameter->first) != Dependencies.end())
3491  {
3492  Infix = (*itComp)->getInitialExpression();
3493  stringReplace(Infix, Old, New);
3494  (*itComp)->setInitialExpression(Infix);
3495  }
3496  }
3497  }
3498 
3499  // Species (expression, initial expression)
3502 
3503  for (; itMetab != endMetab; ++itMetab)
3504  {
3505  const CExpression * pExpression = (*itMetab)->getExpressionPtr();
3506 
3507  if (pExpression != NULL)
3508  {
3509  const CCopasiObject::DataObjectSet & Dependencies = pExpression->getDirectDependencies();
3510 
3511  if (Dependencies.find(itParameter->first) != Dependencies.end())
3512  {
3513  Infix = (*itMetab)->getExpression();
3514  stringReplace(Infix, Old, New);
3515  (*itMetab)->setExpression(Infix);
3516  }
3517  }
3518 
3519  pExpression = (*itMetab)->getInitialExpressionPtr();
3520 
3521  if (pExpression != NULL)
3522  {
3523  const CCopasiObject::DataObjectSet & Dependencies = pExpression->getDirectDependencies();
3524 
3525  if (Dependencies.find(itParameter->first) != Dependencies.end())
3526  {
3527  Infix = (*itMetab)->getInitialExpression();
3528  stringReplace(Infix, Old, New);
3529  (*itMetab)->setInitialExpression(Infix);
3530  }
3531  }
3532  }
3533 
3534  // Model Quantities (expression, initial expression)
3537 
3538  for (; itValue != endValue; ++itValue)
3539  {
3540  const CExpression * pExpression = (*itValue)->getExpressionPtr();
3541 
3542  if (pExpression != NULL)
3543  {
3544  const CCopasiObject::DataObjectSet & Dependencies = pExpression->getDirectDependencies();
3545 
3546  if (Dependencies.find(itParameter->first) != Dependencies.end())
3547  {
3548  Infix = (*itValue)->getExpression();
3549  stringReplace(Infix, Old, New);
3550  (*itValue)->setExpression(Infix);
3551  }
3552  }
3553 
3554  pExpression = (*itValue)->getInitialExpressionPtr();
3555 
3556  if (pExpression != NULL)
3557  {
3558  const CCopasiObject::DataObjectSet & Dependencies = pExpression->getDirectDependencies();
3559 
3560  if (Dependencies.find(itParameter->first) != Dependencies.end())
3561  {
3562  Infix = (*itValue)->getInitialExpression();
3563  stringReplace(Infix, Old, New);
3564  (*itValue)->setInitialExpression(Infix);
3565  }
3566  }
3567  }
3568 
3569  // Events (trigger, delay, assignments)
3572 
3573  for (; itEvent != endEvent; ++itEvent)
3574  {
3575  const CExpression * pExpression = (*itEvent)->getTriggerExpressionPtr();
3576 
3577  if (pExpression != NULL)
3578  {
3579  const CCopasiObject::DataObjectSet & Dependencies = pExpression->getDirectDependencies();
3580 
3581  if (Dependencies.find(itParameter->first) != Dependencies.end())
3582  {
3583  Infix = (*itEvent)->getTriggerExpression();
3584  stringReplace(Infix, Old, New);
3585  (*itEvent)->setTriggerExpression(Infix);
3586  }
3587  }
3588 
3589  pExpression = (*itEvent)->getDelayExpressionPtr();
3590 
3591  if (pExpression != NULL)
3592  {
3593  const CCopasiObject::DataObjectSet & Dependencies = pExpression->getDirectDependencies();
3594 
3595  if (Dependencies.find(itParameter->first) != Dependencies.end())
3596  {
3597  Infix = (*itEvent)->getDelayExpression();
3598  stringReplace(Infix, Old, New);
3599  (*itEvent)->setDelayExpression(Infix);
3600  }
3601  }
3602 
3603  CCopasiVector< CEventAssignment >::iterator itAssignment = (*itEvent)->getAssignments().begin();
3604  CCopasiVector< CEventAssignment >::iterator endAssignment = (*itEvent)->getAssignments().end();
3605 
3606  for (; itAssignment != endAssignment; ++itAssignment)
3607  {
3608  pExpression = (*itAssignment)->getExpressionPtr();
3609 
3610  if (pExpression != NULL)
3611  {
3612  const CCopasiObject::DataObjectSet & Dependencies = pExpression->getDirectDependencies();
3613 
3614  if (Dependencies.find(itParameter->first) != Dependencies.end())
3615  {
3616  Infix = (*itAssignment)->getExpression();
3617  stringReplace(Infix, Old, New);
3618  (*itAssignment)->setExpression(Infix);
3619  }
3620  }
3621  }
3622  }
3623  }
3624 
3625  // BUG 1848. We need to replace all references to the flux and particle flux
3626  // with the difference of the forward and backward reaction fluxes and particle fluxes, i.e,
3627  // flux = forward.flux - backward.flux
3628 
3629  Old = "<" + reac0->getFluxReference()->getCN() + ">";
3630  New = "(<" + reac1->getFluxReference()->getCN() + "> - <" + reac2->getFluxReference()->getCN() + ">)";
3631 
3632  // Find all objects which directly depend on the flux or particle flux.
3633  std::set< const CCopasiObject * > Flux;
3634  Flux.insert(reac0->getFluxReference());
3635  std::set< const CCopasiObject * > FluxDependents;
3636 
3637  // Initial Expression and Expression
3638  appendDependentCompartments(Flux, FluxDependents);
3639  appendDependentModelValues(Flux, FluxDependents);
3640  appendDependentMetabolites(Flux, FluxDependents);
3641 
3642  std::set< const CCopasiObject * >::iterator it = FluxDependents.begin();
3643  std::set< const CCopasiObject * >::iterator end = FluxDependents.end();
3644 
3645  for (; it != end; ++it)
3646  {
3647  CModelEntity * pEntity = static_cast< CModelEntity * >(const_cast< CCopasiObject * >(*it));
3648 
3649  // Expression
3650  Infix = pEntity->getExpression();
3651 
3652  if (stringReplace(Infix, Old, New))
3653  {
3654  pEntity->setExpression(Infix);
3655  }
3656 
3657  // Initial Expression
3658  if (pEntity->getStatus() != CModelEntity::ASSIGNMENT)
3659  {
3660  Infix = pEntity->getInitialExpression();
3661 
3662  if (stringReplace(Infix, Old, New))
3663  {
3664  pEntity->setInitialExpression(Infix);
3665  }
3666  }
3667  }
3668 
3669  FluxDependents.clear();
3670 
3671  // Trigger and Assignments
3672  appendDependentEvents(Flux, FluxDependents);
3673 
3674  it = FluxDependents.begin();
3675  end = FluxDependents.end();
3676 
3677  for (; it != end; ++it)
3678  {
3679  CEvent * pEvent = static_cast< CEvent * >(const_cast< CCopasiObject * >(*it));
3680 
3681  // Trigger Expression
3682  Infix = pEvent->getTriggerExpression();
3683 
3684  if (stringReplace(Infix, Old, New))
3685  {
3686  pEvent->setTriggerExpression(Infix);
3687  }
3688 
3689  // Delay Expression
3690  Infix = pEvent->getDelayExpression();
3691 
3692  if (stringReplace(Infix, Old, New))
3693  {
3694  pEvent->setDelayExpression(Infix);
3695  }
3696 
3697  // Priority Expression
3698  Infix = pEvent->getPriorityExpression();
3699 
3700  if (stringReplace(Infix, Old, New))
3701  {
3702  pEvent->setPriorityExpression(Infix);
3703  }
3704 
3705  // Assignments
3708 
3709  for (; itAssignment != endAssignment; ++itAssignment)
3710  {
3711  Infix = (*itAssignment)->getExpression();
3712 
3713  if (stringReplace(Infix, Old, New))
3714  {
3715  (*itAssignment)->setExpression(Infix);
3716  }
3717  }
3718  }
3719 
3720  FluxDependents.clear();
3721  Flux.clear();
3722 
3723  // particleFlux = forward.particleFlux - backward.particleFlux
3724  Old = "<" + reac0->getParticleFluxReference()->getCN() + ">";
3725  New = "(<" + reac1->getParticleFluxReference()->getCN() + "> - <" + reac2->getParticleFluxReference()->getCN() + ">)";
3726 
3727  Flux.insert(reac0->getParticleFluxReference());
3728 
3729  // Initial Expression and Expression
3730  appendDependentCompartments(Flux, FluxDependents);
3731  appendDependentModelValues(Flux, FluxDependents);
3732  appendDependentMetabolites(Flux, FluxDependents);
3733 
3734  it = FluxDependents.begin();
3735  end = FluxDependents.end();
3736 
3737  for (; it != end; ++it)
3738  {
3739  CModelEntity * pEntity = static_cast< CModelEntity * >(const_cast< CCopasiObject * >(*it));
3740 
3741  // Expression
3742  Infix = pEntity->getExpression();
3743 
3744  if (stringReplace(Infix, Old, New))
3745  {
3746  pEntity->setExpression(Infix);
3747  }
3748 
3749  // Initial Expression
3750  if (pEntity->getStatus() != CModelEntity::ASSIGNMENT)
3751  {
3752  Infix = pEntity->getInitialExpression();
3753 
3754  if (stringReplace(Infix, Old, New))
3755  {
3756  pEntity->setInitialExpression(Infix);
3757  }
3758  }
3759  }
3760 
3761  FluxDependents.clear();
3762 
3763  // Trigger and Assignments
3764  appendDependentEvents(Flux, FluxDependents);
3765 
3766  it = FluxDependents.begin();
3767  end = FluxDependents.end();
3768 
3769  for (; it != end; ++it)
3770  {
3771  CEvent * pEvent = static_cast< CEvent * >(const_cast< CCopasiObject * >(*it));
3772 
3773  // Trigger Expression
3774  Infix = pEvent->getTriggerExpression();
3775 
3776  if (stringReplace(Infix, Old, New))
3777  {
3778  pEvent->setTriggerExpression(Infix);
3779  }
3780 
3781  // Delay Expression
3782  Infix = pEvent->getDelayExpression();
3783 
3784  if (stringReplace(Infix, Old, New))
3785  {
3786  pEvent->setDelayExpression(Infix);
3787  }
3788 
3789  // Priority Expression
3790  Infix = pEvent->getPriorityExpression();
3791 
3792  if (stringReplace(Infix, Old, New))
3793  {
3794  pEvent->setPriorityExpression(Infix);
3795  }
3796 
3797  // Assignments
3800 
3801  for (; itAssignment != endAssignment; ++itAssignment)
3802  {
3803  Infix = (*itAssignment)->getExpression();
3804 
3805  if (stringReplace(Infix, Old, New))
3806  {
3807  (*itAssignment)->setExpression(Infix);
3808  }
3809  }
3810  }
3811 
3812  // Schedule the old reaction for removal.
3813  reactionsToDelete.push_back(reac0);
3814  }
3815 
3816  imax = reactionsToDelete.size();
3817 
3818  for (i = 0; i < imax; ++i)
3819  {
3820  delete reactionsToDelete[i];
3821  }
3822 
3823  if (imax != 0)
3824  {
3825  setCompileFlag(true);
3826  }
3827 
3828  return success;
3829 }
3830 
3831 //**********************************************************************
3832 
3834 {
3835  mKey = CCopasiRootContainer::getKeyFactory()->add("Model", this);
3836 
3837  // The regular CModelEntity mechanism does not work since
3838  // CModel is created before mStateTemplate :(
3839  C_FLOAT64 InitialValue = *mpIValue;
3840  C_FLOAT64 Value = *mpValue;
3841  pdelete(mpIValue);
3842  pdelete(mpValue);
3843  mStateTemplate.add(this);
3844  *mpIValue = InitialValue;
3845  *mpValue = Value;
3846 
3847  mpIValueReference->setObjectName("Initial Time");
3849 
3850  mRate = 1.0;
3851 
3852  addObjectReference("Comments", *const_cast<std::string *>(&getNotes()));
3853 
3854  // These are broken since they contain pointers to values :(
3855  // addVectorReference("Fluxes", mFluxes, CCopasiObject::ValueDbl);
3856  // addVectorReference("Particle Fluxes", mParticleFluxes, CCopasiObject::ValueDbl);
3857 
3859  addMatrixReference("Reduced Model Stoichiometry", mRedStoi, CCopasiObject::ValueDbl);
3860 
3862  addObjectReference("Quantity Unit", mQuantityUnit);
3863  addObjectReference("Quantity Conversion Factor", mQuantity2NumberFactor, CCopasiObject::ValueDbl);
3864  addObjectReference("Avogadro Constant", mAvogadro, CCopasiObject::ValueDbl);
3865 
3866  mpStoiAnnotation = new CArrayAnnotation("Stoichiometry(ann)", this, new CCopasiMatrixInterface<CMatrix<C_FLOAT64> >(&mStoi), true);
3867  mpStoiAnnotation->setDescription("Stoichiometry Matrix");
3869  mpStoiAnnotation->setDimensionDescription(0, "Species that are controlled by reactions");
3871  mpStoiAnnotation->setDimensionDescription(1, "Reactions");
3873 
3874  mpRedStoiAnnotation = new CArrayAnnotation("Reduced stoichiometry(ann)", this, new CCopasiMatrixInterface<CMatrix<C_FLOAT64> >(&mRedStoi), true);
3875  mpRedStoiAnnotation->setDescription("Reduced stoichiometry Matrix");
3877  mpRedStoiAnnotation->setDimensionDescription(0, "Species (reduced system)");
3881 
3882  mpLinkMatrixAnnotation = new CArrayAnnotation("Link matrix(ann)", this, new CCopasiMatrixInterface<CLinkMatrixView>(&mLView), true);
3883  mpLinkMatrixAnnotation->setDescription("Link matrix");
3885  mpLinkMatrixAnnotation->setDimensionDescription(0, "Species that are controlled by reactions (full system)");
3887  mpLinkMatrixAnnotation->setDimensionDescription(1, "Species (reduced system)");
3888 
3889  mpMathModel = new CMathModel(this);
3890 }
3891 
3893 {
3894  size_t i, imax = mSteps.size();
3895 
3896  for (i = 0; i < imax; ++i) if (mSteps[i]->isReversible()) return true;
3897 
3898  return false;
3899 }
3900 
3902 {
3903  size_t i, reactSize = mSteps.size();
3904  C_INT32 multInt;
3905  size_t j;
3906  C_FLOAT64 multFloat;
3907  // C_INT32 metabSize = mMetabolites->size();
3908 
3909  for (i = 0; i < reactSize; i++) // for every reaction
3910  {
3911  // TEST getCompartmentNumber() == 1
3912  //if (mSteps[i]->getCompartmentNumber() != 1) return - 1;
3913 
3914  // TEST isReversible() == 0
3915  if (mSteps[i]->isReversible() != 0)
3916  return "At least one reaction is reversible. That means stochastic simulation is not possible. \nYou can use \"Tools|Convert to irreversible\" which will split the reversible reactions \n into two irreversible reactions. However you should check the kinetics afterwards.";
3917 
3918  // TEST integer stoichiometry
3919  // Iterate through each the metabolites
3920  // Juergen: the number of rows of mStoiInternal equals the number of non-fixed metabs!
3921  // for (j=0; i<metabSize; j++)
3922  for (j = 0; j < mStoiInternal.numRows(); j++)
3923  {
3924  multFloat = mStoiInternal(j, i);
3925  multInt = static_cast<C_INT32>(floor(multFloat + 0.5)); // +0.5 to get a rounding out of the static_cast to int!
3926 
3927  if ((multFloat - multInt) > 0.01)
3928  return "Not all stoichiometries are integer numbers. \nThat means that discrete simulation is not possible.";
3929  }
3930  }
3931 
3932  for (i = 0; i < mMetabolites.size(); ++i)
3933  {
3935  return "At least one particle number in the initial state is too big.";
3936  }
3937 
3938  return ""; // Model is appropriate for hybrid simulation
3939 }
3940 
3941 #ifdef COPASI_DEBUG
3942 void CModel::check() const
3943 {}
3944 #endif
3945 
3947 {
3950  mStoi = mStoiInternal;
3951  mL.doRowPivot(mStoi);
3952 
3953  return;
3954 }
3955 
3956 const bool & CModel::isAutonomous() const
3957 {return mIsAutonomous;}
3958 
3960 {
3961  mIsAutonomous = true;
3962 
3963  // If the model is not empty we check whether anything depends on time
3964  if (mCompartments.size() != 0 ||
3965  mValues.size() != 0)
3966  {
3967  std::set< const CCopasiObject * > TimeDependent;
3968 
3969  appendDependentReactions(getDeletedObjects(), TimeDependent);
3973  appendDependentEvents(getDeletedObjects(), TimeDependent);
3974 
3975  mIsAutonomous = (TimeDependent.begin() == TimeDependent.end());
3976  }
3977 
3978  // An autonomous models always start simulation at T = 0
3979  if (mIsAutonomous)
3980  setInitialValue(0.0);
3981 }
3982 
3983 bool CModel::isStateVariable(const CCopasiObject * pObject) const
3984 {
3985  if (pObject == NULL)
3986  {
3987  return false;
3988  }
3989 
3990  // We check whether the object itself or the parent object is a state variable
3991  // A state variable is an independent model entity, a dependent species, or
3992  // a fixed entity, which is an event target.
3993 
3994  const CModelEntity * pEntity = dynamic_cast< const CModelEntity * >(pObject);
3995 
3996  if (pEntity == NULL)
3997  {
3998  pEntity = dynamic_cast< const CModelEntity * >(pObject->getObjectParent());
3999  }
4000 
4001  if (pEntity == NULL)
4002  {
4003  return false;
4004  }
4005 
4007  CModelEntity * const* end = mStateTemplate.endDependent();
4008 
4009  for (; it != end; ++it)
4010  {
4011  if (*it == pEntity)
4012  {
4013  return true;
4014  }
4015  }
4016 
4017  std::set< const CModelEntity * > EventTargets = CObjectLists::getEventTargets(this);
4018  std::set< const CModelEntity * >::const_iterator itSet = EventTargets.begin();
4019  std::set< const CModelEntity * >::const_iterator endSet = EventTargets.end();
4020 
4021  for (; itSet != endSet; ++itSet)
4022  {
4023  if (*itSet == pEntity)
4024  {
4025  return true;
4026  }
4027  }
4028 
4029  return false;
4030 }
4031 
4033 {
4034  // CModelEntities and derived classes are the only object which have initial and transient values
4035  // Note, for species we have distinguish between particle number and concentration.
4036 
4037  const CModelEntity * pEntity = dynamic_cast< const CModelEntity * >(pObject);
4038 
4039  if (pEntity == NULL)
4040  {
4041  pEntity = dynamic_cast< const CModelEntity * >(pObject->getObjectParent());
4042  }
4043 
4044  if (pEntity == NULL)
4045  {
4046  return const_cast< CCopasiObject * >(pObject);
4047  }
4048 
4049  const CMetab * pMetab = dynamic_cast< const CMetab * >(pEntity);
4050 
4051  if (pMetab != NULL && pMetab->getInitialConcentrationReference() == pObject)
4052  {
4053  return pMetab->getConcentrationReference();
4054  }
4055 
4056  return pEntity->getValueReference();
4057 }
4058 
4059 std::vector< const CEvaluationTree * > CModel::getTreesWithDiscontinuities() const
4060 {
4061  std::vector< const CEvaluationTree * > TreesWithDiscontinuities;
4062 
4063  // Check all expressions for entities of type ASSIGNMENT and ODE
4064  CModelEntity *const* ppEntity = mStateTemplate.getEntities();
4065  CModelEntity *const* ppEntityEnd = ppEntity + mStateTemplate.size();
4066 
4067  for (; ppEntity != ppEntityEnd; ++ppEntity)
4068  {
4069  switch ((*ppEntity)->getStatus())
4070  {
4071  case ODE:
4072  case ASSIGNMENT:
4073 
4074  if ((*ppEntity)->getExpressionPtr() &&
4075  (*ppEntity)->getExpressionPtr()->hasDiscontinuity())
4076  {
4077  TreesWithDiscontinuities.push_back((*ppEntity)->getExpressionPtr());
4078  }
4079 
4080  break;
4081 
4082  default:
4083  break;
4084  }
4085  }
4086 
4087  // Check all kinetic functions.
4090 
4091  for (; itReaction != endReaction; ++itReaction)
4092  {
4093  if ((*itReaction)->getFunction() &&
4094  (*itReaction)->getFunction()->hasDiscontinuity())
4095  {
4096  TreesWithDiscontinuities.push_back((*itReaction)->getFunction());
4097  }
4098  }
4099 
4100  // Check all event triggers
4103 
4104  for (; itEvent != endEvent; ++itEvent)
4105  {
4106  if ((*itEvent)->getTriggerExpressionPtr() &&
4107  (*itEvent)->getTriggerExpressionPtr()->hasDiscontinuity())
4108  {
4109  TreesWithDiscontinuities.push_back((*itEvent)->getTriggerExpressionPtr());
4110  }
4111  }
4112 
4113  return TreesWithDiscontinuities;
4114 }
4115 
4117 {
4118  bool success = true;
4119 
4120  std::vector< CCopasiContainer * > ListOfContainer;
4121 
4124 
4125  for (; it != end; ++ it)
4126  {
4127  success &= (*it)->compile(ListOfContainer);
4128  }
4129 
4130  return success;
4131 }
4132 
4133 const std::vector< Refresh * > & CModel::getListOfInitialRefreshes() const
4134 {return mInitialRefreshes;}
4135 
4136 const std::vector< Refresh * > & CModel::getListOfSimulatedRefreshes() const
4137 {return mSimulatedRefreshes;}
4138 
4139 const std::vector< Refresh * > & CModel::getListOfConstantRefreshes() const
4141 
4142 const std::vector< Refresh * > & CModel::getListOfNonSimulatedRefreshes() const
4143 {return mNonSimulatedRefreshes;}
4144 
4145 void
4146 CModel::updateInitialValues(std::set< const CCopasiObject * > & changedObjects)
4147 {
4148  std::vector< Refresh * > refreshes = buildInitialRefreshSequence(changedObjects);
4149  std::vector< Refresh * >::iterator it = refreshes.begin(), endIt = refreshes.end();
4150 
4151  while (it != endIt)
4152  (**it++)();
4153 }
4154 
4155 void
4157 {
4158  std::set<const CCopasiObject*> changedObjects;
4159  changedObjects.insert(changedObject);
4160  updateInitialValues(changedObjects);
4161 }
4162 
4163 std::vector< Refresh * >
4164 CModel::buildInitialRefreshSequence(std::set< const CCopasiObject * > & changedObjects)
4165 {
4166  // First we remove all objects which are of type assignment from the changed objects
4167  // since this may not be changed as they are under control of the assignment.
4168  std::set< const CCopasiObject * >::iterator itSet;
4169  std::set< const CCopasiObject * >::iterator endSet;
4170  std::set< const CCopasiObject * > Objects;
4171 
4172  CModelEntity **ppEntity;
4173  CModelEntity **ppEntityEnd = mStateTemplate.endFixed();
4174 
4175  const CModelEntity * pEntity;
4176  CMetab * pMetab;
4177 
4178  // If the changed objects are empty we assume that all changeable objects have been changed
4179  if (changedObjects.size() == 0)
4180  {
4181  // The objects which are changed are all initial values of of all model entities including
4182  // fixed and unused once. Additionally, all kinetic parameters are possibly changed.
4183  // This is basically all the parameters in the parameter overview whose value is editable.
4184 
4185  // :TODO: Theoretically, it is possible that also task parameters influence the initial
4186  // state of a model but that is currently not handled.
4187 
4188  // The initial values of the model entities
4189  ppEntity = mStateTemplate.beginIndependent() - 1; // Offset for time
4190 
4191  for (; ppEntity != ppEntityEnd; ++ppEntity)
4192  {
4193  // If we have an initial expression we have no initial values
4194  if (((*ppEntity)->getInitialExpression() != "" ||
4195  (*ppEntity)->getStatus() == ASSIGNMENT) &&
4196  (*ppEntity)->getInitialValueReference()->getDirectDependencies().size() > 0)
4197  continue;
4198 
4199  // Metabolites have two initial values
4200  if ((pMetab = dynamic_cast< CMetab * >(*ppEntity)) != NULL)
4201  {
4202  // The concentration is assumed to be fix accept when this would lead to circular dependencies,
4203  // for the parent's compartment's initial volume.
4205  changedObjects.insert(pMetab->getInitialConcentrationReference());
4206  else
4207  changedObjects.insert(pMetab->getInitialValueReference());
4208  }
4209  else
4210  changedObjects.insert((*ppEntity)->getInitialValueReference());
4211  }
4212 
4213  // The reaction parameters
4216  size_t i, imax;
4217 
4218  for (; itReaction != endReaction; ++itReaction)
4219  {
4220  const CCopasiParameterGroup & Group = (*itReaction)->getParameters();
4221 
4222  for (i = 0, imax = Group.size(); i < imax; i++)
4223  changedObjects.insert(Group.getParameter(i));
4224  }
4225 
4226  // Fix for Issue 1170: We need to add elements of the stoichiometry, reduced stoichiometry,
4227  // and link matrices.
4228  if (mpStoiAnnotation != NULL)
4229  mpStoiAnnotation->appendElementReferences(changedObjects);
4230 
4231  if (mpRedStoiAnnotation != NULL)
4233 
4234  if (mpLinkMatrixAnnotation != NULL)
4236  }
4237  else
4238  {
4239  // Remove all objects with initial assignments
4240  itSet = changedObjects.begin();
4241  endSet = changedObjects.end();
4242 
4243  for (; itSet != endSet; ++itSet)
4244  {
4245  if ((pEntity = dynamic_cast< const CModelEntity * >((*itSet)->getObjectParent())) != NULL &&
4246  (pEntity->getInitialExpression() != "" ||
4247  pEntity->getStatus() == ASSIGNMENT) &&
4248  pEntity->getInitialValueReference()->getDirectDependencies().size() > 0)
4249  Objects.insert(*itSet);
4250  }
4251 
4252  for (itSet = Objects.begin(), endSet = Objects.end(); itSet != endSet; ++itSet)
4253  changedObjects.erase(*itSet);
4254  }
4255 
4256  // We need to add all initial values which are dynamically calculated.
4257  // These are initial assignments and either a metabolite's initial particle number
4258  // or concentration.
4259  ppEntity = mStateTemplate.beginIndependent() - 1; // Offset for time
4260  Objects.clear();
4261 
4262  for (; ppEntity != ppEntityEnd; ++ppEntity)
4263  {
4264  if (((*ppEntity)->getInitialExpression() != "" ||
4265  (*ppEntity)->getStatus() == ASSIGNMENT))
4266  {
4267  Objects.insert((*ppEntity)->getInitialValueReference());
4268  continue;
4269  }
4270 
4271  // For metabolites we need to add the initial concentration or the initial
4272  // particle number..
4273  if ((pMetab = dynamic_cast< CMetab * >(*ppEntity)) != NULL)
4274  {
4275  if (changedObjects.count(pMetab->getInitialConcentrationReference()) != 0)
4276  {
4277  Objects.insert(pMetab->getInitialValueReference());
4278  }
4279  else
4280  {
4281  Objects.insert(pMetab->getInitialConcentrationReference());
4282  }
4283  }
4284  }
4285 
4286  // We need to add the total particle number of moieties.
4289 
4290  for (; itMoiety != endMoiety; ++itMoiety)
4291  Objects.insert((*itMoiety)->getInitialValueReference());
4292 
4293  std::set< const CCopasiObject * > DependencySet;
4294  std::set< const CCopasiObject * > VerifiedSet;
4295  std::pair<std::set< const CCopasiObject * >::iterator, bool> InsertedObject;
4296 
4297  assert(Objects.count(NULL) == 0);
4298 
4299  // Check whether we have any circular dependencies
4300  for (itSet = Objects.begin(), endSet = Objects.end(); itSet != endSet; ++itSet)
4301  if ((*itSet)->hasCircularDependencies(DependencySet, VerifiedSet, changedObjects))
4302  CCopasiMessage(CCopasiMessage::EXCEPTION, MCObject + 1, (*itSet)->getCN().c_str());
4303 
4304  // Build the complete set of dependencies
4305  for (itSet = Objects.begin(); itSet != endSet; ++itSet)
4306  {
4307  // At least the object itself needs to be up to date.
4308  InsertedObject = DependencySet.insert(*itSet);
4309 
4310  // Add all its dependencies
4311  if (InsertedObject.second)
4312  (*itSet)->getAllDependencies(DependencySet, changedObjects);
4313  }
4314 
4315  // Remove all objects which do not depend on the changed objects, or do not have a
4316  // refresh method.
4317  Objects.clear();
4318 
4319  // We now check which objects we need to refresh
4320  for (itSet = DependencySet.begin(), endSet = DependencySet.end(); itSet != endSet; ++itSet)
4321  {
4322  // No refresh method
4323  if ((*itSet)->getRefresh() == NULL)
4324  Objects.insert(*itSet);
4325  // Is a changed object
4326  else if (changedObjects.count(*itSet) != 0)
4327  Objects.insert(*itSet);
4328  // Not dependent on the changed objects.
4329  else if (!(*itSet)->dependsOn(changedObjects, changedObjects))
4330  Objects.insert(*itSet);
4331  }
4332 
4333  for (itSet = Objects.begin(), endSet = Objects.end(); itSet != endSet; ++itSet)
4334  DependencySet.erase(*itSet);
4335 
4336  // Create a properly sorted list.
4337  std::list< const CCopasiObject * > SortedList =
4338  sortObjectsByDependency(DependencySet.begin(), DependencySet.end(), changedObjects);
4339 
4340  std::list< const CCopasiObject * >::iterator itList;
4341  std::list< const CCopasiObject * >::iterator endList;
4342 
4343  // Build the vector of pointers to refresh methods
4344  Refresh * pRefresh;
4345 
4346  std::vector< Refresh * > UpdateVector;
4347  std::vector< Refresh * >::const_iterator itUpdate;
4348  std::vector< Refresh * >::const_iterator endUpdate;
4349 
4350  itList = SortedList.begin();
4351  endList = SortedList.end();
4352 
4353  for (; itList != endList; ++itList)
4354  {
4355  pRefresh = (*itList)->getRefresh();
4356  itUpdate = UpdateVector.begin();
4357  endUpdate = UpdateVector.end();
4358 
4359  while (itUpdate != endUpdate && !(*itUpdate)->isEqual(pRefresh)) ++itUpdate;
4360 
4361  if (itUpdate == endUpdate)
4362  UpdateVector.push_back(pRefresh);
4363  }
4364 
4365  return UpdateVector;
4366 }
4367 
4368 CVector< C_FLOAT64 > CModel::initializeAtolVector(const C_FLOAT64 & atol, const bool & reducedModel) const
4369 {
4370  CVector< C_FLOAT64 > Atol;
4371 
4372  if (reducedModel)
4374  else
4376 
4377  C_FLOAT64 * pAtol = Atol.array();
4378  C_FLOAT64 * pEnd = pAtol + Atol.size();
4379 
4380  C_FLOAT64 InitialValue;
4381  C_FLOAT64 Limit;
4382 
4383  CModelEntity *const* ppEntity = getStateTemplate().beginIndependent();
4384  const CMetab * pMetab;
4385 
4386  for (; pAtol != pEnd; ++pAtol, ++ppEntity)
4387  {
4388  *pAtol = atol;
4389 
4390  InitialValue = fabs((*ppEntity)->getInitialValue());
4391 
4392  if ((pMetab = dynamic_cast< const CMetab * >(*ppEntity)) != NULL)
4393  {
4394  Limit =
4396 
4397  if (InitialValue != 0.0)
4398  *pAtol *= std::min(Limit, InitialValue);
4399  else
4400  *pAtol *= std::max(1.0, Limit);
4401  }
4402  else if (InitialValue != 0.0)
4403  *pAtol *= std::min(1.0, InitialValue);
4404  }
4405 
4406  return Atol;
4407 }
4408 
4409 #include "utilities/CDimension.h"
4410 
4412 {
4413  std::ostringstream oss;
4414  CModel* model = this;
4415 
4416  oss << "Initial time: " << model->getInitialTime() << " " << model->getTimeUnitName() << std::endl;
4417 
4418  oss << std::endl;
4419 
4420  size_t i, imax, j, jmax;
4421 
4422  //Compartments
4423  const CCopasiVector< CCompartment > & comps = model->getCompartments();
4424  imax = comps.size();
4425 
4426  if (imax)
4427  {
4428  oss << "Initial volumes:\n\n";
4429 
4430  for (i = 0; i < imax; ++i)
4431  oss << comps[i]->getObjectName() << " \t" << comps[i]->getInitialValue()
4432  << " " << model->getVolumeUnitsDisplayString() << "\n";
4433 
4434  oss << "\n";
4435  }
4436 
4437  //Species
4438  const CCopasiVector< CMetab > & metabs = model->getMetabolites();
4439  imax = metabs.size();
4440 
4441  if (imax)
4442  {
4443  oss << "Initial concentrations:\n\n";
4444 
4445  for (i = 0; i < imax; ++i)
4446  oss << CMetabNameInterface::getDisplayName(model, *metabs[i], false) << " \t"
4447  << metabs[i]->getInitialConcentration() << " "
4448  << model->getConcentrationUnitsDisplayString() << "\n";
4449 
4450  oss << "\n";
4451  }
4452 
4453  //global Parameters
4454  const CCopasiVector< CModelValue > & params = model->getModelValues();
4455  imax = params.size();
4456 
4457  if (imax)
4458  {
4459  oss << "Initial values of global quantities:\n\n";
4460 
4461  for (i = 0; i < imax; ++i)
4462  oss << params[i]->getObjectName() << " \t"
4463  << params[i]->getInitialValue() << "\n";
4464 
4465  oss << "\n";
4466  }
4467 
4468  //Reactions
4469  const CCopasiVector< CReaction > & reacs = model->getReactions();
4470  imax = reacs.size();
4471 
4472  if (imax)
4473  {
4474  oss << "Reaction parameters:\n\n";
4475  CReaction* reac;
4476 
4477  for (i = 0; i < imax; ++i)
4478  {
4479  reac = reacs[i];
4480  oss << reac->getObjectName() << "\n";
4481 
4482  //calculate units
4488  units.setUseHeuristics(true);
4489  units.setChemicalEquation(&reac->getChemEq());
4490  units.findDimensions(reac->getCompartmentNumber() > 1);
4491 
4492  const CFunctionParameters & params = reac->getFunctionParameters();
4493  jmax = params.size();
4494 
4495  for (j = 0; j < jmax; ++j)
4496  if (params[j]->getUsage() == CFunctionParameter::PARAMETER)
4497  {
4499 
4500  if (!obj) continue;
4501 
4502  if (reac->isLocalParameter(j))
4503  {
4504  CCopasiParameter * par = dynamic_cast<CCopasiParameter*>(obj); //must be a CCopasiParameter
4505 
4506  if (!par) continue; //or rather fatal error?
4507 
4508  oss << " " << params[j]->getObjectName() << " \t"
4509  << *par->getValue().pDOUBLE << " "
4510  << units.getDimensions()[j].getDisplayString(this) << "\n";
4511  }
4512  else
4513  {
4514  CModelValue * par = dynamic_cast<CModelValue*>(obj); //must be a CModelValue
4515 
4516  if (!par) continue; //or rather fatal error?
4517 
4518  oss << " " << params[j]->getObjectName() << " \t"
4519  << "-> " + par->getObjectName()
4520  << " (" << units.getDimensions()[j].getDisplayString(this) << ")\n";
4521  }
4522  }
4523 
4524  oss << "\n";
4525  }
4526  }
4527 
4528  return oss.str();
4529 }
4530 
4532 {
4534  return "";
4535 
4536  return TimeUnitNames[mTimeUnit];
4537 }
4538 
4540 {
4542  return "";
4543 
4544  return std::string("1/") + TimeUnitNames[mTimeUnit];
4545 }
4546 
4548 {
4550  return "";
4551 
4552  return VolumeUnitNames[mVolumeUnit];
4553 }
4554 
4556 {
4558  return "";
4559 
4560  return AreaUnitNames[mAreaUnit];
4561 }
4562 
4564 {
4566  return "";
4567 
4568  return LengthUnitNames[mLengthUnit];
4569 }
4570 
4572 {
4574  {
4576  return "";
4577 
4578  return std::string("1/") + TimeUnitNames[mTimeUnit];
4579  }
4580 
4582  return VolumeUnitNames[mVolumeUnit];
4583 
4584  return std::string(VolumeUnitNames[mVolumeUnit]) + "/" + TimeUnitNames[mTimeUnit];
4585 }
4586 
4588 {
4589  std::string Units;
4590 
4592  {
4594  return "";
4595 
4596  return std::string("1/") + VolumeUnitNames[mVolumeUnit];
4597  }
4598 
4600 
4602  return Units;
4603 
4604  return Units + "/" + VolumeUnitNames[mVolumeUnit];
4605 }
4606 
4608 {
4609  std::string Units;
4610 
4612  {
4613  Units = "1";
4614 
4616  {
4618  return "";
4619 
4620  return Units + "/" + TimeUnitNames[mTimeUnit];
4621  }
4622  else
4623  {
4625  return Units + "/" + VolumeUnitNames[mVolumeUnit];
4626 
4627  return Units + "/(" + VolumeUnitNames[mVolumeUnit] + "*" + TimeUnitNames[mTimeUnit] + ")";
4628  }
4629  }
4630 
4632 
4634  {
4636  return Units;
4637 
4638  return Units + "/" + TimeUnitNames[mTimeUnit];
4639  }
4640 
4642  return Units + "/" + VolumeUnitNames[mVolumeUnit];
4643 
4644  return Units + "/(" + VolumeUnitNames[mVolumeUnit] + "*" + TimeUnitNames[mTimeUnit] + ")";
4645 }
4646 
4648 {
4650  {
4651  return "";
4652  }
4653 
4655 }
4656 
4658 {
4659  std::string Units;
4660 
4662  {
4664  return "";
4665 
4666  return std::string("1/") + TimeUnitNames[mTimeUnit];
4667  }
4668 
4670 
4672  return Units;
4673 
4674  return Units + "/" + TimeUnitNames[mTimeUnit];
4675 }
4676 
4677 /****** Below will be removed when the math model completed ******/
4678 
4680  const bool & ignoreDiscrete)
4681 {
4682  return mpMathModel->evaluateRoots(rootValues, ignoreDiscrete);
4683 }
4684 
4686  const bool & equality,
4687  CProcessQueue::resolveSimultaneousAssignments pResolveSimultaneousAssignments)
4688 {
4689  return mpMathModel->processQueue(time, equality, pResolveSimultaneousAssignments);
4690 }
4691 
4693  const bool & equality,
4694  const bool & correct,
4695  const CVector< C_INT > & roots)
4696 {
4697  mpMathModel->processRoots(time, equality, correct, roots);
4698 
4699  return;
4700 }
4701 
4703 {
4705 }
4706 
4707 size_t CModel::getNumRoots() const
4708 {
4709  return mpMathModel->getNumRoots();
4710 }
4711 
4713 {
4714  return mpMathModel->calculateRootDerivatives(rootDerivatives);
4715 }
4716 
4718 {
4719  return mpMathModel->getRootFinders();
4720 }
4721 
4723 {return mpMathModel;}
4724 
4726 {return mpMathModel;}
4727 
4728 #ifdef USE_MATH_CONTAINER
4729 const CMathContainer* CModel::getMathContainer() const
4730 {return mpMathContainer;}
4731 
4732 CMathContainer* CModel::getMathContainer()
4733 {return mpMathContainer;}
4734 #endif // USE_MATH_CONTAINER
std::string unQuote(const std::string &name)
Definition: utility.cpp:120
bool stringReplace(std::string &str, const std::string &target, const std::string &replacement)
Definition: utility.cpp:165
Definition: CEvent.h:152
CCopasiDataModel * getObjectDataModel()
virtual Refresh * getRefresh() const
static const char * LengthUnitNames[]
Definition: CModel.h:81
void setInitialConcentration(const C_FLOAT64 &initialConcentration)
Definition: CMetab.cpp:257
#define C_INT
Definition: copasi.h:115
bool remove(const std::string &key)
void initObjects()
Definition: CModel.cpp:3833
CArrayAnnotation * mpRedStoiAnnotation
Definition: CModel.h:1399
CModel::AreaUnit getAreaUnitEnum() const
Definition: CModel.cpp:2198
const C_FLOAT64 & getAvogadro() const
Definition: CModel.cpp:2349
void setReversible(bool reversible)
Definition: CReaction.cpp:247
bool build(const CMatrix< C_FLOAT64 > &matrix, size_t maxRank=C_INVALID_INDEX)
Definition: CLinkMatrix.cpp:44
const std::string & getFileName() const
const CCopasiVectorN< CEventAssignment > & getAssignments() const
Definition: CEvent.cpp:678
static const char * ModelTypeNames[]
Definition: CModel.h:117
const CModelParameterSet & getModelParameterSet() const
Definition: CModel.cpp:1072
CFunction * addAndAdaptName(CFunction *pFunction)
const C_FLOAT64 & getProcessQueueExecutionTime() const
Definition: CMathModel.cpp:357
virtual std::string getObjectDisplayName(bool regular=true, bool richtext=false) const
#define pdelete(p)
Definition: copasi.h:215
Header file of class CModelEntity and CModelValue.
bool setAreaUnit(const std::string &name)
Definition: CModel.cpp:2182
const CFunctionParameter::DataType & getType() const
CCopasiObjectReference< C_FLOAT64 > * mpIValueReference
Definition: CModelValue.h:356
CEvent * createEvent(const std::string &name)
Definition: CModel.cpp:2928
void setUserOrder(const CVector< CModelEntity * > &userOrder)
Definition: CState.cpp:189
CCopasiObject * getCorrespondingTransientObject(const CCopasiObject *pObject) const
Definition: CModel.cpp:4032
CCopasiObject * getParticleFluxReference()
Definition: CReaction.cpp:207
virtual void load(CReadConfig &configbuffer, size_t size)
CCopasiVectorN< CEvent > & getEvents()
Definition: CModel.cpp:1110
bool setVolumeUnit(const std::string &name)
Definition: CModel.cpp:2159
CCopasiProblem * getProblem()
std::string getTimeUnitsDisplayString() const
Definition: CModel.cpp:4531
#define K
LengthUnit mLengthUnit
Definition: CModel.h:1304
virtual CCopasiObjectName getCN() const
bool addSubstrate(const std::string &metabKey, const C_FLOAT64 &multiplicity=1.0)
Definition: CReaction.cpp:232
void setContinueSimultaneousEvents(const bool &continueSimultaneousEvents)
CCopasiVectorN< CEvent > mEvents
Definition: CModel.h:1344
const CLinkMatrixView & getL() const
Definition: CModel.cpp:1166
const CMatrix< C_FLOAT64 > & getDependentData() const
const CCopasiVector< CMetab > & getMetabolites() const
Definition: CModel.cpp:1051
bool mIsAutonomous
Definition: CModel.h:1517
std::string getConcentrationRateUnitsDisplayString() const
Definition: CModel.cpp:4607
C_FLOAT64 calculateDivergence() const
Definition: CModel.cpp:2153
const std::string & getObjectName() const
bool buildUserOrder()
Definition: CModel.cpp:1349
size_t getNumIndependent() const
Definition: CState.cpp:222
std::string getFrequencyUnitsDisplayString() const
Definition: CModel.cpp:4539
const CCopasiVectorN< CModelValue > & getModelValues() const
Definition: CModel.cpp:1060
virtual void setStatus(const CModelEntity::Status &status)
Definition: CMetab.cpp:291
virtual CCopasiObjectName getCN() const =0
void add(CModelEntity *entity)
Definition: CState.cpp:63
CModel::QuantityUnit getQuantityUnitEnum() const
Definition: CModel.cpp:2331
size_t mNumMetabolitesReactionIndependent
Definition: CModel.h:1441
virtual size_t size() const
std::string getQuantityUnitOldXMLName() const
Definition: CModel.cpp:2326
const CVector< CMathTrigger::CRootFinder * > & getRootFinders() const
Definition: CModel.cpp:4717
C_FLOAT64 mRate
Definition: CModelValue.h:331
bool removeMetabolite(const std::string &key, const bool &recursive=true)
Definition: CModel.cpp:2667
C_FLOAT64 * mpValue
Definition: CModelValue.h:321
bool setLengthUnit(const std::string &name)
Definition: CModel.cpp:2204
void buildRedStoi()
Definition: CModel.cpp:829
void evaluateRoots(CVectorCore< C_FLOAT64 > &rootValues, const bool &ignoreDiscrete)
Definition: CModel.cpp:4679
CCopasiObject * get(const std::string &key)
void updateNonSimulatedValues(void)
Definition: CModel.cpp:1892
CModelParameterSet mParameterSet
Definition: CModel.h:1359
CMatrix< C_FLOAT64 > mElasticities
Definition: CModel.h:1404
size_t getNumMetabs() const
Definition: CModel.cpp:1118
AreaUnit
Definition: CModel.h:66
const CCopasiMessage & getMessage() const
const CVector< C_FLOAT64 > & getParticleFlux() const
Definition: CModel.cpp:1045
CCopasiObject * getInitialValueReference() const
static CRenameHandler * smpRenameHandler
const std::vector< Refresh * > & getListOfConstantRefreshes() const
Definition: CModel.cpp:4139
virtual bool setName(const std::string &name)
void setInitialState(const CState &state)
Definition: CModel.cpp:1774
bool setTimeUnit(const std::string &name)
Definition: CModel.cpp:2227
void applyActiveParameterSet()
Definition: CModel.cpp:1082
iterator begin()
void updateSimulatedValues(const bool &updateMoieties)
Definition: CModel.cpp:1851
virtual size_t numRows() const
Definition: CMatrix.h:138
void setParameterValue(const std::string &parameterName, const C_FLOAT64 &value, const bool &updateStatus=true)
Definition: CReaction.cpp:303
void calculateJacobianX(CMatrix< C_FLOAT64 > &jacobianX, const C_FLOAT64 &derivationFactor, const C_FLOAT64 &resolution)
Definition: CModel.cpp:2082
bool removeCompartment(const size_t index, const bool &recursive=true)
Definition: CModel.cpp:2720
#define fatalError()
virtual void load(CReadConfig &configbuffer, size_t size)
std::string getTimeUnitName() const
Definition: CModel.cpp:2238
void setNotes(const std::string &notes)
bool buildNonSimulatedSequence()
Definition: CModel.cpp:1678
std::string getQuantityRateUnitsDisplayString() const
Definition: CModel.cpp:4657
CCopasiObject * getValueReference() const
CCopasiVectorNS< CCompartment > mCompartments
Definition: CModel.h:1324
const size_t & getNumber() const
const CExpression * getExpressionPtr() const
Definition: CEvent.cpp:226
void setAnnotationCN(size_t d, size_t i, const std::string cn)
void evaluateRoots(CVectorCore< C_FLOAT64 > &rootValues, const bool &ignoreDiscrete)
Definition: CMathModel.cpp:169
bool setExpression(const std::string &expression)
Definition: CEvent.cpp:167
Definition: CState.h:305
#define MNumMetabolitesReactionDependent
Definition: CModel.cpp:64
bool addProduct(const std::string &metabKey, const C_FLOAT64 &multiplicity=1.0)
Definition: CReaction.cpp:236
bool isLocalParameter(const size_t &index) const
Definition: CReaction.cpp:449
std::string printParameterOverview()
Definition: CModel.cpp:4411
#define CCHECK
Definition: CModel.cpp:61
std::string getVolumeRateUnitsDisplayString() const
Definition: CModel.cpp:4571
const CMathModel * getMathModel() const
Definition: CModel.cpp:4722
#define MCObject
CLinkMatrixView mLView
Definition: CModel.h:1456
void resize(size_t size, const bool &copy=false)
Definition: CVector.h:301
const std::set< const CCopasiObject * > & getUptoDateObjects() const
Definition: CModel.cpp:1178
bool setInitialExpression(const std::string &expression)
#define C_INVALID_INDEX
Definition: copasi.h:222
std::string getExpression() const
void applyInitialValues()
Definition: CMathModel.cpp:372
std::string getTriggerExpression() const
Definition: CEvent.cpp:524
const CExperimentSet & getExperiementSet() const
const std::map< CCopasiObject *, size_t > & getDependentObjects() const
CExperiment * getExperiment(const size_t &index)
void calculateJacobian(CMatrix< C_FLOAT64 > &jacobian, const C_FLOAT64 &derivationFactor, const C_FLOAT64 &resolution)
Definition: CModel.cpp:2006
void setDescription(const std::string &s)
virtual bool compile(std::vector< CCopasiContainer * > listOfContainer=CCopasiContainer::EmptyList)
Definition: CExpression.cpp:97
CModelValue * createModelValue(const std::string &name, const C_FLOAT64 &value=0.0)
Definition: CModel.cpp:2838
virtual size_t getIndex(const std::string &name) const
size_t getNumODEMetabs() const
Definition: CModel.cpp:1124
CModelEntity ** endDependent()
Definition: CState.cpp:211
void clearPivoting()
CMatrix< C_FLOAT64 > mStoi
Definition: CModel.h:1389
#define C_INT32
Definition: copasi.h:90
CMatrix< C_FLOAT64 > mRedStoi
Definition: CModel.h:1394
static std::string getDisplayName(const CModel *model, const std::string &key, const bool &quoted)
void setParameterMappingVector(const std::string &parameterName, const std::vector< std::string > &keys)
Definition: CReaction.cpp:390
bool buildStateTemplate()
Definition: CModel.cpp:1262
CLinkMatrix mL
Definition: CModel.h:1446
void updateMoietyValues()
Definition: CModel.cpp:896
Definition: CMetab.h:178
const size_t & getNumIndependent() const
const size_t & size() const
Definition: CState.cpp:242
void applyInitialValues()
Definition: CModel.cpp:1236
void setUseHeuristics(bool flag)
Definition: CDimension.cpp:330
std::vector< Refresh * > mNonSimulatedRefreshes
Definition: CModel.h:1507
CCopasiVectorN< CModelValue > mValues
Definition: CModel.h:1354
const std::vector< Refresh * > & getListOfSimulatedRefreshes() const
Definition: CModel.cpp:4136
std::vector< Refresh * > DataUpdateSequence
const CMatrix< C_FLOAT64 > & getRedStoi() const
Definition: CModel.cpp:1154
virtual std::set< const CCopasiObject * > getDeletedObjects() const
virtual bool progressItem(const size_t &handle)
CModelEntity ** beginDependent()
Definition: CState.cpp:210
virtual std::set< const CCopasiObject * > getDeletedObjects() const
Definition: CReaction.cpp:959
void setTime(const C_FLOAT64 &time)
Definition: CState.cpp:326
const CCopasiVector< CMoiety > & getMoieties() const
Definition: CModel.cpp:1163
void calculateRootDerivatives(CVector< C_FLOAT64 > &rootDerivatives)
Definition: CModel.cpp:4712
C_FLOAT64 * mpIValue
Definition: CModelValue.h:326
size_t mNumMetabolitesUnused
Definition: CModel.h:1420
std::list< const CCopasiObject * > sortObjectsByDependency(ForwardAccessIterator begin, ForwardAccessIterator end, const CCopasiObject::DataObjectSet &context)
size_t getTotSteps() const
Definition: CModel.cpp:1136
VolumeUnit mVolumeUnit
Definition: CModel.h:1294
VolumeUnit
Definition: CModel.h:56
bool hasReversibleReaction() const
Definition: CModel.cpp:3892
CMathModel * mpMathModel
Definition: CModel.h:1532
ModelType
Definition: CModel.h:112
virtual DataObjectSet getDeletedObjects() const
Definition: CMetab.cpp:806
virtual const std::string & getKey() const
void setInitialTime(const C_FLOAT64 &time)
Definition: CModel.cpp:1181
CCopasiVector< CMetab > mMetabolites
Definition: CModel.h:1329
bool removeModelValue(const CModelValue *pModelValue, const bool &recursive=true)
Definition: CModel.cpp:2903
std::string mKey
Definition: CAnnotation.h:119
virtual bool compile()
bool appendDependentReactions(std::set< const CCopasiObject * > candidates, std::set< const CCopasiObject * > &dependents) const
Definition: CModel.cpp:2472
size_t getNumDataRows() const
std::set< const CCopasiObject * > mSimulatedUpToDateObjects
Definition: CModel.h:1285
const C_FLOAT64 & getInitialValue() const
void compile()
Definition: CReaction.cpp:583
LengthUnit
Definition: CModel.h:76
const CCopasiVector< CChemEqElement > & getProducts() const
Definition: CChemEq.cpp:63
const CFunction * getFunction() const
Definition: CReaction.cpp:252
#define DESTRUCTOR_TRACE
Definition: copasi.h:206
std::vector< Refresh * > mApplyInitialValuesRefreshes
Definition: CModel.h:1501
std::string getLengthUnitsDisplayString() const
Definition: CModel.cpp:4563
iterator end()
void addDirectDependency(const CCopasiObject *pObject)
CStateTemplate mStateTemplate
Definition: CModel.h:1279
CArrayAnnotation * mpLinkMatrixAnnotation
Definition: CModel.h:1451
std::vector< CType * >::const_iterator const_iterator
Definition: CCopasiVector.h:57
bool getInitialUpdateSequence(const CMath::SimulationContextFlag &context, const CCopasiObject::DataObjectSet &changedObjects, const CCopasiObject::DataObjectSet &requestedObjects, CCopasiObject::DataUpdateSequence &updateSequence) const
Definition: CModel.cpp:1792
const C_FLOAT64 & getQuantity2NumberFactor() const
Definition: CModel.cpp:2354
bool appendDependentModelObjects(const std::set< const CCopasiObject * > &candidates, std::set< const CCopasiObject * > &dependentReactions, std::set< const CCopasiObject * > &dependentMetabolites, std::set< const CCopasiObject * > &dependentCompartments, std::set< const CCopasiObject * > &dependentModelValues, std::set< const CCopasiObject * > &dependentEvents) const
Definition: CModel.cpp:2364
CVector< C_FLOAT64 > mParticleFluxes
Definition: CModel.h:1349
C_FLOAT64 mNumber2QuantityFactor
Definition: CModel.h:1473
const CFunctionParameters & getFunctionParameters() const
Definition: CReaction.cpp:576
const C_FLOAT64 & getNumber2QuantityFactor() const
Definition: CModel.cpp:2357
virtual bool add(const CType &src)
bool doRowPivot(CMatrix< C_FLOAT64 > &matrix) const
size_t getNumVariableMetabs() const
Definition: CModel.cpp:1121
const std::vector< Refresh * > & getListOfNonSimulatedRefreshes() const
Definition: CModel.cpp:4142
size_t getNumDependentReactionMetabs() const
Definition: CModel.cpp:1133
size_t findMetabByName(const std::string &Target) const
Definition: CModel.cpp:1198
std::string getVolumeUnitsDisplayString() const
Definition: CModel.cpp:4547
virtual const std::string & getKey() const
bool mReorderNeeded
Definition: CModel.h:1512
CModelEntity ** endFixed()
Definition: CState.cpp:213
void setState(const CState &state)
Definition: CModel.cpp:1785
bool addModifier(const std::string &metabKey, const C_FLOAT64 &multiplicity=1.0)
Definition: CReaction.cpp:240
void removeDependentModelObjects(const std::set< const CCopasiObject * > &deletedObjects)
Definition: CModel.cpp:2859
bool compile()
Definition: CModel.cpp:381
bool buildDependencyGraphs()
Definition: CModel.cpp:547
const bool & isInitialConcentrationChangeAllowed() const
Definition: CMetab.cpp:664
bool applyRowPivot(CVectorCore< CType > &vector) const
Definition: CLinkMatrix.h:115
const std::string & getNotes() const
size_t addItem(const std::string &name, const std::string &value, const std::string *pEndValue=NULL)
void determineIsAutonomous()
Definition: CModel.cpp:3959
bool setFunction(const std::string &functionName)
Definition: CReaction.cpp:255
virtual void setStatus(const CModelEntity::Status &status)
virtual void cleanup()
size_t mNumMetabolitesAssignment
Definition: CModel.h:1435
void clearMoieties()
Definition: CModel.cpp:1257
void processRoots(const C_FLOAT64 &time, const bool &equality, const bool &correct, const CVector< C_INT > &roots)
Definition: CModel.cpp:4692
CReaction * createReaction(const std::string &name)
Definition: CModel.cpp:2760
void refreshActiveParameterSet()
Definition: CModel.cpp:1105
std::pair< CFunction *, CFunction * > splitFunction(const CEvaluationNode *node, const std::string &name1, const std::string &name2) const
Definition: CFunction.cpp:445
const bool & isAutonomous() const
Definition: CModel.cpp:3956
bool setTriggerExpression(const std::string &expression)
Definition: CEvent.cpp:474
CCopasiVector< CMetab > mMetabolitesX
Definition: CModel.h:1334
virtual bool add(const CType &src)
std::string getInitialExpression() const
const Value & getValue() const
const C_FLOAT64 & getInitialTime() const
Definition: CModel.cpp:1184
std::vector< Refresh * > mInitialRefreshes
Definition: CModel.h:1490
bool appendDependentCompartments(std::set< const CCopasiObject * > candidates, std::set< const CCopasiObject * > &dependents) const
Definition: CModel.cpp:2538
void calculateElasticityMatrix(const C_FLOAT64 &factor, const C_FLOAT64 &resolution)
Definition: CModel.cpp:1968
CCopasiVectorN< CCopasiTask > * getTaskList()
void setUsed(const bool &used)
void initializeMetabolites()
Definition: CModel.cpp:942
bool getUpdateSequence(const CMath::SimulationContextFlag &context, const CObjectInterface::ObjectSet &changedObjects, const CObjectInterface::ObjectSet &requestedObjects, CObjectInterface::UpdateSequence &updateSequence)
void appendElementReferences(std::set< const CCopasiObject * > &objects) const
bool mCompileIsNecessary
Definition: CModel.h:1479
void exportDOTFormat(std::ostream &os, const std::string &name) const
std::string getLengthUnitName() const
Definition: CModel.cpp:2215
QuantityUnit mQuantityUnit
Definition: CModel.h:1314
#define MCReadConfig
CModel::LengthUnit getLengthUnitEnum() const
Definition: CModel.cpp:2220
std::vector< Refresh * > mSimulatedRefreshes
Definition: CModel.h:1495
size_t getNumIndependent() const
Definition: CState.cpp:342
bool handleUnusedMetabolites()
Definition: CModel.cpp:723
void setCopasiVector(size_t d, const CCopasiContainer *v)
bool appendDependentMetabolites(std::set< const CCopasiObject * > candidates, std::set< const CCopasiObject * > &dependents) const
Definition: CModel.cpp:2501
CState mInitialState
Definition: CModel.h:1272
virtual bool finishItem(const size_t &handle)
bool getTransientUpdateSequence(const CMath::SimulationContextFlag &context, const CCopasiObject::DataObjectSet &changedObjects, const CCopasiObject::DataObjectSet &requestedObjects, CCopasiObject::DataUpdateSequence &updateSequence) const
Definition: CModel.cpp:1804
virtual void refreshInitialValue()
Definition: CMetab.cpp:264
CModelEntity ** beginIndependent()
Definition: CState.cpp:208
std::vector< const CEvaluationTree * > getTreesWithDiscontinuities() const
Definition: CModel.cpp:4059
bool buildApplyInitialValuesSequence()
Definition: CModel.cpp:1627
CCopasiVectorN< CModelParameterSet > mParameterSets
Definition: CModel.h:1364
const CCopasiVectorN< CModelParameterSet > & getModelParameterSets() const
Definition: CModel.cpp:1066
CCopasiParameter * getParameter(const std::string &name)
bool setPriorityExpression(const std::string &expression)
Definition: CEvent.cpp:611
TimeUnit
Definition: CModel.h:86
int dgemm_(char *transa, char *transb, integer *m, integer *n, integer *k, doublereal *alpha, doublereal *a, integer *lda, doublereal *b, integer *ldb, doublereal *beta, doublereal *c, integer *ldc)
void buildStoi()
Definition: CModel.cpp:642
CCopasiObjectReference< C_FLOAT64 > * mpValueReference
Definition: CModelValue.h:357
std::string add(const std::string &prefix, CCopasiObject *pObject)
const CCopasiVector< CChemEqElement > & getSubstrates() const
Definition: CChemEq.cpp:60
const std::string & getKey() const
Definition: CModel.cpp:1142
size_t findMoiety(const std::string &Target) const
Definition: CModel.cpp:1216
const bool & isUsed() const
CMathDependencyGraph mTransientDependencies
Definition: CModel.h:1288
static CFunctionDB * getFunctionList()
virtual void resize(size_t rows, size_t cols, const bool &copy=false)
Definition: CMatrix.h:151
std::string getDelayExpression() const
Definition: CEvent.cpp:591
virtual std::set< const CCopasiObject * > getDeletedObjects() const
size_t getNumIndependentReactionMetabs() const
Definition: CModel.cpp:1130
const C_FLOAT64 & getTime() const
Definition: CModel.cpp:1190
CType toEnum(const char *attribute, const char **enumNames, const CType &enumDefault)
Definition: utility.h:107
void setCompileFlag(bool flag=true)
Definition: CModel.cpp:607
bool processQueue(const C_FLOAT64 &time, const bool &equality, CProcessQueue::resolveSimultaneousAssignments pResolveSimultaneousAssignments)
Definition: CModel.cpp:4685
virtual void resize(const size_t &newSize)
TimeUnit mTimeUnit
Definition: CModel.h:1309
#define MCReaction
bool mBuildInitialSequence
Definition: CModel.h:1523
bool dependsOn(DataObjectSet candidates, const DataObjectSet &context=DataObjectSet()) const
CProcessReport * mpCompileHandler
Definition: CModel.h:1485
size_t size() const
Definition: CVector.h:100
long int flag
Definition: f2c.h:52
static const char * AreaUnitNames[]
Definition: CModel.h:71
void setTime(const C_FLOAT64 &time)
Definition: CModel.cpp:1187
virtual void remove(const size_t &index)
virtual const DataObjectSet & getDirectDependencies(const DataObjectSet &context=DataObjectSet()) const
void setMode(size_t d, Mode m)
void calculateDerivativesX(C_FLOAT64 *derivativesX)
Definition: CModel.cpp:1941
std::string getAreaUnitName() const
Definition: CModel.cpp:2193
std::string getConcentrationUnitsDisplayString() const
Definition: CModel.cpp:4587
bool getUpdateSequence(CMathDependencyGraph &dependencyGraph, const CMath::SimulationContextFlag &context, const CCopasiObject::DataObjectSet &changedObjects, const CCopasiObject::DataObjectSet &requestedObjects, CCopasiObject::DataUpdateSequence &updateSequence) const
Definition: CModel.cpp:1816
CCopasiVectorNS< CCompartment > & getCompartments()
Definition: CModel.cpp:1145
void setDimensionDescription(size_t d, const std::string &s)
iterator addObject(const CObjectInterface *pObject)
const CExpression * getTriggerExpressionPtr() const
Definition: CEvent.cpp:534
const CCopasiTask::Type & getExperimentType() const
std::string getQuantityUnitsDisplayString() const
Definition: CModel.cpp:4647
void updateMatrixAnnotations()
Definition: CModel.cpp:863
const C_FLOAT64 & getValue() const
void buildLinkZero()
Definition: CModel.cpp:3946
virtual size_t getIndex(const CCopasiObject *pObject) const
CModel::TimeUnit getTimeUnitEnum() const
Definition: CModel.cpp:2243
virtual std::string getChildObjectUnits(const CCopasiObject *pObject) const
Definition: CModel.cpp:246
const CVector< C_FLOAT64 > & getTimeData() const
iterator(* resolveSimultaneousAssignments)(const std::multimap< CKey, CAction > &, const C_FLOAT64 &, const bool &, const size_t &)
static CKeyFactory * getKeyFactory()
Header file of class CArrayAnnotation.
void add(C_FLOAT64 value, CMetab *metabolite)
Definition: CMoiety.cpp:96
bool setDelayExpression(const std::string &expression)
Definition: CEvent.cpp:544
CMathDependencyGraph mInitialDependencies
Definition: CModel.h:1287
size_t mNumMetabolitesReaction
Definition: CModel.h:1430
std::string getQuantityUnitName() const
Definition: CModel.cpp:2321
#define C_FLOAT64
Definition: copasi.h:92
const CState & getInitialState() const
Definition: CModel.cpp:1768
bool compileIfNecessary(CProcessReport *pProcessReport)
Definition: CModel.cpp:612
std::string mActiveParameterSetKey
Definition: CModel.h:1369
const CCopasiParameterGroup & getParameters() const
Definition: CReaction.cpp:333
CMatrix< C_FLOAT64 > mStoiInternal
Definition: CModel.h:1379
CFunction * findFunction(CCopasiVectorN< CFunction > &db, const CFunction *func)
size_t getNumAssignmentMetabs() const
Definition: CModel.cpp:1127
CType * array()
Definition: CVector.h:139
void buildMoieties()
Definition: CModel.cpp:905
const CStateTemplate & getStateTemplate() const
Definition: CModel.cpp:1172
void setParameterMapping(const size_t &index, const std::string &key)
Definition: CReaction.cpp:339
size_t getCompartmentNumber() const
Definition: CReaction.cpp:873
size_t getNumModelValues() const
Definition: CModel.cpp:1139
CCopasiVectorNS< CReaction > mSteps
Definition: CModel.h:1339
ModelType mType
Definition: CModel.h:1319
bool compile(const std::vector< CCopasiContainer * > listOfContainer=CCopasiContainer::EmptyList)
std::vector< CType * >::iterator iterator
Definition: CCopasiVector.h:56
~CModel()
Definition: CModel.cpp:227
const ModelType & getModelType() const
Definition: CModel.cpp:2339
AreaUnit mAreaUnit
Definition: CModel.h:1299
The class for handling a chemical kinetic function.
Definition: CFunction.h:29
CCopasiVectorNS< CReaction > & getReactions()
Definition: CModel.cpp:1039
static std::vector< Refresh * > buildUpdateSequence(const DataObjectSet &objects, const DataObjectSet &uptoDateObjects, const DataObjectSet &context=DataObjectSet())
bool setQuantityUnit(const std::string &name)
Definition: CModel.cpp:2250
bool appendDependentEvents(std::set< const CCopasiObject * > candidates, std::set< const CCopasiObject * > &dependents) const
Definition: CModel.cpp:2593
virtual void clear()
const C_FLOAT64 & getProcessQueueExecutionTime() const
Definition: CModel.cpp:4702
bool updateInitialValues()
Definition: CModel.cpp:1461
void calculateDerivatives(C_FLOAT64 *derivatives)
Definition: CModel.cpp:1903
const C_FLOAT64 & getConcentration() const
Definition: CMetab.cpp:218
virtual void setInitialValue(const C_FLOAT64 &initialValue)
virtual bool refreshFromModel(const bool &modifyExistence)
const CLinkMatrix & getL0() const
Definition: CModel.cpp:1169
bool compileEvents()
Definition: CModel.cpp:4116
std::string getAreaUnitsDisplayString() const
Definition: CModel.cpp:4555
CVector< C_FLOAT64 > initializeAtolVector(const C_FLOAT64 &baseTolerance, const bool &reducedModel) const
Definition: CModel.cpp:4368
Definition: CModel.h:50
C_FLOAT64 mAvogadro
Definition: CModel.h:1461
bool convert2NonReversible()
Definition: CModel.cpp:3205
C_INT32 getVariable(const std::string &name, const std::string &type, void *pout, CReadConfig::Mode mode=CReadConfig::NEXT)
Definition: CReadConfig.cpp:81
const std::vector< Refresh * > & getListOfInitialRefreshes() const
Definition: CModel.cpp:4133
virtual const CObjectInterface * getObject(const CCopasiObjectName &cn) const
const CState & getState() const
Definition: CModel.cpp:1771
static CLocaleString fromUtf8(const std::string &utf8)
void setAvogadro(const C_FLOAT64 &avogadro)
Definition: CModel.cpp:2342
const CModelEntity::Status & getStatus() const
std::string getVolumeUnitName() const
Definition: CModel.cpp:2170
void calculateRootDerivatives(CVector< C_FLOAT64 > &rootDerivatives)
Definition: CMathModel.cpp:417
const CCopasiVector< CMetab > & getMetabolitesX() const
Definition: CModel.cpp:1057
CModel::VolumeUnit getVolumeUnitEnum() const
Definition: CModel.cpp:2175
const CMatrix< C_FLOAT64 > & getStoi() const
Definition: CModel.cpp:1160
const CVector< size_t > & getUserOrder() const
Definition: CState.cpp:201
static const char * TimeUnitNames[]
Definition: CModel.h:91
const std::vector< std::vector< std::string > > & getParameterMappings() const
Definition: CReaction.h:285
CConcentrationReference * getInitialConcentrationReference() const
Definition: CMetab.cpp:861
virtual size_t numCols() const
Definition: CMatrix.h:144
std::vector< Refresh * > buildInitialRefreshSequence(std::set< const CCopasiObject * > &changedObjects)
Definition: CModel.cpp:4164
const unsigned C_INT32 & getNumColumns() const
const bool & getContinueSimultaneousEvents() const
C_FLOAT64 * beginDependent()
Definition: CState.cpp:330
bool setExpression(const std::string &expression)
bool setObjectName(const std::string &name)