Search code examples
plotoctavescilab

Why is there an inconstancy between these two ways of plotting a function in Octave?

Context

I have a function named ibm that I defined with a simple if condition: it should be 0 until a certain value (500), then is linear after that value:

function f = ibm(x)
  if(nargin != 1)
    usage("ibm(x)");
  endif
  if(x <= 500)
    f = 0;
  else
    f = (x-500) *0.02;
  endif
endfunction

Then in a script, I wanted to plot this function with x in range from 0 to 1000, so I tried to plot it with a simple

x=0:1:1000;
plot(x,ibm(x));

But the result did not show what I wanted: it showed a linear function without a flat part (see figure 1, the red graph). It did NOT correspond to the function results, and I checked some values below 500: if I type ibm(34) the result is 0 as expected, but not in the graph.
Why this discrepancy?

So after a while puzzling over this issue, I tried adding the result values "manually" in a for loop, and then the result became what I expected (see figure 2, the blue graph).

A complete code sample (MVE) here:

This produces both figures (wrong and right) with the same function ibm.

1; # start of file should not start with function, so useless statement here
# The function we want to plot
function f = ibm(x)
  if(nargin != 1)
    usage("ibm(x)");
  endif
  if(x <= 500)
    f = 0;
  else
    f = (x-500) *0.02;
  endif
endfunction

# For figure 1, directly computed (wrong plot result, in red)
x=0:1:1000;
figure(1);
plot( x, ibm(x), "xr");
xlabel("x");
ylabel("ibm(x)");
title("ibm(x) directly computed");
legend("ibm(x)");

# for figure 2, in for loop (correct plot result, in blue)
a = [];
for(i=x)
  a = [a, ibm(i)];
end
figure(2);
plot(x,a, "ob");
xlabel("x");
ylabel("ibm(x)");
title("ibm(x) in for loop");
legend("ibm(x)");

Question

Why the 1st (and most straightforward) method doesn't produce the graph result that I want? Should I always plot graphs with a for loop? It's as if the plot function doesn't take into account the if(x <= 500) and goes straight to the linear part.

As a note, I even tried the 1st way of plotting with Scilab, and ended up with the same result (red wrong curve), and I still don't know why.

Figures

Wrong

Function shown is linear and does not correspond to function output

Right

Correctly plotted function (with manual for)

Solution

  • The problem is your definition of ibm. You only coded it for a x of length 1. What does the code do for ibm([0,1000])? How does the if condition react to that x?

    I think this is a confusion of how does vectorization occurs. Octave will not call 1000 ibm each with its own individual x, but will call ibm once with all 1000 x as input. How you deal with vector inputs its on you, and your function simply doesn't.

    for example, if you replace the if with:

     f = (x-500) *0.02;
 f(x<500) = 0 ;

    the function should work (and therefore the plot should work), as both of those lines are vectorized and work for any size of x.

    If you want more info, check octaves manual for if.

    The condition in an if statement is considered true if its value is nonzero, and false if its value is zero. If the value of the conditional expression in an if statement is a vector or a matrix, it is considered true only if it is non-empty and all of the elements are nonzero. The conceptually equivalent code when condition is a matrix is shown below.

    if (matrix) ≡ if (all (matrix(:)))