石 巖
(吉林建筑大學(xué)市政與環(huán)境工程學(xué)院,長春 130118)
?
井筒空間內(nèi)以CO2為載熱流體的多相流流動算法分析*
石巖
(吉林建筑大學(xué)市政與環(huán)境工程學(xué)院,長春130118)
摘要:在地?zé)衢_采井井筒空間內(nèi)流體的流速計算通常采用動量守恒方程,而對于多相流系統(tǒng)直接求解動量方程比較困難,同時也很難將井筒流和儲層流模型耦合在一起.為了解決這一問題,通過介紹引入漂移流模型DFM(Drift Flux Model)的方法,將多相流系統(tǒng)不同相的流速建立關(guān)系,從而將動量守恒方程中多個未知數(shù)轉(zhuǎn)為一個未知數(shù),實現(xiàn)多相流動量守恒方程各參數(shù)的求解.
關(guān)鍵詞:漂移流模型;動量守恒;多相流
深部地?zé)豳Y源開發(fā)過程中,井筒(水平、垂直或傾斜)是地表和地下熱儲層發(fā)生物質(zhì)和能量交換的主要通道.由于井筒流動空間大、表面光滑等特征,流體在其中的流動與在多孔介質(zhì)中截然不同,不再滿足達西定律,而是滿足更為普遍的動量守恒定律,而如何實現(xiàn)動量守恒方程的求解是計算井筒空間內(nèi)流體流速的關(guān)鍵問題.
1質(zhì)量和能量守恒方程
首先假設(shè)井筒具有圓形的一維線性特征,井內(nèi)流體的流態(tài)為紊流,流體的能量運移及轉(zhuǎn)換除了考慮傳導(dǎo)和對流以外,增加考慮了動能和勢能的影響和相互轉(zhuǎn)化,同時考慮井筒內(nèi)流體與周圍介質(zhì)的熱量交換過程.井筒中的流體的流速計算采用的是動量守恒方程,在多相流流動系統(tǒng)中,直接求解動量方程比較困難,同時也很難和儲層的模型耦合在一起.為了解決這個困難,引入漂移流模型DFM(Drift Flux Model)[1-2],將多相系統(tǒng)中不同相的流速建立關(guān)系,進而把動量守恒方程中的多個未知數(shù)轉(zhuǎn)為一個未知數(shù),實現(xiàn)動量守恒方程的求解.在完成動量守恒方程求解后,再利用DFM獲得每個相的流速,統(tǒng)一求解儲層和井筒的多相流動、熱傳導(dǎo)及對流方程.
根據(jù)質(zhì)量和能量守恒原理,廣義的井筒質(zhì)量和能量守恒方程如下:
(1)
式中,κ為組份指數(shù):κ=1表示水,κ=2表示CO2,κ=3表示能量(包含動能和內(nèi)能);Mκ為組份κ的累積項;qκ為質(zhì)量或能量的源(匯)項;Fκ為平流過程井筒內(nèi)質(zhì)量和能量的運移項.
單相或兩相系統(tǒng)中對于質(zhì)量組份(水和CO2)公式(1)的累積相Mκ計算公式如下:
(2)
式中,Xβκ為β相κ組份的質(zhì)量分?jǐn)?shù)(β=G表示氣相,β=L表示液相);ρβ為β相的密度;Sβ為β相的局部飽和度,計算公式如下:
(3)
式中,A為井筒的橫截面積,m2;AG和AL分別為被氣相和液相占據(jù)的給定高度(或井的高程)的橫斷面積,m2.
累積相的能量方程為:
(4)
流體沿著井筒的運移是通過對流、擴散和彌散等方式實現(xiàn).κ組份的對流通量被定義為:
(5)
式中,uβ為井筒內(nèi)的平均速度矢量,m/s;A為井筒的橫截面積,m2;z為沿井筒的坐標(biāo)(可以是垂直的,也可以是傾斜或水平的),m.
井筒內(nèi)流體能量的運移包括對流換熱、動能、勢能以及流體與側(cè)邊界的熱損失.一維能量運移項可以如下表示:
(6)
式中,hβ為流體β相的比焓,J/kg;g為重力加速度,Pa/m;θ為井筒的傾斜角度,°;q″為單位長度井筒的得(失)熱量(當(dāng)井筒及圍巖構(gòu)造不確定時是可選項),W;ρ為氣-液混合密度,kg/m3;λ為井筒截面積上的平均導(dǎo)熱系數(shù),W/(m·℃);T為溫度,℃.
為簡化計算過程,上述公式省略了井筒與周圍構(gòu)造的質(zhì)量或能量的交換.
井筒與圍巖之間的熱量交換通過兩種方式計算:
(1) 如果存在圍巖網(wǎng)格,可以通過TOUGH2中的流體通過孔隙介質(zhì)的導(dǎo)熱及對流公式計算,如將網(wǎng)格節(jié)點到井筒側(cè)面接口的距離設(shè)置為零時此項功能即取消;
(7)
其中,Qi3既包括井筒側(cè)壁面的傳熱量,也包括勢能的轉(zhuǎn)化(摩擦損失被轉(zhuǎn)化成熱,因此不會影響整體的能量平衡).式中,Awi為井筒與圍巖間的側(cè)面積,m2;Kwi為井筒與圍巖間的傳熱系數(shù),W/m2·K;Ti表示井筒第i個網(wǎng)格的溫度,℃;T∞(z)表示周圍環(huán)境的溫度,℃;r表示井筒的半徑,m;f(t)表示Ramey井流熱損失.
(8)
其中,α是圍巖的熱擴散系數(shù).
2漂移流模型(DFM)的動量守恒方程
孔隙介質(zhì)中的流動其流量或流速可以通過達西定律的壓力和重力梯度進行簡化計算,但井筒內(nèi)的流速是通過求解動量守恒方程得到的.直接求解與儲層模擬器耦合的兩相流體動量方程比較困難,因此,我們通過調(diào)用DFM模型描述井筒內(nèi)單相和多相流動,從而獲得對流傳輸項(質(zhì)量或能量流量Fβ與流速uβ).
漂移流模型最早由Zuber,F(xiàn)indlay及Wallis提出[1-2],其基本理論是將氣相速度uG、液相速度uL用體積流量j、系數(shù)C0及漂移速度ud來表達.
(9)
式中,C0為形狀參數(shù),影響井筒橫斷面的局部氣體飽和度和速度.由定義,體積流量j為體積加權(quán)的平均速度:
(10)
因此,液相流速uL可以定義為:
(11)
漂移流模型公式(9)~(11),井筒兩相流的動量方程可以簡化成混合流速um和漂移流速ud的單相方程.漂移流模型中所涉及的所有變量除明確指出,否則均被認為是井筒截面的面積平均值或是常數(shù).
忽略井筒軸向應(yīng)力時,井筒多相流動量守恒方程如下[3]:
2016年企業(yè)研發(fā)投入的參數(shù)β2=7.41,即在其他條件不變的情況下,研發(fā)費用每增加1萬元,主營業(yè)務(wù)利潤增加7.41萬元。P值為0.058,在10%的顯著性水平下,表明2016年的企業(yè)研發(fā)投入對于2017年的主營業(yè)務(wù)利潤有顯著影響。
(12)
混合密度ρm和um定義如下:
(13)
(14)
受形狀影響的平均密度ρm*定義如下:
(15)
因此,通過采用DFM方法,使復(fù)雜的兩相動量方程簡化成兩個步驟:
(1) 通過求解簡化的動量方程(12)獲得混合流速um,以及通過經(jīng)驗關(guān)系獲得漂移流速ud;
(2) 再通過下式計算氣體流速和液體流速:
(16)
不同界面上的兩相流相互作用,形成不同的流動機制,因此如何準(zhǔn)確計算漂移流速ud及形狀參數(shù)C0是DFM的難點.不同流動機制C0和ud關(guān)系函數(shù)是不同的.依據(jù)Shi等提出的形狀參數(shù)和漂移流速的功能表,C0與ud是通過大量廣泛的管道流實驗獲得的[4].這些實驗是Oddie等在多種條件下針對單相、兩相和三相流體流動所做的,并被連續(xù)應(yīng)用到所有流動機制[5].TOUGH2Bore中采用的漂移流速數(shù)學(xué)公式是由Shi等提出[4].首先,通過氣體飽和度和其他流體性質(zhì)計算漂移流速,公式如下:
(17)
式中,m(θ)被描述為m(θ)=m0(cosθ)n1(1+sinθ)n2.其中m0,n1和n2是安裝參數(shù);Ku是庫塔杰拉茲(Kutateladze)數(shù),為邦德數(shù)NB的函數(shù)(井徑的無量綱數(shù)的平方):
(18)
式中,Cw是壁面摩擦因素(在代碼被取為固定數(shù)0.008),邦德數(shù)NB定義如下:
(19)
式中,d是井筒直徑,參照圖1相同文獻的實驗數(shù)據(jù),Cku在TOUGH2Bore中取為142.
其中,特征流速uc是液體沸騰上升流速的測量值,計算公式如下:
(20)
式中,σGL是氣相和液相之間的表面張力.
公式(17)中K的功能是使漂移流速在沸騰上升階段和浸水層流階段形成一個平穩(wěn)過渡.不同于Shi等線性插入方法,我們使用下面的功能:
(21)
這里a1和a2是Shi等[4]提出的氣體飽和度的兩個轉(zhuǎn)換點.當(dāng)SG≥a2和SG≤a1時,K不依賴氣體飽和度;當(dāng)a1≤SG≤a2時,K的功能通過dK/dSG=0搭建起來,使導(dǎo)數(shù)在整個區(qū)間均連續(xù).對于擬合參數(shù)m0,n1,n2,a1和a2的賦值來自Shi等中水/氣實例,其值取決于Cmax(一個用戶指定的最大形狀參數(shù)介于1.0和1.5之間),如表1所示.
表1 TOUGH2Bore 中DFM的經(jīng)驗參數(shù)[4]
形狀參數(shù)Cmax通過下式進行計算[4]:
(22)
這里η是反映形狀參數(shù)對流動狀態(tài)影響的參數(shù),計算公式如下:
(23)
式中,B是高于C0低于Cmax的極限參數(shù).β的計算公式如下:
(24)
與文獻[4]中所述略有不同的是,我們使用的極限參數(shù)B是由Cmax決定:
(25)
式中提供當(dāng)Cmax=1.2時,B=0.6,這與文獻[4]中a1(=0.06)和a2(=0.12)的計算值是一致的.B的變化區(qū)間由Cmax=1.0時的0.933 3到Cmax=1.5時的0.266 6.需注意的是:如果Cmax=1.0,C0將不受氣體飽和度或流速的限制,形狀影響消失.
形狀平的或多或少對氣相流速的影響可以通過調(diào)整公式(24)的FV(默認值=1),因此usgf可以通過下式計算:
(26)
3結(jié)語
本文通過介紹在井筒空間內(nèi)引入漂移流模型DFM(Drift Flux Model)的方法,將多相流動量守恒方程中多個未知數(shù)轉(zhuǎn)化為一個未知數(shù),從而獲得不同相流速的求解方法,實現(xiàn)了井筒空間內(nèi)以CO2為載熱流體的多相流動量守恒方程參數(shù)求解.
參考文獻
[1] Zuber N,Findlay JA.Average volumetric concentration in two-phase flow systems[J].J Heat Transfer ASME,1965,87(4):453-468.
[2] Wallis.One-Dimension Two Phase Flow[M].New York,1969:10-54.
[3] Brennen,CE.Fundamentals of multiphase flows[M].Cambridge:Cambridge University Press,2009:34-35.
[4] Shi H,Holmes JA,Durlofsky LJ,Aziz K,Diaz LR,Alkaya B,Oddie G.Drift-flux modeling of two-phase flow in wellbores[J].Soc Pet Eng J,2005,10(1):24-33.
[5] Oddie G,Shi H,Durofsky LJ,Aziz K,Pfeffer B,Holmes JA.Experimental study of two and three phase flows in large diameter inclined pipes[J].International Journal of Multiphase Flow,2003,29(4):527-558.
Algorithm Analysis of the Flow on Multi-phase Fluid in
Wellbore with CO2as the Heat Transfer Fluid
SHI Yan
(SchoolofMunicipalandEnvironmentalEngineering,JilinJianzhuUniversity,Changchun,China130118)
Abstract:In a multiphase flow system,the flowrate of wellbore fluid is calculated according to the momentum conservation equation.It is difficult to directly solve the momentum equation,and it is also difficult to couple with the reservoir model together.In order to solve this problem,the drift flow model (DFM)is conductied to build a relationship of the flow rate of the different phases in the multiphase system,and then many of the in the momentuconservation equation is converted to one unknown parameter,so it can be to realize the solution of the momentum conservation equation. After completion of the momentum conservation equations,the each phase velocity is can be get by using DFM.
Keywords:drift flux model;conservation of momentum;multi-phase fluid
*基金項目:吉林省科技廳科技計劃項目(20130522070JH);吉林省教育廳“十二五”科學(xué)技術(shù)研究項目(吉教科合字[2014]第239號).
作者簡介:石巖(1977~),女, 吉林省長春市人,副教授,博士.
收稿日期:2014-10-10.
中圖分類號:TB 302.2
文獻標(biāo)志碼:A
文章編號:2095-8919(2015)03-0031-05