成品油管道分輸計劃優化的并行模擬退火算法

陳海宏 ，左麗麗，吳長春，李清平

1 中國石油大學(北京)油氣管道輸送安全國家工程實驗室/城市油氣輸配技術北京市重點實驗室，北京 102249

2 中海油研究總院深水工程國家重點實驗室，北京 100028

成品油管道運行與管理的重要內容之一是制定批次計劃?？尚械呐斡媱澬铦M足站場分輸/注入流量、站場分輸/注入條件、管段運行流量、批量和批次次序等約束[1]。對于中石油現有成品油管道，借助于批次計劃模擬軟件，計劃編制人員手工制定一個月的管道批次計劃通常需耗時半天以上。為提高計劃編制人員的工作效率，許多研究人員提出了批次計劃優化模型與算法，并開發了相應的優化軟件[2-9]，明顯縮短了計劃編制時間。針對成品油管道批次計劃優化問題，文獻[10～13]采用離散時間方程建立了混合整數線性規劃(MILP)模型，而文獻[14～19]采用連續時間方程建立了批次計劃優化的MILP模型。與離散模型相比，連續模型的規模更小，計算時間更短，優化結果質量更高[20-22]。當問題規模較小時，通過CPLEX商業優化軟件一般可在較短時間內獲取MILP模型最優解。隨著問題規模的擴大，CPLEX軟件的計算時間將大幅增加，其往往只能獲取MILP模型較優解，有時甚至無法獲取可行解，故CPLEX不適合于管道規模大、計劃周期長、批次數目多、油品種類多的情況。為保證在較短時間內產生較優解，文獻[11]提出采用模擬退火(SA)算法優化單點注入、多點分輸成品油管道分輸計劃。然而，對于大規模問題，SA算法應用效果依舊不佳。為改善SA算法的優化效果，實現SA算法并行計算，本文提出構造分輸計劃新方案的兩階段法，并改進了目標函數。

以單點注入、多點分輸成品油管道為研究對象，已知管道沿線各站場里程、各管段內徑、各分輸站分輸流量上/下限值、各管段運行流量上/下限值、計劃期內管道起點輸入各批次的時間與流量、計劃期開始時刻管道內各批次的油品種類、體積與位置、計劃期內各中間分輸站每種油品的分輸需求體積，在滿足管道運行約束的前提下，優化管道沿線各中間站分輸計劃，包括站場分輸各批次的開始時刻、結束時刻、流量與體積。

1 并行模擬退火算法

1.1 算法總體框架

成品油管道分輸計劃優化的并行SA算法邏輯框架(圖1)主要包含8個要素：目標函數、初始解構造方法、新解構造方法、新解接收準則、初始溫度Tmax、終了溫度Tmin、降溫速率r、Markov鏈長度S。

(1)分輸計劃初始方案構造方法采用空間遞推法(Space Recursive Method, SRM)[11]。按照SRM法，從上游站場至下游站場，依次基于經驗規則編制各分輸站分輸計劃初始方案。

(2)分輸計劃新方案構造過程分為兩個階段：第一階段利用變量鄰域解構造方法局部調整分輸計劃原方案；第二階段基于分輸作業銜接性規則微調第一階段產生的分輸計劃新方案。

(3)在每次迭代過程中，同時分別構造多個新的分輸計劃，并選擇一個最好的新解作為本次迭代的最終解，從而實現SA算法并行計算。若同時存在多個最好的新解，則從中隨機選擇一個。

(4)分輸計劃新方案接收準則為Metropolis準則：若新解優于舊解，則接收新解；否則，按概率接收新解。概率計算公式為：

其中，prof為差解被接收的概率；fn為新解目標函數值；fo為舊解目標函數值；T為SA算法的當前溫度。

1.2 目標函數

以計劃期內管道沿線所有管段運行流量波動幅度之和最小作為目標函數。

其中，JS為站場總數；KS為時間節點總數；g1為計劃期內各管段運行流量總體波動幅度之和，m3/h；g2為計劃期內各管段運行流量瞬時波動幅度之和，m3/h；cpj為管段j權重系數，用于避免運行流量的波動集中在某些管徑較大、里程較長的管段上，建議取值為管段j體積容量與管道總體積容量的比值；csj,k為0、1變量，csj,k=1 表示時段(tk-1，tk)和(tk，tk+1)內管段j運行流量不相等，csj,k=0則表示運行流量相等；QDj,k為時段(tk，tk+1)內管段j運行流量，m3/h；QAj表示計劃期內管段j運行流量的平均值，m3/h。計劃期劃分方法：選取計劃期內管道沿線各站場執行分輸/注入作業的開始時刻與結束時刻作為時間節點，將計劃期劃分為若干時段。

1.3 分輸計劃初始方案構造方法

為降低成品油管道分輸計劃優化問題的求解難度，假設計劃期內任一中間站對同一批次至多執行一次分輸作業。實際上，目前中石油成品油管道分輸計劃的編制基本上符合該假設?；谠摷僭O，采用SRM法[11]構造分輸計劃初始方案。按照SRM法，從上游站場至下游站場，依次基于經驗規則構造管道沿線各分輸站分輸計劃初始方案。確定某站場是否分輸某批次的原則：計劃期內該批次可以到達該站場、該站場對該批次所含油品存在分輸需求。確定某站場分輸各批次體積的原則：按計劃期內批次過站體積的比例將該站場分輸每種油品的需求體積分配至各個批次。確定某站場分輸某批次流量的原則：優先選取分輸流量允許范圍內的中間值。確定站場分輸某批次時間的原則：為保證初始方案有較大的調整空間，將每次分輸作業安排在批次過站時間窗的中間時段。

1.4 分輸計劃新方案構造方法

分輸計劃新方案構造過程(圖2)可分為兩個階段。第一階段為采用變量鄰域解構造方法調整分輸計劃原方案；第二階段為基于分輸作業銜接性規則微調第一階段產生的分輸計劃新方案。

圖1 成品油管道分輸計劃優化的并行SA算法示意圖Fig. 1 Schematic diagram for optimizing delivery schedules of a multiproduct pipeline by parallel SA algorithm

(1)基于變量鄰域解構造方法調整分輸計劃原方案

A 隨機選擇一個中間分輸站sj。

B 隨機選擇一個批次si，其油品類型為sp。

C 若計劃期內被選中站場未分輸被選中批次，則增加一個新的分輸作業；否則，從分輸時間、分輸流量、分輸體積3個變量中隨機選擇一個變量，利用變量鄰域解構造公式調整被選中分輸作業的執行計劃。變量鄰域解構造公式[23]如下：

其中，γ為決策變量，Xa≤γ≤Xb；γ’表示γ的鄰域解；T表示SA算法的當前溫度；σ為0～1隨機數；sgn為符號函數。

(2)基于分輸作業銜接性規則微調分輸計劃

為提高分輸計劃新方案的質量，第一階段基于變量鄰域解構造法調整被選中分輸作業后，嘗試將該分輸作業與其它中間站執行的分輸作業進行銜接。若第一階段調整的分輸作業或者新增分輸作業的執行開始/結束時刻與相鄰某站場所執行分輸作業的開始/結束時刻的時間間隔小于φ小時，則銜接這兩個分輸作業。φ一般可取5 h。

圖2 分輸計劃新方案構造流程Fig. 2 Flowchart for generating a new delivery schedule

1.5 參數設置

根據SA算法理論可知，算法迭代總次數由外循環與內循環共同控制。初始溫度、終止溫度與降溫函數可用于控制外循環執行次數，Markov鏈長度可用于控制內循環執行次數。根據文獻[11]中SA算法參數設置方法進行參數設置，具體如下：

初始溫度：Tmax= 1000 ℃。

終止溫度：Tmin= 1 ℃。

Markov鏈長度：S =1000。

降溫函數：T’=r×T。其中，r表示降溫速率，0＜r＜1。需根據具體測試情況，綜合比較優化結果的質量與計算耗時確定r最佳取值。

并行線程數PN：需根據具體測試情況，綜合比較優化結果的質量與計算耗時確定PN最佳取值。

2 算例分析

以某成品油管道PA為例，驗證并行SA算法優化管道分輸計劃的有效性。管道PA總長2000 km、內徑544.8 mm、沿線共建有8個站場(編號分別為E0～E7)；管道沿線各管段運行流量允許范圍為400～1200 m3/h，沿線各中間分輸站分輸流量范圍見表1。管道PA測試的計劃周期為0～670 h，已知計劃期內管道起點輸入計劃(表2)、計劃期內管道沿線各中間站每種油品分輸需求體積(表3)和初始時刻(0 h)管道內各批次的油品種類、體積及位置(表4)，利用并行SA算法優化管道PA分輸計劃。

2.1 冷卻系數r的影響

影響SA算法穩定性的最主要因素是冷卻系數r。一般來說，r越大，SA算法優化結果質量越高，優化耗時越長，算法穩定性也越高。根據具體管道測試情況，綜合優化結果質量與耗時可確定r推薦值。不考慮并行計算，設定不同的r值，應用SA算法優化管道PA分輸計劃(圖3)。當0.05≤r≤0.85時，隨著r增加，優化結果質量呈明顯上升趨勢。當0.85≤r≤0.95時，隨著r增加，優化結果質量提升幅度較小，但計算耗時顯著增加。因此，將r推薦值設定為0.85。

表1 管道PA沿線中間分輸站分輸流量允許范圍Table 1 Delivering flow rate limits of intermediate stations along pipeline PA

表2 計劃期內管道PA起點輸入計劃Table 2 Inputting schedules for the initial station of pipeline PA during the scheduling horizon

表3 計劃期內管道PA中間分輸站對各油品的分輸需求體積Table 3 Product demands of intermediate stations along pipeline PA during the scheduling horizon

表4 初始時刻管道PA內各批次油品的種類、體積及位置Table 4 Product type, volume and location of every batch in pipeline PA at the initial moment

2.2 并行線程數PN的影響

令r=0.85，設定不同的并行線程數目PN，應用并行SA算法優化管道PA的分輸計劃(圖4)。

(1)當1≤PN≤8時，隨著PN增加，并行SA算法計算耗時呈線性增長。并行SA算法采用的并行策略為：每次迭代時，同時產生若干個新解，并選擇質量最高的解作為本次迭代的最終解?？偟螖挡蛔?，但每次迭代時由于增加了各線程之間的等待耗時，使得并行SA算法的計算耗時隨PN增加而增大。當PN=8時，計算時間約為6 min，滿足單條管道分輸計劃編制工作允許計算時間的要求。

(2)當1≤PN≤6時，隨著PN增加，優化結果質量呈明顯上升趨勢，目標函數值由5010 m3/h降低至4490 m3/h。當6≤PN≤8時，隨著PN增加，優化結果質量提升幅度較小，目標函數值由4490 m3/h降低至4480 m3/h，但計算耗時卻由325 s增加至386 s，增長約19%。

綜合考慮優化結果質量與計算耗時，將PN推薦值設定為6。與PN=1相比，當PN=6時，優化結果目標函數值降低約10%。

圖3 不同冷卻系數下，SA算法優化管道PA分輸計劃的結果Fig. 3 The optimization results for pipeline PA provided by SA algorithm with different cooling factors

圖4 當r=0.85時，不同并行線程數下，并行SA算法優化管道PA分輸計劃的結果Fig. 4 The optimization results for pipeline PA provided by parallel SA algorithm with different parallel threads when r=0.85

2.3 優化結果展示

令r=0.85、PN=6，20次應用并行SA算法優化管道PA分輸計劃，平均耗時約5.4 min，優化結果目標函數平均值為4490 m3/h。圖5為某次利用并行SA算法產生的管道PA分輸計劃，優化結果目標函數值為4386 m3/h。圖5中P1、P2、P3分別表示0#柴油、92#汽油、95#汽油。計劃期內，中間分輸站E1—E2、E3—E4、E5—E6之間所執行的分輸計劃具有很好的銜接性，保證了管道沿線各管段流量運行平穩。

2.4 并行SA算法與CPLEX對比

文獻[10]采用離散時間方程建立了成品油管道分輸計劃優化的MILP模型，并利用CPLEX求解模型?；谖墨I[10]建立的MILP模型，利用CPLEX優化管道PA的分輸計劃。實際生產中，單條管道分輸計劃編制工作允許的計算時間一般不超過10 min。將計算時間設定為10 min，CPLEX成功產生了管道PA的分輸計劃，優化結果目標函數值為6021 m3/h。CPLEX與并行SA算法的計算過程如圖6所示。與CPLEX相比，并行SA算法優化結果的目標函數值降低約27%，計算耗時降低約46%。從優化結果來看，并行SA算法明顯優于CPLEX。

圖5 當r=0.85、PN=6時，應用并行SA算法產生的管道PA分輸計劃Fig. 5 The delivery schedules of pipeline PA provided by parallel SA algorithm, when r=0.85 and PN=6

圖6 CPLEX、并行SA算法優化管道PA分輸計劃的迭代過程Fig. 6 Computational processes about optimizing delivery schedules of pipeline PA using CPLEX and parallel SA algorithm

3 結束語

通過在SA算法中引入并行計算，建立了適用于成品油管道中間站分輸計劃優化的并行SA算法，其綜合性能不僅優于非并行SA算法，而且優于CPLEX商業優化求解器。不同于非并行SA算法，并行SA算法最大的特點是：從當前迭代點出發同時構造并檢驗多個鄰域解，選擇其中最好者作為本次迭代產生的新迭代點。對于某成品油管道PA，將并行計算線程數取為6個，并行SA算法可以保證在計劃編制工作可接受的計算時間內(不超過10 min)產生較優的成品油管道分輸計劃?；诓⑿蠸A算法產生的管道PA分輸計劃具有很好的銜接性，符合計劃編制人員編制批次計劃的要求，保證了管道沿線各管段流量運行平穩。與非并行SA算法相比，并行SA算法的優化結果目標函數值降低約10%。與CPLEX優化求解器相比，在10 min的計算時間內，并行SA算法的優化結果的目標函數值降低約27%。