摘 要 本文給出了我校《海洋油氣集輸》課程中計(jì)算機(jī)編程項(xiàng)目的詳細(xì)求解方法及參照答案，這些計(jì)算機(jī)編程項(xiàng)目在教學(xué)中應(yīng)用效果良好，能夠有效的幫助學(xué)生掌握油氣集輸課程的重要的理論，同時(shí)發(fā)展學(xué)生解決實(shí)際工程問題的能力?，F(xiàn)將實(shí)踐教學(xué)中使用的編程項(xiàng)目的解法及參考答案整理成文，方便相關(guān)專業(yè)的廣大師生在本科教學(xué)過程中使用。

關(guān)鍵詞 油氣集輸 編程項(xiàng)目 算例 閃蒸計(jì)算 管網(wǎng)計(jì)算

中圖分類號(hào)：G642文獻(xiàn)標(biāo)識(shí)碼：A

在《石油與天然氣工程挑戰(zhàn)性課程教學(xué)實(shí)踐——海洋油氣集輸計(jì)算機(jī)編程項(xiàng)目設(shè)計(jì)》一文中，介紹了中國石油大學(xué)（北京）在海洋油氣集輸全英文教學(xué)改革過程中使用的計(jì)算機(jī)編程項(xiàng)目的內(nèi)容設(shè)計(jì)和教學(xué)效果，本文詳細(xì)給出了三個(gè)計(jì)算機(jī)編程項(xiàng)目的求解方法和算例，以便開設(shè)類似課程的大專院校使用或在此基礎(chǔ)上對(duì)教學(xué)內(nèi)容做進(jìn)一步改進(jìn)。

1計(jì)算機(jī)編程項(xiàng)目1（Computer Project 1， CP1）：采用DAK方程進(jìn)行天然氣壓縮因子計(jì)算

采用DAK方程進(jìn)行天然氣壓縮因子計(jì)算的基本原理請(qǐng)參見文獻(xiàn)[1][2]?，F(xiàn)給出算例如下：

算例1：設(shè)氣體溫度為103.73F（313K），摩爾分?jǐn)?shù)如下表（馮叔初等，2006;例2-3， p.101）。求：壓力為1015.264psia（7.0MPa）時(shí)的氣體壓縮因子？

注：表中氣體摩爾分?jǐn)?shù)取自馮叔初等;臨界溫度、臨界壓力取自Ahmed;R = F + 459.67。

【解】：

第1步：計(jì)算擬臨界溫度（）、擬臨界壓力（）：;;注意在所有的計(jì)算過程中，溫度均需采用絕對(duì)溫標(biāo)：熱力學(xué)溫標(biāo)（K）或朗肯溫標(biāo)（R），本文使用朗肯溫標(biāo)。

第2步：如果天然氣組分中不含CO2和H2S，則直接跳至第3步。本例中，CO2和H2S組分的存在會(huì)使得天然氣壓縮因子偏離Standing-Katz壓縮因子圖（即DAK方程原始參考數(shù)據(jù)），Wichert和Aziz提出了一種針對(duì)酸氣擬臨界參數(shù)（和）的修正方法，使得經(jīng)過修正后的酸氣壓縮因子符合Standing-Katz壓縮因子圖。對(duì)于朗肯溫標(biāo)，Wichert和Aziz對(duì)擬臨界溫度的修正量可表示為：

（1）

此表達(dá)式與溫標(biāo)的選取有關(guān)（對(duì)于開氏溫標(biāo)：）。式（1）中代表CO2和H2S在天然氣中總的摩爾分?jǐn)?shù);B代表H2S在天然氣中的摩爾分?jǐn)?shù)。本例中，

第3步：計(jì)算擬對(duì)比溫度（）、擬對(duì)比壓力（）：;。將所得的與代入目標(biāo)函數(shù)（參文獻(xiàn)[1]中式（6）），形成以為唯一未知量的方程。采用牛頓法對(duì)其進(jìn)行數(shù)值求解（參文獻(xiàn)[1]中式（7）），求得（表2給出了詳細(xì)的迭代計(jì)算結(jié)果）：

表2：DAK狀態(tài)方程求解壓縮因子計(jì)算迭代結(jié)果

意味著初始迭代時(shí)需給出的初始值，由參文獻(xiàn)[1]中式（4）可得， （取壓縮因子;臨界壓縮因子）。注：表2在計(jì)算時(shí)，收斂標(biāo)準(zhǔn)設(shè)為

;

第4步：參文獻(xiàn)[1]中式（4），

。

在進(jìn)行CP1要求的天然氣壓縮因子圖（參文獻(xiàn)[1]，圖1、圖2）計(jì)算時(shí)，在給定的值下，取的取值范圍為[0，15]即可做出不同下的天然氣壓縮因子曲線;此外，將換算為即可得到壓縮因子與的關(guān)系曲線。

2計(jì)算機(jī)編程項(xiàng)目2（Computer Project 2， CP2）：二級(jí)地面油氣分離系統(tǒng)閃蒸計(jì)算

算例2：CP2的主要教學(xué)任務(wù)和教學(xué)目標(biāo)參見文獻(xiàn)[1]。

【解】求解CP2時(shí)，需要對(duì)每個(gè)分離器分別實(shí)施閃蒸計(jì)算。在完成閃蒸計(jì)算程序時(shí)，需要采用模塊化的編程思想將任務(wù)分解：（1）編寫一個(gè)求解立方形狀態(tài)方程的子函數(shù);（2）編寫一個(gè)求解Rachford-Rice目標(biāo)函數(shù)的子函數(shù);（3）編寫一段代碼更新平衡常數(shù)，調(diào)用前兩個(gè)子函數(shù)完成相平衡計(jì)算。具體求解過程如下：

第1步：求解立方形狀態(tài)方程（編寫Zfactor子函數(shù)）。

流體系統(tǒng)的PVT關(guān)系由狀態(tài)方程表述。Peng-Robinson狀態(tài)方程（PR EOS）和Soave-Redlich-Kwong狀態(tài)方程（SRK EOS）是目前工業(yè)界使用最廣泛的立方形狀態(tài)方程，其表達(dá)式分別如下：

（2）

（3）

式（2、3）中，為壓力（psia），為溫度（R），為摩爾容積（ft3/lbmol），為通用氣體常數(shù)10.73 ;和分別為狀態(tài)方程相對(duì)于混合物的引力項(xiàng)和斥力項(xiàng)參數(shù)，可由純組分參數(shù)（）和相應(yīng)的混合規(guī)則計(jì)算得到：

（4）

（5）

（6）

（7）

式（4-7中），分別為狀態(tài)方程對(duì)于第種組分對(duì)應(yīng)的純組分物質(zhì)的引力項(xiàng)和斥力項(xiàng)參數(shù);代表混合物中組分的數(shù)量;和是純組分參數(shù)，對(duì)于，;對(duì)于，;為混合物體系中第種組分的摩爾分?jǐn)?shù);和為混合物體系中第種組分的臨界溫度和臨界壓力;為混合物中第種組分和第種組分間的二元相互作用參數(shù)（Binary interaction parameter，BIP）;為第種組分的對(duì)比溫度（）;為與純組分偏心系數(shù)（）有關(guān)的表達(dá)式，對(duì)于PR EOS和SRK EOS分別如式（8.1）和式（8.2）所示：

（8.1）

（8.2）

在石油和天然氣行業(yè)的應(yīng)用中，立方形狀態(tài)方程常常表示為壓縮因子（z）的形式，根據(jù)實(shí)際氣體狀態(tài)方程，式（2、3）可分別改寫為：

（9）

（10）

其中，

;? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?（11）

Coats進(jìn)一步指出式（9）和（10）可以寫成統(tǒng)一的形式：

（12）

式（12）中，對(duì)于，;對(duì)于。求解式（12）的過程就是求三次方程的根。在教學(xué)實(shí)踐中，理解根的選擇是個(gè)教學(xué)難點(diǎn)。三次方程必然有三個(gè)根，但對(duì)于給定的溫度、壓力、組分的體系，只對(duì)應(yīng)一個(gè)具有物理意義的壓縮因子，需要在求解的三次方程時(shí)進(jìn)行根的選擇（參圖1）。

圖1中，在溫度下，立方型狀態(tài)方程反映的P-v變化關(guān)系如圖中較粗的黑色曲線所示，該曲線在氣液兩相區(qū)內(nèi)分別與泡點(diǎn)線交于A點(diǎn)并與露點(diǎn)線交于E點(diǎn)。AE水平連線對(duì)應(yīng)的壓力（即為溫度下的飽和蒸汽壓。A、E點(diǎn)分別對(duì)應(yīng)“可能的”液相或氣相摩爾容積和，“可能的”意味著在求解時(shí)，我們通常并不預(yù)先知道所給定的混合物是以氣相還是液相存在。圖中C點(diǎn)位于曲線BCD段上，該段曲線上，有，說明在該段曲線上對(duì)于閉口系統(tǒng)容積增大壓力也增大，不符合實(shí)際情況;因此C點(diǎn)處對(duì)應(yīng)的摩爾容積僅僅是一個(gè)數(shù)學(xué)意義上的解，沒有實(shí)際的物理意義，在求解立方形狀態(tài)方程時(shí)應(yīng)當(dāng)舍去。

此時(shí)，在剩下的兩個(gè)具有物理意義的解和中，必須再舍棄一個(gè)，因?yàn)樵谟?jì)算過程中針對(duì)液相（）或氣相（）組分進(jìn)行計(jì)算最終只能有一個(gè)具有物理意義的實(shí)根。通常有如下兩種處理方式：

如果預(yù)先知道計(jì)算過程針對(duì)的是氣相組分（），此時(shí)式（4）和式（5）中：，圖1中和兩個(gè)根選擇較大的根;對(duì)應(yīng)式（12）的解也選擇較大的根，。如果預(yù)先知道計(jì)算過程針對(duì)液相組分（），此時(shí)式（4）和式（5）中：，對(duì)應(yīng)式（12）的解也選擇較小的根，即。

在兩個(gè)具有物理意義的根（和）中，直接比較二者得到的系統(tǒng)Gibbs自由能大小，因?yàn)橄到y(tǒng)總是傾向于處于最小Gibbs自由能的狀態(tài)。Danesh指出：

（13）

式（13）大于零時(shí)，說明使得系統(tǒng)Gibbs自由能最小，取舍;反之，取。

驗(yàn)證壓縮因子子程序計(jì)算的正確性。采用文獻(xiàn)[1]中表3提供的原油組分及物性參數(shù)數(shù)據(jù)，將自行編寫的壓縮因子計(jì)算子程序的計(jì)算結(jié)果與下表驗(yàn)證后再進(jìn)行下一步：

第2步：求解Rachford-Rice目標(biāo)函數(shù)（編寫VLE子函數(shù)）。

由氣液相物料守恒可以推導(dǎo)得出Rachford-Rice目標(biāo)函數(shù)具有如下形式：

（14）

式（14）中，為體系中第種組分在體系中總的摩爾分?jǐn)?shù);為體系達(dá)到平衡狀態(tài)時(shí)的氣化分率（氣相中物質(zhì)的量與體系中總物質(zhì)的量的比值）;為第種組分的平衡常數(shù)，可由Wilson經(jīng)驗(yàn)公式近似計(jì)算得到：

（15）

式中，，。由式（14-15）可見，在各組分平衡常數(shù)已知的情況下，式（14）中唯一的未知量便是氣化分率，可采用數(shù)值解法（如：Newton-Raphson法）求解方程（14）得到氣化分率：

（16）

當(dāng)時(shí)停止迭代（為預(yù)先擬定的收斂標(biāo)準(zhǔn)）。在迭代計(jì)算過程中，需要給賦初值使得迭代計(jì)算開始（通常取），氣化分率的計(jì)算值雖然可能超出實(shí)際取值區(qū)間[0，1]，只要其最終的計(jì)算結(jié)果在區(qū)間[0，1]之間即可，為了在氣化分率超出實(shí)際取值區(qū)間時(shí)迭代計(jì)算仍然能夠繼續(xù)進(jìn)行，在每一步迭代后都進(jìn)行如下判斷：

（17）

式（17）中，和分別為組分平衡常數(shù)的最大值和最小值。當(dāng)?shù)諗亢?，氣液相的平衡組分可由式（18-19）求出：

（18）

（19）

此時(shí)，我們得到了平衡氣相組分（）以及平衡液相組分（）。通過采用氣相或液相組分替代式（4-5）中的（），我們可以得到氣相或液相的壓縮因子?，F(xiàn)采用文獻(xiàn)[1]中表3提供的原油組分及物性參數(shù)數(shù)據(jù)，讀者可將自行編寫的程序運(yùn)行結(jié)果與下表對(duì)照，驗(yàn)證求解Rachford-Rice目標(biāo)函數(shù)的正確性;學(xué)生經(jīng)與此表4內(nèi)的提供的數(shù)據(jù)驗(yàn)證正確后，再進(jìn)行下一步計(jì)算。

第3步：獲得準(zhǔn)確的平衡常數(shù)（編寫MAIN程序）。

命名MAIN主程序，實(shí)施平衡常數(shù)更新。第2步計(jì)算的前提是平衡常數(shù)已知（由Wilson公式計(jì)算）。但是在實(shí)際的工程應(yīng)用中，在特定的壓力、溫度和組分條件下，平衡常數(shù)的值往往并不能在一開始計(jì)算時(shí)就準(zhǔn)確獲得，Wilson經(jīng)驗(yàn)公式只是較粗略的估算。準(zhǔn)確的平衡常數(shù)需要根據(jù)相平衡理論反復(fù)迭代更新后才能得到。

當(dāng)體系中氣液兩相處于平衡狀態(tài)時(shí)，相間各組分逸度相等，可表示為：

（20）

式（20）表示含有種組分的混合物體系處于氣液兩相平衡時(shí)，任意第種組分在氣相和液相中的逸度相等（上標(biāo)L表示液相;V表示氣相）。通過逸度的定義（）以及平衡常數(shù)的定義，可得：

（21）

當(dāng)且僅當(dāng)時(shí)，式（21）表示體系的平衡常數(shù)，因此可用連續(xù)迭代法（式22）對(duì)計(jì)算中求得的平衡常數(shù)進(jìn)行更新：

（22）

在迭代開始時(shí)，需要給平衡常數(shù)賦初值，該初值往往由Wilson公式（式15）估算得到。在迭代過程中，當(dāng)滿足如下收斂標(biāo)準(zhǔn)時(shí)（式23）停止迭代，認(rèn)為平衡常數(shù)已不再變化。

（23）

由式（22）可以看出，平衡常數(shù)由逸度系數(shù)計(jì)算更新得到，而氣液相的逸度系數(shù)（）經(jīng)過嚴(yán)格的熱力學(xué)推導(dǎo)可由狀態(tài)方程得到（式24）：

（24）

式中，

（24.1）

（24.2）

（24.3）

（24.4）

（24.5）

注意在應(yīng)用式（24，24.3，24.5）時(shí)，要注意對(duì)于氣、液相來說，需分別帶入氣相組分（）或液相組分（）進(jìn)行計(jì)算。圖2給出了CP2閃蒸計(jì)算主程序MAIN的算法流程圖。

表5給出了在一定的溫度和壓力條件下，采用文獻(xiàn)[1]中表3提供的原油組分及物性參數(shù)數(shù)據(jù)、依照?qǐng)D2流程計(jì)算獲得的計(jì)算結(jié)果，如果學(xué)生自行編程計(jì)算的結(jié)果與下表各項(xiàng)數(shù)據(jù)吻合，則說明閃蒸計(jì)算流程及程序運(yùn)行均正確無誤。

第4步：將步驟1至3中獲得的閃蒸計(jì)算程序應(yīng)用至2級(jí)分離計(jì)算（參文獻(xiàn)[1]CP2主要內(nèi)容）中，計(jì)算儲(chǔ)罐中原油收率、產(chǎn)出液氣油比（GOR）、產(chǎn)出原油API重度等數(shù)據(jù)。

（1）儲(chǔ)罐中原油收率：

假設(shè)1磅摩爾（1lbmol = 453.5924mol）井流進(jìn)入高壓分離器，由計(jì)算知經(jīng)高壓分離器后得到0.5882磅摩爾液相原油;該部分原油再進(jìn)入中壓分離器，閃蒸分離后得到0.4931磅摩爾原油;該部分原油進(jìn)入儲(chǔ)罐后，在標(biāo)況下閃蒸分離后最終得到0.4821磅摩爾原油。由此可見，1磅摩爾井流最終在儲(chǔ)罐中得到0.4821磅摩爾原油，則儲(chǔ)罐中原油收率為48.2%。各級(jí)分離器閃蒸分離的計(jì)算結(jié)果如圖3所示。

（2）產(chǎn)出液氣油比（GOR）：

氣油比可以定義為產(chǎn)出氣的總體積（SCF， standard cubic feet）與儲(chǔ)罐中得到的原油體積（STB， stock tank barrel）之比。即：

（25）

式（25）中表示氣體在標(biāo)況下的體積;為原油在儲(chǔ)罐中的容積，由實(shí)際氣體狀態(tài)方程可得：

（26）

其中為原油在儲(chǔ)罐中的摩爾數(shù)（lb-mol）;為原油在儲(chǔ)罐條件下的壓縮因子;和分別為儲(chǔ)罐的溫度（R）和壓力（psia）。式（26）計(jì)算出的原油體積單位為SCF，轉(zhuǎn)化為STB需要除以5.615（1STB = 5.615SCF）。此外，假設(shè)產(chǎn)出氣在標(biāo)況（14.7psia， 60F）下可用理想氣體近似，1磅摩爾的氣體具有379.4SCF的體積，得：

（27）

式（27）中，代表經(jīng)各級(jí)分離器分離出氣體得總摩爾數(shù)。將式（26-27）帶入式（25）可得：

（28）

結(jié)合圖3中計(jì)算結(jié)果，，，帶入式（28），得到（SCF/STB）。

（3）儲(chǔ)罐原油API重度：

從實(shí)際氣體狀態(tài)方程，可得：

（29）

（30）

式中，為油水比重;為儲(chǔ)罐中原油的平均分子量;水密度（）取為62.4lbm/SCF，液相壓縮因子由閃蒸計(jì)算得到，則本例得到的結(jié)果為，為輕質(zhì)原油。

3計(jì)算機(jī)編程項(xiàng)目3（Computer Project 3， CP3）：天然氣管網(wǎng)管段流量及節(jié)點(diǎn)壓力計(jì)算

本例中管網(wǎng)系統(tǒng)求解時(shí)，要求學(xué)生分別使用兩種解法：（1）Q-formulation：以管段流量為未知數(shù)求解;（2）P-formulation：以節(jié)點(diǎn)壓力為未知數(shù)求解。如圖4所示（假設(shè)節(jié)點(diǎn)壓力已知）：

圖4：管網(wǎng)樣例（N：節(jié)點(diǎn);B：管段;L：閉環(huán);S：氣源;D：用氣）

對(duì)于任意給定的天然氣管網(wǎng)穩(wěn)態(tài)流動(dòng)過程，（1）其每個(gè)節(jié)點(diǎn)處流量守恒（流入節(jié)點(diǎn)流量等于流出節(jié)點(diǎn)流量）;（2）如氣流形成封閉環(huán)路，規(guī)定環(huán)路方向后（順時(shí)針或逆時(shí)針）則環(huán)路中各管段的壓力差之和為零。數(shù)學(xué)表達(dá)式如下：

（對(duì)所有節(jié)點(diǎn)） ? ? ? ? ? ? ? ? ? ? ? ?（31）

（對(duì)管流環(huán)路） ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? （32）

Q-formulation法：圖4中，當(dāng)采用管段流量（）為未知數(shù)求解管網(wǎng)問題，將管流方程寫成如下形式：

（33）

對(duì)于管網(wǎng)中每個(gè)節(jié)點(diǎn)依照式（31）可列方程如下：

節(jié)點(diǎn)：? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? （34.1）

節(jié)點(diǎn)：? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?（34.2）

節(jié)點(diǎn)：? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?（34.3）

節(jié)點(diǎn)：? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? （34.4）

考慮到，式（34.4）其實(shí)是式（34.1至34.3）的線性組合（式（34.4）=式（34.1）+式（34.2）+式（34.3）），說明式（34.1至34.4）中，只有3個(gè)獨(dú)立方程，而未知量卻有5個(gè)，因此仍需要兩個(gè)獨(dú)立方程求解方程組。觀察圖4中的管路形成了2個(gè)閉環(huán)（假設(shè)氣流在閉環(huán)內(nèi)以順時(shí)針方向流動(dòng)）：L1和L2;對(duì)每個(gè)閉環(huán)依照式（32）可列方程如下：

閉環(huán)：? ? ? ? ? ? ? ? ? ?（35）

閉環(huán)：? ? ? ? ? ? ? ? ? （36）

式（35、36）中采用式（33）將各管段的壓降轉(zhuǎn)換為管段阻力與流量的表達(dá)式;當(dāng)假設(shè)的管段流動(dòng)方向與閉環(huán)的流動(dòng)方向相反時(shí)，表示壓降需添加負(fù)號(hào)。聯(lián)立式（34.1至34.3）和（35-36）得到如下方程組（未知量）：

（37.1）

（37.2）

（37.3）

（37.4）

（37.5）

若記，由多變量牛頓法可得：

（38）

式（38）中，為Jacobian矩陣：

（39）

在用式（38）進(jìn)行計(jì)算時(shí)，需要提供的初始值。各管段流量的初始值一般根據(jù)節(jié)點(diǎn)流量和管段阻力加以估算，例如：圖4中計(jì)算初始值可估算如下：;;;;。在進(jìn)行迭代時(shí)，可將循環(huán)收斂標(biāo)準(zhǔn)定為當(dāng)向量（）中的最大元素小于。管段流量算出后，結(jié)合已知管網(wǎng)中某節(jié)點(diǎn)的壓力，可由式（33）解出官網(wǎng)中各節(jié)點(diǎn)的壓力。

P-formulation法：直接以管道中的節(jié)點(diǎn)壓力作為未知量進(jìn)行求解。將管流方程寫成如下形式：

（40）

由式（31）列出任意三個(gè)節(jié)點(diǎn)的方程如下：

節(jié)點(diǎn)

（41.1）

節(jié)點(diǎn)：? （41.2）

節(jié)點(diǎn)：

（41.3）

由于節(jié)點(diǎn)壓力已知，方程組（41）中含有三個(gè)方程、三個(gè)未知量（，，），可直接用牛頓法求解。該法看似較Q-formulation法簡單，但對(duì)于氣體管網(wǎng)，方程的非線性很強(qiáng)，且對(duì)初值的選取非常敏感，直接求解可能會(huì)遇到不能收斂的問題。對(duì)此法的改進(jìn)本文不做闡述，感興趣的讀者可參閱文獻(xiàn)[15]。

在使用通用性管流方程（Generalized）時(shí)，需要采用Colebrook公式求解Moody管道摩擦系數(shù)：

（42）

式（42）中，為管道絕對(duì)粗糙度，in;為管徑，ft;為雷諾數(shù)：

（43）

式（43）中，為氣體比重;為流量，MMSCFD;為氣體粘度，cp。當(dāng)采用牛頓法求解管道摩擦系數(shù)時(shí)，初始值可以令式（42）中的為無限大計(jì)算得到。結(jié)合上述求解方法，針對(duì)CP3給出的幾個(gè)穩(wěn)態(tài)工況下的氣體管網(wǎng)案例，現(xiàn)給出節(jié)點(diǎn)壓力及管段流量參考解。

算例3：管網(wǎng)系統(tǒng)及假設(shè)的管段流動(dòng)方向如圖5所示（管道及氣體參數(shù)見文獻(xiàn)[1]）

【解】：

算例4：管網(wǎng)系統(tǒng)及假設(shè)的管段流動(dòng)方向如圖6所示（管道及氣體參數(shù)見文獻(xiàn)[1]）

【解】：

算例5：管網(wǎng)系統(tǒng)及假設(shè)的管段流動(dòng)方向如圖7所示（管道及氣體參數(shù)見文獻(xiàn)[1]）

【解】：

算例6：管網(wǎng)系統(tǒng)及假設(shè)的管段流動(dòng)方向如圖8所示（管道及氣體參數(shù)見文獻(xiàn)[1]）

【解】：

4總結(jié)

本文給出了我校海洋油氣工程專業(yè)的本科生《海洋油氣集輸》挑戰(zhàn)性課程計(jì)算機(jī)編程項(xiàng)目的詳細(xì)求解方法和參照答案。該套編程項(xiàng)目對(duì)于中國石油大學(xué)（北京）大三的同學(xué)獨(dú)立完成存在一定困難，但在教學(xué)中取得了大多數(shù)學(xué)生的理解和認(rèn)可。本文所述的求解方法能夠幫助學(xué)生更好的完成相關(guān)教學(xué)環(huán)節(jié)，在獨(dú)立完成這些教學(xué)環(huán)節(jié)之后，學(xué)生不但能夠感受到自己在計(jì)算、分析及應(yīng)用能力上的進(jìn)步，同時(shí)也對(duì)所學(xué)的內(nèi)容理解更加透徹。

基金項(xiàng)目：中國石油大學(xué)（北京）教育教學(xué)改革項(xiàng)目“海洋油氣集輸全英文核心課程建設(shè)”。

參考文獻(xiàn)

[1] 趙海寧.石油與天然氣工程挑戰(zhàn)性課程教學(xué)實(shí)踐——海洋油氣集輸計(jì)算機(jī)編程項(xiàng)目設(shè)計(jì)[J].化工高等教育，2019，36（06）：40-50.

[2] Dranchuk，P.M.and Abou-Kassem，J.H.Calculation of Z factors for natural gas using equations of state[J]. Journal of Canadian Petroleum Technology，1975，14（03）：34-36.

[3] 馮叔初，郭揆常等.油氣集輸與礦場加工（第二版）[M].中國石油大學(xué)出版社，2006.

[4] Ahmed，T.Equations of state and PVT analysis applications for improved reservoir modeling， 2nd edition[M].Gulf Professional Publishing，2016.

[5] Standing，M.B.&Katz，D.L.Density of natural gases[J].Pet. Trans.AIME，1942（146）：140-149.

[6] Wichert， E.&Aziz，K.Calculation of Z's for sour gases[J].Hydrocarbon Processing，1972，51（05）：119-122.

[7] Rachford Jr. H.H.&Rice，J.D.Procedure for use of electronic digital computers in calculating flash vaporization hydrocarbon equilibrium[J].Pet. Trans. AIME，1952（195）：327-328.

[8] Peng， D.-Y.&Robinson，D.B.A new two-constant equation of state[J].Industrial Eng. Chem. Fundam，1976，15（01）：59-64.

[9] Soave，G.Equilibrium constants from a modified Redlich-Kwong equation of state[J].Chemical Engineering Science，1972（27）：1197-1203.

[10] Coats，K.H.Simulation of gas condensate reservoir performance[J].Journal of Petroleum Technology，1985，37（10）：1870-1886.

[11] Danesh，A.PVT and phase behaviour of petroleum reservoir fluids[M].Elsevier，1998.

[12] Wilson，G.A modified Redlich-Kwong Equation of state[A].Application to general physical data calculations [C].American Institute of Chemical Engineers National Meeting，1968.

[13] Whitson，C.H.&Michelsen，M.L.The negative flash[J].Fluid Phase Equilibria，1989（53）：51-71.

[14] Essenfeld，D.Network modeling and tracking of liquid preferential routes in natural gas pipeline systems [D].Master Thesis，The Pennsylvania State University，2010.

[15] Ayala， L.F.&Leong，C.Y.A robust linear-pressure analog for the analysis of natural gas transportation networks[J].Journal of Natural Gas Science and Engineering，2013（14）：174-184.

[16] Colebrook，C.F.Turbulent flow in pipes， with particular reference to the transition region between smooth and rough pipe laws[J].Journal of the Institution of Civil Engineers，1939，11（04）：133-156.

[17] Moody，L.F.Friction factors for pipe flow[J].Trans.ASME，1944（66）：671-684.