pyBaram：高性能Python并行可压缩流求解器

≥软件X 20（2022）101272原始软件出版物pyBaram：用于教学和研究的高性能Python并行可压缩流求解器朴镇锡韩国仁川市延寿区盖别路36号仁荷大学航空宇宙工学系ar t i cl e i nf o文章历史记录：收到2022年8月2日收到修订版2022年11月7日接受2022年关键词：可压缩流有限体积法非结构网格Pythona b st ra ctpyBaram是一个用Python编写的并行可压缩流求解器。它基于非结构网格上的有限体积法，可以很容易地处理二维和三维复杂的几何形状。稳健的数值格式和湍流模型被应用于解决高速和高雷诺数流动。计算可以通过高性能Python包加速和并行化。易于阅读的代码和文档有助于理解和修改它们。建议的求解器作为一个教学工具来描述现代计算流体动力学（CFD）的可压缩流。此外，它的目的是提供一个研究框架，以解决工业和学术问题的空气动力学与并行计算。本文介绍了结构和功能的求解器与说明性的例子。©2022作者由爱思唯尔公司出版这是CC BY许可下的开放获取文章（http://creativecommons.org/licenses/by/4.0/）中找到。代码元数据当前代码版本0.2.4用于此代码版本的代码/存储库的永久链接https://github.com/ElsevierSoftwareX/SOFTX-D-22-00220可复制胶囊的永久链接法律代码许可证BSD 3条款使用git的代码版本控制系统软件代码使用的语言、工具和服务Python编译要求、运行环境依赖性Python≥ 3.7结节≥1.10Numba≥ 0.5h5py≥ 2.6mpi4py≥2.0如果可用，链接到开发人员文档/手册tqdm 4.0inha.gitlab.io/pyBaram/技术支持邮箱jinseok. inha.ac.kr1. 动机和意义计算流体动力学（CFD）已被广泛应用于许多工程应用，如航空航天，汽车，生物医学和推进工业。由于计算技术的快速发展，模拟真实复杂应用周围的流动成为可能。特别是，现代计算流体动力学可以解决整个飞机周围的湍流可压缩流，它可以提供可靠的气动力 和 力 矩 在 稳 态飞 行 条 件 。 在 各 种 数 值 方 法 中 ， 有 限 体 积 法（FVM）是解决可压缩流动问题的主流方法.此外，具有非结构化电子邮件地址：jinseok. inha.ac.kr。https://doi.org/10.1016/j.softx.2022.101272网格可以容易地处理复杂的几何形状。通常，在空间中获得二阶精度，并且鲁棒的激波捕获方案使得能够模拟超音速和高超音速流。采用雷诺平均Navier-Stokes（RANS）湍流模型模拟高雷诺数湍流流动。隐式时间积分方案用于有效地获得稳态解。因此，一些传统的求解器[1-这些求解器是非常成熟的，他们也被扩展到解决多物理问题。然而，这些求解器变得非常巨大的代码库，使得可能不容易理解它们的整体结构和实现。学生和实践者可以使用这些求解器来模拟空气动力学问题。2352-7110/©2022作者。由爱思唯尔公司出版这是CC BY许可下的开放获取文章（http://creativecommons.org/licenses/by/4.0/）。可在ScienceDirect上获得目录列表SoftwareX期刊主页：www.elsevier.com/locate/softxJ. S. 公园软件X 20（2022）1012722∇∇然而，他们感到挑战，以实现新的想法，并扩展这些解决方案，为他们的学术研究。这项研究提出了pyBaram，一个用Python编写的并行可压缩流求解器。它是基于传统的混合非结构网格有限体积法。采用成熟的格式和物理模型求解可压缩湍流。所有代码都是用Python脚本语言编写的。计算可以通过高性能Python包（如Numba）加速和并行化[8]和mpi4py [9]。因此，pyBaram同时支持进程和线程并行。计算速度接近于C/C++和Fortran等编译语言编写的代码该求解器旨在作为描述可压缩流的现代CFD理论的教学工具，也作为解决空气动力学中的工业和学术问题的有许多项目具有类似的目标，以传播CFD技术。其中有些是用脚本语言描述新发展的数值格式。然而，他们专注于用一个简单的模型问题来演示方法[10]。其中一些是用Fortran或C/C++编写的，并支持并行计算。然而，他们专注于低速流[11，12]。本研究的重点是简单的Python脚本语言的可压缩流的CFD方法的直接执行。求解器不使用可能隐藏基本和重要逻辑的专用库。因此，仅掌握Python语言知识的学生或从业者可以很容易地理解这些实现。该求解器还可以用于研究活动，进行大规模的并行模拟现实的空气动力学问题。本文的结构如下。第二节简要介绍pyBaram求解器及其架构和功能的描述。第3节给出了ONERA M6机翼的跨音速绕流和HB-2标准模型的超音速和高超音速绕流的示例。分析了pyBaram的并行性能，并与传统的CFD求解器进行了比较。此外，本求解器的影响进行了讨论，并得出结论，在其余两个部分。2. 软件描述2.1. 数值格式可压缩流的控制系统可以写成如下形式：简体中文t+其中，U对应于状态变量矢量;Fc和Fv分别表示无粘通量和粘性通量。S是源向量。pyBaram在非结构网格上应用了格心有限体积法，可以处理二维和三维空间中的各种混合单元Eq.的半离散化形式（1）对于每个单元格，在等式中呈现（二）、在F面处。它们可以通过MUSCL型重建来获得（三）、Uf+=U<$i+φi<$Ui·xi，f，（ 3）其中，Ui对应于单元i处的状态变量的梯度，并且xi，f表示从单元中心到脸上此外，φi是单元i处的斜率限制器，用于鲁棒地捕捉激波不连续性;Uf−可以在相邻单元处类似地计算计算右手边等式（2）可以概括如下：• 计算每个像元• 计算每个像元• 使用MUSCL类型重建计算Uf+和Uf−• 计算数值通量：Hc和Hv• 计算通量和源项的负散度显式（2）及时。表1总结了pyBaram中实现的数值方案。在pyBaram中应用的标准FVM理论的细节可以在文献[13- 15 ]中找到。2.2. 软件构架PyBaram中的所有代码都是用Python编写的计算密集型循环可以使用Numba包[8]使用共享内存模型加速和并行化。此外，消息传递并行性也可以用mpi4py包实现[9]。因此，pyBaram可以在由多核处理器组成的集群机器上进行大规模并行模拟。所有代码都非常容易阅读，因为遵守Python编码标准。此外，使用面向对象的编程范式，因此，新功能可以很容易地实现继承或重用现有的类。pyBaram由几个模块组成，如图所示。1.一、这些模块负责预处理、模拟流程和后处理。现将其简要说明如下：“reader”和“partitions”模块处理网格软件生成的非结构化网格文件的导入和分区。The ‘‘partitions’’ module uses the METISlibrary [网格文件可以划分为MPI并行仿真。“积分器”模块控制流动模拟的迭代。集成器对象是从这个模块中定义的类实例化的。这个对象调用“求解器”模块中的系统类（二）、获得该项，并在积分器对象中更新解此外，插件被这个集成器对象调用以进行后处理。“求解器”模块计算数值通量及其散度，这是有限体积法的核心部分。首先，系统对象对应的治理系统，U<$i=−1∑（Hc（U+，U−，nf）−H（U<$f，<$Uf，nf））<$Af+S<$i，（2）tems在此模块中实例化系统对象还∂t∆ Vif调用元素、接口和顶点对象。它们能-erate计算内核，如2.1节所述。其中，U<$i和S<$i分别表示单元i处的单元平均状态变量向量和源项向量。 Hc和Hv分别表示数值无粘通量和粘性通量。乌夫 和Uf分别对应于面f 此外，nf和Af表示单位分别是面F的法向量和面积。是单元格i的体积。Uf+和Uf−是左和右状态向量The ‘‘backends" and ‘‘utils" modules accelerate the purePythonThe可以生成以下插件对象：收集残余历史，积分特定边界表面上的力或属性，并生成输出数据。····ff·vJ. S. 公园软件X 20（2022）1012723表1pyBaram中实现的数值方案的总结。Term描述控制系统Euler，几何尺寸2D、3D元素三角形，四边形，四面体，棱柱，金字塔，六面体梯度空间离散限制器MLP-u1/u2 [17，18]无粘通量RoeM [19]，AUSMPW+[20]，AUSM+up [21]，HLLEM [22]，Rusanov [23]时间积分5阶湍流模型SA 1-eq模型[27]，kω-SST 2-eq模型[28]Fig. 1. pyBaram的软件架构图“writers "模块将输出数据转换为可视化文件。源代码和文档在存储库中提供;因此，感兴趣的用户可以跟踪详细的实现。此外，还给出了示例，用户指南描述了运行它的指令。2.3. 软件功能pyBaram提供了一个控制台脚本。它可以执行五种模式：导入，分区，运行，重新启动和导出。导入和分割模式转换和分割由网格化软件生成的非结构化网格文件。虽然pyBaram没有网格化功能，但标准CGNS格式[30]或开源Gmsh [31]网格文件可以轻松转换。运行和重新启动模式通过求解Euler系统、Navier-Stokes系统或RANS系统进行流动模拟。需要配置文件以及转换后的网格文件。该文件以INI文件格式编写，因此方案和参数可以轻松配置。用户手册中描述了详细的选项和参数。导出模式将模拟结果转换为可视化软件，如Tecplot [32]或Paraview [33]。3. 说明性实例3.1. ONERA M6机翼ONERA M6机翼上的跨音速流是可压缩空气动力学的标准基准。在跨音速马赫数下，形成了一个大激波。强大而准确的激波捕捉方法和适当的湍流模型是必要的，以捕捉这种流动物理。良好的网格组织和详细的数值结果的传统求解器表2ONERA M6机翼升阻系数的比较升力系数阻力系数USM3D（最细网格）0.269 0.0170FUN3D（最细网格）0.269 0.0170pyBaram（细网）0.266 0.0170由NASA湍流建模资源提供[34]。这项研究使用了这个网页上第二精细的混合网格，由650万个六面体和220万个棱柱组成。自由飞行马赫数和攻角设定为0.84，3.06度，分别。基于根弦的雷诺数为1460万。空间离散采用二阶精度的MUSCL型MLP-μ 2限制器和RoeM数值无粘通量格式。时间积分采用LU-SGS格式，湍流模型采用kω图图2示出了上表面和下表面处的压力系数的等值线。在上表面上的半圆形激波结构清晰可见。比较了四个展向站的压力系数剖面分布（见图1）。（3）第三章。计算结果与实验数据吻合较好。表2将基于机翼平面形状面积的升力和阻力系数与USM3D和FUN3D解算器进行了比较。这些求解器使用了最精细的网格和不同的湍流模型; pyBaram的空气动力系数与这些结果吻合得很好。3.2. HB-2标准模型HB-2，超高速弹道型模型2，是数值模拟和实验的超音速和高超音速流的标准基准之一[35]。它包括一个·J. S. 公园软件X 20（2022）1012724××图二、O N E R A M6机翼上（左）下（右）表面压力系数的等值线。图三. ONERA M6机翼根部不同翼展位置截面压力系数的比较。x是从前缘沿局部弦线，其长度为c。钝头圆柱体外展强冲击波和膨胀波在这个弹丸周围发展以前的几个实验和数值数据是通过测量大范围超音速和高超音速流的气动力和力矩系数得到的。本研究考虑了马赫数为2和8时的流动基于物体直径的自由雷诺数为1。7 10 6和2. 1106，分别对应马赫数2和8。采用了由60万个四面体、40万个六面体、0.05万个棱锥体和0.03万个棱柱体组成的混合非结构化网格。空间离散采用二阶精度的MUSCL型MLP-μ 2限制器和AUSMPW+数值无粘通量格式时间积分采用LU-SGS格式，湍流模型采用kω图图4显示了对称平面上的马赫数等值线。机头附近的强弓形激波和喇叭口附近的斜向弱激波都得到了很好的分辨，两个马赫数 图图5比较了法向力系数随迎角的变化。pyBaram模拟与实验数据取得了良好的一致性[35]。3.3. 性能pyBaram的强大的可扩展性与线程和进程并行性进行了评估的计算节点与一个AMD EPYC 7513处理器和256 GB的内存。同时多线程（SMT）被禁用，因此最多可使用32个进程或32个线程. 表3总结了pyBaram的Python编译器和库。对HB-2标准模型绕流进行了定常数值模拟测量了超过5，000次迭代的运行时间;但是，排除了第一次迭代的时间，因为它包含JIT编译时间。图图6显示了通过纯MPI并行化和使用OpenMP或英特尔线程构建块（TBB）层的纯多线程实现的加速。MPI并行化实现了近J. S. 公园软件X 20（2022）1012725图四、H B -2在对称平面上马赫数等值线的比较。图五、H B - 2 垂直力系数随迎角变化的比较。表3可伸缩性测试中使用的编译器和库PythonNumbampi4py英特尔TBB版本3.9.100.55.13.0.32021.6最多可扩展至8个进程。当使用所有过程时，观察到24倍的加速但是，这两个多线程层的可伸缩性都略低于进程并行化。目前，Numba没有提供用户参数来控制线程关联或控制手动同步的方法。在其他应用程序中也观察到了Numba类似的多线程性能[36]。如果Numba的多线程支持足够成熟，pyBaram的线程并行化有望得到增强。此外，整个求解器的性能与SU 2求解器（一种由C/C++编译的开源CFD求解器）一致[7]。计算了同样的问题，并比较了运行时间。在数值方法方面，有一个明显的差异。SU 2是基于单元顶点有限体积法（FVM），pyBaram是基于单元中心有限体积法。在相同的非结构化网格上，以单元为中心的FVM的自由度大于单元顶点FVM的自由度;然而，两种方法的总体计算成本和性能相似[37]。除了这一点，两个求解器都使用了类似的数值方案。特别地，对于大型隐式运算矩阵，采用同样的松弛方法进行求解。 此外，一个类似的二阶重建方法与斜率限制器和相同的湍流模型应用于两个求解器。这两个解算器的配置文件作为电气补充材料提供。本基准测试使用了来自存储库的SU2 MPI 7.4.0for Linux二进制文件两个解算器利用纯MPI并行化使用整个32个进程。见图6。对比pyBaram使用线程和进程并行的强大可扩展性。表4中比较了两种求解器每次迭代的平均挂钟时间。两个求解器的计算成本之间的差异考虑到数值方法的差异，可以得出结论，pyBaram的性能与用编译语言编写的传统CFD求解器相同4. 影响pyBaram的主要目的之一是帮助学生和从业者学习航空航天领域的当前CFD方法J. S. 公园软件X 20（2022）1012726表4模拟HB-2标准模型周围超音速流的每次迭代平均挂钟时间的比较。pyBaram SU2时间0.194秒0.180秒工程.它采用有限体积法和高精度的数值格式来求解可压缩湍流流动。这些技术是模拟飞行器和航天器绕流的常用方法。pyBaram是用Python语言编写的，只有8000行代码，非常轻量级。由于Python的简单性和代码的文档化，新手可以很容易地理解整个代码结构。教科书和参考文献中的数值格式的直接实施也有助于他们理解CFD理论。此外，pyBaram还针对学术研究人员进行流动模拟或开发新算法。由于利用了高性能的Python包，pyBaram中的内核可以提供与C/C++或Fortran编写的代码近似等效的性能。此外，采用了消息传递和共享内存的并行性。因此，pyBaram甚至可以对完整的航空和航天飞行器周围的流动进行大规模由于pyBaram的模块化结构， 研究人员可 以轻松地实现 和修改算法并进 行模拟。pyBaram依赖于许多开源软件;它也受到以前开发良好的开源代码的启发。预计pyBaram还可以向其他研究人员传播现代CFD理论及其实施。5. 结论pyBaram是一个用Python编写的可压缩流求解器。实现了在航空航天工程中广泛应用的现代CFD方法pyBaram采用有限体积法;它可以处理二维和三维非结构化网格。鲁棒性和准确的空间离散方案和高效的解决方案的实施。所有函数和例程都是通过纯Python实现的。此外，高性能Python包用于加速和并行化它们。说明性的例子证实，pyBaram可以进行平行模拟的高速和高雷诺数的三维机翼和弹丸周围的流动，即使有强大的冲击波和膨胀波。预计pyBaram可以向学生提供新的功能，如实现新的算法，扩展系统，并利用GPU并行编程，计划作为未来的工作。我们欢迎任何贡献和合作。竞合利益作者声明，他们没有已知的竞争性财务利益或个人关系，可能会影响本文报告的工作数据可用性数据将根据要求提供。确认作者要感谢韩国仁荷大学通过研究基金提供的支持附录A. 补充数据与本文相关的补充材料可以在https://doi.org/10.1016/j.softx.2022.101272上找到。引用[1]弗林克NT。非结构四面体网格上求解欧拉方程的迎风格式。AIAA J 1992;30：70-7. http://dx.doi.org/10.2514/www.example.com[2]Manual F ， Biedron RT ， Derlaga JM ， Gnoffo PA ， Hammond DP ，Jacobson KE，等. FUN3D手册：13.7.2020年。[3]Baurle RA ， White JA ， Drozda TG ， Norri AT. VULCAN-CFD 用 户 手 册 。NASA;2020年。[4]MatssonJE.使用ANSYS Fluent 2021。SDC出版物; 2021.[5]西门子PLM。STAR CCM+用户指南。版本; 2019.[6]格林希尔兹角OpenFOAM V10用户指南。London，UK：The OpenFOAMFoundation; 2022，URLhttps://doc.cfd.direct/openfoam/user-guide-v10.[7]Economon TD ， Palacios F ， Copeland SR ， Lukaczyk TW ， Alonso JJ.SU2 ： 用 于 多 物 理 场 仿 真 和 设 计 的 开 源 套 件 。 AIAA J 2016;54 ： 828-46.http://dx.doi.org/10.2514/1.J053813上下载。[8]放大图片作者：Lam SK，Pitrou A，Seibert S. Numba：基于LLVM的PythonJIT编译器。 在：LLVM-HPC 2015会议录：关于HPC中LLVM编译器基础设施的第二次研讨会-与SC 2015联合举办：高性能计算、网络、存储和分析国际会议，Vol. 2015-Janua。2015年，http://dx.doi.org/10.1145/2833157。2833162[9]Dalcin L，Fang YLL. Mpi4py：经过12年的发展后的状态更新。计算机科学与工 程 2021;23 ： 47-54. http://dx.doi.org/10.1109/MCSE 的 网 站 上 下 载 。2021.3083216。[10] Ching EJ，Bornhoft B，Lasemi A，Ihme M. Quail：一个轻量级的Python开源不连续Galerkin代码，用于教学和原型设计.SoftwareX2022;17：100982. http://dx.doi.org/10.1016/j.softx.2022.100982上下载。[11]Blais B，Barbeau L，Bibeau V，Gauvin S，Geitani TE，Golshan S，et al.Lethe ： An open-source parallel high-order adaptative CFD solver for in-compressableflows.SoftwareX2020;12：100579.http://dx.doi.org/10.1016/j.softx.2020.100579天啊[12]Alinovi E，Guerrero J. FLUBIO-An非结构化，并行，基于有限体积的Navier-Stokes和对流扩散类方程求解器，用于教学和研究目的。SoftwareX 2021;13：100655. http://dx.doi的网站上下载。org/10.1016/j.softx.2020.100655。[13]Blazek J.计算流体动力学：原理与应用。Butterwest-Heinemann;2015.[14]Toro EF.流体动力学的黎曼解算器和数值方法：实用介绍。Springer ScienceBusiness Media; 2013.[15]莱维克双曲型问题的有限体积法，第31卷。北京：清华大学出版社.[16]Shima E，Kitamura K，Haga T.任意多面体非结构网格的格林-高斯/加权最小二乘混合梯度重建。AIAA J 2013;51：2740-7. http://dx.doi.org/10.2514/1.J052095上下载。[17]作者：Park JS，Yoon SH，Kim C非结构网格上双曲守恒律的多维极限过程。J Comput Phys 2010;229：788-812. http://dx.doi.org/10.1016/j.jcp.2009.10.011上下载。[18]作者：Kim C.非结构网格上有限体积法的多维极限处理。Comput Fluids2012;65：8-24. 网址：//dx.doi.org/10.1016/j.compfluid.2012.04.015上下载。[19]张文，等.激波不稳定性的计 算方法：一种激波稳定的Roe格式的发展. JComput Phys 2003;185 ： 342-74. http ： //dx.doi.org/10.1016/S0021-9991（02）00037-2.[20]Kim KH，Kim C，Rho OH.高超音速流的精确计算方法。I. AUSMPW+方案。J Comput Phys 2001;174：38-80. http://dx.doi.org/10.1006/jcph.2001.6873上下载。[21]AUSM 的 续 集 ， 第 二 部 分 ： AUSM+- 适 合 所 有 速 度 。 J Comput Phys2006;214：137-70. http://dx.doi.org/10.1016/j.jcp.2005.09.020上下载。[22] Einfeldt B，Munz CD，Roe PL，Sjögreen B.关于低密度附近的Godunov型 方 法 J Comput Phys 1991;92 ： 273-95. http://dx.doi.org/10.1016/0021-9991（91）90211-3.[23]Rusanov VV.非定常激波与障碍物相互作用的计算。NRC，机械工程部;1962年。[24]放大图片作者：Martinelli L.雷诺平均方程多重网格法的验证。AIAA第26届航空航天科学会议。一九八八年[25]Gottlieb S，Shu C-W.总变差递减Runge-Kutta格式数学比较1998;67：73-85.http://dx.doi.org/10.1090/s0025-5718-98-00913-2.[26]杨S，詹姆森A. Euler和Navier-Stokes方程的上下对称高斯-赛德尔方法。AIAA J 1988;26：1025-6. http://dx.doi的网站上下载。org/10.2514/3.10007。[27]Spalart P，Allmaras S.空气动力流的单方程紊流模型。Rech Aerosp1994;5-21.J. S. 公园软件X 20（2022）1012727[28]Menter FR.工程应用的双方程涡粘性湍流模型。AIAA J 1994;32：1598-605.http://dx.doi.org/10.2514/3的网站上下载。12149[29]Karypis G ， Kumar V. A fast and high quality multilevel scheme forpartitioning unregular graphs.SIAM J Sci Comput 1998;20：359网址：//dx.doi.org/10.1137/S1064827595287997上下载。[30]拉 姆 齐 角 计 算 流 体 动 力 学 通 用 符 号 系 统 CGNS 用 户 指 南 . 2022 ，URLhttps://cgns.github.io。[31] Geuzaine C，Remacle JF. Gmsh：一个内置预处理和后处理设施的三维有限 元 网 格 生 成 器 。 Internat J Numer Methods Engrg 2009;79 ：1309http://dx.doi.org/10.1002/nme.2579上下载。[32] Tecplot Tecplot 360 EX用户2022，URLhttp://www.tecplot.com。[33]张文辉，张文辉，张文辉. Paraview：一个用于大数据可视化的终端用户工具。在：可视化手册，卷。717，2005年。[34]李伟杰，李伟杰.湍流模型验证和确认网站资源的描述。第40届流体动力学会议和展览。2010年，p.四七四二[35]Gray D，lindsay yon Karman EE.标准超高速弹道模型HB-1和HB-2在马赫数1.5至10时的力试验。一九六三年[36]张文辉，张文辉.比较研究用于并行元分析的高生产率高性能编程语言。Swarm Evol Comput2020;57：100720.[37]放大图片Wang G，Schwöppe A，Heinrich R.湍流粘性流非结构τ-代码中单元中心和单元顶点离散化的 比 较 和 评 价 。在：ECCOMAS CFD，vol.2010年2010年