Index: /branches/3073_IA_constraint_splitting/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/HeuristicLab.Problems.DataAnalysis.Symbolic-3.4.csproj
===================================================================
--- /branches/3073_IA_constraint_splitting/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/HeuristicLab.Problems.DataAnalysis.Symbolic-3.4.csproj	(revision 17771)
+++ /branches/3073_IA_constraint_splitting/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/HeuristicLab.Problems.DataAnalysis.Symbolic-3.4.csproj	(revision 17772)
@@ -241,4 +241,5 @@
     <Compile Include="Interpreter\SymbolicDataAnalysisExpressionCompiledTreeInterpreter.cs" />
     <Compile Include="Interpreter\SymbolicDataAnalysisExpressionTreeBatchInterpreter.cs" />
+    <Compile Include="Interpreter\IACompiledExpressionBoundsEstimator.cs" />
     <Compile Include="Interpreter\SymbolicDataAnalysisExpressionTreeNativeInterpreter.cs" />
     <Compile Include="IntervalUtil.cs" />
Index: /branches/3073_IA_constraint_splitting/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/IACompiledExpressionBoundsEstimator.cs
===================================================================
--- /branches/3073_IA_constraint_splitting/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/IACompiledExpressionBoundsEstimator.cs	(revision 17772)
+++ /branches/3073_IA_constraint_splitting/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/IACompiledExpressionBoundsEstimator.cs	(revision 17772)
@@ -0,0 +1,400 @@
+﻿using HEAL.Attic;
+
+using HeuristicLab.Common;
+using HeuristicLab.Core;
+using HeuristicLab.Data;
+using HeuristicLab.Encodings.SymbolicExpressionTreeEncoding;
+using HeuristicLab.Parameters;
+
+using System;
+using System.Collections.Generic;
+using System.Linq;
+using System.Linq.Expressions;
+
+namespace HeuristicLab.Problems.DataAnalysis.Symbolic
+{
+  [StorableType("60015D64-5D8B-408A-90A1-E4111BC114D4")]
+  [Item("IA Compiled Expression Bounds Estimator", "Compile a symbolic model into a lambda and use it to evaluate model bounds.")]
+  public class IACompiledExpressionBoundsEstimator : ParameterizedNamedItem, IBoundsEstimator
+  {
+    // interval method names
+    private static readonly Dictionary<byte, string> methodName = new Dictionary<byte, string>() {
+      { OpCodes.Add, "Add" },
+      { OpCodes.Sub, "Subtract" },
+      { OpCodes.Mul, "Multiply" },
+      { OpCodes.Div, "Divide" },
+      { OpCodes.Sin, "Sine" },
+      { OpCodes.Cos, "Cosine" },
+      { OpCodes.Tan, "Tangens" },
+      { OpCodes.Tanh, "HyperbolicTangent" },
+      { OpCodes.Log, "Logarithm" },
+      { OpCodes.Exp, "Exponential" },
+      { OpCodes.Square, "Square" },
+      { OpCodes.Cube, "Cube" },
+      { OpCodes.SquareRoot, "SquareRoot" },
+      { OpCodes.CubeRoot, "CubicRoot" },
+      { OpCodes.Absolute, "Absolute" },
+      { OpCodes.AnalyticQuotient, "AnalyticalQuotient" },
+    };
+
+    private const string EvaluatedSolutionsParameterName = "EvaluatedSolutions";
+    private const string UseIntervalSplittingParameterName = "Use Interval splitting";
+    private const string MaxSplitParameterName = "MaxSplit";
+    private const string MinWidthParameterName = "MinWidth";
+
+    public IFixedValueParameter<IntValue> EvaluatedSolutionsParameter {
+      get => (IFixedValueParameter<IntValue>)Parameters[EvaluatedSolutionsParameterName];
+    }
+
+    public IFixedValueParameter<BoolValue> UseIntervalSplittingParameter {
+      get => (IFixedValueParameter<BoolValue>)Parameters[UseIntervalSplittingParameterName];
+    }
+
+    public IFixedValueParameter<IntValue> MaxSplitParameter {
+      get => (IFixedValueParameter<IntValue>)Parameters[MaxSplitParameterName];
+    }
+
+    public IFixedValueParameter<DoubleValue> MinWidthParameter {
+      get => (IFixedValueParameter<DoubleValue>)Parameters[MinWidthParameterName];
+    }
+
+    public int MaxSplit {
+      get => MaxSplitParameter.Value.Value;
+      set => MaxSplitParameter.Value.Value = value;
+    }
+
+    public double MinWidth {
+      get => MinWidthParameter.Value.Value;
+      set => MinWidthParameter.Value.Value = value;
+    }
+
+    public int EvaluatedSolutions {
+      get => EvaluatedSolutionsParameter.Value.Value;
+      set => EvaluatedSolutionsParameter.Value.Value = value;
+    }
+
+    public bool UseIntervalSplitting {
+      get => UseIntervalSplittingParameter.Value.Value;
+      set => UseIntervalSplittingParameter.Value.Value = value;
+    }
+
+    private readonly object syncRoot = new object();
+
+    public IACompiledExpressionBoundsEstimator() : base("IA Bounds Estimator",
+      "Estimates the bounds of the model with interval arithmetic, by first compiling the model into a lambda.") {
+      Parameters.Add(new FixedValueParameter<IntValue>(EvaluatedSolutionsParameterName,
+        "A counter for the total number of solutions the estimator has evaluated.", new IntValue(0)));
+      Parameters.Add(new FixedValueParameter<BoolValue>(UseIntervalSplittingParameterName,
+        "Defines whether interval splitting is activated or not.", new BoolValue(false)));
+      Parameters.Add(new FixedValueParameter<IntValue>(MaxSplitParameterName,
+        "Defines the number of iterations of splitting.", new IntValue(200)));
+      Parameters.Add(new FixedValueParameter<DoubleValue>(MinWidthParameterName,
+        "Width of interval, after the splitting should stop.", new DoubleValue(0.0)));
+    }
+
+    [StorableConstructor]
+    private IACompiledExpressionBoundsEstimator(StorableConstructorFlag _) : base(_) { }
+
+    private IACompiledExpressionBoundsEstimator(IACompiledExpressionBoundsEstimator original, Cloner cloner) : base(original, cloner) { }
+
+    public override IDeepCloneable Clone(Cloner cloner) {
+      return new IACompiledExpressionBoundsEstimator(this, cloner);
+    }
+
+
+
+    public double CheckConstraint(ISymbolicExpressionTree tree, IntervalCollection variableRanges, IntervalConstraint constraint) {
+      if (!UseIntervalSplitting) {
+        var modelBound = GetModelBound(tree, variableRanges);
+        if (constraint.Interval.Contains(modelBound)) return 0.0;
+        return Math.Abs(modelBound.LowerBound - constraint.Interval.LowerBound) +
+               Math.Abs(modelBound.UpperBound - constraint.Interval.UpperBound);
+      }
+
+      if (double.IsNegativeInfinity(constraint.Interval.LowerBound) &&
+          double.IsPositiveInfinity(constraint.Interval.UpperBound)) {
+        return 0.0;
+      }
+
+      //ContainsVariableMultipleTimes(tree, out var variables);
+
+      lock (syncRoot) { EvaluatedSolutions++; }
+
+      double upperBound;
+      if (double.IsNegativeInfinity(constraint.Interval.LowerBound)) {
+        upperBound = EstimateUpperBound(tree, variableRanges.GetReadonlyDictionary(), MaxSplit, MinWidth);
+
+        return upperBound <= constraint.Interval.UpperBound
+          ? 0.0
+          : Math.Abs(upperBound - constraint.Interval.UpperBound);
+      }
+
+      double lowerBound;
+      if (double.IsPositiveInfinity(constraint.Interval.UpperBound)) {
+        lowerBound = EstimateLowerBound(tree, variableRanges.GetReadonlyDictionary(), MaxSplit, MinWidth);
+
+        return lowerBound <= constraint.Interval.LowerBound
+          ? 0.0
+          : Math.Abs(lowerBound - constraint.Interval.LowerBound);
+      }
+
+      var ranges = variableRanges.GetReadonlyDictionary();
+      lowerBound = EstimateLowerBound(tree, ranges, MaxSplit, MinWidth);
+      upperBound = EstimateUpperBound(tree, ranges, MaxSplit, MinWidth);
+
+      var res = 0.0;
+
+      res += upperBound <= constraint.Interval.UpperBound ? 0.0 : Math.Abs(upperBound - constraint.Interval.UpperBound);
+      res += lowerBound <= constraint.Interval.LowerBound ? 0.0 : Math.Abs(lowerBound - constraint.Interval.LowerBound);
+
+      return res;
+    }
+
+    public void ClearState() {
+      EvaluatedSolutions = 0;
+    }
+
+    public Interval GetModelBound(ISymbolicExpressionTree tree, IntervalCollection variableRanges) {
+      lock (syncRoot) { EvaluatedSolutions++; }
+      var resultInterval = UseIntervalSplitting
+        ? EstimateBounds(tree, variableRanges.GetReadonlyDictionary(), MaxSplit, MinWidth)
+        : EstimateBounds(tree, variableRanges.GetReadonlyDictionary());
+
+      if (resultInterval.IsInfiniteOrUndefined || resultInterval.LowerBound <= resultInterval.UpperBound)
+        return resultInterval;
+      return new Interval(resultInterval.UpperBound, resultInterval.LowerBound);
+    }
+
+    public IDictionary<ISymbolicExpressionTreeNode, Interval> GetModelNodesBounds(ISymbolicExpressionTree tree, IntervalCollection variableRanges) {
+      throw new NotSupportedException("Model nodes bounds are not supported.");
+    }
+
+    public void InitializeState() {
+      EvaluatedSolutions = 0;
+    }
+
+    public bool IsCompatible(ISymbolicExpressionTree tree) {
+      var containsUnknownSymbols = (
+        from n in tree.Root.GetSubtree(0).IterateNodesPrefix()
+        where
+          !(n.Symbol is Variable) &&
+          !(n.Symbol is Constant) &&
+          !(n.Symbol is StartSymbol) &&
+          !(n.Symbol is Addition) &&
+          !(n.Symbol is Subtraction) &&
+          !(n.Symbol is Multiplication) &&
+          !(n.Symbol is Division) &&
+          !(n.Symbol is Sine) &&
+          !(n.Symbol is Cosine) &&
+          !(n.Symbol is Tangent) &&
+          !(n.Symbol is HyperbolicTangent) &&
+          !(n.Symbol is Logarithm) &&
+          !(n.Symbol is Exponential) &&
+          !(n.Symbol is Square) &&
+          !(n.Symbol is SquareRoot) &&
+          !(n.Symbol is Cube) &&
+          !(n.Symbol is CubeRoot) &&
+          !(n.Symbol is Absolute) &&
+          !(n.Symbol is AnalyticQuotient)
+        select n).Any();
+      return !containsUnknownSymbols;
+    }
+
+    #region compile a tree into a IA arithmetic lambda and estimate bounds
+    static Expression MakeExpr(ISymbolicExpressionTreeNode node, IReadOnlyDictionary<string, Interval> variableRanges, IReadOnlyDictionary<string, int> variableIndices, Expression args) {
+      Expression expr(ISymbolicExpressionTreeNode n) => MakeExpr(n, variableRanges, variableIndices, args);
+      var opCode = OpCodes.MapSymbolToOpCode(node);
+
+      switch (opCode) {
+        case OpCodes.Variable: {
+            var name = (node as VariableTreeNode).VariableName;
+            var weight = (node as VariableTreeNode).Weight;
+            var index = variableIndices[name];
+            return Expression.Multiply(
+              Expression.Constant(weight, typeof(double)),
+              Expression.ArrayIndex(args, Expression.Constant(index, typeof(int)))
+            );
+          }
+        case OpCodes.Constant: {
+            var v = (node as ConstantTreeNode).Value;
+            // we have to make an interval out of the constant because this may be the root of the tree (and we are expected to return an Interval)
+            return Expression.Constant(new Interval(v, v), typeof(Interval));
+          }
+        case OpCodes.Add: {
+            var e = expr(node.GetSubtree(0));
+            foreach (var s in node.Subtrees.Skip(1)) {
+              e = Expression.Add(e, expr(s));
+            }
+            return e;
+          }
+        case OpCodes.Sub: {
+            var e = expr(node.GetSubtree(0));
+            if (node.SubtreeCount == 1) {
+              return Expression.Subtract(Expression.Constant(0.0, typeof(double)), e);
+            }
+            foreach (var s in node.Subtrees.Skip(1)) {
+              e = Expression.Subtract(e, expr(s));
+            }
+            return e;
+          }
+        case OpCodes.Mul: {
+            var e = expr(node.GetSubtree(0));
+            foreach (var s in node.Subtrees.Skip(1)) {
+              e = Expression.Multiply(e, expr(s));
+            }
+            return e;
+          }
+        case OpCodes.Div: {
+            var e1 = expr(node.GetSubtree(0));
+            if (node.SubtreeCount == 1) {
+              return Expression.Divide(Expression.Constant(1.0, typeof(double)), e1);
+            }
+            // division is more expensive than multiplication so we use this construct
+            var e2 = expr(node.GetSubtree(1));
+            foreach (var s in node.Subtrees.Skip(2)) {
+              e2 = Expression.Multiply(e2, expr(s));
+            }
+            return Expression.Divide(e1, e2);
+          }
+        // all these cases share the same code: get method info by name, emit call expression
+        case OpCodes.Exp:
+        case OpCodes.Log:
+        case OpCodes.Sin:
+        case OpCodes.Cos:
+        case OpCodes.Tan:
+        case OpCodes.Tanh:
+        case OpCodes.Square:
+        case OpCodes.Cube:
+        case OpCodes.SquareRoot:
+        case OpCodes.CubeRoot:
+        case OpCodes.Absolute:
+        case OpCodes.AnalyticQuotient: {
+            var arg = expr(node.GetSubtree(0));
+            var fun = typeof(Interval).GetMethod(methodName[opCode], new[] { arg.Type });
+            return Expression.Call(fun, arg);
+          }
+        default: {
+            throw new Exception($"Unsupported OpCode {opCode} encountered.");
+          }
+      }
+    }
+
+    public static IReadOnlyDictionary<string, int> GetVariableIndices(ISymbolicExpressionTree tree, IReadOnlyDictionary<string, Interval> variableIntervals, out Interval[] inputIntervals) {
+      var variableIndices = new Dictionary<string, int>();
+      var root = tree.Root;
+      while (root.Symbol is ProgramRootSymbol || root.Symbol is StartSymbol) {
+        root = root.GetSubtree(0);
+      }
+      inputIntervals = new Interval[variableIntervals.Count];
+      int count = 0;
+      foreach (var node in root.IterateNodesPrefix()) {
+        if (node is VariableTreeNode varNode) {
+          var name = varNode.VariableName;
+          if (!variableIndices.ContainsKey(name)) {
+            variableIndices[name] = count;
+            inputIntervals[count] = variableIntervals[name];
+            ++count;
+          }
+        }
+      }
+      Array.Resize<Interval>(ref inputIntervals, count);
+      return variableIndices;
+    }
+
+    public static Func<Interval[], Interval> Compile(ISymbolicExpressionTree tree, IReadOnlyDictionary<string, Interval> variableRanges, IReadOnlyDictionary<string, int> variableIndices) {
+      var root = tree.Root.GetSubtree(0).GetSubtree(0);
+      var args = Expression.Parameter(typeof(Interval[]));
+      var expr = MakeExpr(root, variableRanges, variableIndices, args);
+      return Expression.Lambda<Func<Interval[], Interval>>(expr, args).Compile();
+    }
+
+    public static Interval EstimateBounds(ISymbolicExpressionTree tree, IReadOnlyDictionary<string, Interval> variableRanges, int n = 0, double w = 1e-5) {
+      var variableIndices = GetVariableIndices(tree, variableRanges, out Interval[] x);
+      var f = Compile(tree, variableRanges, variableIndices);
+      if (n == 0) return f(x);
+      var inf = EstimateBound(x, f, true, n, w);
+      var sup = EstimateBound(x, f, false, n, w);
+      return inf < sup ? new Interval(inf, sup) : new Interval(sup, inf);
+    }
+    public  double EstimateLowerBound(ISymbolicExpressionTree tree, IReadOnlyDictionary<string, Interval> variableRanges, int n = 1000, double w = 1e-5) {
+      var variableIndices = GetVariableIndices(tree, variableRanges, out Interval[] x);
+      var f = Compile(tree, variableRanges, variableIndices);
+      var inf = EstimateBound(x, f, true, n, w);
+      return inf;
+    }
+
+    public double EstimateUpperBound(ISymbolicExpressionTree tree, IReadOnlyDictionary<string, Interval> variableRanges, int n = 1000, double w = 1e-5) {
+      var variableIndices = GetVariableIndices(tree, variableRanges, out Interval[] x);
+      var f = Compile(tree, variableRanges, variableIndices);
+      var sup = EstimateBound(x, f, false, n, w);
+      return sup;
+    }
+
+    static double EstimateBound(Interval[] x, Func<Interval[], Interval> f, bool m = false, int n = 1000, double w = 1e-4) {
+      double getBound(Interval iv) => m ? iv.LowerBound : -iv.UpperBound;
+      double getVolume(Interval[] box) => box.Aggregate(1.0, (acc, iv) => acc * iv.Width);
+
+      var splits = Enumerable.Range(0, x.Length).Select(_ => new List<Interval>()).ToArray();
+      var newbox = new Interval[x.Length];
+
+      int compare(Tuple<double, double, Interval[]> a, Tuple<double, double, Interval[]> b) {
+        var res = a.Item1.CompareTo(b.Item1);
+        if (res == 0) {
+          res = b.Item2.CompareTo(a.Item2);
+        }
+        return res;
+      }
+
+      var q = new SortedSet<Tuple<double, double, Interval[]>>(Comparer<Tuple<double, double, Interval[]>>.Create(compare)) {
+        Tuple.Create(getBound(f(x)), getVolume(x), x)
+      };
+
+
+      var bestBound = double.MaxValue;
+
+      // examine all the ordered pairs in the cartesian product
+      void next_pair(int i) {
+        if (i == splits.Length) {
+          var tmp = newbox.ToArray(); // make a copy to put in the queue
+          q.Add(Tuple.Create(getBound(f(tmp)), getVolume(tmp), tmp));
+          return;
+        }
+
+        foreach (var iv in splits[i]) {
+          newbox[i] = iv;
+          next_pair(i + 1);
+        }
+      }
+
+      while (q.Count > 0 && n-- > 0) {
+        var currentBound = q.Min; q.Remove(currentBound);
+        var bound = currentBound.Item1;
+        var box = currentBound.Item3;
+
+        if (!box.Any(b => b.Width > w)) {
+          bestBound = Math.Min(bestBound, bound);
+          continue;
+        }
+
+        // do the splits
+        for (int i = 0; i < box.Length; ++i) {
+          splits[i].Clear();
+          var iv = box[i];
+          if (iv.Width > w) {
+            var t = iv.Split();
+            splits[i].Add(t.Item1);
+            splits[i].Add(t.Item2);
+          } else {
+            splits[i].Add(iv);
+          }
+        }
+        next_pair(0);
+      }
+      if (q.Count > 0) {
+        bestBound = Math.Min(bestBound, q.First().Item1);
+      }
+      return m ? bestBound : -bestBound;
+    }
+    #endregion
+  }
+}
