, , ,

(華東理工大學(xué)承壓系統(tǒng)與安全教育部重點(diǎn)實(shí)驗(yàn)室,上海 200237)

天然氣物性數(shù)據(jù)庫(kù)的開(kāi)發(fā)及其在閥門(mén)仿真中的應(yīng)用

湯俊,于新海,李松,周邵萍

(華東理工大學(xué)承壓系統(tǒng)與安全教育部重點(diǎn)實(shí)驗(yàn)室,上海200237)

利用目前在天然氣應(yīng)用領(lǐng)域中較為全面和精確的GERG—2004方程建立天然氣數(shù)據(jù)庫(kù),使用FORTRAN程序,VBA宏語(yǔ)言和MATLAB編寫(xiě)了插值算法,在CFX仿真計(jì)算中調(diào)用了數(shù)據(jù)庫(kù)中物性數(shù)據(jù)。以天然氣迷宮式調(diào)節(jié)閥三維模型為研究對(duì)象,基于k-Epsilon流動(dòng)模型,調(diào)用天然氣物性數(shù)據(jù)庫(kù),實(shí)現(xiàn)了以天然氣為介質(zhì)的迷宮式調(diào)節(jié)閥仿真,為以天然氣為介質(zhì)的過(guò)流部件仿真提供了解決思路。

數(shù)據(jù)庫(kù); 迷宮式調(diào)節(jié)閥; 物性模型

迷宮式調(diào)節(jié)閥廣泛應(yīng)用于天然氣領(lǐng)域以及大型電站的鍋爐等非常惡劣的工況下,是流量控制和壓力控制的主要元件之一[1-2]。數(shù)值模擬是研究迷宮式調(diào)節(jié)閥內(nèi)流場(chǎng)設(shè)計(jì)、外特性分析和結(jié)構(gòu)評(píng)估優(yōu)化的主要工具。但天然氣組成的氣體種類(lèi)繁多,且其物理性質(zhì)隨著組分、溫度、壓強(qiáng)的變化而變化[3],不僅給實(shí)驗(yàn)物性參數(shù)研究帶來(lái)困難,也給利用天然氣物性進(jìn)行仿真計(jì)算帶來(lái)巨大的挑戰(zhàn)。國(guó)內(nèi)外流體商業(yè)軟件尚未建立全面、高精度、高穩(wěn)定性的天然氣數(shù)據(jù)庫(kù),并且未提出針對(duì)多種類(lèi)型天然氣的多相流數(shù)值仿真方法。

基于此,本文利用GERG—2004方程[4],結(jié)合VBA宏語(yǔ)言和MATLAB建立天然氣數(shù)據(jù)庫(kù),并以天然氣迷宮式調(diào)節(jié)閥三維模型為研究對(duì)象,基于k-Epsilon流動(dòng)模型,調(diào)用天然氣物性數(shù)據(jù)庫(kù),探究物性對(duì)于閥門(mén)性能的影響。

1 數(shù)據(jù)庫(kù)建立方法

1.1基于GERG—2004方程的數(shù)據(jù)庫(kù)建立

首先,利用GERG—2004方程建立氣相區(qū)和兩相區(qū)天然氣數(shù)據(jù)庫(kù)。然后用FORTRAN程序編寫(xiě)插值算法,數(shù)據(jù)庫(kù)的生成與插值程序流程如圖1所示。最后,在CFX仿真計(jì)算中調(diào)用數(shù)據(jù)庫(kù)中物性數(shù)據(jù),完成以天然氣為介質(zhì)的閥門(mén)內(nèi)流體仿真計(jì)算。為減少內(nèi)存消耗,提高程序運(yùn)算的速度,同時(shí)為避免預(yù)設(shè)范圍之外的點(diǎn)的插值,增強(qiáng)程序的魯棒性和精確度,本數(shù)據(jù)庫(kù)采用了動(dòng)態(tài)表形式。

1.2基于MATLAB程序獲得和采集數(shù)據(jù)

GERG—2004方程可計(jì)算包含甲烷、乙烷、二氧化碳、氮?dú)獾?8種天然氣組分混合的物性參數(shù),計(jì)算區(qū)域覆蓋了氣相區(qū)與液相區(qū)。

表1～表2為采用Excel程序和DOS程序?qū)σ欢ńM分、溫度、壓力和體積下物性的計(jì)算結(jié)果,分別采集單相區(qū)和兩相區(qū)數(shù)據(jù)。

圖1 數(shù)據(jù)庫(kù)生成與插值程序流程圖Fig.1 Program flow chart of database generation and interpolation

表1 GERG—2004方程Excel程序計(jì)算結(jié)果

表2 GERG—2004方程DOS程序計(jì)算結(jié)果

為高效獲取GERG—2004方程中的數(shù)據(jù),運(yùn)用VBA宏語(yǔ)言和MATLAB等工具實(shí)現(xiàn)數(shù)據(jù)的輸入、可視化處理和提取分析。在完成組分、壓力、溫度等數(shù)據(jù)輸入后,DOS程序?qū)⑸筛鹘M分、溫度和壓力工況下的物性數(shù)據(jù),數(shù)據(jù)采集流程如圖2所示。文本采用正則表達(dá)式,提取密度和氣液相體積分?jǐn)?shù)等參數(shù),數(shù)據(jù)以單個(gè)字符串來(lái)描述和匹配文本中字符串。最后,使用MATLAB整理歸納提取的物性數(shù)據(jù),完成數(shù)據(jù)采集。

1.3基于線性插值法的數(shù)據(jù)處理

對(duì)采集的數(shù)據(jù)進(jìn)行插值處理,所選用的二維插值方法為雙線性插值[5]。二維插值方法的實(shí)現(xiàn)過(guò)程如圖3所示。P1、P2、P3、P4是所需求解點(diǎn)P(x,y)四周的插值節(jié)點(diǎn)。

圖2 數(shù)據(jù)采集流程圖Fig.2 Flow chart of data acquisition

圖3 雙線性插值Fig.3 Flow chart of data acquisition

雙線性插值的基本函數(shù)形式如下:

f(x,y)=(ax+b)(cy+d)

(1)

式(1)中含有4個(gè)待定系數(shù)。通過(guò)4個(gè)插值節(jié)點(diǎn)的數(shù)據(jù)值,可求解出這4個(gè)待定系數(shù)。

首先,在x方向上進(jìn)行線性插值,得

(2)

(3)

然后,在y方向上進(jìn)行線性插值,得

(4)

整理得

f(x,y)=f(P3)(1-x)(1-y)+

f(P1)(1-x)y+f(P4)x(1-y)

(5)

線性插值的結(jié)果并不會(huì)因?yàn)椴逯淀樞虻牟煌淖?。本文三維插值算法也采用線性插值,與二維采用的雙線性插值方法相同,原理如圖4所示[6]。

圖4 三維線性插值Fig.4 Interpolation in 3-dimension

1.4CFX調(diào)用FORTRAN程序

將自定義的物理模型和物性特征寫(xiě)入FORTRAN程序中,求解器求解時(shí)將調(diào)用這些插值處理后的數(shù)據(jù)。CFX通過(guò)CEL函數(shù)調(diào)用自定義FORTRAN程序,流程如圖5所示。

編譯平臺(tái)搭建成功后,在命令行程序執(zhí)行cfx5mkext語(yǔ)句,將編寫(xiě)的F格式FORTRAN文件編譯成DLL動(dòng)態(tài)共享庫(kù)。求解器進(jìn)行求解時(shí),將通過(guò)在預(yù)處理器中預(yù)先定義的CEL語(yǔ)句和CEL函數(shù)對(duì)象對(duì)共享庫(kù)進(jìn)行調(diào)用和計(jì)算。

圖5 CFX調(diào)用FORTRAN程序流程Fig.5 Calling of user FORTRAN program by CFX

2 數(shù)據(jù)庫(kù)生成與實(shí)驗(yàn)驗(yàn)證

2.1組分、壓力、溫度范圍確定

在獲得初始數(shù)據(jù)后,需對(duì)天然氣數(shù)據(jù)庫(kù)組分、壓力、溫度的范圍進(jìn)行確定和節(jié)點(diǎn)數(shù)據(jù)進(jìn)行定義。本文數(shù)值模擬所選用天然氣組分?jǐn)?shù)據(jù)來(lái)源于四川盆地新場(chǎng)氣田,其天然氣組成如表3所示。

表3 四川新場(chǎng)氣田某層系天然氣組分及其摩爾分?jǐn)?shù)

本文選取迷宮式氣動(dòng)調(diào)節(jié)閥為研究對(duì)象,結(jié)構(gòu)如圖6所示。實(shí)現(xiàn)天然氣數(shù)據(jù)庫(kù)的調(diào)用,可完成以不同壓差工況下天然氣為介質(zhì)的閥內(nèi)流場(chǎng)仿真,有利于進(jìn)一步研究迷宮盤(pán)片內(nèi)部流動(dòng)在不同工況下的變化規(guī)律,對(duì)天然氣閥門(mén)設(shè)計(jì)性能進(jìn)行評(píng)估。其設(shè)計(jì)入口絕對(duì)壓力為34.9 MPa,設(shè)計(jì)出口絕對(duì)壓力為0.4 MPa,數(shù)據(jù)庫(kù)的壓力范圍需要包含最大、最小壓力,壓力邊界設(shè)定為0.01 MPa≤P≤36.01 MPa。天然氣的入口溫度為70 ℃。選取純CH4為介質(zhì)進(jìn)行模擬,CH4的物性數(shù)據(jù)來(lái)自于CFX中的材料庫(kù)。通過(guò)試算調(diào)節(jié)閥內(nèi)溫度場(chǎng),預(yù)估數(shù)據(jù)庫(kù)溫度范圍。針對(duì)調(diào)節(jié)閥的高壓差、高流速的工況,選用Total Energy熱傳遞模型和k-Epsilon湍流模型。由于高壓差的邊界條件不穩(wěn)健,對(duì)初始條件較為敏感,在模擬中先將預(yù)先估計(jì)的入口流量和出口靜壓作為邊界條件進(jìn)行求解,待收斂后,再將該結(jié)果作為初值代入兩端壓差的邊界條件的算例中模擬。由于盤(pán)片閥結(jié)構(gòu)比較復(fù)雜,采用ICEM軟件對(duì)其流道模型進(jìn)行整體非結(jié)構(gòu)網(wǎng)格劃分,對(duì)盤(pán)片、閥芯等關(guān)鍵位置的流體域進(jìn)行網(wǎng)格加密如圖7所示,網(wǎng)格總計(jì)3.58×106。同時(shí),在對(duì)網(wǎng)格進(jìn)行無(wú)關(guān)性驗(yàn)證的基礎(chǔ)上,分析閥門(mén)內(nèi)的流場(chǎng)情況。

1—閥桿;2—閥蓋;3—閥瓣;4—迷宮盤(pán)片閥芯;5—閥體;6—閥座

圖7 盤(pán)片處網(wǎng)格Fig.7 Mesh of the valve discs

圖8為調(diào)節(jié)閥中心對(duì)稱(chēng)溫度場(chǎng)云圖。從圖中可以看出,由于閥門(mén)的節(jié)流效應(yīng),天然氣在流經(jīng)閥門(mén)的過(guò)程中溫度有所降低,由入口管段的345.7 K降至出口管段的200 K左右。數(shù)據(jù)庫(kù)溫度需包含最大值和最小值,因此將其溫度范圍確定為180 K≤T≤350 K。

圖8 中心對(duì)稱(chēng)面溫度云圖Fig.8 Temperature contour of symmetric plane

2.2節(jié)點(diǎn)數(shù)確定

數(shù)據(jù)庫(kù)中包含的數(shù)據(jù)節(jié)點(diǎn)數(shù)越多,所得到插值的結(jié)果會(huì)越接近實(shí)際值。但數(shù)據(jù)節(jié)點(diǎn)過(guò)多將導(dǎo)致數(shù)據(jù)庫(kù)變得很大,模擬時(shí)也占用大量的RAM容量,使得程序運(yùn)行變慢。

為得到較為合適的節(jié)點(diǎn)數(shù),首先生成不同體積分?jǐn)?shù)的CH4和CO2的二元混合物的數(shù)據(jù)庫(kù),對(duì)其進(jìn)行插值,并將插值結(jié)果與已有的實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比。以CH4和CO2二元混合物的插值結(jié)果作為參考,擬定合適的節(jié)點(diǎn)數(shù)。初步將網(wǎng)格步長(zhǎng)設(shè)為10 K和1 MPa。以體積分?jǐn)?shù)為0.9的CH4和體積分?jǐn)?shù)為0.1的CO2的二元混合物為例。插值結(jié)果與文獻(xiàn)的實(shí)驗(yàn)結(jié)果[7-9]對(duì)比如表4所示。

表4 氣相區(qū)插值結(jié)果與實(shí)驗(yàn)結(jié)果對(duì)比

氣相區(qū)內(nèi)實(shí)驗(yàn)值與插值計(jì)算值的最大誤差為0.097 78%,可滿(mǎn)足模擬的精確度需求。因此,數(shù)據(jù)庫(kù)的溫度步長(zhǎng)和壓力步長(zhǎng)滿(mǎn)足條件。溫度方向網(wǎng)格節(jié)點(diǎn)為18個(gè),壓力方向網(wǎng)格節(jié)點(diǎn)為37個(gè),數(shù)據(jù)庫(kù)總節(jié)點(diǎn)數(shù)為666個(gè)。將這些網(wǎng)格節(jié)點(diǎn)的密度數(shù)據(jù)和線性插值算法編入FORTRAN程序,并在命令行程序中使用cfx5mkext語(yǔ)句將該FORTRAN程序編譯成DLL共享庫(kù)。

3 基于數(shù)據(jù)庫(kù)的閥門(mén)仿真

將所建立的氣相區(qū)數(shù)據(jù)庫(kù)應(yīng)用于天然氣迷宮式調(diào)節(jié)閥的數(shù)值模擬進(jìn)行計(jì)算。介質(zhì)的物性數(shù)據(jù)來(lái)源于數(shù)據(jù)庫(kù),選取體積分?jǐn)?shù)為0.25的CH4和體積分?jǐn)?shù)為0.75的CO2為仿真介質(zhì)。入口總壓設(shè)置為34.9 MPa,出口靜壓設(shè)置為0.4 MPa。選用了總能模型獲取調(diào)節(jié)閥內(nèi)的溫度場(chǎng)等信息,介質(zhì)的入口溫度為343.15 K。模擬結(jié)果如圖9所示。

圖9示出了壓力、溫度和密度在閥門(mén)內(nèi)部隨介質(zhì)流動(dòng)的變化,3個(gè)監(jiān)測(cè)值都在盤(pán)片流道中發(fā)生劇烈的變化。模擬計(jì)算得到入口平均密度為579.9 kg/m3,而由GERG—2004方程計(jì)算得到的密度為16.193 28 mol/L,換算得599.1 kg/m3,入口密度誤差約為3.205%。在不同的溫度和壓力工況下,對(duì)比多個(gè)監(jiān)測(cè)點(diǎn)的仿真結(jié)果與GERG—2004方程計(jì)算結(jié)果,如表5所示,發(fā)現(xiàn)密度誤差最大值為3.432%,插值計(jì)算誤差較小,仿真精確度滿(mǎn)足要求。同時(shí),數(shù)據(jù)庫(kù)在計(jì)算高壓差的工況調(diào)用計(jì)算過(guò)程中,表現(xiàn)出良好的適應(yīng)性、魯棒性和精確性。

圖9 對(duì)稱(chēng)面的壓力、溫度及密度云圖Fig.9 Symmetric plane pressure contour,temperature contour and density contour

表5 計(jì)算結(jié)果對(duì)比

4 結(jié) 論

以天然氣迷宮式調(diào)節(jié)閥為研究對(duì)象,建立天然氣物性數(shù)據(jù)庫(kù),使用CFX中的流動(dòng)模型,實(shí)現(xiàn)了高壓差仿真計(jì)算并得出以下結(jié)論:

(1) 數(shù)據(jù)庫(kù)密度實(shí)驗(yàn)值與插值計(jì)算值的最大誤差為0.098%,仿真結(jié)果與GERG—2004方程密度計(jì)算結(jié)果最大誤差為3.432%。因此,該方法得出的天然數(shù)據(jù)庫(kù)滿(mǎn)足計(jì)算的精確度和穩(wěn)定性要求。

(2) 創(chuàng)新性地調(diào)用GERG—2004方程和FORTRAN插值程序,實(shí)現(xiàn)了以天然氣為介質(zhì)的迷宮式調(diào)節(jié)閥仿真,為以天然氣為介質(zhì)的過(guò)流部件仿真提供了解決思路。

[1] 鄭建農(nóng),孫葉柱.進(jìn)口電站閥門(mén)選型分析[J].熱力發(fā)電,2007,36(2):83-84.

[2] 賈青,陶正良.從對(duì)我國(guó)電站調(diào)節(jié)閥的調(diào)研看其現(xiàn)狀及發(fā)展前景[J].華電技術(shù),2005,27(1):8-10.

[3] 胡見(jiàn)義,鄭俊章.天然氣形成與分布的自身屬性[C]//中國(guó)石油地質(zhì)年會(huì).北京:石油工業(yè)出版社,2005:45-52.

[4] KUNZ O,WAGNER W.The GERG—2008 wide-range equation of state for natural gases and other mixtures:An expansion of GERG—2004[J].Journal of Chemical & Engineering Data,2012,57(11):3032-3091.

[5] 龔昌來(lái),羅聰,楊冬濤.一種基于邊緣方向的雙線性插值方法[J].激光與紅外,2010,40(7):788-791.

[6] 楊代倫.基于質(zhì)心的二維和三維插值算法[J].計(jì)算機(jī)工程與應(yīng)用,2004,40(12):77-78.

[7] LI H,JAKOBSEN J P,WILHELMSEN ?,etal.PVTxyproperties of CO2mixtures relevant for CO2capture,transport and storage:Review of available experimental data and theoretical models[J].Applied Energy,2011,88(11):3567-3579.

[8] KESTIN J,YATA J.Viscosity and diffusion coefficient of six binary mixtures[J].Journal of Chemical Physics,1968,49(11):4780-4791.

[9] KESTIN J,RO S T,WAKEHAM W A.The transport properties of binary mixtures of hydrogen with CO,CO2,and CH4[J].Physica A Statistical Mechanics & Its Applications,1983,119(3):615-638.

DevelopmentandApplicationsinValveSimulationofNaturalGasPhysicalPropertyDatabase

TANGJun,YUXin-hai,LISong,ZHOUShao-ping

(KeyLaboratoryofPressurizedSystemandSafety,MinistryofEducation,EastChinaUniversityofScienceandTechnology,Shanghai200237,China)

Natural gas database is established by GERG—2004equation which is very comprehensive and accurate in the field of natural gas application at present.Interpolation algorithm is defined by user FORTRAN program,VBA macro language and MATLAB.Physical property data in database is called in the numerical simulation of CFX.This paper takes the natural gas labyrinth type valve as an example,calling physical property data in the numerical simulation based onk-Epsilon flow model,completing the simulation of labyrinth control valve using natural gas as medium,which provides a solution for the simulation of over-current components with natural gas as medium.

database; labyrinth control valve; physical property model

TH3

A

1006-3080(2017)05-0712-05

10.14135/j.cnki.1006-3080.2017.05.017

2016-11-22

湯 俊(1992-),男,江蘇淮安人,碩士生,主要從事泵閥的數(shù)值模擬與結(jié)構(gòu)優(yōu)化。

周邵萍,E-mail:shpzhou@ecust.edu.cn