Simpson's 1/3 Rule Implementation Issue (Code Provided)

9 Ansichten (letzte 30 Tage)
N/A
N/A am 2 Apr. 2020
Beantwortet: David Hill am 2 Apr. 2020
I have written the following code to Implement Simpson's 1/3 Rule, but it does not seem to give the correct answer. I am unsure what I have done wrong in my code.
MATLAB CODE:
function Integration(a,b,x,n,true)
a=0; %Lower Integral Bound
b=1; %Upper Integral Bound
x=[a,.25,.5,.75,b];
n=4; %Even Number of Segments
f=zeros(n+1,1);
true=.602298; %True Value of Equation
h=(b-a)/n; %width of each segment
f(1)=(x(1)^0.1)*(1.2-x(1))*(1-exp(20*(x(1)-1))); %Value For Initial Bound
for i=2:n
f(i)=2*(x(i)^0.1)*(1.2-x(i))*(1-exp(20*(x(i)-1))); %Values for In-Between Bounds
end
f(n+1)=(x(n+1)^0.1)*(1.2-x(n+1))*(1-exp(20*(x(n+1)-1))); %Value for Final Bound
q=sum(f,'all');
I=(h/2)*(q); %Width Times Average Height
disp('USING TRAPEZOIDAL')
fprintf('I=')
disp(I)
et=abs((true-I)/true)*100;
fprintf('True Error=')
disp(et)
disp("USING SIMPSON'S 1/3")
f2=zeros(n,1);
f2(1)=(x(1)^0.1)*(1.2-x(1))*(1-exp(20*(x(1)-1)));
for j=2:2:n-2
f2(j)=2*(x(j)^0.1)*(1.2-x(j))*(1-exp(20*(x(j)-1)));
end
for j=3:2:n-1
f2(j)=4*(x(j)^0.1)*(1.2-x(j))*(1-exp(20*(x(j)-1)));
end
f2(n)=(x(n)^0.1)*(1.2-x(n))*(1-exp(20*(x(n)-1)));
v=sum(f2,'all');
I2=(h/3)*v;
fprintf('I=')
disp(I2)
%ea=abs((true-I2)/true)*100;
%printf('True Error=')
%disp(ea)
%Nobody Cares, Work Harder
%Keep Hammering
end

Akzeptierte Antwort

David Hill
David Hill am 2 Apr. 2020
The following function might help.
function o = simp(a,b,f,n)
x=linspace(b,a,n+1);
o=0;
for k=1:n
xx=[x(k);(x(k)+x(k+1))/2;x(k+1)];
y=f(xx);
c=([xx.^2,xx,ones(3,1)]\y)';
p=polyint(c);
o=o+polyval(p,x(k))-polyval(p,x(k+1));
end
Evaluate the function using your equation desiring integration.
f=@(x)x.^.1.*(1.2-x).*(1-exp(20*(x-1)));
output=simp(0,1,f,4);

Weitere Antworten (0)

Kategorien

Mehr zu Numerical Integration and Differential Equations finden Sie in Help Center und File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by