Changeset 17727 for branches/2825-NSGA3

Timestamp:

08/26/20 21:27:42 (4 years ago)

Author:

dleko

Message:

#2825 Refactoring.

Location:

branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3

Files:

• NSGA3.cs (modified) (22 diffs)
• NSGA3Selection.cs (modified) (9 diffs)
• ReferencePoint.cs (modified) (2 diffs)
• SimulatedBinaryCrossover.cs (modified) (3 diffs)
• Solution.cs (modified) (2 diffs)
• Utility.cs (modified) (2 diffs)

branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/NSGA3.cs

@@ -1 +1 @@
 using System;
 using System.Collections.Generic;
+using System.IO;
 using System.Linq;
 using System.Threading;
@@ -11 +12 @@
 using HeuristicLab.Parameters;
 using HeuristicLab.Problems.TestFunctions.MultiObjective;
+using HeuristicLab.Problems.TestFunctions.MultiObjective.TestFunctions;
 using HeuristicLab.Random;
 
@@ -66 +68 @@
         private NSGA3Selection selection;
 
+        [Storable]
+        private double[] allIgds;
+
         #endregion Storable fields
 
@@ -75 +80 @@
         private const string MaximumGenerationsName = "Maximum Generations";
         private const string CrossoverProbabilityName = "Crossover Probability";
+        private const string CrossoverEtaName = "Crossover eta";
         private const string MutationProbabilityName = "Mutation Probability";
         private const string DominateOnEqualQualitiesName = "Dominate On Equal Qualities";
@@ -92 +98 @@
         private const string ScatterPlotResultName = "Scatter Plot";
         private const string CurrentFrontResultName = "Pareto Front"; // Do not touch this
+        private const string RandomSeedUsedName = "Random Seed Used";
 
         #endregion ParameterAndResultsNames
@@ -112 +119 @@
         }
 
+        private IFixedValueParameter<DoubleValue> CrossoverEtaParameter
+        {
+            get { return (IFixedValueParameter<DoubleValue>)Parameters[CrossoverEtaName]; }
+        }
+
         private IFixedValueParameter<PercentValue> MutationProbabilityParameter
         {
@@ -146 +158 @@
 
         public PercentValue CrossoverProbability => CrossoverProbabilityParameter.Value;
+
+        public DoubleValue CrossoverEta => CrossoverEtaParameter.Value;
 
         public PercentValue MutationProbability => MutationProbabilityParameter.Value;
@@ -213 +227 @@
             get { return (DoubleMatrix)Results[CurrentFrontResultName].Value; }
             set { Results[CurrentFrontResultName].Value = value; }
+        }
+
+        public IntValue ResultsRandomSeedUsed
+        {
+            get { return (IntValue)Results[RandomSeedUsedName].Value; }
+            set { Results[RandomSeedUsedName].Value = value; }
         }
 
@@ -224 +244 @@
             Parameters.Add(new FixedValueParameter<IntValue>(MaximumGenerationsName, "The maximum number of generations which should be processed.", new IntValue(1000)));
             Parameters.Add(new FixedValueParameter<PercentValue>(CrossoverProbabilityName, "The probability that the crossover operator is applied on two parents.", new PercentValue(1.0)));
+            Parameters.Add(new FixedValueParameter<DoubleValue>(CrossoverEtaName, "TODO: description of this parameter", new DoubleValue(30)));
             Parameters.Add(new FixedValueParameter<PercentValue>(MutationProbabilityName, "The probability that the mutation operator is applied on a Individual.", new PercentValue(0.05)));
             Parameters.Add(new FixedValueParameter<BoolValue>(DominateOnEqualQualitiesName, "Flag which determines wether Individuals with equal quality values should be treated as dominated.", new BoolValue(true)));
@@ -231 +252 @@
         }
 
-        // Persistence uses this ctor to improve deserialization efficiency. If we would use the
-        // default ctor instead this would completely initialize the object (e.g. creating
-        // parameters) even though the data is later overwritten by the stored data.
         [StorableConstructor]
         public NSGA3(StorableConstructorFlag _) : base(_) { }
@@ -245 +263 @@
             referencePoints = original.referencePoints?.Select(cloner.Clone).ToList();
             selection = cloner.Clone(original.selection);
+            if (original.allIgds != null)
+            {
+                allIgds = new double[original.allIgds.Length];
+                Array.Copy(original.allIgds, allIgds, allIgds.Length);
+            }
         }
 
@@ -262 +285 @@
             InitResults();
             InitFields();
-            Analyze();
+            Analyze(false);
         }
 
@@ -276 +299 @@
             Results.Add(new Result(ScatterPlotResultName, "A scatterplot displaying the evaluated solutions and (if available) the analytically optimal front", new ParetoFrontScatterPlot()));
             Results.Add(new Result(CurrentFrontResultName, "The Pareto Front", new DoubleMatrix()));
+            Results.Add(new Result(RandomSeedUsedName, "The random seed used in the end", new IntValue(-1)));
 
             if (!(Problem is MultiObjectiveTestFunctionProblem problem)) return;
-            // todo: add BestKnownFront parameter
-            ResultsScatterPlot = new ParetoFrontScatterPlot(new double[0][], new double[0][], problem.BestKnownFront.ToJaggedArray(), problem.Objectives, problem.ProblemSize);
+            var bestKnownFront = problem.BestKnownFront.ToJaggedArray();
             if (problem.BestKnownFront == null) return;
+            ResultsScatterPlot = new ParetoFrontScatterPlot(new double[0][], new double[0][], bestKnownFront, problem.Objectives, problem.ProblemSize);
             ResultsBestKnownHypervolume = new DoubleValue(Hypervolume.Calculate(problem.BestKnownFront.ToJaggedArray(), problem.ReferencePoint.CloneAsArray(), problem.Maximization));
         }
@@ -286 +310 @@
         private void InitFields()
         {
-            random = new MersenneTwister();
-            if (!SetSeedRandomly.Value)
-                random.Reset(Seed.Value);
+            if (SetSeedRandomly.Value)
+            {
+                System.Random seedRand = new System.Random();
+                ResultsRandomSeedUsed.Value = seedRand.Next();
+                random = new MersenneTwister(Convert.ToUInt32(ResultsRandomSeedUsed.Value));
+            }
+            else
+            {
+                ResultsRandomSeedUsed.Value = Seed.Value;
+                random = new MersenneTwister(Convert.ToUInt32(Seed.Value));
+            }
 
             solutions = GetInitialPopulation();
@@ -296 +328 @@
 
             selection = new NSGA3Selection(problem.Objectives, problem.Maximization, PopulationSize.Value, random, DominateOnEqualQualities.Value);
+
+            allIgds = new double[MaximumGenerations.Value];
         }
 
@@ -311 +345 @@
                 RealVector randomRealVector = new RealVector(Problem.Encoding.Length, random, minBound, maxBound);
                 var solution = new Solution(randomRealVector);
-                solution.Fitness = Evaluate(solution.Chromosome);
+                solution.Objectives = Evaluate(solution.Chromosome);
                 solutions.Add(solution);
             }
@@ -328 +362 @@
                 while (ResultsCurrentGeneration.Value < MaximumGenerations.Value)
                 {
-                    // todo: make parameter out of this
                     List<Solution> qt = Mutate(Recombine(solutions));
                     foreach (var solution in qt)
-                        solution.Fitness = Evaluate(solution.Chromosome);
+                    {
+                        solution.Objectives = Evaluate(solution.Chromosome);
+                    }
 
                     List<Solution> rt = Utility.Concat(solutions, qt);
@@ -337 +372 @@
                     solutions = selection.SelectSolutionsForNextGeneration(rt, GetCopyOfReferencePoints());
 
-                    if (AnalyzeEveryGeneration.Value)
-                        Analyze();
+                    Analyze(!AnalyzeEveryGeneration.Value);
 
                     ResultsCurrentGeneration.Value++;
@@ -354 +388 @@
             finally
             {
-                Analyze();
-            }
+                Analyze(false);
+                WriteIGDsToFile();
+            }
+        }
+
+        private void WriteIGDsToFile()
+        {
+            string resultsDirectory = "Results";
+            Directory.CreateDirectory(resultsDirectory);
+            var problem = Problem as MultiObjectiveTestFunctionProblem;
+            string fileName = Path.Combine(resultsDirectory, $"{problem.TestFunction.Name}({problem.Objectives})_DOEQ={DominateOnEqualQualities.Value}_Seed={ResultsRandomSeedUsed.Value}.txt");
+
+            File.WriteAllLines(fileName, allIgds.Select(igd => igd.ToString()));
         }
 
@@ -369 +414 @@
         }
 
-        private void Analyze()
-        {
-            ResultsScatterPlot = new ParetoFrontScatterPlot(solutions.Select(x => x.Fitness).ToArray(), solutions.Select(x => x.Chromosome.ToArray()).ToArray(), ResultsScatterPlot.ParetoFront, ResultsScatterPlot.Objectives, ResultsScatterPlot.ProblemSize);
-            ResultsSolutions = solutions.Select(s => s.Chromosome.ToArray()).ToMatrix();
+        private void Analyze(bool calculateOnlyIGD)
+        {
+            if (!calculateOnlyIGD)
+            {
+                ResultsScatterPlot = new ParetoFrontScatterPlot(solutions.Select(x => x.Objectives).ToArray(), solutions.Select(x => x.Chromosome.ToArray()).ToArray(), ResultsScatterPlot.ParetoFront, ResultsScatterPlot.Objectives, ResultsScatterPlot.ProblemSize);
+                ResultsSolutions = solutions.Select(s => s.Chromosome.ToArray()).ToMatrix();
+            }
 
             if (!(Problem is MultiObjectiveTestFunctionProblem problem)) return;
 
-            // Indicators
-            ResultsGenerationalDistance = new DoubleValue(problem.BestKnownFront != null ? GenerationalDistance.Calculate(solutions.Select(s => s.Fitness), problem.BestKnownFront.ToJaggedArray(), 1) : double.NaN);
-            ResultsInvertedGenerationalDistance = new DoubleValue(problem.BestKnownFront != null ? InvertedGenerationalDistance.Calculate(solutions.Select(s => s.Fitness), problem.BestKnownFront.ToJaggedArray(), 1) : double.NaN);
-
-            var front = NonDominatedSelect.GetDominatingVectors(solutions.Select(x => x.Fitness), problem.ReferencePoint.CloneAsArray(), Problem.Maximization, true).ToArray();
-            if (front.Length == 0) return;
-            ResultsHypervolume = new DoubleValue(Hypervolume.Calculate(front, problem.ReferencePoint.CloneAsArray(), problem.Maximization));
-            ResultsDifferenceToBestKnownHypervolume = new DoubleValue(ResultsBestKnownHypervolume.Value - ResultsHypervolume.Value);
-
-            Problem.Analyze(
-                solutions.Select(s => (Individual)new SingleEncodingIndividual(Problem.Encoding, new Scope { Variables = { new Variable(Problem.Encoding.Name, s.Chromosome) } })).ToArray(),
-                solutions.Select(s => s.Fitness).ToArray(),
-                Results,
-                random);
+            CalculateInvertedGenerationalDistance(problem);
+
+            if (!calculateOnlyIGD)
+            {
+                ResultsGenerationalDistance = new DoubleValue(problem.BestKnownFront != null ? GenerationalDistance.Calculate(solutions.Select(s => s.Objectives), problem.BestKnownFront.ToJaggedArray(), 1) : double.NaN);
+                var front = NonDominatedSelect.GetDominatingVectors(solutions.Select(x => x.Objectives), problem.ReferencePoint.CloneAsArray(), Problem.Maximization, true).ToArray();
+                if (front.Length == 0) return;
+                ResultsHypervolume = new DoubleValue(Hypervolume.Calculate(front, problem.ReferencePoint.CloneAsArray(), problem.Maximization));
+                ResultsDifferenceToBestKnownHypervolume = new DoubleValue(ResultsBestKnownHypervolume.Value - ResultsHypervolume.Value);
+
+                Problem.Analyze(
+                    solutions.Select(s => (Individual)new SingleEncodingIndividual(Problem.Encoding, new Scope { Variables = { new Variable(Problem.Encoding.Name, s.Chromosome) } })).ToArray(),
+                    solutions.Select(s => s.Objectives).ToArray(),
+                    Results,
+                    random);
+            }
+        }
+
+        private void CalculateInvertedGenerationalDistance(MultiObjectiveTestFunctionProblem problem)
+        {
+            double igd;
+            if (problem.TestFunction is ScaledDTLZ dtlz)
+            {
+                // Normalize objectives before calculating the igd
+                var scalingFactor = dtlz.GetScalingFactor(problem.Objectives);
+                var front = new double[solutions.Count][];
+
+                for (int i = 0; i < solutions.Count; i++)
+                {
+                    front[i] = new double[problem.Objectives];
+                    Array.Copy(solutions[i].Objectives, front[i], front[i].Length);
+                    for (int objective = 0; objective < problem.Objectives; objective++)
+                    {
+                        // get normalized value
+                        front[i][objective] = front[i][objective] / Math.Pow(scalingFactor, objective);
+                    }
+                }
+
+                // Normalize best known front before calculating the igd
+                var bestKnownFront = problem.BestKnownFront.ToJaggedArray();
+                foreach (var point in bestKnownFront)
+                {
+                    for (int objective = 0; objective < problem.Objectives; objective++)
+                    {
+                        point[objective] = point[objective] / Math.Pow(scalingFactor, objective);
+                    }
+                }
+
+                igd = InvertedGenerationalDistance.Calculate(front, bestKnownFront, 1);
+            }
+            else
+            {
+                igd = problem.BestKnownFront != null ? InvertedGenerationalDistance.Calculate(solutions.Select(s => s.Objectives), problem.BestKnownFront.ToJaggedArray(), 1) : double.NaN;
+            }
+            ResultsInvertedGenerationalDistance = new DoubleValue(igd);
+            allIgds[ResultsCurrentGeneration.Value - 1] = ResultsInvertedGenerationalDistance.Value;
         }
 
@@ -403 +493 @@
         }
 
-        // todo: try both recombination methods
         private List<Solution> Recombine(List<Solution> solutions)
         {
-            //List<Solution> childSolutions = new List<Solution>();
-            //for (int i = 0; i < solutions.Count; i++)
-            //{
-            //    var parent1 = solutions[random.Next(solutions.Count)];
-            //    var parent2 = solutions[random.Next(solutions.Count)];
-            //    var childChromosome = HeuristicLab.Encodings.RealVectorEncoding.SimulatedBinaryCrossover.Apply(random, parent1.Chromosome, parent2.Chromosome, new DoubleValue(5));
-
-            // BoundsChecker.Apply(childChromosome, Problem.Encoding.Bounds);
-
-            //    Solution child = new Solution(childChromosome);
-            //    childSolutions.Add(child);
-            //}
-
-            //return childSolutions;
-
             List<Solution> childSolutions = new List<Solution>();
 
             for (int i = 0; i < solutions.Count; i += 2)
             {
+                if (solutions.Count != PopulationSize.Value) throw new Exception("solutions.Count != PopulationSize");
                 int parentIndex1 = random.Next(solutions.Count);
                 int parentIndex2 = random.Next(solutions.Count);
@@ -434 +509 @@
                 // Do crossover with crossoverProbabilty == 1 in order to guarantee that a crossover happens
                 var children = SimulatedBinaryCrossover.Apply(random,
-                    Problem.Encoding.Bounds, parent1.Chromosome, parent2.Chromosome, CrossoverProbability.Value);
+                    Problem.Encoding.Bounds, parent1.Chromosome, parent2.Chromosome, CrossoverProbability.Value, CrossoverEta.Value);
 
                 childSolutions.Add(new Solution(children.Item1));

branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/NSGA3Selection.cs

@@ -67 +67 @@
 
         /// <summary>
-        /// TODO: Description If the fitnesses of all solutions in <paramref name="rt" /> do not
-        /// have the same number of objectives, undefined behavior happens.
+        /// Select the solutions for the next generation from <see cref="rt" /> according to the
+        /// NSGA-III algorithm.
         /// </summary>
-        /// <param name="rt">The previous generation and its mutated children.</param>
+        /// <param name="rt">The current generation and its mutated children.</param>
         /// <param name="referencePoints"></param>
         /// <returns></returns>
@@ -84 +84 @@
                 .Select(list => new List<Solution>(list.Select(pair => pair.Item1))).ToList();
 
-            int lf = -1; // last front to be included
-            for (int i = 0; i < fronts.Count && st.Count < populationSize; i++)
-            {
-                st = Utility.Concat(st, fronts[i]);
-                lf = i;
-            }
+            int last = 0; // last front to be included
+            while (st.Count < populationSize)
+            {
+                st = Utility.Concat(st, fronts[last++]);
+            }
+
             // remove useless fronts
-            fronts.RemoveRange(lf + 1, fronts.Count - lf - 1);
+            fronts.RemoveRange(last, fronts.Count - last);
 
             List<Solution> nextGeneration = new List<Solution>();
-            if (st.Count == populationSize) // no special selection needs to be done
+            if (st.Count == populationSize)
+            { // no special selection needs to be done
                 nextGeneration = st;
+            }
             else
             {
                 // Add every front to the next generation list except for the one that is added only partially
                 nextGeneration = new List<Solution>();
-                for (int f = 0; f < lf; f++)
+                for (int f = 0; f < last - 1; f++)
                     nextGeneration = Utility.Concat(nextGeneration, fronts[f]);
+
                 int k = populationSize - nextGeneration.Count;
                 Normalize(rt, fronts);
@@ -115 +118 @@
         {
             // Find the ideal point
+            double[] idealPoint = TranslateObjectives(solutions, fronts[0]);
+
+            // Find the extreme points
+            List<Solution> extremePoints = GetExtremePoints(fronts[0]);
+
+            // Compute intercepts
+            double[] intercepts = GetIntercepts(extremePoints, fronts[0]);
+
+            // Normalize objectives
+            NormalizeObjectives(solutions, intercepts, idealPoint);
+        }
+
+        private double[] TranslateObjectives(List<Solution> solutions, List<Solution> firstFront)
+        {
             double[] idealPoint = new double[numberOfObjectives];
-
             foreach (var solution in solutions)
-                solution.ConvertedFitness = new double[numberOfObjectives];
+            {
+                if (solution.NormalizedObjectives == null) solution.NormalizedObjectives = new double[numberOfObjectives];
+            }
 
             for (int j = 0; j < numberOfObjectives; j++)
             {
                 // Compute ideal point
-                // todo: search only in first front
-                idealPoint[j] = Utility.Min(s => s.Fitness[j], solutions);
-                var x = Utility.Min(s => s.Fitness[j], fronts.First());
-                if (Math.Abs(idealPoint[j] - x) > 10e-5) throw new Exception("idealPoint value is not in first front");
-
-                // Translate objectives and save the modified fitness values in ConvertedFitness
+                idealPoint[j] = Utility.Min(s => s.Objectives[j], firstFront);
+
+                // Translate objectives and save the modified fitness values in NormalizedObjectives
                 foreach (var solution in solutions)
-                    solution.ConvertedFitness[j] = solution.Fitness[j] - idealPoint[j];
-            }
-
-            // Find the extreme points
+                    solution.NormalizedObjectives[j] = solution.Objectives[j] - idealPoint[j];
+            }
+
+            return idealPoint;
+        }
+
+        private List<Solution> GetExtremePoints(List<Solution> firstFront)
+        {
             List<Solution> extremePoints = new List<Solution>(numberOfObjectives);
             for (int j = 0; j < numberOfObjectives; j++)
@@ -142 +161 @@
                 weights[j] = 1;
 
-                var extremePoint = Utility.ArgMin(s => ASF(s.ConvertedFitness, weights), solutions);
+                var extremePoint = Utility.ArgMin(s => ASF(s.NormalizedObjectives, weights), firstFront);
                 extremePoints.Add(extremePoint);
             }
 
-            // Compute intercepts
-            double[] intercepts = GetIntercepts(extremePoints, solutions);
-
-            // Normalize objectives
-            NormalizeObjectives(solutions, intercepts, idealPoint);
+            return extremePoints;
+        }
+
+        private double ASF(double[] convertedFitness, double[] weight)
+        {
+            var dims = Enumerable.Range(0, numberOfObjectives);
+
+            return Utility.Max((dim) => convertedFitness[dim] / weight[dim], dims);
+        }
+
+        private double[] GetIntercepts(List<Solution> extremePoints, List<Solution> solutions)
+        {
+            // Check whether there are duplicate extreme points. This might happen but the original
+            // paper does not mention how to deal with it.
+            bool duplicate = false;
+            for (int i = 0; !duplicate && i < extremePoints.Count; i++)
+            {
+                for (int j = i + 1; !duplicate && j < extremePoints.Count; j++)
+                {
+                    duplicate = extremePoints[i] == extremePoints[j];
+                }
+            }
+
+            double[] intercepts = new double[numberOfObjectives];
+
+            bool negative_intercept = false;
+            if (!duplicate)
+            {
+                // Find the equation of the hyperplane
+                List<double> b = new List<double>(Enumerable.Repeat(1.0, numberOfObjectives));
+
+                List<List<double>> a = new List<List<double>>(extremePoints.Count);
+
+                for (int i = 0; i < extremePoints.Count; i++)
+                {
+                    List<double> convertedFitness = new List<double>(extremePoints[i].NormalizedObjectives);
+                    a.Add(convertedFitness);
+                }
+                List<double> x = GaussianElimination(a, b);
+
+                // Find intercepts
+                for (int i = 0; i < intercepts.Length; i++)
+                {
+                    intercepts[i] = 1.0 / x[i];
+
+                    if (x[i] < 0)
+                    {
+                        negative_intercept = true;
+                        break;
+                    }
+                }
+            }
+
+            // follow the method in Yuan et al. (GECCO 2015)
+            if (duplicate || negative_intercept)
+            {
+                double[] maxObjectives = FindMaxObjectives(solutions);
+                for (int i = 0; i < intercepts.Length; i++)
+                    intercepts[i] = maxObjectives[i];
+            }
+
+            return intercepts;
+        }
+
+        private List<double> GaussianElimination(List<List<double>> a, List<double> b)
+        {
+            List<double> x = new List<double>();
+
+            int n = a.Count;
+            for (int i = 0; i < n; i++)
+                a[i].Add(b[i]);
+
+            for (int @base = 0; @base < n - 1; @base++)
+                for (int target = @base + 1; target < n; target++)
+                {
+                    double ratio = a[target][@base] / a[@base][@base];
+                    for (int term = 0; term < a[@base].Count; term++)
+                        a[target][term] = a[target][term] - a[@base][term] * ratio;
+                }
+
+            for (int i = 0; i < n; i++)
+                x.Add(0.0);
+
+            for (int i = n - 1; i >= 0; i--)
+            {
+                for (int known = i + 1; known < n; known++)
+                    a[i][n] = a[i][n] - a[i][known] * x[known];
+                x[i] = a[i][n] / a[i][i];
+            }
+
+            return x;
+        }
+
+        private double[] FindMaxObjectives(List<Solution> solutions)
+        {
+            double[] maxPoint = new double[numberOfObjectives];
+            for (int i = 0; i < maxPoint.Length; i++) maxPoint[i] = double.MinValue;
+
+            for (int f = 0; f < numberOfObjectives; f++)
+            {
+                maxPoint[f] = Utility.Max(s => s.Objectives[f], solutions);
+            }
+
+            return maxPoint;
         }
 
@@ -159 +277 @@
                 for (int i = 0; i < numberOfObjectives; i++)
                 {
-                    // if (Math.Abs(intercepts[i] - idealPoints[i]) > EPSILON)
-                    if (Math.Abs(intercepts[i]) > EPSILON)
-                    {
-                        //solution.ConvertedFitness[i] = solution.ConvertedFitness[i] / (intercepts[i] - idealPoint[i]);
-                        // todo: try this
-                        solution.ConvertedFitness[i] = solution.ConvertedFitness[i] / intercepts[i];
+                    if (Math.Abs(intercepts[i] - idealPoint[i]) > 10e-10)
+                    {
+                        solution.NormalizedObjectives[i] = solution.NormalizedObjectives[i] / (intercepts[i] - idealPoint[i]);
                     }
                     else
                     {
-                        solution.ConvertedFitness[i] = solution.ConvertedFitness[i] / EPSILON;
+                        solution.NormalizedObjectives[i] = solution.NormalizedObjectives[i] / 10e-10;
                     }
                 }
@@ -182 +297 @@
                     // find reference point for which the perpendicular distance to the current
                     // solution is the lowest
-                    var rpAndDist = Utility.MinArgMin(rp => GetPerpendicularDistance(rp.Values, solution.ConvertedFitness), referencePoints);
+                    var rpAndDist = Utility.MinArgMin(rp => GetPerpendicularDistance(rp.Values, solution.NormalizedObjectives), referencePoints);
                     // associated reference point
                     var arp = rpAndDist.Item1;
@@ -190 +305 @@
                     if (f + 1 != fronts.Count)
                     {
-                        // Todo: Add member for reference point on index min_rp
                         arp.NumberOfAssociatedSolutions++;
                     }
                     else
                     {
-                        // Todo: Add potential member for reference point on index min_rp
                         arp.AddPotentialAssociatedSolution(solution, dist);
                     }
                 }
             }
+        }
+
+        private double GetPerpendicularDistance(double[] direction, double[] point)
+        {
+            double numerator = 0;
+            double denominator = 0;
+            for (int i = 0; i < direction.Length; i++)
+            {
+                numerator += direction[i] * point[i];
+                denominator += Math.Pow(direction[i], 2);
+            }
+            double k = numerator / denominator;
+
+            double d = 0;
+            for (int i = 0; i < direction.Length; i++)
+            {
+                d += Math.Pow(k * direction[i] - point[i], 2);
+            }
+            return Math.Sqrt(d);
         }
 
@@ -230 +362 @@
             int minNumber = Utility.Min(rp => rp.NumberOfAssociatedSolutions, referencePoints);
 
-            // todo: try this instead of the current implementation
-            //var min_rps = new List<ReferencePoint>(referencePoints.Where(r => r.NumberOfAssociatedSolutions == minNumber));
-            //if (min_rps.Count == 0) return min_rps.Single();
-            //else return min_rps[random.Next(min_rps.Count)];
-
-            // the reference points that share the number of associated solutions where that number
-            // is equal to minNumber
-            List<ReferencePoint> minAssociatedReferencePoints = new List<ReferencePoint>();
-            foreach (var referencePoint in referencePoints)
-                if (referencePoint.NumberOfAssociatedSolutions == minNumber)
-                    minAssociatedReferencePoints.Add(referencePoint);
-
-            if (minAssociatedReferencePoints.Count > 1)
-                return minAssociatedReferencePoints[random.Next(minAssociatedReferencePoints.Count)];
-            else
-                return minAssociatedReferencePoints.Single();
+            var minReferencePoints = new List<ReferencePoint>(referencePoints.Where(r => r.NumberOfAssociatedSolutions == minNumber));
+            if (minReferencePoints.Count == 1) return minReferencePoints.Single();
+            else return minReferencePoints[random.Next(minReferencePoints.Count)];
         }
 
@@ -260 +379 @@
             return chosen;
         }
-
-        private double GetPerpendicularDistance(double[] values, double[] fitness)
-        {
-            double numerator = 0;
-            double denominator = 0;
-            for (int i = 0; i < values.Length; i++)
-            {
-                numerator += values[i] * fitness[i];
-                denominator += Math.Pow(values[i], 2);
-            }
-            double k = numerator / denominator;
-
-            double d = 0;
-            for (int i = 0; i < values.Length; i++)
-            {
-                d += Math.Pow(k * values[i] - fitness[i], 2);
-            }
-            return Math.Sqrt(d);
-        }
-
-        private double ASF(double[] x, double[] weight)
-        {
-            List<int> dimensions = new List<int>();
-            for (int i = 0; i < numberOfObjectives; i++)
-                dimensions.Add(i);
-
-            return Utility.Max((dim) => x[dim] / weight[dim], dimensions);
-        }
-
-        private double[] GetIntercepts(List<Solution> extremePoints, List<Solution> solutions)
-        {
-            // Check whether there are duplicate extreme points. This might happen but the original
-            // paper does not mention how to deal with it.
-            bool duplicate = false;
-            for (int i = 0; !duplicate && i < extremePoints.Count; i++)
-            {
-                for (int j = i + 1; !duplicate && j < extremePoints.Count; j++)
-                {
-                    duplicate = extremePoints[i] == extremePoints[j];
-                }
-            }
-
-            double[] intercepts = new double[numberOfObjectives];
-
-            bool negative_intercept = false;
-            if (!duplicate)
-            {
-                // Find the equation of the hyperplane
-                List<double> b = new List<double>(numberOfObjectives);
-                for (int i = 0; i < numberOfObjectives; i++)
-                    b.Add(1.0);
-
-                List<List<double>> a = new List<List<double>>(extremePoints.Count);
-
-                for (int i = 0; i < extremePoints.Count; i++)
-                {
-                    List<double> convertedFitness = new List<double>(extremePoints[i].ConvertedFitness);
-                    a.Add(convertedFitness);
-                }
-                List<double> x = GaussianElimination(a, b);
-
-                // Find intercepts
-                for (int i = 0; i < intercepts.Length; i++)
-                {
-                    intercepts[i] = 1.0 / x[i];
-
-                    if (x[i] < 0)
-                    {
-                        negative_intercept = true;
-                        break;
-                    }
-                }
-            }
-
-            // follow the method in Yuan et al. (GECCO 2015)
-            if (duplicate || negative_intercept)
-            {
-                double[] maxObjectives = FindMaxObjectives(solutions);
-                for (int i = 0; i < intercepts.Length; i++)
-                    intercepts[i] = maxObjectives[i];
-            }
-
-            return intercepts;
-        }
-
-        private double[] FindMaxObjectives(List<Solution> solutions)
-        {
-            double[] maxPoint = new double[numberOfObjectives];
-            for (int i = 0; i < maxPoint.Length; i++) maxPoint[i] = double.MinValue;
-
-            foreach (var solution in solutions)
-            {
-                for (int f = 0; f < numberOfObjectives; f++)
-                {
-                    maxPoint[f] = Math.Max(maxPoint[f], solution.ConvertedFitness[f]);
-                }
-            }
-
-            return maxPoint;
-        }
-
-        private List<double> GaussianElimination(List<List<double>> a, List<double> b)
-        {
-            List<double> x = new List<double>();
-
-            int n = a.Count;
-            for (int i = 0; i < n; i++)
-                a[i].Add(b[i]);
-
-            for (int @base = 0; @base < n - 1; @base++)
-                for (int target = @base + 1; target < n; target++)
-                {
-                    double ratio = a[target][@base] / a[@base][@base];
-                    for (int term = 0; term < a[@base].Count; term++)
-                        a[target][term] = a[target][term] - a[@base][term] * ratio;
-                }
-
-            for (int i = 0; i < n; i++)
-                x.Add(0.0);
-
-            for (int i = n - 1; i >= 0; i--)
-            {
-                for (int known = i + 1; known < n; known++)
-                    a[i][n] = a[i][n] - a[i][known] * x[known];
-                x[i] = a[i][n] / a[i][i];
-            }
-
-            return x;
-        }
     }
 }

branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/ReferencePoint.cs

@@ -227 +227 @@
         }
 
-        internal Solution FindClosestMember()
+        public Solution FindClosestMember()
         {
             if (!potentialAssociatedSolutions.Any())
@@ -235 +235 @@
         }
 
-        internal Solution RandomMember()
+        public Solution RandomMember()
         {
             if (!potentialAssociatedSolutions.Any())

branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/SimulatedBinaryCrossover.cs

@@ -11 +11 @@
      *
      * The simulated binary crossover in HeuristicLab.Encodings.RealVectorEncoding/Crossovers/SimulatedBinaryCrossover.cs
-     * is insufficient for NSGA3, because the static Apply-method only returns 1 child.
+     * is not usable, because it takes different parameters than stated in the NSGA-III paper. It takes contiguity as a
+     * parameter instead of an eta value.
      */
 
@@ -18 +19 @@
         private const double EPSILON = 10e-6; // a tiny number that is greater than 0
 
-        public static Tuple<RealVector, RealVector> Apply(IRandom random, DoubleMatrix bounds, RealVector parent1, RealVector parent2, double crossoverProbability = 1.0, double eta = 30)
+        public static Tuple<RealVector, RealVector> Apply(IRandom random, DoubleMatrix bounds, RealVector parent1, RealVector parent2, double crossoverProbability = 1.0, double eta = 20)
         {
             var child1 = new RealVector(parent1);
@@ -33 +34 @@
                 double y2 = Math.Max(parent1[i], parent2[i]);
 
-                // todo: is this correct? test it with Encoding.Length != Number of objectives
                 double lb;
                 double ub;

branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/Solution.cs

@@ -15 +15 @@
         // actual fitness of solution as given by Problem
         [Storable]
-        public double[] Fitness { get; set; }
+        public double[] Objectives { get; set; }
 
         // normalized fitness used in selection process (in order to not overwrite the original Fitness)
         [Storable]
-        public double[] ConvertedFitness { get; set; }
+        public double[] NormalizedObjectives { get; set; }
 
         [StorableConstructor]
@@ -32 +32 @@
         {
             Chromosome = cloner.Clone(solution.Chromosome);
-            if (solution.Fitness != null)
+            if (solution.Objectives != null)
             {
-                Fitness = new double[solution.Fitness.Length];
-                Array.Copy(solution.Fitness, Fitness, solution.Fitness.Length);
+                Objectives = new double[solution.Objectives.Length];
+                Array.Copy(solution.Objectives, Objectives, solution.Objectives.Length);
             }
-            if (solution.ConvertedFitness != null)
+            if (solution.NormalizedObjectives != null)
             {
-                ConvertedFitness = new double[solution.ConvertedFitness.Length];
-                Array.Copy(solution.ConvertedFitness, ConvertedFitness, solution.ConvertedFitness.Length);
+                NormalizedObjectives = new double[solution.NormalizedObjectives.Length];
+                Array.Copy(solution.NormalizedObjectives, NormalizedObjectives, solution.NormalizedObjectives.Length);
             }
         }

branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/Utility.cs

@@ -69 +69 @@
         {
             if (m == null) return null;
-            var i = m.Rows - 1;
-            var a = new double[i][];
-            for (i--; i >= 0; i--)
-            {
-                var j = m.Columns;
-                a[i] = new double[j];
-                for (j--; j >= 0; j--) a[i][j] = m[i, j];
+            var a = new double[m.Rows][];
+            for (int i = m.Rows - 1; i >= 0; i--)
+            {
+                a[i] = new double[m.Columns];
+                for (int j = m.Columns - 1; j >= 0; j--) a[i][j] = m[i, j];
             }
             return a;
@@ -91 +89 @@
             {
                 Solution solution = solutions[i];
-                data[i] = new double[solution.Fitness.Length];
-                for (int j = 0; j < solution.Fitness.Length; j++)
+                data[i] = new double[solution.Objectives.Length];
+                for (int j = 0; j < solution.Objectives.Length; j++)
                 {
-                    data[i][j] = solution.Fitness[j];
+                    data[i][j] = solution.Objectives[j];
                 }
             }