龙格库塔方法的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（1、广东技术师范学院天河学院，广东广州5105402、广东环境保护工程职业学院，广东佛山528200）在许多科学研究和工程测量中，根据实际问题所建立的数学模型大多数是和微分方程有关系的，而并不是所有的微分方程都可以求出其精确解，相当大一部分微分方程要求出其精确解是极其困难，甚至是不可能的。
因此，在求解微分方程时，我们常常用近似的方法求的数值解，而求数值解的方法除去一些特殊类型的方程，其他的大多数采用的是利用计算机的计算能力结合数值分析方法进行求解。
研究和选择好的数值计算方法在某种意义上来说比提高计算机计算速度更为重要，本文就龙格-库塔方法谈一点自己的想法。
1龙格-库塔方法的基本思想根据著名的欧拉公式，很显然这是一个显式的差分方程，可以改写成，用它计算y n+1，需要计算一次f （x ，y ），若假设y （x n ）=y n ，那么y n+1的表达式与y （x n+1）的泰勒展开式的前两项完全相同，也就是说局部截断误差为O （h 2）。
因此，我们可以得到龙格-库塔方法的基本思想：设法计算f （x ，y ）在某些点上的函数值，然后对这些函数值作线性组合，构造近似计算公式，最后在把近似公式与解的泰勒展开式进行比较，使得前面若干项完全一样，就得到了一定精度（由误差决定）的数值计算公式。
一般的显式的龙格-库塔方法的形式具体到构造时，先引进若干参数：其中，参数ωi ，αi ，βij 既与f （x,y ）无关，又与步数无关的常数，选择的原则是要使得最后求出来的近似解与解的泰勒展开式进行比较，有更多项吻合。
通常用的较多的是r 取2，3，4，同样为了计算的方便，可适当选择ωi ，αi ，βij 。
显然，取法有很多，都可以使相应的计算公式的局部截断误差到达一定精度，但是不能再进一步的提高局部截断误差的阶。
也就是说在计算两次函数值的情况下，局部截断误差的最高阶就是O （h 3），要想进一步提高阶，必须通过继续增加计算函数值的次数的途径。
第3章03龙格-库塔方法
这样有计算公式：
yn1 yn h[(1 )k1 k2 )
k1 f (xn , yn ), k2 f (xn ph, yn phk1)
y(h) n1
1 16
y( xn1)
(h)
y2 n1
1 15
(h)
y2 n1
y(h) n1
误差事后估计
可以通过检查步长折半前后两次计算结果的偏差
y y ( h ) 2
(h)
n1
n1
来判断所选取的步长是否合适.
1) 对于给定的精度ε,如果δ> ε,则反复将步长折半 进行计算,直到δ< ε为止,这时取步长折半后的 新值作为结果;
k4
k3 f ( xn 2h, yn 21hk1 22hk2 ) f ( xn 3h, yn 31hk1 32hk2 33hk3)
四阶龙格－库塔公式每步要四次计算函数值，具有四阶精 度，局部截断误差是O(h5) .
四阶龙格－库塔格式
下列经典格式是其中常用的一种:
yn1 k2
yk
1.000000 1.005000 1.019025 1.041218 1.070800 1.107076 1.149404 1.197210 1.249975 1.307228 1.368514
y( xk ) yk
0.0 1.6×10－4 2.9×10－4 4.0×10－4 4.8×10－4 5.5×10－4 5.9×10－4 6.2×10－4 6.5×10－4 6.6×10－4 6.6×10－4
调用 2 次函数。有 2 阶精度。
龙格库塔实验报告
龙格库塔实验报告龙格库塔实验报告引言：龙格库塔法（Runge-Kutta method）是一种常用的数值求解常微分方程（ODE）的方法。
它是由德国数学家卡尔·龙格（Carl Runge）和马丁·威尔海姆·库塔（Martin Wilhelm Kutta）于19世纪末独立发展而来。
龙格库塔法通过将微分方程转化为一系列逼近值，从而实现对微分方程的数值求解。
本实验旨在通过对龙格库塔法的研究和实践，深入了解其原理和应用。
一、龙格库塔法的原理龙格库塔法的基本思想是将微分方程的解分解为一系列逼近值，然后通过逐步迭代的方式计算这些逼近值。
具体而言，龙格库塔法将求解微分方程的过程分为多个步骤，每个步骤都计算一个逼近值，并利用这些逼近值逐步逼近真实解。
二、龙格库塔法的步骤1. 确定微分方程和初始条件：首先，需要明确待求解的微分方程以及初始条件。
微分方程可以是一阶或高阶的，初始条件是方程在某一点上的已知值。
2. 确定步长：步长（或称为时间间隔）决定了逼近值的精度。
较小的步长能够提高逼近值的准确性，但也会增加计算量。
因此，需要根据具体问题的需求选择合适的步长。
3. 迭代计算：根据龙格库塔法的公式，进行逐步的迭代计算。
首先，根据初始条件计算出第一个逼近值。
然后，利用该逼近值和微分方程的导数计算出下一个逼近值。
重复这一过程，直到达到所需的迭代次数或满足精度要求。
4. 计算误差：为了评估逼近值的准确性，需要计算误差。
误差可以通过与解析解的比较得出，或者通过比较不同步长下的逼近值得出。
较小的误差表明逼近值较为准确。
三、龙格库塔法的应用龙格库塔法广泛应用于科学和工程领域，特别是在求解微分方程模型的数值计算中。
以下是一些常见的应用场景：1. 物理学：龙格库塔法可以用于模拟物体在重力场中的运动，如自由落体、抛体运动等。
通过求解微分方程，可以得到物体的位置、速度等随时间的变化规律。
2. 经济学：龙格库塔法可以用于经济学模型的求解，如经济增长模型、投资模型等。
龙格库塔法-原理及程序实现
龙格库塔法一、基本原理：“龙格-库塔(Runge-Kutta)方法”是一种在工程上应用广泛的“高精度单步算法”。
由于此算法精度高，采取措施对误差进行抑制，可以得出四阶龙格－库塔公式，也就是在工程中应用广泛的经典龙格－库塔算法，即：yi＋１＝yi＋h*( K1+ 2*K2 +2*K3+ K4)/6K1=f(xi,yi)K2=f(xi+h/2,yi+h*K1/2)K3=f(xi+h/2,yi+h*K2/2)K4=f(xi+h,yi+h*K3)通常所说的龙格-库塔法就是指四阶——龙格库塔法，我们可以仿二阶、三阶的情形推导出常用的标准四阶龙格-库塔法公式。
（1）计算公式（1）的局部截断误差是。
龙格-库塔法具有精度高，收敛，稳定（在一定条件下），计算过程中可以改变步长，不需要计算高阶导数等优点，但仍需计算在一些点上的值，如四阶龙格-库塔法每计算一步需要计算四次的值，这给实际计算带来一定的复杂性，因此，多用来计算“表头”。
二、小程序#include<stdio.h>#include<math.h>#define f(x,y) (-1*(x)*(y)*(y))void main(void){double a,b,x0,y0,k1,k2,k3,k4,h;int n,i;printf("input a,b,x0,y0,n:");scanf("%lf%lf%lf%lf%d",&a,&b,&x0,&y0,&n); printf("x0\ty0\tk1\tk2\tk3\tk4\n");for(h=(b-a)/n,i=0;i!=n;i++){k1=f(x0,y0);k2=f(x0+h/2,y0+k1*h/2);k3=f(x0+h/2,y0+k2*h/2);k4=f(x0+h,y0+h*k3);printf("%lf\t%lf\t",x0,y0);printf("%lf\t%lf\t",k1,k2);printf("%lf\t%lf\n",k3,k4);y0+=h*(k1+2*k2+2*k3+k4)/6;x0+=h;}printf("xn=%lf\tyn=%lf\n",x0,y0);}运行结果：input a,b,x0,y0,n:0 5 0 2 20x0y0k1k2k3 k40.000000 2.000000-0.000000-0.500000-0.469238-0.8861310.250000 1.882308-0.885771-1.176945-1.129082-1.2800600.500000 1.599896-1.279834-1.295851-1.292250-1.2227280.750000 1.279948-1.228700-1.110102-1.139515-0.9901621.000000 1.000027-1.000054-0.861368-0.895837-0.7528521.2500000.780556-0.761584-0.645858-0.673410-0.5621891.5000000.615459-0.568185-0.481668-0.500993-0.4205371.7500000.492374-0.424257-0.361915-0.374868-0.3178552.0000000.400054-0.320087-0.275466-0.284067-0.2435982.2500000.329940-0.244935-0.212786-0.218538-0.1894822.5000000.275895-0.190295-0.166841-0.170744-0.1495632.7500000.233602-0.150068-0.132704-0.135399-0.1197033.0000000.200020-0.120024-0.106973-0.108868-0.0970483.2500000.172989-0.097256-0.087300-0.088657-0.0796183.5000000.150956-0.079757-0.072054-0.073042-0.0660303.7500000.132790-0.066124-0.060087-0.060818-0.0553054.0000000.117655-0.055371-0.050580-0.051129-0.0467434.2500000.104924-0.046789-0.042945-0.043363-0.0398334.5000000.094123-0.039866-0.036750-0.037072-0.0342024.7500000.084885-0.034226-0.031675-0.031926-0.029571xn=5.000000yn=0.076927。
龙格库塔法的基本原理
龙格库塔法的基本原理嘿，朋友们！今天咱来唠唠龙格库塔法的基本原理。
咱就说啊，这龙格库塔法就像是一个超级厉害的解题高手！它能帮我们在面对那些复杂得让人头疼的数学问题时找到答案。
你想想，有时候那些数学式子就像一团乱麻，让你无从下手。
但龙格库塔法呢，它就有办法一点点地把这团乱麻给解开。
它通过巧妙地计算和推测，一步一步地接近问题的真相。
就好比你要去一个陌生的地方，你不知道该怎么走。
龙格库塔法就像是一个聪明的向导，它能给你指出一条清晰的路。
它不是那种随随便便指个方向就完事儿的，而是经过深思熟虑，考虑了各种因素之后才做出的决定。
龙格库塔法在很多领域都大显身手呢！比如在物理、工程这些地方，那可真是立下了汗马功劳。
没有它，好多难题都没法解决，好多伟大的工程可能都没法实现。
它的原理其实也不复杂，就是通过一系列的计算步骤来逼近真实的答案。
但可别小看这些步骤哦，它们就像一个个小小的拼图，组合起来就能呈现出一幅完整的画面。
比如说，它会在不同的点上进行计算，就像在不同的地方观察一样，然后把这些观察到的信息综合起来，得出一个比较准确的结果。
这是不是很神奇？而且啊，龙格库塔法还特别灵活。
它可以根据不同的问题进行调整和优化，就像一个武林高手能根据对手的情况随时改变招式一样。
你说，这么厉害的龙格库塔法，咱能不好好了解了解吗？咱得知道它是怎么工作的，怎么帮我们解决问题的呀！不然岂不是浪费了这么好的一个工具？总之，龙格库塔法就是数学世界里的一颗璀璨明星，照亮了我们解决难题的道路。
它让那些看似不可能的问题变得有可能，让我们在知识的海洋里能更加自由地航行。
所以啊，大家可别小瞧了它哟！。
龙格库塔方法的Miline-Hamming预测-校正算法实验报告知识讲解
龙格库塔方法的M i l i n e-H a m m i n g预测-校正算法实验报告2011-2012学年第2学期实验报告实验名称：微分方程数值解实验学院：******专业：**************班级：**********班内序号：**学号：********姓名：******任课教师：******北京邮电大学时间：****年**月**日实验目标用多环节Miline-Hamming 预测-校正算法求下列方程的解{y‘=y −2xy ,y (0)=1, 0≤x ≤4 其中解析解为 y (x )=√1+2x实验原理计算龙格库塔显示公式计算预测值，然后用隐式公式进行校正。
Miline-Hamming 预测-校正公式为{p n+1=u n−3+43h(2f n −f n−1+f n−2)m n+1=p n+1+112121(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-1K 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+nh4） 输出(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+hp=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否则停止。
龙格-库塔方法基本原理
c1 c2 c3 1
a2c2
a3c3
1 2
,
a22c2 a32c31 3,b221c2
(b31
b32 )2 c3
1 3
a 2 b32 c 2
1 6
,
b21c2
(b31
b32 )c3
1 2
a2b21c2
a3c3 (b31
b32 )
1 3
1 b21b32c3 6
2021/5/25
6
同理，改进Euler公式可改写成
y
i
1
yi
1 2
K1
1 2
K2
K
1
hf
( xi ,
yi )
K
2
hf
( xi
h, yi
K1)
局部截断误差为O(h3)
上述两组公式在形式上共同点:都是用f(x,y)在某些 点上值的线性组合得出y(xi+1)的近似值yi+1, 且增加计 算的次数f(x,y)的次数,可提高截断误差的阶。如欧拉 法:每步计算一次f(x,y)的值,为一阶方法。改进欧拉法 需计算两次f(x,y)的值，为二阶方法。
17
若取 a2b2 塔公式。
11 2,
c1，就0,是c2另一1种形式的二阶龙格
-
库
yi1 yi K2
K1 hf(xi, yi)
K2
hf(xi
1h, 2
yi
12K1)
此计算公式称为变形的二阶龙格—库塔法。式中
xi
1为h区间
2
x的i , x中i1点。也称中点公式。
Q:为获得更高的精度，应该如何进一步推广？
第四讲龙格-库塔方法
第四讲龙格-库塔⽅法龙格-库塔⽅法3.2 Runge-Kutta法3.2.1 显式Runge-Kutta法的⼀般形式上节已给出与初值问题(1.2.1)等价的积分形式(3.2.1)只要对右端积分⽤不同的数值求积公式近似就可得到不同的求解初值问题(1.2.1)的数值⽅法，若⽤显式单步法(3.2.2)当，即数值求积⽤左矩形公式，它就是Euler法(3.1.2)，⽅法只有⼀阶精度，若取(3.2.3)就是改进Euler法，这时数值求积公式是梯形公式的⼀种近似，计算时要⽤⼆个右端函数f的值，但⽅法是⼆阶精度的.若要得到更⾼阶的公式，则求积分时必须⽤更多的f值，根据数值积分公式，可将(3.2.1)右端积分表⽰为注意，右端f中还不能直接得到，需要像改进Euler法(3.1.11)⼀样，⽤前⾯已算得的f值表⽰为(3.2.3)，⼀般情况可将(3.2.2)的表⽰为(3.2.4)其中这⾥均为待定常数，公式(3.2.2)，(3.2.4)称为r级的显式Runge-Kutta法，简称R-K⽅法.它每步计算r个f值(即)，⽽k由前⾯(i-1)个已算出的表⽰，故公式是显式的.例i如当r=2时，公式可表⽰为(3.2.5) 其中.改进Euler 法(3.1.11)就是⼀个⼆级显式R-K ⽅法.参数取不同的值，可得到不同公式.3.2.2 ⼆、三级显式R-K ⽅法对r=2的显式R-K ⽅法(3.2.5)，要求选择参数，使公式的精度阶p 尽量⾼，由局部截断误差定义11122211()()[(,())(,)]n n n n n n n T y x y x h c f x y x c f x a h y b hk ++=--+++ (3.2.6)令，对(3.2.6)式在处按Taylor 公式展开，由于将上述结果代⼊(3.2.6)得要使公式(3.2.5)具有的阶p=2,即，必须(3.2.7)即由此三式求的解不唯⼀.因r=2，由（3.2.5）式可知，于是有解(3.2.8)它表明使(3.2.5)具有⼆阶的⽅法很多，只要都可得到⼆阶精度R-K⽅法.若取，则，则得改进Euler法(3.1.11)，若取，则得，此时(3.2.5)为(3.2.9)其中称为中点公式.改进的Euler法(3.1.11)及中点公式(3.2.9)是两个常⽤的⼆级R-K⽅法，注意⼆级R-K⽅法只能达到⼆阶，⽽不可能达到三阶.因为r=2只有4个参数，要达到p=3则在(3.2.6)的展开式中要增加3项，即增加三个⽅程,加上(3.2.7)的三个⽅程，共计六个⽅程求4个待定参数，验证得出是⽆解的.当然r=2,p=2的R-K⽅法(3.2.5)当取其他数时，也可得到其他公式，但系数较复杂，⼀般不再给出.对r=3的情形，要计算三个k值，即其中将按⼆元函数在处按Taylor公式展开，然后代⼊局部截断误差表达式，可得可得三阶⽅法，其系数共有8个，所应满⾜的⽅程为(3.2.10)这是8个未知数6个⽅程的⽅程组，解也是不唯⼀的，通常.⼀种常见的三级三阶R-K⽅法是下⾯的三级Kutta⽅法：(3.2.11)附：R-K 的三级Kutta ⽅法程序如下function y = DELGKT3_kuta(f, h,a,b,y0,varvec) format long; N = (b-a)/h;y = zeros(N+1,1); y(1) = y0; x = a:h:b;var = findsym(f); for i=2:N+1K1 = Funval(f,varvec,[x(i-1) y(i-1)]);K2 = Funval(f,varvec,[x(i-1)+h/2 y(i-1)+K1*h/2]); K3 = Funval(f,varvec,[x(i-1)+h y(i-1)-h*K1+K2*2*h]);y(i) = y(i-1)+h*(K1+4*K2+K3)/6; %满⾜c1+c2+c3=1,(1/6 4/6 1/6)endformat short; 3.2.3 四阶R-K ⽅法及步长的⾃动选择利⽤⼆元函数Taylor 展开式可以确定(3.2.4)中r=4,p=4的R-K ⽅法，其迭代公式为111223344()n n y y h c k c k c k c k +=++++其中1(,)n n k f x y =，2221(,(,))n n n n k f x a h y b hf x y =++，⽽33311322(,)n n k f x a h y b hk b hk =+++ 44411422433(,)n n k f x a h y b hk b hk b hk =++++共计13个参数待定，Taylor 展开分析局部截断误差，使得精度达到四阶，即误差为5()O h 。
龙格-库塔方法微积分教学
龙格-库塔方法8.2.1 二阶龙格-库塔方法常微分方程初值问题：做在点的泰勒展开：这里。
取，就有(8.11)截断可得到近似值的计算公式，即欧拉公式：若取，式（8.11）可写成：或（8.12）截断可得到近似值的计算公式：或上式为二阶方法，一般优于一阶的欧拉公式（8.2），但是在计算时，需要计算在点的值，因此，此法不可取。
龙格-库塔设想用在点和值的线性组合逼近式（8.12）的主体，即用（8.13）逼近得到数值公式：（8.14）或更一般地写成对式（8.13）在点泰勒展开得到：将上式与式（8.12）比较，知当满足时有最好的逼近效果，此时式（8.13）－式（8.14）。
这是4个未知数的3个方程，显然方程组有无数组解。
若取，则有二阶龙格-库塔公式，也称为改进欧拉公式：（8.15）若取，则得另一种形式的二阶龙格-库塔公式，也称中点公式：（8.16）从公式建立过程中可看到，二阶龙格-库塔公式的局部截断误差仍为，是二阶精度的计算公式。
类似地，可建立高阶的龙格-库塔公式，同时可知四阶龙格-库塔公式的局部截断误差为，是四阶精度的计算公式。
欧拉法是低精度的方法，适合于方程的解或其导数有间断的情况以及精度要求不高的情况，当解需要高精度时，必须用高阶的龙格-库塔等方法。
四阶龙格-库塔方法应用面较广，具有自动起步和便于改变步长的优点，但计算量比一般方法略大。
为了保证方法的收敛性，有时需要步长取得较小，因此，不适于解病态方程。
8.2.2 四阶龙格-库塔公式下面列出常用的三阶、四阶龙格-库塔计算公式。
三阶龙格-库塔公式（1）（8.17）（8.18）(3)（8.19）四阶龙格-库塔公式（1）（8.20）（8.21）例8.3用四阶龙格-库塔公式（8.20）解初值问题：解：取步长，计算公式为：计算结果列表8.3中。
表8.3 计算结果1 0.2 1.24789 1.24792 0.000032 3 4 0.40.60.81.637622.296183.533891.637782.296963.538020.000160.000780.004138.2.3 步长的自适应欧拉方法和龙格-库塔方法在计算时仅用到前一步的值，我们称这样的方法为单步法。
- 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次迭代的结果会对第五轮求得的数值产生影响。
因此，一旦前四轮轮迭代所得的结果有偏差，下一轮结果的偏差将大于之前的偏差。
因此会导致伴随迭代次数的增加而产生更大的偏差。