source: trunk/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs @ 17963

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

#3121: support integer powers in interval arithmetic (and two simplification rules)

  • Property svn:mergeinfo set to (toggle deleted branches)
    /branches/2971_named_intervals/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.csmergedeligible
    /branches/3073_IA_constraint_splitting/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.csmergedeligible
    /stable/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.csmergedeligible
    /trunk/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.csmergedeligible
    /branches/2435-alglib_3_15/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs16683-17036
    /branches/2520_PersistenceReintegration/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs16451-16564
    /branches/2839_HiveProjectManagement/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs15761-16116
    /branches/2892_LR-prediction-intervals/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs15743-16388
    /branches/2904_CalculateImpacts/3.4/Implementation/Interval/Interval.cs15808-16421
    /branches/2915-AbsoluteSymbol/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs15943-16355
    /branches/2947_ConfigurableIndexedDataTable/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs16148-16526
    /branches/2965_CancelablePersistence/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs16321-16439
    /branches/2966_interval_calculation/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs16320-16406
    /branches/3073_IA_constraint_splitting_reintegration/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs17885-17900
    /branches/Async/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs13329-15286
    /branches/Classification-Extensions/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs11606-11761
    /branches/ClassificationModelComparison/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs9073-13099
    /branches/CloningRefactoring/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs4656-4721
    /branches/DataAnalysis Refactoring/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs5471-5808
    /branches/DataAnalysis SolutionEnsembles/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs5815-6180
    /branches/DataAnalysis/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs4220,​4226,​4236-4238,​4389,​4458-4459,​4462,​4464
    /branches/DataAnalysisCSVImport/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs8713-8875
    /branches/DataPreprocessing/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs10085-11101
    /branches/DatasetFeatureCorrelation/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs8035-8538
    /branches/GP.Grammar.Editor/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs6284-6795
    /branches/GP.Symbols (TimeLag, Diff, Integral)/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs5060
    /branches/HeuristicLab.DatasetRefactor/sources/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs11570-12508
    /branches/HeuristicLab.Problems.Orienteering/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs11130-12721
    /branches/HeuristicLab.RegressionSolutionGradientView/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs13819-14091
    /branches/HeuristicLab.TimeSeries/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs7098-8789
    /branches/LogResidualEvaluator/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs10202-10483
    /branches/NET40/sources/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs5138-5162
    /branches/ParallelEngine/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs5175-5192
    /branches/ProblemInstancesRegressionAndClassification/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs7570-7810
    /branches/QAPAlgorithms/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs6350-6627
    /branches/Restructure trunk solution/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs6828
    /branches/SensitivityEvaluator/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs12212-16784
    /branches/SimplifierViewsProgress/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs15318-15370
    /branches/SpectralKernelForGaussianProcesses/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs10204-10479
    /branches/Trunk/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs6829-6865
    /branches/histogram/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs5959-6341
    /branches/symbreg-factors-2650/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval/Interval.cs14232-14825
    /trunk/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/Interval.cs17301-17367
File size: 14.0 KB
Line 
1#region License Information
2/* HeuristicLab
3 * Copyright (C) Heuristic and Evolutionary Algorithms Laboratory (HEAL)
4 *
5 * This file is part of HeuristicLab.
6 *
7 * HeuristicLab is free software: you can redistribute it and/or modify
8 * it under the terms of the GNU General Public License as published by
9 * the Free Software Foundation, either version 3 of the License, or
10 * (at your option) any later version.
11 *
12 * HeuristicLab is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15 * GNU General Public License for more details.
16 *
17 * You should have received a copy of the GNU General Public License
18 * along with HeuristicLab. If not, see <http://www.gnu.org/licenses/>.
19 */
20#endregion
21
22using System;
23using System.Collections.Generic;
24using System.Linq;
25using HEAL.Attic;
26using HeuristicLab.Common;
27
28namespace HeuristicLab.Problems.DataAnalysis {
29  [StorableType("849e42d3-8934-419d-9aff-64ad81c06b67")]
30  public class Interval : IEquatable<Interval> {
31    [Storable]
32    public double LowerBound { get; private set; }
33    [Storable]
34    public double UpperBound { get; private set; }
35
36    public double Width => UpperBound - LowerBound;
37
38    [StorableConstructor]
39    protected Interval(StorableConstructorFlag _) { }
40
41    /// <summary>
42    /// Creates an interval with given bounds, where lower bound must be smaller than
43    /// the upper bound. Floating point precision errors trough calculations are fixed by,
44    /// checking if the intervals are almost the same (E-12). If this is the case, the bounds
45    /// will be set to the bound closer to zero.
46    /// </summary>
47    /// <param name="lowerBound">Lower bound of the interval</param>
48    /// <param name="upperBound">Upper bound of the interval</param>
49    public Interval(double lowerBound, double upperBound) {
50      if (lowerBound.IsAlmost(upperBound)) {
51        //If the bounds go over zero
52        if (lowerBound <= 0 && upperBound >= 0) {
53          lowerBound = 0.0;
54          upperBound = 0.0;
55          //Interval is negative
56        } else if (upperBound < 0) {
57          lowerBound = upperBound;
58          //Interval is positive
59        } else {
60          upperBound = lowerBound;
61        }
62      }
63
64      if (lowerBound > upperBound)
65        throw new ArgumentException("lowerBound must be smaller than or equal to upperBound.");
66
67      this.LowerBound = lowerBound;
68      this.UpperBound = upperBound;
69    }
70
71    private Interval(double v) : this(v, v) { }
72
73    public bool Contains(double value) {
74      return LowerBound <= value && value <= UpperBound;
75    }
76
77    public bool Contains(Interval other) {
78      if (double.IsNegativeInfinity(LowerBound) && double.IsPositiveInfinity(UpperBound)) return true;
79      if (other.LowerBound >= LowerBound && other.UpperBound <= UpperBound) return true;
80
81      return false;
82    }
83
84    public override string ToString() {
85      return "Interval: [" + LowerBound + ", " + UpperBound + "]";
86    }
87
88    public bool IsInfiniteOrUndefined {
89      get {
90        return double.IsInfinity(LowerBound) || double.IsInfinity(UpperBound) ||
91                double.IsNaN(LowerBound) || double.IsNaN(UpperBound);
92      }
93    }
94
95    /// <summary>
96    /// True if the interval is positive without zero
97    /// </summary>
98    public bool IsPositive {
99      get => LowerBound > 0.0;
100    }
101
102    /// <summary>
103    /// True if the interval is negative without zero
104    /// </summary>
105    public bool IsNegative {
106      get => UpperBound < 0.0;
107    }
108
109    public static Interval GetInterval(IEnumerable<double> values) {
110      if (values == null) throw new ArgumentNullException("values");
111      if (!values.Any()) throw new ArgumentException($"No values are present.");
112
113      var min = double.MaxValue;
114      var max = double.MinValue;
115
116      foreach (var value in values) {
117        //If an value is NaN return an interval [NaN, NaN]
118        if (double.IsNaN(value)) return new Interval(double.NaN, double.NaN);
119
120        if (value < min) min = value;
121        if (value > max) max = value;
122      }
123
124      return new Interval(min, max);
125    }
126
127    #region Equals, GetHashCode, == , !=
128    public bool Equals(Interval other) {
129      if (other == null)
130        return false;
131
132      return (UpperBound == other.UpperBound || (double.IsNaN(UpperBound) && double.IsNaN(other.UpperBound)))
133        && (LowerBound == other.LowerBound || (double.IsNaN(LowerBound) && double.IsNaN(other.LowerBound)));
134    }
135
136    public override bool Equals(object obj) {
137      return Equals(obj as Interval);
138    }
139
140    public override int GetHashCode() {
141      return LowerBound.GetHashCode() ^ UpperBound.GetHashCode();
142    }
143
144    public static bool operator ==(Interval interval1, Interval interval2) {
145      if (ReferenceEquals(interval1, null)) return ReferenceEquals(interval2, null);
146      return interval1.Equals(interval2);
147    }
148    public static bool operator !=(Interval interval1, Interval interval2) {
149      return !(interval1 == interval2);
150    }
151    #endregion
152
153    #region operations
154
155    // [x1,x2] + [y1,y2] = [x1 + y1,x2 + y2]
156    public static Interval Add(Interval a, Interval b) {
157      return new Interval(a.LowerBound + b.LowerBound, a.UpperBound + b.UpperBound);
158    }
159
160    // [x1,x2] − [y1,y2] = [x1 − y2,x2 − y1]
161    public static Interval Subtract(Interval a, Interval b) {
162      return new Interval(a.LowerBound - b.UpperBound, a.UpperBound - b.LowerBound);
163    }
164
165    // [x1,x2] * [y1,y2] = [min(x1*y1,x1*y2,x2*y1,x2*y2),max(x1*y1,x1*y2,x2*y1,x2*y2)]
166    public static Interval Multiply(Interval a, Interval b) {
167      double v1 = a.LowerBound * b.LowerBound;
168      double v2 = a.LowerBound * b.UpperBound;
169      double v3 = a.UpperBound * b.LowerBound;
170      double v4 = a.UpperBound * b.UpperBound;
171
172      double min = Math.Min(Math.Min(v1, v2), Math.Min(v3, v4));
173      double max = Math.Max(Math.Max(v1, v2), Math.Max(v3, v4));
174      return new Interval(min, max);
175    }
176
177    //Division by intervals containing 0 is implemented as defined in
178    //http://en.wikipedia.org/wiki/Interval_arithmetic
179    public static Interval Divide(Interval a, Interval b) {
180      if (b.Contains(0.0)) {
181        if (b.LowerBound.IsAlmost(0.0)) return Interval.Multiply(a, new Interval(1.0 / b.UpperBound, double.PositiveInfinity));
182        else if (b.UpperBound.IsAlmost(0.0)) return Interval.Multiply(a, new Interval(double.NegativeInfinity, 1.0 / b.LowerBound));
183        else return new Interval(double.NegativeInfinity, double.PositiveInfinity);
184      }
185      return Interval.Multiply(a, new Interval(1.0 / b.UpperBound, 1.0 / b.LowerBound));
186    }
187
188    public static Interval Sine(Interval a) {
189      if (Math.Abs(a.UpperBound - a.LowerBound) >= Math.PI * 2) return new Interval(-1, 1);
190
191      //divide the interval by PI/2 so that the optima lie at x element of N (0,1,2,3,4,...)
192      double Pihalf = Math.PI / 2;
193      Interval scaled = Interval.Divide(a, new Interval(Pihalf, Pihalf));
194      //move to positive scale
195      if (scaled.LowerBound < 0) {
196        int periodsToMove = Math.Abs((int)scaled.LowerBound / 4) + 1;
197        scaled = Interval.Add(scaled, new Interval(periodsToMove * 4, periodsToMove * 4));
198      }
199
200      double scaledLowerBound = scaled.LowerBound % 4.0;
201      double scaledUpperBound = scaled.UpperBound % 4.0;
202      if (scaledUpperBound < scaledLowerBound) scaledUpperBound += 4.0;
203      List<double> sinValues = new List<double>();
204      sinValues.Add(Math.Sin(scaledLowerBound * Pihalf));
205      sinValues.Add(Math.Sin(scaledUpperBound * Pihalf));
206
207      int startValue = (int)Math.Ceiling(scaledLowerBound);
208      while (startValue < scaledUpperBound) {
209        sinValues.Add(Math.Sin(startValue * Pihalf));
210        startValue += 1;
211      }
212
213      return new Interval(sinValues.Min(), sinValues.Max());
214    }
215    public static Interval Cosine(Interval a) {
216      return Interval.Sine(Interval.Add(a, new Interval(Math.PI / 2, Math.PI / 2)));
217    }
218    public static Interval Tangens(Interval a) {
219      return Interval.Divide(Interval.Sine(a), Interval.Cosine(a));
220    }
221    public static Interval HyperbolicTangent(Interval a) {
222      return new Interval(Math.Tanh(a.LowerBound), Math.Tanh(a.UpperBound));
223    }
224
225    public static Interval Logarithm(Interval a) {
226      return new Interval(Math.Log(a.LowerBound), Math.Log(a.UpperBound));
227    }
228    public static Interval Exponential(Interval a) {
229      return new Interval(Math.Exp(a.LowerBound), Math.Exp(a.UpperBound));
230    }
231
232    public static Interval Square(Interval a) {
233      if (a.UpperBound <= 0) return new Interval(a.UpperBound * a.UpperBound, a.LowerBound * a.LowerBound);     // interval is negative
234      else if (a.LowerBound >= 0) return new Interval(a.LowerBound * a.LowerBound, a.UpperBound * a.UpperBound); // interval is positive
235      else return new Interval(0, Math.Max(a.LowerBound * a.LowerBound, a.UpperBound * a.UpperBound)); // interval goes over zero
236    }
237
238    public static Interval Cube(Interval a) {
239      return new Interval(Math.Pow(a.LowerBound, 3), Math.Pow(a.UpperBound, 3));
240    }
241
242    public static Interval Power(Interval a, int b) {
243      if (b < 0) return Power(1.0 / a, -b); // a^(-b) = 1/(a^b)
244      if (b == 0 && (a.Contains(0.0) || a.IsInfiniteOrUndefined)) return new Interval(double.NaN, double.NaN);  // 0^0, +/-inf^0 are undefined
245      if (b == 0) return new Interval(1.0, 1.0); // x^0 = 1
246      if (b == 1) return a;
247      if (b % 2 == 0) {
248        // even powers (see x²)
249        if (a.UpperBound <= 0) return new Interval(Math.Pow(a.UpperBound, b), Math.Pow(a.LowerBound, b));     // interval is negative
250        if (a.LowerBound >= 0) return new Interval(Math.Pow(a.LowerBound, b), Math.Pow(a.UpperBound, b)); // interval is positive
251        return new Interval(0, Math.Max(Math.Pow(a.LowerBound, b), Math.Pow(a.UpperBound, b))); // interval goes over zero
252      } else {
253        // odd powers (see x³)
254        return new Interval(Math.Pow(a.LowerBound, b), Math.Pow(a.UpperBound, b));
255      }
256    }
257
258    /// <summary>
259    /// The interval contains both possible results of the calculated square root +-sqrt(x). That results in a wider
260    /// interval, but it contains all possible solutions.
261    /// </summary>
262    /// <param name="a">Interval to build square root from.</param>
263    /// <returns></returns>
264    public static Interval SquareRoot(Interval a) {
265      if (a.LowerBound < 0) return new Interval(double.NaN, double.NaN);
266      return new Interval(-Math.Sqrt(a.UpperBound), Math.Sqrt(a.UpperBound));
267    }
268
269    public static Interval CubicRoot(Interval a) {
270      var lower = (a.LowerBound < 0) ? -Math.Pow(-a.LowerBound, 1d / 3d) : Math.Pow(a.LowerBound, 1d / 3d);
271      var upper = (a.UpperBound < 0) ? -Math.Pow(-a.UpperBound, 1d / 3d) : Math.Pow(a.UpperBound, 1d / 3d);
272
273      return new Interval(lower, upper);
274    }
275
276    public static Interval Absolute(Interval a) {
277      var absLower = Math.Abs(a.LowerBound);
278      var absUpper = Math.Abs(a.UpperBound);
279      var min = Math.Min(absLower, absUpper);
280      var max = Math.Max(absLower, absUpper);
281
282      if (a.Contains(0.0)) {
283        min = 0.0;
284      }
285
286      return new Interval(min, max);
287    }
288
289    public static Interval AnalyticQuotient(Interval a, Interval b) {
290      var dividend = a;
291      var divisor = Add(Square(b), new Interval(1.0, 1.0));
292      divisor = SquareRoot(divisor);
293
294      var quotient = Divide(dividend, divisor);
295      return quotient;
296    }
297    #endregion
298
299    #region arithmetic overloads
300    public static Interval operator +(Interval a, Interval b) => Add(a, b);
301    public static Interval operator +(Interval a, double b) => Add(a, new Interval(b));
302    public static Interval operator +(double a, Interval b) => Add(new Interval(a), b);
303    public static Interval operator -(Interval a, Interval b) => Subtract(a, b);
304    public static Interval operator -(Interval a, double b) => Subtract(a, new Interval(b));
305    public static Interval operator -(double a, Interval b) => Subtract(new Interval(a), b);
306    public static Interval operator -(Interval a) => Subtract(new Interval(0), a);
307    public static Interval operator *(Interval a, Interval b) => Multiply(a, b);
308    public static Interval operator *(Interval a, double b) => Multiply(a, new Interval(b));
309    public static Interval operator *(double a, Interval b) => Multiply(new Interval(a), b);
310    public static Interval operator /(Interval a, Interval b) => Divide(a, b);
311    public static Interval operator /(Interval a, double b) => Divide(a, new Interval(b));
312    public static Interval operator /(double a, Interval b) => Divide(new Interval(a), b);
313    public static Interval Exponential(double a) { return Exponential(new Interval(a)); }
314    public static Interval Logarithm(double a) { return Logarithm(new Interval(a)); }
315    public static Interval Sine(double a) { return Sine(new Interval(a)); }
316    public static Interval Cosine(double a) { return Cosine(new Interval(a)); }
317    public static Interval Tangens(double a) { return Tangens(new Interval(a)); }
318    public static Interval HyperbolicTangent(double a) { return HyperbolicTangent(new Interval(a)); }
319    public static Interval Square(double a) { return Square(new Interval(a)); }
320    public static Interval Cube(double a) { return Cube(new Interval(a)); }
321    public static Interval SquareRoot(double a) { return SquareRoot(new Interval(a)); }
322    public static Interval CubicRoot(double a) { return CubicRoot(new Interval(a)); }
323    public static Interval Absolute(double a) { return Absolute(new Interval(a)); }
324    public static Interval AnalyticQuotient(Interval a, double b) { return AnalyticQuotient(a, new Interval(b)); }
325    public static Interval AnalyticQuotient(double a, Interval b) { return AnalyticQuotient(new Interval(a), b); }
326    public static Interval AnalyticQuotient(double a, double b) { return AnalyticQuotient(new Interval(a), new Interval(b)); }
327    #endregion
328  }
329}
Note: See TracBrowser for help on using the repository browser.