龙贝格算法求数值积分程序说明及流程图

《龙贝格算法求数值积分程序说明及流程图》由会员上传分享，免费在线阅读，更多相关内容在行业资料-天天文库。

1、陆韶琦3110000441bsinx程序说明：本程序计算dx，程序初始要求输出需要得到的精度，最后输出得到ax数值积分值。输入精度eps，积分上下限b,a流程图：定义函数f(x)=1x=0;sin(x)/x其他。计算Ti(0

2、Ri+1-Ri

3、>eps2i+42j−1b−aT2^(i+5)=T2^(i+4)/2+j=1f

4、(2i+5+1+?)Si+4=4Ti+5/3-Ti+4/3Ci+3=16Si+4/15-S3/15Ri+2=16Ci+3/15-Ci+2/15i=i+1输出Ri+1程序代码：%定义函数functionf=sin1(x)if(x==0)f=1;elsef=sin(x)/x;end%龙贝格算法求数值积分clc;clearall;flag=zeros(1,4);In=zeros(25,4);eps1=input('whatprecisiondoyouwanttoget?');i=0;%计算一些初始序列while(flag(4)<2)fork=flag(1)+1:i+4forj=1:2^(k

5、-1)In(k,1)=In(k,1)+sin1(a+(2*j-1)*(b-a)/(2^k));endifk>1In(k,1)=In(k,1)/2^k+In(k-1,1)/2;elseIn(k,1)=(b-a)*(sin1(b)+sin1(a))/2;endendflag(1)=i+4;fork=flag(2)+1:i+3In(k,2)=4*In(k+1,1)/3-In(k,1)/3;endflag(2)=i+3;fork=flag(3)+1:i+2In(k,3)=16*In(k+1,2)/15-In(k,2)/15;endflag(3)=i+2;fork=flag(4)+1:i+1I

6、n(k,4)=64*In(k+1,3)/63-In(k,3)/63;endflag(4)=i+1;i=i+1;end%根据递推式计算龙贝格序列i=1;while(abs(In(i+1,4)-In(i,4))>eps1)forj=1:2^(i+4)In(i+5,1)=In(i+5,1)+sin1(a+(2*j-1)/(2^(i+5)));endIn(i+5,1)=In(i+5,1)/2^(i+5)+In(i+4,1)/2;In(i+4,2)=4*In(i+5,1)/3-In(i+4,1)/3;In(i+3,3)=16*In(i+4,2)/15-In(i+3,2)/15;In(i+2,4

7、)=64*In(i+3,3)/63-In(i+2,3)/63;i=i+1;end%输出高精度计算值vpa(In(i+1,4))程序执行结果：