5
\$\begingroup\$

I've implemented the trapezoidal rule to compute the integral for a function \$x^2\$. I would like to see another style of the same code. It seems Matlab hates for a matrix to be expanded without specifying the size.

clear all 
clc
n = 5; % number  subdivisions
a = 1.0; % lower limit
b = 2.0; % upper limit
sum = 0.0;
dx = (b-a)/(n-1); % step size

% Generate samples
for i = 1:n
    x(i) = a + (i-1)*dx; 
end

% Generate function's values
for i = 1:n
    y(i) = x(i).^2;
end

% Compute area using Trapz. method 
for i = 1:n
    if ( i == 1 || i == n) % for the first and last data
        sum = sum + y(i)./2;
    else
        sum = sum + y(i); % for the rest of data
    end
end
area = sum * dx;


area
\$\endgroup\$
8
  • \$\begingroup\$ Your code does not appear to be correct as the integral of \$x^2\$ from 1 to 2 using 5 areas (4 subdivisions) is supposed to be 2.34 exactly, while your code seems to yield the answer for 4 areas (2.3438) \$\endgroup\$
    – mleyfman
    Aug 21, 2014 at 6:07
  • \$\begingroup\$ @mleyfman, don't forget the error. The result 2.3333 not what you said 2.34 is using the analytical form of the function however trapezoidal method is a method to approximate the analytical function. The error here is 0.45 % which means we have some level of uncertainty in the estimated result. \$\endgroup\$
    – CroCo
    Aug 21, 2014 at 6:14
  • 1
    \$\begingroup\$ I was going by this: nastyaccident.com/calculators/calculus/trapezoidalRule. The integral is indeed 2.333 if done analytically. Trapezoidal rule with n=5 should yield 2.34 exactly, whereas n=4 should yield 2.34375 \$\endgroup\$
    – mleyfman
    Aug 21, 2014 at 6:17
  • \$\begingroup\$ @mleyfman, according to the link you gave Answer: 2.34375 which is same of mine. Even though n=5 in my code but the loop start from 0 to 4. Check the implementation of my code. \$\endgroup\$
    – CroCo
    Aug 21, 2014 at 6:20
  • \$\begingroup\$ Are you sure you ran it with amount=5? Also see this: wolframalpha.com/input/… \$\endgroup\$
    – mleyfman
    Aug 21, 2014 at 6:28

2 Answers 2

6
\$\begingroup\$

Avoid variable names that collide with built-in functions

sum is a built-in function. Redefining it causes unconventional behavior:

>> sum([1 2 3])

ans =

     6

>> sum = 0.0

sum =

     0

>> sum([1 2 3]) 
Index exceeds matrix dimensions.

>> clear all
>> sum([1 2 3])

ans =

     6

Idiomatic MATLAB

The whole point of using MATLAB, rather than C, is that the language is designed to work on many values at once. Therefore, you should avoid looping wherever possible, and let MATLAB do its job on vectors.

% Givens
n = 5                                  % number of subdivisions
a = 1.0                                % x coordinate left endpoint
b = 2.0                                % x coordinate right endpoint

x = linspace(a, b, n)                  % x coordinates of samples
y = x .^ 2                             % y coordinates of samples
sum_y = sum(y) - (y(1) + y(end)) ./ 2  % sum y coordinates, taking just half of endpoints
dx = (b - a) / (n - 1)                 % width of each trapezoid
area = sum_y * dx

Reinventing the wheel

For your information, MATLAB has built-in support for numerical integration:

>> integral(@(x) x .^ 2, 1, 2)

ans =

    2.3333
\$\endgroup\$
0
3
\$\begingroup\$

Setting aside that I don't know MATLAB, I would like to mention some things

Naming

If you need a comment to explain a variable then the variable is poorly named.

n = 5; % number subdivided areas
a = 1.0; % lower limit
b = 2.0; % upper limit
dx = (b-a)/(n-1); % step size  

should become

subdividedAreaCount = 5;
lowerLimit = 1.0;
upperLimit = 2.0;
stepSize = (upperLimit - lowerLimit) / (subdividedAreaCount -1 );  

Refactoring

Now the loops for generating samples and function values would look like this if we also rename x and y to sampleValue and functionValue.

% Generate samples
for i = 1:subdividedAreaCount
    sampleValue(i) = lowerLimit + (i-1)*stepSize; 
end

% Generate function values
for i = 1:subdividedAreaCount
    functionValue(i) = sampleValue(i).^2;
end  

The computing of the area would look like

% Compute area using Trapz. method 
for i = 1:subdividedAreaCount
    if ( i == 1 || i == subdividedAreaCount) % for the first and last data
        sum = sum + functionValue(i)./2;
    else
        sum = sum + functionValue(i); % for the rest of data
    end
end  

but we can refactor this also to exclude the if condition like

% Compute area using Trapz. method 
for i = 2:subdividedAreaCount-1
    sum = sum + functionValue(i); % for the rest of data
end 
sum = sum + functionValue(2)./2;
sum = sum + functionValue(subdividedAreaCount)./2;
\$\endgroup\$
2
  • \$\begingroup\$ +1 for naming variables. Is there any benefit of refactoring if statement? I think mine is more readable than the approach you stated. \$\endgroup\$
    – CroCo
    Aug 21, 2014 at 6:08
  • \$\begingroup\$ This just depends how big subdividedAreaCount can be. If subdividedAreaCount = 100000 it would be at least faster without using the if . \$\endgroup\$
    – Heslacher
    Aug 21, 2014 at 6:18

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.