唐宇峰,史君林，李 濤

(1.四川輕化工大學(xué) 機(jī)械工程學(xué)院,四川 宜賓 644000；2.過程裝備與控制工程四川省高校重點(diǎn)實(shí)驗(yàn)室，四川 宜賓 644000)

0 引言

隨城鎮(zhèn)化快速發(fā)展，諸如道路施工、地基開挖、基礎(chǔ)夯土等第3方作用導(dǎo)致埋地管道損傷或破壞事故屢見不鮮，造成嚴(yán)重經(jīng)濟(jì)損失。其中，挖掘破壞導(dǎo)致管道損傷與破壞事故最為嚴(yán)重，且比例最大。目前，針對第3方挖掘造成管道破壞研究方法主要為有限元法：藍(lán)朝遜等[1]將管道破壞分為接觸、屈服、挖裂和挖穿4個(gè)階段，并采用有限元法分析管道彎頭在不同階段力學(xué)響應(yīng)；李軍等[2]采用有限元分析基于應(yīng)力與應(yīng)變2種失效準(zhǔn)則在管道失效過程中適用性；周立國等[3]采用ABAQUS分析挖掘載荷作用下燃?xì)夤艿懒W(xué)特征及動(dòng)力響應(yīng)；Neacsa[4]模擬分析不同挖掘及影響因素對管道損失影響。在挖掘作用下，管道受土體、斗齒等多方面綜合作用力，出現(xiàn)裂紋擴(kuò)展或破壞大變形，而基于網(wǎng)格的有限元法研究大變形問題時(shí)會(huì)產(chǎn)生網(wǎng)格畸變，即使采用網(wǎng)格重構(gòu)技術(shù)，也存在計(jì)算費(fèi)用昂貴，精度低甚至計(jì)算終止等問題。

光滑粒子流體動(dòng)力學(xué)(Smoothed Particle Hydrodynamics,SPH)是1種新型純拉格朗日無網(wǎng)格計(jì)算方法，在諸如流體、巖土、爆炸與沖擊等計(jì)算力學(xué)領(lǐng)域得到廣泛應(yīng)用[5-11]，粒子間無網(wǎng)格連接，避免因大變形導(dǎo)致的網(wǎng)格畸變，以表征大變形過程實(shí)時(shí)運(yùn)動(dòng)狀態(tài)。但在許多領(lǐng)域SPH仍處于起步階段，尤其在挖齒-土體-管道等多方接觸或沖擊作用等領(lǐng)域。因此，本文基于SPH理論，建立挖齒-土體-管道間相互接觸和沖擊力學(xué)模型，分析土體和管道在挖齒作用下動(dòng)態(tài)響應(yīng)。研究結(jié)果可為挖齒-土體-管道等多方接觸或沖擊領(lǐng)域研究提供新的計(jì)算方法。

1 SPH理論基礎(chǔ)及沖擊方程離散

1.1 SPH理論基礎(chǔ)

SPH法關(guān)鍵步驟為對任意宏觀變量函數(shù)f(x)的積分近似[12]，以及對連續(xù)積分表達(dá)式的粒子近似。任意宏觀變量函數(shù)f(x)的積分近似如式(1)所示：

(1)

(2)

連續(xù)積分表達(dá)式的粒子近似是將相關(guān)連續(xù)積分表達(dá)式轉(zhuǎn)化為支持域內(nèi)所有與粒子疊加求和的離散化方程，如式(3)所示：

(3)

(4)

在式(1)～(4)中，光滑函數(shù)Wij的選擇對離散近似起重要作用。本文選取3次光滑核函數(shù)，如式(5)所示：

(5)

式中：R為點(diǎn)x和x′處2粒子間相對距離，R=r/h=|x-x′|/h，h為光滑長度，m；αd在1，2，3維空間下分別等于1/h，15/7πh2，3/2πh3。

1.2 基于彈性及彈塑性動(dòng)力學(xué)理論的沖擊方程SPH離散

基于SPH基本理論，對質(zhì)量守恒方程、動(dòng)量守恒方程離散得到SPH形式下連續(xù)性方程及動(dòng)量方程，如式(6)～(7)所示：

(6)

(7)

式中：vi與vj分別表示粒子i與j速度，m/s；α，β代表坐標(biāo)方向；σαβ為總應(yīng)力張量，N；Fα表示外力引起加速度，如重力、齒輪齒條相互作用力及固壁邊界作用力等;∏ij為人工黏度[13]，防止粒子接近發(fā)生非物理穿透，如式(8)所示：

(8)

(9)

(10)

(11)

(12)

vij=vi-vj

(13)

xij=xi-xj

(14)

式中：vij代表粒子i與粒子j速度差，m/s；xij代表粒子i與粒子j距離,m。

同時(shí)，引入XSPH速度平均方法，如式(15)所示：

(15)

式中：ε為修正因子，取0.3。

PE管道為具有黏彈性的聚合物，其彈性模量隨時(shí)間增長而減少。由于挖齒與管道沖擊速度極快，可忽略時(shí)間影響，將模型簡化為分段彈性模型。對于土體模型，采用基于D-P屈服準(zhǔn)則的彈塑性本構(gòu)模型，如式(16)～(17)所示：

(16)

(17)

(18)

(19)

2 挖齒-土體-管道接觸模型

定義挖齒、土體、管道為不同部件。其中，挖齒作為剛體不考慮其變形，土體假設(shè)為連續(xù)均勻彈塑性介質(zhì)，管體因受力時(shí)間短而假設(shè)為分段彈性介質(zhì)。在SPH方法中，粒子一般通過核估計(jì)計(jì)算相互作用，如果不額外定義特殊接觸特性，當(dāng)各部件相互分離時(shí)，接觸面處粒子將產(chǎn)生與實(shí)際情況不符的剪切及拉伸力。因此，本文作以下3個(gè)假設(shè)：1)將土體、管體內(nèi)部視為連續(xù)體，挖齒-土體-管道接觸面處視為離散體。2)當(dāng)接觸面處不同部件粒子間距為各部件粒子平均間距1/2時(shí)，部件間粒子視為相互接觸。3)當(dāng)部件間粒子接觸時(shí)，給部件間粒子施加1個(gè)具有彈簧特性的相互作用力，該力與各部件變形及彈性模量有關(guān)，絕對值隨接觸粒子間距離減小而增大，且施加給各部件作用力大小相等，方向相反。當(dāng)部件間粒子間距大于部件粒子平均間距1/2時(shí)，表示部件未接觸，即作用力為0；當(dāng)部件接觸時(shí)，該力在2維條件下如式(20)所示：

(20)

基于SPH采用FORTRAN編制埋地燃?xì)夤艿劳诰蚱茐臄?shù)值模擬程序。其中，密度采用連續(xù)密度法求解，粒子信息更新采用跳蛙法，計(jì)算流程如圖1所示。

圖1 基于SPH的埋地燃?xì)夤艿劳诰蚱茐姆治鲇?jì)算流程

3 數(shù)值算例

本文以挖齒挖掘薄層土覆蓋管道為研究對象，管道外徑0.104 m，壁厚0.008 m，土體長1 m，高0.24 m，挖齒截面為扇形，圓心角35°。挖齒視為剛體，并以1 m/s速度垂直向下運(yùn)動(dòng)；土體彈性模量20 MPa，泊松比0.25，黏聚力12.38 KPa，摩擦角20°；管體彈性模量1.115 GPa，泊松比0.3。模型共離散為1 986個(gè)粒子，其中土粒子1 262個(gè)，初始間距0.017 m；管粒子480個(gè)，初始間距0.005 m；挖齒粒子110個(gè)，初始間距0.01 m；邊界粒子132個(gè)，初始間距0.012 5 m。時(shí)間步長4×10-7s，共計(jì)算8×105個(gè)時(shí)間步，即0.32 s。PE管道塑性階段應(yīng)力應(yīng)變見表1。挖齒-土體-管道模型初始狀態(tài)如圖2所示。模型分別在0.208，0.288，0.32 s發(fā)生接觸、挖裂和挖穿，模型發(fā)生接觸、挖裂、挖穿時(shí)粒子分布如圖3所示。

圖2 挖齒-土體-管道模型粒子初始狀態(tài)

表1 PE管塑性階段應(yīng)力應(yīng)變

由圖3可知，SPH法可呈現(xiàn)挖掘過程挖齒-土體-管道實(shí)時(shí)動(dòng)態(tài)及破壞狀態(tài)。在挖齒直接接觸管道前，管道變形較??；當(dāng)挖齒直接作用管道后，管道以接觸點(diǎn)為中心產(chǎn)生較大變形直至斷裂。

圖3 模型發(fā)生接觸、挖裂、挖穿時(shí)粒子分布

進(jìn)一步采用橢圓度對管道變形定量評價(jià)，如式(21)所示：

(21)

式中：Dmax、Dmin分別為管道最大外徑與最小外徑，m。挖齒挖掘土體及管道全過程，管道橢圓度隨挖掘時(shí)間變化關(guān)系如圖4所示。

圖4 管道橢圓度隨挖掘時(shí)間變化關(guān)系

由圖4可知，盡管在0.208 s前挖齒尚未直接接觸管道，但部分挖齒作用力由土體傳遞給管道，仍會(huì)使管道產(chǎn)生一定橢圓度變化，但此時(shí)橢圓度不足0.1，不會(huì)對管道產(chǎn)生破壞；當(dāng)挖齒與管道直接接觸后，管道橢圓度迅速增大并且逐漸挖裂，挖裂后管道抵抗變形能力大幅降低，直至挖穿時(shí)橢圓度約為0.45。

各挖掘階段管道環(huán)向應(yīng)力及環(huán)向應(yīng)變分布如圖5～6所示。由圖5～6可知，管道環(huán)向應(yīng)力隨挖齒向下作用而逐漸增大，管道由接觸至挖穿過程中，最大應(yīng)力應(yīng)變發(fā)生在挖齒與管道接觸部分，即管道正上方；當(dāng)管道發(fā)生挖裂時(shí)，最大應(yīng)力約35 MPa，等效應(yīng)變約0.045；當(dāng)管道挖穿時(shí)，環(huán)向應(yīng)力約65 MPa，等效應(yīng)變約0.08。此外，當(dāng)管道貫通后，管道上方約45°位置會(huì)出現(xiàn)較大應(yīng)力應(yīng)變，易發(fā)生進(jìn)一步破壞。

圖5 各挖掘階段管道環(huán)向應(yīng)力分布

圖6 各挖掘階段管道環(huán)向應(yīng)變分布

采用有限元法對相同模型模擬，以驗(yàn)證SPH法準(zhǔn)確性。但有限元法挖裂后網(wǎng)格趨于畸變，因此，以挖裂初環(huán)向應(yīng)力為例進(jìn)行對比，如圖7所示。有限元方法在挖裂階段得到環(huán)向最大應(yīng)力約35 MPa，與SPH法結(jié)果一致。

圖7 挖裂階段管道環(huán)向應(yīng)力分布(有限元)

4 結(jié)論

1)當(dāng)挖齒作用于土體而尚未作用于管道時(shí)，管道橢圓度受一定影響，但不會(huì)產(chǎn)生斷裂或破壞；當(dāng)挖齒直接作用于管道且繼續(xù)向管道方向挖掘時(shí)，管道會(huì)發(fā)生局部斷裂直至挖穿。

2)挖齒與管道直接接觸點(diǎn)為管道發(fā)生破壞或斷裂第1危險(xiǎn)點(diǎn)；當(dāng)管道繼續(xù)大變形后，其左右兩側(cè)受拉力作用可能產(chǎn)生斷裂破壞，為斷裂第2危險(xiǎn)點(diǎn)。

3)模型離散粒子數(shù)量及計(jì)算域等相關(guān)因素對計(jì)算精度產(chǎn)生一定影響。