compare paretosearch and gamultiobj -凯发k8网页登录
this example shows how to create a set of points on the pareto front using both paretosearch
and gamultiobj
. the objective function has two objectives and a two-dimensional control variable x
. the objective function mymulti3
is available in your matlab® session when you click the button to edit or try this example. alternatively, copy the mymulti3
code to your session. for speed of calculation, the function is vectorized.
type mymulti3
function f = mymulti3(x) % f(:,1) = x(:,1).^4 x(:,2).^4 x(:,1).*x(:,2) - (x(:,1).^2).*(x(:,2).^2) - 9*x(:,1).^2; f(:,2) = x(:,2).^4 x(:,1).^4 x(:,1).*x(:,2) - (x(:,1).^2).*(x(:,2).^2) 3*x(:,2).^3;
basic example and plots
find pareto sets for the objective functions using paretosearch
and gamultiobj
. set the usevectorized
option to true
for added speed. include a plot function to visualize the pareto set.
rng default nvars = 2; opts = optimoptions(@gamultiobj,'usevectorized',true,'plotfcn','gaplotpareto'); [xga,fvalga,~,gaoutput] = gamultiobj(@(x)mymulti3(x),nvars,[],[],[],[],[],[],[],opts);
optimization terminated: average change in the spread of pareto solutions less than options.functiontolerance.
optsp = optimoptions('paretosearch','usevectorized',true,'plotfcn',{'psplotparetof' 'psplotparetox'}); [xp,fvalp,~,psoutput] = paretosearch(@(x)mymulti3(x),nvars,[],[],[],[],[],[],[],optsp);
pareto set found that satisfies the constraints. optimization completed because the relative change in the volume of the pareto set is less than 'options.paretosetchangetolerance' and constraints are satisfied to within 'options.constrainttolerance'.
compute theoretically exact points on the pareto front by using mymulti4
. the mymulti4
function is available in your matlab session when you click the button to edit or try this example.
type mymulti4
function mout = mymulti4(x) % gg = [4*x(1)^3 x(2)-2*x(1)*(x(2)^2) - 18*x(1); x(1) 4*x(2)^3-2*(x(1)^2)*x(2)]; gf = gg [18*x(1);9*x(2)^2]; mout = gf(1)*gg(2) - gf(2)*gg(1);
the mymulti4
function evaluates the gradients of the two objective functions. next, for a range of values of x(2)
, use fzero
to locate the point where the gradients are exactly parallel, which is where the output mout
= 0.
a = [fzero(@(t)mymulti4([t,-3.15]),[2,3]),-3.15]; for jj = linspace(-3.125,-1.89,50) a = [a;[fzero(@(t)mymulti4([t,jj]),[2,3]),jj]]; end figure plot(fvalp(:,1),fvalp(:,2),'bo'); hold on fs = mymulti3(a); plot(fvalga(:,1),fvalga(:,2),'r*'); plot(fs(:,1),fs(:,2),'k.') legend('paretosearch','gamultiobj','true') xlabel('objective 1') ylabel('objective 2') hold off
gamultiobj
finds points with a slightly wider spread in objective function space. plot the solutions in decision variable space, along with the theoretical optimal pareto curve and a contour plot of the two objective functions.
[x,y] = meshgrid(1.9:.01:3.1,-3.2:.01:-1.8); mydata = mymulti3([x(:),y(:)]); myff = sqrt(mydata(:,1) 39);% spaces the contours better mygg = sqrt(mydata(:,2) 28);% spaces the contours better myff = reshape(myff,size(x)); mygg = reshape(mygg,size(x)); figure; hold on contour(x,y,mygg,50) contour(x,y,myff,50) plot(xp(:,1),xp(:,2),'bo') plot(xga(:,1),xga(:,2),'r*') plot(a(:,1),a(:,2),'-k') xlabel('x(1)') ylabel('x(2)') hold off
unlike the paretosearch
solution, the gamultiobj
solution has points at the extreme ends of the range in objective function space. however, the paretosearch
solution has more points that are closer to the true solution in both objective function space and decision variable space. the number of points on the pareto front is different for each solver when you use the default options.
shifted problem
what happens if the solution to your problem has control variables that are large? examine this case by shifting the problem variables. for an unconstrained problem, gamultiobj
can fail, while paretosearch
is more robust to such shifts.
for easier comparison, specify 35 points on the pareto front for each solver.
shift = [20,-30];
fun = @(x)mymulti3(x shift);
opts.populationsize = 100; % opts.paretofraction = 35
[xgash,fvalgash,~,gashoutput] = gamultiobj(fun,nvars,[],[],[],[],[],[],opts);
optimization terminated: average change in the spread of pareto solutions less than options.functiontolerance.
gamultiobj
fails to find a useful pareto set.
optsp.plotfcn = []; optsp.paretosetsize = 35; [xpsh,fvalpsh,~,pshoutput] = paretosearch(fun,nvars,[],[],[],[],[],[],[],optsp);
pareto set found that satisfies the constraints. optimization completed because the relative change in the volume of the pareto set is less than 'options.paretosetchangetolerance' and constraints are satisfied to within 'options.constrainttolerance'.
figure plot(fvalpsh(:,1),fvalpsh(:,2),'bo'); hold on plot(fs(:,1),fs(:,2),'k.') legend('paretosearch','true') xlabel('objective 1') ylabel('objective 2') hold off
paretosearch
finds solution points spread evenly over nearly the entire possible range.
adding bounds, even fairly loose ones, helps both gamultiobj
and paretosearch
to find appropriate solutions. set lower bounds of -50
in each component, and upper bounds of 50
.
opts.plotfcn = []; optsp.plotfcn = []; lb = [-50,-50]; ub = -lb; [xgash,fvalgash,~,gashoutput] = gamultiobj(fun,nvars,[],[],[],[],lb,ub,opts);
optimization terminated: average change in the spread of pareto solutions less than options.functiontolerance.
[xpsh2,fvalpsh2,~,pshoutput2] = paretosearch(fun,nvars,[],[],[],[],lb,ub,[],optsp);
pareto set found that satisfies the constraints. optimization completed because the relative change in the volume of the pareto set is less than 'options.paretosetchangetolerance' and constraints are satisfied to within 'options.constrainttolerance'.
figure plot(fvalpsh2(:,1),fvalpsh2(:,2),'bo'); hold on plot(fvalgash(:,1),fvalgash(:,2),'r*'); plot(fs(:,1),fs(:,2),'k.') legend('paretosearch','gamultiobj','true') xlabel('objective 1') ylabel('objective 2') hold off
in this case, both solvers find good solutions.
start paretosearch
from gamultiobj
solution
obtain a similar range of solutions from the solvers by starting paretosearch
from the gamultiobj
solution.
optsp.initialpoints = xgash; [xpsh3,fvalpsh3,~,pshoutput3] = paretosearch(fun,nvars,[],[],[],[],lb,ub,[],optsp);
pareto set found that satisfies the constraints. optimization completed because the relative change in the volume of the pareto set is less than 'options.paretosetchangetolerance' and constraints are satisfied to within 'options.constrainttolerance'.
figure plot(fvalpsh3(:,1),fvalpsh3(:,2),'bo'); hold on plot(fvalgash(:,1),fvalgash(:,2),'r*'); plot(fs(:,1),fs(:,2),'k.') legend('paretosearch','gamultiobj','true') xlabel('objective 1') ylabel('objective 2') hold off
now the paretosearch
solution is similar to the gamultiobj
solution, although some of the solution points differ.
start from single-objective solutions
another way of obtaining a good solution is to start from the points that minimize each objective function separately.
from the multiobjective function, create a single-objective function that chooses each objective in turn. use the shifted function from the previous section. because you are giving good start points to the solvers, you do not need to specify bounds.
nobj = 2; % number of objectives x0 = -shift; % initial point for single-objective minimization uncmin = cell(nobj,1); % cell array to hold the single-objective minima allfuns = zeros(nobj,2); % hold the objective function values eflag = zeros(nobj,1); fopts = optimoptions('patternsearch','display','off'); % use an appropriate solver here for i = 1:nobj indi = zeros(nobj,1); % choose the objective to minimize indi(i) = 1; funi = @(x)dot(fun(x),indi); [uncmin{i},~,eflag(i)] = patternsearch(funi,x0,[],[],[],[],[],[],[],fopts); % minimize objective i allfuns(i,:) = fun(uncmin{i}); end uncmin = cell2mat(uncmin); % matrix of start points
start paretosearch
from the single-objective minimum points and note that it has a full range in its solutions. paretosearch
adds random initial points to the supplied ones in order to have a population of at least options.paretosetsize
individuals. similarly, gamultiobj
adds random points to the supplied ones to obtain a population of at least (options.populationsize)*(options.paretofraction)
individuals.
optsp = optimoptions(optsp,'initialpoints',uncmin);
[xpinit,fvalpinit,~,outputpinit] = paretosearch(fun,nvars,[],[],[],[],[],[],[],optsp);
pareto set found that satisfies the constraints. optimization completed because the relative change in the volume of the pareto set is less than 'options.paretosetchangetolerance' and constraints are satisfied to within 'options.constrainttolerance'.
now solve the problem using gamultiobj
starting from the initial points.
opts = optimoptions(opts,'initialpopulationmatrix',uncmin);
[xgash2,fvalgash2,~,gashoutput2] = gamultiobj(fun,nvars,[],[],[],[],[],[],opts);
optimization terminated: average change in the spread of pareto solutions less than options.functiontolerance.
figure plot(fvalpinit(:,1),fvalpinit(:,2),'bo'); hold on plot(fvalgash2(:,1),fvalgash2(:,2),'r*'); plot(fs(:,1),fs(:,2),'k.') plot(allfuns(:,1),allfuns(:,2),'gs','markersize',12) legend('paretosearch','gamultiobj','true','start points') xlabel('objective 1') ylabel('objective 2') hold off
both solvers fill in the pareto front between the extreme points, with reasonably accurate and well-spaced solutions.
view the final points in decision variable space.
figure; hold on xx = x - shift(1); yy = y - shift(2); contour(xx,yy,mygg,50) contour(xx,yy,myff,50) plot(xpinit(:,1),xpinit(:,2),'bo') plot(xgash2(:,1),xgash2(:,2),'r*') ashift = a - shift; plot(ashift(:,1),ashift(:,2),'-k') plot(uncmin(:,1),uncmin(:,2),'gs','markersize',12); xlabel('x(1)') ylabel('x(2)') hold off
see also
|