运筹学与最优化MATLAB编程.docx
- 文档编号:9568381
- 上传时间:2023-02-05
- 格式:DOCX
- 页数:12
- 大小:103.75KB
运筹学与最优化MATLAB编程.docx
《运筹学与最优化MATLAB编程.docx》由会员分享,可在线阅读,更多相关《运筹学与最优化MATLAB编程.docx(12页珍藏版)》请在冰豆网上搜索。
运筹学与最优化MATLAB编程
运筹学与最优化MATLAB编程
实验报告
院系:
专业:
姓名:
学号:
指导老师:
完成日期:
割平面法求解整数规划问题
一、引言:
通过对MATLAB实践设计的学习,学会使用MATLAB解决现实生活中的问题。
该设计是在MATLAB程序设计语言的基础上,对实际问题建立数学模型并设计程序,使用割平面法解决一个整数规划问题。
经实验,该算法可成功运行并求解出最优整数解。
二、算法说明:
割平面法有许多种类型,本次设计的原理是依据Gomory的割平面法。
Gomory割平面法首先求解非整数约束的线性规划,再选择一个不是整数的基变量,定义新的约束,增加到原来的约束中,新的约束缩小了可行域,但是保留了原问题的全部整数可行解。
算法具体设计步骤如下:
1、首先,求解原整数规划对应的线性规划
,设最优解为x*。
2、如果最优解的分量均为整数,则x*为原整数规划的最优解;否则任选一个x*中不为整数的分量,设其对应的基变量为xp,定义包含这个基变量的切割约束方程
,其中xp为非基变量。
3、令
,
,其中[]为高斯函数符号,表示不大于某数的最大整数。
将切割约束方程变换为
,由于0<
<1,0<
<1,所以有
,因为自变量为整数,则
也为整数,所以进一步有
。
4、将切割方程加入约束方程中,用对偶单纯形法求解线性规划
,然后在转入步骤2进行求解,直到求出最优整数解停止迭代。
三、程序实现:
程序设计流程图如图1,具体设计代码(见附录)。
图1
四、算例分析:
已知AM工厂是一个拥有四个车间的玩具生产厂商,该厂商今年新设计出A、B、C、D、E、F六种玩具模型,根据以前的生产情况及市场调查预测,得知生产每单位产品所需的工时、每个车间在一季度的工时上限以及产品的预测价格,如下表所示。
问:
每种设计产品在这个季度各应生产多少,才能使AM工厂这个季度的生产总值达到最大?
产品
车间
A
B
C
D
E
F
每个车间每季度工时上限
甲
乙
丙
丁
0.01
0.02
0.01
0.02
0.01
0.03
0.03
0.05
0.03
0.05
0.03
0.08
850
700
100
900
单价(元)
20
14
16
36
32
30
1、问题分析并建立模型:
由题意可知这是一个求解产量使产值最大的整数规划问题。
根据上述问题和已知数据,可以假设每种产品在这个季度各应生产产量分别为:
x1、x2、x3、x4、x5、x6,则有以下线性方程组
maxZ=20x1+14x2+16x3+36x4+32x5+30x6
2、实验步骤:
首先引入松弛变量x7、x8、x9、x10,使其化为标准型
minZ=-20x1-14x2-16x3-36x4-32x5-30x6
其次从标准型可表示出约束系数矩阵、右端项常数矩阵、目标系数矩阵分别为A、b、c。
然后调用DividePlane函数,使用割平面法进行求解。
在MATLAB的命令窗口输入一下命令:
>>A=[0.010.010.010.030.030.031000;0.02000.05000100;00.02000.0500010;000.03000.080001];
>>c=[-20;-14;-16;-36;-32;-30];
>>b=[850;700;100;900];
>>[intx,intf]=DividePlane(A,c,b,[78910])
3、实验结果及分析:
intx=
35000500030000000
intf=
-1250000
实验结果求出的目标函数值是化为标准型的最小值,则转化为原问题的目标函数值应取相反数,所以从实验结果可知:
生产各种产品的产量分别应为为,生产A35000、生产B5000、生产C30000、生产D、E、F均为0,此时的季度产值为最大即1250000元。
该结果是可信的,故通过该实例说明该程序能够运用于实际,用来解决实际生活中求解整数规划的问题。
五、结束语:
Matlab是个很强大的软件,提供了大量的函数来处理各种数学、工程、运筹等的问题,并且含有处理二维、三维图形的功能,使用matlab能够解决许多实际生活中的问题。
通过这个学期的学习,仅是了解了matlab的部分函数功能和简单的GUI界面设计,掌握了一些基本的程序编写技能,同时,在老师的指导下简单了解了使用LinGo和Excel解决线性和非线性规划问题的求解方法,收获相当丰富,同时认识到要学好matlab仍然需要一个长期的过程。
六、参考文献:
[1]龚纯,王正林.精通MATLAB最优化计算.北京:
电子工业出版社,2009
[2]吴祈宗,郑志勇,邓伟等.运筹学与最优化MATLAB编程.北京:
机械工业出版社,2009
[3]邓成梁.运筹学的原理和方法(第二版).武汉:
华中科技大学出版社,2002
七、附录:
function[intx,intf]=DividePlane(A,c,b,baseVector)
%功能:
用割平面法求解整数规划
%调用格式:
[intx,intf]=DividePlane(A,c,b,baseVector)
%其中,A:
约束矩阵;
%c:
目标函数系数向量;
%b:
约束右端向量;
%baseVector:
初始基向量;
%intx:
目标函数取最小值时的自变量值;
%intf:
目标函数的最小值;
sz=size(A);
nVia=sz
(2);
n=sz
(1);
xx=1:
nVia;
iflength(baseVector)~=n
disp('基变量的个数要与约束矩阵的行数相等!
');
mx=NaN;
mf=NaN;
return;
end
M=0;
sigma=-[transpose(c)zeros(1,(nVia-length(c)))];
xb=b;
%首先用单纯形法求出最优解
while1
[maxs,ind]=max(sigma);
%--------------------用单纯形法求最优解--------------------------------------
ifmaxs<=0%当检验数均小于0时,求得最优解。
vr=find(c~=0,1,'last');
forl=1:
vr
ele=find(baseVector==l,1);
if(isempty(ele))
mx(l)=0;
else
mx(l)=xb(ele);
end
end
ifmax(abs(round(mx)-mx))<1.0e-7%判断最优解是否为整数解,如果是整数解。
intx=mx;
intf=mx*c;
return;
else%如果最优解不是整数解时,构建切割方程
sz=size(A);
sr=sz
(1);
sc=sz
(2);
[max_x,index_x]=max(abs(round(mx)-mx));
[isB,num]=find(index_x==baseVector);
fi=xb(num)-floor(xb(num));
fori=1:
(index_x-1)
Atmp(1,i)=A(num,i)-floor(A(num,i));
end
fori=(index_x+1):
sc
Atmp(1,i)=A(num,i)-floor(A(num,i));
end
Atmp(1,index_x)=0;%构建对偶单纯形法的初始表格
A=[Azeros(sr,1);-Atmp(1,:
)1];
xb=[xb;-fi];
baseVector=[baseVectorsc+1];
sigma=[sigma0];
%-------------------对偶单纯形法的迭代过程----------------------
while1
%----------------------------------------------------------
ifxb>=0%判断如果右端向量均大于0,求得最优解
ifmax(abs(round(xb)-xb))<1.0e-7%如果用对偶单纯形法求得了整数解,则返回最优整数解
vr=find(c~=0,1,'last');
forl=1:
vr
ele=find(baseVector==l,1);
if(isempty(ele))
mx_1(l)=0;
else
mx_1(l)=xb(ele);
end
end
intx=mx_1;
intf=mx_1*c;
return;
else%如果对偶单纯形法求得的最优解不是整数解,继续添加切割方程
sz=size(A);
sr=sz
(1);
sc=sz
(2);
[max_x,index_x]=max(abs(round(mx_1)-mx_1));
[isB,num]=find(index_x==baseVector);
fi=xb(num)-floor(xb(num));
fori=1:
(index_x-1)
Atmp(1,i)=A(num,i)-floor(A(num,i));
end
fori=(index_x+1):
sc
Atmp(1,i)=A(num,i)-floor(A(num,i));
end
Atmp(1,index_x)=0;%下一次对偶单纯形迭代的初始表格
A=[Azeros(sr,1);-Atmp(1,:
)1];
xb=[xb;-fi];
baseVector=[baseVectorsc+1];
sigma=[sigma0];
continue;
end
else%如果右端向量不全大于0,则进行对偶单纯形法的换基变量过程
minb_1=inf;
chagB_1=inf;
sA=size(A);
[br,idb]=min(xb);
forj=1:
sA
(2)
ifA(idb,j)<0
bm=sigma(j)/A(idb,j);
ifbm minb_1=bm; chagB_1=j; end end end sigma=sigma-A(idb,: )*minb_1; xb(idb)=xb(idb)/A(idb,chagB_1); A(idb,: )=A(idb,: )/A(idb,chagB_1); fori=1: sA (1) ifi~=idb xb(i)=xb(i)-A(i,chagB_1)*xb(idb); A(i,: )=A(i,: )-A(i,chagB_1)*A(idb,: ); end end baseVector(idb)=chagB_1; end %------------------------------------------------------------ end %--------------------对偶单纯形法的迭代过程--------------------- end else%如果检验数有不小于0的,则进行单纯形算法的迭代过程 minb=inf; chagB=inf; forj=1: n ifA(j,ind)>0 bz=xb(j)/A(j,ind); ifbz minb=bz; chagB=j; end end end sigma=sigma-A(chagB,: )*maxs/A(chagB,ind); xb(chagB)=xb(chagB)/A(chagB,ind); A(chagB,: )=A(chagB,: )/A(chagB,ind); fori=1: n ifi~=chagB xb(i)=xb(i)-A(i,ind)*xb(chagB); A(i,: )=A(i,: )-A(i,ind)*A(chagB,: ); end end baseVector(chagB)=ind; end M=M+1; if(M==1000000) disp('找不到最优解! '); mx=NaN; minf=NaN; return; end end
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 运筹学 优化 MATLAB 编程