面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法，它属于工业过程的数据驱动故障诊断领域。本发明是为解决现有方法无法对含有未知扰动信号的系统进行准确的故障诊断的问题。本发明的故障诊断方法不依赖于工业系统模型，仅利用采集到的系统数据即可进行故障的诊断；且本发明方法对未知扰动具有强鲁棒性，在计算方面负担小，易于在线实现；本发明研究了含有未知扰动的子空间辨识问题，通过将过程数据投影到不同的子空间中辨识闭环数据驱动稳定核表示，获得对未知扰动鲁棒的残差，从而实现更好的故障诊断。本发明可以应用于对复杂工业系统进行故障诊断。
1.面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法，其特征在于，所述方法具体包括以下步骤：
步骤一、测量工业系统在正常工况下的输入信号u以及输出信号y，根据给定的样本组数N，过去步长sp和未来步长sf构造汉克尔矩阵Up,N、Uf,N、Yp,N和Yf,N；
其中，Up,N是用过去步长sp构建的输入信号汉克尔矩阵，Yp,N是用过去步长sp构建的输出信号汉克尔矩阵，Uf,N是用未来步长sf构建的输入信号汉克尔矩阵，Yf,N是用未来步长sf构建的输出信号汉克尔矩阵；
步骤二、根据步骤一中构造的汉克尔矩阵Up,N、Uf,N、Yp,N和Yf,N，计算中间变量矩阵Mf,N和Zp,N；
所述中间变量矩阵Zp,N为：其中，上角标T代表矩阵的转置；
步骤三、将未知扰动的子空间描述为Db,N，对矩阵 做LQ分解；
步骤四、根据LQ分解结果计算中间变量矩阵 和所述步骤四的具体过程为：
其中，L31、L32、L33、L41、L42和L43均为进行LQ分解得到的下三角矩阵中的子块；
步骤五、根据中间变量矩阵 和 确定数据驱动闭环稳定核表示步骤六、根据 计算残差rd,k；
步骤七、利用卡方分布表确定χα，并基于χα设置阈值Jth；
其中，α为常数，prob{χ＞χα}表示χ＞χα的概率，则设置阈值Jth为：Jth＝χα/2
步骤八、采集工业控制系统在线数据进行监测设置窗口宽度为kw，分别计算每个窗口内的残差和rd,win(w)：并基于各个窗口内的残差和，计算在线统计量：其中，w表示第w个窗口，Σr为协方差矩阵；
步骤九、判断系统是否发生故障若J≤Jth，则说明工业系统无故障，若J＞Jth，则说明工业系统有故障。
2.根据权利要求1所述的面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法，其特征在于，所述步骤一的具体过程为：其中，l表示输入信号的维数，m表示输出信号的维数，R为实数， 表示矩阵Up,N的行数为spl列数为N， 表示矩阵Uf,N的行数为sfl列数为N， 表示矩阵Yp,N的行数为spm列数为N， 表示矩阵Yf,N的行数为sfm列数为N；
us,k和ys,k分别为：其中，uk为输入信号在第k个采样点处的数据， 为输入信号在第k+sf‑1个采样点处的数据，yk为输出信号在第k个采样点处的数据， 为输出信号在第k+sf‑1个采样点处的数据。
3.根据权利要求2所述的面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法，其特征在于，所述步骤三的具体过程为：其中，L11、L21、L22和L44均为进行LQ分解得到的下三角矩阵中的子块，Q1、Q2、Q3和Q4为酉矩阵的子块。
4.根据权利要求3所述的面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法，其特征在于，所述步骤五的具体过程为：其中，
5.根据权利要求4所述的面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法，其特征在于，所述步骤六的具体过程为：
面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法
技术领域
[0001] 本发明属于工业过程的数据驱动故障诊断领域，具体涉及一种面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法。
背景技术
[0002] 为保证复杂工业系统的安全性和可靠性，在系统发生故障后，需要及时精准地诊断出故障。但复杂工业系统的数学模型通常难以精准地建立，且在运行过程中通常受到环
境的影响，系统中往往存在扰动，这种扰动信号通常是未知且不可测量的，因此需要对含有
未知扰动信号的系统进行准确的故障诊断，保证对未知扰动的鲁棒性，这对于进一步研究
复杂工业系统的故障诊断带来了一定的挑战。
发明内容
[0003] 本发明的目的是为解决现有方法无法对含有未知扰动信号的系统进行准确的故障诊断的问题，提出一种面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法。
[0004] 本发明为解决上述技术问题所采取的技术方案是：
[0005] 面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法，所述方法具体包括以下步骤：
[0006] 步骤一、测量工业系统在正常工况下的输入信号u以及输出信号y，根据给定的样本组数N，过去步长sp和未来步长sf构造汉克尔矩阵Up,N、Uf,N、Yp,N和Yf,N；
[0007] 其中，Up,N是用过去步长sp构建的输入信号汉克尔矩阵，Yp,N是用过去步长sp构建的输出信号汉克尔矩阵，Uf,N是用未来步长sf构建的输入信号汉克尔矩阵，Yf,N是用未来步长sf
构建的输出信号汉克尔矩阵；
[0008] 步骤二、根据步骤一中构造的汉克尔矩阵Up,N、Uf,N、Yp,N和Yf,N，计算中间变量矩阵Mf,N和Zp,N；
[0009] 步骤三、将未知扰动的子空间描述为Db,N，对矩阵 做LQ分解；
[0010] 步骤四、根据LQ分解结果计算中间变量矩阵 和
[0011] 步骤五、根据中间变量矩阵 和 确定数据驱动闭环稳定核表示
[0012] 步骤六、根据 计算残差rd,k；
[0013] 步骤七、利用卡方分布表确定χα，并基于χα设置阈值Jth；
[0014] 步骤八、采集工业控制系统在线数据进行监测
[0015] 设置窗口宽度为kw，分别计算每个窗口内的残差和rd,win(w)：
[0016]
[0017] 并基于各个窗口内的残差和，计算在线统计量：
[0018]
[0019] 其中，w表示第w个窗口，Σr为协方差矩阵；
[0020] 步骤九、判断系统是否发生故障
[0021] 若J≤Jth，则说明工业系统无故障，若J＞Jth，则说明工业系统有故障。
[0022] 本发明的有益效果是：本发明提出了一种面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法，本发明设计的故障诊断方法不依赖于工业系统模型，仅利用采集到
的系统数据即可进行故障的诊断，避免因模型不精确影响故障诊断效果；本发明设计的故
障诊断方法对未知扰动具有强鲁棒性，在计算方面负担小，易于在线实现；本发明研究了含
有未知扰动的子空间辨识问题，通过将过程数据投影到不同的子空间中辨识闭环数据驱动
稳定核表示，获得对未知扰动鲁棒的残差，从而实现更好的故障诊断。
附图说明
[0023] 图1是本发明的一种面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法的流程图；
[0024] 图2是本发明实施例的残差曲线图；
[0025] 图3为本发明实施例的评价指标结果图。
具体实施方式
[0026] 具体实施方式一、结合图1说明本实施方式。本实施方式所述的一种面对含有未知扰动的工业系统的数据驱动鲁棒故障诊断方法，所述方法具体包括以下步骤：
[0027] 步骤一、测量工业系统在正常工况下的输入信号u以及输出信号y，根据给定的样本组数N，过去步长sp和未来步长sf构造汉克尔矩阵Up,N、Uf,N、Yp,N和Yf,N；
[0028] 其中，Up,N是用过去步长sp构建的输入信号汉克尔矩阵，Yp,N是用过去步长sp构建的输出信号汉克尔矩阵，Uf,N是用未来步长sf构建的输入信号汉克尔矩阵，Yf,N是用未来步长sf
构建的输出信号汉克尔矩阵；
[0029] 步骤二、根据步骤一中构造的汉克尔矩阵Up,N、Uf,N、Yp,N和Yf,N，计算中间变量矩阵Mf,N和Zp,N；
[0030] 步骤三、将未知扰动的子空间描述为Db,N，对矩阵 做LQ分解；
[0031] 步骤四、根据LQ分解结果计算中间变量矩阵 和
[0032] 步骤五、根据中间变量矩阵 和 确定数据驱动闭环稳定核表示
[0033] 步骤六、根据 计算残差rd,k；
[0034] 步骤七、利用卡方分布表确定χα，并基于χα设置阈值Jth；
[0035] 步骤八、采集工业控制系统在线数据进行监测
[0036] 设置窗口宽度为kw，分别计算每个窗口内的残差和rd,win(w)：
[0037]
[0038] 并基于各个窗口内的残差和，计算在线统计量：
[0039]
[0040] 其中，w表示第w个窗口，Σr为协方差矩阵；
[0041] 步骤九、判断系统是否发生故障
[0042] 若J≤Jth，则说明工业系统无故障，若J＞Jth，则说明工业系统有故障。
[0043] 具体实施方式二：本实施方式与具体实施方式一不同的是，所述步骤一的具体过程为：
[0044]
[0045]
[0046]
[0047]
[0048] 其中，l表示输入信号的维数，m表示输出信号的维数，R为实数， 表示矩阵Up,N的行数为spl列数为N， 表示矩阵Uf,N的行数为sfl列数为N， 表示矩阵Yp,N的
行数为spm列数为N， 表示矩阵Yf,N的行数为sfm列数为N；
[0049] us,k和ys，k分别为：
[0050]
[0051] 其中，uk为输入信号在第k个采样点处的数据， 为输入信号在第k+sf‑1个采样点处的数据，yk为输出信号在第k个采样点处的数据， 为输出信号在第k+sf‑1个采
样点处的数据。
[0052] 参数 和 的定义同理。
[0053] 其它步骤及参数与具体实施方式一相同。
[0054] 具体实施方式三：本实施方式与具体实施方式一至二之一不同的是，所述中间变量矩阵Mf,N为：
[0055]
[0056] 其中，中间变量矩阵
[0057] 其它步骤及参数与具体实施方式一或二相同。
[0058] 具体实施方式四：本实施方式与具体实施方式一至三之一不同的是，所述中间变量矩阵Zp,N为：
[0059]
[0060] 其中，上角标T代表矩阵的转置。
[0061] 其它步骤及参数与具体实施方式一至三之一相同。
[0062] 具体实施方式五：本实施方式与具体实施方式一至四之一不同的是，所述步骤三的具体过程为：
[0063]
[0064] 其中，L11、L21、L22、L31、L32、L33、L41、L42、L43和L44均为进行LQ分解得到的下三角矩阵中的子块，Q1、Q2、Q3和Q4为酉矩阵的子块。
[0065] 其它步骤及参数与具体实施方式一至四之一相同。
[0066] 具体实施方式六：本实施方式与具体实施方式一至五之一不同的是，所述步骤四的具体过程为：
[0067]
[0068] 其它步骤及参数与具体实施方式一至五之一相同。
[0069] 具体实施方式七：本实施方式与具体实施方式一至六之一不同的是，所述步骤五的具体过程为：
[0070]
[0071] 其中，
[0072] 其它步骤及参数与具体实施方式一至七之一相同。
[0073] 具体实施方式八：本实施方式与具体实施方式一至七之一不同的是，所述步骤六的具体过程为：
[0074]
[0075] 其中，∈k为噪 声，且满足分布 Σr为协方差矩阵 ，
[0076] 其它步骤及参数与具体实施方式一至七之一相同。
[0077] 具体实施方式九：本实施方式与具体实施方式一至八之一不同的是，所述步骤七的具体过程为：
[0078] prib{χ＞χα}＝α
[0079] 其中，α为常数(可以根据残差rd,k人为设定)，prob{χ＞χα}表示χ＞χα的概率，则设置阈值Jth为：
[0080] Jth＝χα/2
[0081] 其它步骤及参数与具体实施方式一至八之一相同。
[0082] 实施例
[0083] 以下将结合具体实验结果来说明本发明方法的有效性。
[0084] 本实验中选取的系统为电机，电机的扰动选取周期为10s，幅值为0.1A的矩形波负载。电机转速为3000r/min，采样周期为0.01s，采用比例积分控制器。在0s‑2000s，电机处于
正常运行状态，在第2000s时，电机发生设备执行器故障，脉冲宽度调制(PWM)的增益下降
20％。在阈值设定上，选取α＝2％，窗口宽度kw设置为5。
[0085] 步骤一、测量电机正常工况下的输入信号(即控制器的输出信号)u，系统的输出信号y。并对于给定的样本组数N＝3000，过去步长sp＝1和未来步长sf＝250构造下列几个汉克
尔矩阵，
[0086]
[0087]
[0088]
[0089]
[0090] 其中，k为采样点，Db,N为未知矩形波的子空间描述。us,k和ys,k分别为：
[0091]
[0092] 步骤二、计算矩阵Mf,N和Zp,N，
[0093]
[0094]
[0095] 其中， Ac、Bc、Cc、Dc为控制器的系统矩阵。
[0096] 步骤三、对汉克尔矩阵做LQ分解，
[0097]
[0098] 步骤四、通过下面式子得到矩阵 和
[0099]
[0100] 步骤五、确定闭环数据驱动稳定核表示
[0101]
[0102] 其中，
[0103] 步骤六、计算残差rd,k，
[0104]
[0105] 其中 ，∈k 为噪 声且满 足分 布 Σr为 协方差 矩阵 ，
[0106] 步骤七、设置阈值Jth。自定义α＝2％，利用卡方分布表确定χα，
[0107] prob{χ＞χα}＝2％，
[0108] 其中，prob{χ＞χα}表示χ＞χα的概率。
[0109] 此时，确定阈值Jth为
[0110] Jth＝χα/2
[0111] 步骤八、采集电机在线数据进行监测。设置窗口宽度kw，这里选取kw＝5，计算每个窗口内的残差之和。
[0112]
[0113] 并计算在线统计量，
[0114]
[0115] 其中，w表示第w个窗口。
[0116] 步骤九、判断系统是否发生故障。如果J≤Jth则说明无故障，J＞Jth则说明有故障。
[0117] 本发明提出的含有未知扰动的数据驱动鲁棒故障诊断方法的仿真结果如图2和图3所示，图2为残差曲线图，图3为评价指标结果图。
[0118] 本发明的上述算例仅为详细地说明本发明的计算模型和计算流程，而并非是对本发明的实施方式的限定。对于所属领域的普通技术人员来说，在上述说明的基础上还可以
做出其它不同形式的变化或变动，这里无法对所有的实施方式予以穷举，凡是属于本发明
的技术方案所引伸出的显而易见的变化或变动仍处于本发明的保护范围之列。
