CN106372299A - 一种基于逆高斯过程和em‑pf的管道剩余寿命预测方法 - Google Patents

Abstract

本发明公开了一种基于逆高斯过程和EM‑PF的管道剩余寿命预测方法，包括以下步骤：获取油气管道的腐蚀数据，再根据油气管道的腐蚀数据利用随机过程‑逆高斯过程描述油气管道的腐蚀退化过程，建立油气管道状态空间模型的状态方程及观测方程，再根据油气管道状态空间模型的状态方程及观测方程建立逆高斯‑状态腐蚀退化模型，然后根据逆高斯‑状态腐蚀退化模型预测油气管管道的剩余寿命，该方法能够准确、快速的预测油气管道的剩余寿命。

Description

一种基于逆高斯过程和EM-PF的管道剩余寿命预测方法

技术领域

本发明涉及一种剩余寿命预测方法，具体涉及一种基于逆高斯过程和EM-PF的管道剩余寿命预测方法。

背景技术

石油天然气管道作为油气储运中主要的集输设备，发挥着重要的作用。随着油气管道向海底不断延伸，影响其腐蚀老化的因素和环境条件日趋复杂，潜在安全隐患随之增大。油气管道的剩余寿命预测及维修是其完整性管理中的关键技术，可以降低发生风险的概率进而为优化维修方案提供依据。

国外J.L.Alamilla，A.Valor等研究基于随机过程的管道腐蚀退化模型；Gasperin等研究线性函数模型退化过程模型下的剩余寿命分布。Peng W建立逆高斯过程模型下的退化模型，进而分析其可靠性；X.Wang研究基于逆高斯过程的腐蚀退化模型；张增刚等量化腐蚀油气管道结构性能的衰减程度，将概率论引入其腐蚀评估领域；赵新伟等建立管道点蚀损伤沿壁厚方向的分布模型，提出点蚀弥散损伤管道剩余寿命的预测方法；Wang建立假设模型参数服从维纳分布的随机退化过程；樊红东等在文献基础上建立非线性函数退化模型，进而估计剩余寿命；任淑红等提出基于性能退化下的剩余寿命预测方法，提高热老化设备可靠性；司小胜等建立带测量误差的非线性退化过程模型进而估计剩余寿命；张星辉等研究基于混合高斯输出贝叶斯信念网络模型的设备退化状态识别与剩余使用寿命预测；陈亮基于伽马过程建立腐蚀退化模型，研究测量误差对估计剩余寿命的影响；袁庚采用策略改进算法对马尔科夫过程求解，分析管道维修措施及相应费用，得出最优管线维修策略；张新生等基于伽马过程建立考虑维修深度的剩余寿命预测模型，最终实现油气管道腐蚀系统定期维修经济最优化，然而上述技术中没有能够准确的预测油气管道剩余寿命的方法，因此油气管道的可靠性失效时间数据短时间难以得到。

发明内容

本发明的目的在于克服上述现有技术的缺点，提供了一种基于逆高斯过程和EM-PF的管道剩余寿命预测方法，该方法能够准确、快速的预测油气管道的剩余寿命。

为达到上述目的，本发明所述的基于逆高斯过程和EM-PF的管道剩余寿命预测方法包括以下步骤：

获取油气管道的腐蚀数据，再根据油气管道的腐蚀数据利用随机过程-逆高斯过程描述油气管道的腐蚀退化过程，建立油气管道状态空间模型的状态方程及观测方程，再根据油气管道状态空间模型的状态方程及观测方程建立逆高斯-状态腐蚀退化模型，然后根据逆高斯-状态腐蚀退化模型预测油气管管道的剩余寿命。

设油气管道管壁的腐蚀状态退化过程符合逆高斯过程，y_i为油气管道管壁t_i时刻的腐蚀管段观测变量，x_i为油气管道管壁t_i时刻的腐蚀状态，则建立的油气管道状态空间模型的状态方程及观测方程分别为：

x_i+1-x_i～IG(Λ(t_i+1)-Λ(t_i),λ) (4)

y_i＝H(x_i)+ε (5)。

设当油气管道管壁t_i时刻腐蚀状态x_i达到油气管道管壁的允许最大值X_d时，油气管道管壁的腐蚀管道穿孔泄露，则油气管道管壁的腐蚀状态首次达到管壁的充许最大值X_d所需的时间T_d为：

T_d＝inf{t|X(t)≥X_d,t＞0} (6)

设油气管道的管壁状态退化过程为逆高斯过程，形状参数Λ(t)为单调增长函数，即Λ(t)＝ηt，设y＝c·x_i+ε，则在观测序列下，逆高斯-状态腐蚀退化模型为：

x_i-x_i-1～IG(Λ(t_i)-Λ(t_i-1),λ) (7)

y_i＝c·x_i+ε (8)

其中，ε服从正态分布，y～N(cx,σ)，η,λ,c,σ采用EM-PF算法估计得到。

令N_l为粒子个数，T_d为腐蚀油气管道泄漏管壁破裂的时刻，ω_i为第i个粒子的权重，X_d为油气管道壁厚的最大值，为当前时刻的观测序列，为当前时刻管壁观测值序列，t_i为腐蚀油气管道第i次时的管壁状态检测时刻，x_i为t_i时刻的管壁状态值，为第n次检测管道所获得的状态序列，为第1次检测管道所获得的状态序列，y_i为t_i时刻的退化量观测值，为第n次检测管道所获得的观测序列，得到t_i时刻的剩余寿命τ_i分布函数以及当腐蚀信息值时，腐蚀油气管道的平均剩余寿命为：

$F\left({\mathrm{\τ}}_i\right|y_0^c)=P(t_i+{\mathrm{\τ}}_i\≥T_d|y_0^i)=\∫p({\mathrm{\τ}}_i+t_i\≥T_d|x_i)\·p\left(x_i\right|y_0^i)dx---\left(10\right)$

$p({\mathrm{\τ}}_i+t_i\≥T_d|x_i)=p\left(x\right({\mathrm{\τ}}_i+t_i)\≥X_d|x_i)=\frac{IG\left(\mathrm{\Λ}\right(t_d)-\mathrm{\Λ}\left(t_i\right),\mathrm{\λ})}{IG\left(\mathrm{\Λ}\right(t_i))}---\left(11\right)$

$p\left(x_i\right|y_0^i)=f\left(x_0^i\right|y_0^i)\≈\underset{l=1}{\overset{N_l}{\Σ}}{\mathrm{\ω}}_i^l\mathrm{\δ}(x_i-x_i^l)---\left(12\right)$

预测tk时刻的油气管道的状态概率密度函数为：

$p\left(x_k\right|y_0^i)=\∫p\left(x_k\right|x_i)p\left(x_i\right|y_0^i){\mathrm{dx}}_i=\underset{l=1}{\overset{N_l}{\Σ}}p\left(x_k\right|x_i^l)\·{\mathrm{\ω}}_i^l---\left(13\right)$

通过表示式(12)，则油气管道的剩余寿命分布函数及其概率密度函数分别为：

$\begin{array}{c}F\left({\mathrm{\τ}}_i\right|y_0^c)=p({\mathrm{\τ}}_i+t_i\≥T_d|y_0^i)\\ =\∫p({\mathrm{\τ}}_i+t_i\≥T_d|x_i)\·p\left(x_i\right|y_0^i)dx\\ \≈\underset{l=1}{\overset{N_l}{\Σ}}p({\mathrm{\τ}}_i+t_i\≥T_d|x_i^l){\mathrm{\ω}}_i^l\\ =\underset{l=1}{\overset{N_l}{\Σ}}\frac{IG\left(\mathrm{\Λ}\right(t_d)-\mathrm{\Λ}\left(t_i\right),\mathrm{\λ})}{IG\left(\mathrm{\Λ}\right(t_i))}\·{\mathrm{\ω}}_i^l\end{array}---\left(14\right)$

$f\left({\mathrm{\τ}}_i\right|y_0^c)=\frac{dF\left({\mathrm{\τ}}_i\right|y_0^c)}{{\mathrm{d\τ}}_i}---\left(15\right)$

则当腐蚀信息值为时，油气管道的平均剩余寿命为：

$\stackrel{\‾}{T}=\underset{0}{\overset{\mathrm{\∞}}{\∫}}(1-F({\mathrm{\τ}}_i|y_0^c\left)\right)dt---\left(16\right)$

本发明具有以下有益效果：

本发明所述的基于逆高斯过程和EM-PF的管道剩余寿命预测方法将管道腐蚀退化状态与腐蚀数据结合起来建立退化模型进行预测油气管道的剩余寿命，在具体操作时，以随机过程理论为基础，根据油气管道的腐蚀数据利用随机过程-逆高斯过程描述油气管道的腐蚀退化过程，同时利用油气管道状态空间模型的状态方程及观测方程建立逆高斯-状态腐蚀退化模型，从而准确、快速的预测出油气管管道的剩余寿命，在实际运用时，可以根据预测出来结果获取油气管道腐蚀管段的最佳维修时间和最小维修费用，为设备维护及检修提供了有力的技术支持。

附图说明

图1(a)为本发明中参数η迭代轨迹图；

图1(b)为本发明中参数λ迭代轨迹图；

图2为本发明中参数η迭代历史图；

图3为本发明中参数λ迭代历史图；

图4为本发明中不同时刻的剩余寿命概率密度分布图；

图5为本发明中计算不同检测时刻下的维修、更换时间分布图。

具体实施方式

下面结合附图对本发明做进一步详细描述：

参考图1，本发明所述的基于逆高斯过程和EM-PF的管道剩余寿命预测方法包括以下步骤：获取油气管道的腐蚀数据，再根据油气管道的腐蚀数据利用随机过程-逆高斯过程描述油气管道的腐蚀退化过程，建立油气管道状态空间模型的状态方程及观测方程，再根据油气管道状态空间模型的状态方程及观测方程建立逆高斯-状态腐蚀退化模型，然后根据逆高斯-状态腐蚀退化模型预测油气管管道的剩余寿命。

设油气管道管壁的腐蚀状态退化过程符合逆高斯过程，y_i为油气管道管壁t_i时刻的腐蚀管段观测变量，x_i为油气管道管壁t_i时刻的腐蚀状态，则建立的油气管道状态空间模型的状态方程及观测方程分别为：

x_i+1-x_i～IG(Λ(t_i+1)-Λ(t_i),λ) (4)

y_i＝H(x_i)+ε (5)。

设当油气管道管壁t_i时刻腐蚀状态x_i达到油气管道管壁的允许最大值X_d时，油气管道管壁的腐蚀管道穿孔泄露，则油气管道管壁的腐蚀状态首次达到管壁的充许最大值X_d所需的时间T_d为：

T_d＝inf{t|X(t)≥X_d,t＞0} (6)

设油气管道的管壁状态退化过程为逆高斯过程，形状参数Λ(t)为单调增长函数，即Λ(t)＝ηt，设y＝c·x_i+ε，则在观测序列下，逆高斯-状态腐蚀退化模型为：

x_i-x_i-1～IG(Λ(t_i)-Λ(t_i-1),λ) (7)

y_i＝c·x_i+ε (8)

其中，ε服从正态分布，y～N(cx,σ)，η,λ,c,σ采用EM-PF算法估计得到。

令N_l为粒子个数，T_d为腐蚀油气管道泄漏管壁破裂的时刻，ω_i为第i个粒子的权重，X_d为油气管道壁厚的最大值，为当前时刻的观测序列，为当前时刻管壁观测值序列，t_i为腐蚀油气管道第i次时的管壁状态检测时刻，x_i为t_i时刻的管壁状态值，为第n次检测管道所获得的状态序列，为第1次检测管道所获得的状态序列，y_i为t_i时刻的退化量观测值，为第n次检测管道所获得的观测序列，得到t_i时刻的剩余寿命τ_i分布函数以及当腐蚀信息值时，腐蚀油气管道的平均剩余寿命为：

$F\left({\mathrm{\τ}}_i\right|y_0^c)=P(t_i+{\mathrm{\τ}}_i\≥T_d|y_0^i)=\∫p({\mathrm{\τ}}_i+t_i\≥T_d|x_i)\·p\left(x_i\right|y_0^i)dx---\left(10\right)$

$p({\mathrm{\τ}}_i+t_i\≥T_d|x_i)=p\left(x\right({\mathrm{\τ}}_i+t_i)\≥X_d|x_i)=\frac{IG\left(\mathrm{\Λ}\right(t_d)-\mathrm{\Λ}\left(t_i\right),\mathrm{\λ})}{IG\left(\mathrm{\Λ}\right(t_i))}---\left(11\right)$

$p\left(x_i\right|y_0^i)=f\left(x_0^i\right|y_0^i)\≈\underset{l=1}{\overset{N_l}{\Σ}}{\mathrm{\ω}}_i^l\mathrm{\δ}(x_i-x_i^l)---\left(12\right)$

预测t_k时刻的油气管道的状态概率密度函数为：

$p\left(x_k\right|y_0^i)=\∫p\left(x_k\right|x_i)p\left(x_i\right|y_0^i){\mathrm{dx}}_i=\underset{l=1}{\overset{N_l}{\Σ}}p\left(x_k\right|x_i^l)\·{\mathrm{\ω}}_i^l---\left(13\right)$

通过表示式(12)，则油气管道的剩余寿命分布函数及其概率密度函数分别为：

$\begin{array}{c}F\left({\mathrm{\τ}}_i\right|y_0^c)=p({\mathrm{\τ}}_i+t_i\≥T_d|y_0^i)\\ =\∫p({\mathrm{\τ}}_i+t_i\≥T_d|x_i)\·p\left(x_i\right|y_0^i)dx\\ \≈\underset{l=1}{\overset{N_l}{\Σ}}p({\mathrm{\τ}}_i+t_i\≥T_d|x_i^l){\mathrm{\ω}}_i^l\\ =\underset{l=1}{\overset{N_l}{\Σ}}\frac{IG\left(\mathrm{\Λ}\right(t_d)-\mathrm{\Λ}\left(t_i\right),\mathrm{\λ})}{IG\left(\mathrm{\Λ}\right(t_i))}\·{\mathrm{\ω}}_i^l\end{array}---\left(14\right)$

$f\left({\mathrm{\τ}}_i\right|y_0^c)=\frac{dF\left({\mathrm{\τ}}_i\right|y_0^c)}{{\mathrm{d\τ}}_i}---\left(15\right)$

则当腐蚀信息值为时，油气管道的平均剩余寿命为：

$\stackrel{\‾}{T}=\underset{0}{\overset{\mathrm{\∞}}{\∫}}(1-F({\mathrm{\τ}}_i|y_0^c\left)\right)dt---\left(16\right)$

仿真实验

某油气管道为APILX52级钢制管道，管道内压为10MPa，外径为340.8mm，管壁厚度为9.9mm，最小极限应力为359MPa，流变系数为1.1，采集数据时，针对油气管道腐蚀状态在线监控系统对腐蚀状态所检测的腐蚀状态作为间接状态观测数据y，管壁磨损量为状态数据x，故对其检测的原始腐蚀，经过分析数据见表1。

表1

首先进行参数的估计，具体的，根据EM-PF算法，选取正态分布作为η和λ的先验分布，其中，η～N(0.001,10.0)，λ～N(0.01,5.0)，结合表1所述的数据估计模型参数进行参数估计，获得参数的吉布斯迭代轨迹(trace)和迭代历史(history)，其中，迭代轨迹如图1(a)、图1(b)、图2及图3所示。

根据EM-PF算法，对腐蚀油气管道IG-SSM退化过程模型中的参数进行估计，结果见表2。

表2

根据表3中的参数估计值，确定IG-SSM退化模型，进而计算求得油气管道的剩余寿命概率密度f(τ_i|y_0:i)，可以通过软件实现不同检测时刻下的剩余寿命概率密度函数曲线图以及相对应预测剩余寿命，如图4所示。

表3

由表3可知，本发明对于预测腐蚀油气管道的剩余寿命方面具有较好的实用性。

在实际应用时，假设C_p＝1200元，C_f＝2000元，在得到腐蚀油气管道的剩余寿命概率密度函数f(τ_i|y_0:i)之后，将其代入油气管道的维修决策模型，各维修时刻的单位时间费用以及相对应的检测时费用最低的更换时间计算结果如图5所示，由图5可知，随着油气管道更换时间的增加，其单位时间费用先减少后增加，在极值点处可以取得合理维修时间。

通过编程计算求解油气管道的维修决策模型，可以得出首次满足T_r-t_i≤Δt条件时，T_r＝30a，即当管道服役30a时对其进行更换，C(T_r)将达到最低，约为376元。

本发明将油气管道的腐蚀数据以及腐蚀状态结合起来进行腐蚀油气管道退化模型的建立；在进行退化模型参数估计的过程中，基于超参数的估计问题，当进行极大似然估计时，主要考虑到以下几点：如果数据缺失满足、观测数据不完全或似然函数表达式较为复杂，主要采用EM算法进行参数估计，同时，PF算法是在序贯重要性采样和贝叶斯理论基础上，采用蒙特卡罗模拟方法和递推贝叶斯估计的一种统计滤波方法，其能够处理非线性、逆高斯问题，具体操作为：结合整体状态向量的经验分布，首先在状态空间形成一组粒子，然后根据所观测数据调整粒子的位置和重要程度，根据所进行调整后的粒子信息对最初的先验分布进行修正，该方法的本质在于由粒子位置及重要程度形成概率分布，根据算法递推进行离散随机测度的更新，本发明采用兼有各自优点的EM-PF算法，可以更加有效地对模型中的参数进行求解。

本发明充分考虑油气管道的腐蚀状态对油气管道腐蚀的影响，尽量能更加符合油气管道的腐蚀状况，以此建立符合油气管道实际腐蚀状况的退化模型，同时结合EM-PF算法进行模型参数的估计，进而预测腐蚀油气管道的剩余寿命，以便确定合理的维修时间，极大提高了腐蚀油气管道剩余寿命预测的效率和准确度，进而实现优化油气管道维修决策的目的，同时仿真结果表明本发明在腐蚀油气管道剩余寿命预测及维修决策优化方面比文献中给出的检测方法更为准确、有效。

以上实施方式仅用于说明本发明，而并非对本发明的限制。尽管为说明目的公开了本发明的最佳实施例和附图，但是本领域的技术人员可以理解；在不脱离本发明及所附的权利要求的精神和范围内，各种替换、变化、修改都是可能的。因此，所有等同的技术方案也属于本发明的范畴，本发明的专利保护范围应由权利要求限定，而不应局限于最佳实施例和附图所公开的内容。

Claims (4)

1.一种基于逆高斯过程和EM-PF的管道剩余寿命预测方法，其特征在于，包括以下步骤：

获取油气管道的腐蚀数据，再根据油气管道的腐蚀数据利用随机过程-逆高斯过程描述油气管道的腐蚀退化过程，建立油气管道状态空间模型的状态方程及观测方程，再根据油气管道状态空间模型的状态方程及观测方程建立逆高斯-状态腐蚀退化模型，然后根据逆高斯-状态腐蚀退化模型预测油气管管道的剩余寿命。

2.根据权利要求1所述的基于逆高斯过程和EM-PF的管道剩余寿命预测方法，其特征在于，

设油气管道管壁的腐蚀状态退化过程符合逆高斯过程，y_i为油气管道管壁t_i时刻的腐蚀管段观测变量，x_i为油气管道管壁t_i时刻的腐蚀状态，则建立的油气管道状态空间模型的状态方程及观测方程分别为：

x_i+1-x_i～IG(Λ(t_i+1)-Λ(t_i),λ) (4)

y_i＝H(x_i)+ε (5)。

3.根据权利要求1所述的基于逆高斯过程和EM-PF的管道剩余寿命预测方法，其特征在于，

设当油气管道管壁t_i时刻腐蚀状态x_i达到油气管道管壁的允许最大值X_d时，油气管道管壁的腐蚀管道穿孔泄露，则油气管道管壁的腐蚀状态首次达到管壁的充许最大值X_d所需的时间T_d为：

T_d＝inf{t|X(t)≥X_d,t>0} (6)

设油气管道的管壁状态退化过程为逆高斯过程，形状参数Λ(t)为单调增长函数，即Λ(t)＝ηt，设y＝c·x_i+ε，则在观测序列下，逆高斯-状态腐蚀退化模型为：

x_i-x_i-1～IG(Λ(t_i)-Λ(t_i-1),λ) (7)

y_i＝c·x_i+ε (8)

其中，ε服从正态分布，y～N(cx,σ)，η,λ,c,σ采用EM-PF算法估计得到。

4.根据权利要求1所述的基于逆高斯过程和EM-PF的管道剩余寿命预测方法，其特征在于，

令N_l为粒子个数，T_d为腐蚀油气管道泄漏管壁破裂的时刻，ω_i为第i个粒子的权重，X_d为油气管道壁厚的最大值，为当前时刻的观测序列，为当前时刻管壁观测值序列，t_i为腐蚀油气管道第i次时的管壁状态检测时刻，x_i为t_i时刻的管壁状态值，为第n次检测管道所获得的状态序列，为第1次检测管道所获得的状态序列，y_i为t_i时刻的退化量观测值，为第n次检测管道所获得的观测序列，得到t_i时刻的剩余寿命τ_i分布函数以及当腐蚀信息值时，腐蚀油气管道的平均剩余寿命为：

$F\left({\mathrm{\τ}}_i\right|y_0^c)=P(t_i+{\mathrm{\τ}}_i\≥T_d|y_0^i)=\∫p({\mathrm{\τ}}_i+t_i\≥T_d|x_i)\·p\left(x_i\right|y_0^i)dx---\left(10\right)$

$p({\mathrm{\τ}}_i+t_i\≥T_d|x_i)=p\left(x\right({\mathrm{\τ}}_i+t_i)\≥X_d|x_i)=\frac{IG\left(\mathrm{\Λ}\right(t_d)-\mathrm{\Λ}(t_i),\mathrm{\λ})}{IG\left(\mathrm{\Λ}\right(t_i\left)\right)}---\left(11\right)$

$p\left(x_i\right|y_0^i)=f\left(x_0^i\right|y_0^i)\≈\underset{l=1}{\overset{N_l}{\Σ}}{\mathrm{\ω}}_i^l\mathrm{\δ}(x_i-x_i^l)---\left(12\right)$

预测t_k时刻的油气管道的状态概率密度函数为：

$p\left(x_k\right|y_0^i)=\∫p\left(x_k\right|x_i)p\left(x_i\right|y_0^i){\mathrm{dx}}_i=\underset{l=1}{\overset{N_l}{\mathrm{\Σ}}}p\left(x_k\right|x_i^l)\·{\mathrm{\ω}}_i^l---\left(13\right)$

通过表示式(12)，则油气管道的剩余寿命分布函数及其概率密度函数分别为：

$\begin{array}{c}F\left({\mathrm{\τ}}_i|y_0^c\right)=p\left({\mathrm{\τ}}_i+t_i\≥T_d|y_0^i\right)\\ =\∫p\left({\mathrm{\τ}}_i+t_i\≥T_d|x_i\right)\·p\left(x_i|y_0^i\right)dx\\ \≈\underset{l=1}{\overset{N_l}{\mathrm{\Σ}}}p\left({\mathrm{\τ}}_i+t_i\≥T_d|x_i^l\right){\mathrm{\ω}}_i^l\\ \≈\underset{l=1}{\overset{N_l}{\mathrm{\Σ}}}\frac{IG\left(\mathrm{\Λ}\left(t_d\right)-\mathrm{\Λ}\left(t_i\right),\mathrm{\λ}\right)}{IG\left(\mathrm{\Λ}\left(t_i\right)\right)}\·{\mathrm{\ω}}_i^l\end{array}---\left(14\right)$

$f\left({\mathrm{\τ}}_i\right|y_0^c)=\frac{dF\left({\mathrm{\τ}}_i\right|y_0^c)}{{\mathrm{d\τ}}_i}---\left(15\right)$

则当腐蚀信息值为时，油气管道的平均剩余寿命为：

$\stackrel{\‾}{T}=\underset{0}{\overset{\mathrm{\∞}}{\∫}}(1-F({\mathrm{\τ}}_i|y_0^c\left)\right)dt---\left(16\right).$