蒸汽管網的動力學建模與迭代方法1)

王安然 曾 波 張 勍 楊榮超 張 凱，2)

趙 瑤* 張 明** 趙曉東* 宿 彬* 李祎萍* 張鵬飛?

*(中國計量大學計量測試工程學院，杭州 310018)

?(鄭州煙草研究院中國煙草標準化研究中心，鄭州 450001)

**(國家電投集團河南電力有限公司，鄭州 450018)

目前，熱蒸汽作為熱電聯產生產中重要的傳輸介質，在蒸汽輸送過程中狀態參數變化大，伴隨有相變過程，具體的網絡計算方法涉及到傳熱學、流體力學方面的問題。目前對管網的設計和運行調節計算較為簡單[1-3]，水力計算主要靠人工查取水力計算圖進行近似計算，熱力計算甚至有時假設忽略，未考慮蒸汽的可壓縮性和物性改變，這給蒸汽管網的計算帶來了較大誤差。

對蒸汽管道中的流量進行計量是一個復雜的課題[4-5]，流體管網計算涉及到幾個基本方程：質量守恒、動量守恒和能量守恒[6]，運用這三個方程可以構建基礎流體管網模型[7-8]。目前該領域的研究較少，且集中在地熱蒸汽管道的應用方面[9-10]。計算蒸汽管網首先需要建立蒸汽管網模型和模擬蒸汽管網結構，管網主體是節點與管道[11-13]。

現在，存在軟件程序如PIPEPHASE 和SIMSNET 等，可以對流體管網進行數值模擬[14]，在管道網絡構建及改造過程中能發揮重大作用[15-16]，但多用于石油化工等領域。在處理蒸汽管網過程中，會存在很多問題，例如存在排水損失的問題，Wang 等[17]研究了適用于排水損失的熱力模型，考慮了管網內冷凝水的影響，證明了疏水損失對蒸汽溫度和密度沿管道的分布有顯著影響。王雅峰等[18]研究了二回路的兩相流問題，將混合相認為是兩個單相，提出了兩相流的網絡建模方法，計算結果與Jtopmeret 軟件進行了對比，證明了算法的有效性。Cazarez-Candia等[19]研究了兩相低壓多組分的數學模型，該模型將氣相和液相(油和水) 視為偽單相(混合物)，建立的均相模型可以很容易地應用于計算氣-液兩相流中的壓降、混合物溫度和混合物速度，模型預測結果與現場數據吻合。

算法研究方面，Desideri 等[20]描述了一種用于通過輸入和輸出條件(由壓力流量相關性定義)來解決管網問題的算法，該方法特別適合求解邊界條件由壓力流量二次曲線定義的管網。Huang 等[21-22]開發了用于模擬和分析地熱管網系統性能的計算機模型，該模型中使用網絡節點和回路處的質量流量和壓力平衡建立非線性方程組，具有較好收斂性。高魯鋒等[23]研究了蒸汽管網的水力熱力耦合計算方法，但未對用戶流量分配的準確性做出驗證。Karney等[24]研究從能量守恒角度闡述管網瞬態行為，通過實際案例表明能量方程對某些部件設備水錘計算是有效的，使用能量方法是現有管網計算中一個有效補充。Wang 等[25]提出了一種環狀網絡的分析方法，可以處理環狀蒸汽網絡的非等溫和可壓縮流動特性等問題。

目前，多數對蒸汽管網的算法研究使用假設管網進行驗證，本文結合熱工流體網絡理論，將蒸汽管網抽象為流體網絡模型，與實際的管網運行參數進行對比。該模型的精度取決于模型的機理性、理論假設與特性參數的準確性，模型將管系內流體的交匯與分流問題簡化為對節點壓力和支路流量的瞬態求解問題，采用流體網絡基本方程– 質量守恒方程、動量守恒方程以及能量守恒方程，將壓力與焓值相互耦合并進行非穩態求解，能夠適應各種工況和參數。

1 數學模型

在這項研究中，輸送管道內的供熱蒸汽一般為過熱蒸汽，可假設輸送過程中沒有相變，所以可以認為是單相可壓縮流體。本研究中為構建模型所做的假設可以表述如下：

(1) 蒸汽管網流體網絡內沒有動力源。

(2)節點內流體狀態均勻(內部壓力處處相等)。

(3)流動阻力僅考慮局部節流阻力與管道摩擦沿程阻力，并設流動阻力系數為常數。

(4)同一支路管道截面積不變，工質參數使用相連節點參數加權平均來表示。

構建的管網模型主要有以下兩種組件：

(1) 節點，包括管網中管道的轉折，網絡中重要組件(廠家，用戶等)。

(2) 分支，兩個節點之間的聯通組件。

特別的，管網中不同的對象在網絡中其含義不同，如圖1 所示，節點可以分為網絡流量入口的源節點、連接網絡內部分支的中間節點和用戶節點，用戶節點一般為廠家。分支元素可用作模擬管道、彎頭、閥門等組件。

圖1 蒸汽供熱管網的有向拓撲圖

1.1 單相可壓縮流體的網絡模型

對穩態蒸汽管網，流動過程中的連續性方程(質量守恒) 可以表示為

式中，Ci=Vi?ρi/(?pi) 是節點i處工質的可壓縮能力(kg/MPa)；pi是節點i處工質的壓力；ρi是節點i處工質的密度；t是時間；N是節點總數；Gij是節點i，j之間的質量流量；Vi是節點i的容積；Dij是節點i，j之間的連接方式(i= 1，2，···，N；j=1，2，···，N)，具體數值為

連續性方程表示了一個節點處的流量守恒的關系，該節點的流動狀態受到和它相連的所有節點影響。圖2 顯示了節點i與其相連的所有節點體積元素示意圖，需要注意的是在邏輯上一個節點能與N個節點相連，但實際模型中一個節點至多與三個節點相連，所以在這個情況下該節點就是一個三通類型，其中P和T分別為對應節點處壓力和溫度。

圖2 管道的第i 個節點的拓撲結構示意圖

流體模型中動量守恒方程可以寫為

式中，Mi=Viρi，是節點i內流體質量；hi是節點i的焓；Qi=γ[Tenv-0.5(Tin+Tout)]是換熱量，Tenv是環境溫度，Tin和Tout分別是流進和流出的氣體溫度，γ是管段上的散熱系數；Vidpi/dt是推進功對流體溫度的影響，通常比傳熱的影響要小，所以可以忽略此項作用，改為

式中，pa是環境大氣的壓力(MPa)；Tci是孤立節點i泄露到環境大氣的時間常數(s)。

1.2 方程組的進一步簡化

方程(2)，(4)，(6)是非線性微分、代數等式。對于動態仿真，重要的要求是要保證實時性，并且要求計算流程的魯棒性高。

式(2)表示節點i壓力與流量關系，流量常有遲滯。因此結合式(2) 和式(7) 可以寫作

為保證動態仿真算法的實時性和穩定性，采用穩定的隱式歐拉積分算法，由式(7) 和式(9) 可得

式中，r是離散時間的變量，時間t=t0+rΔt，r=0，1，2，···。

1.3 線性方程的迭代解法

可得到壓力矩陣，首先求解壓力矩陣，得到壓力，再根據壓力計算出流量，然后將流量代入節點能量方程，即可解出節點的焓值。

總體而言，國內一些機構與學者致力于引介漢學研究成果，為我們換個鏡子認識自身做出了卓越貢獻。但漢學或者說新漢學中的政治學研究，對于國內學界而言，仍然是新生事物，零星的介紹成果不少，但系統的梳理與完整的解讀，仍然非常不夠。我們將時間切片設置為1年，時間跨度選擇2008—2018年，Node types選擇為“author”，最終得到該領域高產作者結果。具體地說，在國內政治學界，專注于評介海外新漢學中的中國政治研究成果的作者，根據發文量由高到低，前五位的作者依次是路克利、梁怡、韓強、朱政惠、劉杉。

同樣由式(10) 得到節點的能量守恒方程

其中

上文的公式是計算中的具體矩陣函數，作為聯立方程的基礎，對每個節點建立的方程有兩個，將密度表示為壓力和焓值的函數ρ=ρ(P，h)，對時間求導后聯立的質量方程和能量方程為

對N個節點，能得到一個具有2N個方程的一階微分方程組，共有2N個待求變量，因此，方程組是封閉的，其計算流程如圖3 所示。

2 算例驗證

為說明耦合解法的功能，我們設計了一個假設的管網模型，其中蒸汽從蒸汽入口通入，作為管網的入口，沿途經過兩條用戶支路，最后流向管網出口，如圖4 所示，該管網系統總共有5 條管道，這些管道的幾何參數如表1 所示。在定義邊界條件的過程中，可以給定壓力或流量，我們采用給定壓力邊界的方式，壓力大的一端即為入口。

圖4 示例仿真管網示意圖

表1 仿真測試結構參數表

在圖4 所示的管網中，計算從管網入口開始，從給定的管網入口壓力計算出管道1 的出口，即節點1 的壓力。對于管道2，其入口壓力指定為管道1 的出口壓力，管道3 的入口壓力指定為管道2 的出口壓力，對其他管道的操作是同樣的。在壓力迭代完成后，依次對流量和焓值進行迭代計算，重復該程序直到管網中各狀態參數穩定。

采用Fluent 計算實驗仿真平臺驗證我們的計算結果，Fluent 計算實驗仿真平臺通過對給定的流體模型進行流域建立、網格劃分及求解計算，得出流體模型內的流場特征及流場分布，其求解結果比較精確，可以作為實驗對照組。由于管道較短，溫度降低不是很顯著，所以采用絕熱模型。

表2 給出了兩種方式的模擬結果，在1.3 MPa的入口壓力下，我們的管網計算模型與計算機實驗仿真結果流量十分接近，從結果分析，管道4 相比于管道5 更接近入口，在同樣的結構條件下，分配的流量更多，在兩種計算結果中顯示了這一點。管道4 和管道5 的流量在Fluent 中的流量相差較小，而計算模型的流量相差較大，這是由于在Fluent 中流體流動更接近于實際狀態，在管道三通處會存在壓力分布不均、壁面處會產生邊界層分離等問題，這會導致流阻比我們的計算模型更大，所以管道4 流量比計算模型低。而對于管道5 的流量，計算模型的流量低于Fluent 仿真，這是由于計算模型采取壓力邊界條件，在同樣的供汽壓力下，程序中的流量為275.2 kg/s，低于Fluent 的流量278.5 kg/s，且在管道4 處分流較多，導致管道5 的流量偏低。

表2 Fluent 計算機實驗仿真和計算模型的結果對比

結果表明，計算模型所得出的結果與計算機實驗仿真結果接近，存在的誤差主要由于管網結構在微觀細節上的不同，這一點在計算方法中需要加入局部阻力部件的管損模型，就可以使得我們構建的管網模型更精確。因此，在這個實例中，驗證了我們所構建的蒸汽管網動態仿真算法的有效性。

3 結果和討論

為了測試和驗證耦合算法程序在實際運行規則下的計算結果，在模擬計算中使用了中國臨港某熱電廠蒸汽管網數據，管網現場部分情況如圖5 所示，管道布局長度在3 km 左右，有13 個用戶。在一系列模擬計算中，對比了程序的運行和實際參數的結果。

圖5 臨港蒸汽管網現場圖

結合網絡元與網絡拓撲理論，對現場蒸汽管網結構進行重構，可以得到如圖6 所示蒸汽管網概念模型，其中包含一個熱電廠，十三家用戶，用戶編號是1 到13。每條管道的預期流向都被標注，網絡沒有回路，本文的耦合算法可以應用在這樣的管網結構上。

圖6 臨港熱電廠蒸汽管網結構圖

取一組臨港某熱電廠用戶開閥較多的某工作日時刻蒸汽管網運行狀態作為參照，熱電廠蒸汽入網壓力為1.36 MPa，溫度為243°C，流量為53.239 kg/s，在程序中建立管網模型，供汽入口處的參數取管網入口實際參數，進行動態仿真計算，結果如表 3所示。

表3 動態仿真結果與運行數據對比

計算模型與現場運行數據對比見圖7。從表3 與圖7 可以進行分析，由于動態仿真計算考慮了蒸汽流動過程中的狀態變化，將壓力、焓值在迭代過程中不斷更新與反饋，溫度與流量的計算結果與實際運行數據吻合度較好，溫度誤差多在5%以內，小流量的情況下流量的波動使得誤差結果多在10% 以下。值得關注的是壓力的對比情況，可以發現計算模型的壓力由于沿途的管損不斷降低，而現場壓力計所采得的壓力曲線波動較大，主路管線的壓力并沒有呈現出合理的壓降關系。某些用戶不僅在壓力，且在溫度上明顯與周圍廠家差異較大，這可能是由現場壓力計的讀數波動或者數據讀取延時所造成的，我們有理由認為這部分壓力或溫度可能不準確。

圖7 計算模型與現場運行數據對比圖

結果表明，用戶在各工況下的流量仿真結果均貼近實際測量數據，但計算過程中仍存在一定的偏差，這是由于構建的管網結構較為簡單，現場管網存在多彎管、多閥門，管網結構老化、凝結水等一系列問題，針對這種問題，需要對現場實際情況進行調查，在管網模型中添加流阻變量，并不斷調整相關參數，使得管網結構更趨向于實際結構，仿真結果更準確。

4 結論

本文提出了一種新的蒸汽供熱網絡動態仿真方法，它利用了內部流體參數的守恒關系，用于模擬流體管路的溫度、壓力、流量的非穩態特性，具有計算迅速、參數可控等特點，可解決目前熱電聯產過程中蒸汽管網的供汽計算忽略熱力損失對管網壓力影響所導致計算結果誤差偏大的問題。在計算過程中，由于局部散熱，壓力突變等問題，蒸汽輸送過程中的熱力參數比較多變，因此水蒸氣的熱力性質比單純的供水管網更為復雜，須采用實際氣體狀態方程求解蒸汽流動過程中的狀態參數，本文采用IAPWSIF97 工業標準中過熱蒸汽區參數狀態方程，用于求解蒸汽流動過程中密度、壓縮因子等參數。

與實際蒸汽管網運行數據作對比，結果較為吻合。示例驗證表明，管網的動態仿真計算中壓力、溫度與流量相對誤差在10%以下，計算結果精度較高，可以應用于管網狀態監測、運行故障診斷、管網優化設計等方面。進一步的工作在于管網實際結構的精確化建模，包括阻力部件、三維空間上管線布局帶來的損失等，從而可以使得結果更具參考價值。