Ignore:
Timestamp:
11/29/16 15:46:48 (4 years ago)
Author:
abeham
Message:

#2701, #2708: Made a new branch from ProblemRefactoring and removed ScopedBasicAlgorithm branch (which becomes MemPR branch)

File:
1 moved

Legend:

Unmodified
Added
Removed
  • branches/PerformanceComparison/HeuristicLab.OptimizationExpertSystem.Common/3.3/ProblemCharacteristicAnalysis/QAP/QAPDirectedWalk.cs

    r14428 r14429  
    8585
    8686    public override IEnumerable<IResult> Calculate() {
     87      IRandom random = Seed.HasValue ? new MersenneTwister((uint)Seed.Value) : new MersenneTwister();
     88      var qap = (QuadraticAssignmentProblem)Problem;
    8789      var pathCount = Paths;
    88       var random = Seed.HasValue ? new MersenneTwister((uint)Seed.Value) : new MersenneTwister();
    89       var bestImprovement = BestImprovement;
    90       var qap = (QuadraticAssignmentProblem)Problem;
    91       var N = qap.Weights.Rows;
    92       var start = new Permutation(PermutationTypes.Absolute, N, random);
    93       var end = start;
    94       List<double> ips = new List<double>(), ups = new List<double>(), sd2 = new List<double>();
    95       var pathsD1 = new List<List<Tuple<double, double>>>();
    96       for (var i = 0; i < pathCount; i++) {
    97         var trajD1 = new List<Tuple<double, double>>();
    98         pathsD1.Add(trajD1);
    99 
    100         if ((i + 1) % N == 0) {
    101           start = new Permutation(PermutationTypes.Absolute, N, random);
    102           end = start;
    103         }
    104         end = (Permutation)end.Clone();
    105         Rot1(end);
    106 
    107         var hist = new Tuple<Permutation, double>[5];
    108         var walkDist = 0.0;
    109         var prevVal = double.NaN;
    110         var sumD2 = 0.0;
    111         var inflectionPoints = 0;
    112         var undulationPoints = 0;
    113         var countPoints = 0;
    114         var counter = 0;
    115         var path = bestImprovement ? BestImprovementWalk(qap, start, QAPEvaluator.Apply(start, qap.Weights, qap.Distances), end)
    116                                   : FirstImprovementWalk(qap, start, QAPEvaluator.Apply(start, qap.Weights, qap.Distances), end, random);
    117         foreach (var next in path) {
    118           if (hist[0] != null) {
    119             var dist = Dist(next.Item1, hist[0].Item1);
    120             walkDist += dist;
     90
     91      var perm = new Permutation(PermutationTypes.Absolute, qap.Weights.Rows, random);
     92      var permutations = new List<Permutation> { perm };
     93      while (permutations.Count < pathCount + 1) {
     94        perm = (Permutation)permutations.Last().Clone();
     95        BiasedShuffle(perm, random);
     96        permutations.Add(perm);
     97      }
     98
     99      var trajectories = Run(random, (QuadraticAssignmentProblem)Problem, permutations, BestImprovement).ToList();
     100      var firstDerivatives = trajectories.Select(path => ApproximateDerivative(path).ToList()).ToList();
     101      var secondDerivatives = firstDerivatives.Select(d1 => ApproximateDerivative(d1).ToList()).ToList();
     102     
     103      var props = GetCharacteristics(trajectories, firstDerivatives, secondDerivatives).ToDictionary(x => x.Item1, x => x.Item2);
     104      foreach (var chara in characteristics.CheckedItems.Select(x => x.Value.Value)) {
     105        if (chara == "Swap2.Sharpness") yield return new Result("Swap2.Sharpness", new DoubleValue(props["Sharpness"]));
     106        if (chara == "Swap2.Bumpiness") yield return new Result("Swap2.Bumpiness", new DoubleValue(props["Bumpiness"]));
     107        if (chara == "Swap2.Flatness") yield return new Result("Swap2.Flatness", new DoubleValue(props["Flatness"]));
     108        if (chara == "Swap2.Steadiness") yield return new Result("Swap2.Steadiness", new DoubleValue(props["Steadiness"]));
     109      }
     110    }
     111
     112    public static IEnumerable<List<Tuple<Permutation, double>>> Run(IRandom random, QuadraticAssignmentProblem qap, IEnumerable<Permutation> permutations, bool bestImprovement = true) {
     113      var iter = permutations.GetEnumerator();
     114      if (!iter.MoveNext()) yield break;
     115
     116      var start = iter.Current;
     117      while (iter.MoveNext()) {
     118        var end = iter.Current;
     119
     120        var walk = (bestImprovement ? BestDirectedWalk(qap, start, end) : FirstDirectedWalk(random, qap, start, end)).ToList();
     121        var max = walk.Max(x => x.Item2);
     122        var min = walk.Min(x => x.Item2);
     123        if (max > min)
     124          yield return walk.Select(x => Tuple.Create(x.Item1, (x.Item2 - min) / (max - min))).ToList();
     125        else yield return walk.Select(x => Tuple.Create(x.Item1, 0.0)).ToList();
     126        start = end;
     127      } // end paths
     128    }
     129
     130    private IEnumerable<Tuple<string, double>> GetCharacteristics(List<List<Tuple<Permutation, double>>> f, List<List<Tuple<Permutation, double>>> f1, List<List<Tuple<Permutation, double>>> f2) {
     131      var sharpness = f2.Average(x => Area(x));
     132      var bumpiness = 0.0;
     133      var flatness = 0.0;
     134      var downPointing = f1.Where(x => x.Min(y => y.Item2) < 0).ToList();
     135
     136      var steadiness = 0.0;
     137      foreach (var path in downPointing) {
     138        steadiness += ComBelowZero(path);
     139      }
     140      if (downPointing.Count > 0) steadiness /= downPointing.Count;
     141
     142      for (var p = 0; p < f2.Count; p++) {
     143        var bump = 0;
     144        var flat = 0;
     145        for (var i = 0; i < f2[p].Count - 1; i++) {
     146          if ((f2[p][i].Item2 > 0 && f2[p][i + 1].Item2 < 0) || (f2[p][i].Item2 < 0 && f2[p][i + 1].Item2 > 0)) {
     147            bump++;
     148          } else if (f2[p][i].Item2 == 0) {
     149            flat++;
    121150          }
    122           // from the past 5 values we can calculate the 2nd derivative
    123           // first derivative in point 2 as differential between points 1 and 3
    124           // first derivative in point 4 as differential between points 3 and 5
    125           // second derivative in point 3 as differential between the first derivatives in points 2 and 4
    126           hist[4] = hist[3];
    127           hist[3] = hist[2];
    128           hist[2] = hist[1];
    129           hist[1] = hist[0];
    130           hist[0] = next;
    131           counter++;
    132 
    133           if (counter < 3) continue;
    134           var grad1 = (hist[0].Item2 - hist[2].Item2) / Dist(hist[0].Item1, hist[2].Item1);
    135 
    136           if (!double.IsNaN(grad1) && !double.IsInfinity(grad1))
    137             trajD1.Add(Tuple.Create(walkDist, grad1));
    138 
    139           if (counter < 5) continue;
    140           countPoints++;
    141           var grad2 = (hist[2].Item2 - hist[4].Item2) / Dist(hist[2].Item1, hist[4].Item1);
    142           var dgrad = (grad1 - grad2) / Dist(hist[1].Item1, hist[3].Item1);
    143 
    144           if (double.IsNaN(dgrad) || double.IsInfinity(dgrad)) continue;
    145           if (!double.IsNaN(prevVal)) {
    146             if (prevVal < 0 && dgrad > 0 || prevVal > 0 && dgrad < 0)
    147               inflectionPoints++;
    148             else if (prevVal.IsAlmost(0) && dgrad.IsAlmost(0)
    149                 || prevVal.IsAlmost(0) && !dgrad.IsAlmost(0)
    150                 || !prevVal.IsAlmost(0) && dgrad.IsAlmost(0))
    151               undulationPoints++;
    152           }
    153           sumD2 += Math.Abs(grad1 - grad2);
    154           prevVal = dgrad;
    155         }
    156         start = end;
    157         ips.Add(inflectionPoints / (double)countPoints);
    158         ups.Add(undulationPoints / (double)countPoints);
    159         sd2.Add(sumD2 / walkDist);
    160       } // end paths
    161       var avgZero = pathsD1.Select(path => path.SkipWhile(v => v.Item2 < 0).First().Item1 / path.Last().Item1).Median();
    162      
    163       foreach (var chara in characteristics.CheckedItems.Select(x => x.Value.Value)) {
    164         if (chara == "Swap2.Sharpness") yield return new Result("Swap2.Sharpness", new DoubleValue(sd2.Average()));
    165         if (chara == "Swap2.Bumpiness") yield return new Result("Swap2.Bumpiness", new DoubleValue(ips.Average()));
    166         if (chara == "Swap2.Flatness") yield return new Result("Swap2.Flatness", new DoubleValue(ups.Average()));
    167         if (chara == "Swap2.Steadiness") yield return new Result("Swap2.Steadiness", new DoubleValue(avgZero));
    168       }
    169     }
    170 
    171     public IEnumerable<Tuple<Permutation, double>> BestImprovementWalk(QuadraticAssignmentProblem qap, Permutation start, double fitness, Permutation end) {
     151        }
     152        bumpiness += bump / (f2[p].Count - 1.0);
     153        flatness += flat / (f2[p].Count - 1.0);
     154      }
     155      bumpiness /= f2.Count;
     156      flatness /= f2.Count;
     157      return new[] {
     158      Tuple.Create("Sharpness", sharpness),
     159      Tuple.Create("Bumpiness", bumpiness),
     160      Tuple.Create("Flatness", flatness),
     161      Tuple.Create("Steadiness", steadiness)
     162    };
     163    }
     164
     165    public static IEnumerable<Tuple<Permutation, double>> BestDirectedWalk(QuadraticAssignmentProblem qap, Permutation start, Permutation end) {
    172166      var N = qap.Weights.Rows;
    173167      var sol = start;
     168      var fitness = QAPEvaluator.Apply(sol, qap.Weights, qap.Distances);
     169      yield return Tuple.Create(sol, fitness);
     170
    174171      var invSol = GetInverse(sol);
    175172      // we require at most N-1 steps to move from one permutation to another
     
    197194    }
    198195
    199     public IEnumerable<Tuple<Permutation, double>> FirstImprovementWalk(QuadraticAssignmentProblem qap, Permutation start, double fitness, Permutation end, IRandom random) {
     196    public static IEnumerable<Tuple<Permutation, double>> FirstDirectedWalk(IRandom random, QuadraticAssignmentProblem qap, Permutation start, Permutation end) {
    200197      var N = qap.Weights.Rows;
    201198      var sol = start;
     199      var fitness = QAPEvaluator.Apply(sol, qap.Weights, qap.Distances);
     200      yield return Tuple.Create(sol, fitness);
     201
    202202      var invSol = GetInverse(sol);
    203203      // randomize the order in which improvements are tried
     
    229229    }
    230230
     231    private static double Area(IEnumerable<Tuple<Permutation, double>> path) {
     232      var iter = path.GetEnumerator();
     233      if (!iter.MoveNext()) return 0.0;
     234      var area = 0.0;
     235      var prev = iter.Current;
     236      while (iter.MoveNext()) {
     237        area += TrapezoidArea(prev, iter.Current);
     238        prev = iter.Current;
     239      }
     240      return area;
     241    }
     242
     243    private static double TrapezoidArea(Tuple<Permutation, double> a, Tuple<Permutation, double> b) {
     244      var area = 0.0;
     245      var dist = Dist(a.Item1, b.Item1);
     246      if ((a.Item2 <= 0 && b.Item2 <= 0) || (a.Item2 >= 0 && b.Item2 >= 0))
     247        area += dist * (Math.Abs(a.Item2) + Math.Abs(b.Item2)) / 2.0;
     248      else {
     249        var k = (b.Item2 - a.Item2) / dist;
     250        var d = a.Item2;
     251        var x = -d / k;
     252        area += Math.Abs(x * a.Item2 / 2.0);
     253        area += Math.Abs((dist - x) * b.Item2 / 2.0);
     254      }
     255      return area;
     256    }
     257
     258    private static double ComBelowZero(IEnumerable<Tuple<Permutation, double>> path) {
     259      var area = 0.0;
     260      var com = 0.0;
     261      var nwalkDist = 0.0;
     262      Tuple<Permutation, double> prev = null;
     263      var iter = path.GetEnumerator();
     264      while (iter.MoveNext()) {
     265        var c = iter.Current;
     266        if (prev != null) {
     267          var ndist = Dist(prev.Item1, c.Item1) / (double)c.Item1.Length;
     268          nwalkDist += ndist;
     269          if (prev.Item2 < 0 || c.Item2 < 0) {
     270            var a = TrapezoidArea(prev, c) / (double)c.Item1.Length;
     271            area += a;
     272            com += (nwalkDist - (ndist / 2.0)) * a;
     273          }
     274        }
     275        prev = c;
     276      }
     277      return com / area;
     278    }
     279
     280    private static IEnumerable<Tuple<Permutation, double>> ApproximateDerivative(IEnumerable<Tuple<Permutation, double>> data) {
     281      Tuple<Permutation, double> prev = null, prev2 = null;
     282      foreach (var d in data) {
     283        if (prev == null) {
     284          prev = d;
     285          continue;
     286        }
     287        if (prev2 == null) {
     288          prev2 = prev;
     289          prev = d;
     290          continue;
     291        }
     292        var dist = Dist(prev2.Item1, d.Item1);
     293        yield return Tuple.Create(prev.Item1, (d.Item2 - prev2.Item2) / (double)dist);
     294        prev2 = prev;
     295        prev = d;
     296      }
     297    }
     298
    231299    private static double Dist(Permutation a, Permutation b) {
    232       return a.Where((t, i) => t != b[i]).Count();
     300      var dist = 0;
     301      for (var i = 0; i < a.Length; i++)
     302        if (a[i] != b[i]) dist++;
     303      return dist;
    233304    }
    234305
    235306    private static int[] GetInverse(Permutation p) {
    236307      var inv = new int[p.Length];
    237       for (var i = 0; i < p.Length; i++) inv[p[i]] = i;
     308      for (var i = 0; i < p.Length; i++) {
     309        inv[p[i]] = i;
     310      }
    238311      return inv;
    239312    }
    240313
    241     private static void Rot1(Permutation p) {
    242       var first = p[0];
    243       for (var i = 0; i < p.Length - 1; i++) p[i] = p[i + 1];
    244       p[p.Length - 1] = first;
     314    // permutation must be strictly different in every position
     315    private static void BiasedShuffle(Permutation p, IRandom random) {
     316      for (var i = p.Length - 1; i > 0; i--) {
     317        // Swap element "i" with a random earlier element (excluding itself)
     318        var swapIndex = random.Next(i);
     319        var h = p[swapIndex];
     320        p[swapIndex] = p[i];
     321        p[i] = h;
     322      }
    245323    }
    246324  }
Note: See TracChangeset for help on using the changeset viewer.