Beyond awesome | 越而胜己

本文是对 GAMES201线上课程:高级物理引擎实战指南2020 第二课的部分笔记,主要对讲义中复杂的公式做一些注释,以供后续参考。

平滑粒子流体动力学(Smoothed-Particle Hydrodynamics)

SPH原本是天体物理学中采用的方法。其中的核心公式为:

$$A(\mathbf{x}) = \sum_i A_i \frac{m_i}{\rho_i}W(\left| \mathbf{x} - \mathbf{x}_i \right| _2, h)$$

其中,$A$ 可以是任意量,如流速、密度、压强等。特别地,当 $A=\rho$ 即密度时,公式简化为:

$$\rho(\mathbf{x}) = \sum_i m_i W(\left| \mathbf{x} - \mathbf{x}_i \right| _2, h)$$

注意粒子本身是没有密度的,只有质量。这个公式提供了根据粒子的质量和位置来计算流体任意位置密度的方法,适用于自由表面的流体。其中 $W$ 是核函数(kernel,类似于卷积核的权重函数),是一个中间大周围小的偶函数,宽度由参数 $h$ 控制。

一种简单的实现方法是 Weakly Compressible SPH (WCSPH)。注意液体是 incompressible 的,这里的 weakly compressible 是一种近似。

$$\frac{D \mathbf{v}}{D t} = - \frac{1}{\rho} \nabla p + \mathbf{g}, \quad p = B \left( \left( \frac{\rho}{\rho_0}\right)^\gamma - 1 \right)$$

其中,左边指的是加速度受到压强的梯度以及密度的影响;$\mathbf{g}$ 指的是重力加速度,有时也用来表示所有外力给的加速度。$B$ 是体积模量,$\rho_0$ 是理想密度,$\gamma$ 是常数,约为 $7$。这里导数使用的 $D$ 表示材料积分。

WCSPH 中提供了计算 $A$ 的梯度的公式:

$$\nabla A_i = \rho_i \sum_{j} m_j \left( \frac{A_i}{\rho_i ^2} - \frac{A_j}{\rho_j ^2} \right) \nabla_{\mathbf{x}_i} W( \left| \mathbf{x}_i - \mathbf{x}_j \right|_2, h )$$

这个公式并不十分准确,但是由于对称性,能保持动量守恒。另外还有计算 Laplacian 的公式。由此公式可以计算 $\nabla p$。

SPH 的模拟循环如下:

  1. 对每个粒子 $i$,计算 $\rho_i = \sum_j m_j W(\left| \mathbf{x}_i - \mathbf{x}_j \right| _2, h)$;
  2. 对每个粒子 $i$,使用梯度公式计算 $\nabla p_i$;
  3. 计算 Symplectic Euler step:$\mathbf{v}_{t+1} = \mathbf{v}_t + \Delta t \frac{D \mathbf{v}_t}{Dt}, \quad \mathbf{x}_{t+1} = \mathbf{x}_t + \Delta t \mathbf{v}_{t+1}$。

此外还有很多 SPH 的变种。

Courant-Friedlichs-Lewy 条件(CFL condition)

CFL 条件指出了时间步长 $\Delta t$ 的一个上界:

$$C = \frac{u \Delta t}{\Delta x} \leq C_\text{max} \sim 1$$

$C$ 称为 CFL Number 或 Courant Number,$\Delta x$ 是长度的区间(比如粒子半径、网格大小),$u$ 是最大速度。在不同的图形学应用中,$C_\text{max}$ 有不同的取值:SPH 约为 $0.4$,MPM 约为 $0.3 \sim 1$,FLIP 液体(烟雾)约为 $1 \sim 5+$。

加速 SPH:邻居搜索(Neighborhood Search)

通过一个空间数据结构,对每个粒子只需要遍历它周围的粒子,可以将 $O(n^2)$ 降为 $O(n)$。