Computing determinants of a 3D array
10 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
Yanir Hainick
il 26 Nov 2012
Modificato: Matt J
il 8 Apr 2020
Let’s say I have an NxNxL array. L is typically 10^4-10^5, and N is typically 10^0-10^1.
My goal is to calculate a vector of length L, where the i-th cell contains the determinant of (:,:,i).
I currently use for-loop, as det(A) accepts 3D arrays with the last dimension being a singleton, so this code works:
for i = 1:L
Vec(i) = det(Mat(:,:,i));
end
However, it seems weird that i can't implement this in a vectorial fasion. Can anybody think of any way to get rid of the for loop? Note that the format of an 3D array is pretty stiff, i.e. i can't change the input to cell array and use cellfun.
Thanks!
Yanir.
2 Commenti
Matt J
il 26 Nov 2012
Hopefully, the reason you're asking for this is not for the purpose of solving many linear systems.
Risposta accettata
Sean de Wolski
il 26 Nov 2012
What's wrong with the for-loop?
n=100;
pages = 1e4;
X = rand(n,n,pages);
D = zeros(pages,1);
tic;
for ii = 1:pages
D(ii) = det(X(:,:,ii));
end
toc;
%Elapsed time is 1.805616 seconds.
1 Commento
James Tursa
il 27 Nov 2012
My 2 cents to all:
The performance advantage of vectorized approaches to this is that one wishes to avoid the data copy involved with the X(:,:,ii) slices and the overhead of the loop. That being said, this data copy & loop overhead is in all likelihood swamped by the numerical calculations (and possible data copy) involved in the determinant calculation itself. So even if one were to get a vectorized one-liner to this that did not involve explicit slices, it probably wouldn't run significantly faster (if at all) than straight forward loops (the small size explicit code e.g. that Matt shows excepted).
Più risposte (2)
Matt J
il 26 Nov 2012
Modificato: Matt J
il 8 Apr 2020
For small N, it would be an advantage to vectorize the determinant formula explicitly (example for for N=2 below). For larger N, maybe you could do the same thing recursively.
%fake data
N=2;
L=1e5;
Mat=rand(N,N,L);
tic
Vec=zeros(1,L);
for ii=1:L
Vec(ii)=det(Mat(:,:,ii));
end
toc;
%Elapsed time is 0.194555 seconds.
Mat=reshape(Mat,[],L);
tic;
Vec =Mat(1,:).*Mat(4,:) - Mat(2,:).*Mat(3,:);
toc
%Elapsed time is 0.000378 seconds.
2 Commenti
Pi Ting
il 8 Nov 2017
The line
Vec =Mat(1,:).*Mat(3,:) - Mat(2,:).*Mat(4,:);
should be
Vec =Mat(1,:).*Mat(4,:) - Mat(2,:).*Mat(3,:)?
Al in St. Louis
il 8 Apr 2020
I had to use Pi Ting's expression to get the correct answers. This is exactly what I need to do!
Vedere anche
Prodotti
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!