基于耦合分析的氣井井筒溫度壓力剖面計算方法

李亞輝 ,王 樂 ,張丹丹

(1.榆林學院 化學與化工學院,陜西 榆林 719000;2.中國石油西南油氣田公司工程技術研究院,四川 成都 610000)

溫度和壓力對氣井的生產和管理至關重要,常通過井底或井口數(shù)據采用理論計算獲取,目前常見的計算模型大多都存在一些缺點和適用范圍[1-13]。氣體和井筒周圍環(huán)境之間存在溫差,導致井深減小氣體溫度逐漸降低;井深減小壓力降低,氣體焦耳湯姆遜效應顯著,使得氣體溫度降低;高速流體與井筒間的摩擦使流體溫度升高;氣體物性參數(shù)都是關于溫度和壓力的函數(shù);因此新模型需考慮熱傳遞、焦耳湯姆遜效應、摩擦、動能和重力勢能的影響,且必須與壓力進行耦合計算。

1 理論分析

1.1 模型假設

在進行理論模型推導前,首先作如下假設[1,14]:

(1)井筒中的傳熱為穩(wěn)態(tài)傳熱,地層中傳熱為非穩(wěn)態(tài)傳熱,傳熱過程中不考慮熱損失;

(2)熱量損失在井筒及地層中僅考慮徑向,不考慮沿井深方向的傳熱;

(3)井筒內任意截面上的溫度處處相等;

(4)地層溫度按線性變化,且地溫梯度已知。

1.2 模型推導

在井筒上取長度為dz的微元體(圖1),取井底為坐標原點,垂直向上為正,根據能量守恒定律得出:流入微元體的熱量QIN等于流出微元體的熱量QOT與微元體內損失的熱量ΔQ的總和,即

圖1 井筒微元體能量守恒示意圖

式中,QIN、QOT、ΔQ分別為單位時間內微元段入口、出口處和段內損失的熱量,J/s。

根據熱量與溫度的關系[1],結合式(1)可知單位時間內微元體內損失的熱量為

式中,Wm為流體質量流量,kg/s;Cpm為混合流體定壓比熱,J/(kg·℃);Tf(a)、Tf(b) 分別為微元段入口與出口處流體的溫度,℃。

微元段內損失的熱量ΔQ包含向地層傳熱Qr、焦湯效應引起的熱損失QJT、重力做功W和動能K改變引起的熱損失以及摩擦生熱Qf,即

根據牛頓冷卻定律可知,微元段向水泥環(huán)/地層界面?zhèn)鬟f的熱量為

式中,r為油管外徑,m;U為總傳熱系數(shù),J/(s·m2·℃);Th為微元段位置對應的水泥環(huán)/地層界面溫度,℃。

界面向地層徑向傳遞的熱量為

式中,ke為地層導熱系數(shù),J/(s·m·℃);Te 為微元段位置對應的地層溫度,℃。

不考慮傳熱過程中的熱損失,因此微元段向界面?zhèn)鬟f的熱量Qr1與界面向地層傳遞的熱量Qr2相等,即

將式(4)、(5)代入式(6)并化簡得

由于地層溫度為線性變化,即

式中,Ter為儲層中部地層溫度,℃;gT為地層溫度梯度,℃/m。將式(7)～(8)代入式(4)得

在井筒中沿流向壓力逐漸降低,氣體膨脹導致氣體溫度降低,即產生焦耳湯姆遜效應。微元段內氣體因焦耳湯姆遜效應而損失的熱量為

式中,Wg為氣體質量流量,kg/s;Cpg為氣體定壓比熱,J/(kg·℃);μJT為焦耳湯姆遜系數(shù),℃/Pa;p為壓力,Pa。

重力導致的熱量損失量為

式中,g為重力加速度(=9.81),m/s2;Wm為混合流體質量流量,kg/s。

動能增加導致熱量損失量[15]為

式中,v為氣體流速,m/s。

摩擦產生的熱量為[16]式中,λ為摩阻系數(shù),無因次;D為井筒直徑,m。

將式(9)～(13)代入式(3)得

將式(2)代入上式(14)左邊,然后兩端同時對z求導數(shù)并化簡得

式中,qsc為地面產量,m3/d;Z為氣體偏差因子,無因次;T為氣體平均溫度,K。

氣體沿井筒從微元段的a端流向b端,以a端的溫度和壓力求得的氣體的偏差系數(shù)Z作為微元段ab的平均偏差系數(shù),則偏差系數(shù)對于式(17) 而言即為常數(shù),則由式(17) 可得

將式(20)～(24) 代入式(15) 并化簡得

求解式(24)、(26)的微分方程即可得到井筒溫度與壓力分布的耦合預測模型。

式中,C1和C2為待定系數(shù),可根據溫度和壓力已知條件進行求解。

當井筒流體為氣液兩相流時,氣液混合物的物性參數(shù)采用Hagedorn-Brown 方法計算。

2 溫度、壓力耦合模型求解

2.1 熱力學參數(shù)的計算

BWRS 狀態(tài)方程是目前公認的最精確的經驗狀態(tài)方程,本文采用BWRS 狀態(tài)方程[17-18]求取焦耳湯姆遜系數(shù)μJ_T。由BWRS 狀態(tài)方程求出(?ρ/?T)p與ρ 后[19-20],即可求得焦耳湯姆遜系數(shù)。

式中,μJ_T為焦耳湯姆遜系數(shù),K/Pa;Cp為氣體定壓比熱容,J/(kg·℃);Tg為氣體絕對溫度,K。

定壓比熱容的計算常采用經驗擬合公式和狀態(tài)方程兩種方法進行計算。擬合經驗公式常采用下式計算[18]。

式中,M為氣體摩爾質量,kg/kmol。

天然氣的偏差系數(shù)采用基于BWRS 狀態(tài)方程所建立的新方法[21-22],該方法計算精度比D-A-K方法精度更高。

式中,A1～A16為模型通用參數(shù);Tr為對比溫度。

瞬態(tài)傳熱函數(shù)f(t)的精確求解過程非常復雜,因此多采用能夠滿足工程精度要求的近似模型進行計算,常采用如下模型[4]進行近似計算。

式中,t為時間,s;α為地層熱擴散系數(shù),m2/s。

總傳熱系數(shù)U 的計算也相對比較復雜,目前主要用的是Ramey[5]在1967 年所提出的模型進行計算。由于油管內流體與管壁之間的傳熱系數(shù)非常大,故忽略其熱阻;油管和套管的導熱系數(shù)比水泥環(huán)和地層的導熱系數(shù)大得多,且相對于地層和水泥環(huán)的尺寸,油管和套管的壁厚很小,故在參數(shù)不全的情況下,可以忽略油管和套管造成的溫度變化,則總傳熱系數(shù)可簡化為[1]

式中,U為總傳熱系數(shù),J/(s·m2·℃);hc、hr分別為環(huán)空對流、輻射傳熱系數(shù),J/(s·m2·℃);rh、rco、rto分別為井眼半徑、套管外半徑、油管外半徑,m;kcem為水泥環(huán)導熱系數(shù),J/(s·m·℃)。

2.2 耦合模型求解

假設已知井底條件,從井底向井口逐步求解井筒溫度和壓力分布。

(1)給定井底溫度、壓力,z=Hw、Tf=Tini、Tei=Tini、pf=pini,將整個井深劃分為n 段,每一段作為一個計算單元。

(2)以微元段下端的溫度和壓力數(shù)據,即Ta和pa計算在a點時氣體的v、μg、Z、μJT、Cp、ρ、λ等物性參數(shù)。

(3)用步驟(2)中計算的物性參數(shù)計算微元段上端b 的溫度Tb和壓力pb。

(4)以ab段的平均溫度(Ta+Tb)/2 和平均壓力(pa+pb)/2 計算ab段天然氣的平均物性參數(shù)。

(5)以(4)中計算的平均物性參數(shù)重新計算微元段b點的溫度Tb和壓力pb。

(6)比較步驟(3)和(5)所計算的溫度與壓力,如絕對誤差滿足要求,則轉步驟(2)計算下一微元段,否則用步驟(5)計算的結果更新b點溫度和壓力,再轉步驟(4)重新計算微元體ab段。重復步驟(4)～(6),直到滿足精度要求。

(7)依次類推,即可求出整個井筒的溫度和壓力分布。

3 模型驗證與分析

3.1 模型驗證

采用文獻[21]中A 氣田氣井數(shù)據對模型進行驗證,氣藏溫度為110 ℃、天然氣相對密度為0.598、油管直徑62 mm、管壁相對粗糙度0.004 5,測試井的基本參數(shù)如表1 所示。

表1 氣井基本參數(shù)

不同模型計算的井口壓力結果如下表2 所示。

表2 模型計算井口壓力與實測壓力對比

表2 表明,本文所提模型計算的井口壓力的精度比H-B 法更高,但由于缺少現(xiàn)場井口溫度實測數(shù)據,不能進行溫度誤差分析。根據文獻[23]中的觀點,當溫度和壓力是通過耦合計算所得,只要壓力計算結果準確,則認為溫度計算結果也比較可靠,因此本文所建立的耦合模型是可靠的,能滿足工程計算要求。

3.2 模型分析

對于無水氣井分析采用文獻[2]中某海上氣井X 的數(shù)據(表3)進行計算。

表3 某海上X 氣井基本參數(shù)表

井筒溫度分布計算結果如圖2 所示。圖2 表明,在海床以下不同模型計算的結果趨勢相同,由于考慮了氣體的焦湯效應使得溫度降低速率逐漸加快;而在海水中溫度變化趨勢存在一定差異,主要是在海水中環(huán)境溫度在逐漸升高,導致向海水傳遞的熱量在逐漸減少;最終減緩了溫度降低的速率。

圖2 X 氣井采用不同方法計算的生產40 天的井筒溫度剖面

圖3 表明整個溫度剖面在地層溫度越高的情況下,井筒溫度下降越慢,且隨生產時間變化較小;隨著生產時間的增加,井口溫度先逐漸增加,在生產30 天后,井口溫度已趨于穩(wěn)定,這表明地層熱傳遞過程已由瞬態(tài)傳熱逐漸趨于穩(wěn)態(tài)傳熱;現(xiàn)場上通常認為井筒溫度達到穩(wěn)定一般需要一個月以上,該溫度剖面的結果也驗證了這一點。

圖3 不同時間不同深度對應的溫度分布

產水氣井耦合分析采用如表4 所示的參數(shù)。

表4 耦合分析基本參數(shù)

3.2.1 耦合性及模型對比分析

在進行耦合性分析時,以耦合模型計算的井口溫度與井底溫度的平均值計算壓力剖面,以及用耦合模型計算的井口壓力與井底壓力的平均值計算溫度剖面,然后與耦合模型計算的溫度剖面進行對比,結果如圖4 所示。

圖4 耦合計算對溫度剖面的影響

圖4 表明耦合模型和平均壓力模型計算的結果之間存在一定的誤差,從井深約4/5 處開始誤差開始變得相對比較明顯,在井口誤差達到最大。主要是因為氣體的物性參數(shù)是溫度和壓力的函數(shù),不考慮壓力的變化,則計算的氣體物性參數(shù)存在較大的誤差,從而導致井筒溫度剖面存在差異,這種情況在低壓下更為明顯。

目前常用井筒瞬態(tài)溫度模型有Ramey[6]模型和Hasan&Kabir[4]模型。分別將本文耦合模型與Ramey及Hasan&Kabir 模型計算結果進行對比,結果如圖5所示。

圖5 不同模型在生產1 天和100 天的溫度剖面

圖5 表明隨著生產時間的增加,各模型計算的溫度剖面都在不斷增大。對于短期預測和長期預測,Hasan&Kabir 模型計算的值都偏小,主要是由于該型沒有考慮摩擦力做功對溫度剖面的影響;Ramey 模型未考慮焦耳湯姆遜效應的影響導致計算結果偏大,本模型所計算的井口溫度介于兩者之間,具有較高的精度。

3.2.2 產氣量敏感性分析

分別計算產氣量為5×104、25×104、50×104及60×104m3/d 在生產30 天時的溫度和壓力剖面,結果如圖6 所示。

圖6 不同產氣量對應的壓力、溫度剖面

圖6 表明隨著產氣量的增加,井筒溫度逐漸升高,主要是由于增大產氣量使得井筒流體攜帶的熱量增加而使溫度升高;隨著產氣量的增加,井筒壓力呈逐漸降低的趨勢,這主要是因為產氣量增加,從而導致摩擦阻力增大。

將井底壓力與井口壓力之差定義為舉升壓力,作不同產氣量與舉升壓力的關系曲線圖,如圖7 所示。

圖7 不同產氣量對應的舉升壓力

圖7 表明舉升壓力呈先下降后上升的趨勢,這主要是因為在產氣量較低的情況下,井筒內摩阻較小,使井筒壓降損失主要由流體重力損失造成;在較高產氣量下,井筒混合流體的密度逐漸降低,使得流體重力損失降低,但是產氣量增大使得流體與壁面之間的摩阻增大。即在低產氣量下,流體重力損失占主導作用;而產氣量較高的情況下,摩阻損失占主導作用。這也解釋了為什么圖6 中產氣量為5×104m3/d 和25×104m3/d 兩種情況下壓力剖面比較接近的原因。

3.2.3 產水量敏感性分析

分別計算產水量為1、3、6 和12 m3/d 生產30 天時的溫度和壓力剖面,結果如下所示。

圖8 表明隨著產水量的增加,井筒壓力剖面逐漸降低,主要是由于產水量增加使得井筒流體密度增加。圖9 表明上在產水量較小的情況下,隨著產水量的增加井筒溫度剖面逐漸降低,這主要是由于產水量增加使得井筒壓降損失增加,導致因氣體焦耳湯姆遜效應損失的熱量增多。

圖8 不同產水量對應的壓力、溫度剖面

圖9 焦-湯效應在不同產水量下對溫度剖面的影響

4 結論及認識

(1)綜合考慮焦耳湯姆遜等微熱效應的影響,通過建立的基于耦合計算的溫度分布模型對氣井的井口壓力進行了預測,并與其它模型進行對比,結果表明本文的耦合模型具有更高的預測精度,可用于提高節(jié)點系統(tǒng)分析結果精度及井筒中天然氣水合物預測研究,指導現(xiàn)場實際生產。

(2)氣井井筒溫度與壓力分布預測,必須以耦合分析為條件,并考慮焦耳湯姆遜效應及摩擦生熱的影響,預測結果精度才會進一步提高。

(3)氣井產氣量與產水量對井筒壓力剖面影響較大,且井筒流體舉升所需壓差呈先減小后增大的趨勢,合理的產量有助于減小井底生產壓差。