ACM通信
研究突出了
一种面向异构体系结构的大规模并行自适应快速多极方法
来源:法国
我们描述了一种并行快速多极方法(FMM)的高度非均匀分布的粒子。我们使用分布式内存并行(通过MPI)和共享内存并行(通过OpenMP和GPU加速)在异构高性能计算架构上快速计算三维二体非振荡势。我们在ORNL基于AMD/ crayen的Jaguar系统上，在196608个核上使用多达300亿个粒子进行了可伸缩性测试。在一个gpu支持的系统上(美国国家科学基金会在乔治亚理工学院/ORNL的Keeneland)，我们观察到单核CPU的速度提高了30倍，多核CPU实现的速度提高了7倍。通过将gpu与MPI相结合，我们在192个gpu上使用4800万个非均匀分布的粒子时，实现了小于10纳斯/粒子和6位数的精度。
1.简介
N-体问题是数学物理中各种问题的计算基石，16机器学习,8和近似理论。4从形式上讲，问题是如何有效地进行评估
这本质上是O（N2)之间所有成对相互作用的计算N粒子。我们将f随着“潜力”随着“源密度,“x作为目标粒子的坐标，y为源粒子坐标，K为相互作用核。例如，在静电或引力相互作用中
与R = x y和|r|是欧几里德范数。式(1)表示n体问题．
挑战在于O（N2)的复杂性使得计算对于实际问题的规模来说代价高昂。例如，了解血液流动机制21或者理解宇宙的起源26需要进行数十亿粒子相互作用的模拟。
鉴于…的重要性N人们已经付出了相当大的努力来加速这些问题。的大小和方差K（r)与|一起衰变r|，所以当粒子或粒子组被很好地分离时，相互作用可以用很小的计算代价近似。但细节是棘手的。我们如何进行近似呢?我们能否提供严格的误差和复杂性界限?我们如何有效地实现这些近似呢?
第一批成功的快速计划尝试之一是BarnesHut方法1这使得复杂性降到O（N日志N)．格林加德和Rokhlin3.，10，22解决了这个问题。他们引入了快速多极方法(FMM)，该方法使用了快速收敛的近似方案，具有严格的误差边界，并将总体复杂性(在三维空间)降低到O(日志(1/7
）3.N),是精确和的近似值的期望精度。这个结果令人印象深刻:我们可以求和的值方程(1)线性时间的任意精度．
那么FMM是如何工作的呢?这是一张草图。给定我们想要计算所有成对相互作用的粒子的位置，我们构建了一个八叉树(从今以后我们只考虑三维的FMM)，这样任何叶节点(或八叉树)都有大约规定数量的粒子;然后我们执行树遍历来计算和。我们执行一个海报(自底向上)遍历，然后执行一个预订遍历(自顶向下)。海报遍历由涉及一个八分区及其子分区的计算组成。预订遍历更加复杂，因为它涉及到被访问的八分区周围的邻域中的几个八分区。
在这个抽象级别上，我们开始看到如何并行化FMM。当粒子的分布是均匀的时，分析就相对简单:我们以一层一层的方式执行遍历，在每一层中使用跨分区的数据并行性。最长的依赖链大概是O(日志N)．挑战在于粒子分布不均匀的情况下，在这种情况下，即使是连续的情况也会变得非常复杂。3.
方法总结:我们描述并评估了用于大规模并行分布式内存系统的FMM实现，包括基于多核CPU和GPU协处理器的异构架构。我们使用消息传递接口(MPI)来表示分布式内存并行，使用OpenMP和CUDA来表示共享内存和细粒度并行。根据核和相互作用近似格式的不同，FMM有多种变体。我们的实现基于kernel-independent FMM．28，29然而，这种选择只影响算法中使用的某些运算符的构造。它不会影响算法的整体“流程”或它的分布式内存并行化。(它确实会影响单核性能，因为向量化和细粒度并行依赖于FMM方案。)
我们的并行算法可以总结如下。给定在MPI过程中以任意方式分布的粒子及其相关的源密度，我们试图计算每个点的电位。(为简单起见，本文假设源点和目标点重合。)首先，我们按照莫顿顺序对点进行排序，24一种空间填充曲线，然后重新分配它们，这样每个进程都拥有排序数组中的一个连续块。其次，我们为并行的每个流程创建所谓的局部基本树(LET，在第3节中定义)。第三，我们在MPI进程中使用let和在每个进程中使用GPU和OpenMP加速来计算总和。算法的所有步骤都是并行的，包括树的构造和评估。
相关工作:有相当多的关于并行算法的文献N身体的问题。对该算法进行了深入分析O（N)工作,O(日志N)深度的粒子分布，导致树木的深度被O(日志N)．Greengard和Gropp在粒子分布均匀的情况下，使用并行读、独占写(CREW) PRAM模型分析并实现了该算法。9使用逐级遍历树和使用数据并行性，实现非常简单。该算法的复杂度，省略了与精度相关的项，为O（N / p) +O(日志p),p是线程数。这个结果扩展到消息传递实现。
当处理粒子的非均匀分布时，实现、并行可伸缩性和复杂性分析变得更加困难。卡拉汉和Kosaraju2提出一种专属读专属写(EREW)-PRAM算法p = N处理器是最优工作状态，记录日志p时间，不取决于点的分布。据我们所知，该算法尚未在实践中使用。我们的算法使用了许多最早出现在Warren和Salmon开创性工作中的思想，他们引入了空间填充曲线，使用morton排序在处理器之间划分点，局部基本树(LET)的概念，以及并行树构造算法。27这些思想在其他作品中被广泛采用和分析。7，23
Teng是第一个提供复杂性分析的人，该分析考虑了树深度为log所限制的粒子分布的通信成本N(在有限浮点精度下，所有粒子分布都满足这个假设)。25关键思想是构建一个通信图，其中图节点对应于树的分区，而边对应于分区之间的交互。Teng证明了FMM通信图有一个等分线，它的边缘切割(在三维)是由O（N2／3(日志N）4/3）;为了进行比较，均匀网格的平分线为O（N2／3)．这一结果表明，可扩展的FMM计算在理论上是可能的。Teng概述了一种分治并行算法来构建和划分树。morton顺序排序用于构建树，并使用几何图分区方法进行重新分区，以保证低边切和良好的负载平衡。18假设(1)使用并行基数排序来划分点(2)FMM树的深度以对数增长N，整个算法的规模为O（N / p) +O(日志p)．复杂度估计中的常量涉及触发器率、内存延迟和底层机器架构的带宽参数。
Kurzak和Pettitt等人更关注实际实现和性能分析15和Ogata等人。19有几个团队一直在研究GPU加速。5，11，12，20.，30.最后，对于几乎均匀分布的粒子，可以使用快速傅里叶变换(粒子在单元中的方法)来计算在O（N日志N)时间。
贡献:据我们所知，目前还没有FMM实现能够在数千个MPI进程上演示高度不均匀粒子分布的可伸缩性。我们的主要贡献是产生并演示这种利用所有级别并行性的实现的可伸缩性。在第3节中，我们讨论了分布式内存并行性，并介绍了一种新的基于超立方体路由的特定于fmm的全约简算法。在第4节中，我们将介绍实现的共享内存并行性，在第5节中，我们将讨论我们的实验和可伸缩性结果。
2.FMM轮廓
在本节中，我们将描述FMM数据结构、主要算法步骤以及该方法的并行可伸缩性。在不丧失本文其余部分的一般性的前提下，我们假设源粒子yj目标粒子x我一致的。FMM树结构确保每个叶八分之一包含不超过问粒子。顺序算法很简单:根八分单位被选为包含所有粒子的立方体;我们将粒子一个一个地插入树中;当一个八分位数包含超过时，我们将它细分问粒子。并行结构更加复杂，我们将在第3节中描述它。
在树形结构完成后，我们需要为每一个八分区构造它交互列表．这些列表包含树中的其他分区。每个列表都代表一个精确定义的空间邻域，对于每个列表，我们必须计算该列表和该列表中的分区之间的“相互作用”。通过“交互”，我们指的是对应于矩阵矩阵或矩阵向量乘法的浮点运算。一旦构建了交互列表，树将被遍历两次:首先是自底向上，然后是自顶向下。在对每个八分区的每次遍历中，我们执行涉及其交互列表中的其他八分区的计算。
FMM编码的一个难点是如何有效地用数学表示和实现辛量相互作用矩阵，这样我们就可以在不影响精度的情况下实现算法效率。矩阵的大小各不相同，但根据实现的不同，典型的尺寸为1001000。FMM算法设计者寻求减少矩阵的大小，使用对称性来减少存储需求，使用预计算，并使用快速逼近(例如，快速傅里叶变换或低秩奇异值分解)。我们需要记住的是，一旦我们确定了所需的精度，就可以完成一个八分之一与其交互列表中的八分之一之间的交互计算O(1)时间。
现在让我们给出这些列表的精确定义。对于树中的每个八分区，我们创建四个八分区列表，即所谓的U-、V-、W-和x -列表。3.u列表只定义叶分区。对于叶八分位，的U-list由所有相邻的叶八分位组成，包括它自己。(我们说if和共享一个顶点、一条边或一个面。)对于任何八分之一，它的同事们定义为位于同一树级别的相邻分区。一个八分域(叶或非叶)的V-list由父八分域的同事的孩子组成，P()，它们不相邻。W-list仅为叶八分位创建，并且包含一个八分位，当且仅当它是的同事的后代，不相邻，且其父的相邻。一个八分之一的x表由那些在w表上的八分之一组成。这些定义总结在表1．
对于每一个八分仪T，我们存储它的级别l两个向量,u而且d．的u矢量是由内部源密度产生的电势的近似表示。这种表示是足够准确的，只有当评估粒子外卷由和同事所覆盖．为u向量，我们将使用“向上等效密度”或简单地称为“向上密度”。这个术语是由Ying等人提出的。28，29
的d矢量是由源密度产生的电势的近似表示外卷由和同事所覆盖．这种近似表示是足够准确的，只有当评估粒子包围．为d向量，我们将使用“向下等效密度”或简单地称为“向下密度”。对于叶分区(l)，我们也储存x,年代,f．根据上述定义，式(1)中求和的近似求值需要向上计算一次T，然后向下计算T，在算法1中可以更具体地看到。让我们重申的细节“互动”算法1中的步骤取决于列表类型和底层FMM方案。
2.1.并行化评估阶段
通过使用分布式和共享内存编程模型，可以理解算法1中的依赖关系并设计有效的并行变体。然后可以使用MPI、OpenMP和GPU应用程序编程接口实现这些模型。
在FMM中有多个并发级别:跨步骤(例如，S2U和ULI步骤没有依赖关系)，在步骤内(例如，在VLI步骤中减少一个八分位数)，以及每八分位数计算的向量化。事实上，可以通过管道连接S2U、U2U和VLI步骤来公开更多的并行性。我们没有利用管道，因为它降低了实现的模块化和通用性。
算法1中计算的一般依赖关系如下:APPROXIMATE interaction(步骤(5a)和(5b)除外)和DIRECT interaction部分可以并发执行。对于近似的交互计算，步骤的顺序表示依赖性，例如，步骤(2)必须在步骤(1)完成后开始。步骤(3a)和(3b)可以以任何顺序执行。然而，步骤(3a)和步骤(3b)的并发执行需要有累加的并发写操作。对于步骤(5a)和步骤(5b)也是如此，它们是独立于并发写入的。
我们的总体策略是使用基于mpi的分布式内存并行性来构建全局树，并将其划分为重叠的子树(let)，以便删除步骤(3a)/(5a)和(3b)/(5b)之间的依赖关系。我们显式地处理并发写入，并在每个MPI进程中使用gpu上基于共享内存的并行，以加速算法1中间接交互的直接交互和步骤(1)、(3)、(5)。在下面的小节中，我们将详细介绍我们的方法。
3.分布式内存并行性
我们的方案的主要组成部分是(1)树结构，其中构建全局FMM树，每个进程接收它的本地基本树和一组叶子，它承担所有权和(2)评价，其中每个进程计算它所拥有的叶分区粒子的和。这个总和有两个组成部分:直接交互(精确评估)和间接交互(近似评估)。
输入由粒子和它们的源密度组成。算法的输出是粒子处的势。在描述算法之前，我们需要以下定义:27
局部基本树(LET):给定一个l跨流程，因此lk叶分区集是否分配给流程k(即，这些区域的潜力是通过过程计算的k)，用LET表示过程k定义为所有拥有叶及其祖先的交互列表的并集:
的基本思想27FMM算法的分布式内存实现是跨进程划分FMM树的叶，构造每个进程的LET，然后计算N-平行体和。该方法有两个通信密集的阶段:第一个阶段是LET构建，第二个阶段是评估期间的全减少。接下来，我们讨论这两个阶段的主要组成部分。
3.1.树结构
树构造过程的输入是粒子。输出是每个过程的局部本质树，随后在计算中使用它，以及跨MPI过程的单元立方体的几何域分解。后者贯穿整个算法。树的构造包括(1)构造一个只包含叶的八分之一的分布式线性完整八叉树(2)构造每个过程的let。
我们首先创建分布式的全局morton排序数组，其中包含全局树中的所有叶节点。这个任务的算法是已知的(我们使用Sundar等人描述的算法;24参见Hariharan和S. Aluru13)，它们的主要成分是平行莫顿类粒子。这种排序决定了树结构的总体复杂性。
叶子在进程之间的分布引出了单元立方体的几何划分:每个进程控制它所拥有的叶子覆盖的体积。通过k，我们将表示过程“控制”的区域k．每个进程存储整体的几何分区:我们使用MPI_AllGather交换其区域的第一个和最后一个分区，这是定义分区所需要的全部内容。
接下来，每个进程将所有祖先分区添加到其本地叶中，从而创建本地树。交换“幽灵”分区信息的任务仍有待完成。为此，让我们介绍以下定义:
- “贡献者”八分音符的过程师:
Pc(): =k1……p:与k - “用户”八分音符的过程:
Pu(): =k1……p: p()或C（P())重叠k．
让我乐的是进程所对应的八分之一的集合k贡献和哪个过程k”用途。然后,过程k必须把所有分区都派进去吗我乐的MPI过程k”．图2提供此过程的说明。也就是说，每个八分之一分区的每个进程贡献者将八分之一数据发送给该八分之一分区的所有用户。曾经的交换我乐的列表完成后，所有MPI进程将收到的分区插入到它们的本地树中。这就结束了每个流程let的构造。我们总结了算法2中局部本质树的构造。
构造的正确性建立在LET和FMM之间的直接关系上。考虑由某些八分之一围合的源所产生的电位。为了评估这个势能外的体积所覆盖C（P())，人们不需要关于源或相关的上行密度的信息，因为会使用的某些祖先的上行密度。这个观察结果可以用来形式化LET结构的正确性。参见Lashuk等人。17以获取更详细的讨论。
在算法2的最后一步中，每个进程都独立地为其LET中的分域元素构建U、V、W和X列表，这些分域元素包含局部粒子(在那里要评估势)。LET中已经存在所有必要的分区，因此在此步骤中不需要进一步的通信。
3.2.负载平衡
在对非均匀八叉树进行交互评估时，为每个进程分配相等的叶块可能会导致严重的负载不平衡。为了克服这个困难，我们使用以下负载平衡算法。
在LET设置之后，根据与它的U、V、W和X列表相关的计算工作(这些列表已经在这个粒子上构建)，为每个叶分配一个权重。然后，我们重新划分叶，以确保每个进程拥有的叶的总权重近似相等。我们使用Sundar等人的算法1。24以并行方式执行此重新分区。注意，与上一节类似，每个进程都获得一个连续的全局(分布式)morton叶数组块。在最后一步中，我们在每个进程上重新构建LET和U、V、W和X列表。
注意，我们完全基于工作平衡重新划分叶，忽略了通信成本。这种方法不是最优的，但计算成本不高，在实践中工作得相当好。除了基于叶的分区，还可以考虑在更粗的级别上进行分区;25我们还没有测试过这种方法。
3.3.FMM评价
潜在评估需要三个沟通步骤。第一个是传递U-list计算的源密度。这种通信是“本地的”，从某种意义上说，每个进程通常只与它的空间邻居通信。因此，一个简单的实现(例如，使用MPI_Isend)是可以接受的。
第二个通信步骤是将每个八分区的所有贡献者的向上密度相加，即U2U步骤。在此步骤(自底向上遍历)中，每个流程只构建部分它的LET中八分之一的上升密度。一个八分之一的部分上升密度不包括属于其他MPI过程的八分之一的后代的贡献。因此，我们需要第三步，即通信步骤，来收集每个八分区用户的向上密度。该通信步骤必须发生在U2U步骤之后，VLI和XLI步骤之前。此步骤完成后，每个流程都执行其LET的自顶向下遍历，而不与其他流程通信。
算法3描述了一种将上述第二步和第三步结合起来的通信过程。该算法类似于超立方体上的标准“AllReduce”算法。下面是算法的非正式描述。
开始时，每个处理器r形成一个池年代是“共享的”，也就是说，不被使用的八分之一(和它们向上的密度)r一个人。然后，MPI通信器(所有MPI进程的集合)被分成两半(为了简单起见，我们假设通信器的大小是2的幂)，来自前一半的每个处理器与来自另一半的“对等”处理器配对。
现在考虑这些区域年代它们被通信器的另一半的某些处理器“使用”(不一定是对等的?r)．如果存在这样的区域，r将它们(连同向上的密度)发送给它的“同伴”。根据同样的规则，“peer”将发送一些(可能没有)分区和密度给r．接收到的分区将合并为年代，消除重复的分区求和重复分区的密度。
在这一点上，两个部分之间不需要进一步的通信，算法通过将自己应用到每一个部分来递归地进行。我们最后注意到，在每一轮沟通之后，年代是否从不再需要通信和不在本地使用的“瞬态”分区中清除r．
该算法的时间复杂度不低于O（
)．更具体地说，假设没有进程使用超过米共享分区，没有过程贡献超过米对于一个超立方体互联，算法3的通信复杂度为O（t年代日志p + tw米(32)),在那里t年代而且tw分别是延迟和带宽常数。参见Lashuk等人。17对于一个证明。
经过这三个通信步骤后，算法1的所有剩余步骤都可以进行，无需再进行通信。
3.4.粒子均匀分布的复杂性
我们有所有必要的成分来推导粒子均匀分布的分布式内存算法的整体复杂性。让N是粒子的数量，让p是进程的数量。分域数与粒子数成正比。第一个通信代价与输入粒子的并行排序有关。我们用比较排序代替了装箱，所以它的时间复杂度为
(样本排序和双元排序的组合)。6交换“幽灵”分区的复杂性与3.3节中描述的减少广播算法相同，即:O（米),米是两个进程之间“共享”分区的最大数目。对于一个统一的网格，米可以估计为O（(N / p）2／3) /p．通信还包括源密度的交换。假设系统的带宽相当高，我们可以忽略所有低阶通信项。总而言之(假设带宽足够大，但不考虑延迟)，设置阶段的总体复杂性是
．对于评估，我们有
．
对于非均匀粒子分布，设置和评估阶段的复杂性估计包括一个额外的tw米术语。因为我们没有一个(非平凡的)边界米，这个结果比Teng算法理论上可能得到的结果要差。(他使用工作和通信成本来划分树，而不是只使用叶分区的工作。)
4.Intra-Node并行性
提高整个算法效率的一个关键因素是高度调优的节点内性能。许多当前的超级计算系统具有异构节点，这意味着它们既有传统的通用多核CPU处理器，也有更专门的高速协处理器，如gpu。因此，我们的实现可以利用cpu(通过OpenMP)和gpu(通过CUDA)。
共享内存并行化对于这两种类型的处理器都是常见的。对于S2U、D2T、ULI、WLI、VLI和XLI步骤，我们可以以一种令人尴尬的并行方式访问分区。此外，所有的八分区访问包括八分区之间或八分区与粒子之间的相互作用，表示为矩阵向量乘法。除了VLI计算外，矩阵是密集的，它对应于对角平移。因此，总体上有两个级别的并行性:跨分区和跨对应矩阵的行。相同的并行方法适用于CPU和GPU协同处理器，因为它们具有本质上相同的用于利用并行的架构特征:共享内存地址空间、多层内存层次结构和用于常规数据并行的向量单元。U2U和D2D在我们当前的实现中保持顺序，尽管原则上可以使用rake和压缩方法并行化它们。14
虽然cpu和gpu之间的架构特性是相似的，但是它们之间的差异规模直接影响算法和实现。FMM的主要算法调优参数是每八分单位的最大粒子数，问．增加问使叶分区更大，而整个树更短，这反过来增加了计算绑定的ULI、WLI和XLI步骤执行的触发器数量，同时减少了内存绑定的VLI和其他阶段的触发器成本。图3(一个)显示了此行为的一个示例。处理器的性能和带宽的相对平衡会改变这些曲线，从而改变的最优值问．
事实上，这一现象导致了整体的性能行为图3 (b)其中，我们比较了三种单套接字实现的单套接字执行时间:(i)使用SSE intrinsic显式向量化的单线程CPU;(ii)具有四线程和SSE的多线程CPU;(iii)使用CUDA的GPU。注意，线程化和向量化的组合比单线程代码提高了7.4×。GPU代码比单线程CPU代码快32倍，比多线程CPU代码快4.3倍，甚至包括从GPU发送/接收数据的时间。有趣的是，GPU的胜利不仅仅来自于并行化和调优，还因为更快的处理能够带来更好的效果算法调优问．特别是，我们执行计算密集型ULI步骤的速度越快，我们就能做得越大问从而减少VLI和其他内存限制的步骤。换句话说，在更快的处理和算法调优之间存在协同作用，我们可以用更少的内存通信(如VLI步)来交换更多的浮点数(ULI步)。随着多核处理器使用的增加，这类技术可能会得到更广泛的应用。
5.数值实验
本节评估我们的实现在两种不同体系结构上的整体可伸缩性，以及在强和弱伸缩机制下的可伸缩性。(我们用斯托克斯核代替拉普拉斯核。K（r)定义为1/|r| (1 -r
r) /r2),是两个三维向量的外积。)
平台和体系结构:我们在两个超级计算平台上测试了我们的代码，一个基于传统的CPU，另一个基于混合CPU+GPU节点。CPU系统为捷豹PF在国家计算科学中心(UT/ORNL)， (http://www.nccs.gov/computing-resources/jaguar/) Cray XT5, 224,256核(使用2.6 ghz六核AMD Opteron cpu, 2GB/核)，三维环面网络(9.6 GB/s链路带宽)。捷豹在500强中排名第三(www.top500.org)截至2011年6月。GPU系统是基恩兰初始输送系统，同样在UT/ORNL (http://keeneland.gatech.edu)．Keeneland基于惠普SL390服务器，使用NVIDIA特斯拉M2070 gpu进行加速。Keeneland有120个计算节点，每个节点都有双插座，六核Intel X5660 2.8 ghz Westmere处理器，每个节点有3个gpu，有24GB的DDR3主机内存。节点间采用单轨QDR ib组网互联。
MPI，在Jaguar上的强可伸缩性测试:图4 (a)总结了强尺度实验。这个问题有1亿个未知数。我们使用非均匀(线)分布。在24576个核心上，评估阶段只需要1.7秒，当增加512倍的核心时，这构成了205倍的加速。设置阶段(树形结构)开始在强尺度系统中占主导地位。
MPI, Jaguar上的弱可伸缩性测试:图4 (b)总结了弱尺度变换的结果。每个核心的问题大小被固定为大约672,000个未知数。由于几个因素，FMM评估阶段和设置阶段(树形结构)的次数增加了。首先，由于我们有高度不均匀的分布(对于均匀树，叶子级别的最大值是23，最小值是4，叶子级别将是9)，随着我们增加问题的规模，树变得更深，每个核心的成本增加(即，我们观察到aO（N日志N)缩放)，因为我们还没有达到渐近O（N)阶段。其次，对于非均匀树，很难实现FMM各阶段的负载平衡。解决这些缩放问题的方法是对FMM的所有阶段使用超立方类树广播算法。(目前只在评估阶段的海报沟通阶段使用。)最后，设置阶段不是多线程的，这是未来工作的主题。尽管如此，我们仍然获得了相当好的计算资源利用:评估阶段维持超过1.2 GFlopsper核(双精度)。
在Keeneland上进行GPU强可伸缩性测试:图5在Keeneland系统的多达192个gpu上，显示了4800万个粒子问题的强大缩放能力。每个GPU使用1个MPI进程，每个节点使用3个GPU(总共64个节点)。在最好的情况下，评估仅在0.47秒内完成。作为比较，我们估计等效的纯直接(O（n2))在192个gpu上计算，忽略所有通信，将需要大约1000倍以上的挂钟时间，这证实了在算法优化方法中使用并行的重要性。
6.讨论和结论
我们提出了几种算法，它们共同暴露和利用了快速多极算法的所有阶段的并发性，并采用了几种并行编程范式。我们展示了我们可以有效地扩展树设置，主要代价是并行排序，这反过来展示了教科书式的可伸缩性。我们描述了FMM算法的一种新的约简方案，并演示了整体的可扩展性。我们探索了在GPU加速器上使用流范式的逐核并发，并提供了出色的加速。FMM是一种具有多个不同阶段的高度非平凡算法，它结合了多分辨率数据结构、快速变换和高度不规则的数据访问。然而，我们能够在异构体系结构上实现显著的加速。
在我们在196608个核心上的最大运行中，我们在评估中观察到大约220 TFlops/s，大约15%的Linpack持续性能。假设在同样大小的混合机器上加速7倍，意味着我们的FMM在没有进一步修改的情况下将超过1千万亿次浮点运算。我们的主要结论是，构建高效的可扩展是可能的N-身体解算器，至少高达数十万核心。但如果是100万个或数百万个核呢?
对我们来说，一个重要的教训是混合并行是绝对必要的;一旦我们达到数千个MPI进程，算法通信密集阶段的MPI资源管理就成了一个问题。在套接字中使用共享内存并行并结合加速器，可以得到更可伸缩的代码。然而，可维护性开销可能是显著的，特别是对于GPU情况和像FMM这样的复杂代码。目前，我们正致力于在树结构中引入共享内存并行，在cpu和GPU之间共享工作负载(本质上将GPU视为另一个MPI进程)，并在向上和向下的计算中引入共享内存并行。
致谢
这项工作得到了美国国家科学基金会(NSF)拨款CNS-0929947, CCF-0833136, CAREER-0953100, OCI-0749285, OCI-0749334, OCI-1047980, CNS-0903447和半导体研究公司(SRC) 1981年奖的部分支持，美国能源部(DOE)拨款DEFC02-10ER26006/DE-SC0004915和美国国防高级研究计划局(DARPA)的拨款。本材料中表达的任何意见、发现和结论或建议都是作者的观点，并不一定反映NSF、SRC、DOE或DARPA的观点。TeraGrid系统上的计算资源是根据TeraGrid分配授权ASC-070050N、CCR-090024、ASC-100019和MCA-04N026提供的。我们要感谢TeraGrid的支持人员，也感谢NCSA、TACC和NICS的工作人员和顾问，我们从他们那里得到了重要的帮助。
参考文献
1.一种分层O(N logN)力计算算法。自然324,4(1986年12月)，446449。
2.kalahan, P.B, Kosaraju, S.R.动态最接近对和n体势场算法。在第六届离散算法ACM-SIAM年会论文集，SODA '95(费城，PA, USA, 1995)，工业与应用数学学会，263272。
3.Carrier, J.， Greengard, L.， Rokhlin, V.粒子模拟的快速自适应多极算法。暹罗j .科学。Stat。第一版。9， 4(1988)， 669686。
4.切里，比特森，r.k.，纽萨姆，G.N.径向基函数的快速求值:广义多二次函数的方法rn．暹罗j .科学。第一版。23， 5(2002)， 15491571。
5.Darve, E.， Cecka, C.， Takahashi, T.并行集群、多核处理器和图形处理器上的快速多极方法。Comptes Rendus Mécanique 339(23)(2011), 185193。
6.Grama, Gupta, A.， Karypis, G.， Kumar, V.《并行计算导论:算法的设计与分析》，第二版，Addison Wesley, 2003。
7.Grama, a.y.， Kumar, V.， Sameh, A. BarnesHut方法用于n体模拟的可扩展并行公式。在1994年超级计算会议记录，超级计算94(Los Alamitos, CA, USA, 1994)， IEEE计算机协会出版社，439448。
8.Gray, A.， Moore, A. N-Body的统计学习问题。放置神经通知。过程系统．(2001), 521527。
9.格林加德，L.，格罗普，W.D.快速多极法的并行版本。第一版。数学。: 20。， 7(1990)， 6371。
10.格林加德，L.，罗克林，V.粒子模拟的快速算法。j .第一版。phy。73(1987), 325348。
11.Gumerov, n.a.， Duraiswami, R.图形处理器的快速多极方法。j .第一版。phy。227， 18(2008)， 82908313。
12.Hamada, T.， Narumi, T.， Yokota, R.， Yasuoka, K.， Nitadori, K.， Taiji, M. 42 TFlops在gpu上的分层n体模拟，在天体物理学和湍流中都有应用。在SC09学报》， SCxy会议系列(波特兰，俄勒冈，2009年11月)，ACM/IEEE。
13.Hariharan, B, Aluru, S.压缩八叉树的高效并行算法和软件及其在分层方法中的应用。平行第一版。31(34)(2005), 311331。
14.JáJá， J.并行算法导论，Addison Wesley, 1992。
15.库扎克，J.， Pettitt, B.M.分布式存储机器的快速多极方法的大规模并行实现。65 .[参考译文， 7(2005)， 870881。
16.经典物理的快速算法。科学265， 5174(1994)， 909914。
17.拉舒克，我，等人。异构体系结构上的大规模并行自适应快速多极方法。在高性能计算网络、存储和分析会议论文集，SC '09(纽约，NY, USA, 2009)， ACM, 58:158:12。
18.Miller, G.， Teng, S.， Thurston, W.， Vavasis, S.用于球形填充和最近邻图的分隔器。J. acm， 1(1997)， 129。
19.绪方，S.等。并行计算机上快速多极方法的可扩展和可移植实现。第一版。理论物理。153年通信。， 3(2003)， 445461。
20.Phillips, J.C, Stone, J.E, Schulten, K.使消息驱动的并行应用程序适应gpu加速的集群。在SC '08: 2008 ACM/IEEE超级计算会议论文集(2008), 19。
21.Rahimian, A.， Lashuk, I.， Veerapaneni, S.， Chandramowlishwaran, A.， Malhotra, D.， Moon, L.， Sampath, R.， Shringarpure, A.， Vetter, J.， Vuduc, R.， Zorin, D.， Biros, G. Petascale在200k核和异构体系结构上血流的直接数值模拟。在SC '10: 2010年ACM/IEEE超级计算会议论文集(皮斯卡塔韦，新泽西州，美国，2010)，IEEE出版社，112。
22.经典势理论积分方程的快速解。j .第一版。理论物理60(1985), 187207。
23.Sevilgen, F.， Aluru, S.， Futamura, N.一个可证明的最优，分布无关的并行快速多极方法。在第14届国际并行与分布式处理学术研讨会论文集，2000。IPDPS 2000(2000), IEEE 7784。
24.孙达，H.， Sampath, r.s.， Biros, G.线性八叉树的自下而上结构和2:1平衡的并行化。暹罗j .科学。第一版。30， 5(2008)， 26752708。
25.可证明良好的并行自适应n体仿真的划分和负载平衡算法。暹罗j .科学。第一版。192(1998)。
26.泰西尔，R.等。全天空弱透镜模拟700亿个粒子。阿斯特朗。497年12,54。， 2(2009)， 335341。
27.一种并行哈希八叉树n体算法。在超级计算机学报》， SCxy会议系列(波特兰，俄勒冈，1993年11月)，ACM/IEEE。
28.应玲，Biros, G. Zorin, D.二维和三维的核无关自适应快速多极方法。j .第一版。phy。196， 2(2004)， 591626。
29.Ying, L.， Biros, G.， Zorin, D.， Langston, H.一种新的并行核无关快速多极算法。在SC03学报》， SCxy会议系列(菲尼克斯，亚利桑那州，2003年11月)，ACM/IEEE。
30.Yokota, R.， Bardhan, J.， Knepley, M.， Barba, L.， Hamada, T.在多达512个gpu和10亿个未知数上使用快速多极BEM的生物分子静电学。第一版。理论物理。182年通信。， 6(2011年6月)，12721283。
数据
图1。FMM列表。在这里，我们展示了一个树八分之一的二维U, V, W和X列表的例子B．我们可以验证这一点我(B)在被包围的范围内C（P（B))。
图2。幽灵分区的交流。进程0将绿色分区发送给进程1，将红色分区发送给进程2，将棕色分区发送给1和2。左下角的白色分区是“内部的”，用于处理0，不发送给任何人。该方法适用于叶和非叶两种分区。(a)较细的八叉树层。(b)较粗的八叉树层次。
图3。Intra-node优化。(a)最佳每八分单位颗粒:增大叶片尺寸增加ULI、WLI和XLI的相对成本，同时降低其他阶段的成本。(b)实施比较:实心横条表示每个平台上观察到的挂钟时间。在CPU+GPU的混合实现中，ULI运行在GPU上，而其他步骤运行在CPU上，从重叠中产生整体效益。
图4。(a)捷豹PF上的强伸缩。22M粒子的强伸缩结果。每个条都标有它的绝对挂钟时间(秒)和相对于48核外壳的加速。每个MPI进程有六个核心(和六个OpenMP线程)。八叉树最精细的层次是九，最粗糙的层次是三。当核心计数增加512x时，评估和设置阶段分别显示出205x和48x的加速。(b) Jaguar PF上的弱缩放。模拟的弱可扩展性可达196,608个核，每个核约有672,000个未知数。每个条都由挂钟时间(秒)和并行效率标记。我们在每个套接字中使用一个MPI进程，在每个MPI进程中使用所有六个OpenMP线程。最精细到最粗糙的八叉树级别从24到4。 Our implementation maintains parallel efficiency levels of 60% or more at scale.
图5。GPU强大的扩展。我们展示了在Keeneland系统多达192个gpu上对4800万个粒子进行FMM评估的强大缩放能力，仅在0.47秒内完成。实点线表示测得的挂钟时间，虚线表示完美加速下的理想时间。每个点都用它的并行效率(理想= 1)来标记。
表
表1。符号。这里我们总结了用于定义FMM树中每个八分单位的交互列表的主要表示法。(符号“\”表示标准集排除。)请注意,W()不一定是叶八分之一;相反,自W()只存在于l，()意味着l．我们说if和共享一个顶点、一条边或一个面。看到
©2012 acm 0001-0782/12/0500 $10.00
允许为个人或课堂使用部分或全部作品制作数字或硬拷贝，但不得为盈利或商业利益而复制或分发，且副本在首页上附有本通知和完整的引用。除ACM外，本作品的其他组件的版权必须受到尊重。允许有信用的文摘。以其他方式复制、重新发布、在服务器上发布或重新分发到列表，都需要事先获得特定的许可和/或费用。请求发布的权限permissions@acm.org传真(212)869-0481。
数字图书馆是由计算机协会出版的。版权所有©2012 ACM, Inc.
没有发现记录