Matlab学习系列26 整数规划.docx
《Matlab学习系列26 整数规划.docx》由会员分享，可在线阅读，更多相关《Matlab学习系列26 整数规划.docx（19页珍藏版）》请在冰豆网上搜索。
Matlab学习系列26整数规划
26.整数规划
全部变量限制为整数的规划问题，称为纯整数规划；部分变量限制为整数的规划问题，称为混合整数规划；变量只取0或1的规划问题，称为0-1整数规划。
整数规划问题，建议使用Lingo软件求解。
常用的整数规划问题解法有：
（1）分枝定界法：
可求纯或混合整数线性规划；
（2）割平面法：
可求纯或混合整数线性规划；
（3）隐枚举法：
用于求解0-1整数规划，有过滤法和分枝法；
（4）匈牙利法：
解决指派问题（0-1规划特殊情形）；
（5）蒙特卡罗法：
求解各种类型规划。
一、分枝定界法
分支定界法的基本思想是：
设有最大化的整数规划问题A，先解与之相应的线性规划问题B，若B的最优解不符合A的整数条件，那么B的最优目标函数必是A的最优目标函数z*的上界，记作z2,而A的任意可行解的目标函数值将是z*的一个下界z1,分支定界法就是将B的可行域分成子区域（称为分支）的方法，逐步减小z2和增大z1,最终求到z*。
例1分枝定界法原理示例：
用Lingo软件求解：
代码：
max5x1+8x2
st
x1+x2<=6
5x1+9x2<=45
end
gin2
运行结果：
Globaloptimalsolutionfound.
Objectivevalue:
40.00000
Objectivebound:
40.00000
Infeasibilities:
0.000000
Extendedsolversteps:
0
Totalsolveriterations:
0
VariableValueReducedCost
X10.000000-5.000000
X25.000000-8.000000
RowSlackorSurplusDualPrice
140.000001.000000
21.0000000.000000
30.0000000.000000
二、0-1整数规划
变量xi只能取值0,1，该约束条件可表示为：
0≤xi≤1,xi∈N或xi（1-xi）=0
1.隐枚举法
例2求解下列0-1规划问题：
求解思路：
（1）先试探性地求一个可行解，易看出（x1,x2,x3）=（1,0,0）满足约束条件，故是一个可行解，相应的目标函数值为z=3.
（2）由于是求极大值，故目标值z<3的解，不必检验是否满足约束条件即可删除，于是可增加一个约束条件（称为过滤条件）：
（e）
（3）用全部枚举法，3个变量共23=8种可能的组合，用过滤条件（并计算目标函数值，不断改进过滤条件）筛选每个可能的组合，最终得到问题的最优解。
（x1,x2,x3）
目标值
约束条件
过滤条件
（a）
（b）
（c）
（d）
（e）
（0,0,0）
0
×
（1,0,0）
3
√
√
√
√
√
3x1-2x2+x3≥3
（0,1,0）
-2
×
（0,0,1）
5
√
√
√
√
√
3x1-2x2+x3≥5
（1,1,0）
1
×
（1,0,1）
8
√
√
√
√
√
3x1-2x2+x3≥8
（1,1,1）
6
×
（0,1,1）
3
×
从而得到最优解为（1,0,1）,最优值为z*=8.
用Lingo软件求解：
代码：
max3x1-2x2+5x3
st
x1+2x2-x3<=2
x1+4x2+x3<=4
x1+x2<=3
4x2+x3<=6
end
int3
运行结果：
Globaloptimalsolutionfound.
Objectivevalue:
8.000000
Objectivebound:
8.000000
Infeasibilities:
0.000000
Extendedsolversteps:
0
Totalsolveriterations:
0
VariableValueReducedCost
X11.000000-3.000000
X20.0000002.000000
X31.000000-5.000000
RowSlackorSurplusDualPrice
18.0000001.000000
22.0000000.000000
32.0000000.000000
42.0000000.000000
55.0000000.000000
2.指派问题
某公司准备分派n个工人X1,…,Xn做n件工作Y1,…,Yn,但由于任务性质和各人专长不同，因此各人完成每件工作的效率也就不同，于是产生一个问题：
应指派哪个人去完成哪件工作，使得工作效率最高？
用效率矩阵表示：
其中，元素cij≥0（i,j=1,…,n）表示指派第i人去完成第j项任务时的效率（或时间、成本等）。
令
则指派问题的一般数学模型为：
对于上述问题，其实通俗来说，就是n×n矩阵中，选取n个元素，每行每列各1个元素，使得和最小。
根据指派问题的特点，可以使用匈牙利算法来进行求解。
匈牙利算法原理（匈牙利数学家狄·康尼格证明的两个定理）：
定理1如果从指派问题效率矩阵[cij]的每一行元素中分别减去（或加上）一个常数ui（称为该行的位势），从每一列分别减去（或加上）一个常数vj（称为该列的位势），得到一个新的效率矩阵[bij],若其中bij=cij-ui-vj,则[bij]的最优解的结果等价于[cij]的最优解的结构。
因此，指派问题的最优解具有性质：
若从矩阵的一行（列）各元素中分别减去该行（列）的最小元素，得到规约矩阵，其最优解和原矩阵的最优解相同。
经过减去行列最小元素的转化之后的效率矩阵[bij]中，把位于不同行不同列的0元素，简称为独立0元素。
目的便是找出位于不同行不同列的n个独立的0元素，从而也就得到了原问题的最优解。
定理2若效率矩阵C的元素看可分成“0”和“非0”两部分，则覆盖所有“0”元素的最少直线数=独立的“0”元素的最多个数。
例3.设指派问题的效率矩阵为
Lingo代码：
model:
sets:
var/1..5/;
link（var,var）:
c,x;
endsets
data:
c=382103
87297
64275
84235
9106910;
enddata
min=@sum（link:
c*x）;
@for（var（i）:
@sum（var（j）:
x（i,j））=1）;
@for（var（j）:
@sum（var（i）:
x（i,j））=1）;
@for（link:
@bin（x））;
end
运行结果（部分）：
Globaloptimalsolutionfound.
Objectivevalue:
21.00000
Objectivebound:
21.00000
Infeasibilities:
0.000000
Extendedsolversteps:
0
Totalsolveriterations:
0
VariableValueReducedCost
X（1,1）0.0000003.000000
X（1,2）0.0000008.000000
X（1,3）0.0000002.000000
X（1,4）0.00000010.00000
X（1,5）1.0000003.000000
X（2,1）0.0000008.000000
X（2,2）0.0000007.000000
X（2,3）1.0000002.000000
X（2,4）0.0000009.000000
X（2,5）0.0000007.000000
X（3,1）0.0000006.000000
X（3,2）1.0000004.000000
X（3,3）0.0000002.000000
X（3,4）0.0000007.000000
X（3,5）0.0000005.000000
X（4,1）0.0000008.000000
X（4,2）0.0000004.000000
X（4,3）0.0000002.000000
X（4,4）1.0000003.000000
X（4,5）0.0000005.000000
X（5,1）1.0000009.000000
X（5,2）0.00000010.00000
X（5,3）0.0000006.000000
X（5,4）0.0000009.000000
X（5,5）0.00000010.00000
RowSlackorSurplusDualPrice
121.00000-1.000000
20.0000000.000000
三、蒙特卡罗法（随机取样法）
前面的方法，主要是针对线性整数规划而言，对于非线性整数规划没有通用的有效解法。
整数规划由于限制变量是整数，增加了求解难度，但整数解是有限个，所以可以采用枚举法。
当枚举个数很多时，显性枚举是不现实的，但利用蒙特卡罗随机取样法，在一定的计算量下是可以得到满意解的。
例4.求解如下非线性整数规划问题：
若用显枚举法，共需1005=1010个点，计算量非常大。
但用蒙特卡罗法随机计算106个点，便可找到满意解。
代码：
目标函数：
function[f,g]=fun1（x）
f=x
（1）^2+x
（2）^2+3*x（3）^2+4*x（4）^2+2*x（5）^2-8*x
（1）-2*x
（2）-3*x（3）-x（4）-2*x（5）;
g
（1）=sum（x）-400;
g
（2）=x
（1）+2*x
（2）+2*x（3）+x（4）+6*x（5）-800;
g（3）=2*x
（1）+x
（2）+6*x（3）-200;
g（4）=x（3）+x（4）+5*x（5）-200;
主程序：
rand（'state',sum（clock））;
p0=0;
tic
fori=1:
10^6
x=99*rand（5,1）;
x1=floor（x）;
x2=ceil（x）;
[f,g]=fun1（x1）;
ifsum（g<=0）==4
ifp0<=f
x0=x1;
p0=f;
end
end
[f,g]=fun1（x2）;
ifsum（g<=0）==4
ifp0<=f
x0=x2;
p0=f;
end
end
end
x0,p0
toc
运行结果：
x0=43941995
p0=49298
Elapsedtimeis45.494952seconds.
注意：
由于是随机取样，故每次的运行结果可能不一样。
利用lingo软件可以求得全局最优解：
代码：
model:
sets:
row/1..4/:
b;
col/1..5/:
c1,c2,x;
link（row,col）:
a;
endsets
data:
c1=1,1,3,4,2;
c2=-8,-2,-3,-1,-2;
a=11111
12216
21600
00115;
b=400,800,200,200;
enddata
max=@sum（col:
c1*x^2+c2*x）;
@for（row（i）:
@sum（col（j）:
a（i,j）*x（j））@for（col:
@gin（x））;
@for（col:
@bnd（0,x,99））;
End
运行结果（部分）：
Localoptimalsolutionfound.
Objectivevalue:
51568.00
Objectivebound:
51568.00
Infeasibilities:
0.000000
Extendedsolversteps:
0
Totalsolveriterations:
51
VariableValueReducedCost
X
（1）50.00000-92.00000
X
（2）99.00000-196.0000
X（3）0.0000003.000000
X（4）99.00000-791.0000
X（5）20.00000-78.00000
RowSlackorSurplusDualPrice
151568.001.000000
2132.00000.000000
3333.00000.000000
41.0000000.000000
51.0000000.000000
四、混合整数规划（配送选址问题）
图1三级配送的网络结构
物流配送中心选址问题，是在给定某一地区所有备选点的地址集合中选出一定数目的地址建立配送中心，从而建立一系列的配送区域，以实现选出点建立的配送中心与各需求点和工厂（供货点）形成的配送系统总物流费用最小。
总费用=运输费用+配送费用+仓储费用+固定成本
设有m个供应工厂，n个客户，l个备选配送中心；
pk表示第k个工厂的产量（k=1,…,m）;
dj表示第j个客户的需求量（j=1,…,n）;
gi表示第i个配送中心的单位管理成本（i=1,…,l）;
fi表示建立第i个配送中心的固定成本；
zi指示是否选中第i个配送中心，1表示选中，0表示未选中；
zmax表示配送中心的最大限制数目；
cki表示第k个供应工厂到第i个配送中心的运输单价；
hij表示第i个配送中心到第j个客户的配送单价；
wki表示第k个供应工厂到第i个配送中心的运量；
xij表示第i个配送中心到第j个客户的运量。
则配送选址问题数学模型的一般形式如下：
例5设有4个备选物流配送中心地址，6个工厂为其供货，6个客户需要产品，最多设置3个物流配送中心。
工厂到物流配送中心的运输价格见表1，物流配送中心到客户的运输价格见表2，工厂的总生产能力见表3，物流配送中心的固定成本、单位管理成本，及容量见表4，客户的需求量见表5：
表1工厂到配送中心的运输价格
w1
w2
w3
w4
p1
6
5
4
2
p2
2
3
4
9
p3
6
8
7
5
p4
7
4
2
3
p5
4
2
5
1
p6
3
4
1
7
表2配送中心到客户的运输价格
c1
c2
c3
c4
c5
c6
w1
3
2
7
4
7
5
w2
6
1
4
2
5
3
w3
2
4
5
3
6
8
w4
5
6
3
7
4
6
表3工厂的总生产能力
工厂
p1
p2
p3
p4
p5
p6
总生产能力（p）
40,000
50,000
60,000
70,000
60,000
40,000
表4备选物流配送中心的固定成本，单位管理成本，容量
物流配送中心
w1
w2
w3
w4
固定成本（f）
500,000
300,000
400,000
400,000
单位管理成本（g）
3
2
5
4
仓库容量（a）
10,000
60,000
70,000
50,000
表5客户的需求量
顾客
c1
c2
c3
c4
c5
c6
需求（d）
10,000
20,000
10,000
20,000
30,000
10,000
利用Lingo软件求解以上混合整数规划。
代码：
model:
sets:
factory/p1..p6/:
p;
warhouse/w1..w4/:
a,f,g;
customer/c1..c6/:
d;
tr/tr1..tr4/:
z;
link1（factory,warhouse）:
c,w;
link2（warhouse,customer）:
h,x;
endsets
data:
p=40000,50000,60000,70000,60000,40000;
a=70000,60000,70000,50000;
f=500000,300000,400000,400000;
g=3,2,5,4;
d=10000,20000,10000,20000,30000,10000;
c=6542
2349
6875
7423
4251
3417;
H=327475
614253
245368
563746;
enddata
min=@sum（link1（k,i）:
c（k,i）*w（k,i））+@sum（link2（i,j）:
h（i,j）*x（i,j））
+@sum（link1（k,i）:
g（i）*w（k,i））+@sum（warhouse（i）:
f（i）*z（i））;
@for（factory（k）:
@sum（link1（k,i）:
w（k,i））<=p（k））;
@for（warhouse（i）:
@sum（link2（i,j）:
x（i,j））=@sum（link1（k,i）:
w（k,i）））;
@for（customer（j）:
@sum（link2（i,j）:
x（i,j））>=d（j））;
@for（warhouse（i）:
@sum（link1（k,i）:
w（k,i））<=（a（i）*z（i）））;
@sum（tr（i）:
z（i））<=3;
@for（tr（i）:
@bin（z））;
End
运行结果（部分）：
Globaloptimalsolutionfound.
Objectivevalue:
1480000.
Objectivebound:
1480000.
Infeasibilities:
0.000000
Extendedsolversteps:
4
Totalsolveriterations:
50
VariableValueReducedCost
Z（TR1）0.000000290000.0
Z（TR2）1.000000300000.0
Z（TR3）0.000000190000.0
Z（TR4）1.000000400000.0
C（P1,W1）6.0000000.000000
C（P1,W2）5.0000000.000000
C（P1,W3）4.0000000.000000
C（P1,W4）2.0000000.000000
C（P2,W1）2.0000000.000000
C（P2,W2）3.0000000.000000
C（P2,W3）4.0000000.000000
C（P2,W4）9.0000000.000000
C（P3,W1）6.0000000.000000
C（P3,W2）8.0000000.000000
C（P3,W3）7.0000000.000000
C（P3,W4）5.0000000.000000
C（P4,W1）7.0000000.000000
C（P4,W2）4.0000000.000000
C（P4,W3）2.0000000.000000
C（P4,W4）3.0000000.000000
C（P5,W1）4.0000000.000000
C（P5,W2）2.0000000.000000
C（P5,W3）5.0000000.000000
C（P5,W4）1.0000000.000000
C（P6,W1）3.0000000.000000
C（P6,W2）4.0000000.000000
C（P6,W3）1.0000000.000000
C（P6,W4）7.0000000.000000
W（P1,W1）0.0000004.000000
W（P1,W2）0.0000002.000000
W（P1,W3）0.0000003.000000
W（P1,W4）40000.000.000000
W（P2,W1）0.0000000.000000
W（P2,W2）0.0000000.000000
W（P2,W3）0.0000003.000000
W（P2,W4）0.0000007.000000
W（P3,W1）0.0000004.000000
W（P3,W2）0.0000005.000000
W（P3,W3）0.0000006.000000
W（P3,W4）0.0000003.000000
W（P4,W1）0.0000005.000000
W（P4,W2）0.0000001.000000
W（P4,W3）0.0000001.000000
W（P4,W4）0.0000001.000000
W（P5,W1）0.0000003.000000
W（P5,W2）60000.000.000000
W（P5,W3）0.0000005.000000
W（P5,W4）0.0000000.000000
W（P6,W1）0.0000001.000000
W（P6,W2）0.0000001.000000
W（P6,W3）0.0000000.000000
W（P6,W4）0.0000005.000000
H（W1,C1）3.0000000.000000
H（W1,C2）2.0000000.000000
H（W1,C3）7.0000000.000000
H（W1,C4）4.0000000.000000
H（W1,C5）7.0000000.000000
H（W1,C6）5.0000000.000000
H（W2,C1）6.0000000.000000
H（W2,C2）1.0000000.000000
H（W2,C3）4.0000000.000000
H（W2,C4）2.0000000.000000
H（W2,C5）5.0000000.000000
H（W2,C6）3.0000000.000000
H（W3,C1）2.0000000.000000
升级会员 