郭亞永，馮興梅

（1.北京市勘察設(shè)計(jì)研究院有限公司，北京 100038；2.河北瑞三元環(huán)境科技有限公司，河北 石家莊 050000）

0 引言

深入分析和研究非飽和土體災(zāi)害的成因是一個(gè)非常艱巨的任務(wù)，而且通過實(shí)驗(yàn)往往很難準(zhǔn)確地分析整個(gè)過程，這也是這項(xiàng)工作很難開展的原因。針對(duì)這種情況，本文應(yīng)用基于浸沒邊界法的三相求解器來(lái)描述非飽和土體的滲流情況，同時(shí)應(yīng)用離散元模型對(duì)不同材料非飽和滲流中顆粒的作用力進(jìn)行分析。

目前有限元法作為滲流應(yīng)力場(chǎng)研究的方法被廣泛應(yīng)用［1-3］，然而表層土邊坡土壤壓實(shí)程度小，顆粒之間雖然接觸但仍存在小孔隙，處于離散狀態(tài)，因而采用有限元法很難分析顆粒間的滲流狀態(tài)。離散元法可以模擬與實(shí)際的土體顆粒相似的孔隙模型，反映流體通過顆粒孔隙時(shí)的狀態(tài)，更深入地了解滲流機(jī)理。隨著離散元法在許多方面的成功應(yīng)用，采用CFD-DEM方法模擬滲流導(dǎo)致的滑坡或者泥石流運(yùn)動(dòng)已逐漸成為主導(dǎo)方向［4-8］。目前一般通過CFD的流場(chǎng)模擬和DEM的多孔模型來(lái)對(duì)土體滲流進(jìn)行數(shù)值分析。對(duì)于耦合的CFD-DEM，通過求解CFD中的雷諾平均N-S方程來(lái)求解流體相的行為。Jing等［9］擴(kuò)展了CFD-DEM與體積分?jǐn)?shù)（VOF）的耦合方法，用于液固相之間體積替換，并根據(jù)自由流體動(dòng)力學(xué)提出了多孔球法，使流體網(wǎng)格和顆粒粒徑比的范圍變得更大。

筆者基于對(duì)降雨非飽和滲流的研究，從土體非飽和降雨不同水頭壓力出發(fā)，設(shè)計(jì)可人工電腦組合控制降雨的設(shè)備，在測(cè)試降雨入滲量的同時(shí)，通過不同水頭高度降雨來(lái)測(cè)試土體內(nèi)水流壓力的變化規(guī)律，并結(jié)合數(shù)值模擬來(lái)分析降雨入滲過程定水頭和變水頭階段滲透系數(shù)的變化特點(diǎn)；此外，借助resolve CFD-DEM的耦合來(lái)實(shí)現(xiàn)雨水充滿土體顆粒區(qū)域的流場(chǎng)分析，用離散元軟件的二次開發(fā)功能研究了在滲流過程中顆粒之間的接觸力與雨水滲流的變化狀態(tài)，揭示了滲流過程土體顆粒壓力與速度的變化特征。

1 實(shí)驗(yàn)材料以及裝置設(shè)計(jì)

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

本實(shí)驗(yàn)選用的3種材料分別為鐵沖流域原狀土、鐵沖流域土顆粒與水泥按照一定比例制成的水泥土、按照相似配比組成的突水突泥隧道巖土（以下簡(jiǎn)稱為相似材料土）。

為了研究滑坡區(qū)域在不同水壓下土壤水分的入滲規(guī)律，設(shè)置不同的水頭高度（100、100、200 mm），將鐵沖流域表層土清理干凈之后將模具打進(jìn)土壤中，制作長(zhǎng)度90 mm、直徑100 mm的表層土壤試件，如圖1所示。

圖1 用模具制作的表層土壤試件

在實(shí)驗(yàn)中要驗(yàn)證不同水泥含量對(duì)水泥土滲透率的影響，因此要制作2組水泥土試件，水泥含量分別為10%、20%。水泥土試件的制作方法：先將紅黏土樣品風(fēng)干，然后碾磨壓碎，經(jīng)過2 mm孔徑的篩子，放入（110±2）℃的恒溫箱烘干；再以不同的比例與水泥混合均勻，按比例加入適量的水，攪拌均勻；最后用液壓千斤頂靜壓方法在模具中壓實(shí)水泥土。2組水泥土試件的參數(shù)如表1所示。

表1 2組水泥土試件材料的質(zhì)量占比 %

與水泥土相比，突水突泥的深層巖土含沙量較高，因此在相似材料土試件的配制中，選取山沙為骨料，C32.5礦渣硅酸鹽水泥為膠結(jié)材料，水為輔助材料，黏土為調(diào)節(jié)劑。相似材料土試件的制作方法：用2 mm孔徑的篩子除去山沙中的粗顆粒；用0.5 mm孔徑的篩子過濾黏土；黏土的粒徑比山沙小且具有黏性，在配比時(shí)適當(dāng)增加黏土含量。模型實(shí)驗(yàn)斷層填充材料由黏土、山沙、碎石子按一定比例配制而成，黏土粒徑小于1 mm，碎石子粒徑范圍為4.75～9.50 mm。不同相似材料土試件中各材料的質(zhì)量配比見表2。

表2 不同相似材料土試件中各材料的質(zhì)量配比

1.2 滲透系數(shù)實(shí)驗(yàn)平臺(tái)搭建以及滲透系數(shù)實(shí)驗(yàn)設(shè)計(jì)

實(shí)驗(yàn)平臺(tái)搭建的主要內(nèi)容是實(shí)驗(yàn)設(shè)備的安裝，實(shí)驗(yàn)設(shè)備主要由2個(gè)部分組成，分別為試件安裝的測(cè)試裝置、數(shù)據(jù)傳輸?shù)膫鞲袃x器，如圖2所示。

圖2 滲透系數(shù)測(cè)試裝置以及各組成部件

首先測(cè)試不同水頭（負(fù)壓）下滲透系數(shù)的動(dòng)態(tài)變化情況，選用原狀土試件進(jìn)行實(shí)驗(yàn)，將模擬降雨噴頭的水頭分別固定于離原狀土試件100、150、200 mm的高度處進(jìn)行3組實(shí)驗(yàn)，記錄其動(dòng)態(tài)滲透系數(shù)的變化。

1.3 滲流規(guī)律的數(shù)學(xué)模型

目前模擬土壤水分入滲過程的模型有很多，其中Richard方程能夠反映時(shí)間與含水率、導(dǎo)水率的關(guān)系［9］，如下式所示：

式（1）中：θ為土壤體積含水率（cm3/cm3）；ψm為非飽和土壤總水勢(shì)；D（θ）為非飽和土壤擴(kuò)散率（cm3/cm3）；K（θ）為非飽和土壤導(dǎo)水率（cm3/min）；t為時(shí)間。

含水量是整個(gè)測(cè)試過程中最為主要的測(cè)試變量，所以該數(shù)學(xué)模型能夠反映整個(gè)滲流過程的滲流規(guī)律。

1.4 設(shè)備測(cè)試機(jī)理

滲透系數(shù)測(cè)定及堵塞裝置的工作原理：在實(shí)驗(yàn)裝置上安裝了2個(gè)電子水壓力傳感器和1個(gè)超聲波流速器，所有傳感器通過模數(shù)轉(zhuǎn)換器連接電腦。通過水壓力傳感器可測(cè)得試件上、下表面的水頭損失；同時(shí)通過流速傳感器可測(cè)出水管內(nèi)水的流速。

在滲透實(shí)驗(yàn)過程中應(yīng)注意避免土壤試件側(cè)壁的滲漏，必須在線連續(xù)記錄滲透系數(shù)的變化過程，模擬土壤在一段時(shí)間內(nèi)的滲流過程并且記錄動(dòng)態(tài)滲透率的變化［10］。滲透系數(shù)的表達(dá)式如下：

式（2）中：K為滲透系數(shù)（mm/s）；a為實(shí)驗(yàn)管截面積（mm2）；L為土壤試件的高度（mm）；A為試樣的橫截面積（mm2）；h1和h2為水面高度（m）；t為滲流時(shí)間（s）。

2 實(shí)驗(yàn)結(jié)果與分析

2.1 實(shí)驗(yàn)過程數(shù)據(jù)的記錄

為了測(cè)試非飽和滲透系數(shù)的變化，在降雨模擬開始之前開啟滲流測(cè)試裝置，然后開啟淋雨泵噴水來(lái)做一定水頭高度的模擬降雨，電腦會(huì)自動(dòng)記錄每秒的滲流數(shù)據(jù)，實(shí)驗(yàn)總時(shí)長(zhǎng)為100 s。在實(shí)驗(yàn)完成之后將數(shù)據(jù)轉(zhuǎn)換成為相應(yīng)的滲透系數(shù)，得到如圖3所示的曲線。

圖3 原狀土、水泥土、相似材料土非飽和滲透系數(shù)的變化

2.2 實(shí)驗(yàn)數(shù)據(jù)分析

將轉(zhuǎn)換之后的數(shù)據(jù)通過數(shù)據(jù)處理軟件制作土壤滲透系數(shù)隨時(shí)間變化的曲線（圖4）。從圖4可以看出：原狀土的滲透系數(shù)最大，20%水泥含量的水泥土次之，而相似材料土的最?。幌嗨撇牧贤恋臐B透系數(shù)先快速降低，然后保持在一個(gè)相對(duì)平穩(wěn)的狀態(tài)；原狀土和水泥土的滲透系數(shù)的下降趨勢(shì)基本相同，且長(zhǎng)時(shí)間保持一個(gè)穩(wěn)定的下降速度。這是因?yàn)樵瓲钔梁退嗤林性瓲钔恋暮慷急容^高，而相似材料土的原狀土含量較低。說明不同原狀土含量的土體的滲透系數(shù)有很大的不同，因此在對(duì)有裂縫的土體進(jìn)行防滲處理時(shí)，可以用水泥和河沙一起加固土體來(lái)降低滲流。

圖4 3種材料的滲透系數(shù)隨時(shí)間的變化

3 CFD-DEM土體滲流模型的應(yīng)用分析

目前關(guān)于滲流方面的數(shù)值模擬主要是CFD的數(shù)值分析，而多孔材料的數(shù)值分析主要集中在CT掃描之后的三維建模，但是CT的成本高，且CT掃描對(duì)試件的尺寸有嚴(yán)格要求，很難滿足工程的需要，因此對(duì)多孔材料滲流數(shù)值模擬的研究多采用CFD-DEM耦合模型，常用的軟件類型主要為Fluent-EDEM軟件耦合以及開源耦合；由于工程上的顆粒數(shù)量多，本次實(shí)驗(yàn)采用穩(wěn)定性更高的Fluent-EDEM耦合進(jìn)行滲流模擬研究。

3.1 CFD-DEM土體滲流耦合理論

借助Fluent通過CFD-DEM耦合，同時(shí)基于離散元軟件EDEM的API功能進(jìn)行二次開發(fā)，來(lái)實(shí)現(xiàn)滲流過程中顆粒之間的接觸力與雨水滲流的變化過程，突出滲流過程與失穩(wěn)過程的變化情況。該求解原理是基于未解析的CFD-DEM模型，由于顆粒遠(yuǎn)小于流體網(wǎng)格尺寸的土體顆粒，同時(shí)顆粒的數(shù)量大，導(dǎo)致解析的網(wǎng)格以及DEM的計(jì)算量大，因而無(wú)法完全解析CFD-DEM，因此使用未解析的CFD-DEM，具體的求解原理如下列公式所示：

上式中：σ表示流場(chǎng)的應(yīng)力張量；是顆粒的外法向矢量；是流體作用于顆粒的牽引矢量；粒與壁面接觸產(chǎn)生的力。

式（3）和式（4）分別為不可壓縮流體的運(yùn)動(dòng)方程和連續(xù)方程，即Navier-Stokes方程。此類方程適用于整個(gè)領(lǐng)域。Dirichlet邊界條件（5）和初始條件（6）完成了流動(dòng)方程組。式（9）描述了拉格朗日粒子的運(yùn)動(dòng)。式（7）和式（8）負(fù)責(zé)流體和固相之間的耦合，其中式（7）用于在流體速度場(chǎng)上傳遞物體的速度，式（8）稱為界面條件，描述了流體和固體之間的應(yīng)力，它可以轉(zhuǎn)化為1個(gè)阻力項(xiàng)。

對(duì)充滿土體顆粒的區(qū)域進(jìn)行滲流的流場(chǎng)分析時(shí)，可以采用解析CFD-DEM的耦合方法來(lái)實(shí)現(xiàn)。由于解析CFD-DEM對(duì)計(jì)算資源的消耗較大，因此可以選取局部土體顆粒區(qū)域進(jìn)行滲流的流場(chǎng)研究。

3.2 不同非飽和滲流對(duì)比

通過局部非飽和滲流數(shù)值模擬過程（圖5）來(lái)研究滲流過程中滲流的變化情況，并采用解析CFD-DEM數(shù)值模擬分析，詳細(xì)分析局部小區(qū)域滲流流場(chǎng)的變化，以降低計(jì)算量。土體顆粒模型如圖6～圖9所示。

圖5 整體土體模型和截?cái)嗝孢x取

圖6 在非飽和的初始階段氣液兩相以及壓力分布

圖9 液面上壓力和速度的變化

圖7 滲流過程中不同時(shí)間的相分布云圖

圖8 在土體滲流穩(wěn)定階段壓力分布和速度矢量

隨著滲流時(shí)間的延長(zhǎng)，液面會(huì)發(fā)生一系列的變化，而且內(nèi)部的壓力也會(huì)發(fā)生一系列變化。保持上面的工況，設(shè)置同樣的數(shù)值模型來(lái)計(jì)算。相似材料土體中顆粒的內(nèi)部孔隙被細(xì)小的沙粒和水泥填充，從而形成了水泥土，其滲流特征如圖10～圖11所示。

圖10 顆?？紫兑约肮羌苣Ｐ?/p>

從圖11中可以看到在非飽和流動(dòng)過程中氣液兩相界面在滲流過程中的變化，以及液面附近的流動(dòng)狀況，其壓力速度隨時(shí)間的變化如圖12所示。

圖11 滲流過程中不同時(shí)間的液面以及速度矢量

圖12 液面壓力和速度的變化

隨著滲流時(shí)間的變化，在相似材料初始液面附近壓力和速度的變化與土體顆粒完全相同，只是在數(shù)值上有所不同。為了解不同材料的孔隙內(nèi)部流動(dòng)規(guī)律，監(jiān)測(cè)不同材料流動(dòng)的同一位置，其流速變化如圖13所示。

圖13 土顆粒以及由土顆粒和水泥沙組成的相似材料在不同時(shí)刻的滲流速度

對(duì)比上述2種材料在不同時(shí)間的滲流速度，發(fā)現(xiàn)這2種材料孔隙內(nèi)部的滲流速度隨入滲距離的變化趨勢(shì)基本一致，只是在數(shù)值上略有不同。另外，對(duì)于相似材料來(lái)說，在土體顆粒模型的孔隙中填充很少的顆粒就能夠?qū)⒘鲌?chǎng)改變，使平均滲流速度下降50%以上，這對(duì)于改造材料形成新的材料有重要意義。該數(shù)值模擬結(jié)果與實(shí)驗(yàn)結(jié)果完全吻合，即在土體顆粒之間增加更小的水泥顆粒會(huì)導(dǎo)致孔隙內(nèi)流速減小。

4 小結(jié)與討論

本實(shí)驗(yàn)結(jié)果表明：在不同材料滲流的過程中滲透系數(shù)呈動(dòng)態(tài)下降的趨勢(shì)，原狀土含量越高，滲透系數(shù)的下降趨勢(shì)越明顯；相似材料土的原狀土含量低，其滲透系數(shù)的后期下降趨勢(shì)不明顯；在土體滲流過程中因液相填充在土體的孔隙之中而導(dǎo)致滲透系數(shù)的變化；不同材料初始液面附近的壓力和速度隨入滲深度增加的變化趨勢(shì)一致；在不同時(shí)間點(diǎn)相同材料孔隙內(nèi)部的滲流速度隨入滲距離的變化趨勢(shì)也基本一致；相似材料在土體顆粒模型的孔隙中填充很少的顆粒就能夠?qū)⒘鲌?chǎng)改變，使平均滲流速度下降50%以上，因此向土體顆粒摻入細(xì)小的石灰或者沙石可用于防滲工程。

本實(shí)驗(yàn)并沒有考慮顆粒的遇水膨脹以及顆粒之間從充滿空氣到水分滲入的流動(dòng)狀態(tài)，今后擬繼續(xù)開發(fā)相關(guān)的數(shù)值分析求解器和模型來(lái)解決這些問題，以便能更加深入地了解滲流過程的變化情況。