ppk计算公式复合梯形公式与复合辛普森公式求积分

复合梯形公式与复合辛普森公式求积分
2010-12-26 09:37:23| 分类： 数值分析 | 标签： |字号大中小 订阅
一 实验目的
1. 掌握复合梯形公式与复合辛普森公式的基本思想。
2. 编程实现用复合梯形公式与复合辛普森公式求积分。
3. 熟悉matlab软件的使用。
二 实验内容
1、用复合梯形公式计算积分 I=4(1+x2)dx ,求它0到1的积分。精确度为10-5.（0.00001）
，精确到
●1 计算公式
h=(b-a)n

h=h2[(f(x0)+f(x1))+(f(x 1)+f(x2))+(f(x2)+f(x3)+...+(f(xn-1)+f(xn)]

。

l1 算法分析
En=h212[f'(b)-f'(a)]
将区间[a,b]等分成n个小区间，在小区间上 分别应用低次积分公式来构造公式，通过for循环来实现，
分的越细，越接近实际结果，精确度越高。
l2 源程序
function f1=fun4(x) %原函数
f1=4(1+x^2);
function ff=fun2(x) %函数对x求导
ff=-8*x((1+x^2)^2);
function f=tixing(a,b) %a，b是区间
a=0;b=1;
disp('******复合梯形公式******')
h=0.008; %h表示区间被等分成若干份后，每两个数的间距
m=(a:h:b); %形成一维矩阵，每两个数间的间隔是h
n=length(m); %求上矩阵的长度，即元素个数
for i=1:n-1
D(i)=fun4(m(i))+fun4(m(i+1));
end
R=h2*sum(D); %积分结果
E=-(h^2)*(fun2(b)-fun2(a))12; %余项，即精度
t=pi-R;
[R;E;t]
实验结果讨论和分析
通过 对h的值的改变，发现h值越小，即等分的区间越小，结果越精确，精确度越高。通过手算得到
积分结果 为π，实验结果为3.313，结果正确，可见复合梯形公式的精确度较高，运算次数
为125.
2、用复合辛普森公式计算积分I=4(1+x2)dx ,求它0到1的积分。精确度为10-5.（0.00001）
l5 计算公式
h=(b-a)2n=(xi+1-xi)2 ；（i=0，1，…n-1）
S=h3[f(xi)+4f(xi+12)+f(xi+1)] (i=0，1，…n-1)
l6 算法分析
复合辛普森公式来求积分是将区间等分为2n份，在每两个相邻 的数间再取中间值，利用for循环实现
辛普森公式。该公式等分的份数更多，是的精度也更高。
l7 源程序
function f1=fun4(x)
f1=4(1+x^2); %公式f（x）
function f=xinpusen(a,b) %a，b分别为区间的端点值
a=0;b=1;
disp('******复合辛普森形公式******')
h1=0.25; %h表示区间被等分成若干份后，每两个相邻数的间距
m=(a:h1:b);
h=h12;
n=length(m);
for i=1:n-1
Z(i)=(m(i)+m(i+1))2;
D(i)=fun4(m(i))+fun4(m(i+1))+4*fun4(Z(i));
end
R=h3*sum(D);
t=pi-R; %精度
[R;t]
l9 实验结果讨论和分析
从计算结果可以看到，复合 辛普森你公式结果更接近精确解，精确度更高，而且运算次数只有40次，
大大减少了运算次数，比复合 梯形公式收敛性高。
三 本次实验总结
在本次实验过程中，我掌握了复合梯形公式与复合 辛普森公式的基本算法与思想，通过编程来实现用复合
梯形公式与复合辛普森公式求积分。而且通过上机 实验，可以看到复合辛普森公式得到的结果更加精确，
运算次数比较少。同时对matlab的使用也更 加熟练，对其中常用语句运用更自如。









MATLAB 用辛普森系列公式求积分
function [I,step] = IntSimpson(f,a,b,type,eps)
%type = 1 辛普森公式
%type = 2 辛普森38公式
%type = 3 复合辛普森公式
if(type==3 && nargin==4)
eps=1.0e-4; %缺省精度为0.0001
end
I=0;
switch type
case 1,
I=((b-a)6)*(subs(sym(f),findsym(sym(f)),a)+...
4*subs(sym(f),findsym(sym(f)),(a+b)2)+...
subs(sym(f),findsym(sym(f)),b));
step=1;

case 2,
I=((b-a)8)*(subs(sym(f),findsym(sym(f)),a)+...
3*subs(sym(f),findsym(sym(f)),(2*a+b)3)+ ...
3*subs(sym(f),findsym(sym(f)),(a+2*b)3)+subs(sym(f ),findsy
m(sym(f)),b));
step=1;

case 3,
n=2;
h=(b-a)2;
I1=0;
I2=(subs(s ym(f),findsym(sym(f)),a)+subs(sym(f),findsym(sym(f )
),b))h;
while abs(I2-I1)>eps
n=n+1;
h=(b-a)n;
I1=I2;
I2=0;
for i=0:n-1
x=a+h*i;
x1=x+h;
I2=I2+(h6)*(subs(sym(f),findsym(sym(f)),x)+...
4*subs(sym(f),findsym(sym(f)),(x+x1)2)+...
subs(sym(f),findsym(sym(f)),x1));
end
end
I=I2;
step=n;
end