熊文濤, 雍龍泉
(1.湖北工程學(xué)院 數(shù)學(xué)與統(tǒng)計學(xué)院, 湖北 孝感 432000;2.陜西理工學(xué)院 數(shù)學(xué)與計算機(jī)科學(xué)學(xué)院, 陜西 漢中 723001)
整數(shù)規(guī)劃模型是運籌學(xué)中一類常見的數(shù)學(xué)模型,在現(xiàn)實生活中有著廣泛的應(yīng)用。在許多仿真實驗中,求解整數(shù)規(guī)劃最常見的工具有Matlab和Lingo軟件。然而,Matlab軟件并沒有提供直接求解一般整數(shù)規(guī)劃的內(nèi)部函數(shù)。盡管它提供了求解0-1規(guī)劃的bintprog函數(shù),但其輸入的約束條件要求采用矩陣乘積形式,這對于人們習(xí)慣用求和或任意形式構(gòu)造成的模型,不便于轉(zhuǎn)換成標(biāo)準(zhǔn)的矩陣乘積形式,求解起來非常費時。雖然Lingo軟件能根據(jù)人們計算習(xí)慣將求解一般整數(shù)規(guī)劃問題的模型轉(zhuǎn)換成程序語句,但無法自動多次求解整數(shù)規(guī)劃,僅限于每次求解一個優(yōu)化模型。為方便整數(shù)優(yōu)化模型求解,本文將介紹求解最優(yōu)化模型的Yalmip工具箱。Yalmip是在Matlab軟件的仿真平臺上,采用Matlab的程序設(shè)計語言開發(fā)的一個工具箱。它不僅可以利用Matlab強(qiáng)大的計算能力,調(diào)用Matlab軟件自帶的優(yōu)化工具箱和Matlab語言編寫的其他工具箱中的函數(shù),而且還能將多種優(yōu)化求解器(如Cplex和Gurobi等工具)結(jié)合在一起使用。本文闡述了利用Yalmip工具箱求解整數(shù)規(guī)劃的一般方法,并通過一個具體的實例說明使用Yalmip工具箱在求解整數(shù)規(guī)劃模型的一般方法。
Yalmip是由Lofberg開發(fā)的一種免費的Matlab工具箱[1],其最初目的是為求解控制與系統(tǒng)理論中的半定規(guī)劃和線性矩陣不等式而設(shè)計的。后來經(jīng)過不斷發(fā)展,可以求解更多的優(yōu)化模型,如二階錐規(guī)劃、混合整數(shù)規(guī)劃、正項幾何規(guī)劃、雙線性矩陣不等式以及多參數(shù)線性規(guī)劃和二次規(guī)劃等。Yalmip最大的特色在于能集成許多外部的最優(yōu)化求解器,無論這些求解器是否采用Matlab程序語言編寫的。通常,這些求解器把優(yōu)化問題描述成一個非常緊湊的格式,學(xué)習(xí)、使用起來比較耗時,并且編程時容易出錯。為了克服這一方面的不足,Yalmip提供了一種合適的建模語言和編程接口,以方便不同求解器之間的集成。Yalmip也能調(diào)用Matlab自帶的優(yōu)化工具箱,特別是對于未知變量較多的線性規(guī)劃問題,調(diào)用Linprog函數(shù)不需要轉(zhuǎn)換成矩陣的形式,使用起來非常方便。Yalmip也包含內(nèi)置的整數(shù)規(guī)劃求解器,可以求解一些小型的整數(shù)規(guī)劃問題,并且與外部的整數(shù)規(guī)劃求解器具有較好的兼容性。
Yalmip作為一個Matlab的工具箱,與一般的外部工具箱安裝相同。使用者可以通過該工具箱的網(wǎng)址http://users.isy.liu.se/johanl/yalmip/pmwiki.php?n=Tutorials.Installation下載安裝包并解壓后,可通過如下步驟將對應(yīng)的文件夾放到Matlab路徑里:首先點擊File中的Set Path菜單,然后點擊Add with subfolders菜單,找到解壓好的文件夾Yalmip,最后確定、保存即可。對于外部的一些優(yōu)化求解器,可采用類似的方法放入到Matlab路徑中。
一般的整數(shù)規(guī)劃可表示為:
(1)
其中,f(x)為目標(biāo)函數(shù),gi(x)和hj(x)是定義在Rn上的多元實值函數(shù),決策變量x至少有一個分量為整數(shù)。當(dāng)x全部為整數(shù)時,模型(1)為純整數(shù)規(guī)劃。
特別地,當(dāng)f(x),gi(x),hj(x)均為線性函數(shù)時,此模型(1)稱為整數(shù)線性規(guī)劃,模型(1)可重新改寫為如下的矩陣表示形式:
(2)
其中c是n維列向量,A和Aeq分別是m×n矩陣和l×n矩陣,b是m維列向量,beq是l維列向量。
對于許多實際問題,人們習(xí)慣于用求和形式與任意的符號形式進(jìn)行表示,而不是矩陣的乘積形式,于是模型(2)又可寫為模型(3)的形式,即:
(3)
其中ci,bi,beqs分別向量c,b,beq的分量,aij,aeqsj則分別A,Aeq是矩陣中元素。
為了方便,假設(shè)模型(3)中所有未知變量 全部為整數(shù)。若部分變量不是整數(shù),則只需改變申明變量類型的函數(shù)即可。表1給出了求解整數(shù)線性規(guī)劃的一般程序。
表1 求解整數(shù)線性規(guī)劃的一般程序模式
以文獻(xiàn)[2]中一個電力生產(chǎn)問題作為實例,介紹使用Yalmip工具箱求解整數(shù)規(guī)劃模型的求解方法。為滿足每日電力需求(單位為MW,見表2),有4種不同類型的發(fā)電機(jī)可供選用。所有發(fā)電機(jī)都存在一個啟動成本(在每個時段開始時才允許啟動或關(guān)閉),以及工作于最小功率狀態(tài)時的固定成本,并且如果功率高于最小功率,則超出部分的功率每兆瓦每小時還存在一個邊際成本。另外,每種發(fā)電機(jī)都有一個最大發(fā)電能力,當(dāng)接入電網(wǎng)時,其輸出功率不應(yīng)低于某一最小輸出功率,具體數(shù)據(jù)見表3。與啟動發(fā)電機(jī)不同,關(guān)閉發(fā)電機(jī)不需要付出任何代價。在任意時刻,正在工作的發(fā)電機(jī)組必須留出20%的發(fā)電能力余量,以防用電量突然上升。每個時段應(yīng)分別使用哪些發(fā)電機(jī)才能使每天的總成本最???
表2 每日用電需求(MW)
表3 發(fā)電機(jī)情況
假設(shè)P為發(fā)電機(jī)型號集合,T={1,2,…,NT}為每天的時段集合,Lenj,Demj分別表示第j個時段的長度(小時數(shù))和電網(wǎng)的用電量需求,Pmini,Pmaxi,Availi,Cstarti,Cmini,Caddi為i型發(fā)電機(jī)的最小功率、最大功率、可用的發(fā)電機(jī)數(shù)目、啟動成本、最小功率每小時工作成本、以及每兆每小時的邊際成本。未知變量startij,workij及paddij分別表示在時段j開始工作的型號i的發(fā)電機(jī)數(shù)量,在時段j內(nèi)工作型號i的發(fā)電機(jī)數(shù)量,以及型號i的發(fā)電機(jī)在時段j內(nèi)功率超出其最低功率之上的量。根據(jù)上述問題描述可以得到如下的整數(shù)規(guī)劃模型。
(4)
其中,約束條件(4.1)表示對于每個類型和每個時段的發(fā)電機(jī)(型號為i),其額外功率都不能超過Pmaxi-Pmini;約束條件(4.2)表示每個時段內(nèi)發(fā)電量都應(yīng)滿足用電需求;約束條件(4.3)則表示在不啟動新的發(fā)電機(jī)的前提下,當(dāng)前正在工作的發(fā)電機(jī)的發(fā)電能力必須高于實際用電量的20%;另外,根據(jù)需求,當(dāng)不關(guān)閉發(fā)電機(jī)時,在開始時段時,投入開始工作的發(fā)電機(jī)數(shù)目應(yīng)等于時段j和j-1間工作的發(fā)電機(jī)數(shù)量之差。如果有些發(fā)電機(jī)可能會被關(guān)閉,那么開始工作的發(fā)電機(jī)數(shù)量應(yīng)至少等于此差值,于是得到約束條件(4.4);約束條件 (4.5)表示當(dāng)日午夜和次日凌晨這兩個時段的情形。
針對此模型,可以通過表4的程序語句計算出其最優(yōu)解和最優(yōu)值,得到發(fā)電機(jī)的使用計劃,其結(jié)果見表5。在Windows xp操作系統(tǒng),3.20 GHz Intel Core i5 CPU和 4 GB RAM的PC上求解此整數(shù)規(guī)劃,時間僅用了0.3 s。表5的結(jié)果表明:對于型號1的發(fā)電機(jī),時段1時有3臺在工作,時段2開始時啟動1臺,時段4開始時啟動3臺,時段4之后關(guān)閉4臺,4臺型號2的發(fā)電機(jī)都將連續(xù)工作;對于型號3的發(fā)電機(jī),時段1時有2臺在工作,時段2 開始時啟動6臺,時段6結(jié)束后關(guān)閉4臺,一天結(jié)束時再關(guān)閉2臺;對于型號4的發(fā)電機(jī),時段2開始時啟動3臺,到下一時段時關(guān)閉2臺,然后在下面的幾個時段這2臺反復(fù)開閉,
表4 求解模型(4)的Matlab程序
直到一天結(jié)束時,所有發(fā)電機(jī)都關(guān)閉。每種型號的發(fā)電機(jī)在每個時段的總功率可以通過發(fā)電機(jī)使用數(shù)量、最小輸出功率和額外功率計算,如型號1的發(fā)電機(jī)在第2時段的總功率為4600 MW。
整數(shù)規(guī)劃模型是運籌學(xué)中一類常見的優(yōu)化求解問題,在實際應(yīng)用中有著非常廣泛的應(yīng)用。由于Matlab軟件并沒有提供直接求解一般整數(shù)規(guī)劃的內(nèi)部函數(shù),不便于人們直接對一般整數(shù)規(guī)劃優(yōu)化求解問題。為方便整數(shù)規(guī)劃問題求解,本文利用Matlab平臺,介紹了一種使用Yalmip工具箱求解整數(shù)規(guī)劃問題的一般方法。該方法能利用Matlab的強(qiáng)大計算能力,以Matlab為平臺調(diào)用Yalmip工具箱對優(yōu)化問題進(jìn)行求解,編程方式具有靈活、直觀等優(yōu)點。而且利用Yalmip能夠集成多種優(yōu)化求解器的優(yōu)點,該方法能方便求解其他類型的優(yōu)化問題,如線性規(guī)劃、二次規(guī)劃和一般的非線性規(guī)劃,便于對多個求解器計算的結(jié)果進(jìn)行比較和選擇。以電力生產(chǎn)問題作為實例,說明了使用Yalmip工具箱求解整數(shù)規(guī)劃模型的過程。計算結(jié)果表明了本文方法的有效性。
表5 發(fā)電機(jī)的使用計劃
[參 考 文 獻(xiàn)]
[1] Lofberg J. YALMIP: A toolbox for modeling and optimization in MATLAB[C]//Proceedings of 2004 IEEE International Symposium on Computer Aided Control Systems Design, 2004:284-289.
[2] Gueret C,Prins C,Sevaux M.運籌學(xué)案例[M].北京: 北京林森科技發(fā)展有限公司,2006.