適應鍋爐調峰運行的水冷壁高溫腐蝕預測模型

鄧磊，袁茂博，楊家輝，岳洋，姜家豪，車得福

（西安交通大學能源與動力工程學院，陜西 西安 710049）

2022 年我國風能、水能、太陽能等可再生能源發電占總發電量的28.7%，該比例在“雙碳”目標的驅動下仍會增長[1-2]。但由于風、光、水資源的時空分布不平衡，傳統火電站在調配能源中的地位將顯著增強。同時，大型電站煤粉鍋爐普遍采用分級配風降低NOx排放[3]，這也導致主燃區形成了強還原性氣氛，加劇了水冷壁高溫腐蝕，直接影響到鍋爐的安全運行[4]。Sun等[5]和Xiong等[6]均在實施空氣分級改造后的鍋爐中發現了腐蝕深度超過1.5mm 的水冷壁管。于英利等[7]通過調節雙切圓鍋爐同層二次風以緩解水冷壁高溫腐蝕問題。王新宇等[8]以對沖燃燒鍋爐為研究對象，模擬水冷壁近壁面的H2S 和CO 分布并提出高溫腐蝕防治措施。然而前人研究主要為特定負荷下的爐膛煙氣場模擬和腐蝕現象分析，難以直接應用于鍋爐調峰運行水冷壁高溫腐蝕程度的評估。Kung[9]通過實驗得出了鍋爐鋼年腐蝕速率與壁溫、H2S濃度以及合金中鉻元素含量的擬合式，如式(1)所示。

式中，CR 為年腐蝕深度，mm；T為金屬溫度，K；CH2S為H2S 濃度，μL/L；ωCr為合金中鉻元素的質量分數，%。

式(1)的年腐蝕速率是以實驗腐蝕速率線性外推至全年得到的，但Xu等[10]的相關實驗表明，H2S高溫腐蝕速率與腐蝕時間是非線性。同時，水冷壁的年腐蝕深度是由不同工況下的腐蝕深度疊加組成，而各個工況的運行時長也存在差異，因此時間維度在鍋爐調峰運行下水冷壁高溫腐蝕深度的計算中尤為重要。對此，Yuan 等[11]提出了包含時間維度的H2S高溫腐蝕數學模型，如式(2)所示。

式中，d-為腐蝕深度，μm；t為腐蝕時間，h；E為腐蝕活化能，kJ/mol；R為氣體常數；系數a、b和指數n取決于金屬性質。

式(2)將腐蝕深度與壁面溫度、H2S 濃度以及腐蝕時間聯系起來，根據水冷壁壁溫分布、近壁面H2S 濃度分布和運行時間即可得到水冷壁高溫腐蝕狀態。在前期工作[12-13]的基礎上，本文將水冷壁熱流分布模擬計算與水動力特性計算相耦合，進一步結合管壁溫計算，提出水冷壁壁溫分布計算方法。同時，根據燃料硫釋放以及含硫組分的相互轉化特性，建立一種適用于還原性氣氛的SOx生成模型，可計算爐膛內H2S 濃度分布。綜合爐膛數值模擬、水冷壁壁溫耦合計算以及高溫腐蝕數學模型，構建出適應鍋爐調峰運行的水冷壁高溫腐蝕預測模型并基于Matlab GUI 平臺開發對應預測軟件。選取一臺600MW 四角切圓燃煤鍋爐為研究對象并建立爐膛數值模型，對100%鍋爐最大連續蒸發量工況（BMCR）、75%熱耗率驗收工況（THA）、50%THA 以及35%BMCR 四種典型負荷開展數值模擬計算，獲得壁面熱流密度分布以及近壁面H2S 濃度分布；將數值模擬結果導入高溫腐蝕預測模型可得到單一負荷下或多個負荷下的水冷壁高溫腐蝕狀態，進而為鍋爐調峰運行水冷壁安全監測提供參考。

1 鍋爐數值計算模型

1.1 物理模型與工況設置

本文選取一臺超臨界600MW 直流四角切圓燃煤鍋爐為研究對象。由于實際鍋爐結構復雜，如圖1 所示對其幾何形狀和燃燒系統進行了合理簡化。爐膛截面深度為17.696m、寬度為18.816m、鍋爐總高為63.750m，水冷壁由傾角13.95°的螺旋管圈和垂直管圈組合構成。燃燒系統由煤粉噴嘴（PA）、 輔 助 風 噴 嘴（AUX）、 緊 湊 燃 盡 風（CCOFA）以及可分離燃盡風（SOFA）組成，輔助風噴嘴位于相鄰的煤粉噴嘴之間，并由上下兩只預置水平偏角的輔助風（CFS）噴嘴以及一只直吹風噴嘴組成。燃燒器中心線和爐壁的夾角為51°和48°，旨在爐膛內形成切圓燃燒，CFS 噴嘴和爐壁夾角為22°和19°。燃燒系統將下爐膛分為主燃區、還原區和燃盡區三個部分。爐膛上部從前到后依次布置有分隔屏、后屏、末級再熱器和末級過熱器，均簡化為無厚度平面。100%BMCR、75%THA、50%THA 和35%BMCR 四個負荷的具體參數設置如表1所示。

表1 不同負荷下鍋爐運行參數

圖1 鍋爐布置簡圖

1.2 數值模型選擇

參照圖1所示鍋爐布置簡圖建立爐膛的三維模型，應用Fluent軟件設置計算模型，對爐膛傳熱傳質過程做穩態計算。選取Realizablek-ε湍流模型和DO輻射模型分別對流場和輻射傳熱進行計算[14]；采用灰色氣體加權求和模型計算煙氣的輻射吸收系數；采用基于歐拉-拉格朗日方法的隨機軌道模型模擬煤粉顆粒的運動軌跡，粒徑大小遵循Rosin-Rammler 分 布[15]，范 圍 為1～100μm，平 均 粒 徑61μm，煤樣的元素分析和工業分析見表2。采用Multiple surface reaction 模型計算煤焦與O2/CO2/H2O三種氣體的異相反應速率[16]。爐膛傳熱傳質計算收斂后，加載UDF自定義SOx生成模型，進一步計算H2S、COS 和SO2的濃度場。SOx生成模型定義燃料硫向H2S、COS以及SO2轉化的速率由式(3)～式(5)表示；含硫組分的相互轉化為表3所示的10個總包反應[17]，最終三種氣體的生成速率（rf）與消耗速率（rr）由式(6)～式(11)表示。

表2 煤樣的工業分析與元素分析

表3 含硫組分相互轉化反應方程[17]

式中，α和γ分別為揮發分和灰分中硫的質量分數，分別取值0.55 和0.1；rvol和rc分別為揮發分的揮發速率和焦炭的燃燒速率，kg/s；Mw,S為硫的摩爾質量，0.032kg/mol；Vcell為網格單元體積，m3；pc為焦炭氧化反應速率占焦炭總反應速率的比例，計算方法參見本文作者課題組[13]前期的工作。

1.3 數值模型驗證

為選擇合適的網格密度，對比了網格數量分別為1.34×106、1.75×106、2.07×106和2.48×106的四個模型沿爐高方向截面平均煙溫，結果如圖2所示。網格數量為207 萬和248 萬的溫度曲線十分接近，繼續增加網格數量對計算結果準確度的提高有限。綜合準確性和經濟性兩方面因素，選用207 萬網格數量的爐膛模型開展研究；選取100%BMCR 負荷下實爐測量數據與模擬結果進行對比，驗證數值計算的可靠性。實測數據的獲取方式為：水冷壁吸熱量是根據省煤器出口和汽水分離器入口的工質焓差計算得到；爐膛出口氧濃度和后屏出口煙溫是根據布置在煙道中的若干測點數據取平均得到；H2S濃度值是在水冷壁近壁面采用煙氣分析儀測量得到。如表4所示，水冷壁吸熱量、爐膛出口氧濃度和后屏出口煙溫的相對偏差分別為7.50%、2.76%和4.79%。鍋爐運行過程中壁面附近的H2S濃度會在一定范圍內波動，因此數值模擬計算的H2S濃度與實際測量值之間存在一定的偏差。標高25m、29m 和35m 分別對應于A 層、C 層和F 層燃燒器高度位置，前墻中心線上H2S 濃度的相對偏差分別為7.76%、8.21% 和8.37%，可見各維度的計算偏差均在合理范圍內。

表4 100% BMCR負荷下實爐測量結果與模擬結果對比

圖2 網格無關性驗證

2 水冷壁壁溫耦合計算方法

螺旋管圈水冷壁為大容量電站鍋爐下爐膛水冷壁的主要形式，傳統的水動力回路劃分是對爐膛四面墻進行分割，沒有考慮螺旋管圈水冷壁的特殊布置形式，而工質在螺旋管圈水冷壁內盤旋上升的行程是連續的。本節提出了一種基于坐標變換的熱量再分配方法，建立了水動力回路吸熱量與工質在流動過程吸熱量的映射關系。進一步耦合煙氣側與工質側傳熱，結合管壁溫計算，形成水冷壁壁溫耦合計算方法。算法在Matlab平臺上實現，具有集成度高、靈活性強的特點。

2.1 螺旋管圈水冷壁熱量再分配方法

圖3 為螺旋管圈水冷壁展開與重構的示意圖。首先將水冷壁展開形成由前墻、右墻、后墻、左墻依次排列的平面，原三維螺旋管圈中某連續回路（ABCDEF）在展開平面上被分為ABCDE和EF的兩段。為使回路重新拼合，對AE連線以上區域整體向右平移，直至E點與（E）點重合，移動距離為四墻寬度的總和。由于螺旋管圈普遍匝數在1.5 左右，因此E點落在左墻的高上，此時還需將（E）F連線以上區域向右平移兩個四墻寬度總和。平移過程同樣是壁面熱流密度點坐標變換的過程。對變換后的平面進行網格劃分，得到壁面重構后的熱量分布。

圖3 螺旋管圈展開與重構示意圖

2.2 水動力計算方法

鍋爐水動力特性計算的核心是獲取鍋爐各個管道中工質流量、壓力和焓達到平衡時的各項參數。設第i回路壓降為該回路中各段的重位壓降（Δph）和摩擦壓降（Δpf）的累加，如式(12)所示。

單相重位壓降如式(13)所示。

式中，h為第i回路第j段垂直上升距離，m；ρ為第i回路第j段的工質密度，kg/m3。

汽水混合物重位壓降如式(14)所示。

式中，ρl為第i回路第j段的飽和水密度，kg/m3；ρs為第i回路第j段的飽和蒸汽密度，kg/m3。

光管單相摩擦壓降如式(15)所示[18]。

式中，L為第i回路第j段長度，m；D為管內徑，m；G為單位面積質量流量，kg/(m2·s)；ν為第i回路第j段的比熱容，m3/kg；f為湍流光管單相摩擦系數。

汽水混合物的摩擦壓降如式(16)所示[19]。

式中，x為第i回路第j段含干度；νl為飽和水比熱容，m3/kg；νs為飽和蒸汽比熱容，m3/kg；ψ為摩擦阻力修正系數。

根據并聯管組壓降相等原則以及管路質量守恒，設并聯管壓降和總質量流量分別為ΔP和M，得式(17)壓降方程和式(18)質量守恒方程。

壓降方程和質量守恒方程可組成包含i+1 個方程的非線性方程組，采用離散牛頓法對方程組進行求解。由于摩擦壓降和重位壓降的計算方程中含有工質物性相關量，物性參數的確定也依賴于回路流量，因此方程組的解若不收斂，需用各回路流量更新方程組中的物性參數，直至計算收斂，獲得最終的各管路流量和回路總壓降（計算中水/蒸汽的物性參數從IAPWS-IF97 數據庫[20]中獲取）。

2.3 管壁溫計算方法

亞臨界圓管的單相傳熱系數采用Dittus-Boelter公式進行計算，如式(19)所示。

式中，λ為管的熱導率，W/(m·K)；Re和Pr分別為雷諾數和普朗特數，實驗驗證范圍為Re=104～1.2×105，Pr=0.7～120。

超臨界光管的對流傳熱系數如式(20)所示[21]。

式中，αl為飽和水的對流傳熱系數，可采用式(19)計算得到，W/(m2·K)；μl為第i回路第j段的飽和水動力黏度，N/(m·s)；μs為第i回路第j段的飽和蒸汽動力黏度，N/(m·s)；p為工質壓力，Pa；pcr為臨界壓力，22.115MPa；G為質量流量，kg/(m2·s)。

各回路各段的管內壁溫度如式(22)[22]所示。

式中，Tf為流體溫度，K；τ為管外徑與內徑的比值，De/D；Jn為向火側內壁均流系數；q為壁面熱流密度。

管外壁溫度如式(23)[23]所示。

式中，δ為管壁厚度，m；λp為管壁熱導率，W/(m·K)。

2.4 計算方法驗證

研究選取的超臨界600MW 四角切圓鍋爐在螺旋管圈水冷壁出口高度沿周向均勻布置了70 個溫度測點，為驗證壁溫耦合計算方法準確性，將螺旋管圈水冷壁劃分為70 個回路。圖4 對比了100%BMCR 負荷下螺旋管圈水冷壁出口的實測溫度和計算溫度。結果顯示，100%BMCR 負荷下的計算溫度與測點數據的最大相對偏差（即最大絕對偏差值與對應測點溫度的比值）為2.7%。壁溫耦合計算方法具有較高的計算精度。

圖4 100%BMCR負荷下螺旋管圈水冷壁出口溫度分布對比

3 水冷壁高溫腐蝕預測模型

圖5為水冷壁高溫腐蝕預測模型的示意圖。將鍋爐各工況的運行參數作為邊界條件輸入爐膛數值模擬，結合SOx生成模型和壁溫耦合計算方法，輸出對應工況下的水冷壁壁溫分布和近壁面H2S濃度分布。進一步將各工況的持續時間、壁溫分布和H2S濃度分布代入高溫腐蝕模型，累加即得到一段時間內的水冷壁高溫腐蝕狀態。本文選取的鍋爐鋼H2S 高溫腐蝕數學模型如式(24)所示，各項系數及腐蝕活化能等參數由前期工作中[11]腐蝕實驗數據的回歸分析得到。

圖5 水冷壁高溫腐蝕預測模型示意圖

式中，t為腐蝕時間，h；T為金屬溫度，K；CH2S為H2S濃度，μL/L。

根據圖5的計算思路，基于Matlab GUI 平臺開發了如圖6 所示的水冷壁高溫腐蝕預測軟件，實現預測算法與用戶的交互。軟件界面主要包含三個模塊：計算區域設定模塊、水動力參數模塊以及高溫腐蝕預測模塊。使用軟件時首先確定計算區域與網格數量，如輸入爐膛尺寸、高度區間、高度方向和寬度方向網格數量等；其次提供對應工況下的水冷壁內工質壓力、工質的單位面積質量流量以及進入計算區域的工質焓；最后導入該工況下的數值模擬計算結果，設置工況運行時間。數值模擬數據中存儲有壁面網格中心對應的三維坐標、熱流密度和H2S 濃度，點擊“計算腐蝕深度”按鈕激活后臺程序。后臺程序會抓取計算區域設定模塊內的爐膛邊界和網格密度數據，對導入的數值模擬數據進行處理。同時抓取界面上的水動力計算參數，結合數值模擬結果，為不同工作壓力下的水冷壁壁溫分布耦合計算匹配對應函數，具體的函數及功能如表5 所示。根據壁溫和H2S濃度分布獲得特定工況下的水冷壁高溫腐蝕狀態。通過累加不同工況下的腐蝕深度，實現鍋爐調峰運行下水冷壁高溫腐蝕狀態的預測。積累大量數據后，可利用機器學習算法建立鍋爐運行參數（鍋爐負荷、一次風率、燃燒器豎直擺角等）與水冷壁高溫腐蝕程度的映射關系，從而實現水冷壁高溫腐蝕狀態的實時更新。

表5 壁溫計算函數匯總

圖6 水冷壁高溫腐蝕預測GUI界面

4 模型計算結果

4.1 水冷壁高溫腐蝕與壁溫分布和H2S 濃度分布的關系

相關研究發現切圓鍋爐下爐膛前、后墻高溫腐蝕較為嚴重[24]，故選取100%BMCR 負荷主燃區水冷壁前墻為研究對象。圖7 為該墻壁溫、H2S 濃度以及腐蝕深度的分布云圖。如圖7(a)所示，壁溫處于650～730K 之間，燃燒器區域壁溫沿高度方向快速上升，且右側壁溫明顯低于左側；進入SOFA 燃盡區后壁溫增速大幅降低甚至出現負增長。這是由于順時針切圓燃燒方式下，射入爐膛的低溫二次風對前墻右側產生了冷卻效果，主燃區燃燒溫度高、輻射換熱量大，壁溫沿高度方向增速高；進入SOFA 區域后，輻射換熱量降低，壁溫增速緩慢。可見，熱力-水動力耦合計算模型能夠充分反映水冷壁的壁溫分布，為高溫腐蝕的準確預測提供基礎。

圖7 100%BMCR負荷水冷壁前墻壁溫分布、H2S濃度分布以及年腐蝕深度分布

圖7(b)顯示高濃度H2S 主要位于燃燒器區域下部，隨著爐膛高度的增加，H2S濃度波動下降，并在燃燒器層和SOFA 層之間區域出現局部的極大值。進入SOFA 區域后H2S 濃度迅速下降。通過數值模擬得到的在爐膛內H2S 濃度分布規律與孟繁兵等[25]的實爐測試規律吻合度高。圖7(c)給出了前墻水冷壁高溫腐蝕程度，云圖顯示底層燃燒器區域和SOFA 區域腐蝕程度較低，這兩個區域分別對應于高H2S 濃度和高壁溫。最大年腐蝕深度出現在燃燒器層和SOFA 層之間區域，為276μm。該區域特征是壁溫和H2S 濃度都相對較高，當高壁溫和高腐蝕氣體濃度重疊時，水冷壁易發生嚴重的高溫腐蝕。

4.2 不同負荷下水冷壁高溫腐蝕特征

表1 中的四個負荷下水冷壁前墻的壁溫分布、H2S 濃度分布以及高溫腐蝕狀態分別如圖8～圖10所示，年運行時間設為8760h。100%BMCR負荷下燃燒器與SOFA之間的區域腐蝕最為嚴重，最大年腐蝕深度為276μm。此外，燃燒器區域的高溫腐蝕同樣較為嚴重，75%THA負荷下的高溫腐蝕較為嚴重區域與100%BMCR 負荷基本一致，但整體的腐蝕程度較低；50%THA 和35%BMCR 負荷下，高溫腐蝕深度在F 層燃燒器高度位置迅速達到最大值，其中35% BMCR 負荷下最大年腐蝕深度達到256μm。各個負荷下前墻左側的腐蝕程度均高于左側，這與熊小鶴等[24]觀測的現象一致。

圖8 不同負荷下水冷壁前墻壁溫分布

圖9 不同負荷下水冷壁前墻H2S濃度分布

圖10 不同負荷下水冷壁前墻高溫腐蝕狀態

結合圖8 和圖9 所示的壁溫分布和H2S 濃度，高負荷下鍋爐水冷壁的平均壁溫較高，且燃燒器層和SOFA層之間區域的水冷壁缺少水平偏置二次風的保護，該區域H2S 濃度較高，使得燃燒器層與SOFA 層之間的區域成為高溫腐蝕最為嚴重的區域。項岱軍等[26]同樣發現SOFA 風下方區域高溫腐蝕比主燃區更嚴重。此外，為抑制氮氧化物的生成，高負荷下爐膛空氣分級度高，一、二次風強度減弱，切圓燃燒偏斜程度大，導致燃燒器區域H2S濃度較高，高溫腐蝕也較為嚴重。50%THA 和35%BMCR 負荷下工質處于亞臨界狀態，水冷壁內工質發生相變的高度位置容易出現壁溫激升。如圖8 所示，壁溫在燃燒器區域內出現峰值，35%BMCR 負荷下壁溫峰值位置更高，與壁面高濃度H2S 區域重合，這使得35%BMCR 負荷下燃燒器區域的上部腐蝕深度超過50%THA 工況。通過調整燃燒器擺角、主燃區空氣過量系數等運行參數，錯開高壁溫區與高濃度H2S區，即可有效減緩低負荷下的水冷壁高溫腐蝕程度。

4.3 調峰運行下水冷壁高溫腐蝕狀態

為模擬鍋爐調峰運行，設計了鍋爐在100%BMCR、75%THA、50%THA 和35%BMCR 四個負荷下各運行2190h 共計8760h 的情境。各工況的運行參數依照表1進行設置，使用水冷壁高溫腐蝕預測GUI界面進行計算。該情境下的水冷壁高溫腐蝕狀態如圖11 所示。其中，后墻的高溫腐蝕最為顯著，最大腐蝕深度為364μm。腐蝕嚴重區域位于F 層燃燒器高度位置以及燃燒器層和SOFA 層中間的位置。前墻、左墻和右墻的特征也與之類似，但腐蝕程度較輕。

圖11 調峰運行下水冷壁高溫腐蝕狀態

結合圖10 所示的不同負荷下水冷壁前墻高溫腐蝕狀態，多工況運行下水冷壁高溫腐蝕狀態表現為各個負荷運行時的高溫腐蝕狀態的時空疊加。通過本文的水冷壁高溫腐蝕預測模型，能夠掌握鍋爐調峰運行下水冷壁高溫腐蝕狀態的空間分布和時間分布，為構建鍋爐安全監測系統提供有利參考。同時，鍋爐的運行參數能夠影響水冷壁的壁溫分布以及近壁面的H2S 濃度分布，進而影響水冷壁的高溫腐蝕程度。運用高溫腐蝕預測模型可對鍋爐運行參數進行優化，為減緩水冷壁高溫腐蝕提供支持。

5 結論

本文綜合爐膛數值模擬、水冷壁壁溫耦合計算以及包含時間維度的管壁高溫腐蝕模型，提出一種適應鍋爐調峰運行的水冷壁高溫腐蝕預測模型。為預測驗證模型可靠性，選取一臺600MW 四角切圓燃煤鍋爐為研究對象，計算了100%BMCR、75%THA、50%THA 以及35%BMCR 四種典型負荷下壁溫分布和近壁面H2S濃度分布，進而預測出單一工況以及多工況疊加情境水冷壁的高溫腐蝕狀態，主要結論如下。

（1）根據螺旋管圈水冷壁的特征采用一種基于坐標變換的熱量再分配方法，進一步耦合煙氣側與工質側傳熱，結合管壁溫計算，形成水冷壁壁溫分布耦合計算方法。同時，建立了用于確定還原性氣氛下的燃料硫釋放以及含硫組分的相互轉化過程的SOx生成模型。通過壁溫耦合計算和SOx生成模型分別得到水冷壁的壁溫分布和H2S濃度分布，為高溫腐蝕的準確預測提供基礎。

（2）不同負荷下水冷壁高溫腐蝕特征存在區別，壁面腐蝕程度整體上隨負荷降低而降低。100%BMCR 與75%THA 負荷下水冷壁前墻燃燒器層與SOFA層之間的區域腐蝕最為嚴重，最大年腐蝕 深 度 分 別 為276μm 和233μm；50%THA 與35%BMCR 負荷下高溫腐蝕深度在燃燒器區域的上部迅速增加至最大值，分別為224μm和256μm。

（3）多工況運行水冷壁高溫腐蝕狀態表現為各工況腐蝕狀態的時空疊加。在100%BMCR、75%THA、50%THA 以及35%BMCR 四個負荷下各運行2190h情境下，后墻高溫腐蝕程度最嚴重，頂層燃燒器高度位置腐蝕深度達到364μm，其他三面墻的腐蝕特征與之類似。

（4）采用Matlab GUI 平臺開發了水冷壁高溫腐蝕預測軟件，高度集成了爐膛數值模擬、水冷壁壁溫耦合計算方法以及高溫腐蝕數學模型，實現了通過鍋爐運行參數和運行時間對多工況下水冷壁高溫腐蝕狀態的時空分布進行預測，為靈活調峰需求下鍋爐安全監控系統的構建和開發提供了有力支持。