搜索
计算方法上机报告27附录4龙贝格积分法的matlab程序clc;clear;s=input('请输入被积函数表达式：\nf(x)=','s');a=input('\n请输入积分区间下限：\n');b=input('\n请输入积分区间上限：\n');e=input('\n请输入积分误差：\n');N=-log10(e);f=inline(s);i=1;error=1;symsxf_a=limit(s,x,a);f_b=limit(s,x,b);T(1)=(b-a)/2*(double(f_a)+double(f_b));S(1)=0;C(1:2)=0;R(1:3)=0;whileerrore||i5i=i+1;sigma=0;forj=1:2^(i-2)sigma=sigma+f(a+(2*j-1)*(b-a)/2^(i-1));endT(i)=1/2*T(i-1)+(b-a)/2^(i-1)*sigma;S(i)=T(i)+1/3*(T(i)-T(i-1));ifi=3C(i)=S(i)+1/15*(S(i)-S(i-1));ifi=4R(i)=C(i)+1/63*(C(i)-C(i-1));error=abs(R(i)-R(i-1));endendendfork=1:iT_int=fix(T(k));T_dec=round((T(k)-T_int)*10^N);S_int=fix(S(k));S_dec=round((S(k)-S_int)*10^N);C_int=fix(C(k));C_dec=round((C(k)-C_int)*10^N);R_int=fix(R(k));R_dec=round((R(k)-R_int)*10^N);ifk==1fprintf('\n%0.f.%0.f\n',T_int,T_dec);elseifk==2fprintf('%0.f.%0.f%0.f.%0.f\n',T_int,T_dec,S_int,S_dec);elseifk==3fprintf('%0.f.%0.f%0.f.%0.f%0.f.%0.f\n',T_int,T_dec,S_int,S_dec,C_int,C_dec);elsefprintf('%0.f.%0.f%0.f.%0.f%0.f.%0.f%0.f.%0.f\n',T_int,T_dec,S_int,S_dec,C_int,C_dec,R_int,R_dec);endendfprintf('\n所求积分值为%0.f.%0.f,误差为%.e。\n',R_int,R_dec,e);