(8.1)--数值分析第七章常微分方程数值解.pdf

《(8.1)--数值分析第七章常微分方程数值解.pdf》由会员分享，可在线阅读，更多相关《(8.1)--数值分析第七章常微分方程数值解.pdf（91页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、数值分析课件数值分析课件同学们好！现在开始讲授第七章常微分方程初值问题的数值解法常微分方程初值问题的数值解法首先讲授数值分析课件数值分析课件讲授：讲授：求常微分方程初值问题的公式构造技术和有关知识。求常微分方程初值问题的公式构造技术和有关知识。重点论述重点论述:Euler方法方法、Runge-Kutta方法和线性多步法原理、构方法和线性多步法原理、构造、局部截断误差和稳定性等。造、局部截断误差和稳定性等。数值分析课件数值分析课件本章讲授内容本章讲授内容7.1 引例7.2 基本概念7.3 数值解法的误差、阶与绝对稳定性7.4 Euler 方法的有关问题7.5 Runge-Kutta方法7.6 线

2、性多法数值分析课件数值分析课件7.1 7.1 引例引例为什么需要求微分方程的数值解？北京交通大学北京交通大学 数值分析课件数值分析课件单摆运动可看作工程技术里某些振动问题的简化。单摆运动如右图所示，其中摆线长为l，小球质量为m，图中竖线是小球的平衡位置。若将小球偏离平衡位置一个角度，然后将重力作用于它，小球将会沿圆弧摆动。假设小球摆动只受重力作用，不考虑空气阻力，试考查单摆在=30的运动规律。本问题由牛顿第二定律，可得微分方程：22sindmlmgdt lmg北京交通大学北京交通大学 数值分析课件数值分析课件加入初始条件，本问题归为求解如下微分方程的特解(t)问题：220sin0(0)30,0

3、tdgdtlddt遗憾的是上面微分方程特解问题不能用通常的方法来求解，事实上它没有解析解。怎样求解不能用解析方法求解的微分方程初值问题就是本章要解决的问题。#北京交通大学北京交通大学 数值分析课件数值分析课件北京交通大学北京交通大学 数值分析课件数值分析课件同学们好！下面讲授常微分方程初值问题的数值解法的基本概念,该内容分2集讲授。第1集：北京交通大学北京交通大学 数值分析课件数值分析课件7.2 7.2 基本概念基本概念什么是微分方程数值解？北京交通大学北京交通大学 数值分析课件数值分析课件1、微分方程数值解有关的概念2、初值问题化为差分方程的方法北京交通大学北京交通大学 数值分析课件数值分析

4、课件微分方程数值解有关的概念北京交通大学北京交通大学 数值分析课件数值分析课件1、常微分方程初值问题一般形式0(,)()yf x yaxby ay已知函数初始条件2、数值方法3、数值解(),1,2,kky xy kn kx求y=y(x)在离散点处近似值的方法。ky数值解北京交通大学北京交通大学 数值分析课件数值分析课件01naxxxb设区间a,b上的一组节点为步长步长：1kkkhxx求数值解一般是从初值开始按递推方式求解，即：012nyy ayyyy b初值问题的解法有单步法和多步法：单步法单步法：1kkyy多步法多步法：1,11,k lk lkkkyyyyy 数值解法还有显格式和隐格式之分。

5、北京交通大学北京交通大学 数值分析课件数值分析课件北京交通大学北京交通大学 数值分析课件数值分析课件同学们好！下面继续讲授常微分方程初值问题的数值解法的基本概念。第2集：北京交通大学北京交通大学 数值分析课件数值分析课件初值问题化为差分方程的方法北京交通大学北京交通大学 数值分析课件数值分析课件基本思想基本思想用数值微分法、数值积分法、Taylor展开法等离散化方法将初值问题化为差分方程后再求解。初值问题化为差分方程的方法初值问题化为差分方程的方法1、用离散方法去掉微分方程中的导数得到近似离散化方程；3、在近似离散化方程中将“”换为“=”2、在近似离散化方程中用替换ky()ky x北京交通大学

6、北京交通大学 数值分析课件数值分析课件1 1、用数值微分法构造、用数值微分法构造得近似离散化方程 (,)(,)()(,()kkkyf x yy xf x y xy xf x y x用数值微分的2点前差公式代替导数，有11()()()(,()kkkkkkky xy xy xf xy xxx11()()(,()kkkkkky xy xf x y xxx北京交通大学北京交通大学 数值分析课件数值分析课件用初值问题化差分方程方法11(,),kkkkkkyyf x yhxxh,kkyy x整理得1(,)kkkkyyhf x y0,1,2,1kn11()()(,()kkkkkky xy xf xy xxx

7、是显式单步法。差分方程Euler公式北京交通大学北京交通大学 数值分析课件数值分析课件2 2、用数值积分法构造、用数值积分法构造111()()(,()kkkkxxkkxxy xy xy dxf x y x dx采用梯形公式()()2babaf x dxf af b得近似离散化方程111()()(,()(,()2kkkkkkhy xy xf xy xf xy x111(,)(,)2kkkkkkhyyf x yf xy是隐式单步法(,)yf x y对两边积分可得,kkyy x用数值积分梯形方法北京交通大学北京交通大学 数值分析课件数值分析课件3 3、用、用TaylorTaylor展开法构造展开法构

8、造取上式右端前2项，得近似离散化方程:函数 y(x)的Taylor展式为221()()()()()(,()(,()2!2!kkkkkkkkx xhhdy xy xhy xy xy xhf x y xf x y xdx1()()(,()kkkky xy xhf x y x,kkyy x做1(,)kkkkyyhf x y()(,()kkky xf x y x1kkxxh北京交通大学北京交通大学 数值分析课件数值分析课件北京交通大学北京交通大学 数值分析课件数值分析课件同学们好！下面讲授微分方程数值解法的误差、阶与绝对稳定性。该内容分2集讲授。7.3 7.3 数值解法的误差、阶与绝对稳定性（数值解法

9、的误差、阶与绝对稳定性（1）第1集：北京交通大学北京交通大学 数值分析课件数值分析课件7.3 7.3 数值解法的误差、阶与绝对稳定性数值解法的误差、阶与绝对稳定性怎样衡量求数值解公式的好坏？北京交通大学北京交通大学 数值分析课件数值分析课件2、数值解法的误差误差的区别3、数值解法的阶和稳定性阶和稳定性北京交通大学北京交通大学 数值分析课件数值分析课件的区别北京交通大学北京交通大学 数值分析课件数值分析课件解 y(x)在kx的准确值，:kyky x的近似值，是数值解，:ky计算ky给出的计算解，没有误差；有截断误差；有舍入误差。本章的解有3种：准确解、数值解、计算解。():ky x准确解北京交通

10、大学北京交通大学 数值分析课件数值分析课件数值解法的误差误差北京交通大学北京交通大学 数值分析课件数值分析课件1 1、单步法数学描述、单步法数学描述显式：1(,)kkkkyyhx y h隐式:11(,)kkkkkyyhxyyh增量函数：与 f(x,y)有关。11111(,)(,)(,)(,)22kkkkkkkkkkf x yf xyf x yf xyyyh 1(,)(,)kkkkkkyyh f x yf x y1(,)kkkkyyh f x y111(,)(,)2kkkkkkf x yf xyyyh北京交通大学北京交通大学 数值分析课件数值分析课件2 2、数值解法的误差、数值解法的误差整体截断

11、误差:111()kkkey xy局部截断误差:111()kkkTy xy1()(,(),)kkkkyy xhx y xh1(,)kkkkyyhx y h()kkyy x11()()(,(),)kkkkkTy xy xhxy xh北京交通大学北京交通大学 数值分析课件数值分析课件北京交通大学北京交通大学 数值分析课件数值分析课件同学们好！下面讲授微分方程数值解法的误差、阶与绝对稳定性。该内容分2集讲授。7.3 7.3 数值解法的误差、阶与绝对稳定性（数值解法的误差、阶与绝对稳定性（2）第2集：北京交通大学北京交通大学 数值分析课件数值分析课件数值解法的阶、稳定性阶、稳定性北京交通大学北京交通大学

12、 数值分析课件数值分析课件P阶方法(精度）:11()PkTO h方法的阶P越大，方法越好。主局部截断误差:1(,()Pkkg x y xh1 1、数值解法的阶、数值解法的阶1121()(,()()PPPkkkTO hg x y xhO h若11()()(,(),)kkkkkTy xy xhx y xh北京交通大学北京交通大学 数值分析课件数值分析课件例1、常微分方程初值问题的单步法为111(,)2(,)3kkkkkkhyyf x yf xy试求其局部截断误差主项并确定其精度阶。这里初值问题及步长h为1100(,),()kkkkyf x yhxxxxhy xy解：1111()()(,()2(,(

13、)3kkkkkkkhTy xy xf x y xf xy x()()()2()3kkkkhy xhy xy xy xh(),kkky xf x y x北京交通大学北京交通大学 数值分析课件数值分析课件函数y=f(x)在 x某邻域内有n+1阶导数，有泰勒公式()(1)1000000()()()()()()()()!(1)!nnnnfxff xf xfxxxxxxxnn介于 x 与 x 之间。()(1)2100000()()()()()()2!(1)!nnnnfxfxff xhf xf x hhhhnn()(1)2100000()()()()()()2!(1)!nnnny xyxyy xhy xy

14、 x hhhhnn复习：北京交通大学北京交通大学 数值分析课件数值分析课件23212()()()()()()()()()2!3!332!kkkkkkkkhhhhTy xhy xy xyy xy xh y xhy xy23232121()()()()()()236kkh y xO hh y xO hO h 故局部截断误差主项是21()6kh y x精度是一阶。同阶无穷小O运算规则：21212(),0,pppO ha ha happ()()(),pqpO hO hO hpq()()()2()3kkkkhy xhy xy xy xh北京交通大学北京交通大学 数值分析课件数值分析课件数值方法是绝对稳定

15、的：*kkkyy1kk绝对稳定域：1,kkhyy绝对稳定区间：绝对稳定域与复平面实轴的交。绝对稳定域越大，方法的绝对稳定性越好。绝对稳定域越大，方法的绝对稳定性越好。2 2、数值解法的稳定性、数值解法的稳定性试验方程：,yy为复数.,f x yy北京交通大学北京交通大学 数值分析课件数值分析课件北京交通大学北京交通大学 数值分析课件数值分析课件同学们好！下面讲授Euler方法的有关问题，该内容分2集讲授。第1集：北京交通大学北京交通大学 数值分析课件数值分析课件7.4 7.4 Euler方法的有关问题方法的有关问题怎样理解微分方程数值解的一些概念？北京交通大学北京交通大学 数值分析课件数值分析

16、课件1 1、Euler方法的几何意义及误差方法的几何意义及误差2 2、Euler方法的稳定域方法的稳定域3 3、改进的改进的Euler方法方法北京交通大学北京交通大学 数值分析课件数值分析课件Euler方法的几何意义及误差方法的几何意义及误差北京交通大学北京交通大学 数值分析课件数值分析课件11(,),kkkkkkyyhf x yhxxEuler公式公式00(,)()yf x yy xyEuler方法常称为折线法。1 1、EulerEuler方法的几何意义方法的几何意义北京交通大学北京交通大学 数值分析课件数值分析课件2、EulerEuler方法的误差方法的误差Euler方法的局部截断误差:2

17、()kTO hEuler方法的总体截断误差:()keO h说明Euler方法计算所得数值解可以逼近准确解，从而Euler方法是收敛的。11(,),kkkkkkyyhf xyhxx()0,0kkkey xyO hhkkyy x北京交通大学北京交通大学 数值分析课件数值分析课件北京交通大学北京交通大学 数值分析课件数值分析课件同学们好！下面讲授Euler方法的有关问题，该内容分2集讲授。第2集：北京交通大学北京交通大学 数值分析课件数值分析课件Euler方法的稳定域方法的稳定域北京交通大学北京交通大学 数值分析课件数值分析课件将Euler公式用于试验方程:yy得到1(1)kkkkyyh yh y设

18、计算ky,0,1,2,kk时有误差:11(1)()kkkkyh y1(1)kkh111kkh Euler方法绝对稳定域为：11hEuler方法绝对稳定区间为：2Re0h 1(,)kkkkyyhf x y北京交通大学北京交通大学 数值分析课件数值分析课件Euler方法绝对稳定区间为:2Re0h 若是负实数，则有当步长h满足:可保证20h Euler方法的计算绝对稳定。例2、设初值问题为 100,00.101yyxy 取h=0.025,用Euler方法求其数值解并与其准确解作比较。解:本题准确解为,100 xyeEuler法计算公式为101001 100,1kkkkyyhyh yy北京交通大学北京

19、交通大学 数值分析课件数值分析课件直接计算结果与对应准确值的误差有如下表0.0250.050.0750.1-1.52.25-3.3755.06251.58-2.243.37-5.06kkkkxyy xy计算结果波动大，不稳定，计算结果失真！原因：h的稳定性范围为200.02,100h 本题的 h=0.025 0.02，故计算结果不稳定。北京交通大学北京交通大学 数值分析课件数值分析课件改进的Euler方法北京交通大学北京交通大学 数值分析课件数值分析课件#3、改进的Euler方法预测：校正：1(,)kkkkyyhf x y111(,)(,)2kkkkkkhyyf x yf xy是预测校正公式。

20、上面组合公式称为改进的Euler方法。Euler方法梯形方法其简写形式：11(,)(,(,)2kkkkkkkkhyyf x yf xyhf x y易证改进的Euler方法是二阶方法。北京交通大学北京交通大学 数值分析课件数值分析课件北京交通大学北京交通大学 数值分析课件数值分析课件同学们好！下面讲授求微分方程初值问题数值解的Runge-Kutta方法。该内容分2集讲授。第1集：北京交通大学北京交通大学 数值分析课件数值分析课件7.5 7.5 Runge-Kutta 方法怎样构造求微分方程数值解高阶方法？北京交通大学北京交通大学 数值分析课件数值分析课件基本思想将微分方程初值问题转化为积分方程问

21、题，再对积分方程中的定积分使用待定的m点插值型求积公式构造高阶的函数展开模式以获得高阶求数值解方法。理论上，由Taylor展开法可以构造出解初值问题的高阶求数值解方法，但这涉及到要计算f(x,y)的高阶导数，很不方便。本节的R-K方法从函数本身着手，避开计算 f(x,y)的高阶导数来构造高阶数值方法。RungeRunge-KuttaKutta方法简称R-K方法。北京交通大学北京交通大学 数值分析课件数值分析课件1、R-K方法的一般公式2、R-K方法的构造北京交通大学北京交通大学 数值分析课件数值分析课件R-K方法的一般公式北京交通大学北京交通大学 数值分析课件数值分析课件(,)yf x y11

22、1()()(,()kkkkxxkkxxy xy xy d xf x y x d x用m个点插值型求积公式，得近似离散化方程11()()(,()mkkiiiiy xy xhc fy11,ikkkkx xhxx11111,(2,3,)mkkiiikkrrkrkrjjjyyhc KKf x yKfxa h yhb KrmR-K一般公式对两边积分，有ikixah取，并对取不同Taylor公式，得(,()iify北京交通大学北京交通大学 数值分析课件数值分析课件11111,(2,3,)mkkiiikkrrkrkrjjjyyhc KKf xyKfxa h yhb Krm利用Taylor公式对参数,rrji

23、a b c适当选择就可以构造高阶方法。R-K方法一般公式R-K方法的增量函数：1,miiix y hc K北京交通大学北京交通大学 数值分析课件数值分析课件二元Taylor展开公式221,(,)22!xyxxxyyyf xh ykf x yf hf kf hhkff k ,xxyyff x y xffx y xffx y式中 ,xxxxxyxyyyyyffx y xffx y xffx y21()2!xxxf xhf xf hf h一元Taylor展开公式北京交通大学北京交通大学 数值分析课件数值分析课件北京交通大学北京交通大学 数值分析课件数值分析课件同学们好！下面讲授求微分方程初值问题数值

24、解的Runge-Kutta方法。该内容分2集讲授。第2集：北京交通大学北京交通大学 数值分析课件数值分析课件R-K方法的构造北京交通大学北京交通大学 数值分析课件数值分析课件只讲授 m=2的R-K方法的构造过程：11122kkyyh c Kc Km=2的R-K公式122211,kkkkKf xyKf xa h yhb K其增量函数为 12221,x y xhc f x y xc f xa h y xhb f x y x有4个待定参数12112,a bc c北京交通大学北京交通大学 数值分析课件数值分析课件12221,kkkkkkkkx y xhc f x y xc f xa h y xhb f

25、 x y x,kkkkkTy xhy xhxy xh对增量函数做Taylor展开，有,x y xhy xh ,Ty xhy xhx y x h为展开方便，记 22332!2!xyyxhy xhy xy x hhO hy xhffffO h,xyyf x yyxfff北京交通大学北京交通大学 数值分析课件数值分析课件 12221,x y x hc f x y xc f xa h y xhb f x y x 21222121222221,xyxyx y xhc fcfa hfhb f fO hccfc a fc b f fhO h对增量函数做Taylor展开221,(,)22!xyxxxyyyf

26、xh ykf x yf hf kf hhkff k北京交通大学北京交通大学 数值分析课件数值分析课件2312222 2111122xyTccf hc afc bf f hO h取12222 21111,0,022ccc ac b 3TO h从方程组解出即得到一组二阶R-K公式。12212,c c ba因方程组有3个方程4个参数，有无穷多解。20ct 令122111,2ct abt 12222 2110.500.50ccc ac b12t 1,2kkkkkkkkhyyf x yf xh yhf x y取改进的Euler公式北京交通大学北京交通大学 数值分析课件数值分析课件R-K方法的阶与级的关系

27、计算f(x,y)的次数m:1 2 3 4 5 6 7对应方法的最高阶1 2 3 4 4 5 6经典经典R R-K K公式公式112341213243226,22,22,kkkkkkkkkkhyyKKKKKf x yhhKfxyKhhKfxyKKf xh yhK北京交通大学北京交通大学 数值分析课件数值分析课件北京交通大学北京交通大学 数值分析课件数值分析课件同学们好！下面讲授求微分方程初值问题数值解的线性多步法。该内容分2集讲授。第1集：北京交通大学北京交通大学 数值分析课件数值分析课件7.5 7.5 线性多步法线性多步法线性多步法的线性和多步体现在哪里？北京交通大学北京交通大学 数值分析课件

28、数值分析课件1211,k nk nkkkyyyyy 基本思想用未知函数 y(x)及其导数 y =f(x,y)在若干节点处函数值的线性组合来构造求数值解的公式。,ii 为常数，2211(,),=-,0kkkk iknnff xyxx ihn个y值11100,0,1,2,nnkik iik iiiyyhfk 线性n步法公式：北京交通大学北京交通大学 数值分析课件数值分析课件1、用数值积分法构造线性多步公式2、用Taylor展开法构造线性多步公式北京交通大学北京交通大学 数值分析课件数值分析课件用数值积分法构造线性多步公式北京交通大学北京交通大学 数值分析课件数值分析课件对 y=f(x,y)两边积分

29、，得11()()(,()k mk mxk mk mxy xy xf x y x dx 对右端的定积分选择 f(x,y(x)的不同插值函数代替被积函数 f(x,y(x)、对公式中的m取不同的正整数即可得到不同的线性多步法计算公式。这里只介绍Adams的做法。北京交通大学北京交通大学 数值分析课件数值分析课件AdamsAdams方法给出构造过程方法给出构造过程设初值问题的解y(x)在n个点11,kkk nx xx 的数值解为11,kkk ny yy,iiiff x y构造 f(x,y(x)的n-1次Lagrange插值多项式(1)1()(,()()()!nnnyf x y xLxxn1110()(

30、)nnink iiLxlx f1111kkk nkkk nxxxfff 考虑北京交通大学北京交通大学 数值分析课件数值分析课件积分中取 m=1，有定积分111(1)11()()()(,()()()()!kkkkkknxxxkkknnxxxyy xy xf x y x dxy xLx dxx dxn110()()nkkik iiy xy xhf110nkkik iiyyhfn步Adams显示（外推）公式1(+1)()()!kknxnnxyRx dxn11111100()()kkkknnxxnk iinik ixxiiLx dxflx dxhf1111,()kkxkkiinxhxxlx dxh北京

31、交通大学北京交通大学 数值分析课件数值分析课件采用同样的方式，可得1212kkk nkkk nxxxfff 若考虑1111kkk nkkk nxxxfff 1110nkkikiiyyhf n步Adams隐式（内插）公式n步Adams隐式和隐式公式的局部截断误差都是11nkTO h北京交通大学北京交通大学 数值分析课件数值分析课件例3、构造3阶Adams显示公式。做 f(x,y(x)的二次Lagrange插值多项式122121212112221()()()()()()()()()()()()()kkkkkkkkkkkkkkkkkkkkkxxxxxxxxxxxxL xfffxxxxxxxxxxxx

32、解：取三个相邻的等距节点对应的数据1212kkkkkkxxxfff,0,1,2k ikxxih i北京交通大学北京交通大学 数值分析课件数值分析课件 2,()f x y xL x1112120(1)(2)1()()()(2)(1)22kkxkkkkkxtty xy xL x dxhft tft tfdt12(23165)12kkkhfff得到3阶Adams外推公式112(23165)12kkkkkhyyfffkxxth 1112()()(,)()kkkkxxkkxxy xy xf x y x dxL x dx北京交通大学北京交通大学 数值分析课件数值分析课件其局部截断误差为11(4)4(4)1

33、20414(4)(4)0()()()(1)(2)()()3!3!(1)(2)3()()3!8kkxkkkxxxxxxxt ttTydxhydtt tthh ydty如果把本题改为构造3阶Adams隐式公式，应该怎样做？北京交通大学北京交通大学 数值分析课件数值分析课件北京交通大学北京交通大学 数值分析课件数值分析课件同学们好！下面讲授求微分方程初值问题数值解的线性多步法。该内容分2集讲授。第2集：北京交通大学北京交通大学 数值分析课件数值分析课件用Taylor展开法构造线性多步公式北京交通大学北京交通大学 数值分析课件数值分析课件先给出线性多步法的计算公式模式，再对局部截断误差表达式处作Tay

34、lor展开确定公式的系数。23(3)23(3)(4)()()()()()()()()()2!3!()()()()()()()()()2!3!k ikkkkkk ikkkkkihihy xy xihy xih y xy xyxihihy xy xihy xih y xyxyx基于Taylor展开的构造方法做法：常用的Taylor展开公式：,kkky xf xy x北京交通大学北京交通大学 数值分析课件数值分析课件例4设求初值问题线性3步公式具有如下形式1212()kkkkkyyh afbfcfh为步长，试求系数a,b,c 使该公式的阶数尽可能高，并写出其局部截断误差。解：局部截断误差为11211

35、221212,kkkkkkkkkkkkkkTy xy xh af xy xbf xy xcf xy xy xy xh ay xby xcy x,kkky xf x y x北京交通大学北京交通大学 数值分析课件数值分析课件 4234511123!4!kkkkky xy xhh yxyxhyxhO h 23445111222223!4!kkkkky xy xhhyxyxhyxhO h 423423441123!1122223!kkkkkkkkkhay xbh y xhyxyxhyxhO hch y xhyxyxhyxhO h11212kkkkkkTy xy xh ay xby xcy x 42345331514322222863kkkka b c hy xbc h yxbc h yxbc h yxO h 北京交通大学北京交通大学 数值分析课件数值分析课件因为有3个待定系数，选择截断误差的前3项的系数为零，得线性方程组301.5201.5 0.520a b cbcbc 解得2.25,0,0.75abc514308638bc是3阶方法，局部截断误差为 44538kkTh yxO h 423451331514322222863kkkkkTa b c hy xbc h yxbc h yxbc h yxO h 北京交通大学北京交通大学 数值分析课件数值分析课件