一种基于动态优化的多模型广义预测控制系统及其控制方法

摘要

一种基于动态优化的多模型广义预测控制系统，其包括动态优化层、MPC层以及基础控制层；所述动态优化层位于上层，其计算出关键控制变量的优化值作为MPC层的最优设定值；所述MPC层位于下层在待优化变量满足模型动态行为的条件下采用滚动优化的预测算法对该待优化变量进行调节，使其跟踪S1中得出的最优设定值；所述基础控制层位于底层，将该待优化变量的最终优化值送至执行机构。本发明降低系统成本消耗，提高系统经济效益，并且可以提高系统暂态性能以及系统模型参数跳变的调节能力，同时还可以有效的消除扰动对系统输出的干扰。

主权项

一种基于动态优化的多模型广义预测控制系统的优化方法，其特征在于，包括以下步骤：S1：动态优化层针对系统的经济目标函数及其时变约束采用控制向量参数化与粒子群优化算法的结合的方法获取关键变量的最优值轨迹，并将该轨迹将作为下层MPC层的最优设定值参考轨迹；S2：MPC层采用滚动优化的预测算法对过程中待优化的变量满足模型参数变化和干扰影响动态行为的条件下进行调节，使其跟踪S1中得出的该MPC层的最优设定值参考轨迹，并且采用多个固定模型和自适应模型来并行辨识系统的动态特性；S3：基础层通过PID作用抑制、消除进入到过程中的扰动对输出的影响，并将该变量的最终优化值送到执行结构；所述关键变量的最优值轨迹的获取过程包括以下步骤：S11：首先将时间区间[t₀,t_f]分成相同的多段，每段用分段常数轨迹近似最优轨迹，得到多个待优化控制参数；S12：粒子群初始化：设置粒子数m、维数D、两个学习因子c₁,c₂，位置上下界x_max,x_min及最大迭代次数M，初始位置，初始速度，初始化全局最优解gBest以及局部最优解pBest；S13：计算各粒子适应度，更新局部最优值和全局最优值；S14：计算更新速度与更新位置，如果位置超过上下界则设定为边界值；S15：判断是否达到最大迭代次数，若没有则返回S13继续计算,若满足条件则输出当前最优值；在步骤S2中，被控对象为：A(z^‑1)y(k)＝B(z^‑1)u(k‑1)+ξ(k)/Δ，多模型集表示为：Δy_i(k)＝φ_i(k)^Tθ₀(k)+ξ_i(k)i＝1,2…,m,m+1,m+2其中为实际优化目标对象模型的系统参数，φ(k)为泛指由实际优化目标对象输入输出值按φ(k)＝[‑Δy(k‑1)…‑Δy(k‑n_a)Δu(k‑1)+…Δu(k‑n_b‑1)]形式组成的向量集合，φ_i(k)为第i个模型的输入输出值按φ(k)＝[‑Δy(k‑1)…‑Δy(k‑n_a)Δu(k‑1)+…Δu(k‑n_b‑1)]形式组成的向量集合，其中，Δy(k‑1)为k‑1时刻系统的输出值y(k‑1)与k‑2时刻系统输出值y(k‑2)的差值，Δu(k‑1)为k‑1时刻系统输入值u(k‑1)与k‑2时刻系统输入值u(k‑2)的差值；其中，ξ_i(k)表示在k时刻系统的噪声；n_a表示被控对象；A(z^‑1)y(k)＝B(z^‑1)u(k‑1)+ξ(k)/Δ中的A(z^‑1)的阶次，n_b表示被控对象A(z^‑1)y(k)＝B(z^‑1)u(k‑1)+ξ(k)/Δ的B(z^‑1)的阶次；表示根据递推最小二乘法进行系统动态特性辨识时在k‑1时刻的θ₀估计值；当i＝1,2,…,m时，θ_i(k)为固定模型的恒定值，相对于单个模型，m个固定模型可以提高系统的暂态性能，其中，$A\left(z^{-1}\right)=1+a_1{z}^{-1}+\mathrm{......}+a_{n_a}{z}^{-n_a}$$B\left(z^{-1}\right)=b_0+b_1{z}^{-1}+\mathrm{......}+b_{n_b}{z}^{-n_b}$式中z^‑1为后移算子；y(k)、u(k)、ξ(k)分别为系统的输出，输入和均值为零的白噪声序列；Δ＝1‑z^‑1为差分算子；自适应模型采用如下递推最小二乘法进行系统动态特性的辨识：K(k)＝P(k‑1)φ(k)[φ(k)^TP(k‑1)φ(k)+μ]^‑1$P\left(k\right)=\frac{1}{\mu }[1-K\left(k\right)\phi {\left(k\right)}^T]P(k-1)$式中0<μ<1为遗忘因子，K(k)为权因子，P(k)为正定协方差阵；多模型切换指标表示为：$J_i={\mathrm{\gamma e}}_i^2\left(k\right)+\eta \underset{j=1}{\overset{L}{\Sigma }}{\rho }^j{e}_i^2(k-j),i=1,2,\mathrm{...},m+2$上式为模型i在时刻k的性能指标，其中e_i(k)为第i个模型在k时刻的输出误差，γ和η是当前和过去时刻的误差权重，ρ为误差遗忘因子，L为过去时刻误差长度，在k时刻，会根据性能指标J_i最小切换到相应模型；所述PID的性能优化指标为:$\min J\left(k\right)=E\{\underset{j=1}{\overset{N}{\Sigma }}{[y(k+j|k)-w_r(k+j)]}^2+\underset{j=1}{\overset{N_u}{\Sigma }}\lambda {\left[\Delta u(k+j-1|k)\right]}^2\},$式中E{·}为数学期望，w_r为w_r(k)的缩写且为输出期望值，N_u为控制时域，λ为控制加权系数，k表示时刻，y(k+j/k)为根据k时刻的输入输出预测得到的k+j时刻的输出；u(k+j‑1/k)为根据k时刻的输入输出预测得到的k+j‑1时刻的输入；对象输出的期望值为w_r(k+j)＝w(k+j)j∈{1,2,…,N}，式中w(k+j)为上层动态优化所给的关键变量设定值，N为优化时域终止时刻；对象输出的期望值为：w_r(k+j)＝w(k+j)j∈{1,2,…,N}式中w(k+j)为上层动态优化所给的关键变量设定值，即模型输出设定值，N为优化时域终止时刻；为了计算出第i步后输出预测值，首先引入丢番图方程： 1＝E_j(z^‑1)A(z^‑1)Δ+z^‑jF_j(z^‑1) E_j(z^‑1)B(z^‑1)＝G_j(z^‑1)+z^‑jH_j(z^‑1)式中： E_j(z^‑1)＝e_j,0+e_j,1z^‑1+…+e_j,j‑1z^‑(j‑1)$F_j\left(z^{-1}\right)=f_{j,0}+f_{j,1}{z}^{-1}+\mathrm{...}+f_{j,n_a}{z}^{-n_a}$ G_j(z^‑1)＝g_j,0+g_j,1z^‑1+…+g_j,j‑1z^‑(j‑1)$H_j\left(z^{-1}\right)=h_{j,0}+h_{j,1}{z}^{-1}+\mathrm{...}+h_{j,n_b-1}{z}^{-(n_b-1)}$$\Delta u\left(k\right)=\frac{Q^T(w_r-F(z^{-1}\left)y\right(k)-H\left(z^{-1}\right)\Delta u\left(k-1\right))}{Q^{T}Q+\lambda (1+{\beta }^2+\mathrm{...}{\beta }^{2(N_u-1)})}$式中：$Q=\left|\begin{array}{c}{g}_{1,0}\\ g_{2,1}+{\mathrm{\beta g}}_{1,0}\\ .\\ .\\ .\\ g_{i,N_u-1}+{\mathrm{\beta g}}_{i-1,N_u-2}+\mathrm{...}{\beta }^{N_u-1}{g}_{1,0}\\ .\\ .\\ .\\ g_{j,N-1}+{\mathrm{\beta g}}_{j-1,N-2}+\mathrm{...}{\beta }^{N_u-1}{g}_{j-N_u+1,N-N_u}\end{array}\right|$F(z^‑1)＝[F₁(z^‑1),…,F_N(z^‑1)]^TH(z^‑1)＝[H₁(z^‑1),…H_N(z^‑1)]^T其中，g_j,i为多项式G_j(z^‑1)的各项系数，由此可得实际输出控制量为：u(k)＝u(k‑1)+Δu(k|k)；E_j(z^‑1)与F_j(z^‑1)是在分别满足E_j(z^‑1)＝e_j,0+e_j,1z^‑1+…+e_j,j‑1z^‑(j‑1)与多项式形式的基础上，并在已知A(z^‑1)的前提下根据1＝E_j(z^‑1)A(z^‑1)Δ+z^‑jF_j(z^‑1)求解得到的多项式；其中n_a表示被控对象A(z^‑1)y(k)＝B(z^‑1)u(k‑1)+ξ(k)/Δ中的A(z^‑1)的阶次；j为期望优化步长，j∈{1,2,…,N}且N为优化时域终止时刻；e_j.0……e_j.j‑1表示E_j(z^‑1)对应求出的系数值，表示F_j(z^‑1)对应求出的系数值；然后，在求出E_j(z^‑1)与已知B(z^‑1)的前提下，根据公式E_j(z^‑1)B(z^‑1)＝G_j(z^‑1)+z^‑jH_j(z^‑1)求出G_j(z^‑1)＝g_j,0+g_j,1z^‑1+…+g_j,j‑1z^‑(j‑1)、形式的多项式G_j(z^‑1)、H_j(z^‑1)，其中g_j.0……g_j.j‑1、分别表示G_j(z^‑1)、H_j(z^‑1)求出的对应系数；n_b表示被控对象A(z^‑1)y(k)＝B(z^‑1)u(k‑1)+ξ(k)/Δ的B(z^‑1)的阶次；其中β为影响系数且β∈[0,1]；当j分别取1、2、…、N得到F_j(z^‑1)后，根据公式F(z^‑1)＝[F₁(z^‑1),…,F_N(z^‑1)]^T可获取F(z^‑1)；同理，当j分别取1、2、…、N得到H_j(z^‑1)后，可获取H(z^‑1)；在得出所有变量后可以得出k时刻控制器输出增量Δu(k|k)，Δu(k|k)表示k时刻的输入值与k‑1时刻的输入值之差。