数值分析与算法作业

第七章上机题作业。

1 解：由于步长h要改变，故将其设为变量，则编程实现romberg求积算法7.1的程序为：

romberg 使用romberg积分算法求积分值。

function [i,n,err]=romberg(fun,a,b,h,ep)

fun是被积函数，a,b是积分下、上限，h是步长，ep是误差精度。

i返回积分近似值，n返回计算次数即t数表的列数，err返回结果误差。

n=1;k=1;

t(1,1)=h*(feval(fun,a)+feval(fun,b))/2;

t(n+1,k)=t(n,k)/2+sum(feval(fun,a+h/2:h:b-h/2))*h/2;

根据递推化的复合梯形公式计算t（n+1,k）

t(n+1,k+1)=(4^k*t(n+1,k)-t(n,k))/4^k-1); t(n+1,k+1)的外推公式。

while(abs(t(n+1,n+1)-t(n,n))>ep)

n=n+1;

h=h/2;

t(n+1,1)=t(n,1)/2+sum(feval(fun,a+h/2:h:b-h/2))*h/2;

for k=1:n

t(n+1,k+1)=(4^k*t(n+1,k)-t(n,k))/4^k-1);

endend

i=t(n+1,k+1);

err=i-pi;

要计算的近似值，将matlab中的输出数据格式设为“long”，以便观察精确位数，误差精度取为1e-6，输入：

i,n,err]=romberg(inline('4./(1+x.^2)')0,1,1,1e-6)

结果为：i =3.141592653638244

n =5err =4.845102097306153e-011

1）改变不同的步长，将步长分别设为1/3,1/4,1/5,1/6,得到结果误差分别为：

err=-6.166799795259692e-007，-6.916489905250955e-007，-7.

363078342592644e-007，-7.659459657638479e-007

由上述结果可知，随着步长的不断减小，结果误差的绝对值不断增大，说明结果误差越来越大，故不存在最小的h，使结果准确度不再改善。

2）自适应求积的matlab程序quadtx书上有源程序，程序如下：

function [q,fcount]=quadtx(f,a,b,tol,varargin)

f表示单变量被积函数f(x)，积分区间为[a,b]，tol表示误差阈值，默认值为1e-6，varargin表示传给函数f的额外参数，返回值q为计算出的积分值，fcount表示积分计算过程中计算f(x)的次数，它不是必须的参数。

if ischar(f)&exist(f)~=2

f=inline(f);

elseif isa(f,'sym')

f=inline(char(f));

endif nargin<4|isempty(tol) %设置默认的tol值。

tol=end

c=(a+b)/2; %计算simpson公式所需的三个函数值。

fa=feval(f,a,varargin);

fc=feval(f,c,varargin);

fb=feval(f,b,varargin);

递归调用，k为计算函数值的次数。

q,k]=quadtxstep(f,a,b,tol,fa,fc,fb,varargin);

fcount=k+3;

function [q,fcount]=quadtxstep(f,a,b,tol,fa,fc,fb,varargin)

quadtx使用的递归调用子程序。

h=b-a;

c=(a+b)/2;

fd=feval(f,(a+c)/2,varargin);

fe=feval(f,(c+b)/2,varargin);

q1=h/6*(fa+4*fc+fb);

q2=h/12*(fa+4*fd+2*fc+4*fe+fb);

if abs(q2-q1)<=tol

q=q2+(q2-q1)/15;

fcount=2;

else[qa,ka]=quadtxstep(f,a,c,tol,fa,fd,fc,varargin);

[qb,kb]=quadtxstep(f,c,b,tol,fc,fe,fb,varargin);

q=qa+qb;

fcount=ka+kb+2;

end输入：

q,fcount]=quadtx(inline('4./(1+x.^2)')0,1,1e-6)

结果为：q =3.141592653708037

fcount =33

设置不同的误差阈值，得到自适应积分的计算量随阈值变化的情况表如下：

输入：x=0:0.01:1;

y=4./(1+x.^2);

plot(x,y,'.

得到函数曲线图，积分节点的分布情况示于曲线图中：

3）比较（1）、（2）两种方法的计算量可知，当误差阈值均为1e-6，设置的步长h均为b-a时，方法（1）的计算量为5次，方法（2）的计算量为33次，说明方法（2）的计算量比方法（1）更大，方法（2）的结果误差是随着阈值的减小而逐渐减小的。

数值分析与算法作业

数值分析与算法作业

算法与分析平时作业

数值分析作业

其他用户还读了