CN105373517A - 基于Spark的分布式稠密矩阵求逆并行化运算方法 - Google Patents

Abstract

本发明公开了一种基于Spark的分布式稠密矩阵求逆并行化运算方法，包括以下步骤：首先在迭代过程中对输入矩阵实施并行化的LU分解运算；接下来，以LU分解得到的分布式上三角矩阵以及分布式下三角矩阵为基础，使用递归算法求取分布式上（下）三角矩阵的逆矩阵；最后，以上面两步得到的置换矩阵和三角矩阵的逆矩阵为基础，实施分布式矩阵乘法，得到原任意输入矩阵的逆矩阵。本发明能够处理维数很大的稠密矩阵，同时拥有较高的运算效率、较好的容错性与可扩展性。

Description

基于Spark的分布式稠密矩阵求逆并行化运算方法

技术领域

本发明涉及线性代数（LinearAlgebra）运算技术领域，尤其涉及一种基于一站式大数据处理平台Spark的分布式稠密矩阵求逆运算的分布式并行计算方法。

背景技术

随着大数据时代的来临，数据量得到爆炸性的增长，人们需要计算与分析的数据规模越来越大，对运行效率与精度的要求也越来越高。而在科学计算、数据挖掘、机器学习等诸多领域，很多问题的求解都可以抽象为以矩阵运算为核心的一系列操作。然而作为常见的复杂数据计算与分析任务的中间计算步骤，大规模矩阵的求逆（inversion）运算是一个非常耗时的过程，一个主要原因是逆矩阵（inverse）中的一个元素的计算往往与输入矩阵中的多个元素存在相互依赖关系。高效的大规模矩阵求逆运算方法可以提高整个计算与分析业务的效率，节省经济成本与时间成本。

矩阵运算中，求逆运算是指求解原方块矩阵（squarematrix）的逆矩阵的运算过程。一个方阵的逆矩阵指的是和原矩阵做矩阵乘法得到的结果为单位矩阵(identitymatrix）的矩阵。常见的矩阵求逆方法包括高斯消元法、QR分解法、LU分解法等等。这三种方法中，高斯消元法直接对输入矩阵进行迭代更新；后两种方法首先对输入矩阵进行矩阵分解(factorization)，接着利用分解得到的中间结果继续求解最终结果，这样可以一定程度上减小计算的时间复杂度，提高运算效率。相比其他两种方法，LU分解法的主要优点在于：通过适当的变形，它可以一定程度上消除计算过程中矩阵各个元素之间的相互依赖关系，通过“分而治之”的方法将整个计算任务划分为多个子计算任务，再分布到各个计算机上并行计算，从而极大地提高计算效率；另一方面受益于能够多个节点并行计算，这种方法可以处理规模巨大的矩阵。

传统的单个计算机矩阵求逆串行方法由于受到计算机有限的内存、CPU等计算资源的制约，存在适用矩阵规模偏小，扩展性差等问题，在处理规模巨大的矩阵时几乎无能为力，无法满足人们对数据处理日益增长的需求。目前处理大规模矩阵求逆的主流方法有两类：一类是提高单节点计算机的硬件性能，在“纵向”提高计算机处理数据的能力，例如基于多核/众核、GPU计算等特定架构的算法实现；另一类通过分布式计算平台实现计算任务的并行化，通过增加计算机的个数，在“横向”上获得极大的数据容量以及强大的计算能力。前者对硬件的定制性较强；后者可以使用一定数量的普通计算机搭建集群，更具操作性。

然而，现有的分布式矩阵求逆实现也存在不足，主要体现在计算过程缺少统一的容错机制，当集群上一个计算机发生错误后，整个计算过程难以高效地自动恢复继续执行；运算方法的可扩展性比较差，即集群的计算能力与计算机的个数无法呈现线性或近似线性的关系；对用户不友好，无法作为一个中间步骤与其他运算步骤无缝整合在一起。

发明内容

本发明的目的是：提出一种基于Spark的分布式稠密矩阵求逆操作并行化运算方法，能够处理维数很大的稠密矩阵，同时拥有较高的运算效率、较好的容错性与可扩展性，此外该运算能够与其他运算步骤无缝对接，对用户充分友好。本发明针对的两个主要问题是：现有的矩阵规模很大，传统单个计算机串行运算方法不可行；现有的分布式矩阵并行化方案容错性、可扩展性较差，并且对用户不友好。

本发明的技术解决方案是：基于Spark分布式处理平台，利用LU分解法实现分布式矩阵求逆操作的并行化。采用此方法，用户首先把待处理的矩阵放置于分布式文件系统中，以供集群中所有的计算机共同访问。基于Spark分布式处理平台的分布式矩阵求逆的并行化方案包括下面步骤：

（1）对输入矩阵进行分布式稠密矩阵的LU分解运算，得到结果<L,U,P>，其中L为下三角矩阵，U为上三角矩阵，P为分解过程中产生的置换矩阵；

（2）利用分布式三角矩阵的求逆算法，分别计算得到分布式下三角矩阵L以及分布式上三角矩阵U的逆矩阵；

（3）利用（1）中得到的置换矩阵以及（2）中得到的两个三角矩阵的逆矩阵，实施矩阵乘法操作，得到任意输入分布式矩阵的逆矩阵。

上述第一步“分布式矩阵LU分解”阶段的工作流程是：首先对于存储在分布式文件系统上的行矩阵，按照指定的分块大小，将分布式行矩阵转化为分布式分块矩阵，各个分块均匀地分布在集群中各个计算机上。然后根据矩阵主对角线上分块的个数，确定并行化LU分解迭代计算过程所需要的迭代次数。接着在每一轮迭代过程中，在主对角线上寻找本次迭代过程中的“主分块”，对该分块进行单节点的串行LU分解算法，并把分解得到的结果广播到集群所有其他计算机上。每个计算机再利用广播得到的结果更新节点本地拥有的部分分块。每次迭代完成后，当前工作集也就是仍需要继续计算更新的分块数目越来越少，这部分分块组成的子矩阵作为下一轮迭代的输入矩阵。最后当所有的迭代完全结束，就完成了对原输入矩阵的LU分解，得到分布式下三角矩阵L，分布式上三角矩阵U，以及一个置换矩阵P。

第二步“计算上（下）三角矩阵的逆矩阵”阶段的工作流程是：对于从第一步得到的两个分布式上（下）三角稠密矩阵，首先判断矩阵的维数是否足够的小（一般低于1000），如果规模足够小可以直接在计算机本地进行串行求逆运算求出其逆矩阵；否则，对三角矩阵继续阵进行划分，得到两个子三角矩阵和一个“矩形”的子矩阵。接着使用“分而治之”的思想，对这两个子三角矩阵分别递归求解它们的逆矩阵。然后使用得到的逆矩阵，更新剩下的矩形的子矩阵。最后合并结果得到原输入上（下）三角矩阵的逆矩阵。

第三步“分布式矩阵乘法”阶段的工作流程是：首先对每一个分块根据它的标识号，算出这个分块最终会在乘积矩阵的哪个分块中被使用到；然后把集群上所有的分块全部打乱（shuffle），通过定制的哈希函数控制每个分块在集群节点上的分布，把结果矩阵的一个分块所依赖的所有输入分块全部发布到同一个计算机上；接着每个节点并行地计算乘积矩阵中的一个或多个分块，得出最终的乘积矩阵。

本发明的有益效果是：（1）针对分布式分块矩阵，通过对LU公式进行变形，消除分块矩阵各个分块之间计算过程中的相互依赖关系，在每一轮的迭代过程中计算可以并行地完成，提高运算效率。（2）在迭代计算过程中，通过不断缩小当前工作集的规模来不断减小待处理的数据量，减少不必要的中间结果的数量，以减少计算机内存的使用空间。（3）针对分布式三角矩阵，使用递归的思想将原来规模过于巨大的规模降解为规模较小的子矩阵，在多个计算机的CPU上并行地处理，从而对超大规模的三角矩阵求逆。（4）基于Spark分布式处理平台，整个计算过程拥有良好的容错性。即使集群中一个计算机发生错误，仍可以快速的恢复错误继续进行计算任务。

附图说明

图1是本发明的整体处理流程示意图；

图2是本发明的分布式矩阵LU分解流程示意图；

图3是本发明的分布式三角矩阵求逆流程示意图；

图4（A）是本发明的数据可扩展性示意图，图4（B）是本发明的节点可扩展性示意图。

具体的实施方式

如图1所示，本发明的实施是通过将分布式文件系统中的分布式行矩阵转换为分布式分块矩阵，运行分布式分块矩阵的LU分解并行化算法，接着对分解结果中的分布式上三角矩阵和分布式下三角矩阵分别运行求逆的并行化算法求出它们的逆矩阵，然后利用原来LU分解得到的置换矩阵以及求得的两个三角矩阵的逆矩阵，运行分布式矩阵乘法，得到最终的结果，即原输入矩阵的逆矩阵。

本发明的完整流程包括分布式稠密矩阵的LU分解、分布式三角矩阵（即为权利要求所述的上三角矩阵和下三角矩阵）的求逆和分布式矩阵的乘法3个部分。具体的实施方式分别说明如下：

分布式稠密矩阵的LU分解的具体实施方式为：首先从分布式文件系统中按行读取出矩阵，每一行的格式为：<rowNumber:element1,element2…>，其中rowNumber指的是当前这一行在矩阵中的行数，element指这一行的所有元素。从分布式行矩阵生成分布式分块矩阵的关键是需要判断出每一个矩阵元素属于哪一个分块，其中分块用<blockRow,blockColumn>标识。每行按照设定的分块数将自己分成若干段，并为每一段附加一个该段所属的分块标志<blockRow,blockColumn>，再将分块标志相同所有的分段聚集起来，并将它们组织为键值为分块标志的分块。这样就成功将原行矩阵变为特定结构的分块矩阵。

然后迭代计算矩阵的LU分解。需要迭代的次数由分布式分块矩阵正对角线上的分块个数确定。在进行第i轮迭代时，首先找到正对角线上的第i个分块，我们把它称为“主分块”。因为这个分块的规模足够小，所以直接把它收集到驱动程序所在的计算机，在该机器的CPU上实施行串行的矩阵LU分解运算，得到计算结果<L,U,P>。

紧接着把这三个计算结果广播（broadcast）到集群上所有工作节点上，利用L更新blockRow与主分块相同、但是blockColumn比主分块更大的所有其他分块；利用U更新blockColumn与住分块相同、但是blockRow比主分块大的所有其他分块；利用L、U与P更新blockRow比主分块大并且blockColumn也比主分块的所有其他分块。由于经过广播之后，每个计算机都拥有<L,U,P>的备份，因此整个更新计算过程可以在计算机的本地进行，而不需要通过网络传输获取数据。这一轮迭代完成后，主分块所在的blockRow这一行上所有的分块以及主分块所在的blockColumn这一列所有的分块已经是最终结果，在下面的迭代过程中，它们不属于需要关注的分块集合，将它们放入“已经完成”的分块集合中去。剩下迭代过程中的存储在内存中的“工作集”规模不断的缩小，对集群中计算机CPU资源的消耗也越来越小。也就是说，随着迭代次数的增加，每一次迭代消耗的时间越来越少。

然后选取主对角线上的下一个分块为主分块，继续进行下一轮的迭代。当所有迭代全部完成后，“已经完成”的分块集合就是原分布式矩阵LU分解的结果。同时迭代过程中还产生一个置换矩阵，它反映了分解过程中矩阵哪些行发生了互换（swap）操作。置换矩阵用一个数组表示，数组中的每个元素表示分布式矩阵中哪一行与矩阵中处于该元素下标行发生了行互换操作。这样的表示方法减小了数据的大小，同时也减小了数据在网络上的时间开销。

分布式上（下）三角矩阵递归求逆的具体实施方式为：对于一个上三角矩阵，将它的矩阵规模与一个预先设定的阈值比较。如果它的阶数比较小，可以在单个计算机上高效地完成求逆运算，就直接在单个计算机上进行串行的求逆算法；否则如果矩阵的阶数过于巨大，则可以把矩阵重新划分为三个子矩阵，其中两个子矩阵仍然是上三角矩阵，另一个剩余的分块是一个呈“矩形”的矩阵，分别把两个分布式上三角子矩阵当作新的输入矩阵，递归求得他们的逆矩阵。得到子矩阵的逆矩阵后，再运行分布式矩阵的乘法，更新矩形的子矩阵。由于递归算法的思路就是将复杂的任务分解为多个简单的子任务分别对子任务求解，最后再合并所有子任务的解得到原问题的解，因此需要合并（union）子三角矩阵的逆矩阵，得到原分布式三角矩阵的逆矩阵。分布式上三角矩阵和分布式下三角的求逆方法是非常类似的。

分布式矩阵乘法的具体实施方式是：考察参与乘法的分布式矩阵A(m*k个分块)和分布式矩阵B(k*n个分块)，它们的最终结果矩阵C（m*n个分块）中的每一个分块依赖于A中的一行(k个分块)以及B中的一列(k个分块)。单独来看，A中的每一个分块在C中使用了n次，B中的每一个分块在C中使用了m次。首先将A中每一个分块复制n份，B中每一个分块复制m份。然后以<blockRow,blockColumn>为键值，调用连接操作（join)，这样结果矩阵C中的分块就获得了它所需要的所有的分块，即A中的第i行分块和B中的第j列分块。然后每个计算机将一个结果分块所依赖的所有输入分块做相应的分块之间的乘法，再把这些乘积进行累加，得到最终分块。由于整个过程是并行的，这样就得到了最终的乘积C。实际计算中，三角矩阵与三角矩阵的乘法运算，以及三角矩阵与置换矩阵之间的乘法，都按照这种方法进行计算。

本发明提出的方法的评测结果如图4所示。测试时采用的矩阵每个元素是随机生成的双精度浮点数，矩阵元素的个数达到亿级。从测试结果可以看出，本发明提出的方法具有良好的数据可扩展性，运行时间与矩阵规模呈近似线性关系；同时计算机的个数与运行时间也呈近线性关系，表明本发明的方法具有良好的随计算机硬件规模增长而性能提升的可扩展性。

Claims (5)

1.基于Spark的分布式稠密矩阵求逆并行化运算方法，其特征是通过对分布式矩阵LU分解的并行化以及分布式三角矩阵求逆的并行化对任意输入稠密矩阵进行求逆，包括以下步骤：

1)把待处理的大规模稠密矩阵存储到分布式文件系统中，得到分布式稠密矩阵，以供集群中所有的计算机后续共同访问；

2)对分布式稠密矩阵实施并行化的LU分解运算，得到上三角矩阵、下三角矩阵和置换矩阵；

3)针对分解得到的上三角矩阵和下三角矩阵，分别实施分布式矩阵的并行化求逆运算，得到上三角矩阵的逆矩阵和下三角矩阵的逆矩阵；

4)最后利用分解得到的置换矩阵、上三角矩阵的逆矩阵以及下三角矩阵的逆矩阵，实施分布式矩阵的并行化乘法运算，实现最终的逆矩阵的求解。

2.根据权利要求1所述基于Spark的分布式稠密矩阵求逆并行化运算方法，其特征在于，所述步骤2）包括：首先按照指定的分块大小，将所述分布式稠密矩阵转化为相应的分布式分块矩阵，各个分块均匀地分布在集群中各个计算机上；然后迭代计算矩阵的分解结果，根据矩阵主对角线上分块的个数，确定总共的迭代轮数；接着在每一轮迭代过程中，在主对角线上寻找属于本次迭代过程中的主分块，在单个节点上对所述主分块实施串行LU分解运算，并把分解得到的结果广播到集群所有其他计算机上；接下来每个节点利用广播得到的结果更新节点本地拥有的矩阵分块；每次迭代完成后，当前工作集的分块数目越来越少，这部分分块组成的子矩阵作为下一轮迭代的输入矩阵；最后当所有的迭代完全结束，就完成了对原输入矩阵的LU分解。

3.根据权利要求1所述基于Spark的分布式稠密矩阵求逆并行化运算方法，其特征在于，所述步骤3）包括：对于从所述步骤2）得到的上三角矩阵以及下三角矩阵，首先判断矩阵的维数是否足够的小，如果维数足够小，直接在计算机本地进行串行求逆运算求出其逆矩阵；否则，对上三角矩阵和下三角矩阵继续进行划分，得到两个子三角矩阵和一个矩形的子矩阵；接着对所述两个子三角矩阵分别递归求解它们的逆矩阵；然后使用得到的逆矩阵，更新剩下的矩形的子矩阵；最后合并结果得到原输入上三角矩阵的逆矩阵以及下三角矩阵的逆矩阵。

4.根据权利要求2所述基于Spark的分布式稠密矩阵求逆并行化运算方法，其特征在于，所述步骤4）包括：首先对每一个分块根据其标识号，算出该分块最终会在乘积矩阵的哪个分块中被使用到；然后把集群上所有的分块全部打乱，通过哈希函数控制每个分块在集群节点上的分布，把结果矩阵的一个分块依赖的所有输入分块全部分布到同一个计算机上；接着每个节点并行地计算乘积矩阵中的一个或多个分块，得出最终的乘积矩阵。

5.根据权利要求3所述基于Spark的分布式稠密矩阵求逆并行化运算方法，其特征在于，若所述矩阵的维数低于1000，则认为维数足够的小。