CN107341849B - 一种快速实时烟雾模拟算法 - Google Patents

Abstract

本发明公开了一种快速实时烟雾模拟算法，根据粒子位置，在整个计算空间创建八叉树；从下至上遍历八叉树，计算八叉树中每一节点的总质量ρ和中心质量ρ_c，每一个粒子从上至下遍历八叉树来计算它的受力，更新粒子的速度和位移，获得实时的烟雾效果。本发明通过自适应八叉树算法，成功地将烟雾的涡粒子计算复杂度从O(N²)降至O(N log N)量级，有效地减少了计算量，加快了求解速度。本发明在10⁴量级的粒子数情况下，计算速度比单纯的直接求和计算大约提高了1.85倍。在更多的粒子数的情况，计算速度会比直接求和更快。

Description

一种快速实时烟雾模拟算法

技术领域

本发明属于实时流体模拟技术领域，涉及一种流体模拟的动态建模算法。

背景技术

烟雾模拟在游戏和影视特效等方面应用非常广泛，形象而逼真的烟雾能给用户带来良好的视觉体验。烟雾现象由于其外表并非光滑的表面，十分不规则，而且随着时间的改变而不断地变化，这使得烟雾的模拟具有一定的挑战性。同时，在使用基于物理的方法来模拟烟雾流体时，巨大的计算量往往成为提高模拟速度的瓶颈，因此如何快速有效地进行模拟具有十分重要的意义。

在早期的时候，受到计算机软硬件条件的限制，人们使用简单的数学模型来对各种流体进行模拟，由于这些模型并不能真正的表现出流体的运动规律，只能说是形似而神不似。随后，人们又开始利用粒子系统、纹理映射等方法来实现流体的模拟，但是随着影视特效、游戏效果对真实感的要求的提高以及科学研究的需要，这些传统的方法越来越难以胜任。此时基于物理的模拟由于其真实地反映了流体的运动规律，在模拟效果更加逼真的同时又可以表现出流体运动过程中丰富的细节，从而逐渐成为流体模拟研究的热点。近年来，随着虚拟现实、电脑游戏等行业的兴起，当前的流体模拟在实时性方面已经难以适应它们的要求，如何提高流体模拟的实时性，已成为一个亟需解决的问题。

涡方法作为一类重要的数值计算方法，因其计算中不需要计算网格且适合于数值分析大雷诺数流场，在计算过程中数值扩展很小等优点而受到人们的重视。涡方法以涡量动力学方程为基础，将速度的涡量场离散成涡元，然后把涡元当作粒子，利用拉格朗日法求解涡量动力学方程。涡方法是一种在研究大雷诺数流动中极为有效的数值方法。它在展现湍流丰富的细节上具有出色的表现，而且没有其他方法所带有的数值扩散或是投影有关的耗散等缺陷，同时由于涡旋在视觉上的重要性，许多图形学研究者纷纷采用涡方法进行流体模拟。

在涡方法求解速度的过程中，将其抽象成一个个涡粒子，这些涡粒子之间不断地相互作用，需要求解这些涡粒子之间的相互作用力。由此便可归结为N体问题的求解。N体问题模拟在数值上近似地表示一个多体系统的演化过程，该系统中的每一个体都持续地与所有其他的体相互作用。一个相似的例子是天体物理学仿真，在该仿真中，每个体代表一个星系或者一个独立的行星，各个体之间通过万有引力相互吸引。在很多其他的计算机科学问题中也会出现N体仿真，例如，蛋白质折叠研究就使用N体仿真计算静电和范德华力。其他使用N体仿真的例子还有湍流流场仿真与全局光照计算等计算机图形学中的问题。

N体问题可以准确地描述如下：在一定的物理空间中，分布有一定数量的粒子，每对粒子间都存在相互作用(如万有引力、库仑力等)。他们初始的速度和位置是确定的，每隔一定的时间步，由于粒子间的作用，需要更新粒子的速度和位置。更新给定的粒子的速度和位移的计算很简单，就是把另外N-1个粒子对它作用的结果叠加起来，但这导致了每一时间步O(N²)的计算量。

可见，在涡粒子方法中应用N体问题时，具有两个重要的特征：第一是计算规模大，因为包含了大量的的涡粒子，粒子的规模大到数百万、千万。由于在系统中任意的两个粒子间都存在着作用，因此直接就按粒子间的相互作用的量级就是O(N²)。第二，系统是动态变化的。为了反映系统的实时变化，要求时间步长需要足够小。这两个特征决定了计算机模拟时巨大的计算量，这对于任何高性能的单台计算机来说都是一个很难突破的瓶颈。

以上所提到的N体仿真方法是一种直接方法，它计算N个体之间所有成对的相互作用。这是一种相对简单的方法，但由于这种方法具有O(N²)的计算复杂度，因此在大系统的仿真中通常不能单独使用它。由此，便衍生出了一系列在大系统中的N体求和算法，这种形式的快速N体算法包括Barnes-Hutch方法(BH)^[1]、快速多级方法(FMM)^[2][3]，以及粒子网格(PM)方法^[4][5]。

本发明的参考文献如下：

[1]Barnes,J.,and P.Hut.1986.“A Hierarchical O(n log n)ForceCalculation Algorithm.”Nature 324。

[2]Greengard,L.1987.The Rapid Evaluation of Potential Fields inParticle Systems.ACM Press。

[3]Greengard,Leslie F.,and Jingfang Huang.2002.“ANew Version oftheFast Multipole Method for Screened Coulomb Interactions in Three Dimensions.”Journal ofComputational Physics 180(2),pp.642–658。

[4]Hockney,R.,and J.Eastwood.1981.Computer Simulation UsingParticles.McGraw-Hill。

[5]Darden,T.,D.York,and L.Pederson.1993.“Particle Mesh Ewald:An N log(N)Method for Ewald Sums in Large Systems.”Journal of Chemical Physics98(12),p.10089。

发明内容

为解决现有技术存在的上述问题，本发明要设计一种计算量少且计算速度快的快速实时烟雾模拟算法。

为了实现上述目的，本发明技术方案如下：一种快速实时烟雾模拟算法，包括以下步骤：

A、设共有N个粒子，已知第j个粒子的质量ρ_j、位移x_j和速度v_j；

B、根据粒子位置，在整个计算空间创建八叉树；所述的八叉树为自适应的八叉树；

C、从下至上遍历八叉树，计算八叉树中每一节点的总质量ρ和中心质量ρ_c，公式如下：

其中N为粒子个数，ρ_j为第j个粒子质量，v_j为第j个粒子速度，x为当前位置，x_j为第j个粒子的位置；公式中x≠x_j表示x与x_j不能相等，||x-x_j||表示x与x_j的差恒为正。

D、每一个粒子从上至下遍历八叉树来计算它的受力，假定八叉树三维空间各边长为d，若粒子距离某一节点大于该边长的一半d/2，则利用总质量和质量中心进行计算，否则遍历它的子节点；此时的速度通过以下公式来计算：

其中，v_j为第j个粒子的速度，ρ_c为中心质量，x为当前位置，x_j为第j个粒子的位置，r为粒子半径，ω(x)为粒子所携带的涡量；

E、更新粒子的速度和位移，转到步骤B继续迭代，获得实时的烟雾效果。

本发明基于以下理论：如果粒子相互间作用力的大小随着距离的增大下降很快时，当一组粒子距离被作用的粒子足够远时，它们的作用可以近似地用一个等价粒子来代替。

本发明用树形数据结构来代表对整个物理空间划分形成的粒子集，在给定的精度范围内来近似粒子和粒子集间的相互作用。在本发明中，整个物理空间被递归地分成不同层次的单元格，在三维空间就是八叉树，二维空间中是四叉树。在最细级别的单元(分级树的叶子)包含了粒子，较高层次的单元(节点)代表了以它们为根的子树中粒子的集合。创建过程大致为：首先从根节点开始，根节点代表了整个物理空间，包含了所有的粒子，将物理空间递归的分成8个(三维)子节点和4个(二维)子节点来形成树的每一层次，直到最细的单元为止。

与现有技术比，本发明的有益效果如下：

1、本发明通过自适应八叉树算法，成功地将烟雾的涡粒子计算复杂度从O(N²)降至O(N log N)量级，有效地减少了计算量，加快了求解速度。

2、本发明在10⁴量级的粒子数情况下，计算速度比单纯的直接求和计算大约提高了1.85倍。在更多的粒子数的情况，计算速度会比直接求和更快。

附图说明

图1是本发明所使用的自适应八叉树二维结构图。

图2是本发明所使用的自适应八叉树树型结构图。

图3是本发明的求和计算示意图中的根据树节点自上而下求和。

图4是本发明的求和计算示意图中的求和区域(深色区域为近程区域，可以采用直接求和的方式，空白区域为远程区域，求和是去其总质量和质量中心)。

具体实施方式

本发明中，远程一组粒子的作用通常用它们的总质量和中心质量来近似，并且在八叉树的最细层次上每个单元格最多有一个粒子。在计算粒子的受力时，粒子从八叉树的根节点开始遍历八叉树，如果节点和粒子间的距离大于空间边长的一半d/2的判据，则利用节点的总质量来计算对粒子的作用，否则就需要细分到更细的层次来进一步比较。这样，可以将粒子间作用的计算时间减少到O(N log N)量级。计算粒子间受力的精确度可以通过采用高阶的近似，而不是一阶的质量中心来代替。

一种快速实时烟雾模拟算法，包括以下步骤：

A、初始化粒子的质量、位移和速度。如果节点不包含涡粒子，则把新的涡粒子放入其中。如果节点是一个内部节点，就更新总质量和质心。递归地将物体插入到图1中分叉中的某个。

B、根据粒子位置，在整个计算空间创建八叉树。所述的八叉树为自适应的八叉树，如图2所示，如果节点是一个外部节点，它包含有一个涡粒子，那么也就是说在一个相同的区域里包含有两个涡粒子和。那么就把这个区域进一步划分成四个子区域。然后递归地插入涡粒子的和到对应的分叉中。因为涡粒子的和可能仍然位于相同的子区域中，因此在一个单独的插入操作中可能涉及多个子区域划分。最终更新节点的质心和总质量。

C、从下至上遍历八叉树，如图3所示，计算八叉树中每一节点的总质量和质量中心(或更高阶近似)。

D、每一个涡粒子从上至下遍历八叉树来计算它的受力，如图4所示，假定八叉树三维空间各边长为d，若粒子距离某一节点大于该边长的一半d/2，则为远程区域，在图中为空白区域，利用总质量和质量中心(或更高阶近似)进行计算，否则，则为近程区域，在图中为深色区域，需要遍历它的子节点，并采用直接求和的方式来计算受力。如果当前节点是一个外部节点，计算当前节点施加在物体上的力，并将其加到合力上。否则，计算粒子之间的距离x-x_j与边长d的商(x-x_j)/d的值。

E、如果将这个内部节点看成一个单独的物体，则计算其施加在物体上的力，并将其加到的合力上，否则在当前节点的每个孩子节点上递归地执行上述步骤。该过程结束后，转到步骤B继续进行迭代更新至下一帧，最终获得烟雾的实时效果。

本发明不局限于本实施例，任何在本发明披露的技术范围内的等同构思或者改变，均列为本发明的保护范围。

Claims (1)

1.一种快速实时烟雾模拟算法，其特征在于：包括以下步骤：

A、设共有N个粒子，已知第j个粒子的质量ρ_j、位移x_j和速度v_j；

B、根据粒子位置，在整个计算空间创建八叉树；所述的八叉树为自适应的八叉树；

C、从下至上遍历八叉树，计算八叉树中每一节点的总质量ρ和中心质量ρ_c，公式如下：

其中N为粒子个数，ρ_j为第j个粒子质量，v_j为第j个粒子速度，x为当前位置，x_j为第j个粒子的位置；公式中x≠x_j表示x与x_j不能相等，||x-x_j||表示x与x_j的差恒为正；

D、每一个粒子从上至下遍历八叉树来计算它的受力，假定八叉树三维空间各边长为d，若粒子距离某一节点大于该边长的一半d/2，则利用总质量和质量中心进行计算，否则遍历它的子节点；此时的速度通过以下公式来计算：

其中，v_j为第j个粒子的速度，ρ_c为中心质量，x为当前位置，x_j为第j个粒子的位置，r为粒子半径，ω(x)为粒子所携带的涡量；

E、更新粒子的速度和位移，转到步骤B继续迭代，获得实时的烟雾效果。

Images