线性演算(二)
线性演算求解
1 分支定界法(The Branch and Bound Method)
分支定界法主要用于解决整数线性规划问题(integer linear programs). 就像 Simplex 算法的例子一样, 分支定界法使用来解决优化问题的, 但是我们还是讨论其在判定问题中的应用.
定义:松弛系统(relaxed system) 给定一个线性规划系统(integer linear system) \(S\), 该系统没有整数的要求, 也就是说, 变量没有要求必须是整型.
用 relaxed(\(S\)) 来表示系统 \(S\) 的松弛问题. 假定存在一个程序 \(LP_{feasible}\), 该程序接收一个线性规划系统作为输入, 如果 \(S\) 是可满足的则返回"unsat", 否则返回一组可满足的赋值. 可以通过通用 Simplex 算法的变式来实现 \(LP_{feasible}\). 如下算法, 判定一个整型线性约束规划系统的可满足性问题. \[ \begin{aligned} &\mathbf{Algorithm\ 1:\ }{\rm FEASIBILITY\mbox{-}BRANCH\mbox{-}AND\mbox{-}BOUND}\\ &\mathbf{Input:}\ {\rm An\ integer\ linear\ system}\ S\\ &\mathbf{Output:}\ {\rm ''Satisfiable''\ if}\ S\ {\rm is\ satisfiable,\ and\ ''Unsatisfiable''\ otherwise}\\ \\ &1.\hspace{2mm}\mathbf{procedure}\ {\rm SEARCH\mbox{-}INTEGRAL\mbox{-}SOLUTION}(S) \\ &2.\hspace{10mm}{\rm res}=LP_{feasible}({\rm relaxed}(S)); \\ &3.\hspace{10mm}\mathbf{if\ }{\rm res\ =\ ''Unsatisfiable''}\ \mathbf{then\ return};\hspace{30mm}\triangleright\ {\rm prune\ branch}\\ &4.\hspace{10mm}\mathbf{else\ }\\ &5.\hspace{18mm}\mathbf{if}\ {\rm res\ is\ integral\ }\mathbf{then}\hspace{40mm}\triangleright\ {\rm integer\ solution\ found}\\ &\hspace{30mm}\mathbf{abort\ }{\rm (''Satisfiable'')}; \\ &6.\hspace{18mm}\mathbf{else\ }\\ &7.\hspace{26mm} {\rm Select\ a\ variable\ }v{\rm\ that\ is\ assigned\ a\ nonintegral\ value\ }r;\\ &8.\hspace{26mm} {\rm SEARCH\mbox{-}INTEGRAL\mbox{-}SOLUTION}(S\ \cup\ (v\le\lfloor r\rfloor));\\ &9.\hspace{26mm} {\rm SEARCH\mbox{-}INTEGRAL\mbox{-}SOLUTION}(S\ \cup\ (v\ge\lceil r\rceil));\\ &10.\hspace{78mm}\triangleright\ {\rm no\ integer\ solution\ in\ this\ branch} \\ \\ \\ \\ &11.\hspace{2mm}\mathbf{procedure}\ {\rm FEASIBILITY\mbox{-}BRANCH\mbox{-}AND\mbox{-}BOUND}(S) \\ &12.\hspace{8mm}{\rm SEARCH\mbox{-}INTEGRAL\mbox{-}SOLUTION}(S); \\ &13.\hspace{8mm}\mathbf{return\ }{\rm(''Unsatisfiable'')};\\ \end{aligned} \] 该算法通过程序 \(LP_{feasible}\) 来解决松弛问题; 如果该问题是不可满足的, 那么就回溯, 因为在该分支下没有整数解. 相反, 如果该问题是可满足的, 同时 \(LP_{feasible}\) 返回了一组整数解, 那么该问题终止. 否则, 该问题分解为两个子问题, 并递归调用自身. 下面用个例子说明问题的分解是如果处理的.
例1
设 \(x_1,...,x_4\) 是系统 \(S\) 的变量. 假定 \(LP_{feasible}\) 返回一组解 \[ (1,0.7,2.5,3)\tag{1} \] SEARCH-INTEGRAL-SOLUTION方法在变量 \(x_2\) 和 \(x_3\) 之间进行选择, \(x_2,x_3\) 是非整数赋值的变量. 假定选择了变量 \(x_2\). 在算法的第 8 行, 然后, 将 \(S\) (在当前递归层级求解的线性系统)增加约束, 并发送给更深的递归层级求解. \[ x_2\le0\tag{2} \] 如果在当前分支下, 没有合适的解, 那么对系统 \(S\) 将使用另外的约束条件, 并同时进一步递归求解. \[ x_2\ge1\tag{3} \] 如果这两次调用都返回了, 就说明 \(S\) 没有可满足的解, 判定程序就回溯. 应当注意, 从初始的递归层回溯, 将会导致FEASIBILITY-BRANCH-AND-BOUND算法返回"Unsatisfiable", 即不可满足.
上述算法存在一个问题, 该算法不是完备的(complete): 会有让其一直分支下去的例子, 即进入递归的死循环. 比如说线性系统 \(1\le 3x-3y\le2\), 没有整数解, 但是也没有有界实数解, 所以该系统会让算法一直循环下去. 为了算法的完备性, 必须依靠这些公式具有的小模型属性. 也就是说, 如果存在可满足的解, 那么在有限域(finite bound)中, 也应该存在这样的解, 对于该理论, 这是可计算的. 因此, 一旦我们在每个变量的域上计算了此界限, 就可以在传值之后停止搜索. 同样的界限同样可以应用于可行性问题. 简要地说, 给定一个整数线性规划系统 \(S\), \(M\times N\) 的系数矩阵 \(A\), 如果有一组解, 那么系统 \(S\) 的凸包的极点也是一组解, 并且任一解 \(x^0\) 有如下约束: \[ x_j^0\le((M+N)\cdot N\cdot \theta)^N\ {\rm for}\ j = 1,...,N, \] 其中, \(\theta\) 是该问题中的最大元素. 因此上式给定每一变量设定了界限, 也就是给每一个变量添加了显式的约束, 强制算法终止.
其实, 分支定界法可以直接扩展成一种解决一些变量是整型而其他的是实型的系统. 在优化问题中, 这类问题就是混合线性规划问题(mixed integer programming).
2 切平面(Cutting Planes)
切平面是添加在线性系统中的约束, 用以移除非整数解; 如果系统中的解所有都是可满足的整数解, 那么维持这种可满足性. 这些新的约束提高了求解整数线性规划系统过程中松弛的紧密度, 因此可以使分支定界算法运行得更快(这种组合称作分支剪界, branch-and-cut). 此外, 如果满足某些条件, 以下所述类型的 Simplex 和切平面算法将形成整数线性演算的判定过程.
我们讨论切平面的一类, 称作 Gomory 剪切(Gomory cuts). 先说明例子, 然后推广.
例2
假定问题包含整型变量 \(x_1,x_2,x_3\), 有下界 \(1\le x_1\) 和 \(0.5\le x_2\). 此外, 假定通用 Simplex 算法的最终的表(tableau)包含下面的约束 \[ x_3 = 0.5x_1+2.5x_2, \] 并且有一组解 \(\alpha:\{x_3\mapsto1.75, x_1\mapsto1, x_2\mapsto0.5\}\), 显然, 该赋值是满足上面的式子的. 每个变量减去该赋值得到 \[ x_3-1.75=0.5(x_1-1)+2.5(x_2-0.5). \] 现在将上面式子改写成坐式是整数的形式: \[ x_3-1=0.75+0.5(x_1-1)+2.5(x_2-0.5). \] 最右边的两项必须是正值, 因为 1 和 0.5 分别是 \(x_1\) 和 \(x_2\) 的下界. 因此, 式子的右边累加起来也必须是一个整数, 也就是说 \[ 0.75+0.5(x_1-1)+2.5(x_2-0.5)\ge1. \] 但是, 该约束在赋值 \(\alpha\) 下是不可满足的, 因为通过构造, 除分数 0.75 外, 左侧的所有元素在 \(\alpha\) 下均等于零. 也就是说, 将此约束加入到松弛系统中将会把该解剔除掉. 从另一个角度来说, 它不能移除任何整数解.
让我们将此示例概括为生成此类切割平面的方法. 泛化还涉及为变量分配上限, 负系数和正系数的情况. 为了从约束中导出Gomory 剪切, 该约束必须满足两个条件:首先, 对于基础变量的赋值要是分数; 其次, 对于所有非基础变量的赋值要对应于它们其中一个有限域值. 以下用例子说明,
考虑第 \(i\) 个约束: \[ x_i=\Sigma_{x_j\in\mathcal{N}}a_{ij}x_j,\tag{4} \] 其中 \(x_i\in \mathcal{B}\). 令 \(\alpha\) 是由通用 Simplex 算法返回的赋值, 因此, \[ \alpha(x_i)=\Sigma_{x_j\in\mathcal{N}}a_{ij}\alpha(x_j).\tag{5} \] 现在将非基本变量划分为当前已分配其下限的变量和当前已分配其上限的变量两部分: \[ J=\{j\ |\ x_j\in\mathcal{N}\land \alpha(x_j)=l_j\},\\ K=\{j\ |\ x_j\in\mathcal{N}\land \alpha(x_j)=u_j\}. \] (4) 式子减 (5) 式子得 \[ x_i-\alpha(x_i)=\Sigma_{j\in J}a_{ij}(x_j-l_j)-\Sigma_{j\in K}a_{ij}(u_j-x_j).\tag{6} \] 令 \(f_0=\alpha(x_i)-\lfloor\alpha(x_i)\rfloor\). 因为我们假定 \(\alpha(x_i)\) 不是整数, 所以 \(0\lt f_0 \le1\). 可以将 (6) 重写成 \[ x_i-\lfloor\alpha(x_i)\rfloor=f_0+\Sigma_{j\in J}a_{ij}(x_j-l_j)-\Sigma_{j\in K}a_{ij}(u_j-x_j). \] 注意等式左边是一个整数. 现在考虑两种情形.
如果 \(\Sigma_{j\in J}a_{ij}(x_j-l_j)-\Sigma_{j\in K}a_{ij}(u_j-x_j)\gt0\), 又因等式左边要是整型, 所以, \[ f_0+\Sigma_{j\in J}a_{ij}(x_j-l_j)-\Sigma_{j\in K}a_{ij}(u_j-x_j)\ge1.\tag{7} \] 现在将 \(J\) 和 \(K\) 分别分解成如下两个部分: \[ \begin{aligned} &J^+=\{j\ |\ j\in J\land a_{ij}\gt0\},\\ &J^-=\{j\ |\ j\in J\land a_{ij}\lt0\},\\ &K^+=\{j\ |\ j\in K\land a_{ij}\gt0\},\\ &K^-=\{j\ |\ j\in K\land a_{ij}\lt0\}.\\ \end{aligned} \] 仅收集 (7) 中不等式左边的正元素得 \[ \Sigma_{j\in J^+}a_{ij}(x_j-l_j)-\Sigma_{j\in K^-}a_{ij}(u_j-x_j)\ge1-f_0, \] 或者, 等价地, \[ \Sigma_{j\in J^+}\frac{a_{ij}}{1-f_0}(x_j-l_j)-\Sigma_{j\in K^-}\frac{a_{ij}}{1-f_0}(u_j-x_j)\ge1 \tag{8} \]
如果 \(\Sigma_{j\in J}a_{ij}(x_j-l_j)-\Sigma_{j\in K}a_{ij}(u_j-x_j)\le0\), 又因等式左边要是整型, 所以, \[ f_0+\Sigma_{j\in J}a_{ij}(x_j-l_j)-\Sigma_{j\in K}a_{ij}(u_j-x_j)\le0.\tag{9} \] 由式子 (9) 得出 \[ \Sigma_{j\in J^-}a_{ij}(x_j-l_j)-\Sigma_{j\in K^+}a_{ij}(u_j-x_j)\le-f_0. \] 两边同除以 \(-f_0\) 得 \[ -\Sigma_{j\in J^-}\frac{a_{ij}}{f_0}(x_j-l_j)+\Sigma_{j\in K^+}\frac{a_{ij}}{f_0}(u_j-x_j)\ge1.\tag{10} \]
注意式子 (8) 和 (10) 的左边都是大于 0 的, 因此这两个式子联立可以得到 \[ \begin{aligned} &\Sigma_{j\in J^+}\frac{a_{ij}}{1-f_0}(x_j-l_j)-\Sigma_{j\in J^-}\frac{a_{ij}}{f_0}(x_j-l_j)\\ &+\Sigma_{j\in K^+}\frac{a_{ij}}{f_0}(u_j-x_j)-\Sigma_{j\in K^-}\frac{a_{ij}}{1-f_0}(u_j-x_j)\ge1 \end{aligned} \] 由于在当前赋值 \(\alpha\) 下, 左侧的每个元素都等于 0, 因此该赋值 \(\alpha\) 被新的约束排除. 换句话说, 可以保证对带有约束条件的线性问题的解与先前的解不同.
3 Fourier-Motzkin Variable Elimination
与 Simplex 方法类似, Fourier-Motzkin 变量消除接收实数变量线性约束的合取式, 并判定其可满足性. 效率并不如 Simplex算法, 但是针对一些较小的公式还是有优势的. 在实践中, 它主要用以消除存在量词.
让 \(m\) 代表约束的数量, \(x_1,...,x_n\) 表示这些约束所涉及到的变量. 先从消除等式开始.
3.1 等式约束的处理
首先第一步, 如下形式的等式约束要被消除: \[ \Sigma_{j=1}^na_{i,j}\cdot x_j=b_i. \] 我们在第 \(i\) 个等式约束中, 选择变量 \(x_j\) 的非零系数 \(a_{i,j}\). 不失一般性地, 假定要消除的变量是 \(x_n\). 上述约束可以改写成 \[ x_n=\frac{b_i}{a_{i,n}}-\Sigma_{j=1}^{n-1}\frac{a_{i,j}}{a_{i,n}}\cdot x_j.\tag{11} \] 在其他所有的约束中, 将等式 (11) 右边将 \(x_n\) 替换掉, 同时移除第 \(i\) 个约束. 迭代这个过程, 指导所有的等式都没移除掉. 这样系统中就只剩下不等式的形式了 \[ \bigwedge_{i=1}^m\Sigma_{j=1}^na_{i,j}x_j\le b_i. \]
3.2 变量消除
变量消除的主要思想是, 启发式地选择一个变量, 然后通过将其约束投影到系统的其余部分来消除它, 结果形成新的约束.
一个简单地例子来说明, 考虑如下约束 \[ 0\le x\le 1, 0\le y\le 1, \frac{3}{4}\le z\le 1 \] 其约束如图1所示
上述约束构造了一个长方体, 将这些约束投影到 \(x\) 轴和 \(y\) 轴, 由此就消除了变量 \(z\), 并形成了一个长方形. \[ 0\le x\le 1, 0\le y\le 1. \] 如图2所示,
图2所示的三角形区域由下面的约束构造 \[ x\le y+10, y\le 15, y\ge -x+20. \] 将该三角形投影到 \(x\) 轴就形成了一条线段 \[ 5\le x\le 25. \] 推广, 假定 \(x_n\) 是要被消除的变量. 约束根据x的系数进行划分. 考虑下标是\(i\)的约束: \[ \Sigma_{j=1}^{n}a_{i,j}\cdot x_j\le b_i. \] 拆分累加和, 可将上面的式子重写成 \[ a_{i,n}\cdot x_n \le b_i - \Sigma_{j=1}^{n-1}a_{i,j}\cdot x_j. \] 如果 \(a_{i,n}\) 为 0, 当消除 \(x_n\) 的时候可以忽略该约束. 否则, 将上述不等式两边同除以 \(a_{i,n}\). 如果 \(a_{i,n}\)为正, 得到 \[ x_n \le \frac{b_i}{a_{i,n}}-\Sigma_{j=1}^{n-1}\frac{a_{i,j}}{a_{i,n}}\cdot x_j\tag{1} \] 所以, 如果 \(a_{i,n}\gt 0\), 该约束就是 \(x_n\) 的上界. 如果 \(a_{i,n}\le 0\), 该约束就是下界.
4 The Omega Test
Omega 测试算法是判定整型变量线性约束可满足问题的算法. 可以看作是 Fourier-Motzkin 算法的变式. 这两者都不算是最快的判定过程, 但是它们都可以用作存在量词的消除.
每一个合取式可以看作是如下形式的等式或者非严格等式 \[ \Sigma_{i=1}^na_ix_i=b,\\ \Sigma_{i=1}^na_ix_i\le b. \] 系数 \(a_i\) 假定为整数; 如果不是, 则可以通过使用系数合理的假设, 将约束条件乘以分母的最小公倍数, 将问题转换为整数系数. Omega 测试算法的执行时间取决于系数 \(a_i\) 的大小. 因此要将约束进行转换, 以获得较小的系数. 这可以通过将每个约束的系数 \(a_1,a_2,...,a_n\) 除以它们的最大公约数 \(g\) 来完成. 所得的约束称为正则化(normalized)约束. 如果该约束为等式约束, 则得到 \[ \Sigma_{i=1}^n\frac{a_i}{g}x_i=\frac{b}{g}. \] 如果 \(g\) 不能整除 \(b\), 则该系统是不可满足的. 如果该约束是一个不等式, 可以通过舍入常量来增强约束: \[ \Sigma_{i=1}^n\frac{a_i}{g}x_i\le \lfloor\frac{b}{g}\rfloor. \] 比如说, 等式 \(3x+3y=2\) 可以正则化成 \(x+y=\frac{2}{3}\), 由此是不可满足的. 约束 \(8x+6y\le 0\) 可以正则化为 \(4x+3y\le 0\). 约束 \(1\le 4y\) 可以增强为 \(1\le y\).
在 Fourier-Motzkin 算法中, 等式形式和不等式形式的约束是分开考虑的; 在对不等式约束进行处理之前, 所有的等式约束是被移除了的.