【预测控制】动态矩阵控制(DMC)及Matlab仿真

这篇具有很好参考价值的文章主要介绍了【预测控制】动态矩阵控制(DMC)及Matlab仿真。希望对大家有所帮助。如果存在错误或未考虑完全的地方,请大家不吝赐教,您也可以点击"举报违法"按钮提交疑问。


相关文章:
【预测控制1】模型预测控制概论
【预测控制2】预测控制基本原理
【预测控制4】基于状态空间模型的预测控制


前言

动态矩阵控制(Dynamic Matrix Control,DMC)起源于20世纪70年代,是一种基于数学模型的先进控制算法,在化工、电力、冶金、制药等领域得到了广泛的应用。传统控制算法主要是基于PID控制器,但对于复杂的动态过程,PID控制器并不能提供良好的控制效果。相比之下,DMC控制算法的控制性能更好,可以适应更为复杂的过程控制。

特别说明本文中介绍的是SISO系统且无约束条件的DMC算法


一、基于阶跃响应的动态矩阵控制

DMC是一种基于阶跃响应的预测控制算法,因而适用于渐进稳定的线性系统。对于弱非线性对象,可以先对其线性化,再采用DMC,而对于不稳定对象,可以先用PID控制使其稳定,然后再使用DMC算法。
 简单的说,DMC就是通过阶跃响应模型预测系统未来输出,然后求解优化问题得到控制增量,将控制增量的第一个元素作用于系统。
DMC算法包括预测模型滚动优化反馈校正三个部分组成,接下来分别介绍三个部分。

1.预测模型

预测模型,就是通过系统的阶跃响应模型和线性系统的叠加性原理来预测对象在未来时刻的输出值。
在DMC算法中,是通过系统的单位阶跃响应模型来预测对象的未来输出。 首先需要得到预测对象的单位阶跃响应的采样值 a i = a ( i T ) , i = 1 , 2 , . . . a_i=a(iT),i=1,2,... ai=a(iT),i=1,2,...,其中 T T T为采样时间。对于渐进稳定的系统,阶跃响应会在某一时刻 t N = N T t_N=NT tN=NT趋于稳定平稳,在 T N TN TN时刻后的采样值 a N , a N + 1 , . . . . a_N,a_{N+1},.... aNaN+1,....都趋于稳态值,所以可以认为 a N = a N + 1 = a N + 2 = . . . . a_N=a_{N+1}=a_{N+2}=.... aN=aN+1=aN+2=....。这样,对象的动态信息就可以近似地用有限集合 { a 1 , a 2 , . . . . , a N − 1 , a N } \{a_1,a_2,....,a_{N-1},a_N\} {a1,a2,....,aN1,aN}加以描述,这个集合的参数构成了DMC的模型参数,向量 a = [ a 1 , a 2 , . . . . , a N − 1 , a N ] T \pmb a=[a_1,a_2,....,a_{N-1},a_N]^T aaa=[a1,a2,....,aN1,aN]T称为模型向量,也是系统的阶跃响应模型 N N N则称为建模时域
举个例子,假如100斤的小明想要把体重增长到120斤,他想制定一个计划,那么他肯定需要结合他之前每天吃多少饭涨了多少斤来计划以后吃多少。假如小明之前连续15天吃3碗饭,每一天都记录他的体重。在第10天的时候,小明的体重首次到达了105斤,第15天,小明的体重依然是105斤,我们可以认为小明吃3碗饭的稳态就是105斤,但这远远不符合他希望达到的120斤,小明就可以根据前10天体重变化为依据,来预测他今后吃多少饭在几天过后能达到120斤。而这前10天的体重就相当于模型向量 a \pmb a aaa,10天为建模时域 N N N

dmc控制算法matlab实现,矩阵,matlab,算法,自动化,线性代数

 下面我们就可以根据系统的模型向量 a \pmb a aaa和线性系统的比例和叠加性质来预测对象在未来时刻的输出值。
 首先,在k时刻,假定控制作用保持不变时对未来 N N N个时刻的输出具有初始预测值 y 0 ~ ( k + 1 ∣ k ) , i = 1 , . . . , N \tilde{y_0}(k+1|k),i=1,...,N y0~(k+1k),i=1,...,N,其中 k + i ∣ k k+i|k k+ik表示在k时刻对 k + i k+i k+i时刻的预测。
 当 k k k时刻输入 u u u的增量为 Δ u ( k ) \Delta u(k) Δu(k),在其作用下未来时刻输出预测值为 y ~ 1 ( k + i ∣ k ) = y 0 ~ ( k + 1 ∣ k ) + a i Δ u ( k ) . \tilde{y}_1(k+i|k)=\tilde{y_0}(k+1|k)+a_i\Delta u(k). y~1(k+ik)=y0~(k+1k)+aiΔu(k). k + 1 k+1 k+1时刻输入 u u u再变为 Δ u ( k + 1 ) \Delta u(k+1) Δu(k+1),则未来输入预测值为 y ~ 1 ( k + 1 ∣ k ) = y ~ 0 ( k + 1 ∣ k ) + a 1 Δ u ( k ) y ~ 1 ( k + 2 ∣ k ) = y ~ 0 ( k + 2 ∣ k ) + a 2 Δ u ( k ) + a 1 Δ u ( k + 1 ) ⋮ y ~ 1 ( k + N ∣ k ) = y ~ 0 ( k + N ∣ k ) + a N Δ u ( k ) + a N − 1 Δ u ( k + 1 ) \begin{aligned}\tilde{y}_1(k+1|k)&=\tilde{y}_0(k+1|k)+a_1\Delta u(k) \\\tilde{y}_1(k+2|k)&=\tilde{y}_0(k+2|k)+a_2\Delta u(k)+a_1\Delta u(k+1) \\\vdots \\\tilde{y}_1(k+N|k)&=\tilde{y}_0(k+N|k)+a_N\Delta u(k)+a_{N-1}\Delta u(k+1) \end{aligned} y~1(k+1k)y~1(k+2k)y~1(k+Nk)=y~0(k+1k)+a1Δu(k)=y~0(k+2k)+a2Δu(k)+a1Δu(k+1)=y~0(k+Nk)+aNΔu(k)+aN1Δu(k+1)
以此类推,连续M个输入 Δ u ( k ) , . . . , Δ u ( k + M − 1 ) \Delta u(k),...,\Delta u(k+M-1) Δu(k),...,Δu(k+M1)作用下,未来各时刻的输出值为 y ~ M ( k + i ∣ k ) = y ~ 0 ( k + i ∣ k ) + ∑ j = 1 m i n ( M , i ) a i − j + 1 Δ u ( k + j − 1 ) , i = 1 , . . . N (1-1) \tilde{y}_M(k+i|k)=\tilde{y}_0(k+i|k)+\sum_{j=1}^{min(M,i)}a_{i-j+1}\Delta u(k+j-1),i=1,...N\tag{1-1} y~M(k+ik)=y~0(k+ik)+j=1min(M,i)aij+1Δu(k+j1),i=1,...N(1-1)
其中, y y y的下标表示控制量变化的次数,即 y M y_M yM表示控制量 Δ u \Delta u Δu连续变化了 M M M次。根据上面的式子,只要我们知道对象的初始预测值 y ~ 0 ( k + i ∣ k ) \tilde{y}_0(k+i|k) y~0(k+ik),就可以根据系统的模型向量 a \pmb a aaa和控制增量来预测对象的未来输出。式(1-1)就被称为预测模型
dmc控制算法matlab实现,矩阵,matlab,算法,自动化,线性代数

2.滚动优化

 当我们建立好预测模型后,现在的问题就是我们要怎么计算 Δ u ( k ) , . . . , Δ u ( k + M − 1 ) \Delta u(k),...,\Delta u(k+M-1) Δu(k),...,Δu(k+M1)的最优值,来让系统的输出达到期望值,这就是滚动优化的核心问题。也就是说在每一时刻 k k k,要确定从该时刻起的 M M M个控制增量 Δ u ( k ) , . . . , Δ u ( k + M − 1 ) \Delta u(k),...,\Delta u(k+M-1) Δu(k),...,Δu(k+M1),使我们的系统的输出预测值 y ~ M ( k + i ∣ k ) \tilde{y}_M(k+i|k) y~M(k+ik)在未来P时刻尽可能接近给定的期望值 ω ( k + i ) , i = 1 , . . . , P \omega(k+i),i=1,...,P ω(k+i),i=1,...,P。这里的 M 、 P M、P MP分别被称为控制时域和优化时域( M ⩽ P ⩽ N M\leqslant P\leqslant N MPN
继续小明的增肥计划,小明根据自己的模型向量 a \pmb a aaa,计划5天分布增加饭量 Δ u ( 1 ) , . . . . , Δ u ( 5 ) \Delta u(1),....,\Delta u(5) Δu(1),....,Δu(5),然后预测接下7天的体重,让其尽可能的接近120斤。这里的5天相当于控制时域,7天相当于优化时域。那么如何求出 Δ u ( 1 ) , . . . . , Δ u ( 5 ) \Delta u(1),....,\Delta u(5) Δu(1),....,Δu(5)的最优值呢?
根据最优控制,我们可以写出在 k k k时刻的性能指标 m i n J ( k ) = ∑ i = 1 P q i [ ω ( k + i ) − y ~ M ( k + i ∣ k ) ] 2 minJ(k)=\sum_{i=1}^{P}q_i[\omega(k+i)-\tilde y_M(k+i|k)]^2 minJ(k)=i=1Pqi[ω(k+i)y~M(k+ik)]2除了要求让输出达到期望值以外,我们还要求不要让 Δ u ( k ) \Delta u(k) Δu(k)剧烈的变化,因此我们再在性能指标后加一项来软约束 Δ u ( k ) \Delta u(k) Δu(k),最终性能指标可取 m i n J ( k ) = ∑ i = 1 P q i [ ω ( k + i ) − y ~ M ( k + i ∣ k ) ] 2 + ∑ j = 1 M r j Δ u 2 ( k + j − 1 ) (2-1) minJ(k)=\sum_{i=1}^{P}q_i[\omega(k+i)-\tilde y_M(k+i|k)]^2+\sum_{j=1}^Mr_j\Delta u^2(k+j-1)\tag{2-1} minJ(k)=i=1Pqi[ω(k+i)y~M(k+ik)]2+j=1MrjΔu2(k+j1)(2-1)其中, q i 、 r j q_i、r_j qirj是权系数,分别表示对跟踪误差及控制量变化的抑制程度, q i q_i qi值越大,表明我们期望对应的控制输出越接近期望值; r j r_j rj的值越大,表明我们期望对应的控制动作变化越小。
 上述问题就是动态模型(1-1)约束下,以 Δ u M ( k ) = [ Δ u ( k ) , . . . , Δ u ( k + M − 1 ) ] T \Delta u_M(k)=[\Delta u(k),...,\Delta u(k+M-1)]^T ΔuM(k)=[Δu(k),...,Δu(k+M1)]T为优化变量,使性能指标(2-1)最小的优化问题。
 现在我们回过头来将式(1-1)写为向量的形式,其中要强调的是,优化时域为P,表示我们要使预测输出 y ~ M ( k + 1 ∣ k ) , . . . , y ~ M ( k + P ∣ k ) \tilde{y}_M(k+1|k),...,\tilde{y}_M(k+P|k) y~M(k+1k),...,y~M(k+Pk)尽可能的接近我们的期望输出;控制时域为M,表示我们有连续M个输入 Δ u ( k ) , . . . , Δ u ( k + M − 1 ) \Delta u(k),...,\Delta u(k+M-1) Δu(k),...,Δu(k+M1)作用,因此可以将式(1-1)展开写为 y ~ M ( k + 1 ∣ k ) = y ~ 0 ( k + 1 ∣ k ) + a 1 Δ u ( k ) y ~ M ( k + 2 ∣ k ) = y ~ 0 ( k + 2 ∣ k ) + a 2 Δ u ( k ) + a 1 Δ u ( k + 1 ) y ~ M ( k + 3 ∣ k ) = y ~ 0 ( k + 3 ∣ k ) + a 3 Δ u ( k ) + a 2 Δ u ( k + 1 ) + a 1 Δ u ( k + 2 ) . . . . . . y ~ M ( k + P ∣ k ) = y ~ 0 ( k + P ∣ k ) + a P Δ u ( k ) + a P − 1 Δ u ( k + 1 ) + . . . + a P − M + 1 Δ u ( k + M − 1 ) \begin{aligned}\tilde{y}_M(k+1|k)&=\tilde{y}_0(k+1|k)+a_1\Delta u(k) \\\tilde{y}_M(k+2|k)&=\tilde{y}_0(k+2|k)+a_2\Delta u(k)+a_1\Delta u(k+1) \\\tilde{y}_M(k+3|k)&=\tilde{y}_0(k+3|k)+a_3\Delta u(k)+a_2\Delta u(k+1)+a_1\Delta u(k+2) \\&...... \\\tilde{y}_M(k+P|k)&=\tilde{y}_0(k+P|k)+a_P\Delta u(k)+a_{P-1}\Delta u(k+1)+...+a_{P-M+1}\Delta u(k+M-1) \end{aligned} y~M(k+1k)y~M(k+2k)y~M(k+3k)y~M(k+Pk)=y~0(k+1k)+a1Δu(k)=y~0(k+2k)+a2Δu(k)+a1Δu(k+1)=y~0(k+3k)+a3Δu(k)+a2Δu(k+1)+a1Δu(k+2)......=y~0(k+Pk)+aPΔu(k)+aP1Δu(k+1)+...+aPM+1Δu(k+M1)现在定义 Y ~ P M ( k ) = [ y ~ M ( k + 1 ∣ k ) y ~ M ( k + 2 ∣ k ) . . . y ~ M ( k + P ∣ k ) ] , Y ~ P 0 ( k ) = [ y ~ 0 ( k + 1 ∣ k ) y ~ 0 ( k + 2 ∣ k ) . . . y ~ 0 ( k + P ∣ k ) ] , Δ U M ( k ) = [ Δ u ( k ) Δ u ( k + 1 ) . . . Δ u ( k + M − 1 ) ] , A = [ a 1 0 . . . 0 a 2 a 1 . . . 0 ⋮ ⋮ ⋱ ⋮ a M a M − 1 . . . a 1 ⋮ ⋮ ⋮ a P a P − 1 . . . a P − M + 1 ] . \tilde Y_{PM}(k)=\begin{bmatrix}\tilde{y}_M(k+1|k)\\\tilde{y}_M(k+2|k)\\...\\\tilde{y}_M(k+P|k)\end{bmatrix},\tilde Y_{P0}(k)=\begin{bmatrix}\tilde{y}_0(k+1|k)\\\tilde{y}_0(k+2|k)\\...\\\tilde{y}_0(k+P|k)\end{bmatrix}, \\\Delta U_M(k)=\begin{bmatrix}\Delta u(k)\\\Delta u(k+1)\\...\\\Delta u(k+M-1)\end{bmatrix},A=\begin{bmatrix} a_1&0&...&0\\a_2&a_1&...&0\\\vdots&\vdots&\ddots&\vdots\\a_M&a_{M-1}&...&a_1\\\vdots&\vdots&&\vdots\\a_P&a_{P-1}&...&a_{P-M+1}\end{bmatrix}. Y~PM(k)=y~M(k+1k)y~M(k+2k)...y~M(k+Pk),Y~P0(k)=y~0(k+1k)y~0(k+2k)...y~0(k+Pk),ΔUM(k)=Δu(k)Δu(k+1)...Δu(k+M1),A=a1a2aMaP0a1aM1aP1............00a1aPM+1.其中 A \pmb A AAA被称为动态矩阵
因此,可以将式(1-1)写为向量形式 Y ~ P M = Y ~ P 0 + A Δ U M ( k ) . (2-2) \tilde Y_{PM}=\tilde Y_{P0}+A\Delta U_M(k).\tag{2-2} Y~PM=Y~P0+AΔUM(k).(2-2)同样,性能指标(2-1)也可以写为向量形式 m i n J ( k ) = ∥ ω P ( k ) − Y ~ P M ( k ) ∥ Q 2 + ∥ Δ U M ( k ) ∥ R 2 (2-3) minJ(k)=\begin{Vmatrix}\pmb \omega_P(k)-\tilde{Y}_{PM}(k)\end{Vmatrix}^2_{\pmb Q}+\begin{Vmatrix}\Delta U_M(k)\end{Vmatrix}^2_{\pmb R}\tag{2-3} minJ(k)=ωωωP(k)Y~PM(k)QQQ2+ΔUM(k)RRR2(2-3)其中, ω P ( k ) = [ ω P ( k + 1 ) , . . . , ω P ( k + P ) ] T , Q = d i a g ( q 1 , . . . , q P ) , R = d i a g ( r 1 , . . . , r M ) . \pmb \omega_P(k)=[\omega_P(k+1),...,\omega_P(k+P)]^T, \\\pmb Q=diag(q_1,...,q_P), \\\pmb R=diag(r_1,...,r_M). ωωωP(k)=[ωP(k+1),...,ωP(k+P)]T,QQQ=diag(q1,...,qP),RRR=diag(r1,...,rM).由权系数构成的对角矩阵 Q 、 R Q、R QR分别被称为误差权矩阵控制权矩阵
将式(2-2)带入式(2-3)可得, m i n J ( k ) = ∥ ω P ( k ) − Y ~ P 0 − A Δ U M ( k ) ∥ Q 2 + ∥ Δ U M ( k ) ∥ R 2 minJ(k)=\begin{Vmatrix}\pmb \omega_P(k)-\tilde Y_{P0}-\pmb A\Delta U_M(k)\end{Vmatrix}^2_{\pmb Q}+\begin{Vmatrix}\Delta U_M(k)\end{Vmatrix}^2_{\pmb R} minJ(k)=ωωωP(k)Y~P0AAAΔUM(k)QQQ2+ΔUM(k)RRR2为了求 Δ U M ( k ) \Delta U_M(k) ΔUM(k)满足上式的最小值,对上式关于 Δ U M ( k ) \Delta U_M(k) ΔUM(k)求导并令其等于0可以得到, Δ U M ( k ) = ( A T Q A + R ) − 1 A T Q [ ω P ( k ) − Y ~ P 0 ( k ) ] \Delta U_M(k)=(\pmb{A^TQA}+\pmb R)^{-1}\pmb{A^TQ}[\pmb \omega_P(k)-\tilde Y_{P0}(k)] ΔUM(k)=(ATQAATQAATQA+RRR)1ATQATQATQ[ωωωP(k)Y~P0(k)]
但我们不把它们全部实施,我们只取 Δ U M ( k ) \Delta U_M(k) ΔUM(k)的首元素来作用于对象,可以表示为, Δ u ( k ) = c T Δ U M ( k ) = c T ( A T Q A + R ) − 1 A T Q [ ω P ( k ) − Y ~ P 0 ( k ) ] \Delta u(k)=\pmb c^T\Delta U_M(k)=\pmb c^T(\pmb{A^TQA}+\pmb R)^{-1}\pmb{A^TQ}[\pmb \omega_P(k)-\tilde Y_{P0}(k)] Δu(k)=cccTΔUM(k)=cccT(ATQAATQAATQA+RRR)1ATQATQATQ[ωωωP(k)Y~P0(k)]其中, c T = [ 1 0 . . . 0 ] \pmb c^T=\begin{bmatrix}1&0&...&0\end{bmatrix} cccT=[10...0]。然后定义控制增益 d T = c T ( A T Q A + R ) − 1 A T Q \pmb d^T=\pmb c^T(\pmb{A^TQA}+\pmb R)^{-1}\pmb{A^TQ} dddT=cccT(ATQAATQAATQA+RRR)1ATQATQATQ最终可以得到 Δ u ( k ) = d T [ ω P ( k ) − Y ~ P 0 ( k ) ] (2-4) \Delta u(k)=\pmb d^T[\pmb \omega_P(k)-\tilde Y_{P0}(k)]\tag{2-4} Δu(k)=dddT[ωωωP(k)Y~P0(k)](2-4)在求出控制增量 Δ u ( k ) \Delta u(k) Δu(k)后,当前k时刻的实际控制量为 u ( k ) = u ( k − 1 ) + Δ u ( k ) (2-5) u(k)=u(k-1)+\Delta u(k)\tag{2-5} u(k)=u(k1)+Δu(k)(2-5)然后将 u ( k ) u(k) u(k)作用于对象,到下一时刻,又以 k + 1 k+1 k+1取代 k k k提出同样的优化问题求出 Δ u ( k + 1 ) \Delta u(k+1) Δu(k+1),得到 u ( k + 1 ) u(k+1) u(k+1)作用于对象。如此滚动的进行,这就是滚动优化的含义。

3.反馈校正

 当我们在k时刻将 u ( k ) u(k) u(k)作用于对象之后,相当于我们在 k k k时刻给对象的输入端加上了一个幅值为 Δ u ( k ) \Delta u(k) Δu(k)的阶跃信号,根据预测模型(1-1),我们可以得到在接下来 N N N个时刻的输出 y ~ N 1 ( k + 1 ∣ k ) , . . . , y ~ N 1 ( k + N ∣ k ) \tilde y_{N1}(k+1|k),...,\tilde y_{N1}(k+N|k) y~N1(k+1k),...,y~N1(k+Nk),即 Y ~ N 1 ( k ) = Y ~ N 0 ( k ) + a Δ u ( k ) , (3-1) \tilde Y_{N1}(k)=\tilde Y_{N0}(k)+\pmb a\Delta u(k),\tag{3-1} Y~N1(k)=Y~N0(k)+aaaΔu(k),(3-1)其中, Y ~ N 1 ( k ) = [ y ~ N 1 ( k + 1 ∣ k ) y ~ N 1 ( k + 2 ∣ k ) . . . y ~ N 1 ( k + N ∣ k ) ] , Y ~ N 0 ( k ) = [ y ~ N 0 ( k + 1 ∣ k ) y ~ N 0 ( k + 2 ∣ k ) . . . y ~ N 0 ( k + N ∣ k ) ] , a = [ a 1 , a 2 , . . . . , a N − 1 , a N ] T . \tilde Y_{N1}(k)=\begin{bmatrix}\tilde{y}_{N1}(k+1|k)\\\tilde{y}_{N1}(k+2|k)\\...\\\tilde{y}_{N1}(k+N|k)\end{bmatrix},\tilde Y_{N0}(k)=\begin{bmatrix}\tilde{y}_{N0}(k+1|k)\\\tilde{y}_{N0}(k+2|k)\\...\\\tilde{y}_{N0}(k+N|k)\end{bmatrix}, \\\pmb a=[a_1,a_2,....,a_{N-1},a_N]^T. Y~N1(k)=y~N1(k+1k)y~N1(k+2k)...y~N1(k+Nk),Y~N0(k)=y~N0(k+1k)y~N0(k+2k)...y~N0(k+Nk),aaa=[a1,a2,....,aN1,aN]T. y ~ N 0 \tilde y_{N0} y~N0表示控制作用保持不变时接下来N个时刻的预测输出, y ~ N 1 \tilde y_{N1} y~N1则表示在当前时刻施加了一个控制增量 Δ u ( k ) \Delta u(k) Δu(k)后接下来N个时刻的预测输出。
 但是在实际情况中,由于环境因素等的干扰,系统的实际输出可能会比我们预测的输出偏小或偏大,就像小明他有可能会感冒生病,导致他吃相同的饭量会比以前增长的体重少。这种情况我们必须要考虑,因为我们下一步的优化是建立在这一次预测输出之上,如果我们这一次预测输出出现了偏差,那么我们下一步优化偏差值可能更大。那么我们要如何减小这个偏差呢?在DMC算法中,我们首先将我们预测到的下一时刻的输出 y ~ 1 ( k + 1 ∣ k ) \tilde y_1(k+1|k) y~1(k+1k)和系统的实际输出 y ( k + 1 ) y(k+1) y(k+1)相比较,得到偏差值 e ( k + 1 ) = y ( k + 1 ) − y ~ 1 ( k + 1 ∣ k ) (3-2) e(k+1)=y(k+1)-\tilde y_1(k+1|k)\tag{3-2} e(k+1)=y(k+1)y~1(k+1k)(3-2)在上面的式子中,系统的实际输出 y ( k + 1 ) y(k+1) y(k+1)是已知的,我们的预测输出 y ~ 1 ( k + 1 ∣ k ) \tilde y_1(k+1|k) y~1(k+1k)也是已知的,说明偏差 e ( k + 1 ) e(k+1) e(k+1)也是已知的,DMC就想,既然偏差都已知了,那我们直接把这个偏差 e ( k + 1 ) e(k+1) e(k+1)补偿到我们的接下来的预测输出 Y ~ N 1 ( k ) \tilde Y_{N1}(k) Y~N1(k)不就可以减小偏差了吗,因此可以得到 Y ~ c o r ( k + 1 ) = Y ~ N 1 ( k ) + h e ( k + 1 ) , (3-3) \tilde Y_{cor}(k+1)=\tilde Y_{N1}(k)+\pmb he(k+1)\tag{3-3}, Y~cor(k+1)=Y~N1(k)+hhhe(k+1),(3-3)其中, Y ~ c o r ( k + 1 ) = [ y ~ c o r ( k + 1 ∣ k + 1 ) y ~ c o r ( k + 2 ∣ k + 1 ) . . . y ~ c o r ( k + N ∣ k + 1 ) ] \tilde Y_{cor}(k+1)=\begin{bmatrix}\tilde{y}_{cor}(k+1|k+1)\\\tilde{y}_{cor}(k+2|k+1)\\...\\\tilde{y}_{cor}(k+N|k+1)\end{bmatrix} Y~cor(k+1)=y~cor(k+1k+1)y~cor(k+2k+1)...y~cor(k+Nk+1)为校正后的输出预测量,有权系数组成的 N N N维向量 h = [ h 1 , . . . h N ] T \pmb h=[h_1,...h_N]^T hhh=[h1,...hN]T称为校正量。
 以上就是反馈校正的主要思想。


 现在我们都知道,预测模型是需要我们提前准备好的,而滚动优化反馈校正是每一次循环都要进行的,那我们如何将滚动优化和反馈校正两部分形成"闭环"呢?回顾公式,我们需要把 k k k时刻经过修正的预测输出 Y ~ c o r ( k + 1 ) \tilde Y_{cor}(k+1) Y~cor(k+1)作为 k + 1 k+1 k+1时刻的初始预测值,也就是初值 Y ~ N 0 ( k + 1 ) \tilde Y_{N0}(k+1) Y~N0(k+1)。但我们不能直接令 Y ~ c o r ( k + 1 ) = Y ~ N 0 ( k + 1 ) \tilde Y_{cor}(k+1)=\tilde Y_{N0}(k+1) Y~cor(k+1)=Y~N0(k+1)。因为对于 k k k时刻而言,我们得到的预测输出 Y ~ c o r ( k + 1 ) = [ y ~ c o r ( k + 1 ∣ k + 1 ) , . . . , y ~ c o r ( k + N ∣ k + 1 ) ] \tilde Y_{cor}(k+1)=[\tilde y_{cor}(k+1|k+1),...,\tilde y_{cor}(k+N|k+1)] Y~cor(k+1)=[y~cor(k+1k+1),...,y~cor(k+Nk+1)],而在 k + 1 k+1 k+1时刻的初始预测值为 Y ~ N 0 ( k + 1 ) = [ y ~ N 0 ( k + 2 ∣ k + 1 ) , . . . , y ~ N 0 ( k + N + 1 ∣ k + 1 ) ] \tilde Y_{N0}(k+1)=[\tilde{y}_{N0}(k+2|k+1),...,\tilde{y}_{N0}(k+N+1|k+1)] Y~N0(k+1)=[y~N0(k+2k+1),...,y~N0(k+N+1k+1)],所以我们令 y ~ N 0 ( k + 2 ∣ k + 1 ) = y ~ c o r ( k + 2 ∣ k + 1 ) y ~ N 0 ( k + 3 ∣ k + 1 ) = y ~ c o r ( k + 3 ∣ k + 1 ) . . . y ~ N 0 ( k + N ∣ k + 1 ) = y ~ c o r ( k + N ∣ k + 1 ) \tilde{y}_{N0}(k+2|k+1)=\tilde y_{cor}(k+2|k+1)\\\tilde{y}_{N0}(k+3|k+1)=\tilde y_{cor}(k+3|k+1)\\...\\\tilde{y}_{N0}(k+N|k+1)=\tilde y_{cor}(k+N|k+1) y~N0(k+2k+1)=y~cor(k+2k+1)y~N0(k+3k+1)=y~cor(k+3k+1)...y~N0(k+Nk+1)=y~cor(k+Nk+1)由于 Y ~ c o r ( k + 1 ) \tilde Y_{cor}(k+1) Y~cor(k+1)中没有 y ~ c o r ( k + N + 1 ∣ k + 1 ) \tilde y_{cor}(k+N+1|k+1) y~cor(k+N+1k+1)这一项,那么我们用 y ~ c o r ( k + N ∣ k + 1 ) \tilde y_{cor}(k+N|k+1) y~cor(k+Nk+1)取近似。最终可以表示为 Y ~ N 0 ( k + 1 ) = S Y ~ c o r ( k + 1 ) \tilde Y_{N0}(k+1)=\pmb S\tilde Y_{cor}(k+1) Y~N0(k+1)=SSSY~cor(k+1)其中, S = [ 0 1 0 . . . 0 0 0 1 . . . 0 0 0 0 ⋱ 0 ⋮ ⋮ ⋮ 0 1 0 . . . 0 1 ] S=\begin{bmatrix}0&1&0&...&0\\0&0&1&...&0\\0&0&0&\ddots &0\\\vdots&\vdots&\vdots&0&1\\0&...&&0&1 \end{bmatrix} S=0000100...010......0000011有了 Y ~ N 0 ( k + 1 ) \tilde Y_{N0}(k+1) Y~N0(k+1),就可以取其前P项得到 Y ~ P 0 ( k + 1 ) \tilde Y_{P0}(k+1) Y~P0(k+1),然后可以根据式子(2-4)计算出 Δ u ( k + 1 ) \Delta u(k+1) Δu(k+1),整个控制过程就是滚动优化结合反馈校正反复在线进行的。

二、Matlab仿真示例

1.算法结构

根据上面整个DMC算法的流程,我们可以得到它的算法结构,如下图所示。图中粗箭头表示向量流,细箭头表示数值流。
dmc控制算法matlab实现,矩阵,matlab,算法,自动化,线性代数

2.计算流程

dmc控制算法matlab实现,矩阵,matlab,算法,自动化,线性代数

流程图中的算式依次对应于式(3-2)、式(3-3)式(2-4)、式(2-5)、式(3-1)

3.仿真示例

例: 考虑如下SISO系统,其传递函数为 G ( s ) = 3 4 s + 1 G(s)=\frac{3}{4s+1} G(s)=4s+13

程序代码:

clear all;
close all;

P=6;M=2;N=25;%优化时域、控制时域和建模时域
ts=1;%采样时间
k_step=50;
times=zeros(1,k_step);

den=[4,1];
num=[3];
sys=tf(num,den);%建立传递函数
sysd=c2d(sys,ts);
[numd,dend]=tfdata(sysd,'v');
yk_1=0;uk_1=0;uk=0;
[a,t]=step(num,den,1:ts:N*ts);%获得模型向量a


S=diag(ones(1,N-1),1);S(N,N)=1;%移位矩阵S
yr=3;%期望输出

dum=zeros(M,1);%优化变量
du=0;%控制增量
yn0=zeros(N,1);%前一时刻预测输出
yn1=zeros(N,1);%当前预测输出
yp0=zeros(P,1);%优化输出
us=zeros(1,k_step);%实际控制量
ys=zeros(1,k_step);%实际输出


R=eye(M);Q=eye(P);%控制权矩阵和误差权矩阵
h=0.5*ones(N,1);
h(1)=1;%校正系数矩阵

A=zeros(P,M);
for i=1:M       
    A(i:P,i)=a(1:P-i+1);
end  %建立动态矩阵A
c=zeros(M,1);c(1)=1;
d=c'*inv(A'*Q*A+R)*A'*Q;d=d';%控制增益

%第一步预测
ys(1)=0;
yk_1=ys(1);
ek=ys(1)-yn1(1);
yn0=S*(yn1+h*ek);%预测值校正
yp0=yn0(1:P);
du=d'*(yr-yp0);%计算这一时刻控制增量
us(1)=0+du;%施加控制增量,得到实际控制量
uk=us(1);
yn1=yn0+a(1:N)*du;%得到施加控制增量后的预测输出

%滚动优化
for k=2:k_step
    times(k)=k;
    ys(k)=-dend(2)*yk_1+numd(2)*uk;%得到系统的实际输出
    yk_1=ys(k);
    ek=ys(k)-yn1(1);%计算偏差,公式(3-2)
    yn0=S*(yn1+h*ek);%预测值校正,公式(3-3)
    yp0=yn0(1:P);
    du=d'*(yr-yp0);%计算这一时刻控制增量,公式(2-4)
    us(k)=us(k-1)+du;%施加控制增量,得到实际控制量,公式(2-5)
    uk=us(k);
    yn1=yn0+a(1:N)*du;%得到施加控制增量后的预测输出,公式(3-1)
end

figure(1)
plot(times,ys');
title('输出曲线')
figure(2)
stairs(times,us');
title('控制量变化曲线')

运行结果:
dmc控制算法matlab实现,矩阵,matlab,算法,自动化,线性代数dmc控制算法matlab实现,矩阵,matlab,算法,自动化,线性代数


总结

DMC算法(Dynamic Matrix Control)是一种经典的模型预测控制算法,适用于线性、纯时滞、开环渐进稳定的非最小相位系统,其优缺点如下:
优点:

  1. 较好的控制性能:DMC算法能够很好地解决传统PID控制算法难以处理的复杂控制问题,能够提供更好的控制性能。
  2. 可以处理多变量控制问题:DMC算法可以对多个变量同时进行控制,能够处理多变量控制问题。
  3. 具有鲁棒性:DMC算法对于模型不确定性、测量噪声等影响具有很强的鲁棒性。

缺点:文章来源地址https://www.toymoban.com/news/detail-763515.html

  1. 计算量大:DMC算法计算比较复杂,需要大量的计算资源,处理速度较慢。
  2. 模型依赖性强:DMC算法的控制性能依赖于模型的精度,如果模型不准确,控制效果会受到影响。
  3. 参数调节困难:DMC算法需要调节多个参数,参数调节比较困难,需要一定的经验和技巧。

到了这里,关于【预测控制】动态矩阵控制(DMC)及Matlab仿真的文章就介绍完了。如果您还想了解更多内容,请在右上角搜索TOY模板网以前的文章或继续浏览下面的相关文章,希望大家以后多多支持TOY模板网!

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处: 如若内容造成侵权/违法违规/事实不符,请点击违法举报进行投诉反馈,一经查实,立即删除!

领支付宝红包 赞助服务器费用

相关文章

  • 基于预测控制模型的自适应巡航控制仿真与机器人实现(Matlab代码实现)

         目录 💥1 概述 📚2 运行结果 🎉3 参考文献 👨‍💻4 Matlab代码 自适应巡航控制技术为目前由于汽车保有量不断增长而带来的行车安全、驾驶舒适性及交通拥堵等问题提供了一条有效的解决途径,因此本文通过理论分析、仿真验证及实车实验对自适应巡航控制中的若干

    2024年02月16日
    浏览(47)
  • 基于龙格库塔算法的SIR病毒扩散预测matlab仿真

    目录 1.程序功能描述 2.测试软件版本以及运行结果展示 3.核心程序 4.本算法原理 5.完整程序       基于龙格库塔算法的SIR病毒扩散预测,通过龙格库塔算法求解传染病模型的微分方程。输出易受感染人群数量曲线,感染人群数量曲线,康复人群数量曲线。 MATLAB2022a版本运行  

    2024年01月23日
    浏览(33)
  • 基于HMM隐马尔可夫模型的金融数据预测算法matlab仿真

    目录 1.程序功能描述 2.测试软件版本以及运行结果展示 3.核心程序 4.本算法原理 5.完整程序         基于HMM隐马尔可夫模型的金融数据预测算法.程序实现HMM模型的训练,使用训练后的模型进行预测。 MATLAB2022A版本运行        隐马尔可夫模型(Hidden Markov Model, HMM)是一种概

    2024年04月25日
    浏览(33)
  • 车辆行驶控制运动学模型的matlab建模与仿真,仿真输出车辆动态行驶过程

    目录 1.课题概述 2.系统仿真结果 3.核心程序与模型 4.系统原理简介 4.1 基本假设 4.2 运动学方程 5.完整工程文件 车辆行驶控制运动学模型的matlab建模与仿真,仿真输出车辆动态行驶过程. 版本:MATLAB2022a        车辆运动学模型从几何学的角度研究车辆的运动规律。包括车辆的空

    2024年01月20日
    浏览(64)
  • 实验四、最少拍控制算法matlab仿真实验

    1.掌握最少拍有纹波、无纹波系统的设计方法; 2.学会最少拍控制系统的分析方法; 3.了解输入信号对最少拍控制系统的影响及改进措施 MATLAB 软件( 2022a)   1、写出广义被控对象的脉冲传递函数G(z)。并求出广义被控对象的差分方程形式。 2、画出未加控制器时系统的单

    2024年02月04日
    浏览(37)
  • 毕业设计-基于 PID 控制算法仿真算法研究- Matlab

    目录 前言 课题背景和意义 实现技术思路 一、 基本原理  二、无超调 PID 控制器的设计 三、无超调 PID 设计的验证 代码 实现效果图样例 最后     📅大四是整个大学期间最忙碌的时光,一边要忙着备考或实习为毕业后面临的就业升学做准备,一边要为毕业设计耗费大量精力。

    2024年02月06日
    浏览(51)
  • 水下机器人双机械手系统动态建模与控制仿真(Matlab代码实现)

    ​       目录 💥1 概述 📚2 运行结果 🎉3 参考文献 👨‍💻4 Matlab代码 水下机器人-机械手系统(Underwater vehicle-manipulator systems, UVMS)可以完成除观测之外的水下采样、抓取、操作等任务,在海洋科学考察、海洋工程等领域得到广泛应用。通过对近年来国内外UVMS的研究现状

    2024年02月08日
    浏览(53)
  • 浅谈BP神经网络PID控制算法及matlab仿真

    本文是对BP神经网络PID控制算法的数学描述及仿真实验,若有错误之处,欢迎指正! 老规矩不废话,直接上链接 BP神经网络维基百科 BP神经网络是人工神经网络中的一种常用结构,其由输入层(input)-隐含层(hidding)-输出层三层构成(output)。 上图中, B 1 B1 B 1 是输入层, B 2 B2 B

    2024年02月05日
    浏览(38)
  • 轨迹规划 | 图解动态窗口算法DWA(附ROS C++/Python/Matlab仿真)

    🔥附C++/Python/Matlab全套代码🔥课程设计、毕业设计、创新竞赛必备!详细介绍全局规划(图搜索、采样法、智能算法等);局部规划(DWA、APF等);曲线优化(贝塞尔曲线、B样条曲线等)。 🚀详情:图解自动驾驶中的运动规划(Moti

    2024年02月08日
    浏览(52)
  • 路径规划 | 图解动态A*(D*)算法(附ROS C++/Python/Matlab仿真)

    🔥附C++/Python/Matlab全套代码🔥课程设计、毕业设计、创新竞赛必备!详细介绍全局规划(图搜索、采样法、智能算法等);局部规划(DWA、APF等);曲线优化(贝塞尔曲线、B样条曲线等)。 🚀详情:图解自动驾驶中的运动规划(Motion Planning),附几十种规划算法

    2024年02月03日
    浏览(65)

觉得文章有用就打赏一下文章作者

支付宝扫一扫打赏

博客赞助

微信扫一扫打赏

请作者喝杯咖啡吧~博客赞助

支付宝扫一扫领取红包,优惠每天领

二维码1

领取红包

二维码2

领红包