一种基于RAS对等式交互的有限元程序并行化方法

本发明属于力学工程设计与计算领域，具体涉及一种基于ras对等式交互的有限元程序并行化方法。

背景技术：

有限元方法作为一种应用于解决工程学和数学物理学的数值方法，被广泛应用于结构分析等工程领域典型问题。通过将结构的特定研究区域离散为数个特定形状的单元，能够将最初描述物理现象的偏微分方程转变为某种形式的矩阵方程，利用特定离散点的计算结果来表征结构内任一点的物理量计算结果，从而简化问题。利用有限元方法对工程问题进行分析，能够为工程实际中问题的研究提供指导，缩短结构设计周期、节约试验成本。

然而随着某些工程问题越来越复杂，有限元方法的求解规模也越来越大，传统的开源串行计算方法无论在性能还是求解速度上都难以满足要求，同时也对计算机内存提出了更高的要求。除此之外，随着eigen、mkl等线性代数运算库的飞速发展，如何将有限元过程与现有的功能强大的库函数相结合从而显著减少总运算时间，也成为了一个严峻的问题。

当考虑结构计算的材料非线性、结构非线性和接触非线性时，有限元计算规模迅速增加，如果考虑抗震特性或者可靠度等情况，将进一步增加计算耗时。随着科技水平的发展，计算机硬件技术得到了快速发展，除了大型化的超级计算机面向科研人员开放使用，小型超级计算集群也成为了很多科研单位的计算资源。然而，在硬件技术提高的同时，应该开发对应性能的有限元分析软件。

采用并行计算是提高有限元计算效率的必要前提，有限元离散的基本思想结合并行计算的分而治之，使有限元计算的并行化有了依据。当前实现有限元并行算法的方法较少，仅部分商用软件实现了该功能，且存在很多局限性：1)商用软件不开放其源代码，使得求解问题的局限性增大，对其参数的设置需要不断尝试，收敛性也不易控制；2)不利用使用者自主开发高效的本构模型和算法，针对不断涌现的工程问题不能及时的更新；3)商业软件一般价格比较昂贵，且每次重新升级都需要重新授权。与此同时，针对并行环境计算的各种开源程序不断涌现，使得自主开发并行程序成为了可能。因此亟需一种基于高效、便捷、与其他库函数结合方便的有限元并行分析方法，从而大幅缩短有限元方法中大型矩阵的形成时间，减少试验成本。

技术实现要素：

本发明的目的是提供一种基于ras对等式交互的有限元程序并行化方法，通过分析现有有限元程序求解的特点，分析子空间迭代方法中常用算法的求解流程，在利用mpi函数库实现并行化的基础上，在对现有有限元程序改动最小的情况下，实现高效率有限元并行计算。

为实现上述目的，本发明采用的技术方案是：

一种基于ras对等式交互的有限元程序并行化方法，包括以下步骤：

步骤一，初始化mpi并行环境，获取进程号和总进程个数；

步骤二，读取数据文件；

步骤三，建立局部到整体的映射关系；

步骤四，确定本进程存储刚度矩阵的起始行istart和终止行iend，预分配分布式刚度矩阵；

步骤五，按照未知量eqmap映射关系组装分布式整体刚度矩阵k；

步骤六，使用并行求解器求解有限元方程组；

步骤七，按照未知量eqmap映射获取本进程结果，并散发其他进程需求的解；

步骤八，后处理，终止并行环境。

具体地，步骤一中，所述mpi是一种跨语言的通讯协议，用于编写并行计算机；每个所述进程都对应一个进程号。

具体地，步骤二中，所述数据文件包括有限元模型数据，所述有限元模型数据包括节点、单元、分区、材料、载荷、约束和交互界面；其中，di为第i个分区，为重叠型分区，为限制分区，为非重叠分区；分区di包含限制分区和镜像点区域。

具体地，步骤三中，所述映射关系包括自由度映射dofmap和未知量映射eqmap；其中，自由度映射dofmap是局部自由度编号和整体自由度编号之间的关系；未知量映射eqmap是局部未知量编号和整体未知量编号之间的关系，由不包含约束的节点自由度按照dof映射累加之后填充得到。

具体地，步骤四中，所述预分配刚度矩阵按照未知量映射eqmap对整体刚度矩阵进行遍历，形成整体刚度矩阵的指示矩阵，确定第i个进程的整体刚度矩阵ki中的主对角阵和非主对角阵每一行非零元素个数，然后预分配分布式刚度矩阵。

具体地，步骤五中，所述分布式整体刚度矩阵k按照每个节点分别组装ki；

其中，为di刚度矩阵非重叠部分，为di刚度矩阵重叠部分，i∈1，2，…，n，n为进程个数。

具体地，步骤六中，所述有限元方程组为：

ku＝f

其中，u为有限元控制方程中的位移，f为有限元控制方程中的等效节点荷载；

对于重叠型分区，第i个分区di区域上的解为：

ui＝riu

其中，ri为di的限制矩阵，将整个区域的解集成可以得到：

对于非重叠型分区，上的解为：

其中，为的限制矩阵；将整个区域的解集成可以得到：

其中，整个区域的解由所有的限制分区上的解组成，对于di的重叠区域的解和则由相邻的分区和传递得到，在上不求解；为限制区域的解，为相邻区域传递回来的解，i∈1，2，…，n，n为进程个数。

与现有技术相比，本发明的有益效果是：本发明的有限元程序并行化方法通过分析现有有限元程序求解的特点，分析子空间迭代方法中常用算法的求解流程，在利用mpi函数库实现并行化的基础上，在对现有有限元程序改动最小的情况下，显著提高了有限元并行计算的效率。

附图说明

图1为本发明实施例一种基于ras对等式交互的有限元程序并行化方法的流程图；

图2为本发明实施例中验证算例有限元网格和分区示意图；

图3为本发明实施例中ras区域分解示意图；

图4为本发明实施例中验证算例z向位移图；

图5为本发明实施例中验证算例求解消耗总时间与进程数量关系图；

图6为本发明实施例中验证算例求解加速比与进程数量关系图。

具体实施方式

下面将结合本发明中的附图，对本发明的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动条件下所获得的所有其它实施例，都属于本发明保护的范围。

如图1所示，本实施例提供了一种基于ras对等式交互的有限元程序并行化方法，主要包括以下几个步骤：

第一步，初始化mpi并行环境，获取进程号和总进程个数；

第二步，读取数据文件，尺寸为200m×200m×100m的立方体，网格划分为100×100×10，底部固定约束，受自重荷载作用，有限元网格和分区如附图2所示；材料参数如表1所示；

表1材料参数

如附图3所示，di为第i个分区，为重叠型分区，为限制分区，为非重叠分区；分区di包含限制分区和镜像点区域；

第三步，建立局部到整体映射关系，包括自由度映射dofmap和未知量映射eqmap；其中，自由度映射dofmap是局部自由度编号和整体自由度编号之间的关系；未知量映射eqmap是局部未知量编号和整体未知量编号之间的关系，由不包含约束的节点自由度按照dof映射累加之后填充得到；

第四步，确定本进程存储刚度矩阵的起始行istart和终止行iend，按照未知量映射eqmap对整体刚度矩阵进行遍历，形成整体刚度矩阵的指示矩阵，确定第i个进程的整体刚度矩阵ki中的主对角阵和非主对角阵每一行非零元素个数，然后预分配分布式刚度矩阵；

第五步，按照未知量eqmap映射关系组装分布式整体刚度矩阵k，所述分布式整体刚度矩阵k按照每个节点分别组装ki；

其中，为di刚度矩阵非重叠部分，为di刚度矩阵重叠部分，i∈1，2，…，n，n为进程个数。

第六步，使用并行求解器求解有限元方程组；有限元方程组为：

ku＝f

其中，u为有限元控制方程中的位移，f为有限元控制方程中的等效节点荷载；

对于重叠型分区，第i个分区di区域上的解为：

ui＝riu

其中，为的限制矩阵，将整个区域的解集成可以得到：

对于非重叠型分区，上的解为：

其中，为的限制矩阵；将整个区域的解集成可以得到：

其中，整个区域的解由所有的限制分区上的解组成，对于di的重叠区域的解和则由相邻的分区和传递得到，在上不求解；为限制区域的解，为相邻区域传递回来的解，i∈1，2，…，n，n为进程个数。

第七步，按照未知量eqmap映射获取本进程结果，并散发其他进程需要的解；

第八步，后处理，终止并行环境。

采用本实施例的基于ras对等式交互的有限元程序并行化方法对有限元进行计算，得到的算例z向位移如图4所示，验证算例求解消耗总时间与进程数量关系如图5所示，验证算例求解加速比与进程数量关系如图6所示；从上图4至6可以看出，采用本实施例的基于ras对等式交互的有限元程序并行化方法，相较于现有的有限元计算方法，显著提高了有限元并行计算的效率。

尽管已经示出和描述了本发明的实施例，对于本领域的普通技术人员而言，可以理解在不脱离本发明的原理和精神的情况下可以对这些实施例进行多种变化、修改、替换和变型，本发明的范围由所附权利要求及其等同物限定。