二阶静电粒子网格算法在多核架构上的性能

×可在www.sciencedirect.com在线获取理论计算机科学电子笔记340（2018）67-84www.elsevier.com/locate/entcs二阶静电粒子网格算法在Modern上的性能多核架构多米尼克S.布朗，1史蒂文·A.作者：Stephen A.贾维斯华威大学计算机科学系英国考文垂摘要在本文中，我们提出了一种新的静电，二阶粒子在细胞（PIC）算法的轮廓，它利用“鬼粒子”位于真实粒子的位置，以代表电荷分布。我们实现我们的算法内EMPIRE-PIC，PIC代码在桑迪亚国家实验室开发。我们测试了我们的算法在各种多核架构，包括NVIDIA GPU，传统的CPU，和英特尔的骑士登陆的性能。我们的初步结果表明，与前几代众核硬件相比，这些架构上的PIC应用程序的二阶方法的可行性。具体来说，我们看到Tesla K20和Tesla P100 GPU设备之间的二阶方法的性能提高了一个数量级，尽管这两个设备之间的理论峰值性能只提高了4。虽然这些初步结果显示，一阶方法的运行时间大幅增加，但我们希望能够在未来显示出改进的缩放行为和提高的模拟精度。关键词：Particle-in-Cell; PIC;二阶算法;众核; P100; KNL; K20; GPU; Broadwell1引言等离子体在各种环境中的行为是科学界的极大兴趣，特别是在聚变能研究领域。例如，法国的ITER项目等企业旨在开发世界然而，进行这些实验可能既危险又昂贵。因此，现代高性能系统的一个常见用途是模拟各种条件下的等离子体粒子单元（PIC）方法通常用于进行此类模拟。1电子邮件：Dominic. warwick.ac.ukhttps://doi.org/10.1016/j.entcs.2018.09.0061571-0661/© 2018由Elsevier B. V.这是一篇基于CC BY-NC-ND许可证的开放获取文章（http://creativecommons.org/licenses/by-nc-nd/4.0/）。68D.A.S. Brown等人/理论计算机科学电子笔记340（2018）67现代众核体系结构的出现，超越了一个极端的并行水平，使方法，以前是不可行的，由于计算费用现在可行。PIC代码往往无法充分利用这种增加的性能潜力，由于其高使用内存带宽。使用高阶PIC方法可以通过显著提高模拟精度来解决这个问题，与一阶对应方法相比，可以增加计算强度。这种更大的开销伴随着模拟期间所需的内存吞吐量的增长在本文中，我们将展示二阶PIC算法的性能。我们的实现使用了二阶有限元素和粒子，这些元素和粒子由周围的幽灵粒子集合表示。这些幽灵粒子每个都具有围绕真实粒子位置的相关权重和权重集，因此表示电荷分布。我们还测试了我们的PIC实现对各种现代计算架构上的一阶算法。具体而言，这项工作的贡献如下：– 本文介绍了在Sandia国家实验室的EMPIRE-PIC程序中开发的一种新的二阶静电粒子池算法该方法实现了标准的PIC程序在二阶，并利用– 我们比较我们的二阶算法的性能，其一阶对应的各种计算架构。这些包括传统的英特尔CPU，NVIDIA Tesla GPU（开普勒和Pascal架构）和英特尔的Xeon Phi Knights Landing。在这个分析中，我们考虑了诸如总执行时间、各种关键PIC内核的执行时间以及应用程序实现的平均内存带宽等指标;– 我们强调EMPIRE-PIC中的关键性能问题，以促进未来的研究和性能调整。本文的其余部分结构如下：第2节概述了与PIC代码相关的相关先前工作，以及高阶PIC方法;第3节介绍了EMPIRE-PIC以及PIC算法的基本原理;第4节包含了对我们新算法的描述;第5节和第6节定义了用于检查我们算法性能的测试用例，并介绍了从这些结果中获得的结果;最后，在第7节中，我们总结了本文，并强调了未来研究的潜在途径。2相关工作有许多产生式代码使用PIC算法来建模带电粒子在电场作用下的行为。EPOCH是一种电磁PIC代码，它模拟交错的Yee网格来表示电场和磁场[1，18]。这允许代码使用Yee的有限差分然后应用蛙跳方法以更新粒子位置。此外，EPOCH还使D.A.S. Brown等人/理论计算机科学电子笔记340（2018）6769使用高阶插值方案。通过使用迷你应用程序优化EPOCH的尝试也已经完成，粒子排序和提高编译器执行的自动向量化量导致显着的加速[2]。由劳伦斯伯克利国家实验室（LBNL）和劳伦斯利弗莫尔国家实验室（LLNL）开发的WARP代码实现了显式三维PIC算法，以模拟高强度离子束的行为[6]。与EPOCH类似，WARP采用二阶蛙跳方法来推进粒子，但仅在粒子和网格点之间使用线性插值函数。代替FDTD方法，静电泊松求解器用于进行场求解。WARP有几个版本，其中包括能够自适应网格细化（AMR）的代码，以及二维和三维电磁实现。Wang等人提出了由普林斯顿大学维护的GTC-P，试图开发现代的高度并行的PIC代码[16]。 GTC-P在环形几何中求解5维Vlasov方程，以便在托卡马克装置内进行磁约束聚变实验的模拟。该应用程序采用多层并行，并已移植到各种现代计算架构，如NVIDIAGTC-P的性能已经在许多著名的HPC系统上进行了测试，如红杉，Piz Daint，泰坦和天河2号[14，15，17]。这项工作发现，该应用程序在各种输入大小的强问题和弱问题缩放方面表现良好，最大的问题由8000万个单独的网格点组成。Decyk提出了一种针对图形处理单元优化的流式PIC算法，并使用二维静电PIC代码演示了结果[4]。类似于EPOCH，为了改进流访问和数据局部性，基于位置对粒子进行排序，使得紧密定位的粒子所需的网格点更可能被保持在高速缓存或寄存器中。在粒子永远不会离开细胞的冷等离子体的理想情况下，与CPU版本相比，流式算法在GPU上显示出每个粒子，每个时间步30倍Kong等人详细介绍了另一种GPU PIC实现[11]。这项工作的重要联系是提出了一种方法，以减少线程之间的写冲突的影响，当模拟中的粒子将电荷沉积回网格节点时，这是PIC方法中最重要的性能瓶颈之一。这是通过将问题空间分裂成相等大小的簇，并且只允许非相邻簇内的粒子同时沉积电荷来实现的，本质上是一种“着色”方法。虽然这并没有完全消除对原子的需求，但它确实减少了所使用的原子操作的总量。其他作者也对高阶方法在PIC算法中的应用进行了研究。Pointon提出了一种二阶电荷守恒算法，用于电磁PIC模拟中积累电荷和电流，70D.A.S. Brown等人/理论计算机科学电子笔记340（2018）67（13）。 该算法使用二阶加权函数来将电荷沉积在与粒子相邻的所有单元上。在问题的边界上，电荷权重从二阶平滑地变为一阶，一旦粒子本身在边界上，则完全变为一阶。 通过使用测试问题，在正方形网格中心的导电球的情况下，证明了电荷在舍入误差的范围内是完全守恒的。然而，该算法需要进一步测试，以确定其在具有强电场和磁场的三维问题上的性能。Jacobs和Hesthaven展示了一种PIC算法，该算法能够处理使用非结构化网格的四阶问题[9]。为了在非结构化网格上求解麦克斯韦方程组，使用Galerkin方法，该方法将问题空间划分为非重叠的三角形所用的正确的顺序[10]。 粒子按顺序被转换到欧拉坐标系中把他们的能量重新分配到问题网格上。类似于上面讨论的算法，该方法也是电荷守恒的。通过对各种一维和二维情况的测试，作者表明，该算法能够模拟基本的等离子体现象，同时具有很大的几何可变性。Moon等人还提出了一个高阶的，电荷守恒算法的非结构网格[12]。该方法的独特之处在于，电荷守恒是从第一原理实现的，高斯定律的保存本文提出了一种新的非结构网格的二阶静电PIC方法.与以前的工作相比，我们的算法利用一组鬼粒子包围的粒子，以模仿电荷分布。我们的实现大量使用了Trilinos库[7]，并使用Kokkos作为其并行编程模型[5]。我们还在各种现代众核硬件上测试了我们的算法的性能，包括NVIDIA GPU，Intel CPU和Intel3背景EMPIRE-PIC是桑迪亚国家实验室（SNL）开发的非结构化PIC应用程序，能够在二维和三维空间中进行静电和电磁它使用了Trilinos库的许多软件包[7]，特别是它的线性求解器，并使用Kokkos [5]作为其并行编程模型。这使得代码可以在各种各样的计算架构上使用，而不需要更改原始程序源代码。这些架构包括传统的CPU和英特尔D.A.S. Brown等人/理论计算机科学电子笔记340（2018）6771NE（xi）=j=03.1粒子网格法在EMPIRE-PIC中实现的PIC方法用于通过在非结构化网格上表示电场（E）和磁场（B）来模拟等离子体的行为。这些场的值存储在每个网格节点上，插值用于确定单元内特定点的场值模拟还包含表示被模拟的等离子体的离散粒子。每一个都与网格上的位置、电荷、质量和速度相关联。标准PIC算法基本上由四个主要计算步骤组成。图1中的流程图对这些进行了总结。其中第一个是在当前时间步求解E和B的新值为了做到这一点，我们必须解决麦克斯韦然而，由于静电近似，我们只需要求解高斯#» ρ·E=（一）其中ρ是电荷密度，ε0是自由空间的介电常数，ε 1是发散算子。然后使用有限元方法（FEM），- 刚度矩阵S，其可以通过各种迭代或直接技术来求解，以便获得电势φ的值。然后可以通过以下等式确定电场，其中v表示所使用的有限元素类型的基函数：#»Ej=0φjvj（2）一旦获得了网格节点处的电场值，就有必要在更新它们的速度和位置之前通过插值来确定每个粒子位置处的E的精确值在结构化代码中，这通常通过使用简单的方法来实现，例如面积加权（或三维体积加权）。这尤其适用于使用交错网格的规范，如Yee [18]提出的规范。然而，在非结构化的情况下，我们必须应用不同的方法，因为网格单元可以采取任意形状。在EMPIRE-PIC中，我们将电场聚集到粒子位置，如等式3所示，其中，re是相对基函数，xi是粒子位置。#»广告下一步是更新模拟中粒子的速度，这可以通过确定当前场施加到粒子上的力来找到。这个力由洛伦兹力定律描述：#»q#»#»F=m（E+v×B）（4）N0Ejej（xi）（3）Σ72D.A.S. Brown等人/理论计算机科学电子笔记340（2018）67Fig. 1.总结PIC算法求解洛仑兹力的标准方法，从而更新粒子PIC代码中的值称为BorisM#ethod [3]（有时称为Boris Pusher）。在我们的静电情况下，B= 0，因此计算减少在给定的粒子的电荷质量比的乘积，以及在其当前位置的电场值。一旦确定了速度，我们就可以通过乘以时间步长来更新粒子的位置来计算粒子应该移动多少。在跨越元素边界的情况下，移动被分解为多个段，并且该过程依次应用于这些段中的每一个。核心静电PIC算法的最后一步是让粒子在下一场求解开始之前将其电荷贡献沉积回网格节点。这个过程几乎与将场加权到粒子的过程相同粒子必须对每个节点有贡献，并且相应地更新相关节点。然而，节点更新必须使用原子加法操作来执行，以防止多个线程试图同时更新同一节点的电荷值，这将导致错误的结果。4算法在本节中，我们将介绍我们的二阶算法的概要及其一些实现细节。该算法的伪代码在算法1中示出。在总结EMPIRE-PIC中实现的二阶算法之前，首先有必要引入一阶和二阶有限元素的概念。考虑由四边形单元组成的二维网格。如第3节所述，电场的值存储在该网格的每个元素的节点处。在第一阶情况下，节点位于单元格的四个角中的每一个，并且必须使用基函数来插值单元格内特定位置处的E值在二阶时，我们仍然必须进行插值以确定精确的网格值。使用二阶元素可以提供更多的插值点D.A.S. Brown等人/理论计算机科学电子笔记340（2018）6773算法1我们的二阶PIC算法1：对E进行初始场求解2：对于i←0到时间步长，3：对于每个中心粒子，4：从每个相关子粒子中累积场值5：结束6：对于每个中心粒子，7：使用Boris方法8：结束9：对于每个中心粒子，10：根据新速度更新位置11：结束12：对于每个中心粒子，13：来自每个相关联的亚粒子的14：结束15：对E进行16：结束值，从而得到更精确的解决方案。然而，这种提高的准确性带来了显着更高的计算成本，而所需的内存带宽增加较少。为清楚起见，图2显示了一阶和二阶四边形有限元的示例，其中突出显示的点表示电场的插值点。为了将EMPIRE-PIC的核心算法提高到这种更高的精度水平，传统PIC算法的大部分组件都被修改为使用这些二阶元素。唯一的例外是场解，标准FEM可以像以前一样应用于场解，尽管现在的计算更加复杂。我们的算法的一个独特的功能是使用“鬼粒子”，以代表一个平滑的电荷分布。与传统的PIC不同，每个粒子占据网格上的一个位置，每个物理粒子被一组幽灵粒子取代，中心粒子占据“真实”位置。这些中的每一个都与确定由幽灵粒子携带的总电荷的比例的权重和确定其相对于主要中心粒子的位置的权重集相关联。图3显示了我们算法中粒子的一个简单示例。4.1粒子的场加权当对由场解确定的E值加权到我们的新型粒子时，我们遵循与第3节定义的过程非常相似的过程，将电场乘以粒子位置处元素基函数的值然而，在我们的二阶算法中，我们必须累积电场对与被处理粒子相关的所有子粒子的贡献，然后才能确定主中心粒子位置处的场的正确值。为此，我们首先获得位置处的基函数值74D.A.S. Brown等人/理论计算机科学电子笔记340（2018）67(a) 一阶（b）二阶图二、一阶和二阶四边形有限单元示例图三.我们的二阶算法中的一个粒子的例子。 幽灵粒子以灰色显示，真实粒子位置以白色显示子粒子的（通过将它的o集加到中心粒子位置），然后乘以正常的电场值。然后根据与特定子粒子相关的权重对该结果进行加权。通过这种方式，我们计算每个子粒子感受到的电场值，然后通过简单地将每个子粒子对每个网格维度的贡献求和来获得真实粒子上的值是微不足道的。4.2粒子移动在我们的算法中，使用鲍里斯推来计算粒子的速度是不变的。 这是因为所有的子粒子都是在一个特定的位置，从真正的粒子位置。因此，我们不需要根据子粒子的速度来计算它的新位置然而，我们不能在粒子移动例程中简单地绕过子粒子，因为有必要更新每个子粒子所在的有限元素。这是通过获取中心粒子的速度和子粒子的位置来实现的，然后使用这些来确定子粒子的运动是否会导致它停留在同一个元素中，或者穿过边界进入相邻的元素。然后，子粒子的元素位置将相应更新。虽然这确实会导致一些冗余的数据计算，但执行的计算量很小，因此对应用程序性能的影响非常小中心粒子的位置使用先前计算的速度进行更新，如正常的那样。D.A.S. Brown等人/理论计算机科学电子笔记340（2018）67754.3向电网与传统的PIC算法非常相似，我们的二阶方法中的粒子必须在下一个场求解之前将其电荷贡献沉积回网格节点。这遵循的过程很像我们的场加权过程，只是相反。对于每一个粒子，它的每一个子粒子都单独地沉积它们的电荷。像往常一样，我们从它的位置集合中获得子粒子的位置，以确定在其位置处的元素基函数的值。 然后，我们将粒子电荷乘以其权重，基函数值以及粒子在模拟中表示的物理粒子的数量，然后将其放置在每个网格节点上。应该指出的是，我们的算法仍然支持线程之间的竞争条件的问题时，存放此费用。 目前，这个问题是通过使用原子操作来解决的，在许多其他代码中，尽管这由于减少的并行性而导致性能开销。5实验装置为了检查我们算法的性能，我们将其行为与一阶实现进行比较。为此，我们使用一个二维两个流不稳定性问题作为我们的主要测试用例。 它由两个光束组成由带正电荷的离子组成的静止流将电子向相反的方向运动。这导致等离子体波的能量增加我们还提出了朗道阻尼问题的性能结果，在那里我们模拟了带电等离子体环境中纵波能量指数下降的影响这可以被认为是所描述的两个流场景的逆。我们最大的测试用例在一个网格上使用了总共100万个粒子，该网格由x维的64个单元和y维的2个单元组成，总共500个时间步。这导致每个网格单元中有大量的粒子，便于以后研究原子操作开销对算法性能的影响除此之外，我们使用各种粒子计数的模拟来评估我们的PIC实现的可扩展性。 除了总执行时间和各种关键内核的时间，我们还检查了应用程序在第一和第二阶使用的平均内存带宽，以确定移动到第二阶是否会大幅改变代码的内存限制程度。我们在以下计算架构上测试EMPIRE-PIC：NVIDIA Kepler GPU、NVIDIAPascal GPU、Intel Broadwell CPU和Intel的Knights Landing（KNL）。本工作中使用的硬件信息见表1。为了这项工作的目的，使用的所有系统都由一个节点组成。对于在Broadwell系统上运行的实验，我们使用了28个线程，其中线程绑定到内核，并且禁用了超线程。在KNL上，我们发现使用绑定线程，每个核心两个超线程，可以获得最佳结果，因此我们在KNL实验中总共使用了128个线程76D.A.S. Brown等人/理论计算机科学电子笔记340（2018）67系统峰值TFLOPS内存带宽GB/s英特尔至强E5-2660 v4双插槽（Broadwell）0.448102英特尔至强融核7210（KNL）3.0475NVIDIA Tesla K201.17208NVIDIA Tesla P1004.7732表1用于测试EMPIRE-PIC图四、使用100万个粒子，500步的双流问题的执行时间6结果6.1双流问题为了获得EMPIRE-PIC在两个流问题上的性能的总体概述，我们首先比较表1中描述的所有架构上的代码的总执行时间。图4显示了在一阶和二阶精度水平下运行仿真所获得的时序。很明显，更现代的系统比传统的系统性能好得多。Tesla K20，同时考虑一阶和二阶问题。具体来说，我们看到K20和P100之间的二阶数量级差异，尽管P100的理论峰值性能仅比K20高4倍。这种差异的一个可能的解释是Pascal架构上原子操作的硬件实现，而不是Kepler基于锁的软件实现的原子。因此，当将粒子电荷加权回网格节点时，大量使用原子操作在Tesla K20上将更加明显虽然KNL在二阶时的表现大约是P100的一半，但它在一阶时的表现几乎一样好。尽管KNL和P100之间存在差异，但KNL的性能与Broadwell CPU相当，在一阶和二阶实验中只需要稍微长一点。图5显示了在一阶和二阶实验中执行代码的关键内核所花费的时间，其中每个内核代表PIC算法的四个阶段之一在所有使用的架构上，我们看到D.A.S. Brown等人/理论计算机科学电子笔记340（2018）6777图五、双流实验中关键EMPIRE-PIC核的执行时间内核负责将场加权到粒子，并将粒子电荷加权回网格。由于这些核中的每一个都需要在每个粒子的位置处计算基函数，因此这种差异并不令人惊讶，因为二阶基函数比一阶基函数的计算成本更高。由于幽灵粒子的存在而需要的电荷沉积量的增长也有助于这种增加，因为每个时间步执行的原子操作的数量更大。与其他三个系统相比，P100上的这种开销不太明显，这也是因为原子操作的硬件实现。执行加权操作的内核的高性能成本也是为了评估基函数而将物理点转换为参考点的费用。这是一个涉及矩阵计算的迭代过程，因此即使对于简单的几何形状也是昂贵的。关于场解，我们观察到求解E所需的时间几乎没有变化，所有系统之间的差异约为1秒。这是可以预料的，因为所使用的方法基本上是相同的，相反，我们只是使用了一个稍微复杂的元素，具有更多的节点。78D.A.S. Brown等人/理论计算机科学电子笔记340（2018）67见图6。EMPIRE-PIC的加权电荷核的执行时间，启用和不启用原子为了评估原子操作对电荷加权内核的开销，我们还进行了实验，其中禁用了对电荷阵列的原子访问限制。我们注意到，虽然由于随之而来的数据风险而产生了不正确的模拟结果，但它确实提供了对内核和应用程序性能的洞察。图6显示了在一阶和二阶启用和不启用原子的情况下执行权重电荷内核所花费的时间。在第一阶，我们看到特斯拉K20上最明显的差异，其性能差异大于一个数量级。然而，我们在其余系统上看到的差异要小得多，TeslaP100上的执行时间差异与机器噪音无法区分。在第二个订单的差异仍然可见，与K20继续显示一个数量级的差异在运行时。在布罗德韦尔和KNL上，性能差距也变得更加明显。P100上的差异仍然很小，突出了硬件实现的原子操作的影响。这些发现表明，静电PIC中的快速电荷沉积问题与现代NVIDIA加速器的相关性较小。然而，这个问题的解决方案对于PIC应用程序在所有现代架构上的性能都是有价值的我们还考虑了EMPIRE-PIC的执行时间如何随着模拟中粒子数量的增加而变化图7显示了EMPIRE-PIC在一阶和二阶算法下运行各种粒子计数所需的时间毫不奇怪，K20的扩展性比现代架构差得多，其最大的一阶运行时间比所有现代系统上最大的二阶运行时间都长。首先，我们看到P100的性能略差于Broadwell系统，这是由于将特定数据从主机传输到设备的开销。这种差异随着粒子数量的增加而减小，因为所需计算量的增加降低了数据传输成本对整体性能的影响在二阶时，P100的缩放比Broadwell和KNL好得多，尽管所有三个系统都是可比的。P100性能的改善，D.A.S. Brown等人/理论计算机科学电子笔记340（2018）6779一千100101数千个粒子(a) 一阶数千个粒子(b) 二阶K20 P100 KNL Broadwell图第七章双流实验的定标结果系统/问题内存带宽GB/sK20一阶17.50K20二阶17.32P100第一订单29.12P100二阶45.88KNL First Order20.55KNL二阶47.04布罗德韦尔一阶5.93布罗德韦尔二阶4.73表2表中显示了所有系统上一阶Broadwell是由从DRAM移动的每字节的FLOP的增加引起的，这是由于二阶方法的更高的浮动点强度。从图7b中可以看出，P100在相对较低的颗粒计数下超过了Broadwell系统，这与一阶的类似行为形成了对比英特尔系统显示出非常相似的可扩展性，随着问题规模的增加，两个系统之间的初始差距逐渐缩小。因此，当非加速剂系统相互竞争时，在我们的算法上，P100显然是优越的，这是合理的结论，这将是对更大的问题大小的情况。然而，在第一阶，所有三种现代架构都表现出类似的缩放，三者之间的差异很小。缩放行为确实表明，如果一阶问题的大小继续增加，P100和英特尔系统之间将出现性能差距虽然在所有问题规模下，与一阶运行相比，执行我们的二阶算法所需的时间之间存在显着的差异，但很明显，二阶方法在现代众核硬件上比这些系统的早期版本更可行。表2显示了我们的EMPIRE-PIC中实现的平均内存带宽执行时间31.2562.5125250500100031.2562.5125250500100080D.A.S. Brown等人/理论计算机科学电子笔记340（2018）67图八、使用100万个粒子，500步的Landau阻尼问题的执行时间实验我们看到Tesla K20实现的带宽几乎没有差异，而Tesla P100和KNL上的二阶问题明显增加。但是，由于两个系统的带宽限制都很高，这些增加并不会使应用程序的内存限制显著增加。令人惊讶的是，我们观察到约1 GB/s的减少，在Broadwell系统上实现的带宽在二阶。这一结果值得进一步调查。6.2Landau阻尼问题图8-10示出了在朗道阻尼实验上执行用于双流问题的相同测试程序的结果朗道阻尼问题表现出几乎相同的性能的两个流的问题，除了粒子移动核。在所有的架构上，这个内核的执行时间都有显著的增加。这可以通过问题中表示的物理空间的大小小于两个流输入的大小来解释。因此，由于单元格的数量保持不变，因此减少了一个单元格所表示的物理空间量。这导致粒子在单个步骤内可能跨越更多边界，从而导致每个粒子平均执行更多次移动例程（因为当跨越边界时，移动被分成一系列较小的移动缩放行为在两种类型的问题中也是一致的，在一阶和二阶精度下，Lan-dau阻尼实验在所有架构中的缩放仅略差这种缩放能力的下降是由粒子移动的执行时间增加引起的。增长的效果在具有高粒子计数的二阶运行中最明显，其中粒子移动执行得更频繁。表3 显示了我们的朗道阻尼实验中EMPIRE-PIC实现的平均存储带宽。对于Tesla K20、Broadwell和KNL，我们看到两个流问题的行为几乎相同，而Tesla P100在一阶和二阶运行中显示出类似的结果。D.A.S. Brown等人/理论计算机科学电子笔记340（2018）6781图第九章Landau阻尼实验中关键EMPIRE-PIC核的执行时间系统/问题内存带宽GB/sK20一阶16.29K20二阶17.11P100第一订单43.31P100二阶44.89KNL First Order22.64KNL二阶42.16布罗德韦尔一阶8.19布罗德韦尔二阶6.56表3表中显示了朗道阻尼问题7结论本文介绍了一种新的静电二阶粒子网格算法及其在EMPIRE-PIC中的实现。EMPIRE-PIC是Sandia国家实验室开发的PIC应用程序。我们已经研究了该算法的性能和缩放时相比，其一阶对应的变化，82D.A.S. Brown等人/理论计算机科学电子笔记340（2018）67一千100101数千个粒子(a) 一阶数千个粒子(b) 二阶K20 P100 KNL Broadwell图10. Landau阻尼实验的标度结果包括传统CPU、NVIDIA GPU和英特尔的Knights Landing在内的多核系统我们的结果显示，Kepler和Pascal GPU架构在二阶问题上的性能存在数量级差异然而，其他现代众核系统仍然相互竞争，并与P100相当例如，Intel的Broadwell和KNL在两个流和Lan- dau阻尼实验中在大颗粒计数下显示出几乎相同的我们还表明，由于现代系统上可用的高带宽，虽然二阶方法目前的性能和规模比一阶更差，但我们的研究结果表明，与前几代众核硬件相比，二阶PIC方法在现代架构上变得更加可行。我们还强调，对于我们的非结构化代码，重要的性能问题位于涉及将电场值加权到粒子位置的例程中，反之亦然。这是由于原子操作的使用，以及为了确定基函数值而将物理点转换为参考点。7.1今后工作这是第一篇论文来记录我们的静电二阶算法及其在各种计算架构上的性能行为。在未来的工作中，我们的目标是使用收敛性分析，量化的精度提高，通过使用我们的方法。此外，通过解决性能问题，我们希望展示改进的性能和扩展行为，而不是一阶问题。最后，我们将把我们的算法扩展到电磁模拟，其中粒子受到磁场和电场的影响。执行时间31.2562.5125250500100031.2562.51252505001000D.A.S. Brown等人/理论计算机科学电子笔记340（2018）6783确认这项工作得到了联合王国原子武器研究所（AWE）的支持，资助号为CDK0724 （AWE 技术推广方案）。Stephen Jarvis教授是AWE William PenneyFellow。作者还要感谢马修·贝当古博士感谢桑迪亚国家实验室对该项目的指导和工作，以及提供对EMPIRE-PIC源代码的访问。引用[1] Arber，T.D、K. 贝内特角S. Brady，A.Lawrence-Douglas，M.G. Ramsay，N.J. Sircombe，P.吉利斯，R. G. Evans，H. Schmitz，A. R. Bell和C. P. Ridgers，Contemporary Particle-in-Cell Approach toLaser-Plasma Modelling，Plasma Physics and Controlled Fusion57（2015），p. 十一万三千零一。[2] 伯德河F.、P. Gillies，M. R. Bareford，J. A. Herdman和S. A. Jarvis，性能优化使用微型应用程序的惯性确认融合代码，国际高性能计算应用杂志（2016年）。[3] 鲍里斯，J.，相对论等离子体模拟：混合代码的优化，在：第四届等离子体数值模拟会议论文集，海军研究实验室，华盛顿特区，1971年，第100页。3-68[4] Decyk，V.K.和T. V. Singh，Adaptable Particle-in-Cell Algorithms for Graphical Processing Units，Computer Physics Communications182（2011），pp. 641-648[5] 爱德华兹，H. C.的方法，C. R. Trott和D. Sunderland，Kokkos：通过多态内存访问模式实现众核性能可移植性，并行和分布式计算杂志74（2014），pp。3202-3216，用于高性能计算的特定领域语言和高级框架。[6] Grote ， D.P. ， A. Friedman ， J. 维 和 我 Haber ， The WARP Code ： Modeling High Intensity IonBeams，AIP Conference Proceedings749（2005），pp.55-58.[7] Heroux，M.一、R. A. 巴特利特河谷E. 豪勒河J. Hoekstra，J.J. Hu，T.G. 科尔达河，巴西-地B. 勒胡克，K. R.隆河，西-地P. Pawlowski，E. T. Phipps，A. G. Salinger，H. K.索恩奎斯特河S. Tuminaro，J. M.Willenbring，A. Williams和K. S. Stanley，An Overview of the Trilinos Project，ACM Transactionson Mathematical Software31（2005），pp. 397-423.[8] Hesthaven，J.和T.Warburton，非结构网格上的节点高阶方法：I。Maxwell方程的时域解，计算物理杂志181（2002），pp。186-221。[9] Jacobs，G.和J. Hethhaven，High-order Nodal Discontinuous Galerkin Particle-in-Cell Method onUnstructured Grids，Journal of Computational Physics214（2006），pp.96比121[10] Jacobs，G.，D. Kopriva和F.马沙耶克，多域交错网格谱方法的粒子跟踪算法，载于：第39届航空航天科学会议和展览，航空航天科学会议，美国航空航天研究所，2001年。[11] 孔 ， X. M. C. 黄 角 Ren 和 V.K. Decyk ， Particle-in-Cell Simulations with Charge-conserving CurrentDeposition on Graphic Processing Units ， Journal of Computational Physics230 （ 2011 ） ，pp. 1676-1685年。[12] Moon，H.，F. L. Teixeira和Y. A. Omelchenko，非结构化网格上粒子模拟的精确电荷守恒散射聚集算法：几何视角，计算机物理通信194（2015），pp。43比53[13] Pointon ， T. ， Second-order ， Exact Charge Conservation for Electromagnetic Particle-in-CellSimulation in Complex Geometry，Computer Physics Communications179（2008），pp.535-544[14] 唐，W.，B. Wang，S.Alzheimer，G.克瓦希涅夫斯基，T.Hoe Escherer，K.Z. 易卜拉欣，K。Madduri，S.威廉姆斯，L.奥利克角Rosales-Fernandez和T.极端尺度等离子体湍流模拟全球顶级超级计算机，在：高性能计算，网络，存储和分析国际会议论文集，SC '16（2016），pp。43：1-43：12。84D.A.S. Brown等人/理论计算机科学电子笔记340（2018）67[15] 王，B.，S. W.唐氏T.威廉斯，K. Z.易卜拉欣，K。Madduri，S. Williams和L.奥利克，Leadership-class系统极端尺度下的动力学湍流模拟，在：高性能计算，网络，存储和分析国际会议论文集，SC'13（2013），pp。82：1-82：12。[16] 王，B.，S. W. M.唐，K. Z.易卜拉欣，K。Madduri，S. Williams和L. Oliker，顶级超级计算机上聚变等离子体的现代回旋粒子模拟，国际高性能计算应用杂志（2017）。[17] Wang，E.，S.吴角，澳-地张建，刘杰，W. Zhang，Z. Lin，Y. Lu，Y. Du和X. Zhu，天河二号超级计算机上聚变等离子体的回旋运动粒子模拟，在：第七届研讨会大规模系统可扩展算法（ScalA）的，2016年，pp.25比32[18] Yee，K.美国，各向同性介质中涉及麦克斯韦方程组的初始边值问题的数值解302-307