-
摘要: 运用数据挖掘技术进行铁路事故类型预测及成因分析, 对于建立铁路事故预警机制具有重要意义. 为此, 本文提出一种基于梯度提升决策树(Grandient boosting decision tree, GBDT)的铁路事故类型预测及成因分析算法. 针对铁路事故记录数据缺失的问题, 提出一种基于属性分布概率的补全算法, 最大程度保持原有数据分布, 从而降低数据缺失对事故类型预测造成的影响. 针对铁路事故记录数据类别失衡的问题, 提出一种集成的GBDT模型, 完成对事故类型的鲁棒性预测. 在此基础上, 根据GBDT预测模型中特征重要度排序, 实现事故成因分析. 通过在开放数据库上进行实验, 验证了本文模型的有效性.Abstract: The application of data mining technology in railway accident type prediction and cause analysis is of great significance to establish railway accident early warning mechanism. This paper proposes a gradient boosting decision tree (GBDT) based algorithm for railway accident type prediction and cause analysis. In order to solve the problem of data missing in railway accident record dataset, we propose a new data complement algorithm based on the attribute distribution probability, which can keep the distribution of original data as much as possible, thus reducing the impact of data missing on predicting railway accident type. To reduce the impact of unbalanced categories of data in railway accident dataset, an ensemble GBDT model is proposed to predict the types of accidents effectively and robustly. On these bases, according to the importance of features in GBDT prediction model, we complete the cause analysis of railway accidents. Experimental results on an open database show that our proposed method can predict the types and causes of railway accidents effectively.
-
近年来, 经济模型预测控制(Economic model predictive control, EMPC)在工业界和学术界引起了广泛关注[1-2]. 作为一种新近发展的先进控制技术, EMPC有望成为解决复杂系统节能、降耗和增效优化控制问题的重要手段, 目前已应用于能源、造纸、车辆等系统的能效优化控制[3-9]. 除具有常规模型预测控制 (Model predictive control, MPC)的显式处理约束和多变量控制的优点外, EMPC还能优化“经济”类目标函数, 通常这类函数不是设定值跟踪偏差的正定函数, 而是系统状态和控制变量的非凸或非正定函数[1-9]. 因此, 把以设定值跟踪偏差的正定函数为优化目标的常规MPC称为目标跟踪MPC, 而不以跟踪偏差的正定函数为优化目标的MPC统称为经济MPC[2]. 现有研究表明: 经济最优性目标与闭环系统的稳定性目标具有一定的冲突性[1-2], 因此近年来EMPC的稳定性综合策略得到了广泛研究.
为建立EMPC关于经济平衡点的稳定性, 一种主要方法是构造基于经济优化目标函数的Lyapunov函数[5, 10-17]. 例如, 使用终端等式约束和强对偶性假设, 定义经济目标函数的旋转代价函数并将其作为闭环系统的一个Lyapunov函数[10], 而引入广义终端等式约束[15], 建立了经济性能变化下的递推可行性与闭环系统的有界稳定性[11]. 进一步, 采用严格耗散性条件、终端不等式约束和终端代价函数代替, 降低了EMPC稳定性综合策略的保守性[12-13]. 在无终端约束EMPC策略中, 闭环轨迹在足够长的预测时域情况下收敛到平衡点的邻域[16-17]. 虽然无终端约束增大了闭环系统的吸引域, 但长时域预测将大大增加了在线优化的计算负担. 进一步, EMPC稳定性和经济性是一对存在冲突的控制目标[18-19], 且稳定性和经济性目标无法统一度量, 难以通过权重标定. 对此, 从多目标优化控制角度, 考虑非线性系统强对偶性或耗散性条件难以满足情况, 文献[20-21]通过构造稳定性收缩约束, 建立闭环系统关于最优经济平衡点的渐近稳定性.
实际系统总是存在不确定扰动, 现有EMPC策略通常难以保证受扰系统的可行性和稳定性. 对于目标跟踪MPC, 目前已有较多鲁棒稳定性结果[22-32], 主要包括本质鲁棒MPC[22]、Tube鲁棒MPC[24-25]以及min-max MPC[26-32]等, 其中min-max MPC采用微分对策原理, 在使最坏扰动输入情况下系统的性能指标上界达到最小. 相比于本质鲁棒MPC和Tube鲁棒MPC, min-max MPC能大大降低鲁棒MPC的保守性, 但会增加优化问题的在线计算量[1]. 为降低min-max MPC的在线计算量, 文献[30]采用仿射输入结构, 使MPC含有抑制扰动的闭环成分和易于求解的开环优化. 另一方面, 输入到状态稳定性(Input-to-state stability, ISS)成为分析不确定系统鲁棒稳定性的一个有效工具[23, 27-32], 并应用到了EMPC鲁棒性研究, 如文献[33-34]采用强对偶性假设和约束紧缩方法, 证明了周期性扰动下线性系统EMPC闭环收敛性, 文献[35]获得了非线性系统EMPC的有界稳定性结果, 提高了经济性能优化的灵活性, 文献[36]将稳定性目标和经济性目标相加, 证明EMPC线性系统关于经济目标的最大值是ISS的, 文献[37]施加保证鲁棒稳定性的显式收缩约束, 提出两种非线性系统鲁棒EMPC算法, 文献[38]提出Lipschitz连续非线性系统的隐式收缩鲁棒EMPC策略, 提高了系统的平均经济性能.
本文针对含有未知有界扰动的不确定非线性系统, 提出一种新的具有递推可行性以及ISS保证的鲁棒EMPC策略. 该策略明确考虑经济最优性和鲁棒稳定性控制目标的矛盾特点, 采用微分对策原理在线滚动优化计算这对冲突目标的min-max问题. 离线计算最优经济平衡点, 并利用状态与该平衡点的偏差定义鲁棒稳定性目标函数, 而经济目标函数则由系统的经济性能给定. 通过特殊设计EMPC优化问题的隐式收缩约束, 并在鲁棒稳定性目标优化问题中引入一个新约束, 保证EMPC优化的递推可行性和闭环系统关于不确定扰动输入的ISS. 相比现有鲁棒EMPC策略, 本文首先建立了约束非线性系统具有ISS的鲁棒EMPC策略; 其次, EMPC递推可行性和鲁棒稳定性无需强对偶性或耗散性假设条件, 从而扩大了鲁棒 EMPC的应用范围; 最后, 采用微分对策原理得到了保守性更低的容许扰动上界. 采用一个受扰非线性连续搅拌釜反应器(Continuous stirred tank reactor, CSTR)的仿真实例, 验证本文提出策略的有效性与优越性.
符号说明:
${{{\bf{Z}}} ^ + }$ 表示非负整数集,${{{\bf{I}}} _{a:b}}$ 表示集合$\{ i \in {{\mathbf{Z}}^ + }:a \leq i \leq b,\;a \in {{\mathbf{Z}}^ + },\;b \in {{\mathbf{Z}}^ + }\}$ ,${{{\bf{I\,}}} _{ \geq j}}$ 表示集合$\{ i \in {{\mathbf{Z}}^ + }:i \geq j,\;j \in {{\mathbf{Z}}^ + }\}$ .${{{\boldsymbol{x}}} ^ + }$ 表示${{\boldsymbol{x}}}$ 的下一时刻状态,$\left| {{\boldsymbol{x}}} \right|$ 表示${{\boldsymbol{x}}}$ 的欧几里得范数,$\left\| {{\boldsymbol{x}}} \right\| = {{{\rm{sup}}\{ }}\left| {{{\boldsymbol{x}}} \left( k \right)} \right|{{,}}\; $ $ k \in {{\mathbf{Z}}^ + }{{\} }}$ ,${{\boldsymbol{x}}} \left( {i|k} \right)$ 表示在第$ k $ 时刻对未来第$ k + i $ 时刻的预测变量. 连续函数$ {h_1}:{{\mathbf{R}}^ + } \to {{\mathbf{R}}^ + } $ 称为$ K $ 类函数系指该函数单调递增, 且$ {h_1}\left( 0 \right) = 0 $ ; 函数${h_2}:{{\mathbf{R}}^ + } \to $ $ {{\mathbf{R}}^ + }$ 称为$ {K_\infty } $ 类函数系指该函数是$ K $ 类函数, 且当$s \to $ $ \infty$ 时, 有$ h\left( s \right) \to \infty $ ; 函数${h_3}:{{\mathbf{R}}^ + } \times {{\mathbf{R}}^ + } \to $ $ {{\mathbf{R}}^ + }$ 称为$ KL $ 函数系指对于任意固定的$ t \geq 0 $ ,$ {h_3}\left( { \cdot ,t} \right) $ 是$ K $ 类函数, 而对于任意固定的$ s > 0 $ ,$ {h_3}\left( {s, \cdot } \right) $ 是单调递减, 且当$ t \to \infty $ 时,$ {h_3}\left( {s,t} \right) \to 0 $ .1. 问题描述
考虑不确定离散时间非线性系统
$$ \left\{ {\begin{aligned} & {{\boldsymbol{x}}(k + 1) = {f_1}({\boldsymbol{x}}(k)) + {f_2}({\boldsymbol{x}}(k)){\boldsymbol{u}}(k) + {f_3}({\boldsymbol{x}}(k)){\boldsymbol{w}}(k)} \\ & {{\boldsymbol{z}}(k) = \left[ {\begin{array}{*{20}{c}} {g{{(}}{\boldsymbol{x}}(k))} \\ {{\boldsymbol{u}}(k)} \end{array}} \right],{{ }}\quad \quad k \in {{\mathbf{{\rm {\bf{Z}}}}}^ + }} \end{aligned}} \right. $$ (1) 其中,
${{\boldsymbol{x}}} \left( k \right) \in {{\mathbf{R}}^n}$ ,${{\boldsymbol{u}}} \left( k \right) \in {{\mathbf{R}}^m}$ ,${{\boldsymbol{w}}} \left( k \right) \in {{\mathbf{R}}^q}$ 和${{\boldsymbol{z}}} \left( k \right) \in $ $ {{\mathbf{R}}^s}$ 分别表示系统在k时刻的状态、控制输入、扰动输入和辅助输出. 假设$ {f_1}\left( \cdot \right) $ ,$ {f_2}\left( \cdot \right) $ ,$ {f_3}\left( \cdot \right) $ 和$ g\left( \cdot \right) $ 分别为定义在$ {{\mathbf{R}}^n} $ 上的光滑函数, 满足$ {f_1}\left( 0 \right) = 0 $ ,${f_2}\left( 0 \right) = $ $ 0$ ,$ {f_3}\left( 0 \right) = 0 $ 以及$ g\left( 0 \right) = 0 $ . 系统(1)在$ k $ 时刻的状态${\boldsymbol{x}}\left( {{k}} \right) = \varphi \left( {{{k}};{{\boldsymbol{x}}_0},{\boldsymbol{u}},{\boldsymbol{w}}} \right)$ 由初始状态$ {{\boldsymbol{x}}_0}, $ 控制序列${\boldsymbol{u}} = \{ {\boldsymbol{u}}\left( 0 \right),{\boldsymbol{u}}\left( 1 \right), \cdots \}$ 和扰动序列${{\boldsymbol{w}}} = \{ {\boldsymbol{w}}\left( 0 \right), $ $ {\boldsymbol{w}}\left( 1 \right), \cdots \}$ 表示. 假设系统(1)的状态完全可测, 且系统状态和控制输入满足约束$${\boldsymbol{x}}(k) \in X,\;{\rm{ }}{\boldsymbol{u}} (k) \in U,\;{\rm{ }}k \in {{\bf{Z}}^ + }$$ (2) 其中, 集合
$ X \subset {{\mathbf{R}}^n} $ 和$ U \subset {{\mathbf{R}}^m} $ 均为凸的紧集, 且它们内部包含某些平衡点$ \left( {{{\boldsymbol{x}}_e},\;{{\boldsymbol{u}}_{e} }} \right) $ .假设 1. 扰动输入
$ {\boldsymbol{w}} $ 满足如下约束$${\boldsymbol{w}}(k) \in W,\;\;{\rm{ }}k \in {{\bf{Z}}^ + }$$ (3) 其中,
$ W $ 为包含平衡点$\left( {{{\boldsymbol{x}}_e},\;{{\boldsymbol{u}}_{e} }} \right)$ 的紧集. 满足约束(3)的扰动称为容许扰动. 注意$ {\boldsymbol{w}} $ 表示参数不确定性、模型失配以及持续外部扰动等多种有界不确定性[29], 其上界表示为$\left\| {\boldsymbol{w}} \right\| = {{\sup\{ }}\left| {{{\boldsymbol{w}}_{t} }} \right|{{:}}{{\boldsymbol{w}}_{t} } \in {{\boldsymbol{W}}} {{,}}\;{t} \in $ $ {{\mathbf{Z}}^ + }{{\} }}$ .考虑系统(1)的经济性能函数
$ {{L} _e}{{:}}{X} \times {U} \to {\mathbf{R}} ,$ 基于系统(1)的名义模型${{\boldsymbol{x}}^ + } = {{{f}}_1}\left( {\boldsymbol{x}} \right) + {{{f}}_2}\left( {\boldsymbol{x}} \right){\boldsymbol{u}}$ 离线计算如下最优经济平衡点$$\begin{split} & ({\boldsymbol{x}}_e^ * ,{\boldsymbol{u}}_e^ * ) = \arg \mathop {\min }\limits_{x,u} {L_e}({\boldsymbol{x}},{\boldsymbol{u}}) \\ &{\rm{s}}{\rm{.t}}{\rm{. }}\;\;\;\;{\boldsymbol{x}} = {f_1}{\rm{(}}{\boldsymbol{x}}) + {f_2}{\rm{(}}{\boldsymbol{x}}{\rm{)}}{\boldsymbol{u}} \\ &\,{\rm{ }}\quad \quad \;{\boldsymbol{x}} \in X, {\boldsymbol{u}} \in U \\ \end{split} $$ (4) 不失一般性, 后文假设系统(1)的最优经济平衡点为原点.
下面介绍输入到状态稳定(ISS)相关理论, 首先回顾鲁棒不变集的概念. 考虑一般离散不确定非线性系统
$ {{\boldsymbol{x}}^ + } = {{f}}\left( {{\boldsymbol{x}},{\boldsymbol{w}}} \right) $ .定义 1[27]. 考虑系统
${{\boldsymbol{x}}^ + } = {{f}}\left( {{\boldsymbol{x}},{\boldsymbol{w}}} \right)$ 和集${{{\boldsymbol{\Theta }}}} \subseteq {{{\bf{R}}} ^n}$ , 如果对于任意${{\boldsymbol{x}}} \in {{{\boldsymbol{\Theta }}}}$ 和${\boldsymbol{w }}\in {W}$ , 该系统满足${{\boldsymbol{x}}^ + } = $ $ {{f}}\left( {{\boldsymbol{x}},{\boldsymbol{w}}} \right) \in {{{\boldsymbol{\Theta}} }}$ , 则称${{{\boldsymbol{\Theta}} }}$ 为该系统的一个鲁棒不变集.定义 2[27]. 考虑系统
${{\boldsymbol{x}}^ + } = {{f}}\left( {{\boldsymbol{x}},{\boldsymbol{w}}} \right)$ 及其内含原点的鲁棒不变集${{{\boldsymbol{\Theta}} }} \subseteq X$ , 如果对于任意初始状态${{{\boldsymbol{x}}} _{{0}}} \in {{{\boldsymbol{\Theta}} }}$ 及${\boldsymbol{w}} \in {W}$ , 存在$ {K_\infty } $ 类函数$ \alpha $ 和$ KL $ 类函数$ \beta $ 使系统满足$$|\varphi (k;{{\boldsymbol{x}}_0},{\boldsymbol{w}})| \le \beta (|{{\boldsymbol{x}}_0}|,k) + \alpha (||{\boldsymbol{w}}||), {\rm{ }}\;\;\forall k \in {{\bf{Z}}^ + }$$ (5) 则该系统在
$ {{{\boldsymbol{\Theta}} }} $ 内是ISS的.引理 1[29]. 考虑系统
${{\boldsymbol{x}}^ + } = {{f}}\left( {{\boldsymbol{x}},{\boldsymbol{w}}} \right)$ 及其内含原点的鲁棒不变集${{{\boldsymbol{\Theta }}}} \subseteq X$ , 如果对于任意${{\boldsymbol{x}}} \in {{\boldsymbol{{\Theta}} }}$ 和${\boldsymbol{ w}} \in {W}$ , 存在连续函数$ V:{{\mathbf{R}}^n} \to {{\mathbf{R}}^ + } $ 满足$$ \begin{split} &{{\xi _1}(|{\boldsymbol{x}}|) \le V({\boldsymbol{x}}) \le {\xi _2}(|{\boldsymbol{x}}|) + {\rho _1}(||{\boldsymbol{w}}||)} \\ & {V({{\boldsymbol{x}}^ + }) - V({\boldsymbol{x}}) \le - {\xi _3}(|{\boldsymbol{x}}|) + {\rho _2}(|{\boldsymbol{w}}|)} \end{split} $$ (6) 其中,
${\xi _1} $ ,${\xi _2} $ 和${\xi _3} $ 为K∞类函数,${\rho _1} $ 和${\rho _2} $ 为K类函数, 则该系统在${{{\boldsymbol{\Theta }}}} $ 内ISS, 称V为该系统的ISS-Lyapunov函数.注 1. 由定义2可知, 当系统不受扰动或仅受衰减扰动作用时, 系统最终在原点处渐近稳定; 当受持续有界扰动作用时, 系统有界稳定, 且状态轨迹最终收敛的范围与持续扰动的上界有关.
本文目标是针对不确定非线性系统(1), 通过极小化经济性能函数在线计算鲁棒EMPC控制器, 要求相应闭环系统满足约束条件(2), 且闭环系统的最优经济平衡点相对于容许扰动(3)具有ISS.
2. 鲁棒经济模型预测控制
考虑有限预测步长
$N \in {{{\bf{I}}} _{ \geq 1}}$ , 定义$k \in {{\bf{Z}}^ + }$ 时刻的$ N $ 步控制序列${{\boldsymbol{u}}} ( {{k}} ) = \{{\boldsymbol{u}}( {0|{{k}}} ),{\boldsymbol{u}}( {1|{{k}}} ), \cdots , {\boldsymbol{u}}( {{{N}} - 1|{{k}}} )\}$ 、扰动序列${{\boldsymbol{w}}} ( {{k}} ) = \{{\boldsymbol{w}}( {0|{{k}}} ),\;{\boldsymbol{w}}( {1|{{k}}} ), \cdots ,\;{\boldsymbol{w}}( {{{N}} - 1|{{k}}} )\}$ 以及对应的预测状态序列${{\boldsymbol{x}}} ( {{k}} ) = {{\{ }}{\boldsymbol{x}}( {1|{{k}}} ){{,}}\; {\boldsymbol{x}}( {2|{{k}}} ){{,}} \cdots {{,}} $ $ \;{\boldsymbol{x}}( {{{N}}|{{k}}} ){{\} }}$ , 并且${\boldsymbol{x}}( {{{i|k}}} ) = \varphi ( {{{i}};{\boldsymbol{x}}( {{k}} ),{{{\boldsymbol{u}}}}( {{k}} ),{{{\boldsymbol{w}}}}( {{k}} )} )$ , 其中${\boldsymbol{x}}( {{k}} )$ 为当前时刻系统状态. 考虑系统(1)的经济性能函数$ {{L} _e} $ , 定义$ N $ 步经济目标函数$${J_e}({\boldsymbol{x}}(k),{\boldsymbol{u}}(k)) = \sum\limits_{i = 0}^{N - 1} {{L_e}({\boldsymbol{x}}(i|k),{\boldsymbol{u}}(i|k))} \tag{7}$$ 其中,
${\boldsymbol{x}}\left( {0{{|k}}} \right) = {\boldsymbol{x}}\left( {{k}} \right)$ 为计算鲁棒EMPC控制器, 在每个时刻优化经济目标函数(7). 由于扰动的存在, 在鲁棒EMPC中求解如下min-max经济最优控制问题:
$$\tag{8a}({{\boldsymbol{u}}^ * }(k),{{\boldsymbol{w}}^ * }(k)) = \arg \mathop {\min }\limits_{{\boldsymbol{u}}(k)} \mathop {\max }\limits_{{\boldsymbol{w}}(k)} {J_e}({\boldsymbol{x}}(k),{\boldsymbol{u}}(k))\;\;$$ $$\tag{8b}\begin{split} {\rm{s}}{\rm{.t}}{\rm{. }}\;\;{\boldsymbol{x}}(i + 1|k) =\;& {f_1}({\boldsymbol{x}}(i|k))\; + {f_2}({\boldsymbol{x}}(i|k)){\boldsymbol{u}}(i|k) \;+ \\ & {f_3}({\boldsymbol{x}}(i|k)){\boldsymbol{w}}(i|k) \\[-10pt] \end{split} $$ $$\tag{8c}{\boldsymbol{x}}(i|k) \in X,\quad {\boldsymbol{u}}(i|k) \in U,\quad {\boldsymbol{w}}(i|k) \in W\,\;\;\;\quad$$ $$\tag{8d}{\boldsymbol{x}}(0|k) = {\boldsymbol{x}}(k),\quad {\boldsymbol{x}}(N|k) \in {\boldsymbol{\Omega}} ,\;{\rm{ }}\forall i \in {{\bf{I}}_{0:N - 1}}\quad$$ $$\tag{8e}{V_R}(x(k),{\boldsymbol{u}}(k),{\boldsymbol{w}}(k)) \le \eta ({\boldsymbol{x}}(k),\lambda )\;\,\quad\quad\quad\quad$$ 其中,
$\left( {{{\boldsymbol{u}}^*}\left( {{k}} \right),{{\boldsymbol{w}}^*}\left( {{k}} \right)} \right)$ 为$ k $ 时刻的经济最优解, 对应最优预测状态序列${{\boldsymbol{x}}^*}\left( {{k}} \right)$ ;${\boldsymbol{x}}\left( {0{{|k}}} \right) = {\boldsymbol{x}}\left( {{k}} \right)$ 为初始条件;${\boldsymbol{x}}\left( {{{N|k}}} \right) \in {{{{{\boldsymbol{\Omega}}}} }}$ 为终端状态约束, 终端约束集${{{{\boldsymbol{\Omega}} }}} \subset $ $ {{X}}$ 且内含原点; (8e)为待设计的鲁棒稳定性收缩约束,$ \lambda \in \left[ {0,1} \right) $ 为收缩因子. 进一步定义如下关于最优经济平衡点的鲁棒稳定性目标函数:$$\begin{split} &{V_R}({\boldsymbol{x}}(k),{\boldsymbol{u}}(k),{\boldsymbol{w}}(k)) = {\rm{ }}E({\boldsymbol{x}}(N|k))\; +\\ &\qquad \sum\limits_{i = 0}^{N - 1} {\{ L({\boldsymbol{x}}(i|k),{\boldsymbol{u}}(i|k)) - {L_w}({\boldsymbol{w}}(i|k))\} } \end{split} $$ (9) 其中, 函数
$ {L} {{:}}{X} \times {U} \to {{\mathbf{R}}^ + } $ ,${{L} _{{w}}}{{:W}} \to {{\mathbf{R}}^ + }$ 和${E} {{:}}{X} \to $ $ {{\mathbf{R}}^ + }$ 为连续有界函数. 为书写简洁, 令${{V} _R}\left( {\boldsymbol{x}} \right) = {{V} _R}( {\boldsymbol{x}}, $ $ {\boldsymbol{u}},{\boldsymbol{w}})$ .假设 2. 存在
$ {K_\infty } $ 类函数$ {\alpha _l}, $ $ {\alpha _w} $ 和$ {\beta _w} ,$ 使鲁棒稳定性目标函数${{V} _R}\left( {\boldsymbol{x}} \right)$ 满足:${L} \left( {{\boldsymbol{x}},{\boldsymbol{u}}} \right) \geq {\alpha _l}\left( {|{\boldsymbol{x}}|} \right),\; \forall {\boldsymbol{x}} \in $ $ {X} ,\;\forall {\boldsymbol{u}} \in {U} ;$ ${\alpha _w}\left( {|{\boldsymbol{w}}|} \right) \leq {{L} _w}\left( {\boldsymbol{w}} \right) \leq {\beta _w}\left( {|{\boldsymbol{w}}|} \right),\;\forall {\boldsymbol{w}} \in {W}.$ 假设 3. 在终端约束集
${{{{{\boldsymbol{\Omega}}}} }}$ 内存在局部反馈控制律${\boldsymbol{u}} = {\boldsymbol{\pi}} \left( {\boldsymbol{x}} \right)$ 使得${\boldsymbol{\pi}} \left( {\boldsymbol{x}} \right) \in {U},\;\forall {\boldsymbol{x}} \in {X}$ 和如下不等式成立:$$\begin{split} E({{\boldsymbol{x}}^ + }) - E({\boldsymbol{x}}) \le \;& - L({\boldsymbol{x}},{\boldsymbol{\pi}} ({\boldsymbol{x}})) + {L_w}({\boldsymbol{w}}){\rm{ }} \\ &\forall {\boldsymbol{x}} \in {\boldsymbol{\Omega}} ,\;{\boldsymbol{w}} \in W \end{split} $$ (10) 其中,
${\alpha _f}\left( {|{\boldsymbol{x}}|} \right) \leq {E} \left( {\boldsymbol{x}} \right) \leq {\beta _{f} }\left( {|{\boldsymbol{x}}|} \right)$ ,$ {\alpha _f} $ 和$ {\beta _f} $ 为$ {K_\infty } $ 类函数.注意, 局部控制律
${\boldsymbol{u}} = {\boldsymbol{\pi}} \left( {\boldsymbol{x}} \right)$ 可采用如$ {H_\infty } $ 控制方法求解[26, 30]. 为此, 定义对称矩阵$${\boldsymbol{R}}({\boldsymbol{x}}) = \left[ {\begin{array}{*{20}{c}} {f_2^{\rm{T}}{{({\boldsymbol{x}})}}{\boldsymbol{P}}{f_2}({\boldsymbol{x}}) + {\boldsymbol{I}}}&{f_2^{\rm{T}}{{({\boldsymbol{x}})}}{\boldsymbol{P}}{f_3}({\boldsymbol{x}})} \\ {f_3^{\rm{T}}{{({\boldsymbol{x}})}}{\boldsymbol{P}}{f_2}({\boldsymbol{x}})}&{f_3^{\rm{T}}{{({\boldsymbol{x}})}}{\boldsymbol{P}}{f_3}({\boldsymbol{x}}) - {\gamma ^2}{\boldsymbol{I}}} \end{array}} \right]$$ 其中,
${{\boldsymbol{P}}}$ 为对称正定矩阵,${{\boldsymbol{I}}}$ 为单位矩阵. 令${{\boldsymbol{R}}} = $ $ {\boldsymbol{R}}\left( 0 \right)$ . 不失一般性, 令扰动抑制水平$ \gamma = 1 $ .引理 2[30]. 假设矩阵P满足Riccati不等式方程
$$\left\{ {\begin{aligned} &{{\boldsymbol{F}}_1}^{\rm{T}}{\boldsymbol{P}}{{\boldsymbol{F}}_1}-{\boldsymbol{P}}+{{\boldsymbol{G}}^{\rm{T}}}{\boldsymbol{G}}-\\ &{{\boldsymbol{A}}^{\rm{T}}}{\boldsymbol{P}}[{{F}_2},{{{F}}_3}]{{\boldsymbol{R}}^{ - 1}}{{[{{{F}}_2},{{{F}}_3}]}^{\rm{T}}}{\boldsymbol{P}}{{\boldsymbol{F}}_1} < 0 \\ &{{{\boldsymbol{G}}^{\rm{T}}}{\boldsymbol{P}}{\boldsymbol{G}} - {\boldsymbol{I}}<0} \end{aligned}} \right.$$ (11) 其中,
${{\boldsymbol{F}}} _1 = \dfrac{\partial {{f} _1}} {\partial {\boldsymbol{x}}}| _{{\boldsymbol{x}} \,= \,0}$ ,$ {{F} _2} = {{f} _2}\left( 0 \right) $ ,$ {{F} _3} = {{f} _3}\left( 0 \right) $ 和${{\boldsymbol{G}}} = $ $ {\dfrac{\partial {g} }{\partial{\boldsymbol{x}}}{ \left| {_{{\boldsymbol{x}} \,=\, 0}} \right.}}$ . 则存在集合${{{\boldsymbol{\Omega}} }} = {{\{ }}{\boldsymbol{x}} \in {{\mathbf{R}}^{n} }{{:}}{E} \left( {\boldsymbol{x}} \right) = {{\boldsymbol{x}}^{{{\rm{T}}}}}{{\boldsymbol{P}}} {\boldsymbol{x}} \leq $ $ {c} , \;{c} > 0{{\} }}$ 和局部控制律${\boldsymbol{u}} = {{{\boldsymbol{\pi}}}} \left( {\boldsymbol{x}} \right)$ $$\left[ {\begin{array}{*{20}{c}} {{\boldsymbol{\pi}} ({\boldsymbol{x}})} \\ {{\boldsymbol{w}}({\boldsymbol{x}})} \end{array}} \right] = - {\boldsymbol{R}}^{ - 1}{({\boldsymbol{x}})}{[{f_2}({\boldsymbol{x}}),{f_3}({\boldsymbol{x}})]^{\rm{T}}}{\boldsymbol{P}}{f_1}\left( {\boldsymbol{x}} \right)$$ (12) 使得
${\boldsymbol{\pi}} \left( {\boldsymbol{x }}\right) \in {U} ,\;\forall {\boldsymbol{x}} \in {{{\boldsymbol{\Omega }}}}$ 成立, 且局部闭环系统$${{\boldsymbol{x}}^ + } = {f_1}({\boldsymbol{x}}) + {f_2}({\boldsymbol{x}}){\boldsymbol{\pi}} ({\boldsymbol{x}}) + {f_3}({\boldsymbol{x}}){\boldsymbol{w}}$$ (13) 满足
$$E({{\boldsymbol{x}}^ + }) - E({\boldsymbol{x}}) \le |{\boldsymbol{w}}{|^2} - |{\boldsymbol{z}}{|^2},\quad \forall {\boldsymbol{x}} \in {\boldsymbol{\Omega}} ,\;\;{\boldsymbol{w}} \in W$$ (14) 注 2. 如果假设3成立, 则对于任意
$ {\boldsymbol{w}} \in W $ ,${{{{{\boldsymbol{\Omega}}}} }}$ 为闭环系统(13)的一个鲁棒不变集, 同时$E({\boldsymbol{x}})$ 为该系统在${{{{{\boldsymbol{\Omega}}}} }}$ 内的一个ISS-Lyapunov函数. 文献[30]给出了保证${{{{{\boldsymbol{\Omega}}}} }}$ 鲁棒不变性的一个充分条件, 即容许扰动${\boldsymbol{w}}$ 的上界$$||{\boldsymbol{w}}|| \le \frac{{\sqrt {{\lambda _{\min }}} - \sqrt {{\lambda _{\max }}} {l_x}}}{{{l_w}}}\sqrt {\frac{c}{{{\lambda _{\max }}{\lambda _{\min }}}}} $$ (15) 其中,
$ {l_x} $ 和$ {l_w} $ 分别为闭环系统(13)关于${\boldsymbol{ x}}$ 和${\boldsymbol{w }}$ 的局部Lipschitz常数,$ {\lambda _{\max }} $ 和$ {\lambda _{\min }} $ 分别为矩阵${{\boldsymbol{P}}}$ 的最大和最小特征值. 如扰动${\boldsymbol{w}}$ 满足式(15), 则${{{{{\boldsymbol{\Omega}}}} }}$ 为闭环系统(13)的鲁棒不变集, 且${\boldsymbol{\pi}} \left( {\boldsymbol{x}} \right) \in {U} ,\;\forall {\boldsymbol{x}} \in {{{\boldsymbol{\Omega}} }}$ .下面构造收缩约束函数
$\eta $ . 令x(k)为当前时刻系统状态, 求解如下有限时域鲁棒性最优控制问题$$\tag{16a}({{\boldsymbol{u}}^0}(k),{{\boldsymbol{w}}^0}(k)) = \arg \mathop {\min }\limits_{{\boldsymbol{u}}(k)} \mathop {\max }\limits_{{\boldsymbol{w}}(k)} {V_R}({\boldsymbol{x}}(k))\;\;\;\;\;\;$$ $$\tag{16b}\begin{split} {\rm{s}}{\rm{.t}}{\rm{. }}\;\;{\boldsymbol{x}}(i + 1|k) =\;& {f_1}({\boldsymbol{x}}(i|k)) + {f_2}({\boldsymbol{x}}(i|k)){\boldsymbol{u}}(i|k) + \\ & {f_3}({\boldsymbol{x}}(i|k)){\boldsymbol{w}}(i|k) \\[-10pt] \end{split} $$ $$\tag{16c}{\boldsymbol{x}}(i|k) \in X,\quad {\boldsymbol{u}}(i|k) \in U,\quad {\boldsymbol{w}}(i|k) \in W\;\;\;$$ $$\tag{16d}{\boldsymbol{x}}(0|k) = {\boldsymbol{x}}(k),\quad {\boldsymbol{x}}(N|k) \in {\boldsymbol{\Omega}} ,{\rm{ }}\forall i \in {{\bf{I}}_{0:N - 1}}$$ $$\tag{16e}{V_R}({\boldsymbol{x}}(k)) \le V_R^1({\boldsymbol{x}}(k)) + {\lambda _1}{L_w}(||{\boldsymbol{w}}||)\qquad\;$$ 其中,
$\left( {{{\boldsymbol{u}}^0}\left( {{k}} \right),{{\boldsymbol{w}}^0}\left( {{k}} \right)} \right)$ 为$ k $ 时刻的鲁棒最优解; 系数$ {\lambda _1} \geq 1 $ ;${V} _R^1( {{\boldsymbol{x}}( {k} )} ) = {{V} _R}( {{\boldsymbol{x}}( {k} ),\;{{\boldsymbol{u}}^1}( {k} ),\;{{\boldsymbol{w}}^1}( {k} )} )$ ,$( {{\boldsymbol{u}}^1}( {k} ), $ $ \;{{\boldsymbol{w}}^1}( {k} ) )$ 是$ k $ 时刻的一组序列对$$\begin{split} & ({{\boldsymbol{u}}^1}(k),{{\boldsymbol{w}}^1}(k)) = \\ &\quad\left\{ {\begin{array}{*{20}{l}} & ({{\boldsymbol{u}}^ * }(i + 1|k - 1),{{\boldsymbol{w}}^ * }(i + 1|k - 1)),&i \in {{\bf{I}}_{0:N - 2}} \\ & ({\boldsymbol{\pi}} ({{\boldsymbol{x}}^ * }(N|k)),{{\boldsymbol{w}}^1}(N - 1|k),&i = N - 1 \end{array}} \right. \end{split} $$ (17) 其中,
$ ( {{{\boldsymbol{u}}^*}( {{{i|k}} - 1} ),{{\boldsymbol{w}}^*}( {{{i|k}} - 1} )} ) $ 是优化问题(8)在$ k - $ $ 1 $ 时刻最优解$( {{\boldsymbol{u}}^*}( {{k}} - $ $ 1 ),{{\boldsymbol{w}}^*}( {{{k}} - 1} ) )$ 的分量. 将$( {{\boldsymbol{u}}^*}( {{k}} - $ $ 1 ),{{\boldsymbol{w}}^*}( {{{k}} - 1} ) )$ 和$( {{{\boldsymbol{u}}^0}( {{k}} ),{{\boldsymbol{w}}^0}( {{k}} )} )$ 分别代入${{V} _R}( {\boldsymbol{x}} )$ , 据此定义函数$$\begin{split} \eta ({\boldsymbol{x}}(k){\rm{,}}\lambda ) =\;& V_R^0({\boldsymbol{x}}(k)) + \lambda [V_R^*({\boldsymbol{x}}(k - 1)) - \\ &V_R^0({\boldsymbol{x}}(k)) + (1 + {\lambda _1}){L_w}(||{\boldsymbol{w}}||)] \end{split} $$ (18) 其中, 值函数
${V} _R^*( {{\boldsymbol{x}}( {{k} - 1} )} ) = {{V} _R}( {\boldsymbol{x}}( {{k} - 1} ), {{\boldsymbol{u}}^*}( {k} - 1 ), $ $ {{\boldsymbol{w}}^*}( {{k} - 1} ) )$ 和${V} _R^0( {{\boldsymbol{x}}( {k} )} ) = {{V} _R}( {{\boldsymbol{x}}( {k} ),{{\boldsymbol{u}}^0}( {k} ),{{\boldsymbol{w}}^0}( {k} )} )$ .如果优化问题(8)在k时刻可行, 则根据滚动时域控制原理, 将
${{\boldsymbol{u}}^*}\left( {{k}} \right)$ 的第1个分量定义为鲁棒EMPC控制律, 即${\boldsymbol{u}}( {{k}} ) = {{\boldsymbol{u}}^{{\rm{mpc}}}}( {{\boldsymbol{x}}( {{k}} )} ) = {{\boldsymbol{u}}^*}( {0|{{k}}} )$ , 对应闭环系统为$$\begin{array}{l} {\boldsymbol{x}}(k + 1) =\\ \qquad{f_1}({\boldsymbol{x}}(k)) + {f_2}({\boldsymbol{x}}(k)){u^{{\rm{mpc}}}}({\boldsymbol{x}}(k)) + {f_3}({\boldsymbol{x}}(k)){\boldsymbol{w}}(k) \end{array} $$ (19) 计算鲁棒EMPC控制器的算法总结如下:
算法 1. 鲁棒EMPC算法
步骤 1. 初始化. 设定
$N \in {{{\bf{I}}} _{ \geq 1}},$ $ \lambda \in \left[ {0,1} \right) ,$ $ {\lambda _1} \geq 1 ,$ ${{{L}}_{{e}}}\left( {{\boldsymbol{x}},{\boldsymbol{u}}} \right),$ ${{L}}\left( {x,u} \right),$ ${{{L}}_{{w}}}\left( {\boldsymbol{w}} \right)$ 和${{E}}\left( {\boldsymbol{x}} \right)$ ; 离线计算${\boldsymbol{\pi}} \left( {\boldsymbol{x}} \right)$ 以及${{{\boldsymbol{\Omega}} }}$ ; 令$ k = 0 $ , 考虑初始状态${{\boldsymbol{x}}_0}$ , 令$ \eta $ 充分大; 求解优化问题(8), 得到经济最优解$\left( {{{\boldsymbol{u}}^*}\left( 0 \right),{{\boldsymbol{w}}^*}\left( 0 \right)} \right)$ , 并将${{\boldsymbol{u}}^*}\left( 0 \right)$ 的首个分量${{\boldsymbol{u}}^*}\left( {0|0} \right)$ 作用于系统(1).步骤 2. 在
$ k $ 时刻, 利用$ k - 1 $ 时刻经济最优解$\left( {{\boldsymbol{u}}^*}\left( {{{k}} - 1} \right), {{\boldsymbol{w}}^*}\left( {{{k}} - 1} \right) \right)$ 构造$\left( {{{\boldsymbol{u}}^1}\left( {{k}} \right),{{\boldsymbol{w}}^1}\left( {{k}} \right)} \right)$ , 更新 (16e); 求解优化问题(16), 得到当前时刻鲁棒最优解$( {\boldsymbol{u}}^0( {k}),{{\boldsymbol{w}}^0}( {{k}} ) )$ .步骤 3. 计算
${V} _R^0\left( {{\boldsymbol{x}}\left( {k} \right)} \right)$ 和${V} _R^*\left( {{\boldsymbol{x}}\left( {{k} - 1} \right)} \right)$ , 更新$ \eta $ ; 求解优化问题(8), 得到当前时刻经济最优解$\left( {{{\boldsymbol{u}}^*}\left( {{k}} \right),{{\boldsymbol{w}}^*}\left( {t{k}} \right)} \right)$ .步骤 4. 将
${{\boldsymbol{u}}^*}\left( {{k}} \right)$ 的首个分量${{\boldsymbol{u}}^*}\left( {0|{{k}}} \right)$ 作用于系统(1).步骤 5. 令
$ k = k + 1 $ , 测量系统(1)的状态; 返回步骤2.为更清晰地描述算法1的运行过程, 图1给出了每个时刻两个优化问题的求解顺序.
由算法1和图1可知, 在初始时刻
$ k = 0 $ , 由于约束(8e)不起作用, 故在初始时刻无需求解优化问题(16), 并可令$\eta {{(}}{{{\boldsymbol{x}}} _0},\;\lambda {{)}} \to \infty$ .3. 输入到状态稳定性
闭环系统(19)具有ISS性质的前提是需要保证算法1具有递推可行性, 即两个优化问题在每个时刻都至少存在一组可行解(不一定最优), 使得对于任意容许扰动, 经济优化问题(8)和鲁棒性优化问题(16)的所有约束均满足. 注意, (8b)~(8d)与(16b)~(16d)具有相同约束.
考虑约束不确定系统(1)~(3)和状态
$\tau = $ $ {\boldsymbol{x}}\left( {0{{|k}}} \right) \in {{X}}$ . 对于优化问题(8), 定义$ N $ 步可行控制序列集为$${U_N}({\boldsymbol{\tau}} ) = \left\{ {{\boldsymbol{u}}\;\left| {\begin{array}{*{20}{l}} {\;{\boldsymbol{u}}(i) \in U,{\rm{ }}\varphi (i;{\boldsymbol{\tau}} ,{\boldsymbol{u}},{\boldsymbol{w}}) \in X} \\ {\;\varphi (N;{\boldsymbol{\tau}} ,{\boldsymbol{u}},{\boldsymbol{w}}) \in \Omega } \\ {\;{V_R}({\boldsymbol{\tau}} ,{\boldsymbol{u}},{\boldsymbol{w}}) \le \eta (\tau ,\lambda )} \\ {\;\forall i \in {{\bf{I}}_{0:N - 1}},{\rm{ }}\forall {\boldsymbol{w}} \in {W^N}} \end{array}} \right.} \right\}$$ 其中,
$ N $ 步扰动约束集$ {W^N} = W \times W \times \cdots \times W $ .定义 3. 考虑系统(1)和状态
$\tau = {\boldsymbol{x}}\left( {0{{|k}}} \right) \in {{X}}$ , 如果可行控制序列集$ {{{U}}_{{N}}}\left( \tau \right) $ 非空, 则称$ \tau $ 为该系统的初始可行状态. 所有初始可行状态$ \tau $ 的集合称为系统的初始可行状态集$ {{{X}}_{{N}}} $ .定理 1. 如果假设1~3成立, 且容许扰动满足式(15), 则优化问题(16)具有递推可行性.
证明. 考虑序列对
$( {{{\boldsymbol{u}}^1}( {{k}} ),{{\boldsymbol{w}}^1}( {{k}} )} )$ 并代入系统(1)得到状态序列${{{\boldsymbol{x}}} ^1}( {{k}} ) = {{\{ }}{{\boldsymbol{x}}^*}( {1|{{k}} - 1} ){{,}} \cdots {{,}} {{\boldsymbol{x}}^*}( {{N}}|{{k}} - 1 ){{,}} $ $ {{\boldsymbol{x}}^1}( {{{N}}|{{k}}} ){{\} }}$ , 其中${{\boldsymbol{x}}^*}( {{{N|k}} - 1} ) \in {{{\boldsymbol{\Omega}} }}$ 是对应于${{\boldsymbol{u}}^*}( {{{k}} - 1} )$ 的终端预测状态. 由假设3和注2可知,${{{\boldsymbol{\Omega}} }}$ 是闭环系统(13)的鲁棒不变集, 故对于任意${{\boldsymbol{w}}^1}( {{{N}} - 1{{|k}}} ) \in $ $ {{W}}$ ,${{\boldsymbol{u}}^1}( {{{N}} - 1{{|k}}} ) = {\boldsymbol{\pi}} ( {{{\boldsymbol{x}}^*}( {{{N|k}} - 1} )} ) \in {{U}}$ 和${{\boldsymbol{x}}^1}( {{{N|k}}} ) = $ $ \varphi ( {{{N}};{{\boldsymbol{x}}^*}( {{{N}}|{{k}} - 1} ),{{\boldsymbol{u}}^1}( {{{N}} - 1|{{k}}} ),{{\boldsymbol{w}}^1}( {{{N}} - 1|{{k}}} )} ) \in {{{\boldsymbol{\Omega}} }}$ 成立. 利用序列对$( {{{\boldsymbol{u}}^1}( {{k}} ),{{\boldsymbol{w}}^1}( {{k}} )} )$ 构造优化问题(16)在$ k $ 时刻的备选解$$\begin{split} &({{\boldsymbol{u}}^2}(k),{{\boldsymbol{w}}^2}(k)) =\\ &\quad \left\{ {\begin{array}{*{20}{l}} &({{\boldsymbol{u}}^1}(i + 1|k),{{\boldsymbol{w}}^1}(i + 1|k)),&i \in {{\bf{I}}_{0:N - 2}} \\ &({\boldsymbol{\pi}} ({{\boldsymbol{x}}^1}(N|k)),{{\boldsymbol{w}}^2}(N - 1|k),&i = N - 1 \end{array}} \right. \end{split}$$ (20) 其对应状态序列
${{{\boldsymbol{x}}} ^2}( {{k}} ) \;=\; {{\{ }}{{\boldsymbol{x}}^1}( {1|{{k}}} ){{,}} \cdots {{,}}\;\;{{\boldsymbol{x}}^1}( {{{N}}|{{k}}} ){{,}} $ $ {{\boldsymbol{x}}^2}( {{{N}}|{{k}}} ){{\} }}$ , 其中${{\boldsymbol{x}}^1}( {{{N|k}}} ) \in {{{\boldsymbol{\Omega}} }}.$ 又${{{\boldsymbol{\Omega }}}}$ 是鲁棒不变的, 则对于任意${{\boldsymbol{w}}^2}( {{{N}} - 1{{|k}}} ) \in {{W}},$ ${{\boldsymbol{u}}^2}( {{{N}} - 1{{|k}}} ) = $ $ {\boldsymbol{\pi}} ( {{{\boldsymbol{x}}^1}( {{{N|k}}} )} ) \in {{U}}$ 和${{\boldsymbol{x}}^2}( {{{N|k}}} ) = \varphi ( {{N}};{x^1}( {{{N}}|{{k}}} ),{{\boldsymbol{u}}^2}( {{N}} - $ $ 1|{{k}} ), {{\boldsymbol{w}}^2}({{N}} - 1|{{k}} ) ) \in {{{\boldsymbol{\Omega}} }}$ 成立, 因此约束(16b)~(16d)满足.再分别将
$\left( {{{\boldsymbol{u}}^1}\left( {{k}} \right),{{\boldsymbol{w}}^1}\left( {{k}} \right)} \right)$ 和$\left( {{{\boldsymbol{u}}^2}\left( {{k}} \right),{{\boldsymbol{w}}^2}\left( {{k}} \right)} \right)$ 代入函数$ {{V} _R}\left( x \right) $ 中, 计算得$$\begin{split} & V_R^2({\boldsymbol{x}}(k)) - V_R^1({\boldsymbol{x}}(k)) - {\lambda _1}{L_w}(||{\boldsymbol{w}}||) = \\ & \qquad E({{\boldsymbol{x}}^2}(N|k)) - E({{\boldsymbol{x}}^1}(N|k)) \;+ \\ & \qquad{\rm{ }}L({{\boldsymbol{x}}^1}(N|k),{\boldsymbol{\pi}} ({{\boldsymbol{x}}^1}(N|k))) - {L_w}({{\boldsymbol{w}}^1}(N - 1|k))\; - \\ & \qquad L({{\boldsymbol{x}}^1}(0|k),{{\boldsymbol{u}}^1}(0|k))+{L_w}({{\boldsymbol{w}}^1}(0|k))- {\lambda _1}{L_w}(||{\boldsymbol{w}}||) \end{split} $$ (21) 其中,
${V} _R^2\left( {{\boldsymbol{x}}\left( {k} \right)} \right) = {{V} _R}\left( {{\boldsymbol{x}}\left( {k} \right),{{{{\boldsymbol{u}}}}^2}\left( {k} \right),{{{{\boldsymbol{w}}}}^2}\left( {k} \right)} \right)$ . 再次应用假设3, 将式(10)代入式(21), 则对于任意$ {\lambda _1} \geq 1 $ , 有$$V_R^2({\boldsymbol{x}}(k)) \le V_R^1({\boldsymbol{x}}(k)) + {\lambda _1}{L_w}(||{\boldsymbol{w}}||)$$ (22) 故约束(16e)满足, 从而
$\left( {{{\boldsymbol{u}}^2}\left( {{k}} \right),{{\boldsymbol{w}}^2}\left( {{k}} \right)} \right)$ 为问题(16)在$ k $ 时刻的可行解, 即优化问题(16)具有递推可行性. □定理 2. 如果假设1~3成立, 且容许扰动满足式(15), 则优化问题(8)具有递推可行性.
证明. 令
$( {{{\boldsymbol{u}}^0}( {{k}} ),{{\boldsymbol{w}}^0}( {{k}} )} )$ 为优化问题(16)在$ k $ 时刻的最优解, 由定理1可知${V} _R^0( {{\boldsymbol{x}}( {k} )} ) \leq {V} _R^1( {{\boldsymbol{x}}( {k} )} ) + $ $ {\lambda _1}{{L} _{w} }( {\| {\boldsymbol{w}} \|} )$ . 对${V} _R^0( {{\boldsymbol{x}}( {k} )} )$ 和${V} _R^*( {{\boldsymbol{x}}( {{k} - 1} )} )$ 做如下运算:$$\begin{split} & V_R^0({\boldsymbol{x}}(k)) - V_R^ * ({\boldsymbol{x}}(k - 1)) \le \\ &\qquad{\rm{ }}V_R^1({\boldsymbol{x}}(k)) + {\lambda _1}{L_{\boldsymbol{w}}}(||w||) - V_R^ * ({\boldsymbol{x}}(k - 1)) = \\ & \qquad{\rm{ }}E({{\boldsymbol{x}}^1}(N|k)) - E({{\boldsymbol{x}}^*}(N|k - 1))\; + \\ & \qquad\,L({{\boldsymbol{x}}^*}(N|k - 1),{\boldsymbol{\pi}} ({{\boldsymbol{x}}^*}(N|k - 1))) \;- \\ & \qquad{L_w}({{\boldsymbol{w}}^1}(N - 1|k))-L({{\boldsymbol{x}}^*}(0|k-1),{{\boldsymbol{u}}^*}(0|k-1))\; + \\ & \qquad\,{L_w}({{\boldsymbol{w}}^*}(0|k - 1)) + {\lambda _1}{L_w}(||{\boldsymbol{w}}||) \\[-10pt] \end{split} $$ (23) 应用假设2和假设3, 式(23)改写为
$$\begin{split} &V_R^ * ({\boldsymbol{x}}(k - 1)) - V_R^0({\boldsymbol{x}}(k)) + (1 + {\lambda _1}){L_w}(||{\boldsymbol{w}}||) \ge \\ &\qquad L({\boldsymbol{x}}(k - 1),{\boldsymbol{u}}(k - 1)) \\[-10pt] \end{split} $$ (24) 将式(24)代入函数
$ \eta $ , 则对于任意$ \lambda \in \left[ {0,1} \right) $ , 有$$V_R^0({\boldsymbol{x}}(k)) \le \eta ({\boldsymbol{x}}(k){\rm{,}}\lambda )$$ (25) 故约束(8e)成立. 由定理1可知,
$\left( {{{\boldsymbol{u}}^0}\left( {{k}} \right),{{\boldsymbol{w}}^0}\left( {{k}} \right)} \right)$ 作为优化问题(16)在$ k $ 时刻的最优解满足约束(8b)~(8d). 因此, 在$ k $ 时刻总存在可行控制序列$ {{\boldsymbol{u}}^0}\left( {{k}} \right) $ 满足优化问题(8)的所有约束, 从而优化问题(8)具有递推可行性. □根据定理1和定理2, 本文鲁棒EMPC策略中的双目标优化问题皆具有递推可行性. 下面将基于定理1和定理2给出闭环系统(19)的ISS结果.
定理 3. 如果假设1~3成立, 且容许扰动满足式(15), 则当优化问题(8)在初始时刻存在可行解时, 闭环系统(19)在鲁棒不变集XN内相对于扰动具有ISS.
证明. 由定理1和定理2可知, 当优化问题(8)在初始时刻存在可行解时, 该优化问题在任意k时刻都是可行的, 则对任意
${\boldsymbol{x}}\left( {{k}} \right) \in {{{X}}_{{N}}}$ 和${\boldsymbol{w}}\left( {{k}} \right) \in {{W}},$ 闭环系统(19)满足${\boldsymbol{x}}\left( {{{k + }}1} \right) \in {{{X}}_{{N}}}.$ 由定义3可知,$ {{{X}}_{{N}}} $ 为闭环系统(19)的一个鲁棒不变集.令
${{V} _{N} }\left( {{\boldsymbol{x}}\left( {k} \right)} \right) = {V} _R^*\left( {{\boldsymbol{x}}\left( {k} \right)} \right) + {N} {{L} _{w} }\left( {\left\| {\boldsymbol{w}} \right\|} \right)$ , 并选择${{V} _{N} }\left( {{\boldsymbol{x}}\left( {k} \right)} \right)$ 作为闭环系统(19)的备选ISS-Lyapunov函数. 考虑假设2, 对于任意${\boldsymbol{x}}\left( {{k}} \right) \in {{{X}}_{{N}}},$ 有$$\begin{split} &{V_N}({\boldsymbol{x}}(k)) = V_R^ * ({\boldsymbol{x}}(k)) + N{L_w}(||{\boldsymbol{w}}||) \ge \\ &\qquad L({\boldsymbol{x}}(k),{\boldsymbol{u}}(k)) \ge{\alpha _1}(|{\boldsymbol{x}}(k)|) \end{split} $$ (26) 其中,
$ {\alpha _1} = {\alpha _l} $ 为$ {K_\infty } $ 类函数.假设函数
${V} _R^{}\left( {{\boldsymbol{x}}\left( {k} \right)} \right)$ 存在上界${{v}} < + \infty ,$ 且$ {{v}} $ 为一个足够大的常数. 定义内部包含原点且半径为$ {{r}} $ 的球形区域${{{B}}_{{r}}} = \{ {\boldsymbol{x}} \in {{\mathbf{R}}^{{n}}}:\left| {\boldsymbol{x}} \right| \leq {{r}}\} \subseteq {\boldsymbol{\Omega}}.$ 设常数$ \varepsilon = \max \{ 1,{{v}}/{\beta _{{f}}}\left( {{r}} \right)\} \geq 1, $ 令函数${\alpha _2} = \varepsilon {\beta _f}.$ 下面考虑两种情况:1) 当
${\boldsymbol{x}}\left( {{k}} \right) \in {{{\boldsymbol{\Omega}} }}$ 时, 由假设3可知,$$\begin{split} &E({\boldsymbol{x}}(i + 1|k)) - E({\boldsymbol{x}}(i|k)) \le \\ &\qquad {L_w}({\boldsymbol{w}}(i|k)) - L({\boldsymbol{x}}(i|k),{\boldsymbol{\pi}} ({\boldsymbol{x}}(i|k))),\;\;{\rm{ }}i \in {{\bf{I}}_{0:N - 1}} \end{split} $$ (27) 其中,
${\boldsymbol{x}}( i + 1{{|k}} ) = \varphi ( {i + 1;{\boldsymbol{x}}( {{{i}}|{{k}}} ),{\boldsymbol{\pi}} ( {{\boldsymbol{x}}( {{{i}}|{{k}}} )} ),{\boldsymbol{w}}( {{{i}}|{{k}}} )} )$ 且${\boldsymbol{x}}( {{{i|k}}} ) \in {{{\boldsymbol{\Omega}} }}$ . 由于${{{\boldsymbol{\Omega}} }}$ 是鲁棒不变集, 故对于任意容许扰动${\boldsymbol{w}}( {{{i|k}}} )$ , 满足${\boldsymbol{\pi}} ( {{\boldsymbol{x}}( {{{i|k}}} )} ) \in {{U}}$ 和${\boldsymbol{x}}( {i + 1{{|k}}} ) \in {{{\boldsymbol{\Omega}} }}$ . 从${{i}} = 0,$ 到$ {{N}} - 1 $ , 将不等式(27)累加得$$\begin{split} E({\boldsymbol{x}}(k)) \le \;&\sum\limits_{i = 0}^{N - 1} {(L({\boldsymbol{x}}(i|k),{\boldsymbol{u}}(i|k)) - {L_w}({\boldsymbol{w}}(i|k)))} {\rm{ + }} \\ &E({\boldsymbol{x}}(N|k))\\[-10pt] \end{split} $$ (28) 结合假设2, 存在
$ {K_\infty } $ 类函数$ {\beta _f} $ 和$ {\alpha _2} $ 满足$$V_R^ * ({\boldsymbol{x}}(k)) \le E({\boldsymbol{x}}(k)) \le {\beta _f}(|{\boldsymbol{x}}(k)|) \le {\alpha _2}(|{\boldsymbol{x}}(k)|)$$ (29) 2) 当
${\boldsymbol{x}}\left( {{k}} \right) \notin {{{\boldsymbol{\Omega}} }}$ 时,$ {\boldsymbol{x}}\left( {{k}} \right) \notin {{{B}}_{{r}}} $ , 则${\beta _{{f}}}\left( {|{\boldsymbol{x}}\left( {{k}} \right)|} \right) \geq $ $ {\beta _{{f}}}\left( {{r}} \right)$ , 得$$\begin{split} V_R^ * ({\boldsymbol{x}}(k)) \le\;& v \le v \times \dfrac{{{\beta _f}(|{\boldsymbol{x}}(k)|)}}{{{\beta _f}(r)}} \le \\ &\varepsilon {\beta _f}(|{\boldsymbol{x}}(k)|) = {\alpha _2}(|{\boldsymbol{x}}(k)|) \end{split} $$ (30) 联立上述两种情况, 整理得
$$V_R^ * ({\boldsymbol{x}}(k)) \le {\alpha _2}(|{\boldsymbol{x}}(k)|),\;\;{\rm{ }}\forall {\boldsymbol{x}}(k) \in {X_N}$$ (31) 则对于任意容许扰动
${\boldsymbol{w}}\left( {{k}} \right)$ , 值函数${{V} _{N} }\left( {{\boldsymbol{x}}\left( {k} \right)} \right)$ 满足$${V_N}({\boldsymbol{x}}(k)) \le {\alpha _2}(|{\boldsymbol{x}}(k)|) + {\delta _1}(||{\boldsymbol{w}}||),\;\;{\rm{ }}\forall {\boldsymbol{x}}(k) \in {X_N}$$ (32) 其中,
$ {\delta _1} $ 为$ K $ 类函数.令
$\left( {{{\boldsymbol{u}}^*}\left( {{{k}} - 1} \right),{{\boldsymbol{w}}^*}\left( {{{k}} - 1} \right)} \right)$ 和$\left( {{{\boldsymbol{u}}^*}\left( {{k}} \right),{{\boldsymbol{w}}^*}\left( {{k}} \right)} \right)$ 分别为优化问题(8)在$ {{k}} - 1 $ 和$ {{k}} $ 时刻的经济最优解, 考虑函数(18), 将函数${V} _R^*\left({\boldsymbol{ x}} \right)$ 沿着闭环状态轨迹(19)做差分运算$$\begin{split} &V_R^ * ({\boldsymbol{x}}(k)) - V_R^ * ({\boldsymbol{x}}(k - 1)) \le \\ &\qquad \eta ({\boldsymbol{x}}(k),\lambda ) - V_R^ * ({\boldsymbol{x}}(k - 1)) = \\ & \qquad (1 - \lambda )(V_R^0({\boldsymbol{x}}(k)) - V_R^ * ({\boldsymbol{x}}(k - 1))) \;+ \\ &\qquad \lambda (1 + {\lambda _1}){L_w}(||{\boldsymbol{w}}||) \end{split} $$ (33) 由定理1可知,
${V} _R^0( {{\boldsymbol{x}}( {k} )} ) \leq {V} _R^1( {{\boldsymbol{x}}( {k} )} ) + {\lambda _1}{{L} _{w} }( {\| {\boldsymbol{w}} \|} )$ , 将该不等式代入式(33), 整理可得$$\begin{split} & V_R^ * ({\boldsymbol{x}}(k)) - V_R^ * ({\boldsymbol{x}}(k - 1)) \le \\ &\qquad (1 - \lambda )(V_R^1({\boldsymbol{x}}(k)) - V_R^ * ({\boldsymbol{x}}(k - 1))\; + \\ & \qquad {\lambda _1}{L_w}(||{\boldsymbol{w}}||)) + \lambda (1 + {\lambda _1}){L_w}(||{\boldsymbol{w}}||) \le \\ & \qquad (\lambda -1)L({\boldsymbol{x}}(k-1),{\boldsymbol{u}}(k-1)) + (1+ {\lambda _1}){L_w}(||{\boldsymbol{w}}||) \end{split} $$ (34) 考虑函数
${{V} _{N} }\left( {{\boldsymbol{x}}\left( {k} \right)} \right)$ 定义, 式(34)等价于$$\begin{split} &{V_N}({\boldsymbol{x}}(k)) - {V_N}({\boldsymbol{x}}(k - 1)) = \\ &\qquad V_R^ * ({\boldsymbol{x}}(k)) - V_R^ * ({\boldsymbol{x}}(k - 1)) \le \\ &\qquad (\lambda -1)L({\boldsymbol{x}}(k-1),{\boldsymbol{u}}(k-1)) + (1+{\lambda _1}){L_w}(||{\boldsymbol{w}}||) \end{split} $$ (35) 由于
$ \lambda \in \left[ {0,1} \right) $ 和$ {\lambda _1} \geq 1 $ , 故存在$ {K_\infty } $ 类函数$ {\alpha _3} $ 和$ K $ 类函数$ {\delta _2} $ , 使得如下不等式$$\begin{split} &{V_N}({\boldsymbol{x}}(k)) - {V_N}({\boldsymbol{x}}(k - 1)) \le \\ &\qquad {\delta _2}(||{\boldsymbol{w}}||) - {\alpha _3}(|{\boldsymbol{x}}(k - 1)|),\;\;{\rm{ }}\forall {\boldsymbol{x}}(k - 1) \in {X_N} \end{split} $$ (36) 成立. 则联立不等式(26), (32)和(36), 由引理1可知, 值函数VN(x(k))为闭环系统(19)的一个ISS-Lyapunov函数. 因此, 闭环系统(19)在XN内相对于扰动具有输入到状态稳定性. □
4. 实例仿真
考虑不确定非线性连续搅拌釜反应器(CSTR)
$$\left\{ {\begin{aligned} &{\dfrac{{{\rm{d}}{c_A}(t)}}{{{\rm{d}}t}} = \dfrac{{Q(t)({c_{Af}} - {c_A}(t) + \Delta c(t))}}{V} - {k_0}c_A^{\rm{2}}(t)} \\ &{\dfrac{{{\rm{d}}{c_B}(t)}}{{{\rm{d}}t}} = \dfrac{{Q(t)({c_{Bf}} - {c_B}(t))}}{V} + {k_0}c_A^{\rm{2}}(t)} \end{aligned}} \right.$$ (37) 其中,
$ {{{c}}_{A} } $ 和$ {{{c}}_{B} } $ 分别为组分A和B的浓度,$ {{Q}} $ 为反应器进料流量,${{{c}}_{{A} {{\rm{f}}} }}$ 和${{{c}}_{{B} {{\rm{f}}} }}$ 分别为进料中组分A和B的浓度,$ \Delta c $ 为进料中组分A浓度的不确定波动, 体积$ {{V}} $ 和反应动力学参数$ {{{k}}_0} $ . 当忽略进料中组分浓度的波动, 该模型广泛用于名义稳定EMPC综合策略的验证[10, 13, 21]. 这里假设进料中组分A的浓度波动是有界的, 用于验证本文鲁棒EMPC的有效性. 取模型参数[10]:${{{c}}_{{A} {{\rm{f}}} }}$ = 1.0 mol/l,${{{c}}_{{B} {{\rm{f}}} }}$ = 0,$ {{V}} $ = 10 l和$ {{{k}}_0} $ = 1.2 l/(mol⋅min).令
${\left[ {{{{c}}_A},{{{c}}_B}} \right]^{\rm{T}}}$ 为系统状态${\left[ {{{{x}}_1},{{{x}}_2}} \right]^{\rm{T}}}$ ,$ {{Q}} $ 为控制输入$ {{u}} $ ,$ {{w}} $ 为扰动$ \Delta c $ . 进一步, 定义状态约束和控制约束$${x_1} \in [0,1],{\rm{ }}\;{x_2} \in [0,1],\;{\rm{ }}u \in [0,15]$$ (38) 选择采样周期
$ {{{T}}_{{s}}} = 0.1 $ min, 利用欧拉差分法离散化连续时间模型(37), 得CSTR非线性离散时间模型$$\left\{ {\begin{aligned} & \left[ {\begin{array}{*{20}{c}} {{x_1}(k + 1)} \\ {{x_2}(k + 1)} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {{x_1}(k) - 1.2{T_s}{x_1}^2(k)} \\ {{x_2}(k) + 1.2{T_s}{x_1}^2(k)} \end{array}} \right] + \\ &\qquad \dfrac{{{T_s}}}{{10}}\left[ {\begin{array}{*{20}{c}} {1 - {x_1}(k)} \\ { - {x_2}(k)} \end{array}} \right]u(k) + \dfrac{{{T_s}}}{{10}}\left[ {\begin{array}{*{20}{c}} {u(k)} \\ 0 \end{array}} \right]w(k) \\ & {z(k) = {{\left[ {\begin{array}{*{20}{c}} {{\boldsymbol{x}}^{\rm{T}}{{(k)}}}&{u(k)} \end{array}} \right]}^{\rm{T}}}} \end{aligned}} \right.$$ (39) 其中, k为采样时刻. 定义CSTR经济性能指标[10]
$${L_e}({\boldsymbol{x}},u) = 0.5u(1 - 4{x_2})$$ (40) 由式(4)计算最优经济平衡点为
${\boldsymbol{x}}_{{e}}^* = ( 0.5785, $ $ 0.4215)$ 和${{{{u}}}}_{{e}}^* = 9.5258$ . 通过求解关于平衡点$\left( {{\boldsymbol{x}}_{{e}}^*,{{{{u}}}}_{{e}}^*} \right)$ 的Riccati不等式方程(11), 得$$\begin{split} &{\boldsymbol{P}} = \left[ {\begin{array}{*{20}{l}} {20}&{0.5} \\ {0.5}&{20} \end{array}} \right]\\ & {{\boldsymbol{\pi}} {\rm{(}}{\boldsymbol{x}}{\rm{)}} = [ - 0.0655,\;0.0739]({\boldsymbol{x}} - {\boldsymbol{x}}_e^*) + u_e^*} \\ &{\boldsymbol{\Omega}} =\{ {\boldsymbol{x}} \in {{\bf{R}}^2}{\rm{:}}\;{{{\rm{(}}{\boldsymbol{x}} - {\boldsymbol{x}}_e^*)}^{\rm{T}}}{\boldsymbol{P}}{\rm{(}}{\boldsymbol{x}} - {\boldsymbol{x}}_e^*) \le {\rm{2}}{\rm{.12\} }} \end{split}$$ 进一步, 为保证终端状态集
$\boldsymbol{\Omega} $ 的鲁棒不变性, 由式(15)计算容许扰动上界为0.1436.令预测步长
$ N = 5 $ , 仿真总步长${{{T}}_{{{{\rm{sim}}}}}} = 100$ , 系数$ {\lambda _1} = 1 $ , 并考虑不同收缩因子$ \lambda \in \left[ {0,1} \right) $ 的控制效果. 选择系统初始状态${{\boldsymbol{x}}_0} = \left( {0.15,0.7} \right)^{\rm{T}}$ , 仿真结果如图2和图3所示, 图2对应持续扰动${{w}}( {{k}})= 0.1436\sin ( {{{k/}}2} )$ 的仿真结果, 图3对应衰减扰动${{w}}\left( {{k}} \right)=0.1436 $ $ \exp( - {{k/}}10)$ 的仿真结果. 分析图2可知, 在持续扰动下, 闭环系统最终在最优经济平衡点$ \left( {{\boldsymbol{x}}_{{e}}^*,{{{{u}}}}_{{e}}^*} \right) $ 附近有界稳定; 分析图3可知, 在衰减扰动下, 闭环状态轨迹渐近收敛于$\left( {{\boldsymbol{x}}_{{e}}^*,{{{{u}}}}_{{e}}^*} \right)$ . 根据图2和图3可知, 对于不同收缩因子$ \lambda \in \left[ {0,1} \right) $ , 由算法1得到的鲁棒EMPC控制器使得闭环系统在$\left( {{\boldsymbol{x}}_{{e}}^*,{{{{u}}}}_{{e}}^*} \right)$ 相对于容许扰动总是ISS的, 但不同$\lambda $ 对应的闭环状态轨迹的动态响应不同. 以衰减扰动${{w}}\left( {{k}} \right) = 0.1436\exp( - {{k/}}10)$ 为例, 由图3可知, 因子$\lambda $ 越小, 闭环状态轨迹的收敛过渡时间$ {{{T}}_{{{\rm{{tr}}}}}} $ 越短, 具体数据如表1最右列所示.
图 2 持续扰动$w(k)=0.1436\sin(k/2) $ 的仿真结果Fig. 2 Simulation results under continuous disturbance$w(k)=0.1436\sin(k/2) $ 
图 3 持续扰动$w(k)=0.1436\exp(-k/10) $ 的仿真结果Fig. 3 Simulation results under continuous disturbance$w(k)=0.1436\exp(-k/10) $ 表 1 平均经济性能和收敛过渡时间Table 1 Average economic performance and transient time$\lambda $ $w(k)= 0.1436\sin ( {k/} 2 )$ $w(k)=0.1436\exp( - {{k/} }10)$ ${J_{{\rm{ave}}} } $ ${J_{{\rm{ave}}} } $ ${ {{T} }_{ {\text{tr} } } }$ 0.1 −3.4361 −3.3953 49Ts 0.3 −3.4464 −3.4019 54Ts 0.5 −3.4553 −3.4083 58Ts 0.7 −3.4623 −3.4143 65Ts 0.9 −3.4712 −3.4194 74Ts 定义闭环系统T步平均经济性能指标
$${J_{{\rm{ave}}}} = \frac{{\sum\limits_{k = 0}^T {{L_e}({\boldsymbol{x}}(k),{{\boldsymbol{u}}^{{\rm{mpc}}}}({\boldsymbol{x}}(k)))} }}{T}$$ (41) 取
${T} = 100 ,$ 表1给出了不同$ \lambda $ 值对应的闭环系统平均经济性能. 分析表1可以看出, 因子$ \lambda $ 越大, 闭环系统平均经济性能越好. 这表明, 无论是名义系统还是不确定系统, 闭环系统的经济最优性和稳定性是互相冲突的双控制目标[10, 21]. 在本文鲁棒EMPC策略中, 可以通过调节收缩因子$ \lambda $ 对经济性控制目标和鲁棒稳定性控制目标进行权衡, 从而实现经济性和鲁棒稳定性综合控制效果.进一步, 为了验证本文鲁棒EMPC策略(简记为MM-EMPC)的优越性, 对比文献[38]中鲁棒收缩EMPC策略(简记为RC-EMPC), 研究不同容许扰动上界对系统稳定性的影响. 对约束CSTR系统(37)和(38), 在系统参数相同情况下, RC-EMPC策略容许扰动上界d1 = 0.0136, MM-EMPC容许扰动上界d2 = 0.1436. 选择持续扰动w(k) = 0.04sin(k/5), 令初始状态x0 = (0.5,
$0.5)^{\rm{T}} $ , 仿真步长Tsim = 500. 在相同仿真环境下, 两种鲁棒EMPC策略的控制结果如图4所示, 其中, 实线表示MM-EMPC策略仿真结果, 虚线表示RC-EMPC策略仿真结果.从图4可以看出, MM-EMPC策略对应的闭环状态轨迹和控制输入曲线均收敛于
$( {{\boldsymbol{x}}_{{e}}^*,{{{\boldsymbol{u}}}}_{{e}}^*})$ 附近, 而RC-EMPC策略对应的闭环轨迹虽然收敛, 但存在较大的稳态误差, 且MM-EMPC闭环系统的动态响应更加快速. 进一步, 考虑衰减扰动${{w}}\left( {{k}} \right)= $ $ 0.04\exp \left( { - {{k/}}5} \right)$ , 选择4个不同的初始状态对两种控制策略进行仿真对比分析, 控制结果如图5所示. 综合图4和图5分析可知, 当容许扰动上界$ {{d}} \in \left( {{{{d}}_1},{{{d}}_2}} \right) $ 时, MM-EMPC策略下的闭环系统相对于扰动是ISS的, 而RC-EMPC却无法保证闭环系统相对于扰动具有ISS, 其中一个重要原因是RC-EMPC对大扰动会丢失递推可行性. 因此, 相较于RC-EMPC策略, 本文策略在保证闭环系统ISS同时, 能够获得更大的容许扰动上界, 从而降低鲁棒EMPC控制器的保守性.
图 5 不同初始状态的闭环系统相轨迹Fig. 5 Phase trajectories of the closed-loop system from different initial states5. 结束语
本文针对有界扰动下的约束不确定仿射输入非线性系统, 提出了一种新的鲁棒EMPC策略. 基于微分对策原理分别对经济目标函数和关于最优经济平衡点的鲁棒稳定性目标函数进行优化, 利用得到的鲁棒稳定性目标最优值函数构造隐式收缩约束, 保证了双控制目标优化问题的递推可行性, 并建立了闭环系统在最优经济平衡点处相对于扰动的输入到状态稳定性结果. 通过对不确定CSTR经济优化控制的对比仿真实验, 验证了本文策略的有效性和优越性.
尽管基于微分对策的min-max鲁棒EMPC在理论上能有效提高鲁棒EMPC的性能, 但min-max优化在线计算复杂, 将阻碍鲁棒EMPC在快速响应系统中的应用. 因此, 降低min-max鲁棒EMPC的在线优化计算量和设计更高效的鲁棒EMPC策略(如Tube鲁棒EMPC、本质鲁棒EMPC等)将是后续研究重点.
-
图 1 基于GBDT的铁路事故类型预测及成因分析框架
Fig. 1 The framework of GBDT-based railroad accident type prediction and cause analysis
图 6 两类事故致因中不同因素的比例
Fig. 6 Proportion of different factors in causes of two types of railroad accidents
表 1 原始数据描述
Table 1 Description of original data
Record Accident type Attribute Number 5 434 11 144 
下载: 导出CSV
表 2 事故类型描述
Table 2 Description of accident types
Type Description 1 Derailment 2 Head on collision 3 Rearend collision 4 Side collision 5 Raking collision 6 Broken train collision 7 Hwy-rail crossing 8 RR grade crossing 9 Obstruction 10 Fire 11 Other impacts
下载: 导出CSV
表 3 数据集部分示例
Table 3 Examples of the dataset
Name Description Number Type RAILROAD Railroad code 5 434 Object CARS Num. of cars carrying hazmat 5 434 Int64 TYPSPD Train speed type 5 086 Object TRNDIR Train direction 5 161 Float64 TONS Gross tonnage, excluding power units 5 434 Int64 TYPEQ Type of consist 5 081 Object EQATT Equipment attended 5 074 Object CDTRHR Num. of hours conductors on duty 3 628 Int64 ENGHR Num. of hours engineers on duty 4 201 Int64 TRKNAME Track identification 5 434 Object
下载: 导出CSV
表 4 预处理后数据描述
Table 4 Description of preprocessed data
Record Accident type Attribute dimension Number 5 434 11 119
下载: 导出CSV
表 5 三种方法补全前后特征TRNDIR取值分布
Table 5 Distribution of the attribute TRNDIR values before and after three completion methods
Algorithm $a_j=1$ $a_j=2$ $a_j=3$ $a_j=4$ Before completion 0.22 0.20 0.31 0.27 Interpolation 0.21 0.19 0.30 0.30 Mode 0.21 0.19 0.34 0.26 Our algorithm 0.22 0.20 0.31 0.27
下载: 导出CSV
表 6 不同采样率下集成GBDT分类准确率
Table 6 Accuracy of classifiers with different sampling rates
$\alpha$ 0.6 0.7 0.8 0.9 1.0 Accuracy 0.841 0.846 0.845 0.852 0.848
下载: 导出CSV
表 7 各分类器性能对比
Table 7 Performance comparison of classifiers
Classifier Accuracy Precision Recall F1 DT 0.728 0.73 0.73 0.73 RF 0.773 0.74 0.77 0.75 ET 0.734 0.70 0.73 0.71 GBDT 0.841 0.84 0.84 0.84 Ensemble GBDT 0.852 0.85 0.85 0.85
下载: 导出CSV
表 8 重要度排名前15的特征
Table 8 Features of top 15 in importance
No. Name Description 1 Latitude Latitude in decimal degrees 2 Longitude Longitude in decimal degrees 3 CNTYCD FIPS county code 4 HIGHSPD Maximum speed 5 TRKNAME Track identification 6 RRCAR1 Car initials (fist involved) 7 TEMP Temperature in degrees fahrenheit 8 MILEPOST Milepost 9 STATION Nearest city and town 10 TRNSPD Speed of train in miles per hour 11 RRCAR2 Car initials (causing) 12 SUBDIV Railroad subdivision 13 ENGHR Num. of hours engineers on duty 14 CDTRHR Num. of hours conductors on duty 15 TONS Gross tonnage
下载: 导出CSV
-
[1] Ming L. Data mining: concepts, models, methods, and algorithms. IIE Transaction, 2004, 36(5): 495-496 [2] 冯士雍. 回归分析方法. 北京: 科学出版社, 1974Feng Shi-Yong. Regression Analysis Method. Beijing: Science Press, 1974 [3] Rutkowski L, Jaworski M, Pietruczuk L, Duda P. Decision trees for mining data streams based on the gaussian approximation. IEEE Transactions on Knowledge and Data Engineering, 2014, 26(1): 108−119 doi: 10.1109/TKDE.2013.34 [4] 李定启, 程远平, 王海峰, 王亮, 周红星, 孙建华. 基于决策树ID3改进算法的煤与瓦斯突出预测. 煤炭学报, 2011, 36(4): 619−622Li Ding-Qi, Cheng Yuan-Ping, Wang Hai-Feng, Wang Liang, Zhou Hong-Xing, Sun Jian-Hua. . Coal and gas outburst prediction based on improved decision tree ID3 algorithm. Journal of China Coal Society, 2011, 36(4): 619−622 [5] Breiman L. Random forest. Machine Learning, 2001, 45(1): 5−32 doi: 10.1023/A:1010933404324 [6] Friedman J H. Greedy function approximation: a gradient boosting machine. The Annals of Statistics, 2001, 29(5): 1189−1232 [7] Friedman J H. Stochastic gradient boosting. Computational Statistics and Data Analysis, 2002, 38(4): 367−378 doi: 10.1016/S0167-9473(01)00065-2 [8] 周志华. 机器学习. 北京: 清华大学出版社, 2016.Zhou Zhi-Hua. Machine Learning. Beijing: Tsinghua University Press, 2016. [9] Schonlau M. Boosted regression (boosting): an introductory tutorial and a stata plugin. The Stata Journal, 2005, 5(3): 330−354. doi: 10.1177/1536867X0500500304 [10] 翁小雄, 吕攀龙. 基于 GBDT 算法的地铁 IC 卡通勤人群识别. 重庆交通大学学报 (自然科学版), 2019, 38(5): 8−12Weng Xiao-Xiong, Lv Pan-Long. Subway IC card commuter crowd identification based on GBDT algorithm. Journal of Chongqing Jiaotong University(Natural Science), 2019, 38(5): 8−12 [11] Mursalin M, Zhang Yuan, Chen Yue-Hui, Chawla N V. Automated epileptic seizure detection using improved correlation-based feature selection with random forest classifier. Neurocomputing, 2017, 241: 204−214 doi: 10.1016/j.neucom.2017.02.053 [12] Cheng J, Li G, Chen X H. Research on travel time prediction model of freeway based on gradient boosting decision tree. IEEE Access, 2018, 7: 7466−7480 [13] Ma X, Ding C, Luan S, Wang Y, Wang Y P. Prioritizing influential factors for freeway incident clearance time prediction using the gradient boosting decision trees method. IEEE Transactions on Intelligent Transportation Systems, 2017, 18(9): 2303−2310 doi: 10.1109/TITS.2016.2635719 [14] Su H W, Zhang W J, Li Z H. Analysis and prediction of water traffic accidents in jingtang port based on improved GM(1, 1) model. In: Proceedings of the 37th Chinese Control Conference (CCC). New York, USA: IEEE, 2018.2212−2217 [15] Das S, Sun X D. Investigating the pattern of traffic crashes under rainy weather by association rules in data mining. In: Proceedings of the 93rd Transportation Research Board (TRB) Annual Meeting. Washington D.C., USA: Nation Academy of Sciences, 2014 [16] 金勇进. 缺失数据的统计处理, 北京: 中国统计出版社, 2009.Jin Yong-Jin. Statistical Processing of Missing Data. Beijing: China Statistics Press, 2009. [17] 金勇进. 调查中的数据缺失及处理 (I)-缺失数据及其影响. 数理统计与管理, 2001, 20(4): 58−60 doi: 10.3969/j.issn.1002-1566.2001.04.012Jin Yong-Jin. Data loss and processing in survey(I)) data missing and impact. Journal of Applied Statistics and Management, 2001, 20(4): 58−60 doi: 10.3969/j.issn.1002-1566.2001.04.012 [18] Collell G, Prelec D, Patil K R. A simple plug-in bagging ensemble based on threshold-moving for classifying binary and multiclass imbalanced data. Neurocomputing, 2018, 275: 330−340 doi: 10.1016/j.neucom.2017.08.035 [19] Galar M, Fernandez A, Barrenechea E, Bustince H, Herrera F. A review on ensembles for the class imbalance problem: bagging-, boosting-, and hybrid-based approaches. IEEE Transactions on Systems, Man and Cybernetics, Part C (Applications and Reviews), 2012, 42(4): 463−484 doi: 10.1109/TSMCC.2011.2161285 [20] 朱振峰, 汤静远, 常冬霞, 赵耀. 基于 GBDT 的商品分配层次化预测模型. 北京交通大学学报, 2018, 42(2): 9−13+45 doi: 10.11860/j.issn.1673-0291.2018.02.002Zhu Zhen-Feng, Tang Jing-Yuan, Chang Dong-Xia, Zhao Yao. GBDT based hierarchical model for commodity distribution prediction. Journal of Beijing Jiaotong University, 2018, 42(2): 9−13+45. doi: 10.11860/j.issn.1673-0291.2018.02.002 [21] 杨连报, 李平, 薛蕊, 马小宁, 吴艳华, 邹丹. 基于不平衡文本数据挖掘的铁路信号设备故障智能分类. 铁道学报, 2018, 40(2): 59−66 doi: 10.3969/j.issn.1001-8360.2018.02.009Yang Lian-Bao, Li Ping, Xue Rui, Ma Xiao-Ning, Wu YanHua, Zou Dan. Intelligent classification of faults of railway signal equipment based on imbalancd text data mining. Journal of the China Railway Society, 2018, 40(2): 59−66 doi: 10.3969/j.issn.1001-8360.2018.02.009 [22] Federal Railroad Administration Office of Safety Analysis [Online], available: https://safetydata.fra.dot.gov/OfficeofSafety/Default.aspx, June 1, 2019 期刊类型引用(6)
1. 任凯龙，薛斌强. 优先级分布式模型预测控制的多智能体系统. 电子设计工程. 2024(01): 82-85+90 . 
百度学术2. 闫茂德，金鑫，杨盼盼，刘军浩. 基于MPC的智能网联卡车队列最优能耗控制. 计算机仿真. 2024(03): 128-133 . 百度学术3. 曹强，党晓圆，高俊，南亮生，吴胜磊，曾庆欣，张杰. 模型预测控制理论的研究与应用. 自动化应用. 2024(09): 8-14+18 . 百度学术4. 巢瑞云，刘源. 基于PSO优化双子支持向量机的电商经济预测研究. 贵阳学院学报(自然科学版). 2024(02): 11-15+42 . 百度学术5. 马小涵，刘晓华，高荣. 基于经济模型预测控制的证券投资组合策略. 鲁东大学学报(自然科学版). 2023(02): 153-158+164 . 百度学术6. 穆建彬，冯阳辉，何德峰. 有界扰动下异质车辆队列分布式鲁棒经济预测控制. 控制与决策. 2023(05): 1386-1394 . 百度学术其他类型引用(12)
-
下载:
计量
- 文章访问数: 2019
- HTML全文浏览量: 855
- PDF下载量: 251
- 被引次数: 18
