常微分方程数值解
第九章 常微方程数值解法

第8章 序
许多科学技术问题,例如天文学中的星体运动, 许多科学技术问题,例如天文学中的星体运动,空间 技术中的物体飞行,自动控制中的系统分析, 技术中的物体飞行,自动控制中的系统分析,力学中的振 动,工程问题中的电路分析等,都可归结为常微分方程的 工程问题中的电路分析等, 初值问题。 初值问题。 所谓初值问题, 所谓初值问题,是函数及其必要的导数在积分的起始 点为已知的一类问题,一般形式为: 点为已知的一类问题,一般形式为:
⇒ y n +1 = y n −1 + 2hf ( xn , y n )
第9章 常微分方程数值解法
(8 - 4)
8-10
Euler公式的推导( Euler公式的推导(续5) 公式的推导
上对y )=f 四、利用数值积分公式:在[xn,xn+1]上对y′(x)=f (x,y(x)) 积分 利用数值积分公式:
x0 < x1 < L < xn < L
(i=1,2,…,n)构造插值函数作为近似函数。上述离散点 i=1,2,…,n)构造插值函数作为近似函数。 相 邻两点间的距离h 称为步长, 邻两点间的距离hi=xi-1-xi 称为步长,若hi 都相等为一定数 h, 则称为定步长,否则为变步长。( x, y ( x)) 则称为定步长,否则为变步长。 a≤ x≤b y ′( x) = f 本章重点讨论如下 y (a ) = y0 一阶微分方程: 一阶微分方程: 在此基础上介绍一阶微分方程组与 8-5 第9章 常微分方程数值解法 高阶微分方程的数值解法。 高阶微分方程的数值解法。
⇒ yn +1 = yn + hf ( xn , yn ) + E ( xn , h) ⇒ yn +1 = yn + hf ( xn , yn )
常微分方程的数值解

f ( x, y1 ) f ( x, y2 ) L y1 y2
(其中 L 为 Lipschitz 常数)则初值问题( 1 )存 在唯一的连续解。
求问题(1)的数值解,就是要寻找解函数在一 系列离散节点x1 < x2 <……< xn < xn+1 上的近似 值y1, y 2,…,yn 。 为了计算方便,可取 xn=x0+nh,(n=0,1,2,…), h称为步长。
(1),(2)式称为初值问题,(3)式称为边值问题。 在实际应用中还经常需要求解常微分方程组:
f1 ( x, y1 , y2 ) y1 ( x0 ) y10 y1 (4) f 2 ( x, y1 , y2 ) y2 ( x0 ) y20 y2
本章主要研究问题(1)的数值解法,对(2)~(4)只 作简单介绍。
得 yn1 yn hf ( xn1 , yn1 )
上式称后退的Euler方法,又称隐式Euler方法。 可用迭代法求解
二、梯形方法 由
y( xn1 ) y( xn )
xn1 xn
f ( x, y( x))dx
利用梯形求积公式: x h x f ( x, y( x))dx 2 f ( xn , y( xn )) f ( xn1 , y( xn1 ))
常微分方程的数言 简单的数值方法 Runge-Kutta方法 一阶常微分方程组和高阶方程
引言
在高等数学中我们见过以下常微分方程:
y f ( x, y, y) a x b y f ( x, y ) a x b (2) (1) (1) y ( x ) y , y ( x ) y 0 0 0 0 y ( x0 ) y0 y f ( x, y, y) a x b (3) y(a) y0 , y(b) yn
求常微分方程的数值解

求常微分方程的数值解一、背景介绍常微分方程(Ordinary Differential Equation,ODE)是描述自然界中变化的数学模型。
常微分方程的解析解往往难以求得,因此需要寻找数值解来近似地描述其行为。
求解常微分方程的数值方法主要有欧拉法、改进欧拉法、龙格-库塔法等。
二、数值方法1. 欧拉法欧拉法是最简单的求解常微分方程的数值方法之一。
它基于导数的定义,将微分方程转化为差分方程,通过迭代计算得到近似解。
欧拉法的公式如下:$$y_{n+1}=y_n+f(t_n,y_n)\Delta t$$其中,$y_n$表示第$n$个时间步长处的函数值,$f(t_n,y_n)$表示在$(t_n,y_n)$处的导数,$\Delta t$表示时间步长。
欧拉法具有易于实现和理解的优点,但精度较低。
2. 改进欧拉法(Heun方法)改进欧拉法又称Heun方法或两步龙格-库塔方法,是对欧拉法进行了精度上提升后得到的一种方法。
它利用两个斜率来近似函数值,并通过加权平均来计算下一个时间步长处的函数值。
改进欧拉法的公式如下:$$k_1=f(t_n,y_n)$$$$k_2=f(t_n+\Delta t,y_n+k_1\Delta t)$$$$y_{n+1}=y_n+\frac{1}{2}(k_1+k_2)\Delta t$$改进欧拉法比欧拉法精度更高,但计算量也更大。
3. 龙格-库塔法(RK4方法)龙格-库塔法是求解常微分方程中最常用的数值方法之一。
它通过计算多个斜率来近似函数值,并通过加权平均来计算下一个时间步长处的函数值。
RK4方法是龙格-库塔法中最常用的一种方法,其公式如下:$$k_1=f(t_n,y_n)$$$$k_2=f(t_n+\frac{\Delta t}{2},y_n+\frac{k_1\Delta t}{2})$$ $$k_3=f(t_n+\frac{\Delta t}{2},y_n+\frac{k_2\Delta t}{2})$$ $$k_4=f(t_n+\Delta t,y_n+k_3\Delta t)$$$$y_{n+1}=y_n+\frac{1}{6}(k_1+2k_2+2k_3+k_4)\Delta t$$三、数值求解步骤对于给定的常微分方程,可以通过以下步骤求解其数值解:1. 确定初值条件:确定$t=0$时刻的函数值$y(0)$。
常微分方程的数值解算法

常微分方程的数值解算法常微分方程的数值解算法是一种对常微分方程进行数值计算的方法,这可以帮助我们更好地理解和研究自然现象和工程问题。
在本文中,我们将介绍一些常用的数值解算法,探讨它们的优缺点和适用范围。
常微分方程(ODE)是描述自然现象和工程问题的重要数学工具。
然而,对于许多ODE解析解是无法求出的,因此我们需要通过数值方法对其进行求解。
常微分方程可以写作:y' = f(t, y)其中,y是函数,f是给定的函数,表示y随t的变化率。
这个方程可以写成初始值问题(IVP)的形式:y'(t) = f(t,y(t)),y(t0) = y0其中,y(t0)=y0是方程的初始条件。
解决IVP问题的典型方法是数值方法。
欧拉方法欧拉方法是最简单的一阶数值方法。
在欧拉方法中,我们从初始条件开始,并在t = t0到t = tn的时间内,用以下公式逐步递推求解:y n+1 = y n + hf (t n, y n)其中,f(t n,y n)是点(t n,y n)处的导数, h = tn - tn-1是时间间隔。
欧拉方法的优点是简单易懂,容易实现。
然而,它的缺点是在整个时间段上的精度不一致。
程度取决于使用的时间间隔。
改进的欧拉方法如果我们使用欧拉方法中每个时间段的中间点而不是起始点来估计下一个时间点,精度就会有所提高。
这个方法叫做改进的欧拉方法(或Heun方法)。
公式为:y n+1 = y n + h½[f(t n, y n)+f(tn+1, yn + h f (tn, yn))]这是一个二阶方法,精度比欧拉方法高,但计算量也大一些。
对于易受噪声干扰的问题,改进的欧拉方法是个很好的选择。
Runge-Kutta方法Runge-Kutta方法是ODE计算的最常用的二阶和高阶数值方法之一。
这个方法对定义域内的每个点都计算一个导数。
显式四阶Runge-Kutta方法(RK4)是最常用的Runge-Kutta方法之一,并已得到大量实践的验证。
常微分方程数值解法

介绍常微分方程数值解法常微分方程(ordinary differential equations,ODE)可用于描述许多日常存在的物理系统。
处理ODE问题常常被称为数值求解法,这指的是找到概括ODE或者其他适用于数学模型的解决方案来模括这些ODE。
这种解决方案可能在一系列不同方案中发挥重要作用,以此来提供更好的解释和预测。
常微分方程与几何图形更为相关,它利用二维或者三维空间中曲线的绘制以及分析。
通过引入一些不同的方法,可以对不同的常微分方程中的量进行描述,使得可以通过数值方法的解析来进行研究。
数值解法可能是时间消耗较多的,但有助于验证几何图形中的某些过程,以此帮助揭示数学模型。
四种常见的常微分方程数值解法四种常见的常微分方程数值解法是:前向差分法、向后差分法、中点法和全分方法。
•前向差分法:前向差分法的基本概念是利用ODE的特定解来表达时间步的影响。
这是一种基本的数值法,可以在ODE中确定任意位置的点作为终点。
在这里,任何这样的点都可以表示为ODE右边的时间步。
•向后差分法:它是反过来基于前向差分法。
它要求对ODE中的时间步进行逆向推导,以获得某一特定点的解。
向后差分法要求推导反向解中点,以便可以从每一步中获取该点的解。
•中点法:这是一种非常基本的数值解法,可以用来求解ODE中的某一步的解,但不具有直观的方法解释。
主要的思想是在每一次时间步中通过求出ODE的中点来寻找解。
•全分方法:这是一种更复杂的数值解法,它要求将ODE中的每一步解细分并解决。
与前面提到的三种解法不同,它首先要求将ODE分解成若干离散区间,然后计算每一段区间中的点。
这种解法可以更准确地进行处理,但时间消耗较多,因此比较少被使用。
优化方案在需要解决常微分方程时,为了得到最佳的结果,有必要考虑一些优化措施。
•首先,应考虑将一个复杂的ODE拆分成一些更易解决的问题。
这样做的结果是,预见到解决此ODR的总耗时将会降低。
•其次,为了加快计算速度,可以考虑使用预解算法。
常微分方程的数值解法

常微分方程的数值解法1. 引言常微分方程是自变量只有一个的微分方程,广泛应用于自然科学、工程技术和社会科学等领域。
由于常微分方程的解析解不易得到或难以求得,数值解法成为解决常微分方程问题的重要手段之一。
本文将介绍几种常用的常微分方程的数值解法。
2. 欧拉方法欧拉方法是最简单的一种数值解法,其具体步骤如下:- 将自变量的区间等分为n个子区间;- 在每个子区间上假设解函数为线性函数,即通过给定的初始条件在每个子区间上构造切线;- 使用切线的斜率(即导数)逼近每个子区间上的解函数,并将其作为下一个子区间的初始条件;- 重复上述过程直至达到所需的精度。
3. 改进的欧拉方法改进的欧拉方法是对欧拉方法的一种改进,主要思想是利用两个切线的斜率的平均值来逼近每个子区间上的解函数。
具体步骤如下: - 将自变量的区间等分为n个子区间;- 在每个子区间上构造两个切线,分别通过给定的初始条件和通过欧拉方法得到的下一个初始条件;- 取两个切线的斜率的平均值,将其作为该子区间上解函数的斜率,并计算下一个子区间的初始条件;- 重复上述过程直至达到所需的精度。
4. 二阶龙格-库塔方法二阶龙格-库塔方法是一种更为精确的数值解法,其基本思想是通过近似计算解函数在每个子区间上的平均斜率。
具体步骤如下: - 将自变量的区间等分为n个子区间;- 在每个子区间上计算解函数的斜率,并以该斜率的平均值近似表示该子区间上解函数的斜率;- 利用该斜率近似值计算下一个子区间的初始条件,并进一步逼近解函数;- 重复上述过程直至达到所需的精度。
5. 龙格-库塔法(四阶)龙格-库塔法是目前常用的数值解法之一,其精度较高。
四阶龙格-库塔法是其中较为常用的一种,其具体步骤如下:- 将自变量的区间等分为n个子区间;- 在每个子区间上进行多次迭代计算,得到该子区间上解函数的近似值;- 利用近似值计算每个子区间上的斜率,并以其加权平均值逼近解函数的斜率;- 计算下一个子区间的初始条件,并进一步逼近解函数;- 重复上述过程直至达到所需的精度。
常微分方程数值解

常微分方程数值解常微分方程数值解是数学中的一门重要学科,主要研究如何求解常微分方程,在科学计算中有着重要的应用。
常微分方程模型是自然界中广泛存在的现象描述方法,有着广泛的应用领域。
比如,在物理学中,运动中的物体的位置、速度和加速度随时间的关系就可以通过微分方程描述;在经济学中,经济变化随时间的变化也可以用微分方程来描述。
而常微分方程数值解的求解方法则提供了一种快速、高效的计算手段。
一、常微分方程数值解的基本概念常微分方程就是一个描述自变量(通常是时间)与其导数之间关系的方程。
其一般形式如下:$\frac{dy}{dt} = f(y,t)$其中 $f(y,t)$ 是一个已知的函数。
常微分方程数值解就是对于一个常微分方程,对其进行数字计算求解的方法。
常微分方程数值解常使用数值积分的方法来求解。
由于常微分方程很少有解析解,因此数值解的求解方法显得尤为重要。
二、常微分方程数值解的求解方法常微分方程数值解的求解方法很多,以下介绍其中两种方法。
1.欧拉法欧拉法是最简单的一种数值算法,其思想是通过将一个微分方程转化为一个数值积分方程来求解。
其数值积分方程为:$y_{i+1}=y_i+hf(y_i,t_i)$其中 $h$ 为步长,可以理解为每次计算的间隔。
欧拉法的主要缺点是其精度比较低,收敛速度比较慢。
因此,当需要高精度的数值解时就需要使用其他的算法。
2.级数展开方法级数展开法是通过将一个待求解的微分方程进行Taylor级数展开来求解。
通过对Taylor级数展开的前若干项进行求和,可以得到微分方程与其解的近似解。
由于级数展开法的收敛速度很快,因此可以得到相对较高精度的数值解。
但是,当级数过多时,会出现截断误差。
因此,在实际应用中需要根据所需精度和计算资源的限制来选择适当的级数。
三、常微分方程数值解的应用常微分方程数值解在现代科学技术中有着广泛的应用。
以下介绍其中两个应用领域。
1.物理建模常微分方程的物理建模是常见的应用领域。
常微分方程的数值解法

主要内容
§1、引言 §2、初值问题的数值解法--单步法 §3、龙格-库塔方法 §4、收敛性与稳定性 §5、初值问题的数值解法―多步法 §6、方程组和刚性方程 §7、习题和总结
§1、 引 言 主要内容 ➢研究的问题 ➢数值解法的意义
1.什么是微分方程 ? 现实世界中大多数事物
使得对任意的x [a,b]及y1, y2都成立
则称 f (x,y) 对y 满足李普希兹条件,L 称为 Lipschitz常数.
就可保证方程解的存在唯一性
若 f (x,y) 在区域 G连续,关于y
满足李普希兹 条件
一阶常微分方程的初值问题的解存在且唯一. 我们以下的讨论,都在满足上述条件下进行.
一阶常微分方程组常表述为:
y(x0
)
y0
(1.2)
种 数 值 解
法
其中f (x,y)是已知函数,(1.2)是定解条件也称为 初值条件。
常微分方程的理论指出:
当 f (x,y) 定义在区域 G=(a≤x≤b,|y|<∞)
若存在正的常数 L 使:
(Lipschitz)条件
| f (x, y1) f (x, y2) | L | y1 y2 | (1.3)
节点 xi a ihi,一般取hi h( (b a) / n)即等距
要计算出解函数 y(x) 在一系列节点
a x0 x1 xn b
处的近似值 yi y(xi )
y f (x, y)
y
(
x0
)
y0
a xb
(1.1) (1.2)
对微分方程(1.1)两端从 xn到xn1 进行积分
内部联系非常复杂
其状态随着 时间、地点、条件 的不同而不同
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
第四章常微分方程数值解[课时安排]6学时[教学课型]理论课[教学目的和要求]了解常微分方程初值问题数值解法的一些基本概念,如单步法和多步法,显式和隐式,方法的阶数,整体截断误差和局部截断误差的区别和关系等;掌握一阶常微分方程初值问题的一些常用的数值计算方法,例如欧拉(Euler)方法、改进的欧拉方法、龙贝-库塔(Runge-Kutta)方法、阿达姆斯(Adams)方法等,要注意各方法的特点及有关的理论分析;掌握构造常微分方程数值解的数值积分的构造方法和泰勒展开的构造方法的基本思想,并能具体应用它们导出一些常用的数值计算公式及评估截断误差;熟练掌握龙格-库塔(R-K)方法的基本思想,公式的推导,R-K公式中系数的确定,特别是能应用“标准四阶R-K公式”解题;掌握数值方法的收敛性和稳定性的概念,并能确定给定方法的绝对稳定性区域。
[教学重点与难点]重点:欧拉方法,改进的欧拉方法,龙贝-库塔方法。
难点:R—K方法,预估-校正公式。
[教学内容与过程]4.1 引言本章讨论常微分方程初值问题(4.1.1)的数值解法,这也是科学与工程计算经常遇到的问题,由于只有很特殊的方程能用解析方法求解,而用计算机求解常微分方程的初值问题都要采用数值方法.通常我们假定(4.1.1)中f(x,y)对y满足Lipschitz条件,即存在常数L>0,使对,有(4.1.2) 则初值问题(4.1.1)的解存在唯一.假定(4.1.1)的精确解为,求它的数值解就是要在区间上的一组离散点上求的近似.通常取,h称为步长,求(4.1.1)的数值解是按节点的顺序逐步推进求得.首先,要对方程做离散逼近,求出数值解的公式,再研究公式的局部截断误差,计算稳定性以及数值解的收敛性与整体误差等问题.4.2 简单的单步法及基本概念4.2.1 Euler法、后退Euler法与梯形法求初值问题(4.1.1)的一种最简单方法是将节点的导数用差商代替,于是(4.1.1)的方程可近似写成(4.2.1)从出发,由(4.2.1)求得再将代入(4.2.1)右端,得到的近似,一般写成(4.2.2) 称为解初值问题的Euler法.Euler法的几何意义如图4-1所示.初值问题(4.1.1)的解曲线y=y(x)过点,从出发,以为斜率作一段直线,与直线交点于,显然有,再从出发,以为斜率作直线推进到上一点,其余类推,这样得到解曲线的一条近似曲线,它就是折线.Euler法也可利用的Taylor展开式得到,由(4.2.3) 略去余项,以,就得到近似计算公式(4.2.2).另外,还可对(4.1.1)的方程两端由到积分得(4.2.4)若右端积分用左矩形公式,用,,则得(4.2.2).如果在(4.2.4)的积分中用右矩形公式,则得(4.2.5)称为后退(隐式)Euler法.若在(4.2.4)的积分中用梯形公式,则得(4.2.6)称为梯形方法.上述三个公式(4.2.2),(4.2.5)及(4.2.6)都是由计算,这种只用前一步即可算出的公式称为单步法,其中(4.2.2)可由逐次求出的值,称为显式方法,而(4.2.5)及(4.2.6)右端含有当f对y非线性时它不能直接求出,此时应把它看作一个方程,求解,这类方法称为稳式方法.此时可将(4.2.5)或(4.2.6)写成不动点形式的方程这里对式(4.2.5)有,对(4.2.6)则,g与无关,可构造迭代法(4.2.7)由于对y满足条件(4.1.2),故有当或,迭代法(4.2.4)收敛到,因此只要步长h足够小,就可保证迭代(4.2.4)收敛.对后退Euler法(4.2.5),当时迭代收敛,对梯形法(4.2.6),当时迭代序列收敛.例4.1用Euler法、隐式Euler法、梯形法解取h=0.1,计算到x=0.5,并与精确解比较.解本题可直接用给出公式计算.由于,Euler法的计算公式为n=0时,.其余n=1,2,3,4的计算结果见表4-1.对隐式Euler法,计算公式为解出当n=0时,.其余n=1,2,3,4的计算结果见表4-1.表4-1 例4.1的三种方法及精确解的计算结果对梯形法,计算公式为解得当n=0时,.其余n=1,2,3,4的计算结果见表4-1.本题的精确解为,表4-1列出三种方法及精确解的计算结果.4.2.2 单步法的局部截断误差解初值问题(4.1.1)的单步法可表示为(4.2.8)其中与有关,称为增量函数,当含有时,是隐式单步法,如(4.2.5)及(4.2.6)均为隐式单步法,而当不含时,则为显式单步法,它表示为(4.2.9)如Euler法(4.2.2),.为讨论方便,我们只对显式单步法(4.2.9)给出局部截断误差概念.定义2.1设y(x)是初值问题(4.1.1)的精确解,记(4.2.10)称为显式单步法(4.2.9)在的局部截断误差.之所以称为局部截断误差,可理解为用公式(4.2.9)计算时,前面各步都没有误差,即,只考虑由计算到这一步的误差,此时由(4.2.10)有局部截断误差(4.2.10)实际上是将精确解代入(4.2.9)产生的公式误差,利用Taylor展开式可得到.例如对Euler法(4.2.2)有,故它表明Euler法(4.2.2)的局部截断误差为,称为局部截断误差主项.定义2.2 设是初值问题(4.1.1)的精确解,若显式单步法(4.2.9)的局部截断误差,是展开式的最大整数,称为单步法(4.2.9)的阶,含的项称为局部截断误差主项.根据定义,Euler法(4.2.2)中的=1故此方法为一阶方法.对隐式单步法(4.2.8)也可类似求其局部截断误差和阶,如对后退Euler法(4.2.5)有局部截断误差故此方法的局部截断误差主项为,也是一阶方法.对梯形法(4.2.6)同样有它的局部误差主项为,方法是二阶的.4.2.3 改进Euler法上述三种简单的单步法中,梯形法(4.2.6)为二阶方法,且局部截断误差最小,但方法是隐式的,计算要用迭代法.为避免迭代,可先用Euler法计算出的近似,将(4.2.6)改为(4.2.11)称为改进Euler法,它实际上是显式方法.即(4.2.12)右端已不含.可以证明,=2,故方法仍为二阶的,与梯形法一样,但用(4.2.11)计算不用迭代.例4.2用改进Euler法求例4.1的初值问题并与Euler法和梯形法比较误差的大小.解将改进Euler法用于例4.1的计算公式当n=0时,.其余结果见表4-2.表4-2 改进Euler法及三种方法的误差比较从表4-2中看到改进Euler法的误差数量级与梯形法大致相同,而比Euler法小得多,它优于Euler法.讲解:求初值问题(4.1.1)的数值解就是在假定初值问题解存在唯一的前提下在给定区间上的一组离散点上求解析解的一组近似为此先要建立求数值解的计算公式,通常称为差分公式,简单的单步法就是由计算下一步,构造差分公式有三种方法,一是用均差(即差商)近似,二是用等价的积分方程(4.2.4)用数值积分方法,三是用函数的Taylor展开,其中Taylor展开最有普遍性,可以得到任何数值解的计算公式及其局部截断误差。
计算公式是微分方程的一种近似,局部截断误差的概念就是刻划这种逼迫的好坏。
当为微分方程的解,即,而用,定义局部截断误差,它表示用精确解代入计算公式(4.2.9)产生的公式误差为越大表明公式逼近微分方程的精度越高,因此就定义为公式的阶,通常的公式才能用于计算初值问题(4.1.1)的数值解。
利用Taylor展开时,只要将的表达式在处展开成Taylor公式就可得到不同公式的局部截断误差。
如4.2.2所给出的Euler法。
后退Euler法和梯形法,它们只需用一元函数的Taylor展开,与后面4.5节的多步法完全一致,而通常单步法(4.2.9)的一般情况则需要用二元函数的Taylor展开,才能得到公式的具体形式和局部截断误差。
例如对改进Euler法,其局部截断误差由(4.2.12)可得要求出它的结果就要用到二元函数的Taylor展开,将在4.3节再作介绍。
4.3 Runge-Kutta方法4.3.1 显式 Runge-Kutta法的一般形式上节已给出与初值问题(4.1.1)等价的积分形式(4.3.1)只要对右端积分用不同的数值求积公式近似就可得到不同的求解初值问题(4.1.1)的数值方法,若用显式单步法(4.3.2)当,即数值求积用左矩形公式,它就是Euler法(4.2.2),方法只有一阶,若取(4.3.3)就是改进Euler法,这时数值求积公式是梯形公式的一种近似,计算时要用二个右端函数f的值,但方法是二阶的.若要得到更高阶的公式,则求积分时必须用更多的f值,根据数值积分公式,可将(4.3.1)右端积分表示为注意,右端f中还不能直接得到,需要像改进Euler法(4.2.11)一样,用前面已算得的f值表示为(4.3.3),一般情况可将(4.3.2)的 表示为(4.3.4)其中这里均为待定常数,公式(4.3.2),(4.3.4)称为r级的显式Runge-Kutta法,简称R-K方法.它每步计算r个f值(即),而ki 由前面(i-1)个已算出的表示,故公式是显式的.例如当r=2时,公式可表示为(4.3.5)其中.改进Euler法(4.2.11)就是一个二级显式R-K方法.参数取不同的值,可得到不同公式.4.3.2 二、三级显式R-K方法对r=2的显式R-K方法(4.3.5),要求选择参数,使公式的阶p尽量高,由局部截断误差定义(4.3.6)令,对(4.3.6)式在处按Taylor公式展开,由于将上述结果代入(4.3.6)得要使公式(4.3.5)具有的阶p=2,即,必须(4.3.4)即由此三式求的解不唯一.因r=2,故,于是有解(4.3.8)它表明使(4.3.5)具有二阶的方法很多,只要都可得到二阶R-K方法.若取,则,则得改进Euler法(4.2.11),若取,则得,此时(4.3.5)为(4.3.9)其中称为中点公式.后退Euler法(4.2.11)及中点公式(4.3.9)是两个常用的二级R-K方法,注意二级R-K方法只能达到二阶,而不可能达到三阶.因为r=2只有4个参数,要达到p=3则在(4.3.6)的展开式中要增加3项,即增加三个方程.加上(4.3.4)的三个方程求4个待定参数是无解的.当然r=2,p=2的R-K方法(4.3.5)当取其他数时,也可得到其他公式,但系数较复杂,一般不再给出.对r=3的情形,要计算三个k值,即其中将按二元函数在处按Taylor公式展开,然后代入局部截断误差表达式,可得可得三阶方法,其系数应满足方程(4.3.10)这是8个未知数6个方程的方程组,解也是不唯一的,通常.一种常见的三级三阶R-K方法是下面的Kutta三阶方法:(4.3.11)4.3.3 四阶R-K方法及步长的自动选择利用二元函数Taylor展开式可以确定(4.3.4)中r=4,p=4的R-K方法,经典的四阶R-K 方法是:(4.3.12)它的局部截断误差,故p=4,这是最常用的四阶R-K方法,数学库中都有用此方法求解初值问题的软件.这种方法的优点是精度较高,缺点是每步要算4个右端函数值,计算量较大.例4.3用经典四阶R-K方法解例4.1的初值问题,仍取h=0.1,计算到,并与改进Euler法、梯形法在处比较其误差大小.解用四阶R-K方法公式(4.3.12),此处,于是当n=0时于是,按公式(4.3.12)可算出此方法误差:改进Euler法误差:梯形法误差:可见四阶R-K方法的精度比二阶方法高得多.用四阶R-K方法求解初值问题(4.1.1)精度较高,但要从理论上给出误差的估计式则比较困难.那么应如何判断计算结果的精度以及如何选择合适的步长h?通常是通过不同步长在计算机上的计算结果近似估计.设在处的值,当时,的近似为,于是由四阶R-K方法有若以为步长,计算两步到,则有于是得即或(4.3.13)它给出了误差的近似估计.如果(ε为给定精度),则认为以为步长的计算结果满足精度要求,若,则还可放大步长.因此(4.3.13)提供了自动选择步长的方法.讲解:求初值问题(4.1.1)的单步法主要是指Runge-Kutta法,本节主要讨论显式R-K方法,建立具体的计算公式使用的是Taylor展开,形如(4.3.4)的显式R-K方法,当r=1时就是Euler法,因此只要讨论的计算公式,在r确定后如何推导公式都是一样的,只是r越大计算越复杂,为了掌握了解公式来源,只要以r=2为例推导计算公式即可。