董蘭芳， 章 恒， 陳駿雄

(1. 中國(guó)科學(xué)技術(shù)大學(xué)計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院，安徽 合肥 230027；2. 浙江大學(xué)流體動(dòng)力與機(jī)電系統(tǒng)國(guó)家重點(diǎn)實(shí)驗(yàn)室，浙江 杭州 310027)

基于IIF模型的流體表面張力模擬

董蘭芳1,2， 章 恒1， 陳駿雄1

(1. 中國(guó)科學(xué)技術(shù)大學(xué)計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院，安徽 合肥 230027；2. 浙江大學(xué)流體動(dòng)力與機(jī)電系統(tǒng)國(guó)家重點(diǎn)實(shí)驗(yàn)室，浙江 杭州 310027)

基于光滑粒子動(dòng)力學(xué)(SPH)方法模擬流體時(shí)流體表面張力的作用在固液、氣液交界處不可忽視，其影響模擬的準(zhǔn)確性和視覺真實(shí)感。目前已有的表面張力模型如連續(xù)表面力(CSF)模型、粒子間相互作用力(IIF)模型都存在各自的缺陷。針對(duì) IIF模型模擬表面張力時(shí)所產(chǎn)生的粒子非物理聚集、流體表面形狀不規(guī)則等現(xiàn)象，采用基于類Lennard-Jones勢(shì)函數(shù)的分子間聚斥力對(duì)表面張力建模，并定義了基于法向差的SPH張力修正項(xiàng)以解決IIF模型不能保證流體表面面積最小化問題。實(shí)驗(yàn)結(jié)果表明，該方法能夠穩(wěn)定和正確地模擬兩相交界處的表面張力的效果。

計(jì)算機(jī)圖形學(xué)；流體模擬；光滑粒子動(dòng)力學(xué)方法；表面張力；粒子間相互作用力模型

表面張力作為流體的一種常見和重要的物理性質(zhì)，隨著對(duì)流體細(xì)節(jié)和模擬精度要求的增強(qiáng)，在諸如計(jì)算機(jī)圖形學(xué)領(lǐng)域和工程數(shù)值計(jì)算領(lǐng)域?qū)Ρ砻鎻埩Φ哪M也越來越重視。從微觀層面說，表面張力的產(chǎn)生是由于表面區(qū)域的水分子受到各種方向分子間作用力的合力不為零，指向流體內(nèi)部。從宏觀表現(xiàn)上看，根據(jù)拉普拉斯定律，表面張力的作用就是最小化表面區(qū)域面積，如自然界中液滴的形成、聚合和分裂。

光滑粒子動(dòng)力學(xué)方法(smoothed particle hydrodynamics, SPH)作為一種典型的拉格朗日流體模擬方法，由于具備相較于傳統(tǒng)基于網(wǎng)格方法的易編程、方便處理復(fù)雜邊界和大形變問題等優(yōu)點(diǎn)而被廣泛應(yīng)用于流體模擬領(lǐng)域[1-3]。目前基于SPH的表面張力模擬方法主要分為兩種：基于連續(xù)表面力模型(continuum surface force, CSF)和基于粒子間相互作用力模型(inter-particle interaction force, IIF)。CSF模型是將作用在流體表面上的力轉(zhuǎn)換為周圍流體體積力，最早由Morris[4]在2000年提出，并與基于網(wǎng)格的流體體積函數(shù) (volume of fluid, VOF)方法進(jìn)行了對(duì)比，得到可信的結(jié)果。在Morris的基礎(chǔ)上，Müller等[5]采用顏色域函數(shù)對(duì)其進(jìn)行簡(jiǎn)化，Hu和Adams[6]則通過配置新的粘性項(xiàng)，解決了速度和剪切應(yīng)力在多相交界處的連續(xù)性問題。2011年強(qiáng)洪夫等[7]采用SPH方法中解決邊界核函數(shù)插值問題的校正平滑粒子方法(corrective smoothed particle method, CSPM)對(duì)CSF模型中表面張力計(jì)算進(jìn)行修正，并對(duì)二維溶液中油滴的碰撞破碎過程進(jìn)行了模擬。如文獻(xiàn)[7]所述，CSF模型模擬表面張力時(shí)在尖角、邊界粒子缺失嚴(yán)重的部位曲率計(jì)算誤差很大，需要花費(fèi)額外的計(jì)算代價(jià)來進(jìn)行修正。另外，CSF模型以一種非對(duì)稱的形式將力施加到SPH粒子上，這使得流體的動(dòng)量不再守恒。

IIF模型是一種基于分子間相互作用力，從微觀尺度上模擬表面張力產(chǎn)生的模型。在流體內(nèi)部，分子受到的周圍領(lǐng)域內(nèi)分子在各種方向上作用力的合力為零，與此同時(shí)，處于流體表面的分子僅受到流體內(nèi)部分子的作用力而導(dǎo)致合力非零，宏觀表現(xiàn)為表面張力。SPH方法模擬表面張力時(shí)，將流體用粒子表示，應(yīng)用IIF模型思想則對(duì)每對(duì)粒子施加一對(duì)對(duì)稱的作用力[8-9]，在粒子分布近似均勻的前提下，對(duì)整個(gè)流體系統(tǒng)施加這一作用力可近似模擬表面張力效果。文獻(xiàn)[8]和[9]的不同之處在于采用不同的分子間作用力模型。IIF模型克服了CSF模型的缺點(diǎn)，不再需要計(jì)算曲率，因此編程簡(jiǎn)單，計(jì)算效率高。同時(shí)由于對(duì)稱形式作用力的引入保證流體的動(dòng)量守恒。然而實(shí)驗(yàn)證明，僅僅引入粒子間相互作用力不能保證表面區(qū)域面積趨向最小，因?yàn)榱黧w平衡狀態(tài)可以在非最小表面面積狀態(tài)下達(dá)到。

本文所基于IIF模型的表面張力模型，采用基于距離的類 Lennard-Jones勢(shì)函數(shù)為粒子間作用力建模。同時(shí)，為了解決IIF模型不保證表面面積最小問題，利用文獻(xiàn)[5]中的一些結(jié)論引入了一個(gè)張力修正項(xiàng)，通過對(duì)立方體液滴在去除重力下形變過程的模擬，得到了與文獻(xiàn)[7]中修正后的CSF模型幾乎一致的結(jié)果，證明了該方法的正確性。通過對(duì)碰撞破碎過程對(duì)比和多液滴合并的模擬證明該方法的有效性。

1 流體模擬框架

1.1 SPH基本形式

SPH方法中用粒子來承載諸如密度、質(zhì)量等各個(gè)物理量，連續(xù)的流體場(chǎng)被離散成大量獨(dú)立的粒子，流體連續(xù)的物理量根據(jù)周圍流體粒子相對(duì)應(yīng)的物理量插值得到。SPH方法中的核函數(shù)用于確定空間中周圍粒子對(duì)中心粒子插值的權(quán)重，一般來說，核函數(shù)的選擇使得離中心粒子越近的粒子插值權(quán)重越大。SPH基本離散化的形式[10]為：

其中，W為插值所用的光滑核函數(shù)；h為其核半徑；ri和rj分別為粒子i、j的位置；Aj為粒子j所對(duì)應(yīng)的物理量；Vj為粒子j的體積。

以當(dāng)前粒子為中心，h為半徑的球體范圍稱為當(dāng)前粒子的支持域。直觀上，SPH方法通過枚舉支持域內(nèi)所有粒子以用于插值，從而得到當(dāng)前粒子的相應(yīng)屬性值。

1.2 SPH流體運(yùn)動(dòng)方程

流體運(yùn)動(dòng)用經(jīng)典的Navier-Stokes方程來刻畫，如式(2)所示：

其中，ρ為流體的密度；u為速度矢量；p為壓強(qiáng)；μ為粘力系數(shù)；f為合外力項(xiàng)諸如重力和浮力等。

基于文獻(xiàn)[6]用SPH方法離散化Navier-Stokes方程，先將密度項(xiàng)根據(jù)粒子的分布計(jì)算出來，再用于計(jì)算壓力項(xiàng)和粘力加速度項(xiàng)，如式(3)～(5)所示。

為了從密度求出壓強(qiáng)，采用文獻(xiàn)[9]提出的弱可壓縮法(weakly compressible SPH，WCSPH)，模擬時(shí)間步長(zhǎng)根據(jù)CFL條件[10]來估計(jì)。對(duì)于在流體與空氣交界處表面由于缺少鄰居粒子，導(dǎo)致密度估計(jì)失真而產(chǎn)生的拉伸不穩(wěn)定現(xiàn)象，采用Shepard filter來對(duì)流體密度進(jìn)行平滑。對(duì)于在流體與固體交界處為了避免 sticking particles問題，采用文獻(xiàn)[11]中利用單層固體邊界粒子對(duì)周圍鄰居流體粒子進(jìn)行密度校正的方法解決。

2 修正的IIF表面張力模型

流體的表面張力由分子間的內(nèi)聚力產(chǎn)生，然而，在采用單純的IIF模型模擬SPH流體是在宏觀層面上考慮粒子間的相互作用力，由于對(duì)稱的粒子作用力的引入，整個(gè)流體系統(tǒng)的動(dòng)量是守恒的，流體系統(tǒng)的平衡態(tài)可以在流體的任意形態(tài)達(dá)到，這樣導(dǎo)致的結(jié)果就是流體表面區(qū)域面積的最小化難以保證，需要通過其他手段來進(jìn)行修正。本文中表面張力模型中采用一個(gè)類Lennard-Jones勢(shì)函數(shù)為粒子間的聚斥力建模，同時(shí)考慮使得表面區(qū)域最小化的表面力項(xiàng)。

2.1 粒子間聚斥力模型

IIF模型中為粒子間施加一個(gè)基于粒子間距離的作用力，基本形式如式(6)所示，F(xiàn)cohesioni＜-j描述粒子j對(duì)粒子i的作用力。其中，mi和mj是粒子i和j的質(zhì)量，L(d)是一個(gè)核函數(shù)控制粒子間作用力的正負(fù)(分別對(duì)應(yīng)著排斥力和吸引力)，文獻(xiàn)[8]最早提出基于一個(gè)類余弦函數(shù)的粒子間聚斥力，但該函數(shù)在粒子距離 d趨于0時(shí)，函數(shù)值也趨于0，意味著粒子間的排斥力隨著距離的減小而減小，這會(huì)導(dǎo)致在模擬弱可壓縮流體時(shí)出現(xiàn)非物理聚集現(xiàn)象。文獻(xiàn)[9]采用一個(gè)SPH核函數(shù)來為粒子間吸引力建模，由于缺少排斥力項(xiàng)，該模型也會(huì)導(dǎo)致非物理聚集現(xiàn)象的發(fā)生。本文采用文獻(xiàn)[12]中模擬流體和可變形固體交互時(shí)流固間作用力所使用的類Lennard-Jones勢(shì)函數(shù)作用核函數(shù)：

其中，h為作用范圍；k為控制表面張力影響程度的系數(shù)；d0為靜止距離；d為當(dāng)前兩粒子間距離。

式(7)中不同作用范圍下L(d)的值隨粒子間距離d的變化而變化(圖1)，其中系數(shù)k值取0.001，d0值取4×10-3mm，h分別取1.5 d0, 2 d0, 2.5 d0。

圖1 不同作用范圍下L(d)隨粒子間距的變化

由圖1可見，當(dāng)d=d0時(shí)，粒子間作用力為0，當(dāng) h＞d＞d0時(shí),函數(shù)值為負(fù)，表現(xiàn)出吸引力，當(dāng)0＜d＜d0時(shí)，函數(shù)值為正，變現(xiàn)為排斥力。另外,當(dāng)d趨于 0時(shí)，函數(shù)值趨于 k，這樣避免了文獻(xiàn)[8] 和[9]中的非物理聚集現(xiàn)象發(fā)生。

2.2 表面張力修正項(xiàng)

基于類Lennard-Jones勢(shì)函數(shù)的粒子間聚斥力可以有效地解決粒子非物理聚集的問題，并且保證動(dòng)量守恒，但不能保證表面面積最小化。為了解決這個(gè)問題，本文引入一個(gè)表面張力的修正項(xiàng)：

其中，mi為當(dāng)前粒子i的質(zhì)量；ni和nj分別為粒子i處的法向和粒子j處的法向；λ為修正項(xiàng)系數(shù)。

該修正項(xiàng)的思想受文獻(xiàn)[13]中的分析啟發(fā)，從離散的觀點(diǎn)來說，表面張力趨向于將一個(gè)局部區(qū)域內(nèi)所有離散采樣點(diǎn)的法向變得一致。直觀地說，即撫平該區(qū)域的曲率。對(duì)于粒子 i，的意義在于作用于粒子i，使得i和j之間的法向差變小，推廣到整個(gè)計(jì)算域中所有流體粒子，項(xiàng)的施加將使得曲率大的地方(即流體表面)的法向差變小，即減小曲率，從而最小化表面面積。推廣到整個(gè)計(jì)算域上，粒子i的表面張力修正項(xiàng)用SPH形式表示有以公式：

其中，F(xiàn)fixi為粒子i的表面張力修正項(xiàng)。

對(duì)于一對(duì)粒子i和j而言：

粒子間的表面張力修正項(xiàng)也是對(duì)稱的，因此修正項(xiàng)的引入不破壞動(dòng)量守恒條件。關(guān)于粒子法向的計(jì)算，采用文獻(xiàn)[5]的方法，即先為整個(gè)計(jì)算域定義一個(gè)顏色域函數(shù)，然后計(jì)算顏色函數(shù)的梯度，得到的梯度函數(shù)就是計(jì)算域上每處的未歸一化法向。用SPH的形式表示如下：

在CSF模型中，通過每個(gè)粒子處法向的散度來估計(jì)該處的曲率，然后定義一個(gè)大小與曲率相關(guān)，方向?yàn)榉ㄏ虻牧?，即該粒子處所受到的表面張力。整個(gè)過程需要對(duì)法向進(jìn)行歸一化，而在流體內(nèi)部，對(duì)法向歸一化可能導(dǎo)致流體內(nèi)部的法向計(jì)算不準(zhǔn)確，從而需要額外處理來修正。在本文模型中，由于避免了直接計(jì)算粒子處的曲率，從而也無需對(duì)法向進(jìn)行歸一化，因此也就不需要進(jìn)行額外的處理。

3 算法實(shí)現(xiàn)和參數(shù)選擇

本文模型中計(jì)算每個(gè)粒子的表面張力時(shí)需要考慮周圍鄰居粒子的貢獻(xiàn)，所以對(duì)表面張力的模擬是耦合在流體模擬框架中的。單獨(dú)來看，具體步驟如圖2所示。

圖2 表面張力模擬流程

關(guān)于 SPH流體的模擬，采用文獻(xiàn)[6]和[14]中所設(shè)計(jì)的幾個(gè)多項(xiàng)式核函數(shù)，針對(duì)計(jì)算不同物理量選取不同的核函數(shù)。SPH的光滑核半徑取值為1.5R，R是 SPH粒子的直徑，粒子半徑的選取與采樣流體所使用的粒子數(shù)有關(guān)，壓強(qiáng)的計(jì)算采用文獻(xiàn)[9]中的WCSPH方法，時(shí)間步長(zhǎng)的選擇根據(jù)CFL條件[10]來估計(jì)，在實(shí)驗(yàn)中時(shí)間步長(zhǎng)設(shè)為0.005 s。數(shù)值積分采用簡(jiǎn)單的leap-frog模式[15]。關(guān)于表面張力模型中涉及的2個(gè)主要參數(shù)：類Lennard-Jones勢(shì)函數(shù)中的聚斥力縮放系數(shù) k和表面張力修正項(xiàng)中修正系數(shù)λ，這兩個(gè)模型參數(shù)都是經(jīng)驗(yàn)值。實(shí)驗(yàn)表明，兩個(gè)參數(shù)的取值應(yīng)在同一數(shù)量級(jí)上，在實(shí)驗(yàn)中將其保持與流體的粘度系數(shù)一致(0.001)，能夠在維持流體粒子的1%密度震蕩幅度的同時(shí)取得良好的可視化效果。

4 實(shí)驗(yàn)算例和結(jié)果

實(shí)驗(yàn)配備 2.8 GHz Intel Core i5-2300 四核CPU和GeForce GTS 450顯卡以及4 G內(nèi)存的硬件機(jī)器，流體粒子形式的渲染基于OpenGL平臺(tái)實(shí)時(shí)繪制。

4.1 立方體液體形變?yōu)橐?guī)則球體

算例1. 模擬了0.5×0.5×0.5(m3)的水柱在無重力的情況下由表面張力驅(qū)動(dòng)，形變成規(guī)則球體的過程，其中流體粒子4 913個(gè)。在無重力等外力的情況下，流體穩(wěn)定狀態(tài)下的形狀主要由表面張力決定，表面張力使得流體各個(gè)地方的曲率趨于一致、表面面積趨于最小化，立方體液塊在表面張力的作用下將收縮震蕩，最終由于流體的粘性作用而穩(wěn)定成一個(gè)規(guī)則的球體。該算例的目的在于證明該表面張力算法的正確性。圖 3是液塊形變模擬過程的部分關(guān)鍵幀截圖信息(流體為三維表示，觀察方向?yàn)樗椒较?。

圖3 立方體液塊形變過程

首先，液塊由立方體先是收縮，再迅速形變?yōu)橐粋€(gè)球體(圖 3(a)～(d))；然后，液塊由一個(gè)球體震蕩為一個(gè)具有6個(gè)尖角的極限狀態(tài)(圖3(e)～(h))；最后，流體由于粘力對(duì)動(dòng)能的耗散作用使得液塊形最終穩(wěn)定在球體狀態(tài)下(圖3(i)～(l))。該結(jié)果與文獻(xiàn)[7]中采用修正的CSF模型在二維情況下得到的數(shù)值模擬結(jié)果類似。

4.2 液塊碰撞破碎和液滴合并

算例2. 模擬了在體積為1 m×1 m×1 m的固體容器中，體積為0.25 m×0.25 m×0.25 m的液塊在重力的作用下下落到與容器底碰撞破碎形成多個(gè)液滴，以及液滴在表面張力的作用下合并的過程，并與不施加表面張力情況下做出對(duì)比。其中粒子數(shù)為810個(gè)。圖4給出液塊在碰撞破碎后施加表面張力與不施加表面張力模擬過程對(duì)比(流體為三維表示，觀察方向?yàn)閺娜萜黜敳康降撞康拇怪狈较?。

圖4 碰撞破碎過程比較

圖 4(a)和(b)自左向右分別是施加表面張力和不施加表面張力的液塊自下落開始第60、70、80、90幀的截圖信息。從2個(gè)模擬過程的對(duì)比中可以看出，在施加表面張力的情況下液塊在與容器底部碰撞后破碎形成液滴，如圖3(a)所示，而不施加表面張力時(shí)破碎后成網(wǎng)狀失真，并且不能形成分離的圓形液滴，如圖3(e)所示。另外，隨著流體的運(yùn)動(dòng)，施加表面的張力的模擬過程中還伴有液滴間的碰撞合并現(xiàn)象。

液體運(yùn)動(dòng)趨于穩(wěn)定后的比較如圖5所示(觀察方向?yàn)樗椒较?，圖 5(a)為施加表面張力后形成的一個(gè)單獨(dú)液滴，包含3層SPH粒子，圖5(b)則為沒有施加表面張力的情況下，流體粒子無法聚攏，從而平鋪在容器底部，只有1層SPH粒子。

圖5 流體靜止?fàn)顟B(tài)比較

圖6為8幅流體粒子位置截圖，顯示碰撞破碎產(chǎn)生的多個(gè)獨(dú)立的運(yùn)動(dòng)液滴在表面張力的作用下逐漸合并為一個(gè)單獨(dú)的穩(wěn)定的圓形液滴的模擬過程，該過程發(fā)生液體在與容器壁發(fā)生碰撞后，多個(gè)分散的液滴向容器中央聚攏。

圖6 多液滴合并過程

4.3 計(jì)算性能對(duì)比

為了測(cè)試本文方法的計(jì)算性能，使用方形液滴自然變化和油滴互溶兩個(gè)算例，選取運(yùn)行40幀后添加表面張力的效果與傳統(tǒng)的 IIF模型和 CSF模型進(jìn)行對(duì)比，時(shí)間為每幀的計(jì)算時(shí)間(表1)。

表1 計(jì)算速度對(duì)比(s)

本文方法由于在計(jì)算表面張力修正項(xiàng)過程中直接使用顏色域梯度求解粒子法向，相較于傳統(tǒng)CSF模型省去了法向歸一化以及相應(yīng)修正處理過程，在單幀計(jì)算速度上略有提高，模擬的場(chǎng)景越復(fù)雜，效果越明顯。與傳統(tǒng)IIF模型相比，本文方法在保證模擬速度的基礎(chǔ)上，通過改進(jìn)聚斥力模型和修正表面張力進(jìn)一步提高了模擬精度。

5 結(jié) 論

本文針對(duì)傳統(tǒng)的表面張力模擬所采用的IIF模型中存在的缺陷，定義一種最小化表面面積的表面張力修正項(xiàng)，結(jié)合基于類L-J勢(shì)函數(shù)的粒子間聚斥力，為SPH方法模擬流體時(shí)表面張力的模擬建模。通過模擬立方體液塊僅由表面張力驅(qū)動(dòng)時(shí)形變?yōu)榍蝮w的算例驗(yàn)證了本方法的正確性。最后通過比較液塊在施加表面張力和不施加情況下碰撞破碎過程的對(duì)比以及多液滴合并算例驗(yàn)證了本文方法的有效性。

[1] Solenthaler B, Pajarola R. Predictive-corrective incompressible SPH [J]. Acm Transactions on Graphics, 2009, 28(3): 341-352.

[2] 溫嬋娟, 歐嘉蔚, 賈金原. GPU通用計(jì)算平臺(tái)上的SPH 流體模擬[J]. 計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào), 2010, 22(3): 406-411.

[3] Zhu Y N, Bridson R. Animating sand as a fluid [J]. Acm Transactions on Graphics, 2006, 24(3): 965-972.

[4] Morris J P. Simulating surface tension with smoothed particle hydrodynam ics [J]. International Journal for Numerical Methods in Fluids, 2000, 33: 333-353.

[5] Müller M, Charypar D, Gross M. Particle-based fluid simulation for interactive applications [C]//Preceedings of Symposium on Computer Animation. Aire?La?Ville: Eurographics Association Press, 2003: 154-159.

[6] Hu X Y, Adams N A. A multi-phase SPH method for macroscopic and mesoscopic flow s [J]. Journal of Computational Physics, 2006, 213(2): 844-861.

[7] 強(qiáng)洪夫, 陳福振, 高巍然. 基于 SPH方法的表面張力修正算法及其應(yīng)用[J]. 計(jì)算力學(xué)學(xué)報(bào), 2011, 28(3): 37-42.

[8] Tartakovsky A, Meakin P. Modeling of surface tension and contact angles with smoothed particle hydrodynam ics [J]. Physical Review E, 2005, 72(2): 026301.

[9] Becker M, Teschner M. Weakly compressible SPH for free surface flows [C]//Proceedings of the 2007 ACM SIGGRAPH/ Eurographics symposium on Computer Animation. Aire?La?Ville: Eurographics Association Press, 2007: 209-217.

[10] Monaghan J J. Smooth particle hydrodynamics [J]. Reports on Progress in Physics, 2005, 68(8): 1703-1759.

[11] Akinci N, Ihmsen M, Akinci G, et al. Versatile rigid-fluid coupling for incompressible SPH [J]. ACM Transactions on Graphics (TOG), 2012, 31(4): 62.

[12] Müller M, Schirm S, Teschner M, et al. Interaction of fluids with deformable solids [J]. Computer Animation and Virtual Worlds, 2004, 15(3/4): 159-171.

[13] He X W, Wang H M, Zhang F J, et al. Robust Simulation of Small-Scale Thin Features in SPH-based Free Surface Flows [EB/OL]. [2015-08-16]. http://life.kunzhou.net/ 2013/SPHsurfacetension.pdf, 2013.

[14] Desbrun M, Gascuel M P. Smoothed particles: a new paradigm for animating highly deformable bodies [M]. Vienna: Springer Vienna, 1996: 61-76.

[15] Kelager M. Lagrangian fluid dynam ics using smoothed particle hydrodynamics [D]. Denmark: University of Copenhagen, 2006.

Surface Tension Simulation for SPH Fluid Based on IIF Model

Dong Lanfang1,2, Zhang Heng1, Chen Junxiong1

(1. Department of Computer Science and Technology, University of Science and Technology of China, Hefei Anhui 230027, China; 2. The State Key Laboratory of Fluid Power Transm ission and Control, Zhejiang University, Hangzhou Zhejiang 310027, China)

Surface tension effect at fluid-air and fluid-solid interfaces should not be ignored when simulating fluid based on smoothed particle hydrodynam ics (SPH) method, it affects the accuracy of the simulation and visual realism. The existing surface tension models such as the continuum surface force (CSF) model and the inter-particle interaction force (IIF) model both have their own flaws. This paper focus on the IIF model which can causes non-physical particle clustering and the irregular shape of fluid surface, and surface tension is modeled by an attraction-repulsion force which is formalized by a Lennard-Jones like potential function, and a SPH tension correction term which based on the normal difference between neighboring particles is introduced to m inim ize the surface area and smooth the surface. Experimental results show that the proposed method can simulate the surface tension effect accurately and stably.

computer graphics; fluid simulation; smoothed particle hydrodynam ics method; surface tension; inter-particle interaction force method

TP 39

10.11996/JG.j.2095-302X.2016030302

A

2095-302X(2016)03-0302-06

2015-09-28；定稿日期：2015-11-13

浙江大學(xué)流體動(dòng)力與機(jī)電系統(tǒng)國(guó)家重點(diǎn)實(shí)驗(yàn)室開放基金項(xiàng)目(GZKF-201318)

董蘭芳(1970-)，女，安徽合肥人，副教授，碩士。主要研究方向?yàn)榭茖W(xué)計(jì)算可視化、計(jì)算機(jī)動(dòng)畫、視頻圖像智能分析。E-m ial：lfdong@ustc.edu.cn