How to avoid the round-off errors in the larger calculations?

Question

Now I need to sum up more than one thousands of terms and then make the four-dimmensional integral in my Fortran program. I found that there are some numerical errors. Can you give me some suggestions that how to check and avoid the possible round-off errors for this kind of general problem ?

Answer 1

Below are some tips to reduce the effect of round off errors. A short method is to increment the floating point precision, for example from float to double, but many times this is too expensive or not possible.

Kahan summation

In the Kahan summation the idea is to make up for the mistake made in the previous step.

function KahanSum(input)
var sum = 0.0
var c = 0.0                  // A running compensation for lost low-order bits.
for i = 1 to input.length do
    var y = input[i] - c         // So far, so good: c is zero.
    var t = sum + y              // Alas, sum is big, y small, so low-order digits of y are lost.
    c = (t - sum) - y        // (t - sum) cancels the high-order part of y; subtracting y recovers negative (low part of y)
    sum = t                  // Algebraically, c should always be zero. Beware overly-aggressive optimizing compilers!
next i                       // Next time around, the lost low part will be added to y in a fresh attempt.
return sum

The variable c is the correction to use, step by step, with the new sum element input[i]. The error during c evaluation is more small, because the numbers in game are in the same order of magnitude. The important step is that the line

c = (t - sum) - y 

is explicit evaluated.

The implemented code can not be effective due to compiler, see wiki, i.e. the compiler puts all in a unique and direct command so the previous line are not evaluated alone.

Sort

This is a trick quite simple to use, before the sum sort the array and start from the smallest value. The idea is that small sum grows and became more large so the round off errors are mitigate.

It is not guaranteed to work, but it is a good practice to use.

Subtraction, array with mixed signed values

Floating point subtraction is not a stable operation. In the case the array contains mixed signed values, need some extra attentions to avoid that the two numbers in absolute value are very close.

Reduce Operations

This is a general advise, not strictly related to a generic sum of values.

Try to reduce the operation to obtain the result. An example is the Horner's method where a polynomial can be rewrite in the form $$ P_N(x) = a_N + x ( a_{N-1} + x ( a_{N-2} + \ldots + x (a_1 + a_0 x) \ldots)) $$

In this form the polynomial can be evaluate in $N$ additions and $N$ multiplications, against the $N$ additions and $\frac{N(N+1)}{N}$ multiplications of the normal form. (more details, wiki )

Answer 2

You can try starting from the smaller elements and adding them and then going to bigger values. This reduces the roundoff error.

Also try operator precedence for parentheses. To elaborate, assume that you want to compute $$ y = ax + bx^2 + cx^3 $$ note that you have 10 arithmetic (2 add + 8 mult) with some roundoff error. In your implementation instead if you write like this $$ y = x(a + x(b + cx)) $$ you will have 7 operations (2 add + 3 mult) and you may get less roundoff.

Answer 3

Rounding errors is due to finite representation of real numbers in the computer. Computer systems can not work with real numbers accurate, but only with rational approximations thereof. Consequently, the actual numbers can not be represented in the computer than with a finite number of significant figures.

In what follows, it is useful to represent real numbers in base 10, with a f fractional part and an exponent n. x = f • 10n. (2.8) Moreover, for any number outside of 0 through convenient choice of n, the Fractional satisfy 0.1 ≤ | f | <1. For example 3.14 = 0,314 * 10 ^ 1, -0.007856 = • 10 ^ 2 -0.7856. In terms of the convention forms, figures are called fractional part significant numbers.