贵州师范学院毕业论文(设计)
for j=0:i
x(j+1)=a+j*h; s=1; if j==0 s=s; else for k=1:j s=s*k; end end r=1; if i-j==0 r=r; else
for k=1:(i-j) r=r*k; end end
if mod((i-j),2)==1 q=-(i*s*r); else q=i*s*r; end
17
贵州师范学院毕业论文(设计)
y=1; for k=0:i if k~=j
y=y*(sym('t')-k); end end l=int(y,0,i); C(j+1)=l/q;
z=z+C(j+1)*f1(m,x(j+1)); end
g=(b-a)*z
1)当输入a?0,b?1,n?1,m?2时,即在MATLAB命令窗口输入 >> NCotes(0,1,1,2)即可得用梯形公式的积分值和相应科特斯系数如图1.1
2)当输入a?0,b?1,n?2,m?2时,即在MATLAB命令窗口输入 >> NCotes(0,1,2,2)即可得用辛浦生公式的积分值和相应科特斯系数 如图1.2
3)当输入a?0,b?1,n?4,m?2时,即在MATLAB命令窗口输入 >> NCotes(0,1,4,2)即可得用科特斯公式的积分值和相应科特斯系数 如图1.3
18
贵州师范学院毕业论文(设计)
图 1.1
图 1.2
19
贵州师范学院毕业论文(设计)
图1.3
2 复化梯形求积公式的MATLAB实现
通过f(x)的n?1个等步长节点逼近积分
其中,xk?a?kh,x0?a, 程序二:
function s=trapr1(f,a,b,n) % f是被积函数;
% a,b分别为积分的上下限; % n是子区间的个数; % s是梯形总面积;
20