【CG物理模拟系列】流体模拟--粒子法之SPH（理论）

SPH法介绍

SPH(Smoothed Particle Hydrodynamics)是早年银河系碰撞，天体形成等宇宙物理学模拟所使用的方法[Lucy1977]，近年来被应用到流体，热等其他现象处理中。但是，由于是曾经宇宙物理学中常产生的压缩性流体问题的处理方法，并不适用于像水和流速较慢的空气(流速<<音速)这样的非压缩性流体 (强行转换非压缩性的研究近年来也有很多[Becker2007,Solenthaler2009])。

[Desbrun1996]将SPH法第一次应用于CG领域， 此后[Muller2003]将其改进，使其也能适用于3次元的粘性流体。 以下，主要对[Muller2003]的方法进行详细介绍。

SPH法的离散化方程

SPH法中的物理量

的离散化方程如下所示。





这里，N是临近粒子的集合，m是粒子质量，

是粒子密度，
W是kernel方程。内核方程的有效半径是h，有效半径外值为0(Compact Support
)。 同时，它的积分设定为1(

)

根据上述方程，通过周围粒子的权重相加和来近似计算物理量。 另外，物理量的梯度

可以使用内核函数的导数来表示。 物理量

的梯度公式为，





物理量

的laplacian算子为，





使用上述方程，流体的密度方程如下所示。





↑

SPH法中NS方程（Navier-Stokes equations）的解法

支配方程（非压缩性的navier-stokes方程）如下所示．









这里，

是流体速度，

是动粘性系数，

是流体的密度，p为圧力，

是外力，包含重力，表面张力等。
式一是质量守恒方程，式二为动量守恒方程。

先将支配方程用粒子离散化后，通过SPH法求解。为了将粒子表示为液体，只要粒子质量保持不变就可以保证质量的守恒性，因此，不必求解质量守恒方程（这里仅质量守恒，体积不守恒）。质量保存式パーティクル質量が変化しないかぎり質量保存性が常に保持され，質量保存式を解く必要がない(あくまで質量保存であり，体積保存ではないことに注意)．同时，粒子法不必像网格法一样计算平流项，只需移动粒子即可。为此，不必离散化平流项

。所以，各粒子i的速度

的计算公式更新如下。











的各项随后叙述。

↑

粘性扩散项 

粘性项根据流体的速度构成。速度

的laplacian算子离散化后（ 把离散化方程

用

替换），由于粒子i到j 所施加的力与j到i所施加的力不同(asymmetric) 。 不做处理的话，速度快的粒子对速度慢的粒子有影响，反之，影响力小的粒子的速度则会发散。

粘性力便是由流体的速度差来决定的力，通过下式即可实现symmetric的作用。





↑

圧力项

圧力项是由压力值的梯度组成。通过使用梯度值用的离散化方程，虽然能离散化压力项，但是也会发生同粘性扩散项相同的asymmetric问题，根据[Muller2003]，





可以使用压力的平均值来解决。



可以通过表示压力和密度关系的状态方程计算。理想气体的状态方程为，





这里，k是气体常数。大多数的SPH模拟，使用的是[Desbrun1996]提出的公式。





这里，

是模拟对象的流体密度。例如，水的密度为

，编程时，初期设定使用最大密度即可。这个公式简单易懂，但会导致流体压缩。非压缩性强制转换中最常用的方法是 求解Poisson 方程

，但为此也要求解一个庞大的线形系统，非常耗时。

[Becker2007]使用Tait方程，可以确保一定程度的非压缩性。Tait方程如下．





这里，

．

B是圧力定数．





这里，

是音速。通过 

，从

的关系中计算得出。
密度的変化率

常用0.01表示。

这个方法保有很高的非压缩性（1%以下），为了使计算更安定，timestep的值要设定的很小才可以。timestep的大小可以通过CFL(Courant-Friedrichs-Lewy)计算





这里，

是外力，

是粘性定数，根据[Becker2007]取值范围为0.08~0.5。

↑

参考文献 

[Becker2007] M. Becker and M. Teschner, Weakly Compressible SPH for Free Furface Flows, In Proc. SCA2007, pp.209-217, 2007.
[Desbrun1996] M. Desbrun and M.-P. Cani, Smoothed Particles: A New Paradigm for Animating Highly Deformable Bodies, Eurographics Workshop on Computer Animation and Simulation (EGCAS), pp.61-76. 1996.
[Lucy1977] L. B. Lucy, A Numerical Approach to the Testing of the Fission Hypothesis, The Astronomical Journal, Vol.82, No.12, pp.1013-1024, 1977.
[Muller2003] M. Muller, D. Charypar and M. Gross, Particle-based Fluid Simulation for Interactive Applications, In Proc. SCA2003, pp.154-159, 2003.
[Solenthaler2009] B. Solenthaler and R. Pajarola, Predictive-Corrective Incompressible SPH, In Proc. SIGGRAPH2009, pp.1-6, 2009.