线性代数学习笔记7-3：解微分方程、矩阵的指数函数（特征值的应用）

之前介绍了求解一阶差分方程，本文介绍求解一阶导数常系数微分方程

常系数微分方程的解是指数形式的$e^{\lambda t}$，基于这个事实，我们得以将问题转为线性代数的问题，求解其指数和系数
另外，在微分方程中将会看到特征值的另一应用：除了可以帮助我们求矩阵的幂，还可以求矩阵的指数

举例：一阶微分方程

对于一阶微分方程 { d u 1 d t = − u 1 + 2 u 2 d u 2 d t = u 1 − 2 u 2 \left\{

\right. {dtdu1​​dtdu2​​​=−u1​+2u2​=u1​−2u2​​，初始条件为$u_1=1,u_2=0$

类比之前的差分方程，微分方程可以表示为一个方程组（未知量组成一个整体的列向量$\mathbf{u}$，系数为$\mathbf{A}$）$$\frac{d\mathbf{u}}{dt}=\mathbf{A}\mathbf{u}$$其中 A = [ − 1 2 1 − 2 ] , u ( 0 ) = [ 1 0 ] \boldsymbol{A}=\left[

\right], \quad \mathbf{u}(0)=\left[\right] A=[−11​2−2​],u(0)=[10​]
我们的最终目标是找到$u(t)$，或者说追踪$u$随时间的变化及最终的值

[结论] 方程$\frac{d\mathbf{u}}{dt}=\mathbf{A}\mathbf{u}$的通解形式为$$\mathbf{u}(\mathrm{t})=c_1{e}^{{\lambda }_{1}t}{\mathbf{x}}_1+c_2{e}^{{\lambda }_{2}t}{\mathbf{x}}_2，其中系数c_i由初值条件给出，{\mathbf{x}}_1和{\mathbf{x}}_2为相应的特征向量$$

• 注意，之前讲的差分方程中，出现纯幂形式$c_1{\lambda }_1^k{\mathbf{x}}_1+c_2{\lambda }_2^k{\mathbf{x}}_2+...$；
• 而微分方程中出现纯指数形式$c_1{e}^{{\lambda }_{1}t}{\mathbf{x}}_1+c_2{e}^{{\lambda }_{2}t}{\mathbf{x}}_2+...$，两者作用类似
  例如将特解$c_1{e}^{{\lambda }_{1}t}{\mathbf{x}}_1$代入方程$\frac{d\mathbf{u}}{dt}=\mathbf{A}\mathbf{u}$中验证：得到$c_1{\lambda }_1{e}^{{\lambda }_{1}t}{\mathbf{x}}_1=\mathbf{A}\mathbf{u}=c_1\mathbf{A}{e}^{{\lambda }_{1}t}{\mathbf{x}}_1$，即${\lambda }_1{\mathbf{x}}_1=\mathbf{A}{\mathbf{x}}_1$，可见，这里纯指数主要是用于保证“微分后形式仍保持不变”

证明：方程的通解$\mathbf{u}(\mathrm{t})=c_1{e}^{{\lambda }_{1}t}{\mathbf{x}}_1+c_2{e}^{{\lambda }_{2}t}{\mathbf{x}}_2$ / 探究为何微分方程中出现指数函数
微分方程$\frac{d\mathbf{u}}{dt}=\mathbf{A}\mathbf{u}$中的$\frac{d\mathbf{u}}{dt}$，就相当于差分方程中前后两项之差${\mathbf{u}}_{k+1}-{\mathbf{u}}_k$，只不过两项的距离无限小；
或者说，从$u(0)$到$u(t)$认为经历了$N\to \infty$个差分方程
那么，将某个时刻的值$u(t){|}_{t=k\Delta }$记为$u_k$，则微分方程可以暂时视为$$\frac{d\mathbf{u}}{dt}=\frac{{\mathbf{u}}_{k+1}-{\mathbf{u}}_k}{\Delta t}=\frac{{\mathbf{u}}_{k+1}-{\mathbf{u}}_k}{t/N}=\mathbf{A}{\mathbf{u}}_k$$
改写得到$${\mathbf{u}}_{k+1}=\left(\frac{t}{N}\mathbf{A}+1\right){\mathbf{u}}_k$$这是上一节所学的内容，可以知道最终的值$\mathbf{u}(\mathrm{t})$，等价于差分方程的最终值${\mathbf{u}}_N$$$\mathbf{u}(t)={\mathbf{u}}_N=c_1{\left(\frac{t}{N}{\lambda }_1+1\right)}^N{\mathbf{x}}_1+c_2{\left(\frac{t}{N}{\lambda }_2+1\right)}^N{\mathbf{x}}_2$$，并且考虑到${\lim }_{k\to \infty }{\left(\frac{1}{k}+1\right)}^k\to e$，从差分变为微分，得到其通解形式：$$\mathbf{u}(\mathrm{t})=c_1{e}^{{\lambda }_{1}t}{\mathbf{x}}_1+c_2{e}^{{\lambda }_{2}t}{\mathbf{x}}_2$$

求解$\frac{d\mathbf{u}}{dt}=\mathbf{A}\mathbf{u}$，其中 A = [ − 1 2 1 − 2 ] , u ( 0 ) = [ 1 0 ] \boldsymbol{A}=\left[

\right], \quad \mathbf{u}(0)=\left[\right] A=[−11​2−2​],u(0)=[10​]：

• 首先求$\mathbf{A}$的特征值： λ 1 = 0 , x 1 = [ 2 1 ] \lambda_{1}=0,\mathbf{x} 1=\left[
  21" role="presentation" style="position: relative;">21$$\begin{array}{l}2\\ 1\end{array}$$
  \right] λ1​=0,x1=[21​]和 λ 1 = − 3 , x 1 = [ 1 − 1 ] \lambda_{1}=-3,\mathbf{x} 1=\left[
  1−1" role="presentation" style="position: relative;">1−1$$\begin{array}{l}1\\ -1\end{array}$$
  \right] λ1​=−3,x1=[1−1​]
• 将初始状态表示为特征向量的线性组合，也就是求解$\mathbf{S}\mathbf{c}=\mathbf{u}(0)$： u ( 0 ) = [ 1 0 ] = c 1 [ 2 1 ] + c 2 [ 1 − 1 ] ,其中c  1 = c 2 = 1 / 3 \mathbf{u}(0)=\left[
  10" role="presentation" style="position: relative;">10$$\begin{array}{l}1\\ 0\end{array}$$
  \right]=c_{1}\left[
  21" role="presentation" style="position: relative;">21$$\begin{array}{l}2\\ 1\end{array}$$
  \right]+c_{2}\left[
  1−1" role="presentation" style="position: relative;">1−1$$\begin{array}{c}1\\ -1\end{array}$$
  \right] \text {,其中c } 1=\mathrm{c} 2=1 / 3 u(0)=[10​]=c1​[21​]+c2​[1−1​],其中c 1=c2=1/3
• 根据上述的通解$\mathbf{u}(\mathrm{t})=c_1{e}^{{\lambda }_{1}t}{\mathbf{x}}_1+c_2{e}^{{\lambda }_{2}t}{\mathbf{x}}_2$，代入系数$c$和特征值$\lambda$、特征向量$\mathbf{x}$解得 u ( t ) = c 1 e λ 1 t x 1 + c 2 e λ 2 t x 2 = c 1 e 0 [ 2 1 ] + c 2 e − 3 t [ 1 − 1 ] = 1 3 [ 2 1 ] + 1 3 e − 3 t [ 1 − 1 ]
  \begin{array}{l} \mathbf{u}(\mathrm{t})=c_{1} e^{\lambda_{1} t} \mathbf{x}_{1}+c_{2} e^{\lambda_{2} t} \mathbf{x}_{2} \\=c_{1} e^{0}\left[\begin{array}{l}2 \\1\end{array}" role="presentation" style="position: relative;">\begin{array}{l} \mathbf{u}(\mathrm{t})=c_{1} e^{\lambda_{1} t} \mathbf{x}_{1}+c_{2} e^{\lambda_{2} t} \mathbf{x}_{2} \\=c_{1} e^{0}\left[\begin{array}{l}2 \\1\end{array}$$\text{begin{array}{l} mathbf{u}(mathrm{t})=c\_1 e^{lambda\_1 t} mathbf{x}\_1+c\_2 e^{lambda\_2 t} mathbf{x}\_2 =c\_1 e^{0}left[begin{array}{l}2 1end{array}}$$
  \right] +c_{2} e^{-3 t}\left[
  1−1" role="presentation" style="position: relative;">1−1$$\begin{array}{c}1\\ -1\end{array}$$
  \right] \\=\frac{1}{3}\left[
  21" role="presentation" style="position: relative;">21$$\begin{array}{l}2\\ 1\end{array}$$
  \right]+\frac{1}{3} e^{-3 t}\left[
  1−1" role="presentation" style="position: relative;">1−1$$\begin{array}{c}1\\ -1\end{array}$$
  \right] \end{array} u(t)=c1​eλ1​tx1​+c2​eλ2​tx2​=c1​e0[21​]+c2​e−3t[1−1​]=31​[21​]+31​e−3t[1−1​]​

分析：随着$t$的增大，答案中第二项会消失，而第一项为稳态（当有0特征值时，就出现稳态SteadyState）’
如果初态 u ( 0 ) = [ 1 0 ] \mathbf{u}(0)=\left[

10" role="presentation" style="position: relative;">10$$\begin{array}{l}1\\ 0\end{array}$$

\right] u(0)=[10​]，那么最终的稳态就是 u ( t ) ∣ t = ∞ = 1 3 [ 2 1 ]

\begin{array}{l} \mathbf{u}(\mathrm{t})|_{t=\infty}=\frac{1}{3}\left[\begin{array}{l}2 \\1\end{array}" role="presentation" style="position: relative;">\begin{array}{l} \mathbf{u}(\mathrm{t})|_{t=\infty}=\frac{1}{3}\left[\begin{array}{l}2 \\1\end{array}$$\text{begin{array}{l} mathbf{u}(mathrm{t})|\_{t=infty}=frac{1}{3}left[begin{array}{l}2 1end{array}}$$

\right]\end{array} u(t)∣t=∞​=31​[21​]​

---

关于“稳态”：

对于通解$\mathbf{u}(\mathrm{t})=c_1{e}^{{\lambda }_{1}t}{\mathbf{x}}_1+c_2{e}^{{\lambda }_{2}t}{\mathbf{x}}_2$
矩阵的特征值给出了$\mathbf{u}(\mathrm{t})$的发展趋势（稳定性）：

• 将$e^{{\lambda }_{1}t}$视为$A{e}^{j\varphi }$的形式，则实部 R e { λ } Re\{\lambda\} Re{λ}决定了稳定性（即决定幅值的增长速度，因为$|e^{a+jb}|=|e^a||e^{jb}|=|e^a|$），虚部 I m { λ } Im\{\lambda\} Im{λ}对应了单位圆上的相位旋转

• R e { λ } > 0 Re\{\lambda\}>0 Re{λ}>0，对应项发散； R e { λ } = 0 Re\{\lambda\}=0 Re{λ}=0，对应项幅值稳定不变； R e { λ } < 0 Re\{\lambda\}<0 Re{λ}<0，对应项消失（$t\to \infty 时\mathbf{u}(\mathrm{t})\to 0$）

对比复习：
①之前的矩阵的幂的情况
如果其所有特征值$|{\lambda }_i|<1$，则$k\to \infty 时{\mathbf{A}}^k\to 0$（因为${\mathbf{\Lambda }}^k\to 0$，故${\mathbf{A}}^k=\mathbf{S}{\mathbf{\Lambda }}^k{\mathbf{S}}^{-1}\to 0$）
②之前的差分方程的稳态情况（主要取决于幅值$|{\lambda }_i|$）
对于实数特征值，特征值$|{\lambda }_i|<1$的项最终会消失，特征值$|{\lambda }_i|=1$的项恒定，特征值$|{\lambda }_i|>1$的项最终不断增长
对于复数特征值，虚部引入了复平面上的“旋转”，故特征值的幅值仍然确定稳态，而相位则对应了每次做矩阵乘法时特征向量的旋转角度

• 综合考虑所有特征值对于解的稳态的影响：
  ①若所有 R e { λ } ≤ 0 Re\{\lambda\}\leq0 Re{λ}≤0（对应项消失/幅值稳定不变），则可以进入稳态
  ②一旦存在 R e { λ } > 0 Re\{\lambda\}>0 Re{λ}>0，则发散

另外，我们比较关注二阶系统的稳定性，一个推论是：

• 假如行列式$det(\mathbf{A})>0$而$trace(\mathbf{A})<0$（这等价于说二阶矩阵有两个实部为负的特征值），则微分方程的解可以进入稳态

扩展：从方程解耦（特征值分解）的角度求解微分方程

下文一切讨论的前提：$\mathbf{A}$的特征向量矩阵$\mathbf{S}$可逆/$\mathbf{A}$有n个无关的特征向量，因为此时才能保证可对角化（从而用于解耦的特征向量数量是足够的）
将矩阵对角化（特征值分解）为$\mathbf{A}=\mathbf{S}\mathbf{\Lambda }{\mathbf{S}}^{-1}$

下面从另一角度考虑一般的一阶微分方程的求解原理，探究为什么一阶微分方程的解是指数函数$e^{\lambda t}$的和的形式（之前已经证明，这里从另一角度理解）；

出发角度是：
将问题转换到另一坐标系（基向量为特征向量），从而解耦方程，得到解，再转换回原坐标系；
具体的解耦方法：
对于未知量相互耦合的方程组，用特征值和特征向量来对角化方程的系数矩阵，可以实现解耦（各未知量没有关系）

对于一阶微分方程 { d u 1 d t = − u 1 + 2 u 2 d u 2 d t = u 1 − 2 u 2 \left\{

du1dt=−u1+2u2du2dt=u1−2u2" role="presentation" style="position: relative;">du1dtdu2dt=−u1+2u2=u1−2u2$$\begin{array}{cc}\frac{d{u}_1}{dt}& =-u_1+2u_2\\ \frac{d{u}_2}{dt}& =u_1-2u_2\end{array}$$

\right. {dtdu1​​dtdu2​​​=−u1​+2u2​=u1​−2u2​​，也就是$\frac{d\mathbf{u}}{dt}=\mathbf{A}\mathbf{u}$

• 原方程有两个相互耦合(coupled)的未知函数$u_1,u_2$；
• 找出$\mathbf{A}$的特征值和特征向量（即对角化），可以实现解耦
  这对应了 最后的通解$\mathbf{u}(\mathrm{t})=c_1{e}^{{\lambda }_{1}t}{\mathbf{x}}_1+c_2{e}^{{\lambda }_{2}t}{\mathbf{x}}_2$有多个指数函数项，但他们之间互不相干；
  并且，各个指数函数前面的系数$c_i$完全由初始条件$u(0)$确定，一旦确定后系数恒定不变

---

进一步的，下面研究如何将解表示为特征值矩阵$\mathbf{\Lambda }$和特征向量矩阵$\mathbf{S}$的形式

对于$\frac{d\mathbf{u}}{dt}=\mathbf{A}\mathbf{u}$方程，原方程未知数$u_1,u_2$耦合，即$\mathbf{A}$不是对角阵
希望解耦，就是说希望将$\mathbf{A}$对角化（出现对角的系数矩阵意味着各未知数互不干扰）

• 希望将$\mathbf{u}$写为特征向量的线性组合形式$\mathbf{u}=\mathbf{Sv}$
  其中$\mathbf{v}$为新的未知量（代替原来的未知量$\mathbf{u}$），$\mathbf{S}$矩阵是$\mathbf{A}$的特征向量矩阵

• $\mathbf{u}=\mathbf{Sv}$带入原方程$\frac{d\mathbf{u}}{dt}=\mathbf{A}\mathbf{u}$，得到 S d v d t = A S v ⇒ d v d t = S − 1 A S v = Λ v

  Sdvdt=ASv⇒dvdt=S−1ASv=Λv" role="presentation">Sdvdt=ASv⇒dvdt=S−1ASv=Λv$$\begin{array}{l}\mathit{S}\frac{d\mathbf{v}}{dt}=\mathit{A}\mathit{S}\mathbf{v}\\ \Rightarrow \frac{d\mathbf{v}}{dt}={\mathit{S}}^{-1}\mathit{A}\mathit{S}\mathbf{v}=\mathbf{\Lambda }\mathbf{v}\end{array}$$
  Sdtdv​=ASv⇒dtdv​=S−1ASv=Λv​

可以理解为：

• 以特征向量为基（带入$\mathbf{u}=\mathbf{Sv}$，原未知量$\mathbf{u}$在新的基下的坐标为$\mathbf{v}$，即新的未知量为$\mathbf{v}$）
• 问题变为：求解关于新的未知量$\mathbf{v}$的对角化方程组$\frac{d\mathbf{v}}{dt}=\mathbf{\Lambda }\mathbf{v}$

新方程组的系数矩阵$\mathbf{\Lambda }$为对角阵，即方程组每一行都形如$\frac{d{v}_i}{dt}={\lambda }_i{v}_i$；
此时新方程组不存在耦合，或者说方程组各个未知量之间没有联系（系数矩阵为对角阵导致的）
新的 对角化方程组$\frac{d\mathbf{v}}{dt}=\mathbf{\Lambda }\mathbf{v}$的解为$$\mathbf{v}(t)=e^{\Lambda t}\mathbf{v}(0)$$

• 通过坐标变换/基变换，可以得到原方程$\frac{d\mathbf{u}}{dt}=\mathbf{A}\mathbf{u}$的解为$$\mathbf{u}(t)=\mathbf{S}{e}^{\mathbf{\Lambda }t}{\mathbf{S}}^{-1}\mathbf{u}(0)=e^{\mathbf{A}t}\mathbf{u}(0)$$也就是说方程$\frac{d\mathbf{u}}{dt}=\mathbf{A}\mathbf{u}$的解就是$\mathbf{u}(t)=e^{\mathbf{A}t}\mathbf{u}(0)$，其中$e^{\mathbf{A}t}=\mathbf{S}{e}^{\mathbf{\Lambda }t}{\mathbf{S}}^{-1}$，意义是坐标变换（下面将会证明）

结论：

• （前提：若$\mathbf{A}$有n个无关的特征向量方程）$\frac{d\mathbf{u}}{dt}=\mathbf{A}\mathbf{u}$的解解就是$\mathbf{u}(t)=e^{\mathbf{A}t}\mathbf{u}(0)$，其中$e^{\mathbf{A}t}=\mathbf{S}{e}^{\mathbf{\Lambda }t}{\mathbf{S}}^{-1}$，意义是坐标变换（这里的$\mathbf{\Lambda }$和$\mathbf{S}$来自对角化$\mathbf{A}=\mathbf{S}\mathbf{\Lambda }{\mathbf{S}}^{-1}$）
• 上面的通解$\mathbf{u}(\mathrm{t})=c_1{e}^{{\lambda }_{1}t}{\mathbf{x}}_1+c_2{e}^{{\lambda }_{2}t}{\mathbf{x}}_2$是这里的一般方法的一个特例，可以带入验证，$\mathbf{u}(t)=\mathbf{S}{e}^{\mathbf{\Lambda }t}{\mathbf{S}}^{-1}\mathbf{u}(0)=c_1{e}^{{\lambda }_{1}t}{\mathbf{x}}_1+c_2{e}^{{\lambda }_{2}t}{\mathbf{x}}_2$
  这里我们将解表示成了特征值矩阵$\mathbf{\Lambda }$和特征向量矩阵$\mathbf{S}$的形式，这里的通解更加通用（上面是这里的特例）

证明：矩阵的指数函数$e^{\mathbf{A}t}=\mathbf{S}{e}^{\mathbf{\Lambda }t}{\mathbf{S}}^{-1}$

我们希望研究$e^{\mathbf{A}t}$，这是个矩阵，但是其元素没有显式计算表示，而是需要通过下面的幂级数公式来计算！！！

目标：证明$e^{\mathbf{A}t}$与$\mathbf{A}$的特征值矩阵$\mathbf{\Lambda }$和特征向量矩阵$\mathbf{S}$的关系为$e^{\mathbf{A}t}=\mathbf{S}{e}^{\mathbf{\Lambda }t}{\mathbf{S}}^{-1}$

利用指数函数的幂级数公式$e^x={\sum }_{n=0}^{\infty }\frac{x^n}{n!}=1+x+\frac{x^2}{2}+\frac{x^3}{6}+\cdots$，可以将指数部分$e^x$变为幂次项$x^n$之和的形式，类比得到矩阵的指数函数$e^{\mathbf{A}t}$，可以写为$$e^{\mathbf{A}t}=I+\mathbf{A}t+\frac{(\mathbf{A}t{)}^2}{2}+\frac{(\mathbf{A}t{)}^3}{6}+\cdots$$利用${\mathbf{A}}^k=\mathbf{S}{\mathbf{\Lambda }}^k{\mathbf{S}}^{-1}$，可以得到 e A t = I + A t + ( A t ) 2 2 + ( A t ) 3 6 + ⋯ = S S − 1 + S Λ S − 1 t + S Λ 2 S − 1 2 t 2 + S Λ 3 S − 1 6 t 3 + ⋯ = S ( I + Λ t + Λ 2 2 t 2 + Λ 3 6 t 3 + ⋯   ) S − 1 = S e Λ t S − 1

eAt=I+At+(At)22+(At)36+⋯=SS−1+SΛS−1t+SΛ2S−12t2+SΛ3S−16t3+⋯=S(I+Λt+Λ22t2+Λ36t3+⋯)S−1=SeΛtS−1" role="presentation" style="position: relative;">eAt=I+At+(At)22+(At)36+⋯=SS−1+SΛS−1t+SΛ2S−12t2+SΛ3S−16t3+⋯=S(I+Λt+Λ22t2+Λ36t3+⋯)S−1=SeΛtS−1$$\begin{array}{rl}{e}^{\mathit{A}t}& =I+\mathit{A}t+\frac{(\mathit{A}t{)}^2}{2}+\frac{(\mathit{A}t{)}^3}{6}+\cdots \\ & =\mathit{S}{\mathit{S}}^{-1}+\mathit{S}\mathbf{\Lambda }{\mathit{S}}^{-1}t+\frac{\mathit{S}{\mathbf{\Lambda }}^2{\mathit{S}}^{-1}}{2}{t}^2+\frac{\mathit{S}{\mathbf{\Lambda }}^3{\mathit{S}}^{-1}}{6}{t}^3+\cdots \\ & =\mathit{S}\left(I+\mathbf{\Lambda }t+\frac{{\mathrm{\Lambda }}^2}{2}{t}^2+\frac{{\mathrm{\Lambda }}^3}{6}{t}^3+\cdots \right){\mathit{S}}^{-1}\\ & =\mathit{S}{e}^{\mathbf{\Lambda }t}{\mathit{S}}^{-1}\end{array}$$

eAt​=I+At+2(At)2​+6(At)3​+⋯=SS−1+SΛS−1t+2SΛ2S−1​t2+6SΛ3S−1​t3+⋯=S(I+Λt+2Λ2​t2+6Λ3​t3+⋯)S−1=SeΛtS−1​
注意，这里出现了${\mathbf{S}}^{-1}$，那么下面的一切成立的前提是，矩阵A必须具有n个线性无关的特征向量，从而矩阵才能对角化

最终，我们可以将$e^{\mathbf{A}t}$分解为 e A t = S e Λ t S − 1 , 其中 e Λ t = [ e λ 1 t 0 ⋯ 0 0 e λ 2 t 0 ⋮ ⋱ ⋮ 0 ⋯ 0 e λ n t ] e^{\mathbf At}=\boldsymbol{S} e^{\boldsymbol{\Lambda} t} \boldsymbol{S}^{-1},其中e^{\Lambda t}=\left[

\right] eAt=SeΛtS−1,其中eΛt=⎣ ⎡​eλ1​t0⋮0​0eλ2​t⋯​⋯⋱0​00⋮eλn​t​⎦ ⎤​

可见，若$\mathbf{A}$有n个无关的特征向量方程，则$e^{\mathbf{A}t}$与$\mathbf{A}$的特征值矩阵$\mathbf{\Lambda }$和特征向量矩阵$\mathbf{S}$的关系为$e^{\mathbf{A}t}=\mathbf{S}{e}^{\mathbf{\Lambda }t}{\mathbf{S}}^{-1}$（给出了$e^{\mathbf{A}t}$的显式计算式）

二阶微分方程

像之前的Fibonacci数列的例子，可以从二阶差分方程构造一阶的差分方程
同样的，给出二阶微分方程（同时出现了前后三项）$y^{\prime \prime }+b{y}^{\prime }+ky=0$，我们也可以增加一个方程，得到一个方程组（可表示为矩阵向量乘法），从而将整体列向量变为新的变量，得到一阶微分方程

构造方程组 { y ′ ′ = − b y ′ − k y y ′ = y ′ \left\{

\right. {y′′=−by′−kyy′=y′​并令新的未知量 u = [ y ′ y ] \mathbf{u}=\left[

y′y" role="presentation">y′y$$\begin{array}{l}{y}^{\mathrm{\prime }}\\ y\end{array}$$

\right] u=[y′y​]，则得到新方程组 [ y ′ ′ y ′ ] = [ − b − k 1 0 ] [ y ′ y ] 即 u ′ = [ − b − k 1 0 ] u \left[

y′′y′" role="presentation">y′′y′$$\begin{array}{l}{y}^{\mathrm{\prime }\mathrm{\prime }}\\ y^{\mathrm{\prime }}\end{array}$$

\right] =\left[

−b−k10" role="presentation">−b1−k0$$\begin{array}{cc}-b& -k\\ 1& 0\end{array}$$

\right]\left[

y′y" role="presentation">y′y$$\begin{array}{l}{y}^{\mathrm{\prime }}\\ y\end{array}$$

\right]即\mathbf{u}^{\prime}=\left[

−b−k10" role="presentation">−b1−k0$$\begin{array}{cc}-b& -k\\ 1& 0\end{array}$$

\right]\mathbf{u} [y′′y′​]=[−b1​−k0​][y′y​]即u′=[−b1​−k0​]u

可以推广到k阶微分方程，我们仍然添加(k-1)个额外方程，则${\mathbf{u}}^{\prime }=\mathbf{Au}$的系数矩阵$\mathbf{A}$是一个k x k矩阵，这样k阶微分方程可以转化为一阶微分方程