GAMES103课程笔记12-Smoothed Particle Hydrodynamics
这个系列是GAMES103:基于物理的计算机动画入门(GAMES 103: Intro to Physics Based Animation)的同步课程笔记。本课程是计算机动画的入门课程,着重介绍各种基于物理的动画仿真模拟技术。本节课主要介绍SPH在流体仿真中的应用。
SPH Model
光滑粒子流体动力学(smoothed particle hydrodynamics, SPH)是一种无网格的流体仿真方法。和前面课程中介绍过的方法相比,SPH是一种Lagrangian视角下的模拟方法。在SPH中我们通过大量的粒子和它们携带的物理量(速度、压强、质量)来描述流体的运动。

Smoothed Interpolation
对于流场中的任意位置,SPH通过对它附近粒子求平均的方式来估计该处的物理量。最简单的方法是对半径\(R\)范围内的所有粒子求平均值。

显然这种直接进行平均的方式存在一些缺陷,比如说它没有考虑到粒子的分布。以下图为例,尽管黄颜色的粒子数量比较多,但它们的分布却只集中在一个小范围内,在进行估计时我们不应该过多地考虑这些集中在一起的粒子。

因此更好的方法是使用加权的方式来计算平均。假设每个粒子都具有一定的体积\(V_j\),我们可以利用体积来对不同粒子携带的物理量进行加权。

但是这种加权方法仍然有一些问题,实际上利用体积进行加权后得到的函数可能是不够光滑的。考虑红色和绿色两个位置的函数值,尽管它们在空间上是临近的但由于绿色位置考虑了更多的粒子,它的函数值可能要比红色位置大得多。

因此在SPH中利用核函数(kernel)来做进一步的加权。当粒子离计算位置比较远时它的权重会比较低,而当粒子距离计算位置比较近时则会具有相对较高的权重。

Particle Volume Estimation
那么怎么计算粒子的体积呢?实际上我们可以通过类似的方式来更新粒子的体积。对于每个单独的粒子,它的体积由质量除以密度来表示。当粒子不断运动时空间中的粒子分布会发生变化,因此同样可以通过加权平均的方法来估计当前时刻粒子的密度,进而更新粒子的体积。

将两个更新步骤结合起来就得到了SPH的更新算法:

Smoothing Kernel
最后我们来考虑核函数。在SPH中最常用的核函数为smoothing kernel:

使用核函数的另一个好处在于我们可以把对流场的求导转换为对核函数的求导。这里我们假设每个粒子携带的物理量和它的体积为常量,这样流场的导数为核函数导数的加权和:

而核函数的导数可以显式计算,比如说smoothing kernel的梯度和Laplacian为:


SPH-based Fluids
Fluid Dynamics
在SPH中我们需要更新每个粒子的状态来实现流场的仿真。对于每个粒子需要考虑3个方面的外力:重力(gravity force)、压力(pressure force)以及粘滞力(viscosity force)。
Gravity Force
重力比较简单,只要把重力施加到每个粒子上即可。

Pressure Force
对于压力首先需要考虑每个粒子的压强,它与粒子的分布有关而且可以通过经验公式来估计:

而压力则来自于压强的变化

更严格的表述是压力等于压强的负梯度乘以粒子的体积,我们可以利用核函数的梯度来进行计算:

Viscosity Force
最后一项粘滞力来自于流体自身的粘性,它会使粒子的速度趋于一致:

粘滞力可以通过对速度场求Laplacian来计算,我们同样对核函数来计算Laplacian:

我们把以上几种力施加到粒子上并按照粒子系统进行更新就得到了SPH的模拟算法:

Spatial Partition
不难发现SPH进行仿真时查询每个粒子的邻域会占用大量的计算资源。如果使用暴力计算的话它的计算复杂度为\(O(N^2)\),这对于大规模粒子系统的仿真是不可接受的。

因此在实际进行仿真时需要使用一些数据结构来加速查询的过程,比如说可以利用spatial partion来进行加速:

如果粒子的分布不够均匀还可以考虑八叉树等其它类型的数据结构。

Fluid Display
最后需要说明的是目前市面上的渲染器是无法直接对粒子进行渲染的。因此要对流体进行可视化需要首先将粒子重建为网格,或者是将粒子系统转换为SDF才能进行渲染。
