5

Assume we have a long array of doubles, say, N == 1000000. 

array<double, N> arr;

There are two naive approaches to compute the average. First 

double result = 0;
for (double x : arr) {
    result += x;
}
result /= arr.size();

This may be inaccurate when the sum of values is very big. Floating point numbers lose precision then.

Another approach is:

double result = 0;
for (double x : arr) {
    result += x / arr.size();
}

This may lose precision when the numbers are small.

Is there any fail-safe way to calculate a simple average of floating point numbers? Solutions, which use only the standard library are appreciated.

marmistrz
  • 5,974
  • 10
  • 42
  • 94
  • 3
    The claim is nonsensical, you don't lose precision based on the magnitude of the value. That's what the "float" in floating point takes care of. The number of significant digits that *double* can store does not depend on the value, you'll get this wrong only if you display more. – Hans Passant May 20 '17 at 11:01
  • 3
    @HansPassant There is a limit in the binary representation. When adding two floats, the exponents are matched. That means if you can only represent 4 decimal digits, adding 1.234e0 to 1.234e3 will first covert the first to 0.001e3, discarding the low digits. – Khouri Giordano May 20 '17 at 17:12

3 Answers3

4

If you want to squeeze more accuracy out of doubles, you can use Kahan summation and finally division by number of elements. There is however no standard library implementation of Kahan summation I know of.

An easy, standard way (almost like cheating) would of course be calculation using long doubles, basically using your first implementation and only converting the result back to double precision.

Community
  • 1
  • 1
Peter G.
  • 14,786
  • 7
  • 57
  • 75
2

The so-called naive ways are not naive. What do the data mean, and how accurately can you measure those values? Unless the answer is something very unusual, the simple method with doubles is fine. However floats are a bit under-powered for general use.

If you add the small absolute values first you might get an extra bit or so of precision. That requires a sort. If the data are all above a certain threshold, subtracting the minimum may also give you another bit.

You can also store a partial total, and a partial mean, and check at each stage that partial mean * number processed is within a certain tolerance of the partial total. That won't give you any extra accuracy, but it will tell you if the fpu is too inaccurate for your purposes.

You can also use long double, or even code your own extended-precision floating point library (or use someone else's). However the solutions get increasingly heroic.

Malcolm McLean
  • 6,258
  • 1
  • 17
  • 18
1

One way to reduce loss of precision would be to sort the doubles and then add them together in sorted order, starting with the smallest values and then at the end divide the final sum by the number of doubles.

So the tools you need would be std::sort and std::accumulate and plain old division /.

Jesper Juhl
  • 30,449
  • 3
  • 47
  • 70