延遲焦化加熱爐熱效率的機理建模與實時估計應用

黃健，趙眾

（北京化工大學信息科學與技術學院，北京100029）

引 言

延遲焦化裝置加工的是劣質重質原料油,具有高溫操作的特點,裝置相互關聯嚴重[1]。在延遲焦化裝置中，焦化加熱爐屬于其中的關鍵設備，焦化加熱爐的燃料氣消耗量超過總耗能77%[2]，在衡量加熱爐操作先進性方面，熱效率是其中一項關鍵指標[3]，然而加熱爐熱效率難以在線測量，實時估計更加困難。影響其熱效率的因素有很多[4]，出口排煙的溫度、爐膛煙氣的氧含量、入爐的空氣量都會對熱效率產生影響。實際操作中熱效率主要依靠離線分析[5]得到，通過離線分析數據有正平衡法與反平衡法兩種方式計算[6]。正平衡法通過直接測量加熱爐輸入熱量和輸出熱量計算得到熱效率。而對于反平衡計算方法，則是通過測試和計算加熱爐各項熱損失。國內外主要通過反平衡法計算。粒子濾波要追溯到20 世紀40 年代提出的Monte Carlo 方法，Hammersley 等［7］基于貝葉斯采樣估計提出順序重要性采樣濾波思想，其基本原理是尋找一組在狀態空間傳播的隨機粒子（樣本）描述系統的狀態，通過Monte Carlo 方法處理貝葉斯估計中的積分運算，從而得到系統狀態的最小均方差估計。當粒子數量趨于無窮時可以逼近服從任意概率分布的系統狀態。Gordon 等[8]通過重采樣技術的引入有效地解決了樣本權值退化問題，提出了粒子濾波算法。近年來，在現代信號處理、通信、人工智能、生物信息學、計算機視覺、自主移動機器人、目標跟蹤及統計學等領域不斷發展。粒子濾波的研究在中國最早見于2003 年，袁澤劍等[9]提出了一種性能優于標準粒子濾波的高斯-厄米特粒子濾波器，莫以為等[10]把粒子濾波用于混合系統狀態監測與診斷。隨后，粒子濾波得到了國內學者廣泛的關注并被嘗試應用于目標跟蹤、組合導航和通信等領域，但需要較多的粒子數和計算時間仍然是亟待解決的問題。本文結合管式加熱爐物料，熱量平衡條件，構造了加熱爐的狀態空間預測模型并提出了一種基于粒子濾波的聯合估計方法對模型預測值進行實時估計，并分析了方法的收斂性。在中國石化塔河煉化大型延遲焦化裝置的應用結果證實所提方法的可行性和有效性，為實現加熱爐熱效率實時優化控制奠定了基礎。

1 工藝流程簡介及需求分析

中國石化塔河煉化延遲焦化裝置采用的是三爐六塔結構，工藝主要包括三部分，第一部分是加熱爐焦炭塔、第二部分是分餾塔、第三部分是穩定吸收部分。原料渣油經一系列換熱器預熱至280℃后，進入分餾塔下段的換熱洗滌區，與來自焦炭塔頂的420℃高溫油氣接觸換熱，油氣中的重蠟油以上的組分被冷凝下來作為循環油與原料油一起流入分餾塔底，與原料油匯在一起成為350℃焦化油，塔底焦化油分四路進入加熱爐加熱至500℃。加熱后的高溫油氣，通過四通閥進入焦炭塔進行反應。高溫焦化油在焦炭塔內發生裂解、縮合等一系列反應，生成反應油氣和焦炭。焦炭由下至上聚結在焦炭塔內，反應油氣由焦炭塔頂逸出，進入分餾塔換熱塔板下部，進入分餾塔的焦化油氣與原料進行接觸換熱，循環油流入塔底，換熱后油氣上升進入分餾段，從下往上分餾出蠟油、柴油、粗汽油和富氣[11]。其中加熱爐作為加熱原料的龍頭部分，其在加熱介質的過程中，加熱爐熱效率是加熱爐的一項重要的操作指標[12]。

2 加熱爐熱效率狀態空間模型

管式加熱爐是將爐管中通過的油料(或其他介質)加熱到所需的溫度,以滿足生產工藝的需要[13]。工業管式加熱爐是石化行業中必不可少的裝置。它通過燃料燃燒產生的火焰與高溫煙氣，對爐管內高速流動的原料進行加熱，使其達到工藝規定的溫度。被加熱介質在對流室內與煙氣進行對流換熱，在輻射室內通過火焰或者高溫煙氣進行輻射傳熱。管式加熱爐的結構如圖1所示。

圖1 工業管式加熱爐Fig.1 Industrial tubular heating furnace

根據別洛康方法[14]及微元法，有：

穩態時，煙氣對輻射室爐墻的熱射熱速率與輻射室爐墻對輻射管反射的熱速率相同，有：

根據以上傳熱方程及熱平衡條件，燃料氣的流量調節回路，對加熱爐煙氣理論氧含量[15]計算如下：

而對于加熱爐熱效率，目前有兩種方法計算得到：正平衡法和反平衡法。可根據實際爐子熱負荷和實際燃料消耗量計算，即正平衡方法，也可根據各種損失由反平衡法計算。本文根據反平衡法計算[16]：

對于本文研究的中石化某煉油廠的渣油，Q=39330 kJ/kg，L=13.225，H=202.6 kJ/kg。圓筒爐及四壁排管的加熱爐q3取3，有廢熱回收系統在公式外額外增加0.5，根據以上條件將式(6)簡化：

其中對于化學不完全燃燒熱損失q2，在一氧化碳所占煙氣百分比為1%以下時根據經驗公式（8）計算得到：

式中，CO為煙氣中的一氧化碳濃度。根據化驗室得到的數據得知常年煙氣中的一氧化碳濃度在0.01%左右，可以使用式(8)。

熱效率還與過剩空氣系數α和排煙熱損失q1有關。在工業生產中，燃料在化學平衡中所需空氣量即理論空氣量下是不可能完全燃燒的,因此需要過提供一定量的空氣即過剩空氣,以保證燃料的完全燃燒[17]，這就是過剩空氣系數。各種燃煤的工業爐窯和鍋爐排出的煙溫較高，燃料熱的相當大部分未加利用，此即排煙熱損失[18]。熱效率與過剩空氣系數和排煙熱損失的關系如圖2所示。相關參數具體推導見附錄1[19-21]。

結合式(6)～式(10)，根據中石化某煉油廠一個月的化驗數據回歸得到熱效率擬合輸出為：

圖2 熱效率與過剩空氣系數和排煙熱損失的關系Fig.2 Thermal efficiency versus excess air and flue gas heat loss

其中，O2為煙氣中氧含量，T為出口煙氣溫度。結合式(5)～式(11)選取出口介質溫度Tof，爐膛煙氣溫度Tfg，煙氣氧含量Ofg，燃料氣進料Fg為狀態變量，爐膛煙氣出口溫度Tofg，燃料氣進料設定值Fsg，單位時間內加熱爐空氣流入的體積Fa與進料介質溫度Tif為輸入變量，W為狀態噪聲，V為觀測噪聲。得到如下狀態空間方程：

其中

3 基于粒子濾波聯合估計的加熱爐熱效率實時估計

3.1 Kalman濾波

實際工業過程由于建模偏差、噪聲及擾動的影響，會造成質量指標的模型預測偏差，因此需要利用離線分析數據對模型預測進行實時反饋修正[22]。在線性離散系統中經典Kalman 濾波（KF）是一種最優估計方法[23]。考慮如下的線性離散動態系統：

其中xk是k時刻的系統狀態，uk是k時刻對系統的控制量，wk和vk分別表示過程和測量的噪聲，假設成高斯白噪聲，協方差分別為Q和R。假設：

（1）過程噪聲wk和測量噪聲vk均是高斯白噪聲序列，且二者互不相關或者二者δ相關，即

式中，μw是wk的均值矩陣，μv是vk的均值矩陣，Qk為非負定方差矩陣，是系統過程噪聲wk的方差矩陣，Rk為正定方差矩陣，是系統測量噪聲vk的方差矩陣，Sk是wk和vk的互協方差矩陣，δkj是克羅尼克δ符號。

（2）系統初始狀態X0是某種已知分布或正態分布的隨機向量，其均值、協方差分別是:

（3）系統的過程噪聲wk和測量噪聲vk都與初始狀態X0互不相關，即

對于式(14)所示系統，滿足以上條件，其時間更新方程[24]如下：

3.2 粒子濾波

粒子濾波（PF）[25]是一種基于貝葉斯濾波思想，應用重要性采樣的概念，來獲得狀態向量后驗概率密度的序貫蒙特卡洛方法。考慮式（21）描述系統：

wk-1和vk為獨立同分布噪聲（非白噪聲），根據貝葉斯定理，假設xk服從一階馬爾科夫過程，且狀態序列xk和量測序列yk相互獨立，則狀態預測方程和更新方程分別如式（22）、式（23）所示：

其中p(yk|y1:k-1)為歸一化常數。在對后驗分布進行蒙特卡羅數值模擬計算時，直接利用未知的后驗概率密度是非常困難的。因此，根據貝葉斯重要性采樣定理[26-27]，引入一個容易采樣的參考分布q(x0:k|y1:k)，則后驗概率密度函數可近似為：

歸一化重要性權值，就可由式（24）計算得到狀態的概率分布，然而，隨著序貫重要性采樣算法的進行，一部分粒子會因為權值逐漸變小而失去作用，這對粒子濾波的進行非常不利。為了克服這種狀況，加入重采樣策略[28-29]，其基本思想是剔除權值小的粒子，同時對大權值粒子根據權值大小進行復制。

3.3 聯合估計

粒子濾波可以很好地克服有色噪聲對模型帶來的影響，但在一些過程中PF為了滿足估計精度要求，就需要更多的粒子數和更長的計算時間，這樣就大大削弱了實時估計的優勢[30]。針對這種情況，本文參照狀態與未知參數聯合估計方法[31]，結合KF 和PF，提出了基于粒子濾波的狀態和校正系數聯合估計方法（JPF）來修正預測模型。考慮如下系統：

基于預測模型的狀態估計如下所示：

雖然對于式（27）可以利用KF 進行估計，但是這樣不能很好地處理非高斯噪聲和模型的大誤差。并且KF 穩態濾波也不是常數，而應該是一個變值。故可以將K描述成：

其中K0由KF 得到，ξ3為噪聲。于是，可將系統擴維成：

從這個擴維系統中可以看出，濾波增益的動態方程中由于存在噪聲，整個增益輸出不夠平滑。為此，對其采用一階慣性濾波，即：

其中λ(0 <λ< 1)為濾波常數，整理得到最終的擴維系統如式(31)所示：

證明見附錄2[32-33]。由定理可見，當粒子數和時間趨于無窮時，該方法收斂于最優濾波。

4 工業應用

本文建立了中石化某煉油廠如圖3所示的延遲焦化加熱爐裝置部分的狀態空間模型。熱效率輸出模型參數使用最小二乘法擬合得到，對所得到的模型進行模型預測驗證結果如圖4 所示，紅色點劃線部分為廠里化驗分析數據，藍色實線部分為狀態空間模型預測值。數據分析對比見表1，可以看到通過預測模型計算得到的熱效率數值與化驗分析得到的數值的誤差不大于0.3，誤差相比較分析計算得出的熱效率誤差不大于0.3%，完全可以作為加熱爐熱效率的實時輸出。將化驗分析值分別與聯合粒子濾波、粒子濾波和經典Kalman濾波的估計值進行比較，效果如圖5所示，具體結果見表2。

圖3 延遲焦化流程圖Fig.3 Process flow diagram

圖4 模型計算對比Fig.4 Model prediction comparison with analysis data

表1 模型計算對比Table 1 Model prediction results

圖5 濾波效果圖Fig.5 Filter effect diagram

表2 各種濾波比較Table 2 RMSE comparisons with different filters

通過圖5 與表2 可以看出運用聯合粒子濾波比單獨使用Kalman濾波或單獨使用粒子濾波，標準誤差分別下降61.53%和70.18%。粒子濾波聯合估計的精度要比PF 和KF 高，這說明聯合粒子濾波在整體性能上是優于粒子濾波和經典Kalman 濾波的。但無論粒子濾波還是聯合粒子濾波都存在合適粒子數問題，具體分析結果見表3。

表3 粒子濾波與聯合粒子濾波對比Table 3 Comparison of particle filtering with joint particle filtering（PF/JPF）

通過表3 可以看出隨著粒子數的增加，無論是PF 還是JPF，均方根誤差（RMSE）都在減小，相對應的運行時間均在延長。但二者處于相同的均方根誤差（RMSE）時，使用JPF 比PF 的運行時間縮短20%，很好地克服了粒子濾波的需要更多的粒子數和更長的計算時間的問題。將模型參數由上層工程師站組態進DCS 中，建立軟測量界面，加熱爐熱效率的實時輸出如圖6 所示，為實現裝置的先進控制奠定了基礎。

圖6 熱效率在線輸出曲線Fig.6 Thermal efficiency online output

5 結 論

對于延遲焦化裝置，管式加熱爐消耗了能耗的大部分，但熱效率并不能實時測量，操作人員不能根據熱效率的實時變化進行相應的操作，以達到節能降耗的目的。本文結合管式加熱爐物料，熱量平衡條件，構造了加熱爐的狀態空間方程并提出了一種基于粒子濾波的聯合方法對模型預測值進行實時估計，并對方法收斂性進行分析，工業應用結果證實所提方法的可行性和有效性。

符 號 說 明

A——爐管有效傳熱面積，m2

Ai,Ao,Am——分別為管內、管外、平均表面積，m2

AR——輻射外表面積，m2

Cf,Cfg——分別為被加熱介質焓值與煙氣焓值，kJ/kg

C0——煙氣氧含量的容量系數

Fa,Fif,Fofg——分別為單位時刻內加熱爐空氣流入量、被加熱介質流量以及煙氣出口流量，kg/h

Fg,Fsg——分別為燃料氣流量和流量設定值，kg/h

Hs——絕對黑體表面熱度，J/(K4·s)

hRC,kc——分別為煙氣傳熱系數與總的傳熱系數，W/(m2·K)

L——單位質量燃料完全燃燒所需要的理論空氣量，是其比值

Ofg——爐膛內煙氣氧含量，%

Q,H——分別為燃料低熱值與單位質量空氣在預熱溫度下的熱焓，kJ/kg

QR——輻射室的傳熱速率，J/s

q1,q2,q3——分別為排煙熱損失、化學不完全燃燒熱損失與爐墻散熱損失，%

Ri,Ro——分別為加熱爐管內與管外結構熱阻，K/W

r——混合程度系數

Tg,Tt——分別為離開輻射室煙氣溫度及輻射外表面的平均溫度，K

Tif,Tof,Tfg,Tofg——分別為被加熱介質的入爐溫度、出爐溫度、爐膛煙氣溫度以及爐出口煙氣溫度，K

U——平均表面熱強度，kJ/(m2·h)

Vf,Vfg——分別為進料爐管體積與煙道體積，m3

Vmf——在一定熱負荷下理論上所需要氧氣量，m3

α——過剩空氣系數

εg,εr,εw——分別為煙氣、輻射管表面和輻射室爐墻的黑度

ρf,ρfg——分別為被加熱介質與煙氣密度，kg/m3

τg——燃料氣回路的時間常數