数值逼近数值积分MATLAB实现

问题描述

积分函数如下

积分区间为[0,1]，我们可以很容易地发现这个积分值是π（4*arctan1）。

数值积分

复化积分

首先将[0,1]区间8等分，分别用8 阶复化梯形公式、4 阶复化Simpson 公式、2 阶复化Cotes 公式来计算

x=[0 0.125 0.25 0.375 0.5 0.625 0.75 0.875 1];
%---------------------
fprintf('------------------------------\n');
T=0;
for i=2:8
    y=2*4/(1+x(i)^2);
    T=T+y;
end
T=T+4/(1+x(1)^2)+4/(1+x(9)^2);
T=T/16;
fprintf('利用8阶复化梯形公式\n');
vpa(T,11)
%---------------------
fprintf('------------------------------\n');
T=0;
for i=2:2:8
    T=T+4*4/(1+x(i)^2);
end
for i=3:2:7
    T=T+2*4/(1+x(i)^2);
end
T=T+4/(1+x(1)^2)+4/(1+x(9)^2);
T=T/24;
fprintf('利用4阶复化Simpson公式\n');
vpa(T,11)
%---------------------
fprintf('------------------------------\n');
T=0;
for i=2:2:8
    T=T+32*4/(1+x(i)^2);
end
T=T+12*(4/(1+x(3)^2)+4/(1+x(7)^2));
T=T+14*4/(1+x(5)^2);
T=T+7*(4/(1+x(1)^2)+4/(1+x(9)^2));
T=T/180;
fprintf('利用2阶复化Cotes公式\n');
vpa(T,11)

运行结果如下

------------------------------
利用8阶复化梯形公式
ans =
3.1389884945
------------------------------
利用4阶复化Simpson公式
ans =
3.1415925025
------------------------------
利用2阶复化Cotes公式
ans =
3.1415940941

Gauss-Legendre 求积

在使用这个方法时，要首先将积分区间化为[-1,1]再进行计算。分别计算n=1,3,5时候的积分值。

%函数化为2/(1+(x+1)^2/4),积分区间为[-1,1]
I1=2/(1+(-0.5773502691896250+1)^2/4)+2/(1+(0.5773502691896250+1)^2/4);
I2=0.3478548451374530*2/(1+(-0.8611363115940520+1)^2/4)+...
    0.6521451548625460*2/(1+(-0.3399810435848560+1)^2/4)+...
    0.6521451548625460*2/(1+(0.3399810435848560+1)^2/4)+...
    0.3478548451374530*2/(1+(0.8611363115940520+1)^2/4);
I3=0.1713244923791700*2/(1+(-0.9324695142031520+1)^2/4)+...
    0.3607615730481380*2/(1+(-0.6612093864662640+1)^2/4)+...
    0.4679139345726910*2/(1+(-0.2386191860831960+1)^2/4)+...
    0.4679139345726910*2/(1+(0.2386191860831960+1)^2/4)+...
    0.3607615730481380*2/(1+(0.6612093864662640+1)^2/4)+...
    0.1713244923791700*2/(1+(0.9324695142031520+1)^2/4);
fprintf('------------------------------\n');
fprintf('n=1时利用Gauss-Legendre求积公式\n');
vpa(I1,11)
fprintf('------------------------------\n');
fprintf('n=3时利用Gauss-Legendre求积公式\n');
vpa(I2,11)
fprintf('------------------------------\n');
fprintf('n=5时利用Gauss-Legendre求积公式\n');
vpa(I3,11)

运行结果如下

------------------------------
n=1时利用Gauss-Legendre求积公式
ans =
3.1475409836
------------------------------
n=3时利用Gauss-Legendre求积公式
ans =
3.1416119052
------------------------------
n=5时利用Gauss-Legendre求积公式
ans =
3.1415926112

积分中的Gauss节点可以在网上找到（这里的N应该是2,4,6）

总结

可以对比一下上述两种积分方法的计算结果的误差

其中效果最好的是n=5时的Gauss-Legendre求积公式。
不过这里的代码写的也是很傻……但是看上去是直观的……