張志華,白金鋒,周 洋,李紹軍

(1.遼寧科技大學 化學工程學院,遼寧 鞍山114051;2.遼寧科技大學 電子與信息工程學院,遼寧 鞍山114051;3.華東理工大學 信息科學與工程學院,上海200237)

1 前 言

煤氣化技術是實現煤炭清潔高效利用的重要途徑。對煤氣化過程的優(yōu)化，可以挖掘和提升氣化裝置生產潛能，對工業(yè)實際應用具有重要意義。Aspen Plus[1-2]和Fluent 等仿真軟件[3-4]常用于煤氣化反應過程研究，但其高精度尋優(yōu)過程卻是計算昂貴的問題。目前，代理模型輔助進化算法(surrogate-assisted evolutionary algorithms，SAEAs)已經廣泛用于昂貴問題優(yōu)化中，在化工過程仿真模擬優(yōu)化中逐漸被研究者關注[5]?；诖砟Ｐ偷膬?yōu)化，訓練樣本點選取是一個關鍵因素。因此，在有限計算成本下采樣策略的構建成為SAEAs的研究重點。

最簡單、直接的采樣策略是選擇當前代理模型預測的最好解加入訓練數據集[6]，改善當前最優(yōu)解的精度，也就是局部開發(fā)[7]。局部開發(fā)的采樣準則更多集中在局部區(qū)域尋找最優(yōu)解，如采用信賴域法[8]、序列二次規(guī)劃法[9]等。另一種采樣策略是全局探索，通過搜尋未探索區(qū)域選擇當前代理模型的最不確定點，跳出局部最優(yōu)，有效提高代理模型精度[7]。在最不確定的樣本點選擇中，常以均方誤差(mean squared error,MSE)[10-11]、采樣點和已知點之間距離函數[12-15]最大為目標選擇采樣點。Xu 等[16]基于Voronoi 圖，根據誤差和距離選擇樣本點，對低維度(維度d≤10)問題獲得較好的最優(yōu)解。Wang 等[17]基于委員會的主動學習方法，以不同代理模型之間預測偏差最大的樣本點作為最不確定的點，有效地解決了復雜高維問題(d=20，30)，但在低維度問題(d=10)上不具優(yōu)勢。

為了在有限的精確函數評價(function evaluations,FEs)內獲得更好的最優(yōu)解，本文提出了基于交叉驗證和Voronoi 圖的代理輔助進化算法(cross validation Voronoi surrogate-assisted evolutionary algorithm,CVVSAEA)。CVVSAEA 通過留一交叉驗證(leave-one-out crossvalidation,LOOCV)和Voronoi 圖選擇不確定性較高的樣本點，采用全局集成模型和局部代理模型輔助粒子群優(yōu)化(particle swarm optimization,PSO)算法尋優(yōu)。通過對中等規(guī)模復雜問題的測試說明了本文方法的有效性，并將其應用到煤氣化過程中。

2 代理模型及PSO算法

2.1 代理模型

常用的代理模型有Kriging[18]、徑向基函數(radial basisfunction,RBF)[8]、多項式回歸響應面(polynomial regression surface,PRS)、支持向量回歸以及多個單獨代理模型構成的集成代理模型[17]。Kriging 模型是一種采用插值方法的無偏估計模型，除了可以預測某一特定點處的值外，還可以預測方差，方差可以作為一種不確定性衡量。RBF是一種采用多變量數據插值方法的模型，通過基函數加權求和近似函數，具有較好的局部逼近能力。

不同的代理模型適用于求解不同特性的問題，集成代理模型可以彌補單個代理模型的缺點，更有效解決復雜黑箱優(yōu)化問題[11]。集成代理模型一般采用單個代理模型權重求和的方法建立，其形式如下[7]：

2.2 PSO算法

PSO算法由于其簡單且尋優(yōu)能力強在工程優(yōu)化設計領域中得到了廣泛應用[17]。本文采用一種具有慣性權重因子的改進PSO算法，粒子每代的速度v g和位置x g更新公式如下[17]：

式中：xpbest為單個粒子最好解；xgbest為整個群體中所有粒子的最好解；g為當前迭代代數；r1和r2為[0,1]的隨機數；c1和c2為學習因子，這里參照Wang 等[17]數據，取c1=c2=1.494 45。?g為慣性權重因子，這里?g從0.9到0.4之間隨著迭代代數g線性下降，具體形式如下[20]：

式中：gmax為最大迭代代數。

3 CVVSAEA的采樣策略

以基于委員會主動學習的代理模型輔助粒子群優(yōu)化(committee-based active learning for surrogate-assisted particle swarm optimization,CAL-SAPSO)[17]為基本框架，本文分別從尋找最不確定點、全局最優(yōu)點、局部最優(yōu)點3個方面選擇樣本點，在全局搜索和局部搜索之間進行采樣準則轉換。

假設已知樣本點X=[x1,…,x i,…,x t]T，x i∈d? ，d為每個樣本x i的維度，t為當前集合中樣本個數，對應精確函數值Y=[y(x1),…,y(x i),…y(x t)]T=[y1,…,yi,…,yt]T，構成了訓練數據集(X,Y)t={(x1,y1),…,(x i,yi),…(x t,yt)}T。

3.1 最不確定點x U 選取

方差最大點xσmax具有最大的不確定性，xσmax周圍的點相對其他區(qū)域也應具有較大的不確定性。在xσmax周圍選擇點加入已知樣本中，則會極大地提高模型的精度。采用Kriging 模型預測，其相關函數選擇常用的高斯函數。

(2)空間分割和xσmax區(qū)域的界定

根據Voronoi 圖方法，由設計空間中任意兩個已知樣本點垂直平分線為界把整個設計空間分割成t個單元{C1,C2,… ,Ct}，每個Ci包含一個已知點x i[16]。對于多維問題，如果要精確計算出形狀不規(guī)則的單元Ci的邊界，需要比較復雜的計算過程。采用Monte Carlo方法可以簡單有效地對Voronoi 圖中各個單元進行分割，近似確定單元Ci的邊界范圍[16]。

首先，在整個設計空間內，產生大量均勻分布隨機點xr，數量為N=t×d×ωr，ωr為隨機點生成數量的權值；把xr歸入距離最近已知點x i的單元Ci中，由包含xσmax點的單元Cσmax中所有隨機點近似描述單元邊界。Cσmax中隨機點相對其他單元具有較大方差，在其單元內選取采樣點可以提高模型近似精度。

(3)采樣點xU選擇

單元Cσmax內的隨機點xr在模型構建上具有和xσmax相似的影響作用，直接從這些隨機點中選取采樣點更加簡單有效。為了增加采樣點的多樣性，在單元Cσmax內選擇距離xσmax最遠的隨機點。為了提高模型的近似精度，還要考慮模型預測的偏差大小。因此，建立了以距離和偏差為目標的函數fU(xr)，

式中：xr為單元Cσmax內的隨機點，?r( )y x為代理模型預測xr的函數值，為了減少距離差和函數值差范圍不同的影響，對求得差值分別進行最大最小歸一化處理，norm 表示歸一化處理。

滿足maxfU(x)的解x即為最不確定的采樣點xU。xU描述如下：

xU選取流程如圖1所示。

圖1 x U 選取的流程圖Fig.1 Flowsheet of selecting x U

3.2 當前全局最優(yōu)點x Gb 選取

全局模型采用Kriging、RBF、PRS構成集成代理模型?G( )y x。根據LOOCV的誤差，獲得各個代理模型的權值[21]。權值ωk(x)形式為

3.3 局部最優(yōu)點x Lb 選擇

3.4 搜索區(qū)域轉換

3.5 CVVSAEA 具體流程

圖2 CVVSAEA 流程圖Fig.2 Flowsheet of CVVSAEA

CVVSAEA 流程圖如圖2所示，具體步驟如下：Step1：通過試驗設計方法在設計空間中得到m個樣本的初始數據集(X,Y)m，設定最大函數評價次數tmax；Step2：全局搜索

Step2.1：采用Kriging 模型，通過LOOCV 和Voronoi圖，找不確定性較大區(qū)域，根據式(7)，選擇出最不確定點xU，并計算其精確函數值yU，保存到當前數據集(X,Y)t；

Step2.2：由已知數據集(X,Y)t建立Kriging+RBF+PRS模型?G( )y x。根據式(9)，采用PSO算法找到全局最優(yōu)點xGb，計算精確函數值yGb，保存到(X,Y)t，根據函數性能改善情況進行全局搜索或者局部搜索；

4 數值實例的實驗結果及分析

為了驗證所提出算法的性能，把CVVSAEA 和先進的代理模型輔助優(yōu)化算法CAL-SAPSO[17]、代理模型加權平均法之一(one of the weighted average method，WTA1)[19]、代理模型輔助進化策略之期望改進(metamodel assisted evolution strategieswith expected improvement，MAES-ExI)[22]進行了比較。

4.1 參數設置

本文采用拉丁超立方采樣(latin hypercubesampling，LHS)生成初始樣本點，取初始樣本個數m= 5d。函數評價最大次數tmax=11d，作為整個算法終止條件。其中，隨機點生成權值ωr=100，全局和局部的轉換閾值ε=1.0×10-3，top=20%，局部代理模型采用RBF，其基函數為適合擬合復雜函數的多重二次函數。PSO 算法中初始種群大小為100個，最大迭代次數為100代。為了保證公平性，PSO 在其他的算法中設置相同。對廣泛采用的5個函數維度d=10 和d=20 情況分別進行了測試，函數信息如表1所示。

4.2 數值結果

為了防止隨機性，每個算法獨立運行30次。CVVSAEA、CAL-SAPSO、WTA1、MAES-ExI算法的平均最好函數值及方差和平均運行時間如表2所示，收斂曲線如圖3所示。由于前5d個FEs用于初始訓練數據的評價，因此圖中只從5d個FEs開始顯示。

表1 測試函數Table1 Test functions

表2 CVVSAEA、CAL-SAPSO、WTA1、MAES-ExI平均最好函數值(均值±方差)和平均運行時間的比較Table 2 Averagebest function values(AVG±STD)and average runtime by CVVSAEA,CAL-SAPSO,WTA1 and MAES-ExI

圖3 CVVSAEA、CAL-SAPSO、WTA1和MAES-ExI 對測試函數的收斂曲線(d = 10,d = 20)Fig.3 Convergence curves of CVVSAEA,CAL-SAPSO,WTA1 and MAES-ExI on the test functions of d = 10,d = 20

對表2中的平均最好函數值進行了Friedman 檢驗[23]，p-value=7.9×10-5<0.05，表明4種算法有顯著差異。進一步采用Nemenyi 后續(xù)檢驗[23]比較各算法之間的差異，CVVSAEA 性能明顯好于MAES-ExI算法，同時CVVSAEA 的平均序值也比CAL-SAPSO和WTA1好。

Ellipsoid 函數是一個單模態(tài)問題。CVVSAEA 性能最好，尤其是d= 10，基本得到了真實最優(yōu)值0，遠遠好于其他算法。CVVSAEA 采樣策略更易跳出局部最優(yōu)，比較適合解決單模態(tài)優(yōu)化問題。Griewank函數是一個有很多規(guī)律分布局部極小點的多模態(tài)問題。在d= 10時，CVVSAEA 明顯好于其他算法。但在d= 20時，未能取得好于CAL-SAPSO的最優(yōu)解。由于Griewank 具有大量局部極小點，在非常有限訓練數據下，建立的代理模型可能未能找到最優(yōu)解附近足夠的樣本信息。Schwefel 函數是一個典型的欺騙復雜問題，有很多局部極小值、全局最優(yōu)值與最近的局部最小值距離很遠，很難跳出局部最優(yōu)。在d=10下，CVVSAEA 好于其他3種算法，在d= 20下，略差于WTA1，但好于CAL-SAPSO和MAES-ExI。WTA1只從全局尋找最優(yōu)解，不去探索不確定性和局部區(qū)域，因此具有更多全局搜索的機會。Rosenbrock函數是一個全局最優(yōu)值位于一個狹長平滑拋物線形谷底的多模態(tài)問題。在d=10，d= 20條件下，CVVSAEA 尋優(yōu)性能最好，加入最不確定點和局部最優(yōu)點的CVVSAEA、CAL-SAPSO比僅取全局最優(yōu)點的WTA1更好。Rastrigin 函數與Griewank 函數有些相似，都是有很多規(guī)律分布的局部極小點，但它更復雜。在d=10下，CAL-SAPSO算法在搜索后期略超過CVVSAEA，得到更好的解。Rastrigin 函數由于具有大量局部極小點，在有限次FEs內很難跳出局部最優(yōu)，可能是ωr略小未有效探索不確定點或者局部取top點有些過多造成的。隨著維度增加在d=20下，CVVSAEA 尋優(yōu)性能好于其他3個算法。除了F9外，采用集成模型的CVVSAEA、CAL-SAPSO以及WTA1均取得比僅用Kriging 模型的MAES-ExI 好的解。

總體來說，基于集成模型的采樣策略具有更好的魯棒性，最不確定樣本點加入可以更好地提高模型精度，CVVSAEA 的取點策略要好于其他3個算法。從表2可以看出，在有限的FEs內CVVSAEA 在10種情況下7次求得解最好，其余均排序在第2位；在運行時間方面，MAES-ExI算法雖然最快，但求得的優(yōu)化解最差，CVVSAEA 運行時間在d= 10 時僅次于MAES-ExI 算法，在d=20 總體略差于CAL-SAPSO算法。因此，綜合考慮算法的求解精度和收斂速度，CVVSAEA 具有更好的尋優(yōu)性能。

5 煤炭氣化過程的優(yōu)化

煤氣化過程是通過煤炭與氣化劑(水蒸氣H2O、CO2、空氣、O2等)在高溫下經過一系列復雜的物理化學過程，產生以CO、H2等為主的清潔合成氣的過程。工業(yè)上常見的氣化爐有Lurgi 爐、英國燃氣Lurgi爐(British Gas-Lurgi，BGL)、高溫溫科勒(high temperaturewinkler，HTW)爐和Texaco爐等。氣流床氣化技術是目前煤氣化技術發(fā)展的主要方向[24-25]，本文以氣流床煤氣化爐Texaco為研究對象，CVVSAEA 作為一種優(yōu)化求解策略對其他氣化爐同樣適用。首先經過研磨，滿足粒度要求的煤粉先與水混合成水煤漿，然后水煤漿和O2同時從Texaco氣化爐頂部噴入，依次發(fā)生煤熱解、揮發(fā)分燃燒和焦炭氣化反應，產生以CO和H2為主的合成氣，接著合成氣和熔渣及未反應的碳離開反應區(qū)后進入爐子底部的激冷室，冷卻后的合成氣和飽和水蒸氣從激冷室上部排出氣化爐，熔渣經激冷固化后被分離出來最終進入渣罐，未反應的碳進入黑水處理系統。Texaco爐氣化過程模擬流程如圖4所示。

圖4 Texaco煤氣化過程模擬流程圖Fig.4 Simulation flowsheet of Texaco coal gasification process

在煤炭氣化過程中，有效氣產率是一個重要的評價指標，代表了單位質量的煤產出的有效氣產量，有效氣產率主要受氣化原料組成及性質、進料量、氧煤比、蒸汽煤比、氣化溫度、氣化壓力等工藝參數影響。根據Texaco的Aspen Plus模型[1]，以Illinois No.6煤為研究對象，氣化壓力維持在2 431 800 Pa不變，選取干煤粉進料流量及其溫度、氧煤比和O2溫度、蒸汽煤比和水蒸氣溫度6個參數為操作變量，各個操作變量的取值范圍根據文獻[26]確定，最大化有效氣(CO+H2)產率為優(yōu)化目標。采用LHS方法生成5d個初始樣本點，調用Aspen Plus仿真模型計算其對應的有效氣產率，得到初始訓練數據集。仿真模型調用最大次數tmax=11d。隨機點生成權值ωr=100，轉換閾值ε=1.0×10-3，局部取函數值排序top=30%，局部代理模型采用Kriging，全局模型為Kriging+RBF+PRS。PSO算法初始種群為100個，最大迭代次數為100代。

從圖5中可以看出，大約從第40次FEs開始，CVVSAEA 找到的優(yōu)化解開始好于其他3種算法，具有更好的尋優(yōu)性能。分析表3中數據，與其他3種算法相比，CVVSAEA 使有效氣產率分別提高了0.32%、0.16%、0.48%，對煤炭氣化工業(yè)流程優(yōu)化具有一定的應用價值和指導意義。

圖5 有效氣產率優(yōu)化過程Fig.5 Optimization process of effective syngas yield

表3 CVVSAEA、CAL-SAPSO、WTA1和MAES-ExI優(yōu)化結果比較Table 3 Comparison of optimization results by CVVSAEA,CAL-SAPSO,WTA1 and MAES-ExI

6 結 論

基于LOOCV 和Voronoi圖，本文提出了一種選取最不確定點的采樣準則優(yōu)化代理模型，進而更加有效輔助進化算法解決其優(yōu)化時計算昂貴問題。界定不確定性最高的區(qū)域，綜合考慮多樣性和性能指標選點，且無需進化算法求解，快速精準地改善了代理模型精度。采用集成模型預測具有更穩(wěn)定的近似能力，選擇較好點建立局部模型加快了收斂過程。全局搜索和局部搜索的自適應切換達到了平衡全局探索和局部開發(fā)的目的。通過不同模態(tài)不同維度下10個復雜測試問題和煤氣化過程實例，發(fā)現在有限昂貴函數評價限制內，CVVSAEA 具有更好的尋優(yōu)性能，在煤氣化過程優(yōu)化中應用前景較好。盡管CVVSAEA在昂貴函數評價上有一定優(yōu)勢，但同時也存在一些不足。如問題維度過高時，Kriging 模型計算成本較高，如何在高維問題下提高其計算效率是今后需要進一步深入研究的方向。