自动化学报  2017, Vol. 43 Issue (1): 101-113   PDF    
一类非线性离散时间动态系统的交替辨识算法及应用
张亚军1, 柴天佑1,2, 杨杰1     
1. 东北大学 流程工业综合自动化国家重点实验室 沈阳 110004;
2. 东北大学自动化研究中心 沈阳 110004
摘要: 表征产品在工业过程加工的质量、效率、成本、能耗或物耗等的运行指标与过程控制系统的输出密切相关,它们之间的动态模型往往机理不清,具有强非线性,难以用精确数学模型描述,但运行指标的预报对运行操作具有重要意义.本文利用工业过程在工作点附近工作的特点,将过程控制系统的输出与运行指标之间的动态模型描述成线性模型与高阶非线性项即未建模动态组成,对线性模型以及未建模动态提出了一种由改进的投影算法与未建模动态估计算法组成的交替辨识算法.最后,通过数值仿真实验和电熔镁炉的真实数据进行功率预报实验,实验结果表明了所提方法的有效性.
关键词: 复杂工业过程     监控系统     运行指标     交替辨识     神经网络    
Alternating Identification Algorithm and Its Application to a Class of Nonlinear Discrete-time Dynamical Systems
ZHANG Ya-Jun1, CHAI Tian-You1,2, YANG Jie1     
1. State Key Laboratory of Synthetical Automation for Process Industries, Northeastern University, Shenyang 110004;
2. Research Center of Automation, Northeastern University, Shenyang 110004
Received: 2015-12-10, Accepted: 2016-03-07.
Foundation Item: Supported by National Natural Science Foundation of China (61403071, 61603168), China Postdoctoral Science Foundation (2014M561246), Fundamental Research Funds for the Seed Foundation (N140804001), Natural Science Foundation of Liaon-ing Province (2015020144)
Author brief: ZHANG Ya-Jun Postdoctoral at Northeastern University. His research interest covers nonlinear fuzzy adaptive control theory, gen-eralized predictive control, multiple models and switchingsystems, intelligent decoupling control, data-based driven control, big data-driven modelling theory, method and technology of intelligent control systems, process industries and their applications.E-mail:zhangyajun79@gmail.com;
YANG Jie Ph. D. candidate at the State Key Laboratory of Synthet-ical Automation for Process Industries.His research interest covers data driven modeling technology and application for industrial pro-cess.E-mail:yjercou@126.com
Corresponding author. CHAI Tian-You Academician of Chinese Academy of Engineering, professor at Northeast-ern University, IEEE Fellow, IFAC Fellow, IEAS Fellow.He received his Ph. D. degree from Northeastern Univer-sity in 1985. His research interest covers adaptive control,intelligent decoupling control, and integrated automation theory, method and technology of industrial process. Cor-responding author of this paper.). E-mail:tychai@mail.neu.edu.cn
Abstract: The major operational indexes such as quality, efficiency, cost, energy and material consumptions in industrial process of product processing are closely related to the output of the process control system; their dynamic models between the operational indexes and the output of the process control system are often nonlinear and with unclear structure nature generally. Therefore, it is difficult to obtain an accurate model. However, the prediction of operational index is of great significance to the operational operations. In this paper, based upon the characteristic that complex industrial systems often work near an operating point, the dynamic model between the operational indexes and the output of the process control system is represented by a linear model plus a higher order nonlinear term (unmodeled dynamics). With the above development, an alternating identification algorithm which consists of improved projection algorithm for the linear model and an estimation algorithm for unmodeled dynamic are proposed. Finally, through simulation study and a power forecast experiment using real data of an electric-melting magnesia furnace, the effectiveness of the proposed algorithm is justified.
Key words: Complex industrial process     monitoring system     operational indexes     alternating identification     neural networks    

工业过程的运行指标往往表示产品在某加工过程中的质量、效率、成本或者能耗、物耗等指标[1-2].过程控制不仅要保证被控对象的输出尽可能好地跟踪设定值[3],而且要对运行指标进行监控,保证运行指标在目标值范围内,尽可能地提高质量与效率指标,尽可能地降低成本或能耗、物耗等指标,实现性能指标最优[1-2, 4].例如电熔镁炉的运行指标是用电量,即功率,功率高,说明能耗高,如果功率超过工业的限幅值,就得拉闸断电以保证生产的正常运行,但过多的拉闸断电会影响电熔镁砂的产品品位.由于运行指标监控系统由过程控制系统的输出和运行指标之间的动态模型组成,其中过程控制系统一般由控制器、传感器、执行机构以及被控对象组成.因此,过程控制系统的输出影响运行指标的变化.过去对运行指标的监控主要关注的是所监控的量是否在某个时刻到达工业的限幅值,然后再采取措施,并没有考虑过程控制系统对运行指标下一时刻变化趋势的动态影响.复杂工业过程,如电熔镁砂,虽然功率在某个时刻达到工业限幅值,但很有可能在下一时刻,由于电流闭环控制系统的作用使得功率下降,如果不能及时预测功率的这种变化趋势,则往往造成误操作,从而导致拉闸次数过多. 因此,对运行指标的监控和预测十分重要,并成为复杂工业过程领域的研究热点,引起了工业界和学术界的广泛关注[5-17]. 如文献[7]采用电信网络管理标准设计了电力系统的一种监控方案;文献[7]提出了一种采用混合神经网络与粒子群优化算法组成的风能预测方法;文献[8]通过对某大型火力发电厂的冷却过程设计监控系统,通过采集到的数据提出了一种监控方法;文献[8]采用两个神经网络直接对风力发电场的发电量区间进行短期预测.文献[10-12]针对某地区的电力负荷数据,提出了电力负荷的相似性时间序列预报算法以及基于神经网络和支持向量机的电力负荷预报方法;文献[15]针对电力系统设计了具有监控和分析功能的多目标平台; 文献[14]提出了能量监控系统,并用来扩大磁性能量收集电路中无线传感系统可持续性.文献[13]研究了针对小规模太阳光伏发电系统超前一天的功率预报方法.上述方法存在预报周期较长,并且没有考虑由于过程控制系统的动态特性使得运行指标非线性动态变化的问题,本质上属于一种静态预报方法.如果考虑过程控制系统输出对运行指标的影响,那么整个监控系统是非线性动态变化的. 因此,上述监控方法难以直接应用到类似于电熔镁炉功率监控系统的复杂工业过程.电熔镁炉功率的动态模型与电流和电阻有关,电流是闭环控制系统的输出,是动态变化的. 电阻包括电弧电阻和熔池电阻,并且受原料熔点、电阻率、弧长、电弧半径、电弧温度以及熔池深度等因素的影响,具有强非线性特性,难以建立精确数学模型,从而导致功率的强非线性变化.考虑到电流闭环控制系统是稳定的,输出电流有界,监控的运行指标即功率也是有界的.因此,根据以上特性,可将电流与功率之间的动态模型在工作点附近表示成由线性模型与未建模动态和的形式,从而对功率指标进行预报,根据预报值提前采取措施使得运行指标达到最优.

综上所述,本文把过程控制系统闭环系统输出对运行指标的影响综合考虑,充分利用了过程控制系统的有界输出以及运行指标本身的有界性,在工作点附近将这类运行指标模型用线性模型与未建模动态来描述,建立了运行指标的动态预报模型,并对线性模型的参数以及未建模动态给出了一种交替辨识算法.为了突出所提辨识算法的原理简单,对未知的模型参数采用改进的投影算法进行估计,对未建模动态采用最早最常用的BP神经网络(Back-propagation neural networks)[16]进行估计.采用上述两种估计算法交替辨识的方式对监控系统的运行指标进行预报,总结了所提算法的具体实现步骤,并且通过数值仿真实验以及电熔镁炉真实的电流和功率数据进行功率预报实验,结果表明了所提算法的有效性. 本文的主要贡献和创新工作总结如下:

1) 针对运行指标与过程控制系统的输出之间的动态模型机理不清,具有强非线性,难以用精确数学模型描述的问题,在工作点附近建立了运行指标的动态预报模型.

2) 分别在运行指标动态模型线性模型参数已知,未建模动态未知以及线性模型参数和未建模动态均未知的情况,提出了基于神经网络的未建模动态估计算法以及一种由改进的投影算法与未建模动态估计算法组成的交替辨识算法.

3) 通过数值仿真实验和电熔镁炉的真实数据进行功率预报实验,实验结果表明了所提方法的有效性.

本文的组织结构如下:第1节,建立了过程控制系统的输出和监控系统运行指标之间的动态预报模型;第2节针对运行指标预报模型,提出了由改进的投影算法和BP神经网络估计算法组成的交替辨识算法;第3节给出了数值仿真实验以及在电熔镁炉功率预报的应用实验;第4节总结了本文的主要结果.

1 监控系统运行指标的动态模型

复杂工业过程的运行指标由过程控制闭环系统的输出和运行指标之间的动态模型组成,但过程系统的输出与运行指标模型之间往往机理不清,具有强非线性特性,难以建立精确数学模型. 因此,考虑过程控制系统的输出对运行指标影响,它们之间动态模型可以描述为

$\begin{align} &y(k)=f[y(k-1),\cdots ,y(k-{{n}_{A}}),u(k-d),\cdots , \\ &u(k-d-{{n}_{B}})] \\ \end{align}$ (1)

式中,d为延时,u(k)为运行指标模型的输入,它代表过程控制系统的闭环输出,且$0\leq\left|{u(k)}\right|<\infty$.y(k)表示运行指标,且$0\leq\left|{y(k)} \right|<\infty$.$n_{A}$$n_{B}$ 为运行指标模型的阶次; $f(\cdot)\in{\bf R}$是未知的非线性函数.

$\begin{align} &\varphi (k-d)=[y(k-\text{1}),\cdots ,y(k-{{n}_{A}}), \\ &u(k-d),\cdots ,u(k-d-{{n}_{B}}){{]}^{\text{T}}} \\ \end{align}$ (2)

${\bf\pmb{{\varphi }}}{\rm{(}}k - d{\rm{)}}$是维数为$p={n_A} + {n_B}+ 1$ 的数据向量,则模型(1) 可表示为

$y(k)=f[\varphi (k-d)]$ (3)

由于过程控制系统是稳定的,输出是有界的; 运行指标本身也是有界的,因此,运行指标的动态模型具有输入输出有界性的特点.在一般情况下,由于运行指标在工作点附近工作,从而可将运行指标模型在工作点附近表示成由线性模型与高阶非线性项和的形式.为此,在工作点附近将运行指标模型(3) 线性化就可得到由线性模型与高阶非线性项组成的数学模型.假设工作点为原点(0,0),也就是运行指标模型(1) 的平衡点 (若工作点偏离原点,可用坐标变换将其移至原点),对应于式(1) ,该平衡点可表示为y = 0,u = 0. 将$f[{\bf\pmb{{\varphi }}}{\rm{(}}k - d{\rm{)}}]$ 在附近Taylor 展开,并令其在该点处的一阶Taylor系数分别为y=0,u=0.将y=0,u=0在附近Taylor 展开,并令其在该点处的一阶Taylor系数分别为

${{a}_{i}}=-\frac{\partial f[\varphi (k-d)]}{\partial y(k-i)}{{|}_{y\text{=0},u\text{=0}}},i=1,\cdots ,{{n}_{A}}$ (4)
${{b}_{i}}=-\frac{\partial f[\mathbf{\varphi }(k-d)]}{\partial u(k-d-j)}{{|}_{y\text{=0},u\text{=0}}},j=1,\cdots ,{{n}_{B}}$ (5)

由式(4) 和式(5) 组成${z^{ - 1}}$的多项式$A({z^{ - 1}})$$B({z^{ -1}})$的系数,$A({z^{ - 1}})$$B({z^{ - 1}})$分别为

$A({{z}^{-1}})=1+{{a}_{1}}{{z}^{-1}}+\cdots +{{a}_{{{n}_{A}}}}{{z}^{-{{n}_{A}}}}$ (6)
$B({{z}^{-1}})={{b}_{0}}+{{b}_{1}}{{z}^{-1}}+\cdots +{{b}_{{{n}_{B}}}}{{z}^{-{{n}_{B}}}}$ (7)

于是,运行指标的动态模型(1) 在平衡点(0,0)附近的等价模型如下:

$A({{z}^{-1}})y(k)=B({{z}^{-1}})u(k-d)+v[\mathbf{\varphi }(k-d)]$ (8)

式中,$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k-d{\rm{)]}}$是高阶非线性函数[18],称之为未建模动态.

为了突出方法原理和推导简单起见,设d = 1,于是运行指标的预报模型(8) 可写为

$A({{z}^{-1}})y(k+1)=B({{z}^{-1}})u(k)+v[\varphi (k)]$ (9)

针对模型(9) ,由于$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$未知,要得到运行指标的预报值须首先对k时刻的未建模动态$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$进行估计.因此,如果模型(9)中的$A({z^{ - 1}})$$B({z^{ - 1}})$ 已知,那么,$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$ 的准确估计是影响运行指标的关键.如果模型(9) 中的$A({z^{ - 1}})$$B({z^{ - 1}})$ 未知,则不但要估计未建模动态$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$,而且同时需要估计多项式$A({z^{ - 1}})$$B({z^{ - 1}})$ 的参数,此时,两种估计算法并存,既互相影响,又互为补充. 为了使所提算法具有一般性,下面采用常用的改进投影算法与基于BP神经网络的未建模动态估计算法组成的交替辨识方法对运行指标模型进行辨识.

注 1. 为简单明了,本文只考虑了时滞d=1的简单情况.当复杂非线性系统的时滞d>1时,可将本文交替辨识的思想做相应的推广. 然而,不得不强调的是,当被控对象为具有大时滞的复杂工业过程,如何有效地对系统进行预报仍然是一个相当具有挑战性的研究课题.

2 运行指标模型的辨识算法

首先考虑当运行指标模型(9) 中的$A({z^{ - 1}})$$B({z^{ -1}})$已知时,$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$的估计算法. 为了突出算法的原理,选用最常用的BP神经网络对$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$进行估计. 其次,考虑$A({z^{ - 1}})$$B({z^{ -1}})$$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$均未知时的非线性模型辨识算法.

由式(9) 和式(6) 可知,

$\begin{array}{l} v\left[ {\varphi \left( k \right)} \right] = A\left( {{z^{ - {\rm{1}}}}} \right)y\left( {k + 1} \right) - B\left( {{z^{ - {\rm{1}}}}} \right)u\left( k \right) = \\ y\left( {k + 1} \right) + \left[ {A\left( {{z^{ - {\rm{1}}}}} \right) - 1} \right]y\left( {k + 1} \right) - B\left( {{z^{ - {\rm{1}}}}} \right)u\left( k \right) = \\ y\left( {k + 1} \right) + z\left[ {A\left( {{z^{ - {\rm{1}}}}} \right) - 1} \right]y\left( k \right) - B\left( {{z^{ - {\rm{1}}}}} \right)u\left( k \right) = \\ y\left( {k + 1} \right) + \bar A\left( {{z^{ - {\rm{1}}}}} \right)y\left( k \right) - B\left( {{z^{ - {\rm{1}}}}} \right)u\left( k \right) \end{array}$ (10)

式中

$\bar{A}({{z}^{-1}})=z[A({{z}^{-1}})-\text{1}]$ (11)
$\begin{align} &\varphi (k)=[y(k-1),\cdots ,y(k-{{n}_{A}}), \\ &u(k),\cdots ,u(k-{{n}_{B}}){{]}^{\text{T}}} \\ \end{align}$

虽然$\bar A({z^{ - 1}})$$B({z^{ - 1}})$已知,但$y(k{\rm{ +}}1) $未知,因此,$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$未知. 由于$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$ 为高阶非线性函数,可以采用BP神经网络对$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$进行估计,得到其估计值$\hat v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$. 采用$u{\rm{(}}k{\rm{)}}$作用于运行指标模型(9) 可得预报值$y(k{\rm{ + }}1) $. 由式(10) 可知:

$y(k\text{+}1)=-\bar{A}({{z}^{-1}})y(k)+B({{z}^{-1}})u(k)+\hat{v}[\varphi (k)]$ (12)

因此,神经网络的导师信号为

$\begin{align} &v[\varphi (k)]\text{=}y(k\text{+}1)+\bar{A}({{z}^{-1}})y(k)-B({{z}^{-1}})u(k)= \\ &A({{z}^{-1}})y(k\text{+}1)-B({{z}^{-1}})u(k) \\ \end{align}$ (13)

由数据向量${\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)}}$以及导师信号$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$ 对神经网络进行训练产生下一时刻的估计值$\hat v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$.

注 2. 值得注意的是,式(10) 和(13) 是有区别的,式(10) 中的$y(k{\rm{ + }}1) $ 未知,但式(13) 中的$y(k{\rm{ + }}1) $是已知的,用于获取神经网络的导师信号.

下面介绍采用BP神经网络对$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$ 的估计算法.

BP神经网络是1986年由D. E. Rumelhart和J. L.McClelland提出的一种利用误差反向传播训练算法的神经网络,简称BP网络,它是一种无反馈的前向网络,网络中的神经元分层排列,由输入层、隐含层和输出层组成; 隐含层一般选一层,但隐含层的节点数可采用凑试的方法来选择,也可以采用正交最小二乘法来逐个增加隐藏层的节点数,然后权衡预报误差的指标和神经网络的计算复杂度,使得神经网络的预报误差符合辨识的精度要求,具体详见文献[19-20].本文对$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$的估计采用三层BP神经网络. 其中,输入层的数据向量采用式(2) 所示的,即:

$\begin{align} &\varphi (k)={{[{{\varphi }_{1}},\cdots ,{{\varphi }_{i}},\cdots ,{{\varphi }_{p}}]}^{\text{T}}}= \\ &[y(k-1),\cdots ,y(k-{{n}_{A}}), \\ &u(k),\cdots ,u(k-{{n}_{B}}){{]}^{\text{T}}},i=1,\text{2},\cdots ,p \\ \end{align}$ (14)

式中,$p=n_A + n_B + 1$为输入层神经元的个数.隐含层神经元节点数为q. 神经网络的输出为$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$ 的估计值$\hat v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$. BP神经网络的结构如图 1 所示.

图 1 BP神经网络估计$v[\varphi (k)]$的结构 Figure 1 The structure of the estimation for $v[\varphi (\text{k})]$ by BP neural networks[

图 1中,${W}{\rm{ = [}}{{W}^1},{{W}^2}{\rm{]}}$代表网络的权向量,${{W}^1}$为输入层到隐含层的链接权矩阵,即

${{W}^{1}}=\left[ \begin{matrix} w_{11}^{1}&\cdots &w_{1q}^{1} \\ \vdots &\ddots &\vdots \\ w_{p1}^{1}&\cdots &w_{pq}^{1} \\ \end{matrix} \right]$ (15)

式中,$w_{ij}^1$(i = 1,2,\cdots ,p; j = 1,2,\cdots {\rm{,}} q)$ 表示输入层神经元i与隐含层神经元j之间的连接权;${{W}^2}$ 为隐含层到输出层的链接权矩阵,即

${{W}^{2}}=\left[ \begin{array}{*{35}{l}} w_{11}^{2} \\ \vdots \\ w_{q1}^{2} \\ \end{array} \right]$ (16)

式中,$w_{j{\rm{1}}}^2$表示隐含层神经元j与输出层神经元之间的连接权.

因此,BP神经网络对$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$end{document}的估计可简单地表示为

$\hat{v}[\varphi (k)]=NN[W(k),\varphi (k)]$ (17)

式中,${\pmb N}{\pmb N}[{W},{\bf\pmb{{\varphi}}}]$代表结构如图 1所示的${W}(k)$${\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)}}$ 的神经网络的函数,${\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)}}$为神经网络的输入数据向量,${W}(k){\rm{ =[}}{{W}^1}(k),{{W}^2}(k){\rm{]}}$k时刻对权${W}$ 的估计.在图 1所示的估计的结构中,隐含层神经元的输入${\pmb N}$

$N={{[{{W}^{1}}(k)]}^{\text{T}}}\varphi (k)$ (18)

式中,${{W}^1}{\rm{(}}k{\rm{)}}$${{W}^1}$的估计值; ${\pmb{N}} =[{n_{\rm{1}}},\cdots ,{n_j},\cdots ,{n_q}{{\rm{]}}^{\rm{T}}}$,${n_j}$为隐含层第j个神经元的输入,表示如下:

${{n}_{j}}=\sum\limits_{i\text{=}1}^{p}{w_{ij}^{1}(k){{\varphi }_{i}}(k)},j=1,2,\cdots ,q$ (19)

${w_{ij}}{\rm{(}}k{\rm{)}}$的调整量按照下面的式(27) 和(28) 进行校正.

隐含层神经元的输出为

$M=g(N)$ (20)

式中,${\pmb{M}} = {[{m_1},\cdots ,{m_j},\cdots ,{m_q}]^{\rm{T}}}$; ${m_j}$ 为隐含层第$j(j=1,2,\cdots ,q$个神经元的输出. ${\pmb{g}} = {[{g_1},\cdots ,{g_j},\cdots ,{g_q}]^{\rm{T}}}$,${g_j}$为sigmoid 函数且

${{m}_{j}}={{g}_{j}}({{n}_{j}})=\frac{1}{1+\exp [-({{n}_{j}}+\sigma )]}$ (21)

式中,$\sigma$为阈值或偏置值.$\sigma>0$表示使式(21) 描述的函数曲线沿横坐标左移; 反之,则右移.隐含层的第$j$(j = 1,2,\cdots ,q)$个神经元的输出${m_j}$将通过权系数向前传播到输出层神经元,并作为它的将通过权系数向前传播到输出层神经元并作为它的输入之一,而输出层的总输入$\bar n$

$\begin{align} &\bar{n}={{({{W}^{2}}(k))}^{\text{T}}}M= \\ &\sum\limits_{j\text{=}1}^{q}{w_{j1}^{2}(k){{m}_{j}}},j=1,2,\cdots ,q \\ \end{align}$ (22)

式中,${{W}^2}{\rm{(}}k{\rm{)}}$${{W}^2}$ 的估计值,$w_{j1}^2(k)$的调整量按照式(24) 和(25) 进行校正.

输出层的输出即为未建模动态$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$ 的估计值,即

$\begin{align} &\hat{v}[\varphi (k)]=g(\bar{n})=g\left[ \sum\limits_{j\text{=}1}^{q}{w_{j1}^{2}(k){{m}_{j}}} \right]= \\ &\frac{1}{1+\exp \left[ -\sum\limits_{j\text{=}1}^{q}{w_{j1}^{2}(k){{m}_{j}}}-\sigma \right]} \\ \end{align}$ (23)

式中,g为前面介绍的sigmoid函数.输出层的权$w_{j1}^2{\rm{(}}k{\rm{)}}$的调权律为

$w_{j1}^{2}(k\text{+}1)=w_{j1}^{2}(k)\text{+}\Delta w_{j1}^{2}(k)$ (24)
$\begin{align} &\Delta w_{j1}^{2}(k)=\eta \hat{v}[\varphi (k)]\left\{ \text{1}-\hat{v}[\varphi (k)] \right\}e(k){{m}_{j}}, \\ &j=1,2,\cdots ,q \\ \end{align}$ (25)

式中,权的初值$w_{j1}^2(0) $选为较小的随机数,$\eta$为学习速率,$\eta>0$,${m_j}$ 为隐含层第${j}$个神经元的输出,$e{\rm{(}}k{\rm{)}}$ 为神经网络的估计误差:

$e(k)=v[\varphi (k)]-\hat{v}[\varphi (k)]$ (26)

$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$为神经网络的导师信号,按照式(13) 计算.

隐含层的权$w_{ij}^1{\rm{(}}k{\rm{)}}$的调权律为

$w_{ij}^{1}(k\text{+}1)=w_{ij}^{1}(k)\text{+}\Delta w_{ij}^{1}(k)$ (27)
$\Delta w_{ij}^{1}(k)=\eta {{m}_{j}}(\text{1}-{{m}_{j}})\left( -\frac{\partial E}{\partial {{m}_{j}}} \right)\hat{v}[\varphi (k)]$ (28)

式中,权的初值$w_{ij}^1{\rm{(}}0{\rm{)}}$ 选为较小的随机数,$\eta$为学习速率,$\eta>0$,${m_j}$ 为隐含层第j个神经元的输出,E为总误差函数:

$E=\frac{1}{\text{2}}\sum\limits_{\tau \text{=1}}^{L}{e{{(\tau )}^{\text{2}}}}$ (29)

式中,L为样本总数,$e{\rm{(}}\tau {\rm{)}}$按照式(26) 计算.

基于BP神经网络的未建模动态估计算法步骤总结如下:

1) 设$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$估计值的初始值,用$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k-1{\rm{)]}}$ 作为$\hat v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$ 的初始值; 对BP 神经网络权${W}{\rm{ =[}}{{W}^1},{{W}^2}{\rm{]}}$赋初值,初值一般选较小的随机数;

2) 组成数据向量${\bf\pmb{{\varphi }}}(k)$,由式(13) 计算神经网络的导师信号$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$;

3) 根据式(18) ~(23) 分别计算隐含层、输出层各神经元的输出,得到$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$的估计值;

4) 采用式(26) 计算估计误差;

5) 由式(27) ~(28) 求下一时刻神经网络输出层的权值;由式(24) ~(25) 求下一时刻神经网络隐含层的权值;

6) 返回2) ,估计下一时刻未建模动态的值.

2.1 运行指标线性模型参数未知时的交替辨识算法

由于许多复杂工业过程的运行指标难以建立精确的数学模型,因此$A({z^{- 1}})$$B({z^{ - 1}})$以及未建模动态$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$均未知,此时,须首先采用参数辨识算法辨识系统参数,对未建模动态采用BP神经网络进行估计,从而实现对整个运行指标模型的辨识. 由式(9) 和式(13) 可得:

$v[\varphi (k)]=A({{z}^{-1}})y(k+1)-B({{z}^{-1}})u(k)$ (30)
$y(k+1)-v[\varphi (k)]=-\bar{A}({{z}^{-1}})y(k)\text{+}B({{z}^{-1}})u(k)$ (31)

从式(30) 可以看出,由于$A({z^{ - 1}})$$B({z^{ - 1}})$未知,因此$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$ 未知.故在辨识时,首先对$A({z^{ - 1}})$$B({z^{ - 1}})$赋予初值,从而得到$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$ 的估计值,将得到的$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$的估计值作为神经网络的导师信号,采用神经网络对$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$进行估计,得到其估计值$\hat v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$,将其代入式(31) 并采用改进投影算法辨识$\bar A({z^{ - 1}})$$B({z^{ -1}})$,从而得到$A({z^{ - 1}})$$B({z^{ - 1}})$ 的辨识值$\hat A({z^{ - 1}})$$\hat B({z^{ - 1}})$. 然后将$\hat A({z^{ -1}})$$\hat B({z^{ - 1}})$代入式(30) 获得下一时刻的导师信号 ,即

$\bar{v}[\varphi (k)]=\hat{A}({{z}^{-1}})y(k\text{+}1)-\hat{B}({{z}^{-1}})u(k)$ (32)

再通过下一时刻的数据向量${\bf\pmb{{\varphi }}}(k)$,采用BP神经网络进行估计,获得下一时刻$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$的估计值,如此反复,即采用估计$\left\{{A({z^{ - 1}}),B({z^{ - 1}})} \right\}$的改进投影算法与$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$的神经网络估计算法的交替辨识的方式对非线性模型式(9) 进行辨识.

综上所述,$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$的神经网络估计的导师信号方程和估计$\left\{{A({z^{ - 1}}),B({z^{ - 1}})} \right\}$ 的参数辨识方程如下:

$\bar{v}[\varphi (k)]=\hat{A}({{z}^{-1}})y(k\text{+}1)-\hat{B}({{z}^{-1}})u(k)$ (33)
$y(k\text{+}1)-\hat{v}[\varphi (k)]=-\bar{A}({{z}^{-1}})y(k)\text{+}B({{z}^{-1}})u(k)$ (34)

为了辨识$A({z^{ - 1}})$$B({z^{ - 1}})$,将辨识方程(34) 写成如下形式:

$y(k)={{\varphi }^{\text{T}}}(k-1)\widehat{\theta }(k-1)+\hat{v}[\varphi (k-1)]$ (35)

式中,$\widehat{\theta }(k-1)={{[-{{a}_{\text{1}}}(k-1),-{{a}_{\text{2}}}(k-1),\cdots ,-{{a}_{{{n}_{A}}}}(k-1),{{b}_{\text{0}}}(k-1),{{b}_{1}}(k-1),\cdots ,{{b}_{{{n}_{B}}}}(k-1)]}^{\text{T}}}$$k-1$时刻对参数${\pmb\theta} = [- {a_{\rm{1}}}{\rm{,}} - {a_{\rm{2}}}{\rm{,}} \cdots {\rm{,}} - {a_{{n_A}}},{b_{\rm{0}}},{b_1},\cdots ,b_{{n_B}}{{\rm{]}}^{\rm T}}$ 的估计值,采用如下的改进投影算法进行估计.

$\widehat{\theta }(k)=\widehat{\theta }(k-1)+\frac{\varphi (k-1)e(k)}{1+{{\varphi }^{\text{T}}}(k-1)\varphi (k-1)}$ (36)

式中

$e(k)=y(k)-{{\varphi }^{\text{T}}}(k-1)\widehat{\theta }(k-1)-\hat{v}[\varphi (k-1)]$ (37)

最后,把$\hat A({z^{ - 1}})$$\hat B({z^{ - 1}})$的估计值代入式(33) ,获得$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$神经网络的估计算法的导师信号$\bar v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$.构建下一时刻的数据向量,采用BP 神经网络估计算法式(18) 和式(29) 对$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k-1{\rm{)]}}$进行估计,将获得的估计值$\hat v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k - 1{\rm{)]}}$代入$\left\{ {A({z^{ - 1}}),B({z^{ -1}})} \right\}$ 的参数辨识方程(35) ,采用辨识算法式(36) 和式(37) 对$\left\{ {A({z^{ - 1}}),B({z^{ - 1}})}\right\}$进行辨识.

上述非线性模型辨识算法步骤总结如下:

1) 对$A({z^{ - 1}})$$B({z^{ - 1}})$赋予初值,即给定初值$\hat\theta (0) $;

2) 根据选定$\left\{ {A({z^{ - 1}}),B({z^{ - 1}})} \right\}$ 的初值,通过式(33) 获得与$\left\{ {A({z^{ - 1}}),B({z^{ - 1}})}\right\}$的初值相对应的神经网络的导师信号$\bar v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$;

3) 测取$y(k),\cdots,u(k),\cdots$,构成数据向量${\bf\pmb{{\varphi }}}(k)$,并与导师信号$\bar v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$训练神经网络,获得$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$的估计值$\hat v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$;

4) 将估计值$\hat v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$代入$\left\{ {A({z^{ - 1}}),B({z^{ - 1}})} \right\}$的参数辨识方程(34) ,采用辨识算法式(36) 和式(37) 估计$\left\{ {A({z^{- 1}}),B({z^{ - 1}})} \right\}$;

5) 将$\left\{ {A({z^{ - 1}}),B({z^{ - 1}})} \right\}$的估计值$\left\{ {\hat A({z^{ - 1}}),\hat B({z^{ - 1}})}\right\}$代入式(33) 并获得下一时刻的神经网络估计的导师信号 ;

6) 返回3) ,估计下一时刻的$\left\{ {A({z^{ - 1}}),B({z^{ - 1}})}\right\}$$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}k{\rm{)]}}$.

3 仿真 3.1 数值仿真

为了说明本文所提估计算法的有效性,假设运行指标的非线性模型如下:

$\begin{align} &y(k+1)=0.6y(k)-0.2y(k-1)+ \\ &u(k)+\text{0}.\text{5}u(k-1)+ \\ &\sin [u(k)+u(k-1)+y(k)+y(k-1)]- \\ &\frac{u(k)+u(k-1)+y(k)+y(k-1)}{1+u{{(k)}^{2}}+u{{(k-1)}^{2}}+y{{(k)}^{2}}+y{{(k-1)}^{2}}} \\ \end{align}$

原点是系统的平衡点.

$\begin{align} &A({{z}^{-1}})=1-0.6{{z}^{-1}}+0.2{{z}^{-2}} \\ &B({{z}^{-1}})=1+0.5{{z}^{-1}} \\ &v[\varphi (k)]=\sin [u(k)+u(k-1)+y(k)+y(k-1)]- \\ &\frac{u(k)+u(k-1)+y(k)+y(k-1)}{1+u{{(k)}^{2}}+u{{(k-1)}^{2}}+y{{(k)}^{2}}+y{{(k-1)}^{2}}} \\ &\varphi (k)={{[y(k),y(k-1),u(k),u(k-1)]}^{\text{T}}} \\ \end{align}$

易知该模型的未建模动态全局有界.选择BP神经网络的隐含层节点数$l=18$,激活函数为“S-型”传输函数;输出层的激活函数取线性传输函数. 学习率$lr=1$,动量因子$mc=0.95$.

在仿真中,选择输入信号为

$u(k)=\sin (\frac{2\pi k}{15})+\sin (\frac{2\pi k}{25})$

采用本文所提的辨识算法时,初始状态为$y(0) = 0$,$u(0) = 0$,${\pmb\theta} (0) = [0.5,- 0.15,1.05,0.45]$,$v{\rm{[}}{\bf\pmb{{\varphi }}}{\rm{(}}0{\rm{)] = 0}}$.

首先,在系统线性模型参数已知的情况下,采用BP神经网络开环辨识被控对象的未建模动态,仿真结果如图 2~6所示.

图 2 输入信号u Figure 2 Control input signal u
图 3 未建模动态(星线)及其估计值(圈线) Figure 3 The un-modeled dynamics (star line) and its estimation value (circle line)
图 4 采用BP神经网络估计未建模动态的误差 Figure 4 Estimation error for the un-modeled dynamics which produced by BP neural networks
图 5 采用BP神经网络的估计误差(残差)的概率密度函数 Figure 5 The probability density functions of prediction errors (residuals) by BP neural networks
图 6 BP神经网络的性能、训练状况和相关系数 Figure 6 The performance,training state,and regression of BP neural networks

其次,在系统线性模型参数未知的情况下,采用改进投影算法以及基于~~BP~~神经网络估计算法组成的交替辨识算法进行非线性系统的辨识,仿真结果如图 7~11所示.

图 7 $A({z^{ - 1}})$ 中的参数和$B({z^{ - 1}})$中的参数 Figure 7 Parameters in A(z-1) and parameters in B(z-1)
图 8 未建模动态(星线)及其估计值(圈线) Figure 8 The un-modeled dynamics (star line) and its estimation value (circle line) by BP neural networks
图 9 采用BP神经网络估计未建模动态的误差 Figure 9 Estimation error for the un-modeled dynamics which produced by BP neural networks with he proposed method
图 10 采用BP神经网络的估计误差(残差)的概率密度函数 Figure 10 The probability density functions of prediction errors (residuals) by BP neural networks
图 11 BP 神经网络的性能、训练状况和相关系数 Figure 11 The performance, training state, and regression of BP neural networks

图 2是运行指标模型的输入信号,它是过程控制系统闭环控制系统的有界输出. 图 3 是采用BP神经网络估计未建模动态的效果图. 从图 3中可以看出,基于BP神经网络的未建模动态估计算法能保证大部分时刻估计方向的准确性,说明了估计算法的有效性. 图 4 是估计误差,从图 4中可以看出,估计误差较小,并且在零点附近分布密集,说明BP神经网络能较好地估计未建模动态.图 5是采用BP神经网络估计未建模动态的误差分布的概率密度函数直方图,从图 5中可以看出,相对误差或残差具有高斯型的概率密度函数,说明估计精度高,估计效果好.图 6是BP神经网络训练的性能(图 6(a))、训练状况(图 6(b))、相关系数(图 6(c)).从图 6中可以明显地看出,网络在200次迭代后就达到了最优的性能0.0020787. 从数据的回归曲线可以看出,估计值与真值之间的相关系数为0.99812,说明神经网络的估计值与未建模动态的真值相关性很大,保证了估计方向的正确性以及估计精度.

当非线性系统模型完全未知时,图 7是采用改进的投影算法估计被控对象线性模型参数的情况.从图 7中可以看出,由于受到未建模动态的影响,参数的估计值在真值附近变化,但并不收敛于真值,从而导致参数估计的残差也包含于未建模动态中. 图 8~10是BP神经网络估计系统未建模动态的情况. 由于系统未知,未建模动态实际上包含了参数估计的残差以及高阶非线性部分.从图 8中可以看出,此时,基于BP神经网络的未建模动态估计算法也能保证大部分时刻估计方向的准确性,说明了估计算法的有效性. 图 9是估计误差,从图 9 中可以看出,估计误差较小. 图 10是BP 神经网络估计误差分布的概率密度函数直方图.从图 10中可以看出,估计的相对误差或者残差也具有高斯型的概率密度函数,说明估计精度高,估计效果好. 比较图 10图 5可以看出,由于线性模型参数未知时的未建模动态比系统线性模型参数已知时更加复杂,更难以估计,因此,估计效果稍显逊色. 比较图 5图 10可知,图 5的概率密度函数直方图比图 10的更窄、更高,这说明图 5的估计效果优于图 10.这也体现了线性模型的参数估计算法与未建模动态估计算法确实是相互影响的,只有两个估计算法交替协调进行辨识,才能达到辨识整个非线性模型的目的.图 11是BP神经网络训练的性能(图 11(a))、训练状况(图 11(b))、相关系数(图 11(c)).从图 11中可以看出,网络在200 次迭代后就达到了最优的训练性能0.00436.从回归曲线可以看出,相关系数为0.99612,以上两个关键系数都比系统线性模型参数已知时的稍小,这正是由于两种估计算法互相影响所致.

3.2 在电熔镁炉功率监控系统中的应用

将本文所提的方法应用在电熔镁炉的功率监控系统进行功率预报实验.电熔镁炉的生产过程如图 12所示,其中,电熔镁砂的熔炼过程中发生如下的物理化学反应[21-22]:

图 12 电熔镁炉结构示意图 Figure 12 The structure of electric-melting magnesia furnace
$\text{MgC}{{\text{O}}_{3}}(\text{S})\xrightarrow{{{600}^{{}^\circ }}\text{C}~\sim {{800}^{{}^\circ }}\text{C}}\text{C}{{\text{O}}_{2}}\uparrow +\text{MgO}(\text{S})$ (38)
$\text{MgO}(\text{S})\xrightarrow{{{2800}^{{}^\circ }}\text{C}}\text{MgO}(\text{S})$ (39)

图 12所示的群炉的用电功率$p(k)$

$\begin{align} &p(k)=\sum\limits_{i=1}^{m}{{{p}_{i}}(k)=} \\ &\sqrt{3}\cdot \sum\limits_{i=1}^{m}{\left( {{I}_{i}}^{\text{2}}(k)\cdot {{R}_{i}}(k)\cdot \cos {{\varphi }_{i}} \right)},i=1,2,\cdots ,m \\ \end{align}$ (40)

其中,${p_i}(k)$为第i台电熔镁炉在k时刻的功率,m为由同一主变压器供电的电熔镁炉的台数,${I_i}(k)$为第i台电熔镁炉在k时刻的熔炼电流,${R_i}(k)$ 为第i台电熔镁炉在k时刻的工作电阻,${U_i}$为第i台电熔镁炉的熔炼电压,熔炼过程中${U_i}$ 基本不变,$\cos {\varphi _i}$ 为第i台电熔镁炉的功率因数,${\varphi_i}$表示第i台电熔镁炉的熔炼电压${U_i}$ 与熔炼电流${I_i}(k)$之间的相位差.

根据电路原理可知电熔镁炉的工作电阻${R_i}(k)$ 等于电弧电阻$R_i^a(k)$和熔池电阻$R_i^b(k)$之和[23-24]:

${{R}_{i}}(k)=R_{i}^{a}(k)+R_{i}^{b}(k)$ (41)

其中,电弧电阻$R_i^a(k)$由电弧电阻率${\rho}_i^a(k)$、弧长$L_i^a(k)$、电弧半径$r_i^a(k)$决定,而${\rho}_i^a(k)$ 与电弧温度${T_1}(k)$ 有关.间接影响$R_i^a(k)$的因素主要有三相电机控制量${{u}_{ij}}(k)(i=1,2,\cdots ,m.j=1,2,3)$ 熔体电阻率${\rho}_i^b(k)$,原料熔点${T_{rm}}$; 熔体杂质成分${B_{rm}}$. 因此,工作电阻$R_i^a(k)$的非线性特性较强,可表示成如下形式:

$R_{i}^{a}(k)={{f}_{i}}\left( {{u}_{i1}}(k),{{u}_{i}}_{2}(k),{{u}_{i}}_{3}(k),{{T}_{rm}},\rho _{i}^{b}(k),{{B}_{rm}} \right)$ (42)

其中,${f_i}\left( \cdot \right)$为未知非线性函数.

熔池电阻$R_i^b(k)$ 具有如下关系式[21]:

$R_{i}^{b}(k)=\frac{\rho _{i}^{b}(k)}{\pi d}\left( 1-\frac{d}{2h_{i}^{b}(k)} \right)$ (43)

其中,$\rho_i^b(k)$K 时刻炉内氧化镁熔体的电阻率,d为电机直径,$h_i^b(k)$ 为熔体层深度.

由式(42) 和式(43) 可以看出,电熔镁炉的工作电阻${R_i}(k)$具有非线性特性. 设工作电阻${R_i}(k)$ 的变化量$\Delta {R_i}(k)$

$\Delta {{R}_{i}}(k)={{R}_{i}}(k)-{{R}_{i}}(k-1)$ (44)

由式(40) 和式(44) 可得功率p(k) 的变化量$\Delta p(k)$

$\begin{array}{*{35}{l}} \Delta p(k)=\sqrt{3}\cdot \sum\limits_{i=1}^{m}{\left( \frac{-\Delta {{R}_{i}}(k){{U}_{i}}^{\text{2}}\cos {{\varphi }_{i}}}{{{R}_{i}}(k){{R}_{i}}(k-1)} \right)}= \\ \sqrt{3}\cdot \sum\limits_{i=1}^{m}{\left( -\Delta {{R}_{i}}(k){{I}_{i}}(k){{I}_{i}}(k-1)\cos {{\varphi }_{i}} \right)} \\ \end{array}$ (45)

由式(40) 和式(45) 得到p(k)与p(k-1)、$\Delta{R_i}(k)$${I_i}(k)$ 之间的关系:

$\begin{align} &p(k)=p(k-1)+ \\ &\sqrt{3}\cdot \left( \sum\limits_{i=1}^{m}{\left( -\Delta {{R}_{i}}(k){{I}_{i}}(k){{I}_{i}}(k-1)\cos {{\varphi }_{i}} \right)} \right) \\ \end{align}$ (46)

由式(46) 可以看出,p(k)由p(k-1)、$\Delta{R_i}(k)$${I_i}(k)$${I_i}(k-1) $ 以及$\cos {\varphi _i}$决定.由于$\Delta {R_i}(k)$${I_i}(k)$ 难以建模,因此,p(k)的动态特性预报非常困难. 式(46) 也蕴含着当前时刻的功率p(k)与过去时刻的功率$p(k - 1) ,p(k - 2) ,\cdots,p(k - {n_p})$以及过去时刻工作电阻以及工作电流的变化规律.由于电阻和电流无法建模,但电弧电阻可以通过过去时刻的功率来间接代替,电流是闭环控制系统的输出,可通过量测得到. 因此,以功率与电流作为功率模型的输入,p(k)可写成如下非线性函数:

$\begin{array}{*{35}{l}} p(k)={{f}_{P}}\left( p(k-1),\cdots ,p(k-{{n}_{p}}), \right. \\ ~\left. I(k-1),\cdots ,I(k-{{n}_{I}}) \right) \\ \end{array}$ (47)

${f_p}\left( \cdot \right)$是描述 p(k)与$R_i^a(k)$ 以及I(k)之间的非线性函数. 由式(46) 和式(47) ,可以利用功率p(k) 和I(k)的历史数据建立p(k+1)的低阶线性输入输出模型.该模型的输入是前${n_p}$ 个过去时刻的功率和${m_p}$ 个过去时刻的电流,模型输出为$k+1$ 时刻的功率预报值$\hat p(k{\rm{ + 1}})$; 因此,功率预报模型可以表示成如下形式:

$\begin{align} &p(k+1)={{f}_{P}}[p(k),\cdots ,p(k-{{n}_{p}}+1), \\ &I(k),\cdots ,I(k-{{m}_{I}})] \\ \end{align}$ (48)

式中,I(k)和p(k) 分别为电熔镁炉在k时刻的输入和输出;模型阶次${n_p}$${m_I}$ 未知;${f_P}( \cdot )\in {\bf R}$是未知的非线性函数.

下面首先利用功率p(k) 和电流I(k)的历史数据建立p(k+1)动态特性的低阶线性输入输出模型,则

$\begin{align} &{{{\bar{p}}}^{*}}(k+1)=-{{a}_{\text{1}}}p(k)-{{a}_{\text{2}}}p(k-1)\cdots \\ &-{{a}_{{{n}_{A}}}}p(k-{{n}_{\text{A}}}+1)+{{b}_{0}}I(k)\text{+}\cdots + \\ &{{b}_{{{n}_{B}}}}I(k-{{n}_{B}})=z[1-A({{z}^{-1}})]p(k)+ \\ &B({{z}^{-1}})I(k) \\ \end{align}$ (49)

式中,${\bar p^*}(k + 1) $ 表示该模型的输出,$A({z^{ - 1}})$是阶次为${{n}_{A}}({{n}_{A}}\le {{n}_{p}})$$z^{ - 1}$的多项式,$B({z^{ -1}})$ 是阶次为$\begin{array}{*{35}{l}} {{n}_{B}}({{n}_{B}}\le {{M}_{I}}) \\ \end{array}$$z^{ - 1}$ 的多项式,$A({z^{ - 1}})$$B({z^{ - 1}})$分别为

$A({{z}^{-1}})=1+{{a}_{1}}{{z}^{-1}}+\cdots +{{a}_{{{n}_{A}}}}{{z}^{-{{n}_{A}}}}$ (50)
$B({{z}^{-1}})={{b}_{0}}+{{b}_{1}}{{z}^{-1}}+\cdots +{{b}_{{{n}_{B}}}}{{z}^{-{{n}_{B}}}}$ (51)

由于模型(49) 并没有考虑电熔电阻、熔炼电流的非线性变化所产生的功率,因此,模型(49) 的输出与式(48) 输出之间必然不匹配,所产生的未建模动态为

$v[\varphi (k)]={{f}_{P}}(\cdot )-{{{\bar{p}}}^{*}}(k+1)$ (52)

式中,$v{\rm{[}}{\bf\pmb{{\varphi}}}{\rm{(}}k{\rm{)]}}$为未知的非线性连续向量函数,$\mathrm{ }\!\!\varphi\!\!\text{ }(k)=[p(k),\cdots ,p(k-{{n}_{A}}\text{+1})$,${I(k)},\cdots$,${I{\rm{(}}k -{n_B}{\rm{)]}}^{\rm{T}}}$.

由式(48) 、式(49) 和式(52) 可得:

$\begin{align} &p(k+1)={{{\bar{p}}}^{*}}(k+1)+v[\varphi (k)]= \\ &-{{a}_{\text{1}}}p(k)-\cdots -{{a}_{{{n}_{a}}}}p(k-{{n}_{\text{A}}}+1)+ \\ &{{b}_{\text{0}}}I(k)+\cdots +{{b}_{{{n}_{B}}}}I(k-{{n}_{B}})+v[\varphi (k)] \\ \end{align}$

即:

$A({{z}^{-1}})p(k\text{+}1)=B({{z}^{-1}})I(k)+v[\varphi (k)]$ (53)

针对式(53) ,利用采集到的5657组工业数据,从中选取1000组数据进行仿真实验,其中500组用来预报,500组用来检验.采用本文所提的交替辨识算法对功率进行预报,实验结果如图 13~17所示.首先对一个熔炼批次内的功率的时间序列作出ACF (Auto-correlation function)和PACF (Partial autocorrelation function)分析,结果如图 13所示. 从图中可以看出,自相关函数呈截尾特性,偏相关函数在第3个点之后PACF值很小.

图 13 功率的相关性分析 Figure 13 The correlation analysis of power
图 14 功率(星线) 及其估计值(圈线) Figure 14 The power (star line) and its estimation value (circle line) by the proposed method
图 15 功率的估计误差 Figure 15 Estimation error of power
图 16 估计误差(残差) 的概率密度函数 Figure 16 The probability density functions of prediction errors (residuals)
图 17 BP 神经网络的性能、训练状况和相关系数 Figure 17 The performance, training state, and regression of BP neural networks

采用本文的交替辨识算法得到电熔镁炉功率的预报值和真实值如图 14所示,预报误差如图 15所示. 从图 15中可以看出,大部分功率预报值与真实值之间的误差在[-600,600]内,在该批次中真实需量小于需量限幅值21500kW.进一步计算辨识误差的平均绝对值误差、方差以及标准差,其具体数值如表 1所示.

表 1 辨识的平均绝对值误差、方差以及标准差 Table 1 The average absolute error,variance and the standard deviation which produced by the proposed method

通过表 1可以看出,平均绝对值误差和绝对误差都在工业要求的区间范围内.由此可见本文所提的估计算法是有效的.

图 16是采用交替辨识算法的误差分布的概率密度函数直方图,从图 16中可以看出,相对误差或者残差具有高斯型的概率密度函数,说明估计效果较好. 因此,本文的辨识算法预报功率是有效的.

在实验中,BP神经网络的训练情况如图 17所示.

图 17中可以明显地看出,网络在7次迭代后就达到了最优的性能332666.3414.从数据的回归曲线可以看出,估计值与真值之间的相关系数为0.87425,说明采用本文辨识算法所得到的估计值与功率的真实值相关性很大,保证了估计方向的正确性以及估计精度.

4 结论

本文针对复杂工业过程中存在的一类运行指标监控系统难以准确预报运行指标的问题,把过程控制系统输出对运行指标模型的影响综合考虑,建立了监控系统运行指标的动态预报模型,并在工作点附近将这类预报模型化为线性模型与未建模动态和的形式,对线性模型的参数以及未建模动态提出了一种由改进的投影算法与基于BP神经网络的未建模动态估计算法组成的交替辨识算法,给出了算法的具体实现步骤.最后通过数值仿真实验以及在电熔镁炉功率监控系统的功率预报实验,实验结果说明了所提算法能有效地预报运行指标的动态变化.

参考文献
1 Chai Tian-You, Ding Jin-Liang, Wang Hong, Su Chun-Yi. Hybrid intelligent optimal control method for operation of complex industrial processes. Acta Automatica Sinica, 2008, 34 (5): 505–515.
( 柴天佑, 丁进良, 王宏, 苏春翌. 复杂工业过程运行的混合智能优化控制方法. 自动化学报, 2008, 34 (5): 505–515. )
2 Chai Tian-You. Operational optimization and feedback control for complex industrial. Acta Automatica Sinica, 2013, 39 (11): 1744–1757.
( 柴天佑. 复杂工业过程运行优化与反馈控制. 自动化学报, 2013, 39 (11): 1744–1757. DOI:10.3724/SP.J.1004.2013.01744 )
3 Fu Jun, Chai Tian-You, Su Chun-Yi, Xie Wen-Fang. Adaptive output tracking control of a class of nonlinear systems. Control Engineering of China, 2015, 22 (4): 731–736.
4 Wang Ai-ping, Wang Hong. Performance analysis for operational optimal control for complex industrial processesthe square impact principle. Control Engineering of China, 2013, 20 (6): 991–995.
( 王爱平, 王宏. 复杂工业系统运行优化品质分析——平方影响度原理. 控制工程, 2013, 20 (6): 991–995. )
5 Li H X, Guan S P. Hybrid intelligent control strategy. Supervising a DCS-controlled batch process. IEEE Control System, 2001, 21 (3): 36–48. DOI:10.1109/37.924796
6 Diaz S, Luque J, Romero M C, Escudero J I. Power systems monitoring and control using telecom network management standards. IEEE Transactions on Power Delivery, 2005, 20 (2): 1349–1356. DOI:10.1109/TPWRD.2004.833918
7 Nima A, Farshid K, Hamidreza Z. Wind power prediction by a new forecast engine composed of modified hybrid neural network and enhanced particle swarm optimization. IEEE Transactions on Sustainable Energy, 2011, 2 (3): 265–276. DOI:10.1109/TSTE.2011.2114680
8 Gyan R B, Maheshwari R P, Dewal M L. Modeling, control, and monitoring of S3RS-based hydrogen cooling system in thermal power plant. IEEE Transactions on Industrial Electronics, 2012, 59 (1): 562–570. DOI:10.1109/TIE.2011.2134059
9 Khosravi A, Nahavandi S, Creighton D. Prediction intervals for short-term wind farm power generation forecasts. IEEE Transactions on Sustainable Energy, 2013, 4 (3): 602–610. DOI:10.1109/TSTE.2012.2232944
10 Paparoditis E, Sapatinas T. Short-term load forecasting:the similar shape functional time-series predictor. IEEE Transactions on Power Systems, 2013, 28 (4): 3818–3825. DOI:10.1109/TPWRS.2013.2272326
11 Quan H, Srinivasan D, Khosravi A. Short-term load and wind power forecasting using neural network-based prediction intervals. IEEE Transactions on Neural Networks and Learning Systems, 2014, 25 (2): 303–315. DOI:10.1109/TNNLS.2013.2276053
12 Ceperic E, Ceperic V, Baric A. A strategy for short-term load forecasting by support vector regression machines. IEEE Transactions on Power Systems, 2013, 28 (4): 4356–4364. DOI:10.1109/TPWRS.2013.2269803
13 Atalik T, Cadirci I, Demirci T, Ermis M, Inan T, Kalaycioglu A S, Salor O. Multipurpose platform for power system monitoring and analysis with sample grid applications. IEEE Transactions on Instrumentation and Measurement, 2014, 63 (3): 566–582. DOI:10.1109/TIM.19
14 Huang T C, Du M J, Kang Y C, Peng R H, Chen K H, Lin Y H, Tsai T Y, Lee C C, Chen L D, Chen J L. 120% percent harvesting energy improvement by maximum power extracting control for high sustainability magnetic power monitoring and harvesting system. IEEE Transactions on Power Electronics, 2015, 30 (4): 2262–2274. DOI:10.1109/TPEL.2014.2330868
15 Zhang Y, Beaudin M, Taheri R, Zareipour H, Wood D. Day-ahead power output forecasting for small-scale solar photovoltaic electricity generators. IEEE Transactions on Smart Grid, 2015, 6 (5): 2253–2262. DOI:10.1109/TSG.2015.2397003
16 Rumelhart D E, Hinton G E, Williams R J. Learning internal representations by error propagation. Parallel Distributed Processing. Cambridge:MIT Press, 1986 : 318–362.
17 Chen L J, Narendra K S. Nonlinear adaptive control using neural networks and multiple models. Automatica, 2001, 37 (8): 1245–1255. DOI:10.1016/S0005-1098(01)00072-3
18 Chen S, Cowan C F N, Grant P M. Orthogonal least squares learning algorithm for radial basis function networks. IEEE Transactions On Neural Networks, 1991, 2 (2): 302–309. DOI:10.1109/72.80341
19 Li Shi-Zhe, Wang Yin-Song, Zhang Tao, Tian Jing-Yu. Performance assessment of load control system based on time domain index. Control Engineering of China, 2016, 23 (3): 864–869.
( 李士哲, 王印松, 张涛, 田靖雨. 基于时域指标的负荷控制系统性能评价. 控制工程, 2016, 23 (3): 864–869. )
20 Du K L, Swamy M N S. Radial basis function networks. Neural Networks and Statistical Learning. London:Springer-Verlag, 2014.
21 Wu Zhi-Wei, Chai Tian-You, Wu Yong-Jian. A Hybrid prediction model of energy consumption per ton for fused magnesia. Acta Automatica Sinica, 2013, 39 (12): 2002–2011.
( 吴志伟, 柴天佑, 吴永建. 电熔镁砂产品单吨能耗混合预报模型. 自动化学报, 2013, 39 (12): 2002–2011. )
22 Guo Mao-Xian. Industry Furnace. Beijing: Metallurgical Industry Press, 2002.
( 郭茂先. 工业电炉. 北京: 冶金工业出版社, 2002. )
23 Zhao W X, Chen H F, Zheng W X. Recursive identification for nonlinear ARX systems based on stochastic approximation algorithm. IEEE Transactions on Automatic Control, 2010, 55 (6): 1287–1299. DOI:10.1109/TAC.2010.2042236
24 Novak A, Simon L, Kadlec F, Lotton P. Nonlinear system identification using exponential swept-sine signal. IEEE Transactions on Instrumentation and Measurement, 2010, 59 (8): 2220–2229. DOI:10.1109/TIM.2009.2031836