复化梯形法-复化矩形法-变步长梯形-变步长辛普森

陕西科技大学机械教改班用C++的积分其实积分的思想就是,微分—>求和—>取极限，如果是用纯手工法那就是先对一个函数微分，再求出它的面积，在取极限，因为我们的计算速度和计算量有限，现在有了计算机这个速度很快的机器，我们可以把微分后的每个小的面积加起来，为了满足精度，我们可以加大分区，即使实现不了微分出无限小的极限情况，我们也至少可以用有限次去接近他，下面我分析了四种不同的积分方法，和一个综合通用程序。一．积分的根本思想1、思路：微分—>求和—>取极限。2、Newton—Leibniz公式其中，被积函数的原函数。3、用计算机积分的思路在积分区间内“微分—>求和—>控制精度〞。因为计算机求和不可以取极限，也就是不可以无限次的加下去，所以要控制精度。二．现有的理论1、一阶求积公式---梯形公式他只能精确计算被积函数为0、1次多项式时的积分。2、二阶求积分公式——牛顿、科特斯公式他只能精确计算被积函数为0、1、2、3次多项式时的积分。四种实现方法复化矩形法将积分区间[a,b]等分成n个子区间:那么h=(b-a)/n,区间端点值=a+kh............................源程序：#include<iostream.h>#include<math.h>doublef(doublex)//计算被积函数{ doubley; y=log(1+x)/(1+x*x);//被积函数 returny;}doubleTn(doublea,doubleb,intn)//求Tn{ doublet=0.0; doublexk;//区间端点值 doublet1,t2;//用来判断精度do { doubleh=(b-a)/n; for(intk=1;k<=n-1;k++)//每一小段的矩形叠加 { t1=t; xk=a+k*h; t+=h*f(xk); t2=t; } n++;//如果精度不够就对区间再次细分，直到到达精度要求 } while(fabs(t1-t2)<=1e-7);//判断计算精度 returnt;}voidmain(){ doublea=0.0;//积分下线 doubleb=2.0;//积分上限 intn=1024;//把区间分为1024段 cout<<Tn(a,b,n)<<endl;//输出积分结果}执行结果：复化梯形法方法和复化矩形法类似，只是把原来的矩形小面积变成了梯形小面积，但是精确度明显提高了，也就是说到达同样的精度需要的时间少了。变形一下：源程序：#include<iostream.h>#include<math.h>doublef(doublex)//计算被积函数{ doubley; y=log(1+x)/(1+x*x);//被积函数 returny;}doubleTn(doublea,doubleb,intn)//求Tn{ doublet=0.0; doublexk;//区间端点值 doublet1,t2,h=(b-a)/n;//用来判断精度do { h=(b-a)/n; for(intk=1;k<=n-1;k++)//余项叠加，相当于每一个小梯形相加 { t1=t; xk=a+k*h; t+=f(xk); t2=t; } n++;//如果精度不够就对区间再次细分，直到到达精度要求 } while(fabs(t1-t2)<=1e-7);//判断计算精度t=h*(f(a)+f(b))/2+t*h; //加上初项就是积分结果了 returnt;}voidmain(){ doublea=0.0;//积分下线 doubleb=2.0;//积分上线 intn=1024;//把区间分为1024段 cout<<Tn(a,b,n)<<endl;//输出积分结果}执行结果：变步长梯形法上面我们在应用复化积分法的时候会对区间进行分割，但是在要求精度是我们也不知道事先应该分割成多少段，分的少了精度达不到，分的多了计算量太大，为了解决这个问题我在上面加了一个循环，来不断地对区间进行再次划分，直到到达精度要求，变步长积分法与我用的方法类似，只不过变步长积分法的区间划分时成倍增加的。实现方法；由复化梯形法知道；步长h=(b-a)/n现在将步长分半，即在相邻两节点的中点处增加一个求积节点，那么变形一下：源程序：#include<iostream.h>#include<math.h>doublef(doublex)//计算被积函数的值{ doubley; y=log(1+x)/(1+x*x); returny;}doublet2ntn(doublea,doubleb){ intn=1; doublehn=b-a;//原步长 doubletn=0.0; doublet2n=(f(a)+f(b))*hn/2.0; while(fabs(t2n-tn)>1e-7)//判断精度 { tn=t2n; t2n=0.0; for(intk=0;k<=n-1;k++)//循环叠加 t2n+=f(a+hn/2.0+k*hn); t2n=tn/2.0+t2n*hn/2.0; n=n*2; hn=hn/2.0;//步长分半 } returnt2n;}voidmain(){ doublea=0.0; doubleb=2.0; cout<<t2ntn(a,b)<<endl;}执行结果：变步长辛普森法之前的积分斜边都是直线，如果用抛物线接近就会更准确复化辛普森求积公式然后就只要每次让他的积分区间加倍就行直到到达要求的精度#include<stdio.h>#include<iostream.h>#include<math.h>doublea=0.0,b=2.0,e=1e-7;//积分上下线，和精度要求intn=1024;doublef(doublex){ Doubley; y=log(1+x)/(1+x*x);//被积函数 returny;}floatsimpson(){ inti; doubles,sn,s2n,p,x,h; h=(b-a)/n;//步长 s=f(a)+f(b); p=0; x=a+h;//积分端点 for(i=1;i<=n-1;i++) { if((i%2)!=0) { p=p+4*f(x);//在区间中间时乘4 x=x+h; } else { s=s+2*f(x);//积分端点时乘2 x=x+h; } } s=s+p;//第一次求和 s2n=s*h/3;//积分值 do { sn=s2n; x=a+h/2;//变步长 s=s-p/2; p=0; for(i=1;i<=n;i++)//变步长只需要加就行了 { p=p+4*f(x); x=x+h; } s=s+p; n=n*2; h=h/2; s2n=s*h/3; }while(fabs(s2n-sn)>e);//控制精度 returns2n;}voidmain(){ cout<<simpson()<<endl;}执行结果：用C++写的综合程序#include<stdio.h>#include<iostream.h>#include<math.h>classBjhs//抽象类 { public: virtualdoublef(doublex)=0;//虚计算被积函数 virtualvoidprint()=0;//输出函数 virtualdoubleTn()=0;//虚函数 }; classFhjx:publicBjhs//复化矩形法类 { public: Fhjx(){a=0;b=2;e=1e-7;n=1024;}//构造函数付初值 doublef(doublex);//计算被积函数 doubleTn()//求Tn { doublet=0.0; doublexk;//区间端点值 doublet1,t2;//用来判断精度 do { doubleh=(b-a)/n; for(intk=1;k<=n-1;k++)//每一小段的矩形叠加 { t1=t; xk=a+k*h; t+=h*f(xk); t2=t; } n++;//如果精度不够就对区间再次细分，直到到达精度要求 } while(fabs(t1-t2)<=e);//判断计算精度 returnt; } voidprint()//输出 { cout<<"用复化矩形法计算的结果"<<Tn()<<endl; }private: doublea,b,e; intn;};doubleFhjx::f(doublex)//外联函数{ doubley; y=log(1+x)/(1+x*x);//被积函数 returny;}classFhtx:publicBjhs//复化梯形法{ public: Fhtx(){a=0;b=2;e=1e-7;n=1024;} doublef(doublex);//计算被积函数 doubleTn()//求Tn { doublet=0.0; doublexk;//区间端点值 doublet1,t2,h=(b-a)/n;//用来判断精度 do { h=(b-a)/n; for(intk=1;k<=n-1;k++)//余项叠加，相当于每一个小梯形相加 { t1=t; xk=a+k*h; t+=f(xk); t2=t; } n++;//如果精度不够就对区间再次细分，直到到达精度要求 } while(fabs(t1-t2)<=1e-7);//判断计算精度 t=h*(f(a)+f(b))/2+t*h; //加上初项就是积分结果了 returnt; } voidprint() { cout<<"用复化梯形法计算的结果"<<Tn()<<endl; }private: doublea,b,e; intn;};doubleFhtx::f(doublex){ doubley; y=log(1+x)/(1+x*x);//被积函数 returny;}classBbctx:publicBjhs//变步长梯形法{public:Bbctx(){a=0;b=2;e=1e-7;n=1;tn=0;} doublef(doublex);//计算被积函数的值doubleTn() { doublehn=b-a;//原步长 doublet2n=(f(a)+f(b))*hn/2.0; while(fabs(t2n-tn)>e)//判断精度 { tn=t2n; t2n=0.0; for(intk=0;k<=n-1;k++)//循环叠加 t2n+=f(a+hn/2.0+k*hn); t2n=tn/2.0+t2n*hn/2.0; n=n*2; hn=hn/2.0;//步长分半 } returnt2n; } voidprint() { cout<<"用变步长梯形法计算的结果"<<Tn()<<endl; }private: doublea,b,e;doubletn; intn;};doubleBbctx::f(doublex){ doubley; y=log(1+x)/(1+x*x);//被积函数 returny;}classBbcxps:publicBjhs//变步长辛普森法{public: Bbcxps(){n=1024;a=0;b=2;e=1e-7;}//积分上下线，和精度要求doublef(doublex);doubleTn(){ inti; doubles,sn,s2n,p,x,h; h=(b-a)/n;//步长 s=f(a)+f(b); p=0; x=a+h;//积分端点 for(i=1;i<=n-1;i++) { if((i%2)!=0)//判奇偶，也就是看哪点乘几 { p=p+4*f(x);//在区间中间时乘4 x=x+h; } else { s=s+2*f(x);//积分端点时乘2 x=x+h; } } s=s+p;//第一次求和 s2n=s*h/3;//积分值 do { sn=s2n; x=a+h/2;//变步长 s=s-p/2; p=0; for(i=1;i<=n;i++)//变步长只需要加就行了 { p=p+4*f(x); x=x+h; } s=s+p; n=n*2; h=h/2; s2n=s*h/3; }while(fabs(s2n-sn)<=