企业官网建设流程全解析

简介：面向数学建模竞赛参与者和核应急相关研究者的论文文档，围绕核电站泄漏后放射性气体浓度分布与扩散规律，建立连续源、瞬时源两类泄漏情形下的多种扩散模型，包括高斯烟羽模型、一维抛物型扩散模型、三维空间扩散模型、有限时间内泄漏扩散模型及高斯烟团模型，并利用 MATLAB 求解得到浓度变化规律，同时结合日本福岛核泄漏事件，从空气扩散、食品与工业产品传播等角度分析对我国东海岸与美国西海岸的影响。资源为单篇 doc 格式文档，压缩包大小 1.18MB，可完整查看赛题、模型假设、符号说明、求解过程与结果，适合数学建模备赛、参考论文结构或作为核污染扩散预测课程设计素材。目前已有 204 人浏览学习，可作为建模方法与 MATLAB 数值模拟的实用范例。

1. 这份数学建模竞赛论文到底能拆出什么：核泄漏气体扩散模型全流程拆解

拿到这份《核电站泄漏后放射性气体浓度分布规律和气体扩散模型研究》的数学建模竞赛论文时，我第一反应不是去看摘要里的公式，而是先确认它的定位：这是 2011 年陕西师范大学组织的全国大学生数学建模竞赛模拟赛题答卷，完整覆盖了从泄漏源类型判定、一维到三维扩散方程、有风修正，再到福岛案例应用的全流程。对准备数学建模竞赛的人来说，它等于把「一篇建模论文该怎么搭骨架」示范了一遍；对做气体泄漏扩散评估的人来说，高斯烟羽、高斯烟团、抛物型扩散三套模型正好对应事故性泄漏里最常见的两种源强假设。这篇拆解就按这条主线来，尽量还原每一步推导的取舍逻辑和实际计算中的参数坑。

2. 连续源还是瞬时源：高斯烟羽模型和扩散方程的选型逻辑

2.1 泄漏源判据：上来就套高斯公式是新手最常见的翻车点

拿到核泄漏扩散题，最容易犯的错误是直接把你背过的高斯烟羽公式往里代参数。数模竞赛的阅卷流程里，评委最先看的就是你有没有做泄漏源类型判断。泄漏源就两类：连续源和瞬时源。连续源指容器、管道、阀门损坏后长时间稳定泄放，特点是泄放持续时间长，浓度场最终趋于稳态，浓度只跟离源距离和空间位置有关，不再随时间变化；瞬时源指设备或容器爆炸破裂瞬间，气体在很短时间内形成气云团，特点是动态扩散，浓度既跟位置有关也跟时间强相关。

论文里给的判断方法很实用，核心是比较泄漏持续时间、平均风速和顺风扩散系数三个量的相对大小。泄漏时间相对扩散时间很短，按瞬时源处理；泄漏时间接近或超过扩散时间，按连续源处理；两者时长可比时最难办，需要具体分析。这里有个数模竞赛里的血泪经验：核泄漏类题目，十有八九都落在「泄漏时间远小于扩散时间」这个情形里，因为反应堆安全壳破裂是瞬间事件，后续渗漏虽然存在，但主要影响来自那一两个瞬间的大量释放。所以这篇文章从第 5.3 节开始就基本放弃连续源，全力推瞬时源模型。我的习惯是先画一条时间轴：泄漏总时长 T、气体扩散到关心点的特征时间、顺风扩散的特征尺度，把它们列出来比较，再决定要不要扛着连续源模型算到底。这一步不只是写给评委看的，它直接决定后面模型的复杂度——选错了源类型，后面所有公式都是白搭。

连续源和瞬时源的判断还可以用一组简明的对照表来把握，下表是我在做这类题时固定会列出的东西：

这张表的好处是建模一开始就锁定方向，后面每一步都围绕你选定的源类型展开，不会写到一半发现模型假设和题目设定矛盾。

2.2 连续源的高斯烟羽模型：有效源高与地面浓度

原文第 5.2 节建立的高斯烟羽模型Ⅰ，是连续源情形下的标准做法。这里有个关键概念叫有效源高 H，它不是烟囱的物理高度，而是排放口高度 h 加上放射性气体被热浮力和初始冲力抬升的高度 Δh。事故性泄漏往往伴随高温高压，气体喷出后有明显抬升，这个 Δh 在工程评估里非常重要，因为它直接决定地面最大浓度出现在下风多远的位置。坐标系取法是固定的：以泄漏点在地面的投影为原点，x 轴沿风向，z 轴竖直向上，有效源位于 (0,0,H)。

无风时，连续源扩散只考虑气体浓度随距离点源半径的变化，不随时间改变，形成稳态扩散。这时浓度表达式里不含风速项，也没有时间项，是一个只依赖空间距离的分布。有风时，浓度表达式里出现风速作为分母，物理含义很清楚：相同源强下，风速越大，气团被吹得越远越稀薄，浓度反而被稀释。原文对无风情况做了几层退化讨论：先令 z=0 得到地面浓度公式，再令 H=0 得到地面点源扩散公式，然后讨论地面浓度最大值的出现位置，用导函数为零求极值。考场上这一套「通解—退化—求极值」的写法非常加分，展示的不是你会背公式，而是你理解公式在什么条件下成立。

不过这里要提醒一句：连续源高斯烟羽模型的成立前提是泄漏持续不断，并且扩散已经达到稳态或准稳态。一旦泄漏停止，气体还在继续往四周扩散，这套公式立刻失效，必须切到瞬时源模型。很多选手在论文里写完连续源就往下冲，结果题目实际给的数据根本不满足连续源条件，评委一眼就能看出来你没理解物理过程。

2.3 为什么本题最后几乎都要落到瞬时源

顺着原文的思路走，第 5.2 节末尾那句话其实是全文真正的转折点：题目条件下气体扩散时间远大于泄漏时间，所以按瞬时泄漏处理。为什么现实中的核事故都符合这个判断？因为堆芯熔化、安全壳超压破裂这类事件，放射性物质的主要释放窗口很短，几小时到几十小时，而放射性气体随风扩散到数百公里外需要几天甚至更久。扩散时间比泄漏时间长一个量级以上，后面建模型只需要考虑瞬时性泄漏的动态扩散过程。

化学事故里也同理：储罐爆炸、管道爆裂，主释放阶段往往以分钟计，而气云扩散到厂区外需要更长时间。所以模型Ⅱ到模型Ⅵ，本质上都在围绕「瞬时源 + 扩散方程」这根主线做文章。看懂这一层，你就不会在高斯烟羽里死磕太久，可以直接跳到抛物型扩散方程那部分，那里才是这份论文从「套公式」转进到「建方程」的地方。

3. 从一维推到三维：模型Ⅱ、Ⅲ、Ⅳ的推导主线与 MATLAB 复现

3.1 一维抛物型扩散模型：微元法写出偏微分方程

原文第 5.3.1 节的推导是标准的教材式路径：假设气体只沿一条直线扩散，取核电站所在位置为原点，扩散方向为 x 轴，任意点 x 在 t 时刻的浓度为 C(x,t)。取一个微元段，分析时间间隔 Δt 内微元内放射性物质总浓度的变化。这里用了热传导方程的同款套路：一方面，时间内微元段内物质浓度变化可以用浓度对时间的偏导乘微元长度表示；另一方面，根据扩散定律，单位时间内沿 x 轴正向通过单位法向面积的流量与浓度梯度成正比，负号表示浓度从高处流向低处。

把「浓度变化」和「扩散流入减扩散流出」两个表达式联立，令 Δt 和 Δx 都趋于零，就得到一维抛物型扩散方程：C_t = D · C_xx，其中 D 是放射性气体扩散系数，与气体本身性质有关。这个方程你可以在任何一本传热学或流体力学教材里找到，它描述的是「浓度随时间的变化率 = 扩散系数 × 浓度的二阶空间导数」，意思是浓度分布越弯曲的地方，随时间变化越快。

求解用分离变量法，令 C(x,t) = X(x)·T(t)，代回方程后把时间部分和空间部分分开。关键细节出现在这里：分离过程中，时间方向会出现指数增长的根；如果对应的常数不为零，浓度会随时间无限增大，这显然违背物理事实，所以必须把这个模式的系数取为零。最终解出来是一个高斯型函数：C(x,t) 正比于 exp(-x²/(4Dt))，再除以一个带根号 πDt 的归一化因子。这个解的图像就是原文图 8 那种：以源点为中心，随时间不断展宽、不断变矮的高斯曲线。t 越大，曲线越扁平；t 越小，曲线越尖锐。把 t 固定，浓度随距离 x 呈对称分布，峰值始终在泄漏点。

3.2 三维抛物型扩散模型：高斯公式与质量守恒

一维解只是热身，实际泄漏发生在三维空间。原文第 5.3.2 节用高斯公式把质量守恒写成了漂亮的体积积分形式：任取一个封闭曲面，曲面围成的区域记为 Ω，时间间隔内穿过曲面进入区域的总质量用面积分表示，再用高斯公式把它化成体积积分；同一时间段内，区域内部浓度变化引起的质量增量也写成体积积分。两边相等，区域任意性导致被积函数相等，就得到三维扩散方程。

方程里出现了三个方向的扩散系数 Dx、Dy、Dz，它们不需要相等。真实大气里，水平方向的湍流扩散和垂直方向的湍流扩散强度差别很大，尤其在近地面层，垂直方向受热力层结和地面粗糙度影响，扩散能力远低于水平方向。所以论文在符号说明里专门区分了三个方向的扩散系数，这个细节在答辩时非常容易被问到：为什么不用同一个 D？标准回答是：大气不是各向同性介质，污染物在三个方向的散布速率本来就不一样。

初始条件取点源，边界条件取地面不吸收，用傅里叶变换求解，得到的解是三维高斯分布。分母里的系数 8(πt)^(3/2) 是三维归一化因子，保证任意时刻对整个空间积分得到的浓度总量等于已泄漏的总质量 Q。这里给一个验证公式的笨办法：如果忘记归一化系数，你就对浓度在整个空间做三重积分，令它等于 Q，把系数反推出来。这个办法在考场上很救急，比死记公式可靠得多。

3.3 有限时间泄漏：瞬时源叠加出模型Ⅳ

实际核泄漏都不是理想瞬时点源，泄漏总有一定持续时间。论文的修正是把时段内的连续排放看成是多个瞬时排放在空间某点造成的浓度叠加，这本质上是线性系统里的冲激响应叠加原理。做法是把泄漏持续期 [0,T] 切成很多小段，每小段 Δτ 内的泄漏当做一个瞬时点源，每个瞬时源在 t 时刻对空间任一点贡献一个高斯分布，然后把所有贡献沿时间积分。

这样得到的浓度表达式里必然出现误差函数 erf。误差函数的作用是把「还在泄漏中」和「泄漏已停止」两段的边界效应编码进去：当观测时刻 t 还小于泄漏总时长 T 时，积分上限被截断，浓度还在累积；当 t 大于 T 时，整个泄漏周期都已经算入，后续只随扩散继续衰减。这个分段行为用解析表达式写出来很啰嗦，用 MATLAB 数值积分反而更直观。这也是为什么论文里说用 MATLAB 求解模型Ⅳ。

我复现这类模型时，一般先写一个一维脚本验证思路，参数归一化之后再套三维。下面这个脚本就是干这个的，可以理解为模型Ⅳ的一维演示版：

% 一维有限时间连续泄漏的浓度模拟 % 泄漏速率 Qdot，扩散系数 D，泄漏时长 T_emit Qdot = 100; % 单位时间泄漏量（kg/s） D = 2.0; % 扩散系数（m^2/s） T_emit = 300; % 泄漏总时长（s） x = -500:5:500; % 空间网格（m） t_vec = 50:50:1000; % 观测时刻（s） C = zeros(length(t_vec), length(x)); for k = 1:length(t_vec) t = t_vec(k); dtau = 10; % 叠加步长（s），取 T_emit 的 1/30 tau = 0:dtau:min(T_emit, t); for i = 1:length(tau) dt = t - tau(i); if dt <= 0, continue; end C(k,:) = C(k,:) + Qdot * dtau / (2 * sqrt(pi * D * dt)) ... .* exp(-x.^2 / (4 * D * dt)); end end surf(x, t_vec, C, 'EdgeColor', 'none') xlabel('x (m)'); ylabel('t (s)'); zlabel('C')

逻辑说明：外层循环遍历每个观测时刻，内层循环把 [0,T] 内的连续泄漏切成多个瞬时源，每个瞬时源按一维高斯解贡献浓度，累加得到该时刻整条浓度曲线。关键参数是 dtau，叠加步长。步长越小越接近真实积分，但内层循环次数越多，计算越慢；步长太大，起始段浓度曲线会出现锯齿状跳变，看起来像数值不稳定。我一般取泄漏总时长的 1/50 到 1/100，比如 300 秒的泄漏取 10 秒，曲线已经很平滑。注意 Qdot * dtau 这个乘积的单位是 kg，对应每个切片源的源强 Q，而 Qdot 是单位时间泄漏速率。想套三维时，把分母改成 8(πdt)^(3/2)·sqrt(Dx·Dy·Dz)，指数项改成三个方向的平方和，代码结构完全不变。

4. 有风场景的模型Ⅴ与高斯烟团Ⅵ：复现浓度计算的四个避坑点

4.1 模型Ⅴ：风向确定时的三维扩散与坐标平移

有风之后，最物理的做法不是重新推一套方程，而是把无风解做坐标平移。原文第 5.4.1 节用了微小长方体流量平衡：取一个微小长方体，单位时间内沿 x 方向流入的浓度流量减去流出的浓度流量，等于长方体内部浓度随时间的变化率。这个等式里，x 方向的流量由两部分组成：一部分是风带着气体整体平移的平流通量 u·C，另一部分是浓度梯度造成的扩散通量 -D·(∂C/∂x)。y、z 方向没有平均风，只有扩散通量。

把三个方向的流量平衡式合在一起，就得到带平流项的三维扩散方程。这个方程的物理图像非常清楚：风的作用是把整个烟团往下风方向搬运，扩散的作用是让烟团在搬运过程中不断展宽。求解时做坐标变换，把观察点放到随风移动的气团中心上，无风解里的 x 坐标换成 x - ut，y、z 保持不变。换成大白话：风只是把整个高斯分布整体挪了个位置，形状展宽规律和原来的无风解一模一样。地面浓度只需要取 z=0 代入，得到的就是原文里的公式⑺。这个模型在福岛案例里用来算美国西海岸（下风方向）和我国东海岸（上风方向）的浓度，本质就是一套坐标变换加代入。

4.2 模型Ⅵ：高斯烟团与 Pasquill-Gifford 扩散曲线

高斯烟团模型是工程上最常用的简化：假设浓度在 x、y、z 三个方向上都服从正态分布，直接写出烟团表达式。这个表达式里最关键的不是前面的指数形式，而是三个扩散参数 σx、σy、σz。它们不再是扩散系数 D，而是随着下风距离、泄漏时间和大气稳定度变化的量。原文引入了 Pasquill-Gifford 扩散曲线法，简称 P-G 法，它是气象学和环境影响评价里的标准做法。

P-G 法的步骤是：根据观测到的地面风速、日照强度、云量等气象资料，先把大气稳定度划分成 A 到 F 六个级别，A 为极不稳定，B 不稳定，C 弱不稳定，D 中性，E 弱稳定，F 稳定。然后在图 4、图 5 的 P-G 扩散曲线上，按当前稳定度级别和下风距离查出对应的水平和垂直扩散参数。不同稳定度下，同样的下风距离对应的 σ 值可以差好几倍，尤其在强不稳定（A 级）和稳定（F 级）之间。所以论文里强调扩散参数与泄漏时间、大气稳定度和地面有效粗糙度有关，这三点缺一不可。

这里有一个我经常反复强调的问题：σ 是带单位的量。P-G 曲线横坐标单位是米，纵坐标单位也是米；但题目里「上风/下风公里处」用的却是公里。如果直接拿公里数值去查曲线或者代入幂函数表达式，相当于把扩散参数放大了 1000 倍，算出的浓度量级会彻底崩掉。用模型以前，第一件事就是把所有距离统一成米。

4.3 复现时最容易翻车的四个地方

下面这四条是我按照这份论文的思路复现时踩过的坑，每条都按「现象—原因—解决」写成，方便你排查自己的代码和计算过程。

坑点一：下风距离单位混用。现象：代入公式后浓度小到离谱，或者浓度曲线的形状明显不对称。原因：P-G 曲线的横坐标是米，而题目里说的是公里，直接把公里数值塞进 σ 的幂函数表达式，相当于把扩散参数擅自放大了 1000 倍。解决：进入方程前先把所有距离统一成米，查 P-G 曲线也用米，算完浓度后再讨论公里尺度的问题。

坑点二：上风处风速符号取反。现象：算「上风公里处浓度」时，结果比下风同距离处还大，显然不物理。原因：上风方向的风是逆着扩散方向吹的，气团被吹离观察点，坐标变换时应取风速为负值，而不是把空间坐标的正负号搞错。解决：按原文做法，上风处直接把公式中的风速 u 改成 -u，再看浓度随时间的变化。想偷懒的话，把坐标系整体旋转 180° 也行，但出图前要想清楚 x 轴方向代表什么。

坑点三：地面浓度计算时丢了有效源高 H。现象：z=0 代入后指数项里只有一个距离平方，没看到「源高越高、地面浓度越低」的衰减趋势。原因：高斯烟团的指数项里是 (z-H)，地面 z=0 时仍然保留 H²/(2σz²) 这一项。有人图省事把 H 设成 0，相当于把泄漏源直接放在地面上，结果近地面浓度明显偏高。解决：保留 H 项，先计算有效源高，再代入地面浓度公式。事故场景下有效源高往往是几十米量级，对近地面浓度影响很大，不能省略。

坑点四：灵敏度分析里风速变化方向搞反。现象：风速从 1 m/s 升到 3.5 m/s，算出的浓度反而升高。原因：连续源模型里风速 u 在分母上，风速增大意味着稀释增强；瞬时烟团模型里风速只控制气团平流位置，不改变峰值扩散衰减规律。两套模型混用后，容易把「风速稀释」和「风速搬运」混为一谈。解决：先确认自己用的是哪套模型；用模型Ⅰ或模型Ⅴ时把风速放在分母上，用模型Ⅵ时不要把风速写进峰值公式的分母。原文表 2 的灵敏度数据可以作为参照：风速从 1 到 3.5 m/s 区间，浓度从 1.2428 降到 0.9228，整体下降约 25.7%，这个趋势和「风速越大稀释越强」的物理直觉是一致的。如果自己的结果趋势相反，优先检查公式里的风速位置。

5. 把模型套到福岛案例里：三个校验技巧与风向改变的空间变换

5.1 三个立刻能用的校验技巧

第一个校验技巧是用对称性。无风模型里浓度分布关于源点对称，画出来应该是左右对称的钟形曲面。如果出图后明显偏向一侧，基本可以断定公式写错或者网格定义偏了。这个检查成本极低，却能在答辩前拦住大部分低级错误。

第二个校验技巧是做质量守恒积分。三维空间里把浓度对整个空间积分，任意时刻得到的都应该等于截止该时刻已释放的总质量。这个性质不随扩散系数、风速变化而改变，是模型本身的内禀属性。只要积分误差超过 1%，就回头检查时间步长、空间网格边界是不是截断得太早。我一般把计算域横向扩到浓度的百分之一以下再截断，误差基本可控。

第三个校验技巧是预判浓度量级。从福岛到我国东海岸是上千公里跨海输运，正常稀释之后浓度量级应该在 10^-8 以下。如果算出 10^-1 这种量级，不用怀疑模型，直接去查参数——大概率是扩散参数单位错了，或者把泄漏总量写成了泄漏速率，又或者把归一化因子丢了。量级判断是建模人最基本的工程直觉，也是评委最爱问的「你觉得这个结果合理吗」的底层答案。

5.2 风向改变的坐标旋转法

原文模型改进部分处理了风向可变的情况：风向改变一次，把基准坐标系旋转到新风向，再套用浓度公式；风向改变多次，就按改变时刻分段，每一段里反复做旋转和平移。实际做题时不需要真的把每个采样点都手动旋转，我写 MATLAB 时只做两件事：先把浓度场存成结构体，再写一个旋转矩阵，浓度计算前先把观测点坐标旋转到当前风向坐标系，算完再旋转回基准坐标系。旋转一段的示意代码如下：

% 风向改变 theta 弧度时，把观测坐标旋到风坐标系 theta = deg2rad(30); % 风向变化角，按实际情况改 R = [cos(theta), sin(theta); -sin(theta), cos(theta)]; X_wind = R * [x_obs; y_obs]; % 旋到风坐标系 % 算完浓度后，图像坐标仍用原始 [x_obs; y_obs] 来画

代码里 R 是二维旋转矩阵，作用是把基准坐标系下的观测点坐标旋转到当前风向坐标系。theta 取正值代表风向逆时针变化，实际使用前先确认自己的风向定义。算完浓度后再把结果存回原始坐标，画图时不用额外处理。

这个处理方式在答辩环节很加分，因为它直接回应了「实际风向不可能始终不变」这个最常被质疑的点。从那以后，我只要拿到这类扩散建模题，都会强制走一遍固定流程：先定泄漏源类型，再选基础模型，然后查大气稳定度级别，统一单位算扩散参数，最后用对称性和质量守恒做双重校验。这个习惯帮我躲开了不少翻车现场，希望帮到你。

本文还有配套的精品资源，点击获取