传热设备热传导的三对角异构众核并行求解方法及系统

技术领域

本发明涉及计算机处理技术领域，尤其涉及一种传热设备热传导的三对角异构众核并行求解方法及系统。

背景技术

本部分的陈述仅仅是提供了与本发明相关的背景技术信息，不必然构成在先技术。

现有工业上有许多以热传导为主的传热过程，如橡胶制品的加热硫化、钢锻件的热处理等。在窑炉、传热设备和热绝缘的设计计算及催化剂颗粒的温度分布分析中，热传导规律都占有重要地位。在高温高压设备(如氨合成塔及大型乙烯装置中的废热锅炉等)的设计中，也需用热传导规律来计算设备各传热间壁内的温度分布，以便进行热应力分析。以传热设备为例，传热设备的热传导过程是一个复杂的物理现象，即其传热过程难以直观显示且不易被分析。现有技术中常采用经验法或者仿真模拟的方法进行热传导过程进行分析或预测。经验法依靠操作人员的以往经验，得到的分析结果偏差较大，往往与实际情况不符。所以一般会选择构建仿真模型建立热传导方程并求解的方法从而获得更为准确的热传导情况分析结果。

发明人发现，在热传导方程求解的过程中会产生三对角矩阵方程，而现有技术对于三对角矩阵方程的并行算法已有多种，包括循环递增的并行算法、循环约化的并行算法、分裂法等，但并行规模非常有限。在传热设备热传导过程中随着模拟仿真问题逐渐复杂，数值模型涉及的参数愈多，形成的三对角矩阵方程规模越来越大，矩阵规模高达百万维甚至千万维以上，工作站和小规模单集群处理器性能无法满足大规模计算对于计算速度、精确度的要求，导致热传导过程的分析过程耗时较长，结果不够准确。

发明内容

针对现有技术存在的不足，本发明的目的是提供一种传热设备热传导的三对角异构众核并行求解方法及系统，能够实现在新一代神威架构上对大规模三对角矩阵方程进行快速求解，解决传热设备在热传导模拟仿真过程中，需要求解的三对角矩阵方程规模较大，导致热传导过程的分析过程耗时较长，结果不够准确的问题。

为了实现上述目的，本发明是通过如下的技术方案来实现：

本发明第一方面提供了一种传热设备热传导的三对角异构众核并行求解方法，包括以下步骤：

采集传热设备的热传导参数，建立传热模型；输入热传导参数至传热模型，构建热传导方程；

对热传导方程通过差分和转化得到关于热传导的三对角线性方程组；

根据三对角线性方程组建立进程，进程包括由主核组成的主核组和由从核组成的从核组；通过并行消元的方法对三对角矩阵方程数据消除数据依赖；对并行消元后形成的新的三对角矩阵方程抽取部分组成小型三对角矩阵方程；采用追赶法求解缩减三对角矩阵方程；

将各个进程以任务并行的方式回代输出三对角矩阵方程剩余全部解；并根据求解结果绘制温度变化曲线，获得热传导过程的温度变化。

进一步的，热传导参数包括介质的热传导率、介质的比热、介质的密度和表面热传导系数。

进一步的，通过并行消元的方法对三对角矩阵方程数据消除数据依赖之前，使用条件运算符判断进程和矩阵规模之间是否整除，若整除，则通过并行消元的方法对三对角矩阵方程数据消除数据依赖，若非整除，则在数据传输过程中使用第一个进程和最后一个进程发送和接收非整除数据。

进一步的，通过并行消元的方法对三对角矩阵方程数据消除数据依赖的具体步骤为：

通过分段插桩计时，找到时间占比最长即计算热点部分为三对角矩阵方程并行消元部分，此时将此部分放进从核优化，进一步缩短求解时间。

通过修改求解表达式消除块间数据依赖，将三对角矩阵分块进行消元，使用直接访问内存(DMA)操作将主存中的数据放到从核局部存储空间进行访存优化。

更进一步的，将三对角矩阵分块进行消元的具体步骤为：将三对角矩阵进行分块得到多个子块，从核组中的每个从核负责处理一个子块，对子块中的三对角矩阵利用多个进程并行消去上对角线，再消去下对角线，分别引入非零元。

更进一步的，通过开辟两个数组大小的数据空间，将计算和通信叠加，在计算的同时获取下一次的数据。

进一步的，抽取的部分矩阵为：第一个进程的从核组中，在从核分块后的第一个子块的最后一行数据，其余所有子块的第一行和最后一行数据；最后一个进程的从核组中，在从核分块后的最后一个子块的第一行数据，其余子块的第一行和最后一行数据；剩余进程的从核组中，在从核分块后的每个子块的第一行和最后一行数据。

进一步的，将各个进程以任务并行的方式回代输出三对角矩阵方程剩余全部解的过程具体包括三个部分：第一个进程的求解，最后一个进程的求解以及其余所有进程的求解；第一个进程和最后一个进程各有两个求解表达式，其余所有进程仅有一个求解表达式。

进一步的，还包括验证求解结果正确性及误差，验证求解结果正确性及误差的具体步骤为：通过设置检查函数，在得出求解结果后，用残差来衡量误差程度。

本发明第二方面提供了一种传热设备热传导的三对角异构众核并行求解系统，包括：

数据获取模块，被配置为采集传热设备的热传导参数，建立传热模型；输入热传导参数至传热模型，构建热传导方程；

方程转化模块，被配置为对热传导方程通过差分和转化得到关于热传导的三对角线性方程组；

方程组求解模块，被配置为根据三对角线性方程组建立进程，进程包括由主核组成的主核组和由从核组成的从核组；通过并行消元的方法对三对角矩阵方程数据消除数据依赖；对并行消元后形成的新的三对角矩阵方程抽取部分组成小型三对角矩阵方程；采用追赶法求解缩减三对角矩阵方程；

温度变化分析模块，被配置为将各个进程以任务并行的方式回代输出三对角矩阵方程剩余全部解；并根据求解结果绘制温度变化曲线，获得热传导过程的温度变化。

以上一个或多个技术方案存在以下有益效果：

本发明提供了一种传热设备热传导的三对角异构众核并行求解方法及系统，在温度变化分析过程中，针对现有技术无法解决传热设备热传导现象模拟仿真过程中，需要求解的三对角矩阵方程规模较大，导致热传导过程的分析过程耗时较长，结果不够准确的问题，基于新一代国产“神威”系列的超级计算机，采用申威26010pro芯片、异构众核的体系架构，在原有的算法基础上结合追赶法，增加非整除部分、增加验证部分，使算法更加完善。在并行层次划分上结合三种方式，从核上进行数据分块，采用访存优化、计算和通信叠加优化、SIMD数据向量化，提高在神威架构上求解三对角矩阵方程的计算效率。本发明改进后的算法能够缩短求解方程组阶段的时间，之后得到离散点上的数值解，最后经过后处理，就可以分析得到准确的传热设备在热传导过程中的温度变化。

本发明附加方面的优点将在下面的描述中部分给出，部分将从下面的描述中变得明显，或通过本发明的实践了解到。

附图说明

构成本发明的一部分的说明书附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。

图1为本发明实施例一中传热设备热传导的三对角异构众核并行求解方法的流程图；

图2为本发明实施例一中进程并行和线程级并行方案的框架结构图；

图3为本发明实施例一中矩阵分块并行分裂法主从核工作状态示意图；

图4为本发明实施例一中实验验证过程的实验效果图；

图5为本发明实施例二中方程组求解模块与温度变化分析模块结构示意图。

具体实施方式

应该指出，以下详细说明都是示例性的，旨在对本申请提供进一步的说明。除非另有指明，本文使用的所有技术和科学术语具有与本申请所属技术领域的普通技术人员通常理解的相同含义。

三对角矩阵是稀疏矩阵中较为典型的一类结构。偏微分方程隐式差分格式的计算往往涉及到求解三对角矩阵方程，而偏微分方程的数值求解又出现在很多科学工程计算领域中。比如在模拟热传导、油气资源勘探、数值海洋模型、流体动力学等的计算最终都归结为求解三对角矩阵方程。

新一代国产“神威”系列的超级计算机，采用申威26010pro芯片、异构众核的体系架构，由多个核组组成，每个核组含有1个主核和64个从核，每个从核的局部存储空间为256KB，从核阵列为8*8分布。在新一代申威众核处理器中已经部署了Xmath、SWBlas、SWSparse等基础求解库，但目前尚未有针对大规模三对角矩阵方程的并行分裂法求解器。由于申威众核处理器中受到从核局部存储空间的限制，使得现有并行分裂法无法充分结合神威的硬件架构，求解三对角矩阵方程面临计算时间长，矩阵规模不够大等难题的困扰。

实施例一：

本发明实施例一提供了一种传热设备热传导的三对角异构众核并行求解方法，如图1所示，包括以下步骤：

S1：采集传热设备的热传导参数，建立传热模型；输入热传导参数至传热模型，构建热传导方程。

以换热器为例，通过换热器设备信息得到介质的热传导率λ、介质的比热c、介质的密度ρ以及表面热传导系数h等热传导参数。

一维热传导方程的表达式如下：

其中，温度u＝u(x,t)，

S2：对热传导方程通过差分和转化得到关于热传导的三对角线性方程组。

其具体过程为：

对

对

因此一维热传导表达式可以转化为：

有限差分法具体位置将待求解的区域划分成无数多个微小的网格，网格上承载在位置和温度信息，用网格上的温度近似替代物体所在位置真实的温度。因此对横坐标轴k,纵坐标轴j坐标轴对应的两个方向进行网格拆分，记u

由差分近似方程可以得到

其中，网格傅里叶数

S3：根据三对角线性方程组建立进程，进程包括由主核组成的主核组和由从核组成的从核组；通过并行消元的方法对三对角矩阵方程数据消除数据依赖；对并行消元后形成的新的三对角矩阵方程抽取部分组成小型三对角矩阵方程；采用追赶法求解缩减三对角矩阵方程；将各个进程以任务并行的方式回代输出三对角矩阵方程剩余全部解；并根据求解结果绘制温度变化曲线，获得热传导过程的温度变化。

具体求解过程包括：

步骤1：获取三对角矩阵方程数据，根据三对角矩阵方程数据建立进程，通过并行消元的方法对三对角矩阵方程数据消除数据依赖。

其中，进程包括由主核组成的主核组和由从核组成的从核组，本实施例设置使用的进程数为p，从核数为64，每个进程负责的矩阵大小为m，每个从核负责的数据大小为m/64，进程并行和线程级并行方案具体如图2所示。

步骤2：对并行消元后形成的新的三对角矩阵方程抽取部分组成小型三对角矩阵方程:。

步骤3：采用追赶法求解缩减三对角矩阵方程。

步骤4：将各个进程以任务并行的方式回代输出三对角矩阵方程剩余全部解。

步骤5：验证求解结果正确性及误差。

三对角矩阵方程的公式为：Ax＝r，其中x是需要求解的未知量，系数矩阵A是三对角矩阵，r是已知的右端向量。非零元分布在主对角线及上下两条次对角线上，其他所有元素都为0。本实施例中表达式

如下：

步骤1中，通过并行消元的方法对三对角矩阵方程数据消除数据依赖之前，使用条件运算符判断进程和矩阵规模之间是否整除，若整除，则通过并行消元的方法对三对角矩阵方程数据消除数据依赖，若非整除，则在数据传输过程中使用第一个进程和最后一个进程发送和接收非整除数据。

为了实现上述过程，本实施例首先在并行消元之前进行声明变量操作，声明变量时就通过条件运算符来对各个进程负责的行数进行分工，条件运算符的格式如下所示，表示符合条件执行表达式1，否则2：

条件？表达式1:表达式2

对应代码中的运算是：mem＝(rank＝＝size-1)？(m+final):m

其中rank为进程号，size为进程数，final为矩阵规模和进程不能整除时的余数，m为矩阵规模和进程整除时每个进程负责的行数，mem为各个进程负责的行数。

若进程和矩阵规模不能够整除，则并行消元前后、并行求解前后数据的发送和接收部分，都需要传输非整除部分的数据。此时都需要MPI的发送和接收的相关函数进行处理，本实施例在此采用条件判断语句来实现。

通过并行消元的方法对三对角矩阵方程数据消除数据依赖的具体步骤为：

(1)通过分段插桩计时，找到时间占比最长即计算热点部分为三对角矩阵方程并行消元部分，此时将此部分放进从核优化，进一步缩短求解时间。本实施例中，将耗时占比百分之六十的三对角矩阵放进从核优化。

(2)通过修改求解表达式消除块间数据依赖，将三对角矩阵分块进行消元，使用DMA操作将主存中的数据放到从核局部存储空间进行访存优化。

需要特别说明的是，新一代申威众核处理器众核并行模式共分为四种，主从加速并行、主从协同并行、主从异步并行、主从动态并行。本实施例在并行消元模块的并行消元过程中采用的是主从加速并行模式，如图3所示，主核完成的是从核部分不可并行的计算和通信，从核完成计算任务时，主核处于等待中。

通过数据分块实现线程级并行，在对从核中数据进行分块时，通过修改算法，消除块间依赖。对于要在从核中进行计算的矩阵均分给64个从核，实现线程并行和进程并行相结合，提高计算效率。

具体的，将三对角矩阵分块进行消元的具体步骤为：将三对角矩阵进行分块得到多个子块，从核组中的每个从核负责处理一个子块，对子块中的三对角矩阵利用多个进程并行消去上对角线，再消去下对角线，分别引入非零元。

在新一代申威众核处理器上，LDM指的是每个从核具有的一块高速本地局部数据存储空间，使用DMA操作可以将直接访问主存变为对本地LDM的访问，因此，本实施例采取的DMA数据传输是非阻塞的方式，采用DMA操作来实现访存优化，提高并行求解方法的计算效率。

从核中同时存在计算和访存，在计算时进行下一次的访存能够充分将计算和通信叠加，省去了访存的时间开销，提高计算效率，因此本实施例采用计算和通信叠加对算法的对并行消元过程进行优化。

具体的，通过开辟两个数组大小的数据空间，将计算和通信叠加，在计算的同时获取下一次的数据。

步骤2中，抽取的部分矩阵为：第一个进程的从核组中，在从核分块后的第一个子块的最后一行数据，其余所有子块的第一行和最后一行数据；最后一个进程的从核组中，在从核分块后的最后一个子块的第一行数据，其余子块的第一行和最后一行数据；剩余进程的从核组中，在从核分块后的每个子块的第一行和最后一行数据。

步骤4中，将各个进程以任务并行的方式回代输出三对角矩阵方程剩余全部解的过程具体包括三个部分：第一个进程的求解，最后一个进程的求解以及其余所有进程的求解；第一个进程和最后一个进程各有两个求解表达式，其余所有进程仅有一个求解表达式。

具体的，第一个进程分为两个求解公式，第一个从核求解和其余63个从核求解。

最后一个进程分为两个求解公式，最后一个从核求解和其余63个从核求解。

其余所有进程含有一个求解公式，使用64次for循环实现。

形成小规模三对角矩阵方程后，在第一个进程使用追赶法进行求解缩减三对角矩阵方程。

根据缩减三对角矩阵方程和原方程的解之间的关系，使用任务并行的方式各个进程分别回代求出矩阵方程全部解。

本实施例中，采用数据并行向量化对算法进行优化，新一代申威众核处理器主核支持256位SIMD扩展指令，从核支持512位SIMD扩展指令。为了实现数据处理的向量化，将并行求解和新一代申威众核处理器芯片的向量指令结合，本实施例中为并行求解过程中加入SIMD数据向量化操作。其中，在主核中向量化包括以下依次进行的循环的展开和分裂、变量替换以及操作替换步骤。

具体的：主核一次可处理4个双精度浮点运算，向量化访存操作要求必须是连续的，数组正好适用。

向量化之前使用load映射到扩展向量类型，向量化之后使用store映射到标准类型，主核上双精度扩展变量类型是duoblev4。向量化相当于添加一个for循环，可以降低对指令访问带宽的要求。在主核中向量化的实现如下：

for(i＝0；i＜n；i＝i+4)

{

s_x[i+3：i]＝(s_r[x+3：i]-s_b[x+3：i])/s_a[x+3：i]

}

通过上述方法最终求解得出关于u(x，t)的值，从而绘制温度变化曲线，进行进一步的分析一个区域内温度是如何随时间变化的。

本实施例通过使用网格拆分得出各个点对应温度，从而可以观测到一个区域内温度的变化。在高温高压设备的设计中，也需用热传导规律来计算设备各传热间壁内的温度分布，以便进行热应力分析。因此本发明传热设备热传导的三对角异构众核并行求解方法能够快速得到温度变化和分布，为高温高压设备设计提供了新的技术支持。

步骤5中，验证求解结果正确性及误差的具体步骤为：通过设置检查函数，在得出求解结果后，用残差来衡量误差程度。

在用编程语言求解矩阵方程时，由于只能用有限的二进制位数来表示一个数，所以由

具体的，通过

为了验证本实施例方法的有效性，对其进行了实验验证，测试环境为国产异构新一代神威超级计算机，测试对象为三对角矩阵方程，测试核心为1个主核，64个从核，实验结果如图4所示。

图4中是使用两个核组，矩阵维数为832000的情况下纯主核运行时间和优化加速后时间的对比效果图，包括热点部分优化前后、整个程序运行总时间优化前后的对比，从中可以看出热点优化前后可达到19倍的加速，总时间优化前后可达到2.69倍的加速。因此本发明实现了在新一代神威架构上对大规模三对角矩阵方程进行快速求解，解决模拟仿真问题中，在新一代申威众核处理器上需要求解的三对角矩阵方程规模较大，导致计算时间过长的问题。

针对现有并行规模有限，在传热设备热传导过程中随着模拟仿真问题逐渐复杂，数值模型涉及的参数愈多，形成的三对角矩阵方程规模越来越大，矩阵规模高达百万维甚至千万维以上，工作站和小规模单集群处理器性能无法满足大规模计算对于计算速度、精确度的要求的缺点，本发明通过矩阵分块并行分裂法解决了上述问题。

本发明的矩阵分块并行分裂法增加非整除部分，使得矩阵在求解时不受限于矩阵规模和进程不能够整除的情况，在求解三对角矩阵方程时，适用于任意情况。在求解缩减三对角矩阵方程时，采用最经典的串行追赶法，访存量小，适合求解小规模矩阵方程。在并行求解的过程中，将并行求解过程进行SIMD向量化，降低了对指令访问带宽的要求，提高了计算效率。

本发明通过采用分段插桩计时，将部分三对角矩阵放进从核优化，通过修改表达式，消除并行消元模块内的数据依赖。将数据平均分块，实现线程级并行，达到负载均衡，使并行分裂法和神威特殊的架构紧密结合，提高计算效率。在并行消元过程中使用DMA操作，结合批量回传、计算和通信叠加，提高访存效率，极大的缩短计算时间。

本发明还增加了验证部分，通过设置检查函数用残差验证最终求解的正确性及误差。

实施例二：

本发明实施例二提供了一种传热设备热传导的三对角异构众核并行求解系统，包括：

数据获取模块，被配置为采集传热设备的热传导参数，建立传热模型；输入热传导参数至传热模型，构建热传导方程；

方程转化模块，被配置为对热传导方程通过差分和转化得到关于热传导的三对角线性方程组；

方程组求解模块与温度变化分析模块，如图5所示：

方程组求解模块包括：

并行消元模块，被配置为获取三对角矩阵方程数据，根据三对角矩阵方程数据建立进程，进程包括由主核组成的主核组和由从核组成的从核组；通过并行消元的方法对三对角矩阵方程数据消除数据依赖；

第一缩减模块，被配置为对并行消元后形成的新的三对角矩阵方程抽取部分组成小型三对角矩阵方程；

第二缩减模块，被配置为采用追赶法求解缩减三对角矩阵方程；

并行求解模块，被配置为将各个进程以任务并行的方式回代求出三对角矩阵方程剩余全部解；

温度变化分析模块包括：

并行求解模块，被配置为将各个进程以任务并行的方式回代输出三对角矩阵方程剩余全部解；并根据求解结果绘制温度变化曲线，获得热传导过程的温度变化；

验证模块，被配置为验证求解结果正确性及误差。

在整个求解方法中，并行方式分为任务并行、线程并行、SIMD数据并行三种。因此本实施例在实现并行算法时结合了这三种并行方式合理分配任务和计算量，充分结合神威的硬件架构，提高计算效率。具体体现在以下方面：

(1)并行求解模块使用任务并行，不同的进程执行不同的求解任务。第一个进程、最后一个进程以及其他所有进程。求解x采用不同的公式，但是进程处理矩阵的行数几乎无差别，相比与大规模的矩阵可以忽略不计。不同的进行并行的执行不同的求解公式，提高计算效率。

(2)从核中采用的是线程并行，在本发明中体现在64个从核完成并行消元模块的计算任务。此模块通过数据平均分块给64个从核以达到负载均衡。

(3)并行求解模块使用SIMD数据并行，通过分段插桩计时，发现在矩阵规模不断增大时，除了并行消元模块占比最大之外，其次就是并行求解模块。此模块通过SIMD数据并行以减少指令数，降低对指令访问带宽的要求。

在传热设备热传导的三对角异构众核并行求解系统中，本实施例首先使用分段插桩计时的方式找到算法的热点模块，即并行消元模块的耗时最长，达到整个并行求解算法的百分之六十。

其次在并行消元模块通过第一个进程分块后，块内串行，块间并行。因此将进程级的矩阵块划分到线程级的64个从核时，块间存在数据依赖。此时修改算法，将进程级的矩阵块进行均分到64个从核，通过修改并行求解模块的表达式及缩减模块的矩阵规模，矩阵块数变为原来的64倍，缩减模块的矩阵维数变为原来的64倍，并行求解模块表达式由3种求解公式变为5种求解公式。

再将并行消元模块放进从核，根据矩阵规模进行均匀划分到64个从核，矩阵规模N、每个从核负责的数组大小arr_num、从核个数number之间的关系式为：

N＝arr_num×number

随着矩阵规模成倍的扩大，达到上百万上千万，由于受到从核局部存储空间的限制，不能一次性直接将全部的数据都放进从核优化，本实施例的从核部分采用批量回传的方式，即需要进行多次的读写操作。

具体的操作就是当64个从核读取完数据并计算传回主存之后，再调用64个从核进行读取和回传，不断重复数组的读取和回传直至整个数组都读取计算回传完毕。此时for循环次数和每个从核传数据大小核矩阵规模间之间的关系为：

count(循环次数)＝N/size(进程数)÷64÷arr_num

之后，本实施例将并行消元模块的数据传进从核，并在从核中经过计算将结果回传给主核，使用CRTS_dma_iget,CRTS_dma_iput函数。使用DMA操作，将直接访问主存变为访问LDM，提高计算效率。

需要特别说明的是，本实施例在DMA操作进行访存优化的基础上使用计算和通信叠加的优化策略，在从核中开辟两个数组空间，在并行消元模块完成各个参数的计算时发起第二次DMA访存操作，使计算和通信叠加，提高整体的计算效率。

以上实施例二中涉及的各步骤与方法实施例一相对应，具体实施方式可参见实施例一的相关说明部分。

本领域技术人员应该明白，上述本发明的各模块或各步骤可以用通用的计算机装置来实现，可选地，它们可以用计算装置可执行的程序代码来实现，从而，可以将它们存储在存储装置中由计算装置来执行，或者将它们分别制作成各个集成电路模块，或者将它们中的多个模块或步骤制作成单个集成电路模块来实现。本发明不限制于任何特定的硬件和软件的结合。

上述虽然结合附图对本发明的具体实施方式进行了描述，但并非对本发明保护范围的限制，所属领域技术人员应该明白，在本发明的技术方案的基础上，本领域技术人员不需要付出创造性劳动即可做出的各种修改或变形仍在本发明的保护范围以内。