放射气体模型的预估模型
摘要
本文是以日本福岛核电站遭遇自然灾害发生核泄漏的背景而提出的。且结合了高斯烟羽模型、线性拟合，以及微分方程模型，运用MATLAB软件，分析泄漏源强度、风速、大气稳定度参数等因素对放射性气体扩散的影响，预测了放射性气体浓度在不同时间，不同地区的浓度变化，并且本文模型中数据可以根据不同的实际情况而加以改变，因而是本文的应用范围大大增加，可以适用于具有较强的应用型。
对于问题一，讨论在无风的情况下，放射性气体以s m/s的匀速在大气中向四周扩散。本问中由于不考虑风力的影响，且扩散出来的气体匀速向四周散开，这样经过任意时刻t，扩散的气体围成一个半径为st的球，且距球心位置不同的地方浓度值不同。采用列数列的表现方法，设定相同时间段t，把条件进行整理，并经过简单计算得出每段时间所预测得到的扩散距离r和浓度C。利用MATLAB软件对数据进行线性拟合，采用微分方程模型得到核电站周边放射性气体在不同地区，不同时间段的浓度变化，得出随着离泄漏源距离的延伸，最后放射性物质的浓度越来越小，趋近于零，即当x趋向无穷时，C(x,y,z,t)趋向于零；当时间趋向于无穷时，C(x,y,z,t)也趋于无穷。
对于问题二，要探究风速对放射性物质浓度分布的影响。风速的处理是问题的中心，采用大气污染的经典高斯扩散模型,实现了高斯烟团气体扩散模型的动态预测,分析计算了气体扩散过程中的各关键参数。
对于问题三，本文在问题二的基础上，结合考虑风速和放射性物质扩散速度在空间中的矢量运算，将在上风和在下风不同情况下与传播速度s之间的比较的分析，利用高斯烟羽模型对核电站周边地区的浓度进行预测，然后，利用MATLAB软件，将相关数据代入程序，我们得到核电站周边地区的浓度分布的等高曲线。
对于问题四，本文参阅整理大量气象、地理、新闻资料，选择我国东海岸典型地域---山东半岛作为研究对象，综合考虑对应海域平均风速及风向、地理距离、海水对放射性物质扩散的部分反射系数等因素，集合核电站周边的浓度等高线，可。
关键词：放射性气体扩散浓度变化高斯修正模型预测
1 问题的提出
由于重大的突发性核泄漏紧急灾害事件具有爆发性、空间分布不连续性、对周边地形和气象条件的敏感性的特点，研究核事故所释放的物质的时空分布需要高度精确的技术，但是在对于更好地保护环境有着极其重要的意义。在有一座核电站遇自然灾害发生泄漏，浓度为p0的放射性气体以匀速排出，速度为m在无风的情况下，匀速在大气中向四周扩散，速度为sm/skg/s，
。
问题一，若能建立一个描述核电站周边不同距离地区、不同时段放射性物质浓度的预测模型，这对于研究核污染模式具有重要的意义。
问题二，当风速为km/s时，给出核电站周边放射性物质浓度的变化情况，这对于研究核电站附近浓度的在实际环境下有着重要的作用。
问题三，当风速为km/s时，计算出上风和下风L公里处的放射性物质浓度
的预测模型就显得更加急迫的了。
问题四，将建立的模型应用于福岛核电站的泄漏，计算出福岛核电站的泄漏对我国东海岸的影响，这一实际的有意义预测，可以明确我们实际的污染情况，为我们的核应急决策提供技术支持。
2 问题的分析
对于问题一，在无风的情况下，放射性气体s以sm/s的速度，匀速在大气中向四周扩散。在此条件下，探求一个模型来对核电站周边不同距离地区、不同时段放射性物质浓度进行预测。我们要明确此问题研究的核扩散是点源连续泄露的扩散问题。虽然只是要求考虑在无风情况下放射性物质浓度分布，但为了使模型更贴切实际，需考虑地面反射、泄漏源有效高度等因素对浓度分布的影响。
根据“泄露放射性物质质量守恒定律”和“气体泄漏连续性原理”进行分析，发现要得出核电站周边不同距离地区、不同时段放射性物质浓度的预测模型，最后对于该方程进行分析求解。
对于问题二，为了探究风速对发生核泄漏的核电站周边放射性物质浓度分布的影响，运用概率学知识，通过图解和数学推导得出“连续点源放射性物质高斯扩散模型”。应在“连续点源放射性物质高斯扩散模型”的基础上经多次合理修正后得到更好的“优化高斯模型”。
对于问题三，该问题要求建立泄漏源上风口和下风口处放射性物质浓度的
预测模型，在参考第二问的基础上，主要考虑风速和放射性物质扩散速度在空间中的矢量运算。在对上风口分析时，要分类讨论风速和自然扩散速度之间的大小关系，当风速小于自然扩散速度时，放射性物质是无法到达上风口的。
对于问题四，应参考大量气象、地理、新闻资料，选择我国东海岸典型地域---山东半岛，作为研究对象，综合考虑对应海域平均风速及风向、地理距离、海水对放射性物质扩散的部分反射系数等因素，预测出放射性核物质与实际情况比较。
3 模型的假设
考虑到放射性气体扩散的复杂性，为简单起见，在讨论扩散模型时都作了如下假设；
（1）瞬时泄漏假定瞬时完成，连续泄漏假定泄漏速率恒定； （2）气云在平整、无障碍物的地面上空扩散； （3）气云中不发生化学反应，地面对气云无吸收； （4）为水平风向，风速和风向不随时间变；
（5）气体的传播服从扩散定律，即单位时间通过单位发祥面积的流量与他的浓度梯度成正比;
(6)气体的扩散看作空中末已连续店员向四周等强度瞬时释放气体，放射性气体在无穷空间适房的过程不发生性质变化。
4 符号说明及名词解释
4.1符号说明
符号
s k H
说明
放射性气体的传播速度 风速，单位m/s 泄漏点O距有效地面的高度
任意扩散时刻
空间任意一点的放射性物质浓度
t
C(x,y,z,t)
?i(i?x,y,z) ?i(i?x,y,z)?空间任意一点的放射性物质的扩散系数的方差 空间任意一点的放射性物质的扩散系数
空间域 空间域其体积 一规则的球面面积 在(t,t??t)内通过?的流量
?
V S
Q1 Q2
内放射性物质的增量
说明
符号
Q0
从泄漏源泄漏的放射性物质的总量
附加高度 核泄漏出口处的温度
环境温度 设地面反射系数 源强，单位为kg/s
?h Ts T0
?
Q
?x,?y,?z
分别为用浓度标准差表示的x,y,z轴上的扩散参数
沉降速度，单位为m/s
地面干沉积率 冲洗系数 放射性核素的半衰期
Vs Wd
? T0.5
4.2名词解释
烟羽又称烟云（smoke cloud)、烟流（smoke plume)：从烟囱中连续排放到大气中的烟气流。由于烟羽各部分的运动速度不同，因而其外形也千变万化。不同的烟羽形状表示污染物浓度的空间分布不同。它与大气湍流、大气稳定度、地
形地物、排放参数等有密切的关系。
动力抬升：暖气流受锋面、辐合气流的作用被迫上抬，或者在运行中受地形阻挡产生上升运动，这种空气在运动中由外力（不包括重力和浮力）使一部分空气被抬上升。
湍流扩散：是指湍流运动导致大气或水体中的污染物质或其他物质与周围洁 流体的混合。
5 模型的建立和求解
5.1问题一模型的建立与求解 5.1.1模型一的建立 模型一的建立和求解
核电站源源不断泄漏引起气体扩散船舶可以看做在无穷空间有连续点源导致的扩散过程，能够由二阶抛物型偏微分方程描述放射性气体扩散过程中浓度变化规律。本问中由于不考虑风力的影响，且扩散出来的气体匀速向四周散开，这样经过任意时刻t，扩散的气体围成一个半径为st的球，且距球心位置不同的地方浓度值不同。
讨论在无风的情况下，放射性气体以s m/s的匀速在大气中向四周扩散。核电站源源不断泄漏引起气体扩散船舶可以看做在无穷空间有连续点源导致的扩散过程，能够由二阶抛物型偏微分方程描述放射性气体扩散过程中浓度变化规律。本问中由于不考虑风力的影响，且扩散出来的气体匀速向四周散开，这样经过任意时刻t，扩散的气体围成一个半径为st的球，且距球心位置不同的地方浓度值不同。利用MATLAB软件对数据进行线性拟合，采用微分方程模型得到核电站周边放射性气体在不同地区，不同时间段的浓度变化。
假设核电站泄漏点O距离有效地面的高度为H，以核泄漏点正下方的地面为坐标原点，X轴指向下风向，Y轴水平垂直于风向轴，Z轴为铅直方向，建立空间直角坐标系，则核泄漏点的空间坐标为O（0，0，H）。
图1 空间坐标系示意图
将气体从泄露时刻记作t=0，泄漏点选为坐标原点，时刻t无穷空间中任一点（x,y,z）的气体浓度记为C(x,y,z,t)假设单位时间内通过单位法向量面积的流量与浓度梯度成正比，则有
?q???i?gradC （1）
其中?i(i?x,y,z)是扩散系数，
grad表示浓度梯度，其中的负号代表放射性物质的浓
度是由高到低的地方扩散。
假设空间域?的体积为V，包围空间域?的曲面为一规则的球面，设其表面
面积为S,外法线向量为n??(-xz,-yz,1)，
则在(t,t??t)内流通过空间域?的流量可表示为：
Q1??t??tt??s??q?nd?dt
（2）
空间域 ?所包围的区域内放射性物质的增量可表示为：
Q2????[C(x,y,z,t??t)?C(x,y,z,t)]dVV（3）
而由泄漏源泄漏出的放射性物质的总量可表示为：
Q0??t??tt????p0dVdt （4）
根据质量守恒定律和连续性原理，单位时间内通过所选曲面S的向外扩散的放射性物质与S曲面内放射性物质增量之和，等于泄漏源在单位时间内向外泄漏的放射性物质。则有：
Q0?Q1?Q2即
（5）
：
???[Cxyzt??t?CxyztdV?V?t??tt??s?q?nd?dt??t??tt????p0dVdt（6）
根据曲面积分的Gauss公式得：
??s??q?nd?????V?divqdV
（7）
则式子（6）可以转换成
????[VC(x,y,z,t??t)?C(x,y,z,t)?t??t]dV??t??tt???V?divqdVdt??t??tt????p0dVdt （8）
由于
?C?tC(x,y,z,t??t)?C(x,y,z,t)??lim?t?0
?t
?lim?t?0?
t??ttkdiv(gradC)dt?t （9）
故式子（8）即可转换成
???[?tV?C]dV??t????V?divqdV??t?????p0dV??t
（10）
即：
???[?tV?C]dV????V?divqdV?p0
（11）
根据A.Fick扩散微分方程式中：
?C?t?u?Cx?x?u?Cy?y?u?Cz?z??x?V?x22???Vy2?y2???Vz2?z2
其中：C为气体浓度；t为时间；ux,uy,uz为x,y,z方向风速；?x方向上的扩散系数。
由于无风时ux?uy?uz?0，为简化模型便于求解，假设扩散各向同性且扩散系数为常数，即?X?C?t??x?V?x22??为x,y,z
yz??Y??Z??2,则根据式（1），无风条件下扩散方程为
???Vy2?y2???Vz?z2
结合式子（12），可解出式子（11）得到结果为：
p0(4?t)1.5222?xy(z?H)?exp????4?zt??4?xt4?ytC(x,y,z,t)?(?x?y?z)0.5?????（13）
这结果表明，对于任意时刻t烟雾浓度C的等值是球面x2?y2?z2?R2,并且随着球面半径R的增加C的只是连续减少的；当R??或t??时，C(x,y,z,t)?0
此模型仅能预测在无边界空间环境中，且不考虑放射性物质释放出时的初动
量的理想化情况下不同时间里放射性物质的浓度变化情况，为使模型更贴切实际，下面分别通过考虑有效泄漏源及地方反射这两方面因素对模型进行优化。
考虑有效泄漏源对模型一进行优化
图（2）
如图2所示，核泄漏口的有效高度为h，放射性物质从泄漏口排出时由于受到热力抬升和本身动力抬升，产生的附加高度为Δh,因此有效泄漏源的高度H=h+Δh。
放射性物质受热力抬升的高度主要由其排放时本身的初始动量和因其温度高于环境温度产生的静浮力决定，同时还会受到温度和风速等因素的影响，根据烟气抬升高度的综合分析公式
?h?(0.92VsD?5.25Fbhus
0.40.6)（14）
其中： h 泄漏源的实际高度；
usDk 泄漏源出口处的风速，已知为m/s；
泄漏源出口的有效直径；
s 放射性气体的扩散速度，已知为m/s；
4m浮力通量，
VsFb/s；
根据Briggs公式，浮力通量
Ts?Ta?D?Fb???VsgTs?2?（15）
其中：Ts—核泄漏出口处的温度K；T0—环境温度K。
a：在有风(us?1.0m/s2)且释放气体温度与环境温度差≥35K(Ts0.4?Ta?35K)时，抬
升高度：
?h?(0.92VsD?0.792Qus2h0.6)
（16）
其中：Q?0.275PDVsTs?T0T0，P为大气压强hPa
b：小风(us?1.0m/s)时，且温度差≥35K(Ts?Ta?35K)，抬升高度为：
?h?5.5QdT0dZ0.25?dT0??0.0098??dZ???0.375
/m )。
（17）
其中：为泄漏源的有效高度处上的环境温度梯度(Kc：当温度差≤35K(Ts?T0?35K)时，抬升高度：
?h?2(1.5VsD?0.01Q)us?1.0m/s （18）
其中，当风速us时，取us= 1.0m/s。
由本题题目知风速为k,放射性物质的扩散速度为s,分别代入式子（16）（17），和（18），可得抬升高度
?(0.92sD?0.792Qh),适用于k?1.0m/s且Ts?T0?35K?k??0.375?dT?0.25?0?h??5.5Q?0.0098,适用于k?1.0m/s且Ts?T0?35K???dZ???2(1.5sD?0.01Q)，适用于k不限（k<1.0时，取k=1.0），Ts?T0?35K?k?其中：H为泄漏源的有效高度； K为泄漏源出口处的风速； D为泄漏源出口处的有效直径； S为放射性物质的扩散速度； Ta为泄漏源出口处的温度（K）； T0为环境温度（K）；
Q?0.275PDV20.40.6
Ts?T0sT0(P为标准大气压)
综上所述，泄漏源的有效高度为：H?h??h
考虑地面反射对模型一进行优化
地面会对排放出来的放射性物质有一定的反射作用，而放射性物质因受沉降和重力等因素的影响并不能全部被地面反射回去，因此放射性物质只能部分返回空气中，由此可知空气中的放射性物质由两部分组成，其中包括由泄漏口直接排向空气中的和由地面反射回到空气中的，假设地面反射系数为?。
根据几何物理学知识可得，通过地面反射进入大气的核扩散物质相当于虚拟泄漏源O泄露的核辐射物质在原来空间的一个浓度叠加。如果所设核泄漏点源O在距有效地面的高度为H的地方，本文求空间任一点A(x,y,z)的浓度值，
)则实际泄漏源对A点的影响部分可用exp[?(z?H在考虑反射系数?后，]来表示。22?z'虚泄漏源O对A点的影响部分可用?模型可以优化为：
'?exp[?z?H2?z2]来表示，因此式子（13）所得
C(x,y,z,t)?p0(4?t)1.5(?x?y?z)0.5222?xy(z?H)z?H?exp???????exp(?24?zt2?z?4?xt4?yt???)? ??其中：泄漏源的有效抬升高度为H?h??h。
5.2问题二模型的建立与求解5.2.1模型二的建立
问题二，针对在有风速km/s时，由于地面风速对大气扩散的影响起着至关
重要的作用，不同的风速其浓度计算方法有较大差异，因此，下面的分析从有风 （大于 1. 5 m/s）、小风（小于1. 5m/s大于0. 5m/s）和静风（小于0. 5m/ s），三种情况分别论述。在此模型中，我们同样采用模型一中的空间坐标系。 5.2.1.1 无风时
当风速k<1m/s时，可以看作是无风时的情况，即此时采用模型一的瞬时点源扩散模型求解。即此时的浓度场有：
C(x,y,z,t)?p0(4?t)1.5222?xy(z?H)?exp????4?zt?4?xt4?yt?(?x?y?z)0.5?????（14）
再利用部分反射倾斜烟云模式，可得出该情况下的浓度分布：
Tc(x,y,z,H)??(2?)0p03/2?x?y?exp(?z(x?uT)2?2x2)exp(?y22y2?){exp[?(z?H?VsT)2?2z2]??exp[?(z?H?VsT)2?2z（15）
上式中假定?x??y，T代表积分时间，即静风持续时间。通常按照每小时排放6个烟团进行积分，即可达到需要的精度。 5.2.1.2 有风时
.即：当风速大于1.0m/s时，选取高斯连续点源扩散模型进行分析，这也是在所
有有风速时，放射性物质扩散模型中最常见也最方便的一种泄漏方式。 因此时的风速u不变，y,z方向上的速度为0，则此时的扩散微分方程为：
u?c?x??y??c?y22??z??c?z22?p0（16）
其中
模拟连续点源烟羽扩散的定解条件是
c?0当x,y,z??
??ucdydx??Q
p022?y(z)?exp????4?yt4?zt?则解是：
C(x,y,z,t)?4?tu(?y?z)0.5?????（17）
以上是一三维高斯函数式，表示连续排放时的烟团浓度呈以烟团中心为原点的正态分布，以?x2,?2?z2表示分布的方差，他们和扩散系数的关系是
y,?x?2?xt?22y?2?yt?z2?2?zt（18）
由于高斯连续点源扩散模型有边界点源，这里设泄漏源有效高度为H，取其地面投影为坐标原点，x轴指向风向，继续考虑地面反射作用，
图4 边界点源示意图
则此时的（17）式可以改成为
c(x,y,z,H)?p02?k?y?zexp(?y22y2?){exp[?(z?H)2?2z2]?exp[?(z?H)2?2z2]}（19） 分别为
其中，Q为源强，kg/s；c为污染物质量浓度，kg/m；?x,?y,?z3用浓度标准差表示的x,y,z轴上的扩散参数；k为泄漏高度的平均风速，m/ s；H为泄漏有效高度，m。
考虑到动力抬升的作用，将有效高度（即有效烟云高度）He?H??H。
图5 动力抬升示意图
因此（19）式可以写成为：
c(x,y,z,H)?p02?k?y?zexp(?y22y2?){exp[?(z?He)2?2z2]?exp[?(z?He)2?2z2]}（20）
式（20）是高架连续点源的一般解析式，由（20）式可以推导出各种条件下的常用大气扩散模型。 (1).泄露源为近地面
当为地面连续点源时，即He=0时，由（20）得：
c(x,y,z,0)?p0exp(?zy2?22y?k?y?)exp(?z22z2?)（21）
地面连续点源在地面上任一点产生的浓度为：
c(x,y,0,0)?p0exp(?zy22y?k?y?2?)（22）
（2）．泄漏源为高架点时
当为高架连续点源时，地面任意一点浓度模型。此时令（20中的）z=0 即可得到高架连续点源的地面浓度模型：
c(x,y,0,He)?p0exp(?zy22y?k?y?2??He2?22z)（23）
5.2.2模型的计算求解 5.2.2.1模型参数的确定
为了求解出核电站周边放射性物质浓度的具体的变化情况，则应首先把大气
条件分类，从而确定出每类天气条件下的扩散参数?y,?z的值。在事故过程模型中所需参数的选取及确定十分重要，通常情况下气象参数的选取是利用该地区多年气象资料，采取工业安全与环保统计的方法进行有关参数的确定，而其他扩散参数是以实际测定为基础的。 a:大气稳定度的计算
对于大气条件的分类，根据国家标准(GB/ T 13201 -- 1991) )制定地方大气污染物排放标准的技术方法的规定，划分大气稳定度的级别，共分为6级A ~ F，A为极不稳定；F为极稳定。其具体分类情况如下：
表 1 大气稳定度的级别参考表
地面风速
(m?s?1白天太阳辐射 强 A B B-C C-D C
中 A -B B B-C C-D D
弱 B C C
阴天的白天或夜
间 D D D
有云的夜晚
薄云遮天或低云≧
0.5 -
E D D
云量≦0.4 - F E D D
)
<3 2-3 3-5 5-6 >6
D D D
D
D
b：扩散参数的计算
有风时的扩散参数，(?y,?z)的确定采用Briggs给出一套扩散参数，则每一个大气稳定度所对应的扩散参数如表1和表2所示。
表2 Briggs 扩散参数（开阔平原田野）
大气稳定度
A B C D E
?y
0.22x(1?0.0001x)0.16x(1?0.0001x)0.11x(1?0.0001x)12?z
0.2x 0.12x
0.08x(1?0.0002x)121212
0.08x(1?0.0001x)0.06x(1?0.0001x)120.06x(1?0.0015x)12120.03x(1?0.0003x)
F
0.04x(1?0.0001x)12
0.016x(1?0.0003x)
表3 Briggs 扩散参数（城市）
大气稳定度
A-B C D E-F
?x
0.32x(1?0.0004x)?12?y
0.24x(1?0.0001x)?12 0.22x(1?0.0004x)0.16x(1?0.0004x)0.11x(1?0.004x)?120.20x
0.14x(1?0.0003x)0.08x(1?0.0015x)?12?12
?12
?12则欲求解出核电站周边放射性物质浓度的变化情况，只需将各个大气稳定度所对应的?x,?y带入计算即可。
5.2.2.2模型的求解思路整合
图6 模型二计算流程图
5.2.3对高斯浓度计算模型的修正
利用标准的高斯模型，计算大气中放射性核素云团的扩散，目前已有较多的研究，而实际的核素扩散过程还存在着粒子的重力沉降、雨洗沉积以及核素衰变等的因素对浓度分布的影响。因此，有必要对高斯浓度扩散模式进行修正，方可较为准确、真实地反映实际的核素扩散规律。 a：考虑干沉积时的连续点源扩散
粒径大于10 m 的粒子有明显的重力沉降，粒子的沉降速度取决于空气阻力和重力平衡，可用斯托克斯公式表示：
Vs??gD18?2（24）
式中，?：粒子密度，kgg/m3；
：重力加速度，9. 806 5m/s；
;
可取1.8?10?5kg/(m?s)；
D：粒子直径，m?：空气的动力粘性系数, Vs：沉降速度，m/s，含碘放射性核素的干沉积速度为Vs= 1. 1 cm/ s, 该
值也可由式( ) 计算得到。
因为在扩散过程中同时有重力沉降的位移迭加到羽流中心线上，中心线就会向下倾斜，所有粒子相当于在下倾的中心线上扩散。该类扩散和沉降的迭加可认为是羽流运行过程中，实源以Vs的速度向下移动，在x处向下移动的高度为
Vst?Vsxu，即源高由H降到了H?Vsxu。实际上，由于大气湍动及其他动力作用，
地面不是全吸收表面，应考虑地面的反射作用。但毕竟存在粒子的沉降作用，又不是全反射，因此，需对反射项乘以一个反射系数?(?源高度也变成了H?Vsxu?1)，由于反射项的有效
，故相应的浓度计算公式为：
(z?H?2yc(x,y,z,H)?Q2?k?y?zexp(?y2?2Vsxu2z)2(z?H?]?exp[?2?zVsxu2)2){exp[?2? ]}（25）
式( 20) 可得到经过干沉积作用后的浓度分布结果。式( 20) 中, 反射系数?需外部给定，?(??1)，通常对于放射性核素可取0.5。
b:地面的干沉积量
根据扩散理论和动量传递的普朗特理论，可以导出干沉积的地面沉积量为：
Wd?Vs?c(x,y,0)（27）
式中，c(x,y,0)前面各种情形计算的地面污染物的浓度；
Vs干沉积速度，m/s
Wd地面干沉积率，Bq?m?2?s?1。
c: 雨洗作用
降雨对烟羽中的颗粒物及气溶胶具有清洗作用，可溶性气体与蒸汽亦可溶于雨水中，降雨过程造成的这类湿沉积是导致放射性气溶胶和气体向地面沉积的另一重要机制。通常以冲洗系数?(S?1)，描述降雨对烟羽中污染物清洗作用的大小。
?与雨强的关系可以表达为：??aIb（28）
式中，I为雨强(mm/ h) ；
a，b为经验系数。
式中，按释放物质为含碘、不含碘情况分别取值。对含碘物质，取
a?8?10?5,b?0.6；对于不含碘物质，取a?1.2?10?5,b?0.5。
对于湿沉积导致的烟羽耗减，可采用湿沉积耗减因子对源强Q进行修正，有
Q(x)?Qexp(??xu)（29）
d：放射性衰变的影响
除了干沉积和湿沉积外，放射性物质的衰变也是影响大气中核素浓度分布的主要原因。由于放射性物质服从简单的衰变规律，其浓度随时间的变化可由下式计算：
c?c0e??t（30）
式中，c0为初始浓度；
?衰变常数；
t经过的时间；
c0可由无衰变的浓度公式计算。
参照前面的模式，同样采用衰变耗减因子来对源项进行修正。经过理论推导，可以得出：
Q(x)?Qexp(?0.693x3600T0.5u)（31）
其中，T0.5放射性核素的半衰期。
如以上情况都存在，应同时考虑源项修正，即式(24) 和式( 26) 相乘后代替Q5.3问题三模型的建立与求解
5.3.1模型三的建立和求解
在问题二的基础上，结合风速和放射性物质的扩散速度在空间中的矢量运算，将在上风和下风情况下的风速与放射性物质的传播速度s之间的比较分析，知道必须在满足风速大于放射性物质的传播速度的情况下，放射性物质才能到达上风口。
Ⅰ、上风情况下，满足条件k m/s≥s m/s，将模型二中的求解结果中的风速k用
k-s代替，即可求出在上风情况下放射性物质浓度的变化：
c(x,y,z,H)?Q2?(k?s)?y?zexp(?y22y2?){exp[?(z?He)2?2z2]?exp[?(z?He)2?2z2]}
（31）
Ⅱ、下风情况下，风速的大小对放射性物质能否到达风口均无影响，将模型二中的风速k用k+s代替，即可求出在下风情况下放射性物质浓度的变化：
c(x,y,z,H)?Q2?(k?s)?y?zexp(?y22y2?){exp[?(z?He)2?2z2]?exp[?(z?He)2?2z2]}（32）
模型四的建立和求解
在结合模型二、模型三的条件下，在参阅整理大量气象、地理、新闻资料，选择我国东海岸——山东半岛作为研究对象，综合考虑对应海域平均风速及风向、地理距离、海水对放射性物质扩散的部分反射系数等因素，并通过C++编程模拟计算，预测出放射性核物质将经过6.5天到达我国东海岸——山东半岛，且131I浓度值为：0.100mBq?m?3与，与实际情况比较吻合。
6结果分析与检验
由于在前四问题中均无具体的数据，所以对于结果，在检验上存在一些麻烦，但是在模型三预测上下风的浓度具有较好的实际意义，而在模型四中，对于预测出浓度的变化和对我国的影响也较为合理，故该模型具有可靠性和参考性，可以进行推广和其它方面的应用。
7模型的优缺点与改进方向
优点：该方法可计算核事故中连续点源和瞬时点源在不同气象、地形条件下的浓度分布。该结果在核事故的应急救援过程中，对救援人员划定警戒区和确定周围居民的疏散范围具有重要意义，并可为制定救援方案和应急决策提供科学依据。放射性云团在空中迁移和扩散提供浓度的定量描述，对可能发生的核事件的放射性核素浓度监测及监测时间范围提供相关信息。
缺点：模型在建立时只考虑了两种下垫面情况，即农田和城市，实际的扩散情况可能因地形条件的复杂化而有所改变。由于扩散是一个非常复杂、影响因素众多的过程，上述方法还有一定的局限性，有待于进一步检验。另外，该方法只给出了浓度的计算结果，没有换算成人的吸收剂量，今后将对其计算过程和参数的选取进行多方面验证，使其进一步完善。
参考文献
[1]金勇进,统计学,北京：中国人民出版社,2010。
[2]梅长林,周家良,实用统计方法,北京:科学出版社,2002。
[3]罗艾民,魏利军.有毒重气泄露安全距离数值方法[J].中国安全科学学报, 2005,15（8）: 98--100
[4]周燕等,Matlab在统计与工程数据分析中的应用,北京:电子工业出版社,2010。 [5]丁信伟,王淑兰,徐国庆.可燃及毒性气体泄漏扩散研究综述.化学工业与工程,1999,2(16):118—122.
[6]胡世明.气体释放源的三维瞬态重气扩散研究.劳动保护科学技术,2O02,3(20):28—30．
[7]蒋军成,潘旭海.描述重气扩散的一种新模型.南京工业大学学报,2OO2,I(24):41—46．
附录
附录一：模型二在上下风L公里的浓度求解源代码
#include
voidqiuxy(char wending ,double x,double&dy,double&dz) {
switch(wending) }
void main() {
double c1,c2,x,z,t,p0,H,dx,dy,dz,a,s,k,I,Ts;
chardaqi;
cout<<\请按以下格式输入：\cout<<\下风口L处坐标（米）：\cin>>x;
cout<<\地点的纵坐标（米）：\{
case 'A': dy=0.22*x/(sqrt(1+0.0001*x)), dz=0.2*x; break; case 'B': dy=0.16*x/(sqrt(1+0.0001*x)), dz=0.12*x; break;
case 'C': dy=0.11*x/(sqrt(1+0.0001*x)), dz=0.08*x/(sqrt(1+0.0002*x)); break; case 'D': dy=0.08*x/(sqrt(1+0.0001*x)), dz=0.06*x/(sqrt(1+0.0015*x)); break; case 'E': dy=0.06*x/(sqrt(1+0.0001*x)), dz=0.03*x/(1+0.0003*x); break; case 'F': dy=0.04*x/(sqrt(1+0.0001*x)), dz=0.016*x/(1+0.0003*x); break; }
cout<<\
cout<<\地点的铅直坐标（米）：\
cin>>z;
cout<<\扩散时间（秒）：\cin>>t;
cout<<\地面反射系数：\cin>>a;
cout<<\放射性气体扩散速度（m/s）：\
cin>>s;
cout<<\风速（m/s）：\cin>>k;
cout<<\雨水强度（毫米/小时）：\
cin>>I;
cout<<\核物质衰变周期（天）：\
cin>>Ts;
cout<<\初始浓度（千克/立方米）：\
c1=p0/(2*PI*(s-k)*dy*dz)*pow(e,-(0.693*x)/(3600*(Ts*24)*(k+s))-(8*e-5*pow(I,0.6)*x)/(k+s))*(pow(e,-((z-H)+(0.01*x/(s-k)))*((z-H)+(0.01*x/(s-k)))/2*dz*dz)+a*pow(e,-(z+H+(0.01*x/(s-k)))*(z+H+(0.01*x/(s-k)))/(2*dz*dz)));
cout<<\下风口L处放射性物质浓度= \
c2=p0/(2*PI*(s-k)*dy*dz)*pow(e,-(0.693*(-x))/(3600*(Ts*24)*(k+s))-(8*e-5*pow(I,0.6)*(-cin>>p0;
cout<<\有效高度（米）：\cin>>H;
cout<<\大气稳定度等级：\cin>>daqi; qiuxy(daqi,x,dy,dz); dx=dz;
x))/(k+s))*(pow(e,-((z-H)+(0.01*(-x)/(s-k)))*((z-H)+(0.01*(-x)/(s-k)))/2*dz*dz)+a*pow(e,-(z+H+(0.01*(-x)/(s-k)))*(z+H+(0.01*(-x)/(s-k)))/(2*dz*dz)));
cout<<\上风口-L处放射性物质浓度= \}
//-x*x/(4*dx*t)-y*y/(4*dy*t)
//-(z-H)*(z-H)/(4*dz*t)-a*(z+H)*(z+H)/(4*dz*t)
附录二：核泄漏模拟图的绘制 Q=1;%输入强源 u=2;%输入风速 d=1;%步长
Z=0.4;%地面粗糙参数 x=10:d:100;%下风向距离 y=-100:d:100;%横风向距离 [x,y]=meshgrid(x,y);
by0=0.08*x.*(1+0.0001*x).^(-1/2); bz0=0.06*x.*(1+0.0015*x).^(-1/2); by=by0.*(1+0.38*Z);
fz=(2.53-0.13*log(x)).*(0.55+0.042*log(x)).^(-1).*Z.^(0.35-0.03*log(x)); bz=bz0.*fz;
tempy1=-y.*y./by./by./2; tempy2=2.718282.^(tempy1); c=Q/pi/u*((by.*bz).^(-1)).*tempy2 Cs=20;%输入求解条数 contour(x,y,c,Cs); shadinginterp; colorbar; grid;
title('核泄漏模拟图')