馬飛英, 王 林,, 何勇明, 王洪輝, 劉大偉

(1.廣東石油化工學(xué)院 石油工程學(xué)院, 廣東 茂名 525000； 2.油氣藏地質(zhì)及開發(fā)工程國(guó)家重點(diǎn)實(shí)驗(yàn)室(成都理工大學(xué)),成都 610059)

自發(fā)滲吸(簡(jiǎn)稱“自吸”)是指多孔介質(zhì)中潤(rùn)濕性流體在毛細(xì)管力的作用下自發(fā)驅(qū)替非潤(rùn)濕相的過程。潤(rùn)濕性是影響自吸的一個(gè)關(guān)鍵因素，學(xué)者們很早就開始了巖石的潤(rùn)濕性研究，但大部分的研究集中在均質(zhì)潤(rùn)濕性上[1]，而實(shí)際上地層中巖石的潤(rùn)濕性是非均質(zhì)的[2-3]。儲(chǔ)層巖石中含有多種礦物，有些礦物是親水的，如石英、長(zhǎng)石、伊利石等，有些礦物是親油的，如有機(jī)質(zhì)、黃鐵礦、含鐵綠泥石等[4-7]。甚至有些親水的巖石，在與地下原油長(zhǎng)期接觸后，原油中的極性物質(zhì)會(huì)改變巖石的潤(rùn)濕性，使原來親水的巖石變?yōu)橛H油的巖石[8]。因此儲(chǔ)層巖石中可能既存在親油孔隙，又存在親水孔隙，甚至還存在既親油又親水的混合潤(rùn)濕孔隙。親水的孔隙會(huì)自發(fā)吸水排油，親油的孔隙會(huì)自發(fā)吸油排水，那么混合潤(rùn)濕孔隙的自吸行為又是怎樣的呢？目前關(guān)于自吸的研究，大多將巖石孔隙看作是單一潤(rùn)濕的，很少考慮孔隙的混合潤(rùn)濕特征，這可能會(huì)導(dǎo)致一些認(rèn)識(shí)出現(xiàn)偏差。

隨著常規(guī)油氣資源的枯竭，頁巖油、頁巖氣等非常規(guī)油氣資源很快成為了常規(guī)油氣的接替能源[9]。頁巖的自吸對(duì)油氣采收率、井壁穩(wěn)定性、返排率等有重要的影響[6,10-13]。與常規(guī)砂巖和碳酸鹽巖等儲(chǔ)層不同的是，頁巖油、頁巖氣等非常規(guī)儲(chǔ)層巖石中有機(jī)質(zhì)含量較高，而這些有機(jī)質(zhì)通常是親油(疏水)的，大部分無機(jī)物卻是親水的，因此頁巖油氣儲(chǔ)層巖石中存在親油(疏水)的有機(jī)孔、親水的無機(jī)孔以及由有機(jī)物和無機(jī)物共同構(gòu)成的混合潤(rùn)濕孔[14-16]，圖1展示了這3種類型的孔隙[16]。有些頁巖的混合潤(rùn)濕孔隙比例是非常高的，Ishank Gupta等[17]通過場(chǎng)發(fā)射掃描電鏡圖對(duì)Marcellus、Eagle Ford以及Woodford地區(qū)的頁巖樣品進(jìn)行了分析，發(fā)現(xiàn)這3個(gè)地區(qū)的巖樣混合潤(rùn)濕孔隙的比例分別為57%、69%、68%，具有非常高的混合潤(rùn)濕孔隙比例?？紫兜幕旌蠞?rùn)濕特征使得頁巖儲(chǔ)層的自吸與常規(guī)儲(chǔ)層存在差別。構(gòu)成孔隙的礦物成分、含量等參數(shù)不同，會(huì)產(chǎn)生不同的自吸結(jié)果；然而國(guó)內(nèi)外學(xué)者對(duì)自吸的研究主要集中在孔隙結(jié)構(gòu)、流體物性等對(duì)自吸能力的影響上[18-22]，對(duì)混合潤(rùn)濕孔隙的自吸鮮有研究。本文通過建立平板混合潤(rùn)濕模型，定量分析自吸速度與界面張力、縫寬、接觸角、黏度、界面位置的函數(shù)關(guān)系，為判斷頁巖孔隙是否發(fā)生水的自吸提供理論依據(jù)。

圖1 Marcellus頁巖場(chǎng)發(fā)射掃描電鏡圖[16]Fig.1 SEM image from the Marcellus shale sample

1 平板模型的自吸受力分析

無機(jī)孔(疏水)、有機(jī)孔(親水)屬于相對(duì)簡(jiǎn)單的單一潤(rùn)濕孔隙，目前已有大量學(xué)者對(duì)其自吸進(jìn)行了研究，而混合潤(rùn)濕孔隙中流體的自吸鮮有報(bào)道。由于裂縫形狀與平板模型相似，于是將裂縫系統(tǒng)簡(jiǎn)化為平板模型。下面我們將通過建立平板模型對(duì)混合潤(rùn)濕孔隙中流體的微觀受力進(jìn)行分析。

根據(jù)學(xué)者們對(duì)頁巖的微觀圖像研究，頁巖孔隙中存在一部分壁面親水，另一部分壁面疏水的混合潤(rùn)濕孔隙[16-17]。為了方便定量描述，將一面為親水的無機(jī)物，另一面為疏水的有機(jī)物的裂縫系統(tǒng)簡(jiǎn)化為混合潤(rùn)濕平板模型(圖2)。在建立方程之前，首先需要進(jìn)行幾個(gè)假設(shè)：①假設(shè)裂縫面為光滑的；②不考慮流體重力的影響；③忽略氣體的黏滯力；④平板模型水平放置；⑤流體的流動(dòng)狀態(tài)為層流。

圖2 混合潤(rùn)濕平板模型中氣-水自吸示意圖Fig.2 Schematic diagram of gas-water spontaneous imbibition in a mixed wetting plate model

自吸模型如圖2所示，流體受到界面張力與黏滯阻力作用。m、n點(diǎn)分別為氣、液、固三相接觸點(diǎn)，x軸向右為正方向、向左為負(fù)方向，則氣-水在x軸方向上受到的合力為

Fr=F1+F2+Ff

(1)

式中：Fr為x軸方向上的合力；F1為流體在m點(diǎn)所受到的x軸方向上的力；F2為流體在n點(diǎn)所受到的x軸方向上的力；Ff為流體在x軸方向上受到的摩擦力。

根據(jù)界面張力知識(shí)，界面張力垂直三相接觸線，與氣-水界面相切，則疏水界面上的流體在x軸方向上受到的力為

F1=Lσgwcosθ1

(2)

式中：L為氣、液、固三相接觸線的長(zhǎng)度；σgw為天然氣與水的界面張力；θ1為氣-水在面板M上的接觸角。

親水界面上的流體在x軸方向上受到的力為

F2=Lσgwcosθ2

(3)

式中：θ2為氣-水在面板N上的接觸角。

假設(shè)在平板壁面與兩平板的中心線之間流體流動(dòng)的速度呈線性分布，則根據(jù)牛頓內(nèi)摩擦定律，平板單位面積上受到的切應(yīng)力為

τ=-4vμw/d

(4)

式中：τ為流體的切應(yīng)力；μw水的黏度；v為氣-水界面的移動(dòng)速度；d為平板間的距離。

水與上、下平板接觸的面積為

A=2Lx

(5)

式中：A為水與上、下平板接觸的面積；x為氣-水界面運(yùn)移的距離。

由于氣體的黏度相對(duì)于水的黏度小得多，因此可以忽略氣體的黏度，則流體受到上下兩塊平板壁面的摩擦阻力為

Ff=2Aτ=-8vμwLx/d

(6)

將式(2)、式(3)、式(6)代入式(1)可求得流體在x軸方向上的合力

Fr=Lσgw(cosθ1+cosθ2)-8vμwLx/d

(7)

根據(jù)牛頓第二定律，可得

Lσgw(cosθ1+cosθ2)-8vμwLx/d=mfa

(8)

式中：mf為流體質(zhì)量；a為流體加速度。

由于微小孔隙中流體的雷諾數(shù)較小，慣性力遠(yuǎn)小于黏滯力，因此可以忽略慣性力的影響，則式(8)可轉(zhuǎn)化為

σgw(cosθ1+cosθ2) = 8vμwx/d

(9)

根據(jù)式(9)，可求得氣-水界面的運(yùn)動(dòng)速度為

v=σgwd(cosθ1+cosθ2)/(8μwx)

(10)

式(10)在形式上與Lucas-Washburn方程相似[23-24]。由于式(10)采用的推導(dǎo)方法與Lucas-Washburn方程的推導(dǎo)方法相同，因此將式(10)稱為混合潤(rùn)濕平板模型的Lucas-Washburn方程。

為了更加直觀地表達(dá)混合潤(rùn)濕孔隙的接觸角對(duì)自吸方向的影響，還需對(duì)式(10)做進(jìn)一步的轉(zhuǎn)化。

根據(jù)三角函數(shù)關(guān)系，可得

(11)

由于平板M疏水，平板N親水，根據(jù)親水與疏水接觸角的劃分規(guī)定，可知接觸角的范圍

(π/2)<θ1<π

(12)

0<θ2<(π/2)

(13)

根據(jù)式(12)和式(13)，可知

(14)

(15)

根據(jù)式(14)可知

(16)

由式(15)可知

cos(θ1/2-θ2/2)>0

(17)

結(jié)合式(16)、式(17)、式(11)以及式(10)，可知界面運(yùn)動(dòng)速度v的正、負(fù)方向

(18)

由式(18)可知，當(dāng)θ1+θ2<π時(shí)，v>0，此時(shí)只可自發(fā)吸水；當(dāng)θ1+θ2>π時(shí)，v<0，只可自發(fā)吸氣；當(dāng)θ1+θ2=π時(shí)，v=0，氣、水均無法自吸。

2 數(shù)值模擬

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

水平集方法可以用來追蹤兩相流體界面的運(yùn)動(dòng)[25]，水平集函數(shù)Φ為光滑函數(shù)。氣-水界面由水平集函數(shù)Φ為0.5的等值線表示；在氣體中，Φ=0；在水中，Φ=1。兩相流界面的運(yùn)動(dòng)方程可用下式表示[26]

(19)

式中：Φ為水平集函數(shù)；γ為初始化強(qiáng)度；ε為界面厚度控制參數(shù)。

兩相流體的密度和黏度分別用下式表示[26]

ρ=ρg+(ρw-ρg)Φ

(20)

μ=μg+(μw-μg)Φ

(21)

式中：ρ為兩相流體密度；μ為兩相流體黏度；ρg、ρw分別為氣、水的密度；μg、μw分別為氣、水的黏度。

控制方程包括了連續(xù)方程、結(jié)合了水平集函數(shù)的Navier-Stokes方程。

連續(xù)方程為

▽·v=0

(22)

考慮了界面張力作用下的Navier-Stokes方程為[26]

μ(v+(v)T)]+Fst

(23)

式中：p為壓力；I為單位矩陣；Fst為氣-水界面張力產(chǎn)生的力。

氣-水界面張力產(chǎn)生的力Fst可用下式表示[27]

Fst=σwokδn

(24)

式中：k為界面曲率；δ為Dirac delta 函數(shù)；n為界面單位向量。

δ可近似表示為[27]

δ=6|Φ(1-Φ)||Φ|

(25)

界面單位向量可用下式表示

(26)

界面曲率可用下式表示

(27)

2.2 數(shù)值模擬驗(yàn)證

采用Comsol Multiphysics仿真軟件建立三維兩相流水平集模型，設(shè)置幾何模型的平板長(zhǎng)(l)40 mm，寬(b)5 mm，上下兩平板的間距(d)為0.1 mm，上平板疏水，下平板親水，幾何模型如圖3所示。模型中的流體為天然氣與水，水的黏度為1 mPa·s，氣-水界面張力為0.07 N/m。初始?xì)?水界面(x=0 mm)距離入口端20 mm。平板出口與入口端壓力為10 MPa。滑移長(zhǎng)度為0.001 m。構(gòu)建好三維幾何模型后，采用自由四邊形掃掠方式生成六面體網(wǎng)格。為了分析接觸角對(duì)自吸方向的影響，模擬了不同親水壁和疏水壁接觸角下(表1)的氣-水界面運(yùn)動(dòng)位置(圖4)。從圖4可知，氣-水界面運(yùn)動(dòng)方向與式(18)判斷的結(jié)果一致。

表1 模型上-下壁面的接觸角Table 1 The contact angle between the upper and lower walls of the model

3 結(jié)果與討論

按照表1中的參數(shù)，運(yùn)用兩相流水平集方法和Lucas-Washburn方程法進(jìn)行計(jì)算，得到了不同時(shí)間下的界面運(yùn)動(dòng)速度與接觸角之和的關(guān)系曲線(圖5)，v為正數(shù)時(shí)，表示吸水排氣；v為負(fù)數(shù)時(shí)，表示吸氣排水。對(duì)比圖5-A與圖5-B，可以發(fā)現(xiàn)兩者的曲線變化趨勢(shì)以及自吸方向是相同的，但這兩種方法求得的速度大小是有差異的。造成這種差異的原因是Lucas-Washburn方程采用的是受力平衡時(shí)流體充分發(fā)展下建立起來的[28]，但實(shí)際的流體并未達(dá)到平衡，兩相流水平集數(shù)值模擬方法的結(jié)果更接近實(shí)際。雖然兩種方法求得的速度大小有差異，但兩種方法對(duì)流體自吸運(yùn)動(dòng)方向的判斷是一致的。當(dāng)接觸角之和小于3.14時(shí)，毛細(xì)管吸水排氣，且吸水速度隨接觸角的增大而減?。划?dāng)接觸角之和大于3.14時(shí)，流體自吸方向發(fā)生改變，毛細(xì)管變?yōu)槲鼩馀潘遗潘俣入S接觸角之和增大而增大。兩種方法得到運(yùn)動(dòng)速度為0時(shí)，對(duì)應(yīng)的接觸角之和均接近π， 這與式(18)相吻合。

圖4 不同接觸角之和條件下的氣-水界面位置圖Fig.4 Gas-water interface position under the condition of the sum of different contact angles 箭頭表示水的流動(dòng)方向， t=40 ms

圖5 界面運(yùn)動(dòng)速度與接觸角之和的關(guān)系曲線Fig.5 The relation curve between interface velocity and the sum of contact angle

接觸角影響毛細(xì)管力的大小，在親水和疏水物質(zhì)含量一定的情況下，親水物質(zhì)和疏水物質(zhì)的接觸角是影響毛細(xì)管自吸的關(guān)鍵因素。兩種物質(zhì)的接觸角之和越小，則吸水能力越強(qiáng)；兩種物質(zhì)的接觸角之和越大，則疏水能力越強(qiáng)，越不容易發(fā)生水的自吸。通過該方法，可以快速判斷頁巖地層中的混合潤(rùn)濕孔是否會(huì)發(fā)生水的自吸。因此，對(duì)于有機(jī)質(zhì)含量低、接觸角小的頁巖氣儲(chǔ)層進(jìn)行水力壓裂，應(yīng)考慮水自吸進(jìn)入頁巖后引起水鎖、水敏等問題，從而指導(dǎo)我們采取合理的措施來防止這些對(duì)生產(chǎn)不利的因素。

本文還存在一些局限性，比如考慮了2種不同潤(rùn)濕性的物質(zhì)組成的模型，而實(shí)際巖層中存在多種不同潤(rùn)濕性的物質(zhì)；為了分析簡(jiǎn)便，文中將裂縫系統(tǒng)簡(jiǎn)化為平板模型，而地下巖層孔隙和裂縫形狀比平板模型復(fù)雜得多，但這并不影響該模型的理論指導(dǎo)意義。

4 結(jié) 論

a.頁巖儲(chǔ)層為非均質(zhì)潤(rùn)濕，存在無機(jī)親水孔隙、有機(jī)疏水孔隙和混合潤(rùn)濕孔隙。根據(jù)裂縫系統(tǒng)的形狀和礦物的分布，可以將裂縫系統(tǒng)簡(jiǎn)化為混合潤(rùn)濕平板模型，實(shí)現(xiàn)了裂縫系統(tǒng)潤(rùn)濕性和自吸特征的定量分析。

b.混合潤(rùn)濕平板模型自吸速度的大小受潤(rùn)濕性、界面張力、縫寬、黏度等多種因素影響。界面張力與縫寬越大，自吸速度越快；黏度越小，自吸速度越慢；潤(rùn)濕性越強(qiáng)，自吸速度越快。

c.對(duì)于混合潤(rùn)濕平板模型，氣-水的自吸方向只由親水壁和疏水壁的接觸角之和決定。當(dāng)θ1+θ2<π時(shí)，v>0，此時(shí)只可自發(fā)吸水排氣；當(dāng)θ1+θ2>π時(shí)，v<0，此時(shí)只可自發(fā)吸氣排水；當(dāng)θ1+θ2=π時(shí)，v=0，氣、水均無法自吸。