CN105518681A - 实施各向同性和伽利略不变性的晶格玻尔兹曼碰撞算子 - Google Patents

Abstract

一种方法，包括：在晶格速度集合中，模拟流体体积中的粒子的移动，其中该移动引起粒子之间的碰撞；基于所模拟的移动，确定在该体积内的特定位置处的粒子的相对粒子速度，该相对粒子速度是(i)在该体积内的该特定位置处的并且在该体积的零流动之下测出的粒子的绝对速度与(ii)在该体积内的该特定位置处的粒子中的一个或多个粒子的平均速度之差；及基于相对粒子速度，确定表示碰撞的指定阶的非平衡碰撞后分布函数。

Description

实施各向同性和伽利略不变性的晶格玻尔兹曼碰撞算子

优先权保护

本申请援引35U.S.C.§119(e)要求对2013年7月24日提交的美国临时专利申请序列No.61/858,051的优先权，该申请的全部内容通过引用被结合于此。

背景技术

碰撞过程是多粒子系统中的两个基本动力学过程之一——另一个是平流过程。碰撞过程对于让单个的粒子相互作用并形成集体行为是必不可少的。在碰撞过程中，质量、动量和能量在粒子当中服从守恒定律地进行交换。这些守恒定律确保参与的粒子当中的整体质量和动量(以及有时还有能量)在碰撞前后保持不变。

发明内容

一般而言，本文档描述了用于在晶格速度集合中模拟流体体积中粒子的移动的技术，其中移动引起粒子之中的碰撞；基于所模拟的移动，确定在该体积内的特定位置处的粒子的相对粒子速度，该相对粒子速度是(i)在该体积内的该特定位置处的并且在该体积的零流动之下测出的粒子的绝对速度与(ii)在该体积内的该特定位置处的粒子中的一个或多个粒子的平均速度之差；以及基于该相对粒子速度，确定表示碰撞的指定阶的非平衡碰撞后分布函数。这方面的其它实施例包括对应的计算机系统、装置、机器可读硬件存储设备和记录在一个或多个计算机存储设备上的计算机程序，它们中的每一个被配置为执行该方法的动作和特征。一个或多个计算机的系统可被配置为通过具有安装在系统上的、在操作中使系统执行动作的软件、固件、硬件或者其组合来执行特定的操作或动作。一个或多个计算机程序可被配置为通过包括指令来执行特定的操作或动作，该指令在被数据处理装置执行时使装置执行动作。

前述和其它实施例每个都可以可选地包括单独或组合的一个或多个以下特征。特别地，一种实施例可以组合地包括以下所有特征。这些特征包括：由一个或多个计算机系统提供支持高至粒子速度的阶的流体动力移动的晶格速度集合；其中模拟包括由该一个或多个计算机系统进行模拟。这些特征还包括：对晶格速度集合所支持的阶小于并且不同于非平衡碰撞后分布函数的指定阶；以及其中用于非平衡碰撞后分布函数的指定阶由粒子速度的阶确定。

这些特征还包括：在该体积内的该特定位置处的粒子中的该一个或多个粒子的平均速度包括在该特定位置处的特定类型粒子的平均速度。这些特征还包括：晶格速度集合是与晶格玻尔兹曼方法关联的状态向量的集合。这些特征还包括：非平衡碰撞后分布函数(i)保留用于预定义物理量的非平衡矩，并且(ii)消除用于未定义物理量的高至该指定阶的非平衡矩。这些特征还包括：该指定阶是和流体速度与晶格声速的比相关联的指数值，其中晶格速度集合支持该指数值。这些特征还包括：晶格速度集合包括在被限定到晶格的空间中的动量状态的集合。这些特征还包括：相对粒子速度是从该体积内的特定位置处的并在该体积的零流动之下测出的该粒子的绝对速度中减去的在该体积内的该特定位置处的粒子中的该一个或多个粒子的平均速度。这些特征还包括：非平衡碰撞后分布函数是伽利略不变过滤算子。这些特征还包括：基于非平衡碰撞后分布函数对该流体体积中的粒子的碰撞过程进行建模。这些特征还包括：就用于为流体动力矩提供第一阶支持的晶格速度集合的马赫数(Machnumber)而言，非平衡碰撞后分布函数是第一阶伽利略不变性的碰撞算子并且其中该碰撞算子是根据下式定义的：

$C_i^{\left(1\right)}(x,t)=(1-\frac{1}{T})\frac{w_i}{2T_0}\[(1+\frac{c_i\·u(x,t)}{T_0})(\frac{c_i{c}_i}{T_0}-I)-\frac{c_{i}u(x,t)+u(x,t)c_i}{T_0}\]:{\Π}^{neq}(x,t);$

其中x是该体积内的特定位置；其中t是特定的时间点；其中i是该集合中的晶格速度的索引号；其中T₀是常量晶格温度；其中c_i是粒子在碰撞之前的速度向量；其中u(x,t)是在时间t在特定位置x处的粒子当中的平均速度；其中I是二秩统一张量；其中τ是碰撞关系时间；其中w_i是常量加权因子；并且其中∏^neq是非平衡动量通量。

这些特征还包括：非平衡碰撞后分布函数是用于为流体动力矩提供无限阶支持的晶格速度集合的碰撞算子C_i(x，t)，并且其中该碰撞算子是根据下式定义的：

$C_i(x,t)=(1-\frac{1}{T})\frac{f_i^{eq}(x,t)}{2\mathrm{\ρ}(x,t)T_0}\[\frac{c_i^{\′}(x,t)c_i^{\′}(x,t)}{T_0}-I\]:{\Π}^{neq}(x,t);$

其中x是该体积内的特定位置；其中t是特定的时间点；其中i是该集合中的晶格速度的索引号；其中T₀是常量晶格温度；其中I是二秩统一张量；其中τ是碰撞关系时间；其中c′_i(x，t)是相对粒子速度；其中ρ是流体密度；其中是平衡分布函数；并且其中∏^neq是非平衡动量通量。这些特征还包括：就用于对流体动力矩提供第二阶支持的晶格速度集合的马赫数而言，非平衡碰撞后分布函数是第二阶伽利略不变性的碰撞算子并且其中该碰撞算子是根据下式定义的：

其中x是该体积内的特定位置；其中t是特定的时间点；其中i是该集合中的晶格速度的索引号；其中T₀是常量晶格温度；其中c_i是粒子在碰撞之前的速度向量；其中u(x,t)是在时间t在特定位置x处的粒子当中的平均速度；其中I是二秩统一张量；其中τ是碰撞关系时间；其中wi是常量加权因子；并且其中∏^neq是非平衡动量通量。这些特征还包括：预定义物理量，该预定义物理量包括该特定体积中的流体的质量、该特定体积中的流体的动量以及该特定体积中的流体的能量。

这些特征还包括：非平衡碰撞后分布函数是关于能量通量的碰撞算子C_i(x，t)，并且其中该碰撞算子是根据下式定义的：

其中x是该体积内的特定位置；其中t是特定的时间点；其中i是该集合中的晶格速度的索引号；其中T₀是常量晶格温度；其中I是二秩统一张量；其中τ是碰撞关系时间；其中c′_i(x，t)是相对粒子速度；其中是平衡分布函数；并且其中W^neq是非平衡能量通量。

以上讨论的技术的实现方式可以包括方法或过程、系统或装置，或者计算机可访问的介质上的计算机软件。

这些系统和方法以及技术可以利用各种类型的数值模拟方法来实现，该数值模拟方法诸如用于多相流的Shan-Chen方法和晶格玻尔兹曼公式。本文将描述关于晶格玻尔兹曼公式的进一步信息。但是，本文所述的系统和技术并不限于利用晶格玻尔兹曼公式的模拟，而是可以应用到其它数值模拟方法。

这些系统和技术可以利用采用晶格玻尔兹曼公式的晶格气体模拟来实现。传统的晶格气体模拟假设在每个晶格位点有有限数量的粒子，其中粒子用位的短向量表示。每一位表示在特定方向上移动的一个粒子。例如，向量中的一位可以表示沿特定方向移动的粒子的存在(当设置为1时)或不存在(当设置为0时)。这种向量可以具有六位，例如，值110000指示两个粒子在沿X轴的相反方向移动，并且没有粒子沿Y轴和Z轴移动。一组碰撞规则支配在每个位点处的粒子之间的碰撞行为(例如，110000向量可以变成001100向量，指示沿X轴移动的两个粒子之间的碰撞产生沿Y轴移动的两个粒子)。规则是通过向查找表提供状态向量来实现的，该查找表对位执行置换(例如，将110000变换成001100)。然后，粒子移动到邻接的位点(例如，沿Y轴移动的两个粒子将移至沿Y轴向左和向右的相邻位点)。

在增强的系统中，在每个晶格位点处的状态向量包括多得多的位(例如，用于亚音速流的54位)，以提供粒子能量和移动方向的变化，并且采用涉及满状态向量的子集的碰撞规则。在进一步增强的系统中，多于单个粒子被允许在每个晶格位点或体元(这两个术语贯穿本文档可互换使用)处的每个动量状态中存在。例如，在八位的实现方式中，0-255个粒子可以在特定体元处在特定方向上移动。状态向量是整数的集合(例如，提供在0至255范围内的整数的八位字节的集合)，而不是位的集合，其中每个整数表示在给定状态下的粒子的数量。

在进一步的增强中，晶格玻尔兹曼方法(LBM)使用流体的介观表示在比利用常规计算流体动力(“CFD”)方法可能的层次更深的层次上模拟复杂几何形状中的3D不稳定可压缩紊流过程。在下面提供LBM方法的简要概述。

玻尔兹曼-层次介观表示

在统计物理学中众所周知，流体系统可以通过在所谓“介观”层次上的动力学方程来表示。在这个层次上，不需要确定单个粒子的详细运动。相反，流体的性质由粒子分布函数表示，该粒子分布函数利用单个粒子相空间来进行定义，f＝f(x,v,t)，其中x是空间坐标而v是粒子速度坐标。典型的流体动力量，诸如质量、密度、流体速度和温度，是粒子分布函数的简单矩。粒子分布函数的动力学服从玻尔兹曼方程：

${\∂}_{t}f+v{\▿}_{x}f+F(x,t){\▿}_{v}f=C\left\{f\right\}$ 方程(1)

其中，F(x,t)表示在(x,t)处在外部或自洽生成的主体力。碰撞项C表示各个速度和位置的粒子的相互作用。重要的是要强调，在不指定用于碰撞项C的特定形式的情况下，以上玻尔兹曼方程适用于所有流体系统，而不仅仅适用于众所周知的稀薄气体的情形(如最初由玻尔兹曼构造的)。

一般而言，C包括两点相关函数的复杂多维积分。为了形成具有单独的分布函数f的密闭系统并且为了高效计算的目的，最方便和物理上一致的形式之一是众所周知的BGK算子。BGK算子是根据(不管碰撞的细节)分布函数都经由碰撞接近由{f^eq(x,v,t)}给出的良好定义的局部平衡的物理论点来构造的：

$C=-\frac{1}{\mathrm{\τ}}(f-f^{eq})$ 方程(2)

其中参数τ表示经由碰撞达到平衡的特征弛豫时间。对于粒子(例如，原子或分子)，弛豫时间通常被取为常数。在“混合”(水动力学)表示中，这个弛豫时间是流体动力学变量(像应变率、湍流动能和其它)的函数。因此，湍流可被表示为具有局部确定的特征性质的紊流粒子(“涡旋”)的气体。

玻尔兹曼-BGK方程的数值解比Navier-Stokes方程的解有几个计算优势。首先，可以立即认识到，在该方程中没有复杂的非线性项或高阶空间导数，并且因此几乎不存在关于平流不稳定性的问题。在这个描述层次上，该方程是局部的，因为不需要处理压力，这对算法并行化提供了相当大的优势。连同没有具有第二阶空间导数的扩散算子一起，线性平流算子的另一个期望的特征是它容易以模仿粒子如何真实地与现实中的固体表面相互作用的方式实现诸如无滑动表面或滑动表面的物理边界条件，而不是实现用于流体偏微分方程(“PDE”)的数学条件。其中一个直接好处就是不存在处理固体表面上的界面的移动的问题，这有助于使基于晶格-玻尔兹曼的模拟软件能够成功地模拟复杂湍流空气动力学。此外，来自边界的某些物理性质，诸如有限粗糙度表面，也可以在力中被结合。此外，BGK碰撞算子仅仅是局部的，而自洽主体力的计算可以仅经由近邻信息来完成。因此，玻尔兹曼-BGK方程的计算可以有效地适于并行处理。

晶格玻尔兹曼公式

求解连续玻尔兹曼方程表示在以下方面的重大挑战：它牵涉对在位置和速度相空间中的积分-微分方程的数值评估。当观察到不仅位置而且速度相空间都可以被离散化时，大大的简化发生，这导致了用于玻尔兹曼方程的解的高效数值算法。流体动力量可以以简单求和的形式来写，该简单求和最多依赖于最邻近的信息。即使在历史上晶格玻尔兹曼方程的公式是基于规定粒子在速度的离散集合v(∈{c_i,i＝1,...,b})上的演化的晶格气体模型，这个方程也可以作为连续玻尔兹曼方程的离散化从第一原理系统地得出。其结果是，LBE不遭受与晶格气体方法关联的众所周知的问题的困扰。因此，代替处理相空间中的连续分布函数f(x,v,t),，仅需要跟踪离散分布的有限集合f_i(x,t)，其中下标标示离散速度索引。处理这个动力学方程而不是宏观描述的关键优点是系统的增加的相空间被问题的局部性抵消了。

由于对称性考虑，速度值集合以这样一种方式来选择：该方式使得它们在配置空间中跨越时形成某种晶格结构。这种离散系统的动力学遵循具有形式f_i(x+c_i,t+1)-f_i(x,t)＝C_i(x,t)的LBE，其中碰撞算子通常采取如上所述的BGK形式。通过平衡分布形式的适当选择，可以从理论上示出，晶格玻尔兹曼方程产生正确的流体动力学和热-流体动力学。即，从f_i(x,t)导出的流体动力矩在宏观限制中服从Navier-Stokes方程。这些矩被定义为：

$\mathrm{\ρ}(x,t)=\underset{i}{\Σ}{f}_i(x,t);\mathrm{\ρ}u(x,t)=\underset{i}{\Σ}{c}_i{f}_i(x,t);DT(x,t)=\underset{i}{\Σ}{(c_i-u)}^2{f}_i(x,t)$ 方程(3)

其中，ρ、u和T分别是流体密度、速度和温度，并且D是离散化的速度空间的维度(完全不等于物理空间维度)。

根据包括附图和权利要求的以下描述，其它特征和优点将是显而易见的。

附图说明

图1和图2示出了两个LBM模型的速度分量。

图3是由物理过程模拟系统遵循的规程的流程图。

图4是微观块的透视图。

图5A和图5B是由图3的系统使用的晶格结构的说明。

图6和图7示出了可变分辨率技术。

图8示出了受表面的刻面影响的区域。

图9示出了粒子从体元到表面的移动。

图10示出了粒子从表面到表面的移动。

图11是用于执行表面动力学的规程的流程图。

图12是用于确定指定阶的非平衡碰撞后分布函数的过程的流程图。

图13是用于确定指定阶的非平衡碰撞后分布函数的系统的部件的框图。

具体实施方式

A.保留指定的非平衡移动并消除其它非平衡移动的碰撞算子

在诸如晶格玻尔兹曼模拟的模拟系统中，所模拟的空间被划分成由直线连接的多个离散的点并且因此提供离散的多个点和方向。该模拟也受限于时间步长的离散集合。在这种系统中，为了让模拟近似现实世界流，多个不同的量必须守恒。例如，该系统使质量、动量和能量守恒。因此，该模拟需要被配置为具有适当的质量通量、动量通量和能量通量。这些守恒的量连同它们的通量一起是与该真实物理世界关联的模拟系统中的基本矩。但是，当使这些量守恒时，模拟会由于离散的速度空间(例如，粒子可以按给定的时间步长行进的方向和距离的离散集合)而无意中激发附加的矩量。无意中生成的量(在本文被称为，无意或不期望的不变量，守恒的或非平衡的矩)会不利地影响模拟结果。例如，这种不期望的量会导致错误的流体动力行为和计算结果的数值不稳定性。

为了减少无意中生成的不变量的影响，在本文中描述这样的碰撞算子：该碰撞算子只保留用于守恒的物理量的非平衡矩，同时消除所有其余的非平衡矩直至期望的阶。在这里，碰撞算子的阶是关于晶格马赫数(流体速度和晶格声速的比)的指数定义的。用于确保不期望的非平衡矩将被过滤掉的权衡是增加的计算时间和处理。在本文所述的碰撞算子中，对动量非平衡通量和能量非平衡通量二者都系统地构造高至任意阶的理论形式。因此，这些碰撞算子满足质量、动量、能量守恒，确保高至选定的阶的正确的质量通量、动量通量和能量通量，而同时消除高至选定阶的所有不期望的非平衡矩。在模拟系统中，在处理高速流或者低粘度时，当与较低的阶相比时，对过滤方案选择更高的阶产生更小的非物理效果或影响。

B.对模拟空间建模

在基于LBM的物理过程模拟系统中，流体流可以由按离散速度的集合c_i进行评估的分布函数值f_i表示。分布函数的动力学是由方程4支配的，其中f_i(0)已知为平衡分布函数，其被定义为：

$f_{\mathrm{\α}}^{\left(0\right)}=w_{\mathrm{\α}}\mathrm{\ρ}\[1+u_{\mathrm{\α}}+\frac{u_{\mathrm{\α}}^2-u^2}{2}+\frac{u_{\mathrm{\α}}(u_{\mathrm{\α}}^2-3u^2)}{6}\]$ 方程(4)其中 $u_{\mathrm{\α}}=\frac{c_i\stackrel{\·}{u}}{T}.$

$f_i(\underset{\‾}{x}+{\underset{\‾}{e}}_i,t+1)-f_i(\underset{\‾}{x},t)=\frac{1}{\mathrm{\τ}}\[f_i(\underset{\‾}{x},t)-f_i^{\left(eq\right)}(\underset{\‾}{x},t)\]$ 方程(5)这个方程是描述分布函数f_i的时间演化的众所周知的晶格玻尔兹曼方程。左手侧表示由于所谓的“流化过程”造成的分布的变化。流化过程是一块流体在网格位置出发，然后沿其中一个速度向量移动到下一个网格位置的时候。在那个点，计算“碰撞算子”，即，附近流体块对出发的流体块的影响。流体只能移动到另一个网格位置，因此速度向量的适当选择是必要的，使得所有速度的所有分量都是共同速率(commonspeed)的倍数。

第一个方程的右手侧是上面提到的“碰撞算子”，其表示由于流体块之间的碰撞造成的分布函数的变化。这里使用的碰撞算子的特定形式归功于Bhatnagar、Gross和Krook(BGK)。它迫使分布函数去向由第二个方程给出的规定值，其中第二个方程是“平衡”形式。

根据这个模拟，常规的流体变量，诸如质量ρ和流体速度u，作为方程(3)中的简单求和而被获得。在这里，c_i和w_i的集合值定义LBM模型。LBM模型可以高效地在可扩展计算机平台上实现并且对于时间不恒定流和复杂的边界条件以最大的健壮性运行。

从玻尔兹曼方程获得用于流体系统的运动的宏观方程的标准技术是Chapman-Enskog方法，在该方法中采取全玻尔兹曼方程的逐次逼近。

在流体系统中，密度的小扰动以声速行进。在气体系统中，声速一般由温度确定。流体中可压缩性的效果的重要性是由特征速度与声速的比来测量的，该比被称为马赫数。

参照图1，第一模型(2D-1)100是包括21个速度的二维模型。在这21个速度中，一个速度(105)表示不移动的粒子；三组四个速度表示在沿晶格的x或y轴的正方向或负方向以归一化速率(r)(110-113)、归一化速率的两倍(2r)(120-123)或归一化速率的三倍(3r)(130-133)移动的粒子；以及两组四个速度表示相对于x和y晶格轴二者以归一化速率(r)(140-143)或归一化速率的两倍(2r)(150-153)移动的粒子。

还如图2中所示，第二模型(3D-1)200是包括39个速度的三维模型，其中每个速度由图2的箭头之一表示。在这39个速度中，一个速度表示不移动的粒子；三组六个速度表示在沿晶格的x、y或z轴的正方向或负方向以归一化速率(r)、归一化速率的两倍(2r)或归一化速率的三倍(3r)移动的粒子；八个速度表示相对于x、y、z晶格轴当中的全部三个晶格轴以归一化速率(r)移动的粒子；以及十二个速度表示相对于x、y、z晶格轴当中的两个晶格轴以归一化速率的两倍(2r)移动的粒子。

也可以使用更复杂的模型，诸如包括101个速度的3D-2模型以及包括37个速度的2D-2模型。

对于三维模型3D-2，在101个速度中，一个速度表示不移动的粒子(组1)；三组六个速度表示在沿晶格的x、y或z轴的正方向或负方向以归一化速率(r)、归一化速率的两倍(2r)或归一化速率的三倍(3r)移动的粒子(组2、4和7)；三组八个速度表示相对于x、y、z晶格轴当中的全部三个晶格轴以归一化速率(r)、归一化速率的两倍(2r)或归一化速率的三倍(3r)移动的粒子(组3、8和10)；十二个速度表示相对于x、y、z晶格轴当中的两个晶格轴以归一化速率的两倍(2r)移动的粒子(组6)；二十四个速度表示相对于x、y、z晶格轴的当中两个晶格轴以归一化速率(r)和归一化速率的两倍(2r)移动并且不相对于剩下的轴移动的粒子(组5)；以及二十四个速度表示相对于x、y、z晶格轴当中的两个晶格轴以归一化速率(r)移动并且相对于剩下的轴以归一化速率的三倍(3r)移动的粒子(组9)。

对于二维模式2D-2，在37个速度中，一个速度表示不移动的粒子(组1)；三组四个速度表示在沿晶格的x或y轴的正或负方向以归一化速率(r)、归一化速率的两倍(2r)或归一化速率的三倍(3r)移动的粒子(组2、4和7)；两组四个速度表示相对于x晶格轴和y晶格轴二者以归一化速率(r)或归一化速率的两倍(2r)移动的粒子；八个速度表示相对于x晶格轴和y晶格轴之一以归一化速率(r)移动并且相对于另一个轴以归一化速率的两倍(2r)移动的粒子；以及八个速度表示相对于x晶格轴和y晶格轴之一以归一化速率(r)移动并且相对于另一个轴以归一化速率的三倍(3r)移动的粒子。

上述LBM模型在二维和三维中为流的数值模拟提供了特定类的高效且健壮的离散速度动力学模型。这种类型的模型包括离散速度的特定集合和与那些速度关联的权重。这些速度与速度空间中的笛卡尔坐标的网格点相符，这便于离散速度模型(尤其是被称为晶格玻尔兹曼模型的种类)的准确和高效实现。利用这种模型，流可以以高保真度被模拟。

参照图3，物理过程模拟系统根据规程300来操作，以模拟诸如流体流的物理过程。在模拟之前，模拟空间被建模为体元的集合(步骤302)。通常，利用计算机辅助设计(CAD)程序生成模拟空间。例如，CAD程序可以被用来绘制位于风洞中的微型设备。其后，由CAD程序产生的数据被处理，以添加具有适当分辨率的晶格结构并说明模拟空间中的对象和表面。

晶格的分辨率可以基于正在被模拟的系统的雷诺数(Reynoldsnumber)来选择。雷诺数涉及流的粘度(v)、流中的对象的特征长度(L)和流的特征速度(u)：

Re＝uL/ν方程(6)

对象的特征长度表示对象的大尺度特征。例如，如果微型设备周围的流正在被模拟，则该微型设备的高度可以被认为是特征长度。当对象的小区域(例如，汽车的侧视镜)周围的流是所关心的时，模拟的分辨率可以增加，或者增加分辨率的区域可以在所关心的区域周围被采用。体元的维度随着晶格的分辨率增加而减小。

状态空间被表示为f_i(x,t)，其中f_i表示在时间t在由三维向量x表示的晶格位点处在状态i下的每单位体积的元素或粒子的数量(即，状态i下的粒子的密度)。对于已知的时间增量，粒子的数量被简称为f_i(x)。晶格位点的所有状态的组合被表示为f(x)。

状态的数量由每个能级内可能的速度向量的数量来确定。速度向量由具有三个维度x、y和z的空间中的整数线性速率组成。对于多种属模拟，状态的数量增加。

每个状态i表示处于特定能级(即，能级零、一或二)的不同速度向量。每个状态的速度ci利用其在三个维度当中的每一个维度中的“速率”指示如下：

c_i＝(c_i,x,c_i,y,c_i,z)方程(7)

能级零状态表示在任何维度都不移动的停止的粒子，即，c_stopped＝(0,0,0)。能级一状态表示在三个维度之一中具有±1速率并且在其它两个维度中具有零速率的粒子。能级二状态表示在所有三个维度中都具有±1速率、或者在三个维度之一中具有±2速率并且在其它两个维度中具有零速率的粒子。

生成三个能级的所有可能的排列给出总共39个可能的状态(一个能量零状态、6个能量一状态、8个能量三状态、6个能量四状态、12个能量八状态和6个能量九状态)。

每个体元(即，每个晶格位点)由状态向量f(x)表示。该状态向量完全定义体元的状态并且包括39个条目。这39个条目对应于一个能量零状态、6个能量一状态、8个能量三状态、6个能量四状态、12个能量八状态和6个能量九状态。通过使用这个速度集合，系统可以对实现的平衡状态向量产生麦克斯韦-玻耳兹曼统计。

为了处理效率，体元被分组在被称为微块的2x2x2的体积中。微块被组织成允许体元的并行处理并且最小化与数据结构关联的开销。微块中用于体元的速记符号被定义为N_i(n)，其中n表示微块中晶格位点的相对位置并且n∈{0,1,2,...,7}。微块在图4中示出。

参照图5A和5B，表面S在模拟空间(图5B)中被表示为刻面F_α的集合：

S＝{F_α}方程(8)

其中α是列举特定刻面的索引。刻面不限于体元的边界，但是通常具有量级为与该刻面相邻的体元的尺寸的尺寸，或者具有稍小于与该刻面相邻的体元的尺寸的尺寸，以使得刻面影响相对少量的体元。为了实现表面动力学，向刻面分配属性。特别地，每个刻面F_α具有单位法线(n_α)、表面积(A_α)、中心位置(x_α)和描述刻面的表面动态属性的刻面分布函数(f_i(α))。

参照图6，不同的分辨率水平可以在模拟空间的不同区域中使用，以提高处理效率。通常，对象655周围的区域650是最关心的并且因此利用最高分辨率进行模拟。因为粘度的影响随着离对象的距离而减小，所以采用降低的分辨率水平(即，扩大的体元体积)来模拟在离对象655按增加的距离隔开的区域660，665。类似地，如图7中所示，较低的分辨率水平可以被用来模拟对象775的较不显著特征周围的区域770，而最高的分辨率水平被用来模拟对象775的最显著特征(例如，前沿和后缘表面)周围的区域780。边远区域785利用最低的分辨率水平和最大的体元来模拟。

C.识别受刻面影响的体元

再次参照图3，一旦模拟空间已经被建模(步骤302)，受一个或多个刻面影响的体元就被识别(步骤304)。体元可以以多种方式受刻面影响。首先，被一个或多个刻面相交的体元受影响在于：该体元相对于非相交的体元具有减小的体积。这会发生是因为刻面以及在由该刻面表示的表面下面的材料占据了体元的一部分。分数因子P_f(x)指示体元的不受刻面影响的部分(即，可以被流体或为对其模拟流的其它材料占据的部分)。对于非相交体元，P_f(x)等于1。

通过将粒子传送到刻面或者从刻面接收粒子而与一个或多个刻面相交的体元也被识别为受刻面影响的体元。被刻面相交的所有体元都将包括从刻面接收粒子的至少一个状态以及向刻面传送粒子的至少一个状态。在大多数情况下，附加的体元也将包括这种状态。

参照图8，对于具有非零速度向量c_i的每个状态i，刻面F_α从由平行六面体G_iα定义的区域接收粒子，或向其传送粒子，其中平行六面体G_iα具有由速度向量c_i和刻面(|c_in_i|)的单位法线n_α的向量点积的量值定义的高度以及由刻面的表面积A_α定义的基部，使得平行六面体G_iα的体积V_iα等于：

V_iα＝|c_in_α|A_α方程(9)

当状态的速度向量指向刻面时(|c_in_i|<0)，刻面F_α从体积V_iα接收粒子，并且当状态的速度向量指向远离刻面时(|c_in_i|>0)，向该区域传送粒子。如下面将要讨论的，当另一个刻面占据平行六面体G_iα的一部分时，即，在诸如内角的非凸特征部附近会发生的一种状况，这个表达式必须被修改。

刻面F_α的平行六面体G_iα可以重叠多个体元的部分或全部。体元的数量或其中的部分依赖于相对于体元尺寸的刻面尺寸、状态的能量以及刻面相对于晶格结构的朝向。受影响的体元的数量随着刻面的尺寸而增大。因此，如上面所指出的，刻面的尺寸通常被选择为位于该刻面附近的体元的尺寸的量级或小于位于该刻面附近的体元的尺寸。

被平行六面体G_iα重叠的体元N(x)的部分被定义为V_iα(x)。利用这个术语，在体元V_iα(x)和刻面F_α之间移动的状态i粒子的通量Γ_iα(x)等于该体元中状态i粒子的密度(N_i(x))乘以与该体元重叠的区域的体积(V_iα(x))：

Γ_iα(x)＝N_i(x)V_iα(x)方程(10)

当平行六面体G_iα被一个或多个刻面相交时，以下条件为真：

V_iα＝∑V_α(x)+∑V_iα(β)方程(11)

其中第一个求和说明被G_iα重叠的所有体元并且第二项说明与G_iα相交的所有刻面。当平行六面体G_iα不被另一刻面相交时，这个表达式简化为：

V_iα＝∑V_iα(x)方程(12)

D.执行模拟

一旦受一个或多个刻面影响的体元被识别(步骤304)，定时器就被初始化，以开始模拟(步骤306)。在模拟的每个时间增量期间，粒子从体元到体元的移动由说明粒子与表面刻面的相互作用的平流阶段模拟(步骤308-316)。接下来，碰撞阶段(步骤318)模拟每个体元内粒子的相互作用。其后，定时器递增(步骤320)。如果递增后的定时器不指示模拟完成(步骤322)，则平流阶段和碰撞阶段(步骤308-320)重复。如果递增后的计时器指示模拟已完成(步骤322)，则模拟的结果被存储和/或显示(步骤324)。

1.用于表面的边界条件

为了正确模拟与表面的相互作用，每个刻面必须满足四个边界条件。首先，由刻面接收的粒子的组合质量必须等于由刻面传送的粒子的组合质量(即，到刻面的净质量通量必须等于零)。第二，由刻面接收的粒子的组合能量必须等于由刻面传送的粒子的组合能量(即，到刻面的净能量通量必须等于零)。这两个条件可以通过要求在每个能级(即，能级一和能级二)的净质量通量等于零来满足。

其它两个边界条件与和刻面相互作用的粒子的净动量相关。对于没有皮肤摩擦的表面(在本文被称为滑动表面)，净切向动量通量必须等于零并且净法线动量通量必须等于在刻面处的局部压力。因此，与刻面的法线n_α垂直的组合接收和传送的动量的分量(即，切向分量)必须相等，而与刻面的法线n_α平行的组合接收和传送的动量的分量(即，法线分量)之间的差必须等于在该刻面的局部压力。对于非滑动表面，表面的摩擦相对于由刻面接收的粒子的组合切向动量将由刻面传送的粒子的组合切向动量减小了与摩擦量相关的一因子。

2.从体元到刻面的聚集

作为模拟粒子和表面之间的相互作用的第一步，粒子从体元被聚集并被提供给刻面(步骤308)。如以上所指出的，体元N(x)和刻面F_α之间的状态i粒子的通量为：

Γ_iα(x)＝N_i(x)V_iα(x)方程(13)

由此看来，对于指向刻面F_α的每个状态i(c_in_α<0)，由体元提供给刻面F_α的粒子的数量为：

${\mathrm{\&Gamma;}}_{i\mathrm{\&alpha;}V\&RightArrow;F}=\underset{x}{\&Sigma;}{\mathrm{\&Gamma;}}_{i\mathrm{\&alpha;}}\left(x\right)=\underset{x}{\&Sigma;}{N}_i\left(x\right)V_{i\mathrm{\&alpha;}}\left(x\right)$ 方程(14)

只有其V_iα(x)具有非零值的体元必须被求和。如以上所指出的，刻面的尺寸被选择为使得V_iα(x)仅对于少数体元具有非零值。因为V_iα(x)和P_f(x)可以具有非整数值，所以Γ_α(x)作为实数被存储和处理。

3.从刻面到刻面的移动

接下来，粒子在刻面之间被移动(步骤310)。如果用于刻面F_α的传入状态(c_in_α<0)的平行六面体G_iα被另一刻面F_β相交，则由F_α接收的状态i粒子的一部分将来自刻面F_β。特别地，刻面F_α将接收在前一时间递增期间由刻面F_β产生的状态i粒子的部分。这种关系在图10中示出，其中被刻面F_β相交的平行六面体G_iα的部分1000等于被刻面F_α相交的平行六面体G_iβ的部分1005。如以上所指出的，相交的部分被表示为V_iα(β)。利用这个项，刻面F_β与刻面F_α之间的状态i粒子的通量可以被描述为：

Γ_iα(β,t-1)＝Γ_i(β)V_iα(β)/V_iα方程(15)

其中Γ_i(β,t-1)是在前一时间递增期间由刻面F_β产生的状态i粒子的测量。由此看来，对于指向刻面F_α的每个状态i(c_in_α<0)，由其它刻面提供给刻面F_α的粒子的数量为：

${\mathrm{\&Gamma;}}_{i\mathrm{\&alpha;}F\&RightArrow;F}=\underset{\mathrm{\&beta;}}{\&Sigma;}{\mathrm{\&Gamma;}}_{i\mathrm{\&alpha;}}\left(\mathrm{\&beta;}\right)=\underset{\mathrm{\&beta;}}{\&Sigma;}{\mathrm{\&Gamma;}}_i(\mathrm{\&beta;},t-1)V_{i\mathrm{\&alpha;}}\left(\mathrm{\&beta;}\right)/V_{i\mathrm{\&alpha;}}$ 方程(16)

并且到刻面中的状态i粒子的总通量是：

${\mathrm{\&Gamma;}}_{iIN}\left(\mathrm{\&alpha;}\right)={\mathrm{\&Gamma;}}_{i\mathrm{\&alpha;}V\&RightArrow;F}+{\mathrm{\&Gamma;}}_{i\mathrm{\&alpha;}F\&RightArrow;F}=\underset{x}{\&Sigma;}{N}_i\left(x\right)V_{i\mathrm{\&alpha;}}\left(x\right)\underset{\mathrm{\&beta;}}{\&Sigma;}{\mathrm{\&Gamma;}}_i(\mathrm{\&beta;},t-1)V_{i\mathrm{\&alpha;}}\left(\mathrm{\&beta;}\right)/V_{i\mathrm{\&alpha;}}$ 方程(17)

用于刻面的状态向量N(α)(也被称为刻面分布函数)具有对应于体元状态向量的M个条目的M个条目。M是离散晶格速率的数量。刻面分布函数N(α)的输入状态被设置为等于进入那些状态中的粒子的通量除以体积V_iα：

对于c_in_α<0，N_i(α)＝Γ_iIN(α)/V_iα方程(18)

刻面分布函数是用于从刻面生成输出通量的模拟工具，并且不一定表示实际的粒子。为了生成准确的输出通量，值被分配给分布函数的其它状态。向外的状态是利用上述用于填充向内的状态的技术来填充的：

对于c_in_α≥0，N_i(α)＝Γ_iOTHER(α)/V方程(19)其中Γ_iOTHER(α)是利用上述用于生成Γ_iIN(α)的技术来确定的，但是将该技术应用到除传入状态(c_in_α<0)之外的状态(c_in_α≥0)。在替代的方法中，Γ_iOTHER(α)可以利用来自前一时间步长的Γ_iOTHER(α)的值生成，使得：

Γ_iOTHER(α,t)＝Γ_iOUT(α,t-1)方程(20)

对于平行状态(c_in_α＝0)，V_iα和V_iα(x)都是零。在用于N_i(α)的表达式中，V_iα(x)出现在分子中(根据用于Γ_iOTHER(α)的表达式且V_iα出现在分母中(根据用于N_i(α)的表达式)。因此，当V_iα和V_iα(x)接近零时，用于并行状态的N_i(α)被确定为N_i(α)的极限。

具有零速度的状态(即，静止状态和状态(0,0,0,2)和(0,0,0,-2))的值在模拟的开始基于针对温度和压力的初始条件被初始化。然后，这些值随时间被调整。

4.执行刻面表面动力学

接下来，对每个刻面执行表面动力学，以满足以上讨论的四个边界条件(步骤312)。用于为刻面执行表面动力学的规程在图11中示出。最初，通过确定粒子在刻面处的组合动量P(α)来确定到刻面F_α的组合动量(步骤1105)：

对于所有i，方程(21)

根据这个方程，法线动量P_n(α)被确定为：

P_n(α)＝n_α·P(α)方程(22)

然后，这个法线动量利用推/拉技术被消除(步骤1110)，以产生N_n-(α)。根据这种技术，粒子以只影响法线动量的方式在状态之间被移动。推/拉技术在美国专利No.5,594,671中被描述，其通过引用被结合于此。

其后，N_n-(α)的粒子被碰撞，以产生玻尔兹曼分布N_n-β(α)(步骤1115)。如下面关于执行流体动力学所描述的，玻尔兹曼分布可以通过向N_n-(α)应用一组碰撞规则来实现。

然后，基于传入的通量分布和玻尔兹曼分布，确定用于刻面F_α的传出通量分布(步骤1120)。首先，传入的通量分布Γ_i(α)和玻尔兹曼分布之间的差被确定为：

ΔΓ_i(α)＝Γ_iIN(α)-N_n-βi(α)V_iα方程(23)

使用该差，传出的通量分布是：

对于n_αc_i>0，Γ_iOUT(α)＝N_n-βi(α)V_iα-.Δ.Γ_i*(α)方程(24)并且其中i*是与状态i具有相反方向的状态。例如，如果状态i是(1,1,0,0)，则状态i*是(-1,-1,0,0)。为了说明皮肤摩擦及其它因素，传出通量分布可以进一步细化为：

对于n_αc_i>0，方程(25)

其中C_f是皮肤摩擦的函数，t_iα是与n_α垂直的第一切线向量，t_2α是与n_α和t_1α都垂直的第二切线向量，并且ΔN_j,1和ΔN_j,2是对应于状态i和所指示的切线向量的能量(j)的分布函数。该分布函数根据下式确定：

方程(26)其中j对于能级1状态等于1并且对于能级2状态等于2。

用于Γ_iOUT(α)的方程的每一项的功能如下。第一项和第二项实施法线动量通量边界条件至碰撞已有效产生玻尔兹曼分布但包括切向动量通量异常的程度。第四项和第五项校正这种异常，这种异常可能由于离散效应或因碰撞不足造成的非玻尔兹曼结构引起。最后，第三项添加指定量的皮肤摩擦，以实施在表面上的切线动量通量的期望改变。摩擦系数C_f的生成在下面被描述。需要注意的是，涉及向量操纵的所有项都是可以在开始模拟之前被计算的几何因子。

由此看来，切向速度被确定为：

u_i(α)＝(P(α)-P_n(α)n_α)/ρ，方程(27)

其中ρ是刻面分布的密度：

$\mathrm{\&rho;}=\underset{i}{\&Sigma;}{N}_i\left(\mathrm{\&alpha;}\right)$ 方程(28)

和以前一样，传入的通量分布与玻尔兹曼分布之差被确定为：

ΔΓ_i(α)＝Γ_iIN(α)-N_n-βi(α)V_iα方程(29)

于是，传出通量分布变为：

Γ_iOUT(α)＝N_n-βi(α)V_iα-ΔΓ_i*(α)+C_f(n_αc_i)[N_n-βi*(α)-N_n-βi(α)]V_iα，方程(30)其对应于由之前的技术确定的传出通量分布的前两行，但不需要对异常切向通量的校正。

使用任一种方法，结果产生的通量分布都满足所有的动量通量条件，即：

方程(31)

其中，p_α是在刻面F_α处的平衡压力并且基于向该刻面提供粒子的体元的平均密度和温度值，并且u_α是在该刻面处的平均速度。

为了确保质量和能量边界条件被满足，对每个能级j测量输入能量与输出能量之差：

方程(32)其中索引j表示状态i的能量。然后，这个能量差被用来生成差异项：

对于c_jin_α>0，方程(33)

这个差异项被用来修改传出通量，使得该通量变为：

对于c_jin_α>0，Γ_αjiOUTf＝Γ_αjiOUT+δΓ_αji方程(34)这种操作校正质量和能量通量，同时保留切向动量通量未被更改。如果流在刻面的邻域是大致均匀的并且接近平衡，则这种调整较小。在调整之后，基于邻域平均属性加上由于邻域的非均匀性或非平衡属性造成的校正，结果产生的法线动量通量稍微被更改为是平衡压力的值。

5.从体元到体元的移动

再次参照图3，粒子沿三维直线晶格在体元之间移动(步骤314)。这种体元到体元的移动是对不与刻面相互作用的体元(即，不位于表面附近的体元)执行的唯一移动操作。在典型的模拟中，不位于足够接近表面以便与表面相互作用的体元构成体元的大多数。

每个分开的状态表示在三个维度x、y和z当中的每一个中以整数速率沿晶格移动的粒子。整数速率包括：0，±1和±2。速率的符号指示粒子正在沿对应的轴移动的方向。

对于不与表面相互作用的体元，移动操作在计算上是非常简单的。状态的整个布居(population)在每个时间增量期间从其目前的体元移动到目的地体元。同时，目的地体元的粒子从那个体元移动到其自己的目的地体元。例如，在+1x和+1y方向(1,0,0)上移动的能级1粒子从其当前的体元移动到在x方向是+1并且对其它方向是0的体元。粒子以其在移动之前具有的状态(1,0,0)终止在其目的地体元。基于与其它粒子和表面的局部相互作用，体元内的相互作用将有可能改变用于那个状态的粒子计数。如果不是，则粒子将继续以相同的速率和方向沿晶格移动。

对于与一个或多个表面相互作用的体元，移动操作变得稍微更复杂。这会导致一个或多个分数粒子被传送到刻面。这种分数粒子到刻面的传送导致分数粒子留在体元中。这些分数粒子被传送到被刻面占据的体元中。例如，参照图9，当体元905的状态i粒子的一部分900被移动到刻面910时(步骤308)，剩余的部分915被移动到刻面910在其中定位并且状态i粒子从其指向刻面910的体元920。因此，如果状态布居等于25并且V_iα(x)等于0.25(即，体元的四分之一与平行六面体G_iα相交)，则6.25个粒子将被移动到刻面F_α并且18.75个粒子将被移动到被刻面F_α占据的体元。因为多个刻面会与单个体元相交，所以传送到被一个或多个刻面占据的体元N(f)的状态i粒子的数量为：

$N_i\left(f\right)=N_i\left(x\right)(1-\underset{\mathrm{\&alpha;}}{\&Sigma;}{V}_{i\mathrm{\&alpha;}}(x\left)\right)$ 方程(35)

其中N(x)是源体元。

6.从刻面到体元的分散

接下来，来自每个刻面的传出粒子被分散到体元(步骤316)。本质上，这一步是将粒子从体元移动到刻面的聚集步骤的逆。从F_α刻面移动到体元N(x)的状态i粒子的数量是：

$N_{\mathrm{\&alpha;}iF\&RightArrow;V}=\frac{1}{P_f\left(x\right)}{V}_{\mathrm{\&alpha;}i}\left(x\right){\mathrm{\&Gamma;}}_{{\mathrm{\&alpha;iOUT}}_f}/V_{\mathrm{\&alpha;}i}$ 方程(36)

其中P_f(x)说明部分体元的体积减小。由此看来，对于每个状态i，从刻面指向体元N_(x)的粒子的总数为：

$N_{iF\&RightArrow;V}=\frac{1}{P_f\left(x\right)}\underset{\mathrm{\&alpha;}}{\&Sigma;}{V}_{\mathrm{\&alpha;}i}\left(x\right){\mathrm{\&Gamma;}}_{{\mathrm{\&alpha;iOUT}}_f}/V_{\mathrm{\&alpha;}i}$ 方程(37)

在从刻面到体元分散粒子之后，将粒子与已经从周围的体元平流的粒子组合，并且将结果整数化，有可能某些体元中的某些方向可能下溢(变成负)或者上溢(在八位实现方式中超过255)。这将导致在质量、动量和能量被截断以适应允许的值范围之后这些量的增益或损失。为了防止这种情况的发生，出界的质量、动量和能量在违规状态的截断之前被累加。对于状态所属的能量，等于增益(由于下溢)或损失(由于上溢)的值的质量的量被加回到随机(或顺序)选择的具有相同能量并且本身不受制于上溢或下溢的状态。从这种质量和能量的添加而造成的额外动量被累加并被加到来自截断的动量。通过仅将质量添加到相同的能量状态，当质量计数器达到零时，质量和能量都被校正。最后，动量利用推/拉技术被校正，直到动量累加器返回到零。

7.执行流体动力学

最后，执行流体动力学(步骤318)。这个步骤可以被称为微动力学或体元内操作。类似地，平流规程可以称为体元间操作。以下描述的微动力学操作也可以被用来在刻面碰撞粒子，以产生玻尔兹曼分布。

流体动力学在晶格玻尔兹曼方程模型中通过已知为BGK碰撞模型的特定碰撞算子来确保。这种碰撞模型模仿真实流体系统中的分布的动力学。碰撞过程可以通过方程1和方程2的右手侧很好地得到描述。在平流步骤之后，利用方程3从分布函数获得流体系统的守恒量，具体而言是密度、动量和能量。从这些量，由方程(2)中的f^eq指出的平衡分布函数完全由方程(4)规定。速度向量集合c_i、权重的选择都在表1中列出，连同方程2一起确保宏观行为遵循正确的流体动力学方程。

E.碰撞过程

为了再现相关的流体物理，晶格玻尔兹曼系统中的碰撞过程起着与真实流体系统中相同的基本作用并受制于相同的基本守恒要求。为方便起见，下面的方程将被编号，从方程(1.1)开始。令f_i(x，t)作为“碰撞前”粒子分布函数(即，在时间t每单位体积的在位置x处的粒子数，并且在碰撞之前全都具有速度向量值c_i)，则这个分布在碰撞之后被改变为f′_i(x，t)(即，“碰撞后”分布)。质量守恒如下被满足，

$\underset{i}{\&Sigma;}{f}_i(x,t)=\underset{i}{\&Sigma;}{f}_{\stackrel{\&CenterDot;}{i}}^{\&prime;}(x,t)=\mathrm{\&rho;}(x,t)---\left(1.1\right)$

其中ρ(x，t)是等于在位置x和时间t的所有速度矢量值的粒子数密度的流体密度。(1.1)中(以及后续方程中)的求和是用于晶格玻尔兹曼模型中的所有可能的粒子速度向量值的。动量守恒由下式给出，

$\underset{i}{\&Sigma;}{c}_i{f}_i(x,t)=\underset{i}{\&Sigma;}{c}_i{f}_i^{\&prime;}(x,t)=\mathrm{\&rho;}(x,t)u(x,t)---\left(1.2\right)$

其中u(x,t)是简单地为在位置x和时间t的粒子之中的平均速度的流体速度。

对于其中粒子动能也守恒的某些流体系统(由许多粒子组成)，除了(1.1)和(1.2)，还定义下面的关系，

$\underset{i}{\&Sigma;}{e}_i{f}_i(x,t)=\underset{i}{\&Sigma;}{e}_i{f}_i^{\&prime;}(x,t)=\frac{1}{2}\mathrm{\&rho;}(x,t)\left(u^2\right(x,t)+DT(x,t\left)\right)---\left(1.3\right)$

其中，并且常量D是粒子运动的维度。T(x,t)是流体系统在x和t的温度。

显然，由于守恒定律，ρ(x,t)和u(x,t)(以及对于能量守恒系统还有T(x,t))的值在碰撞过程中不变。对于给定的ρ(x,t)和u(x,t)(以及T(x,t))，还存在特殊类型的粒子分布函数，被称为平衡分布函数。平衡分布函数具有在方程(1.1)和(1.2)(以及(1.3))中定义的相同的质量和动量(以及能量)值，并且它被完全确定为根据ρ(x,t)和u(x,t)(以及T(x,t))的变化而变化。

关于分布函数(对粒子速度向量值)的求和的量一般被称为分布函数的矩。除了对应于质量、动量和能量的守恒宏观量的(1.1)-(1.3)中的三个基本的矩，对应于其通量的矩具有相同的重要性。这些可以关于碰撞后分布函数的求和来定义，

$\mathrm{\&rho;}(x,t)u(x,t)=\underset{i}{\&Sigma;}{c}_i{f}_i^{\&prime;}(x,t)---\left(1.4\right)$

$\&Pi;(x,t)=\underset{i}{\&Sigma;}{c}_i{c}_i{f}_i^{\&prime;}(x,t)---\left(1.5\right)$

$Q(x,t)=\underset{i}{\&Sigma;}{c}_i{c}_i{c}_i{f}_i^{\&prime;}(x,t)---\left(1.6\right)$

在这里，Π(x,t)和Q(x,t)分别是在位置x和时间t的动量和能量通量张量。回忆(1.2)的定义，只有动量通量和能量通量独立于上面的(在(1.1)-(1.3)中定义的)守恒宏观量。

对于平衡中的分布(即，代替地使用(1.5)和(1.6)中的)，平衡动量通量和能量通量分别具有根据(连续)气体的基本动力学理论而众所周知的形式

Π^eq(x，t)＝p(x，t)u(x，t)u(x，t)+p(x，t)I(1.7)

$Tr\&lsqb;Q^{eq}(x,t)\&rsqb;=\frac{1}{2}\mathrm{\&rho;}(x,t)\&lsqb;(D+2)T(x,t)u(x,t)+u^2(x,t)u(x,t)\&rsqb;---\left(1.8\right)$

其中，p(x,t)是压力，p(x,t)＝ρ(x,t)T(x,t)基于理想的气体定律，并且(1.7)中的“I”表示2秩统一张量。Tr[.]是求迹操作。晶格玻尔兹曼方法的中心主题是恢复(1.7)和(1.8)的形式。

(1.5)和(1.6)的矩也可以被表示为具有平衡部分和非平衡部分的项，

Π(x，t)＝Π^eq(x，t)+Π′^neq(x，t)

Q(x，t)＝Q^eq(x，t)+Q′^neq(x，t)

显然，

${\&Pi;}^{\&prime;neq}(x,t)=\underset{i}{\&Sigma;}{c}_i{c}_i{f}_i^{neq}(x,t)=\underset{i}{\&Sigma;}{c}_i{c}_i\&lsqb;f_i^{\&prime;}(x,t)-f_i^{eq}(x,t)\&rsqb;---\left(1.9\right)$

$Q^{\&prime;neq}(x,t)=\underset{i}{\&Sigma;}{c}_i{c}_{\stackrel{\&CenterDot;}{i}}{c}_i{f}_i^{neq}(x,t)=\underset{i}{\&Sigma;}{c}_i{c}_i{c}_i\&lsqb;f_i^{\&prime;}(x,t)-f_i^{eq}(x,t)\&rsqb;---\left(1.10\right)$

非平衡矩通量在确定流体系统中的运输行为时起关键作用。例如，非平衡矩通量(1.9)直接确定牛顿流体系统中的粘度，而(1.10)的非平衡矩通量确定热扩散性。虽然可以有从粒子分布中构造出的无限数量的矩，但是只有由(1.1)-(1.10)给出的那些矩在物理流体系统中是宏观相关的。

该系统将与平衡的偏离度表示为然后物理稳定的碰撞过程总是在减少偏离的方向上工作。即，与平衡的碰撞后偏离比碰撞前的小，

$|f_i^{\&prime;}(x,t)-f_i^{eq}(x,t)|\&le;|f_i(x,t)-f_i^{eq}(x,t)|---\left(1.11\right)$

现实的多粒子系统中的碰撞会相当复杂。最简单的碰撞模型(算子)既满足守恒定律要求又收敛到(1.11)的平衡要求。“BGK”碰撞过程在数学上被表示如下，

$f_i^{\&prime;}(x,t)=f_i^{eq}(x,t)+C_i(x,t)---\left(1.12\right)$

其中碰撞算子由下式给出：

$C_i(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\left(f_i\right(x,t)-L_i^{eq}(x,t\left)\right)---\left(1.13\right)$

因为所有这三个分布函数都给出相同的ρ(x,t)和u(x,t)(以及对于能量守恒系统还有T(x,t))，所以质量和动量(和能量)守恒(即(1.1)-(1.3))被自动满足。此外，根据(1.12)，碰撞后分布与平衡的偏离与碰撞前与平衡的偏离按因子(1-1/τ)成比例。因此，只要参数值τ(被称为碰撞弛豫时间)大于1/2并且碰撞后偏离在τ＝1时消失，收敛条件(1.11)就被满足。

根据BGK关系(1.12)和(1.13)，方程(1.9)和(1.10)可以被改写为，

${\&Pi;}^{\&prime;neq}(x,t)=(1-\frac{1}{\mathrm{\&tau;}}){\&Pi;}^{neq}(x,t)---\left(1.14\right)$

$Q^{\&prime;neq}(x,t)=(1-\frac{1}{\mathrm{\&tau;}})Q^{neq}(x,t)---\left(1.15\right)$

其中Π^neq(x,t)和Q^neq(x,t)分别是下面定义的碰撞前非平衡动量通量和能量通量，

${\&Pi;}^{neq}(x,t)=\underset{i}{\&Sigma;}{c}_i{c}_i\&lsqb;f_i(x,t)-f_i^{eq}(x,t)\&rsqb;---\left(1.16\right)$

$Q^{neq}(x,t)=\underset{i}{\&Sigma;}{c}_i{c}_i{c}_{\stackrel{\&CenterDot;}{i}}\&lsqb;f_i(x,t)-f_i^{eq}(x,t)\&rsqb;---\left(1.17\right)$

一旦τ的值被选定，晶格玻尔兹曼流体中具有BGK碰撞算子(1.12)的移动粘度值就被确定为，

v＝(τ-0.5)T₀(1.18)

在这里，常量T₀是标准晶格温度，并且所谓的晶格单位约定被使用，使得晶格单元中的时间增量是统一的。BGK碰撞算子一直是晶格玻尔兹曼模型中最常用的算子。它已在过去的20多年中显示出各种程度的成功。另一方面，BGK算子有一些固有的局限性。除了统一普朗特数(＝粘度与热扩散性之比)，除了(1.9)-(1.10)的基本问题，还有其中一个问题是所有的碰撞后非平衡矩都在τ不统一的任何时候生成。事实上，对于BGK，前非平衡矩和后非平衡矩都表现出以下关系，

$M_n^{\&prime;}(x,t)=(1-\frac{1}{\mathrm{\&tau;}})*M_n(x,t)---\left(1.19\right)$

其中M′_n(x，t)和M_n(x，t)分别表示任何n阶矩的碰撞后和碰撞前。

F.过滤后的碰撞算子

任何晶格玻尔兹曼模型的最通用特征都是它具有有限且恒定的粒子速度向量值集合。对于给定的晶格玻尔兹曼模型，指定常量向量值的集合。其结果是，只有从离散速度集合中构造出的粒子分布函数的矩的有限集合可以在现实流体中恢复它们的对等部分。用于恢复现实流体的高至任意给定阶的矩的通用框架是严格定义的。不同的晶格玻尔兹曼向量值集合可以支持高至不同的阶的物理矩。例如，所谓的D3Q19和D3Q15仅在低马赫数限制的情况下支持高至平衡动量通量和线性偏离(牛顿)非平衡动量通量的矩。另一方面，所谓的更高阶晶格玻尔兹曼模型，诸如D3Q39，可以在有限马赫数和Knudsen数的情况下支持高至超过线性偏离的平衡能量通量和非平衡动量通量的矩。

由于仅物理相关的矩是平衡和非平衡动量和能量矩以及它们的通量，因此，为了实现正确的物理流体行为和更大的数值稳定性，期望设计碰撞算子，使得：除非平衡动量通量(和可能的能量通量)之外，所有碰撞后非平衡矩都消失。具体而言，非平衡动量通量、能量通量根据下面的方程表示：

${\&Pi;}^{\&prime;neq}(x,t)=(1-\frac{1}{\mathrm{\&tau;}}){\&Pi;}^{neq}(x,t)---\left(2.1\right)$

$Q^{\&prime;neq}(x,t)=(1-\frac{1}{{\mathrm{\&tau;}}_e})Q^{neq}(x,t)---\left(2.2\right)$

但其它碰撞后非平衡矩都消失，

M′_n(x，t)＝0.(2.3)

回忆(1.19)，BGK碰撞生成所有非平衡矩，实现(2.1)-(2.3)的碰撞算子也被称为过滤后的碰撞算子，因为它滤掉了不必要的非平衡矩。在消除所有其它非平衡矩的同时，过滤后的碰撞算子保留如(2.1)的非平衡动量通量。这种过滤后的碰撞算子的明确表达式在下面给出，

$C_i(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{w_i}{2T_0}\&lsqb;\frac{c_i{c}_i}{T_0}-I\&rsqb;:{\&Pi;}^{neq}(x,t)---\left(2.4\right)$

其中，常量w_i是一旦特定的晶格速度集合被选定就被确定的加权因子值。晶格速度集合是在局限到晶格的空间中的微观速度或动量的离散集合。该加权因子值的集合对于不同的晶格玻尔兹曼粒子速度集合是不同的，并且其目的是实现到希望阶的矩各向同性。由于新的(过滤后的)碰撞算子(2.4)产生与BGK相同的碰撞后非平衡动量通量(1.9)，因此它在(1.18)中自动引起与BGK相同的粘性效果。

对于某些晶格系统，获取不消失的非平衡能量通量也是期望的，如在(2.2)示为具有非统一τ_e。导出被示为实现这种目的的特定碰撞算子，

注意，值τ_e可能不需要等于τ，使得普朗特数不限于统一，这与BGK中相反。由于动量通量与能量通量之间的矩正交性，更广义的过滤后的碰撞形式是通过(2.4)和(2.5)中的两种形式的直接求和构造的，使得(2.1)和(2.2)中的非平衡通量被自动满足，而其余的非平衡矩消失。

通过过滤掉不被给定晶格速度集合支持的不想要的非平衡矩，广义的碰撞算子((2.4)、(2.5)或它们的组合)已被证明比BGK有显著改进的流体流各向同性，但它扔保留像BGK那样的期望的动量通量和能量通量。重要的是，广义碰撞算子形式(2.4)和(2.5)(以及它们的组合)的类不仅适用于Navier-Stokes流体的粘度和热扩散性，而且还确保在更宽流体区域内涉及有限Knudsen数的正确流体动力学。

G.用于导出N阶伽利略不变过滤后的算子的规程

与本公开内容一致的系统基于何种量的速率可被模拟支持来生成过滤后的碰撞算子，该过滤后的碰撞算子具有关于晶格马赫数的指定指数。这个过滤后的碰撞算子使用(例如，速率或能量的)相对值，而不是(例如，速率或能量的)绝对值。本文描述的是为不包括不想要的矩的碰撞算子适当地构造正确理论形式的规程。非平衡矩的量对应于可被模拟支持的速率的量。因此，过滤后的碰撞算子表示由给定的晶格玻尔兹曼模型的粒子速率集合支持的非平衡矩，因此非平衡矩对应于在物理世界中实际发生的情况。相对速率的使用使得能够构造具有想要的非平衡贡献的过滤后的碰撞算子，其中想要的非平衡贡献由模拟的速率支持。

即，通过分别关于相对粒子速度和能量确定非平衡动量和能量通量，可以构造过滤掉(例如，排除)在模拟中不被晶格玻尔兹曼模型的晶格速率集合支持的更高阶项的过滤后的碰撞算子。在例子中，碰撞算子被识别为延展形式，因此，对于多个晶格速度，用户或系统有办法知道高至哪一阶的哪些项应当被保留并且超出那一阶的项应当(从碰撞算子中)被消除。为了这么做，该延展(例如，碰撞算子)必须以相对速度来表示。

在例子中，速度模型具有粒子速率的有限集合。因此，粒子移动的真实物理矩只能在模型中准确表示到某一阶。为了在更高的流速下具有伽利略不变碰撞操作和模拟速率，碰撞算子的形式基于相对速度，其中粒子速度是针对其自身的流速测量的。至少由于这个原因，参数c′_i(表示相对速度或能量)在下面的方程3.5中被使用，而不是使用方程3.5中的参数c_i。参数c_i在方程3.5中的使用不允许按流速(例如，u(x,t))的幂的碰撞形式的扩展，并且因此就对局部流速的相对粒子速度而言不提供根据伽利略不变性的适当碰撞操作的构造。但是，参数c′_i在方程3.5中的使用允许按流速(例如，u(x,t))的幂(例如，N阶)的碰撞形式的扩展，并且因此提供N阶伽利略不变过滤后的算子。

在例子中，依据BGK模型的碰撞算子包括无限阶的非平衡矩。高至某一阶的非平衡矩变得非物理，因为某些速度模型不能正确支持高至无限阶的矩。因此，这个系统必须消除高阶非平衡矩，以防止包括物理假象。例如，用于19速率模型的相对矩只高至第一阶。因此，更高阶矩是不相关的。但是，当研究模拟利用BGK算子运行时，超出由模型支持的那些矩的非平衡矩被包括在内。因此，本文所描述的过滤后的碰撞算子更紧密地对应于在物理世界中发生的情况。

上述过滤后的碰撞形式只用于处于非常小流体速度的流的情况。根据伽利略变换不变性的基本原理，物理流体的所有统计性质都应当只随对于平均流体速度的粒子相对速度的变化而变化。具体而言，物理流体系统中的粒子非平衡分布函数和相关的非平衡矩应当仅依赖于相对速度(c_i-u(x,t))，而不是(在晶格玻尔兹曼系统中在零流体参考帧中测出的)绝对值c_i。因此，代替(1.14)和(1.15)，物理上更有意义的非平衡动量通量和能量通量分别是

${\&Pi;}^{neq}(x,t)=\underset{i}{\&Sigma;}{c}_i^{\&prime;}(x,t)c_i^{\&prime;}(x,t)\&lsqb;f_i(x,t)-f_i^{eq}(x,t)\&rsqb;---\left(3.1\right)$

$Q^{neq}(x,t)=\underset{i}{\&Sigma;}{c}_i^{\&prime;}(x,t)c_i^{\&prime;}(x,t)c_i^{\&prime;}(x,t)\&lsqb;f_i(x,t)-f_i^{eq}(x,t)\&rsqb;---\left(3.2\right)$

其中相对粒子速度和能量由下式给出，

c′_i(x，t)＝c_i-u(x，t)(3.3)

如在上面的方程3.3中所示，绝对粒子速度被相对速度代替。有趣的是，由于质量和动量守恒，非平衡动量通量(3.1)被证明与(1.16)相同。另一方面，非平衡能量通量(3.2)不能简化成(1.17)。

与本公开内容一致的系统为非平衡动量通量创建伽利略不变碰撞算子。根据动力学理论的基本物理规律，由于流体速度不均匀性造成的主导阶非平衡粒子分布可以表示为，

$f^{neq}~f^{eq}*\&lsqb;UU-\frac{1}{D}{U}^{2}I\&rsqb;:{\&Pi;}^{neq}---\left(3.4\right)$

其中U是连续动力学理论中相对局部平均流体的粒子相对速度。受这个概念的启发，本文所描述的系统在晶格玻尔兹曼中识别用于根据伽利略不变性的非平衡分布的类似表达式。该系统将以下的明确形式识别为完全对应的非平衡碰撞后分布函数，

$C_i(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{f_i^{eq}(x,t)}{2\mathrm{\&rho;}(x,t)T_0}\&lsqb;\frac{c_i^{\&prime;}(x,t)c_i^{\&prime;}(x,t)}{T_0}-I\&rsqb;:{\&Pi;}^{neq}(x,t)---\left(3.5\right)$

其中Π^neq(x,t)由(1.14)或(3.1)给出。如在上面的方程3.2中所示，分布函数既包括平衡分量又包括非平衡分量。此外，上面的方程需要到无限阶的粒子速度集合。

(3.5)中的平衡分布函数是用于任意马赫数的完整形式，

$f_i^{eq}(x,t)=w_{,i}{e}^{(c_i^2/2T_0)}\mathrm{\&rho;}(x,t)e^{\&lsqb;-c_i^{\&prime;}{(x,t)}^2/2T_0\&rsqb;}---\left(3.6\right)$

(3.5)和(3.6)中分别用于非平衡分布函数和平衡分布函数的完整形式只有在晶格速度集合支持到所有阶的正确流体动力学矩时才可实现。换句话说，对于具有粒子速度值的有限集合的任意给定的晶格速度集合，该完整形式是不能实现的。但是，这些可以利用至对应阶N的给定的有限晶格速度集合来实现，如下面详细描述的。

表达式(3.5)和(3.6)可以以u(x,t)的幂的多项式形式进行扩展。众所周知，平衡分布(3.6)可以用一系列Hermite多项式表示，

$f_i^{eq}(x,t)=w_i\mathrm{\&rho;}(x,t)\underset{n=0}{\overset{\mathrm{\&infin;}}{\&Sigma;}}\frac{H^{\left(n\right)}\left({\mathrm{\&xi;}}_i\right)}{n!}V{(x,t)}^{\&lsqb;n\&rsqb;}---\left(3.7\right)$

在这里，V(x,t)^[n]是向量V(x,t)的n次方直接乘积的短记号。和第n阶Hermite多项式H⁽ⁿ)(ξ_i)是标准(标量)n阶Hermite函数的n秩张量推广。在没有证据的情况下，非平衡贡献(3.5)也可以用一系列Hermite多项式表示，

$C_i(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{w_i}{2T_0}\underset{n=2}{\overset{\mathrm{\&infin;}}{\&Sigma;}}\frac{H^{\left(n\right)}\left({\mathrm{\&xi;}}_i\right)}{n!}V{(x,t)}^{\&lsqb;n-2\&rsqb;}{\&Pi;}^{neq}(x,t)---\left(3.8\right)$

利用Hermite多项式的正交性质，系统通过简单地截断两个无穷级数(3.7)和(3.8)获得对这些形式的N阶近似，从而对于m≤N保留与u^m(x,t)成比例的项。一个速度模型是19速率立方D3Q19晶格，其将每个晶格点连接到其第一邻居和第二邻居。具体而言，对于D3Q19(或对于D3Q15)类型的晶格速度集合，系统扩展(3.6)中的至u³(x,t)，并且扩展(3.5)中的C_i(x,t)至u(x,t)的第一(线性)幂：

$C_i^{\left(1\right)}(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{w_i}{2T_0}\&lsqb;(1+\frac{c_i\&CenterDot;u(x,t)}{T_0})(\frac{c_i{c}_i}{T_0}-I)-\frac{c_{i}u(x,t)+u(x,t)c_i}{T_0}\&rsqb;:{\&Pi;}^{neq}(x,t)---\left(3.9\right)$

如在上面的方程3.9中所示，这个碰撞操作对于可以支持十九个速率的模拟系统是准确的。此外，这个碰撞算子是(3.5)中C_i(x,t)高至第一(线性)幂的延展形式，以便对于可以支持十九个速率的模拟对应于物理世界中发生的情况。这个碰撞算子可以在方程1.12中用来修改分布函数。在上面的方程3.9中，x是体积内的特定位置，t是特定的时间点，i是集合中的晶格速度的索引号；T₀是常量晶格温度；c_i是所述粒子在碰撞之前的速度向量；u(x,t)是在时间t在特定位置x处的粒子当中的平均速度；I是二秩统一张量；τ是碰撞关系时间；w_i是常量加权因子；并且Π^neq是非平衡动量通量。

另一速度模型是D3Q39晶格，其支持高至三十九个速率。对于D3Q39，扩展可以执行到u(x,t)的平方或立方幂。至u²(x,t)的截断在下面明确给出，

$\begin{array}{c}{C}_i^{\left(2\right)}(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{w_i}{2T_0}\&lsqb;(1+\frac{c_i\&CenterDot;u(x,t)}{T_0}+\frac{{(c_i\&CenterDot;u(x,t\left)\right)}^2}{2T_0^2}-\frac{u^2(x,t)}{2T_0})(\frac{c_i{c}_i}{T_0}-I)\\ -(1+\frac{c_i\&CenterDot;U(x,t)}{T_0})\frac{c_{i}u(x,t)+u(x,t)c_i}{T_0}+u(x,t)u(x,t)\&rsqb;:{\&Pi;}^{neq}(x,t)\end{array}---\left(3.10\right)$

如在上述方程3.10中所示，通过在扩展后的形式中不包括高于平方幂并且因此不被模拟支持的更高阶矩，碰撞算子的扩展形式排除了不想要的矩。在这个例子中，x是体积内的特定位置；t是特定的时间点，i是集合中晶格速度的索引号；T₀是常量晶格温度；c_i是所述粒子在碰撞之前的速度向量；u(x,t)是在时间t在特定位置x处的粒子当中的平均速度；I是二秩统一张量；τ是碰撞关系时间；w_i是常量加权因子；并且Π^neq是非平衡动量通量。

因此，系统将(2.4)的碰撞形式重新解释为(3.5)的第0阶近似(即，C_i ⁽⁰⁾(x,t))，因为它仅包括来自(3.5)的独立于u(x,t)的项。对于诸如D3Q39的更高阶晶格速度集合，可以在平衡分布函数(3.6)中保留高至u⁵(x,t)的项，而在碰撞后非平衡分布(3.5)中保留高至u³(x,t)的项。一般而言，当晶格速度集合对流体动力矩提供恰当的对应阶的支持时，这个系统性过程可以被执行到(按u(x,t)的幂的)任何任意阶。虽然，在给定的有限阶，伽利略不变性被确切地满足，但是随着越来越高的晶格速度集合被使用以及(3.5)和(3.6)中越来越高阶的扩展形式被使用，误差朝越来越高的阶移动。

在上面，系统生成并执行系统性规程，以构造对于按流体速度的幂的任何给定阶的广义过滤后的碰撞算子。特别地，对于碰撞后非平衡通量，对预先存在的过滤器算子(2.4)的第一阶校正和第二阶校正分别在(3.9)和(3.10)中明确地表示。除了期望的那些矩之外，过滤后的碰撞算子过滤掉非平衡矩。过滤器算子(2.4)和(3.9)(或(3.10))都起到保留非平衡动量通量的目的。另一方面，虽然(3.9)(和(3.10))给出了与BGK和(2.4)相同的非平衡动量通量(和粘度值)，但这种碰撞算子实现了改进的数值稳定性和伽利略不变性。

相同的规程也可以为关于能源通量的过滤后的碰撞算子制定。通用的完全伽利略不变形式类似于(3.5)中的形式，关于相对速度c′_i(x，t)被表示为：

其中W^neq(x,t)是(3.1)和(3.2)的Q^neq(x,t)和u(x,t)Π^neq(x,t)的适当线性组合。根据相同的规程，可以系统地获得由给定晶格速度集合充分支持的按u(x,t)的幂的任何有限阶形式。由于动量通量与能量通量之间的矩正交性，通用的碰撞后形式被简单地产生为(3.5)和(3.11)的(以适当的扩展形式)相加，其能够独立地实现期望的粘度和热扩散性。在上面的方程3.11中，x是体积内的特定位置；t是在特定的时间点；i是集合中的晶格速度的索引号；T₀是常量晶格温度；I是二秩统一张量；τ是碰撞关系时间；c′_i(x，t)是相对粒子速度；是平衡分布函数；并且W^neq是非平衡能量通量。

参照图12，与本公开内容一致的系统在确定非平衡碰撞后分布函数时执行过程1200。在操作中，系统提供(或以其它方式获得)(1202)支持高至特定阶的粒子速度的流体动力学移动的晶格速度集合。例如，系统获得支持流体动力学移动的高至第一阶的19个速率的D3Q19模型。在这个例子中，对晶格速度集合所支持的阶(例如，19个速率)小于并且不同于非平衡碰撞后分布函数的指定阶(例如，第一阶或线性)，并且非平衡碰撞后分布函数的指定阶是由粒子速度的阶确定的。即，对于支持19个速率的速度模型，非平衡碰撞后分布函数是线性的或者只支持高至第一阶(因为19个速度速率与第一阶流体速度(即，u(x,t))关联)。在这个例子中，系统被配置为访问(在数据存储库中或在另一个系统中)对于特定速度模型(或对于特定数量的速度速率)的流体速度的被支持的阶。在这个例子中，系统访问表或对于这个信息的其它映射。

在例子中，晶格速度集合是与晶格玻尔兹曼方法关联的状态向量的集合。在这个例子中，状态向量是指示在晶格位点的粒子之间的碰撞行为的一系列二进制位。在另一个例子中，晶格速度集合包括在局限到晶格的空间中的动量状态的集合。

系统在晶格速度集合中模拟(1204)流体体积中的粒子移动，其中该移动造成粒子之间的碰撞。模拟过程在前面已描述。基于所模拟的移动，系统确定(1206)在体积内的特定位置处的粒子的相对粒子速度，其中相对粒子速度是(i)在体积内的特定位置处的并且在体积的零流动之下测出的粒子的绝对速度与(ii)在体积内的特定位置处的粒子中的一个或多个粒子的平均速度的差。在例子中，相对粒子速度是从在体积内的特定位置处的并在体积的零流动之下测出的粒子的绝对速度减去的在体积内的特定位置处的粒子中的一个或多个粒子的平均速度。在例子中，在体积内的特定位置处的粒子中的一个或多个粒子的平均速度包括在该特定位置处的特定类型粒子的平均速度。例如，流体体积可以包括各种不同类型的粒子。在这个例子中，系统被配置为确定特定类型的粒子的至少子集的平均速度。如前所述，系统基于方程(3.3)确定相对粒子速度。

系统还基于相对粒子速度确定(1208)表示碰撞的指定阶的非平衡碰撞后分布函数。在例子中，非平衡碰撞后分布函数(i)保留用于预定义物理量的非平衡矩，和(ii)消除用于未定义物理量的高至指定阶的非平衡矩。在这个例子中，通过在非平衡碰撞后分布函数的扩展形式中包括表示预定义物理量的项，非平衡碰撞后分布函数保留用于这些预定义物理量的非平衡矩。通过在扩展形式中不包括用于未定义的物理量的非平衡矩，例如，通过截断在方程3.7、3.8中所示的无穷级数并且保留与指定阶成比例的项，非平衡碰撞后分布函数消除用于未定义物理量的高至指定阶的非平衡矩。在这个例子中，指定阶是和流体速度与晶格声速的比相关联的指数值，其中晶格速度集合支持该指数值。

图13是网络环境1300的部件的框图。网络环境1300还包括系统1302，该系统1302包括存储器1304、总线系统1306和处理器1308。存储器1304可以包括硬盘驱动器和随机存取存储器存储设备，诸如动态随机存取存储器、机器可读硬件存储设备、机器可读介质或其它类型的非瞬时机器可读存储设备。包括例如数据总线和母板的总线系统1306可以被用来建立并控制系统1302的部件之间的数据通信。处理器1308可以包括一个或多个微处理器和/或处理设备。一般而言，处理器1308可以包括能够接收并存储数据并且能够经网络(未示出)通信的任何适当的处理器和/或逻辑。

系统1302可以是各种能够接收数据的计算设备当中的任一种，诸如服务器、分布式计算系统、台式计算机、膝上型计算机、手机、机架式安装的服务器等等。系统1302可以是在相同位置或在不同位置的单个服务器或一组服务器。所示出的系统1302可以经由输入/输出(“I/O”)接口1310接收数据。I/O接口1310可以是能够经由网络接收数据的任何类型的接口，诸如以太网接口、无线联网接口、光纤联网接口、调制解调器等等。系统1302被配置用于与数据存储库1312通信，该数据存储库1312可以被配置为存储速度模型、模拟数据等。

利用本文所述的技术，系统被描述为用于生成非平衡碰撞后分布函数，例如伽利略不变过滤后的算子。利用这些技术，生成各种类型的非平衡碰撞后分布函数，例如，如方程3.5、3.9、3.10和3.11中所示。所生成的非平衡碰撞后分布函数在建模流体体积中粒子的碰撞过程时被使用，如方程1.12中所示。

实施例可以在数字电子电路中或者在计算机硬件、固件、软件或者在它们的组合中实现。本文所述技术的装置可以在有形地体现或存储在机器可读介质(例如，硬件存储设备)中、用于由可编程处理器执行的计算机程序产品中实现；并且方法动作可以由执行指令程序以便通过对输入数据执行操作并生成输出来执行本文所述技术的操作的可编程处理器执行。本文所述的技术可以在一个或多个可在包括至少一个可编程处理器的可编程系统上执行的计算机程序中实现，其中编程处理器被耦合成从数据存储系统、至少一个输入设备和至少一个输出设备接收数据和指令，并向其发送数据和指令。每个计算机程序可以以高级过程化或面向对象的编程语言或者(如果期望的话)以汇编或机器语言实现；并且在任何情况下，该语言都可以是已编译或已解释语言。

作为例子，合适的处理器既包括通用微处理器，又包括专用微处理器。一般而言，处理器将从只读存储器和/或随机存取存储器接收指令和数据。一般而言，计算机将包括用于存储数据文件的一个或多个海量存储设备；这种设备包括磁盘，诸如内部硬盘和可移除盘；磁-光盘；以及光盘。适于有形地体现计算机程序指令和数据的存储设备包括所有形式的非易失性存储器，作为例子包括半导体存储器设备，诸如EPROM、EEPROM和闪速存储器设备；磁盘，诸如内部硬盘和可移除盘；磁-光盘；以及CD-ROM盘。前述中的任一种都可以用ASIC(专用集成电路)补充或结合在其中。

已经描述了许多实现方式。不过，应当理解的是，在不背离权利要求的精神和范围的情况下，可以进行各种修改。因此，其它实现方式在以下权利要求的范围之内。

Claims (48)

1.一种方法，包括：

在晶格速度集合中模拟流体体积中粒子的移动，其中所述移动引起所述粒子之间的碰撞；

基于所模拟的移动，

确定在所述体积内的特定位置处的粒子的相对粒子速度，其中所述相对粒子速度是(i)在所述体积内的所述特定位置处的并且在所述体积的零流动之下测出的所述粒子的绝对速度与(ii)在所述体积内的所述特定位置处的所述粒子中的一个或多个粒子的平均速度之差；及

基于所述相对粒子速度，确定表示所述碰撞的指定阶的非平衡碰撞后分布函数。

2.如权利要求1所述的方法，还包括：

由一个或多个计算机系统提供支持高至粒子速度的阶的流体动力移动的晶格速度集合；

其中模拟包括由所述一个或多个计算机系统进行模拟。

3.如权利要求2所述的方法，其中对晶格速度集合所支持的阶小于并且不同于非平衡碰撞后分布函数的所述指定阶；及

其中用于非平衡碰撞后分布函数的所述指定阶由所述粒子速度的阶来确定。

4.如权利要求1所述的方法，其中在所述体积内的所述特定位置处的粒子中的所述一个或多个粒子的平均速度包括在所述特定位置处的特定类型的粒子的平均速度。

5.如权利要求1所述的方法，其中晶格速度集合是与晶格玻尔兹曼方法关联的状态向量的集合。

6.如权利要求1所述的方法，其中非平衡碰撞后分布函数(i)保留用于预定义物理量的非平衡矩，以及(ii)消除用于未定义物理量的高至所述指定阶的非平衡矩。

7.如权利要求1所述的方法，其中所述指定阶是和流体速度与晶格声速的比相关联的指数值，其中晶格速度集合支持所述指数值。

8.如权利要求1所述的方法，其中晶格速度集合包括在被限定到晶格的空间中的动量状态的集合。

9.如权利要求1所述的方法，其中所述相对粒子速度是从所述体积内的特定位置处的并在所述体积的零流动之下测出的粒子的绝对速度中减去的在所述体积内的所述特定位置处的粒子中的所述一个或多个粒子的平均速度。

10.如权利要求1所述的方法，其中非平衡碰撞后分布函数是伽利略不变过滤后的算子。

11.如权利要求1所述的方法，还包括：

基于非平衡碰撞后分布函数，对所述流体体积中的粒子的碰撞过程进行建模。

12.如权利要求1所述的方法，其中，就用于为流体动力矩提供第一阶支持的晶格速度集合的马赫数而言，非平衡碰撞后分布函数是第一阶伽利略不变性的碰撞算子并且其中所述碰撞算子是根据下式定义的：

$C_i^{\left(i\right)}(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{w_i}{2T_0}\&lsqb;(1+\frac{c_i\&CenterDot;u(x,t)}{T_0})(\frac{c_i{c}_i}{T_0}-I)-\frac{c_{i}u(x,t)+u(x,t)c_i}{T_0}\&rsqb;:{\&Pi;}^{neq}(x,t)$

其中x是所述体积内的特定位置；

其中t是特定的时间点；

其中i是所述集合中的晶格速度的索引号；

其中T₀是常量晶格温度；

其中c_i是所述粒子在碰撞之前的速度向量；

其中u(x,t)是在时间t在特定位置x处的粒子当中的平均速度；

其中I是二秩统一张量；

其中τ是碰撞关系时间；

其中w_i是常量加权因子；及

其中Π^neq是非平衡动量通量。

13.如权利要求1所述的方法，其中非平衡碰撞后分布函数是用于为流体动力矩提供无限阶支持的晶格速度集合的碰撞算子C_i(x，t)，并且其中所述碰撞算子是根据下式定义的：

$C_i(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{f_i^{eq}(x,t)}{2\mathrm{\&rho;}(x,t)T_0}\&lsqb;\frac{c_i^{\&prime;}(x,t)c_i^{\&prime;}(x,t)}{T_0}-I\&rsqb;:{\&Pi;}^{neq}(x,t)$

其中x是所述体积内的特定位置；

其中t是特定的时间点；

其中i是所述集合中的晶格速度的索引号；

其中T₀是常量晶格温度；

其中I是二秩统一张量；

其中τ是碰撞关系时间；

其中c′_i(x，t)是相对粒子速度；

其中ρ是流体密度；

其中是平衡分布函数；及

其中Π^neq是非平衡动量通量。

14.如权利要求1所述的方法，其中，就用于对流体动力矩提供第二阶支持的晶格速度集合的马赫数而言，非平衡碰撞后分布函数是第二阶伽利略不变性的碰撞算子并且其中所述碰撞算子是根据下式定义的：

$\begin{array}{c}{C}_i^{\left(2\right)}(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{w_i}{2T_0}\&lsqb;(1+\frac{c_i\&CenterDot;u(x,t)}{T_0}+\frac{{(c_i\&CenterDot;u(x,t\left)\right)}^2}{2T_0^2}-\&CenterDot;\frac{u^2(x,t)}{2T_0})(\frac{c_i{c}_i}{T_0}-I)\\ -(1+\frac{c_i\&CenterDot;u(x,t)}{T_0})\frac{c_{i}u(x,t)+u(x,t)c_i}{T_0}+u(x,t)u(x,t)\&rsqb;:{\&Pi;}^{nep}(x,t)\end{array};$

其中x是所述体积内的特定位置；

其中t是特定的时间点；

其中i是所述集合中的晶格速度的索引号；

其中T₀是常量晶格温度；

其中c_i是所述粒子在碰撞之前的速度向量；

其中u(x,t)是在时间t在特定位置x处的粒子当中的平均速度；

其中I是二秩统一张量；

其中τ是碰撞关系时间；

其中w_i是常量加权因子；及

其中Π^neq是非平衡动量通量。

15.如权利要求1所述的方法，其中所述预定义物理量包括所述特定体积中的流体的质量、所述特定体积中的流体的动量以及所述特定体积中的流体的能量。

16.如权利要求1所述的方法，其中非平衡碰撞后分布函数是关于能量通量的碰撞算子C_i(x，t)，并且其中所述碰撞算子是根据下式定义的：

其中x是所述体积内的特定位置；

其中t是特定的时间点；

其中i是所述集合中的晶格速度的索引号；

其中T₀是常量晶格温度；

其中I是二秩统一张量；

其中τ是碰撞关系时间；

其中c′_i(x，t)是相对粒子速度；

其中是平衡分布函数；及

其中W^neq是非平衡能量通量。

17.一种或多种机器可读的硬件存储设备，所述硬件存储设备存储能够由一个或多个处理设备执行的指令以执行操作，所述操作包括：

在晶格速度集合中模拟流体体积中的粒子的移动，其中所述移动引起粒子之间的碰撞；

基于所模拟的移动，

确定在所述体积内的特定位置处的粒子的相对粒子速度，其中所述相对粒子速度是(i)在所述体积内的所述特定位置处的并且在所述体积的零流动之下测出的粒子的绝对速度与(ii)在所述体积内的所述特定位置处的粒子中的一个或多个粒子的平均速度之差；及

基于所述相对粒子速度，确定表示所述碰撞的指定阶的非平衡碰撞后分布函数。

18.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中所述操作还包括：

提供支持高至粒子速度的阶的流体动力移动的晶格速度集合。

19.如权利要求18所述的一种或多种机器可读的硬件存储设备，其中对晶格速度集合所支持的阶小于并且不同于非平衡碰撞后分布函数的所述指定阶；及

其中用于非平衡碰撞后分布函数的所述指定阶由所述粒子速度的阶确定。

20.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中在所述体积内的所述特定位置处的粒子中的所述一个或多个粒子的平均速度包括在所述特定位置处的特定类型粒子的平均速度。

21.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中晶格速度集合是与晶格玻尔兹曼方法关联的状态向量的集合。

22.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中非平衡碰撞后分布函数(i)保留用于预定义物理量的非平衡矩，以及(ii)消除用于未定义物理量的高至所述指定阶的非平衡矩。

23.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中所述指定阶是和流体速度与晶格声速的比相关联的指数值，其中晶格速度集合支持所述指数值。

24.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中晶格速度集合包括在被限定到晶格的空间中的动量状态的集合。

25.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中所述相对粒子速度是从所述体积内的特定位置处的并在所述体积的零流动之下测出的粒子的绝对速度中减去的在所述体积内的所述特定位置处的粒子中的所述一个或多个粒子的平均速度。

26.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中非平衡碰撞后分布函数是伽利略不变过滤后的算子。

27.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中所述操作还包括：

基于非平衡碰撞后分布函数，对所述流体体积中的粒子的碰撞过程进行建模。

28.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中，就用于为流体动力矩提供第一阶支持的晶格速度集合的马赫数而言，非平衡碰撞后分布函数是第一阶伽利略不变性的碰撞算子并且其中所述碰撞算子是根据下式定义的：

$C_i^{\left(1\right)}(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{w_i}{2T_0}\&lsqb;(1+\frac{c_i\&CenterDot;u(x,t)}{T_0})(\frac{c_i{c}_i}{T_0}-1)-\frac{c_{i}u(x,t)+u(x,t)c_i}{T_0}\&rsqb;:{\&Pi;}^{neq}(x,t)$

其中x是所述体积内的特定位置；

其中t是特定的时间点；

其中i是所述集合中的晶格速度的索引号；

其中T₀是常量晶格温度；

其中c_i是所述粒子在碰撞之前的速度向量；

其中u(x,t)是在时间t在特定位置x处的粒子当中的平均速度；

其中I是二秩统一张量；

其中τ是碰撞关系时间；

其中w_i是常量加权因子；及

其中Π^neq是非平衡动量通量。

29.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中非平衡碰撞后分布函数是用于为流体动力矩提供无限阶支持的晶格速度集合的碰撞算子C_i(x，t)，并且其中所述碰撞算子是根据下式定义的：

$C_i(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{f_i^{eq}(x,t)}{2\mathrm{\&rho;}(x,t)T_0}\&lsqb;\frac{c_i^{\&prime;}(x,t)c_i^{\&prime;}(x,t)}{T_0}-I\&rsqb;:{\&Pi;}^{neq}(x,t)$

其中x是所述体积内的特定位置；

其中t是特定的时间点；

其中i是所述集合中的晶格速度的索引号；

其中T₀是常量晶格温度；

其中I是二秩统一张量；

其中τ是碰撞关系时间；

其中c′_i(x，t)是相对粒子速度；

其中ρ是流体密度；

其中是平衡分布函数；及

其中Π^neq是非平衡动量通量。

30.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中，就用于对流体动力矩提供第二阶支持的晶格速度集合的马赫数而言，非平衡碰撞后分布函数是第二阶伽利略不变性的碰撞算子并且其中所述碰撞算子是根据下式定义的：

$\begin{array}{c}{C}_i^{\left(2\right)}(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{w_i}{2T_0}\&lsqb;(1+\frac{c_i\&CenterDot;u(x,t)}{T_0}+\frac{{(c_i\&CenterDot;u(x,t\left)\right)}^2}{2T_0^2}-\frac{u^2(x,t)}{2T_0})(\frac{c_i{c}_i}{T_0}-I)\\ -(1+\frac{c_i\&CenterDot;u(x,t)}{T_0})\frac{c_{i}u(x,t)+u(x,t)c_i}{T_0}+u(x,t)u(x,t)\&rsqb;:{\&Pi;}^{nep}(x,t)\end{array};$

其中x是所述体积内的特定位置；

其中t是特定的时间点；

其中i是所述集合中的晶格速度的索引号；

其中T₀是常量晶格温度；

其中c_i是所述粒子在碰撞之前的速度向量；

其中u(x,t)是在时间t在特定位置x处的粒子当中的平均速度；

其中I是二秩统一张量；

其中τ是碰撞关系时间；

其中w_i是常量加权因子；及

其中Π^neq是非平衡动量通量。

31.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中所述预定义物理量包括所述特定体积中的流体的质量、所述特定体积中的流体的动量以及所述特定体积中的流体的能量。

32.如权利要求17所述的一种或多种机器可读的硬件存储设备，其中非平衡碰撞后分布函数是关于能量通量的碰撞算子C_i(x，t)，并且其中所述碰撞算子是根据下式定义的：

其中x是所述体积内的特定位置；

其中t是特定的时间点；

其中i是所述集合中的晶格速度的索引号；

其中T₀是常量晶格温度；

其中I是二秩统一张量；

其中τ是碰撞关系时间；

其中c′_i(x，t)是相对粒子速度；

其中是平衡分布函数；及

其中W^neq是非平衡能量通量。

33.一种系统，包括：

一个或多个处理设备；及

一个或多个机器可读的硬件存储设备，所述硬件存储设备存储能够由一个或多个处理设备执行的指令以执行操作，所述操作包括：

在晶格速度集合中模拟流体体积中的粒子的移动，其中所述移动引起粒子之间的碰撞；

基于所模拟的移动，

确定在所述体积内的特定位置处的粒子的相对粒子速度，所述相对粒子速度是(i)在所述体积内的所述特定位置处的并且在所述体积的零流动之下测出的粒子的绝对速度与(ii)在所述体积内的所述特定位置处的粒子中的一个或多个粒子的平均速度之差；及

基于所述相对粒子速度，确定表示所述碰撞的指定阶的非平衡碰撞后分布函数。

34.如权利要求33所述的系统，其中所述操作还包括：

提供支持高至粒子速度的阶的流体动力移动的晶格速度集合。

35.如权利要求34所述的系统，其中对晶格速度集合所支持的阶小于并且不同于非平衡碰撞后分布函数的所述指定阶；及

其中用于非平衡碰撞后分布函数的所述指定阶由所述粒子速度的阶确定。

36.如权利要求33所述的系统，其中在所述体积内的所述特定位置处的粒子中的所述一个或多个粒子的平均速度包括在所述特定位置处的特定类型粒子的平均速度。

37.如权利要求33所述的系统，其中晶格速度集合是与晶格玻尔兹曼方法关联的状态向量的集合。

38.如权利要求33所述的系统，其中非平衡碰撞后分布函数(i)保留用于预定义物理量的非平衡矩，以及(ii)消除用于未定义物理量的高至所述指定阶的非平衡矩。

39.如权利要求33所述的系统，其中所述指定阶是和流体速度与晶格声速的比相关联的指数值，其中晶格速度集合支持所述指数值。

40.如权利要求33所述的系统，其中晶格速度集合包括在被限定到晶格的空间中的动量状态的集合。

41.如权利要求33所述的系统，其中所述相对粒子速度是从所述体积内的特定位置处的并在所述体积的零流动之下测出的粒子的绝对速度中减去的在所述体积内的所述特定位置处的粒子中的所述一个或多个粒子的平均速度。

42.如权利要求33所述的系统，其中非平衡碰撞后分布函数是伽利略不变过滤后的算子。

43.如权利要求33所述的系统，其中所述操作还包括：

基于非平衡碰撞后分布函数，对所述流体体积中的粒子的碰撞过程进行建模。

44.如权利要求33所述的系统，其中，就用于为流体动力矩提供第一阶支持的晶格速度集合的马赫数而言，非平衡碰撞后分布函数是第一阶伽利略不变性的碰撞算子并且其中所述碰撞算子是根据下式定义的：

$C_i^{\left(1\right)}(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{w_i}{2T_0}\&lsqb;(1+\frac{c_i\&CenterDot;u(x,t)}{T_0})(\frac{c_i{c}_i}{T_0}-1)-\frac{c_{i}u(x,t)+u(x,t)c_i}{T_0}\&rsqb;:{\&Pi;}^{neq}(x,t)$

其中x是所述体积内的特定位置；

其中t是特定的时间点；

其中i是所述集合中的晶格速度的索引号；

其中T₀是常量晶格温度；

其中c_i是所述粒子在碰撞之前的速度向量；

其中u(x,t)是在时间t在特定位置x处的粒子当中的平均速度；

其中I是二秩统一张量；

其中τ是碰撞关系时间；

其中w_i是常量加权因子；及

其中Π^neq是非平衡动量通量。

45.如权利要求33所述的系统，其中非平衡碰撞后分布函数是用于为流体动力矩提供无限阶支持的晶格速度集合的碰撞算子C_i(x，t)，并且其中所述碰撞算子是根据下式定义的：

$C_i(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{f_i^{eq}(x,t)}{2\mathrm{\&rho;}(x,t)T_0}\&lsqb;\frac{c_i^{\&prime;}(x,t)c_i^{\&prime;}(x,t)}{T_0}-I\&rsqb;:{\&Pi;}^{neq}(x,t)$

其中x是所述体积内的特定位置；

其中t是特定的时间点；

其中i是所述集合中的晶格速度的索引号；

其中T₀是常量晶格温度；

其中I是二秩统一张量；

其中τ是碰撞关系时间；

其中c′_i(x，t)是相对粒子速度；

其中ρ是流体密度；

其中是平衡分布函数；及

其中Π^neq是非平衡动量通量。

46.如权利要求33所述的系统，其中，就用于对流体动力矩提供第二阶支持的晶格速度集合的马赫数而言，非平衡碰撞后分布函数是第二阶伽利略不变性的碰撞算子并且其中所述碰撞算子是根据下式定义的：

$\begin{array}{c}{C}_i^{\left(2\right)}(x,t)=(1-\frac{1}{\mathrm{\&tau;}})\frac{w_i}{2T_0}\&lsqb;(1+\frac{c_i\&CenterDot;u(x,t)}{T_0}+\frac{{(c_i\&CenterDot;u(x,t\left)\right)}^2}{2T_0^2}-\frac{u^2(x,t)}{2T_0})(\frac{c_i{c}_i}{T_0}-I)\\ -(1+\frac{c_i\&CenterDot;u(x,t)}{T_0})\frac{c_{i}u(x,t)+u(x,t)c_i}{T_0}+u(x,t)u(x,t)\&rsqb;:{\&Pi;}^{nep}(x,t)\end{array};$

其中x是所述体积内的特定位置；

其中t是特定的时间点；

其中i是所述集合中的晶格速度的索引号；

其中T₀是常量晶格温度；

其中c_i是所述粒子在碰撞之前的速度向量；

其中u(x,t)是在时间t在特定位置x处的粒子当中的平均速度；

其中I是二秩统一张量；

其中τ是碰撞关系时间；

其中w_i是常量加权因子；及

其中Π^neq是非平衡动量通量。

47.如权利要求33所述的系统，其中所述预定义物理量包括所述特定体积中的流体的质量、所述特定体积中的流体的动量以及所述特定体积中的流体的能量。

48.如权利要求33所述的系统，其中非平衡碰撞后分布函数是关于能量通量的碰撞算子C_i(x，t)，并且其中所述碰撞算子是根据下式定义的：

其中x是所述体积内的特定位置；

其中t是特定的时间点；

其中i是所述集合中的晶格速度的索引号；

其中T₀是常量晶格温度；

其中I是二秩统一张量；

其中τ是碰撞关系时间；

其中c′_i(x，t)是相对粒子速度；

其中是平衡分布函数；及

其中W^neq是非平衡能量通量。