error in matrix multiplication
1 view (last 30 days)
Show older comments
Banafshe Shafieizargar
on 9 Sep 2020
Edited: James Tursa
on 9 Sep 2020
Hi
I'm getting a very confusing logic error in a MATLAB code.
two matrices K and D' are equal (isequal(K,D') gives 1)
while B*D' and B*K are not.
I know if I remove the permute and reshaping, then the problem is solved but this code here is a part of a larger project and it's not possible to do so.
The minimal code to reproduce the error is
clear all
close all
GT = rand(7,48,48);
GT_perm = permute(GT,[2 3 1]);
t_GT = rand(122,3);
B = rand(122,7);
param = [t_GT(:)',GT_perm(:)'];
K = reshape(GT,7,[]);
D = reshape(param(3*122+1:end),48^2,size(B,2));
isequal(K,D') %I receive 1 at this point which means these 2 variables are equal
ex_DB = B*D';
Image = B*K;
isequal(Image,ex_DB) %I receive 0 at this point which means the multiplication of the equal variables is not equal
sum(Image-ex_DB,'all')
0 Comments
Accepted Answer
Bruno Luong
on 9 Sep 2020
Edited: Bruno Luong
on 9 Sep 2020
MATLAB uses two different BLAS functions to perform
B*D'
and
B*K
In the first case, it use D and it never explicitly compute D'. Whereas in the second case, D is conjugate-transposed; assigned to K, then direct matrix product between B and K is computed.
That the reason why they are different (and also B*D' is faster than two statements K=D', B*K)
1 Comment
James Tursa
on 9 Sep 2020
Edited: James Tursa
on 9 Sep 2020
@Banafshe: Bottom line is that the individual calculations are done in a different order, so you can't expect the floating point results to be the same. Side Note: The same top level BLAS function is called for both (e.g. DGEMM for real double), but with slightly different inputs (as Bruno points out). In one case the pointer to the D data is passed with the transpose flag set, and in the other case the pointer to the K data is passed with the transpose flag not set.
More Answers (0)
See Also
Categories
Find more on Matrices and Arrays in Help Center and File Exchange
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!