线性代数与编程：特征值、对角化、Jordan 标准型

稳定性

假设有一个「黑盒」，如果值 $u$ 进入这个箱子，就会变成值 $\xi$ 。如果这个「黑盒」每时每刻都在输入 $u$ 和输出 $\xi$ 。但需要强调时间变化时，我们就会写成 $u(t)$ 和 $\xi(t)$ 。但在这里，「黑箱」输出的 $\xi(t)$ 不仅仅与当前输入 $u(t)$ 有关，还与过去输入的 $u$ 有关系。

控制系统模型
- 设 $u$ 为踩油门的力度， $\xi$ 为汽车的速度
- 放开油门的状态， $u = 0$ ，但是这一个瞬间，汽车的速度 $\xi$ 并不是 0，而是在慢慢减少
信号传输模型
- 设 $u$ 为无线通信中发出的信号， $\xi$ 为接收到的信号
- 理想状态下， $\xi(t) = u(t)$ ，但在实际传输过程中会发生衰减、延迟、扭曲甚至反射，最终都会影响信号

我们在这里感兴趣的是模型会不会失控，换句话说，无论从什么状态出发， $\xi(t)$ 都不会超出有限的范围（不会失控）；还是某些情况下， $|\xi(t)|$ 会趋近于无穷大（失控）。

对角矩阵

考虑下面这个模型：

$x(t) = \left( \begin{matrix} 5 & 0 & 0 \\ 0 & -3 & 0 \\ 0 & 0 & 0.8 \end{matrix} \right) x(t-1)$

一计算就变成：

$\left( \begin{matrix} x_1(t) \\ x_2(t) \\ x_3(t) \end{matrix} \right) = \left( \begin{matrix} 5 x_1(t-1) \\ -3 x_2(t-1)\\ 0.8 x_3(t-1) \end{matrix} \right)$

解出每个式子的答案（递推连乘），就是

$\begin{align*} x_1(t) &= 5^{t}x_1(0)\\ x_2(t) &= (-3)^{t}x_2(0)\\ x_3(t) &= 0.8^{t}x_3(0) \end{align*}$

还可以写成

$x(t) = \left( \begin{matrix} 5^t & 0 & 0 \\ 0 & (-3)^t & 0 \\ 0 & 0 & 0.8^t \end{matrix} \right) x(0) = \left( \begin{matrix} 5 & 0 & 0 \\ 0 & -3 & 0 \\ 0 & 0 & 0.8 \end{matrix} \right)^t x(0)$

只要初始值 $x(0)$ 不满足 $x_1(0)=x_2(0)=0$ ，当 $x\to \infty$ ， $x_1(t)$ 和 $x_2(t)$ 都会越走越远，这个系统是会失控的。

能简单判断是否失控的关键在于系数矩阵是 对角矩阵 。只要有一个绝对值大于 1，系统就会失控。

$x(t)=\left( \begin{matrix} a_1^t & & \\ & \ddots & \\ & & a_n^t \end{matrix} \right) x(0) = \left( \begin{matrix} a_1 & & \\ & \ddots & \\ & & a_n \end{matrix} \right)^t x(0)$

可对角化

如果 $A$ 是对角矩阵，问题就可以迎刃而解。那么，对于一般的 $A$ ，是否可以想办法归结为对角矩阵。大部分情况下，这是可行的。

考虑下面这个系统：

$\left( \begin{matrix} x_1(t) \\ x_2(t) \end{matrix} \right) = \left( \begin{matrix} 5 & 1 \\ 1 & 5 \end{matrix} \right) \left( \begin{matrix} x_1(t-1) \\ x_2(t-1) \end{matrix} \right)$

展开成方程：

$\begin{align*} x_1(t) &= 5x_1(t-1) + x_2(t-1)\\ x_2(t) &= x_1(t-1) + 5x_2(t-1) \end{align*}$

这种形式没有办法处理，那么可以考虑这个方程：

$\begin{align*} y_1(t) &= x_1(t) + x_2(t) \\ y_2(t) &= x_1(t) - x_2(t) \end{align*}$

于是，我们就有一个漂亮的结果：

$\begin{align*} y_1(t) &= x_1(t) + x_2(t) \\ &= (5x_1(t-1) + x_2(t-1)) + (x_1(t-1) + 5x_2(t-1)) \\ &= 6x_1(t - 1) + 6x_2(t - 1) \\ &=6y_1(t - 1)\\ \\ y_2(t) &= x_1(t) - x_2(t)\\ &= (5x_1(t-1) + x_2(t-1))-(x_1(t-1) + 5x_2(t-1))\\ &= 4x_1(t - 1) - 4x_2(t - 1)\\ &= 4y_2(t-1) \end{align*}$

这样就变成一个简单的问题了：

$\begin{align*} y_1(t) &= 6^t y_1(0)\\ y_2(t) &= 4^t y_2(0) \end{align*}$

很容易知道 $y1,y2$ 都会失控。

接下去，只需要把 $y1,y_2$ 还原成 $x_1,x_2$ ，就可以了。只需要：

$\begin{align*} x_1(t) &= \frac{y_1(t) + y_2(t)}{2} \\ x_2(t) &= \frac{y_1(t) - y_2(t)}{2} \end{align*}$

代入上面的式子，

$\begin{align*} x_1(t) &= \frac{6^t y_1(0) + 4^t y_2(0)}{2} \\ &=(\frac{6^t+4^t}{2})x_1(0) + (\frac{6^t-4^t}{2})x_2(0) \\ \\ x_2(t) &= \frac{6^t y_1(0) - 4^t y_2(0)}{2} \\ &= (\frac{6^t - 4^t}{2})x_1(0) + (\frac{6^t+4^t}{2})x_2(0) \end{align*}$

当 $t\to \infty$ 时， $(6^t \pm 4^t)/2 \to \infty$ ，所以， $x1,x2$ 的系统一样会失控。

用矩阵的语言表达

$y(t) = Cx(t), C=\left( \begin{matrix} 1 & 1 \\ 1 & -1 \end{matrix}\right)$

原来的问题 $x(t) = Ax(t - 1)$ 就变成：

$y(t) = \Lambda y(t-1), \Lambda=\left(\begin{matrix} 6 & 0 \\ 0 & 4 \end{matrix}\right)$

这里的 $\Lambda$ 是对角矩阵，所以很容易的得出

$y(t)=\Lambda^t y(0)=\left(\begin{matrix} 6^t & 0 \\ 0 & 4^t \end{matrix}\right)\left(\begin{matrix}y_1(0) \\ y_2(0) \end{matrix}\right)$

接下来，把 $y$ 还原成 $x$ 即可。还原时使用变换：

$\left( \begin{matrix} x_1(t) \\ x_2(t) \end{matrix} \right) =\left( \begin{matrix} \frac{y_1(t) + y_2(t)}{2} \\ \frac{y_1(t) - y_2(t)}{2} \end{matrix} \right)$

本质上其实是，

$x(t)=C^{-1}y(t)=\left(\begin{matrix} 1/2 & 1/2 \\ 1/2 & -1/2\end{matrix}\right)\left(\begin{matrix}y_1(t) \\ y_2(t) \end{matrix}\right)$

使用上面的矩阵还原 $x$ ，可得

$\begin{align*} \left( \begin{matrix} x_1(t) \\ x_2(t) \end{matrix} \right) &=\left(\begin{matrix} 1/2 & 1/2 \\ 1/2 & -1/2\end{matrix}\right)\left(\begin{matrix}y_1(t) \\ y_2(t) \end{matrix}\right) \\ &=\left(\begin{matrix} 1/2 & 1/2 \\ 1/2 & -1/2\end{matrix}\right)\left(\begin{matrix} 6^t & 0 \\ 0 & 4^t \end{matrix}\right)\left(\begin{matrix}y_1(0) \\ y_2(0) \end{matrix}\right)\\ &=\left(\begin{matrix} 6^t/2 & 4^t/2 \\ 6^t/2 & -4^t/2\end{matrix}\right)\left(\begin{matrix}y_1(0) \\ y_2(0) \end{matrix}\right) \end{align*}$

其中，初始值也满足 $y(0)=Cx(0)$ ，于是：

$\begin{align*} \left( \begin{matrix} x_1(t) \\ x_2(t) \end{matrix} \right) &=\left(\begin{matrix} 6^t/2 & 4^t/2 \\ 6^t/2 & -4^t/2\end{matrix}\right)\left(\begin{matrix}1 & 1 \\ 1 & -1 \end{matrix}\right)\left(\begin{matrix}x_1(0) \\ x_2(0) \end{matrix}\right)\\ &=\left(\begin{matrix} \frac{6^t+4^t}{2} & \frac{6^t-4^t}{2} \\ \frac{6^t-4^t}{2} & \frac{6^t + 4^t}{2}\end{matrix}\right)\left(\begin{matrix}x_1(0) \\ x_2(0) \end{matrix}\right) \end{align*}$

综上，用矩阵运算得出了相同的结果。

一般化

总结一下这些流程：

选定某矩阵 $C$ ，用其将 $x(t)$ 变成另一组变量 $y(t) = Cx(t)$ ；
把关于 $x(t)$ 的差分方程改写成关于 $y(t)$ 的方程；
新的方程的系数矩阵如果是对角矩阵，那么就很容易求解；
把解出来的 $y(t)$ 还原为 $x(t)$ 。

整个过程的关键在于，变换后的 $y(t)$ 的方程要具有对角矩阵的形式。

并不是任意的 $C$ 都可以，它一定要在 $x$ 和 $y$ 之间建立一一对应（双射）。

我们考虑的问题是：

$\begin{align*} y(t) &= Cx(t)\\ x(t) &= C^{-1}y(t) \end{align*}$

也就是 $x\to y$ 为顺， $y\to x$ 为逆。现在反过来，改写成：

$\begin{align*} x(t) &= Py(t)\\ y(t) &= P^{-1}x(t) \end{align*}$

这样 $y\to x$ 为顺， $x \to y$ 为逆。

很显然，

$\begin{align*} C &= P^{-1}\\ P & = C^{-1} \end{align*}$

于是有：

$\begin{align*} y(t) &= P^{-1}x(t)=P^{-1}Ax(t-1) \\ &=P^{-1}A(Py(t-1))=(P^{-1}AP)y(t-1) \end{align*}$

从 $y$ 的角度， $x(t)=Ax(t-1)$ 呈现出以下形态：

$\begin{align*} y(t) &= \Lambda y(t-1)\\ \Lambda &= P^{-1}AP \end{align*}$

如果 $\Lambda$ 是对角矩阵，则可以用 $y(t)=\Lambda^t y(0)$ 求出 $y(t)$ 。

接下来，根据 $x(t)=Py(t)$ 和 $y(0)=P^{-1}x(0)$ ，可以得到最终的 $x$ ：

$x(t)=Py(t)=P\Lambda^t y(0)=P\Lambda^t P^{-1}x(0)$

问题的重点就在于 如何选择合适的可逆矩阵 $P$ ，使得 $P^{-1}AP$ 是对角矩阵。

特征值、特征向量

我们的任务就是要找到对角化所需要的 $P$

$P^{-1}AP\equiv \Lambda=diag(\lambda_1,\cdots,\lambda_n)$

两边左乘 $P$ ，就变成 $AP=P\Lambda$ ，也就是：

$A(p_1,\cdots,p_n)=(p_1,\cdots,p_n) \begin{align*} \left( \begin{matrix} \lambda_1 \\ &\ddots \\ & & \lambda_n \end{matrix} \right) \end{align*}$

按照分块矩阵乘法，分别计算左右两边，可以得到：

$(Ap_1,\cdots,Ap_n)=(\lambda_1 p_1,\cdots,\lambda_n p_n)$

数 $\lambda$ 和向量 $p$ 分别称为「特征值」和「特征向量」。

所以，求合适的 $P$ 分为两步：

求 $A$ 的特征值 $\lambda_1,\cdots,\lambda_n$ 及其对应的特征向量 $p_1,\cdots,p_n$
将求得的特征向量排列起来，即 $P=(p_1,\cdots,p_n)$

严格来讲，上面并不严谨，得确认 $P$ 不是奇异，才可以求 $P^{-1}$ 。

MATLAB：

1
2
3

A = [ 5 3 -4;6 8 -8; 6 9 -9]
eig(A)
X \ A*X % 特征值

Mathematica:

Eigenvalues[({
   {5, 3, -4},
   {6, 8, -8},
   {6, 9, -9}
  })]

从几何学角度，特征向量乘上 $A$ ，除了长度会有伸缩变化，方向不发生改变。这里的长度变化率，就是特征值。

所以判断是否是特征向量就可以通过这个特性：

例：设 $A=\left[\begin{matrix} 1 & 6 \\ 5 & 2\end{matrix}\right]$ ， $\textbf{u}=\left[\begin{matrix} 6 \\ -5\end{matrix}\right]$ 和 $\textbf{v}=\left[\begin{matrix} 3 \\ -2\end{matrix}\right]$ ， $\textbf{u}$ 和 $\textbf{v}$ 是否是 $A$ 的特征向量？

解

$\begin{align*} A\textbf{u}&=\left[\begin{matrix} 1 & 6 \\ 5 & 2\end{matrix}\right]\left[\begin{matrix} 6 \\ -5\end{matrix}\right]=\left[\begin{matrix} -24 \\ 20\end{matrix}\right]=-4\left[\begin{matrix} 6 \\ -5\end{matrix}\right]=-4\textbf{u} \\ A\textbf{v}&=\left[\begin{matrix} 1 & 6 \\ 5 & 2\end{matrix}\right]\left[\begin{matrix} 3 \\ -2\end{matrix}\right] =\left[\begin{matrix} -9 \\ 11\end{matrix}\right] \ne\lambda\left[\begin{matrix}3 \\ -2\end{matrix}\right] \end{align*}$

因此， $\textbf{u}$ 是特征值 -4 对应的特征向量，但 $A\textbf{v}$ 不是 $\textbf{v}$ 的倍数，所以 $\textbf{v}$ 不是 $A$ 的特征向量。

特征方程

特征方程是一个数量方程，包含了有关特征值的信息。

假设 $\lambda$ 是特征值，那么我们知道 $Ax=\lambda x$ 。那么可以得到一个方程，即

$(A-\lambda I)x=0$

如果 $A=\left[\begin{matrix} 1 & 6 \\ 5 & 2\end{matrix}\right]$ ，可以得到：

$(A-\lambda I)=\left[\begin{matrix} 1 & 6 \\ 5 & 2\end{matrix}\right] - \left[\begin{matrix} \lambda & 0 \\ 0 & \lambda \end{matrix}\right] = \left[\begin{matrix} 1 - \lambda & 6 \\ 5 & 2 - \lambda \end{matrix}\right]$