陸 海， 羅恩博， 孫鑫南， 周 懿， 蘇 適， 鐘 崴、

(1.云南電網(wǎng)有限責(zé)任公司電力科學(xué)研究院，云南昆明650217；2.浙江大學(xué)工程師學(xué)院，浙江杭州310011；3.浙江大學(xué)能源工程學(xué)院熱工與動(dòng)力系統(tǒng)研究所，浙江杭州310027)

1 概述

在蒸汽輸送過程中，蒸汽的狀態(tài)參數(shù)將發(fā)生變化以及可能發(fā)生相變[1]。王威等人[2]綜合熱平衡原理和蒸汽管網(wǎng)的傳熱特點(diǎn)，建立了蒸汽管網(wǎng)散熱損失數(shù)學(xué)模型、沿途溫降計(jì)算模型和凝結(jié)水量計(jì)算模型。張?jiān)鰟俒3]通過對蒸汽熱力性質(zhì)計(jì)算方法、蒸汽管網(wǎng)熱損失分析、水力和熱力計(jì)算理論的系統(tǒng)性研究，建立了蒸汽管網(wǎng)水力-熱力耦合計(jì)算模型。文獻(xiàn)[4]、[5]均利用Excel的函數(shù)編輯和計(jì)算功能，采用IAPWS-IF97相關(guān)計(jì)算公式，耦合水力和熱力計(jì)算，量化蒸汽管道能耗損失，為長距離輸送蒸汽設(shè)計(jì)提供依據(jù)。

此外，針對蒸汽管道流動(dòng)的模擬計(jì)算而涌現(xiàn)出的一批商業(yè)軟件，如Flowmaster、PROSS、SynerGEE gas、TERMIS等大都側(cè)重穩(wěn)態(tài)計(jì)算，在蒸汽管道流動(dòng)的瞬變特性方面研究較少。孫雅慧[6]基于一維仿真軟件Flowmaster，從蒸汽管網(wǎng)的整體和變工況運(yùn)行條件出發(fā)，對蒸汽管網(wǎng)進(jìn)行水力仿真，得到蒸汽管網(wǎng)的水力特性，依此改善蒸汽管網(wǎng)的設(shè)計(jì)布局和流量分配。張紅玉[7]基于SynerGEE gas仿真模擬軟件建立了天津空港加工區(qū)蒸汽供熱系統(tǒng)穩(wěn)態(tài)模型，并利用該模型分析評估了不同工況下蒸汽供熱系統(tǒng)的運(yùn)行情況，從而確定了管網(wǎng)質(zhì)量損失與散熱損失的根源，并給出了相應(yīng)的運(yùn)行管理建議。

Modelica語言由Elmqvist等人在1997年的第一屆世界系統(tǒng)仿真會(huì)議上首次提出[8]。Modelica實(shí)現(xiàn)了不同應(yīng)用領(lǐng)域的復(fù)雜物理系統(tǒng)統(tǒng)一建模，可有效促進(jìn)建模知識的可重用性和可拓展性。文獻(xiàn)[9]-[11]從概念、特征、用法、實(shí)例等方面講述了Modelica語言是如何描述物理系統(tǒng)的行為的。Modelica語言有以下特點(diǎn)：面向?qū)ο蠼?、陳述式物理建模、多領(lǐng)域統(tǒng)一建模，建模過程中不過分關(guān)注方程的求解過程，而重點(diǎn)關(guān)注如何去精準(zhǔn)地描述一個(gè)實(shí)際物理系統(tǒng)。Rickard[12]為模擬區(qū)域供熱網(wǎng)絡(luò)，改進(jìn)了Modelica語言標(biāo)準(zhǔn)庫中的源端、末端以及管道模型，增加了模擬散熱損失等功能，基于Dymola軟件(Dynamic Modeling language，Dymola是Modelica語言的商用化軟件)先在組件級別上進(jìn)行測試，然后用于模擬更大的區(qū)域供熱網(wǎng)絡(luò)，并與現(xiàn)有區(qū)域供熱網(wǎng)絡(luò)數(shù)據(jù)進(jìn)行比較。仿真結(jié)果表明，模擬與實(shí)測數(shù)據(jù)差異很小。

本文綜合考慮蒸汽在管道的流動(dòng)過程中的熱力、水力耦合因素以及動(dòng)態(tài)特性，建立基于Modelica語言的單根蒸汽管道傳輸動(dòng)態(tài)仿真模型。以上海地區(qū)某化工園區(qū)某段架空蒸汽管道的熱源蒸汽運(yùn)行參數(shù)(壓力、溫度、質(zhì)量流量)作為模型仿真輸入，仿真計(jì)算用戶蒸汽運(yùn)行參數(shù)(壓力、溫度)，將仿真值與實(shí)測值進(jìn)行對比，驗(yàn)證仿真模型的計(jì)算精度。利用仿真模型對用戶端的熱力響應(yīng)、水力響應(yīng)進(jìn)行研究。

2 建模機(jī)理

2.1 數(shù)學(xué)模型

考慮一根長為L，流道截面積為A的管道，認(rèn)為管內(nèi)介質(zhì)充分均勻混合，且在同一截面內(nèi)有均勻的流速，介質(zhì)做一維流動(dòng)(流向沿Ox方向)，可由以下方程、計(jì)算式描述。

質(zhì)量方程[3]：

式中ρ——蒸汽密度，kg/m3

A——流道截面積，m2

t——時(shí)間，s

w——蒸汽流速，m/s

x——Ox軸坐標(biāo)，m

能量方程[3]：

K(T-To)πdodx=0

式中h——蒸汽比焓，J/kg

p——蒸汽壓力(絕對壓力)，Pa

g——重力加速度，m/s2，本文取9.8 m/s2

z——管道中心線距地面高度，m

u——蒸汽比熱力學(xué)能，J/kg

K——管道總傳熱系數(shù)，W/(m2·K)

T——蒸汽溫度，K

To——外界環(huán)境溫度，K

do——管道外直徑，m

動(dòng)量方程[3]：

式中λ——管道摩擦阻力系數(shù)

din——管道內(nèi)直徑，m

管道沿程阻力Δp的計(jì)算式為[13]：

式中 Δp——管道沿程阻力，Pa

Re——蒸汽的雷諾數(shù)

L——管道長度，m

μ——蒸汽動(dòng)力黏度，Pa·s

qm——蒸汽質(zhì)量流量，kg/s

單位長度蒸汽管道傳熱方程式見文獻(xiàn)[14]。

2.2 模型建立

采用Modelica語言，在OpenModelica平臺上對蒸汽管道進(jìn)行建模，組成具有分布式的質(zhì)量、能量和動(dòng)量平衡的直管Modelica動(dòng)態(tài)仿真模型。Modelica動(dòng)態(tài)仿真模型由多個(gè)流動(dòng)段組成，每個(gè)流動(dòng)段包括3個(gè)模塊：管道接口、流動(dòng)模型、傳熱模型。動(dòng)態(tài)仿真模型使用有限體積法和動(dòng)量平衡的交錯(cuò)網(wǎng)格方案處理偏微分方程。

管道接口模型由2個(gè)接口連接器(即1個(gè)流動(dòng)段的進(jìn)出口)組成，用于傳遞流體介質(zhì)的流量、壓力與焓等參數(shù)，對節(jié)點(diǎn)連接處使用基爾霍夫定律進(jìn)行計(jì)算。當(dāng)管網(wǎng)結(jié)構(gòu)復(fù)雜呈現(xiàn)樹狀或者環(huán)狀時(shí)，可將若干流動(dòng)段進(jìn)行串聯(lián)，流體介質(zhì)的流通通過相鄰接口模型實(shí)現(xiàn)。

流動(dòng)模型用于描述管道中的壁面摩擦和重力引起的壓力損失，以及動(dòng)量守恒形式，使用默認(rèn)的上游離散化方案(Upstream Discretization Scheme)提供穩(wěn)態(tài)或動(dòng)態(tài)動(dòng)量平衡。上游離散化即兩流動(dòng)段接口處的跨越段邊界的流體介質(zhì)比焓以來流流體介質(zhì)比焓為準(zhǔn)，忽略下游強(qiáng)度量的影響。對于給定的流體介質(zhì)，在流動(dòng)模型中指定流體介質(zhì)的熱力學(xué)狀態(tài)，指定管道幾何參數(shù)，定義管道的結(jié)構(gòu)參數(shù)，包括長度、流道截面積、粗糙度等。

傳熱模型基于傳熱方程計(jì)算管道內(nèi)流體介質(zhì)的流量、流動(dòng)段之間的熱量傳遞、管道散熱量。流動(dòng)段散熱量由給定流體介質(zhì)流動(dòng)段的熱力學(xué)狀態(tài)函數(shù)計(jì)算，包括密度、雷諾數(shù)、比焓、動(dòng)力黏度、努塞爾數(shù)及管道表面積等。

流體介質(zhì)的5個(gè)熱力學(xué)變量分別為壓力、溫度、密度、比熱力學(xué)能、比焓，任意3個(gè)變量都是其余2個(gè)變量的函數(shù)。IAPWS-IF97方程提供了代數(shù)方程式，用于求解未知量。

3 研究對象

研究對象為上海某化工園區(qū)蒸汽管網(wǎng)，園區(qū)以煉化一體化項(xiàng)目為龍頭，發(fā)展以烯烴和芳烴為原料的中下游石油化工裝置以及精細(xì)化工深加工系列，形成乙烯、丙烯、碳四、芳烴為原料的產(chǎn)品鏈[14]。目前，園區(qū)內(nèi)有兩個(gè)集中熱源，有30多家熱用戶單位，蒸汽管道采用架空敷設(shè)，管網(wǎng)總長超過60 km。

為方便管網(wǎng)動(dòng)態(tài)特性研究，研究對象的運(yùn)行調(diào)節(jié)相對系統(tǒng)其余部分應(yīng)當(dāng)保持較高獨(dú)立性，因此本次研究選取管網(wǎng)內(nèi)一段單獨(dú)長直管段進(jìn)行研究。管段包含1個(gè)熱源、1座應(yīng)急鍋爐房以及1個(gè)用戶。研究期間應(yīng)急鍋爐不工作，且分支管被關(guān)斷，管段可視為無支路單熱源單用戶管段。管段長度為7 km，工作管外直徑為530 mm，壁厚為10 mm，鋼材為10號鋼。內(nèi)層保溫材料為硅酸鋁，保溫層厚度為80 mm。外層保溫材料為礦棉，保溫層厚度為60 mm。硅酸鋁、礦棉的熱導(dǎo)率根據(jù)文獻(xiàn)[15]給出的方法以及相關(guān)實(shí)測數(shù)據(jù)計(jì)算確定。在熱源、用戶側(cè)分別安裝溫度、壓力測點(diǎn)，管段中部安裝質(zhì)量流量測點(diǎn)。

4 仿真精度驗(yàn)證

① 輸入條件與目的

在Modelica動(dòng)態(tài)仿真模型中，將熱源蒸汽壓力、溫度及質(zhì)量流量作為輸入條件，模擬計(jì)算用戶蒸汽壓力、溫度。

② 驗(yàn)證結(jié)果

研究階段，對熱源、用戶的蒸汽壓力、溫度及質(zhì)量流量進(jìn)行了24 h實(shí)測。實(shí)測在用戶穩(wěn)定生產(chǎn)期間進(jìn)行，起始時(shí)間為0:00，截止時(shí)間為24:00，采集數(shù)據(jù)間隔為1 800 s。熱源蒸汽始終保持以下參數(shù)：壓力為4.51 MPa，溫度為323.13 ℃，質(zhì)量流量為25.60 kg/s。為避免測量初期的不穩(wěn)定性，我們選取2×104～8.64×104s的仿真值、實(shí)測值進(jìn)行比較，分別見圖1、2。由圖1、2可知，仿真值與實(shí)測值比較吻合。用戶蒸汽溫度仿真值與實(shí)測值的最大相對誤差為3.1%，用戶蒸汽壓力仿真值與實(shí)測值的最大相對誤差為-0.6%。因此，Modelica動(dòng)態(tài)仿真模型的仿真精度可滿足工程要求。

圖1 用戶蒸汽溫度仿真值與實(shí)測值隨時(shí)間的變化

圖2 用戶蒸汽壓力仿真值與實(shí)測值隨時(shí)間的變化

5 仿真結(jié)果與分析

① 熱力響應(yīng)

將熱源蒸汽壓力設(shè)定為4.5 MPa，質(zhì)量流量設(shè)定為30 kg/s。在2×104s時(shí)，將熱源蒸汽溫度在50 s內(nèi)由320 ℃升至345 ℃。在5×104s時(shí)，將熱源蒸汽溫度由345 ℃降至320 ℃。基于以上設(shè)定，仿真用戶蒸汽溫度的響應(yīng)過程，仿真結(jié)果見圖3。由圖3可知，用戶蒸汽溫度的響應(yīng)明顯滯后于熱源蒸汽溫度變化，響應(yīng)完成時(shí)間接近9 000 s。由以上分析可知，蒸汽熱網(wǎng)具有強(qiáng)烈的熱慣性和滯后性。

圖3 用戶蒸汽溫度隨熱源蒸汽溫度的變化

② 水力響應(yīng)

將熱源蒸汽溫度設(shè)定為320 ℃，質(zhì)量流量設(shè)定為30 kg/s。在2×104s時(shí)，將熱源蒸汽壓力在50 s內(nèi)由4.5 MPa降至4.0 MPa，仿真用戶蒸汽壓力的響應(yīng)過程，仿真結(jié)果見圖4。由圖4可知，用戶蒸汽壓力對熱源壓力的變化響應(yīng)迅速，響應(yīng)完成時(shí)間約200 s。

圖4 用戶蒸汽壓力隨熱源蒸汽壓力的變化

將熱源蒸汽溫度設(shè)定為320 ℃，蒸汽壓力設(shè)定為4.5 MPa。在2×104s時(shí)，將熱源蒸汽質(zhì)量流量在100 s內(nèi)由30 kg/s增至35 kg/s，仿真用戶蒸汽質(zhì)量流量的響應(yīng)過程，仿真結(jié)果見圖5。由圖5可知，用戶蒸汽質(zhì)量流量對熱源質(zhì)量流量的變化響應(yīng)迅速，響應(yīng)完成時(shí)間約180 s。

圖5 用戶蒸汽質(zhì)量流量隨熱源質(zhì)量流量的變化

6 結(jié)論

計(jì)算精度驗(yàn)證結(jié)果表明，仿真值與實(shí)測值比較吻合，用戶蒸汽溫度仿真值與實(shí)測值的最大相對誤差為3.1%，用戶蒸汽壓力仿真值與實(shí)測值的最大相對誤差為-0.6%，仿真模型的計(jì)算精度可滿足工程要求。

對于熱力響應(yīng)，用戶蒸汽溫度的響應(yīng)明顯滯后于熱源蒸汽溫度變化，蒸汽熱網(wǎng)具有比較強(qiáng)的熱慣性。對于水力響應(yīng)，用戶蒸汽壓力對熱源壓力的變化響應(yīng)迅速，用戶蒸汽質(zhì)量流量對熱源質(zhì)量流量的變化響應(yīng)迅速。