龙贝格积分-matlab通用程序

数值分析课程的大作业，教材《数值分析》李乃成.梅立泉

clear
clc
format long
f=input('请输入原函数f=','s');
a=input('积分下限a=');
b=input('积分上限b=');
eps1=input('精度eps1=');
T(1)=double((b-a)/2*(limit(sym(f),findsym(sym(f)),a)+limit(sym(f),findsym(sym(f)),b)));
for k=2:4
    sum1=0;
    for i=1:2^(k-2)
    sum1=sum1+subs(sym(f),findsym(sym(f)),(a+(2*i-1)*(b-a)/2^(k-1)));
    end 
    T(k)=1/2*T(k-1)+(b-a)/(2^(k-1))*sum1;
end
for k=1:3
    S(k)=T(k+1)+1/(4-1)*(T(k+1)-T(k));
end
for k=1:2
    C(k)=S(k+1)+1/(4^2-1)*(S(k+1)-S(k));
end