变步长复化辛普森公式计算积分

变步长复化辛普森公式计算积分 matlab编程

2. 编写用变步长复化辛普森公式计算积分 b

af(x)dx 的程序。

1用上面编写的程序计算下列积分并分析计算结果 （1

）

0cosxdx （2

）0xcosxdx （3） 220xdx

程序：

function S=bianfuhuasimpson(fx,a,b,eps,M)

% 变步长复合simpson求积公式

% 调用方式： S=fuhuasimpson(@fx,a,b,epsilon)

% fx -- 求积函数（函数文件）

% a, b -- 求积区间

% eps -- 计算精度

% M--最大允许输出划分数

n=1;

h=(b-a)/n;

T1=h*(feval(fx,a)-feval(fx,b))/2;

Hn=h*feval(fx,(a+b)/2);

S1=(T1+2*Hn)/3;

n=2*n;

% 最好与倒数第三行保持一致（变步长）

while n<=M

T2=(T1+Hn)/2;

Hn=0;

h=(b-a)/n;

for j=1:n

x(j)=a+(j-1/2)*h;

y(j)=feval(fx,x(j));

Hn=Hn+y(j);

end

Hn=h*Hn;

S2=(T2+2*Hn)/3;

fprintf(' n=%2d S2=%-12.9f S2-S1=%-12.9f\n',n,S2,abs(S2-S1)); if abs(S2-S1)<eps

break;

else

T1=T2;

S1=S2;

n=2*n;

end

end

变步长复化辛普森公式计算积分 matlab编程

S=S2;

% 达到下列条件之一，则运算终止:

% (1).abs(S2-S1)<eps

% (2).下一次的n>M

% 输入1：S=bianfuhuasimpson(inline('sqrt(x)*cos(x)'),0,pi,10e-6,2000)

% 输入2：S=bianfuhuasimpson(inline('2*x^2*cos(x^2)'),0,sqrt(pi),10e-6,2000) % 输入3：S=bianfuhuasimpson(inline('sqrt(x)'),0,1,10e-6,2000)

输出结果：

（1）

S=bianfuhuasimpson(inline('sqrt(x)*cos(x)'),0,pi,10e-6,2000)

n= 2 S2=-0.016369112 S2-S1=0.944423778

n= 4 S2=-0.450266122 S2-S1=0.433897010

n= 8 S2=-0.669839370 S2-S1=0.219573248

n=16 S2=-0.781318443 S2-S1=0.111479074

n=32 S2=-0.837710689 S2-S1=0.056392245

n=64 S2=-0.866141900 S2-S1=0.028431211

n=128 S2=-0.880440980 S2-S1=0.014299080

n=256 S2=-0.887620063 S2-S1=0.007179083

n=512 S2=-0.891220052 S2-S1=0.003599989

n=1024 S2=-0.893023740 S2-S1=0.001803689

S =

-0.8930

（2）

S=bianfuhuasimpson(inline('2*x^2*cos(x^2)'),0,sqrt(pi),10e-6,2000)

n= 2 S2=1.076354541 S2-S1=2.092222287

n= 4 S2=0.039359358 S2-S1=1.036995183

n= 8 S2=-0.430456535 S2-S1=0.469815894

n=16 S2=-0.662796649 S2-S1=0.232340113

n=32 S2=-0.778823323 S2-S1=0.116026674

n=64 S2=-0.836827971 S2-S1=0.058004648

n=128 S2=-0.865829756 S2-S1=0.029001785

n=256 S2=-0.880330615 S2-S1=0.014500859

n=512 S2=-0.887581042 S2-S1=0.007250427

n=1024 S2=-0.891206256 S2-S1=0.003625214

S =

变步长复化辛普森公式计算积分 matlab编程

-0.8912

（3）

S=bianfuhuasimpson(inline('sqrt(x)'),0,1,10e-6,2000) n= 2 S2=0.489859598 S2-S1=0.185121744 n= 4 S2=0.579745947 S2-S1=0.089886349 n= 8 S2=0.623731522 S2-S1=0.043985575 n=16 S2=0.645384849 S2-S1=0.021653327 n=32 S2=0.656091436 S2-S1=0.010706587 n=64 S2=0.661402273 S2-S1=0.005310836 n=128 S2=0.664042680 S2-S1=0.002640407 n=256 S2=0.665357576 S2-S1=0.001314896 n=512 S2=0.666013147 S2-S1=0.000655572 n=1024 S2=0.666340270 S2-S1=0.000327122

S =

0.6663