CN109143995B - 质量相关慢特征分解的闭环***精细运行状态监测方法 - Google Patents

质量相关慢特征分解的闭环***精细运行状态监测方法 Download PDF

Info

Publication number
CN109143995B
CN109143995B CN201810771664.6A CN201810771664A CN109143995B CN 109143995 B CN109143995 B CN 109143995B CN 201810771664 A CN201810771664 A CN 201810771664A CN 109143995 B CN109143995 B CN 109143995B
Authority
CN
China
Prior art keywords
quality
matrix
variable
dynamic
space
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201810771664.6A
Other languages
English (en)
Other versions
CN109143995A (zh
Inventor
赵春晖
郑嘉乐
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhejiang University ZJU
Original Assignee
Zhejiang University ZJU
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhejiang University ZJU filed Critical Zhejiang University ZJU
Priority to CN201810771664.6A priority Critical patent/CN109143995B/zh
Publication of CN109143995A publication Critical patent/CN109143995A/zh
Application granted granted Critical
Publication of CN109143995B publication Critical patent/CN109143995B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B19/00—Programme-control systems
    • G05B19/02—Programme-control systems electric
    • G05B19/418—Total factory control, i.e. centrally controlling a plurality of machines, e.g. direct or distributed numerical control [DNC], flexible manufacturing systems [FMS], integrated manufacturing systems [IMS] or computer integrated manufacturing [CIM]
    • G05B19/41875—Total factory control, i.e. centrally controlling a plurality of machines, e.g. direct or distributed numerical control [DNC], flexible manufacturing systems [FMS], integrated manufacturing systems [IMS] or computer integrated manufacturing [CIM] characterised by quality surveillance of production
    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B19/00—Programme-control systems
    • G05B19/02—Programme-control systems electric
    • G05B19/418—Total factory control, i.e. centrally controlling a plurality of machines, e.g. direct or distributed numerical control [DNC], flexible manufacturing systems [FMS], integrated manufacturing systems [IMS] or computer integrated manufacturing [CIM]
    • G05B19/41885—Total factory control, i.e. centrally controlling a plurality of machines, e.g. direct or distributed numerical control [DNC], flexible manufacturing systems [FMS], integrated manufacturing systems [IMS] or computer integrated manufacturing [CIM] characterised by modeling, simulation of the manufacturing system
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02P—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN THE PRODUCTION OR PROCESSING OF GOODS
    • Y02P90/00—Enabling technologies with a potential contribution to greenhouse gas [GHG] emissions mitigation
    • Y02P90/02—Total factory control, e.g. smart factories, flexible manufacturing systems [FMS] or integrated manufacturing systems [IMS]

Landscapes

  • Engineering & Computer Science (AREA)
  • Manufacturing & Machinery (AREA)
  • General Engineering & Computer Science (AREA)
  • Quality & Reliability (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Automation & Control Theory (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)
  • Feedback Control In General (AREA)

Abstract

本发明公开了一种基于质量相关慢特征充分分解的闭环***精细运行状态监测方法。本方法通过提取过程变量与质量变量的静态相关性，将闭环大规模过程变量分解成与质量相关和与质量无关的两个空间，并在每个空间中提取过程的静态和动态信息对其进行建模，分别在两个空间内建立静态指标和动态指标，针对两个空间进行协同的故障检测。该方法不仅可以实现对闭环***动态异常情况和工况变化的监测，有效识别闭环故障***中真正故障的发生，并且通过对变量空间进行质量相关分解，可以进一步实现对闭环***产品质量的相关监测，判断出故障和异常具体发生在哪个空间，是否影响产品质量。

Description

质量相关慢特征分解的闭环***精细运行状态监测方法
技术领域
本发明属于工业***过程监测领域，特别是针对闭环***产品质量相关的精细运行状态监测方法。
背景技术
随着科学技术的发展，现代工业***的规模和复杂程度都在日益提高。复杂大***一旦出现异常，都可能带来重大的财产损失和人员伤亡。通过对***设备日常的运行状态监测，提前检测出设备异常，使设备得到及时的维修，不仅使得检修周期延长，维修费用大大降低，而且检修更有针对性，设备使用寿命得到延长，避免了严重安全事故的发生。因此，为了保证***在运行的安全性与可靠性，及时发现***运行中的异常情况并进行处理，减少生产中的安全隐患，提高设备使用周期，十分有必要采用有效手段对***进行实时监测和故障检测。同时随着传感技术的发展，在工业现场获得数据变得越来越容易，过程数据中蕴含了大量的过程信息，基于数据的状态监测和故障监测逐渐成为研究的热点。
在过去的几十年时间里，过程监测和故障检测技术得到了广泛研究和发展，大量研究成果得到发表，前人对基于数据的故障检测和故障诊断作出了相应的研究。主成分分析(PCA)，偏最小二乘(PLS)和费舍尔判别分析(FDA)等多元统计分析方法已经被广泛应用于基于数据的过程监测领域。然而现有的大部分的过程监测和故障检测方法都是针对开环***设计的，没有考虑到闭环控制规律的影响。但是在真实工业***中，为满足稳定性、快速性、准确性要求，达到实际生产目标保证产品质量维持在合理的水平，需要对***施加闭环反馈控制，大量闭环控制律，如PID控制、最优控制等在实际生产***中有着大量的应用。闭环控制律的引入，使得输入输出关系、***变量之间关系发生改变，并且闭环控制规律使得***对于外部干扰更具鲁棒性，使得故障处于早期阶段或幅值较小时，所带来的影响在***中被淹没，无法通过残差信号及时诊断，造成误报和漏报。这些限制了之前研究的诊断方法在实际中得应用，因此有必要在闭环前提下研究***状态监测和故障检测方法，以解决实际的生产安全性问题。目前基于数据的闭环***监测与故障检测研究成果十分有限，还处于初步的探索阶段，有待进一步深入研究。
本发明针对闭环***的过程监测与故障检测，并从工业***产品质量的角度出发，提出了一种结合慢特征分析和典型相关性分析的闭环***质量相关精细运行状态监测的算法。本方法通过提取过程数据与产品质量变量的相关性，实现对过程数据空间进行充分分解，分成与质量密切相关和与质量无关的两个空间，并建立模型对两个空间的动态和静态信息进行监测，该方法通过对动态和静态信息的分别监测可以实现对闭环***工况变化和真正故障发生的区分，并且通过质量相关的空间分解可以识别出异常具体发生在哪一个空间是否影响产品质量，提高了对***状态监测的精准，到目前为止，尚未见到与本发明相关的研究报道。
发明内容
本发明的目的在于针对工业闭环***过程故障检测，并且关注过程产品质量，提供了一种基于质量相关慢特征充分分解的闭环***精细运行状态监测方法。
本发明的目的是通过以下技术方案实现的：一种基于质量相关慢特征充分分解的闭环***精细运行状态监测方法，包含以下步骤：
(1)选择***正常过程测量数据与所关注的产品质量数据：对于一个工业闭环***，设其运行过程包含m个可测过程变量，在t时刻采样可以得到一个1×m的向量x(t)＝(x1(t),x2(t),…,xm(t))，经过N次采样后得到正常过程下过程测量变量的数据矩阵XM＝(x(t),x(t+1),…,x(t+N))T，根据工业过程选择k个所关注的产品质量变量，经N次采样后，得到质量变量的数据矩阵Yk(∈N*k)。
(2)数据标准化：分别对过程数据矩阵XM和质量变量数据矩阵Yk按列减去该列的均值，并除以该列标准差进行标准化。
(3)对(2)中标准化后的过程数据矩阵X进行慢特征分析建模，获得初始慢特征。
这里考虑线性慢特征分析，每一个慢特征都可以看作是过程变量的线性组合，所以从过程数据矩阵X到慢特征S(t)＝[s1(t),s2(t),…,sm(t)]的映射表示成：
S(t)＝W×X (24)
其中W＝[w1,w2,…,wm]T是SFA需要优化的系数矩阵，m为慢特征的个数，与可测过程变量相同。求解慢特征的问题是使得所获得慢特征变化速率Δ(si)最小，求解慢特征问题转化成求解下列广义特征值问题：
AW＝BWΩ (25)
其中A表示输入的过程数据矩阵X的一阶差分矩阵
Figure GDA0002570881790000031
的协方差矩阵，表示成
Figure GDA0002570881790000032
B表示X矩阵的协方差矩阵<XXT>t，对应的，W＝[w1,w2,…,wm]T为特征向量构成的特征矩阵，m过程输入数据的变量个数；Ω为对应广义特征值构成的对角矩阵。利用两步奇异值分解(SVD)来求解慢特征，具体过程如下：
(3.1)首先对矩阵X进行白化，即对B矩阵的奇异值分解。矩阵X的协方差矩阵可表示成<XXT>t，对其进行奇异值分解(SVD)得：
<XXT>t＝UΛUT (26)
其中U是协方差矩阵<XXT>的特征向量组成的矩阵，Λ是一个对角阵，其每一个对角线上的元素就是一个特征值，则白化后得数据矩阵可表示成：
Z＝Λ-1/2UTX＝QX (27)
白化后数据矩阵Z得协方差满足<ZZT>t＝Q<XXT>tQT。
(3.2)求解线性SFA问题的目标等同于寻找一个矩阵P＝WQ-1使得S＝P*Z，并且使得S满足<SST>t＝P<ZZT>tPT＝I。
矩阵P可以通过对Z的差分矩阵
Figure GDA0002570881790000036
的协方差矩阵
Figure GDA0002570881790000033
进行奇异值分解求得：
Figure GDA0002570881790000034
由此可以计算得到初始慢特征系数矩阵为：
W＝PΛ-1/2UT (29)
初始慢特征为：
S＝WX＝PΛ-1/2UTX (30)
(4)利用式(30)中获得的初始慢特征S与质量变量Yk进行典型相关性分析(CCA)获取典型变量，具体过程如下：
对于初始慢特征S与质量变量Yk，定义u，v两个变量，其中u，v是始慢特征S与质量变量Yk的线性组合如下：
Figure GDA0002570881790000035
其中Ψ1＝[a1,a2,…,am]T，Ψ2＝[b1,b2,…,bk]T。求解典型变量u，v就是求解系数矩阵Ψ1，Ψ2使得u，v间的Pearson系数最大。转化成如下优化问题求解：
Figure GDA0002570881790000041
Subjectto：
Figure GDA0002570881790000042
其中corr(u,v)表示u，v之间的Pearson系数，cov(u,v)为u，v的协方差，Var(u)和Var(v)分别为u，v的方差。此问题通过构造Lagrangian等式来求解，求解对应最大Pearson系数的两个典型变量u，v。
同理，在此基础上继续求解寻找对应第二大Pearson系数的典型变量u2和v2：
Figure GDA0002570881790000043
一共可以求取n对典型变量对，n＝min(m,k)，对应的Pearson系数依次减小，即表示典型变量间的相关性越来越小。根据所求得的Pearson系数选取前p个对应Pearson系数较大的典型相关变量对，舍去剩下对应Pearson系数较小的(n-p)个典型相关变量对。利用前p个典型相关变量对，利用所求得关于慢特征S的线性组合的典型变量ui(i＝1,2,…,p)，组合成新的矩阵Uy＝[u1,u2,…,up]，其中ui(i＝1,2,…,p)对应的Pearson系数依次减小，并且ui之间相互正交。
(5)将(4)中所获得的矩阵Uy利用偏最小二乘公式重构过程数据X空间，重构公式如下：
Figure GDA0002570881790000044
Eo为重构残差，即与Uy无关的部分。由此可将X空间分成与Yk相关的部分Uy和与Yk无关的部分Xo两个部分，实现对X空间的质量充分分解。
X＝Uy+Xo (35)
(6)重新对Xo进行慢特征分析，得到慢特征，并按变化速度从小到大的顺序排列：
So＝WoXo＝[So1,So2,…,Som] (36)
Xo是原始过程数据减去质量变量Yk相关部分数据所获得的，所以获得的慢特征不是全部都有意义，选取满足
Figure GDA0002570881790000045
的g个慢特征，E{}为求期望运算符，I是单位矩阵；
根据所获得慢特征的变化快慢程度，q个变化速度较慢特征组成矩阵Sod＝[so1,so2,…,soq]，剩下(g-q)个变化速度较快的特征组成另一矩阵：Soe＝[so(q+1),so(q+2),…,sog]。
(7)建立静态和动态指标体系。
考虑到对于实际工业***,所关注的产品质量变量可能为一个或者两个以上,对于两种情况,建立不同的指标体系对***进行监测.
(7.1)对于单个质量变量的情况
由于质量变量只有一个,所以所提取到的典型变量也只有一个,所以对于与质量相关的空间直接利用典型变量Uy进行监测。其中静态指标为：
Figure GDA0002570881790000051
动态指标为：
Figure GDA0002570881790000052
对于动静态指标均采用shewhart控制图的方法进行监测，对于静态指标：
Figure GDA0002570881790000053
其中UTH1为
Figure GDA0002570881790000054
的控制上限，LTH1为
Figure GDA0002570881790000055
的控制下限，μ1为
Figure GDA0002570881790000056
的均值，
Figure GDA0002570881790000057
为其标准差，b1为控制图门限宽度。对于动态指标：
Figure GDA0002570881790000058
其中UTH2为
Figure GDA0002570881790000059
的控制上限，LTH2为
Figure GDA00025708817900000510
的控制下限，μ2为
Figure GDA00025708817900000511
的均值，
Figure GDA00025708817900000512
为其标准差，b2为控制图门限宽度。
对于与质量变量无关的空间，静态指标为：
Figure GDA00025708817900000513
Figure GDA00025708817900000514
动态指标为：
Figure GDA00025708817900000515
Figure GDA0002570881790000061
Figure GDA0002570881790000062
表示Sod的协方差矩阵特征值所构成的对角阵的逆，
Figure GDA0002570881790000063
表示Soe的协方差矩阵特征值所构成的对角阵的逆，指标的控制限均通过核密度函数在0.95的置信水平下获得，分别为
Figure GDA0002570881790000064
分别为
Figure GDA0002570881790000065
协方差矩阵的特征值构成的对角阵的逆，两个动态指标所对应的控制限均通过核密度函数在0.95的置信水平下计算得到，分别为
Figure GDA0002570881790000066
(7.2)对于有两个以上质量变量的情况
所获得的典型变量并不唯一，所以需要计算统计量来对分Uy进行监测，静态指标计算如下:
Figure GDA0002570881790000067
其中
Figure GDA0002570881790000068
表示Uy的协方差矩阵对应的特征值构成的对角阵的逆，指标的控制限通过核密度函数在0.95的置信水平下获得分别为
Figure GDA0002570881790000069
对Uy进行一次差分得到
Figure GDA00025708817900000610
动态指标为：
Figure GDA00025708817900000611
Figure GDA00025708817900000612
为
Figure GDA00025708817900000613
协方差矩阵的特征值构成的对角阵的逆，
Figure GDA00025708817900000614
为其三个动态指标所对应的控制限，通过核密度函数在0.95的置信水平下计算得到。与质量变量无关的空间的指标计算方式与步骤7.1中单个质量变量的情况相同。
(8)利用动静态指标可以分别对闭环***的过程动态性和稳态操作情况进行监测。
对于与质量变量相关的空间:
1)静态指标
Figure GDA00025708817900000615
超限,动态指标
Figure GDA00025708817900000616
超限后又恢复到控制限以下,表明与质量变量相关的空间经过控制作用调节到达一个新的工况,质量变量出现短暂波动后又恢复到设定水平；
2)静态指标
Figure GDA00025708817900000617
超限,动态指标
Figure GDA00025708817900000618
也一直超限,表明与质量变量相关的空间出现故障,控制器调节失败,产品质量遭到破坏；
3)静态指标
Figure GDA00025708817900000619
动态指标
Figure GDA00025708817900000620
均未超限,则表明与质量变量相关的空间处于良好的运行状态,产品质量维持在设定值附近。
对于与质量变量无关的空间:
1)质量无关的静态指标
Figure GDA0002570881790000071
与动态指标
Figure GDA0002570881790000072
均短暂超限后回到控制限以下,表明在与质量变量无关的空间出现短暂异常,后通过闭环控制作用的调节,该空间又恢复到与异常发生前相同的稳态工况下运行；
2)与质量无关空间的静态指标
Figure GDA0002570881790000073
一直超限,动态指标
Figure GDA0002570881790000074
短暂超限后回到控制限以下,表明与质量变量无关的空间出现短暂异常,但是通过闭环控制作用的调节该空间恢复到一个异常发生前并不一样的稳态工况下运行,静态性能已经发生改变,但是动态性能仍维持正常；
3)质量无关的静态指标
Figure GDA0002570881790000075
与动态指标
Figure GDA0002570881790000076
均超限,表明该空间出现故障,且该空间相关控制器调节失败,静动态性能都遭到破坏；
4)质量无关的静态指标
Figure GDA0002570881790000077
与动态指标
Figure GDA0002570881790000078
均维持正常,表明该空间保持一个良好的运行状态。
本发明的有益效果为：本发明可以将工业闭环***划分成与产品质量相关和与产品质量无关的两个空间，不仅可以有效区分故障具体发生在哪一个空间，进而判断故障是否影响产品质量，而且通过对***动静态信息的监测，可以有效地识别工况的变化与真正故障的区分，实现闭环***的精细运行状态监测。
附图说明
图1是本发明的流程图；
图2是本方法具体应用的闭环TE过程工艺流程图；
图3是TE过程提供的故障类展示图；
图4是本方法应用在TE过程故障1情况下的静态指标监测结果；
图5是本方法应用在TE过程故障1情况下的动态指标监测结果；
图6是TE过程故障1情况下对应的产品质量A与正常情况下产品质量A对比图；
图7是TE过程故障1情况下对应的产品质量B与正常情况下产品质量B对比图；
图8是TE过程故障1情况下对应的产品质量C与正常情况下产品质量C对比图；
图9是本方法应用在TE过程故障4情况下的静态指标监测结果；
图10是本方法应用在TE过程故障4情况下的动态指标监测结果；
图11是TE过程故障4情况下对应的产品质量A与正常情况下产品质量A对比图；
图12是TE过程故障4情况下对应的产品质量A与正常情况下产品质量A对比图；
图13是TE过程故障4情况下对应的产品质量A与正常情况下产品质量A对比图。
具体实施方式
下面结合附图及具体实例，对本发明作进一步详细说明。
Tenessee Eastman(TE)过程是由Downs等人基于Tenessee Eastman化学公司某实际化工生产过程提出的一个仿真***,在过程***工程领域的研究中,TE过程是一个常用的标准问题(Benchmark problem),其较好地模拟了实际复杂工业过程***的许多典型特征,因此被作为仿真例子广泛应用于控制、优化、过程监控与故障诊断的研究中。本研究中，采用的是Lyman和Georgakis提出的全厂控制策略的闭环TE过程。TE过程包含四种气体原料A、C、D和E,两种液态产物G和H,还包含副产物F和惰性气体B。
TE过程包括五个主要单元：反应器、冷凝器、压缩机、分离器和汽提塔，如图2所示，共包括41个测量变量和12个控制变量。如图3，TE过程提供了21个故障类,可供使用。
如图1所示，本发明包括以下步骤：
(1)数据选择:选择TE过程正常运行时的过程测量数据XMEAS(1-22)和XMV(1-11)共33个变量作为过程输入变量X，XMEAS(29-31)三个变量所关注的产品量变量Y，每次采样可以得到一个向量x(t)＝(x1(t),x2(t),…,xm(t))维度为1×m的，经过N次采样后得到过程测量变量X＝(x(t),x(t+1),…,x(t+N))T，产品质量变量，经N次采样后，得到质量变量的数据矩阵Yk(∈N*k)。
(2)数据标准化：分别对过程数据矩阵X和质量变量数据矩阵Yk按列去均值并除以标准差进行标准化。
(3)对(2)中标准化后的过程数据矩阵X进行慢特征分析建模，获得初始慢特征。
这里考虑线性慢特征分析，每一个慢特征都可以看作是过程变量的线性组合，所以从过程数据矩阵到慢特征S(t)＝[s1(t),s2(t),…,sm(t)]的映射可以表示成：
S(t)＝W×X
其中W＝[w1,w2,…,wm]T是SFA需要优化的系数矩阵，m为慢特征的个数，与过程变量数相同。求解慢特征的问题是使得所获得慢特征变化速率Δ(si)最小，求解慢特征问题可以转化成求解下列广义特征值问题：
AW＝BWΩ
其中A表示输入数据X的一阶差分矩阵
Figure GDA0002570881790000091
的协方差矩阵，可表示成
Figure GDA0002570881790000092
B表示X矩阵的协方差矩阵<XXT>t，W＝[w1,w2,…,wm]T为特征向量构成的特征矩阵，m为过程输入数据的变量个数，Ω为对应广义特征值构成的对角矩阵，同样也是。同样可以利用两步奇异值分解(SVD)来求解慢特征，具体过程如下：
(3.1)首先对矩阵X进行白化，即对B矩阵的奇异值分解。矩阵X的协方差矩阵可表示成<XXT>t，对其进行奇异值分解(SVD)得：
<XXT>t＝UΛUT
其中U是协方差矩阵<XXT>t的特征向量组成的矩阵，Λ是一个对角阵，其每一个对角线上的元素就是一个特征值，则白化后得数据矩阵可表示成：
Z＝Λ-1/2UTX＝QX
白化后数据矩阵Z得协方差满足<ZZT>t＝Q<XXT>tQT。
(4.2)求解线性SFA问题的目标等同于寻找一个矩阵P＝WQ-1使得S＝P*Z，并且使得S满足<SST>t＝P<ZZT>tPT＝I。
矩阵P可以通过对Z的差分矩阵
Figure GDA0002570881790000093
的协方差矩阵
Figure GDA0002570881790000094
进行奇异值分解求得：
Figure GDA0002570881790000095
由此可以计算得到初始慢特征系数矩阵为：
W＝PΛ-1/2UT
初始慢特征为：
S＝WX＝PΛ-1/2UTX
(5)利用步骤(4)中获得的初始慢特征S与质量变量Yk进行典型相关性分析(CCA)获取典型变量，具体过程如下：
(5.1)对于初始慢特征S与质量变量Yk，定义u，v两个变量，其中u，v是始慢特征S与质量变量Yk的线性组合如下：
Figure GDA0002570881790000101
求解典型变量u，v就是求解系数矩阵A，B使得u，v间的Pearson系数最大。
可转化成如下优化问题求解：
Maximizeu2
Subjectto：
Figure GDA0002570881790000102
其中corr(u,v)表示Pearson系数，cov(u,v)为u，v的协方差，Var(u)和Var(v)分别u，v的方差。此问题可以通过构造Lagrangian等式求解，求解对应最大Pearson系数的两个典型变量u，v。
(5.1)同理，可以在此基础上继续求解寻找对应第二大Pearson系数的典型变量u2和v2：
Figure GDA0002570881790000103
一共可以求取n对典型变量对，n＝min(m,k)，对应的Pearson系数依次减小，即典型变量间的相关性越来越小。根据需要一般根据所求得的Pearson系数选取前p个对应最大相关系数的典型相关变量对，舍去剩下对应Pearson系数很小的(n-p)个典型相关变量对。利用前p个典型相关变量对，利用对应慢特征S的线性组合的典型变量组合成新的矩阵Uy＝[u1,u2,…,up]，ui(i＝1,2,…,p)为求得的典型变量，对应的Pearson系数依次减小，并且ui之间互不相关。
(6)将(5)中所获得的典型变量Uy利用偏最小二乘公式重构回Xd空间，重构公式如下：
X＝UyPT+E
PT＝(Uy TUy)-1Uy TXd
并将X空间分成与Yk相关的部分Uy和与Yk无关的部分Xo两个部分，实现对X空间的质量充分分解。
X＝Uy+Xo
(7)重新对Xo进行慢特征分析，得到慢特征，并按变化速度从小到大的顺序排列：
So＝WoXo＝[So1,So2,…,Som]
Xo是原始过程数据减去质量变量Yk相关部分数据所获得的，所以获得的慢特征不是全部都有意义，选取满足
Figure GDA0002570881790000111
的g个慢特征。并根据所获得慢特征的变化快慢程度，选取q个变化速度较慢特征组成Sod＝[so1,so2,…,soq]，剩下(g-q)个变化速度较快的特征Soe＝[so(q+1),so(q+2),…,sog]。
(8)建立静态和动态指标体系。
考虑到对于实际工业***,所关注的产品质量变量可能为一个或者多个,对于两种情况,建立不同的指标体系对***进行监测.
(8.1)对于单个质量变量的情况
由于质量变量只有一个,所以所提取到的典型变量也只有一个,所以对于与质量相关的空间直接利用典型变量Uy进行监测。其中静态指标为：
Figure GDA0002570881790000112
动态指标为：
Figure GDA0002570881790000113
对于动静态指标均采用shewhart控制图的方法进行监测，以静态指标为例：
Figure GDA0002570881790000114
其中UTH为
Figure GDA0002570881790000115
的控制上限，LTH为控制下限，μ为
Figure GDA0002570881790000116
的均值，
Figure GDA0002570881790000117
为其标准差，b为控制图门限宽度。
对于与质量变量无关的空间，静态指标为：
Figure GDA0002570881790000118
Figure GDA0002570881790000119
动态指标为：
Figure GDA00025708817900001110
Figure GDA00025708817900001111
Figure GDA00025708817900001112
表示Sod的协方差矩阵特征值所构成的对角阵的逆，
Figure GDA00025708817900001113
表示Soe的协方差矩阵特征值所构成的对角阵的逆，指标的控制限均通过核密度函数在0.95的置信水平下获得分别为
Figure GDA0002570881790000121
分别为
Figure GDA0002570881790000122
协方差矩阵的特征值构成的对角阵的逆，
Figure GDA0002570881790000123
分别为两个动态指标所对应的控制限，通过核密度函数在0.95的置信水平下计算得到。
(8.2)对于有多个质量变量的情况
所获得的典型变量并不唯一，所以需要计算统计量来对分Uy进行监测，静态指标计算如下:
Figure GDA0002570881790000124
其中
Figure GDA0002570881790000125
表示Uy的协方差矩阵对应的特征值构成的对角阵的逆，指标的控制限均通过核密度函数在0.95的置信水平下获得分别为
Figure GDA0002570881790000126
对Uy进行一次差分得到
Figure GDA0002570881790000127
动态指标为：
Figure GDA0002570881790000128
Figure GDA0002570881790000129
为
Figure GDA00025708817900001210
协方差矩阵的特征值构成的对角阵的逆，
Figure GDA00025708817900001211
为其三个动态指标所对应的控制限，通过核密度函数在0.95的置信水平下计算得到,
与质量变量无关的空间的指标计算方式同上。
利用动静态指标可以分别对闭环***的过程动态性和稳态操作情况进行监测。对于与质量变量相关的空间:
1)静态指标
Figure GDA00025708817900001212
超限,动态指标
Figure GDA00025708817900001213
超限后又恢复到控制限以下,表明与质量变量相关的空间经过控制作用调节到达一个新的工况,质量变量出现短暂波动后又恢复到设定水平；
2)静态指标
Figure GDA00025708817900001214
超限,动态指标
Figure GDA00025708817900001215
也一直超限,表明与质量变量相关的空间出现故障,控制器调节失败,产品质量遭到破坏；
3)静态指标
Figure GDA00025708817900001216
动态指标
Figure GDA00025708817900001217
均未超限,则表明与质量变量相关的空间处于良好的运行状态,产品质量维持再设定值附近。
对于与质量变量无关的空间:
1)质量无关的静态指标
Figure GDA00025708817900001218
与动态指标
Figure GDA00025708817900001219
均短暂超限后回到控制限以下,表明在与质量变量无关的空间出现短暂异常,后通过闭环控制作用的调节,该空间又恢复到与异常发生前相同的稳态工况下运行；
2)与质量无关空间的静态指标
Figure GDA0002570881790000131
一直超限,动态指标
Figure GDA0002570881790000132
短暂超限后回到控制限以下,表明与质量变量无关的空间出现短暂异常,但是通过闭环控制作用的调节该空间恢复到一个异常发生前并不一样的稳态工况下运行,静态性能已经发生改变,但是动态性能仍维持正常；
3)质量无关的静态指标
Figure GDA0002570881790000133
与动态指标
Figure GDA0002570881790000134
均超限,表明该空间出现故障,且该空间相关控制器调节失败,静动态性能都遭到破坏；
4)质量无关的静态指标
Figure GDA0002570881790000135
与动态指标
Figure GDA0002570881790000136
均维持正常,表明该空间保持一个良好的运行状态。
两个空间相互独立,相应的动静态指标分别指示各自空间的动静态信息的变化,可以帮助我们判断故障及异常具体发生在哪个空间,是否会影响产品质量。
图4、5本方法用于TE过程故障1数据的监测结果，图6、7、8是所选择的三个产品质量变量正常情况下变化和在故障1情况下的变化对比。根据图可以看出，与产品质量y相关的空间，静态指标和动态指标超限一段时间后均恢复到控制限以下，说明该空间受到故障的短暂干扰，通过控制器调节一段时间后静动态特性均恢复到初始水平，产品质量也恢复到初始设定值；与产品质量无关的空间，静态指标一直超限，但趋于平稳，动态指标超限一段时间后恢复到控制限以下，说明该空间受到干扰后通过调节恢复到一个新的稳态工况下。图9、10是本方法用于TE过程故障4数据的监测结果，从图中可以看出与产品质量y相关的空间静动态指标均未超限，说明该空间运行正常，产品质量未产生变化，从图11、12、13中各个产品质量正常情况下和故障四情况下的变化对比也可以看出监测结果的判断是正确的；在与产品质量无关的空间，静态指标一直超限，动态指标短暂超限后恢复到控制限以下，说明故障发生在与产品质量无关的闭环中，通过控制器的调节又达到一个新的稳态工况。传统的慢特征分析只能区分工况异常和故障发生，而通过本方法可以进一步区分故障所发生的空间以及是否影响产品质量。

Claims (1)

1.一种基于质量相关慢特征充分分解的闭环***精细运行状态监测方法，其特征在于，包含以下步骤：
(1)选择***正常过程测量数据与所关注的产品质量数据：对于一个工业闭环***，设其运行过程包含m个可测过程变量，在t时刻采样可以得到一个1×m的向量x(t)＝(x1(t),x2(t),…,xm(t))，经过N次采样后得到正常过程下过程变量的数据矩阵XM＝(x(t),x(t+1),…,x(t+N))T，根据工业过程选择k个所关注的产品质量变量，经N次采样后，得到质量变量的数据矩阵Yk∈N*k；
(2)数据标准化：分别对过程数据矩阵XM和质量变量数据矩阵Yk按列减去该列的均值，并除以该列标准差进行标准化；
(3)对步骤(2)中标准化后的过程数据矩阵X进行慢特征分析建模，获得初始慢特征；
这里考虑线性慢特征分析，每一个慢特征都可以看作是过程变量的线性组合，所以从过程数据矩阵X到慢特征S(t)＝[s1(t),s2(t),…,sm(t)]的映射表示成：
S(t)＝W×X (1)
其中W＝[w1,w2,…,wm]T是SFA需要优化的系数矩阵，m为慢特征的个数，与可测过程变量相同；求解慢特征的问题是使得所获得慢特征变化速率Δ(si)最小，求解慢特征问题转化成求解下列广义特征值问题：
AW＝BWΩ (2)
其中A表示输入的过程数据矩阵X的一阶差分矩阵
Figure FDA0002503452140000011
的协方差矩阵，表示成
Figure FDA0002503452140000012
B表示X矩阵的协方差矩阵<XXT>t，对应的，W＝[w1,w2,…,wm]T为特征向量构成的特征矩阵，m过程输入数据的变量个数；Ω为对应广义特征值构成的对角矩阵；利用两步奇异值分解SVD来求解慢特征，具体过程如下：
(3.1)首先对矩阵X进行白化，即对B矩阵的奇异值分解；矩阵X的协方差矩阵可表示成<XXT>，对其进行奇异值分解SVD得：
<XXT>t＝UΛUT (3)
其中U是协方差矩阵<XXT>t的特征向量组成的矩阵，Λ是一个对角阵，其每一个对角线上的元素就是一个特征值，则白化后得数据矩阵可表示成：
Z＝Λ-1/2UTX＝QX (4)
白化后数据矩阵Z得协方差满足<ZZT>t＝Q<XXT>tQT；
(3.2)求解线性SFA问题的目标等同于寻找一个矩阵P＝WQ-1使得S＝P*Z，并且使得S满足<SST>t＝P<ZZT>tPT＝I；
矩阵P可以通过对Z的差分矩阵
Figure FDA0002503452140000021
的协方差矩阵
Figure FDA0002503452140000022
进行奇异值分解求得：
Figure FDA0002503452140000023
由此可以计算得到初始慢特征系数矩阵为：
W＝PΛ-1/2UT (6)
初始慢特征为：
S＝WX＝PΛ-1/2UTX (7)
(4)利用(7)中获得的初始慢特征S与质量变量Yk进行典型相关性分析CCA获取典型变量，具体过程如下：
对于初始慢特征S与质量变量Yk，定义u，v两个变量，其中u，v是初始慢特征S与质量变量Yk的线性组合如下：
Figure FDA0002503452140000024
其中Ψ1＝[a1,a2,…,am]T，Ψ2＝[b1,b2,…,bk]T；求解典型变量u，v就是求解系数矩阵Ψ1，Ψ2使得u，v间的Pearson系数最大；转化成如下优化问题求解：
Figure FDA0002503452140000025
Figure FDA0002503452140000026
其中corr(u,v)表示u，v之间的Pearson系数，cov(u,v)为u，v的协方差，Var(u)和Var(v)分别为u，v的方差；此问题通过构造Lagrangian等式来求解，求解对应最大Pearson系数的两个典型变量u，v；
同理，在此基础上继续求解寻找对应第二大Pearson系数的典型变量u2和v2：
Figure FDA0002503452140000027
一共可以求取n对典型变量对，n＝min(m,k)，对应的Pearson系数依次减小，即表示典型变量间的相关性越来越小；根据所求得的Pearson系数选取前p个对应Pearson系数较大的典型相关变量对，舍去剩下对应Pearson系数较小的n-p个典型相关变量对；利用前p个典型相关变量对，利用所求得关于慢特征S的线性组合的典型变量ui，其中i＝1,2,…,p，组合成新的矩阵Uy＝[u1,u2,…,up]，其中ui(i＝1,2,…,p)对应的Pearson系数依次减小，并且ui之间相互正交；
(5)将步骤(4)中所获得的矩阵Uy利用偏最小二乘公式重构过程数据X空间，重构公式如下：
Figure FDA0002503452140000031
Eo为重构残差，即与Uy无关的部分；由此可将X空间分成与Yk相关的部分Uy和与Yk无关的部分Xo两个部分，实现对X空间的质量充分分解；
X＝Uy+Xo (12)
(6)重新对Xo进行慢特征分析，得到慢特征，并按变化速度从小到大的顺序排列：
So＝WoXo＝[So1,So2,…,Som] (13)
Xo是原始过程数据减去质量变量Yk相关部分数据所获得的，所以获得的慢特征不是全部都有意义，选取满足
Figure FDA0002503452140000032
的g个慢特征，E{}为求期望运算符，I是单位矩阵；
根据所获得慢特征的变化快慢程度，q个变化速度较慢特征组成矩阵Sod＝[so1,so2,…,soq]，剩下g-q个变化速度较快的特征组成另一矩阵：Soe＝[so(q+1),so(q+2),…,sog]；
(7)建立静态和动态指标体系；
考虑到对于实际工业***,所关注的产品质量变量可能为一个或者两个以上,对于两种情况,建立不同的指标体系对***进行监测；
(7.1)对于单个质量变量的情况
由于质量变量只有一个,所以所提取到的典型变量也只有一个,所以对于与质量相关的空间直接利用典型变量Uy进行监测；其中静态指标为：
Figure FDA0002503452140000033
动态指标为：
Figure FDA0002503452140000041
对于动静态指标均采用shewhart控制图的方法进行监测，对于静态指标：
Figure FDA0002503452140000042
其中UTH1为
Figure FDA0002503452140000043
的控制上限，LTH1为
Figure FDA0002503452140000044
的控制下限，μ1为
Figure FDA0002503452140000045
的均值，
Figure FDA0002503452140000046
为其标准差，b1为控制图门限宽度；对于动态指标：
Figure FDA0002503452140000047
其中UTH2为
Figure FDA0002503452140000048
的控制上限，LTH2为
Figure FDA0002503452140000049
的控制下限，μ2为
Figure FDA00025034521400000410
的均值，
Figure FDA00025034521400000411
为其标准差，b2为控制图门限宽度；
对于与质量变量无关的空间，静态指标为：
Figure FDA00025034521400000412
Figure FDA00025034521400000413
动态指标为：
Figure FDA00025034521400000414
Figure FDA00025034521400000415
Figure FDA00025034521400000416
表示Sod的协方差矩阵特征值所构成的对角阵的逆，
Figure FDA00025034521400000417
表示Soe的协方差矩阵特征值所构成的对角阵的逆，指标的控制限均通过核密度函数在0.95的置信水平下获得，分别为
Figure FDA00025034521400000418
Figure FDA00025034521400000419
分别为
Figure FDA00025034521400000420
协方差矩阵的特征值构成的对角阵的逆，两个动态指标所对应的控制限均通过核密度函数在0.95的置信水平下计算得到，分别为
Figure FDA00025034521400000421
(7.2)对于有两个以上质量变量的情况
所获得的典型变量并不唯一，所以需要计算统计量来对分Uy进行监测，静态指标计算如下:
Figure FDA00025034521400000422
其中
Figure FDA00025034521400000423
表示Uy的协方差矩阵对应的特征值构成的对角阵的逆，指标的控制限通过核密度函数在0.95的置信水平下获得分别为
Figure FDA0002503452140000051
对Uy进行一次差分得到
Figure FDA0002503452140000052
动态指标为：
Figure FDA0002503452140000053
Figure FDA0002503452140000054
为
Figure FDA0002503452140000055
协方差矩阵的特征值构成的对角阵的逆，
Figure FDA0002503452140000056
为其三个动态指标所对应的控制限，通过核密度函数在0.95的置信水平下计算得到；与质量变量无关的空间的指标计算方式与步骤7.1中单个质量变量的情况相同；
(8)利用动静态指标可以分别对闭环***的过程动态性和稳态操作情况进行监测；
对于与质量变量相关的空间:
1)静态指标
Figure FDA0002503452140000057
超限,动态指标
Figure FDA0002503452140000058
超限后又恢复到控制限以下,表明与质量变量相关的空间经过控制作用调节到达一个新的工况,质量变量出现短暂波动后又恢复到设定水平；
2)静态指标
Figure FDA0002503452140000059
超限,动态指标
Figure FDA00025034521400000510
也一直超限,表明与质量变量相关的空间出现故障,控制器调节失败,产品质量遭到破坏；
3)静态指标
Figure FDA00025034521400000511
动态指标
Figure FDA00025034521400000512
均未超限,则表明与质量变量相关的空间处于良好的运行状态,产品质量维持在设定值附近；
对于与质量变量无关的空间:
1)质量无关的静态指标
Figure FDA00025034521400000513
与动态指标
Figure FDA00025034521400000514
均短暂超限后回到控制限以下,表明在与质量变量无关的空间出现短暂异常,后通过闭环控制作用的调节,该空间又恢复到与异常发生前相同的稳态工况下运行；
2)与质量无关空间的静态指标
Figure FDA00025034521400000515
一直超限,动态指标
Figure FDA00025034521400000516
短暂超限后回到控制限以下,表明与质量变量无关的空间出现短暂异常,但是通过闭环控制作用的调节该空间恢复到一个异常发生前并不一样的稳态工况下运行,静态性能已经发生改变,但是动态性能仍维持正常；
3)质量无关的静态指标
Figure FDA00025034521400000517
与动态指标
Figure FDA00025034521400000518
均超限,表明该空间出现故障,且该空间相关控制器调节失败,静动态性能都遭到破坏；
4)质量无关的静态指标
Figure FDA00025034521400000519
与动态指标
Figure FDA00025034521400000520
均维持正常,表明该空间保持一个良好的运行状态。
CN201810771664.6A 2018-07-13 2018-07-13 质量相关慢特征分解的闭环***精细运行状态监测方法 Active CN109143995B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810771664.6A CN109143995B (zh) 2018-07-13 2018-07-13 质量相关慢特征分解的闭环***精细运行状态监测方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810771664.6A CN109143995B (zh) 2018-07-13 2018-07-13 质量相关慢特征分解的闭环***精细运行状态监测方法

Publications (2)

Publication Number Publication Date
CN109143995A CN109143995A (zh) 2019-01-04
CN109143995B true CN109143995B (zh) 2020-09-01

Family

ID=64800694

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810771664.6A Active CN109143995B (zh) 2018-07-13 2018-07-13 质量相关慢特征分解的闭环***精细运行状态监测方法

Country Status (1)

Country Link
CN (1) CN109143995B (zh)

Families Citing this family (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110222393B (zh) * 2019-05-28 2020-12-18 浙江大学 基于细粒度风电发电状态划分的风机叶片结冰异常监测方法
CN111401940B (zh) * 2020-03-05 2023-07-04 杭州网易再顾科技有限公司 特征预测方法、装置、电子设备及存储介质
CN111413949B (zh) * 2020-03-30 2023-12-01 南京富岛信息工程有限公司 一种降低工业过程故障预警误报率的方法
CN111624979B (zh) * 2020-05-18 2021-07-06 浙江大学 一种基于慢特征分析的工业闭环控制回路多振荡检测与溯源方法
CN111949003B (zh) * 2020-07-17 2021-09-03 浙江浙能技术研究院有限公司 一种基于SFA与Hellinger距离的闭环控制回路性能评价方法
CN112199409B (zh) * 2020-08-17 2024-06-11 浙江中控技术股份有限公司 一种用于催化重整装置实时工况的监测方法及装置
CN112541017A (zh) * 2020-12-02 2021-03-23 中海油信息科技有限公司 一种工业生产过程状态监测方法、装置、设备及存储介质
CN117048113B (zh) * 2023-10-13 2023-12-22 江苏爱箔乐铝箔制品有限公司 铝箔餐盒冲压的动态精度检测方法及***

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104077608A (zh) * 2014-06-11 2014-10-01 华南理工大学 基于稀疏编码慢特征函数的行为识别方法
CN104537260A (zh) * 2015-01-14 2015-04-22 清华大学 基于缓慢特征回归的动态软测量方法和***
CN105678723A (zh) * 2015-12-29 2016-06-15 内蒙古科技大学 基于稀疏分解和差分图像的多聚焦图像融合方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104077608A (zh) * 2014-06-11 2014-10-01 华南理工大学 基于稀疏编码慢特征函数的行为识别方法
CN104537260A (zh) * 2015-01-14 2015-04-22 清华大学 基于缓慢特征回归的动态软测量方法和***
CN105678723A (zh) * 2015-12-29 2016-06-15 内蒙古科技大学 基于稀疏分解和差分图像的多聚焦图像融合方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
Concurrent monitoring of operating condition deviations and process dynamics anomalies with slow feature analysis;Shang, C.等;《AIChE Journal》;20151231;第3666-3682页 *
基于核慢特征回归与互信息的常压塔软测量建模;蒋昕祎 等;《化工学报》;20170531;第68卷(第5期);第1977-1986页 *

Also Published As

Publication number Publication date
CN109143995A (zh) 2019-01-04

Similar Documents

Publication Publication Date Title
CN109143995B (zh) 质量相关慢特征分解的闭环***精细运行状态监测方法
CN101169623B (zh) 基于核主元分析贡献图的非线性过程故障辨识方法
CN108062565B (zh) 基于化工te过程的双主元-动态核主元分析故障诊断方法
CN104714537B (zh) 一种基于联合相对变化分析和自回归模型的故障预测方法
Jiang et al. Weighted kernel principal component analysis based on probability density estimation and moving window and its application in nonlinear chemical process monitoring
US20140067327A1 (en) Similarity curve-based equipment fault early detection and operation optimization methodology and system
Ralston et al. Computer-based monitoring and fault diagnosis: a chemical process case study
CN111368428B (zh) 一种基于监控二阶统计量的传感器精度下降故障检测方法
CN115994337B (zh) 一种带钢热连轧非平稳过程微小故障检测方法及装置
CN110942258B (zh) 一种性能驱动的工业过程异常监测方法
Xu et al. A novel kernel dynamic inner slow feature analysis method for dynamic nonlinear process concurrent monitoring of operating point deviations and process dynamics anomalies
Linzhe et al. A nonlinear quality-relevant process monitoring method with kernel input-output canonical variate analysis
Li et al. Dynamic reconstruction principal component analysis for process monitoring and fault detection in the cold rolling industry
CN109283912B (zh) 一种面向智能电厂大型燃煤发电机组制粉***的分布式动静协同综合监测方法
Shaikh et al. Data-driven based fault diagnosis using principal component analysis
CN114995338A (zh) 一种基于规范变量分析与js散度融合的工业过程微小故障检测方法
Hanyuan et al. Batch process monitoring based on batch dynamic Kernel slow feature analysis
Galicia et al. A comprehensive evaluation of Statistics Pattern Analysis based process monitoring
Alcala et al. Monitoring of dynamic processes with subspace identification and principal component analysis
Ma et al. Process monitoring of the pneumatic control valve using canonical variate analysis
Ji et al. Orthogonal projection based statistical feature extraction for continuous process monitoring
Qi et al. Data-driven fault diagnosis and prognosis for process faults using principal component analysis and extreme learning machine
CN110928263B (zh) 一种预先考虑动态关系的复杂过程的故障检测方法及***
Li et al. Fault detection and diagnosis based on new ensemble kernel principal component analysis
Zhang et al. A modified PCA-based approach for process monitoring

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant