基于CT的黃土孔隙尺度優(yōu)先流特性

寧瑞浩， 冷艷秋*， 何芝遠(yuǎn)， 李澤坤， 馬哲

(1. 長(zhǎng)安大學(xué)地質(zhì)工程與測(cè)繪學(xué)院，西安 710054； 2.西部礦產(chǎn)資源與地質(zhì)環(huán)境教育部重點(diǎn)實(shí)驗(yàn)室，西安 710054)

黃土在中國(guó)主要分布在山西、陜西、甘肅、寧夏的黃土高原一帶，是一種疏松多孔、豎向節(jié)理、裂隙發(fā)育，具有很強(qiáng)的濕陷性和水敏性的結(jié)構(gòu)性土；因其特殊的濕陷性和水敏性，水便成為了黃土地區(qū)誘發(fā)地質(zhì)災(zāi)害的關(guān)鍵因素。常規(guī)的滲水試驗(yàn)表明降雨或者灌溉的入滲深度只有2～3 m[1]，但是黑方臺(tái)和涇陽(yáng)等地的滑坡實(shí)例卻顯示水運(yùn)移到坡體以下數(shù)十米，優(yōu)先流便是這類地質(zhì)現(xiàn)象的合理解釋，因此探索黃土孔隙、裂隙中優(yōu)先流的滲流特性，對(duì)防控黃土地質(zhì)災(zāi)害具有重大意義。

優(yōu)勢(shì)流又稱為優(yōu)先流，是一種非平衡流，主要特征為部分水流繞過(guò)大部分土體基質(zhì)沿著優(yōu)勢(shì)通道迅速入滲到土體深部[2]。目前對(duì)優(yōu)先流的研究主要集中在兩個(gè)尺度，即宏觀尺度和孔隙尺度[3]，Watanabe等[4]通過(guò)構(gòu)建局部非孔洞模型，研究了碳酸鹽巖的孔隙度與滲透率和內(nèi)部?jī)?yōu)勢(shì)流的關(guān)系，為碳酸鹽巖儲(chǔ)層的開(kāi)發(fā)利用提供了依據(jù)。趙寬耀等[5]通過(guò)物探探測(cè)與數(shù)值模擬相結(jié)合的方法，分析了黃土邊坡在優(yōu)勢(shì)流與基質(zhì)流共同作用下不同灌溉強(qiáng)度水的滲流過(guò)程。陳思婕等[6]將雙重滲透模型與無(wú)限邊坡穩(wěn)定分析法相耦合，根據(jù)土壤含水量及孔隙水壓力觀測(cè)數(shù)據(jù)，模擬了強(qiáng)降雨條件下滑體內(nèi)的水動(dòng)力過(guò)程，以此分析了優(yōu)先流對(duì)滑坡觸發(fā)機(jī)理的影響。然而這些研究都是集中在宏觀層面的。潘網(wǎng)生等[7]認(rèn)為定量化描述黃土微細(xì)觀孔隙、裂隙結(jié)構(gòu)模型,建立從微細(xì)觀到宏觀不同尺度下研究黃土滑坡機(jī)理的系統(tǒng)理論方法，對(duì)認(rèn)識(shí)優(yōu)先流對(duì)滑坡的觸發(fā)機(jī)理尤為重要，僅僅開(kāi)展宏觀尺度的研究并不能揭示黃土中優(yōu)先流的本質(zhì)規(guī)律。

研究孔隙尺度內(nèi)的優(yōu)先流特性首先要獲取真實(shí)的黃土內(nèi)部孔隙結(jié)構(gòu)，在眾多方法中計(jì)算機(jī)斷層掃描(computed tomography,CT)技術(shù)作為一種無(wú)損檢測(cè)技術(shù)正被廣泛用于巖土材料的孔隙定量化表征。Li等[8]利用CT掃描技術(shù)定量化研究了馬蘭黃土的大孔隙結(jié)構(gòu)特征，結(jié)果表明馬蘭黃土在垂直和水平方向上的孔隙形狀、聯(lián)通特征有很大差異，這與室內(nèi)試驗(yàn)與數(shù)值模擬結(jié)果較為吻合。王超等[9]采用微米CT建立了黃土巖的孔隙網(wǎng)絡(luò)模型，研究了黃土巖中礦物成分、形狀對(duì)黃土巖的滲透率各向異性的影響。Wang等[10]通過(guò)CT掃描技術(shù)研究了黃土區(qū)露天煤礦排土場(chǎng)重構(gòu)土壤孔隙三維分布的多重分形表征，為礦區(qū)土壤重構(gòu)提供了依據(jù)。上述成果說(shuō)明利用CT技術(shù)獲取黃土內(nèi)部真實(shí)孔隙結(jié)構(gòu)并進(jìn)行定量化的研究是可行的。

目前已有眾多學(xué)者在真實(shí)孔隙結(jié)構(gòu)的基礎(chǔ)上開(kāi)展優(yōu)先流滲流特性研究，Li等[11]、Larsbo等[12]在孔隙尺度上研究了孔隙與孔隙網(wǎng)絡(luò)特征(平均直徑、配位數(shù)等)對(duì)土體內(nèi)部?jī)?yōu)勢(shì)流滲流規(guī)律的影響。Kim等[13]通過(guò)數(shù)值模擬的方法研究了孔隙尺度內(nèi)河床附近優(yōu)勢(shì)流及湍流對(duì)溶質(zhì)異常輸運(yùn)特性的影響。Zhang等[14]采用低場(chǎng)核磁共振技術(shù)對(duì)納米顆粒在非均質(zhì)多孔介質(zhì)中的運(yùn)移方式實(shí)時(shí)監(jiān)測(cè)，觀察到了明顯的優(yōu)勢(shì)流現(xiàn)象。魯拓等[15]研究了基于多種分形模型，計(jì)算了不同尺度范圍孔隙的分形維數(shù)，并討論了二者之間的聯(lián)系。Sarah等[16]利用高分辨率的微米CT研究了土壤內(nèi)部的細(xì)微觀結(jié)構(gòu)并利用貝葉斯統(tǒng)計(jì)方法探討了土壤的微觀特征與其飽和導(dǎo)水率等參數(shù)的相關(guān)性。Li等[17]利用AVIZO軟件研究了馬蘭黃土的大孔隙特征并基于形狀因子對(duì)其進(jìn)行分類，為研究孔隙尺度內(nèi)優(yōu)勢(shì)流特征提供了依據(jù)。

前人所做研究采用CT掃描的分辨率多數(shù)大于50 μm，雷祥義[18]根據(jù)壓汞實(shí)驗(yàn)的結(jié)果將大孔隙劃分為孔隙半徑為16～250 μm，中孔隙為孔隙半徑為4～16 μm，采用50 μm以上的分辨率顯然并不能滿足精度需求，因此在提取黃土內(nèi)部真實(shí)孔隙結(jié)構(gòu)的基礎(chǔ)上開(kāi)展微細(xì)觀尺度的優(yōu)先流滲流特性研究。采用高精度的微米CT，掃描分辨率為4 μm，以此構(gòu)建真實(shí)孔隙結(jié)構(gòu)，研究延安黃土孔隙尺度內(nèi)的優(yōu)先流滲流特性。

1 實(shí)驗(yàn)材料與方法

1.1 實(shí)驗(yàn)材料

實(shí)驗(yàn)所用土樣取自陜西延安邊坡坡腳，取樣深度約為75 m，所取原狀土樣為圓柱體，直徑和高度分別為25 cm和50 cm，在密封和標(biāo)記后運(yùn)至實(shí)驗(yàn)室，運(yùn)輸過(guò)程中采取合理的避震措施從而避免擾動(dòng)。通過(guò)對(duì)所取試樣隨機(jī)測(cè)量，研究區(qū)黃土試樣的基本物理性質(zhì)見(jiàn)表1。圖1為試驗(yàn)所采用黃土試樣的粒度分布曲線，其中粉粒(2～75 μm)含量81.72%，黏粒(<2 μm)含量為6.65%。結(jié)合粒度分布曲線及塑性指數(shù)，按照文獻(xiàn)[19]對(duì)黃土的分類，研究區(qū)黃土屬于黃土質(zhì)砂黏土。

表1 土樣基本物理性質(zhì)Table 1 Basic physical properties of soil samples

圖1 粒度分布曲線Fig.1 Particle size distribution curve

1.2 微米CT掃描

實(shí)驗(yàn)所采用CT掃描設(shè)備為Phoenix v|tome|x s工業(yè)CT(圖2)，該CT系統(tǒng)組合有180 kV/15 W 高功率nanofocus X射線管和240 kV/320 W 的微焦點(diǎn)管，是一種靈活可靠的掃描分析工具。本次測(cè)試試樣為圓柱體，直徑7 mm，高15 mm(圖3)，掃描分辨率為4 μm，共獲取X、Y、Z方向切片14 200張，每張二維切片像素點(diǎn)的灰度值代表著該點(diǎn)的X射線衰減系數(shù)，密度越大衰減系數(shù)越大，借助這一原理再結(jié)合適當(dāng)?shù)姆指罘椒?，三維重建等手段，便可以把試樣內(nèi)部的孔隙裂隙識(shí)別重建出來(lái)以用于二維、三位參數(shù)的提取，具體技術(shù)路線見(jiàn)圖4。

圖2 CT掃描設(shè)備Fig.2 CT scanning equipment

圖3 掃描土樣Fig.3 Scanned soil sample

圖4 技術(shù)流程圖Fig.4 Technical flow chart

1.3 室內(nèi)滲透實(shí)驗(yàn)

為研究黃土微細(xì)觀與宏細(xì)觀滲流之間的聯(lián)系并驗(yàn)證計(jì)算結(jié)果可靠性，進(jìn)行室內(nèi)滲透實(shí)驗(yàn)與數(shù)值計(jì)算結(jié)果進(jìn)行對(duì)比。本次室內(nèi)變水頭滲透實(shí)驗(yàn)采用土樣為延安黃土原狀土樣與重塑土樣(圖5)，各自設(shè)置5組平行試驗(yàn)，重塑土樣采用分層靜壓法制備，將原狀黃土進(jìn)行碾碎、過(guò)篩后放入涂好凡士林的壓樣筒內(nèi)，控制重塑土樣干密度與原狀樣相等，分五層壓實(shí)，層間接觸部位進(jìn)行刮毛，每層壓實(shí)時(shí)間不少于45 min。試樣為標(biāo)準(zhǔn)尺寸，直徑、高度分別為61.8、40 mm，在真空缸中飽和后進(jìn)行滲透實(shí)驗(yàn)。滲透實(shí)驗(yàn)流程按照《土工試驗(yàn)方法標(biāo)準(zhǔn)》(GB/T 50123—2019)規(guī)定，待水頭管中水頭穩(wěn)定后在起始時(shí)刻標(biāo)定水頭高度，每隔20～40 s測(cè)定水頭和時(shí)間的變化，如此連續(xù)測(cè)記10次后再使水頭回升至適當(dāng)高度，重復(fù)試驗(yàn)5次，所得結(jié)果如表2所示。

圖5 實(shí)驗(yàn)土樣Fig.5 Experimental soil sample

表2 室內(nèi)滲透實(shí)驗(yàn)所得滲透系數(shù)Table 2 Permeability coefficient obtained from indoor permeability test

2 CT圖像預(yù)處理與三維重建

2.1 建模軟件選取

在目前常用的多孔介質(zhì)二、三維孔隙結(jié)構(gòu)及顆粒定量化分析平臺(tái)中，AVIZO軟件可以滿足圖像降噪處理-圖像二值分割-三維體重建-數(shù)據(jù)分析-滲流模擬-絕對(duì)滲透率計(jì)算等各方面的需求，在科學(xué)和工業(yè)領(lǐng)域、生命科學(xué)和生物醫(yī)學(xué)、材料科學(xué)、土壤學(xué)、巖土領(lǐng)域有著廣泛的應(yīng)用，能夠滿足本研究的相關(guān)需求，因此選用AVIZO平臺(tái)進(jìn)行優(yōu)先流的相關(guān)研究。

2.2 圖像預(yù)處理

在CT圖像的采集和獲取過(guò)程中，會(huì)因?yàn)椴杉O(shè)備和儀器的少量干擾和周?chē)h(huán)境因素的影響，使得采集的圖像或多或少的存在部分噪聲，為了進(jìn)一步提高掃描圖像的質(zhì)量，突出圖像中孔隙的邊緣輪廓特征，以便能夠準(zhǔn)確分割目標(biāo)孔隙和背景，提取出感興趣的孔隙部分，通常采用圖像噪聲濾波方法和邊緣增強(qiáng)進(jìn)行預(yù)處理。

目前常用的濾波方法有高斯濾波、中值濾波、非局部均值濾波法等，圖6分別為采用上述濾波方法的預(yù)處理結(jié)果。由圖6可知，非局部均值濾波法幾乎能完全去除圖像中的高斯噪聲、脈沖噪聲等，同時(shí)還能夠很好地保留圖像細(xì)節(jié)，很大程度地減輕邊緣模糊[圖6(d)]；其他方法雖然能夠消除噪聲但是不能很好的保留圖片信息，致使邊緣模糊，因此選用非局部均值濾波法對(duì)所獲取的CT圖像進(jìn)行降噪處理。非局部均值濾波法的基本思想為根據(jù)圖像的自相似性來(lái)計(jì)算鄰域像素權(quán)重，通過(guò)加權(quán)的形式將最相近的幾個(gè)像素塊中的中心點(diǎn)結(jié)合起來(lái)估計(jì)真實(shí)值，缺點(diǎn)是對(duì)計(jì)算機(jī)的計(jì)算能力要求較高；雖然采用的濾波方法在有效去除噪聲的同時(shí)能夠很好地保留圖像細(xì)節(jié)，但仍然存在一定的邊緣模糊現(xiàn)象，即黃土固相和孔隙之間交界處的過(guò)渡帶[圖6(e)]。采用AVIZO軟件中的unsharp mask模塊對(duì)圖像進(jìn)行邊緣增強(qiáng)處理[圖6(f)]，消減固相和孔隙之間交界處的過(guò)渡帶，大大增加分割精度。

2.3 圖像分割

在黃土孔隙結(jié)構(gòu)的三維重建及要素提取研究中，針對(duì)黃土CT掃描圖像的分割技術(shù)是整個(gè)研究的核心與關(guān)鍵點(diǎn)之一，其分割結(jié)果直接影響到后續(xù)孔隙結(jié)構(gòu)三維模型重建的精確性。到目前為止還沒(méi)有任何一種完備的分割方法能夠適用各種應(yīng)用環(huán)境，常用的圖像分割方法可以劃分為:基于閾值的分割方法、基于區(qū)域的分割方法、基于邊緣的分割方法以及基于某些特定理論的分割方法等。本研究采用基于閾值的分割方法。閾值分割法的基本原理為給定一個(gè)恰當(dāng)?shù)幕叶乳撝?，圖像某點(diǎn)處的灰度值大于該閾值認(rèn)為是土顆粒，小于該閾值則認(rèn)為是孔隙，因此閾值的選取尤為關(guān)鍵，過(guò)大會(huì)導(dǎo)致過(guò)分割，過(guò)小則會(huì)導(dǎo)致欠分割；常用的閾值分割方法有最大熵法、Ostu法、矩量法等，對(duì)同一張灰度圖像采用不同的分割方法的局部放大圖如圖7所示，最大熵法與矩量法確定的閾值偏大，存在明顯的過(guò)分割現(xiàn)象，Ostu分割法的分割效果與實(shí)際情況最為接近，然而還存在輕微的過(guò)分割現(xiàn)象，需手動(dòng)調(diào)試閾值?；诖?，分割閾值采用Ostu法確定參考閾值，在AVIZO軟件中通過(guò)手動(dòng)調(diào)整確定最佳閾值[圖7(d)]，對(duì)部分難以識(shí)別小孔隙、微小裂隙通過(guò)AVIZO軟件中的Top-hat算法進(jìn)行識(shí)別[圖7(e)]，最終取二者并集[圖7(f)]作為最終孔隙分割結(jié)果。經(jīng)計(jì)算分割后的微觀孔隙度約為19%，小于實(shí)際孔隙度，這是由于小于4 μm的孔隙無(wú)法被識(shí)別出來(lái)。

圖6 圖像預(yù)處理結(jié)果Fig.6 Image preprocessing results

2.4 孔隙結(jié)構(gòu)三維重建

2.4.1 代表性體元確定

目前，基于CT掃描圖像三維重構(gòu)的基本方法主要包括體繪制和面繪制兩種方法，采用可以避免圖像信息丟失，能夠保留數(shù)據(jù)完整性的體繪制方法對(duì)土樣中的孔隙網(wǎng)絡(luò)進(jìn)行三維重構(gòu)?？紤]到計(jì)算機(jī)的計(jì)算能力有限并且必須兼顧樣品具有夠的代表性，需確定代表性體積單元(representative elementary volume，REM)[20]。本文研究確定代表性體積單元的方法為在土樣任意一點(diǎn)處取一系列不同邊長(zhǎng)的立方體(圖8)，計(jì)算其孔隙度，當(dāng)孔隙度趨于穩(wěn)定時(shí)即認(rèn)為該邊長(zhǎng)的立方體為最小代表性體積單元。如圖9所示，在立方體邊長(zhǎng)達(dá)到1 250 μm時(shí)，代表性立方體孔隙度已達(dá)到穩(wěn)定值，即立方體邊長(zhǎng)達(dá)到1 250 μm以上時(shí)即可認(rèn)為該子體積具有代表性。

圖7 圖像分割結(jié)果Fig.7 Image segmentation results

圖8 REV選取示意圖Fig.8 Schematic diagram of REV selection

目前對(duì)優(yōu)先流的常見(jiàn)分類為：裂隙流、管流、指流、側(cè)向流，指流和側(cè)向流分別發(fā)育在上細(xì)下粗和含有較大阻水塊體的土體中，而本次所取土樣土質(zhì)較均一且不含巖石、樹(shù)根等阻水塊體?；诖耍趻呙柰翗又蟹謩e提取含裂隙(C1、C2，圖10)、含管狀通道(C3、C4，圖10)與不含優(yōu)勢(shì)通道的(C5、C6，圖10)2 000 μm×2 000 μm×2 000 μm的立方體，排除孤立孔隙，在各立方體連通孔隙基礎(chǔ)上對(duì)比分析裂隙流、管流的滲流特性。

圖9 孔隙度與立方體邊長(zhǎng)關(guān)系曲線Fig.9 Relation curve between porosity and cube size s

圖10 各代表性立方體Fig.10 Representative elementary volume

2.4.2 隙網(wǎng)絡(luò)模型的構(gòu)建

孔隙網(wǎng)絡(luò)模型(pore network model,PNM)的基本假定為一個(gè)相互連通的孔隙網(wǎng)絡(luò)有多個(gè)孔隙組成，不同孔隙之間通過(guò)孔道連接，孔喉是孔道中橫截面積最小的地方。在AVIZO中將孔隙網(wǎng)絡(luò)模型用球棍模型來(lái)表示，分別代表孔隙與孔隙之間的聯(lián)通關(guān)系，通過(guò)構(gòu)建孔隙網(wǎng)絡(luò)模型可以獲取孔隙之間的孔道長(zhǎng)度、孔喉面積以及配位數(shù)等關(guān)鍵信息，孔隙網(wǎng)絡(luò)模型如圖11所示。

圖11 孔隙網(wǎng)絡(luò)模型Fig.11 Pore network model

2.5 重建結(jié)果分析

2.5.1 孔隙與孔喉分布特征分析

根據(jù)所建立孔隙網(wǎng)絡(luò)模型所得數(shù)據(jù)，孔喉和孔隙等效直徑的分布范圍基本一致，對(duì)孔隙和孔喉的等效直徑以10 μm為間隔進(jìn)行直方圖統(tǒng)計(jì)，發(fā)現(xiàn)孔隙與孔喉近似服從對(duì)數(shù)正態(tài)分布。根據(jù)圖12，孔隙主要分布在孔徑為40～140 μm的區(qū)間內(nèi)，占總數(shù)目的91.83%，孔喉主要分布在孔徑為10～50 μm的區(qū)間內(nèi)，占總數(shù)目的85.07%。

圖12 孔隙與孔喉分布特征Fig.12 Distribution characteristics of pores and pore throats

2.5.2 孔隙連通特征分析

孔隙的連通性可以用配位數(shù)表征，一個(gè)孔隙的配位數(shù)是指與其他孔隙連接的孔喉的數(shù)量，配位數(shù)越多則代表著與周?chē)目紫哆B通性越好，根據(jù)構(gòu)建的孔隙網(wǎng)絡(luò)模型，對(duì)獲取的孔隙配位數(shù)進(jìn)行統(tǒng)計(jì)，發(fā)現(xiàn)配位數(shù)近似服從對(duì)數(shù)正態(tài)分布(圖13)，其值介于1～28之間；含優(yōu)勢(shì)通道的立方體中平均值為7.17，眾數(shù)為6，不含優(yōu)勢(shì)通道的立方體中平均值為6.15，眾數(shù)為5，這說(shuō)明在含優(yōu)勢(shì)通道的立方體中，孔隙的整體連通性要好于不含優(yōu)勢(shì)通道的立方體。

圖13 孔隙配位數(shù)分布特征Fig.13 Distribution characteristics of pore coordination number

3 滲流模擬

3.1 基本原理

Navier-Stokes方程是用于描述不可壓縮黏性流體的動(dòng)量守恒方程，其矢量形式為

(1)

式(1)中:ρ流體密度，kg/m3；V為速度矢量；P為流體壓力,kPa；f為單位體積流體所受外力，kN/m3；μ為流體動(dòng)力黏度,Pa·s。對(duì)V有

(2)

根據(jù)式(2)將式(1)展開(kāi)，得

(3)

一般條件下很難得到其精確解，考慮以下條件將Stokes方程簡(jiǎn)化：流體為不可壓縮流體；流體為牛頓流體，動(dòng)力黏度μ為常數(shù)；流體為穩(wěn)定流，流速不隨時(shí)間改變；流體流動(dòng)方式為層流；對(duì)黃土孔隙尺度內(nèi)流動(dòng)的水，上述條件均可滿足。

簡(jiǎn)化后的Stokes方程為

(4)

為保證方程在整個(gè)體積的有效性，將其轉(zhuǎn)化為體積平均形式，使方程在體積空間內(nèi)更加平滑，即

(5)

式(5)中:D為速度擾動(dòng)場(chǎng)張量；d壓力擾動(dòng)場(chǎng)向量；I單位張量。

通過(guò)求解系統(tǒng)體積v上速度擾動(dòng)場(chǎng)張量的平均值來(lái)求得滲透率張量k，即

(6)

在AVIZO中用有限體積法求解方程組，方程在交錯(cuò)網(wǎng)格上離散，假定體素為立方體，待求解壓力位于體素中心，而速度在體素面上分解，待求解量的時(shí)間導(dǎo)數(shù)趨于零時(shí)認(rèn)為收斂，誤差計(jì)算方法為

(7)

式(7)中:n為當(dāng)前迭代次數(shù)；?t為時(shí)間增量；c2為虛擬壓縮系數(shù)[21]。

對(duì)于本次計(jì)算，當(dāng)誤差小于10-4時(shí)即認(rèn)為收斂，當(dāng)求得絕對(duì)滲透率后，滲透系數(shù)計(jì)算公式為

(8)

3.2 模擬結(jié)果

模擬滲流方向?yàn)樽陨隙?，上下邊界為定水頭邊界，四周為不透水邊界，上下邊界的水頭差為0.002 m，即代表性立方體的高度，各立方體內(nèi)部流速分布如圖14～圖16所示。根據(jù)圖14(a)(e)～圖16(a)(e)，在含優(yōu)勢(shì)通道的代表性立方體中(C1、C2、C3、C4)高流速的滲流路徑沿著裂隙展布方向集中，裂隙、管狀通道內(nèi)的流線更加順直平滑，滲流壓力大，而其周?chē)拇罂紫秲?nèi)滲流路徑則較為分散，滲流壓力也更小，壓力分布更加不均勻；而在不含優(yōu)勢(shì)通道立方體中(C5、C6)，流線及壓力的分布較為均勻，流速也更慢，表現(xiàn)為均勻入滲,各代表性立方體滲透系數(shù)及絕對(duì)滲透率如表3所示。

表3 滲流模擬所得滲透系數(shù)Table 3 Permeability coefficients obtained from seepage simulations

圖14 含裂隙代表性體元模擬結(jié)果Fig.14 Representative voxel simulation results with fractures

綜合室內(nèi)滲透實(shí)驗(yàn)及數(shù)值計(jì)算結(jié)果可知重塑黃土由于內(nèi)部結(jié)構(gòu)被完全破壞導(dǎo)致孔隙之間連通性極差，滲透系數(shù)最小，要比原狀黃土小近一個(gè)數(shù)量級(jí)。不含優(yōu)勢(shì)通道的代表性立方體滲透系數(shù)與原狀黃土室內(nèi)滲透實(shí)驗(yàn)所得滲透系數(shù)在同一個(gè)數(shù)量級(jí)上，含優(yōu)勢(shì)通道的代表性立方體滲透系數(shù)要比原狀黃土大近一個(gè)數(shù)量級(jí)，這與前人所得結(jié)果一致[22]，說(shuō)明優(yōu)勢(shì)通道對(duì)黃土的滲透性起著主導(dǎo)作用，也間接證明了黃土微細(xì)觀尺度與宏細(xì)觀尺度滲流的統(tǒng)一性。

為詳細(xì)研究黃土內(nèi)部?jī)?yōu)先流的分布及流速特征，借助前述對(duì)黃土內(nèi)部孔隙進(jìn)行概化的孔隙網(wǎng)絡(luò)模型將孔喉按流量大小進(jìn)行展示[圖14(b)(e)～圖16(b)(e)]，通過(guò)孔喉截面面積及流量換算流速并分組統(tǒng)計(jì)，依據(jù)數(shù)值計(jì)算及室內(nèi)滲透實(shí)驗(yàn)所得滲透系數(shù)將孔喉按流速劃分為三部分，不流動(dòng)-低流速孔喉(T1)：流速小于1 μm/s，多位于孔隙連通性較差部位，流速與重塑黃土滲透系數(shù)相當(dāng)，處于同一數(shù)量級(jí)；均勻流速孔喉(T2)：流速介于1～100 μm/s，距離優(yōu)勢(shì)通道較遠(yuǎn)，但孔隙連通性一般部位，與原狀黃土及含優(yōu)勢(shì)通道的代表性立方體滲透系數(shù)處于同一數(shù)量級(jí)；優(yōu)先流發(fā)育孔喉(T3)：流速大于100 μm/s，多位于優(yōu)勢(shì)通道及孔隙連通性較好部位。

圖15 含管狀通道代表性體元模擬結(jié)果Fig.15 Representative voxel simulation results with pipes

對(duì)三個(gè)流速段的孔喉數(shù)量進(jìn)行統(tǒng)計(jì)，統(tǒng)計(jì)結(jié)果見(jiàn)圖17。三種立方體中均以均勻流速孔喉為主，其中無(wú)優(yōu)勢(shì)通道立方體中的T2占比最高，含裂隙立方體與含管狀通道立方體中的T2通道數(shù)量幾乎相同。將含裂隙立方體與含管狀通道立方體進(jìn)行對(duì)比分析發(fā)現(xiàn)前者T3占比為30%，其含量超過(guò)T1，同時(shí)也比后者T3占比高，因此含裂隙立方體的絕對(duì)滲透率比含管狀通道立方體大。將含裂隙與不含優(yōu)勢(shì)通道的代表性立方體孔喉流速分布進(jìn)行對(duì)比，發(fā)現(xiàn)裂隙發(fā)育的立方體中T1所占百分比明顯降低，T2數(shù)量約減少10%，而T3所占比例顯著增加，約增加15%，說(shuō)明裂隙對(duì)整體的滲透特性影響為均勻促進(jìn)型，即裂隙的存在會(huì)整體地提升立方體的滲透性；而管狀通道發(fā)育的立方體與不含優(yōu)勢(shì)通道的立方體相比T1占比約增加5%，T2數(shù)量約減少10%，T3所占比例增加幅度較裂隙發(fā)育的立方體小，但管狀通道內(nèi)的最大流速顯著高于裂隙內(nèi)部，這說(shuō)明管狀優(yōu)勢(shì)通道對(duì)整體的滲透特性影響為滲流集中型，即滲流路徑沿著管狀通道集中，對(duì)非優(yōu)勢(shì)通道在滲流過(guò)程中所發(fā)揮的作用有所削弱。在不含優(yōu)勢(shì)通道的代表性立方體中仍有相當(dāng)一部分孔喉流速大于100 μm/s，這說(shuō)明局部聯(lián)通性較好的大孔隙也可形成優(yōu)勢(shì)通道。對(duì)三種流速段孔喉直徑進(jìn)行統(tǒng)計(jì)發(fā)現(xiàn)T1部分孔喉平均直徑為11.2 μm, T2部分孔喉平均直徑為19.4 μm, T3部分孔喉平均直徑為29.1 μm,由此可見(jiàn)孔喉流速與孔喉直徑關(guān)系為正相關(guān)。

圖16 不含優(yōu)勢(shì)通道代表性體元模擬結(jié)果Fig.16 Representative voxel simulation results without dominant channels

圖17 各流速段孔喉百分比Fig.17 Percentage of orifice throat at each flow rate section

4 結(jié)論與展望

在對(duì)黃土內(nèi)部孔隙進(jìn)行三位重建的基礎(chǔ)上進(jìn)行滲流模擬，結(jié)合室內(nèi)滲透實(shí)驗(yàn)，可得出以下結(jié)論。

(1)基于真實(shí)孔隙結(jié)構(gòu)數(shù)值計(jì)算得到的滲透系數(shù)與室內(nèi)滲透實(shí)驗(yàn)所得原狀黃土的滲透系數(shù)在同一數(shù)量級(jí)上，證明了計(jì)算結(jié)果的可靠性。

(2)從壓力分布圖可知，裂隙中不僅有著更高的流速而且有著更高的壓力，這種壓力差會(huì)促進(jìn)聯(lián)通孔隙向著微裂隙轉(zhuǎn)變，即優(yōu)先流對(duì)黃土內(nèi)部的微觀孔隙結(jié)構(gòu)存在著一定的改造作用，但這種作用的大小依賴于邊界壓差。

(3)根據(jù)原狀黃土、重塑黃土滲透系數(shù)將孔喉按流速分組統(tǒng)計(jì)，發(fā)現(xiàn)裂隙與管狀通道對(duì)整體滲透特性的影響分別為均勻促進(jìn)型與滲流集中型。

更高的分辨率意味著更小的試樣尺寸及更大的計(jì)算機(jī)能力要求，隨著未來(lái)CT掃描技術(shù)與計(jì)算機(jī)計(jì)算能力的提升，可以進(jìn)行更大尺度更高分辨率的研究，從而更好的研究黃土微觀結(jié)構(gòu)和宏、細(xì)觀結(jié)構(gòu)之間的聯(lián)系。