BGK 与 MRT 碰撞算子

碰撞算子是格子玻尔兹曼方法的"物理内核"——全部耗散与非平衡弛豫都由它决定。单松弛的 BGK 算子（源自 Bhatnagar-Gross-Krook 1954 对 Boltzmann 方程的模型化）最简单也最流行；多松弛时间（multiple-relaxation-time, MRT）碰撞把不同物理模态的松弛率解耦，是低粘度模拟保持稳定的标准升级。

图:D2Q9 格式中的碰撞 f'=f+Π(f) 与流送方程(来源:QC069)

BGK：单松弛推导

碰撞的最简模型是把全部非平衡偏离按同一速率松弛回平衡态：

$$f_i^{\mathrm{coll}}(\mathbf{x},t)=f_i(\mathbf{x},t)-\frac{1}{\tau }\left[f_i(\mathbf{x},t)-f_i^{\mathrm{eq}}(\mathbf{x},t)\right],$$

无量纲松弛时间 $\tau ={\tau }_c/\Delta t$。它自动继承碰撞不变量：对守恒矩求和时 ${\sum }_i{f}_i^{\mathrm{eq}}=\rho$、${\sum }_i{\mathbf{c}}_i{f}_i^{\mathrm{eq}}=\rho \mathbf{u}$，故碰撞不改变 $\rho$ 与 $\rho \mathbf{u}$——离散层面精确守恒。

$\tau$ 与运动粘度的关系可以从非平衡动量通量的松弛直接看出。定义动量通量张量（二阶矩）${\Pi }_{\alpha \beta }={\sum }_i{c}_{i\alpha }{c}_{i\beta }{f}_i$。对 LBE 取二阶矩、在局部平衡附近线性化，非平衡部分满足

$${\Pi }_{\alpha \beta }^{\mathrm{neq},\mathrm{coll}}=\left(1-\frac{1}{\tau }\right){\Pi }_{\alpha \beta }^{\mathrm{neq}},{\Pi }_{\alpha \beta }^{\mathrm{neq}}\approx -\tau \Delta t\rho c_s^2\left({\partial }_{\alpha }{u}_{\beta }+{\partial }_{\beta }{u}_{\alpha }\right)\cdot \frac{1}{c_s^2\Delta t}\cdots$$

把非平衡通量与牛顿本构 ${\Pi }_{\alpha \beta }^{\mathrm{neq}}=-2\rho \nu S_{\alpha \beta }$ 对齐（完整多尺度推导见 chapman-enskog），得

$$\boxed{\nu =c_s^2\left(\tau -\frac{1}{2}\right)\Delta t}$$

对 D2Q9（$c_s^2=1/3$）即 $\nu =\frac{1}{3}(\tau -\frac{1}{2})\Delta t$。那个 $-\frac{1}{2}$ 是流送步骤自带的隐式扩散（时间离散 $O(\Delta t)$ 截断），不是物理。

稳定性下限 $\tau >1/2$

由上式，$\nu >0$ 要求 $\tau >1/2$；$\tau \to 1/2^{+}$ 对应”零粘度”，理论上诱人，实践上 BGK 在此崩溃：所有九个模态被同一速率松弛，非物理的”幽灵模态”（checkerboard 型奇偶振荡）得不到额外阻尼，随低粘度一起失控，表现为网格量级的锯齿噪声直至发散。经验上 BGK 的稳定区大约在 $\tau \gtrsim 0.55$，这把可模拟的 $Re$ 压在较低水平。

MRT：松弛谱与模态分离

MRT（d’Humières 1992）的关键一步是换到矩空间做碰撞。D2Q9 的九个分布线性组合成九个矩，例如（一种常见选择）

$$\mathbf{m}=(|e{|}^2,e_x^2+e_y^2,e_x^2-e_y^2,e_x{e}_y,e_x,e_y,|e{|}^2{e}_x,|e{|}^2{e}_y,|e{|}^4{)}^{\mathsf{T}}=\mathbf{M}\mathbf{f},$$

分别对应密度、能量、两条对角/交叉剪切应力、动量、两个能量流（幽灵模态）、能量平方（幽灵模态）。碰撞在矩空间逐模态进行：

$$m_i^{\mathrm{coll}}=m_i-s_i\left(m_i-m_i^{\mathrm{eq}}\right),$$

再变换回速度空间完成流送。这样松弛谱 $\{s_i\}$ 自由：物理模态——动量（$s=0$ 不松弛，守恒）与剪切应力（$s_{{\tau }_{\eta }}$ 定粘度 $\nu =c_s^2(\frac{1}{s_{{\tau }_{\eta }}}-\frac{1}{2})\Delta t$）——按流体力学要求设置；幽灵模态可取 $s\to 1$（一步强阻尼）。低粘度时剪切模态松弛慢、幽灵模态松弛快，两者不再互相拖累，稳定域显著扩宽。

图:方腔流算例中不同松弛矩阵 S=diag(0,0,0,1/τ_a,…) 配置的 MRT 验证(来源:QC069)

Wang 与 Shi 对 D2Q9-MRT 做了系统的线性分析：把碰撞算子按特征值分解，逐模态考察松弛时间变化对剪切波、剪切-剪切耦合波传播与衰减的影响，并给出与 Navier-Stokes 一致性及方腔流数值的双重验证Wang & Shi 2015 (arXiv)。实用上还有折中方案 TRT（双松弛时间，只用奇偶两组速率）与带小普朗特数修正的矩格式，谱系综述见Dzanic et al. 2026 (JNNFM)。

与量子路线的关系

碰撞算子是量子 LBM 最关键的一站：流送天然酉（格点平移，见 qlbm-streaming-circuit），而碰撞是非线性的非酉映射，必须借助 Carleman 线性化、量子线性求解或二次量子化嵌入处理——分别对应 carleman-qlbm、qalb 与 bosonic-encoding-lbm 三条路线。

延伸阅读

• lattice-boltzmann
• d2q9-lattice
• chapman-enskog
• qalb
• carleman-qlbm