ode23s

在对求解器的调用中，可将只有一个解分量的简单 ODE 指定为匿名函数。该匿名函数必须同时接受两个输入 (t,y)，即使在该函数中一个输入未使用也是如此。

解算 ODE

指定时间区间 [0 2] 和初始条件 y0 = 1。

tspan = [0 2];
y0 = 1;
[t,y] = ode23s(@(t,y) -10*t, tspan, y0);

对解绘图。

plot(t,y,'-o')

表示张弛振荡的范德波尔方程组是一个刚性方程组示例。范围环的某些区域中的解分量变化缓慢且问题呈刚性，这些区域与变化显著的非刚性区域交替出现。

方程组为：

初始条件为 和 。函数 vdp1000 随 MATLAB® 一起提供，用于对方程进行编码。

function dydt = vdp1000(t,y)
%VDP1000  Evaluate the van der Pol ODEs for mu = 1000.
%
%   See also ODE15S, ODE23S, ODE23T, ODE23TB.

%   Jacek Kierzenka and Lawrence F. Shampine
%   Copyright 1984-2014 The MathWorks, Inc.

dydt = [y(2); 1000*(1-y(1)^2)*y(2)-y(1)];

使用 ode45 以及默认的相对和绝对误差容限（分别为 1e-3 和 1e-6）解算此方程组特别慢，需要好几分钟才能完成解算和绘制。ode45 需要几百万个时间步才能完成积分，因为它要在刚度区域满足容差要求。

下面是通过 ode45 获得的解的绘图，需要很长时间来计算。请注意，通过刚度区域需要大量时间步。

[t,y] = ode23s(@vdp1000,[0 3000],[2 0]);
plot(t,y(:,1),'-o')

ode23s 仅适用于使用两个输入参量（t 和 y）的函数。但是，通过在函数外部定义参数并在指定函数句柄时传递这些参数，可以传入额外参数。

解算 ODE

将该方程重写为一阶方程组可以得到

odefcn.m 将此方程组表示为接受四个输入参量（t、y、A 和 B）的函数。

function dydt = odefcn(t,y,A,B)
dydt = zeros(2,1);
dydt(1) = y(2);
dydt(2) = (A/B)*t.*y(1);

A = 1;
B = 2;
tspan = [0 5];
y0 = [0 0.01];
[t,y] = ode23s(@(t,y) odefcn(t,y,A,B), tspan, y0);

plot(t,y(:,1),'-o',t,y(:,2),'-.')

ode15s 求解器是适用于大多数刚性问题的首选求解器。但是，对于特定类型的问题，其他刚性求解器可能更高效。本示例使用所有四个刚性 ODE 求解器解算刚性测试方程。

请考虑以下测试方程：

随着 $$\lambda$$ 的量级增加，方程的刚性逐渐增强。使用 $$\lambda =1\times 10^9$$、初始条件 $$y(0)=1$$ 和时间区间 [0 0.5]。这些值使该问题足够刚性，从而使 ode45 和 ode23 需要对该方程进行积分。此外，还使用 odeset 传入常量雅可比$$J=\frac{\partial f}{\partial y}=-\lambda$$ 并打开求解器统计信息的显示。

lambda = 1e9;
y0 = 1;
tspan = [0 0.5];
opts = odeset('Jacobian',-lambda,'Stats','on');

使用 ode15s、ode23s、ode23t 和 ode23tb 对方程求解。生成子图进行比较。

subplot(2,2,1)
tic, ode15s(@(t,y) -lambda*y, tspan, y0, opts), toc

104 successful steps
1 failed attempts
212 function evaluations
0 partial derivatives
21 LU decompositions
210 solutions of linear systems
Elapsed time is 2.131028 seconds.

title('ode15s')
subplot(2,2,2)
tic, ode23s(@(t,y) -lambda*y, tspan, y0, opts), toc

63 successful steps
0 failed attempts
191 function evaluations
0 partial derivatives
63 LU decompositions
189 solutions of linear systems
Elapsed time is 0.302041 seconds.

title('ode23s')
subplot(2,2,3)
tic, ode23t(@(t,y) -lambda*y, tspan, y0, opts), toc

95 successful steps
0 failed attempts
125 function evaluations
0 partial derivatives
28 LU decompositions
123 solutions of linear systems
Elapsed time is 0.241553 seconds.

title('ode23t')
subplot(2,2,4)
tic, ode23tb(@(t,y) -lambda*y, tspan, y0, opts), toc

71 successful steps
0 failed attempts
167 function evaluations
0 partial derivatives
23 LU decompositions
236 solutions of linear systems
Elapsed time is 0.278277 seconds.

title('ode23tb')

这些刚性求解器都表现不错，但对于这个特定问题，ode23s 完成积分所用的步长最少，运行速度最快。由于指定了常量雅可比，任何求解器都不需要计算偏导数即可计算出解。指定雅可比对 ode23s 最有利，因为它通常需要在每个步长中计算雅可比。

对于一般的刚性问题，刚性求解器的性能因问题的形式和指定的选项而异。提供雅可比矩阵或稀疏模式始终可以提高求解器解决刚性问题的效率。但是，由于刚性求解器使用雅可比的方式不同，计算效率的提高程度也有很大差异。实际上，如果方程组非常大，或者要解算很多次，则非常值得研究一下不同求解器的性能，以最大程度地缩短执行时间。