龙格库塔方法的Miline-Hamming预测-校正算法实验报告知识讲解
龙格库塔算法
龙格库塔算法龙格库塔算法（Runge-Kutta method）是一种常用的数值解微分方程的方法，其基本原理是通过逐步逼近的方式，根据初始条件和微分方程的表达式，计算出方程的近似解。
该方法具有较高的精度和稳定性，在科学计算、物理模拟、工程建模等领域得到广泛应用。
龙格库塔算法的核心思想是将微分方程的解按照一定的步长进行离散化，从而将连续的求解问题转化为离散的迭代过程。
具体来说，龙格库塔算法通过计算函数在一定步长内的平均斜率，来估计下一个点的函数值。
这个平均斜率是通过多次计算函数在不同点上的导数得到的，从而提高了计算的精度。
龙格库塔算法的一般形式可以表示为：k1 = f(tn, yn)k2 = f(tn + h/2, yn + h/2 * k1)k3 = f(tn + h/2, yn + h/2 * k2)k4 = f(tn + h, yn + h * k3)yn+1 = yn + h/6 * (k1 + 2k2 + 2k3 + k4)其中，tn是当前时间点，yn是当前函数值，h是步长，f是微分方程的表达式。
通过多次迭代，可以逐渐逼近微分方程的解。
龙格库塔算法的优点在于其精确度较高，可以通过调整步长来控制计算的精度和效率。
此外，该算法具有较好的数值稳定性，可以有效处理非线性、刚性或高阶微分方程等复杂问题。
因此，在科学和工程计算中，龙格库塔算法被广泛应用于各种数值模拟和求解问题。
需要注意的是，龙格库塔算法并非万能的，对于一些特殊的问题，可能存在数值不稳定性或计算精度不够的情况。
此外，算法的步长选择也需要根据具体问题进行调整，过小的步长会增加计算量，而过大的步长可能导致精度下降。
因此，在使用龙格库塔算法时，需要根据具体问题的特点和要求来选择合适的步长和算法参数，以获得满意的计算结果。
总结起来，龙格库塔算法是一种常用的数值解微分方程的方法，具有较高的精度和稳定性。
通过离散化和迭代的方式，可以逐步逼近微分方程的解。
MATLAB中龙格库塔（RUNGEKUTTA）方法原理及实现
MATLAB中龙格库塔（RUNGEKUTTA）方法原理及实现函数功能ode是专门用于解微分方程的功能函数，他有ode23,ode45,ode23s等等，采用的是Runge-Kutta算法。
ode45表示采用四阶，五阶runge-kutta单步算法,截断误差为(Δx)3。
解决的是Nonstiff(非刚性)的常微分方程.是解决数值解问题的首选方法,若长时间没结果,应该就是刚性的,换用ode23来解.使用方法[T,Y]=ode45(odefun,tspan,y0)odefun是函数句柄,可以是函数文件名,匿名函数句柄或内联函数名tspan是区间[t0tf]或者一系列散点[t0,t1,...,tf]y0是初始值向量T返回列向量的时间点Y返回对应T的求解列向量[T,Y]=ode45(odefun,tspan,y0,options)options是求解参数设置,可以用odeset在计算前设定误差,输出参数,事件等[T,Y,TE,YE,IE]=ode45(odefun,tspan,y0,options)在设置了事件参数后的对应输出TE事件发生时间YE事件解决时间IE事件消失时间sol=ode45(odefun,[t0tf],y0...)sol结构体输出结果应用举例1求解一阶常微分方程程序:一阶常微分方程odefun=@(t,y)(y+3*t)/t^2;%定义函数tspan=[14];%求解区间y0=-2;%初值[t,y]=ode45(odefun,tspan,y0);plot(t,y)%作图title('t^2y''=y+3t,y(1)=-2,1<t<4')< p="">legend('t^2y''=y+3t')xlabel('t')ylabel('y')%精确解%dsolve('t^2*Dy=y+3*t','y(1)=-2')%ans=一阶求解结果图%(3*Ei(1)-2*exp(1))/exp(1/t)-(3*Ei(1/t))/exp(1/t)2求解高阶常微分方程关键是将高阶转为一阶,odefun的书写.F(y,y',y''...y(n-1),t)=0用变量替换,y1=y,y2=y'...注意odefun方程定义为列向量dxdy=[y(1),y(2)....]程序:function Testode45tspan=[3.94.0];%求解区间y0=[28];%初值[t,x]=ode45(@odefun,tspan,y0);plot(t,x(:,1),'-o',t,x(:,2),'-*')legend('y1','y2')title('y''''=-t*y+e^t*y''+3sin2t')xlabel('t')ylabel('y')function y=odefun(t,x)y=zeros(2,1);%列向量y(1)=x(2);y(2)=-t*x(1)+exp(t)*x(2)+3*sin(2*t);endend高阶求解结果图相关函数ode23,ode45,ode113,ode15s,ode23s,ode23t,ode23tbMatlab中龙格-库塔(Runge-Kutta)方法原理及实现龙格-库塔(Runge-Kutta)方法是一种在工程上应用广泛的高精度单步算法。
龙格库塔实验报告
一、实验背景常微分方程（ODE）在自然科学、工程技术等领域中具有广泛的应用。
然而，许多微分方程无法得到精确解析解，因此需要借助数值方法进行求解。
龙格-库塔（Runge-Kutta）方法是一种常用的数值求解常微分方程的方法，具有精度高、稳定性好等优点。
本实验旨在通过编写程序，实现四阶龙格-库塔方法，并验证其在求解常微分方程中的有效性和准确性。
二、实验目的1. 理解四阶龙格-库塔方法的基本原理和计算步骤。
2. 编写程序实现四阶龙格-库塔方法。
3. 选取典型常微分方程，验证四阶龙格-库塔方法的求解精度和稳定性。
三、实验原理四阶龙格-库塔方法是一种基于泰勒级数展开的数值方法，其基本思想是将微分方程的解在某个区间内进行近似，并通过迭代计算得到近似解。
具体步骤如下：1. 初始化：给定初始条件y0，步长h，求解区间[a, b]。
2. 迭代计算：对于k=1, 2, ..., n（n为迭代次数），- 计算k1 = f(xk-1, yk-1)（f为微分方程的右端函数）；- 计算k2 = f(xk-1 + h/2, yk-1 + h/2 k1)；- 计算k3 = f(xk-1 + h/2, yk-1 + h/2 k2)；- 计算k4 = f(xk-1 + h, yk-1 + h k3)；- 更新y值：yk = yk-1 + (h/6) (k1 + 2k2 + 2k3 + k4)；- 更新x值：xk = xk-1 + h；3. 输出结果：输出最终的近似解y(n)。
四、实验步骤1. 编写程序实现四阶龙格-库塔方法。
2. 选取典型常微分方程，如：- y' = -y，初始条件y(0) = 1，求解区间[0, 2π]；- y' = y^2，初始条件y(0) = 1，求解区间[0, 1]。
3. 对每个常微分方程，设置不同的步长h和迭代次数n，分别计算近似解y(n)。
4. 将计算得到的近似解与解析解进行比较，分析四阶龙格-库塔方法的精度和稳定性。
龙格库塔法介绍
h 0.005稳定.
2) 改进欧拉法(预测 — 校正，即二阶R K法)：
yn
1
yn
h
k1 2
k2 2
,
k1
f (xn, yn ) yn,
k2 f (xn h, yn hk1) ( yn hyn )，
即yn1 故
yn
当初值准确即e0 0时，整体误差为en O(h p ).
证明
考察单步法的收敛性归结为验证增量函数(x, y,h)是否
满足Lipschitz条件.
欧拉法 : (x, y,h) f (x, y)，L L.
改进的欧拉法：
yn1
yn
h[ 2
f
(xn, yn )
f
( xn 1,
xn
f
( x,
y ( x))dx
r
h ci
i 1
f
( xn
ih,
y ( xn
ih)).
yn1 yn h(xn, yn, h),
(3.4)
其中
r
(xn, yn, h) ciki ,
(3.5)
i 1
k1 f (xn, yn )，
欧拉法r 1, p 1.改进 欧拉法r 2, p 2.
k1)
k3）
k3 f (xn h, yn hk1 2hk2 )
称为库塔三阶方法.
阶数p和段数r（计算函数值次数）的关系
r12 p1 2
3 4 5 6 7 r≥8 3 4 4 5 6 r-2
常用的经典四阶龙格 库塔方法:
龙格库塔方法理解及应用
龙格库塔方法理解及应用龙格库塔方法是一种常用的数值解微分方程的方法，也是许多科学、工程和经济领域中常用的算法之一。
本文将介绍龙格库塔方法的原理及其应用。
一、龙格库塔方法原理在计算微分方程时，往往需要对方程进行离散化，采用数值方法处理。
龙格库塔方法（Runge-Kutta method）就是一种离散化的数值方法，其原理可以概括为：通过相应的递推公式，将微分方程在离散时间点上进行逼近，从而得到近似的解。
具体来说，假设要求解如下形式的一阶常微分方程：$$ y'=f(t,y) $$其中，$f(t,y)$是已知的函数，$y(t)$是未知函数，并且已知初值$y(t_0)=y_0$。
为了离散化这个方程，我们可以采用以下的递推公式：$$ \begin{aligned} y_1 &=y_0 + h\varphi_1 \\ y_2 &=y_0 +h\varphi_2 \\ \cdots &=\cdots \\ y_n &=y_0 + h\varphi_n \end{aligned} $$其中，$h$是离散时间点的时间步长，$t_n=t_0+nh$，$\varphi_i$是与$t_i$有关的递推公式。
根据龙格库塔方法的不同级别，$\varphi_i$也有不同的形式。
二、龙格库塔方法的应用由于龙格库塔方法的较高精度和鲁棒性，以及易于实现等特点，它在各个领域都有着广泛的应用。
1. 数学领域在数学领域，龙格库塔方法可以用于求解常微分方程、偏微分方程、常微分方程组等等，特别是对于复杂的高阶微分方程，龙格库塔方法更是可以发挥其优势。
2. 物理学领域在物理学领域，各种微分方程是研究物理过程的基础。
龙格库塔方法在求解各种物理问题时也得到了广泛的应用，如天体力学、流体力学、电磁场问题等等。
3. 经济学领域在经济学领域，许多经济问题可以通过微分方程的形式进行建模，并采用龙格库塔方法进行数值求解。
龙格库塔法解轨道参数
龙格库塔法解轨道参数龙格库塔法（Runge-Kutta method）是一种常用的数值解微分方程的方法，可以用于解决轨道参数的计算问题。
轨道参数是描述天体在空间中运动的重要参数，包括轨道半径、轨道倾角、轨道偏心率等。
龙格库塔法的基本思想是通过一系列的近似计算，逐步逼近微分方程的解。
该方法通常适用于一阶常微分方程，但也可以通过将高阶微分方程转化为一阶形式来应用。
在解轨道参数的问题中，我们可以将轨道的运动方程转化为一阶微分方程组。
以行星绕太阳运动为例，我们可以将行星的轨道运动描述为：$frac{dx}{dt} = v_x$$frac{dy}{dt} = v_y$$frac{dv_x}{dt} = frac{-GMx}{r^3}$$frac{dv_y}{dt} = frac{-GMy}{r^3}$其中，$x$和$y$分别表示行星在直角坐标系中的位置，$v_x$和$v_y$表示行星在$x$和$y$方向上的速度，$r$表示行星到太阳的距离，$G$是引力常数，$M$是太阳的质量。
接下来，我们可以使用龙格库塔法来逐步计算行星的位置和速度。
首先，我们可以选择一个初始的位置和速度，然后根据上述微分方程组进行迭代计算。
每一步的计算都依赖于前一步的结果，通过不断迭代可以得到轨道参数的数值解。
龙格库塔法的主要步骤包括计算$k_1$、$k_2$、$k_3$和$k_4$：$k_1 = h f(t_n, y_n)$$k_2 = h f(t_n + frac{h}{2}, y_n + frac{k_1}{2})$$k_3 = h f(t_n + frac{h}{2}, y_n + frac{k_2}{2})$$k_4 = h f(t_n + h, y_n + k_3)$其中，$h$表示步长，$t_n$表示当前时间，$y_n$表示当前位置和速度。
$f(t_n, y_n)$表示微分方程组的右端项。
通过上述计算得到的$k_1$、$k_2$、$k_3$和$k_4$，可以计算下一步的位置和速度：$y_{n+1} = y_n + frac{1}{6}(k_1 + 2k_2 + 2k_3 + k_4)$通过不断重复上述步骤，我们可以得到行星轨道的数值解。
滤波算法 龙格库塔算法-概述说明以及解释
滤波算法龙格库塔算法-概述说明以及解释1.引言1.1 概述概述：滤波算法和龙格库塔算法是计算机科学领域中常用的算法之一，它们在数据处理和数值计算中有着重要的应用价值。
滤波算法被广泛应用于信号处理、图像处理、通信系统等领域，用于消除信号中的噪声和提高数据的质量。
而龙格库塔算法则是一种常用的数值求解微分方程的方法，能够有效地对复杂的数学模型进行数值求解，具有较高的准确性和稳定性。
本文将分别介绍滤波算法和龙格库塔算法的原理、优缺点以及应用领域，希望读者通过本文能够对这两种算法有更深入的了解，并在实际应用中能够灵活运用。
1.2 文章结构本文将分为四个部分来探讨滤波算法和龙格库塔算法。
首先在引言部分，对滤波算法和龙格库塔算法进行简要介绍，并说明本文的结构和目的。
接着在第二部分，详细介绍滤波算法的概念、常见算法和应用场景，以便读者对滤波算法有个全面的了解。
然后在第三部分，深入探讨龙格库塔算法的简介、原理和优缺点，帮助读者更好地理解这一种数值计算方法。
最后，在结论部分对两种算法进行总结，并展望未来可能的发展方向，以及得出结论。
通过以上四个部分的内容，读者能够全面了解和掌握滤波算法和龙格库塔算法的相关知识。
1.3 目的本文的主要目的是介绍和探讨滤波算法和龙格库塔算法这两种在计算机科学和工程领域中广泛应用的算法。
通过对这两种算法的概述、原理和应用进行详细分析，能够帮助读者全面了解它们的工作原理和特点。
同时，通过对这两种算法的比较和讨论，可以帮助读者更好地理解它们在不同应用场景下的适用性和优缺点。
此外，本文还旨在为读者提供一个深入学习和掌握这两种算法的基础知识和入门指南。
通过本文的学习，读者可以加深对滤波算法和龙格库塔算法的理解，为进一步的研究和实践打下坚实的基础。
同时，希望本文能够激发读者对算法领域的兴趣，促使他们深入研究和探索更多先进的算法及其应用。
2.滤波算法2.1 滤波算法概述滤波算法是一种用于处理信号或数据的技术，其主要目的是通过去除噪声或不需要的信息，从而提取出所需的信号或数据。
龙格—库塔实验报告
龙格—库塔法解常微分方程实验报告一、实验题目求解初值问题⎪⎩⎪⎨⎧=<<-=1)0()10(2'y x y x y y二、实验引言1、实验目的进一步理解龙格—库塔方法的设计思路和算法流程，培养动手实践能力和分析能力。
2、实验意义龙格—库塔方法的推导基于泰勒展开方法，因而它要求的解具有较好的光滑性质。
反之，如果解得光滑性差，那么，使用四阶龙格—库塔法方法求得的数值解，其精度可能反而不如梯形方法。
实际计算中，应当针对问题的具体特点选择合适的算法。
三、算法设计1. 基本思想由Lagrange 微分中值定理，11()()'()()()(,())n n n n n y x y x y x x y x hf y ξξξ++=+-=+记*(,())k hf y ξξ=，则得到*1()()n n y x y x k +=+这样，给出*k 的一种算法，就得到求解微分方程初值问题的一种计算公式。
四阶龙格_库塔法是用1k ，2k ，3k 和4k 的加权平均值来近似*k 。
最经典的四阶龙格—库塔公式为：121324311234(,)(,)22(,)22(,)y (22)2n n n n n nn n n n K f x y h h K f x y K h h K f x y K K f x h y hK h y K K K K +=⎧⎪⎪=++⎪⎪⎪=++⎨⎪=++⎪⎪⎪=++++⎪⎩四阶龙格—库塔法的误差估计局部截断误差为5()O h 。
2.算法流程图四、程序设计program longgekutaimplicit nonereal,parameter::b=1real::h=0.2integer::nreal::x,K1,K2,K3,K4,yreal,external::fx=0y=1open (unit=10,file='1.txt')do while(x<=b)K1=f(x,y)K2=f(x+h/2,y+K1*h/2)K3=f(x+h/2,y+K2*h/2)K4=F(x+h,y+K3*h)y=y+(k1+2*K2+2*K3+K4)*h/6 x=x+hwrite(10,*) x,yend doendfunction f(x,y)implicit nonereal::f,x,yf=y-2*x/yend function五、结果及讨论1.实验结果0.2000000 1.183229 0.4000000 1.341667 0.6000000 1.4832810.8000000 1.6125141.000000 1.7321421.200000 1.844040六、算法评价1、本次实验实现了常微分方程初值问题数值解法中的四阶龙格—库塔法2、对欧拉法和龙格—库塔法进行比较：在相同步长的情况下，欧拉法每步只计算一个函数值，四阶龙格—库塔法每步需计算四个函数值，就是说，四阶龙格—库塔法的计算量差不多是欧拉法的四倍，为了比较它们的计算精度，可以将欧拉法的步长取为h,将四阶龙格—库塔法的步长取为4h。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
龙格库塔方法的
M i l i n e-H a m m i n g预测-校正算法实验报告
2011-2012学年第2学期实验报告
实验名称：微分方程数值解实验学院：******
专业：**************
班级：**********
班内序号：**
学号：********
姓名：******
任课教师：******
北京邮电大学
时间：****年**月**日
实验目标
用多环节Miline-Hamming 预测-校正算法求下列方程的解
{y‘=y −2x
y ,y (0)=1, 0≤x ≤4 其中解析解为 y (x )=√1+2x
实验原理
计算龙格库塔显示公式计算预测值，然后用隐式公式进行校正。
Miline-Hamming 预测-校正公式为
{
p n+1=u n−3+4
3
h(2f n −f n−1+f n−2)
m n+1=p n+1+112
121(c n
−p n )
c n+1=18
(9u n −u n−2)+38h[f (t n+1,m n+1)+2f n −f n−1]
u n+1=c n+1−9(c n+1
−p n+1)
其对应的算法流程为
1） 输入a ，b ，f(t,u)，N ，u 0 2） 置h=(b-a)/N ，t 0=a ，n=1 3） 计算f n-1=f(t n-1,u n-1)
K 1=hf n-1
K 2=hf(t n-1+h/2, u n-1+K 1/2) K 3=hf(t n-1+h/2, u n-1+K 2/2) K 4=hf(t n-1+h, u n-1+K 3)
u n = u n-1+1/6(K 1 +2K 2 +2K 3 +K 4) t n =a+nh
4） 输出(t n ,u n )
5） 若n<3，置n+1→n ，返回3；
否则，置n+1→n ，0→p 0，0→c 0，转6.
6） 计算
f 3=f(t 3,u 3) t= t 3+h
p=u 0+4/3(2 f 3 –f 2 +2f 1) m=p+112/121(c 0-p 0)
c=1/8(9u 3- u 1)+3/8h[f(t,m)+ 2 f 3 –f 2] u=c-9/121(c-p)
输出(t,u)
7）如n<N，则置n+1→n，t j+1→ t j，u j+1→u j，f j+1→f j(j=0,1,2)，t= t3，u= u3，p= p0，c= c0，转6
否则停止。
实验过程
我们不妨设步长h=0.2，编程实现如下：
clear
clf
clc
%直接求解微分方程
y=dsolve('Dy=y-2*t/y','y(0)=1','t');
%Miline-Hamming预测-校正法
h=0.2;
t=0:h:4;
n=length(t);
u=zeros(1,n);
u(1)=1;
zbu(1,1)=t(1);
zbu(2,1)=u(1);
f=zeros(1,n);
p=zeros(1,n);
c=zeros(1,n);
m=zeros(1,n);
for i=2:n
if i-1<=3
f(i-1)=u(i-1)-2*t(i-1)/u(i-1);
k1=h*f(i-1);
k2=h*((u(i-1)+k1/2)-2*(t(i-1)+h/2)/(u(i-1)+k1/2)); k3=h*((u(i-1)+k2/2)-2*(t(i-1)+h/2)/(u(i-1)+k2/2)); k4=h*((u(i-1)+k3)-2*(t(i-1)+h)/(u(i-1)+k3));
u(i)=u(i-1)+1/6*(k1+2*k2+2*k3+k4);
zbu(1,i)=t(i);
zbu(2,i)=u(i);
else
f(i-1)=u(i-1)-2*t(i-1)/u(i-1);
p(i)=u(i-4)+4/3*h*(2*f(i-1)-f(i-2)+2*f(i-3));
m(i)=p(i)+112/121*(c(i-1)-p(i-1));
c(i)=1/8*(9*u(i-1)-u(i-3))+3/8*h*(m(i)-2*t(i)/m(i)+2*f(i-1)-f(i-2));
u(i)=c(i)-9/121*(c(i)-p(i));
zbu(1,i)=t(i);
zbu(2,i)=u(i);
end
end
zbu
%作图
plot(t,u,'r*','markersize',10)
hold on,
ezplot(y,[0,4])
hold on,
title('Miline-HammingÔ¤²â-УÕýËã·¨')
grid on
legend('Miline-HammingÔ¤²â-УÕýËã·¨','½âÎö½â')
%解真值
h=0.2;
t=0:h:4;
n=length(t);
for i=1:n
y(i)=(1+2*t(i))^(1/2);
zby(1,i)=t(i);
zby(2,i)=y(i);
end
zby
我们可以得到计算后的结果图像如图1所示
图1 Miline-Hamming预测-校正法与解析解比较（h=0.2）同时我们可以得到Miline-Hamming预测-校正法和解析解在各点处的数值分别如下表1所示：
表1 Miline-Hamming预测-校正法与解析解在各点数值比较（h=0.1）
为了评判Miline-Hamming预测-校正法的算法精度，在这里我们利用相对误差的概念进行评判。
对于Miline-Hamming预测-校正法的每个的估计值有：
相对误差=|估计值－真值|
真值
从而我们可以通过计算得到如下的相对误差表：
t坐00.20000.40000.60000.8000 1.0000 1.2000 1.4000 1.6000
t坐 1.8000 2.0000 2.2000 2.4000 2.6000 2.8000 3.0000 3.2000 3.4000
很明显，当我们对各点处的相对误差取平均后，该平均值小于0.01。
因
此，我们可以认为Miline-Hamming预测-校正法的在h=0.2时的算法精度相对较高，所得到的结果与真值较为接近。
接下来我们在对h=0.5和h=0.1的情况进行计算，可以得到结果如下表3和表4所示
00.5000 1.0000 1.5000 2.0000 2.5000 3.0000 3.5000 4.0000 t坐
标
H=.5 1.0000 1.4155 1.7355 2.0084 2.2517 2.4886 2.7454 3.0750 3.5964 真 1.0000 1.4142 1.7321 2.0000 2.2361 2.4495 2.6458 2.8284 3.0000
t坐00.10000.20000.30000.40000.50000.60000.70000.8000
值
表4 Miline-Hamming预测-校正法在各点值及相对误差比较（h=0.1）
其中表3中相对误差的平均值为0.0393。
而表4中的误差值小于是10的-3次方，在此不列出。
可以的到结果图像如下图2和图3所示：
图1 Miline-Hamming预测-校正法与解析解比较（h=0.5）
图1 Miline-Hamming预测-校正法与解析解比较（h=0.1）
实验结果
通过以上计算，我们可以得到如下的结论：
1．Miline-Hamming预测-校正法的计算精度相对较高，，当步长h=0.2时平均相对误差已小于0.01，因此可以认为这种方法可以得到和解析解较为接近的数值解。
2．伴随着步长的增加，我们可以发现相对误差的平均值随之减小。
因此我们认为当步长h越小时计算精度越高。
因此在计算能力允许的范围内，选取步长越小可以得到更加精确的结果。
3．在利用Miline-Hamming预测-校正法的过程中，前4次迭代的结果会对第五轮求得的数值产生影响。
因此，一旦前四轮轮迭代所得的结果有偏差，下一轮结果的偏差将大于之前的偏差。
因此会导致伴随迭代次数的增加而产生更大的偏差。