一、说明

   关于流体物质的仿真和模拟，需要流体理论方面的一般知识。我们这里从基本流体方程入手，详细解释如何实现流体仿真的每一个具体步骤。

二、摘要

   构建类似流体运动的动画工具是一项重要且重要的工作。计算机图形学中的许多应用都具有挑战性的问题。
  使用基于物理的流体流动模型可以极大地帮助在创建此类工具时。与关键帧或基于程序的技术不同，物理模型允许动画师几乎毫不费力地创造有趣的、类似旋转流体的行为。此外，流与物体和虚拟力的相互作用也得到了优雅的处理。
  直到最近，人们还认为物理流体模型太实现实时交互的成本很高。这很大程度上是由于事实上，以前的模型使用不稳定的方案来求解控制流体的物理方程。在这篇论文中，第一次，我们提出了一个无条件稳定的模型，它仍然产生复杂的流体状流动。同样，我们的方法非常容易实现。我们模型的稳定性使我们能够花更多的时间步骤，从而实现更快的模拟。我们已经用我们的模型与平流实体纹理相结合，创造出许多二维和三维交互式流体动画。
  CR 类别：I.3.7 [计算机图形学]：三维图形与写实主义——动画

三、简介

   计算机图形学中最有趣的问题之一是模拟流体行为。一个好的流体求解器是伟大的在许多不同领域的重要性。在特效行业对令人信服地模仿外观有很高的要求烟、水和火等流体的行为。绘画程序还可以受益于流体求解器来模拟传统技术例如水彩和油画。纹理合成是另一种可能的应用。事实上，许多纹理是由类似流体的过程产生的，例如侵蚀。流体的建模和模拟是，当然，在大多数科学学科中也至关重要和工程领域。流体力学被用作这些模拟所基于的标准数学框架。那里科学家们一致认为，纳维-斯托克斯方程是一个非常好的流体流动模型。数以千计的书籍和关于如何计算的文章已在各个领域发表这些方程的数值。在实践中使用哪种求解器很大程度上取决于手头的问题和计算能力可用的。大多数工程任务要求模拟提供所涉及的物理量的准确范围以回答与安全、性能等相关的问题。 视觉外观在这些应用中，流动（形状）是次要的。
   另一方面，在计算机图形学中，流体的形状和行为是主要关注点，而物理准确性则是最重要的。次要的或在某些情况下无关的。流体求解器，用于计算机图形，理想情况下应该为用户提供一个工具，使他们能够实时实现类似流体的效果。这些因素比严格的物理准确性更重要，后者需要太多的时间计算能力。
   事实上，以前计算机图形学中的大多数模型都是由通过视觉外观而不是物理准确性。早期流动模型是根据简单的原语构建的。各种组合这些基元允许粒子系统的动画 [15, 17]或简单的几何形状，例如叶子[23]。的复杂性随着随机湍流的引入，流动得到了极大的改善 [16, 20]。这些湍流是质量守恒的，并且，因此，自动表现出旋转运动。而且湍流在空间和时间上都是周期性的，这对于运动来说是理想的“纹理映射”[19]。由叠加建立的流流原语都有一个缺点，就是不响应动态地适应用户施加的外力。动力学模型基于纳维-斯托克斯方程的流体分析首先在二维中实现。耶格尔、厄普森和加米托等人。使用涡旋方法与泊松解算器相结合来创建流体的二维动画 [24, 8]。后来，陈等人等人。根据纳维-斯托克斯方程 [2] 的二维模拟给出的压力项生成动画水面。他们的与我们不同的方法不仅限于二维而且不稳定。 Kass 和 Miller 将浅水方程线性化为模拟液体[12]。然而，这种简化并没有捕捉到流体有趣的旋转运动特征。更多的最近，Foster 和 Metaxas 清楚地展示了使用完整三维纳维-斯托克斯方程创建方程的优势类似流体的动画[7]。许多难以设置关键帧的效果手动获取诸如旋转运动和流过物体的信息是自动获得的。他们的算法主要基于工作哈洛和韦尔奇在计算流体动力学方面的研究，可追溯到回到1965年[11]。从那时起，福斯特和梅塔克萨斯可以使用的许多其他技术已经被开发出来。然而，他们的模型的优点是编码简单，因为它是基于纳维-斯托克斯方程的有限差分和显式时间求解器。类似的求解器及其源代码是也可从 Griebel 等人的书中找到。 [9]。显式求解器的主要问题是数值格式可以变成对于大时间步长来说不稳定。不稳定会导致数值模拟“爆炸”，因此必须重新启动较小的时间步长。这些显式算法集的不稳定性速度和交互性受到严重限制。理想情况下，用户应该是能够与流体求解器实时交互，而无需担心可能出现的“爆炸”。
   在本文中，我们首次提出了一种稳定的算法求解完整的纳维-斯托克斯方程。我们的算法非常容易实现，并允许用户与实时交互图形工作站上的三维流体。我们实现这是通过使用比 Foster 和 Metaxas 使用的时间步长大得多的时间步长来实现的。为了获得稳定的求解器，我们离开了 Foster以及 Metaxas 的解决方法。而不是他们明确的欧拉方案中，我们使用拉格朗日法和隐式方法来求解纳维-斯托克斯方程。我们的方法在计算流体文献中找不到，因为它是为计算机定制的图形应用程序。该模型不够准确对于大多数工程应用。确实是受了太多的苦“数值耗散”，即流量往往衰减得太快与实际实验相比。另一方面，在计算机图形应用程序中，这并没有那么糟糕，特别是在交互式系统中，动画师应用程序“保持流程处于活动状态”外力。事实上，根本不衰减的流动可能过于混乱且难以控制。正如我们的结果所示尽管数值有限，我们仍然能够产生漂亮的旋流耗散。
   在本文中，我们主要将流程应用于模拟类似气体的现象。我们使用求解器来更新两者流动和流动内的密度运动。进一步为了增加动画的复杂性，我们将纹理坐标与密度一起平流[13]。通过这种方式我们能够即使在低浓度的情况下也能合成高度详细的“纤细”气流分辨率网格。我们相信，基于物理的结合流体求解器和固体纹理是最有前途的方法在计算机图形学中实现高度复杂的流程。
   下一节介绍纳维-斯托克斯方程和推导得出我们的解决方法。该部分包含获得知识所需的所有基本思想和技术稳定的流体求解器。由于它依赖于复杂的数学技术，它是用数学物理术语写的大多数计算机图形学研究人员应该熟悉在基于物理的建模中。以应用为导向的读者仅希望实现我们的求解器可以完全跳过第 2 节，直接进入第 3 节。在那里我们描述了我们的实现求解器，提供足够的信息来编码我们的技术。
   第 4 节专门介绍几个应用程序，这些应用程序演示了我们新求解器的强大功能。最后，在第 5 节中我们总结并讨论未来的研究。为了将其保持在短期内论文中，我们决定不包含“教程类型”的部分流体动力学，因为有许多优秀的教科书提供了必要的背景和直觉。不会的读者具有流体动力学背景并且希望完全理解本文中的方法的人应该参考这样的文本。
   有数学倾向的读者可能希望从 Chorin 和 Marsden 的优秀著作开始 [3]。具有工程学的读者另一方面可以参考 Abbott 的教学书籍 [1]。此外，Foster 和 Metaxas 的论文很好地介绍了从流体动力学到计算机图形学社区的概念。

四、稳定的纳维-斯托克斯

4.1 基本方程

   在本节中，我们将介绍纳维-斯托克斯方程导致我们稳定求解器的操作。密度和温度几乎恒定的流体由速度描述场 u 和压力场 p。这些数量通常会有所不同在空间和时间上取决于周围的边界体液。我们用 x 表示空间坐标，对于二维流体来说是 x = (x; y)，对于三维流体来说是等于（x；y；z）。我们决定不专门研究我们的结果对于特定维度。因此，除非另有说明，所有结果对于二维和三维流均有效。
  假设在某个初始时间 t = 0 内速度和压力已知，则给出这些量随时间的演变由纳维-斯托克斯方程 [3]：

  其中 $\nu$是流体的运动粘度， $\rho$是其密度，$f$是外力。有的读者可能对此比较陌生纳维-斯托克斯方程的紧凑版本。等式2 是三个（二维中的两个）分量的矢量方程速度场。 “$\cdot$”表示向量之间的点积，而符号$▽$是空间偏导数的向量。更准确地说，二维 $▽=[\partial /\partial x,\partial /\partial y]$且三维$▽=[\partial /\partial x,\partial /\partial y,\partial /\partial z]$。我们也用过简写符号${▽}^2=▽\cdot ▽$ 。纳维-斯托克斯方程通过假定流体质量守恒而获得（方程 1）和动量（方程 2）。我们建议读者参考任何标准文本流体力学的实际推导。这些方程也必须补充边界条件。在本文中我们将考虑在实际应用中有用的两种边界条件：周期性边界条件和固定边界条件边界条件。在周期性边界的情况下，流体为定义在 n 维环面 (n = 2; 3) 上。在这种情况下有没有墙壁，只是包裹周围的流体。尽管在实践中不会遇到这种流体，但它们在创建不断演变的纹理贴图。此外，这些边界条件导致非常优雅的实现，使用快速傅里叶变换作为如下所示。我们考虑的第二种边界条件是当流体位于某个有界域 D 中时。在这种情况下，边界条件由定义在上的函数 uD 给出域的边界$\partial D$。请参阅 Foster 和 Metaxas 的工作在热的情况下对这些边界条件进行了很好的讨论液体[7]。无论如何，边界条件应该是这样的速度场的法向分量在边界处为零；无论如何都要穿越墙壁。
  纳维斯托克斯方程中出现的压力场和速度场实际上是相关的。通过组合方程可以得到速度的单个方程-1 和等式-2.我们简单介绍一下概述得出该方程的步骤，因为它是我们算法的基础。我们遵循 Chorin 和 Marsden 的治疗主题（第 36 页，[3]）。数学结果，称为Helmholtz-Hodge 分解，指出任何向量场 w 都可以唯一地分解为：
for。

  其中 u 的散度为零：$▽\cdot u=0$ 并且 q 是标量场。任何矢量场是质量守恒场和梯度之和场地。这个结果允许我们定义一个算子 P 来投影任意向量场 w 到其无散度部分 $u=Pw$。这运算符实际上是通过将两边相乘$▽$来隐式定义的等式-3 得到：

  这是标量场 q 的泊松方程，其中 Neumann $\partial D$ 上的边界条件$\frac{\partial q}{\partial n}=0$。该方程的一个解是用于计算投影 u

  如果我们在方程两边应用这个投影算子： 2 我们得到速度的单个方程

  (5)我们使用了 Pu = u 和 Prp = 0 的事实。这是我们将根据它开发稳定流体的基本方程所以。

图 1：我们的求解器的一个模拟步骤由多个步骤组成。前三个步骤可能会将这个领域带出发散的空间自由字段。最后的投影步骤确保场是发散的整个模拟步骤结束后免费。

图 2：为了求解平流部分，我们追踪了领域在时间上落后。因此 x 处的新速度为t 时间前粒子在旧位置的速度p(x;t)
  

4.2 解决方法

   等式。 5 从初始状态 u0 = u(x; 0) 通过步进求解以时间步长 t 遍历时间。让我们假设该字段有已在时间 t 解决，并且我们希望计算 a 处的场稍后时间 t + t。我们解决方程。 5 在时间跨度 t 的四分之四内脚步。我们从前面的解 w0(x) = u(x; t) 开始时间步长，然后依次求解右侧的每一项等式的一侧5，然后投影到发散的自由场上。
一般流程如图 1 所示。步骤是

X）：
   然后，时间 t + t 处的解由最后一个速度场给出：u(x; t + t) = w4(x)。通过迭代这些获得模拟脚步。我们现在更详细地解释如何计算每个步骤。
   最容易解决的项是添加外力 f 。如果我们假设力在整个过程中变化不大时间步长，

t) 是力对场的影响的良好近似时间步长t。在交互式系统中，这是一个很好的近似值，因为力仅在每个时间步开始时施加。
  下一步考虑流体的平流（或对流）对其自身的影响。流体中某处的扰动按照表达式 (u r)u 传播。本术语使纳维-斯托克斯方程成为非线性的。福斯特和梅塔克萨斯使用有限差分求解该分量。他们的方法仅当时间步足够小时才稳定
  t <=juj，其中是它们计算的间距网格。因此，对于小的分离和/或大的速度，非常必须采取小的时间步骤。另一方面，我们使用完全不同的方法，从而产生无条件稳定的结果求解器。无论时间步长有多大，我们的模拟都会永远不要“爆炸”。我们的方法基于一种称为特征法的偏微分方程求解技术。
  由于这种方法对于获得我们的稳定至关重要求解器，我们在附录 A 中提供了所有数学细节。
但是，方法可以直观地理解。在每个时间步所有流体粒子都通过流体本身的速度移动。
  因此，要获得新时刻 t+t 时 x 点的速度，我们在一段时间内通过速度场 w1 回溯点 xt。这定义了对应于速度场的部分流线的路径p(x;s)。 x 点的新速度为设置为粒子现在位于 x 时的速度地点$\Delta t$ 时间

t））：
   图 2 说明了上述内容。这种方法有几个优点。最重要的是它是无条件稳定的。确实，从上面的方程我们观察到新的最大值字段永远不会大于前一个字段的最大值。
  其次，该方法非常容易实现。实践中所需要的只是一个粒子示踪器和一个线性插值器（参见下一节）。因此该方法既稳定又简单实现任何计算机图形流体解算器的两个非常理想的属性。我们采用了类似的方案来移动密度通过用户定义的速度场[19]。方法的版本其他研究人员也使用了这些特征。应用程序用于可视化流场 [13, 18] 或改进气体模拟的渲染 [21, 5]。我们的应用这项技术根本不同，因为我们用它来更新速度场，以前的研究人员没有动态地动画。第三步解决了粘度的影响，等价于到扩散方程

r2w2:
   这是一个标准方程，许多数值过程都适用于此已经开发了。解决这个问题最直接的方法方程是离散化扩散算子r2然后做正如 Foster 和 Metaxas 所做的那样，是一个明确的时间步长 [7]。然而，这当粘度较大时，该方法不稳定。因此，我们更愿意，使用隐式

   w2(x)；其中 I 是恒等运算符。当扩散算子为离散化，这会导致未知数的稀疏线性系统字段 w3 。然而，解决这样的系统可以有效地完成（见下文）。
  第四步涉及投影步骤，这使得由此产生的场发散自由。正如上一小节所指出的，这涉及到定义的泊松问题的解决.

w3 请求：
   因此，投影步骤需要一个好的泊松解算器。Foster 和 Metaxas 使用松弛法求解了类似的方程方案。然而，放松方案的收敛性较差，并且通常需要多次迭代。福斯特和梅塔克萨斯报道说即使经过很少数量的松弛步骤，他们也获得了良好的结果。然而，由于我们使用不同的方法解决平流步骤，我们必须使用更准确的.
   事实上，当场接近发散自由。更重要的是从视觉角度来看视图中，投影步骤迫使场产生涡流，从而导致更多类似漩涡的运动。由于这些原因，我们使用了
更准确的投影步骤求解器。
   泊松方程在空间上离散化后，变为稀疏线性系统。因此，投影和粘度步骤都涉及大型稀疏方程组的求解。
例如，多重网格方法可以解决稀疏线性系统线性时间[10]。由于我们的平流求解器在时间上也是线性的，我们提出的算法的复杂度是 O(N)。
   Foster 和 Metaxas 的求解器具有相同的复杂性。该性能理论上是最佳的，因为对于复杂的流体，任何算法必须至少参考计算网格的每个单元。