CN104408322A - 综合多源故障概率似然信度的旋转机械设备故障诊断方法 - Google Patents

Abstract

本发明涉及一种综合多源故障似然信度的旋转机械设备故障诊断方法。本发明能够对不同故障模式下获取的故障特征监测数据进行区间划分，构造故障特征参数与故障类型的映射关系表；从该表中统计出各故障模式发生的似然信度向量；结合训练样本集合构造目标函数，并优化该函数获得可靠性因子，用该可靠性因子对似然信度向量进行修正，生成诊断证据；在线获取多种故障特征的取值后，分别计算它们激活的诊断证据，再将这些被激活的诊断证据综合，根据综合结果决策得到该在线故障特征所对应的故障类型；本发明对海量故障特征数据进行处理时，未对数据的变化特性、格式、精度等做出任何限制和约束，便于工程实现。

Description

综合多源故障概率似然信度的旋转机械设备故障诊断方法

技术领域

本发明涉及一种综合多源故障概率似然信度的旋转机械设备故障诊断方法，属于旋转机械设备状态监测与故障诊断技术领域。

背景技术

在线监测与故障诊断技术的使用，能够使维护工程师及时监测旋转机械设备的工作异常，发现危及设备安全运行的故障，并对设备的定期或事情维护提供必要的历史运行资料，这对于降低设备维护成本并提升其工作效率起到了至关重要的作用。对于现代化的大型复杂设备，如大型旋转机械设备，在其各个关键部位可以布置密集的振动传感器采集设备的振动信号，从这些振动信号中提取的故障特征信息可以反映设备的各种故障。此时面临两个问题，一是需要找到一种在工程上简单易行的方式，对海量的故障特征监测数据进行分析处理，使得其提供的诊断证据客观可信；二是如何实现故障特征信息的综合。通常利用单个传感器提供的特征信息不能完全反应设备的故障，需要多个传感器提供的信息进行融合来提升诊断的精准度。

发明内容

本发明的目的是提出一种综合多源故障概率似然信度的旋转机械设备故障诊断方法，从海量故障特征数据中统计各故障发生的似然信度，给出可靠性因子对该似然信度进行修正得到可靠的诊断证据，并给出诊断证据综合方法对各个故障特征提供的诊断证据进行综合，利用综合后的证据做出故障决策，该方法对海量故障特征数据进行处理时，未对数据的变化特性、格式、精度等做出任何限制和约束，便于工程实现，且通过综合后证据可以做出精确的诊断。

本发明提出的综合多源故障概率似然信度的旋转机械设备故障诊断方法，包括以下各步骤：

(1)设定旋转机械设备的故障集合Θ＝{F₁,...F_i,...F_N}，F_i代表故障集合Θ中的第i个故障，i＝1,2,...,N，N为设备含有故障的个数。

(2)设x为能够反映故障集合Θ中每个故障F_i的故障特征参数，该特征参数由某传感器信息源提供，建立故障特征参数x关于每个故障F_i的特征参数取值变化区间其中，和分别代表区间的左右端点，具体获取的步骤如下：

(2-1)当故障集合Θ中的故障F_i发生时，可获取故障特征参数x的δ个测量样本，记这些测量样本构成的集合为一般500≥δ≥300。

(2-2)对于步骤(2-1)中获取的测量样本集合分别求取这些样本的最小值和最大值，令它们分别为和并构成区间

(3)按照步骤(2)可以获取故障集合Θ中N种故障的特征参数x的测量样本变化区间将这N个区间的2N个左右端点按照从大到小的顺序排序，组成含有2N个点的序列S＝{s₁,s₂,…,s_2N}，其中按照s_j的排序，生成故障特征参数x关于故障集合Θ的2N+1个样本变化区间 $I_1^x=[0,s_1),I_2^x=[s_1,s_2),I_3^x=[s_2,s_3),...,$ $I_{2N}^x={[x}_{2N-1},s_{2N}),I_{2N+1}^x=[s_{2N+1},\∞).$

(4)在故障集合Θ中的N种故障分别发生时，通过步骤(2-1)可获得每种故障的测量样本集合 $V_1^x=\{v_{\mathrm{1,1}}^x,v_{\mathrm{1,2}}^x,...,v_{1,\mathrm{\δ}}^x\},V_2^x=\{v_{\mathrm{2,1}}^x,v_{\mathrm{2,2}}^x,...,v_{2,\mathrm{\δ}}^x\},...,$ $V_N^x=\{v_{N,1}^x,v_{N,2}^x,...,v_{N,\mathrm{\δ}}^x\},$ 共计可获得Nδ个样本，用它们可构造故障特征x和故障F₁,F₂,…,F_N之间的映射关系表，如下表1所示，其中l＝1,2，…,2N+1，为步骤(3)中构造的特征参数x测量样本变化区间的个数，a_il表示故障F_i的故障特征参数x的样本落入区间[s_l-1,s_l)中的个数，并有η_l为落入区间[s_l-1,s_l)中的样本个数总和，并有

表1 故障特征参数x与故障类型映射关系表

(5)根据步骤(4)中获取的映射关系表，可获得当故障F_i发生时，故障特征参数x的取值落入区间[s_l-1,s_l)的似然函数为：

$p\left(I_l^x\right|F_i)=\frac{a_{\mathrm{il}}}{\mathrm{\δ}}---\left(1\right)$

将式(1)中的似然函数进行归一化，获得当故障特征x的取值落入区间时，故障F_i发生的似然信度为：

$b_x\left(F_i\right)=\frac{p\left(I_l^x\right|F_i)}{{\mathrm{\Σ}}_{i=1}^{N}p\left(I_l^x\right|F_i)}---\left(2\right)$

并有则可以定义此时获取的似然信度向量为：

B_x＝[b_x(F₁),b_x(F₂),…,b_x(F_N),b_x(Θ)]    (3)

其中b_x(Θ)＝0，表示为对故障集合Θ信度赋值为0；若任何故障的特征参数x的样本都未落入区间中，亦即η_l＝0，则此时b_x(F_i)＝0，b_x(Θ)＝1，似然信度向量的取值为：

B_x＝[0,0,…,0,1]            (4)

给定一个故障特征参数x的取值，其必然落入中的某一个区间，此时该区间所对应的似然信度向量被激活，则可以被激活的似然信度向量的个数为2N+1，它们分别由式(3)和式(4)给出。

(6)定义可靠性因子α_x用其对似然信度向量B_x中的元素进行折扣，生成关于特征参数x的诊断证据

$m_x\left(A\right)=\left\{\begin{array}{cc}(1-{\mathrm{\α}}_x)b_x\left(A\right)& \∀A\⋐\mathrm{\Θ},A\≠\mathrm{\Θ}\\ {\mathrm{\α}}_x+(1-{\mathrm{\α}}_x)b_x\left(A\right)& A=\mathrm{\Θ}\end{array}\right.---\left(5\right)$

求取可靠性因子α_x的步骤如下：

(6-1)在故障集合Θ中的N种故障分别发生时，再一次利用步骤(2-1)获得每种故障在故障特征参数x下的训练样本集合 $T_1^x=\{t_{\mathrm{1,1}}^x,t_{\mathrm{1,2}}^x,...,t_{1,\mathrm{\ω}}^x\},T_2^x=\{t_{\mathrm{2,1}}^x,t_{\mathrm{2,2}}^x,...,t_{2,\mathrm{\ω}}^x\},...,T_N^x=\{t_{N,1}^x,t_{M,2}^x,...,t_{N,\mathrm{\ω}}^x\},$ 一般200≥ω≥50，共计可获得Nω个样本，它们组成的集合定义为 $T^x=T_1^x\∪T_2^x\∪...\∪T_N^x.$

(6-2)构造目标函数

它是关于可靠性因子α_x的函数，其中，b_x(F_i)，b_x(Θ)为根据落入的特征参数x的样本变化区间，由步骤(5)确定的故障类型F_i和故障集合Θ的似然信度赋值，当则否则

(6-3)设定α_x的取值范围为0到1之间，利用数学计算软件Matlab中的fmincom非线性优化函数，最小化Obj(α_x)的取值，此时对应获得的Obj函数的输入值即为寻找的可靠性因子α_x的取值。

(7)除了故障特征参数x，若存在另一个故障特征参数y也能反映Θ中的各个故障，则重复以上步骤(2)至(6)，可以获得关于y的诊断证据：

$m_y\left(A\right)=\left\{\begin{array}{cc}(1-{\mathrm{\α}}_y)b_y\left(A\right)& \∀A\⋐\mathrm{\Θ},A\≠\mathrm{\Θ}\\ {\mathrm{\α}}_y+(1-{\mathrm{\α}}_y)b_y\left(A\right)& A=\mathrm{\Θ}\end{array}\right.---\left(7\right)$

(8)当在线获取故障特征参数x和y的取值后，它们分别会激活各自2N+1个似然信度向量中的某一个，并分别通过式(5)和式(7)得到相应的诊断证据m_x和m_y，对它们进行综合，得到综合后的诊断证据为：

其中，C可以取故障集合Θ中的任何一个故障类型F₁,F₂,...F_N，或者故障集合Θ本身；

(9)利用步骤(8)得到的综合后诊断证据m_xy(C)，对旋转机械设备的故障进行诊断：那个取值最大的m_xy(C)所对应的C即为真实发生的故障类型。

本发明提出的综合多源故障概率似然信度的旋转机械设备故障诊断方法，对不同故障模式下获取的故障特征监测数据进行区间划分，构造故障特征参数与故障类型的映射关系表；从该表中统计出各故障模式发生的似然信度向量；结合训练样本集合构造目标函数计算得到可靠性因子，用获取的可靠性因子对似然信度向量进行修正，生成诊断证据；在线获取多种故障特征的取值后，分别计算它们激活的诊断证据，再将这些被激活的诊断证据综合，根据综合结果决策得到该在线故障特征所对应的故障类型。根据本发明方法编制的程序(编译环境LabVIEW，C++等)可以在监控计算机上运行，并联合传感器、数据采集器等硬件组成在线监测系统，进行实时的旋转机械设备状态监测与故障诊断。

附图说明

图1是本发明方法的流程框图。

图2是电机转子故障诊断系统图。

图3是本发明方法的实施例中电机转子故障诊断系统结构图。

具体实施方式

本发明提出的一种综合多源故障概率似然信度的旋转机械设备故障诊断方法，其流程框图如图1所示，包括以下各步骤：

(1)设定旋转机械设备的故障集合Θ＝{F₁,...F_i,...F_N}，F_i代表故障集合Θ中的第i个故障，i＝1,2,...,N，N为设备含有故障的个数。

(2)设x为能够反映故障集合Θ中每个故障F_i的故障特征参数，该特征参数由某传感器信息源提供，建立故障特征参数x关于每个故障F_i的特征参数取值变化区间其中，和分别代表区间的左右端点，具体获取的步骤如下：

(2-1)当故障集合Θ中的故障F_i发生时，可获取故障特征参数x的δ个测量样本，记这些测量样本构成的集合为一般500≥δ≥300。

(2-2)对于步骤(2-1)中获取的测量样本集合分别求取这些样本的最小值和最大值，令它们分别为和并构成区间

(3)按照步骤(2)可以获取故障集合Θ中N种故障的特征参数x的测量样本变化区间将这N个区间的2N个左右端点按照从小到大的顺序排序，组成含有2N个点的序列S＝{s₁,s₂,…,s_2N}，其中按照s_j的排序，生成故障特征参数x关于故障集合Θ的2N+1个样本变化区间 $I_1^x=[0,s_1),I_2^x=[s_1,s_2),I_3^x=[s_2,s_3),...,$ $I_{2N}^x={[x}_{2N-1},s_{2N}),I_{2N+1}^x=[s_{2N+1},\∞).$

为了加深对故障特征参数的样本变化区间的理解，这里举例说明。设图3所示的电机转子系统有N＝4种典型故障模式F_i：正常运行F₁，转子不平衡F₂，转子不对中F₃，基座松动F₄，则故障集合Θ＝{F₁,F₂,F₃,F₄}，它们共同的故障特征为振动加速度频谱中1X(1倍频)的幅值x，该振动幅值由振动加速度传感器采集的时域振动信号经快速傅里叶变换后得到。通过步骤(2)获得δ＝500时的4个测量样本集合以及4种故障的特征参数x的测量样本变化区间 $[x_1^L,x_1^R]=\left[\mathrm{0.0175,0.1194}\right],[x_2^L,x_2^R]=\left[\mathrm{0.1511,0.2006}\right],[x_3^L,x_3^R]=\left[\mathrm{0.1578,0.2038}\right],$ $[x_4^L,x_4^R]=[0.3006,0.3476],$ 在步骤(3)中，将这4个区间的8个左右端点按照从小到大的顺序排序，组成含有8个点的序列S＝{0.0175,0.1194,0.1511,0.1578,0.2006,0.2038,0.3006,0.3476}，那么就可以得到故障特征参数x关于故障集合Θ的9个样本变化区间：

$I_1^x=\left[\mathrm{0,0.0175}\right),I_2^x=\left[\mathrm{0.0175,0.1194}\right),I_3^x=\left[\mathrm{0.1194,0.1511}\right),$ $I_4^x=\left[\mathrm{0.1511,0.1578}\right),I_5^x=\left[\mathrm{0.1578,0.2006}\right),I_6^x=\left[\mathrm{0.2006,0.238}\right),$ $I_7^x=\left[\mathrm{0.2038,0.3006}\right),I_8^x\left[\mathrm{0.3006,0.3476}\right),I_9^x=[0.3476,\∞);$

(4)在故障集合Θ中的N种故障分别发生时，通过步骤(2-1)可获得每种故障的测量样本集合 $V_1^x=\{v_{\mathrm{1,1}}^x,v_{\mathrm{1,2}}^x,...,v_{1,\mathrm{\δ}}^x\},V_2^x=\{v_{\mathrm{2,1}}^x,v_{\mathrm{2,2}}^x,...,v_{2,\mathrm{\δ}}^x\},...,$ $V_N^x=\{v_{N,1}^x,v_{N,2}^x,...,v_{N,\mathrm{\δ}}^x\},$ 共计可获得Nδ个样本，用它们可构造故障特征x和故障F₁,F₂,…,F_N之间的映射关系表，如下表1所示，其中l＝1,2，…,2N+1，为步骤(3)中构造的特征参数x测量样本变化区间的个数，a_il表示故障F_i的故障特征参数x的样本落入区间[s_l-1,s_l)中的个数，并有η_l为落入区间[s_l-1,s_l)中的样本个数总和，并有

表1 故障特征参数x与故障类型映射关系表

为了理解上表所示的映射关系表，沿用步骤(3)例子中δ＝500的4个测量样本集合共计获得Nδ＝2000个样本，以及故障特征参数x关于4种故障类型的9个样本变化区间l＝1,2，…,9，即可构造出故障特征参数x和故障F₁，F₂，F₃，F₄之间的映射关系表，如下表2所示

表2 故障特征1X的幅值x与4种故障类型映射关系表

(5)根据步骤(4)中获取的映射关系表，可获得当故障F_i发生时，故障特征参数x的取值落入区间[s_l-1,s_l)的似然函数为：

$p\left(I_l^x\right|F_i)=\frac{a_{\mathrm{il}}}{\mathrm{\δ}}---\left(1\right)$

将式(1)中的似然函数进行归一化，获得当故障特征x的取值落入区间时，故障F_i发生的似然信度为：

$b_x\left(F_i\right)=\frac{p\left(I_l^x\right|F_i)}{{\mathrm{\Σ}}_{i=1}^{N}p\left(I_l^x\right|F_i)}---\left(2\right)$

并有则可以定义此时获取的似然信度向量为：

B_x＝[b_x(F₁),b_x(F₂),…,b_x(F_N),b_x(Θ)]    (3)

其中b_x(Θ)＝0，表示为对故障集合Θ信度赋值为0；若任何故障的特征参数x的样本都未落入区间中，亦即η_l＝0，则此时b_x(F_i)＝0，b_x(Θ)＝1，似然信度向量的取值为：

B_x＝[0,0,…,0,1]     (4)

给定一个故障特征参数x的取值，其必然落入中的某一个区间，此时该区间所对应的似然信度向量被激活，则可以被激活的似然信度向量的个数为2N+1，它们分别由式(3)和式(4)给出。

继续沿用步骤(3)中的例子来加深对似然信度向量的理解，设转子不对中故障F₃发生时，故障特征参数x的一个取值落入映射关系表2的区间中，此时该区间所对应的似然信度向量被激活，由式(1)得该取值对应的各似然函数分别为 $p\left(I_5^x\right|F_1)=\frac{0}{500}=0,p\left(I_5^x\right|F_2)=\frac{426}{500}=0.852,$ $p\left(I_5^x\right|F_3)=\frac{489}{500}=0.978,p\left(I_5^x\right|F_4)=\frac{0}{500}=0,$ 然后由式(2)可得故障F₃发生的似然信度为：

$b_x\left(F_3\right)=\frac{p\left(I_5^x\right|F_3)}{{\mathrm{\Σ}}_{i=1}^{4}p\left(I_5^x\right|F_i)}=\frac{0.789}{0+0.852+0+0.978}=0.5344$

同样地，可求取其它故障发生的似然信度b_x(F₁)＝0，b_x(F₂)＝0.4656，b_x(F₄)＝0，并且b_x(Θ)＝0，那么就可获取该取值对应的似然信度向量为：

B_x＝[0,0.4656,0.5344,0,0]

(6)定义可靠性因子α_x用其对似然信度向量B_x中的元素进行折扣，生成关于特征参数x的诊断证据：

$m_x\left(A\right)=\left\{\begin{array}{cc}(1-{\mathrm{\α}}_x)b_x\left(A\right)& \∀A\⋐\mathrm{\Θ},A\≠\mathrm{\Θ}\\ {\mathrm{\α}}_x+(1-{\mathrm{\α}}_x)b_x\left(A\right)& A=\mathrm{\Θ}\end{array}\right.---\left(5\right)$

求取可靠性因子α_x的步骤如下：

(6-1)在故障集合Θ中的N种故障分别发生时，再一次利用步骤(2-1)获得每种故障在故障特征参数x下的训练样本集合 $T_1^x=\{t_{\mathrm{1,1}}^x,t_{\mathrm{1,2}}^x,...,t_{1,\mathrm{\ω}}^x\},T_2^x=\{t_{\mathrm{2,1}}^x,t_{\mathrm{2,2}}^x,...,t_{2,\mathrm{\ω}}^x\},...,T_N^x=\{t_{N,1}^x,t_{M,2}^x,...,t_{N,\mathrm{\ω}}^x\},$ 一般200≥ω≥50，共计可获得Nω个样本，它们组成的集合定义为

(6-2)构造目标函数

它是关于可靠性因子α_x的函数，其中，b_x(Θ)为根据落入的特征参数x的样本变化区间，由步骤(5)确定的故障类型F_i和故障集合Θ的似然信度赋值，当则否则

(6-3)设定α_x的取值范围为0到1之间，利用数学计算软件Matlab中的fmincom非线性优化函数，最小化Obj(α_x)的取值，此时对应获得的Obj函数的输入值即为寻找的可靠性因子α_x的取值。

为了加深对取值的理解，在前例的基础上，假设给定一个在故障转子不平衡F₂发生时获得的训练样本并且落入区间中，得到似然信度向量B_x＝[0,0.4656,0.5344,0,0]，已知该样本是在故障F₂发生时采集的，那么 的取值均为0。

(7)除了故障特征参数x，若存在另一种故障特征参数y也能反映Θ中的各个故障，则重复以上步骤(2)至(6)，可以获得关于y的诊断证据：

$m_y\left(A\right)=\left\{\begin{array}{cc}(1-{\mathrm{\α}}_y)b_y\left(A\right)& \∀A\⋐\mathrm{\Θ},A\≠\mathrm{\Θ}\\ {\mathrm{\α}}_y+(1-{\mathrm{\α}}_y)b_y\left(A\right)& A=\mathrm{\Θ}\end{array}\right.---\left(7\right)$

(8)当在线获取故障特征参数x和y的取值后，它们分别会激活各自2N+1个似然信度向量中的某一个，并分别通过式(5)和式(7)得到相应的诊断证据m_x和m_y，对它们进行综合，得到综合后的诊断证据为：

其中，C都可以取故障集合Θ中的任何一个故障类型F₁,F₂,...F_N，或者故障集合Θ本身。

为加深对式(8)的理解，这里举例说明，若有振动加速度频谱中1X(1倍频)的幅值x和2X(2倍频)的幅值y两个故障特征参数通过式(5)和式(7)获得的诊断证据为：

m_x＝[m_x(F₁),m_x(F₂),m_x(F₃),m_x(F₄),m_x(Θ)]＝[0.3,0.2,0.5,0,0]

m_y＝[m_y(F₁),m_y(F₂),m_y(F₃),m_y(F₄),m_y(Θ)]＝[0,0.5,0.5,0,0]

诊断证据m_x和m_y之间的冲突程度则为：

那么利用(8)式综合后对各故障模式以及全集的诊断结果为：

这里以上两证据综合的公式也适用于多个故障特征参数提供的多个诊断证据的融合，只需将某两个证据综合的结果再与另一证据用式(8)综合，这三个证据的综合结果再与第四条证据用式(8)综合，以此类推综合所有证据即可。

(9)利用步骤(8)得到的综合后诊断证据m_xy(C)，对旋转机械设备的故障进行诊断：那个取值最大的m_xy(C)所对应的C即为真实发生的故障类型。

以下结合附图，详细介绍本发明方法的实施例：

本发明方法的流程图如图1所示，核心部分是：确定故障特征参数关于各故障类型的样本变化区间；构造故障特征参数与故障类型的映射关系表；从该表中统计出各故障模式发生的似然信度向量；结合训练样本集合构造目标函数，并优化该函数获得可靠性因子，用该可靠性因子对似然信度向量进行修正，生成诊断证据；在线获取多种故障特征的取值后，分别计算它们激活的诊断证据，再将这些被激活的诊断证据综合，根据综合结果决策得到该在线故障特征所对应的故障类型。

以下结合图2中电机转子故障诊断系统的最佳实施例，详细介绍本发明方法的各个步骤。

1、电机转子故障诊断系统设置实例

实验设备如图3中的ZHS-2型多功能电机柔性转子系统，振动位移传感器和振动加速度传感器分别安置在转子支撑座的水平和垂直方向采集转子振动信号，两个传感器采集到的振动信号传入HG-8902数据采集箱，经过信号调理电路处理后，最终经过A/D转换器输出至监控计算机，然后利用Labview环境下的HG-8902数据分析软件得到转子振动加速度频谱以及时域振动位移平均幅值作为故障特征信号。

2、电机转子故障设置及故障特征参数的选取

根据试验台的具体特性，分别在试验台上设置了以下4种典型故障模式：正常运行，转子不平衡，转子不对中，基座松动。通过对大量实验数据的分析可知，引发异常振动的故障源都会产生一定频率成分的振动幅值增加或减少。因此，这里选取1～3X倍频以及时域振动位移平均幅值作为故障特征量。设定转子转速为1500r/m，则基频1X为25Hz，n倍频nX，n＝1,2,3,…,为(n×25)Hz。将频域的1X～3X的振动幅值以及时域振动位移4种特征信息进行综合做出综合决策。

3、确定振动加速度1X～3X的幅值和时域振动位移平均幅值这4种故障特征参数分别关于4种故障模式变化的特征区间

利用本发明方法步骤(2)，分别在电机转子上设置“F₁”、“F₂”、“F₃”、“F₄”这N＝4种故障，并利用振动加速度传感器和振动位移传感器分别获取振动加速度1X、2X、3X的幅值和时域振动位移平均幅值这4种故障特征参数。对于每种故障特征参数，分别在4种故障模式下，以时间间隔Δt＝16s内连续采集δ＝500次测量，即可获取2000个测量样本，用于建立故障特征参数关于4种故障的2N+1＝9个样本变化区间，根据本发明方法步骤(2)得到4种故障特征参数对应4种故障，共4组这样的样本变化区间。

4、构造1X、2X、3X振动幅值和时域振动位移平均幅值(简称“位移”)这4种故障特征参数与故障F₁-F₄之间的映射关系表

以故障特征参数1X为例，利用本发明方法步骤(2)-(3)获得传感器测量故障特征参数1X的2000个测量数据和故障特征参数1X关于4种故障的9个样本变化区间，构造如本发明方法步骤(4)中的表1所示的映射关系表，如下表3所示。按照以上过程可以构造其他3种故障特征参数与故障类型的映射关系表，分别如下表4、表5和表6所示。

表3 故障特征参数1X与故障类型映射关系表

表4 故障特征参数2X与故障类型映射关系表

表5 故障特征参数3X与故障类型映射关系表

表6 故障特征参数时域振动位移平均幅值与故障类型映射关系表

表3-表6中样本变化区间的取值如下：

$I_1^{1X}=\left[\mathrm{0,0.0175}\right),I_2^{1X}=\left[\mathrm{0.0175,0.1194}\right),I_3^{1X}\left[\mathrm{0.1194,0.1511}\right),$ $I_4^{1X}=\left[\mathrm{0.1511,0.1578}\right),I_5^{1X}=\left[\mathrm{0.1578,0.2006}\right),I_6^{1X}=\left[\mathrm{0.2006,0.2038}\right),$ $I_7^{1X}=\left[\mathrm{0.2038,0.3006}\right),I_8^{1X}=\left[\mathrm{0.3006,0.476}\right),I_9^{1X}=[0.3476,\∞);I_1^{2X}=\left[\mathrm{0,0.0169}\right),$ $I_2^{2X}=\left[\mathrm{0.0169,0.1222}\right),I_3^{2X}=\left[\mathrm{0.1222,0.1321}\right),I_4^{2X}=\left[\mathrm{0.1321,0.1652}\right),$ $I_5^{2X}=\left[\mathrm{0.1652,0.2801}\right),I_6^{2X}=\left[\mathrm{0.2801,0.3079}\right),I_7^{2X}=\left[\mathrm{0.3079,0.3510}\right),$ $I_8^{2X}=\left[\mathrm{0.3510,0.3647}\right),I_9^{2X}=[0.3647,\∞);I_1^{3X}=\left[\mathrm{0,0.0446}\right),I_2^{3X}=\left[\mathrm{0.0446,0.0934}\right),$ $I_3^{3X}=\left[\mathrm{0.0934,0.1179}\right),I_4^{3X}=\left[\mathrm{0.1179,0.1296}\right),I_5^{3X}=\left[\mathrm{0.1296,0.1384}\right),$ $I_6^{3X}=\left[\mathrm{0.1384,0.864}\right),I_7^{3X}=\left[\mathrm{0.1864,0.1869}\right),I_8^{3X}=\left[\mathrm{0.1869,0.3040}\right),$

5、根据本发明方法(步骤5)求取4种故障特征参数的取值落入样本变化区间时的似然信度向量

在根据本发明方法步骤(4)获取各故障特征参数与故障类型的映射关系表之后，依照本发明方法的步骤(5)获得振动加速度1X、2X、3X的幅值和时域振动位移平均幅值这4种特征参数各区间的似然信度向量，如下表7至表10所示。

表7 故障特征参数1X各区间对应的似然信度向量

表8 故障特征参数2X各区间对应的似然信度向量

表9 故障特征参数3X各区间对应的似然信度向量

表10 故障特征参数时域振动位移平均幅值各区间对应的似然信度向量

6、根据本发明方法步骤(6)获取可靠性因子，并用其对4种故障特征参数对应的似然信度向量中的元素进行修正，具体过程如下：

对于每种故障，同时在时间间隔Δt＝4s内连续采集故障特征参数1X的ω＝50次测量值，获得总共200个样本构成训练样本集合T^x，用其构造形如式(6)的目标函数，最小化该函数可获取1X的似然信度向量的可靠性因子α_1X＝0.0071

按照以上过程，可以得到其它3种故障特征参数2X、3X和时域振动位移平均幅值的可靠性因子的取值分别为：

α_2X＝0.0245，α_3X＝0.0143，α_位移＝0.5868

7、当在线获取故障特征参数1X、2X、3X和时域振动位移平均幅值的取值之后，它分别会激活各自7个似然信度向量中的某一个，并分别通过本发明方法步骤(6)中的式(5)得到相应的诊断证据m_1X、m_2X，m_3X和m_位移。例如当1X＝0.1618，2X＝0.1529、3X＝0.1094和时域振动位移平均幅值＝4.3645时，获取的诊断证据分别为m_1X(F₁)＝0、m_1X(F₂)＝0.4623、m_1X(F₃)＝0.5306，m_1X(F₄)＝0，m_1X(Θ)＝0.0071；

m_2X(F₁)＝0、m_2X(F₂)＝0.9755、m_2X(F₃)＝0，m_2X(F₄)＝0，m_1X(Θ)＝0.0245；

m_3X(F₁)＝0.2639、m_3X(F₂)＝0.7218、m_3X(F₃)＝0，m_3X(F₄)＝0，m_3X(Θ)＝0.0143；

m_位移(F₁)＝0、m_位移(F₂)＝0.2800、m_位移(F₃)＝0.1333，m_位移(F₄)＝0，m_位移(Θ)＝0.5868

针对该例中的这组“1X＝0.1618，2X＝0.1529、3X＝0.1094和时域振动位移平均幅值＝4.3645”采样值，其分别激活的样本变化区间为和经各自特征的可靠性因子修正后得到如上的4个诊断证据。

根据本发明方法步骤(8)将4种故障特征参数提供的诊断证据综合，得到如表11所示的综合结果

表11 测试样本的综合诊断证据

诊断证据	m(F1)	m(F2)	m(F3)	m(F4)	m(Θ)
						综合结果	0.0001	0.9995	0.0004	0	0

按照本发明方法步骤(9)中的决策规则，可判断故障F₂发生，这与这组样本“1X＝0.1618，2X＝0.1529、3X＝0.1094和时域振动位移平均幅值＝4.3645”采集时设定的故障一致，说明决策结果正确。

对于4种故障，分别取4种特征参数的150组在线样本进行测试，测试样本的故障诊断结果如表12所示，总故障确准率达到了99.7％，达到一般诊断系统的确准率要求。

表12 测试样本的故障诊断结果

故障模式	F1	F2	F3	F4
					确准率	100％	100％	99.33％	99.33％

Claims (1)

1.综合多源故障概率似然信度的旋转机械设备故障诊断方法，其特征在于该方法包括以下各步骤：

(1)设定旋转机械设备的故障集合Θ＝{F₁,...F_i,...F_N}，F_i代表故障集合Θ中的第i个故障，i＝1,2,...,N，N为设备含有故障的个数；

(2)设x为能够反映故障集合Θ中每个故障F_i的故障特征参数，该特征参数由某传感器信息源提供，建立故障特征参数x关于每个故障F_i的特征参数取值变化区间其中，和分别代表区间的左右端点，具体获取的步骤如下：

(2-1)当故障集合Θ中的故障F_i发生时，可获取故障特征参数x的δ个测量样本，记这些测量样本构成的集合为

(2-2)对于步骤(2-1)中获取的测量样本集合分别求取这些样本的最小值和最大值，令它们分别为和并构成区间

(3)按照步骤(2)可以获取故障集合Θ中N种故障的特征参数x的测量样本变化区间将这N个区间的2N个左右端点按照从大到小的顺序排序，组成含有2N个点的序列S＝{s₁,s₂,…,s_2N}，其中按照s_j的排序，生成故障特征参数x关于故障集合Θ的2N+1个样本变化区间 $I_1^x=[0,s_1),I_2^x=[s_1,s_2),I_3^x=[s_2,s_3),...,I_{2N}^x=[s_{2N-1},s_{2N}),I_{2N+1}^x=[s_{2N+1},\∞);$

(4)在故障集合Θ中的N种故障分别发生时，通过步骤(2-1)可获得每种故障的测量样本集合 $V_1^x=\{v_{\mathrm{1,1}}^x,v_{\mathrm{1,2}}^x,...,v_{1,\mathrm{\δ}}^x\},V_2^x=\{v_{\mathrm{2,1}}^x,v_{\mathrm{2,2}}^x,...,v_{2,\mathrm{\δ}}^x\},...,V_N^x=\{v_{N,1}^x,v_{N,2}^x,...,v_{N,\mathrm{\δ}}^x\},$ 共计可获得Nδ个样本，用它们可构造故障特征x和故障F₁,F₂,…,F_N之间的映射关系表，如下表1所示，其中l＝1,2，…,2N+1，为步骤(3)中构造的特征参数x测量样本变化区间的个数，a_il表示故障F_i的故障特征参数x的样本落入区间[s_l-1,s_l)中的个数，并有η_l为落入区间[s_l-1,s_l)中的样本个数总和，并有

表1 故障特征参数x与故障类型映射关系表

(5)根据步骤(4)中获取的映射关系表，可获得当故障F_i发生时，故障特征参数x的取值落入区间[s_l-1,s_l)的似然函数为

$p\left(I_l^x\right|F_i)=\frac{a_{\mathrm{il}}}{\mathrm{\δ}}---\left(1\right)$

将式(1)中的似然函数进行归一化，获得当故障特征x的取值落入区间时，故障F_i发生的似然信度为

$b_x\left(F_i\right)=\frac{p\left(I_l^x\right|F_i)}{{\mathrm{\Σ}}_{i=1}^{N}p\left(I_l^x\right|F_i)}---\left(2\right)$

并有则可以定义此时获取的似然信度向量为

B_x＝[b_x(F₁),b_x(F₂),…,b_x(F_N),b_x(Θ)]   (3)

其中b_x(Θ)＝0，表示为对故障集合Θ信度赋值为0；若任何故障的特征参数x的样本都未落入区间中，亦即η_l＝0，则此时b_x(F_i)＝0，b_x(Θ)＝1，似然信度向量的取值为

B_x＝[0,0,…,0,1]           (4)

给定一个故障特征参数x的取值，其必然落入中的某一个区间，此时该区间所对应的似然信度向量被激活，则可以被激活的似然信度向量的个数为2N+1，它们分别由式(3)和式(4)给出；

(6)定义可靠性因子α_x用其对似然信度向量B_x中的元素进行折扣，生成关于特征参数x的诊断证据

$m_x\left(A\right)=\left\{\begin{array}{cc}(1-{\mathrm{\α}}_x)b_x\left(A\right)& \∀A\⋐\mathrm{\Θ},A\≠\mathrm{\Θ}\\ {\mathrm{\α}}_x+(1-{\mathrm{\α}}_x)b_x\left(A\right)& A=\mathrm{\Θ}\end{array}\right.---\left(5\right)$

求取可靠性因子α_x的步骤如下

(6-1)在故障集合Θ中的N种故障分别发生时，再一次利用步骤(2-1)获得每种故障在故障特征参数x下的训练样本集合 $T_1^x=\{t_{\mathrm{1,1}}^x,t_{\mathrm{1,2}}^x,...,t_{1,\mathrm{\ω}}^x\},T_2^x=\{t_{\mathrm{2,1}}^x,t_{\mathrm{2,2}}^x,...,t_{2,\mathrm{\ω}}^x\},{...,T}_N^x=\{t_{N,1}^x,t_{N,2}^x,...,t_{N,\mathrm{\ω}}^x\},$ 共计可获得Nω个样本，它们组成的集合定义为

(6-2)构造目标函数

它是关于可靠性因子α_x的函数，其中，b_x(F_i)，b_x(Θ)为根据落入的特征参数x的样本变化区间，由步骤(5)确定的故障类型F_i和故障集合Θ的似然信度赋值，当 $t_k\∈T_i^x,$ 则否则

(6-3)设定α_x的取值范围为0到1之间，利用数学计算软件Matlab中的fmincom非线性优化函数，最小化Obj(α_x)的取值，此时对应获得的Obj函数的输入值即为寻找的可靠性因子α_x的取值；

(7)除了故障特征参数x，若存在另一个故障特征参数y也能反映Θ中的各个故障，则重复以上步骤(2)至(6)，获得关于y的诊断证据

$m_y\left(A\right)=\left\{\begin{array}{cc}(1-{\mathrm{\α}}_y)b_y\left(A\right)& \∀A\⋐\mathrm{\Θ},A\≠\mathrm{\Θ}\\ {\mathrm{\α}}_y+(1-{\mathrm{\α}}_y)b_y\left(A\right)& A=\mathrm{\Θ}\end{array}\right.---\left(7\right)$

(8)当在线获取故障特征参数x和y的取值后，它们分别会激活各自2N+1个似然信度向量中的某一个，并分别通过式(5)和式(7)得到相应的诊断证据m_x和m_y，对它们进行综合，得到综合后的诊断证据为

其中，C可以取故障集合Θ中的任何一个故障类型F₁,F₂,...F_N，或者故障集合Θ本身；

(9)利用步骤(8)得到的综合后诊断证据m_xy(C)，对旋转机械设备的故障进行诊断：那个取值最大的m_xy(C)所对应的C即为真实发生的故障类型。