矩阵微分方程

微分方程是变量的未知函数的数学方程，将函数值与不同阶导数联系起来。矩阵微分方程包含多个函数，以向量形式堆叠在一起，并由一个矩阵将它们与导数联系起来。

例如，一阶矩阵常微分方程为

$\dot{𝐱}(t)=𝐀(t)𝐱(t)$

其中$𝐱(t)$是基变量$t$的函数的$n\times 1$向量，$\dot{𝐱}(t)$是函数的一阶导，$𝐀(t)$是$n\times n$系数矩阵。

若$𝐀$为常数且有n个线性无关的特征向量，微分方程有如下一般解：

$𝐱(t)=c_1{e}^{{\lambda }_{1}t}{𝐮}_1+c_2{e}^{{\lambda }_{2}t}{𝐮}_2+\cdots +c_n{e}^{{\lambda }_{n}t}{𝐮}_n,$

其中Template:Nowrap是A的特征值，Template:Nowrap是A相应的特征向量；Template:Nowrap为常数。

更一般地说，若$𝐀(t)$等于其积分$\underset{a}{\overset{t}{\int }}𝐀(s)ds$，则马格努斯展开降为前导阶，微分方程的一般解是

$𝐱(t)=e^{\underset{a}{\overset{t}{\int }}𝐀(s)ds}𝐜,$

其中$𝐜$是$n\times 1$常向量。

通过使用哈密尔顿–凯莱定理和类范德蒙矩阵，这种形式化的矩阵指数解可简化为一种简单的形式。下面，我们将用普策算法（Putzer's algorithm）来展示这一方法。

矩阵系统的稳定性与稳态

矩阵方程

$\dot{𝐱}(t)=𝐀𝐱(t)+𝐛$

当且仅当常数矩阵A的所有特征值的实部都为负，n×1参数常数向量b才稳定。

稳定时收敛到的稳态x*可置

${\dot{𝐱}}^{*}(t)=\mathrm{𝟎},$

找到，因此有

${𝐱}^{*}=-{𝐀}^{-1}𝐛,$

假设A可逆。

因此，原方程可用偏离稳态的齐次形式来写

$\dot{𝐱}(t)=𝐀[𝐱(t)-{𝐱}^{*}].$

一种等效的表达是，x*是非齐次方程的一个特解，而所有解的形式都是

${𝐱}_h+{𝐱}^{*},$

其中${𝐱}_h$是齐次方程（b=0）的解。

双状态变量情形的稳定性

n = 2（2个状态变量）时，稳定条件为：过渡矩阵A的两个特征值均有负实部，等价于A的迹为负、行列式为正。

矩阵形式的解

$\dot{𝐱}(t)=𝐀[𝐱(t)-{𝐱}^{*}]$的形式解为矩阵指数形式

$𝐱(t)={𝐱}^{*}+e^{𝐀t}[𝐱(0)-{𝐱}^{*}],$

使用多种技术中的任何一种进行评估。

计算Template:Math的普策算法

给定特征值为${\lambda }_1,{\lambda }_2,\dots ,{\lambda }_n$的矩阵A，

$e^{𝐀t}=\sum _{j=0}^{n-1}{r}_{j+1}\left(t\right){𝐏}_j$

其中

${𝐏}_0=𝐈$

${𝐏}_j=\prod _{k=1}^j(𝐀-{\lambda }_k𝐈)={𝐏}_{j-1}(𝐀-{\lambda }_j𝐈),j=1,2,\dots ,n-1$

${\dot{r}}_1={\lambda }_1{r}_1$

$r_1\left(0\right)=1$

${\dot{r}}_j={\lambda }_j{r}_j+r_{j-1},j=2,3,\dots ,n$

$r_j\left(0\right)=0,j=2,3,\dots ,n$

$r_i(t)$的方程是简单的一阶非齐次常微分方程。

注意该算法并不要求矩阵A可对角化，并绕过了通常使用的若尔当标准形的计算。

矩阵常微分方程解构示例

一阶齐次矩阵常微分方程包含两个函数x(t)、y(t)，从矩阵形式解出后有如下形式：

$\frac{dx}{dt}=a_{1}x+b_{1}y,\frac{dy}{dt}=a_{2}x+b_{2}y$

其中$a_1$、$a_2$、$b_1$、$b_2$可为任意标量。 高阶矩阵ODE的形式可能复杂得多。

解分解后的矩阵常微分方程

求解上述方程并找到这种特定阶次和形式的所需函数的过程大概分3步。每个步骤的简要说明如下：

• 找到特征值
• 找到特征向量
• 找到所需函数

第三部通常是把前两步的结果代入专门形式的一般方程中，下详。

矩阵ODE已解示例

Template:See also 要按上述3步解矩阵ODE，并在过程中使用简单矩阵，具体来说，现在下面的一阶齐次线性ODE中找到函数Template:Mvar、函数Template:Mvar，都用单一自变量Template:Mvar表示：

$\frac{dx}{dt}=3x-4y,\frac{dy}{dt}=4x-7y.$

要解这个常微分方程系统，在过程中的某时刻需要一组两个初始条件（对应起点的两个状态变量）。这时先取Template:Math。

第一步

第一步即找到A的特征值

$\left[\begin{array}{c}{x}^{\prime }\\ y^{\prime }\end{array}\right]=\left[\begin{array}{cc}3& -4\\ 4& -7\end{array}\right]\left[\begin{array}{c}x\\ y\end{array}\right].$

上面的导数记号x′等称为拉格朗日记法（由约瑟夫·拉格朗日提出，等同于前面方程里的dx/dt，这是莱布尼兹记法，得名于戈特弗里德·莱布尼茨）。

一旦两个变量的系数被写为上述矩阵形式A，就可估计特征值了。为此，可求矩阵行列式，即从上述系数矩阵中减去单位矩阵$I_n$乘常数Template:Mvar，再得到特征多项式

$\det (\left[\begin{array}{cc}3& -4\\ 4& -7\end{array}\right]-\lambda \left[\begin{array}{cc}1& 0\\ 0& 1\end{array}\right]),$

再解得其零点。

进一步简化、应用矩阵加法的基本规则，得出

$\det \left[\begin{array}{cc}3-\lambda & -4\\ 4& -7-\lambda \end{array}\right].$

应用求单一2×2矩阵行列式的规则，可得下列一元二次方程

$\det \left[\begin{array}{cc}3-\lambda & -4\\ 4& -7-\lambda \end{array}\right]=0$

$-21-3\lambda +7\lambda +{\lambda }^2+16=0$

可以进一步简化

${\lambda }^2+4\lambda -5=0.$

应用因式分解得到给定一元二次方程的两个根${\lambda }_1$、${\lambda }_2$

${\lambda }^2+5\lambda -\lambda -5=0$

$\lambda (\lambda +5)-1(\lambda +5)=0$

$(\lambda -1)(\lambda +5)=0$

$\lambda =1,-5.$

上面算出的${\lambda }_1=1$、${\lambda }_2=-5$即所求A的特征值。 矩阵ODE的特征值可能是复数，求解过程的下一步及最终形式和解法可能会有巨大变化。

第二步

第二步即找到A的特征向量。

对算出的每个特征值，都有单独的特征向量。例如对第一个特征值即${\lambda }_1=1$，有

$\left[\begin{array}{cc}3& -4\\ 4& -7\end{array}\right]\left[\begin{array}{c}\alpha \\ \beta \end{array}\right]=1\left[\begin{array}{c}\alpha \\ \beta \end{array}\right].$

应用矩阵乘法规则简化上式，得到

$3\alpha -4\beta =\alpha$

$\alpha =2\beta .$

所有计算都是为了得到最后一个式子，本例中就是Template:Math。现在任取一个无关紧要的小值（这样更容易处理），代入Template:Math中的Template:Mvar或Template:Mvar（选哪个并不重要），这样就得到了一个简单的向量，就是这个特定特征值所需的特征向量。在本例中，我们取Template:Math，得Template:Math。用标准的向量符号来写，向量是这样的

${\widehat{𝐯}}_1=\left[\begin{array}{c}2\\ 1\end{array}\right].$

对第二个特征值$\lambda =-5$进行相同的计算，得到第二个特征向量，结果为

${\widehat{𝐯}}_2=\left[\begin{array}{c}1\\ 2\end{array}\right].$

第三步

最后一步是找到“隐藏”在导数背后的所求函数。有两个函数，因为微分方程涉及两个变量。

方程包含之前得到的所有信息，形式如下：

$\left[\begin{array}{c}x\\ y\end{array}\right]=A{e}^{{\lambda }_{1}t}{\widehat{𝐯}}_1+B{e}^{{\lambda }_{2}t}{\widehat{𝐯}}_2.$

代入特征值和特征向量，得到

$\left[\begin{array}{c}x\\ y\end{array}\right]=A{e}^t\left[\begin{array}{c}2\\ 1\end{array}\right]+B{e}^{-5t}\left[\begin{array}{c}1\\ 2\end{array}\right].$

简化

$\left[\begin{array}{c}x\\ y\end{array}\right]=\left[\begin{array}{cc}2& 1\\ 1& 2\end{array}\right]\left[\begin{array}{c}A{e}^t\\ B{e}^{-5t}\end{array}\right].$

再简化，分别写出函数Template:Mvar、Template:Mvar的方程

$x=2A{e}^t+B{e}^{-5t}$

$y=A{e}^t+2B{e}^{-5t}.$

上述方程就是所求的一般函数，但只是一般形式（Template:Mvar、Template:Mvar的值未指定），但我们想找到它们的精确形式和解。因此现在，考虑问题的给定初始条件（即所谓初值问题）。假设给定了$x(0)=y(0)=1$，是ODE的起点；条件的应用指定了常数Template:Mvar、Template:Mvar。从条件$x(0)=y(0)=1$可以看出，Template:Math时，上述方程的左式等于1，由此可构造下列线性方程组

$1=2A+B$

$1=A+2B.$

求解这些等式，发现常数Template:Mvar、Template:Mvar都等于1/3。因此将这些值代入这两个函数的一般形式，就可以得到它们的精确形式 $$x=\frac{2}{3}{e}^t+\frac{1}{3}{e}^{-5t}$$ $$y=\frac{1}{3}{e}^t+\frac{2}{3}{e}^{-5t},$$ 所求的两个函数。

使用矩阵指数

上述问题可以直接应用矩阵指数法解决。也就是说，可以说

$\left[\begin{array}{c}x(t)\\ y(t)\end{array}\right]=\exp \left(\left[\begin{array}{cc}3& -4\\ 4& -7\end{array}\right]t\right)\left[\begin{array}{c}{x}_0(t)\\ y_0(t)\end{array}\right]$

给出了（可用MATLAB的expm工具包之类，或通过对角化，并利用对角矩阵的矩阵指数与元素的指数化相等这一特性来计算）

$\exp \left(\left[\begin{array}{cc}3& -4\\ 4& -7\end{array}\right]t\right)=\left[\begin{array}{cc}4e^t/3-e^{-5t}/3& 2e^{-5t}/3-2e^t/3\\ 2e^t/3-2e^{-5t}/3& 4e^{-5t}/3-e^t/3\end{array}\right]$

得到最终解

$\left[\begin{array}{c}x(t)\\ y(t)\end{array}\right]=\left[\begin{array}{cc}4e^t/3-e^{-5t}/3& 2e^{-5t}/3-2e^t/3\\ 2e^t/3-2e^{-5t}/3& 4e^{-5t}/3-e^t/3\end{array}\right]\left[\begin{array}{c}1\\ 1\end{array}\right]$

$\left[\begin{array}{c}x(t)\\ y(t)\end{array}\right]=\left[\begin{array}{c}{e}^{-5t}/3+2e^t/3\\ e^t/3+2e^{-5t}/3\end{array}\right]$

这与之前展示的特征向量方法相同。

另见

• 齐次微分方程
• 矩阵差分方程
• 冷却定律
• 斐波那契数列
• 差分方程
• 波动方程
• 自治系统 (数学)