一种机电控制系统的初值预估迭代学故障诊断算法

摘要

本发明公开了一种机电控制系统的初值预估迭代学故障诊断算法。首先建立执行器单元发生故障的机电控制系统模型；构建执行器故障机电控制系统的离散状态空间方程；然后设计迭代学故障诊断算法；进一步分析故障诊断算法的收敛性和阈值选择条件；最后在迭代学故障诊断算法中加入初值预估算法，实现对机电控制系统执行器单元的实时故障诊断。其优点是：初值预估迭代学故障诊断算法可以减小诊断过程中前面离散采样点故障的估计误差对后续离散采样点故障诊断的影响，从而有效减少迭代次数,提高故障诊断效率；故障算法结构简单，不仅可以检测和重构机电控制系统的各类型执行器故障，而且易于工程实现，便于实时在线诊断。

主权项

一种机电控制系统的初值预估迭代学习故障诊断算法，其特征包括：建立执行器单元发生故障的机电控制系统模型；构建执行器故障机电控制系统的离散状态空间方程；设计迭代学习故障诊断算法；分析迭代学习故障诊断算法的收敛性和阈值选择条件；在迭代学习故障诊断算法中加入初值预估算法，实现机电控制系统的实时故障诊断；第一步：建立执行器单元发生故障的机电控制系统模型考虑到机电控制系统在实际应用过程中会受到扰动的影响，因此执行器单元发生故障的机电控制系统模型可以描述为：$\left\{\begin{array}{c}{R}_a{i}_a\left(t\right)+L_a\frac{{\mathrm{di}}_a\left(t\right)}{dt}+C_e\omega \left(t\right)=e\left(t\right)+w_1\left(t\right)\\ J\frac{d\omega \left(t\right)}{dt}+C_f\omega \left(t\right)+f\left(t\right)=C_M{i}_a\left(t\right)+w_2\left(t\right)\end{array}\right.---\left(1\right)$其中R_a为电枢电阻，J为转子转动惯量，L_a为电枢电感，C_e为反电动势系数，C_f为电机轴机械阻尼系数，C_M为转矩系数，ω(t)为电动机转速，i_a(t)为电枢电流，e(t)为输入电压，f(t)为执行器单元发生的故障，故障信号在电动机运行过程中通常表现为一类非线性函数，w₁(t)和w₂(t)分别为电动机运行过程中受到的扰动；第二步：构建执行器故障机电控制系统的离散状态空间方程利用电动机的电枢电流和转速定义状态变量x(t)＝[i_a(t) ω(t)]^T，利用输入电压定义系统输入变量u(t)＝e(t)，定义扰动变量w(t)＝[w₁(t) w₂(t)]^T，则式(1)所示的一种机电控制系统可描述为：$\stackrel{·}{x}\left(t\right)=\left[\begin{array}{cc}-\frac{R_a}{L_a}& -\frac{C_e}{L_a}\\ \frac{C_M}{J}& -\frac{C_f}{J}\end{array}\right]x\left(t\right)+\left[\begin{array}{c}\frac{1}{L_a}\\ 0\end{array}\right]u\left(t\right)+\left[\begin{array}{c}0\\ -\frac{1}{J}\end{array}\right]f\left(t\right)+\left[\begin{array}{cc}\frac{1}{L_a}& 0\\ 0& \frac{1}{J}\end{array}\right]w\left(t\right)---\left(2\right)$显然式(2)为连续系统模型，因此需要对式(2)进行离散化，选取满足香农采样定理的采样周期T_s，并将电动机转速作为输出，同时考虑传感器检测电动机转速过程中受到的扰动影响，进一步可以得到执行器故障机电控制系统的离散状态空间模型：$\begin{array}{c}x(i+1)+Ax\left(i\right)+B_{f}f\left(i\right)+Bu\left(i\right)+B_{w}w\left(i\right)\\ y\left(i\right)=Cx\left(i\right)+D_{v}v\left(i\right)\end{array}---\left(3\right)$式中i为对机电控制系统(2)离散化后的采样时刻；u(i)，y(i)和x(i)分别是机电控制系统在各个离散采样时刻的输入，输出和状态向量；f(i)为待诊断的执行器故障在各个采样时刻的信号；w(i)，v(i)分别表示式(3)离散采样系统的状态和传感器输出部分在采样时刻受到的有限能量扰动，即||v(i)||＜d_v，||w(i)||＜d_w；A，B，C，B_f，B_w和D_v为离散采样系统(3)中的矩阵；矩阵对(C,A)可观测且系统运行的初始状态x(0)＝x₀；第三步：设计迭代学习故障诊断算法针对式(3)设计如下形式的迭代学习故障跟踪估计器：$\begin{array}{c}{\hat{x}}_k(i+1)=A{\hat{x}}_k\left(i\right)+B_f{\hat{f}}_k\left(i\right)+Bu\left(i\right)+L\left(y\right(i)-\hat{y}(i\left)\right)\\ {\hat{y}}_k\left(i\right)=C\hat{x}\left(i\right)\end{array}---\left(4\right)$其中分别为迭代学习故障跟踪估计器的离散输出和状态向量，下标k表示故障跟踪估计器在每个离散采样点的迭代估计次数，故障跟踪估计器的初始状态矩阵L是预先设计的增益矩阵，保证矩阵(A‑LC)的特征根位于单位圆内；为故障跟踪估计器在第i个离散采样点经第k次迭代学习调节后的虚拟故障信号；结合式(3)和式(4)的输出信号定义离散采样系统的残差：$r_k\left(i\right)=y\left(i\right)-{\hat{y}}_k\left(i\right)---\left(5\right)$迭代学习故障诊断算法中的虚拟故障估计算法设计为：$\{\begin{array}{cc}{\hat{f}}_{k+1}\left(i\right)={\hat{f}}_k\left(i\right)+{\mathrm{\Gamma r}}_k(i+1,)& ||r_k(i+1)||>\delta \left(i\right)\\ {\hat{f}}_{k+1}\left(i\right)={\hat{f}}_k\left(i\right)& ||r_k(i+1)||\le \delta \left(i\right)\end{array}---\left(6\right)$式(6)中虚拟故障估计算法的启动条件：$||y\left(i\right)-{\hat{y}}_0\left(i\right)||>\delta \left(i\right)---\left(7\right)$式中r_k(i)是在第i个采样点进行的第k次故障迭代估计后系统(3)的实际输出和故障跟踪估计器(4)输出之间的残差信号；虚拟故障估计的过程就是不断利用残差r_k(i)来更新虚拟故障Γ为所选取的迭代学习参数，δ(i)是判定式(3)的采样点是否发生故障的阈值；本发明的故障诊断算法利用故障跟踪估计器(4)对离散采样系统(3)的输出进行跟踪，反复更新故障跟踪估计器(4)中的虚拟故障(6)使其逐渐接近系统实际故障从而达到故障估计的目的；首先根据顺序依次将离散系统输出和故障跟踪估计器输出作比较得到残差信号，若其残差信号在阈值限定范围内则判定该点无故障发生，否则判定为发生故障并启动迭代学习故障估计算法，然后利用残差信号调节虚拟故障的值并再次比较，反复执行此步骤直到残差值收敛到阈值范围内，则故障诊断过程结束并将此时得到的虚拟故障作为实际故障的估计值，然后再利用同样方法进行下一个采样点的故障诊断；第四步：分析迭代学习故障诊断算法的收敛性和阈值选择条件定义采样点经第k次迭代学习诊断后的状态误差为：$e_k\left(i\right)=x\left(i\right)-{\hat{x}}_k\left(i\right)---\left(8\right)$虚拟故障的估计误差为：${\mathrm{\Delta f}}_k\left(i\right)=f\left(i\right)-\hat{f}\left(i\right)---\left(9\right)$由式(8)、(9)和式(3)、(4)可得误差系统的状态方程为：$\begin{array}{c}{e}_k(i+1)=(A-LC)e_k\left(i\right)+B_f{\mathrm{\Delta f}}_k\left(i\right)+B_{w}w\left(i\right)-{\mathrm{LD}}_{v}v\left(i\right)\\ r_k\left(i\right)={\mathrm{Ce}}_k\left(i\right)+D_{v}v\left(i\right)\end{array}---\left(10\right)$由式(10)可得$e_k(i+1)=\Phi (i,0)e_k\left(0\right)+\underset{s=0}{\overset{i}{\Sigma }}\Phi (i,s+1)\left(B_f\Delta f\right(s)+B_{w}w(s)-{\mathrm{LD}}_{v}v\left(s\right))---\left(11\right)$其中Φ(·)是式(10)的状态转移矩阵，对式(11)两边取范数后可得$||e_k(i+1)||\le ||\Phi (i,0)e_k\left(0\right)||+\underset{s=0}{\overset{i}{\Sigma }}||\Phi (i,s+1)B_f\Delta f\left(s\right)||+\underset{s=0}{\overset{i}{\Sigma }}||\Phi (i,s+1)\left(B_{w}w\right(s)-{\mathrm{LD}}_{v}v(s\left)\right)||---\left(12\right)$又由式(10)可得r_k(i)≤||C||||e_k(i)||+||D_vv(i)||＝c||e_k(i)||+d₁    (13)其中c＝||C||，结合式(12)可得r_k(i)≤cβ+c(a₁+a₂)+d₁＝δ(i)    (14)其中为系统状态干扰的累积量,代表系统输出部分传感器受到的干扰累积量,是前面虚拟故障估计误差的累积量；因为算法是按照采样顺序逐点诊断故障，所以每个采样点都要进行故障判定，即每个点都需要有其对应的阈值，阈值大小由系统干扰和故障误差确定，而故障估计误差的范围又是通过干扰确定的，所以阈值实际上是由干扰和系统参数决定，只要获取系统参数值和干扰范围就可以计算出各采样点的阈值大小；针对式(3)设计的迭代学习故障跟踪估计器(4)，若其参数满足||I‑ΓCB_f||＜1的条件，则当迭代学习次数k→∞时,在范数意义下迭代学习故障跟踪估计器的最终输出逼近离散采样系统的实际采样输出y(i)，虚拟故障逼近采样时刻的实际故障f(i)；本发明的迭代学习故障诊断算法交换了故障诊断过程中时域和迭代域的次序，因而能够根据采样顺序逐点诊断故障，进而减小前面采样点故障在诊断过程中的估计误差对后续采样点的影响，更重要的是为后续加入故障初值预估算法提供了条件；由于需要对每个采样点进行收敛性分析，常规的λ范数分析法不再适用，所里利用递推法进行收敛性分析；由式(10)可得r_k(i+1)＝C(A‑LC)e_k(i)+CB_fΔf_k(i)+CB_ww(i)‑CLD_vv(i)+D_vv(i)    (15)进一步由式(6)和式(9)可得Δf_k+1(i)＝Δf_k(i)‑Γr_k(i+1)    (16)将式(15)带入式(16)后可以得到$\begin{array}{c}{\mathrm{\Delta f}}_{k+1}\left(i\right)={\mathrm{\Delta f}}_k\left(i\right)-{\mathrm{\Gamma CB}}_f{\mathrm{\Delta f}}_k\left(i\right)-\Gamma C(A-LC)e_k\left(i\right)-{\mathrm{\Gamma CB}}_{w}w\left(i\right)-\Gamma (CL-I)D_{v}v\left(i\right)\\ =(I-{\mathrm{\Gamma CB}}_f){\mathrm{\Delta f}}_k\left(i\right)-\Gamma C(A-LC)e_k\left(i\right)-{\mathrm{\Gamma CB}}_{w}w\left(i\right)-\Gamma (CL-I)D_{v}v\left(i\right)\end{array}---\left(17\right)$再将式(10)带入式(17)后可得$\begin{array}{c}{\mathrm{\Delta f}}_{k+1}\left(i\right)=(I-{\mathrm{\Gamma CB}}_f){\mathrm{\Delta f}}_k\left(i\right)-{\mathrm{\Gamma CB}}_{w}w\left(i\right)-\Gamma (CL-I)D_{v}v\left(i\right)-\\ \Gamma C(A-LC)\underset{s=0}{\overset{i-1}{\Sigma }}\Phi (i,s+1)\left(B_f\Delta f\right(s)+B_{w}w(s)-{\mathrm{LD}}_{v}v(s\left)\right)\end{array}---\left(18\right)$当i＝0时Δf_k+1(0)＝(I‑ΓCB_f)Δf_k(0)‑ΓCB_ww(0)‑Γ(CL‑I)D_vv(0)    (19)对式(19)两边取范数可得$\begin{array}{c}||{\mathrm{\Delta f}}_{k+1}\left(0\right)||\le ||I-{\mathrm{\Gamma CB}}_f||||{\mathrm{\Delta f}}_k\left(0\right)||+||{\mathrm{\Gamma CB}}_{w}w\left(0\right)||+||\Gamma (CL-I)D_{v}v\left(0\right)||\\ \le \rho ||{\mathrm{\Delta f}}_k\left(0\right)||+a_3\end{array}---\left(20\right)$其中ρ＝||I‑ΓCB_f||,进一步可得$\begin{array}{c}||{\mathrm{\Delta f}}_{k+1}\left(1\right)||\le \rho ||{\mathrm{\Delta f}}_k\left(0\right)||+a_3\le {\rho }^2||{\mathrm{\Delta f}}_{k-1}\left(0\right)||+{\mathrm{\rho a}}_3+a_3\\ \le {\rho }^3||{\mathrm{\Delta f}}_{k-2}\left(0\right)||+{\rho }^2{a}_3+{\mathrm{\rho a}}_3+a_3\le \mathrm{...}\le {\rho }^k||{\mathrm{\Delta f}}_1\left(0\right)||+{\rho }^{k-1}{a}_3+\mathrm{...}+{\rho }^2{a}_3+{\mathrm{\rho a}}_3+a_3\\ ={\rho }^k||{\mathrm{\Delta f}}_1\left(0\right)||+\frac{1-{\rho }^k}{1-\rho }{a}_3\end{array}---\left(21\right)$由于0＜ρ＜1,则$||{\mathrm{\Delta f}}_k\left(0\right)||\to \frac{1-{\rho }^k}{1-\rho }{a}_3,k\to \infty ---\left(22\right)$必然存在0＜K(0)＜∞,ε＞0使得k＝K(0)，同时$||{\mathrm{\Delta f}}_k\left(0\right)||\le \frac{1-{\rho }^k}{1-\rho }{a}_3+ϵ=\alpha \left(0\right)+\eta \left(0\right)ϵ=\Delta F\left(0\right)---\left(23\right)$其中η(0)＝1，由式(10)可得e_K(0)(1)＝Φ(1,1)(B_fΔf_K(0)(0)+B_ww(0)‑LD_vv(0))    (24)对式(24)两边取范数后进一步得到其中θ(1)＝a₄+bα(0)，b＝||B_f||，而且其中γ(1)＝cθ(1)+d₁,R(1)为此时的残差最大值；当i＝1时$\begin{array}{c}{\mathrm{\Delta f}}_{k+1}\left(1\right)=(I-{\mathrm{\Gamma CB}}_f){\mathrm{\Delta f}}_k\left(1\right)-\Gamma C(A-LC)e_k\left(1\right)\\ -{\mathrm{\Gamma CB}}_{w}w\left(i\right)-\Gamma (CL-I)D_{v}v\left(i\right)\end{array}---\left(27\right)$对式(27)两边取范数后得到$\begin{array}{c}||{\mathrm{\Delta f}}_{k+1}\left(1\right)||\le ||I-{\mathrm{\Gamma CB}}_f||||{\mathrm{\Delta f}}_k\left(1\right)||+||\Gamma C(A-LC)e_{K\left(0\right)}\left(1\right)||\\ +||{\mathrm{\Gamma CB}}_{w}w\left(1\right)||+||\Gamma (CL-I)D_{v}v\left(1\right)||\\ \le \rho ||{\mathrm{\Delta f}}_k\left(1\right)||+a_3+lE\left(1\right)\end{array}---\left(28\right)$其中l＝||Γ(i)C(A‑LC)||,则$\begin{array}{c}||{\mathrm{\Delta f}}_{k+1}\left(1\right)||\le \rho ||{\mathrm{\Delta f}}_k\left(1\right)||+a_1+lE\left(1\right)\le {\rho }^2||{\mathrm{\Delta f}}_{k-1}\left(1\right)||+\rho (a_1+lE(1\left)\right)+a_1+lE\left(1\right)\\ \le {\rho }^3||{\mathrm{\Delta f}}_{k-2}\left(1\right)||+({\rho }^2+\rho +1)(a_1+lE(1\left)\right)\le \mathrm{...}\le {\rho }^k||{\mathrm{\Delta f}}_1\left(1\right)||+({\rho }^{k-1}+\mathrm{...}+{\rho }^2+\rho +1)(a_1+lE(1\left)\right)\\ ={\rho }^k||{\mathrm{\Delta f}}_1\left(1\right)||+\frac{1-{\rho }^k}{1-\rho }{a}_1+\frac{1-{\rho }^k}{1-\rho }lE\left(1\right)\end{array}---\left(29\right)$由于ρ＜1,因此$||{\mathrm{\Delta f}}_k\left(1\right)||\to \frac{1-{\rho }^k}{1-\rho }{a}_1+\frac{1-{\rho }^k}{1-\rho }lE\left(1\right),k\to \infty ---\left(30\right)$此时存在0＜K(1)＜∞，ε＞0使得k＝K(1)，并且其中由式(15)可得$e_{K\left(1\right)}\left(2\right)=\underset{s=0}{\overset{1}{\Sigma }}\Phi (2,s+1)\left(B_f{\mathrm{\Delta f}}_{K\left(s\right)}\right(s)+B_{w}w(s)-{\mathrm{LD}}_{v}v(s\left)\right)---\left(32\right)$对式(32)两边取范数后得到其中则其中γ(2)＝cθ(2)+d₁，R(2)为此时的残差最大值；以此递推，当i＝n‑1时$\begin{array}{c}{\mathrm{\Delta f}}_{k+1}(n-1)=(I-{\mathrm{\Gamma CB}}_f){\mathrm{\Delta f}}_k(n-1)-{\mathrm{\Gamma CB}}_{w}w(n-1)\\ -\Gamma C(A-LC)e_{K\left(n\right)}(n-1)-\Gamma (CL-I)D_{v}v(n-1)\end{array}---\left(35\right)$同样对式(35)两边取范数后可得$\begin{array}{c}||{\mathrm{\Delta f}}_{k+1}(n-1)||\le ||I-{\mathrm{\Gamma CB}}_f||||{\mathrm{\Delta f}}_k(n-1)||+||\Gamma C(A-LC)||||e_{K(n-1)}(n-1)||\\ +||{\mathrm{\Gamma CB}}_{w}w(n-1)||+||\Gamma (CL-I)D_{v}v(n-1)||\le \rho ||{\mathrm{\Delta f}}_k(n-1)||+a_3+lE(n-1)\end{array}---\left(36\right)$进一步展开后得到$\begin{array}{c}||{\mathrm{\Delta f}}_{k+1}(n-1)||\le \rho ||{\mathrm{\Delta f}}_k(n-1)||+a_3+lE(n-1)\\ \le {\rho }^2||{\mathrm{\Delta f}}_{k-1}(n-1)||+\rho (a_3+lE(n-1\left)\right)+a_3+lE(n-1)\\ \le {\rho }^3||{\mathrm{\Delta f}}_{k-2}(n-1)||+({\rho }^2+\rho +1)(a_3+lE(n-1\left)\right)\\ \le \mathrm{...}\\ \le {\rho }^k||{\mathrm{\Delta f}}_1(n-1)||+({\rho }^{k+1}+\mathrm{...}+{\rho }^2+\rho +1)(a_3+lE(n-1\left)\right)\\ ={\rho }^k||{\mathrm{\Delta f}}_1(n-1)||+\frac{1-{\rho }^k}{1-\rho }(a_3+lE(n-1\left)\right)\end{array}---\left(37\right)$由于0＜ρ＜1,则$||{\mathrm{\Delta f}}_k(n-1)||\to \frac{1-{\rho }^k}{1-\rho }(a_3+lE(n-1\left)\right),k\to \infty ---\left(38\right)$此时必然存在0＜K(n‑1)＜∞，ε＞0使得k＝K(n‑1)，而且$||{\mathrm{\Delta f}}_k(n-1)||\le \frac{1-{\rho }^k}{1-\rho }(a_3+lE(n-1\left)\right)+ϵ=\alpha (n-1)+\eta (n-1)ϵ=\Delta F(n-1)---\left(39\right)$其中结合式(15)可得$e_{K(n-1)}\left(n\right)=\underset{s=0}{\overset{n-1}{\Sigma }}\Phi (n,s+1)\left(B_f{\mathrm{\Delta f}}_{K\left(s\right)}\right(s)+B_{w}w(s)-{\mathrm{LD}}_{v}v(s\left)\right)---\left(40\right)$对式(40)两边取范数后得到其中此时系统的残差为其中γ(n)＝cθ(n)+d₁，R(n)是最终的残差值；因此，最终系统的故障估计误差收敛于序列ΔF(0)，ΔF(1)，...，ΔF(n‑1)，输出误差收敛于序列R(1)，R(2)，...R(n)，若ε→0，则ΔF(i)→{α(n)}，R(i)→{γ(n)}，i＝0,1,2…n；在收敛性分析里加入ε项，是为了给故障诊断过程设定一个平衡诊断时间和精度的指标，从而使诊断算法更加灵活；通常阈值范围若选取得比较小，在诊断过程中就会出现迭代次数过多导致诊断时间较长的情况，此时如果更侧重诊断的快速性则可以将ε值适当增大以减小故障诊断的时间代价；第五步：在迭代学习故障诊断算法中加入初值预估算法本发明基于离散采样点设计迭代学习故障诊断算法，可以根据采样顺序逐点诊断出故障，因此在诊断某采样点故障时不必将虚拟故障初值设置为零，可以充分利用该采样点前面已估计出的故障信息对该点故障进行预估，并将估计值作为其虚拟故障的初值；这样选取虚拟故障初值使得故障预估值会尽量接近实际故障值，最好的情况是预估值和故障实际值误差在所要求范围内，表明此预估值足以反映故障实际值，则此采样点的故障诊断结束，所需迭代次数为0，当然这种情况可能不会经常发生，多数情况下故障预估值和实际故障值有较大偏差；同时由于故障模型未知，所以在预估算法选择上，不能直接利用所有已估计出的故障点；本发明根据滑动平均滤波原理提出一种滑动故障初值预估方法，只选取采样点之前距离最近一个时段内的故障信息对该点的故障进行预测估计；设计初值预估算法如下：$Pn\left(i\right)=polyfit\left(p_H\right(i),{\hat{f}}_H(i),n)---\left(43\right)$${\hat{f}}_0\left(i\right)=polyval\left(Pn\right(i),i)---\left(44\right)$式中p_H(i)＝[i‑H,…,i‑1]是第i个采样点之前的H个采样时间点向量，H是滑动窗口长度，为相对应的H个已估计出的故障值，Pn(i)＝[p_n(i),…,p₀(i)]为第i个采样点之前的由H个采样点故障估计值线性拟合得到的多项式系数向量，n是所用拟合多项式的次数，是第i个采样点的故障初值预估值；因此，本发明进一步将线性拟合与滑动平均滤波方法中滑动窗口的思想相结合，在诊断某一采样点的故障之前，首先对其之前已估计出的H个故障点进行线性拟合，利用拟合出的多项式函数求得第i个故障点的故障值即为第i个采样点的故障初值预估值，然后采用前面算法对该采样点进行故障诊断；此采样点诊断结束后，诊断结果又会作为已知值来预估下一个采样点的故障值，此时滑动窗口加入同时丢掉使窗口长度仍为H；由此可知，窗口将始终保持固定长度，并沿着离散采样序列滑动，当向前滑动一步时，窗口前面将进入一个新的数据，后面会丢弃一个旧数据，滑动窗口中始终存放着与待诊断采样点距离最近的故障估计值；本发明的优点是：以机电控制系统这类在工业中广泛应用的工业设备为研究对象，利用其在实际应用过程中的离散采样特性进行故障诊断；本发明的初值预估迭代学习故障诊断算法可以减小诊断过程中前面离散采样点故障的估计误差对后续离散采样点故障诊断的影响，从而有效减少迭代次数,提高故障诊断效率；故障算法结构简单，不仅可以检测和重构机电控制系统的各类型执行器故障，而且易于工程实现，便于实时在线诊断故障；可进一步推广应用于机械臂，磁盘驱动读取系统等实际工程对象。