LuosBlog
0%

线性状态观测器

发表于 更新于 分类于
本文字数： 5.5k 阅读时长 ≈ 5 分钟
线性状态观测器

前言

系统的可观测性

对于一个开环系统

X˙=AXy=CX\dot{X}=AX\\y=CX

可以得到系统的状态转移方程

X=eAtX(0)X=e^{At}X(0)

根据状态方程的解，可以知道

y=CX=CeAtX(0)y=CX=Ce^{At}X(0)

根据矩阵质数定义的状态转移矩阵有

eAt=I+At+12!A2t2+...e^{At}=I+At+\frac{1}{2!}A^2t^2+...

取前 n 项

eAt=I+At+12!A2t2+...+1(n−1)!An−1tn−1e^{At}=I+At+\frac{1}{2!}A^2t^2+...+\frac{1}{(n-1)!}A^{n-1}t^{n-1}

带入到上面 y 的式子中

y=(CI+CAt+...+1(n−1)!CAn−1tn−1)X(0)=[It...1(n−1)!tn−1][CCA...CAn−1]X(0)y=(CI+CAt+...+\frac{1}{(n-1)!}CA^{n-1}t^{n-1})X(0)\\=\begin{bmatrix}I&t&...&\frac{1}{(n-1)!}t^{n-1}\end{bmatrix}\begin{bmatrix}C\\CA\\...\\CA^{n-1}\end{bmatrix}X(0)

这是一个非齐次线性方程组，要从上面的式子中唯一的解出 n 阶初始状态 x(0) 的充要条件就是下面这个矩阵为满秩。实际的意义就是保证从任意状态来说都能推导出一个唯一的初始状态

定义一个矩阵 O

O=[CCA...CAn−1]O=\begin{bmatrix}C\\CA\\...\\CA^{n-1}\end{bmatrix}

当 Rank(O)=nRank(O)=n 就称该系统是可观测的

连续系统线性状态观测器

对于系统的状态，一般都是可以直接测定的，但是有些情况，系统的状态不能够直接测定，所以要根据系统的输入和输出来估计系统的状态

对于一个系统

X˙=AX+Buy=CX+Du\dot{X}=AX+Bu\\y=CX+Du

设定 X^\hat{X} 为系统状态的估计值， y^\hat{y} 为输出的估计值

X^˙=AX^+Bu+L(y−y^)y^=CX^+Du\dot{\hat{X}}=A\hat{X}+Bu+L(y-\hat{y})\\\hat{y}=C\hat{X}+Du

两式合并得

X^˙=(A−LC)X^+(B−LD)u+Ly\dot{\hat{X}}=(A-LC)\hat{X}+(B-LD)u+Ly

这就是观测器形式。此时将估计值与真实值相减，可以得到

X˙−X^˙=AX+Bu−((A−LC)X^+(B−LD)u+Ly)⇓X˙−X^˙=(A−LC)X−(A−LC)X^=(A−LC)(X−X^)\dot{X}-\dot{\hat{X}}=AX+Bu-((A-LC)\hat{X}+(B-LD)u+Ly)\\\Downarrow\\\dot{X}-\dot{\hat{X}}=(A-LC)X-(A-LC)\hat{X}=(A-LC)(X-\hat{X})

此时，令

ex=X−X^e_x={X}-{\hat{X}}

是估计值与实际值之间的误差。由于我们的目标肯定是使得实际值与估计值是一致的，所以误差值的期望是变小的，也就是 ex→0e_x→0 ，所以根据上式可以得到一个这样的式子

ex˙=(A−LC)ex\dot{e_x}=(A-LC)e_x

所以就需要 A−LCA-LC 的特征值小于 0 ，才能使得该误差值逐渐趋近于 0。所以需要建立一个新的反馈系统使得 ex→0e_x→0 ，所以此时令

L=[l1l2]L=\begin{bmatrix}l_1\\l_2\end{bmatrix}

经过计算之后得到 LL 的值，带入之后选取适当的数值可以使得系统稳定

例子

一个倒立摆模型

Untitled

建立系统状态空间方程，令

X=[x1x2]x1=ϕx2=ϕ˙（不可测）u=1Lx¨blockX=\begin{bmatrix}x_1\\x_2\end{bmatrix}\\x_1=\phi\\x_2=\dot{\phi}（不可测）\\u=\frac{1}{L}\ddot{x}_{block}

则从上式可以得到状态空间方程

X˙=[01gL0]X+[0−1]uy=[10]X+[0]u\dot{X}=\begin{bmatrix}0&1\\\frac{g}{L}&0\end{bmatrix}X+\begin{bmatrix}0\\-1\end{bmatrix}u\\y=\begin{bmatrix}1&0\end{bmatrix}X+[0]u

可以得到

A−LC=[01gL0]−[l1l2][10]=[−l11gL−l20]A-LC=\begin{bmatrix}0&1\\\frac{g}{L}&0\end{bmatrix}-\begin{bmatrix}l_1\\l_2\end{bmatrix}\begin{bmatrix}1&0\end{bmatrix}=\begin{bmatrix}-l_1&1\\\frac{g}{L}-l_2&0\end{bmatrix}

可以假设

g=10 N/m2L=1 mg=10~N/m^2\\L=1~m

带入后得到

A−LC=[−l1110−l20]A-LC=\begin{bmatrix}-l_1&1\\10-l_2&0\end{bmatrix}

求其特征方程为

∣λI−(A−LC)∣=0↓λ2+l1λ+l2−10=0|\lambda I-(A-LC)|=0\\\downarrow\\\lambda^2+l_1\lambda+l_2-10=0

所以，令该特征方程的两个根都为 -1 ，可以得到

l1=2l2=11l_1=2\\l_2=11

所以带入到估计值的方程中可以得到

X^˙=(A−LC)X^+(B−LD)u+Ly⇓X^˙=[−21−10]X^+[0−1]u+[211]y\dot{\hat{X}}=(A-LC)\hat{X}+(B-LD)u+Ly\\\Downarrow\\\dot{\hat{X}}=\begin{bmatrix}-2&1\\-1&0\end{bmatrix}\hat{X}+\begin{bmatrix}0\\-1\end{bmatrix}u+\begin{bmatrix}2\\11\end{bmatrix}y

离散系统线性状态观测器

对于一个离散系统

X(k)=AX(k−1)+Buy(k)=CX(k)+DuX(k)=AX(k-1)+Bu\\y(k)=CX(k)+Du

设定 X^\hat{X} 为系统状态的估计值， y^\hat{y} 为输出的估计值

X^(k)=AX^(k−1)+Bu+L(y(k)−y^(k−1))y^(k)=CX^(k)+Du\hat{X}(k)=A\hat{X}(k-1)+Bu+L\Big(y(k)-\hat{y}(k-1)\Big)\\\hat{y}(k)=C\hat{X}(k)+Du

两式合并得

X^(k)=(A−LC)X^(k−1)+(B−LD)u+Ly(k)\hat{X}(k)=(A-LC)\hat{X}(k-1)+(B-LD)u+Ly(k)

这就是观测器形式。此时将估计值与真实值相减，可以得到

X(k)−X^(k)=AX(k−1)+Bu−((A−LC)X^(k−1)+(B−LD)u+Ly(k))⇓X(k)−X^(k)=(A−LC)(X(k−1)−X^(k−1))X(k)-\hat{X}(k)=AX(k-1)+Bu-((A-LC)\hat{X}(k-1)+(B-LD)u+Ly(k))\\\Downarrow\\X(k)-\hat{X}(k)=(A-LC)(X(k-1)-\hat{X}(k-1))

此时要使估计的误差 e(k)=X(k)−X^(k)e(k)=X(k)-\hat{X}(k) 渐近收敛到 0，矩阵 A−LCA-LC 的特征值必须严格位于复平面单位圆内部

例子

对于上述的倒立摆的例子

X˙=[01gL0]X+[0−1]uy=[10]X+[0]u\dot{X}=\begin{bmatrix}0&1\\\frac{g}{L}&0\end{bmatrix}X+\begin{bmatrix}0\\-1\end{bmatrix}u\\y=\begin{bmatrix}1&0\end{bmatrix}X+[0]u

利用泰勒展开近似计算其离散系统传递函数

X(k)=[1TgLT1]X(k−1)+[0−T]uy(k)=[10]X+[0]uX(k)=\begin{bmatrix}1&T\\\frac{g}{L}T&1\end{bmatrix}X(k-1)+\begin{bmatrix}0\\-T\end{bmatrix}u\\y(k)=\begin{bmatrix}1&0\end{bmatrix}X+\begin{bmatrix}0\end{bmatrix}u

得到

A−LC=[1TgLT1]−[l1l2][10]=[1−l1TgLT−l21]A-LC=\begin{bmatrix}1&T\\\frac{g}{L}T&1\end{bmatrix}-\begin{bmatrix}l_1\\l_2\end{bmatrix}\begin{bmatrix}1&0\end{bmatrix}=\begin{bmatrix}1-l_1&T\\\frac{g}{L}T-l_2&1\end{bmatrix}

可以假设

g=10 N/m2L=1 mg=10~N/m^2\\L=1~m

计算 A−LCA-LC 的特征值，且设定特征值的范围位于单位圆内，可以计算设置

l1=0.3l2=20.01l_1=0.3\\l_2=20.01

所以带入到估计值的方程中可以得到

X^(k)=(A−LC)X^(k−1)+(B−LD)u+Ly(k)⇓X^(k)=[0.70.001−201]X^(k−1)+[0−0.0010−0.001]u+[0.320.01]y(k)\hat{X}(k)=(A-LC)\hat{X}(k-1)+(B-LD)u+Ly(k)\\\Downarrow\\\hat{X}(k)=\begin{bmatrix}0.7&0.001\\-20&1\end{bmatrix}\hat{X}(k-1)+\begin{bmatrix}0&-0.001\\0&-0.001\end{bmatrix}u+\begin{bmatrix}0.3\\20.01\end{bmatrix}y(k)