source: trunk/sources/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/SymbolicDataAnalysisExpressionTreeSimplifier.cs @ 14826

Last change on this file since 14826 was 14826, checked in by gkronber, 6 months ago

#2650: merged the factors branch into trunk

File size: 57.7 KB
Line 
1#region License Information
2
3/* HeuristicLab
4 * Copyright (C) 2002-2016 Heuristic and Evolutionary Algorithms Laboratory (HEAL)
5 *
6 * This file is part of HeuristicLab.
7 *
8 * HeuristicLab is free software: you can redistribute it and/or modify
9 * it under the terms of the GNU General Public License as published by
10 * the Free Software Foundation, either version 3 of the License, or
11 * (at your option) any later version.
12 *
13 * HeuristicLab is distributed in the hope that it will be useful,
14 * but WITHOUT ANY WARRANTY; without even the implied warranty of
15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
16 * GNU General Public License for more details.
17 *
18 * You should have received a copy of the GNU General Public License
19 * along with HeuristicLab. If not, see <http://www.gnu.org/licenses/>.
20 */
21
22#endregion
23
24using System;
25using System.Collections.Generic;
26using System.Diagnostics;
27using System.Linq;
28using HeuristicLab.Common;
29using HeuristicLab.Core;
30using HeuristicLab.Encodings.SymbolicExpressionTreeEncoding;
31
32namespace HeuristicLab.Problems.DataAnalysis.Symbolic {
33  /// <summary>
34  /// Simplifier for symbolic expressions
35  /// </summary>
36  public class SymbolicDataAnalysisExpressionTreeSimplifier {
37    private Addition addSymbol = new Addition();
38    private Multiplication mulSymbol = new Multiplication();
39    private Division divSymbol = new Division();
40    private Constant constSymbol = new Constant();
41    private Variable varSymbol = new Variable();
42    private Logarithm logSymbol = new Logarithm();
43    private Exponential expSymbol = new Exponential();
44    private Root rootSymbol = new Root();
45    private Square sqrSymbol = new Square();
46    private SquareRoot sqrtSymbol = new SquareRoot();
47    private Power powSymbol = new Power();
48    private Sine sineSymbol = new Sine();
49    private Cosine cosineSymbol = new Cosine();
50    private Tangent tanSymbol = new Tangent();
51    private IfThenElse ifThenElseSymbol = new IfThenElse();
52    private And andSymbol = new And();
53    private Or orSymbol = new Or();
54    private Not notSymbol = new Not();
55    private GreaterThan gtSymbol = new GreaterThan();
56    private LessThan ltSymbol = new LessThan();
57    private Integral integralSymbol = new Integral();
58    private LaggedVariable laggedVariableSymbol = new LaggedVariable();
59    private TimeLag timeLagSymbol = new TimeLag();
60
61    public ISymbolicExpressionTree Simplify(ISymbolicExpressionTree originalTree) {
62      var clone = (ISymbolicExpressionTreeNode)originalTree.Root.Clone();
63      // macro expand (initially no argument trees)
64      var macroExpandedTree = MacroExpand(clone, clone.GetSubtree(0), new List<ISymbolicExpressionTreeNode>());
65      ISymbolicExpressionTreeNode rootNode = (new ProgramRootSymbol()).CreateTreeNode();
66      rootNode.AddSubtree(GetSimplifiedTree(macroExpandedTree));
67
68#if DEBUG
69      // check that each node is only referenced once
70      var nodes = rootNode.IterateNodesPrefix().ToArray();
71      foreach(var n in nodes) if(nodes.Count(ni => ni == n) > 1) throw new InvalidOperationException();
72#endif
73      return new SymbolicExpressionTree(rootNode);
74    }
75
76    // the argumentTrees list contains already expanded trees used as arguments for invocations
77    private ISymbolicExpressionTreeNode MacroExpand(ISymbolicExpressionTreeNode root, ISymbolicExpressionTreeNode node,
78      IList<ISymbolicExpressionTreeNode> argumentTrees) {
79      List<ISymbolicExpressionTreeNode> subtrees = new List<ISymbolicExpressionTreeNode>(node.Subtrees);
80      while (node.SubtreeCount > 0) node.RemoveSubtree(0);
81      if (node.Symbol is InvokeFunction) {
82        var invokeSym = node.Symbol as InvokeFunction;
83        var defunNode = FindFunctionDefinition(root, invokeSym.FunctionName);
84        var macroExpandedArguments = new List<ISymbolicExpressionTreeNode>();
85        foreach (var subtree in subtrees) {
86          macroExpandedArguments.Add(MacroExpand(root, subtree, argumentTrees));
87        }
88        return MacroExpand(root, defunNode, macroExpandedArguments);
89      } else if (node.Symbol is Argument) {
90        var argSym = node.Symbol as Argument;
91        // return the correct argument sub-tree (already macro-expanded)
92        return (SymbolicExpressionTreeNode)argumentTrees[argSym.ArgumentIndex].Clone();
93      } else {
94        // recursive application
95        foreach (var subtree in subtrees) {
96          node.AddSubtree(MacroExpand(root, subtree, argumentTrees));
97        }
98        return node;
99      }
100    }
101
102    private ISymbolicExpressionTreeNode FindFunctionDefinition(ISymbolicExpressionTreeNode root, string functionName) {
103      foreach (var subtree in root.Subtrees.OfType<DefunTreeNode>()) {
104        if (subtree.FunctionName == functionName) return subtree.GetSubtree(0);
105      }
106
107      throw new ArgumentException("Definition of function " + functionName + " not found.");
108    }
109
110    #region symbol predicates
111
112    // arithmetic
113    private bool IsDivision(ISymbolicExpressionTreeNode node) {
114      return node.Symbol is Division;
115    }
116
117    private bool IsMultiplication(ISymbolicExpressionTreeNode node) {
118      return node.Symbol is Multiplication;
119    }
120
121    private bool IsSubtraction(ISymbolicExpressionTreeNode node) {
122      return node.Symbol is Subtraction;
123    }
124
125    private bool IsAddition(ISymbolicExpressionTreeNode node) {
126      return node.Symbol is Addition;
127    }
128
129    private bool IsAverage(ISymbolicExpressionTreeNode node) {
130      return node.Symbol is Average;
131    }
132
133    // exponential
134    private bool IsLog(ISymbolicExpressionTreeNode node) {
135      return node.Symbol is Logarithm;
136    }
137
138    private bool IsExp(ISymbolicExpressionTreeNode node) {
139      return node.Symbol is Exponential;
140    }
141
142    private bool IsRoot(ISymbolicExpressionTreeNode node) {
143      return node.Symbol is Root;
144    }
145
146    private bool IsSquare(ISymbolicExpressionTreeNode node) {
147      return node.Symbol is Square;
148    }
149
150    private bool IsSquareRoot(ISymbolicExpressionTreeNode node) {
151      return node.Symbol is SquareRoot;
152    }
153
154    private bool IsPower(ISymbolicExpressionTreeNode node) {
155      return node.Symbol is Power;
156    }
157
158    // trigonometric
159    private bool IsSine(ISymbolicExpressionTreeNode node) {
160      return node.Symbol is Sine;
161    }
162
163    private bool IsCosine(ISymbolicExpressionTreeNode node) {
164      return node.Symbol is Cosine;
165    }
166
167    private bool IsTangent(ISymbolicExpressionTreeNode node) {
168      return node.Symbol is Tangent;
169    }
170
171    // boolean
172    private bool IsIfThenElse(ISymbolicExpressionTreeNode node) {
173      return node.Symbol is IfThenElse;
174    }
175
176    private bool IsAnd(ISymbolicExpressionTreeNode node) {
177      return node.Symbol is And;
178    }
179
180    private bool IsOr(ISymbolicExpressionTreeNode node) {
181      return node.Symbol is Or;
182    }
183
184    private bool IsNot(ISymbolicExpressionTreeNode node) {
185      return node.Symbol is Not;
186    }
187
188    // comparison
189    private bool IsGreaterThan(ISymbolicExpressionTreeNode node) {
190      return node.Symbol is GreaterThan;
191    }
192
193    private bool IsLessThan(ISymbolicExpressionTreeNode node) {
194      return node.Symbol is LessThan;
195    }
196
197    private bool IsBoolean(ISymbolicExpressionTreeNode node) {
198      return
199        node.Symbol is GreaterThan ||
200        node.Symbol is LessThan ||
201        node.Symbol is And ||
202        node.Symbol is Or;
203    }
204
205    // terminals
206    private bool IsVariable(ISymbolicExpressionTreeNode node) {
207      return node.Symbol is Variable;
208    }
209
210    private bool IsVariableBase(ISymbolicExpressionTreeNode node) {
211      return node is VariableTreeNodeBase;
212    }
213
214    private bool IsFactor(ISymbolicExpressionTreeNode node) {
215      return node is FactorVariableTreeNode;
216    }
217
218    private bool IsBinFactor(ISymbolicExpressionTreeNode node) {
219      return node is BinaryFactorVariableTreeNode;
220    }
221
222    private bool IsConstant(ISymbolicExpressionTreeNode node) {
223      return node.Symbol is Constant;
224    }
225
226    // dynamic
227    private bool IsTimeLag(ISymbolicExpressionTreeNode node) {
228      return node.Symbol is TimeLag;
229    }
230
231    private bool IsIntegral(ISymbolicExpressionTreeNode node) {
232      return node.Symbol is Integral;
233    }
234
235    #endregion
236
237    /// <summary>
238    /// Creates a new simplified tree
239    /// </summary>
240    /// <param name="original"></param>
241    /// <returns></returns>
242    public ISymbolicExpressionTreeNode GetSimplifiedTree(ISymbolicExpressionTreeNode original) {
243      if (IsConstant(original) || IsVariableBase(original)) {
244        return (ISymbolicExpressionTreeNode)original.Clone();
245      } else if (IsAddition(original)) {
246        return SimplifyAddition(original);
247      } else if (IsSubtraction(original)) {
248        return SimplifySubtraction(original);
249      } else if (IsMultiplication(original)) {
250        return SimplifyMultiplication(original);
251      } else if (IsDivision(original)) {
252        return SimplifyDivision(original);
253      } else if (IsAverage(original)) {
254        return SimplifyAverage(original);
255      } else if (IsLog(original)) {
256        return SimplifyLog(original);
257      } else if (IsExp(original)) {
258        return SimplifyExp(original);
259      } else if (IsSquare(original)) {
260        return SimplifySquare(original);
261      } else if (IsSquareRoot(original)) {
262        return SimplifySquareRoot(original);
263      } else if (IsPower(original)) {
264        return SimplifyPower(original);
265      } else if (IsRoot(original)) {
266        return SimplifyRoot(original);
267      } else if (IsSine(original)) {
268        return SimplifySine(original);
269      } else if (IsCosine(original)) {
270        return SimplifyCosine(original);
271      } else if (IsTangent(original)) {
272        return SimplifyTangent(original);
273      } else if (IsIfThenElse(original)) {
274        return SimplifyIfThenElse(original);
275      } else if (IsGreaterThan(original)) {
276        return SimplifyGreaterThan(original);
277      } else if (IsLessThan(original)) {
278        return SimplifyLessThan(original);
279      } else if (IsAnd(original)) {
280        return SimplifyAnd(original);
281      } else if (IsOr(original)) {
282        return SimplifyOr(original);
283      } else if (IsNot(original)) {
284        return SimplifyNot(original);
285      } else if (IsTimeLag(original)) {
286        return SimplifyTimeLag(original);
287      } else if (IsIntegral(original)) {
288        return SimplifyIntegral(original);
289      } else {
290        return SimplifyAny(original);
291      }
292    }
293
294    #region specific simplification routines
295
296    private ISymbolicExpressionTreeNode SimplifyAny(ISymbolicExpressionTreeNode original) {
297      // can't simplify this function but simplify all subtrees
298      List<ISymbolicExpressionTreeNode> subtrees = new List<ISymbolicExpressionTreeNode>(original.Subtrees);
299      while (original.Subtrees.Count() > 0) original.RemoveSubtree(0);
300      var clone = (SymbolicExpressionTreeNode)original.Clone();
301      List<ISymbolicExpressionTreeNode> simplifiedSubtrees = new List<ISymbolicExpressionTreeNode>();
302      foreach (var subtree in subtrees) {
303        simplifiedSubtrees.Add(GetSimplifiedTree(subtree));
304        original.AddSubtree(subtree);
305      }
306      foreach (var simplifiedSubtree in simplifiedSubtrees) {
307        clone.AddSubtree(simplifiedSubtree);
308      }
309      if (simplifiedSubtrees.TrueForAll(t => IsConstant(t))) {
310        SimplifyConstantExpression(clone);
311      }
312      return clone;
313    }
314
315    private ISymbolicExpressionTreeNode SimplifyConstantExpression(ISymbolicExpressionTreeNode original) {
316      // not yet implemented
317      return original;
318    }
319
320    private ISymbolicExpressionTreeNode SimplifyAverage(ISymbolicExpressionTreeNode original) {
321      if (original.Subtrees.Count() == 1) {
322        return GetSimplifiedTree(original.GetSubtree(0));
323      } else {
324        // simplify expressions x0..xn
325        // make sum(x0..xn) / n
326        var sum = original.Subtrees
327          .Select(GetSimplifiedTree)
328          .Aggregate(MakeSum);
329        return MakeFraction(sum, MakeConstant(original.Subtrees.Count()));
330      }
331    }
332
333    private ISymbolicExpressionTreeNode SimplifyDivision(ISymbolicExpressionTreeNode original) {
334      if (original.Subtrees.Count() == 1) {
335        return Invert(GetSimplifiedTree(original.GetSubtree(0)));
336      } else {
337        // simplify expressions x0..xn
338        // make multiplication (x0 * 1/(x1 * x1 * .. * xn))
339        var first = original.GetSubtree(0);
340        var second = original.GetSubtree(1);
341        var remaining = original.Subtrees.Skip(2);
342        return
343          MakeProduct(GetSimplifiedTree(first),
344            Invert(remaining.Aggregate(GetSimplifiedTree(second), (a, b) => MakeProduct(a, GetSimplifiedTree(b)))));
345      }
346    }
347
348    private ISymbolicExpressionTreeNode SimplifyMultiplication(ISymbolicExpressionTreeNode original) {
349      if (original.Subtrees.Count() == 1) {
350        return GetSimplifiedTree(original.GetSubtree(0));
351      } else {
352        return original.Subtrees
353          .Select(GetSimplifiedTree)
354          .Aggregate(MakeProduct);
355      }
356    }
357
358    private ISymbolicExpressionTreeNode SimplifySubtraction(ISymbolicExpressionTreeNode original) {
359      if (original.Subtrees.Count() == 1) {
360        return Negate(GetSimplifiedTree(original.GetSubtree(0)));
361      } else {
362        // simplify expressions x0..xn
363        // make addition (x0,-x1..-xn)
364        var first = original.Subtrees.First();
365        var remaining = original.Subtrees.Skip(1);
366        return remaining.Aggregate(GetSimplifiedTree(first), (a, b) => MakeSum(a, Negate(GetSimplifiedTree(b))));
367      }
368    }
369
370    private ISymbolicExpressionTreeNode SimplifyAddition(ISymbolicExpressionTreeNode original) {
371      if (original.Subtrees.Count() == 1) {
372        return GetSimplifiedTree(original.GetSubtree(0));
373      } else {
374        // simplify expression x0..xn
375        // make addition (x0..xn)
376        return original.Subtrees
377          .Select(GetSimplifiedTree)
378          .Aggregate(MakeSum);
379      }
380    }
381
382    private ISymbolicExpressionTreeNode SimplifyNot(ISymbolicExpressionTreeNode original) {
383      return MakeNot(GetSimplifiedTree(original.GetSubtree(0)));
384    }
385
386    private ISymbolicExpressionTreeNode SimplifyOr(ISymbolicExpressionTreeNode original) {
387      return original.Subtrees
388        .Select(GetSimplifiedTree)
389        .Aggregate(MakeOr);
390    }
391
392    private ISymbolicExpressionTreeNode SimplifyAnd(ISymbolicExpressionTreeNode original) {
393      return original.Subtrees
394        .Select(GetSimplifiedTree)
395        .Aggregate(MakeAnd);
396    }
397
398    private ISymbolicExpressionTreeNode SimplifyLessThan(ISymbolicExpressionTreeNode original) {
399      return MakeLessThan(GetSimplifiedTree(original.GetSubtree(0)), GetSimplifiedTree(original.GetSubtree(1)));
400    }
401
402    private ISymbolicExpressionTreeNode SimplifyGreaterThan(ISymbolicExpressionTreeNode original) {
403      return MakeGreaterThan(GetSimplifiedTree(original.GetSubtree(0)), GetSimplifiedTree(original.GetSubtree(1)));
404    }
405
406    private ISymbolicExpressionTreeNode SimplifyIfThenElse(ISymbolicExpressionTreeNode original) {
407      return MakeIfThenElse(GetSimplifiedTree(original.GetSubtree(0)), GetSimplifiedTree(original.GetSubtree(1)),
408        GetSimplifiedTree(original.GetSubtree(2)));
409    }
410
411    private ISymbolicExpressionTreeNode SimplifyTangent(ISymbolicExpressionTreeNode original) {
412      return MakeTangent(GetSimplifiedTree(original.GetSubtree(0)));
413    }
414
415    private ISymbolicExpressionTreeNode SimplifyCosine(ISymbolicExpressionTreeNode original) {
416      return MakeCosine(GetSimplifiedTree(original.GetSubtree(0)));
417    }
418
419    private ISymbolicExpressionTreeNode SimplifySine(ISymbolicExpressionTreeNode original) {
420      return MakeSine(GetSimplifiedTree(original.GetSubtree(0)));
421    }
422
423    private ISymbolicExpressionTreeNode SimplifyExp(ISymbolicExpressionTreeNode original) {
424      return MakeExp(GetSimplifiedTree(original.GetSubtree(0)));
425    }
426
427    private ISymbolicExpressionTreeNode SimplifySquare(ISymbolicExpressionTreeNode original) {
428      return MakeSquare(GetSimplifiedTree(original.GetSubtree(0)));
429    }
430
431    private ISymbolicExpressionTreeNode SimplifySquareRoot(ISymbolicExpressionTreeNode original) {
432      return MakeSquareRoot(GetSimplifiedTree(original.GetSubtree(0)));
433    }
434
435    private ISymbolicExpressionTreeNode SimplifyLog(ISymbolicExpressionTreeNode original) {
436      return MakeLog(GetSimplifiedTree(original.GetSubtree(0)));
437    }
438
439    private ISymbolicExpressionTreeNode SimplifyRoot(ISymbolicExpressionTreeNode original) {
440      return MakeRoot(GetSimplifiedTree(original.GetSubtree(0)), GetSimplifiedTree(original.GetSubtree(1)));
441    }
442
443    private ISymbolicExpressionTreeNode SimplifyPower(ISymbolicExpressionTreeNode original) {
444      return MakePower(GetSimplifiedTree(original.GetSubtree(0)), GetSimplifiedTree(original.GetSubtree(1)));
445    }
446
447    private ISymbolicExpressionTreeNode SimplifyTimeLag(ISymbolicExpressionTreeNode original) {
448      var laggedTreeNode = original as ILaggedTreeNode;
449      var simplifiedSubtree = GetSimplifiedTree(original.GetSubtree(0));
450      if (!ContainsVariableCondition(simplifiedSubtree)) {
451        return AddLagToDynamicNodes(simplifiedSubtree, laggedTreeNode.Lag);
452      } else {
453        return MakeTimeLag(simplifiedSubtree, laggedTreeNode.Lag);
454      }
455    }
456
457    private ISymbolicExpressionTreeNode SimplifyIntegral(ISymbolicExpressionTreeNode original) {
458      var laggedTreeNode = original as ILaggedTreeNode;
459      var simplifiedSubtree = GetSimplifiedTree(original.GetSubtree(0));
460      if (IsConstant(simplifiedSubtree)) {
461        return GetSimplifiedTree(MakeProduct(simplifiedSubtree, MakeConstant(-laggedTreeNode.Lag)));
462      } else {
463        return MakeIntegral(simplifiedSubtree, laggedTreeNode.Lag);
464      }
465    }
466
467    #endregion
468
469    #region low level tree restructuring
470
471    private ISymbolicExpressionTreeNode MakeTimeLag(ISymbolicExpressionTreeNode subtree, int lag) {
472      if (lag == 0) return subtree;
473      if (IsConstant(subtree)) return subtree;
474      var lagNode = (LaggedTreeNode)timeLagSymbol.CreateTreeNode();
475      lagNode.Lag = lag;
476      lagNode.AddSubtree(subtree);
477      return lagNode;
478    }
479
480    private ISymbolicExpressionTreeNode MakeIntegral(ISymbolicExpressionTreeNode subtree, int lag) {
481      if (lag == 0) return subtree;
482      else if (lag == -1 || lag == 1) {
483        return MakeSum(subtree, AddLagToDynamicNodes((ISymbolicExpressionTreeNode)subtree.Clone(), lag));
484      } else {
485        var node = (LaggedTreeNode)integralSymbol.CreateTreeNode();
486        node.Lag = lag;
487        node.AddSubtree(subtree);
488        return node;
489      }
490    }
491
492    private ISymbolicExpressionTreeNode MakeNot(ISymbolicExpressionTreeNode t) {
493      if (IsConstant(t)) {
494        var constNode = t as ConstantTreeNode;
495        if (constNode.Value > 0) return MakeConstant(-1.0);
496        else return MakeConstant(1.0);
497      } else if (IsNot(t)) {
498        return t.GetSubtree(0);
499      } else if (!IsBoolean(t)) {
500        var gtNode = gtSymbol.CreateTreeNode();
501        gtNode.AddSubtree(t);
502        gtNode.AddSubtree(MakeConstant(0.0));
503        var notNode = notSymbol.CreateTreeNode();
504        notNode.AddSubtree(gtNode);
505        return notNode;
506      } else {
507        var notNode = notSymbol.CreateTreeNode();
508        notNode.AddSubtree(t);
509        return notNode;
510      }
511    }
512
513    private ISymbolicExpressionTreeNode MakeOr(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
514      if (IsConstant(a) && IsConstant(b)) {
515        var constA = a as ConstantTreeNode;
516        var constB = b as ConstantTreeNode;
517        if (constA.Value > 0.0 || constB.Value > 0.0) {
518          return MakeConstant(1.0);
519        } else {
520          return MakeConstant(-1.0);
521        }
522      } else if (IsConstant(a)) {
523        return MakeOr(b, a);
524      } else if (IsConstant(b)) {
525        var constT = b as ConstantTreeNode;
526        if (constT.Value > 0.0) {
527          // boolean expression is necessarily true
528          return MakeConstant(1.0);
529        } else {
530          // the constant value has no effect on the result of the boolean condition so we can drop the constant term
531          var orNode = orSymbol.CreateTreeNode();
532          orNode.AddSubtree(a);
533          return orNode;
534        }
535      } else {
536        var orNode = orSymbol.CreateTreeNode();
537        orNode.AddSubtree(a);
538        orNode.AddSubtree(b);
539        return orNode;
540      }
541    }
542
543    private ISymbolicExpressionTreeNode MakeAnd(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
544      if (IsConstant(a) && IsConstant(b)) {
545        var constA = a as ConstantTreeNode;
546        var constB = b as ConstantTreeNode;
547        if (constA.Value > 0.0 && constB.Value > 0.0) {
548          return MakeConstant(1.0);
549        } else {
550          return MakeConstant(-1.0);
551        }
552      } else if (IsConstant(a)) {
553        return MakeAnd(b, a);
554      } else if (IsConstant(b)) {
555        var constB = b as ConstantTreeNode;
556        if (constB.Value > 0.0) {
557          // the constant value has no effect on the result of the boolean condition so we can drop the constant term
558          var andNode = andSymbol.CreateTreeNode();
559          andNode.AddSubtree(a);
560          return andNode;
561        } else {
562          // boolean expression is necessarily false
563          return MakeConstant(-1.0);
564        }
565      } else {
566        var andNode = andSymbol.CreateTreeNode();
567        andNode.AddSubtree(a);
568        andNode.AddSubtree(b);
569        return andNode;
570      }
571    }
572
573    private ISymbolicExpressionTreeNode MakeLessThan(ISymbolicExpressionTreeNode leftSide,
574      ISymbolicExpressionTreeNode rightSide) {
575      if (IsConstant(leftSide) && IsConstant(rightSide)) {
576        var lsConst = leftSide as ConstantTreeNode;
577        var rsConst = rightSide as ConstantTreeNode;
578        if (lsConst.Value < rsConst.Value) return MakeConstant(1.0);
579        else return MakeConstant(-1.0);
580      } else {
581        var ltNode = ltSymbol.CreateTreeNode();
582        ltNode.AddSubtree(leftSide);
583        ltNode.AddSubtree(rightSide);
584        return ltNode;
585      }
586    }
587
588    private ISymbolicExpressionTreeNode MakeGreaterThan(ISymbolicExpressionTreeNode leftSide,
589      ISymbolicExpressionTreeNode rightSide) {
590      if (IsConstant(leftSide) && IsConstant(rightSide)) {
591        var lsConst = leftSide as ConstantTreeNode;
592        var rsConst = rightSide as ConstantTreeNode;
593        if (lsConst.Value > rsConst.Value) return MakeConstant(1.0);
594        else return MakeConstant(-1.0);
595      } else {
596        var gtNode = gtSymbol.CreateTreeNode();
597        gtNode.AddSubtree(leftSide);
598        gtNode.AddSubtree(rightSide);
599        return gtNode;
600      }
601    }
602
603    private ISymbolicExpressionTreeNode MakeIfThenElse(ISymbolicExpressionTreeNode condition,
604      ISymbolicExpressionTreeNode trueBranch, ISymbolicExpressionTreeNode falseBranch) {
605      if (IsConstant(condition)) {
606        var constT = condition as ConstantTreeNode;
607        if (constT.Value > 0.0) return trueBranch;
608        else return falseBranch;
609      } else {
610        var ifNode = ifThenElseSymbol.CreateTreeNode();
611        if (IsBoolean(condition)) {
612          ifNode.AddSubtree(condition);
613        } else {
614          var gtNode = gtSymbol.CreateTreeNode();
615          gtNode.AddSubtree(condition);
616          gtNode.AddSubtree(MakeConstant(0.0));
617          ifNode.AddSubtree(gtNode);
618        }
619        ifNode.AddSubtree(trueBranch);
620        ifNode.AddSubtree(falseBranch);
621        return ifNode;
622      }
623    }
624
625    private ISymbolicExpressionTreeNode MakeSine(ISymbolicExpressionTreeNode node) {
626      if (IsConstant(node)) {
627        var constT = node as ConstantTreeNode;
628        return MakeConstant(Math.Sin(constT.Value));
629      } else if (IsFactor(node)) {
630        var factor = node as FactorVariableTreeNode;
631        return MakeFactor(factor.Symbol, factor.VariableName, factor.Weights.Select(Math.Sin));
632      } else if (IsBinFactor(node)) {
633        var binFactor = node as BinaryFactorVariableTreeNode;
634        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Sin(binFactor.Weight));
635      } else {
636        var sineNode = sineSymbol.CreateTreeNode();
637        sineNode.AddSubtree(node);
638        return sineNode;
639      }
640    }
641
642    private ISymbolicExpressionTreeNode MakeTangent(ISymbolicExpressionTreeNode node) {
643      if (IsConstant(node)) {
644        var constT = node as ConstantTreeNode;
645        return MakeConstant(Math.Tan(constT.Value));
646      } else if (IsFactor(node)) {
647        var factor = node as FactorVariableTreeNode;
648        return MakeFactor(factor.Symbol, factor.VariableName, factor.Weights.Select(Math.Tan));
649      } else if (IsBinFactor(node)) {
650        var binFactor = node as BinaryFactorVariableTreeNode;
651        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Tan(binFactor.Weight));
652      } else {
653        var tanNode = tanSymbol.CreateTreeNode();
654        tanNode.AddSubtree(node);
655        return tanNode;
656      }
657    }
658
659    private ISymbolicExpressionTreeNode MakeCosine(ISymbolicExpressionTreeNode node) {
660      if (IsConstant(node)) {
661        var constT = node as ConstantTreeNode;
662        return MakeConstant(Math.Cos(constT.Value));
663      } else if (IsFactor(node)) {
664        var factor = node as FactorVariableTreeNode;
665        return MakeFactor(factor.Symbol, factor.VariableName, factor.Weights.Select(Math.Cos));
666      } else if (IsBinFactor(node)) {
667        var binFactor = node as BinaryFactorVariableTreeNode;
668        // cos(0) = 1 see similar case for Exp(binfactor)
669        return MakeSum(MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Cos(binFactor.Weight) - 1),
670          MakeConstant(1.0));
671      } else {
672        var cosNode = cosineSymbol.CreateTreeNode();
673        cosNode.AddSubtree(node);
674        return cosNode;
675      }
676    }
677
678    private ISymbolicExpressionTreeNode MakeExp(ISymbolicExpressionTreeNode node) {
679      if (IsConstant(node)) {
680        var constT = node as ConstantTreeNode;
681        return MakeConstant(Math.Exp(constT.Value));
682      } else if (IsFactor(node)) {
683        var factNode = node as FactorVariableTreeNode;
684        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Exp(w)));
685      } else if (IsBinFactor(node)) {
686        // exp( binfactor w val=a) = if(val=a) exp(w) else exp(0) = binfactor( (exp(w) - 1) val a) + 1
687        var binFactor = node as BinaryFactorVariableTreeNode;
688        return
689          MakeSum(MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Exp(binFactor.Weight) - 1), MakeConstant(1.0));
690      } else if (IsLog(node)) {
691        return node.GetSubtree(0);
692      } else if (IsAddition(node)) {
693        return node.Subtrees.Select(s => MakeExp(s)).Aggregate((s, t) => MakeProduct(s, t));
694      } else if (IsSubtraction(node)) {
695        return node.Subtrees.Select(s => MakeExp(s)).Aggregate((s, t) => MakeProduct(s, Negate(t)));
696      } else {
697        var expNode = expSymbol.CreateTreeNode();
698        expNode.AddSubtree(node);
699        return expNode;
700      }
701    }
702    private ISymbolicExpressionTreeNode MakeLog(ISymbolicExpressionTreeNode node) {
703      if (IsConstant(node)) {
704        var constT = node as ConstantTreeNode;
705        return MakeConstant(Math.Log(constT.Value));
706      } else if (IsFactor(node)) {
707        var factNode = node as FactorVariableTreeNode;
708        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Log(w)));
709      } else if (IsExp(node)) {
710        return node.GetSubtree(0);
711      } else if (IsSquareRoot(node)) {
712        return MakeFraction(MakeLog(node.GetSubtree(0)), MakeConstant(2.0));
713      } else {
714        var logNode = logSymbol.CreateTreeNode();
715        logNode.AddSubtree(node);
716        return logNode;
717      }
718    }
719
720    private ISymbolicExpressionTreeNode MakeSquare(ISymbolicExpressionTreeNode node) {
721      if (IsConstant(node)) {
722        var constT = node as ConstantTreeNode;
723        return MakeConstant(constT.Value * constT.Value);
724      } else if (IsFactor(node)) {
725        var factNode = node as FactorVariableTreeNode;
726        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => w * w));
727      } else if (IsBinFactor(node)) {
728        var binFactor = node as BinaryFactorVariableTreeNode;
729        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, binFactor.Weight * binFactor.Weight);
730      } else if (IsSquareRoot(node)) {
731        return node.GetSubtree(0);
732      } else {
733        var sqrNode = sqrSymbol.CreateTreeNode();
734        sqrNode.AddSubtree(node);
735        return sqrNode;
736      }
737    }
738
739    private ISymbolicExpressionTreeNode MakeSquareRoot(ISymbolicExpressionTreeNode node) {
740      if (IsConstant(node)) {
741        var constT = node as ConstantTreeNode;
742        return MakeConstant(Math.Sqrt(constT.Value));
743      } else if (IsFactor(node)) {
744        var factNode = node as FactorVariableTreeNode;
745        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Sqrt(w)));
746      } else if (IsBinFactor(node)) {
747        var binFactor = node as BinaryFactorVariableTreeNode;
748        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Sqrt(binFactor.Weight));
749      } else if (IsSquare(node)) {
750        return node.GetSubtree(0);
751      } else {
752        var sqrtNode = sqrtSymbol.CreateTreeNode();
753        sqrtNode.AddSubtree(node);
754        return sqrtNode;
755      }
756    }
757
758    private ISymbolicExpressionTreeNode MakeRoot(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
759      if (IsConstant(a) && IsConstant(b)) {
760        var constA = a as ConstantTreeNode;
761        var constB = b as ConstantTreeNode;
762        return MakeConstant(Math.Pow(constA.Value, 1.0 / Math.Round(constB.Value)));
763      } else if (IsFactor(a) && IsConstant(b)) {
764        var factNode = a as FactorVariableTreeNode;
765        var constNode = b as ConstantTreeNode;
766        return MakeFactor(factNode.Symbol, factNode.VariableName,
767          factNode.Weights.Select(w => Math.Pow(w, 1.0 / Math.Round(constNode.Value))));
768      } else if (IsBinFactor(a) && IsConstant(b)) {
769        var binFactor = a as BinaryFactorVariableTreeNode;
770        var constNode = b as ConstantTreeNode;
771        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Pow(binFactor.Weight, 1.0 / Math.Round(constNode.Value)));
772      } else if (IsConstant(a) && IsFactor(b)) {
773        var constNode = a as ConstantTreeNode;
774        var factNode = b as FactorVariableTreeNode;
775        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Pow(constNode.Value, 1.0 / Math.Round(w))));
776      } else if (IsConstant(a) && IsBinFactor(b)) {
777        var constNode = a as ConstantTreeNode;
778        var factNode = b as BinaryFactorVariableTreeNode;
779        return MakeBinFactor(factNode.Symbol, factNode.VariableName, factNode.VariableValue, Math.Pow(constNode.Value, 1.0 / Math.Round(factNode.Weight)));
780      } else if (IsFactor(a) && IsFactor(b) && AreSameTypeAndVariable(a, b)) {
781        var node0 = a as FactorVariableTreeNode;
782        var node1 = b as FactorVariableTreeNode;
783        return MakeFactor(node0.Symbol, node0.VariableName, node0.Weights.Zip(node1.Weights, (u, v) => Math.Pow(u, 1.0 / Math.Round(v))));
784      } else if (IsConstant(b)) {
785        var constB = b as ConstantTreeNode;
786        var constBValue = Math.Round(constB.Value);
787        if (constBValue.IsAlmost(1.0)) {
788          return a;
789        } else if (constBValue.IsAlmost(0.0)) {
790          return MakeConstant(1.0);
791        } else if (constBValue.IsAlmost(-1.0)) {
792          return MakeFraction(MakeConstant(1.0), a);
793        } else if (constBValue < 0) {
794          var rootNode = rootSymbol.CreateTreeNode();
795          rootNode.AddSubtree(a);
796          rootNode.AddSubtree(MakeConstant(-1.0 * constBValue));
797          return MakeFraction(MakeConstant(1.0), rootNode);
798        } else {
799          var rootNode = rootSymbol.CreateTreeNode();
800          rootNode.AddSubtree(a);
801          rootNode.AddSubtree(MakeConstant(constBValue));
802          return rootNode;
803        }
804      } else {
805        var rootNode = rootSymbol.CreateTreeNode();
806        rootNode.AddSubtree(a);
807        rootNode.AddSubtree(b);
808        return rootNode;
809      }
810    }
811
812
813    private ISymbolicExpressionTreeNode MakePower(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
814      if (IsConstant(a) && IsConstant(b)) {
815        var constA = a as ConstantTreeNode;
816        var constB = b as ConstantTreeNode;
817        return MakeConstant(Math.Pow(constA.Value, Math.Round(constB.Value)));
818      } else if (IsFactor(a) && IsConstant(b)) {
819        var factNode = a as FactorVariableTreeNode;
820        var constNode = b as ConstantTreeNode;
821        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Pow(w, Math.Round(constNode.Value))));
822      } else if (IsBinFactor(a) && IsConstant(b)) {
823        var binFactor = a as BinaryFactorVariableTreeNode;
824        var constNode = b as ConstantTreeNode;
825        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Pow(binFactor.Weight, Math.Round(constNode.Value)));
826      } else if (IsConstant(a) && IsFactor(b)) {
827        var constNode = a as ConstantTreeNode;
828        var factNode = b as FactorVariableTreeNode;
829        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Pow(constNode.Value, Math.Round(w))));
830      } else if (IsConstant(a) && IsBinFactor(b)) {
831        var constNode = a as ConstantTreeNode;
832        var factNode = b as BinaryFactorVariableTreeNode;
833        return MakeBinFactor(factNode.Symbol, factNode.VariableName, factNode.VariableValue, Math.Pow(constNode.Value, Math.Round(factNode.Weight)));
834      } else if (IsFactor(a) && IsFactor(b) && AreSameTypeAndVariable(a, b)) {
835        var node0 = a as FactorVariableTreeNode;
836        var node1 = b as FactorVariableTreeNode;
837        return MakeFactor(node0.Symbol, node0.VariableName, node0.Weights.Zip(node1.Weights, (u, v) => Math.Pow(u, Math.Round(v))));
838      } else if (IsConstant(b)) {
839        var constB = b as ConstantTreeNode;
840        double exponent = Math.Round(constB.Value);
841        if (exponent.IsAlmost(0.0)) {
842          return MakeConstant(1.0);
843        } else if (exponent.IsAlmost(1.0)) {
844          return a;
845        } else if (exponent.IsAlmost(-1.0)) {
846          return MakeFraction(MakeConstant(1.0), a);
847        } else if (exponent < 0) {
848          var powNode = powSymbol.CreateTreeNode();
849          powNode.AddSubtree(a);
850          powNode.AddSubtree(MakeConstant(-1.0 * exponent));
851          return MakeFraction(MakeConstant(1.0), powNode);
852        } else {
853          var powNode = powSymbol.CreateTreeNode();
854          powNode.AddSubtree(a);
855          powNode.AddSubtree(MakeConstant(exponent));
856          return powNode;
857        }
858      } else {
859        var powNode = powSymbol.CreateTreeNode();
860        powNode.AddSubtree(a);
861        powNode.AddSubtree(b);
862        return powNode;
863      }
864    }
865
866
867    // MakeFraction, MakeProduct and MakeSum take two already simplified trees and create a new simplified tree
868    private ISymbolicExpressionTreeNode MakeFraction(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
869      if (IsConstant(a) && IsConstant(b)) {
870        // fold constants
871        return MakeConstant(((ConstantTreeNode)a).Value / ((ConstantTreeNode)b).Value);
872      } else if ((IsConstant(a) && !((ConstantTreeNode)a).Value.IsAlmost(1.0))) {
873        return MakeFraction(MakeConstant(1.0), MakeProduct(b, Invert(a)));
874      } else if (IsVariableBase(a) && IsConstant(b)) {
875        // merge constant values into variable weights
876        var constB = ((ConstantTreeNode)b).Value;
877        ((VariableTreeNodeBase)a).Weight /= constB;
878        return a;
879      } else if (IsFactor(a) && IsConstant(b)) {
880        var factNode = a as FactorVariableTreeNode;
881        var constNode = b as ConstantTreeNode;
882        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => w / constNode.Value));
883      } else if (IsBinFactor(a) && IsConstant(b)) {
884        var factNode = a as BinaryFactorVariableTreeNode;
885        var constNode = b as ConstantTreeNode;
886        return MakeBinFactor(factNode.Symbol, factNode.VariableName, factNode.VariableValue, factNode.Weight / constNode.Value);
887      } else if (IsFactor(a) && IsFactor(b) && AreSameTypeAndVariable(a, b)) {
888        var node0 = a as FactorVariableTreeNode;
889        var node1 = b as FactorVariableTreeNode;
890        return MakeFactor(node0.Symbol, node0.VariableName, node0.Weights.Zip(node1.Weights, (u, v) => u / v));
891      } else if (IsFactor(a) && IsBinFactor(b) && ((IVariableTreeNode)a).VariableName == ((IVariableTreeNode)b).VariableName) {
892        var node0 = a as FactorVariableTreeNode;
893        var node1 = b as BinaryFactorVariableTreeNode;
894        var varValues = node0.Symbol.GetVariableValues(node0.VariableName).ToArray();
895        var wi = Array.IndexOf(varValues, node1.VariableValue);
896        if (wi < 0) throw new ArgumentException();
897        var newWeighs = new double[varValues.Length];
898        node0.Weights.CopyTo(newWeighs, 0);
899        for (int i = 0; i < newWeighs.Length; i++)
900          if (wi == i) newWeighs[i] /= node1.Weight;
901          else newWeighs[i] /= 0.0;
902        return MakeFactor(node0.Symbol, node0.VariableName, newWeighs);
903      } else if (IsFactor(a)) {
904        return MakeFraction(MakeConstant(1.0), MakeProduct(b, Invert(a)));
905      } else if (IsVariableBase(a) && IsVariableBase(b) && AreSameTypeAndVariable(a, b) && !IsBinFactor(b)) {
906        // cancel variables (not allowed for bin factors because of division by zero)
907        var aVar = a as VariableTreeNode;
908        var bVar = b as VariableTreeNode;
909        return MakeConstant(aVar.Weight / bVar.Weight);
910      } else if (IsAddition(a) && IsConstant(b)) {
911        return a.Subtrees
912          .Select(x => GetSimplifiedTree(x))
913          .Select(x => MakeFraction(x, GetSimplifiedTree(b)))
914          .Aggregate((c, d) => MakeSum(c, d));
915      } else if (IsMultiplication(a) && IsConstant(b)) {
916        return MakeProduct(a, Invert(b));
917      } else if (IsDivision(a) && IsConstant(b)) {
918        // (a1 / a2) / c => (a1 / (a2 * c))
919        return MakeFraction(a.GetSubtree(0), MakeProduct(a.GetSubtree(1), b));
920      } else if (IsDivision(a) && IsDivision(b)) {
921        // (a1 / a2) / (b1 / b2) =>
922        return MakeFraction(MakeProduct(a.GetSubtree(0), b.GetSubtree(1)), MakeProduct(a.GetSubtree(1), b.GetSubtree(0)));
923      } else if (IsDivision(a)) {
924        // (a1 / a2) / b => (a1 / (a2 * b))
925        return MakeFraction(a.GetSubtree(0), MakeProduct(a.GetSubtree(1), b));
926      } else if (IsDivision(b)) {
927        // a / (b1 / b2) => (a * b2) / b1
928        return MakeFraction(MakeProduct(a, b.GetSubtree(1)), b.GetSubtree(0));
929      } else {
930        var div = divSymbol.CreateTreeNode();
931        div.AddSubtree(a);
932        div.AddSubtree(b);
933        return div;
934      }
935    }
936
937    private ISymbolicExpressionTreeNode MakeSum(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
938      if (IsConstant(a) && IsConstant(b)) {
939        // fold constants
940        ((ConstantTreeNode)a).Value += ((ConstantTreeNode)b).Value;
941        return a;
942      } else if (IsConstant(a)) {
943        // c + x => x + c
944        // b is not constant => make sure constant is on the right
945        return MakeSum(b, a);
946      } else if (IsConstant(b) && ((ConstantTreeNode)b).Value.IsAlmost(0.0)) {
947        // x + 0 => x
948        return a;
949      } else if (IsFactor(a) && IsConstant(b)) {
950        var factNode = a as FactorVariableTreeNode;
951        var constNode = b as ConstantTreeNode;
952        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select((w) => w + constNode.Value));
953      } else if (IsFactor(a) && IsFactor(b) && AreSameTypeAndVariable(a, b)) {
954        var node0 = a as FactorVariableTreeNode;
955        var node1 = b as FactorVariableTreeNode;
956        return MakeFactor(node0.Symbol, node0.VariableName, node0.Weights.Zip(node1.Weights, (u, v) => u + v));
957      } else if (IsBinFactor(a) && IsFactor(b)) {
958        return MakeSum(b, a);
959      } else if (IsFactor(a) && IsBinFactor(b) &&
960        ((IVariableTreeNode)a).VariableName == ((IVariableTreeNode)b).VariableName) {
961        var node0 = a as FactorVariableTreeNode;
962        var node1 = b as BinaryFactorVariableTreeNode;
963        var varValues = node0.Symbol.GetVariableValues(node0.VariableName).ToArray();
964        var wi = Array.IndexOf(varValues, node1.VariableValue);
965        if (wi < 0) throw new ArgumentException();
966        var newWeighs = new double[varValues.Length];
967        node0.Weights.CopyTo(newWeighs, 0);
968        newWeighs[wi] += node1.Weight;
969        return MakeFactor(node0.Symbol, node0.VariableName, newWeighs);
970      } else if (IsAddition(a) && IsAddition(b)) {
971        // merge additions
972        var add = addSymbol.CreateTreeNode();
973        // add all sub trees except for the last
974        for (int i = 0; i < a.Subtrees.Count() - 1; i++) add.AddSubtree(a.GetSubtree(i));
975        for (int i = 0; i < b.Subtrees.Count() - 1; i++) add.AddSubtree(b.GetSubtree(i));
976        if (IsConstant(a.Subtrees.Last()) && IsConstant(b.Subtrees.Last())) {
977          add.AddSubtree(MakeSum(a.Subtrees.Last(), b.Subtrees.Last()));
978        } else if (IsConstant(a.Subtrees.Last())) {
979          add.AddSubtree(b.Subtrees.Last());
980          add.AddSubtree(a.Subtrees.Last());
981        } else {
982          add.AddSubtree(a.Subtrees.Last());
983          add.AddSubtree(b.Subtrees.Last());
984        }
985        MergeVariablesInSum(add);
986        if (add.Subtrees.Count() == 1) {
987          return add.GetSubtree(0);
988        } else {
989          return add;
990        }
991      } else if (IsAddition(b)) {
992        return MakeSum(b, a);
993      } else if (IsAddition(a) && IsConstant(b)) {
994        // a is an addition and b is a constant => append b to a and make sure the constants are merged
995        var add = addSymbol.CreateTreeNode();
996        // add all sub trees except for the last
997        for (int i = 0; i < a.Subtrees.Count() - 1; i++) add.AddSubtree(a.GetSubtree(i));
998        if (IsConstant(a.Subtrees.Last()))
999          add.AddSubtree(MakeSum(a.Subtrees.Last(), b));
1000        else {
1001          add.AddSubtree(a.Subtrees.Last());
1002          add.AddSubtree(b);
1003        }
1004        return add;
1005      } else if (IsAddition(a)) {
1006        // a is already an addition => append b
1007        var add = addSymbol.CreateTreeNode();
1008        add.AddSubtree(b);
1009        foreach (var subtree in a.Subtrees) {
1010          add.AddSubtree(subtree);
1011        }
1012        MergeVariablesInSum(add);
1013        if (add.Subtrees.Count() == 1) {
1014          return add.GetSubtree(0);
1015        } else {
1016          return add;
1017        }
1018      } else {
1019        var add = addSymbol.CreateTreeNode();
1020        add.AddSubtree(a);
1021        add.AddSubtree(b);
1022        MergeVariablesInSum(add);
1023        if (add.Subtrees.Count() == 1) {
1024          return add.GetSubtree(0);
1025        } else {
1026          return add;
1027        }
1028      }
1029    }
1030
1031    // makes sure variable symbols in sums are combined
1032    private void MergeVariablesInSum(ISymbolicExpressionTreeNode sum) {
1033      var subtrees = new List<ISymbolicExpressionTreeNode>(sum.Subtrees);
1034      while (sum.Subtrees.Any()) sum.RemoveSubtree(0);
1035      var groupedVarNodes = from node in subtrees.OfType<IVariableTreeNode>()
1036                            where node.SubtreeCount == 0
1037                            group node by GroupId(node) into g
1038                            select g;
1039      var constant = (from node in subtrees.OfType<ConstantTreeNode>()
1040                      select node.Value).DefaultIfEmpty(0.0).Sum();
1041      var unchangedSubtrees = subtrees.Where(t => t.SubtreeCount > 0 || !(t is IVariableTreeNode) && !(t is ConstantTreeNode));
1042
1043      foreach (var variableNodeGroup in groupedVarNodes) {
1044        var firstNode = variableNodeGroup.First();
1045        if (firstNode is VariableTreeNodeBase) {
1046          var representative = firstNode as VariableTreeNodeBase;
1047          var weightSum = variableNodeGroup.Cast<VariableTreeNodeBase>().Select(t => t.Weight).Sum();
1048          representative.Weight = weightSum;
1049          sum.AddSubtree(representative);
1050        } else if (firstNode is FactorVariableTreeNode) {
1051          var representative = firstNode as FactorVariableTreeNode;
1052          foreach (var node in variableNodeGroup.Skip(1).Cast<FactorVariableTreeNode>()) {
1053            for (int j = 0; j < representative.Weights.Length; j++) {
1054              representative.Weights[j] += node.Weights[j];
1055            }
1056          }
1057          for (int j = 0; j < representative.Weights.Length; j++) {
1058            representative.Weights[j] += constant;
1059          }
1060          sum.AddSubtree(representative);
1061        }
1062      }
1063      foreach (var unchangedSubtree in unchangedSubtrees)
1064        sum.AddSubtree(unchangedSubtree);
1065      if (!constant.IsAlmost(0.0)) {
1066        sum.AddSubtree(MakeConstant(constant));
1067      }
1068    }
1069
1070    // nodes referencing variables can be grouped if they have
1071    private string GroupId(IVariableTreeNode node) {
1072      var binaryFactorNode = node as BinaryFactorVariableTreeNode;
1073      var factorNode = node as FactorVariableTreeNode;
1074      var variableNode = node as VariableTreeNode;
1075      var laggedVarNode = node as LaggedVariableTreeNode;
1076      if (variableNode != null) {
1077        return "var " + variableNode.VariableName;
1078      } else if (binaryFactorNode != null) {
1079        return "binfactor " + binaryFactorNode.VariableName + " " + binaryFactorNode.VariableValue;
1080      } else if (factorNode != null) {
1081        return "factor " + factorNode.VariableName;
1082      } else if (laggedVarNode != null) {
1083        return "lagged " + laggedVarNode.VariableName + " " + laggedVarNode.Lag;
1084      } else {
1085        throw new NotSupportedException();
1086      }
1087    }
1088
1089
1090    private ISymbolicExpressionTreeNode MakeProduct(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
1091      if (IsConstant(a) && IsConstant(b)) {
1092        // fold constants
1093        return MakeConstant(((ConstantTreeNode)a).Value * ((ConstantTreeNode)b).Value);
1094      } else if (IsConstant(a)) {
1095        // a * $ => $ * a
1096        return MakeProduct(b, a);
1097      } else if (IsFactor(a) && IsFactor(b) && AreSameTypeAndVariable(a, b)) {
1098        var node0 = a as FactorVariableTreeNode;
1099        var node1 = b as FactorVariableTreeNode;
1100        return MakeFactor(node0.Symbol, node0.VariableName, node0.Weights.Zip(node1.Weights, (u, v) => u * v));
1101      } else if (IsBinFactor(a) && IsBinFactor(b) && AreSameTypeAndVariable(a, b)) {
1102        var node0 = a as BinaryFactorVariableTreeNode;
1103        var node1 = b as BinaryFactorVariableTreeNode;
1104        return MakeBinFactor(node0.Symbol, node0.VariableName, node0.VariableValue, node0.Weight * node1.Weight);
1105      } else if (IsFactor(a) && IsConstant(b)) {
1106        var node0 = a as FactorVariableTreeNode;
1107        var node1 = b as ConstantTreeNode;
1108        return MakeFactor(node0.Symbol, node0.VariableName, node0.Weights.Select(w => w * node1.Value));
1109      } else if (IsBinFactor(a) && IsConstant(b)) {
1110        var node0 = a as BinaryFactorVariableTreeNode;
1111        var node1 = b as ConstantTreeNode;
1112        return MakeBinFactor(node0.Symbol, node0.VariableName, node0.VariableValue, node0.Weight * node1.Value);
1113      } else if (IsBinFactor(a) && IsFactor(b)) {
1114        return MakeProduct(b, a);
1115      } else if (IsFactor(a) && IsBinFactor(b) &&
1116        ((IVariableTreeNode)a).VariableName == ((IVariableTreeNode)b).VariableName) {
1117        var node0 = a as FactorVariableTreeNode;
1118        var node1 = b as BinaryFactorVariableTreeNode;
1119        var varValues = node0.Symbol.GetVariableValues(node0.VariableName).ToArray();
1120        var wi = Array.IndexOf(varValues, node1.VariableValue);
1121        if (wi < 0) throw new ArgumentException();
1122        return MakeBinFactor(node1.Symbol, node1.VariableName, node1.VariableValue, node1.Weight * node0.Weights[wi]);
1123      } else if (IsConstant(b) && ((ConstantTreeNode)b).Value.IsAlmost(1.0)) {
1124        // $ * 1.0 => $
1125        return a;
1126      } else if (IsConstant(b) && IsVariableBase(a)) {
1127        // multiply constants into variables weights
1128        ((VariableTreeNodeBase)a).Weight *= ((ConstantTreeNode)b).Value;
1129        return a;
1130      } else if (IsConstant(b) && IsAddition(a) ||
1131          IsFactor(b) && IsAddition(a) ||
1132          IsBinFactor(b) && IsAddition(a)) {
1133        // multiply constants into additions
1134        return a.Subtrees.Select(x => MakeProduct(GetSimplifiedTree(x), GetSimplifiedTree(b))).Aggregate((c, d) => MakeSum(c, d));
1135      } else if (IsDivision(a) && IsDivision(b)) {
1136        // (a1 / a2) * (b1 / b2) => (a1 * b1) / (a2 * b2)
1137        return MakeFraction(MakeProduct(a.GetSubtree(0), b.GetSubtree(0)), MakeProduct(a.GetSubtree(1), b.GetSubtree(1)));
1138      } else if (IsDivision(a)) {
1139        // (a1 / a2) * b => (a1 * b) / a2
1140        return MakeFraction(MakeProduct(a.GetSubtree(0), b), a.GetSubtree(1));
1141      } else if (IsDivision(b)) {
1142        // a * (b1 / b2) => (b1 * a) / b2
1143        return MakeFraction(MakeProduct(b.GetSubtree(0), a), b.GetSubtree(1));
1144      } else if (IsMultiplication(a) && IsMultiplication(b)) {
1145        // merge multiplications (make sure constants are merged)
1146        var mul = mulSymbol.CreateTreeNode();
1147        for (int i = 0; i < a.Subtrees.Count(); i++) mul.AddSubtree(a.GetSubtree(i));
1148        for (int i = 0; i < b.Subtrees.Count(); i++) mul.AddSubtree(b.GetSubtree(i));
1149        MergeVariablesAndConstantsInProduct(mul);
1150        return mul;
1151      } else if (IsMultiplication(b)) {
1152        return MakeProduct(b, a);
1153      } else if (IsMultiplication(a)) {
1154        // a is already an multiplication => append b
1155        a.AddSubtree(GetSimplifiedTree(b));
1156        MergeVariablesAndConstantsInProduct(a);
1157        return a;
1158      } else {
1159        var mul = mulSymbol.CreateTreeNode();
1160        mul.AddSubtree(a);
1161        mul.AddSubtree(b);
1162        MergeVariablesAndConstantsInProduct(mul);
1163        return mul;
1164      }
1165    }
1166
1167    #endregion
1168
1169    #region helper functions
1170
1171    private bool ContainsVariableCondition(ISymbolicExpressionTreeNode node) {
1172      if (node.Symbol is VariableCondition) return true;
1173      foreach (var subtree in node.Subtrees)
1174        if (ContainsVariableCondition(subtree)) return true;
1175      return false;
1176    }
1177
1178    private ISymbolicExpressionTreeNode AddLagToDynamicNodes(ISymbolicExpressionTreeNode node, int lag) {
1179      var laggedTreeNode = node as ILaggedTreeNode;
1180      var variableNode = node as VariableTreeNode;
1181      var variableConditionNode = node as VariableConditionTreeNode;
1182      if (laggedTreeNode != null)
1183        laggedTreeNode.Lag += lag;
1184      else if (variableNode != null) {
1185        var laggedVariableNode = (LaggedVariableTreeNode)laggedVariableSymbol.CreateTreeNode();
1186        laggedVariableNode.Lag = lag;
1187        laggedVariableNode.VariableName = variableNode.VariableName;
1188        return laggedVariableNode;
1189      } else if (variableConditionNode != null) {
1190        throw new NotSupportedException("Removal of time lags around variable condition symbols is not allowed.");
1191      }
1192      var subtrees = new List<ISymbolicExpressionTreeNode>(node.Subtrees);
1193      while (node.SubtreeCount > 0) node.RemoveSubtree(0);
1194      foreach (var subtree in subtrees) {
1195        node.AddSubtree(AddLagToDynamicNodes(subtree, lag));
1196      }
1197      return node;
1198    }
1199
1200    private bool AreSameTypeAndVariable(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
1201      return GroupId((IVariableTreeNode)a) == GroupId((IVariableTreeNode)b);
1202    }
1203
1204    // helper to combine the constant factors in products and to combine variables (powers of 2, 3...)
1205    private void MergeVariablesAndConstantsInProduct(ISymbolicExpressionTreeNode prod) {
1206      var subtrees = new List<ISymbolicExpressionTreeNode>(prod.Subtrees);
1207      while (prod.Subtrees.Any()) prod.RemoveSubtree(0);
1208      var groupedVarNodes = from node in subtrees.OfType<IVariableTreeNode>()
1209                            where node.SubtreeCount == 0
1210                            group node by GroupId(node) into g
1211                            orderby g.Count()
1212                            select g;
1213      var constantProduct = (from node in subtrees.OfType<VariableTreeNodeBase>()
1214                             select node.Weight)
1215        .Concat(from node in subtrees.OfType<ConstantTreeNode>()
1216                select node.Value)
1217        .DefaultIfEmpty(1.0)
1218        .Aggregate((c1, c2) => c1 * c2);
1219
1220      var unchangedSubtrees = from tree in subtrees
1221                              where tree.SubtreeCount > 0 || !(tree is IVariableTreeNode) && !(tree is ConstantTreeNode)
1222                              select tree;
1223
1224      foreach (var variableNodeGroup in groupedVarNodes) {
1225        var firstNode = variableNodeGroup.First();
1226        if (firstNode is VariableTreeNodeBase) {
1227          var representative = (VariableTreeNodeBase)firstNode;
1228          representative.Weight = 1.0;
1229          if (variableNodeGroup.Count() > 1) {
1230            var poly = mulSymbol.CreateTreeNode();
1231            for (int p = 0; p < variableNodeGroup.Count(); p++) {
1232              poly.AddSubtree((ISymbolicExpressionTreeNode)representative.Clone());
1233            }
1234            prod.AddSubtree(poly);
1235          } else {
1236            prod.AddSubtree(representative);
1237          }
1238        } else if (firstNode is FactorVariableTreeNode) {
1239          var representative = (FactorVariableTreeNode)firstNode;
1240          foreach (var node in variableNodeGroup.Skip(1).Cast<FactorVariableTreeNode>()) {
1241            for (int j = 0; j < representative.Weights.Length; j++) {
1242              representative.Weights[j] *= node.Weights[j];
1243            }
1244          }
1245          for (int j = 0; j < representative.Weights.Length; j++) {
1246            representative.Weights[j] *= constantProduct;
1247          }
1248          constantProduct = 1.0;
1249          // if the product already contains a factor it is not necessary to multiply a constant below
1250          prod.AddSubtree(representative);
1251        }
1252      }
1253
1254      foreach (var unchangedSubtree in unchangedSubtrees)
1255        prod.AddSubtree(unchangedSubtree);
1256
1257      if (!constantProduct.IsAlmost(1.0)) {
1258        prod.AddSubtree(MakeConstant(constantProduct));
1259      }
1260    }
1261
1262
1263    /// <summary>
1264    /// x => x * -1
1265    /// Is only used in cases where it is not necessary to create new tree nodes. Manipulates x directly.
1266    /// </summary>
1267    /// <param name="x"></param>
1268    /// <returns>-x</returns>
1269    private ISymbolicExpressionTreeNode Negate(ISymbolicExpressionTreeNode x) {
1270      if (IsConstant(x)) {
1271        ((ConstantTreeNode)x).Value *= -1;
1272      } else if (IsVariableBase(x)) {
1273        var variableTree = (VariableTreeNodeBase)x;
1274        variableTree.Weight *= -1.0;
1275      } else if (IsFactor(x)) {
1276        var factorNode = (FactorVariableTreeNode)x;
1277        for (int i = 0; i < factorNode.Weights.Length; i++) factorNode.Weights[i] *= -1;
1278      } else if (IsBinFactor(x)) {
1279        var factorNode = (BinaryFactorVariableTreeNode)x;
1280        factorNode.Weight *= -1;
1281      } else if (IsAddition(x)) {
1282        // (x0 + x1 + .. + xn) * -1 => (-x0 + -x1 + .. + -xn)       
1283        var subtrees = new List<ISymbolicExpressionTreeNode>(x.Subtrees);
1284        while (x.Subtrees.Any()) x.RemoveSubtree(0);
1285        foreach (var subtree in subtrees) {
1286          x.AddSubtree(Negate(subtree));
1287        }
1288      } else if (IsMultiplication(x) || IsDivision(x)) {
1289        // x0 * x1 * .. * xn * -1 => x0 * x1 * .. * -xn
1290        var lastSubTree = x.Subtrees.Last();
1291        x.RemoveSubtree(x.SubtreeCount - 1);
1292        x.AddSubtree(Negate(lastSubTree)); // last is maybe a constant, prefer to negate the constant
1293      } else {
1294        // any other function
1295        return MakeProduct(x, MakeConstant(-1));
1296      }
1297      return x;
1298    }
1299
1300    /// <summary>
1301    /// x => 1/x
1302    /// Must create new tree nodes
1303    /// </summary>
1304    /// <param name="x"></param>
1305    /// <returns></returns>
1306    private ISymbolicExpressionTreeNode Invert(ISymbolicExpressionTreeNode x) {
1307      if (IsConstant(x)) {
1308        return MakeConstant(1.0 / ((ConstantTreeNode)x).Value);
1309      } else if (IsFactor(x)) {
1310        var factorNode = (FactorVariableTreeNode)x;
1311        return MakeFactor(factorNode.Symbol, factorNode.VariableName, factorNode.Weights.Select(w => 1.0 / w));
1312      } else if (IsDivision(x)) {
1313        return MakeFraction(x.GetSubtree(1), x.GetSubtree(0));
1314      } else {
1315        // any other function
1316        return MakeFraction(MakeConstant(1), x);
1317      }
1318    }
1319
1320    private ISymbolicExpressionTreeNode MakeConstant(double value) {
1321      ConstantTreeNode constantTreeNode = (ConstantTreeNode)(constSymbol.CreateTreeNode());
1322      constantTreeNode.Value = value;
1323      return constantTreeNode;
1324    }
1325
1326    private ISymbolicExpressionTreeNode MakeFactor(FactorVariable sy, string variableName, IEnumerable<double> weights) {
1327      var tree = (FactorVariableTreeNode)sy.CreateTreeNode();
1328      tree.VariableName = variableName;
1329      tree.Weights = weights.ToArray();
1330      return tree;
1331    }
1332    private ISymbolicExpressionTreeNode MakeBinFactor(BinaryFactorVariable sy, string variableName, string variableValue, double weight) {
1333      var tree = (BinaryFactorVariableTreeNode)sy.CreateTreeNode();
1334      tree.VariableName = variableName;
1335      tree.VariableValue = variableValue;
1336      tree.Weight = weight;
1337      return tree;
1338    }
1339
1340
1341    #endregion
1342  }
1343}
Note: See TracBrowser for help on using the repository browser.