最新复化梯形法复化矩形法变步长梯形变步长辛普森.docx
《最新复化梯形法复化矩形法变步长梯形变步长辛普森.docx》由会员分享,可在线阅读,更多相关《最新复化梯形法复化矩形法变步长梯形变步长辛普森.docx(16页珍藏版)》请在冰豆网上搜索。
![最新复化梯形法复化矩形法变步长梯形变步长辛普森.docx](https://file1.bdocx.com/fileroot1/2023-2/1/8a460af2-2299-48b6-929c-d67456fa5036/8a460af2-2299-48b6-929c-d67456fa50361.gif)
最新复化梯形法复化矩形法变步长梯形变步长辛普森
复化梯形法复化矩形法变步长梯形变步长辛普森
陕西科技大学
机械教改班
用C++的积分
其实积分的思想就是,微分—>求和—>取极限,如果是用纯手工法那就是先对一个函数微分,再求出它的面积,在取极限,因为我们的计算速度和计算量有限,现在有了计算机这个速度很快的机器,我们可以把微分后的每个小的面积加起来,为了满足精度,我们可以加大分区,即使实现不了微分出无限小的极限情况,我们也至少可以用有限次去接近他,下面我分析了四种不同的积分方法,和一个综合通用程序。
一.积分的基本思想
1、思路:
微分—>求和—>取极限。
2、Newton—Leibniz公式
«SkipRecordIf...»
其中,«SkipRecordIf...»被积函数«SkipRecordIf...»的原函数。
3、用计算机积分的思路
在积分区间内“微分—>求和—>控制精度”。
因为计算机求和不可以取极限,也就是不可以无限次的加下去,所以要控制精度。
二.现有的理论
1、一阶求积公式---梯形公式
«SkipRecordIf...»
他只能精确计算被积函数为0、1次多项式时的积分。
2、二阶求积分公式——牛顿、科特斯公式
«SkipRecordIf...»
他只能精确计算被积函数为0、1、2、3次多项式时的积分。
3.四种实现方法
1.复化矩形法
将积分区间[a,b]等分成n个子区间:
«SkipRecordIf...»
则h=(b-a)/n,区间端点值«SkipRecordIf...»=a+kh
«SkipRecordIf...»
«SkipRecordIf...»
............................
«SkipRecordIf...»
«SkipRecordIf...»
源程序:
#include
#include
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<}
执行结果:
2.复化梯形法
方法和复化矩形法类似,只是把原来的矩形小面积变成了梯形小面积,但是精确度明显提高了,也就是说达到同样的精度需要的时间少了。
«SkipRecordIf...»
变形一下:
«SkipRecordIf...»
源程序:
#include
#include
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<}
执行结果:
3.变步长梯形法
上面我们在应用复化积分法的时候会对区间进行分割,但是在要求精度是我们也不知道事先应该分割成多少段,分的少了精度达不到,分的多了计算量太大,为了解决这个问题我在上面加了一个循环,来不断地对区间进行再次划分,直到达到精度要求,变步长积分法与我用的方法类似,只不过变步长积分法的区间划分时成倍增加的。
实现方法;
由复化梯形法知道;
«SkipRecordIf...»
步长h=(b-a)/n
现在将步长分半,即在相邻两节点«SkipRecordIf...»的中点«SkipRecordIf...»处增加一个求积节点,那么
«SkipRecordIf...»
变形一下:
«SkipRecordIf...»
源程序:
#include
#include
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<}
执行结果:
4.变步长辛普森法
之前的积分斜边都是直线,如果用抛物线接近就会更准确
复化辛普森求积公式
«SkipRecordIf...»
然后就只要每次让他的积分区间加倍就行直到达到要求的精度
#include
#include
#include
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<}
执行结果:
4.用C++写的综合程序
#include
#include
#include
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<<"用复化矩形法计算的结果"<}
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<<"用复化梯形法计算的结果"<}
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<<"用变步长梯形法计算的结果"<}
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)<=e);//控制精度
returns2n;
}
voidprint()
{
cout<<"用变步长辛普森法计算的结果"<}
private:
doublea,b,e,p;
intn;
};
doubleBbcxps:
:
f(doublex)
{
doubley;
y=log(1+x)/(1+x*x);//被积函数
returny;
}
voiddisplay(Bjhs&q)//输出函数
{
q.Tn();
q.print();
};
voidmain()//选择你想要的方法,用switch语句
{
chara;
Fhjxb;
Fhtxc;
Bbctxd;
Bbcxpse;
cout<<"求函数y=log(1+x)/(1+x*x)的积分值"<cout<<"请输入你要选用的方法的编号"<cout<<"A.复化矩形法B.复化梯形法"<cout<<"C.变步长梯形法D.变步长辛普森法"<cin>>a;
switch(a){
case'A':
display(b);break;
case'B':
display(c);break;
case'C':
display(d);break;
case'D':
display(e);break;
default:
cout<<"输入代码错误"<}
执行结果: