Free cookie consent management tool by TermsFeed Policy Generator

source: branches/3040_VectorBasedGP/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Converters/VectorTreeSimplifier.cs @ 17596

Last change on this file since 17596 was 17596, checked in by pfleck, 4 years ago

#3040: added subtraction/division simplification for sum and mean symbols by converting them to sums/products.

File size: 74.1 KB
Line 
1#region License Information
2
3/* HeuristicLab
4 * Copyright (C) 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.Linq;
27using HeuristicLab.Common;
28using HeuristicLab.Encodings.SymbolicExpressionTreeEncoding;
29
30namespace HeuristicLab.Problems.DataAnalysis.Symbolic {
31  /// <summary>
32  /// Simplifier for symbolic expressions
33  /// </summary>
34  public class VectorTreeSimplifier {
35    private static readonly Addition addSymbol = new Addition();
36    private static readonly Multiplication mulSymbol = new Multiplication();
37    private static readonly Division divSymbol = new Division();
38    private static readonly Constant constSymbol = new Constant();
39    private static readonly Absolute absSymbol = new Absolute();
40    private static readonly Logarithm logSymbol = new Logarithm();
41    private static readonly Exponential expSymbol = new Exponential();
42    private static readonly Root rootSymbol = new Root();
43    private static readonly Square sqrSymbol = new Square();
44    private static readonly SquareRoot sqrtSymbol = new SquareRoot();
45    private static readonly AnalyticQuotient aqSymbol = new AnalyticQuotient();
46    private static readonly Cube cubeSymbol = new Cube();
47    private static readonly CubeRoot cubeRootSymbol = new CubeRoot();
48    private static readonly Power powSymbol = new Power();
49    private static readonly Sine sineSymbol = new Sine();
50    private static readonly Cosine cosineSymbol = new Cosine();
51    private static readonly Tangent tanSymbol = new Tangent();
52    private static readonly IfThenElse ifThenElseSymbol = new IfThenElse();
53    private static readonly And andSymbol = new And();
54    private static readonly Or orSymbol = new Or();
55    private static readonly Not notSymbol = new Not();
56    private static readonly GreaterThan gtSymbol = new GreaterThan();
57    private static readonly LessThan ltSymbol = new LessThan();
58    private static readonly Integral integralSymbol = new Integral();
59    private static readonly LaggedVariable laggedVariableSymbol = new LaggedVariable();
60    private static readonly TimeLag timeLagSymbol = new TimeLag();
61    private static readonly Sum sumSymbol = new Sum();
62    private static readonly Mean meanSymbol = new Mean();
63    private static readonly Length lengthSymbol = new Length();
64
65    private readonly SymbolicDataAnalysisExpressionTreeVectorInterpreter interpreter;
66    private readonly IDataAnalysisProblemData problemData;
67
68    public VectorTreeSimplifier(SymbolicDataAnalysisExpressionTreeVectorInterpreter interpreter, IDataAnalysisProblemData problemData) {
69      this.interpreter = interpreter;
70      this.problemData = problemData;
71    }
72
73    public ISymbolicExpressionTree Simplify(ISymbolicExpressionTree originalTree) {
74      var clone = (ISymbolicExpressionTreeNode)originalTree.Root.Clone();
75      // macro expand (initially no argument trees)
76      var macroExpandedTree = MacroExpand(clone, clone.GetSubtree(0), new List<ISymbolicExpressionTreeNode>());
77      ISymbolicExpressionTreeNode rootNode = (new ProgramRootSymbol()).CreateTreeNode();
78      rootNode.AddSubtree(GetSimplifiedTree(macroExpandedTree));
79
80#if DEBUG
81      // check that each node is only referenced once
82      var nodes = rootNode.IterateNodesPrefix().ToArray();
83      foreach (var n in nodes) if (nodes.Count(ni => ni == n) > 1) throw new InvalidOperationException();
84#endif
85      return new SymbolicExpressionTree(rootNode);
86    }
87
88    // the argumentTrees list contains already expanded trees used as arguments for invocations
89    private static ISymbolicExpressionTreeNode MacroExpand(ISymbolicExpressionTreeNode root, ISymbolicExpressionTreeNode node,
90      IList<ISymbolicExpressionTreeNode> argumentTrees) {
91      List<ISymbolicExpressionTreeNode> subtrees = new List<ISymbolicExpressionTreeNode>(node.Subtrees);
92      while (node.SubtreeCount > 0) node.RemoveSubtree(0);
93      if (node.Symbol is InvokeFunction) {
94        var invokeSym = node.Symbol as InvokeFunction;
95        var defunNode = FindFunctionDefinition(root, invokeSym.FunctionName);
96        var macroExpandedArguments = new List<ISymbolicExpressionTreeNode>();
97        foreach (var subtree in subtrees) {
98          macroExpandedArguments.Add(MacroExpand(root, subtree, argumentTrees));
99        }
100        return MacroExpand(root, defunNode, macroExpandedArguments);
101      } else if (node.Symbol is Argument) {
102        var argSym = node.Symbol as Argument;
103        // return the correct argument sub-tree (already macro-expanded)
104        return (SymbolicExpressionTreeNode)argumentTrees[argSym.ArgumentIndex].Clone();
105      } else {
106        // recursive application
107        foreach (var subtree in subtrees) {
108          node.AddSubtree(MacroExpand(root, subtree, argumentTrees));
109        }
110        return node;
111      }
112    }
113
114    private static ISymbolicExpressionTreeNode FindFunctionDefinition(ISymbolicExpressionTreeNode root, string functionName) {
115      foreach (var subtree in root.Subtrees.OfType<DefunTreeNode>()) {
116        if (subtree.FunctionName == functionName) return subtree.GetSubtree(0);
117      }
118
119      throw new ArgumentException("Definition of function " + functionName + " not found.");
120    }
121
122    #region symbol predicates
123
124    // arithmetic
125    private static bool IsDivision(ISymbolicExpressionTreeNode node) {
126      return node.Symbol is Division;
127    }
128
129    private static bool IsMultiplication(ISymbolicExpressionTreeNode node) {
130      return node.Symbol is Multiplication;
131    }
132
133    private static bool IsSubtraction(ISymbolicExpressionTreeNode node) {
134      return node.Symbol is Subtraction;
135    }
136
137    private static bool IsAddition(ISymbolicExpressionTreeNode node) {
138      return node.Symbol is Addition;
139    }
140
141    private static bool IsAverage(ISymbolicExpressionTreeNode node) {
142      return node.Symbol is Average;
143    }
144
145    private static bool IsAbsolute(ISymbolicExpressionTreeNode node) {
146      return node.Symbol is Absolute;
147    }
148
149    // exponential
150    private static bool IsLog(ISymbolicExpressionTreeNode node) {
151      return node.Symbol is Logarithm;
152    }
153
154    private static bool IsExp(ISymbolicExpressionTreeNode node) {
155      return node.Symbol is Exponential;
156    }
157
158    private static bool IsRoot(ISymbolicExpressionTreeNode node) {
159      return node.Symbol is Root;
160    }
161
162    private static bool IsSquare(ISymbolicExpressionTreeNode node) {
163      return node.Symbol is Square;
164    }
165
166    private static bool IsSquareRoot(ISymbolicExpressionTreeNode node) {
167      return node.Symbol is SquareRoot;
168    }
169
170    private static bool IsCube(ISymbolicExpressionTreeNode node) {
171      return node.Symbol is Cube;
172    }
173
174    private static bool IsCubeRoot(ISymbolicExpressionTreeNode node) {
175      return node.Symbol is CubeRoot;
176    }
177
178    private static bool IsPower(ISymbolicExpressionTreeNode node) {
179      return node.Symbol is Power;
180    }
181
182    // trigonometric
183    private static bool IsSine(ISymbolicExpressionTreeNode node) {
184      return node.Symbol is Sine;
185    }
186
187    private static bool IsCosine(ISymbolicExpressionTreeNode node) {
188      return node.Symbol is Cosine;
189    }
190
191    private static bool IsTangent(ISymbolicExpressionTreeNode node) {
192      return node.Symbol is Tangent;
193    }
194
195    private static bool IsAnalyticalQuotient(ISymbolicExpressionTreeNode node) {
196      return node.Symbol is AnalyticQuotient;
197    }
198
199    // boolean
200    private static bool IsIfThenElse(ISymbolicExpressionTreeNode node) {
201      return node.Symbol is IfThenElse;
202    }
203
204    private static bool IsAnd(ISymbolicExpressionTreeNode node) {
205      return node.Symbol is And;
206    }
207
208    private static bool IsOr(ISymbolicExpressionTreeNode node) {
209      return node.Symbol is Or;
210    }
211
212    private static bool IsNot(ISymbolicExpressionTreeNode node) {
213      return node.Symbol is Not;
214    }
215
216    // comparison
217    private static bool IsGreaterThan(ISymbolicExpressionTreeNode node) {
218      return node.Symbol is GreaterThan;
219    }
220
221    private static bool IsLessThan(ISymbolicExpressionTreeNode node) {
222      return node.Symbol is LessThan;
223    }
224
225    private static bool IsBoolean(ISymbolicExpressionTreeNode node) {
226      return
227        node.Symbol is GreaterThan ||
228        node.Symbol is LessThan ||
229        node.Symbol is And ||
230        node.Symbol is Or;
231    }
232
233    // terminals
234    private static bool IsVariable(ISymbolicExpressionTreeNode node) {
235      return node.Symbol is Variable;
236    }
237
238    private static bool IsVariableBase(ISymbolicExpressionTreeNode node) {
239      return node is VariableTreeNodeBase;
240    }
241
242    private static bool IsFactor(ISymbolicExpressionTreeNode node) {
243      return node is FactorVariableTreeNode;
244    }
245
246    private static bool IsBinFactor(ISymbolicExpressionTreeNode node) {
247      return node is BinaryFactorVariableTreeNode;
248    }
249
250    private static bool IsConstant(ISymbolicExpressionTreeNode node) {
251      return node.Symbol is Constant;
252    }
253
254    // dynamic
255    private static bool IsTimeLag(ISymbolicExpressionTreeNode node) {
256      return node.Symbol is TimeLag;
257    }
258
259    private static bool IsIntegral(ISymbolicExpressionTreeNode node) {
260      return node.Symbol is Integral;
261    }
262
263    private static bool IsSum(ISymbolicExpressionTreeNode node) {
264      return node.Symbol is Sum;
265    }
266
267    private static bool IsMean(ISymbolicExpressionTreeNode node) {
268      return node.Symbol is Mean;
269    }
270    #endregion
271
272    #region type predicates
273    public bool IsScalarNode(ISymbolicExpressionTreeNode node) {
274      var results = interpreter.EvaluateNode(node, problemData.Dataset, problemData.TrainingIndices);
275      return results.All(r => r.IsScalar);
276    }
277    public bool IsVectorNode(ISymbolicExpressionTreeNode node) {
278      var results = interpreter.EvaluateNode(node, problemData.Dataset, problemData.TrainingIndices);
279      return results.All(r => r.IsVector);
280    }
281    #endregion
282
283    /// <summary>
284    /// Creates a new simplified tree
285    /// </summary>
286    /// <param name="original"></param>
287    /// <returns></returns>
288    public ISymbolicExpressionTreeNode GetSimplifiedTree(ISymbolicExpressionTreeNode original) {
289      if (IsConstant(original) || IsVariableBase(original)) {
290        return (ISymbolicExpressionTreeNode)original.Clone();
291      } else if (IsAbsolute(original)) {
292        return SimplifyAbsolute(original);
293      } else if (IsAddition(original)) {
294        return SimplifyAddition(original);
295      } else if (IsSubtraction(original)) {
296        return SimplifySubtraction(original);
297      } else if (IsMultiplication(original)) {
298        return SimplifyMultiplication(original);
299      } else if (IsDivision(original)) {
300        return SimplifyDivision(original);
301      } else if (IsAverage(original)) {
302        return SimplifyAverage(original);
303      } else if (IsLog(original)) {
304        return SimplifyLog(original);
305      } else if (IsExp(original)) {
306        return SimplifyExp(original);
307      } else if (IsSquare(original)) {
308        return SimplifySquare(original);
309      } else if (IsSquareRoot(original)) {
310        return SimplifySquareRoot(original);
311      } else if (IsCube(original)) {
312        return SimplifyCube(original);
313      } else if (IsCubeRoot(original)) {
314        return SimplifyCubeRoot(original);
315      } else if (IsPower(original)) {
316        return SimplifyPower(original);
317      } else if (IsRoot(original)) {
318        return SimplifyRoot(original);
319      } else if (IsSine(original)) {
320        return SimplifySine(original);
321      } else if (IsCosine(original)) {
322        return SimplifyCosine(original);
323      } else if (IsTangent(original)) {
324        return SimplifyTangent(original);
325      } else if (IsAnalyticalQuotient(original)) {
326        return SimplifyAnalyticalQuotient(original);
327      } else if (IsIfThenElse(original)) {
328        return SimplifyIfThenElse(original);
329      } else if (IsGreaterThan(original)) {
330        return SimplifyGreaterThan(original);
331      } else if (IsLessThan(original)) {
332        return SimplifyLessThan(original);
333      } else if (IsAnd(original)) {
334        return SimplifyAnd(original);
335      } else if (IsOr(original)) {
336        return SimplifyOr(original);
337      } else if (IsNot(original)) {
338        return SimplifyNot(original);
339      } else if (IsTimeLag(original)) {
340        return SimplifyTimeLag(original);
341      } else if (IsIntegral(original)) {
342        return SimplifyIntegral(original);
343      } else if (IsSum(original)) {
344        return SimplifySumAggregation(original);
345      } else if (IsMean(original)) {
346        return SimplifyMeanAggregation(original);
347      } else {
348        return SimplifyAny(original);
349      }
350    }
351
352    #region specific simplification routines
353
354    private ISymbolicExpressionTreeNode SimplifyAny(ISymbolicExpressionTreeNode original) {
355      // can't simplify this function but simplify all subtrees
356      List<ISymbolicExpressionTreeNode> subtrees = new List<ISymbolicExpressionTreeNode>(original.Subtrees);
357      while (original.Subtrees.Count() > 0) original.RemoveSubtree(0);
358      var clone = (SymbolicExpressionTreeNode)original.Clone();
359      List<ISymbolicExpressionTreeNode> simplifiedSubtrees = new List<ISymbolicExpressionTreeNode>();
360      foreach (var subtree in subtrees) {
361        simplifiedSubtrees.Add(GetSimplifiedTree(subtree));
362        original.AddSubtree(subtree);
363      }
364      foreach (var simplifiedSubtree in simplifiedSubtrees) {
365        clone.AddSubtree(simplifiedSubtree);
366      }
367      if (simplifiedSubtrees.TrueForAll(t => IsConstant(t))) {
368        SimplifyConstantExpression(clone);
369      }
370      return clone;
371    }
372
373    private ISymbolicExpressionTreeNode SimplifyConstantExpression(ISymbolicExpressionTreeNode original) {
374      // not yet implemented
375      return original;
376    }
377
378    private ISymbolicExpressionTreeNode SimplifyAverage(ISymbolicExpressionTreeNode original) {
379      if (original.Subtrees.Count() == 1) {
380        return GetSimplifiedTree(original.GetSubtree(0));
381      } else {
382        // simplify expressions x0..xn
383        // make sum(x0..xn) / n
384        var sum = original.Subtrees
385          .Select(GetSimplifiedTree)
386          .Aggregate(MakeSum);
387        return MakeFraction(sum, MakeConstant(original.Subtrees.Count()));
388      }
389    }
390
391    private ISymbolicExpressionTreeNode SimplifyDivision(ISymbolicExpressionTreeNode original) {
392      if (original.Subtrees.Count() == 1) {
393        return Invert(GetSimplifiedTree(original.GetSubtree(0)));
394      } else {
395        // simplify expressions x0..xn
396        // make multiplication (x0 * 1/(x1 * x1 * .. * xn))
397        var first = original.GetSubtree(0);
398        var second = original.GetSubtree(1);
399        var remaining = original.Subtrees.Skip(2);
400        return
401          MakeProduct(GetSimplifiedTree(first),
402            Invert(remaining.Aggregate(GetSimplifiedTree(second), (a, b) => MakeProduct(a, GetSimplifiedTree(b)))));
403      }
404    }
405
406    private ISymbolicExpressionTreeNode SimplifyMultiplication(ISymbolicExpressionTreeNode original) {
407      if (original.Subtrees.Count() == 1) {
408        return GetSimplifiedTree(original.GetSubtree(0));
409      } else {
410        return original.Subtrees
411          .Select(GetSimplifiedTree)
412          .Aggregate(MakeProduct);
413      }
414    }
415
416    private ISymbolicExpressionTreeNode SimplifySubtraction(ISymbolicExpressionTreeNode original) {
417      if (original.Subtrees.Count() == 1) {
418        return Negate(GetSimplifiedTree(original.GetSubtree(0)));
419      } else {
420        // simplify expressions x0..xn
421        // make addition (x0,-x1..-xn)
422        var first = original.Subtrees.First();
423        var remaining = original.Subtrees.Skip(1);
424        return remaining.Aggregate(GetSimplifiedTree(first), (a, b) => MakeSum(a, Negate(GetSimplifiedTree(b))));
425      }
426    }
427
428    private ISymbolicExpressionTreeNode SimplifyAddition(ISymbolicExpressionTreeNode original) {
429      if (original.Subtrees.Count() == 1) {
430        return GetSimplifiedTree(original.GetSubtree(0));
431      } else {
432        // simplify expression x0..xn
433        // make addition (x0..xn)
434        return original.Subtrees
435          .Select(GetSimplifiedTree)
436          .Aggregate(MakeSum);
437      }
438    }
439
440    private ISymbolicExpressionTreeNode SimplifyAbsolute(ISymbolicExpressionTreeNode original) {
441      return MakeAbs(GetSimplifiedTree(original.GetSubtree(0)));
442    }
443
444    private ISymbolicExpressionTreeNode SimplifyNot(ISymbolicExpressionTreeNode original) {
445      return MakeNot(GetSimplifiedTree(original.GetSubtree(0)));
446    }
447
448    private ISymbolicExpressionTreeNode SimplifyOr(ISymbolicExpressionTreeNode original) {
449      return original.Subtrees
450        .Select(GetSimplifiedTree)
451        .Aggregate(MakeOr);
452    }
453
454    private ISymbolicExpressionTreeNode SimplifyAnd(ISymbolicExpressionTreeNode original) {
455      return original.Subtrees
456        .Select(GetSimplifiedTree)
457        .Aggregate(MakeAnd);
458    }
459
460    private ISymbolicExpressionTreeNode SimplifyLessThan(ISymbolicExpressionTreeNode original) {
461      return MakeLessThan(GetSimplifiedTree(original.GetSubtree(0)), GetSimplifiedTree(original.GetSubtree(1)));
462    }
463
464    private ISymbolicExpressionTreeNode SimplifyGreaterThan(ISymbolicExpressionTreeNode original) {
465      return MakeGreaterThan(GetSimplifiedTree(original.GetSubtree(0)), GetSimplifiedTree(original.GetSubtree(1)));
466    }
467
468    private ISymbolicExpressionTreeNode SimplifyIfThenElse(ISymbolicExpressionTreeNode original) {
469      return MakeIfThenElse(GetSimplifiedTree(original.GetSubtree(0)), GetSimplifiedTree(original.GetSubtree(1)),
470        GetSimplifiedTree(original.GetSubtree(2)));
471    }
472
473    private ISymbolicExpressionTreeNode SimplifyTangent(ISymbolicExpressionTreeNode original) {
474      return MakeTangent(GetSimplifiedTree(original.GetSubtree(0)));
475    }
476
477    private ISymbolicExpressionTreeNode SimplifyCosine(ISymbolicExpressionTreeNode original) {
478      return MakeCosine(GetSimplifiedTree(original.GetSubtree(0)));
479    }
480
481    private ISymbolicExpressionTreeNode SimplifySine(ISymbolicExpressionTreeNode original) {
482      return MakeSine(GetSimplifiedTree(original.GetSubtree(0)));
483    }
484
485    private ISymbolicExpressionTreeNode SimplifyExp(ISymbolicExpressionTreeNode original) {
486      return MakeExp(GetSimplifiedTree(original.GetSubtree(0)));
487    }
488
489    private ISymbolicExpressionTreeNode SimplifySquare(ISymbolicExpressionTreeNode original) {
490      return MakeSquare(GetSimplifiedTree(original.GetSubtree(0)));
491    }
492
493    private ISymbolicExpressionTreeNode SimplifySquareRoot(ISymbolicExpressionTreeNode original) {
494      return MakeSquareRoot(GetSimplifiedTree(original.GetSubtree(0)));
495    }
496    private ISymbolicExpressionTreeNode SimplifyCube(ISymbolicExpressionTreeNode original) {
497      return MakeCube(GetSimplifiedTree(original.GetSubtree(0)));
498    }
499
500    private ISymbolicExpressionTreeNode SimplifyCubeRoot(ISymbolicExpressionTreeNode original) {
501      return MakeCubeRoot(GetSimplifiedTree(original.GetSubtree(0)));
502    }
503
504    private ISymbolicExpressionTreeNode SimplifyLog(ISymbolicExpressionTreeNode original) {
505      return MakeLog(GetSimplifiedTree(original.GetSubtree(0)));
506    }
507
508    private ISymbolicExpressionTreeNode SimplifyRoot(ISymbolicExpressionTreeNode original) {
509      return MakeRoot(GetSimplifiedTree(original.GetSubtree(0)), GetSimplifiedTree(original.GetSubtree(1)));
510    }
511
512    private ISymbolicExpressionTreeNode SimplifyPower(ISymbolicExpressionTreeNode original) {
513      return MakePower(GetSimplifiedTree(original.GetSubtree(0)), GetSimplifiedTree(original.GetSubtree(1)));
514    }
515
516    private ISymbolicExpressionTreeNode SimplifyAnalyticalQuotient(ISymbolicExpressionTreeNode original) {
517      return MakeAnalyticalQuotient(GetSimplifiedTree(original.GetSubtree(0)), GetSimplifiedTree(original.GetSubtree(1)));
518    }
519
520    private ISymbolicExpressionTreeNode SimplifyTimeLag(ISymbolicExpressionTreeNode original) {
521      var laggedTreeNode = original as ILaggedTreeNode;
522      var simplifiedSubtree = GetSimplifiedTree(original.GetSubtree(0));
523      if (!ContainsVariableCondition(simplifiedSubtree)) {
524        return AddLagToDynamicNodes(simplifiedSubtree, laggedTreeNode.Lag);
525      } else {
526        return MakeTimeLag(simplifiedSubtree, laggedTreeNode.Lag);
527      }
528    }
529
530    private ISymbolicExpressionTreeNode SimplifyIntegral(ISymbolicExpressionTreeNode original) {
531      var laggedTreeNode = original as ILaggedTreeNode;
532      var simplifiedSubtree = GetSimplifiedTree(original.GetSubtree(0));
533      if (IsConstant(simplifiedSubtree)) {
534        return GetSimplifiedTree(MakeProduct(simplifiedSubtree, MakeConstant(-laggedTreeNode.Lag)));
535      } else {
536        return MakeIntegral(simplifiedSubtree, laggedTreeNode.Lag);
537      }
538    }
539
540    private ISymbolicExpressionTreeNode SimplifySumAggregation(ISymbolicExpressionTreeNode original) {
541      return MakeSumAggregation(GetSimplifiedTree(original.GetSubtree(0)));
542    }
543
544    private ISymbolicExpressionTreeNode SimplifyMeanAggregation(ISymbolicExpressionTreeNode original) {
545      return MakeMeanAggregation(GetSimplifiedTree(original.GetSubtree(0)));
546    }
547
548    #endregion
549
550    #region low level tree restructuring
551
552    private ISymbolicExpressionTreeNode MakeTimeLag(ISymbolicExpressionTreeNode subtree, int lag) {
553      if (lag == 0) return subtree;
554      if (IsConstant(subtree)) return subtree;
555      var lagNode = (LaggedTreeNode)timeLagSymbol.CreateTreeNode();
556      lagNode.Lag = lag;
557      lagNode.AddSubtree(subtree);
558      return lagNode;
559    }
560
561    private ISymbolicExpressionTreeNode MakeIntegral(ISymbolicExpressionTreeNode subtree, int lag) {
562      if (lag == 0) return subtree;
563      else if (lag == -1 || lag == 1) {
564        return MakeSum(subtree, AddLagToDynamicNodes((ISymbolicExpressionTreeNode)subtree.Clone(), lag));
565      } else {
566        var node = (LaggedTreeNode)integralSymbol.CreateTreeNode();
567        node.Lag = lag;
568        node.AddSubtree(subtree);
569        return node;
570      }
571    }
572
573    private ISymbolicExpressionTreeNode MakeNot(ISymbolicExpressionTreeNode t) {
574      if (IsConstant(t)) {
575        var constNode = t as ConstantTreeNode;
576        if (constNode.Value > 0) return MakeConstant(-1.0);
577        else return MakeConstant(1.0);
578      } else if (IsNot(t)) {
579        return t.GetSubtree(0);
580      } else if (!IsBoolean(t)) {
581        var gtNode = gtSymbol.CreateTreeNode();
582        gtNode.AddSubtree(t);
583        gtNode.AddSubtree(MakeConstant(0.0));
584        var notNode = notSymbol.CreateTreeNode();
585        notNode.AddSubtree(gtNode);
586        return notNode;
587      } else {
588        var notNode = notSymbol.CreateTreeNode();
589        notNode.AddSubtree(t);
590        return notNode;
591      }
592    }
593
594    private ISymbolicExpressionTreeNode MakeOr(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
595      if (IsConstant(a) && IsConstant(b)) {
596        var constA = a as ConstantTreeNode;
597        var constB = b as ConstantTreeNode;
598        if (constA.Value > 0.0 || constB.Value > 0.0) {
599          return MakeConstant(1.0);
600        } else {
601          return MakeConstant(-1.0);
602        }
603      } else if (IsConstant(a)) {
604        return MakeOr(b, a);
605      } else if (IsConstant(b)) {
606        var constT = b as ConstantTreeNode;
607        if (constT.Value > 0.0) {
608          // boolean expression is necessarily true
609          return MakeConstant(1.0);
610        } else {
611          // the constant value has no effect on the result of the boolean condition so we can drop the constant term
612          var orNode = orSymbol.CreateTreeNode();
613          orNode.AddSubtree(a);
614          return orNode;
615        }
616      } else {
617        var orNode = orSymbol.CreateTreeNode();
618        orNode.AddSubtree(a);
619        orNode.AddSubtree(b);
620        return orNode;
621      }
622    }
623
624    private ISymbolicExpressionTreeNode MakeAnd(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
625      if (IsConstant(a) && IsConstant(b)) {
626        var constA = a as ConstantTreeNode;
627        var constB = b as ConstantTreeNode;
628        if (constA.Value > 0.0 && constB.Value > 0.0) {
629          return MakeConstant(1.0);
630        } else {
631          return MakeConstant(-1.0);
632        }
633      } else if (IsConstant(a)) {
634        return MakeAnd(b, a);
635      } else if (IsConstant(b)) {
636        var constB = b as ConstantTreeNode;
637        if (constB.Value > 0.0) {
638          // the constant value has no effect on the result of the boolean condition so we can drop the constant term
639          var andNode = andSymbol.CreateTreeNode();
640          andNode.AddSubtree(a);
641          return andNode;
642        } else {
643          // boolean expression is necessarily false
644          return MakeConstant(-1.0);
645        }
646      } else {
647        var andNode = andSymbol.CreateTreeNode();
648        andNode.AddSubtree(a);
649        andNode.AddSubtree(b);
650        return andNode;
651      }
652    }
653
654    private ISymbolicExpressionTreeNode MakeLessThan(ISymbolicExpressionTreeNode leftSide,
655      ISymbolicExpressionTreeNode rightSide) {
656      if (IsConstant(leftSide) && IsConstant(rightSide)) {
657        var lsConst = leftSide as ConstantTreeNode;
658        var rsConst = rightSide as ConstantTreeNode;
659        if (lsConst.Value < rsConst.Value) return MakeConstant(1.0);
660        else return MakeConstant(-1.0);
661      } else {
662        var ltNode = ltSymbol.CreateTreeNode();
663        ltNode.AddSubtree(leftSide);
664        ltNode.AddSubtree(rightSide);
665        return ltNode;
666      }
667    }
668
669    private ISymbolicExpressionTreeNode MakeGreaterThan(ISymbolicExpressionTreeNode leftSide,
670      ISymbolicExpressionTreeNode rightSide) {
671      if (IsConstant(leftSide) && IsConstant(rightSide)) {
672        var lsConst = leftSide as ConstantTreeNode;
673        var rsConst = rightSide as ConstantTreeNode;
674        if (lsConst.Value > rsConst.Value) return MakeConstant(1.0);
675        else return MakeConstant(-1.0);
676      } else {
677        var gtNode = gtSymbol.CreateTreeNode();
678        gtNode.AddSubtree(leftSide);
679        gtNode.AddSubtree(rightSide);
680        return gtNode;
681      }
682    }
683
684    private ISymbolicExpressionTreeNode MakeIfThenElse(ISymbolicExpressionTreeNode condition,
685      ISymbolicExpressionTreeNode trueBranch, ISymbolicExpressionTreeNode falseBranch) {
686      if (IsConstant(condition)) {
687        var constT = condition as ConstantTreeNode;
688        if (constT.Value > 0.0) return trueBranch;
689        else return falseBranch;
690      } else {
691        var ifNode = ifThenElseSymbol.CreateTreeNode();
692        if (IsBoolean(condition)) {
693          ifNode.AddSubtree(condition);
694        } else {
695          var gtNode = gtSymbol.CreateTreeNode();
696          gtNode.AddSubtree(condition);
697          gtNode.AddSubtree(MakeConstant(0.0));
698          ifNode.AddSubtree(gtNode);
699        }
700        ifNode.AddSubtree(trueBranch);
701        ifNode.AddSubtree(falseBranch);
702        return ifNode;
703      }
704    }
705
706    private ISymbolicExpressionTreeNode MakeSine(ISymbolicExpressionTreeNode node) {
707      if (IsConstant(node)) {
708        var constT = node as ConstantTreeNode;
709        return MakeConstant(Math.Sin(constT.Value));
710      } else if (IsFactor(node)) {
711        var factor = node as FactorVariableTreeNode;
712        return MakeFactor(factor.Symbol, factor.VariableName, factor.Weights.Select(Math.Sin));
713      } else if (IsBinFactor(node)) {
714        var binFactor = node as BinaryFactorVariableTreeNode;
715        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Sin(binFactor.Weight));
716      } else {
717        var sineNode = sineSymbol.CreateTreeNode();
718        sineNode.AddSubtree(node);
719        return sineNode;
720      }
721    }
722
723    private ISymbolicExpressionTreeNode MakeTangent(ISymbolicExpressionTreeNode node) {
724      if (IsConstant(node)) {
725        var constT = node as ConstantTreeNode;
726        return MakeConstant(Math.Tan(constT.Value));
727      } else if (IsFactor(node)) {
728        var factor = node as FactorVariableTreeNode;
729        return MakeFactor(factor.Symbol, factor.VariableName, factor.Weights.Select(Math.Tan));
730      } else if (IsBinFactor(node)) {
731        var binFactor = node as BinaryFactorVariableTreeNode;
732        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Tan(binFactor.Weight));
733      } else {
734        var tanNode = tanSymbol.CreateTreeNode();
735        tanNode.AddSubtree(node);
736        return tanNode;
737      }
738    }
739
740    private ISymbolicExpressionTreeNode MakeCosine(ISymbolicExpressionTreeNode node) {
741      if (IsConstant(node)) {
742        var constT = node as ConstantTreeNode;
743        return MakeConstant(Math.Cos(constT.Value));
744      } else if (IsFactor(node)) {
745        var factor = node as FactorVariableTreeNode;
746        return MakeFactor(factor.Symbol, factor.VariableName, factor.Weights.Select(Math.Cos));
747      } else if (IsBinFactor(node)) {
748        var binFactor = node as BinaryFactorVariableTreeNode;
749        // cos(0) = 1 see similar case for Exp(binfactor)
750        return MakeSum(MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Cos(binFactor.Weight) - 1),
751          MakeConstant(1.0));
752      } else {
753        var cosNode = cosineSymbol.CreateTreeNode();
754        cosNode.AddSubtree(node);
755        return cosNode;
756      }
757    }
758
759    private ISymbolicExpressionTreeNode MakeExp(ISymbolicExpressionTreeNode node) {
760      if (IsConstant(node)) {
761        var constT = node as ConstantTreeNode;
762        return MakeConstant(Math.Exp(constT.Value));
763      } else if (IsFactor(node)) {
764        var factNode = node as FactorVariableTreeNode;
765        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Exp(w)));
766      } else if (IsBinFactor(node)) {
767        // exp( binfactor w val=a) = if(val=a) exp(w) else exp(0) = binfactor( (exp(w) - 1) val a) + 1
768        var binFactor = node as BinaryFactorVariableTreeNode;
769        return
770          MakeSum(MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Exp(binFactor.Weight) - 1), MakeConstant(1.0));
771      } else if (IsLog(node)) {
772        return node.GetSubtree(0);
773      } else if (IsAddition(node)) {
774        return node.Subtrees.Select(s => MakeExp(s)).Aggregate((s, t) => MakeProduct(s, t));
775      } else if (IsSubtraction(node)) {
776        return node.Subtrees.Select(s => MakeExp(s)).Aggregate((s, t) => MakeProduct(s, Negate(t)));
777      } else {
778        var expNode = expSymbol.CreateTreeNode();
779        expNode.AddSubtree(node);
780        return expNode;
781      }
782    }
783    private ISymbolicExpressionTreeNode MakeLog(ISymbolicExpressionTreeNode node) {
784      if (IsConstant(node)) {
785        var constT = node as ConstantTreeNode;
786        return MakeConstant(Math.Log(constT.Value));
787      } else if (IsFactor(node)) {
788        var factNode = node as FactorVariableTreeNode;
789        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Log(w)));
790      } else if (IsExp(node)) {
791        return node.GetSubtree(0);
792      } else if (IsSquareRoot(node)) {
793        return MakeFraction(MakeLog(node.GetSubtree(0)), MakeConstant(2.0));
794      } else {
795        var logNode = logSymbol.CreateTreeNode();
796        logNode.AddSubtree(node);
797        return logNode;
798      }
799    }
800
801    private ISymbolicExpressionTreeNode MakeSquare(ISymbolicExpressionTreeNode node) {
802      if (IsConstant(node)) {
803        var constT = node as ConstantTreeNode;
804        return MakeConstant(constT.Value * constT.Value);
805      } else if (IsFactor(node)) {
806        var factNode = node as FactorVariableTreeNode;
807        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => w * w));
808      } else if (IsBinFactor(node)) {
809        var binFactor = node as BinaryFactorVariableTreeNode;
810        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, binFactor.Weight * binFactor.Weight);
811      } else if (IsSquareRoot(node)) {
812        return node.GetSubtree(0);
813      } else if (IsMultiplication(node)) {
814        // sqr( x * y ) = sqr(x) * sqr(y)
815        var mulNode = mulSymbol.CreateTreeNode();
816        foreach (var subtree in node.Subtrees) {
817          mulNode.AddSubtree(MakeSquare(subtree));
818        }
819        return mulNode;
820      } else if (IsAbsolute(node)) {
821        return MakeSquare(node.GetSubtree(0)); // sqr(abs(x)) = sqr(x)
822      } else if (IsExp(node)) {
823        return MakeExp(MakeProduct(node.GetSubtree(0), MakeConstant(2.0))); // sqr(exp(x)) = exp(2x)
824      } else if (IsCube(node)) {
825        return MakePower(node.GetSubtree(0), MakeConstant(6));
826      } else {
827        var sqrNode = sqrSymbol.CreateTreeNode();
828        sqrNode.AddSubtree(node);
829        return sqrNode;
830      }
831    }
832
833    private ISymbolicExpressionTreeNode MakeCube(ISymbolicExpressionTreeNode node) {
834      if (IsConstant(node)) {
835        var constT = node as ConstantTreeNode;
836        return MakeConstant(constT.Value * constT.Value * constT.Value);
837      } else if (IsFactor(node)) {
838        var factNode = node as FactorVariableTreeNode;
839        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => w * w * w));
840      } else if (IsBinFactor(node)) {
841        var binFactor = node as BinaryFactorVariableTreeNode;
842        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, binFactor.Weight * binFactor.Weight * binFactor.Weight);
843      } else if (IsCubeRoot(node)) {
844        return node.GetSubtree(0); // NOTE: not really accurate because cuberoot(x) with negative x is evaluated to NaN and after this simplification we evaluate as x
845      } else if (IsExp(node)) {
846        return MakeExp(MakeProduct(node.GetSubtree(0), MakeConstant(3)));
847      } else if (IsSquare(node)) {
848        return MakePower(node.GetSubtree(0), MakeConstant(6));
849      } else {
850        var cubeNode = cubeSymbol.CreateTreeNode();
851        cubeNode.AddSubtree(node);
852        return cubeNode;
853      }
854    }
855
856    private ISymbolicExpressionTreeNode MakeAbs(ISymbolicExpressionTreeNode node) {
857      if (IsConstant(node)) {
858        var constT = node as ConstantTreeNode;
859        return MakeConstant(Math.Abs(constT.Value));
860      } else if (IsFactor(node)) {
861        var factNode = node as FactorVariableTreeNode;
862        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Abs(w)));
863      } else if (IsBinFactor(node)) {
864        var binFactor = node as BinaryFactorVariableTreeNode;
865        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Abs(binFactor.Weight));
866      } else if (IsSquare(node) || IsExp(node) || IsSquareRoot(node) || IsCubeRoot(node)) {
867        return node; // abs(sqr(x)) = sqr(x), abs(exp(x)) = exp(x) ...
868      } else if (IsMultiplication(node)) {
869        var mul = mulSymbol.CreateTreeNode();
870        foreach (var st in node.Subtrees) {
871          mul.AddSubtree(MakeAbs(st));
872        }
873        return mul;
874      } else if (IsDivision(node)) {
875        var div = divSymbol.CreateTreeNode();
876        foreach (var st in node.Subtrees) {
877          div.AddSubtree(MakeAbs(st));
878        }
879        return div;
880      } else {
881        var absNode = absSymbol.CreateTreeNode();
882        absNode.AddSubtree(node);
883        return absNode;
884      }
885    }
886
887    // constant folding only
888    private ISymbolicExpressionTreeNode MakeAnalyticalQuotient(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
889      if (IsConstant(b)) {
890        var c = b as ConstantTreeNode;
891        return MakeFraction(a, MakeConstant(Math.Sqrt(1.0 + c.Value * c.Value)));
892      } else if (IsFactor(b)) {
893        var factNode = b as FactorVariableTreeNode;
894        return MakeFraction(a, MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Sqrt(1.0 + w * w))));
895      } else if (IsBinFactor(b)) {
896        var binFactor = b as BinaryFactorVariableTreeNode;
897        return MakeFraction(a, MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Sqrt(1.0 + binFactor.Weight * binFactor.Weight)));
898      } else {
899        var aqNode = aqSymbol.CreateTreeNode();
900        aqNode.AddSubtree(a);
901        aqNode.AddSubtree(b);
902        return aqNode;
903      }
904    }
905
906    private ISymbolicExpressionTreeNode MakeSquareRoot(ISymbolicExpressionTreeNode node) {
907      if (IsConstant(node)) {
908        var constT = node as ConstantTreeNode;
909        return MakeConstant(Math.Sqrt(constT.Value));
910      } else if (IsFactor(node)) {
911        var factNode = node as FactorVariableTreeNode;
912        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Sqrt(w)));
913      } else if (IsBinFactor(node)) {
914        var binFactor = node as BinaryFactorVariableTreeNode;
915        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Sqrt(binFactor.Weight));
916      } else if (IsSquare(node)) {
917        return node.GetSubtree(0); // NOTE: not really accurate because sqrt(x) with negative x is evaluated to NaN and after this simplification we evaluate as x
918      } else {
919        var sqrtNode = sqrtSymbol.CreateTreeNode();
920        sqrtNode.AddSubtree(node);
921        return sqrtNode;
922      }
923    }
924
925    private ISymbolicExpressionTreeNode MakeCubeRoot(ISymbolicExpressionTreeNode node) {
926      if (IsConstant(node)) {
927        var constT = node as ConstantTreeNode;
928        return MakeConstant(Math.Pow(constT.Value, 1.0 / 3.0));
929      } else if (IsFactor(node)) {
930        var factNode = node as FactorVariableTreeNode;
931        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Pow(w, 1.0 / 3.0)));
932      } else if (IsBinFactor(node)) {
933        var binFactor = node as BinaryFactorVariableTreeNode;
934        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Sqrt(Math.Pow(binFactor.Weight, 1.0 / 3.0)));
935      } else if (IsCube(node)) {
936        return node.GetSubtree(0);
937      } else {
938        var cubeRootNode = cubeRootSymbol.CreateTreeNode();
939        cubeRootNode.AddSubtree(node);
940        return cubeRootNode;
941      }
942    }
943
944    private ISymbolicExpressionTreeNode MakeRoot(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
945      if (IsConstant(a) && IsConstant(b)) {
946        var constA = a as ConstantTreeNode;
947        var constB = b as ConstantTreeNode;
948        return MakeConstant(Math.Pow(constA.Value, 1.0 / Math.Round(constB.Value)));
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,
953          factNode.Weights.Select(w => Math.Pow(w, 1.0 / Math.Round(constNode.Value))));
954      } else if (IsBinFactor(a) && IsConstant(b)) {
955        var binFactor = a as BinaryFactorVariableTreeNode;
956        var constNode = b as ConstantTreeNode;
957        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Pow(binFactor.Weight, 1.0 / Math.Round(constNode.Value)));
958      } else if (IsConstant(a) && IsFactor(b)) {
959        var constNode = a as ConstantTreeNode;
960        var factNode = b as FactorVariableTreeNode;
961        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Pow(constNode.Value, 1.0 / Math.Round(w))));
962      } else if (IsConstant(a) && IsBinFactor(b)) {
963        var constNode = a as ConstantTreeNode;
964        var factNode = b as BinaryFactorVariableTreeNode;
965        return MakeBinFactor(factNode.Symbol, factNode.VariableName, factNode.VariableValue, Math.Pow(constNode.Value, 1.0 / Math.Round(factNode.Weight)));
966      } else if (IsFactor(a) && IsFactor(b) && AreSameTypeAndVariable(a, b)) {
967        var node0 = a as FactorVariableTreeNode;
968        var node1 = b as FactorVariableTreeNode;
969        return MakeFactor(node0.Symbol, node0.VariableName, node0.Weights.Zip(node1.Weights, (u, v) => Math.Pow(u, 1.0 / Math.Round(v))));
970      } else if (IsConstant(b)) {
971        var constB = b as ConstantTreeNode;
972        var constBValue = Math.Round(constB.Value);
973        if (constBValue.IsAlmost(1.0)) {
974          return a;
975        } else if (constBValue.IsAlmost(0.0)) {
976          return MakeConstant(1.0);
977        } else if (constBValue.IsAlmost(-1.0)) {
978          return MakeFraction(MakeConstant(1.0), a);
979        } else if (constBValue < 0) {
980          var rootNode = rootSymbol.CreateTreeNode();
981          rootNode.AddSubtree(a);
982          rootNode.AddSubtree(MakeConstant(-1.0 * constBValue));
983          return MakeFraction(MakeConstant(1.0), rootNode);
984        } else {
985          var rootNode = rootSymbol.CreateTreeNode();
986          rootNode.AddSubtree(a);
987          rootNode.AddSubtree(MakeConstant(constBValue));
988          return rootNode;
989        }
990      } else {
991        var rootNode = rootSymbol.CreateTreeNode();
992        rootNode.AddSubtree(a);
993        rootNode.AddSubtree(b);
994        return rootNode;
995      }
996    }
997
998
999    private ISymbolicExpressionTreeNode MakePower(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
1000      if (IsConstant(a) && IsConstant(b)) {
1001        var constA = a as ConstantTreeNode;
1002        var constB = b as ConstantTreeNode;
1003        return MakeConstant(Math.Pow(constA.Value, Math.Round(constB.Value)));
1004      } else if (IsFactor(a) && IsConstant(b)) {
1005        var factNode = a as FactorVariableTreeNode;
1006        var constNode = b as ConstantTreeNode;
1007        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Pow(w, Math.Round(constNode.Value))));
1008      } else if (IsBinFactor(a) && IsConstant(b)) {
1009        var binFactor = a as BinaryFactorVariableTreeNode;
1010        var constNode = b as ConstantTreeNode;
1011        return MakeBinFactor(binFactor.Symbol, binFactor.VariableName, binFactor.VariableValue, Math.Pow(binFactor.Weight, Math.Round(constNode.Value)));
1012      } else if (IsConstant(a) && IsFactor(b)) {
1013        var constNode = a as ConstantTreeNode;
1014        var factNode = b as FactorVariableTreeNode;
1015        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => Math.Pow(constNode.Value, Math.Round(w))));
1016      } else if (IsConstant(a) && IsBinFactor(b)) {
1017        var constNode = a as ConstantTreeNode;
1018        var factNode = b as BinaryFactorVariableTreeNode;
1019        return MakeBinFactor(factNode.Symbol, factNode.VariableName, factNode.VariableValue, Math.Pow(constNode.Value, Math.Round(factNode.Weight)));
1020      } else if (IsFactor(a) && IsFactor(b) && AreSameTypeAndVariable(a, b)) {
1021        var node0 = a as FactorVariableTreeNode;
1022        var node1 = b as FactorVariableTreeNode;
1023        return MakeFactor(node0.Symbol, node0.VariableName, node0.Weights.Zip(node1.Weights, (u, v) => Math.Pow(u, Math.Round(v))));
1024      } else if (IsConstant(b)) {
1025        var constB = b as ConstantTreeNode;
1026        double exponent = Math.Round(constB.Value);
1027        if (exponent.IsAlmost(0.0)) {
1028          return MakeConstant(1.0);
1029        } else if (exponent.IsAlmost(1.0)) {
1030          return a;
1031        } else if (exponent.IsAlmost(-1.0)) {
1032          return MakeFraction(MakeConstant(1.0), a);
1033        } else if (exponent < 0) {
1034          var powNode = powSymbol.CreateTreeNode();
1035          powNode.AddSubtree(a);
1036          powNode.AddSubtree(MakeConstant(-1.0 * exponent));
1037          return MakeFraction(MakeConstant(1.0), powNode);
1038        } else {
1039          var powNode = powSymbol.CreateTreeNode();
1040          powNode.AddSubtree(a);
1041          powNode.AddSubtree(MakeConstant(exponent));
1042          return powNode;
1043        }
1044      } else {
1045        var powNode = powSymbol.CreateTreeNode();
1046        powNode.AddSubtree(a);
1047        powNode.AddSubtree(b);
1048        return powNode;
1049      }
1050    }
1051
1052
1053    // MakeFraction, MakeProduct and MakeSum take two already simplified trees and create a new simplified tree
1054    private ISymbolicExpressionTreeNode MakeFraction(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
1055      if (IsConstant(a) && IsConstant(b)) {
1056        // fold constants
1057        return MakeConstant(((ConstantTreeNode)a).Value / ((ConstantTreeNode)b).Value);
1058      } else if ((IsConstant(a) && !((ConstantTreeNode)a).Value.IsAlmost(1.0))) {
1059        return MakeFraction(MakeConstant(1.0), MakeProduct(b, Invert(a)));
1060      } else if (IsVariableBase(a) && IsConstant(b)) {
1061        // merge constant values into variable weights
1062        var constB = ((ConstantTreeNode)b).Value;
1063        ((VariableTreeNodeBase)a).Weight /= constB;
1064        return a;
1065      } else if (IsFactor(a) && IsConstant(b)) {
1066        var factNode = a as FactorVariableTreeNode;
1067        var constNode = b as ConstantTreeNode;
1068        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select(w => w / constNode.Value));
1069      } else if (IsBinFactor(a) && IsConstant(b)) {
1070        var factNode = a as BinaryFactorVariableTreeNode;
1071        var constNode = b as ConstantTreeNode;
1072        return MakeBinFactor(factNode.Symbol, factNode.VariableName, factNode.VariableValue, factNode.Weight / constNode.Value);
1073      } else if (IsFactor(a) && IsFactor(b) && AreSameTypeAndVariable(a, b)) {
1074        var node0 = a as FactorVariableTreeNode;
1075        var node1 = b as FactorVariableTreeNode;
1076        return MakeFactor(node0.Symbol, node0.VariableName, node0.Weights.Zip(node1.Weights, (u, v) => u / v));
1077      } else if (IsFactor(a) && IsBinFactor(b) && ((IVariableTreeNode)a).VariableName == ((IVariableTreeNode)b).VariableName) {
1078        var node0 = a as FactorVariableTreeNode;
1079        var node1 = b as BinaryFactorVariableTreeNode;
1080        var varValues = node0.Symbol.GetVariableValues(node0.VariableName).ToArray();
1081        var wi = Array.IndexOf(varValues, node1.VariableValue);
1082        if (wi < 0) throw new ArgumentException();
1083        var newWeighs = new double[varValues.Length];
1084        node0.Weights.CopyTo(newWeighs, 0);
1085        for (int i = 0; i < newWeighs.Length; i++)
1086          if (wi == i) newWeighs[i] /= node1.Weight;
1087          else newWeighs[i] /= 0.0;
1088        return MakeFactor(node0.Symbol, node0.VariableName, newWeighs);
1089      } else if (IsFactor(a)) {
1090        return MakeFraction(MakeConstant(1.0), MakeProduct(b, Invert(a)));
1091      } else if (IsVariableBase(a) && IsVariableBase(b) && AreSameTypeAndVariable(a, b) && !IsBinFactor(b)) {
1092        // cancel variables (not allowed for bin factors because of division by zero)
1093        var aVar = a as VariableTreeNode;
1094        var bVar = b as VariableTreeNode;
1095        return MakeConstant(aVar.Weight / bVar.Weight);
1096      } else if (IsAddition(a) && IsConstant(b)) {
1097        return a.Subtrees
1098          .Select(x => GetSimplifiedTree(x))
1099          .Select(x => MakeFraction(x, GetSimplifiedTree(b)))
1100          .Aggregate((c, d) => MakeSum(c, d));
1101      } else if (IsMultiplication(a) && IsConstant(b)) {
1102        return MakeProduct(a, Invert(b));
1103      } else if (IsDivision(a) && IsConstant(b)) {
1104        // (a1 / a2) / c => (a1 / (a2 * c))
1105        return MakeFraction(a.GetSubtree(0), MakeProduct(a.GetSubtree(1), b));
1106      } else if (IsDivision(a) && IsDivision(b)) {
1107        // (a1 / a2) / (b1 / b2) =>
1108        return MakeFraction(MakeProduct(a.GetSubtree(0), b.GetSubtree(1)), MakeProduct(a.GetSubtree(1), b.GetSubtree(0)));
1109      } else if (IsDivision(a)) {
1110        // (a1 / a2) / b => (a1 / (a2 * b))
1111        return MakeFraction(a.GetSubtree(0), MakeProduct(a.GetSubtree(1), b));
1112      } else if (IsDivision(b)) {
1113        // a / (b1 / b2) => (a * b2) / b1
1114        return MakeFraction(MakeProduct(a, b.GetSubtree(1)), b.GetSubtree(0));
1115      } else if (IsAnalyticalQuotient(a)) {
1116        return MakeAnalyticalQuotient(a.GetSubtree(0), MakeProduct(a.GetSubtree(1), b));
1117      } else {
1118        var div = divSymbol.CreateTreeNode();
1119        div.AddSubtree(a);
1120        div.AddSubtree(b);
1121        return div;
1122      }
1123    }
1124
1125    private ISymbolicExpressionTreeNode MakeSum(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
1126      if (IsConstant(a) && IsConstant(b)) {
1127        // fold constants
1128        ((ConstantTreeNode)a).Value += ((ConstantTreeNode)b).Value;
1129        return a;
1130      } else if (IsConstant(a)) {
1131        // c + x => x + c
1132        // b is not constant => make sure constant is on the right
1133        return MakeSum(b, a);
1134      } else if (IsConstant(b) && ((ConstantTreeNode)b).Value.IsAlmost(0.0)) {
1135        // x + 0 => x
1136        return a;
1137      } else if (IsFactor(a) && IsConstant(b)) {
1138        var factNode = a as FactorVariableTreeNode;
1139        var constNode = b as ConstantTreeNode;
1140        return MakeFactor(factNode.Symbol, factNode.VariableName, factNode.Weights.Select((w) => w + constNode.Value));
1141      } else if (IsFactor(a) && IsFactor(b) && AreSameTypeAndVariable(a, b)) {
1142        var node0 = a as FactorVariableTreeNode;
1143        var node1 = b as FactorVariableTreeNode;
1144        return MakeFactor(node0.Symbol, node0.VariableName, node0.Weights.Zip(node1.Weights, (u, v) => u + v));
1145      } else if (IsBinFactor(a) && IsFactor(b)) {
1146        return MakeSum(b, a);
1147      } else if (IsFactor(a) && IsBinFactor(b) &&
1148        ((IVariableTreeNode)a).VariableName == ((IVariableTreeNode)b).VariableName) {
1149        var node0 = a as FactorVariableTreeNode;
1150        var node1 = b as BinaryFactorVariableTreeNode;
1151        var varValues = node0.Symbol.GetVariableValues(node0.VariableName).ToArray();
1152        var wi = Array.IndexOf(varValues, node1.VariableValue);
1153        if (wi < 0) throw new ArgumentException();
1154        var newWeighs = new double[varValues.Length];
1155        node0.Weights.CopyTo(newWeighs, 0);
1156        newWeighs[wi] += node1.Weight;
1157        return MakeFactor(node0.Symbol, node0.VariableName, newWeighs);
1158      } else if (IsAddition(a) && IsAddition(b)) {
1159        // merge additions
1160        var add = addSymbol.CreateTreeNode();
1161        // add all sub trees except for the last
1162        for (int i = 0; i < a.Subtrees.Count() - 1; i++) add.AddSubtree(a.GetSubtree(i));
1163        for (int i = 0; i < b.Subtrees.Count() - 1; i++) add.AddSubtree(b.GetSubtree(i));
1164        if (IsConstant(a.Subtrees.Last()) && IsConstant(b.Subtrees.Last())) {
1165          add.AddSubtree(MakeSum(a.Subtrees.Last(), b.Subtrees.Last()));
1166        } else if (IsConstant(a.Subtrees.Last())) {
1167          add.AddSubtree(b.Subtrees.Last());
1168          add.AddSubtree(a.Subtrees.Last());
1169        } else {
1170          add.AddSubtree(a.Subtrees.Last());
1171          add.AddSubtree(b.Subtrees.Last());
1172        }
1173        MergeVariablesInSum(add);
1174        if (add.Subtrees.Count() == 1) {
1175          return add.GetSubtree(0);
1176        } else {
1177          return add;
1178        }
1179      } else if (IsAddition(b)) {
1180        return MakeSum(b, a);
1181      } else if (IsAddition(a) && IsConstant(b)) {
1182        // a is an addition and b is a constant => append b to a and make sure the constants are merged
1183        var add = addSymbol.CreateTreeNode();
1184        // add all sub trees except for the last
1185        for (int i = 0; i < a.Subtrees.Count() - 1; i++) add.AddSubtree(a.GetSubtree(i));
1186        if (IsConstant(a.Subtrees.Last()))
1187          add.AddSubtree(MakeSum(a.Subtrees.Last(), b));
1188        else {
1189          add.AddSubtree(a.Subtrees.Last());
1190          add.AddSubtree(b);
1191        }
1192        return add;
1193      } else if (IsAddition(a)) {
1194        // a is already an addition => append b
1195        var add = addSymbol.CreateTreeNode();
1196        add.AddSubtree(b);
1197        foreach (var subtree in a.Subtrees) {
1198          add.AddSubtree(subtree);
1199        }
1200        MergeVariablesInSum(add);
1201        if (add.Subtrees.Count() == 1) {
1202          return add.GetSubtree(0);
1203        } else {
1204          return add;
1205        }
1206      } else {
1207        var add = addSymbol.CreateTreeNode();
1208        add.AddSubtree(a);
1209        add.AddSubtree(b);
1210        MergeVariablesInSum(add);
1211        if (add.Subtrees.Count() == 1) {
1212          return add.GetSubtree(0);
1213        } else {
1214          return add;
1215        }
1216      }
1217    }
1218
1219    // makes sure variable symbols in sums are combined
1220    private void MergeVariablesInSum(ISymbolicExpressionTreeNode sum) {
1221      var subtrees = new List<ISymbolicExpressionTreeNode>(sum.Subtrees);
1222      while (sum.Subtrees.Any()) sum.RemoveSubtree(0);
1223      var groupedVarNodes = from node in subtrees.OfType<IVariableTreeNode>()
1224                            where node.SubtreeCount == 0
1225                            group node by GroupId(node) into g
1226                            select g;
1227      var constant = (from node in subtrees.OfType<ConstantTreeNode>()
1228                      select node.Value).DefaultIfEmpty(0.0).Sum();
1229      var unchangedSubtrees = subtrees.Where(t => t.SubtreeCount > 0 || !(t is IVariableTreeNode) && !(t is ConstantTreeNode));
1230
1231      foreach (var variableNodeGroup in groupedVarNodes) {
1232        var firstNode = variableNodeGroup.First();
1233        if (firstNode is VariableTreeNodeBase) {
1234          var representative = firstNode as VariableTreeNodeBase;
1235          var weightSum = variableNodeGroup.Cast<VariableTreeNodeBase>().Select(t => t.Weight).Sum();
1236          representative.Weight = weightSum;
1237          sum.AddSubtree(representative);
1238        } else if (firstNode is FactorVariableTreeNode) {
1239          var representative = firstNode as FactorVariableTreeNode;
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          sum.AddSubtree(representative);
1246        }
1247      }
1248      foreach (var unchangedSubtree in unchangedSubtrees)
1249        sum.AddSubtree(unchangedSubtree);
1250      if (!constant.IsAlmost(0.0)) {
1251        sum.AddSubtree(MakeConstant(constant));
1252      }
1253    }
1254
1255    // nodes referencing variables can be grouped if they have
1256    private string GroupId(IVariableTreeNode node) {
1257      var binaryFactorNode = node as BinaryFactorVariableTreeNode;
1258      var factorNode = node as FactorVariableTreeNode;
1259      var variableNode = node as VariableTreeNode;
1260      var laggedVarNode = node as LaggedVariableTreeNode;
1261      if (variableNode != null) {
1262        return "var " + variableNode.VariableName;
1263      } else if (binaryFactorNode != null) {
1264        return "binfactor " + binaryFactorNode.VariableName + " " + binaryFactorNode.VariableValue;
1265      } else if (factorNode != null) {
1266        return "factor " + factorNode.VariableName;
1267      } else if (laggedVarNode != null) {
1268        return "lagged " + laggedVarNode.VariableName + " " + laggedVarNode.Lag;
1269      } else {
1270        throw new NotSupportedException();
1271      }
1272    }
1273
1274
1275    private ISymbolicExpressionTreeNode MakeProduct(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
1276      if (IsConstant(a) && IsConstant(b)) {
1277        // fold constants
1278        return MakeConstant(((ConstantTreeNode)a).Value * ((ConstantTreeNode)b).Value);
1279      } else if (IsConstant(a)) {
1280        // a * $ => $ * a
1281        return MakeProduct(b, a);
1282      } else if (IsFactor(a) && IsFactor(b) && AreSameTypeAndVariable(a, b)) {
1283        var node0 = a as FactorVariableTreeNode;
1284        var node1 = b as FactorVariableTreeNode;
1285        return MakeFactor(node0.Symbol, node0.VariableName, node0.Weights.Zip(node1.Weights, (u, v) => u * v));
1286      } else if (IsBinFactor(a) && IsBinFactor(b) && AreSameTypeAndVariable(a, b)) {
1287        var node0 = a as BinaryFactorVariableTreeNode;
1288        var node1 = b as BinaryFactorVariableTreeNode;
1289        return MakeBinFactor(node0.Symbol, node0.VariableName, node0.VariableValue, node0.Weight * node1.Weight);
1290      } else if (IsFactor(a) && IsConstant(b)) {
1291        var node0 = a as FactorVariableTreeNode;
1292        var node1 = b as ConstantTreeNode;
1293        return MakeFactor(node0.Symbol, node0.VariableName, node0.Weights.Select(w => w * node1.Value));
1294      } else if (IsBinFactor(a) && IsConstant(b)) {
1295        var node0 = a as BinaryFactorVariableTreeNode;
1296        var node1 = b as ConstantTreeNode;
1297        return MakeBinFactor(node0.Symbol, node0.VariableName, node0.VariableValue, node0.Weight * node1.Value);
1298      } else if (IsBinFactor(a) && IsFactor(b)) {
1299        return MakeProduct(b, a);
1300      } else if (IsFactor(a) && IsBinFactor(b) &&
1301        ((IVariableTreeNode)a).VariableName == ((IVariableTreeNode)b).VariableName) {
1302        var node0 = a as FactorVariableTreeNode;
1303        var node1 = b as BinaryFactorVariableTreeNode;
1304        var varValues = node0.Symbol.GetVariableValues(node0.VariableName).ToArray();
1305        var wi = Array.IndexOf(varValues, node1.VariableValue);
1306        if (wi < 0) throw new ArgumentException();
1307        return MakeBinFactor(node1.Symbol, node1.VariableName, node1.VariableValue, node1.Weight * node0.Weights[wi]);
1308      } else if (IsConstant(b) && ((ConstantTreeNode)b).Value.IsAlmost(1.0)) {
1309        // $ * 1.0 => $
1310        return a;
1311      } else if (IsConstant(b) && ((ConstantTreeNode)b).Value.IsAlmost(0.0)) {
1312        return MakeConstant(0);
1313      } else if (IsConstant(b) && IsVariableBase(a)) {
1314        // multiply constants into variables weights
1315        ((VariableTreeNodeBase)a).Weight *= ((ConstantTreeNode)b).Value;
1316        return a;
1317      } else if (IsConstant(b) && IsAddition(a) ||
1318          IsFactor(b) && IsAddition(a) ||
1319          IsBinFactor(b) && IsAddition(a)) {
1320        // multiply constants into additions
1321        return a.Subtrees.Select(x => MakeProduct(GetSimplifiedTree(x), GetSimplifiedTree(b))).Aggregate((c, d) => MakeSum(c, d));
1322      } else if (IsDivision(a) && IsDivision(b)) {
1323        // (a1 / a2) * (b1 / b2) => (a1 * b1) / (a2 * b2)
1324        return MakeFraction(MakeProduct(a.GetSubtree(0), b.GetSubtree(0)), MakeProduct(a.GetSubtree(1), b.GetSubtree(1)));
1325      } else if (IsDivision(a)) {
1326        // (a1 / a2) * b => (a1 * b) / a2
1327        return MakeFraction(MakeProduct(a.GetSubtree(0), b), a.GetSubtree(1));
1328      } else if (IsDivision(b)) {
1329        // a * (b1 / b2) => (b1 * a) / b2
1330        return MakeFraction(MakeProduct(b.GetSubtree(0), a), b.GetSubtree(1));
1331      } else if (IsMultiplication(a) && IsMultiplication(b)) {
1332        // merge multiplications (make sure constants are merged)
1333        var mul = mulSymbol.CreateTreeNode();
1334        for (int i = 0; i < a.Subtrees.Count(); i++) mul.AddSubtree(a.GetSubtree(i));
1335        for (int i = 0; i < b.Subtrees.Count(); i++) mul.AddSubtree(b.GetSubtree(i));
1336        MergeVariablesAndConstantsInProduct(mul);
1337        return mul;
1338      } else if (IsMultiplication(b)) {
1339        return MakeProduct(b, a);
1340      } else if (IsMultiplication(a)) {
1341        // a is already an multiplication => append b
1342        a.AddSubtree(GetSimplifiedTree(b));
1343        MergeVariablesAndConstantsInProduct(a);
1344        return a;
1345      } else if (IsAbsolute(a) && IsAbsolute(b)) {
1346        return MakeAbs(MakeProduct(a.GetSubtree(0), b.GetSubtree(0)));
1347      } else if (IsAbsolute(a) && IsConstant(b)) {
1348        var constNode = b as ConstantTreeNode;
1349        var posF = Math.Abs(constNode.Value);
1350        if (constNode.Value > 0) {
1351          return MakeAbs(MakeProduct(a.GetSubtree(0), MakeConstant(posF)));
1352        } else {
1353          var mul = mulSymbol.CreateTreeNode();
1354          mul.AddSubtree(MakeAbs(MakeProduct(a.GetSubtree(0), MakeConstant(posF))));
1355          mul.AddSubtree(MakeConstant(-1.0));
1356          return mul;
1357        }
1358      } else if (IsAnalyticalQuotient(a)) {
1359        return MakeAnalyticalQuotient(MakeProduct(a.GetSubtree(0), b), a.GetSubtree(1));
1360      } else {
1361        var mul = mulSymbol.CreateTreeNode();
1362        mul.AddSubtree(a);
1363        mul.AddSubtree(b);
1364        MergeVariablesAndConstantsInProduct(mul);
1365        return mul;
1366      }
1367    }
1368
1369    ISymbolicExpressionTreeNode MakeSumAggregation(ISymbolicExpressionTreeNode node) {
1370      if (IsConstant(node)) { // assumes scalar constant
1371        return node;
1372      } else if (IsAddition(node) || IsSubtraction(node)) {
1373        var terms = node.Subtrees;
1374        if (IsSubtraction(node)) {
1375          if (terms.Count() == 1) {
1376            terms = terms.Select(Negate);
1377          } else {
1378            var first = terms.First().ToEnumerable();
1379            var remaining = terms.Skip(1).Select(Negate);
1380            terms = first.Concat(remaining).ToList();
1381          }
1382        }
1383
1384        var scalarTerms = terms.Where(IsScalarNode).ToList();
1385        var remainingTerms = terms.Except(scalarTerms).ToList();
1386
1387        if (scalarTerms.Any() && remainingTerms.Any()) {
1388          var scalarNode = scalarTerms.Aggregate(MakeSum);
1389          var vectorNode = remainingTerms.Aggregate(MakeSum);
1390
1391          var lengthNode = lengthSymbol.CreateTreeNode();
1392          lengthNode.AddSubtree((ISymbolicExpressionTreeNode)vectorNode.Clone());
1393          var scalarMulNode = mulSymbol.CreateTreeNode();
1394          scalarMulNode.AddSubtree(scalarNode);
1395          scalarMulNode.AddSubtree(lengthNode);
1396
1397          var sumNode = sumSymbol.CreateTreeNode();
1398          sumNode.AddSubtree(vectorNode);
1399
1400          return MakeSum(scalarMulNode, sumNode);
1401        } else if (scalarTerms.Any()) {
1402          var scalarNode = scalarTerms.Aggregate(MakeSum);
1403          return scalarNode;
1404        } else if (remainingTerms.Any()) {
1405          var vectorNode = remainingTerms.Aggregate(MakeSum);
1406          var sumNode = sumSymbol.CreateTreeNode();
1407          sumNode.AddSubtree(vectorNode);
1408          return sumNode;
1409        } else
1410          throw new InvalidOperationException("Addition does not contain any terms to simplify.");
1411      } else if (IsMultiplication(node) || IsDivision(node)) {
1412        var factors = node.Subtrees;
1413        if (IsDivision(node)) {
1414          if (factors.Count() == 1) {
1415            factors = factors.Select(Invert);
1416          } else {
1417            var first = factors.First().ToEnumerable();
1418            var remaining = factors.Skip(1).Select(Invert);
1419            factors = first.Concat(remaining).ToList();
1420          }
1421        }
1422
1423        var scalarFactors = factors.Where(IsScalarNode).ToList();
1424        var remainingFactors = factors.Except(scalarFactors).ToList();
1425
1426        if (scalarFactors.Any() && remainingFactors.Any()) {
1427          var scalarNode = scalarFactors.Aggregate(MakeProduct);
1428          var vectorNode = remainingFactors.Aggregate(MakeProduct);
1429
1430          var sumNode = sumSymbol.CreateTreeNode();
1431          sumNode.AddSubtree(vectorNode);
1432
1433          return MakeProduct(scalarNode, sumNode);
1434        } else if (scalarFactors.Any()) {
1435          var scalarNode = scalarFactors.Aggregate(MakeProduct);
1436          return scalarNode;
1437        } else if (remainingFactors.Any()) {
1438          var vectorNode = remainingFactors.Aggregate(MakeProduct);
1439          var sumNode = sumSymbol.CreateTreeNode();
1440          sumNode.AddSubtree(vectorNode);
1441          return sumNode;
1442        } else
1443          throw new InvalidOperationException("Multiplication does not contain any terms to simplify.");
1444      } else {
1445        var sumNode = sumSymbol.CreateTreeNode();
1446        sumNode.AddSubtree(node);
1447        return sumNode;
1448      }
1449    }
1450
1451    ISymbolicExpressionTreeNode MakeMeanAggregation(ISymbolicExpressionTreeNode node) {
1452      if (IsConstant(node)) { // assumes scalar constant
1453        return node;
1454      } else if (IsAddition(node) || IsSubtraction(node)) {
1455        var terms = node.Subtrees;
1456        if (IsSubtraction(node)) {
1457          if (terms.Count() == 1) {
1458            terms = terms.Select(Negate);
1459          } else {
1460            var first = terms.First().ToEnumerable();
1461            var remaining = terms.Skip(1).Select(Negate);
1462            terms = first.Concat(remaining).ToList();
1463          }
1464        }
1465
1466        var scalarTerms = terms.Where(IsScalarNode).ToList();
1467        var remainingTerms = terms.Except(scalarTerms).ToList();
1468
1469        if (scalarTerms.Any() && remainingTerms.Any()) {
1470          var scalarNode = scalarTerms.Aggregate(MakeSum);
1471          var vectorNode = remainingTerms.Aggregate(MakeSum);
1472
1473          var meanNode = meanSymbol.CreateTreeNode();
1474          meanNode.AddSubtree(vectorNode);
1475
1476          return MakeSum(scalarNode, meanNode);
1477        } else if (scalarTerms.Any()) {
1478          var scalarNode = scalarTerms.Aggregate(MakeSum);
1479          return scalarNode;
1480        } else if (remainingTerms.Any()) {
1481          var vectorNode = remainingTerms.Aggregate(MakeSum);
1482          var meanNode = meanSymbol.CreateTreeNode();
1483          meanNode.AddSubtree(vectorNode);
1484          return meanNode;
1485        } else
1486          throw new InvalidOperationException("Addition does not contain any terms to simplify.");
1487      } else if (IsMultiplication(node) || IsDivision(node)) {
1488        var factors = node.Subtrees;
1489        if (IsDivision(node)) {
1490          if (factors.Count() == 1) {
1491            factors = factors.Select(Invert);
1492          } else {
1493            var first = factors.First().ToEnumerable();
1494            var remaining = factors.Skip(1).Select(Invert);
1495            factors = first.Concat(remaining).ToList();
1496          }
1497        }
1498
1499        var scalarFactors = factors.Where(IsScalarNode).ToList();
1500        var remainingFactors = factors.Except(scalarFactors).ToList();
1501
1502        if (scalarFactors.Any() && remainingFactors.Any()) {
1503          var scalarNode = scalarFactors.Aggregate(MakeProduct);
1504          var vectorNode = remainingFactors.Aggregate(MakeProduct);
1505
1506          var meanNode = meanSymbol.CreateTreeNode();
1507          meanNode.AddSubtree(vectorNode);
1508
1509          return MakeProduct(scalarNode, meanNode);
1510        } else if (scalarFactors.Any()) {
1511          var scalarNode = scalarFactors.Aggregate(MakeProduct);
1512          return scalarNode;
1513        } else if (remainingFactors.Any()) {
1514          var vectorNode = remainingFactors.Aggregate(MakeProduct);
1515          var meanNode = meanSymbol.CreateTreeNode();
1516          meanNode.AddSubtree(vectorNode);
1517          return meanNode;
1518        } else
1519          throw new InvalidOperationException("Multiplication does not contain any terms to simplify.");
1520      } else {
1521        var meanNode = meanSymbol.CreateTreeNode();
1522        meanNode.AddSubtree(node);
1523        return meanNode;
1524      }
1525    }
1526
1527    #endregion
1528
1529    #region helper functions
1530
1531    private bool ContainsVariableCondition(ISymbolicExpressionTreeNode node) {
1532      if (node.Symbol is VariableCondition) return true;
1533      foreach (var subtree in node.Subtrees)
1534        if (ContainsVariableCondition(subtree)) return true;
1535      return false;
1536    }
1537
1538    private ISymbolicExpressionTreeNode AddLagToDynamicNodes(ISymbolicExpressionTreeNode node, int lag) {
1539      var laggedTreeNode = node as ILaggedTreeNode;
1540      var variableNode = node as VariableTreeNode;
1541      var variableConditionNode = node as VariableConditionTreeNode;
1542      if (laggedTreeNode != null)
1543        laggedTreeNode.Lag += lag;
1544      else if (variableNode != null) {
1545        var laggedVariableNode = (LaggedVariableTreeNode)laggedVariableSymbol.CreateTreeNode();
1546        laggedVariableNode.Lag = lag;
1547        laggedVariableNode.VariableName = variableNode.VariableName;
1548        return laggedVariableNode;
1549      } else if (variableConditionNode != null) {
1550        throw new NotSupportedException("Removal of time lags around variable condition symbols is not allowed.");
1551      }
1552      var subtrees = new List<ISymbolicExpressionTreeNode>(node.Subtrees);
1553      while (node.SubtreeCount > 0) node.RemoveSubtree(0);
1554      foreach (var subtree in subtrees) {
1555        node.AddSubtree(AddLagToDynamicNodes(subtree, lag));
1556      }
1557      return node;
1558    }
1559
1560    private bool AreSameTypeAndVariable(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b) {
1561      return GroupId((IVariableTreeNode)a) == GroupId((IVariableTreeNode)b);
1562    }
1563
1564    // helper to combine the constant factors in products and to combine variables (powers of 2, 3...)
1565    private void MergeVariablesAndConstantsInProduct(ISymbolicExpressionTreeNode prod) {
1566      var subtrees = new List<ISymbolicExpressionTreeNode>(prod.Subtrees);
1567      while (prod.Subtrees.Any()) prod.RemoveSubtree(0);
1568      var groupedVarNodes = from node in subtrees.OfType<IVariableTreeNode>()
1569                            where node.SubtreeCount == 0
1570                            group node by GroupId(node) into g
1571                            orderby g.Count()
1572                            select g;
1573      var constantProduct = (from node in subtrees.OfType<VariableTreeNodeBase>()
1574                             select node.Weight)
1575        .Concat(from node in subtrees.OfType<ConstantTreeNode>()
1576                select node.Value)
1577        .DefaultIfEmpty(1.0)
1578        .Aggregate((c1, c2) => c1 * c2);
1579
1580      var unchangedSubtrees = from tree in subtrees
1581                              where tree.SubtreeCount > 0 || !(tree is IVariableTreeNode) && !(tree is ConstantTreeNode)
1582                              select tree;
1583
1584      foreach (var variableNodeGroup in groupedVarNodes) {
1585        var firstNode = variableNodeGroup.First();
1586        if (firstNode is VariableTreeNodeBase) {
1587          var representative = (VariableTreeNodeBase)firstNode;
1588          representative.Weight = 1.0;
1589          if (variableNodeGroup.Count() > 1) {
1590            var poly = mulSymbol.CreateTreeNode();
1591            for (int p = 0; p < variableNodeGroup.Count(); p++) {
1592              poly.AddSubtree((ISymbolicExpressionTreeNode)representative.Clone());
1593            }
1594            prod.AddSubtree(poly);
1595          } else {
1596            prod.AddSubtree(representative);
1597          }
1598        } else if (firstNode is FactorVariableTreeNode) {
1599          var representative = (FactorVariableTreeNode)firstNode;
1600          foreach (var node in variableNodeGroup.Skip(1).Cast<FactorVariableTreeNode>()) {
1601            for (int j = 0; j < representative.Weights.Length; j++) {
1602              representative.Weights[j] *= node.Weights[j];
1603            }
1604          }
1605          for (int j = 0; j < representative.Weights.Length; j++) {
1606            representative.Weights[j] *= constantProduct;
1607          }
1608          constantProduct = 1.0;
1609          // if the product already contains a factor it is not necessary to multiply a constant below
1610          prod.AddSubtree(representative);
1611        }
1612      }
1613
1614      foreach (var unchangedSubtree in unchangedSubtrees)
1615        prod.AddSubtree(unchangedSubtree);
1616
1617      if (!constantProduct.IsAlmost(1.0)) {
1618        prod.AddSubtree(MakeConstant(constantProduct));
1619      }
1620    }
1621
1622
1623    /// <summary>
1624    /// x => x * -1
1625    /// Is only used in cases where it is not necessary to create new tree nodes. Manipulates x directly.
1626    /// </summary>
1627    /// <param name="x"></param>
1628    /// <returns>-x</returns>
1629    private ISymbolicExpressionTreeNode Negate(ISymbolicExpressionTreeNode x) {
1630      if (IsConstant(x)) {
1631        ((ConstantTreeNode)x).Value *= -1;
1632      } else if (IsVariableBase(x)) {
1633        var variableTree = (VariableTreeNodeBase)x;
1634        variableTree.Weight *= -1.0;
1635      } else if (IsFactor(x)) {
1636        var factorNode = (FactorVariableTreeNode)x;
1637        for (int i = 0; i < factorNode.Weights.Length; i++) factorNode.Weights[i] *= -1;
1638      } else if (IsBinFactor(x)) {
1639        var factorNode = (BinaryFactorVariableTreeNode)x;
1640        factorNode.Weight *= -1;
1641      } else if (IsAddition(x)) {
1642        // (x0 + x1 + .. + xn) * -1 => (-x0 + -x1 + .. + -xn)       
1643        var subtrees = new List<ISymbolicExpressionTreeNode>(x.Subtrees);
1644        while (x.Subtrees.Any()) x.RemoveSubtree(0);
1645        foreach (var subtree in subtrees) {
1646          x.AddSubtree(Negate(subtree));
1647        }
1648      } else if (IsMultiplication(x) || IsDivision(x)) {
1649        // x0 * x1 * .. * xn * -1 => x0 * x1 * .. * -xn
1650        var lastSubTree = x.Subtrees.Last();
1651        x.RemoveSubtree(x.SubtreeCount - 1);
1652        x.AddSubtree(Negate(lastSubTree)); // last is maybe a constant, prefer to negate the constant
1653      } else {
1654        // any other function
1655        return MakeProduct(x, MakeConstant(-1));
1656      }
1657      return x;
1658    }
1659
1660    /// <summary>
1661    /// x => 1/x
1662    /// Must create new tree nodes
1663    /// </summary>
1664    /// <param name="x"></param>
1665    /// <returns></returns>
1666    private ISymbolicExpressionTreeNode Invert(ISymbolicExpressionTreeNode x) {
1667      if (IsConstant(x)) {
1668        return MakeConstant(1.0 / ((ConstantTreeNode)x).Value);
1669      } else if (IsFactor(x)) {
1670        var factorNode = (FactorVariableTreeNode)x;
1671        return MakeFactor(factorNode.Symbol, factorNode.VariableName, factorNode.Weights.Select(w => 1.0 / w));
1672      } else if (IsDivision(x)) {
1673        return MakeFraction(x.GetSubtree(1), x.GetSubtree(0));
1674      } else {
1675        // any other function
1676        return MakeFraction(MakeConstant(1), x);
1677      }
1678    }
1679
1680    private ISymbolicExpressionTreeNode MakeConstant(double value) {
1681      ConstantTreeNode constantTreeNode = (ConstantTreeNode)(constSymbol.CreateTreeNode());
1682      constantTreeNode.Value = value;
1683      return constantTreeNode;
1684    }
1685
1686    private ISymbolicExpressionTreeNode MakeFactor(FactorVariable sy, string variableName, IEnumerable<double> weights) {
1687      var tree = (FactorVariableTreeNode)sy.CreateTreeNode();
1688      tree.VariableName = variableName;
1689      tree.Weights = weights.ToArray();
1690      return tree;
1691    }
1692    private ISymbolicExpressionTreeNode MakeBinFactor(BinaryFactorVariable sy, string variableName, string variableValue, double weight) {
1693      var tree = (BinaryFactorVariableTreeNode)sy.CreateTreeNode();
1694      tree.VariableName = variableName;
1695      tree.VariableValue = variableValue;
1696      tree.Weight = weight;
1697      return tree;
1698    }
1699
1700
1701    #endregion
1702  }
1703}
Note: See TracBrowser for help on using the repository browser.