用图形处理器仿真流体流动

将图形处理单元用作通用处理器(可称作通用GPU,缩写为GPGPU),已成为计算机图形领域中日益流行的研究课题。当处理过程从CPU转移到GPU时,随着GPU处理能力的增长,几项并行技术正在加速发展。

不过,为了说明某些通用算法如何在GPU上实现,先了解现代图形硬件架构是很重要的。本文着重介绍一个采用Navier-Stokes方程进行不可压缩流体仿真的应用示例,将证明GPU如何很好地满足流体仿真的要求。

图形处理器已经被高度优化以绘制三维几何图形,而且具有几条平行的顶点管道来同时处理多个顶点。同样地,一个三角扫描转换器从处理过的顶点上产生像素,然后大量的平行像素管线被用来生成最终图像。

总之,图形管线的设计使不同顶点及像素之间只有很少或没有相关处理。因此,单个的顶点和像素可以被单独处理,而在其后的管线阶段中合并处理结果。

图1:描述流体流动行为的Navier-Stokes方程。

为了充分利用GPU的处理能力,所采用的算法应具有与3D几何绘制相符的特征。这些特征包括:

统一的网格数据输出。绘制算法的输出是二维图像。在统一采样网格上运行的算法可最大限度地利用GPU的处理优势。

局部计算。对于每个被处理值,GPU都会执行一系列独立的数据采集与处理操作。而全局操作(例如计算一个大型二维阵列的最大值),则要求有多个数据处理Pass,这在GPU上并非特别有效率。

无散布数据。GPU不支持散布操作。这就是说,GPU不能指定其输出位置,因为这是由光栅化(rasterization)过程预先确定的。

为了充分利用硬件功能,GPGPU算法经常将输入数据存储在一个或多个片外浮点纹理图像中,并执行几个绘制pass(其中执行了计算)。你可以将这些纹理图像看成是二维浮点阵列。大多数GPGPU算法都采用可绘制四边形的处理pass,该四边形覆盖了一个表现纹理的整个矩形范围。

在这些绘制pass期间,可编程的像素着色器程序读取输入纹理、处理数据并写入结果。诸如流体仿真这样的复杂算法,可通过采用几个这样的pass处理数据来实现。

近来,GPU已被用于各种算法,例如具有上述特征的图像与声音处理。接下来我们讨论如何采用GPU仿真和绘制流体流动。

图1给出著名的Navier-Stokes方程式,它描述了流体流动的行为。这些方程式提供了一种方法来计算在一小段时间内流体密度与速度的变化。

前两个方程式计算流体速度与密度的变化。特别指出,这两个场都基于速度场(第一项)平流输送,然后慢慢地随时间散布(第二项)并由外力驱动(第三项)。最后一个方程式应用在被称为发射操作的最后步骤中。发射步骤用于更新速度场以使其维持大部分保留(mass-conserving)。

为了在图形处理器的离散领域中对这些连续方程式建模,我们沿用Jos Stam在1999年提出的“semi-La-grangian”技术。由于该技术能执行所有局部计算而不要求数据散布,因此它能被扩展至GPU。

我们首先将二维空间分成由小单元(即样本)组成的正方形网格,并执行所有局部计算,根据其自身及前一次迭代中的相邻值来更新采样的密度与速度。然后,通过相邻样本间的有限差分执行所有的微分计算。

为了将流体密度与速度存储到二维网格中,我们使用两个浮点纹理图像。我们还采用额外的发散、压力草稿图像,以及多个绘制pass来计算每个仿真时间段的稳定结果。该算法步骤如下:

外力。首先执行两个绘制pass,根据密度与速度纹理的增量来计算注入流体的外力。

扩散。流速场的扩散是通过局部松驰来执行的。实质上,每个样本都提取其相邻样本的流速值,并相应更新其自身值。

图2:该图显示了一个着色器。
它提取四个纹理,执行计算并将结
果写入其中一个表现纹理。

平流。对于给定样本,通过提取该样本的流速值来执行平流输送。基于这个流速,我们对那些现在应该被输送到当前单元位置上的数据进行采样。流速场用于平流输送密度场与流速场本身。

大部分保留。在最后的步骤,必须确保已被前面各步骤所改变的流速场能被大部分保留。计算流速场中的大部分保留成分涉及到求解泊松方程,这可采用局部松弛来近似,并使用几个绘制pass。

由于大部分保留步骤相当昂贵,因此我们想将计算集中在高压区域以便更有效地分配GPU资源。下面将介绍通过使用初阶Z精选硬件来完成此项工作。初阶Z精选是一种GPU功能,用于有效图形处理,但我们可以灵活地利用该特性来优化通用计算。

在进行像素处理之前,图形硬件要执行内插深度或Z值检查,并与Z缓冲器中的Z值进行比较。如果此项测试失败,则不处理该像素。实时绘制中所使用的初阶Z精选可用于优化流体仿真算法。

昂贵的大部分保留步骤可被作为一系列绘制pass来执行,以便求解采用松弛方法的线性系统。当结果近似于此线性系统时,迭代次数(绘制pass)越多,结果越精确。我们不是在所有单元上都执行相同次数的pass,而是用初阶Z精选对高压区执行更多的pass,对低压区到无压区执行较少的pass。

这可通过执行其他廉价的绘制pass来实现,这些绘制pass根据压力缓冲器中的单元及其相邻单元的最大值来设置仿真中每个单元的Z值。当执行昂贵的压力计算pass时,我们将深度比较状态设置为“小于/等于”,并线性地增加每个发射pass的Z值。

随着这些pass被绘制出来,被处理的单元数目逐渐减少,从而提高了仿真流速。

作者:Pedro Sander

psander@

应用研究部成员

ATI研究公司