周 浩，徐志宏，鄒 順，湯文輝

(1.國(guó)防科技大學(xué) 理學(xué)院， 湖南 長(zhǎng)沙 410073； 2.國(guó)防科技大學(xué) 計(jì)算機(jī)學(xué)院， 湖南 長(zhǎng)沙 410073)

SPH并行計(jì)算方案及其在自由表面流動(dòng)模擬中的應(yīng)用*

周 浩1，徐志宏1，鄒 順2，湯文輝1

(1.國(guó)防科技大學(xué) 理學(xué)院， 湖南 長(zhǎng)沙 410073； 2.國(guó)防科技大學(xué) 計(jì)算機(jī)學(xué)院， 湖南 長(zhǎng)沙 410073)

針對(duì)光滑粒子動(dòng)力學(xué)主要計(jì)算量是近鄰粒子搜索這一特點(diǎn)，提出了一種基于粒子分解的光滑粒子動(dòng)力學(xué)并行計(jì)算方案。利用該方案可以方便地將任意串行光滑粒子動(dòng)力學(xué)代碼并行計(jì)算，而且每一個(gè)時(shí)間步內(nèi)的信息傳遞量只和粒子總數(shù)有關(guān)，而和粒子的分布無(wú)關(guān)，因而特別適合于自由表面流動(dòng)等大變形問(wèn)題的并行數(shù)值模擬。對(duì)一個(gè)粒子總數(shù)為40萬(wàn)的三維潰壩問(wèn)題的模擬結(jié)果表明：此方案能達(dá)到的最大加速比約為16，這一結(jié)果可能比空間分解方案(不考慮動(dòng)態(tài)負(fù)載均衡)更優(yōu)。

光滑粒子動(dòng)力學(xué)；并行方案；粒子分解；空間分解

(1.CollegeofScience,NationalUniversityofDefenseTechnology,Changsha410073,China；

2.CollegeofComputer,NationalUniversityofDefenseTechnology,Changsha410073,China)

多數(shù)數(shù)值方法在處理大變形時(shí)存在一定困難。光滑粒子動(dòng)力學(xué)(SmoothedParticleHydrodynamics,SPH)[1-2]方法作為一種無(wú)網(wǎng)格的拉格朗日方法，在處理自由表面流動(dòng)等大變形問(wèn)題時(shí)具有很大的優(yōu)勢(shì)，近年來(lái)其應(yīng)用已經(jīng)擴(kuò)展到很多領(lǐng)域。又因?yàn)樗睦窭嗜仗匦?，界面附近物質(zhì)的運(yùn)動(dòng)即代表了界面的運(yùn)動(dòng)，因此處理多相流也比較自然，且不需要額外的界面追蹤算法。但是，SPH方法比有限差分、有限元等網(wǎng)格方法的計(jì)算量要大，因此在模擬規(guī)模較大的問(wèn)題時(shí)經(jīng)常需要并行計(jì)算。

SPH的計(jì)算過(guò)程和分子動(dòng)力學(xué)(MolecularDynamics,MD)類(lèi)似，所以其并行也和MD類(lèi)似。Plimpton[3]總結(jié)了3種不同的MD并行方案，即空間分解，粒子分解和力分解。目前，文獻(xiàn)中絕大多數(shù)SPH方法并行采用的都是空間分解方案，因?yàn)槠湫畔鬟f量小，其缺點(diǎn)是編程實(shí)現(xiàn)復(fù)雜。

空間分解方案又可以分為靜態(tài)和動(dòng)態(tài)兩種。靜態(tài)空間分解時(shí)，每個(gè)處理器分配的空間始終固定。當(dāng)所有粒子始終在計(jì)算域近似均勻分布時(shí)，這種方法能夠得到極高的加速比。但是，在自由表面流動(dòng)等實(shí)際問(wèn)題中，粒子往往在空間中分布非常不均勻，這會(huì)使得負(fù)載嚴(yán)重不平衡，導(dǎo)致加速比較小。動(dòng)態(tài)空間分解時(shí)，需要根據(jù)負(fù)載平衡原則每隔一段時(shí)間動(dòng)態(tài)調(diào)整各個(gè)處理器空間的大小和位置，來(lái)保證所有處理器的計(jì)算量大致相等。調(diào)整處理器空間需要額外的計(jì)算時(shí)間，因此動(dòng)態(tài)調(diào)整的頻率要適當(dāng)，而且調(diào)整算法也必須高效。

動(dòng)態(tài)空間分解雖然能夠高效地模擬大變形問(wèn)題，但是實(shí)現(xiàn)起來(lái)非常復(fù)雜。如Ferrari等[4]在采用并行SPH方法模擬三維潰壩問(wèn)題時(shí)，采用第三方軟件包METIS來(lái)動(dòng)態(tài)調(diào)整處理器空間。在METIS中，首先需要在空間布滿(mǎn)背景網(wǎng)格，再根據(jù)每個(gè)小網(wǎng)格中粒子數(shù)的多少來(lái)給每個(gè)小網(wǎng)格一個(gè)權(quán)重，來(lái)實(shí)現(xiàn)每個(gè)處理器空間大小不等但是計(jì)算量近似相等，用以實(shí)現(xiàn)負(fù)載平衡。Marrone等[5]采用并行SPH方法模擬行船產(chǎn)生的波浪及其破碎時(shí)，只在主要流動(dòng)方向上動(dòng)態(tài)調(diào)整，因?yàn)檫@能極大地簡(jiǎn)化編程，并且指出存在主要流動(dòng)方向的特征越明顯，并行模擬的加速比越高。并行SPH程序JOSEPHINE[6]和開(kāi)源軟件parallelSPHYSICSv2.0在模擬二維自由表面流動(dòng)時(shí)，也都只在一個(gè)主要流動(dòng)方向上實(shí)現(xiàn)了動(dòng)態(tài)調(diào)整，因此，當(dāng)流動(dòng)有多個(gè)主要方向時(shí)，加速比會(huì)下降?？偟膩?lái)說(shuō)，動(dòng)態(tài)空間分解的實(shí)現(xiàn)復(fù)雜度嚴(yán)重限制了其在SPH并行模擬中的廣泛使用。

Plimpton描述了MD的粒子分解并行方案，其基本思想是每個(gè)處理器負(fù)責(zé)維護(hù)固定的某一部分粒子。這種方案的最大優(yōu)點(diǎn)在于并行很容易在串行代碼的基礎(chǔ)上實(shí)現(xiàn)[3]。但是，由于MD中的分子由若干原子組成，分子具有一定內(nèi)部結(jié)構(gòu)，需要考慮其旋轉(zhuǎn)等因素，而且與每個(gè)分子有相互作用的分子列表可以若干時(shí)間步更新一次，所以其主要計(jì)算量在于分子間相互作用力的計(jì)算。而SPH中粒子間的相互作用力相對(duì)簡(jiǎn)單，而且一般每個(gè)時(shí)間步都需要更新每個(gè)粒子的近鄰粒子，因此其主要計(jì)算量在于近鄰粒子搜索[7]。

針對(duì)SPH方法的上述特點(diǎn)，提出了一種簡(jiǎn)單的粒子并行方案。這種方案能夠達(dá)到可觀加速比，有效節(jié)約計(jì)算時(shí)間，而且方案實(shí)現(xiàn)非常簡(jiǎn)單，負(fù)載均衡，每個(gè)時(shí)間步的信息傳遞量與粒子分布無(wú)關(guān)，可以十分方便地采取任意奇數(shù)個(gè)處理器并行計(jì)算。但是，這種方案的信息傳遞量較大，不適合上百個(gè)處理器的大規(guī)模并行計(jì)算。

以二維圓柱繞流和三維潰壩問(wèn)題為例，比較了這種粒子并行方案和空間分解方案(不考慮動(dòng)態(tài)負(fù)載均衡)的計(jì)算時(shí)間或加速比。在這兩個(gè)算例中，前者粒子分布始終非常均勻，而后者則是典型的大變形問(wèn)題，粒子分布嚴(yán)重不均勻。這兩個(gè)算例的對(duì)比可以明顯地展示出粒子并行方案只適合于粒子分布嚴(yán)重不均勻問(wèn)題的并行模擬。SPH靜態(tài)空間分解方案在開(kāi)源軟件平臺(tái)OpenFOAM的基礎(chǔ)上實(shí)現(xiàn)，粒子分解方案采用Fortran90語(yǔ)言加上訊息傳遞接口(MessagePassingInterface，MPI)消息傳遞庫(kù)實(shí)現(xiàn)。本文所有模擬都是在廣州超算中心的天河二號(hào)計(jì)算機(jī)上完成。

1 控制方程及其SPH離散

假設(shè)系統(tǒng)中沒(méi)有溫度梯度，黏性流動(dòng)的控制方程為如下N-S方程

(1)

(2)

其中:ρ，u,p,t分別為密度、速度、壓力和時(shí)間；υ2u為黏性力；υ為運(yùn)動(dòng)黏性系數(shù)；fe為外力。控制方程式(1)和式(2)右邊各導(dǎo)數(shù)項(xiàng)的SPH離散形式如下：

(-ρ

(3)

(4)

(5)

其中，式(5)為Morris基于有限差分和SPH離散提出的黏性力計(jì)算模型[8]。本文模擬中時(shí)間積分采用蛙跳法。

對(duì)不可壓縮流體采用弱可壓來(lái)近似，通過(guò)將流體的密度變化控制在1%以?xún)?nèi)來(lái)近似實(shí)現(xiàn)不可壓[9]。本文中采用如下最常用的物態(tài)方程其中，γ=7，c0取粒子最大速度的10倍。

(6)

2 SPH并行方案

現(xiàn)階段文獻(xiàn)中SPH并行絕大多數(shù)采用的是空間分解方案。本文實(shí)現(xiàn)并比較了空間分解和粒子分解這兩種并行方案各自的優(yōu)缺點(diǎn)。

2.1 空間分解方案

空間分解方案是目前使用最普遍的方案，其基本思想就是首先將整個(gè)計(jì)算域劃分為若干不重疊的子域，然后每個(gè)處理器負(fù)責(zé)一個(gè)子域，如圖1所示。

圖1 空間分解示意圖Fig.1 Domain decomposition

當(dāng)計(jì)算滿(mǎn)足局域性(任何一點(diǎn)在下一時(shí)刻的值只和包括這點(diǎn)的某局部空間相關(guān))時(shí)，只有在子域之間的交界處才需要傳遞信息，如圖1中的實(shí)心粒子部分。因此，劃分子域的一個(gè)原則就是讓交界面積盡量小，這代表了信息傳遞少，即并行帶來(lái)的額外時(shí)間開(kāi)銷(xiāo)小。另一個(gè)原則是盡量保持每個(gè)子域的計(jì)算量相等，即負(fù)載均衡，這樣才能達(dá)到較高加速比。

在歐拉類(lèi)網(wǎng)格方法中，往往是網(wǎng)格的數(shù)量代表了計(jì)算量。由于網(wǎng)格是不動(dòng)的，所以負(fù)載均衡比較容易實(shí)現(xiàn)。但是SPH方法屬于拉格朗日方法，粒子數(shù)量大致代表了計(jì)算量，大變形條件下，粒子占據(jù)的空間可能發(fā)生很大變化，如模擬潰壩、氣體擴(kuò)散、超高速碰撞等問(wèn)題。這會(huì)導(dǎo)致每個(gè)子域的粒子數(shù)量嚴(yán)重不平衡，導(dǎo)致并行效率低下。一種做法是采取動(dòng)態(tài)負(fù)載均衡，即模擬過(guò)程中如果發(fā)現(xiàn)負(fù)載不均衡，就按照一定策略將空間重新分解，如開(kāi)源軟件parallelSPHYSICSv2.0的使用手冊(cè)中給出的二維潰壩問(wèn)題的模擬結(jié)果，其動(dòng)態(tài)負(fù)載均衡如圖2所示。

圖2 動(dòng)態(tài)負(fù)載均衡Fig. 2 Dynamic load balancing

圖2中有3個(gè)垂直方向的處理器邊界，它們將整個(gè)計(jì)算域劃分為4個(gè)子域，每個(gè)子域用一個(gè)處理器計(jì)算。當(dāng)粒子整體向右移動(dòng)時(shí)，每個(gè)處理器邊界也適當(dāng)向右移動(dòng)，來(lái)盡量保證每個(gè)子域的粒子數(shù)相等。

動(dòng)態(tài)空間分解雖然能夠解決負(fù)載不均衡問(wèn)題，但是實(shí)現(xiàn)起來(lái)復(fù)雜，所以本文不考慮動(dòng)態(tài)空間分解。

2.2 粒子分解方案

本文針對(duì)SPH方法的主要計(jì)算量在于近鄰粒子搜索這一特點(diǎn)，提出并實(shí)現(xiàn)了一種新的粒子分解方案，其基本思想是將找出所有近鄰粒子對(duì)這樣一個(gè)大任務(wù)劃分為np(處理器數(shù))個(gè)小任務(wù)，每個(gè)處理器完成一個(gè)。

首先將所有粒子依據(jù)編號(hào)從小到大分為np個(gè)組Gm(m=1,2,…,np)，每個(gè)組的粒子數(shù)近似相等，約為N/np(其中N為粒子總數(shù))，處理器m維護(hù)組Gm。當(dāng)粒子i和粒子j是近鄰粒子時(shí)，可以用(i,j)來(lái)表示這一近鄰粒子對(duì)，為了避免重復(fù)，因此要求i

由于m≤n，可知子集(Gm,Gn)有np(np+1)/2個(gè)，且相互不重疊。每個(gè)處理器負(fù)責(zé)(np+1)/2個(gè)。顯然，當(dāng)np為奇數(shù)時(shí)，(np+1)/2為整數(shù)，此時(shí)，每個(gè)處理器的負(fù)載非常均衡。下面以np=5為例來(lái)說(shuō)明這些子集在處理器上的分配策略，如表1所示。

表1 近鄰粒子對(duì)在不同處理器之間分配策略

由表1可以看到，處理器P1負(fù)責(zé)維護(hù)G1這組粒子。在每個(gè)時(shí)間步長(zhǎng)內(nèi)，P1首先需要將G1組粒子的所有信息發(fā)送給P4和P5，因?yàn)镻4需要找出近鄰粒子對(duì)子集(G1,G4)，P5需要找出(G1,G5)；然后P1搜索被分配的近鄰粒子對(duì)(本文采用樹(shù)搜索算法)并同時(shí)計(jì)算相互作用力；最后P4和P5將所計(jì)算的G1組粒子受到的部分作用力發(fā)送回P1，由P1更新G1這組粒子的狀態(tài)。

處理器P2的計(jì)算過(guò)程類(lèi)似，即P2首先將G2組粒子的所有信息發(fā)送給P1和P5；然后尋找被分配的近鄰粒子對(duì)并同時(shí)計(jì)算相互作用力；最后P1和P5將所計(jì)算的G2組粒子受到的部分作用力發(fā)送回P2，由P2更新G2這組粒子的狀態(tài)。其他處理器的計(jì)算過(guò)程也類(lèi)似。

于是，每個(gè)時(shí)間步可以分為如下4步：粒子信息發(fā)送、近鄰粒子搜索并同時(shí)計(jì)算粒子間相互作用力、相互作用力信息交換和更新粒子狀態(tài)。分別統(tǒng)計(jì)每步的計(jì)算時(shí)間，可以得到：第1步和第3步為信息傳遞，時(shí)間隨著處理器的增加而緩慢增加；第2步是主要計(jì)算量，隨著處理器的增加而單調(diào)減??；第4步也是計(jì)算，時(shí)間幾乎可以忽略不計(jì)[10]。

可以看到，本粒子分解方案有如下特點(diǎn)：①一個(gè)時(shí)間步長(zhǎng)內(nèi)的所有通信為每個(gè)處理器調(diào)用一次廣播函數(shù)和一次求和歸約函數(shù)，實(shí)現(xiàn)非常簡(jiǎn)單；②每個(gè)處理器的計(jì)算量和通信量都大致相等，負(fù)載非常均衡；③每步的信息傳遞量為O(N)，與粒子的分布無(wú)關(guān)；④可以十分方便地采取任意奇數(shù)個(gè)處理器并行計(jì)算。

3 兩種并行方案的比較

3.1 二維圓柱繞流

模型參數(shù)與文獻(xiàn)[8]中的圓柱繞流算例一致，雷諾數(shù)為0.03。初始時(shí)，X方向和Y方向各均勻分布1000個(gè)粒子，總粒子數(shù)為100萬(wàn)。X方向?yàn)橹芷谛赃吔鐥l件，對(duì)應(yīng)的實(shí)際物理問(wèn)題為無(wú)限長(zhǎng)管道中均勻分布無(wú)數(shù)個(gè)圓柱。圖3給出了SPH模擬結(jié)果和有限元模擬結(jié)果的對(duì)比，其中圖3(a)為文獻(xiàn)[8]中的有限元模擬結(jié)果，圖3(b)為本文模擬結(jié)果。

(a)有限單元法(a)Finte Elmemt Modeling

(b) 光滑粒子動(dòng)力學(xué)(b)Smoothed particle hydrodynamics圖3 SPH和有限元模擬的速度大小等值線(xiàn)對(duì)比Fig.3 Contour lines of velocity magnitude using FEM and SPH for Re=0.03

分別采取兩種并行方案，測(cè)得的加速比如圖4所示。

圖4 空間分解和粒子分解加速比隨處理器的變化Fig. 4 Speed-up ratio of domain decomposition and particle decomposition with the change of processors

從圖4可以看到，空間分解時(shí)，當(dāng)處理器數(shù)目在小于100時(shí)，加速比反常地略高于理想值，其原因可能和緩存等硬件有關(guān)，文獻(xiàn)[4-5]以及開(kāi)源軟件parallelSPHYSICSv2.0使用手冊(cè)中都報(bào)告了類(lèi)似的情況。當(dāng)處理器數(shù)目為320時(shí)，加速比約為234.69，并行效率為73.34%。如果繼續(xù)增加處理器數(shù)至640，加速比迅速降低為100.13，對(duì)應(yīng)的并行效率為15.65%，因此結(jié)果沒(méi)有在圖4中給出。

粒子分解時(shí)，最多采用了93個(gè)處理器，此時(shí)的加速比為23.33，并行效率為25.09%。如果繼續(xù)增加處理器數(shù)量，加速比則緩慢下降。

因此，對(duì)于粒子分布均勻的情況，空間分解的加速比遠(yuǎn)大于粒子分解。而且可以發(fā)現(xiàn)，粒子分解不適合上百個(gè)處理器的大規(guī)模并行計(jì)算。

3.2 三維潰壩問(wèn)題

模型參數(shù)與文獻(xiàn)[11]中的潰壩算例一致，初始水柱高0.292m，寬0.146m，容器寬度為0.584m。初始粒子總數(shù)大約為40萬(wàn)，時(shí)間步長(zhǎng)固定為20μs，共計(jì)算了50 000個(gè)時(shí)間步長(zhǎng)。圖5為文獻(xiàn)[11]給出的二維流體體積法(VolumeOfFluid，VOF)模擬結(jié)果，從左到右的時(shí)刻依次為0.2s，0.4s，0.6s和0.8s。

圖5 潰壩問(wèn)題的二維模擬結(jié)果[10]Fig. 5 2D VOF simulation of dam break

本文采用靜態(tài)空間分解和粒子分解兩種不同的并行方案對(duì)同樣的問(wèn)題進(jìn)行了三維數(shù)值模擬。模擬的粒子分布如圖6所示。

(a) t=0.2s

(b) t=0.4s

(c) t=0.6s

(d) t=0.8s圖6 三維潰壩問(wèn)題的SPH模擬Fig. 6 3D SPH simulation of dam break

可以看到，三維模擬結(jié)果和文獻(xiàn)[11]中的實(shí)驗(yàn)及VOF模擬結(jié)果基本一致。為了和二維模擬結(jié)果對(duì)比，圖7給出了0.8s時(shí)刻粒子分布的一個(gè)二維視圖。

圖7 0.8s時(shí)刻的二維視圖Fig.7 2D view point at t=0.8s

從圖7可以看到，流體在0.8s時(shí)刻形成了一個(gè)空腔，這和VOF計(jì)算結(jié)果吻合。

因?yàn)槿萜髟赬，Y，Z方向上的尺寸約為1∶4∶2，因此在空間分解時(shí)，X，Y，Z方向上的處理器個(gè)數(shù)取為n，4n，2n，從而讓每個(gè)處理器分配的空間大小近似相等，總處理器個(gè)數(shù)為8n3。但是由于粒子在整個(gè)容器中的分布極不均勻，每個(gè)處理器上分配的粒子數(shù)實(shí)際上極不均勻，甚至多數(shù)處理器上沒(méi)有粒子。本算例中，空間分解時(shí)處理器個(gè)數(shù)分別取為8，64，216，分別對(duì)應(yīng)于n=1，2，3。

兩種并行方案測(cè)得的計(jì)算時(shí)間如圖8所示。

圖8 空間分解和粒子分解計(jì)算時(shí)間隨處理器的變化Fig.8 Total computing time of domain decomposition and particle decomposition with the change of processors

可以看到，空間分解的計(jì)算時(shí)間比粒子分解要長(zhǎng)，而且計(jì)算時(shí)間隨處理器數(shù)的增加緩慢下降。當(dāng)處理器數(shù)為216時(shí)，空間分解需要約7.7h，加速比為9.0。此處加速比定義為speedup=T1/Tn，其中Tn為采用n個(gè)處理器并行計(jì)算的時(shí)間，T1為串行計(jì)算時(shí)間。由于空間分解程序的計(jì)算時(shí)間較長(zhǎng)，沒(méi)有測(cè)試其串行計(jì)算時(shí)間，所以其串行計(jì)算時(shí)間采用粒子分解程序的串行時(shí)間。計(jì)算中還發(fā)現(xiàn)，當(dāng)處理器數(shù)為216時(shí)，絕大多數(shù)處理器中沒(méi)有粒子，這是由潰壩問(wèn)題中粒子分布不均勻的特點(diǎn)以及沒(méi)有動(dòng)態(tài)負(fù)載均衡造成的。

粒子分解時(shí)，當(dāng)處理器數(shù)為23時(shí)，計(jì)算時(shí)間為5.7h，幾乎已經(jīng)降到最低值，且計(jì)算時(shí)間比空間分解方案的最優(yōu)值少。當(dāng)處理器數(shù)為95時(shí)，需要的時(shí)間最少，為4.3h，最大加速比為16.1。如果繼續(xù)增加處理器數(shù)，計(jì)算時(shí)間緩慢增加。

4 結(jié)論

本文提出了一種非常簡(jiǎn)單的基于粒子分解的SPH并行方案，并以二維圓柱繞流和三維潰壩問(wèn)題的SPH并行模擬為例，比較了粒子分解和空間分解的計(jì)算時(shí)間或加速比，得到了以下主要結(jié)論：

1)當(dāng)粒子分布均勻時(shí)，SPH并行應(yīng)該采用空間分解方案。此時(shí)，只有處理器邊界處有信息傳遞，信息傳遞量小，能夠達(dá)到非常大的加速比。而粒子分解的加速比相對(duì)較小。

2)當(dāng)粒子分布非常不均勻時(shí)，SPH并行可以采用粒子分解方案。此時(shí)，靜態(tài)空間分解的負(fù)載嚴(yán)重不平衡，導(dǎo)致加速比較小。而粒子分解采用較少處理器就能達(dá)到可觀加速比。

3)從文中兩個(gè)算例可以看出，當(dāng)問(wèn)題規(guī)模在40萬(wàn)～100萬(wàn)左右時(shí)，粒子分解方案的加速比固定在16～23左右。在粒子分布非常不均勻的情況下，這一結(jié)果可能比靜態(tài)空間分解更優(yōu)。因而，本文提出的粒子分解方案在自由表面流動(dòng)模擬中具有一定優(yōu)勢(shì)。

References)

[1]LucyLB.Anumericalapproachtothetestingofthefissionhypothesis[J].AstronomicalJournal, 1977, 82: 1013-1024.[2]MonaghanJJ,GingoldRA.ShocksimulationbytheparticlemethodSPH[J].JournalofComputationalPhysics, 1983, 52(2): 374-389.

[3]PlimtonS.Fastparallelalgorithmsforshort-rangemoleculardynamics[J].JournalofComputationalPhysics, 1995, 117(1): 1-19.

[4]FerrariA,DumbserM,ToroEF,etal.Anew3DparallelSPHschemeforfreesurfaceflows[J].Computers&Fluids, 2009, 38(6): 1203-1217.

[5]MarroneS,BouscasseB,ColagrossiA,etal.Studyofshipwavebreakingpatternsusing3DparallelSPHsimulations[J].Computers&Fluids, 2012, 69: 54-66.

[6]CherfilsJM,PinonG,RivoalenE.JOSEPHINE:aparallelSPHcodeforfree-surfaceflows[J].ComputerPhysicsCommunication, 2012, 183(7): 1468-1480.

[7]MoulinecC,IssaR,LatinoD,etal.High-performancecomputingandsmoothedparticlehydrodynamics[C]//ProceedingsofParallelComputationalFluidDynamics,Lyon, 2008.

[8]MorrisJP,FoxPJ,ZhuY.ModelinglowreynoldsnumberincompressibleflowsusingSPH[J].JournalofComputationalPhysics, 1997, 136(1): 214-226.

[9]MonahanJJ.SimulatingfreesurfaceflowswithSPH[J].JournalofComputationalPhysics, 1994, 110(2): 399-406.[10] 周浩，湯文輝，冉憲文，等. 光滑粒子流體動(dòng)力學(xué)的一種并行數(shù)值計(jì)算方案[J]. 航天器環(huán)境工程， 2012， 29(1): 23-26.ZHOUHao,TANGWenhui,RANXianwen,etal.Spacecraftenvironmentengineering[J].SpacecraftEnvironmentEngineering, 2012, 29(1): 23-26.(inChinese)[11]DoringM,AndrillonA,AlessandriniB,etal.SPHfreesurfaceflowsimulation[C]//Proceedingsofthe18thInternationalWorkshoponWaterWavesandFloatingBodies,LeCroisic, 2003.

SPH parallel schemes and its application in free surface flow simulation

ZHOU Hao1, XU Zhihong1, ZOU Shun2, TANG Wenhui1

Forthemaincomputationofsmoothedparticlehydrodynamicswasfinitenearestparticlesearch,anovelschemetoparallelizethesmoothedparticlehydrodynamicsmethodbasedontheconceptofparticledecompositionwasproposed.Anyserialsmoothedparticlehydrodynamicscodecouldbeeasilyparallelizedbyusingtheproposedscheme.Theamountofinformation,whichwastransformedineachtimestep,dependedonlyonthetotalnumberofparticles,butnotonthespatialdistributionofparticles.Therefore,theproposedschemewasparticularlyusefulfortheparallelsimulationofcasesinvolvedviolentfreesurfacemovements.Fromthesimulationresultsofa3Ddambreakcasewithatotalnumberof0.4million,theproposedschemecanachieveaspeedupratioabout16,whichprovesthattheproposedschememaybeisbetterthanthedomaindecompositionscheme(withoutconsideringdynamicloadbalance).

smoothedparticlehydrodynamics;parallelscheme;particledecomposition;domaindecomposition

2014-07-08

國(guó)家自然科學(xué)基金資助項(xiàng)目(61221491，61303071)，廣州科信局基金資助項(xiàng)目(134200026)

周浩(1986—)，男，湖南常德人，博士研究生， E-mail: zhouhao12029014@163.com；湯文輝(通信作者)，男，教授，博士， 博士生導(dǎo)師，E-mail：wenhuitang@163.com

10.11887/j.cn.201504012

http://journal.nudt.edu.cn

TB

A