Changeset 17294


Ignore:
Timestamp:
10/02/19 17:06:47 (12 days ago)
Author:
gkronber
Message:

#2994: min() and max() are problematic in AutoDiff. We have used max(a,b) = (a+b + abs(b-a)) / 2 but this is problematic because of catastrophic cancellation. Min() and Max() are necessary for interval calculation. However, there we only need to determine min and max for double values (or duals). We do not need min() and max() in general.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2994-AutoDiffForIntervals/HeuristicLab.Problems.DataAnalysis.Symbolic/3.4/Interpreter/Interpreter.cs

    r17292 r17294  
    422422      // we assume that the for variable nodes ( v(x,w) = w * x ) the gradient is returned for parameter w
    423423      instruction.value = new AlgebraicInterval(
    424         low: new MultivariateDual<AlgebraicDouble>(v: lower, key: variable, dv: intervals[variable.VariableName].LowerBound), 
     424        low: new MultivariateDual<AlgebraicDouble>(v: lower, key: variable, dv: intervals[variable.VariableName].LowerBound),
    425425        high: new MultivariateDual<AlgebraicDouble>(v: upper, key: variable, dv: intervals[variable.VariableName].UpperBound)
    426426        );
     
    433433
    434434  public interface IAlgebraicType<T> {
    435     T Zero { get; }
     435    T Zero { get; } // Zero and One must create new objects
    436436    T One { get; }
    437437
     
    453453    T AssignIntRoot(T a, int r);
    454454    T AssignSgn(T a); // set this to sign(a)
     455    T AssignMin(T other); // set this min(this, other)
     456    T AssignMax(T other); // set this max(this, other)
    455457    T Clone();
    456458  }
     
    468470    public static T IntRoot<T>(this T a, int r) where T : IAlgebraicType<T> { a.AssignIntRoot(a.Clone(), r); return a; }
    469471
    470     public static T Max<T>(T a, T b) where T : IAlgebraicType<T> {
    471       // ((a + b) + abs(b - a)) / 2
    472       return a.Clone().Add(b).Add(b.Clone().Sub(a).Abs()).Scale(1.0 / 2.0);
    473     }
    474     public static T Min<T>(T a, T b) where T : IAlgebraicType<T> {
    475       // ((a + b) - abs(a - b)) / 2
    476       return a.Clone().Add(b).Sub(a.Clone().Sub(b).Abs()).Scale(1.0 / 2.0);
    477     }
     472    internal static T Min<T>(T a, T b) where T : IAlgebraicType<T> { return a.Clone().AssignMin(b); }
     473    internal static T Max<T>(T a, T b) where T : IAlgebraicType<T> { return a.Clone().AssignMax(b); }
     474
     475    // public static T Max<T>(T a, T b) where T : IAlgebraicType<T> {
     476    //   // ((a + b) + abs(b - a)) / 2
     477    //   return a.Clone().Add(b).Add(b.Clone().Sub(a).Abs()).Scale(1.0 / 2.0);
     478    // }
     479    // public static T Min<T>(T a, T b) where T : IAlgebraicType<T> {
     480    //   // ((a + b) - abs(a - b)) / 2
     481    //   return a.Clone().Add(b).Sub(a.Clone().Sub(b).Abs()).Scale(1.0 / 2.0);
     482    // }
    478483  }
    479484
     
    511516    public AlgebraicDouble AssignIntPower(AlgebraicDouble a, int p) { Value = Math.Pow(a.Value, p); return this; }
    512517    public AlgebraicDouble AssignIntRoot(AlgebraicDouble a, int r) { Value = IntRoot(a.Value, r); return this; }
     518    public AlgebraicDouble AssignMin(AlgebraicDouble other) { Value = Math.Min(Value, other.Value); return this; }
     519    public AlgebraicDouble AssignMax(AlgebraicDouble other) { Value = Math.Max(Value, other.Value); return this; }
    513520
    514521    // helper
     
    563570    public AlgebraicDoubleVector AssignIntPower(AlgebraicDoubleVector a, int p) { for (int i = 0; i < arr.Length; ++i) { arr[i] = Math.Pow(a.arr[i], p); } return this; }
    564571    public AlgebraicDoubleVector AssignIntRoot(AlgebraicDoubleVector a, int r) { for (int i = 0; i < arr.Length; ++i) { arr[i] = IntRoot(a.arr[i], r); } return this; }
     572    public AlgebraicDoubleVector AssignMin(AlgebraicDoubleVector other) { for (int i = 0; i < arr.Length; ++i) { arr[i] = Math.Min(arr[i], other.arr[i]); } return this; }
     573    public AlgebraicDoubleVector AssignMax(AlgebraicDoubleVector other) { for (int i = 0; i < arr.Length; ++i) { arr[i] = Math.Max(arr[i], other.arr[i]); } return this; }
    565574
    566575    // helper
     
    778787    }
    779788
     789    public AlgebraicSparseVector<K, T> AssignMin(AlgebraicSparseVector<K, T> other) {
     790      // assumes that keys without a matching key in other are zero and vice versa
     791      foreach (var kvp in elems) if (!other.elems.ContainsKey(kvp.Key)) kvp.Value.AssignMin(kvp.Value.Zero); // min(v, 0)
     792      foreach (var kvp in other.elems) {
     793        if (elems.TryGetValue(kvp.Key, out T value))
     794          value.AssignMin(kvp.Value);
     795        else
     796          elems.Add(kvp.Key, kvp.Value.Zero.AssignMin(kvp.Value));
     797      }
     798      return this;
     799    }
     800
     801    public AlgebraicSparseVector<K, T> AssignMax(AlgebraicSparseVector<K, T> other) {
     802      // assumes that keys without a matching key in other are zero and vice versa
     803      foreach (var kvp in elems) if (!other.elems.ContainsKey(kvp.Key)) kvp.Value.AssignMax(kvp.Value.Zero); // max(v, 0)
     804      foreach (var kvp in other.elems) {
     805        if (elems.TryGetValue(kvp.Key, out T value))
     806          value.AssignMax(kvp.Value);
     807        else
     808          elems.Add(kvp.Key, kvp.Value.Zero.AssignMax(kvp.Value));
     809      }
     810      return this;
     811    }
     812
     813
    780814    public AlgebraicSparseVector<K, T> Clone() {
    781815      return new AlgebraicSparseVector<K, T>(this);
     
    836870    }
    837871
    838     // algebraic min() / max() do not work for infinities
    839     // detect and handle infinite values explicitly
     872
    840873    private static MultivariateDual<AlgebraicDouble> Min(MultivariateDual<AlgebraicDouble> a, MultivariateDual<AlgebraicDouble> b) {
    841       if (a.Value.IsInfinity || b.Value.IsInfinity) return a.Value < b.Value ? a : b; // NOTE: this is not differentiable but for infinite values we do not care
    842       else return Algebraic.Min(a, b); // differentiable statement
     874      return a.Value < b.Value ? a : b;
    843875    }
    844876    private static MultivariateDual<AlgebraicDouble> Max(MultivariateDual<AlgebraicDouble> a, MultivariateDual<AlgebraicDouble> b) {
    845       if (a.Value.IsInfinity || b.Value.IsInfinity) return a.Value >= b.Value ? a : b; // NOTE: this is not differentiable but for infinite values we do not care
    846       else return Algebraic.Max(a, b); // differentiable statement
     877      return a.Value > b.Value ? a : b;
    847878    }
    848879
     
    919950        if (a.Contains(0.0)) {
    920951          low = new MultivariateDual<AlgebraicDouble>(0.0);
    921           high = Algebraic.Max(a.low.IntPower(p), a.high.IntPower(p));
     952          high = a.low.IntPower(p).AssignMax(a.high.IntPower(p));
    922953        } else {
    923954          var lowPower = a.low.IntPower(p);
    924955          var highPower = a.high.IntPower(p);
    925           low = Algebraic.Min(lowPower, highPower);
    926           high = Algebraic.Max(lowPower, highPower);
     956          low = lowPower.AssignMin(highPower);
     957          high = lowPower.AssignMax(highPower);
    927958        }
    928959      } else {
     
    934965          var lowPower = a.low.IntPower(p);
    935966          var highPower = a.high.IntPower(p);
    936           low = Algebraic.Min(lowPower, highPower);
    937           high = Algebraic.Max(lowPower, highPower);
     967          low = lowPower.AssignMin(highPower);
     968          high = lowPower.AssignMax(highPower);
    938969        }
    939970      }
     
    10641095      return this;
    10651096    }
     1097
     1098    public AlgebraicInterval AssignMin(AlgebraicInterval other) {
     1099      low.AssignMin(other.low);
     1100      high.AssignMin(other.high);
     1101      return this;
     1102    }
     1103
     1104    public AlgebraicInterval AssignMax(AlgebraicInterval other) {
     1105      low.AssignMax(other.low);
     1106      high.AssignMax(other.high);
     1107      return this;
     1108    }
    10661109  }
    10671110
     
    11161159    public Dual<V> Clone() { return new Dual<V>(v.Clone(), dv.Clone()); }
    11171160
     1161    public Dual<V> AssignMin(Dual<V> other) {
     1162      throw new NotImplementedException();
     1163    }
     1164
     1165    public Dual<V> AssignMax(Dual<V> other) {
     1166      throw new NotImplementedException();
     1167    }
    11181168  }
    11191169
     
    11991249    public MultivariateDual<V> AssignSgn(MultivariateDual<V> a) { v.AssignSgn(a.v); dv = a.dv.Zero; return this; } // sign(f(x))' = 0;     
    12001250
     1251    public MultivariateDual<V> AssignMin(MultivariateDual<V> other) {
     1252      XXX
     1253    }
     1254
     1255    public MultivariateDual<V> AssignMax(MultivariateDual<V> other) {
     1256      XXX
     1257    }
    12011258  }
    12021259}
Note: See TracChangeset for help on using the changeset viewer.