１．２NＧS方程

流体的控制方程采用描述流体运动规律的NavierＧStokes(NＧS)方程，其形式为：∇·u＝０(２) ∂ u ∂ t＋u·∇ u＋１ ρ∇ P＝ g＋ ν∇２ u(３)其中，u为速度场；P为压力场；ν为流体的动力学黏度；ρ为密度；∇为空间梯度算符；g为外力；∇２为拉普拉斯运算符。(２)式、(３)式分别代表流体运动的质量与动量守恒定律。使用SPH方法求解NＧS方程的基本过程如下：ρi＝∑ N j＝１ mjW(xi－xj, h)(４) pi＝ k( ρi－ ρ０)(５) f p i＝－∑ N j＝１ mj pi＋ pj ２ ρj ∇ W(xi－xj, h)(６) f ν i＝ μ∑ N j＝１ mj uj－ ui ρj ∇２ W(xi－xj, h)(７)其中，k为气体状态方程常数；ρ０为流体的静止密度；f p i、f ν i分别为粒子所受的压力和黏性力；pi、pj分别为第i个、第j个粒子的压强；ui、 uj分别为第i个、第j个粒子的速度。

２改进的CSF流体表面张力模型

２．１CSF流体表面张力模型

CSF模型流体表面张力的实施首先是对流场定义色标函数，位于流体粒子位置的色标值为１，位于流体粒子之外位置的色标值为０。通过已定义的色标函数场，可以得到表面法向以及表面曲率，通过表面法向和表面曲率的信息得到表面张力，计算公式为：fs＝ σk( x) ^n(８)其中，σ为表面张力系数；k( x)为表面 x处的曲率；^n为界面的单位法向。法向 n的计算公式为：n＝∇ c( x)(９)其中，c( x)为人为引入的色标函数。曲率计算公式为：k＝－(∇· ^n)(１０)

２．２改进的CSF流体表面张力模型

２．２．１流体表面的切向力

传统的流体表面张力的求解只是对色标函数求梯度，表面张力的方向就是色标函数的梯度方向，也就是说传统的流体表面张力的方向和表面法向相同；传统的流体表面张力模拟将表面的法向力作为表面张力，由于近似计算导致的计算误差使得流体表面张力的模拟很难到达自然界中的平滑状态，会出现凹凸不平，甚至尖角的模拟效果。本文在传统的CSF流体表面张力的基础上，考虑流体表面的法向力和切向力的综合作用，使得流体表面达到足够光滑。流体表面张力的切向力f切向i计算公式为：f切向i＝－∑j mimj ρiρj λ( ni－ nj) W( xi－ xj, h) (１１)其中，λ为一个可调节的参数。对于任意一对相邻的流体粒子i和j，有f切向i← j＝－ mimj ρiρj λ( ni－ nj) W( xi－ xj, h), f切向j← i＝－ f切向 i← j ì î í ïï ïï (１２)其中，f切向i← j、f切向 j← i分别为粒子j对粒子i的切向力、粒子i对粒子j的切向力。因为所施加的表面张力切向力是对称的，所以切向力的施加能够保持系统的动量守恒条件。

２．２．２密度修正

针对CSF流体表面张力模型中流体表面的粒子数目较少，流体表面的流体粒子密度求和的边界截断误差较大，进而导致流体表面张力模拟不稳定的问题，本文对边界粒子的密度进行修正，计算公式为：ρi＝∑ N j＝１ mjWij∑ N j＝１ mj ρj Wij(１３)其中，Wij＝ W( xi－ xj, h)。

２．３本文算法流程

（１）流体粒子初始布局。（２）利用修正的密度求和公式计算流体粒子密度。（３）对所有粒子j，计算其压强pj。（４）对所有粒子j，计算其压力f p j和黏性力f ν j。（５）计算流体粒子j的表面张力。（６）更新流体粒子j的速度与位置。（７）返回步骤（２）。

３仿真实验结果

分别用文献方法、文献方法及本文方法模拟在表面张力作用下水珠的形成过程，对比它们的模拟效果，如图１所示。图１方形液滴在表面张力作用下的变形３４１第１期汪欢欢，等：改进的CSF流体表面张力模型实验在Window７平台上进行，采用Visual Studio２０１３。实验环境为Intel(R)Core(TM) i５Ｇ２４００，３．１０GHzCPU，４G内存，NVIDIA GeＧ forceGT６２０显卡，实验所采用的流体粒子数目为１０２４。模拟总帧数为５００，文献方法、文献方法及本文方法分别用时９．１２、９．８７、９．３４s。从实验结果可以看出，本文方法与文献方法、文献方法的模拟时间差异很小，但本文方法模拟的水珠达到稳定状态时表面更加圆滑，这说明本文方法可以更准确、有效地模拟流体表面张力。

４结论

本文提出一种改进的流体表面张力处理算法，该算法在传统CSF模型流体表面张力的基础上，施加切向力，并且对于流体的边界粒子进行密度修正。从实验结果可以看出，本文提出的改进的CSF流体表面张力算法能够很好地模拟流体表面张力效果；与其他方法模拟结果相比，本文方法改善了流体表面的光滑性。