CN101350106B - 用于邻近粒子搜索的数据结构构建方法 - Google Patents

Abstract

提供一种邻近粒子搜索的数据结构的构建方法。该数据结构的构建方法，通过GPU等流处理器执行，将3维数字图像数据变换为2维数字图像数据，由此来对物理对象物的模拟进行运算。在使用GPU的情况下使用桶纹理，该桶纹理由将对应于出现物理对象物的空间而生成的存储器空间进行分割的桶构成。在GPU的情况下，具有使用顶点着色器，读出在桶中容纳的粒子的所述粒子号码和粒子坐标，并写入2维图像空间的步骤。

Description

用于邻近粒子搜索的数据结构构建方法 

技术领域

本发明涉及通过用GPU等流处理器(Streaming processer)完全执行粒子法模拟，来实现将邻近粒子搜索高速化的方法，例如使实时地进行许多刚体的撞击计算等成为可能的用于邻近粒子搜索的数据结构的构建方法、其程序以及存储该程序的存储介质。

背景技术

所谓“粒子法模拟”，是将对象物体作为许多粒子的集团(集合)来获取，关注许多粒子中各个粒子的行为来进行计算的模拟方法。在计算机的计算方面，在上述粒子中通常使用直径固定的“球”。当通过球的集合体来表现物体的形状时，根据该球的大小(直径)的决定方法，物体的形状的近似度或空间析像度发生变化。另外，通过使球的直径变化、使空间析像度变化，可以控制基于计算机的模拟计算中的计算精度和计算速度。

近年来，利用上述粒子法模拟进行许多刚体的撞击计算、流体或粉体等的行为计算，在计算机的显示装置的画面上进行了撞击或流动等的模拟的研究。该研究结果被记载在例如田中正幸等，“使用了粒子法的刚体计算方法的开发和向计算机图形学的应用”，日本机械学会，第19次计算力学讲演会讲演论文集，701-702页，2006。这些模拟技术的研究作为计算机图形学(CG)的基础技术是重要的。

在上述粒子法模拟的研究中，以往大多进行利用CPU(Central ProcessingUnit)生成计算算法的方法。但是，不太进行积极利用GPU(Graphics ProcessingUnit)等流处理器(Streaming processor)的粒子法模拟的研究。例如，根据利用GPU的现有的某个研究，通过CPU执行邻近粒子搜索，将通过CPU得到的搜索数据转发给GPU来进行图像处理。根据这种技术以如下方式构成：将CPU用作主要计算机执行部，通过该CPU执行邻近粒子搜索的计算算法，GPU作为辅助计算机执行部仅执行图像显示处理的计算算法。因此，根据现有的计 算技术，无法构成仅在GPU的内部完全执行粒子法模拟中的邻近粒子搜索以及图像处理的计算算法。

由于在现有的粒子法模拟的计算中，无法仅通过GPU来完全执行全部的计算(用于邻近粒子搜索的数据结构的构建、图像处理等)，因此计算效率不高，GPU无法充分发挥其本来具有的图像处理能力，无法进行许多刚体的撞击计算等实时的刚体模拟、流体模拟等。所以，要求开发出完全利用GPU或其它流处理器，在GPU等上执行用于邻近粒子搜索的数据结构的构建，最大限度地发挥其图像处理能力的方法。

发明内容

本发明的目的在于，鉴于上述问题，提供用于邻近粒子搜索的数据结构的构建方法、其程序以及存储了该程序的存储介质。所述方法可以在GPU等流处理器上完成粒子法模拟的邻近粒子搜索中所使用的数据结构的构建，高效地利用GPU等，使其发挥最大限度的(图像)处理能力，使高速计算成为可能，例如实现实时的刚体模拟或流体模拟等的计算处理/显示处理。

为了实现上述目的，本发明的用于邻近粒子搜索的数据结构的构建方法、其程序以及存储了该程序的存储介质，如下构成。

数据结构的构建方法，通过GPU等图形运算处理器来执行，被用于通过将3维数字图像数据变换为2维数字图像数据来对物理对象物的模拟进行运算的方法中，所述3维数字图像数据是根据在视频存储器的内部生成的图像数据而构建的，在该数据结构的构建方法中，使用了以许多粒子表现物理对象物，根据视频存储器内的3维图像数据结构中的许多的粒子的相互的容纳位置关系来决定物理对象物的存在状态的方法。3维图像数据结构中包含的桶纹理(bucket texture)，由将对应于出现物理对象物的空间而生成的存储器空间进行分割的桶(bucket)构成，将分配给许多粒子的粒子号码存储在许多桶的每一个中。所述数据结构的构建方法具有：通过图形运算处理器使用顶点着色器，读出在桶中存储的粒子的所述粒子号码和粒子坐标，写入2维图像空间的写入步骤。在该写入步骤中，其特征在于，当桶中存在多个粒子时，在模拟的一个时步(time step)中分为多次写入(用于多次处理的渲染)，一个一个地写入所述多个粒子的各个的粒子号码和粒子坐标。

上述数据结构的构建方法中，在作为图形运算处理器而使用GPU的情况下，使用由GPU准备的深度测试、颜色蒙板、模板测试的各个功能，来进行用于在多次处理的渲染的各个渲染中写入的粒子的粒子号码和粒子坐标的选择、以及写入的像素中的通道的选择。

另外，在使用桶含GPU等在内的流处理器的情况下，也可以使用由流处理器准备的多个计数器的各功能，来进行在多次处理周期的各步骤中写入的粒子的粒子号码和存储器位置(Memory location)的选择、以及为了体素(Voxel)而分配的存储器中的存储器位置的选择。

即，用于邻近粒子搜索的数据结构的构建方法，其由流处理器执行，用于通过将3维数字图像数据变换为2维数字图像数据来对物理对象物的模拟进行运算的方法，所述3维数字图像数据是根据在存储器内部生成的数据阵列而构建的。在所述数据结构的构建方法中，使用了下述方法，即通过许多粒子来表现物理对象物，根据存储器内的3维图像数据结构中的许多粒子的相互的容纳位置关系，决定物理对象物的存在状态，3维图像数据结构中包含的桶阵列，由将对应于出现物理对象物的空间而生成的存储器空间进行分割的桶构成，将分配给许多粒子的粒子号码存储在许多所述桶的每一个中。具有写入步骤：其通过流处理器，读出在桶中存储的粒子的所述粒子号码和粒子坐标，使用流处理器的分散操作(scattering operation)来写入存储器中的阵列中。在写入步骤中，当桶中存在多个粒子时，在模拟的一个时步中划分为用于多次处理的内核(kernel)，一个一个地写入多个粒子的各个粒子的粒子号码和粒子坐标。

数据结构的构建方法，在上述方法中，当桶中存在最多4个粒子时，在由写入粒子号码的像素等的RGBA组成的4个通道的每一个通道中，独立地从按照升序或者降序的粒子号码的粒子中顺次存储该粒子的号码。

在上述数据结构的构建方法中，在作为图形运算处理器使用GPU的情况下，更具体的说，具有以下步骤：当假定4个粒子的粒子号码为i_b0、i_b1、i_b2、i_b3(升序的情况下i_b0＜i_b1＜i_b2＜i_b3，降序的情况下i_b0＞i_b1＞i_b2＞i_b3)时，在向第一通路的像素的写入中，以最大值对深度缓冲区进行初始化，使用深度测试使具有较小值的深度缓冲区合格，由此，在所述R通道中写入粒子号码i_b0，进行用于第一次处理的渲染的步骤；在向第二通路的像素的写入中，使用颜色蒙板 来禁止向所述R通道的写入并且使用了模板测试，由此，在所述G通道中写入粒子号码i_b1，进行用于第二次处理的渲染的步骤；在向第三通路的像素的写入中，使用颜色蒙板来禁止向所述R和G的各通道的写入并且使用了模板测试，由此，在所述B通道中写入粒子号码i_b2，进行用于第三次处理的渲染的步骤；在向第四通路的像素的写入中，使用颜色蒙板来禁止向所述R、G和B的各通道的写入并且使用了模板测试，由此，在所述A通道中写入粒子号码i_b3，进行用于第四次处理的渲染的步骤。

一般，在使用流处理器的情况下，代替上述像素(written pixel)而使用分配给体素的存储器，代替上述通道而使用为了体素而分配的存储器中的存储器位置(Memory Location)，而且，代替由GPU准备的上述深度测试、颜色蒙板、模板测试的各功能而使用与其分别对应的计数器的功能。

即，具有以下步骤：当假定4个粒子的粒子号码为i_b0、i_b1、i_b2、i_b3(升序的情况下i_b0＜i_b1＜i_b2＜i_b3，降序的情况下i_b0＞i_b1＞i_b2＞i_b3)时，在向第一通路中的为了体素而分配的存储器进行的写入中，通过将计数器初始化为0并且接受向第一存储器位置的最初的写入，在第一存储器位置中写入粒子号码i_b0，执行用于第一次的例程(routine)的内核的步骤；在向第二通路中的为了体素而分配的存储器进行的写入中，通过使用计数器来禁止向第一存储器位置的写入并且向第二存储器位置写入下一值，由此，在第二存储器位置中写入粒子号码i_b1，执行用于第二次的例程的内核的步骤；在向第三通路中的为了体素而分配的存储器进行的写入中，通过使用计数器来禁止向第一和第二存储器位置的写入并且向第三存储器位置写入下一值，由此，在第三存储器位置中写入粒子号码i_b2，执行用于第三次的例程的内核的步骤；以及在向第四通路中的为了体素而分配的存储器进行的写入中，通过使用计数器来禁止向第一、第二和第三存储器位置的写入并且向第四存储器位置写入下一值，由此，在第四存储器位置中写入粒子号码i_b3，执行用于第四次的例程的内核的步骤。

数据结构的构建方法，在上述方法中，理想的是上述物理对象物是许多刚体，物理对象物的模拟图像是基于许多刚体的撞击计算的模拟图像。

关于数据结构的构建方法，作为流处理器的一例的GPU是图形运算处理器。

另外，本发明的用于邻近粒子搜索的数据结构的构建程序，是用于通过包含GPU等在内的流处理器执行构成上述邻近粒子搜索方法的各步骤等的程序。

根据本发明，收到到如下效果。

应用粒子法模拟来近似刚体、流体、粉体等，通过实施用于邻近粒子搜索的数据结构的构建方法，在计算机上对物理现象进行模拟的方法中，在GPU等流处理器上完成用于执行用于邻近粒子搜索的数据结构的构建计算的各步骤，高效地利用GPU等，可以发挥其最大限度的图像运算处理能力，可以进行高速计算，由此可以实现例如实时的刚体模拟或流体模拟等。

另外，通过将刚体模拟或流体模拟等作为模拟CG而用于图像显示，作为最佳的游戏等的重要要素技术，可以用于水的流动等自然物的表现。

附图说明

图1是表示实施本发明的邻近粒子搜索中所使用的数据结构的构建方法的计算机系统的基本结构的系统结构图。

图2是表示在计算机系统中安装的GPU的内部基本结构的框图。

图3是基于与10922个圆环落下·撞击的物理现象相关的刚体模拟的结果图像(CG)的图。

图4是基于与16384个国际象棋中的棋子落下·撞击的物理现象相关的刚体模拟的结果图像(CG)的图。

图5是表示使用了GPU的刚体模拟的高速化的计算算法的流程图。

图6是表示桶构建的处理步骤的详细内容的流程图。

图7是表示使用了GPU的刚体模拟的高速化的计算算法的变形例的流程图。

图8是基于与多个球体形状的流体落入流体时的物理现象相关的流体模拟的结果图像(CG)的图。

图9是基于漏斗状容器中容纳的粉体向下方流动，从容器的下端开口部落到外侧，被配置在下侧的容器容纳时的粉体模拟的结果图像(CG)的图。

具体实施方式

以下，根据附图说明本发明的优选实施例。图1表示在粒子法模拟中，实 施了在邻近粒子搜索中使用的数据结构的构建方法的计算机系统的基本结构。

在该实施例的说明中，作为主要实施例，说明作为图形运算处理器而使用例如GPU的例子。

计算机系统10例如利用PC(Personal Computer)而构成，具备作为运算处理部的CPU11和GPU12。另外，计算机系统10作为存储装置而具有主存储器13A和视频存储器(VRAM)13B，作为外围装置而至少具有输入装置14以及显示装置15。CPU11、GPU12、主存储器13A、输入装置14和显示装置15经由总线16相互连接。输入装置14经由输入接口17与总线16连接，显示装置15经由输出接口18与总线16连接。另外，视频存储器13B与GPU12直接相连。对于视频存储器13B，通过GPU12进行数据交换。

CPU(中央运算处理装置)11例如是PC等所配备的通常的运算处理部。另外，GPU12是“Graphics Processing Unit(图形运算处理器)”，执行图像运算处理。该GPU12内置多个处理器而形成，被用作并行计算机。作为GPU12，例如使用NVIDIA公司生产的“GeForce 7900GTX”。此外，在GPU12中所使用的产品不限于此，可以使用具有类似的结构·功能的图像处理用处理器或其它流处理器。该GPU12执行存储在主存储器13A中的邻近粒子搜索程序13A-1，实施粒子法模拟中的邻近粒子搜索运算。邻近粒子搜索程序13A-1执行在邻近粒子搜索中使用的数据结构的构建。另外，GPU12执行存储在主存储器13A中的通常的图像处理程序13A-2，使用在邻近粒子搜索的计算中获得、并且存储在像素中的图像数据来实施图像处理，使用所生成的图像数据在显示装置15的画面上进行描绘(用于显示的渲染(Rendering))。在显示装置15的画面上，根据粒子法模拟的计算来显示CG图像。

在本实施例的以下说明中，对实时的刚体模拟的运算进行说明，所述运算基于执行依据粒子法模拟的许多刚体的撞击计算。在该模拟的运算中，仅通过GPU12执行用于邻近粒子搜索的数据结构的构建方法。

接下来，图2表示所述GPU12内的重要部分的基本结构。GPU12具有并行配置的多个处理器21，基于这些处理器进行的并行计算处理，来执行用于显示(描绘)的渲染。特别是根据本实施例的邻近粒子搜索中使用的数据结构的构建(桶构建)的方法，如后所述，当1个桶中存在4个粒子时，将用于处 理的渲染重复4次，将4个粒子的粒子号码和粒子坐标写入1个像素的RGBA通道的每一个通道中。通过GPU12中的该桶构建的处理，实现图像运算处理的高速化。更具体而言，当GPU12是“GeForce 7900GTX”时，其具有在前级并行配置的8个顶点着色器(VS)22和在后级并行配置的24个片元着色器(FS)23。通过顶点(Vertex)处理器(21)执行的程序是顶点着色器(VertexShader)22，通过片元处理器(21)执行的程序是片元着色器(Fragment Shader)23。

如上所述，顶点处理器是并行配置的处理器，是专门用于一齐处理许多顶点的坐标变换(将输入的顶点坐标变换为被渲染的图像的空间中的坐标)的处理器。GPU12与CPU相比，可以高速地计算顶点的坐标变换。由顶点着色器指定顶点处理器的动作。

而且，在GPU12中具有深度测试功能部24、颜色蒙板(color mask)功能部25、模板(stencil)测试功能部26。另外，方框27是着色器指令分配部。

用于实现本发明的用于邻近粒子搜索的数据结构的构建方法的邻近粒子搜索程序13A-1，利用上述GPU12的内部构造，仅通过GPU12可完全执行。

接着说明许多刚体撞击的刚体模拟的实施例。通过该刚体模拟的说明，可以明白基于邻近粒子搜索程序13A-1的执行而实施的、用于邻近粒子搜索的数据结构的构建方法的内容。

作为许多刚体的例子，例如是许多圆环(torus)材料、或者是许多国际象棋中的棋子。图3和图4表示作为本实施例的刚体模拟结果而得到的CG的例子。图3表示例如基于与10922个圆环落下·撞击的物理现象相关的刚体模拟的撞击计算的结果图像(CG)，图4表示例如基于与16384个国际象棋中的棋子落下·撞击的物理现象相关的刚体模拟的撞击计算的结果图像(CG)。

在上述的许多刚体的撞击计算中，通过应用粒子方法(粒子法模拟)，将许多刚体的各自的形状全部作为相同大小的球(粒子)的集合体来近似表现地进行处理。基于球的刚体的形状近似，根据该刚体的多边形模型(polygonmodel)来进行，例如使用由本发明人等提出的应用了GPU的体素(voxel)化方法(“Fast solid voxelization using hardware”，Transactions of JSCES，page No.20060023，2006)而生成。此外，考虑计算精度和计算速度而适当地决定所使用的球的大小(空间析像度)。

在各个刚体作为球的集合体来近似表现而形成的许多刚体的撞击计算中，执行(1)各个刚体的计算、(2)刚体间的撞击的检测处理(构成刚体的球之间的距离的计算)、以及(3)基于撞击的响应处理(施加于球的力和转矩的计算)。以下，分别说明这些计算的思路。

(1)刚体的计算：

刚体的计算所需要的物理量如下述所示，包括重心的坐标、平移速度、表示旋转量的四元数(quaternion)、角速度、惯性张量。

X：重心的坐标

V：平移速度

q：表示旋转量的四元数

w：角速度

I：惯性张量

刚体的运动的计算，分为平移运动和旋转运动来计算。平移运动的计算内容由(1)-(3)式表示，旋转运动的计算内容由(4)-(9)式表示。

施加于刚体的力F如(1)式所示，可分解为由于撞击而产生的Fc和除此以外的外力Fe。

F＝Fc+Fe    (1)

使用施加于刚体的力F，如(2)式所示那样计算刚体的平移运动量P的时间微分。

$\frac{\mathrm{dP}}{\mathrm{dt}}=F---\left(2\right)$

使用平移运动量P，由P＝Mv(M：刚体的质量)计算重心的速度v。当使用重心的速度v时，刚体的重心位置x的时间微分如(3)式所示。

$\frac{\mathrm{dx}}{\mathrm{dt}}=v---\left(3\right)$

施加于刚体的F引起刚体的旋转运动，使角运动量L变化。角运动量L的时间微分以(4)式表现。r是从刚体的重心到力的作用点的矢量。

$\frac{\mathrm{dL}}{\mathrm{dt}}=r\×F---\left(4\right)$

使用上述角运动量L，通过(5)式计算角速度w。

w＝I(t)^-1L  (5)

在(5)式中，I(t)^-1是在时刻t的惯性张量的逆矩阵。

使用初始状态下的刚体的惯性张量I(0)的逆矩阵I(0)^-1，通过(6)式计算I(t)^-1。

I(t)^-1＝R(t)I(0)^-1R(t)^T  (6)

在(6)式中，R(t)和R(t)^T是在时刻t的旋转矩阵及其转置矩阵。使用(5)式中计算出的角速度w来更新四元数q。根据角速度w，通过(7)式计算四元数的变化量dq。

$\mathrm{dq}=[\cos \left(\frac{\mathrm{\θ}}{2}\right),a\sin \left(\frac{\mathrm{\θ}}{2}\right)]---\left(7\right)$

在此，旋转轴a和旋转角θ通过(8)式和(9)式来计算。

$a=\frac{w}{\left|w\right|}---\left(8\right)$

θ＝|wdt|(9)

使用(7)式的dq和在时刻t的q(t)，通过(10)式计算dt后的四元数q(t+dt)。

q(t+dt)＝dq×q(t)(10)

(2)撞击的检测处理：

为了检测刚体的撞击，计算构成刚体的球间的距离。若球间的距离小于球的直径则视为两球已撞击。在许多刚体的撞击检测中，需要针对构成各个刚体的全部球来进行上述距离计算。计算量与球数的平方成比例，刚体的数量增多时，球的数量也增多，计算量也增大。因此，在本实施例的撞击检测处理的结构中，为了减少计算量，将许多刚体所存在的虚拟空间(存储器中的3维数据存储空间)分割成被称为“桶(或者也称为格栅(Grid))”的格子(3维空间格子)，将表现许多刚体中的每一个刚体的许多球容纳在桶中。1个桶构成1个立方体，假定在使用了使桶的一边的长度与球直径相同的3维数据结构时，可能与某个球i撞击的球被限定为在与容纳球i的桶相邻的3³个各桶的内部所容纳的球。根据如上所述地构建桶来进行撞击检测，由此可以减少计算量，并且可以大幅度地降低计算成本。在后面描述本实施例的桶构建的详细情况。

(3)基于撞击的响应：

在撞击响应中计算施加在构成1个刚体的球上的力。在撞击力的计算中，关于接触的两个球，使用了在离散单元法(distinct element method)中使用的线性弹性常数和阻尼常数。关于该撞击力的计算的内容表示如下。

设弹性常数为κ、阻尼常数为η，对于正在撞击的球施加与陷入量成比例的回弹力以及与相对速度成比例的阻尼力。有球i和球j两个球，若设它们的直径为d，则当两球距离|rij|小于d时发生撞击，施加了(11)式的基于弹性的力fs和(12)式的基于阻尼(damper)的阻尼力fd。

$\mathrm{fs}=-\mathrm{\κ}(d-|\mathrm{rij}\left|\right)\frac{\mathrm{rij}}{\left|\mathrm{rij}\right|}---\left(11\right)$

fd＝-η(vj-vi)      (12)

在(11)式和(12)式中，rij＝rj-ri，ri、rj是球i、j的位置矢量。另外，剪切方向的摩擦力ft作为与剪切方向的速度vt成比例的力，通过(13)式进行计算。

ft＝κtvt        (13)

通过(14)式计算施加于球i的力fi。

fi＝fs+fd+ft     (14)

使用施加于球i的力fi，分别通过(15)式和(16)式计算由于撞击而施加于刚体的力Fc和转矩Tc。

(16)式中，r’i是球相对于刚体重心的相对位置矢量。

接下来，对于仅使用GPU12的情况下的邻近粒子搜索所使用的数据结构的构建方法的高速化进行说明。根据本实施方式，与许多刚体的撞击相关的上述各计算全部在GPU12上执行。由此可以实现上述计算的高速化。下面，说明GPU12中基于邻近粒子搜索程序13A-1的计算算法。

最初，说明为了上述计算算法而准备的数据结构。

当通过GPU12进行所述各计算时，将所述各物理量作为视频存储器13B中桶含的多块纹理13B-1(图1所示)而进行存储、保持。

具体而言，为了进行刚体计算，关于各刚体的重心坐标、四元数、平移运动量、角运动量，分别使用了存储在视频存储器13B中的两块纹理13B-1。为了进行撞击计算，关于各球的中心坐标、速度、桶、施加于各球的力、施加于各球的转矩，分别使用存储在存储器13B中的一块纹理13B-1。

在本实施例的情况下，在基于粒子法模拟的许多刚体的撞击计算中设想3维空间，而且准备如前所述那样分割该3维空间的3维空间格子、即桶。但是，通过当前可以使用的GPU12无法向3维纹理进行写入。因此，在本实施例中，利用嵌入(D×D)个所需个数的2维格子而组成的1块2维纹理来表现由(L×L×L)个格子数组成的3维桶。在这种情况下，当设3维桶上的下标为(i，j，k)，并且设2维纹理上的下标为(s，t)时，在所述两个下标之间，通过下面的(17)式和(18)式来建立关联。

s＝i+L(k-D[k/D])(17)

t＝j+Lf[k/D](18)

此外，i、j、k和s、t分别取桶含在[0，L-1]和[0，D×L-1]的范围中的值。

接下来，参照图5所示的流程图，说明使用了GPU12的高速化的计算算法。

当将进行用于在显示装置15的显示画面中进行描绘的1次计算(模拟)的单位时间设为1个时步时，该1个时步中的计算通过以下8个阶段来进行计算。

阶段1：球的物理量的计算         …(步骤S11)

阶段2：桶构建                …(步骤S12)

阶段3：撞击检测                 …(步骤S13)

阶段4：撞击响应(施加于球的力和转矩的计算)…(步骤S14)

阶段5：刚体平移运动量的计算          …(步骤S15)

阶段6：刚体角运动量的计算            …(步骤S16)

阶段7：刚体重心坐标的计算            …(步骤S17)

阶段8：刚体四元数的计算              …(步骤S18)

使用在上述计算(模拟)中最终得到的刚体的重心坐标和四元数进行用于显示的渲染处理，进行模拟结果图像的描绘(步骤S19)。在接下来的判断步骤S20中判断是否继续模拟，在判断结果为“是”的情况下，返回上述的步骤S11，重复步骤S11～S19；在判断结果为“否”的情况下，结束模拟。如上所述，针对用于一次显示的渲染处理(步骤S19)，在一个时步的时间内执行上述8阶段的计算步骤S11～S18。

如上所述，在步骤S12中执行的“桶构建”的计算中，进行邻近搜索的处理计算。该计算仅通过GPU12完全进行，由此实现计算的高速化。以下说明桶构建的计算的详细内容。

以下，说明上述步骤S11～S19的计算内容。

步骤S11：球的物理量的计算

为了通过粒子法计算许多刚体的撞击计算，计算构成各刚体的球的中心坐标和速度。与球的中心坐标和速度相关的计算内容如下述所示。

设刚体j的重心坐标、四元数、速度、角速度分别为Xj、Qj、Vj、Wj。

当设未旋转状态下，球i的中心对于刚体重心的相对位置矢量为ri时，通过(19)式计算以四元数Qj旋转状态下的相对位置矢量r’i。

r′i＝QjriQj^T    (19)

当设构成刚体j的球i的中心坐标和速度分别为xi、vi时，以(20)式和(21)式进行计算。

xi＝Xj+r′i    (20)

vi＝Vj+Wj×r′i    (21)

步骤S12：桶构建

所谓“桶构建”，是将为了表现许多刚体中的各个刚体而准备的许多球(粒子)的全部的号码，如前所述作为视频存储器13B的内部的数据结构，存储在作为分割空间的要素而准备的3维空间格子、即桶中的处理。具体而言，是将分配给许多球中的各个球的号码存储在许多桶的各个的内部的处理，是使与许多球中的各个球相关的数据分散(scattering)的处理。

在与桶的构建相关的上述处理中，使用了顶点(Vertex)着色器的写入功能(用于处理的渲染)。该写入功能，是可以将顶点写入作为输入的坐标中的功能。换言之，在顶点着色器中具有用于访问视频存储器13B的纹理的顶点纹理获取(VTE)功能，通过使用顶点着色器的该功能可以参照球的坐标。通过将顶点配置在与该球的坐标相对应的桶的坐标中，在GPU12中可以进行上述数据分散处理。

当将球的数量设为n(i₀，i₁，...，i_n-1)时，仅准备与球的数量相同的数量n(p₀，p₁，...，P_n-1)的上述顶点。这些顶点的各个对应于1个球。将顶点p₀，p₁，...，P_n-1各自对应的球的号码设为i₀，i₁，...，i_n-1。假设球的号码例如按照i₀＜i₁＜...＜i_n-1的顺序(升序)排列。此外，球的排列方法也可以是降序。

以上述对应关系为前提，顶点着色器使用上述顶点纹理获取(VTF)，参照前一次的时步中的全部球的坐标，根据该坐标计算桶纹理上的坐标，向该坐标输出球的号码，由此可以向桶纹理写入球的号码。

然而，在上述处理的情况下，仅在桶的内部存在1个球时可以正确地进行桶构建，反之，在1个桶的内部存在多个球时无法正确地构建。当有可能在1个桶内部存在多个球时，为了实现正确的桶构建，必须在参照该桶中容纳的球的数量的号码的同时，顺次将球逐个容纳在桶中。通过GPU12中所使用的并列处理，在用于一次处理的渲染中，无法参照已写入的球的号码，因此在该并行处理后，仅在各桶中写入1个号码。例如，当1个桶中存在4个球时，将4个顶点写入同一像素，因此最终在该像素中写入关于1个球的号码。

因此，在本实施例的桶构建中，构成为可以在1个桶b₁的内部容纳多个球。参照图6的流程图说明该方法。

在此，考虑在桶b_j中放入最多4个球的情况。将这4个球的号码作为i_b0、i_b1、i_b2、i_b3，写入1个像素的RGBA通道的各个的值。假设这4个球的号码按照i_b0＜i_b1＜i_b2＜i_b3的顺序(升序)排列。此外，4个球的号码的排列方法也可以是降序。并且，在桶构建中，通过将针对1个像素的写入划分为4次球号码的写入，即使在1个桶的内部存在最多4个球的情况下，也可以使用GPU12来进行桶构建。换言之，在桶b₁中存在4个球的情况下，通过在1时步中划分为4次写入(用于处理的渲染)，基于GPU12进行桶构建。

在所述顶点着色器(VS)22中，使用与顶点对应的球的坐标来计算该坐标的桶内的坐标，将球的号码传送到片元着色器(FS)23。在片元着色器(FS)23中将球的号码作为“颜色”和“深度”而写出。

参照图6，具体描述桶构建S12的处理内容。桶构建S12的工序由顺次进行的4条通路(pass)的步骤S101～S104组成。

向第一通路中的1个像素进行的最初的写入(步骤S101)，是在R通道中写入球号码i_b0的步骤。球号码i_b0是4个球中最小的号码。即，作为具有最小深度的像素而被写入，进行用于最初的处理的渲染。在该球号码i_b0的写入的情况下，首先设定用于选择R通道的颜色蒙板(步骤S31)，设定深度测试(步骤S32)，接着以最大值对深度缓冲区进行初始化(步骤S33)，使用所设定的上述深度测试，使具有较小值的深度缓冲区为合格，使该写入(顶点的渲染S34)成为可能。

向第二通路中的上述同一像素进行的写入(步骤S102)，是为了不覆盖已写入R通道的值，而在G通道中写入球号码i_b1的步骤。在该写入中，通过使用用于选择G通道的颜色蒙板(步骤S41)，禁止对RBA的各通道的写入。深度缓冲区仍然使用在第一通路中所使用的深度缓冲区，并且对深度测试进行设定，以便使具有较大值的深度缓冲区为合格(步骤S42)。但是，这样一来，在上述像素的G通道中写入最大的球号码i_b3，因此为了禁止该写入而使用模板测试(步骤S43)。在模板测试中以0对模板缓冲区进行初始化(步骤S44)，在模板测试中进行设定以使值增加。而且，在模板测试中进行设定，以便当值在1以上时使写入失败(步骤S43)。由此，在像素中进行了一次写入处理后，可以不对写入失败的球号码进行写入。于是，在向同一像素的第二通路的写入中，在G通道中写入球号码i_b1(步骤S45)。

此后，在向同一像素中的第三通路的B通道的球号码i_b2的写入(步骤S103)、第四通路的球号码i_b3的写入(步骤S104)的情况下，也进行与上述相同的处理，将球号码i_b2和球号码i_b3分别写入B通道和A通道。第三通路的步骤S103具有用于选择B通道的颜色蒙板的设定、模板缓冲区的清除、顶点的渲染(写入)的步骤S51～S53，第四通路的步骤S104具有用于选择A通道的颜色蒙板的设定、模板缓冲区的清除、顶点的渲染(写入)的步骤S61～S63。

如上所述，在使用GPU12的邻近粒子搜索中所使用的数据结构的构建、即桶构建中，在桶中存在多个(例如4个)球的情况下，当将这些球的号码存储在1个桶中并且与像素对应起来时，将向桶的存储处理(向像素的RGBA通道写入球号码的处理)分为用于多次(例如4次)处理的渲染来进行。在用于1次处理的渲染中，按照前文所述的规则，将1个球的号码写入像素的1个通道中。利用GPU12具有的深度测试、颜色蒙板、模板测试的功能，来决定与各渲染中的球号码的写入顺序以及写入通道的选择相关的所述规则。

步骤S13：撞击计算

在该撞击计算中，进行可能与某球i撞击的、存在于附近的球的搜索。在该搜索中，以作为球的坐标纹理而生成的桶纹理为输入，使用片元着色器进行处理。根据该处理，通过参照桶而得到桶围容纳了某球i的桶的3³个桶中容纳的球的号码，通过使用该球的号码来参照存储了球的坐标的纹理，得到该球坐标(球的中心坐标)。使用各球的球坐标来进行这些球和球i的距离计算，进行撞击计算。

步骤S14：撞击响应(施加于球的力和力矩的计算)

关于判断为已撞击的两个球，可以使用两个球的坐标(位置矢量)、它们的速度(速度矢量)和对于刚体重心的相对位置矢量，根据所述(12)式计算施加于球的力、和由于该力而在刚体上产生的转矩。

步骤S15、S16：刚体的平移运动量和角运动量的计算

在步骤S14中计算了施加于球的力和由于该力在刚体上产生的转矩，因此，使用所述(15)式和(16)式将它们相加，计算施加于许多刚体中的各个刚体的力。接着，使用所述(2)式和(4)式，更新刚体的平移运动量和角运动量。

步骤S17、S18：刚体的重心坐标和四元数的计算

使用所述(3)式来计算刚体的重心坐标。在本实施例中，以刚体的重心坐标和重心速度纹理作为输入，对另1块刚体重心坐标纹理进行输出。另外，使用所述(10)式进行刚体的四元数的更新。在本实施例中，进一步以刚体的四元数和角速度纹理作为输入，对另1块刚体四元数纹理写出更新后的值。

步骤S19：渲染处理(用于显示的渲染)

使用如上所述计算出的刚体的重心坐标和四元数，根据在下述中所示的计算内容，进行用于在显示装置15的画面上进行描绘的显示的渲染。

若最初使用四元数Q，使顶点i的距离未旋转时的刚体的中心的相对位置矢量ri旋转当前的刚体旋转量，然后使其移动刚体的重心坐标X的量，则通过(22)式计算顶点i的当前的相对位置矢量r’。

r′i＝X+QriQ^T    (22)

接着，使用模型视图投影矩阵M，根据(23)式，变换为当前正在渲染的坐标系中的坐标r”。

r″i＝Mr′i    (23)

在通过GPU12进行的计算中，基本上通过重复以上8个阶段的算法，可以在显示装置15的显示画面上实时地显示许多刚体的撞击模拟的图像。

此外，图7表示图5所示的刚体模拟的高速化的计算算法的变形例。在图7的处理的流程中，对于与图5中说明的步骤实质上相同的步骤赋予相同符号。根据图7所示的计算算法，可以构成为通过并行的处理，同时执行刚体平移运动量的计算S15和刚体重心坐标的计算S16、刚体角运动量的计算S17和刚体四元数的计算S18。这利用了GPU12中的并行处理。其它处理的流程与图5中说明的内容相同。由此，可以进一步缩短刚体模拟的计算算法中的1个时步的时间，可以实现进一步高速化。

在图3所示的10922个圆环的落下·撞击的刚体模拟图像中，1个时步所需要的时间是16.6毫秒，渲染时的帧频是23FPS。另外，在图4所示的16384个国际象棋的棋子的落下·撞击的刚体模拟图像中，1个时步所需要的计算时间是12.8毫秒，渲染时的帧频是21.2FPS。

在上述的实施例的说明中，关于刚体的方向，使用了“四元数”这样的物理量，但在计算上也可以代替它而使用“旋转矩阵”。“四元数”和“旋转矩阵”都是表示“物体姿势”的物理量。

所述实施例中的粒子法模拟的邻近粒子搜索中所使用的数据结构的构建方法，被应用于与许多刚体的撞击相关的刚体模拟，但本发明的该数据结构的构建方法不限于此，例如也可以如图8和图9所示，应用于流体或粉体的流动行动等自然物表现。

图8是表示多个球体形状的流体(液滴)落入到水等流体中时的流体的变化状态的模拟的结果图像。图9是表示漏斗状容器中容纳的粉体向下方流动，从容器的下端开口部落到外侧，被配置于下侧的容器容纳的变化状态的模拟的结果图像。

通过GPU12执行图1所示的邻近粒子搜索程序13A-1，由此实现所述邻近粒子搜索中所使用的数据结构的构建方法。邻近粒子搜索中所使用的数据结构的构建方法涉及到的邻近粒子搜索程序13A-1，其自身可以被存储在CDROM等存储介质上，可以构成独立的程序产品。

在所述的实施例的说明中，说明了在通过GPU的程序实现中使用着色器(Shader)的情况。这是使用了专门用于GPU的图形功能的程序模型的例子。但本发明不限于此。作为图形运算处理器，可以使用包含上述GPU在内的一般的流处理器，而且，作为编程模型，可以使用一般的编程模型来实现。

在包含上述GPU12在内的流处理器(Streaming Processor)中，一般除了上述着色器(Shader)以外，还可以使用统一计算设备架构(CUDA)、CTM、Brook等编程模型来实现。在所述实施例的GPU中的着色器(Shader)的实现中，将用于GPU12上的计算的计算空间中的切片(slice)、体素(voxel)以及与粒子相关的数据作为纹理13B-1而存储在视频存储器13B中。在作为程序实现而使用了统一计算设备架构(CUDA)等的情况下，代替纹理而使用一般的“阵列(Array)”。

如前所述，在数据结构的构建方法中作为图形运算处理器而使用GPU的情况下，使用由GPU准备的深度测试、颜色蒙板、模板测试的各功能，来进行用于多次的处理器的渲染的每个渲染中写入的粒子号码和粒子坐标的选择、以及被写入的像素中的通道的选择。

另一方面，在作为图形运算处理器而使用流处理器的情况下，一般代替上述像素(written pixel)而使用存储器位置(Memory location)，代替上述通道而使用用于体素(Voxel)的存储器中的存储器位置，而且，代替由GPU准备的上述深度测试、颜色蒙板、模板测试的各功能而使用由流处理器准备的计数器的功能。该计数器代替上述深度测试、颜色蒙板、模板测试而工作。

因此，在作为图形运算处理器而使用了流处理器的情况下，使用由流处理器准备的计数器的功能，进行多次处理器周期的各步骤中写入的粒子的粒子号 码和存储器位置(Memory Location)的选择、以及用于体素(Voxel)的存储器中的存储器位置的选择。

如前所述，数据结构的构建方法、即图6所示的桶构建S12的工序，在作为图形运算处理器而使用GPU12进行程序实现的情况下，具有如下步骤：当假定4个粒子的粒子号码是i_b0、i_b1、i_b2、i_b3(i_b0＜i_b1＜i_b2＜i_b3)时，在向第一通路的像素进行的写入中，以最大值对深度缓冲区进行初始化，使用深度测试使具有较小值的深度缓冲区合格，由此，在R通道中写入粒子号码i_b0，进行用于第一次处理的渲染的步骤(S101)；在向第二通路的像素进行的写入中，使用颜色蒙板来禁止向R通道的写入并且使用了模板测试，由此，在G通道中写入粒子号码i_b1，进行用于第二次处理的渲染的步骤(S102)；在向第三通路的像素进行的写入中，使用颜色蒙板来禁止向R和G的各通道的写入并且使用了模板测试，由此，在B通道中写入粒子号码i_b2，进行用于第三次处理的渲染的步骤(S103)；以及在向第四通路的像素进行的写入中，使用颜色蒙板来禁止向R、G和B的各通道的写入并且使用了模板测试，由此，在A通道中写入粒子号码i_b3，进行用于第四次处理的渲染的步骤(S104)。

与上述桶构建S12的工序相关的处理的意思，一般而言，是这样的意思，亦即，在通过流处理器写入最大4个值时分为4次的步骤(S101～S104)来进行写入。在4次的步骤的每一个步骤中，写入1个此前未写入的值。

在安装了使用上述GPU12的程序的情况下，使用了GPU12的图形特有功能“深度测试”、“颜色蒙板”、“模板测试”。

称之为“深度测试”的功能，进行“排除此前写入的值”这样的处理。称之为“颜色蒙板”的功能，进行“选择写入目的地”这样的处理。称之为“模板测试”的功能，进行“一个步骤中仅许可写入一个值”这样的处理。

此外，除了“在一个步骤中，在数据的写入目的地的存储器中写入粒子号码”这样的处理以外，还更新“深度缓冲区”和“模板缓冲区”。

在使用包含上述GPU12在内的流处理器的情况下，一般在程序实现中可以使用统一计算设备架构(CUDA)这种程序模型。在使用统一计算设备架构(CUDA)的情况下，无法使用GPU12的图形功能、即“深度测试”、“颜色蒙板”、“模板测试”。因此，在使用统一计算设备架构(CUDA)的情况下， 为了作成与这些图形功能等同的功能，要准备写入目的地存储器和计数器。

根据上述计数器，构成为：通过计数器进行“排除此前所写入的值”的处理(“深度测试”的功能)，而且进行“选择写入目的地”的处理(“颜色蒙板”的功能)和“一个步骤中仅许可写入一个值”的处理(“模板测试”的功能)。

而且，除了“在一个步骤中，在数据的写入目的地的存储器中写入粒子号码”的处理以外，还更新上述“计数器”。在此，“计数器的更新”意味着将计数器增加1。

如上所述，当使用GPU等流处理器且作为编程模型而实现统一计算设备架构时，代替纹理而使用了在视频存储器中生成的阵列。此外，纹理是被包含在“在视频存储器中生成的阵列”的概念。“在视频存储器中生成的阵列”一般可以解释为“在专用存储器中生成的阵列”。

在GPU等流处理器中，在作为编程模型代替着色器使用统一计算设备架构来实现的情况下，代替像素(written pixel)而使用存储器位置(MemoryLocation)，代替通道而使用用于体素(Voxel)的存储器中的存储器位置，而且代替深度测试、颜色蒙板、模板测试的各功能而使用与之分别对应的计数器的功能。

Claims (10)

1.一种使用邻近粒子搜索的物理对象物的模拟方法，其由图形运算处理器GPU执行，被用于通过将3维数字图像数据变换为2维数字图像数据来对物理对象物的模拟进行运算的方法，所述3维数字图像数据是根据在视频存储器内部生成的图像数据而构建的，所述物理对象物的模拟方法的特征在于，

使用了下述方法，即通过许多粒子来表现所述物理对象物，根据所述视频存储器内的3维数字图像数据结构中的所述许多粒子的相互的容纳位置关系，决定所述物理对象物的存在状态，

所述3维数字图像数据结构中包含的桶纹理，由将对应于出现所述物理对象物的空间而生成的存储器空间进行分割的桶构成，将分配给所述许多粒子的粒子号码存储在许多所述桶的每一个中，

具有写入步骤：其通过所述GPU使用顶点着色器，读出在所述桶中存储的粒子的粒子号码和粒子坐标，写入2维图像空间，

在所述写入步骤中，当所述桶中存在多个粒子时，在模拟的一个时步中划分为用于多次处理的渲染，一个一个地写入所述多个粒子的各个粒子的粒子号码和粒子坐标。

2.根据权利要求1所述的方法，其中，

使用由所述GPU准备的深度测试、颜色蒙板、模板测试的各功能，来进行在用于所述多次处理的渲染的各个渲染中写入的粒子的粒子号码和写入的像素的选择、以及被写入的像素中的通道的选择。

3.根据权利要求2所述的方法，其中，

当所述桶中存在4个粒子时，在由写入粒子号码的像素的RGBA组成的4个通道的每一个中，独立地从按照升序或降序的粒子号码的粒子中顺次存储粒子号码。

4.根据权利要求3所述的方法，其中，

具有如下步骤：

当假定4个粒子的粒子号码是i_b0、i_b1、i_b2、i_b3时，

在向第一通路的所述像素进行的写入中，以最大值对深度缓冲区进行初始化，使用深度测试使具有较小值的深度缓冲区合格，由此，在R通道中写入粒子号码i_b0，进行用于第一次处理的渲染；

在向第二通路的所述像素进行的写入中，使用颜色蒙板来禁止向所述R通道的写入并且使用模板测试，由此，在G通道中写入粒子号码i_b1，进行用于第二次处理的渲染；

在向第三通路的所述像素进行的写入中，使用颜色蒙板来禁止向所述R和G的各通道的写入并且使用模板测试，由此，在B通道中写入粒子号码i_b2，进行用于第三次处理的渲染；以及

在向第四通路的所述像素进行的写入中，使用颜色蒙板来禁止向所述R、G和B的各通道的写入并且使用模板测试，由此，在A通道中写入粒子号码i_b3，进行用于第四次处理的渲染。

5.根据权利要求1所述的方法，其中，

所述物理对象物是许多刚体，

所述物理对象物的模拟方法，执行基于所述许多刚体的撞击计算的模拟的运算。

6.一种使用邻近粒子搜索的物理对象物的模拟方法，其由流处理器执行，被用于通过将3维数字图像数据变换为2维数字图像数据来对物理对象物的模拟进行运算的方法，所述3维数字图像数据是根据在存储器内部生成的数据阵列而构建的，所述物理对象物的模拟方法的特征在于，

使用了下述方法，即通过许多粒子来表现所述物理对象物，根据所述存储器内的3维数字图像数据结构中的所述许多粒子的相互的容纳位置关系，决定所述物理对象物的存在状态，

所述3维数字图像数据结构中包含的桶阵列，由将对应于出现所述物理对象物的空间而生成的存储器空间进行分割的桶构成，将分配给所述许多粒子的粒子号码存储在许多所述桶的每一个中，

具有写入步骤：其通过所述流处理器，使用该流处理器的分散操作，读出在所述桶中存储的粒子的粒子号码和粒子坐标，并写入所述存储器中的阵列中，

在所述写入步骤中，当所述桶中存在多个粒子时，在模拟的一个时步中划分为用于多次处理的内核，一个一个地写入所述多个粒子的各个粒子的粒子号码和粒子坐标。

7.根据权利要求6所述的方法，其中，

使用由所述流处理器准备的计数器的功能，来进行在用于所述多次处理的各步骤中写入的粒子的粒子号码和写入数据的选择、以及为了体素而分配的存储器中的多个存储器位置的选择。

8.根据权利要求7所述的方法，其中，

当所述桶中存在4个粒子时，在写入粒子号码的第1至第4所述存储器位置的每一个中，独立地从按照升序或降序的粒子号码的粒子中顺次存储粒子号码。

9.根据权利要求7所述的方法，其中，

具有如下步骤：

当假定4个粒子的粒子号码是i_b0、i_b1、i_b2、i_b3时，

在向第一通路中的为了所述体素而分配的所述存储器的写入中，通过初始化成所述计数器为0并且接受向第一存储器位置的写入，在所述第一存储器位置中写入粒子号码i_b0，执行用于第一次的例程的内核；

在向第二通路中的为了所述体素而分配的所述存储器的写入中，通过使用所述计数器来禁止向所述第一存储器位置的写入并且向第二存储器位置写入下一值，由此，在所述第二存储器位置中写入粒子号码i_b1，执行用于第二次的例程的内核；

在向第三通路中的为了所述体素而分配的所述存储器的写入中，通过使用所述计数器来禁止向所述第一和第二存储器位置的写入并且向第三存储器位置写入下一值，由此，在所述第三存储器位置中写入粒子号码i_b2，执行用于第三次的例程的内核；以及

在向第四通路中的为了所述体素而分配的所述存储器的写入中，通过使用所述计数器来禁止向所述第一、第二和第三存储器位置的写入并且向第四存储器位置写入下一值，由此，在所述第四存储器位置中写入粒子号码i_b3，执行用于第四次的例程的内核。

10.根据权利要求6所述的方法，其中，

所述物理对象物是许多刚体，

所述物理对象物的模拟方法，执行基于所述许多刚体的撞击计算的模拟的运算。

Images