力學(xué)大會(huì)2015集滑坡變形破壞的三維SPH模擬方法研究

2），，（中國(guó)力學(xué)流固耦合系統(tǒng)力學(xué)，滑坡是一種常見的嚴(yán)重地質(zhì)失穩(wěn)問題和變形破壞問題最受工程關(guān)注。而由于現(xiàn)有數(shù)值實(shí)驗(yàn)技術(shù)的限制，目前對(duì)邊坡失穩(wěn)后的變形特性，尤其是坡體的大變形和滑帶的三維擴(kuò)究相對(duì)較少。無網(wǎng)格方法，如方法，在描述滑坡發(fā)生后的變形特性上具有獨(dú)特優(yōu)勢(shì)P計(jì)算框架的三維土體大變形模型，采用了非關(guān)聯(lián)流動(dòng)法則下的-土體本構(gòu)來描述土體破壞的發(fā)生、發(fā)展到最終穩(wěn)定的全過程，并得到了實(shí)驗(yàn)數(shù)據(jù)的有效驗(yàn)證?；谠撃Ｐ?，本文對(duì)兩類帶特殊坡面的三維土坡失穩(wěn)案例進(jìn)行了模擬研究，分析了坡體最終形態(tài)帶形狀和滑移距離等特征，發(fā)現(xiàn)三維土坡變形破壞過程與二維條件有較大的差別實(shí)際工程中應(yīng)慮三維效應(yīng)。引的數(shù)值而導(dǎo)致的。為克服這些數(shù)值，越來越多的學(xué)者開始嘗試用基于連續(xù)介質(zhì)力學(xué)的無網(wǎng)格類節(jié)點(diǎn)關(guān)系，在此離散基礎(chǔ)上通過求解連續(xù)介質(zhì)的控制方程并施加各種邊界條件來描述粒子行蹤。這樣這類方法在處理大變形和自由表面問題時(shí)就具有明顯的優(yōu)勢(shì)[65無網(wǎng)格類方法中，光滑粒子動(dòng)力學(xué)(SPH)由于發(fā)展的相對(duì)比較成熟而應(yīng)用較廣。學(xué)者們利用方法模擬土體大變形問題種處理方法：一種是從流體力學(xué)出發(fā)，將ngh流動(dòng)模型引入控制方程（t等[1g等，等[2，ani等），擬流動(dòng)形式的滑坡，在二維和三維情況都展示出對(duì)滑坡運(yùn)動(dòng)分析的能力，但是它不能合理的模擬出滑帶的發(fā)育和裂紋的擴(kuò)展過程；另一種是從巖土力學(xué)出發(fā)i–Stokes方程結(jié)合（Bui等[3,4]，Chen等[5]），可以描述坡體從穩(wěn)定到失穩(wěn)后崩塌流動(dòng)再到最1資助：國(guó)家自然科學(xué)基金 );國(guó)家重點(diǎn)基礎(chǔ)研究發(fā)展計(jì)劃

ρσ1 σ fα ρ其中x,y,z是空間坐標(biāo)；α,β滿足約定；ρ,v分別代表土體密度、速度，fα是外力引起的加速度(本文指重力加速度g)；σ為應(yīng)力張量，可由下式彈塑性應(yīng)αγ γβ αγ γβ γγ αβ σ σ i

iK

，G 3(1 其是塑性因子變化率，它決定著塑性應(yīng)變的大小；Sαβ為偏應(yīng)力率張量；KG分別是體積模量和剪切模量；υ為E為楊mn為啞指標(biāo)；δ為Kronecker’sdelta函數(shù)；為總應(yīng)變率張量。這里采用非流動(dòng)法則下的Prager本構(gòu)模型，屈服函數(shù)f(I1,J2)和塑性位勢(shì)函數(shù)g(I1,J2)分別可以寫f(I1,J2) g(I1,J2) 其中I1和J2分別是第一和第二應(yīng)力不變量；αφ和kc為Prager常量，它們跟Mohr–Coulomb破壞準(zhǔn)則的粘聚力c和內(nèi)摩擦角φ有關(guān)；αψ采用與SPH模擬基于連續(xù)介質(zhì)力學(xué)，它將整個(gè)計(jì)算域離散成有限的但足夠多的粒子來表

Nm(vαv

j

j ii

σ im( -

n

xβ i

j i αγ γβ

αβ σ σ

Kε

3αKδ i jj

f(I1,J2) f(I1,J2)N N

α

N j

β

(vv) (vv) x j j

N

α

β ω (vv) (vv) x j j Wij為核函數(shù)，為抑制計(jì)算過程中的張力不穩(wěn)定，這里采用構(gòu)造的二階導(dǎo)數(shù)負(fù)的核函數(shù)。歸一化系數(shù)αD在一維、二維和三的取值分別為1/(7h),1/(3πh2)F15/(62πh3) F F ，指數(shù)nW(0,

用來控制人工應(yīng)力的大小，本文取n2.55 W(x, W(x, 00ε為0~1之間的常數(shù)，本文取ε=0.3

通過二維典型算例來檢驗(yàn)所建立的土體大變形模型的有效性。ui等[41以圓柱形的鋁棒代替散粒體，做了近似二維的散粒體傾倒試驗(yàn)為方便觀測(cè)散粒體的變形，試驗(yàn)中對(duì)散粒體布置了0.02m0.02m的柵格，土槽的尺寸為0.6m1.2m們對(duì)此過程進(jìn)擬為一矩形土體，土體尺寸0.2m.1m圖1展示了散粒體的最終形態(tài)們將模擬的結(jié)果與實(shí)驗(yàn)觀測(cè)的散粒體外形輪廓和柵格變形對(duì)好。 （a）試驗(yàn)觀 （b）SPH模擬結(jié)1外形輪廓與內(nèi)部變形的對(duì)本文模型同樣可用于粘性土的研究，Bui[97]模擬了二維均質(zhì)土坡由強(qiáng)度折減誘發(fā)穩(wěn)的過程。這里我們用本文的模型重復(fù)了他的工作。本文采用c、φ同比例折減的方法折減后的剪切強(qiáng)度參數(shù)cr和φrc

和φarctan(

SRFa.滑坡的初始形態(tài) b.模擬滑帶位置與Bui以及圓弧法結(jié)果2二維粘性土滑帶的模擬結(jié)果本例邊坡的幾何尺寸如圖2a所示。圖2b展示了邊坡滑動(dòng)后的最終狀態(tài)，這與Bui表1剪脹實(shí)際邊坡的坡面大多都不是平整的，總是呈現(xiàn)一定的弧度。根據(jù)二維平面應(yīng)變理土坡的幾何尺寸：(a)為俯視圖，(b)為中間斷面M，M1在坡頂，M2在坡腳。黑實(shí)線、紅3凹坡、平坡和凸坡三種圖5曲面三維土坡變形后的最終形態(tài)(a)b)(c)分別針對(duì)凹坡、平坡和凸坡三種①②③④分別為位移圖、累積塑性應(yīng)變圖、滑坡體的形態(tài)和斷面N的塑性區(qū)5(④)中可以看出，在橫斷面N處，凹坡的滑帶比凸坡薄很多，而平坡介于兩者之間。如6M1的垂向位移和坡腳M2的水平位移，我們做了定量的對(duì)比，如圖76最終狀態(tài)下二維情況(a)與三維凹坡(b)、平坡(c)和凸坡(d)在斷面M上的累積塑性應(yīng)變對(duì)比3D3D3D3DM2M2位位位位位位76543210--1012345678位位

63D3D3D3D3DM1位位位位位43210-1012345678圖臺(tái)兩種典型狀況（圖8），拐角是90°，坡度是45°。它們的幾何尺寸由圖9給出：(a)凹臺(tái)；(b)凸臺(tái)。P1和P2分別在斷面P的坡頂和坡腳；Q1和Q2分別在斷面Q的坡圖9直角三維土坡的俯視圖。(a)凹臺(tái)；(b)P1和P2分別在斷面P的坡頂和坡腳；Q1和Q2分別在斷面Q的坡頂和坡腳10直角三維土坡變形后的最終形態(tài)(a)(b)分別針對(duì)凹臺(tái)和凸臺(tái)。①②③分別為位移圖、累積塑性應(yīng)變圖和滑坡體的形態(tài)圖我們分別截取斷面P和Q的數(shù)據(jù)并同二維情況作了對(duì)比，如位位位位A1位位位位位43210-(a)坡頂

3D3D2D2D-1012345678

位位P位位P2位位位位位位765 3D3D 2D 2D凸0--1012345678位位(b)坡腳12斷面P上坡頂P1和坡腳P2的位移。(a)坡頂P1(b)坡腳位位Q1位位位Q1位位位位位位B2位位B2位位位位位 - -

3D3D2D2D-1012345678位位

-1012345678位位(a)坡頂 (b)坡腳圖13斷面Q上坡頂Q1和坡腳Q2的位移。(a)坡頂 (b)坡腳對(duì)比。從中可以看出，在坡頂P1處，凹臺(tái)幾乎沒有動(dòng)，凸臺(tái)的垂向位移正好介于截面P圖13是針對(duì)斷面上的坡頂1的垂向位移和坡腳2的水平位移同二維1處的垂向位移和坡腳2處的水位位移各自都介于截面和截面兩個(gè)相對(duì)應(yīng)的二維情況之間；而凸臺(tái)在坡頂1處的垂向位移和坡腳2處的水位位移都比截面相對(duì)應(yīng)的二維坡腳處都比較安全，而凸臺(tái)拐角處的坡頂和凹臺(tái)的側(cè)坡更容易發(fā)生變形。3 下的Drucker-帶形狀和滑移距離方面做了分析研究，并同二維情況作了對(duì)比。結(jié)有明顯的三維效應(yīng)；坡面形態(tài)影響坡體最終形態(tài)帶形狀和滑移距離；二維和三維模擬的結(jié)果存在差異，差異與坡面形態(tài)有關(guān)針對(duì)凹、平、凸三種典型坡面，隨著坡面凹陷程度加，滑帶的凹陷程度也越明顯；坡面越凸越容易發(fā)生失穩(wěn)破壞，坡頂?shù)某两稻嚯x越大，坡腳的平移距離也越大。③針對(duì)帶拐角的凸臺(tái)和凹臺(tái)，三維效應(yīng)更加明顯PastorM,ManzD,FernándezMerodoJA.FromSolidstoFluidizedSoils:DiffuseFailureMechanismsinwithApplicationstoFastCatastrophicLandslides.GranularMatter,2010,12(3): HuM,LiuMB,XieMW.Three-DimensionalRun-OutysisandPredictionofFlow-LikeLandslidesUsingSmoothedParticleHydrodynamics.EnvironmentalEarthSciences,2014BuiHH,FukagawaR,SakoK.SlopeStabilityysisandDiscontinuousSlopeFailureSimulationbyElasto-sticSmoothedParticleHydrodynamics(SPH).Geotechnique,2011,61(7):565-574BuiHH,SakoK.SPH-BasedNumericalSimulationsforLargeDeformationofGeomaterialConsideringSoil-Structure ction.2008ChenW,QiuT.SimulationofEarthquake-InducedSlopeDeformationUsingSPHMethod.InternationalJournalforNumericalandyticalMethodsinGeomechanics,2014,38(3):297-330RESEARCHONTHREE-DIMENSIONALSIMULATIONMETHODOFSLOPEFAILUEWITHSPHMETHODANYiWUQiangLIU(1KeyLaboratoryforMechanicsinFluidSolidCouplingSystems,InstituteofMechanics,CAS,No.15BeisihuanxiRoad,Beijing100190,)ASPHmodelforlargedeformationofsoilwasdevelopedinthispaper,and3Dsimulationofsoildeformationandfailureprocesswasimplemented.TheDrucker-PragermodelundertheunassociatedflowrulewasimplementedintotheSPH-codetodescribetheelastic-sticsoilbehaviorwhiletherocksweresimulatedasrigidbodiesbyusingclassicalrigidmotionequations.Firstly,thepost-failureoftheearth-rockaggregatelandslidewassimulatedandtheinfluenceontheslidedistancewerediscussedintermsofslopegradientandrocksize.Secondly,thismodelwasappliedtotwospecial3Dslopeswithdifferentgeometricconfigurations,includingcurvingslope