GPU上稀疏矩阵的LU分解:符号执行和动态并行算法优化及实现

288→××GPU上大型矩阵的端到端LU分解杨霞俄亥俄州立大学计算机科学与工程系关闭WI，USAosu.edu加甘·阿格拉瓦尔奥古斯塔大学计算机与网络科学学院关闭IN，USAgagrawal@augusta.edu摘要稀疏矩阵的LU分解是电路仿真等许多工程和科学问题的重要计算步骤。 已经有许多努力来并行化和扩展该算法，其中包括最近针对GPU的努力。然而，由于高存储器要求和数据依赖性，在GPU上部署完整的稀疏LU分解工作流仍然是困难的。在本文中，我们提出了稀疏LU分解的第一个完整的GPU解决方案 为了实现这一目标，我们提出了一个符号执行阶段的核外实现，从而消除了由于大型中间数据结构的瓶颈。接下来，我们提出了一个动态并行实现的Kahn算法的拓扑排序的GPUs上。最后，对于数值因式分解阶段，我们通过消除大矩阵的内存限制来增加并行度，与现有的实现方法相比。实验结果表明，与GLU3.0修改后的实现相比，我们的核外版本实现了1.13- 32.65倍的加速比。此外，我们的核外实现实现了1.2-2.2的加速比，在GPU上的优化的统一内存实现。最后，我们表明，我们介绍的数值因式分解的优化结果是有效的。允许免费制作本作品的全部或部分的数字或硬拷贝，以供个人或课堂使用，前提是制作或分发副本的目的不是为了盈利或商业利益，并且副本的第一页上有本声明和完整的引用必须尊重作者以外的其他人拥有的本作品组件的版权允许使用学分进行摘要以其他方式复制、重新发布、在服务器上发布或重新分发到列表中，需要事先获得特定许可和/或付费。请求权限请发邮件至permissions@acm.org。PPoPP©2023版权归所有者/作者所有。授权给ACM的出版权。ACM ISBN 979-8-4007-0015-6/23/02。. . 15美元https://doi.org/10.1145/3572848.3577486彭江爱荷华大学计算机科学系美国爱荷华州爱荷华市peng-jiang@uiowa.edu拉吉夫·拉姆纳特俄亥俄州立大学计算机科学与工程系关闭WI，USAramnath@cse.ohio-state.eduCCS概念：·计算方法并行算法;线性代数算法。关键词：GPU加速，LU分解，内存限制1介绍广泛的工程和科学计算应用需要求解大型线性代数方程组，即，解决这个问题的直接方法涉及将矩阵变换为两个矩阵：下三角矩阵和上三角矩阵，使得 =。 这种转换是用一种通常称为LU分解的方法完成的。 在调用该过程之后， 可以通过求解涉及这两个三角矩阵的方程来容易地获得解，这在计算上要容易得多。对于许多应用（例如，电路仿真[16，24，30]），矩阵 可以非常大且稀疏。因此，大型稀疏矩阵的LU分解成为这些应用的关键核心。稀疏LU分解通常引入新的填充（即：非零元素的非零元素）。这些新填充的位置在运行时之前是未知的。因此，LU因式分解涉及符号因式分解阶段，其中确定这些填充的数量（和可能的位置）。此步骤之后是数值因子分解阶段，现在计算实际的非零值。由于LU分解的重要性，加速稀疏LU分解的研究已经持续了几十年[3- 5，8 - 11，15，19，27，32 ]。随着高性能计算（HPC）领域逐渐被异构和基于加速器的计算所主导，最近的几项研究工作专门考虑了GPU上的稀疏LU因子化[6，13，29]。然而，这些努力中的每一个都只加速了使用GPU的计算的一部分-具体地说，要么只是数值因式分解[ 19，28，32，36 ]。289××PPoPP或者仅仅是符号因式分解（并且进一步限制于计算非零的数量而不是它们的位置）[11]。因此，这些工作不能充分利用GPU的并行性在本文中，我们提供了第一个高效的实现，在GPU上进行稀疏LU求解器的所有步骤，同时还可以处理中间内存需求可能超过GPU上可用内存的大型矩阵。为了实现这一目标，我们要解决三个主要问题。第一个问题是在sym期间使用大量内存曲分解步骤为了将实现扩展到大矩阵，我们提出了一个核外GPU实现，迭代地执行符号因式分解。第二个问题是调度程序的并行化。作为背景，数值因式分解的实现是逐列的，并且需要检测列之间的依赖性并确定它们的顺序以进行数值因式分解。 以前的工作要么使用消除树[10，38]，要么使用分层算法[32，33]，但在CPU上部署该过程。在这项工作中，(a) 示例矩阵A（b）图形表示矩阵A（d）不同的列ID我们观察到这个调度步骤本质上是-实际上是一个拓扑排序，并提出了一个有效的GPU实现与动态并行功能的卡恩（c）矩阵A的依赖图水平算法[20]。最后，以前的GPU实现[19，32]使用矩阵的密集数据格式，以在数值因式分解步骤期间实现快速数据访问。 在这种情况下，每个列需要O（ ）内存空间，并且可以存储在GPU上（因此可以并行处理）的列数受到限制。在处理大型矩阵（即，随着温度的升高）。为了解决这个问题，我们建议在矩阵大小增加时切换到稀疏数据格式虽然数据访问（即对于给定的列ID搜索行ID）在稀疏数据格式上的效率会较低，但我们表明，我们仍然可以用基于二分搜索的访问算法实现更好的性能，因为消除了对并行列我们使用SuiteSparse矩阵集合[7]中的一组矩阵来评估我们的实现 我们将我们的观察总结如下。首先，我们的核外GPU符号因式分解可以支持大型矩阵，并在最近出版的高效CPU实现上提供1.13至32.65的加速比（即，GLU 3.0[32]）。其次，与符号因式分解的优化统一存储解决方案相比，我们的核外GPU解决方案实现了1.2-2.2范围内的加速比，因为它可以更有效地访问GPU设备内存然后，我们证明了我们对数值因式分解的优化进一步提高了大矩阵的性能高达3.33倍，这是由于内存效率导致并行性增加。图1. 示例矩阵。（a）稀疏矩阵A上的符号分解，其中第9行在分析中，红色圆圈表示新的填充（b）矩阵A的图形表示，读取线表示新的填充（c）矩阵A的依赖图，虚线将列从不同级别分隔开（d）示出不同级别的列id的表。2背景的LU分解矩阵， 具有以下形式其中是下三角矩阵， 是上三角矩阵。LU分解的实现一般涉及预处理、符号分解和数值分解的主要步骤 在预处理过程中，为了提高数值稳定性并减少 和 矩阵[9，10，27]中的填充数量，执行行和列置换。由于我们的实现主要涉及数字和符号因式分解阶段，我们在这里回顾它们。具体来说，我们首先详细介绍符号因式分解，然后概述数字因式分解的GPU实现：GLU3.0[19，32]。最后，我们简要说明了我们工作的动机。2.1符号分解与GPU实现LU因式分解涉及一系列基本行操作，其中旋转行乘以非零标量并更新为下面的另一行该另一行在旋转列中具有非零。因此，对于稀疏的290（）下一页（）下一页（）下一页（）下一页（）（）（）（）下一页（）下一页（）下一页（）下一页（）下一页（）下一页（）下一页（）下一页（）下一页（）下一页（）下一页∈（）∈（）（）下一页（）下一页（）（）（）（）ΣΣ ΣΣΣΣ（−）×（−）（−）× ×（−）（）下一页GPUs上大型矩阵的端到端LU分解PPoPP'23，2023年2月25日至3月1日，加拿大魁北克省矩阵，LU分解经常引入新的非零条目，这被称为填充。图1显示了一个示例矩阵 ，将在本文中使用。 考虑图1（a）中的行9，由于行9在列5处具有非零，所以行5需要在行9上执行行操作。结果，它产生新的填充（9，8）。如前所述，需要符号因式分解步骤来识别和矩阵的非零结构。新填充的位置由以下定理确定：定理1. 当且仅当存在从到的有向路径时，引入索引处的填充，其中中间顶点小于 和 [34]。基于定理1，已经提出了几种实现符号因式分解的算法[14，15，34]. 在本节中，我们简要总结了fill2算法，该算法具有高度并行性，因此适用于GPU。详细的过程在算法1中示出。它使用一个数组来���������表示一个已经访问过的顶点，设置为=。���������������������������������在算法开始时，它执行初始化，然后对于被视为边界的每个阈值，它检查其邻居，更新邻居的状态，并向,:or,:,as well as to the:. 随后，它将前进到第11行中的下一个阈值顶点。最近，提出了一种fill2算法的GPU实现，以加速符号因子分解[11]。 如算法1所示，来自每个源行的符号因式分解是独立的。因此，我们可以从所有列执行并行遍历然而，在这种独立遍历的情况下，每个节点都需要100万内存，这导致了100万2 的总内存需求。虽然使用分布式资源集合可以增加总可用内存（之前的研究[11]部署了多达44个节点和264个GPU），但更好的内存效率是可取的。还应该注意的是，他们的解决方案只是计算每行中新填充的数量，这对于随后的数字因式分解步骤（他们不在GPU上实现）来说是不够的。2.2数字因子分解和GPU实现传统的数字分解方法是右视LU分解：算法1基于定理1的填充2算法。此算法显示的程序为src列中的行.���输入：src-矩阵的第src行， *：，：-原始矩阵输出：填充矩阵L和U1： 0=0;2：=;3：forv in，：do4：=;5：如果v src，则6：将v加到，：;7：其他8：将v加到，：;9：如果结束10：结束11：对于阈值= 0：src-1 s. t. ���������������ℎ��������� ℎ���������==���������do十二：将阈值添加到一个新的阈值：13： 对 于 每 个 ��������������������������������������� 边 界���������������������������<���������������������������������16：=������������������������������������������17：如果邻居>阈值，则18：将邻居添加到邻居，：或邻居，：;19：其他20：将邻居添加到邻居列表中：;21：如果结束22：如果结束23：结束24：结束二十五日：swap（交换，交换，交换）26：转到13号线;27：结束方法先求解矩阵的一行 ，然后求解矩阵的一列 ，并通过迭代递归求解矩阵 。然而，这种方法具有顺序数据依赖性，这将限制并行性的量。因此，为了克服这个问题，以前的努力[19，23，32]提出了一种混合的基于列的右视算法，它可以利用列级并行。 该过程在算法2中示出。对于每一列，第一步是计算当前列的第二部分，如第2-6行所示 然后，它看起来是正确的，以查找所有列���������������������（这样的列称为列的子列���。然后，它可以并行分解的所有子列���，如第8-14行所示。1121 2211 122211 12=121 122为了安排不同列的因子分解顺序算法需要检查依赖关系在不同的列之间。例如，对于任意的n（n，n）n0，在此，X11、 X11是标量并且X11 = 1，X21是大小为X11 的列向量，X12是大小为1X11的行向量，并且X22和X22是X11X11子矩阵。为了计算矩阵和，我们首先可以发现 11= 11，12=���然后，我们可以解出22×22=���22− ���21 × 12。可以看出，传统的右视我们可以得出这样的结论，即列的ε取决于列- 是的具体地说，这是因为该列 是该列的子列 ，并且算法在对该列进行因子分解时将读取，如算法2中的第8-14行所示。还有其他的依赖关系，但为了简洁起见，我们建议读者参考早期的出版物[32]。GLU3.0然后派生291（）下一页（）下一页（）（）/（）12：如果结束（）下一页PPoPP算法2基于列的混合右视算法。输入：矩阵-符号化后A的非零填充矩阵分析1：对于=1;=;++do2：对于=j+1;= ;++do3：如果是������，则4：， =，你好，5：如果结束6：结束锻造7：对于=j+1;= ;++do8：如果是，则���为0，9：对于=j+1;=;++do10：如果（���，）0，则������11：（， ）=���（���， ）−（，）×（，）������������������������13：结束14：如果结束15：结束16：结束列的所有依赖信息，并构造依赖图。 图1（b）显示了矩阵的依赖图 。图中的边表示该列 依赖于该列 。 基于依赖图，该算法将列分组为多个级别，使得级别内的列彼此独立，从而可以并行地进行因子分解。 图1（c）显示了矩阵的级别信息 。例如，列1、2、3、6和7彼此独立，并且它们的处理可以并行。 用于确定这样的级别的过程本质上是拓扑排序，但也称为级别化。GLU3.0工作观察到潜在的并行性在各个级别上不断变化。总的来说，他们将这些级别分为三类。 在因素化的开始阶段，水平是“A型”水平。 这样的级别通常具有大量的可并行化列，而每个列具有非常少的相关联的子列。因此，它们使用一个线程块来分解一列，并将一个线程束分配给一个子列。相比之下，C型水平位于因子分解过程的末尾。在此阶段，级别具有有限数量的列，而每个列通常具有大量的子列。为了利用子列的并行性，将线程块分配给每个子列，并将内核调用而不是块分配给每个列。处于过渡阶段的B类层次，栏目数量多，同时栏目也有很多子栏目。2.3当前工作的局限性和挑战先前的研究工作证明GPU用于数值因式分解阶段[19，32]或部分符号因式分解阶段[11]。然而，还没有提出一个完整的GPU解决方案来解决LU分解。更图2. 大型稀疏矩阵的GPU LU分解的总体框架。具体来说，最密切相关的先前工作[11]提供了基于GPU的符号执行的部分解决方案，即调度阶段仅在CPU上执行此外，这项工作不集成符号和数值执行阶段。基于GPU的设备内存越来越大，越来越多的应用程序被移植到GPU上的事实，我们建议开发第一个实现这一目标的主要挑战是：1）中间阶段的大量内存需求，以及2）涉及数据依赖性的计算。3端到端GPU LU分解在本节中，我们首先介绍我们的端到端GPU解决方案的框架然后，我们展示了我们的核外GPU实现，以执行符号因式分解和GPU实现的拓扑排序的目的调度。最后，我们提出了优化，以增加parallelism大型矩阵在数字因式分解步骤。3.1总体框架稀疏LU分解的总体框架如图2所示。按照惯例，我们首先执行某些预处理步骤，即：行和列置换，目的是减少填充和提高数值稳定性。然后，我们分两个阶段执行符号因子分解。 之后，在GPU上进行并行版本的水平化-该步骤的输出用于在数值因子化阶段期间调度并行计算。最后是数值因式分解292×/（）下一页（）下一页GPUs上大型矩阵的端到端LU分解PPoPP'23，2023年2月25日至3月1日，加拿大魁北克省实现，其中我们的新贡献是当行数变得很大时切换到使用稀疏数据格式3.2基于动态矩阵分配的核外符号分解我们首先说明了内存限制的象征fac-torization。如算法1所示，与每个源行的处理相关联的是分配几个阵列 的要求，并且每个阵列需要O（N）存储器空间（其中N是行数）。因此，总共需要O（102）的内存空间，即使对于相对较小的矩阵大小，也超过了内存限制注意，原始矩阵是一个稀疏矩阵，���(a) pre2（b）audikw_1图3.每次迭代的边界大小（y轴）（x轴）你好，我很抱歉。为了保证中间数据结构驻留在GPU上，它将使用保守值������������������������并相应地限制并行度存储器需求远低于1002。解决内存限制的一个解决方案是使用最近可用的功能（至少在NVIDIAGPU的上下文中是最近的），称为统一内存[17，31]。 此功能允许应用程序透明地访问主机端的内存，并将数据加载到物理GPU内存中，同时处理页面错误。 与典型CPU操作系统上的虚拟内存概念类似，此功能可以显著简化编程，实际上，此选项最近已用于几个核心外GPU实现[1，2，12，25，26，40]。然而，事实证明，基于这种方法的实现将有显着的额外数据移动成本，特别是考虑到不规则的访问与符号分解 作为替代方案，我们关注的是一个数据移动被显式控制的版本。 我们首先提出了一种用于符号因子化的朴素核外GPU实现，如算法3所示。 在第一步中，它计算迭代次数，这是基于GPU的内存大小。假设GPU的设备内存大小为。每个源行最多需要1000���个存储空间，用于图遍历以存储中间顶点之类的值–然后，������������_������������= /（ × ）。相应地，迭代次数（表示为_）为 ���n���������_������������。我们还注意到，在我们的端到端核外GPU实现中，基于GPU的符号因式分解[11]（仅执行部分符号因式分解）的唯一其他先前工作中存在两个问题：首先，它只计算LU因式分解期间产生的新填充的总数。这对于随后的数值因式分解来说是不足够的信息，即，从算法2中可以看出，数字因式分解算法需要：1）每行的新填充的数量2）每个新填充的确切位置。第二，[11]使用固定值算法3用于符号因式分解的朴素核外GPU实现1：num_iter = n/chunk_size2：计算每一行中的填充数3： foriter=0;iternum_iter;iter++do4 ： symbolic_1<<<> 结 束������������符号_字符串5：end for6：记录每一行中非零的个数������������7：前缀_和（前缀_和_前缀_和_和）8：为分解矩阵分配内存。9：计算内核symbolic_2中填充的位置10：foriter = 0; iter num_iter; iter++do11：symbolic_2<<<>象征符号���_象征符号12：结束为了解决[11]中的第一个问题，我们的实现包括两个阶段。 对于第一阶段，我们只计算分解矩阵每行中非零的数量。在每次迭代期间，我们启动内核计算器，以并行计算 非 零 行 的 非 零 计 数 ， 如 算 法 3 中 的 第 3-5 行 所 示 。������������������������根据最近的出版物[11]修改了密码锁密码锁密码锁_1的实现每行中非零的计数存储在数组������������-- 由于我们使用压缩稀疏行（CSR）数据格式来存储因子分解矩阵，因此我们在array_上应用前缀和的GPU实现������������，以获得每行的起始位置和非零的总数（第7行）。有了这些信息，我们就可以为分解矩阵分配设备内存（第8行）。最后，对于第二阶段（第10-12行），在每次迭代中，我们都会������������为每个行启动内核与_1相比，_2实现的主要区别在于，一旦我们找到一个填充，我们也会存储它的位置。这也是为什么需要提前执行数据库管理器1以找到空间需求的原因。对于第二个问题，我们观察到每个源行的内存需求随着源行标识符值的增加而增加。更确切地说，这是293（）下一页（）下一页PPoPP定理1，即，由于中间顶点需要具有小的值，因此随着源行标识符值的增加，将有更多数量的中间顶点需要考虑 我们使用示例矩阵pre 2和audikw_1验证了这一点，结果如图3所示-具有比源行标识符更小的标识符并且需要主动考虑的中间顶点被表示为边界。从图中可以看出，边界的数量通常在最后几次迭代中很大，否则很小。因此，我们提出了一个动态的并行分配实现，它被概括为算法4。我们首先将行划分为两部分：第一部分包含1001行，第二部分包含剩余的（1002）行 关键的区别在于，我们为这两个部分计算并使用不同的���k���������_������������。具体来说，对于第一部分，对应于 每 一 行 的内存需求较 小 ， 因 此 我 们 分 配 一 个 大 的���k���������_������������来增加并行性。然后，我们计算每个部分的迭代次数，如算法中的第1-2行所示最后，对于每个部分，我们迭代地启动内核_1来计数非零，这分别在第4-6行和第8-9行中显示。第二阶段的过程是类似的，在算法中省略。请注意，在执行此优化时，可以探索使用2个以上的阶段，但这也意味着更多的内核启动。算法4用于具有动态并行分配的符号因式分解的1：num_iter_1= n1 /chunk_size_12：num_iter_2 = n2/chunk_size_23：计算每一行中的填充数4： foriter=0;iternum_iter_1;iter++do5：symbolic_1<<<>结束������������字符串_16：end for七：……8： foriter=0;iternum_iter_2;iter++do9 ：symbolic_1<<<>结束������������字符串_210：end for11：第二阶段是类似的。12：...3.3GPU上的并行调度过程基于依赖图，以前的工作计算每列的水平数如下：���������������( )=���������(−1,���������������( 1),���������������( 2),)+ 1其中 1， 2，. . . 是节点的子节点 。该过程本质上是串行的，因为不同列之间存在相关性：列E1的级别编号取决于列E1、 E2的级别编号，依此类推。因此，在本发明中，先前对LU因式分解的努力都在CPU上执行层次化，因此没有在GPU上实现端到端LU因式分解。为了在GPU上并行化层次化，这本质上是一种拓扑排序，我们首先注意到以前将该内核映射到GPU的努力的局限性已经有一些通用的GPU图形处理系统[21，41]报告它们可以支持拓扑排序。然而，它们都没有明确优化这个算法。一些出版物[37]明确报告了GPU拓扑排序的实现，但它们使用CPU来启动内核，因此无法充分利用GPU提供的并行性。在算法5GPU上的并行层次化实现。一曰：全球public void run（）{2：level_num = 03：队列 _队列_队列是所有没有传入边的节点的队列4：cons_queue<<<>（...）5：level_num++6：whileqsize > 0do7：update<<<>（...）;return 0;9：cons_<<（.）10：level_num++;11：结束时第十二章：十三日： 程序级别IzA tI on14：cons_graph<<<>（...）15：cnt_indegree<<<>（.）16：拓扑<<<>（...）;17：结束程序在这项工作中，我们提出了一个纯粹的GPU实现与动态并行功能的CUDA。 与[37]相比，我们的动态并行实现具有以下优点：首先，它避免了CPU和GPU之间的同步和数据传输。其次，通过在GPU内调用函数，内核启动开销大大减少。详细的过程在算法5中示出。该方法假设���已经基于第2节（第14行）中提到的方法构建了依赖图。然后，我们使用内核cnt_indegree计算每个节点的入度（第15行）。实际的拓扑排序过程在我们启动内核Topo_Sort时开始（第16行）。在内核内部，我们首先创建d_queue数据结构，它表示要处理的节点。最初，这个集合包括所有没有传入节点的节点-更具体地说，该内核检查每个节点的入度，如果节点的入度为零，则该方法将该节点放入d_queue并将级别数设置为0。在循环中的每次迭代期间，启动子内核更新以更新d_queue中节点的邻居的入度值。然后，我们为这样的邻居构造一个新的d_queue（cons_queue过程，第9行）。我们提高水平294（+）×（）下一页（）下一页GPUs上大型矩阵的端到端LU分解PPoPP'23，2023年2月25日至3月1日，加拿大魁北克省第10行中每次迭代结束时的编号我们重复这个过程，直到没有节点没有传入节点，即_与我们意识到的这个问题的现有工作相比[37]，我们的改进在于在GPU内调用函数，而不是使用CPU来启动内核。虽然由于基线代码不可用而无法进行直接比较，但我们可以预期，由于内核启动开销被删除，因此会有显著的改进。 序贯拓扑排序的计算复杂度是，是节点数， 是边数。并行拓扑排序的跨度（最长执行步骤）是层次的数量。算法6二叉搜索以访问As（i，j）。输入：列偏移_列偏移-CSC格式的列偏移CSC格式的行ID表1. Nvidia Tesla V100个gpuTesla V100#SM80FP 32 CUDA核心/GPU5120存储器接口4096位HBM2寄存器文件大小/SM（KB）65536最大寄存器/线程数255共享内存大小/SM（KB）可配置高达96 KB最大线程块大小1024越来越大最后，as 变得非常大， 可能小于并发线程的最大数量（表示为_）。这将导致实现无法利用GPU的足够并行性。为了解决这个问题，我们建议采用压缩稀疏列（CSC）数据格式来处理数据，其中我们确定-CSC格式的值证明大于_×. 查尔-一：……2：fs = col_offset[j]3：fe = col_offset[j+1]4：当fe >= fs时，5：mid =（fs + fe）/26：如果row_id[mid]== i），则7：As（i，k）= As（i，k）-val[mid] As（j，k）8：休息9：else ifrow_id[mid]> ithen10：fe = mid-111：其他12：fs = mid +113：如果结束14：结束时3.4通过消除数值分解的内存限制来提高并行性先前在GPU上实现的数值因式分解- 特别是GLU实现[19，32，33]������������为了进一步说明，在算法2中，我们需要搜索一个行idn，它大于列IDn。因此，当我们使用密集格式时，我们可以有效地访问数据，因为位置是直接的 。然而，我们观察到这增加了总的内存需求，限制了可以存储在块中的行数，从而减少了并行度。 假设总可用设备数为 ，最大可并行列数 可计算为：=���������������������×，其中n是顶点的数量由于我们使用线程块来执行一列的数值因子分解然而，这种变化的挑战在于，对于给定的列IDn，在这种情况下，我们不能直接获得大于列IDn的行IDn因此，我们利用CSC格式中的升序属性，并执行二进制搜索以找到行ID大于列ID1的位置。详细的过程在算法6中示出。 在该算法中， 表示列id， 表示行id，并且_[j]和_[j+1]之间的索引被排序。表示要搜索的最小可能索引，其被初始化为_[j]，并且表示要搜索的最大位置索引，其被初始化为_[j+1]。 在每次迭代中，我们将索引的中间值，，与 。 如果中的行id与相同 ，则我们已经找到了索引，即。否则，如果数据库中的行id大于 ，则数据库将更新为���������1. 或者，如果行id在“行”中小于“行”，������is updated as���������− 1,and we continue the search.4评估和绩效研究在本节中，我们将展示我们在一组大型矩阵上的实验结果，以证明我们的端到端GPU实现的有效性。我们首先展示了我们的实验环境和所选输入矩阵的特征然后，为了显示我们的符号因式分解的核外GPU实现的有效性，我们将其与从GLU修改的并行实现进行比较。3.0 [32]和优化的统一内存实现。接下来，我们评估从已经引入的优化中获得的好处，包括在最大矩阵的数值因式分解过程中使用稀疏数据表示。表示最大可能并发线程块。 作为从这个表达式中可以看出， 当1时会更小注意，我们的CSC格式是排序的。295//PPoPP表2.符号因式分解的内存需求超过GPU内存大小的输入矩阵。矩阵abbrnnnznnz/ng7jac200scG75931083793614.1rma10RM46835237400150.7pre2PR65903359592829.0inline_1在50371218660027 37.0曲轴段_2CR2638387106348111.3bmwcra_1BMC148770539638636.3曲轴段_1CR1528045333507101.0bmw7st_1Bm7141347374050726.5Apache2AP71517627665233.9s3dkq4m2S3490449245567027.1s3dkt3m2S3390449192195521.2一元2OT2360572276286.3rajat15R153726144357311.9bbmatBB38744177172245.7mixtank_newMI29957199504166.6古德温_054去32510103087831.7一元1OT1360573410889.5windtunnel_3dWi40816273060066.94.1实验设计环境：我们在Nvidia Tesla V100上进行实验。GPU的规格如表1所示。GPU连接到运行在2.4 GHz的Intel（R）Xeon（R）CPU E5-2680（2013 Ivy Bridge）-CPU包含14个物理内核，并提供每个内核2个线程的超线程，用于我们的在我们的实验中，主机内存的大小是128 GB。我们实验的主机操作系统是CentOS Linux版本7.4.1708（核心）。我们的GPU实现基于CUDA 11.2工具包，NVCCV11.2.152用于编译我们的程序。输入矩阵：我们从SuiteSparse Matrix Collection [7]中选择了18个矩阵进行详细的研究和分析。之所以选择这些矩阵，是因为可以对这些矩阵进行LU分解，并且正如我们所验证的那样，中间数据结构的内存要求超过了Nvidia Tesla V100的设备内存大小。换句话说，对于这些矩阵中的每一个，如果没有显式的数据移动或使用统一内存，就不能在GPU上执行符号因式分解。基质的质量标准见表2。我们的实验使用float作为数据类型。基线的选择：我们主要比较了我们的核外GPU实现与从GLU3.0[ 32 ]修改的并行实现的性能，GLU3.0是最近的高效实现。 对于符号部分，因为我们的实现直接从Gaihre等人的代码开始。 [11]并在功能上改进它（即，计算数值相位的非零位置和防止较大数据的存储器不足）而不是性能，则性能的比较将没有意义。其他数值阶段优化的努力并不总是使代码可用，此外，数值阶段是完整代码的执行时间的一小部分最后我们已经广泛地与另一种设计方案进行了比较这种设计方案是使用统一内存。4.2与修改的GLU3.0实现的比较比较结果如图4所示。执行时间分别按符号因式分解和数值因式分解所花费的时间进行分解。 我们注意到加速比（对于整个执行）在1.13-32.65之间。还可以看出，我们的核外GPU实现与GLU3.0实现之间的差异主要来自符号因子分解阶段。虽然多核CPU和GPU之间的相对性能差异很大，但GPU加速似乎取决于每行非零的数量 。 当比率较大时，加速比往往较大。例如，矩阵WI和MI都具有最高的比值Δ P/ΔT，并且在最高的加速比中，其中AP和OT 2在相反的频谱上这与一般观察结果一致，即GPU随着计算变得（相对）密集而变得更高效。4.3与统一内存解决方案的我们进一步比较我们的核心GPU实现与统一内存实现。首先，我们注意到，即使是统一内存解决方案也会受到CPU主内存大小的限制。因此，对于这个实验，我们选择了18个矩阵中的7个，对于这些矩阵，中间数据大小可以适合CPU主内存，但不适合GPU设备内存。具体地，这些是表2中具有7个最小值的矩阵，所有都具有少于41，000行。对于我们的实验，我们进一步调优了统一内存实现，并尝试了不同的优化。我们发现预取中间数据结构可以提高效率。比较我们的实现与统一内存实现（启用预取）的结果如图5所示。如前图所示，执行时间分别按符号和数字因式分解阶段所花费的时间进行分解 我们看到，我们的核外GPU实现比优化的统一内存实现效率高1.06-2.22倍。更仔细的分析表明，对于密度相对较高的矩阵，即，WI和MI，统一内存实现是相当有竞争力的。另一方面，对于具有最低密度的R15和OT2，统一存储器开销较大。 这符合我们的预期，即由于计算较少，页面错误的影响会更大。我们进一步创建了一个版本的统一内存imple-没有任何预取的mentation我们比较了我们的核外GPU实现与统一的基于内存的实现的符号因式分解阶段的执行时间，如图6所示。如图所示，如果没有预取，统一内存实现的性能会更差。他们的296GPUs上大型矩阵的端到端LU分解PPoPP'23，2023年2月25日至3月1日，加拿大魁北克省图4.核心外GPU实现和修改后的GLU3.0基线的标准化端到端执行时间（符号和数值阶段分开的时间）。图5.针对核外GPU实施和统一内存GPU基线的标准化端到端执行时间（符号和数值阶段分开的时间）图6.我们的核外GPU实现的标准化符号执行阶段时间，有和没有预取的统一内存实现。对于密度较低的矩阵，如R15和OT2，相对性能变差关于三个版本之间的性能差异的进一步阐述见表3。主要统一存储器实现中的按需分页的性能缺点是GPU页面错误的开销。如表中所示，统一内存版本的大量时间用于处理页面错误，而297（）下一页PPoPP表3. 比较GPU页面错误组的数量以及在不使用预取和使用预取的情况下处理GPU页面错误的时间百分比。WP表示具有预取，而WOP表示不具有预取。矩阵GPU故障数量故障wpPC. wo p（%）PC. 可湿性粉剂（%）PC. ooc（%）OT216734463878.3756.600.06R1517322439286.2165.460.15BB19753579846.9826.940.09MI12803437736.1521.610.10去13670384878.1957.390.33OT116884471769.5845.180.06Wi24977856933.1119.540.01表4. 原始版本的大矩阵规格和最大并行线程块数矩阵秩序nnz最大块数hugetrace-0002016,002,41347,997,626124Delaunay_n2416,777,216100,663,202119hugebubbles-00000 18,318,14354,940,162109产品编号：0001019,458,08758,359,528102图8.我们的二叉搜索实现和原始实现的归一化数字因式分解时间。我们在其他实验中使用过的。具体而言，这些矩阵如表4所示，���大于_×.由于这些矩阵图7.我们的动态并行分配实现和原始符号因子分解实现的执行时间。核外GPU实现在数据移动上花费非常少量的时间。回头参考图6，我们还观察到，对于具有显著计算开销的矩阵，例如MI和WI（它们恰好是更密集的矩阵），用于服务GPU页面错误的时间百分比较小，并且相应地，我们的实现的益处也较小。4.4优化评估动态并行分配：我们进一步证明了动态并行分配的有效性和局限性。我们比较我们的动态并行分配实现与本机的核心外实现符号因式分解两个大矩阵。选择这些矩阵是因为它们很大，迭代次数也很大。这种比较如图7所示。我们观察到，动态实现比朴素实现的性能提高了10%我们还注意到，性能改进受到实现的限制。在符号因式分解的某些步骤中，并行度由边界的个数决定因此，当边界的数量非常大时，这些步骤的性能改进将受到限制。数值分解的内存优化：我们接下来证明了我们的优化的有效性，以增加并行的数值因式分解。这些好处只在非常大的矩阵中才能看到，碰巧不是LU可因子分解的（它们不是满秩的）。在我们的实验中，我们用一个非零的数字（1000）替换了它们的0对角元素，使它们可分解。表4还显示了这些矩阵的并行线程块的最大数量因为我们的GPU的线程块的最大数量是160，原始的数值因子化实现不能利用我们的GPU2上的完全并行性。接下来，我们将二进制搜索实现的执行时间与图8中的原始实现进行了比较。从图中可以看出，由于并行列的数量增加，我们的二进制5相关工作在本节中，我们将讨论加速LU分解的相关工作，以及在图处理和线性代数的背景下提出的某些相关的核外实现。加 速 LU 分 解 ： 受 超 节 点 方 法 在 加 速 对 称 正 定 矩 阵（SPD）的方面的成功启发，提出了用于非对称矩阵的超节点LU方法[9，10，27]。 沿着这些路线，Demmelet al. [9]介绍了五种