source: branches/GeneralizedQAP/HeuristicLab.Problems.GeneralizedQuadraticAssignment/3.3/Operators/Crossovers/CordeauCrossover.cs @ 15504

#region License Information
/* HeuristicLab
 * Copyright (C) 2002-2017 Heuristic and Evolutionary Algorithms Laboratory (HEAL)
 *
 * This file is part of HeuristicLab.
 *
 * HeuristicLab is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * HeuristicLab is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with HeuristicLab. If not, see <http://www.gnu.org/licenses/>.
 */
#endregion

using System;
using System.Linq;
using HeuristicLab.Common;
using HeuristicLab.Core;
using HeuristicLab.Data;
using HeuristicLab.Encodings.IntegerVectorEncoding;
using HeuristicLab.Parameters;
using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;
using HeuristicLab.Random;

namespace HeuristicLab.Problems.GeneralizedQuadraticAssignment {
  [Item("CordeauCrossover", "The merge procedure that is described in Cordeau, J.-F., Gaudioso, M., Laporte, G., Moccia, L. 2006. A memetic heuristic for the generalized quadratic assignment problem. INFORMS Journal on Computing, 18, pp. 433–443.")]
  [StorableClass]
  public class CordeauCrossover : GQAPCrossover,
    IQualitiesAwareGQAPOperator, IProblemInstanceAwareGQAPOperator {

    public ILookupParameter<BoolValue> MaximizationParameter {
      get { return (ILookupParameter<BoolValue>)Parameters["Maximization"]; }
    }
    public IScopeTreeLookupParameter<DoubleValue> QualityParameter {
      get { return (IScopeTreeLookupParameter<DoubleValue>)Parameters["Quality"]; }
    }
    public IScopeTreeLookupParameter<Evaluation> EvaluationParameter {
      get { return (IScopeTreeLookupParameter<Evaluation>)Parameters["Evaluation"]; }
    }
    public ILookupParameter<IntValue> EvaluatedSolutionsParameter {
      get { return (ILookupParameter<IntValue>)Parameters["EvaluatedSolutions"]; }
    }

    [StorableConstructor]
    protected CordeauCrossover(bool deserializing) : base(deserializing) { }
    protected CordeauCrossover(CordeauCrossover original, Cloner cloner)
      : base(original, cloner) {
    }
    public CordeauCrossover()
      : base() {
      Parameters.Add(new LookupParameter<BoolValue>("Maximization", ""));
      Parameters.Add(new ScopeTreeLookupParameter<DoubleValue>("Quality", "The quality of the parents", 1));
      Parameters.Add(new ScopeTreeLookupParameter<Evaluation>("Evaluation", GQAP.EvaluationDescription, 1));
      Parameters.Add(new LookupParameter<IntValue>("EvaluatedSolutions", "The number of evaluated solutions."));
    }

    public override IDeepCloneable Clone(Cloner cloner) {
      return new CordeauCrossover(this, cloner);
    }

    public static IntegerVector Apply(IRandom random, bool maximization,
      IntegerVector parent1, DoubleValue quality1,
      IntegerVector parent2, DoubleValue quality2,
      GQAPInstance problemInstance, IntValue evaluatedSolutions) {
      var distances = problemInstance.Distances;
      var capacities = problemInstance.Capacities;
      var demands = problemInstance.Demands;

      var medianDistances = Enumerable.Range(0, distances.Rows).Select(x => distances.GetRow(x).Median()).ToArray();

      int m = capacities.Length;
      int n = demands.Length;
      
      bool onefound = false;
      double fbest, fbestAttempt = maximization ? double.MinValue : double.MaxValue;
      IntegerVector bestAttempt = null;
      IntegerVector result = null;

      fbest = quality1.Value;
      if (maximization && quality1.Value < quality2.Value
        || !maximization && quality1.Value > quality2.Value) {
        var temp = parent1;
        parent1 = parent2;
        parent2 = temp;
        fbest = quality2.Value;
      }
      var cap = new double[m];
      for (var i = 0; i < m; i++) {
        int unassigned;
        Array.Clear(cap, 0, m);
        var child = Merge(parent1, parent2, distances, demands, medianDistances, m, n, i, cap, out unassigned);
        if (unassigned > 0)
          TryRandomAssignment(random, demands, capacities, m, n, cap, child, ref unassigned);
        if (unassigned == 0) {
          var childFit = problemInstance.ToSingleObjective(problemInstance.Evaluate(child));
          evaluatedSolutions.Value++;
          if (maximization && childFit >= fbest
            || !maximization && childFit <= fbest) {
            fbest = childFit;
            result = child;
            onefound = true;
          }
          if (!onefound && (maximization && fbestAttempt < childFit || !maximization && fbestAttempt > childFit)) {
            bestAttempt = child;
            fbestAttempt = childFit;
          }
        }
      }

      if (!onefound) {
        var i = random.Next(m);
        int unassigned;
        Array.Clear(cap, 0, m);
        var child = Merge(parent1, parent2, distances, demands, medianDistances, m, n, i, cap, out unassigned);
        RandomAssignment(random, demands, capacities, m, n, cap, child, ref unassigned);

        var childFit = problemInstance.ToSingleObjective(problemInstance.Evaluate(child));
        evaluatedSolutions.Value++;
        if (childFit < fbest) {
          fbest = childFit;
          result = child;
          onefound = true;
        }

        if (!onefound && (maximization && fbestAttempt < childFit || !maximization && fbestAttempt > childFit)) {
          bestAttempt = child;
          fbestAttempt = childFit;
        }
      }
        /*if (tabufix(&son, 0.5 * sqrt(n * m), round(n * m * log10(n)), &tabufix_it)) {
          solution_cost(&son);
          if (son.cost < fbest) {
            fbest = son.cost;
            *sptr = son;
            onefound = TRUE;
            merge_fixed++;
          }*/
      return result ?? bestAttempt;
    }

    private static IntegerVector Merge(IntegerVector p1, IntegerVector p2,
      DoubleMatrix distances, DoubleArray demands, double[] mediana,
      int m, int n, int i, double[] cap, out int unassigned) {
      unassigned = n;
      var child = new IntegerVector(n);
      for (var k = 0; k < n; k++) {
        child[k] = -1;
        var ik1 = p1[k];
        var ik2 = p2[k];
        if (distances[i, ik1] < mediana[i]) {
          child[k] = ik1;
          cap[ik1] += demands[k];
          unassigned--;
        } else if (distances[i, ik2] > mediana[i]) {
          child[k] = ik2;
          cap[ik2] += demands[k];
          unassigned--;
        };
      }
      return child;
    }

    private static bool TryRandomAssignment(IRandom random, DoubleArray demands, DoubleArray capacities, int m, int n, double[] cap, IntegerVector child, ref int unassigned) {
      var unbiasedOrder = Enumerable.Range(0, n).Shuffle(random).ToList();
      for (var idx = 0; idx < n; idx++) {
        var k = unbiasedOrder[idx];
        if (child[k] < 0) {
          var feasibleInserts = Enumerable.Range(0, m)
            .Select((v, i) => new { Pos = i, Slack = capacities[i] - cap[i] })
            .Where(x => x.Slack >= demands[k]).ToList();
          if (feasibleInserts.Count == 0) return false;
          var j = feasibleInserts.SampleRandom(random).Pos;
          child[k] = j;
          cap[j] += demands[k];
          unassigned--;
        }
      }
      return true;
    }

    private static void RandomAssignment(IRandom random, DoubleArray demands, DoubleArray capacities, int m, int n, double[] cap, IntegerVector child, ref int unassigned) {
      for (var k = 0; k < n; k++) {
        if (child[k] < 0) {
          var j = random.Next(m);
          child[k] = j;
          cap[j] += demands[k];
          unassigned--;
        }
      }
    }

    protected override IntegerVector Cross(IRandom random, ItemArray<IntegerVector> parents,
      GQAPInstance problemInstance) {
      if (parents == null) throw new ArgumentNullException("parents");
      if (parents.Length != 2) throw new ArgumentException(Name + " works only with exactly two parents.");

      var qualities = QualityParameter.ActualValue;
      return Apply(random, MaximizationParameter.ActualValue.Value,
        parents[0], qualities[0],
        parents[1], qualities[1],
        problemInstance,
        EvaluatedSolutionsParameter.ActualValue);
    }
  }
}