【建模算法】零基础学优化—多目标规划

硬件电子工程师 · 发表于 2023-9-6 00:01:00

点击上方蓝字和“好玩的MATLAB”一起快乐玩耍吧！

ib4ivre4n2w6405197350.jpg

好玩的matlab
带你解锁不一样的matlab新玩法

喜欢此推文的小伙伴们记得点赞+关注+分享！
01
算法用途

ejnjna1yms46405197450.png

优化算法是指在满足一定条件下,在众多方案中或者参数中最优方案,或者参数值,以使得某个或者多个功能指标达到最优,或使得系统的某些性能指标达到最大值或者最小值。
多目标规划的概念是 1961年由美国数学家查尔斯和库柏首先提出的。多目标最优化思想，最早是在1896年由法国经济学家V.帕雷托提出来的。他从政治经济学的角度考虑把本质上是不可比较的许多目标化成单个目标的最优化问题，从而涉及了多目标规划问题和多目标的概念。

02
实例分析

多目标规划用得比较多的主要有以下几种形式
第一种是通过加权的方式将多目标优化转化成单目标优化，不过加权权重比较难确定

lsjgn2loyaq6405197550.png

第二种是goal attain 的形式，将原来的多目标转化成新的优化问题，和单个目标的最优化进行比较，目标是与原来单个最优化目标最小

q2qw0l4mlay6405197650.png

最后就是一些多目标优化的智能优化算法，智能优化算法能在一次运行过程中找到Pareto最优解集中的多个解，且不限于Parato前沿解集的形状和连续性，易于解决不连续的，凹的Pareto前沿。智能优化算法解决多目标优化问题等到下一期出奥

03
MATLAB求解

案例一加权方式求解

fmincon是MATLAB的非线性规划求解函数
[x,fval]=fmincon(fun,x0,A,b,Aeq,Beq,LB,UB,nonlcon)
x：求得最优情况下变量的解
fval：求得最优目标值
x0:初始解
fun：目标函数（符号按最小值标准，若目标是求解机大值可以通过添加负号改成求极小值）
A：不等式约束的变量系数（符合按小于标准，如果是大于约束可通过加负号变成小于）
b：不等式约束的常量
Aeq：等式约束的变量系数Beq：等式约束的常量LB：变量的下限UB：变量的上限 nonlcon :非线性约束函数表达式

优化问题

pffnmefeycx6405197750.png

目标函数

nz4vgzsxo2j6405197850.jpg

求极大值的问题要转化成极小值，因此加个负号变成求最小值的问题。通过附权将两种相加，权重可以自定可以根据因变量的范围而定不一定要限制在一定的范围

w1=0.5;w2=0.5; %fun = @(x)(-x(1)^2+x(2)^2-x(2)*x(3))*w1+(2*x(1)^2-x(2)^3+2*x(2)*x(3))*w2;
不等式线性约束

%不等约束A=[2,1,3];%左边特征矩阵b=[6]; %右边
等式线性约束

无

Aeq=[];Beq=[];
变量约束

%变量约束,上限，下限LB=zeros(3,1);UB=ones(3,1);
非线性约束

nonlcon = @unitdisk;function [c,ceq] = unitdisk(x)%c为不等式非线性约束%ceq为等式非线性约束c=x(1)^2+x(1)*x(2)+x(2)*x(3)-x(2)-6;ceq = [];end

全部代码

clc;clear;close all;%%%初始解，随意给个初始解x0=zeros(3,1);
%不等约束A=[2,1,3];%左边特征矩阵b=[6]; %右边
%等式约束Aeq=[];Beq=[];
%变量约束,上限，下限LB=zeros(3,1);UB=1*ones(3,1);%优化求解%%w1=0.5;w2=0.5; %fun = @(x)(-x(1)^2+x(2)^2-x(2)*x(3))*w1+(2*x(1)^2-x(2)^3+2*x(2)*x(3))*w2;
nonlcon = @unitdisk;[x1,fval1]=fmincon(fun,x0,A,b,Aeq,Beq,LB,UB,nonlcon);objstr=['目标函数最优值：',num2str(fval1)];disp(objstr)for i=1:length(x1) xstr=['x',num2str(i),'的值为：',num2str(x1(i))]; disp(xstr)end
function [c,ceq] = unitdisk(x)%c为不等式非线性约束%ceq为等式非线性约束c=x(1)^2+x(1)*x(2)+x(2)*x(3)-x(2)-6;ceq = [];end

zjrnlu20b256405197950.png

goal attain方式求解1

5af5bnu5m2s6405198051.png

优化目标
例如有以下两个目标，两个目标的非线性约束有一点差别，其他约束都一样

kzg2mfu5g3j6405198151.png

分别得到单个目标最优解做了两次单独的非线性优化，和前面的非线性优化篇方法是一样的，做好初始解，做好线性约束和非线性约束，写好目标带入函数求解即可

%初始解，随意给个初始解x0=zeros(3,1);
%不等约束A=[2,1,3];%左边特征矩阵b=[6]; %右边
%等式约束Aeq=[];Beq=[];
%变量约束,上限，下限LB=zeros(3,1);UB=1*ones(3,1);%优化求解%%fun1 = @(x)-x(1)^2+x(2)^2-x(2)*x(3);fun2 = @(x)2*x(1)^2-x(2)^3+2*x(2)*x(3);%%nonlcon = @unitdisk;[x1,fval1]=fmincon(fun1,x0,A,b,Aeq,Beq,LB,UB,nonlcon);objstr=['目标函数最优值：',num2str(fval1)];disp(objstr)for i=1:length(x1) xstr=['x',num2str(i),'的值为：',num2str(x1(i))]; disp(xstr)endnonlcon1 = @unitdisk1;[x2,fval2]=fmincon(fun2,x0,A,b,Aeq,Beq,LB,UB,nonlcon1);objstr=['目标函数最优值：',num2str(fval2)];disp(objstr)for i=1:length(x2) xstr=['x',num2str(i),'的值为：',num2str(x2(i))]; disp(xstr)end
function [c,ceq] = unitdisk(x)%c为不等式非线性约束%ceq为等式非线性约束c=[x(1)^2+x(1)*x(2)+x(2)*x(3)-x(2)-6];ceq = [];endfunction [c,ceq] = unitdisk1(x)%c为不等式非线性约束%ceq为等式非线性约束c=x(1)^2+x(1)*x(2)+x(2)*x(3)-6;ceq = [];end 可以求解得到两个单目标函数的最优分别如下，他们取得单个最优的x值不同，因此两者结合肯定会相互牺牲

目标函数最优值：-1.25
x1的值为：1
x2的值为：0.5
x3的值为：1

目标函数最优值：-1
x1的值为：0.00056595
x2的值为：1
x3的值为：6.4e-07

确定两个约束的单个最优

goal=[fval1,fval2];

将两个单目标优化结合
约束条件是两个优化目标均需要满足的，因此可以把所有的约束条件直接融合在一起

015uonvlzsj6405198251.png

目标函数 两个目标函数都比较简单，可以用数学表达式表达，因此目标函数可以写成以下的形式，将目标函数分号隔开

func = @(x)[-x(1)^2+x(2)^2-x(2)*x(3);2*x(1)^2-x(2)^3+2*x(2)*x(3)];

不等约束
直接用分号将两个非线性约束结合

function [c,ceq] = unitdisk3(x)%c为不等式非线性约束%ceq为等式非线性约束c=[x(1)^2+x(1)*x(2)+x(2)*x(3)-6;x(1)^2+x(1)*x(2)+x(2)*x(3)-6];ceq = [];end

权重
离两个分别最优的函数差值的权重，可以用于调节两个优化目标量纲不对等的情况

weight=[1,1];

完整MATLAB代码

clc;clear;close all;%%%初始解，随意给个初始解x0=zeros(3,1);
%不等约束A=[2,1,3];%左边特征矩阵b=[6]; %右边
%等式约束Aeq=[];Beq=[];
%变量约束,上限，下限LB=zeros(3,1);UB=1*ones(3,1);%优化求解%%fun1 = @(x)-x(1)^2+x(2)^2-x(2)*x(3);fun2 = @(x)2*x(1)^2-x(2)^3+2*x(2)*x(3);%%nonlcon = @unitdisk;[x1,fval1]=fmincon(fun1,x0,A,b,Aeq,Beq,LB,UB,nonlcon);objstr=['目标函数最优值：',num2str(fval1)];disp(objstr)for i=1:length(x1) xstr=['x',num2str(i),'的值为：',num2str(x1(i))]; disp(xstr)endnonlcon1 = @unitdisk1;[x2,fval2]=fmincon(fun2,x0,A,b,Aeq,Beq,LB,UB,nonlcon1);objstr=['目标函数最优值：',num2str(fval2)];disp(objstr)for i=1:length(x2) xstr=['x',num2str(i),'的值为：',num2str(x2(i))]; disp(xstr)end%% 多目标规划goal=[fval1,fval2];%%nonlcon3=@unitdisk3;func = @(x)[-x(1)^2+x(2)^2-x(2)*x(3);2*x(1)^2-x(2)^3+2*x(2)*x(3)];weight=[1,1];[x,fival]=fgoalattain(func,x0,goal,weight,A,b,Aeq,Beq,LB,UB,nonlcon3);disp('在两个目标的优化结果为')disp(func(x))for i=1:length(x) xstr=['x',num2str(i),'的值为：',num2str(x(i))]; disp(xstr)end%%function [c,ceq] = unitdisk(x)%c为不等式非线性约束%ceq为等式非线性约束c=[x(1)^2+x(1)*x(2)+x(2)*x(3)-x(2)-6];ceq = [];endfunction [c,ceq] = unitdisk1(x)%c为不等式非线性约束%ceq为等式非线性约束c=x(1)^2+x(1)*x(2)+x(2)*x(3)-6;ceq = [];endfunction [c,ceq] = unitdisk3(x)%c为不等式非线性约束%ceq为等式非线性约束c=[x(1)^2+x(1)*x(2)+x(2)*x(3)-6;x(1)^2+x(1)*x(2)+x(2)*x(3)-6];ceq = [];end
可以得到以下的优化结果：
在两个目标的优化结果为
1.0e-15 *

-0.2220
0.4441

x1的值为：1.4901e-08
x2的值为：-3.6734e-40
x3的值为：0

goal attain方式求解2

优化问题
接下来看一个稍微复杂的优化问题，是一个10维的变量，有以下两个目标

h3uzw3isyex6405198351.png

目标函数
整体思路和前面差不多，主要是目标函数变复杂了很多，目标函数不好直接用数学表达式写，需要用函数写，因此把目标做成了函数

单个目标函数

function y1=obj1(x) [dim, num] = size(x); tmp = zeros(dim,num); tmp(2:dim,:)= (x(2:dim,:) - sin(6.0*pi*repmat(x(1,:),[dim-1,1]) + pi/dim*repmat((2:dim)',[1,num]))).^2; tmp1 = sum(tmp(3:2:dim,:)); % odd index tmp2 = sum(tmp(2:2:dim,:)); % even index y1 = x(1,:) + 2.0*tmp1/size(3:2:dim,2);% y(2,:) = 1.0 - sqrt(x(1,:)) + 2.0*tmp2/size(2:2:dim,2);endfunction y2=obj2(x) [dim, num] = size(x); tmp = zeros(dim,num); tmp(2:dim,:)= (x(2:dim,:) - sin(6.0*pi*repmat(x(1,:),[dim-1,1]) + pi/dim*repmat((2:dim)',[1,num]))).^2; tmp1 = sum(tmp(3:2:dim,:)); % odd index tmp2 = sum(tmp(2:2:dim,:)); % even index% y2 = x(1,:) + 2.0*tmp1/size(3:2:dim,2); y2 = 1.0 - sqrt(x(1,:)) + 2.0*tmp2/size(2:2:dim,2);end

整合两个目标函数
将两个目标函数写在一起，分行写

function y=obj3(x) [dim, num] = size(x); tmp = zeros(dim,num); tmp(2:dim,:)= (x(2:dim,:) - sin(6.0*pi*repmat(x(1,:),[dim-1,1]) + pi/dim*repmat((2:dim)',[1,num]))).^2; tmp1 = sum(tmp(3:2:dim,:)); % odd index tmp2 = sum(tmp(2:2:dim,:)); % even index y(1,:) = x(1,:) + 2.0*tmp1/size(3:2:dim,2); y(2,:) = 1.0 - sqrt(x(1,:)) + 2.0*tmp2/size(2:2:dim,2);end
约束条件基本没什么，设置为空就行

全部求解代码

%% 非线性规划1clc;clear;close all;%%%初始解，随意给个初始解x0=zeros(10,1);%%%不等约束A=[];%左边特征矩阵b=[]; %右边%等式约束Aeq=[];Beq=[];%%%变量约束,上限，下限LB=-1*ones(10,1);LB(1)=0;UB=1*ones(10,1);%优化求解%%fun1 = @obj1;fun2 = @obj2;%%nonlcon = [];[x1,fval1]=fmincon(fun1,x0,A,b,Aeq,Beq,LB,UB,nonlcon);objstr=['目标函数最优值：',num2str(fval1)];disp(objstr)for i=1:length(x1) xstr=['x',num2str(i),'的值为：',num2str(x1(i))]; disp(xstr)end[x2,fval2]=fmincon(fun2,x0,A,b,Aeq,Beq,LB,UB,nonlcon);objstr=['目标函数最优值：',num2str(fval2)];disp(objstr)for i=1:length(x2) xstr=['x',num2str(i),'的值为：',num2str(x2(i))]; disp(xstr)end%% 多目标规划goal=[fval1,fval2];func = @obj3;weight=[1,1];[x,fival]=fgoalattain(func,x0,goal,weight,A,b,Aeq,Beq,LB,UB,nonlcon);disp('在两个目标的优化结果为')disp(func(x))for i=1:length(x) xstr=['x',num2str(i),'的值为：',num2str(x(i))]; disp(xstr)end%%function y1=obj1(x) [dim, num] = size(x); tmp = zeros(dim,num); tmp(2:dim,:)= (x(2:dim,:) - sin(6.0*pi*repmat(x(1,:),[dim-1,1]) + pi/dim*repmat((2:dim)',[1,num]))).^2; tmp1 = sum(tmp(3:2:dim,:)); % odd index tmp2 = sum(tmp(2:2:dim,:)); % even index y1 = x(1,:) + 2.0*tmp1/size(3:2:dim,2);% y(2,:) = 1.0 - sqrt(x(1,:)) + 2.0*tmp2/size(2:2:dim,2);endfunction y2=obj2(x) [dim, num] = size(x); tmp = zeros(dim,num); tmp(2:dim,:)= (x(2:dim,:) - sin(6.0*pi*repmat(x(1,:),[dim-1,1]) + pi/dim*repmat((2:dim)',[1,num]))).^2; tmp1 = sum(tmp(3:2:dim,:)); % odd index tmp2 = sum(tmp(2:2:dim,:)); % even index% y2 = x(1,:) + 2.0*tmp1/size(3:2:dim,2); y2 = 1.0 - sqrt(x(1,:)) + 2.0*tmp2/size(2:2:dim,2);endfunction y=obj3(x) [dim, num] = size(x); tmp = zeros(dim,num); tmp(2:dim,:)= (x(2:dim,:) - sin(6.0*pi*repmat(x(1,:),[dim-1,1]) + pi/dim*repmat((2:dim)',[1,num]))).^2; tmp1 = sum(tmp(3:2:dim,:)); % odd index tmp2 = sum(tmp(2:2:dim,:)); % even index y(1,:) = x(1,:) + 2.0*tmp1/size(3:2:dim,2); y(2,:) = 1.0 - sqrt(x(1,:)) + 2.0*tmp2/size(2:2:dim,2);end

求解得到结果如下所示：
在两个目标的优化结果为
0.3837
0.5008

x1的值为：0.3773
x2的值为：1
x3的值为：1
x4的值为：0.45933
x5的值为：0.72692
x6的值为：0.41477
x7的值为：0.13379
x8的值为：0.098991
x9的值为：-0.58879
x10的值为：-0.90934

参考资料：
【1】https://item.jd.com/13422442.html MATLAB智能优化算法：从写代码到算法思想
【2】https://ww2.mathworks.cn/help/optim/ug/intlinprog.html
往期精彩回顾

ihyzhcjvttl6405198451.png

推荐 | 【建模算法】BP神经网络推荐 | 【建模算法】K-means聚类推荐 | 【建模算法】多层复杂BP神经网络推荐 | MATLAB｜聚类散点图｜边际图｜核密度填充图

bwqrk00qgox6405198551.png

an3fpiepk5h6405198651.jpg