李 棟,鄭雙金，任春梅，趙雪峰，濮 御，4，常澤輝

(1.東北石油大學(xué) 土木建筑工程學(xué)院，黑龍江 大慶 163318； 2.大慶石油管理局 水務(wù)公司制水第三分公司，黑龍江 大慶 163000； 3.大慶油田工程有限公司，黑龍江 大慶 163453； 4.東北石油大學(xué) 提高油氣采收率教育部重點(diǎn)實驗室，黑龍江 大慶 163318； 5.內(nèi)蒙古工業(yè)大學(xué) 能源與動力工程學(xué)院，內(nèi)蒙古 呼和浩特 010051)

0 引言

頁巖油在開采過程中孔隙尺度流動特性復(fù)雜，容易發(fā)生運(yùn)移殘留現(xiàn)象，對頁巖油高效開發(fā)產(chǎn)生不利影響[1-4]。因此，研究頁巖油孔隙尺度運(yùn)移殘留規(guī)律對提高頁巖油采油效率具有重要意義[5-7]。人們對多孔介質(zhì)石油類運(yùn)移進(jìn)行研究[8-10]，通過將孔隙尺度儲層與油類介質(zhì)流動特性進(jìn)行耦合，分析頁巖油在儲層中滲流特征[11-13]。雷浩等[14]建立頁巖油封閉式高精度滲流評價系統(tǒng)，開展頁巖油巖心滲流實驗，分析頁巖油在儲層中的滲流特征。劉向君等[15]建立三維數(shù)字巖心模型，研究微觀孔隙結(jié)構(gòu)氣水兩相運(yùn)移規(guī)律。宋文輝等[16]建立不同尺度下數(shù)字巖心模型，研究頁巖油儲層多尺度孔隙結(jié)構(gòu)特征與油相流動能力。

由于頁巖油在孔隙尺度儲層運(yùn)移復(fù)雜，在耦合研究時很難準(zhǔn)確分析流動特性，有必要研究油類和非油類介質(zhì)兩相流孔隙尺度運(yùn)移特性[17-18]。高亞軍等[19]建立油水兩相微觀滲流數(shù)學(xué)模型，采用Level Set數(shù)學(xué)方法和N-S控制方程，分析兩相微觀孔隙運(yùn)移規(guī)律，發(fā)現(xiàn)優(yōu)勢滲流通道特征。王寶等[20]采用油水互驅(qū)微觀實驗，分析油水兩相微觀介質(zhì)滲流規(guī)律。油類介質(zhì)孔隙尺度運(yùn)移影響因素較多，如注入速度、流體黏度比和流體密度等。ZHANG Faqiang等[21]建立孔隙網(wǎng)絡(luò)模型，研究流體黏度比、界面張力和流體密度對石油運(yùn)移的影響。張莉等[22]采用數(shù)值模擬和壓強(qiáng)分析法，分析微觀剩余油驅(qū)替動力和驅(qū)替阻力的影響因素，發(fā)現(xiàn)增加驅(qū)替流度比、降低界面張力可以提高驅(qū)替倍數(shù)。DAI Caili等[23]采用巖心動態(tài)滲吸實驗，研究注入速度和界面張力等因素對致密油運(yùn)移的影響。人們對石油類介質(zhì)運(yùn)移特性和運(yùn)移影響進(jìn)行研究，但有關(guān)頁巖油孔隙尺度運(yùn)移殘留特性的研究較少。

基于水平集方法和N-S控制方程，筆者建立頁巖油多孔介質(zhì)孔隙尺度運(yùn)移模型，采用COMSOL Multiphysics軟件和向后差分法等求解，研究頁巖油在多孔介質(zhì)中的運(yùn)移過程，考慮多孔介質(zhì)中頁巖油初始運(yùn)移速度、密度和動力黏度等因素對運(yùn)移過程中頁巖油體積分?jǐn)?shù)的影響，分析頁巖油在多孔介質(zhì)中的運(yùn)移規(guī)律和運(yùn)移殘留特征。

1 模型建立

1.1 物理模型

根據(jù)微米級CT掃描成像構(gòu)建多孔介質(zhì)結(jié)構(gòu)[24]，選取部分具有代表性區(qū)域建立物理模型，兩相流區(qū)域多孔介質(zhì)幾何模型見圖1。其中多孔介質(zhì)區(qū)域尺寸為90.0 μm×80.0 μm，填充空氣介質(zhì)，左側(cè)區(qū)域設(shè)置尺寸為15.0 μm×80.0 μm的初始區(qū)域且填充油類介質(zhì)。左側(cè)區(qū)域的左側(cè)邊設(shè)置速度入口邊界，多孔介質(zhì)區(qū)域的右側(cè)邊設(shè)置壓力出口邊界。整個幾何模型的上下側(cè)設(shè)置對稱邊界條件。多孔介質(zhì)區(qū)域內(nèi)孔隙度為39.1%。

假設(shè)條件：多孔介質(zhì)骨架內(nèi)初始狀況填充氣相，油相從左側(cè)區(qū)域注入多孔介質(zhì)骨架；頁巖油在多孔介質(zhì)孔隙尺度運(yùn)移過程中不受重力的影響；頁巖油在多孔介質(zhì)孔隙尺度下運(yùn)移時，不考慮溫度產(chǎn)生的影響。

1.2 數(shù)學(xué)模型

頁巖油多孔介質(zhì)孔隙尺度運(yùn)移由水平集方法和N-S控制方程，多孔介質(zhì)骨架內(nèi)充滿空氣，頁巖油在多孔介質(zhì)內(nèi)運(yùn)移的本質(zhì)是氣液兩相驅(qū)替過程。N-S控制方程用于描述基于連續(xù)介質(zhì)的假定黏性流體動量守恒?？紤]表面張力，建立兩相流動不可壓縮N-S控制方程[25]為

ρ(u·)u=·[-pI+μ(u+(u)T)]+Fst，

(1)

ρ·u=0,

(2)

式中：ρ為流體密度；μ為流體動力黏度；u為流體速度；p為流體壓力；I為單位矩陣；Fst為兩相表面張力。

水平集方法是歐拉法中經(jīng)典的界面追蹤方法，主要用于描述運(yùn)動界面動態(tài)變化。在兩相流驅(qū)替運(yùn)動中引入水平集函數(shù)φ，表示水平集變量中兩相體積分?jǐn)?shù)，0≤φ≤1。當(dāng)φ=0時，流體為空氣介質(zhì)；當(dāng)φ=1時，流體為頁巖油介質(zhì)。兩相流運(yùn)動水平集方程為

(3)

式中：t為時間；γ為重新初始化參數(shù)；εIS為界面厚度控制參數(shù)。

頁巖油多孔介質(zhì)孔隙尺度運(yùn)移模型中流體物性[25]表示為

(4)

式中：ρg和ρw分別為空氣和頁巖油密度；μg和μw分別為空氣和頁巖油動力黏度。

初始狀態(tài)多孔介質(zhì)骨架內(nèi)為飽和氣相，左側(cè)矩形內(nèi)為油相，油氣兩相的接觸面為兩相的初始界面，左右兩側(cè)的進(jìn)出口流動方程表示為

u=u0n,

(5)

[-pI+μ(u+(u)T)]n=-p0n，

(6)

式(5—6)中：u0為入口頁巖油流速；p0為邊界處相對壓強(qiáng)；n為界面法向量，n=φ/︳φ∣。

在出口邊界設(shè)置抑制回流，調(diào)整出口壓力并防止流體通過邊界進(jìn)入多孔介質(zhì)骨架。

在水平集方法中，通過滑移長度β和接觸角θw實現(xiàn)不同潤濕性的設(shè)定，為防止界面移動速度過大，引入表面摩擦力Fsr約束滑移長度[26]：

(7)

接觸角θw為壁面與頁巖油界面處的夾角，由液相指向氣相，由楊式方程限定的接觸角表示為

σcosθw+γsw=γsg，

(8)

式中：σ為空氣與頁巖油表面的張力系數(shù)；γsw為頁巖油面能量密度；γsg為空氣面能量密度。

分子吸引力在液相與氣相的分界面處產(chǎn)生極其微小的拉力，即表面張力Fst[26]表示為

Fsr=δ(φ)[-σn(·n)+(1-nnT)σ]，

(9)

式中：δ(φ)為Dirac函數(shù)，在相界面上為0時，Dirac函數(shù)可表示為

δ(φ)=σ︳φ∣φ(1-φ)。

(10)

模型中多孔介質(zhì)骨架孔隙和喉道復(fù)雜區(qū)域網(wǎng)格細(xì)密，孔隙和喉道單一區(qū)域網(wǎng)格稀疏，采用COMSOL Multiphysics軟件的三角形網(wǎng)格對幾何區(qū)域進(jìn)行剖分，其中最大單元尺寸設(shè)置分別為5.36、3.60和2.80 μm，分別生成 51 184、75 546、98 051個網(wǎng)格。當(dāng)網(wǎng)格數(shù)超過51 184時，頁巖油在一次運(yùn)移后體積分?jǐn)?shù)誤差不超過5%，以達(dá)到網(wǎng)格獨(dú)立解。模型采用層流兩相流模塊和N-S控制方程耦合水平集模塊，設(shè)定頁巖油瞬態(tài)運(yùn)移為0.005 s進(jìn)行求解。時間步進(jìn)方法采用向后差分法，求解器計算步長設(shè)置為自由，瞬態(tài)求解器的容差因子為0.05，采用非線性控制器控制求解。

2 數(shù)值模擬結(jié)果

通過數(shù)值模擬分析頁巖油在多孔介質(zhì)中的運(yùn)移殘留規(guī)律，以及頁巖油初始運(yùn)移速度、密度和動力黏度等因素對頁巖油體積分?jǐn)?shù)及運(yùn)移周期的影響，研究頁巖油多孔介質(zhì)孔隙尺度運(yùn)移特征。

2.1 運(yùn)移殘留規(guī)律

選取界面入口速度為4.00 m·s-1，對氣油兩相流模擬結(jié)果進(jìn)行分析，隨時間的增加，頁巖油多孔介質(zhì)運(yùn)移過程見圖2，其中藍(lán)色代表油相，即φ=1；綠色代表氣相，即φ=0；介于藍(lán)色與綠色之間的部分同時包含兩相，即0<φ<1。

由圖2可知，頁巖油初次運(yùn)移過程符合兩相流驅(qū)替模擬基本特征[25]。頁巖油在具有孔喉和孔徑結(jié)構(gòu)的多孔介質(zhì)中運(yùn)移，隨運(yùn)移的進(jìn)行，頁巖油在多孔介質(zhì)中體積占比增加，空氣體積占比減小。頁巖油形成小尺度的指進(jìn)現(xiàn)象(見圖2(a))，在多孔介質(zhì)中優(yōu)先大孔隙通道運(yùn)移(見圖2(b))，原因是在毛細(xì)指進(jìn)過程中毛細(xì)管力對流動過程起主導(dǎo)作用，當(dāng)其中一條流動路徑到達(dá)出口端，就形成優(yōu)先流動，頁巖油在一定的多孔介質(zhì)區(qū)域內(nèi)形成優(yōu)勢通道。在多孔介質(zhì)內(nèi)，部分孔道在初次運(yùn)移時無法被油相波及而成為殘余氣的主要賦存空間(見圖2(c))。在多孔介質(zhì)通道運(yùn)移過程中，頁巖油主流通道旁孔隙內(nèi)賦存殘余氣，原因是當(dāng)黏性力占主導(dǎo)地位時，相界面移動受毛細(xì)管力的影響越來越小，相界面在黏性力作用下快速突破相鄰的喉道，導(dǎo)致部分喉道內(nèi)的氣體被圈閉。左側(cè)區(qū)域的頁巖油容積初始時是固定的，隨頁巖油逐步向多孔介質(zhì)運(yùn)移，可以看出比較明顯的綠色和青色區(qū)域殘余氣(見圖2(d))，殘余氣分為單孔道和多孔道殘余氣，大多分布在流動路徑分支附近，原因是流動路徑在運(yùn)移過程中自主通道向支路橫向發(fā)展，導(dǎo)致位于前緣后部的殘余氣被截斷，殘余氣滯留在多孔介質(zhì)區(qū)域內(nèi)。頁巖油在多孔介質(zhì)內(nèi)的初次運(yùn)移基本結(jié)束，二次運(yùn)移開始(見圖2(f))。頁巖油多次運(yùn)移是以優(yōu)勢通道為主(見圖2)，初次運(yùn)移時殘余氣區(qū)域由頁巖油占據(jù)而形成殘余油，是由多孔介質(zhì)前緣橫向發(fā)展和多孔介質(zhì)區(qū)域空間限制等因素引起的。

頁巖油在多孔介質(zhì)運(yùn)移過程中體積分?jǐn)?shù)隨時間變化曲線見圖3，其中字母標(biāo)注具有代表性的波峰和波谷。A-H點(diǎn)對應(yīng)圖2(e-l)的頁巖油體積分?jǐn)?shù)分別為0.862、0.372、0.821、0.313、0.431、0.356、0.631和0.269。A點(diǎn)時，頁巖油初次運(yùn)移且?guī)缀醭錆M多孔介質(zhì)骨架，頁巖油體積分?jǐn)?shù)在整個運(yùn)移過程中達(dá)到最大，只有部分殘余氣留存；B點(diǎn)時，頁巖油初次運(yùn)移結(jié)束；C點(diǎn)時，頁巖油二次運(yùn)移充滿多孔介質(zhì)骨架；D點(diǎn)時，頁巖油二次運(yùn)移結(jié)束，油相將殘余氣排出，殘余油留存；E點(diǎn)時，在壓差作用下，左側(cè)區(qū)域下方的頁巖油運(yùn)移到多孔介質(zhì)區(qū)域；F點(diǎn)時，頁巖油運(yùn)移出多孔介質(zhì)區(qū)域；G點(diǎn)時，頁巖油再次運(yùn)移到多孔介質(zhì)區(qū)域，受主流優(yōu)勢作用影響，頁巖油在運(yùn)移時以大通道為主；H點(diǎn)時，頁巖油再次運(yùn)移結(jié)束，殘余油留存。在相界面和表面張力作用下，隨時間的增加，頁巖油在多孔介質(zhì)內(nèi)體積分?jǐn)?shù)的波動基本呈周期性變化。

模型穩(wěn)定運(yùn)移是指在多孔介質(zhì)內(nèi)頁巖油體積分?jǐn)?shù)呈明顯周期性變化的運(yùn)移狀態(tài)(見圖3中N點(diǎn)后頁巖油的運(yùn)移狀態(tài))。在多孔介質(zhì)運(yùn)移穩(wěn)定后，頁巖油在多孔介質(zhì)內(nèi)的運(yùn)移過程呈周期性變化，穩(wěn)定運(yùn)移時波峰為n-r，波谷為N-R。模型穩(wěn)定運(yùn)移后，頁巖油波峰、波谷體積分?jǐn)?shù)及相鄰波峰時間差隨運(yùn)移次數(shù)變化曲線見圖4。由圖4可知，隨頁巖油運(yùn)移次數(shù)的增加，多孔介質(zhì)內(nèi)頁巖油在波峰和波谷處的體積分?jǐn)?shù)逐漸增加，頁巖油運(yùn)移時間逐漸減小。在表面張力作用下，隨頁巖油在多孔介質(zhì)內(nèi)運(yùn)移次數(shù)的增加，多孔介質(zhì)內(nèi)殘余油增加，頁巖油運(yùn)移周期變短。

2.2 影響因素

2.2.1 初始運(yùn)移速度

選取頁巖油初始運(yùn)移速度為0.08、0.10、0.20、0.40和0.80 m·s-1開展運(yùn)移數(shù)值模擬，在多孔介質(zhì)內(nèi)頁巖油體積分?jǐn)?shù)變化曲線見圖5。由圖5可知，當(dāng)初始運(yùn)移速度為0.08、0.10 m·s-1時，頁巖油體積分?jǐn)?shù)不呈周期性變化，在多孔介質(zhì)中沒有達(dá)到穩(wěn)定運(yùn)移的條件。在多孔介質(zhì)運(yùn)移穩(wěn)定后，隨時間的推移，頁巖油在波峰和波谷處的體積分?jǐn)?shù)隨初始運(yùn)移速度的變化較大，在多孔介質(zhì)孔隙尺度運(yùn)移時，初始運(yùn)移速度的變化對殘余油含量的影響較大。頁巖油在多孔介質(zhì)內(nèi)運(yùn)移時受表面張力作用控制，初始運(yùn)移速度為0.40 m·s-1時，在波谷處的體積分?jǐn)?shù)最小，在多孔介質(zhì)內(nèi)運(yùn)移留下的殘余油體積分?jǐn)?shù)最少。在多孔介質(zhì)運(yùn)移穩(wěn)定后，在相同時間范圍內(nèi)，根據(jù)頁巖油體積分?jǐn)?shù)波峰數(shù)可以判斷運(yùn)移周期數(shù)。初始運(yùn)移速度為0.20 m·s-1時，頁巖油具有9個波峰，在多孔介質(zhì)內(nèi)運(yùn)移周期為9個。初始運(yùn)移速度分別為0.40、0.80 m·s-1時，頁巖油在多孔介質(zhì)內(nèi)運(yùn)移周期分別為11、12個。這表明隨初始運(yùn)移速度的增加，頁巖油運(yùn)移周期縮短。

2.2.2 頁巖油密度

選取頁巖油密度為728.1、733.6和754.4 kg·m-3開展運(yùn)移數(shù)值模擬，在多孔介質(zhì)內(nèi)頁巖油體積分?jǐn)?shù)變化曲線見圖6。由圖6可知，在多孔介質(zhì)運(yùn)移穩(wěn)定后，隨時間的推移，頁巖油在波峰處的體積分?jǐn)?shù)為0.870左右，波谷處的體積分?jǐn)?shù)為0.520左右，頁巖油在波峰和波谷處的體積分?jǐn)?shù)隨密度的變化較小。在多孔介質(zhì)內(nèi)運(yùn)移時，密度變化對殘余油體積分?jǐn)?shù)的影響較小。在多孔介質(zhì)運(yùn)移穩(wěn)定后，在相同時間范圍內(nèi)，頁巖油物性參數(shù)改變多孔介質(zhì)孔隙的接觸角和分子間作用力，不同頁巖油密度使運(yùn)移周期發(fā)生改變。密度分別為728.1、733.6和754.4 kg·m-3時，頁巖油在多孔介質(zhì)內(nèi)運(yùn)移周期分別為9、11和11.25個。這表明隨密度的增加，頁巖油運(yùn)移周期縮短。

2.2.3 動力黏度

選取頁巖油動力黏度為4.36×10-4、4.48×10-4和5.53×10-4Pa·s開展運(yùn)移數(shù)值模擬，在多孔介質(zhì)內(nèi)頁巖油體積分?jǐn)?shù)變化曲線見圖7。由圖7可知，在多孔介質(zhì)運(yùn)移穩(wěn)定后，隨時間的推移，頁巖油在波峰和波谷處的體積分?jǐn)?shù)隨動力黏度的變化較小，在多孔介質(zhì)中運(yùn)移時，動力黏度變化對殘余油體積分?jǐn)?shù)影響較小。在多孔介質(zhì)運(yùn)移穩(wěn)定后，頁巖油物性參數(shù)改變多孔介質(zhì)孔隙的接觸角和分子間作用力，不同頁巖油動力黏度使運(yùn)移周期發(fā)生改變，在相同時間范圍內(nèi)，動力黏度分別為4.36×10-4、4.48×10-4和5.53×10-4Pa·s 時，頁巖油在多孔介質(zhì)內(nèi)運(yùn)移周期分別為11、12和8個。動力黏度為4.48×10-4Pa·s時，頁巖油運(yùn)移周期最短。

3 結(jié)論

(1)基于水平集方法和N-S控制方程，建立多孔介質(zhì)孔隙尺度頁巖油運(yùn)移模型，研究頁巖油在多孔介質(zhì)中運(yùn)移的體積分?jǐn)?shù)變化，分析頁巖油進(jìn)入多孔介質(zhì)的初始運(yùn)移速度、密度和動力黏度對頁巖油在多孔介質(zhì)中運(yùn)移的影響，以及頁巖油多孔介質(zhì)孔隙尺度運(yùn)移殘留規(guī)律。

(2)在孔隙尺度下初次運(yùn)移時，頁巖油油相優(yōu)先從大通道中突破并形成優(yōu)勢流動通道，運(yùn)移路徑以優(yōu)勢流動通道為主。在多孔介質(zhì)穩(wěn)定運(yùn)移后，隨運(yùn)移次數(shù)的增加，多孔介質(zhì)殘余油越多，頁巖油運(yùn)移周期越短。

(3)頁巖油初始運(yùn)移速度對殘余油影響較大，初始運(yùn)移速度為0.40 m·s-1時，孔隙內(nèi)殘余油體積分?jǐn)?shù)最小，頁巖油密度和動力黏度對殘余油體積分?jǐn)?shù)影響較??；隨頁巖油初始運(yùn)移速度和密度的增加，運(yùn)移周期縮短；動力黏度為4.48×10-4Pa·s時，頁巖油運(yùn)移周期最短。