This is expected. Floating-point addition is non-associative and these codes add values in different orders. Julia's built-in sum function uses recursive pairwise summation, which is the most accurate of these three approaches but nearly as fast as linear scanning. A slower but more accurate algorithm is compensated summation <http://en.wikipedia.org/wiki/Kahan_summation_algorithm>. Even compensated summation is not completely accurate; completely accurate floating-point summation is tricky.
On Sat, Nov 15, 2014 at 11:52 AM, <[email protected]> wrote: > Dear all, > > Some (serious?) issue was raised on a mailing list concerning the "sum" > functions involved in computing the sum of elements of an array. > > It seems that we obtain different results for the three following > procedures: > > - using sum over the whole matrix > - using sum over the columns/rows of the matrix and summing these results > within a loop > - computing the sum with a double loop over all the elements of the matrix > > Attached is a file with the different procedures and their results for a > particular matrix. > > What is wrong with the implementation of the sum function? > This could have important consequences with iterative procedures. > > Many thanks, >
