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

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

0 引言

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

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

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

1 模型建立

1.1 物理模型

根據(jù)微米級(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è)置壓力出口邊界。整個(gè)幾何模型的上下側(cè)設(shè)置對(duì)稱邊界條件。多孔介質(zhì)區(qū)域內(nèi)孔隙度為39.1%。

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

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

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

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

(1)

ρ·u=0,

(2)

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

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

(3)

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

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

(4)

式中：ρg和ρw分別為空氣和頁(yè)巖油密度；μg和μw分別為空氣和頁(yè)巖油動(dòng)力黏度。

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

u=u0n,

(5)

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

(6)

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

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

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

(7)

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

σcosθw+γsw=γsg，

(8)

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

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

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

(9)

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

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

(10)

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

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

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

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

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

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

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

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

2.2 影響因素

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

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

2.2.2 頁(yè)巖油密度

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

2.2.3 動(dòng)力黏度

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

3 結(jié)論

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

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

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