王乾 王聪

非线性系统特征提取是智能系统、模式识别、故障诊断和生物医学工程等领域的基础和关键问题.在实际对象中,系统产生的信号本质上都是非线性非平稳的,如何从这些非线性非平稳信号中得到有效的特征描述是研究者关心的主要研究方向之一[1−2].早期的信号处理方法如时域分析(均值、方差等)、频域分析(傅里叶分析),虽然简单方便,但是只对平稳和线性的信号有作用[3−4].随着技术的发展,一些适用于非线性非平稳信号的时频分析方法被提出[5−7],例如短时傅里叶变换、Wigner-Ville分布、小波变换、Hilbert-Huang变换、局部均值分解等.一些研究者同时提出了能用于刻画系统非线性程度的动力学不变量特征指标[8−10],例如Lyapunov指数、熵及LZ复杂度等.

尽管上述非线性系统特征提取方法取得了一定进展,但它们是基于系统的状态轨迹[11]对系统进行特征提取和分析,提取的特征容易丢失系统内在的动态信息,导致很难识别出系统潜在的、早期的微弱信号特征,不能对系统产生的非线性非平稳信号进行充分的特征表达[12−14].因此,如何从原始系统信号中得到既包含系统状态信息又包含系统动态信息的数据,并对这些数据进行特征表达是一个新的研究重点和难点.对此,确定学习理论提供了一个有前途的方向.

确定学习理论可用于对非平稳环境下的时变或动力学模式进行准确建模和快速识别[15−16].确定学习运用自适应控制和非线性动力学系统的概念与方法,研究未知动态环境下的知识获取、表达、存储和再利用等问题.针对产生周期或回归轨迹的非线性动力学系统,通过选择局部径向基函数(Radial basis function,RBF)神经网络作为参数化的模型结构,确定学习理论可以对其未知的系统动态进行局部准确建模/辨识,将系统状态轨迹带入动力学建模结果中,可得到系统的动力学轨迹,以时不变且空间分布的方式进行表达,并以常值存储在RBF神经网络中.这种表达方式包含了系统状态轨迹和系统内在动力学的全部信息[17−20].

本文从动力学轨迹这一新角度对非线性系统进行特征提取和分析.1)上述非线性系统特征提取方法中,由于LZ复杂度是一种非概率的复杂性测量方法且较容易实施[21−24],本文使用LZ复杂度对上述系统动力学轨迹进行特征表达.LZ复杂度反映了一个数据序列随着序列长度的增加出现新模式的速率,可定量分析复杂数据序列的有序性.通常,复杂度值越大代表越复杂、无序、不规则的动态系统,复杂度值越小代表动态系统越规则、有序.在LZ复杂度算法的基础上,本文提出时空LZ复杂度.时空LZ复杂度包含时间LZ复杂度(Temporal-LZ complexity,TLZC)和空间LZ复杂度(Saptio-LZ complexity,SLZC)两个指标,指标值的大小和数据序列的复杂程度成正比例关系[25].因此,时空LZ复杂度可以在时域和空间域上对系统动力学轨迹的复杂程度进行分析.2)将系统的动力学轨迹用时空LZ复杂度特征指标定量表达后,可对其敏感性进行分析.敏感性分析能反映出非线性动力学系统的状态或行为相对于系统参数变化的敏感程度.Savageau提出参数敏感度可作为评价和比较生化系统性能的标准[26].Wu等阐述了局部敏感度[27]、全局敏感度[28]和目标敏感度[29]的概念,并给出了具体的计算过程,其中目标敏感度表示系统的某一性能指标对于系统参数变化的敏感程度[30].由于可将时空LZ复杂度特征指标作为系统不同参数/状态的目标函数,本文选择目标敏感度对动态特征提取方法进行敏感性分析,定量评价该方法相对系统状态变化的敏感程度.3)结合经典的非线性动力学系统(Rossler系统)进行数值仿真,结果表明,相比于状态轨迹,从动力学轨迹中提取的特征可以更加敏感地反映出非线性动力学系统的周期、倍周期及混沌状态.4)使用时空LZ复杂度对低速轴流压气机旋转失速过程进行动态特征表达,动力学轨迹的TLZC和SLZC指标能更加明显地显示出系统从失速前进入旋转失速初始扰动的过程.本文提出的动态特征提取方法的优点是从系统内在动态的角度对原系统进行更好的表达.

本文安排如下:第1节为预备知识,简要介绍确定学习理论和Lempel-Ziv复杂度算法;第2节阐述了动态特征提取的过程;第3节对动态特征提取方法进行敏感度分析;第4节和第5节分别通过数值仿真和实验分析验证本文所提方法的有效性;第6节对本文工作进行总结.

1 预备知识

1.1 确定学习理论

确定学习理论运用自适应控制和非线性动力学系统的概念与方法,对产生周期或回归轨迹的非线性动力学系统,其未知系统动态可被局部准确建模/辨识和快速识别[15−16].其基本要素包括:1)使用 RBF神经网络;2)周期或回归的系统轨迹可满足部分持续激励条件;3)在系统轨迹的邻域内实现对非线性系统动态的局部准确神经网络建模/辨识;4)所学的知识以时不变且空间分布的方式进行表达,以常值神经网络权值的方式进行存储.

考虑如下的非线性动态系统

其中,x=[x1,···,xN]T∈RN是可测量的系统状态,p是系统的常值参数向量(通常不同的p可以产生不同的系统动态行为),F(x;p)=[f1(x;p),···,fN(x;p)]T表示光滑未知的系统动态,fi(x;p)是未知的连续非线性函数.始于初值x0的系统轨迹记为φζ(x0),假设系统(1)状态x保持一致有界,如x(t)∈Ω∈RN,其中Ω是一个紧集,且系统轨迹φζ(x0)是回归轨迹[21].

确定学习理论采用如下动态RBF神经网络对非线性系统(1)的未知系统动态F(x;p)进行辨识:

文献[16]指出,对于非线性动力学系统(1)产生的周期轨迹或更一般的回归轨迹,结合辨识模型(2)和权值更新律(3),RBF神经网络中沿着回归轨迹的神经元函数构成的子向量可以满足部分PE条件,即靠近轨迹φζ(x0)的神经网络回归向量Sζ满足PE条件.这个部分PE条件可以使得辨识误差系统满足指数稳定,进而在沿周期或回归轨迹的局部区域实现对非线性系统动态的准确神经网络逼近.

基于确定学习机制[16],通过式(4)对非线性系统的未知动态进行局部准确辨识后,可将随时间变化的动态模式以时不变且空间分布的方式有效地表达.这种表达方式是一种包含全部状态和动态信息的全息表达方法[26−27].

1.2 Lempel-Ziv复杂度算法

复杂度算法最初由Kolmogorov提出[31],但是并没有给出具体的计算过程.随后 Lempel和 Ziv在数学上给出了严谨的证明和可行的计算过程,使得复杂度算法得以在实际中应用,在以后的文献中把这种复杂度算法统一称为 Lempel-Ziv复杂度(LZ复杂度)[27].LZ复杂度反映了一个数据序列随着序列长度的增加出现新模式的速率,其值越大代表数据序列越复杂.在计算LZ复杂度之前,任何一个数据序列y(i)(i=1,···,n,n=length(y(i)))需要转换成0-1二值符号序列u(i)(i=1,···,n).

其中,yave表示数据序列y(i)的平均值.

在复杂度计算过程中,从左到右对符号序列u(i)进行扫描,每次有新子串出现时,其复杂度值c(n)增加1,计算过程参见文献[32−33].

为了使得到的复杂度值适用于不同数据长度,c(n)需要进行归一化处理.如果序列的长度为n,在0-1符号集中不同符号的个数为2,已有文献表明复杂度c(n)的上界为[10,34]

那么,复杂度c(n)可以通过b(n)进行归一化

通常,将归一化的复杂度C记为LZ复杂度的值.

注1.文献[21]证明了当n足够大时,归一化式(7)成立.然而,文献[35−36]表明,如果数据序列比较小时,归一化的复杂度值C可能大于1.为了削弱数据序列长度的限制,文献[25]和文献[34]分别给出了一个经验值n≥3600及n≥6000时,归一化的复杂度值C可在[0,1]区间内.如文献[37−38]所述,本文选取数据序列的长度为n=104,详见附录A.

如文献[21]所述的Lempel-Ziv复杂度算法中粗粒化过程,数据序列y(i)经过二值粗粒化过程被转换为0-1符号序列(Binary coarse-graining LZC,BLZC).然而,文献[39]表明BLZC可能丢失非线性系统一些有用的信息.为了解决这个问题,多种多值粗粒化(Multi-valued coarse-graining LZC,MLZC)过程被陆续提出[40−41].然而,多值粗粒化机制中的粗粒化程度需要外部人工干预,凭借经验进行选择.本文提出一种新的多值粗粒化过程— 极值粗粒化(Extremum coarse-graining LZC,ELZC):找到数据序列y(i)所有的极值点(包括极大值和极小值),进行升序排列记为E(j)(j=1,···,l,l=length(E)),并记β={0,1,···,l−1}代表符号集.那么,非线性动力学数据序列y(i)可通过如下极值粗粒化过程转化为符号序列u(i)(u(i)∈{β(1),β(2),···,β(l)}).

基于上述极值粗粒化机制,复杂度计算过程中的粗粒化过程不仅解决了二值粗粒化可能丢失一些有用信息的问题,其粗粒化程度还可以根据信号自身进行选择,避免了凭借经验进行选择的外部干预影响.二值粗粒化、多值粗粒化和极值粗粒化过程的对比效果见附录A.

2 动态特征提取

本节基于确定学习理论和时空Lempel-Ziv复杂度算法提出动态特征提取新方法.

首先,通过确定学习理论,非线性系统(1)的未知动态F(x,p)可以通过式(4)进行准确地建模/辨识.由于辨识得到的系统动力学以常值的形式存储在RBF神经网络中,因此可将其作为有限的数据序列进行分析.

其次,在LZ复杂度算法的基础上,提出了时空LZ复杂度对非线性动力学数据序列进行分析.时空LZ复杂度包含时间复杂度TLZC和空间复杂度SLZC两个指标,可以从时间域和空间域上刻画系统的复杂程度.

时间复杂度TLZC的计算过程如下[33,35]:

步骤1.非线性动力学数据序列Y=[y1,···,yN]T需要转换为新的符号数据序列U=[u1,···,uN]T.根据第1.2节中的极值粗粒化机制,对非线性动力学数据序列yk(i)(k=1,···,N,i=1,···,n,n=length(yk(i)))找到其所有的极值点(包括极大值和极小值),进行升序排列记为Ek(j)(j=1,···,l,l=length(Ek)),并记βk={0,1,···,l−1}代表符号集.那么,通过式(8)可将数据序列yk(i)转化为符号序列uk(i)(uk(i)∈{βk(1),βk(2),···,βk(l)}).

步骤2.对每个符号序列uk(i),其复杂度ck(n)初始值置为1,并令P={uk(1)},Q={uk(2)},PQ={uk(1),uk(2)}及PQπ=uk(1)(P和Q代表符号序列uk的两个子数据序列,PQ是P和Q的结合,PQπ表示PQ删除最后一个字符,v(PQπ)代表PQπ所有不同的子数据序列).

步骤3.在计算过程中,当P={uk(1),uk(2),···,uk(r)}和Q={uk(r+1)}时,PQπ={uk(1),uk(2),···,uk(r)}. 如果Q属于v(PQπ),更新Q为{uk(r+1),uk(r+2)}.如果不属于,则更新P和Q为P={uk(1),uk(2),···,uk(r+1)},Q=Qk(r+2),同时,其复杂度值加1:ck(n)=ck(n)+1.

步骤4.重复步骤2和步骤3,直到Q为最后一个字符.

步骤5.每个数据序列uk的归一化复杂度可以通过下式获得.

其中,k=1,···,N,n是数据序列yk的长度,lk是第k个符号集中不同字符的个数.

基于上述步骤1～5,非线性系统动力学轨迹的时间复杂度TLZC为

其中,RMS(·)是均方根值函数.TLZC值可以从时域的角度对非线性动力学系统的复杂度进行度量.

空间复杂度SLZC是从系统动力学轨迹数据序列方向导数的角度进行计算.

步骤1.非线性系统动力学数据序列Y的方向导数可以通过下式进行近似计算[42],并记为数据序列Z.

步骤2.每个方向导数序列Zk的归一化复杂度值SCk的计算过程可类似于时间复杂度的步骤1和步骤5.

其中,m=n−1是每个方向导数序列的长度,lk是第k个符号集中不同字符的个数.那么,其相应的空间复杂度SLZC为

方向导数数据序列反映了系统动力学轨迹在空间上的变化速率[42−43],相应的空间复杂度值SLZC可以作为对非线性动力学系统在空间域上的复杂度度量.

3 动态特征提取方法的敏感度分析

使用时空LZ复杂度对系统的动力学轨迹进行特征表达后,可通过敏感度分析评价提出方法的性能.目标敏感度可以反映出系统的某一性能指标相对于系统参数变化的敏感程度[30].由于时间复杂度TLZC和空间复杂度SLZC是系统的复杂性表达,因此可用目标敏感度对时空LZ复杂度性能指标进行分析.数学上,相对目标敏感度系数为[26]

其中,C和p分别是状态轨迹/动力学轨迹复杂度和系统参数,η是系统的目标敏感度系数.相对的意思是对复杂度和系统参数分别进行归一化.由于系统的性能指标和参数值可能在一个很大的区间内变化,所以经常使用归一化的复杂度值对系统的性能进行比较.

在实际计算中,目标敏感度系数可使用有限差分法进行近似[26].

其中,∆C和∆p分别是状态轨迹/动力学轨迹复杂度和系统参数在一个小区间内的差值.为了较合适地比较系统的状态轨迹和动力学轨迹复杂度指标的敏感性,需使用相同的归一化尺度.同时,在比较时使用敏感度系数η的绝对值.|η|的值越大,相应的性能指标对系统参数的变化越敏感.

4 数值仿真

基于上述敏感度分析,使用经典的非线性动力学系统—Rossler系统[44]验证提出方法的有效性.

其中,x=[x1,x2,x3]T∈R3是可测量的系统状态向量,p=[p1,p2,p3]T是系统常值参数向量,f1(x;p)=−x2−x3,f2(x;p)=x1+p1x3,f3(x;p)=p2+x3(x1−p3)是系统的未知动态.第1.1节已表明,系统的未知动态f1(x;p),f2(x;p),f3(x;p)可以被确定学习理论准确地建模/辨识,过程参见文献[15−16].

根据文献[45],固定参数p1=p2=0.2,变化参数p3,Rossler系统(17)可以产生不同的状态.例如,选择参数p3的变化范围为2.6～4.38,步长为0.01,系统将产生倍周期分岔过程:单周期(p3=2.6～2.8),2倍周期(p3=2.81～3.82),4倍周期(p3=3.83～4.12),8 倍周期 (p3=4.13～4.18),···,混沌状态(p3=4.22～4.38).图1显示了Rossler系统状态x1的倍周期分岔过程.

微分方程组(17)可以通过四阶Runge-Kutta进行求解(时间步长设置为0.01,系统状态序列的大小取为104,与注1一致).系统的状态轨迹和相应的动力学轨迹如图2所示.本文选取单周期(p3=2.7),4倍周期(p3=3.96)和混沌状态(p3=4.28)进行对比图示.

图1 Rossler系统状态x1的倍周期分岔过程Fig.1 The period-doubling bifurcation diagram of the statex1of the Rossler system

从图2可以看出,系统的动力学轨迹和状态轨迹的形状相一致.下面使用时空LZ复杂度分别对系统的状态轨迹和动力学轨迹进行分析.其时间复杂度指标TLZC如图3所示,空间复杂度指标SLZC如图4所示.图3和图4中M代表状态轨迹的复杂度指标曲线,代表动力学轨迹的复杂度指标曲线.

在对时间复杂度和空间复杂度两个性能指标进行敏感度分析之前,计算系统不同状态区间的复杂度值算术平均,用该算术平均值表示系统不同状态的复杂度值.系统状态轨迹和动力学轨迹的时空LZ复杂度(TLZC和SLZC)的敏感度系数可根据敏感度系数式(16)得到.为了对系统不同状态变化进行详细的敏感度对比分析,本文分别计算了单周期～2倍周期、2倍周期～4倍周期、4倍周期～8倍周期和8倍周期～混沌状态的敏感度系数,结果如表1所示.η1表示系统状态轨迹的敏感度系数,η2表示系统动力学轨迹的敏感度系数.

表1 Rossler系统的敏感度系数Table 1 The sensitivity coefficients of the Rossler system

图2 Rossler系统的状态轨迹和动力学轨迹图Fig.2 The state trajectory and dynamics trajectory of the Rossler system

图3 系统的状态轨迹和动力学轨迹时间复杂度指标图Fig.3 The TLZC indices of state trajectory and dynamics trajectory of the Rossler system

图4 系统的状态轨迹和动力学轨迹空间复杂度指标图Fig.4 The SLZC indices of state trajectory and dynamics trajectory of the Rossler system

从表1可以看出,在时域和空间域上,动力学轨迹的复杂度性能指标相比于状态轨迹的复杂度指标,对系统状态变化时的敏感度系数整体上都要大一些.当系统从单周期到2倍周期变化时,动力学轨迹的TLZC和SLZC指标敏感度系数分别为0.0229和0.0268,相比于状态轨迹的TLZC和SLZC指标敏感度系数0.0018和0.0247有所增加,虽然不很明显,但也能区分出系统的不同状态,与图3和图4中的复杂度指标曲线相一致.随着系统状态的不断变化,系统动力学轨迹的复杂度指标敏感度系数越来越大,特别是系统从8倍周期状态变化到混沌状态时,动力学轨迹的TLZC和SLZC指标敏感度系数(2.2767,1.3192)明显大于状态轨迹的TLZC和SLZC指标敏感度系数(1.4589,0.9603).从上述分析可以得出,基于动力学轨迹的复杂度特征表达可以更敏感地反映出系统的周期、倍周期和混沌状态.

除了使用时空LZ复杂度在时域和空间域上对非线性动力学系统进行复杂程度分析,本文使用时频分析方法(Hilbert-Huang变换)对上述Rossler系统(17)进行分析,以对比系统状态和系统动态在频谱上的分布情况.Hilbert-Huang变换方法于1998年提出,用于分析非线性非平稳的信号,主要包含两个过程:经验模态分解(Empirical mode decomposition,EMD)和Hilbert谱分析[1].对信号进行Hilbert-Huang变换后,得到的结果包含模态分解图、Hilbert谱图和Hilbert边际谱图.模态分解图表示原始信号经 EMD得到的本征模态函数(Intrinsic mode function,IMF),即各个频率成份;Hilbert谱图表示信号能量在频率和时间轴上的分布情况;Hilbert边际谱图表示信号的幅值在整个频率段内随频率的变化情况.选取上述Rossler系统(17)的2倍周期为例,使用Hilbert-Huang变换对其系统状态x1和系统动态f1(x,3.3)进行分析,得到的模态分解图、Hilbert谱图和Hilbert边际谱图分别如图5～7所示.

图5 Rossler系统2倍周期模态分解图Fig.5 The EMD of period-2 of Rossler system

图6 Rossler系统2倍周期Hilbert谱图Fig.6 The Hilbert spectrum of period-2 of Rossler system

图7 Rossler系统2倍周期Hilbert边际谱图Fig.7 The Hilbert marginal spectrum of period-2 of Rossler system

从图5可以看出,Hilbert-Huang变换可以对系统状态和动态进行很好地分解,都包含一个高频成份、一个低频成份和一个均值非零的低强度剩余信号,即分解误差.

从图6的Hilbert谱图可以看出,系统状态的谱图(图6(a))和系统动态的谱图(图6(b))频率成份是一致的,都包含两个主要的频率成份.高频成份从0.09Hz～0.15Hz之间变化,低频成份主要集中在0.055Hz.但是,对照系统动态的模态分解图5(b)和Hilbert谱图6(b),分解得到的系统动态的高频成份和低频成份相比于系统状态的模态分解图5(a)和谱图6(a),其周期性要弱一些,表明系统内在动态变化会较状态更明显,这与系统的Hilbert边际谱图相对应,如图7所示.系统状态和动态的边际谱图可以很清晰显示出系统的能量或频率成份主要分布在0.055Hz的低频成份和0.1Hz附近的高频成份上.但是,相比于系统状态的边际谱图7(a),系统动态的边际谱图7(b)的高频成份分布有稍微大一些的波动,能量分布平坦一些,说明系统内在动态的变化较明显.

通过上述分析,从信号的模态分解图、Hilbert谱图和边际谱图可以看出,系统状态和系统动态在频谱上的分布虽然有微小差别,但是基本相一致.频谱分析方法和复杂度表征方法是从不同的角度对系统的非线性非平稳信号进行分析.频谱分析方法给出的是信号的频谱分布,频谱分布主要是将信号从时域转换到频域,对其频率特性和能量分布进行分析;而LZ复杂度特征分析给出的是非线性系统确定的复杂程度,并且可以从时间和空间的角度上,对系统状态和动态的复杂度特征进行明显区分.因此,结合图3和图4可得,LZ复杂度分析方法可以敏感地反映出系统的不同状态.相对于频谱分析,复杂度量化指标对系统的特征表达更为直观,并且可以在更大程度上对系统的不同状态进行区分.

5 实验分析

基于文中提出的动态特征提取方法,对涡轮风扇发动机中的旋转失速过程进行动态特征表达,从实验分析上验证提出方法的可行性和有效性.

涡轮风扇发动机是目前世界上军用和大型民用飞机最常用的动力装置,最主要的特点是在高亚音速/超音速飞行条件下具有很高的效率.轴流压气机是涡扇发动机的核心部件之一,当代大型航空涡扇发动机的发展方向是追求更高的单级压比和更少的级数.单级压比的提高,引起流动分离,导致发动机内部产生不稳定流动问题,例如旋转失速和喘振等.喘振和旋转失速都是发动机内流的系统性失稳,它们限定了发动机的稳定工作区域,在系统喘振发生的最初,总是伴有压气机的旋转失速,旋转失速被认为是喘振先兆.旋转失速检测已经成为压气机研究领域中重要且困难的问题.轴流压气机旋转失速和喘振发生过程可以分为四个阶段,1)失速前阶段;2)旋转失速初始扰动阶段;3)完全失速阶段;4)喘振阶段.各阶段详细描述可参见文献[46].

关于旋转失速的建模,Moore和Greitzer从压气机转子和静子的整体性质出发,在流体拟定常和压气机半激盘等假设下,推导出一组描述压气机整体流场特性的偏微分方程(Partial differential equation,PDE),并在此基础上得到一个由三阶常微分方程(Ordinary differential equation,ODE)组描述的Moore-Greitzer模型.由于旋转失速是由偏微分方程描述的无限维分布参数系统产生的复杂动态现象,其无限维特性意味着任何基于有限状态测量的建模都是有限维近似建模.在Moore和Greitzer推导的偏微分方程基础上,Paduano和Mansoux等利用离散傅里叶变换及其逆变换推导出一个高阶ODE周向离散化Mansoux模型.

其中,φ=[φ1,φ2,···,φM]T是压气机周向上的M个测量点的流量,是压气机的平均压力增长,其他参数参见文献[19].

该Mansoux模型是一个描述旋转失速过程的有限维ODE系统,能够定量描述多种压气机旋转失速的发展过程,并可以在一定精度内产生与压气机试验台失速初始扰动相似的仿真结果,其系统状态可以看作是在压气机周向均匀布置2N+1(N为8的倍数)个流量传感器和在压气机入口和出口布置2个压力传感器获取的流量和压力信号.基于此,本节对高阶Mansoux模型进行研究和分析,对其系统动态进行建模和辨识,并提取出系统的动态时空复杂度特征.

通过确定学习理论,可对压气机失速前和旋转失速初始扰动的内在系统动态进行准确地建模和辨识,并将辨识得到的系统动力学轨迹以常值RBF神经网络的形式进行保存,可用于对压气机旋转失速前和失速初始扰动的状态和动力学轨迹进行动态特征提取和对比分析，详细过程见文献[19].轴流压气机失速前和失速初始扰动过程如图8所示.

当系统从失速前阶段(825转～875转)进入到旋转失速初始扰动阶段(876转～925转)时,对系统的状态轨迹和通过确定学习理论得到的动力学轨迹使用时空LZ复杂度指标进行特征表达,在计算具体复杂度时,采用步长为5转,在整个间隔内计算出21个时间点的复杂度,系统各个阶段状态轨迹和动力学轨迹的时间复杂度指标TLZC如图9所示,空间复杂度指标SLZC如图10所示.

图8 系统从失速前进入到旋转失速初始扰动阶段的过程Fig.8 Time evolution of the first flow state before rotating stall

从图9(a)和图9(b)可以看出,当系统从失速前进入到旋转失速初始扰动阶段时,系统动力学轨迹的时间复杂度指标TLZC变化的幅度为0.1921−0.1681=0.024,比系统状态轨迹的变化幅度0.2472−0.2284=0.019要大一些,说明从动力学轨迹更加能敏感反映出系统的状态变化.图10(a)和图10(b)的系统空间复杂度指标中也有相同的结果.

结合图9(a)、图9(b)和图10(a)、图10(b)系统状态轨迹和动力学轨迹的时空LZ复杂度,根据文中的敏感度系数式(16),系统状态轨迹和动力学轨迹的时间复杂度TLZC指标和空间复杂度SLZC指标相对系统参数变化的敏感度系数见表2,其中η1表示系统状态轨迹的敏感度系数,η2表示系统动力学轨迹的敏感度系数.

从表2中敏感度系数结果可以看出,通过确定学习理论对轴流压气机旋转失速过程进行准确地建模和辨识,提取出系统的动态特征复杂度指标,相比于状态轨迹复杂度指标,动力学轨迹可以更加敏感地反映系统从失速前到初始扰动过程的变化.此外,从时空复杂度特征指标图可以看出,在刚进入旋转失速初始扰动时(875转),系统动力学轨迹的时空复杂度指标有明显的改变,这与文献[19]的检测时间(当发生旋转失速时,可在885转检测出)相对应,从而验证了本文提出方法的可行性和优越性.

图9 失速前到旋转失速初始扰动阶段状态、动力学轨迹的时间复杂度TLZCFig.9 The TLZC index of the system state and dynamics trajectory before rotating stall

表2 失速前到初始扰动过程的时空复杂度指标敏感度系数Table 2 The sensitivity coefficients of the normal system to stall precursors

6 结论

本文提出了一个新的非线性系统动态特征提取方法.通过确定学习理论和时空LZ复杂度算法,实现了对非线性动力学系统从时域和空间域的角度进行动态特征提取.从Rossler系统状态轨迹和动力学轨迹的时间复杂度TLZC及空间复杂度SLZC指标对比结果以及频谱分布可以看出,与传统的时频分析方法和基于状态轨迹的复杂度特征提取方法相比,本文提出的新的动态特征提取方法可以更加敏感地反映出非线性动力学系统的周期、倍周期及混沌状态,其优点是从系统内在动力学的角度对原系统进行更好的表达.数值仿真和实验分析虽然验证了本方法的有效性,未来还需要对实际应用(例如航空发动机微小故障诊断等领域)进行进一步验证.

图10 失速前到旋转失速初始扰动阶段状态、动力学轨迹的空间复杂度SLZCFig.10 The SLZC index of the system state and dynamics trajectory before rotating stall

附录A

为了验证极值粗粒化过程的有效性以及数据序列长度对归一化复杂度的影响,选取一个高斯分布的随机序列作为例子进行说明[34,36].随机数据序列代表了大部分的混沌过程,其归一化复杂度接近1.在计算复杂度的过程中,采用50组不同的随机数据序列,然后对复杂度值求取算数平均.经过二值粗粒化BLZC、多值粗粒化MLZC(多值粗粒化中l选取6[36,47])和极值粗粒化ELZC过程后,归一化复杂度结果如图A1所示,实线代表经过BLZC的复杂度值,虚线代表经过MLZC的复杂度值,点线代表经过ELZC的复杂度值.

图A1 随机数据序列的归一化复杂度Fig.A1 The normalized Lempel-Ziv complexity of random sequence

从图A1可以看出,随机数据序列经过二值粗粒化BLZC过程后,其归一化的复杂度值都大于1,与文献[34−36]表述一致.同样,在数据序列的长度较小时,经过多值粗粒化MLZC过程的归一化复杂度值也大于1.而经过极值粗粒化ELZC过程得到的归一化复杂度值都小于1.此外,三个指标曲线中,当数据序列的长度为大于104时,归一化的复杂度值将趋于稳定.因此,本文选取的数据序列长度104.这个例子显示MLZC和ELZC都可以对数据序列进行较好处理,不同的是ELZC的粗粒化程度可以根据信号自身进行选择,而MLZC需要外部干预,凭借经验进行选择[40−42].