常微分方程初值问题数值解法的比较

常微分方程初值问题数值解法的比较
常微分方程初值问题数值解法的比较

常微分方程初值问题数值解法的比较

数值计算实践—课程设计报告

指导教师:徐红敏2013年1 月17 日

常微分方程边值问题的数值解法

第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 αβ''-=<

常微分方程初值问题数值解法.

常微分方程初值问题数值解法 朱欲辉 (浙江海洋学院数理信息学院, 浙江舟山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. 让学生掌握用Euler 法, Runge-Kutta 法求解常微分方程初值问题. 2. 培养Matlab 编程与上机调试能力. 三、实验环境 每人一台计算机,要求安装Windows XP 操作系统,Microsoft office2003、MATLAB6.5(或7.0). 四、实验内容 1. 取步长h=0.1,0.05,0.01, ,用Euler 法及经典4阶Runge-Kutta 法求解初值 问题 ?? ?=≤≤++-=1 )0() 10(2222'y t t t y y 要求: 1) 画出准确解(准确解22t e y t +=-)的曲线,近似解折线; 2) 把节点0.1和0.5上的精确解与近似解比较,观察误差变化情况. 2. 用 Euler 法,隐式Euler 法和经典4阶R-K 法取不同步长解初值问题 ?? ? ??= ∈-=21 )0(],1,0[,50'y x y y 并画出曲线观察稳定性. 注:题1必须写实验报告 五、算法描述及实验步骤 Euler 法: 输入 000),(,,,),,(y a x x h b a y x f = 输出 Euler 解y 步1 ),,2,1(;m n h n a x h a b m n =?+=-? 步2 对1,,2,1,0-=m n 执行),(1n n n n y x f h y y ?+?+

步3 输出T m y y y y ),,,(21 = 经典4阶R-K 法: 输入 000),(,,,),,(y a x x h b a y x f = 输出 4阶R-K 解y 步1 ),,2,1(;m n h n a x h a b m n =?+=-? 步2 对1,,2,1,0-=m n 执行),(1n n y x f K ?,)5.0,(15.02hK y x f K n n +?+, )5.0,(25.03hK y x f K n n +?+,),(314hK y x f K n n +?+ )22(6 43211K K K K h y y n n ++++?+ 步3 输出T m y y y y ),,,(21 = 六、调试过程及实验结果 >> shiyan6 Y1 = 0.8000 0.6620 0.5776 0.5401 0.5441 0.5853 0.6602 0.7662 0.9009 1.0627 Y2 = 0.8287 0.7103 0.6388 0.6093 0.6179 0.6612 0.7366 0.8419 0.9753 1.1353

一阶常微分方程解法总结

第 一 章 一阶微分方程的解法的小结 ⑴、可分离变量的方程: ①、形如 )()(y g x f dx dy = 当0)(≠y g 时,得到 dx x f y g dy )() (=,两边积分即可得到结果; 当0)(0=ηg 时,则0)(η=x y 也是方程的解。 例1.1、 xy dx dy = 解:当0≠y 时,有 xdx y dy =,两边积分得到)(2ln 2为常数C C x y += 所以)(112 12 C x e C C e C y ±==为非零常数且 0=y 显然是原方程的解; 综上所述,原方程的解为)(12 12 为常数C e C y x = ②、形如0)()()()(=+dy y Q x P dx y N x M 当0)()(≠y N x P 时,可有 dy y N y Q dx x P x M ) () ()()(=,两边积分可得结果; 当0)(0=y N 时,0y y =为原方程的解,当0(0=) x P 时,0x x =为原方程的解。 例1.2、0)1()1(2 2 =-+-dy x y dx y x 解:当0)1)(1(2 2 ≠--y x 时,有 dx x x dy y y 1 122-=-两边积分得到 )0(ln 1ln 1ln 22≠=-+-C C y x ,所以有)0()1)(1(22≠=--C C y x ; 当0)1)(1(2 2 =--y x 时,也是原方程的解; 综上所述,原方程的解为)()1)(1(2 2 为常数C C y x =--。 ⑵可化为变量可分离方程的方程: ①、形如 )(x y g dx dy = 解法:令x y u =,则udx xdu dy +=,代入得到)(u g u dx du x =+为变量可分离方程,得到

常微分方程解题方法总结.doc

常微分方程解题方法总结 来源:文都教育 复习过半, 课本上的知识点相信大部分考生已经学习过一遍 . 接下来, 如何将零散的知 识点有机地结合起来, 而不容易遗忘是大多数考生面临的问题 . 为了加强记忆, 使知识自成 体系,建议将知识点进行分类系统总结 . 著名数学家华罗庚的读书方法值得借鉴, 他强调读 书要“由薄到厚、由厚到薄”,对同学们的复习尤为重要 . 以常微分方程为例, 本部分内容涉及可分离变量、 一阶齐次、 一阶非齐次、 全微分方程、 高阶线性微分方程等内容, 在看完这部分内容会发现要掌握的解题方法太多, 遇到具体的题 目不知该如何下手, 这种情况往往是因为没有很好地总结和归纳解题方法 . 下面以表格的形 式将常微分方程中的解题方法加以总结,一目了然,便于记忆和查询 . 常微分方程 通解公式或解法 ( 名称、形式 ) 当 g( y) 0 时,得到 dy f (x)dx , g( y) 可分离变量的方程 dy f ( x) g( y) 两边积分即可得到结果; dx 当 g( 0 ) 0 时,则 y( x) 0 也是方程的 解 . 解法:令 u y xdu udx ,代入 ,则 dy 齐次微分方程 dy g( y ) x dx x u g (u) 化为可分离变量方程 得到 x du dx 一 阶 线 性 微 分 方 程 P ( x)dx P ( x) dx dy Q(x) y ( e Q( x)dx C )e P( x) y dx

伯努利方程 解法:令 u y1 n,有 du (1 n) y n dy , dy P( x) y Q( x) y n(n≠0,1)代入得到du (1 n) P(x)u (1 n)Q(x) dx dx 求解特征方程:2 pq 三种情况: 二阶常系数齐次线性微分方程 y p x y q x y0 二阶常系数非齐次线性微分方程 y p x y q x y f ( x) (1)两个不等实根:1, 2 通解: y c1 e 1x c2 e 2x (2) 两个相等实根:1 2 通解: y c1 c2 x e x (3) 一对共轭复根:i , 通解: y e x c1 cos x c2 sin x 通解为 y p x y q x y 0 的通解与 y p x y q x y f ( x) 的特解之和. 常见的 f (x) 有两种情况: x ( 1)f ( x)e P m ( x) 若不是特征方程的根,令特解 y Q m ( x)e x;若是特征方程的单根,令特 解 y xQ m ( x)e x;若是特征方程的重根, 令特解 y*x2Q m (x)e x; (2)f (x) e x[ P m ( x) cos x p n ( x)sin x]

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

第七章 常微分方程初值问题的数值解法 --------学习小结 一、本章学习体会 通过本章的学习,我了解了常微分方程初值问题的计算方法,对于解决那些很难求解出解析表达式的,甚至有解析表达式但是解不出具体的值的常微分方程非常有用。在这一章里求解常微分方程的基本思想是将初值问题进行离散化,然后进行迭代求解。在这里将初值问题离散化的方法有三种,分别是差商代替导数的方法、Taylor 级数法和数值积分法。常微分方程初值问题的数值解法的分类有显示方法和隐式方法,或者可以分为单步法和多步法。在这里单步法是指计算第n+1个y 的值时,只用到前一步的值,而多步法则是指计算第n+1个y 的值时,用到了前几步的值。通过对本章的学习,已经能熟练掌握如何用Taylor 级数法去求解单步法中各方法的公式和截断误差,但是对线性多步法的求解理解不怎么透切,特别是计算过程较复杂的推理。 在本章的学习过程中还遇到不少问题,比如本章知识点多,公式多,在做题时容易混淆,其次对几种R-K 公式的理解不够透彻,处理一个实际问题时,不知道选取哪一种公式,通过课本里面几种方法的计算比较得知其误差并不一样,,这个还需要自己在往后的实际应用中多多实践留意并总结。 二、本章知识梳理 常微分方程初值问题的数值解法一般概念 步长h ,取节点0,(0,1,...,)n t t nh n M =+=,且M t T ≤,则初值问题000 '(,),()y f t y t t T y t y =≤≤?? =?的数值解法的一般形式是 1(,,,...,,)0,(0,1,...,)n n n n k F t y y y h n M k ++==-

常微分方程初值问题

常微分方程初值问题 12.1引言 在数学模型中经常出现的常微分方程在科学的许多分支中同样出现,例如工程和经济学。不幸的是却很少出现这些方程可得到表示在封闭的形式的解的情况,所以通常采用数值方法来寻找近似解。如今,这通常可以非常方便的达到高精度和在解析解和数值逼近之间可靠的误差界。在本节我们将关注一阶微分方程(12.1)形式关于实值函数y的实变 量x的结构和数值分析方法,其中和f是一个给定的实值函数的两个变量。为了从解曲线的无限族选择一个特定的积分构成(12.1)的通解,微分方程将与初始条件一起考虑:给定两个实数和,我们寻求一个(12.1)的解决方案,对于有 (12.2) 微分方程(12.1)与初始条件(12.2)被称为一个初值问题。如果你认为任何(12.1),(12.2)形式的初始值问题具有一个唯一解,看看以下例子。 例12.1考虑微分方程,初始条件,其中α是一个固定的实数,α∈(0,1)。 这是一个关于上述想法的简单验证,对于任何非负实数C, 是初值问题在区间[ 0,∞)上的一个解。因此解的存在性是肯定的,但解不一定唯一;事实上,初始值问题的解有一个无限族,当参数。 我们注意到,在与α∈(0,1)相反的情况下,当α≥1,初值问题,具有唯一解y(x)≡0。 例12.1表明函数f必须遵循相对于它的第二个参数的一定的增长性条件,以保证(12.1),(12.2)有唯一解。精确的保证初始值问题(12.1),(12.2)假设f解的存在惟一基于下面的定理。 定理12.1(Picard theorem)假定实值函数是连续的矩形区域D定义 ;当时;且f 满足Lipschitz条件:存在L>0则 。

常微分方程数值解法

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 。

常微分方程初值问题答案

1.(10分)对常微分方程初值问题(0)1(01) dy y dx y x ?=-???=≤≤? 取步长0.1,h = 分别用改进的Euler 法和标准的四阶Runge-Kutta 法作数值计算,写出公式和简要推导过程,并把结果填入表内。 解:(1) 改进的Euler 方法: 代入公式得10.905n n y y +=,即0.905n n y = …2分 (2)标准的四阶Runge-Kutta 方法: 1 12341213 2430.1(22)0.90483756 (0.05)0.95(0.05)0.9525(0.1)0.90475n n n n n n n n n n y y k k k k y k y k y k y k y k y k y k y +?=++++=?? =-?? =-+=-??=-+=-??=-+=-?? 即0.9048375n n y = ……(4分) 2. 对常微分方程初值问题12 (0)1(01) dy y dx y x ?=-???=≤≤? 取步长0.1,h = 分别用改进的Euler 法和标准的四阶Runge-Kutta 法作数值计算,写出公式和推导过程,并把结果填入表内。

解:(1) 改进的Euler 方法: 代入公式得10.95125n n y y +=,即0.95125n n y = ……………….(2分) (2)标准的四阶Runge-Kutta 方法: 1 12341213 2430.1(22)0.9512196 /2(0.05)/20.4875(0.05)/20.4878125(0.1)/20.47622n n n n n n n n n n y y k k k k y k y k y k y k y k y k y k y +?=++++=?? =-?? =-+=-??=-+=-??=-+=-?? 即0.95145314n n y =……(4分) 《数值分析》复习题 一、填空题 1.绝对误差限=末位的一半+单位,相对误差限=绝对误差限/原值*100% 1. 度量一根杆子长250厘米,则其绝对误差限为 ,相对误差限是 。 2. 测量一支铅笔长是16cm , 那么测量的绝对误差限是 ,测量的相对误差限是 。 3. 称量一件商品的质量为50千克,则其绝对误差限为 ,相对误差限是 。 2.利用平方差的方法 4. 在数值计算中,当a _____________

各类微分方程的解法大全

各类微分方程的解法 1.可分离变量的微分方程解法 一般形式:g(y)dy=f(x)dx 直接解得∫g(y)dy=∫f(x)dx 设g(y)及f(x)的原函数依次为G(y)及F(x),则G(y)=F(x)+C为微分方程的隐式通解 2.齐次方程解法 一般形式:dy/dx=φ(y/x) 令u=y/x则y=xu,dy/dx=u+xdu/dx,所以u+xdu/dx=φ(u),即du/[φ(u)-u]=dx/x 两端积分,得∫du/[φ(u)-u]=∫dx/x 最后用y/x代替u,便得所给齐次方程的通解 3.一阶线性微分方程解法 一般形式:dy/dx+P(x)y=Q(x) 先令Q(x)=0则dy/dx+P(x)y=0解得y=Ce- ∫P(x)dx,再令y=u e-∫P(x)dx代入原方程解得u=∫Q(x) e∫P(x)dx dx+C,所以y=e-∫P(x)dx[∫Q(x)e∫P(x)dx dx+C] 即y=Ce-∫P(x)dx +e- ∫P(x)dx∫Q(x)e∫P(x)dx dx为一阶线性微分方程的通解 4.可降阶的高阶微分方程解法 ①y(n)=f(x)型的微分方程 y(n)=f(x) y(n-1)= ∫f(x)dx+C1 y(n-2)= ∫[∫f(x)dx+C1]dx+C2 依次类推,接连积分n次,便得方程y(n)=f(x)的含有n个任意常数的通解②y”=f(x,y’) 型的微分方程 令y’=p则y”=p’,所以p’=f(x,p),再求解得p=φ(x,C1) 即dy/dx=φ(x,C1),所以y=∫φ(x,C1)dx+C2 ③y”=f(y,y’) 型的微分方程

令y ’=p 则y ”=pdp/dy,所以pdp/dy=f(y,p),再求解得p=φ(y,C 1) 即dy/dx=φ(y,C 1),即dy/φ(y,C 1)=dx,所以∫dy/φ(y,C 1)=x+C 2 5.二阶常系数齐次线性微分方程解法 一般形式:y ”+py ’+qy=0,特征方程r 2+pr+q=0 6.二阶常系数非齐次线性微分方程解法 一般形式: y ”+py ’+qy=f(x) 先求y ”+py ’+qy=0的通解y 0(x),再求y ”+py ’+qy=f(x)的一个特解y*(x) 则y(x)=y 0(x)+y*(x)即为微分方程y ”+py ’+qy=f(x)的通解 求y ”+py ’+qy=f(x)特解的方法: ① f(x)=P m (x)e λx 型 令y*=x k Q m (x)e λx [k 按λ不是特征方程的根,是特征方程的单根或特征方程的重根依次取0,1或2]再代入原方程,确定Q m (x)的m+1个系数 ② f(x)=e λx [P l(x)cos ωx+P n (x)sin ωx ]型 令y*=x k e λx [Q m (x)cos ωx+R m (x)sin ωx ][m=max ﹛l,n ﹜,k 按λ+i ω不是特征方程的根或是特征方程的单根依次取0或1]再代入原方程,分别确定Q m (x)和R m (x)的m+1个系数

常微分方程数值解法

第八章 常微分方程数值解法 考核知识点: 欧拉法,改进欧拉法,龙格-库塔法,单步法的收敛性与稳定性。 考核要求: 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≤

常微分方程数值解法

第八章 常微分方程的数值解法 一.内容要点 考虑一阶常微分方程初值问题:?????==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 阶精度。

常微分方程数值解法

常微分方程数值解法 【作用】微分方程建模是数学建模的重要方法,因为许多实际问题的数学描述将导致求解微分方程的定解问题。把形形色色的实际问题化成微分方程的定解问题,大体上可以按以下几步: 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(非刚性)常微分方程。

常微分方程的数值解

实验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

(完整版)专题一(二阶常微分方程解法)

二阶微分方程: 时为非齐次 时为齐次,0)(0)()()()(22≠≡=++x f x f x f y x Q dx dy x P dx y d 二阶常系数齐次线性微分方程及其解法: 2 122,)(2,,(*)0)(1,0(*)r r y y y r r q pr r q p qy y p y 式的两个根、求出的系数; 式中的系数及常数项恰好是,,其中、写出特征方程:求解步骤: 为常数; ,其中?'''=++?=+'+''式的通解:出的不同情况,按下表写、根据(*),321r r 二阶常系数非齐次线性微分方程 型 为常数; 型,为常数 ,]sin )(cos )([)()()(,)(x x P x x P e x f x P e x f q p x f qy y p y n l x m x ωωλλλ+===+'+'' 二阶常系数非齐次线性微分方程的一般形式是 ''+'+=y py qy f x () (1) 其中p q ,是常数。 方程(1)的通解为对应的齐次方程 0=+'+''qy y p y (2) 的通解Y 和方程(1)的一个特解*y 之和。即 *y Y y +=.我们已解决了求二阶常系数齐 次线性方程通解的问题,所以,我们只需讨论求二阶常系数非齐次线性微分方程的特解* y 的方法。 下面我们只介绍当方程(1)中的)(x f 为如下两种常见形式时求其特解*y 的方法。 一、 f x e P x x m ()()=?λ型 由于方程(1)右端函数f x ()是指数函数e x λ?与m 次多项式P x m ()的乘积,而指数

函数与多项式的乘积的导数仍是这类函数,因此,我们推测: 方程(1)的特解应为 y e Q x x *?=λ()( Q x ()是某个次数待定的多项式 ) y e Q x e Q x x x *??'=+'λλλ()() y e Q x Q x Q x x *?"=?+'+''λλλ[()()()]22 代入方程(1),得 e Q x p Q x p q Q x e P x x x m λλλλλ???''++'+++≡?[()()()()()]()22 消去e x λ?,得 ''++'+++≡Q x p Q x p q Q x P x m ()()()()()()22λλλ (3) 讨论 01、如果λ不是特征方程 r pr q 20++=的根。 即 02≠++q p λλ 由于P x m ()是一个m 次的多项式,欲使(3)的两端恒等,那未Q x ()必为一个m 次多项式,设为 Q x b x b x b x b m m m m m ()=++++--0111Λ 将之代入(3),比较恒等式两端x 的同次幂的系数,就得到以b b b b m m 01 1,,,,Λ-为未知数的m +1个线性方程的联立方程组,解此方程组可得到这m +1个待定的系数,并得到特解 y e Q x x m *?=λ() 02、如果λ是特征方程 r pr q 20++=的单根。 即 λλ20++=p q ,但 20λ+≠p 欲使(3)式的两端恒等,那么'Q x ()必是一个m 次多项式。 因此,可令 Q x x Q x m ()()=? 并且用同样的方法来确定)(x Q 的系数b b b b m m 0 11,,,,Λ-。 03、如果λ是特征方程 r pr q 20++=的二重根。 即 λλ20++=p q ,且 20λ+=p 。 欲使(3)式的两端恒等,那么''Q x ()必是一个m 次多项式 因此, 可令 Q x x Q x m ()()=?2 并且用同样的方法来确定)(x Q 的系数b b b b m m 011,,,,Λ-。

(整理)常微分方程数值解法

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 引言 常微分方程是现代数学的一个重要分支,是人们解决各种实际问题的有效工具。微分方程的理论和方法从17世纪末开始发展起来,很快成了研究自然现象的强有力工具,在17到18世纪,在力学、天文、科学技术、物理中,就已借助微分方程取得了巨大的成就。1864年Leverrer根据这个方程预见了海王星的存在,并确定出海王星在天空中的位置。现在,常微分方程在许多方面获得了日新月异的应用。这些应用也为常微分方程的进一步发展提供了新的问题,促使人们对微分方程进行更深入的研究,以便适应科学技术飞速发展的需要。 研究常微分方程常用数值解是数学工作者的一项基本的且重要的工作。在国内外众多数学家的不懈努力,使此学科基本上形成了一套完美的体系。微分方程的首要问题是如何求一个给定方程的通解或特解。到目前为止,人们已经对许多微分方程得出了求解的一般方法。由于在生产实际和科学研究中所遇到的微分方程问题比较复杂,使这些问题的解即使能求出解析表达式,也往往因计算量太大而难于求出,而对于一些典型的微分方程则可以运用基本方法求出其解析解,并可以根据初值问题的条件把其中的任意常数确定下来。 由于求通解存在许多困难,人们就开始研究带某种定解条件的特解。首先是Cauchy对微分方程初始解的存在惟一性进行了研究。目前解的存在惟一性、延拓性、大范围的存在性以及解对初始解和参数的延续性和可微性等理论问题都已发展成熟。与此同时,人们开始采取各种近似方法来求微分方程的特解,例如求微分方程数值解的Euler折线法、Runge-Kutta法等,可以求得若干个点上微分方程的近似解。最后,由于当代高科技的发展为数学的广泛应用和深入研究提供了更好的手段。用计算机结合Matlab软件求方程的精确解、近似解,对解的性态进行图示和定性、稳定性研究都十分方便有效。 本章先介绍常微分的一般概念、导出微分方程的一些典型例子及求解微分方程的思路分析。从而得到常微分方程的常用数值解法。

相关文档
最新文档