Index: /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/NSGA3.cs
===================================================================
--- /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/NSGA3.cs	(revision 17723)
+++ /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/NSGA3.cs	(revision 17724)
@@ -63,4 +63,7 @@
         private List<ReferencePoint> referencePoints;
 
+        [Storable]
+        private NSGA3Selection selection;
+
         #endregion Storable fields
 
@@ -222,5 +225,5 @@
             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<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(false)));
+            Parameters.Add(new FixedValueParameter<BoolValue>(DominateOnEqualQualitiesName, "Flag which determines wether Individuals with equal quality values should be treated as dominated.", new BoolValue(true)));
             Parameters.Add(new FixedValueParameter<BoolValue>(SetSeedRandomlyName, "True if the random seed should be set to a random value, otherwise false.", new BoolValue(true)));
             Parameters.Add(new FixedValueParameter<IntValue>(SeedName, "The random seed used to initialize the new pseudo random number generator.", new IntValue(0)));
@@ -241,4 +244,5 @@
             solutions = original.solutions?.Select(cloner.Clone).ToList();
             referencePoints = original.referencePoints?.Select(cloner.Clone).ToList();
+            selection = cloner.Clone(original.selection);
         }
 
@@ -289,4 +293,7 @@
             referencePoints = ReferencePoint.GenerateReferencePoints(random, Objectives);
             ResultsGeneratedReferencePoints = Utility.ConvertToDoubleMatrix(referencePoints);
+            var problem = Problem as MultiObjectiveTestFunctionProblem;
+
+            selection = new NSGA3Selection(problem.Objectives, problem.Maximization, PopulationSize.Value, random, DominateOnEqualQualities.Value);
         }
 
@@ -328,5 +335,5 @@
                     List<Solution> rt = Utility.Concat(solutions, qt);
 
-                    solutions = NSGA3Selection.SelectSolutionsForNextGeneration(rt, GetCopyOfReferencePoints(), Problem.Maximization, PopulationSize.Value, random);
+                    solutions = selection.SelectSolutionsForNextGeneration(rt, GetCopyOfReferencePoints());
 
                     if (AnalyzeEveryGeneration.Value)
@@ -359,5 +366,5 @@
             if (referencePoints == null) return null;
 
-            return referencePoints.Select(rp => (ReferencePoint)rp.Clone()).ToList();
+            return referencePoints.Select(rp => new ReferencePoint(rp)).ToList();
         }
 
@@ -377,4 +384,10 @@
             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);
         }
 
@@ -393,39 +406,39 @@
         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++)
-            {
-                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);
+
+            for (int i = 0; i < solutions.Count; i += 2)
+            {
+                int parentIndex1 = random.Next(solutions.Count);
+                int parentIndex2 = random.Next(solutions.Count);
+                // ensure that the parents are not the same object
+                if (parentIndex1 == parentIndex2) parentIndex2 = (parentIndex2 + 1) % solutions.Count;
+                var parent1 = solutions[parentIndex1];
+                var parent2 = solutions[parentIndex2];
+
+                // 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);
+
+                childSolutions.Add(new Solution(children.Item1));
+                childSolutions.Add(new Solution(children.Item2));
             }
 
             return childSolutions;
-
-            //List<Solution> childSolutions = new List<Solution>();
-
-            //for (int i = 0; i < solutions.Count; i += 2)
-            //{
-            //    int parentIndex1 = random.Next(solutions.Count);
-            //    int parentIndex2 = random.Next(solutions.Count);
-            //    // ensure that the parents are not the same object
-            //    if (parentIndex1 == parentIndex2) parentIndex2 = (parentIndex2 + 1) % solutions.Count;
-            //    var parent1 = solutions[parentIndex1];
-            //    var parent2 = solutions[parentIndex2];
-
-            // // 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);
-
-            //    childSolutions.Add(new Solution(children.Item1));
-            //    childSolutions.Add(new Solution(children.Item2));
-            //}
-
-            //return childSolutions;
         }
 
Index: /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/NSGA3Selection.cs
===================================================================
--- /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/NSGA3Selection.cs	(revision 17723)
+++ /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/NSGA3Selection.cs	(revision 17724)
@@ -2,4 +2,6 @@
 using System.Collections.Generic;
 using System.Linq;
+using HEAL.Attic;
+using HeuristicLab.Common;
 using HeuristicLab.Core;
 using HeuristicLab.Optimization;
@@ -7,7 +9,60 @@
 namespace HeuristicLab.Algorithms.NSGA3
 {
-    public static class NSGA3Selection
+    [StorableType("53158EA5-4D1C-42DB-90F5-EC555EC0A450")]
+    public class NSGA3Selection : IDeepCloneable
     {
         private const double EPSILON = 10e-6; // a tiny number that is greater than 0
+
+        [Storable]
+        private int numberOfObjectives;
+
+        [Storable]
+        private bool[] maximization;
+
+        [Storable]
+        private int populationSize;
+
+        [Storable]
+        private IRandom random;
+
+        [Storable]
+        private bool dominateOnEqualQualities;
+
+        [StorableConstructor]
+        public NSGA3Selection(StorableConstructorFlag _) { }
+
+        public NSGA3Selection(int numberOfObjectives, bool[] maximization, int populationSize, IRandom random, bool dominateOnEqualQualities)
+        {
+            if (maximization == null) throw new ArgumentNullException(nameof(maximization));
+            if (random == null) throw new ArgumentNullException(nameof(random));
+            if (numberOfObjectives != maximization.Length) throw new ArgumentException($"{nameof(numberOfObjectives)} and the length of {nameof(maximization)} do not match.");
+            if (populationSize < 1) throw new ArgumentOutOfRangeException($"{nameof(populationSize)} has to be >= 1.");
+
+            this.numberOfObjectives = numberOfObjectives;
+            this.maximization = maximization;
+            this.populationSize = populationSize;
+            this.random = random;
+            this.dominateOnEqualQualities = dominateOnEqualQualities;
+        }
+
+        public NSGA3Selection(NSGA3Selection other, Cloner cloner)
+        {
+            numberOfObjectives = other.numberOfObjectives;
+            maximization = new bool[other.maximization.Length];
+            Array.Copy(other.maximization, maximization, maximization.Length);
+            populationSize = other.populationSize;
+            random = cloner.Clone(other.random);
+            dominateOnEqualQualities = other.dominateOnEqualQualities;
+        }
+
+        public IDeepCloneable Clone(Cloner cloner)
+        {
+            return new NSGA3Selection(this, cloner);
+        }
+
+        public object Clone()
+        {
+            return new Cloner().Clone(this);
+        }
 
         /// <summary>
@@ -17,11 +72,7 @@
         /// <param name="rt">The previous generation and its mutated children.</param>
         /// <param name="referencePoints"></param>
-        /// <param name="maximization"></param>
-        /// <param name="N">Number of solutions to select.</param>
-        /// <param name="random">The <see cref="IRandom" /> to use for randomness.</param>
         /// <returns></returns>
-        public static List<Solution> SelectSolutionsForNextGeneration(List<Solution> rt, List<ReferencePoint> referencePoints, bool[] maximization, int N, IRandom random)
-        {
-            int numberOfObjectives = rt.First().Fitness.Length;
+        public List<Solution> SelectSolutionsForNextGeneration(List<Solution> rt, List<ReferencePoint> referencePoints)
+        {
             List<Solution> st = new List<Solution>();
 
@@ -30,9 +81,9 @@
 
             // compute the pareto fronts using the DominationCalculator
-            var fronts = DominationCalculator<Solution>.CalculateAllParetoFronts(rt.ToArray(), qualities, maximization, out int[] rank, true)
+            List<List<Solution>> fronts = DominationCalculator<Solution>.CalculateAllParetoFronts(rt.ToArray(), qualities, maximization, out int[] rank, dominateOnEqualQualities)
                 .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 < N; i++)
+            for (int i = 0; i < fronts.Count && st.Count < populationSize; i++)
             {
                 st = Utility.Concat(st, fronts[i]);
@@ -43,5 +94,5 @@
 
             List<Solution> nextGeneration = new List<Solution>();
-            if (st.Count == N) // no special selection needs to be done
+            if (st.Count == populationSize) // no special selection needs to be done
                 nextGeneration = st;
             else
@@ -51,6 +102,6 @@
                 for (int f = 0; f < lf; f++)
                     nextGeneration = Utility.Concat(nextGeneration, fronts[f]);
-                int k = N - nextGeneration.Count;
-                Normalize(st, numberOfObjectives);
+                int k = populationSize - nextGeneration.Count;
+                Normalize(rt, fronts);
                 Associate(fronts, referencePoints);
                 List<Solution> solutionsToAdd = Niching(k, referencePoints, random);
@@ -61,10 +112,10 @@
         }
 
-        private static void Normalize(List<Solution> st, int numberOfObjectives)
+        private void Normalize(List<Solution> solutions, List<List<Solution>> fronts)
         {
             // Find the ideal point
             double[] idealPoint = new double[numberOfObjectives];
 
-            foreach (var solution in st)
+            foreach (var solution in solutions)
                 solution.ConvertedFitness = new double[numberOfObjectives];
 
@@ -72,13 +123,16 @@
             {
                 // Compute ideal point
-                idealPoint[j] = Utility.Min(s => s.Fitness[j], st);
+                // 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
-                foreach (var solution in st)
+                foreach (var solution in solutions)
                     solution.ConvertedFitness[j] = solution.Fitness[j] - idealPoint[j];
             }
 
             // Find the extreme points
-            Solution[] extremePoints = new Solution[numberOfObjectives];
+            List<Solution> extremePoints = new List<Solution>(numberOfObjectives);
             for (int j = 0; j < numberOfObjectives; j++)
             {
@@ -88,27 +142,27 @@
                 weights[j] = 1;
 
-                extremePoints[j] = Utility.ArgMin(s => ASF(s.ConvertedFitness, weights), st);
+                var extremePoint = Utility.ArgMin(s => ASF(s.ConvertedFitness, weights), solutions);
+                extremePoints.Add(extremePoint);
             }
 
             // Compute intercepts
-            List<double> intercepts = GetIntercepts(extremePoints.ToList());
+            double[] intercepts = GetIntercepts(extremePoints, solutions);
 
             // Normalize objectives
-            NormalizeObjectives(st, intercepts, idealPoint);
-        }
-
-        private static void NormalizeObjectives(List<Solution> st, List<double> intercepts, double[] idealPoint)
-        {
-            int numberOfObjectives = st.First().Fitness.Length;
-
-            foreach (var solution in st)
+            NormalizeObjectives(solutions, intercepts, idealPoint);
+        }
+
+        private void NormalizeObjectives(List<Solution> solutions, double[] intercepts, double[] idealPoint)
+        {
+            foreach (var solution in solutions)
             {
                 for (int i = 0; i < numberOfObjectives; i++)
                 {
-                    if (Math.Abs(intercepts[i] - idealPoint[i]) > EPSILON)
-                    {
-                        solution.ConvertedFitness[i] = solution.ConvertedFitness[i] / (intercepts[i] - idealPoint[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];
+                        solution.ConvertedFitness[i] = solution.ConvertedFitness[i] / intercepts[i];
                     }
                     else
@@ -120,5 +174,5 @@
         }
 
-        private static void Associate(List<List<Solution>> fronts, List<ReferencePoint> referencePoints)
+        private void Associate(List<List<Solution>> fronts, List<ReferencePoint> referencePoints)
         {
             for (int f = 0; f < fronts.Count; f++)
@@ -148,5 +202,5 @@
         }
 
-        private static List<Solution> Niching(int k, List<ReferencePoint> referencePoints, IRandom random)
+        private List<Solution> Niching(int k, List<ReferencePoint> referencePoints, IRandom random)
         {
             List<Solution> solutions = new List<Solution>();
@@ -171,8 +225,13 @@
         }
 
-        private static ReferencePoint FindNicheReferencePoint(List<ReferencePoint> referencePoints, IRandom random)
+        private ReferencePoint FindNicheReferencePoint(List<ReferencePoint> referencePoints, IRandom random)
         {
             // the minimum number of associated solutions for a reference point over all reference points
             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
@@ -189,5 +248,5 @@
         }
 
-        private static Solution SelectClusterMember(ReferencePoint referencePoint)
+        private Solution SelectClusterMember(ReferencePoint referencePoint)
         {
             Solution chosen = null;
@@ -202,5 +261,5 @@
         }
 
-        private static double GetPerpendicularDistance(double[] values, double[] fitness)
+        private double GetPerpendicularDistance(double[] values, double[] fitness)
         {
             double numerator = 0;
@@ -221,7 +280,6 @@
         }
 
-        private static double ASF(double[] x, double[] weight)
-        {
-            int numberOfObjectives = x.Length;
+        private double ASF(double[] x, double[] weight)
+        {
             List<int> dimensions = new List<int>();
             for (int i = 0; i < numberOfObjectives; i++)
@@ -231,8 +289,6 @@
         }
 
-        private static List<double> GetIntercepts(List<Solution> extremePoints)
-        {
-            int numberOfObjectives = extremePoints.First().Fitness.Length;
-
+        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.
@@ -246,40 +302,64 @@
             }
 
-            List<double> intercepts = new List<double>();
-
-            // todo: do handling of this case as in Tsung Che Chiang's implementation
-            if (duplicate)
-            { // cannot construct the unique hyperplane (this is a casual method to deal with the condition)
+            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++)
                 {
-                    // extreme_points[f] stands for the individual with the largest value of
-                    // objective f
-                    intercepts.Add(extremePoints[f].ConvertedFitness[f]);
-                }
-            }
-            else
-            {
-                // Find the equation of the hyperplane
-                List<double> b = new List<double>(); //(pop[0].objs().size(), 1.0);
-                for (int i = 0; i < numberOfObjectives; i++)
-                {
-                    b.Add(1.0);
-                }
-
-                List<List<double>> a = extremePoints.Select(p => p.ConvertedFitness.ToList()).ToList();
-
-                List<double> x = GaussianElimination(a, b);
-
-                // Find intercepts
-                for (int f = 0; f < numberOfObjectives; f++)
-                {
-                    intercepts.Add(1.0 / x[f]);
-                }
-            }
-
-            return intercepts;
-        }
-
-        private static List<double> GaussianElimination(List<List<double>> a, List<double> b)
+                    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>();
Index: /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/ReferencePoint.cs
===================================================================
--- /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/ReferencePoint.cs	(revision 17723)
+++ /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/ReferencePoint.cs	(revision 17724)
@@ -25,6 +25,4 @@
         public int NumberOfAssociatedSolutions { get; set; } = 0;
 
-        public int Objectives => Values.Length;
-
         #endregion Properties
 
@@ -32,7 +30,5 @@
 
         [StorableConstructor]
-        public ReferencePoint(StorableConstructorFlag _)
-        {
-        }
+        public ReferencePoint(StorableConstructorFlag _) { }
 
         public ReferencePoint(IRandom random, int obj)
@@ -45,6 +41,13 @@
         {
             random = cloner.Clone(other.random);
-            Values = new double[other.Objectives];
-            other.Values.CopyTo(Values, 0);
+            Values = new double[other.Values.Length];
+            Array.Copy(other.Values, Values, Values.Length);
+        }
+
+        public ReferencePoint(ReferencePoint other)
+        {
+            random = other.random;
+            Values = new double[other.Values.Length];
+            Array.Copy(other.Values, Values, Values.Length);
         }
 
Index: /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/Solution.cs
===================================================================
--- /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/Solution.cs	(revision 17723)
+++ /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/Solution.cs	(revision 17724)
@@ -22,7 +22,5 @@
 
         [StorableConstructor]
-        public Solution(StorableConstructorFlag _)
-        {
-        }
+        public Solution(StorableConstructorFlag _) { }
 
         public Solution(RealVector chromosome)
Index: /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/Utility.cs
===================================================================
--- /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/Utility.cs	(revision 17723)
+++ /branches/2825-NSGA3/HeuristicLab.Algorithms.NSGA3/3.3/Utility.cs	(revision 17724)
@@ -114,5 +114,5 @@
         {
             if (referencePoints == null || !referencePoints.Any()) throw new ArgumentException($"{nameof(referencePoints)} is null or empty");
-            int obj = referencePoints.First().Objectives;
+            int obj = referencePoints.First().Values.Length;
             int pointCount = referencePoints.Count;
             double[,] array = new double[pointCount, obj];
