夏 霄，甘鵬路

（中國(guó)電建集團(tuán)華東勘察設(shè)計(jì)研究院有限公司，浙江 杭州 311122）

1 概述

管涌作為一種土體內(nèi)部侵蝕現(xiàn)象是基坑工程中的常見災(zāi)害，它引起的基坑擋墻倒塌事故時(shí)有發(fā)生[1]。例如，2015 年杭州某濱江基坑開挖中發(fā)生了管涌誘發(fā)的擋墻倒塌事故，承受較高水力梯度的砂質(zhì)土體從擋墻裂縫處涌出，造成坑外土體大幅沉降和擋墻的嚴(yán)重變形[2]。類似基坑擋墻薄弱處（如樁或地連墻連接處）破裂引發(fā)的管涌災(zāi)害屢有發(fā)生[3-4]。

正是由于管涌細(xì)觀機(jī)制的復(fù)雜性，相應(yīng)的數(shù)值分析研究也逐步展開。一類是利用有限元將土體視為連續(xù)性介質(zhì)并與滲透水流進(jìn)行耦合的唯象分析[5-7]，它很難捕捉管涌中土顆粒流失，重新排列、沉積等重要特征。因此，另一類基于流體力學(xué)（CFD）和散體介質(zhì)（DEM）的多相耦合仿真分析方法（CFD-DEM）在管涌研究中越發(fā)常見[8]。這種方法將土體視為不同相的組合，如土骨架以及土骨架孔隙中的水流與細(xì)顆粒，它可考慮內(nèi)部侵蝕或管涌過(guò)程中細(xì)顆粒從土骨架中的流失。CFD-DEM耦合仿真分析方法可捕捉基坑管涌過(guò)程中土體的主要力學(xué)特征，也可用于解決實(shí)際基坑管涌問(wèn)題[9]。

本文利用CFD-DEM 耦合仿真分析方法研究細(xì)粒含量（Fc）與水力梯度（i）對(duì)間斷級(jí)配土基坑中管涌災(zāi)害的影響。首先介紹了CFD-DEM 流固耦合仿真模擬方法，然后通過(guò)建立基坑管涌模型，考慮不同的Fc和i進(jìn)行一系列參數(shù)分析，得到不同工況下基坑土體顆粒速度場(chǎng)和流場(chǎng)特征，驗(yàn)證了該方法在基坑管涌模擬中的優(yōu)越性。

2 CFD-DEM流固耦合仿真分析方法

CFD-DEM 方法由三個(gè)部分構(gòu)成，分別是離散元法（DEM），計(jì)算流體力學(xué)（CFD）及CFD-DEM耦合算法。DEM主要用于對(duì)土顆粒的仿真，CFD則是仿真土中的滲透水流。CFD-DEM耦合算法主要計(jì)算流體—顆粒間相互作用力，以及實(shí)現(xiàn)CFD 與DEM 間信息交互。開源 DEM 軟件 LIGGGHTS 和 CFD 軟件 Open?FOAM 是本研究中具體使用的多相耦合仿真分析工具。下面簡(jiǎn)要介紹這三個(gè)仿真分析方法的基本原理和控制方程。

2.1 DEM方法控制方程

DEM 中顆粒運(yùn)動(dòng)遵循牛頓第二定律,顆粒i的平動(dòng)和轉(zhuǎn)動(dòng)由下式計(jì)算：

式中：mi、It——顆粒i的質(zhì)量和轉(zhuǎn)動(dòng)慣量；

Ui、ωi——顆粒i的平動(dòng)速度與轉(zhuǎn)動(dòng)角速度；

本研究中顆粒間接觸模型采用Hertz 非線性接觸定律[10]。對(duì)于切向上的滑動(dòng)阻力，利用Mohr-Coulomb準(zhǔn)則計(jì)算。

2.2 CFD方法控制方程

CFD 域中，流體運(yùn)動(dòng)主要由連續(xù)性方程和動(dòng)量方程控制，如下所示：

式中：Uf、p——流體網(wǎng)格的平均速度和壓強(qiáng)；

εf——流體網(wǎng)格中顆粒體的孔隙率；

fp——單位體積顆粒與流體間平均相互作用力，即

ΔV——流體網(wǎng)絡(luò)體積；

k——顆粒數(shù)；

ρf、μf——流體密度和粘度；

g——重力加速度，本文中為50倍的常重力速度。

2.3 流體—顆粒相互作用力計(jì)算

在管涌問(wèn)題中，流體—顆粒相互作用力Ff主要包括拖曳力、粘滯力和動(dòng)壓差力。本研究中也主要考慮上述三種流固相互作用力，忽略其它相互作用力對(duì)管涌特性的影響。Ff的計(jì)算如公式（3）所示：

流體對(duì)土顆粒產(chǎn)生的拖曳力采用Buijtenen 等（2011）[11]方法計(jì)算，如公式（4）所示，它適用于顆粒體密度較大的情況，在流體雷諾數(shù)變化較大時(shí)也成立：

式中：dp、Up——顆粒的直徑與速度；

Vp——單個(gè)顆粒的體積；

εp——流體網(wǎng)格中顆粒所占體積百分比，εp+εf=1；

Rep——單個(gè)顆粒的雷諾數(shù)。

作用于顆粒上的動(dòng)壓差力和粘滯力計(jì)算見公式（5）和（6）[12]：

式中：τ——流體剪切力。

3 基坑管涌仿真分析工況與參數(shù)設(shè)置

3.1 仿真分析工況

本文主要通過(guò)不同的土體細(xì)粒含量Fc與水力梯度i工況來(lái)研究其對(duì)基坑管涌的影響。一般認(rèn)為Fc直接決定土體組構(gòu)，采用Fc=15%和35%的試樣研究?jī)煞N不同的典型土體組構(gòu)對(duì)基坑管涌的影響。

為研究基坑管涌啟動(dòng)前后土體力學(xué)性質(zhì)演變，在仿真過(guò)程中，水力梯度由i=0.8 增至1.8（i= Δp/ρgL，Δp是基坑內(nèi)處水壓差，L是模型中流體最短流動(dòng)路徑，本研究中是0.5m，g是模型中采用的重力加速度值）。當(dāng)i=0.8時(shí)，管涌并未啟動(dòng)，而i增至1.8時(shí)，管涌開始發(fā)生并持續(xù)至仿真結(jié)束。表1總結(jié)了本研究考慮的四個(gè)管涌仿真工況。圖1是本研究中采用的土體的級(jí)配曲線，粗粒組和細(xì)粒組均為單一直徑顆粒，分別為2mm和0.5mm?？紤]到CFD-DEM 耦合仿真計(jì)算效率，每一種工況的仿真時(shí)長(zhǎng)均為1.0s。結(jié)果表明，在此較短的仿真時(shí)長(zhǎng)內(nèi)，基坑變形也可達(dá)到穩(wěn)定狀態(tài)。

圖1 試樣級(jí)配曲線

表1 仿真分析工況

3.2 模型尺寸與參數(shù)設(shè)置

圖2 是基坑DEM 模型的尺寸示意圖。由于本研究主要針對(duì)管涌過(guò)程中顆粒與流體力學(xué)特性演化，因此擋墻假定為鋼性體，在仿真過(guò)程中不發(fā)生變形。為研究擋墻裂縫引起的基坑管涌災(zāi)害，在擋墻10cm高度處開有2.5cm 高度的孔洞。擋墻右側(cè)和左側(cè)分別是基坑的開挖側(cè)與支護(hù)側(cè)，高度分別為50cm 和25cm，長(zhǎng)度分別均為25cm和50cm。整個(gè)模型厚度為5cm，將實(shí)際的基坑管涌簡(jiǎn)化成平面應(yīng)變問(wèn)題。結(jié)果表明這一尺寸設(shè)置可以減小邊界效應(yīng)的影響。仿真模型的CFD域與DEM域范圍大小相同，保證所有顆粒均受滲透水流的作用。表2列出了模型中顆粒與流體的力學(xué)性質(zhì)參數(shù)。

圖2 基坑管涌仿真DEM模型

表2 CFD與DEM參數(shù)匯總

3.3 邊界條件設(shè)置

為模擬實(shí)際中基坑擋墻兩側(cè)水壓差，在CFD 域的開挖側(cè)頂面施加一定水壓，使得基坑內(nèi)外水力梯度在整個(gè)仿真過(guò)程中均維持在0.8或1.8。CFD域其余面均為不透水面。對(duì)于顆粒與DEM邊界的接觸，假定它們之間的摩擦系數(shù)等于0，法向接觸剛度是顆粒彈性模量的10倍。

3.4 仿真分析步驟

首先在基坑區(qū)域上方生成一定數(shù)量的顆粒并使其在重力作用下自由下落。在填至模型1/5高度時(shí)，將顆粒體豎向縮至一定的相對(duì)密實(shí)度，本研究中將土體相對(duì)密實(shí)度設(shè)定為Re=65%～67%。重復(fù)上述兩個(gè)步驟直至顆粒填滿基坑支護(hù)側(cè)與開挖側(cè)。在DEM模型準(zhǔn)備完成后，立即在CFD域的支護(hù)與開挖側(cè)頂面上施加固定的水壓力差以開始基坑管涌的仿真模擬。所有顆粒及它們間的接觸信息在整個(gè)仿真過(guò)程中每隔0.01s輸出。

4 基坑管涌仿真分析主要結(jié)果

圖3是細(xì)粒含量分別為Fc=15%和35%情況下基坑管涌啟動(dòng)時(shí)刻（t=0.2s，i=1.8）地層變形與顆粒速度場(chǎng)。在管涌啟動(dòng)時(shí)刻，擋墻裂縫處和基坑開挖側(cè)坑底土體首先開始流動(dòng)，具有較大的初始速度，而支護(hù)側(cè)土體相對(duì)流動(dòng)速度較小。圖4 是在Fc=15%和35%時(shí)，基坑管涌仿真結(jié)束時(shí)刻地層變形與顆粒速度場(chǎng)。此時(shí)，支護(hù)側(cè)土體出現(xiàn)明顯沉降，開挖側(cè)坑底土體明顯隆起。在擋墻裂縫前方出現(xiàn)一條管涌通道，支護(hù)側(cè)土顆粒經(jīng)此流入開挖側(cè)并逐漸趨于穩(wěn)定。當(dāng)土體中含較多細(xì)粒時(shí)，支護(hù)側(cè)土體沉降較大。這主要由于裂縫寬度與顆粒粒徑之比較大時(shí)，顆粒不易堵塞，更易通過(guò)裂縫，因而在相同水力條件下流動(dòng)性更強(qiáng)，導(dǎo)致最終支護(hù)側(cè)顆粒流失量增大。

圖4 最終時(shí)刻(t=1.0s)基坑變形與顆粒速度場(chǎng)

圖5 是管涌過(guò)程基坑流場(chǎng)的CFD-DEM 仿真結(jié)果?；又ёo(hù)側(cè)水體在水頭差作用下通過(guò)擋墻裂縫流入開挖側(cè)，裂縫周圍水流速度達(dá)到最大。

圖5 基坑管涌過(guò)程流場(chǎng)仿真結(jié)果

此外，通過(guò)比較不同工況下支護(hù)側(cè)土體流失量，可以發(fā)現(xiàn)，對(duì)于土體細(xì)粒含量較低（Fc=15%）情況，土體流失量較多。如前分析，這主要由于裂縫寬度與顆粒粒徑之比較大時(shí)，顆粒在裂縫處不易堵塞，更易流動(dòng)導(dǎo)致的。

進(jìn)一步對(duì)比Fc=15%和35%兩種工況下的支護(hù)側(cè)地表沉降演化，可以發(fā)現(xiàn)，地表沉降都在靠近擋墻處最大，且在管涌啟動(dòng)時(shí)刻（t=0.2s）不明顯，而在t=0.4～0.8s 時(shí)刻內(nèi)發(fā)展較快，隨后沉降增量趨于減小。對(duì)于Fc=15%的情況，由于其支護(hù)側(cè)土體流失量較小，地表沉降范圍主要局限在0

然而，開挖側(cè)坑底土體隆起演化結(jié)果與支護(hù)側(cè)地表沉降演化不同。在管涌啟動(dòng)時(shí)刻（t=0.2s），坑底隆起即達(dá)到最大值，隨后下降并由于支護(hù)側(cè)土體流入而最終略有上升。這一現(xiàn)象與圖3中開挖側(cè)土體在管涌啟動(dòng)時(shí)刻即具有較大速度相吻合。這說(shuō)明坑底隆起主要是管涌流對(duì)坑底土顆粒向上的拖曳力，使得土體懸浮于管涌流中引起。

5 結(jié)論

本文介紹了一種可反映流體與砂土顆粒間相互作用的CFD-DEM 方法，并利用該法對(duì)深基坑開挖過(guò)程擋墻裂縫導(dǎo)致的管涌問(wèn)題進(jìn)行了仿真分析。通過(guò)砂土基坑管涌的宏細(xì)觀仿真結(jié)果，可以發(fā)現(xiàn)CFD-DEM 方法可以較好重現(xiàn)基坑管涌的主要特征，如坑底隆起、地表沉降和管涌通道中流體流動(dòng)特征。管涌過(guò)程中地表沉降逐步增大，而坑底隆起在管涌啟動(dòng)時(shí)即達(dá)到最大值。管涌災(zāi)害主要由滲流力對(duì)基坑開挖側(cè)土體力學(xué)性質(zhì)的削弱引起，在實(shí)際工程中為了防止管涌災(zāi)害的發(fā)生，除提高擋墻或止水帷幕施工質(zhì)量外，仍需對(duì)開挖側(cè)土體實(shí)施良好的加固措施。