multi array vectorization problem (reformulated)
1 view (last 30 days)
Show older comments
How to vectorize (for-loop elimination) the following code?
% init
A = [ -0.8013 -0.4981; -0.2278 -0.9009];
t = 0:0.01:1e2;
% eigen space
[V,D]=eig(A,'vector');
D = reshape(D,1,[]);
% computing
expmAt = zeros([size(A),length(t)]);
for i = 1:length(t)
expmAt(:,:,i) = V.*(exp(D*t(i)))/V;
end
8 Comments
Accepted Answer
Matt J
on 19 May 2022
Edited: Matt J
on 19 May 2022
runtest(1e2)
function runtest(N)
A = [ -0.8013 -0.4981; -0.2278 -0.9009];
t = 0:0.01:N;
b = [2;3];
[V,d]=eig(A,'vector');
tic;
E=exp(d*t);
expmAtb=repmat(eye(size(V)),1,1,numel(t));
expmAtb(logical(expmAtb))=E(:);
expmAtb3=squeeze(pagemtimes( pagemtimes(V,expmAtb), V\b));
toc
tic;
expmAtb4=V*(exp(d*t).*(V\b));
toc
end
More Answers (2)
Bruno Luong
on 19 May 2022
A = [ -0.8013 -0.4981; -0.2278 -0.9009];
t = 0:3;
% eigen space
[V,D]=eig(A,'vector');
D = reshape(D,1,[]);
% computing
expmAt = zeros([size(A),length(t)]);
for i = 1:length(t)
expmAt(:,:,i) = V.*(exp(D*t(i)))/V;
end
expmAt
expmAt2 = pagemrdivide(V.*exp(D.*reshape(t,1,1,[])),V)
5 Comments
Catalytic
on 19 May 2022
Edited: Catalytic
on 19 May 2022
% init
A = [ -0.8013 -0.4981; -0.2278 -0.9009];
t = 0:0.01:1e2;
tic;
% eigen space
[V,D]=eig(A,'vector');
D = reshape(D,1,[]);
% computing
expmAt = zeros([size(A),length(t)]);
for i = 1:length(t)
expmAt(:,:,i) = V.*(exp(D*t(i)))/V;
end
toc
tic
out=theTask(A,t);
toc
[lb,ub]=bounds(expmAt-out,'all')
function expmAt=theTask(A,t)
[V,d]=eig(A,'vector');
E=exp(d*t);
expmAt=repmat(eye(size(A)),1,1,numel(t));
expmAt(logical(expmAt))=E(:);
expmAt=pagemtimes( pagemtimes(V,expmAt), inv(V));
end
3 Comments
Matt J
on 19 May 2022
I wouldn't expect the addition of squeeze() to impact timing significantly.
Also, I find it strange that pagemrdivide(A,B) is not optimized for the case where B has only one page.
See Also
Categories
Find more on Loops and Conditional Statements in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!