OUTLINE
- 前言
- 預備知識預警
- 什么是column generation
- 相關概念科普
- Cutting Stock Problem
- CG求解Cutting Stock Problem
- 列生成代碼
- reference
00 前言
這幾天勤奮的小編一直在精確算法的快樂學習之中不能自拔。到列生成算法這一塊，看了好幾天總算把這塊硬骨頭給啃下來了。然后發現網上關於列生成的教學資料也不是很多，大部分講的不是那么通俗易懂。所以今天就打算寫一寫這個算法，盡可能寫得通俗易懂。
01 預備知識預警
由於列生成算法涉及的知識點非常多，所以在開始之前希望讀者必須要具備以下的基礎知識，不然就沒法往下玩了：
- 線性規划以及線性規划對偶問題
- 單純形法原理
- 原問題的影子價格（shadow price）以及對偶變量
- 單純形法非基變量進基時非基變量檢驗數（reduce cost）的計算
以上內容我就不展開科普了。如果對這些概念還有不熟悉的小伙伴，一定要回去搞清楚再往下看哦。
02 什么是column generation？
2.1 相關背景
Column generation 是一種用於求解大規模線性優化問題的非常高效的算法。[3]其理論基礎是由Danzig等於1960年提出。本質上而言，列生成算法就是單純形法的一種形式，是用來求解線性規划問題的。列生成算法已被應用於求解如下著名的NP-hard優化問題：機組人員調度問題(Crew Assignment Problem)、切割問題(Cutting Stock Problem)、車輛路徑問題(Vehicle Routing Problem)、單資源工廠選址問題(The single facility location problem )等。
2.2 larger linear programs
在某些線性優化問題的模型中,約束的數目有限，但是變量的數目隨着問題規模的增長會爆炸式的增長，因此不能把所有的變量都顯性的在模型中表達出來。比如下面一個線性規划問題：
從中可以看出，約束條件只有三個，但是當n=10000時，其變量數就達到了10000個。這類問題就是大規模的線性優化問題了。
2.3 column generation
單純型法雖然能保證在數次迭代后找到最優解，但是其面對變量很多的線性規划問題就顯得很弱了。因為它需要去在眾多變量里進行基變換，就上面的問題而言，你想想你要在近10000個變量里面找個能進基的，活着不好嗎？
再有，在用單純形法求解這類線性規划問題時，基變量(basic variable)只與約束的個數相關，每次迭代只會有一個新的非基變量(non-basic variable)進基，因此，在整個求解過程中其實只有很少一部分變量會被涉及到。
因此，有人基於單純型法提出了列生成算法。其思路大概如下：[1]
-
先把原問題(master problem)restrict到一個規模更小（即變量數比原問題少的）的restricted master problem，在restricted master problem上用單純型法求最優解，但是此時求得的最優解只是restricted master problem上的，並不是master problem的最優解。
-
此時，就需要通過一個subproblem去check在那些未被考慮的變量中是否有使得reduced cost小於零的？如果有，那么就把這個變量的相關系數列加入到restricted master problem的系數矩陣中，回到第1步。
經過反復的迭代，直到subproblem中的reduced cost rate大於等於零，那么master problem就求到了最優解。
看算法流程圖會更加直觀哦：[2]
03 相關概念科普
剛剛講的內容涉及到了幾個概念，master problem，restricted master problem，subproblem等，這一節來把這幾個概念給講清楚。基於一個問題線性規划問題：
3.1 master problem(MP)
master problem可以認為是原問題。即
3.2 restricted master problem(RMP)
前面我們說過，把原問題(master problem)restrict到一個規模更小（即變量數比原問題少的）的就是restricted master problem了。比如可以用啟發式算法，在上面的master problem找出滿足條件（也就是形成的restricted master problem必須要有可行解）的k個列，得到如下的restricted master problem：
可以看到，相比原來的master problem，restricted master problem相當於把\(y_{k+1}...y_m\)強制限制為非基變量了。[4]
3.3 subproblem
核能預警，如果這部分看不懂，請確保預備知識過關。如果預備知識不過關，請在運籌學老師的陪同下觀看，謝謝合作！
上面的限制主問題求解完成后，我們想使用單純型法進行基變量的轉換，看看\(y_{k+1}...y_m\)中，是否有可以轉入基變量的列。還記得怎么找進基的非基變量嗎？（不記得就去問你們的運籌學老師）。當然是通過非基變量的檢驗數辣，通過\(\sigma_j = c_j - c_BB^{-1}a_j\)，在\(y_{k+1}...y_m\)中尋找檢驗數最小並且為負數的變量，將變量對應的那一列添加到RMP中。
那么，在檢驗數的計算公式中，大家還記得\(c_BB^{-1}\)是什么嗎？\(c_BB^{-1}\)有兩重含義：
- 通過求解RMP問題得到的影子價格（shadow price）。
- 通過求解RMP對偶問題得到的對偶變量(dual variable)。
所以在開始之前小編一直強調預備知識一定要過關。這兩個含義意味着我們有上面兩種方式得到\(c_BB^{-1}\)，不過我們一般傾向於使用第二種，WHY？
雖然通過單純型法直接求解restricted master problem能得到\(c_BB^{-1}\)。但是restricted master problem也可能是一個變量很多的線性規划。前面也說過了，單純型法對變量很多的問題是無能為力的。因此通過單純型法求restricted master problemde的對偶問題（將restricted master problem對偶一下，就能使得變量數大幅減小，因為這些變量轉換成了對偶問題中的限制條件了），能更快地得到子問題想要的\(c_BB^{-1}\)。[1]
所以我們總結一下：
通過求解RMP問題或者RMP對偶問題，得到我們想要的\(c_BB^{-1}\)以后，subproblem就是通過\(\sigma_j = c_j - c_BB^{-1}a_j\)這條公式，在\(y_{k+1}...y_m\)中尋找檢驗數為負並且最小的變量，將變量對應的那一列添加到RMP中。
3.4 算法流程圖
通過上面講了這么多以后，這里在給出一個更詳細的流程圖：[5]
04 Cutting Stock Problem[1]
講column generation怎么可能少得了Cutting Stock Problem這個經典的問題呢！
我們有以下問題，原紙卷每個長為L=17m，顧客們分別需要25個3m長，20個5m長，18個7m長的紙卷。那么需要怎樣切割才能使得浪費最小呢？
Master Problem
Column Generation Formulation：
- \(P\)是所有可行的裁剪方案的集合，里面方案的總數為n(我們並不需要確切的知道這個值是多少，只需要知道它很大)。
- \(a_{ij}\) 表示第j種方案里類別i的個數。
- \(y_{j}\)表示第 j 種方案的選擇個數。
於是，我們得到如下模型：
這樣，我們得到了Cutting Stock Problem的Master Problem。
05 CG求解Cutting Stock Problem
通過上面的問題分析和建模以后，我們這一步一步一步來求解該問題，讓大家徹底理解column generation這個過程。該過程模擬需要用到一個線性求解器，大家還記得小編以前講過的lpsolve的教程嗎?趕緊去翻一下以前的教程，把lpsolveIDE裝上，然后跟着小編的腳步一步一步往下走。
5.1 restricted master problem(RMP)
前面我們完成了問題的建模，得到了Cutting Stock Problem的Master Problem。現在，我們可以用啟發式算法找到一個滿足客戶需要的初始解：
首先，一個卷筒有三種切割方案：
方案1：切成5個3m
方案2：切成2個6m
方案3：切成2個7m
很容易得出，5個方案1、10個方案2、8個方案3，是能滿足所有客戶需求的。即得到MP的一個RMP如下：
其中，
這三列分別對應着5個方案1、10個方案2、8個方案3。還有一點需要注意的，對於每一列，都需要滿足：
\(3a_{1j} + 6a_{2j}+ 7a_{3j} \le 16\)，也就是每一卷紙只有16的長度，不能超出這個長度。這個叫列生成規則，不同問題有不同的規則約束。subproblem在尋找某些列或者生成某些列時，就是受到列生成規則的約束的。
5.2 開始列生成過程
iteration 1
RMP：
將該模型輸入lpsolve，得到對偶變量如下：
得到\(c_BB^{-1} = [0.2, 0.5, 0.5]\)。現在要找一列加入RMP，是哪一列呢？現在還不知道，我們暫記為\(\alpha_4 = [a_{14},a_{24},a_{34}]^T\)。非基變量檢驗數\(\sigma_4 = c_4 - c_BB^{-1}\alpha_4 = 1 - 0.2a_{14}-0.5a_{24}-0.5a_{34}\)。
subproblem：
求解結果得$ \alpha_4 = [1,2,0]^T, \sigma_4= -0.2 < 0\(,reduced cost 為負數，因此將\) \alpha_4$加入RMP，開始第二輪迭代。
iteration 2
RMP：
將該模型輸入lpsolve，得到對偶變量如下：
得到\(c_BB^{-1} = [0.2, 0.4, 0.5]\)。現在要找一列加入RMP，是哪一列呢？現在還不知道，我們暫記為\(\alpha_5 = [a_{15},a_{25},a_{35}]^T\)。非基變量檢驗數\(\sigma_5 = c_5 - c_BB^{-1}\alpha_5 = 1 - 0.2a_{15}-0.4a_{25}-0.5a_{35}\)。
subproblem：
求解結果得$ \alpha_5 = [1,1,1]^T, \sigma_5= -0.1 < 0\(,reduced cost 為負數，因此將\) \alpha_5$加入RMP，開始第三輪迭代。
iteration 3
RMP：
將該模型輸入lpsolve，得到對偶變量如下：
得到\(c_BB^{-1} = [0.2, 0.4, 0.4]\)。現在要找一列加入RMP，是哪一列呢？現在還不知道，我們暫記為\(\alpha_6 = [a_{16},a_{26},a_{36}]^T\)。非基變量檢驗數\(\sigma_6 = c_6 - c_BB^{-1}\alpha_6 = 1 - 0.2a_{16}-0.4a_{26}-0.5a_{36}\)。
subproblem：
求解結果得$ \alpha_6 = [5,0,0]^T, \sigma_6 = 0\(,reduced cost 不為負數，因此不用將\) \alpha_6$加入RMP，列生成算法結束。
最終，我們求解最后一次迭代的RMP：
得到RM的最優解\(y = [1.2, 0,0,1, 18]\),聰明的同學已經主要到了，\(y_1=1.2\)怎么出現了小數呢，按理說y應該是整數才對啊。回到原問題RM：
我們並沒有加上\(y_i \in Z\)這個約束，這是因為我們在用列生成的時候把這個模型給松弛為了線性模型，畢竟列生成是用於求解linear program的。如果要求解大規模整數規划問題，列生成是無法辦到的，后面我們會介紹結合column generation的branch and price方法。
至此，我們已經完完整整把列生成算法給走了一遍。相信列生成算法的原理已經深入各位讀者的心里啦。
06 列生成代碼
關於Cutting Stock Problem的列生成java代碼，可以關注我們的公眾號：
請關注公眾號【程序猿聲】，后台回復【cgcsp】，不包括【】即可下載！
