MATLAB: Having issues with Simpson’s Rule ln(x)

integralMATLABnatural lognatural logarithmplottingplotting functionssimpsons rule

Hi,

Before I found out MATLAB had an integral function, I created my own integral function by using the composite Simpson's Rule.

The code not only finds the integral but also plots the curve of the integral. I have tried many functions and so far they have all excecuted properly with the exception of the function ln(x). When I include this integral, the function is solved and I obtain the correct "area under the curve". However, when I try plotting the function, it will only plot the last two points of the array which correspond to the last x value (length of the step size of the integral). I have included two codes, the first one using MATLAB's function so you know how the integral should look like. The second one is the script I wrote that solves the Simpson's Rule.

% Function I am integrating
function fval = myFunInt(x)
fval = log(x);
end

Using MATLAB's function

clear
a=1;            % lower limit

b=30;           % upper limit

n=1000;         % subintervals

h = (b-a)/n;    % Spacing

x = 0:30;
z = zeros(1,n+1);
int = zeros(1,n+1);
for j = 0:n    
    x_j=a+j*h;      % x values are being allocated in the empty array

    x(:,j+1)=x_j;
    fun = @myFunInt;
    y=integral(fun,a,x(:,j+1));
    int(:,j+1)=y;
end
plot(x,int);

Using Simpson's Rule

clear
disp('******************************************************************');
a=1;            % lower limit
b=30;           % upper limit
n=1000;         % subintervals
h = (b-a)/n;    % Spacing
% Creating an empty array where all the arguments will be generated. 
% These arguments are the x values that will be inputted into myFunInt(x).
% The x values will be created in the for loop below where x_j is defined
x = zeros(1,n+1);
z = zeros(1,n+1);
for j=0:n
    x_j=a+j*h;      % x values are being allocated in the empty array
    x(:,j+1)=x_j;
%***************************************************************************
% Creating an empty array where the even values will be allocated.
% This term is enabled when n>2.
even_term = zeros(1,n/2-1);
if  j>=4
    i = 1:j/2-1;
    even_term(:,i) = 2.*myFunInt(x(:,2*i+1));
end
%**************************************************************************
% Creating an empty array where the odd values will be allocated. 
odd_term = zeros(1,n/2);
if j>=2
    k = 1:j/2;
    odd_term(:,k) = 4.*myFunInt(x(:,2*k));
end
%All the terms in the function are now added to obtain the final
% integral value (area under the curve) of the function.
first_term = myFunInt(x(:,a+1));
last_term = myFunInt(x(:,n));
final_integral = h./3*(first_term+sum(even_term)+sum(odd_term)+last_term);
z_j = final_integral;
z(:,j+1)=z_j;
end
plot(x,z);

Any recommendation is appreciated.

In the end I know I can end up just using MATLAB's integration function but I would rather know why my code is not properly plotting the integral of ln(x).

Thank you!

Guillermo Naranjo

MATLAB: Having issues with Simpson’s Rule ln(x)

Best Answer

Related Question