数值积分 (论文)
for(j=1;j<=sum_num;j++)
{
add_T=add_T+2*f_x(s_point+(j-1)*d_point)-4*f_x(s_point-sd_point+(j-1)*d_point)-4*f_x(s_point+sd_point+(j-1)*d_point);
}
add_T=add_T*xishu;
else
{
T0=T1;
T1=0;
add_T=0;
err_T=0;
}
}
在这个函数中我们将复化simpose公式和积分过程都用计算机语言表示出来。首先我们给出复化simpose公式,进行迭代,直到精确度达到设定பைடு நூலகம்求,算出最后结果。
3.3 测试结果
用复化simpose迭代取有效数字四位求得的结果如下:
用复化simpose迭代取有效数字七位求得的结果如下:
double add_T=0;
T0=(a-b)*(f_x(a)+f_x(b))/2;//n=1时的cotes公式即梯形公式。
for(i=1;i<=100;i++)//the first base number
{
sum_num=pow(2,i-1);
xishu=double(a-b)/sum_num;
在微积分理论中,我们知道了牛顿-莱布尼茨(Newton-Leibniz)公式
其中 是被积函数 的某个原函数。但是随着学习的深入,我们发现一个问题:对很多实际问题,上述公式却无能为力。这主要是因为:它们或是被积函数没有解析形式的原函数,或是只知道被积函数在一些点上的值,而不知道函数的形式,对此,牛顿—莱布尼茨(Newton-Leibniz)公式就无能为力了。此外,即使被积函数存在原函数,但因找原函数很复杂,人们也不愿花费太多的时间在求原函数上,这些都促使人们寻找定积分近似计算方法的研究,特别是有了计算机后,人们希望这种定积分近似计算方法能在计算机上实现,并保证计算结果的精度,具有这种特性的定积分近似计算方法称为数值积分。由定积分知识,定积分只与被积函数和积分区间有关,而在对被积函数做插值逼近时,多项式的次数越高,对被积函数的光滑程度要求也越高,且会出现Runge现象。如 时,Newton-Cotes公式就是不稳定的。因而,人们把目标转向积分区间,类似分段插值,把积分区间分割成若干小区间,在每个小区间上使用次数较低的Newton-Cotes公式,然后把每个小区间上的结果加起来作为函数在整个区间上积分的近似,这就是复化的基本思想。本文主要研究的公式有:复化梯形公式﹑复化Simpson公式﹑复化Cotes公式﹑Romberg积分法。
把[a,b] 2k等分,分点xi=a+(b-a)/ 2k ·i(i =0,1,2 · · · 2k)每个小区间长度为(b-a)/2k,由归纳法可得面所说的的第一个公式.
(二)计算公式如下:
整个程序就是循着这四个公式进行计算的。Sn,Cn, Rn分别代表特例梯形积分,抛物线积分,龙贝格积分.当然,编程的时候统一处理即可。
由以上结果可以看出两次不同精度要求的计算可以看出不同精度计算计算次数相差较多。精度为四和七间计算次数相差了三次。
第四章 复化cotes公式
4.1 复化cotes公式的算法描述
复化求积公式
当L=4时可得复化梯形公式:
=
复化cotes公式=
4.2 复化cotes公式在C语言中的实现
复化cotes公式运用的程序如下:
sum_num=pow(2,i-1); //the same as T
xishu=double(a-b)/sum_num/6;
s_point=double(b)+double(a-b)/pow(2,i);
d_point=double(a-b)/pow(2,i-1);
sd_point=double(a-b)/sum_num/4;
复化求积公式
当L=1时可得复化梯形公式:
=
复化梯形公式=
2.2 复化梯形公式在C语言中的实现
复化梯形公式运用的程序如下:
T0=(a-b)*(f_x(a)+f_x(b))/2;//n=1时的cotes公式即梯形公式
for(i=1;i<=100;i++)
{
//计算sum_num、xishu、s_point(start point)、d_point
T0=(a-b)*(7*f_x(1)+32*f_x(2)+12*f_x(3)+32*f_x(4)+7*f_x(5))/90;//四阶(n=4)cotes公式。
for(i=1;i<=100;i++)
{
sum_num=pow(2,i-1); //the same as T
xishu=double(a-b)/sum_num/90;
s_point=double(b)+double(a-b)/pow(2,i);
d_point=double(a-b)/pow(2,i-1);
for(j=1;j<=sum_num;j++)
{
add_T=add_T+f_x(s_point+(j-1)*d_point);
}
add_T=add_T*xishu;
T0=(a-b)*(f_x(a)+4*f_x((a+b)/2)+f_x(b))/6;//n=2的cotes公式即simpson公式
for(i=1;i<=100;i++)
{
//计算sum_num、xishu、s_point(start point)、d_point
//long powl (long double x, long double y)
T1=(T0+add_T)/2;
add_T=0;
//计算S1
S1=4*T1/3-T0/3;
if(i>=2)
C1=16*S1/15-S0/15;
if(i>=3)
s_point=double(b)+double(a-b)/pow(2,i);
d_point=double(a-b)/pow(2,i-1);
sd_point=double(a-b)/sum_num/8;
for(j=1;j<=sum_num;j++)
{
add_T=add_T-2*f_x(s_point+(j-1)*d_point)-32*f_x(s_point-sd_point+(j-1)*d_point)+20*f_x(s_point-2*sd_point+(j-1)*d_point)-32*f_x(s_point-3*sd_point+(j-1)*d_point)-32*f_x(s_point+sd_point+(j-1)*d_point)+20*f_x(s_point+2*sd_point+(j-1)*d_point)-32*f_x(s_point+3*sd_point+(j-1)*d_point);
1.2 问题重述
本文主要介绍微积分方程的复化解法。通过运用复化梯形公式、复化Simpose公式、复化cotes公式和Romberg积分法这四种积分法方法,解出微分方程的近似解。并进行误差分析和结果比较。
当积分区间[a,b]的长度较大,而节点个数n + 1固定时,直接使用Newton-Cotes公式的余项将会较大,而如果增加节点个数,即n + 1增加时,公式的舍入误差又很难得到控制,为提高公式的精度,又使算法简单易行,往往使用复化方法。即将积分区间[a,b]分成若干个子区间,然后在每个小区间上使用低阶Newton-Cotes公式,最后将每个小区间上的积分的近似值相加。
sum_num=pow(2,i-1);
xishu=double(a-b)/sum_num;
s_point=double(b)+double(a-b)/pow(2,i);
d_point=double(a-b)/pow(2,i-1);
for(j=1;j<=sum_num;j++)
{
add_T=add_T+f_x(s_point+(j-1)*d_point);
}
add_T=add_T*xishu;
T1=(T0-add_T)/2;
err_T=(T1-T0)/63;
//output
printf("%d %d %10.8f %10.8f",i,pow(2,i),T1,err_T);
printf("\n");
if(err_T<=0)
err_T=(-1)*err_T;
对区间[a,b],令h=b-a构造梯形值序列{T2K}。
T1=h[f(a)+f(b)]/2
把区间二等分,每个小区间长度为h/2=(b-a)/2,于是
T2 =T1/2+[h/22]f(a+h/2)
把区间四(22)等分,每个小区间长度为h/22 =(b-a)/4,于是
T4 =T2/2+[h/2][f(a+h/4)+f(a+3h/4).....................
由以上结果可以看出取两个不同的精度相对误差比较小,但计算次数大大的增加,复化梯形公式计算次数多。
第三章 复化simpson公式
3.1 复化simpson公式的算法描述
复化求积公式
当L=2时可得复化梯形公式:
=
复化simpson公式=
3.2 复化simpson公式在C语言中的实现