CN102819454A - 基于gpu的有限元显式并行求解仿真方法 - Google Patents

Abstract

本发明公开了一种基于GPU的有限元显式并行求解仿真方法，实施步骤如下：1）将数据读入主机内存中，通过主机CPU识别有限元计算模型数据中的节点和单元，查找并保存节点和单元之间的归属关系，将所述有限元计算模型数据读入GPU全局存储器；2）设定GPU执行参数；3）通过GPU读取GPU全局存储器中存储的数据，根据GPU执行参数通过多线程的方法计算各个节点的节点作用力；4）选择性输出计算结果，CPU根据收敛判断情况决定退出或者返回执行步骤3）。本发明具有运算速度快、计算成本低、硬件成本低，使用简单方便，尤其适用于板料成形仿真优化、汽车碰撞仿真等力学领域的动态非线性计算。

Description

基于GPU的有限元显式并行求解仿真方法

技术领域

本发明涉及计算机辅助工程中的有限元仿真领域，具体涉及一种基于GPU的有限元显式并行求解仿真方法。

背景技术

显式有限元计算方法不需要进行任何方程组的求解，具有很好的稳定性和收敛性，目前已经被广泛应用于汽车碰撞和金属薄板成形等强非线性问题的仿真计算。随着现代工程问题复杂度的提升，有限元计算模型规模越来越大，对于条件稳定的显式有限元算法，计算时间步长一般在微秒级，而实际仿真时间一般达到秒级，这意味着巨大的计算量。

随着计算技术的发展，并行计算方法是解决计算效率问题的有效手段，并且已经在多个领域等到了成功的应用。现有技术的并行有限元计算方法，在技术层面上主要采用分布式计算、并行机或多线程等并行处理技术，所使用的计算硬件主要是以CPU为计算核心的并行计算平台，在取得较好加速比的同时还存在以下的下足：一方面，对于大规模的有限元模型，庞大的单元、节点规模导致分布式计算机之间或多线程之间的通讯消耗提升，导致加速比难以提升；另一方面，并行机价格昂贵，使用和维护相对复杂，对研究人员算法改写和编程能力要求也较高，难于普及。此外，目前基于GPU和CUDA的通用计算平台，硬件成本低，软件实现容易，计算效率提升明显，目前已成为国内外一个新的研究热点，目前在力学领域，显式有限元计算仍然不能利用基于GPU和CUDA的通用计算平台来实现并行求解仿真。

发明内容

本发明要解决的技术问题是提供一种运算速度快、计算成本低、硬件成本低，使用简单方便，尤其适用于板料成形优化、汽车碰撞等力学领域的非线性计算的基于GPU的有限元显式并行求解仿真方法。

为了解决上述技术问题，本发明采用的技术方案为：

一种基于GPU的有限元显式并行求解仿真方法，其实施步骤如下：

1）将有限元计算模型数据读入主机内存中，通过主机CPU识别所述有限元计算模型数据中的节点和单元，找到每一个节点归属的单元并保存节点和单元之间的归属关系，将所述有限元计算模型数据读入GPU全局存储器；

2）设定GPU执行的线程块数量和每一个线程块所包含的线程数量；

3）通过GPU读取所述GPU全局存储器中存储的有限元计算模型数据，并根据线程块数量和每一个线程块所包含的线程数量，通过多线程的方法计算单元的最小时间步长、各个节点的节点加速度以及各个节点的节点作用力；

4）判断是否需要输出计算结果，如果需要则采用异步的方式将GPU全局存储器上的计算结果复制到主机的内存中，并写入到输出文件；通过CPU判断总仿真时间或者总迭代次数是否达到设定值，如果达到设定值则退出；否则返回执行步骤3）。

所述步骤3）的详细步骤包括：

3.1）通过GPU读取GPU全局存储器中存储的所有单元数据，然后采用线程与单元一一对应的执行方式计算出每个单元对应的时间步长，并通过并行缩减的方法在所述时间步长中查找最小值作为最小时间步长；

3.2）通过GPU读取GPU全局存储器中存储的所有节点数据，然后采用线程与节点一一对应的执行方式计算每个所述节点所对应的节点加速度；

3.3）通过GPU根据约束边界条件以及所述最小时间步长、节点加速度采用线程与节点对应的方法计算出每个节点新的速度和位移；

3.4）通过GPU读取GPU全局存储器中存储的所有节点数据和单元数据，然后GPU根据有限元计算模型的力边界条件计算得到各个节点的节点外力，通过并行运算的方法分别计算各个单元的单元内力，采用线程与节点一一对应的执行方式将单元内力离散到整体节点内力中，采用线程与节点一一对应的执行方式将节点外力和整体节点内力做差值运算得到各个节点的节点作用力。

本发明具有下述优点：

1、本发明的并行计算基于GPU通用计算平台实现，具有计算成本低，通用性强的优点。

2、本发明在显式有限元计算的全部迭代过程均移植到GPU上进行计算，数据输出采用异步执行的方式，从而有效避免了由于主机与GPU的数据交换所产生的额外时间。

3、本发明能够充分利用GPU计算特性，保证了程序的计算效率，能够利用现有的NVIDIA公司所提供的免费数值计算库CUBLAS，保证了程序代码的简洁性和可操作性。

4、本实施例通过在程序开始阶段查找到每个节点所归属的单元，并记录到相关数据中，从而避免了单元作用力与离散节点作用力操作在GPU上并行执行时所产生的并行竞写错误。

5、本发明进一步在采用GPU进行计算时，根据所计算对象的不同，分别采用不同的线程映射模式，主要包括两种映射模式：单元与线程的一一映射，节点与线程的一一映射。

6、本发明通过建立有限元模型和GPU线程之间的映射关系，针对非线性问题，在保证求解精度的前提下，能够大幅度提高显式有限元算法的计算效率。同传统的显式有限元并行计算方法不同，本发明通过建立有限元计算模型中的单元或节点与GPU线程间的一一对应关系，利用GPU所提供大规模计算线程，在计算流程中由一个线程负责一个单元或节点的计算，从而能够在包含成百上千个处理核心的GPU上并行、高效的实现显式有限元计算以及在工程中应用，具有高效、计算成本低的特点，能够应用于任何力学领域有限元问题的求解，尤其适用于类似于板料成形优化、汽车碰撞等的非线性计算领域。

附图说明

图1为本发明实施例的实施流程示意图。

图2为用于验证本发明实施例的Spherical cap模型示意图。

图3为本发明实施例与传统CPU并行计算结果的精度对比示意图。

具体实施方式

如图1所示，本实施例基于GPU的有限元显式并行求解仿真方法的实施步骤如下：

1）将有限元计算模型数据读入主机内存中；通过主机CPU识别有限元计算模型数据中的节点和单元，找到每一个节点归属的单元并保存节点和单元之间的归属关系；将有限元计算模型数据读入GPU全局存储器；

2）设定GPU执行的线程块数量和每一个线程块所包含的线程数量；

3）通过GPU读取GPU全局存储器中存储的有限元计算模型数据，并根据线程块数量和每一个线程块所包含的线程数量，通过多线程的方法计算单元的最小时间步长、各个节点的节点加速度以及各个节点的节点作用力；

4）判断是否需要输出计算结果，如果需要则采用异步的方式将GPU全局存储器上的计算结果复制到主机的内存中，并写入到输出文件；通过CPU判断总仿真时间或者总迭代次数是否达到设定值，如果达到设定值则退出；否则返回执行步骤3）。

其中，步骤3）的详细步骤包括：

3.1）通过GPU读取GPU全局存储器中存储的所有单元数据，然后采用线程与单元一一对应的执行方式计算出每个单元对应的时间步长，并通过并行缩减的方法在时间步长中查找最小值作为最小时间步长；

3.2）通过GPU读取GPU全局存储器中存储的所有节点数据，然后采用线程与节点一一对应的执行方式计算每个节点所对应的节点加速度；

3.3）通过GPU根据约束边界条件以及最小时间步长、节点加速度采用线程与节点对应的方法计算出每个节点新的速度和位移；

3.4）通过GPU读取GPU全局存储器中存储的所有节点数据和单元数据，然后GPU根据有限元计算模型的力边界条件计算得到各个节点的节点外力，通过并行运算的方法分别计算各个单元的单元内力，采用线程与节点一一对应的执行方式将单元内力离散到整体节点内力中，采用线程与节点一一对应的执行方式将节点外力和整体节点内力做差值运算得到各个节点的节点作用力。

本实施例中使用的有限元计算模型数据包括1600个四边形单元和1681个节点。本实施例中针对该有限元计算模型数据的具体实施步骤如下：

第一步：将有限元计算模型数据（包括模型网格信息、单元信息、节点信息、边界信息、材料参数等）存储到输入文件input.inp中。在程序运行界面输入input.inp，input.inp中的数据将被程序读取到主机内存中，本实施例中以A_c表示主机内存中的有限元计算模型数据。一个标准的input.inp文件格式如下：

第二步：依据A_c的大小在GPU的全局存储器上分配与内存中数量和容量大小一致的存储空间A_d，然后将内存上的数据A_c复制到全局存储器的存储空间A_d中。本实施例的存储空间A_d中，用于存储单元信息的存储空间大小为1600*sizeof(int)*4个字节，用于存储节点信息的存储空间大小为1600*sizeof(double)*3个字节。

第三步：为了方便应力组装过程在GPU上顺利进行，在主机CPU上找到每个节点归属于哪些单元，记录归属的单元总数到数组P₁中，记录归属的单元编号到数组P₂中。在GPU的全局存储器上分配与P₁和P₂空间一样的存储空间，记为D₁、D₂。最后，复制内存中P₁和P₂到显存中的D₁、D₂中。

第四步：根据有限元模型的单元数和节点数，针对基于单元的线程和单元一一对应和基于节点的线程和节点一一对应的两种计算模式设定不同的GPU执行参数。本实施例中，基于单元的GPU执行参数为<<<7,256>>>,基于节点的GPU执行参数为<<<7,256>>>，其中7代表线程块的数目，256代表每一个线程块中所含线程的数目，两者的乘积则为线程总数。

第五步：在GPU上进行最小时间步长Δt的计算。首先由每一个线程计算到每个单元所对应的时间步长，存储到GPU全局存储器的临时数组T中，再在GPU上采用并行缩减的方法在数组T中找到最小值，即为最小时间步长Δt，并将最小时间步长的值复制到CPU中。

第六步：在GPU上采用线程与节点对应的方法，计算出每个节点所对应的加速度。

第七步：由节点加速度、时间步长、约束边界条件，在GPU上采用线程与节点对应的方法计算出节点的新的速度和位移。

第八步：在GPU上根据计算模型的力边界条件计算得到节点的外力，记为f^ext。

第九步：在GPU上并行计算各个单元的内力

在每一个线程中执行如下相同的指令：（1）从全局存储器中集合单元节点位移和速度；（2）在单元积分点中由上本构方程计算得到

第十步：根据第三步中的D1、D2，由每一个线程对应一个节点的方式并行将单元内力

离散到整体节点内力f^int中。

第十一步：在GPU上采用线程与节点对应的方法，由公式f＝f^ext-f^int计算节点作用力f。

第十二步：由当前时间步长累加，更新计算时间节点。

第十三步：判断是否需要输出计算结果，如果需要则采用异步的方式将GPU全局存储器上的计算结果复制到主机的内存中，并写入到输出文件。

第十四步：由CPU判断收敛条件，CPU判断总仿真时间或者总迭代次数是否达到设定值，如果达到设定值（符合收敛条件）则退出；否则（不符合收敛条件）则返回执行第五步。

本实施例能够充分利用GPU在并行计算方面的优势，采用GPU进行中心差分法格式显式有限元并行化求解，实现了在GPU上进行平面非线性动态问题的快速求解，为相关领域的数值模拟提供了一个高效、简单的新途径，而且通过GPU全局存储器建立合理的存储方式，能够减少计算过程中GPU与主机间的数据交换，程序优化以达到最佳计算加速度，具有运算速度快、计算成本低、硬件成本低，使用简单方便的优点。通过在显式有限元计算方法中的迭代部分引入了GPU计算，相比传统CPU计算方法，采用计算处理器不同，可以极大的提高计算效率。另一方面，与传统并行计算方法相比，本方法的硬件成本较低，易于在工程中推广使用。本实施例基于GPU通用计算平台的显式有限元并行计算，具有很强的易用性和经济性，同时本实施例中通过在计算流程中多项适应于GPU的单元和节点并行处理方法中，将线程分配和节点、单元进行有机结合起来，能够与其它细粒度并行计算方法一样避免并行计算过程中有可能出现的并行竞写等错误，能够充分利用GPU在并行计算方面的优势，具有高效、计算成本低的特点，能够应用于任何力学领域有限元问题的求解，尤其适用于类似于板料成形优化、汽车碰撞等的非线性计算领域。

为了证明本实施例的可行性，本实施例采用名为Spherical cap的计算模型进行数值验证。并且，对不同的计算模型采用不同的单元规模分别在CPU和GPU计算平台上进行计算，在对比两者计算结果精度的同时，将CPU执行时间和GPU执行时间相除得到计算加速比。本实施例中采用的CPU为Intel Core Quad Q8200，GPU为NVIDIA公司的GeForce GTX 460。如图2所示，该计算模型的半径r为22.27mm、厚度h为0.41mm，图中x、y、z为三个坐标轴，夹角α为26.67°。根据计算模型的网格对称性，采用四分之一的计算模型进行计算，并在计算模型表面上加均布的加载p，大小为600Pa。该计算模型的材料采用弹塑性本构，其杨氏模量E=10.5E6Mpa、泊松比v=0.3、弹性模量Ep=0.22E6Mpa、密度为2.45E-6kg/mm³、屈服应力σy=24000。本实施例中为了验证本实施例的计算时间，分别将有限元计算模型数据输入的单元数量规模分为4种，分别针对本实施例的方法和现有技术采用CPU计算的方法，结合上述计算模型进行实验，模拟加载时间为0.001/s，最终得到不同单元数量所对应的计算时间和计算加速比如下表1所示。

表1

单元数目	CPU计算时间/s	GPU计算时间/s	加速比
				7498	197.286	10.097	19.5
29991	1364.282	56.342	24.2
				67500	6041.722	149.829	40.3
119983	19426.470	450.003	43.2

其中，CPU计算时间代表现有技术采用CPU计算方法的计算时间，GPU计算时间代表采用本实施例基于GPU的有限元显式并行求解仿真方法的计算时间。由上表可知，对于同一个计算模型，本实施例的计算效率要远高于传统CPU串行计算方法。同时，通过对比不同计算规模下的加速比，可知本发明的计算时间加速比在材料及边界信息一致的前提下随着单元数目的增加而增加。

为了比较GPU和CPU计算结果的精度，对比计算模型的中点的位移随时间的变化趋势如图3所示，图中“GPU”关联的曲线代表本实施例并行计算结果的精度，“CPU”关联的曲线代表现有技术通过CPU并行计算结果的精度，由图3可见两者计算结果曲线完全重合，因此可知本实施例并行计算结果的精度和现有技术通过CPU并行计算结果的精度完全一致，从而证实了本实施例方法的可行性和稳定性。由以上验证算例可知本实施例可以有效解决有限元计算效率的问题，在保证计算精度的同时，明显的减少了计算时间。

以上所述仅是本发明的优选实施方式，本发明的保护范围并不仅局限于上述实施例，凡属于本发明思路下的技术方案均属于本发明的保护范围。应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理前提下的若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (2)

1.一种基于GPU的有限元显式并行求解仿真方法，其特征在于其实施步骤如下：

1）将有限元计算模型数据读入主机内存中，通过主机CPU识别所述有限元计算模型数据中的节点和单元，找到每一个节点归属的单元并保存节点和单元之间的归属关系，将所述有限元计算模型数据读入GPU全局存储器；

2）设定GPU执行的线程块数量和每一个线程块所包含的线程数量；

3）通过GPU读取所述GPU全局存储器中存储的有限元计算模型数据，并根据线程块数量和每一个线程块所包含的线程数量，通过多线程的方法计算单元的最小时间步长、各个节点的节点加速度以及各个节点的节点作用力；

4）判断是否需要输出计算结果，如果需要则采用异步的方式将GPU全局存储器上的计算结果复制到主机的内存中，并写入到输出文件；通过CPU判断总仿真时间或者总迭代次数是否达到设定值，如果达到设定值则退出；否则返回执行步骤3）。

2.根据权利要求1所述的基于GPU的有限元显式并行求解仿真方法，其特征在于，所述步骤3）的详细步骤包括：

3.1）通过GPU读取GPU全局存储器中存储的所有单元数据，然后采用线程与单元一一对应的执行方式计算出每个单元对应的时间步长，并通过并行缩减的方法在所述时间步长中查找最小值作为最小时间步长；

3.2）通过GPU读取GPU全局存储器中存储的所有节点数据，然后采用线程与节点一一对应的执行方式计算每个所述节点所对应的节点加速度；

3.3）通过GPU根据约束边界条件以及所述最小时间步长、节点加速度采用线程与节点对应的方法计算出每个节点新的速度和位移；

3.4）通过GPU读取GPU全局存储器中存储的所有节点数据和单元数据，然后GPU根据有限元计算模型的力边界条件计算得到各个节点的节点外力，通过并行运算的方法分别计算各个单元的单元内力，采用线程与节点一一对应的执行方式将单元内力离散到整体节点内力中，采用线程与节点一一对应的执行方式将节点外力和整体节点内力做差值运算得到各个节点的节点作用力。

Images