目录 目录 摘要 一、前言 1 二、相关工作 1 21研究背景 1 S22高斯扩散模型.·········。·。·.···· 1 2.2.1无界空间连续点源扩散的高斯模型········。······· 2 2.22高架连续点源扩散的地表浓度。··········· 2 2.3拉格朗日扩散模型。.························ 3 2.3.1拉格朗日粒子模型 ………”……… 3 232拉格朗日烟团模型。·······。··············· 4 $2.4常见模型的分析 。。。。。。。。。。。,。。。，，。,，。。，，，。, 4 三、问题分析 4 31建模假设···· 5 311灰尘性质… 5 312迎风扩散··。··……········ 3.13有效性验证。·········…········……·· 5 S3.2符号说明·········· 5 四、模型建立 6 S41对流-扩散方程·····。······。············ 6 42有限差分法·。。。···。。···。······ 7 421基本介绍… 7 422FTCS格式..·......·.....··..···.·. 423CN格式··········…··· 8 S4.3稳定性分析···· 8 4.3.1 von Neumman方法。，。。·。。：··.。··············· 4.3.2TCS格式分析·····…······ 9 433CN格式分析1
目录 目录 I 摘要 III 一、 前言 1 二、 相关工作 1 §2.1 研究背景 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 §2.2 高斯扩散模型 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 2.2.1 无界空间连续点源扩散的高斯模型 . . . . . . . . . . . . . . . . . . . . 2 2.2.2 高架连续点源扩散的地表浓度 . . . . . . . . . . . . . . . . . . . . . . . 2 §2.3 拉格朗日扩散模型 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 2.3.1 拉格朗日粒子模型 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 2.3.2 拉格朗日烟团模型 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 §2.4 常见模型的分析 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 三、 问题分析 4 §3.1 建模假设 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 3.1.1 灰尘性质 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 3.1.2 迎风扩散 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 3.1.3 有效性验证 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 §3.2 符号说明 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 四、 模型建立 6 §4.1 对流-扩散方程 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6 §4.2 有限差分法 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 4.2.1 基本介绍 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 4.2.2 FTCS 格式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 4.2.3 C-N 格式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 §4.3 稳定性分析 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 4.3.1 von Neumman 方法 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 4.3.2 FTCS 格式分析 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 4.3.3 C-N 格式分析 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 I
五、求解方案选择 11 S5.1伪边界处理 511一维问题· 11 5.1.2边界方案 。 12 5.1.3效果对比 $5.2稳定性 15 5.2.1 二维问题与稳态解 。。，。。。，，，。，，，，，，，，，，，，，，，，， 15 5.2.2预估-校正实现C-N 5.2.3效果对比 六、三维扩散展示 18 S6.1基本结果. 19 6.1.1简化处理 19 6.1.2扩散过程 19 6.1.3风速与扩散系数 6.1.4高度影响 $6.2重力项. 子 6.2.1重力项影响.. 62.2自适应重力更新 $6.3地面边界..... 6.3.1吸收系数的影响 28 6.3.2蒙版法.··· 6.3.3 复杂情况 。。。。。。。。。。。。，。。。。。，，，。,，，，，，，，， 七、总结与讨论 污 S71总结 33 7.1.1效果实现 $7.2讨论· 34 7.2.1应用举例 34 722改进方向····，············· 34 参考文献 36 附录-文件列表 37
五、 求解方案选择 11 §5.1 伪边界处理 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 5.1.1 一维问题 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 5.1.2 边界方案 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 5.1.3 效果对比 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 §5.2 稳定性 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 5.2.1 二维问题与稳态解 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 5.2.2 预估-校正实现 C-N . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 5.2.3 效果对比 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 六、 三维扩散展示 18 §6.1 基本结果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 6.1.1 简化处理 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 6.1.2 扩散过程 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 6.1.3 风速与扩散系数 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23 6.1.4 高度影响 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 §6.2 重力项 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 6.2.1 重力项影响 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 6.2.2 自适应重力更新 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27 §6.3 地面边界 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28 6.3.1 吸收系数的影响 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28 6.3.2 蒙版法 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 6.3.3 复杂情况 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31 七、 总结与讨论 33 §7.1 总结 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 7.1.1 效果实现 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 7.1.2 模型比较 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 §7.2 讨论 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 7.2.1 应用举例 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 7.2.2 改进方向 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 参考文献 36 附录-文件列表 37 II
摘要 针对工厂烟肉排放的灰尘扩散问题，本文基于微分方程的大气扩散数学模型，建立一个 能够模拟灰尘浓度分布特征的模型，用于确定工厂周围地表灰尘浓度。该模型考虑了烟肉排 放的灰尘量、风速和地形等因素，并使用了数值方法进行求解。考虑到简捷性，在模型建立 中对现实条件进行适当简化： (1)灰尘排放：设定为烟肉点源排放，单位时间内排放的灰尘质量Q保持不变： (2)灰尘扩散：灰尘本身不发生分解、衰变等变化，满足质量守恒： (3)外界条件：气温、大气稳定，风的平均流场稳定、风速U恒定、风向固定且平直：地 势不对风场产生影响： 通过各种数值方法对模型的计算求解后，我们得到了不同条件下的灰尘浓度分布情况， 为工厂的环境保护提供了参考。 Abstract Aiming at the problem of dust dispersion emitted by factory chimneys,this paper estab- lishes a model that can simulate the distribution characteristics of dust concentration based on the mathematical model of atmospheric dispersion based on differential equations,and is used to determine the dust concentration on the surface around the factory.The model takes into account factors such as the amount of dust emitted from the chimney,wind speed,and terrain.and uses numerical methods to solve for it.Taking into account the simplicity.the realistic conditions are appropriately simplified in the model building: (1)Dust emission:set as chimney point source emission,the mass of dust emitted pe unit time Q remains unchanged; (2)Dust diffusion:the dust itself does not undergo decomposition,decay and other changes,which meets the conservation of quality (3)External conditions:air temperature,stable atmosphere,stable average flow field of wind,constant wind speed U constant,fixed and straight wind direction;The terrain is flat and does not affect the wind field: Calculating the model by several numerical methods,we can obtain the dust concentra tion distribution under different conditions,which provides a reference for the environmental protection of the factory
摘要 针对工厂烟囱排放的灰尘扩散问题，本文基于微分方程的大气扩散数学模型，建立一个 能够模拟灰尘浓度分布特征的模型，用于确定工厂周围地表灰尘浓度。该模型考虑了烟囱排 放的灰尘量、风速和地形等因素，并使用了数值方法进行求解。考虑到简捷性，在模型建立 中对现实条件进行适当简化： (1) 灰尘排放：设定为烟囱点源排放，单位时间内排放的灰尘质量 Q 保持不变； (2) 灰尘扩散：灰尘本身不发生分解、衰变等变化，满足质量守恒； (3) 外界条件：气温、大气稳定，风的平均流场稳定、风速 U 恒定、风向固定且平直；地 势不对风场产生影响； 通过各种数值方法对模型的计算求解后，我们得到了不同条件下的灰尘浓度分布情况， 为工厂的环境保护提供了参考。 Abstract Aiming at the problem of dust dispersion emitted by factory chimneys, this paper establishes a model that can simulate the distribution characteristics of dust concentration based on the mathematical model of atmospheric dispersion based on differential equations, and is used to determine the dust concentration on the surface around the factory. The model takes into account factors such as the amount of dust emitted from the chimney, wind speed, and terrain, and uses numerical methods to solve for it. Taking into account the simplicity, the realistic conditions are appropriately simplified in the model building: (1) Dust emission: set as chimney point source emission, the mass of dust emitted per unit time Q remains unchanged; (2) Dust diffusion: the dust itself does not undergo decomposition, decay and other changes, which meets the conservation of quality; (3) External conditions: air temperature, stable atmosphere, stable average flow field of wind, constant wind speed U constant, fixed and straight wind direction; The terrain is flat and does not affect the wind field; Calculating the model by several numerical methods, we can obtain the dust concentration distribution under different conditions, which provides a reference for the environmental protection of the factory. III
期末作业 Final Assignment 2023年5月28日 数学建模课程 May 28th,2023 一、前言 随着工业化的发展，工厂生产对环境污染成为了人们生活上的困扰。我国许多城市的生 活区内大烟囱林立，排放的灰尘浓烟对人体健康和环境造成了潜在的危害。因此，对于工厂 烟肉周边地表灰尘扩散浓度的研究具有重要的现实意义。 本文旨在针对目前的排放问题，通过建立模型，对工厂烟囱排放的灰尘在周围环境扩散 分布浓度等进行模拟评估，研究工厂排放污染的扩散规律，估计工厂更利于居民健康的建设 地点，为相关部门提供科学的决策依据，以保障人民健康和环境的安全。我们将通过以下几 个方面来实现研究目的： 1.分析工厂烟囱周边地表灰尘扩散的影响因素，包括气象条件、地形地貌、工厂烟肉高 度、烟肉排放速率等： 2。建立数学模型，描述灰尘扩散的过程，并考虑上述影响因素的作用： 3.考察对模型求解的各种数值条件与方式，确定求解方法 4。通过数值模拟，模拟不同气象条件、地形地貌、工厂烟囱高度、烟囱排放速率等情况 下的灰尘扩散浓度分布情况： 5.对模拟结果进行分析和评价，探讨工厂烟囱周边地表灰尘扩散的规律，并提出相应的 控制措施和建议。 二、相关工作 S2.1研究背景 1932年，微气象学家Sutton提出描述大气扩散过程的高斯方程，首次通过野外试验参 数推算了大气污染物的扩散状况。 0世纪70年代，随着有限元法、有限差分法、有限体积法等近似计算方法的发展，基于 流体力学、数值计算方法以及计算机图形学建立的流体动力学方法形成，基于Navier-Stokes 方程的三维流体力学模型预测方法得以建立，大气扩散模拟的精度得到大幅度提升，模拟过 程也更为便捷。 目前国内常用的大气污染扩散模型主要有：以高斯理论为基础、以拉格朗日方法为基础 的模型，采用Navier--Stokes方程和扩散方程的偏微分方程模型，等。 $2.2高斯扩散模型 以高斯理论为基础的扩散模式-以烟囱排放点作为原点，平均风向为x轴正方向，烟 囱垂直于地表水平面朝上为z轴正方向，建立右手坐标系Oy2,并做如下假设： 1
期末作业 2023 年 5 月 28 日 数学建模课程 Final Assignment May 28th, 2023 一、 前言 随着工业化的发展，工厂生产对环境污染成为了人们生活上的困扰。我国许多城市的生 活区内大烟囱林立，排放的灰尘浓烟对人体健康和环境造成了潜在的危害。因此，对于工厂 烟囱周边地表灰尘扩散浓度的研究具有重要的现实意义。 本文旨在针对目前的排放问题，通过建立模型，对工厂烟囱排放的灰尘在周围环境扩散 分布浓度等进行模拟评估，研究工厂排放污染的扩散规律，估计工厂更利于居民健康的建设 地点，为相关部门提供科学的决策依据，以保障人民健康和环境的安全。我们将通过以下几 个方面来实现研究目的： 1. 分析工厂烟囱周边地表灰尘扩散的影响因素，包括气象条件、地形地貌、工厂烟囱高 度、烟囱排放速率等； 2. 建立数学模型，描述灰尘扩散的过程，并考虑上述影响因素的作用； 3. 考察对模型求解的各种数值条件与方式，确定求解方法； 4. 通过数值模拟，模拟不同气象条件、地形地貌、工厂烟囱高度、烟囱排放速率等情况 下的灰尘扩散浓度分布情况； 5. 对模拟结果进行分析和评价，探讨工厂烟囱周边地表灰尘扩散的规律，并提出相应的 控制措施和建议。 二、 相关工作 §2.1 研究背景 1932 年，微气象学家 Sutton 提出描述大气扩散过程的高斯方程，首次通过野外试验参 数推算了大气污染物的扩散状况。 20 世纪 70 年代，随着有限元法、有限差分法、有限体积法等近似计算方法的发展，基于 流体力学、数值计算方法以及计算机图形学建立的流体动力学方法形成，基于 Navier-Stokes 方程的三维流体力学模型预测方法得以建立，大气扩散模拟的精度得到大幅度提升，模拟过 程也更为便捷。 目前国内常用的大气污染扩散模型主要有：以高斯理论为基础、以拉格朗日方法为基础 的模型，采用 Navier-Stokes 方程和扩散方程的偏微分方程模型，等。 §2.2 高斯扩散模型 以高斯理论为基础的扩散模式 [1–3] 以烟囱排放点作为原点，平均风向为 x 轴正方向，烟 囱垂直于地表水平面朝上为 z 轴正方向，建立右手坐标系 Oxyz，并做如下假设： 1
期末作业 Final assignment 2023年5月28日 数学建模课程 May 28th,2023 1.烟肉排放的灰尘浓度在4，之轴向上的分布符合正态分布： 2.在全空间内风速均匀： 3.排放源强连续均匀： 4.在扩散过程中灰尘符合质量守恒。 2.21无界空间连续点源扩散的高斯模型 对于下风向任意一点（，，），灰尘浓度为： C(,y,2)=A()e-av'ebz (2.1) 定义在，之方向的方差为扩散参数σg,:,它们与大气稳定程度有关，由于风的存在 它们随x增加而增加，且由等式： =d山y (2.2) 反解得到： =玩6=编 (2.3) 根据质量守恒，源强： Q-ucdvdz (2.4) 其中，U为风速，Q为单位时间的排放质量。 将式(2.1)、式(2.3)带入式(2.4)得： A()=2nUova: (2.5) 最后整理得到无界空间连续点源扩散的高斯模型： cu=2t8。即（盖+】 (2.6) 2.2.2高架连续点源扩散的地表浓度 这里考虑烟囱（高架）高度为，并假设地面对灰尘起全反射作用。 那么类似221分析，对于下风向的点(c,,2),灰尘浓度： C(z,y,2)=C1+C2 (2.7） 其中
期末作业 2023 年 5 月 28 日 数学建模课程 Final Assignment May 28th, 2023 1. 烟囱排放的灰尘浓度在 y, z 轴向上的分布符合正态分布； 2. 在全空间内风速均匀； 3. 排放源强连续均匀； 4. 在扩散过程中灰尘符合质量守恒。 2.2.1 无界空间连续点源扩散的高斯模型 对于下风向任意一点 (x, y, z) ，灰尘浓度为： C(x, y, z) = A(x)e −ay2 e −bz2 (2.1) 定义在 y, z 方向的方差为扩散参数 σy, σz，它们与大气稳定程度有关，由于风的存在， 它们随 x 增加而增加，且由等式： σ 2 y = ∫ ∞ 0 y 2 dy ∫ ∞ 0 C dy , σ2 z = ∫ ∞ 0 z 2 dz ∫ ∞ 0 C dz (2.2) 反解得到： a = 1 2σ 2 y , b = 1 2σ 2 z (2.3) 根据质量守恒，源强： Q = ∫ ∞ −∞ ∫ ∞ −∞ UC dy dz (2.4) 其中，U 为风速，Q 为单位时间的排放质量。 将式 (2.1)、式 (2.3) 带入式 (2.4) 得： A(x) = Q 2πUσyσz (2.5) 最后整理得到无界空间连续点源扩散的高斯模型： C(x, y, z) = Q 2πUσyσz exp [ − ( y 2 2σ 2 y + z 2 2σ 2 z )] (2.6) 2.2.2 高架连续点源扩散的地表浓度 这里考虑烟囱 (高架) 高度为 h ，并假设地面对灰尘起全反射作用。 那么类似 2.2.1 分析，对于下风向的点 (x, y, z) ，灰尘浓度： C(x, y, z) = C1 + C2 (2.7) 其中 C1 = Q 2πUσyσz exp [ − ( y 2 2σ 2 y + (z − h) 2 2σ 2 z )] 2
期末作业 Final Assignment 2023年5月28日 数学建模课程 May 28th,2023 是排放点源(0,0，)提供的浓度： =8-（品+】 是全反射像源(0,0，-h)提供的浓度。 计算地面灰尘浓度即令之=0，得到高架连续点源扩散的地表浓度模型： Q =8.（+）】 2 若要计算地面最大浓度，即求： Q mg C(,0.0)-Udy.exp2) (2.9) 为简化计算不妨令0，/o:为常数S,在式(2.8)中令： dcc,0,0=0 (2.10) da 求解式2.9)得到： 20S Cnar-URe (2.11) S2.3拉格朗日扩散模型 拉格朗日大气扩散模型可以细分为两类：拉格朗日粒子模型与拉格朗日烟团模型。 2.3.1拉格朗日粒子模型 该模型是基于湍流统计理论得到的扩散模型，可以细致追踪到每一个粒子在大气中的运 动状态，计算其时间、空间概率分布情况，以估算出灰尘浓度变化同。 粒子模型将粒子速度视为风场平均速度与湍流速度之和： V=V+V (2.12) 而湍流速度矢量满足： V'(t+△)=V'()R(△)+p (2.13) 其中，t为时间：△t为时间步长；R为拉格朗日自相关系数：P为蒙特卡洛分量。 在实际计算中，粒子模型将带预测区域划分为若干个网络，通过记录第i个粒子在第 个网络中的停留时间T,算出第n个网络中灰尘的浓度： C.(..)-NATAUA (2.14) 其中，Q为源强：N为灰尘粒子总数：△工，△g,△2分别为,,z方向上的网格分辨率。 3
期末作业 2023 年 5 月 28 日 数学建模课程 Final Assignment May 28th, 2023 是排放点源 (0, 0, h) 提供的浓度； C2 = Q 2πUσyσz exp [ − ( y 2 2σ 2 y + (z + h) 2 2σ 2 z )] 是全反射像源 (0, 0, −h) 提供的浓度。 计算地面灰尘浓度即令 z = 0 ，得到高架连续点源扩散的地表浓度模型： C(x, y, 0) = Q 2πUσyσz exp [ − ( y 2 2σ 2 y + (−h) 2 2σ 2 z )] + Q 2πUσyσz exp [ − ( y 2 2σ 2 y + h 2 2σ 2 z )] = Q πUσyσz exp [ − ( y 2 2σ 2 y + h 2 2σ 2 z )] (2.8) 若要计算地面最大浓度，即求: max x C(x, 0, 0) = Q πUσyσz exp ( − h 2 2σ 2 z ) (2.9) 为简化计算不妨令 σy/σz 为常数 S ，在式 (2.8) 中令： dC(x, 0, 0) dx = 0 (2.10) 求解式 (2.9) 得到： Cmax = 2QS πUh2e (2.11) §2.3 拉格朗日扩散模型 拉格朗日大气扩散模型可以细分为两类 [4]：拉格朗日粒子模型与拉格朗日烟团模型。 2.3.1 拉格朗日粒子模型 该模型是基于湍流统计理论得到的扩散模型，可以细致追踪到每一个粒子在大气中的运 动状态，计算其时间、空间概率分布情况，以估算出灰尘浓度变化 [5]。 粒子模型将粒子速度视为风场平均速度与湍流速度之和： V = V + V ′ (2.12) 而湍流速度矢量满足： V ′ (t + ∆t) = V ′ (t)RL(∆t) + ρ (2.13) 其中，t 为时间；∆t 为时间步长；RL 为拉格朗日自相关系数；ρ 为蒙特卡洛分量。 在实际计算中，粒子模型将带预测区域划分为若干个网络，通过记录第 i 个粒子在第 n 个网络中的停留时间 Tni ，算出第 n 个网络中灰尘的浓度： Cn(x, y, z) = Q ∑N i=1 Tni N∆x∆y∆z (2.14) 其中，Q 为源强；N 为灰尘粒子总数；∆x, ∆y, ∆z 分别为 x, y, z 方向上的网格分辨率。 3
期末作业 Final assignment 2023年5月28日 数学建模课程 May 28th,2023 2.3.2拉格朗日烟团模型 该模型将同一种类的不同烟尘粒子合并为烟团，将烟尘视为有限个烟团的叠加。假设烟 团的浓度分布在水平与垂直方向上均满足高斯分布，利用22中的分析，第i个烟团对空间 内任一点(r,)的浓度作用： c=9a即(-)m(-) (2.15) 其中，Q:为改烟团的质量：？为距离点源的水平距离：z为距离点源的垂直距离。 注意到空间任一点会受多个烟团的影响，且空间内排放源是连续不断的，故该点的浓度 为所有的浓度贡献之和： C(r,2)=〉C(r,) (2.16) 对各个烟团的浓度分布进行模拟或经验求解后，通过简单叠加即可实现对整体的求解。 S2.4常见模型的分析 高斯模型是基于简化假设而建立的统计学模型，其优点在于其简单易用，适用于大气污 染物浓度分布较为均匀的情况。但是，该模型忽略了地形、建筑物等因素对大气扩散的影响 在某些场景下，如建筑密集区域、狭长山谷等地方，衍射和折射现象十分重要，而高斯模型 缺乏对这些现象的解释，因此得到的浓度分布精度很低：同时高斯模型通常不适用于时间尺 度较短、风速变化的情况，在这种情况下，必须使用微分方程模型来计算污染物扩散过程。 而拉格朗日粒子模型可模拟非稳态、复杂大气流场下污染物的扩散过程，并能有效考虑 污染物在大气中传输、湍流等因素对扩散的影响，对污染物浓度分布结果准确性更高。但是 由于其模拟了每个粒子的运动过程，计算中需耗费大量的计算资源和时间，参数比较繁多 使用不太方便，适用范围较窄，主要应用于小尺寸场景的污染扩散预测与评估。 拉格朗日烟尘模型将两者的优点结合，但仍然避免不了先验假设浓度在特定方向形成了 正态分布稳态，只要考虑地形因素，这个假设往往无法成立。 于是，我们希望能得到一个介于这两者间的模型，既不预先假设整体浓度分布，也不细 致到每个粒子，而是将浓度作为整体，考虑烟尘浓度的传输。 三、问题分析 本文考虑一个高度为h的烟肉单位时间排放质量Q的灰尘在风速为U的西风作用下 在周围地表的扩散过程。我们需要建立一个数学模型来描述灰尘的扩散过程，并且通过实验 验证模型的有效性。具体来说，除了己知条件以外，我们需要在求解浓度变化前先回答以下 问题：
期末作业 2023 年 5 月 28 日 数学建模课程 Final Assignment May 28th, 2023 2.3.2 拉格朗日烟团模型 该模型将同一种类的不同烟尘粒子合并为烟团，将烟尘视为有限个烟团的叠加。假设烟 团的浓度分布在水平与垂直方向上均满足高斯分布，利用 2.2 中的分析，第 i 个烟团对空间 内任一点 (r, z) 的浓度作用： Ci(r, z) = Qi (2π) 3/2σrσz exp ( − r 2 2σ 2 r ) exp ( − z 2 2σ 2 z ) (2.15) 其中，Qi 为改烟团的质量；r 为距离点源的水平距离；z 为距离点源的垂直距离。 注意到空间任一点会受多个烟团的影响，且空间内排放源是连续不断的，故该点的浓度 为所有的浓度贡献之和： C(r, z) = ∑ M i=1 Ci(r, z) (2.16) 对各个烟团的浓度分布进行模拟或经验求解后，通过简单叠加即可实现对整体的求解。 §2.4 常见模型的分析 高斯模型是基于简化假设而建立的统计学模型，其优点在于其简单易用，适用于大气污 染物浓度分布较为均匀的情况。但是，该模型忽略了地形、建筑物等因素对大气扩散的影响， 在某些场景下，如建筑密集区域、狭长山谷等地方，衍射和折射现象十分重要，而高斯模型 缺乏对这些现象的解释，因此得到的浓度分布精度很低；同时高斯模型通常不适用于时间尺 度较短、风速变化的情况，在这种情况下，必须使用微分方程模型来计算污染物扩散过程。 而拉格朗日粒子模型可模拟非稳态、复杂大气流场下污染物的扩散过程，并能有效考虑 污染物在大气中传输、湍流等因素对扩散的影响，对污染物浓度分布结果准确性更高。但是 由于其模拟了每个粒子的运动过程，计算中需耗费大量的计算资源和时间，参数比较繁多， 使用不太方便，适用范围较窄，主要应用于小尺寸场景的污染扩散预测与评估。 拉格朗日烟尘模型将两者的优点结合，但仍然避免不了先验假设浓度在特定方向形成了 正态分布稳态，只要考虑地形因素，这个假设往往无法成立。 于是，我们希望能得到一个介于这两者间的模型，既不预先假设整体浓度分布，也不细 致到每个粒子，而是将浓度作为整体，考虑烟尘浓度的传输。 三、 问题分析 本文考虑一个高度为 h 的烟囱单位时间排放质量 Q 的灰尘在风速为 U 的西风作用下 在周围地表的扩散过程。我们需要建立一个数学模型来描述灰尘的扩散过程，并且通过实验 验证模型的有效性。具体来说，除了已知条件以外，我们需要在求解浓度变化前先回答以下 问题： 4
期末作业 Final Assignment 2023年5月28日 数学建模课程 May 28th,2023 1.烟肉排放的灰尘有怎样的性质？ 2.烟囱排放的灰尘如何在风的作用下扩散？ 3.如何验证模型的有效性？ $3.1建模假设 3.1.1灰尘性质 现实中，工厂排放的烟尘可以看作很多种不同颗粒，每种具有各异的扩散系数、质量与 表面积。扩散系数会影响扩散的速度，质量与表面积则是会影响重力与风力的作用。 为方便计算，假设烟尘颗粒相对空气足够稀薄，它们的碰撞与相互作用可以忽略不计 于是，最终的结果也可以看作每种颗粒分别独立叠加，进而只需要对特定的某种粒子考察即 可。此后不妨假设灰尘颗粒是均匀的，且其需要满足在扩散过程中质量守恒的基本要求。 此外，由于无法确定工厂具体的排放情况，假定每个时刻排出的烟尘是完全一致的，且 排出时的速度可以忽略。 3.1.2迎风扩散 由已知条件，我们假定风速V恒定，地势对风场的影响较小。在忽略了灰尘颗粒的相 互作用后，灰尘在空气中的浓度分布由四个条件决定：沿x轴正方向的风力作用、z轴方向 的重力与空气浮力复合作用、颗粒在空气中的扩散作用、边界处的情况。 这其中，风力U是给定的常数，扩散系数D也可通过实验测量，由于粒子均匀，视为 常数。因此，后文对模型应用的讨论主要侧重于重力的影响和边界的不同情况。 3.1.3有效性验证 由于我们的模型只是对单一烟尘考虑，而实际情况有各种复杂的经验公式，我们无法直 接通过测量的数据进行结果验证。 不过，我们可以将方程模拟结果与前述的高斯模型结果进行对比，并考察不同地形对结 果影响的适用性，进而验证模型的有效性。 S3.2符号说明 由于实际问题中的物理量都是存在量纲的，为了方便数学求解，需要进行无量纲化。例 如，排放强度Q的大小差别仅在于每点浓度乘常数倍，不影响分布，可设为1或其他值：而 对z坐标进行乏=2/h的无量纲化后，则可设h=1(在对不同高度的效果进行对比时可以 在1附近调整)。其他量以各自单位基准进行无量纲化后，即可以进行数学上的建模。 为方便讨论，我们的模型中采用以下的符号： 5
期末作业 2023 年 5 月 28 日 数学建模课程 Final Assignment May 28th, 2023 1. 烟囱排放的灰尘有怎样的性质？ 2. 烟囱排放的灰尘如何在风的作用下扩散？ 3. 如何验证模型的有效性？ §3.1 建模假设 3.1.1 灰尘性质 现实中，工厂排放的烟尘可以看作很多种不同颗粒，每种具有各异的扩散系数、质量与 表面积。扩散系数会影响扩散的速度，质量与表面积则是会影响重力与风力的作用。 为方便计算，假设烟尘颗粒相对空气足够稀薄，它们的碰撞与相互作用可以忽略不计。 于是，最终的结果也可以看作每种颗粒分别独立叠加，进而只需要对特定的某种粒子考察即 可。此后不妨假设灰尘颗粒是均匀的，且其需要满足在扩散过程中质量守恒的基本要求。 此外，由于无法确定工厂具体的排放情况，假定每个时刻排出的烟尘是完全一致的，且 排出时的速度可以忽略。 3.1.2 迎风扩散 由已知条件，我们假定风速 U 恒定，地势对风场的影响较小。在忽略了灰尘颗粒的相 互作用后，灰尘在空气中的浓度分布由四个条件决定：沿 x 轴正方向的风力作用、z 轴方向 的重力与空气浮力复合作用、颗粒在空气中的扩散作用、边界处的情况。 这其中，风力 U 是给定的常数，扩散系数 D 也可通过实验测量，由于粒子均匀，视为 常数。因此，后文对模型应用的讨论主要侧重于重力的影响和边界的不同情况。 3.1.3 有效性验证 由于我们的模型只是对单一烟尘考虑，而实际情况有各种复杂的经验公式，我们无法直 接通过测量的数据进行结果验证。 不过，我们可以将方程模拟结果与前述的高斯模型结果进行对比，并考察不同地形对结 果影响的适用性，进而验证模型的有效性。 §3.2 符号说明 由于实际问题中的物理量都是存在量纲的，为了方便数学求解，需要进行无量纲化。例 如，排放强度 Q 的大小差别仅在于每点浓度乘常数倍，不影响分布，可设为 1 或其他值；而 对 z 坐标进行 z˜ = z/h 的无量纲化后，则可设 h = 1（在对不同高度的效果进行对比时可以 在 1 附近调整）。其他量以各自单位基准进行无量纲化后，即可以进行数学上的建模。 为方便讨论，我们的模型中采用以下的符号： 5
期末作业 Final assignment 2023年5月28日 数学建模课程 May 28th,2023 表1：符号说明 符号 说明 单位 h 烟囱（排放点）高度 m Q 单位时间内排放的灰尘质量 kg/s U 风速 m/s 刀 扩散系数 m2/s (红，)地表位置坐标，烟肉为原点(0,0)，西风吹向为x正方向 (m,m） 垂直地表高度 m t 扩散时间 C(A,t) 灰尘浓度，表示空间点A,时间为t时的灰尘浓度 kg/(m3.s) s(r,y) 地面高度，表示(x,)处地面的之坐标 m 四、模型建立 S4.1对流-扩散方程 由流体力学，扩散作用项的表达是D△C,其中 22 △=+亦+0 是拉普拉斯算子。 而由空气运动所产生的对流项的表达是1VC,其中u(红，，云，)为流速， -(品”) 为Nabla算子。 由于x方向风速恒定，其余只有之方向速度，可化简为得到以下偏微分方程： 架=架++p(祭+器+9)+6 (4.1) 式(4.1)是一个对流扩散方程，对其有多种不同的数值求解策略同。 根据流体力学推导，扩散系数D可以表示为红，其中k是玻尔兹曼常数，T是温度， 是空气粘度，r是灰尘粒子半径，这里可不妨设其为己知常数。 U与m的系数为正，意味着均指向对应轴的负方向。U指向x轴的负方向是为方便之 后作图，而心由于主要代表重力作用，也应指向之轴负方向。 这里Co表示新生成的浓度，由单位时间生成浓度为Q,C%应满足∫Codz dydz=Q, 无量纲化后可任意指定Q为方便计算的值。由于已经假定了每个时刻排出烟尘的一致，C 与t无关。 6
期末作业 2023 年 5 月 28 日 数学建模课程 Final Assignment May 28th, 2023 表 1: 符号说明 符号 说明 单位 h 烟囱 (排放点) 高度 m Q 单位时间内排放的灰尘质量 kg/s U 风速 m/s D 扩散系数 m2/s (x, y) 地表位置坐标，烟囱为原点 (0, 0), 西风吹向为 x 正方向 (m,m) z 垂直地表高度 m t 扩散时间 s C(A, t) 灰尘浓度，表示空间点 A，时间为 t 时的灰尘浓度 kg/(m3 · s) s(x, y) 地面高度，表示 (x, y) 处地面的 z 坐标 m 四、 模型建立 §4.1 对流-扩散方程 由流体力学，扩散作用项的表达是 D∆C，其中 ∆ = ∂ 2 ∂x2 + ∂ 2 ∂y2 + ∂ 2 ∂z2 是拉普拉斯算子。 而由空气运动所产生的对流项的表达是 u · ∇C，其中 u(x, y, z, t) 为流速， ∇ = ( ∂ ∂x, ∂ ∂y , ∂ ∂z) 为 Nabla 算子。 由于 x 方向风速恒定，其余只有 z 方向速度，可化简为得到以下偏微分方程： ∂C ∂t = U ∂C ∂x + w(z, t) ∂C ∂z + D ( ∂ 2C ∂x2 + ∂ 2C ∂y2 + ∂ 2C ∂z2 ) + C0(x, y, z) (4.1) 式 (4.1) 是一个对流-扩散方程，对其有多种不同的数值求解策略 [6]。 根据流体力学推导，扩散系数 D 可以表示为 kT 6πηr，其中 k 是玻尔兹曼常数，T 是温度， η 是空气粘度，r 是灰尘粒子半径，这里可不妨设其为已知常数。 U 与 w 的系数为正，意味着均指向对应轴的负方向。U 指向 x 轴的负方向是为方便之 后作图，而 w 由于主要代表重力作用，也应指向 z 轴负方向。 这里 C0 表示新生成的浓度，由单位时间生成浓度为 Q，C0 应满足 ∫∫∫ C0 dx dy dz = Q， 无量纲化后可任意指定 Q 为方便计算的值。由于已经假定了每个时刻排出烟尘的一致，C0 与 t 无关。 6
期末作业 Final Assignment 2023年5月28日 数学建模课程 May 28th,2023 此外，假定地面边界有形状s(红，)，可能具有一定的吸收与反射比例，将在之后边界条 件的部分进行详细讨论。 S4.2有限差分法 4.2.1基本介绍 为了求解上述偏微分方程，我们采用有限差分法进行数值模拟。具体地，我们将三维空 间离散化，得到一个网格状的空间，并在每个网格点上计算灰尘浓度，根据扩散方程进行更 新。 对偏微分方程的另一个常用的离散方法为控制容积法。由于选取一定的型线可以得到和 有限差分法相同的形式，且控制容积的写法不易于从数学上进行分析，我们采取较直接的差 分方式。 假设工，y,z的模拟范围分别是1，x,h,h,[a1,剑，选取的空间步长为九，时间步长 为，则问题变为，已知此刻的三维矩阵C经*XX,如何得到下一刻的矩阵C女x’其 中xs=2+1是x方向以h为步长离散出的格点数（含边界），其他方向同理。只要能 够正确迭代，通过取之=0即可得到每一刻的地表浓度分布。 接下来介绍两个常用的迭代求解方法，也是本文实现的做法。值得注意的是，本节中介 绍的方法均不考虑边界条件，也即迭代规律只对内部的点有效。 4.2.2FTCS格式 若Ck表示h:时刻(h,jh,kh)位置的浓度，=w(h,h:),Ck表示此点对a 求n阶导数的值，由泰勒展开可推得： r=f+f-+o 2h f=f+月-f@+O h g=+)+g-月-2@+0u 2 因此有 acC+o0h) h acA-c盛t-C达+o时 2h gC4-c出+C-2C+o的 h2 7
期末作业 2023 年 5 月 28 日 数学建模课程 Final Assignment May 28th, 2023 此外，假定地面边界有形状 s(x, y)，可能具有一定的吸收与反射比例，将在之后边界条 件的部分进行详细讨论。 §4.2 有限差分法 4.2.1 基本介绍 为了求解上述偏微分方程，我们采用有限差分法进行数值模拟。具体地，我们将三维空 间离散化，得到一个网格状的空间，并在每个网格点上计算灰尘浓度，根据扩散方程进行更 新。 对偏微分方程的另一个常用的离散方法为控制容积法。由于选取一定的型线可以得到和 有限差分法相同的形式，且控制容积的写法不易于从数学上进行分析，我们采取较直接的差 分方式。 假设 x, y, z 的模拟范围分别是 [x1, x2], [y1, y2], [z1, z2]，选取的空间步长为 h，时间步长 为 ht，则问题变为，已知此刻的三维矩阵 C i xs×ys×zs，如何得到下一刻的矩阵 C i+1 xs×ys×zs，其 中 xs = x2−x1 h + 1 是 x 方向以 h 为步长离散出的格点数（含边界），其他方向同理。只要能 够正确迭代，通过取 z = 0 即可得到每一刻的地表浓度分布。 接下来介绍两个常用的迭代求解方法，也是本文实现的做法。值得注意的是，本节中介 绍的方法均不考虑边界条件，也即迭代规律只对内部的点有效。 4.2.2 FTCS 格式 若 C l i,j,k 表示 lht 时刻 (ih, jh, kh) 位置的浓度，w l k = w(kh, lht)，∂ n a C l i,j,k 表示此点对 a 求 n 阶导数的值，由泰勒展开可推得： f ′ (x) = f(x + h) − f(x − h) 2h + O(h 2 ) f ′ (x) = f(x + h) − f(x) h + O(h) f ′′(x) = f(x + h) + f(x − h) − 2f(x) h 2 + O(h 2 ) 因此有 ∂tC l i,j,k = C l+1 i,j,k − C l i,j,k ht + O(ht) ∂xC l i,j,k = C l+1 i+1,j,k − C l i−1,j,k 2h + O(h 2 ) ∂ 2 xC l i,j,k = C l+1 i+1,j,k + C l i−1,j,k − 2C l i,j,k h 2 + O(h 2 ) 7
