组合辛普生公式求解积分.m

来自「拉格朗日插值多项式拟合,牛顿插值多项式,欧拉方程解偏微分方程,使用极限微分求解导」· M 代码 · 共 20 行

%____________________________________________________
%组合辛普生公式
%组合梯形公式求解积分
%说明：积分区间(a,b) ,被积函数f(x)，结果为h/3(f(a)+f(b))+2h/3∑f(x2k)+4h/3∑f(x2k-1),k=1到M-1,k=0,1,2,...M
%xk=a+kh,注意x0=a,x2M=b
%_____________________________________________________
function s=simpr1(f,a,b,M)
h=(b-a)/(2*M);
s1=0;
s2=0;
for k=1:M
    x=a+h*(2*k-1);
    s1=s1+feval(f,x);
end
for k=1:(M-1)
    x=a+h*2*k;
    s2=s2+feval(f,x);
end
s=h*(feval(f,a)+feval(f,b)+4*s1+2*s2)/3;