数值双曲型方程的无粘通量分裂方法
好的,我们来详细讲解这个词条。无粘通量分裂方法是计算流体力学和双曲型守恒律方程数值求解中的核心基础技术。它的核心思想是将复杂的、具有方向性的物理通量,分解为几个更简单的、具有明确方向性的分量,从而构造出稳定、高效的数值格式。下面我们循序渐进地展开。
第一步:核心动机——迎风格式的局限与通量的方向性
我们从最简单的一维标量双曲守恒律方程开始:
\[\frac{\partial u}{\partial t} + \frac{\partial f(u)}{\partial x} = 0 \]
其中 \(f(u)\) 是通量函数,其特征速度为 \(a(u) = f'(u)\)。
- 迎风格式的本质:迎风格式之所以稳定,是因为它在空间离散时,始终从信息传播过来的“上游”方向取数据。例如:
- 如果 \(a(u) > 0\)(信息从左向右传播),则应使用向后差分:\(\frac{\partial f}{\partial x} \approx \frac{f_i - f_{i-1}}{\Delta x}\)。
- 如果 \(a(u) < 0\)(信息从右向左传播),则应使用向前差分:\(\frac{\partial f}{\partial x} \approx \frac{f_{i+1} - f_i}{\Delta x}\)。
- 这保证了数值格式与物理信息的传播方向一致,从而具有稳定性。
- 问题所在:对于非线性方程(如欧拉方程),特征速度 \(a(u)\) 不是常数,甚至在同一位置,不同波的特征速度可能方向相反(例如,在跨音速流中,既有亚音速特征也有超音速特征)。此时,我们无法简单地为整个通量 \(f\) 判断一个单一的“迎风方向”。
第二步:核心思想——通量分解
为了克服上述困难,无粘通量分裂方法的核心思想应运而生:将通量函数 \(f(u)\) 分裂为两个部分:
\[f(u) = f^+(u) + f^-(u) \]
并要求这两个部分满足特定的“方向性”条件:
- \(f^+(u)\) 对应正向(或右行)特征:其雅可比矩阵 \(\frac{\partial f^+}{\partial u}\) 的特征值(即分裂后的“正向特征速度”)全部 非负(≥ 0)。
- \(f^-(u)\) 对应负向(或左行)特征:其雅可比矩阵 \(\frac{\partial f^-}{\partial u}\) 的特征值全部 非正(≤ 0)。
几何/物理图像:想象一下在某个网格界面 \(i+1/2\) 处。来自左侧单元 \(i\) 的信息主要由 \(f^+\) 携带向右传播,而来自右侧单元 \(i+1\) 的信息主要由 \(f^-\) 携带向左传播。分裂的目的,就是要将通量中“本该向右走”的部分和“本该向左走”的部分分离开来。
第三步:利用分裂通量构造数值通量
一旦完成了通量分裂,构造数值通量就变得非常自然和稳定。在网格界面 \(i+1/2\) 处,定义数值通量 \(F_{i+1/2}\) 为:
\[F_{i+1/2} = f^+(u_L) + f^-(u_R) \]
其中 \(u_L\) 和 \(u_R\) 分别是界面左侧和右侧的物理量重构值(可以是简单的 \(u_i\) 和 \(u_{i+1}\),也可以是高阶重构结果)。
- 为什么这样构造是迎风的?
- \(f^+(u_L)\) 只使用左侧信息 \(u_L\)。因为 \(f^+\) 的特征速度非负,其信息传播方向是向右的,所以从左侧取值是合理的“迎风”选择。
- \(f^-(u_R)\) 只使用右侧信息 \(u_R\)。因为 \(f^-\) 的特征速度非正,其信息传播方向是向左的,所以从右侧取值是合理的“迎风”选择。
这种构造方法自动保证了格式的迎风特性,是构造高分辨率、无振荡的TVD、ENO/WENO格式的基石之一。
第四步:经典的通量分裂方法
具体如何实现分裂 \(f(u) = f^+(u) + f^-(u)\) 呢?以下是几种经典方法,以标量方程和欧拉方程为例:
- Lax-Friedrichs分裂(全局分裂):
\[ f^{\pm}(u) = \frac{1}{2} (f(u) \pm \alpha u) \]
其中 \(\alpha = \max_u |a(u)|\) 是全局最大特征速度。这是最简单的分裂,\(f^\pm\) 的雅可比矩阵特征值分别为 \((a \pm \alpha)/2\),显然满足非负/非正的条件。它非常鲁棒,但数值耗散较大,因为 \(\alpha\) 是全局最大值。
- Roe分裂(基于局部线性化):
对于方程组,在界面处定义一个平均状态 \(\tilde{u}\)(称为Roe平均状态),使得 \(f(u_R) - f(u_L) = A(\tilde{u})(u_R - u_L)\),其中 \(A = \partial f / \partial u\) 是通量雅可比矩阵。然后将矩阵 \(A(\tilde{u})\) 进行特征分解 \(A = R \Lambda R^{-1}\),其中 \(\Lambda\) 是对角特征值矩阵。接着将特征值矩阵分裂为:
\[ \Lambda = \Lambda^+ + \Lambda^-, \quad \Lambda^\pm = \frac{1}{2}(\Lambda \pm |\Lambda|) \]
最后定义分裂通量为:
\[ f^{\pm}(u) = \frac{1}{2} (f(u) \pm R|\Lambda| R^{-1} u) \]
Roe分裂在界面处是精确的,能很好捕捉激波等间断,但可能在某些情况下违反熵条件(产生非物理的膨胀激波),需要“熵修正”。
- Steger-Warming分裂:
直接将通量雅可比矩阵 \(A(u)\) 在局部状态 \(u\) 处进行特征分解 \(A = R \Lambda R^{-1}\),然后定义:
\[ f^{\pm}(u) = R \Lambda^\pm R^{-1} u \]
其中 \(\Lambda^\pm\) 由正负特征值构成的对角阵。它在单点处是精确的,但在界面处用左右状态分别计算 \(f^+\) 和 \(f^-\) 时,可能不满足连续性,在低马赫数流中可能带来问题。
- Van Leer分裂:
针对欧拉方程,Van Leer设计了一种在声速点(马赫数 \(M=\pm1\))处光滑(可微)的分裂函数。其核心是构造 \(f^\pm\) 的解析表达式,使得分裂后的马赫数多项式函数及其一阶导数在 \(M=\pm1\) 处连续。这种分裂在低马赫数流动和边界层计算中表现更优。
第五步:总结与意义
- 核心作用:无粘通量分裂方法将复杂的、多方向信息传播的物理问题,分解为几个具有明确方向性的子问题,为构建迎风型数值格式提供了系统、坚实的框架。
- 与格式构造的关系:它是构造Godunov型格式(如Roe、HLLC格式)和通量差分裂格式的理论基础之一,同时也是高阶ENO/WENO格式中数值通量构建的关键组成部分。
- 主要考量:选择哪种分裂方法,需要在计算精度(如Roe分裂的尖锐激波分辨率)、数值耗散(如Lax-Friedrichs分裂较耗散)、计算效率和鲁棒性(如对低马赫数流、复杂方程的稳定性)之间进行权衡。