Peculiar behavior regarding matrix operations
2 views (last 30 days)
Show older comments
Hi,
I just discovered something that I found really weird when I was working with some element-wise multiplication and division for large matrices. I found that the order in which I carry out the operations gives me different results, the line of code giving this behavior is the following:
efrog2=efrog./abs(efrog).*spectrogram;
where efrog and spectrogram are 1325x1325 (efrog contains complex values) matrices. Now if I change the order to:
efrog3=spectrogram.*efrog./abs(efrog);
I get a different result, I if look at the maximum difference I get:
max(max(abs(efrog3-efrog2)))=6.0024e+08
However, if I change the order to:
efrog4=efrog.*spectrogram./abs(efrog);
I get:
max(max(abs(efrog3-efrog4)))=0
I tried the same thing with small 4x4 matrices, but for them the order (as one would suspect) didn't matter. Does anybody have any idea about what is going on here?
Thanks in advance
Cheers
Robert
2 Comments
Walter Roberson
on 1 Feb 2016
what is max(abs(efrog(:)), max(spectrogram(:)) ?
Is it possible that 6.0024e+08 is on the order of max(eps(efrog(:)) or max(eps(spectrogram)) ?
Accepted Answer
John D'Errico
on 1 Feb 2016
Edited: John D'Errico
on 1 Feb 2016
Congratulations. You are the 10 millionth person to have discovered that floating point operations can have subtly different results depending on the order in which those operations were performed. Of course, you award for having done so will be your name appearing in USA today, the Washington Post, and Time magazine, to name just a few. No money though.
Essentially, the basic properties of arithmetic that you learned in grade school no longer EXACTLY apply in floating point arithmetic. In fact, that which you do on a computer only resembles/approximates true mathematics when done in floating point arithmetic.
2 Comments
Walter Roberson
on 1 Feb 2016
There is no "correct" order, really. But division has more potential to lose precision than multiplication does, so best do that last provided that you don't need to divide in the middle to avoid overflow.
More Answers (0)
See Also
Categories
Find more on Matrix Indexing 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!