基于不确定理论的退役铀尾矿库环境稳定性分析及预测方法

摘要

本文公开了一种基于不确定理论的退役铀尾矿库环境稳定性分析及预测方法，稳定性分析包括铀尾矿库环境稳定性指标体系的建立、各指标稳定区间的计算、铀尾矿库各坝段环境稳定率的分析；环境稳定性预测包括由稳定状态转化为不稳定状态所需时间预测、由不稳定状态转化为稳定状态所需时间预测。得出了铀尾矿库主要环境污染指标的稳定区间，并计算了铀尾矿库在不同时间点的环境稳定率，以精确的数学语言定义了铀尾矿库环境稳定性的模糊性概念，实现了对铀尾矿库环境稳定性定性分析和定量计算的有机结合，为退役铀尾矿库安全管理的及决策分析提供了理论依据。

说明书

技术领域

本发明属于能源及信息处理技术领域，涉及一种基于不确定理论的退役铀尾矿库环境稳定性分析及预测方法。

背景技术

随着我国对能源的需求，将产生大量的铀尾矿和废石，用于储存铀尾矿和废石的铀尾矿库对环境和公众的潜在安全问题非常突出(铀尾矿库是指筑坝拦截谷口或围地构成的，用以堆存铀矿废石、尾渣的工业场所)。虽然经过退役治理后其周边环境质量和生态状况得到了明显的改善，但由于其放射性核素及其他化学有毒元素在覆盖治理及闭坑后的相当长时间内依然通过多种途径扩散、迁移，溃坝、塌方等事故时有发生，严重威胁公众生命财产安全；另一方面由于土地资源匮乏，需尽快恢复到采冶前的生态环境水平，确保退役铀矿冶设施所在区域的可持续发展，因此，加速铀尾矿库的安全稳定化进程成为铀矿冶设施退役治理的重要研究课题。

无论是从保护生态环境还是保障公众生命财产安全的角度来讲，对退役铀尾矿库稳定化进程进行研究和分析都具有十分重要的意义。然而，目前对铀尾矿库稳定性的分析仅仅只考虑了尾矿坝的力学稳定性，仅以坝体的力学稳定性指标作为退役铀尾矿库安全与否的评价指标已经无法满足退役铀尾矿库生态环境保护和可持续发展的需要。针对上述问题和铀尾矿库的环境现状，本研究基于不确定理论提出了一种退役铀尾矿库环境稳定性的分析及预测方法，进一步丰富和完善退役铀尾矿库综合评价技术，对铀尾矿库的退役治理和环境效应有着积极的意义和影响。其中，不确定理论是从测度论观点出发，具有规范性、自对偶性、单调性、次可加性和乘积测度公理的数学系统。现有技术采用如下方法：

1.铀尾矿库环境稳定性的分析。

利用概率论中期望与方差的理论得出各指标渐近稳定区间，运用概率统计中的切比雪夫不等式，找出影响各环境指标的稳定概率大于或等于α时的稳定区间。

定义1：设随机变量X为某个指标的观察值，X的取值为x₁,x₂,L,x_n，Y_k(k＝1,2,Λ,n)为某个指标在时间[t₁,t₂]内连续k次观察值构成的随机变量序列，若Y_k的数学期望EY_k和方差DY_k都存在，且对任意的α∈[0,1]，存在ε>0，使得不等式

P(|Y_k-EY_k|<ε)≥α成立，

则称该指标在[t₁,t₂]上(α,ε)稳定，称[EY_k-ε,EY_k+ε]是该指标在时间[t₁,t₂]内对应概率为α的稳定区间。结合切比雪夫不等式和指标v_i稳定的定义，可以得出α、ε之间的关系为$>\alpha =1-\frac{D\left(Y_k\right)}{{ϵ}^2},$>即$>ϵ=\sqrt{\frac{D\left(Y_k\right)}{1-\alpha }}.$>

为了对铀尾矿库的环境稳定性进行定量分析，给出了环境稳定率的定义及计算方法。

定义2：设x_i为各指标在某一时刻的观察值，w_i为各环境指标的权重，为各指标稳定区间的示性函数，为指标A_i稳定区间的下限值，为指标A_i稳定区间的上限值，环境稳定率ρ为某一时刻尾矿库环境的稳定性程度，且

$>\rho =w_1{\chi }_{A_1}\left(x_1\right)+w_2{\chi }_{A_2}\left(x_2\right)+\Lambda +w_m{\chi }_{A_m}\left(x_m\right)=\underset{i=1}{\overset{m}{\Sigma }}{w}_i·{\chi }_{A_i}\left(x_i\right)$>

其中，示性函数可表示为

定义3：设时间[t₁,t_k]内某环境指标随机变量序列Y_k的所对应的稳定区间$>[a,b]=\underset{i=1}{\overset{k}{I}}[{\mathrm{EY}}_i-{ϵ}_i,{\mathrm{EY}}_i+{ϵ}_i],$>若区间$>[a^{*},b^{*}]=[a-\frac{b-a}{l},b+\frac{b-a}{l}]$>使得

$>\underset{i=1}{\overset{k}{I}}[{\mathrm{EY}}_i-{ϵ}_i,{\mathrm{EY}}_i+{ϵ}_i]⋐[a^{*},b^{*}]⋐\underset{i=1}{\overset{k}{Y}}[{\mathrm{EY}}_i-{ϵ}_i,{\mathrm{EY}}_i+{ϵ}_i]$>成立，

则称区间[a^*,b^*]为该指标稳定性的弹性区间。

2.铀尾矿库环境稳定性的预测

对退役铀尾矿库的环境稳定化进程进行预测，就是对各环境指标的稳定化趋势进行预测，通过构造随机变量和随机变量序列得到各指标的稳定区间，构造预测函数计算新的稳定区间与原有的稳定区间作比较，预测环境指标由稳定状态转化为不稳定状态和由不稳定状态转化为稳定状态所需的时间。由稳定状态转化为不稳定状态的预测函数仅考虑直线增长型、指数增长型和周期变化型，若新的稳定区间包含于原来的稳定区间，则该指标仍处于稳定状态，若新的稳定区间突破了原来的稳定区间，则该指标已转化为不稳定状态；由不稳定状态转化为稳定状态的预测函数只考虑指数衰减型，若原来的稳定区间包含于新的稳定区间，则该指标仍处于不稳定状态，若新的稳定区间进入了原来的稳定区间，则该指标已转化为稳定状态。两种情况下指标突破或进入稳定区间所需的时间s需满足的关系式分别为：

$>\underset{i=1}{\overset{k}{I}}[{\mathrm{EY}}_i-{ϵ}_i-\frac{b-a}{l},{\mathrm{EY}}_i+{ϵ}_i+\frac{b-a}{l}]⋐\underset{i=1}{\overset{k+s}{I}}[{\mathrm{EY}}_i-{ϵ}_i,{\mathrm{EY}}_i+{ϵ}_i]$>

或$>\underset{i=1}{\overset{k+s}{I}}[{\mathrm{EY}}_i-{ϵ}_i,{\mathrm{EY}}_i+{ϵ}_i]⋐\underset{i=1}{\overset{k}{I}}[{\mathrm{EY}}_i-{ϵ}_i-\frac{b-a}{l},{\mathrm{EY}}_i+{ϵ}_i+\frac{b-a}{l}]$>

其中：EY_i表示环境指标前i个月监测值的数学期望，k表示环境指标监测的第k个月份，同时也是环境稳定性预测的起点时刻，表示指标稳定区间的弹性范围，对应的区间$>\underset{i=1}{\overset{k}{I}}[{\mathrm{EY}}_i-{ϵ}_i-\frac{b-a}{l},{\mathrm{EY}}_i+{ϵ}_1+\frac{b-a}{l}]$>称为指标的稳定弹性区间。

现有针对铀尾矿库环境方面相关的研究大多仅考虑了某些环境因素，尚未从一个大的角度出发进行全面而系统的分析，而对铀尾矿库环境稳定性的分析和研究还存在相当大的一片空白，目前的一些评价方法和理论还无法解决这一技术问题。

不确定理论在尾矿库领域的应用还只是局限在尾矿坝坝体边坡稳定性，面对铀尾矿库日益严峻的环境问题及复杂多变的影响因素，显然这已经无法满足铀尾矿库退役治理的生态环保和可持续发展需求，如何对其稳定状态进行有效的分析和预测也是目前的研究瓶颈。

发明内容

本发明的目的是提供一种基于不确定理论的退役铀尾矿库环境稳定性分析及预测方法，解决了现有技术中存在的问题，实现了对铀尾矿库环境稳定性定性分析和定量计算的有机结合，对铀尾矿库环境指标由稳定状态转化为不稳定状态以及由不稳定状态转化为稳定状态所需时间进行了预测。

本发明所采用的技术方案是，一种基于不确定理论的退役铀尾矿库环境稳定性分析及预测方法，其特征在于，按照以下步骤进行：

步骤一、环境稳定性分析；

1)铀尾矿库环境稳定性指标体系的建立；

2)各指标稳定区间的计算；

3)铀尾矿库各坝段环境稳定率的分析；

步骤二、环境稳定性预测；

1)由稳定状态转化为不稳定状态所需时间预测；

2)由不稳定状态转化为稳定状态所需时间预测。

本发明的特征还在于，进一步的，所述步骤一中，铀尾矿库环境稳定性指标体系的建立，利用层次分析法求各环境指标的权重，具体按照以下步骤进行：

步骤1：构造判断矩阵；

判断矩阵中的元素值是各元素相对重要性判断的定量化指标，构造的判断矩阵C＝(C_ij)_n×n如表a所示；

表a

准则B_kC₁C₂ΛC_nC₁C₁₁C₁₂ΛC_1nC₂C₂₁C₂₂ΛC_2nΜΜΜΛΜC_nC_n1C_n2ΛC_nn

判断矩阵C具有如下性质：(1)C_ij>0；(2)C_ij＝1/C_ji(i≠j)；(3)C_ii＝1(i,j＝1,2,Λ,n)；判断矩阵中各元素的数值是对各因素相对重要程度作出判断，然后根据比率标度将判断定量化而获得的，采用1～9尺度法，判断矩阵标度及其含义见表b；

表b

序号重要性等级C_ij赋值1i,j两元素同等重要12i元素比j元素稍重要33i元素比j元素明显重要54i元素比j元素强烈重要75i元素比j元素极端重要96i元素比j元素稍不重要1/37i元素比j元素明显不重要1/58i元素比j元素强烈不重要1/79i元素比j元素极端不重要1/9

M_i为每一行元素的乘积，设a_ij为判断矩阵中第i行第j列的元素值，则

$>M_i=\underset{j=1}{\overset{n}{\Pi }}{a}_{ij}$>

其中，i＝1，2，…n；

步骤2：计算M_i的n次方根

$>{\bar{w}}_i=\sqrt[n]{M_i}$>

步骤3：对向量归一化处理，归一化公式为

$>w_i=\frac{{\bar{w}}_i}{\underset{j=1}{\overset{n}{\Sigma }}{\bar{w}}_j}$>

w_i为各指标的权重，则w＝[w₁,w₂,Λ,w_n]^T即为所求的特征向量；

步骤4：计算判断矩阵的最大特征根λ_max，设(Cw)_i表示向量Cw的第i个元素，则

$>{\lambda }_{max}=\underset{i=1}{\overset{n}{\Sigma }}\frac{{\left(Cw\right)}_i}{{\mathrm{nw}}_i}.$>

进一步的，对判断矩阵进行一致性检验步骤如下：

步骤1：计算一致性检验指标CI；

设λ_max判断矩阵最大特征值，则

$>CI=\frac{{\lambda }_{max}-n}{n-1}$>

步骤2：根据表1查找相应的平均随机一致性指标RI，其中n表示判断矩阵的阶数；

表1

n1234567891011RI000.580.901.121.241.321.411.451.491.51

步骤3：计算随机一致性比率CR，计算公式为

$>CR=\frac{CI}{RI}$>

当CR<0.10时，即认为判断矩阵具有满意的一致性，否则就需要调整判断矩阵的元素取值，使之具有满意的一致性。

进一步的，所述步骤一中，各指标稳定区间的计算按照以下步骤进行：

算法1.1：稳定区间的算法；

步骤1：输入X＝[x₁,x₂,Λ,x_k]，随机变量X为某指标的监测值，x_k≥0，x_k表示某指标第k个月份的监测值，α＝0.95；

步骤2：计算EY_k和DY_k，其中$>{\mathrm{EY}}_k=\frac{1}{k}\underset{i=1}{\overset{k}{\Sigma }}{x}_i,{\mathrm{DY}}_k=\frac{1}{k}\underset{i=1}{\overset{k}{\Sigma }}{(L-{\mathrm{EY}}_k)}^2;$>

其中：EY_k为环境指标随机变量序列Y_k的数学期望，DY_k为环境指标随机变量序列Y_k的方差，k表示环境指标监测的第k个月份，同时也是环境稳定性预测的起点时刻；

步骤3：计算ε_k，ε_k是计算指标稳定区间的中间变量，反应了不同监测时刻稳定区间偏离随机变量序列期望值的程度，其中

步骤4：计算a_k，b_k；a_k，b_k分别表示前k个月份指标所对应的稳定区间的下限和上限，其中a_k＝EY_k-ε_k，b_k＝EY_k+ε_k；

步骤5：输出[a_k,b_k]；计算$>\underset{i=1}{\overset{k}{I}}[a_k,b_k]=[a,b];$>

则计算结果[a,b]即为所求的前k个月的稳定区间。

进一步的，铀尾矿库各坝段环境稳定率的分析采用以下步骤进行：

将各坝段多个月份的环境监测数据代入稳定区间的算法，得到各指标的稳定区间，再结合环境稳定率的计算公式，得到各坝段环境稳定率ρ随时间的变化曲线；

其中，环境稳定率的定义及计算公式如下：

设x_i为各指标在某一时刻的观察值，w_i为各环境指标的权重，为各指标稳定区间的示性函数，为指标A_i稳定区间的下限值，为指标A_i稳定区间的上限值，环境稳定率ρ为某一时刻尾矿库环境的稳定性程度，且

$>\rho =w_1{\chi }_{A_1}\left(x_1\right)+w_2{\chi }_{A_2}\left(x_2\right)+\Lambda +w_m{\chi }_{A_m}\left(x_m\right)=\underset{i=1}{\overset{m}{\Sigma }}{w}_i·{\chi }_{A_i}\left(x_i\right)$>

其中，示性函数可表示为

进一步的，所述步骤二，由稳定状态转化为不稳定状态所需时间预测具体按照以下步骤进行：

假设t_k时刻某指标的初始值a₀为稳定区间的均值，即a₀＝EY_k，分别预测监测值按直线变化，指数变化及周期变化时突破稳定区间所需的时间；

由稳定状态转化为不稳定状态，新的稳定区间突破原有的稳定弹性区间需满足的关系式为：

$>\frac{\underset{i=1}{\overset{k+s}{\Sigma }}{y}_i}{k+s}+\sqrt{20\times \frac{\underset{i=1}{\overset{k+s}{\Sigma }}{(y_i-\frac{\underset{i=1}{\overset{k+s}{\Sigma }}{y}_i}{k+s})}^2}{k+s}}\ge \frac{\underset{i=1}{\overset{k}{\Sigma }}{y}_i}{k}+\sqrt{20\times \frac{\underset{i=1}{\overset{k}{\Sigma }}{(y_i-\frac{\underset{i=1}{\overset{k}{\Sigma }}{y}_i}{k})}^2}{k}}+\frac{b-a}{l}$>

其中：yi表示该指标各月份的监测值，当1≤i≤k时，yi等于前k个月份的监测值，当i>k时，y_i分别对应三种增长函数，且

设计算法1.2计算指标的弹性区间，算法1.3计算指标由稳定状态转化为不稳定状态所需的时间；

算法1.2：指标弹性区间的算法；

步骤1-1：输入ε_k，a，b，l；其中，a、b分别表示指标稳定区间的下限值和上限值，l为弹性区间的弹性系数；

步骤1-2：计算a^*，b^*，其中i＝1,2,Λ,k；a^*，b^*分别表示指标弹性区间的下限值和上限值；

步骤1-3：输出[a^*,b^*]；

其中[a^*,b^*]即为前k个月某指标稳定区间所对应的弹性区间；

算法1.3：稳定状态到不稳定状态所需时间的算法；

步骤2-1：输入X＝[x₁,x₂,Λ,x_k]以及前k个月的稳定区间[a_k,b_k]；X表示前k个月份某指标监测值构成的集合；

步骤2-2：计算x(n+1)＝f(t_n+1)，x＝f(t)为构造的浓度变化函数；更新X＝[x₁,x₂,Λ,x_k,x_k+1]；

步骤2-3：利用算法1.1计算出相应的稳定区间[a_k+1,b_k+1]，若则计算结束并输出第k+1个开始不稳定；反之，则输出第k+1个月稳定，并更新k＝k+1，回到步骤2-2循环计算，直到输出不稳定时间为止；

所述步骤二，由不稳定状态转化为稳定状态所需时间预测具体按照以下步骤进行：

假设某指标的初始值a₀为历史监测数据的最大值y_max，在此基础上预测指标监测值按指数衰减时进入稳定区间所需的时间；

由不稳定状态转化为稳定状态，新的稳定区间重新进入原有的稳定弹性区间需满足的关系式为：

$>\frac{\underset{i=1}{\overset{k+s}{\Sigma }}{y}_i}{k+s}+\sqrt{20\times \frac{\underset{i=1}{\overset{k+s}{\Sigma }}{(y_i-\frac{\underset{i=1}{\overset{k+s}{\Sigma }}{y}_i}{k+s})}^2}{k+s}}\le \frac{\underset{i=1}{\overset{k}{\Sigma }}{y}_i}{k}+\sqrt{20\times \frac{\underset{i=1}{\overset{k}{\Sigma }}{(y_i-\frac{\underset{i=1}{\overset{k}{\Sigma }}{y}_i}{k})}^2}{k}}+\frac{b-a}{l}$>

其中，y_i表示该指标各月份的观测值；当1≤i≤k时，y_i等于前k个月的监测值；当i>k时，y_i＝EY_k+(a₀-EY_k)·e^-mi，其中m≥0；i＝0，1，2，…s；e表示自然对数；a₀表示某指标的初始监测值；m为预测函数的衰减指数。

本发明的有益效果是，将不确定理论运用到铀尾矿库的环境稳定性分析和预测，一方面实现了对铀尾矿库环境稳定性的定量分析，得到了各个时刻尾矿库的环境稳定性动态变化规律，克服传统评价分析方法的缺陷，使铀尾矿库的稳定性计算分析更加合理；另一方面，通过不确定性理论对铀尾矿库的环境稳定性进行了预测，并充分考虑了从稳定状态到不稳定状态和从不稳定的状态到稳定状态的多种变化情况，并编程得出了环境指标在不同情况下突破稳定状态和重新回到稳定状态所需的时间，从而实现了对铀尾矿库环境稳定性的预测，更加科学和准确的分析了尾矿库的稳定性。

附图说明

图1是铀尾矿库环境稳定性指标体系图。

图2是跃进坝环境稳定率变化规律图。

图3是南坡坝环境稳定率变化规律图。

图4是松林坝环境稳定率变化规律图。

图5是西眉坝环境稳定率变化规律图。

图6是直线增长型突破稳定区间所需的时间图。

图7是指数增长型突破稳定区间所需的时间图。

图8是周期变化型突破稳定区间所需的时间图。

图9是指数衰减型进入稳定区间所需的时间图。

具体实施方式

下面结合附图和具体实施方式对本发明进行详细说明。

一种基于不确定理论的退役铀尾矿库环境稳定性分析及预测方法，包括环境稳定性分析与环境稳定性预测。

步骤一，环境稳定性分析。

1)铀尾矿库环境稳定性指标体系的建立；

铀尾矿库的环境污染物主要通过水的流动和大气的扩散迁移，影响周边的水塘、湘江水系以及土壤和大气，根据主成分分析和相关方法得出了影响环境的主要指标，并结合铀尾矿的实际情况，从尾矿库渗水、尾矿库大气环境、放射性污染三个方面选取了十二个环境影响指标，如图1所示，经过前期的退役治理大气环境指标和放射性污染指标已经达到了国家限定标准并基本趋于稳定，因此铀尾矿库渗水的指标是影响环境稳定性的主要因素。

由于退役铀尾矿库稳定化评价数据有限，无法采用客观的数理统计的方法分析各个因素指标权重的分布规律，只能根据人为主观判断来确定各个指标的权重，而层次分析法正(AHP)好克服了这些不足，可以适合各种缺乏数据支持、指标体系复杂的情况，并具有定性和定量分析相结合的优点，能把人为的评断结果通过数量的形式表达出来并进行科学处理分析，能够全面地反映问题，因此决定采用层次分析法来进行权重赋值。

利用层次分析法求各环境指标的权重，其具体步骤如下：

步骤1：构造判断矩阵。(判断矩阵的定义：把系统中每一层次各因素的相对重要性的判断结果用数值表示出来，写成矩阵形式就是判断矩阵)。判断矩阵中的元素值是各元素相对重要性判断的定量化指标，构造的判断矩阵C＝(C_ij)_n×n如表a所示。

表a判断矩阵一般形式

准则B_kC₁C₂ΛC_nC₁C₁₁C₁₂ΛC_1nC₂C₂₁C₂₂ΛC_2nΜΜΜΛΜC_nC_n1C_n2ΛC_nn

判断矩阵C具有如下性质：(1)C_ij>0；(2)C_ij＝1/C_ji(i≠j)；(3)C_ii＝1(i,j＝1,2,Λ,n)。

判断矩阵中各元素的数值是通过专家组对各因素相对重要程度作出判断，然后根据一定的比率标度将判断定量化而获得的。一般采用1～9尺度法，判断矩阵标度及其含义见表b。

表b判断矩阵标度及其含义

序号重要性等级C_ij赋值1i,j两元素同等重要12i元素比j元素稍重要33i元素比j元素明显重要54i元素比j元素强烈重要75i元素比j元素极端重要96i元素比j元素稍不重要1/37i元素比j元素明显不重要1/58i元素比j元素强烈不重要1/79i元素比j元素极端不重要1/9

每一行元素的乘积M_i，设a_ij为判断矩阵中第i行第j列的元素值，则

$>M_i=\underset{j=1}{\overset{n}{\Pi }}{a}_{ij}$>

其中，i＝1，2，…n；

步骤2：计算M_i的n次方根

$>{\bar{w}}_i=\sqrt[n]{M_i}$>

步骤3：对向量归一化处理，归一化公式为

$>w_i=\frac{{\bar{w}}_i}{\underset{j=1}{\overset{n}{\Sigma }}{\bar{w}}_j}$>

w_i为各指标的权重，则w＝[w₁,w₂,Λ,w_n]^T即为所求的特征向量；

步骤4：计算判断矩阵的最大特征根λ_max，设(Cw)_i表示向量Cw的第i个元素，则

$>{\lambda }_{max}=\underset{i=1}{\overset{n}{\Sigma }}\frac{{\left(Cw\right)}_i}{{\mathrm{nw}}_i}$>

一致性检验：为了保证结论的合理性，需要对判断矩阵进行一致性检验。对判断矩阵进行一致性检验步骤如下：

步骤1：计算一致性检验指标CI；

设λ_max判断矩阵最大特征值，则

$>CI=\frac{{\lambda }_{max}-n}{n-1}$>

步骤2：根据表1(注：表1是现有的公开技术，且所有RI都是统一的)查找相应的平均随机一致性指标RI，其中n表示判断矩阵的阶数。

表1平均随机一致性指标

n1234567891011RI000.580.901.121.241.321.411.451.491.51

步骤3：计算随机一致性比率CR，计算公式为

$>CR=\frac{CI}{RI}$>

当CR<0.10时，即认为判断矩阵具有满意的一致性，否则就需要调整判断矩阵的元素取值(即需要重新考虑模型或重新构造那些一致性比率CR较大的成对比较矩阵)，使之具有满意的一致性。

经专家组(由该领域的科研人员和工作人员组成)对铀尾矿库的环境指标进行两两比较，最终构造的判断矩阵A如下所示：

$>\left(\begin{array}{cccccc}1& 3/4& 3/2& 3& 5/2& 2\\ 4/3& 1& 2& 4& 7/2& 8/3\\ 2/3& 1/2& 1& 2& 5/3& 4/3\\ 1/3& 1/4& 1/2& 1& 5/6& 2/3\\ 2/5& 1/3& 3/5& 6/5& 1& 4/5\\ 1/2& 3/8& 3/4& 3/2& 5/4& 1\end{array}\right)$>

根据以上步骤进行一致性检验：CR＝0.0044＜0.1满足要求。

最终得到各指标的权重向量：w＝(0.23520.31620.15680.07840.09580.1176)。

2)各指标稳定区间的计算；

根据稳定区间的定义和算法1.1，取指标的稳定概率为0.95，即α＝0.95，则可算出跃进坝、南陂坝、松林坝及西眉坝各指标的稳定区间，见表2所示。

算法1.1：稳定区间的算法

步骤1：输入X＝[x₁,x₂,Λ,x_k]，(注：随机变量X为某指标的监测值，同前文)x_k≥0(x_k表示某指标第k个月份的监测值)，α＝0.95；

步骤2：计算EY_k和DY_k，其中$>{\mathrm{EY}}_k=\frac{1}{k}\underset{i=1}{\overset{k}{\Sigma }}{x}_i,{\mathrm{DY}}_k=\frac{1}{k}\underset{i=1}{\overset{k}{\Sigma }}{(x_i-{\mathrm{EY}}_k)}^2;$>

其中：EY_k为环境指标随机变量序列Y_k的数学期望，DY_k为环境指标随机变量序列Y_k的方差，k表示环境指标监测的第k个月份，同时也是环境稳定性预测的起点时刻；

步骤3：计算ε_k(ε_k是计算指标稳定区间的中间变量，反应了不同监测时刻稳定区间偏离随机变量序列期望值的程度)，其中

步骤4：计算a_k，b_k(a_k，b_k分别表示前k个月份指标所对应的稳定区间的下限和上限)，其中a_k＝EY_k-ε_k，b_k＝EY_k+ε_k；

步骤5：输出[a_k,b_k]；计算$>\underset{i=1}{\overset{k}{I}}[a_k,b_k]=[a,b].$>

则计算结果[a,b]即为所求的前k个月的稳定区间。

表2各坝段主要环境指标的稳定区间

3)铀尾矿库各坝段环境稳定率的分析；

将2010-2012年跃进坝、南坡坝、松林坝、西眉坝36个月的环境监测数据代入算法1.1得到各指标的稳定区间，再结合定义2中环境稳定率的计算公式即可得到各坝段环境稳定率ρ随时间的变化曲线，计算结果如图2-图5所示。

由图2可知跃进坝总体环境稳定性情况较好，除了第12个月，第30个月和第33个月以外，其余月份的环境稳定率均为1，由此说明跃进坝大部分月份的主要环境指标监测值都处于稳定状态，环境稳定性较好。

由图3可知南陂坝总体环境稳定性情况较差，36个月中仅有13个月环境稳定率等于1，其余月份环境稳定率均小于1，且稳定率曲线波动较大，由此说明南坡坝大部分月份的主要环境指标监测值都超出了稳定状态，环境稳定性较差。

由图4可知松林坝总体环境稳定性情况较好，36个月中有29个月环境稳定率等于1，其余月份的环境稳定率小于1，稳定率曲线在某些月份存在小幅度波动，但总体情况较好，由此说明松林坝环境基本上处于稳定状态。

由图5可知西眉坝总体环境稳定性情况很差，36个月中仅有6个月环境稳定率等于1，其余月份环境稳定率均小于1，某些月份的稳定率甚至小于0.5，且稳定率曲线波动的幅度和范围都很大，说明跃进坝环境处于极不稳定状态。

步骤二、环境稳定性预测；

1)由稳定状态转化为不稳定状态所需时间预测；

由稳定状态变为不稳定状态主要是因为铀尾矿库未来一段时间的监测值累积贡献突破原来的稳定区间，假设t_k时刻某指标的初始值a₀为稳定区间的均值，即a₀＝EY_k，在此基础上分别预测监测值按直线变化，指数变化及周期变化时突破稳定区间所需的时间。实际中直线增长型对应稳定排放，指数增长对应加速排放，周期变化型对应有规律排放等，其增长函数y_i分别如图6-8所示。

由稳定状态转化为不稳定状态，新的稳定区间突破原有的稳定弹性区间需满足的关系式为：

$>\frac{\underset{i=1}{\overset{k+s}{\Sigma }}{y}_i}{k+s}+\sqrt{20\times \frac{\underset{i=1}{\overset{k+s}{\Sigma }}{(y_i-\frac{\underset{i=1}{\overset{k+s}{\Sigma }}{y}_i}{k+s})}^2}{k+s}}\ge \frac{\underset{i=1}{\overset{k}{\Sigma }}{y}_i}{k}+\sqrt{20\times \frac{\underset{i=1}{\overset{k}{\Sigma }}{(y_i-\frac{\underset{i=1}{\overset{k}{\Sigma }}{y}_i}{k})}^2}{k}}+\frac{b-a}{l}$>

其中：y_i表示该指标各月份的监测值，当1≤i≤k时，y_i等于前k个月份的监测值，当i>k时，y_i分别对应三种增长函数，且本测试实例中k＝36。

根据计算及分析过程，设计了算法1.2计算指标的弹性区间，算法1.3计算指标由稳定状态转化为不稳定状态所需的时间。

算法1.2指标弹性区间的算法；

步骤1：输入ε_k，a，b，l；(a、b分别表示指标稳定区间的下限值和上限值，l为弹性区间的弹性系数)

步骤2：计算a^*，b^*，其中i＝1,2,Λ,k；(a^*，b^*分别表示指标弹性区间的下限值和上限值)

步骤3：输出[a^*,b^*]。

其中[a^*,b^*]即为前k个月某指标稳定区间所对应的弹性区间。

算法1.3稳定状态到不稳定状态所需时间的算法；

步骤1：输入X＝[x₁,x₂,Λ,x_k]以及前k个月的稳定区间[a_k,b_k]；

(此处X与算法1.1中的X含义一致，均表示前k个月份某指标监测值构成的集合)

步骤2：计算x(n+1)＝f(t_n+1)，x＝f(t)为构造的浓度变化函数；更新X＝[x₁,x₂,Λ,x_k,x_k+1]；

步骤3：利用算法1.1计算出相应的稳定区间[a_k+1,b_k+1]，若则计算结束并输出第k+1个开始不稳定；反之，则输出第k+1个月稳定，并更新k＝k+1，回到第二步循环计算，直到输出不稳定时间为止。

现以2010-2012年跃进坝F离子为例，详细数据已由表3列出，根据表3中36个月的监测数据预测F离子监测值按图6-8变化时由稳定状态变为不稳定状态所需的时间，则初始时刻t_k＝t₃₆，算出[t₁,t₃₆]时间内F离子的稳定区间根据定义3中的条件确定F离子l的范围，因为稳定区间确定后，可知前36个月的平均区间接近区间[a₃₆,b₃₆]，所以选取l使得$>a^{*}=a_{36}-\frac{b_{36}-a_{36}}{20}$>或$>b^{*}=b_{36}-\frac{b_{36}-a_{36}}{20},$>利用算法1.2计算稳定区间所对应的弹性区间[a^*,b^*]。

表32010-2012年跃进坝F离子监测数据

由表3中F离子36次监测数据算出随机变量序列Y₃₆的数学期望EY₃₆和方差DY₃₆，以及a₀的初值分别为：

EY₃₆＝1.6295，DY₃₆＝1.0759，a₀＝1.6295.

假设s个月后F离子新的稳定区间恰好突破了原来36个月监测值的稳定区间，则有

$>\frac{\underset{i=1}{\overset{36+s}{\Sigma }}{y}_i}{36+s}+\sqrt{20\times \frac{\underset{i=1}{\overset{36+s}{\Sigma }}{(y_i-\frac{\underset{i=1}{\overset{36+s}{\Sigma }}{y}_i}{36+s})}^2}{36+s}}\ge \frac{\underset{i=1}{\overset{36}{\Sigma }}{y}_i}{36}+\sqrt{20\times \frac{\underset{i=1}{\overset{36}{\Sigma }}{(y_i-\frac{\underset{i=1}{\overset{36}{\Sigma }}{y}_i}{36})}^2}{36}}+\frac{b-a}{20}$>

其中，y_i表示F离子各月份的观测值，当1≤i≤36时，y_i等于前36个月的监测值，当i>36时，y_i可由增长函数求得，且

根据算法1.3及表3中的监测数据，利用C⁺⁺语言进行编程，以跃进坝为实例，对三种情况下F离子由稳定状态变为不稳定状态所需要的时间进行预测，详细结果见表4-6所示，同理可预测跃进坝其他指标按直线增长时突破稳定区间所需的时间，此处就不一一列出。

表4跃进坝F离子按直线增长时突破稳定区间的时间预测

m取值预测时间(月)m取值预测时间(月)≤0.011570.13～0.14120.02750.14～0.15110.03470.16～0.17100.04370.18～0.2090.05290.21～0.2380.06250.24～0.2870.07210.29～0.3560.08190.36～0.4750.09170.48～0.6740.10150.68～1.0930.11141.10～2.4320.1213≥2.441

表5跃进坝F离子按指数增长时突破稳定区间的时间预测

表6跃进坝F离子按周期变化时突破稳定区间的时间预测

m取值预测时间(月)m取值预测时间(月)0.01170.88～1.3080.02151.31～1.4370.03～0.0414＝1.4460.05～0.06131.45～1.5050.07～0.11121.51～1.7440.12～0.21111.75～2.4430.22～0.44102.45～4.8720.45～0.879＝4.881

2)由不稳定状态转化为稳定状态所需时间预测；

由不稳定状态变为稳定状态需要监测值不断累积减少才能进入原来的稳定区间，假设某指标的初始值a₀为历史监测数据的最大值y_max，在此基础上预测指标监测值按指数衰减时进入稳定区间所需的时间，实际中指数衰减型对应自然净化和人工干预等，其衰减函数y_i如图9所示。

由不稳定状态转化为稳定状态，新的稳定区间重新进入原有的稳定弹性区间需满足的关系式为：

$>\frac{\underset{i=1}{\overset{k+s}{\Sigma }}{y}_i}{k+s}+\sqrt{20\times \frac{\underset{i=1}{\overset{k+s}{\Sigma }}{(y_i-\frac{\underset{i=1}{\overset{k+s}{\Sigma }}{y}_i}{k+s})}^2}{k+s}}\le \frac{\underset{i=1}{\overset{k}{\Sigma }}{y}_i}{k}+\sqrt{20\times \frac{\underset{i=1}{\overset{k}{\Sigma }}{(y_i-\frac{\underset{i=1}{\overset{k}{\Sigma }}{y}_i}{k})}^2}{k}}+\frac{b-a}{l}$>

其中，y_i表示该指标各月份的观测值。当1≤i≤k时，y_i等于前k个月的监测值；当i>k时，y_i＝EY_k+(a₀-EY_k)·e^-mi，其中m≥0；i＝0，1，2，…s；(e表示自然对数，约等于2.72；a₀表示某指标的初始监测值；m为预测函数的衰减指数)，测试中k＝36，在实例分析中y_i＝EY_k+(y_max-EY_k)·e^-mi，且其中，1≤i≤k。y_max为某指标历史监测数据的最大值。

以2010-2012年西眉坝F离子为例，根据表7中36个月的监测数据预测F离子按指数衰减时由不稳定状态变为稳定状态所需的时间，则初始时刻t_k＝t₃₆，算出[t₁,t₃₆]时间内F离子的稳定区间根据定义确定l的范围，类似直线

增长型，利用算法1.2算出了F离子稳定性的弹性区间[a^*,b^*]。

表72010-2012年西眉坝F离子监测数据

由表7中F离子的36次监测数据算出随机变量序列Y₃₆的数学期望EY₃₆和方差DY₃₆，以及a₀的初值分别为

EY₃₆＝3.64，DY₃₆＝16.89，a₀＝25.25.

假设s个月后F离子新的稳定区间恰好进入原来的36个月监测值的稳定区间，根据前文分析可得到如下不等式

$>\frac{\underset{i=1}{\overset{36+s}{\Sigma }}{y}_i}{36+s}+\sqrt{20\times \frac{\underset{i=1}{\overset{36+s}{\Sigma }}{(y_i-\frac{\underset{i=1}{\overset{36+s}{\Sigma }}{y}_i}{36+s})}^2}{36+s}}\le \frac{\underset{i=1}{\overset{36}{\Sigma }}{y}_i}{36}+\sqrt{20\times \frac{\underset{i=1}{\overset{36}{\Sigma }}{(y_i-\frac{\underset{i=1}{\overset{36}{\Sigma }}{y}_i}{36})}^2}{36}}+\frac{b-a}{20}$>

其中，y_i表示F离子各月份的观测值。当1≤i≤36时，y_i等于前36个月份的监测值，当i>36时，y_i＝3.64+21.61·e^-mi，且

根据算法1.3及表7中的原始监测数据，利用C⁺⁺语言进行编程，以西眉坝为例，对F离子由不稳定状态变为稳定状态所需要的时间进行预测，详细结果见表8，同理可预测西眉坝其他指标按指数衰减时进入稳定区间所需的时间，此处就不一一列出。

表8西眉坝F离子按指数衰减时进入稳定区间的时间预测

m取值预测时间(月)m取值预测时间(月)0.15149＝0.80470.20117＝0.85460.2598＝0.90450.3085＝0.95440.3576＝1.00430.40701.05～1.10420.4564＝1.15410.50601.20～1.30400.55571.35～1.45390.60541.50～1.65380.65521.70～1.95370.70502.00～2.75360.7548＝2.8035

本发明的关键点是不确定理论对退役铀尾矿库环境稳定性的定量分析，以及退役铀尾矿库环境指标稳定与不稳定状态时间的预测。对退役铀尾矿库的环境稳定化进程预测主要有环境指标由稳定状态转化为不稳定状态和不稳定状态转化为稳定状态的两种预测理论方法，可以预测出稳定状态下的指标进入不稳定状态的时间以及不稳定状态下的指标进入稳定状态的时间，为尾矿库环境治理提供理论依据。

本发明利用概率分析等相关方法，得出了铀尾矿库主要环境污染指标的稳定区间，并计算了铀尾矿库在不同时间点的环境稳定率，以精确的数学语言定义了铀尾矿库环境稳定性的模糊性概念，实现了对铀尾矿库环境稳定性定性分析和定量计算的有机结合。对铀尾矿库环境指标由稳定状态转化为不稳定状态所需时间进行了预测，并考虑了环境指标监测值按直线增长、指数增长及周期变化时突破稳定区间所需的时间；同时，还对铀尾矿库环境指标由不稳定状态转化为稳定状态所需时间进行了预测，并考虑了环境指标监测值按指数衰减变化时进入稳定区间所需的时间，为退役铀尾矿库安全管理的及决策分析提供了理论依据。