面向SPH流体模拟的GPGPU并行计算方法技术

本发明专利技术公开了一种面向SPH流体模拟的GPGPU并行计算方法。该方法将SPH方法中粒子的计算任务进行了分组来迎合GPGPU中的线程层次结构，同一个协作线程队列（CTA）中的线程能够处理逻辑上相似的工作。这个分组方法针对了GPGPU中的CTA调度规则开发并在GPGPU上实现，从而能快速分配任务的同时增加系统的cache命中率。以此分配方法为基础，片内共享内存能够被用来缓存来自全局内存中读入的数据。并且本发明专利技术能够根据全局内存的存储访问特性充分存储器带宽。同时，通过对 GPGPU的SIMD特性的利用，由于使用共享内存所带来的线程同步开销也能被有效规避。

【技术实现步骤摘要】

本专利技术属于计算机图形学与高性能计算领域，具体地说是一个新的基于通用计算图形处理器(GPGPU)并行方法，用于以SPH算法为基础的流体模拟计算。涉及到计算机图形学以及仿真，计算流体动力学，单指令流多数据流(SIMD)架构，GPGPU等。
技术介绍
SPH方法是一种无网格拉格朗日数值方法。该方法在近年来被广泛运用于流体模拟、天体动力学、分子动力学以及生物模拟中。作为一个粒子方法，SPH需要大量的粒子来模拟高精度大规模的场景。而大量的高密度粒子会导致海量的计算需求。所以对于用户来说，他们总是被迫于在模拟效果和系统性能之间做出权衡。因此，提升SPH方法的性能总是非常有意义的。如今，非常多的研究致力于提升SPH方法的性能。其中一个主题是从系统体系结构的角度来做性能的提升。由于SPH方法本身非常适合并行化，非常多的研究得益于近年来并行硬件以及开发环境的发展并成功地将SPH方法和它们结合起来。然而，大多数这些研究，尤其是那些利用GPU的研究工作，没有能够充分地利用这些新硬件和开发环境的能力。比如说，主流的传统SPH的GPGPU实现中，尽管这种实现非常易于实现，并且相比起CPU代码来说，有1到2个数量级的性能提升，但是GPU在这种实现中只是被单纯地看作是一个SIMT处理器，而没有能够运用到很多GPU中的高级特性。而这种做法导致了传统的实现让SPH成为了内存受限型的应用。于是，也有一些研究试图去利用GPU中的各种内存空间来改进其在内存访问方面的高额代价，特别是共享内存。然而，已有的共享内存解决方案并没有带来理想中的性能提升，甚至要比传统方法慢上20％～70％。这个解决方案使用了显式的CTA中线程同步的方法来保证数据的正确性，但是线程同步带来的非常高昂的性能开销；同时该方法过多地占用共享内存使得其在大多数设备上运行时有着较低的设备占用率；其次，其低效的全局内存访问模式没有很好地利用全局内存的访问特性；而其随机的任务分配策略会导致较低的cache命中率。总而言之，尽管利用GPGPU设备的高并行处理能力，SPH的性能得到了显著的提升，但是其潜力仍然远远没有被充分挖掘。
技术实现思路
本专利技术的目的是提供一个新的GPU并行计算架构来加速以SPH方法为基础的流体模拟在GPGPU设备上的性能。该方法包含了一个新的面向SPH流体模拟的GPGPU并行计算方法。本专利技术的目的是这样实现的：一种面向SPH流体模拟的GPGPU并行计算方法，该方法包括以下具体步骤：步骤1：设SPH方法中的光滑核半径为h，将整个模拟空间划分为以h为边长的互不相交的正方体单元格；步骤2：使用基于GPGPU的计数排序算法，将空间内的粒子根据它们所在的单元格来进行排序，在计数排序算法中，首先计算每个单元格中所包含的粒子数量，然后通过并行前缀求和算法计算每个单元格中的粒子在内存中的偏移量，之后根据每个粒子所在单元格在内存中的偏移量以及其在该单元格中的位置进行复制；步骤3：根据每个单元格中粒子的数量，为接下来的粒子交互阶段中每一个CTA分配任务，此过程完全在GPGPU上实现：i)计算每个单元格所包含的粒子在内存中的范围，其结果保存在两个数组中，数组中每个元素描述一个单元格的相关信息：(1)单元格中的粒子在内存中的偏移量；(2)单元格所包含的粒子数量B；这步由个GPU线程并行完成，为空间中所包含的单元格的数量；ii)计算每一个单元格将在步骤4中所需要用来处理其内部粒子的交互的CTA数量，并将结果保存在数组B’中；所需的CTA数量为其中，Scta表示一个CTA中包含了线程数量，表示单元格cj中所包含的粒子的数量，这些数据保存在数组B中；这步由个GPU线程并行完成；iii)对数组B’中的元素做前缀求和并将结果保存在数组O’中；O’中的元素表示每个单元格所需的任务在任务队列T中的起始位置；这步由个GPU线程并行完成；iv)T中每个元素ti＝{ci，ki本文档来自技高网...

【技术保护点】
一种面向SPH流体模拟的GPGPU并行计算方法，其特征在于该方法包括以下具体步骤：步骤1：设SPH方法中的光滑核半径为h，将整个模拟空间划分为以h为边长的互不相交的正方体单元格；步骤2：使用基于GPGPU的计数排序算法，将空间内的粒子根据它们所在的单元格来进行排序，在计数排序算法中，首先计算每个单元格中所包含的粒子数量，然后通过并行前缀求和算法计算每个单元格中的粒子在内存中的偏移量，之后根据每个粒子所在单元格在内存中的偏移量以及其在该单元格中的位置进行复制；步骤3：根据每个单元格中粒子的数量，为接下来的粒子交互阶段中每一个CTA分配任务，此过程完全在GPGPU上实现：i)计算每个单元格所包含的粒子在内存中的范围，其结果保存在两个数组中，数组中每个元素描述一个单元格的相关信息：(1)单元格中的粒子在内存中的偏移量；(2)单元格所包含的粒子数量B；这步由个GPU线程并行完成，为空间中所包含的单元格的数量；ii)计算每一个单元格将在步骤4中所需要用来处理其内部粒子的交互的CTA数量，并将结果保存在数组B’中；所需的CTA数量为其中，Scta表示一个CTA中包含了线程数量，表示单元格cj中所包含的粒子的数量，这些数据保存在数组B中；这步由个GPU线程并行完成；iii)对数组B’中的元素做前缀求和并将结果保存在数组O’中；O’中的元素表示每个单元格所需的任务在任务队列T中的起始位置；这步由个GPU线程并行完成；iv)T中每个元素ti＝{ci，ki}，其中ci表示当前这个任务需要处理哪个单元格中的粒子，而ki代表自己是第几个处理该单元格的任务；这步由个GPU线程并行完成；v)获取总共需要的CTA数量Sgrid；该值通过将O’中最后一个值与最后一个非空单元格所需的CTA数量相加得到；步骤4：进行粒子交互；粒子交互包括两个独立的遍历邻居粒子阶段：读取邻居粒子的位置与质量来计算每个离散粒子表示的液体密度；通过邻居粒子的密度信息来计算每个离散粒子的受力；这两个阶段遍历邻居粒子的实现流程是相同的，具体步骤如下：i)每个线程将自己所负责的粒子的必要数据从全局内存中读取到寄存器中；如果该线程不负责任何粒子，则跳过此步骤；ii)一个CTA中前9个线程计算它们所负责的那个单元格周围的单元格的粒子在内存中的保存范围，结果保存在共享内存中；周围的单元格的粒子在内存中的范围是9个连续区域即3D场景下；其余线程跳过此步骤；iii)读取一个邻居粒子区域中最多Scta个粒子，保存到共享内存中；CTA中每个线程最多负责读取一个；iv)凡是有粒子负责的线程通过在寄存器中自己负责的粒子数据以及共享内存中最多Scta个邻居粒子的数据来计算自己负责的粒子的密度或是受力；v)Scta设为32，也就是一个线程束中所包含的线程的数量；因此内存访问与计算之间不需要显式的线程同步；步骤5：根据步骤4中计算出的每个粒子的受力情况，来移动粒子；移动粒子的时候要考虑外力影响与边界情况；当粒子在边界上时给予一个反向的力将粒子弹回模拟空间内；步骤6：一个时间步中所需要的计算就完成了；之后，输出粒子的位置信息，然后计算下一个时间步中粒子的变化，回到步骤1。...

【技术特征摘要】
1.一种面向SPH流体模拟的GPGPU并行计算方法，其特征在于该方法包括以下具体步骤：步骤1：设SPH方法中的光滑核半径为h，将整个模拟空间划分为以h为边长的互不相交的正方体单元格；步骤2：使用基于GPGPU的计数排序算法，将空间内的粒子根据它们所在的单元格来进行排序，在计数排序算法中，首先计算每个单元格中所包含的粒子数量，然后通过并行前缀求和算法计算每个单元格中的粒子在内存中的偏移量，之后根据每个粒子所在单元格在内存中的偏移量以及其在该单元格中的位置进行复制；步骤3：根据每个单元格中粒子的数量，为接下来的粒子交互阶段中每一个CTA分配任务，此过程完全在GPGPU上实现：i)计算每个单元格所包含的粒子在内存中的范围，其结果...

【专利技术属性】
技术研发人员：阮骥鸣，王长波，秦洪，
申请(专利权)人：华东师范大学，
类型：发明
国别省市：上海;31