基于FICK定律和高斯烟羽模型的放射性气体扩散研究.doc

题目：基于FICK定律和高斯烟羽模型的放射性气体扩散研究摘要日本核污染扩散问题不仅对该国公众健康造成巨大危害，还对其对外政治关系、全球环境，乃至世界经济格局产生了深远的影响；因此，建立有效的模型模拟放射性气体的扩散，并预测放射性气体在不同地区的浓度变化情况可为决策者提供及时准确的信息，从而尽可能的减少核泄漏带来的损失。对于问题一，我们对放射性气体的扩散过程进行合理的简化和抽象，在不考虑气体受到的重力、浮力和风速的影响时，气体呈放射状向四周扩散。我们首先考虑质量守恒定律，再由Fick定律求出扩散系数D，进而得出扩散的粒子流量与其浓度梯度的正比关系，得出描述扩散情况的偏微分方程。随后，用傅里叶变换求解得到扩散方程的解，并在MATLAB中绘制此微分方程的图形，发现预测图形与东京市测得的实际数据的图形基本吻合，即离泄漏源越远浓度越低。对于问题二，要探究风速对放射性物质浓度分布的影响。风速的处理是此问题的核心，因此我们采用大气污染的经典高斯扩散模型。通过查阅相关资料，我们发现连续点源的平均气体流，其浓度分布符合正态分布规律，因此污染物浓度在y、z轴上的分布为正态分布。取烟云轴线为x方向（平均风向），得出无界情况下下风向空间某一点的浓度函数。同时考虑到泄漏点的高度，对模型进行了修正，得出的浓度分布情况与实际情况相似 。对于问题三，由问题二得到的浓度分布函数，可以分上风和下风两种情况，根据当时的实际情况，假设自然风速大于泄露的自身扩散速度，则可将第二问中的风速替换成与的线性组合即可，即下风向的速度为，上风向的速度为。将平均风速分别代入浓度函数，得出了上风向和下风向L处浓度分布函数。结果显示两个地点气体浓度变化情况与第二问得出的结果一致。对于问题四，本文参考了大量地理、气象、专题报道等资料。假设风向不随时间变化而变化，即我国东海岸一直处于上风，美国一直处于下风，并且取大气稳定度均在D时的扩散参数。中国取东海岸的上海为例，取纵向位移，美国取核电站到西海岸的直线距离为，带入第三问的模型可以得出我国东海岸山东半岛辐射浓度的预测值，由于题中的各个参数无法准确设置，所以我们将预测的浓度与东京的值相比得到一个相对值，得出了3月23日上海相对浓度为 ，美国西海岸的相对浓度为。预测结果要比实际结果大，但是与实际值相比可以接受，所以这个结果是合理的。关键词：气体扩散，Fick定律，高斯烟羽模型，扩散系数，浓度函数一、问题重述2011年3月11日, 日本近海发生9.0级地震并引发了大海啸，沿海的核电站受到破坏，开始释放出大量具有放射性的物质。由于放射性气体泄露事故的发生和发展都具有很大的不确定性，且不随人们的意志为转移，因此事故后果通常都比较严重，会造成人员伤亡并给生态环境带来毁灭性的后果。放射性气体与常规的大气污染物发生污染事故相比具有很多不同之处，正确地分析放射性气体的大气扩散过程，研究其扩散后造成的伤害时期和区域对实际的抢险救灾有很强的指导作用。 题目以日本福岛核泄漏为背景，假设有一座核电站遇自然灾害发生泄漏，给出的已知条件有：放射性气体的浓度为，扩散速度为，在无风的情况下，匀速在大气中向四周扩散, 速度为。1）请你建立一个描述核电站周边不同距离地区、不同时段放射性物质浓度的预测模型。2）当风速为 m/s时，给出核电站周边放射性物质浓度的变化情况。3）当风速为 m/s时，分别给出上风和下风公里处，放射性物质浓度的预测模型。4）将你建立的模型应用于福岛核电站的泄漏，计算出福岛核电站的泄漏对日本本土、我国东海岸、及美国西海岸的影响。二、问题分析2.1问题的重要性分析（社会背景）3月11日发生在日本的特大地震成为全世界关注的焦点，大地震发生后，福岛第一核电站机组相继发生爆炸，使人们产生了对“核爆炸”的恐惧。对此，各国为预防日本核泄漏对环境和人类的影响采取积极措施：韩国100所学校因担心“辐射鱼”临时停课；印度全面禁止进口日本食品；巴西、俄罗斯、新加坡等国相继对部分日本食品采取了进口限制措施。由此可见，日本核已对各国的经济社会和人民生活均受到了不同程度的影响，为了尽量减轻核危机造成的不便和损失，及时准确的预测核辐射对不同地区的污染程度已成亟待解决的问题。2.2有关方面在这个问题上做过的研究张子波、李自立等人曾研究了一种适合计算天然气扩散的模型，考虑了初始喷射和浮力对扩散的影响，用流体计算软件FLUENT计算一个约700m*200m*300m的扩散区域的浓度场，8核的计算机计算完毕用时24h。因此，用此类方法计算天然气泄露后扩散的浓度场时，达不到“实时”计算的效果，很容易厌恶应急指挥决策的时机，增大产生爆炸事故的。而我们建立的模型不仅具有实时性，更具有前瞻性，能较为准确的预测放射性气体在不同时间，不同地区的浓度。有学者采用有限元方法数值模拟气体的泄露扩散过程，但此方法只能针对特定装置，缺乏。我们组运用的高斯烟羽模型能推广到一般，用于描述不同背景下的气体扩散问题。2.3问题分析对于问题一，首先，确定此问题实际上是一个点源连续泄露的扩散问题，需要建立一个简单直观的模型对核电站周边不同地区、不同时段放射性气体浓度进行描述与预测。在无风的情况下，放射性气体以恒定的速率，在大气中匀速向各个方向呈放射状扩散。考虑到模型的一般性和简洁性，我们忽略了放射性气体在空气中受到的重力和浮力、核电站的高度、地面物体的高度以及地形的变化等因素，把放射性气体的影响范围假想成一个从球心不断向外各个方向均匀扩大的半球体，建立扩散模型，根据Fick定律和傅里叶变换求出放射性气体的浓度方程。对于问题二，题目引入了风速变量对放射性气体浓度的影响，选取高斯连续点源扩散模型进行分析，这也是在考虑风速时，气体扩散模型中最为经典的模型之一。在不考虑垂直风速，假设空间放射性气体的浓度服从高斯分布的前提下，运用高斯模型可以合理的计算出核电站周边地区的放射性气体浓度。对于问题三， 本文直接用第二问的模型来预测上风和下风 公里处，放射性物的浓度。对于问题四，我们将第二问的模型应用于福岛核电站的泄漏，通过参考世界地图的距离，计算出福岛核电站的泄漏对日本本土、我国东海岸、及美国西海岸的影响。三、基本假设（1） 连续泄露时放射性气体泄露的速率恒定；（2） 放射性气体在平整、无障碍的地面上空扩散；（3） 气体在空气中不发生化学反应，地面及地标地物对放射性气体无吸收；（4） 风向水平，风速和风向恒定。（5） 污染源的源强是连续且均匀的，初始时刻放射性气体内部的浓度、温度呈均匀分布。（6） 风速大于无风情况下放射性气体扩散的速度。四、符号说明4.1模型一符号说明t从释放污染物开始的时间（天）u放射性气体边缘经过的距离（公里）x放射性气体和核电站之间的距离s在时间t放射性气体的范围P核电站附近上空放射性气体的浓度C(x,t)时刻t在位置x污染物的相对浓度q(x,t)粒子流量在时间区间内物质的净减少量C放射性气体浓度Q源强水平方向扩散参数垂直方向扩散参数平均风速泄漏源抬升的高度气云出口速度出口直径环境风速五、模型的建立与求解5.1问题一5.1.1模型的建立首先，我们使用扩散模型模拟放射性气体泄露后在空气中的扩散过程，用偏微分方程来描述不同时刻不同位置放射性气体的浓度。然后，考虑到模型的一般性和简洁性，我们忽略了放射性气体在空气中受到的重力和浮力、核电站的高度、地物的高程以及地形的变化，把放射性气体的影响范围简化为一个从球心不断向外各个方向均匀扩大的半球体（如图1所示）。最后，根据Fick定律列出扩散方程，再用傅里叶变换求解放射性气体浓度的表达式。图1 放射性气体模拟扩散图扩散是由于微小的随机游动使粒子产生的散布现象。时刻t在位置x污染物的相对浓度C(x,t)由扩散方程的偏微分方差来描述。这里所说的相对浓度是指已经规范化的浓度函数，所以粒子的总质量。这将有利于强调它与概率论的联系。扩散方程是基于如下两个考虑得到的：首先，由质量守恒定律有，在位于x点宽度为的小区间内浓度的变化是由在小区间的左端点进入的粒子流量（即单位时间通过点x的粒子的个数），减去在其右端点流出的粒子流量。 如果，则在右端流出的物质快于在左端进入的物质，浓度将是减少的。在时间区间内物质的净减少量为，因此，浓度的改变量为，由此可得。令到取极限就得到方程（5-1-1）。(5-1-1)其次，考虑到Fick定律，观测数据表明，扩散的粒子流量正比于其浓度的梯度（粒子一般是从高浓度区向低浓度区扩散），于是有：(5-1-2)其中D0是个常数，称为扩散系数，其大小取决于温度、流体黏度与分子大小，并与扩散分子流动的平均速度成正比。在很多情况下，为了方便起见，用下面经验表达式来代替上面的扩散系数计算式：(5-1-3)将（5-1-1）和（5-1-2）结合就得到扩散方程(5-1-4)可以用以上方程预报放射性气体的扩散。由于偏微分方程不容易直接求解，因此先对扩散方程（5-1-4）作傅里叶变换。傅立叶变换的定义如下：(5-1-5)有以上定义得到： (5-1-6)这是一个非常简单的常微分方程（）。对于所有的k,给定初始条件，这个微分方程的解为。这个初始条件意味着在时刻t=0，污染物的云雾聚集在点x=0（的极限情形）。然后做傅里叶逆变换，得到扩散方程（5-1-4）的解，也就是放射性气体浓度随距离和时间的表达式：(5-1-7)5.1.2模型的求解利用MATLAB软件绘制出（5-1-7）式，选定某一观测点，即x值不变，作出放射性气体浓度随时间变化的图像（如图2），竖轴表示放射性气体的浓度C，横轴表示核电站泄露后的时间。图2 放射性元素预测浓度变化图然后，选定距福岛核电站约280公里的东京市作为观测点，我们访问了日本东京市气象网站，1搜集到了在东京市检测到的三种主要放射性元素（碘131、铯134和铯137）的浓度在核泄漏后一个月内的变化数据，并作出放射性元素实际浓度变化图（如图3-5）。图3 东京放射性碘131实际浓度变化图 图4 东京放射性铯134实际浓度变化图图5 东京放射性铯137实际浓度变化图将图2分别与图3、图4、图5比较，不难发现我们预测的放射性气体浓度随时间推移的变化曲线与实际情况近似吻合，有一定的关联性。由于本模型是忽略了诸多因素（尤其是风速）后的简化模型，所以预测结果与实际结果有一定的差异。5.2问题二5.2.1模型的建立5.2.1.1假设处理通过查阅相关资料，我们发现在实验和理论中，连续点源的平均气体流，其浓度分布符合正态分布规律，因此先由之前的假设得出以下结论：1）污染物浓度在y、z轴上的分布为正态分布，即在y轴和z轴上分别有 2）风只在一个方面作稳定的水平运动，即=常数。3）放射性气体在扩散中没有发生化学变化和衰减，即4）在X轴方向上，风速大于扩散速度，即平流输送作用远大于扩散作用，在湍流扩散方程中：5）一定时间后，放射性气体覆盖的区域中气体浓度分布不随时间改变，即6）地表面足够平坦，不考虑地面物体对气体的扰动。5.2.1.2浓度公式的确定首先讨论无界空间（无地面影响）的情形，取烟云轴线为x方向（平均风向）。由于湍流场在水平和垂直方向是均匀的，因此，可导出无界情况下的扩散模式。由假设1可写出污染物浓度的分布函数(5-2-1)由统计理论写出方差的表达式 (5-2-2)由假设3的连续条件可以写出(5-2-3)由上述四个方程组组成一个方程组，其中可以测量及计算的量为源强Q，平均风速及， 。未知量为：浓度C，待定函数A（x）,待定系数a,b。故次方程组可以求解。将（5-2-1）式代入（5-2-2）式 中，解之得 将（5-2-1）式和（5-2-4）式代入（5-2-3） 式中，解得将 （5-2-3）式，（5-2-4）式代入（5-2-1）式中，得(5-2-4)这个公式是无界情况下风向空间某一点的浓度，原点与放射性气体源重合。考虑地面的影响时，原点取在核电站在地面的铅直投影点上，这时高斯扩散公式如下： (5-2-5)式中H为烟云有效高度，这个公式是无界情况高架连续点源下风向空间某一点的浓度，原点与放射性气体源重合。如果地面对放射性气体完全吸收，公式中无反射项2。5.2.1.3相关参数的确定下面给出扩散系数，和有效高度H的确定过程：1）扩散系数，的确定过程：中国国家标准规定的方法：根据时间、地理位置确定日倾角、太阳高度角，利用天气条件确定辐射等级，然后利用辐射等级和风速确定大气稳定度，最后查扩散参数幂函数表，确定扩散参数。根据国家标准(GB/ T 13201 - 1991) 制定地方大气污染物排放标准的技术方法的规定，划分大气稳定度的级别，共分为6级A - F，A为极不稳定；F为极稳定。首先，根据释放源所在地的经度和纬度以及泄漏的日期和时间计算当时的太阳高度角h0；然后，由太阳高度角h0和云量查出太阳辐射等级；最后，再根据地面风速确定当时的大气稳定度，计算细节可参考文献。表1：大气稳定度的级别参考表地面风速白天太阳辐射阴天的白天或夜间有云的夜晚强中弱薄云遮天或低云0.5云量0.46CDDDDD有风时的扩散参数根据计算公式（5-2-2），并联合采用Briggs给出一套扩散参数，如表2和表3所示。表2：Briggs 扩散参数（开阔平原田野）大气稳定度ABCDEF表3：Briggs 扩散参数（城市）大气稳定度A-BCDE-F2）有效高度H的确定过程：影响泄漏源抬升高度的因素很多，主要包括：扩散气体的初始速度和方向、初始温度、泄漏口直径、环境风速及风速随高度的变化率、环境温度以及大气温度。目前大多采用半经验公式计算抬升高度。此处采用20 世纪80 年代初W ilson 根据管道破裂泄漏实验所得到的经验公式:(5-2-6)此公式适合于出口喷射方向竖直向上, 喷射路径上无障碍物的情况，这与我们的假设相同，所以在本题中可以用此方法确定放射源的有效高度。5.2.3模型的求解根据高斯扩散公式：我们在MATLAB中绘制出核电站周边放射性物质浓度的变化情况，如图6所示，我们发现放射性气体的浓度从泄露中心开始向周围逐渐降低，沿顺风方向浓度降低的较慢，而沿逆风方向浓度迅速降低。图6核电站周边放射性物质浓度的分布图5.3问题三本问可以直接套用问题二所建立的模型进行求解，由于假设中已经给出自然风速大于泄露的自身扩散速度，则可将第二问中的风速替换成与的线性组合即可，即下风向的速度为，上风向的速度为。由此可以得到第三问的浓度计算表达式：(5-2-7)研究上风和下风公里处，放射性物质浓度的预测模型。预测结果如下图7和图8所示：图7上风公里处浓度的预测结果 图 8上风公里处浓度的预测结果由上图可知，上风向处，在距核电站的距离从10到22千米的过程中，放射性气体的浓度从0.02降低到0.006；下风向处，放射系气体浓度同样从0.02降低到0.006，但范围为距核电站50到270千米的区域。5.4问题四5.4.1问题四概述对于此问题，分两步进行考虑，即核辐射对我国东海岸的影响以及对美国西海岸的影响。对于这个问题，我们假设风向不随时间变化而变化，即我国东海岸一直处于上风，美国一直处于下风，并且取大气稳定度均在D时的扩散参数。1） 对我国东海岸的影响首先，我们找出了福岛核电站到我国东海岸各个城市的距离，如下图所示：图9 福岛核电站到我国东海岸各个城市的距离取东海岸的上海为例，取纵向位移进行求解，带入第三问的模型可以得出我国东海岸山东半岛辐射浓度的预测值，由于题中的各个参数无法准确设置，所以模型得到的浓度值与准确值没有可比性，所以我们又求解出了日本岛上的一个地区（东京，距离福岛大约280）的浓度值，将这两者相比得到 ，再由网上查得江苏23日检测到的放射性碘浓度为1，而东京当天的浓度为13000 ，两者相比得，这个数据与实际值相差较大是因为我们在建模过程中将模型相对的理想化了，也就是说有很多因素没有考虑进去。2） 对美国西海岸的影响对于美国，设定传播期间的风速为，再由google地球大致算出福岛到美国西海岸的距离为，忽略纵向位移，从而可以得到预测的美国西海岸的浓度值与东京的浓度值之比，为。而根据美国气象杂志的分析，到达美国境内的放射性物质只有福岛核电站方圆80范围内的放射性强度的，由于本模型忽略了其他因素的影响，我们得到的结果要比实际结果大，但是与实际值相比可以接受，所以这个结果是合理的。六、模型的分析6.1假设的合理性分析 假设（1）由题目条件给出，在此不做和理性分析，下面对假设（2）、（3）、（4）、（5）进行分析：日本第一福岛核电站高20米，其周围建筑物大都比它低，因此不会对放射性气体的扩散造成较大的阻碍；因此，我们假设放射性气体在平整、无障碍的地面上空扩散是十分合理的，而且还有利于模型的一般性，可推广到其它地方的气体扩散预测问题。众所周知，放射性气体多为惰性气体，一般情况下化学性质稳定；因此，放射性气体在空气中不发生化学反应的假设是十分科学可靠的。同时，气体在空气中扩散的过程受气象影响明显，任何风向风速的变化都会对模型的建立和求解造成很到的困难；由此可知，假设风向水平，风速和风向恒定对于构建简单直观的模型是必不可少的。除此之外，放射性气体的浓度一般不可能发生突变，故我们可以认为污染源的源强是连续且均匀的。初始时刻放射性气体聚集在一个范围极小的区域，没有必要研究其内部的浓度、温度的分布，于是我们假设它们呈均匀分布。6.2原理误差分析本文建立的模型没有考虑风向、风速、地形的变化以及周围建筑物对气体的阻滞作用，忽略了气体在空气中收到的重力和浮力；然而，常识告诉我们，气体在空气中扩散的速度和路径对于气象变化、地物的阻碍作用十分敏感。因此，根据本文的模型预测出福岛核电站的泄漏对日本本土、我国东海岸、及美国西海岸的影响与实际的扩散情况差距较大。由此可得，理想情况下的预测模型由于误差太大，并不能运用与实际的有害气体扩散的预测中，更不能根据预测的结果进行灾后抢险决策。七、模型的检验7.1模型一的检验将模型一预测的浓度变化图形与根据实际数据在MATLAB绘出的图形相比较，我们发现两者近似吻合。说明我们建立的模型精度较高，具有极强的现实作用，可运用于有害气体泄露扩散问题的研究。八，模型的推广本文所建立的模型不仅针对放射性气体扩散模型，还可运用于其他气体的扩散问题。有如下几种不同的情况：第一，气体扩散的原因是多样的，不仅包括盛放气体的压力容器出现气象小孔，从而出现比爆炸相对缓慢的泄露事故，例如毒气罐破裂导致毒气的泄露、天然气管道老化出现的气体泄露的问题。第二，本题研究的气体扩散问题是气体的不断向外扩散，其扩散到空气中的总量不断增大；而实际的气体扩散问题可能出现气体一旦扩散后总量不变的情况九、模型的评价与优化9.1模型的优缺点分析9.1.1模型的优点 本模型为理想模型，未考虑重力和地面的阻力以及建筑物等对扩散的阻滞，由于本模型主要用来预测事故发生后可能波及的危险范围，不考虑这些因素会使估计的范围扩大，相对安全性会提高。因此，使用此模型对核电站发生泄露事故后撤离范围的估计是合适的。基于污染物的浓度分布符合正态分布这一假设的高斯模式是目前应用最广也最受欢迎的气体扩散模型之一，具有物理概念清晰、实验基础坚实、数学计算简单等优点3。我们通过一些数学方法对经典的高斯模型加以修正，以扩展它的应用范围。这一模型能够在风向变化的情况下计算扩散区域任一点的地面有害气体（包括放射性气体）的浓度3，具有通用性。不仅可以作为突发事件应急管理模拟系统的一个组成部分对各种气体瞬间扩散的范围进行估计和模拟，而且可以根据事故时的风力和风向判断危险范围的实时变化和维持时间，确定疏散范围的大小，在实际的事故救援中发挥重要的协助作用。9.1.2模型的缺点我们构建模型的时候没有注意到云内密度分层，大气湍流作用下的空气夹带明显减少，云团在上风向的距离扩散。另外，需要考虑某些放射性物质在常态条件下与水蒸气之间的物理和化学反应4。高斯烟团气体扩散模型气体覆盖范围在二维平面上的投影为椭圆形，在城市建筑物密度较高的地区，这显然是不太符合实际的。所以高斯烟团模型比较适用于建筑物普遍较低或建筑物较少的平原地区或郊区。9.2模型的优化9.2.1模型的优化方案本文运用的高斯烟团模式，成立的前提条件是在扩散过程中风向保持不变。但在实际应用中风向却是常常变化的。于是我们在已建立的模型的基础上引入坐标变换的数学方法，建立了一个能够应对风向变化条件的有害气体扩散模型。该模型包括网格化预处理、高斯模式选择和坐标变化三个主要环节。9.2.2优化模型的建立、求解与分析我们对经典的高斯模式进行改进，建立新的数学模型对有害气体扩散事故的地面浓度进行预测，通过坐标变换的思路，在经典的高斯模式上加入了对风向的考虑5。预测工作基于对有害气体扩散区域的网格化预处理，预测过程中的数据资料借助历史数据库和空间信息智能系统提供服务，预测结果将以数字和模型相结合的方式实时直观地显示在网格化的地图上。十、参考文献1吴笑，龙长江，安全评价中的气体扩散模型及应用，工业安全与环保，第31卷：5051,20052张启平，麻德贤，危险物泄露扩散过程的重气效应，北京化工大学学报，第25卷：8990，19983张梦瑶，崔晋川，适用于风向变化条件的毒气扩散模型研究，运筹与管理，第17卷：79，20084张子波，李自立，等，用修正的高斯模型计算天然气稳态扩散，石油天然气学报，第32卷：369，20105Mark M.Meerschaert,Mathematical Modeling,184186,America:China Machine Press,2009十一、附录附录一：Q=1;%