运筹学与最优化MATLAB编程pctools5.0
实验报告
一、 引言：
医用钛通过对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.01 0.01 0.01 0.03 0.03 0.03 1 0 0 0;0.02 0 0 0.05 0 0 0 1 0 0;0 0.02 0 0 0.05 0 0 0 1 0;0 0 0.03 0 0 0.08 0 0 0 1];呋喃树脂砂
>> c=[-20;-14;-16;-36;-32;-30];
>> b=[850;700;100;900];
>> [intx,intf]=DividePlane(A,c,b,[7 8 9 10])
3、实验结果及分析：
intx =
35000 5000 30000 0 0 0
intf =
-1250000
实验结果求出的目标函数值是化为标准型的最小值，则转化为原问题的目标函数值应取相反数，所以从实验结果可知：生产各种产品的产量分别应为为，生产A 35000、生产B 5000、生产C 30000、生产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;
if length(baseVector) ~= n
disp('基变量的个数要与约束矩阵的行数相等！');
mx = NaN;
mf = NaN;
return;
end
M = 0;
sigma = -[transpose(c) zeros(1,(nVia-length(c)))];
xb = b;
%首先用单纯形法求出最优解
while 1
[maxs,ind] = max(sigma);
%--------------------用单纯形法求最优解--------------------------------------
if maxs <= 0 %当检验数均小于0时，求得最优解。
vr = find(c~=0 ,1,'last');
for l=1:vr
ele = find(baseVector == l,1);
if(isempty(ele))
mx(l) = 0;rbd-573
else
mx(l)=xb(ele);
end
end
if max(abs(round(mx) - mx))<1.0e-7 %判断最优解是否为整数解，如果是整数解。
intx = mx;松下vs6
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));
for i=1:(index_x-1)
Atmp(1,i) = A(num,i) - floor(A(num,i));
end
for i=(index_x+1):sc
Atmp(1,i) = A(num,i) - floor(A(num,i));
end
Atmp(1,index_x) = 0; %构建对偶单纯形法的初始表格
A = [A zeros(sr,1);-Atmp(1,:) 1];
豫公网安备41160202000603
站长QQ:729038198
关于我们
投诉建议