different results for expf on x86_64

different results for expf on x86_64

[paul zimmermann]

       Hi,

sorry if I ask a dummy question, I'm new to the list.

While performing my exhaustive tests on binary32 math functions on x86_64,
I've noticed that for some functions, I get different numbers of incorrectly
rounded results on different processors (all with the brand new glibc-2.31).

For example for expf, on a i5-4590, I get 170648 incorrectly rounded results,
while on a E7540, I get only 170646, i.e., 2 less.

The two inputs are 3.256463242e+01 and -6.309946060e+01. For example
on the i5-4590 I get exp(3.256463242e+01) = 1.388801499e+14.

After investigating more, if I simply add printf ("z=%.16e\n", z) after
z = InvLn2N * xd (line 80) in e_expf.c, I get 1.388801415e+14 instead,
like on the E7540.

Is that issue known? How many "flavours" of x86_64 can we have? To which
flavour does correspond the "Wrong count: 170635" at the top of e_expf.c?

Best regards,
Paul

Re: different results for expf on x86_64

[Joseph Myers]

On Tue, 4 Feb 2020, Wilco Dijkstra wrote:

> In this particular case x86_64 uses the generic math code, but compiles expf
> with and without FMA. FMA causes minor rounding differences, but this is fine
> since the new math functions have been designed to be used with FMA (AArch64
> uses FMA by default so there is only one version).

Furthermore, you can get different results with different GCC versions (or 
optimzation options) used to build glibc that is then run on the same 
processor, because the compilers may make different choices about whether 
to contract operations that were written as a*b+c to use FMA.

-- 
Joseph S. Myers
joseph@codesourcery.com

Re: different results for expf on x86_64

[Szabolcs Nagy]

On 04/02/2020 11:29, Wilco Dijkstra wrote:
> Hi Paul,
> 
>> sorry if I ask a dummy question, I'm new to the list.
>>
>> While performing my exhaustive tests on binary32 math functions on x86_64,
>> I've noticed that for some functions, I get different numbers of incorrectly
>> rounded results on different processors (all with the brand new glibc-2.31).
>>
>> For example for expf, on a i5-4590, I get 170648 incorrectly rounded results,
>> while on a E7540, I get only 170646, i.e., 2 less.
>>
>> The two inputs are 3.256463242e+01 and -6.309946060e+01. For example
>> on the i5-4590 I get exp(3.256463242e+01) = 1.388801499e+14.
>>
>> After investigating more, if I simply add printf ("z=%.16e\n", z) after
>> z = InvLn2N * xd (line 80) in e_expf.c, I get 1.388801415e+14 instead,
>> like on the E7540.
>>
>> Is that issue known? How many "flavours" of x86_64 can we have? To which
>> flavour does correspond the "Wrong count: 170635" at the top of e_expf.c?
> 
> x86 has a large number of target specific implementations, including for math
> functions. So you can expect different results depending on the CPU. Even on
> the same CPU there are huge differences between 32-bit and 64-bit mode.

now 32bit x86 uses the generic expf code, instead
of asm based on the f2xm1 instruction, so i would
not expect different behaviour on i686.

> 
> In this particular case x86_64 uses the generic math code, but compiles expf
> with and without FMA. FMA causes minor rounding differences, but this is fine
> since the new math functions have been designed to be used with FMA (AArch64
> uses FMA by default so there is only one version).

this is the likely reason for the difference.
(because of double prec arithmetics for single
prec results, the rounding change caused by
fma should have very minor effect.. 2 cases
out of all binary32 inputs sounds about right)

> 
> One of the goals of rewriting math functions has been to remove the many target
> specific math functions since they are often inaccurate as well as slow. The generic
> math code in C outperforms existing assembler implementations by a large factor!
> 
> Cheers,
> Wilco
> 

Re: different results for expf on x86_64

[Wilco Dijkstra]

Hi Paul,

> sorry if I ask a dummy question, I'm new to the list.
>
> While performing my exhaustive tests on binary32 math functions on x86_64,
> I've noticed that for some functions, I get different numbers of incorrectly
> rounded results on different processors (all with the brand new glibc-2.31).
>
> For example for expf, on a i5-4590, I get 170648 incorrectly rounded results,
> while on a E7540, I get only 170646, i.e., 2 less.
>
> The two inputs are 3.256463242e+01 and -6.309946060e+01. For example
> on the i5-4590 I get exp(3.256463242e+01) = 1.388801499e+14.
>
> After investigating more, if I simply add printf ("z=%.16e\n", z) after
> z = InvLn2N * xd (line 80) in e_expf.c, I get 1.388801415e+14 instead,
> like on the E7540.
>
> Is that issue known? How many "flavours" of x86_64 can we have? To which
> flavour does correspond the "Wrong count: 170635" at the top of e_expf.c?

x86 has a large number of target specific implementations, including for math
functions. So you can expect different results depending on the CPU. Even on
the same CPU there are huge differences between 32-bit and 64-bit mode.

In this particular case x86_64 uses the generic math code, but compiles expf
with and without FMA. FMA causes minor rounding differences, but this is fine
since the new math functions have been designed to be used with FMA (AArch64
uses FMA by default so there is only one version).

One of the goals of rewriting math functions has been to remove the many target
specific math functions since they are often inaccurate as well as slow. The generic
math code in C outperforms existing assembler implementations by a large factor!

Cheers,
Wilco