Is it possible to make this vectorized MATLAB code faster?

Question

In my program I need to calculate the sum:

.

I calculate this sum 2500 times with new values of C and z.

Argument z may be a vector. I wrote straightforward for-loop and vectorized version code as follows:

K = 200;
n_z = 40000;
C = ones(K,1); % an example, in real life the arey some coefficients (at next call will be new)
k = 0:K-1;
z = linspace(0, 2*pi, n_z); % at next call will be new

tic;
    my_sum_for = zeros(1, K);
    for i=1:n_z
       my_sum_for(i) = C' * tan(k' * z(i));
    end
toc; % Elapsed time is 1.820485 seconds.

tic;
     my_sum = C' * tan(k' * z);
toc; % Elapsed time is 0.160924 seconds.

Vectorized version is faster, but not enough. Is it possible to improve vectorized version?

After Dominique Jacquel's answer I have this vectorized version, it's faster:

    K = 200;
    n_z = 40000;
    C = ones(K,1)';  % an example, in real life they are some coefficients (at next call will be new)
    k = (0:K-1)';
    z = linspace(0, 2*pi, n_z);  % at next call will be new

    tic;
        my_sum_for = zeros(1, K);
        for i=1:n_z
           my_sum_for(i) = C * tan(k * z(i));
        end
    toc; % Elapsed time is 1.521587 seconds.

    tic;
         my_sum = C * tan(k * z);
    toc; % Elapsed time is 0.125468 seconds.

Is it possible to improve vectorized version even more (bsxfun, arrayfun or something)? The time of 250 seconds is still slow for me (it is 75% of all compuations).

Answer 1

I think you're pretty close to the hardware limits here. Matrix multiplications in Matlab are done with the BLAS libraries, which have proven tough to beat when it comes to performance.

AFAIK, the tangent function has actual dedicated hardware to compute its value. Also, Matlab automatically distributes trig functions of large matrices over multiple cores, so there's basically little to improve there.

Also, correct me if I'm wrong, but given the probable data overheads and memory issues, I think this computation would actually be slower on GPUs.

Now if you compare these:

tic;
for ii = 1:10
    my_sum = C * tan(k * z);
end
toc 

tic;
for ii = 1:10
    my_sum_notan = C * k * z;
end
toc

You'll see that all of the pain comes from the tangent function, so you had best focus on this. As you can read here, speeding trig functions up is basically only possible if you sacrifice some accuracy.

In all you should ask yourself these questions:

1. Are you willing to give up full double precision, or are, say, 6 digits "close enough"?

2. Can't you re-formulate the problem so that the tangent is computed afterwards? or before? Or anyway on a significantly smaller amount of elements? In the problem as stated this is obviously not possible, but I don't know the full code -- there might be some nice trig-identities that could apply to your problem.

3. Given all of the above, does the amount of effort needed to optimize this further really outweigh the longer runtime? 250 seconds does not sound too bad compared to writing custom, hardly portable MEX functions that implement crippled-but-fast trig functions.

Answer 2

Do as many matrix operations as possible upfront (transposition in this case) to save some time in the loop

K = 200;
n_z = 40000;
C = ones(K,1)';
k = (0:K-1)';
z = linspace(0, 2*pi, n_z);

tic;
    my_sum_for = zeros(1, K);
    for i=1:n_z
       my_sum_for(i) = C * tan(k * z(i));
    end
toc

tic;
     my_sum = C * tan(k * z);
toc;

my execution times before

Elapsed time is 1.266158 seconds.
Elapsed time is 0.531173 seconds.

and after

Elapsed time is 0.496803 seconds.
Elapsed time is 0.185396 seconds.

Answer 3

Here is a version that runs slightly faster on my computer:

k = repmat((0:K-1)', 1, n_z);
z = repmat(linspace(0, 2*pi, n_z), K, 1);
C = ones(1, K);
tic
my_sum = C*tan(k.*z);
toc

Essentially, instead of an outer product of k and z I operate directly on matrices.

First version

Elapsed time is 0.652923 seconds.
Elapsed time is 0.240300 seconds.

After Dominique Jacquel's answer

Elapsed time is 0.376218 seconds.
Elapsed time is 0.214047 seconds.

My version

Elapsed time is 0.168535 seconds.

You may have to add the cost of repmats, but maybe you can do that once only, I don't know the rest of the code.

I fully agree with Rody Oldenhuis. The majority of the work lies in the tangent function. I can tell you more. k.*z computation is very efficient and can not be improved much. If you calculate the memory bandwidth, it gets around 10GB/s on my computer. The peak I can get is around 16GB/s, so its close. Not many possibilities there. Same with C*T. That is a simple BLAS2 matrix-vector multiplication, which is memory bounded. For the system sizes you are showing the MATLAB overhead is not too big.

Edit: as Rody mentioned, new versions of MATLAB already do parallelize tan(). So not much here either.

You can only hope to improve the tan() - possibly by running it in parallel. After all, it is a trivially parallelizable task... Consider exporting just this to a MEX file, which will use OpenMP. Very simple work, lots of speedup if you have a few spare cores.