水平集方法

定义

设计算区域位于 n 维空间中，其中有一张把区域分成两部分的 n-1 维运动界面 Γ(t)。水平集方法引入一个连续函数 φ(x,t)，用它的零等值面表示界面，即 $\Gamma(t)=\{x \mid \phi(x,t)=0\}$，该函数称为水平集函数^[2]。约定界面两侧分别取负值与正值，符号改变之处即界面所在位置^[4]。

多数计算进一步要求 φ 是符号距离函数，即处处满足 $|\nabla\phi|=1$，此时 |φ(x)| 就等于点 x 到界面的最短距离^[4]。这一性质使界面附近的梯度既不过陡也不过平，便于保证精度，也便于按与界面的距离加密网格^[5]。界面的单位法向量与曲率可以直接由 φ 的导数算出：$n=\nabla\phi/|\nabla\phi|$，$\kappa=\nabla\cdot n$^[4]。

背景

在水平集方法出现之前，处理运动界面的常见做法是显式地给出曲线或曲面：用一组节点及连接它们的线或样条表示界面，令节点按 $dx^{(i)}/dt=v(x^{(i)},t)$ 运动，再用差商近似法向量、曲率等几何量^[3]。这类显式方法的不足包括：需要不断重新分布节点并因此引入误差；无法给出尖角处的正确解；界面分裂或合并时必须专门处理；当速度依赖曲率时还受到稳定性条件的限制^[3]。

另外一些方法走的是不同路线。以拉格朗日标记点为基础的界面追踪法在曲率变化平缓时效率较高，但拓扑改变的结果依赖一个影响很大的相互作用时间参数；流体体积法用每个网格单元内各相的体积分数来描述界面，主要困难在于界面重构，二维已相当繁琐，三维代价更高，并因此带来曲率与表面张力计算的不确定性^[6]。

水平集方法属于界面捕捉类方法。它不逐点追踪界面，而是演化覆盖整个区域的函数，界面位置由零等值面隐式给出；因此分裂与融合等拓扑变化是方法本身固有的能力，不需要附加手续^[2]。

内容

当界面随速度场 u 运动时，对 φ(x(t),t)=0 求时间导数可得移流方程 $\phi_t+\mathbf{u}\cdot\nabla\phi=0$，它在水平集方法中被称为水平集方程^[4][2]。若界面只沿法向以速度 F 前进，代入法向量的表达式后方程化为 $\phi_t+F|\nabla\phi|=0$，这是非线性的哈密顿-雅可比方程，其解可能出现不连续，必须附加熵条件才能确定物理上正确的解^[3][2]。

数值求解时，空间导数采用迎风型差分，高阶情形可用 ENO、WENO 等格式，与曲率有关的部分用中心差分；时间推进则常用 TVD Runge-Kutta 格式^[3][4]。由于结果精度主要取决于零等值面附近，可采用窄带法，每个时间步只更新界面邻域内的网格点^[7]。重新初始化时，可先显式求出界面附近节点的距离，其余节点用快速行进法求解程函方程 $|\nabla T|F=1$ 得到，从而在全区构造距离函数^[3][7]。

移流之后 φ 一般不再满足 $|\nabla\phi|=1$，梯度在界面附近会变得过陡或过平，既降低精度又要求更小的时间步长，因此需要周期性重新初始化：保持零等值面不动，另外算出一个具有距离函数性质的新 φ^[3]。常用做法是对重新初始化方程 $\phi_\tau+\mathrm{sign}(\phi_0)(|\nabla\phi|-1)=0$ 在虚构时间上迭代若干步^[3]。但重新初始化会使零等值面发生少量位移，这类偏差会随时间累积，可能使计算出的界面逐渐偏离真实界面^[8]。

影响与争议

方法提出后进入了许多领域。1994 年 Sussman、Smereka 与 Osher 给出了把它用于不可压缩两相流的做法，将水平集函数与 Navier-Stokes 方程耦合求解^[9]。在图像处理中，1997 年 Caselles 等人提出测地线活动轮廓模型，2001 年 Chan 与 Vese 以简化的 Mumford-Shah 泛函为基础建立了分段常值的水平集分割模型（CV 模型），2007 年 Li 等人又提出利用局部灰度信息的 LBF 模型^[10]。此外，它在计算几何与最优化等问题中也有应用^[1]。

标准水平集方法不保证质量守恒，移流格式的离散误差会使某一相的质量随时间增减^[11][12]。为弥补这一点，Chang、Hou、Merriman 与 Osher 曾引入第二个重新初始化方程来抵消质量损失，但该方程把曲率当作非线性速度，当界面曲率变化剧烈时数值解会失稳并出现振荡^[13]。2005 年 Olsson 与 Kreiss 提出守恒水平集方法，改用光滑 Heaviside 函数作为水平集函数，令其在一个相内取 0、另一相内取 1，以 φ=0.5 的等值面代表界面，配合守恒型离散格式显著改善了质量守恒性^[14]。

与其它界面表示耦合是另一条改进途径。粒子水平集方法在界面两侧布置无质量粒子，由穿过界面的粒子指示局部误差并据此修正水平集函数；水平集与流体体积法耦合的 CLSVOF 等格式则试图同时保留两种方法的长处^[5][11]。在图像分割方面，CV 模型存在结果依赖初始化的位置与大小、容易收敛到局部极小以及计算效率偏低等缺陷，后续研究通过凸松弛与全局优化模型加以改进^[15]。

参见

• 流体体积法 —— 另一种界面捕捉方法，用网格单元内的体积分数描述界面，常与水平集方法比较或耦合使用。

• 快速行进法 —— 求解程函方程的快速算法，常用于水平集函数的重新初始化与距离函数构造。

• 符号距离函数 —— 水平集函数通常采用的取值形式，其绝对值表示点到界面的距离。

• 哈密顿-雅可比方程 —— 界面沿法向运动时水平集函数所满足的方程类型。

• 斯坦利·奥舍 —— 与塞西安共同提出水平集方法的数学家。

• 詹姆斯·塞西安 —— 水平集方法的共同提出者，并发展了快速行进法。