自适应Simpson积分算法(MATLAB及C++实现代码)

自适应Simpson积分算法（MATLAB及C++实现代码）(计算数学课用)在CSDN论坛中找到了却要金币，无奈之下自己写了一份。对于类似问题，改一下积分函数和区间即可。针对问题：数学上已经证明了∫41+𝑥2𝑑𝑥=𝜋10成立，所以可以通过数值积分来求π的近似值。试利用自适应Simpson算法计算积分近似值。C++版：(直接复制粘贴在VC++6.0即可运行)/*用自适应Simpson积分方法计算积分值*/#includeiostream#includecmathintn=0;//设置全局变量n,用来记录最高迭代次数，避免递归一直进行下去。doublepi=3.141592653589793238462643;//设置近似精确值，用以比较doublee1=0.00001;//设置误差容限为10^-5doublef(double);//要积分的函数doubleSimpson(double,double,double,double);//迭代函数usingnamespacestd;//主函数intmain(){doublea=0,b=1,t,h,S;//积分区间h=(b-a)/2;S=h/3*(f(a)+f(b)+4*f((a+b)/2));//第一次Simpson公式积分值t=Simpson(a,b,e1,S);cout积分值为：tendl;cout最大迭代次数为：nendl;cout设置误差容限为e1\n误差为:pi-tendl;return0;}//子函数1（积分函数）doublef(doublex){return4/(1+x*x);}//子函数2（迭代函数）doubleSimpson(doubleA,doubleB,doublee,doubleS){doubleh,S1,S2;h=(B-A)/2;n++;//统计迭代次数if(n500){cout方法有误，跳出递归endl;return0;}S1=h/6*(f(A)+f(A+h)+4*f(A+h/2));//在[A,(A+B)/2]区间上计算Simpson积分值S2=h/6*(f(A+h)+f(B)+4*f(A+3/2*h));//在[(A+B)/2，B]区间上计算Simpson积分值if(fabs(S-S1-S2)15*e)returnS1+S2;//如果满足误差容限要求，就以S1+S2作为此时对应区间上的函数的近似值elsereturnSimpson(A,(A+B)/2,e/2,S1)+Simpson((A+B)/2,B,e/2,S2);//递归调用}MATLAB版：（两个函数文件加一个脚本文件）1.编写积分函数文件：functiony=f(x)y=4./(1+x.^2);end2.编写Simpson迭代函数文件functiony=Simpson(A,B,e,S)h=(B-A)/2;S1=h/6*(f(A)+f(A+h/2)+4*f(A+h/2));S2=h/6*(f(A+h)+4*af(A+3/2*h)+f(B));ifabs(S-S1-S2)10*ey=S1+S2;elsey=Simpson(A,(A+B)/2,e/2,S1)+Simpson((A+B)/2,B,e/2,S2);endend3.编写脚本调用文件ticclear;a=0;b=1;%积分区间e=0.0000001;%误差容限h=(b-a)/2;S=h/3*(f(a)+f(b)+4*f(1/2*(a+b)));%第一次Simpson积分值t=Simpson(a,b,e,S)%最终自适应方法积分值abs(pi-t)%实际误差e%设置的误差容限toc%返回所用时间亲测可用。这两个代码本质上是一样的。我先用C++语言写好，然后又换用成MATLAB语言。MATLAB好像可以把误差容限调到10^-7以下，而C++则只能到10^-5左右。原因不甚了解，猜测可能是由于C++计算时字节长度不够，导致精度不够，要递归调用很多次才能达到所需精度。