Fixed small bug in ND multiobj optimization and added benchmark from …

…paper with known analytical solution

docs/make.jl

@@ -34,7 +34,8 @@ makedocs(
         "Branin function" => "BraninFunction.md",
         "Ackley function" => "ackley.md",
         "Gramacy & Lee Function" => "gramacylee.md",
-        "Salustowicz Benchmark" => "Salustowicz.md"
+        "Salustowicz Benchmark" => "Salustowicz.md",
+        "Multi objective optimization" => "multi_objective_opt.md"
         ]
     "Contributing" => "contributing.md"
     ]

docs/src/multi_objective_opt.md

@@ -0,0 +1,46 @@
+# Multi objective optimization benchmarks
+
+## Case 1: Non colliding objective functions
+
+```@example multi_obj
+using Surrogates
+#EGO
+m = 10
+f  = x -> [x^i for i = 1:m]
+lb = 1.0
+ub = 10.0
+x  = sample(100, lb, ub, SobolSample())
+y  = f.(x)
+my_radial_basis_ego = RadialBasis(x, y, lb, ub, rad = linearRadial)
+pareto_set, pareto_front = surrogate_optimize(f,EGO(),lb,ub,my_radial_basis_ego,SobolSample())
+
+m = 5
+f  = x -> [x^i for i =1:m]
+lb = 1.0
+ub = 10.0
+x  = sample(100, lb, ub, SobolSample())
+y  = f.(x)
+my_radial_basis_rtea = RadialBasis(x, y, lb, ub, rad = linearRadial)
+Z = 0.8
+K = 2
+p_cross = 0.5
+n_c = 1.0
+sigma = 1.5
+surrogate_optimize(f,RTEA(Z,K,p_cross,n_c,sigma),lb,ub,my_radial_basis_rtea,SobolSample())
+```
+
+## Case 2: objective functions with conflicting minima
+
+```@example multi_obj
+#EGO
+f  = x -> [sqrt((x[1] - 4)^2 + 25*(x[2])^2),
+           sqrt((x[1]+4)^2 + 25*(x[2])^2),
+           sqrt((x[1]-3)^2 + (x[2]-1)^2)]
+lb = [2.5,-0.5]
+ub = [3.5,0.5]
+x  = sample(100, lb, ub, SobolSample())
+y  = f.(x)
+my_radial_basis_ego = RadialBasis(x, y, lb, ub, rad = linearRadial)
+#I can find my pareto set and pareto front by calling again the surrogate_optimize function:
+pareto_set, pareto_front = surrogate_optimize(f,EGO(),lb,ub,my_radial_basis_ego,SobolSample(),maxiters=30);
+```

src/Optimization.jl

@@ -1485,15 +1485,15 @@ function surrogate_optimize(obj::Function,ego::EGO,lb,ub,surrEGOND::AbstractSurr
         add_point!(surrEGOND,x_new,y_new)
     end
     #Find and return Pareto
-    y = surrEGO.y
+    y = surrEGOND.y
     y = permutedims(reshape(hcat(y...),(length(y[1]), length(y)))) #2d matrix
     Fronts = _nonDominatedSorting(y) #this returns the indexes
     pareto_front_index = Fronts[1]
     pareto_set = []
     pareto_front = []
     for i = 1:length(pareto_front_index)
-        push!(pareto_set,surrEGO.x[pareto_front_index[i]])
-        push!(pareto_front,surrEGO.y[pareto_front_index[i]])
+        push!(pareto_set,surrEGOND.x[pareto_front_index[i]])
+        push!(pareto_front,surrEGOND.y[pareto_front_index[i]])
     end
     return pareto_set,pareto_front
 end
@@ -1606,29 +1606,27 @@ end
      return pareto_set,pareto_front
  end
 
-
-
  function surrogate_optimize(obj,rtea::RTEA,lb,ub,surrRTEAND::AbstractSurrogate,sample_type::SamplingAlgorithm;maxiters=100, n_new_look = 1000)
      Z = rtea.z
      K = rtea.k
      p_cross = rtea.p
      n_c = rtea.n_c
      sigma = rtea.sigma
      #find pareto set of the first evaluations: (estimated pareto)
-     y = surrRTEA.y
+     y = surrRTEAND.y
      y = permutedims(reshape(hcat(y...),(length(y[1]), length(y)))) #2d matrix
      Fronts = _nonDominatedSorting(y) #this returns the indexes
      pareto_front_index = Fronts[1]
      pareto_set = []
      pareto_front = []
      for i = 1:length(pareto_front_index)
-         push!(pareto_set,surrRTEA.x[pareto_front_index[i]])
-         push!(pareto_front,surrRTEA.y[pareto_front_index[i]])
+         push!(pareto_set,surrRTEAND.x[pareto_front_index[i]])
+         push!(pareto_front,surrRTEAND.y[pareto_front_index[i]])
      end
      number_of_revaluations = zeros(Int,length(pareto_set))
      iter = 1
      d = length(lb)
-     dim_out = length(surrRTEA.y[1])
+     dim_out = length(surrRTEAND.y[1])
      while iter < maxiters
 
          if iter < (1-Z)*maxiters