基于格子玻爾茲曼方法的致密砂巖驅替模擬

唐明明, 盧雙舫, 辛 盈, 檀洪坤, 劉雪萍, 時 賢, 張洪亮

(1.中國石油大學(華東)深層油氣重點實驗室,山東青島 266580; 2.中國石油國際勘探開發(fā)有限公司,北京100621;3.中國石油大學(華東)石油工程學院,山東青島266580; 4.中國石油大慶油田第九采油廠地質大隊,黑龍江大慶 163853)

隨著能源需求的日益攀升以及常規(guī)油氣資源的日漸枯竭,具有較大資源潛力的非常規(guī)油氣,特別是致密油氣,已經(jīng)開始重塑全球能源格局[1-3]。然而目前致密油藏衰竭式開發(fā)的平均采收率仍然較低(4%～10%)[4]。水驅作為油藏能量補充的重要手段,其在致密油藏中的應用研究受到人們的廣泛關注[5-6]。目前,致密砂巖巖心驅替物理試驗研究,可以分析水驅采收率及殘余油分布特征,但是巖心物理試驗成本高、周期長、計量困難,并且難以分析流體三維動態(tài)變化規(guī)律。通過油藏尺度水驅模擬研究可以得到致密油藏水驅產(chǎn)能特征,然而主要針對宏觀尺度分析,無法揭示微觀滲流機制[7-13]。筆者基于三維數(shù)字巖心建模方法和格子玻爾茲曼多相流方法[14-15],構建致密砂巖三維數(shù)字巖心模型,通過恒壓條件與恒速條件下的模擬,分析致密砂巖水驅的三維微觀過程,研究不同條件下的驅替效率和換油率。

1 模擬方法

格子玻爾茲曼方法(LBM)的核心是建立微觀尺度與宏觀尺度的橋梁,它不考慮單個粒子的運動而是將所有粒子的運動視為一個整體,粒子的整體運動特性由分布函數(shù)表示,粒子存在的空間用格子排列表示[16]。目前,對于三維格子排列模型主要有3種D3Q15、D3Q19和D3Q27,其中“3”代表格子排列模型的維度,“15”、“19”和“27”代表速度模型中的離散的速度分量。速度分量越多,數(shù)值模擬的準確性越高。考慮到模型孔喉結構復雜,本次模擬采用D3Q27模型(圖1)，以達到最高計算精度。

D3Q27的LMB的模型如圖1所示,其中數(shù)字標注的為離散速度矢量(共27個離散速度矢量)。第0個離散速度矢量位于格子單元中心,離散速度為0,其他的離散速度矢量分為以下幾種:①從模型中心指向面的中心點的離散速度矢量(圖1中1～6),其離散速度為1;②從模型中心指向各邊的中點的離散速度矢量(圖1中7～18),其離散速度為2;③從模型中心指向模型頂點的離散速度矢量(圖1中19～26),其離散速度為3。按照LBM理論,所有離散速度值均為無量綱形式。

在D3Q27模型中,分布函數(shù)滿足以下格子玻爾茲曼輸運方程[14-20]:

fα,k(x+ekΔt,t+Δt)-fα,k(x,t)=

(1)

圖1 格子玻爾茲曼D3Q27模型Fig.1 D3Q27 model of lattice Boltzmann method

方程式(1)可由BGK方程導出,BGK方程[14-20]為

(2)

其中弛豫時間λα與流體黏度μα存在以下相關關系:

(3)

(4)

(5)

(6)

在D3Q27模型的基礎上,采用Shan-Chen模型開展多相流模擬研究。Shan-Chen模型的基本思想是利用一個依賴于局部密度和相互作用強度的勢函數(shù)來表示不同相態(tài)分子間的相互作用[23-24],已有的Shan-Chen模型研究表明,在D3Q27模型中加入合適的勢函數(shù),不同相流體會自動分離[25-26]。Shan-Chen模型已成功地應用于多相流的廣泛研究中,特別是多孔介質中的多相流[25-26]。

Shan-Chen模型中,作用于流體相α的總力Fα,ff(x,t)[27]可以表示為

(7)

式中,ψα(x,t)為偽勢;α和β表示α相流體和β相流體;Gα,β為相α相與β相之間的相互作用強度;S為流體相的數(shù)目;B為離散速度矢量的數(shù)目。

基于相互作用力,對于多相流LBM提供了模擬方法。多相流的LBM模擬通過協(xié)同速度v(x,t)′實現(xiàn),其中

采用格子量綱,而被Fα,ff(x,t)改變的瞬時速度vα(x,t)[24]可以表示為

(8)

由于本文中只考慮兩相流,因此Gα,β在之后模擬中始終為常數(shù)G,通過改變相互作用強度矩陣中的Gα,β來改變流體相的狀態(tài)[25-26]。

接觸角也由固相和流體相之間的相互作用力控制,可將式(7)中的ψβ(x+ekΔt,t)改為π(x+ekΔt,t),由此流體相和固體表面的作用力[18-19]可以表示為

(9)

式中,Uα為流體相與固體相之間的相互作用強度,采用格子量綱;π(x+ekΔt,t)等價于1代表固體,等價于0代表流體相。因此作用在流體節(jié)點上的總作用力為Fα(x,t)=Fα,ff(x,t)+Fα,fw(x,t)。

2 試驗數(shù)據(jù)

本次研究樣品Y采自松遼盆地龍26致密油區(qū)塊W井高臺子油層。Y樣品直徑DY為2.5 cm, 軸長LY為10 cm,孔隙度φe為9.7%,滲透率Ke為0.95×10-3μm2。Y樣品采回后,在實驗室內分為2段,其中一段為用于進行CT掃描試驗的樣品A,另外一段為用于進行室內驅替試驗的樣品B,其中樣品B長度為9 cm。再將B樣品分割為B1和B2兩個樣品,其長度均為4.2 cm。B1樣品應用于恒壓驅試驗,B2樣品應用于恒速驅替試驗。其中物理試驗選用二氯甲烷和丙酮體積比為3∶1的混合有機溶劑作為洗油溶劑。模擬地層水為NaHCO3型,礦化度為5 497 mg/L。驅替使用的模擬油為淡黃色煤油,密度為0.8 g/cm3,黏度為1.0～2.0 mm2/s。每次試驗前首先利用洗油溶劑對致密砂巖樣品進行洗油,時間為72 h,之后進行24 h的高溫烘干,烘干溫度為110 ℃。烘干完成后,將樣品放入飽和儀內抽真空,飽和模擬地層水,飽和8 h。為了保證數(shù)值模擬試驗與物理試驗的可對比性,恒壓驅替條件下,物理試驗和數(shù)值試驗按照等壓力梯度原則設計試驗壓差,恒速驅替條件下,按照流速相同原則,設定樣品端口流量。

內部控制評價主體由內部評價主體和外部評價主體兩部分組成，兩個評價主體分別有不同的評價目標，因此就使得兩個主體之間必然會產(chǎn)生矛盾。兩個評價主體會在考慮自身利益的基礎上進行內部控制評價，這就使得在內部控制評價過程中兩個評價主體的效率有一定的差異。公司內部評價主體在內部控制評價時更注重公司的經(jīng)營戰(zhàn)略目標，這從很大程度上保證了公司戰(zhàn)略目標的順利實現(xiàn)；而外部評價主體關注的是報告的合法合規(guī)性，外部評價主體在保證自身獨立、權威的同時協(xié)助政府監(jiān)督機構對所評價公司的內部控制進行監(jiān)督審查，提高所評價公司的內部控制質量，從而促進會計穩(wěn)健性的提高。根據(jù)分析提出第一個假設：

圖2 三維數(shù)字巖心構建Fig.2 Three dimensional digital rock construction

數(shù)值模擬中使用的三維數(shù)字巖心采用A樣品,使用ZEISS CT 3D X射線顯微鏡進行樣品掃描,其掃描分辨率為1.0 μm,掃描視場為1 747 μm,掃描電壓為50 kV,電流為50 mA,曝光時間0.7 s,掃描總時長2 h。圖2(a)為樣品的一張CT掃描圖像局部(400 μm×400 μm),圖片分辨率為1 μm/像素。為了提高計算精度,對CT圖像數(shù)據(jù)進行了10倍加密,加密后像素大小為0.1加密。應用中值濾波處理,中值濾波處理使骨架與孔隙之間的過渡變得更加平滑自然,并且去除了礦物表面的孤立噪點。中值濾波利用MATLAB中medfilter2濾波函數(shù),采用最鄰近中值濾波模板參數(shù)。濾波之后通過設定分割閾值,對圖像組分進行分割,將灰度圖轉化為二值的黑白圖像,如圖2(b),其中黑色代表孔隙,白色代表骨架。經(jīng)過二值化分割后的圖像沿著切片縱向疊加并通過閾值分割提取孔隙,得到致密砂巖樣品的三維數(shù)字巖心。然后基于二值化模型,計算數(shù)字巖心的孔隙度和滲透率。然后將數(shù)字巖心得到的物性數(shù)據(jù)與試驗數(shù)據(jù)進行對比,不斷優(yōu)化分割閾值,最終優(yōu)化分割閾值為Th=51,數(shù)字巖心孔隙度和滲透率分別為9.9%和1.02×10-3μm2?；谏鲜鰯?shù)字巖心模型,利用格子玻爾茲模型,開展驅替模擬。格子玻爾茲曼模型中,固體邊界采用全反彈模型。模擬物理空間步長為Δxp=1×10-7m,物理時間步長為Δtp=1×10-9s,格子特征速度為uLB=1×10-5,物理特征速度為up=1×10-3m/s。

為了縮減模擬時間,從重構的三維數(shù)字巖心中截取尺寸為210 μm×80 μm×40 μm的模型A(圖2(c))作為LBM水驅模擬平臺,模擬像素精度為0.1 μm,模擬結果輸出間隔為0.1 s。模型參數(shù)設置如下:模擬溫度T=100 ℃,壓力條件p=19 MPa,模擬地層水黏度μw=1.0 mPa·s,模擬地層原油黏度μo=1.45 mPa·s,模擬水相密度ρw=1.0 g·cm-3,模擬原油密度ρo=0.78 g·cm-3。

3 模擬結果

結合研究區(qū)實際參數(shù)設定模型系統(tǒng)參數(shù)后,可應用給予Shan-Chen模型的D3Q27模型開展水驅模擬過程模擬。模擬一共分為2個階段:階段Ⅰ為充注油階段,階段Ⅱ為水驅油階段。階段Ⅰ充注油的目的是為階段II水驅油模擬生成初始油水分布。模型計算中時間步長度10-7s,為了加快輸出和后處理的速度,每間隔0.1 s設定一個輸出步(記s為輸出步數(shù)量)。

階段Ⅰ充注油階段,該階段油相充注速度為100 μm/s,從模型左側注入油相,恒速驅替到含水飽和趨向不變?yōu)橹?。圖3為油充注過程模擬得到的油水分布圖(其中固體礦物透明化處理)。圖3(a)為初始飽和水狀態(tài),模擬樣品烘干后飽和水的狀態(tài);圖3(b)～(d)為充注過程中含油飽和度分布。圖3表明:①隨著充注的進行油相充填的孔隙不斷增多;②油相前緣為凸狀形態(tài)(圖3(b));③模型右側靠近出口端的水相分布大于模型入口端附近的水相分布;④充注過程中出現(xiàn)了黏性指進現(xiàn)象(圖3(c));⑤在被油填充的孔隙壁面上仍然發(fā)現(xiàn)有潤濕相的水存在(圖3(d)),有一部分吼部較小的孔隙在油充注過程中未被填充,油相優(yōu)先填充了后續(xù)孔吼半徑較大的孔隙。

圖3 油驅水過程模擬及巖心中油相初始分布Fig.3 Simulations of oil displacing water and initial oil distribution in digital rock

階段II水驅油階段,該階段在階段Ⅰ之后進行,從模型左側注入水相。階段Ⅱ水驅油又根據(jù)注入條件又分成恒壓條件下水驅油和恒速條件下水驅油兩類。其中恒壓條件下水驅油為研究不同注入壓力對驅替效率的影響,本次研究共設定了4個注入壓差,4個模擬方案分別為pH2O_p1、pH2O_p2、pH2O_p3和pH2O_p4(表1)。另外恒速條件下水驅油為研究不同注入速度對驅替效率的影響,本次研究共設定了4個注入速度,4個模擬方案分別為vH2O_v1、vH2O_v2、vH2O_v3和vH2O_v4(表1)。

表1 水驅數(shù)值模擬方案

3.1 恒壓條件下三維動態(tài)驅油過程模擬

恒壓條件下水驅油過程中油水分布變化與恒速驅替過程中相似,對比圖4(a)和(b)可以看出高壓p2條件下模型pH2O_p2的驅替速度明顯高于p1條件下pH2O_p1模型的驅替速度。同時恒壓驅替過程中,最終殘余油分布區(qū)域也與恒速驅替下的殘余油分布區(qū)域相似。

圖4 恒壓條件下水驅油過程模擬Fig.4 Numerical simulation of water displacing oil with constant injecting pressure

3.2 恒壓驅替效率和換油率

基于LBM水驅模擬結果,可以計算恒壓條件下含油飽和度隨著注入水量的變化,以及驅替效率和換油率隨著注入量的變化。圖5為不同注入壓力水驅飽和度隨時間變化關系。從圖5可以看出:①恒壓注入條件下含油飽和度隨著水驅的進行開始逐漸下降;②4個模型的含油飽和度下降速度關系為pH2O_p1

圖5 不同注入壓力水驅飽和度隨時間變化關系Fig.5 Curves of oil saturation versus time with different injecting pressure

圖6為不同注入壓力水驅換油率和驅替效率與注入體積關系(圖中橫坐標為注入液量Vi與數(shù)字巖心模型孔隙總體積V的比值)。從圖6可以看出:恒壓條件下水驅換油率同樣隨著注入量的增加先增加后減少;4種注入壓力下,最優(yōu)換油率注入體積為1.1V;恒壓水驅驅替效率隨著注入量增加首先增加,然后漸進平穩(wěn);當注入量大于1.5V以后,驅替效率逐漸趨向于穩(wěn)定;當累積注入量為2.4V時,不同注入壓力下的驅替效率基本不變。對比不同驅替壓差模擬條件下的含油飽和度-時間變化曲線以及驅替效率-注入體積曲線可以發(fā)現(xiàn),致密油恒壓驅替壓差-驅替效率曲線存在最優(yōu)拐點。當驅替壓差從3.0×104Pa增加到3.5×104Pa后,模擬驅替效率僅增加了4.4%,而當驅替壓差從3.5×104Pa增加到4.0×104Pa時,模擬驅替效率僅增加了0.6%。

圖6 不同注入壓力水驅換油率和驅替效率與累積注入體積關系Fig.6 Curves of oil exchange ratio versus cumulative injection volume and displacement efficiency versus cumulative injection volume with different injection pressure

表2為最終驅替效率-注入壓力統(tǒng)計。驅替試驗所用壓力梯度與模擬所用壓力梯度相同。對比試驗驅替結果和數(shù)值模擬驅替結果表明,數(shù)值模擬得到的恒壓水驅的平均驅替效率為47.3%,物理驅替試驗得到恒壓水驅平均驅替效率為44.7%,4種條件下數(shù)值模擬驅替結果與室內試驗結果吻合度均大于90%。

表2 不同注入壓力下最終驅替效率Table 2 Ultimate displacement efficiency with different injection pressure

3.3 恒速條件下三維動態(tài)驅油過程模擬

恒速條件下共開展了4種不同注入速度下的水驅油模擬(圖7)。vH2O_v1、vH2O_v2、vH2O_v3和vH2O_v4的水驅油模型,均基于相同的模型,以增加模擬的可對比性,同時4個模型中也使用了相同的油充注結果,消除初始模型的影響。由圖7可以看出:①隨著水驅油過程的進行,油相分布逐漸減少,水相分布逐漸增加;②水相分布也出現(xiàn)黏性指進現(xiàn)象,在水相掃略過程中部分油相滯留下來;③vH2O_v2模型水相驅替速度大于vH2O_v1(圖7(a)和(b));④vH2O_v1和vH2O_v2模型,最終殘余油相分布區(qū)域相似(圖7(a)和(b)),都集中在模型的中部區(qū)域。圖7(c)和(d)為vH2O_v3和vH2O_v4模型的驅替過程模擬結果,與圖7(a)和(b)對比表明,最終殘余油分布與vH2O_v1和vH2O_v2模型相似,可見對于相同的孔隙空間分布,采用水驅油方式時,最終殘余油分布的有利區(qū)域隨注入速度改變較小。

3.4 恒速條件下驅替效率和換油率

根據(jù)LBM水驅模擬結果,通過綜合計算可以給出含油飽和度隨著注入水量的變化,并據(jù)此計算了驅替效率和換油率隨著注入量的變化。圖8為不同注入速度水驅飽和度-時間變化關系。從圖8可以看出:①含油飽和度隨著水驅的進行開始逐漸下降;② 4個模型的含油飽和度下降速度關系為vH2O_v1

圖9為不同注入速度水驅換油率和驅替效率與注入體積關系。從圖9可以看出:水驅換油率隨著注入量的增加先增加后減少;4種注入速度下,最優(yōu)換油率注入體積為2.2V;水驅驅替效率隨著注入量增加首先增加,然后漸進平穩(wěn);當注入量大于2V以后,驅替效率逐漸趨向于穩(wěn)定;當累積注入量為3V時,不同注入速度下的驅替效率關系為vH2O_v1

圖7 恒速條件下水驅油過程模擬Fig.7 Numerical simulation of water displacing oil with constant injecting velocity

圖8 不同注入速度水驅飽和度隨時間變化關系Fig.8 Curves of oil saturation versus time with different injecting velocity

圖9 不同注入速度水驅換油率和驅替效率與累積注入體積關系Fig.9 Curves of oil exchange ratio versus cumulative injection volume and displacement efficiency versus cumulative injection volume with different injection velocity

表3為不同注入速度水驅最終驅替效率-注入速度統(tǒng)計結果。試驗所用注入速度與模擬注入速度相同。表3表明隨著注入速度的增加,最終驅替效率逐漸增加。這與圖8和9的分析結果一致,總體來看,數(shù)值模擬得到的恒速水驅的平均驅替效率為49.6%,驅替物理試驗得到的恒速水驅的平均驅替效率為47.1%,4種注入速度條件下,試驗驅替效率與模擬驅替效率的吻合度均大于90%。

表3 不同注入速度下最終驅替效率Table 3 Ultimate displacement efficiency with different injection velocity

3.5 分割閾值對模擬結果的影響討論

分割閾值是影響數(shù)字巖心物性的關鍵參數(shù)之一,在將灰度圖轉化為二值的黑白圖像過程,分割閾值控制著孔隙的連通性,也隨之影響模擬結果。為了進一步分析分割閾值對于模擬結果的影響。選取4個分割閾值,分別為47、49、51和53,構建相應的數(shù)字巖心,并以恒速驅替為例,分析分割閾值對于驅替效果的影響。表4為分割閾值對數(shù)字巖心模型總孔隙度和連通孔隙度的影響。結果表明,數(shù)字巖心總孔隙度隨著分割閾值增大而升高,連通孔隙度也相應升高。同時表4表明,隨著分割閾值的減少,不連通孔隙在總孔隙度中所占的比例也不斷提高,當分割閾值小于等于47時,孔隙無法建立連通網(wǎng)絡。

表4 數(shù)字巖心模型分割閾值對總孔隙度和連通孔隙度的影響Table 4 Effects of threshold of digital rock model on the total porosity and connective porosity

圖10為不同分割閾值、相同注入速度(25 μm/s)條件下,驅替數(shù)值模擬得到的飽和度-驅替時間變化曲線(圖10(a))和驅替效率-注入體積變化曲線(圖10(b))。結果表明分割閾值對于初始含油飽和度和驅替效率存在較強的控制作用。分割閾值與初始含油飽和度和驅替效率存在正相關關系,分割閾值越小,有效的連通孔隙度越小,通過階段I充注試驗得到的初始含油飽和度越小,最終驅替效率也較小。對比不同分割閾值得到的模擬驅替效率與試驗得到的驅替效率(表4),表明當分割閾值Th=51時得到的連通孔隙度為9.65%,與試驗氣測孔隙度(9.7%)符合較好,同時模擬得到的驅替效率(46.1%)也與試驗得到的驅替效率(45.9%)符合較好;而當分割閾值變?yōu)?9時,數(shù)字巖心的連通孔隙度(6.63%)與試驗氣測孔隙度(9.7%)符合度變差,模擬得到的驅替效率與試驗驅替效率符合度也相應變差。

圖10 數(shù)字巖心模型分割閾值對水驅飽和度和驅替效率影響關系Fig.10 Effects of threshold of digital rock model on oil saturation and displacement efficiency

4 結 論

(1)通過聯(lián)合三維數(shù)字巖心模型和Shan-Chen格子玻爾茲曼多相流模擬方法,可以定量分析致密油充注及水驅微觀過程,以及含油飽和度、驅替效率和換油率動態(tài)變化特征。

(2)恒速狀態(tài)下水驅換油率隨著注入量的增加先增加后減少。4種注入速度下,最優(yōu)換油率注入體積為2.2V;恒速狀態(tài)下水驅驅替效率隨著注入量增加首先增加,然后漸進平穩(wěn)。當注入量大于2V以后,驅替效率逐漸趨向于穩(wěn)定;當累積注入量為3V時,不同注入速度下的驅替效率關系為vH2O_v1

(3)恒壓注入條件下含油飽和度隨著水驅的進行開始逐漸下降。4個模型的含油飽和度下降速度關系為pH2O_P1

(4)恒壓條件下水驅換油率同樣隨著注入量的增加先增加后減少;4種注入壓力下,最優(yōu)換油率注入體積為1.1V;恒壓水驅驅替效率隨著注入量增加首先增加,然后漸進平穩(wěn);當注入量大于1.5V后,驅替效率逐漸趨向于穩(wěn)定;當累積注入量為2.4V時,不同注入壓力下的驅替效率基本不變。

(5)數(shù)字巖心分割閾值對于模擬初始含油飽和度和驅替效率存在較強的控制作用,分割閾值越小,有效的連通孔隙度越小,充注模擬得到的初始含油飽和度越小,最終驅替效率也較小。