上海交通大学学报(自然版), 2021, 55(11): 1438-1444 doi: 10.16183/j.cnki.jsjtu.2020.256

基于MSOP的蒸汽动力系统单参数运行稳定性评估方法

郑奕扬1, 倪何,2, 金家善3

1.海军92118部队,浙江 舟山 316000

2.海军工程大学 动力工程学院,武汉 430033

3.海军工程大学 船舶与海洋学院,武汉 430033

An Operation Stability Assessment Method of a Single-Parameter in team Power System Based on MSOP

ZHENG Yiyang1, NI He,2, JIN Jiashan3

1. No.92118 Troop of PLA Navy, Zhoushan 316000, Zhejiang, China

2. College of Power ngineering, Naval University of Engineering, Wuhan 430033, China

3. College of Naval rchitecture and Marine Engineering, Naval University of Engineering, Wuhan 430033, China

通讯作者: 倪 何,男,副教授,电话(Tel.):13871124621;E-mail:elegance28@sina.com.

责任编辑: 石易文

收稿日期: 2020-08-11  

基金资助: 国家自然科学基金资助项目(51909254)

Received: 2020-08-11  

作者简介 About authors

郑奕扬(1996-),男,浙江省台州市人,硕士生,主要从事稳定性评估方法研究.

摘要

针对当前蒸汽动力系统缺乏有效稳定性评估方法的问题,提出一种适用于单个运行参数的稳定性评估方法.该方法是一种复合评估方法,首先综合运用中值回归经验模态分解和奇异值分解,对运行参数的时间序列进行分解并提取其隐藏的趋势项,然后依据各分量的最佳算法参数排列熵选取分量进行重构,最后通过非平稳时间序列分析中常用的整合滑动平均自回归模型预测趋势项和扰动项的走势并提取两者的分布特征,进而计算得到运行参数在预测趋势上各点的失稳概率并由此对其稳定性进行定量评估.经实际案例验证,证明该方法能够有效评估蒸汽动力系统单参数的运行稳定性,具有一定的理论创新性和工程应用价值.

关键词: 中值回归经验模态分解; 排列熵; 失稳概率; 稳定性评估

Abstract

Aimed at the lack of effective stability evaluation methods for the current steam power system, an operation stability assessment method suitable for single parameter is proposed. This method is a composite method, which first applied the midpoint and regression based empirical mode decomposition (MREMD) and singular value decomposition (SVD) to decompose the time series of operation parameters and extract their hidden trend terms. Then, the components are selected for reconstruction according to the optimal algorithm parameter permutation entropy (OAPPE) of each component. Finally, the auto-regressive integrated moving average (ARIMA ) model commonly used in the non-stationary time series analysis is utilized to predict the trend and the disturbance of parameters, and their distribution characteristics are also extracted in this process, based on which, the probability of instability (PI) of operation parameters at each point on the predicted trend are calculated, and their stabilities are quantitatively evaluated. The actual case proves that this method can effectively assess the operation stability of a single parameter of the steam power system, which has a certain theoretical innovation and engineering application value.

Keywords: midpoint and regression based empirical mode decomposition (MREMD); permutation entropy; probability of instability (PI); stability assessment

PDF (972KB) 元数据 多维度评价 相关文章 导出 EndNote| Ris| Bibtex  收藏本文

本文引用格式

郑奕扬, 倪何, 金家善. 基于MSOP的蒸汽动力系统单参数运行稳定性评估方法[J]. 上海交通大学学报(自然版), 2021, 55(11): 1438-1444 doi:10.16183/j.cnki.jsjtu.2020.256

ZHENG Yiyang, NI He, JIN Jiashan. An Operation Stability Assessment Method of a Single-Parameter in team Power System Based on MSOP[J]. Journal of shanghai Jiaotong University, 2021, 55(11): 1438-1444 doi:10.16183/j.cnki.jsjtu.2020.256

蒸汽动力系统因其组成设备多、系统耦合关系复杂的特点,在不同的运行工况下都保持稳定运行是相当困难的,因此蒸汽动力系统的运行稳定性评估相当重要.目前,稳定性评估在多系统耦合的复杂大系统中的应用停留在理论研究的层面上[1],不能全面评估系统的真实运行稳定性.而系统在一段运行时间内参数的时间序列则真实反应了系统的运行状态,分解参数的时间序列,对重构后得到的趋势项和扰动项进行预测分析可以从中提取反应系统运行稳定性的评估指标.

在非平稳时间序列趋势项提取方面,基于经验模态分解(EMD)及基本模式分量(IMF)重构是一种常用的算法[2,3,4].由于传统的经验模态分解不能完全提取所有IMF分量中包含的趋势信息,出现了大量的改进算法[5,6,7].文献[8]提出了集合经验模态分解(EEMD)、奇异值分解(SVD)和排列熵(PE)的趋势提取方法,并验证了该方法相对于传统分解算法的优越性,为本文研究提供了思路.但是,EEMD算法在原信号的基础上添加了不同的白噪声并通过重复集合平均进行抵消,分解效果取决于添加白噪声的去除,重构误差较大.同时,在排列熵的计算上较为粗糙,无法保证以此为依据进行选取分量的合理性.在时间序列分析中,非平稳时间序列的常用预测模型为整合滑动平均自回归模型(ARIMA)[9];由于该模型比较成熟且具有良好的通用性,所以本文在计算失稳概率(PI)时,采用ARIMA模型对运行参数时间序列的趋势项和扰动项进行预测.

综上所述,本文提出了一种名为MSOP的单参数运行稳定性评估方法,该方法是一种按照特定流程进行的复合评估方法,由中值回归经验模态分解(MREMD)、奇异值分解、基于最佳算法参数排列熵(OAPPE)重组和失稳概率(PI)计算等4个要素组成,其名称MSOP即为这4个要素的英文首字母组合.该方法经实际案例验证,可以根据运行参数的时间序列来定量评估蒸汽动力系统单参数的运行稳定性,可为蒸汽动力系统的整体运行稳定性评估提供底层评估指标输入,对蒸汽动力系统的在线运行稳定性管理和决策研究提供技术基础.

1 基于MSOP的运行稳定性评估方法

1.1 MREMD分解

MREMD[10]是经验模态分解的改进算法,通过优化均值包络的生成方式来抑制端点效应.对于输入长度为N的时间序列s(t),MREMD分解的具体步骤如下.

步骤1s(t)的左右端点处采用自回归(AR)模型进行延拓处理,找出时间序列所有极值点构成的序列X={x1,x2,…,xk},计算极值点序列的统计有义值.

P xi|xixσ;i=1,2,,k≥0.68

式中:k为极值点的个数;xσ为极值点序列的统计有义值;P{·}为概率.式(1)的数学意义为“全部极值点的绝对值小于统计有义值的概率不小于0.68”.

步骤2 取相邻极值点的均值构成时间序列的均值点序列 X-={ x-1, x-2, …, x-k+1},计算局部限制参数序列θ0:

θ0= x-1xσ,x-2xσ,,x-k+1xσ

采用三次样条曲线拟合得到均值点序列 X-的包络线 s-(t),计算得到时间序列s(t)的第1个分量s1(t)如下:

s1(t)=s(t)- s-(t)

步骤3 输出时序s1(t)为1阶IMF分量,并判断下列终止条件是否成立.

σ0-σ1/xσ0.25Pθ1|θ1θ00.95

式中:σ0和σ1分别为时间序列s(t)和s1(t)的标准差;θ1为时间序列s1(t)的局部限制参数序列.P θ1|θ1θ0的数学意义为“θ1中各参数不大于θ0中对应参数的统计概率不小于0.95”.

若终止条件式(4)成立,则直接进入步骤4,否则将s1(t)作为新的输入时序重复步骤1~2,计算得到2阶IMF分量;如此反复,直至满足终止条件式(4),记录此时IMF分量的最大阶数为n,并输出各阶IMF分量sg(t)(g=1,2,…,n).

步骤4 计算第n阶残余分量rn(t)=sn-1(t)-sn(t),若rn(t)为单调函数或常数,则结束分解过程;否则将rn(t)作为新的输入时序重复步骤1~3,最终得到的分解结果为

s(t)= g=0nsg(t)+rn(t)

将每一个IMF分量看作一个行向量,得到IMF分量矩阵如下:

SN×(n+1)= s1(t)Ts2(t)Tsn(t)Trn(t)T

1.2 SVD分解

由于趋势项具有低复杂度和低频等特点,舍去最为复杂的高频1阶IMF分量,得到矩阵SN×n=[s2(t)Ts3(t)Tsn(t)Trn(t)T],计算SN×n的协方差矩阵

Cn×n:Cn×n=ESN×n-E(SN×n)TSN×n-ESN×n

式中:E(·)为取均值计算.E(SN×n)为矩阵SN×n的均值矩阵,形式如下:

E(SN×n)= E(s2(t)T)E(s3(t)T)E(sn(t)T)E(rn(t)T)

对矩阵Cn×n进行SVD分解[11],重构后的奇异值分量矩阵QN×K

QN×K=SN×KUK×K= q1q2qK

式中:q1,q2,…,qK为分解得到的奇异值分量;U为奇异矩阵;K为非零奇异值个数.

1.3 基于OAPPE的趋势项重构

排列熵反映了时序的复杂程度,根据排列熵的大小可以对时间序列的所有极值点进行排序,从分解得到的奇异值分量矩阵QN×K中得到时间序列的预测趋势项,具体步骤如下.

步骤1 对所计算的极值点时间序列X={x1,x2,…,xk}进行相空间重构.

x1x1+τx1+(m-1)τx2x2+τx2+(m-1)τxkxk+τxk+(m-1)τ

式中:m为嵌入维数;τ为延迟时间.

步骤2 将每一行升序排列获得元素初始位置的索引值序列,计算m!种排列方式在索引值序列中出现的概率Pz,并由此计算排列熵Hpe:

Hpe=- z=1m!PzlnPz

步骤3 根据文献[12]提出的互信息算法确定延迟时间τ.首先,确定分段数d=1.87(n-τ-1)0.4,将一个二维坐标轴空间均匀划分成d2个网格,并将每个网格记为Rf(f=1,2,…,d2);然后,取时序X的前k-τ项{x1,x2,…,xk-τ}和后k-τ项{xτ+1,xτ+2,…,xk},构成A、B两个与延迟时间τ相关的时序空间,计算两个空间的互信息值I(τ),则有:

I(τ)=- f=1d2PAB(Rf)lgPAB(Rf)PA(Rf)PB(Rf)

式中:PA(Rf)和PB(Rf)为时序空间A和B中的点全部落在网格Rf中的概率;PAB(Rf)为时序空间A和B中的点同时落在网格Rf中的概率.

互信息值I(τ)实际上是联合分布PAB(Rf)与边缘分布PA(Rf)、PB(Rf)的相对熵,反映了两个空间系统的统计相关性,互信息值越小, 两个系统间的相关性就越小,本文取I(τ)的第1个极小值点对应的时间τ为最佳延迟时间.

步骤4 在获得最佳延迟时间的基础上,采用伪近邻法确定最小嵌入维数.对于时序X={x1,x2,…,xk},当嵌入维数为m时,重构时序X的集合Ω可以表示为

Ω={Xi|Xi={xi,xi+τ,…,xi+(m-1)τ};i=1,2,…,n-(m-1)τ}

由于集合Ω为有限集合,所以对于集合Ω中任意Xi,存在一个最近邻点 XMi满足:

Dm,i= Xi-XMi= minij(Xi-Xj) = miniju=0m-1[xi+-xj+]2

式中:Dm,i为XiXMi之间的距离.

假设嵌入维数为m和m+1时,Xi与最近邻点 XMi的距离分别为Dm,i和Dm+1,i,如果两者满足如下关系,判断 XMi为Xi的伪近邻点,则有:

Dm+1,i2-Dm,i212Dm,i>athDm+1,i2-Dm,i2121nj=1n(xj-1nj=1nxj)2>bth

式中:athbth为根据计算精度要求自定的阈值.

将集合Ω中所有元素的近邻点计算出来,组成集合Ψ可以表示为

Ψ={ XMi| XMi={ xMi, xMi+τ,…, xMi+(m-1)τ}; i=1,2,…,n-(m-1)τ}

根据式(14)筛选集合Ψ中所有伪近邻点,并计算伪近邻率(FNNP)FΨ,可以表示为

FΨ=NΨ/k

式中:NΨ近邻点集合Ψ中伪近邻点的个数.

当伪近邻率小于5%或者伪近邻率不再随着嵌入维数m的增大而减小时,则该点对应的嵌入维数即为最小嵌入维数[13].

步骤5 对于奇异值分量矩阵QN×K中所有的奇异值分量,按照步骤1~4计算其排列熵,然后采用K-means++算法[14]对计算得到的排列熵进行聚类分析,选择排列熵值小的奇异值分量进行叠加重构,获得时间序列的预测趋势项.

1.4 基于ARIMA模型的失稳概率计算

根据所获得的预测趋势项,预测运行参数时间序列在未来一小段时间内的趋势走向,计算这段时间内参数超出稳定范围的失稳概率,并以此作为其运行稳定性的评估指标.考虑到复杂热力系统的运行参数都具有趋势项非平稳的特点,所以选择在非平稳时间序列分析中较为常用的ARIMA模型进行失稳概率计算.ARIMA模型通过多次差分将非平稳时间序列变为平稳时间序列,在计算时有3个关键参数:自回归项数、滑动平均项数和差分次数,其计算流程和计算方法在文献[9]中已有详细描述,本节不再赘述.

对于预测趋势线上任一时间点的参数y,采用正态分布模型计算其失稳概率

Py:Py=P yymaxyymin

y~N(y-2)

式中:ymaxymin分别为稳定区域的上界和下界; y-为正态分布的均值,其为预测值;σ为标准差,由原时序与趋势项之间的偏差计算获得.

1.5 稳定性评估流程图

综上所述,基于MOSP的蒸汽动力系统单参数运行稳定性评估流程,如图1所示.

图1

图1   蒸汽动力系统单参数运行稳定性评估流程图

Fig.1   Flow chart of single parameter operation stability assessment of steam power system


对于输入的时间序列依次进行MREMD分解与SVD分解;通过最佳算法参数独立确定方法计算奇异值分量矩阵的排列熵;以排列熵为依据选取最低排列熵值的奇异值分量进行叠加重构,获得原时间序列的趋势项与扰动项;通过ARIMA模型预测未来一段时间内参数的变化趋势,并计算预测趋势上各点的失稳概率.

2 蒸汽动力系统单参数运行稳定性评估

为了验证所提方法的有效性,取某型蒸汽动力系统在一段时间内的微过热蒸汽压力为验证数据,在MATLAB 2016环境下编写代码完成评估全过程,各步骤的运行结果如下.

(1) 对输入的运行参数(微过热蒸汽压力)进行MREMD分解,获得的7个IMF分量(s1~s7)和残余分量(re)如图2所示,其中:t为时间.

图2

图2   由MREMD分解得到的IMF分量

Fig.2   IMF components obtained by MREMD decomposition


(2) 在去除1阶IMF分量后继续进行SVD分解,获得的7个奇异值分量(q1~q7)如图3所示.

图3

图3   由SVD分解得到奇异值分量

Fig.3   Singular value components obtained by SVD decomposition


(3) 为计算各奇异值分量的排列熵,通过互信息法确定最佳延迟时间,计算当τ0,100s时,各阶奇异值分量的互信息值(I1~I7)随τ的变化情况, 如图4所示.

图4

图4   互信息值随延迟时间的变化

Fig.4   Mutual information versus delay time


(4) 以最佳延迟时间为输入,计算当m0,100时,各阶奇异值分量的伪近邻率(F1~ F7)随m的变化情况,如图5所示.

图5

图5   伪近邻率随嵌入维数的变化

Fig.5   False nearest neighbor percents versus embedding dimension


(5) 根据1.3节给出的依据选取最小嵌入维数,并计算各奇异值分量的排列熵,计算结果如表1所示.

表1   各奇异值分量的最优排列熵算法参数及排列熵

Tab.1  Optimal permutation entropy algorithm parameters and permutation entropy of each singular value component

奇异值分量最佳延迟时间/s最小嵌入维数排列熵
q1644.5252
q2455.9989
q3844.5294
q42331.6453
q51631.9715
q61732.6032
q71932.7345

新窗口打开| 下载CSV


(6) 根据排列熵对各奇异值分量进行K-means聚类,并选择排列熵值小的分量进行重组.根据分类结果,选取q4q5q6q7进行叠加重构,得到的微过热蒸汽压力的预测趋势项与扰动项如图6所示.其中:pst为微过热蒸汽的压力值;psr为微过热蒸汽压力的扰动值.

图6

图6   微过热蒸汽压力数据的预测趋势项与扰动项

Fig.6   Predicted trends and disturbance items of slightly superheated steam pressure data


(7) 选择自回归项数为3、滑动平均项数为2、差分次数为1,建立ARIMA(3,1,2)模型对微过热蒸汽压力进行预测.为检验模型的预测效果,取图6(a)中原信号数据的前80%作为训练集,后20%作为验证集,对模型预测能力进行检验,对比结果如图7所示.

图7

图7   ARIMA模型预测能力检验

Fig.7   Predictive ability test of ARIMA model


图7可见,模型输出的预测趋势项均落在验证数据的95%置信区间内,说明模型参数选择合理、预测精度较高.

(8) 该型蒸汽动力系统的微过热蒸汽压力整定值为2.45 MPa,取其±3%的波动区间为稳定区间,应用ARIMA(3,1,2)模型进行预测,计算得到的微过热蒸汽压力在未来60 s的变化趋势及失稳概率如图8所示,其中:P'为失稳概率.

图8

图8   预测趋势及失稳概率

Fig.8   Predicted trends and PI of each point


图8可知,微过热蒸汽压力在未来60 s内的失稳概率较小,均在0.2以下,说明其运行稳定性较好,失稳风险较低.

3 结语

针对目前蒸汽动力系统缺少定量运行稳定性评估方法的问题,提出一种基于MSOP的单参数稳定性评估方法,并通过对某型蒸汽动力系统微过热蒸汽压力的评估分析,验证了该方法的科学性.本文的研究成果具有一定的理论创新性和工程应用价值,通过增加数据接口读取实际系统的实时运行数据,实现蒸汽动力系统的在线运行稳定性评估,为操作人员的使用管理提供参考.同时,本研究也可为蒸汽动力系统的全系统运行稳定性提供底层评估指标,通过建立蒸汽动力系统运行稳定性评估指标体系,实现整个系统的运行稳定性评估.

参考文献

唐韧博, 杨建东.

水泵水轮机不同工况点的稳定性分析

[J]. 水力发电学报, 2016, 35(5):117-122.

[本文引用: 1]

TANG Renbo, YANG Jiandong.

Operation stability analysis of pump-turbine in different working conditions

[J]. Journal of Hydroelectric Engineering, 2016, 35(5):117-122.

[本文引用: 1]

LU S Z, WANG X H, YU H L, et al.

Trend extraction and identification method of cement burning zone flame temperature based on EMD and least square

[J]. Measurement, 2017, 111:208-215.

DOI:10.1016/j.measurement.2017.07.047      URL     [本文引用: 1]

CHEN Y Q, ZHENG W, WEI X B.

Application of EMD to integrated signal trend extraction

[J]. Advanced Materials Research, 2012, 591/592/593:2072-2076.

[本文引用: 1]

DONG Y, SUN Y.

Trend extraction by difference method for least square combining with EMD

[J]. International Journal of Advancements in Computing Technology, 2012, 4(15):512-519.

DOI:10.4156/ijact      URL     [本文引用: 1]

BHUIYAN S M A, CAREY S, KHAN J F.

Analysis of signal trend by a pseudo-EMD method with applications to weather and speech data

[J]. Advances in Adaptive Data Analysis, 2013, 5(2):1350009.

DOI:10.1142/S179353691350009X      URL     [本文引用: 1]

BREAKER L C.

Energy production trend extraction using ensemble empirical mode decomposition

[J]. International Journal of Energy and Statistics, 2013, 1(3):195-204.

DOI:10.1142/S2335680413500130      URL     [本文引用: 1]

李精明, 魏海军, 魏立队, .

摩擦振动信号的EEMD和多重分形去趋势波动分析

[J]. 哈尔滨工程大学学报, 2016, 37(9):1204-1208.

[本文引用: 1]

LI Jingming, WEI Haijun, WEI Lidui, et al.

Ensemble empirical mode decomposition and multifractal detrended fluctuation analysis of frictional vibration signals

[J]. Journal of Harbin Engineering University, 2016, 37(9):1204-1208.

[本文引用: 1]

陈亮, 刘宏立, 郑倩, .

基于EEMD-SVD-PE的轨道波磨趋势项提取

[J]. 哈尔滨工业大学学报, 2019, 51(5):171-177.

[本文引用: 1]

CHEN Liang, LIU Hongli, ZHENG Qian, et al.

An EEMD-SVD-PE approach to extract the trend of track irregularity

[J]. Journal of Harbin Institute of Technology, 2019, 51(5):171-177.

[本文引用: 1]

KIRBAŞ İ, SÖZEN A, TUNCER A D, et al.

Comparative analysis and forecasting of COVID-19 cases in various European countries with ARIMA, NARNN and LSTM approaches

[J]. Chaos, Solitons & Fractals, 2020, 138:110015.

DOI:10.1016/j.chaos.2020.110015      URL     [本文引用: 2]

黄礼敏.

海浪中非平稳非线性舰船运动在线预报研究

[D]. 哈尔滨: 哈尔滨工程大学, 2016.

[本文引用: 1]

HUANG Limin.

On-line prediction of non-stationary and nonlinear ship motions at sea

[D]. Harbin: Harbin Engineering University, 2016.

[本文引用: 1]

UDDIN J, KIM J M, ISLAM R.

Texture analysis based feature extraction using Gabor filter and SVD for reliable fault diagnosis of an induction motor

[J]. International Journal of Information Technology and Management, 2018, 17(1/2):20.

DOI:10.1504/IJITM.2018.089452      URL     [本文引用: 1]

FRASER A M, SWINNEY H L.

Independent coordinates for strange attractors from mutual information

[J]. Physical Review A, General Physics, 1986, 33(2):1134-1140.

DOI:10.1103/PhysRevA.33.1134      URL     [本文引用: 1]

KENNEL, BROWN, ABARBANEL.

Determining embedding dimension for phase-space reconstruction using a geometrical construction

[J]. Physical Review A, Atomic, Molecular, and Optical Physics, 1992, 45(6):3403-3411.

DOI:10.1103/PhysRevA.45.3403      URL     [本文引用: 1]

VASSILVITSKII S, ARTHUR D. K-means++: The advantages of careful seeding[EB/OL]. [2020-04-11]. https://maheshgoyani.in/resource/workshop/kmeans++.pdf.

URL     [本文引用: 1]

/