Improve Speed of Piecewise Function in MATLAB

Question

I have a small piecewise function that profiling reveals is taking 60% of the runtime of the program. It is called very often because it goes within some integrals that I perform quite a lot in my code.

According to profiling, it is called 213560 times, taking 47.786 s in total, corresponding to ~220 microseconds per call.

I want to pass it an array, and it should return an array, operating element wise.

I know that using loops in Matlab is very slow and should be avoided, but I'm not sure how to vectorise this sort of function.

function bottleradius = aux_bottle_radius(z_array)
%AUXBOTTLERADIUS Radius of aux bottle
%   This returns the radius of the aux bottle as a function of z. It works for all
%   heights of aux bottle, just remember to integrate over the right height
%   range
    bottleradius = zeros(size(z_array));
    for i = 1 : max(size(z_array))
        if z_array(i)<-30e-3
            %door cavity
            bottleradius(i) = 34e-3;
        elseif z_array(i)>=-30e-3 && z_array(i)<-20e-3
            %radiussing door cavity
            bottleradius(i) = 34e-3 + 10e-3 - sqrt((10e-3).^2 - (z_array(i)+30e-3).^2);
        elseif z_array(i)>=-20e-3 && z_array(i)<-10e-3
            %aluminium plate
            bottleradius(i) = 46e-3;
        elseif z_array(i)>=-10e-3 && z_array(i)<0e-3
            %radiussing aluminium plate to main bottle
            bottleradius(i) = 46e-3 + 10e-3 - sqrt((10e-3).^2 - (z_array(i)+10e-3).^2);
        elseif z_array(i)>=0e-3
            %top of Al plate, bottom of main bottle
            bottleradius(i) = 185e-3;
        else
            bottleradius(i) = 0;
        end
    end
end

Answer 1

You can do that completely vectorized with logical operators. You can essentially replace that code with:

function bottleradius = aux_bottle_radius(z_array)

%// Declare initial output array of all zeroes
bottleradius = zeros(size(z_array));

%// Condition #1 - Check for all values < -30e-3 and set accordingly
bottleradius(z_array < -30e-3) = 34e-3;

%// Condition #2 - Check for all values >= -30e-3 and < -20e-3 and set accordingly
ind = z_array >= -30e-3 & z_array < -20e-3;
bottleradius(ind) = 34e-3 + 10e-3 - sqrt((10e-3).^2 - (z_array(ind)+30e-3).^2);

%// Condition #3 - Check for all values >= -20e-3 and < -10e-3 and set accordingly
bottleradius(z_array >=-20e-3 & z_array < -10e-3) = 46e-3;

%// Condition #4 - Check for all values >= -10e-3 and < 0 and set accordingly
ind = z_array >=-10e-3 & z_array < 0;
bottleradius(ind) = 46e-3 + 10e-3 - sqrt((10e-3).^2 - (z_array(ind)+10e-3).^2);

%// Condition #5 - Check for all values >= 0 and set accordingly
bottleradius(z_array >= 0) = 185e-3;
end

Minor comments

1. 0e-3 doesn't make any sense precision wise. This is essentially the same as 0 and I've changed that in your code.
2. Note that for conditions #2 and #4, I precompute a logical array that indicates where we would need to access the corresponding values in z_array to make things easier and set those same locations in bottleradius to be the desired computed outputs. I don't do this for the other conditions because you're just setting them to a single constant.
3. Thankfully, you use element-wise operators for conditions #2 and #4 so there wasn't a need to change the expressions for those conditions.