基于SPH-DEM的泥石流沖擊剛性防護(hù)結(jié)構(gòu)的動(dòng)力過程研究

趙德博,郝夢(mèng)潔,熊 昊

(1. 深圳大學(xué) 土木與交通工程學(xué)院,廣東 深圳 518060; 2. 深圳大學(xué) 濱海城市韌性基礎(chǔ)設(shè)施教育部重點(diǎn)實(shí)驗(yàn)室,廣東 深圳 518060)

0 引言

泥石流作為一種山區(qū)地帶多發(fā)的地質(zhì)災(zāi)害,包含泥砂、塊石的固液兩相流體,呈黏性層流或稀性紊流狀態(tài)[1]。由于其往往突然暴發(fā)、來勢(shì)兇猛、破壞力強(qiáng),泥石流對(duì)人類的生命財(cái)產(chǎn)安全構(gòu)成巨大威脅。我國(guó)泥石流幾乎遍布每個(gè)省、市和自治區(qū),據(jù)自然資源部統(tǒng)計(jì),僅2022年上半年,特別是5月底以來,南方多省遭遇多輪強(qiáng)降雨,地質(zhì)災(zāi)害多發(fā)群發(fā),特別是浙江、福建、江西、湖南、廣東、廣西6省(區(qū))發(fā)生4000多處崩塌滑坡泥石流,已超去年全年的數(shù)量[2]。2019年6月27日20時(shí)30分,四川省金川縣曾達(dá)鄉(xiāng)曾達(dá)溝(流域面積125.53 km2)暴發(fā)特大型泥石流災(zāi)害,造成房屋、基礎(chǔ)設(shè)施、農(nóng)業(yè)、公益設(shè)施損失約1.51億元[3]。泥石流災(zāi)害后果極為嚴(yán)重,因此對(duì)于泥石流易發(fā)區(qū)域迫切需要采取治理措施。沿著泥石流路徑建造防護(hù)結(jié)構(gòu)是目前泥石流治理采用的主要工程措施。防護(hù)結(jié)構(gòu)主要包括剛性[4]和柔性[5-6]2種,它的設(shè)計(jì)可以改變泥石流的運(yùn)動(dòng)狀態(tài),減弱甚至避免泥石流對(duì)下游造成破壞。然而目前在工程實(shí)踐中,主要還是采用經(jīng)驗(yàn)方法對(duì)防護(hù)結(jié)構(gòu)進(jìn)行設(shè)計(jì)與施工。究其原因,主要是缺乏對(duì)泥石流運(yùn)動(dòng)特點(diǎn)以及泥石流與防護(hù)結(jié)構(gòu)相互作用機(jī)理的研究。

許多研究人員都通過現(xiàn)場(chǎng)實(shí)測(cè)、室內(nèi)水槽試驗(yàn)和數(shù)值模擬[7-9]等手段對(duì)泥石流對(duì)防護(hù)結(jié)構(gòu)的沖擊開展過探究。劉道川等[10]通過小型水槽試驗(yàn)開展了系列泥石流的沖擊試驗(yàn)研究了兩相泥石流沖擊力的時(shí)空分布特征。而NG等[11]將泥石流簡(jiǎn)化為均一的流體,借助水槽試驗(yàn)研究了玻璃球和砂子分別作用于剛性防護(hù)結(jié)構(gòu)的沖擊特性和機(jī)理,揭示了玻璃球作用下的爬升機(jī)制和砂顆粒作用下的堆積機(jī)制。此外,其他學(xué)者利用水槽試驗(yàn)裝置研究了不同組成成分的泥石流與防護(hù)結(jié)構(gòu)的作用機(jī)制,發(fā)現(xiàn)泥石流的堆積形態(tài)、爬升高度和泥石流的組成密切相關(guān)[12-13]。相比于原型場(chǎng)試驗(yàn),小尺度水槽試驗(yàn)存在明顯的尺寸效應(yīng),無法完全還原泥石流的動(dòng)力學(xué)特征。但大尺度的原型試驗(yàn)成本較高,工作量大,不確定因素多,不適合泥石流對(duì)防護(hù)結(jié)構(gòu)沖擊機(jī)理的探究。因此,近年來數(shù)值模擬在巖土工程方面[13-15],尤其是泥石流運(yùn)動(dòng)機(jī)理方面得到廣泛的應(yīng)用[16]。甘建軍等[8]采用RAMMAS3D數(shù)值軟件準(zhǔn)確模擬了泥石流的運(yùn)動(dòng)變形過程并確定了其影響區(qū)域。HE等[17]利用光滑粒子流體動(dòng)力學(xué)(SPH)方法模擬了干顆粒流和剛性防護(hù)結(jié)構(gòu)的作用,而CUI等[18]利用離散元法研究了不同粒徑的顆粒流撞擊防護(hù)系統(tǒng)的過程,并與試驗(yàn)結(jié)果進(jìn)行比對(duì),取得了較好的效果。CALVETTI等[19]則從微觀角度解釋了干顆粒流對(duì)防護(hù)結(jié)構(gòu)的沖擊。但以上數(shù)值模擬方法都傾向于將泥石流模擬為單一相的流體,這種簡(jiǎn)化限制了對(duì)泥石流這種多相流體的深入理解。

因此,本研究采用SPH-DEM耦合方法開展水槽試驗(yàn)的數(shù)值模擬,其中SPH模擬泥石流中的黏性流體,DEM模擬泥石流中的石塊,探究含石塊的泥石流的運(yùn)動(dòng)過程及其對(duì)剛性防護(hù)結(jié)構(gòu)的沖擊效應(yīng),并研究了不同水槽傾角對(duì)泥石流的運(yùn)動(dòng)形態(tài)以及對(duì)剛性防護(hù)結(jié)構(gòu)的沖擊力的影響,為剛性防護(hù)結(jié)構(gòu)強(qiáng)度的工程設(shè)計(jì)標(biāo)準(zhǔn)提供參考。

1 SPH-DEM數(shù)值模擬方法

1.1 SPH數(shù)值方法

光滑粒子流體動(dòng)力學(xué)(smooth particle hydrodynamics,SPH)是一種純拉格朗日無網(wǎng)格粒子方法,是當(dāng)前數(shù)值計(jì)算的熱點(diǎn)話題之一。該方法基于連續(xù)介質(zhì)假設(shè)以一系列任意分布的粒子來模擬流體運(yùn)動(dòng)規(guī)律,其中每一個(gè)粒子都有獨(dú)立的影響域和插值域,計(jì)算過程中粒子的速度、壓強(qiáng)和位置等性質(zhì)以及它們的梯度分布都是根據(jù)插值域內(nèi)粒子構(gòu)造的核近似函數(shù)插值得到的。

在SPH方法中,位于r處的流體粒子特性的函數(shù)可由以下的積分表達(dá)式定義:

(1)

式中:h為光滑長(zhǎng)度,控制積分域的大小;Ω為支持域;r′為支持域內(nèi)其他粒子的位置;W為核函數(shù),它不僅決定了SPH近似的內(nèi)值精度,而且與數(shù)值穩(wěn)定性有關(guān)。本文采用Wendland核函數(shù)[20]進(jìn)行平滑,其表達(dá)式為:

(2)

式中:αD在二維問題中等于7/(4πh2),在三維問題中等于21/(16πh3);q=|r-r′|/h。

1.2 DEM數(shù)值方法

離散單元法(discrete element method,DEM)最早是由CUNDALL[21]于1971年提出的,已被廣泛應(yīng)用于顆粒材料[22]、巖土工程[23]等領(lǐng)域。DEM的基本思想是利用一系列離散的粒子來模擬固體介質(zhì),其中每個(gè)粒子作為一個(gè)單元,具有固定的大小和形狀,單元之間可以在小范圍內(nèi)相互擠壓重疊從而產(chǎn)生相互作用力。當(dāng)單元i與單元j接觸時(shí),根據(jù)牛頓第二定律,單元i的運(yùn)動(dòng)方程可以寫為:

(3)

(4)

式中:mi為流體粒子i的質(zhì)量;u為單元i的運(yùn)動(dòng)速度;g為重力加速度;Fn,ij和Ft,ij分別是單元j對(duì)單元i的法向和切向作用力;Ii和ωi分別是單元i的轉(zhuǎn)動(dòng)慣量和角速度;而Tt,ij為Fn,ij和Ft,ij對(duì)單元i的力矩。

為了簡(jiǎn)化接觸模型,單元均用球形粒子來表示,那么對(duì)于發(fā)生碰撞的球形粒子,粒子間的法向接觸力如式(5)所示:

Fn=knδn+cnvn

(5)

式中:kn為法向剛度系數(shù);δn為法向疊合量;cn為法向黏性阻尼系數(shù);vn為法向相對(duì)速度。

粒子間的切向力大小與相接觸的粒子間的相對(duì)速度有關(guān),滿足庫(kù)侖定律,如式(6)所示:

Ft=min(ktvt,μFn)

(6)

式中:kt為切向剛度系數(shù);vt為切向相對(duì)速度。

1.3 SPH-DEM耦合方法

SPH和DEM的耦合計(jì)算實(shí)際上是將流體作用于固體顆粒表面上的作用力進(jìn)行積分得到流體對(duì)固體的作用力,而固體顆粒表面單元內(nèi)插到與其相鄰的流體粒子的SPH方程中,從而可以計(jì)算固體對(duì)流體的作用力。

在SPH-DEM耦合方法中,流體粒子的動(dòng)量方程式為:

(7)

固體DEM顆粒所受的作用力主要包括固體顆粒間的接觸力、潤(rùn)滑力、流體施加的力(浮力和阻力)和重力等,其控制方程式為:

(8)

(9)

2 SPH-DEM方法驗(yàn)證

為了驗(yàn)證本文采用的SPH-DEM耦合方法在模擬流體與固體之間復(fù)雜相互作用的可靠性,模擬了單個(gè)非旋轉(zhuǎn)球體的入水測(cè)試,并結(jié)合ARISTOFF等[24]開展的室內(nèi)試驗(yàn)和XU等[25]進(jìn)行的模擬校驗(yàn)結(jié)果。

試驗(yàn)?zāi)Ｐ偷难b置如圖1所示,水箱的長(zhǎng)為0.20 m,寬為0.20 m,高為0.14 m,水深0.11 m。直徑為25.4 mm的球體以2.17 m/s的初始速度與水面相切釋放。數(shù)值模擬中使用的球體和水的參數(shù)如表1所示。模擬時(shí)考慮重力,重力加速度為9.81 m/s2。SPH流體粒子初始間距dp為球體半徑的1/10,采用Verlet時(shí)間積分算法,模擬1 s內(nèi)小球進(jìn)入水中的運(yùn)動(dòng)狀態(tài)。

圖1 球體入水的數(shù)值模型Fig. 1 Numerical model for water-entry of a sphere表1 球體入水模擬的關(guān)鍵參數(shù)Table 1 Parameters used in the simulation for water-entry of a sphere材料參數(shù)取值球體密度/(kg/m3)860楊氏模量/(MPa)100泊松比0.3水密度/(kg/m3)1000黏度/(×10-6 Pa?s)1

不同時(shí)刻小球進(jìn)入水中的運(yùn)動(dòng)過程如圖2所示,與ARISTOFF等[24]的試驗(yàn)在下降趨勢(shì)、水面形狀等方面一致。球體與水面接觸并開始向下運(yùn)動(dòng),沖擊并擠壓水,致使球體周圍的水流向側(cè)向運(yùn)動(dòng),并在球體上方形成一個(gè)空腔。隨著小球繼續(xù)向下運(yùn)動(dòng),空腔塌陷,形成向上的水射流。而球體由于密度較小,速度逐漸減小為0,隨后在浮力的作用下被抬升。球體的整個(gè)運(yùn)動(dòng)軌跡都是沿一條垂直直線,沒有任何的旋轉(zhuǎn)。這是因?yàn)榍蝮w在耦合過程中只受來自垂直方向的力,這與理論情況也是相符的。

圖2 小球入水過程模擬圖Fig. 2 Simulation snapshots of the sphere entering water

球體入水深度隨時(shí)間的變化情況如圖3所示,并與試驗(yàn)數(shù)據(jù)、理論數(shù)據(jù)和XU等[25]的模擬結(jié)果做了對(duì)比,發(fā)現(xiàn)在0.07 s之前,隨著時(shí)間的不斷增大,球體在水中的深度逐漸變大,但深度的增大率,即球體的速度卻逐漸變小。模擬效果與理論結(jié)果較為接近,且由于所采用的粒子間距大于XU等[25]的模擬但結(jié)果相差不大,模擬效率大大提升。由此可見,本文采用的SPH-DEM耦合方法具有可行性,能夠較好地模擬流體與固體之間的相互作用。

圖3 球體入水的深度隨時(shí)間的演變Fig. 3 Evolution of the depth of sphere into the water with time

考慮到多數(shù)泥石流中的流體均為黏性流體,本文進(jìn)一步研究了球體進(jìn)入黏性流體的運(yùn)動(dòng),并提取其深度,如圖4所示?？梢钥吹皆陴ば粤黧w(黏度=1)中,球體的深度隨時(shí)間變化得慢,速度減小得快,在0.8 s左右,速度基本為0。因此,耦合模型模擬黏性流體與固體的相互作用依舊有效且可靠。

3 泥石流沖擊剛性防護(hù)結(jié)構(gòu)

長(zhǎng)期以來,泥石流頻發(fā)區(qū)防護(hù)結(jié)構(gòu)的設(shè)計(jì)一直是人們關(guān)注的問題。目前許多學(xué)者已經(jīng)針對(duì)泥石流與防護(hù)結(jié)構(gòu)的作用過程開展了大量試驗(yàn)研究。其中CHOI等[26]利用小型水槽裝置研究了不同傾角下干性泥石流沖擊下,防護(hù)結(jié)構(gòu)的動(dòng)態(tài)受荷過程,揭示了水槽傾角與沖擊力對(duì)防護(hù)結(jié)構(gòu)設(shè)計(jì)的重要性。本文采用SPH-DEM的耦合方法模擬水槽試驗(yàn)中包含石塊的泥石流在水槽內(nèi)的流動(dòng)以及對(duì)剛性防護(hù)結(jié)構(gòu)的沖擊,在此基礎(chǔ)上改變水槽傾角、石塊的粒徑、泥石流中流體的體積,從而分析了不同條件下防護(hù)結(jié)構(gòu)所受沖擊力以及其上水面高程的演變。

3.1 計(jì)算模型

水槽試驗(yàn)?zāi)M的數(shù)值模型尺寸如圖5所示。水槽的長(zhǎng)、寬和高分別為2.3、0.6、0.7 m,防護(hù)結(jié)構(gòu)距離料箱門1.6 m,料箱的容量為0.168 m3。防護(hù)結(jié)構(gòu)的寬度、高度和厚度分別為600、700、5 mm。水槽底部與水平面的傾角為15°。

圖5 水槽試驗(yàn)?zāi)M的模型示意圖Fig. 5 Schematic diagram of the model for the simulation of the flume experiments表2 水槽試驗(yàn)?zāi)M的關(guān)鍵參數(shù)Table 2 Parameters used in the simulation of the flume experiment材料參數(shù)取值流體密度/(kg/m3)1000黏度/(Pa·s)1玻璃珠密度/(kg/m3)2000楊氏模量/MPa200泊松比0.2動(dòng)態(tài)摩擦系數(shù)0.2恢復(fù)系數(shù)0.6

在初始狀態(tài)時(shí),泥石流位于水槽頂部,由固液兩相組成。其中由SPH模擬泥石流中的黏性流體,流體的體積為0.072 m3,流體的黏度為1 Pa·s,流體粒子之間的間距為0.005 m,流體總粒子數(shù)為573934,而采用直徑為4 cm 的玻璃珠模擬泥石流中的石塊,本模擬中包含27個(gè)玻璃珠,即模擬的泥石流的固體體積分?jǐn)?shù)為1.26%。整個(gè)模擬過程中考慮重力,重力加速度被設(shè)置為9.81 m/s2。模擬需要用到的其它參數(shù)如表2所示。

3.2 泥石流運(yùn)動(dòng)過程分析

泥石流在水槽中運(yùn)動(dòng)過程的模擬圖如圖6所示,SPH粒子的顏色代表流體的流速。在模擬之初,泥石流被生成在水槽頂部的料箱中,料箱門可以有效地阻止泥石流的運(yùn)動(dòng)。當(dāng)料箱中的泥石流達(dá)到穩(wěn)定狀態(tài)(T=1.5 s)后,料箱門向上抽出,模擬泥石流的啟動(dòng)。在重力的作用下,兩相流沿著水槽底部運(yùn)動(dòng)。隨著運(yùn)動(dòng)速度的不斷增大,泥石流中的流體前部逐漸靠近剛性防護(hù)結(jié)構(gòu),當(dāng)其一接觸到防護(hù)結(jié)構(gòu),流體的速度迅速減小,并沿著防護(hù)結(jié)構(gòu)向上爬升(T=2.5 s),在其動(dòng)能減少為0后,流體達(dá)到最大爬升高度,開始向后回滾(T=2.8 s),一部分堆積在防護(hù)結(jié)構(gòu)的底部,形成一個(gè)穩(wěn)定的堆積區(qū),一部分與后續(xù)泥石流相撞,使其速度減少(T=3.1 s),而后續(xù)的泥石流不會(huì)直接撞擊到防護(hù)結(jié)構(gòu),而是直接匯入堆積體,堆積區(qū)的大小逐漸增大,直至達(dá)到穩(wěn)定,而泥石流中的玻璃珠由于與水槽底部的摩擦以及與流體的相互作用的影響,運(yùn)動(dòng)速度要小于流體,故其到達(dá)防護(hù)結(jié)構(gòu)的時(shí)間更晚,在流體已經(jīng)形成堆積體時(shí),玻璃珠會(huì)沿著堆積體向上運(yùn)動(dòng)(T=5.0 s)。

圖6 泥石流的運(yùn)動(dòng)過程Fig. 6 Runout process of the debris flow

為了進(jìn)一步量化分析泥石流對(duì)防護(hù)結(jié)構(gòu)的沖擊效應(yīng),模擬中還提取了泥石流中的黏性流體對(duì)防護(hù)結(jié)構(gòu)的沖擊力隨時(shí)間的演變,如圖8中的綠色曲線所示。在t=2.36 s時(shí),流體前端接觸到剛性防護(hù)結(jié)構(gòu)底部,動(dòng)態(tài)沖擊力由0迅速增大,流體沿著防護(hù)結(jié)構(gòu)向上爬升時(shí),沖擊力增加速度放緩,隨后隨著流體爬升到最大高度后開始回流以及后續(xù)泥石流的沖擊,沖擊力存在波動(dòng)并在2.8 s左右達(dá)到峰值。之后的泥石流由于被堆積在防護(hù)結(jié)構(gòu)后部的流體阻擋,并不能直接接觸到防護(hù)結(jié)構(gòu),故對(duì)防護(hù)結(jié)構(gòu)的沖擊力逐漸減弱并趨于平緩,造成這種現(xiàn)象的原因主要是泥石流體積是一定的,防護(hù)結(jié)構(gòu)后的沖擊力流量不再隨著時(shí)間繼續(xù)增加。而穩(wěn)定前曲線仍存在的微小波動(dòng),是由于防護(hù)結(jié)構(gòu)后堆積體的回流以及后續(xù)流體匯入堆積體對(duì)防護(hù)結(jié)構(gòu)存在擠壓間接造成的。

泥石流與防護(hù)結(jié)構(gòu)的相互作用受多種因素的共同影響,包括泥石流的屬性(如泥石流的體積、固體體積分?jǐn)?shù)、密度和流體黏度等)、防護(hù)結(jié)構(gòu)的位置、坡度以及水槽傾角等,這些都是設(shè)計(jì)防護(hù)結(jié)構(gòu)需要考慮的重要因素。因此,本研究基于以上水槽試驗(yàn)?zāi)Ｐ瓦M(jìn)行多次模擬,探究不同影響因素下泥石流對(duì)防護(hù)結(jié)構(gòu)的沖擊效應(yīng)。

3.3 水槽傾角的影響

為了比較不同水槽傾角對(duì)泥石流與防護(hù)結(jié)構(gòu)相互作用的影響,分別設(shè)置θ=0°、8°、15°、24°、30°進(jìn)行模擬。圖7為在模擬最終時(shí)刻t=6.0 s時(shí),不同水槽傾角下泥石流的最終形態(tài),對(duì)比發(fā)現(xiàn)水槽傾角越大,泥石流形成的堆積區(qū)體積就越大(θ=30°時(shí)除外),且只有θ=30°時(shí),泥石流撞擊防護(hù)結(jié)構(gòu)發(fā)生了溢流。這是因?yàn)樗蹆A角越大,泥石流的動(dòng)能越大,前端接觸到防護(hù)結(jié)構(gòu)的泥石流越不易發(fā)生回流。同時(shí),可以觀察到傾角越大,懸浮在流體上的玻璃珠數(shù)量越少,這是因?yàn)榱黧w對(duì)顆粒向上的浮力減小導(dǎo)致的。

圖7 不同水槽傾角下泥石流的最終形態(tài)Fig. 7 Final shape of debris flow at different inclinations of the flume

通過比較不同水槽傾角下,泥石流對(duì)防護(hù)結(jié)構(gòu)的沖擊力隨時(shí)間的變化情況,進(jìn)一步分析了水平傾角對(duì)于泥石流對(duì)防護(hù)結(jié)構(gòu)的沖擊效應(yīng)的影響,如圖8所示。顯然,沖擊力的變化趨勢(shì)基本相似,先迅速增加,達(dá)到峰值沖擊力后再慢慢下降,并逐漸趨于穩(wěn)定。但θ>15°時(shí),沖擊力曲線會(huì)出現(xiàn)2個(gè)明顯的波峰,且后者為峰值沖擊力,這種現(xiàn)象的發(fā)生是因?yàn)榈谝粋€(gè)波峰主要是由流體前部對(duì)防護(hù)結(jié)構(gòu)的撞擊導(dǎo)致,隨著流體的回流沖擊力減弱,后部流體的動(dòng)能仍較大,超過前部流體對(duì)其的阻礙,故出現(xiàn)了第二個(gè)波峰。而θ<15°時(shí),后部流體的動(dòng)能小于前部流體回流的阻礙,所以只存在一個(gè)波峰。同時(shí),注意到若以θ=0°時(shí)的峰值沖擊力50.26 N為基準(zhǔn),則θ=8°時(shí)峰值沖擊力為其的2.59倍,隨著傾角的增大,其倍數(shù)依次為4.75、8.88和11.45,傾角越大其峰值沖擊力的增長(zhǎng)率越大。且通過圖像可以看出,在趨于穩(wěn)定前,傾角越大,沖擊力波動(dòng)越大,這是因?yàn)榱黧w的反復(fù)沖擊與回流。穩(wěn)定后的沖擊力差異也很大,這是由于防護(hù)結(jié)構(gòu)底部堆積體的體積不一致導(dǎo)致的。提取防護(hù)結(jié)構(gòu)中心軸處水面高程隨時(shí)間的變化情況如圖9所示。由圖9可知,傾角越大,到達(dá)防護(hù)結(jié)構(gòu)的時(shí)間越短,并且沿防護(hù)結(jié)構(gòu)爬升的高度越大,最后穩(wěn)定時(shí)的水面高程也越大,但穩(wěn)定后水面高程的增加率是隨傾角的增大而逐漸減小的。同時(shí)水面高程達(dá)到穩(wěn)定的時(shí)間隨傾角的增大而變長(zhǎng),即泥石流與防護(hù)結(jié)構(gòu)的動(dòng)態(tài)作用時(shí)間隨傾角的增大而增大。

圖8 不同水槽傾角下沖擊力隨時(shí)間的變化Fig. 8 Relationship of impact force-time at different inclinations of the flume

圖9 不同水槽傾角下水面高程隨時(shí)間的變化Fig. 9 Relationship of water surface elevation-time at different inclinations of the flume

綜上所述,改變水槽傾角對(duì)于泥石流對(duì)防護(hù)結(jié)構(gòu)的沖擊效應(yīng)的影響是顯著的,不僅沖擊力會(huì)成倍地變化,當(dāng)防護(hù)結(jié)構(gòu)的高度低于泥石流沿防護(hù)結(jié)構(gòu)爬升的高度時(shí),還會(huì)有溢流現(xiàn)象發(fā)生。因此在工程實(shí)踐中,設(shè)計(jì)防護(hù)結(jié)構(gòu)前考慮泥石流易發(fā)處的地形傾角,模擬泥石流對(duì)防護(hù)結(jié)構(gòu)的沖擊特性是非常必要的。

3.4 泥石流體積的影響

在自然災(zāi)害中,泥石流的暴發(fā)規(guī)模不同,其造成的危害程度是不一致的,那么需要采取的預(yù)防措施便也不同。本文為了探究不同規(guī)模泥石流對(duì)防護(hù)結(jié)構(gòu)沖擊的影響,將泥石流中流體部分的體積V分別設(shè)置為0.024、0.048、0.072、0.096、0.120 m3進(jìn)行模擬。圖10比較t=6.0 s時(shí)刻不同流體體積的泥石流運(yùn)動(dòng)形態(tài),可見流體體積的增多導(dǎo)致固體體積分?jǐn)?shù)的減小,會(huì)促進(jìn)流體與顆粒間的相互作用,從而一定程度上增加玻璃珠的動(dòng)能。圖11進(jìn)一步比較不同流體體積下沖擊力隨時(shí)間的變化,可以看出,流體體積越大其接觸到防護(hù)結(jié)構(gòu)所需的時(shí)間會(huì)越短,達(dá)到穩(wěn)定沖擊力的時(shí)間越長(zhǎng)。同時(shí),流體規(guī)模越大,其對(duì)防護(hù)結(jié)構(gòu)的峰值沖擊力和靜態(tài)沖擊力都越大。隨著流體體積的增大,其峰值沖擊力增長(zhǎng)率也逐漸變大,相較V=0.024 m3分別增長(zhǎng)了0.79倍、1.64倍、3.20倍和4.56倍,但是穩(wěn)定下來的靜態(tài)沖擊力基本與體積的變化成正比。結(jié)合圖12水面高程隨時(shí)間的變化可以看到,雖然水面高程隨流體體積的增大而增高,但流體達(dá)到最大水面高程所需時(shí)間基本相同。而且最大水面高程的增長(zhǎng)率隨體積的增大表現(xiàn)為先增加后減小的趨勢(shì),分別為27%、47%、33%和21%;穩(wěn)定后水面高程的增長(zhǎng)率隨體積的增大則依次遞減的,分別為76%、21%、15%和13%。也就是說,在實(shí)際工程中為能更加經(jīng)濟(jì)環(huán)保地設(shè)計(jì)防護(hù)結(jié)構(gòu),可以考慮在保持一定高度時(shí),適當(dāng)增強(qiáng)防護(hù)結(jié)構(gòu)的強(qiáng)度以阻擋較大規(guī)模的泥石流。

圖10 不同流體體積下泥石流的最終形態(tài)Fig. 10 Final morphology of debris flow with different fluid volumes

圖11 不同流體體積下沖擊力隨時(shí)間的變化Fig. 11 Relationship of impact force-time with different fluid volumes

圖12 不同流體體積下水面高程隨時(shí)間的變化Fig. 12 Relationship of water surface elevation-time with different fluid volumes

3.5 石塊粒徑的影響

此外,本研究還考慮了其他條件一定的情況下,泥石流中石塊粒徑的改變對(duì)泥石流與防護(hù)結(jié)構(gòu)相互作用的影響。在確保石塊總體積不變時(shí),模擬球形石塊半徑R分別為10、20、30 mm時(shí),泥石流對(duì)防護(hù)結(jié)構(gòu)的沖擊。由圖13中t=6.0 s時(shí)刻不同石塊粒徑下泥石流的運(yùn)動(dòng)形態(tài)可以看出,相同體積下,石塊粒徑越小,其與流體的接觸面積越大,對(duì)流體的拖曳力的反作用力越強(qiáng),則越易阻礙流體的運(yùn)動(dòng)。而圖14展現(xiàn)的防護(hù)結(jié)構(gòu)所受沖擊力的變化則表明,由于固體體積分?jǐn)?shù)較小,石塊粒徑的改變對(duì)峰值沖擊力的影響不大,但是石塊粒徑越大,泥石流對(duì)防護(hù)結(jié)構(gòu)的靜態(tài)沖擊力越易波動(dòng),這可能是由于大塊石對(duì)流體回流的阻礙更弱導(dǎo)致的。同時(shí),由圖15可知,石塊粒徑越大,流體的最大水面高程越大,但石塊粒徑對(duì)穩(wěn)定后的水面高程基本沒有影響。因此,對(duì)于含有大塊石的泥石流,在設(shè)計(jì)防護(hù)結(jié)構(gòu)時(shí)可以適當(dāng)提高防護(hù)結(jié)構(gòu)的高度來對(duì)泥石流起到更好的阻擋效果。

圖13 不同石塊粒徑下泥石流的最終形態(tài)Fig. 13 Final morphology of the debris flow with different rock radii

圖14 不同石塊粒徑下沖擊力隨時(shí)間的變化Fig. 14 Relationship of impact force-time with different rock radii

圖15 不同石塊粒徑下水面高程隨時(shí)間的變化Fig. 15 Relationship of water surface elevation-time with different rock radii

4 結(jié)論

本文利用球體入水試驗(yàn)驗(yàn)證了SPH-DEM耦合方法在模擬流固兩相相互作用時(shí)的可靠性?；隍?yàn)證方法進(jìn)一步建立了水槽試驗(yàn)中含石塊的泥石流沖擊剛性防護(hù)結(jié)構(gòu)的模型,并研究了泥石流的運(yùn)動(dòng)過程及對(duì)剛性防護(hù)結(jié)構(gòu)的沖擊效應(yīng),分析了不同水槽傾角、不同泥石流流體體積和不同石塊粒徑的情況下沖擊力和水面高程的變化,得出以下結(jié)論:

1)通過與試驗(yàn)、理論結(jié)果的對(duì)比分析表明,采用耦合的SPH-DEM方法可以較好地模擬固液兩相的相互作用,其中SPH模擬流體,而DEM模擬固體。而且該方法對(duì)模擬含石塊泥石流的運(yùn)動(dòng)過程依舊可靠。

2)含石塊泥石流與防護(hù)結(jié)構(gòu)的相互作用可以分為3個(gè)階段:沖擊、爬升和堆積。當(dāng)水槽傾角較大時(shí),也可能會(huì)有溢流發(fā)生。泥石流中的黏性流體對(duì)防護(hù)結(jié)構(gòu)的沖擊力先迅速上升達(dá)到峰值沖擊力后逐漸下降并趨于穩(wěn)定,穩(wěn)定后的沖擊力不為0。而其中的石塊由于與流體的相互作用以及水槽底部的摩擦的影響,運(yùn)動(dòng)速度小于流體且會(huì)沿著堆積體向上攀爬。

3)水槽傾角、流體體積和石塊粒徑對(duì)于泥石流與防護(hù)結(jié)構(gòu)的相互作用影響顯著。隨著傾角由0°依次增長(zhǎng)為8°、15°、24°和30°,峰值沖擊力增長(zhǎng)倍數(shù)依次為2.59、4.75、8.88和11.45,傾角越大其峰值沖擊力的增長(zhǎng)率越大。但穩(wěn)定后水面高程的增加率是隨傾角的增大而逐漸減小的;隨著流體體積的增大,其峰值沖擊力增長(zhǎng)率也逐漸變大,相較V=0.024 m3分別增長(zhǎng)了0.79倍、1.64倍、3.20倍和4.56倍,但是穩(wěn)定下來的靜態(tài)沖擊力基本與體積的增長(zhǎng)成正比;石塊粒徑的改變對(duì)沖擊力的影響不大,但會(huì)影響峰值水面高程的形成,即石塊粒徑越大,峰值水面高程越高。