Add columns to an uninitialized matrix

Dear all,

if I run the following code

clear all
n = 200000;
rep = 100;
x = randn(n,1);
counter = zeros(1,rep);
%x_new = zeros(n,rep);
for i = 1:rep
  tic
  x_new(:,i) = x;
  counter(i) = toc;
end
sum(counter)
plot(1:rep,counter)

I see a linear growth of the elapsed times at each iteration in octave while I see flat plot with four peaks in matlab. The overall cost is much smaller in matlab (0.16 seconds vs 3.16). Is it expected behavior? I know that I should initialize x_new, but is some circumstances the final size is not known (think to integration of ODEs with variable stepsize).

Thanks,

Marco

That is not surprising, but there are two very easy things you can do to speed it up. The first is to grow your array in reverse direction: for i = rep:-1:1 instead of for i = 1:rep, which reduces it for me from 2.02 seconds to only 50 milliseconds (i.e., 40 times faster). This works because it allocates memory only once, instead of in every iteration.

In your code you know the value of rep so preallocation should not be a problem, but if you don’t know it at all, you can use the second approach, which is to preallocate every 50 or so iterations, somewhat like this:

x_new = zeros(n, 50);
pos = 0;
while (<condition>)
  x_new (:, ++pos) = x;
  if (pos >= columns (x_new))
    x_new(:, end+50) = 0;   ### this reduces resizing to 1/50 of the original
  end
end
x_new = x_new (:, 1:pos);
2 Likes

Yes this is expected behavior. Years ago you would have seen similar linear growth in Matlab. They have made certain optimizations including a Just-in-Time compiler and other Accellerator features that make a significant difference in certain problems, including array growth.
see:

Octave has had several attempts to introduce similar enhancements, but none that have yet proven effective and stable, and it’s still on the developer wish-list. It is expected the next developer meeting will discuss some accellerator options, as per VM for Octave interpreter

2 Likes