最小方差自校正调节器设计

1、前言

自校正控制适用于结构已知但参数未知而恒定的随机系统，也适用于结构已知但参数缓慢变化的随机系统。其能完成调节任务，也能完成伺服跟踪任务。

自校正控制涉及到以下几个方面：

（1）被控对象参数估计：由于系统参数未知或者缓慢变化，无法直接进行控制，因此需要根据系统的输入输出，估计系统的参数；
（2）自校正控制系统性能指标：根据控制系统的性质和对控制系统的要求，有不同的性能指标形式，最普遍的是最小方差控制和极点配置控制。对于最小方差控制的目标函数为误差的二次型；极点配置控制的性能指标不是用目标函数的形式表示的，而是把预期的闭环系统的行为用一组期望传递函数的零极点的位置加以规定的，其目标就是想要实际的闭环系统的零极点收敛于这一组期望的零极点。

最小方差控制包含最小方差自校正调节器和最小方差自校正控制器，本文先从最小方差自校正调节器进行介绍。

2、系统介绍

最小方差自校正调节器结构框图如下：

其包含3个模块：控制对象、参数估计器和调节器。参数估计器根据控制对象的输出$y$和输入$u$的观测序列来估计被控对象的参数$\theta$，将得到的参数估计值$\hat{\theta }$送到调节器，调节器按照设计好的性能指标以及参数估计值$\hat{\theta }$对控制器的参数进行调整，实现系统运行的性能指标达到最优或者接近最优状态。

在最小方差自校正调节器中，用在线辨识方法估计系统的参数，以输出$y(t)$的方差最小作为调节指标，在自校正调节过程中，使系统的输出方差达到最小。

3、调节器设计

3.1 最优控制

被控对象的差分方程为
$$\begin{array}{c}A(q^{-1})y(t)=q^{-d}B(q^{-1})u(t)+C(q^{-1})w(t)\end{array}\text{\hspace{1em}(1)}$$
其中，
$$\begin{array}{c}\begin{array}{rl}A(q^{-1})& =1+a_1{q}^{-1}+a_2{q}^{-2}+\cdots +a_{n_a}{q}^{-n_a}\\ B(q^{-1})& =b_0+b_1{q}^{-1}+\cdots +b_{n_b}{q}^{-n_b}\\ C(q^{-1})& =c_0+c_1{q}^{-1}+\cdots +c_{n_c}{q}^{-n_c}\end{array}\end{array}\text{\hspace{1em}(2)}$$
$y(t)$为控制对象的输出；$u(t)$为其输入；$w(t)$为均值为0，方差为${\sigma }^2$的白噪声序列；$d$为系统的时延。

由于系统存在$d$个周期的时延，即当前的控制作用$u(t)$要滞后$d$的采样周期才能对系统的输出产生影响，因此要使系统输出方差最小，需要对输出量提前$d$步进行预测，根据预测值来调整控制作用$u(t)$，以补偿随机扰动在$(t+d)$时刻对输出的影响。

最小方差自校正调节器二次型指标函数为
$$\begin{array}{c}J=E[y^2(t+d)]=J_{min}\end{array}\text{\hspace{1em}(3)}$$
根据$eq(1)$可得
$$\begin{array}{c}A(q^{-1})y(t)=q^{-d}B(q^{-1})u(t)+C(q^{-1})w(t)\\ \Rightarrow A(q^{-1})y(t+d)=B(q^{-1})u(t)+C(q^{-1})w(t+d)\\ \Rightarrow y(t+d)=\frac{B(q^{-1})}{A(q^{-1})}u(t)+\frac{C(q^{-1})}{A(q^{-1})}w(t+d)\end{array}\text{\hspace{1em}(4)}$$

要使最小方差自校正调节器的解存在，需要满足下列假设：

如果 $B(q^{-1}))$的零点全在 $q^{-1}$复平面单位圆外，则称该系统为最小相位系统，否则为非最小相位系统，有时称 $B(q^{-1}))$的零点全在 $q^{-1}$复平面单位圆外的系统为逆稳定系统，否则为逆不稳定系统。之所以要求 $B(q^{-1}))$和$C(q^{-1}))$的零点全在单位圆外, 与闭环系统的稳定性有关，后面在稳定性这一节进行分析。

由$eq(4)$可得
$$\begin{array}{c}y(t+d)=\frac{B(q^{-1})}{A(q^{-1})}u(t)+\frac{C(q^{-1})}{A(q^{-1})}w(t+d)\end{array}$$
将$\frac{C(q^{-1})}{A(q^{-1})}$分解为
$$\begin{array}{c}\frac{C(q^{-1})}{A(q^{-1})}=F(q^{-1})+\frac{q^{-d}G(q^{-1})}{A(q^{-1})}\end{array}\text{\hspace{1em}(5)}$$
其中，
$$\begin{array}{c}\begin{array}{rl}F(q^{-1})& =1+f_1{q}^{-1}+\cdots +f_{d-1}{q}^{-d+1}\\ G(q^{-1})& =g_0+g_1{q}^{-1}+\cdots +g_{n_a-1}{q}^{-n_a+1}\end{array}\end{array}\text{\hspace{1em}(6)}$$
将$eq(6)$代入$eq(4)$可得
$$\begin{array}{c}y(t+d)=\frac{B(q^{-1})}{A(q^{-1})}u(t)+F(q^{-1})w(t+d)+\frac{G(q^{-1})}{A{(}^{-1})}w(t)\end{array}\text{\hspace{1em}(7)}$$
由$eq(1)$可得
$$\begin{array}{c}w(t)=\frac{A(q^{-1})}{C(q^{-1})}y(t)-\frac{B(q^{-1})}{C(q^{-1})}{q}^{-d}u(t)\end{array}\text{\hspace{1em}(8)}$$
将$eq(8)$代入$eq(7)$可得
$$\begin{array}{c}y(t+d)=\frac{B(q^{-1})}{A(q^{-1})}u(t)+F(q^{-1})w(t+d)+\frac{G(q^{-1})}{A{(}^{-1})}[\frac{A(q^{-1})}{C(q^{-1})}y(t)-\frac{B(q^{-1})}{C(q^{-1})}{q}^{-d}u(t)]\\ \\ =\frac{B(q^{-1})}{A(q^{-1})}u(t)+F(q^{-1})w(t+d)+\frac{G(q^{-1})}{C{(}^{-1})}y(t)-\frac{B(q^{-1})G(q^{-1})}{A{(}^{-1})C(q^{-1})}{q}^{-d}u(t)\\ \\ =F(q^{-1})w(t+d)+\left[\frac{B(q^{-1})}{A(q^{-1})}-\frac{q^{-d}B(q^{-1})G(q^{-1})}{A{(}^{-1})C(q^{-1})}\right]u(t)+\frac{G(q^{-1})}{C{(}^{-1})}y(t)\\ \\ =F(q^{-1})w(t+d)+\frac{B(q^{-1})C(q^{-1})-q^{-d}B(q^{-1})G(q^{-1})}{A{(}^{-1})C(q^{-1})}u(t)+\frac{G(q^{-1})}{C{(}^{-1})}y(t)\end{array}\text{\hspace{1em}(9)}$$
由$eq(5)$可得
$$\begin{array}{c}C(q^{-1})=A(q^{-1})F(q^{-1})+q^{-d}G(q^{-1})\\ \\ \Rightarrow B(q^{-1})C(q^{-1})=A(q^{-1})B(q^{-1})F(q^{-1})+q^{-d}B(q^{-1})G(q^{-1})\end{array}$$
代入$eq(9)$可得
$$\begin{array}{c}y(t+d)=F(q^{-1})w(t+d)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)+\frac{G(q^{-1})}{C{(}^{-1})}y(t)\end{array}\text{\hspace{1em}(10)}$$
$eq(10)$也称为控制对象的预测模型。

代入$eq(3)$可得，方差为
$$\begin{array}{c}\begin{array}{rl}J& =E[y^2(t+d)]=\\ & =E\left\{{\left[F(q^{-1})w(t+d)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)+\frac{G(q^{-1})}{C(q^{-1})}y(t)\right]}^2\right\}=\\ & =E\left\{[F(q^{-1})w(t+d){]}^2\right\}+E\left\{{\left[\frac{G(q^{-1})}{C(q^{-1})}y(t)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)\right]}^2\right\}+\\ & +2E\left\{F(q^{-1})w(t+d)\left[\frac{G(q^{-1})}{C(q^{-1})}y(t)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)\right]\right\}\end{array}\end{array}\text{\hspace{1em}(11)}$$
由于$w(t+d)$与$y(t),u(t)$互不相关，所以$eq(11)$最后一项为0，即
$$\begin{array}{c}\begin{array}{r}E\left\{F(q^{-1})w(t+d)\left[\frac{G(q^{-1})}{C(q^{-1})}y(t)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)\right]\right\}=0\end{array}\end{array}$$
此时，
$$\begin{array}{c}\begin{array}{r}E[y^2(t+d)]=E\left\{[F(q^{-1})w(t+d){]}^2\right\}+E\left\{{\left[\frac{G(q^{-1})}{C(q^{-1})}y(t)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)\right]}^2\right\}\end{array}\end{array}$$
又$E\left\{{\left[\frac{G(q^{-1})}{C(q^{-1})}y(t)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)\right]}^2\right\}\ge 0$，所以
$$\begin{array}{c}\begin{array}{r}E[y^2(t+d)]\ge E\left\{[F(q^{-1})w(t+d){]}^2\right\}\end{array}\end{array}\text{\hspace{1em}(12)}$$
如果在任意时刻，有
$$\begin{array}{c}\begin{array}{r}\frac{G(q^{-1})}{C(q^{-1})}y(t)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)=0\end{array}\end{array}$$
那么
$$\begin{array}{c}\begin{array}{r}E[y^2(t+d)]=E\left\{[F(q^{-1})w(t+d){]}^2\right\}=J_{min}\end{array}\end{array}\text{\hspace{1em}(13)}$$
可得，最优控制为
$$\begin{array}{c}\begin{array}{r}u(t)=-\frac{G(q^{-1})}{B(q^{-1})F(q^{-1})}y(t)\end{array}\end{array}\text{\hspace{1em}(14)}$$
对应的最小方差为
$$\begin{array}{c}\begin{array}{r}{J}_{\min }=E\left\{(F(q^{-1})w(t+d){)}^2\right\}=(1+f_1^2+\cdots +f_{d-1}^2){\sigma }^2\end{array}\end{array}\text{\hspace{1em}(15)}$$

对应的自校正框图为

3.2 最优预测估值

最优预测估值$\hat{y}(t+d\mid t)$为$y(t),y(t-1),\cdots ,u(t),u(t-1),\cdots$的线性组合，其估计误差为
$$\begin{array}{c}\begin{array}{r}\tilde{y}(t+d\mid t)=y(t+d)-\hat{y}(t+d\mid t)\end{array}\end{array}\text{\hspace{1em}(16)}$$
要求估计误差的方差最小，即
$$\begin{array}{c}{J}^{\prime }=E[{\tilde{y}}^2(t+d\mid t)]=E\left\{[y(t+d)-\hat{y}(t+d\mid t){]}^2\right\}=J_{\min }^{\prime }\end{array}\text{\hspace{1em}(17)}$$
将$eq(4)$代入可得
$$\begin{array}{rl}{J}^{\prime }& =E\left\{\right[\frac{G(q^{-1})}{C(q^{-1})}y(t)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)+F(q^{-1})w(t+d)-\hat{y}(t+d\mid t){]}^2\}=\\ & E\left\{\right[\frac{G(q^{-1})}{C(q^{-1})}y(t)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)-\hat{y}(t+d\mid t){]}^2+\\ & E[F(q^{-1})w(t+d){]}^2\}\end{array}$$
其中，第一项为非负项，如果有
$$\begin{array}{c}\begin{array}{r}\frac{G(q^{-1})}{C(q^{-1})}y(t)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)-\hat{y}(t+d\mid t)=0\end{array}\end{array}$$
则，
$$\begin{array}{rl}{J}^{\prime }& =E\{[F(q^{-1})w(t+d){]}^2\}=J_{min}^{\prime }\end{array}$$
可得，最优预测估值为
$$\begin{array}{c}\begin{array}{r}\hat{y}(t+d\mid t)=\frac{G(q^{-1})}{C(q^{-1})}y(t)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)=0\end{array}\end{array}$$
输出误差为
$$\begin{array}{c}\begin{array}{rl}\tilde{y}(t+d\mid t)& =y(t+d)-\hat{y}(t+d\mid t)=y(t+d)\\ & =F(q^{-1})w(t+d)\end{array}\end{array}$$

根据前面最优控制推导，由$eq(14)$可得，在最优控制$u(t)$下，最优预测估值$\hat{y}(t+d\mid t)=0$。

4、算法步骤

4.1 显式最小方差调节器

（1）采样读取新的观测值$y(t)$和$u(t)$；

（2）按照递推最小二乘法估计控制对象的参数$a_1,a_2,\cdots ,a_{n_a},b_0,b_1,\cdots ,b_{n_b},c_1,c_2,\cdots ,c_{n_c}$。这里假定$c_0=1$。由最小二乘递推公式得：
$$\begin{gathered} {\theta }^{\mathrm{T}}=\left[\begin{array}{cccccccccccc}{a}_1& a_2& \cdots & a_{n_a}& b_0& b_1& \cdots & b_{n_b}& c_1& c_2& \cdots & c_{n_c}\end{array}\right] \\ {\mathbf{\phi }}^{\mathrm{T}}(t-1)=\left[\begin{array}{ccccccccc}-y(t-1)& \cdots & -y(t-n_a)& u(t-d)& \cdots & u(t-n_b)& \hat{w}(t-1)& \cdots & \hat{w}(t-n_c)\end{array}\right] \end{gathered}$$
$$\begin{array}{c}\begin{array}{rl}& \hat{\epsilon }(t)=y(t)-{\mathbf{\phi }}^{\mathrm{T}}(t-1)\hat{\mathbf{\theta }}(t-1)\\ & \hat{\mathbf{\theta }}(t)=\hat{\mathbf{\theta }}(t-1)+\mathbf{K}(t)[y(t)-{\mathbf{\phi }}^{\mathrm{T}}(t-1)\hat{\mathbf{\theta }}(t-1)]\\ & K(t)=P(t-1)\mathbf{\phi }(t-1)[1+{\mathbf{\phi }}^{\mathrm{T}}(t-1)\mathbf{P}(t-1)\mathbf{\phi }(t-1){]}^{-1}\\ & P(t)=P(t-1)-P(t-1)\mathbf{\phi }(t-1){\left[1+{\mathbf{\phi }}^{\mathrm{T}}(t-1)\mathbf{P}(t-1)\mathbf{\phi }(t-1)\right]}^{-1}{\mathbf{\phi }}^{\mathrm{T}}(t-1)\mathbf{P}(t-1)\end{array}\end{array}$$
（3）按$eq(5)$计算$F(q^{-1})$和$G(q^{-1})$的系数；

（4）按$eq(14)$计算最优控制$u(t)$；

（5）采样时间由$t$变为$t+1$，重复上述步骤。

显式自校正调节器直接对控制对象的参数进行辨识，其计算比较复杂；其次，由于自校正调节器是闭环控制，还需要检验反馈通道的阶数是否符合闭环辨识的要求。

4.2 隐式最小方差调节器

隐式最小方差调节器不对控制对象的参数进行辨识，而是直接辨识调节器的参数。

根据等式
$$\begin{array}{c}\begin{array}{r}\frac{G(q^{-1})}{C(q^{-1})}y(t)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)=0\end{array}\end{array}\text{\hspace{1em}(18)}$$
设
$$\begin{array}{c}\begin{array}{rl}& \frac{G(q^{-1})}{C(q^{-1})}={\alpha }_0+{\alpha }_1{q}^{-1}+\cdots +{\alpha }_{n_a-1}{q}^{-n_a+1}=A_{\alpha }(q^{-1})\\ & \frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}={\beta }_0+{\beta }_1{q}^{-1}+\cdots +{\beta }_{n_b+d-1}{q}^{-n_b-d+1}=B_{\beta }(q^{-1})\end{array}\end{array}\text{\hspace{1em}(19)}$$
代入$eq(18)$，可得
$$\begin{array}{c}\begin{array}{r}{\alpha }_{0}y(t)+{\alpha }_{1}y(t-1)+\cdots +{\alpha }_{n_a-1}y(t-n_a+1)+\\ {\beta }_{0}u(t)+{\beta }_{1}u(t-1)+\cdots +{\beta }_{n_b+d-1}u(t-n_b-d+1)=0\end{array}\end{array}\text{\hspace{1em}(20)}$$
其中，待估参数为${\alpha }_0,{\alpha }_1,\cdots ,{\alpha }_{n_a-1},{\beta }_0,{\beta }_1,\cdots ,{\beta }_{n_b+d-1}$。按照$eq(20)$不能进行辨识，利用$eq(10)$进行参数辨识。将$eq(19)$代入$eq(10)$得，
$$\begin{array}{c}\begin{array}{r}y(t+d)=F(q^{-1})w(t+d)+\frac{B(q^{-1})F(q^{-1})}{C(q^{-1})}u(t)+\frac{G(q^{-1})}{C{(}^{-1})}y(t)\\ ={\alpha }_{0}y(t)+{\alpha }_{1}y(t-1)+\cdots +{\alpha }_{n_a-1}y(t-n_a+1)+\\ {\beta }_{0}u(t)+{\beta }_{1}u(t-1)+\cdots +{\beta }_{n_b+d-1}u(t-n_b-d+1)+e(t+d)\end{array}\end{array}\text{\hspace{1em}(21)}$$
其中，$e(t+d)=F(q^{-1})w(t+d)$。为了使调节器参数便于辨识，应尽量减少被估参数数目，可凭经验给出${\hat{\beta }}_0$，则
$$\begin{gathered} {\theta }^{\mathrm{T}}=\left[\begin{array}{cccccccc}{\alpha }_0& {\alpha }_1& \cdots & {\alpha }_{n_a-1}& {\beta }_1& {\beta }_2& \cdots & {\beta }_{n_b+d-1}\end{array}\right] \\ {\mathbf{\phi }}^{\mathrm{T}}(t)=\left[\begin{array}{ccccccc}y(t)& y(t-1)& \cdots & y(t-n_a+1)& u(t-1)& \cdots & u(t-n_b-d+1)\end{array}\right] \end{gathered}$$
$eq(21)$可写成
$$\begin{array}{c}\begin{array}{r}y(t+d)-{\hat{\beta }}_{0}u(t)={\mathbf{\phi }}^{⊺}(t)\hat{\mathbf{\theta }}(t)+e(t+d)\end{array}\end{array}\text{\hspace{1em}(22)}$$
递推最小二乘参数估计公式如下：
$$\begin{array}{c}\begin{array}{rl}& \hat{\theta }(t+1)=\hat{\theta }(t)+\mathbf{K}(t)[y(t+d)-{\hat{\beta }}_{0}u(t)-{\phi }^{\mathrm{T}}(t)\hat{\theta }(t)]\\ & \mathbf{K}(t)=\mathbf{P}(t-1)\phi (t-1)[1+{\phi }^{\mathrm{T}}(t-1)\mathbf{P}(t-1)\phi (t-1){]}^{-1}\\ & \mathbf{P}(t)=\mathbf{P}(t-1)-\mathbf{P}(t-1)\phi (t-1)[1+{\phi }^{\mathrm{T}}(t-1)\mathbf{P}(t-1)\phi (t-1){]}^{-1}{\phi }^{\mathrm{T}}(t-1)\mathbf{P}(t-1)\end{array}\end{array}\text{\hspace{1em}(23)}$$
最优控制为
$$\begin{array}{c}\begin{array}{r}u(t)=-\frac{1}{{\hat{\beta }}_0}[{\hat{\alpha }}_{0}y(t)+{\hat{\alpha }}_{1}y(t-1)+\cdots +{\hat{\alpha }}_{n_a-1}y(t-n_a+1)\\ +{\hat{\beta }}_{1}u(t-1)+\cdots +{\hat{\beta }}_{n_b+d-1}u(t-n_b-d+1)]\end{array}\end{array}\text{\hspace{1em}(24)}$$
隐式最小方差自校正调节器计算步骤如下：
（1）采样读取新的观测值$y(t+d)$；

（2）按$eq(23)$估计$\hat{\theta }(t)$；

（3）按$eq(24)$计算最佳控制$\hat{u}(t)$；

（4）采样时间$t+d$变为$t+d+1$，重复上述步骤。

5、稳定性分析

前面在进行最优控制推到时，需要假设$B(q^{-1}),C(q^{-1})$的所有零点都位于$q^{-1}$复平面单位圆外。这一节将对其进行解释：

由$eq(14)$可得最优控制，将其代入被控对象的微分方程$eq(1)$，可得
$$\begin{array}{c}\begin{array}{r}A(q^{-1})y(t)=q^{-d}B(q^{-1})[-\frac{G(q^{-1})}{B(q^{-1})F(q^{-1})}y(t)]+C(q^{-1})w(t)\\ B(q^{-1})\left[A(q^{-1})F(q^{-1})+q^{-d}G(q^{-1})\right]y(t)=B(q^{-1})F(q^{-1})C(q^{-1})w(t)\end{array}\end{array}\text{\hspace{1em}(25)}$$
由$eq(5)$得
$$\begin{array}{c}\begin{array}{r}C(q^{-1})=A(q^{-1})F(q^{-1})+q^{-d}G(q^{-1})\end{array}\end{array}$$
将其代入$eq(25)$可得，
$$\begin{array}{c}\begin{array}{r}B(q^{-1})C(q^{-1})y(t)=B(q^{-1})F(q^{-1})C(q^{-1})w(t)\end{array}\end{array}$$
可得闭环系统的特征方程为
$$\begin{array}{c}\begin{array}{r}B(q^{-1})C(q^{-1})=0\end{array}\end{array}$$
为了使闭环系统稳定，$B(q^{-1})$和$C(q^{-1})$的零点必须在单位圆外。如果$B(q^{-1})$的零点在单位圆外，则由
$$\begin{array}{c}\begin{array}{r}A(q^{-1})y(t+d)=B(q^{-1})u(t)+C(q^{-1})w(t+d)\end{array}\end{array}$$
可得，在$y(t+d)$有界时，能保证序列$u(t)$有界。如果$B(q^{-1})$的零点不在单位圆外，则$y(t+d)$有界时，不能保证序列$u(t)$有界。因此最小方差自校正调节器只适用于$B(q^{-1})$的零点在单位圆外的系统。

6、例题

（1）当$d=1$时，由$eq(5)$可得，
$$\begin{array}{c}\frac{C(q^{-1})}{A(q^{-1})}=F(q^{-1})+\frac{q^{-d}G(q^{-1})}{A(q^{-1})}\end{array}$$
求得
$$\begin{array}{c}F(q^{-1})=1\\ G(q^{-1})=3.2+0.2q^{-1}\end{array}$$
代入$eq(14)$可得，
$$\begin{array}{c}\begin{array}{rl}u(t)& =-\frac{G(q^{-1})}{B(q^{-1})F(q^{-1})}y(t)\\ & =-\frac{3.2+0.2q^{-1}}{1+0.5q^{-1}}y(t)\\ \Rightarrow u(t)& =-0.5u(t-1)-3.2y(t)-0.2y(t-1)\end{array}\end{array}$$
输出误差为
$$\tilde{y}(t)=F(q^{-1})w(t)=w(t)$$
输出误差方差为
$$E({\tilde{y}}^2(t))={\sigma }^2$$
（2）当$d=2$时，表明系统延迟增大，性能变差。由$eq(5)$可得，
$$\begin{array}{c}\frac{C(q^{-1})}{A(q^{-1})}=F(q^{-1})+\frac{q^{-d}G(q^{-1})}{A(q^{-1})}\end{array}$$
求得
$$\begin{array}{c}F(q^{-1})=1+3.2q^{-1}\\ G(q^{-1})=5.64-2.24q^{-1}\end{array}$$
代入$eq(14)$可得，
$$\begin{array}{c}\begin{array}{rl}u(t)& =-\frac{G(q^{-1})}{B(q^{-1})F(q^{-1})}y(t)\\ & =-\frac{5.64-2.24q^{-1}}{(1+0.5q^{-1})(1+3.2q^{-1})}y(t)\\ & =-\frac{5.64-2.24q^{-1}}{1+3.7q^{-1}+1.6q^{-2}}y(t)\\ \Rightarrow u(t)& =-3.7u(t-1)-1.6u(t-2)-5.64y(t)+2.24y(t-1)\end{array}\end{array}$$
输出误差为
$$\tilde{y}(t)=F(q^{-1})w(t)=w(t)+3.2w(t-1)$$
输出误差方差为
$$E({\tilde{y}}^2(t))=(1+3.2^2){\sigma }^2=11.24{\sigma }^2$$

对应的matlab代码如下：

% 隐式最小方差调节器% y_hat : 观测数据
% na : A(q^{-1}) = 1 + a_1*q^{-1} + a_2*q^{-2} + ... + a_{na}*q_{-na}
% nb : B(q^{-1}) = b_0 + b_1*q^{-1} + ... + b_{nb}q_{-nb}
% nc : C(q^{-1}) = c_0 + c_1*q^{-1} + ... + c_{nc}q_{-nc}
% d : 延时步数clc
clearA = [1, -1.7, -0.7];
B = [1, 0.5];
C = [1, 1.5, 0.9];
d = 2;na = length(A) - 1;
nb = length(B) - 1;
nc = length(C) - 1;% 初始化
beta0 = 1; % 经验值
theta = zeros(na + nb + d, 1); % 参数向量
P = eye(na + nb + d); % 协方差矩阵
N = 1000;% 产生初始观测
y_hat = zeros(na + nb + d, 1);
u = zeros(na + nb + d, 1);for t=1:N% 构造数据矩阵phi = [y_hat(end-na+1:end); u(end-nb-d+1:end)];% 计算PP = P - P*phi* inv((1 + phi' * P * phi)) * phi' * P;% 计算KK = P * phi * inv((1 + phi' * P * phi));% 计算thetatheta = theta + K * (y_hat(end) - beta0 * u(end) - phi' * theta);% 计算uu_hat = -1/beta0 * phi' * theta;%u_hat = 0;u = [u; u_hat] ;y = phi' * theta + beta0 * u(end) + 0.1*sin(0.01*t);%+ (2*rand - 1); % 根据u产生观测数据y_hat = [y_hat; y]; % 将新观测到的数据放入数组
end% 绘图
plot(y_hat)
hold on
plot(u)
legend("y_{hat}","u")disp("y_hat 方差为：")
disp(var(y_hat))

在计算$u$时，当输入最优控制时，结果如下：

方差为： 0.0044

当注释掉最优控制，令其为0时，

%u_hat = -1/beta0 * phi' * theta;
u_hat = 0;

方差为：37.5501

可以看到在没有控制时，输出的方差很大。

7、总结

从以上推导可以看出最小方差调节器在得到最优控制$u(t)$时，对应的最优预测$\tilde{y}(t+d\mid t)$为0。在运用该算法时，使用隐式最小方差调节器实现起来较为容易。

下一篇将介绍最小方差自校正控制器，用于跟踪参考信号。