简 介：本文将研究微分方程的解析解算法，介绍在MATLAB 环境中如何用微分方程求解函数直接得出线性微分方程组的解析解，并对一阶简单的非线性微分方程的解析解求解进行探讨，从而得出结论，一般非线性微分方程是没有解析解的。

关键词：微分方程，解析解，MATLAB

§01 总述

---

假设已知常系数线性微分方程的一般描述方法为
$$\begin{gathered} \frac{{\mathrm{d}}^{n}y(t)}{\mathrm{d}{t}^n}+a_1\frac{{\mathrm{d}}^{n-1}y(t)}{\mathrm{d}{t}^{n-1}}+a_2\frac{{\mathrm{d}}^{n-2}y(t)}{\mathrm{d}{t}^{n-2}}+\cdots +a_{n-1}\frac{\mathrm{d}y(t)}{\mathrm{d}t}+a_{n}y(t) \\ =b_1\frac{{\mathrm{d}}^{m}u(t)}{\mathrm{d}{t}^m}+b_2\frac{{\mathrm{d}}^{m-1}u(t)}{\mathrm{d}{t}^{m-1}}+\cdots +b_m\frac{\mathrm{d}u(t)}{\mathrm{d}t}+b_{m+1}u(t) \end{gathered}$$

其中，$a_i$, $b_i$均为常数，对零初值问题有 $\mathcal{L}\left[{\mathrm{d}}^{m}y(t)/\mathrm{d}{t}^m\right]=s^m\mathcal{L}[y(t)]$，可以对应得出下面的多项式代数方程

$$s^n+a_1{s}^{n-1}+a_2{s}^{n-2}+\cdots +a_{n-1}s+a_n=0$$

假设代数方程的特征根$s_i$均可以求出，且假设它们均相异，则可以得出原微分方程的解析解一般形式为

$$y(t)=C_1{\mathrm{e}}^{r_{1}t}+C_2{\mathrm{e}}^{r_{2}t}+\cdots +C_n{\mathrm{e}}^{r_{n}t}+\gamma (t)$$

其中，$C_i$为待定系数，而 $\gamma (t)$是满足 $u(t)$输入的一个特解。$s_i$有重根的情况也有相应的解析解形式。

从得出的代数方程看，由著名的Abel-Ruffini定理可知，4次及以下的多项式代数方程是能求出根的解析解的，故可以得出结论，低阶线性微分方程有一般意义下的解析解，结合多项式方程的数值解法可以得出一般高次多项式代数方程的数值解法，即高阶线性微分方程的准解析解方法。本文将介绍用 MATLAB 语言及其符号运算工具箱求解线性常系数微分方程解析解的方法。

§02 函数调用格式

---

y = dsolve(fun1, fun2, ... , funm)  % 默认自变量t
y = dsolve(fun1, fun2, ... , funm, 'x')   % 指明自变量x

• 可以同时求解多个方程、已知条件
• 注意自变量设置，否则可能得出无用的结果

§03 应用举例

---

例1：微分方程解析解

题目： 假设输入信号为$u(t)={\mathrm{e}}^{-5t}\cos (2t+1)+5$，试求下面微分方程的通解：

$$y^{(4)}(t)+10y^{(3)}(t)+35\ddot{y}(t)+50\dot{y}(t)+24y(t)=5\ddot{u}(t)+4\dot{u}(t)+2u(t)$$

求解：

syms y(t);
u = exp(-5*t)*cos(2*t+1)+5;
uu = 5*diff(u,t,2) + 4*diff(u,t) + 2*u;
f(t) = diff(y,t,4) + 10*diff(y,t,3) + 35*diff(y,t,2) + 50*diff(y,t) + 24*y;
ySol = dsolve(f(t)==uu);
ySol = simplify(ySol)

运行程序可得：

sy =
 
C1*exp(-4*t) - (547*exp(-5*t)*sin(2*t + 1))/520 - (343*exp(-5*t)*cos(2*t + 1))/520 + C2*exp(-3*t) + C3*exp(-2*t) + C4*exp(-t) + 5/12

化简可得 解的数学形式：

$$y(t)=\frac{5}{12}-\frac{343}{520}{\mathrm{e}}^{-5t}\cos (2t+1)-\frac{547}{520}{\mathrm{e}}^{-5t}\sin (2t+1)+C_1{\mathrm{e}}^{-4t}+C_2{\mathrm{e}}^{-3t}+C_3{\mathrm{e}}^{-2t}+C_4{\mathrm{e}}^{-t}$$

解析解的检验：

diff(ySol,4) + 10*diff(ySol,3) + 35*diff(ySol,2) + 50*diff(ySol) + 24*ySol - uu;
simplify(ans)

运行程序可得：

ans =
 
0

已知初值条件：
$$y(0)=3,\dot{y}(0)=2,\ddot{y}(0)=y^{(3)}(0)=0$$

syms y(t);
u = exp(-5*t)*cos(2*t+1)+5;
uu = 5*diff(u,t,2) + 4*diff(u,t) + 2*u;
f(t) = diff(y,t,4) + 10*diff(y,t,3) + 35*diff(y,t,2) + 50*diff(y,t) + 24*y;
Dy = diff(y,t); D2y = diff(y,t,2); D3y = diff(y,t,3);
cond = [y(0)==3, Dy(0)==2, D2y(0)==0, D3y(0)==0];
ySol= dsolve(f(t)==uu, cond);
ySol = simplify(ySol) 
ezplot(ySol,[0,5])

▲ 图1 已知初值条件微分方程解的图像

复杂边界条件：

$$y(0)=1/2,\dot{y}(\pi )=1,\ddot{y}(2\pi )=0,\dot{y}(2\pi )=1/5$$

syms y(t);
u = exp(-5*t)*cos(2*t+1)+5;
uu = 5*diff(u,t,2) + 4*diff(u,t) + 2*u;
f(t) = diff(y,t,4) + 10*diff(y,t,3) + 35*diff(y,t,2) + 50*diff(y,t) + 24*y;
Dy = diff(y,t); D2y = diff(y,t,2);
cond = [y(0)==1/2, Dy(pi)==1, D2y(2*pi)==0, Dy(2*pi)==1/5];
ySol= dsolve(f(t)==uu, cond);
ySol = vpa(ySol);
ySol = simplify(ySol) 
ezplot(ySol,[0,5])

▲ 图2 复杂边界条件微分方程解的图像

例2：求解微分方程组

题目： 试求解线性微分方程组的解析解

$$\{\begin{array}{l}{x}^{\prime \prime }(t)+2x^{\prime }(t)=x(t)+2y(t)-{\mathrm{e}}^{-t}\\ \\ y^{\prime }(t)=4x(t)+3y(t)+4{\mathrm{e}}^{-t}\end{array}$$

syms x(t) y(t)
eqn1 = diff(x,t,2) + 2*diff(x,t) == x + 2*y(t) - exp(-t);
eqn2 = diff(y,t) == 4*x + 3*y(t) + 4*exp(-t);
[xSol, ySol] = dsolve([eqn1, eqn2])

例3：变系数微分方程

题目： 变系数微分方程

$$(2x+3{)}^3{y}^{\prime \prime \prime }+3(2x+3)y^{\prime }-6y=0$$

syms y(x)
eqn = (2*x+3)^3*diff(y,x,3) + 3*(2*x+3)*diff(y,x) - 6*y == 0;
ySol = dsolve(eqn)

例4：多维微分方程组

题目：

$$\{\begin{array}{l}{x}^{\prime \prime }-x+y+z=0\\ \\ x+y^{\prime \prime }-y+z=0\\ \\ x+y+z^{\prime \prime }-z=0\end{array}$$

初值：$x(0)=1,y(0)=z(0)=x^{\prime }(0)=y^{\prime }(0)=z^{\prime }(0)=0$

syms x(t) y(t) z(t)
eqn1 = diff(x,t,2) - x + y + z == 0;
eqn2 = x + diff(y,t,2) - y + z == 0;
eqn3 = x + y + diff(z,t,2) - z == 0;
Dx = diff(x, t); Dy = diff(y, t); Dz = diff(z, t); 
cond = [x(0)==1, y(0)==0, z(0)==0, Dx(0)==0, Dy(0)==0, Dz(0)==0];
[xSol, ySol, zSol] = dsolve([eqn1, eqn2, eqn3], cond)