MATLAB runs significantly slower after storing value in a global search

Question

Suppose that I have a function defined on 0<= x<= pi, 0<=y, z<=2pi. I have written a matlab code to global search the maximum point, with step length 0.01. I also restored the function value. The code is (of course the function is more complicated than x+y+z. This is simply an example.)

x=0;
y=0;
z=0;
i=1;

Maxdeviation= ComputeDeviationSum(0,0,0,mufT, mufD, mufC);
MaximizingArg= [x, y,z] ;
while x<= pi
    y=0;
    z=0;
    while y<= 2*pi
        z=0;
        while z<=2 *pi
            deviation= x+y+z; % Just an example. The function is more complex.
            if deviation> Maxdeviation
                MaximizingArg= [x, y,z];
                Maxdeviation= deviation;
            end
            a(i,:)=[x, y,z, deviation];
            z=z+0.01;
            i=i+1;
        end
        y=y+0.01;
    end
    x=x+ 0.01;
end

And I found this really slow. After I deleted the three lines restoring the value, i.e. deleted

i=1;
a(i,:)=[x, y,z, deviation];
i=i+1

these three lines, it became much faster. (I stopped the execution after 20 seconds and found the latter did 10 times more iteration than the former.)

My question is: Why after stopping restoring value, the code is much faster? What should I do to let the code be faster if I want to restore the value for each point?

Answer 1

The comments on your question ask about whether you see a warning in your editor about the array size increasing every loop, since that's what makes those lines slow. This is what they're referring to:

You can work out up-front how big a is going to be and pre-allocate it. Using for loops makes this a bit clearer, so I've made both changes in the code below:

Maxdeviation= ComputeDeviationSum(0,0,0,mufT, mufD, mufC);
MaximizingArg= [0,0,0];

% Predefine sweep points
X = 0:0.01:pi;
Y = 0:0.01:2*pi;
Z = 0:0.01:2*pi;
i=1;

% Pre-allocate output since we know how big it's going to be
a = NaN(numel(X)*numel(Y)*numel(Z),4);

% Loop over each point array
for ix = 1:numel(X)
    x = X(ix);
    for iy = 1:numel(Y)
        y = Y(iy);
        for iz = 1:numel(Z)
            z = Z(iz);

            deviation= x+y+z; % Just an example. The function is more complex.
            if deviation> Maxdeviation
                MaximizingArg= [x, y,z];
                Maxdeviation= deviation;
            end

            % Assign output at this index to preallocated array
            a(i,:) = [x,y,z,deviation];
            i = i+1;
        end
    end
end

Taking this further, it might be worth noting that storing x,y,z in the output array might be redundant, you could instead pre-allocate a=NaN(numel(X),numel(Y),numel(Z)) and then store the function value as a(ix,iy,iz)=deviation;. Then you can easily index into the resulting array for given indices of x, y and z.

Conversely if you don't want to use the indices ix,iy,iz for anything you could simplify the loops to be something like for x = X and remove the assignments x = X(ix)