选择 ODE 求解器

常微分方程

常微分方程 (ODE) 包含与一个自变量 t（通常称为时间）相关的因变量 y 的一个或多个导数。此处用于表示 y 关于 t 的导数的表示法对于一阶导数为 $y'$ ，对于二阶导数为 $y''$ ，依此类推。ODE 的阶数等于 y 在方程中出现的最高阶导数。

例如，这是一个二阶 ODE：

$y'' = 9 y$

在初始值问题中，从初始状态开始解算 ODE。利用初始条件 $y_{0}$ 以及要在其中求得答案的时间段 $(t_{0}, t_{f})$ ，以迭代方式获取解。在每一步，求解器都对之前各步的结果应用一个特定算法。在第一个这样的时间步，初始条件将提供继续积分所需的必要信息。最终结果是，ODE 求解器返回一个时间步向量 $t = [t_{0}, t_{1}, t_{2}, ..., t_{f}]$ 以及在每一步对应的解 $y = [y_{0}, y_{1}, y_{2}, ..., y_{f}]$ 。

ODE 的类型

MATLAB^® 中的 ODE 求解器可解算以下类型的一阶 ODE：

$y' = f (t, y)$ 形式的显式 ODE。
$M (t, y) y' = f (t, y)$ 形式的线性隐式 ODE，其中 $M (t, y)$ 为非奇异质量矩阵。该质量矩阵可以是时间或状态依赖的矩阵，也可以是常量矩阵。线性隐式 ODE 涉及在质量矩阵中编码的一阶 y 导数的线性组合。
线性隐式 ODE 可随时变换为显式形式 $y' = M^{- 1} (t, y) f (t, y)$ 。不过，将质量矩阵直接指定给 ODE 求解器可避免这种既不方便还可能带来大量计算开销的变换操作。
如果 $y'$ 的某些分量缺失，则这些方程称为微分代数方程或 DAE，并且 DAE 方程组会包含一些代数变量。代数变量是导数未出现在方程中的因变量。可通过对方程求导来将 DAE 方程组重写为等效的一阶 ODE 方程组，以消除代数变量。将 DAE 重写为 ODE 所需的求导次数称为微分指数。ode15s 和 ode23t 求解器可解算微分指数为 1 的 DAE。
$f (t, y, y') = 0$ 形式的完全隐式 ODE。完全隐式 ODE 不能重写为显式形式，还可能包含一些代数变量。ode15i 求解器专为完全隐式问题（包括微分指数为 1 的 DAE）而设计。

ODE 方程组

您可以指定需要解算的任意数量的 ODE 耦合方程，原则上，方程的数量仅受计算机可用内存的限制。如果方程组包含 n 个方程，

$(\begin{matrix} y'_{1} \\ y'_{2} \\ ⋮ \\ y'_{n} \end{matrix}) = (\begin{matrix} f_{1} (t, y_{1}, y_{2}, ..., y_{n}) \\ f_{2} (t, y_{1}, y_{2}, ..., y_{n}) \\ ⋮ \\ f_{n} (t, y_{1}, y_{2}, ..., y_{n}) \end{matrix}),$

则用于编写该方程组代码的函数将返回一个向量，其中包含 n 个元素，对应于 $y'_{1}, y'_{2}, \dots, y'_{n}$ 值。例如，考虑以下包含两个方程的方程组

${\begin{cases} y'_{1} = y_{2} \\ y'_{2} = y_{1} y_{2} - 2 . \end{cases}$

用于编写该方程组代码的函数为

function dy = myODE(t,y)
  dy(1) = y(2);
  dy(2) = y(1)*y(2)-2;
end

高阶 ODE

MATLAB ODE 求解器仅可解算一阶方程。您必须使用常规代换法，将高阶 ODE 重写为等效的一阶方程组

$\begin{array}{l} y_{1} = y \\ y_{2} = y' \\ y_{3} = y'' \\ ⋮ \\ y_{n} = y^{(n - 1)} . \end{array}$

这些代换将生成一个包含 n 个一阶方程的方程组

${\begin{cases} y'_{1} = y_{2} \\ y'_{2} = y_{3} \\ ⋮ \\ y'_{n} = f (t, y_{1}, y_{2}, ..., y_{n}) . \end{cases}$

例如，考虑三阶 ODE

$y''' - y'' y + 1 = 0.$

使用代换法

$\begin{array}{l} y_{1} = y \\ y_{2} = y' \\ y_{3} = y'' \end{array}$

生成等效的一阶方程组

${\begin{cases} y'_{1} = y_{2} \\ y'_{2} = y_{3} \\ y'_{3} = y_{1} y_{3} - 1. \end{cases}$

此方程组的代码则为

function dydt = f(t,y)
  dydt(1) = y(2);
  dydt(2) = y(3);
  dydt(3) = y(1)*y(3)-1;
end

复数 ODE

考虑复数 ODE 方程

$y' = f (t, y),$

其中 $y = y_{1} + i y_{2}$ 。为解算该方程，需要将实部和虚部分解为不同的解分量，最后重新组合相应的结果。从概念上讲，这类似于

$\begin{array}{l} y v = [Real (y) Imag (y)] \\ f v = [Real (f (t, y)) Imag (f (t, y))] . \end{array}$

例如，如果 ODE 为 $y' = y t + 2 i$ ，则可以使用函数文件来表示该方程：

function f = complexf(t,y)
  f = y.*t + 2*i;
end

然后，分解实部和虚部的代码为

function fv = imaginaryODE(t,yv)
  % Construct y from the real and imaginary components
  y = yv(1) + i*yv(2);            

  % Evaluate the function
  yp = complexf(t,y);             

  % Return real and imaginary in separate components
  fv = [real(yp); imag(yp)]; 
end

在运行求解器以获取解时，初始条件 y0 也会分解为实部和虚部，以提供每个解分量的初始条件。

y0 = 1+i;
yv0 = [real(y0); imag(y0)];
tspan = [0 2];
[t,yv] = ode45(@imaginaryODE, tspan, yv0);

获得解后，将实部和虚部分量组合到一起可获得最终结果。

y = yv(:,1) + i*yv(:,2);

基本求解器选择

ode45 适用于大多数 ODE 问题，一般情况下应作为您的首选求解器。但对于精度要求更宽松或更严格的问题而言，ode23、ode78、ode89 和 ode113 可能比 ode45 更加高效。

一些 ODE 问题具有较高的计算刚度或难度。术语“刚度”无法精确定义，但一般而言，当问题的某个位置存在标度差异时，就会出现刚度。例如，如果 ODE 包含的两个解分量在时间标度上差异极大，则该方程可能是刚性方程。如果非刚性求解器（例如 ode45）无法解算某个问题或解算速度极慢，则可以将该问题视为刚性问题。如果您观察到非刚性求解器的速度很慢，请尝试改用 ode15s 等刚性求解器。在使用刚性求解器时，可以通过提供雅可比矩阵或其稀疏模式来提高可靠性和效率。

您可以使用 ode 对象根据问题的属性自动选择求解器。如果您不确定要使用哪个求解器，请参考下表，其中提供关于每个求解器的适用情形的一般指导。

求解器	问题类型	精度	何时使用
`ode45`	非刚性	中	大多数情况下，您应当首先尝试求解器 `ode45`。
`ode23`		低	对于容差较宽松的问题或在刚度适中的情况下，`ode23` 可能比 `ode45` 更加高效。
`ode113`		低到高	对于具有严格误差容限的问题或在 ODE 函数需要大量计算开销的情况下，`ode113` 可能比 `ode45` 更加高效。
`ode78`		高	对于具有高准确度要求的平滑解的问题，`ode78` 可能比 `ode45` 更加高效。
`ode89`		高	对于非常平滑的问题，当在较长的时间区间内进行积分时，或当容差特别严格时，`ode89` 可能比 `ode78` 更加高效。
`ode15s`	刚性	低到中	若 `ode45` 失败或效率低下并且您怀疑面临刚性问题，请尝试 `ode15s`。此外，当解算微分代数方程 (DAE) 时，请使用 `ode15s`。
`ode23s`		低	对于误差容限较宽松的问题，`ode23s` 可能比 `ode15s` 更加高效。它可以解算一些刚性问题，而使用 `ode15s` 解算这些问题的效率不高。 `ode23s` 会在每一步计算雅可比，因此通过 `odeset` 提供雅可比有利于最大限度地提高效率和精度。如果存在质量矩阵，则它必须为常量矩阵。
`ode23t`		低	对于仅仅是刚度适中的问题，并且您需要没有数值阻尼的解，请使用 `ode23t`。 `ode23t` 可解算微分代数方程 (DAE)。
`ode23tb`		低	与 `ode23s` 一样，对于误差容限较宽松的问题，`ode23tb` 求解器可能比 `ode15s` 更加高效。
`ode15i`	完全隐式	低	对于完全隐式问题 f(t,y,y’) = 0 和微分指数为 1 的微分代数方程 (DAE)，请使用 `ode15i`。

有关何时使用每种求解器的详细信息和更多建议，请参阅 [5]。

ODE 示例和文件摘要

有几个示例文件可用作大多数 ODE 问题的有用起点。要运行微分方程示例应用，以便轻松浏览和运行示例，请键入

odeexamples

要打开单独的示例文件进行编辑，请键入

edit exampleFileName.m

要运行示例，请键入

exampleFileName

此表包含可用的 ODE 和 DAE 示例文件及其使用的求解器和选项的列表。其中包含示例子集的链接，这些示例也已直接发布在文档中。

示例文件	使用的求解器	指定的选项	描述	文档链接
`amp1dae`	`ode23t`	`'Mass'`	刚性 DAE - 包含常量奇异质量矩阵的电路	求解刚性晶体管微分代数方程
`ballode`	`ode23`	`'Events'` `'OutputFcn'` `'OutputSel'` `'Refine'` `'InitialStep'` `'MaxStep'`	简单事件位置 - 弹球	ODE 事件位置
`batonode`	`ode45`	`'Mass'`	时间依赖和状态依赖的质量矩阵的 ODE - 短棒的移动	求解抛向空中的短棒的运动方程
`brussode`	`ode15s`	`'JPattern'` `'Vectorized'`	刚性大问题 - 化学反应中的扩散 (Brusselator)	解算刚性 ODE
`burgersode`	`ode15s`	`'Mass'` `'MStateDependence'` `'JPattern'` `'MvPattern'` `'RelTol'` `'AbsTol'`	强烈依赖于状态的质量矩阵的 ODE - 使用移动网格方法解算伯格斯方程	求解具有强状态依赖质量矩阵的 ODE
`fem1ode`	`ode15s`	`'Mass'` `'MStateDependence'` `'Jacobian'`	时间依赖的质量矩阵的刚性问题 - 有限元方法	—
`fem2ode`	`ode23s`	`'Mass'`	常量质量矩阵的刚性问题 - 有限元方法	—
`hb1ode`	`ode15s`	—	在非常长的区间内解算的刚性 ODE 问题 - 罗伯逊化学反应	—
`hb1dae`	`ode15s`	`'Mass'` `'RelTol'` `'AbsTol'` `'Vectorized'`	基于守恒定律的刚性线性隐式 DAE - 罗伯逊化学反应	将罗伯逊问题作为半显式微分代数方程 (DAE) 求解
`ihb1dae`	`ode15i`	`'RelTol'` `'AbsTol'` `'Jacobian'`	刚性完全隐式 DAE - 罗伯逊化学反应	将罗伯逊问题作为隐式微分代数方程 (DAE) 求解
`iburgersode`	`ode15i`	`'RelTol'` `'AbsTol'` `'Jacobian'` `'JPattern'`	隐式 ODE 方程组 - 伯格斯方程	—
`kneeode`	`ode15s`	`'NonNegative'`	具有非负约束的“膝盖问题”	非负 ODE 解
`orbitode`	`ode45`	`'RelTol'` `'AbsTol'` `'Events'` `'OutputFcn'`	高级事件位置 - 限制性三体问题	ODE 事件位置
`rigidode`	`ode45`	—	非刚性问题 - 刚体在不受外力作用时的欧拉方程	求解非刚性 ODE
`vdpode`	`ode15s`	`'Jacobian'`	可参数化的 van der Pol 方程（对于较大 μ，为刚性问题）	解算刚性 ODE

参考

[1] Shampine, L. F. and M. K. Gordon, Computer Solution of Ordinary Differential Equations: the Initial Value Problem, W. H. Freeman, San Francisco, 1975.

[2] Forsythe, G., M. Malcolm, and C. Moler, Computer Methods for Mathematical Computations, Prentice-Hall, New Jersey, 1977.

[3] Kahaner, D., C. Moler, and S. Nash, Numerical Methods and Software, Prentice-Hall, New Jersey, 1989.

[4] Shampine, L. F., Numerical Solution of Ordinary Differential Equations, Chapman & Hall, New York, 1994.

[5] Shampine, L. F. and M. W. Reichelt, “The MATLAB ODE Suite,” SIAM Journal on Scientific Computing, Vol. 18, 1997, pp. 1–22.

[6] Shampine, L. F., Gladwell, I. and S. Thompson, Solving ODEs with MATLAB, Cambridge University Press, Cambridge UK, 2003.

另请参阅

ode | odeset | odextend

外部网站

常微分方程