线性代数-微分方程和矩阵指数-23

Linear Algebra-微分方程和矩阵指数-23

一、知识概要

今天我们来学习线性代数在微分方程中的应用。本节会介绍一阶线性常微分方程的矩阵解法——也就是如何把微分方程用矩阵抽象出来，通过”解耦”思想计算出对应系数，最终得到方程的解。

这里会牵涉到计算 $e^{At}$ 的问题（其中 $A$ 是矩阵），所以我们也会讨论当幂指数是矩阵时，这个算式该如何计算。最后我们还会扩展介绍高阶微分方程的降阶求解方法，把整个知识链串起来。

二、解微分方程

其实解决微分方程问题，重点在于掌握整个解题流程。我们通过一道具体的例题来一步步讲解。

【例 1】：

$$\frac{d u_{1}}{dt}=-u_{1}+2 u_{2};\frac{d u_{2}}{dt}=u_{1}-2 u_{2};u(0)=\begin{bmatrix}1 \\ 0\end{bmatrix}$$

求 $u(t)$。

在开始解题之前，我们首先要搞清楚这个微分方程的实际意义。我们看到初值 $u(0)=\begin{bmatrix}1 \\ 0\end{bmatrix}$，这说明在初始 0 时刻，$u_{1}=1$，所有的值都在$u_{1}$中；而此时$u_{2}=0$。但是随着时间流逝，当 $t$ 增加时，我们可以看到$\frac{d u_{2}}{dt}>0$，说明$u_{2}$的导数大于 0，因此$u_{2}$会慢慢增加，而$u_{1}$会慢慢减少（你可以理解为”值”从$u_{1}$流向了$u_{2}$）。最终会达到一个稳定状态，具体结果需要我们计算得到。

现在我们列出方程$\frac{d u}{d t}=A u$，其中系数矩阵 $A$ 综合了$\frac{d u_{1}}{dt}$与$\frac{d u_{2}}{d t}$的系数，可以写做$\begin{bmatrix}-1 & 2 \\ 1 & -2\end{bmatrix}$。我们先直接给出通解形式：

$$\sum_{i=1}^{n} C_{i} e^{\lambda_{i} t} x_{i}$$

这个通解形式是怎么来的，我们这里不深入研究。你可以简单理解为：通解就是几个纯指数解的线性组合，随便挑一个代入原方程验证一下就能成立。

接下来我们解 $A$ 矩阵的特征值与特征向量。不难得到这个矩阵有两个特征值：$\lambda_{1}=0$，$\lambda_{2}=-3$。对应的特征向量为：$x_{1}=\begin{bmatrix}2 \\ 1\end{bmatrix}$，$x_{2}=\begin{bmatrix}1 \\ -1\end{bmatrix}$。把这些结果代入通解，得到形式如下：

$$u(t)=C_{1} e^{0 t}\begin{bmatrix} 2 \\ 1 \end{bmatrix}+C_{2} e^{-3 t}\begin{bmatrix} 1 \\ -1 \end{bmatrix}$$

再代入初值$u(0)=\begin{bmatrix}1 \\ 0\end{bmatrix}$，就可以确定$C_{1}$与$C_{2}$，最终解为：

$$u(t)=\frac{1}{3}\begin{bmatrix} 2 \\ 1 \end{bmatrix}+\frac{1}{3} e^{-3 t}\begin{bmatrix} 1 \\ -1 \end{bmatrix}$$

现在我们来分析一下这个通解。你会发现，随着时间 $t$ 的增加，后一项$\frac{1}{3} e^{-3 t}\begin{bmatrix}1 \\ -1\end{bmatrix}$逐渐衰减，最后趋于 0；而前一项$\frac{1}{3}\begin{bmatrix}2 \\ 1\end{bmatrix}$不随时间改变。这正好符合我们一开始分析微分方程意义的时候对 $u$ 走势的预判。

通过这道题，我们可以总结出解决微分方程过程中观察到的几个重要特点：

(1) 特征值是负数时，$u(t)$趋于 0

这个特点很直观，但要注意一种特殊情况——当特征值为复数时，比如写成 $a + bi$，我们该怎么判断$u(t)$的趋势呢？答案是：只有实数部分决定$u(t)$的趋势。我们可以画一个复数坐标系来说明：

不难看出，投影到实数轴，只有实数部分$a$决定正负性，而虚部$b$只是在另一条轴上指明方向，所以不影响我们对趋势的判断。

(2) 稳态存在的条件

当稳态存在时（比如我们【例 1】中最后 $t$ 趋于无穷时，$u$ 趋于一个确定值），必然有一个特征值等于 0，其余特征值的实部全部小于 0。

(3) 解的收敛性判断

如果有任何特征值的实数部分大于 0，则解无法收敛。

三、解耦与 $e^{A t}$

3.1 解耦

回到【例 1】，我们有方程 $\frac{d u}{d t}=A u$。在矩阵 $A$ 中，$u_{1}$ 和 $u_{2}$ 是耦合在一起的——也就是说，每个导数都同时包含两个未知函数，它们互相影响。我们处理的方法，就是通过计算特征值与特征向量把 $A$ 解耦。我们设 $S$ 是特征向量矩阵。

令 $u = SV$，那么：

$$\frac{d u}{d t}=S \frac{d v}{d t}=A u=A S V$$

两边提取出 $S$ 后，我们再同时左乘 $S^{-1}$，就得到：

$$\frac{dv}{dt}=S^{-1} ASV=\Lambda V$$

这样就得到了关于 $V$ 的对角化方程组。新方程组不存在耦合，$\frac{d v_{1}}{dt}=\lambda_{1} v_{1}$，$\frac{d v_{2}}{dt}=\lambda_{2} v_{2}$，以此类推。最终可得到：

$$v(t)=e^{\Lambda t} v(0)$$

同理，我们也有：$u(t) = e^{At}u(0) = S e^{\Lambda t} S^{-1}u(0)$。这里就牵扯到了一个新问题：$e^{At}$ 和 $e^{\Lambda t}$ 到底是什么？表面来看，$e^{At}$ 就是 $u(t)$ 的解，那么为什么 $e^{At}$ 等于 $S e^{\Lambda t} S^{-1}$ 呢？它们到底表示什么含义呢？接下来我们就来详细介绍这些式子。

3.2 $e^{At}$ 与 $S(e^{\Lambda t}) S^{-1}$

我们首先回忆一下熟悉的幂级数公式：

$$e^{x}=1+x+\frac{x^{2}}{2}+\frac{x^{3}}{6}+\cdots \cdots+\frac{x^{n}}{n!}+\cdots \cdots$$

这个级数对全体实数都收敛，所以我们不用考虑特征值带来的收敛域问题。

现在我们把它扩展到矩阵的计算中，很自然地，用单位矩阵 $I$ 代替 1，用矩阵 $A$ 代替 $x$，就得到了矩阵指数的定义：

$$e^{A t}=I+A t+\frac{(A t)^{2}}{2}+\frac{(A t)^{3}}{6}+\cdots \cdots+\frac{(A t)^{n}}{n!}+\cdots \cdots$$

接下来我们利用对角化形式化简 $A$，我们知道对于可对角化矩阵 $A = S\Lambda S^{-1}$，代入上式每一项：

$$e^{A t}=S S^{-1}+S \Lambda t S^{-1}+\frac{S(\Lambda t)^{2} S^{-1}}{2}+\frac{S(\Lambda t)^{3} S^{-1}}{6}+\cdots \cdots$$

把 $S$ 和 $S^{-1}$ 从每一项中提取出来，得到：

$$e^{A t}=S\left(I+\Lambda t+\frac{(\Lambda t)^{2}}{2}+\frac{(\Lambda t)^{3}}{6}+\cdots \cdots\right) S^{-1}$$

再根据幂级数公式，括号里的部分正好就是 $e^{\Lambda t}$，因此我们得到了重要结论：

$$e^{A t}=S(e^{\Lambda t}) S^{-1}$$

注意：这步的化简是有条件的，即 $A$ 必须可以对角化——也就是说，$A$ 要有 $n$ 个独立的特征向量，此时 $S^{-1}$ 才存在。

3.3 矩阵指数

上面的等式将对角矩阵与一般矩阵联系了起来，那么其中的$e^{\Lambda t}$ 具体是什么样子呢？

我们知道对角矩阵 $\Lambda$ 具有如下形式：

$$\Lambda=\begin{bmatrix}\lambda_{1} & 0 & \cdots & 0 \\ & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & \lambda_{n}\end{bmatrix}$$

那么 $e^{\Lambda t}$ 就可以简单表示为：

$$\begin{bmatrix} e^{\lambda_{1} t} & 0 & 0 & 0 \\ 0 & e^{\lambda_{2} t} & \cdots & 0 \\ \vdots & 0 & \ddots & \vdots \\ 0 & 0 & 0 & e^{\lambda_{n} t} \end{bmatrix}$$

也就是说，矩阵指数 $e^{\Lambda t}$ 仍然是对角矩阵，其主对角线上的元素都是 $e^{\lambda_{i} t}$，其余位置都为 0。

同样，这里判断是否收敛的方法和微分方程中的判别思路差不多，即通过比较 $\lambda$ 实部的符号来判断，这个我们前面已经讲过了。

四、二阶微分方程的解

讲完了一阶微分方程，我们再来看看二阶微分方程该怎么处理。实际上，我们可以通过降阶处理，把它转化为我们已经会解的一阶形式。具体步骤如下：

对于二阶微分方程：

$$y^{\prime \prime}+b y'+k y=0$$

我们设$u=\begin{bmatrix}y' \\ y\end{bmatrix}$，就可以利用 $u$ 将上面的方程化简为：

$$u'=\begin{bmatrix} y'' \\ y' \end{bmatrix}=\begin{bmatrix} -b & -k \\ 1 & 0 \end{bmatrix}\begin{bmatrix} y' \\ y \end{bmatrix}$$

这样我们就把二阶微分方程化简为了一阶微分方程乘以一个矩阵，问题就转化为了我们刚才学会的解法。

同样的思路可以推广到更高阶，如果是求解一个五阶微分方程，我们只需要像上面那样化简，得到一个 5×5 的矩阵，形式如下：

$$\begin{bmatrix}a & b & c & d & e \\ 1 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0\end{bmatrix}$$

其中 $a$、$b$、$c$、$d$、$e$ 都是原方程中的系数，次对角线上的元素都是 1，其余位置为 0。这样的矩阵就能把五阶微分方程转化为一阶向量方程，接下来只要使用一阶微分方程的方法正常求解就可以了。

这个降阶思路非常巧妙，也非常通用，不管多少阶都可以用这种方法转化。

五、学习感悟

本节内容较多，主要目的是在实际情况下使用矩阵对角化、特征值等方法求解微分方程，给出了一种使用矩阵求解微分方程的通用规律，即高阶降阶，一阶用特征值和特征向量将原系数矩阵 A 解耦，最后得到结果。并介绍了在我们解耦 A 时使用矩阵对角化将其与特征向量联系起来运算的方法。另外介绍了判断收敛性的方法，即看特征值实部绝对值与 1 的大小关系。这些内容都是特征值与特征向量的实际应用，较为重要。