Re: [LAD] Is -ffast-math safe for audio?

From: Robin Gareus <robin@email-addr-hidden>
Date: Fri Nov 23 2018 - 00:29:11 EET

Hi Will,

I just ran your code and -ffast-math does not make any difference.

With or without --ffast-math I get "int: 5 rem: 0.049994"

However optimizing the code with `-O2 --ffast-math` does make a
difference because SSE is used.

Do you also use -O2, or -O3 along with --fast-math?

On 11/22/2018 06:30 PM, Will Godfrey wrote:
[...]
> My suspicion is that the difference is due to accumulated rounding
> errors.
>
> Curiously without the decrements the behavior with and without
> -ffast-math seems to be identical well into the millions.
Yep. you do loose precision. In IEEE 32bit float, there is no 0.1.

0.1 would be approximated as (13421773) * 2^(-27), and 0.05 as
(13421773) * 2^(-28) and truncated.

Note that this is an odd number. The last bit of the mantissa
(0x004ccccd) is lost when the exponent changes.

A simpler example to show this is

#include <stdio.h>
#include <math.h>
int main (int argc, char** argv) {
  float a = 0;
  for (int i = 0; i < 100; ++i) {
    a += 0.1f;
    a -= 0.05f;
    a = fmodf (a, 1.f);
  }
  printf ("%f\n", a);
  return 0;
}

using gcc 6.3.0-18+deb9u1, x86_64, this
prints 1.000000 (when compiled with -O0)
and 0.000001 (when compiled with -O2 --fast-math)

The difference here is that the former calls fmodf(), while in the
latter, optimized, case the compiler uses cvtss2sd SSE instead. The
internal limits of float precision differ for those cases.

--ffast-math issues in audio-dsp are usually due to re-ordering and
reciprocal approximations. e.g. instead of (x / 2) the compiler calls
(0.5 * x), which can lead to different results if the inverse value does
not a precise floating point representation. e.g. 1.0 / 0.1 != 10.0

A good read on the subject is
http://www.itu.dk/~sestoft/bachelor/IEEE754_article.pdf (Goldberg,
David: "What Every Computer Scientist Should Know About Floating-Point
Arithmetic").

If you were change the two constants to
  a += 0.5f;
  a -= 0.25f;
the issue vanishes.

Cheers!
robin
_______________________________________________
Linux-audio-dev mailing list
Linux-audio-dev@lists.linuxaudio.org
https://lists.linuxaudio.org/listinfo/linux-audio-dev
Received on Fri Nov 23 04:15:02 2018

This archive was generated by hypermail 2.1.8 : Fri Nov 23 2018 - 04:15:02 EET