# Pareto set for 2D functions

2 views (last 30 days)
laura bagnale on 20 May 2021
Commented: laura bagnale on 26 May 2021
Hello everyone!
Can someone help me with this problem please?
I represented two functions like these:
fcontour(@(x,y) 0.9097 + (-0.2332*x) + 0.2645*y + (-1.07)*x^.2 + (-0.02072)*x*y + 0.5786*y^.2;
fcontour(@(x,y) 0.8906 -9.46e+12*x +9.46e12*y + 5.20e+13*x.^2 + 0.03807*y*x + 5.208e+13*y.^2;
Then I wrote an optimization problem in order to find their minimum.
Now I would like to find the Pareto set.
I was trying with the following code:
fitness = @simpleMultiObjective;
nvar = 2;
[x, fval] = gamultiobj(fitness, nvar);
function z = simpleMultiObjective(xy)
x = xy(1); y = xy(2);
z(1) = 0.9097 + (-0.2332*x) + 0.2645*y + (-1.07)*x^.2 + (-0.02072)*x*y + 0.5786*y^.2;
z(2) = 0.8906 -9.46e+12*x +9.46e12*y + 5.20e+13*x.^2 + 0.03807*y*x + 5.208e+13*y.^2;
end
But I don't know if I'm on the right way.
Can someone suggest a way to solve this problem, please?
Thank you very much!

Andreas Apostolatos on 20 May 2021
Hi Laura,
You are on the right path, but I would have some recommendations for you accordingly:
1.) When using function 'fcontour' to plot the contours of a function, MATLAB tries to provide the inputs (x,y) of the function as vectors. For this to be possible, the underlying expressions within the function definition should support vector operations. In your definitions above, you used expression 'x.^2' which supports the base to be a vector, thus raising each component of the vector 'x' to 2 while returning a vector containing the results. However, expression 'x*y' cannot be applied on two vectors of the same form (both row or both column vectors). To allow operation 'x*y' to be applied on two vectors elementwise, please consider adding a full stop '.' right before the corresponding operation, namely 'x.*y'. This will render the underlying computations more efficient. In your case, you can change your function definitions as follows,
fcontour(@(x,y) 0.9097 + (-0.2332.*x) + 0.2645.*y + (-1.07).*x.^.2 + (-0.02072).*x.*y + 0.5786.*y.^.2);
fcontour(@(x,y) 0.8906 - 9.46e+12.*x +9.46e12.*y + 5.20e+13.*x.^2 + 0.03807.*y.*x + 5.208e+13.*y.^2);
2.) Your objective function is unconstrained and has expressions such as 'x^.2'. The latter expression is in fact equal to 'x^0.2', which leads to imaginary results if 'x' is negative. Therefore, it may be advisable to define 0 as lower bound to your design variables such that the objective function does not attain imaginary values, of course depending on your application and your needs. I have adapted your script to account also for 0 as lower bound, namely,
%% clear output window and workspace
clc, clear;
%% plot the contour of the objective functions
fcontour(@(x, y) 0.9097 + (-0.2332.*x) + 0.2645.*y + (-1.07).*x.^.2 + (-0.02072).*x.*y + 0.5786.*y.^.2);
fcontour(@(x, y) 0.8906 - 9.46e+12.*x +9.46e12.*y + 5.20e+13.*x.^2 + 0.03807.*y.*x + 5.208e+13.*y.^2);
%% define optimization options and solve the Pareto front problem
% opts = optimoptions('gamultiobj')';
opts = optimoptions('gamultiobj', 'PopulationSize', 61, 'ParetoFraction', 0.71);
fitness = @simpleMultiObjective;
nvar = 2;
lb = zeros(2, 1);
[x, fval] = gamultiobj(fitness, nvar, [], [], [], [], lb, [], [], opts);
figure
plot(fval(:, 1), fval(:, 2),'o')
%% define the multiobjective function
function z = simpleMultiObjective(xy)
x = xy(1); y = xy(2);
z(1) = 0.9097 + (-0.2332.*x) + 0.2645.*y + (-1.07).*x.^.2 + (-0.02072).*x.*y + 0.5786.*y.^.2;
z(2) = 0.8906 - 9.46e+12.*x +9.46e12.*y + 5.20e+13.*x.^2 + 0.03807.*y.*x + 5.208e+13.*y.^2;
end
Please note that I have increased the values for the 'PopulationSize' and 'ParetoFraction' optimization options when using the 'gamultiobj' function in order to obtain a denser Pareto front. For more information, please visit also the following documentation page,
In this way, the following Pareto front is produced when running the latter script, I hope this helps.
Kind regards,
Andreas
Andreas Apostolatos on 21 May 2021
Hi Laura,
Function 'gamultiobj' utilizes the genetic algorithm for finding the Pareto front which does not employ any starting point. Therefore, the optimization options of function 'gamultiobj' do not provide an API for supplying the algorithm with a starting point.
However, function 'paretosearch' uses the pattern search algorithm for finding the Pareto front, which provides an API for supplying a starting point to the algorithm, see the following documentation pages for more information to the subject,
The corresponding optimization option is called 'InitialPoints' and can accommodate a matrix with 'nvars' columns, where each row represents one initial point while 'nvars' stands for the number of variables. You can then supply the corresponding minima of each objective function that comprise your multiobjective function to the latter optimization option of the 'paretosearch' function, and in this way, you can also skip the specification of lower bounds, since providing the minima of the individual objective functions as starting point suffices for finding the desirable Pareto front.
You can find the adaptation of my previous code snippet that utilizes function 'paretosearch' and optimization option 'InitialPoints' for your convenience below,
%% clear output window and workspace
clc, clear;
%% plot the contour of the objective functions
fcontour(@(x,y) 0.9097 + (-0.2332.*x) + 0.2645.*y + (-1.07).*x.^.2 + (-0.02072).*x.*y + 0.5786.*y.^.2);
fcontour(@(x,y) 0.8906 - 9.46e+12.*x +9.46e12.*y + 5.20e+13.*x.^2 + 0.03807.*y.*x + 5.208e+13.*y.^2);
%% find the minimum of objective function 1
p1 = optimproblem;
p1.Description = "minimization";
x1 = optimvar("x", "LowerBound", 0.5, "UpperBound",2);
y1 = optimvar("y", "LowerBound", 0.5, "UpperBound",2);
p1.Objective = 0.9097 + (-0.2332.*x1) + 0.2645.*y1 + (-1.07).*x1.^.2 + (-0.02072).*x1.*y1 + 0.5786.*y1.^.2;
initialPt1.x = 0.3;
initialPt1.y = 0.5;
[sol1, fval1, exitflag1, outout1] = solve(p1, initialPt1);
%% find the minimum of objective function 2
p2 = optimproblem;
p2.Description = "minimization";
x2 = optimvar("x", "LowerBound", 0.5, "UpperBound",2);
y2 = optimvar("y", "LowerBound", 0.5, "UpperBound",2);
p2.Objective = 0.8906 - 9.46e+12.*x2 +9.46e12.*y2 + 5.20e+13.*x2.^2 + 0.03807.*y2.*x2 + 5.208e+13.*y2.^2;
initialPt2.x = 0.3;
initialPt2.y = 0.5;
[sol2, fval2, exitflag2, outout2] = solve(p2, initialPt2);
%% find the Pareto front using the paretosearch algorithm
x0(1, :) = [sol1.x sol1.y];
x0(2, :) = [sol2.x sol2.y];
opts = optimoptions('paretosearch', 'InitialPoints', x0);
fitness = @simpleMultiObjective;
nvar = 2;
[x, fval] = paretosearch(fitness, nvar,[],[],[],[],[],[],[],opts);
figure
plot(fval(:,1), fval(:,2),'o')
xlabel('Objective 1')
ylabel('Objective 2')
title('Pareto front')
%% define the multiobjective function
function z = simpleMultiObjective(xy)
x = xy(1); y = xy(2);
z(1) = 0.9097 + (-0.2332.*x) + 0.2645.*y + (-1.07).*x.^.2 + (-0.02072).*x.*y + 0.5786.*y.^.2;
z(2) = 0.8906 - 9.46e+12.*x +9.46e12.*y + 5.20e+13.*x.^2 + 0.03807.*y.*x + 5.208e+13.*y.^2;
end
There is also a short video that introduces the latter topic, see the link below,
I hope that this information helps you further.
Kind regards,
Andreas

laura bagnale on 24 May 2021
Hi Andreas,
it is perfect for me, thank you very much for your support!
I really want to thank you a lot!
Kind Regards,
Laura!
laura bagnale on 26 May 2021
Dear Andreas,
sorry if I am disturbing you. It seems that I solved the problem by replacing z(1) and z(2) with corrected form in the syntax. In particular I replaced x.^2 with x.^.2 and for y.^2 too. But I don't understand what does it mean exactly.
I try to looking for it in the matlab documentation but if you could clarify this, you would be doing me a great favour.
If you can't, no problem I will find it! :)
Thank you very much!
Laura