[<prev] [next>] [<thread-prev] [thread-next>] [day] [month] [year] [list]
Message-ID: <20260810155545.GX27423@brightrain.aerifal.cx>
Date: Mon, 10 Aug 2026 11:55:45 -0400
From: Rich Felker <dalias@...c.org>
To: Sergey Davidoff <shnatsel@...il.com>
Cc: musl@...ts.openwall.com
Subject: Re: Bug in fmaf: subnormals are rounded incorrectly
On Mon, Aug 10, 2026 at 04:11:19PM +0100, Sergey Davidoff wrote:
> Hi all,
>
> The algorithm in src/math/fmaf.c (software FMA implementation for systems
> without hardware support) rounds subnormal results incorrectly.
>
> For example, these bit patterns reinterpreted (not cast) into floats:
>
> x = 0x97000800
> y = 0x1cfff001
> z = 0x00010002
>
> evaluate to 0x00010002 in musl fmaf, but hardware FMA instructions evaluate
> to 0x00010001.
>
> This was tested on x86_64 so I'm sure it's not running into x87 float
> precision issues. The algorithm is simply incorrect. The root cause is this
> line:
>
> if ((u.i & 0x1fffffff) != 0x10000000 /* not a halfway case */
>
> This only checks halfway cases for normal floats, but fails to consider
> subnormals.
If this test is the only thing causing it to be wrong, it could easily
be adjusted to work with subnormals that are halfway cases, or to just
handle all subnormals in this code path if that works.
Rich
Powered by blists - more mailing lists
Confused about mailing lists and their use?
Read about mailing lists on Wikipedia
and check out these
guidelines on proper formatting of your messages.
