力學大會2015集交通船舶海洋與建筑工程學院800號

水下潰沙問題的SPH交通大學,船舶海洋與建筑,東川路800號,,:水下滑坡問題廣泛存在于海洋中，給對海底結(jié)構(gòu)物帶來極大危害。由于水環(huán)境的存在，水下滑坡呈現(xiàn)出與陸上滑坡顯著不同的特性。比如由于水的阻滯作用，水下滑坡過程要緩慢得多。土的含水量、內(nèi)摩擦角、滲透系數(shù)、強度等對水下滑坡過程均有重要影響。本文介紹了一種水土混合物流動的兩相流動模型。在這個模型中，水、土各占一定的體積分數(shù)，滿足各自的質(zhì)量和動量守恒方程。水、土的相互作用通過孔隙水壓力和阻力來模擬。土被看作為理想彈塑性物質(zhì)，服從Drur-ragr屈服準則和關(guān)聯(lián)流動律。由于土的大變形及自由液面的存在，基于網(wǎng)格的傳統(tǒng)數(shù)值方法難以模擬水下滑坡過程。本文建立了基于moothedrileyrodnmis（H）方法的數(shù)值方法來求解這種混合物流動。將其應用于水下直立沙柱的垮塌過程的數(shù)值模擬，得到了水、土各自的流場、應力場以及體積分數(shù)的變化規(guī)律，這些結(jié)論與實驗觀察吻合。研究工作表明，本文所提理論模型和數(shù)值方法，可以模擬水下滑坡過程，并可應用于其他復雜的水下飽和水、土混合流動問題。:水下,潰沙,海底滑坡廣泛存在于海洋中,對海底結(jié)構(gòu)物造成極大危害，長期以來一直受到人們的重視。海底滑坡過存在復雜的水/土耦合作用,這為研究海底滑坡帶來了極大的。由于水環(huán)境的存在，水下滑坡呈現(xiàn)出與陸上滑坡顯著不同的特性。比如，Meruaneetal.(2010)通過實驗和數(shù)值模擬發(fā)現(xiàn)水壓力可以促進固體從密實堆不可忽略。Rondonetal.(2011)實驗研究了水下沙柱的垮塌過程，他們發(fā)現(xiàn)初始體積分數(shù)對沙柱的垮塌、遷移和堆積過程有重要影響，而與沙柱的關(guān)系不大，始密實堆積的沙柱，其潰沙過程則要緩慢得多。Topinetal.(2012)應用接觸動力學最近，Sgeet.(2014)利用混合物模型和水平集（LlSt）方法，數(shù)值模擬研究了水下潰沙過程。對松散堆積的沙柱，他們的數(shù)值模擬結(jié)果與實驗觀測較為吻合，但對于初始密實堆積的沙柱，數(shù)值模擬結(jié)果卻不太理想。在潰沙過，沙體經(jīng)歷大變形，傳統(tǒng)基于網(wǎng)格的Euler方法難以處理這種問題，因此，??通訊作者:近期基于無網(wǎng)格的Laranin粒子法的數(shù)值方法在顆粒流動的模擬中流行起來。由GingolddMonahan(1977)和Lucy(1977)PH方法，是一種無網(wǎng)格的完全Laranin方法，具有處理大變形問題的能力。Buitl.(2007)率先將這種方法應用于水射流與土的相互作用問題上。他們提出了一個兩相模型，將土看作理想彈塑性材料，水/土耦合作用則通過阻力和孔隙水壓力來模擬。Buitl.(2008)進一步改進了土力學模型，實現(xiàn)了陸上滑坡的數(shù)值模擬。本文應用兩相混合物模型以及SPH方法研究水下潰沙過的水/土作用過程，主要目的是研究周圍水的存在對潰沙過程的影響。不同于Buietal.(2007)的模型，本文所采用的兩相混合物模型可以考慮潰沙過的兩相體積分數(shù)的變化，因水/土水/土耦合流動可以兩相混合物理論來描述。混合物理論(可參考Drew1983;Truesdell1984;Morland1992)假定,空間中的任一點在同一時刻總是被兩相的質(zhì)點以一定的體積分數(shù)占據(jù)。在這個假定下，可定義每相的部分密度（partialdensi-ty）ρα,部分速度vα,和部分應力σα，其中α=l,s分別表示液體和固體。各相均滿 +?·(ρvα)= ?(ραvα)+?·(ραvα?vα)=?·σα+ραg+f 這里，表示并矢積，g為重力加速度，fα是相互作用力。由定義,相互作用力的和為零,即fl+fs=0。利用質(zhì)量守恒方程(1),動量方程(2)可改寫為α α+ραg+f 其中

=?· 在混合物理論中，須對部分（partial）變量作出假設，以便將其與各相的真實（true）變量聯(lián)系起來。在標準的混合物理論中,一般假定部分速度與真實 vα=?α 這里帶的量均為真實量，或稱為內(nèi)部量（intrinsicvariable）；?α是α的體積分數(shù)且有?l?s=1?s σl=?pI+?lτ?l ,,用力fs(i.e.?fl)假設為fs=??s?p+Cd(vl? 真密度ρ?s假設為常數(shù).為使方程組封閉，我們還需給出本構(gòu)關(guān)系和相間阻力的表達在本文中，水被看作弱可壓縮的流體，其應力與應變率張量為線性關(guān)τ?αβ= 0p= ?l 000參考Buietal.(2008),將土看作理想彈塑性材料，無粘聚性，且服從Drucker-Prager屈服準則及關(guān)聯(lián)流動率（associate?owrule），則土的本構(gòu)關(guān)系可以寫為

s其中e˙αβ是偏應變率張量,G為剪切模量,K為體積模量。對平面應變問題,αθsαθ=

tan9+12tan2

率張量。塑性乘子變化率 通過下式計3αθKε˙γγ+

√θλ˙= s. 9α2K+Gθ本文中,阻力系數(shù)Cdin(6)可由定律推導Cd=?f f這里k為水力傳導系數(shù)，其量綱為速度的量綱[LT-1],γw=ρg為以部分密度表達的水的重度。k不但與固相有關(guān)，也與其中的液相有關(guān)，其值可由經(jīng)驗、物理實驗等得到(NieldandBejan,2006)。fSPH關(guān)于SPH方法的基本理論，可參考LiuandLiu(2003)，這里不再贅述。為簡化論述，我們約定下標a、b表示水粒子,下標i、j表示土粒子。在SPHDf Df a m +δh ab bab faab

abρbaijDρi=∑mvα?Wij+δhc∑ψαmj?Wij aiji i

jij si

ijρj這里vmn=vm—vn；ha、hi分別為粒子a和i的光滑長度（本文中取為相同；c、cs考聲速δ、δs為常數(shù)，本文設為0.1；ψijψij=2(ρi—ρj)|r|2 其中rij=xixj。ψαβ類推。式(13)、(14)τ了所謂δ-SPH方法，見Antuonoetal(2010)和Marroneetal(2011a,b)。τDf ??aDf

pa

pb+

+

aa

?x

a a miρ

+ i

(ip p

i =

i2i+j2j+

ij i i+

N

α∑

?Wia+∑mfiaW+ aρiρa aρ i這里ΠijandΠab為人工粘性項(參見Monaghan 由Df?f=Df?s可推導得水粒子的體積分數(shù)? Df

m

?Waj j =?sj

jajaa1?a=1? mjWaj s本文采用如下Wendlandquintic核函數(shù)作為光滑函數(shù){(1?R/2)4(2R+ 0≤R<W(R,h)=αd R≥ 其中R=|x?x′|/h；對二，αd=(π2壁面邊界采用所謂動態(tài)粒子方法處理(見Gomez-Gesteiraetal.2012)，即壁面在對固體應力進行顯式時間步進時，固體應力可能會躍出屈服面外，這是不允許的。本文采用了所謂所“處理”和“比例拉回”的辦法將屈服面外的應力狀態(tài)影射回屈服面上。在時間步長很小的情況下，這種處理不致引起太大誤差。詳細過程可參考ChndMizuo(1990)和Buitl.(2008)。L如圖1所示,長Li高Hi的沙柱浸沒在長方形的箱子中。在LFigure1:水下潰沙計算示意圖Table1:計算所需的物質(zhì)參數(shù)。沙柱初始形狀為6cm×8cm(LiEνκμ1.0×10?3Pa·?f計算中總共用了個流體粒子和.0025m。數(shù)值參數(shù)作如下設置：在狀態(tài)方程(8)中,B=1.4285×104;水中參考聲速為cf=10m/s;固體聲速為cs=215m/s;時間步長?t=5.0×10?6s。圖2給出了不同時刻水壓強的分布情況。可以看到在潰沙起始階段，由于剪逐漸增高，并最終大于周圍水壓力。還可以看到，在潰沙過，水面上也形成了波浪，在里來回晃蕩。上述數(shù)值模擬結(jié)果與Rondonetal.(2011)的實驗觀察相PressurePressuret=t=0t=0.1LowpressurePressurePressuret=t=0.2t=0.3PressurePressuret=t=0.4t=0.5Figure2:不同時刻的壓強場Rondonetal(2011)認為水下潰沙過程主要由初始體積分數(shù)決定。利用我們的界面存在一積分數(shù)的過渡區(qū)。需要的是，由于我們采用的體積分數(shù)重新初本文建立了兩相混合物PH計算模型，并應用其研究了水下潰沙問題。在該模型中，水、土分別滿足各自的質(zhì)量守恒方程和動量方程，體積分數(shù)也作為變量進入到方。土作為理想彈塑性物質(zhì)，且服從Drur-ragr屈服準則及關(guān)聯(lián)流動率。水/土作用通過阻力和孔隙水壓力來模擬。利用體積分數(shù)重初始化方法保持水土界面處的體積分數(shù)的突然變化。數(shù)值模擬結(jié)果與實驗觀察結(jié)果吻合，表明本文所提方法能夠用于研究復雜的水/土耦合作用問題。

t=t=0t=0.1

t=t=0.2t=0.3

t=t=0.4t=0.5

:solvedbymeansofSPHschemeswithnumericaldi?usiveterms.”ComputerPhysicsCommunications,181(3),532–549.Bui,H.H.,Fukagawa,R.,Sako,K.,andOhno,S.(2008).“LagrangianmeshBui,H.H.,Sako,K.,andFukagawa,R.(2007).“Numericalsimulationofsoil-waterinctionusingsmoothedparticlehydrodynamics(SPH)method.”JournalofTerramechanics,44(5),339–346.Chen,W.F.andMizuno,E.(1990).Nonlinearysisinsoilmechanics:theoryandimplementation.Elsevier:Amsterdam.Drew,D.A.(1983).“Mathematicalmodelingoftwo-phase?ow.”Annualreviewoffluidmechanics,15(1),261–291.Gingold,R.A.andMonaghan,J.J.(1977).“Smoothedparticlehydrodynamics:tronomicalSociety,181(3),375–389.1:Theoryandformulations.”Computers&Geosciences,48(0),289–299.Liu,G.R.andLiu,M.B.(2003).Smoothedparticlehydrodynamics:ameshparticlemethod.WorldScienti?c.Lucy,L.B.(1977).“Anumericalapproachtothetestingofthe?ssionTheastronomicaljournal,82,ouze,Marrone,S.,Colagrossi,A.,Antuono,M.,Lugni,C.,andTulin,M.(2011b).“A2D+tSPHmodeltostudythebreakingwavepatterngeneratedbyfastships.”JournalofFluidsandStructures,27(8),1199–1215.Meruane,C.,Tamburrino,A.,andRoche,O.(2010).“OntheroleoftheambientMonaghan,J.J.(1994).“Simulatingsurface?owswithSPH.”JournalofCom-putationalPhysics,110(2),399–406.Morland,L.(1992).“Flowofviscous?uidsthroughaporousdeformableSurveysinGeophysics,13(3),Rondon,L.,Pouliquen,O.,andAussillous,P.(2011).“Granular