散文網(wǎng) » 科技 »學(xué)習(xí) » Benders分解及其python實(shí)現(xiàn)

Benders分解及其python實(shí)現(xiàn)

2023-03-21 00:59 作者:ZJU前校長(zhǎng)-吳朝暉院士 0人讀過 | 我要投稿

本文的思路和代碼部分參考了這篇知乎專欄：【整數(shù)規(guī)劃(九)】Benders 分解（理論分析+Python代碼實(shí)現(xiàn)） - 知乎 (zhihu.com): https://zhuanlan.zhihu.com/p/428706477。本文對(duì)思路和代碼進(jìn)行了改進(jìn)，添加了自己對(duì)理論的理解，同時(shí)修改代碼，刪改了一些我認(rèn)為不必要的變量，并使得y由實(shí)數(shù)域擴(kuò)展到n維實(shí)數(shù)空間。

benders分解的出發(fā)點(diǎn)及原理說明

考慮如下的規(guī)劃問題： $max_{x,y}c^Tx+d^Ty \tag{1}$ $s.t. Ax+By \leq b \tag{2}$ $x \in R^m_+,y\in Y \subset Z^n_+ \tag{3}$

在求解此問題時(shí)，y作為一個(gè)整數(shù)變量是不好處理的（實(shí)際上，Y可以是任何“難以處理的變量”）。如果我們能夠”消掉“y而保留x，那么這個(gè)問題就變成了一個(gè)線性規(guī)劃問題，我們可以用單純形法去求解。自然而然地，我們有了這么一個(gè)想法：我們先固定一組y，此時(shí)問題變?yōu)橐粋€(gè)關(guān)于x的線性規(guī)劃問題；求解這個(gè)線性規(guī)劃問題，我們可以獲得對(duì)應(yīng)于這組y的最優(yōu)解；然后我們遍歷所有可能的y，找到最優(yōu)解最大的那組y，就能獲得全局最優(yōu)解。

基于上述討論，我們可以將問題（1）-（3）分解為兩個(gè)優(yōu)化問題： $W=max_y d^Ty+\phi(y)\\ s.t. y\in Y \subset Z^n_+ \tag{4}$ $\phi(y) = max_x c^Tx\\ s.t.Ax \leq b-By\\ x\in R^m_+ \tag{5}$

有了上述的優(yōu)化模型，我們就可以實(shí)現(xiàn)我們之前的想法：固定y求解x，再遍歷所有的y，最終找到全局最優(yōu)點(diǎn)。然而這個(gè)想法有一個(gè)嚴(yán)重的問題：當(dāng)y的可行域非常大的時(shí)候，我們想要遍歷y是不現(xiàn)實(shí)的。

那么有沒有一種方法可以讓我們只遍歷少數(shù)y就能找到最優(yōu)解呢？實(shí)際上，我們可以發(fā)現(xiàn)兩個(gè)事實(shí)，它們可以減少我們需要遍歷的y的數(shù)量：

很多y會(huì)使得問題（5）無解，對(duì)于這樣使問題無解的y，我們是不用遍歷的。
在遍歷過程中，我們可以記錄當(dāng)前遍歷過程中已經(jīng)找的使W最大的y；對(duì)于那些必然不可能獲得更大W的y，我們是不用遍歷的。

那么，我們?cè)鯓硬拍芑谏鲜鱿敕ǎ懦暨@些不需要遍歷的y呢？我們只需要取問題（5）的對(duì)偶問題即可。 $min_u(b-By)^Tu\\ s.t.A^Tu \geq c \\ u\in R^p_+ \tag{6}$ 問題（6）是問題（5）的對(duì)偶問題。由于問題（5）是一個(gè)線性規(guī)劃問題，根據(jù)對(duì)偶理論，對(duì)偶問題（6）的最優(yōu)解等于原問題（5）的最優(yōu)解；同時(shí)，問題（6）的解為負(fù)無窮，則問題（5）不可行。

我們希望找到使得問題（5）無解的y，也就是找到那些使得問題（6）的解為負(fù)無窮的y。為什么要做這個(gè)轉(zhuǎn)換呢？因?yàn)閱栴}（6）具備一個(gè)優(yōu)良的性質(zhì)：它的可行域是與y無關(guān)的，y僅出現(xiàn)在目標(biāo)函數(shù)值中。對(duì)線性規(guī)劃有了解的同學(xué)或許知道，任何一個(gè)線性規(guī)劃的可行解都可以表示為可行域的極點(diǎn)和極射線的線性組合。為了防止大家不熟悉什么是極點(diǎn)和極射線，下面舉一個(gè)簡(jiǎn)單的例子：

在上面這個(gè)圖中，假設(shè)紫色區(qū)域是可行域，則紅色線就是可行域的極射線（無限延伸的射線），綠色的點(diǎn)就是可行域的極點(diǎn)（頂點(diǎn)）；可行域內(nèi)任何一個(gè)點(diǎn)沿著極射線方向無限延伸后仍然處于可行域內(nèi)。

顯然，如果我們希望問題（6）的解不是負(fù)無窮，就要求對(duì)于可行域內(nèi)的任意一點(diǎn)沿著任意極射線方向延伸，其目標(biāo)函數(shù)值都不會(huì)降低；否則只要無限延伸，就能趨近于負(fù)無窮。從數(shù)學(xué)的語言上來說，假設(shè)極射線的的方向向量為v，問題（6）的解不是負(fù)無窮的充要條件為： $v_n^T(b-By) \geq 0, \forall n \in {1,2,3,...,N} \tag{7}$ 其中N表示極射線的數(shù)目。根據(jù)式（7），假如我們知道問題（6）可行域中所有極射線的方向，那么對(duì)于任意的y，我們不用去求解問題（6）就能知道它會(huì)不會(huì)趨近于負(fù)無窮，也就知道了問題（5）可不可行。如果我們把式（7）對(duì)應(yīng)的約束加到問題（4）中，那么就可以排除掉相當(dāng)大一部分的y，極大地減少計(jì)算量。將問題（5）轉(zhuǎn)換為問題（6）給我們帶來的好處還不止上面這些。學(xué)過線性規(guī)劃的同學(xué)們肯定知道，根據(jù)線性規(guī)劃的基本定理，若問題（6）有解且非無窮大，則解必然出現(xiàn)在可行域的極點(diǎn)上（也就是多面體的頂點(diǎn)上，上圖的綠色點(diǎn)上）。也就是說，如果我們給定一個(gè)y并且求解問題（6），得到的最優(yōu)解u*必然是可行域的一個(gè)頂點(diǎn)。我們假設(shè)多面體極點(diǎn)的集合為：

? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?? ? $(u_1%2Cu_2%2C...%2Cu_M)$

由于問題（6）有解且非無窮大時(shí)解必然出現(xiàn)在可行域的極點(diǎn)上，因此優(yōu)化問題（6）等價(jià)于如下問題（其思想是遍歷所有極點(diǎn)，找到目標(biāo)函數(shù)值最小的那個(gè)極點(diǎn)）：

$%5Cphi(y)%3Dmin_%7Bm%20%5Cin%20%5C%7B1%2C2%2C...%2CM%5C%7D%7D(b-By)%5ETu_m%20%5Ctag%7B8%7D$

極小化函數(shù)等價(jià)于最大化函數(shù)的下界，我們?cè)O(shè)函數(shù)的下界為η，于是問題（8）可以進(jìn)一步等價(jià)為：

$%5Cphi(y)%3Dmax%20%5Cquad%20%5Ceta%20%5C%5C%0A(b-By)%5ETu_m%20%5Cgeq%20%5Ceta%20%2C%20%5Cforall%20m%20%5Cin%20%5C%7B1%2C2%2C...%2CM%5C%7D%20%5C%5C%0A%5Ceta%20%5Cin%20R%5Ctag%7B9%7D$

我們把式（9）中φ(y)的表達(dá)式帶入到式（4）中：

$W%3Dmax_y%20%5CBigg(%20d%5ETy%2Bmax_%5Ceta%20%5Cbigg(%5Ceta%3A%20(b-By)%5ETu_m%20%5Cgeq%20%5Ceta%2C%20%5Cforall%20m%20%5Cin%20%5C%7B1%2C2%2C...%2CM%5C%7D%0A%5Cbigg)%20%2Cy%5Cin%20Y%20%5Csubset%20Z%5En_%2B%20%5CBigg)%20%5Ctag%7B10%7D$

問題（10）是連續(xù)求兩個(gè)最大值，所以中間的括號(hào)可以拆開。我們拆掉括號(hào)，整理得：

$W%3Dmax_%7By%2C%5Ceta%7D%20d%5ETy%2B%5Ceta%20%5C%5C%0As.t.%20(b-By)%5ETu_m%20%5Cgeq%20%5Ceta%20%2C%5Cforall%20m%20%5Cin%20%5C%7B1%2C2%2C...%2CM%5C%7D%5C%5C%0Ay%5Cin%20Y%20%5Csubset%20Z%5En_%2B%20%5Ctag%7B11%7D$

基于之前的討論，我們把式（7）加入到（11）中，保證問題存在可行解：

$W%3Dmax_%7By%2C%5Ceta%7D%20d%5ETy%2B%5Ceta%20%5C%5C%0As.t.%20(b-By)%5ETu_m%20%5Cgeq%20%5Ceta%20%2C%5Cforall%20m%20%5Cin%20%5C%7B1%2C2%2C...%2CM%5C%7D%5C%5C%0Av_n%5ET(b-By)%20%5Cgeq%200%2C%20%5Cforall%20n%20%5Cin%20%5C%7B1%2C2%2C3%2C...%2CN%5C%7D%5C%5C%0Ay%5Cin%20Y%20%5Csubset%20Z%5En_%2B%20%5Ctag%7B12%7D$

理論上來說，我們只要解上面這個(gè)整數(shù)規(guī)劃問題，就能獲得問題的最優(yōu)解。

值得注意的是，問題（12）中間“兩項(xiàng)”限制條件遠(yuǎn)遠(yuǎn)不止兩條——實(shí)際上，當(dāng)u的維數(shù)p很高的時(shí)候，它的可行域的極點(diǎn)與極射線可能會(huì)非常非常多（它們的數(shù)量會(huì)呈現(xiàn)指數(shù)級(jí)上升，這是很容易想象的），所以我們?cè)谇蠼鈫栴}的時(shí)候不可能一口氣把所有限制條件都加進(jìn)去，而是要在求解的過程中不斷添加——實(shí)際上，每給定一組y，求解問題（6)，我們得到的u或者是可行域的極點(diǎn)，或者是可行域的極射線（這里解不作證明了）；我們只要把對(duì)應(yīng)的u或者v帶入到（12）中，就能添加一個(gè)限制條件，減少y的可行域。

有同學(xué)可能會(huì)好奇，既然極點(diǎn)和極射線非常多，而我們求解一次只能添加一條極射線或者極點(diǎn)，會(huì)不會(huì)效率很低？經(jīng)驗(yàn)表明，很多時(shí)候我們只要添加少量的限制條件，就能極大地減少y的可行域，并輕易地找到最優(yōu)解，所以這個(gè)方法是足夠有效的。

Benders算法的流程

Step 1:初始化 $\text{[math]}$ ，初始化上界UB=+∞，下界LB=-∞。

Step 2: 若UB-LB=δ（提前給定的算法精度），則進(jìn)入下一步，否則直接輸出最優(yōu)解

Step 3: 求解對(duì)偶問題（6），得到最優(yōu)解 u*

Step 4: 若u*是無窮大，說明u*是極射線的方向。根據(jù)式（7)，添加約束 $u^{*T}(b-By) \geq 0$ 到主問題(12)中

Step 5: 若u*不是無窮大，說明u*是極點(diǎn)。根據(jù)式（11)，添加約束 $u%5E%7B*T%7D(b-By)%20%5Cgeq%20%5Ceta$ 到問題（12）中，并更新 $LB%3Dmax%5C%7BLB%2Cd%5ETy%5E*%2Bu%5E%7B*T%7D(b-By)%5C%7D$ 。（LB記錄了到目前為止我們找到的最優(yōu)的解）

Step 6: 求解主問題（12）得到最優(yōu)解 y*，η*，然后更新 $UB%3Dmax%5C%7BLB%2Cd%5ETy%5E*%2B%20%5Ceta%5E*%5C%7D$ （因?yàn)榇藭r(shí)主問題（12）中的約束還沒有全部被找到，所以這個(gè)解不一定是可行的，但是一定大于等于所有約束都被添加進(jìn)去以后再被求解出的最優(yōu)解，因此它是真實(shí)最優(yōu)值的上界）

Step 7: 返回step2

總結(jié)與拓展（基于邏輯的Benders分解算法）

回顧我們上述的討論，我們可以發(fā)現(xiàn)，Benders分解算法的核心思想就是我們?cè)陂_始時(shí)提出的兩個(gè)想法：

很多y會(huì)使得問題無解——我們利用原問題無解的充要條件，使用對(duì)偶問題的極射線添加約束。
我們已經(jīng)遍歷過的點(diǎn)為問題的最優(yōu)解提供了一個(gè)下界——我們使用對(duì)偶問題的極點(diǎn)添加約束。

綜上所述，我們的目標(biāo)就是通過添加約束的方式不斷縮小y的可行域，并不斷更新原問題的上下界。

在很多現(xiàn)實(shí)的優(yōu)化問題中，我們可以利用實(shí)際問題的獨(dú)特特征，來添加約束，從而達(dá)到更快地縮小y可行域的目的。這個(gè)方法被稱為基于邏輯的Benders分解。

對(duì)基于邏輯的Benders分解算法以及其他Benders分解算法感興趣的同學(xué)，可以參考這篇綜述文章：Ragheb Rahmaniani, Teodor Gabriel Crainic, Michel Gendreau, Walter Rei, The Benders decomposition algorithm: A literature review

Benders分解算法的python實(shí)現(xiàn)

標(biāo)簽：