课堂记录:常微分方程的数值解法
- 格式:docx
- 大小:44.30 KB
- 文档页数: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) 该格式是显式,也是二阶方法。
常微分方程的数值解法在自然科学的许多领域中,都会遇到常微分方程的求解问题。
然而,我们知道,只有少数十分简单的微分方程能够用初等方法求得它们的解,多数情形只能利用近似方法求解。
在常微分方程课中已经讲过的级数解法,逐步逼近法等就是近似解法。
这些方法可以给出解的近似表达式,通常称为近似解析方法。
还有一类近似方法称为数值方法,它可以给出解在一些离散点上的近似值。
利用计算机解微分方程主要使用数值方法。
我们考虑一阶常微分方程初值问题⎪⎩⎪⎨⎧==00)(),(yx y y x f dx dy在区间[a, b]上的解,其中f (x, y )为x, y 的已知函数,y 0为给定的初始值,将上述问题的精确解记为y (x )。
数值方法的基本思想是:在解的存在区间上取n + 1个节点b x x x x a n =<<<<= 210这里差i i i x x h -=+1,i = 0,1, …, n 称为由x i 到x i +1的步长。
这些h i 可以不相等,但一般取成相等的,这时na b h -=。
在这些节点上采用离散化方法,(通常用数值积分、微分。
泰勒展开等)将上述初值问题化成关于离散变量的相应问题。
把这个相应问题的解y n 作为y (x n )的近似值。
这样求得的y n 就是上述初值问题在节点x n 上的数值解。
一般说来,不同的离散化导致不同的方法。
§1 欧拉法与改进欧拉法 1.欧拉法1.对常微分方程初始问题(9.2))((9.1) ),(00⎪⎩⎪⎨⎧==y x y y x f dx dy用数值方法求解时,我们总是认为(9.1)、(9.2)的解存在且唯一。
欧拉法是解初值问题的最简单的数值方法。
从(9.2)式由于y (x 0) = y 0已给定,因而可以算出),()('000y x f x y =设x 1 = h 充分小,则近似地有:),()(')()(00001y x f x y hx y x y =≈-(9.3)记 ,n ,,i x y y i i 10 )(== 从而我们可以取),(0001y x hf y y ==作为y (x 1)的近似值。
求常微分方程的数值解一、背景介绍常微分方程(Ordinary Differential Equation,ODE)是描述自然界中变化的数学模型。
常微分方程的解析解往往难以求得,因此需要寻找数值解来近似地描述其行为。
求解常微分方程的数值方法主要有欧拉法、改进欧拉法、龙格-库塔法等。
二、数值方法1. 欧拉法欧拉法是最简单的求解常微分方程的数值方法之一。
它基于导数的定义,将微分方程转化为差分方程,通过迭代计算得到近似解。
欧拉法的公式如下:$$y_{n+1}=y_n+f(t_n,y_n)\Delta t$$其中,$y_n$表示第$n$个时间步长处的函数值,$f(t_n,y_n)$表示在$(t_n,y_n)$处的导数,$\Delta t$表示时间步长。
欧拉法具有易于实现和理解的优点,但精度较低。
2. 改进欧拉法(Heun方法)改进欧拉法又称Heun方法或两步龙格-库塔方法,是对欧拉法进行了精度上提升后得到的一种方法。
它利用两个斜率来近似函数值,并通过加权平均来计算下一个时间步长处的函数值。
改进欧拉法的公式如下:$$k_1=f(t_n,y_n)$$$$k_2=f(t_n+\Delta t,y_n+k_1\Delta t)$$$$y_{n+1}=y_n+\frac{1}{2}(k_1+k_2)\Delta t$$改进欧拉法比欧拉法精度更高,但计算量也更大。
3. 龙格-库塔法(RK4方法)龙格-库塔法是求解常微分方程中最常用的数值方法之一。
它通过计算多个斜率来近似函数值,并通过加权平均来计算下一个时间步长处的函数值。
RK4方法是龙格-库塔法中最常用的一种方法,其公式如下:$$k_1=f(t_n,y_n)$$$$k_2=f(t_n+\frac{\Delta t}{2},y_n+\frac{k_1\Delta t}{2})$$ $$k_3=f(t_n+\frac{\Delta t}{2},y_n+\frac{k_2\Delta t}{2})$$ $$k_4=f(t_n+\Delta t,y_n+k_3\Delta t)$$$$y_{n+1}=y_n+\frac{1}{6}(k_1+2k_2+2k_3+k_4)\Delta t$$三、数值求解步骤对于给定的常微分方程,可以通过以下步骤求解其数值解:1. 确定初值条件:确定$t=0$时刻的函数值$y(0)$。
常微分方程的数值解算法常微分方程的数值解算法是一种对常微分方程进行数值计算的方法,这可以帮助我们更好地理解和研究自然现象和工程问题。
在本文中,我们将介绍一些常用的数值解算法,探讨它们的优缺点和适用范围。
常微分方程(ODE)是描述自然现象和工程问题的重要数学工具。
然而,对于许多ODE解析解是无法求出的,因此我们需要通过数值方法对其进行求解。
常微分方程可以写作:y' = f(t, y)其中,y是函数,f是给定的函数,表示y随t的变化率。
这个方程可以写成初始值问题(IVP)的形式:y'(t) = f(t,y(t)),y(t0) = y0其中,y(t0)=y0是方程的初始条件。
解决IVP问题的典型方法是数值方法。
欧拉方法欧拉方法是最简单的一阶数值方法。
在欧拉方法中,我们从初始条件开始,并在t = t0到t = tn的时间内,用以下公式逐步递推求解:y n+1 = y n + hf (t n, y n)其中,f(t n,y n)是点(t n,y n)处的导数, h = tn - tn-1是时间间隔。
欧拉方法的优点是简单易懂,容易实现。
然而,它的缺点是在整个时间段上的精度不一致。
程度取决于使用的时间间隔。
改进的欧拉方法如果我们使用欧拉方法中每个时间段的中间点而不是起始点来估计下一个时间点,精度就会有所提高。
这个方法叫做改进的欧拉方法(或Heun方法)。
公式为:y n+1 = y n + h½[f(t n, y n)+f(tn+1, yn + h f (tn, yn))]这是一个二阶方法,精度比欧拉方法高,但计算量也大一些。
对于易受噪声干扰的问题,改进的欧拉方法是个很好的选择。
Runge-Kutta方法Runge-Kutta方法是ODE计算的最常用的二阶和高阶数值方法之一。
这个方法对定义域内的每个点都计算一个导数。
显式四阶Runge-Kutta方法(RK4)是最常用的Runge-Kutta方法之一,并已得到大量实践的验证。
常微分方程的数值解法1. 引言常微分方程是自变量只有一个的微分方程,广泛应用于自然科学、工程技术和社会科学等领域。
由于常微分方程的解析解不易得到或难以求得,数值解法成为解决常微分方程问题的重要手段之一。
本文将介绍几种常用的常微分方程的数值解法。
2. 欧拉方法欧拉方法是最简单的一种数值解法,其具体步骤如下:- 将自变量的区间等分为n个子区间;- 在每个子区间上假设解函数为线性函数,即通过给定的初始条件在每个子区间上构造切线;- 使用切线的斜率(即导数)逼近每个子区间上的解函数,并将其作为下一个子区间的初始条件;- 重复上述过程直至达到所需的精度。
3. 改进的欧拉方法改进的欧拉方法是对欧拉方法的一种改进,主要思想是利用两个切线的斜率的平均值来逼近每个子区间上的解函数。
具体步骤如下: - 将自变量的区间等分为n个子区间;- 在每个子区间上构造两个切线,分别通过给定的初始条件和通过欧拉方法得到的下一个初始条件;- 取两个切线的斜率的平均值,将其作为该子区间上解函数的斜率,并计算下一个子区间的初始条件;- 重复上述过程直至达到所需的精度。
4. 二阶龙格-库塔方法二阶龙格-库塔方法是一种更为精确的数值解法,其基本思想是通过近似计算解函数在每个子区间上的平均斜率。
具体步骤如下: - 将自变量的区间等分为n个子区间;- 在每个子区间上计算解函数的斜率,并以该斜率的平均值近似表示该子区间上解函数的斜率;- 利用该斜率近似值计算下一个子区间的初始条件,并进一步逼近解函数;- 重复上述过程直至达到所需的精度。
5. 龙格-库塔法(四阶)龙格-库塔法是目前常用的数值解法之一,其精度较高。
四阶龙格-库塔法是其中较为常用的一种,其具体步骤如下:- 将自变量的区间等分为n个子区间;- 在每个子区间上进行多次迭代计算,得到该子区间上解函数的近似值;- 利用近似值计算每个子区间上的斜率,并以其加权平均值逼近解函数的斜率;- 计算下一个子区间的初始条件,并进一步逼近解函数;- 重复上述过程直至达到所需的精度。
常微分方程的数值求解在数学中,常微分方程是一类重要的数学模型,通常用来描述物理、化学、生物等自然现象中的变化规律。
对于一些复杂的微分方程,无法通过解析方法进行求解,这时候就需要借助数值方法来近似求解。
本文将介绍常微分方程的数值求解方法及其应用。
一、数值求解方法常微分方程的数值求解方法主要包括欧拉法、改进的欧拉法、龙格-库塔法等。
欧拉法是最简单也是最常用的数值求解方法,其基本思想是根据微分方程的导数近似求解下一个时间步上的解,并通过逐步迭代来得到整个解的数值近似。
改进的欧拉法在欧拉法的基础上做出了一定的修正,提高了数值求解的精度。
而龙格-库塔法则是一种更加精确的数值求解方法,通过考虑多个点的斜率来进行求解,从而减小误差。
二、应用领域常微分方程的数值求解方法在科学研究和工程实践中有着广泛的应用。
在物理学中,通过数值求解微分方程可以模拟天体运动、粒子运动等现象;在生物学领域,可以模拟生物种群的增长和变化规律;在工程领域,可以通过数值求解微分方程来设计控制系统、优化结构等。
三、实例分析以一个简单的一阶常微分方程为例:dy/dx = -y,初始条件为y(0) = 1。
我们可以用欧拉法来进行数值求解。
将时间间隔取为0.1,通过迭代计算可以得到y(1)的近似值为0.367。
而利用改进的欧拉法或者龙格-库塔法可以得到更加精确的数值近似。
这个例子展示了数值方法在解决微分方程问题上的有效性。
四、总结常微分方程是求解自然界中变化规律的重要数学工具,而数值方法则是解决一些难以解析求解的微分方程的有效途径。
通过本文的介绍,读者可以了解常微分方程的数值求解方法及其应用,希望可以对相关领域的研究和实践有所帮助。
至此,关于常微分方程的数值求解的文章正文部分结束。
浙江大学城市学院实验报告课程名称数值计算方法实验项目名称常微分方程初值问题的数值解法 实验成绩指导老师签名日期2015/12/16 一.实验目的和要求1. 用Matlab 软件掌握求微分方程数值解的欧拉方法和龙格-库塔方法; 2. 通过实例学习用微分方程模型解决简化的实际问题;二.实验内容和原理编程题2-1要求写出Matlab 源程序m 文件,并有适当的注释语句;分析应用题2-2,2-3,2-4,2-5要求将问题的分析过程、Matlab 源程序和运行结果和结果的解释、算法的分析写在实验报告上; 2-1 编程编写用向前欧拉公式和改进欧拉公式求微分方程数值解的Matlab 程序,问题如下:在区间[],a b 内(1)N +个等距点处,逼近下列初值问题的解,并对程序的每一句添上注释语句; Euler 法y=eulera,b,n,y0,f,f1,b1改进Euler 法y=eulerproa,b,n,y0,f,f1,b1 2-2 分析应用题假设等分区间数100n =,用欧拉法和改进欧拉法在区间[0,10]t ∈内求解初值问题()()20(0)10y t y t y '=-⎧⎨=⎩并作出解的曲线图形,同时将方程的解析解也画在同一张图上,并作比较,分析这两种方法的精度; 2-3 分析应用题用以下三种不同的方法求下述微分方程的数值解,取10h = 画出解的图形,与精确值比较并进行分析; 1欧拉法; 2改进欧拉法; 3龙格-库塔方法;2-4 分析应用题考虑一个涉及到社会上与众不同的人的繁衍问题模型;假设在时刻t 单位为年,社会上有人口()x t 人,又假设所有与众不同的人与别的与众不同的人结婚后所生后代也是与众不同的人;而固定比例为r 的所有其他的后代也是与众不同的人;如果对所有人来说出生率假定为常数b ,又如果普通的人和与众不同的人的婚配是任意的,则此问题可以用微分方程表示为:其中变量()()()i p t x t x t =表示在时刻t 社会上与众不同的人的比例,()i x t 表示在时刻t 人口中与众不同的人的数量;1假定(0)0.01,0.02p b ==和0.1r =,当步长为1h =年时,求从0t =到50t =解()p t 的近似值,并作出近似解的曲线图形;2精确求出微分方程的解()p t ,并将你当50t =时在分题b 中得到的结果与此时的精确值进行比较; MATLAB 相关函数求微分方程的解析解及其数值的代入dsolve‘egn1’,‘egn2’,‘x ’ subsexpr,{x,y,…},{x1,y1,…}其中‘egn i ’表示第i 个方程,‘x ’表示微分方程中的自变量,默认时自变量为t ; subs 命令中的expr 、x 、y 为符合型表达式,x 、y 分别用数值x1、x2代入; >>symsxyz>>subs'x+y+z',{x,y,z},{1,2,3} ans= 6>>symsx>>subs'x^2',x,2 ans= 4>>s=dsolve‘12Dy y ∧=+’,‘(0)1y =’,‘x ’ ans= >>symsx >>subss,x,2 ans=右端函数(,)f x y 的自动生成f=inline ‘expr ’,’var1’,‘var2’,……其中’expr ’表示函数的表达式,’var1’,‘var2’表示函数表达式中的变量,运行该函数,生成一个新的函数表达式为fvar1,var2,……; >>f=inline'x+3y','x','y' f=Inlinefunction: fx,y=x+3y >>f2,3 ans= 114,5阶龙格-库塔方法求解微分方程数值解t,x=ode45f,ts,x0,options其中f 是由待解方程写成的m 文件名;x0为函数的初值;t,x 分别为输出的自变量和函数值列向量,t的步长是程序根据误差限自动选定的;若ts=t0,t1,t2,…,tf,则输出在自变量指定值,等步长时用ts=t0:k:tf,输出在等分点;options 用于设定误差限可以缺省,缺省时设定为相对误差310-,绝对误差610-,程序为:options=odeset ‘reltol ’,rt,’abstol ’,at,这里rt,at 分别为设定的相对误差和绝对误差;常用选项见下表;选项名 功能 可选值 省缺值 AbsTol 设定绝对误差正数 RelTol 设定相对误差 正数InitialStep 设定初始步长 正数 自动 MaxStep设定步长上界正数MaxOrder 设定ode15s 的最高阶数 1,2,3,4,5 5 Stats 显示计算成本统计 on,off off BDF 设定ode15s 是否用反向差分on,offoff例:在命令窗口执行>>odefun =inline ‘2*y t y -’,‘t ’,‘y ’;>>[],45(,[0,4],1)t y ode odefun =;ans=>>t y ‘o-’,%解函数图形表示>>45(,[0,4],1)ode odefun %不用输出变量,则直接输出图形 >>[],45(,0:4,1)t y ode odefun =;[],t yans=三.操作方法与实验步骤包括实验数据记录和处理2-1编程编写用向前欧拉公式和改进欧拉公式求微分方程数值解的Matlab 程序,问题如下:在区间[],a b 内(1)N +个等距点处,逼近下列初值问题的解,并对程序的每一句添上注释语句; Euler 法y=eulera,b,n,y0,f,f1,b1改进Euler 法y=eulerproa,b,n,y0,f,f1,b1Euler 法y=eulera,b,n,y0,f,f1,b1 y=zeros1,n+1; y1=y0; h=b-a/n; x=a:h:b; fori=1:n; yi+1=yi+hfxi,yi; end plotx,y holdon%求微分方程的精确解 x1=linspacea,b,100; '精确解为' s=dsolvef1,b1,'x' symsxy1=zeros1,100; for i=1:100y1i=subss,x,x1i; endplotx1,y1,'r'title'红色代表精确解'改进Euler 法y=eulerproa,b,n,y0,f,f1,b1 %求微分方程的数值解 y=zeros1,n+1; y1=y0; h=b-a/n; x=a:h:b; fori=1:n; T1=fxi,yi; T2=fxi+1,yi+hT1; yi+1=yi+h/2T1+T2; end plotx,y holdon%求微分方程的精确解 x1=linspacea,b,100; '精确解为' s=dsolvef1,b1,'x' symsxy1=zeros1,100; fori=1:100 y1i=subss,x,x1i; endplotx1,y1,'r'title'红色代表精确解' 2-2分析应用题假设等分区间数100n =,用欧拉法和改进欧拉法在区间[0,10]t ∈内求解初值问题()()20(0)10y t y t y '=-⎧⎨=⎩并作出解的曲线图形,同时将方程的解析解也画在同一张图上,并作比较,分析这两种方法的精度;1向前欧拉法>>euler0,10,100,10,inline'y-20','x','y','Dy=y-20','y0=10' ans= 精确解为 s= 20-10expx ans= +005Columns1through8(2)改进欧拉法>>eulerpro0,10,100,10,inline'y-20','x','y','Dy=y-20','y0=10' ans= 精确解为 s= 20-10expx ans= +005Columns1through8改进欧拉法的精度比向前欧拉法更高; 2-3分析应用题用以下三种不同的方法求下述微分方程的数值解,取10h = 画出解的图形,与精确值比较并进行分析; 1欧拉法; 2改进欧拉法;2-4分析应用题考虑一个涉及到社会上与众不同的人的繁衍问题模型;假设在时刻t 单位为年,社会上有人口()x t 人,又假设所有与众不同的人与别的与众不同的人结婚后所生后代也是与众不同的人;而固定比例为r 的所有其他的后代也是与众不同的人;如果对所有人来说出生率假定为常数b ,又如果普通的人和与众不同的人的婚配是任意的,则此问题可以用微分方程表示为:其中变量()()()i p t x t x t =表示在时刻t 社会上与众不同的人的比例,()i x t 表示在时刻t 人口中与众不同的人的数量;1假定(0)0.01,0.02p b ==和0.1r =,当步长为1h =年时,求从0t =到50t =解()p t 的近似值,并作出近似解的曲线图形;2精确求出微分方程的解()p t ,并将你当50t =时在分题b 中得到的结果与此时的精确值进行比较;1>>euler0,50,50,,inline'','t','p','Dp=','p0= 1' ans= 精确解为 s=1-99/100expx/500 ans=Columns1through82>>dsolve'Dp=','p0=','t' ans=1-99/100expt/500 >>1-99/100exp ans=与欧拉法求得的精确值差0,0001四.实验结果与分析。
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 。
下面我们用数值积分法重新导出 Euler 法以及其它几种方法。
为此,在区间1[,]n n t t +上积分常微分方程(i.1a ),得11()()(,())n n t n n t u t u t f t u t dt ++=+⎰ (i.5)用各种数值积分公式计算(i.5)中的积分,便导致各种不同的差分法。
常微分方程组数值解法一、引言常微分方程组是数学中的一个重要分支,它在物理、工程、生物等领域都有广泛应用。
对于一些复杂的常微分方程组,往往难以通过解析方法求解,这时候数值解法就显得尤为重要。
本文将介绍常微分方程组数值解法的相关内容。
二、数值解法的基本思想1.欧拉法欧拉法是最基础的数值解法之一,它的思想是将时间连续化,将微分方程转化为差分方程。
对于一个一阶常微分方程y'=f(x,y),其欧拉公式为:y_{n+1}=y_n+hf(x_n,y_n)其中h为步长,x_n和y_n为第n个时间点上x和y的取值。
2.改进欧拉法改进欧拉法是对欧拉法的改良,其公式如下:y_{n+1}=y_n+\frac{h}{2}[f(x_n,y_n)+f(x_{n+1},y_n+hf(x_n,y_n))] 3.四阶龙格-库塔方法四阶龙格-库塔方法是目前最常用的数值解法之一。
其公式如下:k_1=f(x_n,y_n)k_2=f(x_n+\frac{h}{2},y_n+\frac{h}{2}k_1)k_3=f(x_n+\frac{h}{2},y_n+\frac{h}{2}k_2)k_4=f(x_n+h,y_n+hk_3)y_{n+1}=y_n+\frac{h}{6}(k_1+2k_2+2k_3+k_4)其中,k_i为中间变量。
三、常微分方程组的数值解法1.欧拉法对于一个二阶常微分方程组:\begin{cases} y'_1=f_1(x,y_1,y_2) \\ y'_2=f_2(x,y_1,y_2)\end{cases}其欧拉公式为:\begin{cases} y_{n+1,1}=y_{n,1}+hf_1(x_n,y_{n,1},y_{n,2}) \\y_{n+1,2}=y_{n,2}+hf_2(x_n,y_{n,1},y_{n,2}) \end{cases}其中,x_n和y_{n,i}(i=1, 2)为第n个时间点上x和y_i的取值。
1常微分方程的数值解法习题课
1.1 常/偏微分方程
1)电磁计算
2)人口增长
1.2 欧拉法实现
欧拉法算法设计
Input :f(x,y),x0,y0
Output: x0,x1+h,x2+2h,…
y0,y1,y2,…
(或二维数组—矩阵W——2行)
第1行存储xi
第2行存储yi
讲到:
3班:例3 证明修正的Euler法具有2阶精度
% y'=y-xy^2,0<x<2
% y(0)=1,
% y(x+h)-y(x)=h(y(x)-x^2)
% y(x+h)=y(x)+h(y(x)-x^2)
h=0.05;
y0=1;
y(1)=y0;
% code math
%y(1)--y(0)
%y(2)--y(h)
%y(3)--y(2h)
xx=linspace(0,2,51); %0:2/50:2
% for i=2:length(xx)
% x = xx(i-1);%important
% y(i)=y(i-1)+h*(y(i-1)-x^2);
% end
for i=1:length(xx)-1
x = xx(i);%important
y(i+1)=y(i)+h*(y(i)-x^2);
end
plot(xx,y,'r.')
1.3 推导欧拉法
请给出欧拉方法的三种推导过程.
),(y x f y ='
1 )导数定义/差商代替微商
h
x y h x y y x f y h )()(lim ),(0-+=='→ )
,()())(,()()
,())(,()()(,),()()(11n n n n n n n n n n n n n n y x hf x y x y x hf x y y x hf x y x hf x y x y nh
a x x x y x f h
x y h x y +====-+==≈-+++
2 )泰勒公式 h
y x f y y h o h x y h x y h x y x f x y x y y x f y h o h x y h x y h x y x y x y n n n n n n n n n n n n n n n ),()(!
3)(!2)(!1))(,()()(),()(!3)(!2)(!1)()()(133213321+=+'''+''++=='+'''+''+'+
=+++
3) 数值积分
1.4 高阶常微分方程化为一阶常微分方程组 高阶微分方程化为一阶常微分方程组:
例1.),(y x f y y ='+''--高阶
降阶:
令)(1x u y =',则方程),(y x f y y ='+''可化为
⎩
⎨⎧='=+')(),()(111x u y y x f x u u 即⎩⎨⎧-='=')(),()(11
1x u y x f u x u y . 例2.
),(2y x f y y y ='-''+'''--高阶
⎪⎩⎪⎨⎧++-='='='⎪⎩⎪⎨⎧++'-='='='='-''+'''),()(2)()(),()(2)()()
,(2122
211112
211y x f x u u u x u u x u y y x f x u u u x u u x u y y x f y y y
1.5 练习
⎩⎨⎧='=')
,,(),,(y x t g x y x t f y --多个未知函数 请给出Euler 方法求解该方程的迭代格式。
⎩⎨⎧+=+=⎩⎨⎧='='++)
,,(),,(),,(),,(11n n n n n n n n n n y x t hg x x y x t hf y y y x t g x y x t f y
请写出解下列微分方程组的Euler 法迭代公式: ),,(),,(212
211y y x h y y y x f y ='=' xn=x0+nh
))(),(,()()())(),(,()()(2122221111x y x y x f h
x y h x y x y x y x f h x y h x y ≈-+≈-+ ⎩⎨⎧+≈++≈+))(),(,()()())(),(,()()(2122221111x y x y x hf x y h x y x y x y x hf x y h x y 求解该问题的Euler 法迭代公式如下:
⎪⎩⎪⎨⎧+=+=++)
,,(),,()(2)(12)(2)1(2)(2)(11)(1)1(1n n n n n n n n n n y y x hf y y y y x hf y y 求解该问题的改进的Euler 法迭代公式如下:
)],,(),,([2
)],,(),,([2
)
,,()
,,()1(2)1(112)(2)(12)(2)1(2)1(2)1(111)(2)(11)(1)1(1)(2)(12)(2)1(2)(2)(11)(1)1(1++++++++++++=++=+=+=n n n n n n n n n n n n n n n n n n n n n n n n n n y y x f y y x f h y y y y x f y y x f h y y y y x hf y y y y x hf y y
练习:设计求解下列微分方程的Euler 迭代格式: ⎩
⎨⎧='='=''.)0(,)0(),,(100y y y y y y x f y 令y z '=,
⎪⎩
⎪⎨⎧==='='.)0(,)0()
,,(100y z y y z y x f z z y Euler:
算法描述:
⎩
⎨⎧='='=''.)0(,)0(),,(100y y y y y y x f y 算法:求解二阶微分方程的Euler 法
Input :
y0 (即y(0))
y10 (即y ’(0)=z(0))
N (即求[0,N*h]区间上函数的近似解) h 步长
Output:
y_k(k=0,1,2,…,N)
步骤:
x_k=0;
For k=0:N-1
x_k=k*h;
y_(k+1)=y_k+h*z_k
z_(k+1)=z_k + h*f(x_k, y_k, z_k)
if y_(k+1) < 1 then
break
end if
End for。