实验七、常微分方程数值解法.
实验报告——常微分方程的数值解法
实验报告实验项目名称常微分方程的数值解法实验室数学实验室所属课程名称微分方程数值解实验类型上机实验实验日期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 法的精度更高。
第7章 常微分方程初值问题的数值解法
例1 函数f ( t , y ) = t y 在区域D0 = {( t , y ) | 1 ≤ t ≤ 2, −3 ≤ y ≤ 4}
关于y满足Lipschitz条件，相应的Lipschitz常数可取为L = 2
3 存在性定理 定理1 设函数f ( t , y )在凸集D ⊂ R 2中有定义，若存在常数
(7.2.7)
称为显式Runge-Kutta(龙格-库塔 )方法，简称R-K方法，
其中正整数N 称为R-K方法的级，所有ci , ai , bij 都是待定 常数。
根据定义(7.2.7)，N 级R-K方法(7.9)的局部截断误差为
Rn+1 = y( t n+1 ) − y( t n ) − h∑ ci ki
dy 其斜率为 = f ( t0 , y0 ) dt ( t0 , y0 ) 由 点 斜 式 写 出 切线 方 程 dy y = y0 + ( t − t0 ) = y0 + ( t − t0 ) f ( t0 , y0 ) dt ( x0 , y )
0
等步长为h,则t1 - t0 = h, 可由切线算出 y1 ： 则 y1 = y0 + hf ( t0 , y0 ) 按此逐步计算y( tn ), 在tn +1处的值 ： yn+1 = yn + hf ( tn , yn ) y 注意： 这是“ 注意 ： 这是 “ 折 yN 线法” 而非“ 线法 ” 而非 “ 切 线法” 线法 ” 除第一个 点是曲线切线外， 点是曲线切线外 ， 其他点不是切线 y2 而是折线(如右 y1 y0 图所示)。 图所示 。
பைடு நூலகம்
则称数值解法(7.5)为显式方法。否则，称数值解法(7.3) 为隐式方法。
实验七 用Python解常微分方程
实验七用Python解常微分方程摘要本实验使用Python编程语言解决常微分方程（Ordinary Differential Equations, ODEs）问题。
常微分方程是描述自然界的变化和动态过程的数学工具，具有广泛的应用价值。
Python拥有强大的数值计算和绘图功能，因此是解决常微分方程问题的理想工具。
本实验将介绍常见的ODE求解方法，并给出相应的Python代码示例。
引言常微分方程描述了函数及其导数之间的关系，可以用来解决一系列实际问题，如物理、工程、经济等各个领域。
通过数值方法求解ODE，可以得到函数的近似解，从而对实际问题进行定性和定量分析。
ODE的数值解法常见的ODE求解方法包括欧拉法、改进欧拉法、龙格-库塔法等。
这些方法基于数值逼近的原理，通过迭代计算逼近函数的解。
以下是这些方法的简要介绍：1. 欧拉法（Euler's method）：将函数的导数展开为差商的形式，利用差商逼近导数，从而得到函数的近似解。
2. 改进欧拉法（Improved Euler's method）：在欧拉法的基础上引入中间点的逼近，进一步提高近似解的准确性。
3. 龙格-库塔法（Runge-Kutta method）：通过计算不同权重的斜率加权平均来逼近函数的导数，并利用逼近的导数进行迭代计算。
使用Python解常微分方程的示例代码以下是使用Python解常微分方程的示例代码，具体说明见代码注释：import numpy as npimport matplotlib.pyplot as plt定义常微分方程的函数形式def ode_func(t, y):return -y + np.sin(t)欧拉法求解常微分方程def euler_method(ode_func, t_start, t_end, y_init, step_size): num_steps = int((t_end - t_start) / step_size) + 1t = np.linspace(t_start, t_end, num_steps)y = np.zeros(num_steps)y[0] = y_initfor i in range(1, num_steps):y[i] = y[i-1] + step_size * ode_func(t[i-1], y[i-1])return t, y解常微分方程并绘制结果t_start, t_end = 0, 10y_init = 0step_size = 0.1t, y = euler_method(ode_func, t_start, t_end, y_init, step_size)plt.plot(t, y)plt.xlabel('t')plt.ylabel('y')plt.title('Solution of ODE using Euler\'s method')plt.grid(True)plt.show()结论本实验介绍了用Python解常微分方程的基本方法和示例代码，涉及了欧拉法。
第七章常微分方程数值解法
h2 h3 y ( xi 1 ) y ( xi h) y ( xi ) hy '( xi ) y ''( xi ) y '''( xi ) 2! 3!
丢掉高阶项，有
y( xi 1 ) y( xi h) y( xi ) hy '( xi ) yi hf ( xi , yi )
| f ( x, y1 ) f ( x, y2 ) | L | y1 y2 | ,
那么模型问题在 [ a, b] 存在唯一解。
Lipschitz 连续： | f ( x, y1 ) f ( x, y2 ) | L | y1 y2 | .
(1) 比连续性强: y1 y2 可推出 f ( x, y1 ) f ( x, y2 ) ; (2) 比连续的 1 阶导弱：具有连续的 1 阶导，则
f | f ( x, y1 ) f ( x, y2 ) || ( ) || y1 y2 | L | y1 y2 | . y
常微分方程数值解法
目标：计算出解析解 y ( x) 在一系列节点 a x0 x1 xn1 xn b 处的近似值 yi y( xi ) ，即所谓的数值解。节点间距 hi xi 1 xi ，一般 取为等距节点。
常微分方程初值问题的数值解法一般分为两大类： (1)单步法：在计算 yn 1 时，只用到前一步的值，即用到 xn1 , xn , yn ，则给定初
值之后，就可逐步计算。例如 Euler 法、向后欧拉法、梯形公式、龙格-库塔法；
(2) 多步法： 这 类 方 法 在 计算 yn 1 时 ， 除 了 用 到 xn1 , xn , yn 外 ， 还 要 用到
常微分方程的数值解法
常微分方程的数值解法1. 引言常微分方程是自变量只有一个的微分方程，广泛应用于自然科学、工程技术和社会科学等领域。
由于常微分方程的解析解不易得到或难以求得，数值解法成为解决常微分方程问题的重要手段之一。
本文将介绍几种常用的常微分方程的数值解法。
2. 欧拉方法欧拉方法是最简单的一种数值解法，其具体步骤如下：- 将自变量的区间等分为n个子区间；- 在每个子区间上假设解函数为线性函数，即通过给定的初始条件在每个子区间上构造切线；- 使用切线的斜率（即导数）逼近每个子区间上的解函数，并将其作为下一个子区间的初始条件；- 重复上述过程直至达到所需的精度。
3. 改进的欧拉方法改进的欧拉方法是对欧拉方法的一种改进，主要思想是利用两个切线的斜率的平均值来逼近每个子区间上的解函数。
具体步骤如下： - 将自变量的区间等分为n个子区间；- 在每个子区间上构造两个切线，分别通过给定的初始条件和通过欧拉方法得到的下一个初始条件；- 取两个切线的斜率的平均值，将其作为该子区间上解函数的斜率，并计算下一个子区间的初始条件；- 重复上述过程直至达到所需的精度。
4. 二阶龙格-库塔方法二阶龙格-库塔方法是一种更为精确的数值解法，其基本思想是通过近似计算解函数在每个子区间上的平均斜率。
具体步骤如下： - 将自变量的区间等分为n个子区间；- 在每个子区间上计算解函数的斜率，并以该斜率的平均值近似表示该子区间上解函数的斜率；- 利用该斜率近似值计算下一个子区间的初始条件，并进一步逼近解函数；- 重复上述过程直至达到所需的精度。
5. 龙格-库塔法（四阶）龙格-库塔法是目前常用的数值解法之一，其精度较高。
四阶龙格-库塔法是其中较为常用的一种，其具体步骤如下：- 将自变量的区间等分为n个子区间；- 在每个子区间上进行多次迭代计算，得到该子区间上解函数的近似值；- 利用近似值计算每个子区间上的斜率，并以其加权平均值逼近解函数的斜率；- 计算下一个子区间的初始条件，并进一步逼近解函数；- 重复上述过程直至达到所需的精度。
实验七 用Maple解常微分方程
实验七用Maple解常微分方程1. 实验目的本实验旨在通过使用数学建模软件Maple来解常微分方程，加深对常微分方程解法的认识和理解。
通过实际操作和观察结果，提高对Maple软件的运用能力。
2. 实验原理常微分方程是描述物理、化学、工程等领域中的连续变化过程的常见数学工具。
解常微分方程可以帮助我们理解系统的演化规律，从而进行预测和控制。
Maple是一款强大的数学软件，其中包含了丰富的求解常微分方程的函数。
通过输入常微分方程的表达式，Maple可以直接给出解析解或数值解。
在本实验中，我们将使用Maple来解常微分方程。
3. 实验步骤3.1 安装Maple软件3.2 打开Maple软件双击桌面上的Maple图标，打开软件。
3.3 输入常微分方程点击菜单栏中的"输入"，选择"数学输入"，在弹出的对话框中输入常微分方程的表达式。
例如，我们要解的方程是一阶线性常微分方程`dy/dx + y = 0`，则输入表达式为：diff(y(x),x) + y(x) = 03.4 求解方程点击菜单栏中的"执行"，选择"执行工作表"，Maple将根据输入的方程进行求解。
3.5 查看解析解或数值解Maple会给出方程的解析解或数值解。
根据实验需求，可以选择相应的解进行查看和分析。
3.6 导出结果点击菜单栏中的"文件"，选择"导出为"，选择导出格式和保存路径，点击"保存"，将结果导出为文档或图像文件。
4. 实验结果根据实验中输入的常微分方程，Maple求解得到如下解析解：y(x) = C exp(-x)其中C为任意常数。
5. 实验总结通过本次实验，我们研究了使用Maple软件求解常微分方程的方法。
Maple的强大功能和简便操作使得解常微分方程变得更加容易。
通过实际操作，我们可以深入理解常微分方程的解法和物理意义。
第7章 常微分方程数值解法
代入(6―3)式得
h yi 1 yi [ f ( xi , yi ) f ( xi 1 , yi 1 )] 2 i 0,1, 2, , n 1
(6―5)
这样得到的点列仍为一折线，只是用平均斜率 来代替原来一点处的斜率。式(6―5)称为改进的欧拉 公式。
不难发现，欧拉公式(6―3)是关于yi+1 的显式，只
h y xi 1 yi 1 f xi 1 , y xi 1 f xi 1 , yi 1 2 (6―15) h 3 '' f 12
因此
hL h3 y ( xi 1 ) yi 1 y ( xi 1 ) yi 1 f ( ) 2 12 h3 (1 q) y ( xi 1 ) yi 1 f ( ) 12 y ( xi 1 ) yi 1 O ( h 3 )
c 并取 yi 1 yi(1)
(6―7)
虽然式(6―7)仅迭代一次，但因进行了预先估计，
故精度却有较大的提高。 在实际计算时，还常常将式(6―7)写成下列形式：
k1 f ( xi , yi ) k f ( x h, y hk ) i i 1 2 h yi 1 yi 2 (k1 k2 ) i 0,1, 2,
在进行误差分析时，我们假设yi=y(xi)，考虑用
yi+1 代替y(x
i+1)而产生截断误差，确定欧拉公式和改
进的欧拉公式的精确度。 设初值问题(6―1)的准确解为y=y(x)，则利用泰 勒公式
y ( xi 1 ) y ( xi h ) h2 y ( xi ) hy ( xi ) y ( xi ) 2! h3 y ( xi ) 3!
常微分方程数值解法
第七章 常微分方程数值解法常微分方程中只有一些典型方程能求出初等解(用初等函数表示的解)，大部分的方程是求不出初等解的。
另外，有些初值问题虽然有初等解，但由于形式太复杂不便于应用。
因此，有必要探讨常微分方程初值问题的数值解法。
本章主要介绍一阶常微分方程初值问题的欧拉法、龙格-库塔法、阿达姆斯方法，在此基础上推出一阶微分方程组与高阶方程初值问题的 数值解法；此外，还将简要介绍求解二阶常微分方程值问题的差分方法、试射法。
第一节 欧拉法求解常微分方程初值问题⎪⎩⎪⎨⎧==00)(),(y x y y x f dxdy（1）的数值解，就是寻求准确解)(x y 在一系列离散节点 <<<<<n x x x x 210 上的近似值 ,,,,,210n y y y y{}n y 称为问题的数值解，数值解所满足的离散方程统称为差分格式，1--=i i ix x h 称为步长，实用中常取定步长。
显然，只有当初值问题(1)的解存在且唯一时，使用数值解法才有意义，这一前提条件由下 面定理保证。
定理 设函数()y x f ,在区域+∞≤≤-∞≤≤y b x a D ,:上连续，且在区域D 内满足李普希兹(Lipschitz)条件，即存在正数L ，使得对于R 内任意两点()1,y x 与()2,y x ,恒有()()2121,,y y L y x f y x f -≤-则初值问题(1)的解()x y 存在并且唯一。
一、欧拉法(欧拉折线法)若将函数)xy 在点nx处的导数()n x y '用两点式代替， 即()hx y x y x y n n n )()(1-≈'+，再用n y 近似地代替()n x y ,则初值问题(1)变为⎩⎨⎧==++=+ ,2,1,0),()(001n x y y y x hf y y n n n n(2)(2)式就是著名的欧拉(Euler)公式。
以上方法称 为欧 拉法或欧拉折线法。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
,试试编程利用欧拉公式画出 的速度曲线(取).
实验算例
1. 欧拉法 利用matlab中的一个循环语句即可实现欧拉法中的利用第一个节点计
算随后所有节点的工作. function Euler %输出节点的x值和y值, clc %------输入(x0,xn)为求解区间,y0=y(x0)为初始条件,n为区间的 等分个数------x0=0; y0=1; xn=0.5; n=50; %-------------------------------------------------------------------y(1)=y0; h=(xn-x0)/n; x=x0:h:xn; disp('欧拉法结果如下：') disp(['y','(',num2str(x0),')','=',num2str(y0)]); for i=2:n+1 %利用循环语句实现欧拉法利用第一个节 点的函数值计算 y(i)=y(i-1)+h*f(x0,y(i-1)); %随后所有节点函数值的过程. x0=x0+h; disp(['y','(',num2str(x0),')','=',num2str(y(i))]); end plot(x,y,'ro-'); %利用求出的x坐标和y坐标画出解的 近似图形. xlabel('x') ylabel('y') title('欧拉法求出的折线') hold on %-------若原微分方程的理论解能求出，则可画出积分曲线比 较--------------x=0:0.01:0.5; y=x+exp(-x); plot(x,y,'b-') hold off %-----------------------------------------------------------------function y=f(x,z) y=-z+x+1;
欧拉法和改进欧拉法的效果比较图为
实验习题
1. 利用欧拉法程序,求解初值问题 （注意：若理论解未知，解曲线画不出来，应删掉相关程序语句） (1) (2) 并将欧拉法求出的折线和改进欧拉法求出的折线放在同一个图中，比较 效果。 2. 某跳伞者在 时刻从飞机上跳出,假设初始时刻的垂直速度为0,且跳伞者垂直下落.已 知空气阻力为 ,其中c为常数,v为垂直速度,向下方向为正.可求得跳伞者的速度满足微分 方程方程 ,g为重力加速度.若此跳伞者质量为 ,且已知
实验七 常微分方程数值解法
实验目的
1、 掌握常微分方程数值解法中欧拉法和改进欧拉法的基本原理 和算法； 2、 培养Matlab数学软件的编程与上机调试能力.
实验课时
4课时
实验准备
初值问题 的解析解法在常微分方程中有介绍,但是实际中遇到的常微分方程一般 很难得到解析解,所以我们一般用数值方法求出它的数值解. 数值解法的 结果是自变量在一系列离散点处的函数的近似值,一般情况下取步长为 定数. 1.欧拉法 利用一阶微分公式 可得差分方程 其中, , ,则得到初值问题的数值解. 2.后退欧拉法 利用一阶微分公式 可得差分方程 此为后退欧拉法. 4. 梯形求积公式 利用数值积分 可得梯形公式 5. 改进欧拉公式 一般来讲,隐式格式要比显示格式具有较好的数值稳定性,所以常常被 采用,但是为了避免隐式格式的计算量,一般采用预估-校正的方法.对于梯 形公式,用欧拉公式预估再用梯形公式校正,可得改进欧拉公式:
y(1)=y0; h=(xn-x0)/n; x=x0:h:xn; disp('改进欧拉法结果如下：') disp(['y','(',num2str(x0),')','=',num2str(y0)]); for i=2:n+1 %利用循环语句实现欧拉法利用第一个节点的 函数值计算 k(1)=f(x0,y(i-1)); k(2)=f(x0,y(i-1)+h*k(1)); y(i)=y(i-1)+h*(1/2*k(1)+1/2*k(2)); %随后所有节点函数值的过 程. x0=x0+h; disp(['y','(',num2str(x0),')','=',num2str(y(i))]); end plot(x,y,'k*-'); %利用求出的x坐标和y坐标画出解的近 似图形. xlabel('x') ylabel('y') title('改进欧拉法求出的折线') hold on %-------若原微分方程的理论解能求出，则可画出积分曲线比较-------------x=0:0.01:0.5; y=x+exp(-x); plot(x,y,'b-') %-----------------------------------------------------------------legend('改进欧拉法求出的折线','解曲线') hold off function y=f(x,z) y=-z+x+1; 对于例1的结果为 输出图像为
例1. 用欧拉法求解初值问题, 并画出解图形. >> 输出图像为
另外,还可以改变输入的参数n,加密节点,得到更精确的数值解. 2. 改进欧拉法 function Gaijin_Euler %输出节点的x值和y值, clc close all %------输入(x0,xn)为求解区间,y0=y(x0)为初始条件,n为区间的等 分个数------x0=0; y0=1; xn=0.5; n=50; %--------------------------------------------------------------------