【GAMES201学习笔记】物质点法入门

1. 流体力学中的两种描述

1.1 拉格朗日描述

1.1.1 什么是拉格朗日描述

拉格朗日描述 的观察对象是组成物体的 质点 ， 粒子（particle） ，或者说 流体微元 ：

• 随着物体变化，对应 “粒子” 的位置也会随之变化，跟踪的观察对象是单个 “粒子” 的属性变化；
• 拉格朗日描述最显著的特征就是观察对象会发生位置或者形状的变换。

1.1.2 常见的使用拉格朗日描述的方法

• 弹簧质点模型（观察质点的变换）
• 有限元法（观察有限微元的变换）

1.1.3 拉格朗日描述擅长的地方

• 拉格朗日描述擅长对独立粒子的操作；
• 物体的变换操作：这时对粒子施加相应的变换即可；
• 粒子的质量和动量守恒是简单的；
• 定义物体的材质：定义单个粒子相关参数。

1.1.4 拉格朗日描述不擅长的地方

由于粒子经过变换之后，对应的位置和形变都是不确定的，因此虽然可以查询之前已经建立了 邻接信息 的粒子之间的关系（当然，在拉格朗日中找邻接粒子本身是一个困难的问题）；但查询未确立 邻接信息 的粒子之间的关系（例如任意两个粒子之间的距离），会异常困难。

这导致当有查询 未邻接的两个粒子 之间关系的需求时，会遇到一定的麻烦：例如在使用有限元进行仿真时，如果没有合理处理自相交问题，会发生模型穿透的情况。

弹性方块的 FEM 模拟 ：

弹性平面的 FEM 模拟 ：

• 由于仅做了邻接粒子之间、粒子与球之间的碰撞检测，所以在平面折叠时发生了穿透。

不过目前在有限元领域，对于碰撞问题，也已经有了对应的解决方案：

• [SIGGRAPH 2020] Incremental Potential Contact (IPC)

1.2 欧拉描述

1.2.1 什么是欧拉描述

欧拉描述 的观察对象是物体所在 空间的场 ， 网格（Grid） ，或者说 流体元胞（Cell） ：

• 网格是在空间中的绝对参照系的微元，随着物体变化，网格的位置和形状不会随之变化，仅物体映射到网格上的信息发生改变；
• 欧拉描述最显著的特征就是观察对象不会发生位置或者形状的变换。

1.2.2 欧拉描述擅长的地方

• 拉格朗日描述擅长处理粒子之间的相对关系；
• 物体性质的变化：质量密度、速度、温度、熵、焓，甚至单位流体中的磁通量；
• 物体内部的压力压强计算，可以很快的获得能量密度函数；
• 边界碰撞检测：欧拉描述可以很快的获得任意粒子之间的相对关系信息，或者说在欧拉视角下物体自然而然地在进行着碰撞检测。

1.2.3 欧拉描述不擅长的地方

• 物体的变换：显然，拉格朗日描述擅长的地方就是欧拉描述不擅长的地方，在物体发生变换的过程中，需要更新网格和相邻网格的信息，这个过程是较为繁琐的。

2. Particle-In-Cell methods(PIC)

Note ：

• 在混合拉格朗日 - 欧拉方法中，“粒子”是一等公民，“网格”是用来计算场信息的中间量。

2.1 粒子信息映射到网格（P2G）

2.1.1 思路

将粒子信息映射到网格上，或者说从一组粒子重建整个场的信息，这里采用了一个 SPH 中处理类似问题的方法：

• 使用一个径向对称的 核函数 （有时称为 基函数 ），来确定一个粒子对周边空间的 “影响（贡献）”，即一种加权计算 ；
• 其中粒子的影响是局部化的，核函数通常选择 有限支撑 ，即存在一个最大支撑距离；
• 核函数是归一化的： ${\oint }_R\hat{N}(\mathbf{r})d\mathbf{r}=1$ ，其中 $R$ 是整个支撑集。

B 样条核函数 ：

在二维空间下，仅考虑粒子周围的九个网格节点。

2.1.2 代码实现

mpm88.py Line 30,42 （除去了 APIC 的部分）：

for p in x:
        Xp = x[p] / dx
        base = int(Xp - 0.5)
        fx = Xp - base
        # Quadratic B-spline
        w = [0.5 * (1.5 - fx)**2, 0.75 - (fx - 1)**2, 0.5 * (fx - 0.5)**2]

        for i, j in ti.static(ti.ndrange(3, 3)):
            offset = ti.Vector([i, j])
            weight = w[i].x * w[j].y
            grid_v[base + offset] += weight * (p_mass * v[p])
            grid_m[base + offset] += weight * p_mass

Note ：

• 粒子在中间网格里；
• base ：粒子所在周边网格的网格点 ${\mathbf{p}}_{[0,0]}$ ，坐标 (0.5 * dx, 0.5 * dx) ；
• fx ：粒子距离网格点 ${\mathbf{p}}_{[0,0]}$ 的距离；
   
• 距离中间的网格点：$\omega =\frac{3}{4}-|({\mathbf{x}}_p-{\mathbf{x}}_{base})-1{|}^2$
   
• 距离靠左下的网格点：$\omega =\frac{1}{2}(\frac{3}{2}-|{\mathbf{x}}_p-{\mathbf{x}}_{base}|{)}^2$
   
• 距离靠右上的网格点：$\omega =\frac{1}{2}(2-(\frac{3}{2}-|{\mathbf{x}}_p-{\mathbf{x}}_{base}|){)}^2=\frac{1}{2}$