【免费下载】常微分方程的数值解法实验报告

合集下载

常微分方程数值解实验报告

常微分方程数值解实验报告

常微分方程数值解实验报告学院:数学与信息科学专业:信息与计算科学姓名:郑思义学号:201216524课程:常微分方程数值解实验一:常微分方程的数值解法1、分别用Euler 法、改进的Euler 法(预报校正格式)和S —K 法求解初值问题。

(h=0.1)并与真解作比较。

⎩⎨⎧=++-=10(1y')y x y 1.1实验代码:%欧拉法function [x,y]=naeuler(dyfun,xspan,y0,h)%dyfun 是常微分方程,xspan 是x 的取值范围,y0是初值,h 是步长 x=xspan(1):h:xspan(2); y(1)=y0;for n=1:length(x)-1y(n+1)=y(n)+h*feval(dyfun,x(n),y(n)); end%改进的欧拉法function [x,m,y]=naeuler2(dyfun,xspan,y0,h)%dyfun 是常微分方程,xspan 是x 的取值范围,y0是初值,h 是步长。

%返回值x 为x 取值,m 为预报解,y 为校正解 x=xspan(1):h:xspan(2); y(1)=y0;m=zeros(length(x)-1,1); for n=1:length(x)-1 k1=feval(dyfun,x(n),y(n)); y(n+1)=y(n)+h*k1; m(n)=y(n+1);k2=feval(dyfun,x(n+1),y(n+1)); y(n+1)=y(n)+h*(k1+k2)/2; end%四阶S —K 法function [x,y]=rk(dyfun,xspan,y0,h)%dyfun 是常微分方程,xspan 是x 的取值范围,y0是初值,h 是步长。

x=xspan(1):h:xspan(2); y(1)=y0;for n=1:length(x)-1 k1=feval(dyfun,x(n),y(n));k2=feval(dyfun,x(n)+h/2,y(n)+(h*k1)/2); k3=feval(dyfun,x(n)+h/2,y(n)+(h*k2)/2); k4=feval(dyfun,x(n)+h,y(n)+h*k3);y(n+1)=y(n)+(h/6)*(k1+2*k2+2*k3+k4);end%主程序x=[0:0.1:1];y=exp(-x)+x;dyfun=inline('-y+x+1');[x1,y1]=naeuler(dyfun,[0,1],1,0.1);[x2,m,y2]=naeuler2(dyfun,[0,1],1,0.1);[x3,y3]=rk(dyfun,[0,1],1,0.1);plot(x,y,'r',x1,y1,'+',x2,y2,'*',x3,y3,'o');xlabel('x');ylabel('y');legend('y为真解','y1为欧拉解','y2为改进欧拉解','y3为S—K解','Location','NorthWest');1.2实验结果:x 真解y 欧拉解y1 预报值m 校正值y2 S—K解y30.0 1.0000 1.0000 1.0000 1.00000.1 1.0048 1.0000 1.0000 1.0050 1.00480.2 1.0187 1.0100 1.0145 1.0190 1.01870.3 1.0408 1.0290 1.0371 1.0412 1.04080.4 1.0703 1.0561 1.0671 1.0708 1.07030.5 1.1065 1.0905 1.1037 1.1071 1.10650.6 1.1488 1.1314 1.1464 1.1494 1.14880.7 1.1966 1.1783 1.1945 1.1972 1.19660.8 1.2493 1.2305 1.2475 1.2500 1.24930.9 1.3066 1.2874 1.3050 1.3072 1.30661.0 1.3679 1.3487 1.3665 1.3685 1.36792、选取一种理论上收敛但是不稳定的算法对问题1进行计算,并与真解作比较。

实验报告——常微分方程的数值解法

实验报告——常微分方程的数值解法

实验报告实验项目名称常微分方程的数值解法实验室数学实验室所属课程名称微分方程数值解实验类型上机实验实验日期2013年3月11日班级10信息与计算科学学号2010119421姓名叶达伟成绩实验概述:【实验目的及要求】运用不同的数值解法来求解具体问题,并通过具体实例来分析比较各种常微分方程的数值解法的精度,为以后求解一般的常微分方程起到借鉴意义。

【实验原理】各种常微分方程的数值解法的原理,包括Euler法,改进Euler法,梯形法,Runge-Kutta方法,线性多步方法等。

【实验环境】(使用的软硬件)Matlab软件实验内容:【实验方案设计】我们分别运用Euler法,改进Euler法,RK方法和Adams隐式方法对同一问题进行求解,将数值解和解析解画在同一图像中,比较数值解的精度大小,得出结论。

【实验过程】(实验步骤、记录、数据、分析)我们首先来回顾一下原题:对于给定初值问题:1. 求出其解析解并用Matlab画出其图形;2. 采用Euler法取步长为0.5和0.25数值求解(2.16),并将结果画在同一幅图中,比较两者精度;3. 采用改进Euler法求解(2.16),步长取为0.5;4. 采用四级Runge-Kutta法求解(2.16),步长取为0.5;5. 采用Adams四阶隐格式计算(2.16),初值可由四级Runge-Kutta格式确定。

下面,我们分五个步骤来完成这个问题:步骤一,求出(2.16)式的解析解并用Matlab 画出其图形; ,用Matlab 做出函数在上的图像,见下图:00.51 1.52 2.53 3.54 4.550.511.522.533.5x 1015y=exp(1/3 t 3-1.2t)exact solution图一 初值问题的解析解的图像步骤二,采用Euler 法取步长为0.5和0.25数值求解(2.16),并将结果画在同一幅图中,比较两者精度;我们采用Euler 法取步长为0.5和0.25数值求解,并且将数值解与解析解在一个图中呈现,见下图:00.51 1.52 2.53 3.54 4.550.511.522.533.5x 1015Numerical solution of Euler and exact solutionexact solution h=0.5h=0.25图二 Euler 方法的计算结果与解析解的比较从图像中不难看出,采用Euler 法取步长为0.5和0.25数值求解的误差不尽相同,也就是两种方法的计算精度不同,不妨将两者的绝对误差作图,可以使两种方法的精度更加直观化,见下图:00.51 1.52 2.53 3.54 4.550.511.522.533.5x 1015Absolute error of numerical solution and exact solutionh=0.5h=0.25图三 不同步长的Euler 法的计算结果与解析解的绝对误差的比较 从图像中我们不难看出,步长为0.25的Euler 法比步长为0.5的Euler 法的精度更高。

常微分方程数值解实验报告

常微分方程数值解实验报告

常微分方程数值解实验报告实验报告:常微分方程数值解1.引言常微分方程(Ordinary Differential Equations, ODEs)是数学领域中一个重要的研究对象,涉及到许多自然科学和工程技术领域的问题。

解常微分方程的数值方法是一种求解差分方程的方法,通过计算机找到方程的近似解,对于模拟和预测连续过程非常有用。

本实验旨在通过数值解法,验证和应用常微分方程的解,并比较不同数值方法的精度和效率。

2.实验目的2.1理解常微分方程的基本概念和数值解法;2.2掌握将常微分方程转化为数值求解问题的基本方法;2.3运用数值解法求解常微分方程;2.4比较不同数值解法的精度和效率。

3.实验内容3.1 欧拉方法(Euler Method)给定一个一阶常微分方程dy/dx=f(x,y),通过将其离散为差分形式,欧拉方法可以通过以下递推公式来求解:y_{n+1}=y_n+h*f(x_n,y_n)其中,h为步长,x_n和y_n为当前的x和y值。

3.2 改进的欧拉方法(Improved Euler Method)改进的欧拉方法使用欧拉方法的斜率的平均值来估计每一步中的斜率。

具体公式如下:k1=f(x_n,y_n)k2=f(x_n+h,y_n+h*k1)y_{n+1}=y_n+h*((k1+k2)/2)3.3 二阶龙格-库塔法(Second-order Runge-Kutta Method)二阶龙格-库塔法通过计算每个步骤中的两个斜率来估计每个步长中的斜率。

具体公式如下:k1=f(x_n,y_n)k2=f(x_n+h/2,y_n+(h/2)*k1)y_{n+1}=y_n+h*k24.实验步骤4.1选取常微分方程,并将其转化为数值求解问题的形式;4.2根据给定的初始条件和步长,使用欧拉方法、改进的欧拉方法和二阶龙格-库塔法求解该方程;4.3比较三种方法的数值解与理论解的差异,并分析其精度和效率;4.4尝试不同的步长,观察相应的数值解的变化。

数值计算实验报告-欧拉法常微分方程

数值计算实验报告-欧拉法常微分方程

数学与计算科学学院实验报告实验项目名称欧拉法解常微分方程所属课程名称数值计算实验类型验证型实验日期2012-6- 4班级隧道1002班学号201008020233姓名李彬彬成绩一、实验概述:【实验目的】 通过运用相关的数值计算软件,解决最基本的常微分方程的数值计算,并且能够熟练的运用这种方法。

【实验原理】 欧拉法1.对常微分方程初始问题(9.2))((9.1)),(00⎪⎩⎪⎨⎧==y x y y x f dxdy用数值方法求解时,我们总是认为(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)的近似值。

利用y 1及f (x 1, y 1)又可以算出y (x 2)的近似值:),(1112y x hf y y +=一般地,在任意点x n +1 = (n + 1)h 处y (x )的近似值由下式给出),(1n n n n y x hf y y +=+(9.4)这就是欧拉法的计算公式,h 称为步长。

不难看出,近似解的误差首先是由差商近似代替微商(见(9.3))引起的,这种近似代替所产生的误差称为截断误差。

还有一种误差称为舍入误差,这种误差是由于利用(9.4)进行计算时数值舍入引起的。

【实验环境】Windows XP 环境下运行 NumericalAnalyse 软件二、实验内容:【实验方案】在区间[0,1]上以h=0.1为步长,分别用欧拉法与预估-校正法求初值问题y’=y-2x/y且 y|x=0 =1的数值解。

将上述方程输入到软件NumericalAnalyse中步骤如图选择常微分方程的数值解法。

大学数学实验四_常微分方程数值解

大学数学实验四_常微分方程数值解
【实验内容】
2 小型火箭初始质量为 1400kg,其中包括 1080kg 燃料。火箭竖直向上发射时,燃料燃烧率 为 18kg/s,由此产生 32000N 推力,火箭引擎在燃料用尽时关闭。设火箭上升时空气阻力正 比于速度的平方,比例系数为 0.4kg/m,求引擎关闭瞬间火箭的高度、速度、加速度,及火 箭到达最高点时的高度和加速度,并画出高度、速度、加速度随时间变化的图形。
6
输出结果为
故,当 t =60s 时,火箭的速度、高度、加速度为
v (m/s)
向前欧拉公式
266.94
经典龙格-库塔公式
266.94
ode45(@rocket,t,v)
266.93
h (m) 12131 12131 12131
从上表可以看出,用三种不同方法求得的结果相差不大。
60s 后,火箭作竖直上抛运动,加速度为
4
输出的结果为
输出的火箭的 h – t 图线如下:
经典龙格-库塔公式: n=100000; h=60/n; t=0:h:60; v(1)=0; z(1)=0; g(1)=90/7; for i=1:n
a = v(i); s = z(i); x1 = h*i; k1 = (32000-0.4*a.^2)/(1400-18*x1)-10; x1 = x1 + h/2; b = a + h/2*k1; k2 = (32000-0.4*b.^2)/(1400-18*x1)-10; c = a + h/2*k2; k3 = (32000-0.4*c.^2)/(1400-18*x1)-10; x1 = x1 + h/2;
n = 500 267.0539 n = 100000 266.9432

微分方程数值解实验报告

微分方程数值解实验报告

微分方程数值解实验报告实验目的:掌握微分方程数值解的基本方法,能够利用计算机编程求解微分方程。

实验原理:微分方程是自然科学与工程技术中常见的数学模型,它描述了变量之间的关系及其随时间、空间的变化规律。

解微分方程是研究和应用微分方程的基础,但有很多微分方程无法找到解析解,只能通过数值方法进行求解。

本实验采用欧拉方法和改进的欧拉方法求解微分方程的初值问题:$$\begin{cases}\frac{dy}{dt}=f(t,y)\\y(t_0)=y_0\end{cases}$$其中,$f(t,y)$是给定的函数,$y(t_0)=y_0$是已知的初值条件。

欧拉方法是最基本的数值解法,其步骤如下:1.给定$t_0$和$y_0$2.计算$t_{i+1}=t_i+h$,其中$h$是步长3. 计算$y_{i+1}=y_i+hf(t_i,y_i)$4.重复步骤2、3直到达到终止条件改进的欧拉方法是对欧拉方法进行改进,通过利用函数$y(t)$在$t+\frac{1}{2}h$处的斜率来更准确地估计$y_{i+1}$,其步骤如下:1.给定$t_0$和$y_0$2.计算$t_{i+1}=t_i+h$,其中$h$是步长3. 计算$y_*=y_i+\frac{1}{2}hf(t_i,y_i)$4. 计算$y_{i+1}=y_i+hf(t_i+\frac{1}{2}h,y_*)$5.重复步骤2、3、4直到达到终止条件实验步骤:1.编写程序实现欧拉方法和改进的欧拉方法2.给定微分方程和初值条件3.设置步长和终止条件4.利用欧拉方法和改进的欧拉方法求解微分方程5.比较不同步长下的数值解与解析解的误差6.绘制误差-步长曲线,分析数值解的精度和收敛性实验结果:以一阶常微分方程$y'=3ty+t$为例,给定初值$y(0)=1$,取步长$h=0.1$进行数值求解。

利用欧拉方法求解微分方程得到的数值解如下:\begin{array}{cccc}t & y_{\text{exact}} & y_{\text{Euler}} & \text{误差} \\ \hline0.0&1.000&1.000&0.000\\0.1&1.035&1.030&0.005\\0.2&1.104&1.108&0.004\\0.3&1.212&1.217&0.005\\0.4&1.360&1.364&0.004\\0.5&1.554&1.559&0.005\\0.6&1.805&1.810&0.005\\0.7&2.131&2.136&0.005\\0.8&2.554&2.560&0.006\\0.9&3.102&3.107&0.006\\1.0&3.807&3.812&0.005\\\end{array}利用改进的欧拉方法求解微分方程得到的数值解如下:\begin{array}{cccc}t & y_{\text{exact}} & y_{\text{Improved Euler}} & \text{误差} \\\hline0.0&1.000&1.000&0.000\\0.1&1.035&1.035&0.000\\0.2&1.104&1.103&0.001\\0.3&1.212&1.211&0.001\\0.4&1.360&1.358&0.002\\0.5&1.554&1.552&0.002\\0.6&1.805&1.802&0.003\\0.7&2.131&2.126&0.005\\0.8&2.554&2.545&0.009\\0.9&3.102&3.086&0.015\\1.0&3.807&3.774&0.032\\\end{array}误差-步长曲线如下:实验结论:通过对比欧拉方法和改进的欧拉方法的数值解与解析解的误差,可以发现改进的欧拉方法具有更高的精度和收敛性。

实验报告七常微分方程初值问题的数值解法

实验报告七常微分方程初值问题的数值解法

浙江大学城市学院实验报告课程名称数值计算方法实验项目名称常微分方程初值问题的数值解法 实验成绩指导老师签名日期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四.实验结果与分析。

第5次实验报告(常微分方程初值问题的数值解法)

第5次实验报告(常微分方程初值问题的数值解法)

班级: 学号: 姓名: 成绩:实验5 常微分方程初值问题的数值解法实验1实验目的1)熟悉欧拉法、改进欧拉法和龙格-库塔法的原理。

2)根据以上方法,编程求解常微分方程初值问题。

2 实验内容(1)编写程序,用以上各种方法求解教材P232例7-1、习题6、11的初值问题。

(2) 使用系统自带的函数dsolve 和ode45求例7-1的符号解析解和数值解。

3实验原理求解微分方程初值问题00(,)()y f x y y x y '=⎧⎨=⎩ (1) 欧拉法(显式):10(,)n n n n n y y hf x y x x nh +=+⎧⎨=+⎩(2) 改进欧拉法(0)1(0)111(,)(,)(,)2n n n n n n n n n n y y hf x y h y y f x y f x y ++++⎧=+⎪⎨⎡⎤=++⎪⎣⎦⎩ (3) 经典龙格-库塔法(四阶)11234121324300(22)6(,)(,)22(,)22(,)()i i i i i i i i i i h y y K K K K K f x y h h K f x y K h h K f x y K K f x h y hK y y x +⎧=++++⎪⎪=⎪⎪=++⎪⎨⎪=++⎪⎪⎪=++⎪=⎪⎩4实验步骤1)建立函数文件,根据各公式编写程序;2)上机调试程序,运行程序进行计算,记录计算结果;3)分析各公式计算结果,比较各公式的优缺点。

5 程序设计欧拉法改进欧拉法function Euler1(x0,y0,h,n)%(x0,y0):方程的初值%h:步长%n:计算的步数for i=1:nx=x0+h;yp=y0+h*f(x0,y0);yc=y0+h*f(x,yp);y=(yp+yc)/2;x %在屏幕显示每一步的x值y %在屏幕显示每一步计算的方程的数值解 x0=x;y0=y;end经典龙格-库塔法1)函数function f=f(x,y) f=y-2*x/y;end6总结注:若要更改matlab计算的数值类型,可以通过在matlab中设置实现:File -> Preferences ->Array Editor窗口中,Format 下方将Default array format设置为:long解微分方程的MATLAB命令MATLAB中主要用dsolve求微分方程的符号解析解,ode45求数值解。

微分方程数值解实验报告

微分方程数值解实验报告

微分方程数值解法课程设计报告班级:_______姓名: ___学号:__________成绩:2017年 6月 21 日摘要自然界与工程技术中的很多现象,可以归结为微分方程定解问题。

其中,常微分方程求解是微分方程的重要基础内容。

但是,对于许多的微分方程,往往很难得到甚至不存在精确的解析表达式,这时候,数值解提供了一个很好的解决思路。

,针对于此,本文对常微分方程数值解法进行了简单研究,主要讨论了一些常用的数值解法,如欧拉法、改进的欧拉法、Runge—Kutta方法、Adams法以及椭圆型方程、抛物型方程的有限差分方法等,通过具体的算例,结合MATLAB求解画图,初步给出了一般常微分方程数值解法的求解过程。

同时,通过对各种方法的误差分析,让大家对各种方法的特点和适用范围有一个直观的感受。

关键词:微分方程数值解、MATLAB目录摘要 (2)目录 (3)第一章常微分方程数值解法的基本思想与原理 (4)1.1 常微分方程数值解法的基本思路 (4)1.2用matlab编写源程序 (4)1.3 常微分方程数值解法应用举例及结果 (5)第二章常系数扩散方程的经典差分格式的基本思想与原理 (6)2.1 常系数扩散方程的经典差分格式的基本思路 (6)2.2 用matlab编写源程序 (7)2.3 常系数扩散方程的经典差分格式的应用举例及结果 (8)第三章椭圆型方程的五点差分格式的基本思想与原理 (10)3.1 椭圆型方程的五点差分格式的基本思路 (10)3.2 用matlab编写源程序 (10)3.3 椭圆型方程的五点差分格式的应用举例及结果 (12)第四章总结 (12)参考文献 (12)第一章常微分方程数值解法的基本思想与原理1.1常微分方程数值解法的基本思路常微分方程数值解法(numerical methods forordinary differential equations)计算数学的一个分支.是解常微分方程各类定解问题的数值方法.现有的解析方法只能用于求解一些特殊类型的定解问题,实用上许多很有价值的常微分方程的解不能用初等函数来表示,常常需要求其数值解.所谓数值解,是指在求解区间内一系列离散点处给出真解的近似值.这就促成了数值方法的产生与发展.1.2用matlab编写源程序龙格库塔法:M文件:function dx=Lorenz(t,x)%r=28,sigma=10,b=8/3dx=[-10*(x(1)-x(2));-x(1)*x(3)+28*x(1)-x(2);x(1)*x(2)-8*x(3)/3];运行程序:x0=[1,1,1];[t,y]=ode45('Lorenz',[0,100],x0);subplot(2,1,1) %两行一列的图第一个plot(t,y(:,3))xlabel('time');ylabel('z');%画z-t图像subplot(2,2,3) %两行两列的图第三个plot(y(:,1),y(:,2))xlabel('x');ylabel('y'); %画x-y图像subplot(2,2,4)plot3(y(:,1),y(:,2),y(:,3))xlabel('x');ylabel('y');zlabel('z');%画xyz图像欧拉法:h=0.010;a=16;b=4;c=49.52;x=5;y=10;z=10;Y=[];for i=1:800x1=x+h*a*(y-x);y1=y+h*(c*x-x*z-y);z1=z+h*(x*y-b*z);x=x1;y=y1;z=z1;Y(i,:)=[x y z];endplot3(Y(:,1),Y(:,2),Y(:,3));1.3常微分方程数值解法的应用举例及结果应用举例:⎪⎪⎪⎩⎪⎪⎪⎨⎧-=--=--=)()()()()()()()()())()(()(t bz t y t x dt t dz t z t x t y t rx dt t dy t y t x a dt t dx a=10,b=8/3,0<r<+∞,当1<r<24.74时,Lorenz 方程有两个稳定的不动点c()1(-r b ,)1(-r b ,r-1)和c '(-)1(-r b ,-)1(-r b ,r-1),一个稳定的不动点0=(0,0,0),当r>24.74时,c 和c '都变成不稳定的,此时存在混沌和奇怪吸引子。

【清华】实验4-常微分方程数值解(011813)

【清华】实验4-常微分方程数值解(011813)

实验4 常微分方程数值解化学工程系分9班焦阳2009011813 【实验目的】1. 掌握用MATLAB软件求微分方程初值问题数值解的方法;2. 通过实例学习用微分方程模型解决简化的实际问题;3. 了解欧拉方法和龙格-库塔方法的基本思想和计算公式,及稳定性等概念。

【实验内容】1.题目3:小型火箭初始重量为1400kg,其中包括1080kg燃料。

火箭竖直向上发射时燃料燃烧率为18kg/s,由此产生32000N的推力,火箭引擎在燃烧用尽时关闭。

设火箭上升时空气阻力正比于速度的平方,比例系数为0.4kg/m,求引擎关闭瞬间火箭的高度、速度、加速度,及火箭到达最高点时的高度和加速度,并画出高度、速度、加速度随时间变化的图形。

建立模型并进行分析:假设火箭在上升过程中,重力加速度g不随高度而变化,即固定g = 9.8m/s^2。

、(1)从火箭开始上升到引擎关闭:设火箭质量为m,高度为h,速度为v,加速度为a,阻力为f:,,ﭸ由牛顿第二定律可得:总ﭸ综上可得:;ﭸ;初值条件为:,;定义域为:。

根据常微分方程组的初值问题,在MATLAB中计算数值解,记,,, 。

通过解出微分方程的数值解,并进行绘图得到高度-时间曲线,速度-时间曲线,加速度-时间曲线如下:由MA度、速度t(s MATLAB 计算速度、加速度如下t(s) 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 计算得到的火箭度如下表:h(m)06.57326.4459.76106.5166.7240.2326.7425.7536.9659.8793.6937.81091.1254.1425.1604.1790.1983.2181.2384.的火箭从开始上h(m) 0 6.5737 26.444 59.762 106.57 166.79 240.27 326.72 425.79 536.99 659.8 793.63 937.85 1091.8 1254.7 1425.9 1604.8 1790.8 1983.1 2181.2 2384.5 开始上升到关闭引v(m/s)0 13.18926.57740.06253.53566.89 80.02192.829105.22117.11128.43139.14149.18158.55167.23175.22182.55189.22195.27200.75205.7 关闭引擎这段时间/s) 189 577 062 535 021 829 .22 .11 .43 .14 .18 .55 .23 .22 .55 .22 .27 .75 段时间内各时刻a(m/s^2)13.0571413.304513.4532813.4971913.4331313.2613112.9853412.6121912.1519511.6169311.0212710.38 9.7083329.0209048.3309057.6502496.9900526.3593395.7646125.2094884.694626各时刻的高^2)71404532871931313153421919569312733290490524905233961248862621 2592.4 210.18 4.22201222 2804.5 214.19 3.79432623 3020.6 217.79 3.41201724 3240.1 221.01 3.07303925 3462.7 223.92 2.7726326 3687.9 226.56 2.50440927 3915.6 228.97 2.26774828 4145.6 231.14 2.06332529 4377.8 233.11 1.88975930 4611.9 234.91 1.74334931 4847.7 236.57 1.6177932 5085 238.14 1.50617933 5323.8 239.61 1.40954434 5564.1 240.99 1.32933135 5805.8 242.28 1.2650336 6048.7 243.5 1.21393737 6292.9 244.68 1.17083938 6538.1 245.83 1.1302639 6784.5 246.96 1.09469840 7032 248.05 1.06634841 7280.5 249.1 1.04558942 7530.2 250.12 1.03075943 7780.9 251.14 1.01781844 8032.5 252.15 1.00237545 8285.1 253.16 0.98756846 8538.8 254.15 0.97632347 8793.4 255.12 0.96963948 9049 256.07 0.96629149 9305.6 257.03 0.96239250 9563.1 257.99 0.95272351 9821.5 258.95 0.9411852 10081 259.9 0.93372753 10341 260.83 0.93277654 10603 261.75 0.93633455 10865 262.67 0.9369556 11128 263.61 0.92584757 11392 264.54 0.91379358 11657 265.46 0.91062859 11923 266.35 0.9161260 12190 267.26 0.917011根据表格可以很容易得到:关闭引擎的瞬间,h=12190m,v=267.26m/s,a=0.917011m/s^2。

计算方法实验报告-常微分方程的数值解法

计算方法实验报告-常微分方程的数值解法
15.0000 50.1112 50.1112
15.1000 50.1224 50.1224
15.2000 50.1332 50.1332
15.3000 50.1436 50.1436
15.4000 50.1535 50.1535
disp([T',V',df']);
%画图观察效果
figure;
plot(T,df,'k*',T,V,'--r');
grid on;
title('速度v与时间t的函数曲线');
legend('准确值','四阶经典R-K法');
结果如下所示:
步长四阶经典R-K法准确值
0 0 0
0.1000 0.9799 0.9799
5.9000 41.1708 41.1708
6.0000 41.4919 41.4919
6.1000 41.8028 41.8028
6.2000 42.1039 42.1039
6.3000 42.3954 42.3954
6.4000 42.6775 42.6775
6.5000 42.9504 42.9504
11.5000 49.2667 49.2667
11.6000 49.3097 49.3097
11.7000 49.3511 49.3511
11.8000 49.3909 49.3909
11.9000 49.4292 49.4292
12.0000 49.4661 49.4661
12.1000 49.5016 49.5016
(2)某跳伞者在t=0时刻从飞机上跳出,假设初始时刻的垂直速度为0,且跳伞者垂直下落。已知空气阻力为F=cv2,其中c为常数,v为垂直速度,向下方方向为正。写出此跳伞者的速度满足的微分方程;若此跳伞者的质量为M=70kg,且已知c=0.27kg/m,利用四阶龙格-库塔公式计算t<=20s的速度(取h=0.1s)

数值分析实验报告之常微分方程数值解

数值分析实验报告之常微分方程数值解

数学与计算科学学院实验报告
实验项目名称常微分方程数值解
所属课程名称数值方法B
实验类型验证
实验日期2013.11.11
班级
学号
姓名
成绩
图1 h=0.1时三个方法走势图
图2 h=0.05时三个方法走势图
图4 h=0.05时三个方法走势图
附录1:源程序
附录2:实验报告填写说明
1.实验项目名称:要求与实验教学大纲一致。

2.实验目的:目的要明确,要抓住重点,符合实验教学大纲要求。

3.实验原理:简要说明本实验项目所涉及的理论知识。

4.实验环境:实验用的软、硬件环境。

5.实验方案(思路、步骤和方法等):这是实验报告极其重要的内容。

概括整个实验过程。

对于验证性实验,要写明依据何种原理、操作方法进行实验,要写明需要经过哪几个步骤来实现其操作。

对于设计性和综合性实验,在上述内容基础上还应该画出流程图、设
计思路和设计方法,再配以相应的文字说明。

对于创新性实验,还应注明其创新点、特色。

6.实验过程(实验中涉及的记录、数据、分析):写明具体实验方案的具体实施步骤,包括实验过程中的记录、数据和相应的分析。

7.实验结论(结果):根据实验过程中得到的结果,做出结论。

8.实验小结:本次实验心得体会、思考和建议。

9.指导教师评语及成绩:指导教师依据学生的实际报告内容,给出本次实验报告的评价。

实验2常微分方程数值解法

实验2常微分方程数值解法

实验2 常微分方程数值解法实验目的与要求1)熟悉求解常微分方程初值问题的有关方法和理论,主要是改进 欧拉法和四阶龙格库塔法2)会编制改进欧拉法的计算程序;3)针对实验内容编制的程序能够正确调试并运行所需结果实验内容:请用改进欧拉法和四阶龙格库塔算法解算法介绍 1、改进欧拉算法概要 20011112111205(0)2(,)()(,)[(,)(,)](,)[(,)(b an i i i i hi i i i i i i i i i h i i i i i y xy x y y f x y a x b y x y y y hf x y y y f x y f x y y y hf x y y y f x y f x 解一阶常微分方程初值问题将区间[a,b]作n 等分,取步长h=欧拉公式为梯形公式为改进欧拉法的公式为-+++++++ì¢ï=-#ïíï=ïîì¢=#ïïíï=ïî=+=++=+=++111,)](,)(,)1()2i p i i i c i i i i p c y y y hf x y y y hf x y y y y 或表示为+++ìïïíïïîìïïï=+ïïï=+íïïïï=+ïïî2、龙格库塔算法概要算法函数示例:void ModEuler(float(*f)(float,float),float x0,float y0,float xn,int n){int i;float yp,yc,x=x0,y=y0,h=(xn-x0)/n;printf("x[0]=%f\ty[0]=%f\n",x,y);for(i=1;i<=n;i++){yp=y+h*(*f)(x,y);x=x0+i*h;yc=y+h*(*f)(x,yp);y=(yp+yc)/2;printf("x[%d]=%f\ty[%d]=%f\n",i,x,i,y);}}void Runge_Kutta(float(*f)(float x,float y),float a,float b,float y0,int n){float x=a,y=y0,K1,K2,K3,K4;float h=(b-a)/N;int i;printf("x[0]=%f\ty[0]=%f\n",x,y);for(i=1;i<=n;i++){K1=(*f)(x,y);K2=(*f)(x+h/2,y+h*K1/2);K3=(*f)(x+h/2,y+h*K2/2);K4=(*f)(x+h,y+h*K3);上各个节点初的近似值,在区间出发,可得未知函数,由初值取步长为库塔公式常用的四阶龙格对于初值问题:][a )(),()2,2()2,2(),()22(6)(...........).........,(034231214321100'b x y y h hk 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 h y y y x y b x a y x f y i i i i i i i i i i ⎪⎪⎪⎪⎩⎪⎪⎪⎪⎨⎧++=++=++==++++=-⎩⎨⎧=≤≤=+y=y+h*(K1+2*K2+2K3+K4)/6;x=a+i*h;printf("x[%d]=%f\ty[%d]=%f\n",i,x,i,y); }}。

常微分方程的数值解法实验报告

常微分方程的数值解法实验报告

常微分方程的数值解法实验报告实验报告:常微分方程的数值解法摘要:常微分方程(ODE)是描述动力学系统中物理量随时间变化的数学方程,广泛应用于自然科学和工程领域。

然而,对于一些复杂的非线性ODE,很难找到解析解。

因此,我们需要数值解法来求解这些方程。

本实验报告将介绍四种常见的常微分方程数值解法,分别是欧拉法、改进的欧拉法、四阶龙格-库塔法和自适应步长的龙格-库塔法,并通过数值实验比较它们的精度和效率。

1.引言在实际问题中,许多物理量的变化规律可以由常微分方程描述。

然而,对于复杂的非线性ODE,很难找到解析解。

因此,为了解决这类问题,我们需要借助数值方法来求解。

2.方法本实验采用四种常见的常微分方程数值解法:欧拉法、改进的欧拉法、四阶龙格-库塔法和自适应步长的龙格-库塔法。

(1)欧拉法是最简单的数值解法,通过将微分方程转化为差分方程,使用离散的步长来近似微分方程。

(2)改进的欧拉法在欧拉法的基础上进行了改进,使用预估-校正的方法来提高精度。

(3)四阶龙格-库塔法是一种经典的数值解法,通过利用不同步长处的斜率来近似微分方程,具有较高的精度。

(4)自适应步长的龙格-库塔法是在四阶龙格-库塔法的基础上改进而来的,根据步长的大小自适应地选择不同的步长,同时保证精度和效率。

3.实验设计为了比较这四种数值解法的精度和效率,我们设计了两个实验。

实验一是求解一阶常微分方程:dy/dx = -2x,初始条件y(0) = 1,解析解为y = 1 - x^2、实验二是求解二阶常微分方程:d^2y/dx^2 + y = 0,初始条件y(0) = 0,dy/dx(0) = 1,解析解为y = sin(x)。

4.结果与分析实验一中,比较四种数值解法在不同步长下的近似解和解析解,计算其误差。

实验结果表明,四阶龙格-库塔法和自适应步长的龙格-库塔法具有较高的精度,而欧拉法和改进的欧拉法的精度较低。

实验二中,我们比较四种数值解法在不同步长下的近似解和解析解,并计算其误差。

常微分方程数值解---实验报告

常微分方程数值解---实验报告

琼州学院实验报告课程名称:__________ ___开课学期:____ ______院(部):__理工学院_________开课实验室:________ __学生姓名:__梁小叶_______ __专业班级:________ __学号:__________常微分方程数值解---实验报告一 实验目的: 1.掌握用MATLAB 求微分方程初值问题数值解的方法;2.通过实例学习微分方程模型解决简化的实际问题;3.了解欧拉方法和龙格库塔方法的基本思想。

二 实验内容:用欧拉方法和龙格库塔方法求下列微分方程初值问题的数值解,画出解的图形,对结果进行分析比较22'2,(1)(01),322;(0)1',(2)(010).(0)0(0)1x y y x x y e x y y x y x y y =+⎧≤≤=--⎨=⎩⎧=-≤≤⎨==⎩精确解或三 问题分析:怎样设计程序?流程图?变量说明?能否将某算法设计成具有形式参数的函数形式?【程序如下】:function f=f(x,y)f=y+2*x;clc;clear;a=0;b=1; %求解区间[x1,y_r]=ode45('f',[a b],1); %调用龙格库塔求解函数求解数值解;%% 以下利用Euler 方法求解y(1)=1;N=100;h=(b-a)/N;x=a:h:b;for i=1:Ny(i+1)=y(i)+h*f(x(i),y(i));endfigure(1)plot(x1,y_r,'r*',x,y,'b+',x,3*exp(x)-2*x-2,'k-');%数值解与真解图title('数值解与真解图');legend('RK4','Euler','真解');xlabel('x');ylabel('y');figure(2)plot(x1,abs(y_r-(3*exp(x1)-2*x1-2)),'k-');%龙格库塔方法的误差title('龙格库塔方法的误差')xlabel('x');ylabel('Error');figure(3)plot(x,abs(y-(3*exp(x)-2*x-2)),'r-')%Euler 方法的误差title('Euler 方法的误差')xlabel('x');ylabel('Error');【运行结果如下】:5总结体会:自己写了哦。

数学实验报告——利用MALTAB计算常微分方程数值解

数学实验报告——利用MALTAB计算常微分方程数值解

实验二 常微分方程数值解一、火箭飞行器㈠问题描述小型火箭初始质量为1400kg ,其中包括1080kg 燃料,火箭竖直向上发射时燃料燃烧率为18kg/s ,由此产生32000N 的推力,火箭引擎在燃料用尽时关闭,设火箭上升时空气阻力正比于速度的平方,比例系数为0.4kg/m ,求引擎关闭瞬间火箭的高度,速度,加速度及火箭达到最高点时的高度和加速度,并画出高度、速度、加速度随时间变化的图形。

㈡方法与公式1、简要分析本题的求解需要用到常微分方程,而整个过程又被分为两个阶段:火箭加速上升阶段和燃料燃尽后减速的阶段。

由题目易知第一个阶段持续时间T 1=108018=60s 列出第一阶段的方程组:设M0为火箭本身质量,m 为燃料质量,g 为重力加速度 = 9.8m/s 2,燃料燃烧率为a ,空气阻力的比例系数为k ,F 为推进力。

M0 = 1400-1080 = 320kg; v̇=F−kv 2−(M0+m )g M0+mm =−a初值v̇=0,m =1080。

由以上各式可以求出t=T 1时火箭的速度。

再求解第二阶段:v̇=−kv 2−M0g M0m =0可以求出火箭速度降为0的时刻。

将整个过程中的时间向量以及速度向量联合起来,利用第三章所学插值与数值积分的方法可以求得任意时刻火箭的近似高度。

2、方法求解常微分方程时,我分别采用了自己编写的欧拉公式、改进欧拉公式、4级4阶龙格-库塔公式,以及MATLAB自带的龙格-库塔方法,求解数值积分时采用辛普森公式。

由于Matlab自带的Simpson公式是自适应的,因此需要使用自己在上一次实验时所编的Simpson公式。

㈢结果与分析1、各种公式的对比首先,我作出了各种不同公式计算得到的火箭速度随时间变化的图像,图如下:从图中可以看出,各种公式计算得到的结果基本一致,为确定其区别,将图像放大,放大约2000 倍后,得到下图:分析:从图中可以看出,自编欧拉公式距离MATLAB 自带龙格-库塔公式最远,精度最差;自编的改进欧拉公式和自编的龙格-库塔公式结果基本一致,两者中自编龙格-库塔公式距MATLAB 自带龙格-库塔公式的结果稍近。

常微分实习报告-常微分方程数值求解问题

常微分实习报告-常微分方程数值求解问题

(数学类课程)课程实习报告课程名称:常微分方程课程实习实习题目:常微分方程数值求解问题的实习姓名:系:专业:年级:学号:指导教师:职称:年月日课程实习报告结果评定目录1. 实习的目的和任务 (1)2. 实习要求 (1)3. 实习地点 (1)4. 主要仪器设备 (1)5. 实习内容 (1)5.1 用不同格式对同一个初值问题的数值求解及其分析 (1)5.1.1求精确解 (1)5.1.2用欧拉法求解 (6)5.1.3用改进欧拉法求解 (9)5.1.4用4级4阶龙格—库塔法求解 (12)5.1.5 问题讨论与分析 (15)5.2 一个算法不同不长求解同一个初值问题及其分析 (18)6. 结束语 (29)参考文献 (29)常微分方程课程实习1. 实习的目的和任务目的:通过课程实习能够应用MATLAB 软来计算微分方程(组)的数值解;了解常微分方程数值解。

任务:通过具体的问题,利用MATLAB 软件来计算问题的结果,分析问题的结论。

2. 实习要求能够从案例的自然语言描述中,抽象出其中的数学模型;能够熟练应用所学的数值解计算方法;能够熟练使用MATLAB 软件;对常微分方程数值解有所认识,包括对不同算法有所认识和对步长有所认识。

3. 实习地点数学实验室、学生宿舍 4. 主要仪器设备计算机、Microsoft Windows 7 Matlab 7.0 5. 实习内容5.1 用欧拉方法,改进欧拉方法,4阶龙格—库塔方法分别求下面微分方程的初值dy/dx=y*cos(x+2) y(-2)=1 x ∈[-2,0] 5.1.1求精确解 ①变量分离方程情形:形如()*()dyf xg y dx=的方程,这里(),()f x g y 分别是,x y 的连续函数.如果()0g y ≠,我们可将方程改写成()()dyf x dxg y =,这样,变量就”分离”开来了,两边同时积分即可:(),()dyf x dx c cg y =+⎰⎰为任意常数. ②常数变易法:一阶线性微分方程()*()dyf x yg x dx =+,其中(),()f x g x 在考虑区 间上是的连续函数.可先解出方程()*dyf x y dx=的解,这是属于变量分离方程情形,可解得:*exp(())y c f x dx =⎰,这里c 是任意常数.然后将变c易为x 的待定函数()c x ,令()*exp(())y c x f x dx =⎰,将其代入原方程可得:()*exp(())()*()*exp(())()*()*exp(())()dy dc x f x dx c x f x f x dx f x c x f x dx g x dx dx=+=+⎰⎰⎰所以可解得()()*exp(())*1c x f x f x dx dx c =-⎰⎰,这里1c 是任意常数.将()()*exp(())1c x f x f x dx dx c =-+⎰⎰代入()*exp(())y c x f x dx =⎰可得原方程的通解:(()*exp(())1)*exp(()),1y f x f x dx dx c f x dx c =-+⎰⎰⎰为任意常数.③恰当微分方程情形:形如(,)(,)0m x y dx n x y dy +=的一阶微分方程,这里 假设(,),(,)m x y n x y 在某矩形域内是的连续函数,且具有连续的一阶偏导数. 若m ny x∂∂=∂∂,则为恰当微分方程.判断为恰当微分方程后,则可用如下解法: 设(,)u x y c =是原方程的解,则um x∂=∂,所以设(,)()u m x y dx v y =+⎰, 则((,)())((,))()m x y dx v y m x y dx udv y n y y y dy ∂+∂∂==+∂∂∂⎰⎰=,所以((,))()m x y dx dv y n ydy ∂-=∂⎰,由此()[(,)]v y n m x y dx dy y∂=-∂⎰⎰,由此可解得(,)[(,)]u m x y dx n m x y dx dy y ∂=+-∂⎰⎰⎰,所以原方程的通解为(,)[(,)],m x y dx n m x y dx dy c c y∂+-=∂⎰⎰⎰为任意常数。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

改进的欧拉公式可以改写为:
f x, y的值两次,截断误差为 o h3 。
yi1
k1

yi
hf xi
yi1
k1 , yi

,它每一步计算
yi

k1 hf xi , yi
k2 hf xi h, yi k1

改进的欧拉方法之所以比欧拉方法具有更高的精度,是因为在每一步它都
三、实验原理与理论基础
(一) 欧拉法算法设计
对常微分方程初始问题
d y f (x, y) d x
y(x0 ) y0
用数值方法求解时,我们总是认为(6-1)、(6-2)的解存在且唯一。
欧拉法是解初值问题的最简单的数值方法。从(6-2)式由于 y (x0) = y0 已给定, 因而可以算出 y'(x0 ) f (x0 , y0 ) 。
(4)利用四阶标准 R- K 方法求二阶方程初值问题的数值解
(I)
(II)
y 3y 2 y 0

y0

0,
y0
y 0.1(1 y 2 ) y y 0

y0
y


1,
y ex 1
(III) y0 1, y0 0
xi
1 6
k4 hf xi h, yi k3
k1
yi


h, 2
h, 2
2k2
yi
yi


1 2
1 2
2k3
k1
k2


k4
此公式称为标准四阶龙格-库塔(Runge-Kutta)公式,它的截断误差为 o h5 。
虽然用龙格-库塔方法每一步需要四次调用 f ,计算量较改进的欧拉方法大一倍,

f (x0 , y0 )
学号:
(6-1)
(6-2)
(6-3)
一般地,在任意点 xn1 n 1h 处 y(x) 的近似值由下式给出
yn1 yn hf (xn , yn ) 这就是欧拉法的计算公式,h 称为步长。
(2)四阶龙格-库塔法算法设计:
欧拉公式可以改写为:
截断误差为 o h2 。

y'(x0 )
对全部高中资料试卷电气设备,在安装过程中以及安装结束后进行高中资料试卷调整试验;通电检查所有设备高中资料电试力卷保相护互装作置用调与试相技互术关,系电,力根通保据过护生管高产线中工敷资艺设料高技试中术卷资,配料不置试仅技卷可术要以是求解指,决机对吊组电顶在气层进设配行备置继进不电行规保空范护载高与中带资负料荷试下卷高问总中题体资,配料而置试且时卷可,调保需控障要试各在验类最;管大对路限设习度备题内进到来行位确调。保整在机使管组其路高在敷中正设资常过料工程试况中卷下,安与要全过加,度强并工看且作护尽下关可都于能可管地以路缩正高小常中故工资障作料高;试中对卷资于连料继接试电管卷保口破护处坏进理范行高围整中,核资或对料者定试对值卷某,弯些审扁异核度常与固高校定中对盒资图位料纸置试,.卷保编工护写况层复进防杂行腐设自跨备动接与处地装理线置,弯高尤曲中其半资要径料避标试免高卷错等调误,试高要方中求案资技,料术编试交写5、卷底重电保。要气护管设设装线备备置敷4高、调动设中电试作技资气高,术料课中并3中试、件资且包卷管中料拒含试路调试绝线验敷试卷动槽方设技作、案技术,管以术来架及避等系免多统不项启必方动要式方高,案中为;资解对料决整试高套卷中启突语动然文过停电程机气中。课高因件中此中资,管料电壁试力薄卷高、电中接气资口设料不备试严进卷等行保问调护题试装,工置合作调理并试利且技用进术管行,线过要敷关求设运电技行力术高保。中护线资装缆料置敷试做设卷到原技准则术确:指灵在导活分。。线对对盒于于处调差,试动当过保不程护同中装电高置压中高回资中路料资交试料叉卷试时技卷,术调应问试采题技用,术金作是属为指隔调发板试电进人机行员一隔,变开需压处要器理在组;事在同前发一掌生线握内槽图部内 纸故,资障强料时电、,回设需路备要须制进同造行时厂外切家部断出电习具源题高高电中中源资资,料料线试试缆卷卷敷试切设验除完报从毕告而,与采要相用进关高行技中检术资查资料和料试检,卷测并主处且要理了保。解护现装场置设。备高中资料试卷布置情况与有关高中资料试卷电气系统接线等情况,然后根据规范与规程规定,制定设备调试高中资料试卷方案。
比欧拉方法多计算了一次 f x, y的值。因此,要进一步提高精度,可以考虑在 每一步增加计算 f x, y的次数。
如果考虑在每一步计算 f x, y的值四次,则可以推得如下公式:
yi1



k1
k2
k3 f
hf
yi

xi ,


xi
专业班级:信息软件
120108010002
常微分方程的数值解法
姓名:吴中原
一、实验目的
1、熟悉各种初值问题的算法,编出算法程序; 2、明确各种算法的精度与所选步长有密切关系;通过计算更加了解各 种
算法的优越性。
二、实验题目
1、根据初值问题数值算法,分别选择二个初值问题编程计算; 2、试分别取不同步长,考察某节点 x j 处数值解的误差变化情况; 3、试用不同算法求解某初值问题,结果有何异常; 4、分析各个算法的优缺点。
设 x1 = h 充分小,则近似地有:
y(x1 ) h
y(x0 )
记 yi y(xi ) i 0,1,,n 从而我们可以取
y1 y0 hf (x0 , y0 ) 作为 y(x1 ) 的近似值。利用 y 1 及 f (x1, y1)又可以算出 y(x2 ) 的近似值:
y2 y1 hf (x1, y1 )
y(1)=y0;
for k=1:M
end
y(k+1)=y(k)+h*t(k)*y(k)^(1/3);
1,
y0 0
y0
1
0
0 x 1
h 0.02
0 x 1, h 0.1
0 x 1
0 x 2, h 0.1 0 x 4, h 0.2
对全部高中资料试卷电气设备,在安装过程中以及安装结束后进行高中资料试卷调整试验;通电检查所有设备高中资料电试力卷保相护互装作置用调与试相技互术关,系电,力根通保据过护生管高产线中工敷资艺设料高技试中术卷资,配料不置试仅技卷可术要以是求解指,决机对吊组电顶在气层进设配行备置继进不电行规保空范护载高与中带资负料荷试下卷高问总中题体资,配料而置试且时卷可,调保需控障要试各在验类最;管大对路限设习度备题内进到来行位确调。保整在机使管组其路高在敷中正设资常过料工程试况中卷下,安与要全过加,度强并工看且作护尽下关可都于能可管地以路缩正高小常中故工资障作料高;试中对卷资于连料继接试电管卷保口破护处坏进理范行高围整中,核资或对料者定试对值卷某,弯些审扁异核度常与固高校定中对盒资图位料纸置试,.卷保编工护写况层复进防杂行腐设自跨备动接与处地装理线置,弯高尤曲中其半资要径料避标试免高卷错等调误,试高要方中求案资技,料术编试交写5、卷底重电保。要气护管设设装线备备置敷4高、调动设中电试作技资气高,术料课中并3中试、件资且包卷管中料拒含试路调试绝线验敷试卷动槽方设技作、案技术,管以术来架及避等系免多统不项启必方动要式方高,案中为;资解对料决整试高套卷中启突语动然文过停电程机气中。课高因件中此中资,管料电壁试力薄卷高、电中接气资口设料不备试严进卷等行保问调护题试装,工置合作调理并试利且技用进术管行,线过要敷关求设运电技行力术高保。中护线资装缆料置敷试做设卷到原技准则术确:指灵在导活分。。线对对盒于于处调差,试动当过保不程护同中装电高置压中高回资中路料资交试料叉卷试时技卷,术调应问试采题技用,术金作是属为指隔调发板试电进人机行员一隔,变开需压处要器理在组;事在同前发一掌生线握内槽图部内 纸故,资障强料时电、,回设需路备要须制进同造行时厂外切家部断出电习具源题高高电中中源资资,料料线试试缆卷卷敷试切设验除完报从毕告而,与采要相用进关高行技中检术资查资料和料试检,卷测并主处且要理了保。解护现装场置设。备高中资料试卷布置情况与有关高中资料试卷电气系统接线等情况,然后根据规范与规程规定,制定设备调试高中资料试卷方案。
y1 y2 y2 y1
y3 y3
y1 0 1 y2 0 0
y3 0 1
取步长 h 0.01,计算 y0.05, y0.10, y0.15数值解,参考结果
y10.15 0.9880787, y2 0.15 0.1493359, y3 0.15 0.8613125
yb=y(M+1); yy=((t.^2+2)./3).^1.5; det=yy-y; plot(t,y,'r-',t,yy,'b:',t,det);
这里由于龙格-库塔方法的步长增大了一倍,因而两种方法总的计算量相同,但
龙格-库塔方法精确度更高。所以龙格-库塔公式兼顾了精度和计算工作量的较
为理想的公式,在实际计算中最为常用。
四、实验内容
(一)问题重述:
对全部高中资料试卷电气设备,在安装过程中以及安装结束后进行高中资料试卷调整试验;通电检查所有设备高中资料电试力卷保相护互装作置用调与试相技互术关,系电,力根通保据过护生管高产线中工敷资艺设料高技试中术卷资,配料不置试仅技卷可术要以是求解指,决机对吊组电顶在气层进设配行备置继进不电行规保空范护载高与中带资负料荷试下卷高问总中题体资,配料而置试且时卷可,调保需控障要试各在验类最;管大对路限设习度备题内进到来行位确调。保整在机使管组其路高在敷中正设资常过料工程试况中卷下,安与要全过加,度强并工看且作护尽下关可都于能可管地以路缩正高小常中故工资障作料高;试中对卷资于连料继接试电管卷保口破护处坏进理范行高围整中,核资或对料者定试对值卷某,弯些审扁异核度常与固高校定中对盒资图位料纸置试,.卷保编工护写况层复进防杂行腐设自跨备动接与处地装理线置,弯高尤曲中其半资要径料避标试免高卷错等调误,试高要方中求案资技,料术编试交写5、卷底重电保。要气护管设设装线备备置敷4高、调动设中电试作技资气高,术料课中并3中试、件资且包卷管中料拒含试路调试绝线验敷试卷动槽方设技作、案技术,管以术来架及避等系免多统不项启必方动要式方高,案中为;资解对料决整试高套卷中启突语动然文过停电程机气中。课高因件中此中资,管料电壁试力薄卷高、电中接气资口设料不备试严进卷等行保问调护题试装,工置合作调理并试利且技用进术管行,线过要敷关求设运电技行力术高保。中护线资装缆料置敷试做设卷到原技准则术确:指灵在导活分。。线对对盒于于处调差,试动当过保不程护同中装电高置压中高回资中路料资交试料叉卷试时技卷,术调应问试采题技用,术金作是属为指隔调发板试电进人机行员一隔,变开需压处要器理在组;事在同前发一掌生线握内槽图部内 纸故,资障强料时电、,回设需路备要须制进同造行时厂外切家部断出电习具源题高高电中中源资资,料料线试试缆卷卷敷试切设验除完报从毕告而,与采要相用进关高行技中检术资查资料和料试检,卷测并主处且要理了保。解护现装场置设。备高中资料试卷布置情况与有关高中资料试卷电气系统接线等情况,然后根据规范与规程规定,制定设备调试高中资料试卷方案。
相关文档
最新文档