本文为转载文章,原文地址:https://blog.csdn.net/wangzhenyang2/article/details/109203962
(本文内容主要参考《an introduction to management science(13e)》中的章节)
分支界定算法的主要思想是把问题的可行解划分为多个子集,对子集进行系统化的评估直到最优解找到。分支界定只是一种算法框架,因为如何构造子集、在每个子集如何搜索最优解可以按不同情况有不同的具体实现,当用于求解线性整数优化问题时,就需要结合一般的线性优化算法或者求解器来实现。下面先用一个实例来说明算法的思路,然后用代码实现一下。
一家工厂打算购买两种设备,印刷机和车床。老板估计每台印刷机每台可以增加$100的利润,每台车床每天可以增加$150的利润;但同时购买数量也会收到预算和场地面积的约束,每台印刷机或车床的花费和占地如下表所示:
Machine | Required Floor Space | Purchase Price |
---|---|---|
印刷机 | 15 | $8000 |
车床 | 30 | $4000 |
老板的预算是$40000,并且现在有200平方的空地;工厂老板应该怎样购买来让每天的利润最大?
如果我们设x 1 x_1x1是购买多少台印刷机,x 2 x_2x2是购买多少台车床,那么可以写出优化问题的模型
相比线性优化模型,上面的数学模型就多出必须为整数的约束。如果我们把整数约束去掉,那么可以得到最优解:
如果我们把这个解四舍五入,可以得到一个可行的整数解
但是这个整数解并不是实际的最优解,因为一旦加上整数约束,线性优化模型往往是非凸的,线性优化算法往往无能为力。而分支界定法则是把原问题进行松弛(即去掉整数约束),以松弛线性模型为根节点不断地划分子问题和更新问题上下界,通过评判上下界来忽略无效子问题和停止算法。首先我们构造根节点,该节点就是上面模型去掉整数约束后的松弛模型,对于这个松弛线性最大化问题,我们构造一下上下界,用它的最优解(Z = 1055.56)作为问题的上界,用四舍五入后的整数解作为下界;上下界的意义在于,原始整数优化模型的解一定是在上下界之间的,所以对于所有最优解比下界更小的子问题,我们可以直接忽略;而子问题的上界一定不会大于父问题,通过不断分支,这个上下界间的范围将不断缩减。
因为根节点的最优解不是整数解,所以我们选择一个非整数解进行分支。这里我们采取选择离四舍五入值最远的变量,对于根节点,即选择了变量x 2 = 5.56。我们以变量x 2 为划分点向下分支出两个子问题,分别在根节点的松弛优化模型的基础上新增约束x 2 ≤ 5和x 2 ≥ 6 ,也就是说,对于左侧的子节点2,松弛优化模型为:
右侧子节点3的问题模型为:
对于节点2,可以求得最优解为x 1 = 2.5 , x 2 = 5 , Z = 1000 ,而节点3的最优解为x 1 = 1.33 , x 2 = 6 , Z = 1033.33 ,这两个解可分别作为这两个节点的上界。可以看到这两个上界都小于根节点的上界,因为每一次分支,都是在原松弛优化模型的基础上新增约束,它的最优解肯定不会好过原问题。再来分析下界,因为这两个子节点的解都不是完整整数解,所以这两个节点下界保持不变,仍是950
我们还没有找到最优可行解,所以要从结点2或者3继续分支。选择哪个结点呢,这里我们按上界更大优先的原则来选取结点,也就是从结点3进行分支。采用这种原则其实是一种间接地剪枝,如果某个结点的最优解是非整数解并且小于当前的全局下界(即上界比下界更小),那么需要把这个结点及其所有子节点都从搜索数里去除,它们的解肯定不会是最优解了,通过上界越大越优先的原则来选取,这些结点实际上根本不会被选中。结点3的x 2已经是整数了,所以用x 1进行划分:x 1 ≤ 1和x 1 ≥ 2 ,分出了结点4和结点5,结点4的问题模型是在3的基础上加上x 1 ≤ 1
结点5则是在3的基础上加上x 1 ≥ 2
对于结点4,其最优解是(x 1 = 1 , x 2 = 6.17 , Z = 1025 );对于结点5,则无解。那就只需要对结点4进行定界,上界就是最优解1025.5,同时最优解也不是完整整数解,那下界仍保持不变。
我们还是没有找到最优整数解,需要继续分支。现在有3个叶节点:2、4和5,结点5已经无解,所以不需要考虑,而在2和4中,4的上界更大,所以从结点4开始分支。在结点4,x 1已经是整数,所以用x 2 进行划分:x 2 ≤ 6 和x 2 ≥ 7 。结点6的问题模型如下,其最优解为(x 1 = 1 , x 2 = 6 , Z = 1000)
结点7的问题模型如下,它同样无解
现在搜索树如下图所示
在结点6,因为出现了新的整数解,并且该整数解的目标值大于之前的下界950,所以结点6的下界更新为1000,这时结点6的上下界已经相等。现在可比较叶节点是2和6(排除掉无解的5和7),可以发现结点6是完整的整数解并且上界是叶节点中最大的,那就说明结点6的解就是原始整数规划的最优解。为什么呢,结点6的上界比其他叶节点的上界大,而从叶节点继续向下分支的上界总是比叶节点的小,所以结点6的上界肯定是比其他所有子节点的上界更大。由此我们找到了原问题的解,即购买1台印刷机和6台车床可以达到最大日利润$1000
我们来归纳一下分支界定算法的一般步骤。我注意到不同的文章里对算法步骤的描述都有所区别,不过核心的几个要素是一致的。
- Step 1: 计算原始整数规划松弛化后的最优值,如果无解,则原问题就是无解,退出算法;如果最优解就是整数解,那么直接得到最终解,退出算法
- Step 2: 建立根节点,在根节点上,上界是原始整数规划松弛后的最优解,下界则需要预先给出一个可行解来得到
- Step 3: 评估当前所有的叶节点,不过不需要加入无解的结点。如果叶节点某个结点的最优解是完整整数解,并且其上界比其他叶节点都要大,则该结点的解就是全局最优解,可以退出算法;否则继续Step 4;另外如果目前没有任何可行的叶节点,那么就表明除了初始给出的可行解,不存在任何其他可行解,也直接退出算法
- Step 4: 从最大上界的可行叶节点继续分支。选择该结点的一个非整数解,产生两个新的≤ \leq≤和≥ \geq≥约束,分别加到该节点的问题模型上产生两个子节点。计算更新两个子节点的上下界,上界值等于子节点的最优解,下界等于全局的最大下界。
- Step 5: 回到Step 3
最后我用C#代码实现一下上面问题求解。首先需要说明的是,求解松弛问题,我使用了Google OR-Tools的线性求解器。定义了一个LooseLinearConstraintProblem类,来表示基于上面例子的松弛优化问题模型,在内部定义了最基本的约束,并且定义了additionalConstraints来表示相对父问题的额外约束,如果要从一个父松弛问题生成一个子松弛问题,可以通过拷贝构造函数来生成
public class LooseLinearConstraintProblem
{
private Solver solver;
private List<Variable> variables;
public List<Variable> Variables
{
get
{
return variables;
}
}
//private List<Constraint> constraints;
private Objective objective;
//item1: 第几个变量
//item2: 大于等于多少
//item3: 小于等于多少
private List<Tuple<int, double, double>> additionalConstraints;
public List<Tuple<int, double, double>> AdditionalConstraints
{
set
{
this.additionalConstraints = value;
}
get
{
return this.additionalConstraints;
}
}
private void BasicInit()
{
this.solver = Solver.CreateSolver("GLOP");
this.variables = new List<Variable>();
Variable x1 = solver.MakeNumVar(0.0, double.PositiveInfinity, "x1");
Variable x2 = solver.MakeNumVar(0.0, double.PositiveInfinity, "x2");
this.variables.Add(x1);
this.variables.Add(x2);
Constraint c0 = solver.MakeConstraint(double.NegativeInfinity, 40000);
c0.SetCoefficient(x1, 8000);
c0.SetCoefficient(x2, 4000);
Constraint c1 = solver.MakeConstraint(double.NegativeInfinity, 200);
c1.SetCoefficient(x1, 15);
c1.SetCoefficient(x2, 30);
this.objective = solver.Objective();
this.objective.SetCoefficient(x1, 100);
this.objective.SetCoefficient(x2, 150);
this.objective.SetMaximization();
}
public LooseLinearConstraintProblem(List<Tuple<int, double, double>> additionalConstraints)
{
BasicInit();
if (additionalConstraints != null)
{
this.additionalConstraints = additionalConstraints;
foreach(var additionalConstraint in additionalConstraints)
{
Constraint newConstraint = this.solver.MakeConstraint(additionalConstraint.Item2, additionalConstraint.Item3);
for (int j = 0; j < this.variables.Count; j++)
{
if (j != additionalConstraint.Item1)
{
newConstraint.SetCoefficient(this.variables[j], 0);
}
else
{
newConstraint.SetCoefficient(this.variables[j], 1);
}
}
}
}
}
public LooseLinearConstraintProblem(LooseLinearConstraintProblem parentProblem,List<Tuple<int, double, double>> additionalConstraints)
{
BasicInit();
this.additionalConstraints = new List<Tuple<int, double, double>>();
if (parentProblem.AdditionalConstraints != null)
{
foreach (var additionalConstraint in parentProblem.AdditionalConstraints)
{
this.additionalConstraints.Add(additionalConstraint);
Constraint newConstraint = this.solver.MakeConstraint(additionalConstraint.Item2, additionalConstraint.Item3);
for (int j = 0; j < this.variables.Count; j++)
{
if (j != additionalConstraint.Item1)
{
newConstraint.SetCoefficient(this.variables[j], 0);
}
else
{
newConstraint.SetCoefficient(this.variables[j], 1);
}
}
}
}
if (additionalConstraints != null)
{
foreach (var additionalConstraint in additionalConstraints)
{
this.additionalConstraints.Add(additionalConstraint);
Constraint newConstraint = this.solver.MakeConstraint(additionalConstraint.Item2, additionalConstraint.Item3);
for (int j = 0; j < this.variables.Count; j++)
{
if (j != additionalConstraint.Item1)
{
newConstraint.SetCoefficient(this.variables[j], 0);
}
else
{
newConstraint.SetCoefficient(this.variables[j], 1);
}
}
}
}
}
public Tuple<List<double>, double> Solve()
{
var status = solver.Solve();
if (status != Solver.ResultStatus.OPTIMAL)
{
return null;
}
else
{
List<double> vResults = new List<double>();
foreach(var v in this.variables)
{
vResults.Add(v.SolutionValue());
}
double oResult = solver.Objective().Value();
return new Tuple<List<double>, double>(vResults, oResult);
}
}
}