2.845

2023影响因子

(CJCR)

  • 中文核心
  • EI
  • 中国科技核心
  • Scopus
  • CSCD
  • 英国科学文摘

留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

基于群决策的道岔控制电路故障诊断方法

董炜 刘明明 王良顺 赵辉 辜勋

杨杰, 柴天佑, 张亚军, 吴志伟. 数据与模型驱动的电熔镁群炉需量预报方法. 自动化学报, 2018, 44(8): 1460-1474. doi: 10.16383/j.aas.2017.c160597
引用本文: 董炜, 刘明明, 王良顺, 赵辉, 辜勋. 基于群决策的道岔控制电路故障诊断方法. 自动化学报, 2018, 44(6): 1005-1014. doi: 10.16383/j.aas.2017.c160715
YANG Jie, CHAI Tian-You, ZHANG Ya-Jun, WU Zhi-Wei. Data and Model Driven Demand Forecasting Method for Fused Magnesium Furnace Group. ACTA AUTOMATICA SINICA, 2018, 44(8): 1460-1474. doi: 10.16383/j.aas.2017.c160597
Citation: DONG Wei, LIU Ming-Ming, WANG Liang-Shun, ZHAO Hui, GU Xun. Fault Diagnosis for Railway Turnout Control Circuit Based on Group Decision Making. ACTA AUTOMATICA SINICA, 2018, 44(6): 1005-1014. doi: 10.16383/j.aas.2017.c160715

基于群决策的道岔控制电路故障诊断方法

doi: 10.16383/j.aas.2017.c160715
基金项目: 

国家重点研发计划 2017YFB1200700

国家自然科学基金 61490701

苏州-清华创新引领行动专项 2016SZ0202

详细信息
    作者简介:

    董炜 北京信息科学与技术国家研究中心副研究员.主要研究方向为复杂工程系统的建模与仿真, 故障诊断与预测维护, 自动测试与安全评估.E-mail:weidong@mail.tsinghua.edu.cn

    王良顺 北京信息科学与技术国家研究中心博士后.主要研究方向为故障诊断, 智能控制.E-mail:wangliangshun340@163.com

    赵辉北京信息科学与技术国家研究中心博士后.主要研究方向为人工智能, 数据挖掘.E-mail:taiyuanjifeng@tsinghua.edu.cn

    辜勋 中国人民解放军63956部队高级工程师.主要研究方向为仿真建模, 机械装备, 知识产权.E-mail:longouzi@163.com

    通讯作者:

    刘明明  清华大学自动化系硕士研究生.主要研究方向为铁路信号故障诊断系统.本文通信作者.E-mail:13466608442@163.com

Fault Diagnosis for Railway Turnout Control Circuit Based on Group Decision Making

Funds: 

National Key Research and Development Program of China 2017YFB1200700

National Natural Science Foundation of China 61490701

Special Fund of SuzhouTsinghua Innovation Leading Action 2016SZ0202

More Information
    Author Bio:

    Associate professor at Beijing National Research Center for Information Science and Technology. His research interest covers modeling and simulation of complex engineering systems, fault diagnosis and prediction maintenance, automatic testing and safety assessment

    Postdoctoral at Beijing National Research Center for Information Science and Technology. His research interest covers fault diagnosis and intelligent control

    Postdoctoral at Beijing National Research Center for Information Science and Technology. His research interest covers artificial intelligence and data mining

    Senior engineer at PLA 63956 Troops. His research interest covers simulation modeling, mechanized equipment, and intellectual property right

    Corresponding author: LIU Ming-Ming Master student in the Department of Automation, Tsinghua University. His main research interest is fault diagnosis of railway signal. Corresponding author of this paper
  • 摘要: 高速铁路道岔是与高速列车直接接触的重要信号设备,其控制电路的故障检测手段目前仍停留在简单仪器与人的经验相结合的方式.为了实现道岔控制电路故障的智能诊断,提高故障诊断的准确率并降低单一诊断方法带来的不确定性,本文提出一种基于群决策的诊断方法:首先根据道岔控制电路的特点,总结了典型的11个故障模式和对应的8个故障特征;其次,分别采用模糊理论、神经网络和支持向量机(Support vector machine,SVM)对道岔控制电路进行故障诊断;然后引入群决策理论将三种方法视为决策专家,通过群基数效应集结方式实现决策级上的信息融合从而得到群专家综合评判的诊断结果.从仿真数据的验证来看,该方法比单一方法的故障诊断的准确率要高,表明了本文所提方法能够实现三种方法的互补融合,也提高了故障诊断的准确率,在该领域有着良好的应用前景.
  • 电熔镁炉是一种以菱镁矿为原料, 由电流控制器控制熔炼电流来生产电熔镁砂的重要设备.产品电熔镁砂是一种应用于冶金、化工、航天等领域的重要高级耐火材料.电熔镁群炉需量指当前时刻和当前时刻之前一定时间内群炉功率的平均值, 用于度量高耗能电熔镁群炉的用电量.在生产过程中需量不得超过规定的最大需量即需量峰值, 以限制电熔镁群炉的用电量.需量监控系统对群炉需量进行实时监控, 当需量超过需量峰值的限幅值时会切断某台炉供电, 以保证群炉需量不超过需量峰值; 当需量低于限幅值时再恢复该台炉供电, 使该炉继续生产.为了保证电熔镁砂质量, 需要电流控制器将电流控制在工艺设定值附近.在电流控制器的调节作用下, 熔炼过程中原料杂质成分含量增大和颗粒长度变大可能会导致需量尖峰, 即需量先增大超过限幅值而后下降低于限幅值.而需量尖峰会造成切断电熔镁炉的供电.然而, 切断供电会破坏炉内温度场吸热和放热平衡, 降低电熔镁砂质量, 因此对需量进行准确的预报对于避免尖峰时刻的错误拉闸显得十分重要.

    近年来, 针对电力系统的功率预报问题相关学者开展了一系列研究, 多采用时间序列方法[1]、支持向量回归机[2]、神经网络[3]和混合方法[4]等.例如, 文献[1]为了制定某地区的发电量的日计划, 基于该地区过去14天用电总功率的数据(采样周期为15分钟), 采用相似形时间序列法先得到参考基准曲线, 然后预报未来一天用电总功率的曲线.文献[2]为了提高某地电力系统的可靠性, 使用该地过去7天用电总功率的数据(采样周期为1小时), 采用支持向量回归机的方法预报未来一天用电总功率曲线.文献[3]为了合理调度某地电力系统, 用该地过去7天风力发电总功率的数据(采样周期为0.5小时), 采用前馈神经网络和上下限估计的方法预报未来一天风力发电总功率的上下限.文献[4]为了降低某地电力市场的经济损失风险, 基于昨天和一周前的一天的用电功率数据和温度数据(采样周期为1小时), 采用遗传算法-径向基函数神经网络(Genetic algorithm --- radial basis function neural network, GA-RBFNN)方法预报未来一天的用电功率.文献[5]为了准确预报群炉需量变化趋势减少不必要的切断供电, 使用某电熔镁砂厂过去时刻的群炉功率数据(采样周期为7秒), 采用RBFNN方法预报下一时刻的群炉需量.

    在电熔镁炉生产过程中, 控制系统通过调节电极位置, 改变电弧弧长, 进而控制熔炼电流稳定在电流设定值附近, 从而实现在满足产量约束的条件下尽量降低产品单吨能耗的控制目标[6].电熔镁炉的熔炼电流通常在15 000 A左右, 电能消耗巨大.熔炼过程中, 当原料杂质成分含量增大和颗粒长度变大时, 工作电阻减小, 电流变大, 需量上升.此时电流控制器会调节弧长使工作电阻变大, 降低电流, 需量又随之下降, 这样就会出现需量先升高后下降的尖峰现象.因此根据功率变化特性建立需量的动态模型才能更准确地预报需量.而文献[1-4]预报的对象采样周期时间尺度较大, 与电熔镁群炉的功率变化特性不同, 而且只单纯依据过去功率的数据进行预报, 没有研究对象的动态特性分析.因此上述文献难以直接适用预报电熔镁群炉需量.文献[5]虽然对电熔镁群炉需量进行了预报, 但只是将历史功率数据和功率变化率整体视为非线性函数关系进行处理, 预报精度有待提高.

    本文首先建立和分析需量动态模型, 提出由功率变化率的线性项、基于PACF输入变量决策的RBFNN未知非线性项估计和需量计算模型组成的需量预报方法, 通过某电熔镁砂厂实际数据的仿真实验和工业实验表明所提方法能够准确预报需量的变化趋势.

    为了准确预报$(k+1)$时刻需量, 本文首先根据需量定义建立$(k+1)$时刻的需量模型, 然后建立功率变化率与电流控制系统输出电流之间的由线性项与未知非线性项组成的动态模型.

    图 1所示, 在熔炼电压$U$作用下, 电流控制器调节升降电机转速$u_i$, 使电流实际值$y_i$跟踪电流设定值$y^{*}$而产生用电功率.炉内原料吸收电弧释放的热量熔化形成逐渐上涨的MgO熔池, 熔炼结束后熔池经过冷却、结晶、破碎等工序形成产品电熔镁砂.熔炼过程中功率变送器测量电力变压器得到群炉功率数据$p(k)$, $p(k - 1)$, $\cdots$, $p(k - n + 1)$, 由需量计算装置基于定义1可得到当前时刻的群炉需量$\bar P(k)$.

    图 1  电熔镁群炉需量监控原理图
    Fig. 1  Schematic diagram of demand monitoring process for FMFG

    定义1. 群炉需量$\bar P(k)$为$k$时刻和$(k-1)$, $\cdots$, $(k-n+1)$时刻群炉功率$p(k)$的平均值

    $ \begin{align} \label{eq1} \bar P(k) = \frac{1}{n}\sum\limits_{j = 0}^{n - 1} {p(k - j)} \end{align} $

    (1)

    其中, 群炉功率$p(k)$为

    $ \begin{align} \label{eq2} p(k) = \sum\limits_{i = 1}^m {\sqrt 3 U{y_i}(k)\cos \varphi } \end{align} $

    (2)

    其中, $m$为电熔镁群炉台数, 各炉熔炼电压为$U$ (常量), ${y_i}(k)$为第$i$台电熔镁炉的熔炼电流, $\rm cos\varphi$为功率因数.

    由定义1可递推得到$(k+1)$时刻的需量$\bar{P}(k+1)$.

    $ \begin{align} \label{eq3} \bar P(k + 1) = &\ \frac{1}{n}\sum\limits_{j = 0}^{n - 1} {p(k + 1 - j)} =\notag\\ &\ \bar P(k) + \frac{{ p(k) - p(k + 1 - n) + \Delta p(k)}}{n} \end{align} $

    (3)

    其中, $\bar P(k)$, $p(k)$, $p(k + 1 - n)$都已知, 可见预报$\bar P(k + 1)$的关键在于估计未知的功率变化率$\Delta p(k)$.

    式(3)中功率变化率$\Delta p(k)$的定义为

    $ \begin{align} \label{eq4} \Delta p(k) =&\ p(k + 1) - p(k)=\nonumber\\ &\ \sum\limits_{i = 1}^m {\sqrt 3 U\left( {{y_i}(k + 1) - {y_i}(k)} \right)\cos \varphi } \end{align} $

    (4)

    其中, $y_i(k+1)$和$y_i(k)$是电流闭环控制系统的输出.

    首先建立以电极升降电机转速$u_i$为输入, 以电流$y_i$为输出的被控对象模型.电流$y_i$与工作电阻$R_i$之间关系有

    $ \begin{align} {y_i} = \frac{U}{\sqrt{3} {R_i}} \end{align} $

    (5)

    其中, $R_i$[7]为电弧电阻$R_{i, {\rm arc}}$和熔池电阻$R_{i, {\rm pool}}$之和, 即

    $ \begin{align} \label{eq6} R_i = {R_{i, {\rm arc}}} + {R_{i, {\rm pool}}} \end{align} $

    (6)

    其中, 电弧电阻$R_{i, {\rm arc}}$[8]

    $ \begin{align} \label{eq7} {R_{i, {{\rm arc}}}} =&\ {G_0}{L_{i, {{\rm arc}}}} = {G_0}\left( {{h_{i, {{\rm elec}}}} - {h_{i, {{\rm pool}}}}} \right)= \nonumber\\ &\ {G_0}\left( {\int_0^{T} {{u_i}(t)r{\rm d}t} - {h_{i, {{\rm pool}}}}({B_{i, 1}}, {B_{i, 2}}, {y_i})} \right) \end{align} $

    (7)

    其中, 参数${G_0} = {( {{g_0}\pi r_{{{\rm arc}}}^2\exp ( {{{ - {T_0}} / {{T_1}}}} )} )^{-1}}$, $g_{0}$为电弧电导率常数, ${r_{{{\rm arc}}}}$为电弧弧柱半径, $T_{0}$为气体电离温度常数, $T_{1}$为电弧间隙温度, $L_{i, {\rm arc}}$为电弧长度, $h_{i, {\rm elec}}$为电极末端位置, $T$为运行时间, $u_i$为升降电机转速, $r$为升降机构等效齿轮半径[9], $h_{i, {\rm pool}}(B_{i, 1}$, $B_{i, 2}, y_i)$为熔池高度, $B_{i, 1}$为原料杂质成分含量, $B_{i, 2}$为原料颗粒长度.

    熔池电阻$R_{i, {\rm pool}}$[10]

    $ \begin{align} \label{eq8} R_{i, {{\rm pool}}} =&\ \frac{\rho _{i, {\rm pool}}(B_{i, 1}, B_{i, 2})}{{\rm \pi }D} \times \nonumber\\ &\ \left( {1 - \frac{D}{{2{h_{i, {{\rm pool}}}}({B_{i, 1}}, {B_{i, 2}}, {y_i})}}} \right) \end{align} $

    (8)

    其中, $\rho_{i, {{\rm pool}}}$为熔池电阻率, $D$为熔池直径.由式(7)和式(8), 得

    $ \begin{align} \label{eq9} \frac{{{\rm d}{R_{i, {\rm arc}}}}}{{{\rm d}t}} = &\ {G_0}\left( {u_ir - \frac{{{\rm d}{h_{i, {\rm pool}}}({B_{i, 1}}, {B_{i, 2}}, y_i)}}{{{\rm d}t}}} \right) \nonumber\\[2mm] \frac{{{\rm d}{R_{i, {\rm pool}}}}}{{{\rm d}t}} = &\ \frac{{{\rho _{i, {\rm pool}}({B_{i, 1}}, {B_{i, 2}})}}}{{2\pi h_{i, {\rm pool}}^2}} \times\nonumber\\ &\ \frac{{{\rm d}{h_{i, {\rm pool}}}({B_{i, 1}}, {B_{i, 2}}, y_i)}}{{{\rm d}t}} \end{align} $

    (9)

    对式(5)两边求导, 代入式(9)得到以升降电机转速$u_i$为输入, 电流$y_i$为输出的动态模型, 即

    $ \begin{align} \label{eq10} \frac{{{\rm d}y_i}}{{{\rm d}t}} =& - \frac{{\sqrt 3 y_i^2}}{U}\times\frac{{{\rm d}{R_i}}} {{{\rm d}t}} =\notag \\ & -\frac{{\sqrt 3 {y_i^2}}}{U}\left( {\frac{{{\rm d}{R_{i, {\rm arc}}}}}{{{\rm d}t}} + \frac{{{\rm d}{R_{i, {\rm pool}}}}}{{{\rm d}t}}} \right)=\notag\\ & -\Bigg[ \left( {\frac{{{\rho_{i, {\rm pool}}(B_{i, 1}, B_{i, 2})}}} {{2\pi h_{i, {\rm pool}}^2}} - {G_0}} \right)\times\notag \\ &\ \frac{{{\rm d}{h_{i, {\rm pool}}(B_1, B_2, y_i)}}}{{{\rm d}t}} \vphantom{\left( {\frac{{{\rho_{pool}(B_1, B_2)}}} {{2\pi h_{i, {\rm pool}}^2}} - {G_0}} \right)\frac{{\sqrt 3 }}{U} } \Bigg]\frac{{\sqrt 3 }}{U} {y_i^2} - \frac{{\sqrt 3 {G_0}r}} {U}u_i{y_i^2} \end{align} $

    (10)

    显然, 式(10)中${u_i}y_i^2$表明电流$y_i$与升降电机转速$u_i$之间为非线性关系, 且$y_i^2$系数中包含未知的熔池电阻率${\rho _{i, {{\rm pool}}}}({B_{i, 1}}, {B_{i, 2}})$和熔池高度非线性变化$\frac{{\rm d}{h_{i, {\rm pool}}}(B_{i, 1}, B_{i, 2}, y_i)}{{\rm d}t}$.由于工艺要求将电流控制在电流设定值附近, 电熔镁炉在工作点附近运行.将其在工作点附近线性化后, 被控对象(10)可由一阶线性模型与未建模动态的形式表示, 即

    $ \begin{align} \label{eq11} & A({z^{ - 1}})y_i(k + 1) =\notag\\ &\qquad B({z^{ - 1}})u_i(k) + v_i[y_i(k), u_i(k), {\mathit{\boldsymbol{ h}}_{i, {\rm pool}}}] \end{align} $

    (11)

    其中, $A(z^{-1})$, $B(z^{-1})$是关于$z^{-1}$的多项式, $A({z^{ - 1}})$ $=$ $1+{a_1}{z^{ - 1}}$, $B({z^{ - 1}}) = {b_{{\rm 0}}}$, ${\mathit{\boldsymbol{ h}}_{i, {\rm pool}}} \in {{\bf R} ^{{n_c}}}$, $n_c$未知. $v_i\left[ \cdot \right]$为高阶非线性项, 称未建模动态.

    设计电流PID控制器为

    $ \begin{align} \label{eq12} H({z^{ - 1}}){u_i}(k) = G({z^{ - 1}})[{y^*} - {y_i}(k)] \end{align} $

    (12)

    其中, $H({z^{-1}}) = 1 - {z^{-1}}$, $G({z^{-1}}) = {g_0} + {g_1}{z^{-1}} + {g_2}{z^{-2}}$.将控制器(12)作用于被控对象(11), 得电流控制系统的闭环方程, 即

    $ \begin{align} \label{eq13} &T({z^{ - 1}}) y_i(k + 1) =B({z^{-1}})G({z^{-1}}){y^*} +\nonumber\\ &\qquad (1- {z^{-1}})v_i[y_i(k), u_i(k), {\mathit{\boldsymbol{ h}}_{i, \rm {pool}}}] \end{align} $

    (13)

    其中, $T({z^{-1}}) = (1-{z^{- 1}})A({z^{-1}}) + {z^{-1}}B({z^{-1}})$ $\times$ $G({z^{-1}})$, 由式(13)得到电流控制系统闭环方程的输出, 即

    $ \begin{align} \label{eq14} {y_i}(k + 1) =& - {T^*}({z^{ - 1}}){y_i}(k + 1) +\nonumber\\ &\ B({z^{ - 1}})G({z^{ - 1}}){y^*} +\nonumber\\ &\ (1 - {z^{ - 1}}){v_i}[{y_i}(k), {u_i}(k), {\mathit{\boldsymbol{ h}}_{i, \rm {pool}}}] \end{align} $

    (14)

    其中, ${T^*}({z^{ - 1}}) = T({z^{ - 1}}) - 1$.从式(14)可以看出, $y_i(k+1)$与过去时刻的电流$ - {T^*}({z^{ - 1}})y_i(k + 1)$、电流设定值$y^*$、未建模动态${v_i}[{y_i}(k), {u_i}(k), {\mathit{\boldsymbol{ h}}_{i, \rm {pool}}}]$相关.

    将式(14)代入式(4)中可得功率变化率$\Delta p(k)$的动态模型.

    $ \begin{align} \label{eq15} \Delta p(k) = &\ \sum\limits_{i = 1}^m {\sqrt 3 U\left( {{y_i}(k + 1) - {y_i}(k)} \right)\cos \varphi } =\nonumber\\ &\ \left[ - z{T^*}({z^{ - 1}}) - 1\right]C\sum\limits_{i = 1}^m {{y_i}(k)} + \nonumber\\ &\ B({z^{ - 1}})G({z^{ - 1}})mC{y^*} + \nonumber\\ &\ C\sum\limits_{i = 1}^m {(1 - {z^{ - 1}}){v_i}[{y_i}(k), {u_i}(k), {\mathit{\boldsymbol{ h}}_{i, \rm {pool}}}]} \end{align} $

    (15)

    其中, $C = \sqrt 3 U\cos \varphi$, $C\sum_{i = 1}^m (1 - {z^{ - 1}}){v_i}[{y_i}(k), $ ${u_i}(k), $ $\mathit{\boldsymbol{ h}}_{i, \rm {pool}}]$为高阶非线性项.式(15)表明了功率变化率$\Delta p(k)$与电流${y_i}(k)$之间的关系.

    由式(5)可知影响电流${y_i}(k)$的直接因素是工作电阻${R_i}(k)$.将式(5)代入式(15)可得功率变化率$\Delta p(k)$与工作电阻${R_i}(k)$之间的关系.

    $ \begin{align} \label{eq16} & \Delta p(k) = \left[ - z{T^*}({z^{ - 1}}) - 1\right]C\sum\limits_{i = 1}^m {\left[ {\frac{U}{{\sqrt 3 {R_i}(k)}}} \right]} + \nonumber\\ &\qquad B({z^{ - 1}})G({z^{ - 1}})mC{y^*} + \nonumber\\ &\qquad C\sum\limits_{i = 1}^m {(1 - {z^{ - 1}}){v_i}\left[ {\frac{U}{{\sqrt 3 {R_i}(k)}}, {u_i}(k), {\mathit{\boldsymbol{ h}}_{i, \rm {pool}}}} \right]} \end{align} $

    (16)

    再由式(7)和式(8)进一步得到功率变化率$\Delta p(k)$和电流影响因素(原料杂质成分含量$ B_{i, 1} $、颗粒长度$ B_{i, 2} $)之间的关系

    $ \begin{align} &\Delta p(k) = \left[ - z{T^*}({z^{ - 1}}) - 1\right] \times\nonumber\\ &\qquad C\sum\limits_{i = 1}^m {\left[ {\frac{U}{{\sqrt 3 {f_{{R_i}}} \left( {{B_{i, 1}}(k), {B_{i, 2}}(k)} \right)}}} \right]} + \nonumber\\ \label{eq17} &\qquad B({z^{ - 1}})G({z^{ - 1}})mC{y^*} + C\sum\limits_{i = 1}^m (1 - {z^{ - 1}}) \times\nonumber\\ &\qquad {v_i}\left[ {\frac{U}{{\sqrt 3 {f_{{R_i}}}\left( {{B_{i, 1}}(k), {B_{i, 2}}(k)} \right)}}, {u_i}(k), \mathit{\boldsymbol{ h}}_{i, \rm {pool}}} \right] \end{align} $

    (17)

    其中, $ {f_{{R_i}}}\left( {{B_{i, 1}}(k), {B_{i, 2}}(k)} \right) $为关于原料杂质成分含量$ B_{i, 1} $、颗粒长度$ B_{i, 2} $的函数.

    $ \begin{align} \label{eq18} &{f_{{R_i}}}\left( {{B_{i, 1}}(k), {B_{i, 2}}(k)} \right)=\nonumber\\ &\qquad {G_0}\left( {{h_{i, {{\rm elec}}}}(k) - {h_{i, {{\rm pool}}}}[{B_{i, 1}} (k), {B_{i, 2}}(k)]} \right)+ \nonumber\\ &\qquad \frac{{{\rho _{i, {{\rm pool}}}}[{B_{i, 1}}(k), {B_{i, 2}}(k)]}}{{ {{\rm \pi }}D}} \times\nonumber\\ &\qquad \left( {1 - \frac{D}{{2{h_{i, {{\rm pool}}}}[{B_{i, 1}}(k), {B_{i, 2}}(k)]}}} \right) \end{align} $

    (18)

    由式(3), (15)~(18)可以看出群炉需量与电流及电流因素之间的关系.当多台电熔镁炉出现原料杂质成分含量$ B_{i, 1}$增大和颗粒长度$ B_{i, 2} $变大的工况时, 会使熔池高度$ h_{i, {\rm pool}} $升高[11], 电弧长度${L_{i, {{\rm arc}}}}$减小, 进而电阻$ R_i $减小, 电流$ y_i $变大, 功率变化率$ \Delta p$ $>$ $0 $, 需量$ \bar{P} $增大.由于电流控制系统调节电机转速$ u_i $使电弧长度$ L_{i, {\rm arc}} $增大, 又使电流$ y_i $减小到设定值$ y^* $附近, 需量$ \bar{P} $又减小, 导致需量$ \bar{P} $先增大然后减小, 出现需量尖峰.

    $ u_i(k) $与过去时刻的电流$y_i(k)$相关, 可表示为

    $ \begin{align}\label{eq19} {u_i} = {u_i}\left( {{y_i}(k)} \right) \end{align} $

    (19)

    熔炼过程中熔池高度$ {\mathit{\boldsymbol{ h}}_{i, {\rm pool}}} $增大, 电流$y_i(k)$会随之增大, 因此$ {\mathit{\boldsymbol{ h}}_{i, {\rm pool}}} $与$y_i(k)$也相关, 即

    $ \begin{align} \label{eq20} {\mathit{\boldsymbol{ h}}_{i, {\rm pool}}} = {\mathit{\boldsymbol{ h}}_{i, {\rm pool}}} \left(y_i(k)\right) \end{align} $

    (20)

    所以$ {v_i}[{y_i}(k), {u_i}(k), {\mathit{\boldsymbol{ h}}_{i, {\rm pool}}}] $可简记为电流$ \mathit{\boldsymbol{ y}}_i $的非线性函数$ {V_i}({\mathit{\boldsymbol{ y}}_i}) $, $ {\mathit{\boldsymbol{ y}}_i} = {[{y_i}(k), {y_i}(k - 1), \cdots]^{{\rm T}}} $.将式(15)改写为

    $ \begin{align} \label{eq21} \Delta p(k) =& \left[ - z{T^*}({z^{ - 1}}) - 1\right]p(k) + \nonumber\\ &\ B({z^{ - 1}})G({z^{ - 1}}){p^*} + \nonumber\\ &\ C\sum\limits_{i = 1}^m {(1 - {z^{ - 1}}){V_i}\left( {{\mathit{\boldsymbol{ y}}_i}} \right)} \end{align} $

    (21)

    其中, $ C\sum_{i = 1}^m {(1 - {z^{ - 1}}){V_i}\left( {{\mathit{\boldsymbol{ y}}_i}} \right)} $为各炉电流的函数乘以熔炼电压$ U $和功率因数$ \cos\varphi $再累加, 与群炉功率有关.故可将线性组合$[ - z{T^*}({z^{ - 1}}) - 1]p(k) + B({z^{ - 1}})G({z^{ - 1}}){p^*} $之外的部分记为非线性函数$ \bar V(\mathit{\boldsymbol{ p}}) $, 则$ \Delta p(k) $的动态模型进一步改写为

    $ \begin{align} \label{eq22} \Delta p(k) = &\ [ - z{T^*}({z^{ - 1}}) - 1]p(k) + \nonumber\\ &\ B({z^{ - 1}})G({z^{ - 1}}){p^*} + \bar V(\mathit{\boldsymbol{ p}}) \end{align} $

    (22)

    其中, $ \bar V(\mathit{\boldsymbol{ p}}) $为关于$ \mathit{\boldsymbol{ p}} = {\left[ {p(k), p(k - 1), \cdots } \right]^{{\rm T}}} $和过去时刻的$ \bar V(\mathit{\boldsymbol{ p}}) $的未知非线性函数.至此, $ (k+1) $时刻的需量动态模型由式(3)和式(22)组成.

    在需量的动态模型中, 影响需量的主要因素为原料杂质成分含量$ B_{i, 1} $和颗粒长度$ B_{i, 2} $.而熔池内的气泡、电弧闪变抖动等因素对需量的干扰本文视为均值为零, 方差有界的白噪声.在生产过程中电流控制器工作在稳定状态, $ B_{i, 1} $和$ B_{i, 2} $变化引起的电流变化是有界的, 即在一个闭集中.由于电流变化是有界的, 功率也有界, 可以采用神经网络的方法对功率变化率的非线性项进行估计.在此条件下, 本文提出一种数据与模型驱动的需量预报方法.

    根据式(3)和式(22)组成的需量动态模型, 本文提出由线性模型、基于PACF输入变量决策的RBFNN未知非线性函数估计和需量计算模型组成的需量预报方法.

    需量动态模型中, 功率变化率$ \Delta p(k) $包含的线性部分$[ - z{T^*}({z^{ - 1}}) - 1]p(k) + B({z^{ - 1}})G({z^{ - 1}}){p^*} $已知, 将需量动态模型的式(22)改写为

    $ \begin{align} \label{eq23} \Delta p(k) = \Delta {p_1}(k) + \bar V(k) \end{align} $

    (23)

    其中, $ \Delta {p_1}(k) $和$ \bar V(k) $分别为$ \Delta p(k) $的线性部分和未建模动态, 即

    $ \begin{align} \Delta {p_1}(k) = &\ \left[ - z{T^*}({z^{ - 1}}) - 1\right]p(k) + \nonumber\\ &\ B({z^{ - 1}})G({z^{ - 1}}){p^*}\nonumber\\[1mm] \bar V(k) = &\ {f_v}( p(k), p(k - 1), \cdots , \nonumber\\ &\ \bar V(k - 1), \bar V(k - 2), \cdots ) \end{align} $

    (24)

    其中, 参数$ a_1$, $b_0 $由式(11)结合实验确定, 控制器参数已知, 因此可计算求得线性部分$ \Delta p_1(k) $.因为$ \bar{V} (k) $的非线性函数$ f_v(\cdot) $的输入和输出都有界, 是一个闭集, 因此可基于文献[12-13]采用RBFNN估计出, 结合$\Delta p_1(k)$得到$ \Delta \hat{p}(k) $, 再代入式(3)中即可计算出$ (k+1) $时刻需量预报值$ \hat{\bar P}(k+1)$.

    因此, 本文提出由线性模型、基于PACF输入变量决策的RBFNN未知非线性项估计和$(k+1)$时刻需量计算模型组成的预报模型结构, 如图 2所示.

    图 2  电熔镁群炉需量预报方法结构框图
    Fig. 2  The structure diagram of demand forecasting method for FMFG

    1) 线性模型.线性部分$ \Delta p_1(k) $根据工艺实验确定被控对象的线性模型参数$ a_1$, $b_0 $和控制器参数$g_0$, $g_1$, $g_2 $计算.

    2) 非线性项函数$ \bar V(k) $估计.根据采集的功率序列样本$ \{ p(i), i = 1, 2, \cdots , N\} $, 采用PACF方法分析, 决策$ \hat {\bar V}(k) $估计模型的功率历史数据输入个数$ n_f $和未建模动态历史数据输入个数$ n_v $.以过去的功率和非线性项为输入, 采用RBFNN估计出$\hat {\bar V}(k)$.

    3) $ (k+1) $时刻需量计算模型.根据检测的$ p(k) $, $ p(k-n+1) $、非线性函数估计值$\hat {\bar V}(k)$和计算的线性部分$ \Delta p_1(k) $, 求出最终的需量预报值$\hat {\bar P}(k + 1)$.

    基于上述策略, 提出需量预报算法如下:

    $ \begin{align} \label{eq25} \hat {\bar P}(k + 1) =\bar P(k) +\frac{{p(k) - p(k - n{{\rm + 1}}){{\rm + }}\Delta \hat p(k)}}{n} \end{align} $

    (25)

    其中, 功率变化率的估计$ \Delta \hat p(k) $为线性部分$ \Delta p_1(k) $和非线性函数估计值$\hat {\bar V}(k)$之和, 即

    $ \begin{align} \label{eq26} \Delta \hat p(k) = \Delta {p_1}(k) + \hat {\bar V}(k) \end{align} $

    (26)
    2.2.1   线性模型

    由式(22)可得$ \Delta p_1(k) $, 即

    $ \begin{align} \label{eq27} &\Delta {p_1}(k) =\nonumber\\ &\qquad \left[ - z{T^*}({z^{ - 1}}) - 1\right]p(k) + B({z^{ - 1}})G({z^{ - 1}}){p^*}=\nonumber\\ &\qquad\ {a_1}\left[ { - p(k) + p(k - 1)} \right] + {b_0}({g_0} + {g_1} + {g_2}){p^*} - \nonumber\\ &\qquad\ {b_0}\left[ {{g_0}p(k) + {g_1}p(k - 1) + {g_2}p(k - 2)} \right] \end{align} $

    (27)
    2.2.2   nfnv决策算法

    采用文献[3, 5]中的PACF决策$ n_f $和$ n_v $值.首先采集功率样本序列$ \{ p(i), i = 1, 2, \cdots , N\} $, 然后按文献[5]的计算步骤求得PACF函数$ \psi ({n_f}) $.当

    $ \begin{align} \label{eq28} &\psi ({n_f} - 1) < \frac{- 2}{\sqrt {N + n_f}}~\mbox{或}~\psi (n_f - 1) > \frac{2}{{\sqrt {N + n_f} }}\nonumber\\[2mm] & \frac{- 2}{\sqrt {N + n_f} } < \psi (n_f) < \frac{2}{\sqrt {N + n_f} } \end{align} $

    (28)

    则选择此时的$n_f$值.

    决策$n_v$的方法类似, 具体决策条件如下:

    $ \begin{align} \label{eq29} &\psi ({n_v} - 1) < \frac{- 2}{\sqrt {N + n_v} }~\mbox{或}~\psi (n_v - 1) > \frac{2}{{\sqrt {N + n_v} }}\nonumber\\[2mm] & \frac{- 2}{\sqrt {N + n_v} } < \psi (n_v) < \frac{2}{\sqrt {N + n_v} } \end{align} $

    (29)
    2.2.3   V(k)的估计算法

    由于$\bar V (k)$和过去时刻的$ \bar V(k - 1), \cdots, \bar V(k - {n_v}) $相关, 也和$ p(k), \cdots, p(k - {n_f} + 1) $相关.因此, $ \bar V(k) $表示以$ n_f $个过去时刻功率和$ n_v $个过去时刻非线性项为输入的未知非线性函数${f_v}( \cdot ) $, 即

    $ \begin{align} \label{eq30}\bar V(k) =&\ {f_v}\big( p(k), p(k - 1), \cdots, p(k - {n_f}{{\rm + 1}}), \notag\\ &\ \bar V(k - 1), \bar V(k - 2), \cdots, \bar V(k - {n_v}) \big) \end{align} $

    (30)

    根据文献[5], RBFNN输出当前的$ \hat{\bar V}((k)$为

    $ \begin{align} \label{eq31} &\hat {\bar V}(k) = \sum\limits_{j = 1}^H {{\phi _j}{{\rm (}} {\mathit{\boldsymbol{ d}}_j}\mathit{\boldsymbol{ x}}{{\rm )}}} \omega _j^b + \beta\nonumber\\ & {\phi _j}{{\rm (}}{\mathit{\boldsymbol{ d}}_j}\mathit{\boldsymbol{ x}}{{\rm )}} = \exp \left( { - \frac{1}{{2{\sigma ^2}}}{{\left\| {{\mathit{\boldsymbol{ d}}_j}\mathit{\boldsymbol{ x}} - {\mathit{\boldsymbol{ C}}_j}} \right\|}^2}} \right), \nonumber\\ &\qquad\qquad\qquad\qquad\qquad\qquad j = 1, 2, \cdots , H \end{align} $

    (31)

    其中, $H$为RBFNN的隐含层节点的个数, ${\mathit{\boldsymbol{ d}}_j} =$ ${\rm diag}\{\omega ^a_{1j}, \omega ^a_{2j}, \cdots , \omega ^a_{({n_f} + {n_v})j}\}$, $\omega^a _{ij}$表示输入层神经元$i$与隐含层神经元$j$之间的连接权, 全部元素$\omega ^a_{ij}$组成连接权矩阵${{{\omega}}^a}$, $n_f+n_v$行$H$列即${{{\omega}}^a}$ $\in$ ${{\bf R} ^{({n_f} + {n_v}) \times H}}$, $\mathit{\boldsymbol{ x}} = [p(k), \cdots, p(k - {n_f} + 1), \bar V(k - 1)$, $\cdots$, $\bar V(k - {n_v})]^{\rm T}$向量的维数为${n_f} + {n_v}$, $\mathit{\boldsymbol{ x}}\in{{\bf R} ^{{n_f} + {n_v}}}$, $\sigma$为高斯函数的宽度, $\left\| \cdot \right\|$为欧几里得范数, $\mathit{\boldsymbol{ C}}_j$是第$j$个高斯函数中心点, ${{\mathit{\boldsymbol{\omega}}}^{b}}$为输出层权值向量, 维数为$H$, ${\mathit{\boldsymbol{\omega}} ^b} \in {{\bf R} ^H}$, $\beta$为偏置.

    1) 输入层到隐含层的连接权矩阵${{{\omega}}^a}$的设定

    由文献[14]知${{{\omega}}^a}$对RBFNN拟合精度的作用很小, 所以将${{{\omega}}^a}$所有元素设为1, 即

    $ \begin{align} \label{eq32} \omega _{i, j}^a = 1, \quad i = 1, 2, \cdots, {n_f} + {n_v}, ~ j = 1, 2, \cdots, H \end{align} $

    (32)

    则${\phi _j}{{\rm (}}{\mathit{\boldsymbol{ d}}_j}\mathit{\boldsymbol{ x}}{{\rm ) = }}{\phi _j}{{\rm (}}\mathit{\boldsymbol{ x}}{{\rm )}}$.

    2) 隐含层节点数$H$及对应中心点$\mathit{\boldsymbol{ C}}_j$和高斯函数宽度$\sigma$的选择

    对训练样本$\{( {\bar V(k), {\mathit{\boldsymbol{ x}}_k}} ), k = 1, \cdots , N \}$进行5折交叉验证实验, 将5折后每4份数据样本数记为$M$, 即$M = \frac{4}{5}N$, 并作为隐含层中心点的候选集, 剩下$\frac{1}{5}N$的训练样本作为验证集. RBFNN对训练样本集合$\{ ( {\bar V(k), {\mathit{\boldsymbol{ x}}_k}} ), k = 1, \cdots , N\}$的回归映射表示如下[15]:

    $ \begin{align} { \mathit{\boldsymbol{\bar V}}}= \Phi {{\mathit{\boldsymbol{\omega}}}^b} + {{\bf \Xi }} \end{align} $

    (33)

    其中,

    $ \begin{align} \label{eq34} &{\mathit{\boldsymbol{\bar V}}} = {{\rm [}}\bar V(1), \bar V(2), \cdots , \bar V(N){{{\rm ]}}^{\rm T}}\nonumber\\ &\Phi = [{{\mathit{\boldsymbol{\phi}} _{1}}}, {{\mathit{\boldsymbol{\phi}} _{2}}}, \cdots , {{\mathit{\boldsymbol{\phi}} _{M}}}]\nonumber\\ &{\mathit{\boldsymbol{\phi}} _{i}} = {[{\phi _i}({\mathit{\boldsymbol{ x}}_1}), \cdots, {\phi _i} ({\mathit{\boldsymbol{ x}}_N})]^{\rm T}}, \quad 1 \le i \le M\nonumber\\ &{{\bf \Xi }} ={[\varepsilon (1), \varepsilon (2), \cdots , \varepsilon (N)]^{\rm T}} \end{align} $

    (34)

    其中, ${\mathit{\boldsymbol{\bar V}}} \in {{\bf R} ^N}$为导师信号, 即期望输出向量, ${\Phi}$ $\in$ ${{\bf R} ^{N \times {H}}}$为RBFNN隐含层输出矩阵, ${{\bf \Xi }} \in {{\bf R} ^N}$为RBFNN对训练样本建模的误差向量. $ {\phi _i}({\mathit{\boldsymbol{ x}}_j}) $为

    $ \begin{align} \label{eq35} &{\phi _i}({\mathit{\boldsymbol{ x}}_j}) = \exp \left( - \dfrac{1}{{2{\sigma ^2}}}{\left\| {{\mathit{\boldsymbol{ x}}_j} - {\mathit{\boldsymbol{ C}}_i}} \right\|^2}\right), \notag\\ &\qquad\qquad\qquad i = 1, 2, \cdots, M, ~~j = 1, 2, \cdots, N \end{align} $

    (35)

    ${\mathit{\boldsymbol{ x}}_j} = [ p(j), \cdots , p(j - {n_f} + 1), \bar V(j - 1), \cdots, \bar V(j-$ ${n_v})]^{\rm T} $, ${\mathit{\boldsymbol{ C}}_i} = [ p(i), p(i - 1), \cdots , p(i - {n_f} + 1)$, $\bar V(i-$ $1)$, $\cdots$, $\bar V(i - {n_v})]^{\rm T}$.对矩阵${\Phi}$用Gram-Schmidt方法[16]进行正交分解来逐一确定中心点, 步骤如下:

    步骤1. 对于所有$1\le i\le M$, 计算如下公式:

    $ \begin{align} \label{eq36} \begin{cases} \mathit{\boldsymbol{ q}}_1^{(i)} = {\mathit{\boldsymbol{\phi}} _{i}}\\ g_1^{(i)} = \dfrac{\left\langle {\mathit{\boldsymbol{ q}}_1^{(i)}, {{\mathit{\boldsymbol{\bar V}}}}} \right\rangle }{\left\langle {\mathit{\boldsymbol{ q}}_1^{(i)}, \mathit{\boldsymbol{ q}}_1^{(i)}} \right\rangle} \\[6mm] {[err]}_1^{(i)} = g_1^{(i)}\dfrac{\left\langle {\mathit{\boldsymbol{ q}}_1^{(i)}, {\mathit{\boldsymbol{\bar V}}}} \right\rangle }{\left\langle {{{\mathit{\boldsymbol{\bar V}}}}, {{\mathit{\boldsymbol{\bar V}}}}} \right\rangle}\\ \end{cases} \end{align} $

    (36)

    得到$M$个候选样本作为中心点的贡献程度指标$[{ err}]_1^{(i)}$, 最大的贡献度指标为

    $ \begin{align} \label{eq37} [{ err}]_1^{({i_1})} = \max \{ [{ err}]_1^{(i)}, 1 \le i \le M\} \end{align} $

    (37)

    选择

    $ \begin{align} \label{eq38} \mathit{\boldsymbol{ q}}_1 = \mathit{\boldsymbol{ q}}_1^{({i_1})} = {\mathit{\boldsymbol{\phi}} _{{i_1}}} \end{align} $

    (38)

    以上标$i_1$对应的样本${\mathit{\boldsymbol{ x}}_{{i_1}}}$作为第1个中心点

    $ \begin{align} \label{eq39} {\mathit{\boldsymbol{ C}}_1} = {\mathit{\boldsymbol{ x}}_{{i_1}}} \end{align} $

    (39)

    步骤H (2 ≤ HHmax). 对于所有的$1\le i$ $\le$ $M$, $i \ne {i_1}, \cdots, i \ne {i_{H - 1}}$, 其中, ${i_1}, \cdots, {i_{H - 1}}$为前面$H - 1$步已经选择作为中心点的样本下标.计算如下公式:

    $ \begin{align} \label{eq40} \begin{cases} l_{jH}^{(i)} = \dfrac{\left\langle {{\mathit{\boldsymbol{ q}}_j}, {\mathit{\boldsymbol{\phi}} _{i}}} \right\rangle } {\left\langle {{\mathit{\boldsymbol{ q}}_j}, {\mathit{\boldsymbol{ q}}_j}} \right\rangle }, \quad 1 \le j \le H\\[4mm] \mathit{\boldsymbol{ q}}_H^{(i)} = {\mathit{\boldsymbol{\phi}} _{i}} - \sum\limits_{j = 1}^{H - 1} {l_{jH}^{(i)}} {\mathit{\boldsymbol{ q}}_j}\\[4mm] g_H^{(i)} = \dfrac{\left\langle {\mathit{\boldsymbol{ q}}_H^{(i)}, {\mathit{\boldsymbol{\bar V}}}} \right\rangle } {\left\langle {q_H^{(i)}, q_H^{(i)}} \right\rangle }\\[6mm] {[{ err}]}_H^{(i)} = g_H^{(i)}\dfrac{\left\langle {\mathit{\boldsymbol{ q}}_H^{(i)}, {\mathit{\boldsymbol{\bar V}}}} \right\rangle }{\left\langle {{\mathit{\boldsymbol{\bar V}}}, {\mathit{\boldsymbol{\bar V}}}} \right\rangle} \end{cases} \end{align} $

    (40)

    $M$个候选样本中的贡献度指标第$H$大的值为

    $ \begin{align} \label{eq41}& [{ err}]_H^{({i_H})} = \nonumber\\ &\quad\max \left\{ [{ err}]_H^{(i)}, 1 \le i \le M, i \ne {i_1}, \cdots, {i_{H - 1}}\right\} \end{align} $

    (41)

    选择

    $ \begin{align} \label{eq42}{\mathit{\boldsymbol{ q}}_H} = \mathit{\boldsymbol{ q}}_H^{({i_H})} = {\mathit{\boldsymbol{\phi}} _{\;{i_H}}} - \sum\limits_{j = 1}^{H - 1} {l_{jH}^{({i_H})}} {\mathit{\boldsymbol{ q}}_j} \end{align} $

    (42)

    以上标${i_H}$对应的样本${\mathit{\boldsymbol{ x}}_{{i_H}}}$作为第$H$个中心点

    $ \begin{align} \label{eq43} {\mathit{\boldsymbol{ C}}_H} = {\mathit{\boldsymbol{ x}}_{{i_H}}} \end{align} $

    (43)

    当$H = {H_{\max }}$时, ${H_{\max }}$个中心点选择完毕.

    节点数$H$逐一增加, 对候选集拟合的误差越来越小, 但当$H$增至某临界值$H_c$值时, 验证集拟合的误差开始增大, 出现过拟合现象.因此权衡预报精度和避免过拟合的要求, 选择节点数$H$.

    $ \begin{align} \label{eq44} H = {H_c}, \quad 1 < {H_c} < {H_{\max }} \end{align} $

    (44)

    当节点数$H$确定之后, 对完整训练样本进行正交分解得到贡献度较大的前$H$个中心点.

    $ \begin{align} \label{eq45} \{ {\mathit{\boldsymbol{ C}}_j}, j = 1, \cdots, H\} \end{align} $

    (45)

    隐层节点的高斯函数采用相同的$\sigma$, 采用文献[17]的方法, 在一个合理区间$[{\sigma _{\min }}, {\sigma _{\max }}]$内进行5折交叉验证实验, 选择验证集预报误差收敛速度变化不明显的临界值${\sigma _c}$.

    $ \begin{align} \label{eq46} \sigma = {\sigma _c}, \quad {\sigma _{\min }} < {\sigma _c} < {\sigma _{\max }} \end{align} $

    (46)

    3) 输出层权重向量${{\mathit{\boldsymbol{\omega}}}^b}$和偏置$\beta$的求解

    采用最小二乘法求如下方程:

    $ \begin{align} \label{eq47} \left[ {\Phi \;\mathit{\boldsymbol{{s}}}} \right]\left[ {\begin{array}{*{20}{c}} {{{\mathit{\boldsymbol{\omega}}}^b}}\\ \beta \end{array}} \right] = {\mathit{\boldsymbol{\bar V}}} \end{align} $

    (47)

    其中, $\mathit{\boldsymbol{ s}} = {[1, \cdots , 1]^{\rm T}}$, $\mathit{\boldsymbol{ s}} \in {{\bf R} ^N}$, 求得${{\mathit{\boldsymbol{\omega}}}^b}$和$\beta$.

    $ \begin{align} \label{eq48} \left[ {\begin{array}{*{20}{c}} {{{\mathit{\boldsymbol{\omega}}}^b}}\\ \beta \end{array}} \right] = {\left[ {\Phi \;\mathit{\boldsymbol{ s}}} \right]^† }{\mathit{\boldsymbol{\bar V}}} \end{align} $

    (48)

    其中, ${\left[ \cdot \right]^† }$为广义逆矩阵, 由奇异值分解方法解得.至此, $ \hat{\bar V}(k) $由RBFNN的估计算法如下:

    $ \begin{align} \label{eq49} \hat {\bar V}(k) = &\ \sum\limits_{j = 1}^H {\exp \left( { - \frac{1}{{2{\sigma ^2}}} {{\left\| {\mathit{\boldsymbol{ x}} - {\mathit{\boldsymbol{ C}}_j}} \right\|}^2}} \right)} \omega ^b_j + \beta\nonumber\\ \mathit{\boldsymbol{ x}} = &\ \left[p(k), \cdots , p(k - {n_f} + 1), \bar V(k - 1), \cdots, \right. \nonumber\\ &\ \left. \bar V(k - {n_v}) \vphantom{p(k), \cdots , p(k - {n_f} + 1), \bar V(k - 1), } \right]^{{\rm T}} \end{align} $

    (49)

    其中, $n_f$和$n_v$由式(28)和式(29)给出, $H$, $\mathit{\boldsymbol{ C}}_j$和$ \sigma $由式(44)~(46)给出, ${{\mathit{\boldsymbol{\omega}}}^b}$和$\beta$由式(48)给出.

    4) $\bar V(k)$预报模型的在线更新

    由于群炉生产状态会发生迁移, 需根据当前预报误差对预报模型的参数${\mathit{\boldsymbol{\omega}}^b}$进行更新.设${{\mathit{\boldsymbol{\omega}}}^b}$在$k$时刻为${\mathit{\boldsymbol{\omega}}}^b_k$, 求$(k+1)$时刻的权值向量$\mathit{\boldsymbol{\omega }}^b_{k + 1}$.

    $ \begin{align} \label{eq50} \mathit{\boldsymbol{\omega }}^b_{k + 1} = \mathit{\boldsymbol{\omega }}^b_k + \Delta \mathit{\boldsymbol{\omega }}^b_k \end{align} $

    (50)

    采用文献[18]递推正交最小二乘法的更新步骤如下:

    步骤1. $k$时刻${\mathit{\boldsymbol{\omega}}^b}$的值为$\mathit{\boldsymbol{\omega}}^b_k$, 隐含层输出矩阵为${{\Phi} _k}$, 对${{\Phi} _k}$正交分解

    $ \begin{align} \label{eq51} {\Phi _k} = {{Q}_k}\left[ {\begin{array}{*{20}{c}} {{{R}_k}}\\ {O} \end{array}} \right] \end{align} $

    (51)

    其中, ${{Q}_k} \in {{\bf R} ^{N \times N}}$是正交矩阵, ${{R}_k} \in {{\bf R} ^{{H} \times {H}}}$是上三角矩阵, ${O} \in {{\bf R} ^{(N - H) \times {H}}}$是零矩阵.

    步骤2. $k+1$时刻递推正交分解求得${{R}_{k + 1}}$.

    $ \begin{align} \label{eq52} \begin{array}{l} \left[ {\begin{array}{*{20}{c}} {{{R}_k}}\\ {\mathit{\boldsymbol{\vartheta}}_{\;k + 1}^{\rm T}} \end{array}} \right] ={{Q}_{k + 1}}\left[ {\begin{array}{*{20}{c}} {{{R}_{k + 1}}}\\ 0 \end{array}} \right] \end{array} \end{align} $

    (52)

    其中, $\mathit{\boldsymbol{\vartheta}} _{\;k + 1}^{\rm T} = [{\phi _1}({\mathit{\boldsymbol{ x}}_{k{{\rm + }}1}}), \cdots, {\phi _H}({\mathit{\boldsymbol{ x}}_{k{{\rm + }}1}})]$.再由${Q}_{k + 1}$计算${ \Omega} _{k + 1}$

    $ \begin{align} \label{eq53} \begin{array}{l} \left[ {\begin{array}{*{20}{c}} {{ \Omega} _{k + 1}}\\ {\Delta \tilde {\bar V}(k)} \end{array}} \right] ={Q}_{k + 1}^{\rm T}\left[ {\begin{array}{*{20}{c}} {{\rm 0}}\\ {\varepsilon (k + 1)} \end{array}} \right] \end{array} \end{align} $

    (53)

    其中, $\varepsilon (k{{\rm + }}1) =\bar V(k + 1{{\rm ) - }}\mathit{\boldsymbol{\vartheta}} _{\;k + 1}^{\rm T}\mathit{\boldsymbol{\omega }}^b_k$.

    步骤3. 由最小二乘法求得$\Delta \mathit{\boldsymbol{\omega }}^b_k$.

    $ \begin{align} \label{eq54} \Delta \mathit{\boldsymbol{\omega }}^b_k = {\left[ {{R_{k + 1}}} \right]^† }{{ \Omega} _{k + 1}} \end{align} $

    (54)

    代入式(50)得到$\mathit{\boldsymbol{ \omega }}^b_{k + 1}$, 更新$k=k+1$, 返回步骤2.

    综上所述, 本文预报方法由式(25)给出, 其中功率变化率预报值$ \Delta \hat p (k) $由式(26)给出.式(26)中的线性部分$ \Delta p_1 (k)$由式(27)给出. $\hat {\bar V} (k)$为$\Delta p (k)$的未建模动态的神经网络估计, 由式(49)给出.

    将本文所提预报方法对某电熔镁砂厂1号电力变压器的1个炉次的群炉需量数据进行仿真实验, 以验证方法的有效性.实验数据中功率的采样周期为7秒, 需量参数$ n=30 $, 群炉生产台数为4台, 当天设定的需量限幅值为22 100 kVA.从该炉次4 620组数据中选取4 000组数据$p(i)$, $i = 1, \cdots, 4 000 $进行仿真实验, 其中$p(i)$, $i = 1, \cdots, 2 000$为训练集, $p(i)$, $i = 2 001, \cdots, 4 000$为测试集.根据下式求得建模所用的$ \Delta p(i)$

    $ \begin{align} \label{eq55} \Delta p(i) = p(i + 1) - p(i) \end{align} $

    (55)
    3.1.1   线性模型∆p1(i)

    根据现场工艺实验在工作点$ y^*=15 000 $ A, $ u^*$ $=$ $6.1918\times10^{-4} $ rad $\cdot $ s$^{-1}$处线性化被控对象模型, 计算得到参数$ a_1$, $b_0 $如下:

    $ \begin{align} \label{eq56} {a_1} = - 1.01, \quad{b_0} = 0.1 \end{align} $

    (56)

    控制器参数$ g_0 = 6.300 035 $, $g_1 = -11.9$, $g_2 = 5.6$, 则功率变化率线性模型$ \Delta p_1(i) $计算公式为

    $ \begin{align} \label{eq57} \Delta {p_1}(i) =&\ {a_1}\left[ { - p(i) + p(i - 1)} \right] + \nonumber\\ &\ {b_0}({g_0} + {g_1} + {g_2}){p^*} - \nonumber\\ &\ {b_0}\left[ {{g_0}p(i) + {g_1}p(i - 1) + {g_2}p(i - 2)} \right] \end{align} $

    (57)

    其中, ${p^*} = \sqrt 3 mU\cos \varphi {y^*}$, $m = 4$, $U = 190$, $\cos \varphi$ $=$ $0.92$.

    3.1.2   nfnv的决策

    对2 000组训练数据组成的功率时间序列求PACF序列值为

    $ \begin{align} \label{eq58} &\psi (1) = 0.9402, \quad~~ \psi (2) = - 0.2357\nonumber\\ &\psi (3) = 0.1486, \quad~~ \psi (4) = - 0.0050\nonumber\\&\psi (5) = - 0.0078, \quad\psi (6) = 0.0068\nonumber &\vdots \notag\\ &\psi (10) = 0.0070 \end{align} $

    (58)

    PACF的95 %的置信区间为$[ - {{\rm 0}}{{\rm .0448}}, {{\rm 0}}{{\rm .0448}}]$, 当$n_f$ $=4$时, PACF值落入95 %置信区间内, 因此设定

    $ \begin{align} \label{eq59} n_f = 4 \end{align} $

    (59)

    类似地, 对2 000组训练数据组成的未建模动态求PACF序列值为

    $ \begin{align} \label{eq60} &\psi (1) = 0.9269, \quad\psi (2) = - 0.9471\nonumber\\ & \psi (3) = 0.8522, \quad\psi (4) = - 0.6768\nonumber\\ & \psi (5) = 0.2228, \quad\psi (6) = - 0.2071\nonumber\\ & \psi (7) = 0.0239, \quad\quad\cdots\notag \\ &\qquad\qquad\vdots\notag\\ & \psi (10) = 0.0278 \end{align} $

    (60)

    PACF的95 %的置信区间为$[ - {{\rm 0}}{{\rm .0448}}, {{\rm 0}}{{\rm .0448}}]$, 设定

    $ \begin{align} \label{eq61} n_v = 7 \end{align} $

    (61)
    3.1.3   基于RBFNN的V(i)估计

    根据下式整理RBFNN训练所需的输入和输出数据:

    $ \begin{align} \label{eq62} \hat {\bar V}(i) = \sum\limits_{j = 1}^H {{\phi _j}{{\rm (}}\mathit{\boldsymbol{ x}}{{\rm )}}} \omega^b _j + \beta \end{align} $

    (62)

    其中, $ {\bf x} = [p(i), \cdots , p(i - {n_f} + 1), \bar V(i -1), \cdots$, $\bar V(i$ - ${n_v})]^{{\rm T}}$.先将$ \sigma $由经验选为0.8, 将节点数$ H $由1逐一增加至$ H_{\max} = 200 $.

    图 3所示, 进行5折交叉验证实验, 可以看出一开始随着节点个数增加, 候选集和验证集的预报误差指标(MAPE, RMSE)都在减小, 但是当节点数$ H >60$时, 虽然候选集的预报误差指标依然在减小, 但是验证集的预报误差指标开始增大.因此根据实验结果选择隐含层节点数

    图 3  隐层节点数的交叉验证
    Fig. 3  Cross-validation of the number of hidden nodes

    $ \begin{align} \label{eq63} H = 60 \end{align} $

    (63)

    图 4所示, 在区间范围$[{\sigma _{\min}}, {\sigma _{\max }}] =$ $[0.01$, $2.00]$内, 通过5折交叉验证实验, 随着$ \sigma $增大候选集和验证集的预报误差指标都相应减小, 当$ \sigma >1.2$时预报误差指标减小的幅度很小.因此根据实验结果选择高斯函数宽度

    图 4  高斯函数宽度的交叉验证
    Fig. 4  Cross-validation of the width of Gaussian function

    $ \begin{align} \label{eq64} \sigma = 1.2 \end{align} $

    (64)

    确定$H$和$\sigma$之后, 对整个训练集数据进行正交分解, 得到$H$个中心点$ \mathit{\boldsymbol{ C}}_j $的初始值${\mathit{\boldsymbol{ C}}_j}(0)$, $j =1$, $2$, $\cdots$, $60 $, 即

    $ \begin{align} \label{eq65} &\mathit{\boldsymbol{ C}}_1(0) = [0.3231, 0.4961, \cdots , 0.8833, 0.9910]^{\rm T}\notag\\ &\qquad\qquad\qquad\vdots \notag \\ &\mathit{\boldsymbol{ C}}_{60}(0) = [0.8934, 0.8902, \cdots , 0.5330, 0.5468]^{\rm T} \end{align} $

    (65)

    由式(48)求出权值向量${{\mathit{\boldsymbol{\omega}}}^b}$的初始值${{\mathit{\boldsymbol{\omega }}}^b(0)}$和偏置$\beta$初始值$\beta(0)$为

    $ {\mathit{\boldsymbol{\omega}}}^b(0) = [2.1270, \;7.1224, \cdots, - 1 295.4]^{\rm T} $

    (66)

    $ \beta(0) = {{\rm - 0}}{{\rm .0057}} $

    (67)
    3.1.4   ∆p(k)预报算法验证

    将离线训练后的预报初模型对剩余的2 000组数据进行验证, 如图 5所示.可以看出功率变化率的预报值基本包括功率变化率中有规律的信息.功率变化率预报误差为$\Delta {e_p}(k) = \Delta p(k) - \Delta \hat p(k)$, 可以看出误差$ \Delta {e_p}(k) $序列明显呈现一个白噪声的特性.统计的$\Delta p(k)$预报误差指标如表 1所示, 预报误差序列的方差为1.1481E+6, 均方根误差为1 071.3.

    图 5  p(k)预报验证曲线
    Fig. 5  Forecast validation curves of ∆p(k)
    表 1  p(k)预报误差指标
    Table 1  Forecast error indicators of ∆p(k)
    方差 RMSE
    1.1481E+6 1 071.3
    下载: 导出CSV 
    | 显示表格
    3.1.5   P(k + 1)预报算法验证

    验证实验中需量有三次较明显的先上升后下降的趋势, 分别是: $k=171$~$220$时段, 如图 6所示; $k$ $=571$~$620$时段, 如图 7所示; $k=1 356$~$1 405$时段, 如图 8所示.

    图 6  时段1需量曲线
    Fig. 6  Demand curve for the 1st time period
    图 7  时段2需量曲线
    Fig. 7  Demand curve for the 2nd time period
    图 8  时段3需量曲线
    Fig. 8  Demand curve for the 3rd time period

    图 6中, $k=171$时, 需量实际值为$\bar P(k)=21 427$ kW, 对应需量预报值$\hat{\bar P}(k)=21 421$ kW; 当$k$ $=187$时, 需量上升至该时段内最大值$\bar P(k)=21 833$ kW, 对应需量预报值$\hat{\bar P}(k)=21 842$ kW; 之后由于电流控制系统的调节作用, $\bar P(k)$开始下降, 当$k=194$时, 需量降低至$\bar P(k)=21 598$ kW, 对应需量预报值为$\hat{\bar P}(k)=21 620$ kW.

    图 7中, $k=571$时, 需量实际值为$\bar P(k)=21 678$ kW, 对应需量预报值$\bar P(k)=21 687$ kW; $k$ $=$ $592$时, 需量为该时段最大值$\bar P(k)=22 067$ kW, 需量预报值$\hat{\bar P}(k)=22 080$ kW; 之后$\bar P(k)$开始下降, $k=620$时, 需量$\bar P(k)=21 730$ kW, 对应需量预报值为$\hat{\bar P}(k)=21 704$ kW.

    图 8中, $k=1 356$时, 需量实际值为$\bar P(k)=21 450$ kW, 对应需量预报值$\hat{\bar P}(k)=21 445$ kW; $k$ $=1 385$时, 需量为该时段最大值$\bar P(k)=22 011$ kW, 需量预报值$\hat{\bar P}(k)=22 006$ kW; 之后$\bar P(k)$开始下降, $k=1 399$时, 需量$\bar P(k)=21 667$ kW, 对应需量预报值为$\hat{\bar P}(k)=21 633$ kW.

    图 6~8可以看出, 在三个时段中需量实际值和需量预报值的上升和下降趋势基本类似.为定量分析预报方法的性能, 本文采用误差方差、预报精度百分比(Percent better, PB)[19]、均方根误差(Root mean square error, RMSE)[20-23]、平均绝对百分误差(Mean absolute percentage error, MAPE)[3]作为评估预报性能的指标, 计算公式如下:

    $ \begin{align} \label{eq68} \begin{cases} {\rm PB} = \dfrac{{\sum\limits_{k = 1}^N {{j_k}} }}{N} \times 100 \%, \\[2mm] \qquad {j_k} ={\small \begin{cases} 1, &\mbox{若}\left| {\hat {\bar P}(k) - \bar P(k)} \right| < \left| {{{\bar P}_{0.5}}(k) - \bar P(k)} \right|\\[0.4em] 0, &\mbox{若}\left| {\hat {\bar P}(k) - \bar P(k)} \right| \ge \left| {{{\bar P}_{0.5}}(k) - \bar P(k)} \right| \end{cases}}\\[7mm] {\rm RMSE} = \sqrt {\dfrac{1}{N}\sum\limits_{k = 1}^N {{{\left(\hat {\bar P}(k) - \bar P(k)\right)}^2}} } \\[7mm] {\rm MAPE} = \dfrac{1}{N}\sum\limits_{k = 1}^N \left| {\dfrac{{\hat {\bar P}(k) - \bar P(k)}}{{\bar P(k)}}} \right|\times 100 \% \end{cases} \end{align} $

    (68)

    其中, $ {\bar P_{0.5}}(k) = (1 + 0.5 \% ) \times \bar P(k) $, $N$为误差序列的样本数量.对包含图 6~8时间段测试集2 000个采样点的预报实验误差序列由式(68)进行计算, 对应预报性能指标见表 2, 误差序列的方差为1 275.7, 预报精度百分比为97.55 %, RMSE为35.7104, MAPE为0.1054 %, 与文献[5]中预报方法对比, 本文的方法明显减小了预报误差的方差和RMSE, 提高了预报精度百分比.由式(25)可知, 需量预报算法的复杂度和精度取决于$ \Delta p(k) $的估计.由于利用了所建立的$\Delta p(k)$动态模型式(22), 将$\Delta p(k)$的估计分为可精确计算的线性部分$ \Delta p_1(k) $和未建模动态$ \bar V(k) $.其中, $ \Delta p_1(k) $由式(27)计算, 未知的$ \bar V(k) $由式(49)采用神经网络估计.因此, 与采用神经网络整体估计$\Delta p(k)$的文献[5]相比, 本文所提方法对线性部分利用已知信息进行精确四则运算, 没有增加神经网络估计算法的复杂度, 同时使预报精度明显提高, 满足了现场预报的要求.

    表 2  需量预报误差指标
    Table 2  Forecast error indicators of demand
    方法方差PB (%) RMSEMAPE (%)
    文献[5]1 533.097.0539.19210.0979
    本文 1 275.797.5535.71040.1054
    下载: 导出CSV 
    | 显示表格

    将提出的预报方法应用于某电熔镁砂厂1号电力变压器负载的电熔镁群炉需量监控过程, 如图 1所示.该过程主要由1号电力变压器(型号: SF9-22500/66, 额定容量: 22 500 kVA), 高压断路器(10 kV, DJS-10机械闭锁), 5台电炉变压器(型号: HKS-4500/10, 额定容量: 4 500 kVA), 5台电熔镁炉(直径2.7 m, 高3 m), 升降电机(型号: YVP160M-4, 标称功率11 kW), PID电流控制器(CPU型号: 313-6CF03-0AB0).当天设定需量限幅值为21 800 kW, 群炉生产台数为4台.该厂对需量的管控动作分为切断和恢复两种.当需量实际值超过需量限幅值且超过时间大于4个采样周期(28秒)时, 进行切断操作.当需量实际值低于需量限幅值时, 对该台电熔镁炉进行恢复供电的操作.

    3.2.1   需量预报系统简介

    采用本文所提的预报方法研制了需量预报系统, 硬件平台为图 1中的研华IPC-7120需量监控计算机以及Siemens CP5621通讯板卡.

    需量预报软件平台包括: STEP7-Micro/WIN编程软件, PC Access OPC服务器软件, SIMATIC WinCC Explorer过程监视软件, jdk1.7版本的JAVA软件开发工具包, JAVA集成开发环境Eclipse各1套以及Windows 7操作系统.基于上述软件平台和本文所提算法研制的需量预报软件的界面如图 9所示, 其中实线表示需量实际值, 虚线表示需量预报值.

    图 9  需量预报软件界面
    Fig. 9  The interface of demand forecasting software
    3.2.2   预报模型的参数选择

    预报模型的参数为仿真实验结束时的预报模型参数, 即模型输入变量个数${n_f} = 4$, ${n_v} = 7$, 隐含层节点数$H=60$, 高斯函数宽度$\sigma = 1.2$, 其他参数值具体如下, 中心点${ {\mathit{\boldsymbol{ C}}}_j}$, $j = 1, 2, \cdots, 60$为

    $ \begin{align} \label{eq69} &\mathit{\boldsymbol{ C}}_1(0) = [0.3231, 0.4961, \cdots, 0.9910]^{\rm T}\notag\\ &\qquad\qquad \qquad \vdots\notag\\ &\mathit{\boldsymbol{ C}}_{60}(0) = [0.8934, 0.8902, \cdots, 0.5468]^{\rm T} \end{align} $

    (69)

    输出层的权值向量${{\mathit{\boldsymbol{\omega}}}^b}$为

    $ \begin{align} \label{eq70} {{\mathit{\boldsymbol{\omega}}}^b} = [ - 2.5886, \;7.6589, \cdots, - 879.1062]^{\rm T} \end{align} $

    (70)

    偏置$\beta$为

    $ \begin{align} \label{formula62} \beta = 0.2957 \end{align} $

    (71)
    3.2.3   实验结果

    在时段21 : 00 : 00~05 : 22 : 29进行的工业实验中, 正常情况下需量实际值和预报值曲线如图 10~12所示.

    图 10  时段1需量曲线
    Fig. 10  Demand curve for the 1st time period
    图 11  时段2需量曲线
    Fig. 11  Demand curve for the 2nd time period
    图 12  时段3需量曲线
    Fig. 12  Demand curve for the 3rd time period

    在时段1 (22 : 48 : 23~22 : 54 : 13)中, $\bar P(k)$由22 : 48 : 23的20 801 kW上升至22 : 51 : 11的21 689 kW, $\hat{\bar P}(k)$由22 : 48 : 23的20 810 kW上升至22 : 51 : 11的21 690 kW; 之后$\bar P(k)$下降至22 : 54 : 13的21 104 kW, $\hat{\bar P}(k)$降至21 123 kW, 如图 10所示.

    在时段2 (23 : 18 : 43~23 : 24 : 33)中, $\bar P(k)$由23 : 18 : 43的21 253 kW上升至23 : 21 : 10的21 783 kW, $\hat{\bar P}(k)$由23 : 18 : 43的21 247 kW上升至23 : 21 : 10的21 786 kW; 之后$\bar P(k)$下降至23 : 24 : 33的21 041 kW, $\hat{\bar P}(k)$降至21 049 kW, 如图 11所示.

    在时段3 (00 : 46 : 13~00 : 52 : 03)中, $\bar P(k)$由00 : 46 : 13的21 122 kW上升至00 : 49 : 08的21 780 kW, $\hat{\bar P}(k)$由00 : 46 : 13的21 127 kW上升至00 : 49 : 08的21 813 kW; 之后$\bar P(k)$下降至00 : 52 : 03的21 375 kW, $\hat{\bar P}(k)$降至21 369 kW, 如图 12所示.

    超限拉闸情况下, 需量实际值和预报值曲线如图 13所示, 在22 : 38 : 00为21 826 kW时出现了拉闸操作.拉闸前后一段时间(10个采样周期)需量实际值与需量预报值及误差见表 3.

    图 13  超限拉闸时段需量曲线
    Fig. 13  Demand curve for cut off time period
    表 3  超限拉闸时段需量预报误差
    Table 3  Demand forecast errors during cut off time period
    时间需量实际值(kW)需量预报值(kW)误差(kW)
    22 : 36 : 5021 62621 635-9
    22 : 36 : 5721 65421 6504
    22 : 37 : 0421 69221 6857
    22 : 37 : 1121 72821 7253
    22 : 37 : 1821 75421 7514
    22 : 37 : 2521 78821 7871
    22 : 37 : 3221 81221 829-17
    22 : 37 : 3921 83421 849-15
    22 : 37 : 4621 83521 839-4
    22 : 37 : 5321 83321 8294
    22 : 38 : 0021 82621 8224
    22 : 38 : 0721 69121 819-128
    22 : 38 : 1421 51021 46347
    22 : 38 : 2121 33521 31817
    22 : 38 : 2821 19121 17615
    22 : 38 : 3521 04921 056-7
    22 : 38 : 4220 90820 9071
    22 : 38 : 4920 75120 769-18
    22 : 38 : 5620 58120 5774
    22 : 39 : 0320 40720 38522
    下载: 导出CSV 
    | 显示表格

    表 3看出, 在22 : 36 : 50~22 : 37 : 53拉闸前半段时间内需量实际值为上升趋势, 对应的需量预报值同样为上升趋势, 趋势相同.在22 : 37 : 32时$\bar P(k)$为21 812 kW, 开始超出需量限幅值, 对应$\hat {\bar P}(k)$为21 829 kW; 22 : 37 : 32之后连续4次采样时刻需量为21 834 kW, 21 835 kW, 21 833 kW, 21 826 kW, 对应的预报值为21 849 kW, 21 839 kW, 21 829 kW, 21 822 kW, 拉闸前4次采样时刻的需量实际值和需量预报值的趋势相同, 工厂在22 : 38 : 00进行拉闸操作, 拉闸时的$\bar P(k)$为21 826 kW, $\hat {\bar P}(k)$为21 822 kW; 从22 : 38 : 07~22 : 39 : 03拉闸后半段时间$\bar P(k)$下降迅速, 对应的$\hat {\bar P}(k)$也开始下降.

    上述需量变化过程中, 虽然需量实际值超过了需量限幅值, 但没有超过需量峰值.从图 13表 3可以看出, 需量预报值在拉闸时刻附近趋于缓慢下降趋势, 如果不拉闸需量实际值有可能会下降到限幅值以下, 那么这次拉闸就可能是不必要的拉闸.因此, 将本文所提预报方法与现有需量监控系统结合将有助于减少尖峰引起的不必要拉闸, 对提高需量限幅值的设定值, 从而提高生产过程的电能利用率有一定的指导作用.

    超限拉闸时段的需量预报误差变化曲线如图 14所示, 由于预报方法是基于拉闸时和之前的功率数据来预报, 预报的是不拉闸情况的下一时刻需量, 所以在拉闸后下一时刻需量预报误差与之前预报误差相比发生了较大的波动变化; 随着时间推移, 预报模型的输入变量开始包含拉闸之后的功率数据, 在经过2个采样时刻之后, 需量预报误差又逐渐减小, 说明预报模型可以及时地随着需量动态变化而在线调整.

    图 14  超限拉闸时段需量预报误差变化曲线
    Fig. 14  Demand forecast error curve during cut off time period

    根据现场实验的数据, 对恢复动作时间段的需量预报进行了仿真实验, 如图 15所示.图中$a$处22 : 39 : 17时刻对断电的电熔镁炉进行恢复供电的动作; 图中$b$处22 : 40 : 48时刻该台电熔镁炉由于刚恢复供电炉内工况不稳, 电弧闪灭出现了跳闸.图中$c$处22 : 43 : 22该台电熔镁炉再次恢复供电.

    图 15  恢复供电动作下的需量预报
    Fig. 15  Demand forecast curve during restore operations

    在$a$, $c$两处恢复供电的动作之后的需量预报误差变化如图 16所示. $a$处首次恢复供电, 由于该炉停止供电一段时间造成炉内熔池温度下降, 所以恢复供电后的熔炼电流值要低于切断之前的电流值.所以需量会有上升但不会达到断电之前的需量值.当3个采样周期之后需量预报误差减小至正常范围内.在$b$处跳闸之后需量又开始突然下降, 需量预报误差经过2个采样周期之后减小至正常范围.在$c$处再次恢复供电, 需量开始减缓下降并逐渐开始上升, 需量预报误差经过3个采样周期减小至正常范围.说明对于恢复供电的动作, 本文提出的预报方法的预报误差也会产生较大波动, 但随着时间推移预报误差又减小到正常范围内.

    图 16  恢复供电动作下的需量预报误差
    Fig. 16  Demand forecast error curve during restore operations

    对上述21 : 00 : 00~05 : 22 : 29时间段的实验结果进行性能分析, 实验的需量预报误差的概率分布和白度分析如图 17所示, 可以看出需量预报误差大部分在$ [-70, 70] $之内, 误差序列的自相关系数大部分在白度测试[24]的95 %置信区间$ [-0.0299$, $0.0299] $之内, 因此可将预报误差视为白噪声序列, 说明了需量预报值的可靠性.

    图 17  工业实验需量预报误差白度分析
    Fig. 17  The whiteness analysis of demand forecast error in industrial experiment

    根据式(68)计算得到工业实验的需量预报方法性能指标见表 4, 方差为1 049.8, 预报精度百分比为97.68 %, RMSE为32.3981, MAPE为0.0996 %.从表 4可以看出, 采用本文提出的需量预报方法能够较为准确地预报出下一时刻的需量值, 对提高需量限幅值的设定值具有指导意义.

    表 4  工业实验需量预报误差指标
    Table 4  Demand forecast error indicators of industrial experiment
    方差PB (%) RMSEMAPE
    1 049.897.6832.39810.0996
    下载: 导出CSV 
    | 显示表格

    本文建立了$(k+1)$时刻电熔镁群炉需量模型, 利用$(k+1)$时刻需量$ \bar P(k) $取决于功率变化率$ \Delta p(k) $, 以及$\Delta p(k)$取决于电流控制系统输出电流的特点, 提出了电熔镁群炉需量预报方法.该方法由线性模型、基于PACF输入变量个数决策的RBFNN未知非线性函数估计和$(k+1)$时刻需量$\bar P(k+1)$计算模型组成.通过某电熔镁砂厂实际数据的仿真实验和工业实验表明, 该方法可准确预报需量变化趋势, 不仅对预报需量尖峰防止不必要拉闸有实际意义, 而且对于工业过程控制系统的运行指标的预报具有一定参考价值.


  • 本文责任编委 钟麦英
  • 图  1  ZDJ9型转辙机的控制电路

    Fig.  1  Control circuit of ZDJ9 type switch machine

    图  2  基于模糊字典法的故障诊断流程图

    Fig.  2  Flow chart of fault diagnosis based on fuzzy dictionary

    图  3  8个特征量的隶属度函数分布图

    Fig.  3  Distribution of membership functions of eight characteristic quantities

    图  4  基于BP神经网络的故障诊断模型

    Fig.  4  Fault diagnosis based on BP neural network

    图  5  基于模糊字典法的故障诊断测试结果

    Fig.  5  Fault diagnosis test results based on fuzzy dictionary method

    图  6  基于BP神经网络的故障诊断测试结果

    Fig.  6  Fault diagnosis test results based on BP neural network

    图  7  基于支持向量机的故障诊断测试结果

    Fig.  7  Fault diagnosis test results based on SVM

    图  8  基于群决策的故障诊断测试结果

    Fig.  8  Fault diagnosis test results based on group decision making

    图  9  四种方法的诊断成功率对比

    Fig.  9  Comparison of diagnostic success rates between the four methods

    表  1  道岔控制电路故障

    Table  1  Fault of turnout control circuit

    ID描述
    $A0$无故障
    $A1$室外X1支路断线
    $A2$室内1DQJ断线
    $A3$室内1DQJF断线
    $A4$ $R_1$开路
    $A5$室内表示继电器断线
    $A6$室外继电器支路开路
    $A7$室外二极管支路击穿
    $A8$室外二极管支路开路
    $A9$整流匣短路
    $A10$V线圈开路
    下载: 导出CSV

    表  2  ZDJ9型道岔控制电路故障字典

    Table  2  Fault dictionary for ZDJ9 turnout control circuit

    类型 $B1$ $B2$ $B3$ $ B4$ $ B5$ $ B6$ $ B7$ $ B8$
    $A0$50215722572200
    $A1$000011001100
    $A2$00000000
    $A3$15010000000
    $A4$001100110000
    $A5$402000697500
    $A6$80025025020
    $A7$2501050105000
    $A8$40200070757075
    $A9$1040306033
    $A10$663800730730
    下载: 导出CSV

    表  3  道岔控制电路故障模糊集中心点

    Table  3  Fuzzy focus point of fault in turnout control circuit

    $B1$ $ B2 $ $ B3$ $ B4$ $ B5$ $ B6$ $ B7$ $B8$
    000000-0.750
    14.819.5322110232.42.4
    25
    108
    下载: 导出CSV
  • [1] Supavatanakul P, Lunze J, Puig V, Quevedo J. Diagnosis of timed automata:theory and application to the DAMADICS actuator benchmark problem. Control Engineering Practice, 2006, 14(6):609-619 doi: 10.1016/j.conengprac.2005.03.028
    [2] Previdi F, Parisini T. Model-free actuator fault detection using a spectral estimation approach:the case of the DAMADICS benchmark problem. Control Engineering Practice, 2006, 14(6):635-644 doi: 10.1016/j.conengprac.2005.04.001
    [3] Pedregal D J, García F P, Schmid F. RCM, predictive maintenance of railway systems based on unobserved components models. Reliability Engineering and System Safety, 2004, 83(1):103-110 doi: 10.1016/j.ress.2003.09.020
    [4] Zattoni E. Detection of incipient failures by using an H2-norm criterion:application to railway switching points. Control Engineering Practice, 2006, 14(8):885-895 doi: 10.1016/j.conengprac.2005.05.004
    [5] Puig V, Stancu A, Escobet T, Nejjari F, Quevedo J, Patton R J. Passive robust fault detection using interval observers:application to the DAMADICS benchmark problem. Control Engineering Practice, 2006, 14(6):621-633 doi: 10.1016/j.conengprac.2005.03.016
    [6] García Márquez F P, Schmid F, Collado J C. Wear assessment employing remote condition monitoring:a case study. Wear, 2003, 255(7-12):1209-1220 doi: 10.1016/S0043-1648(03)00214-X
    [7] Calado J M F, Sáda Costa J M G, Bartys M, Korbicz J. FDI approach to the DAMADICS benchmark problem based on qualitative reasoning coupled with fuzzy neural networks. Control Engineering Practice, 2006, 14(6):685-698 doi: 10.1016/j.conengprac.2005.03.025
    [8] 杨阳, 陶彩霞, 张睿兴.遗传算法优化支持向量机的道岔控制电路故障诊断.计算机测量与控制, 2013, 21(1):48-50 http://cdmd.cnki.com.cn/Article/CDMD-10287-1011292003.htm

    Yang Yang, Tao Cai-Xia, Zhang Rui-Xing. Fault diagnosis of switch control circuit using support vector machine optimized by genetic algorithm. Computer Measurement and Control, 2013, 21(1):48-50 http://cdmd.cnki.com.cn/Article/CDMD-10287-1011292003.htm
    [9] García Márquez F P, Pedregal Tercero D J, Schmid F. Unobserved component models applied to the assessment of wear in railway points:a case study. European Journal of Operational Research, 2007, 176(3):1703-1712 doi: 10.1016/j.ejor.2005.10.037
    [10] 王思明, 雷烨.一种基于LS-SVM的道岔控制电路故障诊断.兰州交通大学学报, 2010, 29(4):1-5 http://mall.cnki.net/magazine/Article/LZTX201004002.htm

    Wang Si-Ming, Lei Ye. Fault diagnosis for railway switch control circuit based on ARPSO least squares support vector machine. Journal of Lanzhou Jiaotong University, 2010, 29(4):1-5 http://mall.cnki.net/magazine/Article/LZTX201004002.htm
    [11] Huang Z J, Wang Z S, Zhang H G. Multilevel feature moving average ratio method for fault diagnosis of the microgrid inverter switch. IEEE/CAA Journal of Automatica Sinica, 2017, 4(2):177-185 doi: 10.1109/JAS.2017.7510496
    [12] 董炜, 陈卫征, 徐晓滨, 吉吟东.基于可分性测度的模糊隶属函数确定方法.控制与决策, 2014, 29(11):2089-2091 http://www.cnki.com.cn/Article/CJFDTotal-KZYC201411030.htm

    Dong Wei, Chen Wei-Zheng, Xu Xiao-Bin, Ji Yin-Dong. Determination method of fuzzy membership function based on separability measure. Control and Decision, 2014, 29(11):2089-2093 http://www.cnki.com.cn/Article/CJFDTotal-KZYC201411030.htm
    [13] 骆志明, 冯庚斌.机车车辆滚动轴承故障BP网络诊断方法.中国铁道科学, 1998, 19(4):26-32 http://mall.cnki.net/magazine/Article/TDYY200704015.htm

    Luo Zhi-Ming, Feng Geng-Bin. BP network fault diagnosis of ball bearing in locomotive and cars. China Railway Science, 1998, 19(4):26-32 http://mall.cnki.net/magazine/Article/TDYY200704015.htm
    [14] 徐晓滨, 张镇, 李世宝, 文成林.基于诊断证据静态融合与动态更新的故障诊断方法.自动化学报, 2016, 42(1):107-121 http://www.aas.net.cn/CN/abstract/abstract18800.shtml

    Xu Xiao-Bin, Zhang Zhen, Li Shi-Bao, Wen Cheng-Lin. Fault diagnosis based on fusion and updating of diagnosis evidence. Acta Automatica Sinica, 2016, 42(1):107-121 http://www.aas.net.cn/CN/abstract/abstract18800.shtml
    [15] 郭欣, 王蕾, 宣伯凯, 李彩萍.基于有监督Kohonen神经网络的步态识别.自动化学报, 2017, 43(3):430-438 http://www.aas.net.cn/CN/abstract/abstract19021.shtml

    Guo Xin, Wang Lei, Xuan Bo-Kai, Li Cai-Ping. Gait recognition based on supervised Kohonen neural network. Acta Automatica Sinica, 2017, 43(3):430-438 http://www.aas.net.cn/CN/abstract/abstract19021.shtml
    [16] Zhang P Y, Shu S, Zhou M C. An online fault detection model and strategies based on SVM-grid in clouds. IEEE/CAA Journal of Automatica Sinica, 2018, 5(2):445-456 doi: 10.1109/JAS.2017.7510817
    [17] 万俊, 邢焕革, 张晓晖.基于熵理论的多属性群决策专家权重的调整算法.控制与决策, 2010, 25(6):907-910 http://www.oalib.com/paper/4221554

    Wan Jun, Xing Huan-Ge, Zhang Xiao-Hui. Algorithm of adjusting weights of decision-makers in multi-attribute group decision-making based on entropy theory. Control and Decision, 2010, 25(6):907-910 http://www.oalib.com/paper/4221554
    [18] 董春玲, 张勤.用于不确定性故障诊断的权重逻辑推理算法研究.自动化学报, 2014, 40(12):2766-2781 http://www.aas.net.cn/CN/abstract/abstract18556.shtml

    Dong Chun-Ling, Zhang Qin. Research on weighted logical inference for uncertain fault diagnosis. Acta Automatica Sinica, 2014, 40(12):2766-2781 http://www.aas.net.cn/CN/abstract/abstract18556.shtml
  • 期刊类型引用(21)

    1. 李伟,邓志翔. 列车运行控制系统运营中的安全分析方法. 武汉冶金管理干部学院学报. 2023(02): 17-20 . 百度学术
    2. 张友鹏,魏智健,杨妮,张迪. 基于KPCA-SVM的S700K转辙机故障诊断方法. 安全与环境学报. 2023(09): 3089-3097 . 百度学术
    3. 欧阳鑫锋,孔令刚. 基于改进动态时间规整的道岔故障诊断方法. 现代信息科技. 2023(20): 136-139+143 . 百度学术
    4. Yong Chen,Christian Buerger,Miao Lin,Xudong Li,Volker Labenski,Haixia Jin,Hai Wang,Yang Liu,Tsuyoshi Ino,Harald Feifel,Tian Tan,Fangrong Chang. Left-turn-across-path-from-opposite-direction accidents in China:CIDAS accident study. Transportation Safety and Environment. 2023(04): 358-370 . 必应学术
    5. Yunting Zheng,Shaohua Chen,Zhiyong Tan,Yongkui Sun. Research on fault diagnosis of a railway point machine based on a multi-entropy feature extraction method and support vector machine. Transportation Safety and Environment. 2023(04): 338-346 . 必应学术
    6. 陈蕊. 基于电流曲线的道岔卡阻识别算法及实现. 自动化与仪表. 2022(04): 21-27 . 百度学术
    7. 池毅,陈光武. 基于一维卷积神经网络的实时道岔故障诊断. 计算机工程与应用. 2022(20): 293-299 . 百度学术
    8. 吴小雪,丁大伟,任莹莹,刘贺平. 二维FM系统的同时故障检测与控制. 自动化学报. 2021(01): 224-234 . 本站查看
    9. 李婉婉,李国宁. 基于GMM聚类和PNN的道岔故障诊断研究. 控制工程. 2021(03): 429-434 . 百度学术
    10. 郑云水,白邓宇,王妍. 基于相似度的道岔健康状态评估及故障检测方法研究. 铁道科学与工程学报. 2021(04): 877-884 . 百度学术
    11. 刘美容,刘津涛,何怡刚. 基于EMD复合多尺度熵的模拟电路故障诊断方法. 电子测量技术. 2021(04): 51-56 . 百度学术
    12. 李林,于颖. 智能继电保护回路故障监测全数字仿真研究. 计算机仿真. 2021(12): 460-464 . 百度学术
    13. 阮莹,梁利娟. 数字集成电路老化故障高精度预测方法仿真. 计算机仿真. 2020(02): 434-437 . 百度学术
    14. 高亚丽,陈光武. 基于改进FNN的道岔电路故障诊断方法. 科技创新与应用. 2020(15): 125-127 . 百度学术
    15. 孔令刚,焦相萌,陈光武,范多旺. 基于Mallat小波分解与改进GWO-SVM的道岔故障诊断. 铁道科学与工程学报. 2020(05): 1070-1079 . 百度学术
    16. 孔令刚,焦相萌,陈光武,范多旺. 基于多域特征提取与改进PSO-PNN的道岔故障诊断. 铁道科学与工程学报. 2020(06): 1327-1336 . 百度学术
    17. 杨菊花,于苡健,陈光武,司涌波,邢东峰. 基于CNN-GRU模型的道岔故障诊断算法研究. 铁道学报. 2020(07): 102-109 . 百度学术
    18. 姬文江,左元,黑新宏,高橋聖,中村英夫. 基于FastDTW的道岔故障智能诊断方法. 模式识别与人工智能. 2020(11): 1013-1022 . 百度学术
    19. Huidong Wang,Shifan He,Chengdong Li,Xiaohong Pan. Pythagorean Uncertain Linguistic Variable Hamy Mean Operator and Its Application to Multi-attribute Group Decision Making. IEEE/CAA Journal of Automatica Sinica. 2019(02): 527-539 . 必应学术
    20. 张友鹏,江雪莹,赵斌. 融合粗糙集与灰色模型的道岔故障预测. 铁道科学与工程学报. 2019(09): 2331-2338 . 百度学术
    21. 吴永成,阳长琼,何涛. 基于Fretchet距离与TWSVM的多机牵引道岔故障诊断研究. 铁道科学与工程学报. 2019(11): 2866-2872 . 百度学术

    其他类型引用(27)

  • 加载中
  • 图(9) / 表(3)
    计量
    • 文章访问数:  3143
    • HTML全文浏览量:  242
    • PDF下载量:  807
    • 被引次数: 48
    出版历程
    • 收稿日期:  2016-10-12
    • 录用日期:  2017-05-04
    • 刊出日期:  2018-06-20

    目录

    /

    返回文章
    返回