COFFEE-用于微分方程数值演化的MPI并行Python包（2019）

软件X 10（2019）100283COFFEE-用于微分方程数值演化的MPI并行Python包Georgios Doulisa，Jörg Frauendienerb，Chris Stevensc， Ben Whaleda波兰华沙华沙大学物理系理论物理研究所b新西兰达尼丁9054奥塔哥大学数学与统计系cRhodes University，Grahamstown 6140，South Africad澳大利亚伍伦贡大学数学与应用统计学院，伍伦贡，新南威尔士州2522ar t i cl e i nf o文章历史记录：收到2019年2019年7月3日收到修订版接受日期2019年7月3日保留字：Python微分方程解算器a b st ra ctCOFFEE（Conformal Field Equation Evolver）是一个Python软件包，主要用于使用线方法随时间数值地演化偏微分方程系统它包括各种时间积分器和有限差分格式与求和的部分财产，以及作为自旋加权函数的角导数的伪谱函数。一些额外的功能包括在各种不同的几何图形上进行MPI并行化，HDF数据输出和后处理脚本以可视化数据，以及允许用户在每个时间步后创建代码进行分析的操作类。©2019由Elsevier B.V.发布这是一个在CC BY-NC-ND许可下的开放获取文章（http://creativecommons.org/licenses/by-nc-nd/4.0/）中找到。代码元数据当前代码版本v1用于此代码版本的代码/存储库的永久链接https://github.com/ElsevierSoftwareX/SOFTX_2019_93法律代码许可证GNU通用公共许可证（GPL）使用git的代码版本控制系统软件代码语言、工具和服务使用Python 2.7、C、MPI、HDF5、numpy、spinsfast、libfftw3编译要求，操作环境兼容性：mpi4py与兼容的MPI实现，h5py与兼容的HDF实现，numpy，scipy，Gnuplot，matplotlib和各种标准Python模块（例如ABC、日志和数学）。如果可用，代码和存储库中提供的开发人员文档/手册文档的链接可通过gitlab存储库与1. 动机和意义我们提出了一个软件包，共形场方程进化或COFFEE的短，实现技术适用于数值解的时间依赖的微分方程系统（DE）通过线的方法。COFFEE主要用Python实现。它对用户的要求很少，并且是以PEP 8 [1]作为指导思想编写的。虽然COFFEE在速度上无法与一些现有的数值积分器竞争，但它提供了较低的使用门槛和很大的灵活性。*通讯作者。电子邮件地址：c.stevens@ru.ac.za（C. Stevens）。https://doi.org/10.1016/j.softx.2019.100283COFFEE是专门为计算代表弗里德里希共形场方程的双曲型偏微分方程（PDE）系统的解而开发的它已被用于八个研究项目，以数值研究广义相对论的共形性，[3作为说明在[10]中，COFFEE的能力被用于解决初始边界值问题（IBVP）形式的偏微分方程系统，该问题包含46个变量和45个约束，涉及两个不同的高性能集群，使用多达200个进程。它使系统及时发展为一系列分辨率，以多种方式近似空间导数，稳定地施加用户给定的边界条件，并将数据存储在HDF文件中。后处理脚本证明了计算解决方案的收敛性和稳定性，并生成了2352-7110/©2019由Elsevier B. V.发布。这是CC BY-NC-ND许可证下的开放获取文章（http://creativecommons.org/licenses/by-nc-nd/4.0/）。可在ScienceDirect上获得目录列表SoftwareX期刊主页：www.elsevier.com/locate/softx2G. Doulis，J.弗劳恩迪纳角Stevens et al.粤ICP备10028830号-1输出. COFFEE包含严格研究时变偏微分方程系统数值演化所需的工具COFFEE的独特之处在于，没有一个计算PDE软件是以用户友好和灵活的理念设计的，用于解决复杂和具有挑战性的方程系统，例如从爱因斯坦场方程导出的方程这将COFFEE与Cactus[11]，Chombo [12]或PETSc [13]等现有软件进行了2. 软件描述COFFEE旨在显著减少编写代码求解方程组所需的工作量。因此，尽管COFFEE具有数值性质，但它是用Python实现的，并且严重依赖于numpy、mpi4py、hdf5和自定义C代码。用Python实现有明显的缺点。例如，由于Python是一种解释型语言，语法可能对执行速度有很大的影响（比较循环和列表压缩），解释型语言存在结构问题（在Python的情况下，这是GIL的原因），将Python脚本转换为Python字节码然后转换为机器指令的额外开销，以及缺乏强类型语言中的编译时检查。尽管如此，论文[3用Python实现也有优势。特别值得注意的是鸭子类型，动态内省和代码注入，这减少了用户遵守严格的编程模式和理解“COFFEE做事方式”的需要符合Python的哲学，COFFEE评估给定的代码并快速失败，即。停止模拟而不是继续潜在的缺陷。在模拟的每次迭代每个动作都是一段任意的用户代码，可以完全访问当前时间步的所有数据和执行模拟的几乎所有对象。这使用户能够完全控制模拟。COFFEE代码库的编写是为了可读性而不是性能（除了处理MPI和HDF API的代码）。代码中包含大量注释，突出显示了代码中比较复杂的部分，为什么选择了某些算法，以及更改后容易出现错误的代码部分。一般来说，我们认为代码不应该被视为不可变的工作体，而应该像在白板上工作一样：根据需要添加和修改。因此，我们希望用户直接改变COFFEE运行时控制可以在第2节中描述的动作或系统对象中进行。然而，由于使用MPI和HDF APICOFFEE的核心功能是支持MPI的实现，在计算机集群上支持空间导数的谱和有限差分技术的因此，它是一个实现的数值方法所需的发展依赖于时间的系统的常微分方程和偏微分方程，例如抛物线和双曲型系统。特别值得注意的是，COFFEE包括用于施加稳定边界条件的同时近似（SAT）方法[14]的代码，用于部分求和有限差分算子的论文[15-[18] 对 于 快 速 自 旋 加 权 球 谐 变 换 的 数 值 实 现 的 微 积 分 （ eth-calculus），参见例如[19]。IO使用HDF5进行数据存储。它已经运行在桌面工作站，在奥塔哥大学的计算机集群和新西兰eScience基础设施对于本节的其余部分，读者参考COFFEE的UML序列图可能会有所帮助，如图11所示。1.一、IBVP类表示初始边值问题。若要初始化此类对象，则必须提供行为类似于求解器描述了如何沿着模拟线执行步骤。系统计算被模拟函数的时间导数网格对象描述函数的域并管理MPI API。初始化时具有合理默认值的其他选项包括模拟期间要执行的操作列表（例如，数据简化、可视化、误差计算等）、最大迭代次数和最小时间步长。求解器是知道如何集成一维常微分方程的对象。提供了一个抽象基类以及Euler（显式和隐式）、4阶Runge Kutta方法和4阶Runge Kutta方法的变体（包含中间步骤的边界数据）的实现。自适应四阶龙格库塔法已经实现，但没有充分测试，为这个版本的COFFEE。这方面的代码最终将包含在存储库中。系统对象表示要使用直线法求解的微分方程组。抽象基类是作为自定义类提供的API的文档的一种形式提供的，否则将隐式地假定API是由自定义类提供的，填充系统对象的角色。系统类必须有一个返回下一个时间步的方法（这允许自适应模拟），一个计算特定时间点的时间导数的方法（evaluate（）），以及给出初始值和边界值的方法。为了使用户的生活更轻松，已经实现了许多数值（空间）导数运算。这些方法中的10%已经在已发表的论文中用于evaluate（）提供了以下运算符的代码11有限差分格式，使用快速傅立叶变换的9种不同实现，The Geroch–Held–Penrose operators7个分部求和有限差分算子和3个支持耗散算子。The Geroch–Held–Penrose operators are supported by a mod-ule,Grid对象表示正在计算的函数的域。要实例化网格对象，必须提 供 每 个 维 度 中 的 数 据 点 数 量 、 数 据 点 坐 标 值 的 边 界 、MPIInterface对象和描述“边界数据”的对象。MPIInterface对象包装了MPI_COMM本身的一个实例，该实例由mpi4py模块包装（ COFFEE 已 经 针 对 MPICH 和 Open MPI 进 行 了 测 试 ） 。COFFEEMPIInterface对象理解如何将足够的数据传递给相邻进程，以便在通过MPI在多个进程上执行模拟时继续进行模拟。边 界 数 据 的 抽 象 基 类 ABCBouncement 表 示 确 定 通 过MPIInterface的子类在不同进程上的子网格之间发送什么数据所需的信息。 ABCBouncement对象在子栅格“边”之间进行区分，子栅格“边”是由Grid对象表示的完整栅格的内部和外部。原则上，每个内部和外部边缘可以具有不同数量的鬼点和将被传送的“边界”区域上的点。 这些数据通过boundary_slices（）方法传递。一旦模拟开始，COFFEE就进入ibvp类中包含的首先，确定下一个时间步长····G. Doulis，J.弗劳恩迪纳角Stevens et al.粤ICP备10028830号-13Fig. 1. COFFEE的UML序列图。函数调用和变量的名称已在图中保留，t是当前时间，tslice是函数数据的时间片（在advance（）和evaluate（）方法中传递，并在initial_data（）调用中创建）或函数数据的派生（在diffop（）调用）或函数在中间时间的值（从evaluate（）调用返回），i是当前迭代的编号，tslice.data是函数在当前时间的值。请注意，严格来说，tslice对象的生存期是不正确的，因为在每个内部片段循环。求解器负责在每个片段计算之间创建新的时间片。我们把它从图中去掉了，因为它只会使它复杂化。其次，执行动作列表中的每个动作。第三，计算函数在当前时间加上时间步长的值。第四，重复该过程，直到出现异常，例如：发生溢出或到达最后时间主循环每次迭代中的所有数据都存储在时间片对象中。每个时间片包含函数数据、网格和当前时间。咖啡最后一个值得讨论的主要部分是行动。动作是一段代码，它用时间片做一些事情动作是用户可定义的，可以包含对运行时可用的当一个动作被调用时，它会传递当前的时间片。动作具有完全的自由度，因此能够动态地影响模拟。一个动作应 该 继 承 actions.Prototype 类 ， 并 且 至 少 实 现 函 数 _doit（）。COFFEE由附加脚本支持，这些脚本处理生成的hdf文件，以生成二级派生信息的正常数组，例如计算收敛率、错误、可视化和hdf文件中的数据操作在使用COFFEE计算解决方案之前，只需要满足几个要求，参见第3节的示例或参考COFFEE存储库中提供的更详细的示例之一。用户必须提供：1. 具有返回下一个时间步长的方法的对象一种计算系统在给定时间点的时间导数的方法（我们称之为“系统”），2. 选择一个可以求解常微分方程（ODE）或提供自己的（“求解器”）的对象3. 指定要计算的解的域以及要使用的离散化一旦从COFFEE中的库中选择了这些组件，或者自定义编写了这些组件，它们就会在对象初始化期间作为参数传递给ibvp类。现在可以通过调用ibvp.run（）方法来运行模拟。还有一系列额外的选项可以指定，包括IO、MPI拓扑和数据通信、计算空间差异的方法、特定时间的模拟和强制评估的长度以及实时可视化。要运行模拟，用户必须实例化对象，这些对象是COFFEE提供的System 、 Grid 、 Solver 和IBVP 类 。 COFFEE 已 经 实 现 了Grid、Solver和IBVP类，这些类对于大多数模拟来说已经足够了。由于System类表示要求解的DE的系统，因此将其留给用户来实现。这些类共同包含执行模拟所需的所有信息。在上面引用的研究项目中，这些类的初始化和模拟的开始都被收集在我们所谓的设置文件中。清单2是一个最小的安装文件示例。COFFEE存储库包含更详细的示例。为了在研究中使用，奥塔哥数值相对论小组在安装文件中包含了设置文件中的典型选项包括：日志记录、微分运算符、求解器、输出设置、2.1. 软件构架为了与Python的哲学保持一致因此，显式的代码依赖性要么是显而易见的，例如，IO系统依赖于h5py，或者明确声明，例如，Grid对象使用与Grid类一起位于grid.grid模块中的ABCBouncement类的实例。为了记录Python依赖于duck类型而产生的隐式依赖关系，提供了许多抽象基类。这些类记录了关于类API的其他隐式解释。我们鼓励用户继承抽象基类，尽管这不是必需的。4G. Doulis，J.弗劳恩迪纳角Stevens et al.粤ICP备10028830号-1∗×= − =+表1强标度测试，使用[ 8 ]中描述的对称双曲PDE的一维系统，具有：z 1和z 1之间的12801个等距空间网格点，空间步长为0.5，CFL为0.5，时间步长为CFL为0.5，空间导数近似为具有部分求和性质的四阶有限差分算子[17]，边界条件用SAT方法[14]实现。模拟在新西兰eScience Infrastructure的Mahuika集群上运行结进程数电话：+86-21 - 6666666传真：+86-21 - 6666666图二. 由清单1和清单2计算的函数的图形。组件0是求解清单1中指定的波动方程的函数的值。1，而组分1是其衍生物。选择时间是为了显示在两个边界的第一 次 反 射 之 后 ， 但 在 单 独 的 波 包 第 一 次 合 并 之 前 的 函 数 。 您 可 以 通 过 使 用gpl_partial操作仅在time 3.0运行来复制此图。由于隐式依赖关系最有可能给新用户带来问题，因此我们简要介绍一下它们。图1描述了在数据初始化和一次模拟迭代期间的预期数据流，因此提供了COFFEE内部隐式依赖关系的3. 说明性实例我们给出了一个代码示例，该代码解决了代码列表1和2中的一维波动方程（也参见图2）。这些清单中的代码可以在COFFEE存储库中找到。清单1中给出的第一个文件定义了系统对象。它指定要使用的空间微分运算符、如何计算时间步长、初始数据以及如何计算时间导数的解决方案。清单2中给出的第二个文件将模拟所需的对象打包，并将它们交给管理主模拟循环的IBVP4. 影响目前只有少数其他软件包包含与COFFEE类似的功能集。每一种都有不同的哲学、结构和目标。例如，Cactus代码[11，20]是一个非常大的，社区驱动的项目，开发了多年。然而，Cactus和COFFEE的哲学有很大的不同;Cactus是非常优化的，主要用Fortran和C编写，并且对如何完成以及用户可以做什么非常严格。另一方面，COFFEE是用Python编写的，旨在方便使用和灵活性。这些特性使咖啡和仙人掌自然互补，并在自己的权利有用特别是，COFFEE，虽然强大到足以满足高级程序员，也是针对用户从各种不同的领域，可能有有限的编程技能。这些用户想要一种简单的方法来数值演化常微分方程或偏微分方程的系统，而不必担心速度。为了帮助减少运行时间，COFFEE是MPI并行的：假定域被分割成更小的域，每个域具有其自己的存储器和专用核心。使用MPI的性能提高通过表1中给出的强缩放测试结果来举例说明。可以看出，增加MPI进程的数量将运行时间减少到大约32个核心，之后增加MPI进程又开始增加计算时间。这是预期的，因为总网格大小是固定的，并且增加MPI进程的数量停止减少运行时间将随着总网格大小的增加而增加。对于足够大数量的MPI进程，每个子网格之间的通信已经增加到通信本身成为瓶颈的程度。在COFFEE中，用户与MPI的交互是最小的，在设置文件中只需要几行Python代码，详细说明：网格的尺寸，周期性（如果有的话）和拓扑结构（例如笛卡尔）。在系统文件中的进程之间可能还需要多几行通信数据，例如用于空间导数近似。COFFEE的存储库中给出了如何执行上述操作的示例此外，COFFEE包含使用[18]的优化变换算法的自旋加权球谐函数的实现目前涉及COFFEE的研究涉及[18]的修改版本，该版本针对轴对称进行了优化。此代码未在此版本的COFFEE中发布，但将在完成足够的测试这种特殊的自旋加权球谐码在[10，21]中使用COFFEE包是非常通用的，因为它对用户输入的系统几乎没有因此，如果一个过程可以通过微分方程或微分方程组来模拟，那么它就可以在COFFEE中进行数值演化。当然，是否能找到数值解将取决于方程组的具体情况和所选择的数值方法。5. 结论COFFEE是一个用户友好的Python包，用于ODE和PDE（系统）它包含了各种各样的数值算法，用于在时间上进行，近似空间导数和稳定地施加边界条件，以及通过分裂计算域进行MPI并行化。它已经在多个研究项目中进行了严格的测试，并通过actions类在演化过程中执行用户定义的任务COFFEE非常适合那些将用户友好性置于绝对速度之上并希望灵活地将代码用于特定问题的用户。信用作者身份贡献声明Georgios Doulis：软件。JörgFrauendiener：软件。 克里斯史蒂文斯：软件。 Ben Whale：软件。确认作者希望感谢使用新西兰电子科学基础设施（NeSI）高性能计算设 施 和 咨 询 支 持 作 为 本 研 究 的 一 部 分 。 新 西 兰 网 址https://www.nesi.org.nz。政务司司长获G. Doulis，J.弗劳恩迪纳角Stevens et al.粤ICP备10028830号-15奥塔哥大学博士奖学金。这项工作得到了新西兰皇家学会Marsdengrant UOO0922的支持。竞合利益作者声明，他们与其作者身份或本文的发表没有利益冲突引用[1]https://www.python.org/dev/peps/pep-0008/网站。[2]弗里德里希·H爱因斯坦方程和共形结构：反德西特型时空的存在性。J GeomPhys1995;17（2）：125-84.[3]Beyer F，Doulis G，Frauendiener J，Whale B.类空和零无穷大附近的数值时空。闵可夫斯基空间中的自旋2系统。经 典 量子引力2012;29（24）：245013.[4]Doulis G，Frauendiener J.中的二阶自旋-2系统平坦类空空间和零无限空间。Gen Relativity Gravitation2013;45（7）：1365-85.6G. Doulis，J.弗劳恩迪纳角Stevens et al.粤ICP备10028830号-1[5]Beyer F，Doulis G，Frauendiener J，Whale B.类空和零无穷远附近的线性引力波。在：数学相对论，引力和宇宙学的进展。Springer; 2014，p.3比17[6]Beyer F，Doulis G，Frauendiener J，Whale B.闵可夫斯基背景下的自旋2方程。在：数学相对论，引力和宇宙学的进展。Springer; 2014，p.465-8[7]Beyer F，Daszuta B，Frauendiener J，Whale B.用基于自旋加权球谐函数的谱方法计算2-球上的场的数值演化。经典量子引力2014;31（7）：075019.[8]放大图片Frauendiener J，Stevens C，Whale B弗里德里希-纳吉规范中平面引力波的数值演化。Phys Rev D 2014;89（10）：104026。[9]Doulis G，Frauendiener J.包含类空无穷大的闵可夫斯基时空的全局模拟。物理修订版D2017;95（2）：024035。[10]杨志华，李志华，李志华.广义共形场方程的数值初边值问题。物理修订版D2017;96（8）：084020。[11]Goodale T，Allen G，Lanfermann G，Massó J，Radke T，Seidel E，Shalf J.仙人掌框架和工具包：设计和应用。在：计算科学的高性能计算国际会议。Springer; 2002，p. 197-227[12]Colella P，Graves D，Ligocki T，Martin D，ModianoD，SerafiniD，Van Straalen B. Chombo AMR软件包应用程序设计文档。2009年，可在Chombo网站：https：//commons.lbl.gov/display/chombo/Chombo+-+Software+for+Adaptive+解+偏+微分+方程。[13][10]张晓刚，张晓刚. PETSC用户手册。技术报告ANL-95/11-修订版2.1。5，阿贡国家实验室; 2004年。[14]Carpenter MH，Gottlieb D，Abarbanel S.求解双曲型方程的时间稳定边界条件系 统：甲 羟吲 哚ogy 及 其在 高 阶 紧致 格 式中 的 应用 。 J ComputPhys1994;111（2）：220-36.[15]Carpenter MH，Nordström J，GottliebD. 一稳定和 保守的，保守的任意空间精度的有效界面处理。J Comput Phys1999;148（2）：341-65.[16]Diener P，Dorband EN，Schnetter E，Tiglio M.满足分部求和的优化高阶导数和耗散算子及其在三维多块演化中的应用。J Sci Comput2007;32（1）：109-45.[17]斯特兰德湾d/dx有限差分近似的分部求和。J Comput Phys1994;110（1）：47-67.[18]Huffenberger KM ， Wandelt BD. 快 速 精 确 的 自 旋 球 谐 变 换 。 Astrophys JSuppl Ser2010;189（2）：255.[19]彭罗斯R，林德勒W。旋量与时空：第1卷。剑桥：剑桥大学出版社，1984年.[20]http://www.cactuscode.org/网站。[21]Beyer F，Escobar L，Frauendiener J.具有空间拓扑S3和对称群U的宇宙时空爱因斯坦方程的数值解（1）. Phys Rev D 2016;93（4）：043009。