基于AGA8和AGA10的天然氣物性參數(shù)計(jì)算的軟件設(shè)計(jì)

陳玉潔 沈昱明

摘? 要： 根據(jù)美國(guó)燃?xì)鈪f(xié)會(huì)（AGA）發(fā)表的AGA No.8及AGA No.10報(bào)告中提出的計(jì)算方程，采用向量法、機(jī)械求積法等對(duì)AGA推薦算法進(jìn)行了優(yōu)化;并通過MATLAB編程，計(jì)算了天然氣物性參數(shù)，包括：壓縮性系數(shù)、比熱、聲速、等熵指數(shù)、焓和熵;分析了AGA8中天然氣壓縮性系數(shù)隨壓力與溫度的變化關(guān)系。測(cè)試結(jié)果表明，計(jì)算誤差低于0.0003‰，優(yōu)化后程序運(yùn)行平均耗時(shí)僅150 ms左右。

關(guān)鍵詞： MATLAB;天然氣;物性參數(shù);AGA8-92DC方程;AGA10;優(yōu)化算法

中圖分類號(hào)： TP391.75? ? 文獻(xiàn)標(biāo)識(shí)碼： A? ? DOI：10.3969/j.issn.1003-6970.2020.09.027

本文著錄格式：陳玉潔，沈昱明. 基于AGA8和AGA10的天然氣物性參數(shù)計(jì)算的軟件設(shè)計(jì)[J]. 軟件，2020，41（09）：9699+132

【Abstract】： According to the calculation equations proposed in the AGA No. 8 and AGA No. 10 reports issued by the American Gas Association （AGA）， the AGA recommendation algorithm was optimized using the vector method and the mechanical quadrature method; and the physical properties of natural gas was calculated through MATLAB， including compressibility factor， heat capacity， speed of sound， isentropic index， enthalpy and entropy; the relationship between the compressibility factor of natural gas in AGA8 and the pressure and temperature was analyzed. The test results show that the calculation error is less than 0.0003‰， and the average running time of the optimized program is only about 150 ms.

【Key words】： MATLAB; Natural gas; Physical properties; AGA8-92DC equation; AGA10; Optimization

0? 引言

天然氣物性參數(shù)，特別是壓縮性系數(shù)的精確度，對(duì)天然氣流量計(jì)量影響巨大[1]。對(duì)于實(shí)際氣體，常采用PV=ZnRT來(lái)代替理想氣體狀態(tài)方程。通過引入壓縮性系數(shù)Z，以衡量實(shí)際氣體相對(duì)理想氣體的偏離程度。

目前，壓縮性系數(shù)的計(jì)算方法有RK公式、Redlieh- Kwong對(duì)比態(tài)方程[2]、SGERC-88方程、AGA8-92DC方程等。已有大量實(shí)驗(yàn)證實(shí)了AGA8-92DC方程的精確性，因此在行業(yè)中被廣泛采用。

文獻(xiàn)3、4、5實(shí)現(xiàn)了對(duì)壓縮性系數(shù)的編程計(jì)算，但計(jì)算結(jié)果與AGA8所給結(jié)果存在較大誤差;文獻(xiàn)6探討了AGA8-92DC方程計(jì)算值與實(shí)驗(yàn)值的差異;文獻(xiàn)7基于AGA8-92DC的狀態(tài)方程，提出了一種簡(jiǎn)便、準(zhǔn)確的計(jì)算天然氣臨界流質(zhì)量流量的方法。

AGA10雖然給出了計(jì)算其他物性參數(shù)的算法及源代碼，但其算法復(fù)雜、代碼量大，不利于芯片化運(yùn)行。本文在其基礎(chǔ)上進(jìn)行了拓展，增加求解了天然氣其他重要的物性參數(shù);優(yōu)化了AGA8、AGA10中的幾個(gè)偏微分的算法。計(jì)算結(jié)果與AGA8、AGA10中的結(jié)果完全吻合，也大大縮短了計(jì)算時(shí)間。經(jīng)測(cè)試，平均計(jì)算時(shí)間約為150 ms。

1? AGA8-92DC方程及其計(jì)算精度

在1992年6月召開的國(guó)際標(biāo)準(zhǔn)化組織（ISO）上，美國(guó)燃?xì)鈪f(xié)會(huì)（AGA）提出了計(jì)算天然氣壓縮性系數(shù)的AGA8-92DC方程（以下簡(jiǎn)稱AGA8方程）。AGA8方程是參照天然氣組分摩爾分?jǐn)?shù)來(lái)計(jì)算的擴(kuò)展維里方程。在甲烷含量不低于70%，溫度263 K～338 K，絕對(duì)壓力低于20 MPa的工作條件下，計(jì)算不確定度優(yōu)于± 0.2%，能滿足絕大部分的工業(yè)條件。

2? 優(yōu)化計(jì)算

輸入變量為21種天然氣組分的摩爾分?jǐn)?shù)、溫度及絕對(duì)壓力，MATLAB軟件編程計(jì)算后，輸出天然氣物性參數(shù)值，包括壓縮因子、比熱、聲速、等熵指數(shù)、焓和熵。其中，比熱、聲速、等熵指數(shù)、焓、熵的計(jì)算公式參照AGA No.10報(bào)告[8]。

計(jì)算輸入常量說明如下，詳見GB/T17747.2-1999[9]表B.1、表B.2、表B.3和表4：

——狀態(tài)方程參數(shù)，n=1～58

——天然氣組分特征參數(shù)

——二元交互作用參數(shù)

——比熱、焓和熵的計(jì)算參數(shù)

圖1給出了天然氣物性參數(shù)的計(jì)算流程，下文將說明計(jì)算過程采用的優(yōu)化方法。

2.1? 密度檢索法迭代求解壓縮因子

AGA8-92DC方程的計(jì)算關(guān)系如式（1）所示：

式（1）～（3）表明壓縮因子、摩爾密度、對(duì)比密度、絕對(duì)壓力一一對(duì)應(yīng)，計(jì)算壓縮因子即為求解方程。該方程為非線性方程，無(wú)法使用直接求解法求解。為方便說明，令：

本文采用對(duì)密度參量進(jìn)行迭代對(duì)壓縮因子求解。圖2為描點(diǎn)曲線圖，可見曲線單調(diào)，滿足二分法使用條件。一般工況條件下，Z處于0.4至1.2之間[3]，故迭代起始條件為：

當(dāng)滿足時(shí)，迭代結(jié)束，對(duì)應(yīng)的即為所求解。

2.2? 向量化求解偏導(dǎo)

天然氣物性參數(shù)計(jì)算公式較為復(fù)雜，偏導(dǎo)計(jì)算更為冗雜。向量化計(jì)算是一種特殊的并行化計(jì)算方式。由于計(jì)算機(jī)并行的硬件結(jié)構(gòu)特點(diǎn)，本文發(fā)現(xiàn)，對(duì)于偏導(dǎo)計(jì)算采用向量化方法，不但可以減少“for循環(huán)”，還可避免大量重復(fù)運(yùn)算，極大提高了計(jì)算效率。以下以第二維里系數(shù)偏導(dǎo)計(jì)算予以說明。

第二維利系數(shù)B，僅與天然氣溫度與組分有關(guān)，由式（6）計(jì)算。

由此可見，僅需對(duì)（n=1～18）構(gòu)成的三個(gè)向量進(jìn)行點(diǎn)積運(yùn)算即可。

第二維里系數(shù)偏導(dǎo)計(jì)算循環(huán)次數(shù)不多，優(yōu)化效果并不明顯。但后續(xù)壓縮性系數(shù)偏導(dǎo)計(jì)算以及積分計(jì)算循環(huán)量大，引入向量化思想則計(jì)算速度提升較為明顯。

2.3? 機(jī)械求積法計(jì)算積分

比熱、焓和熵的計(jì)算涉及求解關(guān)于摩爾密度的定值積分。被積函數(shù)復(fù)雜，包含累加和偏導(dǎo)計(jì)算，本文使用機(jī)械求積法求解。與定義法相比，機(jī)械求積法能夠劃分更少的區(qū)間而達(dá)到相同的計(jì)算精度，算法時(shí)間復(fù)雜度得以降低。物性參數(shù)計(jì)算方程為式（10）～（12）所示：

將被積函數(shù)參照AGA No.10報(bào)告[8]中式（3.12）和式（3.14）展開，其中常數(shù)項(xiàng)積分直接按照積分表求解，并記被積函數(shù)中關(guān)于摩爾密度的非常數(shù)積分項(xiàng)分別為Q1、Q2、Q3：

積分計(jì)算時(shí)，本文選用梯形公式代入。表1顯示當(dāng)選則不同區(qū)間長(zhǎng)度，計(jì)算相對(duì)誤差與平均耗時(shí)的變化情況。根據(jù)測(cè)試結(jié)果，取區(qū)間長(zhǎng)度為0.000002最為合適。

3? 軟件設(shè)計(jì)

GUI（Graphical User Interface），指利用圖形控件進(jìn)行人機(jī)交互操作。為了方便用戶的使用，本文引入MATLAB的GUI技術(shù)[10]，只需用戶在人機(jī)交互界面上輸入天然氣溫度、壓力和組分摩爾數(shù)，即可自行計(jì)算出物性參數(shù)。

圖形界面設(shè)計(jì)完成后，通過Matlab的自動(dòng)導(dǎo)出功能可將計(jì)算程序*.m文件與GUI圖形界面*.fig文件合并為一個(gè).m文件。隨后C/C++編譯器生成*.exe可執(zhí)行文件。未安裝MATLAB的機(jī)器僅需下載編譯環(huán)境庫(kù)MATLAB Runtime（MCR）或者將MATLAB安裝路徑下的/toolbox/compiler/Deploy/win32的MCRInstaller.exe與編譯的執(zhí)行文件打包發(fā)送即可，執(zhí)行文件運(yùn)行時(shí)會(huì)自動(dòng)調(diào)用下載的動(dòng)態(tài)鏈接庫(kù)DLL。

但MCRInstaller.exe所占內(nèi)存較大，可移植性不高，對(duì)應(yīng)用范圍有一定限制。既需MATLAB強(qiáng)大的計(jì)算能力，又需生成便捷的計(jì)算機(jī)軟件，也可采用其他編程語(yǔ)言調(diào)用MATLAB程序。如Mideva工具編譯*.m文件后可生成能被Visual C++[12]、Visual Basic以及Delphi調(diào)用的DLL動(dòng)態(tài)鏈接庫(kù)。Microsoft VB在軟件開發(fā)界面上具有獨(dú)到之處，但在計(jì)算方面略顯不足，這方面恰可與MATLAB互補(bǔ)[13]。圖3是本文軟件界面設(shè)計(jì)。

4? 測(cè)試

本文將軟件計(jì)算結(jié)果與GB/T 17747-1999與AGA No.10報(bào)告中算例進(jìn)行比較，計(jì)算結(jié)果中壓縮性系數(shù)完全一致，其他參數(shù)如表1所示，計(jì)算相對(duì)誤差小于0.0003‰，說明軟件計(jì)算正確，滿足精度要求。

MATLAB源程序中，本文利用tic， toc語(yǔ)句對(duì)程序進(jìn)行時(shí)間測(cè)試。不采用任何優(yōu)化方法計(jì)算（默認(rèn)使用循環(huán)結(jié)構(gòu)和MATLAB自帶積分函數(shù)），平均耗時(shí)為1.5 s;優(yōu)化后平均耗時(shí)為150 ms，降為原來(lái)的10%，優(yōu)化效果明顯。

表2中顯示的是AGA10更寬的應(yīng)用范圍，本文在AGA8原有數(shù)據(jù)的基礎(chǔ)上，計(jì)算了更高溫度、更大壓力條件的壓縮性系數(shù)。對(duì)于AGA8所展示的6種天然氣，其壓縮性系數(shù)與溫度、壓力的關(guān)系如圖4所示。

可以看出：對(duì)同一種類天然氣，壓力為6 MPa的壓縮性系數(shù)變化曲線與12 MPa情況下相差較大，而12 MPa與20 MPa曲線幾乎重合。

5? 結(jié)語(yǔ)

（1）本文采用向量化方法、機(jī)械求積法、密度檢索法迭代等優(yōu)化算法，極大提高了計(jì)算效率。

（2）計(jì)算結(jié)果與AGA報(bào)告數(shù)據(jù)相比，二者完全吻合。

（3）對(duì)同一種類天然氣，壓力為6 MPa的壓縮性系數(shù)變化曲線與12 MPa情況下相差較大;而12 MPa以上曲線幾乎重合;且隨溫度升高，不同壓力下天然氣壓縮性系數(shù)差值減小。

參考文獻(xiàn)

[1]孟祥適， 姜印平， 劉玉杰， 等. 基于天然氣壓縮系數(shù)Z提高天然氣計(jì)量準(zhǔn)確度的方法[J]. 測(cè)控技術(shù)， 2004， 23（6）： 16-17， 28.

[2]黃希揚(yáng)， 沈昱明. 基于Redlieh-Kwong對(duì)比態(tài)方程的天然氣壓縮因子計(jì)算[J]. 化工自動(dòng)化及儀表， 2019， 46（02）： 107-110.

[3]李丹華， 姜東琪. AGA8-92DC計(jì)算方法天然氣壓縮因子計(jì)算[J]. 煤氣與熱力， 2011， 31（03）： 43-46.

[4]梁光川， 左果， 文雯， 等. 應(yīng)用VB實(shí)現(xiàn)AGA8-92DC氣體壓縮因子計(jì)算方法[J]. 遼寧化工， 2014， 43（05）： 609-612.

[5]陳旻灝， 黃希揚(yáng)， 沈昱明. Matlab與LabVIEW混合編程的天然氣壓縮因子軟件設(shè)計(jì)[J]. 化工自動(dòng)化及儀表， 2019， 46（10）： 844-848.

[6]蘇中良， 劉曰武， 張俊清. AGA8-92DC方法計(jì)算天然氣偏差因子的研究[J]. 油氣井測(cè)試， 2010， 19（06）： 29-36+83.

[7]D.G. Stewart， J.T.R. Watson， A.M. Vaidya. A new correlation for the critical mass flux of natural gas mixtures[J]. Flow Measurement and Instrumentation， 2000， 11（4）.

[8]AGA Report No. 10 Speed of Sound in Natural Gas and Other Related Hydrocarbon Gases[S]. Transmission Measurement Committee. 2002.11.15.

[9]GB/T 17747.2-2011， 天然氣壓縮因子計(jì)算—第2部分： 用摩爾組成進(jìn)行計(jì)算（ISO 12213-2：2006）[S]. 2011， 12， 30.

[10]姚楠， 耿奇. 基于MATLAB GUI的車牌自動(dòng)識(shí)別系統(tǒng)設(shè)計(jì)[J]. 軟件， 2016， 37（3）： 44-46.

[11]曹妍， 陳偉. 基于MATLAB GUI的指紋圖像增強(qiáng)技術(shù)研究[J]. 軟件， 2018， 39（11）： 30-32.

[12]徐治. Visual C++調(diào)用MATLAB函數(shù)庫(kù)的混合編程技術(shù)[J]. 軟件， 2015， 36（2）： 55-58.

[13]曹妍， 陳偉. 基于MATLAB圖像處理的藥片計(jì)數(shù)方法研究[J]. 軟件， 2018， 39（9）： 13-15.

[14]劉艷華. 數(shù)字圖像置亂加密解密技術(shù)的研究與實(shí)現(xiàn)[J]. 軟件， 2018， 39（3）： 141-144.