2

Okay so I'm doing numerical integration using the trapezoid rule. I wrote the m file as such:

function y = trap(f,a,b,N)  
h= (b-a)/N; y=0;  
for i=1: N-1  
    x=(a+i*h);y=y+eval(f);  
end  
y=2*y; x=a; y=y+eval(f); x=b; y=y+eval(f);  
y= (h/2)*y;

It says that the value assigned to x might be unused and I don't know why. If someone could tell me what I am doing wrong that would be awesome.

costrom
  • 443
blargh
  • 21

2 Answers2

1

You need to pass $x$ into the function $f$ for evaluation; MATLAB is giving you the "X may not be used" warning because $x$ is changing in every loop without ever being passed into $f$

I have made a couple changes to the code (for readability, changing eval to feval...)

function y = trap(f,a,b,N)  
h = (b-a)/(N);
y=0;
for i=1:N-1
    x=(a+i*h);
    y=y+feval(f,x);
end
y = y+feval(f,a)+feval(f,b)
y = h*y;
costrom
  • 443
1

I'm assuming everything else is OK, just considering performance in Matlab. If (1:N-1) is many values you will be much better off storing the vector than using a for loop.

function y = trap(f,a,b,N)  
h = (b-a)/(N);
y=0;
% Instead of a for loop we make a vector v storing the values (1,2,...N-1):
v = 1:N-1;
x = (a+v*h);
y = sum(feval(f,x));

y = y+feval(f,a)+feval(f,b); 
y = h*y;

Avoiding for-loops can easily give a speed up of 100-1000 times or more. Maybe it doesn't matter right now, but it will start mattering once you build more advanced and useful stuff.

Edit As mentioned by costrom in the comments : The above assumes that function f can take an array of inputs. If it can not do so, we can still replace the line

y = sum(feval(f,x));

with

y = sum(arrayfun(f,x));

for example the function

f = @(x) x^2;

feval(f,x);

will give error if x is a vector , but the following will work:

arrayfun(f,x);

Community
  • 1
mathreadler
  • 25,824
  • you can combine the endpoints of the integration by changing v to 0:N. the implementation in the answer will work much faster than OP, but with the caveat that the function f has to be able to accept a vector of $x$ – costrom Dec 09 '15 at 17:10
  • yes and if it does not accept a vector we can do an arrayfun. – mathreadler Dec 09 '15 at 17:50