薛向東，宋永紅，張 博，朱 峰，閻 濤

1.西安交通大學(xué)軟件學(xué)院，陜西西安 710049

2.國(guó)家管網(wǎng)集團(tuán)科學(xué)技術(shù)研究總院分公司，河北廊坊 065000

3.國(guó)家管網(wǎng)集團(tuán)建設(shè)項(xiàng)目管理公司，河北廊坊 065000

天然氣管網(wǎng)系統(tǒng)的穩(wěn)態(tài)組分占比計(jì)算是指，在穩(wěn)態(tài)工況下，以水力仿真所得氣體體積流量與流向參數(shù)為基礎(chǔ)，通過(guò)穩(wěn)態(tài)組分跟蹤模型求解管網(wǎng)內(nèi)各節(jié)點(diǎn)與任意管道處的每種組分的占比。天然氣管網(wǎng)的組分跟蹤在管網(wǎng)調(diào)控中具有重要意義。多氣源混輸管網(wǎng)是未來(lái)天然氣管網(wǎng)的發(fā)展趨勢(shì)之一[1-2]，由于不同氣源的氣體組成不同，這意味著管網(wǎng)內(nèi)將發(fā)生頻繁的組分摻混，管網(wǎng)中組分占比也難以確定。由于天然氣中組分的變化會(huì)導(dǎo)致氣體熱值的變化，因此需要進(jìn)行組分占比計(jì)算以確保用戶用氣質(zhì)量的穩(wěn)定。此外，某些組分如H2與H2S等會(huì)對(duì)管道造成氫脆或腐蝕等危害[3-4]，同樣需要計(jì)算這些組分在管網(wǎng)中所有管道內(nèi)的占比以確保管道安全。因此，天然氣管網(wǎng)穩(wěn)態(tài)組分占比計(jì)算對(duì)于管網(wǎng)調(diào)控具有重要意義。

目前，國(guó)際上關(guān)于多氣源管網(wǎng)的組分與熱值占比仿真已經(jīng)取得一些成果。2001 年，德國(guó)在PSIG 年會(huì)上報(bào)告了一種在VNG-Verbundnetz 天然氣股份公司和天然氣進(jìn)口商中應(yīng)用的管道熱值模擬方案，作為管道熱力計(jì)費(fèi)的依據(jù)[5]。FLUENT 也可以實(shí)現(xiàn)對(duì)多氣源混輸管網(wǎng)的組分分布的穩(wěn)態(tài)仿真[6]。除此之外，還有一些國(guó)外學(xué)者進(jìn)行了多氣源混輸管網(wǎng)組分分布的研究。比如，Pellegrino[7]研究了穩(wěn)態(tài)工況下向天然氣管網(wǎng)中混輸綠色燃料（比如氫氣或生物甲烷等）的管網(wǎng)仿真模型；該模型的一個(gè)特點(diǎn)是定義了一個(gè)名為“偽電導(dǎo)”的參數(shù)，將管道內(nèi)流體流動(dòng)的動(dòng)力學(xué)問(wèn)題用電流來(lái)模擬；但是，該方法需要分別對(duì)管道、壓氣站、調(diào)壓站建模，因此模型較為復(fù)雜。同時(shí)，該模型在求解時(shí)需要不斷迭代各節(jié)點(diǎn)壓力、溫度與各節(jié)點(diǎn)組分，當(dāng)初始值選取不當(dāng)時(shí)極易出現(xiàn)不收斂的情況，因此數(shù)值穩(wěn)定性差。Cavana[8]等采用Pellegrino 提出的框架研究了穩(wěn)態(tài)工況下小型燃?xì)夤芫W(wǎng)中混輸生物甲烷的氣體質(zhì)量跟蹤方法；與Pellegrino 相比，其簡(jiǎn)化了非管道元件的建模，但是，求解流程并沒(méi)有做出簡(jiǎn)化，依舊需要依次迭代管網(wǎng)內(nèi)壓力、組分與溫度參數(shù)并最終達(dá)到平衡。Abeysekera[9]通過(guò)兩個(gè)步驟完成分布式注入燃?xì)夤芫W(wǎng)的穩(wěn)態(tài)模擬，即首先建立一個(gè)名為“節(jié)點(diǎn)流分析順序”的算法模塊用以存儲(chǔ)管網(wǎng)中各節(jié)點(diǎn)的分支數(shù)量與流向，并對(duì)每個(gè)分支編號(hào)；在計(jì)算程序處，調(diào)用“節(jié)點(diǎn)流分析順序”模塊進(jìn)行氣體組分計(jì)算、熱值計(jì)算、氣源占比計(jì)算與壓力分布計(jì)算，并最終輸出結(jié)果；這種方法可以保證對(duì)復(fù)雜管網(wǎng)連接信息的提取，因此適用于復(fù)雜管網(wǎng)的計(jì)算，但存儲(chǔ)空間較大。Chaczykowski 等[10-12]研究了不同氣質(zhì)、氣源天然氣注入管網(wǎng)時(shí)氣體的組分跟蹤方法，可以跟蹤管道出口熱值、某些氣體組分（如C2H6、CO2等）的變化。該方法可以用于管網(wǎng)瞬態(tài)組分仿真，但是對(duì)于大型管網(wǎng)來(lái)說(shuō)，仍需要一種可以精確、簡(jiǎn)便、快速的求解方法。相比之下，國(guó)內(nèi)管網(wǎng)仿真軟件對(duì)于氣體成分、熱值仿真與氣源跟蹤領(lǐng)域存在不足。在穩(wěn)態(tài)管網(wǎng)組分仿真模型研究方面，童?？礫13]在傳統(tǒng)天然氣管網(wǎng)水力、熱力模型的基礎(chǔ)上，通過(guò)補(bǔ)充穩(wěn)態(tài)工況下氣體組分方程、熱值方程與氣源比例方程構(gòu)建了多氣源混輸天然氣管網(wǎng)的數(shù)學(xué)模型，并采用節(jié)點(diǎn)-節(jié)點(diǎn)法來(lái)創(chuàng)建管網(wǎng)結(jié)構(gòu)數(shù)組，從而描述管網(wǎng)的拓?fù)浣Y(jié)構(gòu)；然而，通過(guò)這種方法構(gòu)建的管網(wǎng)組分分布計(jì)算模型邏輯上較為復(fù)雜，且對(duì)于大型復(fù)雜管網(wǎng)來(lái)說(shuō)計(jì)算時(shí)間長(zhǎng)，效率不高。

本文提出了一種穩(wěn)態(tài)天然氣管網(wǎng)組分占比計(jì)算模型，該模型基于管網(wǎng)節(jié)點(diǎn)摻混方程與管道組分約束方程，對(duì)管網(wǎng)整體建立求解矩陣，并采用直接LU 分解的方法進(jìn)行求解。采用此方法建立的模型不僅可以求解任意結(jié)構(gòu)管網(wǎng)的組分占比，且由于是對(duì)管網(wǎng)整體進(jìn)行直接求解，因此在應(yīng)用于大型復(fù)雜管網(wǎng)時(shí)具有模型簡(jiǎn)明、計(jì)算效率高、求解速度快的特點(diǎn)。采用國(guó)外某商業(yè)軟件進(jìn)行檢驗(yàn)，驗(yàn)證了模型的準(zhǔn)確性。模型可以用于大型穩(wěn)態(tài)天然氣管網(wǎng)的組分占比計(jì)算。

1 穩(wěn)態(tài)組分分布計(jì)算模型

1.1 基本假設(shè)

對(duì)于穩(wěn)態(tài)工況，管網(wǎng)中所有流體處于平衡狀態(tài)。同時(shí)，可以將管網(wǎng)內(nèi)每條管道中的流體看作一整個(gè)含有恒定組分的大批次。因此，對(duì)于穩(wěn)態(tài)管網(wǎng)，可以做出3 個(gè)基本假設(shè)[14-19]：其一，各個(gè)氣源的組分不隨時(shí)間發(fā)生變化；其二，每條管道內(nèi)的天然氣組分保持一致；其三，管道內(nèi)的天然氣在節(jié)點(diǎn)處均勻混合后從節(jié)點(diǎn)流出。

1.2 約束方程

在穩(wěn)態(tài)管網(wǎng)組分占比計(jì)算模型中，約束方程共有兩類：一類是內(nèi)部節(jié)點(diǎn)（除氣源點(diǎn)外的管網(wǎng)節(jié)點(diǎn)）的組分摻混方程，用于計(jì)算管網(wǎng)內(nèi)部節(jié)點(diǎn)的組分分布；另一類是管道元件的組分約束方程，用于計(jì)算管道元件中的組分分布。

1.2.1 管道元件組分約束方程

對(duì)于穩(wěn)態(tài)模型，可認(rèn)為管網(wǎng)中每條管道內(nèi)組分恒定不變，因此對(duì)于每條管道只需要一個(gè)約束方程來(lái)約束其組分。由于管道與節(jié)點(diǎn)連接，且流體是連續(xù)的，因此認(rèn)為與節(jié)點(diǎn)相連的所有下游管道內(nèi)的組分分布一致?？梢缘玫皆M分約束方程，針對(duì)P個(gè)管道元件可以建立P個(gè)約束方程：

式中：αP,λ為第P個(gè)管道元件中λ組分的摩爾分?jǐn)?shù)，αupP,N(E),λ為第P個(gè)管道元件中的上游節(jié)點(diǎn)N或上游非管道元件E處λ組分的摩爾分?jǐn)?shù)。

1.2.2 節(jié)點(diǎn)摻混方程

對(duì)每一個(gè)內(nèi)部節(jié)點(diǎn)采用組分完全摻混假設(shè)，可以得到組分摻混方程：

式中：αN,λ為節(jié)點(diǎn)N中λ組分的摩爾分?jǐn)?shù)；QP為流入管道P的體積流量，m3/s。

1.3 閉合關(guān)系

除了管道中的組分約束方程和節(jié)點(diǎn)中的摻混方程外，還需要補(bǔ)充邊界條件方程和非管道元件的約束方程，以使方程組閉合。一方面，各氣源點(diǎn)中各組分的摩爾分?jǐn)?shù)為已知量，該條件可以作為邊界條件方程；另一方面，還需要引入第E個(gè)非管道元件的約束方程，針對(duì)第E個(gè)非管道元件可以建立E個(gè)約束方程，即：

式中：、分別為第E個(gè)非管道元件的上游、下游節(jié)點(diǎn)處λ組分的摩爾分?jǐn)?shù)。

1.4 矩陣構(gòu)建

對(duì)于一個(gè)管網(wǎng)系統(tǒng)，假設(shè)其具有n個(gè)節(jié)點(diǎn)（其中有s個(gè)氣源點(diǎn)），e個(gè)非管道元件，p個(gè)管道元件，則對(duì)于組分i，基于管網(wǎng)整體可以聯(lián)立出s個(gè)邊界條件方程、（n-s）個(gè)內(nèi)部節(jié)點(diǎn)摻混方程、p個(gè)管道組分約束方程與e個(gè)非管道元件約束方程，并表示為如下矩陣形式：

Q為流量系數(shù)矩陣，其階數(shù)為（n+p）×（n+p），Q矩陣組成如下：

式中：Qj,k表示第j行第k列節(jié)點(diǎn)的流量。

如果節(jié)點(diǎn)m為氣源點(diǎn)，則：

如果節(jié)點(diǎn)m為非氣源點(diǎn)時(shí)，假設(shè)為非管道元件E的下節(jié)點(diǎn)，則：

式中：下標(biāo)f為以節(jié)點(diǎn)m為下節(jié)點(diǎn)的非管道元件E的上節(jié)點(diǎn)。

如果節(jié)點(diǎn)m為非氣源點(diǎn)時(shí)，假設(shè)k管道為節(jié)點(diǎn)m的上游管道，則：

式中：qm,k為與m節(jié)點(diǎn)相連的k管道中的體積流量，m3/s。

當(dāng)n+1≤i≤p+n時(shí)，Qi,j表示管道的流量矩陣元素，則：

α為待求稀疏矩陣，為一個(gè)（n+p）階列矩陣，其組成為：

式中：αiNn表示節(jié)點(diǎn)n中組分i的摩爾分?jǐn)?shù)，αiPp表示管道p中組分i的摩爾分?jǐn)?shù)。

b為（n+p）階列矩陣，其組成為：

式中：下標(biāo)m表示矩陣中的列。

當(dāng)1≤m≤n時(shí)，表示的是節(jié)點(diǎn)中組分i的b元素值；當(dāng)n+1≤m≤n+p時(shí)，表示的是管道中組分i的b元素值。當(dāng)1≤m≤n且m為氣源點(diǎn)時(shí)，若m節(jié)點(diǎn)即為氣源點(diǎn)j，則等于aij，aij為氣源j中組分i的摩爾分?jǐn)?shù)。

1.5 矩陣求解

對(duì)于上述管網(wǎng)組分分布矩陣，可以采用直接LU 分解的方法對(duì)Q矩陣進(jìn)行求解，分解求解公式如下：

式中：αij為原矩陣中的元素，這里為流量系數(shù)矩陣Q中的元素。因此，可以將原流量系數(shù)矩陣分解為：

式中：L和U分別為一個(gè)下三角矩陣和一個(gè)上三角矩陣，即：

再通過(guò)如下步驟進(jìn)行求解：

式中：y為采用直接LU 分解法求解時(shí)的一個(gè)中間矩陣。

式（17）與式（18）為矩陣向量表示的形式，因此將其改寫(xiě)為：

式（13）、式（19）與式（20）即為該穩(wěn)態(tài)組分分布計(jì)算模型的求解公式。

2 案例研究

2.1 案例分析

以圖1 所示的8 節(jié)點(diǎn)管網(wǎng)為例。假設(shè)該管網(wǎng)的水力工況如表1 所示；假設(shè)有節(jié)點(diǎn)1 與節(jié)點(diǎn)6 兩個(gè)氣源點(diǎn)，其組分構(gòu)成如表2所示。

表1 案例管網(wǎng)水力工況

表2 案例管網(wǎng)氣源組分構(gòu)成

圖1 案例管網(wǎng)示意

2.2 模型驗(yàn)證

基于C#語(yǔ)言邏輯，以本文提出的穩(wěn)態(tài)管網(wǎng)組分占比計(jì)算模型為理論基礎(chǔ)，采用微軟公司開(kāi)發(fā)的Visio_Studio 編程軟件編寫(xiě)8 節(jié)點(diǎn)任意結(jié)構(gòu)管網(wǎng)組分占比計(jì)算程序。

采用某商業(yè)軟件建立模型，分別利用如上8節(jié)點(diǎn)任意結(jié)構(gòu)管網(wǎng)組分占比計(jì)算程序與某商業(yè)軟件按照表1 與表2 的工況進(jìn)行求解，對(duì)本文提出的穩(wěn)態(tài)天然氣管網(wǎng)組分占比計(jì)算模型進(jìn)行驗(yàn)證，某商業(yè)軟件模型建立結(jié)果如圖2 所示。由于節(jié)點(diǎn)2 與節(jié)點(diǎn)7 分別和氣源點(diǎn)1 與氣源點(diǎn)6 通過(guò)單根管道相連，因此節(jié)點(diǎn)2 與氣源點(diǎn)1 的組分一致，節(jié)點(diǎn)7 與氣源點(diǎn)6 的組分一致。因此，在采用某商業(yè)軟件建立模型的過(guò)程中，省略氣源點(diǎn)1 與氣源點(diǎn)6，將節(jié)點(diǎn)2 與節(jié)點(diǎn)7 當(dāng)作氣源點(diǎn)進(jìn)行初始組分參數(shù)的設(shè)置。

圖2 案例管網(wǎng)采用某商業(yè)軟件建立的模型

除節(jié)點(diǎn)1、節(jié)點(diǎn)2、節(jié)點(diǎn)6、節(jié)點(diǎn)7 屬于單一輸入源的節(jié)點(diǎn)，其組分分布與設(shè)定組分一致外，其余節(jié)點(diǎn)的組分驗(yàn)證結(jié)果如表3所示。

表3 組分跟蹤模型驗(yàn)證結(jié)果

從表3 中可以看出，采用本文提出的穩(wěn)態(tài)天然氣管網(wǎng)組分占比計(jì)算模型分別計(jì)算包括節(jié)點(diǎn)3、節(jié)點(diǎn)4、節(jié)點(diǎn)5、節(jié)點(diǎn)8 在內(nèi)的四個(gè)摻混節(jié)點(diǎn)的組分分布，與某商業(yè)軟件模擬的誤差均在1%以內(nèi)。其中，四個(gè)摻混節(jié)點(diǎn)CH4組分的平均計(jì)算誤差為0.041%，四個(gè)摻混節(jié)點(diǎn)C2H6組分的平均計(jì)算誤差為-0.397%，四個(gè)摻混節(jié)點(diǎn)C3H8組分的平均計(jì)算誤差為-0.504%，誤差均在可接受范圍內(nèi)，故此認(rèn)為本文提出的模型具有可行性。

2.3 結(jié)果分析

基于本文提出的穩(wěn)態(tài)管網(wǎng)組分占比計(jì)算模型，采用8節(jié)點(diǎn)任意結(jié)構(gòu)管網(wǎng)組分占比計(jì)算程序?qū)Π咐芫W(wǎng)進(jìn)行求解，得到天然氣管網(wǎng)穩(wěn)態(tài)工況下各節(jié)點(diǎn)與管道中的組分占比（摩爾分?jǐn)?shù)），求解結(jié)果如表4所示。

表4 組分占比模擬結(jié)果

從表4 可以看出，除氣源點(diǎn)外，節(jié)點(diǎn)5 的甲烷占比最高，節(jié)點(diǎn)8 與節(jié)點(diǎn)4 中的甲烷占比次之。這是因?yàn)楣?jié)點(diǎn)5 中氣源6 的比例最高，而節(jié)點(diǎn)8 與節(jié)點(diǎn)4 的氣源6 占比次之，而在案例工況中氣源6 的甲烷組分濃度高于氣源1中的甲烷組分濃度，因此當(dāng)內(nèi)部節(jié)點(diǎn)的氣源占比構(gòu)成中氣源點(diǎn)6的占比較高時(shí)，其甲烷組分的濃度就較高。

從表4 還可以看出，管道3-4 的組分構(gòu)成與節(jié)點(diǎn)3一致；管道5-4、管道5-3的組分構(gòu)成與節(jié)點(diǎn)5一致；管道4-8與節(jié)點(diǎn)4的組分占比構(gòu)成一致。

從圖1 中可以看出，節(jié)點(diǎn)3 是管道3-4 的上游節(jié)點(diǎn)；節(jié)點(diǎn)5 是管道5-4、管道5-3 的上游節(jié)點(diǎn)；節(jié)點(diǎn)4 是管道4-8 的上游節(jié)點(diǎn)，這符合管道元件的組分約束。

通過(guò)如上分析，證明該模型可以正確地表示任意結(jié)構(gòu)管網(wǎng)中的拓?fù)浣Y(jié)構(gòu)。

3 結(jié)論

本文提出了一種天然氣管網(wǎng)穩(wěn)態(tài)組分占比計(jì)算模型，基于上述模型編寫(xiě)任意結(jié)構(gòu)8節(jié)點(diǎn)管網(wǎng)組分占比計(jì)算程序，并采用國(guó)外某商業(yè)軟件針對(duì)一個(gè)管網(wǎng)案例進(jìn)行模型驗(yàn)證。本文得出以下結(jié)論。

1）提出了一種穩(wěn)態(tài)天然氣管網(wǎng)組分占比計(jì)算模型。模型在基于水力計(jì)算結(jié)果的基礎(chǔ)上，以管網(wǎng)整體為對(duì)象，構(gòu)建管網(wǎng)組分分布矩陣，并采用直接LU 分解的方法進(jìn)行求解，得到任意結(jié)構(gòu)穩(wěn)態(tài)天然氣管網(wǎng)中各節(jié)點(diǎn)與管道元件內(nèi)每種天然氣組分的摩爾分?jǐn)?shù)。

2）利用國(guó)外某商業(yè)軟件通過(guò)一個(gè)管網(wǎng)案例進(jìn)行了模型驗(yàn)證，模型求解結(jié)果與某商業(yè)軟件求解結(jié)果間相對(duì)誤差均低于0.6%，具有可行性。

3）本文提出的穩(wěn)態(tài)管網(wǎng)組分占比計(jì)算模型可以用于任意結(jié)構(gòu)、任意規(guī)模的穩(wěn)態(tài)天然氣管網(wǎng)組分占比計(jì)算。在應(yīng)用于大規(guī)模穩(wěn)態(tài)天然氣管網(wǎng)的組分分布計(jì)算時(shí)，該模型具有求解速度快、泛用性強(qiáng)、邏輯清晰、內(nèi)存占用低的特點(diǎn)。

建議在后續(xù)研究中，進(jìn)一步開(kāi)展對(duì)瞬態(tài)工況下管網(wǎng)組分跟蹤模型的研究。