微分方程的分类及其数值解法
微分方程的分类及其数值解法
微分方程的分类：
含有未知函数的导数，如dy/dx=2x 、ds/dt=0.4都是微分方程。 一般的凡是表示未知函数、未知函数的导数与自变量之间的关系的方程，叫做微分方程。未知函数是一元函数的，叫常微分方程；未知函数是多元函数的叫做偏微分方程。微分方程有时也简称方程。
一、常微分方程的数值解法：
1、Euler 法：
00d (,), (1.1)d (), (1.2)
y f x y x y x y ?=???=? 001
(),(,),0,1,,1n n n n y y x y y hf x y n N +=??=+=-? (1.4) 其中0,n b a x x nh h N
-=+=. 用(1.4)求解(1.1)的方法称为Euler 方法。 后退Euler 公式???+==+++),,(),(111
00n n n n y x hf y y x y y 梯形方法公式
)].,(),([2
111+++++=n n n n n n y x f y x f h y y 改进的Euler 方法11(,),(,),1().2p
n n n c n n p n p
c y y hf x y y y hf x y y y y ++?=+??=+???=+???
2、Runge-Kutta 方法：
p 阶方法 ： 1()O h -=?总体截断误差局部截断误差
二阶Runge-Kutta 方法 ???
????++==++=+),,(),,(,2212
1211hk y h x f k y x f k k h k h y y n n n n n n
三阶及四阶Runge-Kutta 方法 ?????????+-+=++==+++=+).2,(),2,2(),,(),4(6213
1213211hk hk y h x f k k h y h x f k y x f k k k k h y y n n n n n n n n
标准的(或经典的)四阶R-K 方法?????????????++=++=++==++++=+).
,(),2,21(),2,2(),
,(),22(61342312143211hk y h x f k k h y h x f k k h y h x f k y x f k k k k k y y n n n n n n n n n n 3、线性多步法：
∑∑-=--=--++=101
11,k i k i i n i i n i n f h y y βα其中(,),(0,1,,1)i j j i f f x y i k α≡=- ，
(1,0,,1)j j k β=-- 为与n 无关的常数。由于按公式(3.1)计算1+n y 时，需要知道11,,,n n n k y y y --+ 这k 个值，所以称为k 步法，又)1,,1,0(-=-k i y i n 及
)1,,1,0,1(--=-k i f i n 都是线性的，所以称(3.1)为线性多步
四阶隐式Adams （外推）公式：
)].,(),(5),(19),(9[24
1d )1()1(),(6
1d )2()1(),(2
1d )2)(1)(1(),(2
1d )2)(1(),(6
1221111102210111010111----++----++++-++=+--+--++--+++=????n n n n n n n n n n n n n n n n n n n y x f y x f y x f y x f h y u u u u y x hf u u u u y x hf u u u u y x hf u u u u y x hf y y 四阶显式Adams （内插）公式：
11122331[55(,)59(,)37(,)9(,)]24
n n n n n n n n n n y y h f x y f x y f x y f x y +------=+
-+-
二、偏微分方程数值解法：
差分法：
差分法是求偏微分方程数值解的重要方法之一，它的主要做法是把偏微分方程中所有偏导数分别用差商代替，从而得到一代数方程组——差分方程，然后对差分方程求解，并以所得的解作为偏微分方程数值解。
为此，必须对区域进行剖分，用网格点来 u=0
代替连续区域，因此差分法亦称“网格法”。 1
我们用一个简单例子来说明差分法的
基本思想和具体要求。 U=0
取一边长为1的正方形均匀薄板， y u πsin = 上下侧面绝热，四周保持恒温（如图10 .1），
求板内各点的稳定温定分布。
这个总是如在数学物理方程中所知，它可
以化为拉普拉斯方程第一边值问题： 0 u=0 1
?????????====<<<<Ω=??+??=?====y
u u u u y x y u x u u x y y x πsin 0)10,10:(011002222 （1）、椭圆型方程的差分解法
椭圆型方程最简单的典型问题就是拉普拉斯方程
02222=??+??=?y
u x u u 和泊松方程
),(2222y x f y
u x u u =??+??=? （2）、抛物型方程的差分解法
抛物方程的最简单的是一维热传导方程：
T t a t x f x u a t u u ≤<>=??-??≡?00),(22 （10.35） 它的定解条件主要有以下两类：
（ⅰ）初值问题：（或称柯西Caucy 问题）
+∞<<∞-==x x u t )(0? (10.36) （ⅱ）边值问题（或称混合问题）
?????====βα?),(,),()(00t a u t a u x u n t （10.37）
求（10.35）满足（ⅰ）或（ⅱ）的解.
（3）、双曲型方程的差分解法
一阶双曲型方程的初值问题为
???=>+∞<<-∞=+)2()()0,()1()0,(0x x u t x au u x t ? a 为常数，亦称（1）为对流方程。称ξ=-at x 为（1）的特征线，ξ 为常数，沿特征线 u (x , t )的方向导数
0d ),(d d =+=+=t x u au t
t at du t u ξ 即u (x , t )沿特征线为常数，再由 （2.），得初值问题（1），（2）的解
)(),(at x t x u -=?
二阶波动方程
02=-xx tt u a u
若令v = u ，w = au x 则得一阶双曲型方程组
???=-=-0
0x t x t av w aw v 再令v w v v w u
-=+=~,~，则得
???=+=-0~~0~~x t x t v a v u a u
常微分方程边值问题的数值解法
第8章 常微分方程边值问题的数值解法 引 言 第7章介绍了求解常微分方程初值问题的常用的数值方法；本章将介绍常微分方程的边值问题的数值方法。 只含边界条件(boundary-value condition)作为定解条件的常微分方程求解问题称为常微分方程的边值问题(boundary-value problem). 为简明起见，我们以二阶边值问题为 则边值问题(8.1.1)有唯一解。 推论 若线性边值问题 ()()()()()(),, (),()y x p x y x q x y x f x a x b y a y b αβ'''=++≤≤?? ==? (8.1.2) 满足 (1) (),()p x q x 和()f x 在[,]a b 上连续； (2) 在[,]a b 上, ()0q x >, 则边值问题(8.1.1)有唯一解。 求边值问题的近似解，有三类基本方法： (1) 差分法(difference method)，也就是用差商代替微分方程及边界条件中的导数，最终化为代数方程求解； (2) 有限元法(finite element method)；
(3) 把边值问题转化为初值问题，然后用求初值问题的方法求解。 差分法 8.2.1 一类特殊类型二阶线性常微分方程的边值问题的差分法 设二阶线性常微分方程的边值问题为 (8.2.1)(8.2.2) ()()()(),,(),(), y x q x y x f x a x b y a y b αβ''-=<
微分方程数值解法
《微分方程数值解法》 【摘要】自然界与工程技术中的很多现象，可以归结为微分方程定解问题。其中，常微分方程求解是微分方程的重要基础内容。但是，对于许多的微分方程，往往很难得到甚至不存在精确的解析表达式，这时候，数值解提供了一个很好的解决思路。，针对于此，本文对常微分方程数值解法进行了简单研究，主要讨论了一些常用的数值解法，如欧拉法、改进的欧拉法、Runge —Kutta 方法、Adams 预估校正法以及勒让德谱方法等，通过具体的算例，结合MA TLAB 求解画图，初步给出了一般常微分方程数值解法的求解过程。同时，通过对各种方法的误差分析，让大家对各种方法的特点和适用范围有一个直观的感受。 【关键词】 常微分方程 数值解法 MA TLAB 误差分析 引言 在我国高校，《微分方程数值解法》作为对数学基础知识要求较高且应用非常广泛的一门课程，不仅 在数学专业，其他的理工科专业的本科及研究生教育中开设这门课程．近四十年来，《微分方程数值解法》不论在理论上还是在方法上都获得了很大的发展．同时，由于微分方程是描述物理、化学和生物现象的数学模型基础，且它的一些最新应用已经扩展到经济、金融预测、图像处理及其他领域 在实际应用中，通过相应的微分方程模型解决具体问题，采用数值方法求得方程的近似解，使具体问题迎刃而解。 2 欧拉法和改进的欧拉法 2.1 欧拉法 2.1.1 欧拉法介绍 首先，我们考虑如下的一阶常微分方程初值问题 ???==0 0)() ,('y x y y x f y （2--1） 事实上，对于更复杂的常微分方程组或者高阶常微分方程，只需要将x 看做向量，（2--1）就成了一个一阶常微分方程组，而高阶常微分方程也可以通过降阶化成一个一阶常微分方程组。 欧拉方法是解常微分方程初值问题最简单最古老的一种数值方法，其基本思路就是把（2--1）中的导数项'y 用差商逼近，从而将一个微分方程转化为一个代数方程，以便求解。 设在[]b a ,中取等距节点h ，因为在节点n x 点上，由（2--1）可得：
各种类型的微分方程及其相应解法教程文件
各种类型的微分方程及其相应解法 专业班级：交土01班 姓名：高云 学号：1201110102 微分方程的类型有很多种，解题时先判断微分方程是哪种类型，可以帮助我们更快解题，所以我们有必要归纳整理一下各类型（主要是一阶和二阶）的微分方程及其相应解法。 一、一阶微分方程的解法 1.可分离变量的方程 dx x f dy y g )()(=,或)()(y g x f dx dy = 其特点是可以把变量x 和y 只分别在等式的两边，解法关键是把变量分离后两边积分。 例1.求微分方程ydy dx y xydy dx +=+2的通解. 解 先合并dx 及dy 的各项,得dx y dy x y )1()1(2-=- 设,01,012≠-≠-x y 分离变量得 dx x dy y y 1112-=- 两端积分??-=-dx x dy y y 1112得 ||ln |1|ln |1|ln 2 112C x y +-=- 于是 2212)1(1-±=-x C y 记,21C C ±=则得到题设方程的通解 .)1(122-=-x C y 2.齐次方程 （1）)(x y f dx dy = （2） )(c by ax f dx dy ++=（a ，b 均不等于0） 例2求解微分方程.2222xy y dy y xy x dx -=+- 解 原方程变形为=+--=2222y xy x xy y dx dy ,1222?? ? ??+--??? ??x y x y x y x y 令,x y u =则,dx du x u dx dy +=方程化为,1222u u u u dx du x u +--=+ 分离变量得?? ????-+--??? ??--112212121u u u u ,x dx du = 两边积分得 ,ln ln ln 2 1)2ln(23)1ln(C x u u u +=----
偏微分方程数值解法答案
1. 课本2p 有证明 2. 课本812,p p 有说明 3. 课本1520,p p 有说明 4. Rit2法，设n u 是u 的n 维子空间，12,...n ???是n u 的一组基底，n u 中的任一元素n u 可 表为1n n i i i u c ?==∑ ，则,11 11()(,)(,)(,)(,)22j n n n n n n i j i j j i j j J u a u u f u a c c c f ???=== -=-∑∑是12,...n c c c 的二次函数，(,)(,)i j j i a a ????=，令 () 0n j J u c ?=?，从而得到12,...n c c c 满足1 (,)(,),1,2...n i j i j i a c f j n ???===∑，通过解线性方程组，求的i c ，代入1 n n i i i u c ?==∑， 从而得到近似解n u 的过程称为Rit2法 简而言之，Rit2法：为得到偏微分方程的有穷维解，构造了一个近似解，1 n n i i i u c ?== ∑， 利用,11 11()(,)(,)(,)(,)22j n n n n n n i j i j j i j j J u a u u f u a c c c f ???===-=-∑∑确定i c ，求得近似解n u 的过程 Galerkin 法：为求得1 n n i i i u c ? == ∑形式的近似解，在系数i c 使n u 关于n V u ∈，满足(,)(,) n a u V f V =，对任 意 n V u ∈或（取 ,1j V j n ?=≤≤） 1 (,)(,),1,2...n i j i j i a c f j n ???===∑的情况下确定i c ，从而得到近似解1 n n i i i u c ?==∑的过程称 Galerkin 法为 Rit2-Galerkin 法方程： 1 (,)(,)n i j i j i a c f ???==∑ 5. 有限元法：将偏微分方程转化为变分形式，选定单元的形状，对求解域作剖分，进而构 造基函数或单元形状函数，形成有限元空间，将偏微分方程转化成了有限元方程，利用 有效的有限元方程的解法，给出偏微分方程近似解的过程称为有限元法。 6. 解：对求解区间进行网格剖分，节点01......i n a x x x x b =<<<<=得到相邻节点1,i i x x -
微分方程的分类及其数值解法
微分方程的分类及其数值解法 微分方程的分类： 含有未知函数的导数，如dy/dx=2x 、ds/dt=0.4都是微分方程。 一般的凡是表示未知函数、未知函数的导数与自变量之间的关系的方程，叫做微分方程。未知函数是一元函数的，叫常微分方程；未知函数是多元函数的叫做偏微分方程。微分方程有时也简称方程。 一、常微分方程的数值解法： 1、Euler 法： 00d (,), (1.1)d (), (1.2) y f x y x y x y ?=???=? 001 (),(,),0,1,,1n n n n y y x y y hf x y n N +=??=+=-? (1.4) 其中0,n b a x x nh h N -=+=. 用(1.4)求解(1.1)的方法称为Euler 方法。 后退Euler 公式???+==+++),,(),(111 00n n n n y x hf y y x y y 梯形方法公式 )].,(),([2 111+++++=n n n n n n y x f y x f h y y 改进的Euler 方法11(,),(,),1().2p n n n c n n p n p c y y hf x y y y hf x y y y y ++?=+??=+???=+??? 2、Runge-Kutta 方法： p 阶方法 ： 1()O h -=?总体截断误差局部截断误差 二阶Runge-Kutta 方法 ??? ????++==++=+),,(),,(,2212 1211hk y h x f k y x f k k h k h y y n n n n n n
常微分方程数值解法
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 。
各种类型的微分方程及其相应解法
各种类型的微分方程及其相应解法 专业班级：交土01班 姓名：高云 学号：1201110102 微分方程的类型有很多种，解题时先判断微分方程是哪种类型，可以帮助我们更快解题，所以我们有必要归纳整理一下各类型（主要是一阶和二阶）的微分方程及其相应解法。 一、一阶微分方程的解法 1.可分离变量的方程 dx x f dy y g )()(=,或)()(y g x f dx dy = 其特点是可以把变量x 和y 只分别在等式的两边，解法关键是把变量分离后两边积分。 例1.求微分方程ydy dx y xydy dx +=+2的通解. 解 先合并dx 及dy 的各项,得dx y dy x y )1()1(2-=- 设,01,012≠-≠-x y 分离变量得 dx x dy y y 1112-=- 两端积分??-=-dx x dy y y 1112得 ||ln |1|ln |1|ln 2 112C x y +-=- 于是 2212)1(1-±=-x C y 记,21C C ±=则得到题设方程的通解 .)1(122-=-x C y 2.齐次方程 （1）)(x y f dx dy = （2） )(c by ax f dx dy ++=（a ，b 均不等于0） 例2求解微分方程.2222xy y dy y xy x dx -=+- 解 原方程变形为=+--=2222y xy x xy y dx dy ,1222?? ? ??+--??? ??x y x y x y x y 令,x y u =则,dx du x u dx dy +=方程化为,1222u u u u dx du x u +--=+ 分离变量得?? ????-+--??? ??--112212121u u u u ,x dx du = 两边积分得 ,ln ln ln 2 1)2ln(23)1ln(C x u u u +=----
常微分方程初值问题数值解法
常微分方程初值问题数值解法 朱欲辉 （浙江海洋学院数理信息学院, 浙江舟山316004） [摘要]：在常微分方程的课程中讨论的都是对一些典型方程求解析解的方法.然而在生产实 际和科学研究中所遇到的问题往往很复杂, 在很多情况下都不可能给出解的解析表达式. 本篇文章详细介绍了常微分方程初值问题的一些数值方法, 导出了若干种数值方法, 如Euler法、改进的Euler法、Runge－Kutta法以及线性多步法中的Adams显隐式公式和预测校正 公式, 并且对其稳定性及收敛性作了理论分析. 最后给出了数值例子, 分别用不同的方法计算出近似解, 从得出的结果对比各种方法的优缺点. [关键词]：常微分方程;初值问题; 数值方法; 收敛性; 稳定性; 误差估计 Numerical Method for Initial-Value Problems Zhu Yuhui (School of Mathematics, Physics, and Information Science, Zhejiang Ocean University, Zhoushan, Zhejiang 316004) [Abstract]:In the course about ordinary differential equations, the methods for analytic solutions of some typical equations are often discussed. However, in scientific research, the problems are very complex and the analytic solutions about these problems can’t be e xpressed explicitly. In this paper, some numerical methods for the initial-value problems are introduced. these methods include Euler method, improved Euler method, Runge-Kutta method and some linear multistep method (e.g. Adams formula and predicted-corrected formula). The stability and convergence about the methods are presented. Some numerical examples are give to demonstrate the effectiveness and accuracy of theoretical analysis. [Keywords]:Ordinary differential equation; Initial-value problem; Numerical method; Convergence; Stability;Error estimate
常微分方程数值解法
第八章 常微分方程数值解法 考核知识点： 欧拉法，改进欧拉法，龙格-库塔法，单步法的收敛性与稳定性。 考核要求： 1． 解欧拉法，改进欧拉法的基本思想；熟练掌握用欧拉法，改进欧拉法、求微 分方程近似解的方法。 2． 了解龙格-库塔法的基本思想；掌握用龙格-库塔法求微分方程近似解的方 法。 3． 了解单步法的收敛性、稳定性与绝对稳定性。 例1 用欧拉法，预估——校正法求一阶微分方程初值问题 ? ??=-='1)0(y y x y ，在0=x （0，1）0.2近似解 解 （1）用1.0=h 欧拉法计算公式 n n n n n n x y y x y y 1.09.0)(1.01+=-+=+，1.0=n 计算得 9.01=y 82.01.01.09.09.02=?+?=y （2）用预估——校正法计算公式 1,0)(05.01.09.0)0(111)0(1=???-+-+=+=++++n y x y x y y x y y n n n n n n n n n 计算得 91.01=y ，83805.02=y 例2 已知一阶初值问题 ???=-='1 )0(5y y y 求使欧拉法绝对稳定的步长n 值。 解 由欧拉法公式 n n n n y h y h y y )51(51-=-=+ n n y h y ~)51(~1-=+
相减得01)51()51(e h e h e n n n -==-=-Λ 当 151≤-h 时，4.00≤ 第六章 二阶线性偏微分方程的分类与小结 一 两个自变量的二阶线性方程 1 方程变换与特征方程 两个自变量的二阶线性偏微分方程总表示成 f cu u b u b u a u a u a y x yy xy xx =+++++212212112 ① 它关于未知函数u 及其一、二阶偏导数都是线性的，其中f u c b b a a a ，，，，，，，21221211都是自变量y x ,的已知函数，假设它们的一阶偏 导数在某平面区域D 内都连续，而且221211a a a ，，不全为0 。 设),(000y x M 是D 内给定的一点，考虑在0M 的领域内对方程进行简化。取自变量变换 ),(y x ξξ=,),(y x ηη= 其中它们具有二连续偏导数，而且在0M 处的雅可比行列式。 = ??),(),(y x ηξy x y x ηηξξ =x y y x ηξηξ- 根据隐函数存在定理，在0M 领域内存在逆变换, ),(ηξx x =,),(ηξy y = 因为 x x x u u u ηξξξ+=，y y y u u u ηξξξ+= xx xx x x x x xx u u u u u u ηξηηξξηξηηξηξξ++++=222 yy yy y y y y yy u u u u u u ηξηηξξηξηηξηξξ++++=222 xy xy y x x y y x x x xy u u u u u u ηξηηηξηξξξηξηηξηξξ+++++=)( 将代入①使其变为 F Cu u B u B u A u A u A =+++++ηξηηξηξξ212212112 经过变换后，方程的阶数不会升高，由变换的可逆性，方程的阶数也不会降低，所以221211,,A A A 不全为0。并可验证 222112122211212))((x y y x a a a A A A ηξηξ--=- 这表明，在可逆变换下2 22112 12A A A -与22112 12 a a a -保持相同的正负号。 定理 在0M 的领域内，不为常数的函数),(y x ?是偏微分方程022*******=++y y x x a a a ????之解的充分必要条件是： C y x ≡),(?是常微分方程的 0)(2)(22212211=++dx a dxdy a dy a 通解。 2 方程的类型及其标准形式 根据以上结论简化方程的问题归结为寻求其特征曲线。为此将特征方程分解成两个方程： 11 22 11 2 12 12 a a a a a dx dy -+=，11 22 11 2 12 12 a a a a a dz dy --= (1) 若在0M 的邻域内022112 12>-a a a 时，方程可以化为 包括基本概念，差分格式的构造、截断误差和稳定性，这些内容是贯穿整个教材的主线。解答问题关键在过程，能够显示出你已经掌握了书上的内容，知道了解题方法。这次考试题目的类型：20分的选择题，主要是基本概念的理解，后面有五个大题，包括差分格式的构造、截断误差和稳定性。 习题一 1． 略 2． y y x f -=),(，梯形公式：n n n n n n y h h y y y h y y )121(),(2111+-+=+- =+++，所以0122)1(01])121[()121()121(y h h y h h y h h y h h n h h n n n +--+--+-+=+-+==+-+= ，当0→h 时， x n e y -→。 同理可以证明预报-校正法收敛到微分方程的解. 3． 局部截断误差的推导同欧拉公式; 整体截断误差: ? ++++++-++≤1 ),())(,(11111n n x x n n n n n n n dx y x f x y x f R εε 11)(++-++≤n n n y x y Lh R ε，这里R R n ≤ 而111)(+++-=n n n y x y ε，所以 R Lh n n += -+εε1)1(，不妨设1 微分方程数值解―― 第二章 习题 1． 设)('x f 为)(x f 的一阶广义导数，试用类似办法定义)(x f 的k 阶广义导数) () (x f k （ ,2,1=k ）。 解：对一维情形，函数的广义导数是通过分部积分来定义的。 我们知，)(x f 的一阶广义导数位)(x g ，如果满足 dx x x f dx x x g b a b a )()()()('?? -=?? 类似的，)(x f 的k 阶广义导数为)()() (x f x g k =，如果有 dx x x f dx x x g b a k k b a )()()1()()()(?? -=?? 2． 试建立与边值问题 ?????====<<=+=) 2.1(0)()(,0)()() 1.1(,''44b u b u a u a u b x a f u dx u d Lu 等价的变分问题。 证明： 设}0)()(,0)()(),(|{' '2====∈=b v a v b v a v I H v v V 对方程)1.1(两边同乘以v ，再关于x 在),(b a 上积分)(V v ∈，得 ??=+b a b a fvdx vdx u dx u d )(44 其中 dx dx dv dx u d dx dx dv dx u d dx u d v dx u d d v vdx dx u d b a b a b a b a b a ???? -=-==33 33333344|)( dx dx v d dx u d dx dv dx u d dx u d d dx dv b a b a b a ??+-=-=22222222|)( dx dx v d dx u d b a ? = 2 222 (*) 记dx uv dx v d dx u d v u a b a ?+=)(),(2 222，?=b a fvdx v f ),(。于是我们得到以下等价变分问题的提法： 实验4 常微分方程的数值解 【实验目的】 1．掌握用MATLAB软件求微分方程初值问题数值解的方法； 2．通过实例用微分方程模型解决简化的实际问题； 3．了解欧拉方法和龙格-库塔方法的基本思想和计算公式，及稳定性等概念。 【实验内容】 题3 小型火箭初始重量为1400kg，其中包括1080kg燃料。火箭竖直向上发射时燃料燃烧率为18kg/s，由此产生32000N的推力，火箭引擎在燃料用尽时关闭。设火箭上升时空气阻力正比于速度的平方，比例系数为0.4kg/m，求引擎关闭瞬间火箭的高度、速度、加速度，及火箭到达最高点的时的高度和加速度，并画出高度、速度、加速度随时间变化的图形。 模型及其求解 火箭在上升的过程可分为两个阶段，在全过程中假设重力加速度始终保持不变，g=9.8m/s2。 在第一个过程中，火箭通过燃烧燃料产生向上的推力，同时它还受到自身重力（包括自重和该时刻剩余燃料的重量）以及与速度平方成正比的空气阻力的作用，根据牛顿第二定律，三个力的合力产生加速度，方向竖直向上。因此有如下二式： a=dv/dt=(F-mg-0.4v2)/m=(32000-0.4v2)/(1400-18t)-9.8 dh/dt=v 又知初始时刻t=0，v=0，h=0。记x(1)=h，x(2)=v，根据MATLAB 可以求出0到60秒内火箭的速度、高度、加速度随时间的变化情况。程序如下： function [ dx ] = rocket( t,x ) a=[(32000-0.4*x(2)^2)/(1400-18*t)]-9.8; dx=[x(2);a]; end ts=0:1:60; x0=[0,0]; [t,x]=ode45(@rocket,ts,x0); h=x(:,1); v=x(:,2); a=[(32000-0.4*(v.^2))./(1400-18*t)]-9.8; [t,h,v,a]; 数据如下： t h v a 0 0 0 13.06 1.00 6.57 13.19 13.30 2.00 26.44 26.58 1 3.45 3.00 59.76 40.06 13.50 4.00 106.57 53.54 13.43 5.00 16 6.79 66.89 13.26 6.00 240.27 80.02 12.99 7.00 326.72 92.83 12.61 8.00 425.79 105.22 12.15 9.00 536.99 117.11 11.62 10.00 659.80 128.43 11.02 11.00 793.63 139.14 10.38 12.00 937.85 149.18 9.71 13.00 1091.79 158.55 9.02 14.00 1254.71 167.23 8.33 15.00 1425.93 175.22 7.65 16.00 1604.83 182.55 6.99 17.00 1790.78 189.22 6.36 18.00 1983.13 195.27 5.76 19.00 2181.24 200.75 5.21 20.00 2384.47 205.70 4.69 21.00 2592.36 210.18 4.22 22.00 2804.52 214.19 3.79 23.00 3020.56 217.79 3.41 24.00 3240.08 221.01 3.07 25.00 3462.65 223.92 2.77 26.00 3687.88 226.56 2.50 27.00 3915.58 228.97 2.27 《微分方程数值解法》 【摘要】自然界与工程技术中得很多现象,可以归结为微分方程定解问题。其中,常微分方程求解就是微分方程得重要基础内容。但就是,对于许多得微分方程,往往很难得到甚至不存在精确得解析表达式,这时候,数值解提供了一个很好得解决思路。,针对于此,本文对常微分方程数值解法进行了简单研究,主要讨论了一些常用得数值解法,如欧拉法、改进得欧拉法、Runge—Kutta方法、Adams预估校正法以及勒让德谱方法等,通过具体得算例,结合MA TLAB求解画图,初步给出了一般常微分方程数值解法得求解过程。同时,通过对各种方法得误差分析,让大家对各种方法得特点与适用范围有一个直观得感受。 【关键词】常微分方程数值解法MA TLAB 误差分析 引言 在我国高校,《微分方程数值解法》作为对数学基础知识要求较高且应用非常广泛得一门课程,不仅在数学专业,其她得理工科专业得本科及研究生教育中开设这门课程.近四十年来,《微分方程数值解法》不论在理论上还就是在方法上都获得了很大得发展.同时,由于微分方程就是描述物理、化学与生物现象得数学模型基础,且它得一些最新应用已经扩展到经济、金融预测、图像处理及其她领域在实际应用中,通过相应得微分方程模型解决具体问题,采用数值方法求得方程得近似解,使具体问题迎刃而解。 2 欧拉法与改进得欧拉法 2、1 欧拉法 2、1、1 欧拉法介绍 首先,我们考虑如下得一阶常微分方程初值问题 (21) 事实上,对于更复杂得常微分方程组或者高阶常微分方程,只需要将瞧做向量,(21)就成了一个一阶常微分方程组,而高阶常微分方程也可以通过降阶化成一个一阶常微分方程组。 欧拉方法就是解常微分方程初值问题最简单最古老得一种数值方法,其基本思路就就是把(21)中得导数项用差商逼近,从而将一个微分方程转化为一个代数方程,以便求解。 设在中取等距节点,因为在节点点上,由(21)可得: , (22) 又由差商得定义可得: (23) 所以有 (24) 用得近似值代入(24),则有计算得欧拉公式 (25) 2、1、2欧拉法误差分析 第八章 常微分方程的数值解法 一．内容要点 考虑一阶常微分方程初值问题:?????==0 0)() ,(y x y y x f dx dy 微分方程的数值解：设微分方程的解y (x )的存在区间是[a,b ]，在[a,b ]内取一系列节 点a= x 0< x 1<…< x n =b ,其中h k =x k+1-x k ；(一般采用等距节点，h=(b-a)/n 称为步长)。在每个节点x k 求解函数y(x)的近似值:y k ≈y(x k )，这样y 0 , y 1 ,...,y n 称为微分方程的数值解。 用数值方法，求得f(x k )的近似值y k ，再用插值或拟合方法就求得y(x)的近似函数。 (一)常微分方程处置问题解得存在唯一性定理 对于常微分方程初值问题:?????==0 0)() ,(y x y y x f dx dy 如果： (1) 在B y y A x x 00≤-≤≤,的矩形内),(y x f 是一个二元连续函数。 (2) ),(y x f 对于y 满足利普希茨条件，即 2121y y L y x f y x f -≤-),(),(则在C x x 0≤≤上方程?????==0 0)() ,(y x y y x f dx dy 的解存在且唯一，这里C=min((A-x 0),x 0+B/L),L 是利普希茨常数。 定义：任何一个一步方法可以写为),,(h y x h y y k k k 1k Φ+=+，其中),,(h y x k k Φ称为算法的增量函数。 收敛性定理：若一步方法满足： (1)是p 解的. (2) 增量函数),,(h y x k k Φ对于y 满足利普希茨条件. (3) 初始值y 0是精确的。则),()()(p h O x y kh y =-kh =x -x 0,也就是有 0x y y lim k x x kh 0h 0 =--=→)( (一)、主要算法 1．局部截断误差 局部截断误差：当y(x k )是精确解时，由y(x k )按照数值方法计算出来的1~ +k y 的误差y (x k+1)- 1~ +k y 称为局部截断误差。 注意：y k+1和1~ +k y 的区别。因而局部截断误差与误差e k +1=y (x k +1) -y k +1不同。 如果局部截断误差是O (h p+1),我们就说该数值方法具有p 阶精度。 《微分方程数值解》教学大纲 一、课程的基本信息 课程名称：《微分方程数值解》 英文名称：Numerical solution for differential equaiton 课程性质：专业方向选修课 课程编号：1623313002 周学时：3学时 总学时：48学时（理论40+实验8） 学分：3学分 适用专业： 适用于信息与计算科学专业 预备知识：数值计算、常微分方程、数值逼近、数理方程 课程教材： 李立康主编，《微分方程数值解法》，复旦大学出版社出版、1999年 参考书目： [1] 戴嘉尊主编，《微分方程数值解法》，东南大学出版社、2008年． [2] 李荣华主编，《微分方程数值解法》(第四版)，高等教育出版社、2009年. 考核方式：考查 制定时间：2013年10月制定 二、课程的目的与任务 《微分方程数值解》是高等院校信息与计算科学专业的专业选修课之一。本课程主要内容为常微分方程和偏微分方程的数值求解问题，包括各种差分方法，有限元方法等的基本理论。通过微分方程数值解的教学，使学生了解和掌握微分方程数值解这一学科的基本概念、理论，培养学生的理论思维能力，为从事信息与计算科学学科的教学和研究打下一定的理论基础。 通过本课程的学习，学生应熟练掌握常微分方程和偏微分方程的常用数值求解方法和分析手段，从能力方面，应使学生初步认识如何从实际问题出发，建立微分方程数学模型，将连续问题离散化，由微分方程转化为差分方程，利用计算机实现数值方法求解一个微分方程的定解问题，并对结果给以几何解释。从教学方法上，着重体现思维方式，注重解决实际问题的方法以及利用计算机进行科学计算的能力 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 -=<<<<=L （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 +=+=-L 方法： （i.4） 从(i.1b) 给出的初始值0u 出发，由上式可以依次算出1,,N t t L 上的差分解1,,N u u L 。 第二章 二阶线性偏微分方程的分类 1．把下列方程化为标准形式： （1）02=+++++u cu bu au au au y x yy xy xx 解：因为 02 22112 12=?-=-a a a a a a 所以该方程是抛物型方程，其特征方程为 12 2 =-± =a a a a dx dy 。 它只有一族实的特征线 c x y =- 在这种情况下，我们设x y -=ξ，x =η（或令y =η，总之，此处η是与ξ无关的任一函数，当然宜取最简单的函数形式x =η或y =η）。 方法一：用抛物型方程的标准形式 ][12122 F Cu u B u B A +++- =ηξηηη 先算出： ? ??? ? ? ?? ? ? ?-====?+?+?+?+?=++++=?+-+?+?+?=++++==?+?+=++=b c C b c b a a a b b a a a B c b a a a b b a a a B a a a a a a a A y x yy xy xx y x yy xy xx y y x x 0F ,1010020 2 1)1(0020 2 002 2212212112 2122121112 221221122ηηηηηξξξξξηηηη ∴])[(1 u bu u c b a u +++--=ηξηη 即 01=+ + -+ u a u a b u a b c u ηξηη 方法二：应用特征方程，作自变量变换，求出 ??? ??=+-=+-=+--==+-= ,2 ,ξξηξξξηηξηξξηηηξξηξξξηξu u u u u u u u u u u u u u u u u u yy xy xx y x 代入原方程得，0)(=++-+u bu u b c au ηξξη 常微分方程数值解法 【作用】微分方程建模是数学建模的重要方法，因为许多实际问题的数学描述将导致求解微分方程的定解问题。把形形色色的实际问题化成微分方程的定解问题，大体上可以按以下几步： 1. 根据实际要求确定要研究的量(自变量、未知函数、必要的参数等)并确定坐标系。 2. 找出这些量所满足的基本规律(物理的、几何的、化学的或生物学的等等)。 3. 运用这些规律列出方程和定解条件。基本模型 1. 发射卫星为什么用三级火箭 2. 人口模型 3. 战争模型 4. 放射性废料的处理通常需要求出方程的解来说明实际现象，并加以检验。如果能得到解析形式的解固然是便于分析和应用的，但是我们知道，只有线性常系数微分方程，并且自由项是某些特殊类型的函数时，才可以得到这样的解，而绝大多数变系数方程、非线性方程都是所谓“解不出来” 的于是对于用微分方程解决实际问题来说，数值解法就是一个十分重要的手段。 1. 改进Euler 法： 2. 龙格—库塔( Runge—Kutta )方法： 【源程序】 1. 改进Euler 法： function [x,y]=eulerpro(fun,x0,x1,y0,n);%fun 为函数，(xO, x1)为x 区间，yO 为初始值，n 为子 区间个数 if nargin<5,n=5O;end h=(x1-xO)/n; x(1)=xO;y(1)=yO; for i=1:n x(i+1)=x(i)+h; y1=y(i)+h*feval(fun,x(i),y(i)); y2=y(i)+h*feval(fun,x(i+1),y1); y(i+1)=(y1+y2)/2; end 调用command 窗口 f=i nlin e('-2*y+2*x A2+2*x') [x,y]=eulerpro(f,O,,1,1O) 2 x +2x , (0 < x < , y(0) = 1 求解函数y'=-2y+2 2. 龙格—库塔( Runge—Kutta )方法： [t,y]=solver('F',tspan ，y0) 这里solver为ode45, ode23, ode113，输入参数F是用M文件定义的微分方程y'= f (x, y)右端的函数。tspan=[t0,tfinal]是求解区间，y0是初值。 注：ode45和ode23变步长的，采用Runge-Kutta算法。 ode45表示采用四阶-五阶Runge-Kutta算法，它用4阶方法提供候选解，5阶方法控制误差，是一种自适应步长(变步长)的常微分方程数值解法，其整体截断误差为(△ 口人5解 决的是Nonstiff(非刚性)常微分方程。二阶线性偏微分方程的分类与小结
微分方程数值解法答案
微分方程数值解――
常微分方程的数值解
微分方程数值解法
常微分方程数值解法
微分方程数值解教学大纲
(整理)常微分方程数值解法
第二章 二阶线性偏微分方程的分类
常微分方程数值解法