Interesting: say you have a 100x100 matrix of doubles. Then you put values in each row from -1 to 1, equally spaced, such that all 100 rows for each column i contain the value:

-1 + i * ((1 - -1) / (100 - 1)) == -1 + i * 2/99

Where the columns are numbered from 0 to 99. By definition, the sum of all elements in this matrix should be 0. For integer numbers, this is true. For floating point numbers, this is tricky.

As a natural born hax0r, you should know that there is something called floating point errors. As the entropy of these values isn’t particularly high, operations on them might result in an error that will be large. But I never thought it could be that different. Watch this [1]:

• Summing all elements individually gives a result of: 3.108624468950438e-15. A large error, this is what we would expect.

However:

• Taking the sum of an individual row (doesn’t matter which one, of course) results in: 3.108624468950438e-15. The same number as summing all elements individually!
• Taking the sum of these rows gives: 3.108624468950438e-13. This is what we would expect, 100 times the sum of a row. The funny thing of course, is the fact that this number is a 100 times different from the sum of all elements.
• Taking the sum of the columns gives: -2.415845301584341e-13. A completely different number than in the previous cases.

This issue surfaced in a debugging session where we were comparing the result of an algorithm written in Matlab and the same algorithm converted to C++. Although I have a copy of “What Every Computer Scientist Should Know About Floating Point Arithmetic”, which tells you that the addition operator + is not necessarily associative for floating point values, I still had to blink a couple of times when I saw the above results for the first time.

Lesson learned: be careful when you take the sum of (a matrix of) floating point numbers as a comparison measure —- or maybe even more general: beware: floating point number calculations ahead.