实验4常微分方程数值解
常微分方程初值问题数值解法

常微分方程初值问题的数值解法在自然科学、工程技术、经济和医学等领域中,常常会遇到一阶常微分方程初值问题:(,),,(),y f x y a x b y a y '=≤≤⎧⎨=⎩ (1) 此处f 为,x y 的已知函数,0y 是给定的初始值。
本章讨论该问题的数值解法,要求f 在区域{(,)|,}G x y a x b y =≤≤<∞内连续,并对y 满足Lipschitz 条件,从而初值问题(1)有唯一的连续可微解()y y x =,且它是适定的。
1 几个简单的数值积分法1.1 Euler 方法(1)向前Euler 公式(显式Euler 公式)10(,),0,1,2,,(),n n n n y y hf x y n y y a +=+=⎧⎨=⎩(2) 其中h 为步长。
由此便可由初值0y 逐步算出一阶常微分方程初值问题(1)的解()y y x =在节点12,,x x 处的近似值12,,y y 。
该公式的局部截断误差为2()O h ,是一阶方法。
(2)向后Euler 公式(隐式Euler 公式)1110(,),0,1,2,,(),n n n n y y hf x y n y y a +++=+=⎧⎨=⎩(3) 这是一个隐格式,也是一阶方法。
这类隐格式的计算比显格式困难,一般采用迭代法求解。
首先用向前Euler 公式提供迭代初值,然后迭代计算:(0)1(1)()111(,),(,),0,1,2,n n n n k k n n n n y y hf x y y y hf x y k +++++⎧=+⎨=+=⎩ (4)1.2 梯形方法1110[(,)(,)],2(),(0,1,2,)n n n n n n h y y f x y f x y y y a n +++⎧=++⎪⎨⎪=⎩= (5) 这也是一个隐格式,是二阶方法。
一般也采用迭代法求解。
迭代公式如下:(0)1(1)()111(,),[(,)(,)],0,1,2,2n n n n k k n n n n n n y y hf x y h y y f x y f x y k +++++⎧=+⎪⎨=++=⎪⎩ (6)1.3 改进的Euler 方法11110(,),[(,)(,)],0,1,2,,2(),n n n n n n n n n n y y hf x y h y y f x y f x y n y y a ++++⎧=+⎪⎪=++=⎨⎪=⎪⎩(7) 为了便于上机编程计算,(7)可改写为110(,),(,),0,1,2,,1(),2(),p n n n cn n p n p c y y hf x y y y hf x y n y y y y y a ++=+⎧⎪=+⎪⎪=⎨=+⎪⎪=⎪⎩(8) 该格式是显式,也是二阶方法。
数值分析常微分方程求解实验

实验报告
一、实验目的
解初值问题各种方法比较。
二、实验题目
给定初值问题
⎪⎩⎪⎨⎧=≤<+=,
0)1(,21y x xe x y dx dy x , 其精确解为)(e e x y x -=,按
(1)改进欧拉法,步长01.0,1.0==h h ;
(2)四阶标准龙格-库塔法,步长1.0=h ;
求在节点)10,...,2,1(1.01=+=k k x k 处的数值解及其误差,比较各个方法的优缺点。
三、实验原理
改进欧拉法程序,四阶标准龙格-库塔法程序。
四、实验内容及结果
五、实验结果分析
实验2中改进欧拉法和四阶标准龙格-库塔法的比较:
结果的第一个ans是x的值与对应的y的值,第二个ans是精确解x的对应值y,第三个ans 是与精确值的误差百分数。
通过误差百分数的比较,可以明显的发现改进欧拉法比四阶标准龙格-库塔法更精确。
数值计算方法实验报告

数值计算方法实验报告一、实验介绍本次实验是关于数值计算方法的实验,旨在通过计算机模拟的方法,实现对于数值计算方法的掌握。
本次实验主要涉及到的内容包括数值微积分、线性方程组的求解、插值与拟合、常微分方程的数值解等。
二、实验内容1. 数值微积分数值微积分是通过计算机模拟的方法,实现对于微积分中的积分运算的近似求解。
本次实验中,我们将会使用梯形公式和辛普森公式对于一定区间上的函数进行积分求解,并比较不同公式的计算误差。
2. 线性方程组的求解线性方程组求解是数值计算领域中的重要内容。
本次实验中,我们将会使用高斯消元法、LU分解法等方法对于给定的线性方程组进行求解,并通过比较不同方法的计算效率和精度,进一步了解不同方法的优缺点。
3. 插值与拟合插值与拟合是数值计算中的另一个重要内容。
本次实验中,我们将会使用拉格朗日插值法和牛顿插值法对于给定的数据进行插值求解,并使用最小二乘法对于给定的函数进行拟合求解。
4. 常微分方程的数值解常微分方程的数值解是数值计算中的难点之一。
本次实验中,我们将会使用欧拉法和龙格-库塔法等方法对于给定的常微分方程进行数值解的求解,并比较不同方法的计算精度和效率。
三、实验结果通过本次实验,我们进一步加深了对于数值计算方法的理解和掌握。
在数值微积分方面,我们发现梯形公式和辛普森公式都能够有效地求解积分,但是辛普森公式的计算精度更高。
在线性方程组求解方面,我们发现LU分解法相对于高斯消元法具有更高的计算效率和更好的数值精度。
在插值与拟合方面,我们发现拉格朗日插值法和牛顿插值法都能够有效地进行插值求解,而最小二乘法则可以更好地进行函数拟合求解。
在常微分方程的数值解方面,我们发现欧拉法和龙格-库塔法都能够有效地进行数值解的求解,但是龙格-库塔法的数值精度更高。
四、实验总结本次实验通过对于数值计算方法的模拟实现,进一步加深了我们对于数值计算方法的理解和掌握。
在实验过程中,我们了解了数值微积分、线性方程组的求解、插值与拟合、常微分方程的数值解等多个方面的内容,在实践中进一步明确了不同方法的特点和优缺点,并可以通过比较不同方法的计算效率和数值精度来选择合适的数值计算方法。
常微分方程(组)的数值解法

刚性常微分方程组求解
function demo1 figure ode23s(@fun,[0,100],[0;1]) hold on, pause ode45(@fun,[0,100],[0;1]) %-------------------------------------------------------------------------function f=fun(x,y) dy1dx = 0.04*(1-y(1))-(1-y(2)).*y(1)+0.0001*(1-y(2)).^2; dy2dx = -1e4*dy1dx + 3000*(1-y(2)).^2; f = [dy1dx; dy2dx];
解算指令的使用方法
[T,Y]=ode45(@fun, TSPAN,Y0) 输出变量T为返回时间列向量;解矩阵Y的每一行对应于T的 一个元素,列数与求解变量数相等。
@fun为函数句柄,为根据待求解的ODE方程所编写的ode文
件(odefile); TSPAN=[T0 TFINAL]是微分系统y'=F(t,y)的积分区间;Y0 为初始条件
2.3 常微分方程(组)的数值解法
知识要点
常微分方程初值问题---ode45,0de23
微分方程在化工模型中的应用
•间歇反应器的计算 •活塞流反应器的计算
•全混流反应器的动态模拟
•定态一维热传导问题
•逆流壁冷式固定床反应器一维模型
•固定床反应器的分散模型
Matlab常微分方程求解问题分类
边值问题:
ode解算指令的选择(2)
2.根据常微分方程组是否为刚性方程
y ' Ay b( x) y ( x0 ) y0
微分方程初值问题的数值解法

积分法:
yk 1 yk h f ( xk , yk ) y ( x0 ) y0
积分项利用矩形公式计算
(1) y( xk 1 ) y( xk )
xk 1
xk
f (t , y(t ))dt
(★)
xk 1
xk
f (t , y(t ))dt h f ( xk , yk ) y( xk 1 ) y( xk ) h f ( xk , yk )
引言
初值问题的数值解法:求初值问题的解在一系列节点的值 y ( xn )的近似值 yn 的方法.本章数值解法的特点:都是采用“步进 式”,即求解过程顺着节点排列的次序一步步向前推进. 常微分方程初值问题: dy f ( x, y ), x [a, b] dx y ( x0 ) y0
替 f (x , y)关于 y 满足Lipschitz条件. 除了要保证(1)有唯一解外,还需保证微分方程本身是稳定的,即 (1)的解连续依赖于初始值和函数 f (x , y). 也就是说, 当初始值 y0 及函数 f (x , y)有微小变化时, 只能引起解的微小变化.
注: 如无特别说明,总假设(1)的解存在唯一且足够光滑. 在 f 连续有界, 则 f (x , y)对变量 y 可微的情形下, 若偏导数 y 可取L为
也称折线法 x
2. 梯形法
若采用梯形公式计算(★)中的积分项,则有 h y ( xk 1 ) y ( xk ) [ f ( xk , y ( xk )) f ( xk 1 , y ( xk 1 ))] 2 h yk 1 yk [ f ( xk , yk ) f ( xk 1 , yk 1 )] 2 称之为梯形公式.这是一个隐式公式,通常用迭代法求解.具体做 法: (0) (0) 先用Euler法求出初值 yk ,1 即 ,将其代入梯形公式 yk 1 yk h f ( xk , yk ) 的右端,使之转化为显式公式,即 h ( l 1) (l ) yk 1 yk [ f ( xk , yk ) f ( xk 1 , yk (☆ ) 1 )] 2
常微分方程中的数值方法

常微分方程中的数值方法常微分方程是数学中的一个重要分支。
它主要研究的对象是随时间变化的函数。
在实际应用中,我们需要求解这些函数的解析解,但通常情况下,解析解并不容易得到,甚至是不可能得到。
因此,我们需要使用数值方法来求解这些函数的数值近似解。
在本文中,我们将介绍常微分方程中的数值方法。
一、欧拉法欧拉法是常微分方程数值解法中最基本的一种方法。
它是根据欧拉公式推导而来的。
具体地,我们可以将一阶常微分方程dy/dt=f(t,y)写成如下形式:y(t+h)=y(t)+hf(t,y(t))其中,h是步长,f(t,y)是t时刻y的导数。
欧拉法就是通过上面的公式进行逐步逼近,然后得到最终的数值解。
欧拉法的计算过程非常简单,但所得到的解可能会出现误差。
这是因为欧拉法忽略了f(t+h,y(t+h))和f(t,y(t))之间的变化。
因此,我们需要使用更为精确的数值方法来解决这个问题。
二、改进欧拉法为了解决欧拉法中的误差问题,我们可以使用改进欧拉法。
改进欧拉法又称作四阶龙格-库塔法。
它的基本思想是对欧拉法公式进行改进,以提高计算精度。
具体地,根据龙格-库塔公式,可将改进欧拉法表示为:y(t+h)=y(t)+1/6(k1+2k2+2k3+k4)其中,k1=h*f(t,y)k2=h*f(t+h/2,y+k1/2)k3=h*f(t+h/2,y+k2/2)k4=h*f(t+h,y+k3)改进欧拉法的计算过程比欧拉法要复杂些,但所得到的数值解比欧拉法更精确。
这种方法适用于一些特殊的问题,但在求解一些更为复杂的问题时,还需要使用其他的数值方法。
三、龙格-库塔法龙格-库塔法是求解常微分方程中数值解的常用方法之一。
它最常用的是四阶龙格-库塔法。
这种方法的基本思想是使用四个不同的斜率来计算数值解。
具体地,我们可以将四阶龙格-库塔法表示为:y(t+h)=y(t)+1/6(k1+2k2+2k3+k4)其中,k1=h*f(t,y)k2=h*f(t+h/2,y+k1/2)k3=h*f(t+h/2,y+k2/2)k4=h*f(t+h,y+k3)与改进欧拉法相比,龙格-库塔法的计算复杂度更高,但所得到的数值解更为精确。
常微分方程的数值解法
主要内容
§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 进行积分
内部联系非常复杂
其状态随着 时间、地点、条件 的不同而不同
常微分方程的数值解与解析解
一、 常微分方程的解析解常微分方程的解析解也就是常微分方程的精确解,也称为常微分方程的符号解;一般可理解为求微分方程的通解或者特解的解析式或表达式;但只有少数的微分方程存在解析解。
在MA TLAB 中,由函数dsolve()求解常微分方程(组)的解析解,其具体格式如下: X=dsolve(‘方程1’,‘方程2’,…‘方程n ’,‘初始条件’,‘自变量’)函数dsolve 用来解符号常微分方程、方程组,如果没有初始条件,则求出通解,如果有初始条件,则求出特解。
例1:求解常微分方程1dy dx x y =+的MA TLAB 程序为:dsolve('Dy=1/(x+y)','x'),注意,系统缺省的自变量为t ,因此这里要把自变量写明。
结果为:-lambertw(-C1*exp(-x-1))-x-1其中:Y=lambertw(X)表示函数关系Y*exp(Y)=X 。
例2:求解常微分方程2'''0yy y -=的MA TLAB 程序为:Y2=dsolve('y*D2y-Dy^2=0’,’x’) 结果为:Y2 =[ exp((x+C2)/C1)][ C2]我们看到有两个解,其中一个是常数。
例3:求常微分方程组253t tdx x y e dt dy x y e dt ⎧++=⎪⎪⎨⎪--=⎪⎩通解的MA TLAB 程序为:[X,Y]=dsolve('Dx+5*x+y=exp(t),Dy-x-3*y=exp(2*t)','t')例4:求常微分方程组020210cos ,224,0t t t dx dy x t x dt dt dx dy y e y dt dt =-=⎧+-==⎪⎪⎨⎪++==⎪⎩通解的MA TLAB 程序为:[X,Y]=dsolve('Dx+2*x-Dy=10*cos(t),Dx+Dy+2*y=4*exp(-2*t)','x(0)=2','y(0)=0')二、 常微分方程的数值解在生产和科研中所处理的微分方程往往很复杂且大多得不出一般解。
常微分方程的数值解法全文
第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 <<.。
常微分方程初值问题的数值解法
常微分方程初值问题数值解法初值问题:即满足初值条件的常微分方程的解y′=f(x,y),x∈[x0,b]y(x0)=y0.定理1(利普希茨条件)若存在正数L,使得对任意,y1,y2,有|f(x,y1)−f(x,y2)|≤L|(y1−y2)|定理2(解存在性)①若函数f在方区域x∈[a,b],y∈R连续,②函数f关于y 满足利普希茨条件,则对任意x∈[a,b],常微分方程存在唯一的连续可微数值解.两类问题:①单步法---计算下一个点的值yn+1只需要用到前面一个点的值yn②多步法---计算下一个点的值yn+1需要用到前面l个点的值yl1、欧拉法---下一个点的计算值等于前一个点的计算值加上步长乘以前一个点的函数值•具体过程一些批注:显式欧拉方程指下一步要计算的值,不在迭代方程中;隐式欧拉方程指下一步要计算的值,在迭代方程中。
怎么计算隐式欧拉方程----要借助显示欧拉迭代计算---一般用迭代法-----迭代---将微分方程在区间[xn,xn+1]进行积分,然后函数f进行近似,即可得到迭代方程-----迭代方程收敛性?由函数关于y满足利普希茨条件,可以推出迭代公式收敛。
•局部截断误差:假设前n步误差为0,我们计算第n+1步的误差,将次误差称为局部截断误差,且局部误差为O(hp+1)•p阶精度:由理论证明:若局部误差阶的时间复杂度为O(hp+1),则整体误差阶为O(hp)我们称公式精度为p。
•显示欧拉法与隐式欧拉法•梯形方法----将显式欧拉迭代方程与隐式欧拉迭代方程做一下加权平均,构造的计算公式.•改进的欧拉方法---思想:因为梯形公式是隐式公式,将显式欧拉公式对下一步的计算值进行预估,用梯形公式对下一步的计算值进行校正.2、龙格-库塔方法思想:根据Lagrange中值定理,下一次的计算值可以用前一次的计算值加上h乘以前一个点的斜率;而这个斜率用该区间上的多个点的斜率的算数平均来逼近。
注意:怎么计算任意斜率Ki?第i个点的斜率Ki有微分方程可以算出f′=f(xn,yn)所以要算的f(xn,yn)值,由欧拉法即可算出, yn+1=yn+hf′•2阶-龙格-库塔方法----类似改进的欧拉法根据Lagrange中值定理,下一次的计算值可以用前一次的计算值加上h乘以斜率;而这个斜率用区间上的端点和中点的斜率的算数平均来逼近。