龙贝格积分实验报告
- 格式:doc
- 大小:110.95 KB
- 文档页数:5
3 龙贝格积分3.1 算法原理及程序框图龙贝格积分法是在复化梯形求积公式、复化辛普森求积公式和复化科茨求积公式关系的基础上,构造出的一种精度更高的数值积分方法。
对于复化梯形求积公式而言,近似积分为()2221[]41n n n n I f T T T T ≈+-=-.(11) 对于复化辛普森求积公式和复化科茨求积公式而言,也有类似的关系,如公式(12)和公式(13)。
()22221[]41n n n n I f S S S S ≈+-=- (12)()22231[]41n n n n I f C C C C ≈+-=- (13)通过对公式(11)~(13)做进一步分析,可得到公式(14)和公式(15)。
()22141n n n n S T T T =+--(14)()222141n n n n C S S S =+-- (15)根据公式(14)和公式(15)表现出来的规律,令龙贝格积分为()223141n n n n R C C C =+-- (16)其截断误差为c R h 8f (8)(η),已经具有很高的精度。
龙贝格积分法是将区间[a , b ]逐次分半进行计算,因此,对已知函数f (x )在区间[a , b ]上的龙贝格积分法的计算公式的算法如下,程序框图如图13所示。
(1) 计算T 1:[]1()()2b aT f a f b -=+;(2) 逐次计算T 2k +1:()1211221121,0,1,2222kk k k k i b a b a T T f a i k +++=--⎛⎫=++-= ⎪⎝⎭∑;(3) 逐次计算S 2k 、C 2k 和R 2k :()()()11111122222222232222141141141kk k k k k k k k k k k S T T T C S S S R C C C ++++++⎧=+-⎪-⎪⎪=+-⎨-⎪⎪=+-⎪-⎩;(4) 若122k k R R ε+-<,则取[]12k I f R +≈;否则,继续计算,直到满足精度为止。
计算方法上机实习题目(四)——龙贝格算法计算椭圆周长一、 题目用龙贝格算法计算椭圆110040022=+y x 的周长,使误差不超过410-。
二、解题方法由题目中椭圆的标准方程 110040022=+y x 知,该椭圆的参数方程为)2,0[ sin 10cos 20π∈⎩⎨⎧==t ty t x 参考微积分中曲线弧长计算公式dt t y t x s ⎰+=βα)()(2'2' 知该椭圆的周长计算公式为 dt t t s ⎰+=π2022cos 100sin 400 故该问题为利用龙贝格算法计算数值积分。
参考教材式(5.42),若以一个二元数组T[i][j]代表式中的)(j i T ,则式(5.42)可化为⎪⎪⎪⎩⎪⎪⎪⎨⎧=⋯=--=--+-+=+-=-+-=-∑-)3,2,1;,2,1,0(144]2)12([221)]()([2][]1[]1[]1[][][21]1[]0[][]0[]0[]0[1m k T T T a b i a f a b T T b f a f a b T m k m k m m k m i l l l l l 其中t t t f 22cos 100sin 400)(+=,a ,b 为区间端点值,π2,0==b a 。
相应的计算程序段为T[0][0]=(b-a)*(f(a)+f(b))/2;k=1;T[0][1]=T[0][0]/2+(b-a)*Sum(1)/pow(2,1);T[1][0]=(4*T[0][1]-T[0][0])/(4-1);k=2;T[0][2]=T[0][1]/2+(b-a)*Sum(2)/pow(2,2);T[1][1]=(4*T[0][2]-T[0][1])/(4-1);T[2][0]=(pow(4,2)*T[1][1]-T[1][0])/(pow(4,2)-1);k=3;T[0][k]=T[0][k-1]/2+(b-a)*Sum(k)/pow(2,k);T[1][k-1]=(4*T[0][k]-T[0][k-1])/(4-1);T[2][k-2]=(pow(4,2)*T[1][k-1]-T[1][k-2])/(pow(4,2)-1);T[3][k-3]=(pow(4,3)*T[2][k-2]-T[2][k-3])/(pow(4,3)-1);for(k=4;fabs(T[3][k-4]-T[3][k-5])>=pow(10,-4);k++){T[0][k]=T[0][k-1]/2+(b-a)*Sum(k)/pow(2,k);T[1][k-1]=(4*T[0][k]-T[0][k-1])/(4-1);T[2][k-2]=(pow(4,2)*T[1][k-1]-T[1][k-2])/(pow(4,2)-1);T[3][k-3]=(pow(4,3)*T[2][k-2]-T[2][k-3])/(pow(4,3)-1);}其中f(x)和Sum(l)为调用的子函数,子函数程序如下:double Sum(int l){i;intdoublesum,a,b;double f(double x);a=0.0;b=2*PI;sum=0;for(i=1;i<=pow(2,l-1);i++)sum=sum+f(a+(2.0*i-1.0)*(b-a)/pow(2,l));sum;return}double f(double x){y;doubley=sqrt(pow(20*sin(x),2)+pow(10*cos(x),2));y;return}因为并不是一开始T[3][j]就同步出现,所以需将k=3之前的各值先单独算出,又因为要10-,所以控制循环结束的条件是fabs(T[3][k-4]-T[3][k-5])>=pow(10,-4),求最后误差不超过4其中fabs(x)为求绝对值的函数,这里需要注意的是,每次判断时已经执行了k++这条指令,所以在判断是应该是T[3][k-4]-T[3][k-5],而不是T[3][k-3]-T[3][k-4],我在最初编写程序时就忽略了这个问题,导致花费很久调试程序。
龙贝格实验报告龙贝格实验报告龙贝格实验是一项经典的心理学实验,旨在探究人类对于延迟满足的选择行为。
实验的设计非常简单,但却能揭示出人类的心理特点和决策模式。
本文将对龙贝格实验进行介绍,并探讨实验结果对于我们日常生活中的决策行为的启示。
龙贝格实验最早由德国心理学家沃尔夫冈·龙贝格于1971年提出,他将实验对象置于一个选择情境中。
在这个情境中,实验对象需要在两个选项之间进行选择,一方面是能够立即获得一定的奖励,另一方面是需要等待一段时间才能获得更大的奖励。
实验对象可以根据自己的意愿选择不同的选项。
实验结果显示,大多数人更倾向于选择立即获得的小额奖励,而不是等待更长时间以获得更大的奖励。
这种行为模式被称为“即时满足偏好”。
这一结果揭示了人类的一种心理特点:我们更倾向于追求即时的满足感,而忽视了长期利益。
那么,为什么人们会有这种即时满足偏好呢?心理学家提出了一些解释。
首先,人们往往更容易预测和评估即时奖励,而对于未来的奖励往往难以准确评估。
这导致了人们对于即时奖励的价值更加敏感,而对于未来奖励的价值更加模糊。
其次,人们对于即时奖励的需求更加迫切,因为它们能够立即满足我们的欲望和需求,而未来奖励则需要等待和付出更多的努力。
然而,即时满足偏好并不总是明智的选择。
在现实生活中,我们常常需要做出一系列决策,这些决策涉及到长期利益和短期满足之间的权衡。
例如,我们是否应该立即购买心仪已久的商品,还是将钱储蓄起来以应对未来的不确定性?我们是否应该立即享受美食,还是节制自己以保持健康的身体?这些决策都需要我们权衡即时满足和长期利益之间的关系。
为了克服即时满足偏好的影响,我们可以采取一些策略。
首先,我们可以通过增强对未来奖励的预期和评估来提高对其价值的认识。
这意味着我们需要更加清晰地设定未来奖励的目标,并了解其对我们个人发展和幸福的重要性。
其次,我们可以通过设定明确的规划和目标来引导自己的行为。
例如,我们可以设定一个长期的储蓄计划,以确保我们能够在未来获得更大的回报。
第1篇一、实验目的本次实验旨在通过数值分析的方法,研究几种常见的数值积分方法,包括梯形法、辛普森法、复化梯形法和龙贝格法,并比较它们在计算精度和效率上的差异。
通过实验,加深对数值积分理论和方法的理解,提高编程能力和实际问题解决能力。
二、实验内容1. 梯形法梯形法是一种基本的数值积分方法,通过将积分区间分割成若干个梯形,计算梯形面积之和来近似积分值。
实验中,我们选取了几个不同的函数,对积分区间进行划分,计算积分近似值,并与实际积分值进行比较。
2. 辛普森法辛普森法是另一种常见的数值积分方法,它通过将积分区间分割成若干个等距的区间,在每个区间上使用二次多项式进行插值,然后计算多项式与x轴围成的面积之和来近似积分值。
实验中,我们对比了辛普森法和梯形法的计算结果,分析了它们的精度差异。
3. 复化梯形法复化梯形法是对梯形法的一种改进,通过将积分区间分割成多个小区间,在每个小区间上使用梯形法进行积分,然后计算所有小区间积分值的和来近似积分值。
实验中,我们对比了复化梯形法和辛普森法的计算结果,分析了它们的精度和效率。
4. 龙贝格法龙贝格法是一种通过外推加速提高计算精度的数值积分方法。
它通过比较使用不同点数(n和2n)的积分结果,得到更高精度的积分结果。
实验中,我们使用龙贝格法对几个函数进行积分,并与其他方法进行了比较。
三、实验步骤1. 编写程序实现梯形法、辛普森法、复化梯形法和龙贝格法。
2. 选取几个不同的函数,对积分区间进行划分。
3. 使用不同方法计算积分近似值,并与实际积分值进行比较。
4. 分析不同方法的精度和效率。
四、实验结果与分析1. 梯形法梯形法在计算精度上相对较低,但当积分区间划分足够细时,其计算结果可以接近实际积分值。
2. 辛普森法辛普森法在计算精度上优于梯形法,但当积分区间划分较细时,计算量较大。
3. 复化梯形法复化梯形法在计算精度上与辛普森法相当,但计算量较小。
4. 龙贝格法龙贝格法在计算精度上优于复化梯形法,且计算量相对较小。
姓名 班级 学号【实验目的】:1通过本实验加深对“龙贝格公式”的认识2通过c++编写程序求解下列积分(E=1E-6)746824.02659330.0sin 11212≈≈⎰⎰-dx e dx xxx 、、3预期达到题中结果4注意实验给出的程序有两个子程序为:f f1使用时只需该程序名即可。
【实验题目】:复合求积公式计算定积分【实验原理及理论基础】:用梯形公式得到的积分近似值n T 2的误差是)(312n n T T -,因此,人们希望用这个误差作为对n T 2的一种补偿,则得到求积公式n n n n n baT T T T T dx x f 3134)(31)(222-=-+≈⎰的代数精度会有所提高。
通过直接验证可知n n n T T S 31342-=梯形公式二分前后的两个积分值n T 与n T 2按如上公式做线性组合其结果正好是抛物线公式得到的积分值n S同理可得到科茨公式积分值n Cn n n S S C 15115162-=和龙贝格公式n n n C C R 63163642-=计算⋅⋅⋅421,,T T T()()()[]312411212212212)(2)()(x f x f ab T T x f ab T T b f a f a b T +-+=-+=⎥⎦⎤⎢⎣⎡+-=【实验内容】:龙贝格求积的计算步奏如下: (1)计算)(a f ,)(b f ,算出1T ; (2)把[a,b]2等分,计算⎪⎭⎫⎝⎛+2b a f ,算出2T 与1S ; (3)把[a,b]4等分,计算⎪⎭⎫ ⎝⎛-+4a b a f ,⎪⎭⎫ ⎝⎛-⨯+43a b a f ,算出4T 、2S 与1C ; (4)把[a,b]8等分,计算7,5,3,1,8=⎪⎭⎫⎝⎛-+i a b ia f ,算出8T 、4S 、2C 与1R ;(5)把[a,b]16等分,计算15,,5,3,1,8⋅⋅⋅=⎪⎭⎫⎝⎛-+i a b ia f ,算出16T 、8S 、4C 与2R ,继续重复进行,直到ε<-n n R R 2(允许误差)时停止运算,2R 就是所求积分值。
数值积分算法实验报告实验目的:验证数值分析理论实验工具:matlab摘要:本实验通过对比龙贝格积分算法和三点,五点勒让德高斯求积公式,对数值分析理论进行验证。
文章中提供了用matlab写的小程序即运行结果。
我们可以看到,龙贝格积分算法方便快捷且效率较高。
三点积分算法不足以满足精度要求,五点法基本可以满足,但是如果追求更高精度,则五点法改进困难,而龙贝格则可适用于任何精度的计算。
高斯求积的效率是比较高的,只需五点就可以达到非常高的精度。
一matlab源程序1.%龙贝格积分算法,jd表示相对精度function z=lbgjf(a,b,jd)h=(b-a);TT(1,1)=h.*(f(b)+f(a))/2;k=2;TT(1,2)=TT(1,1)./2+h/2.*f(a+h/2);TT(2,1)=TT(1,2).*4/3-TT(1,1)./3;z=TT(2,1);while abs((TT(k,1)-TT(k-1,1))./TT(k,1))>=jdk=k+1;h=h./2;for j=1:2.^(k-2)ff(1,j)=f(a+h*(j-1/2));endfff=sum(ff).*h/2;TT(1,k)=TT(1,k-1)./2+fff;for j=2:kTT(j,k-j+1)=4^(j-1).*TT(j-1,k-j+2)./(4^(j-1)-1)-TT(j-1,k-j+1)/(4^(j-1 )-1);z=TT(j,k-j+1);endend2 %五点法·¨function z=fivedlrd(a,b)Ak=[0.2369269 0.4786287 0.5688889 0.4786287 0.2369269];xk=[-0.9061798 -0.5384693 0 0.5384693 0.9061798];for i=1:5ff(i)=Ak(i).*f((b-a).*xk(i)./2+(b+a)./2);endz=(b-a)./2.*sum(ff)3%三点法¨function z=threedlrd(a,b)Ak=[0.5555556 0.8888889 0.5555556];xk=[-0.7745967 0 0.7745967];for i=1:3ff(i)=Ak(i).*f((b-a).*xk(i)./2+(b+a)./2); endz=(b-a)./2.*sum(ff)4%另用一matlab文件来表示函数function y=f(x)y=sin(x).^10;end二实验过程1 f(x)=x.^2,2 f(x)=x.^3可以看到,对于这种简单函数三种方法基本没有区别3 f(x)=sin(x).^10.lbgjf(2,3,10e-6)ans =0.0608>> fivedlrd(2,3)z =0.0608ans =0.0608>> threedlrd(2,3)z =0.0618ans =0.06184f(x)= exp(exp(sin(x).^10))threedlrd(2,3)z =2.9278ans =2.9278>> fivedlrd(2,3)z =2.9304ans =2.9304>> lbgjf(2,3,10e-6)ans =2.9304>>可见三点法已经不再满足要求。
实验8 数值积分实验实验类型:●验证性实验○综合性实验○设计性实验实验目的:进一步熟练掌握变步长数值积分算法,提高编程能力和解决定积分问题的实践技能。
实验内容:用龙贝格积分算法计算⎰102dxx实验原理龙贝格积分算法:生成j>=k 的近似积分结果逼近表,并以R(j+1,j+1)为最终解来逼近积分。
R(j,0)=T(j), j>=0,T(j)为区间逐次减半递推梯形求积分公式算出的结果;R(j,1)=S(j), j>=1,S(j)为区间逐次减半递推辛普森求积分公式算出的结果;R(j,2)=B(j), j>=2,B(j)为递推布尔求积分公式算出的结果;1 4)1 ,1()1,()1,(),(-----+-=kkjRkjRkjRkjR实验步骤1 要求上机实验前先编写出程序代码2 编辑录入程序3 调试程序并记录调试过程中出现的问题及修改程序的过程4 经反复调试后,运行程序并验证程序运行是否正确。
5 记录运行时的输入和输出。
实验总结:练习对龙贝格积分的应用,熟悉算法的编写;实验报告:根据实验情况和结果撰写并递交实验报告。
>> [R,quad,err,h]=roomber('f1',0,1,10,0.001)R =0.5000 0 0 0 00.3750 0.3333 0 0 00.3438 0.3333 0.3333 0 00.3359 0.3333 0.3333 0.3333 00.3340 0.3333 0.3333 0.3333 0.3333quad =0.3333err =h =0.0625参考程序function [R,quad,err,h]=roomber(f,a,b,n,tol)M=1;h=b-a;err=1;J=0;R=zeros(4,4);R(1,1)=h*(feval(f,a)+feval(f,b))/2;while((err>tol)&(J<n))|(J<4)J=J+1;h=h/2;s=0;for p=1:Mx=a+h*(2*p-1);s=s+feval(f,x);endR(J+1,1)=R(J,1)/2+h*s;M=2*M;for k=1:JR(J+1,k+1)=R(J+1,k)+(R(J+1,k)-R(J,k))/(4^k-1); enderr=abs(R(J,J)-R(J+1,k+1));endquad=R(J+1,J+1);。
实验二:龙贝格算法一、实验目的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 →∞=这是一个收敛速度更快的一个数值求积公式,我们称为龙贝格积分法。
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窗口中编辑要调用的函数名与指定的函数名字不同导致出现错误,通过改正与函数名相同即可。
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)。
实验名称实验4实验地点6A-XXX 实验类型设计实验学时 2 实验日期20 /X/X ★撰写注意:版面格式已设置好(不得更改),填入内容即可。
一、实验目的
1.龙贝格积分
二、实验内容
1.实验任务
1.龙贝格积分
2.程序设计
1)数据输入(输入哪些数据、个数、类型、来源、输入方式)
const double eps=0.000005;
const int MaxD=100;
double f(double), double a, double b
2)数据存储(输入数据在内存中的存储)
函数
double f1(double x)
double f2(double x)
double f3(double x)
void Romberg(double f(double), double a, double b)
3)数据处理(说明处理步骤。
若不是非常简单,需要绘制流程图)
1.输入要处理的数据进入变量中
2.进行函数处理
3.输出函数处理结果
4)数据输出(贴图:程序运行结果截图。
图幅大小适当,不能太大)。
目录设计总说明 (1)关键词 (1)绪论 (2)第1章课程设计正文 (3)1.1 需求分析 (3)1.2 界面设计 (3)1.3 功能实现 (4)1.4系统算法和程序代码的实现 (5)第2章总结 (7)参考文献 (8)设计总说明学习了《数值分析》和相关开发软件课程并能够熟练掌握各种算法的原理及其求解方法和特点;1、独立进行面向对象软件的分析、设计及在计算机中进行编写代码和调试工作。
培养自己独立进行面向对象软件的分析、面向对象软件的设计及面向对象软件的实现与测试的能力;2、通过课题的理论设计和在计算机中实验调试代码,加深对数值分析知识的理解,培养软件开发的实践技能,提高分析解决具体问题的能力;3、充分利用面向对象程序设计,用精要的语句设计程序,培养独立创新能力;4、熟练运用MATLAB7开发界面美观、清晰,简单,好学,好用的计算软件;说明书中一般包括任务的提出、方案论证、设计与计算说明、试验调试及结果的分析、结束语等内容。
要求理论依据充分,数据准确,公式推导及计算结果正确。
涉及到计算机软件:学生要独立完成一个软件或较大软件中的一个模块,要有足够的工作量;要写出相关的设计文件;能够进行计算机演示和给出运行结果。
在此软件中,我们运用龙贝格方法进行积分。
此外,只需要在编辑框中输入上下限和含x的被积函数,点击计算按钮,得到计算结果。
关键词:MATLAB 数值分析龙贝格积分绪论数值分析是计算数学的一个主要部分,它研究用计算机求解各种数学问题的数值计算方法及其理论与软件实现。
数值分析一个很强的特点就是面向计算机,要根据计算机特点提供切实可行的有效算法。
因此,本次课程设计是在学习了《数值分析》和相关开发软件课程之后,为了能掌握数值分析的基本概念,结合实际的操作和设计,巩固课堂教学内容,使自己能掌握数值分析的基本概念、原理和技术,将理论与实际相结合,应用现有的数据建模工具和数值分析理论,规范、科学地完成一个小型可视化软件的设计与实现,把理论课与实验课所学内容做一综合,通过实际项目的设计、开发,培养学生独立进行计算软件的建模、在计算机中进行计算、设计、通过相关软件开发系统,并在此基础上强化学生的实践意识、提高其实际动手能力和创新能力。
二、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 。
龙贝格算法报告 中国地质大学姓名:马立华 班号: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 积分法
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 +个点处的函数值在计算出N T 后,在计算2N T 时,需将每个子区间再做二等分,共新增N 个节点。
为了避免重复计算,计算2N T 时,将已计算的1N +个点的数值保留下来,只计算新增N 个节点处的值。
为此,把2N T 表示成两部分之和,即
21
2221
122211
122211
1
[()()2()]21
[()()2()2((21))]21[()()2()]((21))]22N N
N N k N N N N N k k N N
N N N N k k T h f a f b f a kh h f a f b f a kh f a k h h f a f b f a kh h f a k h -=-==-===+++=+++++-=+++++-∑∑∑∑∑ 由此得到梯形值递推公式
2221
1
((21))12N
N
N N N k T T h f a k h ==++-∑
因此
1
11212122,[()()],211
,()
22h h b a T f a f b h h T T h f a h =-=
+==++
由复化梯形公式的截断误差有
2''
112''2222(),12(),12
N N N
N b a I T h f a b b a I T h f a b
ηηηη--=-≤≤--=-≤≤ 若''()f x 变化不大时,即''''12()()f f ηη≈,则有
22241
()
2413
N N N N N T T I T T T -≈
=+--
式(2)表明,用2N T 作为定积分I 的近似值,其误差大致为21
()3
N N T T -,
因此其终止条件为
2N 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 。
由其可构造一个序列(){}k N T ,次序列称为Romberg 序列,并满足如下递推关系:
(0)
(0)(0)1
211[()()],((21)),2222N N N k b a b a b a T
f a f b T T f a k N N
=---=+=++-∑ (1)(1)
24,1,2,41k k k k N N
N
k
T T T k ---==
-
以上递推公式就是Romberg 积分递推公式。
3.Romberg 积分程序
1. 置1N =,精度要求ε,1h b a =-;
2. 计算(0)1[()()]2
b a
T f a f b -=
+;
3. 置22
N N
h h =,并计算(0)
(0)211((21))222N N N k b a b a T T f a k N N =--=++-∑; 4. 置,2,1;M N N N K ===
5. 计算(1)(1)2441
k k k k M M
M
k
T T T ---=-; 6. 若 1M =,则转(7);否则置2
M
M =
,1k k =+转(5); 7. 若()(1)11k k T T ε--≤,则停止计算(输出()1k T ),否则转(3)。
4.Romberg 积分法的应用
function [T,n] = romb(f,a,b,eps) double R ;
if nargin<4,eps=1e-8; end
h=b-a;R(1,1)=(h/2)*(feval(f,a)+feval(f,b)); n=1;J=0;err=1; while (err>eps)
J=J+1;h=h/2;S=0; for i=1:n
x=a+h*(2*i-1); S=S+feval(f,x); end
R(J+1,1)=R(J,1)/2+h*S; for k=1:J
R(J+1,k+1)=(4^k*R(J+1,k)-R(J,k))/(4^k-1); end
err=abs(R(J+1,J+1)-R(J+1,J)); n=2*n; end R;
T=R(J+1,J+1) End
其中输入项:f 为被积函数,ab 为积分区间的端点值,ep 为积分精度;输出项:T 是逐次积分表值,n 是迭代次数,R 是最后积分值。
4.1程序调用
可以将被积分函数编成函数文件,也可以直接使用联函数来表示被积分函数,示例如下:
>>f=inline('1/(1+x.^2)','x');
>> [T,n,R]=romb(f,2,9,1e-9)
运行后得出其迭代次数,最终积分结果以及龙贝格积分矩阵如表2-1所示, 迭代次数N=64,最终的积分值R=0.3530.
表2-1 龙贝格积分矩阵
3.课本例题求解
2
111
2
0000
1ln(1)ln(1)sin()(1),(2),(3),(4)1+1x x x dx dx dx dx x x x x π
+++⎰⎰⎰⎰
1 当迭代精度ep=1e-9的条件下,迭代次数N=32,迭代结果R=0.6931
表2-2 式1对应的龙贝格积分矩阵
2 当迭代精度ep=1e-9的条件下,迭代次数N=32,迭代结果R=0.2722.
表2-3 式2对应的龙贝格积分矩阵
3对于积分1
ln(1)
x dx x +⎰
,由于积分下限0为其奇点,理论上无法进行数值积分,本题中近似取下限为1*10-9来进行计算。
当迭代精度ep=1e-9的条件下,迭代次
数N=16,迭代结果R=0.2722.
表2-4 式3对应的龙贝格积分矩阵
4.对于积分2
0sin()x
dx x
π
⎰,同样积分下限0为积分函数的奇点,理论上无法进行数
值积分运算,本题中仍取积分下限近似为1*10-9进行计算。
当迭代精度ep=1e-9的条件下,迭代次数N=16,迭代结果R=1.3708.
表2-5 式4对应的龙贝格积分矩阵。