基于DantzigWolfe分解和列生成算法的自存儲收益管理大規(guī)模優(yōu)化模型

基于Dantzig-Wolfe分解和列生成算法的自存儲收益管理大規(guī)模優(yōu)化模型摘要：自存儲行業(yè)正處于高速發(fā)展階段，每年10%的美國家庭租有一個自助存化Dantzig-Wolfe分解技術(shù)、分支定界技術(shù)提出了針對該模型的列生成算法，并用MATLAB程序?qū)崿F(xiàn)了此算法。求解和比較表明，本文提出的列生成算法能準(zhǔn)確解突出。關(guān)鍵詞：自存儲、收益管理、整數(shù)規(guī)劃、列生成算法、D-W分解AbstractSelf-storageisaboomingindustryhasreceivedlittleacademicattentionfar.Inabout1outof10householdscurrentlyrentaself-storageunit.Thisconsidersaself-storagecompany,facingasetofreservationsforstorageunitsahorizonwithrevenuerewards.Thecompanyhastostoragerequeststobeacceptedtomaximizerevenue.Wemodelproblemasabinaryintegerprogramming,anditbyDantzig-Wolfedecompositionandcolumngeneration.WeimplementthisalgorithmbyTheexperimentsthatthisalgorithmcanefficientlyproblem,andcanmuchcomparedwithimprovedsimplexalgorithmanddualsimplexalgorithm.Keywords:Self-storageindustry,revenuemanagement,binaryintegerprogramming,columngenerationalgorithm,Dantzig-Wolfedecomposition1.引言20世紀(jì)60年代中期美國德克薩斯州，隨后迅速蔓延至全美及歐洲。至2008年，美國已有30604家自助存儲公51500大型自助存儲倉庫（self-storage國公民擁有6.89平方英尺自助存儲倉庫面積，將近10%美國家庭擁有一個租賃1995年為65%的增長[1]，也意味著自助存儲已被越來越多的家庭認(rèn)識、接受、并使用。上世紀(jì)80年代末，自助存儲行業(yè)在英國出現(xiàn)，隨后開始了其在歐洲的迅速2007869個。歐洲自助存儲協(xié)會聯(lián)盟（FederationofEuropeanSelfStorageAssociation,FEDESSA）和英國自助存儲協(xié)會（SelfStorageAssociationoftheUnitedKingdom）的官方網(wǎng)站上都提供了該行業(yè)在歐洲的發(fā)展概況。不僅在Yellow，正在迪拜興建一個全世界最大的自助存儲倉庫[1]。從以上數(shù)據(jù)及描述可以看出，自助存儲行業(yè)雖然只出現(xiàn)了40多年時間，可是其市低迷及金融危機更賦予此行業(yè)發(fā)展的一個有利時機。本文研究的問題來源于Shurgard，好地控制了成本，因此，收入是該公司的主要關(guān)注點。Shuragrd公司的倉庫都是收據(jù)公司的經(jīng)驗，大約90%的未來三個月內(nèi)倉庫預(yù)訂需求可以提前知曉。有些客戶大公司目前空置的倉庫類別及各類別倉庫數(shù)目限制。有收入。在現(xiàn)實生活中，該問題有五種情況：（1）只有一種倉庫類別。每份訂單只預(yù)訂該倉庫類別中的一個存儲單元。這是最為簡單的一種情況。（2）只有一種倉庫類別。每份訂單可能預(yù)訂不只一個該倉庫類別的存儲單元。這種情況在一些中小型自助存儲公司中很典型。（3）有不只一種倉庫類別。每份訂單總共只預(yù)訂某種倉庫類別的一個存儲單元。這是最普遍的一種情況。在Shurgard，大部分家庭客戶和一些商業(yè)客戶的訂單都屬于這種情況。（4）有不只一種倉庫類別。每份訂單可能預(yù)訂多個存儲單元，并且這些存儲單元同屬于一類倉庫類別。許多商業(yè)客戶的訂單屬于這種情況。（5）有不只一種倉庫類別。每份訂單可能預(yù)訂多種倉庫類別中的多個存儲單元。這種情況中訂單預(yù)訂沒有任何限制，是最一般的情況。由于第3種情況和第4種情況中的訂單只能預(yù)訂一類倉庫，因此，可以將情況3分解為多個屬于情況14分解為多個屬于情況21和情況23和情況4所屬的問題即可相1和情況2所屬的問題可用多項式算法解決，將其轉(zhuǎn)化為最小成本流問題，算法復(fù)雜度為O((nn)2)，從而解決了問題1和問題也解決了問題3和問題4。由于這四種情況并非本文研究問題，這里對其不再贅述。本文主要研究情況55，算法復(fù)雜度為O(nm+1)表示所有類別倉庫數(shù)量總和?？梢钥闯?，當(dāng)倉庫數(shù)量較多的時候，該問題的動態(tài)規(guī)劃算法將會耗費大量的求解時間[2]。本文提出的列生成算法將從另外一個角度對問題5進行建模并對其作更有效地求解。2.模型2.1模型符號定義在建立問題的模型之前我們首先用數(shù)學(xué)標(biāo)號對問題進行數(shù)學(xué)定義問題的一些符號如下表所示：t時t=間單位。j訂單編號。j=2,…,ni倉庫類別總數(shù)。i=2,…,KM所有倉庫總數(shù)。第i類倉庫數(shù)量。第j份訂單預(yù)訂的第i類倉庫的數(shù)量。Sj第j份訂單貨物的開始儲存時間。Cj第j份訂單貨物的結(jié)束儲存時間。每個第i類倉庫每時間單位的收費價格。vj接受第j份訂單產(chǎn)生的收益。值得注意的有以下幾點：（1）（2）

K∑i=mi=1Kj≤因而∑j≤,i=1

j=2,…,n（3）Sj和Cj分別表示第j份訂單的儲存的開始時間和結(jié)束時間，因此，對于訂單1要求其貨物于第2個時間單位到第4個時間單位之間被儲存，則S1=2,C1=4,因此，其被儲存時間為（4-2+1）=3個時間單位。K K（4）

vj=∑iSj?(Cj?Sj+),i=1

j=1,2,…,n.令Wj=∑i=1vj=Wj?(Cj?Sj+2.2原模型及D-W分解模型(1)問題的整數(shù)規(guī)劃模型及松弛模型為了建立此問題的整數(shù)規(guī)劃模型，本文假設(shè)了兩類變量：第一類變量：xj，其為0-1變量。若xj=1，則表示訂單j最終被接受；否則xj=0。第二類變量：xjt，其也為0-1變量。對于任意

j2,…,n}，t∈{Sj,Sj+,Sj+2,…,Cj}

。若xjt=1，則表示t時刻j訂單

t∈{Sj,Sj+Sj+2,…,Cj}=0。顯然，當(dāng)xj=1時，對任意xjtxj=0時此，

Cj∑xjt=(Cj?Sj+?xjt=Sj則問題的整數(shù)規(guī)劃模型（IP）如下：n n CjZIP=max

∑j=1

vjxj

=max

∑∑j=1t=Sj

WjxjtCj∑xjt?(Cj?Sj+?xj=0,t=Sjn

j=2,…,n

（1）∑xjt≤,j=1

i=K;

t=T

（2）

（3）將約束（3）松弛到[0,1]，則原問題可松弛為如下線性規(guī)劃問題(LP)：n n CjZLP=max

∑

vjxj

=max

∑∑j=1t=Sj

WjxjtCj∑xjt?(Cj?Sj+?xj=0,t=Sjn

j=2,…,n

（4）∑xjt≤,j=1

i=K;

t=T

（5）0≤≤1,0≤≤1

（6）其中約束（4）表示若訂單j被接受，則其對應(yīng)的所有變量均為1；否則，其對應(yīng)的所有變量均為0。約束（5）表示倉庫資源約束，即任意一時刻任意一類倉庫的使用數(shù)不能超過這一類倉庫的數(shù)量。(2)D-W分解原理及問題的D-W分解模型令約束（5）的系數(shù)矩陣為S，X為xjt組成的向量，注意到，S為塊對角矩陣，因此，約束（5）可拆分為如下T組約束條件：?Xj11≤m1≤m,Xj2 i2?2?SSX ≤m,i

titi?t=T,i?T其中，Xjt,t2,…,T}表示所有滿足條件：Sj≤t≤Cj的xjt組成的向tt時刻可以被加工的訂單所組成的變量集 Sij 合。表示對應(yīng)變量的系t數(shù)。另外，原模型的約束條件個數(shù)為n+K*T個，因此，若訂單數(shù)和時間范圍量的時間。結(jié)合以上兩點，本文采用分解（Dantzig-Wolfe分解）對原模型進行分解變形，使之相對于原模型有以下幾點改進：（1）變形后的分解模型將會使約束條件大大減少，為n+T個，雖然分解后變量個數(shù)將會比原模型多出許多，然而，本文采用的列生成算法可以很好地解決變量數(shù)目巨大這一問題；（2）由于可將原問題中的系數(shù)矩陣S分解為T組約束條件（每組由K分解后，列生成算法的子問題為T個獨立的具有K個約束條件的線性規(guī)劃問題，這使得子問題的求解難度大大降低。（3）經(jīng)過分解后的模型下界與原模型相比要更為收緊，因此，只需進行有限步迭代即可獲得原問題松弛模型的最優(yōu)解，節(jié)省了求解松弛模型的時間。分解原理的描述如下：考慮任意LP問題，稱為原問題：x∈S

(7)其中A為m*n矩陣假設(shè)S是有界多面體令1,…,xP是S的極點集合則有：∑p ∑p p∈

λxp,

λλ

≥

（8）problemorMP）：

P∑

p(cxppP∑pPpp)λ=b

（9）∑λp=1λp≥n問題（7）和（9）是等價的，它們有同樣的最優(yōu)解。問題(4.9)的任意可行解λ相應(yīng)問題（7）的可行解x,反之亦然。這就是Dantzig-Wolfe分解原理，簡稱分解。在分解中，用單一加權(quán)約束代替了多面體集約束減少了約束數(shù)目，但因這個多面體的頂點可能非常多（Cm種多列的LP分解的優(yōu)勢體現(xiàn)在它收緊了原問題的下界，而且只需有限步迭代即可獲得原問題的最優(yōu)解。n根據(jù)以上對D-W分解原理的描述，約束（5.1）至（5.T）中的每組約束與約束（6）的可行域分別構(gòu)成一個多面體，設(shè)為,P2,…,PT。它們分別代表著tt t時刻倉庫資源約束的可行域多面體。令Xl(l2,…,L})表示多面體Pt的第t tl lt個頂點，其中表示相對應(yīng)于t時刻的多面體的頂點數(shù)。則Pl lt

內(nèi)的點Xt可表示Xt=

Ltl l ll=1

λtXt,其中

Lt∑l=1

λt=1，λt≥0因此，問題LP可被轉(zhuǎn)化為以下模型，稱為DWLP：T ?n ?llZDWLP=maxll

∑∑?∑

Wjxjt?

λtT

t

l?

j?∑∑jt t j j j

（10）tl

λl?(C?S

+?x

=0,

j=2,…,nl∑l=

t=T

（11）0≤λl≤1,0≤

x≤1

（12）t j至此，原問題被轉(zhuǎn)化為分解后的線性規(guī)劃問題?？梢钥吹?，原問題有tn+K*T個約束條件，而經(jīng)過分解后的問題約束條件為n+T中變量為λl），t不過用列生成算法來求解這個問題，變量的個數(shù)將不再成為困擾。3.算法基具化可以為求解本文的問題模型提供更多的啟發(fā)。(1)列生成算法基本思想與原理列生成算法（ColumnGenerationAlgorithm）的基本出發(fā)點是研究這樣一類大規(guī)模線性規(guī)劃問題：x≥0

（13）nx,A,b是有合適維數(shù)的矩陣，在這些問題里變量數(shù)n遠(yuǎn)遠(yuǎn)超過約束數(shù)m，以至于計算機不能對約束矩陣進行有效的存儲和計算。這種例子如分解的主問題里，每一列對應(yīng)一個基可行解，即可行域的頂點，因為有Cm個頂點，當(dāng)n很大時，即使m較小，總的極點數(shù)即主問題列數(shù)（變量）也是一個非常巨大的數(shù)，處理這類問題，可以使用列生成算法。它只是用主問題列的一些子集所構(gòu)成的問題，稱為限制主問題（Restrictedmasterproblem,RMP）去求解原問題，這就是列生成算法的基本思想。事實上，m，這樣就可以使用較少的列獲得原問題的最優(yōu)解，而其它的列只有當(dāng)需要的時候才去產(chǎn)生。n列生成算法的出發(fā)點是這樣的，由于所研究問題的大規(guī)模特性，使得不可法獲得這種大規(guī)模問題的解。列生成算法框架如圖（4-1）所示。假設(shè)在求解過程的某一步，一個限制主問題（RMP）A'x=b，它是原問題（主問題）的一個縮小問題，僅包含原問題系數(shù)矩陣A的部分列。通過解這個限制主問題，獲得單純形乘子（simplexπ=(π1,…,πm)。然后利用限制主問題的解信數(shù)：moj=cj?∑πiaij≥0,i=1

,n則滿足對偶可行性，表示已獲得原問題的最優(yōu)解。但假使存在某一列s，使σj進行直到獲得原問題的最優(yōu)解?？墒菃栴}在于原問題的列非常多，不可能用枚舉方法找到這樣的列，因此需構(gòu)造一個子問題（subproblem,（pricingproblem）來獲得這樣的列或判斷原問題解的最優(yōu)性。典型地，通過使用這種交互方法，只需有限步迭代即可獲得原問題的最優(yōu)解。綜上，列生成算法一般步聚如下：1：用啟發(fā)式算法（或其它方法）生成初始解，每一個子問題至少生成2：解限制主問題，求得單純形乘子π。3：解相應(yīng)價格問題（對偶問題），判斷最優(yōu)性，如果得到最優(yōu)解，終止；否則，轉(zhuǎn)4。4：選擇具有最小負(fù)判別數(shù)，生成添加列。5：返回2。獲得原問題最優(yōu)解獲得原問題最優(yōu)解解限制主問題：min∑c'jxjst. A'x=bx≥0'為單純形乘子將列s加到限制主問題：A'=A+s解子問題：m mβ=min(c'j?∑πiaij)=c's?∑πiaisi=1 i=1β≥0?NY圖1：列生成算法框架(2)列生成算法在塊對角狀大規(guī)模線性規(guī)劃上的應(yīng)用列生成算法最早是被用來求解塊對角狀大規(guī)模線性規(guī)劃問題（DantzigandWolfe,1961;Gilmoreand式如下：X1+C2X2+?+CpXp+X2+?+Xp=

（14）

?

Xp

===Xi≥i=p 10由分解原理，上述LP問題對應(yīng)的主問題（MP）為：p ∑∑ i i ijz=Xjp j∑∑(=j

（15）

j

≥

j則提出問題列生成算法步驟如下：i i1：利用兩階段算法，對每個i=1,…,p求得兩個初始極點X2i iπs=|π0

]到初始的RMP。解得：1，1，λ

p1，λp2

及單純形乘子πs

，分解πs為兩部分：πs=|π0]，其中π對應(yīng)組合約束，π0對應(yīng)p個凸組合約束，k=k0

+1，這里k0=2是初始極點數(shù)。2：解子問題：=?π)Xi?π0iXi=Xi≥0對所有i=p，獲得極點

Xk，其中πi是π

的第i個分量。如果i 0 0ui≥0,i=1,…,p，那么RMP得最優(yōu)解，也就是原問題的最優(yōu)解，算法停止；否則，令us=:ui<0}，轉(zhuǎn)3。3：產(chǎn)生添加列：s s?AXks s=? ?? ?其中，es是p維單位向量，并且它的s分量為1。4：解新的RMP，得,…,λ1,…,λ,

及π'，返回2。p plp s(3)自助存儲優(yōu)化問題的列生成算法根據(jù)以上兩部分對列生成算法的描述，回到自助存儲優(yōu)化問題的模型，對DWLP中的T個多面體，首先對每個多面體取兩個頂點，X1,X2,X1,X2,…,X1,Xλ1,λ2,λ1,λ2,…,λ1,λ1 1 2 2 T T

1 1 2 2 T TDWLP構(gòu)成初始限制主問題（n+T個約束條件，因此，用單純形法即能求解。得到單純形乘子π=j|πt]，其中，πj表示對應(yīng)約束（4.10）第j個約束條件的單純形乘子；πt表示對應(yīng)約束（4.11）第t個約束條件的單純形乘子。由于原問題的資源約束矩陣可以分解成T個獨立的T個獨立的子問題。由于每個子問題是一個具有K個約束條件的線性規(guī)劃，而解T個具有K個約束條件的線性規(guī)劃的難易程度及耗費時間要大大小于解一個具有SP(t)如下：(SP(t)):n nxσt=Wjxt?∑πjxt?πtn

j=1. ∑Sjxt≤i,j=1x

i=1,2,…,Kt為常數(shù)，因此實際上子問題是求目標(biāo)函數(shù)右邊部份的前兩個表達式之差。由于子問題只有K個約束條件，而在現(xiàn)實生活中，K往往是一個不大的數(shù)。這T個子問題可以很容易地由單純形法求解。取στ

=t}，若στ≤0，則所求的解即為DWLP的最優(yōu)解；否則，生成相應(yīng)的列加入RMP中繼續(xù)求解，直到滿足以上的條件，則所求解為最優(yōu)解。由于DWLP的最優(yōu)解不一定是IP問題的最優(yōu)解，甚至不是IP問題的可行解，因此，接下來需要對非整數(shù)變量進行分支定界。本文采用的分支策略如下：ss j對求得的DP的最優(yōu)解*的第一類變量進行判斷取其中一變量x*使其滿足以下條件：|x*?05=max(|x*?05|),j=,…,n，對其進行分支定界，之所以這樣取的原因是因為此變量最不“整數(shù)性（羅守成20。將s分為兩支s=；s=，分別加進原問題約束條件中。事實上，當(dāng)s=0時，與s相對應(yīng)的第二類變量xt均為0；當(dāng)s=1時，與s相對應(yīng)的第二類變量均為ss j由以上分析，本文得出所研究問題的列生成算法步驟：Stepvj按降序排列，利用貪婪算法求得一原問題的可行解，作為原問題的初始下界。Step2：對每個多面體各取兩個頂點X1,X2,X1,X2,…,X1,X1 1 2 2 T T由于變量均小于等于1，且Sij≤，因此，初始兩個頂點可以選取如下：1X1=11X2={0,1,0,…,0}1?TX1=TTX2={0,1,0,…,0}T以上向量的維數(shù)與頂點向量維數(shù)相同（事實上，每個多面體各取一個初始頂將頂點及相對應(yīng)的λ變量代入DWLP形成初始限制主問題如下：T 2?n ?llZRMP=maxll

∑∑?∑

Wjxjt?

?λttT 2

l?

j?∑∑jt t j j jtll2l

λl?(C?S

+?x

=0,

j=2,…,n∑l=1=

t=T0≤λl≤1,0≤

x≤1t j求得初始解及單純形乘子π=j|πt]。Step3：令πj表示對應(yīng)約束（4.10）第j個約束條件的單純形乘子；πt表示對應(yīng)約束（4.11）第t個約束條件的單純形乘子。同時求解T個子問題，得到στ=t}。其對應(yīng)的子問題的解為SPXτ。Step4：判斷解的最優(yōu)性。若στ≤0，則所求的解即為DWLP的最優(yōu)解；否τ則，生成對應(yīng)的列如下加入RMP中重新求解，其對應(yīng)的變量為λ3：ττ(0,…,0,τ

τ,0,…,0,eτ

)T，其中，e

是一個T單位向量，其第τ個元素為1；之前的部分是一個n維向量，SPXτ對應(yīng)于τ時刻的訂單所在的行。jStepRMP，直至求得DWLP的最優(yōu)解λ*和X*j

（第一類變量j最優(yōu)值）。若X*為整數(shù)，則其亦為原IP問題第一類變量的最優(yōu)解，將λ*還原成jtX*，原問題解畢，程序結(jié)束；否則，轉(zhuǎn)6tStep6：對

X*進行判斷，取其中一變量

x*，使其滿足以下條件：js|x*?0.5max(|x*?0.5|),j=2,…,n，對其進行分支定界：x*=0；x*=1。jss j s sStep7：將兩分支分別加入原問題IP，轉(zhuǎn)2。4.應(yīng)用規(guī)模較小的模型為例子，并詳細(xì)說明如何使用上節(jié)所述算法求解該模型。此模型背景如下：有一自助存儲公司，其總共有2類倉庫，每類倉庫數(shù)量分別有2個和35個訂單需要公司來決定作決策，這5個訂單所包含的時間段為第1個時間單位至第5個時間單位。因此，在此模型中，T==

K===3

，并且本文假設(shè)==5，Sj,Cj,Sij都根據(jù)相應(yīng)的約束在excel里隨機生成，如下表：訂單jSjCjjS2jWjvj1231224510483450154141251202由上表可以得到此問題整數(shù)規(guī)劃模型如下：Z

=+8x2+10x3+56x4+20x5=14(x12+x13)+24+x25)+34+x35)41+x42+x43+x44)51+x52)x1213?1=0x2425?2=0x34+x35?2x3=0x41424344?4=0x5152?5=0x41≤2x41+51≤3

tx1242≤2x12+2x42+2x52≤3x1343≤13+43≤3x2444≤x34+44≤3x25≤x35≤3

t=t=t=t=xj

j2,…,5},t對上述模型進行松弛后，對每個多面體各取初始頂點，為了計算簡便，本例對每個多面體各取一個初始頂點：

X1=(x,x

)T=，1 41 51X1=(x,x,x=，

X1=(x,x=，2

3 X1=(x

,x,x=，=(x

,x分別對應(yīng)著λ1,λ1,λ1,λ1,λ14

5

1 2 3 4 5五個變量，代入以上模型松弛后可得該問題的初始RMP為：1 2 3 4 5maxZRMP=14(λ1+λ1+λ1)+1 2 3 4 52 3 1λ1+λ1?2x2 3 14 5 2λ1+λ1?2x4 5 2t jλ1=1,0≤t j

?2=01 4λ1?2x=1 4?2=0≤2,…,5},j用單純形法求解，得最優(yōu)值：Z*RMP=50,單純形乘子π=(0,0,…,0)。分別解5個子問題：SP1=max14x

+10x

SP2=max14x

+14x

+10x41 51s.t.

s.t.

12 42 52≤2+≤22x41+2≤3

2+2x42+2≤30≤x41,≤1

0≤,,≤1SP3=max14(x

+x)

SP4=max4x

+5x

+14x13 43s.t.

s.t.

24 34 44+≤2+≤22(+)≤3+2≤30≤,≤1

0≤x24,x34,x44≤1SP5=max4xs.t.

25+5x35≤2≤30≤x25,x35≤1解得=SP2*=21,SP3*=21,SP4*=23,SP5*=9,其中最大值為SP4*=2