大连理工大学 高等数值分析 常微分方程数值解法-2017

合集下载

常微分方程数值解法

常微分方程数值解法

用分段的折线逼近函数,此为 “折线法”而非“切线法”, 除第一个点是曲线上的切线,
其它都不是。
2、Euler方法的误差估计
1)局部截断误差。 在一步中产生的误差而非累积误差:
~
T x y y

n1
n1
n1
其中
~
y
是当
y
n

y(
x
)
n
(精确解!)时
n1
由Euler法求出的值,即y 无误差! n
y y x , y x y 则得: h f
f ,
n1
n2
n
n
n1
n1
同样与Euler法结合,形成迭代算法,对n 0,1,2,
y y x , y 0 hf
n1
n
n
n
y y x , y x y

k 1
推出总体误差与步长的关系。
由微分方程解的存在唯一性,自然假定 ( f x,y)
充分光滑,或满足 Lipschitz条件:
f
x,ny源自xn

f
x
,
n
y
n

L
yxn
y n
第 n 步 的 总 体 截 断 误 差 记 为
en y
xn
y n
则 对 n 1 步:
e x y x y y y T y y ~ ~
用yn1, yn代替y(xn1), y(xn ), 对右端积分采用 取左端点的矩形公式
则有
xn1 xn
f
(x,
y)dx

hf
(xn ,
yn )

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

数值分析第九章常微分方程数值解法
高斯-赛德尔迭代法
松弛法
通过迭代更新函数值并逐步放松约束 条件来逼近解,适用于刚性和非刚性 问题。
利用线性组合迭代函数值来逼近解, 具有更高的收敛速度和稳定性。
03
数值解法的稳定性分析
数值解法的稳定性定义
数值解法的稳定性是指当微分方程的初值有微小的扰动时, 其数值解的近似值的变化情况。如果数值解在微小扰动下变 化较小,则称该数值方法是稳定的。
更高的精度和稳定性。
数值逼近法
泰勒级数法
将微分方程的解展开为泰勒级数,通过截断级数来逼 近解。
多项式逼近法
利用多项式来逼近微分方程的解,通过选取合适的基 函数和系数来提高逼近精度。
样条插值法
利用样条函数来逼近微分方程的解,具有更好的光滑 性和连续性。
迭代法
雅可比迭代法
通过迭代更新函数值来逼近微分方程 的解,具有简单易行的优点。
初值和边界条件的处理
根据实际问题,合理设定初值和边界 条件,以获得更准确的数值解。
收敛性和误差分析
对数值解进行收敛性和误差分析,评 估解的精度和稳定性。
数值解法的应用案例分析
人口增长模型
通过数值解法求解人口增长模型,预测未来人口数量,为政策制 定提供依据。
化学反应动力学
利用数值解法研究化学反应的动力学过程,模拟反应过程和结果。
数值分析第九章常微分方 程数值解法
• 引言 • 常微分方程数值解法的基本思想 • 数值解法的稳定性分析 • 数值解法的收敛性和误差分析 • 数值解法的实现和应用案例
01
引言
常微分方程的应用背景
自然科学
描述物理、化学、生物等自然 现象的变化规律。
工程领域
控制系统设计、航天器轨道计 算等。

数值分析常微分方程数值解

数值分析常微分方程数值解

许多实际问题的数学模型是微分方程或微分方程的定解问题。

如物体运动、电路振荡、化学反映及生物群体的变化等。

常微分方程可分为线性、非线性、高阶方程与方程组等类;线性方程包含于非线性类中,高阶方程可化为一阶方程组。

若方程组中的所有未知量视作一个向量,则方程组可写成向量形式的单个方程。

因此研究一阶微分方程的初值问题⎪⎩⎪⎨⎧=≤≤=0)(),(y a y bx a y x f dxdy, (9-1) 的数值解法具有典型性。

常微分方程的解能用初等函数、特殊函数或它们的级数与积分表达的很少。

用解析方法只能求出线性常系数等特殊类型的方程的解。

对非线性方程来说,解析方法一般是无能为力的,即使某些解具有解析表达式,这个表达式也可能非常复杂而不便计算。

因此研究微分方程的数值解法是非常必要的。

只有保证问题(9-1)的解存在唯一的前提下,研究其数值解法或者说寻求其数值解才有意义。

由常微分方程的理论知,如果(9-1)中的),(y x f 满足条件(1)),(y x f 在区域} ),({+∞<<∞-≤≤=y b x a y x D ,上连续; (2)),(y x f 在上关于满足Lipschitz 条件,即存在常数,使得y y L y x f y x f -≤-),(),(则初值问题(9-1)在区间],[b a 上存在惟一的连续解)(x y y =。

在下面的讨论中,我们总假定方程满足以上两个条件。

所谓数值解法,就是求问题(9-1)的解)(x y y =在若干点b x x x x a N =<<<<= 210处的近似值),,2,1(N n y n =的方法。

),,2,1(N n y n =称为问题(9-1)的数值解,n n x x h -=+1称为由到1+n x 的步长。

今后如无特别说明,我们总假定步长为常量。

建立数值解法,首先要将微分方程离散化,一般采用以下几种方法: (1) 用差商近似导数在问题(9-1)中,若用向前差商hx y x y n n )()(1-+代替)(n x y ',则得)1,,1,0( ))(,()()(1-=≈-+N n x y x f hx y x y n n n n n)(n x y 用其近似值代替,所得结果作为)(1+n x y 的近似值,记为1+n y ,则有 1(,) (0,1,,1)n n n n y y hf x y n N +=+=-这样,问题(9-1)的近似解可通过求解下述问题100(,) (0,1,,1)()n n n n y y hf x y n N y y x +=+=-⎧⎨=⎩(9-2)得到,按式(9-2)由初值经过步迭代,可逐次算出N y y y ,,21。

数值分析常微分方程数值解法

数值分析常微分方程数值解法
7
第8页/共105页
➢ 数值积分方法(Euler公式)
设将方程 y=f (x, y)的两端从 xn 到xn+1 求积分, 得
y( xn1) y( xn )
xn1 f ( x, y( x))dx :
xn
xn1 F ( x)dx
xn
用不同的数值积分方法近似上式右端积分, 可以得到计算 y(xn+1)的不同的差分格 式.
h2 2
y''( )
Rn1
:
y( xn1)
yn1
h2 2
y''( )
h2 2
y''( xn ) O(h3 ).
局部截断误差主项
19
第20页/共105页
➢ 向后Euler法的局部截断误差
向后Euler法的计算公式
yn1 yn hf ( xn1, yn1 ), n 0, 1, 2,
定义其局部截断误差为
y 计算 的n递1 推公式,此类计算格式统称为差分格式.
3
第4页/共105页
数值求解一阶常微分方程初值问题
y' f ( x, y), a x b,
y(a)
y0
难点: 如何离散 y ?
➢ 常见离散方法
差商近似导数 数值积分方法 Taylor展开方法
4
第5页/共105页
➢ 差商近似导数(Euler公式)
(0 x 1)
y(0) 1.
解 计算公式为
yn1
yn
hfn
yn
h( yn
2xn ), yn
y0 1.0
n 0, 1, 2,
取步长h=0.1, 计算结果见下表
13

(完整版)大连理工大学高等数值分析偏微分方程数值解(双曲方程书稿)

(完整版)大连理工大学高等数值分析偏微分方程数值解(双曲方程书稿)

双曲型方程的有限差分法线性双曲型方程定解问题: (a )一阶线性双曲型方程()0=∂∂+∂∂xux a t u (b )一阶常系数线性双曲型方程组0=∂∂+∂∂xt uA u 其中A ,s 阶常数方程方阵,u 为未知向量函数。

(c )二阶线性双曲型方程(波动方程)()022=⎪⎭⎫⎝⎛∂∂∂∂-∂∂x u x a x t u()x a 为非负函数(d )二维,三维空间变量的波动方程0222222=⎪⎪⎭⎫⎝⎛∂∂+∂∂-∂∂y u x u t u 022222222=⎪⎪⎭⎫ ⎝⎛∂∂+∂∂+∂∂-∂∂z u y u xu t u §1 波动方程的差分逼近 1.1 波动方程及其特征线性双曲型偏微方程的最简单模型是一维波动方程:(1.1) 22222xu a t u ∂∂=∂∂ 其中0>a 是常数。

(1.1)可表示为:022222=∂∂-∂∂xu a t u ,进一步有0=⎪⎭⎫ ⎝⎛∂∂+∂∂⋅⎪⎭⎫ ⎝⎛∂∂-∂∂u x a t x a t由于xat ∂∂±∂∂当a dt dx ±=时为()t x u ,的全导数(=dtdu dt dx x u t u ⋅∂∂+∂∂x ua t u ∂∂±∂∂=),故由此定出两个方向(1.3) adx dt 1±=解常微分方程(1.3)得到两族直线 (1.4) 1C t a x =⋅+ 和 2C t a x =⋅- 称其为特征。

特征在研究波动方程的各种定解问题时,起着非常重要的作用。

比如,我们可通过特征给出(1.1)的通解。

(行波法、特征线法) 将(1.4)视为),(t x 与),(21C C 之间的变量替换。

由复合函数的微分法则212211C uC u x C C u x C C u x u ∂∂+∂∂=∂∂⋅∂∂+∂∂⋅∂∂=∂∂ x C C u C u C x C C u C u C x u ∂∂⎪⎪⎭⎫ ⎝⎛∂∂+∂∂∂∂+∂∂⎪⎪⎭⎫ ⎝⎛∂∂+∂∂∂∂=∂∂2212121122 222122212212C u C C u C C u C u ∂∂+∂∂∂+∂∂∂+∂∂= 2222122122C uC C u C u ∂∂+∂∂∂+∂∂= 同理可得a t t a t C -=∂∂-=∂∂1,a tC=∂∂2 ⎪⎪⎭⎫⎝⎛∂∂-∂∂=∂∂⋅∂∂+∂∂⋅∂∂=∂∂212211C u C u a t C C u t C C u t utC C u C u a C u t C C u C u a C t u ∂∂⎪⎪⎭⎫ ⎝⎛⎪⎪⎭⎫ ⎝⎛∂∂-∂∂⋅∂∂+∂∂⎪⎪⎭⎫ ⎝⎛⎪⎪⎭⎫ ⎝⎛∂∂-∂∂⋅∂∂=∂∂2122112122 ⎥⎦⎤⎢⎣⎡∂∂∂-∂∂+⎥⎦⎤⎢⎣⎡∂∂-∂∂∂-=21222222221222C C u C u a C u C C u a ⎥⎦⎤⎢⎣⎡∂∂+∂∂∂-∂∂=22221221222C u C C u C u a 将22x u ∂∂和22tu∂∂代入(1.1)可得: ⎥⎦⎤⎢⎣⎡∂∂+∂∂∂-∂∂22221221222C u C C u C u a ⎥⎦⎤⎢⎣⎡∂∂+∂∂∂+∂∂=22221221222C u C C u C u a 即有0212=∂∂∂C C u求其对2C 的积分得:()11C f C u=∂∂ 其中()1C f 是1C 的任意可微函数。

第九节常微分方程的数值解法

第九节常微分方程的数值解法
N
卷积
• 在求拉氏逆变换的过程中,卷积往往 有着重要的应用价值。 • 定义 f * f f ( ) f (t )d
t 1 2 0 1 2
为函数f1(t)与f2(t)的卷积。 • 注意:当t<0时,f1(t)=f2(t) 0
• 交换律
f1*f2=f2*f1
• 例1
求t*sint f1*(f2+f3)=f1*f2+f1*f3
• (每一步计算f(x,y)四次,截断误差为O(h5))
• 例 分别用改进的欧拉格式和四阶龙格— 库塔格式解初值问题(取步长h=0.2):
2x y y y y (0) 1
(0 x 1)

表7—4
准确解
• 节点 改进欧拉法 四阶龙格—库塔法

• • •

1 0.2 1.186667 1.183216 0.4 1.348312 1.341641 0.6 1.493704 1.483240 0.8 1.627861
第九节
常微分方程的数值解法
dy f ( x, y ) dx y ( x0 ) y0
• 一阶常微分方程的初值问题:
• 节点:x1<x2< <xn • 步长 h x x 为常数 i i 1
• • • • •

欧拉方法(折线法) yi+1=yi+hf(xi,yi)(i=0,1, , n 优点:计算简单。 缺点:精度不高。 二 改进的欧拉方法
• • • • •


第三步 如果 1 ,则将步长h折 半, xi 1/ 2 从xi出发以区间[xi,xi+1]的中点(记为 ) ( h / 2) (h) 1/ 2 yi 1/ 2 yi 1/ 2 为目标,判别 y ( xi 1/ 2 ) 1 / 2 如果 yi(h1//22) ,则得 xi 1 / 2 的满足精度要 求 的近似值 ,然后从 出发,以 xi+1 为 目标,重复上述步骤,否则继续下一步

常微分方程数值解法5262115页PPT文档

常微分方程数值解法5262115页PPT文档
x 1 ( t ) 表示时刻 t 食饵的密度,x 2 ( t ) 表示捕食者的密度;
r 表示食饵独立生存时的增长率;
d 表示捕食者独立生存时的死亡率;
a 表示捕食者的存在对食饵增长的影响系数,反映捕
食者对食饵的捕获能力;
b 表示食饵的存在对捕食者增长的促进系数,反映食
饵对捕食者的喂养能力
150 100
令 y 1 y ,y 2 y ',y 3 y '', ,y n y ( n 1 )
可以将以上高阶微分方程化为如下一阶常微分方程组
y1 ' y2 y2 ' y3 yn ' an(x)y1
a1(x)yn f (x)
例:P120,1(a),Bessel方程
常微分方程的数值解
一般地,凡表示未知函数,未知函数的导 数与自变量之间的关系的方程叫做微分方 程.未知函数是一元函数的,叫常微分方 程;未知函数是多元函数的,叫做偏微分方 程.

y ' x y'x2y2 y''y'xy
Matlab实现 [t,x]=ode45(f,ts,x0,options,p1,p2,......)
50 0 0
30 20 10
0 0
10
20
50
30
20
10

0
30
0
10
8
6
4
2
100
0
50
100
150
50
100
高阶常微分方程的解法
高阶常微分方程
y ( n ) a 1 ( x ) y ( n 1 ) a ( n 1 ) ( x ) y ' a n ( x ) y f( x )

数值分析课件第9章常微分方程初值问题数值解法

数值分析课件第9章常微分方程初值问题数值解法
隐式方程通常用迭代法求解,而迭代过程的实质 是逐步显式化.
上页 下页
设用欧拉公式
y(0) n1
yn
hf
( xn ,
yn )
给出迭代初值 yn(0)1,用它代入(2.5)式的右端,使之转
化为显式,直接计算得
y(1) n1
yn
hf
( xn1 ,
y(0) n1
),
然后再用 yn(代1)1 入(2.5)式,又有
在yn=y(xn)的前提下,f(xn,yn )=f(xn,y(xn))=y(xn).于是
可得欧拉法(2.1)的公式误差为
y(xn1)
yn1
h2 2
y(n )
h2 2
y(xn ),
(2.3)
称为此方法的局部截断误差.
上页 下页
如果对方程(1.1)从xn到xn+1积分,得
y( xn1) y( xn )
基于上述几何解释,我们从初始点P0(x0, y0)出发, 先依方向场在该点的方向推进到x=x1上一点P1,然后 再从P1点依方向场在该点的方向推进到 x=x2 上一点 P2 , 循环前进做出一条折线P0 P1 P2.
上页 下页
一般地,设已做出该折线的顶点Pn,过Pn(xn, yn)依
方向场的方向再推进到Pn+1(xn+1, yn+1),显然两个顶
(2.13)
校正 yn1 yn 2 [ f ( xn , yn ) f ( xn1, yn1 )] .
或表为下列平均化形式
yp yc
y(2) n1
yn
hf
( xn1 ,
y(1) n1
).
如此反复进行,得
y(k1) n1

数值分析常微分方程的数值解法

数值分析常微分方程的数值解法

《计算机数学基础》数值部分第五单元辅导14 常微分方程的数值解法一、重点内容 1. 欧拉公式:),...,,,(),()(1-210=⎩⎨⎧+=+=≈01+1+n k kh x x y x hf y y x y kk k k k k局部截断误差是O (h 2)。

2. 改进欧拉公式:预报-校正公式:⎪⎩⎪⎨⎧++=+=++++)],(),([2),(1111k k k k k k k k k k y x f y x f hy y y x hf y y 校正值预报值即 ))],(,(),([211k k k k k k k k y x hf y x f y x f hy y +++=++ 或表成平均的形式:⎪⎪⎪⎩⎪⎪⎪⎨⎧+21=+=+=1+1+)(),(),(c p k p k k c k k k p y y y y x hf y y y x hf y y改进欧拉法的局部截断误差是O (h 3)3. 龙格-库塔法二阶龙格-库塔法的局部截断误差是O (h 3) 三阶龙格-库塔法的局部截断误差是O (h 4) 四阶龙格−库塔法公式: )22(643211κκκκ++++=+hy y k k其中 κ1=f (x k ,y k );κ2=f (x n +12h ,y k +21h κ1);κ3=f (x k +12h ,y n +21h κ2);κ4=f (x k +h ,y k +h κ3)四阶龙格-库塔法的局部截断误差是O (h 5)。

二、实例例1 用欧拉法解初值问题⎩⎨⎧1=060≤≤0--='2)().(y x xy y y ,取步长h =0.2。

计算过程保留4位小数。

解h =0.2, f (x )=-y -xy 2。

首先建立欧拉迭代格式),,)((.),(210=-420=--=+=21+k y x y y hx hy y y x hf y y k k k kk k k k k k k当k =0,x 1=0.2时,已知x 0=0,y 0=1,有y (0.2)≈y 1=0.2×1(4-0×1)=0.8000当k =1,x 2=0.4时,已知x 1=0.2, y 1=0.8,有 y (0.4)≈y 2=0.2×0.8×(4-0.2×0.8)=0.614 4 当k =2,x 3=0.6时,已知x 2=0.4,y 2=0.6144,有 y (0.6)≈y 3=0.2×0.6144×(4-0.4×0.4613)=0.8000例2 用欧拉预报-校正公式求解初值问题⎩⎨⎧1=10=++'2)(sin y x y y y ,取步长h =0.2,计算y (0.2),y (0.4)的近似值,计算过程保留5位小数。

常微分方程数值解法2线性多步法

常微分方程数值解法2线性多步法
对于线性多步法,其收敛性取决于微分方程的解的性质和方法的阶数。一般来说,高阶方法具有更好 的收敛性。
03
常见的线性多步法
欧拉方法
总结词
欧拉方法是常微分方程数值解法中最简单的一种方法,它基于线性近似,通过已知的函 数值来估计新的函数值。
详细描述
欧拉方法的基本思想是利用已知的函数值来估计下一个点的函数值。具体来说,假设我 们有一个函数 (y = f(x)),在已知 (x_0) 处的函数值 (y_0 = f(x_0)) 的情况下,欧拉方法 通过线性插值来估计 (x_1) 处的函数值 (y_1),即 (y_1 = y_0 + h f(x_0)),其中 (h) 是
05
线性多步法的优缺点
优点
稳定性好
线性多步法在处理常微分方程时具有较好的数值稳定性, 能够有效地抑制数值振荡,提高计算结果的精度。
01
易于实现
线性多步法的计算过程相对简单,易于 编程实现,适合于大规模数值计算。
02
03
精度可调
通过选择不同的步长和线性多步法公 式,可以灵活地调整计算结果的精度, 满足不同的数值模拟需求。
改进方法的收敛性
研究收敛性条件
深入研究线性多步法的收敛性条件,了解哪些情况下方法可能不收 敛,并寻找改进措施。
优化迭代算法
通过优化迭代算法,提高方法的收敛速度和精度,减少迭代次数, 提高计算效率。
引入预处理技术
利用预处理技术对线性系统进行预处理,改善系统的条件数,提高方 法的收敛性。
拓展应用领域
在工程问题中的应用
控制系统设计
在工程领域中,线性多步法可以用于控制系统设计,通过 建立控制系统的数学模型,设计控制算法和控制器,实现 系统的稳定性和性能优化。

常微分方程的数值解法全文

常微分方程的数值解法全文

第8章常微分方程的数值解法8.4单步法的收敛性与稳定性8.4.1相容性与收敛性上面所介绍的方法都是用离散化的方法,将微分方程初值问题化为差分方程初值问题求解的.这些转化是否合理?即当h →∞时,差分方程是否能无限逼近微分方程,差分方程的解n y 是否能无限逼近微分方程初值问题的准确解()n y x ,这就是相容性与收敛性问题.用单步法(8.3.14)求解初值问题(8.1.1),即用差分方程初值问题100(,,)()n n n n y y h x y h y x y ϕ+=+⎧⎨=⎩(8.4.1)的解作为问题(8.1.1)的近似解,如果近似是合理的,则应有()()(,(),)0 (0)y x h y x x y x h h hϕ+--→→(8.4.2)其中()y x 为问题(8.1.1)的精确解.因为0()()lim ()(,)h y x h y x y x f x y h→+-'==故由(8.4.2)得lim (,,)(,)h x y h f x y ϕ→=如果增量函数(,(),)x y x h ϕ关于h 连续,则有(,,0)(,)x y f x y ϕ=(8.4.3)定义8.3如果单步法的增量函数(,,)x y h ϕ满足条件(8.4.3),则称单步法(8.3.14)与初值问题(8.1.1)相容.通常称(8.4.3)为单步法的相容条件.满足相容条件(8.4.3)是可以用单步法求解初值问题(8.1.1)的必要条件.容易验证欧拉法和改进欧拉法均满足相容性条件.一般地,如果单步法有p 阶精度(1p ≥),则其局部截断误差为[]1()()(,(),)()p y x h y x h x y x h O h ϕ++-+=上式两端同除以h ,得()()(,,)()p y x h y x x y h O h hϕ+--=令0h →,如果(,(),)x y x h ϕ连续,则有()(,,0)0y x x y ϕ'-=所以1p ≥的单步法均与问题(8.1.1)相容.由此即得各阶龙格-库塔法与初值问题(8.1.1)相容.定义8.4一种数值方法称为是收敛的,如果对于任意初值0y 及任意固定的(,]x a b ∈,都有lim () ()n h y y x x a nh →==+其中()y x 为初值问题(8.1.1)的精确解.如果我们取消局部化假定,使用某单步法公式,从0x 出发,一步一步地推算到1n x +处的近似值1n y +.若不计各步的舍入误差,而每一步都有局部截断误差,这些局部截断误差的积累就是整体截断误差.定义8.5称111()n n n e y x y +++=-为某数值方法的整体截断误差.其中()y x 为初值问题(8.1.1)的精确解,1n y +为不计舍入误差时用某数值方法从0x 开始,逐步得到的在1n x +处的近似值(不考虑舍入误差的情况下,局部截断误差的积累).定理8.1设单步法(8.3.14)具有p 阶精度,其增量函数(,,)x y h ϕ关于y 满足利普希茨条件,问题(8.1.1)的初值是精确的,即00()y x y =,则单步法的整体截断误差为111()()p n n n e y x y O h +++=-=证明由已知,(,,)x y h ϕ关于y 满足利普希茨条件,故存在0L >,使得对任意的12,y y 及[,]x a b ∈,00h h <≤,都有1212(,,)(,,)x y h x y h L y y ϕϕ-≤-记1()(,(),)n n n n y y x h x y x h ϕ+=+,因为单步法具有p 阶精度,故存在0M >,使得1111()p n n n R y x y Mh ++++=-≤从而有111111111()()()(,(),)(,,)()(,(),)(,,)n n n n n n n p n n n n n n p n n n n n n e y x y y x y y y Mh y x h x y x h y h x y h Mh y x y h x y x h x y h ϕϕϕϕ+++++++++=-≤-+-≤++--≤+-+-1(1)p nMh hL e +≤++反复递推得11111101110(1)(1)1(1)(1)(1)(1)1(1)p p n n n p n n p n e Mh hL Mh hL e hL hL Mh hL e hL Mh hL e hL+++-+++++⎡⎤≤++++⎣⎦⎡⎤≤+++++++⎣⎦+-≤++因为00()y x y =,即00e =,又(1)n h b a +≤-,于是ln(1)1()(1)(1)b a b a hL n L b a h h hL hL e e --++-+≤+=≤所以()11()p L b a p n M e h e O h L -+⎡⎤≤-=⎣⎦推论设单步法具有p (1p ≥)阶精度,增量函数(,,)x y h ϕ在区域G :, , 0a x b y h h ≤≤-∞<<+∞≤≤上连续,且关于y 满足利普希茨条件,则单步法是收敛的.当(,)f x y 在区域:,D a x b y ≤≤-∞<<+∞上连续,且关于y 满足利普希茨条件时,改进欧拉法,各阶龙格-库塔法的增量函数(,,)x y h ϕ在区域G 上连续,且关于y 满足利普希茨条件,因而它们都是收敛的.关于单步法收敛的一般结果是:定理8.2设增量函数(,,)x y h ϕ在区域G 上连续,且关于y 满足利普希茨条件,则单步法收敛的充分必要条件是相容性条件(8.4.3).8.4.2稳定性稳定性与收敛性是两个不同的概念,收敛性是在假定每一步计算都准确的前提下,讨论当步长0h →时,方法的整体截断误差是否趋于零的问题.而稳定性则是讨论舍入误差的积累能否对计算结果有严重影响的问题.定义8.6若一种数值方法在节点值n y 上有一个大小为δ的扰动,于以后各节点()m y m n >上产生的偏差均不超过δ,则称该方法是稳定的.我们以欧拉法为例进行讨论.假设由于舍入误差,实际得到的不是n y 而是n n n y y δ=+,其中n δ是误差.由此再计算一步,得到1(,)n n n n y y hf x y +=+把它与不考虑舍入误差的欧拉公式相减,并记111n n n y y δ+++=-,就有[]1(,)(,)1(,)n n n n n n y n nh f x y f x y hf x δδηδ+⎡⎤=+-=+⎣⎦其中y f f y∂=∂.如果满足条件1(,)1y n hf x η+≤,(8.4.4)则从n y 到1n y +的计算,误差是不增的,可以认为计算是稳定的.如果条件(8.4.4)不满足,则每步误差将增大.当0y f >时,显然条件(8.4.4)不可能满足,我们认为问题本身具有先天的不稳定性.当0y f <时,为了满足稳定性要求(8.4.4),有时h 要很小.一般的,稳定性与方法有关,也与步长h 的大小有关,当然也与方程中的(,)f x y 有关.为简单起见,通常只考虑数值方法用于求解模型方程的稳定性,模型方程为y y λ'=(8.4.5)其中λ为复数.一般的方程可以通过局部线性化转化为模型方程,例如在(,)x y 的邻域(,)(,)(,)()(,)()x y y f x y f x y f x y x x f x y y y '==+-+-+略去高阶项,再作变量替换就得到u u λ'=的形式.对于模型方程(8.4.5),若Re 0λ>,类似以上分析,可以认为方程是不稳定的.所以我们只考虑Re 0λ<的情形,这时不同的数值方法可能是数值稳定的或者是数值不稳定的.当一个单步法用于试验方程y y λ'=,从n y 计算一步得到1()n n y E h y λ+=(8.4.6)其中()E h λ依赖于所选的方法.因为通过点(,)n n x y 试验方程的解曲线(它满足,()n n y y y x y λ'==)为[]exp ()n n y y x x λ=-,而一个p 阶单步法的局部截断误差在()n n y x y =时有1111()()p n n n T y x y O h ++++=-=,所以有1exp()()()p n n y h E h y O h λλ+-=(8.4.7)这样可以看出()E h λ是h e λ的一个近似值.由(8.4.6)可以看到,若n y 计算中有误差ε,则计算1n y +时将产生误差()E h λε,所以有下面定义.定义8.7如果(8.4.6)式中,()1E h λ<,则称单步法(8.3.14)是绝对稳定的.在复平面上复变量h λ满足()1E h λ<的区域,称为方法(8.3.14)的绝对稳定区域,它与实轴的交称为绝对稳定区间.在上述定义中,规定严格不等式成立,是为了和线性多步法的绝对稳定性定义一致.事实上,()1E h λ=时也可以认为误差不增长.(1)欧拉法的稳定性欧拉法用于模型方程(8.4.5),得1(1)n n y h y λ+=+,所以有()1E h h λλ=+.所以绝对稳定条件是11h λ+<,它的绝对稳定区域是h λ复平面上以(1,0)-为中心的单位圆,见图8.3.而λ为实数时,绝对稳定区间是(2,0)-.Im()h λRe()h λ2-1-O 图8.3欧拉法的绝对稳定区域(2)梯形公式的稳定性对模型方程,梯形公式的具体表达式为11()2n n n n h y y y y λλ++=++,即11212n nh y y h λλ++=-,所以梯形公式的绝对稳定区域为12112h h λλ+<-.化简得Re()0h λ<,因此梯形公式的绝对稳定区域为h λ平面的左半平面,见图8.4.特别地,当λ为负实数时,对任意的0h >,梯形公式都是稳定的.Im()h λRe()h λO 图8.4梯形公式的绝对稳定区域(3)龙格-库塔法的稳定性与前面的讨论相仿,将龙格-库塔法用于模型方程(8.4.5),可得二、三、四阶龙格-库塔法的绝对稳定区域分别为211()12h h λλ++<23111()()126h h h λλλ+++<2341111()()()12624h h h h λλλλ++++<当λ为实数时,二、三、四阶显式龙格-库塔法的绝对稳定区域分别为20h λ-<<、2.510h λ-<<、 2.780h λ-<<.例8.5设有初值问题21010101(0)0xy y x x y ⎧'=-≤≤⎪+⎨⎪=⎩用四阶经典龙格-库塔公式求解时,从绝对稳定性考虑,对步长h 有何限制?解对于所给的微分方程有2100,(010)1f x x y xλ∂==-<≤≤∂+在区间[0,10]上,有201010max ||max51t x x λ<<==+由于四阶经典龙格-库塔公式的绝对稳定区间为 2.7850h λ-<<,则步长h 应满足00.557h <<.。

第九章常微分方程初值问题数值解法数值分析

第九章常微分方程初值问题数值解法数值分析

yn1
由此只要 hL<1 上述迭代就可收敛到 y n+1
又如把积分区间改为[xn-1, xn+1],并对右端用中矩形公式代替
积分,可得到:
y(xn1)
y ( xn1 )
2hf
(xn ,
y(xn ))
h3 3
y( )
(xn , xn1)
推出: yn1 yn1 2hf (xn , yn ) n 1,2,
对任意的 x ∈[a, b] 和任意的y1, y2 成立,则称 f (x , y)对 y 满足利普希茨(Lipschitz)条件, L称利普希茨常数. 若f ( x, y) 在 区域G连续,关于y 满足利普希茨条件,则上述一阶常微分方 程的初值问题的解存在且唯一. 我们以下的讨论,都在满足 上述条件下进行.
决定其精度的是函数 f (xn , yn). 改进成上述增量函数, 就可 提高公式的精度.
前面所得的改进的欧拉法可表示为:
yn1
yn
h[ 2
f
(xn ,
yn )
f
( xn
h, hf
(xn ,
yn ))]
增量函数
(xn ,
yn , h)
1[f 2
(xn ,
yn )
f
( xn
h,
yn
hf
(xn ,
yn ))]
由欧拉法的公式:
yn+1= yn+h f (xn ,yn)
n=0,1,2,… 20
线法. 8
欧拉隐式公式和欧拉中点公式
对(*)式的右端用右矩形公式代替积分,可得到
y(xn1)
y(xn ) hf
(xn1, y(xn1))

数值分析中的常微分方程数值求解

数值分析中的常微分方程数值求解

数值分析中的常微分方程数值求解常微分方程是自然科学中一类最为普遍的数学模型,涉及到热力学、物理、化工等多个领域。

然而,解常微分方程并非易事。

尤其是当我们面对一些复杂、非线性、多维的方程组时,常微分方程数值求解成为了一个十分关键的问题。

因此,数值求解方法成为了常微分方程研究中的重要组成部分。

本文将介绍一些数值解常微分方程的常见方法和应用。

1. 一般线性方法一般线性方法(general linear methods)是经典的常微分方程数值解法之一。

它以一种特殊的形式给出步进公式:$$ y_{n+1}=\sum_{i=0}^{s-1}\alpha_i y_{n-i}+h\sum_{i=0}^{s-1}\beta_i f(t_{n-i},y_{n-i}) $$ 其中,$y_{n}$为第$n$步的项值,$f(t_n,y_n)$为时间$t_n$处函数$y(t)$的导数。

$\alpha_i$和$\beta_i$是常数,可以通过确定如下特征方程来选择:$$ \sum_{i=0}^{s-1}\alpha_i\ lambda^{i}=0,~(\lambda\in C) $$ 与此同时,也可以通过选择$\beta_i$来使方法达到一定的准确性和稳定性。

2. Runge-Kutta方法比一般线性方法更为流行的方法是Runge-Kutta方法。

通常附加一个或多个修正以获得更好的数值稳定性和误差控制。

第1阶Runge-Kutta方法仅使用导数$f(t_n,y_n)$估算下一个项的值:$$y_{n+1}=y_n+hf(t_n,y_n)$$ 许多高阶方法可以使用中间的“插值”来更准确地估计下一个步骤:$$y_{n+1}=y_n+h\sum_{i=1}^kb_ik_i$$$$k_i=f(t_n+c_ih,y_n+h\sum _{j=1}^{i-1}a_{ij}k_j)$$ $k_i$是第$i$台车的估计值,$a_{ij}$和$b_i$在经典Runge-Kutta方法和其他变体中具有不同的取值。

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

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

数值分析第九章常微分方程数值解法常微分方程数值解法是数值分析中非常重要的一部分内容。

常微分方程是描述自然现象中动态变化规律的数学模型,解常微分方程可以揭示系统的变化趋势和规律。

然而,大多数常微分方程是无法通过解析方法求出解的,因此需要借助计算机进行数值计算。

数值解常微分方程方法主要包括:Euler方法、改进的Euler方法、四阶Runge-Kutta方法和龙格-库塔方法。

Euler方法是最简单的方法之一,它采用的是一阶Taylor展开式。

将待求的函数值与函数的一阶导数值代入Taylor展开式中,可以得到函数值在下一个时间步长上的近似值。

Euler方法的优点是简单易于实现,但其精度不够高,容易积累误差。

改进的Euler方法是对Euler方法的改进,它通过使用中间点上的导数值来减小误差。

改进的Euler方法的精度相比Euler方法要高一些,但仍然不够高。

四阶Runge-Kutta方法是目前使用较为广泛的数值解常微分方程的方法之一、它通过计算不同时间点上的斜率来估计函数值,在多个时间点上计算斜率的平均值来提高精度。

四阶Runge-Kutta方法的精度比Euler方法和改进的Euler方法要高,但计算量也相对较大。

龙格-库塔方法是数值解常微分方程中最常用的方法之一、它是四阶Runge-Kutta方法的延伸,通过计算不同时间点上的斜率来估计函数值,然后利用这些估计值计算更准确的斜率,在不同步长上进行迭代计算,直到满足所需精度。

龙格-库塔方法的精度比四阶Runge-Kutta方法要高,但计算量也相对较大。

除了以上几种方法外,还有一些其他数值解常微分方程的方法,如Adams法、Gear法等。

这些方法在不同场景下有着不同的适用性和优劣势。

总结起来,数值解常微分方程方法是研究常微分方程数值计算中的重要内容。

不同的方法有着不同的精度和计算量,可以根据具体问题的特点选择合适的方法进行数值计算。

然而,需要注意的是,数值解只是在给定的步长下对函数的近似值,可能会引入误差。

数值分析常微分方程数值解法

数值分析常微分方程数值解法
10定义若某算法在计算过程中任一步产生的误差在以后的计算中都逐步衰减则称该算法是绝对稳定的absolutelystable常数可以一般分析时为简单起见只考虑试验方程test是复数equation时将某算法应用于上式幵假设只在初值产生误差
(Numerical Methods for Ordinary Differential Equations )
中心差商近似导数
y( x1 )
y( x2 ) y( x0 ) 2h
y( x2 ) y( x0 ) 2h f ( x1 , y( x1 ))
x0
x1
x2
yi1 yi1 2h f ( xi , yi ) i 1, ... , n 1
假设 yi1 y( xi1 ), yi y( xi ) ,则可以导出 Ri y( xi1 ) yi1 O(h3 ) 即中点公式具有 2 阶精度。
只要f (x, y)在a,b R1上连续, 且关于 y 满足
Lipschitz 条件,即存在与 x, y无关的常数 L 使 | f (x, y1) f (x, y2 ) | L | y1 y2 |
对任意定义在 a,b上的 y1 x, y2 x都成立,
则上述IVP存在唯一解。
求函数 y(x) 在一系列节点 a = x0< x1<…< xn= b 处的近似值
yi y(xi ) (i 1, ... , n)
的方法称为微分方程的数值解法。 y1,L , yn 称为微分方程的数值解。
称节点间距 hi xi1 xi (i 0, ... , n 1)为步长, 通常采用等距节点,即取 hi = h (常数)。
三、初值问题的离散化方法
离散化方法的基本特点是依照某一递推公式, 按节点从左至右的顺序依次求出y(xi )的近似 yi 值(i 1, ... , n),取 y0。

高阶常微分方程初值问题的数值解法

高阶常微分方程初值问题的数值解法

高阶常微分方程初值问题的数值解法一、概述在科学和工程领域中,常微分方程是一种常见的数学模型,用于描述各种自然现象和物理过程。

而对于高阶常微分方程初值问题的数值解法,一直是数学与计算领域的研究热点。

本文将介绍关于高阶常微分方程初值问题的数值解法的一些基本概念和方法。

二、高阶常微分方程初值问题的数值解法1. 常微分方程的基本概念我们来了解一下常微分方程的基本概念。

常微分方程是描述一个函数的导数与自身之间的关系的微分方程。

具体来说,对于一个函数y(x),其满足的方程为F(x, y, y', y'', ... ,y^(n)) = 0,其中y^(n)表示y的n阶导数。

2. 高阶常微分方程初值问题在前面的基本概念中,我们已经介绍了常微分方程的基本概念,接下来我们来看一下高阶常微分方程初值问题。

高阶常微分方程初值问题是指对于一个n阶常微分方程y^(n) = f(x, y, y', ... ,y^(n-1)),给定初始条件y(x0) = y0, y'(x0) = y1, ... , y^(n-1)(x0) = yn-1时,求解y(x)在给定区间上的解。

3. 数值解法对于高阶常微分方程初值问题,要求解其解析解通常是非常困难甚至不可能的。

我们通常会采用数值解法来逼近其解。

常见的数值解法包括欧拉方法、改进的欧拉方法、龙格-库塔方法等。

4. 欧拉方法欧拉方法是一种基本的数值解法,用于解决常微分方程初值问题。

其思想是通过离散化微分方程,从初始点开始,逐步迭代得到下一个点的近似解。

然而,欧拉方法的局限性在于其精度较低,容易积累误差。

5. 改进的欧拉方法针对欧拉方法精度较低的问题,人们提出了改进的欧拉方法,其主要思想是通过对步长h上进行两次迭代,从而得到更为精确的近似解。

改进的欧拉方法相对于普通欧拉方法,能够减小误差,提高精度。

6. 龙格-库塔方法龙格-库塔方法是一种高阶的数值解法,通过逐步计算区间上的近似解,其精度相对较高。

大连理工大学 高等数值分析 常微分方程数值解法-2017

大连理工大学 高等数值分析 常微分方程数值解法-2017

i.常微分方程初值问题数值解法i.1 常微分方程差分法考虑常微分方程初值问题:求函数()u t 满足(,), 0du f t u t T dt=<≤ (i.1a ) 0(0)u u = (i.1b)其中(,)f t u 是定义在区域G : 0t T ≤≤, u <∞上的连续函数,0u 和T 是给定的常数。

我们假设(,)f t u 对u 满足Lipschitz 条件,即存在常数L 使得121212(,)(,), [0,]; ,(,)f t u f t u L u u t T u u -≤-∀∈∈-∞∞ (i.2) 这一条件保证了(i.1)的解是适定的,即存在,唯一,而且连续依赖于初值0u 。

通常情况下,(i.1)的精确解不可能用简单的解析表达式给出,只能求近似解。

本章讨论常微分方程最常用的近似数值解法-差分方法。

先来讨论最简单的Euler 法。

为此,首先将求解区域[0,]T 离散化为若干个离散点:0110N N t t t t T -=<<<<= (i.3) 其中n t hn =,0h >称为步长。

在微积分课程中我们熟知,微商(即导数)是差商的极限。

反过来,差商就是微商的近似。

在0t t =处,在(i.1a )中用向前差商10()()u t u t h -代替微商du dt ,便得 10000()()(,())u t u t hf t u t ε=++如果忽略误差项0ε,再换个记号,用i u 代替()i u t 便得到1000(,)u u hf t u -=一般地,我们有1Euler (,), 0,1,,1n n n n u u hf t u n N +=+=-方法: (i.4) 从(i.1b) 给出的初始值0u 出发,由上式可以依次算出1,,N t t 上的差分解1,,N u u 。

下面我们用数值积分法重新导出 Euler 法以及其它几种方法。

为此,在区间1[,]n n t t +上积分常微分方程(i.1a ),得11()()(,())n n t n n t u t u t f t u t dt ++=+⎰ (i.5)用各种数值积分公式计算(i.5)中的积分,便导致各种不同的差分法。

(完整版)大连理工大学高等数值分析有限元简述-2017

(完整版)大连理工大学高等数值分析有限元简述-2017

椭圆与抛物微分方程的有限元法空间m H 作为例子,我们将考虑区间()0,1I =上的微分方程。

用2()L I 表示在I 上勒贝格平方可积函数的集合,()m H I 表示本身以及直到m 阶的导数都属于2()L I 的函数的集合。

我们下面用到的主要是1()H I 。

这里所说的导数准确地说是应该是广义导数,对此我们不予详细说明,只需知道比如说,连续的分片线性函数(折线函数)就属于1()H I ,其广义导数是分片常数函数。

另外,我们还用到空间}0)0(),({)(11=∈=v I H v I H E。

(空间=函数集合。

)微分方程 考虑两点边值问题(), (0,1)pu qu f x ''-+=∈(1) (0)0u = (2)0)1(='u(3)其中, , p q f 都是区间)1,0(上的光滑函数,0≥q , 并且0p p ≥,0p 是一个正常数。

用)(1I H E中任一函数v 乘(1)式两端,并在]1,0[上积分,得10[()]0pu v quv fv dx ''-+-=⎰ (4) 利用分部积分,并注意0)1(='u 和0)0(=v ,得⎰⎰''+'-=''-10110|)(dx v u p v u p vdx u p ⎰''=10dx v u p以此代入到(4)得到⎰=-+''100)(dv fv quv v u p (5) 为了方便,定义()10,w v w vdx =⋅⎰ (7)),(),(),(v qw v w p v w a +''=(8)则相应于微分方程(1)-(3)的变分方程为:求)(1I H u E∈满足),(),(v f v u a = )(1I H v E∈∀ (9)注意在(9)中不出现二阶导数。

可以证明,满足微分方程(1)-(3)的光滑解一定满足变分方程(9)。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

i.常微分方程初值问题数值解法i.1 常微分方程差分法考虑常微分方程初值问题:求函数()u t 满足(,), 0du f t u t T dt=<≤ (i.1a ) 0(0)u u = (i.1b)其中(,)f t u 是定义在区域G : 0t T ≤≤, u <∞上的连续函数,0u 和T 是给定的常数。

我们假设(,)f t u 对u 满足Lipschitz 条件,即存在常数L 使得121212(,)(,), [0,]; ,(,)f t u f t u L u u t T u u -≤-∀∈∈-∞∞ (i.2) 这一条件保证了(i.1)的解是适定的,即存在,唯一,而且连续依赖于初值0u 。

通常情况下,(i.1)的精确解不可能用简单的解析表达式给出,只能求近似解。

本章讨论常微分方程最常用的近似数值解法-差分方法。

先来讨论最简单的Euler 法。

为此,首先将求解区域[0,]T 离散化为若干个离散点:0110N N t t t t T -=<<<<= (i.3) 其中n t hn =,0h >称为步长。

在微积分课程中我们熟知,微商(即导数)是差商的极限。

反过来,差商就是微商的近似。

在0t t =处,在(i.1a )中用向前差商10()()u t u t h -代替微商du dt ,便得 10000()()(,())u t u t hf t u t ε=++如果忽略误差项0ε,再换个记号,用i u 代替()i u t 便得到1000(,)u u hf t u -=一般地,我们有1Euler (,), 0,1,,1n n n n u u hf t u n N +=+=-方法: (i.4) 从(i.1b) 给出的初始值0u 出发,由上式可以依次算出1,,N t t 上的差分解1,,N u u 。

下面我们用数值积分法重新导出 Euler 法以及其它几种方法。

为此,在区间1[,]n n t t +上积分常微分方程(i.1a ),得11()()(,())n n t n n t u t u t f t u t dt ++=+⎰ (i.5)用各种数值积分公式计算(i.5)中的积分,便导致各种不同的差分法。

例如,若用左矩形公式就得到 Euler 法(i.4)。

如果用右矩形公式,便得到下面的:111Euler (,), 0,1,,1n n n n u u h f t u n N +++=+=-隐式方法: (i.6) 类似地,如果用梯形公式,就得到111Euler [(,)(,)], 0,1,,12n n n n n n h u u f t u f t u n N +++=++=-改进的方法 (i.7) 当(,)f t u 关于u 是非线性函数的时候,不能由(i.6)或 (i.7) 从n u 直接算出1n u +,称这一类方法为隐式,通常采用某种迭代法求解。

例如,将一般的隐式方法写成11(,,)n n n n u F t u u ++= (i.8) 则可以利用如下的迭代法由n u 算出1n u +:11101(,,), 0,1, k k n n n n n n u F t u u k u u ++++⎧==⎨=⎩ (i.9)关于k 的迭代通常只需进行很少几步就可以满足精度要求了。

为了避免对隐式方法进行迭代的麻烦,比如说对于改进的Euler 方法(i .7),可以采用某种预估法近似算出11(,)n n f t u ++,然后再用(i .7)作校正,这就导致所谓预估校正法。

下面给出一个例子:1111111(,)2(,)()2n n n n n n n n n n n n n u f t u u u hu u f t u h u u u u +-+++++'⎧=⎧⎪⎪'=+⎨⎪⎪⎪'=⎨⎩⎪⎪''=++⎪⎩预估: 校正: (i.10) 这是一个多步法,即计算节点1n t +上的近似值1n u +时,除了用到前一点的近似值n u 之外,还要用到1n u -,甚至可能用到2,n u -。

而用前面的各种Euler 法计算节点1n t +上的近似值1n u +时,只用到n u ,这样的方法称之为单步法。

下面给出另一个多步法的例子。

在区间2[,]n n t t +上积分(i.1a ),得22()()(,())n n t n n t u t u t f t u t dt ++=+⎰用Simpson 公式(即把被积函数看作二次函数)近似计算积分,便得到212Milne (4), 0,1,,23n n n n n h u u f f f n N +++=+++=-法: (i.11) 用多步法(i.10)或(i.11)计算时,必须先用某种单步法由0u 计算出1u ,称为造表头。

然后再逐次算出23,,,N u u u 。

一般说来,多步法比Euler 法等简单的单步法精度要高一些。

下面我们讨论一类所谓Runge-Kutta 法。

他们是单步法,但是其精度可以与多步法比美。

最常用的是下面的标准Runge-Kutta 法:标准Lunge-Kutta 法⎪⎪⎪⎪⎩⎪⎪⎪⎪⎨⎧++++=++=++=++==+)22(6),()2,2()2,2(),(432113423121K K K K h u u hK u h t f K hK u h t f K hK u h t f K u t f K n n n n n n n n n n (i.12)从几何上, Runge-Kutta 法可以粗略地解释为:在区间1[,]n n t t +中选取若干个点(可以重复)1211n k k n t t ττττ-+=<≤≤≤=,仅仅利用在区间1[,]n n t t +内可以得到的所有信息,依次给出函数(,())f t u t 在这些点上尽可能精确的的近似值12,,,k K K K ,然后把它们组合起来,尽可能精确地近似计算(i.5)中的积分。

下面介绍Runge-Kutta 法的一般构造方式。

选定常数k ,令1(,)n n K f t u =22211(,)n n K f t h u K αβ=++33311322(,)n n K f t h u K K αββ=+++1122,11(,)k n k n k k k k k K f t h u K K K αβββ--=+++++ 11122()n n k k u u h K K K ωωω+=++++ (i.14)选取这些待定常数,,i ij m αβω的原则是:将(i.14)在(,)n n t u 作Taylor 展开,然后按照h 的幂重新整理,使得231123112!31n n u u h h h γγγ+=++++ (i.15) 与微分方程(i.1a )的解()u u t =在n t t =处的Taylor 展式23111()()2!3!n n n n n u t u t f h f h f h +'''=++++ (i.16) 有尽可能多的项重合,即要求123, , , n n n f f f γγγ'''=== 这里(,)(,)(,)(,), n n n n n n t u t u df t u f f t u f dt ='==等等。

按照(i.14)构造出的都是显式Runge-Kutta 方法,每一个i K 可以依次显式地算出。

如果在某一个i K 的表达式中出现j K ,其中j i >,则导致隐式Runge-Kutta 方法,可以迭代求解。

一般来说,隐式Runge-Kutta 方法的稳定性更好一些。

i.2 常微分方程组与高阶常微分方程先来考虑下面的常微分方程组初值问题11122212121100(,,,,)(,,,,)(,,,,)(0),,(0)m m mm m m m du f t u u u dt du f t u u u dt du f t u u u dtu u u u ⎧=⎪⎪⎪=⎪⎪⎨⎪⎪=⎪⎪⎪==⎩ (i.17)利用向量记号,上式可以改写为()(,)(0)t t u '=⎧⎨=⎩u f u u (i.18) 上节中各方法都可以直接应用到常微分方程组(i.18)。

例如,Euler 方法成为 1(,), 0,1,,1n n n n h t n N +=+=-u u f u 再来考虑高阶常微分方程111211(,,,,), 0(0),,(0),(0)m m m m m m m d u du d u f t u t dt dt dt d u du v v u v dt dt ----⎧=>⎪⎪⎨⎪===⎪⎩ (i.19) 这时,可以令1121(,,,)(,,)m m m du d u u u u u dt dt --=≡u (i.20) 231212(,,,)(,,,,(,,,,))m m m f f f u u u f t u u u ==f (i.21) 012(,,,)m v v v =u (i.22)于是可以把高阶常微分方程(i.19)化成一阶常微分方程组(i.18)。

i.3 收敛性与稳定性截断误差 粗略地说,截断误差可以定义为将微分方程解带入到差分方程后得到的误差,代表了微分方程与差分方程之间的误差。

例如,由Taylor 展式和微分方程(i.1a)得到221(,())()()()()()(,())2!2!nn n n n n n n t h h df t u t u t u t hu t u u t hf t u t dt ττ+='''=++=++ 其中n τ是区间1[,]n n t t +上某个常数。

与Euler 法1(,)n n n n u u hf t u +=+相比较,定义余项2(,())2!nt h df t u t dt τ=为Euler 法的截断误差,它关于h 是2阶的,记为2()h O 。

一般地,将上节中讨论的单步法写成1(,,)0n n F u u h +=,则可以定义截断误差为1((),(),)n n F u t u t h +。

对于每一个差分法在适当的点(例如n t t =,1n t +或11()2n n t t ++)作Taylor 展开,就可以得到截断误差的阶。

对多步法可以类似处理。

上节中各方法的截断误差阶分别为:相容性 一个差分方法称为相容的,如果其截断误差至少是2阶的。

收敛性 实用中我们更关心的是近似解的收敛性,即0h →时,是否有()n n u u t →。

在适当的条件下,例如步长h 足够小,右端函数f 和解u 足够光滑等等,可以证明以上讨论的各方法都是收敛的,并且有估计式()pn n u t u Ch -≤ (i.23) 其中常数C 与u 和f 有关,与h 和n 无关;阶数p 等于截断误差的阶数减去1。

相关文档
最新文档