基于SPH 方法一維顆粒仿真分析

程 麗 張 超， 肖道林， 張趙威

（1、沈陽大學機械工程學院，遼寧 沈陽110044 2、中國科學院沈陽自動化研究所，遼寧 沈陽110016）

由于科技的發(fā)展和人類對未知領域探索的需要，力學仿真廣泛存在于人類科學研究中。但是傳統(tǒng)力學仿真所用網格法具有一定局限性。[3-4]

近年，新興的光滑無網格的方法經過已經成熟，光滑粒子流體力學(SPH)是由核近似和粒子近似兩步組成的無網格拉格朗日粒子方法[5]。由于SPH 是在不同時刻計算該區(qū)域任意分布粒子上承載的物理量，所以具有極強的自適應性，并且各個物理量均為矢量。基于SPH 方法此種特性，對其一維進行研究可以為高維的工程問題提供參考。

1 SPH 模型

1.1 SPH 基本形式

在SPH 中，通過光滑核函數引進積分表示式[1]。W(xi-xj,h)稱為光滑核函數，該函數具有狄拉克函數的性質。其中h 表示光滑長度，? 為積分域，在SPH 中也叫支持域。A（x）j表示物理量A 在坐標xj處的值。x 是空間坐標向量，i 與j 表示粒子序號。

1.2 邊界條件的加入

實際工程中流體粒子運動會遇到無法穿透的邊界，模擬時為防止粒子穿透邊界，采用一種Lennard-Jones 邊界力模型

式中，i 和j 分別指邊界粒子和內部流體粒子，n1n2為常數，通常取12 和4。D 為最大速度的平方，r0為截止距離，通常取初始距離一半，粒子運動到邊界粒子r0范圍內將受到邊界力的作用。

1.3 時間積分的方法

為了求得粒子在每一個時間步運動狀態(tài)，引入蛙跳法。根據蛙跳積分的方案，在時間步n 有

2 對此次一維算例進行仿真

圖1 0.005s 粒子運動分布狀態(tài)

圖2 0.01975s 粒子運動分布狀態(tài)

圖3 0.025s 粒子分布狀態(tài)

圖4 速度隨時間變化曲線

仿真采用顆粒密度1.5kg/m，將每一個粒子看為半徑為0.01m。每一個時間步取得到廣泛應用的5*10-5s，共計算500 步。仿真結果如下，粒子呈現(xiàn)三角形為向右移動，呈現(xiàn)圓形為向左移動，正方形為邊界。

經過1、圖2、圖3 可以觀測到，0.005s 時刻所有粒子還是保持向右移動，0.01975s 開始出現(xiàn)由于邊界力作用而向左移動的粒子，最終500 步時即為0.025s 接近一半粒子向左移動，另一半粒子還未運動到邊界繼續(xù)向右移動，圖4 對以上加以驗證。

3 結論

最終得出SPH 可以模擬一維顆粒狀樣本變化，面對高維的工程問題只需提升維度如加上環(huán)境所需要的外力項與粘性項即可。