XX大
学
综
合设
计
报
告
综合设计五
多方法求解数值积分
学生姓名:
学
号:
年级专业:
指导老师:
学
院:
评阅成绩:
评阅意见:
成绩评定教师签名:
时间:
提交日期:2014年X月
多方法求解数值积分
具体题目要求:用不同数值方法计算积分
(1)
取不同的步长,分别用复合梯形及复合辛普森公式计算积分,给出误差中关于的函数,并与积分精确值比较两个公式的精度,是否存在一个最小的,使得精度不能再被改善?
(2)
用龙贝格求积计算完成问题(1);
(3)
用自适应辛普森积分,使其精度达到。
1设计目的、要求
由积分学基本理论,定积分可由公式计算,但是对于一些无法找到原函数的函数(如等)不能通过牛顿—莱布尼兹公式计算,就必须得另寻它法。因此需要我们能够熟练地应用常用的数值积分计算方法(如机械求积、公式等)并掌握结合数值计算软件(、等)及计算机高级语言进行对应算法实现的技能。
熟练数学软件求解数学问题,掌握各种数学问题的求解方法。本设计主要是通过多种复合求积公式求解积分,主要包括复化梯度法、复化辛普森法、龙贝格以及自适应辛普森法等求解方法,利用软件编写相对应的算法进行求解,大大地提高了解题的速度。
2设计原理
由积分中值定理我们可以知道在积分区间内存在一点,使得式子
成立。这个式子在于对于点的具体位置一般是不知道的,因此难以准确算出的值。也就是不同算法求得平均高度,对应的就是一种不同的数值求积方法。更一般地,我们可以在区间上适当选取某些节点,然后用的加权平均得到平均高度的近似值,这样构造出的求积公式具有下列形式:
称为机械求积公式。
复合梯形公式、复合辛普森公式、龙贝格求积公式以及自适应辛普森公式都以此公式的基础,对积分区间进行变步长的划分求得近似的平均高度值,得到积分函数的近似值。也由于牛顿—柯特斯公式在时不具有稳定性,所以不可能通过提高阶的方法来提高求积精度。为了我提高精度通常可把积分区间等分成为若干个子区间,再在每个子区间上用低价求积公式,这就是复合求积方法。但是这样的积分求解方法也是存在不容忽视的误差。因此需要在设计算法时考虑到算法存在的误差(舍入误差、截断误差等),并对误差作出分析。
3采用软件、设备
软件
4设计内容
第一步:复合梯形公式、复合辛普森公式算法
(一)、复合梯形公式计算积分
复化梯形公式的主要思想是利用若干小梯形的面积代替原方程的积分,利用微元法,可以求出坐标面上由函数与坐标轴围城的图像的面积的近似值,符合了计算机计算存储的思想。下面,我们在探讨复化梯形公式的计算规律:
设将求积区间分成等份,则一共有个分点,按梯形公式
计算积分值,需要提供个函数值。
这里代表步长,分点为,其中
(二)、复合辛普森公式计算积分
算法的基本思想是:把积分区间等分成若干个子区间,而在每一个子区间上用辛普森
求积公式:
得到复合辛普森求积公式:
并且用软件来求解。
第二步:龙贝格算法
考虑积分,欲求其近似值,通常有复化的梯形公式、公式和公式。但是给定一个精度,这些公式达到要求的速度很缓慢。如何提高收敛速度,自然是人们极为关心的课题。为此,记,为将区间进行等分的复化的梯形公式计算结果,记,为将区间进行等分的复化的公式计算结果,记,为将区间进行等分的复化的公式计算结果。根据外推加速方法,可以得到收敛速度较快的积分法。其具体的计算公式为:
1、准备初值,计算
2、按梯形公式的递推关系,计算
3、按Romberg积分公式计算加速值,4、精度控制。对给定的精度,若
则终止计算,并取为所求结果;否则返回2重复计算,直至满足要求的精度为止。
第三步:自适应辛普森算法
复合求积方法通常适用于被积函数变化不太大的积分,如果在积分区间被积函数变化很大有的部分函数值变化剧烈而有的部分则是变化平缓,如果此时还是将积分区间等分后用复合求积公式的话计算量很大。而采用针对被积函数在区间上的不同情形而采用不同的步长,使得在满足精度的前提下积分计算量减少,这就是自适应积分方法,能自动地在被积函数变化剧烈的区域增多节点,而在被积函数变化平缓的地方减少节点。因此它是一种不均匀区间的积分方法。题目要求使相邻两个区间的误差达到一定的要求,即用自适应辛普森公式来求积分,先算出积分区间的左右端点函数值,求出区间中点函数值与左右端点的函数差值,再与所要求的精度比较,不满足的对所在区间二等分,接着算出每个子区间端点的函数值判断时否符合精度要求,直到积分每个子区间内都满足精度要求,最后所得各个区间端点的函数值之和即为积分的近似值。
第四步:误差余项以及精度分析:
由插值型的求积公式我们得到求积公式误差余项的表达式:
其中表示求积公式的代数精度,为不依赖于的待定系数,这个结果表明当是次数小于等于的多项式时由于,故此时,即前面的求积公式精确成立。而当时,故可求得:
因为梯形公式的代数精度为1,可以得到的值为
于是得到梯形公式的余项为:
又因为复合梯形公式要满足
综上所述,就得到了复合梯形公式的余项表达式:
同理可得复合辛普森公式的余项表达式:
结果分析:从以上余项的表达式可以看出复合辛普森公式的代数精度为3,而复合梯形公式的代数精度为1,所以复合辛普森比复合梯形精确度更高。对于算法的精度,是通过对设计所得值与准确值之间的误差值来评判。将变步长的复合求积方法每次求得计算结果与准确值进行比较求出误差值,通过画出误差值的变化趋势图比较复合梯形公式与复合辛普森公式这两种算法的精度。经过实验中验证,也表明自己的初步推理是正确的,无论是复合梯形公式还是复合辛普森公式它们最终结果都会随着步长
值的减小而更加精确。复合梯形公式和复合辛普森公式计算出的结果进行比较,发现复合辛普森公式计算出的结果更加的精确。
5原始程序、数据
文件:f.m
function
y=f(x);
y=sqrt(x)*log(x);
1、复合梯形公式求解算法:
文件:trapezoid.m
clc
a=0;
%积分下限
b=1;
%积分上限
T=[];
%用来装不同n值所计算出的结果
R=[];
G=[];
m=120;
%等分数
true=-(4/9);
for
n=2:m;
h=(b-a)/n;
%步长
x=zeros(1,n+1);
%给节点定初值
y=zeros(1,n+1);
for
i=1:n+1
x(i)=a+(i-1)*h;
end
x(1,1)=0.000000001;
for
i=1:n+1
y(i)=x(i).^(1/2)*log(x(i))
%
g=(-(b-a)/12*h^
2)*(-log(x(i))/(4*x(i)*x(i)^(1/2)))
%准确的积分余项(计算误差)
end
%
G=[G,g];
t=0;
r=0;
for
i=1:n
format
long
t=t+h/2*(y(i)+y(i+1))
;
%利用复化梯形公式求值
err=t-floor(t);
digits(7);
%
此处为需要的小数位+1
t=floor(t)+vpa(err,6)
;
%
此处控制显示的小数点位数,更改显示的小数位数
r=t-true;
%计算的值与真实值之差(实际误差)
end
T=[T,t]
;
%把不同n值所计算出的结果装入
T中
R=[R,r];
end
x=linspace(0,1,m-1);
plot(x,R,'*')
%将计算误差与实际误差用图像画出来
2、复合辛普森积分求解算法:
simpon.m
clc
clear
a=0;
%积分下限
b=1;
%积分上限
T=[];
%用来装不同n值所计算出的结果
R=[];
true=-(4/9);
m=20;
%等分数
for
n=2:m
h=(b-a)/(2*n);
%步长
x=zeros(1,2*n+1);
%给节点定初值
y=zeros(1,2*n+1);
for
i=1:2*n+1
x(i)=a+(i-1)*h;
%给节点赋值
end
x(1,1)=0.000000001;
for
i=1:2*n+1
y(i)=x(i).^(1/2)*log(x(i));
%给相应节点处的函数值赋值
end
t=0;
r=0;
for
i=1:n
format
long
t=t+h/3*(y(2*i-1)+4*y(2*i)+y(2*i+1));
%利用复化simpson公式求值
err=t-floor(t)
digits(7);
%
此处为需要的小数位+1
t=floor(t)+vpa(err,6);
r=t-true;
end
T=[T,t]
%把不同n值所计算出的结果装入
T中
R=[T,r];
end
%
R=(-(b-a)/180*((b-a)/2).^4*24)
%积分余项(计算误差)
%
true=quad(@fx1,0,1);
%积分的真实值
x=linspace(0,1,2*m-1);
plot(x,R,'*')
3、龙贝格算法
rebeg.m
%%龙贝格
clear
clc
a=0;
b=1;
%确定积分上下限
eps=10^(-4);
err=1;
k=1;
a=0.0000001;
T(1,1)=(b-a)/2*(f(a)+f(b));
while(err>eps)
h=(b-a)/2^(k-1);
S=0;
for
x=a:h:b-h
S=S+f(x+h/2);
end
T(k+1,1)=1/2*T(k,1)+h/2*S;
k=k+1;
for
i=2:k
T(k,i)=(4^(k-1)*T(k,i-1)+T(k-1,i-1))/(4^(k-1)-1);
end
err=abs(T(i,i)+4/9);
end
fprintf('龙贝格求积算法积分值为%.10f\n',T(k,k));
disp(T)
T
龙贝格求积算法积分值为-0.44438207534、自适应辛普森算法:
%自适应辛普森算法
Self_Adaptive_integral.m
function
s=Self_Adaptive_integral(a,b,tol)
k=0;
w=0;
x=a;
y=b;
t=0;
h=(b-a)/2;
s=0;
i=0;
to=abs(simpson_integral(x,y,2)-simpson_integral(x,y,1));
while
to>=tol
i=i+1;
while
to>=tol
t=x;
if
k==0
x=t;
y=t+h;
to=(abs(simpson_integral(x,y,2)-simpson_integral(x,y,1)))*2^i;
k=1;
w=0;
end
if
w==0
x=t+h;
y=t+2*h;
to=(abs(simpson_integral(x,y,2)-simpson_integral(x,y,1)))*2^i;
k=0;
w=1;
end
end
s=s+simpson_integral(x,y,2);
if
k==0
x=t;
y=t+h;
h=h/2;
to=(abs(simpson_integral(x,y,2)-simpson_integral(x,y,1)))*2^i;
end
if
w==0
x=t+h;
y=t+2*h;
h=h/2;
to=(abs(simpson_integral(x,y,2)-simpson_integral(x,y,1)))*2^i;
end
if
to s=s+simpson_integral(x,y,2); to=tol-10; end end %自适应辛普森算法 simpson_integral.m function s=simpson_integral(a,b,m) h=(b-a)/(2*m); s1=0; s2=0; s=0; if m>1 for i=1:(m-1) x=a+2*i*h; s1=s1+f(x); end for i=1:m x=a+(2*i-1)*h; s2=s2+f(x); end s=h/3*(f(a)+f(b)+2*s1+4*s2); else s=s+h/3*(f(a)+f(b)++4*f((a+b)/2)); end 6结果分析与设计总结 结果分析: 初步分析:通过对步长h的值的改变,只要h值越小(等分数n的值越大),即等分的区间越小,结果应该更加精确,精确度越高。实验结果分析: 1、复合梯度算法: 通过算法的运行结果可得:当等分数n从2开始变化到50时实验计算结果以及与准确值之间的误差可以达到-0.4410968和0.003347665; 当等分数更改为80时实验计算结果以及与准确值之间的误差可以达到-0.4426581和0.001786344; 结果分析:复合梯形求积公式随着区间数不断增加,积分的误差不断减小。运算所得结果如下图所示: 2、复合辛普森算法: 当等分区间数目达到50时,实验计算出来的结果以及与准确值(-4/9)之间的误差值分别为:-0.443796和0.0006484617; 而当n为80时实验结果分别为-0.4441055和0.0003389765; 结果分析:复合辛普森求积公式也是随着等分数不断增加,积分的误差不断减小。算法计算结果如图所示: 将这两种插值形的积分算法所得结果与准确值比较,得出两种方法产生的误差趋势图如下: 经过算法设计验证,也表明自己的初步推理是正确的,无论是复合梯形公式还是复合辛普森公式计算的最终结果都会随着步长值的减小而更加精确,更加趋近于准确值。复合梯形公式和复合辛普森公式计算出的结果进行比较,发现复合辛普森公式计算出的结果更加的精确。 3、龙贝格算法: 通过龙贝格算法,当要求积分计算结果达到精度为时:实验所得结果为-0.***,而当要求的计算精度达到时,计算所得结果为-0.***。 计算结果为: k h …… 0 …… 0 0 …… …… 0 …… …… -0.444***4 0 …… …… -0.*** -0.*** 0 -0.*** -0.*** -0.*** 4、自适应辛普森算法: 输入s1=Self_Adaptive_integral(0.0000001,1,0.01)时 运算所得结果为-0.***; 输入s2=Self_Adaptive_integral(0.0000001,1,0.0001)时,运算结果也是-0.***; 设计总结: 虽然此次课程设计时间不是很长,但是还是让我学会了不少。不仅是运筹学知识的应用还是对于数值分析中多种数值计算方法的回顾,都让我对于专业知识得到进一步地加深理解。本课程设计也让我更加熟练地掌握了应用MATLAB编写相应的算法求解相应的数学问题,将理论知识与实际应用想结合,提高了自身的算法设计能力以及编程程序的技能。这一过程也得利于到老师、同学以及组员的帮助,才能如期完成自己的任务。这一次的课程设计更让学会对问题的分析以及思考,以及查找算法中的不足并作出改进。 参考文献 [1] 李庆扬,王能超,易大义.数值分析第四版[M].北京:清华大学出版社.2001.[2] 董霖.MATLAB使用详解第一版[M].科学出版社.2008 [3] 龚纯,王正林.MATLAB语言常用算法程序集[M].电子工业出版社.2008 [4] 李庆杨.科学计算方法基础[M].北京:清华大学出版社,2006 [5] 白峰杉.数值计算引论[M].北京:高等教育出版社,2004