龙格库塔方法的Miline-Hamming预测-校正算法实验报告
四阶龙格库塔实验报告
三、四阶Runge-Kutta 法求解常微分方程一、龙格库塔法的思想根据第九章的知识可知道，Euler 方法的局部截断误差是2()O h ，而当用Euler 方法估计出1,()(1)n n n n y y hf x y +=+ 再用梯形公式111[(,)(,)](2)2n n n n n n h y y f x y f x y +++=++进行校正，即采用改进Euler 方法得出数值解的截断误差为3()O h 。
由Lagrange 微分中值定理'11()()()()()(,())(3)n n n n n y x y x y x x y x hf y ξξξ++=+-=+ 记*(,())k hf y ξξ=，得到*1()()(4)n n y x y x k +=+这样只要给出一种计算*k 的算法，就能得到相应的计算公式。
用这种观点的来分析Euler 方法和改进Euler 方法，Euler 方法的迭代公式可改写为111(,)n n n n y y k k hf x y +=+=改进Euler 方法的预报-校正公式可改写为 1121211()2(,),(,)n n n n n n y y k k k hf x y k hf x h y k +=++==++ Euler 方法实际上是用一个点处的值1k 近似*k ，而改进Euler 方法是用两个点处的值1k ，和2k ，做算术平均值近似*k 自然改进Euler 方法要优于Euler 方法。
因此，可以想到假如在1[,]n n x x +内多预报几个点值i k ，并用他们的加权平均值作为*k 的近似值，则有可能构造出具有更高精度的计算公式，这就是Runge-Kutta 法的基本思想。
二、四阶龙格库塔法由Runge-Kutta 的基本思想，构造四阶Runge-Kutta 法是利用1234,,k k k k 和的加权平均值来近似*k ，因此令1112233441211132211243312213(,)(,)(5)(,)(,)n n n n n n n n n n y y w K w K w K w K K hf x y K hf x h y K K hf x h y K K K hf x h y K K K αβαβγαβγη+=++++⎧⎪=⎪⎪=++⎨⎪=+++⎪⎪=++++⎩使得511()()n n y x y O h ++-=即其总体截断误差为4()O h 。
龙格库塔方法的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否则停止。
hamming实验报告
hamming实验报告Hamming实验报告引言：Hamming实验是一项重要的计算机科学实验，旨在研究和验证Hamming码的纠错能力。
Hamming码是一种用于纠正单一比特错误的错误检测和纠正编码方式，被广泛应用于数据传输和存储中。
本实验将通过模拟数据传输过程，并使用Hamming码进行纠错，来验证其在实际应用中的有效性。
实验目的：本实验的目的是通过模拟数据传输过程，验证Hamming码的纠错能力。
具体而言，我们将通过引入人为制造的错误，检测和纠正这些错误，以验证Hamming码的可靠性和有效性。
实验步骤：1. 设计Hamming码生成矩阵和校验矩阵。
2. 生成待发送的数据，并使用Hamming码进行编码。
3. 引入人为制造的错误，模拟数据传输过程中的错误。
4. 使用Hamming码进行错误检测和纠正。
5. 比较纠错前后的数据，验证Hamming码的纠错能力。
实验结果：在本次实验中，我们成功设计并实现了Hamming码的纠错过程。
通过引入人为制造的错误，我们模拟了数据传输过程中的错误情况。
使用Hamming码进行错误检测和纠正后，我们成功恢复了原始数据，并验证了Hamming码的纠错能力。
讨论：Hamming码作为一种常用的纠错编码方式，具有较高的纠错能力和可靠性。
通过本次实验，我们进一步验证了Hamming码的有效性。
然而，Hamming码并不能纠正所有错误，它只能纠正单一比特错误。
对于多比特错误或连续错误，Hamming码的纠错能力将受到限制。
因此，在实际应用中，我们需要综合考虑数据传输的可靠性需求，并选择适当的纠错编码方式。
结论：通过本次实验，我们验证了Hamming码的纠错能力。
Hamming码作为一种常用的纠错编码方式，在数据传输和存储中具有重要的应用价值。
然而，我们也需要认识到Hamming码的局限性，它只能纠正单一比特错误。
在实际应用中，我们需要根据具体需求选择适当的纠错编码方式，以确保数据的可靠性和完整性。
四阶龙格库塔实验报告
三、四阶Runge-Kutta 法求解常微分方程一、龙格库塔法的思想根据第九章的知识可知道，Euler 方法的局部截断误差是2()O h ，而当用Euler 方法估计出1,()(1)n n n n y y hf x y +=+ 再用梯形公式111[(,)(,)](2)2n n n n n n h y y f x y f x y +++=++进行校正，即采用改进Euler 方法得出数值解的截断误差为3()O h 。
由Lagrange 微分中值定理'11()()()()()(,())(3)n n n n n y x y x y x x y x hf y ξξξ++=+-=+ 记*(,())k hf y ξξ=，得到*1()()(4)n n y x y x k +=+这样只要给出一种计算*k 的算法，就能得到相应的计算公式。
用这种观点的来分析Euler 方法和改进Euler 方法，Euler 方法的迭代公式可改写为111(,)n n n n y y k k hf x y +=+=改进Euler 方法的预报-校正公式可改写为 1121211()2(,),(,)n n n n n n y y k k k hf x y k hf x h y k +=++==++ Euler 方法实际上是用一个点处的值1k 近似*k ，而改进Euler 方法是用两个点处的值1k ，和2k ，做算术平均值近似*k 自然改进Euler 方法要优于Euler 方法。
因此，可以想到假如在1[,]n n x x +内多预报几个点值i k ，并用他们的加权平均值作为*k 的近似值，则有可能构造出具有更高精度的计算公式，这就是Runge-Kutta 法的基本思想。
二、四阶龙格库塔法由Runge-Kutta 的基本思想，构造四阶Runge-Kutta 法是利用1234,,k k k k 和的加权平均值来近似*k ，因此令1112233441211132211243312213(,)(,)(5)(,)(,)n n n n n n n n n n y y w K w K w K w K K hf x y K hf x h y K K hf x h y K K K hf x h y K K K αβαβγαβγη+=++++⎧⎪=⎪⎪=++⎨⎪=+++⎪⎪=++++⎩使得511()()n n y x y O h ++-=即其总体截断误差为4()O h 。
龙格库塔实验报告
一、实验背景常微分方程（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. 将计算得到的近似解与解析解进行比较，分析四阶龙格-库塔方法的精度和稳定性。
第四讲龙格-库塔方法
龙格-库塔方法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 。
龙格-库塔方法-文档资料
c3a 2b32
c3a3 1
6
; 2
O (h4)
常见的2种三阶方法：
库塔三阶方法
h yn1yn6(k14k2k3)
k1
f(xn,yn);
k2
hh f(xn2,yn2k1)
k 3 f(x n h ,y n h k 1 2 h k 2 ) •5
四级方法：N = 4
局部截断误差 O ( h 5 )
可见误差随着 x n 的增加呈指数函数递减
当 f y 0 时，微分方程是不稳定的； 而 f y 0 时，微分方程是稳定的。
上面讨论的稳定性，与数值方法和方程中 f 有关
•21
实验方程： y y C ,R e () 0
D e f 3 对单步法 yn 1ynh(xn,yn,h )应用实验方程，
e n 1 e n h [ ( x n ,y ( x n ) , h ) ( x n ,y n , h ) ] T n 1
•15
因为单步法是 p 阶的：h0,0hh0满足|Tn1|Chp1
|e n 1| |e n| h L |e n| C h p 1|en |
其中 1hL,C hp1
•18
三、绝对稳定性 /*Absolute Stibility*/ 计算过程中产生的舍入误差对计算结果的影响
首先以Euler公式为例，来讨论一下舍入误差的传播:
yn1ynhf(xn,yn)
设实际计算得到的点 x n 的近似函数值为 yn yn n，
其中 y n 为精确值， n 为误差
yn1ynhf(xn,yn)
通过适当选取参数 1,2和 p 的值，使得公式具有 2阶精度！！
•3
由泰勒公式展开，要使公式具有 2 阶精度，只需
hamming实验报告
hamming实验报告Hamming实验报告引言：Hamming实验是一种用于错误检测和纠正的编码技术。
它由美国科学家理查德·哈明（Richard Hamming）于20世纪50年代提出，并被广泛应用于计算机科学和通信领域。
本文将介绍Hamming实验的背景、原理、实验过程和结果，并探讨其在现代技术中的应用。
一、背景：在信息传输和存储过程中，出现错误是不可避免的。
特别是在电信和计算机网络中，数据传输的准确性至关重要。
传统的奇偶校验方法只能检测错误，而无法纠正。
为了解决这个问题，Hamming提出了一种能够检测和纠正错误的编码技术。
二、原理：Hamming编码是一种基于二进制的线性块编码。
它通过在数据中插入冗余位来实现错误检测和纠正。
冗余位的数量取决于数据位的数量，并且通过一系列算法计算得出。
Hamming编码的关键思想是将数据位与冗余位进行异或运算，以便在接收端进行错误检测和纠正。
三、实验过程：1. 准备数据：选择一组二进制数据作为实验对象，例如1011。
2. 计算冗余位：根据数据位的数量，计算所需的冗余位数量。
在本例中，需要添加两个冗余位。
3. 插入冗余位：将冗余位插入数据中，形成编码后的数据。
例如，将冗余位插入1011变为1101011。
4. 传输数据：将编码后的数据传输到接收端。
5. 接收数据：接收端接收到编码后的数据。
6. 检测错误：接收端使用Hamming编码的算法检测错误位。
7. 纠正错误：如果检测到错误位，接收端可以使用Hamming编码的算法进行错误纠正。
四、实验结果：在本次实验中，我们选择了1011作为实验数据。
根据Hamming编码的原理，我们计算出了两个冗余位，并将其插入到数据中，形成了1101011的编码数据。
在传输过程中，我们模拟了一个错误位，并将其插入到编码数据中。
接收端成功检测到错误位，并进行了错误纠正。
最终，接收端得到了正确的数据1011。
五、应用：Hamming编码在现代技术中有着广泛的应用。
龙格-库塔方法
h 按照数学学 习一般方法， yn 1 yn （k1 2k 2 2k3 +k 4） 6 继续上面的 做法，经过 k1 f ( xn , yn ) 较复杂的数 h 学推理和计 k （16） 2 f ( xn 1 , yn k1 ) 算，可以导 2 2 出四阶龙格 h －库塔格式， k f ( x 1 , yn k 2 ) 这里有一个 3 n 2 2 较有用的四 k f ( x , y hk ) 阶格式。 n 1 n 3 4
？
3、欧拉格式
yn 1 yn hf ( xn , yn ),
取点 xn 的斜率 k
*
(3)
n 0,1, 2,
f ( xn , yn )作为区间 [ xn , xn1 ] 上的平均斜率，精度低。
4、改进的欧拉格式
于是：
h yn 1 yn [ f ( xn , yn ) f ( xn 1 , yn hf ( xn , yn ))] 2 h yn 1 yn [k1 k2 ] 2
龙格－库塔方法
考虑一维经典初值问题
dy f ( x , y ) , y( x0 ) y0 , x [a, b] dx
(1)
龙格－库塔方法
一、设计思想
加权平均斜率。
1、微分中值定理
如果 f (x) 在 [a, b]上连续，在 ( a, b) 内可导，则在( a, b) 内至 少存在一点 ，有 f (b) f (a) f ' ( ) ba
xnq xn qh, 0 p q 1
令
yn1 yn h[(1 )k1 k2 k3 ] (14)
k1 , k 2
龙格—库塔实验报告
龙格—库塔法解常微分方程实验报告一、实验题目求解初值问题⎪⎩⎪⎨⎧=<<-=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。
线性多步法
(6)计算 f3 f ( x3 , y3 ) x x3 h 4 112 p y0 (2 f3 f 2 2 f1 ) m p (c0 p0 ) 3 121 1 3 c (9 y3 y1 ) h[ f ( x, m) 2 f 3 f 2 ] 8 8 9 y c (c p ) 输出(x, y ) 121 (7)若n N，置n 1 n, x j 1 x j , y j 1 y j , f j 1 f j ( j 0,1, 2), x x3 , y y3 , p p0 , c c0 , 转6； 否则停机。
算法 ba （1）输入 a, b, f ( x, y ), N , y0 (2)置h , x0 a, n 1 N (3)计算 f n 1 f ( xn 1 , yn 1 ) K1 hf n 1 h K1 K 2 hf ( xn 1 , yn 1 ) 2 2
(i 1, 2, , N ) (i 1, 2, , N )
若对 y 和 f 采用向量的记号 y ' f ( x, y); y( x0 ) y 0 ;
f ( f1 , f 2 , f N ) ;
T
或表达为 h yi ,n 1 yin ( K i1 2 K i 2 2 K i 3 K i 4 ) 6 (i 1, 2, , N ), 其中 K i1 fi ( xn , y1n , y2n , , yNn ); h h h h K i 2 fi ( xn , y1n K11 , y2 n K 21 , , y Nn K N 1 ); 2 2 2 2 h h h h K i 3 fi ( xn , y1n K12 , y2 n K 22 , , yNn K N 2 ); 2 2 2 2 K i 4 fi ( xn h, y1n hK13 , y2 n hK 23 , , y Nn hK N 3 );
计算方法 实验报告 拉格朗日 龙贝格 龙格库塔
主界面：
/*lagrange.c*/
float real_value(float x) /*由被插值函数计算真实值*/
c=getchar();
if(c=='N'||c=='n') break;
}
}
/*romberg.c*/
double function(double x) /*被积函数*/
{
return 4.0/(1+(x)*(x));
}
double t(double a,double b,int m) /*计算T1*/
实验二（龙贝格公式）
§公式
§算法描述
§流程图
§运行结果
§结果分析：Romberg积分法是在积分步长逐步折半的过程中，用低精度求积公式的组合得到更高精度求积公式的一种方法，它算法简单，且收敛加速效果极其显著。
实验三（四阶龙格库塔）
§公式
k1=h*f(xn,yn);
k2=h*f(xn+h/2,yn+k1/2);
T1=t(a,b,0);
T2=T1/2.0+t(a,b,1);
S1=(4*T2-T1)/3.0;
T1=T2;
T2=T1/2.0+t(a,b,2);
S2=(4*T2-T1)/3.0;
C1=(16*S2-S1)/15.0;
T1=T2;
T2=T1/2.0+t(a,b,3);
S1=S2;
S2=(4*T2-T1)/3.0;
数学实验“欧拉数值算法,Runge-Kutta数值算法”实验报告(内含matlab程序)
西京学院数学软件实验任务书实验二十四实验报告一、实验名称：欧拉数值算法（显式，隐式，欧拉预估-校正法），Runge-Kutta 数值算法。
二、实验目的：进一步熟悉欧拉数值算法（显式，隐式，欧拉预估-校正法），Runge-Kutta 数值算法。
三、实验要求：运用Matlab/C/C++/Java/Maple/Mathematica 等其中一种语言完成程序设计。
四、实验原理：1．欧拉数值算法（显式）：微分方程里最简单的方程形式莫过于一阶常微分方程的初值问题，即：0(,)()dy f x y a x b dx y a y ⎧=≤≤⎪⎨⎪=⎩其中a ，b 为常数。
因为其简单但又是求解其他方程的基础，所以发展了许多典型的解法。
所有算法中的f 就是代表上式中(,)f x y ，而y f 表示(,)f x y y ∂∂，x f 表示(,)f x y x∂∂。
简单欧拉法是一种单步递推算法。
简单欧拉法的公式如下所示：1(,)n n n n y y hf x y +=+简单欧拉法的算法过程介绍如下:给出自变量x 的定义域[,]a b ，初始值0y 及步长h 。
对0,1,()/k b a h =-，计算1(,)k k k k y y hf x y +=+2．欧拉数值算法（隐式）：隐式欧拉法也叫退欧拉法，隐式欧拉法的公式如下所示：111(,)n n n n y y hf x y +++=+隐式欧拉法是一阶精度的方法，比它精度高的公式是：111[(,)(,)]2n n n n n n h y y f x y f x y +++=++ 隐式欧拉的算法过程介绍如下。
给出自变量x 的定义域[,]a b ，初始值0y 及步长h 。
对0,1,()/k b a h =-，用牛顿法或其他方法求解方111(,)k k k k y y hf x y +++=+得出1k y +。
3．欧拉预估-校正法：改进欧拉法是一种二阶显式求解法，其计算公式如下所示：1[,(,)]22n n n n n n h h y y hf x y f x y +=+++11(,)[(,)(,)]2n n n n n n n n t y hf x y h y y f x y f x t ++=+⎧⎪⎨=++⎪⎩ 4．Runge-Kutta 数值算法：二阶龙格-库塔法有多种形式，除了改进的欧拉法外，还有中点法。
龙格库塔方法的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否则停止。
第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)方法”是一种在工程上应用广泛的“高精度单步算法”。
由于此算法精度高，采取措施对误差进行抑制，可以得出四阶龙格－库塔公式，也就是在工程中应用广泛的经典龙格－库塔算法，即：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。
龙格-库塔方法
8.2 龙格-库塔方法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）（2）（8.18）(3)（8.19）四阶龙格-库塔公式（1）（8.20）（2）（8.21）例8.3用四阶龙格-库塔公式（8.20）解初值问题：解：取步长，计算公式为：计算结果列表8.3中。
表8.3 计算结果8.2.3 步长的自适应欧拉方法和龙格-库塔方法在计算时仅用到前一步的值，我们称这样的方法为单步法。
龙格库塔实验报告
龙格库塔实验报告龙格库塔实验报告引言：龙格库塔法（Runge-Kutta method）是一种常用的数值求解常微分方程（ODE）的方法。
它是由德国数学家卡尔·龙格（Carl Runge）和马丁·威尔海姆·库塔（Martin Wilhelm Kutta）于19世纪末独立发展而来。
龙格库塔法通过将微分方程转化为一系列逼近值，从而实现对微分方程的数值求解。
本实验旨在通过对龙格库塔法的研究和实践，深入了解其原理和应用。
一、龙格库塔法的原理龙格库塔法的基本思想是将微分方程的解分解为一系列逼近值，然后通过逐步迭代的方式计算这些逼近值。
具体而言，龙格库塔法将求解微分方程的过程分为多个步骤，每个步骤都计算一个逼近值，并利用这些逼近值逐步逼近真实解。
二、龙格库塔法的步骤1. 确定微分方程和初始条件：首先，需要明确待求解的微分方程以及初始条件。
微分方程可以是一阶或高阶的，初始条件是方程在某一点上的已知值。
2. 确定步长：步长（或称为时间间隔）决定了逼近值的精度。
较小的步长能够提高逼近值的准确性，但也会增加计算量。
因此，需要根据具体问题的需求选择合适的步长。
3. 迭代计算：根据龙格库塔法的公式，进行逐步的迭代计算。
首先，根据初始条件计算出第一个逼近值。
然后，利用该逼近值和微分方程的导数计算出下一个逼近值。
重复这一过程，直到达到所需的迭代次数或满足精度要求。
4. 计算误差：为了评估逼近值的准确性，需要计算误差。
误差可以通过与解析解的比较得出，或者通过比较不同步长下的逼近值得出。
较小的误差表明逼近值较为准确。
三、龙格库塔法的应用龙格库塔法广泛应用于科学和工程领域，特别是在求解微分方程模型的数值计算中。
以下是一些常见的应用场景：1. 物理学：龙格库塔法可以用于模拟物体在重力场中的运动，如自由落体、抛体运动等。
通过求解微分方程，可以得到物体的位置、速度等随时间的变化规律。
2. 经济学：龙格库塔法可以用于经济学模型的求解，如经济增长模型、投资模型等。
Adams与Milline-Hamming预测校正公式例题
例：y(0)=1 (0x1)本题采用adams与Miline-Haming预测校正公式求解的一种解法。
f = sym('y-2*x/y'); %sym为定义符号变量函数。
例f=sym(‘sin2’) 输出结果为f=sin2x = 0; y = 1; %初始条件h = 0.1; %步长a = 0;b = 1; %上下线n = (b-a)/h+1; %步数sjs=[x,y]; %先定义一个数值解；sjs符号自定义即可，没有要求adx=0; %adams方法的x初始值；%下面的是龙格库塔方法求解，这里求得的是[0.1]中的所有解，为和后面的方法进行比较for i = 1:n-1RK(i)=y;L1=subs(f,{'x','y'},{x,y});L2=subs(f,{'x','y'},{x+h/2,y+h*L1/2}); %subs是替代函数，也就是用后面的代替前面的%例如：L2中，用(x+h/2,y+h*L1/2)代替原(x,y);此时函数应为f(x+h/2,y+h*L1/2);%下面同理L3=subs(f,{'x','y'},{x+h/2,y+h*L2/2});L4=subs(f,{'x','y'},{x+h,y+h*L3});y=y+h*(L1+2*L2+2*L3+L4)/6;x=x+h;adx(i+1)=x;sjs=[sjs;x,y];endady(1)=RK(1); mh(1)=RK(1); %ady为adams方法的y值，而mh为另一种方法的y值ady(2)=RK(2); mh(2)=RK(2); %因为这两种方法均要用龙格库塔方法求出前几项ady(3)=RK(3); mh(3)=RK(3);ady(4)=RK(4); mh(4)=RK(4);%下面是adams预测-校正方法for j=5:npre1=ady(j-1)+h*(55*subs(f,{'x','y'},{adx(j-1),ady(j-1)})-59*subs(f,{'x','y'},{adx(j-2),ady(j-2)})+37*subs(f,{'x','y'},{adx(j-3),ady(j-3)})-9*subs(f,{'x','y'},{adx(j-4),ady(j-4)}))/24;ady(j)=ady(j-1)+h*(9*subs(f,{'x','y'},{adx(j),pre1})+19*subs(f,{'x','y'},{adx(j-1),ady(j-1)})-5*subs(f,{'x','y'},{adx(j-2),ady(j-2)})+subs(f, {'x','y'},{adx(j-3),ady(j-3)}))/24;end%下面是Miline-Hamming预测校正方法for k=5:11pre2=mh(k-4)+4*h*(2*subs(f,{'x','y'},{adx(k-1),mh(k-1)})-subs(f,{'x','y'},{adx(k-2),mh(k-2)})+2*subs(f,{'x','y'},{adx(k-3),mh(k-3)}))/3;mh(k)=(9*mh(k-1)-mh(k-3))/8+3*h/8*(subs(f,{'x','y'},{adx(k),pre2})+2*subs(f,{'x','y'},{adx(k-1),mh(k-1)})-subs(f,{'x','y'},{adx(k-2),mh(k-2)})); end%输出结果并打印图像sjs %RK公式结果plot(sjs(:,1),sjs(:,2),'r') %RK公式图像hold onady %adams公式结果plot(sjs(:,1),ady,'ob') %adams公式图像hold onexact=sqrt(1+2*sjs(:,1)) %精确解plot(sjs(:,1),exact,'*g'); %精确解图像mh %mh公式结果hold onplot(sjs(:,1),mh,'k+'); %mh公式图像。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
2011-2012学年第2学期实验报告
实验名称：微分方程数值解实验学院：******
专业：**************班级：**********
班内序号：**
学号：********
姓名：******
任课教师：******
北京邮电大学
时间：****年**月**日
实验目标
用多环节Miline-Hamming预测-校正算法求下列方程的解
y‘=y−2x
,
y0=1,
0≤x≤4
其中解析解为 y x=1+2x
实验原理
计算龙格库塔显示公式计算预测值，然后用隐式公式进行校正。
Miline-Hamming预测-校正公式为
p n+1=u n−3+4
h(2f n−f n−1+f n−2)
m n+1=p n+1+112
121
(c n−p n)
c n+1=1
9u n−u n−2+
3
h[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，u0
2）置h=(b-a)/N，t0=a，n=1
3）计算f n-1=f(t n-1,u n-1)
K1=hf n-1
K2=hf(t n-1+h/2, u n-1+K1/2)
K3=hf(t n-1+h/2, u n-1+K2/2)
K4=hf(t n-1+h, u n-1+K3)
u n= u n-1+1/6(K1 +2K2 +2K3 +K4)
t n=a+nh
4）输出(t n ,u n)
5）若n<3，置n+1→n，返回3；
否则，置n+1→n，0→p0，0→c0，转6. 6）计算
f3=f(t3,u3)
t= t3+h
p=u0+4/3(2 f3–f2 +2f1)
m=p+112/121(c0-p0)
c=1/8(9u3- u1)+3/8h[f(t,m)+ 2 f3–f2]
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所示：
t坐标 0 0.2000 0.4000 0.6000 0.8000 1.0000 1.2000 1.4000 1.6000
真值 1.0000 1.1832 1.3416 1.4832 1.6125 1.7321 1.8439 1.9494 2.0494
M-H法 2.1445 2.2357 2.3233 2.4077 2.4890 2.5677 2.6438 2.7174 2.7887
t坐标
真值
表1 Miline-Hamming预测-校正法与解析解在各点数值比较（h=0.1）
为了评判Miline-Hamming预测-校正法的算法精度，在这里我们利用相对误差
的概念进行评判。
对于Miline-Hamming预测-校正法的每个的估计值有：
相对误差=估计值－真值
真值
从而我们可以通过计算得到如下的相对误差表：
t坐标
t坐标
很明显，当我们对各点处的相对误差取平均后，该平均值小于0.01。
因此，我们可以认为Miline-Hamming预测-校正法的在h=0.2时的算法精度相对较高，所得到的结果与真值较为接近。
接下来我们在对h=0.5和h=0.1的情况进行计算，可以得到结果如下表3和表4所示
t坐标0 0.5000 1.0000 1.5000 2.0000 2.5000 3.0000 3.5000 4.0000 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
1.0000 1.0954 1.1832 1.2649 1.3416 1.4142 1.4832 1.5492 1.6125
真值
真值 2.8636 2.8983 2.9326 2.9665 3.0000
表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次迭代的结果会对第五轮求得的数值产生影响。
因此，一旦前四轮轮迭代所得的结果有偏差，下一轮结果的偏差将大于之前的偏差。
因此会导致伴随迭代次数的增加而产生更大的偏差。