在这篇文章中，我将介绍一系列常用的预处理技术，以进一步提升对偶单纯形法的求解效率。
问题形式转换 Transformation
在本系列之前的介绍中，我们假定原始问题具有以下形式
$$
\begin{aligned}
\min \ & c^T x \\
s.t. \ & Ax = b, \\
& l \leq x \leq u.
\end{aligned}
\label{eq:primal}\tag{1.1}
$$ 相应的对偶问题为
$$
\begin{aligned}
\max \ & b^T \lambda + u^T s_u + l^T s_l \\
s.t. \ & A^T \lambda + s_u + s_l = c, \\
& s_l \geq 0, \ s_u \leq 0.
\end{aligned}
\label{eq:dual}\tag{1.2}
$$ 但在实际的求解场景中，原始问题一般是按以下形式存储的：
$$
\begin{aligned}
\min \ & c^T x \\
s.t. \ & b_l \leq Ax \leq b_u, \\
& l \leq x \leq u. \\
\end{aligned}
\label{eq:primal-in}\tag{1.3}
$$ 该形式的好处在于线性约束的表示能力更强，在极端情况下可以节省一半的表示空间。但是，我们无法对这种问题形式应用DS算法；为此，需要将问题\eqref{eq:primal-in}转换成问题\eqref{eq:primal}的形式。一种简单的方法是引入新变量$x^b$来代替$Ax$，从而将线性约束转为上下界约束：
$$
\begin{aligned}
\min \ & c^T x + 0^T x^b \\
s.t. \ & Ax - x^b = 0, \\
& l \leq x \leq u, \ b_l \leq x^b \leq b_u.
\end{aligned}
\label{eq:primal-in-reform}\tag{1.4}
$$ $x_b$一般被称为逻辑变量/行变量。上述转换的另一个好处在于转换后的矩阵$\bar{A} = [A;-I]$是满秩的，因此降低了迭代过程中出现高条件数的基矩阵$A_B$的几率。然而，当线性约束为等式$Ax=b$时，上述转换方式会引入多余的固定值变量$x^b=b$（虽然这并不会影响求解效率）。
在通过DS流程完成对问题\eqref{eq:primal-in-reform}的求解后，需要将结果转换回原问题\eqref{eq:primal-in}的解。通过比较两者的形式易知，原始解$x$的转换是trivial的；而对于对偶解$(\lambda,s)$，注意到问题\eqref{eq:primal-in-reform}对应的对偶问题为
$$
\begin{aligned}
\max \ & u^T s_u + l^T s_l + b_u^T s_u^b + b_l^T s_l^b \\
s.t. \ & A^T \lambda + s_u + s_l = c, \\
& -\lambda + s_u^b + s_l^b = 0, \\
& s_l,s_l^b \geq 0, \ s_u,s_u^b \leq 0.
\end{aligned}
\label{eq:dual-in}\tag{1.5}
$$ 简单分析可知新增的对偶松弛变量$s_u^b$和$s_l^b$正好分别对应原问题\eqref{eq:primal-in}中$\leq$和$\geq$方向的约束的对偶变量，因此对偶解的转换也是trivial的。
虽然在DS求解流程中上述问题转换步骤是必需的，但是直接对原问题\eqref{eq:primal-in}进行转换会影响求解效率：原问题形式中一些易于化简的结构可能会在形式转换后变得难以化简。因此，一般在完成其他预处理步骤后再进行问题转换。在下面的讨论中，我都将假定原始问题的形式为问题\eqref{eq:primal-in}。
系数缩放 Scaling
在本系列之前的介绍中，我曾提到：DS迭代过程的数值误差很大程度受到矩阵$A_B$的条件数的影响。由于一个矩阵$A$的条件数与$A^TA$的最大最小特征值的比例有关，而特征值的大小与矩阵系数之间的量级差异有关（一个简单的例子是对角阵$diag([10^5,10^{-5}])$）；因此，可以通过对$A$的系数进行缩放来平衡系数的数量级，而这在一定程度上可降低$A_B$的条件数以及问题的数值计算难度。
一种简单的系数缩放思路是按行和列进行缩放。具体来说，假设给定行权重$w_b\in \mathbb{R}^n$和列权重$w_c \in \mathbb{R}^m$，问题\eqref{eq:primal-in}可等价转化为
$$
\begin{aligned}
\min \ & (w_c \odot c)^T x^{(w)} \\
s.t. \ & w_b \odot b_l \leq (diag(w_b) \cdot A \cdot diag(w_c) )x^{(w)} \leq w_b \odot b_u, \\
& (1/w_c) \odot l \leq x^{(w)} \leq (1/w_c) \odot u. \\
\end{aligned}
\label{eq:primal-in-weight}\tag{2.1}
$$ 其中$\odot$表示向量的元素相乘(Hadamard Product)运算，$diag(\cdot)$表示将向量转化成对角矩阵。形式\eqref{eq:primal-in-weight}中同时考虑了行维度和列维度的缩放；通过令$w_b = \mathbf{1}$或$w_c = \mathbf{1}$，可以实现只对行维度或列维度进行缩放。
由于我们的目标是平衡矩阵$A$的系数的量级，因此，行/列权重一般被选为行/列系数的量级的倒数。以下是几种常用的权重选取方式：
- $l_2$范数：对于行$i$，令$w_i=1/\sqrt{\sum_j A_{ij}^2}$；对于列$j$，令$w_j=1/\sqrt{\sum_i A_{ij}^2}$；
- $l_{\infty}$范数（最大值）：对于行$i$，令$w_i=1/\max_j |A_{ij}|$；对于列$j$，令$w_j=1/\max_i |A_{ij}|$；
- 几何平均(geomean)：对于行$i$，令$w_i=1/\sqrt{\max_{\{j| A_{ij} \neq 0\}} |A_{ij}|\cdot \min_{\{j| A_{ij} \neq 0\}} |A_{ij}|}$；
对于列$j$，令$w_j=1/\sqrt{\max_{\{i| A_{ij} \neq 0\}} |A_{ij}|\cdot \min_{\{i| A_{ij} \neq 0\}} |A_{ij}|}$。
在实际操作中，一般会多次调用上述缩放方法并配合不同权重选取方式，以获得更优的系数整体平衡效果。例如，在本系列的代码实现中，我将首先使用基于$l_2$的行+列权重，然后使用基于$l_{\infty}$的行权重，最后再使用基于$l_{\infty}$的列权重。
问题简化 Reductions
根据此前的理论介绍和数值实验可知，线性代数相关计算是DS迭代过程中的主要计算开销来源（之一）；因此，如果能在进入DS迭代过程之前对线性系统$A$做一些简化，则可以一定程度上降低迭代中的计算量。简化思路主要有以下三点：
- 减少行（线性约束）的个数。首先，可以根据变量的上下界约束（和其他线性约束）来判断一个线性约束是否是多余的。其次，在DS迭代中，上下界约束相比线性约束更容易处理（计算量更少），因此可以尝试将一些简单的线性约束（等价地）转换为上下界约束。
- 减少列（决策变量）的个数。首先，可以根据约束和最优性条件，判断变量的取值是否可以固定。其次，可以基于等式约束将变量表示为其他变量的线性组合；不过这种做法可能会导致$A$的非零元素数量增加。
- 对变量的上下界进行调整。通过放松(relax)上下界约束，可以降低可行域的复杂度（多面体的顶点个数），有利于原始单纯形法(PS)迭代。而通过加强(tighten)上下界约束，特别是将无界转换成有界，可以提高问题的对偶可行程度，有利于DS迭代。另外，在调整上下界约束后，可能可以进一步对行和列进行简化。
需要注意的是，简化后问题和变量的定义发生了变化，因此在求解结束后还需要根据简化问题的结果来恢复原始问题的解；该步骤一般被称为后处理(postprocessing/post-solve)。
接下来，我将介绍几种常用的简化手段（及其后处理方法）。
空行/列 Empty Row/Column
空行和空列中没有非零元素，即与其他行和列之间没有相互作用，故在简化步骤中只需要检查相应的（原始/对偶）可行性条件。具体来说，对于一个空行$i$，相应的约束为$b_{l,i} \leq 0 \leq b_{u,i}$，因此根据$b_{l,i}$和$b_{u,i}$的具体取值，该约束要么是多余的，要么会导致问题原始不可行。对于一个空列$j$，相应的变量$j$只受到上下界约束$l_j \leq x_j \leq u_j$的影响，因此可以根据$c_j$的符号直接取上界或下界；如果取到的值为$\pm \infty$（即无界），则问题对偶不可行。
在后处理阶段，对于空列$j$，对应的原始变量$x_j$的取值已知，而对偶松弛变量的取值为
$$s_j = c_j - A_j^T \lambda = c_j - 0^T \lambda = c_j; \tag{3.1}
$$对于空行$i$，相应的对偶变量$\lambda_i$和$s_{i}^b$都可以直接置0。
固定值变量 Fixed-value Column
对于一个变量$j$，其应满足的上下界约束为$l_j \leq x_j \leq u_j$。如果$l_j > u_j$，那么问题原始不可行。如果$l_j = u_j$，那么可以固定$x_j$的取值为$l_j$并去除该变量。进行该简化操作时，$b$的取值需要进行相应改动：
$$b_l \gets b_l - A_j x_j,b_u \gets b_u - A_j x_j. \tag{3.2}
$$ 在后处理阶段，原始变量$x_j$的取值已知，而对偶松弛变量的取值为$s_j = c_j - A_j^T \lambda$。
单变量约束 Singleton Row
如果行$i$只包含一个非零元素$A_{ij}$，那么该线性约束相当于对变量$j$的上下界约束。具体来说，假设该线性约束为$b_{l,i} \leq A_{ij}x_j \leq b_{u,i}$，其中$A_{ij} > 0$，那么可以直接将其与已有的上下界$l_j \leq x_j \leq u_j$进行整合：
$$
\max(l_j,b_{l,i}/A_{ij}) \leq x_j \leq \min(u_j,b_{u,i}/A_{ij}), \tag{3.3}
$$ 并进行相应的可行性判断。
在后处理阶段，需要给出约束$i$对应的对偶变量$\lambda_i,s_{l,i}^b,s_{u,i}^b$的取值。由于$A_j^T\lambda + s_j = c_j$，而$\lambda_{-i}$的取值已知，因此只需要确定$s_j$的取值。如果$x_j$的取值为$l_j$或$u_j$，则恢复约束$i$后$x_j$和$s_j$的取值不变，有$\lambda_i = s_{i}^b = 0$；否则，恢复约束$i$后列$j \in B$，进而有$s_j = 0, \ \lambda_i = (c_j - A_{(-i)j}^T\lambda_{-i})/A_{ij}$。
可行性隐含约束 Feasibly-trivial Row
根据变量的上下界约束$u_j,l_j$，可以计算每个线性约束的隐含上下界：
$$
\begin{aligned}
b_i^+ & = \sum_{\{j|A_{ij} > 0\}} A_{ij}u_j + \sum_{\{j|A_{ij} < 0\}} A_{ij}l_j; \\
b_i^- & = \sum_{\{j|A_{ij} > 0\}} A_{ij}l_j + \sum_{\{j|A_{ij} < 0\}} A_{ij}u_j. \\
\end{aligned}\tag{3.4}
$$ 然后可以根据$b_i^+,b_i^-,b_{u,i},b_{l,i}$的相对大小进行简化。具体来说，如果$b_i^- > b_{u,i}$或者$b_i^+ < b_{l,i}$，那么该线性约束原始不可行。如果$b_i^- = b_{u,i}$或者$b_i^+ = b_{l,i}$，那么可以根据式(3.4)固定相关原始变量的取值。如果$b_{l,i} \leq b_i^- \leq b_i^+ \leq b_{u,i}$，则约束是多余的，可以被去除。
后处理阶段的处理方法需要根据不同的简化情况进行讨论。
- 如果约束$i$是作为多余约束被去除的，那么恢复该约束不影响当前解的可行性，故只需令$i^b\in B$以及对偶变量$\lambda_i=s_i^b=0$。
- 如果约束$i$还固定了变量取值，那么该约束是有效的，需要对对偶变量$\lambda_i$和$s_j$进行计算。不妨假设$b_i^- = b_{u,i}$，则在原问题中约束$i$被固定在上界、值$A_{ij}x_j$被固定在下界。因此，对偶变量需要满足以下条件
$$
\begin{aligned}
\lambda_i & = s_i^b \leq 0; \\
s_j &
\begin{cases}
\geq 0, & A_{ij} > 0; \\
\leq 0, & A_{ij} < 0. \\
\end{cases}
\end{aligned}\tag{3.5}
$$ 又因为新增了一个约束，因此在各恢复的列$j$和$i^b$中间需要有一个进入$B$，其对应的对偶松弛变量必须为0。现在，根据对偶松弛变量的计算公式$s_j = (c_j - A_{(-i)j}^T \lambda_{-i}) - A_{ij} \lambda_i$，可以确定唯一的$\lambda_i \leq 0$使得(3.5)中所有约束都能满足，且各$s_j$和$s_i^b$中恰好有一个为0：$\lambda_i = \min(\min_{j}\{(c_j - A_{(-i)j}^T \lambda_{-i})) / A_{ij} \},0)$，并进而确定所有对偶变量的取值。
最优性隐含变量 Optimally-trivial Column
上面关于空列的讨论指出，在没有线性约束$b_l \leq Ax \leq b_u$的情况下，问题\eqref{eq:primal-in}的最优解是平凡(trivial)的：
$$
x_j = \begin{cases}
u_j, & c_j < 0; \\
l_j, & c_j > 0. \\
\end{cases}\tag{3.5}
$$ 如果线性约束$b_l \leq Ax \leq b_u$对一个变量$x_j$取到上述最优解没有明显约束作用，那么可以将该变量当做空列/固定值处理。具体来说，如果$c_j < 0$且系数$A_{ij}$满足
$$ A_{ij} \begin{cases}
\geq 0, & b_{u,i} = \infty; \\
\leq 0, & b_{l,i} = -\infty, \\
\end{cases}\tag{3.6}
$$ 那么可以在保持其他变量$x_{-j}$不变的情况下令$x_j=u_j$来最优化目标；因此最优解必然满足$x_j \equiv u_j$。同理，我们可以给出$x_j \equiv l_j$的简单检查条件。如果在上述简化步骤中列$j$被去除了，那么在后处理阶段可以按以下方式计算和恢复相应的对偶松弛变量：$s_j = c_j - A_j^T \lambda$。
上述简化思路也可以从对偶可行性的角度来理解。根据对偶可行性，如果对偶松弛变量$s_j$非零，则其符号确定了列$j$的属性$\{L,U\}$及$x_j$的取值。又因为$s_j = s_{u,j} + s_{l,j} = c_j - A_j^T \lambda$，因此可以根据对偶变量$\lambda$的上下界$[\lambda^-,\lambda^+]$来检查$s_j$的符号是否固定。特别地，在没有额外信息的情况下，可以根据等式$\lambda_i = s_{u,i}^b + s_{l,i}^b = s_i^b$以及$s^b$的平凡上下界$\{0,\pm\infty\}$来获取$\lambda$的上下界；这种方式恰好与上述简化思路相对应。也可以进一步通过上下界加强(bound tightening)手段处理$\lambda$的上下界，以期固定更多变量。
双变量等式约束 Doubleton Row
如果约束$i$是一个只包含两个变量的等式约束$A_{ij_1} x_{j_1} + A_{ij_2} x_{j_2} = b_i$，那么可以用其中一个变量$x_{j_2}$来表示另一个变量$x_{j_1}$：$x_{j_1} = (b_i - A_{ij_2} x_{j_2}) / A_{ij_1}$，然后去除变量$x_{j_1}$和约束$i$。注意保留变量$x_{j_2}$的上下界也需要根据$x_{j_1}$的上下界进行调整，可能会变得更紧。
在后处理阶段，原始变量$x_{j_1}$可以直接通过上述表示来获取，而对偶变量$\lambda_i$和$s_{j_1},s_{j_2}$需要根据列$j_2$的属性和$x_{j_2}$的取值进行相应处理。具体来说，
- 如果$x_{j_2}$在原始边界$\{l_{j_2},u_{j_2}\}$上或$j_2$在$B$中，则恢复变量$x_{j_1}$和约束$i$后，$j_2$的属性和$x_{j_2}$的取值不需要进行改变；故可令$j_1 \in B$，进而$s_{j_1} = 0, \ \lambda_i = (c_{j_1} - A_{(-i)j_1}^T \lambda_{-i})/A_{ij_1}$。由于恢复前（简化后）的松弛变量$\tilde{s}_{j_2}$满足$\tilde{s}_{j_2} = s_{j_2} - (A_{ij_2}/A_{ij_1}) s_{j_1}$，因此$s_{j_2}$的取值没有变化。
- 如果$x_{j_2}$在新生成的边界上，那么恢复变量$x_{j_1}$和约束$i$后，$j_2$应进入$B$，故$s_{j_2} = 0$，进而可以算得
$$\lambda_i = (c_{j_2} - A_{(-i)j_2}^T \lambda_{-i})/A_{ij_2}, \ s_{j_2} = - (A_{ij_1}/A_{ij_2})\tilde{s}_{j_2}.$$
接下来，我将在这个notebook中尝试实现上述各预处理方法，并通过数值实验分析效果。