一种基于涡度损失的自适应涡度限制力方法与流程

本发明涉及流体模拟技术领域，更具体的，涉及一种基于涡度损失的自适应涡度限制力方法。

背景技术：

流体模拟技术在电影、游戏等行业有着丰富的应用。电影与游戏中的烟雾、海洋等流体画面常常由于拍摄难度大、需要可控地流动或与游戏者交互等原因，无法以直接拍摄或提前建模的形式实现，必须通过计算机模拟的方式来得到物理真实的流体，亦即计算机图形学中的流体模拟技术。

在流体模拟研究领域中，往往需要均衡两个根本性的问题：更低的性能开销和更丰富的流体视觉效果。流体模拟技术由于其求解精度不足，存在一定的数值耗散，在视觉上表现为湍流细节的损失，而丰富的湍流细节对提升流体模拟的视觉效果、提高其视觉真实感起到重要的作用。为了获得更丰富的细节，最直接的方法是提高模拟场的分辨率，但这将带来指数级的运算量提升。因此，如何合理地补偿求解过程中的精度损失，以较小的计算量获得尽可能多的湍流细节，是流体模拟研究的一个重要问题。

基于欧拉网格求解流体控制方程是最常见的流体模拟方法，前人在此方法上已做了大量工作。比如通过semi-lagrangian法求解控制方程中速度的自对流，该方法是无条件稳定的，但数值耗散导致对流项的求解始终存在着精度不足的问题，在视觉上则表现为被弱化的流体细节。为了减弱对流项中数值耗散的影响，提出了一种涡旋限制力方法，该方法以较小的计算开销有效提升了模拟流体的细节丰富度，得到了广泛的应用。但限于精度等因素的约束，其发挥的效果依然是有限的。另一方面，由于该方法在计算涡旋限制力时采用了全局统一的细节增强系数，使其在模拟复杂涡流场时存在难以调和的矛盾：对于数值耗散严重的区域，必须采用足够大的系数才能弥补细节的损失，然而较大的系数又往往容易造成流体的不稳定和失真感。

考虑到流体的速度场在执行上述自对流步骤时忽视了流体的不可压缩性，其后必须执行投影步骤以消除，因此造成的速度散度非零。尽管这一方法保证了流体模拟的无条件稳定，但消除散度必然带来速度的损失。因此，如何补偿这一速度损失，减弱其造成的数值耗散，提升流体细节丰富度，具有重要研究意义。

技术实现要素：

本发明主要研究如何通过速度场的涡度来评估速度场在自对流过程中的精度损失的问题，提供了一种基于涡度损失的自适应涡度限制力方法，其能准确地弥补对流过程中的数值耗散。

为实现上述本发明目的，采用的技术方案如下：一种基于涡度损失的自适应涡度限制力方法，所述该方法包括以下步骤：

s1：基于欧拉网格法求解流体控制方程，解算每一帧中的流体控制方程；

s2：在每一帧计算对流步之前，对当前速度场uⁿ求旋度，得到其涡度场ωⁿ；通过对流求得其下一状态ωⁿ⁺¹，并记此涡度场为ω′

s3：对当前帧的速度场uⁿ进行对流，得到对流后的速度场uⁿ⁺¹，针对对流后的速度场uⁿ⁺¹求旋度，得到其涡度场ω^*；

s4：计算涡度损失因子，其计算公式如下：

其中，

式中，δω表示涡度损失因子，ω表示流体涡度；

s5：根据涡度损失因子构建自适应涡旋限制力模型，其模型表达式如下：

favc＝ε′·δω(n×ω)

其中，

式中，favc表示自适应涡旋限制力，ε′表示涡旋限制力的细节增强系数；n表示涡度幅值的梯度方向，ω表示流体涡度。

优先地，步骤s1，所述流体控制方程的表达式如下：

式中：u表示流体的速度，t表示流体的流过的时间，ρ表示流体的密度，p表示流体的压强，f表示流体受到的外力。

进一步地，所述步骤s2，当前速度场uⁿ求涡度场ωⁿ的计算公式如下：

式中，为旋度算子，为梯度算子。

再进一步地，步骤s2，根据涡度场ωⁿ，通过对流求其下一状态ωⁿ⁺¹，具体如下：

记对流后的涡度场ωⁿ⁺¹为ω′。

再进一步地，所述步骤s3，其中，对当前帧的速度场un进行对流，得到对流后的速度场uⁿ⁺¹，即：

uⁿ⁺¹＝advect(uⁿ)

根据对流后的速度场uⁿ⁺¹求旋度，得到其涡度场ω^*，其计算公式如下：

将得到的涡度场ωⁿ⁺¹记为ω^*。

再进一步地，步骤s6，将自适应涡旋限制力加入速度场中，并继续进行逐帧模拟，其计算公式如下：

uⁿ⁺¹＝uⁿ⁺¹+δtfavc

式中，δt表示相邻两帧之间的间隔时间；

更新模拟场，得到随时间变化的模拟流体。

本发明的有益效果如下：本发明通过引入涡度损失因子δω，该涡度损失因子的大小随速度损失大小动态地变化，通过将涡度损失因子与涡旋限制力结合，构建自适应涡旋限制力，其具有了动态弥补数值耗散影响的理想特性。

附图说明

图1是本实施例所述方法的步骤流程图。

图2是本实施例所述方法取得的效果图。

图3是本实施例现有技术的方法取得的效果图。

具体实施方式

下面结合附图和具体实施方式对本发明做详细描述。

实施例1

如图1所示，一种基于涡度损失的自适应涡度限制力方法，所述该方法包括以下步骤：

步骤s1：创建流体模型，初始化参数，进行解算每一帧中的流体控制方程：首先建造网格系统以便进行流体模拟，初始化参数包括且不限于流体的初始速度u、密度ρ、温度、所受外力f、涡旋限制力的细节增强系数ε′。根据各流体参数，逐帧求解流体控制方程，不断得到流体的下一帧状态并更新。每一帧执行的运算是相同的。本例中，设定仅受固定浮力作用垂直上升的烟雾，取其中一帧如图2由子图中所示。

其中，所述流体控制方程的表达式如下：

式中：u表示流体的速度，t表示相邻两帧流体间隔的时间，ρ表示流体的密度，p表示流体的压强，f表示流体受到的外力。

本实施例所述的流体控制方程即ns方程描述了流体运动的基本规律。给定一系列初始的流体状态(速度、密度、温度、外力等)，代入ns方程，可以求解得到经过固定时间间隔δt后的下一流体状态。不断重复这一过程，就可以得到随时间变化的一个流体序列，也就是流体随时间运动的过程。

而上述流体的状态需要以一定形式储存，所谓欧拉网格法就是把模拟空间均匀切分成大量细小的网格单元。每一个网格单元储存着该位置的速度、涡度等物理信息。

步骤s2：同步对流速度场和涡度场。首先，在每一帧计算对流步前，对当前速度场uⁿ求旋度，得到其涡度场ωⁿ。其中，涡度场ωⁿ的计算公式如下：

式中，为旋度算子，为梯度算子。

然后，不同于传统的单一速度场对流方法，本实施例同步对流当前速度场uⁿ和涡度场ωⁿ，通过对流求得其下一状态ωⁿ⁺¹，记对流后的涡度场为ω′。

具体如下：

同步对流此速度场和涡度场，记对流后的涡度场ωⁿ⁺¹为ω′。

步骤s3：接着，对当前帧的速度场uⁿ进行对流，得到对流后的速度场uⁿ⁺¹，

即：

uⁿ⁺¹＝advect(uⁿ)

针对对流后的速度场uⁿ⁺¹求旋度，得到其涡度场ω^*，其计算公式如下：

考虑到涡度场始终是无散的，因而在执行对流步后，流体的不可压性仅造成了速度场损失，而不会对涡度场造成损失，因此，可认为涡度场ω′的精度优于ω^*

步骤s4：计算涡度损失因子δω，其计算公式如下：

式中，δω表示涡度损失因子。

理想情况下，涡度场ω^*与ω′应是完全相等的。但由于求解速度自对流的过程中，显然没有考虑到压力因素的影响，即忽视了流体的不可压性，故上述两种方法解得的涡度场存在差异。这一差异实际上就是在速度自对流时产生的精度误差，且当模拟的时间步长越大时，其差异越明显。本实施例通过计算涡度损失因子δω作为衡量精度此误差的标准。

步骤s5：使用自适应涡旋限制力弥补精度损失。涡度损失因子δω提供了速度场在对流过程中的精度损失情况，因而可以作为精度损失的衡量指标，指导各种精度补偿方法更准确、自适应地弥补精度损失。

根据涡度损失因子构建自适应涡旋限制力模型，其模型表达式如下：

favc＝ε′·δω(n×ω)

其中，

式中，favc表示自适应涡旋限制力，ε′表示涡旋限制力的细节增强系数；n表示涡度幅值的梯度方向。，ω表示流体涡度。

根据上式得到自适应涡旋限制力favc，通过将上述的损失因子δω引入到原始的涡旋限制力方法中。不同于常量细节增强因子ε，对流涡度损失δω在整个模拟场中并非固定值，而是随当地的物理特性而存在大小各异的损失。因而将δω与原始的涡旋限制力方法结合，就可以得到具有自适应性的涡旋限制力。

根据损失量的大小，δω越大的区域，补偿的自适应涡旋限制力favc则越大，反之favc越小。本实施例所述的方法具有了对精度损失的自适应性，能更合理地对流体细节进行补偿。本实施例计算涡度损失因子δω的方法独立于求解流体控制方程的各个步骤，因此可以与现有的各种求解方法良好地结合，在发挥其自适应细节补偿的同时不会对求解过程造成影响。

步骤s6：将自适应涡旋限制力加入速度场中，并继续进行逐帧模拟。将自适应涡旋限制力加入速度场的计算公式如下：

uⁿ⁺¹＝uⁿ⁺¹+δtfavc

更新模拟场，得到随时间变化的模拟流体。

本实施例引入自适应因子δω使得涡旋限制力可以根据每个网格单元对流前后的速度损失情况进行动态变化，其能更准确地弥补对流过程中的数值耗散。

本实施例所述的自适应涡旋限制力方法相比原始方法，如图2所示，图3中表示现有技术的方法。从图2、图3中可以得到，本实施例所述的方法更准确地弥补了精度损失，因而取得了更丰富的视觉细节。从整体来看，本实施例所述的自适应涡旋限制力方法产生了更多的涡旋细节；从局部来看，采用本实施例所述的自适应涡旋限制力方法模拟的烟柱顶部保留了更丰富的涡旋层次，而现有技术的方法的相同区域受数值耗散影响较严重。这是因为现有技术的涡旋限制力方法不能自适应地进行补偿，因而视觉细节相对模糊。

显然，本发明的上述实施例仅仅是为清楚地说明本发明所作的举例，而并非是对本发明的实施方式的限定。凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明权利要求的保护范围之内。