4 views (last 30 days)

Show older comments

Essentially I want a fast way of doing: c = (A^k)*b0. But I want the result for multiple values of k (I don't need it for all values of k, just some).

At the moment, I am just doing this in a normal for loop (b1 = A*b, b2 = A*b1, b3 = A*b2, etc.) for all k. But I am wondering if there is a faster way (maybe using GPUs).

Doing loops in GPUs doesn't seem like the way forward. I was thinking I could just request c = (A^k)*b0 (which is very fast on the GPU) for only the k that I want, but if I want many (for example for k = [1:5:1000]) this still ends up being slower than just doing it on a loop on the CPU.

Any suggestions? Thanks -

clear; rng(1);

N = 301;

k = 1000; A = randn(N)/17; b = rand(N,1);

f = @() r(A,b,k);

t = timeit(f);disp(t)

function b = r(A,b,k)

for ix = 1:k

% currently not pulling out b for the k of interest (but is doable here)

b = A*b;

end

end

James Tursa
on 17 Jun 2018

Edited: James Tursa
on 17 Jun 2018

Note that as the power value gets higher, the numerical stability of successive matrix multiplies will degrade and you will not get as accurate an answer as calling the MATLAB mpower function directly (which uses the matrix exponential function). Also, I am not sure all of those successive matrix multiplies will be faster than calling mpower directly anyway. E.g., what kind of timing do you get with this compared to your looping?

X = your matrix

b = your vector

k = 1:5:1000;

tic

result = arrayfun(@(p)X^p*b,k,'uni',false);

toc

Walter Roberson
on 17 Jun 2018

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!