薛東杰趙艾博劉奎昌侯孟冬付艷艷辛翠徐顏卓

中國礦業(yè)大學(北京)力學與建筑工程學院,北京 100083

“滲”這個字意味著液相或氣相在固相的空隙中發(fā)生空間移動,尤其是多孔介質(zhì)中,固相并非連續(xù)分布,不同尺度的孔裂隙互相交錯并連通,為流相的滲透提供了幾何空間。滲透、滲流和突滲是流相發(fā)生轉(zhuǎn)移的不同層次,對應著安全和災害的轉(zhuǎn)換過程。自然界中滲透現(xiàn)象廣泛存在,但當科學家開始研究其規(guī)律和內(nèi)在本質(zhì)時才發(fā)現(xiàn),問題的層次性依然存在。如規(guī)律性的實驗展示了在層流特定條件下,某些物理量之間存在著可基于類似Darcy 定律描述的線性關系。但是當滲透向滲流甚至突滲行為轉(zhuǎn)變時,非線性帶來的挑戰(zhàn)困難得多,如Navier-Stokes 方程的理論求解成為歷史難題;突變時的臨界特征和臨界力學理論仍處于建立之中。這些難點的解決需要觀察視角的轉(zhuǎn)換?！皾B”意味著一種相互作用,液相在固相的空隙空間中發(fā)生轉(zhuǎn)移,若在幾何模型上將固相剔除,可獲得空隙幾何結(jié)構,理論上將空隙結(jié)構作為邊界條件?；贜avier-Stokes 方程建立流體力學模型即可獲得流動規(guī)律的求解,這種觀察視角是不考慮固體力的相互作用,而是僅僅從幾何約束的條件下考慮的。可見,相互作用尤其是力的相互作用在固相-液相間的描述常常是被忽略,而從幾何約束的角度來體現(xiàn)相互作用,本質(zhì)上是一種幾何約束。

多孔介質(zhì)的空隙結(jié)構幾何復雜性遠超人們想象,從幾何角度來破解理論難題是科學家關注的焦點;而工程師往往關注的是宏觀滲流規(guī)律,對于滲透、滲流和突滲行為轉(zhuǎn)換過程中的強非線性只能依賴幾何角度求解,裂隙網(wǎng)絡幾何復雜且在應力作用下會發(fā)生連通,形成新的幾何網(wǎng)絡,問題變得更加復雜。迂曲度是描述裂隙網(wǎng)絡彎折效應的一個結(jié)構變量,其并不直接描述裂隙幾何,而是描述其復雜性。彎折性反映了流體通過固相路程的難易程度,如在土壤中發(fā)生滲透遠比巖石中容易[1],相應的迂曲度也要小很多。在深地煤炭開采工程[2-3]、地下鹽穴儲氣庫工程[4]及花崗巖封存高放核廢料工程[5]中滲透性等評價也是工程重中之重。

在裂隙網(wǎng)絡幾何描述問題解決之前,力學建模和試驗也是重要的研究手段。例如,基于分數(shù)階理論建立滲流模型用以描述非飽和滲流問題[6];模擬煤樣在三軸加載過程中的流固耦合[7];揭示巖石破碎后的非線性滲流行為[8]。這些模型的建立的依據(jù)和試驗過程的準確解析都需要對裂隙幾何進行精準描述。但事實是,精準測量低滲巖體內(nèi)部的空隙結(jié)構是十分困難的,設備精度和視野總是存在不可協(xié)調(diào)的矛盾,高精度設備無法測量大尺度樣品;大尺度樣品又難以探測各種尺度的裂隙[9-10]。若結(jié)構不明或不準確,在定量計算迂曲度時就無法直接依賴試驗所測的裂隙幾何數(shù)據(jù),因此迂曲度的理論計算仍存在著巨大挑戰(zhàn)。通過試驗也可間接得知迂曲度的一些特性,根據(jù)Carman-Kozeny 經(jīng)驗公式,多孔介質(zhì)滲透率與毛細管迂曲度近似呈反比關系[11];孔隙度迂曲度存在正向關系[12-16]。因此,采用滲透率相關原理間接測試迂曲度也是主要方法,即采用理論與實驗結(jié)合的方法可有效避免直接采用幾何測量方法的誤差。Purcell 較早提出用汞開展毛管壓力與滲透率測試[17],但其忽略了迂曲度的影響,導致滲透率誤差較大[18-20];Boundreau 對有限顆粒構成的迂曲度開展了建模研究[21];Sen 等認為孔隙結(jié)構空間分布具有自相似性,即迂曲度可利用分形理論來研究[22]。

對于流體,迂曲度和渦流等現(xiàn)象的關聯(lián)也是研究熱點[23-24],但迂曲度的定義來源及演化十分復雜。早期對迂曲度的認識僅是對滲透率的誤差修正,隨著認識的深入,賦予了迂曲度明確的物理意義,即反映孔隙結(jié)構流程的彎曲程度,才稟賦了其內(nèi)在物理意義[25],國內(nèi)也有文獻提出直接用孔隙度計算迂曲度[26],認為迂曲度曲線具有分維特征[27],提出了分形毛細管力模型[28]。

筆者基于毛細管模型,根據(jù)Hagen-Poiseuille 公式與達西定律建立普適的迂曲度表達式,并針對低滲介質(zhì)提出了一種適用于毛細管力特征明顯的迂曲度表達式,建立了迂曲度分形維數(shù)。針對分叉現(xiàn)象,在滿足質(zhì)量守恒前提下,分析了分叉孔道母孔子孔幾何關系,并指出分叉現(xiàn)象的產(chǎn)生并非完全隨機,而是在滿足優(yōu)化能量前提下的理想分配,驗證了Murray 定律,也適用于巖石孔裂紋系統(tǒng)。本文以鹽巖為案例,驗證了表達式的適用性與可靠性。研究結(jié)果可為求解多孔介質(zhì)材料迂曲度提供理論與實驗參考。

1 迂曲度建模及分形維數(shù)計算

孔隙通道包括非有效連通部分(盲孔)及有效連通部分。連通性描述是揭示流體流通的重要基礎。通道方向分布具有不確定性,截面形狀與尺寸大小位置分布具有隨機性,在水力、溫度等外界因素影響下還可能存在著開閉狀態(tài)的轉(zhuǎn)換,這些復雜因素都給認識孔隙結(jié)構空間分布帶來挑戰(zhàn),因此從統(tǒng)計理論上分析更具有現(xiàn)實意義?；诖?本文將孔隙通道假設為毛細管模型[25](圖1)。

圖1 毛細管模型Fig.1 Capillary model

如圖1所示,假設毛細管道在截面xz上均勻分布,毛細管孔道形狀相同,孔道中流體單相、飽和,其運動規(guī)律符合達西定律。微元體截面xz面積為A,單位面積分布毛細管根數(shù)為n,每根毛細管彎曲度為τ,孔道截面形狀為圓形,半徑為r,單根孔道體積為V,單位橫截面密度為N。迂曲度定義較多[29-31],Le表示流動路徑長度,定義迂曲度[29]

一種方法是根據(jù)Hagen-Poiseuille 公式[32],則通道流量q為

式中,Δp為壓差;μ為黏滯系數(shù)。

根據(jù)達西定律:

式中,K為滲透率。

定義孔隙度[33]

聯(lián)立式(2)、式(3)與式(4)得

引入水力學半徑平均值Rh,則

式中,S為流體流過總截面積;Z為總濕潤周長。

根據(jù)經(jīng)驗[34],存在如下關系:

式中,Dp為粒子平均半徑。

結(jié)合式(5)、式(6)與式(7),得迂曲度

假設系統(tǒng)由N個球形毛細管粒子構成,則

對應的,迂曲度可定義為

另一種方法是利用毛細管壓力公式[25],較大孔徑孔隙通道并不存在典型毛細現(xiàn)象,其適用范圍需要界定,則單根毛管壓力[31]

式中,σ為界面張力;θ為接觸角。

結(jié)合式(2)得

聯(lián)立達式(3)得

考慮飽和度s,則式(13)可改為

則迂曲度重新定義為

假設所有毛細管壓力p相同,則

利用式(15)和式(16)即可測定迂曲度。Wheatcraft 和Tyler 建立了多孔介質(zhì)迂曲流管特性的分形標度關系[34]:

則

式中,λmin為毛細管孔徑的最小值;DT為迂曲度分形維數(shù)。

對于線性孔道滲流,1≤DT≤2。

聯(lián)立式(15)與式(18)得

式(19)說明利用壓汞實驗可以確定迂曲孔道分形維數(shù)。定義分形影響系數(shù):

當迂曲度分形維數(shù)DT確定時[圖2(a)],即毛細管道彎折程度一定時,分形影響系數(shù)與長徑比近似成線性關系。如圖2(b)所示,當毛細管孔道長徑比確定時,分形影響系數(shù)與迂曲度分形維數(shù)呈非線性關系,DT與L/r越大影響越明顯。如當L/r≤5 時,變化不大,此時可以忽略DT的影響,由此可知當考慮DT對的影響時,需要考慮DT是否在長流程環(huán)境下,當L/r≤5 的短流程環(huán)境時,孔道長度和直徑在一個數(shù)量級上,孔道彎折性不顯著,彎折性對流體的阻滯作用不明顯;另一方面說明,粗糙性特征至少精度上要高于流程定義,并且應考慮粗糙在流程上的非線性積分。

圖2 長徑比和迂曲度分維與分形影響系數(shù)的關系Fig.2 Fractal coefficient influenced by the aspect ratio and fractal dimension of tortuosity

2 分叉孔道迂曲度特征分析

通常在描述孔、裂隙通道時會隱含假設,即通道始終是1 條,不存在1 條通道連接2 條孔道甚至更多。分叉現(xiàn)象總是被刻意忽略或者簡單地統(tǒng)計長度。這是因為分叉會引起更多的非線性特征,如分流現(xiàn)象;分叉孔道給分維的計算帶來挑戰(zhàn)。

裂紋分叉是普遍存在的,應引起足夠的重視。自然界中材料分布與裂紋生長過程并非完全隨機現(xiàn)象,最優(yōu)化能量的分配與傳播控制著裂紋生長[35-36]。合理的能量分配最終導致裂紋的半徑與長度生長相互影響,形成具有規(guī)律性的網(wǎng)絡系統(tǒng)。天然生成的孔道具有同樣的規(guī)律,雖然能量形式釋放機制迥異,但都遵循最小化能量的原理。

利用科林斯定義僅統(tǒng)計分叉毛細管長度意義不大,有較多的局限性,首先無法表征長程單條孔道粗糙度及迂曲度的特征,其次忽略了迂曲度對孔道滲透特性的界定。基于迂曲度定義,為了更好地反映滲透特征,有必要建立考慮分叉孔道特性的迂曲度模型。

母孔設為分叉前的單一孔道,子孔為由母孔分叉而來的多條孔道?？紤]最簡單分叉(1 分2)(圖3),考慮流體從總孔道(標號為0,下同)分叉分別流向叉道1 與叉道2。總孔道長度為L0,迂曲度為τ0,等效半徑為r0,流量q0,對應的分叉1 物理意義:L1、τ1、r1與q1;分叉2 物理意義:L2、τ2、r2與q2。

圖3 分叉示意圖Fig.3 Bifurcate path of capillary

根據(jù)式(2)可得流體持續(xù)流動所需要的能量:

孔道表面粗糙性引起的流體能量耗散:

式中,η為表面粗糙性系數(shù)。

總的能量:

根據(jù)最小能量原理,對式(23)求微分,得

式中,ξ為表面粗糙性系數(shù)。

可見,流量與孔徑的立方成正比,這與平行板立方體定律相似。平行板定律是由實驗總結(jié)的,而式(24)是基于Hagen-Poiseuille 公式理論推導的結(jié)果。

根據(jù)質(zhì)量守恒q0=q1+q2,對母孔0、子孔1、2同時使用式(24)得

母孔孔徑的立方等于子孔孔徑的立方和,即Murray 定律[37-39],說明Murray 定律適用于巖石孔裂隙系統(tǒng)。

定義流阻:

母孔與子孔串聯(lián),兩子孔并聯(lián),則總流阻為

考慮分叉結(jié)構的體積,即流量的總體積:

定義Lagrangian 函數(shù):

分別對r0、r1與r2求極值:

式(30)化簡后同式(25)。再次驗證了巖石裂紋分叉孔道滿足Murray 定律。

將分叉孔道整體考慮,則總能量:

等效迂曲度:

忽略粗糙性,則

3 低滲鹽巖毛細孔道迂曲度計算

鹽巖是天然氣儲藏、核廢料等有害物質(zhì)封存的理想地下場所,低滲現(xiàn)象是其致密的孔隙結(jié)構導致的。樣品源于湖北應城鹽礦,主要成分為石鹽,含少量石英。加工成幾何尺寸為直徑25 mm、高30 mm 試樣,進行高壓壓汞實驗。實驗前利用工業(yè)CT 對孔隙掃描,精度為30 μm,難見孔隙分布,可見鹽巖孔隙尺寸十分微小(圖4)。設計壓力范圍0.1～226 MPa,對應孔隙尺寸為7.5 μm～3.2 nm,連續(xù)加壓,進汞速度控制7 MPa/min,設備精度為0.01%。

圖4 鹽巖取樣、加工及CT 掃描Fig.4 Drilling and processed sample for CT scanning

圖5 為所得到的毛細管壓力曲線,初始壓汞從0.008 MPa 開始,迅速增長至30 MPa,此時對應的汞飽和度從0% 到接近10% ;然后壓汞曲線變得平緩,壓力從30 MPa 到最大值200 MPa,對應汞飽和度從10% 到70%。本樣品30 MPa 是汞飽和度出現(xiàn)急速變化的拐點,將壓汞曲線分為兩個部分。30 MPa 之前,汞進入的都是大孔隙,但這種空隙結(jié)構并不多對應的汞飽和度不高。30 MPa 之后,壓力突破了拐點,毛細管孔開始連通,汞飽和度迅速增長,但高壓下鹽巖可能會出現(xiàn)損傷,所測飽和度并非原始狀態(tài)的孔隙特征。之后從200 MPa 開始急速回汞到40 MPa,但汞飽和度維持在70% 幾乎沒有變化;之后回汞到0.015 MPa 對應的汞飽和度下降至45%?？梢?在不同的壓力下孔道的有效連通對汞飽和度影響較大。側(cè)面驗證了鹽巖結(jié)構的復雜性和致密性。根據(jù)實驗測定本樣品鹽巖孔隙度φ=6.5%,滲透率K=0.001md,取汞的界面張力σ=0.48 N/m,接觸角θ=140°,利用式(15)結(jié)合圖5 數(shù)據(jù)積分,可計算得到鹽巖迂曲度τ=58.32。常規(guī)計算方法得到的迂曲度普遍小于5,本方法測定的迂曲度是其十幾倍,更符合真實情況。鹽巖迂曲度較大,說明復雜分布孔道形成的滲透路徑曲折,最終造成流體通過困難,即鹽巖低滲的主要原因是低孔隙度和高迂曲度。

圖5 毛細管壓力曲線Fig.5 Curve of capillary pressure

4 結(jié)論

低滲孔隙介質(zhì)如鹽巖內(nèi)部的空隙結(jié)構幾何和拓撲精準表征在數(shù)學上是很大的挑戰(zhàn),多數(shù)工作關注了幾何信息,而忽略了拓撲連通等核心變量。用于描述空隙結(jié)構的幾何變量較多,但普適性的關鍵幾何或拓撲變量并未最終確定。滲透的本質(zhì)是連通,連通性描述的數(shù)學基礎是拓撲幾何理論,因此,有必要在傳統(tǒng)迂曲度幾何定義的基礎上開展拓撲幾何描述研究。本文嘗試從最簡單的分叉模型推導出迂曲度的理論表達式,但與圓形毛細管模型等相比,表達式仍差異較大。可見,迂曲度的理論推導實現(xiàn)普適性應用仍任重道遠。

(1)考慮流通孔道的分形特性,引入迂曲度分維,并推導了迂曲度分維確定的實驗公式,這是一種基于壓汞實驗求解分維的新方法。定義了分形影響系數(shù),當迂曲度分形維數(shù)DT確定時,分形影響系數(shù)ˉλ與長徑比近似成線性關系。當長徑比確定時,分形影響系數(shù)ˉλ與迂曲度分形維數(shù)DT呈非線性關系。

(2)針對分叉孔道,流量與孔徑立方成正比關系,在滿足質(zhì)量守恒前提下,分叉孔道母孔孔徑的立方等于子孔孔徑的立方和,從兩個方面驗證了巖石裂紋分叉孔道滿足Murray 定律,推導了分叉孔道等效迂曲度數(shù)學表達式。以鹽巖滲透為案例,低孔隙度、高迂曲度是鹽巖低滲的主要原因。