COPASI API  4.16.103
CMathObject.cpp
Go to the documentation of this file.
1 // Copyright (C) 2011 - 2013 by Pedro Mendes, Virginia Tech Intellectual
2 // Properties, Inc., University of Heidelberg, and The University
3 // of Manchester.
4 // All rights reserved.
5 
6 #include <cmath>
7 
8 #include "copasi.h"
9 
10 #include "CMathObject.h"
11 #include "CMathExpression.h"
12 #include "CMathContainer.h"
13 
14 #include "model/CMetab.h"
15 #include "model/CCompartment.h"
16 #include "model/CModel.h"
17 #include "function/CExpression.h"
18 
19 // static
20 C_FLOAT64 CMathObject::InvalidValue = std::numeric_limits< C_FLOAT64 >::quiet_NaN();
21 
22 // static
24  C_FLOAT64 *& pValue,
25  const CMath::ValueType & valueType,
26  const CMath::EntityType & entityType,
27  const CMath::SimulationType & simulationType,
28  const bool & isIntensiveProperty,
29  const bool & isInitialValue,
30  const CCopasiObject * pDataObject)
31 {
32  pObject->mpValue = pValue;
33  pObject->mValueType = valueType;
34  pObject->mEntityType = entityType;
35  pObject->mSimulationType = simulationType;
38  pObject->mpDataObject = pDataObject;
39 
40  pObject++;
41  pValue++;
42 }
43 
46  mpExpression(NULL),
47  mpValue(&InvalidValue),
48  mPrerequisites(),
49  mValueType(CMath::ValueTypeUndefined),
50  mEntityType(CMath::EntityTypeUndefined),
51  mSimulationType(CMath::SimulationTypeUndefined),
52  mIsIntensiveProperty(false),
53  mIsInitialValue(false),
54  mpIntensiveProperty(NULL),
55  mpDataObject(NULL)
56 {}
57 
58 // virtual
60 {
62 }
63 
64 void CMathObject::copy(const CMathObject & src, CMathContainer & container, const size_t & valueOffset, const size_t & objectOffset)
65 {
66  mpValue = (C_FLOAT64 *)((size_t) src.mpValue + valueOffset);
67  mValueType = src.mValueType;
73 
74  if (src.mpIntensiveProperty != NULL)
75  {
76  mpIntensiveProperty = (CMathObject *)((size_t) src.mpIntensiveProperty + objectOffset);
77  }
78  else
79  {
80  mpIntensiveProperty = NULL;
81  }
82 
83  if (src.mpExpression != NULL)
84  {
85  mpExpression = CMathExpression::copy(*src.mpExpression, container, valueOffset, objectOffset);
86  }
87  else
88  {
89  mpExpression = NULL;
90  }
91 
92  ObjectSet::const_iterator it = src.getPrerequisites().begin();
93  ObjectSet::const_iterator end = src.getPrerequisites().end();
94 
95  for (; it != end; ++it)
96  {
97  mPrerequisites.insert((CMathObject *)((size_t) *it + objectOffset));
98  }
99 
100  if (mpExpression != NULL &&
101  mPrerequisites.empty())
102  {
103  calculate();
104  }
105 }
106 
107 // virtual
109 {
110  if (mpDataObject == NULL)
111  return CCopasiObjectName("");
112 
113  return mpDataObject->getCN();
114 }
115 
116 //virtual
118 {
119  if (mpDataObject == NULL)
120  return NULL;
121 
122  const CObjectInterface * pObject = mpDataObject->getObject(cn);
123 
124  if (pObject != mpDataObject)
125  return pObject;
126 
127  return this;
128 }
129 
130 // virtual
132 {
133  return mPrerequisites;
134 }
135 
136 // virtual
138  const CMath::SimulationContextFlag & context,
139  const CObjectInterface::ObjectSet & changedObjects) const
140 {
141  // This method should only be called for objects which are prerequisites.
142  // We check for this only in debug mode.
143  assert(mPrerequisites.find(pObject) != mPrerequisites.end());
144 
145  switch (mEntityType)
146  {
147  case CMath::Moiety:
148 
149  if ((context & CMath::UpdateMoieties) &&
151  {
152  return true;
153  }
154 
155  if ((context & CMath::UseMoieties) &&
157  {
158  return true;
159  }
160 
161  return false;
162  break;
163 
164  case CMath::Species:
165 
166  // For species we need to account for the duality of the intensive and extensive value
167  if (mValueType != CMath::Value)
168  return true;
169 
170  if ((context & CMath::UseMoieties) &&
173  {
174  if (mpIntensiveProperty != pObject)
175  {
176  return true;
177  }
178 
179  return false;
180  }
181 
182  // If the value is in the context, it does not depend on the object.
183  if (changedObjects.find(this) != changedObjects.end())
184  return false;
185 
187  {
188  // Densities which are not in the context have to be recalculated.
189  return true;
190  }
191  else
192  {
193  // Amount which are determine by assignment need to be recalculated.
195  return true;
196 
197  // If the concentration was changed in the context we need to recalculate.
198  if (changedObjects.find(mpIntensiveProperty) != changedObjects.end())
199  return true;
200 
201  return false;
202  }
203 
204  break;
205 
206  case CMath::Event:
207 
208  if ((context & CMath::EventHandling) &&
210  {
211  switch ((int) mpExpression->getRoot()->getType())
212  {
214  {
215  const CMathObject * pMathObject = dynamic_cast< const CMathObject * >(pObject);
216 
217  if (pMathObject != NULL &&
218  pMathObject->mValueType == CMath::EventTrigger)
219  {
220  return false;
221  }
222 
223  return true;
224  }
225  break;
226 
228  return false;
229  break;
230 
232  return false;
233  break;
234 
235  default:
236  return true;
237  }
238  }
239 
240  return true;
241  break;
242 
243  default:
244  return true;
245  }
246 
247  // This should never be reached.
248  return true;
249 }
250 
251 // virtual
252 void CMathObject::print(std::ostream * ostream) const
253 {
254  (*ostream) << *mpValue;
255 }
256 
257 // virtual
259 {
260  return mpValue;
261 }
262 
264 {
265  // This method should only be called if there is something to calculate, i.e.
266  // mpExpression != NULL
267  assert(mpExpression != NULL);
268 
269  *mpValue = mpExpression->value();
270 
271 #ifdef COPASI_DEBUG
272 
273  // Check for NaN
274  if (isnan(*mpValue))
275  {
276  std::cout << "NaN Value for: " << getCN() << std::endl;
277  }
278 
279 #endif // COPASI_DEBUG
280 
281  // For an extensive transient value of a dependent species we have 2
282  // possible assignments depending on the context.
283  // 1) Conversion from the intensive property
284  // 2) Dependent mass off a moiety
285  //
286  // The solution is that the moiety automatically updates the value in conjunction
287  // with the dependency graph omitting the value in the update sequence if the context
288  // is CMath::UseMoieties.
289 }
290 
292 {
293  return mpDataObject;
294 }
295 
297 {
298  return mValueType;
299 }
300 
302 {
303  return mEntityType;
304 }
305 
307 {
308  return mSimulationType;
309 }
310 
312 {
313  mSimulationType = simulationType;
314 }
315 
317 {
318  return mIsIntensiveProperty;
319 }
320 
321 const bool & CMathObject::isInitialValue() const
322 {
323  return mIsInitialValue;
324 }
325 
326 bool CMathObject::setExpression(const std::string & infix,
327  const bool & isBoolean,
328  CMathContainer & container)
329 {
330  bool success = true;
331  CExpression Expression;
332  Expression.setIsBoolean(isBoolean);
333 
334  success &= Expression.setInfix(infix);
335  std::vector< CCopasiContainer * > ListOfContainer;
336  ListOfContainer.push_back(&container);
337  success &= Expression.compile(ListOfContainer);
338  success &= setExpression(Expression, container);
339 
340  return success;
341 }
342 
343 bool CMathObject::setExpression(const CExpression & expression,
344  CMathContainer & container)
345 {
346  bool success = true;
347 
348  success &= createConvertedExpression(&expression, container);
349 
350  return success;
351 }
352 
354 {
355  bool success = true;
356 
358  mPrerequisites.clear();
359 
360  mpExpression = pMathExpression;
361 
362  if (mpExpression != NULL)
363  {
364  success &= mpExpression->compile();
366  }
367  else
368  {
369  success = false;
370  }
371 
372  return success;
373 }
374 
376 {
377  return mpExpression;
378 }
379 
381 {
382 
383  bool success = true;
384 
385  switch (mValueType)
386  {
388  success = false;
389  break;
390 
391  case CMath::Value:
392 
393  if (mIsInitialValue)
394  {
395  success = compileInitialValue(container);
396  }
397  else
398  {
399  success = compileValue(container);
400  }
401 
402  break;
403 
404  case CMath::Rate:
405  success = compileRate(container);
406  break;
407 
408  case CMath::ParticleFlux:
409  success = compileParticleFlux(container);
410  break;
411 
412  case CMath::Flux:
413  success = compileFlux(container);
414  break;
415 
416  case CMath::Propensity:
417  success = compilePropensity(container);
418  break;
419 
420  case CMath::TotalMass:
421  success = compileTotalMass(container);
422  break;
423 
425  success = compileDependentMass(container);
426  break;
427 
429  case CMath::EventDelay:
432  case CMath::EventTrigger:
433  case CMath::EventRoot:
435  // These objects are compiled through the event compile,
436  // which is executed after the object compile. It is therefore
437  // correct to leave the object in its default state.
438  break;
439  }
440 
441  return success;
442 }
443 
445 {
446  bool success = true;
447 
448  // The default value is NaN
450 
451  // Initial values are taken from the data model
452  if (mpDataObject != NULL)
453  {
455  }
456 
457  // Reset the prerequisites
458  mPrerequisites.clear();
459 
460  const CModelEntity * pEntity = dynamic_cast< const CModelEntity * >(mpDataObject->getObjectParent());
461 
463  {
464  // Only species have intensive properties (concentration and concentration rate).
465  const CMetab * pSpecies = static_cast< const CMetab * >(pEntity);
466 
467  switch (mSimulationType)
468  {
469  case CMath::EventTarget:
470  case CMath::Fixed:
471  case CMath::ODE:
472  case CMath::Independent:
473  case CMath::Dependent:
474  case CMath::Conversion:
475  success &= createIntensiveValueExpression(pSpecies, container);
476  break;
477 
478  case CMath::Assignment:
479  // Extensive Property * Conversion / Compartment Size
480  success &= createConvertedExpression(pSpecies->getInitialExpressionPtr(), container);
481 
482  break;
483 
484  case CMath::Time:
486  success = false;
487  break;
488  }
489  }
490  else
491  {
492  switch (mSimulationType)
493  {
494  case CMath::Fixed:
495  break;
496 
497  case CMath::Assignment:
498  success &= createConvertedExpression(pEntity->getInitialExpressionPtr(), container);
499  break;
500 
501  case CMath::Conversion:
502  {
503  const CMetab * pSpecies = static_cast< const CMetab * >(pEntity);
505  success &= createExtensiveValueExpression(pSpecies, container);
506  }
507  break;
508 
510  case CMath::EventTarget:
511  case CMath::Time:
512  case CMath::ODE:
513  case CMath::Independent:
514  case CMath::Dependent:
515  success = false;
516  break;
517  }
518  }
519 
520  return success;
521 }
522 
524 {
525  bool success = true;
526 
527  // The default value is NaN
529 
530  // Reset the prerequisites
531  mPrerequisites.clear();
532 
534  {
535  // Only species have intensive properties.
536  const CMetab * pSpecies = static_cast< const CMetab * >(mpDataObject->getObjectParent());
537 
538  switch (mSimulationType)
539  {
540  case CMath::Assignment:
541  success &= createConvertedExpression(pSpecies->getExpressionPtr(), container);
542  break;
543 
544  case CMath::EventTarget:
545  case CMath::Conversion:
546  success &= createIntensiveValueExpression(pSpecies, container);
547  break;
548 
550  case CMath::Fixed:
551  case CMath::Time:
552  case CMath::ODE:
553  case CMath::Independent:
554  case CMath::Dependent:
555  success = false;
556  break;
557  }
558  }
559  else
560  {
561  // Species need an additional conversion since the event targets the
562  // intensive property.
563  const CMetab * pSpecies = NULL;
564 
566  {
567  pSpecies = static_cast< const CMetab * >(mpDataObject->getObjectParent());
569  success &= createExtensiveValueExpression(pSpecies, container);
570  }
571 
572  switch (mSimulationType)
573  {
574  case CMath::Fixed:
575  case CMath::EventTarget:
576  case CMath::Time:
577  case CMath::ODE:
578  case CMath::Independent:
579  case CMath::Conversion:
580  break;
581 
582  case CMath::Dependent:
583  {
584  // We need to add the dependent number of the moiety as a possible
585  // prerequisite.
586  const CMoiety * pMoiety = pSpecies->getMoiety();
587  const CMathObject * pDependentNumber =
588  container.getMathObject(pMoiety->getDependentNumberReference());
589 
590  mPrerequisites.insert(pDependentNumber);
591  }
592  break;
593 
594  case CMath::Assignment:
595  {
596  const CModelEntity * pEntity = static_cast< const CModelEntity * >(mpDataObject->getObjectParent());
597  success &= createConvertedExpression(pEntity->getExpressionPtr(), container);
598  }
599  break;
600 
602  success = false;
603  break;
604  }
605  }
606 
607  return success;
608 }
609 
611 {
612  bool success = true;
613 
614  // The default value is NaN
616 
617  // Reset the prerequisites
618  mPrerequisites.clear();
619 
621  {
622  // Only species have intensive properties.
623  const CMetab * pSpecies = static_cast< const CMetab * >(mpDataObject->getObjectParent());
624 
625  switch (mSimulationType)
626  {
627  case CMath::Assignment:
628  success &= createIntensiveRateExpression(pSpecies, container);
629  break;
630 
632  case CMath::Fixed:
633  case CMath::EventTarget:
634  case CMath::Time:
635  case CMath::ODE:
636  case CMath::Independent:
637  case CMath::Dependent:
638  case CMath::Conversion:
639  success = false;
640  break;
641  }
642  }
643  else
644  {
645 
646  switch (mSimulationType)
647  {
648  case CMath::Fixed:
649  *mpValue = 0;
650  break;
651 
652  case CMath::Time:
653  *mpValue = 1;
654  break;
655 
656  case CMath::ODE:
657 
659  {
660  const CMetab * pSpecies = static_cast< const CMetab * >(mpDataObject->getObjectParent());
661  success &= createExtensiveODERateExpression(pSpecies, container);
662  }
663  else
664  {
665  const CModelEntity * pEntity = static_cast< const CModelEntity * >(mpDataObject->getObjectParent());
666  success &= createConvertedExpression(pEntity->getExpressionPtr(), container);
667  }
668 
669  break;
670 
671  case CMath::Independent:
672  case CMath::Dependent:
673  {
674  const CMetab * pSpecies = static_cast< const CMetab * >(mpDataObject->getObjectParent());
675  success &= createExtensiveReactionRateExpression(pSpecies, container);
676  }
677  break;
678 
679  case CMath::Assignment:
680  // TODO When we have symbolic differentiation we can deal with this.
681  break;
682 
684  case CMath::EventTarget:
685  case CMath::Conversion:
686  success = false;
687  break;
688  }
689  }
690 
691  return success;
692 }
693 
695 {
696  bool success = true;
697 
698  // The default value is NaN
700 
701  // Reset the prerequisites
702  mPrerequisites.clear();
703 
704  const CReaction * pReaction = static_cast< const CReaction * >(mpDataObject->getObjectParent());
705 
706  // We need to check whether this reaction is a single compartment reaction and scale
707  // it if true.
708  // mParticleFlux = *mUnitScalingFactor * mFlux;
709  // mUnitScalingFactor = & pModel->getQuantity2NumberFactor();
710 
711  std::ostringstream Infix;
712  Infix.imbue(std::locale::classic());
713  Infix.precision(16);
714 
715  Infix << container.getModel().getQuantity2NumberFactor();
716  Infix << "*<";
717  Infix << pReaction->getFluxReference()->getCN();
718  Infix << ">";
719 
720  CExpression E("ParticleExpression", &container);
721 
722  success &= E.setInfix(Infix.str());
723 
724  mpExpression = new CMathExpression(E, container, !mIsInitialValue);
726 
727  return success;
728 }
729 
731 {
732  bool success = true;
733 
734  // The default value is NaN
736 
737  // Reset the prerequisites
738  mPrerequisites.clear();
739 
740  const CReaction * pReaction = static_cast< const CReaction * >(mpDataObject->getObjectParent());
741 
742  // We need to check whether this reaction is a single compartment reaction and scale it if true.
743  // mFlux = *mScalingFactor * mpFunction->calcValue(mMap.getPointers());
744  // mScalingFactor = compartment volume or 1
745 
746  mpExpression = new CMathExpression(*pReaction->getFunction(),
747  pReaction->getCallParameters(),
748  container,
749  !mIsInitialValue);
750 
751  std::set< const CCompartment * > Compartments = pReaction->getChemEq().getCompartments();
752 
753  if (Compartments.size() == 1)
754  {
755  CExpression Tmp(mpExpression->getObjectName(), &container);
756 
757  std::string Infix = "<" + (*Compartments.begin())->getValueReference()->getCN() + ">*(" + mpExpression->getInfix() + ")";
758  success &= Tmp.setInfix(Infix);
759  success &= Tmp.compile();
760 
762  mpExpression = new CMathExpression(Tmp, container, false);
763  }
764 
766 
767  return success;
768 }
769 
771 {
772  bool success = true;
773 
774  // The default value is NaN
776 
777  // Reset the prerequisites
778  mPrerequisites.clear();
779 
780  const CReaction * pReaction = static_cast< const CReaction * >(mpDataObject->getObjectParent());
781 
782  std::ostringstream Infix;
783  Infix.imbue(std::locale::classic());
784  Infix.precision(16);
785 
786  // Propensity for reversible reactions must be NaN
787  if (pReaction->isReversible())
788  {
789  Infix << "NAN";
790  }
791  else if (container.getModel().getModelType() == CModel::deterministic)
792  {
793  Infix << "<" << pReaction->getParticleFluxReference()->getCN() << ">";
794 
795  std::ostringstream Divisor;
796  Divisor.imbue(std::locale::classic());
797  Divisor.precision(16);
798 
799  const CCopasiVector<CChemEqElement> & Substrates = pReaction->getChemEq().getSubstrates();
800  CCopasiVector< CChemEqElement >::const_iterator itSubstrate = Substrates.begin();
801  CCopasiVector< CChemEqElement >::const_iterator endSubstrate = Substrates.end();
802  bool first = true;
803 
804  for (; itSubstrate != endSubstrate; ++itSubstrate)
805  {
806  const std::string NumberCN = (*itSubstrate)->getMetabolite()->getValueReference()->getCN();
807  C_FLOAT64 Multiplicity = (*itSubstrate)->getMultiplicity();
808 
809  Multiplicity -= 1.0; // Nothing to correct if the multiplicity is 1.
810 
811  if (Multiplicity > 2.0 - 100.0 * std::numeric_limits< C_FLOAT64 >::epsilon())
812  {
813  if (!first)
814  {
815  Divisor << "*";
816  }
817 
818  first = false;
819  Divisor << "<" << NumberCN << ">^" << Multiplicity;
820  }
821  else if (Multiplicity > 1.0 - 100.0 * std::numeric_limits< C_FLOAT64 >::epsilon())
822  {
823  if (!first)
824  {
825  Divisor << "*";
826  }
827 
828  first = false;
829  Divisor << "<" << NumberCN << ">";
830  }
831 
832  while (Multiplicity > 1.0 - 100.0 * std::numeric_limits< C_FLOAT64 >::epsilon())
833  {
834  Infix << "*(<" << NumberCN << ">-" << Multiplicity << ")";
835  Multiplicity -= 1.0;
836  }
837  }
838 
839  if (Divisor.str() != "")
840  {
841  Infix << "/(" << Divisor.str() << ")";
842  }
843  }
844  else
845  {
846  // Propensity is the same as the flux.
847  Infix << "<" << pReaction->getParticleFluxReference()->getCN() << ">";
848  }
849 
850  CExpression E("PropensityExpression", &container);
851 
852  success &= E.setInfix(Infix.str());
853 
854  mpExpression = new CMathExpression(E, container, !mIsInitialValue);
856 
857  return success;
858 }
859 
861 {
862  bool success = true;
863 
864  // The default value is NaN
866 
867  // Reset the prerequisites
868  mPrerequisites.clear();
869 
870  const CMoiety * pMoiety = static_cast< const CMoiety *>(mpDataObject->getObjectParent());
871 
872  std::ostringstream Infix;
873  Infix.imbue(std::locale::classic());
874  Infix.precision(16);
875 
876  std::vector< std::pair< C_FLOAT64, CMetab * > >::const_iterator it = pMoiety->getEquation().begin();
877  std::vector< std::pair< C_FLOAT64, CMetab * > >::const_iterator end = pMoiety->getEquation().end();
878  bool First = true;
879 
880  for (; it != end; ++it)
881  {
882  const C_FLOAT64 & Multiplicity = it->first;
883 
884  if (First || Multiplicity < 0.0)
885  {
886  Infix << Multiplicity;
887  }
888  else
889  {
890  Infix << "+" << Multiplicity;
891  }
892 
893  First = false;
894 
895  Infix << "*<";
896  Infix << it->second->getValueReference()->getCN();
897  Infix << ">";
898  }
899 
900  CExpression E("TotalMass", &container);
901 
902  success &= E.setInfix(Infix.str());
903 
904  mpExpression = new CMathExpression(E, container, !mIsInitialValue);
906 
907  return success;
908 }
909 
911 {
912  bool success = true;
913 
914  // The default value is NaN
916 
917  // Reset the prerequisites
918  mPrerequisites.clear();
919 
920  const CMoiety * pMoiety = static_cast< const CMoiety *>(mpDataObject->getObjectParent());
921 
922  std::ostringstream Infix;
923  Infix.imbue(std::locale::classic());
924  Infix.precision(16);
925 
926  Infix << "<" << pMoiety->getTotalNumberReference()->getCN() << ">";
927 
928  std::vector< std::pair< C_FLOAT64, CMetab * > >::const_iterator it = pMoiety->getEquation().begin();
929  std::vector< std::pair< C_FLOAT64, CMetab * > >::const_iterator end = pMoiety->getEquation().end();
930  bool First = true;
931 
932  // The first element in the equation is always the dependent species. We can directly update
933  // its value and therefore point mpValue to it.
934  mpValue = (C_FLOAT64 *) container.getMathObject(it->second->getValueReference())->getValuePointer();
935 
936  ++it;
937 
938  for (; it != end; ++it)
939  {
940  const C_FLOAT64 & Multiplicity = it->first;
941 
942  if (First || Multiplicity >= 0.0)
943  {
944  Infix << "-" << Multiplicity;
945  }
946  else
947  {
948  Infix << "+" << fabs(Multiplicity);
949  }
950 
951  First = false;
952 
953  Infix << "*<";
954  Infix << it->second->getValueReference()->getCN();
955  Infix << ">";
956  }
957 
958  CExpression E("DependentMass", &container);
959 
960  success &= E.setInfix(Infix.str());
961 
962  mpExpression = new CMathExpression(E, container, !mIsInitialValue);
964 
965  return success;
966 }
967 
969 {
970  assert(mpExpression);
971 
972  if (mIsInitialValue)
973  {
975  }
976 
977  mPrerequisites.insert(mpExpression->getPrerequisites().begin(),
978  mpExpression->getPrerequisites().end());
979 
980  if (mPrerequisites.empty())
981  {
982  calculate();
983  }
984 }
985 
987  CMathContainer & container)
988 {
989  assert(pExpression != NULL);
990 
991  bool success = true;
992 
993  bool ReplaceDiscontinousNodes =
994  !mIsInitialValue &&
997 
998  mpExpression = new CMathExpression(*pExpression, container, ReplaceDiscontinousNodes);
1000 
1001  return success;
1002 }
1003 
1005  CMathContainer & container)
1006 {
1007  bool success = true;
1008 
1009  // mConc = *mpValue / mpCompartment->getValue() * mpModel->getNumber2QuantityFactor();
1010  CObjectInterface * pNumber = NULL;
1011  CObjectInterface * pCompartment = NULL;
1012 
1013  if (mIsInitialValue)
1014  {
1015  pNumber = pSpecies->getInitialValueReference();
1016  pCompartment = pSpecies->getCompartment()->getInitialValueReference();
1017  }
1018  else
1019  {
1020  pNumber = pSpecies->getValueReference();
1021  pCompartment = pSpecies->getCompartment()->getValueReference();
1022  }
1023 
1024  std::ostringstream Infix;
1025  Infix.imbue(std::locale::classic());
1026  Infix.precision(16);
1027 
1028  Infix << container.getModel().getNumber2QuantityFactor();
1029  Infix << "*<";
1030  Infix << pNumber->getCN();
1031  Infix << ">/<";
1032  Infix << pCompartment->getCN();
1033  Infix << ">";
1034 
1035  CExpression E("IntensiveValueExpression", &container);
1036 
1037  success &= E.setInfix(Infix.str());
1038 
1039  mpExpression = new CMathExpression(E, container, !mIsInitialValue);
1041 
1042  return success;
1043 }
1044 
1046  CMathContainer & container)
1047 {
1048  bool success = true;
1049 
1050  // mConc * mpCompartment->getValue() * mpModel->getQuantity2NumberFactor();
1051 
1052  CObjectInterface * pDensity = NULL;
1053  CObjectInterface * pCompartment = NULL;
1054 
1055  if (mIsInitialValue)
1056  {
1057  pDensity = pSpecies->getInitialConcentrationReference();
1058  pCompartment = pSpecies->getCompartment()->getInitialValueReference();
1059  }
1060  else
1061  {
1062  pDensity = pSpecies->getConcentrationReference();
1063  pCompartment = pSpecies->getCompartment()->getValueReference();
1064  }
1065 
1066  std::ostringstream Infix;
1067  Infix.imbue(std::locale::classic());
1068  Infix.precision(16);
1069 
1070  Infix << container.getModel().getQuantity2NumberFactor();
1071  Infix << "*<";
1072  Infix << pDensity->getCN();
1073  Infix << ">*<";
1074  Infix << pCompartment->getCN();
1075  Infix << ">";
1076 
1077  CExpression E("ExtensiveValueExpression", &container);
1078 
1079  success &= E.setInfix(Infix.str());
1080 
1081  mpExpression = new CMathExpression(E, container, !mIsInitialValue);
1083 
1084  return success;
1085 }
1086 
1088  CMathContainer & container)
1089 {
1090  bool success = true;
1091 
1092  /*
1093  mConcRate =
1094  (mRate * mpModel->getNumber2QuantityFactor() - mConc * mpCompartment->getRate())
1095  / mpCompartment->getValue();
1096  */
1097 
1098  std::ostringstream Infix;
1099  Infix.imbue(std::locale::classic());
1100  Infix.precision(16);
1101 
1102  Infix << "(<";
1103  Infix << pSpecies->getRateReference()->getCN();
1104  Infix << ">*";
1105  Infix << container.getModel().getNumber2QuantityFactor();
1106  Infix << "-<";
1107  Infix << pSpecies->getCompartment()->getValueReference()->getCN();
1108  Infix << ">*<";
1109  Infix << pSpecies->getCompartment()->getRateReference()->getCN();
1110  Infix << ">)/<";
1111  Infix << pSpecies->getCompartment()->getValueReference()->getCN();
1112  Infix << ">";
1113 
1114  CExpression E("IntensiveRateExpression", &container);
1115 
1116  success &= E.setInfix(Infix.str());
1117 
1118  mpExpression = new CMathExpression(E, container, !mIsInitialValue);
1120 
1121  return success;
1122 }
1123 
1125  CMathContainer & container)
1126 {
1127  bool success = true;
1128 
1129  /*
1130  mRate = mpModel->getQuantity2NumberFactor() *
1131  mpCompartment->getValue() * mpExpression->calcValue();
1132  */
1133 
1134  std::ostringstream Infix;
1135  Infix.imbue(std::locale::classic());
1136  Infix.precision(16);
1137 
1138  Infix << container.getModel().getQuantity2NumberFactor();
1139  Infix << "*<";
1140  Infix << pSpecies->getCompartment()->getValueReference()->getCN();
1141  Infix << ">*";
1142  Infix << pSpecies->getExpression();
1143 
1144  CExpression E("ExtensiveODERateExpression", &container);
1145 
1146  success &= E.setInfix(Infix.str());
1147 
1148  mpExpression = new CMathExpression(E, container, !mIsInitialValue);
1150 
1151  return success;
1152 }
1153 
1155  CMathContainer & container)
1156 {
1157  bool success = true;
1158 
1159  std::ostringstream Infix;
1160  Infix.imbue(std::locale::classic());
1161  Infix.precision(16);
1162 
1163  std::string Key = pSpecies->getKey();
1164  bool First = true;
1165 
1168 
1169  for (; it != end; ++it)
1170  {
1171  const CCopasiVector< CChemEqElement > &Balances =
1172  (*it)->getChemEq().getBalances();
1175 
1176  for (; itChem != endChem; ++itChem)
1177  if ((*itChem)->getMetaboliteKey() == Key)
1178  break;
1179 
1180  if (itChem != endChem)
1181  {
1182  const C_FLOAT64 & Multiplicity = (*itChem)->getMultiplicity();
1183 
1184  if (First || Multiplicity < 0.0)
1185  {
1186  if (Multiplicity == std::numeric_limits< C_FLOAT64 >::infinity())
1187  {
1188  Infix << "infinity";
1189  }
1190  else if (Multiplicity == -std::numeric_limits< C_FLOAT64 >::infinity())
1191  {
1192  Infix << "-infinity";
1193  }
1194  else
1195  {
1196  Infix << Multiplicity;
1197  }
1198  }
1199  else
1200  {
1201  if (Multiplicity == std::numeric_limits< C_FLOAT64 >::infinity())
1202  {
1203  Infix << "+infinity";
1204  }
1205  else
1206  {
1207  Infix << "+" << Multiplicity;
1208  }
1209  }
1210 
1211  First = false;
1212 
1213  Infix << "*<";
1214  Infix << (*it)->getParticleFluxReference()->getCN();
1215  Infix << ">";
1216  }
1217  }
1218 
1219  CExpression E("ExtensiveReactionExpression", &container);
1220 
1221  success &= E.setInfix(Infix.str());
1222 
1223  mpExpression = new CMathExpression(E, container, !mIsInitialValue);
1225 
1226  return success;
1227 }
1228 
1229 std::ostream &operator<<(std::ostream &os, const CMathObject & o)
1230 {
1231  if (o.mpDataObject != NULL)
1232  {
1233  os << o.mpDataObject->getCN() << std::endl;
1234  }
1235  else
1236  {
1237  os << "Data Object = NULL" << std::endl;
1238  }
1239 
1240  os << " Pointer: " << &o << std::endl;
1241  os << " Value Type: ";
1242 
1243  switch (o.mValueType)
1244  {
1246  os << "ValueTypeUndefined" << std::endl;
1247  break;
1248 
1249  case CMath::Value:
1250  os << "Value" << std::endl;
1251  break;
1252 
1253  case CMath::Rate:
1254  os << "ValueRate" << std::endl;
1255  break;
1256 
1257  case CMath::ParticleFlux:
1258  os << "ParticleFlux" << std::endl;
1259  break;
1260 
1261  case CMath::Flux:
1262  os << "Flux" << std::endl;
1263  break;
1264 
1265  case CMath::Propensity:
1266  os << "Propensity" << std::endl;
1267  break;
1268 
1269  case CMath::TotalMass:
1270  os << "TotalMass" << std::endl;
1271  break;
1272 
1273  case CMath::DependentMass:
1274  os << "DependentMass" << std::endl;
1275  break;
1276 
1277  case CMath::Discontinuous:
1278  os << "Discontinuous" << std::endl;
1279  break;
1280 
1281  case CMath::EventDelay:
1282  os << "EventDelay" << std::endl;
1283  break;
1284 
1285  case CMath::EventPriority:
1286  os << "EventPriority" << std::endl;
1287  break;
1288 
1290  os << "EventAssignment" << std::endl;
1291  break;
1292 
1293  case CMath::EventTrigger:
1294  os << "EventTrigger" << std::endl;
1295  break;
1296 
1297  case CMath::EventRoot:
1298  os << "EventRoot" << std::endl;
1299  break;
1300 
1301  case CMath::EventRootState:
1302  os << "EventRootState" << std::endl;
1303  break;
1304  }
1305 
1306  os << " Simulation Type: ";
1307 
1308  switch (o.mSimulationType)
1309  {
1311  os << "SimulationTypeUndefined" << std::endl;
1312  break;
1313 
1314  case CMath::Fixed:
1315  os << "Fixed" << std::endl;
1316  break;
1317 
1318  case CMath::EventTarget:
1319  os << "EventTarget" << std::endl;
1320  break;
1321 
1322  case CMath::Time:
1323  os << "Time" << std::endl;
1324  break;
1325 
1326  case CMath::ODE:
1327  os << "ODE" << std::endl;
1328  break;
1329 
1330  case CMath::Independent:
1331  os << "Independent" << std::endl;
1332  break;
1333 
1334  case CMath::Dependent:
1335  os << "Dependent" << std::endl;
1336  break;
1337 
1338  case CMath::Assignment:
1339  os << "Assignment" << std::endl;
1340  break;
1341 
1342  case CMath::Conversion:
1343  os << "Conversion" << std::endl;
1344  break;
1345  };
1346 
1347  os << " Entity Type: ";
1348 
1349  switch (o.mEntityType)
1350  {
1352  os << "EntityTypeUndefined" << std::endl;
1353  break;
1354 
1355  case CMath::Model:
1356  os << "Model" << std::endl;
1357  break;
1358 
1359  case CMath::Analysis:
1360  os << "Analysis" << std::endl;
1361  break;
1362 
1363  case CMath::GlobalQuantity:
1364  os << "GlobalQuantity" << std::endl;
1365  break;
1366 
1367  case CMath::Compartment:
1368  os << "Compartment" << std::endl;
1369  break;
1370 
1371  case CMath::Species:
1372  os << "Species" << std::endl;
1373  break;
1374 
1376  os << "LocalReactionParameter" << std::endl;
1377  break;
1378 
1380  os << "StoichiometricCoefficients" << std::endl;
1381  break;
1382 
1383  case CMath::Reaction:
1384  os << "Reaction" << std::endl;
1385  break;
1386 
1387  case CMath::Moiety:
1388  os << "Moiety" << std::endl;
1389  break;
1390 
1391  case CMath::Event:
1392  os << "Event" << std::endl;
1393  break;
1394  };
1395 
1396  os << " Is Intensive Property: " << (o.mIsIntensiveProperty ? "true" : "false") << std::endl;
1397 
1398  os << " Is Initial Value: " << (o.mIsInitialValue ? "true" : "false") << std::endl;
1399 
1400  os << " Intensive Property: ";
1401 
1402  if (o.mpIntensiveProperty != NULL)
1403  {
1404  os << o.mpIntensiveProperty->getCN() << std::endl;
1405  }
1406  else
1407  {
1408  os << "NULL" << std::endl;
1409  }
1410 
1411  os << " Value: " << *o.mpValue << " (" << o.mpValue << ")" << std::endl;
1412 
1413  os << " Expression: ";
1414 
1415  if (o.mpExpression != NULL)
1416  {
1417  os << o.mpExpression->getRoot()->buildInfix() << std::endl;
1418  }
1419  else
1420  {
1421  os << "NULL" << std::endl;
1422  }
1423 
1424  // CObjectInterface::ObjectSet mPrerequisites;
1425 
1426  return os;
1427 }
const CExpression * getExpressionPtr() const
Header file of class CExpression.
const CMathExpression * getExpressionPtr() const
CObjectInterface::ObjectSet mPrerequisites
Definition: CMathObject.h:331
void setSimulationType(const CMath::SimulationType &simulationType)
const CMath::EntityType & getEntityType() const
virtual bool setInfix(const std::string &infix)
Definition: CExpression.cpp:63
#define pdelete(p)
Definition: copasi.h:215
const CMoiety * getMoiety() const
Definition: CMetab.cpp:992
CCopasiObject * getParticleFluxReference()
Definition: CReaction.cpp:207
bool createExtensiveReactionRateExpression(const CMetab *pSpecies, CMathContainer &container)
virtual const CObjectInterface::ObjectSet & getPrerequisites() const
virtual CCopasiObjectName getCN() const
CMathObject * getMathObject(const CObjectInterface *pObject) const
SimulationType
Definition: CMathEnum.h:182
const std::string & getObjectName() const
virtual const CObjectInterface * getObject(const CCopasiObjectName &cn) const
virtual CCopasiObjectName getCN() const =0
void copy(const CMathObject &src, CMathContainer &container, const size_t &valueOffset, const size_t &objectOffset)
Definition: CMathObject.cpp:64
CCopasiObject * getInitialValueReference() const
const CMath::SimulationType & getSimulationType() const
bool compilePropensity(CMathContainer &container)
iterator begin()
const Type & getType() const
const CMathObject * mpIntensiveProperty
Definition: CMathObject.h:362
const CCallParameters< C_FLOAT64 > & getCallParameters() const
Definition: CReaction.cpp:216
CCopasiObject * getTotalNumberReference() const
Definition: CMoiety.cpp:148
std::string buildInfix() const
virtual CCopasiObjectName getCN() const
bool mIsIntensiveProperty
Definition: CMathObject.h:351
std::string getExpression() const
void compileExpression()
virtual bool compile(std::vector< CCopasiContainer * > listOfContainer=CCopasiContainer::EmptyList)
Definition: CExpression.cpp:97
const std::vector< std::pair< C_FLOAT64, CMetab * > > & getEquation() const
Definition: CMoiety.cpp:267
C_FLOAT64 * mpValue
Definition: CMathObject.h:326
bool compile(CMathContainer &container)
CMath::EntityType mEntityType
Definition: CMathObject.h:341
void setIsBoolean(const bool &booleanRequired)
Definition: CExpression.cpp:58
Definition: CMetab.h:178
virtual const CObjectInterface::ObjectSet & getPrerequisites() const
virtual void * getValuePointer() const
static C_FLOAT64 InvalidValue
Definition: CMathObject.h:316
bool compileRate(CMathContainer &container)
std::set< const CCompartment * > getCompartments() const
Definition: CChemEq.cpp:115
const CFunction * getFunction() const
Definition: CReaction.cpp:252
EntityType
Definition: CMathEnum.h:195
bool compileTotalMass(CMathContainer &container)
iterator end()
const C_FLOAT64 & getQuantity2NumberFactor() const
Definition: CModel.cpp:2354
const C_FLOAT64 & getNumber2QuantityFactor() const
Definition: CModel.cpp:2357
ValueType
Definition: CMathEnum.h:163
virtual const std::string & getKey() const
const bool & isInitialValue() const
CMath::ValueType mValueType
Definition: CMathObject.h:336
bool mIsInitialValue
Definition: CMathObject.h:357
bool setExpression(const std::string &infix, const bool &isBoolean, CMathContainer &container)
const C_FLOAT64 & value()
bool compileInitialValue(CMathContainer &container)
const CCopasiVector< CChemEqElement > & getSubstrates() const
Definition: CChemEq.cpp:60
bool convertToInitialExpression()
bool isReversible() const
Definition: CReaction.cpp:229
virtual bool compile()
bool createIntensiveValueExpression(const CMetab *pSpecies, CMathContainer &container)
virtual void print(std::ostream *ostream) const
bool compileValue(CMathContainer &container)
std::set< const CObjectInterface * > ObjectSet
static CMathExpression * copy(const CMathExpression &src, CMathContainer &container, const size_t &valueOffset, const size_t &objectOffset)
const CMath::ValueType & getValueType() const
bool compileParticleFlux(CMathContainer &container)
#define C_FLOAT64
Definition: copasi.h:92
static void initialize(CMathObject *&pObject, C_FLOAT64 *&pValue, const CMath::ValueType &valueType, const CMath::EntityType &entityType, const CMath::SimulationType &simulationType, const bool &isIntensiveProperty, const bool &isInitialValue, const CCopasiObject *pDataObject)
Definition: CMathObject.cpp:23
virtual const CObjectInterface * getObject(const CCopasiObjectName &cn) const
const CModel & getModel() const
bool compileDependentMass(CMathContainer &container)
const CCopasiObject * mpDataObject
Definition: CMathObject.h:367
const bool & isIntensiveProperty() const
virtual ~CMathObject()
Definition: CMathObject.cpp:59
const ModelType & getModelType() const
Definition: CModel.cpp:2339
CCopasiVectorNS< CReaction > & getReactions()
Definition: CModel.cpp:1039
bool createIntensiveRateExpression(const CMetab *pSpecies, CMathContainer &container)
CMathExpression * mpExpression
Definition: CMathObject.h:321
bool compileFlux(CMathContainer &container)
virtual void * getValuePointer() const
void calculate()
bool createExtensiveValueExpression(const CMetab *pSpecies, CMathContainer &container)
CConcentrationReference * getInitialConcentrationReference() const
Definition: CMetab.cpp:861
bool createExtensiveODERateExpression(const CMetab *pSpecies, CMathContainer &container)
CCopasiObject * getRateReference() const
const CCompartment * getCompartment() const
Definition: CMetab.cpp:222
const CChemEq & getChemEq() const
Definition: CReaction.cpp:223
CEvaluationNode * getRoot()
bool createConvertedExpression(const CExpression *pExpression, CMathContainer &container)
const std::string & getInfix() const
CCopasiObject * getDependentNumberReference() const
Definition: CMoiety.cpp:151
CCopasiContainer * getObjectParent() const
const CExpression * getInitialExpressionPtr() const
bool setExpressionPtr(CMathExpression *pMathExpression)
std::ostream & operator<<(std::ostream &os, const CMathObject &o)
const CCopasiObject * getDataObject() const
virtual bool isPrerequisiteForContext(const CObjectInterface *pObject, const CMath::SimulationContextFlag &context, const CObjectInterface::ObjectSet &changedObjects) const
CCopasiObject * getFluxReference()
Definition: CReaction.cpp:198
CConcentrationReference * getConcentrationReference() const
Definition: CMetab.cpp:864
CMath::SimulationType mSimulationType
Definition: CMathObject.h:346
CCopasiObject * getValueReference() const