單純形法

<image>

單純形法（simplex algorithm）在數學優化領域中常用於線性規劃問題的數值求解，由喬治·伯納德·丹齊格發明。

下山單純形法（Nelder-Mead method）與單純形法名稱相似，但二者關聯不大。該方法由Nelder和Mead於1965年發明，是用於優化多維無約束問題的一種數值方法，屬於更普遍的搜索算法的類別。這兩種方法都使用了單純形的概念。單純形是 $N$ 維中的 $N+1$ 個頂點的凸包，是一個多胞體：直線上的一個線段，平面上的一個三角形，三維空間中的一個四面體等等，都是單純形。

## 標準形式

假設有n個變量和m個約束。線性規劃的標準形式如下：

$$\begin{align} <pre> & & \max \sum\limits_{1\le k\le n}{{{c}_{k}}{{x}_{k}}} \\ & s.t. & \sum\limits_{1\le k\le n}{{{A}_{1,k}}{{x}_{k}}}\le {{b}_{1}}, \\ & & \sum\limits_{1\le k\le n}{{{A}_{2,k}}{{x}_{k}}\le {{b}_{2}},} \\ & ... \\ & & \sum\limits_{1\le k\le n}{{{A}_{m,k}}{{x}_{k}}}\le {{b}_{m}} \\ & & {{x}_{1}},{{x}_{2,}}...,{{x}_{n}}\ge 0 </pre> \end{align}$$

所有其他形式的線性規劃方程組都可以按照下列方式轉化成標準形式：

- 目標函數並非最大化：將所有${{c}_{k}}$取負。
- 約束條件中存在大於或等於約束：將約束兩邊取負。
- 約束條件中存在等式：將其拆分為兩個不等式（一個大於等於，一個小於等於）
- 有的變量沒有非負約束：加入新變量$x',x''(x',x''>=0)$並用$x'-x''$替換原來的變量$x$

## 鬆弛形式

可以將標準形式的線性規劃轉化為鬆弛形式，以方便運算。

在原來n個變量，m個約束的線性規劃中，加入m個新的變量，將原來的不等式化為等式：

$${{x}_{n+j}}={{b}_{j}}-\sum\limits_{1\le k\le n}{{{A}_{j,k}}{{x}_{k}}}$$

當然，此時${{x}_{n+j}}\ge 0$依然成立。

我們將${{x}_{1}},{{x}_{2}},...,{{x}_{n}}$這些變量稱為非基變量，它們構成的集合記為N。將${{x}_{n+1}},{{x}_{n+2}},...,{{x}_{n+m}}$這些變量稱為基變量，它們構成的集合記為B。簡單地理解，非基變量能夠由基變量唯一確定。

在這樣的定義下，線性規劃的鬆弛形式可以寫為如下形式：

$$\begin{align} <pre> & \max \sum\limits_{k\in N}{{{c}_{k}}{{x}_{k}}} \\ & s.t. \\ & \forall 1\le i\le n+m,{{x}_{i}}\ge 0 \\ & \forall j\in B,{{x}_{j}}={{b}_{j}}-\sum\limits_{k\in N}{{{A}_{j,k}}{{x}_{k}}} \\ </pre> \end{align}$$

因此，線性規劃的鬆弛形式可以由$c, A, b, N, B$唯一確定，$c$是長度為n的向量，$b$是長度為m的向量，$A$是m*(n+m)的矩陣。$N, B$是整數集合，分別表示非基變量集合以及基變量集合。

## 轉軸操作

轉軸操作是單純形法中的核心操作，其作用是將一個基變量與一個非基變量進行互換。可以將轉軸操作理解為從單純形上的一個頂點走向另一個頂點。

設變量${{x}_{n+d}}$屬於B（基變量），變量${{x}_{e}}$屬於N（非基變量），執行轉軸操作pivot(d,e)之後，${{x}_{n+d}}$將變為非基變量，相應地${{x}_{e}}$將變為基變量。

具體地說，一開始我們有

$${{x}_{n+d}}={{b}_{d}}-\sum\limits_{k\in N}{{{A}_{d,k}}{{x}_{k}}}$$

移項，得

$$A_{d,e}x_e = b_d-\sum\limits_{k\in N,k\ne e}A_{d,k}x_k-{x}_{n+d}$$

如果${{A}_{d,e}}\ne 0$，我們有

$${{x}_{e}}=\frac{{{b}_{d}}}{{{A}_{d,e}}}-(\sum\limits_{k\in N,k\ne e}{\frac{{{A}_{d,k}}}{{{A}_{d,e}}}{{x}_{k}})}-\frac{1}{{{A}_{d,e}}}{{x}_{n+d}}$$

將此式代入其他的約束等式以及目標函數，我們就實現了${{x}_{n+d}}$與${{x}_{e}}$在基變量和非基變量上的互換。

## 方法步驟

單純形法的一般解題步驟可歸納如下：

- 把線性規劃問題的約束方程組表達成典範型方程組，找出基本可行解作為初始基本可行解。
- 若基本可行解不存在，即約束條件有矛盾，則問題無解。
- 若基本可行解存在，從初始基本可行解作為起點，根據最優性條件和可行性條件，引入非基變量取代某一基變量，找出目標函數值更優的另一基本可行解。
- 按步驟3進行迭代，直到對應檢驗數滿足最優性條件（這時目標函數值不能再改善），即得到問題的最優解。
- 若迭代過程中發現問題的目標函數值無界，則終止迭代。

## 最優化過程

如果b向量所有元素非負，則顯然我們只需要令所有的變量等於0，就可以得到一個可行解。在這種情況下，通過下述最優化過程，我們可以得到該線性規劃的最優解，或者指出該線性規劃的最優解為無窮大（不存在）。

- 任取一個非基變量${{x}_{e}}$，使得${{c}_{e}}>0$。
- 選取一個基變量${{x}_{d}}$，使得${{A}_{d,e}}>0$，且最小化${{{b}_{d}}}/{{{A}_{d,e}}}\;$
- 執行轉軸操作pivot(d, e)，並轉到第一步繼續演算法。

根據${{{b}_{d}}}/{{{A}_{d,e}}}\;$的最小性不難證明pivot(d, e)不會破壞b的非負性。因此將所有變量取0值仍然是可行解。同時，根據$\Delta v={{c}_{e}}\frac{{{b}_{d}}}{{{A}_{d,e}}}\ge 0$，我們發現v一定是不降的。這就達到了更新解的目的。

不難發現，演算法終止有兩種情況：

- 對於所有的非基變量，c均非正。
- 對於某一個e，所有的${{A}_{d,e}}$均非正。

可以證明，對於第一種情況，我們已經得到了該線性規劃的最優解。當前的v即為答案。嚴格證明比較複雜，但是直觀上是很容易理解的。因為所有的非基變量都是非負的，而所有的c都是非正的，因此只要某個非基變量不為0，就會使得目標函數更小。

對於第二種情況來說，很容易證明此時線性規劃的最優解是無窮大。只要讓其他所有變量均為0，變量${{x}_{e}}$為正無窮。由於所有的${{A}_{d,e}}$都非正，因此非基變量的非負性得到保證。同時由於${{c}_{e}}>0$，目標函數值為正無窮。

### 實例

例：解最優化問題：

$$\min \quad Z = - x_1 - x_2$$

$$s.t. \quad 2x_1 + x_2 + x_3 = 12,$$

$$\quad \quad \quad x_1 + 2x_2+ x_4 = 9,$$

$$\quad \quad \quad x_i \geq 0, i = 1, 2, 3, 4.$$

列單純形表（即矩陣）：

$x_1$｜$x_2$｜$x_3$｜$x_4$｜b

$x_3$｜2｜1｜1｜0｜12

$x_4$｜1｜2｜0｜1｜9

c｜1｜1｜0｜0｜0

從c所在行的正數中最大的一個所對應的變量作為基變量，因為這裡兩者相等，不妨選為$x_1$。

用$x_1$所在列的正係數除b所在列的數並比較大小，有$\frac {12}2 = 6 < \frac 91 = 9$，故取$x_3$離開基變量。

然後對該矩陣進行轉軸操作，使$x_1$所在列變為單位向量：

$x_1$｜$x_2$｜$x_3$｜$x_4$｜b

$x_1$｜1｜1/2｜1/2｜0｜6

$x_4$｜0｜3/2｜-1/2｜1｜3

c｜0｜1/2｜-1/2｜0｜-6

令c所在行其餘最大的正數所在列的變量$x_2$進入基變量，並且根據$\frac 6{1/2} = 12 > \frac 3{3/2} = 2$，令$x_4$離開基變量。

繼續進行行變換，得到

$x_1$｜$x_2$｜$x_3$｜$x_4$｜b

$x_1$｜1｜0｜2/3｜-1/3｜5

$x_2$｜0｜1｜-1/3｜2/3｜2

c｜0｜0｜-1/3｜-1/3｜-7

由於c所在行的所有數均非正，問題結束。最優解為$x_1 = 5, x_2 = 2$，最優值為$Z = - x_1 - x_2 = -7$。

## 初始化過程

如果b向量並不全為非負，則我們需要通過初始化過程來找到一個可行解，然後才可以使用最優化過程進行優化。當然，此時原線性規劃不一定存在可行解。

具體的方法是，加入一個新的非基變量${{x}_{0}}$，並在原線性規劃的基礎上構造一個新的輔助的線性規劃。

$$\begin{align} <pre> & \max -{{x}_{0}} \\ & s.t. \\ & \forall 0\le i\le n+m,{{x}_{i}}\ge 0 \\ & \forall j\in B,{{x}_{j}}={{b}_{j}}-(\sum\limits_{k\in N}{{{A}_{j,k}}{{x}_{k}}})+{{x}_{0}} \\ </pre> \end{align}$$

注意這裡N集合併不包含${{x}_{0}}$。

然後，選擇一個基變量${{x}_{d}}$使得${{b}_{d}}$最小，執行轉軸操作pivot(d, 0)。不難證明該操作過後所有的b值全部非負。然後，使用前文中所述的最優化過程求解該輔助線性規劃。

由於${{x}_{0}}$非負，因此該線性規劃的答案非正。如果答案為負數，則說明原線性規劃不可能讓所有的基變量都非負，因此原線性規劃無可行解。否則，只要令所有變量為0，並去掉${{x}_{0}}$變量，就可以得到可行解。

在從輔助線性規劃轉化到原來的線性規劃的過程中，如果${{x}_{0}}$已經是非基變量，則可以將其從約束條件和目標函數中直接去掉。否則，需要任取一個非基變量${{x}_{e}}$執行pivot(0, e)，將${{x}_{0}}$變為非基變量。由於此時${{x}_{0}}$是基變量且${{x}_{0}}=0$，故${{b}_{0}}=0$一定成立，因此這個轉軸操作不會破壞b向量的非負性。

## 效率分析

在採用Bland's法則選擇用於轉軸操作的d和e（相同值的情況下取字典序最小）之後，可以證明單純形法一定能夠在有限步之後終止，但是最壞情況演算法的時間複雜度為指數函數級別的，而且可以構造出讓單純形法的時間複雜度達到指數級別的具體實例。不過實踐證明在絕大多數情況下單純形法的效率非常令人滿意。

單純形法的最壞時間複雜度為指數級別，並不意味著線性規劃不存在多項式級別的演算法。橢球演算法和內點算法均為解決線性規劃的多項式時間演算法。

## 參考

- Greenberg, Harvey J., Klee-Minty Polytope Shows Exponential Time Complexity of Simplex Method University of Colorado at Denver (1997) PDF download
- Frederick S. Hillier and Gerald J. Lieberman: Introduction to Operations Research, 8th edition. McGraw-Hill. ISBN 0-07-123828-X
- Thomas H. Cormen, Charles E. Leiserson, Ronald L. Rivest, and Clifford Stein. Introduction to Algorithms, Second Edition. MIT Press and McGraw-Hill, 2001. ISBN 0-262-03293-7. Section 29.3: The simplex algorithm, pp.790–804.
- IOI2007國家集訓隊論文，《淺談資訊學競賽中的線性規劃——簡潔高效的單純形法實現與應用》，作者：李宇騫

## 外部連結

- 線性規劃和單純形法簡介 作者Spyros Reveliotis（喬治亞理工）。
- 分步單純形法求解器 可以求解線性規劃問題。
- 基於Java的 交互式單純形工具 由Argonne國家實驗室提供。
- 單純形法 作者Elmer G. Wiens。演示演算法細節，使用了「單純形表」。
- 單純形法教程 作者Stefan Waner，hofstra.edu。
- 單純形法分解 Mazoo學習部落格.