力学大会2015集交通船舶海洋与建筑工程学院800号

水下溃沙问题的SPH交通大学,船舶海洋与建筑,东川路800号,,:水下滑坡问题广泛存在于海洋中，给对海底结构物带来极大危害。由于水环境的存在，水下滑坡呈现出与陆上滑坡显著不同的特性。比如由于水的阻滞作用，水下滑坡过程要缓慢得多。土的含水量、内摩擦角、渗透系数、强度等对水下滑坡过程均有重要影响。本文介绍了一种水土混合物流动的两相流动模型。在这个模型中，水、土各占一定的体积分数，满足各自的质量和动量守恒方程。水、土的相互作用通过孔隙水压力和阻力来模拟。土被看作为理想弹塑性物质，服从Drur-ragr屈服准则和关联流动律。由于土的大变形及自由液面的存在，基于网格的传统数值方法难以模拟水下滑坡过程。本文建立了基于moothedrileyrodnmis（H）方法的数值方法来求解这种混合物流动。将其应用于水下直立沙柱的垮塌过程的数值模拟，得到了水、土各自的流场、应力场以及体积分数的变化规律，这些结论与实验观察吻合。研究工作表明，本文所提理论模型和数值方法，可以模拟水下滑坡过程，并可应用于其他复杂的水下饱和水、土混合流动问题。:水下,溃沙,海底滑坡广泛存在于海洋中,对海底结构物造成极大危害，长期以来一直受到人们的重视。海底滑坡过存在复杂的水/土耦合作用,这为研究海底滑坡带来了极大的。由于水环境的存在，水下滑坡呈现出与陆上滑坡显著不同的特性。比如，Meruaneetal.(2010)通过实验和数值模拟发现水压力可以促进固体从密实堆不可忽略。Rondonetal.(2011)实验研究了水下沙柱的垮塌过程，他们发现初始体积分数对沙柱的垮塌、迁移和堆积过程有重要影响，而与沙柱的关系不大，始密实堆积的沙柱，其溃沙过程则要缓慢得多。Topinetal.(2012)应用接触动力学最近，Sgeet.(2014)利用混合物模型和水平集（LlSt）方法，数值模拟研究了水下溃沙过程。对松散堆积的沙柱，他们的数值模拟结果与实验观测较为吻合，但对于初始密实堆积的沙柱，数值模拟结果却不太理想。在溃沙过，沙体经历大变形，传统基于网格的Euler方法难以处理这种问题，因此，∗∗通讯作者:近期基于无网格的Laranin粒子法的数值方法在颗粒流动的模拟中流行起来。由GingolddMonahan(1977)和Lucy(1977)PH方法，是一种无网格的完全Laranin方法，具有处理大变形问题的能力。Buitl.(2007)率先将这种方法应用于水射流与土的相互作用问题上。他们提出了一个两相模型，将土看作理想弹塑性材料，水/土耦合作用则通过阻力和孔隙水压力来模拟。Buitl.(2008)进一步改进了土力学模型，实现了陆上滑坡的数值模拟。本文应用两相混合物模型以及SPH方法研究水下溃沙过的水/土作用过程，主要目的是研究周围水的存在对溃沙过程的影响。不同于Buietal.(2007)的模型，本文所采用的两相混合物模型可以考虑溃沙过的两相体积分数的变化，因水/土水/土耦合流动可以两相混合物理论来描述。混合物理论(可参考Drew1983;Truesdell1984;Morland1992)假定,空间中的任一点在同一时刻总是被两相的质点以一定的体积分数占据。在这个假定下，可定义每相的部分密度（partialdensi-ty）ρα,部分速度vα,和部分应力σα，其中α=l,s分别表示液体和固体。各相均满 +∇·(ρvα)= ∂(ραvα)+∇·(ραvα⊗vα)=∇·σα+ραg+f 这里，表示并矢积，g为重力加速度，fα是相互作用力。由定义,相互作用力的和为零,即fl+fs=0。利用质量守恒方程(1),动量方程(2)可改写为α α+ραg+f 其中

=∇· 在混合物理论中，须对部分（partial）变量作出假设，以便将其与各相的真实（true）变量联系起来。在标准的混合物理论中,一般假定部分速度与真实 vα=˜α 这里带的量均为真实量，或称为内部量（intrinsicvariable）；ϕα是α的体积分数且有ϕlϕs=1˜s σl=−pI+ϕlτ˜l ,,用力fs(i.e.−fl)假设为fs=−ϕs∇p+Cd(vl− 真密度ρ˜s假设为常数.为使方程组封闭，我们还需给出本构关系和相间阻力的表达在本文中，水被看作弱可压缩的流体，其应力与应变率张量为线性关τ˜αβ= 0p= ˜l 000参考Buietal.(2008),将土看作理想弹塑性材料，无粘聚性，且服从Drucker-Prager屈服准则及关联流动率（associateflowrule），则土的本构关系可以写为

s其中e˙αβ是偏应变率张量,G为剪切模量,K为体积模量。对平面应变问题,αθsαθ=

tan9+12tan2

率张量。塑性乘子变化率 通过下式计3αθKε˙γγ+

√θλ˙= s. 9α2K+Gθ本文中,阻力系数Cdin(6)可由定律推导Cd=ϕf f这里k为水力传导系数，其量纲为速度的量纲[LT-1],γw=ρg为以部分密度表达的水的重度。k不但与固相有关，也与其中的液相有关，其值可由经验、物理实验等得到(NieldandBejan,2006)。fSPH关于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为常数，本文设为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可推导得水粒子的体积分数ϕ Df

m

∂Waj j =−sj

jajaa1ϕa=1− mjWaj s本文采用如下Wendlandquintic核函数作为光滑函数{(1−R/2)4(2R+ 0≤R<W(R,h)=αd R≥ 其中R=|x−x′|/h；对二，αd=(π2壁面边界采用所谓动态粒子方法处理(见Gomez-Gesteiraetal.2012)，即壁面在对固体应力进行显式时间步进时，固体应力可能会跃出屈服面外，这是不允许的。本文采用了所谓所“处理”和“比例拉回”的办法将屈服面外的应力状态影射回屈服面上。在时间步长很小的情况下，这种处理不致引起太大误差。详细过程可参考ChndMizuo(1990)和Buitl.(2008)。L如图1所示,长Li高Hi的沙柱浸没在长方形的箱子中。在LFigure1:水下溃沙计算示意图Table1:计算所需的物质参数。沙柱初始形状为6cm×8cm(LiEνκµ1.0×10−3Pa·˜f计算中总共用了个流体粒子和.0025m。数值参数作如下设置：在状态方程(8)中,B=1.4285×104;水中参考声速为cf=10m/s;固体声速为cs=215m/s;时间步长∆t=5.0×10−6s。图2给出了不同时刻水压强的分布情况。可以看到在溃沙起始阶段，由于剪逐渐增高，并最终大于周围水压力。还可以看到，在溃沙过，水面上也形成了波浪，在里来回晃荡。上述数值模拟结果与Rondonetal.(2011)的实验观察相PressurePressuret=t=0t=0.1LowpressurePressurePressuret=t=0.2t=0.3PressurePressuret=t=0.4t=0.5Figure2:不同时刻的压强场Rondonetal(2011)认为水下溃沙过程主要由初始体积分数决定。利用我们的界面存在一积分数的过渡区。需要的是，由于我们采用的体积分数重新初本文建立了两相混合物PH计算模型，并应用其研究了水下溃沙问题。在该模型中，水、土分别满足各自的质量守恒方程和动量方程，体积分数也作为变量进入到方。土作为理想弹塑性物质，且服从Drur-ragr屈服准则及关联流动率。水/土作用通过阻力和孔隙水压力来模拟。利用体积分数重初始化方法保持水土界面处的体积分数的突然变化。数值模拟结果与实验观察结果吻合，表明本文所提方法能够用于研究复杂的水/土耦合作用问题。

t=t=0t=0.1

t=t=0.2t=0.3

t=t=0.4t=0.5

:solvedbymeansofSPHschemeswithnumericaldiffusiveterms.”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-phaseflow.”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.WorldScientific.Lucy,L.B.(1977).“AnumericalapproachtothetestingofthefissionTheastronomicaljournal,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).“SimulatingsurfaceflowswithSPH.”JournalofCom-putationalPhysics,110(2),399–406.Morland,L.(1992).“FlowofviscousfluidsthroughaporousdeformableSurveysinGeophysics,13(3),Rondon,L.,Pouliquen,O.,andAussillous,P.(2011).“Granular