中心差分法在燃?xì)夤艿绖?dòng)態(tài)模擬中的應(yīng)用*

高蘭周 侯根富 王 東

0 引言

天然氣長(zhǎng)輸管線不穩(wěn)定流動(dòng)的基本方程組可以利用分離變量的方法求出解析解,但是要通過(guò)簡(jiǎn)化和復(fù)雜的數(shù)學(xué)變換,求解比較復(fù)雜。而且求解的結(jié)果與實(shí)際運(yùn)行中的管道有一些差別。如果利用數(shù)值求解,就不存在上述問(wèn)題。

在對(duì)本問(wèn)題的數(shù)值解法中,有特征線法、隱式差分法和顯示法等。對(duì)于特征線法,可將管道上的偏微分方程化為特征差分方程,不用求解龐大的非線性方程組,占用較少的計(jì)算機(jī)內(nèi)存,易于求解,計(jì)算速度快,但為了滿足求解的穩(wěn)定性,時(shí)間層次往往取得很小。然而,對(duì)于天然氣管網(wǎng)系統(tǒng)而言,由于氣體的可壓縮性,其出現(xiàn)不穩(wěn)定程度要比液體管道出現(xiàn)的不穩(wěn)定程度小得多,因而時(shí)間層次取得過(guò)小沒(méi)多大的意義,還會(huì)浪費(fèi)計(jì)算時(shí)間。對(duì)于隱式差分法,由于要將待求時(shí)間層次上所有的未知量聯(lián)立起來(lái)進(jìn)行求解,需要求解的非線性方程組非常龐大,求解需要較長(zhǎng)的時(shí)間。但是這種方法能夠保證求解的穩(wěn)定性,時(shí)間步長(zhǎng)可取得較大[1]。

1 數(shù)學(xué)模型的建立

天然氣長(zhǎng)輸管線的管徑在很長(zhǎng)的距離上是不變的,因此垂直于流線的氣體特性變化率同沿著流線的變化率相比可忽略不計(jì),而且管道的曲率半徑比管道的直徑大得多,因此流動(dòng)可以看作是一元流動(dòng),在管道任意一截面上所有流體參數(shù)都可以看作是均一的,為該斷面的平均值。本文采用以下數(shù)學(xué)模型[1,2]。

其中,ρ為燃?xì)饷芏?kg/m3;t為時(shí)間,s;A為管段的流通橫截面面積,m2;qm為質(zhì)量流量,kg/s;x為管段軸向長(zhǎng)度,m;p為氣體壓力,Pa;λ為摩阻系數(shù);D為管道內(nèi)徑,m;g為重力加速度,m/s2;θ為燃?xì)夤艿缹?duì)水平面的傾斜角,rad;h為比焓,J/kg;K為傳熱系數(shù),W/(m2?K);T為燃?xì)鉁囟?K;T0為土壤溫度,K。

2 差分格式的建立和方程組封閉的分析

非穩(wěn)態(tài)氣體狀態(tài)方程(1)～方程(3)及實(shí)際氣體狀態(tài)方程p=p(ρ,T)和比焓方程 h=h(ρ,T)可以用更一般的形式來(lái)表示,即為=0。此方程為一階雙曲型方程。對(duì)于此形式的偏微分方程,考慮選用有時(shí)步和距步二階精度的Wendoff差分格式[3,4]。

為了求解所形成的差分方程,還需確定輸氣的邊界條件和初始條件。邊界是對(duì)系統(tǒng)而言的,現(xiàn)定義一段管子為研究對(duì)象,管子兩端的連接情況及相應(yīng)的運(yùn)行控制情況為邊界條件。對(duì)一大型輸氣管網(wǎng),邊界條件可能是多種多樣的,主要有:氣源、一個(gè)管子的終點(diǎn)且有分氣量、管子與管子的相交節(jié)點(diǎn)且有分氣量,如忽略節(jié)點(diǎn)的局部阻力和速度的差別則各支管節(jié)點(diǎn)壓力和溫度相等、節(jié)點(diǎn)控制最小壓力[5]。

結(jié)合一典型的輸氣工藝流程進(jìn)行封閉性分析[6]。

該輸氣干線中由4支管段,5個(gè)節(jié)點(diǎn)組成,其中節(jié)點(diǎn)1是氣源,節(jié)點(diǎn)3,4,5是分氣點(diǎn),節(jié)點(diǎn)2是普通節(jié)點(diǎn),將每一段管線進(jìn)行剖分,得截面為n1,n2,n3,n4,因?yàn)槊總€(gè)截面有3個(gè)未知變量p,qm,T,所以共有未知變量3?(n1+n2+n3+n4)。

可建立的差分方程數(shù)為:3?(n1-1)+3?(n2-1)+3?(n3-1)+3?(n4-1)。

按上面的邊界條件確定方法,上游邊界條件有兩個(gè),即氣源的壓力與溫度,下游分氣點(diǎn)的邊界條件有四個(gè),即例中的 qm(t)或p≥pmin,不能兩者兼有,對(duì)節(jié)點(diǎn)2和節(jié)點(diǎn)3可建立節(jié)點(diǎn)壓力、溫度或密度和流量3個(gè)平衡方程,共計(jì)6個(gè)方程。

因此共有邊界條件:2+4+6=12;

因?yàn)?3(n1-1)+3(n2-1)+3(n3-1)+3(n4-1)+12-3×(n1+n2+n3+n4)=0。

方程的個(gè)數(shù)與未知數(shù)的個(gè)數(shù)相等,因此方程是封閉的。

長(zhǎng)輸管線不穩(wěn)定流動(dòng)數(shù)學(xué)模型的中心有限差分法計(jì)算格式及相應(yīng)的初始條件,所形成的是一個(gè)非線性、非齊次的方程組。對(duì)一個(gè)剖面就有5個(gè)方程,那么對(duì)一個(gè)大型的輸氣管網(wǎng)將有多個(gè)方程,因此在選擇計(jì)算方程組的方法上,應(yīng)該能夠快速逼近真值,考慮到該方程組易求導(dǎo)和迭代初值離真值的偏差較小,故考慮選用Newton-Raphson迭代法求解[3]。

3 計(jì)算實(shí)例分析

某市高壓管網(wǎng)半環(huán)(門(mén)站—調(diào)壓站 B—調(diào)壓站 C)管道長(zhǎng)40 km,管道內(nèi)徑為 700 mm,土壤溫度為278.15 K,氣體常數(shù)為8.314 5 J/(mol?K),已知各個(gè)管段步長(zhǎng)節(jié)點(diǎn)的初始時(shí)刻質(zhì)量流量、密度和溫度以及調(diào)壓站B和調(diào)壓站C 24 h的流量數(shù)據(jù),氣源的密度和溫度恒定。

模擬時(shí)間為1 d,模擬程序用Visual C++語(yǔ)言編制,時(shí)間步長(zhǎng)取3 600 s,管段步長(zhǎng)計(jì)算網(wǎng)格為40個(gè),即空間步長(zhǎng)為1 km。選取部分時(shí)刻,其計(jì)算結(jié)果如下:

計(jì)算結(jié)果分析:圖1～圖3分別為天然氣的密度、溫度和壓力在同一時(shí)刻隨距離的分布。從圖1～圖3中可以看出天然氣的密度、溫度和壓力并不是線性分布的。它們?cè)?0 km處產(chǎn)生明顯而劇烈的變化,這是因?yàn)樵?0 km處存在分氣點(diǎn),從而氣體壓力變化劇烈。同時(shí)氣體在此處有明顯的節(jié)流效應(yīng)產(chǎn)生,因而溫度變化明顯。圖4,圖5分別反映了天然氣在同一點(diǎn)處的壓力隨時(shí)間變化的情況。由于用戶(hù)用氣量周期性變化,使得管網(wǎng)各處的壓力也是周期性變化的,該變化趨勢(shì)與用氣量的變化趨勢(shì)相反。由于在12:00及19:00,用戶(hù)用氣量最大,在此時(shí)刻,管道中的壓力卻最小;相反在用氣量的低峰時(shí),管道內(nèi)的壓力相對(duì)最大。另外,在分氣點(diǎn)前,其壓力變化的幅度小,而在分氣點(diǎn)后,氣體壓力變化的幅度大。圖6反映的是分氣點(diǎn)處的流量變化。從圖6中可以看出,在8:00,12:00及19:00等用戶(hù)用氣量較大,在3:00,15:00等用戶(hù)的用氣量較小,與實(shí)際用氣規(guī)律相符。

4 結(jié)語(yǔ)

采用了隱式差分法,對(duì)不穩(wěn)定流動(dòng)模型進(jìn)行求數(shù)值解。通過(guò)研究可以得出如下結(jié)論:

1)數(shù)值解法對(duì)長(zhǎng)輸管線不穩(wěn)定流動(dòng)進(jìn)行研究。不用對(duì)模型方程組進(jìn)行簡(jiǎn)化和復(fù)雜的數(shù)學(xué)變換,因而求解出來(lái)的結(jié)果與實(shí)際運(yùn)行中的數(shù)據(jù)更為接近;

2)求解時(shí)要將待求的時(shí)間層次上所有的未知量聯(lián)立起來(lái)同時(shí)進(jìn)行求解,求解的非線性方程組非常龐大,方程組過(guò)大時(shí)求解需要較長(zhǎng)的時(shí)間。為了保證數(shù)值解的絕對(duì)穩(wěn)定性,因而時(shí)間、空間步長(zhǎng)不能夠取的過(guò)小;

3)數(shù)值解考慮了氣體與外界環(huán)境的熱交換,能夠計(jì)算出天然氣沿線的溫度分布;

4)由于供氣與用氣的波動(dòng),導(dǎo)致管網(wǎng)中氣體的壓力、流量與溫度也是隨時(shí)間波動(dòng)的。如果管網(wǎng)存在有分氣點(diǎn)時(shí),在分氣點(diǎn)處氣體的壓力、溫度、密度的變化最為劇烈。

[1] 李長(zhǎng)俊,曾自強(qiáng),江茂澤.天然氣在管道系統(tǒng)中不穩(wěn)定流動(dòng)的分析[J].天然氣工業(yè),1994,14(6):54-59.

[2] 隋元春.不穩(wěn)定工況下燃?xì)夤芫W(wǎng)的模擬計(jì)算[J].煤氣與熱力,1992(1):36-41.

[3] 翟瑞彩,謝偉松.數(shù)值分析[M].天津:天津大學(xué)出版社,2000.

[4] 李長(zhǎng)俊.天然氣管道輸送[M].北京:石油工業(yè)出版社,2006.

[5] 周 游.天然氣長(zhǎng)輸管線末段與城市高壓外環(huán)不穩(wěn)定流動(dòng)及末段儲(chǔ)氣的研究[D].濟(jì)南:山東建筑大學(xué)碩士學(xué)位論文,2004.

[6] 唐建峰,段長(zhǎng)貴,呂文哲,等.特征線法在燃?xì)夤艿绖?dòng)態(tài)模擬中的應(yīng)用[J].油氣儲(chǔ)運(yùn),2001,20(8):12-17.