CN104950667A - 一种应用在列车主动悬挂系统上的多速率预测控制方法 - Google Patents

Abstract

一种应用在列车主动悬挂系统上的多速率预测控制方法，本发明在高速列车主动悬挂系统模型的基础上，开展基于一种标准正交基函数-Laguerre函数的预测控制算法研究，将悬挂系统的状态空间模型作为预测控制模型，结合预测控制中的滚动优化等基本原理，考虑到列车在运行时受到的不确定轨道干扰激励等影响因素，为高速列车主动悬挂多速率及含时滞的多速率系统设计相应的预测控制器，根据系统模型，不断调整优化控制器参数，使系统达到最佳的控制效果，以优化高速列车在运行时的振动主动控制性能，保证列车运行的操作稳定性和舒适性。

Description

一种应用在列车主动悬挂系统上的多速率预测控制方法

技术领域

本发明涉及一种多速率预测控制方法，尤其涉及一种应用在列车主动悬挂系统上的多速率预测控制方法。

背景技术

目前随着列车的速度不断提升，列车车体的振动量也随之增加，这不但影响了列车行驶平稳性，而且影响了乘坐舒适性，列车的主动悬挂包括四大子系统：外界能源输入系统、传感系统、作动系统以及测量反馈调节系统，其是一个典型的复杂多环节和多变量控制系统，目前的高速列车在实际运行时存在信号采集速率不一致的问题，所以高速列车主动悬挂振动控制系统实质上是一个离散振动控制系统，在复杂的高速列车主动悬挂多变量计算机振动控制系统设计中，保持系统各处采用单一的采样策略是比较理想状态，但是也不利于振动精度的控制，多速率数字控制系统是一个周期时变系统，可以实现许多单速率数字控制系统所不具备的或难以实现的控制功能，如强镇定、同时镇定、分散控制和改善系统鲁棒性等。

近些年来基于状态空间模型的预测控制算法的研究已日趋成熟，针对连续时域系统设计的预测控制器早在上个世纪70年代就已经成功应用于各项工程控制之中，但是其是使用标准正交基函数，特别是Laguerre函数和Kautz函数等的预测控制算法在控制领域中仍是一中较新的方法，其实从模型控制轨迹的角度很容易将系统辨识和预测控制相联系起来，因为预测控制问题实质上就是将模型的未来控制轨迹公式化，这样便可以使用同样的结构来设计连续时域系统和离散系统的预测控制器，当设定的控制时域足够大时，模型也可以和经典的线性二次调节器(Linear Quadratic Regulator，LQR)相连接，而预测控制和LQR的最大区别就是预测控制是使用滚动优化窗口，即反复在线优化来获取最优性能指标，LQR则是全局使用同一优化性能指标，显然滚动优化能获得实时的最优性能指标，从而达到对系统性能更优质的控制效果，众所周知，PID控制因为其带有的积分作用常被应用在工业控制系统中，其实在传统的预测控制系统中，例如GPC和DMC也都具备这样的积分功能，因此为了在MPC模型预测控制算法中加入积分函数，就需要在建立的状态空间模型中加入积分器，和GPC、DMC算法类似，在MPC算法中最优控制轨迹既不是连续时域系统中控制信号的衍生量也不是离散系统中控制信号的增量，这种算法不需要很多的系统信息，对于系统的内部机理及模型结构也没有过多要求，只需要简单的程序就可以实现得到系统模型的状态反馈控制增益，在确定目标函数及目标函数最小性能指标后，不断调整控制器参数直到选取选取合适的预测时域Np和控制时域Nc，设未来输出的参考轨迹包含在预测时域之内，那么根据控制变量u(k)和状态变量x(k)的关系，即可计算出系统的状态反馈增益矩阵K_mpc。

发明内容

本发明的目的在于提供一种应用在列车主动悬挂系统上的多速率预测控制方法，解决了目前的高速列车在实际运行时主动悬挂系统的传感器与作动器频率特性的差异，从而导致主动悬挂系统的控制性能不佳的问题。

本发明是这样实现的，其技术方案为：

1、首先参照图1以高速列车车辆1/4车体为研究对象，建立包含车体以及前后转向架的侧滚、沉浮和点头等九自由度的横向-垂向耦合的列车主动悬挂系统动力学模型；

2、根据步骤1建立的动力学模型得到悬挂系统动力学方程：取状态变量为：

$X=\begin{array}{c}\left[\begin{array}{ccccccccccccc}{Z}_b& \mathrm{\φ}& \mathrm{\θ}& Z_2& {\mathrm{\φ}}_2& {\mathrm{\θ}}_1& Z_1& {\mathrm{\φ}}_1& {\mathrm{\θ}}_2& {\stackrel{.}{Z}}_b& \stackrel{.}{\mathrm{\φ}}& \stackrel{.}{\mathrm{\θ}}& {\stackrel{.}{Z}}_2\end{array}\right.\\ {\left.\begin{array}{ccccc}{\stackrel{.}{\mathrm{\φ}}}_2& {\stackrel{.}{\mathrm{\θ}}}_1& {\stackrel{.}{Z}}_1& {\stackrel{.}{\mathrm{\φ}}}_1& {\stackrel{.}{\mathrm{\θ}}}_2\end{array}\right]}^T\end{array}$

控制变量为：u_c＝[u₁ u′₁ u₂ u′₂ u₃ u′₃ u₄ u′₄ u₅ u′₅ u₆ u′₆]^T输出变量为： $y={\left[\begin{array}{ccc}{\stackrel{..}{Z}}_b& \stackrel{..}{\mathrm{\φ}}& \stackrel{..}{\mathrm{\θ}}\end{array}\right]}^T$ 轨道干扰输入变量为w(t)，设干扰输入矩阵为E_c，可得系统的状态空间方程式

$\left\{\begin{array}{c}\stackrel{.}{x}\left(t\right)=A_{c}x\left(t\right){+B}_{c}u\left(t\right)+E_{c}w\left(t\right)\\ y\left(t\right)=C_{c}x\left(t\right)+D_{c}u\left(t\right)\end{array}\right.$

那么离散化得到单速率系统模型状态方程式

$\left\{\begin{array}{c}x(k+1)=\mathrm{Ax}\left(k\right)+\mathrm{Bu}\left(k\right)+E_{C}w\left(k\right)\\ y\left(k\right)=C_{c}x\left(k\right)+D_{c}u\left(k\right)\end{array}\right.$

然后离散化得到单速率系统模型状态方程式，再利用提升技术得到多速率系统的状态方程 $\left\{\begin{array}{c}x(k+1)=A_{T}x\left(k\right)+B_{T}u\left(k\right)\\ y\left(k\right)=C_{T}x\left(k\right)+D_{T}u\left(k\right)\end{array}\right.$

其中：

$A_T=A_p^q;B_T=\left[\begin{array}{ccccc}{A}_p^{q-1}{B}_p& A_p^{q-2}{B}_p& ...& A_p{B}_p& B_p\end{array}\right]$

$C_T=\left[\begin{array}{c}C\\ {\mathrm{CA}}_q\\ .\\ .\\ .\\ {\mathrm{CA}}^{p-1}\end{array}\right]D_T=\left[\begin{array}{cccccccc}D& 0& ...& ...& ...& ...& ...& 0\\ D_{10}& D_{11}& ...& D_{{1c}_1}& 0& ...& ...& 0\\ D_{20}& D_{21}& ...& ...& D_{{2c}_2}& 0& ...& 0\\ .& .& & & & & & .\\ .& .& & & & & & .\\ .& .& & & & & & .\\ D_{(p-1),0}& D_{(p-1),1}& ...& ...& ...& ...& D_{{(p-1),c}_{p-1}}& 0\end{array}\right];$

3、步骤2中由于系统采用提升技术后，周期时变的多数率控制系统已经转化成为等效的时不变单速率系统，因此，多速率系统模型的预测控制器可根据单速率系统模型，设系统状态变量的维数为n，通过分析可知关于状态变量和控制变量有：

Δx(k)＝x(k)-x(k-1)，Δu(k)＝u(k)-u(k-1)；

4、依据步骤3计算可得：Δx(k+1)＝A_tΔx(k)+B_tΔu(k)，同理可得，Δy(k+1)＝y(k+1)-y(k)，那么将输出y(k)以状态变量Δx(k)来描述，可以推出下式：

Δy(k+1)＝C_tΔx(k+1)+D_tΔu(k+1)

                                     ；

＝C_tA_tΔx(k)+C_tB_tΔu(k)+D_tΔu(k+1)

5、同步骤4的操作可得新的状态变向量x(k)＝[Δx(k)^T y(k)^T]^T和控制变量u(k)＝[Δu(k)^TΔu(k+1)^T]^T后，可到得主动悬挂多速率系统改进的状态空间模型如下：

x(k+1)＝A_ttx(k)+B_ttu(k)

y(k)＝C_ttx(k)

其中： $A_{\mathrm{tt}}=\left[\begin{array}{cc}{A}_t& o_{q\×n}^T\\ C_t{A}_t& I_{q\×q}\end{array}\right]{,B}_{\mathrm{tt}}=\left[\begin{array}{cc}{B}_t& o_{p^2\×n}^T\\ C_t{B}_t& D_t\end{array}\right],$ C_tt＝[o_q×n I_q×q]，o_q×n为q×n维零矩阵，I_q×q为q×q维单位矩阵，o_p×n为p×n维零矩阵，由于状态空间维数的扩展，上式中A_tt，B_tt和C_tt分别扩展成为(n+q)×(n+q)，(n+q)(p*p+p*p)和q×(n+q)维矩阵；

6、系统采样时刻为k_i(k_i＞0)，N_c和N_p分别为控制时域和预测时域，且N_c≤N_p，那么预测状态变量如下：(A_dd，B_dd变为A_tt，B_tt)

$\begin{array}{c}x(k_i+1|k_i)=A_{\mathrm{dd}}x\left(k_i\right)+B_{\mathrm{dd}}u\left(k_i\right)\\ x(k_i+2|k_i)=A_{\mathrm{dd}}x(k_i+1|k_i)+B_{\mathrm{dd}}u(k_i+1)\\ =A_{\mathrm{dd}}^{2}x\left(k_i\right)+A_{\mathrm{dd}}{B}_{\mathrm{dd}}u\left(k_i\right)+B_{\mathrm{dd}}u(k_i+1)\\ .\\ .\\ .\\ x(k_i+N_p|k_i)=A_{\mathrm{dd}}^{N_p}x\left(k_i\right)+A_{\mathrm{dd}}^{N_d-1}{B}_{\mathrm{dd}}u\left(k_i\right)+A_{\mathrm{dd}}^{N_p-2}{B}_{\mathrm{dd}}u(k_i+1)\\ +...+A_{\mathrm{dd}}^{N_p-N_c}{B}_{\mathrm{dd}}u(k_i+N_c-1)\end{array}$

同理，根据预测状态变量可以得出预测输出变量：

$\begin{array}{c}y(k_i+1{|k}_i)=C_{\mathrm{dd}}{A}_{\mathrm{dd}}x\left(k_i\right)+C_{\mathrm{dd}}{B}_{\mathrm{dd}}u\left(k_i\right)\\ y(k_i+2{|k}_i)=C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{2}x\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}{B}_{\mathrm{dd}}u\left(k_i\right)+C_{\mathrm{dd}}{B}_{\mathrm{dd}}u(k_i+1)\\ y(k_i+3{|k}_i)=C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{3}x\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}^2{B}_{\mathrm{dd}}u\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}{B}_{\mathrm{dd}}u(k_i+1)\\ {+C}_{\mathrm{dd}}{B}_{\mathrm{dd}}u(k_i+2)\\ .\\ .\\ .\\ y(k_i+N_p|k_i)=C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{N_p}x\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{N_p-1}{B}_{\mathrm{dd}}u\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{N_p-1}{B}_{\mathrm{dd}}u(k_i+1)\\ +...+C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{N_p-N_c}{B}_{\mathrm{dd}}u(k_i+N_c-1)\end{array}$

根据当前的状态变量信息得出未来的控制变量等所有预测变量，分析可得控制变量U和输出变量Y表示如下：

Y＝Fx(k_i)+ΦU

其中：

$F=\left[\begin{array}{c}{C}_{\mathrm{tt}}{A}_{\mathrm{tt}}\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^2\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^3\\ .\\ .\\ .\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^{N_p}\end{array}\right],\mathrm{\Φ}=\left[\begin{array}{ccccc}{C}_{\mathrm{tt}}{B}_{\mathrm{tt}}& 0& 0& ...& 0\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{B}_{\mathrm{tt}}& 0& ...& 0\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^2{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{A}_{\mathrm{tt}}{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{B}_{\mathrm{tt}}& ...& 0\\ .& & & & \\ .& & & & \\ .& & & & \\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^{N_p-1}{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{A}_{\mathrm{tt}}^{N_p-2}{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{A}_{\mathrm{tt}}^{N_p-3}{B}_{\mathrm{tt}}& ...& {\mathrm{CA}}_{\mathrm{tt}}^{N_p-N_c}{B}_{\mathrm{tt}}\end{array}\right];$

7、对于扩展后状态空间模型的预测控制算法，注意到该多输入-多输出系统中状态变量、控制变量和输出变量的维数，根据当前的状态变量信息便可以得出未来的控制变量等所有预测变量，分析可得控制变量U和输出变量Y表示如下：

U＝[u(k_i)^Tu(k_i+1)^Tu(k_i+2)^T…u(k_i+N_c-1)^T]^T

Y＝[y(k_i+1|k_i)^Ty(k_i+2|k_i)^Ty(k_i+3|k_i)^T…y(k_i+N_p|k_i)^T]^T

采用上述离散标准正交基函数(Laguerre函数)，近似求得控制变量U中包含的u(k_i),u(k_i+1),...,u(k_i+N_c-1)等序列；

8、主动悬挂多速率系统模型属于多输入-多输出系统，其目标函数如下：

$J={\mathrm{\η}}^T\mathrm{\Ω\η}+{2\mathrm{\η}}^T\mathrm{\Ψx}\left(k_i\right)+\underset{p=1}{\overset{N_p}{\mathrm{\Σ}}}x{\left(k_i\right)}^T{\left(A_{\mathrm{tt}}^T\right)}^{p}Q{A}_{\mathrm{tt}}^{p}x\left(k_i\right)$

式中： $\mathrm{\Ω}=\underset{p=1}{\overset{N_p}{\mathrm{\Σ}}}\mathrm{\φ}\left(p\right)\mathrm{Q\φ}{\left(p\right)}^T+R_L,\mathrm{\Ψ}=\underset{p=1}{\overset{N_p}{\mathrm{\Σ}}}\mathrm{\φ}\left(p\right){\mathrm{QA}}_{\mathrm{tt}}^p$

当目标函数满足η＝-Ω^-1Ψx(k_i)时，J为最小值，得到最优解，那么，根据滚动优化控制原理，设未来的参考轨迹在预测时域之内，那么系统的状态反馈控制变量为：

u(k)＝-K_tmpcx(k)

由此可得，高速列车主动悬挂多速率系统模型预测控制的状态反馈控制增益矩阵K_dmpc为：

本发明的技术效果是：本发明在高速列车主动悬挂系统模型的基础上，开展基于一种标准正交基函数-Laguerre函数的预测控制算法研究，将悬挂系统的状态空间模型作为预测控制模型，结合预测控制中的滚动优化等基本原理，考虑到列车在运行时受到的不确定轨道干扰激励等影响因素，为高速列车主动悬挂多速率及含时滞的多速率系统设计相应的预测控制器，根据系统模型，不断调整优化控制器参数，使系统达到最佳的控制效果，以优化高速列车在运行时的振动主动控制性能，保证列车运行的操作稳定性和舒适性。

附图说明

图1为本发明的主动悬挂系统动力学模型示意图。

图2为本发明的多速率数字系统原理图。

图3为本发明的离散Laguerre网络的示意图。

图4为本发明的实施主动悬挂多速率预测控制的示意图。

在图中，M为车体质量；m为前、后转向架质量；I_Φ为车体点头转动惯量；I_θ为车体侧滚转动惯量；I_Φ1、I_Φ2为前后转向架点头转动惯量；I_θ1、I_θ2为前后转向架侧滚转动惯量；k_l、k₂为一、二系悬挂弹簧垂向刚度；c₁、c₂为垂向减振器阻尼；l、l₁为悬挂到构架质心的纵向距离；z_b为车体沉浮位移；z₁为前转向架沉浮位移；z₂为后转向架沉浮位移；Φ为车体点头角位移；θ为车体侧滚角位移；θ₁、θ₂为转向架侧滚角位移；a、b为一、二系悬挂弹簧的横向距离之半；v为列车运行速度；z_v1，z′_v1，z_v2，z′_v2，z_v3，z′_v3，z_v4，z′_v4为轮对处轨道面输入位移；u₁、u′₁，u₂、u′₂，u₃、u′₃，u₄、u′₄为前转向架一系作动器控制力；u₅、u′₅，u₆、u′₆为后转向架一系作动器控制力。

具体实施方式

结合图1、2、3、4来具体说明本发明，其具体有以下步骤:

1、建立列车主动悬挂动力学模型，求得状态空间方程式

$\left\{\begin{array}{c}\stackrel{.}{x}\left(t\right)=A_{c}x\left(t\right){+B}_{c}u\left(t\right)+E_{c}w\left(t\right)\\ y\left(t\right)=C_{c}x\left(t\right)+D_{c}u\left(t\right)\end{array}\right.$

其中： $A_C=\left[\begin{array}{cc}{0}_{9*9}& I_{9*9}\\ D_1& F\end{array}\right],B_C=\left[\begin{array}{cc}{0}_{9*12}& 0_{12*8}\\ R& G\end{array}\right],$

其中：0_9*9为9*9阶零矩阵，I_9*9为9*9阶单位矩阵，

$R=\left[\begin{array}{cccccccccccc}1/M& 1/M& 1/M& 1/M& 0& 0& 0& 0& 0& 0& 0& 0\\ -1/I_{\mathrm{\φ}}& -1/I_{\mathrm{\φ}}& 1/I_{\mathrm{\φ}}& 1/I_{\mathrm{\φ}}& 0& 0& 0& 0& 0& 0& 0& 0\\ -b/I_{\mathrm{\θ}}& b/I_{\mathrm{\θ}}& -b/I_{\mathrm{\θ}}& b/I_{\mathrm{\θ}}& 0& 0& 0& 0& 0& 0& 0& 0\\ 0& 0& -1/m& -1/m& 0& 0& 0& 0& 1/m& 1/m& 1/m& 1/m\\ 0& 0& 0& 0& 0& 0& 0& 0& -1_1/I_{\mathrm{\φ}2}& -1_1/I_{\mathrm{\φ}2}& 1_1/I_{\mathrm{\φ}2}& 1_1/I_{\mathrm{\φ}2}\\ 0& 0& -a/I_{\mathrm{\θ}2}& a/I_{\mathrm{\θ}2}& 0& 0& 0& 0& -a/I_{\mathrm{\θ}2}& a/I_{\mathrm{\θ}2}& -a/I_{\mathrm{\θ}2}& a/I_{\mathrm{\θ}2}\\ -1/m& -1/m& 0& 0& 1/m& 1/m& 1/m& 1/m& 0& 0& 0& 0\\ 0& 0& 0& 0& {-1}_1/I_{\mathrm{\φ}1}& {-1}_1/I_{\mathrm{\φ}1}& 1_1/I_{\mathrm{\φ}1}& 1_1/I_{\mathrm{\φ}1}& 0& 0& 0& 0\\ -a/I_{\mathrm{\θ}2}& a/I_{\mathrm{\θ}2}& 0& 0& -a/I_{\mathrm{\θ}2}& a/I_{\mathrm{\θ}2}& -a/I_{\mathrm{\θ}2}& a/I_{\mathrm{\θ}2}& 0& 0& 0& 0\end{array}\right]$

$G=\left[\begin{array}{cccccccc}0& 0& 0& 0& k_1/m& k_1/m& k_1/m& k_1/m\\ 0& 0& 0& 0& -k_1{1}_1/I_{\mathrm{\φ}2}& {-k}_1{1}_1/I_{\mathrm{\φ}2}& k_1{1}_1/I_{\mathrm{\φ}2}& k_1{1}_1/I_{\mathrm{\φ}2}\\ 0& 0& 0& 0& -k_{1}a/I_{\mathrm{\θ}2}& k_{1}a/I_{\mathrm{\θ}2}& -k_{1}a/I_{\mathrm{\θ}2}& k_{1}a/I_{\mathrm{\θ}2}\\ k_1/m& k_1/m& k_1/m& k_1/m& 0& 0& 0& 0\\ {-k}_1{1}_1/I_{\mathrm{\φ}2}& {-k}_1{1}_1/I_{\mathrm{\φ}2}& k_1{1}_1/I_{\mathrm{\φ}2}& k_1{1}_1/I_{\mathrm{\φ}2}& 0& 0& 0& 0\\ {-k}_{1}a/I_{\mathrm{\θ}2}& k_{1}a/I_{\mathrm{\θ}2}& {-k}_{1}a/I_{\mathrm{\θ}2}& k_{1}a/I_{\mathrm{\θ}2}& 0& 0& 0& 0\end{array}\right]$

；

2、离散化得到单速率系统模型状态方程式

$\left\{\begin{array}{c}x(k+1)=\mathrm{Ax}\left(k\right)+\mathrm{Bu}\left(k\right)+E_{C}w\left(k\right)\\ y\left(k\right)=C_{c}x\left(k\right)+D_{c}u\left(k\right)\end{array}\right.$

其中A＝exp(A_ch)∈R^18*18；

3、多速率系统的状态空间模型：

$\left\{\begin{array}{c}x(k+1)=A_{T}x\left(k\right)+B_{T}u\left(k\right)\\ y\left(k\right)=C_{T}x\left(k\right)+D_{T}u\left(k\right)\end{array}\right.;$

4、多速率系统模型的预测控制器可根据单速率系统模型，设系统状态变量的维数为n，通过分析可知关于状态变量和控制变量有：

Δx(k)＝x(k)-x(k-1)，Δu(k)＝u(k)-u(k-1)

因此计算可得：Δx(k+1)＝A_tΔx(k)+B_tΔu(k)，同理可得，Δy(k+1)＝y(k+1)-y(k)，那么将输出y(k)以状态变量Δx(k)来描述，可以推出下式：

Δy(k+1)＝C_tΔx(k+1)+D_tΔu(k+1)

＝C_tA_tΔx(k)+C_tB_tΔu(k)+D_tΔu(k+1)

同理取新的状态变向量x(k)＝[Δx(k)^T y(k)^T]^T和控制变量u(k)＝[Δu(k)^TΔu(k+1)^T]^T后，可到得主动悬挂多速率系统改进的状态空间模型如下：

x(k+1)＝A_ttx(k)+B_ttu(k)

y(k)＝C_ttx(k)

其中： $A_{\mathrm{tt}}=\left[\begin{array}{cc}{A}_t& o_{q\×n}^T\\ C_t{A}_t& I_{q\×q}\end{array}\right]{,B}_{\mathrm{tt}}=\left[\begin{array}{cc}{B}_t& o_{p^2\×n}^T\\ C_t{B}_t& D_t\end{array}\right],$ C_tt＝[o_q×n I_q×q]，o_q×n为q×n维零矩阵，I_q×q为q×q维单位矩阵，o_p×n为p×n维零矩阵。由于状态空间维数的扩展，上式中A_tt，B_tt和C_tt分别扩展成为(n+q)×(n+q)，(n+q)(p*p+p*p)和q×(n+q)维矩阵，

设系统采样时刻为k_i(k_i＞0)，N_c和N_p分别为控制时域和预测时域，且N_c≤N_p。那么预测状态变量如下：(A_dd，B_dd变为A_tt，B_tt)

$\begin{array}{c}x(k_i+1|k_i)=A_{\mathrm{dd}}x\left(k_i\right)+B_{\mathrm{dd}}u\left(k_i\right)\\ x(k_i+2|k_i)=A_{\mathrm{dd}}x(k_i+1|k_i)+B_{\mathrm{dd}}u(k_i+1)\\ =A_{\mathrm{dd}}^{2}x\left(k_i\right)+A_{\mathrm{dd}}{B}_{\mathrm{dd}}u\left(k_i\right)+B_{\mathrm{dd}}u(k_i+1)\\ .\\ .\\ .\\ x(k_i+N_p|k_i)=A_{\mathrm{dd}}^{N_p}x\left(k_i\right)+A_{\mathrm{dd}}^{N_d-1}{B}_{\mathrm{dd}}u\left(k_i\right)+A_{\mathrm{dd}}^{N_p-2}{B}_{\mathrm{dd}}u(k_i+1)\\ +...+A_{\mathrm{dd}}^{N_p-N_c}{B}_{\mathrm{dd}}u(k_i+N_c-1)\end{array}$

同理，根据预测状态变量可以得出预测输出变量：

$\begin{array}{c}y(k_i+1{|k}_i)=C_{\mathrm{dd}}{A}_{\mathrm{dd}}x\left(k_i\right)+C_{\mathrm{dd}}{B}_{\mathrm{dd}}u\left(k_i\right)\\ y(k_i+2{|k}_i)=C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{2}x\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}{B}_{\mathrm{dd}}u\left(k_i\right)+C_{\mathrm{dd}}{B}_{\mathrm{dd}}u(k_i+1)\\ y(k_i+3{|k}_i)=C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{3}x\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}^2{B}_{\mathrm{dd}}u\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}{B}_{\mathrm{dd}}u(k_i+1)\\ {+C}_{\mathrm{dd}}{B}_{\mathrm{dd}}u(k_i+2)\\ .\\ .\\ .\\ y(k_i+N_p|k_i)=C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{N_p}x\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{N_p-1}{B}_{\mathrm{dd}}u\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{N_p-1}{B}_{\mathrm{dd}}u(k_i+1)\\ +...+C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{N_p-N_c}{B}_{\mathrm{dd}}u(k_i+N_c-1)\end{array}$

根据当前的状态变量信息得出未来的控制变量等所有预测变量，分析可得控制变量U和输出变量Y表示如下：

Y＝Fx(k_i)+ΦU

其中：

$F=\left[\begin{array}{c}{C}_{\mathrm{tt}}{A}_{\mathrm{tt}}\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^2\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^3\\ .\\ .\\ .\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^{N_p}\end{array}\right],\mathrm{\Φ}=\left[\begin{array}{ccccc}{C}_{\mathrm{tt}}{B}_{\mathrm{tt}}& 0& 0& ...& 0\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{B}_{\mathrm{tt}}& 0& ...& 0\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^2{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{A}_{\mathrm{tt}}{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{B}_{\mathrm{tt}}& ...& 0\\ .& & & & \\ .& & & & \\ .& & & & \\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^{N_p-1}{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{A}_{\mathrm{tt}}^{N_p-2}{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{A}_{\mathrm{tt}}^{N_p-3}{B}_{\mathrm{tt}}& ...& {\mathrm{CA}}_{\mathrm{tt}}^{N_p-N_c}{B}_{\mathrm{tt}}\end{array}\right]$

对于扩展后状态空间模型的预测控制算法，注意到该多输入-多输出系统中状态变量、控制变量和输出变量的维数，根据当前的状态变量信息便可以得出未来的控制变量等所有预测变量，分析可得控制变量U和输出变量Y表示如下：

U＝[u(k_i)^Tu(k_i+1)^Tu(k_i+2)^T…u(k_i+N_c-1)^T]^T

Y＝[y(k_i+1|k_i)^Ty(k_i+2|k_i)^Ty(k_i+3|k_i)^T…y(k_i+N_p|k_i)^T]^T

采用上述离散标准正交基函数(Laguerre函数)，近似求得控制变量U中包含的u(k_i),u(k_i+1),...,u(k_i+N_c-1)等序列，

离散时间Laguerre网络的z变换形式如下：

$\begin{array}{c}{\mathrm{\Γ}}_1\left(z\right)=\frac{\sqrt{{1-a}^2}}{1-a{z}^{-1}}\\ {\mathrm{\Γ}}_2\left(z\right)=\frac{\sqrt{1-a^2}}{1-{\mathrm{az}}^{-1}}\frac{z^{-1}-a}{1-{\mathrm{az}}^{-1}}\\ .\\ .\\ .\\ {\mathrm{\Γ}}_N\left(z\right)=\frac{\sqrt{{1-a}^2}}{1-{\mathrm{az}}^{-1}}{\left(\frac{z^{-1}-a}{{1-\mathrm{az}}^{-1}}\right)}^{N-1}\end{array}$

其中：a是离散时间Laguerre网络的极点，又叫比例因子，是自行设定的参数，为使其稳定设0≤a＜1，N是一个给定的项数，

可以得出：

${\mathrm{\Γ}}_k\left(z\right)={\mathrm{\Γ}}_{k-1}\left(z\right)\frac{z^{-1}-a}{1-{\mathrm{az}}^{-1}}$

根据这种关系，Laguerre网络可以概括成如图(3)所示，

用l₁(k)表示Γ₁(z,a)的逆z变换，类似的，用l_N(k)表示Γ_N(z,a)的逆z变换。那么，多输入-多输出系统离散时间Laguerre函数则可以表示为如下形式：

$L_i\left(k\right)={\left[\begin{array}{cccc}{1}_1^i\left(k\right)& 1_2^i\left(k\right)& ...& 1_{N_i}^i\left(k\right)\end{array}\right]}^T$

可以得出该函数满足：

L_i(k+1)＝A_lL_i(k)

式中：A_l是一个包含参数a的N×N矩阵，β＝(1-a²)，初始条件表示如下：

$L{\left(0\right)}^T=\sqrt{\mathrm{\β}}\left[\begin{array}{ccccccc}1& -a& a^2& -a^3& ...& {(-1)}^{N-1}& a^{N-1}\end{array}\right]$

在未来任意采样时刻k，u_i(k)＝L_i(k)^Tη_i

式中：η是一个参数向量，包含N个Laguerre系数， ${\mathrm{\η}}^T=\left[\begin{array}{cccc}{\mathrm{\η}}_1^T& {\mathrm{\η}}_2^T& ...& {\mathrm{\η}}_p^T\end{array}\right],$

多输入-多输出系统中，u(k)＝[u₁(k) u₂(k) ... u_p(k)]^T，B＝[B₁ B₂ ... B_p]^T，p为系统的输入个数，系统中每个输入都有一个对应的N和a，因此N和a的个数和输入的个数相等，在预测控制的实际应用中，灵活的选择参数a和N都可对控制器效果产生影响，根据上文已知的状态变量信息，可以预测未来在p时刻的状态：

$x(k_i+{p|k}_i)=A_{\mathrm{dd}}^{p}x\left(k_i\right)+\mathrm{\φ}{\left(p\right)}^T\mathrm{\η}$

式中： $\mathrm{\φ}{\left(p\right)}^T=\underset{j=0}{\overset{p-1}{\mathrm{\Σ}}}{A}^{p-j-1}\left[\begin{array}{cccc}{B}_1{L}_1{\left(j\right)}^T& B_2{L}_2{\left(j\right)}^T& ...& B_p{L}_p{\left(j\right)}^T\end{array}\right]$

主动悬挂多速率系统模型属于多输入-多输出系统，其目标函数如下：

$J={\mathrm{\η}}^T\mathrm{\Ω\η}+{2\mathrm{\η}}^T\mathrm{\Ψx}\left(k_i\right)+\underset{p=1}{\overset{N_p}{\mathrm{\Σ}}}x{\left(k_i\right)}^T{\left(A_{\mathrm{tt}}^T\right)}^{p}Q{A}_{\mathrm{tt}}^{p}x\left(k_i\right)$

式中： $\mathrm{\Ω}=\underset{p=1}{\overset{N_p}{\mathrm{\Σ}}}\mathrm{\φ}\left(p\right)\mathrm{Q\φ}{\left(p\right)}^T+R_L,\mathrm{\Ψ}=\underset{p=1}{\overset{N_p}{\mathrm{\Σ}}}\mathrm{\φ}\left(p\right){\mathrm{QA}}_{\mathrm{tt}}^p$

当目标函数满足η＝-Ω^-1Ψx(k_i)时，J为最小值，得到最优解。那么，根据滚动优化控制原理，设未来的参考轨迹在预测时域之内，那么系统的状态反馈控制变量为：

u(k)＝-K_tmpcx(k)

由此可得，高速列车主动悬挂多速率系统模型预测控制的状态反馈控制增益矩阵K_dmpc为：

以上所述的实施例仅仅是对本发明的优选实施方式进行描述，并非对本发明的范围进行限定，在不脱离本发明设计精神的前提下，本领域普通技术人员对本发明的技术方案作出的各种变形和改进，均应落入本发明权利要求书确定的保护范围内。

Claims (1)

1.一种应用在列车主动悬挂系统上的多速率预测控制方法，其技术方案步骤为：

1)以高速列车车辆1/4车体为研究对象，建立包含车体以及前后转向架的侧滚、沉浮和点头等九自由度的横向-垂向耦合的列车主动悬挂系统动力学模型；

2)根据步骤一建立的动力学模型得到悬挂系统动力学方程：取状态变量为：

$X=\begin{array}{ccccccccccccc}[Z_b& \mathrm{\φ}& \mathrm{\θ}& Z_2& {\mathrm{\φ}}_2& {\mathrm{\θ}}_1& Z_1& {\mathrm{\φ}}_1& {\mathrm{\θ}}_2& {\stackrel{\·}{Z}}_b& \stackrel{\·}{\mathrm{\φ}}& \stackrel{\·}{\mathrm{\θ}}& {\stackrel{\·}{Z}}_2\\ {\stackrel{\·}{\mathrm{\φ}}}_2& {\stackrel{\·}{\mathrm{\θ}}}_1& {\stackrel{\·}{Z}}_1& {\stackrel{\·}{\mathrm{\φ}}}_1& {\stackrel{\·}{\mathrm{\θ}}}_2{]}^T& & & & & & & & \end{array}$

控制变量为：u_c＝[u₁ u′₁ u₂ u′₂ u₃ u′₃ u₄ u′₄ u₅ u′₅ u₆ u′₆]^T输出变量为： $y={\left[\begin{array}{ccc}{\stackrel{\·\·}{Z}}_b& \stackrel{\·\·}{\mathrm{\φ}}& \stackrel{\·\·}{\mathrm{\θ}}\end{array}\right]}^T$ 轨道干扰输入变量为w(t)，设干扰输入矩阵为E_c，可得系统的状态空间方程式

$\left\{\begin{array}{c}\stackrel{\·}{x}\left(t\right)=A_{c}x\left(t\right)+B_{c}u\left(t\right)+E_{c}w\left(t\right)\\ y\left(t\right)=C_{c}x\left(t\right)+D_{c}u\left(t\right)\end{array}\right.$

那么离散化得到单速率系统模型状态方程式

$\left\{\begin{array}{c}x(k+1)=\mathrm{Ax}\left(k\right)+\mathrm{Bu}\left(k\right)+E_{C}w\left(k\right)\\ y\left(k\right)=C_{c}x\left(k\right)+C_{c}u\left(k\right)\end{array}\right.$

然后离散化得到单速率系统模型状态方程式，再利用提升技术得到多速率系统的状态方程 $\left\{\begin{array}{c}x(k+1)=A_{T}x\left(k\right)+B_{T}u\left(k\right)\\ y\left(k\right)=C_{T}x\left(k\right)+D_{T}u\left(k\right)\end{array}\right.$

其中：

$A_T=A_p^q;B_T=\left[\begin{array}{ccccc}{A}_p^{q-1}{B}_p& A_p^{q-2}{B}_p& ...& A_p{B}_p& B_p\end{array}\right]$

$C_T=\left[\begin{array}{c}C\\ C{A}_q\\ .\\ .\\ .\\ C{A}^{p-1}\end{array}\right]D_T=\left[\begin{array}{cccccccc}D& 0& ...& ...& ...& ...& ...& 0\\ D_{10}& D_{11}& ...& D_{1c_1}& 0& ...& ...& 0\\ D_{20}& D_{21}& ...& ...& D_{2c_2}& 0& ...& 0\\ .& .& & & & & & .\\ .& .& & & & & & .\\ .& .& & & & & & .\\ D_{(p-1),0}& D_{(p-1),1}& ...& ...& ...& ...& D_{(p-1),c_{p-1}}& 0\end{array}\right];$

3)步骤二中由于系统采用提升技术后，周期时变的多数率控制系统已经转化成为等效的时不变单速率系统，因此，多速率系统模型的预测控制器可根据单速率系统模型，设系统状态变量的维数为n，通过分析可知关于状态变量和控制变量有：

Δx(k)＝x(k)-x(k-1)，Δu(k)＝u(k)-u(k-1)；

4)依据步骤三计算可得：Δx(k+1)＝A_tΔx(k)+B_tΔu(k)，同理可得，Δy(k+1)＝y(k+1)-y(k)，那么将输出y(k)以状态变量Δx(k)来描述，可以推出下式：

Δy(k+1)＝C_tΔx(k+1)+D_tΔu(k+1)

＝C_tA_tΔx(k)+C_tB_tΔu(k)+D_tΔu(k+1)；

5)同步骤四的操作可得新的状态变向量x(k)＝[Δx(k)^T y(k)^T]^T和控制变量u(k)＝[Δu(k)^T Δu(k+1)^T]^T后，可到得主动悬挂多速率系统改进的状态空间模型如下：

x(k+1)＝A_ttx(k)+B_ttu(k)

y(k)＝C_ttx(k)

其中： $A_{\mathrm{tt}}=\left[\begin{array}{cc}{A}_t& o_{q\×n}^T\\ C_t{A}_t& I_{q\×q}\end{array}\right],B_{\mathrm{tt}}=\left[\begin{array}{cc}{B}_t& o_{p^2\×n}^T\\ C_t{B}_t& D_t\end{array}\right],$ C_tt＝[o_q×n I_q×q]，o_q×n为q×n维零矩阵，I_q×q为q×q维单位矩阵，o_p×n为p×n维零矩阵，由于状态空间维数的扩展，上式中A_tt，B_tt和C_tt分别扩展成为(n+q)×(n+q)，(n+q)(p*p+p*p)和q×(n+q)维矩阵；

6)系统采样时刻为k_i(k_i＞0)，N_c和N_p分别为控制时域和预测时域，且N_c≤N_p。那么预测状态变量如下：(A_dd，B_dd变为A_tt，B_tt)

x(k_i+1|k_i)＝A_ddx(k_i)+B_ddu(k_i)

$\begin{array}{c}x(k_i+2|k_i)=A_{\mathrm{dd}}x(k_i+1|k_i)+B_{\mathrm{dd}}u(k_i+1)\\ =A_{\mathrm{dd}}^{2}x\left(k_i\right)+A_{\mathrm{dd}}{B}_{\mathrm{dd}}u\left(k_i\right)+B_{\mathrm{dd}}u(k_i+1)\end{array}$

.

.

.

$\begin{array}{c}x(k_i+N_p|k_i)=A_{\mathrm{dd}}^{N_p}x\left(k_i\right)+A_{\mathrm{dd}}^{N_p-1}{B}_{\mathrm{dd}}u\left(k_i\right)+A_{\mathrm{dd}}^{N_p-2}{B}_{\mathrm{dd}}u(k_i+1)\\ +...+A_{\mathrm{dd}}^{N_p-N_c}{B}_{\mathrm{dd}}u(k_i+N_c-1)\end{array}$

同理，根据预测状态变量可以得出预测输出变量：

y(k_i+1|k_i)＝C_ddA_ddx(k_i)+C_ddB_ddu(k_i)

$y(k_i+2|k_i)=C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{2}x\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}{B}_{\mathrm{dd}}u\left(k_i\right)+C_{\mathrm{dd}}{B}_{\mathrm{dd}}u(k_i+1)$

$\begin{array}{c}y(k_i+3|k_i)=C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{3}x\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}^2{B}_{\mathrm{dd}}u\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}{B}_{\mathrm{dd}}u(k_i+1)\\ +C_{\mathrm{dd}}{B}_{\mathrm{dd}}u(k_i+2)\end{array}$

.

.

.

$\begin{array}{c}y(k_i+N_p|k_i)=C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{N_p}x\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{N_p-1}{B}_{\mathrm{dd}}u\left(k_i\right)+C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{N_p-2}{B}_{\mathrm{dd}}u(k_i+1)\\ +...+C_{\mathrm{dd}}{A}_{\mathrm{dd}}^{N_p-N_c}{B}_{\mathrm{dd}}u(k_i+N_c-1)\end{array}$

根据当前的状态变量信息得出未来的控制变量等所有预测变量，分析可得控制变量U和输出变量Y表示如下：

Y＝Fx(k_i)+ΦU

其中：

$F=\left[\begin{array}{c}{C}_{\mathrm{tt}}{A}_{\mathrm{tt}}\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^2\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^3\\ .\\ .\\ .\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^{N_p}\end{array}\right],\mathrm{\Φ}=\left[\begin{array}{ccccc}{C}_{\mathrm{tt}}{B}_{\mathrm{tt}}& 0& 0& ...& 0\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{B}_{\mathrm{tt}}& 0& ...& 0\\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^2{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{A}_{\mathrm{tt}}{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{B}_{\mathrm{tt}}& ...& 0\\ .& & & & \\ .& & & & \\ .& & & & \\ C_{\mathrm{tt}}{A}_{\mathrm{tt}}^{N_p-1}{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{A}_{\mathrm{tt}}^{N_p-2}{B}_{\mathrm{tt}}& C_{\mathrm{tt}}{A}_{\mathrm{tt}}^{N_p-3}{B}_{\mathrm{tt}}& ...& C{A}_{\mathrm{tt}}^{N_p-N_c}{B}_{\mathrm{tt}}\end{array}\right];$

7)对于扩展后状态空间模型的预测控制算法，注意到该多输入-多输出系统中状态变量、控制变量和输出变量的维数，根据当前的状态变量信息便可以得出未来的控制变量等所有预测变量，分析可得控制变量U和输出变量Y表示如下：

U＝[u(k_i)^Tu(k_i+1)^Tu(k_i+2)^T…u(k_i+N_c-1)^T]^T

Y＝[y(k_i+1|k_i)^Ty(k_i+2|k_i)^Ty(k_i+3|k_i)^T…y(k_i+N_p|k_i)^T]^T

采用上述离散标准正交基函数(Laguerre函数)，近似求得控制变量U中包含的u(k_i),u(k_i+1),...,u(k_i+N_c-1)等序列；

8)主动悬挂多速率系统模型属于多输入-多输出系统，其目标函数如下：

$J={\mathrm{\η}}^T\mathrm{\Ω\η}+2{\mathrm{\η}}^T\mathrm{\Ψx}\left(k_i\right)+\underset{p=1}{\overset{N_p}{\mathrm{\Σ}}}x{\left(k_i\right)}^T{\left(A_{\mathrm{tt}}^T\right)}^{p}Q{A}_{\mathrm{tt}}^{p}x\left(k_i\right)$

式中： $\mathrm{\Ω}=\underset{p=1}{\overset{N_p}{\mathrm{\Σ}}}\mathrm{\φ}\left(p\right)\mathrm{Q\φ}{\left(p\right)}^T+R_L,\mathrm{\Ψ}=\underset{p=1}{\overset{N_p}{\mathrm{\Σ}}}\mathrm{\φ}\left(p\right)Q{A}_{\mathrm{tt}}^p$

当目标函数满足η＝-Ω^-1Ψx(k_i)时，J为最小值，得到最优解，那么根据滚动优化控制原理，设未来的参考轨迹在预测时域之内，那么系统的状态反馈控制变量为：

u(k)＝-K_tmpcx(k)

由此可得，高速列车主动悬挂多速率系统模型预测控制的状态反馈控制增益矩阵K_dmpc为：