一种间断伽辽金法求解欧拉方程的GPU加速方法技术

本发明专利技术属于计算流体力学、高性能计算领域，涉及一种间断伽辽金方法(DGM)的GPU并行加速技术，具体为一种间断伽辽金法求解欧拉方程的GPU加速方法。本发明专利技术采用四面体网格对求解区域进行剖分，以基函数、高斯积分、数值通量为基础，GPU为主要计算硬件，CUDA为编程模型建立间断伽辽金方法GPU并行框架。通过CUDA并行框架实现了GPU多线程的管理，通过设计的数据结构和线程访问方式来实现高效的内存访问。在解决面积分需要单元间数据交换而不独立的问题时，采用了按照面网格并行、每个面网格的计算线程处理两个单元的方式，既避开了单元不独立的问题，实现了大规模并行，还减少了计算量。

GPU acceleration method for solving Euler equations by discontinuous Galerkin method

The invention belongs to the field of computational fluid dynamics and high-performance computing, and relates to a GPU parallel acceleration technology of discontinuous Galerkin method (DGM), in particular to a GPU acceleration method of discontinuous Galerkin method for solving Euler equation. The method adopts tetrahedral grid to divide the solution area, and takes the basis function, Gauss integral and numerical flux as the basis, GPU as the main calculation hardware, CUDA as the programming model to establish the GPU parallel framework of discontinuous Galerkin method. GPU multithread management is realized by CUDA parallel framework, and efficient memory access is realized by the designed data structure and thread access mode. In order to solve the problem that area division needs data exchange between cells and is not independent, the method of processing two cells according to the parallel computing thread of face grid and each face grid is adopted, which not only avoids the problem of cell independence, but also realizes large-scale parallel, and reduces the calculation amount.

【技术实现步骤摘要】
一种间断伽辽金法求解欧拉方程的GPU加速方法
本专利技术属于计算流体力学、高性能计算领域，涉及一种间断伽辽金方法(DGM)的GPU并行加速技术，具体为一种间断伽辽金法求解欧拉方程的GPU加速方法。
技术介绍
欧拉方程是流体力学中描述无粘流体的方程组，其形式如下：Ut+▽·F＝0(1)其中U代表守恒量、Ut代表守恒量对时间t的偏导数，F代表守恒通量，▽·F代表守恒通量的散度，且在三维情况下，有其中ρ为气体密度，u、v、w为气体的三个速度分量，e为完全气体的单位体积总能量，p为气体压强。对于上述欧拉方程的数值求解，通常采用以下几种方法：有限差分法、有限体积法、有限元法。其中有限差分法需要采用结构网格，且计算量小，常用于处理结构网格划分的简单几何区域上的求解，而对于复杂几何区域的求解则相对困难。有限体积法可以求解结构或非结构网格，因此可以处理复杂集合区域，应用范围相对较广，但其难以构造高阶格式(需要扩展模板)，且构造的高精度格式要么求解复杂、要么不够紧致。而间断伽辽金方法(DGM)则结合了有限元和有限体积方法，能够处理任意网格和复杂几何区域，而且DGM可以通过简单地增加单元内的解多项式的次数进而增加单元自由度(DOFs)来获得更高的空间精度，是一种高精度的流场求解方法。间断伽辽金法的数值求解过程如下：将欧拉方程的两端乘以试探函数Φ并在体积Ω上积分，经积分变换可以得到伽辽金弱形式将积分区域Ω划分为四面体网格，并取其中的一个单元Ωk来考察上述方程，为了保证单元之间场的连续性，上式的面积分项的被积函数需要使用“数值通量”F*来代替，数值通量由单元两侧的数值计算得到，于是可以将上式重写为其中u+,u-分别代表积分面两侧的数值。若将守恒量u用基函数φi展开，并且试探函数也为φi，可得到其中左端项仅与基函数有关，使用正交的基函数可以得到一个对角的质量矩阵，且仅与本单元相关。右端第一项为体积分项，仅与本单元项相关。以上两项是DG中最直接具有并行性的部分。右端第二项为面积分项，采用数值通量之后，这一项与积分面两侧的场都有关。虽然DGM有更高的精度，也存在的一个计算量较大的问题，往往需要比其他方法更长的计算时间。但单元DOFs在单元间相对独立的特点使得该方法具有天然的并行性，非常适合于大规模并行计算，利用这一点可以弥补其计算量大的问题。在并行计算技术方面，当前主要有多核中央处理器(CPU)并行、图形处理器(GPU)与CPU的异构并行两类。其中，多核CPU并行出现最早，当前超算采用最多的架构就是多核CPU架构，通过增加核心数来增加并行任务的数量。这一方法的缺点也很明显，当前单个的CPU线程数只有最多64个，而再需求更多核心则需要增加CPU个数，这会大幅度提高成本，且对于小型工作站而言CPU数目限制比较大，难以实现较大规模的计算。对于GPU并行，以英伟达TITANV为例，其包含80个流处理器簇(SM)，每个流处理器簇又包含64个流处理器(SP)，在英伟达CUDA编程模型下，可以同时处理海量的线程以实现大规模并行计算。当前的英伟达GPU设备的每秒浮点运算次数已达到15.7TFLOPS，是同时期英特尔CPU的10倍，使用一个GPU工作站将能够取代十台CPU工作站，且成本更低。GPU-CPU异构架构则是一种高效的并行计算架构，以当前世界超算天梯榜首Summit为代表的一系列超算就是采用使用这样的架构。异构计算平台的CPU和GPU的数量比较灵活，最简单的应用仅需要一个CPU与一个GPU，即可实现大规模的并行计算。如果能将GPU并行计算应用于DGM可以很好地弥补DGM计算量偏大的问题，但其中存在若干待解决的问题：在技术上，传统的CPU运行的DGM通常依赖于现成且优化过的线性代数及通信原语，而这些在GPU上都是不可用的；在实施上，虽然半离散格式(5)是基于一个单元展开的，但还是通过右端第二项与其它相邻单元进行了数据交换，也就是这一部分并非是单元独立的，需要通过一定的实施手段来保证其能够并行且效率足够高。
技术实现思路
针对上述存在问题或不足，为解决现有欧拉方程的数值求解DGM结合多核CPU的方法成本相对较高，并且相对效率低下的问题，本专利技术提供了一种间断伽辽金法求解欧拉方程的GPU加速方法，基于间断伽辽金法的特点和GPU并行计算的优势，将GPU并行计算应用于间断伽辽金法求解欧拉方程,能最大程度够发挥间断伽辽金法的优势，加速计算流体力学高精度格式的计算。该间断伽辽金法求解欧拉方程的GPU加速方法的具体技术方案包括以下几个步骤：步骤1、读取计算网格的信息，并将面网格按照边界类型排序；步骤2、在CPU端进行几何参数的预计算，并将结果拷贝到GPU显存；预计算并存储几何参数，是因为几何参数在计算全程不发生改变，只需计算一次即可，且几何参数的计算需要访问非连续的内存，无法实现对齐、合并的访问，若在GPU端的DGM计算过程中频繁计算，将大幅降低计算效率。步骤3、在GPU端完成流场的初始化，并且将时间步初始化为0；初始化时分配与单元数量一致的线程数，每个线程负责初始化单元内的所有场。步骤4、判断计算时间步是否达到预定的终止时间步，若达到，结束计算，否则继续；步骤5、在GPU端进入时间推进，具体为：依次启动面上场计算核函数、面积分核函数、体积分核函数、以及时间更新核函数。步骤6、再次判断当前计算时间步是否达到预定的后处理时间步，若达到，同步GPU线程，并将计算结果拷贝到CPU端，将当前计算得到的流场数据输出；否则不进行任何操作。步骤7、时间步自增1，转至步骤4。上述步骤中GPU显存不存储原始的网格信息(如网格结点坐标)而只存储：计算直接需要的且由原始数据计算得到的数据(法向量、体积和面积)，网格拓扑关系(如一个单元的每一个面的全局编号等)，以及流场数据(场、右端项)，并使用数组结构体的方式组织并存储在GPU全局内存中，便于管理与使用，而CPU端仅需在预计算时保留网格原始信息和流场数据，在预计算完成后即可释放，仅需要保留计算的场数据。对于高斯积分点的基函数值、高斯积分权重、常数质量矩阵三类所有线程都会使用到的常数，存储于GPU的常量内存中。进一步的，所述步骤1中对面网格进行的排序时，以网格编号为关键字，采用桶排序算法进行排序。进一步的，所述步骤3中使用由文件中导入的流场进行初始化。进一步的，所述步骤5中，对于所有计算核函数的线程分配，为了保证对齐与合并访问来提高内存访问效率，需要将4个单元或者面划分到一个线程块内，并在同一个线程束内，且对内存中数组的访问正好实现128字节的对齐以及合并访问，线程结构及对内存访问方式如图2所示。进一步的，所述步骤5中面积分核函数，即式(5)中右端第二项的数值通量，采用计算过程中无需判断的Roe通量格式，以减少线程束分化。进一步的，所述步骤5中，面积分核函数与体积分核函数的计算相互独立，使用CUDA流，将面上场计算核函数、面积分核函数发本文档来自技高网...

【技术保护点】
1.一种间断伽辽金法求解欧拉方程的GPU加速方法，包括以下几个步骤：/n步骤1、读取计算网格的信息，并将面网格按照边界类型排序；/n步骤2、在CPU端进行几何参数的预计算，并将结果拷贝到GPU显存；/n步骤3、在GPU端完成流场的初始化，并且将时间步初始化为0；初始化时分配与单元数量一致的线程数，每个线程负责初始化单元内的所有场；/n步骤4、判断计算时间步是否达到预定的终止时间步，若达到，结束计算，否则继续；/n步骤5、在GPU端进入时间推进，具体为：依次启动面上场计算核函数、面积分核函数、体积分核函数、以及时间更新核函数。/n步骤6、再次判断当前计算时间步是否达到预定的后处理时间步，若达到，同步GPU线程，并将计算结果拷贝到CPU端，将当前计算得到的流场数据输出；否则不进行任何操作。/n步骤7、时间步自增1，转至步骤4。/n上述步骤中GPU显存不存储原始的网格信息而只存储：计算直接需要的且由原始数据计算得到的数据、网格拓扑关系以及流场数据，使用数组结构体的方式组织并存储在GPU全局内存中，而CPU端仅在预计算时保留网格原始信息和流场数据，在预计算完成后释放，仅保留计算的场数据；对于高斯积分点的基函数值、高斯积分权重、常数质量矩阵三类所有线程都会使用到的常数，存储于GPU的常量内存中。/n所述计算直接需要的且由原始数据计算得到的数据为法向量、体积和面积。/n...

【技术特征摘要】
1.一种间断伽辽金法求解欧拉方程的GPU加速方法，包括以下几个步骤：
步骤1、读取计算网格的信息，并将面网格按照边界类型排序；
步骤2、在CPU端进行几何参数的预计算，并将结果拷贝到GPU显存；
步骤3、在GPU端完成流场的初始化，并且将时间步初始化为0；初始化时分配与单元数量一致的线程数，每个线程负责初始化单元内的所有场；
步骤4、判断计算时间步是否达到预定的终止时间步，若达到，结束计算，否则继续；
步骤5、在GPU端进入时间推进，具体为：依次启动面上场计算核函数、面积分核函数、体积分核函数、以及时间更新核函数。
步骤6、再次判断当前计算时间步是否达到预定的后处理时间步，若达到，同步GPU线程，并将计算结果拷贝到CPU端，将当前计算得到的流场数据输出；否则不进行任何操作。
步骤7、时间步自增1，转至步骤4。
上述步骤中GPU显存不存储原始的网格信息而只存储：计算直接需要的且由原始数据计算得到的数据、网格拓扑关系以及流场数据，使用数组结构体的方式组织并存储在GPU全局内存中，而CPU端仅在预计算时保留网格原始信息和流场数据，在预计算完成后释放，仅保留计算的场数据；对于高斯积分点的基函数值、高斯积分权重、常数质量矩阵三类所有线程都会使用到的常数，存储于GPU的常量...

【专利技术属性】
技术研发人员：徐立，王晗，杨中海，李斌，
申请(专利权)人：电子科技大学，
类型：发明
国别省市：四川;51