格子玻尔兹曼机(Lattice Boltzmann Method)系列1：LBM用于流体模拟的基本原理

1.从粒子碰撞模型到玻尔兹曼方程

首先引入一个分布函数f：

LBM在推理过程中的思想实际上跟分子动力学模拟(MD)或是SPH算法很类似，只是推导的过程会相对地十分繁琐，且方向不同导致模型在尺度上的适用性也相差十分大。这些模型中的粒子都规定了一个速度空间和位置空间，以及相应的时间。区别在于，SPH或是MD中，规定的粒子都是可以移动的。而格子玻尔兹曼算法中，粒子被束缚在相应的网格点上，粒子的移动依赖分布函数、密度等参数之间的信息传递来实现。

上述的公式中的右项，可以被划分为两部分：分布函数的变化一部分来源于外力场的影响(下标为d的项)，另一部分则来源于粒子碰撞(下标为c的项)。

如果有力学系背景的话，在《弹性力学》课程中，第一课就会提到一个十分基础的定理：刘维尔定理(Liouville’s theorem)：

刘维尔定理是《弹性力学》的几个基础定理之一，无法从任何一个已知定理中推理得到，但是也不具有“公理”的显然性和一般性。

接下来再看一个粒子的碰撞模型：

在这个粒子碰撞的模型中，可以总结出一个碰撞的公式，这个公式来源于牛顿的力学碰撞体系。因此我们也可以很容易地发现为什么格子玻尔兹曼机不会适用于纳米级的体系。无论是SPH算法还是MD，在进入纳米级时，都会考虑原子之间的作用势函数，但是LBM并不会将其纳入考量。


以上的两个公式描述了分布函数f控制下，速度空间+位置空间下碰撞粒子数N的微分表达式。

结合我们一开始就提到的分布函数f的公式：粒子碰撞在分布函数f上的体现，得到一个新的公式：

考虑已经提及的分布函数被划分为外力场+粒子碰撞两部分(分布函数的外界影响划分)，并将其中的外力场部分使用Liouville’s theorem进行替代，可得：

其中，

最终得到的玻尔兹曼方程显然是一个很难解的偏微分方程，在数值模拟中，会考虑使用差分算法进行近似：

得到近似后的玻尔兹曼方程：

在LBM中，速度空间和位置空间是独立解耦的，因此上述方程可以等价地写为：

至此，我们已经完整地将玻尔兹曼方程的形式进行了推导。

2.Boltzmann’s H-theorem

我们第二部分要讨论的是玻尔兹曼的H定理，目的是讨论格子玻尔兹曼的最终平衡状态。首先定义一个时间相关的H函数：

以及其关于时间的导数：

结合已知的格子玻尔兹曼方程，可以整合得到一个新的形式：

得到了一个十分复杂的形式，因此第一步会先寻找一个方法尽量做到形式的简化。等号右边一共有3项，首先可以考虑其中的第一项。

很显然的是，由于法向速度在边界的无穿刺基本条件，第一项可以简单地被当作零来消去。接下来考虑等号右边的第二项。

第二项可以直接从数学角度验证其一定为0，因此等号右边现在只剩下了第三项：

再考虑碰撞与反碰撞之间的对偶关系：

将以上两式进行之间的加和平均：

同样也存在一个反碰撞的对偶形式：

再一次进行加和平均：

这就是我们需要的最终的简化形式。显然的是，在任意情况下，都有：

且等号成立当且仅当：

这就是在这一部分我们需要的最终结果，即格子玻尔兹曼机的收敛性验证，它被命名为Boltzmann’s H-theorem。由于这个H函数对于时间的导数永远都是小于或等于零，而这是一个十分类似与《热力学与统计物理》中的熵的性质，因此，只要体系的H达到了它的最小值，我们就会称之为：体系达到了细致平衡(Detailed Balance)。

3.BGK approximation

第一部分中我们推导出来了玻尔兹曼方程，第二部分中借助于H定理讨论了该方程的细致平衡。在本次第三部分中，会讨论LBM的BGK假定，单最终目的还是为了讨论分布函数的具体形式。在一个稳定的体系中，分布函数f的表达将不显含时间项：

回顾第一部分中玻尔兹曼方程的推导过程，我们把分布函数的影响分成了两部分(分布函数的外界影响划分)：一部分是外力项的影响，一部分是粒子碰撞的影响。在格子玻尔兹曼方程中，如果我们假定体系已经达到了稳定状态，并且进行粒子碰撞与外场力之间影响的简单解耦，可得：


在collision部分的公式中，进行公式的等价表述：

如果有《流体力学》的基础知识的话，我们知道，流体的控制方程实际上只涉及到了：质量守恒，动量守恒，以及能量守恒。其中能量守恒部分在不涉及传热流体动力学，而仅仅是普通的多相流体力学模拟的话，据我在各个纯流体组里面学(摸)习(鱼)的经验，实际上几乎不会去考虑能量守恒的部分。

总之由以上的讨论，为了满足三大守恒定律，假设了分布函数的参数拟合形式：

实际上使用的一个满足该三大定律的分布函数形式是：

学过《热力学与统计物理》的话，应该会对玻尔兹曼分布有影响。

接下来终于要开始讨论BGK假定了(逐渐忘记标题)，BGK 假定的形式是：


这个BGK假定是由Bhatnagar, Gross & Krook在1954年提出来的。

由BGK假定可以将上述的分布函数形式进行具体化：

其中u是系统的平均速度。

在上图中使用了张量的形式写公式，有力学/物理背景的同学应该都不难看懂。左边一栏的三项积分由上到下分别对应了质量守恒，动量守恒和能量守恒。而右边一栏的1，3，5项积分则可以进行数值上的对应。

4.N-S方程

流体物理/数学物理背景的同学应该对N-S方程已经十分熟悉了，N-S方程作为流体力学最著名也是最重要的一个公式，在偏微分方程的数值模拟解领域，也是一个学者们经常研究的问题。这一节作为本篇文章的最后一部分，主要目的是为了验证格子玻尔兹曼方法是否真的是一个N-S方程的求解器。

仍然是回顾玻尔兹曼方程：

直接进行速度空间上的积分：

在上一节中的结论：三大守恒定律的分布函数积分描述，可以直接调用这些积分进入公式，得到质量守恒定律：

对玻尔兹曼方程的各项乘以一个速度变量，再进行速度空间的积分：

调用上一节中的结论：三大守恒定律的分布函数积分描述，可以直接调用这些积分进入公式，得到：

考虑积分：

则形式更新为动量守恒定律的形式：

正如前文提到的，由于我是主要做多相流体力学的数值模拟，平时没用做过传热，所以对于能量守恒定律的推导，我自己就没有实践过了，各位有兴趣的话可以仿照刚才的方法进行能量守恒定律的推导。

这篇文章到这里就结束啦，主要就是推导出来了玻尔兹曼方程的形式，玻尔兹曼机的细致平衡条件，分布函数的形式，和验证了格子玻尔兹曼机对于流体物理控制方程的契合性。在BGK假定那里写得有点简略(因为看得实在是太烦了)，有兴趣的话可以去看我提到的那篇原论文。

这个系列会持续更新，后面会再写一写详细的算法，边界条件，代码和多相流模型，码字很浪费时间，希望能坚持下来吧。