微分方程常用的两种数值解法：欧拉方法与龙格—库塔法

四川师范大学本科毕业论文微分方程常用的两种数值解法：欧拉方法与龙格—库塔法学生姓名XXX院系名称数学与软件科学学院专业名称信息与计算科学班级2006级4班学号20060640XX指导教师Xxx四川师范大学教务处二○一○年五月微分方程常用的两种数值解法：欧拉方法与龙格—库塔法学生姓名：xxx指导教师：xx【内容摘要】微分方程是最有生命力的数学分支，在自然科学的许多领域中，都会遇到常微分方程的求解问题。当前计算机的发展为常微分方程的应用及理论研究提供了非常有力的工具，利用计算机解微分方程主要使用数值方法，欧拉方法和龙格——库塔方法是求解微分方程最典型常用的数值方法。本文详细研究了这两类数值计算方法的构造过程，分析了它们的优缺点，以及它们的收敛性，相容性，及稳定性。讨论了步长的变化对数值方法的影响和系数不同的同阶龙格—库塔方法的差别。通过编制C程序在计算机上实现这两类方法及对一些典型算例的结果分析比较，能更深切体会它们的功能，优缺点及适用场合，从而在实际应用中能对不同类型和不同要求的常微分方程会选取适当的求解方法。关键词：显式单步法欧拉（Euler）方法龙格—库塔（Runge—Kutta）方法截断误差收敛性Twocommonlyusednumericalsolutionofdifferentialequations：EulermethodandRunge-KuttamethodStudentName:XiongShiyingTutor：ZhangLi【Abstract】Thedifferentialequationisthemostvitalitybranchinmathematics.Inmanydomainsofnaturalscience,wecanmeettheordinarydifferentialequationsolutionquestion.Currently,thedevelopmentofcomputerhasprovidedtheextremelypowerfultoolfortheordinarydifferentialequationapplicationandthefundamentalresearch,thecomputersolvingdifferentialequationmainlyusesvaluemethod.TheEulermethodandtheRunge—Kuttamethodarethemosttypicalcommonlyvaluemethodtosolvethedifferentialequation.Thisarticledissectsthestructureprocessofthesetwokindsofvaluescommonlyvaluemethodtosolvetheanalysestheirgoodandbadpoints,totheirastringency,thecompatibility,andthestabilityhasmadetheproof.Atthesametime,thearticlediscussthelengthofstridetothenumericalmethodchanginginfluenceandthedifferenceofthecoefficientdifferentsamestepRunge—kuttamethod.ThroughestablishingCprogramonthecomputercanrealizethesetwokindofmethods,Anglicizingsomemodelsofcalculateexampleresultcansincerelyrealizetheirfunction,theadvantageanddisadvantagepointsandthesuitablesituation,thusthesuitablesolutionmethodcanbeselectedtosolvethedifferenttypeandthedifferentrequestordinarydifferentialequationinthepracticalapplication.Keywords:Explicitsingle-stepprocessEulermethodRunge—Kuttamethodtruncationerrorconvergence目录1微分方程常用的两种数值解法：欧拉方法与龙格—库塔法前言常微分方程的形成与发展是和力学、天文学、物理学以及其他科学技术的发展密切相关的。数学其他分支的新发展，如复变函数、群、组合拓扑学等，都对常微分方程的发展产生了深刻的影响，当前计算机的发展更是为常微分方程的应用及理论研究提供了非常有力的工具。牛顿在研究天体力学和机械力学的时候，利用了微分方程这个工具，从理论上得到了行星运动规律。后来，法国天文学家勒维烈和英国天文学家亚当斯使用微分方程各自计算出那时尚未发现的海王星的位置。这些都使数学家更加深信微分方程在认识自然、改造自然方面的巨大力量，微分方程也就成了最有生命力的数学分支。然而，我们知道，只有少数十分简单的微分方程能够用初等方法求得它们的解，多数情形只能利用近似方法求解。在常微分方程课中的级数解法，逐步逼近法等就是近似解法。这些方法可以给出解的近似表达式，通常称为近似解析方法。还有一类近似方法称为数值方法，它可以给出解在一些离散点上的近似值，利用计算机解微分方程主要使用数值方法。本文主要讨论一阶常微分方程初值问题00)(),(yxyyxfdxdy（1.1）在区间],[ba上的数值解法，其中),(yxf为关于x，y的已知函数，0y为给定的初始值，将上述问题的精确解记为)(xy。该问题常用的数值解法有：欧拉（Euler）方法、龙格—库塔（Runge—Kutta）方法及一些常用的线性多步法。本文重点介绍欧拉（Euler）方法和龙格—库塔（Runge—Kutta）方法。并对这两种方法编制程序，体会它们的功能、优缺点及适用场合，对不同类型常微分方程会选取适当的求解方法。１基本概念和准备知识一阶常微分方程初值问题是：(1.1.2))((1.1.1)),(00yxyyxfdxdy其中),(yxf是平面上某个区域D上的连续函数，式（1.1.1）的微分方程一般有无穷多个解，式（1.1.2）是确定解的初始条件，如果一元函数)(xy对一切bxa满2足（1）Dxyx))(,(；（2）00)(yxy；（3）'y存在而且))(,()('xyxfxy；则称)(xy是初值问题（1.1）在区间],[ba上的解。误差：假定在计算1ny时，用到的前一步的值ny是准确的(即))(nnxyy，把用)(nxy计算得到的近似值记为1ny，估计误差：1nE=y(xn+1)1ny，这种误差称为局部截断误差。如果不作这一假定，在每一步计算时除局部截断误差以外，还有由于前一步不准确而引起的误差，称为总体截断误差。收敛性：对于解初值问题的数值方法，我们希望它产生的数值解收敛于初值问题的准确解，“收敛性”的一般定义为：对于所有满足引理1.1条件的初值问题（1.1），如果有一种显式单步法：),,(1hyxhyynnnn产生的近似解，对于任意固定的nhxxn0，均有0limh)(nnxyy，则称该显式单步法是收敛的。相容性：显式单步法（1.2.1）称为与原微分方程相容，如果))(,()0),(,(xyxfxyx（1.2.3）成立。并称式（1.2.3）为相容性条件。稳定性：在实际计算中，一方面初始值0y不一定精确，往往带有一定误差，同时由于计算机字长有限，在计算过程中有舍入误差，而且这种误差在式（1.2.1）的递推过程中会传递下去，对以后的结果产生影响。因此要考虑舍入误差的积累是否会得到控制，也即要考虑数值方法的稳定性。当n时，若舍入误差引起的后果是有限的，则可以认为该方法是数值稳定的。2欧拉方法2.1欧拉方法简介对常微分方程初值问题（1.1）用数值方法求解时，我们总是认为（1.1）的解存在且唯一。欧拉方法是解初值问题的最简单、最原始的数值方法，它是显式单步法。下面介绍几种导出欧拉法的途径，每个途径皆可以推导出更为有效的数值方法。（1）Taylor展开在nx点将)()(1hxyxynn作Taylor展开得：3)(!2)(')()()(21nnnnnyhxhyxyhxyxy),(1nnnxx当h充分小时略去误差项)(22nyh，并注意到))(,()('nnnxyxfxy，得))(,()()(1nnnnxyxhfxyxy，以ny近似代替)(nxy，以1ny近似代替)(1nxy，且用“=”代替“”得差分方程初值问题：)(),(001xyyyxhfyynnnnnhxxn0，1,,2,1,0Nn，Nhab(2.1.1)用式(2.1.1)求解初值问题（1.1）的方法称为欧拉方法。（2）数值微分由导数的定义知，对于充分小的h))(,()(')()(1nnnnnxyxfxyhxyxy整理得))(,()()(1nnnnxyxhfxyxy，对此作相似的处理也可以得到欧拉方法(2.1.1)。（3）数值积分在区间],[1nnxx对))(,('xyxfy积分得dxxyxfxyxynnxxnn1))(,()()(1（2.1.2）用数值积分的左矩形公式计算式（2.1.2）右端的积分，得))(,()()(1nnnnxyxhfxyxy，于是同样可以得到欧拉方法(2.1.1)。（4）多项式插值利用解和其导数),(yxf在nx点的值ny，),(nnyxf作一次埃尔米特插值，得到关于h的插值多项式：),()(nnnyxhfyhp，用)(hp近似代替1ny就得到欧拉方法。（5）待定系数法在第n步，已知ny和),('nnnyxfy，利用这两个值估计出下一步的1ny，将已知的值与估计值作线性组合：'1nnnyyy，其中，为待定系数。为确定这两个参数，要求这个估计值对cxy)(和cxxy)(（c为常数）精确成立。如果4cxy)(，则0)('xy，得到方程：0cc，得1。如果cxxy)(，则cxy)('，这样有：)(1nnnxcccxcx，说明hxxnn1，这样估计值为：),('1nnnnnnyxhfyhyyy，即为欧拉方法。欧拉方法的几何意义：由点斜式得切线方程),()()(0000,0000yxfxxydxdyxxyyyx等步长为h，则hxx01，可由切线方程算出1y：),(0001yxhfyy，逐步计算出)(xyy在1nx点的值：),(1nnnnyxhfyy，,2,1,0n，用分段的折线逼近函数为“折线法”而非“切线法”，除第一个点是曲线上的切线，其它都不是。2.2欧拉方法的截断误差，收敛性，相容性，稳定性设)(nnxyy，把)(1nxy在nx处展开成泰勒级数，即)(!2)(')()(21nnnnyhxhyxyxy1,nnnxx再由欧拉方法：)(')())(,()(),(1nnnnnnnnnxhyxyxyxhfxyyxhfyy两式相减得欧拉方法的局部截断误差为：)(2121nnyhE，若)(xy在],[ba上充分)(xyy),(00yx00xix1ixnx),(iiyx),(11iiyx),(nnyx图1欧拉方法的几何意义5光滑，且令Mbax,max)(xy，则1nE)(2122hMh，故欧拉方法是一阶方法或具有一阶精度。欧拉方法的增量函数),,(hyx就是f，由引理1.3、引理1.4知当f满足Lipschitz条件时欧拉方法是收敛的而且是相容的。用欧拉方法求解典型方程