COPASI API  4.16.103
CEvaluationNodeOperator.cpp
Go to the documentation of this file.
1 // Copyright (C) 2010 - 2014 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) 2005 - 2007 by Pedro Mendes, Virginia Tech Intellectual
12 // Properties, Inc. and EML Research, gGmbH.
13 // All rights reserved.
14 
15 #include "copasi.h"
16 #include "CEvaluationNode.h"
17 #include <sstream>
18 #include "CEvaluationTree.h"
19 #include "sbml/math/ASTNode.h"
20 
21 CEvaluationNodeOperator::CEvaluationNodeOperator():
22  CEvaluationNode(CEvaluationNode::INVALID, ""),
23  mpLeft(NULL),
24  mpRight(NULL)
25 {}
26 
27 CEvaluationNodeOperator::CEvaluationNodeOperator(const SubType & subType,
28  const Data & data):
29  CEvaluationNode((Type)(CEvaluationNode::OPERATOR | subType), data),
30  mpLeft(NULL),
31  mpRight(NULL)
32 {
33  switch (mType & 0x00FFFFFF)
34  {
35  case POWER:
36  mPrecedence = PRECEDENCE_OPERATOR_POWER;
37  break;
38 
39  case MULTIPLY:
40  mPrecedence = PRECEDENCE_OPERATOR_MULTIPLY;
41  break;
42 
43  case DIVIDE:
44  mPrecedence = PRECEDENCE_OPERATOR_DIVIDE;
45  break;
46 
47  case MODULUS:
48  mPrecedence = PRECEDENCE_OPERATOR_MODULUS;
49  break;
50 
51  case PLUS:
52  mPrecedence = PRECEDENCE_OPERATOR_PLUS;
53  break;
54 
55  case MINUS:
56  mPrecedence = PRECEDENCE_OPERATOR_MINUS;
57  break;
58 
59  default:
60  break;
61  }
62 }
63 
64 CEvaluationNodeOperator::CEvaluationNodeOperator(const CEvaluationNodeOperator & src):
65  CEvaluationNode(src),
66  mpLeft(NULL),
67  mpRight(NULL)
68 {}
69 
70 CEvaluationNodeOperator::~CEvaluationNodeOperator() {}
71 
72 bool CEvaluationNodeOperator::compile(const CEvaluationTree * /* pTree */)
73 {
74  mpLeft = static_cast<CEvaluationNode *>(getChild());
75 
76  if (mpLeft == NULL) return false;
77 
78  mpRight = static_cast<CEvaluationNode *>(mpLeft->getSibling());
79 
80  if (mpRight == NULL) return false;
81 
82  return (mpRight->getSibling() == NULL); // We must have only two children
83 }
84 
85 // virtual
86 std::string CEvaluationNodeOperator::getInfix(const std::vector< std::string > & children) const
87 {
88  if (const_cast<CEvaluationNodeOperator *>(this)->compile(NULL))
89  {
90  Data Infix;
91 
92  if (*mpLeft < * (CEvaluationNode *)this)
93  Infix = "(" + children[0] + ")";
94  else
95  Infix = children[0];
96 
97  Infix += mData;
98 
99  if (!(*(CEvaluationNode *)this < *mpRight))
100  Infix += "(" + children[1] + ")";
101  else
102  Infix += children[1];
103 
104  return Infix;
105  }
106  else
107  return "@";
108 }
109 
110 // virtual
111 std::string CEvaluationNodeOperator::getDisplayString(const std::vector< std::string > & children) const
112 {
113  if (const_cast<CEvaluationNodeOperator *>(this)->compile(NULL))
114  {
115  Data DisplayString;
116 
117  if (*mpLeft < * (CEvaluationNode *)this)
118  DisplayString = "(" + children[0] + ")";
119  else
120  DisplayString = children[0];
121 
122  DisplayString += mData;
123 
124  if (!(*(CEvaluationNode *)this < *mpRight))
125  DisplayString += "(" + children[1] + ")";
126  else
127  DisplayString += children[1];
128 
129  return DisplayString;
130  }
131  else
132  return "@";
133 }
134 
135 // virtual
136 std::string CEvaluationNodeOperator::getCCodeString(const std::vector< std::string > & children) const
137 {
138  if (const_cast<CEvaluationNodeOperator *>(this)->compile(NULL))
139  {
140  Data DisplayString;
141  SubType subType = (SubType)CEvaluationNode::subType(this->getType());
142 
143  if (subType == POWER)
144  DisplayString = "pow(";
145 
146  if (subType == MODULUS)
147  DisplayString = "(int)";
148 
149  if (*mpLeft < * (CEvaluationNode *)this)
150  DisplayString += "(" + children[0] + ")";
151  else
152  DisplayString += children[0];
153 
154  switch (subType)
155  {
156  case POWER:
157  DisplayString += ",";
158  break;
159 
160  case MODULUS:
161  DisplayString += "%(int)";
162  break;
163 
164  default:
165  DisplayString += mData;
166  break;
167  }
168 
169  if (!(*(CEvaluationNode *)this < *mpRight))
170  DisplayString += "(" + children[1] + ")";
171  else
172  DisplayString += children[1];
173 
174  if (subType == POWER)
175  DisplayString += ")";
176 
177  return DisplayString;
178  }
179  else
180  return "@";
181 }
182 
183 // virtual
184 std::string CEvaluationNodeOperator::getBerkeleyMadonnaString(const std::vector< std::string > & children) const
185 {
186  if (const_cast<CEvaluationNodeOperator *>(this)->compile(NULL))
187  {
188  std::string mdata = "";
189 
190  /* if ((SubType)CEvaluationNode::subType(this->getType()) == MODULUS)
191  mdata = "@";
192  else */
193  mdata = mData;
194 
195  Data DisplayString;
196 
197  if (*mpLeft < * (CEvaluationNode *)this)
198  DisplayString = "(" + children[0] + ")";
199  else
200  DisplayString = children[0];
201 
202  DisplayString += mdata;
203 
204  if (!(*(CEvaluationNode *)this < *mpRight))
205  DisplayString += "(" + children[1] + ")";
206  else
207  DisplayString += children[1];
208 
209  return DisplayString;
210  }
211  else
212  return "@";
213 }
214 
215 // virtual
216 std::string CEvaluationNodeOperator::getXPPString(const std::vector< std::string > & children) const
217 {
218  if (const_cast<CEvaluationNodeOperator *>(this)->compile(NULL))
219  {
220  Data DisplayString;
221  SubType subType = (SubType)CEvaluationNode::subType(this->getType());
222 
223  if (subType == MODULUS)
224  DisplayString = "mod(";
225 
226  if (*mpLeft < * (CEvaluationNode *)this)
227  DisplayString += "(" + children[0] + ")";
228  else
229  DisplayString += children[0];
230 
231  switch (subType)
232  {
233  case MODULUS:
234  DisplayString += ",";
235  break;
236 
237  default:
238  DisplayString += mData;
239  break;
240  }
241 
242  if (!(*(CEvaluationNode *)this < *mpRight))
243  DisplayString += "(" + children[1] + ")";
244  else
245  DisplayString += children[1];
246 
247  if (subType == MODULUS)
248  DisplayString += ")";
249 
250  return DisplayString;
251  }
252  else
253  return "@"; //TODO
254 }
255 
256 // static
257 CEvaluationNode * CEvaluationNodeOperator::fromAST(const ASTNode * pASTNode, const std::vector< CEvaluationNode * > & children)
258 {
259  assert(pASTNode->getNumChildren() == children.size());
260 
261  size_t i = 0, iMax = children.size();
262 
263  int type = (int)pASTNode->getType();
264 
265  SubType subType;
266  std::string data = "";
267 
268  switch (type)
269  {
270  case AST_PLUS:
271  subType = PLUS;
272  data = "+";
273  break;
274 
275  case AST_MINUS:
276  subType = MINUS;
277  data = "-";
278  break;
279 
280  case AST_TIMES:
281  subType = MULTIPLY;
282  data = "*";
283  break;
284 
285  case AST_DIVIDE:
286  subType = DIVIDE;
287  data = "/";
288  break;
289 
290  case AST_POWER:
291  case AST_FUNCTION_POWER:
292  subType = POWER;
293  data = "^";
294  break;
295 
296  default:
297  subType = INVALID;
298  fatalError();
299  break;
300  }
301 
302  CEvaluationNode* pNode = NULL;
303 
304  // handle "-" since it can be unary or binary
305  if (type == AST_MINUS)
306  {
307  switch (iMax)
308  {
309  case 1:
310  delete pNode;
311  pNode = new CEvaluationNodeFunction(CEvaluationNodeFunction::MINUS, data);
312  pNode->addChild(children[0]);
313  break;
314 
315  case 2:
316  pNode = new CEvaluationNodeOperator(subType, data);
317  pNode->addChild(children[0]);
318  pNode->addChild(children[1]);
319  break;
320 
321  default:
322  // error
323  CCopasiMessage(CCopasiMessage::EXCEPTION, MCMathML + 3);
324  break;
325  }
326  }
327  // handle binary operators (POWER,/)
328  else if (type == AST_DIVIDE || type == AST_POWER || type == AST_FUNCTION_POWER)
329  {
330  switch (pASTNode->getNumChildren())
331  {
332  case 2:
333  pNode = new CEvaluationNodeOperator(subType, data);
334  pNode->addChild(children[0]);
335  pNode->addChild(children[1]);
336  break;
337 
338  default:
339  // error
340  CCopasiMessage(CCopasiMessage::EXCEPTION, MCMathML + 4);
341  break;
342  }
343  }
344  // handle n-ary operators (+,*)
345  else if (type == AST_PLUS || type == AST_TIMES)
346  {
347  switch (iMax)
348  {
349  case 0:
350 
351  if (type == AST_PLUS)
352  {
353  pNode = new CEvaluationNodeNumber(CEvaluationNodeNumber::DOUBLE, "0.0");
354  }
355  else
356  {
357  pNode = new CEvaluationNodeNumber(CEvaluationNodeNumber::DOUBLE, "1.0");
358  }
359 
360  break;
361 
362  case 1:
363  // replace the current node with its only child
364  pNode = children[0];
365  break;
366 
367  case 2:
368  {
369  pNode = new CEvaluationNodeOperator(subType, data);
370  CEvaluationNode * pCurrent = pNode;
371 
372  // We have at least 2 children
373  while (i < iMax - 1)
374  {
375  // add the first value
376  pCurrent->addChild(children[i++]);
377 
378  switch (iMax - i)
379  {
380  case 1:
381  // We have only 1 more child
382  pCurrent->addChild(children[i++]);
383  break;
384 
385  default:
386  // We have at least 2 more children
387  {
388  // create a new node with the same operator
389  CEvaluationNode * pTmp = new CEvaluationNodeOperator(subType, data);
390  pCurrent->addChild(pTmp);
391  pCurrent = pTmp;
392  }
393  break;
394  }
395  }
396  }
397  break;
398  }
399  }
400 
401  return pNode;
402 }
403 
404 ASTNode* CEvaluationNodeOperator::toAST(const CCopasiDataModel* pDataModel) const
405 {
406  SubType subType = (SubType)CEvaluationNode::subType(this->getType());
407  ASTNode* node = new ASTNode();
408 
409  switch (subType)
410  {
411  case POWER:
412  node->setType(AST_POWER);
413  break;
414 
415  case MULTIPLY:
416  node->setType(AST_TIMES);
417  break;
418 
419  case DIVIDE:
420  node->setType(AST_DIVIDE);
421  break;
422 
423  case MODULUS:
424  // replace this with a more complex subtree
425  CEvaluationNodeOperator::createModuloTree(this, node, pDataModel);
426  break;
427 
428  case PLUS:
429  node->setType(AST_PLUS);
430  break;
431 
432  case MINUS:
433  node->setType(AST_MINUS);
434  break;
435 
436  case INVALID:
437  break;
438  }
439 
440  // for all but INVALID and MODULUS two children have to be converted
441  if (subType != INVALID && subType != MODULUS)
442  {
443  const CEvaluationNode* child1 = dynamic_cast<const CEvaluationNode*>(this->getChild());
444  const CEvaluationNode* child2 = dynamic_cast<const CEvaluationNode*>(child1->getSibling());
445  node->addChild(child1->toAST(pDataModel));
446  node->addChild(child2->toAST(pDataModel));
447  }
448 
449  return node;
450 }
451 
452 /*
453 CEvaluationNode* CEvaluationNodeOperator::simplifyNode(CEvaluationNode *child1, CEvaluationNode *child2) const
454 {
455  std::vector<CEvaluationNode*> children;
456  children.push_back(child1);
457  children.push_back(child1);
458  return this->simplifyNode(children);
459 }
460  */
461 
462 CEvaluationNode* CEvaluationNodeOperator::simplifyNode(const std::vector<CEvaluationNode*>& children) const
463 {
464  assert(children.size() > 0);
465  CEvaluationNode* child1 = children[0];
466  CEvaluationNode* child2 = NULL;
467 
468  if (children.size() > 1)
469  {
470  child2 = children[1];
471  }
472 
473  switch ((SubType) CEvaluationNode::subType(mType))
474  {
475  case POWER:
476  {
477  if (CEvaluationNode::type(child2->getType()) == NUMBER)
478  {
479  if (CEvaluationNode::type(child1->getType()) == NUMBER)
480  {
481  // both children numbers ->calculate
482  std::stringstream tmp;
483  tmp.precision(18);
484  tmp << pow(child1->getValue(), child2->getValue());
485  CEvaluationNode *newnode = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), tmp.str());
486  delete child1;
487  delete child2;
488  return newnode;
489  }
490 
491  if (fabs(child2->getValue() - 1.0) < 1.0E-100)
492  {
493  // a^1 -> a
494  delete child2;
495  return child1;
496  }
497 
498  if (fabs(child2->getValue()) < 1.0E-100)
499  {
500  // a^0 -> 1
501  CEvaluationNode* newnode = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), "1.0");
502  delete child1;
503  delete child2;
504  return newnode;
505  }
506 
507  if (child2->getValue() < 0.0)
508  {
509  //negative float exponent ->write as fraction
510  CEvaluationNode* newnode = CEvaluationNode::create((Type)(OPERATOR | DIVIDE), "/");
511  CEvaluationNode* newchild1 = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), "1.0");
512  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
513  CEvaluationNode* grandchild1 = child1;
514  std::stringstream tmp;
515  tmp.precision(18);
516  tmp << fabs(child2->getValue());
517  CEvaluationNode* grandchild2 = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), tmp.str());
518  newnode->addChild(newchild1, NULL);
519  newnode->addChild(newchild2, newchild1);
520  newchild2->addChild(grandchild1, NULL);
521  newchild2->addChild(grandchild2, grandchild1);
522  delete child2;
523  return newnode;
524  }
525 
526  if ((child1->getType() == (CEvaluationNode::Type)(FUNCTION | CEvaluationNodeFunction::MINUS))
527  && ((fabs(child2->getValue() - floor(child2->getValue())) < 1.0E-100) || (fabs(child2->getValue() - floor(child2->getValue()) - 1.0) < 1.0E-100)))
528  {
529  // (-a)^n -> (-1)^n * a^n for n int
530  // make negativity a property of product,
531  // easier to recognize a negative exponent.
532  int exp;
533 
534  if (fabs(child2->getValue() - floor(child2->getValue())) < 1.0E-100)
535  exp = (int) floor(child2->getValue());
536  else
537  exp = (int) floor(child2->getValue()) + 1;
538 
539  CEvaluationNode* newpower = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
540  CEvaluationNode* newchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild())->copyBranch();
541  newpower->addChild(newchild1, NULL);
542  newpower->addChild(child2, newchild1);
543  delete child1;
544 
545  if (pow(-1.0, exp) == 1.0)
546  return newpower;
547 
548  CEvaluationNode* newnode = CEvaluationNode::create((Type)(FUNCTION | CEvaluationNodeFunction::MINUS), "-");
549  newnode->addChild(newpower, NULL);
550  return newnode;
551  }
552  }
553 
554  if (child1->getData() == "^")
555  {
556  // (a^b)^c -> a^(b*c)
557  CEvaluationNode* newnode = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
558  CEvaluationNode* newchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild())->copyBranch();
559  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
560  CEvaluationNode* grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild()->getSibling())->copyBranch();
561  CEvaluationNode* grandchild2 = child2;
562  newnode->addChild(newchild1, NULL);
563  newnode->addChild(newchild2, newchild1);
564  newchild2->addChild(grandchild1, NULL);
565  newchild2->addChild(grandchild2, grandchild1);
566  delete child1;
567  return newnode;
568  }
569 
570  /* disable this since negative exponents are actually generated
571  during the simplification. the normalization should get rid of
572  those.
573  if (child2->getType() == (Type)(FUNCTION | CEvaluationNodeFunction::MINUS))
574  {// a^(-b) -> 1/(a^b)
575  // don't want negative exponents
576  CEvaluationNode* newnode = CEvaluationNode::create((Type)(OPERATOR | DIVIDE), "/");
577  CEvaluationNode* newchild1 = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), "1.0");
578  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
579  CEvaluationNode* grandchild1 = child1;
580  CEvaluationNode* grandchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild())->copyBranch();
581  newnode->addChild(newchild1, NULL);
582  newnode->addChild(newchild2, newchild1);
583  newchild2->addChild(grandchild1, NULL);
584  newchild2->addChild(grandchild2, grandchild1);
585  delete child2;
586  return newnode;
587  }
588  */
589  if (child1->getData() == "*")
590  {
591  // (a*b)^c -> a^c * b^c
592  CEvaluationNode* newnode = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
593  CEvaluationNode* newchild1 = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
594  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
595  CEvaluationNode * grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild())->copyBranch();
596  CEvaluationNode * grandchild2 = child2->copyBranch();
597  CEvaluationNode * grandchild3 = dynamic_cast<CEvaluationNode*>(child1->getChild()->getSibling())->copyBranch();
598  CEvaluationNode * grandchild4 = child2;
599  newnode->addChild(newchild1, NULL);
600  newnode->addChild(newchild2, newchild1);
601  newchild1->addChild(grandchild1, NULL);
602  newchild1->addChild(grandchild2, grandchild1);
603  newchild2->addChild(grandchild3, NULL);
604  newchild2->addChild(grandchild4, grandchild3);
605  delete child1;
606  return newnode;
607  }
608 
609  if (child1->getData() == "/")
610  {
611  // (a/b)^c -> a^c/b^c
612  CEvaluationNode* newnode = CEvaluationNode::create((Type)(OPERATOR | DIVIDE), "/");
613  CEvaluationNode* newchild1 = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
614  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
615  CEvaluationNode * grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild())->copyBranch();
616  CEvaluationNode * grandchild2 = child2->copyBranch();
617  CEvaluationNode * grandchild3 = dynamic_cast<CEvaluationNode*>(child1->getChild()->getSibling())->copyBranch();
618  CEvaluationNode * grandchild4 = child2;
619  newnode->addChild(newchild1, NULL);
620  newnode->addChild(newchild2, newchild1);
621  newchild1->addChild(grandchild1, NULL);
622  newchild1->addChild(grandchild2, grandchild1);
623  newchild2->addChild(grandchild3, NULL);
624  newchild2->addChild(grandchild4, grandchild3);
625  delete child1;
626  return newnode;
627  }
628 
629  /*
630  if ((child1->getData() == "+") && (CEvaluationNode::type(child2->getType()) == NUMBER) && (child2->getValue() >= 1.0 + 1.0E-100))
631  {// (a+b)^x -> (a+b) * (a+b)^(x-1) for real x > 1
632  // this is expanded step by step
633  CEvaluationNode* newnode = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
634  CEvaluationNode* newchild1 = child1->copyBranch();
635  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
636  CEvaluationNode * grandchild1 = child1;
637  std::stringstream tmp;
638  tmp.precision(18);
639  tmp << child2->getValue() - 1.0;
640  CEvaluationNode * grandchild2 = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), tmp.str());
641  newnode->addChild(newchild1, NULL);
642  newnode->addChild(newchild2, newchild1);
643  newchild2->addChild(grandchild1, NULL);
644  newchild2->addChild(grandchild2, grandchild1);
645  delete child2;
646  return newnode;
647  }
648  // don't do this during the simplification. This has to be done
649  // during the normalization
650  if (child2->getData() == "+")
651  {// a^(b+c) -> a^b*a^c minimize exponent
652  CEvaluationNode* newnode = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
653  CEvaluationNode* newchild1 = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
654  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
655  CEvaluationNode * grandchild1 = child1->copyBranch();
656  CEvaluationNode * grandchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild())->copyBranch();
657  CEvaluationNode * grandchild3 = child1;
658  CEvaluationNode * grandchild4 = dynamic_cast<CEvaluationNode*>(child2->getChild()->getSibling())->copyBranch();
659  newnode->addChild(newchild1, NULL);
660  newnode->addChild(newchild2, newchild1);
661  newchild1->addChild(grandchild1, NULL);
662  newchild1->addChild(grandchild2, grandchild1);
663  newchild2->addChild(grandchild3, NULL);
664  newchild2->addChild(grandchild4, grandchild3);
665  delete child2;
666  return newnode;
667  }
668  */
669  CEvaluationNode *newnode = copyNode(child1, child2);
670  return newnode;
671  }
672 
673  case MULTIPLY:
674  {
675  if (CEvaluationNode::type(child1->getType()) == NUMBER)
676  {
677  if (CEvaluationNode::type(child2->getType()) == NUMBER)
678  {
679  // both children numbers ->calculate
680  std::stringstream tmp;
681  tmp.precision(18);
682  tmp << child1->getValue() * child2->getValue();
683  CEvaluationNode *newnode = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), tmp.str());
684  delete child1;
685  delete child2;
686  return newnode;
687  }
688 
689  if (fabs(child1->getValue()) < 1.0E-100)
690  {
691  // 0*a -> 0
692  CEvaluationNode *newnode = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), "0.0");
693  delete child1;
694  delete child2;
695  return newnode;
696  }
697 
698  if (fabs(child1->getValue() - 1.0) < 1.0E-100)
699  {
700  // 1*a -> a
701  delete child1;
702  return child2;
703  }
704 
705  if (child1->getValue() < 0.0)
706  {
707  // x*a -> -(|x|*a) x < 0
708  // make negativity a property of product,
709  // easier to recognize a negative exponent.
710  CEvaluationNode *newnode = CEvaluationNode::create((Type)(FUNCTION | CEvaluationNodeFunction::MINUS), "-");
711  CEvaluationNode *newchild1 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
712  std::stringstream tmp;
713  tmp.precision(18);
714  tmp << fabs(child1->getValue());
715  CEvaluationNode *grandchild1 = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), tmp.str());
716  CEvaluationNode *grandchild2 = child2;
717  newnode->addChild(newchild1, NULL);
718  newchild1->addChild(grandchild1, NULL);
719  newchild1->addChild(grandchild2, grandchild1);
720  delete child1;
721  return newnode;
722  }
723  }
724 
725  if (CEvaluationNode::type(child2->getType()) == NUMBER)
726  {
727  // because of commutativity the same as before for child1.
728  if (fabs(child2->getValue()) < 1.0E-100)
729  {
730  // a*0 -> 0
731  CEvaluationNode *newnode = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), "0.0");
732  delete child1;
733  delete child2;
734  return newnode;
735  }
736 
737  if (fabs(child2->getValue() - 1.0) < 1.0E-100)
738  {
739  // a*1 -> a
740  delete child2;
741  return child1;
742  }
743 
744  if (child2->getValue() < 0.0)
745  {
746  // a*x -> -(a*|x|) for x < 0
747  CEvaluationNode *newnode = CEvaluationNode::create((Type)(FUNCTION | CEvaluationNodeFunction::MINUS), "-");
748  CEvaluationNode *newchild1 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
749  std::stringstream tmp;
750  tmp.precision(18);
751  tmp << fabs(child2->getValue());
752  CEvaluationNode *grandchild1 = child1;
753  CEvaluationNode *grandchild2 = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), tmp.str());
754  newnode->addChild(newchild1, NULL);
755  newchild1->addChild(grandchild1, NULL);
756  newchild1->addChild(grandchild2, grandchild1);
757  delete child2;
758  return newnode;
759  }
760  }
761 
762  if (child1->getData() == "/")
763  {
764  if (child2->getData() == "/")
765  {
766  // (a/b)*(c/d) -> (a*d)/(b*c)
767  // this rule is actually already executed by the following two rules, but more efficiently this way
768  CEvaluationNode * newnode = CEvaluationNode::create((Type)(OPERATOR | DIVIDE), "/");
769  CEvaluationNode * newchild1 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
770  CEvaluationNode * newchild2 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
771  CEvaluationNode * grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild())->copyBranch();
772  CEvaluationNode * grandchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild())->copyBranch();
773  CEvaluationNode * grandchild3 = dynamic_cast<CEvaluationNode*>(child1->getChild()->getSibling())->copyBranch();
774  CEvaluationNode * grandchild4 = dynamic_cast<CEvaluationNode*>(child2->getChild()->getSibling())->copyBranch();
775  newnode->addChild(newchild1, NULL);
776  newnode->addChild(newchild2, newchild1);
777  newchild1->addChild(grandchild1, NULL);
778  newchild1->addChild(grandchild2, grandchild1);
779  newchild2->addChild(grandchild3, NULL);
780  newchild2->addChild(grandchild4, grandchild3);
781  delete child1;
782  delete child2;
783  return newnode;
784  }
785 
786  // (a/b) * c -> (a*c)/b
787  CEvaluationNode * newnode = CEvaluationNode::create((Type)(OPERATOR | DIVIDE), "/");
788  CEvaluationNode * newchild1 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
789  CEvaluationNode * newchild2 = dynamic_cast<CEvaluationNode*>(child1->getChild()->getSibling())->copyBranch();
790  CEvaluationNode * grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild())->copyBranch();
791  CEvaluationNode * grandchild2 = child2;
792  newnode->addChild(newchild1, NULL);
793  newnode->addChild(newchild2, newchild1);
794  newchild1->addChild(grandchild1, NULL);
795  newchild1->addChild(grandchild2, grandchild1);
796  delete child1;
797  return newnode;
798  }
799 
800  if (child2->getData() == "/")
801  {
802  // a * (b/c) -> (a*b)/c
803  CEvaluationNode * newnode = CEvaluationNode::create((Type)(OPERATOR | DIVIDE), "/");
804  CEvaluationNode * newchild1 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
805  CEvaluationNode * newchild2 = dynamic_cast<CEvaluationNode*>((child2->getChild())->getSibling())->copyBranch();
806  CEvaluationNode * grandchild1 = child1;
807  CEvaluationNode * grandchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild())->copyBranch();
808  newnode->addChild(newchild1, NULL);
809  newnode->addChild(newchild2, newchild1);
810  newchild1->addChild(grandchild1, NULL);
811  newchild1->addChild(grandchild2, grandchild1);
812  delete child2;
813  return newnode;
814  }
815 
816  if (CEvaluationNode::type(child1->getType()) == (Type)(FUNCTION | CEvaluationNodeFunction::MINUS))
817  {
818  // (-a) * b -> -(a*b)
819  // make negativity a property of product
820  CEvaluationNode* newnode = CEvaluationNode::create((Type)(FUNCTION | CEvaluationNodeFunction::MINUS), "-");
821  CEvaluationNode* newchild1 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
822  CEvaluationNode* grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild())->copyBranch();
823  CEvaluationNode* grandchild2 = child2;
824  newnode->addChild(newchild1, NULL);
825  newchild1->addChild(grandchild1, NULL);
826  newchild1->addChild(grandchild2, grandchild1);
827  delete child1;
828  return newnode;
829  }
830 
831  if (CEvaluationNode::type(child2->getType()) == (Type)(FUNCTION | CEvaluationNodeFunction::MINUS))
832  {
833  // a*(-b) -> -(a*b)
834  // make negativity a property of product,
835  // easier to recognize a negative exponent.
836  CEvaluationNode* newnode = CEvaluationNode::create((Type)(FUNCTION | CEvaluationNodeFunction::MINUS), "-");
837  CEvaluationNode* newchild1 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
838  CEvaluationNode* grandchild1 = child1;
839  CEvaluationNode* grandchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild())->copyBranch();
840  newnode->addChild(newchild1, NULL);
841  newchild1->addChild(grandchild1, NULL);
842  newchild1->addChild(grandchild2, grandchild1);
843  delete child2;
844  return newnode;
845  }
846 
847  if (child1->getData() == "+")
848  {
849  // expand (a+b)*c -> a*c + b*c
850  CEvaluationNode * newnode = CEvaluationNode::create((Type)(OPERATOR | PLUS), "+");
851  CEvaluationNode * newchild1 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
852  CEvaluationNode * newchild2 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
853  CEvaluationNode * grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild())->copyBranch();
854  CEvaluationNode * grandchild2 = child2->copyBranch();
855  CEvaluationNode * grandchild3 = dynamic_cast<CEvaluationNode*>(child1->getChild()->getSibling())->copyBranch();
856  CEvaluationNode * grandchild4 = child2;
857  newnode->addChild(newchild1, NULL);
858  newnode->addChild(newchild2, newchild1);
859  newchild1->addChild(grandchild1, NULL);
860  newchild1->addChild(grandchild2, grandchild1);
861  newchild2->addChild(grandchild3, NULL);
862  newchild2->addChild(grandchild4, grandchild3);
863  delete child1;
864  return newnode;
865  }
866 
867  if (child2->getData() == "+")
868  {
869  // expand a*(b+c) -> a*b + a*c
870  CEvaluationNode * newnode = CEvaluationNode::create((Type)(OPERATOR | PLUS), "+");
871  CEvaluationNode * newchild1 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
872  CEvaluationNode * newchild2 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
873  CEvaluationNode * grandchild1 = child1->copyBranch();
874  CEvaluationNode * grandchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild())->copyBranch();
875  CEvaluationNode * grandchild3 = child1;
876  CEvaluationNode * grandchild4 = dynamic_cast<CEvaluationNode*>(child2->getChild()->getSibling())->copyBranch();
877  newnode->addChild(newchild1, NULL);
878  newnode->addChild(newchild2, newchild1);
879  newchild1->addChild(grandchild1, NULL);
880  newchild1->addChild(grandchild2, grandchild1);
881  newchild2->addChild(grandchild3, NULL);
882  newchild2->addChild(grandchild4, grandchild3);
883  delete child2;
884  return newnode;
885  }
886 
887  if (CEvaluationNode::type(child1->getType()) == CEvaluationNode::OPERATOR && ((CEvaluationNodeOperator::SubType)CEvaluationNode::subType(child1->getType())) == CEvaluationNodeOperator::POWER)
888  {
889  // A^n*A^l -> A^(n+l) this way exponents can be simplified
890  // check if the second child is also a power item
891  if (CEvaluationNode::type(child2->getType()) == CEvaluationNode::OPERATOR && ((CEvaluationNodeOperator::SubType)CEvaluationNode::subType(child2->getType())) == CEvaluationNodeOperator::POWER)
892  {
893  if (dynamic_cast<const CEvaluationNode*>(child2->getChild())->buildInfix() == dynamic_cast<const CEvaluationNode*>(child1->getChild())->buildInfix())
894  {
895  CEvaluationNode* newNode = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
896  CEvaluationNode* newchild1 = dynamic_cast<CEvaluationNode*>(child2->getChild())->copyBranch();
897  newNode->addChild(newchild1);
898  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | PLUS), "+");
899  CEvaluationNode* grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild()->getSibling())->copyBranch();
900  CEvaluationNode* grandchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild()->getSibling())->copyBranch();
901  // simplify newchild2
902  std::vector<CEvaluationNode*> children;
903  children.push_back(grandchild1);
904  children.push_back(grandchild2);
905  newNode->addChild(newchild2->simplifyNode(children));
906  // simplify the result again since a power node with an
907  // exponent of 1 or 0 could have been created.
908  delete child1;
909  children.clear();
910  children.push_back(dynamic_cast<CEvaluationNode*>(newNode->getChild())->copyBranch());
911  children.push_back(dynamic_cast<CEvaluationNode*>(newNode->getChild()->getSibling())->copyBranch());
912  child1 = newNode->simplifyNode(children);
913  delete newNode;
914  newNode = child1;
915  delete newchild2;
916  delete child2;
917  return newNode;
918  }
919  }
920  else
921  {
922  // A^n * A -> A^(n+1)
923  // check if the second child is the same as the first child to
924  // the power operator
925  if (child2->buildInfix() == dynamic_cast<const CEvaluationNode*>(child1->getChild())->buildInfix())
926  {
927  CEvaluationNode* newNode = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
928  CEvaluationNode* newchild1 = dynamic_cast<CEvaluationNode*>(child2->getChild())->copyBranch();
929  newNode->addChild(newchild1);
930  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | PLUS), "+");
931  CEvaluationNode* grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild()->getSibling())->copyBranch();
932  CEvaluationNode* grandchild2 = new CEvaluationNodeNumber(CEvaluationNodeNumber::DOUBLE, "1.0");
933  // simplify newchild2
934  std::vector<CEvaluationNode*> children;
935  children.push_back(grandchild1);
936  children.push_back(grandchild2);
937  newNode->addChild(newchild2->simplifyNode(children));
938  // simplify the result again since a power node with an
939  // exponent of 1 or 0 could have been created.
940  delete child1;
941  children.clear();
942  children.push_back(dynamic_cast<CEvaluationNode*>(newNode->getChild())->copyBranch());
943  children.push_back(dynamic_cast<CEvaluationNode*>(newNode->getChild()->getSibling())->copyBranch());
944  child1 = newNode->simplifyNode(children);
945  delete newNode;
946  newNode = child1;
947  delete newchild2;
948  delete child2;
949  return newNode;
950  }
951  }
952  }
953  else if (CEvaluationNode::type(child2->getType()) == CEvaluationNode::OPERATOR && ((CEvaluationNodeOperator::SubType)CEvaluationNode::subType(child2->getType())) == CEvaluationNodeOperator::POWER)
954  {
955  // A*A^n -> A^(n+1)
956  // check if child 1 is the same as the first child to the power
957  // operator
958  if (child1->buildInfix() == dynamic_cast<const CEvaluationNode*>(child2->getChild())->buildInfix())
959  {
960  CEvaluationNode* newNode = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
961  newNode->addChild(child1);
962  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | PLUS), "+");
963  CEvaluationNode* grandchild1 = new CEvaluationNodeNumber(CEvaluationNodeNumber::DOUBLE, "1.0");
964  CEvaluationNode* grandchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild()->getSibling())->copyBranch();
965  // simplify newchild2
966  std::vector<CEvaluationNode*> children;
967  children.push_back(grandchild1);
968  children.push_back(grandchild2);
969  newNode->addChild(newchild2->simplifyNode(children));
970  // simplify the result again since a power node with an
971  // exponent of 1 or 0 could have been created.
972  children.clear();
973  children.push_back(dynamic_cast<CEvaluationNode*>(newNode->getChild())->copyBranch());
974  children.push_back(dynamic_cast<CEvaluationNode*>(newNode->getChild()->getSibling())->copyBranch());
975  child1 = newNode->simplifyNode(children);
976  delete newNode;
977  newNode = child1;
978  delete newchild2;
979  delete child2;
980  return newNode;
981  }
982  }
983 
984  // default: copy
985  CEvaluationNode *newnode = copyNode(child1, child2);
986  return newnode;
987  }
988 
989  case DIVIDE:
990  {
991  if (CEvaluationNode::type(child1->getType()) == NUMBER)
992  {
993  if (CEvaluationNode::type(child2->getType()) == NUMBER)
994  {
995  // both children numbers ->calculate
996  std::stringstream tmp;
997  tmp.precision(18);
998  tmp << child1->getValue() / child2->getValue();
999  CEvaluationNode *newnode = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), tmp.str());
1000  delete child1;
1001  delete child2;
1002  return newnode;
1003  }
1004 
1005  if (fabs(child1->getValue()) < 1.0E-100)
1006  {
1007  // 0/a -> a
1008  CEvaluationNode *newnode = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), "0.0");
1009  delete child1;
1010  delete child2;
1011  return newnode;
1012  }
1013  }
1014 
1015  if (CEvaluationNode::type(child2->getType()) == NUMBER)
1016  {
1017  // a/0 -> false
1018  if (fabs(child2->getValue()) < 1.0E-100)
1019  return false;
1020 
1021  if (fabs(child2->getValue() - 1.0) < 1.0E-100)
1022  {
1023  // a/1 -> a
1024  delete child2;
1025  return child1;
1026  }
1027  }
1028 
1029  if (child1->buildInfix() == child2->buildInfix())
1030  {
1031  // a/a -> 1
1032  CEvaluationNode *newnode = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), "1.0");
1033  delete child1;
1034  delete child2;
1035  return newnode;
1036  }
1037 
1038  if (child1->getData() == "/")
1039  {
1040  if (child2->getData() == "/")
1041  {
1042  // (a/b)/(c/d) -> (a*d)/(b*c)
1043  // this rule is actually already executed by the following two rules, but more efficiently this way
1044  CEvaluationNode* newnode = CEvaluationNode::create((Type)(OPERATOR | DIVIDE), "/");
1045  CEvaluationNode* newchild1 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
1046  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
1047  CEvaluationNode* grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild())->copyBranch();
1048  CEvaluationNode* grandchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild()->getSibling())->copyBranch();
1049  CEvaluationNode* grandchild3 = dynamic_cast<CEvaluationNode*>(child1->getChild()->getSibling())->copyBranch();
1050  CEvaluationNode* grandchild4 = dynamic_cast<CEvaluationNode*>(child2->getChild())->copyBranch();
1051  newnode->addChild(newchild1, NULL);
1052  newnode->addChild(newchild2, newchild1);
1053  newchild1->addChild(grandchild1, NULL);
1054  newchild1->addChild(grandchild2, grandchild1);
1055  newchild2->addChild(grandchild3, NULL);
1056  newchild2->addChild(grandchild4, grandchild3);
1057  delete child1;
1058  delete child2;
1059  return newnode;
1060  }
1061 
1062  // (a/b)/c -> a/(b*c)
1063  CEvaluationNode* newnode = CEvaluationNode::create((Type)(OPERATOR | DIVIDE), "/");
1064  CEvaluationNode* newchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild())->copyBranch();
1065  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
1066  CEvaluationNode* grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild()->getSibling())->copyBranch();
1067  CEvaluationNode* grandchild2 = child2;
1068  newnode->addChild(newchild1, NULL);
1069  newnode->addChild(newchild2, newchild1);
1070  newchild2->addChild(grandchild1, NULL);
1071  newchild2->addChild(grandchild2, grandchild1);
1072  delete child1;
1073  return newnode;
1074  }
1075 
1076  if (child2->getData() == "/")
1077  {
1078  // a/(b/c) -> (a*c)/b
1079  CEvaluationNode* newnode = CEvaluationNode::create((Type)(OPERATOR | DIVIDE), "/");
1080  CEvaluationNode* newchild1 = CEvaluationNode::create((Type)(OPERATOR | MULTIPLY), "*");
1081  CEvaluationNode* newchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild())->copyBranch();
1082  CEvaluationNode* grandchild1 = child1;
1083  CEvaluationNode* grandchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild()->getSibling())->copyBranch();
1084  newnode->addChild(newchild1, NULL);
1085  newnode->addChild(newchild2, newchild1);
1086  newchild1->addChild(grandchild1, NULL);
1087  newchild1->addChild(grandchild2, grandchild1);
1088  delete child2;
1089  return newnode;
1090  }
1091 
1092  if (CEvaluationNode::type(child1->getType()) == CEvaluationNode::OPERATOR && ((CEvaluationNodeOperator::SubType)CEvaluationNode::subType(child1->getType())) == CEvaluationNodeOperator::POWER)
1093  {
1094  // A^n / A^l -> A^(n-l)
1095  // check if the second child is also a power item
1096  if (CEvaluationNode::type(child2->getType()) == CEvaluationNode::OPERATOR && ((CEvaluationNodeOperator::SubType)CEvaluationNode::subType(child2->getType())) == CEvaluationNodeOperator::POWER)
1097  {
1098  if (dynamic_cast<const CEvaluationNode*>(child2->getChild())->buildInfix() == dynamic_cast<const CEvaluationNode*>(child1->getChild())->buildInfix())
1099  {
1100  CEvaluationNode* newNode = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
1101  CEvaluationNode* newchild1 = dynamic_cast<CEvaluationNode*>(child2->getChild())->copyBranch();
1102  newNode->addChild(newchild1);
1103  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | MINUS), "-");
1104  CEvaluationNode* grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild()->getSibling())->copyBranch();
1105  CEvaluationNode* grandchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild()->getSibling())->copyBranch();
1106  // simplify newchild2
1107  std::vector<CEvaluationNode*> children;
1108  children.push_back(grandchild1);
1109  children.push_back(grandchild2);
1110  newNode->addChild(newchild2->simplifyNode(children));
1111  // simplify the result again since a power node with an
1112  // exponent of 1 or 0 could have been created.
1113  delete child1;
1114  children.clear();
1115  children.push_back(dynamic_cast<CEvaluationNode*>(newNode->getChild())->copyBranch());
1116  children.push_back(dynamic_cast<CEvaluationNode*>(newNode->getChild()->getSibling())->copyBranch());
1117  child1 = newNode->simplifyNode(children);
1118  delete newNode;
1119  newNode = child1;
1120  delete newchild2;
1121  delete child2;
1122  return newNode;
1123  }
1124  }
1125  else
1126  {
1127  // A^n/A -> A^(n-1)
1128  // check if the second child is the same as the first child to
1129  // the power operator
1130  if (child2->buildInfix() == dynamic_cast<const CEvaluationNode*>(child1->getChild())->buildInfix())
1131  {
1132  CEvaluationNode* newNode = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
1133  CEvaluationNode* newchild1 = dynamic_cast<CEvaluationNode*>(child2->getChild())->copyBranch();
1134  newNode->addChild(newchild1);
1135  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | MINUS), "-");
1136  CEvaluationNode* grandchild1 = dynamic_cast<CEvaluationNode*>(child1->getChild()->getSibling())->copyBranch();
1137  CEvaluationNode* grandchild2 = new CEvaluationNodeNumber(CEvaluationNodeNumber::DOUBLE, "1.0");
1138  // simplify newchild2
1139  std::vector<CEvaluationNode*> children;
1140  children.push_back(grandchild1);
1141  children.push_back(grandchild2);
1142  newNode->addChild(newchild2->simplifyNode(children));
1143  // simplify the result again since a power node with an
1144  // exponent of 1 or 0 could have been created.
1145  delete child1;
1146  children.clear();
1147  children.push_back(dynamic_cast<CEvaluationNode*>(newNode->getChild())->copyBranch());
1148  children.push_back(dynamic_cast<CEvaluationNode*>(newNode->getChild()->getSibling())->copyBranch());
1149  child1 = newNode->simplifyNode(children);
1150  delete newNode;
1151  newNode = child1;
1152  delete newchild2;
1153  delete child2;
1154  return newNode;
1155  }
1156  }
1157  }
1158  else if (CEvaluationNode::type(child2->getType()) == CEvaluationNode::OPERATOR && ((CEvaluationNodeOperator::SubType)CEvaluationNode::subType(child2->getType())) == CEvaluationNodeOperator::POWER)
1159  {
1160  // A / A^n -> A^(1-n)
1161  // check if child 1 is the same as the first child to the power
1162  // operator
1163  if (child1->buildInfix() == dynamic_cast<const CEvaluationNode*>(child2->getChild())->buildInfix())
1164  {
1165  CEvaluationNode* newNode = CEvaluationNode::create((Type)(OPERATOR | POWER), "^");
1166  newNode->addChild(child1);
1167  CEvaluationNode* newchild2 = CEvaluationNode::create((Type)(OPERATOR | MINUS), "-");
1168  CEvaluationNode* grandchild1 = new CEvaluationNodeNumber(CEvaluationNodeNumber::DOUBLE, "1.0");
1169  CEvaluationNode* grandchild2 = dynamic_cast<CEvaluationNode*>(child2->getChild()->getSibling())->copyBranch();
1170  // simplify newchild2
1171  std::vector<CEvaluationNode*> children;
1172  children.push_back(grandchild1);
1173  children.push_back(grandchild2);
1174  newNode->addChild(newchild2->simplifyNode(children));
1175  // simplify the result again since a power node with an
1176  // exponent of 1 or 0 could have been created.
1177  children.clear();
1178  children.push_back(dynamic_cast<CEvaluationNode*>(newNode->getChild())->copyBranch());
1179  children.push_back(dynamic_cast<CEvaluationNode*>(newNode->getChild()->getSibling())->copyBranch());
1180  child1 = newNode->simplifyNode(children);
1181  delete newNode;
1182  newNode = child1;
1183  delete newchild2;
1184  delete child2;
1185  return newNode;
1186  }
1187  }
1188 
1189  // default: copy
1190  CEvaluationNode *newnode = copyNode(child1, child2);
1191  return newnode;
1192  }
1193 
1194  case PLUS:
1195  {
1196  if (CEvaluationNode::type(child1->getType()) == NUMBER)
1197  {
1198  if (CEvaluationNode::type(child2->getType()) == NUMBER)
1199  {
1200  // both children numbers ->calculate
1201  std::stringstream tmp;
1202  tmp.precision(18);
1203  tmp << child1->getValue() + child2->getValue();
1204  CEvaluationNode *newnode = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), tmp.str());
1205  delete child1;
1206  delete child2;
1207  return newnode;
1208  }
1209 
1210  if (fabs(child1->getValue()) < 1.0E-100)
1211  {
1212  // 0 + a -> a
1213  delete child1;
1214  return child2;
1215  }
1216  }
1217 
1218  if ((CEvaluationNode::type(child2->getType()) == NUMBER) && (fabs(child2->getValue()) < 1.0E-100))
1219  {
1220  // a + 0 -> a
1221  delete child2;
1222  return child1;
1223  }
1224 
1225  if (((child1->getType() == (Type)(FUNCTION | CEvaluationNodeFunction::MINUS))
1226  && (dynamic_cast<CEvaluationNode*>(child1->getChild())->buildInfix() == child2->buildInfix()))
1227  ||
1228  ((child2->getType() == (Type)(FUNCTION | CEvaluationNodeFunction::MINUS))
1229  && (dynamic_cast<CEvaluationNode*>(child2->getChild())->buildInfix() == child1->buildInfix())))
1230  {
1231  // -(a) + a and a + (-a) -> 0
1232  CEvaluationNode *newnode = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), "0.0");
1233  delete child1;
1234  delete child2;
1235  return newnode;
1236  }
1237 
1238  // default: copy
1239  CEvaluationNode *newnode = copyNode(child1, child2);
1240  return newnode;
1241  }
1242 
1243  case MINUS:
1244  {
1245  if (CEvaluationNode::type(child1->getType()) == NUMBER)
1246  {
1247  if (CEvaluationNode::type(child2->getType()) == NUMBER)
1248  {
1249  // both children numbers ->calculate
1250  std::stringstream tmp;
1251  tmp.precision(18);
1252  tmp << child1->getValue() - child2->getValue();
1253  CEvaluationNode *newnode = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), tmp.str());
1254  delete child1;
1255  delete child2;
1256  return newnode;
1257  }
1258 
1259  if (fabs(child1->getValue()) < 1.0E-100)
1260  {
1261  // 0-a -> -(a)
1262  CEvaluationNode *newnode = CEvaluationNode::create((Type)(FUNCTION | CEvaluationNodeFunction::MINUS), "-");
1263  newnode->addChild(child2, NULL);
1264  delete child1;
1265  return newnode;
1266  }
1267  }
1268 
1269  if ((CEvaluationNode::type(child2->getType()) == NUMBER) && (fabs(child2->getValue()) < 1.0E-100))
1270  {
1271  // a-0 -> a
1272  delete child2;
1273  return child1;
1274  }
1275 
1276  if (child1->buildInfix() == child2->buildInfix())
1277  {
1278  // a-a -> 0
1279  CEvaluationNode *newnode = CEvaluationNode::create((Type)(NUMBER | CEvaluationNodeNumber::DOUBLE), "0.0");
1280  delete child1;
1281  delete child2;
1282  return newnode;
1283  }
1284 
1285  // default: a - b -> a + (-b)
1286  CEvaluationNode *newnode = CEvaluationNode::create((Type)(OPERATOR | PLUS), "+");
1287  CEvaluationNode *newchild2 = CEvaluationNode::create((Type)(FUNCTION | CEvaluationNodeFunction::MINUS), "-");
1288  newnode->addChild(child1, NULL);
1289  newnode->addChild(newchild2, child1);
1290  newchild2->addChild(child2, NULL);
1291  return newnode;
1292  }
1293 
1294  default: //case MODULUS
1295  {
1296  CEvaluationNode *newnode = copyNode(child1, child2);
1297  return newnode;
1298  }
1299  }
1300 }
1301 
1302 bool CEvaluationNodeOperator::createModuloTree(const CEvaluationNodeOperator* pNode, ASTNode* pASTNode, const CCopasiDataModel* pDataModel) const
1303 {
1304  bool result = false;
1305 
1306  if ((SubType)CEvaluationNode::subType(pNode->getType()) == MODULUS)
1307  {
1308  // the node has two children x and y
1309  const CEvaluationNode* x = dynamic_cast<const CEvaluationNode*>(pNode->getChild());
1310 
1311  if (x != NULL)
1312  {
1313  const CEvaluationNode* y = dynamic_cast<const CEvaluationNode*>(x->getSibling());
1314 
1315  if (y != NULL)
1316  {
1317  // Frank noticed that this should actually be implemented as a
1318  // piecewise function because if one of the arguments is
1319  // negative, the definition is different
1320  pASTNode->setType(AST_FUNCTION_PIECEWISE);
1321  // the first child is the true branch
1322  // x%y -> x-ceil(x/y)*y
1323  ASTNode* pASTNodeTrue = new ASTNode();
1324  pASTNodeTrue->setType(AST_MINUS);
1325  ASTNode* tmpASTNode = new ASTNode(AST_DIVIDE);
1326  tmpASTNode->addChild(x->toAST(pDataModel));
1327  tmpASTNode->addChild(y->toAST(pDataModel));
1328  ASTNode* tmpASTNode2 = new ASTNode(AST_FUNCTION_CEILING);
1329  tmpASTNode2->addChild(tmpASTNode);
1330  tmpASTNode = new ASTNode(AST_TIMES);
1331  tmpASTNode->addChild(y->toAST(pDataModel));
1332  tmpASTNode->addChild(tmpASTNode2);
1333  pASTNodeTrue->addChild(x->toAST(pDataModel));
1334  pASTNodeTrue->addChild(tmpASTNode);
1335  pASTNode->addChild(pASTNodeTrue);
1336  // now comes the condition
1337  // if exactly one of the arguments to modulo is a negative number
1338  // we use the ceil branch, else we use the floor branch
1339  // x < 0 xor y < 0
1340  // xor
1341  ASTNode* pASTNodeCondition = new ASTNode();
1342  pASTNodeCondition->setType(AST_LOGICAL_XOR);
1343  // x < 0
1344  // <
1345  tmpASTNode = new ASTNode(AST_RELATIONAL_LT);
1346  // x
1347  tmpASTNode->addChild(x->toAST(pDataModel));
1348  // 0
1349  tmpASTNode2 = new ASTNode(AST_INTEGER);
1350  tmpASTNode2->setValue(0);
1351  tmpASTNode->addChild(tmpASTNode2);
1352  pASTNodeCondition->addChild(tmpASTNode);
1353  // y < 0
1354  // <
1355  tmpASTNode = new ASTNode(AST_RELATIONAL_LT);
1356  // y
1357  tmpASTNode->addChild(y->toAST(pDataModel));
1358  // 0
1359  tmpASTNode2 = new ASTNode(AST_INTEGER);
1360  tmpASTNode2->setValue(0);
1361  tmpASTNode->addChild(tmpASTNode2);
1362  pASTNodeCondition->addChild(tmpASTNode);
1363  pASTNode->addChild(pASTNodeCondition);
1364  // last is the false branch
1365  // x%y -> x-floor(x/y)*y
1366  ASTNode* pASTNodeFalse = new ASTNode();
1367  pASTNodeFalse->setType(AST_MINUS);
1368  tmpASTNode = new ASTNode(AST_DIVIDE);
1369  tmpASTNode->addChild(x->toAST(pDataModel));
1370  tmpASTNode->addChild(y->toAST(pDataModel));
1371  tmpASTNode2 = new ASTNode(AST_FUNCTION_FLOOR);
1372  tmpASTNode2->addChild(tmpASTNode);
1373  tmpASTNode = new ASTNode(AST_TIMES);
1374  tmpASTNode->addChild(y->toAST(pDataModel));
1375  tmpASTNode->addChild(tmpASTNode2);
1376  pASTNodeFalse->addChild(x->toAST(pDataModel));
1377  pASTNodeFalse->addChild(tmpASTNode);
1378  pASTNode->addChild(pASTNodeFalse);
1379  result = true;
1380  }
1381  }
1382  }
1383 
1384  return result;
1385 }
1386 
1387 CEvaluationNode * CEvaluationNodeOperator::getLeft()
1388 {return mpLeft;}
1389 const CEvaluationNode * CEvaluationNodeOperator::getLeft() const
1390 {return mpLeft;}
1391 CEvaluationNode * CEvaluationNodeOperator::getRight()
1392 {return mpRight;}
1393 const CEvaluationNode * CEvaluationNodeOperator::getRight() const
1394 {return mpRight;}
1395 
1396 #include "utilities/copasimathml.h"
1397 
1398 // virtual
1399 std::string CEvaluationNodeOperator::getMMLString(const std::vector< std::string > & children,
1400  bool expand,
1401  const std::vector< std::vector< std::string > > & /* variables */) const
1402 {
1403  std::ostringstream out;
1404 
1405  bool flag;
1406 
1407  switch (mType & 0x00FFFFFF)
1408  {
1409  case PLUS:
1410 
1411  out << "<mrow>" << std::endl;
1412  out << children[0];
1413  out << "<mo>" << "+" << "</mo>" << std::endl;
1414  out << children[1];
1415  out << "</mrow>" << std::endl;
1416  break;
1417 
1418  case MINUS:
1419  out << "<mrow>" << std::endl;
1420  out << children[0];
1421  out << "<mo>" << "-" << "</mo>" << std::endl;
1422 
1423  flag = ((mpRight->getType() == (CEvaluationNode::OPERATOR | PLUS))
1424  || (mpRight->getType() == (CEvaluationNode::OPERATOR | MINUS))
1425  || (((mpRight->getType() & 0xFF000000) == CEvaluationNode::CALL) && expand)
1426  );
1427 
1428  if (flag) out << "<mfenced>" << std::endl;
1429 
1430  out << children[1];
1431 
1432  if (flag) out << "</mfenced>" << std::endl; // ???
1433 
1434  out << "</mrow>" << std::endl;
1435  break;
1436 
1437  case MULTIPLY:
1438  out << "<mrow>" << std::endl;
1439 
1440  //do we need "()" ?
1441  flag = ((mpLeft->getType() == (CEvaluationNode::OPERATOR | PLUS))
1442  || (mpLeft->getType() == (CEvaluationNode::OPERATOR | MINUS))
1443  || (((mpLeft->getType() & 0xFF000000) == CEvaluationNode::CALL) && expand)
1444  );
1445 
1446  if (flag) out << "<mfenced>" << std::endl;
1447 
1448  out << children[0];
1449 
1450  if (flag) out << "</mfenced>" << std::endl;
1451 
1452  out << "<mo>" << "&CenterDot;" << "</mo>" << std::endl;
1453 
1454  flag = ((mpRight->getType() == (CEvaluationNode::OPERATOR | PLUS))
1455  || (mpRight->getType() == (CEvaluationNode::OPERATOR | MINUS))
1456  || (((mpRight->getType() & 0xFF000000) == CEvaluationNode::CALL) && expand)
1457  );
1458 
1459  if (flag) out << "<mfenced>" << std::endl;
1460 
1461  out << children[1];
1462 
1463  if (flag) out << "</mfenced>" << std::endl;
1464 
1465  out << "</mrow>" << std::endl;
1466 
1467  break;
1468 
1469  case DIVIDE:
1470  out << "<mfrac>" << std::endl;
1471 
1472  //out << "<mrow>" << std::endl;
1473  out << children[0];
1474  //out << "</mrow>" << std::endl;
1475 
1476  //out << "<mrow>" << std::endl;
1477  out << children[1];
1478  //out << "</mrow>" << std::endl;
1479 
1480  out << "</mfrac>" << std::endl;
1481  break;
1482 
1483  case POWER:
1484  out << "<msup>" << std::endl;
1485 
1486  //do we need "()" ?
1487  flag = ((mpLeft->getType() == (CEvaluationNode::OPERATOR | PLUS))
1488  || (mpLeft->getType() == (CEvaluationNode::OPERATOR | MINUS))
1489  || (mpLeft->getType() == (CEvaluationNode::OPERATOR | MULTIPLY))
1490  || (mpLeft->getType() == (CEvaluationNode::OPERATOR | DIVIDE))
1491  || (mpLeft->getType() == (CEvaluationNode::OPERATOR | POWER))
1492  || (((mpLeft->getType() & 0xFF000000) == CEvaluationNode::CALL) && expand)
1493  );
1494 
1495  if (flag) out << "<mfenced>" << std::endl;
1496 
1497  out << children[0];
1498 
1499  if (flag) out << "</mfenced>" << std::endl;
1500 
1501  out << "<mrow>" << std::endl;
1502  out << children[1];
1503  out << "</mrow>" << std::endl;
1504 
1505  out << "</msup>" << std::endl;
1506  break;
1507 
1508  case MODULUS:
1509  out << "<mrow>" << std::endl;
1510 
1511  //do we need "()" ?
1512  flag = true;
1513 
1514  if (flag) out << "<mfenced>" << std::endl;
1515 
1516  out << children[0];
1517 
1518  if (flag) out << "</mfenced>" << std::endl;
1519 
1520  out << "<mo>" << "%" << "</mo>" << std::endl;
1521 
1522  flag = true;
1523 
1524  if (flag) out << "<mfenced>" << std::endl;
1525 
1526  out << children[1];
1527 
1528  if (flag) out << "</mfenced>" << std::endl;
1529 
1530  out << "</mrow>" << std::endl;
1531  break;
1532  }
1533 
1534  return out.str();
1535 }
CEvaluationNodeOperator::getDisplayString
virtual std::string getDisplayString(const std::vector< std::string > &children) const
Definition: CEvaluationNodeOperator.cpp:111
CEvaluationNodeOperator::getBerkeleyMadonnaString
virtual std::string getBerkeleyMadonnaString(const std::vector< std::string > &children) const
Definition: CEvaluationNodeOperator.cpp:184
CEvaluationNode::toAST
virtual ASTNode * toAST(const CCopasiDataModel *pDataModel) const
Definition: CEvaluationNode.cpp:436
CEvaluationNode::copyBranch
CEvaluationNode * copyBranch() const
Definition: CEvaluationNode.cpp:404
MCMathML
#define MCMathML
Definition: CCopasiMessage.h:50
CCopasiNode::getSibling
CCopasiNode< Data > * getSibling()
Definition: CCopasiNode.h:353
CEvaluationNodeOperator::createModuloTree
bool createModuloTree(const CEvaluationNodeOperator *pNode, ASTNode *pASTNode, const CCopasiDataModel *pDataModel) const
Definition: CEvaluationNodeOperator.cpp:1302
CEvaluationNodeOperator::getMMLString
virtual std::string getMMLString(const std::vector< std::string > &children, bool expand, const std::vector< std::vector< std::string > > &variables) const
Definition: CEvaluationNodeOperator.cpp:1399
CEvaluationNode::copyNode
CEvaluationNode * copyNode(CEvaluationNode *child1, CEvaluationNode *child2) const
Definition: CEvaluationNode.cpp:378
CEvaluationNodeOperator::compile
virtual bool compile(const CEvaluationTree *pTree)
Definition: CEvaluationNodeOperator.cpp:72
CEvaluationNode::getValue
const C_FLOAT64 & getValue() const
Definition: CEvaluationNode.h:144
fatalError
#define fatalError()
Definition: CCopasiMessage.h:74
CEvaluationNode::getType
const Type & getType() const
Definition: CEvaluationNode.cpp:357
CEvaluationNodeOperator::toAST
virtual ASTNode * toAST(const CCopasiDataModel *pDataModel) const
Definition: CEvaluationNodeOperator.cpp:404
CEvaluationNode::create
static CEvaluationNode * create(const Type &type, const Data &data)
Definition: CEvaluationNode.cpp:37
CEvaluationNode::buildInfix
std::string buildInfix() const
Definition: CEvaluationNode.cpp:223
CEvaluationNodeOperator::getXPPString
virtual std::string getXPPString(const std::vector< std::string > &children) const
Definition: CEvaluationNodeOperator.cpp:216
PRECEDENCE_OPERATOR_MODULUS
#define PRECEDENCE_OPERATOR_MODULUS
Definition: CEvaluationNode.h:442
PRECEDENCE_OPERATOR_PLUS
#define PRECEDENCE_OPERATOR_PLUS
Definition: CEvaluationNode.h:441
CEvaluationNode::type
static Type type(const Type &type)
Definition: CEvaluationNode.cpp:123
CEvaluationNodeOperator::getCCodeString
virtual std::string getCCodeString(const std::vector< std::string > &children) const
Definition: CEvaluationNodeOperator.cpp:136
PRECEDENCE_OPERATOR_DIVIDE
#define PRECEDENCE_OPERATOR_DIVIDE
Definition: CEvaluationNode.h:443
CCopasiNode::addChild
virtual bool addChild(CCopasiNode< Data > *pChild, CCopasiNode< Data > *pAfter=NULL)
Definition: CCopasiNode.h:156
CEvaluationNodeOperator::simplifyNode
virtual CEvaluationNode * simplifyNode(const std::vector< CEvaluationNode * > &children) const
Definition: CEvaluationNodeOperator.cpp:462
flag
long int flag
Definition: f2c.h:52
CEvaluationNode::subType
static Type subType(const Type &type)
Definition: CEvaluationNode.cpp:120
PRECEDENCE_OPERATOR_POWER
#define PRECEDENCE_OPERATOR_POWER
Definition: CEvaluationNode.h:445
PRECEDENCE_OPERATOR_MINUS
#define PRECEDENCE_OPERATOR_MINUS
Definition: CEvaluationNode.h:440
CEvaluationNode::mPrecedence
class CEvaluationNode::CPrecedence mPrecedence
CEvaluationNodeOperator::getInfix
virtual std::string getInfix(const std::vector< std::string > &children) const
Definition: CEvaluationNodeOperator.cpp:86
PRECEDENCE_OPERATOR_MULTIPLY
#define PRECEDENCE_OPERATOR_MULTIPLY
Definition: CEvaluationNode.h:444
CEvaluationNode::simplifyNode
virtual CEvaluationNode * simplifyNode(const std::vector< CEvaluationNode * > &children) const
Definition: CEvaluationNode.cpp:430
CCopasiNode::getData
virtual const Data & getData() const
Definition: CCopasiNode.h:118
CCopasiNode< std::string >::getChild
CCopasiNode< Data > * getChild()
Definition: CCopasiNode.h:210
CEvaluationNodeOperator::fromAST
static CEvaluationNode * fromAST(const ASTNode *pASTNode, const std::vector< CEvaluationNode * > &children)
Definition: CEvaluationNodeOperator.cpp:257
Generated on Thu Jul 2 2015 14:50:12 for COPASI API by   doxygen 1.8.6