Romberg龙贝格算法实验报告.

合集下载

龙贝格求 积分

龙贝格求 积分

龙贝格(Romberg )求积法1.算法理论Romberg 求积方法是以复化梯形公式为基础,应用Richardson 外推法导出的数值求积方法。

由复化梯形公式 )]()(2)([2222b f h a f a f h T +++=可以化为)]()]()([2[212112h a f h b f a f hT +++==)]([21211h a f h T ++一般地,把区间[a,b ]逐次分半k -1次,(k =1,2,……,n)区间长度(步长)为kk m a b h -=,其中mk =2k -1。

记k T =)1(k T由)1(k T =]))12(([21211)1(1∑=---++km j k k k h j a f h T 从而⎰badxx f )(=)1(kT-)(''122k f h a b ξ- (1)按Richardson 外推思想,可将(1)看成关于k h ,误差为)(2k h O 的一个近似公式,因而,复化梯形公式的误差公式为⎰badxx f )(-)1(k T =......4221++kkh K h K =∑∞=12i i k i h K (2)取1+k h =k h 21有 ⎰ba dx x f )(-)1(1+k T =∑∞=+121221i ik ii hK (3)误差为)(2jh O 的误差公式 )(j kT=)1(-j kT+141)1(1)1(------j j k j k T T2。

误差及收敛性分析(1)误差,对复化梯形公式误差估计时,是估计出每个子区间上的误差,然后将n 个子区间上的误差相加作为整个积分区间上的误差。

(2)收敛性,记h x i =∆,由于∑=++=ni i i n x f x f h f T 01))]()([2)(=))()((21101∑∑-==∆+∆n i ni i i i i x x f x x f上面两个累加式都是积分和,由于)(x f 在区间],[b a 上可积可知,只要],[b a 的分划的最大子区间的长度0→λ时,也即∞→n 时,它们的极限都等于积分值)(f I 。

龙贝格实验报告

龙贝格实验报告

龙贝格实验报告龙贝格实验报告龙贝格实验是一项经典的心理学实验,旨在探究人类对于延迟满足的选择行为。

实验的设计非常简单,但却能揭示出人类的心理特点和决策模式。

本文将对龙贝格实验进行介绍,并探讨实验结果对于我们日常生活中的决策行为的启示。

龙贝格实验最早由德国心理学家沃尔夫冈·龙贝格于1971年提出,他将实验对象置于一个选择情境中。

在这个情境中,实验对象需要在两个选项之间进行选择,一方面是能够立即获得一定的奖励,另一方面是需要等待一段时间才能获得更大的奖励。

实验对象可以根据自己的意愿选择不同的选项。

实验结果显示,大多数人更倾向于选择立即获得的小额奖励,而不是等待更长时间以获得更大的奖励。

这种行为模式被称为“即时满足偏好”。

这一结果揭示了人类的一种心理特点:我们更倾向于追求即时的满足感,而忽视了长期利益。

那么,为什么人们会有这种即时满足偏好呢?心理学家提出了一些解释。

首先,人们往往更容易预测和评估即时奖励,而对于未来的奖励往往难以准确评估。

这导致了人们对于即时奖励的价值更加敏感,而对于未来奖励的价值更加模糊。

其次,人们对于即时奖励的需求更加迫切,因为它们能够立即满足我们的欲望和需求,而未来奖励则需要等待和付出更多的努力。

然而,即时满足偏好并不总是明智的选择。

在现实生活中,我们常常需要做出一系列决策,这些决策涉及到长期利益和短期满足之间的权衡。

例如,我们是否应该立即购买心仪已久的商品,还是将钱储蓄起来以应对未来的不确定性?我们是否应该立即享受美食,还是节制自己以保持健康的身体?这些决策都需要我们权衡即时满足和长期利益之间的关系。

为了克服即时满足偏好的影响,我们可以采取一些策略。

首先,我们可以通过增强对未来奖励的预期和评估来提高对其价值的认识。

这意味着我们需要更加清晰地设定未来奖励的目标,并了解其对我们个人发展和幸福的重要性。

其次,我们可以通过设定明确的规划和目标来引导自己的行为。

例如,我们可以设定一个长期的储蓄计划,以确保我们能够在未来获得更大的回报。

龙贝格算法-数值分析-实验报告

龙贝格算法-数值分析-实验报告

3 / 13
一、任务综述于分析
1.1 任务
1、 根据所掌握龙贝格算法,独立编写程序,实现具体问题的求解等; 2、 自己设计一个多项式,根据设定的代数多项式进行测试,调式程序; 3、 对所设计的程序解决实际问题并分析,撰写分析报告
1.2 分析
1、 首先先明确什么是龙贝格算法,了解他的基本原理,画出流程图,编写程序代码。 2、 自己设计一个多项式,对他利用龙贝格算法求出他的积分。并与真实值比较,检测程序是否有错误。 并进一步修改完善程序。可以多试几个多项式,来判断程序的正确性。 3、 任务三实际上是求平面曲线长的问题。利用弧长微分公式可求得结果。
四、测试 ................................................................................................................................................. 8
主要参考资料
《数值分析简明教程》 《数值分析课程设计》
王能超 陈越、童若锋
高等教育出版社 浙江大学出版社
教研室签字: 年 月 日
院签字: 年 月 日
2 / 13
目录
目录 ......................................................................................................................................................... 3 一、任务综述于分析 ............................................................................................................................. 4

Romberg算法的实验报告

Romberg算法的实验报告

Romberg 算法一、实验目的:学会数值求积的Romberg 算法,并应用该算法于实际问题.二、实验内容:求定积分 ⎰15.0dxx三、实验要求:(1)要求程序不断加密对积分区间的等分,自动地控制Romberg 算法中的加速收敛过程,直到定积分近似值的误差不超过610-为止,输出求得的定积分近似值。

(2)可用MATLAB 中的内部函数int 求得此定积分的准确值与Romberg 算法计算的近似值进行比较。

四、实验基本原理Romberg 方法是使用行很强的一种数值积分方法,其收敛速度很快,这里直接给出Romberg 积分的计算方法。

(1)计算)]()()[(21)0,0(b f a f a b R +-=(2)计算∑-=--⎪⎪⎭⎫ ⎝⎛⎪⎭⎫ ⎝⎛-++-=12111212)0,1(21)0,(i k i i h k a f h i R r R (3)14)1,1()1,(4),(11-----=--j j j m R j m R j m R这样就构成了Romberg 积分的基本步骤,其计算步骤可以表1.1来表示:表1.1 Romberg 积分R(1,1)R(2,1) R(2,2)R(3,1) R(3,2) R(3,3)R(4,1) R(4,2) R(4,3) R(4,4)R(5,1) R(5,2) R(5,3) R(5,4) R(5,5)可以证明Romberg 方法是数值稳定的。

五、实验过程:1、编写主函数。

打开Editor编辑器,输入romberg法主程序语句:function [R,wugu,h]=romberg(fun,a,b, wucha,m)n=1;h=b-a; wugu=1; x=a;k=0; RT=zeros(4,4);RT(1,1)=h*(feval(fun,a)+feval(fun,b))/2;while((wugu>wucha)&(k<m)|(k<4))k=k+1; h=h/2; s=0;for j=1:nx=a+h*(2*j-1); s=s+feval(fun,x);endRT(k+1,1)= RT(k,1)/2+h*s; n=2*n;for i=1:kRT(k+1,i+1)=((4^i)*RT(k+1,i)-RT(k,i))/(4^i-1);endwugu=abs(RT(k+1,k)-RT(k+1,k+1));endR=RT(k+1,k+1);以文件名romberg.m保存。

龙贝格算法

龙贝格算法
2
2)理查森外推加速 从梯形公式出发, 将区间[a,b]逐次二分可以提高求积公式精度, 当[a,b] 分为 n 等份时,若记T������ = T ℎ ,当区间[a,b]分为 2n 等份时,则有 T2������ = T
ℎ 2
。再有泰勒公式展开为: T ℎ = ������ + ������1 ℎ2 + ������2 ℎ+ ⋯ +
建立一个命名为Romberg.m的function文件:
function[T]=Romberg(f,a,b,e) T=zeros(10,10); T(1,1)=(b-a)/2*(f(a)+f(b)); for k=2:10, sum=0; for i=1:2^(k-2), x=a+(2*i-1)*(b-a)/2^(k-1); sum=sum+f(x); end T(k,1)=T(k-1,1)/2+[(b-a)/2^(k-1)]*sum; for j=2:k, end %第一列用递推梯形公式 %定义龙贝格函数 %定义10阶的零元矩阵
1 3/2 ������ 0
������������。
算,并取������ ������, ������ ≈ ������ ;否则令 k=k+1,转(3)继续计算。 6)下图为我按照自己的算法所设计的示意表: 算法设计表:
k 1 2 3 4 …
h b-a (b-a)/2 (b-a)/4 (b-a)/8 …
������ ������, 1 T(1,1) T(2,1) T(3,1) T(4,1)
������ ������, 2
������ ������, 3
������ ������, 4
������ ������, 5

龙贝格(Romberg)求积法

龙贝格(Romberg)求积法
数值计算方法
龙贝格(Romberg)求积法
复化求积方法对于提高计算精度是行之有效的 方法,但复化公式的一个主要缺点在于要先估计出 步长。若步长太大,则难以保证计算精度,若步长 太小,则计算量太大,并且积累误差也会增大。在 实际计算中通常采用变步长的方法,即把步长逐次 分半,直至达到某种精度为止。
1.1变步长的梯形公式 变步长复化求积法的基本思想是在求积过程中,
)
0.9397933
进一步二分求积区间,并计算新分点上的函数值
f
(
1 4
)
0.9896158
,
f
(
3 4
)
0.9088516

T4
1 2 T2
1 4
f
(
1 4
)
f
(
3 4
)
0.9445135
这样不断二分下去,计算结果如P139列表所 示。积分的准确值为0.9460831,从表中可
看出用变步长二分10次可得此结果。
的考积察分T值与nS等n 。份辛卜生公式 S n之间的关系。将
复化梯形公式
Tn
h 2
f
(a)
2
n1 k 1
f
(xk ) f (b)
梯形变步长公式
T2 n
Tn 2
h n1
2 k 0
f (xk1 ) 2
代入(6.9) T 表达式得
h
n1
n1
T
6 f (a) 4k0
f
(
x k
1
)
2

2
输入 a,b,ε


b-ah,
h 2
f(a)+f(b)T1

数值积分龙贝格matlab

数值积分龙贝格matlab

《数值分析》课程实验报告一、实验目的1、进一步熟悉向量矩阵的运算;2、掌握龙贝格(Romberg )算法,并能用高级程序语言MATLAB 编写实现此算法的程序;3、进而加深对龙贝格(Romberg )算法的理解。

二、实验内容1. 使用Romberg 积分,对于计算下列⎰+4802)cos (1dx x 各近似值a.确定1,51,41,31,21,1,,,,R R R R Rb.确定5,54,43,32,2,,,R R R Rc.6,65,64,63,62,61,6,,,,,R R R R R Rd.确定10,109,98,87,7,,,R R R R三、实验步骤1. 编写程序龙贝格积分方法如下:n=5;a=0;b=48;h(1,1)=b-a;fa=sqrt(1+(cos(a))^2);fb=sqrt(1+(cos(b))^2);r(1,1)=h(1,1)/2*(fa+fb);disp('R11,R21,R31,R41,R51分别为');disp(r(1,1));for i=2:nh(i,1)=(b-a)/(2^(i-1));sum=0;for k=1:2^(i-2)x=a+(2*k-1)*h(i,1);sum=sum+sqrt(1+(cos(x)).^2);endr(i,1)=0.5*(r(i-1,1)+h(i-1,1)*sum);disp(r(i,1));enddisp('R22,R33,R44,R55分别为');for k=2:nfor j=2:kr(k,j)=r(k,j-1)+(r(k,j-1)-r(k-1,j-1))/(4^(j-1)-1);enddisp(r(k,k));enddisp('R61,R62,R63,R64,R65,R66分别为');n=6;for i=2:nh(i,1)=(b-a)/(2^(i-1));sum=0;for k=1:2^(i-2)x=a+(2*k-1)*h(i,1);sum=sum+sqrt(1+(cos(x)).^2);endr(i,1)=0.5*(r(i-1,1)+h(i-1,1)*sum);endfor k=2:nfor j=2:kr(k,j)=r(k,j-1)+(r(k,j-1)-r(k-1,j-1))/(4^(j-1)-1);endendfor i=1:ndisp(r(6,i));enddisp('R77,R88,R99,R10,10分别为');n=10;for i=2:nh(i,1)=(b-a)/(2^(i-1));sum=0;for k=1:2^(i-2)x=a+(2*k-1)*h(i,1);sum=sum+sqrt(1+(cos(x)).^2);endr(i,1)=0.5*(r(i-1,1)+h(i-1,1)*sum);endfor k=2:nfor j=2:kr(k,j)=r(k,j-1)+(r(k,j-1)-r(k-1,j-1))/(4^(j-1)-1);endendfor i=7:10disp(r(i,i));end运行结果如下:R11,R21,R31,R41,R51分别为62.437457.288656.443856.263156.2188R22,R33,R44,R55分别为55.572356.201556.205656.2041R61,R62,R63,R64,R65,R66分别为58.362759.077359.268959.317559.329759.3328R77,R88,R99,R10,10分别为58.422158.470758.470558.4705四、实验小结在这次编程中我学到了很多东西,把程序写进软件中也出现了很多错误,细节问题使我们必须注意的,自己有了很多的收获,自己进一步理解和学习了Matlab软件。

龙贝格算法实验报告

龙贝格算法实验报告

实验二:龙贝格算法一、实验目的1、通过本实验理解数值积分与微分的基本原理2、掌握数值积分中常见的复合求积公式的编程实现3、掌握龙贝格算法的基本思路和迭代步骤二、实验原理三、运行结果三、代码using System;using System.Collections.Generic;using System.Linq;using System.Text;namespace ConsoleApplication4{public delegate double F(double x);class Program{const double Precision = 0.00000000001;const int MAXRepeat = 10;static double f1(double x){double s=4/(1+x*x );return s;}static double Romberg(double a,double b, F f){int m,n,k;double[] y = new double[MAXRepeat];double h,ep,p,xk,s,q=0;h=b-a;y[0]=h*(f(a)+f(b))/2.0;//计算T`1`(h)=1/2(b-a)(f(a)+f(b));m=1;n=1;ep=Precision+1;while((ep>=Precision)&&(m<MAXRepeat)){p=0.0;for(k=0;k<n;k++){xk = a + (k + 0.5) * h; // n-1p = p + f(xk); //计算∑f(xk+h/2),T} // k=0p = (y[0] + h * p) / 2.0; //T`m`(h/2),变步长梯形求积公式s = 1.0;for (k = 1; k <= m; k++){s = 4.0 * s;// pow(4,m)q = (s * p - y[k - 1]) / (s - 1.0);//[pow(4,m)T`m`(h/2)-T`m`(h)]/[pow(4,m)-1],2m阶牛顿柯斯特公式,即龙贝格公式y[k - 1] = p;p = q;}ep = Math.Abs(q - y[m - 1]);//前后两步计算结果比较求精度m = m + 1;y[m - 1] = q;n = n + n; // 2 4 8 16h = h / 2.0;//二倍分割区间}return q;}static void Main(string[] args){double a, b, result;Console.WriteLine("请输入积分下限:");a = Convert.ToDouble(Console.ReadLine());Console.WriteLine("请输入积分上限:");b = Convert.ToDouble(Console.ReadLine());result = Romberg(a, b, new F(f1));Console.Write("定积分计算结果为:{0}:", result);Console.ReadLine();}}}四、分析本次试验使我认识到了计算机计算能力的强大,通过本次实验对数值积分与微分的基本原理有了深刻理解。

龙贝格积分算法实验

龙贝格积分算法实验

实验题目2 Romberg 积分法摘要考虑积分()()b aI f f x dx =⎰欲求其近似值,可以采用如下公式: (复化)梯形公式 11[()()]2n ii i hT f x f x-+==+∑2()12b a E h f η-''=-[,]a b η∈ (复化)辛卜生公式 11102[()4()()]6n i i i i hS f x f x f x -++==++∑4(4)()1802b a h E f η-⎛⎫=- ⎪⎝⎭ [,]a b η∈(复化)柯特斯公式 111042[7()32()12()90n i i i i hC f x f x f x -++==+++∑31432()7()]i i f xf x +++6(6)2()()9454b a h E f η-⎛⎫=- ⎪⎝⎭[,]a b η∈这里,梯形公式显得算法简单,具有如下递推关系121021()22n n n i i h T T f x -+==+∑因此,很容易实现从低阶的计算结果推算出高阶的近似值,而只需要花费较少的附加函数计算。

但是,由于梯形公式收敛阶较低,收敛速度缓慢。

所以,如何提高收敛速度,自然是人们极为关心的课题。

为此,记0,k T 为将区间[,]a b 进行2k等份的复化梯形积分结果,1,k T 为将区间[,]a b 进行2k等份的复化辛卜生积分结果,2,k T 为将区间[,]a b 进行2k等份的复化柯特斯积分结果。

根据李查逊(Richardson )外推加速方法,可得到1,11,,0,1,2,40,1,2,41m m k m km k m k T T T m -+-=-⎛⎫=⎪=-⎝⎭可以证明,如果()f x 充分光滑,则有,lim ()m k k T I f →∞= (m 固定),0lim ()m m T I f →∞=这是一个收敛速度更快的一个数值求积公式,我们称为龙贝格积分法。

该方法的计算可按下表进行0,0T 0,1T 0,2T … 0,m T 1,0T 1,1T … 1,1m T - 2,0T … 2,2m T - … … ,0m T很明显,龙贝格计算过程在计算机上实现时,只需开辟一个一维数组,即每次计算的结果,m k T ,可存放在0,k T 位置上,其最终结果,0m T 是存放在0,0T 位置上。

Romberg数值积分实验报告

Romberg数值积分实验报告

Romberg数值积分实验报告班级:__ 姓名:学号:日期:_一、实验目的1、加深外推法的原理理解, 掌握Romberg外推法的计算方法。

2、用matlab软件实现Romberg数值积分来计算题目的运算。

二、基本理论及背景1、理论推导:对区间[a,b],令h=b-a构造梯形值序列{T2K}。

T1=h[f(a)+f(b)]/2把区间二等分,每个小区间长度为h/2=(b-a)/2,于是T2 =T1/2+[h/2]f(a+h/2)把区间四(2)等分,每个小区间长度为h/2 =(b-a)/4,于是T4 =T2/2+[h/2][f(a+h/4)+f(a+3h/4).....................把[a,b] 2等分,分点xi=a+(b-a)/ 2 ·i (i =0,1,2 · · · 2k)每个小区间长度为(b-a)/ 2 .2、参考Romberg数值积分,实现积分的数值求解,完成下列题目:三、算法设计及实现1、算法设计(a)function y=fun1(x)y=sqrt(4-(sin(x)).^2);(b)function y=fun2 (x)y=log(1+2*x)/3;(c)function y=fun3(x)y=exp(-x).*sin(x.^2)四、实验步骤1、○1打开matlab软件,新建Romberg.m文件,在窗口中编辑Romber数值积分函数程序代码,并保存在指定的文件夹下,在Current Directory窗口右边点击《Browse For Folder》按钮指向Romberg.m文件;○2在Command Window中编辑相应要计算的题目的数值函数及相应的题目的表达式。

2、输出结果和初步分析说明(见附件一)。

五、使用说明实验结果分析1、在Command Window窗口中编辑要调用的函数名与指定的函数名字不同导致出现错误,通过改正与函数名相同即可。

龙贝格(Romberg)算法的应用实验报告

龙贝格(Romberg)算法的应用实验报告

Lab4 龙贝格(Romberg)算法的应用下面图1中的塑料雨蓬材料是由图2中所示的长方形平板塑料材料压制而成。

图1 图2已知图1的横截面曲线形状满足函数,则给定了雨蓬的长度后,要求需要平板原材料的长度。

函数接口定义:double Integral(double a, double b, double (*f)(double x, double y, double z), double TOL, double l, double t)在接口定义中:a、b分别为定积分的上、下界,f是积分核函数,其中x是积分哑元,y、z是本题目定义的特殊参数,分别对应中的l和t;TOL是要求积分达到的精度;l和t传入裁判输入的参数l和t的值。

另注意:的单位是厘米,输出的长度值要求以米为单位。

裁判程序样例如下:#include<stdio.h>#include<math.h>double f0( double x, double l, double t ){ /* 弧长积分的核函数*/return sqrt(1.0+l*l*t*t*cos(t*x)*cos(t*x));}double Integral(double a, double b, double (*f)(double x, double y, double z), double TOL, double l, double t);int main(){double a=0.0, b, TOL=0.005, l, t;while (scanf("%lf %lf %lf", &l, &b, &t) != EOF)printf("%.2f\n", Integral(a, b, f0, TOL, l, t));return 0;}裁判输入样例:2 100 1标准输出样例:1.68实验报告:1.求解步骤参照书上的龙贝格求积算法2.步骤①利用k=0;h=b-a;T[0][0]=(h/2)*(f(a,l,t)+f(b,l,t));求解T0(0)3.求解T0(0)到T0(k)的值由公式h=b−an 及h=b−a2k可得n=2k其中h为步长,n为二分次数又由递推公式:T2n=12T n+ℎ2∑f(xk+12)n−1k=0得T2k+1=12T2k+b−a2k+1∑f[a+b−a2k+1(2i−1)]2k−1i=1,k=0,1,2,3~其中xk+12= a+ℎ2(2i−1)。

计算方法 实验报告 拉格朗日 龙贝格 龙格库塔

计算方法 实验报告 拉格朗日 龙贝格  龙格库塔
agrange.c,romberg.c和runge.c三个子文件包含(#include)在了main.c中,所以,要查看运行结果,老师只需直接用win-tc打开main.c,运行即可。感谢使用!
主界面:
/*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;

romberg龙贝格算法实验报告

romberg龙贝格算法实验报告

Romberg龙贝格算法实验报告课程实验报告课程名称:专业班级:CS1306班学号:U20**14967姓名:段沛云指导教师:报告日期:计算机科学与技术学院目录1 实验目的 .. 13 算法设计与流程框图 .. 76 结果分析 .. 71 实验目的掌握Romberg公式的用法,适用范围及精度,熟悉Romberg算法的流程,并能够设计算法计算积分31得到结果并输出,Romberg龙贝格算法实验报告。

1x2 实验原理1 取k=0,h=b-a,求T0=数)。

2 求梯形值T03.1 先求T3.2再由公式T 求T 4-14-13.3 在求出的同时比较T与T的大小,如果二者之差的绝对值小于1e-5,就停止求T;此时的k就是所求的二分次数,而此时的T就是最终的结果3.4 打印出所有的T;程序流程图4 源程序#inlude #inlude #inlude #inlude int main{float f{float ;=1x;return ; }float a,b,e,h,s,k,x,T1=0,T2=0,S1=0,S2=0,C1=0,C2=0,R1=0,R2=0;inti=0; printf{T1=T2; S1=S2; h=h2; k=k+1; }elseif{C2=S2+{R2=C2+{printf;break;} else{R1=R2;C1=C2;T1=T2;S1=S2;h=h2;k=k+1; } }i++;printf;}hile;sstem;return 0; }5 程序运行6 结果分析如上所示的结果与课本中求得的结果完全一样,表明程序编写正确,且符合要求,事实上,只要再将所求值的精度设置得更小,则所求的结果将更加准确,最终将无限接近于标准值,由上表也可以看出用龙贝格积分法求函数的积分值在精度比较低的情况下就能求到很准确的值!7 实验体会本次实验较为简单,主要时间是耗费在循环判断上面,因为书上已经给了流程图,都是基本的C语言,难度不大,工作报告《Romberg龙贝格算法实验报告》。

龙贝格算法-数值分析-实验报告

龙贝格算法-数值分析-实验报告

主要参考资料
《数值分析简明教程》 《数值分析课程设计》
王能超 陈越、童若锋
高等教育出版社 浙江大学出版社
教研室签字: 年 月 日
院签字: 年 月 日
2 / 13
目录
目录 ......................................................................................................................................................... 3 一、任务综述于分析 ............................................................................................................................. 4
3 / 13
一、任务综述于分析
1.1 任务
1、 根据所掌握龙贝格算法,独立编写程序,实现具体问题的求解等; 2、 自己设计一个多项式,根据设定的代数多项式进行测试,调式程序; 3、 对所设计的程序解决实际问题并分析,撰写分析报告
1.2 分析
1、 首先先明确什么是龙贝格算法,了解他的基本原理,画出流程图,编写程序代码。 2、 自己设计一个多项式,对他利用龙贝格算法求出他的积分。并与真实值比较,检测程序是否有错误。 并进一步修改完善程序。可以多试几个多项式,来判断程序的正确性。 3、 任务三实际上是求平面曲线长的问题。利用弧长微分公式可求得结果。
)
(2)
4 / 13
直接验证易知 (3) 这就是说,用梯形法二分前后两个积分值 法的积分值 。 与 按式(3)做线性组合,结果得到辛普生

龙贝格积分实验报告

龙贝格积分实验报告

二、Romberg 积分法1.变步长Romberg 积分法的原理复化求积方法对于提高精度是行之有效的方法,但复化公式的一个主要缺点在于要事先估计出部长。

若步长过大,则精度难于保证;若步长过小,则计算量又不会太大。

而用复化公式的截断误差来估计步长,其结果是步长往往过小,而且''()f x 和(4)()f x 在区间[,]a b 上的上界M 的估计是较为困难的。

在实际计算中通常采用变步长的方法,即把步长逐次分半(也就是把步长二等分),直到达到某种精度为止,这种方法就是Romberg 积分法的思想。

在步长的逐步分半过程中,要解决两个问题: 1. 在计算出N T 后,如何计算2N T ,即导出2N T 和N T 之间的递推公式; 2. 在计算出N T 后,如何估计其误差,即算法的终止的准则是什么。

首先推导梯形值的递推公式,在计算N T 时,需要计算1N +个点处的函数值在计算出NT 后,在计算2N T 时,需将每个子区间再做二等分,共新增N 个节点。

为了避免重复计算,计算2N T 时,将已计算的1N +个点的数值保留下来,只计算新增N 个节点处的值。

为此,把2N T 表示成两部分之和,即由此得到梯形值递推公式 因此由复化梯形公式的截断误差有若''()f x 变化不大时,即''''12()()f f ηη≈,则有式(2)表明,用2N T 作为定积分I 的近似值,其误差大致为21()3N N T T -,因此其终止条件为其中ε是预先给定的精度。

2.Romberg 积分公式将上述方法不断推广下去,可以得到一个求积分的序列,而且这个序列很快收敛到所求的定积分。

记(0)N N T T =,将区间N 等分的梯形值。

(1)N N T S =,将区间N 等分的Simpson (2)N N T C =,将区间N 等分的Cotes 。

(3)N N T R =,将区间N 等分的Romberg 。

实验报告-龙贝格算法

实验报告-龙贝格算法

ans =
0.9207 0.9398 0.9445 0.9457 0.9460 0.9461 0.9461 0 0 0 0 0 0
0 0.9461 0.9461 0.9461 0.9461 0.9461 0.9461 0 0 0 0 0 0
0 0 0.9461 0.9461 0.9461 0.9461 0.9461 0 0 0 0 0 0
0 0 0 0.9461 0.9461 0.9461 0.9461 0 0 0 0 0 0
0 0
0 0
0 0
0 0
(2) 、计算从 f(x) = 4/(1+x^2)0 到 1 的积分的结果为: >> f=inline('4/(1+x^2)','x')
f=
Inline function: f(x) = 4/(1+x^2)
2. 解决问题: (1)计算 sin(x)/x 从 10^(-100)到 1 的积分: >> f=inline('sin(x)/x','x')
f=
Inline function: f(x) = sin(x)/x
>> romberg(f,10^(-100),1,5*10^(-7)) 答案 0.94608
0 0 3.1421 3.1416 3.1416 3.1416 3.1416 0 0 0 0 0 0 0
0 0 0 3.1416 3.1416 3.1416 3.1416 0 0 0 0 0 0 0
0 另一种代码:
0
0

0
计算从 f(x) = 4/(1+x^2)0 到 1 的积分的结果为:

龙贝格算法报告

龙贝格算法报告

龙贝格算法报告 中国地质大学姓名:马立华 班号:121112 学号:20111000536 (一)、龙贝格算法的推导工程的形成。

()()baI f f x dx =⎰n 等份,分点1、梯形公式是将区间[,a b ⎤⎦ 有为kx =a+kh ,h=b a n- ,k=0,1,2,3....在每个小区间上采用梯形公式为110[()()]2n i i i h T f x f x -+==+∑.通过计算复合梯形公式的余项为. (梯形公式是把大区间分成一些小区间,通过梯形面积近似曲边提醒的体积)2、梯形公式的递推法,复合梯形公式在区间[,a b ⎤⎦ 上,共有n 等份,则共有n+1个分点,将求积区间在二分一次,则分点增加至2n+1个,由梯形公式()()()⎥⎦⎤⎢⎣⎡+-=∑-=1122n k k n b f x f a f h T 得式()∑-=++=121222n k k n nx f h T T 求出n T 2;(梯形公式是通过增加区间每一次都是将上一次区间加倍)3,外推技巧:由1,2,可以看出梯形公式收敛地速度较慢,并且精确地阶数较低由余项得 h=b a n -,。

(外推技巧是通过将梯形公式的递推法得到的公式进行加减乘除法提高其精确度。

) 定理四 设](),f x c a b ∞⎡∈⎣,则有2462123()........l l T h I h h h h αααα=++++++ (1)其中2"()12n b a R h f η-=-]2''(),,12n b a I T h f a b ηη-⎡-=-∈⎣3201x x dx+⎰系数l α(l=1,2,3....)与h 无关。

通过定理四 则有24212().......() (24162)l l h h h h T I ααα=+++++。

(2),通过(1)(2)得46124(/2)()() (3)T h T h S h I h h ββ-==+++发现我们只是通过简单地加减法,就将误差的阶数提高了两阶。

Romberg算法实验

Romberg算法实验

数值计算实验—实验报告3一、实验项目:Romberg 算法 二、实验目的和要求通过本实验可以使学生理解如何在计算机上使用数值方法计算定积分()ba f x dx ⎰的近似值。

三、实验内容h=b-a; k=1;T(1,1)=(b-a)*(f(a)+f(b))/2; err=1;while(err>eps) new=0;for j=1:2^(k-1); x=a+(2*j-1)*h/2; new=new+f(x); endT(k+1,1)=(T(k,1)+h*new)/2; for m=1:kT(k+1,m+1)=(4^m*T(k+1,m)-T(k,m))/(4^m-1); enderr=abs(T(k+1,m+1)-T(k,m)); k=k+1; h=h/2; endk,T,m(1)将上述程序录入计算机,并进行调试(2)用调试好的程序计算10sin()/x xdx ⎰的近似值,取精度为eps=10^-8。

步骤一:先编制计算函数值的程序function z=f(x) if(x!=0)z=sin(x)/x; else z=1; end步骤二:执行上述Romberg 算法,输出的T —数表为(参考)(3)用下列方法计算积分2(1(^2)*(cos())^2))^(0.5)b t dx π-⎰的近似值,其中b分别取为0.1,0.3,0.5.0.7,0.9。

(a) Rombergfunction z=f(x)b=0.1;z=sqrt(1-b*b*cos(x)*cos(x));end(b)复合的Simpson公式(n=16)h=(b-a)/16;X0=f(a)+f(b);X1=0;X2=0;for i=1:15x=a+i*h;if(rem(i,2)==0)X2=X2+f(x);elseX1=X1+f(x);endendX=h*(X0+4*X1+2*X2)/3;b=0.1,X= 6.26744776809318b=0.3,X= 6.13933390904750b=0.5,X= 5.86985288646431b=0.7,X= 5.42275185873374b=0.9,X= 4.68936847260762四、实验总结Romberg方法也称为逐次分半加速法。

龙贝格方法

龙贝格方法

实验三 龙贝格方法【实验内容】1. 理解龙贝格方法的基本思路2. 用龙贝格方法设计算法,编程求解一个数值积分的问题。

【实验方法与步骤】(对于必须编写计算机程序的实验,要附上学生自己编写的程序)1.龙贝格方法的基本思路龙贝格方法是在积分区间逐次二分的过程中,通过对梯形之值进行加速处理,从而获得高精度的积分值。

2.龙贝格方法的算法步骤1 准备初值 和,用梯形计算公式计算出积分近似值()f a ()f b ()()12b a T f a f b -=+⎡⎤⎣⎦ 步骤2 按区间逐次分半计算梯形公式的积分近似值令,2i b a h -=0,1,2,...i = 计算,12102122n n n i i h T T f x -+=⎛⎫=+ ⎪⎝⎭∑2in =步骤3 按下面的公式积分(为便于编程,写下列形式)梯形公式:()223n n n n T T S T -=+ 辛普生公式:()2215n n n n S S C S -=+ 龙贝格公式:()2263n n n n C C R C -=+ 步骤4 精度控制当,(为精度)时,终止计算,并取为近似值,否则,将2n n R R ε-<ε2n R 步长折半,转步骤2。

#include<stdio.h>#include<math.h>float f(float ){return 4/(1+x*x);}float Romberg(float a,float b,float( * f)(float),float eps) {int n=1,k;float h=b-a,y,temp;float T1,T2,S1,S2,C1,C2,R1,R2;T1=(b-a)/2*(( * f)(a)+( * f)(b));while(1){temp=0;for(k=0;k<=n-1;k++){ y=a+k*h+h/2;temp+=( * f)(y);}T2=(T1+temp*h)/2;if(fabs(T2-T1)<eps) return T2;S2=T2+(T2-T1)/3;if(n==1){T1=T2;S1=S2;h/=2;n*=2;continue;}C2=S2+(S2-S1)/15;if(n==2){C1=C2;T1=T2;S1=S2;h/=2;n*=2;continue;}R2=C2+(C2-C1)/63;if(n==4){R1=R2;C1=C2;T1=T2;S1=S2;h/=2;n*=2;continue;} if(fabs(R2-R1)<eps) return R2;R1=R2;C1=C2;T1=T2;S1=S2;h/=2;n*=2;}}void main(){ double eps=5e-6;printf("R=%f",Romberg(0.0,1.0,f,eps));}【实验结果】。

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

Romberg龙贝格算法实验报告
2017-08-09
课程实验报告
课程名称:
专业班级: CS1306班学号: U201314967 姓名:段沛云指导教师:报
告日期:
计算机科学与技术学院
目录
1 实验目的 (1)
2 实验原理 (1)
3 算法设计与流程框图 (2)
4 源程序 (4)
5 程序运行 (7)
6 结果分析 (7)
7 实验体会 (7)
1 实验目的
掌握Romberg公式的用法,适用范围及精度,熟悉Romberg算法的流程,并能够设计算法计算积分
31
得到结果并输出。

1x
2 实验原理
2.1 取k=0,h=b-a,求T0=
数)。

2.2 求梯形值T0(
b-a
),即按递推公式(4.1)计算T0。

k
2
h
[f(a)+f(b)],令1→k,(k记区间[a,b]的二分次2
2.3 求加速值,按公式(4.12)逐个求出T表的第k行其余各元素Tj(k-j)
(j=1,2,….k)。

2.4 若|Tk+1-Tk|
n-1
11T2n=[Tn+hn∑f(xi+)]
22i=0
1
Sn=T2n+(T2n-Tn)
31
Cn=S2n+(S2n-Sn)
151
Rn=C2n+(C2n-Cn)
63
3 算法设计与流程框图
算法设计:(先假定所求积分二分最大次数次数为20) 3.1 先求T[k][0] 3.2 再由公式T
(k)m
4m(k+1)1)=mTm-1-mTm(k-1(k=1,2,) 求T[i][j] 4-14-1
3.3 在求出的同时比较T[k][k]与T[k-1][k-1]的大小,如果二者之差的绝对
值小于1e-5,就停止求T[k][k];此时的k就是所求的二分次数,而此时的T[k][k]就是最终的结果 3.4 打印出所有的T[i][j];程序流程图
#include #include #include #include int main(void) {
float f(float(x)) {
float y; y=1/x; return y; }
float a,b,e,h,s,k,x,T1=0,T2=0,S1=0,S2=0,C1=0,C2=0,R1=0,R2=0; int i=0;
printf("请输入积分下限 : "); scanf("%f",&a);
printf("\n请输入积分上限 :"); scanf("%f",&b);
printf("\n请输入允许误差 :"); scanf("%f",&e); k大学网=1; h=b-a;
T1=h*(f(a)+f(b))/2;
printf("____________________________________________\n"); printf("计算结果如下 : \n");
printf("\nk T2 S2 C2 R2\n");
printf("%d %10.7f %10.7f %10.7f %10.7f\n",i,T1,S1,C1,R1); do {
x=a+h/2; s=0; while(x
{ s=s+f(x); x=x+h; }
T2=(T1+s*h)/2; S2=T2+(T2-T1)/3; if(k==1) {
T1=T2; S1=S2; h=h/2; k=k+1; }
else if(k==2) {
C2=S2+(S2-S1)/15; C1=C2; T1=T2; S1=S2; h=h/2; k=k+1; }
else if(k==3) {
R2=C2+(C2-C1)/63; C2=S2+(S2-S1)/15; C1=C2; T1=T2; S1=S2; h=h/2;
k=k+1; } else {
C2=S2+(S2-S1)/15;
R2=C2+(C2-C1)/63; if(fabs(R2-R1)
printf("%d %10.7f %10.7f %10.7f %10.7f\n",i+1,T2,S2,C2,R2);
} else { R1=R2; C1=C2; T1=T2; S1=S2; h=h/2; k=k+1; } } i++;
printf("%d %10.7f %10.7f %10.7f %10.7f\n",i,T2,S2,C2,R2); }
while(1); system("pause"); return 0; }
5 程序运行
6 结果分析
如上所示的结果与课本中求得的结果完全一样,表明程序编写正确,且符合要求,事实上,只要再将所求值的精度设置得更小,则所求的结果将更加准确,最终将无限接近于标准值,由上表也可以看出用龙贝格积分法求函数的积分值在精度比较低的.情况下就能求到很准确的值!
7 实验体会
本次实验较为简单,主要时间是耗费在循环判断上面,因为书上已经给了流程图,都是基本的C语言,难度不大。

过程中唯一遇到的一点障碍就是在写循环判断时由于多重判断多重循环导致混乱,幸好最后改
正了,最后得到的结果经检验与给定的结果相同。

通过这次实验上机,
使我更进一步了解了龙贝格法的计算思想,其在精度上很有保证,收敛较快,是解积分问题的有效方法。

相关文档
最新文档