視頻推薦：B站_數(shù)學(xué)建模老哥
目錄
一、整數(shù)規(guī)劃基本原理
1.1 整數(shù)規(guī)劃的分類
1.2 整數(shù)規(guī)劃的特點
1.3 案例
1.4? 整數(shù)規(guī)劃的數(shù)學(xué)模型一般形式
二、整數(shù)線性規(guī)劃的求解方法
2.1 分枝定界法
2.1.1 分枝定界法的求解過程
2.1.2 案例
2.1.3 代碼實現(xiàn)
2.2 割平面法
2.2.1 割平面法的基本思想
2.2.2 割平面法的基本步驟
2.2.3 案例和代碼實現(xiàn)
一、整數(shù)規(guī)劃基本原理
數(shù)學(xué)規(guī)劃中的變量（部分或全部）限制為整數(shù)時，稱為整數(shù)規(guī)劃。若在線性規(guī)劃模型中，變量限制為整數(shù)，則稱為整數(shù)線性規(guī)劃。目前所流行的求解整數(shù)規(guī)劃的方法，往往只適用于整數(shù)線性規(guī)劃。目前還沒有一種方法能有效地求解一切整數(shù)規(guī)劃。
1.1 整數(shù)規(guī)劃的分類
(1)變量全限制為整數(shù)時，稱為純（完全）整數(shù)規(guī)劃。
(2)變量部分限制為整數(shù)的，稱為混合整數(shù)規(guī)劃。
(3)變量取值要么為0要么為1，稱為0-1規(guī)劃。
1.2 整數(shù)規(guī)劃的特點
1. 原線性規(guī)劃有最優(yōu)解，當自變量限制為整數(shù)后，其整數(shù)規(guī)劃解出現(xiàn)下述情況：
? ? (1)原線性規(guī)劃最優(yōu)解全是整數(shù)，則整數(shù)規(guī)劃最優(yōu)解與線性規(guī)劃最優(yōu)解一致。
??? (2)整數(shù)規(guī)劃可能不存在可行解。
??? (3)有可行解(當然就存在最優(yōu)解)，但最優(yōu)解值變差。
2. 整數(shù)規(guī)劃最優(yōu)解不能按照實數(shù)最優(yōu)解簡單取整而獲得。
1.3 案例
合理下料問題：
設(shè)用某型號的圓鋼下零件A1,A2...,Am的毛坯。在一根圓鋼上下料的方式有B1,B2,... Bn種，每種下料方式可以得到各種零件的毛坯數(shù)以及每種零件的需要量，如表所示。問怎樣安排下料方式，使得即滿足需要，所用的原材料又最少?
如表所示：
?模型：
其中，xj表示用Bj種方式下料根數(shù)，aij為每種下料方式可以得到各種零件的毛坯數(shù)，bi每種零件的需要量。
1.4? 整數(shù)規(guī)劃的數(shù)學(xué)模型一般形式
另外補充： 整數(shù)規(guī)劃與松弛的線性規(guī)劃之間的關(guān)系。
不難看出兩者之間的關(guān)系。
二、整數(shù)線性規(guī)劃的求解方法
????????從數(shù)學(xué)模型上看整數(shù)規(guī)劃似乎是線性規(guī)劃的一種特殊形式，求解只需在線性規(guī)劃的基礎(chǔ)上，通過舍入取整，尋求滿足整數(shù)要求的解即可。但實際上兩者卻有很大的不同，通過舍入得到的解（整數(shù)）也不一定就是最優(yōu)解，有時甚至不能保證所得到的解是整數(shù)可行解。
2.1 分枝定界法
分枝定界法（branch and bound）是一種求解整數(shù)規(guī)劃問題的最常用算法。分支定界法是一種搜索與迭代的方法，選擇不同的分支變量和子問題進行分支。
2.1.1 分枝定界法的求解過程
- 1. 先求出相應(yīng)松弛問題的最優(yōu)解。
- 2. 若松弛問題無可行解，則整數(shù)線性規(guī)劃問題無可行解。
- 3. 若求得的松弛問題最優(yōu)解符合整數(shù)要求，則是整數(shù)線性規(guī)劃問題的最優(yōu)解。若不滿足整數(shù)條件，則任選一個不滿足整數(shù)條件的變量 來構(gòu)造新的約束，分別添加到松弛問題中形成兩個子問題：
- ????????????????????????????????????????????????????? 和
- 依次在縮小的可行域中求解新構(gòu)造的線性規(guī)劃的最優(yōu)解，并重復(fù)上述過程，直到子問題無解或有整數(shù)最優(yōu)解（被查清）。
Q：為什么在步驟3中不滿足整數(shù)條件時要添加上述兩個約束條件？
A：因為變量 是一個不滿足整數(shù)條件的變量，它注定無法構(gòu)成整數(shù)規(guī)劃的最優(yōu)解，因此我們要向變量 的兩邊去尋找合適的整數(shù)最優(yōu)解。
注意：構(gòu)造新的約束條件是分別到整數(shù)規(guī)劃問題對應(yīng)的松弛問題中，獨立構(gòu)成兩個不同的子問題再求解。
詳情參考：分枝定界法
2.1.2 案例
分枝定界法總體上有點像二叉樹，能看懂下面的案例基本就能理解分枝定界法。
2.1.3 代碼實現(xiàn)
數(shù)學(xué)模型如下：
1. 先創(chuàng)建intprog.m文件：
function [x,fval,status] = intprog(f,A,B,I,Aeq,Beq,lb,ub,e)
%整數(shù)規(guī)劃求解函數(shù) intprog()
% 其中 f為目標函數(shù)向量
% A和B為不等式約束 Aeq與Beq為等式約束
% I為整數(shù)約束
% lb與ub分別為變量下界與上界
% x為最優(yōu)解，fval為最優(yōu)值
%例子:
% maximize 20 x1 + 10 x2
% S.T.
% 5 x1 + 4 x2 <=24
% 2 x1 + 5 x2 <=13
% x1, x2 >=0
% x1, x2是整數(shù)
% f=[-20, -10];
% A=[ 5 4; 2 5];
% B=[24; 13];
% lb=[0 0];
% ub=[inf inf];
% I=[1,2];
% e=0.000001;
% [x v s]= IP(f,A,B,I,[],[],lb,ub,,e)
% x = 4 1 v = -90.0000 s = 1
% 控制輸入?yún)?shù)
if nargin < 9, e = 0.00001;
if nargin < 8, ub = [];
if nargin < 7, lb = [];
if nargin < 6, Beq = [];
if nargin < 5, Aeq = [];
if nargin < 4, I = [1:length(f)];
end, end, end, end, end, end
%求解整數(shù)規(guī)劃對應(yīng)的線性規(guī)劃,判斷是否有解
options = optimset('display','off'); %不需要每次優(yōu)化結(jié)果都輸出
[x0,fval0,exitflag] = linprog(f,A,B,Aeq,Beq,lb,ub,options); %決策向量，最優(yōu)解，是否有解
if exitflag < 0
disp('沒有合適整數(shù)解');
x = x0;
fval = fval0;
status = exitflag;
return;
else
%采用分支定界法求解
bound = inf;
[x,fval,status] = branchbound(f,A,B,I,x0,fval0,bound,Aeq,Beq,lb,ub,e);
end
2. 再創(chuàng)建branchbound.m文件
function [newx,newfval,status,newbound] = branchbound(f,A,B,I,x,fval,bound,Aeq,Beq,lb,ub,e)
% 分支定界法求解整數(shù)規(guī)劃
% f,A,B,Aeq,Beq,lb,ub與線性規(guī)劃相同
% I為整數(shù)限制變量的向量
% x為初始解，fval為初始值
options = optimset('display','off');
[x0,fval0,status0]=linprog(f,A,B,Aeq,Beq,lb,ub,options);
%遞歸中的最終退出條件
%無解或者解比現(xiàn)有上界大則返回原解
if status0 <= 0 || fval0 >= bound
newx = x;
newfval = fval;
newbound = bound;
status = status0;
return;
end
%是否為整數(shù)解,如果是整數(shù)解則返回
intindex = find(abs(x0(I) - round(x0(I))) > e);
if isempty(intindex) %判斷是否為空值
newx(I) = round(x0(I));
newfval = fval0;
newbound = fval0;
status = 1;
return;
end
%當有非整可行解時，則進行分支求解
%此時必定會有整數(shù)解或空解
%找到第一個不滿足整數(shù)要求的變量
n = I(intindex(1));
addA = zeros(1,length(f));
addA(n) = 1;
%構(gòu)造第一個分支 x<=floor(x(n))
A = [A;addA];
B = [B,floor(x(n))];%向下取整
[x1,fval1,status1,bound1] = branchbound(f,A,B,I,x0,fval0,bound,Aeq,Beq,lb,ub,e);
A(end,:) = [];
B(:,end) = [];
%解得第一個分支，若為更優(yōu)解則替換，若不是則保持原狀
status = status1;
if status1 > 0 && bound1 < bound
newx = x1;
newfval = fval1;
bound = fval1;
newbound = bound1;
else
newx = x0;
newfval = fval0;
newbound = bound;
end
%構(gòu)造第二分支
A = [A;-addA];
B = [B,-ceil(x(n))];%向上取整
[x2,fval2,status2,bound2] = branchbound(f,A,B,I,x0,fval0,bound,Aeq,Beq,lb,ub,e);
A(end,:) = [];
B(:,end) = [];
%解得第二分支，并與第一分支做比較，如果更優(yōu)則替換
if status2 > 0 && bound2 < bound
status = status2;
newx = x2;
newfval = fval2;
newbound = bound2;
end3. 測試
function[x. fval.stetus] = intprog(f,A,B,I,Aeq,Beq,lb,ub,e)
f = [-20, -10]
A = [5,4;2,5]
B = [24,13]
lb = [0,0]
[x, fval, status] = intprog(f, A, B, [1,2],[],[],lb)2.2 割平面法
2.2.1 割平面法的基本思想
- 1. 如果松弛問題()無解，則()無解;
- 2. 如果()的最優(yōu)解為整數(shù)向量，則也是()的最優(yōu)解;
- 3. 如果()的解含有非整數(shù)分量，則對()增加割平面條件：即對()增加一個線性約束，將()的可行區(qū)域割掉一塊，使得非整數(shù)解恰好在割掉的一塊中，但又沒有割掉原問題()的可行解，得到問題()，重復(fù)上述的過程。
2.2.2 割平面法的基本步驟
注：這里的第二、三步是推理過程，簡略來看，只需要看第一步、第二步中1和第四步即可。?
第一步：求解整數(shù)規(guī)劃對應(yīng)松弛問題是否有整數(shù)最優(yōu)解。若有整數(shù)最優(yōu)解，則它也是整數(shù)規(guī)劃的最優(yōu)解，否則，轉(zhuǎn)到第二步。
第二步：1. 引入松弛變量轉(zhuǎn)化為等式約束：
其中，為決策變量，為松弛變量，為松弛變量的系數(shù)，與原約束條件一致。
2. 然后，將松弛變量的系數(shù)分別劃分為整數(shù)部分和小數(shù)部分，如下：
其中，表示向下取整，也表示向下取整，表示小數(shù)部分（減去它整數(shù)部分得到的小數(shù)，小于0）。
同樣，將求和結(jié)果劃分為整數(shù)部分和小數(shù)部分，如下：
其中，表示向下取整，也表示向下取整，表示小數(shù)部分（同上，小于0）。
第三步：將劃分后的和代入所構(gòu)建的方程中，可得：
我們將所有的整數(shù)部分放入公式左側(cè)，所有的小數(shù)部分放入公式右側(cè)：
因為小于0，所以
，則可推出：
????????????????????????????????????????????????????????????????
這樣就對決策變量進行了一次割平面（增加約束條件）。
第四步：將約束條件 添加到原數(shù)學(xué)模型中，重復(fù)第一步。
2.2.3 案例和代碼實現(xiàn)
已知AM工廠是一個擁有四個車間的玩具生產(chǎn)廠商，該廠商今年新設(shè)計出A、B、C、D、E、F六種玩具模型，根據(jù)以前的生產(chǎn)情況及市場調(diào)查預(yù)測，得知生產(chǎn)每單位產(chǎn)品所需的工時、每個車間在一季度的工時上限以及產(chǎn)品的預(yù)測價格，如下表所示。問:每種設(shè)計產(chǎn)品在這個季度各應(yīng)生產(chǎn)多少，才能使AM工廠這個季度的生產(chǎn)總值達到最大?
1. 建立模型
?其中，x1，x2，x3，x4，x5，x6是A、B、C、D、E、F六種玩具模型的生產(chǎn)數(shù)量。注意，是車間之間是合作生產(chǎn)，而不是獨立生產(chǎn)。
2. 引入松弛變量轉(zhuǎn)化為等式約束
3. 代碼實現(xiàn)
DividePlane.m文章來源：http://www.zghlxwxcb.cn/news/detail-651281.html
function [intx,intf] = DividePlane(A,c,b,baseVector)
%功能：用割平面法求解整數(shù)規(guī)劃
%調(diào)用格式：[intx,intf]=DividePlane(A,c,b,baseVector)
%其中, A：約束矩陣；
% c：目標函數(shù)系數(shù)向量；
% b：約束右端向量；
% baseVector：初始基向量；
% intx：目標函數(shù)取最值時的自變量值；
% intf：目標函數(shù)的最值；
sz = size(A);
nVia = sz(2);%獲取有多少決策變量
n = sz(1);%獲取有多少約束條件
xx = 1:nVia;
if length(baseVector) ~= n
disp('基變量的個數(shù)要與約束矩陣的行數(shù)相等！');
mx = NaN;
mf = NaN;
return;
end
M = 0;
sigma = -[transpose(c) zeros(1,(nVia-length(c)))];
xb = b;
%首先用單純形法求出最優(yōu)解
while 1
[maxs,ind] = max(sigma);
%--------------------用單純形法求最優(yōu)解--------------------------------------
if maxs <= 0 %當檢驗數(shù)均小于0時，求得最優(yōu)解。
vr = find(c~=0 ,1,'last');
for l=1:vr
ele = find(baseVector == l,1);
if(isempty(ele))
mx(l) = 0;
else
mx(l)=xb(ele);
end
end
if max(abs(round(mx) - mx))<1.0e-7 %判斷最優(yōu)解是否為整數(shù)解，如果是整數(shù)解。
intx = mx;
intf = mx*c;
return;
else %如果最優(yōu)解不是整數(shù)解時，構(gòu)建切割方程
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; %構(gòu)建對偶單純形法的初始表格
A = [A zeros(sr,1);-Atmp(1,:) 1];
xb = [xb;-fi];
baseVector = [baseVector sc+1];
sigma = [sigma 0];
%-------------------對偶單純形法的迭代過程----------------------
while 1
%----------------------------------------------------------
if xb >= 0 %判斷如果右端向量均大于0，求得最優(yōu)解
if max(abs(round(xb) - xb))<1.0e-7 %如果用對偶單純形法求得了整數(shù)解，則返回最優(yōu)整數(shù)解
vr = find(c~=0 ,1,'last');
for l=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 %如果對偶單純形法求得的最優(yōu)解不是整數(shù)解，繼續(xù)添加切割方程
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));
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];
xb = [xb;-fi];
baseVector = [baseVector sc+1];
sigma = [sigma 0];
continue;
end
else %如果右端向量不全大于0，則進行對偶單純形法的換基變量過程
minb_1 = inf;
chagB_1 = inf;
sA = size(A);
[br,idb] = min(xb);
for j=1:sA(2)
if A(idb,j)<0
bm = sigma(j)/A(idb,j);
if bm<minb_1
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);
for i =1:sA(1)
if i ~= 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 %如果檢驗數(shù)有不小于0的，則進行單純形算法的迭代過程
minb = inf;
chagB = inf;
for j=1:n
if A(j,ind)>0
bz = xb(j)/A(j,ind);
if bz<minb
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);
for i =1:n
if i ~= 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('找不到最優(yōu)解！');
mx = NaN;
minf = NaN;
return;
end
end
test.m文章來源地址http://www.zghlxwxcb.cn/news/detail-651281.html
c = [-20 ; -14 ; -16 ; -36 ; -32 ; -30]; % 目標函數(shù)系數(shù)向量
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]; % 加上松弛變量后約束條件的系數(shù)矩陣
b = [850; 700; 100; 900;]; % 約束右端向量
[intx , intf] = DividePlane(A,c,b,[7 8 9 10]); % 初始基向量的下標 7 8 9 10
intx
intf = -intf % 取反求最大到了這里，關(guān)于數(shù)學(xué)建模（三）整數(shù)規(guī)劃的文章就介紹完了。如果您還想了解更多內(nèi)容，請在右上角搜索TOY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關(guān)文章，希望大家以后多多支持TOY模板網(wǎng)！
