用matlab求解优化问题.docx
《用matlab求解优化问题.docx》由会员分享，可在线阅读，更多相关《用matlab求解优化问题.docx（38页珍藏版）》请在冰豆网上搜索。
用matlab求解优化问题
§8.1.1线性规划问题的MATLAB求解方法
与一般线性规划理论一样，在MATLAB中有线性规划的标准型。
在调用MATLAB线性规划函数linprog时，要遵循MATLAB中对标准性的要求。
线性规划问题的MATLAB标准形为：
在上述模型中，有一个需要极小化的目标函数
，以及需要满足的约束条件假设
为
维设计变量，且线性规划问题具有不等式约束
个，等式约束
个，那么：
和
均为
维列向量，
为
维列向量，
为m2维列向量，
为
维矩阵，
为
维矩阵
需要注意的是：
MATLAB标准型是对目标函数求极小，如果遇到是对目标函数求极大的问题，在使用MATLAB求解时，需要在函数前面加一个负号转化为对目标函数求极小的问题；
MATLAB标准型中的不等式约束形式为
，如果在线性规划问题中出现
形式的不等式约束，则我们需要在两边乘以（-1）使其转化为MATLAB中的
形式。
如果在线性规划问题中出现了“<”或者“>”的约束形式，则我们需要通过添加松弛变量使得不等式约束变为等式约束之后，我们只需要将所有的约束（包括不等式约束和等式约束）转化为矩阵形式的即可。
例如，对于如下线性规划模型：
要转化为MATLAB标准形，则要经过：
（1）原问题是对目标函数求极大，故添加负号使目标变为：
；
（2）原问题中存在“≥”的约束条件，故添加负号使其变为：
用MATLAB表达则为
c=[-4;2;-1];%将目标函数转化为求极小
A=[2-11;8-22];b=[12;-8];%不等式约束系数矩阵
Aeq=[-201;110];beq=[3;7];%等式约束系数矩阵
lb=[0;0;0];ub=[Inf;Inf;Inf]%对设计变量的边界约束
MATLAB优化工具箱中求解线性规划问题的命令为linprog，其函数调用方法有多种形式如下所示：
x=linprog（c,A,b）
x=linprog（c,A,b,Aeq,beq）
x=linprog（c,A,b,Aeq,beq,lb,ub）
x=linprog（c,A,b,Aeq,beq,lb,ub,x0）
x=linprog（c,A,b,Aeq,beq,lb,ub,x0,options）
x=linprog（problem）
[x,fval]=linprog（...）
[x,fval,exitflag]=linprog（...）
[x,fval,exitflag,output]=linprog（...）
[x,fval,exitflag,output,lambda]=linprog（...）
输入参数
MATLAB工具箱中的linprog函数在求解线性规划问题时，提供的参数为：
模型参数、初始解参数和算法控制参数。
模型参数x、c、lb、ub、b、beq、A和Aeq在MATLAB标准型中已经介绍了其具体物理意义和获得方法
x0为线性规划问题的初始解，该设置仅在中型规模算法中有效，而在默认的大型规模算法和单纯形算法中，MATLAB将忽略一切初始解。
options为包含算法控制参数的结构变量，我们可以通过optimset命令对这些具体的控制参数进行设置，例如下述格式
options=optimset（'param1',value1,'param2',value2,...）
该命令格式创建一组控制参数结构变量options，将参数的具体值赋给单引号之间的参数，任何未被指定的参数将被赋值为[]，参数值为[]的具体的含义是将该组控制参数传递给优化函数时将使用MATLAB提供的默认值
在线性规划问题中可以用到的设置参数如下表所示
表1linprog中的控制参数设置
参数名称
参数设置
Diagnostics
设置是否显示函数优化中的诊断信息，可以选择on或者off（默认值），该功能主要显示一些退出信息，即linprog函数运算终止的原因
Display
设置显示信息的级别，当该参数值为off时，不显示任何输出信息；当参数值为iter时，将显示每一步迭代的输出信息，iter参数值仅对大型规模算法和中型规模的单纯形算法有效；当参数值为final时，仅显示最终的输出信息
Simplex
当该参数值为on时，函数采用单纯形算法
LargeScale
设置是否采用大型规模算法，当参数值为on（默认值）时，使用大型规模算法；当参数值为off时，使用中型规模算法
MaxIter
算法运行中的最大迭代次数，对于大型规模算法，默认值为85，对于单纯形算法，其默认值为10*设计变量的个数，对于中型有效集算法为10*max（设计变量的个数，
不等式约束的个数+边界约束的个数）
TolFun
函数计算终止的误差限，对于大型规模算法其默认值为1e-8，该控制参数对于中型规模的有效集算法无效。
输出参数
linprog函数返回的输出参数有x、fval、exitflag、lambda和output。
输出参数x为线性规划问题的最优解
输出参数fval为线性规划问题在最优解x处的函数值
输出参数exitflag返回的是优化函数计算终止时的状态指示，说明算法终止的原因，其取值和其代表的具体原因如下表所示
表2参数exitflag的取值和物理意义
值
物理意义
1
已经收敛到解x
0
已经达到最大迭代次数限制options.MaxIter
-2
没有找到问题的可行点
-3
问题无有限最优解
-4
在算法执行过程中遇到了NaN值
-5
原线性规划问题和其对偶问题均不可行
-7
搜索方向变化太小，无法进一步获得更优解，说明原线性规划问题或者约束条件是病态的
输出参数output是一个返回优化过程中相关信息的结构变量，其属性如下表所示
表3参数output所包含信息
属性名称
属性含义
output.iterations
优化过程的实际迭代次数
output.algorithm
优化过程中所采用的具体算法
output.cgiterations
0（仅用于大型规模算法，为了后向兼容性而设置的参数）
output.message
退出信息
输出参数lambda是返回线性规划问题最优解x处的拉格朗日乘子的一个结构变量,其总维数等于约束条件的个数，其非零分量对应于起作用的约束条件，其属性如表所示。
表4参数lambda所包含信息
属性名称
属性含义
ineqlin
线性不等式约束的拉格朗日乘子
eqlin
线性等式约束的拉格朗日乘子
upper
上界约束的拉格朗日乘子
lower
下界约束的拉格朗日乘子
命令详解
1）x=linprog（c,A,b）；该函数调用格式求解线性规划问题：
2）x=linprog（c,A,b,Aeq,beq）；该函数调用格式求解线性规划问题：
即该函数调用格式解决的是既含有线性等式约束，又含有线性不等式约束的线性规划问题，如果在线性规划问题中无线性不等式约束，则可以设A=[]以及b=[]
3）x=linprog（c,A,b,Aeq,beq,lb,ub）；该函数调用格式求解线性规划问题：
即在线性规划问题的求解过程中进一步考虑了对设计变量的约束，其中lb和ub均是和设计变量维数相同的列向量，如果对设计变量没有上界约束，可以设置ub（i）=Inf，如果没有下界约束则可以设置lb（i）=-Inf，和
（2）类似，如果问题中没有等式约束，则可以设Aeq=[]以及beq=[]
4）x=linprog（c,A,b,Aeq,beq,lb,ub,x0）
在前面调用方法的基础上设置线性规划问题求解的初始解为x0，该参数仅在使用有效集算法时生效，否则当使用默认的内点算法时，将忽略任何初始点，即参数无效。
5）x=linprog（c,A,b,Aeq,beq,lb,ub,x0,options）
用options指定的优化参数进行最小化。
可以使用optimset来设置这些参数
上面的函数调用格式仅设置了最优解这一输出参数，如果需要更多的输出参数，则可以参照下面的调用格式：
1）[x,fval]=linprog（...）
在优化计算结束之时返回线性规划问题在解x处的目标函数值fval。
2）[x,fval,exitflag]=linprog（...）
在优化计算结束之时返回exitflag值，描述函数计算的退出条件。
3）[x,fval,exitflag,output]=linprog（...）
在优化计算结束之时返回返回结构变量output
4）[x,fval,exitflag,output,lambda]=linprog（...）
在优化计算结束之时返回线性规划问题最优解x处的拉格朗日乘子lambda
例1用MATLAB求解线性规划问题
算法：
c=[-1;-1];%目标函数，为转化为极小，故取目标函数中设计变量的相反数
A=[1-2;12];%线性不等式约束
b=[4;8];
lb=[0;0];%设计变量的边界约束，由于无上界，故设置ub=[Inf;Inf];
[x,fval]=linprog（c,A,b,[],[],lb,ub）
结果：
Optimizationterminated.
x=
6.0000
1.0000
fval=
-7.0000
例2用MATLAB求解线性规划问题
算法：
c=[-4;-3];%目标函数，为转化为极小，故取目标函数中设计变量的相反数
A=[34;33;42];%线性不等式约束
b=[12;10;8];
lb=[0;0];%设计变量的边界约束，由于无上界，故设置ub=[Inf;Inf]
ub=[Inf;Inf];
[x,fval,exitflag]=linprog（c,A,b,[],[],lb,ub）
结果：
x=
0.8000
2.4000
fval=
-10.4000
exitflag=
1
例3用MATLAB求解线性规划问题
算法：
c=[-1;-3;1];%目标函数，为转化为极小，故取目标函数中设计变量的相反数
Aeq=[112;-121];%线性等式约束
beq=[4;4];
lb=[0;0;0];%设计变量的边界约束，由于无上界，故设置ub=[Inf;Inf;Inf]
[x,fval,exitflag]=linprog（c,[],[],Aeq,beq,lb,ub）
结果：
x=
1.3333
2.6667
0.0000
fval=
-9.3333
exitflag=
1
例4用MATLAB求解线性规划问题
算法：
c=[-3;1;1];%目标函数，为转化为极小，故取目标函数中设计变量的相反数
A=[1-21;4-1-2];%线性不等式约束
b=[11;-3];
Aeq=[-201];%线性等式约束
beq=[1];
lb=[0;0;0];%设计变量的边界约束，由于无上界，故设置ub=[Inf;Inf;Inf]
ub=[Inf;Inf;Inf];
[x,fval,exitflag,output,lambda]=linprog（c,A,b,Aeq,beq,lb,ub）
结果：
x=
4.0000
1.0000
9.0000
fval=
-2.0000
exitflag=
1
output=%输出信息结构变量
iterations:
6%说明优化算法迭代6次
algorithm:
'large-scale:
interiorpoint'%说明采用的是大型规模的内点算法
cgiterations:
0%取值为0，为了向后兼容而设定
message:
'Optimizationterminated.'%退出信息
lambda=%最优解x处的拉格朗日乘子结构变量
ineqlin:
[2x1double]
eqlin:
-0.6667
upper:
[3x1double]
lower:
[3x1double]
§8.1.2整数规划的MATLAB求解方法
（一）用MATLAB求解一般混合整数规划问题
由于MATLAB优化工具箱中并未提供求解纯整数规划和混合整数规划的函数，因而需要自行根据需要和设定相关的算法来实现。
现在有许多用户发布的工具箱可以解决该类问题，例如比较著名的YALMIP，读者可以自行到网上下载相关的工具包并进行学习。
这里我们给出开罗大学的Sherif和Tawfik在MATLABCentral上发布的一个用于求解一般混合整数规划的程序，在此命名为intprog，笔者在原程序的基础上做了简单的修改，将其选择分枝变量的算法由自然序改造成分枝变量选择原则中的一种，即：
选择与整数值相差最大的非整数变量首先进行分枝。
intprog函数的调用格式如下：
[x,fval,exitflag]=intprog（c,A,b,Aeq,beq,lb,ub,M,TolXInteger）
该函数解决的整数规划问题为：
在上述标准问题中，假设
为
维设计变量，且问题具有不等式约束
个，等式约束
个，那么：
、
均为
维列向量，
为
维列向量，
为
维列向量，
为
维矩阵，
为
维矩阵。
在该函数中，输入参数有c,A,b,Aeq,beq,lb,ub,M和TolXInteger。
其中c为目标函数所对应设计变量的系数，A为不等式约束条件方程组构成的系数矩阵，b为不等式约束条件方程组右边的值构成的向量。
Aeq为等式约束方程组构成的系数矩阵，beq为等式约束条件方程组右边的值构成的向量。
lb和ub为设计变量对应的上界和下界。
M为具有整数约束条件限制的设计变量的序号，例如问题中设计变量为
，要求
和
为整数，则M=[2;3;6]；若要求全为整数，则M=1:
6，或者M=[1;2;3;4;5;6]。
TolXInteger为判定整数的误差限，即若某数x和最邻近整数相差小于该误差限，则认为x即为该整数。
在该函数中，输出参数有x,fval和exitflag。
其中x为整数规划问题的最优解向量，fval为整数规划问题的目标函数在最优解向量x处的函数值，exitflag为函数计算终止时的状态指示变量。
在MATLAB中实现intprog的代码和分析如下：
%整数规划的MATLAB实现
%OriginallyDesignedBySherifA.Tawfik,RevisedByLiMing,2009-12-29
%函数调用形式[x,fval,exitflag]=intprog（f,A,b,Aeq,beq,lb,ub,M,TolXInteger）
%函数求解如下形式的整数规划问题
%minf'*x
%subjectto
%A*x<=b
%Aeq*x=beq
%lb<=x<=ub
%M存储有整数约束的变量编号的向量
%TolXInteger是判定整数的误差限
%
%函数返回变量
%x:
整数规划的最优解
%fval:
目标函数在最优解处的函数值
%exitflag=1收敛到解x
%0达到线性规划的最大迭代次数
%-1无解
%
function[x,fval,exitflag]=intprog（f,A,b,Aeq,beq,lb,ub,M,TolXInteger）
%设置不显示求解线性规划过程中的提示信息
options=optimset（'display','off'）;
%上界的初始值
bound=inf;
%求解原问题P0的松弛线性规划Q0，首先获得问题的初始解
[x0,fval0]=linprog（f,A,b,Aeq,beq,lb,ub,[],options）;
%利用递归法进行二叉树的遍历，实现分枝定界法对整数规划的求解
[x,fval,exitflag,b]=rec_BranchBound（f,A,b,Aeq,beq,lb,ub,x0,fval0,M,TolXInteger,bound）;
%分枝定界法的递归算法,x为问题的初始解，v是目标函数在x处的取值
function[xx,fval,exitflag,bb]=rec_BranchBound（f,A,b,Aeq,beq,lb,ub,x,v,M,TolXInteger,bound）
options=optimset（'display','off'）;
%求解不考虑整数约束的松弛线性规划
[x0,fval0,exitflag0]=linprog（f,A,b,Aeq,beq,lb,ub,[],options）;
%如果算法结束状态指示变量为负值，即表示无可行解，返回初始输入
%或者所目标函数值大于已经获得的上界，返回初始输入
ifexitflag0<=0|fval0>bound
xx=x;
fval=v;
exitflag=exitflag0;
bb=bound;
return;
end
%确定所有变量是否均为整数，如是，则返回
ind=find（abs（x0（M）-round（x0（M）））>TolXInteger）;%该条件表示x0（M）不是整数
%如果都是整数
ifisempty（ind）exitflag=1;
%如果当前的解优于已知的最优解，则将当前解作为最优解
iffval0%否则，返回原来的解
else
xx=x;
fval=v;
bb=bound;
end
return;
end
%程序运行至此，说明松弛线性规划的解是一个可行解且目标函数值比当前记录的上界要小，只是某些变量的值并非整数，于是在此选择合适的变量进行分枝形成两个子问题，分别进行递归求解，该处选择与整数值相差最大的非整数变量首先进行分枝形成两个子问题,第一个非整数变量的序号，且记录该变量与其最邻近的整数之差的绝对值
[rowcol]=size（ind）;
br_var=M（ind
（1））;
br_value=x（br_var）;
flag=abs（br_value-floor（br_value）-0.5）;
%用于查找非整数设计变量中整数值相差最大的设计变量，即每当遇到与其最邻近的整数差别更大的非整数设计变量之时，即记录下该设计变量的序号，直至遍历完所有非整数变量
fori=2:
col
tempbr_var=M（br_var）;
tempbr_value=x（br_var）
temp_flag=abs（tempbr_value-floor（tempbr_value）-0.5）;
iftemp_flag>flag
br_var=tempbr_var;
br_value=tempbr_value;
flag=temp_flag;
end
end
ifisempty（A）[rc]=size（Aeq）;
else[rc]=size（A）;
end
%分枝后第1个子问题的设置,
%添加约束条件Xi<=floor（Xi），i即为上面找到的最接近整数的非整数设计变量的序号
A1=[A;zeros（1,c）];
A1（end,br_var）=1;
b1=[b;floor（br_value）];
%第2个子问题的设置,添加约束条件Xi>=ceil（Xi），i即为上面找到的最接近整数的非整数设计变量的序号
A2=[A;zeros（1,c）];
A2（end,br_var）=-1;
b2=[b;-ceil（br_value）];
%分枝后的第一个子问题的递归求解
[x1,fval1,exitflag1,bound1]=rec_BranchBound（f,A1,b1,Aeq,beq,lb,ub,x0,fval0,M,TolXInteger,bound）;
exitflag=exitflag1;
ifexitflag1>0&bound1xx=x1;
fval=fval1;
bound=bound1;
bb=bound1;
else
xx=x0;
fval=fval0;
bb=bound;
end
%分枝后的第二个子问题的递归求解
[x2,fval2,exitflag2,bound2]=rec_BranchBound（f,A2,b2,Aeq,beq,lb,ub,x0,fval0,M,TolXInteger,bound）;
ifexitflag2>0&bound2exitflag=exitflag2;
xx=x2;
fval=fval2;
bb=bound2;
end
例1求解整数规划问题
算法：
c=[-5;-8];
A=[11;59];
b=[6;45];
lb=[0;0];
M=[1;2];%x1,x2均要求为整数变量
Tol=1e-8;%判断是否整数的误差限
[x,fval]=linprog（c,A,b,[],[],lb,[]）%松弛问题的解
[x1,fval1]=intprog（c,A,b,[],[],lb,[],M,Tol）%整数规划的解
结果：
Optimizationterminated.
x=
2.2500
3.7500
fval=
-41.2500
x1=
0
5
fval1=
-40.0000
上例中x和fval为整数规划问题对应松弛线性规划的最优解和最优值，x1和fval1为整数规划的最优解和最优值，可见本例中无论如何对x进行取整操作都无法得到整数规划问题的最优解。
例2求解整数规划问题：
算法：
c=[-1;-1];
A=[-42;42;0-2];
b=[-1;11;-1];
lb=[0;0];
M=[1;2];%均要求为整数变量
Tol=1e-8;%判断是否整数的误差限
[x,fval]=linprog（c,A,b,[],[],lb,[]）%求解原问题松弛线性规划
[x1,fval1]=intprog（c,A,b,[],[],lb,[],M,Tol）%求最优解整数解
结果：
x=%松弛线性规划问题的最优解
1.5000
2.5000
fval=
-4.0000
x1=%整数规划的最优解
2
1
fval2=
-3.0000
（二）用MATLAB求解0-1规划问题
在MATLAB优化工具箱中，提供了专门用于求解0-1规划问题的函数bintprog，其算法基础即为分枝界定法，在MATLAB中调用bintprog函数求解0-1规划时，需要遵循MATLAB中对0-1规划标准性的要求。
0-1规划问题的MATLAB标准型
在上述模型中，目
升级会员 