我们提出了一个并行实现局部尺度空气质量模型的计算框架，该模型由对流-扩散-反应偏微分方程描述，即所谓的反应扩散方程。 使用正向欧拉格式对模型进行时间离散。 使用有限元方法实现了空间离散化。 用于并行实现的策略基于使用消息传递库MPI的分布式内存方法。 模拟的重点是两种与道路交通相关的空气污染物，即颗粒物PM2.5和PM10。 并行实现的效率和可扩展性通过使用集群计算系统的多达128个处理器内核进行的数值实验进行了说明。
1.简介 空气污染是当今对人类健康和环境日益严重的全球性问题，具有重大的经济后果。 由于人口密度高和污染源高度集中，在城市环境中尤为严重[ 1 ]. 空气污染物是指空气中的任何物质，如果浓度足够高，可能会危害人类健康和动物，并破坏环境。 它们可能由固体颗粒、液滴、气体或这些形式的组合组成。 空气污染物分为两大类：一级污染物和二级污染物。 主要污染物直接从可识别的来源排放到大气中，这些来源可以是自然的，例如野火或火山爆发，也可以是人为的，例如工业活动或交通排放。 它们包括一氧化碳（CO）、二氧化硫（SO2）、氮氧化物（NOx）、颗粒物（PM）、氨（NH3）、挥发性有机化合物（VOCs）和有毒金属。 至于二次污染物，它们不是直接排放的，而是大气中主要污染物之间化学相互作用的结果。 二次污染物包括臭氧（O3）、二氧化氮（NO2）、三氧化硫（SO3）、硫酸（H2SO4）、硝酸（HNO3）、过氧化氢（H2O2）和颗粒物（PM）。
空气污染与一系列急性和慢性健康影响有关[ 2 ]例如中风、心脏病、慢性阻塞性肺病、肺癌、急性呼吸道感染、皮肤刺激、眼睛和耳朵炎症以及鼻咽喉（ENT）疾病。 世界卫生组织（WHO）报告称，全世界每年有420万人因接触环境空气污染而过早死亡[ 三 ]. 国际癌症研究机构（IARC）是世卫组织的专门癌症机构，已将室外空气污染归类为对人类致癌[ 4 ]，第1组。 面对健康和环境问题以及空气污染的经济代价，正在逐步实施在地方、区域和全球各级控制和调节空气质量的政策。 在较先进的国家，大量资源用于使用不断发展的技术，包括空气质量传感器，在目标地区的离散站测量和监测环境空气污染物浓度。 由于这些监测器往往很昂贵，有时很难使用，尤其是对于发展中国家而言，因此很难与它们在目标地区建立足够的网络，以便更好地评估空气质量指数，并作出减轻其影响的决定。
空气污染物运动和扩散的建模和模拟是可提供良好结果的替代方案[ 5 ]. 这些是传统上空气质量研究和预测所依据的常用工具。 空气污染建模是最重要和最具挑战性的科学问题之一[ 6 ]通常用于支持空气质量评估和管理中的决策。 它涵盖了污染物在大气中的传输和扩散、干沉降和湿沉降以及化学反应。 它还取决于污染物特性、气象条件、排放数据和地形参数。
本工作的重点是对西非内陆国家马里的政治和经济首都巴马科的道路交通相关空气污染进行建模和模拟（见图 1 ). 该市占国家领土的0.02%，集中了马里12.46%的人口，每公里人口密度为9062人 2 [ 7 ]. 马里70%以上的经济活动都在这里进行，这使得它成为马里主要的工业和商业十字路口。 经济活动的高度集中，伴随着人口的快速增长和日益城市化，是巴马科空气污染的主要驱动力。 此外，巴马科在地理上建在一个被群山环绕的盆地中，这有利于污染物的禁锢，使城市更加脆弱。 尽管发动机技术有了显著改进，但交通流量仍是城市地区的主要污染源[ 1 ]尤其是在巴马科，那里的道路基础设施很差，道路运输车队也过时了。 大量的柴油车和劣质的燃料造成了交通拥堵和大量污染物的排放。
上述因素意味着，与大多数非洲国家首都一样，巴马科是一个污染严重的城市，PM10是一种直径为10的可吸入颗粒物，其年浓度为 μ m或更小，可达到333 μ 克 / 米 三 [ 8 ]日峰值超过600 μ 克 / 米 三 而空气质量指南[ 9 ]世界卫生组织建议每日最大平均浓度为45 μ 克 / 米 三 另一种被广泛用作城市空气污染水平指标的有害污染物是PM2.5，这是一种直径为2.5的细颗粒物 μ m或更小。 根据安装在巴马科美国大使馆大院的监管级空气监测器收集的每小时PM2.5数据，得出2020年10月PM2.5的日平均浓度，如图所示 2 从图中可以看出，这些浓度的峰值达到165 μ 克 / 米 三 ，远远超过15 μ 克 / 米 三 [ 10 ]，WHO空气质量指南以及美国环境保护署（EPA）每日标准水平建议的最大值固定为35 μ 克 / 米 三 .
在本文中，我们提出了巴马科局部尺度确定性和基于欧拉的空气质量模型的并行实现框架和可扩展性结果。 我们将重点关注用于生成以下数值结果的计算模型的并行实现方面和可伸缩性行为[ 12 ]. 在节中 2 描述了控制大气污染物浓度时空演变的数学模型。 章节 三 详细介绍了模型的时间和空间离散化。 章节 4 重点介绍了模型的并行实现和有限元分析。 在节中 5 ，我们给出了该框架的数值实验和计算可扩展性结果。主要结论总结在第节 6 .
2.型号说明 让 是的有界计算域 和 的边界 . 我们考虑的是 形式的 , 哪里 是地面， 表示域的上限边界 和 表示横向边界。 让 是浓度场的矢量，其中每个元素 对应于标记的化学物质（污染物）浓度的标量场 在空中。 浓度的时空演变 在域中 在时间间隔内 由以下对流扩散反应描述[ 13 ]型号： 其中（i） 是标记的污染物浓度 , （ii） 是扩散系数，（iii） 是风速（或风速），（iv） 是描述化学反应的化学源项，（v） 是清除系数，（vi） 是源项，（vii） 是道路交通的地表排放量，（viii） 是干沉积速度，（ix） 是到边界的单位向外法向量 , （x） 是空间坐标，其中 和 是平面的变量 是海拔高度，和（xi） 是物理时间。 表中报告了这些物理场和参数的国际单位制（SI） 1 扩散系数参数化、干沉积、清除过程和化学动力学在[ 12 ]. 型号( 1 )基于化学物质和流场（包括风速和扩散系数）之间没有反馈的假设。 我们假设流体是不可压缩的。 这意味着假定空气密度的变化是恒定的。 我们还假设城市地形是均匀的。 气象参数受到周围大气的强烈影响。 污染物的输送主要发生在大气边界层（ABL）[ 14 ]，对流层的最低部分受到地面影响，因为大多数污染源位于地面。 ABL对城市空气污染建模很感兴趣，因为它是人类居住的大气层[ 15 ]. 我们将重点关注海拔高度，我们表示 , 包含在ABL层中。
在本文的其余部分中，我们将考虑具有常数的一级化学反应动力学。 这意味着模型中物理化学转换的平衡( 1 )可以写为 , 哪里 是反应速率。 特征时间[ 16 ]标记的化学物质 , 我们指的是 , 根据反应速率定义为以下关系：
该特征时间标度，也称为停留时间，表示代表性分子在干沉积、湿清除和化学反应等损失过程中被清除之前在大气中所花费的平均时间。 物种运输的特征距离（空间尺度）与停留时间相关。 如果排放的物种在大气中的停留时间很短（分钟-小时），它们就无法进行远距离传输，因此它们的影响将在局部范围内非常重要。 然而，其他停留时间较长（小时-天）的化学物种可能具有更广的影响区域，需要区域或大陆尺度。 特征时间刻度 与高度大气层中的干沉降有关 对于标记的物种 可根据其干沉积速度计算如下：
3.模型离散化 模型的时间集成( 1 )使用正向欧拉方法进行。 这种经典方法的选择取决于感兴趣的时间尺度，大约10分钟，以及我们感兴趣的空气污染物的物理和化学特性。 考虑第节中列出的假设 2 ，型号( 1 )时间离散化如下： 哪里 和 , . 表示离散化的时间步长。
使用有限元方法实现模型的空间离散化[ 17 , 18 ]. 它被认为是计算复杂科学和工程问题近似解的最合适方法之一，这些问题用复杂几何结构区域上的偏微分方程（PDE）表示。 该方法基于能量最小化原理，在有限维子空间中逼近无穷维空间中以变分形式表示的问题。 近似解最终被计算为离散值的集合，离散值是由离散化产生的线性系统的解向量的分量的形式。 已经熟悉经典有限元方法的读者可以直接跳到第节 4 .
有限元分析中涉及的主要步骤，已在中重述[ 19 ]是（i）问题定义，（ii）离散化（三角剖分和近似空间的定义），（iii）变分（或弱）公式的定义，（iv）系统集合（局部（每个网格元素上）和全局（整个网格上）），以及（v）所得线性系统的解。
让 是一个网格 , 哪里 指网格特征长度。 我们表示 定义在上的分段线性有限元空间 闭合边界处的值为零 . 其定义为 哪里 是满足闭边界上齐次Dirichlet条件的一阶Hilbert-Sobolev空间 . 是在上定义的连续函数集 , 结束 . 问题的变分形式( 1 )如下所示：find 这样的话
哪里 指 . 积分方程( 6 )并考虑问题中定义的边界条件( 1 )，以下弱形式成立：find 这样的话 哪里
Lax-Milgram定理[ 17 ]保证弱问题解的存在唯一性( 7 )，也可以写为：
让 是的一组基本函数 . 近似解 可以表示为 哪里 是节点值。 通过使用方程式( 9 )和( 10 )，问题的以下离散形式成立：
哪里 和 . 矩阵的条目 计算为 . 右侧矢量的分量 是 , 和矢量解 由未知系数组成 . 这里提出的有限元方法是指所谓的Galerkin近似[ 20 ]由于试验空间，其中的解决方案 所属，与测试函数的空间相同。 线性系统( 11 )这种离散化的结果可以使用LU分解等直接求解器求解[ 21 ]或广义最小残差等迭代方法[ 22 ]. 使用预处理技术通常很重要[ 22 ]，包括域分解方法[ 23 ]，用于求解线性系统( 11 )因为它通常是病态的。
4.并行实现 数学家和工程师通常使用偏微分方程问题的有限元解库。 这些软件包通常为实现高级数值方法提供了强大的工具，并支持高性能计算功能。 文献中有许多开源和流行的库提供了这些功能：感觉++[ 24 ]，FreeFem++[ 25 ]，FEniCS项目[ 26 ]，获取FEM[ 27 ]和交易。 二[ 28 ]. 对于这个项目，我们开发了一个新的健壮框架，而不是使用现有的包，因为这些包大多非常通用，因此不足以支持与此工作相关的一些特定功能。 因此，这是一个机会，可以利用灵活的方法使该框架在独立项目中发展，从而实现预期目标。
过去几十年来，并行计算机设计的巨大进步以及多核处理器和硬件加速器技术的最新进展，为解决基于仿真的应用科学和工程中的高要求复杂问题开辟了高性能计算的新时代。 我们的目标是通过求解空气质量模型来再现一些大气污染物的动力学行为( 1 )覆盖巴马科整个地区的计算域（面积约为267 公里 2 )使用足够精细的空间分辨率。 在合理的时间内实际执行此类计算需要使用并行计算[ 29 , 30 ].
在这种情况下，我们使用现代C++17编程特性开发了一个名为AirQMB的C++框架。 该模型代码的主要核心结构和并行化策略受到了中引入的并行计算框架的启发[ 19 ]. 我们考虑使用消息传递接口（MPI）的基于CPU的分布式内存编程方法[ 31 ]. 处理器内核之间的通信使用Boost进行处理。 MPI，一个有用的Boost库[ 32 ]它抽象了标准MPI上的层，以简化用户界面。
对于并行有限元分析，我们感兴趣的第一个基本步骤是网格处理。 基于从OpenStreetMap提取的地图的计算网格[ 33 ]，通过使用名为Gmsh的开源网格生成器在感兴趣区域的预处理中生成[ 34 ]. 软件包METIS[ 35 ]用于对生成的串行网格进行分区，以便网格分区的数量与将部署用于数值模拟的处理器内核的数量相对应。 这个包是从Gmsh接口调用的。 已经构建的并行（分区）有限元网格由所有处理器内核独立和并发加载。 每个处理器核心仅在本地提取与其本地秩对应的分区，再加上进程间（虚）元素，以便有限元基函数所需的所有元素在本地可用。 因此，并行有限元装配不需要进程间通信。
网格处理之后的关键步骤是有限元离散化。 它包括构造自由度（DOF）的集合（或表格），并在每个处理器内核内的局部网格上独立实现。 然后使用最小带宽编号程序在每个处理器中独立引入DOF的本地编号。 这种局部编号将通过称为局部到全局映射（LtGM）的关系与与整个问题相关的全局编号相连接。 构建这些映射需要涉及所有部署处理器的全局通信。 使用LtGM从处理器内核上本地构建的基本矩阵实现并行全局组装。 涉及矩阵和向量的串行和并行代数运算，包括求解线性系统的直接和迭代解算器，由PETSc处理[ 36 ]包装纸。
空气质量模拟所需的输入数据，尤其是风场和交通排放，使用NetCDF从每个处理器核心中的本地文件上传到计算框架中[ 37 ]库，一个支持创建、访问和共享面向阵列的科学数据的界面。 导入后，这些数据将使用健壮高效的插值例程在计算网格上进行插值。 通过使用名为Boost Program Options的有用库初始化模拟环境。 此工具使解析控制台应用程序的命令行选项变得容易。
5.数值实验 这里的数值模拟主要针对空气质量调查中常用的两种主要空气污染物。 第一种是空气动力学直径小于或等于2.5的细颗粒物（PM2.5） μ m，第二个是指空气动力学直径小于或等于10的可吸入颗粒物（PM10） μ m.这些污染物由悬浮在空气中的固体和液体颗粒物的混合物组成，是一个重要的空气质量指标，因为它们是影响短期和长期健康的最广泛分布的与交通有关的空气污染物。
我们假设沉积和清除造成的损失可以忽略不计，这意味着系数 . 最大高度固定为 . PM2.5对应的物理化学参数为 米 2 秒 –1 , 秒 –1 、和 毫秒 –1 PM10的物理化学参数由下式给出 米 2 秒 –1 , 秒 –1 、和 毫秒 –1 计算网格覆盖巴马科市，平均空间分辨率约为10 m如图所示 三 .网格元素（四面体）的数量大约等于 . 网格节点的数量约为 . 物理时间步长 固定为60秒，模拟持续时间设置为24小时。 空间离散化采用一阶有限元近似，风数据来自欧洲中期天气预报中心（ECMWF）。 作为模型输入的道路交通排放数据文件是在预处理过程中根据图中所示的路网进行的交通模拟生成的 4 。在该网络中，我们关注的是在每日高峰时段发生大规模交通拥堵的主要和最频繁的道路。 使用名为SUMO的开源道路交通仿真套件进行交通仿真[ 38 ]. SUMO输入所需的参数，如交通流量、交通密度、车辆类型、车辆发动机和平均速度，来自马里交通和基础设施部门提供的国家数据库。
模拟已在“ 计算、建模与仿真中心 “巴马科科学技术学院（FST）CCMS。 CCMS是一个集群，由十个计算节点（服务器）通过InfiniBand QDR网络互连。这些节点分为两组，每组有两种不同的机器类型。 第一组由五台Dell PowerEdge服务器组成，每个节点有两个Intel Xeon Silver 4110处理器，8核（每个核2个线程），运行速度为2.10 GHz和64 GB RAM。在第二组中，由五台HPE ProLiant服务器组成，每个节点有两个Intel Xeon E5-2623 v4处理器，4核（每个核2个线程），速率为2.60 GHz和16 GB内存。操作系统是Ubuntu服务器。
在并行计算环境中，计算可伸缩性（或伸缩性）是一种强大的工具，广泛用于表示给定并行算法的能力，以最佳地利用并行计算平台。 它与两个常见指标有关，即强伸缩性（或加速）和弱伸缩性（或者效率）。 加速是指在一个处理器内核上运行一个问题所花费的时间除以在上运行相同问题（固定大小）所需的时间 ( ) 处理器内核。 换句话说，让 是单个处理器内核上问题的执行时间 是相同问题的并行执行时间 ( ) 处理器内核。 加速 处理器内核定义为 . 效率 处理器内核计算如下 . 可能的最佳效率是 . 当加速比为线性时，即。， . 对大规模计算的这些度量标准的概括导致了相对加速和相对效率的引入。 让 是问题的并行执行时间 处理器内核和 是相同问题的并行执行时间 ( ) 处理器内核。 相对于 使用时的处理器内核 核心是 . 最佳可能的相对加速比是线性加速比，计算如下 . 相对于以下各项的效率 使用时的处理器内核 核心定义为 .
当处理器内核的数量从8个到128个不等时，此并行计算框架的强缩放和弱缩放结果如图所示 5 .
关于图中所示的强标度分析 5（a） ，使用128个处理器内核时，相对于8个内核的加速比为 而线性相对加速比，即最佳可能值，为 . 该结果对应于约93.88%的性能增益。 根据图中所示的弱标度分析 5（b） ，当使用128个处理器内核时，相对于8个处理器内核的效率为 . 正如我们所看到的，这些弱扩展和强扩展结果清楚地表明，该框架可以扩展到集群计算体系结构的128个处理器内核。
6.结论 我们提出了一个并行实施局部尺度空气质量模型的计算框架，并将其应用于巴马科市的城市部分。 我们简要描述了控制空气污染物浓度时空演变的对流-扩散-反应偏微分方程。 给出了模型的时间积分和空间离散化的数值方法。 特别强调了基于分布式内存方法的并行实现，该方法使用消息传递接口（MPI）和现代C++17编程特性。 模拟的重点是两种广泛分布的道路交通相关空气污染物，即颗粒物PM2.5和PM10。
这里给出的可伸缩性分析显示，高达128个处理器内核的速度和效率都很高。 这些结果很有希望，使我们能够重现指定大气污染物的最新动力学行为。
数据可用性 用于支持本研究结果的数据可向相应作者索取。
利益冲突 提交人声明他们没有利益冲突。
致谢 这项工作由巴马科科技大学校长通过项目资助提供支持。 作者还想感谢“ 计算、建模与仿真中心 “巴马科科学技术学院（FST）（CCMS）通过提供科学支持和访问集群上的计算资源来支持这项工作。
