Skip to content

Reducing rounding error with compensated summation

A method known as compensated summation1 is implemented for all integrators in grav_sim except WHFast:

When we advance our system by \(\Delta t\), we have

\(x_{n+1} = x_n + \delta x\)

Since \(\delta x\) is very small compared to \(x_n\), many digits of precision will be lost. By compensated summation, we keep track of the losing digits using another variable, which allows us to effectively eliminates round off error with very little cost.

The algorithm is as follows:

  1. Calculate \(\delta x\)
  2. \(x_0 \leftarrow x\)
  3. \(e \leftarrow e + \delta x\)
  4. \(x \leftarrow x_0 + e\)
  5. \(e \leftarrow e + (x_0 - x)\)

Example

Consider

\(x_n\) = 1.3124125125124122

\(\delta x\) = 0.000000012412512412512

\(x_{n+1} = x_n + \delta x\) = 1.312412524924924612512

Let us try to add them together in Python: 

>>> x = 1.3124125125124122
>>> dx = 0.000000012412512412512
>>> x = x + dx
>>> print(f"{x:.19f}")
1.3124125249249245506

The result is 1.3124125249249245506.

Round off error \(\approx 6.1912 \times 10^{-17}\)

Now we try compensated summation: 

>>> x = 1.3124125125124122
>>> dx = 0.000000012412512412512
>>> e = 0.0
>>> x_0 = x
>>> e = e + dx
>>> x = x_0 + e
>>> e = e + (x_0 - x)
>>> e
6.12222940391771e-17
As we can see, the value of \(e\) is very close to the round off error! Therefore, the new variable \(e\) effectively keeps track of the lost digits.

  1. (Ernst Hairer), Gerhard Wanner, and Christian Lubich. Geometric Numerical Integration, chapter VIII, pages 272–274. Springer, Berlin, Heidelberg, 2002. ↩

Comments