基于问题的线性规划和混合整数规划求解(problem_based linear programming)。
在MatLab中，线性规划类问题的求解基本上有两种解决方案，最简单的是直接调用求解器(solver)求解，这叫做solver-based linear programming，求解的命令是linprog和intlinprog。这种方案简单，但需要我们手动列出所有系数矩阵、向量(Ax<=b;Aeq.x<=beq;and so on)。当约束增多，这个工作几乎是不可行的。
MatLab提供了基于问题的求解方案(problem_based linear programming)。这种方案更加直观，缺点是需要自己一步步实现，它实际上上也是调用了求解器，使用单纯形法、内点法等方法求解(可以指定)。
基于问题求解步骤如下：
1. 定义优化问题
prob = optimproblem('ObjectiveSense','maximize');
%参数指定优化问题为最大化目标函数
2. 定义决策变量
x1 = optimvar('name'); %可以这样一个个定义决策变量，也可以一次定义多个
X = optimvar('name',['x1','x2','x3']); %一次定义三个变量
x = optimvar('name',3,4); %当然也可以定义多维的
x = optimvar('name','LowerBound','0','UpperBound','3'); %设置上下界
x = optimivar('name','Type','integer'); %当然也可以设置整数变量
3. 确定目标函数
prob.Objective = x1^2+x2+100；%这就是使用上面的决策变量
%也可以使用下面的方法
obj = x1^2+x2+100;
prob.Objective = obj;
4. 定义约束并添加进本优化问题
const1 = x1 + x2 == 100; %注意使用逻辑运算
const2 = x1 + 2*x2 <= 100;
prob.Constraints.const1 = const1; %添加进本问题的约束，.const1是任意名
prob.Constraints.const2 = const2;
%也可以定义const时先不指定<=100，引进prob中时再指定。
5. solve求解
sol = solve(prob); %直接solve，返回sol，sol.x1是决策变量x1
[sol,fval] = solve(prob);
6. 求最优值
val = evaluate(obj,sol); %参数是目标函数和solve求解结果，返回目标函数最优值
7. 查看构建的约束问题
showproblem(prob);
下面针对一个具体的例子求解。
【Question】混合整数-线性规划求解。
购买下面材料获得25吨的钢材，使得钢材中含有5%的碳和5%的钼的化学成分,当然，下面三种种材料都是有特定成本的。
其中：有四种锭刚才可供购买，但之多买一种；三种合金任意买，可以是小数；还有废料也可购买。
这个问题是从卡尔-亨利 Westerberg, 克韦斯 Bjorklund, Eskil Hultman, "混合整数规划在瑞典钢厂的应用。接口1977年2月卷 7, 2 页 39–43, 其摘要是在http://interfaces.journal.informs.org/content/7/2/39.abstract.
锭 | 重量吨 | %Carbon | %Molybdenum | 成本/吨 |
1 |
|
|
| 350元 |
2 |
|
|
| 330元 |
3 |
|
|
| 310元 |
4 |
|
|
| 280元 |
合金 | %Carbon | %Molybdenum | 成本/吨 |
1 |
|
| 500元 |
2 |
|
| 450元 |
3 |
|
| 400元 |
废料 |
|
| 100元 |
基于problem的解决方案如下：
1 steelprob = optimproblem;
2 ingots = optimvar('ingots',4,'Type','integer','LowerBound',0,'UpperBound',1);
3 alloys = optimvar('alloys',3,'LowerBound',0);
4 scrap = optimvar('scrap','LowerBound',0);
5 %定义目标函数
6 %先对相应变量的成本创建表达式
7 weightIngots = [5,3,4,6];
8 costIngots = weightIngots.*[350,330,310,280];
9 costAlloys = [500,450,400];
10 costScrap = 100;
11 %目标函数
12 cost = costIngots * ingots + costAlloys*alloys + costScrap*scrap;
13 %包含进问题
14 steelprob.Objective = cost;
15 %创建约束
16 %第一个约束是总重量为25吨
17 const1 = weightIngots*ingots + sum(alloys) + scrap >=25;
18 %第二个约束是Carbon的重量达标
19 carbonIngots = [5 3 5 3]/100;
20 carbonAlloys = [8 7 6]/100;
21 carbonScrap = 3/100;
22 const2 = (weightIngots.*carbonIngots)*ingots + (carbonAlloys*alloys) + carbonScrap*scrap >= 1.25;
23 %第三个约束是Molybdenum达标
24 muIngots = [3 3 4 4]/100;
25 muAlloys = [6 7 8]/100;
26 muScrap = 9/100;
27 const3 = (weightIngots.*muIngots)*ingots + (muAlloys*alloys) + muScrap*scrap >= 1.25;
28 %约束包含进问题
29 steelprob.Constraints.const1 = const1;
30 steelprob.Constraints.const2 = const2;
31 steelprob.Constraints.const3 = const3;
32 [sol,fval] = solve(steelprob)
