3.6.2

龙格-

库塔方法

改进的欧拉法比欧拉法精度高的原因在于，它在确定平均斜率时，多取了一个点的斜

率值。这样，如果我们在[Xi,X(i+1)]上多取几个点的斜率值，然后对它们作线性组合得到平均

斜率，则有可能构造出精度更高的计算方法。这就是龙格-库塔法的基本思想。龙格-库塔

法可看作是欧拉法思想的提高，属于精度较高的单步法。

龙格-库塔法是求解常微分方程初值问题的最重要的方法之一。MATLAB中提供了几

个采用龙格-库塔法来求解常微分方程的函数，即ode23，ode45，ode113

，ode23s

，ode15s

等，其中最常用的函数是

ode23(

二三阶龙格-库塔函数)和ode45(

四五阶龙格-库塔函数)，

下面分别对它们进行介绍。

1

．二三阶龙格-

库塔函数(ode23)

函数

ode23

的调用格式如下：

(1)

[T,Y]=ODE23('F',TSPAN,Y0)

输入参数中的'F'

是一个字符串，表示微分方程的形

式，也可以是

f

(x

,

y

)的M

文件。TSPAN=[T0

TFINAL]表示积分区间，Y0表示初始条件。

函数

ode23

表示在初始条件

Y0下从

T0到TFINAL

对微分方程

'(,)

yFty

=

进行积分。函数

F(T,

Y)

必须返回一列向量，两个输出参数是列向量

T

与矩阵

Y，其中向量

T

包含估计响应

的积分点，而矩阵

Y

的行数与向量

T

的长度相等。向量

T

中的积分点不是等间距的，这是

为了保持所需的相对精度，而改变了积分算法的步长。为了获得在确定点T0,T1,

"的解，

TSPAN=[T0

T1

TFINAL]

。需要注意的是：TSPAN中的点必须是单调递增或单调递减的。

(2)

[T,Y]=ODE23('F',TSPAN,Y0,OPTIONS)

其中，参数

options

为积分参数，它可由函

数ODESET

来设置。Options参数最常用的是相对误差‘RelTol’(

默认值是

1e-3)和绝对误差

‘AbsTol’(默认值是

1e-6)，其他参数同上。

(3)

[T,Y]=ODE23('F',TSPAN,Y0,OPTIONS,P1,P2,…)

参数P1,P2,

…可直接输入到函数

F

中去．如

F(T,Y,FLAG,P1,P2,…)。如果参数

OPTIONS为空，则输入

OPTIONS=[

]。也可

以在

ODE文件中(可参阅

ODEFILE函数)指明参数

TSPAN、Y0和OPTIONS的值。如果参

数TSPAN

或Y0

是空，则ODE23函数通过调用ODE文件[TSPAN,

Y0,

OPTIONS]

=

F([

],[

],

'init

')来获得

ODE23函数没有被提供的自变量值。如果获得的自变量表示空，则函

数ODE23会忽略，此时为

ODE23('F')。

(4)

[T,Y,TE,YE,IE]=ODE23('F',TSPAN,Y0,OPTIONS)

此时要求在参数

options

中的事

件属性设为'on'

，ODE文件必须被标记，以便

P(T,Y,'events')

能返回合适的信息，详细可参

阅函数

ODEFILE。输出参数中的

TE是一个列向量，矩阵

YE的行与列向量

TE中元素相

对应，向量

IE

表示解的索引。

2

．四五阶龙格-

库塔函数(ode45)

函数

ode45

的调用格式同

ode23

相同，其差别在于内部算法不同。如果'F'

为向量函数，

则ode23

和ode45

也可用来解微分方程组。

【例3.47

】

分别用二三阶龙格-库塔法和四五阶龙格-库塔法解常微分方程的初值问题：

解：先将微分方程写成自定义函数

exam2fun.m

function

f=exam2fun

(x,y)

f=-y-x*y.^2;

f=f(:);

然后在命令窗口输入以下语句：

>>

[x1,y1]=ode23('exam2fun',[0:0.1:1],1)

x1

=

0

0.1000

0.2000

0.3000

0.4000

0.5000

0.6000

0.7000

0.8000

0.9000

1.0000

y1

=

1.0000

0.9006

0.8046

0.7144

0.6314

0.5563

0.4892

0.4296

0.3772

0.3312

0.2910

>>

[x2,y2]=ode45('exam2fun',[0:0.1:1],1)

x2

=

0

0.1000

0.2000

0.3000

0.4000

0.5000

0.6000

0.7000

0.8000

0.9000

1.0000

y2

=

1.0000

0.9006

0.8046

0.7144

0.6315

0.5563

0.4892

0.4296

0.3772

0.3312

0.2910

温馨提示：答案为网友推荐，仅供参考