Simpson 1/3 rule

Question

The solution of the problem it is 1.732400451459101 for Simpson 1/3 Rule. Instead the solution that the program give me is 1.73239801

Can anyone help me out? Thanks in advance.

clc
clear
close all

f = @(x) sin(x);
a = 0.1; 
g = a;
b = 2.4;
k = 19;   

n = 2*k;
S = 0;
h = (b-a)/n;

for i=1:k
    S=S+(h/3)*(f(a)+4*f(a+h)+f(a+2*h));
    a=a+2*h;

end

fprintf('La integral se aproxima a: %0.8f \n',S)

syms x

y = sin(x);

InT = int(y,g,b);

InT = double(InT)

Simpson's rule computes the integral numerically, that is, it's an _approximation_ to the true value of the integral. To improve accuracy you can increase the number of steps `k` — Luis Mendo, Mar 10 '18 at 00:42
"The solution of the problem it is 1.732400451459101 for Simpson 1/3 Rule." How? For how many steps? With what numerical precision in between steps? When the two answers are this close, frequently they are not really different. — Argyll, Mar 10 '18 at 01:06

Hunter Jiang · Answer 1 · 2018-03-12T16:28:42.093

The error in approximating an integral by Composite Simpson's rule is:

which is about 1.7149e-07 in f(x)=sin(x) case, and that means the absolute error bound is 9.8990e-08, which is palatable to me.

Besides, here is an alternative to the code above:

f = @(x) sin(x);
[a,b,g,k]=deal(0.1,2.4,0.1,19);
[S,n,h]=deal(0,2*k,(b-a)/(2*k));
for i=1:k
    S=S+(h/3)*(f(g)+4*f(g+h)+f(g+2*h));
    g=g+2*h;
end

Or, we could just call:

f = @(x) sin(x);
[a,b,k]=deal(0.1,2.4,19);
Int = a:(b-a)/(2*k):b;
S=(b-a)/(6*k) * ( f(Int(1)) + 4*sum(f(Int(2:2:end-1))) ...
                + 2*sum(f(Int(3:2:end-2))) + f(Int(end)));

Simpson 1/3 rule

1 Answers1